5 #include <dune/common/unused.hh> 6 #include <dune/common/version.hh> 8 #include <dune/istl/novlpschwarz.hh> 9 #include <dune/istl/schwarz.hh> 10 #include <dune/istl/paamg/amg.hh> 11 #include <dune/istl/paamg/fastamg.hh> 12 #include <dune/istl/paamg/kamg.hh> 14 #include <amdis/Initfile.hpp> 15 #include <amdis/common/ConceptsBase.hpp> 16 #include <amdis/linearalgebra/istl/Communication.hpp> 17 #include <amdis/linearalgebra/istl/CreatorInterfaces.hpp> 18 #include <amdis/linearalgebra/istl/precompiled/Preconditioners.hpp> 19 #include <amdis/linearalgebra/istl/precompiled/Solvers.hpp> 26 template <
template <
class...>
class AMGSolver,
class Traits>
29 template <
class Smoother>
30 using SmootherArgs =
typename Dune::Amg::SmootherTraits<Smoother>::Arguments;
34 template <
class Traits>
37 using PrecBase =
typename Traits::Prec;
39 template <
class Smoother,
class LinOp,
class Criterion,
class Comm>
40 static std::unique_ptr<PrecBase>
41 create([[maybe_unused]] std::string prefix, LinOp
const& linOp, Criterion
const& criterion, SmootherArgs<Smoother>
const& smootherArgs, Comm
const& comm)
43 using Solver = Dune::Amg::AMG<LinOp, typename Traits::X, Smoother, Comm>;
44 return std::make_unique<Solver>(linOp, criterion, smootherArgs, comm);
54 template <
class Traits>
57 using PrecBase =
typename Traits::Prec;
59 template <
class Smoother,
class LinOp,
class Criterion>
60 static std::unique_ptr<PrecBase>
61 create(std::string prefix, LinOp
const& linOp, Criterion
const& criterion, [[maybe_unused]] SmootherArgs<Smoother>
const& smootherArgs,
62 Dune::Amg::SequentialInformation
const& comm)
64 bool symmetric = Parameters::get<bool>(prefix +
"->symmetric").value_or(
true);
66 using Solver = Dune::Amg::FastAMG<LinOp, typename Traits::X, Dune::Amg::SequentialInformation>;
67 return std::make_unique<Solver>(linOp, criterion, criterion, symmetric, comm);
70 template <
class Smoother,
class LinOp,
class Criterion,
class SmootherArgs,
class Comm,
71 REQUIRES(not std::is_same_v<Comm, Dune::Amg::SequentialInformation>)>
72 static std::unique_ptr<PrecBase>
73 create(std::string, LinOp
const&, Criterion
const&, SmootherArgs
const&, Comm
const&)
75 error_exit(
"Currently only sequential FastAMG supported.");
91 template <
class Traits>
94 using X =
typename Traits::X;
95 using PrecBase =
typename Traits::Prec;
97 template <
class Smoother,
class LinOp,
class Criterion,
class Comm>
98 static std::unique_ptr<PrecBase>
99 create(std::string prefix, LinOp
const& linOp, Criterion
const& criterion, [[maybe_unused]] SmootherArgs<Smoother>
const& smootherArgs, Comm
const& comm)
101 std::string solver = Parameters::get<std::string>(prefix +
"->krylov solver").value_or(
"default");
103 std::size_t maxLevelKrylovSteps = 3;
104 Parameters::get(prefix +
"->krylov solver->maxLevelKrylovSteps", maxLevelKrylovSteps);
105 double minDefectReduction = 1.e-1;
106 Parameters::get(prefix +
"->krylov solver->minDefectReduction", minDefectReduction);
108 if (solver ==
"pcg" ||
111 using Solver = Dune::Amg::KAMG<LinOp, X, Smoother, Comm, Dune::GeneralizedPCGSolver<X>>;
112 return std::make_unique<Solver>(linOp, criterion, smootherArgs, maxLevelKrylovSteps, minDefectReduction, comm);
114 #if DUNE_VERSION_GTE(DUNE_ISTL,2,7) 115 else if (solver ==
"fcg")
117 using Solver = Dune::Amg::KAMG<LinOp, X, Smoother, Comm, Dune::RestartedFCGSolver<X>>;
118 return std::make_unique<Solver>(linOp, criterion, smootherArgs, maxLevelKrylovSteps, minDefectReduction, comm);
120 else if (solver ==
"cfcg")
122 using Solver = Dune::Amg::KAMG<LinOp, X, Smoother, Comm, Dune::CompleteFCGSolver<X>>;
123 return std::make_unique<Solver>(linOp, criterion, smootherArgs, maxLevelKrylovSteps, minDefectReduction, comm);
126 else if (solver ==
"bicgstab" ||
129 using Solver = Dune::Amg::KAMG<LinOp, X, Smoother, Comm, Dune::BiCGSTABSolver<X>>;
130 return std::make_unique<Solver>(linOp, criterion, smootherArgs, maxLevelKrylovSteps, minDefectReduction, comm);
132 else if (solver ==
"minres")
134 using Solver = Dune::Amg::KAMG<LinOp, X, Smoother, Comm, Dune::MINRESSolver<X>>;
135 return std::make_unique<Solver>(linOp, criterion, smootherArgs, maxLevelKrylovSteps, minDefectReduction, comm);
140 else if (solver ==
"gmres")
142 using Solver = Dune::Amg::KAMG<LinOp, X, Smoother, Comm, Dune::RestartedGMResSolver<X>>;
143 return std::make_unique<Solver>(linOp, criterion, smootherArgs, maxLevelKrylovSteps, minDefectReduction, comm);
148 error_exit(
"Unknown coarse space solver {} given for parameter `{}`.", solver, prefix +
"->coarse space solver");
155 template <
class Creator,
class Smoother,
class Traits>
165 template <
template <
class...>
class AMGSolver,
class Traits>
169 using M =
typename Traits::M;
170 using X =
typename Traits::X;
171 using Y =
typename Traits::Y;
176 template <
class Smoother>
183 std::string smoother =
"default";
186 if (smoother ==
"diag" ||
187 smoother ==
"jacobi" ||
188 smoother ==
"default")
191 return creator.createWithString(prefix);
196 else if (smoother ==
"sor")
199 return creator.createWithString(prefix);
204 else if (smoother ==
"gs" ||
205 smoother ==
"gauss_seidel")
208 return creator.createWithString(prefix);
213 else if (smoother ==
"richardson") {
215 return creator.createWithString(prefix);
218 else if (smoother ==
"ssor") {
220 return creator.createWithString(prefix);
223 error_exit(
"Unknown smoother type '{}' given for parameter '{}'", smoother, prefix +
"->smoother");
252 template <
class AMGCreator,
class Smoother,
class Traits>
256 using M =
typename Traits::M;
257 using X =
typename Traits::X;
258 using Y =
typename Traits::Y;
259 using Comm =
typename Traits::Comm;
264 using Interface = Dune::Preconditioner<X,Y>;
265 using LinOperator =
typename Traits::LinOpCreator::Interface;
267 using SolverCategory =
typename Dune::SolverCategory::Category;
269 static const bool is_arithmetic = std::is_arithmetic_v<typename Dune::FieldTraits<M>::field_type>;
270 using Norm = std::conditional_t<is_arithmetic, Dune::Amg::FirstDiagonal, Dune::Amg::RowSum>;
272 using SymCriterion = Dune::Amg::CoarsenCriterion<Dune::Amg::SymmetricCriterion<M,Norm>>;
273 using UnSymCriterion = Dune::Amg::CoarsenCriterion<Dune::Amg::UnSymmetricCriterion<M,Norm>>;
280 return std::make_unique<Self>(prefix);
287 initParameters(prefix);
291 std::unique_ptr<Interface>
create(M
const& mat,
typename Traits::Comm
const& comm)
const override 293 return createImpl1(Dune::SolverCategory::category(comm),mat,comm.get());
297 template <
class Comm,
298 REQUIRES(not std::is_same_v<Comm, Dune::Amg::SequentialInformation>)>
299 std::unique_ptr<Interface>
300 createImpl1(SolverCategory cat, M
const& mat, Comm
const& comm)
const 303 case SolverCategory::sequential:
305 return createImpl1(cat, mat, Dune::Amg::SequentialInformation{});
308 case SolverCategory::overlapping:
310 using LinOp = Dune::OverlappingSchwarzOperator<M,X,Y,Comm>;
311 using ParSmoother = Dune::BlockPreconditioner<X,Y,Comm,Smoother>;
312 LinOp* linOpPtr =
new LinOp(mat, comm);
313 linOperator_.reset(linOpPtr);
314 return createImpl2<ParSmoother>(*linOpPtr, comm);
316 case SolverCategory::nonoverlapping:
318 using LinOp = Dune::NonoverlappingSchwarzOperator<M,X,Y,Comm>;
319 using ParSmoother = Dune::NonoverlappingBlockPreconditioner<Comm,Smoother>;
320 LinOp* linOpPtr =
new LinOp(mat, comm);
321 linOperator_.reset(linOpPtr);
322 return createImpl2<ParSmoother>(*linOpPtr, comm);
326 error_exit(
"Invalid solver category for AMGSolver\n");
331 std::unique_ptr<Interface>
332 createImpl1([[maybe_unused]] SolverCategory cat, M
const& mat, Dune::Amg::SequentialInformation
const& comm)
const 334 assert(cat == SolverCategory::sequential);
335 using LinOp = Dune::MatrixAdapter<M,X,Y>;
336 LinOp* linOpPtr =
new LinOp(mat);
337 linOperator_.reset(linOpPtr);
338 return createImpl2<Smoother>(*linOpPtr, comm);
342 template <
class S,
class LinOp,
class Comm>
343 std::unique_ptr<Interface> createImpl2(LinOp
const& linOp, Comm
const& comm)
const 345 typename Dune::Amg::SmootherTraits<S>::Arguments smootherArgs;
346 initSmootherParameters(prefix_ +
"->smoother", smootherArgs);
348 return symmetricAggregation_
349 ? AMGCreator::template create<S>(prefix_, linOp, SymCriterion(parameters_), smootherArgs, comm)
350 : AMGCreator::template create<S>(prefix_, linOp, UnSymCriterion(parameters_), smootherArgs, comm);
354 void initParameters(std::string
const& prefix)
357 parameters_.setDebugLevel(Parameters::get<int>(prefix +
"->debugLevel").value_or(2));
359 parameters_.setNoPreSmoothSteps(Parameters::get<std::size_t>(prefix +
"->preSmoothSteps").value_or(2));
361 parameters_.setNoPostSmoothSteps(Parameters::get<std::size_t>(prefix +
"->postSmoothSteps").value_or(2));
363 parameters_.setGamma(Parameters::get<std::size_t>(prefix +
"->gamma").value_or(1));
365 parameters_.setAdditive(Parameters::get<bool>(prefix +
"->additive").value_or(
false));
367 symmetricAggregation_ = Parameters::get<bool>(prefix +
"->symmetric aggregation").value_or(
true);
369 initCoarseningParameters(prefix +
"->coarsening");
370 initAggregationParameters(prefix +
"->aggregation");
371 initDependencyParameters(prefix +
"->dependency");
374 void initCoarseningParameters(std::string
const& prefix)
377 parameters_.setMaxLevel(Parameters::get<int>(prefix +
"->maxLevel").value_or(100));
379 parameters_.setCoarsenTarget(Parameters::get<int>(prefix +
"->coarsenTarget").value_or(1000));
381 parameters_.setMinCoarsenRate(Parameters::get<double>(prefix +
"->minCoarsenRate").value_or(1.2));
383 parameters_.setProlongationDampingFactor(Parameters::get<double>(prefix +
"->dampingFactor").value_or(1.6));
386 void initAggregationParameters(std::string
const& prefix)
389 parameters_.setMaxDistance(Parameters::get<std::size_t>(prefix +
"->maxDistance").value_or(2));
391 parameters_.setMinAggregateSize(Parameters::get<std::size_t>(prefix +
"->minAggregateSize").value_or(4));
393 parameters_.setMaxAggregateSize(Parameters::get<std::size_t>(prefix +
"->maxAggregateSize").value_or(6));
395 parameters_.setMaxConnectivity(Parameters::get<std::size_t>(prefix +
"->maxConnectivity").value_or(15));
397 parameters_.setSkipIsolated(Parameters::get<bool>(prefix +
"->skipIsolated").value_or(
false));
400 void initDependencyParameters(std::string
const& prefix)
403 parameters_.setAlpha(Parameters::get<double>(prefix +
"->alpha").value_or(1.0/3.0));
405 parameters_.setBeta(Parameters::get<double>(prefix +
"->beta").value_or(1.e-5));
409 void initSmootherParameters(std::string
const& prefix, SA& smootherArgs)
const 412 Parameters::get(prefix +
"->relaxationFactor", smootherArgs.relaxationFactor);
417 Dune::Amg::Parameters parameters_;
418 bool symmetricAggregation_ =
true;
420 mutable std::shared_ptr<LinOperator> linOperator_;
std::unique_ptr< Super > createWithString(std::string prefix) override
Must be implemented by sub classes of CreatorInterfaceName. Creates a new instance of the sub class o...
Definition: AMGPrecon.hpp:278
Definition: CreatorInterfaces.hpp:11
Definition: AdaptiveGrid.hpp:373
void error_exit(std::string const &str, Args &&... args)
print a message and exit
Definition: Output.hpp:142
Interface for creators with name.
Definition: CreatorInterface.hpp:42
Contains all classes needed for solving linear and non linear equation systems.
Definition: AdaptBase.hpp:6
Definition: AMGPrecon.hpp:27
std::unique_ptr< PreconCreatorBase > createWithString(std::string prefix) override
Must be implemented by sub classes of CreatorInterfaceName. Creates a new instance of the sub class o...
Definition: AMGPrecon.hpp:180
static std::optional< T > get(std::string const &key)
Get parameter-values from parameter-tree.
Definition: Initfile.hpp:25
Definition: AMGPrecon.hpp:156
A creator for AMGPrecon, reads the smoother type from initfile:
Definition: AMGPrecon.hpp:166
std::unique_ptr< Interface > create(M const &mat, typename Traits::Comm const &comm) const override
Implements ISTLPreconCreatorInterface::create.
Definition: AMGPrecon.hpp:291
Definition: AMGPrecon.hpp:276
auto X()
Generator for CoordsFunction.
Definition: CoordsGridFunction.hpp:123