21#ifndef OPM_NON_LINEAR_SOLVER_HPP
22#define OPM_NON_LINEAR_SOLVER_HPP
24#include <opm/simulators/timestepping/SimulatorReport.hpp>
25#include <opm/common/ErrorMacros.hpp>
26#include <opm/simulators/timestepping/SimulatorTimerInterface.hpp>
28#include <opm/models/utils/parametersystem.hh>
29#include <opm/models/utils/propertysystem.hh>
30#include <opm/models/utils/basicproperties.hh>
31#include <opm/models/nonlinear/newtonmethodproperties.hh>
32#include <opm/common/Exceptions.hpp>
34#include <dune/common/fmatrix.hh>
35#include <dune/istl/bcrsmatrix.hh>
38namespace Opm::Properties {
44template<
class TypeTag,
class MyTypeTag>
46 using type = UndefinedProperty;
52template<
class TypeTag,
class MyTypeTag>
54 using type = UndefinedProperty;
56template<
class TypeTag,
class MyTypeTag>
58 using type = UndefinedProperty;
61template<
class TypeTag>
63 using type = GetPropType<TypeTag, Scalar>;
64 static constexpr type value = 0.5;
66template<
class TypeTag>
67struct NewtonMaxIterations<TypeTag, TTag::FlowNonLinearSolver> {
68 static constexpr int value = 20;
70template<
class TypeTag>
72 static constexpr int value = 2;
74template<
class TypeTag>
76 static constexpr auto value =
"dampen";
86enum class NonlinearRelaxType {
94void detectOscillations(
const std::vector<std::vector<double>>& residualHistory,
95 const int it,
const int numPhases,
const double relaxRelTol,
96 bool& oscillate,
bool& stagnate);
100template <
class BVector>
101void stabilizeNonlinearUpdate(BVector& dx, BVector& dxOld,
102 const double omega, NonlinearRelaxType relaxType);
108 template <
class TypeTag,
class PhysicalModel>
111 using Scalar = GetPropType<TypeTag, Properties::Scalar>;
117 NonlinearRelaxType relaxType_;
119 double relaxIncrement_;
130 relaxMax_ = Parameters::get<TypeTag, Properties::NewtonMaxRelax>();
131 maxIter_ = Parameters::get<TypeTag, Properties::NewtonMaxIterations>();
132 minIter_ = Parameters::get<TypeTag, Properties::NewtonMinIterations>();
134 const auto& relaxationTypeString = Parameters::get<TypeTag, Properties::NewtonRelaxationType>();
135 if (relaxationTypeString ==
"dampen") {
136 relaxType_ = NonlinearRelaxType::Dampen;
137 }
else if (relaxationTypeString ==
"sor") {
138 relaxType_ = NonlinearRelaxType::SOR;
140 OPM_THROW(std::runtime_error,
141 "Unknown Relaxtion Type " + relaxationTypeString);
145 static void registerParameters()
147 Parameters::registerParam<TypeTag, Properties::NewtonMaxRelax>
148 (
"The maximum relaxation factor of a Newton iteration");
149 Parameters::registerParam<TypeTag, Properties::NewtonMaxIterations>
150 (
"The maximum number of Newton iterations per time step");
151 Parameters::registerParam<TypeTag, Properties::NewtonMinIterations>
152 (
"The minimum number of Newton iterations per time step");
153 Parameters::registerParam<TypeTag, Properties::NewtonRelaxationType>
154 (
"The type of relaxation used by Newton method");
160 relaxType_ = NonlinearRelaxType::Dampen;
162 relaxIncrement_ = 0.1;
183 std::unique_ptr<PhysicalModel>
model)
185 , model_(std::move(
model))
187 , nonlinearIterations_(0)
188 , linearIterations_(0)
190 , nonlinearIterationsLast_(0)
191 , linearIterationsLast_(0)
192 , wellIterationsLast_(0)
195 OPM_THROW(std::logic_error,
"Must provide a non-null model argument for NonlinearSolver.");
207 report += model_->prepareStep(timer);
214 bool converged =
false;
222 auto iterReport = model_->nonlinearIteration(iteration, timer, *
this);
224 report += iterReport;
225 report.converged = iterReport.converged;
227 converged = report.converged;
233 failureReport_ = report;
234 failureReport_ += model_->failureReport();
238 while ( (!converged && (iteration <=
maxIter())) || (iteration <=
minIter()));
241 failureReport_ = report;
243 std::string msg =
"Solver convergence failure - Failed to complete a time step within " + std::to_string(
maxIter()) +
" iterations.";
244 OPM_THROW_NOLOG(TooManyIterations, msg);
248 report += model_->afterStep(timer);
249 report.converged =
true;
255 {
return failureReport_; }
259 {
return linearizations_; }
263 {
return nonlinearIterations_; }
267 {
return linearIterations_; }
271 {
return wellIterations_; }
275 {
return nonlinearIterationsLast_; }
279 {
return linearIterationsLast_; }
283 {
return wellIterationsLast_; }
285 std::vector<std::vector<double> >
286 computeFluidInPlace(
const std::vector<int>& fipnum)
const
287 {
return model_->computeFluidInPlace(fipnum); }
299 const int it,
bool& oscillate,
bool& stagnate)
const
301 detail::detectOscillations(residualHistory, it, model_->numPhases(),
302 this->relaxRelTol(), oscillate, stagnate);
308 template <
class BVector>
311 detail::stabilizeNonlinearUpdate(dx, dxOld, omega, this->
relaxType());
316 {
return param_.relaxMax_; }
320 {
return param_.relaxIncrement_; }
324 {
return param_.relaxType_; }
328 {
return param_.relaxRelTol_; }
332 {
return param_.maxIter_; }
336 {
return param_.minIter_; }
345 SolverParameters param_;
346 std::unique_ptr<PhysicalModel> model_;
348 int nonlinearIterations_;
349 int linearIterations_;
351 int nonlinearIterationsLast_;
352 int linearIterationsLast_;
353 int wellIterationsLast_;
A nonlinear solver class suitable for general fully-implicit models, as well as pressure,...
Definition NonlinearSolver.hpp:110
int linearIterationsLastStep() const
Number of linear solver iterations used in the last call to step().
Definition NonlinearSolver.hpp:278
int linearizations() const
Number of linearizations used in all calls to step().
Definition NonlinearSolver.hpp:258
int wellIterationsLastStep() const
Number of well iterations used in all calls to step().
Definition NonlinearSolver.hpp:282
int nonlinearIterations() const
Number of full nonlinear solver iterations used in all calls to step().
Definition NonlinearSolver.hpp:262
void setParameters(const SolverParameters ¶m)
Set parameters to override those given at construction time.
Definition NonlinearSolver.hpp:339
double relaxMax() const
The greatest relaxation factor (i.e. smallest factor) allowed.
Definition NonlinearSolver.hpp:315
void detectOscillations(const std::vector< std::vector< double > > &residualHistory, const int it, bool &oscillate, bool &stagnate) const
Detect oscillation or stagnation in a given residual history.
Definition NonlinearSolver.hpp:298
NonlinearRelaxType relaxType() const
The relaxation type (Dampen or SOR).
Definition NonlinearSolver.hpp:323
int minIter() const
The minimum number of nonlinear iterations allowed.
Definition NonlinearSolver.hpp:335
int linearIterations() const
Number of linear solver iterations used in all calls to step().
Definition NonlinearSolver.hpp:266
int nonlinearIterationsLastStep() const
Number of nonlinear solver iterations used in the last call to step().
Definition NonlinearSolver.hpp:274
int maxIter() const
The maximum number of nonlinear iterations allowed.
Definition NonlinearSolver.hpp:331
double relaxRelTol() const
The relaxation relative tolerance.
Definition NonlinearSolver.hpp:327
double relaxIncrement() const
The step-change size for the relaxation factor.
Definition NonlinearSolver.hpp:319
const SimulatorReportSingle & failureReport() const
return the statistics if the step() method failed
Definition NonlinearSolver.hpp:254
int wellIterations() const
Number of well iterations used in all calls to step().
Definition NonlinearSolver.hpp:270
const PhysicalModel & model() const
Reference to physical model.
Definition NonlinearSolver.hpp:290
void stabilizeNonlinearUpdate(BVector &dx, BVector &dxOld, const double omega) const
Apply a stabilization to dx, depending on dxOld and relaxation parameters.
Definition NonlinearSolver.hpp:309
NonlinearSolver(const SolverParameters ¶m, std::unique_ptr< PhysicalModel > model)
Construct solver for a given model.
Definition NonlinearSolver.hpp:182
PhysicalModel & model()
Mutable reference to physical model.
Definition NonlinearSolver.hpp:294
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].
The state of a set of wells, tailored for use by the fully implicit blackoil simulator.
Definition WellState.hpp:61
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition BlackoilPhases.hpp:27
Definition NonlinearSolver.hpp:116
Definition NonlinearSolver.hpp:45
Definition NonlinearSolver.hpp:53
Definition NonlinearSolver.hpp:57
Definition NonlinearSolver.hpp:41
A struct for returning timing data from a simulator to its caller.
Definition SimulatorReport.hpp:34