Removed old amg/cpr code that can be replaced.

The replacement consists of using the FlexibleSolver code.

CMakeLists_files.cmake

@@ -170,10 +170,8 @@ list (APPEND PUBLIC_HEADER_FILES
   opm/simulators/linalg/bda/openclSolverBackend.hpp
   opm/simulators/linalg/bda/MultisegmentWellContribution.hpp
   opm/simulators/linalg/bda/WellContributions.hpp
-  opm/simulators/linalg/BlackoilAmg.hpp
   opm/simulators/linalg/amgcpr.hh
   opm/simulators/linalg/twolevelmethodcpr.hh
-  opm/simulators/linalg/CPRPreconditioner.hpp
   opm/simulators/linalg/ExtractParallelGridInformationToISTL.hpp
   opm/simulators/linalg/FlexibleSolver.hpp
   opm/simulators/linalg/FlexibleSolver_impl.hpp

flow/flow_blackoil_dunecpr.cpp

@@ -50,14 +50,6 @@ namespace Opm {
         static constexpr int value = 100;
     };
     template<class TypeTag>
-    struct UseAmg<TypeTag, TTag::EclFlowProblemSimple> { // probably not used
-        static constexpr bool value = true;
-    };
-    template<class TypeTag>
-    struct UseCpr<TypeTag, TTag::EclFlowProblemSimple> {
-        static constexpr bool value = true;
-    };
-    template<class TypeTag>
     struct CprMaxEllIter<TypeTag,TTag::EclFlowProblemSimple> {
         static constexpr int value = 1;
     };
@@ -70,17 +62,9 @@ namespace Opm {
         static constexpr int value = 3;
     };
     template<class TypeTag>
-    struct CprSolverVerbose<TypeTag, TTag::EclFlowProblemSimple> {
-        static constexpr int value = 0;
-    };
-    template<class TypeTag>
     struct LinearSolverConfiguration<TypeTag, TTag::EclFlowProblemSimple> {
         static constexpr auto value = "ilu0";
     };
-    template<class TypeTag>
-    struct SystemStrategy<TypeTag, TTag::EclFlowProblemSimple> {
-        static constexpr auto value = "quasiimpes";
-    };
 
     template<class TypeTag>
     struct FluidSystem<TypeTag, TTag::EclFlowProblemSimple>

flow/flow_ebos_blackoil.cpp

@@ -16,9 +16,6 @@
 */
 #include "config.h"
 
-// Define making clear that the simulator supports AMG
-#define FLOW_SUPPORT_AMG 1
-
 #include <flow/flow_ebos_blackoil.hpp>
 
 #include <opm/material/common/ResetLocale.hpp>

flow/flow_ebos_gasoil.cpp

@@ -16,9 +16,6 @@
 */
 #include "config.h"
 
-// Define making clear that the simulator supports AMG
-#define FLOW_SUPPORT_AMG 1
-
 #include <flow/flow_ebos_gasoil.hpp>
 
 #include <opm/material/common/ResetLocale.hpp>

flow/flow_ebos_oilwater.cpp

@@ -16,9 +16,6 @@
 */
 #include "config.h"
 
-// Define making clear that the simulator supports AMG
-#define FLOW_SUPPORT_AMG 1
-
 #include <flow/flow_ebos_oilwater.hpp>
 
 #include <opm/material/common/ResetLocale.hpp>

flow/flow_ebos_oilwater_brine.cpp

@@ -16,9 +16,6 @@
 */
 #include "config.h"
 
-// Define making clear that the simulator supports AMG
-#define FLOW_SUPPORT_AMG 1
-
 #include <flow/flow_ebos_oilwater_brine.hpp>
 
 #include <opm/material/common/ResetLocale.hpp>

flow/flow_ebos_oilwater_polymer.cpp

@@ -16,9 +16,6 @@
 */
 #include "config.h"
 
-// Define making clear that the simulator supports AMG
-#define FLOW_SUPPORT_AMG 1
-
 #include <flow/flow_ebos_oilwater_polymer.hpp>
 
 #include <opm/material/common/ResetLocale.hpp>

flow/flow_ebos_oilwater_polymer_injectivity.cpp

@@ -16,9 +16,6 @@
 */
 #include "config.h"
 
-// Define making clear that the simulator supports AMG
-#define FLOW_SUPPORT_AMG 1
-
 #include <flow/flow_ebos_oilwater_polymer_injectivity.hpp>
 
 #include <opm/material/common/ResetLocale.hpp>

opm/simulators/linalg/CPRPreconditioner.hpp

@@ -1,528 +0,0 @@
-/*
-  Copyright 2014 SINTEF ICT, Applied Mathematics.
-  Copyright 2014 IRIS AS.
-  Copyright 2015 Dr. Blatt - HPC-Simulation-Software & Services
-  Copyright 2015 NTNU
-  Copyright 2015 Statoil AS
-
-  This file is part of the Open Porous Media project (OPM).
-
-  OPM is free software: you can redistribute it and/or modify
-  it under the terms of the GNU General Public License as published by
-  the Free Software Foundation, either version 3 of the License, or
-  (at your option) any later version.
-
-  OPM is distributed in the hope that it will be useful,
-  but WITHOUT ANY WARRANTY; without even the implied warranty of
-  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
-  GNU General Public License for more details.
-
-  You should have received a copy of the GNU General Public License
-  along with OPM.  If not, see <http://www.gnu.org/licenses/>.
-*/
-
-#ifndef OPM_CPRPRECONDITIONER_HEADER_INCLUDED
-#define OPM_CPRPRECONDITIONER_HEADER_INCLUDED
-
-#include <memory>
-#include <type_traits>
-
-#include <opm/common/utility/platform_dependent/disable_warnings.h>
-#include <opm/common/utility/parameters/ParameterGroup.hpp>
-#include <opm/simulators/linalg/ParallelOverlappingILU0.hpp>
-#include <opm/simulators/linalg/FlowLinearSolverParameters.hpp>
-
-#include <dune/istl/bvector.hh>
-#include <dune/istl/bcrsmatrix.hh>
-#include <dune/istl/operators.hh>
-#include <dune/istl/io.hh>
-#include <dune/istl/owneroverlapcopy.hh>
-#include <dune/istl/preconditioners.hh>
-#include <dune/istl/schwarz.hh>
-#include <dune/istl/solvers.hh>
-#include <dune/istl/paamg/amg.hh>
-#include <dune/istl/paamg/kamg.hh>
-#include <dune/istl/paamg/pinfo.hh>
-
-#include <opm/common/utility/platform_dependent/reenable_warnings.h>
-
-#include <opm/common/ErrorMacros.hpp>
-#include <opm/common/Exceptions.hpp>
-
-#include <opm/simulators/linalg/ParallelRestrictedAdditiveSchwarz.hpp>
-#include <opm/simulators/linalg/ParallelOverlappingILU0.hpp>
-namespace Opm
-{
-
-template<typename O, typename S, typename C,
-         typename P, std::size_t COMPONENT_INDEX, std::size_t VARIABLE_INDEX>
-class BlackoilAmg;
-
-namespace Amg
-{
-    template<int Row, int Column>
-    class Element;
-}
-
-namespace ISTLUtility
-{
-
-template<class T>
-void setILUParameters(Opm::ParallelOverlappingILU0Args<T>& args,
-                      const CPRParameter& params)
-{
-    args.setN(params.cpr_ilu_n_);
-    args.setMilu(params.cpr_ilu_milu_);
-}
-
-template<class T>
-void setILUParameters(Opm::ParallelOverlappingILU0Args<T>& args,
-                      MILU_VARIANT milu, int n=0)
-{
-    args.setN(n);
-    args.setMilu(milu);
-}
-
-template<class S, class P>
-void setILUParameters(S&, const P&)
-{}
-
-template<class S, class P>
-void setILUParameters(S&, bool, int)
-{}
-
-///
-/// \brief A traits class for selecting the types of the preconditioner.
-///
-/// \tparam M The type of the matrix.
-/// \tparam X The type of the domain of the linear problem.
-/// \tparam Y The type of the range of the linear problem.
-/// \tparam P The type of the parallel information.
-////
-template<class M, class X, class Y, class P>
-struct CPRSelector
-{
-    /// \brief The information about the parallelization and communication
-    using ParallelInformation = P;
-    /// \brief The operator type;
-    typedef Dune::OverlappingSchwarzOperator<M,X,X,ParallelInformation> Operator;
-    typedef ParallelOverlappingILU0<M,X, X, ParallelInformation>
-    EllipticPreconditioner;
-    /// \brief The type of the unique pointer to the preconditioner of the elliptic part.
-    typedef std::unique_ptr<EllipticPreconditioner>
-    EllipticPreconditionerPointer;
-
-    typedef EllipticPreconditioner Smoother;
-    typedef Dune::Amg::AMG<Operator, X, Smoother, ParallelInformation> AMG;
-
-    /// \brief creates an Operator from the matrix
-    /// \param M The matrix to use.
-    /// \param p The parallel information to use.
-    static Operator* makeOperator(const M& m, const ParallelInformation& p)
-    {
-        return new Operator(m, p);
-    }
-};
-
-template<class M, class X, class Y>
-struct CPRSelector<M,X,Y,Dune::Amg::SequentialInformation>
-{
-    /// \brief The information about the parallelization and communication
-    typedef Dune::Amg::SequentialInformation ParallelInformation;
-    /// \brief The operator type;
-    typedef Dune::MatrixAdapter<M, X, Y> Operator;
-    /// \brief The type of the preconditioner used for the elliptic part.
-    typedef ParallelOverlappingILU0<M,X, X, ParallelInformation>
-    EllipticPreconditioner;
-    /// \brief The type of the unique pointer to the preconditioner of the elliptic part.
-    typedef std::unique_ptr<EllipticPreconditioner> EllipticPreconditionerPointer;
-
-    /// \brief type of AMG used to precondition the elliptic system.
-    typedef EllipticPreconditioner Smoother;
-    typedef Dune::Amg::AMG<Operator, X, Smoother, ParallelInformation> AMG;
-
-    /// \brief creates an Operator from the matrix
-    /// \param M The matrix to use.
-    /// \param p The parallel information to use.
-    static Operator* makeOperator(const M& m, const ParallelInformation&)
-    {
-        return new Operator(m);
-    }
-
-};
-//! \brief Creates and initializes a unique pointer to an sequential ILU0 preconditioner.
-//! \param A     The matrix of the linear system to solve.
-//! \param relax The relaxation factor to use.
-template<class M, class X, class C>
-std::shared_ptr<ParallelOverlappingILU0<M,X,X,C> >
-createILU0Ptr(const M& A, const C& comm, double relax, MILU_VARIANT milu)
-{
-    return std::make_shared<ParallelOverlappingILU0<M,X,X,C> >(A, comm, relax, milu);
-}
-//! \brief Creates and initializes a shared pointer to an ILUn preconditioner.
-//! \param A     The matrix of the linear system to solve.
-//! \param ilu_n The n parameter for the extension of the nonzero pattern.
-//! \param relax The relaxation factor to use.
-template<class M, class X, class C>
-std::shared_ptr<ParallelOverlappingILU0<M,X,X,C> >
-createILUnPtr(const M& A, const C& comm, int ilu_n, double relax, MILU_VARIANT milu)
-{
-    return std::make_shared<ParallelOverlappingILU0<M,X,X,C> >( A, comm, ilu_n, relax, milu );
-}
-/// \brief Creates the elliptic preconditioner (ILU0)
-/// \param Ae    The matrix of the elliptic system.
-/// \param relax The relaxation parameter for ILU0.
-/// \param milu  If true, the modified ilu approach is used. Dropped elements
-///              will get added to the diagonal of U to preserve the row sum
-///              for constant vectors (Ae = LUe).
-/// \param comm  The object describing the parallelization information and communication.
-template<class M, class X=typename M::range_type, class P>
-typename CPRSelector<M,X,X,P>::EllipticPreconditionerPointer
-createEllipticPreconditionerPointer(const M& Ae, double relax,
-                                    MILU_VARIANT milu, const P& comm)
-{
-    typedef typename CPRSelector<M,X,X,P >
-        ::EllipticPreconditioner ParallelPreconditioner;
-
-    typedef typename CPRSelector<M,X,X,P>
-        ::EllipticPreconditionerPointer EllipticPreconditionerPointer;
-    return EllipticPreconditionerPointer(new ParallelPreconditioner(Ae, comm, relax, milu));
-}
-
-template < class C, class Op, class P, class S, std::size_t PressureEqnIndex, std::size_t PressureVarIndex, class Vector>
-inline void
-createAMGPreconditionerPointer(Op& opA, const double relax, const P& comm,
-                               std::unique_ptr< BlackoilAmg<Op,S,C,P,PressureEqnIndex,PressureVarIndex> >& amgPtr,
-                               const CPRParameter& params,
-                               const Vector& weights)
-{
-    using AMG = BlackoilAmg<Op,S,C,P,PressureEqnIndex,PressureVarIndex>;
-    int verbosity = 0;
-    if (comm.communicator().rank() == 0) {
-        verbosity = params.cpr_solver_verbose_;
-    }
-    // TODO: revise choice of parameters
-    int coarsenTarget=1200;
-    using Criterion = C;
-    Criterion criterion(15, coarsenTarget);
-    criterion.setDebugLevel( verbosity ); // no debug information, 1 for printing hierarchy information
-    criterion.setDefaultValuesIsotropic(2);
-    criterion.setNoPostSmoothSteps( 1 );
-    criterion.setNoPreSmoothSteps( 1 );
-
-    // Since DUNE 2.2 we also need to pass the smoother args instead of steps directly
-    typedef typename AMG::Smoother Smoother;
-    typedef typename Dune::Amg::SmootherTraits<Smoother>::Arguments  SmootherArgs;
-    SmootherArgs  smootherArgs;
-    smootherArgs.iterations = 1;
-    smootherArgs.relaxationFactor = relax;
-    setILUParameters(smootherArgs, params);
-
-    amgPtr.reset( new AMG( params, weights, opA, criterion, smootherArgs, comm ) );
-}
-
-template < class C, class Op, class P, class AMG >
-inline void
-createAMGPreconditionerPointer(Op& opA, const double relax, const MILU_VARIANT milu, const P& comm, std::unique_ptr< AMG >& amgPtr)
-{
-    // TODO: revise choice of parameters
-    int coarsenTarget=1200;
-    using Criterion = C;
-    Criterion criterion(15, coarsenTarget);
-    criterion.setDebugLevel( 0 ); // no debug information, 1 for printing hierarchy information
-    criterion.setDefaultValuesIsotropic(2);
-    criterion.setNoPostSmoothSteps( 1 );
-    criterion.setNoPreSmoothSteps( 1 );
-
-    // for DUNE 2.2 we also need to pass the smoother args
-    typedef typename AMG::Smoother Smoother;
-    typedef typename Dune::Amg::SmootherTraits<Smoother>::Arguments  SmootherArgs;
-    SmootherArgs  smootherArgs;
-    smootherArgs.iterations = 1;
-    smootherArgs.relaxationFactor = relax;
-    setILUParameters(smootherArgs, milu);
-
-    amgPtr.reset( new AMG(opA, criterion, smootherArgs, comm ) );
-}
-
-/// \brief Creates the elliptic preconditioner (ILU0)
-/// \param opA     The operator representing the matrix of the system.
-/// \param relax   The relaxation parameter for ILU0.
-/// \param comm    The object describing the parallelization information and communication.
-//  \param amgPtr  The unique_ptr to be filled (return)
-template < int PressureEqnIndex, int PressureVarIndex, class Op, class P, class AMG >
-inline void
-createAMGPreconditionerPointer( Op& opA, const double relax, const MILU_VARIANT milu, const P& comm, std::unique_ptr< AMG >& amgPtr )
-{
-    // type of matrix
-    typedef typename Op::matrix_type  M;
-
-    // The coupling metric used in the AMG
-    typedef Opm::Amg::Element<PressureEqnIndex, PressureVarIndex> CouplingMetric;
-
-    // The coupling criterion used in the AMG
-    typedef Dune::Amg::SymmetricCriterion<M, CouplingMetric> CritBase;
-
-    // The coarsening criterion used in the AMG
-    typedef Dune::Amg::CoarsenCriterion<CritBase> Criterion;
-
-    createAMGPreconditionerPointer<Criterion>(opA, relax, milu, comm, amgPtr);
-}
-
-} // end namespace ISTLUtility
-
-    /*!
-      \brief CPR preconditioner.
-
-      This is a two-stage preconditioner, combining an elliptic-type
-      partial solution with ILU0 for the whole system.
-
-      \tparam M The matrix type to operate on
-      \tparam X Type of the update
-      \tparam Y Type of the defect
-      \tparam P Type of the parallel information. If not provided
-                this will be Dune::Amg::SequentialInformation.
-                The preconditioner is parallel if this is
-                Dune::OwnerOverlapCopyCommunication<int,int>
-    */
-    template<class M, class X, class Y,
-             class P=Dune::Amg::SequentialInformation>
-    class CPRPreconditioner : public Dune::Preconditioner<X,Y>
-    {
-        // prohibit copying for now
-        CPRPreconditioner( const CPRPreconditioner& );
-
-    public:
-        //! \brief The type describing the parallel information
-        typedef P ParallelInformation;
-        //! \brief The matrix type the preconditioner is for.
-        typedef typename std::remove_const<M>::type matrix_type;
-        //! \brief The domain type of the preconditioner.
-        typedef X domain_type;
-        //! \brief The range type of the preconditioner.
-        typedef Y range_type;
-        //! \brief The field type of the preconditioner.
-        typedef typename X::field_type field_type;
-
-        // define the category
-        Dune::SolverCategory::Category category() const override
-        {
-          return std::is_same<P,Dune::Amg::SequentialInformation>::value ?
-                 Dune::SolverCategory::sequential : Dune::SolverCategory::overlapping;
-        }
-
-        typedef ISTLUtility::CPRSelector<M,X,X,P>  CPRSelectorType ;
-
-        //! \brief Elliptic Operator
-        typedef typename CPRSelectorType::Operator Operator;
-
-        //! \brief preconditioner for the whole system (here either ILU(0) or ILU(n)
-        typedef Dune::Preconditioner<X,X> WholeSystemPreconditioner;
-
-        //! \brief type of the unique pointer to the ilu-0 preconditioner
-        //! used the for the elliptic system
-        typedef typename CPRSelectorType::EllipticPreconditionerPointer
-        EllipticPreconditionerPointer;
-
-        //! \brief amg preconditioner for the elliptic system
-        typedef typename CPRSelectorType::AMG  AMG;
-
-        /*! \brief Constructor.
-
-          Constructor gets all parameters to operate the prec.
-          \param A       The matrix to operate on.
-          \param Ae      The top-left elliptic part of A.
-          \param relax   The ILU0 relaxation factor.
-          \param useAMG  if true, AMG is used as a preconditioner for the elliptic sub-system, otherwise ilu-0 (default)
-          \param useBiCG if true, BiCG solver is used (default), otherwise CG solver
-          \param paralleInformation The information about the parallelization, if this is a
-                                    parallel run
-        */
-        CPRPreconditioner (const CPRParameter& param, const M& A, const M& Ae,
-                           const ParallelInformation& comm=ParallelInformation(),
-                           const ParallelInformation& commAe=ParallelInformation())
-            : param_( param ),
-              A_(A),
-              Ae_(Ae),
-              de_( Ae_.N() ),
-              ve_( Ae_.M() ),
-              dmodified_( A_.N() ),
-              opAe_(CPRSelectorType::makeOperator(Ae_, commAe)),
-              precond_(), // ilu0 preconditioner for elliptic system
-              amg_(),     // amg  preconditioner for elliptic system
-              pre_(), // copy A will be made be the preconditioner
-              vilu_( A_.N() ),
-              comm_(comm),
-              commAe_(commAe)
-        {
-            // create appropriate preconditioner for elliptic system
-            createEllipticPreconditioner( param_.cpr_use_amg_, commAe_ );
-
-            // create the preconditioner for the whole system.
-            if( param_.cpr_ilu_n_ == 0 ) {
-                pre_ = ISTLUtility::createILU0Ptr<M,X>( A_, comm, param_.cpr_relax_, param_.cpr_ilu_milu_ );
-            }
-            else {
-                pre_ = ISTLUtility::createILUnPtr<M,X>( A_, comm, param_.cpr_ilu_n_, param_.cpr_relax_, param_.cpr_ilu_milu_);
-            }
-        }
-
-        /*!
-          \brief Prepare the preconditioner.
-
-          \copydoc Preconditioner::pre(X&,Y&)
-        */
-        virtual void pre (X& /*x*/, Y& /*b*/)
-        {
-        }
-
-        /*!
-          \brief Apply the preconditoner.
-
-          \copydoc Preconditioner::apply(X&,const Y&)
-        */
-        virtual void apply (X& v, const Y& d)
-        {
-            // Extract part of d corresponding to elliptic part.
-            // Note: Assumes that the elliptic part comes first.
-            std::copy_n(d.begin(), de_.size(), de_.begin());
-
-            // Solve elliptic part, extend solution to full.
-            // reset result
-            ve_ = 0;
-            solveElliptic( ve_, de_ );
-
-            //reset return value
-            v = 0.0;
-            // Again assuming that the elliptic part comes first.
-            std::copy(ve_.begin(), ve_.end(), v.begin());
-
-            // Subtract elliptic residual from initial residual.
-            // dmodified = d - A * vfull
-            dmodified_ = d;
-            A_.mmv(v, dmodified_);
-            // A is not parallel, do communication manually.
-            comm_.copyOwnerToAll(dmodified_, dmodified_);
-
-            // Apply Preconditioner for whole system (relax will be applied already)
-            pre_->apply( vilu_, dmodified_);
-
-            // don't apply relaxation if relax_ == 1
-            if( std::abs( param_.cpr_relax_ - 1.0 ) < 1e-12 ) {
-                v += vilu_;
-            }
-            else {
-                v *= param_.cpr_relax_;
-                v += vilu_;
-            }
-        }
-
-        /*!
-          \brief Clean up.
-
-          \copydoc Preconditioner::post(X&)
-        */
-        virtual void post (X& /*x*/)
-        {
-        }
-
-     protected:
-        void solveElliptic(Y& x, Y& de)
-        {
-            // Linear solver parameters
-            const double tolerance = param_.cpr_solver_tol_;
-            const int maxit        = param_.cpr_max_ell_iter_;
-            const int verbosity    = ( param_.cpr_solver_verbose_ &&
-                                       comm_.communicator().rank()==0 ) ? 1 : 0;
-
-            // operator result containing iterations etc.
-            Dune::InverseOperatorResult result;
-
-            // the scalar product chooser
-            auto sp = Dune::createScalarProduct<X,ParallelInformation>(commAe_, category());
-
-            if( amg_ )
-            {
-                // Solve system with AMG
-                if( param_.cpr_use_bicgstab_ ) {
-                    Dune::BiCGSTABSolver<X> linsolve(*opAe_, *sp, (*amg_), tolerance, maxit, verbosity);
-                    linsolve.apply(x, de, result);
-                }
-                else {
-                    Dune::CGSolver<X> linsolve(*opAe_, *sp, (*amg_), tolerance, maxit, verbosity);
-                    linsolve.apply(x, de, result);
-                }
-            }
-            else
-            {
-                assert( precond_ );
-                // Solve system with ILU-0
-                if( param_.cpr_use_bicgstab_ ) {
-                    Dune::BiCGSTABSolver<X> linsolve(*opAe_, *sp, (*precond_), tolerance, maxit, verbosity);
-                    linsolve.apply(x, de, result);
-                }
-                else {
-                    Dune::CGSolver<X> linsolve(*opAe_, *sp, (*precond_), tolerance, maxit, verbosity);
-                    linsolve.apply(x, de, result);
-                }
-
-            }
-
-            if (!result.converged) {
-                OPM_THROW(LinearSolverProblem, "CPRPreconditioner failed to solve elliptic subsystem.");
-            }
-        }
-
-        //! \brief Parameter collection for CPR
-        const CPRParameter& param_;
-
-        //! \brief The matrix for the full linear problem.
-        const matrix_type& A_;
-        //! \brief The elliptic part of the matrix.
-        const matrix_type& Ae_;
-
-        //! \brief temporary variables for elliptic solve
-        Y de_, ve_, dmodified_;
-
-        //! \brief elliptic operator
-        std::unique_ptr<Operator> opAe_;
-
-        //! \brief ILU0 preconditioner for the elliptic system
-        EllipticPreconditionerPointer precond_;
-        //! \brief AMG preconditioner with ILU0 smoother
-        std::unique_ptr< AMG > amg_;
-
-        //! \brief The preconditioner for the whole system
-        //!
-        //! We have to use a shared_ptr instead of a unique_ptr
-        //! as we need to use a custom allocator based on dynamic
-        //! information. But for unique_ptr the type of this deleter
-        //! has to be available at coompile time.
-        std::shared_ptr< WholeSystemPreconditioner > pre_;
-
-        //! \brief temporary variables for ILU solve
-        Y vilu_;
-
-        //! \brief The information about the parallelization of the whole system.
-        const P& comm_;
-        //! \brief The information about the parallelization of the elliptic part
-        //! of the system
-        const P& commAe_;
-     protected:
-        void createEllipticPreconditioner( const bool amg, const P& comm )
-        {
-            if( amg )
-            {
-                ISTLUtility::createAMGPreconditionerPointer( *opAe_ , param_.cpr_relax_, param_.cpr_ilu_milu_, comm, amg_ );
-            }
-            else
-            {
-                precond_ = ISTLUtility::createEllipticPreconditionerPointer<M,X>( Ae_, param_.cpr_relax_, param_.cpr_ilu_milu_, comm);
-            }
-       }
-    };
-
-
-} // namespace Opm
-
-#endif // OPM_CPRPRECONDITIONER_HEADER_INCLUDED

opm/simulators/linalg/FlowLinearSolverParameters.hpp

@@ -94,34 +94,14 @@ struct LinearSolverIgnoreConvergenceFailure{
     using type = UndefinedProperty;
 };
 template<class TypeTag, class MyTypeTag>
-struct UseAmg {
-    using type = UndefinedProperty;
-};
-template<class TypeTag, class MyTypeTag>
-struct UseCpr {
-    using type = UndefinedProperty;
-};
-template<class TypeTag, class MyTypeTag>
 struct PreconditionerAddWellContributions {
     using type = UndefinedProperty;
 };
 template<class TypeTag, class MyTypeTag>
-struct SystemStrategy {
-    using type = UndefinedProperty;
-};
-template<class TypeTag, class MyTypeTag>
 struct ScaleLinearSystem {
     using type = UndefinedProperty;
 };
 template<class TypeTag, class MyTypeTag>
-struct CprSolverVerbose {
-    using type = UndefinedProperty;
-};
-template<class TypeTag, class MyTypeTag>
-struct CprUseDrs {
-    using type = UndefinedProperty;
-};
-template<class TypeTag, class MyTypeTag>
 struct CprMaxEllIter {
     using type = UndefinedProperty;
 };
@@ -205,14 +185,6 @@ struct LinearSolverIgnoreConvergenceFailure<TypeTag, TTag::FlowIstlSolverParams>
     static constexpr bool value = false;
 };
 template<class TypeTag>
