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>
37#include <opm/simulators/timestepping/SimulatorReport.hpp>
38#include <opm/simulators/timestepping/SimulatorTimerInterface.hpp>
42namespace Opm::Parameters {
55enum class NonlinearRelaxType {
64void detectOscillations(
const std::vector<std::vector<Scalar>>&
residualHistory,
65 const int it,
const int numPhases,
const Scalar relaxRelTol,
71template <
class BVector,
class Scalar>
72void stabilizeNonlinearUpdate(BVector&
dx, BVector&
dxOld,
73 const Scalar
omega, NonlinearRelaxType relaxType);
76void registerNonlinearParameters();
82 template <
class TypeTag,
class PhysicalModel>
91 NonlinearRelaxType relaxType_;
93 Scalar relaxIncrement_;
104 relaxMax_ = Parameters::Get<Parameters::NewtonMaxRelax<Scalar>>();
105 maxIter_ = Parameters::Get<Parameters::NewtonMaxIterations>();
106 minIter_ = Parameters::Get<Parameters::NewtonMinIterations>();
110 relaxType_ = NonlinearRelaxType::Dampen;
112 relaxType_ = NonlinearRelaxType::SOR;
119 static void registerParameters()
121 detail::registerNonlinearParameters<Scalar>();
127 relaxType_ = NonlinearRelaxType::Dampen;
129 relaxIncrement_ = 0.1;
147 std::unique_ptr<PhysicalModel>
model)
151 , nonlinearIterations_(0)
152 , linearIterations_(0)
154 , nonlinearIterationsLast_(0)
155 , linearIterationsLast_(0)
156 , wellIterationsLast_(0)
159 OPM_THROW(std::logic_error,
"Must provide a non-null model argument for NonlinearSolver.");
171 report += model_->prepareStep(timer);
178 bool converged =
false;
191 converged = report.converged;
197 failureReport_ = report;
198 failureReport_ += model_->failureReport();
205 failureReport_ = report;
207 std::string
msg =
"Solver convergence failure - Failed to complete a time step within " + std::to_string(
maxIter()) +
" iterations.";
212 report += model_->afterStep(timer);
213 report.converged =
true;
219 {
return failureReport_; }
223 {
return linearizations_; }
227 {
return nonlinearIterations_; }
231 {
return linearIterations_; }
235 {
return wellIterations_; }
239 {
return nonlinearIterationsLast_; }
243 {
return linearIterationsLast_; }
247 {
return wellIterationsLast_; }
249 std::vector<std::vector<Scalar> >
250 computeFluidInPlace(
const std::vector<int>&
fipnum)
const
251 {
return model_->computeFluidInPlace(
fipnum); }
272 template <
class BVector>
280 {
return param_.relaxMax_; }
284 {
return param_.relaxIncrement_; }
288 {
return param_.relaxType_; }
292 {
return param_.relaxRelTol_; }
296 {
return param_.maxIter_; }
300 {
return param_.minIter_; }
309 SolverParameters param_;
310 std::unique_ptr<PhysicalModel> model_;
312 int nonlinearIterations_;
313 int linearIterations_;
315 int nonlinearIterationsLast_;
316 int linearIterationsLast_;
317 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:84
int linearIterationsLastStep() const
Number of linear solver iterations used in the last call to step().
Definition NonlinearSolver.hpp:242
void stabilizeNonlinearUpdate(BVector &dx, BVector &dxOld, const Scalar omega) const
Apply a stabilization to dx, depending on dxOld and relaxation parameters.
Definition NonlinearSolver.hpp:273
Scalar relaxRelTol() const
The relaxation relative tolerance.
Definition NonlinearSolver.hpp:291
Scalar relaxMax() const
The greatest relaxation factor (i.e. smallest factor) allowed.
Definition NonlinearSolver.hpp:279
int linearizations() const
Number of linearizations used in all calls to step().
Definition NonlinearSolver.hpp:222
int wellIterationsLastStep() const
Number of well iterations used in all calls to step().
Definition NonlinearSolver.hpp:246
int nonlinearIterations() const
Number of full nonlinear solver iterations used in all calls to step().
Definition NonlinearSolver.hpp:226
void setParameters(const SolverParameters ¶m)
Set parameters to override those given at construction time.
Definition NonlinearSolver.hpp:303
NonlinearRelaxType relaxType() const
The relaxation type (Dampen or SOR).
Definition NonlinearSolver.hpp:287
int minIter() const
The minimum number of nonlinear iterations allowed.
Definition NonlinearSolver.hpp:299
int linearIterations() const
Number of linear solver iterations used in all calls to step().
Definition NonlinearSolver.hpp:230
int nonlinearIterationsLastStep() const
Number of nonlinear solver iterations used in the last call to step().
Definition NonlinearSolver.hpp:238
int maxIter() const
The maximum number of nonlinear iterations allowed.
Definition NonlinearSolver.hpp:295
const SimulatorReportSingle & failureReport() const
return the statistics if the step() method failed
Definition NonlinearSolver.hpp:218
int wellIterations() const
Number of well iterations used in all calls to step().
Definition NonlinearSolver.hpp:234
const PhysicalModel & model() const
Reference to physical model.
Definition NonlinearSolver.hpp:254
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:262
NonlinearSolver(const SolverParameters ¶m, std::unique_ptr< PhysicalModel > model)
Construct solver for a given model.
Definition NonlinearSolver.hpp:146
Scalar relaxIncrement() const
The step-change size for the relaxation factor.
Definition NonlinearSolver.hpp:283
PhysicalModel & model()
Mutable reference to physical model.
Definition NonlinearSolver.hpp:258
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
This file provides the infrastructure to retrieve run-time parameters.
The Opm property system, traits with inheritance.
Definition NonlinearSolver.hpp:90
Definition NonlinearSolver.hpp:45
Definition NonlinearSolver.hpp:47
Definition NonlinearSolver.hpp:48
A struct for returning timing data from a simulator to its caller.
Definition SimulatorReport.hpp:34