24 #ifndef OPM_NEWTONITERATIONBLACKOILCPR_HEADER_INCLUDED 25 #define OPM_NEWTONITERATIONBLACKOILCPR_HEADER_INCLUDED 27 #include <opm/autodiff/DuneMatrix.hpp> 28 #include <opm/autodiff/NewtonIterationBlackoilInterface.hpp> 29 #include <opm/autodiff/CPRPreconditioner.hpp> 30 #include <opm/core/utility/parameters/ParameterGroup.hpp> 31 #include <opm/core/linalg/LinearSolverInterface.hpp> 32 #include <dune/istl/scalarproducts.hh> 33 #include <dune/istl/operators.hh> 34 #include <dune/istl/bvector.hh> 48 typedef Dune::FieldVector<double, 1 > VectorBlockType;
49 typedef Dune::FieldMatrix<double, 1, 1> MatrixBlockType;
50 typedef Dune::BCRSMatrix <MatrixBlockType> Mat;
51 typedef Dune::BlockVector<VectorBlockType> Vector;
87 template<
int category=Dune::SolverCategory::sequential,
class O,
class P>
88 void constructPreconditionerAndSolve(O& opA,
DuneMatrix& istlAe,
89 Vector& x, Vector& istlb,
90 const P& parallelInformation_arg,
91 const P& parallelInformationAe,
92 Dune::InverseOperatorResult& result)
const 94 typedef Dune::ScalarProductChooser<Vector,P,category> ScalarProductChooser;
95 std::unique_ptr<typename ScalarProductChooser::ScalarProduct>
96 sp(ScalarProductChooser::construct(parallelInformation_arg));
100 parallelInformation_arg.copyOwnerToAll(istlb, istlb);
101 Preconditioner precond(cpr_param_, opA.getmat(), istlAe, parallelInformation_arg,
102 parallelInformationAe);
107 if ( newton_use_gmres_ ) {
108 Dune::RestartedGMResSolver<Vector> linsolve(opA, *sp, precond,
109 linear_solver_reduction_, linear_solver_restart_, linear_solver_maxiter_, linear_solver_verbosity_);
111 linsolve.apply(x, istlb, result);
114 Dune::BiCGSTABSolver<Vector> linsolve(opA, *sp, precond,
115 linear_solver_reduction_, linear_solver_maxiter_, linear_solver_verbosity_);
117 linsolve.apply(x, istlb, result);
121 CPRParameter cpr_param_;
123 mutable int iterations_;
124 boost::any parallelInformation_;
126 const bool newton_use_gmres_;
127 const double linear_solver_reduction_;
128 const int linear_solver_maxiter_;
129 const int linear_solver_restart_;
130 const int linear_solver_verbosity_;
131 const bool linear_solver_ignoreconvergencefailure_;
136 #endif // OPM_NEWTONITERATIONBLACKOILCPR_HEADER_INCLUDED LinearisedBlackoilResidual::ADB::V SolutionVector
Return type for linearSolve(). A simple, non-ad vector type.
Definition: NewtonIterationBlackoilInterface.hpp:35
NewtonIterationBlackoilCPR(const ParameterGroup ¶m, const boost::any ¶llelInformation=boost::any())
Construct a system solver.
Definition: NewtonIterationBlackoilCPR.cpp:56
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: AdditionalObjectDeleter.hpp:22
Definition: DuneMatrix.hpp:43
CPR preconditioner.
Definition: CPRPreconditioner.hpp:370
This class solves the fully implicit black-oil system by applying a Constrained Pressure Residual pre...
Definition: NewtonIterationBlackoilCPR.hpp:46
Residual structure of the fully implicit solver.
Definition: LinearisedBlackoilResidual.hpp:47
virtual SolutionVector computeNewtonIncrement(const LinearisedBlackoilResidual &residual) const
Solve the system of linear equations Ax = b, with A being the combined derivative matrix of the resid...
Definition: NewtonIterationBlackoilCPR.cpp:80
Interface class for (linear) solvers for the fully implicit black-oil system.
Definition: NewtonIterationBlackoilInterface.hpp:31
virtual const boost::any & parallelInformation() const
Get the information about the parallelization of the grid.
Definition: NewtonIterationBlackoilCPR.cpp:191
virtual int iterations() const
Definition: NewtonIterationBlackoilCPR.hpp:77