-struct UseAmg<TypeTag, TTag::FlowIstlSolverParams> {
-    static constexpr bool value = false;
-};
-template<class TypeTag>
-struct UseCpr<TypeTag, TTag::FlowIstlSolverParams> {
-    static constexpr bool value = false;
-};
-template<class TypeTag>
 struct LinearSolverBackend<TypeTag, TTag::FlowIstlSolverParams> {
     using type = Opm::ISTLSolverEbos<TypeTag>;
 };
@@ -221,22 +193,10 @@ struct PreconditionerAddWellContributions<TypeTag, TTag::FlowIstlSolverParams> {
     static constexpr bool value = false;
 };
 template<class TypeTag>
-struct SystemStrategy<TypeTag, TTag::FlowIstlSolverParams> {
-    static constexpr auto value = "none";
-};
-template<class TypeTag>
 struct ScaleLinearSystem<TypeTag, TTag::FlowIstlSolverParams> {
     static constexpr bool value = false;
 };
 template<class TypeTag>
-struct CprSolverVerbose<TypeTag, TTag::FlowIstlSolverParams> {
-    static constexpr int value = 0;
-};
-template<class TypeTag>
-struct CprUseDrs<TypeTag, TTag::FlowIstlSolverParams> {
-    static constexpr bool value = false;
-};
-template<class TypeTag>
 struct CprMaxEllIter<TypeTag, TTag::FlowIstlSolverParams> {
     static constexpr int value = 20;
 };
