29#ifndef OPM_EQUILIBRATION_HELPERS_HPP
30#define OPM_EQUILIBRATION_HELPERS_HPP
32#include <opm/material/common/Tabulated1DFunction.hpp>
34#include <opm/input/eclipse/EclipseState/InitConfig/Equil.hpp>
94namespace Miscibility {
123 const Scalar sat = 0.0)
const = 0;
130template<
class Scalar>
156 const Scalar = 0.0)
const override
168template <
class Flu
idSystem>
172 using Scalar =
typename FluidSystem::Scalar;
180 RsVD(
const int pvtRegionIdx,
181 const std::vector<Scalar>& depth,
182 const std::vector<Scalar>& rs);
184 virtual ~RsVD() =
default;
204 const Scalar satGas = 0.0)
const override;
207 using RsVsDepthFunc = Tabulated1DFunction<Scalar>;
209 const int pvtRegionIdx_;
210 RsVsDepthFunc rsVsDepth_;
212 Scalar satRs(
const Scalar press,
const Scalar temp)
const;
221template <
class Flu
idSystem>
225 using Scalar =
typename FluidSystem::Scalar;
233 PBVD(
const int pvtRegionIdx,
234 const std::vector<Scalar>& depth,
235 const std::vector<Scalar>& pbub);
237 virtual ~PBVD() =
default;
254 const Scalar cellPress,
256 const Scalar satGas = 0.0)
const override;
259 using PbubVsDepthFunc = Tabulated1DFunction<Scalar>;
261 const int pvtRegionIdx_;
262 PbubVsDepthFunc pbubVsDepth_;
264 Scalar satRs(
const Scalar press,
const Scalar temp)
const;
273template <
class Flu
idSystem>
276 using Scalar =
typename FluidSystem::Scalar;
285 PDVD(
const int pvtRegionIdx,
286 const std::vector<Scalar>& depth,
287 const std::vector<Scalar>& pdew);
289 virtual ~PDVD() =
default;
306 const Scalar cellPress,
308 const Scalar satOil = 0.0)
const override;
311 using PdewVsDepthFunc = Tabulated1DFunction<Scalar>;
313 const int pvtRegionIdx_;
314 PdewVsDepthFunc pdewVsDepth_;
316 Scalar satRv(
const Scalar press,
const Scalar temp)
const;
325template <
class Flu
idSystem>
328 using Scalar =
typename FluidSystem::Scalar;
337 RvVD(
const int pvtRegionIdx,
338 const std::vector<Scalar>& depth,
339 const std::vector<Scalar>& rv);
359 const Scalar satOil = 0.0)
const override;
362 using RvVsDepthFunc = Tabulated1DFunction<Scalar>;
364 const int pvtRegionIdx_;
365 RvVsDepthFunc rvVsDepth_;
367 Scalar satRv(
const Scalar press,
const Scalar temp)
const;
376template <
class Flu
idSystem>
379 using Scalar =
typename FluidSystem::Scalar;
388 RvwVD(
const int pvtRegionIdx,
389 const std::vector<Scalar>& depth,
390 const std::vector<Scalar>& rvw);
410 const Scalar satWat = 0.0)
const override;
413 using RvwVsDepthFunc = Tabulated1DFunction<Scalar>;
415 const int pvtRegionIdx_;
416 RvwVsDepthFunc rvwVsDepth_;
418 Scalar satRvw(
const Scalar press,
const Scalar temp)
const;
436template <
class Flu
idSystem>
439 using Scalar =
typename FluidSystem::Scalar;
449 const Scalar pContact,
450 const Scalar T_contact);
470 const Scalar satGas = 0.0)
const override;
473 const int pvtRegionIdx_;
474 Scalar rsSatContact_;
476 Scalar satRs(
const Scalar press,
const Scalar temp)
const;
494template <
class Flu
idSystem>
497 using Scalar =
typename FluidSystem::Scalar;
507 const Scalar pContact,
508 const Scalar T_contact);
528 const Scalar satOil = 0.0)
const override;
531 const int pvtRegionIdx_;
532 Scalar rvSatContact_;
534 Scalar satRv(
const Scalar press,
const Scalar temp)
const;
551template <
class Flu
idSystem>
554 using Scalar =
typename FluidSystem::Scalar;
564 const Scalar pContact,
565 const Scalar T_contact);
585 const Scalar satWat = 0.0)
const override;
588 const int pvtRegionIdx_;
589 Scalar rvwSatContact_;
591 Scalar satRvw(
const Scalar press,
const Scalar temp)
const;
615template<
class Scalar>
618 using TabulatedFunction = Tabulated1DFunction<Scalar>;
634 const TabulatedFunction& tempVdTable,
635 const TabulatedFunction& saltVdTable,
657 Scalar
datum()
const;
715 const TabulatedFunction& saltVdTable()
const;
716 const TabulatedFunction& tempVdTable()
const;
725 std::shared_ptr<Miscibility::RsFunction<Scalar>> rs_;
726 std::shared_ptr<Miscibility::RsFunction<Scalar>> rv_;
727 std::shared_ptr<Miscibility::RsFunction<Scalar>> rvw_;
728 const TabulatedFunction& tempVdTable_;
729 const TabulatedFunction& saltVdTable_;
738template <
class Flu
idSystem,
class MaterialLawManager>
741 using Scalar =
typename FluidSystem::Scalar;
742 PcEq(
const MaterialLawManager& materialLawManager,
745 const Scalar targetPc);
747 Scalar operator()(Scalar s)
const;
750 const MaterialLawManager& materialLawManager_;
753 const Scalar targetPc_;
756template <
class Flu
idSystem,
class MaterialLawManager>
757typename FluidSystem::Scalar
758minSaturations(
const MaterialLawManager& materialLawManager,
759 const int phase,
const int cell);
761template <
class Flu
idSystem,
class MaterialLawManager>
762typename FluidSystem::Scalar
763maxSaturations(
const MaterialLawManager& materialLawManager,
764 const int phase,
const int cell);
768template <
class Flu
idSystem,
class MaterialLawManager>
769typename FluidSystem::Scalar
770satFromPc(
const MaterialLawManager& materialLawManager,
773 const typename FluidSystem::Scalar targetPc,
774 const bool increasing =
false);
779template <
class Flu
idSystem,
class MaterialLawManager>
782 using Scalar =
typename FluidSystem::Scalar;
783 PcEqSum(
const MaterialLawManager& materialLawManager,
787 const Scalar targetPc);
789 Scalar operator()(Scalar s)
const;
792 const MaterialLawManager& materialLawManager_;
796 const Scalar targetPc_;
802template <
class Flu
idSystem,
class MaterialLawManager>
803typename FluidSystem::Scalar
808 const typename FluidSystem::Scalar targetPc);
811template <
class Flu
idSystem,
class MaterialLawManager>
812typename FluidSystem::Scalar
813satFromDepth(
const MaterialLawManager& materialLawManager,
814 const typename FluidSystem::Scalar cellDepth,
815 const typename FluidSystem::Scalar contactDepth,
818 const bool increasing =
false);
821template <
class Flu
idSystem,
class MaterialLawManager>
822bool isConstPc(
const MaterialLawManager& materialLawManager,
Aggregate information base of an equilibration region.
Definition InitStateEquil.hpp:61
Scalar datum() const
Datum depth in current region.
Definition EquilibrationHelpers_impl.hpp:373
int pvtIdx() const
Retrieve pvtIdx of the region.
Definition EquilibrationHelpers_impl.hpp:450
EquilReg(const EquilRecord &rec, std::shared_ptr< Miscibility::RsFunction< Scalar > > rs, std::shared_ptr< Miscibility::RsFunction< Scalar > > rv, std::shared_ptr< Miscibility::RsFunction< Scalar > > rvw, const TabulatedFunction &tempVdTable, const TabulatedFunction &saltVdTable, const int pvtIdx)
Constructor.
Definition EquilibrationHelpers_impl.hpp:355
Scalar pcgoGoc() const
Gas-oil capillary pressure at gas-oil contact.
Definition EquilibrationHelpers_impl.hpp:403
int equilibrationAccuracy() const
Accuracy/strategy for initial fluid-in-place calculation.
Definition EquilibrationHelpers_impl.hpp:409
const CalcEvaporation & evaporationCalculator() const
Retrieve vapourised oil-gas ratio calculator of current region.
Definition EquilibrationHelpers_impl.hpp:423
Scalar zgoc() const
Depth of gas-oil contact.
Definition EquilibrationHelpers_impl.hpp:397
Scalar pressure() const
Pressure at datum depth in current region.
Definition EquilibrationHelpers_impl.hpp:379
Scalar pcowWoc() const
water-oil capillary pressure at water-oil contact.
Definition EquilibrationHelpers_impl.hpp:391
Scalar zwoc() const
Depth of water-oil contact.
Definition EquilibrationHelpers_impl.hpp:385
const CalcDissolution & dissolutionCalculator() const
Retrieve dissolved gas-oil ratio calculator of current region.
Definition EquilibrationHelpers_impl.hpp:416
const CalcWaterEvaporation & waterEvaporationCalculator() const
Retrieve vapourised water-gas ratio calculator of current region.
Definition EquilibrationHelpers_impl.hpp:430
Type that implements "no phase mixing" policy.
Definition EquilibrationHelpers.hpp:132
Scalar operator()(const Scalar, const Scalar, const Scalar, const Scalar=0.0) const override
Function call.
Definition EquilibrationHelpers.hpp:153
Type that implements "dissolved gas-oil ratio" tabulated as a function of depth policy.
Definition EquilibrationHelpers.hpp:223
Scalar operator()(const Scalar depth, const Scalar cellPress, const Scalar temp, const Scalar satGas=0.0) const override
Function call.
Definition EquilibrationHelpers_impl.hpp:107
PBVD(const int pvtRegionIdx, const std::vector< Scalar > &depth, const std::vector< Scalar > &pbub)
Constructor.
Definition EquilibrationHelpers_impl.hpp:96
Type that implements "vaporized oil-gas ratio" tabulated as a function of depth policy.
Definition EquilibrationHelpers.hpp:275
Scalar operator()(const Scalar depth, const Scalar cellPress, const Scalar temp, const Scalar satOil=0.0) const override
Function call.
Definition EquilibrationHelpers_impl.hpp:144
PDVD(const int pvtRegionIdx, const std::vector< Scalar > &depth, const std::vector< Scalar > &pdew)
Constructor.
Definition EquilibrationHelpers_impl.hpp:133
Base class for phase mixing functions.
Definition InitStateEquil.hpp:62
virtual Scalar operator()(const Scalar depth, const Scalar press, const Scalar temp, const Scalar sat=0.0) const =0
Function call operator.
Type that implements "dissolved gas-oil ratio" tabulated as a function of depth policy.
Definition EquilibrationHelpers.hpp:170
Scalar operator()(const Scalar depth, const Scalar press, const Scalar temp, const Scalar satGas=0.0) const override
Function call.
Definition EquilibrationHelpers_impl.hpp:69
RsVD(const int pvtRegionIdx, const std::vector< Scalar > &depth, const std::vector< Scalar > &rs)
Constructor.
Definition EquilibrationHelpers_impl.hpp:58
Type that implements "vaporized oil-gas ratio" tabulated as a function of depth policy.
Definition EquilibrationHelpers.hpp:327
Scalar operator()(const Scalar depth, const Scalar press, const Scalar temp, const Scalar satOil=0.0) const override
Function call.
Definition EquilibrationHelpers_impl.hpp:181
RvVD(const int pvtRegionIdx, const std::vector< Scalar > &depth, const std::vector< Scalar > &rv)
Constructor.
Definition EquilibrationHelpers_impl.hpp:170
Type that implements "vaporized water-gas ratio" tabulated as a function of depth policy.
Definition EquilibrationHelpers.hpp:378
Scalar operator()(const Scalar depth, const Scalar press, const Scalar temp, const Scalar satWat=0.0) const override
Function call.
Definition EquilibrationHelpers_impl.hpp:224
RvwVD(const int pvtRegionIdx, const std::vector< Scalar > &depth, const std::vector< Scalar > &rvw)
Constructor.
Definition EquilibrationHelpers_impl.hpp:213
FluidSystem::Scalar satFromDepth(const MaterialLawManager &materialLawManager, const typename FluidSystem::Scalar cellDepth, const typename FluidSystem::Scalar contactDepth, const int phase, const int cell, const bool increasing=false)
Compute saturation from depth. Used for constant capillary pressure function.
Definition EquilibrationHelpers_impl.hpp:650
FluidSystem::Scalar satFromPc(const MaterialLawManager &materialLawManager, const int phase, const int cell, const typename FluidSystem::Scalar targetPc, const bool increasing=false)
Compute saturation of some phase corresponding to a given capillary pressure.
Definition EquilibrationHelpers_impl.hpp:586
bool isConstPc(const MaterialLawManager &materialLawManager, const int phase, const int cell)
Return true if capillary pressure function is constant.
Definition EquilibrationHelpers_impl.hpp:670
FluidSystem::Scalar satFromSumOfPcs(const MaterialLawManager &materialLawManager, const int phase1, const int phase2, const int cell, const typename FluidSystem::Scalar targetPc)
Compute saturation of some phase corresponding to a given capillary pressure, where the capillary pre...
Definition EquilibrationHelpers_impl.hpp:619
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
Functor for inverting a sum of capillary pressure functions.
Definition EquilibrationHelpers.hpp:781
Functor for inverting capillary pressure function.
Definition EquilibrationHelpers.hpp:740