diff --git a/dune/istl/paamg/amg.hh b/dune/istl/paamg/amg.hh index 0d311967e6961d514709cdff2467973649e52790..a901c6d26c492de12fd6c3ebf48b1f2395b78659 100644 --- a/dune/istl/paamg/amg.hh +++ b/dune/istl/paamg/amg.hh @@ -123,6 +123,32 @@ namespace Dune const SmootherArgs& smootherArgs=SmootherArgs(), const ParallelInformation& pinfo=ParallelInformation()); + /*! + \brief Constructor an AMG via ParameterTree. + + \param fineOperator The operator on the fine level. + \param configuration ParameterTree containing AMG parameters. + \param pinfo Optionally, specify ParallelInformation + + ParameterTree Key | Meaning + --------------------------|------------ + smootherIterations | The number of iterations to perform. + smootherRelaxation | Common parameter defined [here](@ref ISTL_Factory_Common_Params). + maxLevel | Maximum number of levels allowed in the hierarchy. + coarsenTarget | Maximum number of unknowns on the coarsest level. + prolongationDampingFactor | Damping factor for the prolongation. + alpha | Scaling avlue for marking connections as strong. + beta | Treshold for marking nodes as isolated. + additive | Whether to use additive multigrid. + gamma | 1 for V-cycle, 2 for W-cycle. + preSteps | Number of presmoothing steps. + postSteps | Number of postsmoothing steps. + verbosity | Output verbosity. + + See \ref ISTL_Factory for the ParameterTree layout and examples. + */ + AMG(std::shared_ptr<const Operator> fineOperator, const ParameterTree& configuration, const ParallelInformation& pinfo=ParallelInformation()); + /** * @brief Copy constructor. */ @@ -354,6 +380,60 @@ namespace Dune createHierarchies(criterion, matrixptr, pinfo); } + template<class M, class X, class S, class PI, class A> + AMG<M,X,S,PI,A>::AMG(std::shared_ptr<const Operator> matrixptr, + const ParameterTree& configuration, + const ParallelInformation& pinfo) : + smoothers_(new Hierarchy<Smoother,A>), + solver_(), rhs_(), lhs_(), update_(), scalarProduct_(), buildHierarchy_(true), + coarsesolverconverged(true), coarseSmoother_(), + category_(SolverCategory::category(pinfo)) + { + + if (configuration.hasKey ("smootherIterations")) + smootherArgs_.iterations = configuration.get<int>("smootherIterations"); + + if (configuration.hasKey ("smootherRelaxation")) + smootherArgs_.relaxationFactor = configuration.get<typename SmootherArgs::RelaxationFactor>("smootherRelaxation"); + + typedef Amg::RowSum Norm; + typedef Dune::Amg::CoarsenCriterion<Dune::Amg::UnSymmetricCriterion<typename M::matrix_type,Norm> > Criterion; + + Criterion criterion (configuration.get<int>("maxLevel")); + + if (configuration.hasKey ("coarsenTarget")) + criterion.setCoarsenTarget (configuration.get<int>("coarsenTarget")); + if (configuration.hasKey ("prolongationDampingFactor")) + criterion.setProlongationDampingFactor (configuration.get<double>("prolongationDampingFactor")); + + if (configuration.hasKey ("alpha")) + criterion.setAlpha (configuration.get<double> ("alpha")); + + if (configuration.hasKey ("beta")) + criterion.setBeta (configuration.get<double> ("beta")); + + if (configuration.hasKey ("gamma")) + criterion.setGamma (configuration.get<std::size_t> ("gamma")); + gamma_ = criterion.getGamma(); + + if (configuration.hasKey ("additive")) + criterion.setAdditive (configuration.get<bool>("additive")); + additive = criterion.getAdditive(); + + if (configuration.hasKey ("preSteps")) + criterion.setNoPreSmoothSteps (configuration.get<std::size_t> ("preSteps")); + preSteps_ = criterion.getNoPreSmoothSteps (); + + if (configuration.hasKey ("postSteps")) + criterion.setNoPostSmoothSteps (configuration.get<std::size_t> ("preSteps")); + postSteps_ = criterion.getNoPostSmoothSteps (); + + verbosity_ = configuration.get("verbosity",2); + criterion.setDebugLevel (verbosity_); + + createHierarchies(criterion, matrixptr, pinfo); + } + template <class Matrix, class Vector> struct DirectSolverSelector @@ -423,7 +503,7 @@ namespace Dune template<class C> void AMG<M,X,S,PI,A>::createHierarchies(C& criterion, const std::shared_ptr<const Operator>& matrixptr, - const PI& pinfo) + const PI& pinfo) { Timer watch; matrices_ = std::make_shared<OperatorHierarchy>( diff --git a/dune/istl/preconditioners.hh b/dune/istl/preconditioners.hh index cc9f7a85234412d9f670c6d8ae3caf554f0ef10e..0901e2423200c301ac0969472f2bb403e036104c 100644 --- a/dune/istl/preconditioners.hh +++ b/dune/istl/preconditioners.hh @@ -12,6 +12,7 @@ #include <dune/common/simd/simd.hh> #include <dune/common/unused.hh> +#include <dune/common/parametertree.hh> #include "preconditioner.hh" #include "solver.hh" @@ -159,6 +160,27 @@ namespace Dune { CheckIfDiagonalPresent<M,l>::check(_A_); } + /*! + \brief Constructor. + + \param A The matrix to operate on. + \param configuration ParameterTree containing preconditioner parameters. + + ParameterTree Key | Meaning + ------------------|------------ + iterations | The number of iterations to perform. + relaxation | Common parameter defined [here](@ref ISTL_Factory_Common_Params). + + See \ref ISTL_Factory for the ParameterTree layout and examples. + */ + SeqSSOR (const M& A, const ParameterTree& configuration) + : _A_(A) + { + _n = configuration.get<field_type>("iterations"); + _w = configuration.get<field_type>("relaxation"); + CheckIfDiagonalPresent<M,l>::check(_A_); + } + /*! \brief Prepare the preconditioner. @@ -249,6 +271,17 @@ namespace Dune { CheckIfDiagonalPresent<M,l>::check(_A_); } + /*! + \copydoc SeqSSOR::SeqSSOR(const M&,const ParameterTree&) + */ + SeqSOR (const M& A, const ParameterTree& configuration) + : _A_(A) + { + _n = configuration.get<field_type>("iterations"); + _w = configuration.get<field_type>("relaxation"); + CheckIfDiagonalPresent<M,l>::check(_A_); + } + /*! \brief Prepare the preconditioner. @@ -354,6 +387,17 @@ namespace Dune { CheckIfDiagonalPresent<M,l>::check(_A_); } + /*! + \copydoc SeqSSOR::SeqSSOR(const M&,const ParameterTree&) + */ + SeqGS (const M& A, const ParameterTree& configuration) + : _A_(A) + { + _n = configuration.get<field_type>("iterations"); + _w = configuration.get<field_type>("relaxation"); + CheckIfDiagonalPresent<M,l>::check(_A_); + } + /*! \brief Prepare the preconditioner. @@ -440,6 +484,17 @@ namespace Dune { CheckIfDiagonalPresent<M,l>::check(_A_); } + /*! + \copydoc SeqSSOR::SeqSSOR(const M&,const ParameterTree&) + */ + SeqJac (const M& A, const ParameterTree& configuration) + : _A_(A) + { + _n = configuration.get<field_type>("iterations"); + _w = configuration.get<field_type>("relaxation"); + CheckIfDiagonalPresent<M,l>::check(_A_); + } + /*! \brief Prepare the preconditioner. @@ -679,6 +734,25 @@ namespace Dune { bilu0_decomposition(ILU); } + /*! + \brief Constructor. + + \param A The matrix to operate on. + \param configuration ParameterTree containing preconditioner parameters. + + ParameterTree Key | Meaning + ------------------|------------ + relaxation | Common parameter defined [here](@ref ISTL_Factory_Common_Params). + + See \ref ISTL_Factory for the ParameterTree layout and examples. + */ + SeqILU0 (const M& A, const ParameterTree& configuration) + : ILU(A) // copy A + { + _w = configuration.get<field_type>("relaxation"); + bilu0_decomposition(ILU); + } + /*! \brief Prepare the preconditioner. @@ -769,6 +843,17 @@ namespace Dune { bilu_decomposition(A,n,ILU); } + /*! + \copydoc SeqSSOR::SeqSSOR(const M&,const ParameterTree&) + */ + SeqILUn (const M& A, const ParameterTree& configuration) + : ILU(A.N(),A.M(),M::row_wise) + { + _n = configuration.get<field_type>("iterations"); + _w = configuration.get<field_type>("relaxation"); + bilu_decomposition(A,_n,ILU); + } + /*! \brief Prepare the preconditioner. @@ -847,6 +932,14 @@ namespace Dune { _w(w) {} + /*! + \copydoc SeqILU0::SeqILU0(const M&,const ParameterTree&) + */ + Richardson (const ParameterTree& configuration) + { + _w = configuration.get<field_type>("relaxation"); + } + /*! \brief Prepare the preconditioner. diff --git a/dune/istl/solver.hh b/dune/istl/solver.hh index 7d718cec5f2f25f5cdc6cba5de48e0995b3979fd..9803a5b93eb86522b5a63d44c8ae4725e0ca77d7 100644 --- a/dune/istl/solver.hh +++ b/dune/istl/solver.hh @@ -13,6 +13,7 @@ #include <dune/common/shared_ptr.hh> #include <dune/common/simd/io.hh> #include <dune/common/simd/simd.hh> +#include <dune/common/parametertree.hh> #include <dune/common/timer.hh> #include "solvertype.hh" @@ -269,7 +270,21 @@ namespace Dune DUNE_THROW(InvalidSolverCategory, "LinearOperator and ScalarProduct must have the same SolverCategory!"); } - /** + IterativeSolver (std::shared_ptr<LinearOperator<X,Y> > op, std::shared_ptr<Preconditioner<X,X> > prec, const ParameterTree& configuration) : + IterativeSolver(op,prec, + configuration.get<real_type>("reduction"), + configuration.get<int>("maxit"), + configuration.get<int>("verbose")) + {} + + IterativeSolver (std::shared_ptr<LinearOperator<X,Y> > op, std::shared_ptr<ScalarProduct<X> > sp, std::shared_ptr<Preconditioner<X,X> > prec, const ParameterTree& configuration) : + IterativeSolver(op,sp,prec, + configuration.get<real_type>("reduction"), + configuration.get<int>("maxit"), + configuration.get<int>("verbose")) + {} + + /** \brief General constructor to initialize an iterative solver \param op The operator we solve. @@ -305,7 +320,6 @@ namespace Dune DUNE_THROW(InvalidSolverCategory, "LinearOperator and ScalarProduct must have the same SolverCategory!"); } - // #warning actually we want to have this as the default and just implement the second one // //! \copydoc InverseOperator::apply(X&,Y&,InverseOperatorResult&) // virtual void apply (X& x, Y& b, InverseOperatorResult& res) diff --git a/dune/istl/solvers.hh b/dune/istl/solvers.hh index 322653c5a31e3fc705829ae487c978642f1fba61..eaf95a0f3c7d600c41ee6854d7bd2d86d25a747d 100644 --- a/dune/istl/solvers.hh +++ b/dune/istl/solvers.hh @@ -844,6 +844,16 @@ namespace Dune { _restart(restart) {} + RestartedGMResSolver (std::shared_ptr<LinearOperator<X,Y> > op, std::shared_ptr<Preconditioner<X,X> > prec, const ParameterTree& configuration) : + IterativeSolver<X,Y>::IterativeSolver(op,prec,configuration), + _restart(configuration.get<int>("restart")) + {} + + RestartedGMResSolver (std::shared_ptr<LinearOperator<X,Y> > op, std::shared_ptr<ScalarProduct<X> > sp, std::shared_ptr<Preconditioner<X,X> > prec, const ParameterTree& configuration) : + IterativeSolver<X,Y>::IterativeSolver(op,sp,prec,configuration), + _restart(configuration.get<int>("restart")) + {} + /*! \brief Set up RestartedGMResSolver solver. @@ -1297,6 +1307,16 @@ private: _restart(restart) {} + + GeneralizedPCGSolver (std::shared_ptr<LinearOperator<X,X> > op, std::shared_ptr<Preconditioner<X,X> > prec, const ParameterTree& configuration) : + IterativeSolver<X,X>::IterativeSolver(op,prec,configuration), + _restart(configuration.get<int>("restart")) + {} + + GeneralizedPCGSolver (std::shared_ptr<LinearOperator<X,X> > op, std::shared_ptr<ScalarProduct<X> > sp, std::shared_ptr<Preconditioner<X,X> > prec, const ParameterTree& configuration) : + IterativeSolver<X,X>::IterativeSolver(op,sp,prec,configuration), + _restart(configuration.get<int>("restart")) + {} /*! \brief Set up nonlinear preconditioned conjugate gradient solver.