NonlinearSolver.hpp
Go to the documentation of this file.
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
32
35
39
40#include <memory>
41
42namespace Opm::Parameters {
43
44template<class Scalar>
45struct NewtonMaxRelax { static constexpr Scalar value = 0.5; };
46
47struct NewtonMinIterations { static constexpr int value = 2; };
48struct NewtonRelaxationType { static constexpr auto value = "dampen"; };
49
50} // namespace Opm::Parameters
51
52namespace Opm {
53
54// Available relaxation scheme types.
56 Dampen,
57 SOR
58};
59
60namespace detail {
61
63template<class Scalar>
64void detectOscillations(const std::vector<std::vector<Scalar>>& residualHistory,
65 const int it, const int numPhases, const Scalar relaxRelTol,
66 const int minimumOscillatingPhases,
67 bool& oscillate, bool& stagnate);
68
71template <class BVector, class Scalar>
72void stabilizeNonlinearUpdate(BVector& dx, BVector& dxOld,
73 const Scalar omega, NonlinearRelaxType relaxType);
74
75}
76
77// Solver parameters controlling nonlinear process.
78template<class Scalar>
80{
82 Scalar relaxMax_;
85 int maxIter_; // max nonlinear iterations
86 int minIter_; // min nonlinear iterations
87
89
90 static void registerParameters();
91
92 void reset();
93};
94
97 template <class TypeTag, class PhysicalModel>
99 {
101
102 public:
104
105 // --------- Public methods ---------
106
115 std::unique_ptr<PhysicalModel> model)
116 : param_(param)
117 , model_(std::move(model))
118 , linearizations_(0)
119 , nonlinearIterations_(0)
120 , linearIterations_(0)
121 , wellIterations_(0)
122 , nonlinearIterationsLast_(0)
123 , linearIterationsLast_(0)
124 , wellIterationsLast_(0)
125 {
126 if (!model_) {
127 OPM_THROW(std::logic_error, "Must provide a non-null model argument for NonlinearSolver.");
128 }
129 }
130
131
133 {
135 report.global_time = timer.simulationTimeElapsed();
136 report.timestep_length = timer.currentStepLength();
137
138 // Do model-specific once-per-step calculations.
139 report += model_->prepareStep(timer);
140
141 int iteration = 0;
142
143 // Let the model do one nonlinear iteration.
144
145 // Set up for main solver loop.
146 bool converged = false;
147
148 // ---------- Main nonlinear solver loop ----------
149 do {
150 try {
151 // Do the nonlinear step. If we are in a converged state, the
152 // model will usually do an early return without an expensive
153 // solve, unless the minIter() count has not been reached yet.
154 auto iterReport = model_->nonlinearIteration(iteration, timer, *this);
155 iterReport.global_time = timer.simulationTimeElapsed();
156 report += iterReport;
157 report.converged = iterReport.converged;
158
159 converged = report.converged;
160 iteration += 1;
161 }
162 catch (...) {
163 // if an iteration fails during a time step, all previous iterations
164 // count as a failure as well
165 failureReport_ = report;
166 failureReport_ += model_->failureReport();
167 throw;
168 }
169 }
170 while ( (!converged && (iteration <= maxIter())) || (iteration <= minIter()));
171
172 if (!converged) {
173 failureReport_ = report;
174
175 std::string msg = "Solver convergence failure - Failed to complete a time step within " + std::to_string(maxIter()) + " iterations.";
176 OPM_THROW_NOLOG(TooManyIterations, msg);
177 }
178 auto relativeChange = model_->relativeChange();
179 if (timeStepControl && !timeStepControl->timeStepAccepted(relativeChange, timer.currentStepLength())) {
180 report.converged = false;
181 report.time_step_rejected = true;
182 failureReport_ = report;
183
184 std::string msg = "Relative change in solution for time step was " + std::to_string(relativeChange) + ", which is larger than the tolerance accepted by the timestepping algorithm.";
185 OPM_THROW_NOLOG(TimeSteppingBreakdown, msg);
186 }
187
188 // Do model-specific post-step actions.
189 report += model_->afterStep(timer);
190 report.converged = true;
191 return report;
192 }
193
196 { return failureReport_; }
197
199 int linearizations() const
200 { return linearizations_; }
201
204 { return nonlinearIterations_; }
205
208 { return linearIterations_; }
209
211 int wellIterations() const
212 { return wellIterations_; }
213
216 { return nonlinearIterationsLast_; }
217
220 { return linearIterationsLast_; }
221
224 { return wellIterationsLast_; }
225
226 std::vector<std::vector<Scalar> >
227 computeFluidInPlace(const std::vector<int>& fipnum) const
228 { return model_->computeFluidInPlace(fipnum); }
229
231 const PhysicalModel& model() const
232 { return *model_; }
233
235 PhysicalModel& model()
236 { return *model_; }
237
239 void detectOscillations(const std::vector<std::vector<Scalar>>& residualHistory,
240 const int it, bool& oscillate, bool& stagnate) const
241 {
242 detail::detectOscillations(residualHistory, it, model_->numPhases(),
243 this->relaxRelTol(), 2, oscillate, stagnate);
244 }
245
246
249 template <class BVector>
250 void stabilizeNonlinearUpdate(BVector& dx, BVector& dxOld, const Scalar omega) const
251 {
252 detail::stabilizeNonlinearUpdate(dx, dxOld, omega, this->relaxType());
253 }
254
256 Scalar relaxMax() const
257 { return param_.relaxMax_; }
258
260 Scalar relaxIncrement() const
261 { return param_.relaxIncrement_; }
262
265 { return param_.relaxType_; }
266
268 Scalar relaxRelTol() const
269 { return param_.relaxRelTol_; }
270
272 int maxIter() const
273 { return param_.maxIter_; }
274
276 int minIter() const
277 { return param_.minIter_; }
278
281 { param_ = param; }
282
283 private:
284 // --------- Data members ---------
285 SimulatorReportSingle failureReport_;
286 SolverParameters param_;
287 std::unique_ptr<PhysicalModel> model_;
288 int linearizations_;
289 int nonlinearIterations_;
290 int linearIterations_;
291 int wellIterations_;
292 int nonlinearIterationsLast_;
293 int linearIterationsLast_;
294 int wellIterationsLast_;
295 };
296
297} // namespace Opm
298
299#endif // OPM_NON_LINEAR_SOLVER_HPP
Defines a type tags and some fundamental properties all models.
Definition: NonlinearSolver.hpp:99
int linearIterationsLastStep() const
Number of linear solver iterations used in the last call to step().
Definition: NonlinearSolver.hpp:219
void stabilizeNonlinearUpdate(BVector &dx, BVector &dxOld, const Scalar omega) const
Definition: NonlinearSolver.hpp:250
Scalar relaxRelTol() const
The relaxation relative tolerance.
Definition: NonlinearSolver.hpp:268
Scalar relaxMax() const
The greatest relaxation factor (i.e. smallest factor) allowed.
Definition: NonlinearSolver.hpp:256
int linearizations() const
Number of linearizations used in all calls to step().
Definition: NonlinearSolver.hpp:199
int wellIterationsLastStep() const
Number of well iterations used in all calls to step().
Definition: NonlinearSolver.hpp:223
int nonlinearIterations() const
Number of full nonlinear solver iterations used in all calls to step().
Definition: NonlinearSolver.hpp:203
void setParameters(const SolverParameters &param)
Set parameters to override those given at construction time.
Definition: NonlinearSolver.hpp:280
NonlinearRelaxType relaxType() const
The relaxation type (Dampen or SOR).
Definition: NonlinearSolver.hpp:264
int minIter() const
The minimum number of nonlinear iterations allowed.
Definition: NonlinearSolver.hpp:276
int linearIterations() const
Number of linear solver iterations used in all calls to step().
Definition: NonlinearSolver.hpp:207
SimulatorReportSingle step(const SimulatorTimerInterface &timer, const TimeStepControlInterface *timeStepControl)
Definition: NonlinearSolver.hpp:132
int nonlinearIterationsLastStep() const
Number of nonlinear solver iterations used in the last call to step().
Definition: NonlinearSolver.hpp:215
int maxIter() const
The maximum number of nonlinear iterations allowed.
Definition: NonlinearSolver.hpp:272
const SimulatorReportSingle & failureReport() const
return the statistics if the step() method failed
Definition: NonlinearSolver.hpp:195
int wellIterations() const
Number of well iterations used in all calls to step().
Definition: NonlinearSolver.hpp:211
const PhysicalModel & model() const
Reference to physical model.
Definition: NonlinearSolver.hpp:231
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:239
NonlinearSolverParameters< Scalar > SolverParameters
Definition: NonlinearSolver.hpp:103
NonlinearSolver(const SolverParameters &param, std::unique_ptr< PhysicalModel > model)
Definition: NonlinearSolver.hpp:114
Scalar relaxIncrement() const
The step-change size for the relaxation factor.
Definition: NonlinearSolver.hpp:260
std::vector< std::vector< Scalar > > computeFluidInPlace(const std::vector< int > &fipnum) const
Definition: NonlinearSolver.hpp:227
PhysicalModel & model()
Mutable reference to physical model.
Definition: NonlinearSolver.hpp:235
Interface class for SimulatorTimer objects, to be improved.
Definition: SimulatorTimerInterface.hpp:34
virtual double currentStepLength() const =0
virtual double simulationTimeElapsed() const =0
Definition: TimeStepControlInterface.hpp:51
virtual bool timeStepAccepted(const double error, const double timeStepJustCompleted) const =0
Definition: blackoilnewtonmethodparams.hpp:31
void stabilizeNonlinearUpdate(BVector &dx, BVector &dxOld, const Scalar omega, NonlinearRelaxType relaxType)
void detectOscillations(const std::vector< std::vector< Scalar > > &residualHistory, const int it, const int numPhases, const Scalar relaxRelTol, const int minimumOscillatingPhases, bool &oscillate, bool &stagnate)
Detect oscillation or stagnation in a given residual history.
Definition: blackoilboundaryratevector.hh:39
NonlinearRelaxType
Definition: NonlinearSolver.hpp:55
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:233
std::string to_string(const ConvergenceReport::ReservoirFailure::Type t)
The Opm property system, traits with inheritance.
Definition: NonlinearSolver.hpp:80
int minIter_
Definition: NonlinearSolver.hpp:86
Scalar relaxRelTol_
Definition: NonlinearSolver.hpp:84
int maxIter_
Definition: NonlinearSolver.hpp:85
Scalar relaxMax_
Definition: NonlinearSolver.hpp:82
NonlinearRelaxType relaxType_
Definition: NonlinearSolver.hpp:81
Scalar relaxIncrement_
Definition: NonlinearSolver.hpp:83
Definition: NonlinearSolver.hpp:45
static constexpr Scalar value
Definition: NonlinearSolver.hpp:45
Definition: NonlinearSolver.hpp:47
static constexpr int value
Definition: NonlinearSolver.hpp:47
Definition: NonlinearSolver.hpp:48
static constexpr auto value
Definition: NonlinearSolver.hpp:48
A struct for returning timing data from a simulator to its caller.
Definition: SimulatorReport.hpp:34
bool converged
Definition: SimulatorReport.hpp:55
double global_time
Definition: SimulatorReport.hpp:60
double timestep_length
Definition: SimulatorReport.hpp:61
bool time_step_rejected
Definition: SimulatorReport.hpp:56