BITL_Solver.hpp 2.74 KB
Newer Older
1
2
#pragma once

3
4
#include <amdis/CreatorMap.hpp>
#include <amdis/linear_algebra/mtl/ITL_Solver.hpp>
5
6
7

namespace AMDiS
{
8
9
10
11
12
13
14
15
16
17
18
19
20
21
  /// Adds default creators for linear solvers based on `BlockMTLMatrix`.
  /**
   * Adds creators for block matrix aware solvers.
   * - *cg*: conjugate gradient method, \see cg_solver_type
   * - *cgs*: stabilized conjugate gradient mtheod, \see cgs_solver_type
   * - *bcgs*: stabilized bi-conjugate gradient method, \see bicgstab_type
   * - *tfqmr*: Transposed-Free Quasi-Minimal Residual method, \see tfqmr_solver_type
   * - *bcgsl*: stabilized BiCG(ell) method, \see bicgstab_ell_type
   * - *gmres*: Generalized minimal residula method, \see gmres_type
   * - *fgmres*: Flexible GMRES, \see fgmres_type
   * - *minres*: Minimal residul method, \see minres_solver_type
   * - *gcr*: Generalized conjugate residual method, \see gcr_type
   * - *umfpack*: external UMFPACK solver, \see UmfpackRunner
   **/
22
  template <class SubMatrix, std::size_t _N, std::size_t _M, class Vector>
23
  class DefaultCreators<LinearSolverInterface<BlockMTLMatrix<SubMatrix, _N, _M>, Vector, Vector>>
24
  {
25
    using Matrix = BlockMTLMatrix<SubMatrix, _N, _M>;
26
    using SolverBase = LinearSolverInterface<Matrix, Vector, Vector>;
27

28
    template <class ITLSolver>
29
    using SolverCreator
30
      = typename LinearSolver<Matrix, Vector,
31
          KrylovRunner<ITLSolver, Matrix, BaseVector_t<Vector>>>::Creator;
32

33
#ifdef HAVE_UMFPACK
34
    using UmfpackSolverCreator
35
      = typename LinearSolver<Matrix, Vector,
36
          UmfpackRunner<Matrix, BaseVector_t<Vector>>>::Creator;
37
#endif
38

39
    using Map = CreatorMap<SolverBase>;
40

41
  public:
42
    static void init()
43
    {
44
      auto cg = new SolverCreator<cg_solver_type>;
45
      Map::addCreator("cg", cg);
46

47
      auto cgs = new SolverCreator<cgs_solver_type>;
48
      Map::addCreator("cgs", cgs);
49

50
      auto bcgs = new SolverCreator<bicgstab_type>;
51
52
      Map::addCreator("bicgstab", bcgs);
      Map::addCreator("bcgs", bcgs);
53

54
      auto tfqmr = new SolverCreator<tfqmr_solver_type>;
55
      Map::addCreator("tfqmr", tfqmr);
56

57
      auto bcgsl = new SolverCreator<bicgstab_ell_type>;
58
59
      Map::addCreator("bicgstab_ell", bcgsl);
      Map::addCreator("bcgsl", bcgsl);
60

61
      auto gmres = new SolverCreator<gmres_type>;
62
      Map::addCreator("gmres", gmres);
63

64
      auto fgmres = new SolverCreator<fgmres_type>;
65
      Map::addCreator("fgmres", fgmres);
66

67
      auto minres = new SolverCreator<minres_solver_type>;
68
      Map::addCreator("minres", minres);
69

70
      auto gcr = new SolverCreator<gcr_type>;
71
      Map::addCreator("gcr", gcr);
72

73
#ifdef HAVE_UMFPACK
74
      auto umfpack = new UmfpackSolverCreator;
75
76
      Map::addCreator("umfpack", umfpack);
      Map::addCreator("direct", umfpack);
77
#endif
78

79
      Map::addCreator("default", gmres);
80
81
82
83
    }
  };

} // end namespace AMDiS