Updated after review comments.

opm/autodiff/BlackoilAmg.hpp

@@ -426,7 +426,7 @@ private:
             }
             // Linear solver parameters
             const double tolerance = param_->cpr_solver_tol_;
-            const int maxit        = param_->cpr_max_iter_;
+            const int maxit        = param_->cpr_max_ell_iter_;
             int verbosity = 0;
             if (comm_.communicator().rank() == 0) {
                 verbosity = param_->cpr_solver_verbose_;
@@ -507,6 +507,12 @@ private:
 #endif
             }
 
+            // Warn if unknown options.
+            if (param_->cpr_ell_solvetype_ > 2 && comm_.communicator().rank() == 0) {
+                OpmLog::warning("cpr_ell_solver_type_unknown", "Unknown CPR elliptic solver type specification, using LoopSolver.");
+            }
+
+
 #if ! DUNE_VERSION_NEWER(DUNE_ISTL, 2, 6)
             delete sp;
 #endif
@@ -926,16 +932,17 @@ private:
  * \brief An algebraic twolevel or multigrid approach for solving blackoil (supports CPR with and without AMG)
  *
  * This preconditioner first decouples the component used for coarsening using a simple scaling
- * approach (e.g. Scheichl, Masson 2013,\see scaleMatrixDRS). Then it constructs the first
+ * approach (e.g. Scheichl, Masson 2013,\see scaleMatrixDRS). Then it constructs the
- * coarse level system, either by simply extracting the coupling between the components at COMPONENT_INDEX
+ * coarse level system. The coupling is defined by the weights corresponding to the element located at
- * in the matrix blocks or by extracting them and applying aggregation to them directly. This coarse level
+ * (COMPONENT_INDEX, VARIABLE_INDEX) in the block matrix. Then the coarse level system is constructed
+ * either by extracting these elements, or by applying aggregation to them directly. This coarse level
  * can be solved either by AMG or by ILU. The preconditioner is configured using CPRParameter.
  * \tparam O The type of the operator (encapsulating a BCRSMatrix).
  * \tparam S The type of the smoother.
  * \tparam C The type of coarsening criterion to use.
  * \tparam P The type of the class describing the parallelization.
- * \tparam COMPONENT_INDEX The index of the component to use for coarsening (usually the pressure).
+ * \tparam COMPONENT_INDEX The index of the component to use for coarsening (usually water).
- * \tparam VARIABLE_INDEX The index of the variable to use for coarsening (usually the pressure).
+ * \tparam VARIABLE_INDEX The index of the variable to use for coarsening (usually pressure).
  */
 template<typename O, typename S, typename C,
          typename P, std::size_t COMPONENT_INDEX, std::size_t VARIABLE_INDEX>
@@ -1005,9 +1012,9 @@ public:
         : param_(param),
           weights_(weights),
           scaledMatrixOperator_(Detail::scaleMatrixDRS(fineOperator, comm,
-                                COMPONENT_INDEX, weights, param)),
+                                                       COMPONENT_INDEX, weights, param)),
           smoother_(Detail::constructSmoother<Smoother>(std::get<1>(scaledMatrixOperator_),
-                    smargs, comm)),
+                                                        smargs, comm)),
           levelTransferPolicy_(criterion, comm, param.cpr_pressure_aggregation_),
           coarseSolverPolicy_(&param, smargs, criterion),
           twoLevelMethod_(std::get<1>(scaledMatrixOperator_), smoother_,

opm/autodiff/FlowLinearSolverParameters.hpp

@@ -63,7 +63,7 @@ NEW_PROP_TAG(SystemStrategy);
 NEW_PROP_TAG(ScaleLinearSystem);
 NEW_PROP_TAG(CprSolverVerbose);
 NEW_PROP_TAG(CprUseDrs);
-NEW_PROP_TAG(CprMaxIter);
+NEW_PROP_TAG(CprMaxEllIter);
 NEW_PROP_TAG(CprEllSolvetype);
 NEW_PROP_TAG(CprReuseSetup);
 
@@ -83,11 +83,11 @@ SET_BOOL_PROP(FlowIstlSolverParams, UseAmg, false);
 SET_BOOL_PROP(FlowIstlSolverParams, UseCpr, false);
 SET_TYPE_PROP(FlowIstlSolverParams, LinearSolverBackend, Opm::ISTLSolverEbos<TypeTag>);
 SET_BOOL_PROP(FlowIstlSolverParams, PreconditionerAddWellContributions, false);
-SET_STRING_PROP(FlowIstlSolverParams, SystemStrategy, "original");
+SET_STRING_PROP(FlowIstlSolverParams, SystemStrategy, "none");
 SET_BOOL_PROP(FlowIstlSolverParams, ScaleLinearSystem, false);
 SET_INT_PROP(FlowIstlSolverParams, CprSolverVerbose, 0);
 SET_BOOL_PROP(FlowIstlSolverParams, CprUseDrs, false);
-SET_INT_PROP(FlowIstlSolverParams, CprMaxIter, 20);
+SET_INT_PROP(FlowIstlSolverParams, CprMaxEllIter, 20);
 SET_INT_PROP(FlowIstlSolverParams, CprEllSolvetype, 0);
 SET_INT_PROP(FlowIstlSolverParams, CprReuseSetup, 0);
 
@@ -111,7 +111,6 @@ namespace Opm
         bool cpr_ilu_reorder_sphere_;
         bool cpr_use_drs_;
         int cpr_max_ell_iter_;
-        int cpr_max_iter_;
         int cpr_ell_solvetype_;
         bool cpr_use_amg_;
         bool cpr_use_bicgstab_;
@@ -130,7 +129,6 @@ namespace Opm
             cpr_max_ell_iter_         = 25;
             cpr_ell_solvetype_        = 0;
             cpr_use_drs_              = false;
-            cpr_max_iter_             = 25;
             cpr_use_amg_              = true;
             cpr_use_bicgstab_         = true;
             cpr_solver_verbose_       = 0;
@@ -185,7 +183,7 @@ namespace Opm
             scale_linear_system_ = EWOMS_GET_PARAM(TypeTag, bool, ScaleLinearSystem);
             cpr_solver_verbose_  =  EWOMS_GET_PARAM(TypeTag, int, CprSolverVerbose);
             cpr_use_drs_  =  EWOMS_GET_PARAM(TypeTag, bool, CprUseDrs);
-            cpr_max_iter_  =  EWOMS_GET_PARAM(TypeTag, int, CprMaxIter);
+            cpr_max_ell_iter_  =  EWOMS_GET_PARAM(TypeTag, int, CprMaxEllIter);
             cpr_ell_solvetype_  =  EWOMS_GET_PARAM(TypeTag, int, CprEllSolvetype);
             cpr_reuse_setup_  =  EWOMS_GET_PARAM(TypeTag, int, CprReuseSetup);
         }
