00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021 #ifndef OPM_NONLINEARSOLVER_HEADER_INCLUDED
00022 #define OPM_NONLINEARSOLVER_HEADER_INCLUDED
00023
00024 #include <opm/core/simulator/SimulatorReport.hpp>
00025 #include <opm/core/utility/parameters/ParameterGroup.hpp>
00026 #include <opm/simulators/timestepping/SimulatorTimerInterface.hpp>
00027 #include <dune/common/fmatrix.hh>
00028 #include <dune/istl/bcrsmatrix.hh>
00029 #include <memory>
00030
00031 namespace Opm {
00032
00033
00036 template <class PhysicalModel>
00037 class NonlinearSolver
00038 {
00039 public:
00040
00041 enum RelaxType { DAMPEN, SOR };
00042
00043
00044 struct SolverParameters
00045 {
00046 enum RelaxType relax_type_;
00047 double relax_max_;
00048 double relax_increment_;
00049 double relax_rel_tol_;
00050 int max_iter_;
00051 int min_iter_;
00052
00053 explicit SolverParameters( const ParameterGroup& param );
00054 SolverParameters();
00055
00056 void reset();
00057 };
00058
00059
00060 typedef typename PhysicalModel::ReservoirState ReservoirState;
00061 typedef typename PhysicalModel::WellState WellState;
00062
00063
00064
00072 explicit NonlinearSolver(const SolverParameters& param,
00073 std::unique_ptr<PhysicalModel> model);
00074
00080 SimulatorReport
00081 step(const SimulatorTimerInterface& timer,
00082 ReservoirState& reservoir_state,
00083 WellState& well_state);
00084
00095 SimulatorReport
00096 step(const SimulatorTimerInterface& timer,
00097 const ReservoirState& initial_reservoir_state,
00098 const WellState& initial_well_state,
00099 ReservoirState& reservoir_state,
00100 WellState& well_state);
00101
00103 const SimulatorReport& failureReport() const
00104 { return failureReport_; }
00105
00107 int linearizations() const;
00108
00110 int nonlinearIterations() const;
00111
00113 int linearIterations() const;
00114
00116 int wellIterations() const;
00117
00119 int nonlinearIterationsLastStep() const;
00120
00122 int linearIterationsLastStep() const;
00123
00125 int wellIterationsLastStep() const;
00126
00131 std::vector<std::vector<double> >
00132 computeFluidInPlace(const ReservoirState& x, const std::vector<int>& fipnum) const
00133 {
00134 return model_->computeFluidInPlace(x, fipnum);
00135 }
00136
00137 std::vector<std::vector<double> >
00138 computeFluidInPlace(const std::vector<int>& fipnum) const
00139 {
00140 return model_->computeFluidInPlace(fipnum);
00141 }
00142
00143
00145 const PhysicalModel& model() const;
00146
00148 PhysicalModel& model();
00149
00151 void detectOscillations(const std::vector<std::vector<double>>& residual_history,
00152 const int it, bool& oscillate, bool& stagnate) const;
00153
00156 template <class BVector>
00157 void stabilizeNonlinearUpdate(BVector& dx, BVector& dxOld, const double omega) const;
00158
00160 double relaxMax() const { return param_.relax_max_; }
00161
00163 double relaxIncrement() const { return param_.relax_increment_; }
00164
00166 enum RelaxType relaxType() const { return param_.relax_type_; }
00167
00169 double relaxRelTol() const { return param_.relax_rel_tol_; }
00170
00172 int maxIter() const { return param_.max_iter_; }
00173
00175 int minIter() const { return param_.min_iter_; }
00176
00178 void setParameters(const SolverParameters& param) { param_ = param; }
00179
00180 private:
00181
00182 SimulatorReport failureReport_;
00183 SolverParameters param_;
00184 std::unique_ptr<PhysicalModel> model_;
00185 int linearizations_;
00186 int nonlinearIterations_;
00187 int linearIterations_;
00188 int wellIterations_;
00189 int nonlinearIterationsLast_;
00190 int linearIterationsLast_;
00191 int wellIterationsLast_;
00192 };
00193 }
00194
00195 #include "NonlinearSolver_impl.hpp"
00196
00197 #endif // OPM_NONLINEARSOLVER_HEADER_INCLUDED