23 #ifndef OPM_NONLINEARSOLVER_IMPL_HEADER_INCLUDED 24 #define OPM_NONLINEARSOLVER_IMPL_HEADER_INCLUDED 26 #include <opm/autodiff/NonlinearSolver.hpp> 27 #include <opm/common/Exceptions.hpp> 28 #include <opm/common/ErrorMacros.hpp> 32 template <
class PhysicalModel>
34 std::unique_ptr<PhysicalModel> model_arg)
36 model_(
std::move(model_arg)),
38 nonlinearIterations_(0),
41 nonlinearIterationsLast_(0),
42 linearIterationsLast_(0),
43 wellIterationsLast_(0)
46 OPM_THROW(std::logic_error,
"Must provide a non-null model argument for NonlinearSolver.");
50 template <
class PhysicalModel>
53 return linearizations_;
56 template <
class PhysicalModel>
59 return nonlinearIterations_;
62 template <
class PhysicalModel>
65 return linearIterations_;
68 template <
class PhysicalModel>
71 return wellIterations_;
74 template <
class PhysicalModel>
80 template <
class PhysicalModel>
86 template <
class PhysicalModel>
89 return nonlinearIterationsLast_;
92 template <
class PhysicalModel>
95 return linearIterationsLast_;
98 template <
class PhysicalModel>
101 return wellIterationsLast_;
104 template <
class PhysicalModel>
108 ReservoirState& reservoir_state,
109 WellState& well_state)
111 return step(timer, reservoir_state, well_state, reservoir_state, well_state);
116 template <
class PhysicalModel>
120 const ReservoirState& initial_reservoir_state,
121 const WellState& initial_well_state,
122 ReservoirState& reservoir_state,
123 WellState& well_state)
125 SimulatorReport iterReport;
126 SimulatorReport report;
127 failureReport_ = SimulatorReport();
130 model_->prepareStep(timer, initial_reservoir_state, initial_well_state);
137 bool converged =
false;
145 iterReport = model_->nonlinearIteration(iteration, timer, *
this, reservoir_state, well_state);
147 report += iterReport;
148 report.converged = iterReport.converged;
150 converged = report.converged;
156 failureReport_ += report;
157 failureReport_ += model_->failureReport();
160 }
while ( (!converged && (iteration <= maxIter())) || (iteration <= minIter()));
163 failureReport_ += report;
165 std::string msg =
"Solver convergence failure - Failed to complete a time step within " + std::to_string(maxIter()) +
" iterations.";
166 OPM_THROW_NOLOG(Opm::TooManyIterations, msg);
170 model_->afterStep(timer, reservoir_state, well_state);
171 report.converged =
true;
178 template <
class PhysicalModel>
183 relax_type_ = DAMPEN;
185 relax_increment_ = 0.1;
186 relax_rel_tol_ = 0.2;
191 template <
class PhysicalModel>
192 NonlinearSolver<PhysicalModel>::SolverParameters::
199 template <
class PhysicalModel>
200 NonlinearSolver<PhysicalModel>::SolverParameters::
201 SolverParameters(
const ParameterGroup& param )
207 relax_max_ = param.getDefault(
"relax_max", relax_max_);
208 max_iter_ = param.getDefault(
"max_iter", max_iter_);
209 min_iter_ = param.getDefault(
"min_iter", min_iter_);
211 std::string relaxation_type = param.getDefault(
"relax_type", std::string(
"dampen"));
212 if (relaxation_type ==
"dampen") {
213 relax_type_ = DAMPEN;
214 }
else if (relaxation_type ==
"sor") {
217 OPM_THROW(std::runtime_error,
"Unknown Relaxtion Type " << relaxation_type);
221 template <
class PhysicalModel>
225 bool& oscillate,
bool& stagnate)
const 239 int oscillatePhase = 0;
240 const std::vector<double>& F0 = residual_history[it];
241 const std::vector<double>& F1 = residual_history[it - 1];
242 const std::vector<double>& F2 = residual_history[it - 2];
243 for (
int p= 0; p < model_->numPhases(); ++p){
244 const double d1 = std::abs((F0[p] - F2[p]) / F0[p]);
245 const double d2 = std::abs((F0[p] - F1[p]) / F0[p]);
247 oscillatePhase += (d1 < relaxRelTol()) && (relaxRelTol() < d2);
251 stagnate = (stagnate && !(std::abs((F1[p] - F2[p]) / F2[p]) > 1.0e-3));
254 oscillate = (oscillatePhase > 1);
258 template <
class PhysicalModel>
259 template <
class BVector>
266 BVector tempDxOld = dxOld;
269 switch (relaxType()) {
275 for (i = 0; i < dx.size(); ++i) {
285 for (i = 0; i < dx.size(); ++i) {
287 tempDxOld[i] *= (1.-omega);
288 dx[i] += tempDxOld[i];
293 OPM_THROW(std::runtime_error,
"Can only handle DAMPEN and SOR relaxation type.");
301 #endif // OPM_FULLYIMPLICITSOLVER_IMPL_HEADER_INCLUDED A nonlinear solver class suitable for general fully-implicit models, as well as pressure, transport and sequential models.
Definition: NonlinearSolver.hpp:37
NonlinearSolver(const SolverParameters ¶m, std::unique_ptr< PhysicalModel > model)
Construct solver for a given model.
Definition: NonlinearSolver_impl.hpp:33
Definition: AutoDiff.hpp:297
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: AdditionalObjectDeleter.hpp:22
Interface class for SimulatorTimer objects, to be improved.
Definition: SimulatorTimerInterface.hpp:35