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
43
45
46#include <memory>
47#include <tuple>
48#include <vector>
49
50#include <fmt/format.h>
51
52namespace Opm {
53
60template <class TypeTag>
62{
63public:
64 // --------- Types and enums ---------
78
79 static constexpr int numEq = Indices::numEq;
80 static constexpr int contiSolventEqIdx = Indices::contiSolventEqIdx;
81 static constexpr int contiZfracEqIdx = Indices::contiZfracEqIdx;
82 static constexpr int contiPolymerEqIdx = Indices::contiPolymerEqIdx;
83 static constexpr int contiEnergyEqIdx = Indices::contiEnergyEqIdx;
84 static constexpr int contiPolymerMWEqIdx = Indices::contiPolymerMWEqIdx;
85 static constexpr int contiFoamEqIdx = Indices::contiFoamEqIdx;
86 static constexpr int contiBrineEqIdx = Indices::contiBrineEqIdx;
87 static constexpr int contiMicrobialEqIdx = Indices::contiMicrobialEqIdx;
88 static constexpr int contiOxygenEqIdx = Indices::contiOxygenEqIdx;
89 static constexpr int contiUreaEqIdx = Indices::contiUreaEqIdx;
90 static constexpr int contiBiofilmEqIdx = Indices::contiBiofilmEqIdx;
91 static constexpr int contiCalciteEqIdx = Indices::contiCalciteEqIdx;
92 static constexpr int solventSaturationIdx = Indices::solventSaturationIdx;
93 static constexpr int zFractionIdx = Indices::zFractionIdx;
94 static constexpr int polymerConcentrationIdx = Indices::polymerConcentrationIdx;
95 static constexpr int polymerMoleWeightIdx = Indices::polymerMoleWeightIdx;
96 static constexpr int temperatureIdx = Indices::temperatureIdx;
97 static constexpr int foamConcentrationIdx = Indices::foamConcentrationIdx;
98 static constexpr int saltConcentrationIdx = Indices::saltConcentrationIdx;
99 static constexpr int microbialConcentrationIdx = Indices::microbialConcentrationIdx;
100 static constexpr int oxygenConcentrationIdx = Indices::oxygenConcentrationIdx;
101 static constexpr int ureaConcentrationIdx = Indices::ureaConcentrationIdx;
102 static constexpr int biofilmConcentrationIdx = Indices::biofilmConcentrationIdx;
103 static constexpr int calciteConcentrationIdx = Indices::calciteConcentrationIdx;
104
105 using VectorBlockType = Dune::FieldVector<Scalar, numEq>;
106 using MatrixBlockType = typename SparseMatrixAdapter::MatrixBlock;
107 using Mat = typename SparseMatrixAdapter::IstlMatrix;
108 using BVector = Dune::BlockVector<VectorBlockType>;
109
111
112 // --------- Public methods ---------
113
125 const ModelParameters& param,
126 BlackoilWellModel<TypeTag>& well_model,
127 const bool terminal_output);
128
129 bool isParallel() const
130 { return grid_.comm().size() > 1; }
131
132 const EclipseState& eclState() const
133 { return simulator_.vanguard().eclState(); }
134
138
140 const int iteration,
141 const int minIter,
142 const int maxIter,
143 const SimulatorTimerInterface& timer);
144
152 template <class NonlinearSolverType>
153 SimulatorReportSingle nonlinearIteration(const int iteration,
154 const SimulatorTimerInterface& timer,
155 NonlinearSolverType& nonlinear_solver);
156
157 template <class NonlinearSolverType>
159 const SimulatorTimerInterface& timer,
160 NonlinearSolverType& nonlinear_solver);
161
166
169 const int iterationIdx);
170
171 // compute the "relative" change of the solution between time steps
172 Scalar relativeChange() const;
173
176 { return simulator_.model().newtonMethod().linearSolver().iterations (); }
177
178 // Obtain reference to linear solver setup time
180 { return linear_solve_setup_time_; }
181
185
187 void updateSolution(const BVector& dx);
188
191 { return terminal_output_; }
192
193 std::tuple<Scalar,Scalar>
195 const Scalar pvSumLocal,
196 const Scalar numAquiferPvSumLocal,
197 std::vector<Scalar>& R_sum,
198 std::vector<Scalar>& maxCoeff,
199 std::vector<Scalar>& B_avg);
200
204 std::pair<Scalar,Scalar>
205 localConvergenceData(std::vector<Scalar>& R_sum,
206 std::vector<Scalar>& maxCoeff,
207 std::vector<Scalar>& B_avg,
208 std::vector<int>& maxCoeffCell);
209
213 std::pair<std::vector<double>, std::vector<int>>
214 characteriseCnvPvSplit(const std::vector<Scalar>& B_avg, const double dt);
215
216 void updateTUNING(const Tuning& tuning);
217
219 getReservoirConvergence(const double reportTime,
220 const double dt,
221 const int iteration,
222 const int maxIter,
223 std::vector<Scalar>& B_avg,
224 std::vector<Scalar>& residual_norms);
225
234 const int iteration,
235 const int maxIter,
236 std::vector<Scalar>& residual_norms);
237
239 int numPhases() const
240 { return phaseUsage_.num_phases; }
241
243 template<class T>
244 std::vector<std::vector<Scalar> >
245 computeFluidInPlace(const T&, const std::vector<int>& fipnum) const
246 { return this->computeFluidInPlace(fipnum); }
247
249 std::vector<std::vector<Scalar> >
250 computeFluidInPlace(const std::vector<int>& /*fipnum*/) const;
251
252 const Simulator& simulator() const
253 { return simulator_; }
254
256 { return simulator_; }
257
260 { return failureReport_; }
261
264
266 const std::vector<SimulatorReport>& domainAccumulatedReports() const;
267
269 void writeNonlinearIterationsPerCell(const std::filesystem::path& odir) const;
270
271 const std::vector<StepReport>& stepReports() const
272 { return convergence_reports_; }
273
274 void writePartitions(const std::filesystem::path& odir) const;
275
279 { return well_model_; }
280
282 wellModel() const
283 { return well_model_; }
284
286 { simulator_.problem().beginEpisode(); }
287
289 { simulator_.problem().endEpisode(); }
290
291 template<class FluidState, class Residual>
292 void getMaxCoeff(const unsigned cell_idx,
293 const IntensiveQuantities& intQuants,
294 const FluidState& fs,
295 const Residual& modelResid,
296 const Scalar pvValue,
297 std::vector<Scalar>& B_avg,
298 std::vector<Scalar>& R_sum,
299 std::vector<Scalar>& maxCoeff,
300 std::vector<int>& maxCoeffCell);
301
303 const ModelParameters& param() const
304 { return param_; }
305
308 { return compNames_; }
309
311 bool hasNlddSolver() const
312 { return nlddSolver_ != nullptr; }
313
314protected:
315 // --------- Data members ---------
317 const Grid& grid_;
319 static constexpr bool has_solvent_ = getPropValue<TypeTag, Properties::EnableSolvent>();
320 static constexpr bool has_extbo_ = getPropValue<TypeTag, Properties::EnableExtbo>();
321 static constexpr bool has_polymer_ = getPropValue<TypeTag, Properties::EnablePolymer>();
322 static constexpr bool has_polymermw_ = getPropValue<TypeTag, Properties::EnablePolymerMW>();
323 static constexpr bool has_energy_ = getPropValue<TypeTag, Properties::EnableEnergy>();
324 static constexpr bool has_foam_ = getPropValue<TypeTag, Properties::EnableFoam>();
325 static constexpr bool has_brine_ = getPropValue<TypeTag, Properties::EnableBrine>();
326 static constexpr bool has_micp_ = getPropValue<TypeTag, Properties::EnableMICP>();
327
330
331 // Well Model
333
337 long int global_nc_;
338
339 std::vector<std::vector<Scalar>> residual_norms_history_;
342
343 std::vector<StepReport> convergence_reports_;
345
346 std::unique_ptr<BlackoilModelNldd<TypeTag>> nlddSolver_;
348
349private:
350 Scalar dpMaxRel() const { return param_.dp_max_rel_; }
351 Scalar dsMax() const { return param_.ds_max_; }
352 Scalar drMaxRel() const { return param_.dr_max_rel_; }
353 Scalar maxResidualAllowed() const { return param_.max_residual_allowed_; }
354 double linear_solve_setup_time_;
355 std::vector<bool> wasSwitched_;
356};
357
358} // namespace Opm
359
361
362#endif // OPM_BLACKOILMODEL_HEADER_INCLUDED
Implementation of penalty cards for three-phase black oil.
Definition: BlackoilModelConvergenceMonitor.hpp:35
Definition: BlackoilModel.hpp:62
static constexpr int contiEnergyEqIdx
Definition: BlackoilModel.hpp:83
SimulatorReportSingle failureReport_
Definition: BlackoilModel.hpp:329
BlackoilModel(Simulator &simulator, const ModelParameters &param, BlackoilWellModel< TypeTag > &well_model, const bool terminal_output)
Definition: BlackoilModel_impl.hpp:53
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:245
static constexpr bool has_energy_
Definition: BlackoilModel.hpp:323
ModelParameters param_
Definition: BlackoilModel.hpp:328
int linearIterationsLastSolve() const
Number of linear iterations used in last call to solveJacobianSystem().
Definition: BlackoilModel.hpp:175
SimulatorReportSingle nonlinearIteration(const int iteration, const SimulatorTimerInterface &timer, NonlinearSolverType &nonlinear_solver)
Definition: BlackoilModel_impl.hpp:217
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:1000
static constexpr bool has_extbo_
Definition: BlackoilModel.hpp:320
SimulatorReportSingle assembleReservoir(const SimulatorTimerInterface &, const int iterationIdx)
Assemble the residual and Jacobian of the nonlinear system.
Definition: BlackoilModel_impl.hpp:360
GetPropType< TypeTag, Properties::SparseMatrixAdapter > SparseMatrixAdapter
Definition: BlackoilModel.hpp:69
ConvergenceReport getConvergence(const SimulatorTimerInterface &timer, const int iteration, const int maxIter, std::vector< Scalar > &residual_norms)
Definition: BlackoilModel_impl.hpp:942
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:744
GetPropType< TypeTag, Properties::IntensiveQuantities > IntensiveQuantities
Definition: BlackoilModel.hpp:68
static constexpr bool has_foam_
Definition: BlackoilModel.hpp:324
static constexpr int foamConcentrationIdx
Definition: BlackoilModel.hpp:97
Scalar current_relaxation_
Definition: BlackoilModel.hpp:340
GetPropType< TypeTag, Properties::FluidSystem > FluidSystem
Definition: BlackoilModel.hpp:72
static constexpr int microbialConcentrationIdx
Definition: BlackoilModel.hpp:99
static constexpr int numEq
Definition: BlackoilModel.hpp:79
const std::vector< SimulatorReport > & domainAccumulatedReports() const
return the statistics of local solves accumulated for each domain on this rank
Definition: BlackoilModel_impl.hpp:990
static constexpr bool has_polymermw_
Definition: BlackoilModel.hpp:322
static constexpr int biofilmConcentrationIdx
Definition: BlackoilModel.hpp:102
static constexpr int contiFoamEqIdx
Definition: BlackoilModel.hpp:85
static constexpr int polymerMoleWeightIdx
Definition: BlackoilModel.hpp:95
void endReportStep()
Definition: BlackoilModel.hpp:288
BlackoilModelConvergenceMonitor< Scalar > conv_monitor_
Definition: BlackoilModel.hpp:347
static constexpr int contiCalciteEqIdx
Definition: BlackoilModel.hpp:91
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:607
GetPropType< TypeTag, Properties::ElementContext > ElementContext
Definition: BlackoilModel.hpp:67
bool isParallel() const
Definition: BlackoilModel.hpp:129
const std::vector< StepReport > & stepReports() const
Definition: BlackoilModel.hpp:271
std::vector< StepReport > convergence_reports_
Definition: BlackoilModel.hpp:343
static constexpr int contiMicrobialEqIdx
Definition: BlackoilModel.hpp:87
bool terminalOutputEnabled() const
Return true if output to cout is wanted.
Definition: BlackoilModel.hpp:190
double & linearSolveSetupTime()
Definition: BlackoilModel.hpp:179
int numPhases() const
The number of active fluid phases in the model.
Definition: BlackoilModel.hpp:239
const SimulatorReport & localAccumulatedReports() const
return the statistics of local solves accumulated for this rank
Definition: BlackoilModel_impl.hpp:979
SimulatorReportSingle afterStep(const SimulatorTimerInterface &)
Definition: BlackoilModel_impl.hpp:347
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:547
static constexpr int contiBiofilmEqIdx
Definition: BlackoilModel.hpp:90
const Grid & grid_
Definition: BlackoilModel.hpp:317
static constexpr int contiUreaEqIdx
Definition: BlackoilModel.hpp:89
static constexpr int calciteConcentrationIdx
Definition: BlackoilModel.hpp:103
static constexpr int temperatureIdx
Definition: BlackoilModel.hpp:96
GetPropType< TypeTag, Properties::SolutionVector > SolutionVector
Definition: BlackoilModel.hpp:70
const ComponentName & compNames() const
Returns const reference to component names.
Definition: BlackoilModel.hpp:307
GetPropType< TypeTag, Properties::PrimaryVariables > PrimaryVariables
Definition: BlackoilModel.hpp:71
void initialLinearization(SimulatorReportSingle &report, const int iteration, const int minIter, const int maxIter, const SimulatorTimerInterface &timer)
Definition: BlackoilModel_impl.hpp:158
void updateSolution(const BVector &dx)
Apply an update to the primary variables.
Definition: BlackoilModel_impl.hpp:523
long int global_nc_
The number of cells of the global grid.
Definition: BlackoilModel.hpp:337
GetPropType< TypeTag, Properties::Indices > Indices
Definition: BlackoilModel.hpp:73
void updateTUNING(const Tuning &tuning)
Definition: BlackoilModel_impl.hpp:730
typename SparseMatrixAdapter::MatrixBlock MatrixBlockType
Definition: BlackoilModel.hpp:106
const ModelParameters & param() const
Returns const reference to model parameters.
Definition: BlackoilModel.hpp:303
static constexpr int contiOxygenEqIdx
Definition: BlackoilModel.hpp:88
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:1037
std::unique_ptr< BlackoilModelNldd< TypeTag > > nlddSolver_
Non-linear DD solver.
Definition: BlackoilModel.hpp:346
SimulatorReportSingle nonlinearIterationNewton(const int iteration, const SimulatorTimerInterface &timer, NonlinearSolverType &nonlinear_solver)
Definition: BlackoilModel_impl.hpp:255
void writePartitions(const std::filesystem::path &odir) const
Definition: BlackoilModel_impl.hpp:1010
static constexpr int solventSaturationIdx
Definition: BlackoilModel.hpp:92
bool hasNlddSolver() const
Returns true if an NLDD solver exists.
Definition: BlackoilModel.hpp:311
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:658
void beginReportStep()
Definition: BlackoilModel.hpp:285
static constexpr int polymerConcentrationIdx
Definition: BlackoilModel.hpp:94
GetPropType< TypeTag, Properties::Simulator > Simulator
Definition: BlackoilModel.hpp:65
const SimulatorReportSingle & failureReport() const
return the statistics if the nonlinearIteration() method failed
Definition: BlackoilModel.hpp:259
Simulator & simulator()
Definition: BlackoilModel.hpp:255
GetPropType< TypeTag, Properties::MaterialLaw > MaterialLaw
Definition: BlackoilModel.hpp:74
static constexpr int contiPolymerMWEqIdx
Definition: BlackoilModel.hpp:84
static constexpr bool has_micp_
Definition: BlackoilModel.hpp:326
const Simulator & simulator() const
Definition: BlackoilModel.hpp:252
BlackoilWellModel< TypeTag > & well_model_
Definition: BlackoilModel.hpp:332
BlackoilWellModel< TypeTag > & wellModel()
return the StandardWells object
Definition: BlackoilModel.hpp:278
std::vector< std::vector< Scalar > > residual_norms_history_
Definition: BlackoilModel.hpp:339
typename SparseMatrixAdapter::IstlMatrix Mat
Definition: BlackoilModel.hpp:107
Dune::FieldVector< Scalar, numEq > VectorBlockType
Definition: BlackoilModel.hpp:105
static constexpr bool has_solvent_
Definition: BlackoilModel.hpp:319
const BlackoilWellModel< TypeTag > & wellModel() const
Definition: BlackoilModel.hpp:282
BVector dx_old_
Definition: BlackoilModel.hpp:341
Simulator & simulator_
Definition: BlackoilModel.hpp:316
static constexpr bool has_polymer_
Definition: BlackoilModel.hpp:321
const PhaseUsage phaseUsage_
Definition: BlackoilModel.hpp:318
static constexpr bool has_brine_
Definition: BlackoilModel.hpp:325
void solveJacobianSystem(BVector &x)
Definition: BlackoilModel_impl.hpp:461
static constexpr int oxygenConcentrationIdx
Definition: BlackoilModel.hpp:100
ComponentName compNames_
Definition: BlackoilModel.hpp:344
SimulatorReportSingle prepareStep(const SimulatorTimerInterface &timer)
Definition: BlackoilModel_impl.hpp:89
static constexpr int ureaConcentrationIdx
Definition: BlackoilModel.hpp:101
GetPropType< TypeTag, Properties::Grid > Grid
Definition: BlackoilModel.hpp:66
Dune::BlockVector< VectorBlockType > BVector
Definition: BlackoilModel.hpp:108
GetPropType< TypeTag, Properties::Scalar > Scalar
Definition: BlackoilModel.hpp:76
static constexpr int saltConcentrationIdx
Definition: BlackoilModel.hpp:98
static constexpr int zFractionIdx
Definition: BlackoilModel.hpp:93
static constexpr int contiBrineEqIdx
Definition: BlackoilModel.hpp:86
Scalar relativeChange() const
Definition: BlackoilModel_impl.hpp:374
static constexpr int contiZfracEqIdx
Definition: BlackoilModel.hpp:81
static constexpr int contiSolventEqIdx
Definition: BlackoilModel.hpp:80
GetPropType< TypeTag, Properties::MaterialLawParams > MaterialLawParams
Definition: BlackoilModel.hpp:75
static constexpr int contiPolymerEqIdx
Definition: BlackoilModel.hpp:82
const EclipseState & eclState() const
Definition: BlackoilModel.hpp:132
bool terminal_output_
Whether we print something to std::cout.
Definition: BlackoilModel.hpp:335
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:177
Scalar max_residual_allowed_
Absolute max limit for residuals.
Definition: BlackoilModelParameters.hpp:190
Definition: BlackoilPhases.hpp:46
int num_phases
Definition: BlackoilPhases.hpp:54
Definition: SimulatorReport.hpp:122
A struct for returning timing data from a simulator to its caller.
Definition: SimulatorReport.hpp:34