My Project
elasticity_preconditioners.hpp
Go to the documentation of this file.
1 //==============================================================================
11 //==============================================================================
12 #ifndef ELASTICITY_PRECONDITIONERS_HPP_
13 #define ELASTICITY_PRECONDITIONERS_HPP_
14 
15 #include <opm/common/utility/platform_dependent/disable_warnings.h>
16 
17 #include <dune/common/fmatrix.hh>
18 #include <dune/istl/bcrsmatrix.hh>
19 #include <dune/istl/matrixmatrix.hh>
20 #include <dune/istl/ilu.hh>
21 #include <dune/istl/solvers.hh>
22 #include <dune/istl/preconditioners.hh>
23 #include <dune/istl/superlu.hh>
24 #include <dune/istl/umfpack.hh>
25 #include <dune/istl/paamg/amg.hh>
26 #include <dune/istl/paamg/fastamg.hh>
27 #include <dune/istl/paamg/twolevelmethod.hh>
28 #include <dune/istl/overlappingschwarz.hh>
29 
30 #include <opm/common/utility/platform_dependent/reenable_warnings.h>
31 
32 #include <opm/grid/CpGrid.hpp>
35 
36 
37 namespace Opm {
38 namespace Elasticity {
39 
40 #if defined(HAVE_UMFPACK)
41 typedef Dune::UMFPack<Matrix> LUSolver;
42 #elif defined(HAVE_SUPERLU)
43 typedef Dune::SuperLU<Matrix> LUSolver;
44 #else
45 static_assert(false, "Enable either SuperLU or UMFPACK");
46 #endif
47 
49 typedef Dune::MatrixAdapter<Matrix,Vector,Vector> Operator;
50 
52 typedef Dune::SeqSSOR<Matrix, Vector, Vector> SSORSmoother;
53 
55 typedef Dune::SeqJac<Matrix, Vector, Vector> JACSmoother;
56 
58 #if DUNE_VERSION_NEWER(DUNE_ISTL, 2, 7)
59 typedef Dune::SeqILU<Matrix, Vector, Vector> ILUSmoother;
60 #else
61 typedef Dune::SeqILU0<Matrix, Vector, Vector> ILUSmoother;
62 #endif
63 
65 typedef Dune::SeqOverlappingSchwarz<Matrix,Vector,
66  Dune::SymmetricMultiplicativeSchwarzMode, LUSolver> SchwarzSmoother;
67 
69 struct Schwarz {
70  typedef Dune::SeqOverlappingSchwarz<Matrix, Vector,
71  Dune::SymmetricMultiplicativeSchwarzMode,
72  LUSolver> type;
81  static std::shared_ptr<type>
82  setup(int /* pre */, int /* post */, int /* target */, int /* zcells */,
83  std::shared_ptr<Operator>& op, const Dune::CpGrid& gv,
84  ASMHandler<Dune::CpGrid>& A, bool& copy)
85  {
86  return std::shared_ptr<type>(setup2(op, gv, A, copy));
87  }
88 
94  static type* setup2(std::shared_ptr<Operator>& op, const Dune::CpGrid& gv,
95  ASMHandler<Dune::CpGrid>& A, bool& copy);
96 };
97 
99 template<class Smoother>
100 struct AMG1 {
102  typedef Dune::Amg::FirstDiagonal CouplingMetric;
103 
105  typedef Dune::Amg::SymmetricCriterion<Matrix, CouplingMetric> CritBase;
106 
108  typedef Dune::Amg::CoarsenCriterion<CritBase> Criterion;
109 
110  typedef Dune::Amg::AMG<Operator, Vector, Smoother> type;
111 
120  static std::shared_ptr<type>
121  setup(int pre, int post, int target, int zcells,
122  std::shared_ptr<Operator>& op, const Dune::CpGrid& /* gv */,
123  ASMHandler<Dune::CpGrid>& /* A */, bool& copy)
124  {
125  Criterion crit;
127  args.relaxationFactor = 1.0;
128  crit.setCoarsenTarget(target);
129  crit.setGamma(1);
130  crit.setNoPreSmoothSteps(pre);
131  crit.setNoPostSmoothSteps(post);
132  crit.setDefaultValuesIsotropic(3, zcells);
133 
134  std::cout << "\t collapsing 2x2x" << zcells << " cells per level" << std::endl;
135  copy = true;
136  return std::shared_ptr<type>(new type(*op, crit, args));
137  }
138 };
139 
141 struct FastAMG {
142  typedef Dune::Amg::FastAMG<Operator, Vector> type;
143 
152  static std::shared_ptr<type>
153  setup(int pre, int post, int target, int zcells,
154  std::shared_ptr<Operator>& op, const Dune::CpGrid& gv,
155  ASMHandler<Dune::CpGrid>& A, bool& copy);
156 };
157 
158 
160  template<class Smoother>
161 struct AMG2Level {
163  typedef Dune::Amg::AggregationLevelTransferPolicy<Operator,
165 
166  typedef Dune::Amg::LevelTransferPolicy<Operator, Operator> LevelTransferPolicy;
167 
168  typedef Dune::Amg::OneStepAMGCoarseSolverPolicy<Operator, Smoother,
169  typename AMG1<Smoother>::Criterion> CoarsePolicy;
170 
171  typedef typename Dune::Amg::SmootherTraits<Smoother>::Arguments SmootherArgs;
172 
173  typedef Dune::Amg::TwoLevelMethod<Operator, CoarsePolicy, Schwarz::type> type;
174 
183  static std::shared_ptr<type>
184  setup(int pre, int post, int target, int zcells,
185  std::shared_ptr<Operator>& op, const Dune::CpGrid& gv,
186  ASMHandler<Dune::CpGrid>& A, bool& copy)
187  {
188  typename AMG1<Smoother>::Criterion crit;
189  SmootherArgs args;
190  args.relaxationFactor = 1.0;
191  crit.setCoarsenTarget(target);
192  crit.setGamma(1);
193  crit.setNoPreSmoothSteps(pre);
194  crit.setNoPostSmoothSteps(post);
195  crit.setDefaultValuesIsotropic(3, zcells);
196  CoarsePolicy coarsePolicy(args, crit);
197  TransferPolicy policy(crit);
198  std::shared_ptr<Schwarz::type> fsp(Schwarz::setup2(op, gv, A, copy));
199  copy = true;
200  return std::shared_ptr<type>(new type(*op, fsp, policy, coarsePolicy, pre, post));
201  }
202 };
203 
204 }
205 }
206 
207 #endif
Class handling finite element assembly.
Class handling finite element assembly.
Definition: asmhandler.hpp:35
Dune::MatrixAdapter< Matrix, Vector, Vector > Operator
A linear operator.
Definition: elasticity_preconditioners.hpp:45
Dune::SeqJac< Matrix, Vector, Vector > JACSmoother
GJ AMG smoother.
Definition: elasticity_preconditioners.hpp:55
Dune::SeqILU0< Matrix, Vector, Vector > ILUSmoother
ILU0 AMG smoother.
Definition: elasticity_preconditioners.hpp:61
Dune::SeqOverlappingSchwarz< Matrix, Vector, Dune::SymmetricMultiplicativeSchwarzMode, LUSolver > SchwarzSmoother
Schwarz + ILU0 AMG smoother.
Definition: elasticity_preconditioners.hpp:66
Dune::SeqSSOR< Matrix, Vector, Vector > SSORSmoother
SSOR AMG smoother.
Definition: elasticity_preconditioners.hpp:52
Smoother
Smoother used in the AMG.
Definition: elasticity_upscale.hpp:75
Helper class with some matrix operations.
Dune::BCRSMatrix< Dune::FieldMatrix< double, 1, 1 > > Matrix
A sparse matrix holding our operator.
Definition: matrixops.hpp:27
Dune::BlockVector< Dune::FieldVector< double, 1 > > Vector
A vector holding our RHS.
Definition: matrixops.hpp:33
Inverting small matrices.
Definition: ImplicitAssembly.hpp:43
An AMG.
Definition: elasticity_preconditioners.hpp:100
Dune::Amg::CoarsenCriterion< CritBase > Criterion
The coarsening criterion used in the AMG.
Definition: elasticity_preconditioners.hpp:108
Dune::Amg::FirstDiagonal CouplingMetric
The coupling metric used in the AMG.
Definition: elasticity_preconditioners.hpp:102
static std::shared_ptr< type > setup(int pre, int post, int target, int zcells, std::shared_ptr< Operator > &op, const Dune::CpGrid &, ASMHandler< Dune::CpGrid > &, bool &copy)
Setup preconditioner.
Definition: elasticity_preconditioners.hpp:121
Dune::Amg::SymmetricCriterion< Matrix, CouplingMetric > CritBase
The coupling criterion used in the AMG.
Definition: elasticity_preconditioners.hpp:105
A two-level method with a coarse AMG solver.
Definition: elasticity_preconditioners.hpp:161
static std::shared_ptr< type > setup(int pre, int post, int target, int zcells, std::shared_ptr< Operator > &op, const Dune::CpGrid &gv, ASMHandler< Dune::CpGrid > &A, bool &copy)
Setup preconditioner.
Definition: elasticity_preconditioners.hpp:184
Dune::Amg::AggregationLevelTransferPolicy< Operator, typename AMG1< Smoother >::Criterion > TransferPolicy
AMG transfer policy.
Definition: elasticity_preconditioners.hpp:164
A FastAMG.
Definition: elasticity_preconditioners.hpp:141
static std::shared_ptr< type > setup(int pre, int post, int target, int zcells, std::shared_ptr< Operator > &op, const Dune::CpGrid &gv, ASMHandler< Dune::CpGrid > &A, bool &copy)
Setup preconditioner.
Definition: elasticity_preconditioners.cpp:20
Overlapping Schwarz preconditioner.
Definition: elasticity_preconditioners.hpp:69
static type * setup2(std::shared_ptr< Operator > &op, const Dune::CpGrid &gv, ASMHandler< Dune::CpGrid > &A, bool &copy)
Setup preconditioner.
Definition: elasticity_preconditioners.cpp:37
static std::shared_ptr< type > setup(int, int, int, int, std::shared_ptr< Operator > &op, const Dune::CpGrid &gv, ASMHandler< Dune::CpGrid > &A, bool &copy)
Setup preconditioner.
Definition: elasticity_preconditioners.hpp:82