AMDiS  0.3
The Adaptive Multi-Dimensional Simulation Toolbox
SolverConfig.hpp
1 #pragma once
2 
3 #include <string>
4 
5 #include <Eigen/UmfPackSupport>
6 #include <unsupported/Eigen/IterativeSolvers>
7 
8 #include <amdis/Initfile.hpp>
9 
10 namespace AMDiS
11 {
12  template <class Solver>
13  struct SolverConfig
14  {
15  template <class D>
16  static void init(std::string const& prefix, Eigen::IterativeSolverBase<D>& solver)
17  {
18  using RealScalar = typename Eigen::IterativeSolverBase<D>::RealScalar;
19  auto rtol = Parameters::get<RealScalar>(prefix + "->relative tolerance");
20  if (rtol)
21  solver.setTolerance(rtol.value());
22 
23  int maxIter = 500;
24  Parameters::get(prefix + "->max iteration", maxIter);
25  solver.setMaxIterations(maxIter);
26  }
27 
28  // fall-back default implementation
29  template <class D>
30  static void init(std::string const&, Eigen::SparseSolverBase<D>&) {}
31  };
32 
33  template <class M, class P>
34  struct SolverConfig<Eigen::GMRES<M,P>>
35  {
36  using Base = Eigen::IterativeSolverBase<Eigen::GMRES<M,P>>;
37 
38  static void init(std::string const& prefix, Eigen::GMRES<M,P>& solver)
39  {
40  SolverConfig<Base>::init(prefix, solver);
41 
42  auto restart = Parameters::get<int>(prefix + "->restart");
43  if (restart)
44  solver.set_restart(restart.value());
45  }
46  };
47 
48  template <class M, class P>
49  struct SolverConfig<Eigen::DGMRES<M,P>>
50  {
51  using Base = Eigen::IterativeSolverBase<Eigen::DGMRES<M,P>>;
52 
53  static void init(std::string const& prefix, Eigen::DGMRES<M,P>& solver)
54  {
55  SolverConfig<Base>::init(prefix, solver);
56 
57  auto restart = Parameters::get<int>(prefix + "->restart");
58  if (restart)
59  solver.set_restart(restart.value());
60 
61  // number of eigenvalues to deflate at each restart
62  auto neigs = Parameters::get<int>(prefix + "->num eigenvalues");
63  if (neigs)
64  solver.setEigenv(neigs.value());
65 
66  // maximum size of the deflation subspace
67  auto maxNeig = Parameters::get<int>(prefix + "->max num eigenvalues");
68  if (maxNeig)
69  solver.setMaxEigenv(maxNeig.value());
70  }
71  };
72 
73 #if HAVE_SUITESPARSE_UMFPACK
74  template <class M>
75  struct SolverConfig<Eigen::UmfPackLU<M>>
76  {
77  static void init(std::string const& prefix, Eigen::UmfPackLU<M>& solver)
78  {
79  [[maybe_unused]] auto& control = solver.umfpackControl();
80  // TODO: initialized umfpack parameters
81  }
82  };
83 #endif
84 
85 } // end namespace AMDiS
Contains all classes needed for solving linear and non linear equation systems.
Definition: AdaptBase.hpp:6
static std::optional< T > get(std::string const &key)
Get parameter-values from parameter-tree.
Definition: Initfile.hpp:25
Definition: SolverConfig.hpp:13
void init(int &argc, char **&argv, std::string const &initFileName="")
Initialized the Environment for MPI.
Definition: AMDiS.hpp:29