24#ifndef OPM_BLACKOILMODEL_HEADER_INCLUDED
25#define OPM_BLACKOILMODEL_HEADER_INCLUDED
27#include <opm/simulators/aquifers/BlackoilAquiferModel.hpp>
29#include <opm/simulators/flow/BlackoilModelConvergenceMonitor.hpp>
30#include <opm/simulators/flow/BlackoilModelNldd.hpp>
31#include <opm/simulators/flow/BlackoilModelProperties.hpp>
33#include <opm/simulators/flow/RSTConv.hpp>
35#include <opm/simulators/linalg/ISTLSolver.hpp>
37#include <opm/simulators/timestepping/ConvergenceReport.hpp>
38#include <opm/simulators/timestepping/SimulatorReport.hpp>
39#include <opm/simulators/timestepping/SimulatorTimer.hpp>
41#include <opm/simulators/utils/BlackoilPhases.hpp>
42#include <opm/simulators/utils/ComponentName.hpp>
44#include <opm/simulators/wells/BlackoilWellModel.hpp>
50#include <fmt/format.h>
60template <
class TypeTag>
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;
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>;
129 bool isParallel()
const
130 {
return grid_.comm().size() > 1; }
132 const EclipseState& eclState()
const
133 {
return simulator_.vanguard().eclState(); }
152 template <
class NonlinearSolverType>
157 template <
class NonlinearSolverType>
172 Scalar relativeChange()
const;
176 {
return simulator_.model().newtonMethod().linearSolver().iterations (); }
179 double& linearSolveSetupTime()
180 {
return linear_solve_setup_time_; }
193 std::tuple<Scalar,Scalar>
194 convergenceReduction(Parallel::Communication comm,
197 std::vector<Scalar>&
R_sum,
199 std::vector<Scalar>&
B_avg);
204 std::pair<Scalar,Scalar>
207 std::vector<Scalar>&
B_avg,
213 std::pair<std::vector<double>, std::vector<int>>
219 getReservoirConvergence(
const double reportTime,
223 std::vector<Scalar>&
B_avg,
240 {
return phaseUsage_.num_phases; }
244 std::vector<std::vector<Scalar> >
249 std::vector<std::vector<Scalar> >
253 {
return simulator_; }
255 Simulator& simulator()
256 {
return simulator_; }
260 {
return failureReport_; }
271 const std::vector<StepReport>& stepReports()
const
272 {
return convergence_reports_; }
274 void writePartitions(
const std::filesystem::path&
odir)
const;
279 {
return well_model_; }
283 {
return well_model_; }
285 void beginReportStep()
286 { simulator_.problem().beginEpisode(); }
289 { simulator_.problem().endEpisode(); }
291 template<
class Flu
idState,
class Res
idual>
292 void getMaxCoeff(
const unsigned cell_idx,
293 const IntensiveQuantities& intQuants,
294 const FluidState& fs,
297 std::vector<Scalar>&
B_avg,
298 std::vector<Scalar>&
R_sum,
308 {
return compNames_; }
328 ModelParameters param_;
339 std::vector<std::vector<Scalar>> residual_norms_history_;
340 Scalar current_relaxation_;
343 std::vector<StepReport> convergence_reports_;
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_;
360#include <opm/simulators/flow/BlackoilModel_impl.hpp>
Implementation of penalty cards for three-phase black oil.
Definition BlackoilModelConvergenceMonitor.hpp:35
A model implementation for three-phase black oil.
Definition BlackoilModel.hpp:62
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
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)
Called once per nonlinear iteration.
Definition BlackoilModel_impl.hpp:219
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:1003
SimulatorReportSingle assembleReservoir(const SimulatorTimerInterface &, const int iterationIdx)
Assemble the residual and Jacobian of the nonlinear system.
Definition BlackoilModel_impl.hpp:363
ConvergenceReport getConvergence(const SimulatorTimerInterface &timer, const int iteration, const int maxIter, std::vector< Scalar > &residual_norms)
Compute convergence based on total mass balance (tol_mb) and maximum residual mass balance (tol_cnv).
Definition BlackoilModel_impl.hpp:945
const std::vector< SimulatorReport > & domainAccumulatedReports() const
return the statistics of local solves accumulated for each domain on this rank
Definition BlackoilModel_impl.hpp:993
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:610
bool terminalOutputEnabled() const
Return true if output to cout is wanted.
Definition BlackoilModel.hpp:190
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:982
SimulatorReportSingle afterStep(const SimulatorTimerInterface &)
Called once after each time step.
Definition BlackoilModel_impl.hpp:350
const ComponentName & compNames() const
Returns const reference to component names.
Definition BlackoilModel.hpp:307
void updateSolution(const BVector &dx)
Apply an update to the primary variables.
Definition BlackoilModel_impl.hpp:526
long int global_nc_
The number of cells of the global grid.
Definition BlackoilModel.hpp:337
const ModelParameters & param() const
Returns const reference to model parameters.
Definition BlackoilModel.hpp:303
std::unique_ptr< BlackoilModelNldd< TypeTag > > nlddSolver_
Non-linear DD solver.
Definition BlackoilModel.hpp:346
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:661
const SimulatorReportSingle & failureReport() const
return the statistics if the nonlinearIteration() method failed
Definition BlackoilModel.hpp:259
BlackoilWellModel< TypeTag > & wellModel()
return the StandardWells object
Definition BlackoilModel.hpp:278
void solveJacobianSystem(BVector &x)
Solve the Jacobian system Jx = r where J is the Jacobian and r is the residual.
Definition BlackoilModel_impl.hpp:464
SimulatorReportSingle prepareStep(const SimulatorTimerInterface &timer)
Called once before each time step.
Definition BlackoilModel_impl.hpp:89
bool terminal_output_
Whether we print something to std::cout.
Definition BlackoilModel.hpp:335
Class for handling the blackoil well model.
Definition BlackoilWellModel.hpp:94
Definition ComponentName.hpp:34
Represents the convergence status of the whole simulator, to make it possible to query and store the ...
Definition ConvergenceReport.hpp:38
Interface class for SimulatorTimer objects, to be improved.
Definition SimulatorTimerInterface.hpp:34
Manages the initializing and running of time dependent problems.
Definition simulator.hh:97
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242
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:235
Solver parameters for the BlackoilModel.
Definition BlackoilModelParameters.hpp:174
Definition BlackoilPhases.hpp:46
A struct for returning timing data from a simulator to its caller.
Definition SimulatorReport.hpp:34
Definition SimulatorReport.hpp:119