29#ifndef OPM_INIT_STATE_EQUIL_HPP
30#define OPM_INIT_STATE_EQUIL_HPP
34#include <opm/material/common/Tabulated1DFunction.hpp>
35#include <opm/material/fluidstates/SimpleModularFluidState.hpp>
37#include <opm/simulators/utils/ParallelCommunication.hpp>
50class NumericalAquifers;
62namespace Miscibility {
template<
class Scalar>
class RsFunction; }
65template <
class Scalar,
class RHS>
70 const std::array<Scalar,2>& span,
74 Scalar operator()(
const Scalar x)
const;
78 std::array<Scalar,2> span_;
79 std::vector<Scalar> y_;
80 std::vector<Scalar> f_;
82 Scalar stepsize()
const;
85namespace PhasePressODE {
86template <
class Flu
idSystem>
89 using Scalar =
typename FluidSystem::Scalar;
90 using TabulatedFunction = Tabulated1DFunction<Scalar>;
93 Water(
const TabulatedFunction& tempVdTable,
94 const TabulatedFunction& saltVdTable,
95 const int pvtRegionIdx,
96 const Scalar normGrav);
98 Scalar operator()(
const Scalar depth,
99 const Scalar press)
const;
102 const TabulatedFunction& tempVdTable_;
103 const TabulatedFunction& saltVdTable_;
104 const int pvtRegionIdx_;
107 Scalar density(
const Scalar depth,
108 const Scalar press)
const;
111template <
class Flu
idSystem,
class RS>
114 using Scalar =
typename FluidSystem::Scalar;
115 using TabulatedFunction = Tabulated1DFunction<Scalar>;
118 Oil(
const TabulatedFunction& tempVdTable,
120 const int pvtRegionIdx,
121 const Scalar normGrav);
123 Scalar operator()(
const Scalar depth,
124 const Scalar press)
const;
127 const TabulatedFunction& tempVdTable_;
129 const int pvtRegionIdx_;
132 Scalar density(
const Scalar depth,
133 const Scalar press)
const;
136template <
class Flu
idSystem,
class RV,
class RVW>
139 using Scalar =
typename FluidSystem::Scalar;
140 using TabulatedFunction = Tabulated1DFunction<Scalar>;
143 Gas(
const TabulatedFunction& tempVdTable,
146 const int pvtRegionIdx,
147 const Scalar normGrav);
149 Scalar operator()(
const Scalar depth,
150 const Scalar press)
const;
153 const TabulatedFunction& tempVdTable_;
156 const int pvtRegionIdx_;
159 Scalar density(
const Scalar depth,
160 const Scalar press)
const;
165template <
class Flu
idSystem,
class Region>
169 using Scalar =
typename FluidSystem::Scalar;
170 using VSpan = std::array<Scalar, 2>;
181 const int samplePoints = 2000);
211 void equilibrate(
const Region& reg,
230 Scalar
oil(
const Scalar depth)
const;
239 Scalar
gas(
const Scalar depth)
const;
248 Scalar
water(
const Scalar depth)
const;
252 class PressureFunction
260 explicit PressureFunction(
const ODE& ode,
265 PressureFunction(
const PressureFunction& rhs);
267 PressureFunction(PressureFunction&& rhs) =
default;
269 PressureFunction& operator=(
const PressureFunction& rhs);
271 PressureFunction& operator=(PressureFunction&& rhs);
273 Scalar value(
const Scalar depth)
const;
276 enum Direction : std::size_t { Up, Down, NumDir };
279 using DistrPtr = std::unique_ptr<Distribution>;
282 std::array<DistrPtr, Direction::NumDir> value_;
286 FluidSystem,
typename Region::CalcDissolution
290 FluidSystem,
typename Region::CalcEvaporation,
typename Region::CalcWaterEvaporation
295 using OPress = PressureFunction<OilPressODE>;
296 using GPress = PressureFunction<GasPressODE>;
297 using WPress = PressureFunction<WatPressODE>;
300 (
const Region&,
const VSpan&);
305 std::unique_ptr<OPress> oil_{};
306 std::unique_ptr<GPress> gas_{};
307 std::unique_ptr<WPress> wat_{};
309 template <
typename PressFunc>
310 void checkPtr(
const PressFunc* phasePress,
311 const std::string& phaseName)
const;
313 Strategy selectEquilibrationStrategy(
const Region& reg)
const;
317 void equil_WOG(
const Region& reg,
const VSpan& span);
318 void equil_GOW(
const Region& reg,
const VSpan& span);
319 void equil_OWG(
const Region& reg,
const VSpan& span);
321 void makeOilPressure(
const typename OPress::InitCond& ic,
325 void makeGasPressure(
const typename GPress::InitCond& ic,
329 void makeWatPressure(
const typename WPress::InitCond& ic,
337template<
class Scalar>
345 this->oil += a * rhs.oil;
346 this->gas += a * rhs.gas;
347 this->water += a * rhs.water;
363 this->oil = this->gas = this->water = 0.0;
384template <
class MaterialLawManager,
class Flu
idSystem,
class Region,
typename CellID>
388 using Scalar =
typename FluidSystem::Scalar;
408 const std::vector<Scalar>& swatInit);
450 struct EvaluationPoint {
451 const Position* position{
nullptr};
452 const Region* region {
nullptr};
453 const PTable* ptable {
nullptr};
459 using FluidState = ::Opm::
460 SimpleModularFluidState<Scalar, 3, 3,
472 using MaterialLaw =
typename MaterialLawManager::MaterialLaw;
475 using PhaseIdx = std::remove_cv_t<
476 std::remove_reference_t<
decltype(FluidSystem::oilPhaseIdx)>
480 MaterialLawManager& matLawMgr_;
483 const std::vector<Scalar>& swatInit_;
486 PhaseQuantityValue<Scalar> sat_;
489 PhaseQuantityValue<Scalar> press_;
492 EvaluationPoint evalPt_;
495 FluidState fluidState_;
498 std::array<Scalar, FluidSystem::numPhases> matLawCapPress_;
509 void setEvaluationPoint(
const Position& x,
516 void initializePhaseQuantities();
535 void deriveWaterSat();
539 void fixUnphysicalTransition();
543 void accountForScaledSaturations();
559 std::pair<Scalar, bool> applySwatInit(
const Scalar pcow);
573 std::pair<Scalar, bool> applySwatInit(
const Scalar pc,
const Scalar sw);
577 void computeMaterialLawCapPress();
581 Scalar materialLawCapPressGasOil()
const;
585 Scalar materialLawCapPressOilWater()
const;
589 Scalar materialLawCapPressGasWater()
const;
598 bool isConstCapPress(
const PhaseIdx phaseIdx)
const;
605 bool isOverlappingTransition()
const;
626 Scalar fromDepthTable(
const Scalar contactdepth,
627 const PhaseIdx phasePos,
628 const bool isincr)
const;
647 Scalar invertCapPress(
const Scalar pc,
648 const PhaseIdx phasePos,
649 const bool isincr)
const;
652 PhaseIdx oilPos()
const
654 return FluidSystem::oilPhaseIdx;
658 PhaseIdx gasPos()
const
660 return FluidSystem::gasPhaseIdx;
664 PhaseIdx waterPos()
const
666 return FluidSystem::waterPhaseIdx;
672template <
typename CellRange,
class Scalar>
673void verticalExtent(
const CellRange& cells,
674 const std::vector<std::pair<Scalar, Scalar>>& cellZMinMax,
675 const Parallel::Communication& comm,
676 std::array<Scalar,2>& span);
678template <
class Scalar,
class Element>
679std::pair<Scalar,Scalar> cellZMinMax(
const Element& element);
683namespace DeckDependent {
685template<
class FluidSystem,
689 class CartesianIndexMapper>
692 using Element =
typename GridView::template Codim<0>::Entity;
693 using Scalar =
typename FluidSystem::Scalar;
695 template<
class MaterialLawManager>
697 const EclipseState& eclipseState,
699 const GridView& gridView,
700 const CartesianIndexMapper& cartMapper,
702 const int num_pressure_points = 2000,
703 const bool applySwatInit =
true);
705 using Vec = std::vector<Scalar>;
706 using PVec = std::vector<Vec>;
708 const Vec& temperature()
const {
return temperature_; }
709 const Vec& saltConcentration()
const {
return saltConcentration_; }
710 const Vec& saltSaturation()
const {
return saltSaturation_; }
711 const PVec& press()
const {
return pp_; }
712 const PVec& saturation()
const {
return sat_; }
713 const Vec& rs()
const {
return rs_; }
714 const Vec& rv()
const {
return rv_; }
715 const Vec& rvw()
const {
return rvw_; }
718 template <
class RMap>
719 void updateInitialTemperature_(
const EclipseState& eclState,
const RMap& reg);
721 template <
class RMap>
722 void updateInitialSaltConcentration_(
const EclipseState& eclState,
const RMap& reg);
724 template <
class RMap>
725 void updateInitialSaltSaturation_(
const EclipseState& eclState,
const RMap& reg);
727 void updateCellProps_(
const GridView& gridView,
728 const NumericalAquifers& aquifer);
730 void applyNumericalAquifers_(
const GridView& gridView,
731 const NumericalAquifers& aquifer,
732 const bool co2store_or_h2store);
735 void setRegionPvtIdx(
const EclipseState& eclState,
const RMap& reg);
737 template <
class RMap,
class MaterialLawManager,
class Comm>
738 void calcPressSatRsRv(
const RMap& reg,
739 const std::vector<EquilRecord>& rec,
740 MaterialLawManager& materialLawManager,
744 template <
class CellRange,
class EquilibrationMethod>
745 void cellLoop(
const CellRange& cells,
746 EquilibrationMethod&& eqmethod);
748 template <
class CellRange,
class PressTable,
class PhaseSat>
749 void equilibrateCellCentres(
const CellRange& cells,
751 const PressTable& ptable,
754 template <
class CellRange,
class PressTable,
class PhaseSat>
755 void equilibrateHorizontal(
const CellRange& cells,
758 const PressTable& ptable,
761 std::vector< std::shared_ptr<Miscibility::RsFunction<Scalar>> > rsFunc_;
762 std::vector< std::shared_ptr<Miscibility::RsFunction<Scalar>> > rvFunc_;
763 std::vector< std::shared_ptr<Miscibility::RsFunction<Scalar>> > rvwFunc_;
764 using TabulatedFunction = Tabulated1DFunction<Scalar>;
765 std::vector<TabulatedFunction> tempVdTable_;
766 std::vector<TabulatedFunction> saltVdTable_;
767 std::vector<TabulatedFunction> saltpVdTable_;
768 std::vector<int> regionPvtIdx_;
770 Vec saltConcentration_;
777 const CartesianIndexMapper& cartesianIndexMapper_;
779 Vec cellCenterDepth_;
780 std::vector<std::pair<Scalar,Scalar>> cellZSpan_;
781 std::vector<std::pair<Scalar,Scalar>> cellZMinMax_;
782 int num_pressure_points_;
Definition InitStateEquil.hpp:691
Definition InitStateEquil.hpp:138
Definition InitStateEquil.hpp:113
Definition InitStateEquil.hpp:88
Calculator for phase saturations.
Definition InitStateEquil.hpp:386
const PhaseQuantityValue< Scalar > & deriveSaturations(const Position &x, const Region ®, const PTable &ptable)
Calculate phase saturations at particular point of the simulation model geometry.
Definition InitStateEquil_impl.hpp:596
PhaseSaturations(MaterialLawManager &matLawMgr, const std::vector< Scalar > &swatInit)
Constructor.
Definition InitStateEquil_impl.hpp:572
PressureTable< FluidSystem, Region > PTable
Convenience type alias.
Definition InitStateEquil.hpp:398
PhaseSaturations & operator=(const PhaseSaturations &)=delete
Disabled assignment operator.
PhaseSaturations & operator=(PhaseSaturations &&)=delete
Disabled move-assignment operator.
const PhaseQuantityValue< Scalar > & correctedPhasePressures() const
Retrieve saturation-corrected phase pressures.
Definition InitStateEquil.hpp:442
Definition InitStateEquil.hpp:167
PressureTable & operator=(const PressureTable &rhs)
Assignment operator.
Definition InitStateEquil_impl.hpp:1017
Scalar water(const Scalar depth) const
Evaluate water phase pressure at specified depth.
Definition InitStateEquil_impl.hpp:1097
Scalar gas(const Scalar depth) const
Evaluate gas phase pressure at specified depth.
Definition InitStateEquil_impl.hpp:1086
bool waterActive() const
Predicate for whether or not water is an active phase.
Definition InitStateEquil_impl.hpp:1068
bool gasActive() const
Predicate for whether or not gas is an active phase.
Definition InitStateEquil_impl.hpp:1061
Scalar oil(const Scalar depth) const
Evaluate oil phase pressure at specified depth.
Definition InitStateEquil_impl.hpp:1076
bool oilActive() const
Predicate for whether or not oil is an active phase.
Definition InitStateEquil_impl.hpp:1054
PressureTable(const Scalar gravity, const int samplePoints=2000)
Constructor.
Definition InitStateEquil_impl.hpp:987
Definition InitStateEquil.hpp:67
Aggregate information base of an equilibration region.
Definition InitStateEquil.hpp:61
Base class for phase mixing functions.
Definition InitStateEquil.hpp:62
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
The Opm property system, traits with inheritance.
Simple set of per-phase (named by primary component) quantities.
Definition InitStateEquil.hpp:338
Evaluation point within a model geometry.
Definition InitStateEquil.hpp:392
Definition InitStateEquil.hpp:255