AMDiS 2.10
The Adaptive Multi-Dimensional Simulation Toolbox
Solvers.hpp
1#pragma once
2
3#include <memory>
4
5#include <Eigen/IterativeLinearSolvers>
6#if HAVE_METIS && !HAVE_SCOTCH_METIS
7#include <Eigen/MetisSupport>
8#endif
9#include <Eigen/SparseLU>
10#if HAVE_SUPERLU
11#include <Eigen/SuperLUSupport>
12#endif
13#if HAVE_SUITESPARSE_UMFPACK
14#include <Eigen/UmfPackSupport>
15#endif
16#include <unsupported/Eigen/IterativeSolvers>
17
18#include <amdis/CreatorMap.hpp>
19#include <amdis/Initfile.hpp>
20#include <amdis/Output.hpp>
21
22#include <amdis/linearalgebra/LinearSolverInterface.hpp>
23#include <amdis/linearalgebra/eigen/DirectRunner.hpp>
24#include <amdis/linearalgebra/eigen/IterativeRunner.hpp>
25#include <amdis/linearalgebra/eigen/Traits.hpp>
26
27namespace AMDiS
28{
29 template <class M, class X, class Y, template <class> class IterativeSolver>
31 : public CreatorInterfaceName<LinearSolverInterface<M,X,Y>>
32 {
34 using T = typename M::Scalar;
35
36 template <class Precon>
38
39 template <class Ordering>
40 using IncompleteCholesky =
42
43 std::unique_ptr<SolverBase> createWithString(std::string prefix) override
44 {
45 // get creator string for preconditioner
46 std::string precon = "no";
47 Parameters::get(prefix + "->precon->name", precon);
48
49 if (precon == "diag" || precon == "jacobi") {
50 return std::make_unique<Solver<Eigen::DiagonalPreconditioner<T>>>(prefix);
51 }
52 else if (precon == "ilu") {
53 return std::make_unique<Solver<Eigen::IncompleteLUT<T>>>(prefix);
54 }
55 else if (precon == "ic") {
56 std::string ordering = "amd";
57 Parameters::get(prefix + "->precon->ordering", ordering);
58
59 if (ordering == "amd") {
60 using AMD = Eigen::AMDOrdering<typename M::StorageIndex>;
61 return std::make_unique<IncompleteCholesky<AMD>>(prefix);
62 }
63 else {
64 using NATURAL = Eigen::NaturalOrdering<typename M::StorageIndex>;
65 return std::make_unique<IncompleteCholesky<NATURAL>>(prefix);
66 }
67 }
68 else {
69 return std::make_unique<Solver<Eigen::IdentityPreconditioner>>(prefix);
70 }
71 }
72 };
73
75
85 template <class M, class X, class Y>
87 {
88 template <template <class> class IterativeSolver>
89 using SolverCreator
91
92 template <template <class> class DirectSolver>
93 using DirectSolverCreator
95
97 template <class Precon>
98 using CG = Eigen::ConjugateGradient<M, Eigen::Lower|Eigen::Upper, Precon>;
99
100 template <class Precon>
101 using BCGS = Eigen::BiCGSTAB<M, Precon>;
102
103 template <class Precon>
104 using MINRES = Eigen::MINRES<M, Eigen::Lower|Eigen::Upper, Precon>;
105
106 template <class Precon>
107 using GMRES = Eigen::GMRES<M, Precon>;
108
109 template <class Precon>
110 using DGMRES = Eigen::DGMRES<M, Precon>;
111 // @}
112
115
116 public:
117 static void init()
118 {
119 // conjugate gradient solver
120 auto cg = new SolverCreator<CG>;
121 Map::addCreator("cg", cg);
122
123 // bi-conjugate gradient stabilized solver
124 auto bicgstab = new SolverCreator<BCGS>;
125 Map::addCreator("bicgstab", bicgstab);
126 Map::addCreator("bcgs", bicgstab);
127
128 // Minimal Residual Algorithm (for symmetric matrices)
129 auto minres = new SolverCreator<MINRES>;
130 Map::addCreator("minres", minres);
131
132 // Generalized Minimal Residual Algorithm
133 auto gmres = new SolverCreator<GMRES>;
134 Map::addCreator("gmres", gmres);
135
136 // Restarted GMRES with deflation.
137 auto dgmres = new SolverCreator<DGMRES>;
138 Map::addCreator("dgmres", dgmres);
139
140 // default iterative solver
141 Map::addCreator("default", gmres);
142
143 init_direct(std::is_same<typename Dune::FieldTraits<typename M::Scalar>::real_type, double>{});
144 }
145
146 static void init_direct(std::false_type) {}
147 static void init_direct(std::true_type)
148 {
149#if HAVE_SUITESPARSE_UMFPACK
150 // sparse LU factorization and solver based on UmfPack
151 auto umfpack = new DirectSolverCreator<Eigen::UmfPackLU>;
152 Map::addCreator("umfpack", umfpack);
153#endif
154
155#if HAVE_SUPERLU
156 // sparse direct LU factorization and solver based on the SuperLU library
157 auto superlu = new DirectSolverCreator<Eigen::SuperLU>;
158 Map::addCreator("superlu", superlu);
159#endif
160
161 // default direct solvers
162#if HAVE_SUITESPARSE_UMFPACK
163 Map::addCreator("direct", umfpack);
164#elif HAVE_SUPERLU
165 Map::addCreator("direct", superlu);
166#endif
167 }
168 };
169
170} // end namespace AMDiS
Interface for creators with name.
Definition: CreatorInterface.hpp:44
A CreatorMap is used to construct objects, which types depends on key words determined at run time....
Definition: CreatorMap.hpp:30
Definition: CreatorMap.hpp:16
static std::optional< T > get(std::string const &key)
Get parameter-values from parameter-tree.
Definition: Initfile.hpp:25
Definition: IterativeRunner.hpp:16
Definition: DirectRunner.hpp:27
Default solver creator for iterative solvers.
Definition: Solvers.hpp:32
std::unique_ptr< SolverBase > createWithString(std::string prefix) override
Must be implemented by sub classes of CreatorInterfaceName. Creates a new instance of the sub class o...
Definition: Solvers.hpp:43