21 #ifndef OPM_NONLINEARSOLVER_HEADER_INCLUDED
22 #define OPM_NONLINEARSOLVER_HEADER_INCLUDED
24 #include <opm/core/simulator/SimulatorReport.hpp>
25 #include <opm/core/utility/parameters/ParameterGroup.hpp>
26 #include <opm/simulators/timestepping/SimulatorTimerInterface.hpp>
27 #include <dune/common/fmatrix.hh>
28 #include <dune/istl/bcrsmatrix.hh>
36 template <
class PhysicalModel>
41 enum RelaxType { DAMPEN, SOR };
46 enum RelaxType relax_type_;
48 double relax_increment_;
49 double relax_rel_tol_;
60 typedef typename PhysicalModel::ReservoirState ReservoirState;
61 typedef typename PhysicalModel::WellState WellState;
73 std::unique_ptr<PhysicalModel>
model);
82 ReservoirState& reservoir_state,
83 WellState& well_state);
97 const ReservoirState& initial_reservoir_state,
98 const WellState& initial_well_state,
99 ReservoirState& reservoir_state,
100 WellState& well_state);
104 {
return failureReport_; }
131 std::vector<std::vector<double> >
134 return model_->computeFluidInPlace(x, fipnum);
137 std::vector<std::vector<double> >
140 return model_->computeFluidInPlace(fipnum);
145 const PhysicalModel&
model()
const;
148 PhysicalModel&
model();
152 const int it,
bool& oscillate,
bool& stagnate)
const;
156 template <
class BVector>
160 double relaxMax()
const {
return param_.relax_max_; }
166 enum RelaxType
relaxType()
const {
return param_.relax_type_; }
172 int maxIter()
const {
return param_.max_iter_; }
175 int minIter()
const {
return param_.min_iter_; }
182 SimulatorReport failureReport_;
183 SolverParameters param_;
184 std::unique_ptr<PhysicalModel> model_;
186 int nonlinearIterations_;
187 int linearIterations_;
189 int nonlinearIterationsLast_;
190 int linearIterationsLast_;
191 int wellIterationsLast_;
195 #include "NonlinearSolver_impl.hpp"
197 #endif // OPM_NONLINEARSOLVER_HEADER_INCLUDED
Definition: NonlinearSolver.hpp:44
int maxIter() const
The maximum number of nonlinear iterations allowed.
Definition: NonlinearSolver.hpp:172
int linearIterationsLastStep() const
Number of linear solver iterations used in the last call to step().
Definition: NonlinearSolver_impl.hpp:93
A nonlinear solver class suitable for general fully-implicit models, as well as pressure, transport and sequential models.
Definition: NonlinearSolver.hpp:37
double relaxMax() const
The greatest relaxation factor (i.e. smallest factor) allowed.
Definition: NonlinearSolver.hpp:160
int linearIterations() const
Number of linear solver iterations used in all calls to step().
Definition: NonlinearSolver_impl.hpp:63
void stabilizeNonlinearUpdate(BVector &dx, BVector &dxOld, const double omega) const
Apply a stabilization to dx, depending on dxOld and relaxation parameters.
Definition: NonlinearSolver_impl.hpp:261
NonlinearSolver(const SolverParameters ¶m, std::unique_ptr< PhysicalModel > model)
Construct solver for a given model.
Definition: NonlinearSolver_impl.hpp:33
int nonlinearIterations() const
Number of full nonlinear solver iterations used in all calls to step().
Definition: NonlinearSolver_impl.hpp:57
const PhysicalModel & model() const
Reference to physical model.
Definition: NonlinearSolver_impl.hpp:75
const SimulatorReport & failureReport() const
return the statistics if the step() method failed
Definition: NonlinearSolver.hpp:103
int wellIterationsLastStep() const
Number of well iterations used in all calls to step().
Definition: NonlinearSolver_impl.hpp:99
int nonlinearIterationsLastStep() const
Number of nonlinear solver iterations used in the last call to step().
Definition: NonlinearSolver_impl.hpp:87
int linearizations() const
Number of linearizations used in all calls to step().
Definition: NonlinearSolver_impl.hpp:51
SimulatorReport step(const SimulatorTimerInterface &timer, ReservoirState &reservoir_state, WellState &well_state)
Take a single forward step, after which the states will be modified according to the physical model...
Definition: NonlinearSolver_impl.hpp:107
Interface class for SimulatorTimer objects, to be improved.
Definition: SimulatorTimerInterface.hpp:35
std::vector< std::vector< double > > computeFluidInPlace(const ReservoirState &x, const std::vector< int > &fipnum) const
Compute fluid in place.
Definition: NonlinearSolver.hpp:132
int wellIterations() const
Number of well iterations used in all calls to step().
Definition: NonlinearSolver_impl.hpp:69
void setParameters(const SolverParameters ¶m)
Set parameters to override those given at construction time.
Definition: NonlinearSolver.hpp:178
double relaxRelTol() const
The relaxation relative tolerance.
Definition: NonlinearSolver.hpp:169
int minIter() const
The minimum number of nonlinear iterations allowed.
Definition: NonlinearSolver.hpp:175
void detectOscillations(const std::vector< std::vector< double >> &residual_history, const int it, bool &oscillate, bool &stagnate) const
Detect oscillation or stagnation in a given residual history.
Definition: NonlinearSolver_impl.hpp:223
double relaxIncrement() const
The step-change size for the relaxation factor.
Definition: NonlinearSolver.hpp:163
enum RelaxType relaxType() const
The relaxation type (DAMPEN or SOR).
Definition: NonlinearSolver.hpp:166