27#ifndef OPM_BLACK_OIL_FLUID_SYSTEM_HPP
28#define OPM_BLACK_OIL_FLUID_SYSTEM_HPP
37#include <opm/common/TimingMacros.hpp>
68template <class FluidSystem, class FluidState, class LhsEval>
69LhsEval getRs_(typename std::enable_if<!HasMember_Rs<FluidState>::value, const FluidState&>::type fluidState,
73 decay<LhsEval>(fluidState.massFraction(FluidSystem::oilPhaseIdx, FluidSystem::gasCompIdx));
74 return FluidSystem::convertXoGToRs(XoG, regionIdx);
77template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
78auto getRs_(
typename std::enable_if<HasMember_Rs<FluidState>::value,
const FluidState&>::type fluidState,
80 ->
decltype(decay<LhsEval>(fluidState.Rs()))
81{
return decay<LhsEval>(fluidState.Rs()); }
83template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
84LhsEval getRv_(
typename std::enable_if<!HasMember_Rv<FluidState>::value,
const FluidState&>::type fluidState,
88 decay<LhsEval>(fluidState.massFraction(FluidSystem::gasPhaseIdx, FluidSystem::oilCompIdx));
89 return FluidSystem::convertXgOToRv(XgO, regionIdx);
92template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
93auto getRv_(
typename std::enable_if<HasMember_Rv<FluidState>::value,
const FluidState&>::type fluidState,
95 ->
decltype(decay<LhsEval>(fluidState.Rv()))
96{
return decay<LhsEval>(fluidState.Rv()); }
98template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
99LhsEval getRvw_(
typename std::enable_if<!HasMember_Rvw<FluidState>::value,
const FluidState&>::type fluidState,
103 decay<LhsEval>(fluidState.massFraction(FluidSystem::gasPhaseIdx, FluidSystem::waterCompIdx));
104 return FluidSystem::convertXgWToRvw(XgW, regionIdx);
107template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
108auto getRvw_(
typename std::enable_if<HasMember_Rvw<FluidState>::value,
const FluidState&>::type fluidState,
110 ->
decltype(decay<LhsEval>(fluidState.Rvw()))
111{
return decay<LhsEval>(fluidState.Rvw()); }
113template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
114LhsEval getRsw_(
typename std::enable_if<!HasMember_Rsw<FluidState>::value,
const FluidState&>::type fluidState,
118 decay<LhsEval>(fluidState.massFraction(FluidSystem::waterPhaseIdx, FluidSystem::gasCompIdx));
119 return FluidSystem::convertXwGToRsw(XwG, regionIdx);
122template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
123auto getRsw_(
typename std::enable_if<HasMember_Rsw<FluidState>::value,
const FluidState&>::type fluidState,
125 ->
decltype(decay<LhsEval>(fluidState.Rsw()))
126{
return decay<LhsEval>(fluidState.Rsw()); }
128template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
129LhsEval getSaltConcentration_(
typename std::enable_if<!HasMember_saltConcentration<FluidState>::value,
130 const FluidState&>::type,
134template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
135auto getSaltConcentration_(
typename std::enable_if<HasMember_saltConcentration<FluidState>::value,
const FluidState&>::type fluidState,
137 ->
decltype(decay<LhsEval>(fluidState.saltConcentration()))
138{
return decay<LhsEval>(fluidState.saltConcentration()); }
140template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
141LhsEval getSaltSaturation_(
typename std::enable_if<!HasMember_saltSaturation<FluidState>::value,
142 const FluidState&>::type,
146template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
147auto getSaltSaturation_(
typename std::enable_if<HasMember_saltSaturation<FluidState>::value,
const FluidState&>::type fluidState,
149 ->
decltype(decay<LhsEval>(fluidState.saltSaturation()))
150{
return decay<LhsEval>(fluidState.saltSaturation()); }
160template <
class Scalar,
class IndexTraits = BlackOilDefaultIndexTraits>
171 template <
class EvaluationT>
174 using Evaluation = EvaluationT;
179 maxOilSat_ = maxOilSat;
180 regionIdx_ = regionIdx;
190 template <
class OtherCache>
193 regionIdx_ = other.regionIndex();
194 maxOilSat_ = other.maxOilSat();
205 {
return regionIdx_; }
215 { regionIdx_ = val; }
217 const Evaluation& maxOilSat()
const
218 {
return maxOilSat_; }
220 void setMaxOilSat(
const Evaluation& val)
221 { maxOilSat_ = val; }
224 Evaluation maxOilSat_;
235 static void initFromState(
const EclipseState& eclState,
const Schedule& schedule);
246 static void initBegin(std::size_t numPvtRegions);
255 { enableDissolvedGas_ = yesno; }
264 { enableVaporizedOil_ = yesno; }
273 { enableVaporizedWater_ = yesno; }
282 { enableDissolvedGasInWater_ = yesno; }
289 { enableDiffusion_ = yesno; }
296 { gasPvt_ = pvtObj; }
302 { oilPvt_ = pvtObj; }
308 { waterPvt_ = pvtObj; }
310 static void setVapPars(
const Scalar par1,
const Scalar par2)
313 gasPvt_->setVapPars(par1, par2);
316 oilPvt_->setVapPars(par1, par2);
319 waterPvt_->setVapPars(par1, par2);
340 static bool isInitialized()
341 {
return isInitialized_; }
364 static std::string_view
phaseName(
unsigned phaseIdx);
381 static constexpr unsigned oilCompIdx = IndexTraits::oilCompIdx;
385 static constexpr unsigned gasCompIdx = IndexTraits::gasCompIdx;
388 static unsigned char numActivePhases_;
389 static std::array<bool,numPhases> phaseIsActive_;
394 {
return numActivePhases_; }
400 return phaseIsActive_[phaseIdx];
414 {
return molarMass_[regionIdx][compIdx]; }
442 {
return molarMass_.size(); }
451 {
return enableDissolvedGas_; }
461 {
return enableDissolvedGasInWater_; }
470 {
return enableVaporizedOil_; }
479 {
return enableVaporizedWater_; }
487 {
return enableDiffusion_; }
495 {
return referenceDensity_[regionIdx][phaseIdx]; }
501 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
502 static LhsEval
density(
const FluidState& fluidState,
505 {
return density<FluidState, LhsEval>(fluidState, phaseIdx, paramCache.
regionIndex()); }
508 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
514 return fugacityCoefficient<FluidState, LhsEval>(fluidState,
521 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
525 {
return viscosity<FluidState, LhsEval>(fluidState, phaseIdx, paramCache.
regionIndex()); }
528 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
529 static LhsEval
enthalpy(
const FluidState& fluidState,
532 {
return enthalpy<FluidState, LhsEval>(fluidState, phaseIdx, paramCache.
regionIndex()); }
539 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
540 static LhsEval
density(
const FluidState& fluidState,
547 const LhsEval& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
548 const LhsEval& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
549 const LhsEval& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
555 const LhsEval& Rs = BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
556 const LhsEval& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
564 const LhsEval Rs(0.0);
565 const auto& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
573 const LhsEval& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
574 const LhsEval& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
575 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
584 const LhsEval Rvw(0.0);
585 const LhsEval& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
586 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
594 const LhsEval Rv(0.0);
595 const LhsEval& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
596 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
604 const LhsEval Rv(0.0);
605 const LhsEval Rvw(0.0);
606 const auto& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
613 const LhsEval& Rsw = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState,
waterPhaseIdx, regionIdx);
614 const LhsEval& bw = waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
619 const LhsEval Rsw(0.0);
622 * waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
625 throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
637 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
645 const auto& p = fluidState.pressure(phaseIdx);
646 const auto& T = fluidState.temperature(phaseIdx);
652 const LhsEval& Rs = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState,
oilPhaseIdx, regionIdx);
653 const LhsEval& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
661 const LhsEval Rs(0.0);
662 const LhsEval& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
669 const LhsEval& Rv = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState,
gasPhaseIdx, regionIdx);
670 const LhsEval& Rvw = saturatedVaporizationFactor<FluidState, LhsEval>(fluidState,
gasPhaseIdx, regionIdx);
671 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
681 const LhsEval Rvw(0.0);
682 const LhsEval& Rv = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState,
gasPhaseIdx, regionIdx);
683 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
692 const LhsEval Rv(0.0);
693 const LhsEval& Rvw = saturatedVaporizationFactor<FluidState, LhsEval>(fluidState,
gasPhaseIdx, regionIdx);
694 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
702 const LhsEval Rv(0.0);
703 const LhsEval Rvw(0.0);
704 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
714 const auto& saltConcentration = decay<LhsEval>(fluidState.saltConcentration());
715 const LhsEval& Rsw = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState,
waterPhaseIdx, regionIdx);
716 const LhsEval& bw = waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
723 *inverseFormationVolumeFactor<FluidState, LhsEval>(fluidState,
waterPhaseIdx, regionIdx);
727 throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
738 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
747 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
748 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
753 const auto& Rs = BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
755 && Rs >= (1.0 - 1e-10)*oilPvt_->saturatedGasDissolutionFactor(regionIdx, scalarValue(T), scalarValue(p)))
757 return oilPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
759 return oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
763 const LhsEval Rs(0.0);
764 return oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
768 const auto& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
769 const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
771 && Rvw >= (1.0 - 1e-10)*gasPvt_->saturatedWaterVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p))
773 && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
775 return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
777 return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
782 const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
784 && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
786 return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
788 const LhsEval Rvw(0.0);
789 return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
794 const auto& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
796 && Rvw >= (1.0 - 1e-10)*gasPvt_->saturatedWaterVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
798 return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
800 const LhsEval Rv(0.0);
801 return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
805 const LhsEval Rv(0.0);
806 const LhsEval Rvw(0.0);
807 return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
811 const auto& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
813 const auto& Rsw = BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
815 && Rsw >= (1.0 - 1e-10)*waterPvt_->saturatedGasDissolutionFactor(regionIdx, scalarValue(T), scalarValue(p), scalarValue(saltConcentration)))
817 return waterPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p, saltConcentration);
819 return waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
822 const LhsEval Rsw(0.0);
823 return waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
825 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
838 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
847 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
848 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
849 const auto& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
852 case oilPhaseIdx:
return oilPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
853 case gasPhaseIdx:
return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
854 case waterPhaseIdx:
return waterPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p, saltConcentration);
855 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
860 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
870 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
871 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
876 const LhsEval phi_oO = 20e3/p;
879 constexpr const Scalar phi_gG = 1.0;
883 const LhsEval phi_wW = 30e3/p;
901 const auto& R_vSat = gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p);
905 const auto& R_sSat = oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p);
908 const auto& x_oOSat = 1.0 - x_oGSat;
910 const auto& p_o = decay<LhsEval>(fluidState.pressure(
oilPhaseIdx));
911 const auto& p_g = decay<LhsEval>(fluidState.pressure(
gasPhaseIdx));
913 return phi_oO*p_o*x_oOSat / (p_g*x_gOSat);
921 throw std::logic_error(
"Invalid component index "+std::to_string(compIdx));
937 const auto& R_vSat = gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p);
940 const auto& x_gGSat = 1.0 - x_gOSat;
942 const auto& R_sSat = oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p);
946 const auto& p_o = decay<LhsEval>(fluidState.pressure(
oilPhaseIdx));
947 const auto& p_g = decay<LhsEval>(fluidState.pressure(
gasPhaseIdx));
949 return phi_gG*p_g*x_gGSat / (p_o*x_oGSat);
956 throw std::logic_error(
"Invalid component index "+std::to_string(compIdx));
971 throw std::logic_error(
"Invalid component index "+std::to_string(compIdx));
975 throw std::logic_error(
"Invalid phase index "+std::to_string(phaseIdx));
978 throw std::logic_error(
"Unhandled phase or component index");
982 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
991 const LhsEval& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
992 const LhsEval& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
997 const auto& Rs = BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
999 && Rs >= (1.0 - 1e-10)*oilPvt_->saturatedGasDissolutionFactor(regionIdx, scalarValue(T), scalarValue(p)))
1001 return oilPvt_->saturatedViscosity(regionIdx, T, p);
1003 return oilPvt_->viscosity(regionIdx, T, p, Rs);
1007 const LhsEval Rs(0.0);
1008 return oilPvt_->viscosity(regionIdx, T, p, Rs);
1013 const auto& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1014 const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1016 && Rvw >= (1.0 - 1e-10)*gasPvt_->saturatedWaterVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p))
1018 && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
1020 return gasPvt_->saturatedViscosity(regionIdx, T, p);
1022 return gasPvt_->viscosity(regionIdx, T, p, Rv, Rvw);
1026 const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1028 && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
1030 return gasPvt_->saturatedViscosity(regionIdx, T, p);
1032 const LhsEval Rvw(0.0);
1033 return gasPvt_->viscosity(regionIdx, T, p, Rv, Rvw);
1037 const auto& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1039 && Rvw >= (1.0 - 1e-10)*gasPvt_->saturatedWaterVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
1041 return gasPvt_->saturatedViscosity(regionIdx, T, p);
1043 const LhsEval Rv(0.0);
1044 return gasPvt_->viscosity(regionIdx, T, p, Rv, Rvw);
1048 const LhsEval Rv(0.0);
1049 const LhsEval Rvw(0.0);
1050 return gasPvt_->viscosity(regionIdx, T, p, Rv, Rvw);
1055 const LhsEval& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1057 const auto& Rsw = BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1059 && Rsw >= (1.0 - 1e-10)*waterPvt_->saturatedGasDissolutionFactor(regionIdx, scalarValue(T), scalarValue(p), scalarValue(saltConcentration)))
1061 return waterPvt_->saturatedViscosity(regionIdx, T, p, saltConcentration);
1063 return waterPvt_->viscosity(regionIdx, T, p, Rsw, saltConcentration);
1066 const LhsEval Rsw(0.0);
1067 return waterPvt_->viscosity(regionIdx, T, p, Rsw, saltConcentration);
1071 throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1075 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1083 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1084 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1089 oilPvt_->internalEnergy(regionIdx, T, p, BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx))
1090 + p/density<FluidState, LhsEval>(fluidState, phaseIdx, regionIdx);
1094 gasPvt_->internalEnergy(regionIdx, T, p,
1095 BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1096 BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx))
1097 + p/density<FluidState, LhsEval>(fluidState, phaseIdx, regionIdx);
1101 waterPvt_->internalEnergy(regionIdx, T, p,
1102 BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1103 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx))
1104 + p/density<FluidState, LhsEval>(fluidState, phaseIdx, regionIdx);
1106 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1109 throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1118 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1126 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1127 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1128 const auto& saltConcentration = decay<LhsEval>(fluidState.saltConcentration());
1132 case gasPhaseIdx:
return gasPvt_->saturatedWaterVaporizationFactor(regionIdx, T, p, saltConcentration);
1134 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1144 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1148 const LhsEval& maxOilSaturation)
1154 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1155 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1156 const auto& So = decay<LhsEval>(fluidState.saturation(
oilPhaseIdx));
1159 case oilPhaseIdx:
return oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p, So, maxOilSaturation);
1160 case gasPhaseIdx:
return gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p, So, maxOilSaturation);
1161 case waterPhaseIdx:
return waterPvt_->saturatedGasDissolutionFactor(regionIdx, T, p,
1162 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1163 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1175 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1184 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1185 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1188 case oilPhaseIdx:
return oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p);
1189 case gasPhaseIdx:
return gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p);
1190 case waterPhaseIdx:
return waterPvt_->saturatedGasDissolutionFactor(regionIdx, T, p,
1191 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1192 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1199 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1210 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1227 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1235 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1238 case oilPhaseIdx:
return oilPvt_->saturationPressure(regionIdx, T, BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1239 case gasPhaseIdx:
return gasPvt_->saturationPressure(regionIdx, T, BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1240 case waterPhaseIdx:
return waterPvt_->saturationPressure(regionIdx, T,
1241 BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1242 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1243 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1254 template <
class LhsEval>
1260 return XoG/(1.0 - XoG)*(rho_oRef/rho_gRef);
1267 template <
class LhsEval>
1273 return XwG/(1.0 - XwG)*(rho_wRef/rho_gRef);
1280 template <
class LhsEval>
1286 return XgO/(1.0 - XgO)*(rho_gRef/rho_oRef);
1293 template <
class LhsEval>
1299 return XgW/(1.0 - XgW)*(rho_gRef/rho_wRef);
1307 template <
class LhsEval>
1313 const LhsEval& rho_oG = Rs*rho_gRef;
1314 return rho_oG/(rho_oRef + rho_oG);
1321 template <
class LhsEval>
1327 const LhsEval& rho_wG = Rsw*rho_gRef;
1328 return rho_wG/(rho_wRef + rho_wG);
1335 template <
class LhsEval>
1341 const LhsEval& rho_gO = Rv*rho_oRef;
1342 return rho_gO/(rho_gRef + rho_gO);
1349 template <
class LhsEval>
1355 const LhsEval& rho_gW = Rvw*rho_wRef;
1356 return rho_gW/(rho_gRef + rho_gW);
1362 template <
class LhsEval>
1368 return XgW*MG / (MW*(1 - XgW) + XgW*MG);
1374 template <
class LhsEval>
1380 return XwG*MW / (MG*(1 - XwG) + XwG*MW);
1386 template <
class LhsEval>
1392 return XoG*MO / (MG*(1 - XoG) + XoG*MO);
1398 template <
class LhsEval>
1404 return xoG*MG / (xoG*(MG - MO) + MO);
1410 template <
class LhsEval>
1416 return XgO*MG / (MO*(1 - XgO) + XgO*MG);
1422 template <
class LhsEval>
1428 return xgO*MO / (xgO*(MO - MG) + MG);
1439 {
return *gasPvt_; }
1449 {
return *oilPvt_; }
1459 {
return *waterPvt_; }
1467 {
return reservoirTemperature_; }
1475 { reservoirTemperature_ = value; }
1477 static short activeToCanonicalPhaseIdx(
unsigned activePhaseIdx);
1479 static short canonicalToActivePhaseIdx(
unsigned phaseIdx);
1483 {
return diffusionCoefficients_[regionIdx][
numPhases*compIdx + phaseIdx]; }
1487 { diffusionCoefficients_[regionIdx][
numPhases*compIdx + phaseIdx] = coefficient ; }
1492 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
1503 if(!diffusionCoefficients_.empty()) {
1507 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1508 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1514 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1519 static void resizeArrays_(std::size_t
numRegions);
1521 static Scalar reservoirTemperature_;
1523 static std::shared_ptr<GasPvt> gasPvt_;
1524 static std::shared_ptr<OilPvt> oilPvt_;
1525 static std::shared_ptr<WaterPvt> waterPvt_;
1527 static bool enableDissolvedGas_;
1528 static bool enableDissolvedGasInWater_;
1529 static bool enableVaporizedOil_;
1530 static bool enableVaporizedWater_;
1531 static bool enableDiffusion_;
1536 static std::vector<std::array<
Scalar, 3> > referenceDensity_;
1537 static std::vector<std::array<
Scalar, 3> > molarMass_;
1538 static std::vector<std::array<
Scalar, 3 * 3> > diffusionCoefficients_;
1540 static std::array<short, numPhases> activeToCanonicalPhaseIdx_;
1541 static std::array<short, numPhases> canonicalToActivePhaseIdx_;
1543 static bool isInitialized_;
1546template <
class Scalar,
class IndexTraits>
1547unsigned char BlackOilFluidSystem<Scalar, IndexTraits>::numActivePhases_;
1549template <
class Scalar,
class IndexTraits>
1550std::array<bool, BlackOilFluidSystem<Scalar, IndexTraits>::numPhases> BlackOilFluidSystem<Scalar, IndexTraits>::phaseIsActive_;
1552template <
class Scalar,
class IndexTraits>
1553std::array<short, BlackOilFluidSystem<Scalar, IndexTraits>::numPhases> BlackOilFluidSystem<Scalar, IndexTraits>::activeToCanonicalPhaseIdx_;
1555template <
class Scalar,
class IndexTraits>
1556std::array<short, BlackOilFluidSystem<Scalar, IndexTraits>::numPhases> BlackOilFluidSystem<Scalar, IndexTraits>::canonicalToActivePhaseIdx_;
1558template <
class Scalar,
class IndexTraits>
1562template <
class Scalar,
class IndexTraits>
1566template <
class Scalar,
class IndexTraits>
1568BlackOilFluidSystem<Scalar, IndexTraits>::reservoirTemperature_;
1570template <
class Scalar,
class IndexTraits>
1571bool BlackOilFluidSystem<Scalar, IndexTraits>::enableDissolvedGas_;
1573template <
class Scalar,
class IndexTraits>
1574bool BlackOilFluidSystem<Scalar, IndexTraits>::enableDissolvedGasInWater_;
1577template <
class Scalar,
class IndexTraits>
1578bool BlackOilFluidSystem<Scalar, IndexTraits>::enableVaporizedOil_;
1580template <
class Scalar,
class IndexTraits>
1581bool BlackOilFluidSystem<Scalar, IndexTraits>::enableVaporizedWater_;
1583template <
class Scalar,
class IndexTraits>
1584bool BlackOilFluidSystem<Scalar, IndexTraits>::enableDiffusion_;
1586template <
class Scalar,
class IndexTraits>
1587std::shared_ptr<OilPvtMultiplexer<Scalar> >
1588BlackOilFluidSystem<Scalar, IndexTraits>::oilPvt_;
1590template <
class Scalar,
class IndexTraits>
1591std::shared_ptr<GasPvtMultiplexer<Scalar> >
1592BlackOilFluidSystem<Scalar, IndexTraits>::gasPvt_;
1594template <
class Scalar,
class IndexTraits>
1595std::shared_ptr<WaterPvtMultiplexer<Scalar> >
1596BlackOilFluidSystem<Scalar, IndexTraits>::waterPvt_;
1598template <
class Scalar,
class IndexTraits>
1599std::vector<std::array<Scalar, 3> >
1600BlackOilFluidSystem<Scalar, IndexTraits>::referenceDensity_;
1602template <
class Scalar,
class IndexTraits>
1603std::vector<std::array<Scalar, 3> >
1604BlackOilFluidSystem<Scalar, IndexTraits>::molarMass_;
1606template <
class Scalar,
class IndexTraits>
1607std::vector<std::array<Scalar, 9> >
1608BlackOilFluidSystem<Scalar, IndexTraits>::diffusionCoefficients_;
1610template <
class Scalar,
class IndexTraits>
1611bool BlackOilFluidSystem<Scalar, IndexTraits>::isInitialized_ =
false;
The base class for all fluid systems.
This class represents the Pressure-Volume-Temperature relations of the liquid phase for a CO2-Brine s...
This class represents the Pressure-Volume-Temperature relations of the liquid phase for a H2-Brine sy...
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
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:47
ScalarT Scalar
The type used for scalar quantities.
Definition BaseFluidSystem.hpp:52
A fluid system which uses the black-oil model assumptions to calculate termodynamically meaningful qu...
Definition BlackOilFluidSystem.hpp:162
static constexpr unsigned waterCompIdx
Index of the water component.
Definition BlackOilFluidSystem.hpp:383
static unsigned numActivePhases()
Returns the number of active fluid phases (i.e., usually three)
Definition BlackOilFluidSystem.hpp:393
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:1482
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:1255
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:1308
static constexpr unsigned numComponents
Number of chemical species in the fluid system.
Definition BlackOilFluidSystem.hpp:378
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:1350
static LhsEval convertXwGToxwG(const LhsEval &XwG, unsigned regionIdx)
Convert a gas mass fraction in the water phase the corresponding mole fraction.
Definition BlackOilFluidSystem.hpp:1375
static LhsEval saturatedInverseFormationVolumeFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the formation volume factor of a "saturated" fluid phase.
Definition BlackOilFluidSystem.hpp:839
static LhsEval convertXoGToxoG(const LhsEval &XoG, unsigned regionIdx)
Convert a gas mass fraction in the oil phase the corresponding mole fraction.
Definition BlackOilFluidSystem.hpp:1387
static LhsEval viscosity(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Calculate the dynamic viscosity of a fluid phase [Pa*s].
Definition BlackOilFluidSystem.hpp:983
static void initEnd()
Finish initializing the black oil fluid system.
Definition BlackOilFluidSystem.cpp:260
static Scalar reservoirTemperature(unsigned=0)
Set the temperature of the reservoir.
Definition BlackOilFluidSystem.hpp:1466
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:1322
static bool phaseIsActive(unsigned phaseIdx)
Returns whether a fluid phase is active.
Definition BlackOilFluidSystem.hpp:397
static void setEnableDiffusion(bool yesno)
Specify whether the fluid system should consider diffusion.
Definition BlackOilFluidSystem.hpp:288
static LhsEval bubblePointPressure(const FluidState &fluidState, unsigned regionIdx)
Returns the bubble point pressure $P_b$ using the current Rs.
Definition BlackOilFluidSystem.hpp:1200
static LhsEval density(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Calculate the density [kg/m^3] of a fluid phase.
Definition BlackOilFluidSystem.hpp:540
static LhsEval saturatedDensity(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Compute the density of a saturated fluid phase.
Definition BlackOilFluidSystem.hpp:638
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:522
static bool isIdealMixture(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal mixture.
Definition BlackOilFluidSystem.hpp:417
static bool enableVaporizedWater()
Returns whether the fluid system should consider that the water component can dissolve in the gas pha...
Definition BlackOilFluidSystem.hpp:478
static constexpr unsigned numPhases
Number of fluid phases in the fluid system.
Definition BlackOilFluidSystem.hpp:348
static void setWaterPvt(std::shared_ptr< WaterPvt > pvtObj)
Set the pressure-volume-saturation (PVT) relations for the water phase.
Definition BlackOilFluidSystem.hpp:307
static void initBegin(std::size_t numPvtRegions)
Begin the initialization of the black oil fluid system.
Definition BlackOilFluidSystem.cpp:223
static void setReservoirTemperature(Scalar value)
Return the temperature of the reservoir.
Definition BlackOilFluidSystem.hpp:1474
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:281
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:1145
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:1076
static constexpr unsigned oilCompIdx
Index of the oil component.
Definition BlackOilFluidSystem.hpp:381
static bool isLiquid(unsigned phaseIdx)
Return whether a phase is liquid.
Definition BlackOilFluidSystem.hpp:367
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:529
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:861
static LhsEval convertxoGToXoG(const LhsEval &xoG, unsigned regionIdx)
Convert a gas mole fraction in the oil phase the corresponding mass fraction.
Definition BlackOilFluidSystem.hpp:1399
static Scalar molarMass(unsigned compIdx, unsigned regionIdx=0)
Return the molar mass of a component in [kg/mol].
Definition BlackOilFluidSystem.hpp:413
static LhsEval saturatedDissolutionFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the dissolution factor of a saturated fluid phase.
Definition BlackOilFluidSystem.hpp:1176
static LhsEval dewPointPressure(const FluidState &fluidState, unsigned regionIdx)
Returns the dew point pressure $P_d$ using the current Rv.
Definition BlackOilFluidSystem.hpp:1211
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:1228
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:1268
static std::string_view componentName(unsigned compIdx)
Return the human readable name of a component.
Definition BlackOilFluidSystem.cpp:355
static LhsEval convertXgOToxgO(const LhsEval &XgO, unsigned regionIdx)
Convert a oil mass fraction in the gas phase the corresponding mole fraction.
Definition BlackOilFluidSystem.hpp:1411
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:1493
static LhsEval saturatedVaporizationFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the water vaporization factor of saturated phase.
Definition BlackOilFluidSystem.hpp:1119
static bool enableDissolvedGas()
Returns whether the fluid system should consider that the gas component can dissolve in the oil phase...
Definition BlackOilFluidSystem.hpp:450
static Scalar referenceDensity(unsigned phaseIdx, unsigned regionIdx)
Returns the density of a fluid phase at surface pressure [kg/m^3].
Definition BlackOilFluidSystem.hpp:494
static Scalar surfaceTemperature
The temperature at the surface.
Definition BlackOilFluidSystem.hpp:361
static void setOilPvt(std::shared_ptr< OilPvt > pvtObj)
Set the pressure-volume-saturation (PVT) relations for the oil phase.
Definition BlackOilFluidSystem.hpp:301
static constexpr unsigned waterPhaseIdx
Index of the water phase.
Definition BlackOilFluidSystem.hpp:351
static LhsEval convertxgOToXgO(const LhsEval &xgO, unsigned regionIdx)
Convert a oil mole fraction in the gas phase the corresponding mass fraction.
Definition BlackOilFluidSystem.hpp:1423
static LhsEval inverseFormationVolumeFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the formation volume factor of an "undersaturated" fluid phase.
Definition BlackOilFluidSystem.hpp:739
static LhsEval convertXgWToxgW(const LhsEval &XgW, unsigned regionIdx)
Convert a water mass fraction in the gas phase the corresponding mole fraction.
Definition BlackOilFluidSystem.hpp:1363
static constexpr unsigned gasPhaseIdx
Index of the gas phase.
Definition BlackOilFluidSystem.hpp:355
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:254
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:1281
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:509
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:502
static constexpr unsigned gasCompIdx
Index of the gas component.
Definition BlackOilFluidSystem.hpp:385
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:1336
static const WaterPvt & waterPvt()
Return a reference to the low-level object which calculates the water phase quantities.
Definition BlackOilFluidSystem.hpp:1458
static bool enableVaporizedOil()
Returns whether the fluid system should consider that the oil component can dissolve in the gas phase...
Definition BlackOilFluidSystem.hpp:469
static bool isIdealGas(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal gas.
Definition BlackOilFluidSystem.hpp:429
static unsigned solventComponentIndex(unsigned phaseIdx)
returns the index of "primary" component of a phase (solvent)
Definition BlackOilFluidSystem.cpp:316
static Scalar surfacePressure
The pressure at the surface.
Definition BlackOilFluidSystem.hpp:358
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:272
static constexpr unsigned oilPhaseIdx
Index of the oil phase.
Definition BlackOilFluidSystem.hpp:353
static const OilPvt & oilPvt()
Return a reference to the low-level object which calculates the oil phase quantities.
Definition BlackOilFluidSystem.hpp:1448
static bool enableDissolvedGasInWater()
Returns whether the fluid system should consider that the gas component can dissolve in the water pha...
Definition BlackOilFluidSystem.hpp:460
static unsigned soluteComponentIndex(unsigned phaseIdx)
returns the index of "secondary" component of a phase (solute)
Definition BlackOilFluidSystem.cpp:333
static bool enableDiffusion()
Returns whether the fluid system should consider diffusion.
Definition BlackOilFluidSystem.hpp:486
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:1294
static void setReferenceDensities(Scalar rhoOil, Scalar rhoWater, Scalar rhoGas, unsigned regionIdx)
Initialize the values of the reference densities.
Definition BlackOilFluidSystem.cpp:249
static std::size_t numRegions()
Returns the number of PVT regions which are considered.
Definition BlackOilFluidSystem.hpp:441
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:263
static std::string_view phaseName(unsigned phaseIdx)
Return the human readable name of a fluid phase.
Definition BlackOilFluidSystem.cpp:299
static void setDiffusionCoefficient(Scalar coefficient, unsigned compIdx, unsigned phaseIdx, unsigned regionIdx=0)
Definition BlackOilFluidSystem.hpp:1486
static const GasPvt & gasPvt()
Return a reference to the low-level object which calculates the gas phase quantities.
Definition BlackOilFluidSystem.hpp:1438
static bool isCompressible(unsigned)
Returns true if and only if a fluid phase is assumed to be compressible.
Definition BlackOilFluidSystem.hpp:425
static void setGasPvt(std::shared_ptr< GasPvt > pvtObj)
Set the pressure-volume-saturation (PVT) relations for the gas phase.
Definition BlackOilFluidSystem.hpp:295
This class represents the Pressure-Volume-Temperature relations of the gas phase in the black-oil mod...
Definition GasPvtMultiplexer.hpp:109
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:336
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:104
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:269
This class represents the Pressure-Volume-Temperature relations of the water phase in the black-oil m...
Definition WaterPvtMultiplexer.hpp:89
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:261
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:173
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:191
void setRegionIndex(unsigned val)
Set the index of the region which should be used to determine the thermodynamic properties.
Definition BlackOilFluidSystem.hpp:214
unsigned regionIndex() const
Return the index of the region which should be used to determine the thermodynamic properties.
Definition BlackOilFluidSystem.hpp:204