BlackoilModel.hpp
Go to the documentation of this file.
1/*
2 Copyright 2013, 2015 SINTEF ICT, Applied Mathematics.
3 Copyright 2014, 2015 Dr. Blatt - HPC-Simulation-Software & Services
4 Copyright 2014, 2015 Statoil ASA.
5 Copyright 2015 NTNU
6 Copyright 2015, 2016, 2017 IRIS AS
7
8 This file is part of the Open Porous Media project (OPM).
9
10 OPM is free software: you can redistribute it and/or modify
11 it under the terms of the GNU General Public License as published by
12 the Free Software Foundation, either version 3 of the License, or
13 (at your option) any later version.
14
15 OPM is distributed in the hope that it will be useful,
16 but WITHOUT ANY WARRANTY; without even the implied warranty of
17 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
18 GNU General Public License for more details.
19
20 You should have received a copy of the GNU General Public License
21 along with OPM. If not, see <http://www.gnu.org/licenses/>.
22*/
23
24#ifndef OPM_BLACKOILMODEL_HEADER_INCLUDED
25#define OPM_BLACKOILMODEL_HEADER_INCLUDED
26
28
34
36
40
42
44
45#include <memory>
46#include <tuple>
47#include <vector>
48
49#include <fmt/format.h>
50
51namespace Opm {
52
59template <class TypeTag>
61{
62public:
63 // --------- Types and enums ---------
77
78 static constexpr int numEq = Indices::numEq;
79 static constexpr int contiSolventEqIdx = Indices::contiSolventEqIdx;
80 static constexpr int contiZfracEqIdx = Indices::contiZfracEqIdx;
81 static constexpr int contiPolymerEqIdx = Indices::contiPolymerEqIdx;
82 static constexpr int contiEnergyEqIdx = Indices::contiEnergyEqIdx;
83 static constexpr int contiPolymerMWEqIdx = Indices::contiPolymerMWEqIdx;
84 static constexpr int contiFoamEqIdx = Indices::contiFoamEqIdx;
85 static constexpr int contiBrineEqIdx = Indices::contiBrineEqIdx;
86 static constexpr int contiMicrobialEqIdx = Indices::contiMicrobialEqIdx;
87 static constexpr int contiOxygenEqIdx = Indices::contiOxygenEqIdx;
88 static constexpr int contiUreaEqIdx = Indices::contiUreaEqIdx;
89 static constexpr int contiBiofilmEqIdx = Indices::contiBiofilmEqIdx;
90 static constexpr int contiCalciteEqIdx = Indices::contiCalciteEqIdx;
91 static constexpr int solventSaturationIdx = Indices::solventSaturationIdx;
92 static constexpr int zFractionIdx = Indices::zFractionIdx;
93 static constexpr int polymerConcentrationIdx = Indices::polymerConcentrationIdx;
94 static constexpr int polymerMoleWeightIdx = Indices::polymerMoleWeightIdx;
95 static constexpr int temperatureIdx = Indices::temperatureIdx;
96 static constexpr int foamConcentrationIdx = Indices::foamConcentrationIdx;
97 static constexpr int saltConcentrationIdx = Indices::saltConcentrationIdx;
98 static constexpr int microbialConcentrationIdx = Indices::microbialConcentrationIdx;
99 static constexpr int oxygenConcentrationIdx = Indices::oxygenConcentrationIdx;
100 static constexpr int ureaConcentrationIdx = Indices::ureaConcentrationIdx;
101 static constexpr int biofilmConcentrationIdx = Indices::biofilmConcentrationIdx;
102 static constexpr int calciteConcentrationIdx = Indices::calciteConcentrationIdx;
103
104 using VectorBlockType = Dune::FieldVector<Scalar, numEq>;
105 using MatrixBlockType = typename SparseMatrixAdapter::MatrixBlock;
106 using Mat = typename SparseMatrixAdapter::IstlMatrix;
107 using BVector = Dune::BlockVector<VectorBlockType>;
108
110
111 // --------- Public methods ---------
112
121 const ModelParameters& param,
122 BlackoilWellModel<TypeTag>& well_model,
123 const bool terminal_output);
124
125 bool isParallel() const
126 { return grid_.comm().size() > 1; }
127
128 const EclipseState& eclState() const
129 { return simulator_.vanguard().eclState(); }
130
134
136 const int iteration,
137 const int minIter,
138 const int maxIter,
139 const SimulatorTimerInterface& timer);
140
148 template <class NonlinearSolverType>
149 SimulatorReportSingle nonlinearIteration(const int iteration,
150 const SimulatorTimerInterface& timer,
151 NonlinearSolverType& nonlinear_solver);
152
153 template <class NonlinearSolverType>
155 const SimulatorTimerInterface& timer,
156 NonlinearSolverType& nonlinear_solver);
157
162
165 const int iterationIdx);
166
167 // compute the "relative" change of the solution between time steps
168 Scalar relativeChange() const;
169
172 { return simulator_.model().newtonMethod().linearSolver().iterations (); }
173
174 // Obtain reference to linear solver setup time
176 { return linear_solve_setup_time_; }
177
181
183 void updateSolution(const BVector& dx);
184
187 { return terminal_output_; }
188
189 std::tuple<Scalar,Scalar>
191 const Scalar pvSumLocal,
192 const Scalar numAquiferPvSumLocal,
193 std::vector<Scalar>& R_sum,
194 std::vector<Scalar>& maxCoeff,
195 std::vector<Scalar>& B_avg);
196
200 std::pair<Scalar,Scalar>
201 localConvergenceData(std::vector<Scalar>& R_sum,
202 std::vector<Scalar>& maxCoeff,
203 std::vector<Scalar>& B_avg,
204 std::vector<int>& maxCoeffCell);
205
209 std::pair<std::vector<double>, std::vector<int>>
210 characteriseCnvPvSplit(const std::vector<Scalar>& B_avg, const double dt);
211
214 void convergencePerCell(const std::vector<Scalar>& B_avg,
215 const double dt,
216 const double tol_cnv,
217 const double tol_cnv_energy,
218 const int iteration);
219
220 void updateTUNING(const Tuning& tuning);
221
223 getReservoirConvergence(const double reportTime,
224 const double dt,
225 const int iteration,
226 const int maxIter,
227 std::vector<Scalar>& B_avg,
228 std::vector<Scalar>& residual_norms);
229
238 const int iteration,
239 const int maxIter,
240 std::vector<Scalar>& residual_norms);
241
243 int numPhases() const
244 { return Indices::numPhases; }
245
247 template<class T>
248 std::vector<std::vector<Scalar> >
249 computeFluidInPlace(const T&, const std::vector<int>& fipnum) const
250 { return this->computeFluidInPlace(fipnum); }
251
253 std::vector<std::vector<Scalar> >
254 computeFluidInPlace(const std::vector<int>& /*fipnum*/) const;
255
256 const Simulator& simulator() const
257 { return simulator_; }
258
260 { return simulator_; }
261
264 { return failureReport_; }
265
268
270 const std::vector<SimulatorReport>& domainAccumulatedReports() const;
271
273 void writeNonlinearIterationsPerCell(const std::filesystem::path& odir) const;
274
275 const std::vector<StepReport>& stepReports() const
276 { return convergence_reports_; }
277
278 void writePartitions(const std::filesystem::path& odir) const;
279
283 { return well_model_; }
284
286 wellModel() const
287 { return well_model_; }
288
290 { simulator_.problem().beginEpisode(); }
291
293 { simulator_.problem().endEpisode(); }
294
295 template<class FluidState, class Residual>
296 void getMaxCoeff(const unsigned cell_idx,
297 const IntensiveQuantities& intQuants,
298 const FluidState& fs,
299 const Residual& modelResid,
300 const Scalar pvValue,
301 std::vector<Scalar>& B_avg,
302 std::vector<Scalar>& R_sum,
303 std::vector<Scalar>& maxCoeff,
304 std::vector<int>& maxCoeffCell);
305
307 const ModelParameters& param() const
308 { return param_; }
309
312 { return compNames_; }
313
315 bool hasNlddSolver() const
316 { return nlddSolver_ != nullptr; }
317
318protected:
319 // --------- Data members ---------
321 const Grid& grid_;
322 static constexpr bool has_solvent_ = getPropValue<TypeTag, Properties::EnableSolvent>();
323 static constexpr bool has_extbo_ = getPropValue<TypeTag, Properties::EnableExtbo>();
324 static constexpr bool has_polymer_ = getPropValue<TypeTag, Properties::EnablePolymer>();
325 static constexpr bool has_polymermw_ = getPropValue<TypeTag, Properties::EnablePolymerMW>();
326 static constexpr bool has_energy_ = getPropValue<TypeTag, Properties::EnableEnergy>();
327 static constexpr bool has_foam_ = getPropValue<TypeTag, Properties::EnableFoam>();
328 static constexpr bool has_brine_ = getPropValue<TypeTag, Properties::EnableBrine>();
329 static constexpr bool has_micp_ = getPropValue<TypeTag, Properties::EnableMICP>();
330
333
334 // Well Model
336
340 long int global_nc_;
341
342 std::vector<std::vector<Scalar>> residual_norms_history_;
345
346 std::vector<StepReport> convergence_reports_;
348
349 std::unique_ptr<BlackoilModelNldd<TypeTag>> nlddSolver_;
351
352private:
353 Scalar dpMaxRel() const { return param_.dp_max_rel_; }
354 Scalar dsMax() const { return param_.ds_max_; }
355 Scalar drMaxRel() const { return param_.dr_max_rel_; }
356 Scalar maxResidualAllowed() const { return param_.max_residual_allowed_; }
357 double linear_solve_setup_time_;
358 std::vector<bool> wasSwitched_;
359};
360
361} // namespace Opm
362
364
365#endif // OPM_BLACKOILMODEL_HEADER_INCLUDED
Implementation of penalty cards for three-phase black oil.
Definition: BlackoilModelConvergenceMonitor.hpp:35
Definition: BlackoilModel.hpp:61
static constexpr int contiEnergyEqIdx
Definition: BlackoilModel.hpp:82
SimulatorReportSingle failureReport_
Definition: BlackoilModel.hpp:332
BlackoilModel(Simulator &simulator, const ModelParameters &param, BlackoilWellModel< TypeTag > &well_model, const bool terminal_output)
Definition: BlackoilModel_impl.hpp:51
std::vector< std::vector< Scalar > > computeFluidInPlace(const T &, const std::vector< int > &fipnum) const
Wrapper required due to not following generic API.
Definition: BlackoilModel.hpp:249
static constexpr bool has_energy_
Definition: BlackoilModel.hpp:326
ModelParameters param_
Definition: BlackoilModel.hpp:331
int linearIterationsLastSolve() const
Number of linear iterations used in last call to solveJacobianSystem().
Definition: BlackoilModel.hpp:171
SimulatorReportSingle nonlinearIteration(const int iteration, const SimulatorTimerInterface &timer, NonlinearSolverType &nonlinear_solver)
Definition: BlackoilModel_impl.hpp:214
void writeNonlinearIterationsPerCell(const std::filesystem::path &odir) const
Write the number of nonlinear iterations per cell to a file in ResInsight compatible format.
Definition: BlackoilModel_impl.hpp:1046
static constexpr bool has_extbo_
Definition: BlackoilModel.hpp:323
SimulatorReportSingle assembleReservoir(const SimulatorTimerInterface &, const int iterationIdx)
Assemble the residual and Jacobian of the nonlinear system.
Definition: BlackoilModel_impl.hpp:357
GetPropType< TypeTag, Properties::SparseMatrixAdapter > SparseMatrixAdapter
Definition: BlackoilModel.hpp:68
ConvergenceReport getConvergence(const SimulatorTimerInterface &timer, const int iteration, const int maxIter, std::vector< Scalar > &residual_norms)
Definition: BlackoilModel_impl.hpp:988
ConvergenceReport getReservoirConvergence(const double reportTime, const double dt, const int iteration, const int maxIter, std::vector< Scalar > &B_avg, std::vector< Scalar > &residual_norms)
Definition: BlackoilModel_impl.hpp:740
GetPropType< TypeTag, Properties::IntensiveQuantities > IntensiveQuantities
Definition: BlackoilModel.hpp:67
static constexpr bool has_foam_
Definition: BlackoilModel.hpp:327
static constexpr int foamConcentrationIdx
Definition: BlackoilModel.hpp:96
Scalar current_relaxation_
Definition: BlackoilModel.hpp:343
GetPropType< TypeTag, Properties::FluidSystem > FluidSystem
Definition: BlackoilModel.hpp:71
static constexpr int microbialConcentrationIdx
Definition: BlackoilModel.hpp:98
static constexpr int numEq
Definition: BlackoilModel.hpp:78
const std::vector< SimulatorReport > & domainAccumulatedReports() const
return the statistics of local solves accumulated for each domain on this rank
Definition: BlackoilModel_impl.hpp:1036
static constexpr bool has_polymermw_
Definition: BlackoilModel.hpp:325
static constexpr int biofilmConcentrationIdx
Definition: BlackoilModel.hpp:101
static constexpr int contiFoamEqIdx
Definition: BlackoilModel.hpp:84
static constexpr int polymerMoleWeightIdx
Definition: BlackoilModel.hpp:94
void endReportStep()
Definition: BlackoilModel.hpp:292
BlackoilModelConvergenceMonitor< Scalar > conv_monitor_
Definition: BlackoilModel.hpp:350
static constexpr int contiCalciteEqIdx
Definition: BlackoilModel.hpp:90
std::pair< Scalar, Scalar > localConvergenceData(std::vector< Scalar > &R_sum, std::vector< Scalar > &maxCoeff, std::vector< Scalar > &B_avg, std::vector< int > &maxCoeffCell)
Get reservoir quantities on this process needed for convergence calculations.
Definition: BlackoilModel_impl.hpp:604
GetPropType< TypeTag, Properties::ElementContext > ElementContext
Definition: BlackoilModel.hpp:66
bool isParallel() const
Definition: BlackoilModel.hpp:125
const std::vector< StepReport > & stepReports() const
Definition: BlackoilModel.hpp:275
std::vector< StepReport > convergence_reports_
Definition: BlackoilModel.hpp:346
static constexpr int contiMicrobialEqIdx
Definition: BlackoilModel.hpp:86
bool terminalOutputEnabled() const
Return true if output to cout is wanted.
Definition: BlackoilModel.hpp:186
double & linearSolveSetupTime()
Definition: BlackoilModel.hpp:175
int numPhases() const
The number of active fluid phases in the model.
Definition: BlackoilModel.hpp:243
const SimulatorReport & localAccumulatedReports() const
return the statistics of local solves accumulated for this rank
Definition: BlackoilModel_impl.hpp:1025
SimulatorReportSingle afterStep(const SimulatorTimerInterface &)
Definition: BlackoilModel_impl.hpp:344
std::tuple< Scalar, Scalar > convergenceReduction(Parallel::Communication comm, const Scalar pvSumLocal, const Scalar numAquiferPvSumLocal, std::vector< Scalar > &R_sum, std::vector< Scalar > &maxCoeff, std::vector< Scalar > &B_avg)
Definition: BlackoilModel_impl.hpp:544
static constexpr int contiBiofilmEqIdx
Definition: BlackoilModel.hpp:89
const Grid & grid_
Definition: BlackoilModel.hpp:321
static constexpr int contiUreaEqIdx
Definition: BlackoilModel.hpp:88
static constexpr int calciteConcentrationIdx
Definition: BlackoilModel.hpp:102
static constexpr int temperatureIdx
Definition: BlackoilModel.hpp:95
GetPropType< TypeTag, Properties::SolutionVector > SolutionVector
Definition: BlackoilModel.hpp:69
const ComponentName & compNames() const
Returns const reference to component names.
Definition: BlackoilModel.hpp:311
GetPropType< TypeTag, Properties::PrimaryVariables > PrimaryVariables
Definition: BlackoilModel.hpp:70
void initialLinearization(SimulatorReportSingle &report, const int iteration, const int minIter, const int maxIter, const SimulatorTimerInterface &timer)
Definition: BlackoilModel_impl.hpp:155
void updateSolution(const BVector &dx)
Apply an update to the primary variables.
Definition: BlackoilModel_impl.hpp:520
long int global_nc_
The number of cells of the global grid.
Definition: BlackoilModel.hpp:340
GetPropType< TypeTag, Properties::Indices > Indices
Definition: BlackoilModel.hpp:72
void updateTUNING(const Tuning &tuning)
Definition: BlackoilModel_impl.hpp:727
typename SparseMatrixAdapter::MatrixBlock MatrixBlockType
Definition: BlackoilModel.hpp:105
const ModelParameters & param() const
Returns const reference to model parameters.
Definition: BlackoilModel.hpp:307
static constexpr int contiOxygenEqIdx
Definition: BlackoilModel.hpp:87
void getMaxCoeff(const unsigned cell_idx, const IntensiveQuantities &intQuants, const FluidState &fs, const Residual &modelResid, const Scalar pvValue, std::vector< Scalar > &B_avg, std::vector< Scalar > &R_sum, std::vector< Scalar > &maxCoeff, std::vector< int > &maxCoeffCell)
Definition: BlackoilModel_impl.hpp:1083
std::unique_ptr< BlackoilModelNldd< TypeTag > > nlddSolver_
Non-linear DD solver.
Definition: BlackoilModel.hpp:349
SimulatorReportSingle nonlinearIterationNewton(const int iteration, const SimulatorTimerInterface &timer, NonlinearSolverType &nonlinear_solver)
Definition: BlackoilModel_impl.hpp:252
void writePartitions(const std::filesystem::path &odir) const
Definition: BlackoilModel_impl.hpp:1056
static constexpr int solventSaturationIdx
Definition: BlackoilModel.hpp:91
bool hasNlddSolver() const
Returns true if an NLDD solver exists.
Definition: BlackoilModel.hpp:315
std::pair< std::vector< double >, std::vector< int > > characteriseCnvPvSplit(const std::vector< Scalar > &B_avg, const double dt)
Compute pore-volume/cell count split among "converged", "relaxed converged", "unconverged" cells base...
Definition: BlackoilModel_impl.hpp:655
void beginReportStep()
Definition: BlackoilModel.hpp:289
static constexpr int polymerConcentrationIdx
Definition: BlackoilModel.hpp:93
GetPropType< TypeTag, Properties::Simulator > Simulator
Definition: BlackoilModel.hpp:64
const SimulatorReportSingle & failureReport() const
return the statistics if the nonlinearIteration() method failed
Definition: BlackoilModel.hpp:263
Simulator & simulator()
Definition: BlackoilModel.hpp:259
GetPropType< TypeTag, Properties::MaterialLaw > MaterialLaw
Definition: BlackoilModel.hpp:73
static constexpr int contiPolymerMWEqIdx
Definition: BlackoilModel.hpp:83
static constexpr bool has_micp_
Definition: BlackoilModel.hpp:329
const Simulator & simulator() const
Definition: BlackoilModel.hpp:256
BlackoilWellModel< TypeTag > & well_model_
Definition: BlackoilModel.hpp:335
BlackoilWellModel< TypeTag > & wellModel()
return the StandardWells object
Definition: BlackoilModel.hpp:282
std::vector< std::vector< Scalar > > residual_norms_history_
Definition: BlackoilModel.hpp:342
typename SparseMatrixAdapter::IstlMatrix Mat
Definition: BlackoilModel.hpp:106
Dune::FieldVector< Scalar, numEq > VectorBlockType
Definition: BlackoilModel.hpp:104
static constexpr bool has_solvent_
Definition: BlackoilModel.hpp:322
const BlackoilWellModel< TypeTag > & wellModel() const
Definition: BlackoilModel.hpp:286
BVector dx_old_
Definition: BlackoilModel.hpp:344
Simulator & simulator_
Definition: BlackoilModel.hpp:320
static constexpr bool has_polymer_
Definition: BlackoilModel.hpp:324
static constexpr bool has_brine_
Definition: BlackoilModel.hpp:328
void solveJacobianSystem(BVector &x)
Definition: BlackoilModel_impl.hpp:458
static constexpr int oxygenConcentrationIdx
Definition: BlackoilModel.hpp:99
ComponentName compNames_
Definition: BlackoilModel.hpp:347
SimulatorReportSingle prepareStep(const SimulatorTimerInterface &timer)
Definition: BlackoilModel_impl.hpp:86
static constexpr int ureaConcentrationIdx
Definition: BlackoilModel.hpp:100
GetPropType< TypeTag, Properties::Grid > Grid
Definition: BlackoilModel.hpp:65
void convergencePerCell(const std::vector< Scalar > &B_avg, const double dt, const double tol_cnv, const double tol_cnv_energy, const int iteration)
Compute the number of Newtons required by each cell in order to satisfy the solution change convergen...
Definition: BlackoilModel_impl.hpp:941
Dune::BlockVector< VectorBlockType > BVector
Definition: BlackoilModel.hpp:107
GetPropType< TypeTag, Properties::Scalar > Scalar
Definition: BlackoilModel.hpp:75
static constexpr int saltConcentrationIdx
Definition: BlackoilModel.hpp:97
static constexpr int zFractionIdx
Definition: BlackoilModel.hpp:92
static constexpr int contiBrineEqIdx
Definition: BlackoilModel.hpp:85
Scalar relativeChange() const
Definition: BlackoilModel_impl.hpp:371
static constexpr int contiZfracEqIdx
Definition: BlackoilModel.hpp:80
static constexpr int contiSolventEqIdx
Definition: BlackoilModel.hpp:79
GetPropType< TypeTag, Properties::MaterialLawParams > MaterialLawParams
Definition: BlackoilModel.hpp:74
static constexpr int contiPolymerEqIdx
Definition: BlackoilModel.hpp:81
const EclipseState & eclState() const
Definition: BlackoilModel.hpp:128
bool terminal_output_
Whether we print something to std::cout.
Definition: BlackoilModel.hpp:338
Class for handling the blackoil well model.
Definition: BlackoilWellModel.hpp:102
Definition: ComponentName.hpp:34
Definition: ConvergenceReport.hpp:38
Interface class for SimulatorTimer objects, to be improved.
Definition: SimulatorTimerInterface.hpp:34
Dune::Communication< MPIComm > Communication
Definition: ParallelCommunication.hpp:30
Definition: blackoilboundaryratevector.hh:39
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
Solver parameters for the BlackoilModel.
Definition: BlackoilModelParameters.hpp:179
Scalar max_residual_allowed_
Absolute max limit for residuals.
Definition: BlackoilModelParameters.hpp:192
Definition: SimulatorReport.hpp:122
A struct for returning timing data from a simulator to its caller.
Definition: SimulatorReport.hpp:34