AMDiS  0.3
The Adaptive Multi-Dimensional Simulation Toolbox
ITL_Preconditioner.hpp
1 #pragma once
2 
3 // MTL4 headers
4 #include <boost/numeric/itl/itl.hpp>
5 #include <boost/numeric/itl/pc/ilu_0.hpp>
6 #include <boost/numeric/itl/pc/ic_0.hpp>
7 #include <boost/numeric/mtl/vector/assigner.hpp>
8 
9 #include <amdis/CreatorMap.hpp>
10 #include <amdis/linearalgebra/mtl/HyprePrecon.hpp>
11 #include <amdis/linearalgebra/mtl/Preconditioner.hpp>
12 #include <amdis/linearalgebra/mtl/SolverPrecon.hpp>
13 #include <amdis/linearalgebra/mtl/Traits.hpp>
14 #include <amdis/linearalgebra/mtl/itl/masslumping.hpp>
15 
16 namespace AMDiS
17 {
26  template <class Matrix>
27  using DiagonalPreconditioner = itl::pc::diagonal<Matrix>;
28 
37  template <class Matrix>
38  using MassLumpingPreconditioner = itl::pc::masslumping<Matrix>;
39 
40 
49  template <class Matrix>
50  using IdentityPreconditioner = itl::pc::identity<Matrix>;
51 
52 
63  template <class Matrix>
64  using ILUPreconditioner = itl::pc::ilu_0<Matrix>;
65 
66 
75  template <class Matrix>
76  using ICPreconditioner = itl::pc::ic_0<Matrix>;
77 
78 
79  template <class Mat, class Vec>
81  {
83 
84  template <template <class> class ITLPrecon>
85  using PreconCreator
87 
89 
90  public:
91  static void init()
92  {
93  auto pc_diag = new PreconCreator<DiagonalPreconditioner>;
94  Map::addCreator("diag", pc_diag);
95  Map::addCreator("jacobi", pc_diag);
96 
97  auto pc_mass = new PreconCreator<MassLumpingPreconditioner>;
98  Map::addCreator("masslumping", pc_mass);
99 
100  auto pc_ilu = new PreconCreator<ILUPreconditioner>;
101  Map::addCreator("ilu", pc_ilu);
102  Map::addCreator("ilu0", pc_ilu);
103 
104  auto pc_ic = new PreconCreator<ICPreconditioner>;
105  Map::addCreator("ic", pc_ic);
106  Map::addCreator("ic0", pc_ic);
107 
108  auto pc_id = new PreconCreator<IdentityPreconditioner>;
109  Map::addCreator("identity", pc_id);
110  Map::addCreator("no", pc_id);
111 
112  auto pc_solver = new typename SolverPrecon<Mat,Vec>::Creator;
113  Map::addCreator("solver", pc_solver);
114 
115  Map::addCreator("default", pc_id);
116 
117  init_hypre(std::is_same<typename Dune::FieldTraits<typename Mat::value_type>::real_type, double>{});
118  }
119 
120  static void init_hypre(std::false_type) {}
121  static void init_hypre(std::true_type)
122  {
123 #if AMDIS_HAS_HYPRE && HAVE_MPI
124  auto pc_hypre = new typename HyprePrecon<Mat,Vec>::Creator;
125  Map::addCreator("hypre", pc_hypre);
126 #endif
127  }
128  };
129 
130  template <class Mat, class Vec, class X>
131  itl::pc::solver<PreconditionerInterface<Mat,Vec>, X, false>
132  solve(PreconditionerInterface<Mat,Vec> const& P, X const& vin)
133  {
134  return {P, vin};
135  }
136 
137  template <class Mat, class Vec, class X>
138  itl::pc::solver<PreconditionerInterface<Mat,Vec>, X, true>
139  adjoint_solve(PreconditionerInterface<Mat,Vec> const& P, X const& vin)
140  {
141  return {P, vin};
142  }
143 
144 } // namespace AMDiS
A creator to be used instead of the constructor.
Definition: SolverPrecon.hpp:31
A creator to be used instead of the constructor.
Definition: Preconditioner.hpp:26
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
Definition: CreatorMap.hpp:16
Interface for Preconditioner y = M*x.
Definition: PreconditionerInterface.hpp:9
void init(int &argc, char **&argv, std::string const &initFileName="")
Initialized the Environment for MPI.
Definition: AMDiS.hpp:29
auto X()
Generator for CoordsFunction.
Definition: CoordsGridFunction.hpp:123