21#ifndef OPM_NON_LINEAR_SOLVER_EBOS_HPP
22#define OPM_NON_LINEAR_SOLVER_EBOS_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/common/Exceptions.hpp>
33#include <dune/common/fmatrix.hh>
34#include <dune/istl/bcrsmatrix.hh>
37namespace Opm::Properties {
43template<
class TypeTag,
class MyTypeTag>
45 using type = UndefinedProperty;
47template<
class TypeTag,
class MyTypeTag>
49 using type = UndefinedProperty;
51template<
class TypeTag,
class MyTypeTag>
53 using type = UndefinedProperty;
55template<
class TypeTag,
class MyTypeTag>
57 using type = UndefinedProperty;
60template<
class TypeTag>
62 using type = GetPropType<TypeTag, Scalar>;
63 static constexpr type value = 0.5;
65template<
class TypeTag>
67 static constexpr int value = 20;
69template<
class TypeTag>
71 static constexpr int value = 1;
73template<
class TypeTag>
75 static constexpr auto value =
"dampen";
86 template <
class TypeTag,
class PhysicalModel>
89 using Scalar = GetPropType<TypeTag, Properties::Scalar>;
101 RelaxType relaxType_;
103 double relaxIncrement_;
114 relaxMax_ = EWOMS_GET_PARAM(TypeTag, Scalar, NewtonMaxRelax);
115 maxIter_ = EWOMS_GET_PARAM(TypeTag,
int, FlowNewtonMaxIterations);
116 minIter_ = EWOMS_GET_PARAM(TypeTag,
int, NewtonMinIterations);
118 const auto& relaxationTypeString = EWOMS_GET_PARAM(TypeTag, std::string, NewtonRelaxationType);
119 if (relaxationTypeString ==
"dampen") {
121 }
else if (relaxationTypeString ==
"sor") {
124 OPM_THROW(std::runtime_error,
"Unknown Relaxtion Type " << relaxationTypeString);
128 static void registerParameters()
130 EWOMS_REGISTER_PARAM(TypeTag, Scalar, NewtonMaxRelax,
"The maximum relaxation factor of a Newton iteration used by flow");
131 EWOMS_REGISTER_PARAM(TypeTag,
int, FlowNewtonMaxIterations,
"The maximum number of Newton iterations per time step used by flow");
132 EWOMS_REGISTER_PARAM(TypeTag,
int, NewtonMinIterations,
"The minimum number of Newton iterations per time step used by flow");
133 EWOMS_REGISTER_PARAM(TypeTag, std::string, NewtonRelaxationType,
"The type of relaxation used by flow's Newton method");
141 relaxIncrement_ = 0.1;
162 std::unique_ptr<PhysicalModel>
model)
164 , model_(std::move(
model))
166 , nonlinearIterations_(0)
167 , linearIterations_(0)
169 , nonlinearIterationsLast_(0)
170 , linearIterationsLast_(0)
171 , wellIterationsLast_(0)
174 OPM_THROW(std::logic_error,
"Must provide a non-null model argument for NonlinearSolver.");
186 report += model_->prepareStep(timer);
193 bool converged =
false;
201 auto iterReport = model_->nonlinearIteration(iteration, timer, *
this);
203 report += iterReport;
204 report.converged = iterReport.converged;
206 converged = report.converged;
212 failureReport_ = report;
213 failureReport_ += model_->failureReport();
217 while ( (!converged && (iteration <=
maxIter())) || (iteration <=
minIter()));
220 failureReport_ = report;
222 std::string msg =
"Solver convergence failure - Failed to complete a time step within " + std::to_string(
maxIter()) +
" iterations.";
223 OPM_THROW_NOLOG(TooManyIterations, msg);
227 report += model_->afterStep(timer);
228 report.converged =
true;
234 {
return failureReport_; }
238 {
return linearizations_; }
242 {
return nonlinearIterations_; }
246 {
return linearIterations_; }
250 {
return wellIterations_; }
254 {
return nonlinearIterationsLast_; }
258 {
return linearIterationsLast_; }
262 {
return wellIterationsLast_; }
264 std::vector<std::vector<double> >
265 computeFluidInPlace(
const std::vector<int>& fipnum)
const
266 {
return model_->computeFluidInPlace(fipnum); }
278 const int it,
bool& oscillate,
bool& stagnate)
const
292 int oscillatePhase = 0;
293 const std::vector<double>& F0 = residualHistory[it];
294 const std::vector<double>& F1 = residualHistory[it - 1];
295 const std::vector<double>& F2 = residualHistory[it - 2];
296 for (
int p= 0; p < model_->numPhases(); ++p){
297 const double d1 = std::abs((F0[p] - F2[p]) / F0[p]);
298 const double d2 = std::abs((F0[p] - F1[p]) / F0[p]);
304 stagnate = (stagnate && !(std::abs((F1[p] - F2[p]) / F2[p]) > 1.0e-3));
307 oscillate = (oscillatePhase > 1);
313 template <
class BVector>
319 BVector tempDxOld = dxOld;
328 for (i = 0; i < dx.size(); ++i) {
338 for (i = 0; i < dx.size(); ++i) {
340 tempDxOld[i] *= (1.-omega);
341 dx[i] += tempDxOld[i];
346 OPM_THROW(std::runtime_error,
"Can only handle Dampen and SOR relaxation type.");
354 {
return param_.relaxMax_; }
358 {
return param_.relaxIncrement_; }
362 {
return param_.relaxType_; }
366 {
return param_.relaxRelTol_; }
370 {
return param_.maxIter_; }
374 {
return param_.minIter_; }
383 SolverParameters param_;
384 std::unique_ptr<PhysicalModel> model_;
386 int nonlinearIterations_;
387 int linearIterations_;
389 int nonlinearIterationsLast_;
390 int linearIterationsLast_;
391 int wellIterationsLast_;
A nonlinear solver class suitable for general fully-implicit models, as well as pressure,...
Definition: NonlinearSolverEbos.hpp:88
double relaxMax() const
The greatest relaxation factor (i.e. smallest factor) allowed.
Definition: NonlinearSolverEbos.hpp:353
int maxIter() const
The maximum number of nonlinear iterations allowed.
Definition: NonlinearSolverEbos.hpp:369
const SimulatorReportSingle & failureReport() const
return the statistics if the step() method failed
Definition: NonlinearSolverEbos.hpp:233
PhysicalModel & model()
Mutable reference to physical model.
Definition: NonlinearSolverEbos.hpp:273
const PhysicalModel & model() const
Reference to physical model.
Definition: NonlinearSolverEbos.hpp:269
int wellIterationsLastStep() const
Number of well iterations used in all calls to step().
Definition: NonlinearSolverEbos.hpp:261
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: NonlinearSolverEbos.hpp:277
double relaxIncrement() const
The step-change size for the relaxation factor.
Definition: NonlinearSolverEbos.hpp:357
int linearIterationsLastStep() const
Number of linear solver iterations used in the last call to step().
Definition: NonlinearSolverEbos.hpp:257
int linearIterations() const
Number of linear solver iterations used in all calls to step().
Definition: NonlinearSolverEbos.hpp:245
enum RelaxType relaxType() const
The relaxation type (Dampen or SOR).
Definition: NonlinearSolverEbos.hpp:361
int linearizations() const
Number of linearizations used in all calls to step().
Definition: NonlinearSolverEbos.hpp:237
void stabilizeNonlinearUpdate(BVector &dx, BVector &dxOld, const double omega) const
Apply a stabilization to dx, depending on dxOld and relaxation parameters.
Definition: NonlinearSolverEbos.hpp:314
int minIter() const
The minimum number of nonlinear iterations allowed.
Definition: NonlinearSolverEbos.hpp:373
int nonlinearIterationsLastStep() const
Number of nonlinear solver iterations used in the last call to step().
Definition: NonlinearSolverEbos.hpp:253
double relaxRelTol() const
The relaxation relative tolerance.
Definition: NonlinearSolverEbos.hpp:365
int wellIterations() const
Number of well iterations used in all calls to step().
Definition: NonlinearSolverEbos.hpp:249
NonlinearSolverEbos(const SolverParameters ¶m, std::unique_ptr< PhysicalModel > model)
Construct solver for a given model.
Definition: NonlinearSolverEbos.hpp:161
void setParameters(const SolverParameters ¶m)
Set parameters to override those given at construction time.
Definition: NonlinearSolverEbos.hpp:377
int nonlinearIterations() const
Number of full nonlinear solver iterations used in all calls to step().
Definition: NonlinearSolverEbos.hpp:241
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:56
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: BlackoilPhases.hpp:27
Definition: NonlinearSolverEbos.hpp:100
Definition: NonlinearSolverEbos.hpp:48
Definition: NonlinearSolverEbos.hpp:44
Definition: NonlinearSolverEbos.hpp:52
Definition: NonlinearSolverEbos.hpp:56
Definition: NonlinearSolverEbos.hpp:40
A struct for returning timing data from a simulator to its caller.
Definition: SimulatorReport.hpp:31