21#ifndef OPM_NON_LINEAR_SOLVER_HPP
22#define OPM_NON_LINEAR_SOLVER_HPP
24#include <dune/common/fmatrix.hh>
25#include <dune/istl/bcrsmatrix.hh>
27#include <opm/common/ErrorMacros.hpp>
28#include <opm/common/Exceptions.hpp>
30#include <opm/models/nonlinear/newtonmethodparams.hpp>
31#include <opm/models/nonlinear/newtonmethodproperties.hh>
36#include <opm/simulators/timestepping/SimulatorReport.hpp>
37#include <opm/simulators/timestepping/SimulatorTimerInterface.hpp>
41namespace Opm::Parameters {
54enum class NonlinearRelaxType {
63void detectOscillations(
const std::vector<std::vector<Scalar>>&
residualHistory,
64 const int it,
const int numPhases,
const Scalar relaxRelTol,
70template <
class BVector,
class Scalar>
71void stabilizeNonlinearUpdate(BVector&
dx, BVector&
dxOld,
72 const Scalar
omega, NonlinearRelaxType relaxType);
80 NonlinearRelaxType relaxType_;
82 Scalar relaxIncrement_;
89 static void registerParameters();
96 template <
class TypeTag,
class PhysicalModel>
114 std::unique_ptr<PhysicalModel>
model)
118 , nonlinearIterations_(0)
119 , linearIterations_(0)
121 , nonlinearIterationsLast_(0)
122 , linearIterationsLast_(0)
123 , wellIterationsLast_(0)
126 OPM_THROW(std::logic_error,
"Must provide a non-null model argument for NonlinearSolver.");
138 report += model_->prepareStep(timer);
145 bool converged =
false;
158 converged = report.converged;
164 failureReport_ = report;
165 failureReport_ += model_->failureReport();
172 failureReport_ = report;
174 std::string
msg =
"Solver convergence failure - Failed to complete a time step within " + std::to_string(
maxIter()) +
" iterations.";
179 report += model_->afterStep(timer);
180 report.converged =
true;
186 {
return failureReport_; }
190 {
return linearizations_; }
194 {
return nonlinearIterations_; }
198 {
return linearIterations_; }
202 {
return wellIterations_; }
206 {
return nonlinearIterationsLast_; }
210 {
return linearIterationsLast_; }
214 {
return wellIterationsLast_; }
216 std::vector<std::vector<Scalar> >
217 computeFluidInPlace(
const std::vector<int>&
fipnum)
const
218 {
return model_->computeFluidInPlace(
fipnum); }
239 template <
class BVector>
247 {
return param_.relaxMax_; }
251 {
return param_.relaxIncrement_; }
255 {
return param_.relaxType_; }
259 {
return param_.relaxRelTol_; }
263 {
return param_.maxIter_; }
267 {
return param_.minIter_; }
276 SolverParameters param_;
277 std::unique_ptr<PhysicalModel> model_;
279 int nonlinearIterations_;
280 int linearIterations_;
282 int nonlinearIterationsLast_;
283 int linearIterationsLast_;
284 int wellIterationsLast_;
Defines a type tags and some fundamental properties all models.
A nonlinear solver class suitable for general fully-implicit models, as well as pressure,...
Definition NonlinearSolver.hpp:98
int linearIterationsLastStep() const
Number of linear solver iterations used in the last call to step().
Definition NonlinearSolver.hpp:209
void stabilizeNonlinearUpdate(BVector &dx, BVector &dxOld, const Scalar omega) const
Apply a stabilization to dx, depending on dxOld and relaxation parameters.
Definition NonlinearSolver.hpp:240
Scalar relaxRelTol() const
The relaxation relative tolerance.
Definition NonlinearSolver.hpp:258
Scalar relaxMax() const
The greatest relaxation factor (i.e. smallest factor) allowed.
Definition NonlinearSolver.hpp:246
int linearizations() const
Number of linearizations used in all calls to step().
Definition NonlinearSolver.hpp:189
int wellIterationsLastStep() const
Number of well iterations used in all calls to step().
Definition NonlinearSolver.hpp:213
int nonlinearIterations() const
Number of full nonlinear solver iterations used in all calls to step().
Definition NonlinearSolver.hpp:193
void setParameters(const SolverParameters ¶m)
Set parameters to override those given at construction time.
Definition NonlinearSolver.hpp:270
NonlinearRelaxType relaxType() const
The relaxation type (Dampen or SOR).
Definition NonlinearSolver.hpp:254
int minIter() const
The minimum number of nonlinear iterations allowed.
Definition NonlinearSolver.hpp:266
int linearIterations() const
Number of linear solver iterations used in all calls to step().
Definition NonlinearSolver.hpp:197
int nonlinearIterationsLastStep() const
Number of nonlinear solver iterations used in the last call to step().
Definition NonlinearSolver.hpp:205
int maxIter() const
The maximum number of nonlinear iterations allowed.
Definition NonlinearSolver.hpp:262
const SimulatorReportSingle & failureReport() const
return the statistics if the step() method failed
Definition NonlinearSolver.hpp:185
int wellIterations() const
Number of well iterations used in all calls to step().
Definition NonlinearSolver.hpp:201
const PhysicalModel & model() const
Reference to physical model.
Definition NonlinearSolver.hpp:221
void detectOscillations(const std::vector< std::vector< Scalar > > &residualHistory, const int it, bool &oscillate, bool &stagnate) const
Detect oscillation or stagnation in a given residual history.
Definition NonlinearSolver.hpp:229
NonlinearSolver(const SolverParameters ¶m, std::unique_ptr< PhysicalModel > model)
Construct solver for a given model.
Definition NonlinearSolver.hpp:113
Scalar relaxIncrement() const
The step-change size for the relaxation factor.
Definition NonlinearSolver.hpp:250
PhysicalModel & model()
Mutable reference to physical model.
Definition NonlinearSolver.hpp:225
Interface class for SimulatorTimer objects, to be improved.
Definition SimulatorTimerInterface.hpp:34
virtual double currentStepLength() const =0
Current step length.
virtual double simulationTimeElapsed() const =0
Time elapsed since the start of the simulation until the beginning of the current time step [s].
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242
typename Properties::Detail::GetPropImpl< TypeTag, Property >::type::type GetPropType
get the type alias defined in the property (equivalent to old macro GET_PROP_TYPE(....
Definition propertysystem.hh:235
The Opm property system, traits with inheritance.
Definition NonlinearSolver.hpp:79
Definition NonlinearSolver.hpp:44
Definition NonlinearSolver.hpp:46
Definition NonlinearSolver.hpp:47
A struct for returning timing data from a simulator to its caller.
Definition SimulatorReport.hpp:34