33#include <opm/common/TimingMacros.hpp>
34#include <opm/common/utility/VectorWithDefaultAllocator.hpp>
69template <
class Scalar,
class IndexTraits = BlackOilDefaultIndexTraits,
70 template<
typename>
typename Storage = VectorWithDefaultAllocator,
71 template<
typename>
typename SmartPointer = std::shared_ptr>
81 #ifdef COMPILING_STATIC_FLUID_SYSTEM
83 template <
class EvaluationT>
86 using Evaluation = EvaluationT;
89 explicit ParameterCache(
Scalar maxOilSat = 1.0,
unsigned regionIdx = 0)
90 : maxOilSat_(maxOilSat)
91 , regionIdx_(regionIdx)
102 template <
class OtherCache>
103 void assignPersistentData(
const OtherCache& other)
105 regionIdx_ = other.regionIndex();
106 maxOilSat_ = other.maxOilSat();
116 unsigned regionIndex()
const
117 {
return regionIdx_; }
126 void setRegionIndex(
unsigned val)
127 { regionIdx_ = val; }
129 const Evaluation& maxOilSat()
const
130 {
return maxOilSat_; }
132 void setMaxOilSat(
const Evaluation& val)
133 { maxOilSat_ = val; }
136 Evaluation maxOilSat_;
143 template<
class EvaluationT>
144 using ParameterCache =
148 SmartPointer>::template ParameterCache<EvaluationT>;
151 #ifndef COMPILING_STATIC_FLUID_SYSTEM
157 template<
template<
typename>
typename StorageT,
template<
typename>
typename SmartPointerT>
161 , numActivePhases_(other.numActivePhases_)
162 , phaseIsActive_(other.phaseIsActive_)
163 , reservoirTemperature_(other.reservoirTemperature_)
164 , gasPvt_(SmartPointerT<typename decltype(gasPvt_)::element_type>(other.gasPvt_))
165 , oilPvt_(SmartPointerT<typename decltype(oilPvt_)::element_type>(other.oilPvt_))
166 , waterPvt_(SmartPointerT<typename decltype(waterPvt_)::element_type>(other.waterPvt_))
167 , enableDissolvedGas_(other.enableDissolvedGas_)
168 , enableDissolvedGasInWater_(other.enableDissolvedGasInWater_)
169 , enableVaporizedOil_(other.enableVaporizedOil_)
170 , enableVaporizedWater_(other.enableVaporizedWater_)
171 , enableDiffusion_(other.enableDiffusion_)
172 , referenceDensity_(StorageT<typename decltype(referenceDensity_)::value_type>(other.referenceDensity_))
173 , molarMass_(StorageT<typename decltype(molarMass_)::value_type>(other.molarMass_))
174 , diffusionCoefficients_(StorageT<typename decltype(diffusionCoefficients_)::value_type>(other.diffusionCoefficients_))
175 , activeToCanonicalPhaseIdx_(other.activeToCanonicalPhaseIdx_)
176 , canonicalToActivePhaseIdx_(other.canonicalToActivePhaseIdx_)
177 , isInitialized_(other.isInitialized_)
178 , useSaturatedTables_(other.useSaturatedTables_)
179 , enthalpy_eq_energy_(other.enthalpy_eq_energy_)
184 #ifdef COMPILING_STATIC_FLUID_SYSTEM
192 template<
template<
typename>
typename StorageT = VectorWithDefaultAllocator,
template<
typename>
typename SmartPointerT = std::shared_ptr>
208 STATIC_OR_DEVICE
void initFromState(
const EclipseState& eclState,
const Schedule& schedule);
219 STATIC_OR_DEVICE
void initBegin(std::size_t numPvtRegions);
228 { enableDissolvedGas_ = yesno; }
237 { enableVaporizedOil_ = yesno; }
246 { enableVaporizedWater_ = yesno; }
255 { enableDissolvedGasInWater_ = yesno; }
262 { enableDiffusion_ = yesno; }
270 { useSaturatedTables_ = yesno; }
275 STATIC_OR_DEVICE
void setGasPvt(SmartPointer<GasPvt> pvtObj)
276 { gasPvt_ = pvtObj; }
281 STATIC_OR_DEVICE
void setOilPvt(SmartPointer<OilPvt> pvtObj)
282 { oilPvt_ = pvtObj; }
288 { waterPvt_ = pvtObj; }
290 STATIC_OR_DEVICE
void setVapPars(
const Scalar par1,
const Scalar par2)
293 gasPvt_->setVapPars(par1, par2);
296 oilPvt_->setVapPars(par1, par2);
299 waterPvt_->setVapPars(par1, par2);
318 STATIC_OR_DEVICE
void initEnd();
320 STATIC_OR_DEVICE
bool isInitialized()
321 {
return isInitialized_; }
344 STATIC_OR_DEVICE std::string_view
phaseName(
unsigned phaseIdx);
368 STATIC_OR_NOTHING
unsigned char numActivePhases_;
369 STATIC_OR_NOTHING std::array<bool,numPhases> phaseIsActive_;
374 {
return numActivePhases_; }
380 return phaseIsActive_[phaseIdx];
390 STATIC_OR_DEVICE std::string_view
componentName(
unsigned compIdx);
394 {
return molarMass_[regionIdx][compIdx]; }
422 {
return molarMass_.size(); }
431 {
return enableDissolvedGas_; }
441 {
return enableDissolvedGasInWater_; }
450 {
return enableVaporizedOil_; }
459 {
return enableVaporizedWater_; }
467 {
return enableDiffusion_; }
475 {
return useSaturatedTables_; }
483 {
return referenceDensity_[regionIdx][phaseIdx]; }
489 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
490 STATIC_OR_DEVICE LhsEval
density(
const FluidState& fluidState,
491 const ParameterCache<ParamCacheEval>& paramCache,
493 {
return density<FluidState, LhsEval>(fluidState, phaseIdx, paramCache.regionIndex()); }
496 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
498 const ParameterCache<ParamCacheEval>& paramCache,
502 return fugacityCoefficient<FluidState, LhsEval>(fluidState,
505 paramCache.regionIndex());
509 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
510 STATIC_OR_DEVICE LhsEval
viscosity(
const FluidState& fluidState,
511 const ParameterCache<ParamCacheEval>& paramCache,
513 {
return viscosity<FluidState, LhsEval>(fluidState, phaseIdx, paramCache.regionIndex()); }
516 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
517 STATIC_OR_DEVICE LhsEval
enthalpy(
const FluidState& fluidState,
518 const ParameterCache<ParamCacheEval>& paramCache,
520 {
return enthalpy<FluidState, LhsEval>(fluidState, phaseIdx, paramCache.regionIndex()); }
522 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
523 STATIC_OR_DEVICE LhsEval internalEnergy(
const FluidState& fluidState,
524 const ParameterCache<ParamCacheEval>& paramCache,
526 {
return internalEnergy<FluidState, LhsEval>(fluidState, phaseIdx, paramCache.regionIndex()); }
533 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
534 STATIC_OR_DEVICE LhsEval
density(
const FluidState& fluidState,
541 const LhsEval& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
542 const LhsEval& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
543 const LhsEval& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
549 const LhsEval& Rs = BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
550 const LhsEval& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
558 const LhsEval Rs(0.0);
559 const auto& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
567 const LhsEval& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
568 const LhsEval& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
569 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
578 const LhsEval Rvw(0.0);
579 const LhsEval& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
580 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
588 const LhsEval Rv(0.0);
589 const LhsEval& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
590 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
598 const LhsEval Rv(0.0);
599 const LhsEval Rvw(0.0);
600 const auto& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
607 const LhsEval& Rsw =BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
608 const LhsEval& bw = waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
613 const LhsEval Rsw(0.0);
616 * waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
619 throw std::logic_error(
"Unhandled phase index " + std::to_string(phaseIdx));
631 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
639 const auto& p = fluidState.pressure(phaseIdx);
640 const auto& T = fluidState.temperature(phaseIdx);
646 const LhsEval& Rs = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState,
oilPhaseIdx, regionIdx);
647 const LhsEval& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
655 const LhsEval Rs(0.0);
656 const LhsEval& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
663 const LhsEval& Rv = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState,
gasPhaseIdx, regionIdx);
664 const LhsEval& Rvw = saturatedVaporizationFactor<FluidState, LhsEval>(fluidState,
gasPhaseIdx, regionIdx);
665 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
675 const LhsEval Rvw(0.0);
676 const LhsEval& Rv = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState,
gasPhaseIdx, regionIdx);
677 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
686 const LhsEval Rv(0.0);
687 const LhsEval& Rvw = saturatedVaporizationFactor<FluidState, LhsEval>(fluidState,
gasPhaseIdx, regionIdx);
688 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
696 const LhsEval Rv(0.0);
697 const LhsEval Rvw(0.0);
698 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
708 const auto& saltConcentration = decay<LhsEval>(fluidState.saltConcentration());
709 const LhsEval& Rsw = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState,
waterPhaseIdx, regionIdx);
710 const LhsEval& bw = waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
717 *inverseFormationVolumeFactor<FluidState, LhsEval>(fluidState,
waterPhaseIdx, regionIdx);
721 throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
732 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
741 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
742 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
747 const auto& Rs = BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
749 && Rs >= (1.0 - 1e-10)*oilPvt_->saturatedGasDissolutionFactor(regionIdx, scalarValue(T), scalarValue(p)))
751 return oilPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
753 return oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
757 const LhsEval Rs(0.0);
758 return oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
762 const auto& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
763 const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
765 && Rvw >= (1.0 - 1e-10)*gasPvt_->saturatedWaterVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p))
767 && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
769 return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
771 return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
776 const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
778 && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
780 return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
782 const LhsEval Rvw(0.0);
783 return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
788 const auto& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
790 && Rvw >= (1.0 - 1e-10)*gasPvt_->saturatedWaterVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
792 return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
794 const LhsEval Rv(0.0);
795 return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
799 const LhsEval Rv(0.0);
800 const LhsEval Rvw(0.0);
801 return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
805 const auto& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
807 const auto& Rsw = BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
809 && Rsw >= (1.0 - 1e-10)*waterPvt_->saturatedGasDissolutionFactor(regionIdx, scalarValue(T), scalarValue(p), scalarValue(saltConcentration)))
811 return waterPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p, saltConcentration);
813 return waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
816 const LhsEval Rsw(0.0);
817 return waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
819 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
832 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
841 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
842 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
843 const auto& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
846 case oilPhaseIdx:
return oilPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
847 case gasPhaseIdx:
return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
848 case waterPhaseIdx:
return waterPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p, saltConcentration);
849 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
854 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
864 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
865 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
870 const LhsEval phi_oO = 20e3/p;
873 constexpr const Scalar phi_gG = 1.0;
877 const LhsEval phi_wW = 30e3/p;
895 const auto& R_vSat = gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p);
899 const auto& R_sSat = oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p);
902 const auto& x_oOSat = 1.0 - x_oGSat;
904 const auto& p_o = decay<LhsEval>(fluidState.pressure(
oilPhaseIdx));
905 const auto& p_g = decay<LhsEval>(fluidState.pressure(
gasPhaseIdx));
907 return phi_oO*p_o*x_oOSat / (p_g*x_gOSat);
915 throw std::logic_error(
"Invalid component index "+std::to_string(compIdx));
931 const auto& R_vSat = gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p);
934 const auto& x_gGSat = 1.0 - x_gOSat;
936 const auto& R_sSat = oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p);
940 const auto& p_o = decay<LhsEval>(fluidState.pressure(
oilPhaseIdx));
941 const auto& p_g = decay<LhsEval>(fluidState.pressure(
gasPhaseIdx));
943 return phi_gG*p_g*x_gGSat / (p_o*x_oGSat);
950 throw std::logic_error(
"Invalid component index "+std::to_string(compIdx));
965 throw std::logic_error(
"Invalid component index "+std::to_string(compIdx));
969 throw std::logic_error(
"Invalid phase index "+std::to_string(phaseIdx));
972 throw std::logic_error(
"Unhandled phase or component index");
976 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
977 STATIC_OR_DEVICE LhsEval
viscosity(
const FluidState& fluidState,
985 const LhsEval& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
986 const LhsEval& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
991 const auto& Rs = BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
993 && Rs >= (1.0 - 1e-10)*oilPvt_->saturatedGasDissolutionFactor(regionIdx, scalarValue(T), scalarValue(p)))
995 return oilPvt_->saturatedViscosity(regionIdx, T, p);
997 return oilPvt_->viscosity(regionIdx, T, p, Rs);
1001 const LhsEval Rs(0.0);
1002 return oilPvt_->viscosity(regionIdx, T, p, Rs);
1007 const auto& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1008 const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1010 && Rvw >= (1.0 - 1e-10)*gasPvt_->saturatedWaterVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p))
1012 && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
1014 return gasPvt_->saturatedViscosity(regionIdx, T, p);
1016 return gasPvt_->viscosity(regionIdx, T, p, Rv, Rvw);
1020 const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1022 && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
1024 return gasPvt_->saturatedViscosity(regionIdx, T, p);
1026 const LhsEval Rvw(0.0);
1027 return gasPvt_->viscosity(regionIdx, T, p, Rv, Rvw);
1031 const auto& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1033 && Rvw >= (1.0 - 1e-10)*gasPvt_->saturatedWaterVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
1035 return gasPvt_->saturatedViscosity(regionIdx, T, p);
1037 const LhsEval Rv(0.0);
1038 return gasPvt_->viscosity(regionIdx, T, p, Rv, Rvw);
1042 const LhsEval Rv(0.0);
1043 const LhsEval Rvw(0.0);
1044 return gasPvt_->viscosity(regionIdx, T, p, Rv, Rvw);
1049 const LhsEval& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1051 const auto& Rsw = BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1053 && Rsw >= (1.0 - 1e-10)*waterPvt_->saturatedGasDissolutionFactor(regionIdx, scalarValue(T), scalarValue(p), scalarValue(saltConcentration)))
1055 return waterPvt_->saturatedViscosity(regionIdx, T, p, saltConcentration);
1057 return waterPvt_->viscosity(regionIdx, T, p, Rsw, saltConcentration);
1060 const LhsEval Rsw(0.0);
1061 return waterPvt_->viscosity(regionIdx, T, p, Rsw, saltConcentration);
1065 throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1068 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1069 STATIC_OR_DEVICE LhsEval internalEnergy(
const FluidState& fluidState,
1070 const unsigned phaseIdx,
1071 const unsigned regionIdx)
1073 const auto p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1074 const auto T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1078 if (!oilPvt_->mixingEnergy()) {
1079 return oilPvt_->internalEnergy
1081 BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1086 if (!waterPvt_->mixingEnergy()) {
1087 return waterPvt_->internalEnergy
1089 BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1090 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1095 if (!gasPvt_->mixingEnergy()) {
1096 return gasPvt_->internalEnergy
1098 BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1099 BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1104 throw std::logic_error {
1105 "Phase index " + std::to_string(phaseIdx) +
" does not support internal energy"
1109 return internalMixingTotalEnergy<FluidState,LhsEval>(fluidState, phaseIdx, regionIdx)
1110 / density<FluidState,LhsEval>(fluidState, phaseIdx, regionIdx);
1114 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1115 STATIC_OR_DEVICE LhsEval internalMixingTotalEnergy(
const FluidState& fluidState,
1121 const LhsEval& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1122 const LhsEval& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1123 const LhsEval& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1127 auto oilEnergy = oilPvt_->internalEnergy(regionIdx, T, p,
1128 BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1129 assert(oilPvt_->mixingEnergy());
1133 const LhsEval& Rs = BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1134 const LhsEval& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
1135 const auto& gasEnergy =
1136 gasPvt_->internalEnergy(regionIdx, T, p,
1137 BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1138 BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1139 const auto hVapG = gasPvt_->hVap(regionIdx);
1146 const LhsEval Rs(0.0);
1147 const auto& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
1153 const auto& gasEnergy =
1154 gasPvt_->internalEnergy(regionIdx, T, p,
1155 BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1156 BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1157 assert(gasPvt_->mixingEnergy());
1159 const auto& oilEnergy =
1160 oilPvt_->internalEnergy(regionIdx, T, p,
1161 BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1162 const auto waterEnergy =
1163 waterPvt_->internalEnergy(regionIdx, T, p,
1164 BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1165 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1167 const LhsEval& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1168 const LhsEval& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1169 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
1170 const auto hVapO = oilPvt_->hVap(regionIdx);
1171 const auto hVapW = waterPvt_->hVap(regionIdx);
1178 const auto& oilEnergy =
1179 oilPvt_->internalEnergy(regionIdx, T, p,
1180 BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1182 const LhsEval Rvw(0.0);
1183 const LhsEval& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1184 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
1185 const auto hVapO = oilPvt_->hVap(regionIdx);
1192 const LhsEval Rv(0.0);
1193 const LhsEval& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1194 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
1195 const auto waterEnergy =
1196 waterPvt_->internalEnergy(regionIdx, T, p,
1197 BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1198 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1199 const auto hVapW = waterPvt_->hVap(regionIdx);
1206 const LhsEval Rv(0.0);
1207 const LhsEval Rvw(0.0);
1208 const auto& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
1213 const auto waterEnergy =
1214 waterPvt_->internalEnergy(regionIdx, T, p,
1215 BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1216 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1217 assert(waterPvt_->mixingEnergy());
1219 const auto& gasEnergy =
1220 gasPvt_->internalEnergy(regionIdx, T, p,
1221 BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1222 BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1224 const LhsEval& Rsw = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState,
waterPhaseIdx, regionIdx);
1225 const LhsEval& bw = waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
1230 const LhsEval Rsw(0.0);
1233 * waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
1235 throw std::logic_error(
"Unhandled phase index " + std::to_string(phaseIdx));
1241 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1242 STATIC_OR_DEVICE LhsEval
enthalpy(
const FluidState& fluidState,
1247 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1248 auto energy = internalEnergy<FluidState, LhsEval>(fluidState, phaseIdx, regionIdx);
1249 if(!enthalpy_eq_energy_){
1251 energy += p/density<FluidState, LhsEval>(fluidState, phaseIdx, regionIdx);
1262 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1270 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1271 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1272 const auto& saltConcentration = decay<LhsEval>(fluidState.saltConcentration());
1276 case gasPhaseIdx:
return gasPvt_->saturatedWaterVaporizationFactor(regionIdx, T, p, saltConcentration);
1278 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1288 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1292 const LhsEval& maxOilSaturation)
1298 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1299 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1303 case oilPhaseIdx:
return oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p, So, maxOilSaturation);
1304 case gasPhaseIdx:
return gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p, So, maxOilSaturation);
1305 case waterPhaseIdx:
return waterPvt_->saturatedGasDissolutionFactor(regionIdx, T, p,
1306 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1307 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1319 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1328 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1329 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1332 case oilPhaseIdx:
return oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p);
1333 case gasPhaseIdx:
return gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p);
1334 case waterPhaseIdx:
return waterPvt_->saturatedGasDissolutionFactor(regionIdx, T, p,
1335 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1336 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1343 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1354 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1371 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1379 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1382 case oilPhaseIdx:
return oilPvt_->saturationPressure(regionIdx, T, BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1383 case gasPhaseIdx:
return gasPvt_->saturationPressure(regionIdx, T, BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1384 case waterPhaseIdx:
return waterPvt_->saturationPressure(regionIdx, T,
1385 BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1386 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1387 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1398 template <
class LhsEval>
1404 return XoG/(1.0 - XoG)*(rho_oRef/rho_gRef);
1411 template <
class LhsEval>
1417 return XwG/(1.0 - XwG)*(rho_wRef/rho_gRef);
1424 template <
class LhsEval>
1430 return XgO/(1.0 - XgO)*(rho_gRef/rho_oRef);
1437 template <
class LhsEval>
1443 return XgW/(1.0 - XgW)*(rho_gRef/rho_wRef);
1451 template <
class LhsEval>
1457 const LhsEval& rho_oG = Rs*rho_gRef;
1458 return rho_oG/(rho_oRef + rho_oG);
1465 template <
class LhsEval>
1471 const LhsEval& rho_wG = Rsw*rho_gRef;
1472 return rho_wG/(rho_wRef + rho_wG);
1479 template <
class LhsEval>
1485 const LhsEval& rho_gO = Rv*rho_oRef;
1486 return rho_gO/(rho_gRef + rho_gO);
1493 template <
class LhsEval>
1499 const LhsEval& rho_gW = Rvw*rho_wRef;
1500 return rho_gW/(rho_gRef + rho_gW);
1506 template <
class LhsEval>
1512 return XgW*MG / (MW*(1 - XgW) + XgW*MG);
1518 template <
class LhsEval>
1524 return XwG*MW / (MG*(1 - XwG) + XwG*MW);
1530 template <
class LhsEval>
1536 return XoG*MO / (MG*(1 - XoG) + XoG*MO);
1542 template <
class LhsEval>
1548 return xoG*MG / (xoG*(MG - MO) + MO);
1554 template <
class LhsEval>
1560 return XgO*MG / (MO*(1 - XgO) + XgO*MG);
1566 template <
class LhsEval>
1572 return xgO*MO / (xgO*(MO - MG) + MG);
1583 {
return *gasPvt_; }
1593 {
return *oilPvt_; }
1603 {
return *waterPvt_; }
1611 {
return reservoirTemperature_; }
1619 { reservoirTemperature_ = value; }
1621 STATIC_OR_DEVICE
short activeToCanonicalPhaseIdx(
unsigned activePhaseIdx);
1623 STATIC_OR_DEVICE
short canonicalToActivePhaseIdx(
unsigned phaseIdx);
1627 {
return diffusionCoefficients_[regionIdx][
numPhases*compIdx + phaseIdx]; }
1631 { diffusionCoefficients_[regionIdx][
numPhases*compIdx + phaseIdx] = coefficient ; }
1636 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
1638 const ParameterCache<ParamCacheEval>& paramCache,
1647 if(!diffusionCoefficients_.empty()) {
1651 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1652 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1658 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1661 STATIC_OR_DEVICE
void setEnergyEqualEnthalpy(
bool enthalpy_eq_energy){
1662 enthalpy_eq_energy_ = enthalpy_eq_energy;
1665 STATIC_OR_DEVICE
bool enthalpyEqualEnergy(){
1666 return enthalpy_eq_energy_;
1670 STATIC_OR_NOTHING
void resizeArrays_(std::size_t
numRegions);
1672 STATIC_OR_NOTHING
Scalar reservoirTemperature_;
1674 STATIC_OR_NOTHING SmartPointer<GasPvt> gasPvt_;
1675 STATIC_OR_NOTHING SmartPointer<OilPvt> oilPvt_;
1676 STATIC_OR_NOTHING SmartPointer<WaterPvt> waterPvt_;
1678 STATIC_OR_NOTHING
bool enableDissolvedGas_;
1679 STATIC_OR_NOTHING
bool enableDissolvedGasInWater_;
1680 STATIC_OR_NOTHING
bool enableVaporizedOil_;
1681 STATIC_OR_NOTHING
bool enableVaporizedWater_;
1682 STATIC_OR_NOTHING
bool enableDiffusion_;
1687 STATIC_OR_NOTHING Storage<std::array<
Scalar, 3> > referenceDensity_;
1688 STATIC_OR_NOTHING Storage<std::array<
Scalar, 3> > molarMass_;
1689 STATIC_OR_NOTHING Storage<std::array<
Scalar, 3 * 3> > diffusionCoefficients_;
1691 STATIC_OR_NOTHING std::array<short, numPhases> activeToCanonicalPhaseIdx_;
1692 STATIC_OR_NOTHING std::array<short, numPhases> canonicalToActivePhaseIdx_;
1694 STATIC_OR_NOTHING
bool isInitialized_;
1695 STATIC_OR_NOTHING
bool useSaturatedTables_;
1696 STATIC_OR_NOTHING
bool enthalpy_eq_energy_;
1698 #ifndef COMPILING_STATIC_FLUID_SYSTEM
1699 template<
template<
typename>
typename StorageT,
template<
typename>
typename SmartPointerT>
1700 explicit FLUIDSYSTEM_CLASSNAME(
const FLUIDSYSTEM_CLASSNAME_STATIC<Scalar, IndexTraits, StorageT, SmartPointerT>& other)
1703 , numActivePhases_(other.numActivePhases_)
1704 , phaseIsActive_(other.phaseIsActive_)
1705 , reservoirTemperature_(other.reservoirTemperature_)
1706 , gasPvt_(SmartPointerT<typename decltype(gasPvt_)::element_type>(other.gasPvt_))
1707 , oilPvt_(SmartPointerT<typename decltype(oilPvt_)::element_type>(other.oilPvt_))
1708 , waterPvt_(SmartPointerT<typename decltype(waterPvt_)::element_type>(other.waterPvt_))
1709 , enableDissolvedGas_(other.enableDissolvedGas_)
1710 , enableDissolvedGasInWater_(other.enableDissolvedGasInWater_)
1711 , enableVaporizedOil_(other.enableVaporizedOil_)
1712 , enableVaporizedWater_(other.enableVaporizedWater_)
1713 , enableDiffusion_(other.enableDiffusion_)
1714 , referenceDensity_(StorageT<typename decltype(referenceDensity_)::value_type>(other.referenceDensity_))
1715 , molarMass_(StorageT<typename decltype(molarMass_)::value_type>(other.molarMass_))
1716 , diffusionCoefficients_(StorageT<typename decltype(diffusionCoefficients_)::value_type>(other.diffusionCoefficients_))
1717 , activeToCanonicalPhaseIdx_(other.activeToCanonicalPhaseIdx_)
1718 , canonicalToActivePhaseIdx_(other.canonicalToActivePhaseIdx_)
1719 , isInitialized_(other.isInitialized_)
1720 , useSaturatedTables_(other.useSaturatedTables_)
1721 , enthalpy_eq_energy_(other.enthalpy_eq_energy_)
1723 OPM_ERROR_IF(!other.isInitialized(),
"The fluid system must be initialized before it can be copied.");
1726 template<
class ScalarT,
class IndexTraitsT,
template<
typename>
typename StorageT,
template<
typename>
typename SmartPointerT>
1727 friend class FLUIDSYSTEM_CLASSNAME_STATIC;
1729 template<
class ScalarT,
class IndexTraitsT,
template<
typename>
typename StorageT,
template<
typename>
typename SmartPointerT>
1730 friend class FLUIDSYSTEM_CLASSNAME_NONSTATIC;
1734template <
class Scalar,
class IndexTraits,
template<
typename>
typename Storage,
template<
typename>
typename SmartPointer>
1738 isInitialized_ =
false;
1739 useSaturatedTables_ =
true;
1741 enableDissolvedGas_ =
true;
1742 enableDissolvedGasInWater_ =
false;
1743 enableVaporizedOil_ =
false;
1744 enableVaporizedWater_ =
false;
1745 enableDiffusion_ =
false;
1749 waterPvt_ =
nullptr;
1751 surfaceTemperature = 273.15 + 15.56;
1752 surfacePressure = 1.01325e5;
1753 setReservoirTemperature(surfaceTemperature);
1755 numActivePhases_ = numPhases;
1756 std::fill_n(&phaseIsActive_[0], numPhases,
true);
1758 resizeArrays_(numPvtRegions);
1761template <
class Scalar,
class IndexTraits,
template<
typename>
typename Storage,
template<
typename>
typename SmartPointer>
1768 referenceDensity_[regionIdx][oilPhaseIdx] = rhoOil;
1769 referenceDensity_[regionIdx][waterPhaseIdx] = rhoWater;
1770 referenceDensity_[regionIdx][gasPhaseIdx] = rhoGas;
1773template <
class Scalar,
class IndexTraits,
template<
typename>
typename Storage,
template<
typename>
typename SmartPointer>
1777 const std::size_t num_regions = molarMass_.size();
1778 for (
unsigned regionIdx = 0; regionIdx < num_regions; ++regionIdx) {
1782 molarMass_[regionIdx][waterCompIdx] = 18e-3;
1784 if (phaseIsActive(gasPhaseIdx)) {
1786 Scalar p = surfacePressure;
1787 Scalar T = surfaceTemperature;
1788 Scalar rho_g = referenceDensity_[0][gasPhaseIdx];
1793 molarMass_[regionIdx][gasCompIdx] = 2e-3;
1796 molarMass_[regionIdx][oilCompIdx] = 175e-3;
1800 int activePhaseIdx = 0;
1801 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
1802 if(phaseIsActive(phaseIdx)){
1803 canonicalToActivePhaseIdx_[phaseIdx] = activePhaseIdx;
1804 activeToCanonicalPhaseIdx_[activePhaseIdx] = phaseIdx;
1808 isInitialized_ =
true;
1811template <
class Scalar,
class IndexTraits,
template<
typename>
typename Storage,
template<
typename>
typename SmartPointer>
1824 throw std::logic_error(std::string(
"Phase index ") + std::to_string(phaseIdx) +
" is unknown");
1828template <
class Scalar,
class IndexTraits,
template<
typename>
typename Storage,
template<
typename>
typename SmartPointer>
1834 return waterCompIdx;
1841 throw std::logic_error(std::string(
"Phase index ") + std::to_string(phaseIdx) +
" is unknown");
1845template <
class Scalar,
class IndexTraits,
template<
typename>
typename Storage,
template<
typename>
typename SmartPointer>
1851 if (enableDissolvedGasInWater())
1853 throw std::logic_error(
"The water phase does not have any solutes in the black oil model!");
1857 if (enableVaporizedWater()) {
1858 return waterCompIdx;
1863 throw std::logic_error(std::string(
"Phase index ") + std::to_string(phaseIdx) +
" is unknown");
1867template <
class Scalar,
class IndexTraits,
template<
typename>
typename Storage,
template<
typename>
typename SmartPointer>
1880 throw std::logic_error(std::string(
"Component index ") + std::to_string(compIdx) +
" is unknown");
1884template <
class Scalar,
class IndexTraits,
template<
typename>
typename Storage,
template<
typename>
typename SmartPointer>
1888 assert(activePhaseIdx<numActivePhases());
1889 return activeToCanonicalPhaseIdx_[activePhaseIdx];
1892template <
class Scalar,
class IndexTraits,
template<
typename>
typename Storage,
template<
typename>
typename SmartPointer>
1893short FLUIDSYSTEM_CLASSNAME<Scalar,IndexTraits, Storage, SmartPointer>::
1894canonicalToActivePhaseIdx(
unsigned phaseIdx)
1896 assert(phaseIdx<numPhases);
1897 assert(phaseIsActive(phaseIdx));
1898 return canonicalToActivePhaseIdx_[phaseIdx];
1901template <
class Scalar,
class IndexTraits,
template<
typename>
typename Storage,
template<
typename>
typename SmartPointer>
1902void FLUIDSYSTEM_CLASSNAME<Scalar,IndexTraits,Storage, SmartPointer>::
1903resizeArrays_(std::size_t numRegions)
1905 molarMass_.resize(numRegions);
1906 referenceDensity_.resize(numRegions);
1909#ifdef COMPILING_STATIC_FLUID_SYSTEM
1910template <
typename T>
using BOFS = FLUIDSYSTEM_CLASSNAME<T, BlackOilDefaultIndexTraits, VectorWithDefaultAllocator, std::shared_ptr>;
1912#define DECLARE_INSTANCE(T) \
1913template<> unsigned char BOFS<T>::numActivePhases_; \
1914template<> std::array<bool, BOFS<T>::numPhases> BOFS<T>::phaseIsActive_; \
1915template<> std::array<short, BOFS<T>::numPhases> BOFS<T>::activeToCanonicalPhaseIdx_; \
1916template<> std::array<short, BOFS<T>::numPhases> BOFS<T>::canonicalToActivePhaseIdx_; \
1917template<> T BOFS<T>::surfaceTemperature; \
1918template<> T BOFS<T>::surfacePressure; \
1919template<> T BOFS<T>::reservoirTemperature_; \
1920template<> bool BOFS<T>::enableDissolvedGas_; \
1921template<> bool BOFS<T>::enableDissolvedGasInWater_; \
1922template<> bool BOFS<T>::enableVaporizedOil_; \
1923template<> bool BOFS<T>::enableVaporizedWater_; \
1924template<> bool BOFS<T>::enableDiffusion_; \
1925template<> std::shared_ptr<OilPvtMultiplexer<T>> BOFS<T>::oilPvt_; \
1926template<> std::shared_ptr<GasPvtMultiplexer<T>> BOFS<T>::gasPvt_; \
1927template<> std::shared_ptr<WaterPvtMultiplexer<T>> BOFS<T>::waterPvt_; \
1928template<> std::vector<std::array<T, 3>> BOFS<T>::referenceDensity_; \
1929template<> std::vector<std::array<T, 3>> BOFS<T>::molarMass_; \
1930template<> std::vector<std::array<T, 9>> BOFS<T>::diffusionCoefficients_; \
1931template<> bool BOFS<T>::isInitialized_; \
1932template<> bool BOFS<T>::useSaturatedTables_; \
1933template<> bool BOFS<T>::enthalpy_eq_energy_;
1935DECLARE_INSTANCE(
float)
1936DECLARE_INSTANCE(
double)
1938#undef DECLARE_INSTANCE
The base class for all fluid systems.
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.
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
Definition Constants.hpp:40
Definition BlackOilFluidSystem.hpp:68
Definition BlackOilFluidSystemNonStatic.hpp:70
A fluid system which uses the black-oil model assumptions to calculate termodynamically meaningful qu...
Definition BlackOilFluidSystem_macrotemplate.hpp:73
STATIC_OR_DEVICE 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_macrotemplate.hpp:1242
STATIC_OR_DEVICE void setEnableVaporizedOil(bool yesno)
Specify whether the fluid system should consider that the oil component can dissolve in the gas phase...
Definition BlackOilFluidSystem_macrotemplate.hpp:236
STATIC_OR_DEVICE LhsEval convertxgOToXgO(const LhsEval &xgO, unsigned regionIdx)
Convert a oil mole fraction in the gas phase the corresponding mass fraction.
Definition BlackOilFluidSystem_macrotemplate.hpp:1567
STATIC_OR_DEVICE bool isCompressible(unsigned)
Returns true if and only if a fluid phase is assumed to be compressible.
Definition BlackOilFluidSystem_macrotemplate.hpp:405
STATIC_OR_DEVICE LhsEval bubblePointPressure(const FluidState &fluidState, unsigned regionIdx)
Returns the bubble point pressure $P_b$ using the current Rs.
Definition BlackOilFluidSystem_macrotemplate.hpp:1344
STATIC_OR_DEVICE std::string_view componentName(unsigned compIdx)
Return the human readable name of a component.
Definition BlackOilFluidSystem_macrotemplate.hpp:1869
STATIC_OR_DEVICE unsigned numActivePhases()
Returns the number of active fluid phases (i.e., usually three)
Definition BlackOilFluidSystem_macrotemplate.hpp:373
static constexpr int waterCompIdx
Index of the water component.
Definition BlackOilFluidSystem_macrotemplate.hpp:363
STATIC_OR_DEVICE 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_macrotemplate.hpp:1480
STATIC_OR_DEVICE unsigned solventComponentIndex(unsigned phaseIdx)
returns the index of "primary" component of a phase (solvent)
Definition BlackOilFluidSystem_macrotemplate.hpp:1830
STATIC_OR_DEVICE 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_macrotemplate.hpp:517
STATIC_OR_DEVICE LhsEval density(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Calculate the density [kg/m^3] of a fluid phase.
Definition BlackOilFluidSystem_macrotemplate.hpp:534
STATIC_OR_DEVICE LhsEval convertXoGToxoG(const LhsEval &XoG, unsigned regionIdx)
Convert a gas mass fraction in the oil phase the corresponding mole fraction.
Definition BlackOilFluidSystem_macrotemplate.hpp:1531
STATIC_OR_DEVICE bool enableVaporizedOil()
Returns whether the fluid system should consider that the oil component can dissolve in the gas phase...
Definition BlackOilFluidSystem_macrotemplate.hpp:449
STATIC_OR_DEVICE bool phaseIsActive(unsigned phaseIdx)
Returns whether a fluid phase is active.
Definition BlackOilFluidSystem_macrotemplate.hpp:377
STATIC_OR_DEVICE Scalar reservoirTemperature(unsigned=0)
Set the temperature of the reservoir.
Definition BlackOilFluidSystem_macrotemplate.hpp:1610
static constexpr unsigned gasPhaseIdx
Index of the gas phase.
Definition BlackOilFluidSystem_macrotemplate.hpp:335
STATIC_OR_DEVICE bool enableDissolvedGas()
Returns whether the fluid system should consider that the gas component can dissolve in the oil phase...
Definition BlackOilFluidSystem_macrotemplate.hpp:430
STATIC_OR_DEVICE LhsEval saturationPressure(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the saturation pressure of a given phase [Pa] depending on its composition.
Definition BlackOilFluidSystem_macrotemplate.hpp:1372
STATIC_OR_DEVICE 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_macrotemplate.hpp:1637
STATIC_OR_DEVICE void setEnableVaporizedWater(bool yesno)
Specify whether the fluid system should consider that the water component can dissolve in the gas pha...
Definition BlackOilFluidSystem_macrotemplate.hpp:245
STATIC_OR_DEVICE LhsEval convertXgOToxgO(const LhsEval &XgO, unsigned regionIdx)
Convert a oil mass fraction in the gas phase the corresponding mole fraction.
Definition BlackOilFluidSystem_macrotemplate.hpp:1555
STATIC_OR_NOTHING Scalar surfacePressure
The pressure at the surface.
Definition BlackOilFluidSystem_macrotemplate.hpp:338
STATIC_OR_DEVICE bool useSaturatedTables()
Returns whether the saturated tables should be used.
Definition BlackOilFluidSystem_macrotemplate.hpp:474
STATIC_OR_DEVICE LhsEval viscosity(const FluidState &fluidState, const ParameterCache< ParamCacheEval > ¶mCache, unsigned phaseIdx)
Calculate the dynamic viscosity of a fluid phase [Pa*s].
Definition BlackOilFluidSystem_macrotemplate.hpp:510
STATIC_OR_DEVICE bool isLiquid(unsigned phaseIdx)
Return whether a phase is liquid.
Definition BlackOilFluidSystem_macrotemplate.hpp:347
STATIC_OR_DEVICE 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_macrotemplate.hpp:1452
STATIC_OR_DEVICE LhsEval saturatedInverseFormationVolumeFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the formation volume factor of a "saturated" fluid phase.
Definition BlackOilFluidSystem_macrotemplate.hpp:833
STATIC_OR_DEVICE LhsEval convertXgWToxgW(const LhsEval &XgW, unsigned regionIdx)
Convert a water mass fraction in the gas phase the corresponding mole fraction.
Definition BlackOilFluidSystem_macrotemplate.hpp:1507
STATIC_OR_DEVICE std::string_view phaseName(unsigned phaseIdx)
Return the human readable name of a fluid phase.
Definition BlackOilFluidSystem_macrotemplate.hpp:1813
STATIC_OR_DEVICE 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_macrotemplate.hpp:855
STATIC_OR_DEVICE bool isIdealMixture(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal mixture.
Definition BlackOilFluidSystem_macrotemplate.hpp:397
STATIC_OR_DEVICE LhsEval convertXwGToxwG(const LhsEval &XwG, unsigned regionIdx)
Convert a gas mass fraction in the water phase the corresponding mole fraction.
Definition BlackOilFluidSystem_macrotemplate.hpp:1519
STATIC_OR_DEVICE void setDiffusionCoefficient(Scalar coefficient, unsigned compIdx, unsigned phaseIdx, unsigned regionIdx=0)
Definition BlackOilFluidSystem_macrotemplate.hpp:1630
STATIC_OR_DEVICE 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_macrotemplate.hpp:1412
static constexpr int gasCompIdx
Index of the gas component.
Definition BlackOilFluidSystem_macrotemplate.hpp:365
STATIC_OR_DEVICE void setUseSaturatedTables(bool yesno)
Specify whether the saturated tables should be used.
Definition BlackOilFluidSystem_macrotemplate.hpp:269
STATIC_OR_DEVICE LhsEval dewPointPressure(const FluidState &fluidState, unsigned regionIdx)
Returns the dew point pressure $P_d$ using the current Rv.
Definition BlackOilFluidSystem_macrotemplate.hpp:1355
static constexpr unsigned numComponents
Number of chemical species in the fluid system.
Definition BlackOilFluidSystem_macrotemplate.hpp:358
STATIC_OR_NOTHING Scalar surfaceTemperature
The temperature at the surface.
Definition BlackOilFluidSystem_macrotemplate.hpp:341
STATIC_OR_DEVICE 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_macrotemplate.hpp:1466
STATIC_OR_DEVICE unsigned soluteComponentIndex(unsigned phaseIdx)
returns the index of "secondary" component of a phase (solute)
Definition BlackOilFluidSystem_macrotemplate.hpp:1847
STATIC_OR_DEVICE 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_macrotemplate.hpp:1626
STATIC_OR_DEVICE bool enableVaporizedWater()
Returns whether the fluid system should consider that the water component can dissolve in the gas pha...
Definition BlackOilFluidSystem_macrotemplate.hpp:458
STATIC_OR_DEVICE Scalar referenceDensity(unsigned phaseIdx, unsigned regionIdx)
Returns the density of a fluid phase at surface pressure [kg/m^3].
Definition BlackOilFluidSystem_macrotemplate.hpp:482
STATIC_OR_DEVICE void initEnd()
Finish initializing the black oil fluid system.
Definition BlackOilFluidSystem_macrotemplate.hpp:1774
STATIC_OR_DEVICE void setEnableDiffusion(bool yesno)
Specify whether the fluid system should consider diffusion.
Definition BlackOilFluidSystem_macrotemplate.hpp:261
STATIC_OR_DEVICE LhsEval viscosity(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Calculate the dynamic viscosity of a fluid phase [Pa*s].
Definition BlackOilFluidSystem_macrotemplate.hpp:977
STATIC_OR_DEVICE LhsEval convertxoGToXoG(const LhsEval &xoG, unsigned regionIdx)
Convert a gas mole fraction in the oil phase the corresponding mass fraction.
Definition BlackOilFluidSystem_macrotemplate.hpp:1543
STATIC_OR_DEVICE bool isIdealGas(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal gas.
Definition BlackOilFluidSystem_macrotemplate.hpp:409
STATIC_OR_DEVICE LhsEval saturatedVaporizationFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the water vaporization factor of saturated phase.
Definition BlackOilFluidSystem_macrotemplate.hpp:1263
STATIC_OR_DEVICE std::size_t numRegions()
Returns the number of PVT regions which are considered.
Definition BlackOilFluidSystem_macrotemplate.hpp:421
static constexpr unsigned numPhases
Number of fluid phases in the fluid system.
Definition BlackOilFluidSystem_macrotemplate.hpp:328
STATIC_OR_DEVICE LhsEval density(const FluidState &fluidState, const ParameterCache< ParamCacheEval > ¶mCache, unsigned phaseIdx)
Calculate the density [kg/m^3] of a fluid phase.
Definition BlackOilFluidSystem_macrotemplate.hpp:490
static constexpr unsigned waterPhaseIdx
Index of the water phase.
Definition BlackOilFluidSystem_macrotemplate.hpp:331
STATIC_OR_DEVICE bool enableDissolvedGasInWater()
Returns whether the fluid system should consider that the gas component can dissolve in the water pha...
Definition BlackOilFluidSystem_macrotemplate.hpp:440
static constexpr unsigned oilPhaseIdx
Index of the oil phase.
Definition BlackOilFluidSystem_macrotemplate.hpp:333
STATIC_OR_DEVICE 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_macrotemplate.hpp:1399
STATIC_OR_DEVICE LhsEval saturatedDissolutionFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx, const LhsEval &maxOilSaturation)
Returns the dissolution factor of a saturated fluid phase.
Definition BlackOilFluidSystem_macrotemplate.hpp:1289
FLUIDSYSTEM_CLASSNAME(const FLUIDSYSTEM_CLASSNAME< Scalar, IndexTraits, StorageT, SmartPointerT > &other)
Default copy constructor.
Definition BlackOilFluidSystem_macrotemplate.hpp:158
STATIC_OR_DEVICE void setEnableDissolvedGasInWater(bool yesno)
Specify whether the fluid system should consider that the gas component can dissolve in the water pha...
Definition BlackOilFluidSystem_macrotemplate.hpp:254
STATIC_OR_DEVICE 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_macrotemplate.hpp:1425
STATIC_OR_DEVICE void setReferenceDensities(Scalar rhoOil, Scalar rhoWater, Scalar rhoGas, unsigned regionIdx)
Initialize the values of the reference densities.
Definition BlackOilFluidSystem_macrotemplate.hpp:1763
STATIC_OR_DEVICE LhsEval saturatedDensity(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Compute the density of a saturated fluid phase.
Definition BlackOilFluidSystem_macrotemplate.hpp:632
STATIC_OR_DEVICE const GasPvt & gasPvt()
Return a reference to the low-level object which calculates the gas phase quantities.
Definition BlackOilFluidSystem_macrotemplate.hpp:1582
STATIC_OR_DEVICE LhsEval saturatedDissolutionFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the dissolution factor of a saturated fluid phase.
Definition BlackOilFluidSystem_macrotemplate.hpp:1320
STATIC_OR_DEVICE const OilPvt & oilPvt()
Return a reference to the low-level object which calculates the oil phase quantities.
Definition BlackOilFluidSystem_macrotemplate.hpp:1592
STATIC_OR_DEVICE 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_macrotemplate.hpp:1494
STATIC_OR_DEVICE const WaterPvt & waterPvt()
Return a reference to the low-level object which calculates the water phase quantities.
Definition BlackOilFluidSystem_macrotemplate.hpp:1602
STATIC_OR_DEVICE 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_macrotemplate.hpp:497
STATIC_OR_DEVICE void initBegin(std::size_t numPvtRegions)
Begin the initialization of the black oil fluid system.
Definition BlackOilFluidSystem_macrotemplate.hpp:1736
STATIC_OR_DEVICE void setWaterPvt(SmartPointer< WaterPvt > pvtObj)
Set the pressure-volume-saturation (PVT) relations for the water phase.
Definition BlackOilFluidSystem_macrotemplate.hpp:287
STATIC_OR_DEVICE 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_macrotemplate.hpp:1438
static constexpr int oilCompIdx
Index of the oil component.
Definition BlackOilFluidSystem_macrotemplate.hpp:361
STATIC_OR_DEVICE bool enableDiffusion()
Returns whether the fluid system should consider diffusion.
Definition BlackOilFluidSystem_macrotemplate.hpp:466
STATIC_OR_DEVICE void setReservoirTemperature(Scalar value)
Return the temperature of the reservoir.
Definition BlackOilFluidSystem_macrotemplate.hpp:1618
STATIC_OR_DEVICE void setEnableDissolvedGas(bool yesno)
Specify whether the fluid system should consider that the gas component can dissolve in the oil phase...
Definition BlackOilFluidSystem_macrotemplate.hpp:227
STATIC_OR_DEVICE void setGasPvt(SmartPointer< GasPvt > pvtObj)
Set the pressure-volume-saturation (PVT) relations for the gas phase.
Definition BlackOilFluidSystem_macrotemplate.hpp:275
STATIC_OR_DEVICE LhsEval inverseFormationVolumeFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the formation volume factor of an "undersaturated" fluid phase.
Definition BlackOilFluidSystem_macrotemplate.hpp:733
STATIC_OR_DEVICE Scalar molarMass(unsigned compIdx, unsigned regionIdx=0)
Return the molar mass of a component in [kg/mol].
Definition BlackOilFluidSystem_macrotemplate.hpp:393
STATIC_OR_DEVICE void setOilPvt(SmartPointer< OilPvt > pvtObj)
Set the pressure-volume-saturation (PVT) relations for the oil phase.
Definition BlackOilFluidSystem_macrotemplate.hpp:281
This class represents the Pressure-Volume-Temperature relations of the gas phase in the black-oil mod...
Definition GasPvtMultiplexer.hpp:111
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:269
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