All Classes Namespaces Files Functions Variables Typedefs Enumerator Pages
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_NONLINEARSOLVER_HEADER_INCLUDED
22 #define OPM_NONLINEARSOLVER_HEADER_INCLUDED
23 
24 #include <opm/core/simulator/SimulatorReport.hpp>
25 #include <opm/core/utility/parameters/ParameterGroup.hpp>
26 #include <opm/simulators/timestepping/SimulatorTimerInterface.hpp>
27 #include <dune/common/fmatrix.hh>
28 #include <dune/istl/bcrsmatrix.hh>
29 #include <memory>
30 
31 namespace Opm {
32 
33 
36  template <class PhysicalModel>
38  {
39  public:
40  // Available relaxation scheme types.
41  enum RelaxType { DAMPEN, SOR };
42 
43  // Solver parameters controlling nonlinear process.
45  {
46  enum RelaxType relax_type_;
47  double relax_max_;
48  double relax_increment_;
49  double relax_rel_tol_;
50  int max_iter_; // max nonlinear iterations
51  int min_iter_; // min nonlinear iterations
52 
53  explicit SolverParameters( const ParameterGroup& param );
55 
56  void reset();
57  };
58 
59  // Forwarding types from PhysicalModel.
60  typedef typename PhysicalModel::ReservoirState ReservoirState;
61  typedef typename PhysicalModel::WellState WellState;
62 
63  // --------- Public methods ---------
64 
72  explicit NonlinearSolver(const SolverParameters& param,
73  std::unique_ptr<PhysicalModel> model);
74 
80  SimulatorReport
81  step(const SimulatorTimerInterface& timer,
82  ReservoirState& reservoir_state,
83  WellState& well_state);
84 
95  SimulatorReport
96  step(const SimulatorTimerInterface& timer,
97  const ReservoirState& initial_reservoir_state,
98  const WellState& initial_well_state,
99  ReservoirState& reservoir_state,
100  WellState& well_state);
101 
103  const SimulatorReport& failureReport() const
104  { return failureReport_; }
105 
107  int linearizations() const;
108 
110  int nonlinearIterations() const;
111 
113  int linearIterations() const;
114 
116  int wellIterations() const;
117 
119  int nonlinearIterationsLastStep() const;
120 
122  int linearIterationsLastStep() const;
123 
125  int wellIterationsLastStep() const;
126 
131  std::vector<std::vector<double> >
132  computeFluidInPlace(const ReservoirState& x, const std::vector<int>& fipnum) const
133  {
134  return model_->computeFluidInPlace(x, fipnum);
135  }
136 
137  std::vector<std::vector<double> >
138  computeFluidInPlace(const std::vector<int>& fipnum) const
139  {
140  return model_->computeFluidInPlace(fipnum);
141  }
142 
143 
145  const PhysicalModel& model() const;
146 
148  PhysicalModel& model();
149 
151  void detectOscillations(const std::vector<std::vector<double>>& residual_history,
152  const int it, bool& oscillate, bool& stagnate) const;
153 
156  template <class BVector>
157  void stabilizeNonlinearUpdate(BVector& dx, BVector& dxOld, const double omega) const;
158 
160  double relaxMax() const { return param_.relax_max_; }
161 
163  double relaxIncrement() const { return param_.relax_increment_; }
164 
166  enum RelaxType relaxType() const { return param_.relax_type_; }
167 
169  double relaxRelTol() const { return param_.relax_rel_tol_; }
170 
172  int maxIter() const { return param_.max_iter_; }
173 
175  int minIter() const { return param_.min_iter_; }
176 
178  void setParameters(const SolverParameters& param) { param_ = param; }
179 
180  private:
181  // --------- Data members ---------
182  SimulatorReport failureReport_;
183  SolverParameters param_;
184  std::unique_ptr<PhysicalModel> model_;
185  int linearizations_;
186  int nonlinearIterations_;
187  int linearIterations_;
188  int wellIterations_;
189  int nonlinearIterationsLast_;
190  int linearIterationsLast_;
191  int wellIterationsLast_;
192  };
193 } // namespace Opm
194 
195 #include "NonlinearSolver_impl.hpp"
196 
197 #endif // OPM_NONLINEARSOLVER_HEADER_INCLUDED
Definition: NonlinearSolver.hpp:44
int maxIter() const
The maximum number of nonlinear iterations allowed.
Definition: NonlinearSolver.hpp:172
int linearIterationsLastStep() const
Number of linear solver iterations used in the last call to step().
Definition: NonlinearSolver_impl.hpp:93
A nonlinear solver class suitable for general fully-implicit models, as well as pressure, transport and sequential models.
Definition: NonlinearSolver.hpp:37
double relaxMax() const
The greatest relaxation factor (i.e. smallest factor) allowed.
Definition: NonlinearSolver.hpp:160
int linearIterations() const
Number of linear solver iterations used in all calls to step().
Definition: NonlinearSolver_impl.hpp:63
void stabilizeNonlinearUpdate(BVector &dx, BVector &dxOld, const double omega) const
Apply a stabilization to dx, depending on dxOld and relaxation parameters.
Definition: NonlinearSolver_impl.hpp:261
NonlinearSolver(const SolverParameters &param, std::unique_ptr< PhysicalModel > model)
Construct solver for a given model.
Definition: NonlinearSolver_impl.hpp:33
int nonlinearIterations() const
Number of full nonlinear solver iterations used in all calls to step().
Definition: NonlinearSolver_impl.hpp:57
const PhysicalModel & model() const
Reference to physical model.
Definition: NonlinearSolver_impl.hpp:75
const SimulatorReport & failureReport() const
return the statistics if the step() method failed
Definition: NonlinearSolver.hpp:103
int wellIterationsLastStep() const
Number of well iterations used in all calls to step().
Definition: NonlinearSolver_impl.hpp:99
int nonlinearIterationsLastStep() const
Number of nonlinear solver iterations used in the last call to step().
Definition: NonlinearSolver_impl.hpp:87
int linearizations() const
Number of linearizations used in all calls to step().
Definition: NonlinearSolver_impl.hpp:51
SimulatorReport step(const SimulatorTimerInterface &timer, ReservoirState &reservoir_state, WellState &well_state)
Take a single forward step, after which the states will be modified according to the physical model...
Definition: NonlinearSolver_impl.hpp:107
Interface class for SimulatorTimer objects, to be improved.
Definition: SimulatorTimerInterface.hpp:35
std::vector< std::vector< double > > computeFluidInPlace(const ReservoirState &x, const std::vector< int > &fipnum) const
Compute fluid in place.
Definition: NonlinearSolver.hpp:132
int wellIterations() const
Number of well iterations used in all calls to step().
Definition: NonlinearSolver_impl.hpp:69
void setParameters(const SolverParameters &param)
Set parameters to override those given at construction time.
Definition: NonlinearSolver.hpp:178
double relaxRelTol() const
The relaxation relative tolerance.
Definition: NonlinearSolver.hpp:169
int minIter() const
The minimum number of nonlinear iterations allowed.
Definition: NonlinearSolver.hpp:175
void detectOscillations(const std::vector< std::vector< double >> &residual_history, const int it, bool &oscillate, bool &stagnate) const
Detect oscillation or stagnation in a given residual history.
Definition: NonlinearSolver_impl.hpp:223
double relaxIncrement() const
The step-change size for the relaxation factor.
Definition: NonlinearSolver.hpp:163
enum RelaxType relaxType() const
The relaxation type (DAMPEN or SOR).
Definition: NonlinearSolver.hpp:166