@@ -207,12 +205,12 @@ namespace Opm
             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, UseCpr, "Use CPR as the linear solver's preconditioner");
-            EWOMS_REGISTER_PARAM(TypeTag, std::string, SystemStrategy, "Strategy for reformulating and scale linear system");
+            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, "Verbose for cpr solver");
+            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 weighs");
+            EWOMS_REGISTER_PARAM(TypeTag, bool, CprUseDrs, "Use dynamic row sum using weights");
-            EWOMS_REGISTER_PARAM(TypeTag, int, CprMaxIter, "MaxIterations of the pressure amg solver");
+            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 solve 0 bicgstab 1 cg other 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");
         }
 

opm/autodiff/ISTLSolverEbos.hpp

@@ -249,13 +249,16 @@ protected:
                     bvec[pressureEqnIndex] = 1;
                     weights_ = getSimpleWeights(bvec);
                 } else {
+                    if (parameters_.system_strategy_ != "none") {
+                        OpmLog::warning("unknown_system_strategy", "Unknown linear solver system strategy: '" + parameters_.system_strategy_ + "', applying 'none' strategy.");
+                    }
                     form_cpr = false;
                 }
                 if (parameters_.scale_linear_system_) {
                     // also scale weights
                     this->scaleEquationsAndVariables(weights_);
                 }
-                if (form_cpr && not(parameters_.cpr_use_drs_)) {
+                if (form_cpr && !(parameters_.cpr_use_drs_)) {
                     scaleMatrixAndRhs(weights_);
                 }
                 if (weights_.size() == 0) {
@@ -302,7 +305,7 @@ protected:
                 //Not sure what actual_mat_for_prec is, so put ebosJacIgnoreOverlap as both variables
                 //to be certain that correct matrix is used for preconditioning.
                 Operator opA(ebosJacIgnoreOverlap, ebosJacIgnoreOverlap, wellModel,
-                        parallelInformation_ );
+                             parallelInformation_ );
                 assert( opA.comm() );
                 solve( opA, x, *rhs_, *(opA.comm()) );
             }
@@ -642,6 +645,8 @@ protected:
         }
 
         // Weights to make approximate pressure equations.
+        // Calculated from the storage terms (only) of the
+        // conservation equations, ignoring all other terms.
         Vector getStorageWeights() const
         {
             Vector weights(rhs_->size()); 
@@ -682,6 +687,13 @@ protected:
             return weights;
         }
 
+        // Interaction between the CPR weights (the function argument 'weights')
+        // and the variable and equation weights from
+        // simulator_.model().primaryVarWeight() and
+        // simulator_.model().eqWeight() is nontrivial and does not work
+        // at the moment. Possibly refactoring of ewoms weight treatment
+        // is needed. In the meantime this function shows what needs to be
+        // done to integrate the weights properly.
         void scaleEquationsAndVariables(Vector& weights)
         {
             // loop over primary variables
@@ -702,7 +714,7 @@ protected:
                 for (std::size_t ii = 0; ii < brhs.size(); ii++) {
                     brhs[ii] *= simulator_.model().eqWeight(i.index(), ii);
                 }
-                if (weights_.size() == matrix_->N()) {
+                if (weights.size() == matrix_->N()) {
                     BlockVector& bw = weights[i.index()];
                     for (std::size_t ii = 0; ii < brhs.size(); ii++) {
                         bw[ii] /= simulator_.model().eqWeight(i.index(), ii);