ISTL_Preconditioner.hpp 3.8 KB
Newer Older
1
2
3
4
#pragma once

#include <dune/istl/preconditioners.hh>

5
#include <amdis/CreatorInterface.hpp>
6
7
8
#include <amdis/linearalgebra/istl/AMGPrecon.hpp>
#include <amdis/linearalgebra/istl/Fwd.hpp>
#include <amdis/linearalgebra/istl/ISTLPreconCreatorInterface.hpp>
9

10
11
namespace AMDiS
{
12
  template <class> struct Type {};
13
14
15

  // creators for preconditioners, depending on matrix-type
  // ---------------------------------------------------------------------------
16

17
  template <class Precon, class Matrix>
18
19
  struct ISTLPreconCreator
      : public ISTLPreconCreatorInterface<Matrix, typename Precon::domain_type, typename Precon::range_type>
20
  {
21
22
23
24
    using Sol = typename Precon::domain_type;
    using Rhs = typename Precon::range_type;
    using Super = ISTLPreconCreatorInterface<Matrix, Sol, Rhs>;
    using Self = ISTLPreconCreator;
25

26
27
    struct Creator : CreatorInterfaceName<Super>
    {
28
      std::unique_ptr<Super> createWithString(std::string prefix) override
29
      {
30
        return std::make_unique<Self>(prefix);
31
32
      }
    };
33

34
    ISTLPreconCreator(std::string const& prefix)
35
    {
36
      Parameters::get(prefix + "->relaxation", w_);
37
      Parameters::get(prefix + "->iterations", iter_);
38
    }
39

40
    using PreconBase = Dune::Preconditioner<Sol, Rhs>;
41
    std::unique_ptr<PreconBase> create(Matrix const& A) const override
42
43
    {
      return createImpl(A, Type<Precon>{});
44
    }
45

46
  private:
47
48
49
    template <class P>
    std::unique_ptr<P> createImpl(Matrix const& A, Type<P>) const
    {
50
      return std::make_unique<P>(A, iter_, w_);
51
52
    }

53
54
    std::unique_ptr<Dune::SeqILU<Matrix, Sol, Rhs>>
    createImpl(Matrix const& A, Type<Dune::SeqILU<Matrix, Sol, Rhs>>) const
55
    {
56
      return std::make_unique<Dune::SeqILU<Matrix, Sol, Rhs>>(A, iter_, w_);
57
58
    }

59
60
    std::unique_ptr<Dune::Richardson<Sol, Rhs>>
    createImpl(Matrix const& /*A*/, Type<Dune::Richardson<Sol, Rhs>>) const
61
    {
62
      return std::make_unique<Dune::Richardson<Sol, Rhs>>(w_);
63
64
65
    }

    double w_ = 1.0;
66
    int iter_ = 0;
67
  };
68

69
70
71
72
73
74
75
76
77
78

  /// Adds default creators for linear solvers based on `Dune::BCRSMatrix`.
  /**
   * Adds creators for full-matrix aware solvers.
   * - *cg*: conjugate gradient method, \see Dune::CGSolver
   * - *bcgs*: stabilized bi-conjugate gradient method, \see Dune::BiCGSTABSolver
   * - *minres*: Minimal residul method, \see Dune::MINRESSolver
   * - *gmres*: Generalized minimal residula method, \see Dune::RestartedGMResSolver
   * - *umfpack*: external UMFPACK solver, \see Dune::UMFPack
   **/
79
80
  template <class Mat, class Sol, class Rhs>
  class DefaultCreators<ISTLPreconCreatorInterface<Mat, Sol, Rhs>>
81
  {
82
83
    template <class Precon>
    using PreconCreator
84
      = typename ISTLPreconCreator<Precon, Mat>::Creator;
85

86
    using Map = CreatorMap<ISTLPreconCreatorInterface<Mat, Sol, Rhs>>;
87

88
89
  public:
    static void init()
90
    {
91
      auto jacobi = new PreconCreator<Dune::SeqJac<Mat, Sol, Rhs>>;
92
93
      Map::addCreator("diag", jacobi);
      Map::addCreator("jacobi", jacobi);
94

95
      auto gs = new PreconCreator<Dune::SeqGS<Mat, Sol, Rhs>>;
96
97
98
      Map::addCreator("gs", gs);
      Map::addCreator("gauss_seidel", gs);

99
      auto sor = new PreconCreator<Dune::SeqSOR<Mat, Sol, Rhs>>;
100
101
      Map::addCreator("sor", sor);

102
      auto ssor = new PreconCreator<Dune::SeqSSOR<Mat, Sol, Rhs>>;
103
104
      Map::addCreator("ssor", ssor);

105
      auto ilu = new PreconCreator<Dune::SeqILU<Mat, Sol, Rhs>>;
106
107
108
      Map::addCreator("ilu", ilu);
      Map::addCreator("ilu0", ilu);

109
110
111
112
113
114
      auto amg = new AMGPreconCreator<Dune::Amg::AMG, Mat, Sol, Rhs>;
      Map::addCreator("amg", amg);
      auto fastamg = new AMGPreconCreator<Dune::Amg::FastAMG, Mat, Sol, Rhs>;
      Map::addCreator("fastamg", fastamg);

      auto richardson = new PreconCreator<Dune::Richardson<Sol, Rhs>>;
115
116
117
      Map::addCreator("richardson", richardson);
      Map::addCreator("default", richardson);
    }
118
119
120
  };

} // end namespace AMDiS