27#ifndef OPM_BLACK_OIL_FLUID_SYSTEM_HPP
28#define OPM_BLACK_OIL_FLUID_SYSTEM_HPP
35#include <opm/common/TimingMacros.hpp>
67template <class FluidSystem, class FluidState, class LhsEval>
68LhsEval getRs_(typename std::enable_if<!HasMember_Rs<FluidState>::value, const FluidState&>::type fluidState,
72 decay<LhsEval>(fluidState.massFraction(FluidSystem::oilPhaseIdx, FluidSystem::gasCompIdx));
73 return FluidSystem::convertXoGToRs(XoG, regionIdx);
76template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
77auto getRs_(
typename std::enable_if<HasMember_Rs<FluidState>::value,
const FluidState&>::type fluidState,
79 ->
decltype(decay<LhsEval>(fluidState.Rs()))
80{
return decay<LhsEval>(fluidState.Rs()); }
82template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
83LhsEval getRv_(
typename std::enable_if<!HasMember_Rv<FluidState>::value,
const FluidState&>::type fluidState,
87 decay<LhsEval>(fluidState.massFraction(FluidSystem::gasPhaseIdx, FluidSystem::oilCompIdx));
88 return FluidSystem::convertXgOToRv(XgO, regionIdx);
91template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
92auto getRv_(
typename std::enable_if<HasMember_Rv<FluidState>::value,
const FluidState&>::type fluidState,
94 ->
decltype(decay<LhsEval>(fluidState.Rv()))
95{
return decay<LhsEval>(fluidState.Rv()); }
97template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
98LhsEval getRvw_(
typename std::enable_if<!HasMember_Rvw<FluidState>::value,
const FluidState&>::type fluidState,
102 decay<LhsEval>(fluidState.massFraction(FluidSystem::gasPhaseIdx, FluidSystem::waterCompIdx));
103 return FluidSystem::convertXgWToRvw(XgW, regionIdx);
106template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
107auto getRvw_(
typename std::enable_if<HasMember_Rvw<FluidState>::value,
const FluidState&>::type fluidState,
109 ->
decltype(decay<LhsEval>(fluidState.Rvw()))
110{
return decay<LhsEval>(fluidState.Rvw()); }
112template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
113LhsEval getRsw_(
typename std::enable_if<!HasMember_Rsw<FluidState>::value,
const FluidState&>::type fluidState,
117 decay<LhsEval>(fluidState.massFraction(FluidSystem::waterPhaseIdx, FluidSystem::gasCompIdx));
118 return FluidSystem::convertXwGToRsw(XwG, regionIdx);
121template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
122auto getRsw_(
typename std::enable_if<HasMember_Rsw<FluidState>::value,
const FluidState&>::type fluidState,
124 ->
decltype(decay<LhsEval>(fluidState.Rsw()))
125{
return decay<LhsEval>(fluidState.Rsw()); }
127template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
128LhsEval getSaltConcentration_(
typename std::enable_if<!HasMember_saltConcentration<FluidState>::value,
129 const FluidState&>::type,
133template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
134auto getSaltConcentration_(
typename std::enable_if<HasMember_saltConcentration<FluidState>::value,
const FluidState&>::type fluidState,
136 ->
decltype(decay<LhsEval>(fluidState.saltConcentration()))
137{
return decay<LhsEval>(fluidState.saltConcentration()); }
139template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
140LhsEval getSaltSaturation_(
typename std::enable_if<!HasMember_saltSaturation<FluidState>::value,
141 const FluidState&>::type,
145template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
146auto getSaltSaturation_(
typename std::enable_if<HasMember_saltSaturation<FluidState>::value,
const FluidState&>::type fluidState,
148 ->
decltype(decay<LhsEval>(fluidState.saltSaturation()))
149{
return decay<LhsEval>(fluidState.saltSaturation()); }
159template <
class Scalar,
class IndexTraits = BlackOilDefaultIndexTraits>
170 template <
class EvaluationT>
173 using Evaluation = EvaluationT;
177 : maxOilSat_(maxOilSat)
178 , regionIdx_(regionIdx)
189 template <
class OtherCache>
192 regionIdx_ = other.regionIndex();
193 maxOilSat_ = other.maxOilSat();
204 {
return regionIdx_; }
214 { regionIdx_ = val; }
216 const Evaluation& maxOilSat()
const
217 {
return maxOilSat_; }
219 void setMaxOilSat(
const Evaluation& val)
220 { maxOilSat_ = val; }
223 Evaluation maxOilSat_;
234 static void initFromState(
const EclipseState& eclState,
const Schedule& schedule);
245 static void initBegin(std::size_t numPvtRegions);
254 { enableDissolvedGas_ = yesno; }
263 { enableVaporizedOil_ = yesno; }
272 { enableVaporizedWater_ = yesno; }
281 { enableDissolvedGasInWater_ = yesno; }
288 { enableDiffusion_ = yesno; }
296 { useSaturatedTables_ = yesno; }
302 { gasPvt_ = pvtObj; }
308 { oilPvt_ = pvtObj; }
314 { waterPvt_ = pvtObj; }
316 static void setVapPars(
const Scalar par1,
const Scalar par2)
319 gasPvt_->setVapPars(par1, par2);
322 oilPvt_->setVapPars(par1, par2);
325 waterPvt_->setVapPars(par1, par2);
346 static bool isInitialized()
347 {
return isInitialized_; }
370 static std::string_view
phaseName(
unsigned phaseIdx);
387 static constexpr unsigned oilCompIdx = IndexTraits::oilCompIdx;
391 static constexpr unsigned gasCompIdx = IndexTraits::gasCompIdx;
394 static unsigned char numActivePhases_;
395 static std::array<bool,numPhases> phaseIsActive_;
400 {
return numActivePhases_; }
406 return phaseIsActive_[phaseIdx];
420 {
return molarMass_[regionIdx][compIdx]; }
448 {
return molarMass_.size(); }
457 {
return enableDissolvedGas_; }
467 {
return enableDissolvedGasInWater_; }
476 {
return enableVaporizedOil_; }
485 {
return enableVaporizedWater_; }
493 {
return enableDiffusion_; }
501 {
return useSaturatedTables_; }
509 {
return referenceDensity_[regionIdx][phaseIdx]; }
515 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
516 static LhsEval
density(
const FluidState& fluidState,
519 {
return density<FluidState, LhsEval>(fluidState, phaseIdx, paramCache.
regionIndex()); }
522 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
528 return fugacityCoefficient<FluidState, LhsEval>(fluidState,
535 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
539 {
return viscosity<FluidState, LhsEval>(fluidState, phaseIdx, paramCache.
regionIndex()); }
542 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
543 static LhsEval
enthalpy(
const FluidState& fluidState,
546 {
return enthalpy<FluidState, LhsEval>(fluidState, phaseIdx, paramCache.
regionIndex()); }
548 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
549 static LhsEval internalEnergy(
const FluidState& fluidState,
550 const ParameterCache<ParamCacheEval>& paramCache,
552 {
return internalEnergy<FluidState, LhsEval>(fluidState, phaseIdx, paramCache.regionIndex()); }
559 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
560 static LhsEval
density(
const FluidState& fluidState,
567 const LhsEval& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
568 const LhsEval& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
569 const LhsEval& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
575 const LhsEval& Rs = BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
576 const LhsEval& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
584 const LhsEval Rs(0.0);
585 const auto& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
593 const LhsEval& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
594 const LhsEval& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
595 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
604 const LhsEval Rvw(0.0);
605 const LhsEval& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
606 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
614 const LhsEval Rv(0.0);
615 const LhsEval& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
616 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
624 const LhsEval Rv(0.0);
625 const LhsEval Rvw(0.0);
626 const auto& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
633 const LhsEval& Rsw =BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
634 const LhsEval& bw = waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
639 const LhsEval Rsw(0.0);
642 * waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
645 throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
657 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
665 const auto& p = fluidState.pressure(phaseIdx);
666 const auto& T = fluidState.temperature(phaseIdx);
672 const LhsEval& Rs = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState,
oilPhaseIdx, regionIdx);
673 const LhsEval& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
681 const LhsEval Rs(0.0);
682 const LhsEval& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
689 const LhsEval& Rv = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState,
gasPhaseIdx, regionIdx);
690 const LhsEval& Rvw = saturatedVaporizationFactor<FluidState, LhsEval>(fluidState,
gasPhaseIdx, regionIdx);
691 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
701 const LhsEval Rvw(0.0);
702 const LhsEval& Rv = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState,
gasPhaseIdx, regionIdx);
703 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
712 const LhsEval Rv(0.0);
713 const LhsEval& Rvw = saturatedVaporizationFactor<FluidState, LhsEval>(fluidState,
gasPhaseIdx, regionIdx);
714 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
722 const LhsEval Rv(0.0);
723 const LhsEval Rvw(0.0);
724 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
734 const auto& saltConcentration = decay<LhsEval>(fluidState.saltConcentration());
735 const LhsEval& Rsw = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState,
waterPhaseIdx, regionIdx);
736 const LhsEval& bw = waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
743 *inverseFormationVolumeFactor<FluidState, LhsEval>(fluidState,
waterPhaseIdx, regionIdx);
747 throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
758 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
767 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
768 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
773 const auto& Rs = BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
775 && Rs >= (1.0 - 1e-10)*oilPvt_->saturatedGasDissolutionFactor(regionIdx, scalarValue(T), scalarValue(p)))
777 return oilPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
779 return oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
783 const LhsEval Rs(0.0);
784 return oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
788 const auto& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
789 const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
791 && Rvw >= (1.0 - 1e-10)*gasPvt_->saturatedWaterVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p))
793 && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
795 return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
797 return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
802 const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
804 && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
806 return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
808 const LhsEval Rvw(0.0);
809 return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
814 const auto& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
816 && Rvw >= (1.0 - 1e-10)*gasPvt_->saturatedWaterVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
818 return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
820 const LhsEval Rv(0.0);
821 return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
825 const LhsEval Rv(0.0);
826 const LhsEval Rvw(0.0);
827 return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
831 const auto& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
833 const auto& Rsw = BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
835 && Rsw >= (1.0 - 1e-10)*waterPvt_->saturatedGasDissolutionFactor(regionIdx, scalarValue(T), scalarValue(p), scalarValue(saltConcentration)))
837 return waterPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p, saltConcentration);
839 return waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
842 const LhsEval Rsw(0.0);
843 return waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
845 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
858 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
867 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
868 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
869 const auto& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
872 case oilPhaseIdx:
return oilPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
873 case gasPhaseIdx:
return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
874 case waterPhaseIdx:
return waterPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p, saltConcentration);
875 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
880 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
890 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
891 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
896 const LhsEval phi_oO = 20e3/p;
899 constexpr const Scalar phi_gG = 1.0;
903 const LhsEval phi_wW = 30e3/p;
921 const auto& R_vSat = gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p);
925 const auto& R_sSat = oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p);
928 const auto& x_oOSat = 1.0 - x_oGSat;
930 const auto& p_o = decay<LhsEval>(fluidState.pressure(
oilPhaseIdx));
931 const auto& p_g = decay<LhsEval>(fluidState.pressure(
gasPhaseIdx));
933 return phi_oO*p_o*x_oOSat / (p_g*x_gOSat);
941 throw std::logic_error(
"Invalid component index "+std::to_string(compIdx));
957 const auto& R_vSat = gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p);
960 const auto& x_gGSat = 1.0 - x_gOSat;
962 const auto& R_sSat = oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p);
966 const auto& p_o = decay<LhsEval>(fluidState.pressure(
oilPhaseIdx));
967 const auto& p_g = decay<LhsEval>(fluidState.pressure(
gasPhaseIdx));
969 return phi_gG*p_g*x_gGSat / (p_o*x_oGSat);
976 throw std::logic_error(
"Invalid component index "+std::to_string(compIdx));
991 throw std::logic_error(
"Invalid component index "+std::to_string(compIdx));
995 throw std::logic_error(
"Invalid phase index "+std::to_string(phaseIdx));
998 throw std::logic_error(
"Unhandled phase or component index");
1002 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1011 const LhsEval& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1012 const LhsEval& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1017 const auto& Rs = BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1019 && Rs >= (1.0 - 1e-10)*oilPvt_->saturatedGasDissolutionFactor(regionIdx, scalarValue(T), scalarValue(p)))
1021 return oilPvt_->saturatedViscosity(regionIdx, T, p);
1023 return oilPvt_->viscosity(regionIdx, T, p, Rs);
1027 const LhsEval Rs(0.0);
1028 return oilPvt_->viscosity(regionIdx, T, p, Rs);
1033 const auto& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1034 const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1036 && Rvw >= (1.0 - 1e-10)*gasPvt_->saturatedWaterVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p))
1038 && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
1040 return gasPvt_->saturatedViscosity(regionIdx, T, p);
1042 return gasPvt_->viscosity(regionIdx, T, p, Rv, Rvw);
1046 const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1048 && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
1050 return gasPvt_->saturatedViscosity(regionIdx, T, p);
1052 const LhsEval Rvw(0.0);
1053 return gasPvt_->viscosity(regionIdx, T, p, Rv, Rvw);
1057 const auto& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1059 && Rvw >= (1.0 - 1e-10)*gasPvt_->saturatedWaterVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
1061 return gasPvt_->saturatedViscosity(regionIdx, T, p);
1063 const LhsEval Rv(0.0);
1064 return gasPvt_->viscosity(regionIdx, T, p, Rv, Rvw);
1068 const LhsEval Rv(0.0);
1069 const LhsEval Rvw(0.0);
1070 return gasPvt_->viscosity(regionIdx, T, p, Rv, Rvw);
1075 const LhsEval& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1077 const auto& Rsw = BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1079 && Rsw >= (1.0 - 1e-10)*waterPvt_->saturatedGasDissolutionFactor(regionIdx, scalarValue(T), scalarValue(p), scalarValue(saltConcentration)))
1081 return waterPvt_->saturatedViscosity(regionIdx, T, p, saltConcentration);
1083 return waterPvt_->viscosity(regionIdx, T, p, Rsw, saltConcentration);
1086 const LhsEval Rsw(0.0);
1087 return waterPvt_->viscosity(regionIdx, T, p, Rsw, saltConcentration);
1091 throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1094 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1095 static LhsEval internalEnergy(
const FluidState& fluidState,
1097 unsigned regionIdx){
1098 bool is_mixing =
false;
1099 const LhsEval& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1100 const LhsEval& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1104 if(!oilPvt_->mixingEnergy()){
1105 const auto& oilEnergy =
1106 oilPvt_->internalEnergy(regionIdx, T, p,
1107 BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1115 if(!waterPvt_->mixingEnergy()){
1116 const auto waterEnergy =
1117 waterPvt_->internalEnergy(regionIdx, T, p,
1118 BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1119 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1127 if(!gasPvt_->mixingEnergy()){
1128 const auto gasEnergy =
1129 gasPvt_->internalEnergy(regionIdx, T, p,
1130 BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1131 BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1139 assert(is_mixing==
true);
1140 const auto& energy = internalMixingTotalEnergy<FluidState,LhsEval>(fluidState, phaseIdx, regionIdx);
1141 const auto& phase_density = density<FluidState,LhsEval>(fluidState, phaseIdx, regionIdx);
1142 return energy/phase_density;
1146 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1147 static LhsEval internalMixingTotalEnergy(
const FluidState& fluidState,
1153 const LhsEval& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1154 const LhsEval& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1155 const LhsEval& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1159 auto oilEnergy = oilPvt_->internalEnergy(regionIdx, T, p,
1160 BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1161 assert(oilPvt_->mixingEnergy());
1165 const LhsEval& Rs = BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1166 const LhsEval& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
1167 const auto& gasEnergy =
1168 gasPvt_->internalEnergy(regionIdx, T, p,
1169 BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1170 BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1171 const auto hVapG = gasPvt_->hVap(regionIdx);
1178 const LhsEval Rs(0.0);
1179 const auto& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
1185 const auto& gasEnergy =
1186 gasPvt_->internalEnergy(regionIdx, T, p,
1187 BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1188 BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1189 assert(gasPvt_->mixingEnergy());
1191 const auto& oilEnergy =
1192 oilPvt_->internalEnergy(regionIdx, T, p,
1193 BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1194 const auto waterEnergy =
1195 waterPvt_->internalEnergy(regionIdx, T, p,
1196 BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1197 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1199 const LhsEval& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1200 const LhsEval& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1201 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
1202 const auto hVapO = oilPvt_->hVap(regionIdx);
1203 const auto hVapW = waterPvt_->hVap(regionIdx);
1210 const auto& oilEnergy =
1211 oilPvt_->internalEnergy(regionIdx, T, p,
1212 BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1214 const LhsEval Rvw(0.0);
1215 const LhsEval& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1216 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
1217 const auto hVapO = oilPvt_->hVap(regionIdx);
1224 const LhsEval Rv(0.0);
1225 const LhsEval& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1226 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
1227 const auto waterEnergy =
1228 waterPvt_->internalEnergy(regionIdx, T, p,
1229 BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1230 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1231 const auto hVapW = waterPvt_->hVap(regionIdx);
1238 const LhsEval Rv(0.0);
1239 const LhsEval Rvw(0.0);
1240 const auto& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
1245 const auto waterEnergy =
1246 waterPvt_->internalEnergy(regionIdx, T, p,
1247 BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1248 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1249 assert(waterPvt_->mixingEnergy());
1251 const auto& gasEnergy =
1252 gasPvt_->internalEnergy(regionIdx, T, p,
1253 BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1254 BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1256 const LhsEval& Rsw = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState,
waterPhaseIdx, regionIdx);
1257 const LhsEval& bw = waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
1262 const LhsEval Rsw(0.0);
1265 * waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
1267 throw std::logic_error(
"Unhandled phase index " + std::to_string(phaseIdx));
1273 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1279 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1280 auto energy = internalEnergy<FluidState, LhsEval>(fluidState, phaseIdx, regionIdx);
1281 if(!enthalpy_eq_energy_){
1283 energy += p/density<FluidState, LhsEval>(fluidState, phaseIdx, regionIdx);
1294 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1302 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1303 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1304 const auto& saltConcentration = decay<LhsEval>(fluidState.saltConcentration());
1308 case gasPhaseIdx:
return gasPvt_->saturatedWaterVaporizationFactor(regionIdx, T, p, saltConcentration);
1310 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1320 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1324 const LhsEval& maxOilSaturation)
1330 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1331 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1335 case oilPhaseIdx:
return oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p, So, maxOilSaturation);
1336 case gasPhaseIdx:
return gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p, So, maxOilSaturation);
1337 case waterPhaseIdx:
return waterPvt_->saturatedGasDissolutionFactor(regionIdx, T, p,
1338 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1339 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1351 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1360 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1361 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1364 case oilPhaseIdx:
return oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p);
1365 case gasPhaseIdx:
return gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p);
1366 case waterPhaseIdx:
return waterPvt_->saturatedGasDissolutionFactor(regionIdx, T, p,
1367 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1368 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1375 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1386 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1403 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1411 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1414 case oilPhaseIdx:
return oilPvt_->saturationPressure(regionIdx, T, BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1415 case gasPhaseIdx:
return gasPvt_->saturationPressure(regionIdx, T, BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1416 case waterPhaseIdx:
return waterPvt_->saturationPressure(regionIdx, T,
1417 BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1418 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1419 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1430 template <
class LhsEval>
1436 return XoG/(1.0 - XoG)*(rho_oRef/rho_gRef);
1443 template <
class LhsEval>
1449 return XwG/(1.0 - XwG)*(rho_wRef/rho_gRef);
1456 template <
class LhsEval>
1462 return XgO/(1.0 - XgO)*(rho_gRef/rho_oRef);
1469 template <
class LhsEval>
1475 return XgW/(1.0 - XgW)*(rho_gRef/rho_wRef);
1483 template <
class LhsEval>
1489 const LhsEval& rho_oG = Rs*rho_gRef;
1490 return rho_oG/(rho_oRef + rho_oG);
1497 template <
class LhsEval>
1503 const LhsEval& rho_wG = Rsw*rho_gRef;
1504 return rho_wG/(rho_wRef + rho_wG);
1511 template <
class LhsEval>
1517 const LhsEval& rho_gO = Rv*rho_oRef;
1518 return rho_gO/(rho_gRef + rho_gO);
1525 template <
class LhsEval>
1531 const LhsEval& rho_gW = Rvw*rho_wRef;
1532 return rho_gW/(rho_gRef + rho_gW);
1538 template <
class LhsEval>
1544 return XgW*MG / (MW*(1 - XgW) + XgW*MG);
1550 template <
class LhsEval>
1556 return XwG*MW / (MG*(1 - XwG) + XwG*MW);
1562 template <
class LhsEval>
1568 return XoG*MO / (MG*(1 - XoG) + XoG*MO);
1574 template <
class LhsEval>
1580 return xoG*MG / (xoG*(MG - MO) + MO);
1586 template <
class LhsEval>
1592 return XgO*MG / (MO*(1 - XgO) + XgO*MG);
1598 template <
class LhsEval>
1604 return xgO*MO / (xgO*(MO - MG) + MG);
1615 {
return *gasPvt_; }
1625 {
return *oilPvt_; }
1635 {
return *waterPvt_; }
1643 {
return reservoirTemperature_; }
1651 { reservoirTemperature_ = value; }
1653 static short activeToCanonicalPhaseIdx(
unsigned activePhaseIdx);
1655 static short canonicalToActivePhaseIdx(
unsigned phaseIdx);
1659 {
return diffusionCoefficients_[regionIdx][
numPhases*compIdx + phaseIdx]; }
1663 { diffusionCoefficients_[regionIdx][
numPhases*compIdx + phaseIdx] = coefficient ; }
1668 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
1679 if(!diffusionCoefficients_.empty()) {
1683 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1684 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1690 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1693 static void setEnergyEqualEnthalpy(
bool enthalpy_eq_energy){
1694 enthalpy_eq_energy_ = enthalpy_eq_energy;
1697 static bool enthalpyEqualEnergy(){
1698 return enthalpy_eq_energy_;
1702 static void resizeArrays_(std::size_t
numRegions);
1704 static Scalar reservoirTemperature_;
1706 static std::shared_ptr<GasPvt> gasPvt_;
1707 static std::shared_ptr<OilPvt> oilPvt_;
1708 static std::shared_ptr<WaterPvt> waterPvt_;
1710 static bool enableDissolvedGas_;
1711 static bool enableDissolvedGasInWater_;
1712 static bool enableVaporizedOil_;
1713 static bool enableVaporizedWater_;
1714 static bool enableDiffusion_;
1719 static std::vector<std::array<
Scalar, 3> > referenceDensity_;
1720 static std::vector<std::array<
Scalar, 3> > molarMass_;
1721 static std::vector<std::array<
Scalar, 3 * 3> > diffusionCoefficients_;
1723 static std::array<short, numPhases> activeToCanonicalPhaseIdx_;
1724 static std::array<short, numPhases> canonicalToActivePhaseIdx_;
1726 static bool isInitialized_;
1727 static bool useSaturatedTables_;
1728 inline static bool enthalpy_eq_energy_ =
false;
1731template <
typename T>
using BOFS = BlackOilFluidSystem<T, BlackOilDefaultIndexTraits>;
1733#define DECLARE_INSTANCE(T) \
1734template<> unsigned char BOFS<T>::numActivePhases_; \
1735template<> std::array<bool, BOFS<T>::numPhases> BOFS<T>::phaseIsActive_; \
1736template<> std::array<short, BOFS<T>::numPhases> BOFS<T>::activeToCanonicalPhaseIdx_; \
1737template<> std::array<short, BOFS<T>::numPhases> BOFS<T>::canonicalToActivePhaseIdx_; \
1738template<> T BOFS<T>::surfaceTemperature; \
1739template<> T BOFS<T>::surfacePressure; \
1740template<> T BOFS<T>::reservoirTemperature_; \
1741template<> bool BOFS<T>::enableDissolvedGas_; \
1742template<> bool BOFS<T>::enableDissolvedGasInWater_; \
1743template<> bool BOFS<T>::enableVaporizedOil_; \
1744template<> bool BOFS<T>::enableVaporizedWater_; \
1745template<> bool BOFS<T>::enableDiffusion_; \
1746template<> std::shared_ptr<OilPvtMultiplexer<T>> BOFS<T>::oilPvt_; \
1747template<> std::shared_ptr<GasPvtMultiplexer<T>> BOFS<T>::gasPvt_; \
1748template<> std::shared_ptr<WaterPvtMultiplexer<T>> BOFS<T>::waterPvt_; \
1749template<> std::vector<std::array<T, 3>> BOFS<T>::referenceDensity_; \
1750template<> std::vector<std::array<T, 3>> BOFS<T>::molarMass_; \
1751template<> std::vector<std::array<T, 9>> BOFS<T>::diffusionCoefficients_; \
1752template<> bool BOFS<T>::isInitialized_; \
1753template<> bool BOFS<T>::useSaturatedTables_;
1755DECLARE_INSTANCE(
float)
1756DECLARE_INSTANCE(
double)
1758#undef DECLARE_INSTANCE
The base class for all fluid systems.
A central place for various physical constants occuring in some equations.
This class represents the Pressure-Volume-Temperature relations of the gas phase in the black-oil mod...
This macro generates a class HasMember_${MEMBER_NAME} which can be used for template specialization.
#define OPM_GENERATE_HAS_MEMBER(MEMBER_NAME,...)
This macro generates a class HasMember_${MEMBER_NAME} which can be used for template specialization.
Definition HasMemberGeneratorMacros.hpp:49
A parameter cache which does nothing.
This class represents the Pressure-Volume-Temperature relations of the oil phase in the black-oil mod...
Some templates to wrap the valgrind client request macros.
This class represents the Pressure-Volume-Temperature relations of the water phase in the black-oil m...
The base class for all fluid systems.
Definition BaseFluidSystem.hpp:43
ScalarT Scalar
The type used for scalar quantities.
Definition BaseFluidSystem.hpp:48
A fluid system which uses the black-oil model assumptions to calculate termodynamically meaningful qu...
Definition BlackOilFluidSystem.hpp:161
static constexpr unsigned waterCompIdx
Index of the water component.
Definition BlackOilFluidSystem.hpp:389
static unsigned numActivePhases()
Returns the number of active fluid phases (i.e., usually three)
Definition BlackOilFluidSystem.hpp:399
static Scalar diffusionCoefficient(unsigned compIdx, unsigned phaseIdx, unsigned regionIdx=0)
Calculate the binary molecular diffusion coefficient for a component in a fluid phase [mol^2 * s / (k...
Definition BlackOilFluidSystem.hpp:1658
static LhsEval convertXoGToRs(const LhsEval &XoG, unsigned regionIdx)
Convert the mass fraction of the gas component in the oil phase to the corresponding gas dissolution ...
Definition BlackOilFluidSystem.hpp:1431
static LhsEval convertRsToXoG(const LhsEval &Rs, unsigned regionIdx)
Convert a gas dissolution factor to the the corresponding mass fraction of the gas component in the o...
Definition BlackOilFluidSystem.hpp:1484
static constexpr unsigned numComponents
Number of chemical species in the fluid system.
Definition BlackOilFluidSystem.hpp:384
static LhsEval convertRvwToXgW(const LhsEval &Rvw, unsigned regionIdx)
Convert an water vaporization factor to the corresponding mass fraction of the water component in the...
Definition BlackOilFluidSystem.hpp:1526
static LhsEval convertXwGToxwG(const LhsEval &XwG, unsigned regionIdx)
Convert a gas mass fraction in the water phase the corresponding mole fraction.
Definition BlackOilFluidSystem.hpp:1551
static void setUseSaturatedTables(bool yesno)
Specify whether the saturated tables should be used.
Definition BlackOilFluidSystem.hpp:295
static LhsEval saturatedInverseFormationVolumeFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the formation volume factor of a "saturated" fluid phase.
Definition BlackOilFluidSystem.hpp:859
static LhsEval convertXoGToxoG(const LhsEval &XoG, unsigned regionIdx)
Convert a gas mass fraction in the oil phase the corresponding mole fraction.
Definition BlackOilFluidSystem.hpp:1563
static LhsEval viscosity(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Calculate the dynamic viscosity of a fluid phase [Pa*s].
Definition BlackOilFluidSystem.hpp:1003
static void initEnd()
Finish initializing the black oil fluid system.
Definition BlackOilFluidSystem.cpp:267
static Scalar reservoirTemperature(unsigned=0)
Set the temperature of the reservoir.
Definition BlackOilFluidSystem.hpp:1642
static LhsEval convertRswToXwG(const LhsEval &Rsw, unsigned regionIdx)
Convert a gas dissolution factor to the the corresponding mass fraction of the gas component in the w...
Definition BlackOilFluidSystem.hpp:1498
static bool phaseIsActive(unsigned phaseIdx)
Returns whether a fluid phase is active.
Definition BlackOilFluidSystem.hpp:403
static void setEnableDiffusion(bool yesno)
Specify whether the fluid system should consider diffusion.
Definition BlackOilFluidSystem.hpp:287
static LhsEval bubblePointPressure(const FluidState &fluidState, unsigned regionIdx)
Returns the bubble point pressure $P_b$ using the current Rs.
Definition BlackOilFluidSystem.hpp:1376
static LhsEval density(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Calculate the density [kg/m^3] of a fluid phase.
Definition BlackOilFluidSystem.hpp:560
static LhsEval saturatedDensity(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Compute the density of a saturated fluid phase.
Definition BlackOilFluidSystem.hpp:658
static LhsEval viscosity(const FluidState &fluidState, const ParameterCache< ParamCacheEval > ¶mCache, unsigned phaseIdx)
Calculate the dynamic viscosity of a fluid phase [Pa*s].
Definition BlackOilFluidSystem.hpp:536
static bool isIdealMixture(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal mixture.
Definition BlackOilFluidSystem.hpp:423
static bool enableVaporizedWater()
Returns whether the fluid system should consider that the water component can dissolve in the gas pha...
Definition BlackOilFluidSystem.hpp:484
static constexpr unsigned numPhases
Number of fluid phases in the fluid system.
Definition BlackOilFluidSystem.hpp:354
static void setWaterPvt(std::shared_ptr< WaterPvt > pvtObj)
Set the pressure-volume-saturation (PVT) relations for the water phase.
Definition BlackOilFluidSystem.hpp:313
static void initBegin(std::size_t numPvtRegions)
Begin the initialization of the black oil fluid system.
Definition BlackOilFluidSystem.cpp:229
static void setReservoirTemperature(Scalar value)
Return the temperature of the reservoir.
Definition BlackOilFluidSystem.hpp:1650
static void setEnableDissolvedGasInWater(bool yesno)
Specify whether the fluid system should consider that the gas component can dissolve in the water pha...
Definition BlackOilFluidSystem.hpp:280
static LhsEval saturatedDissolutionFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx, const LhsEval &maxOilSaturation)
Returns the dissolution factor of a saturated fluid phase.
Definition BlackOilFluidSystem.hpp:1321
static LhsEval enthalpy(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Given a phase's composition, temperature, pressure and density, calculate its specific enthalpy [J/kg...
Definition BlackOilFluidSystem.hpp:1274
static constexpr unsigned oilCompIdx
Index of the oil component.
Definition BlackOilFluidSystem.hpp:387
static bool isLiquid(unsigned phaseIdx)
Return whether a phase is liquid.
Definition BlackOilFluidSystem.hpp:373
static LhsEval enthalpy(const FluidState &fluidState, const ParameterCache< ParamCacheEval > ¶mCache, unsigned phaseIdx)
Given a phase's composition, temperature, pressure and density, calculate its specific enthalpy [J/kg...
Definition BlackOilFluidSystem.hpp:543
static LhsEval fugacityCoefficient(const FluidState &fluidState, unsigned phaseIdx, unsigned compIdx, unsigned regionIdx)
Calculate the fugacity coefficient [Pa] of an individual component in a fluid phase.
Definition BlackOilFluidSystem.hpp:881
static LhsEval convertxoGToXoG(const LhsEval &xoG, unsigned regionIdx)
Convert a gas mole fraction in the oil phase the corresponding mass fraction.
Definition BlackOilFluidSystem.hpp:1575
static Scalar molarMass(unsigned compIdx, unsigned regionIdx=0)
Return the molar mass of a component in [kg/mol].
Definition BlackOilFluidSystem.hpp:419
static LhsEval saturatedDissolutionFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the dissolution factor of a saturated fluid phase.
Definition BlackOilFluidSystem.hpp:1352
static LhsEval dewPointPressure(const FluidState &fluidState, unsigned regionIdx)
Returns the dew point pressure $P_d$ using the current Rv.
Definition BlackOilFluidSystem.hpp:1387
static LhsEval saturationPressure(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the saturation pressure of a given phase [Pa] depending on its composition.
Definition BlackOilFluidSystem.hpp:1404
static LhsEval convertXwGToRsw(const LhsEval &XwG, unsigned regionIdx)
Convert the mass fraction of the gas component in the water phase to the corresponding gas dissolutio...
Definition BlackOilFluidSystem.hpp:1444
static std::string_view componentName(unsigned compIdx)
Return the human readable name of a component.
Definition BlackOilFluidSystem.cpp:362
static LhsEval convertXgOToxgO(const LhsEval &XgO, unsigned regionIdx)
Convert a oil mass fraction in the gas phase the corresponding mole fraction.
Definition BlackOilFluidSystem.hpp:1587
static LhsEval diffusionCoefficient(const FluidState &fluidState, const ParameterCache< ParamCacheEval > ¶mCache, unsigned phaseIdx, unsigned compIdx)
Calculate the binary molecular diffusion coefficient for a component in a fluid phase [mol^2 * s / (k...
Definition BlackOilFluidSystem.hpp:1669
static LhsEval saturatedVaporizationFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the water vaporization factor of saturated phase.
Definition BlackOilFluidSystem.hpp:1295
static bool enableDissolvedGas()
Returns whether the fluid system should consider that the gas component can dissolve in the oil phase...
Definition BlackOilFluidSystem.hpp:456
static Scalar referenceDensity(unsigned phaseIdx, unsigned regionIdx)
Returns the density of a fluid phase at surface pressure [kg/m^3].
Definition BlackOilFluidSystem.hpp:508
static Scalar surfaceTemperature
The temperature at the surface.
Definition BlackOilFluidSystem.hpp:367
static void setOilPvt(std::shared_ptr< OilPvt > pvtObj)
Set the pressure-volume-saturation (PVT) relations for the oil phase.
Definition BlackOilFluidSystem.hpp:307
static constexpr unsigned waterPhaseIdx
Index of the water phase.
Definition BlackOilFluidSystem.hpp:357
static bool useSaturatedTables()
Returns whether the saturated tables should be used.
Definition BlackOilFluidSystem.hpp:500
static LhsEval convertxgOToXgO(const LhsEval &xgO, unsigned regionIdx)
Convert a oil mole fraction in the gas phase the corresponding mass fraction.
Definition BlackOilFluidSystem.hpp:1599
static LhsEval inverseFormationVolumeFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the formation volume factor of an "undersaturated" fluid phase.
Definition BlackOilFluidSystem.hpp:759
static LhsEval convertXgWToxgW(const LhsEval &XgW, unsigned regionIdx)
Convert a water mass fraction in the gas phase the corresponding mole fraction.
Definition BlackOilFluidSystem.hpp:1539
static constexpr unsigned gasPhaseIdx
Index of the gas phase.
Definition BlackOilFluidSystem.hpp:361
static void setEnableDissolvedGas(bool yesno)
Specify whether the fluid system should consider that the gas component can dissolve in the oil phase...
Definition BlackOilFluidSystem.hpp:253
static LhsEval convertXgOToRv(const LhsEval &XgO, unsigned regionIdx)
Convert the mass fraction of the oil component in the gas phase to the corresponding oil vaporization...
Definition BlackOilFluidSystem.hpp:1457
static LhsEval fugacityCoefficient(const FluidState &fluidState, const ParameterCache< ParamCacheEval > ¶mCache, unsigned phaseIdx, unsigned compIdx)
Calculate the fugacity coefficient [Pa] of an individual component in a fluid phase.
Definition BlackOilFluidSystem.hpp:523
static LhsEval density(const FluidState &fluidState, const ParameterCache< ParamCacheEval > ¶mCache, unsigned phaseIdx)
Calculate the density [kg/m^3] of a fluid phase.
Definition BlackOilFluidSystem.hpp:516
static constexpr unsigned gasCompIdx
Index of the gas component.
Definition BlackOilFluidSystem.hpp:391
static LhsEval convertRvToXgO(const LhsEval &Rv, unsigned regionIdx)
Convert an oil vaporization factor to the corresponding mass fraction of the oil component in the gas...
Definition BlackOilFluidSystem.hpp:1512
static const WaterPvt & waterPvt()
Return a reference to the low-level object which calculates the water phase quantities.
Definition BlackOilFluidSystem.hpp:1634
static bool enableVaporizedOil()
Returns whether the fluid system should consider that the oil component can dissolve in the gas phase...
Definition BlackOilFluidSystem.hpp:475
static bool isIdealGas(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal gas.
Definition BlackOilFluidSystem.hpp:435
static unsigned solventComponentIndex(unsigned phaseIdx)
returns the index of "primary" component of a phase (solvent)
Definition BlackOilFluidSystem.cpp:323
static Scalar surfacePressure
The pressure at the surface.
Definition BlackOilFluidSystem.hpp:364
static void setEnableVaporizedWater(bool yesno)
Specify whether the fluid system should consider that the water component can dissolve in the gas pha...
Definition BlackOilFluidSystem.hpp:271
static constexpr unsigned oilPhaseIdx
Index of the oil phase.
Definition BlackOilFluidSystem.hpp:359
static const OilPvt & oilPvt()
Return a reference to the low-level object which calculates the oil phase quantities.
Definition BlackOilFluidSystem.hpp:1624
static bool enableDissolvedGasInWater()
Returns whether the fluid system should consider that the gas component can dissolve in the water pha...
Definition BlackOilFluidSystem.hpp:466
static unsigned soluteComponentIndex(unsigned phaseIdx)
returns the index of "secondary" component of a phase (solute)
Definition BlackOilFluidSystem.cpp:340
static bool enableDiffusion()
Returns whether the fluid system should consider diffusion.
Definition BlackOilFluidSystem.hpp:492
static LhsEval convertXgWToRvw(const LhsEval &XgW, unsigned regionIdx)
Convert the mass fraction of the water component in the gas phase to the corresponding water vaporiza...
Definition BlackOilFluidSystem.hpp:1470
static void setReferenceDensities(Scalar rhoOil, Scalar rhoWater, Scalar rhoGas, unsigned regionIdx)
Initialize the values of the reference densities.
Definition BlackOilFluidSystem.cpp:256
static std::size_t numRegions()
Returns the number of PVT regions which are considered.
Definition BlackOilFluidSystem.hpp:447
static void setEnableVaporizedOil(bool yesno)
Specify whether the fluid system should consider that the oil component can dissolve in the gas phase...
Definition BlackOilFluidSystem.hpp:262
static std::string_view phaseName(unsigned phaseIdx)
Return the human readable name of a fluid phase.
Definition BlackOilFluidSystem.cpp:306
static void setDiffusionCoefficient(Scalar coefficient, unsigned compIdx, unsigned phaseIdx, unsigned regionIdx=0)
Definition BlackOilFluidSystem.hpp:1662
static const GasPvt & gasPvt()
Return a reference to the low-level object which calculates the gas phase quantities.
Definition BlackOilFluidSystem.hpp:1614
static bool isCompressible(unsigned)
Returns true if and only if a fluid phase is assumed to be compressible.
Definition BlackOilFluidSystem.hpp:431
static void setGasPvt(std::shared_ptr< GasPvt > pvtObj)
Set the pressure-volume-saturation (PVT) relations for the gas phase.
Definition BlackOilFluidSystem.hpp:301
This class represents the Pressure-Volume-Temperature relations of the gas phase in the black-oil mod...
Definition GasPvtMultiplexer.hpp:110
Evaluation diffusionCoefficient(const Evaluation &temperature, const Evaluation &pressure, unsigned compIdx) const
Calculate the binary molecular diffusion coefficient for a component in a fluid phase [mol^2 * s / (k...
Definition GasPvtMultiplexer.hpp:268
A parameter cache which does nothing.
Definition NullParameterCache.hpp:40
This class represents the Pressure-Volume-Temperature relations of the oil phase in the black-oil mod...
Definition OilPvtMultiplexer.hpp:105
Evaluation diffusionCoefficient(const Evaluation &temperature, const Evaluation &pressure, unsigned compIdx) const
Calculate the binary molecular diffusion coefficient for a component in a fluid phase [mol^2 * s / (k...
Definition OilPvtMultiplexer.hpp:240
This class represents the Pressure-Volume-Temperature relations of the water phase in the black-oil m...
Definition WaterPvtMultiplexer.hpp:90
Evaluation diffusionCoefficient(const Evaluation &temperature, const Evaluation &pressure, unsigned compIdx) const
Calculate the binary molecular diffusion coefficient for a component in a fluid phase [mol^2 * s / (k...
Definition WaterPvtMultiplexer.hpp:231
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30
The type of the fluid system's parameter cache.
Definition BlackOilFluidSystem.hpp:172
void assignPersistentData(const OtherCache &other)
Copy the data which is not dependent on the type of the Scalars from another parameter cache.
Definition BlackOilFluidSystem.hpp:190
void setRegionIndex(unsigned val)
Set the index of the region which should be used to determine the thermodynamic properties.
Definition BlackOilFluidSystem.hpp:213
unsigned regionIndex() const
Return the index of the region which should be used to determine the thermodynamic properties.
Definition BlackOilFluidSystem.hpp:203