@@ -277,48 +237,48 @@ namespace Opm
     /**
      * \brief Parameters used to configure the CPRPreconditioner.
      */
-    struct CPRParameter
+    // struct CPRParameter
-    {
+    // {
-        double cpr_relax_;
+    //     double cpr_relax_;
-        double cpr_solver_tol_;
+    //     double cpr_solver_tol_;
-        int cpr_ilu_n_;
+    //     int cpr_ilu_n_;
-        MILU_VARIANT cpr_ilu_milu_;
+    //     MILU_VARIANT cpr_ilu_milu_;
-        bool cpr_ilu_redblack_;
+    //     bool cpr_ilu_redblack_;
-        bool cpr_ilu_reorder_sphere_;
+    //     bool cpr_ilu_reorder_sphere_;
-        bool cpr_use_drs_;
+    //     bool cpr_use_drs_;
-        int cpr_max_ell_iter_;
+    //     int cpr_max_ell_iter_;
-        int cpr_ell_solvetype_;
+    //     int cpr_ell_solvetype_;
-        bool cpr_use_amg_;
+    //     bool cpr_use_amg_;
-        bool cpr_use_bicgstab_;
+    //     bool cpr_use_bicgstab_;
-        int cpr_solver_verbose_;
+    //     int cpr_solver_verbose_;
-        bool cpr_pressure_aggregation_;
+    //     bool cpr_pressure_aggregation_;
-        int cpr_reuse_setup_;
+    //     int cpr_reuse_setup_;
-        CPRParameter() { reset(); }
+    //     CPRParameter() { reset(); }
 
-        void reset()
+    //     void reset()
-        {
+    //     {
-            cpr_solver_tol_           = 1e-2;
+    //         cpr_solver_tol_           = 1e-2;
-            cpr_ilu_n_                = 0;
+    //         cpr_ilu_n_                = 0;
-            cpr_ilu_milu_             = MILU_VARIANT::ILU;
+    //         cpr_ilu_milu_             = MILU_VARIANT::ILU;
-            cpr_ilu_redblack_         = false;
+    //         cpr_ilu_redblack_         = false;
-            cpr_ilu_reorder_sphere_   = true;
+    //         cpr_ilu_reorder_sphere_   = true;
-            cpr_max_ell_iter_         = 25;
+    //         cpr_max_ell_iter_         = 25;
-            cpr_ell_solvetype_        = 0;
+    //         cpr_ell_solvetype_        = 0;
-            cpr_use_drs_              = false;
+    //         cpr_use_drs_              = false;
-            cpr_use_amg_              = true;
+    //         cpr_use_amg_              = true;
-            cpr_use_bicgstab_         = true;
+    //         cpr_use_bicgstab_         = true;
-            cpr_solver_verbose_       = 0;
+    //         cpr_solver_verbose_       = 0;
-            cpr_pressure_aggregation_ = false;
+    //         cpr_pressure_aggregation_ = false;
-            cpr_reuse_setup_          = 0;
+    //         cpr_reuse_setup_          = 0;
-        }
+    //     }
-    };
+    // };
 
 
 
 
     /// This class carries all parameters for the NewtonIterationBlackoilInterleaved class.
     struct FlowLinearSolverParameters
