12 #ifndef ELASTICITY_PRECONDITIONERS_HPP_
13 #define ELASTICITY_PRECONDITIONERS_HPP_
15 #include <dune/common/fmatrix.hh>
16 #include <dune/istl/bcrsmatrix.hh>
17 #include <dune/istl/matrixmatrix.hh>
18 #include <dune/istl/ilu.hh>
19 #include <dune/istl/solvers.hh>
20 #include <dune/istl/preconditioners.hh>
21 #include <dune/grid/CpGrid.hpp>
26 #include <dune/istl/superlu.hh>
27 #include <dune/istl/umfpack.hh>
28 #include <dune/istl/paamg/amg.hh>
29 #include <dune/istl/paamg/fastamg.hh>
30 #include <dune/istl/paamg/twolevelmethod.hh>
31 #include <dune/istl/overlappingschwarz.hh>
34 namespace Elasticity {
36 #if defined(HAVE_UMFPACK)
37 typedef Dune::UMFPack<Matrix> LUSolver;
38 #elif defined(HAVE_SUPERLU)
39 typedef Dune::SuperLU<Matrix> LUSolver;
41 static_assert(
"Enable either SuperLU or UMFPACK");
45 typedef Dune::MatrixAdapter<Matrix,Vector,Vector> Operator;
48 typedef Dune::SeqSSOR<Matrix, Vector, Vector> SSORSmoother;
51 typedef Dune::SeqJac<Matrix, Vector, Vector> JACSmoother;
54 typedef Dune::SeqILU0<Matrix, Vector, Vector> ILUSmoother;
58 Dune::SymmetricMultiplicativeSchwarzMode, LUSolver> SchwarzSmoother;
63 Dune::SymmetricMultiplicativeSchwarzMode,
73 static std::shared_ptr<type>
74 setup(
int pre,
int post,
int target,
int zcells,
75 std::shared_ptr<Operator>& op,
const Dune::CpGrid& gv,
76 ASMHandler<Dune::CpGrid>& A,
bool& copy)
78 return std::shared_ptr<type>(setup2(op, gv, A, copy));
86 static type* setup2(std::shared_ptr<Operator>& op,
const Dune::CpGrid& gv,
87 ASMHandler<Dune::CpGrid>& A,
bool& copy);
91 template<
class Smoother>
94 typedef Dune::Amg::FirstDiagonal CouplingMetric;
97 typedef Dune::Amg::SymmetricCriterion<Matrix, CouplingMetric> CritBase;
100 typedef Dune::Amg::CoarsenCriterion<CritBase> Criterion;
102 typedef Dune::Amg::AMG<Operator, Vector, Smoother> type;
112 static std::shared_ptr<type>
113 setup(
int pre,
int post,
int target,
int zcells,
114 std::shared_ptr<Operator>& op,
const Dune::CpGrid& gv,
115 ASMHandler<Dune::CpGrid>& A,
bool& copy)
118 typename AMG1<Smoother>::type::SmootherArgs args;
119 args.relaxationFactor = 1.0;
120 crit.setCoarsenTarget(target);
122 crit.setNoPreSmoothSteps(pre);
123 crit.setNoPostSmoothSteps(post);
124 crit.setDefaultValuesIsotropic(3, zcells);
126 std::cout <<
"\t collapsing 2x2x" << zcells <<
" cells per level" << std::endl;
128 return std::shared_ptr<type>(
new type(*op, crit, args));
134 typedef Dune::Amg::FastAMG<Operator, Vector> type;
144 static std::shared_ptr<type>
145 setup(
int pre,
int post,
int target,
int zcells,
146 std::shared_ptr<Operator>& op,
const Dune::CpGrid& gv,
147 ASMHandler<Dune::CpGrid>& A,
bool& copy);
152 template<
class Smoother>
155 typedef Dune::Amg::AggregationLevelTransferPolicy<Operator,
156 typename AMG1<Smoother>::Criterion> TransferPolicy;
158 typedef Dune::Amg::LevelTransferPolicy<Operator, Operator> LevelTransferPolicy;
160 typedef Dune::Amg::OneStepAMGCoarseSolverPolicy<Operator,
Smoother,
161 typename AMG1<Smoother>::Criterion> CoarsePolicy;
163 typedef typename Dune::Amg::SmootherTraits<Smoother>::Arguments SmootherArgs;
165 typedef Dune::Amg::TwoLevelMethod<Operator, CoarsePolicy, Schwarz::type> type;
175 static std::shared_ptr<type>
176 setup(
int pre,
int post,
int target,
int zcells,
177 std::shared_ptr<Operator>& op,
const Dune::CpGrid& gv,
178 ASMHandler<Dune::CpGrid>& A,
bool& copy)
180 typename AMG1<Smoother>::Criterion crit;
182 args.relaxationFactor = 1.0;
183 crit.setCoarsenTarget(target);
185 crit.setNoPreSmoothSteps(pre);
186 crit.setNoPostSmoothSteps(post);
187 crit.setDefaultValuesIsotropic(3, zcells);
188 CoarsePolicy coarsePolicy(args, crit);
189 TransferPolicy policy(crit);
190 Dune::shared_ptr<Schwarz::type> fsp(Schwarz::setup2(op, gv, A, copy));
192 return std::shared_ptr<type>(
new type(*op, fsp, policy, coarsePolicy, pre, post));
Definition: applier.hpp:18
Smoother
Smoother used in the AMG.
Definition: elasticity_upscale.hpp:70
Dune::BCRSMatrix< Dune::FieldMatrix< double, 1, 1 > > Matrix
A sparse matrix holding our operator.
Definition: matrixops.hpp:23
Helper class with some matrix operations.
Class handling finite element assembly.
Dune::BlockVector< Dune::FieldVector< double, 1 > > Vector
A vector holding our RHS.
Definition: matrixops.hpp:29