6 #include <boost/numeric/itl/itl.hpp> 9 #include <amdis/linearalgebra/mtl/itl/minres.hpp> 10 #include <amdis/linearalgebra/mtl/itl/gcr.hpp> 11 #include <amdis/linearalgebra/mtl/itl/fgmres.hpp> 12 #include <amdis/linearalgebra/mtl/itl/fgmres_householder.hpp> 13 #include <amdis/linearalgebra/mtl/itl/gmres2.hpp> 14 #include <amdis/linearalgebra/mtl/itl/gmres_householder.hpp> 15 #include <amdis/linearalgebra/mtl/itl/preonly.hpp> 17 #include <amdis/CreatorMap.hpp> 18 #include <amdis/Initfile.hpp> 19 #include <amdis/linearalgebra/LinearSolver.hpp> 20 #include <amdis/linearalgebra/mtl/KrylovRunner.hpp> 21 #include <amdis/linearalgebra/mtl/Traits.hpp> 22 #include <amdis/linearalgebra/mtl/UmfpackRunner.hpp> 39 template <
class LinOp,
class X,
class B,
class P,
class I>
40 int operator()(LinOp
const& A,
X& x, B
const& b, P
const& precon, I& iter)
42 return itl::cg(A, x, b, precon, iter);
59 template <
class LinOp,
class X,
class B,
class P,
class I>
60 int operator()(LinOp
const& A,
X& x, B
const& b, P
const& precon, I& iter)
62 return itl::cgs(A, x, b, precon, iter);
79 template <
class LinOp,
class X,
class B,
class P,
class I>
80 int operator()(LinOp
const& A,
X& x, B
const& b, P
const& precon, I& iter)
82 return itl::bicg(A, x, b, precon, iter);
99 template <
class LinOp,
class X,
class B,
class P,
class I>
100 int operator()(LinOp
const& A,
X& x, B
const& b, P
const& precon, I& iter)
102 return itl::bicgstab(A, x, b, precon, iter);
119 template <
class LinOp,
class X,
class B,
class P,
class I>
120 int operator()(LinOp
const& A,
X& x, B
const& b, P
const& precon, I& iter)
122 return itl::bicgstab_2(A, x, b, precon, iter);
138 template <
class LinOp,
class X,
class B,
class P,
class I>
139 int operator()(LinOp
const& A,
X& x, B
const& b, P
const& precon, I& iter)
141 itl::pc::identity<LinOp> id(A);
142 return itl::qmr(A, x, b, precon,
id, iter);
159 template <
class LinOp,
class X,
class B,
class P,
class I>
160 int operator()(LinOp
const& A,
X& x, B
const& b, P
const& precon, I& iter)
162 itl::pc::identity<LinOp> id(A);
163 return itl::tfqmr(A, x, b, precon,
id, iter);
184 template <
class LinOp,
class X,
class B,
class P,
class I>
185 int operator()(LinOp
const& A,
X& x, B
const& b, P
const& precon, I& iter)
187 itl::pc::identity<LinOp> id(A);
188 return itl::bicgstab_ell(A, x, b, precon,
id, iter, ell);
206 explicit gmres_type(std::string
const& name) : restart(30), ortho(GRAM_SCHMIDT)
211 template <
class LinOp,
class X,
class B,
class P,
class I>
212 int operator()(LinOp
const& A,
X& x, B
const& b, P
const& precon, I& iter)
214 itl::pc::identity<LinOp> id(A);
219 return itl::gmres2(A, x, b,
id, precon, iter, restart);
222 return itl::gmres_householder(A, x, b, precon, iter, restart);
228 static constexpr
int GRAM_SCHMIDT = 1;
229 static constexpr
int HOUSEHOLDER = 2;
254 explicit idr_s_type(std::string
const& name) : s(30)
258 template <
class LinOp,
class X,
class B,
class P,
class I>
259 int operator()(LinOp
const& A,
X& x, B
const& b, P
const&, I& iter)
261 itl::pc::identity<LinOp> id(A);
262 return itl::idr_s(A, x, b,
id,
id, iter, s);
279 template <
class LinOp,
class X,
class B,
class P,
class I>
280 int operator()(LinOp
const& A,
X& x, B
const& b, P
const& precon, I& iter)
282 return itl::minres(A, x, b, precon, iter);
301 explicit gcr_type(std::string
const& name) : restart(30)
305 template <
class LinOp,
class X,
class B,
class P,
class I>
306 int operator()(LinOp
const& A,
X& x, B
const& b, P
const& precon, I& iter)
308 itl::pc::identity<LinOp> id(A);
309 return itl::gcr(A, x, b,
id, precon, iter, restart);
328 : restart(30), orthogonalization(GRAM_SCHMIDT)
333 template <
class LinOp,
class X,
class B,
class P,
class I>
334 int operator()(LinOp
const& A,
X& x, B
const& b, P
const& precon, I& iter)
336 itl::pc::identity<LinOp> id(A);
337 switch (orthogonalization)
341 return itl::fgmres(A, x, b,
id, precon, iter, restart);
344 return itl::fgmres_householder(A, x, b, precon, iter, restart);
350 static constexpr
int GRAM_SCHMIDT = 1;
351 static constexpr
int HOUSEHOLDER = 2;
354 int orthogonalization;
369 template <
class LinOp,
class X,
class B,
class P,
class I>
370 int operator()(LinOp
const& A,
X& x, B
const& b, P
const& precon, I& iter)
372 return itl::preonly(A, x, b, precon, iter);
398 template <
class Mat,
class Vec>
403 template <
class ITLSolver>
408 using UmfpackSolverCreator
417 auto cg =
new SolverCreator<cg_solver_type>;
418 Map::addCreator(
"cg", cg);
420 auto cgs =
new SolverCreator<cgs_solver_type>;
421 Map::addCreator(
"cgs", cgs);
423 auto bicg =
new SolverCreator<bicg_solver_type>;
424 Map::addCreator(
"bicg", bicg);
426 auto bcgs =
new SolverCreator<bicgstab_type>;
427 Map::addCreator(
"bicgstab", bcgs);
428 Map::addCreator(
"bcgs", bcgs);
430 auto bcgs2 =
new SolverCreator<bicgstab2_type>;
431 Map::addCreator(
"bicgstab2", bcgs2);
432 Map::addCreator(
"bcgs2", bcgs2);
434 auto qmr =
new SolverCreator<qmr_solver_type>;
435 Map::addCreator(
"qmr", qmr);
437 auto tfqmr =
new SolverCreator<tfqmr_solver_type>;
438 Map::addCreator(
"tfqmr", tfqmr);
440 auto bcgsl =
new SolverCreator<bicgstab_ell_type>;
441 Map::addCreator(
"bicgstab_ell", bcgsl);
442 Map::addCreator(
"bcgsl", bcgsl);
444 auto gmres =
new SolverCreator<gmres_type>;
445 Map::addCreator(
"gmres", gmres);
447 auto fgmres =
new SolverCreator<fgmres_type>;
448 Map::addCreator(
"fgmres", fgmres);
450 auto minres =
new SolverCreator<minres_solver_type>;
451 Map::addCreator(
"minres", minres);
453 auto idrs =
new SolverCreator<idr_s_type>;
454 Map::addCreator(
"idrs", idrs);
456 auto gcr =
new SolverCreator<gcr_type>;
457 Map::addCreator(
"gcr", gcr);
459 auto preonly =
new SolverCreator<preonly_type>;
460 Map::addCreator(
"preonly", preonly);
463 Map::addCreator(
"default", gmres);
465 init_direct(std::is_same<
typename Dune::FieldTraits<typename Mat::value_type>::real_type,
double>{});
468 static void init_direct(std::false_type) {}
469 static void init_direct(std::true_type)
472 auto umfpack =
new UmfpackSolverCreator;
473 Map::addCreator(
"umfpack", umfpack);
476 Map::addCreator(
"direct", umfpack);
ITL_Solver <fgmres_type> implementation of flexible GMRes method.
Definition: ITL_Solver.hpp:324
ITL_Solver <qmr_solver_type> implementation of Quasi-Minimal Residual method.
Definition: ITL_Solver.hpp:134
ITL_Solver <bicgstab_type> implementation of stabilized bi-conjugate gradient method.
Definition: ITL_Solver.hpp:95
ITL_Solver <tfqmr_solver_type> implementation of Transposed-Free Quasi-Minimal Residual method...
Definition: ITL_Solver.hpp:155
ITL_Solver <gmres_type> implementation of generalized minimal residual method.
Definition: ITL_Solver.hpp:202
ITL_Solver <bicgstab2_type> implementation of BiCGStab(l) method with l=2.
Definition: ITL_Solver.hpp:115
A CreatorMap is used to construct objects, which types depends on key words determined at run time...
Definition: CreatorMap.hpp:29
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
Abstract base class for linear solvers.
Definition: LinearSolverInterface.hpp:27
ITL_Solver <gcr_type> implementation of generalized conjugate residual method.
Definition: ITL_Solver.hpp:296
Wrapper class for various solvers. These solvers are parametrized by MatrixType and VectorType...
Definition: LinearSolver.hpp:25
Definition: CreatorMap.hpp:16
ITL_Solver <cg_solver_type> implementation of conjugate gradient method.
Definition: ITL_Solver.hpp:35
ITL_Solver <minres_solver_type> implementation of minimal residual method.
Definition: ITL_Solver.hpp:275
ITL_Solver <bicgstab_ell_type> implementation of stabilized BiCG(ell) method.
Definition: ITL_Solver.hpp:176
ITL_Solver <idr_s_type> implementation of Induced Dimension Reduction method.
Definition: ITL_Solver.hpp:250
ITL_Solver <bicg_solver_type> implementation of bi-conjugate gradient method.
Definition: ITL_Solver.hpp:75
void init(int &argc, char **&argv, std::string const &initFileName="")
Initialized the Environment for MPI.
Definition: AMDiS.hpp:29
ITL_Solver <preonly_type> implementation of preconditioner as.
Definition: ITL_Solver.hpp:365
ITL_Solver <cgs_solver_type> implementation of squared conjugate gradient method. ...
Definition: ITL_Solver.hpp:55
auto X()
Generator for CoordsFunction.
Definition: CoordsGridFunction.hpp:123