00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021 #ifndef OPM_BLACKOILSEQUENTIALMODEL_HEADER_INCLUDED
00022 #define OPM_BLACKOILSEQUENTIALMODEL_HEADER_INCLUDED
00023
00024
00025 #include <opm/autodiff/BlackoilModelBase.hpp>
00026 #include <opm/core/simulator/BlackoilState.hpp>
00027 #include <opm/autodiff/WellStateFullyImplicitBlackoil.hpp>
00028 #include <opm/autodiff/BlackoilModelParameters.hpp>
00029 #include <opm/simulators/timestepping/SimulatorTimerInterface.hpp>
00030
00031 namespace Opm {
00032
00033 struct BlackoilSequentialModelParameters : public BlackoilModelParameters
00034 {
00035 bool iterate_to_fully_implicit;
00036 explicit BlackoilSequentialModelParameters( const ParameterGroup& param )
00037 : BlackoilModelParameters(param),
00038 iterate_to_fully_implicit(param.getDefault("iterate_to_fully_implicit", false))
00039 {
00040 }
00041 };
00042
00043
00045 template<class Grid, class WellModel,
00046 template <class G, class W> class PressureModelT,
00047 template <class G, class W> class TransportModelT>
00048 class BlackoilSequentialModel
00049 {
00050 public:
00051 typedef BlackoilState ReservoirState;
00052 typedef WellStateFullyImplicitBlackoil WellState;
00053 typedef BlackoilSequentialModelParameters ModelParameters;
00054 typedef DefaultBlackoilSolutionState SolutionState;
00055
00056 typedef PressureModelT<Grid, WellModel> PressureModel;
00057 typedef TransportModelT<Grid, WellModel> TransportModel;
00058 typedef NonlinearSolver<PressureModel> PressureSolver;
00059 typedef NonlinearSolver<TransportModel> TransportSolver;
00060
00061 typedef typename TransportModel::SimulatorData SimulatorData;
00062 typedef typename TransportModel::FIPDataType FIPDataType;
00063
00079 BlackoilSequentialModel(const ModelParameters& param,
00080 const Grid& grid ,
00081 const BlackoilPropsAdFromDeck& fluid,
00082 const DerivedGeology& geo ,
00083 const RockCompressibility* rock_comp_props,
00084 const WellModel well_model,
00085 const NewtonIterationBlackoilInterface& linsolver,
00086 std::shared_ptr< const EclipseState > eclState,
00087 const bool has_disgas,
00088 const bool has_vapoil,
00089 const bool terminal_output)
00090 : pressure_model_(new PressureModel(param, grid, fluid, geo, rock_comp_props, well_model,
00091 linsolver, eclState, has_disgas, has_vapoil, terminal_output)),
00092 transport_model_(new TransportModel(param, grid, fluid, geo, rock_comp_props, well_model,
00093 linsolver, eclState, has_disgas, has_vapoil, terminal_output)),
00094
00095 pressure_solver_(typename PressureSolver::SolverParameters(), std::move(pressure_model_)),
00096 transport_solver_(typename TransportSolver::SolverParameters(), std::move(transport_model_)),
00097 initial_reservoir_state_(0, 0, 0),
00098 iterate_to_fully_implicit_(param.iterate_to_fully_implicit)
00099 {
00100 typename PressureSolver::SolverParameters pp;
00101 pp.min_iter_ = 0;
00102 pressure_solver_.setParameters(pp);
00103 typename TransportSolver::SolverParameters tp;
00104 tp.min_iter_ = 0;
00105 transport_solver_.setParameters(tp);
00106 }
00107
00108
00109
00110
00111
00116 void prepareStep(const SimulatorTimerInterface& ,
00117 const ReservoirState& reservoir_state,
00118 const WellState& well_state)
00119 {
00120 initial_reservoir_state_ = reservoir_state;
00121 initial_well_state_ = well_state;
00122 }
00123
00124
00125
00126
00127
00136 template <class NonlinearSolverType>
00137 SimulatorReport nonlinearIteration(const int iteration,
00138 const SimulatorTimerInterface& timer,
00139 NonlinearSolverType& ,
00140 ReservoirState& reservoir_state,
00141 WellState& well_state)
00142 {
00143 if (!iterate_to_fully_implicit_) {
00144
00145 if (terminalOutputEnabled()) {
00146 OpmLog::info("Using sequential model.");
00147 }
00148
00149
00150 if (terminalOutputEnabled()) {
00151 OpmLog::info("Solving the pressure equation.");
00152 }
00153 ReservoirState initial_state = reservoir_state;
00154 const SimulatorReport pressure_report = pressure_solver_.step(timer, reservoir_state, well_state);
00155 const int pressure_liniter = pressure_report.total_linear_iterations;
00156 if (pressure_liniter == -1) {
00157 OPM_THROW(std::runtime_error, "Pressure solver failed to converge.");
00158 }
00159
00160
00161 if (terminalOutputEnabled()) {
00162 OpmLog::info("Solving the transport equations.");
00163 }
00164 const SimulatorReport transport_report = transport_solver_.step(timer, initial_state, well_state, reservoir_state, well_state);
00165 const int transport_liniter = transport_report.total_linear_iterations;
00166 if (transport_liniter == -1) {
00167 OPM_THROW(std::runtime_error, "Transport solver failed to converge.");
00168 }
00169
00170
00171 SimulatorReport report;
00172 report.converged = true;
00173 report.total_linear_iterations = pressure_liniter + transport_liniter;
00174 return report;
00175 } else {
00176
00177
00178
00179 if (terminalOutputEnabled()) {
00180 OpmLog::info("Using sequential model in iterative mode, outer iteration " + std::to_string(iteration));
00181 }
00182
00183
00184 if (terminalOutputEnabled()) {
00185 OpmLog::info("Solving the pressure equation.");
00186 }
00187
00188 const SimulatorReport pressure_report = pressure_solver_.step(timer, initial_reservoir_state_, initial_well_state_, reservoir_state, well_state);
00189 const int pressure_liniter = pressure_report.total_linear_iterations;
00190 if (pressure_liniter == -1) {
00191 OPM_THROW(std::runtime_error, "Pressure solver failed to converge.");
00192 }
00193
00194
00195 if (terminalOutputEnabled()) {
00196 OpmLog::info("Solving the transport equations.");
00197 }
00198 const SimulatorReport transport_report = transport_solver_.step(timer, initial_reservoir_state_, initial_well_state_, reservoir_state, well_state);
00199 const int transport_liniter = transport_report.total_linear_iterations;
00200 if (transport_liniter == -1) {
00201 OPM_THROW(std::runtime_error, "Transport solver failed to converge.");
00202 }
00203
00204
00205 bool done = false;
00206 {
00207 auto rstate = reservoir_state;
00208 auto wstate = well_state;
00209 pressure_solver_.model().prepareStep(timer, initial_reservoir_state_, initial_well_state_);
00210 SimulatorReport rep = pressure_solver_.model().nonlinearIteration(0, timer, pressure_solver_, rstate, wstate);
00211 if (rep.converged && rep.total_newton_iterations == 0) {
00212 done = true;
00213 }
00214 }
00215
00216 SimulatorReport report;
00217 report.converged = done;
00218 report.total_linear_iterations = pressure_liniter + transport_liniter;
00219 return report;
00220 }
00221 }
00222
00223
00224
00225
00226
00232 void afterStep(const SimulatorTimerInterface& ,
00233 ReservoirState& ,
00234 WellState& )
00235 {
00236 }
00237
00238
00239
00240
00241
00250 void setThresholdPressures(const std::vector<double>& threshold_pressures_by_face)
00251 {
00252 pressure_solver_.model().setThresholdPressures(threshold_pressures_by_face);
00253 transport_solver_.model().setThresholdPressures(threshold_pressures_by_face);
00254 }
00255
00256
00257
00258
00259
00261 bool terminalOutputEnabled() const
00262 {
00263 return pressure_solver_.model().terminalOutputEnabled();
00264 }
00265
00266
00267
00268
00269
00271 double relativeChange(const SimulationDataContainer& previous,
00272 const SimulationDataContainer& current ) const
00273 {
00274
00275 return std::max(pressure_solver_.model().relativeChange(previous, current),
00276 transport_solver_.model().relativeChange(previous, current));
00277 }
00278
00280 const WellModel& wellModel() const
00281 {
00282 return pressure_solver_.model().wellModel();
00283 }
00284
00285
00290 std::vector<std::vector<double> >
00291 computeFluidInPlace(const ReservoirState& x,
00292 const std::vector<int>& fipnum) const
00293 {
00294 return transport_solver_.computeFluidInPlace(x, fipnum);
00295 }
00296
00297
00299 const SimulatorData& getSimulatorData(const SimulationDataContainer& localState) const {
00300 return transport_solver_.model().getSimulatorData(localState);
00301 }
00302
00304 FIPDataType getFIPData() const {
00305 return transport_solver_.model().getFIPData();
00306 }
00307
00312 const SimulatorReport& failureReport() const
00313 { return failureReport_; }
00314
00315 protected:
00316 SimulatorReport failureReport_;
00317
00318 std::unique_ptr<PressureModel> pressure_model_;
00319 std::unique_ptr<TransportModel> transport_model_;
00320 PressureSolver pressure_solver_;
00321 TransportSolver transport_solver_;
00322
00323 ReservoirState initial_reservoir_state_;
00324 WellState initial_well_state_;
00325
00326 bool iterate_to_fully_implicit_;
00327 };
00328
00329 }
00330
00331
00332
00333 #endif // OPM_BLACKOILSEQUENTIALMODEL_HEADER_INCLUDED