My Project
Loading...
Searching...
No Matches
NonlinearSolver.hpp
1/*
2 Copyright 2015 SINTEF ICT, Applied Mathematics.
3 Copyright 2015 Statoil ASA.
4
5 This file is part of the Open Porous Media project (OPM).
6
7 OPM is free software: you can redistribute it and/or modify
8 it under the terms of the GNU General Public License as published by
9 the Free Software Foundation, either version 3 of the License, or
10 (at your option) any later version.
11
12 OPM is distributed in the hope that it will be useful,
13 but WITHOUT ANY WARRANTY; without even the implied warranty of
14 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 GNU General Public License for more details.
16
17 You should have received a copy of the GNU General Public License
18 along with OPM. If not, see <http://www.gnu.org/licenses/>.
19*/
20
21#ifndef OPM_NON_LINEAR_SOLVER_HPP
22#define OPM_NON_LINEAR_SOLVER_HPP
23
24#include <dune/common/fmatrix.hh>
25#include <dune/istl/bcrsmatrix.hh>
26
27#include <opm/common/ErrorMacros.hpp>
28#include <opm/common/Exceptions.hpp>
29
30#include <opm/models/nonlinear/newtonmethodparams.hpp>
31#include <opm/models/nonlinear/newtonmethodproperties.hh>
32
35
36#include <opm/simulators/timestepping/SimulatorReport.hpp>
37#include <opm/simulators/timestepping/SimulatorTimerInterface.hpp>
38
39#include <memory>
40
41namespace Opm::Parameters {
42
43template<class Scalar>
44struct NewtonMaxRelax { static constexpr Scalar value = 0.5; };
45
46struct NewtonMinIterations { static constexpr int value = 2; };
47struct NewtonRelaxationType { static constexpr auto value = "dampen"; };
48
49} // namespace Opm::Parameters
50
51namespace Opm {
52
53// Available relaxation scheme types.
54enum class NonlinearRelaxType {
55 Dampen,
56 SOR
57};
58
59namespace detail {
60
62template<class Scalar>
63void detectOscillations(const std::vector<std::vector<Scalar>>& residualHistory,
64 const int it, const int numPhases, const Scalar relaxRelTol,
66 bool& oscillate, bool& stagnate);
67
70template <class BVector, class Scalar>
71void stabilizeNonlinearUpdate(BVector& dx, BVector& dxOld,
72 const Scalar omega, NonlinearRelaxType relaxType);
73
74}
75
76// Solver parameters controlling nonlinear process.
77template<class Scalar>
79{
80 NonlinearRelaxType relaxType_;
81 Scalar relaxMax_;
82 Scalar relaxIncrement_;
83 Scalar relaxRelTol_;
84 int maxIter_; // max nonlinear iterations
85 int minIter_; // min nonlinear iterations
86
88
89 static void registerParameters();
90
91 void reset();
92};
93
96 template <class TypeTag, class PhysicalModel>
98 {
100
101 public:
103
104 // --------- Public methods ---------
105
114 std::unique_ptr<PhysicalModel> model)
115 : param_(param)
116 , model_(std::move(model))
117 , linearizations_(0)
118 , nonlinearIterations_(0)
119 , linearIterations_(0)
120 , wellIterations_(0)
121 , nonlinearIterationsLast_(0)
122 , linearIterationsLast_(0)
123 , wellIterationsLast_(0)
124 {
125 if (!model_) {
126 OPM_THROW(std::logic_error, "Must provide a non-null model argument for NonlinearSolver.");
127 }
128 }
129
130
132 {
134 report.global_time = timer.simulationTimeElapsed();
135 report.timestep_length = timer.currentStepLength();
136
137 // Do model-specific once-per-step calculations.
138 report += model_->prepareStep(timer);
139
140 int iteration = 0;
141
142 // Let the model do one nonlinear iteration.
143
144 // Set up for main solver loop.
145 bool converged = false;
146
147 // ---------- Main nonlinear solver loop ----------
148 do {
149 try {
150 // Do the nonlinear step. If we are in a converged state, the
151 // model will usually do an early return without an expensive
152 // solve, unless the minIter() count has not been reached yet.
153 auto iterReport = model_->nonlinearIteration(iteration, timer, *this);
154 iterReport.global_time = timer.simulationTimeElapsed();
155 report += iterReport;
156 report.converged = iterReport.converged;
157
158 converged = report.converged;
159 iteration += 1;
160 }
161 catch (...) {
162 // if an iteration fails during a time step, all previous iterations
163 // count as a failure as well
164 failureReport_ = report;
165 failureReport_ += model_->failureReport();
166 throw;
167 }
168 }
169 while ( (!converged && (iteration <= maxIter())) || (iteration <= minIter()));
170
171 if (!converged) {
172 failureReport_ = report;
173
174 std::string msg = "Solver convergence failure - Failed to complete a time step within " + std::to_string(maxIter()) + " iterations.";
176 }
177
178 // Do model-specific post-step actions.
179 report += model_->afterStep(timer);
180 report.converged = true;
181 return report;
182 }
183
186 { return failureReport_; }
187
189 int linearizations() const
190 { return linearizations_; }
191
194 { return nonlinearIterations_; }
195
198 { return linearIterations_; }
199
201 int wellIterations() const
202 { return wellIterations_; }
203
206 { return nonlinearIterationsLast_; }
207
210 { return linearIterationsLast_; }
211
214 { return wellIterationsLast_; }
215
216 std::vector<std::vector<Scalar> >
217 computeFluidInPlace(const std::vector<int>& fipnum) const
218 { return model_->computeFluidInPlace(fipnum); }
219
221 const PhysicalModel& model() const
222 { return *model_; }
223
226 { return *model_; }
227
229 void detectOscillations(const std::vector<std::vector<Scalar>>& residualHistory,
230 const int it, bool& oscillate, bool& stagnate) const
231 {
232 detail::detectOscillations(residualHistory, it, model_->numPhases(),
233 this->relaxRelTol(), 2, oscillate, stagnate);
234 }
235
236
239 template <class BVector>
240 void stabilizeNonlinearUpdate(BVector& dx, BVector& dxOld, const Scalar omega) const
241 {
242 detail::stabilizeNonlinearUpdate(dx, dxOld, omega, this->relaxType());
243 }
244
246 Scalar relaxMax() const
247 { return param_.relaxMax_; }
248
250 Scalar relaxIncrement() const
251 { return param_.relaxIncrement_; }
252
254 NonlinearRelaxType relaxType() const
255 { return param_.relaxType_; }
256
258 Scalar relaxRelTol() const
259 { return param_.relaxRelTol_; }
260
262 int maxIter() const
263 { return param_.maxIter_; }
264
266 int minIter() const
267 { return param_.minIter_; }
268
271 { param_ = param; }
272
273 private:
274 // --------- Data members ---------
275 SimulatorReportSingle failureReport_;
276 SolverParameters param_;
277 std::unique_ptr<PhysicalModel> model_;
278 int linearizations_;
279 int nonlinearIterations_;
280 int linearIterations_;
281 int wellIterations_;
282 int nonlinearIterationsLast_;
283 int linearIterationsLast_;
284 int wellIterationsLast_;
285 };
286
287} // namespace Opm
288
289#endif // OPM_NON_LINEAR_SOLVER_HPP
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 &param)
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 &param, 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