All Classes Namespaces Files Functions Variables Typedefs Enumerator Pages
BlackoilSequentialModel.hpp
1 /*
2  Copyright 2015, 2016 SINTEF ICT, Applied Mathematics.
3  Copyright 2016 Statoil AS.
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_BLACKOILSEQUENTIALMODEL_HEADER_INCLUDED
22 #define OPM_BLACKOILSEQUENTIALMODEL_HEADER_INCLUDED
23 
24 
25 #include <opm/autodiff/BlackoilModelBase.hpp>
26 #include <opm/core/simulator/BlackoilState.hpp>
27 #include <opm/autodiff/WellStateFullyImplicitBlackoil.hpp>
28 #include <opm/autodiff/BlackoilModelParameters.hpp>
29 #include <opm/simulators/timestepping/SimulatorTimerInterface.hpp>
30 
31 namespace Opm {
32 
34  {
35  bool iterate_to_fully_implicit;
36  explicit BlackoilSequentialModelParameters( const ParameterGroup& param )
37  : BlackoilModelParameters(param),
38  iterate_to_fully_implicit(param.getDefault("iterate_to_fully_implicit", false))
39  {
40  }
41  };
42 
43 
45  template<class Grid, class WellModel,
46  template <class G, class W> class PressureModelT,
47  template <class G, class W> class TransportModelT>
49  {
50  public:
51  typedef BlackoilState ReservoirState;
55 
56  typedef PressureModelT<Grid, WellModel> PressureModel;
57  typedef TransportModelT<Grid, WellModel> TransportModel;
60 
61  typedef typename TransportModel::SimulatorData SimulatorData;
62  typedef typename TransportModel::FIPDataType FIPDataType;
63 
80  const Grid& grid ,
81  const BlackoilPropsAdFromDeck& fluid,
82  const DerivedGeology& geo ,
83  const RockCompressibility* rock_comp_props,
84  const WellModel well_model,
85  const NewtonIterationBlackoilInterface& linsolver,
86  std::shared_ptr< const EclipseState > eclState,
87  const bool has_disgas,
88  const bool has_vapoil,
89  const bool terminal_output)
90  : pressure_model_(new PressureModel(param, grid, fluid, geo, rock_comp_props, well_model,
91  linsolver, eclState, has_disgas, has_vapoil, terminal_output)),
92  transport_model_(new TransportModel(param, grid, fluid, geo, rock_comp_props, well_model,
93  linsolver, eclState, has_disgas, has_vapoil, terminal_output)),
94  // TODO: fix solver parameters for pressure and transport solver.
95  pressure_solver_(typename PressureSolver::SolverParameters(), std::move(pressure_model_)),
96  transport_solver_(typename TransportSolver::SolverParameters(), std::move(transport_model_)),
97  initial_reservoir_state_(0, 0, 0), // will be overwritten
98  iterate_to_fully_implicit_(param.iterate_to_fully_implicit)
99  {
100  typename PressureSolver::SolverParameters pp;
101  pp.min_iter_ = 0;
102  pressure_solver_.setParameters(pp);
103  typename TransportSolver::SolverParameters tp;
104  tp.min_iter_ = 0;
105  transport_solver_.setParameters(tp);
106  }
107 
108 
109 
110 
111 
116  void prepareStep(const SimulatorTimerInterface& /*timer*/,
117  const ReservoirState& reservoir_state,
118  const WellState& well_state)
119  {
120  initial_reservoir_state_ = reservoir_state;
121  initial_well_state_ = well_state;
122  }
123 
124 
125 
126 
127 
136  template <class NonlinearSolverType>
137  SimulatorReport nonlinearIteration(const int iteration,
138  const SimulatorTimerInterface& timer,
139  NonlinearSolverType& /* nonlinear_solver */,
140  ReservoirState& reservoir_state,
141  WellState& well_state)
142  {
143  if (!iterate_to_fully_implicit_) {
144  // Do a single pressure solve, followed by a single transport solve.
145  if (terminalOutputEnabled()) {
146  OpmLog::info("Using sequential model.");
147  }
148 
149  // Pressure solve.
150  if (terminalOutputEnabled()) {
151  OpmLog::info("Solving the pressure equation.");
152  }
153  ReservoirState initial_state = reservoir_state;
154  const SimulatorReport pressure_report = pressure_solver_.step(timer, reservoir_state, well_state);
155  const int pressure_liniter = pressure_report.total_linear_iterations;
156  if (pressure_liniter == -1) {
157  OPM_THROW(std::runtime_error, "Pressure solver failed to converge.");
158  }
159 
160  // Transport solve.
161  if (terminalOutputEnabled()) {
162  OpmLog::info("Solving the transport equations.");
163  }
164  const SimulatorReport transport_report = transport_solver_.step(timer, initial_state, well_state, reservoir_state, well_state);
165  const int transport_liniter = transport_report.total_linear_iterations;
166  if (transport_liniter == -1) {
167  OPM_THROW(std::runtime_error, "Transport solver failed to converge.");
168  }
169 
170  // Report and return.
171  SimulatorReport report;
172  report.converged = true;
173  report.total_linear_iterations = pressure_liniter + transport_liniter;
174  return report;
175  } else {
176  // Iterate to fully implicit solution.
177  // This call is just for a single iteration (one pressure and one transport solve),
178  // we return a 'false' converged status if more are needed
179  if (terminalOutputEnabled()) {
180  OpmLog::info("Using sequential model in iterative mode, outer iteration " + std::to_string(iteration));
181  }
182 
183  // Pressure solve.
184  if (terminalOutputEnabled()) {
185  OpmLog::info("Solving the pressure equation.");
186  }
187 
188  const SimulatorReport pressure_report = pressure_solver_.step(timer, initial_reservoir_state_, initial_well_state_, reservoir_state, well_state);
189  const int pressure_liniter = pressure_report.total_linear_iterations;
190  if (pressure_liniter == -1) {
191  OPM_THROW(std::runtime_error, "Pressure solver failed to converge.");
192  }
193 
194  // Transport solve.
195  if (terminalOutputEnabled()) {
196  OpmLog::info("Solving the transport equations.");
197  }
198  const SimulatorReport transport_report = transport_solver_.step(timer, initial_reservoir_state_, initial_well_state_, reservoir_state, well_state);
199  const int transport_liniter = transport_report.total_linear_iterations;
200  if (transport_liniter == -1) {
201  OPM_THROW(std::runtime_error, "Transport solver failed to converge.");
202  }
203 
204  // Revisit pressure equation to check if it is still converged.
205  bool done = false;
206  {
207  auto rstate = reservoir_state;
208  auto wstate = well_state;
209  pressure_solver_.model().prepareStep(timer, initial_reservoir_state_, initial_well_state_);
210  SimulatorReport rep = pressure_solver_.model().nonlinearIteration(0, timer, pressure_solver_, rstate, wstate);
211  if (rep.converged && rep.total_newton_iterations == 0) {
212  done = true;
213  }
214  }
215 
216  SimulatorReport report;
217  report.converged = done;
218  report.total_linear_iterations = pressure_liniter + transport_liniter;
219  return report;
220  }
221  }
222 
223 
224 
225 
226 
232  void afterStep(const SimulatorTimerInterface& /* timer */,
233  ReservoirState& /* reservoir_state */,
234  WellState& /* well_state */)
235  {
236  }
237 
238 
239 
240 
241 
250  void setThresholdPressures(const std::vector<double>& threshold_pressures_by_face)
251  {
252  pressure_solver_.model().setThresholdPressures(threshold_pressures_by_face);
253  transport_solver_.model().setThresholdPressures(threshold_pressures_by_face);
254  }
255 
256 
257 
258 
259 
262  {
263  return pressure_solver_.model().terminalOutputEnabled();
264  }
265 
266 
267 
268 
269 
271  double relativeChange(const SimulationDataContainer& previous,
272  const SimulationDataContainer& current ) const
273  {
274  // TODO: this is a quick stopgap implementation, and should be evaluated more carefully.
275  return std::max(pressure_solver_.model().relativeChange(previous, current),
276  transport_solver_.model().relativeChange(previous, current));
277  }
278 
280  const WellModel& wellModel() const
281  {
282  return pressure_solver_.model().wellModel();
283  }
284 
285 
290  std::vector<std::vector<double> >
291  computeFluidInPlace(const ReservoirState& x,
292  const std::vector<int>& fipnum) const
293  {
294  return transport_solver_.computeFluidInPlace(x, fipnum);
295  }
296 
297 
299  const SimulatorData& getSimulatorData(const SimulationDataContainer& localState) const {
300  return transport_solver_.model().getSimulatorData(localState);
301  }
302 
304  FIPDataType getFIPData() const {
305  return transport_solver_.model().getFIPData();
306  }
307 
312  const SimulatorReport& failureReport() const
313  { return failureReport_; }
314 
315  protected:
316  SimulatorReport failureReport_;
317 
318  std::unique_ptr<PressureModel> pressure_model_;
319  std::unique_ptr<TransportModel> transport_model_;
320  PressureSolver pressure_solver_;
321  TransportSolver transport_solver_;
322 
323  ReservoirState initial_reservoir_state_;
324  WellState initial_well_state_;
325 
326  bool iterate_to_fully_implicit_;
327  };
328 
329 } // namespace Opm
330 
331 
332 
333 #endif // OPM_BLACKOILSEQUENTIALMODEL_HEADER_INCLUDED
void prepareStep(const SimulatorTimerInterface &, const ReservoirState &reservoir_state, const WellState &well_state)
Called once before each time step.
Definition: BlackoilSequentialModel.hpp:116
const PhysicalModel & model() const
Reference to physical model.
Definition: NonlinearSolver_impl.hpp:75
Solver parameters for the BlackoilModel.
Definition: BlackoilModelParameters.hpp:31
const WellModel & wellModel() const
Return the well model.
Definition: BlackoilSequentialModel.hpp:280
Definition: BlackoilSequentialModel.hpp:33
std::vector< std::vector< double > > computeFluidInPlace(const ReservoirState &x, const std::vector< int > &fipnum) const
Compute fluid in place.
Definition: BlackoilSequentialModel.hpp:291
void setThresholdPressures(const std::vector< double > &threshold_pressures_by_face)
Set threshold pressures that prevent or reduce flow.
Definition: BlackoilSequentialModel.hpp:250
BlackoilModelParameters()
Construct with default parameters.
Definition: BlackoilModelParameters.cpp:28
const SimulatorData & getSimulatorData(const SimulationDataContainer &localState) const
Return reservoir simulation data (for output functionality)
Definition: BlackoilSequentialModel.hpp:299
bool terminalOutputEnabled() const
Return true if output to cout is wanted.
Definition: BlackoilSequentialModel.hpp:261
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
BlackoilSequentialModel(const ModelParameters &param, const Grid &grid, const BlackoilPropsAdFromDeck &fluid, const DerivedGeology &geo, const RockCompressibility *rock_comp_props, const WellModel well_model, const NewtonIterationBlackoilInterface &linsolver, std::shared_ptr< const EclipseState > eclState, const bool has_disgas, const bool has_vapoil, const bool terminal_output)
Construct the model.
Definition: BlackoilSequentialModel.hpp:79
Interface class for SimulatorTimer objects, to be improved.
Definition: SimulatorTimerInterface.hpp:35
FIPDataType getFIPData() const
Return fluid-in-place data (for output functionality)
Definition: BlackoilSequentialModel.hpp:304
The state of a set of wells, tailored for use by the fully implicit blackoil simulator.
Definition: WellStateFullyImplicitBlackoil.hpp:44
std::vector< std::vector< double > > computeFluidInPlace(const ReservoirState &x, const std::vector< int > &fipnum) const
Compute fluid in place.
Definition: NonlinearSolver.hpp:132
Struct for containing iteration variables.
Definition: DefaultBlackoilSolutionState.hpp:29
void setParameters(const SolverParameters &param)
Set parameters to override those given at construction time.
Definition: NonlinearSolver.hpp:178
double relativeChange(const SimulationDataContainer &previous, const SimulationDataContainer &current) const
Return the relative change in variables relevant to this model.
Definition: BlackoilSequentialModel.hpp:271
void afterStep(const SimulatorTimerInterface &, ReservoirState &, WellState &)
Called once after each time step.
Definition: BlackoilSequentialModel.hpp:232
Class containing static geological properties that are derived from grid and petrophysical properties...
Definition: GeoProps.hpp:59
const SimulatorReport & failureReport() const
return the statistics if the nonlinearIteration() method failed.
Definition: BlackoilSequentialModel.hpp:312
Interface class for (linear) solvers for the fully implicit black-oil system.
Definition: NewtonIterationBlackoilInterface.hpp:31
This class implements the AD-adapted fluid interface for three-phase black-oil.
Definition: BlackoilPropsAdFromDeck.hpp:61
SimulatorReport nonlinearIteration(const int iteration, const SimulatorTimerInterface &timer, NonlinearSolverType &, ReservoirState &reservoir_state, WellState &well_state)
Called once per nonlinear iteration.
Definition: BlackoilSequentialModel.hpp:137
A sequential splitting model implementation for three-phase black oil.
Definition: BlackoilSequentialModel.hpp:48