-        : public CPRParameter
+        //: public CPRParameter
     {
         double linear_solver_reduction_;
         double ilu_relaxation_;
@@ -332,15 +292,18 @@ namespace Opm
         bool   newton_use_gmres_;
         bool   require_full_sparsity_pattern_;
         bool   ignoreConvergenceFailure_;
-        bool   linear_solver_use_amg_;
+        // bool   linear_solver_use_amg_;
-        bool   use_cpr_;
+        // bool   use_cpr_;
-        std::string system_strategy_;
+        // std::string system_strategy_;
         bool scale_linear_system_;
         std::string linear_solver_configuration_;
         std::string linear_solver_configuration_json_file_;
         std::string gpu_mode_;
         int bda_device_id_;
         int opencl_platform_id_;
+        int cpr_max_ell_iter_ = 20;
+        int cpr_reuse_setup_ = 0;
+        bool use_gpu_;
 
         template <class TypeTag>
         void init()
@@ -358,14 +321,14 @@ namespace Opm
             newton_use_gmres_ = EWOMS_GET_PARAM(TypeTag, bool, UseGmres);
             require_full_sparsity_pattern_ = EWOMS_GET_PARAM(TypeTag, bool, LinearSolverRequireFullSparsityPattern);
             ignoreConvergenceFailure_ = EWOMS_GET_PARAM(TypeTag, bool, LinearSolverIgnoreConvergenceFailure);
-            linear_solver_use_amg_ = EWOMS_GET_PARAM(TypeTag, bool, UseAmg);
+            // linear_solver_use_amg_ = EWOMS_GET_PARAM(TypeTag, bool, UseAmg);
-            use_cpr_ = EWOMS_GET_PARAM(TypeTag, bool, UseCpr);
+            // use_cpr_ = EWOMS_GET_PARAM(TypeTag, bool, UseCpr);
-            system_strategy_ = EWOMS_GET_PARAM(TypeTag, std::string, SystemStrategy);
+            // system_strategy_ = EWOMS_GET_PARAM(TypeTag, std::string, SystemStrategy);
             scale_linear_system_ = EWOMS_GET_PARAM(TypeTag, bool, ScaleLinearSystem);
-            cpr_solver_verbose_  =  EWOMS_GET_PARAM(TypeTag, int, CprSolverVerbose);
+            // cpr_solver_verbose_  =  EWOMS_GET_PARAM(TypeTag, int, CprSolverVerbose);
-            cpr_use_drs_  =  EWOMS_GET_PARAM(TypeTag, bool, CprUseDrs);
+            // cpr_use_drs_  =  EWOMS_GET_PARAM(TypeTag, bool, CprUseDrs);
             cpr_max_ell_iter_  =  EWOMS_GET_PARAM(TypeTag, int, CprMaxEllIter);
-            cpr_ell_solvetype_  =  EWOMS_GET_PARAM(TypeTag, int, CprEllSolvetype);
+            // cpr_ell_solvetype_  =  EWOMS_GET_PARAM(TypeTag, int, CprEllSolvetype);
             cpr_reuse_setup_  =  EWOMS_GET_PARAM(TypeTag, int, CprReuseSetup);
             linear_solver_configuration_ = EWOMS_GET_PARAM(TypeTag, std::string, LinearSolverConfiguration);
             linear_solver_configuration_json_file_ = EWOMS_GET_PARAM(TypeTag, std::string, LinearSolverConfigurationJsonFile);
@@ -389,14 +352,14 @@ namespace Opm
             EWOMS_REGISTER_PARAM(TypeTag, bool, UseGmres, "Use GMRES as the linear solver");
             EWOMS_REGISTER_PARAM(TypeTag, bool, LinearSolverRequireFullSparsityPattern, "Produce the full sparsity pattern for the linear solver");
             EWOMS_REGISTER_PARAM(TypeTag, bool, LinearSolverIgnoreConvergenceFailure, "Continue with the simulation like nothing happened after the linear solver did not converge");
-            EWOMS_REGISTER_PARAM(TypeTag, bool, UseAmg, "Use AMG as the linear solver's preconditioner");
+            // EWOMS_REGISTER_PARAM(TypeTag, bool, UseAmg, "Use AMG as the linear solver's preconditioner");
-            EWOMS_REGISTER_PARAM(TypeTag, bool, UseCpr, "Use CPR as the linear solver's preconditioner");
+            // EWOMS_REGISTER_PARAM(TypeTag, bool, UseCpr, "Use CPR as the linear solver's preconditioner");
-            EWOMS_REGISTER_PARAM(TypeTag, std::string, SystemStrategy, "Strategy for reformulating and scaling linear system (none: no scaling -- should not be used with CPR, original: use weights that are equivalent to no scaling -- should not be used with CPR, simple: form pressure equation as simple sum of conservation equations, quasiimpes: form pressure equation based on diagonal block, trueimpes: form pressure equation based on linearization of accumulation term)");
+            // EWOMS_REGISTER_PARAM(TypeTag, std::string, SystemStrategy, "Strategy for reformulating and scaling linear system (none: no scaling -- should not be used with CPR, original: use weights that are equivalent to no scaling -- should not be used with CPR, simple: form pressure equation as simple sum of conservation equations, quasiimpes: form pressure equation based on diagonal block, trueimpes: form pressure equation based on linearization of accumulation term)");
             EWOMS_REGISTER_PARAM(TypeTag, bool, ScaleLinearSystem, "Scale linear system according to equation scale and primary variable types");
-            EWOMS_REGISTER_PARAM(TypeTag, int, CprSolverVerbose, "Verbosity of cpr solver (0: silent, 1: print summary of inner linear solver, 2: print extensive information about inner linear solve, including setup information)");
+            // EWOMS_REGISTER_PARAM(TypeTag, int, CprSolverVerbose, "Verbosity of cpr solver (0: silent, 1: print summary of inner linear solver, 2: print extensive information about inner linear solve, including setup information)");
-            EWOMS_REGISTER_PARAM(TypeTag, bool, CprUseDrs, "Use dynamic row sum using weights");
+            // EWOMS_REGISTER_PARAM(TypeTag, bool, CprUseDrs, "Use dynamic row sum using weights");
             EWOMS_REGISTER_PARAM(TypeTag, int, CprMaxEllIter, "MaxIterations of the elliptic pressure part of the cpr solver");
-            EWOMS_REGISTER_PARAM(TypeTag, int, CprEllSolvetype, "Solver type of elliptic pressure solve (0: bicgstab, 1: cg, 2: only amg preconditioner)");
+            // EWOMS_REGISTER_PARAM(TypeTag, int, CprEllSolvetype, "Solver type of elliptic pressure solve (0: bicgstab, 1: cg, 2: only amg preconditioner)");
             EWOMS_REGISTER_PARAM(TypeTag, int, CprReuseSetup, "Reuse Amg Setup");
             EWOMS_REGISTER_PARAM(TypeTag, std::string, LinearSolverConfiguration, "Configuration of solver valid is: ilu0 (default), cpr_quasiimpes, cpr_trueimpes or file (specified in LinearSolverConfigurationJsonFile) ");
             EWOMS_REGISTER_PARAM(TypeTag, std::string, LinearSolverConfigurationJsonFile, "Filename of JSON configuration for flexible linear solver system.");
@@ -410,7 +373,7 @@ namespace Opm
         // set default values
         void reset()
         {
-            use_cpr_     = false;
+            // use_cpr_     = false;
             newton_use_gmres_        = false;
             linear_solver_reduction_ = 1e-2;
             linear_solver_maxiter_   = 150;
@@ -418,7 +381,7 @@ namespace Opm
             linear_solver_verbosity_ = 0;
             require_full_sparsity_pattern_ = false;
             ignoreConvergenceFailure_ = false;
-            linear_solver_use_amg_    = false;
+            // linear_solver_use_amg_    = false;
             ilu_fillin_level_         = 0;
             ilu_relaxation_           = 0.9;
             ilu_milu_                 = MILU_VARIANT::ILU;

opm/simulators/linalg/setupPropertyTree_impl.hpp

@@ -50,18 +50,13 @@ setupPropertyTree(const FlowLinearSolverParameters& p)
 #else
        OPM_THROW(std::invalid_argument,
                  "--linear-solver-configuration=file not supported with "
-                 << "boost version. Needs versoin > 1.48.");
+                 << "boost version. Needs version > 1.48.");
 #endif
     }
     else
     {
         std::string conf =  p.linear_solver_configuration_;
         // Support old UseCpr if not configuration was set
-        if (!EWOMS_PARAM_IS_SET(TypeTag, std::string, LinearSolverConfiguration) && p.use_cpr_)
-        {
-            conf = "cpr_trueimpes";
-        }
-
         if((conf == "cpr_trueimpes") || (conf == "cpr_quasiimpes")){
             prm.put("tol", p.linear_solver_reduction_);
             if (EWOMS_PARAM_IS_SET(TypeTag, int, LinearSolverMaxIter))

tests/test_blackoil_amg.cpp

@@ -21,13 +21,15 @@
 #define BOOST_TEST_MODULE BlackoilAmgTest
 #define BOOST_TEST_NO_MAIN
 #include <boost/test/unit_test.hpp>
-#include <opm/simulators/linalg/BlackoilAmg.hpp>
+#include <opm/simulators/linalg/PreconditionerFactory.hpp>
+#include <opm/simulators/linalg/FlexibleSolver.hpp>
 
 #include <dune/common/parallel/mpihelper.hh>
 #include <dune/common/fmatrix.hh>
 #include <dune/common/unused.hh>
 #include <dune/common/parallel/indexset.hh>
 #include <dune/common/parallel/plocalindex.hh>
+#include <dune/common/version.hh>
 #if DUNE_VERSION_NEWER(DUNE_COMMON, 2, 7)
 #include <dune/common/parallel/communication.hh>
 #else
@@ -35,6 +37,7 @@
 #endif
 #include <dune/istl/bcrsmatrix.hh>
 #include <dune/istl/owneroverlapcopy.hh>
+#include <dune/istl/schwarz.hh>
 
 class MPIError {
 public:
@@ -260,6 +263,7 @@ M setupAnisotropic2d(int N, Dune::ParallelIndexSet<G,L,s>& indices, const Dune::
     return mat;
 }
 
+
 //BOOST_AUTO_TEST_CASE(runBlackoilAmgLaplace)
 void runBlackoilAmgLaplace()
 {
@@ -289,30 +293,18 @@ void runBlackoilAmgLaplace()
     setBoundary(x, b, N, comm.indexSet());
 
     Operator fop(mat, comm);
-    typedef Dune::Amg::CoarsenCriterion<Dune::Amg::SymmetricCriterion<BCRSMat,Dune::Amg::FirstDiagonal> >
-        Criterion;
-    typedef Dune::SeqSSOR<BCRSMat,Vector,Vector> Smoother;
-    //typedef Dune::SeqJac<BCRSMat,Vector,Vector> Smoother;
-    //typedef Dune::SeqILU0<BCRSMat,Vector,Vector> Smoother;
-    //typedef Dune::SeqILUn<BCRSMat,Vector,Vector> Smoother;
-    typedef Dune::BlockPreconditioner<Vector,Vector,Communication,Smoother> ParSmoother;
-    typedef typename Dune::Amg::SmootherTraits<ParSmoother>::Arguments SmootherArgs;
-
     Dune::OverlappingSchwarzScalarProduct<Vector,Communication> sp(comm);
-
     Dune::InverseOperatorResult r;
-    SmootherArgs smootherArgs;
-    Criterion criterion;
 
-    smootherArgs.iterations = 1;
+    boost::property_tree::ptree prm;
-    Opm::CPRParameter param;
+    prm.put("type", "amg");
+    std::function<Vector()> weights = [&mat]() {
+        return Opm::Amg::getQuasiImpesWeights<BCRSMat, Vector>(mat, 0, false);
+    };
+    auto amg = Opm::PreconditionerFactory<Operator, Communication>::create(fop, prm, weights, comm);
+
+    Dune::CGSolver<Vector> amgCG(fop, sp, *amg, 10e-8, 300, (ccomm.rank()==0) ? 2 : 0);
 
-    Opm::BlackoilAmg<Operator,ParSmoother,Criterion,Communication,0,0> amg(param,
-                                                                         {},
-                                                                         fop, criterion,
-                                                                         smootherArgs,
-                                                                         comm);
-    Dune::CGSolver<Vector> amgCG(fop, sp, amg, 10e-8, 300, (ccomm.rank()==0) ? 2 : 0);
     amgCG.apply(x,b,r);
 
 }