27 #ifndef OPM_BLACK_OIL_FLUID_SYSTEM_HPP
28 #define OPM_BLACK_OIL_FLUID_SYSTEM_HPP
45 #include <opm/input/eclipse/EclipseState/EclipseState.hpp>
46 #include <opm/input/eclipse/EclipseState/Tables/FlatTable.hpp>
47 #include <opm/input/eclipse/EclipseState/Tables/TableManager.hpp>
62 template <class FluidSystem, class FluidState, class LhsEval>
63 LhsEval getRs_(typename std::enable_if<!HasMember_Rs<FluidState>::value, const FluidState&>::type fluidState,
67 decay<LhsEval>(fluidState.massFraction(FluidSystem::oilPhaseIdx, FluidSystem::gasCompIdx));
68 return FluidSystem::convertXoGToRs(XoG, regionIdx);
71 template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
72 auto getRs_(
typename std::enable_if<HasMember_Rs<FluidState>::value,
const FluidState&>::type fluidState,
74 -> decltype(decay<LhsEval>(fluidState.Rs()))
75 {
return decay<LhsEval>(fluidState.Rs()); }
77 template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
78 LhsEval getRv_(
typename std::enable_if<!HasMember_Rv<FluidState>::value,
const FluidState&>::type fluidState,
82 decay<LhsEval>(fluidState.massFraction(FluidSystem::gasPhaseIdx, FluidSystem::oilCompIdx));
83 return FluidSystem::convertXgOToRv(XgO, regionIdx);
86 template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
87 auto getRv_(
typename std::enable_if<HasMember_Rv<FluidState>::value,
const FluidState&>::type fluidState,
89 -> decltype(decay<LhsEval>(fluidState.Rv()))
90 {
return decay<LhsEval>(fluidState.Rv()); }
92 template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
93 LhsEval getRvw_(
typename std::enable_if<!HasMember_Rvw<FluidState>::value,
const FluidState&>::type fluidState,
97 decay<LhsEval>(fluidState.massFraction(FluidSystem::gasPhaseIdx, FluidSystem::waterCompIdx));
98 return FluidSystem::convertXgWToRvw(XgW, regionIdx);
101 template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
102 auto getRvw_(
typename std::enable_if<HasMember_Rvw<FluidState>::value,
const FluidState&>::type fluidState,
104 -> decltype(decay<LhsEval>(fluidState.Rvw()))
105 {
return decay<LhsEval>(fluidState.Rvw()); }
107 template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
108 LhsEval getSaltConcentration_(
typename std::enable_if<!HasMember_saltConcentration<FluidState>::value,
109 const FluidState&>::type,
113 template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
114 auto getSaltConcentration_(
typename std::enable_if<HasMember_saltConcentration<FluidState>::value,
const FluidState&>::type fluidState,
116 -> decltype(decay<LhsEval>(fluidState.saltConcentration()))
117 {
return decay<LhsEval>(fluidState.saltConcentration()); }
119 template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
120 LhsEval getSaltSaturation_(
typename std::enable_if<!HasMember_saltSaturation<FluidState>::value,
121 const FluidState&>::type,
125 template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
126 auto getSaltSaturation_(
typename std::enable_if<HasMember_saltSaturation<FluidState>::value,
const FluidState&>::type fluidState,
128 -> decltype(decay<LhsEval>(fluidState.saltSaturation()))
129 {
return decay<LhsEval>(fluidState.saltSaturation()); }
139 template <
class Scalar,
class IndexTraits = BlackOilDefaultIndexTraits>
150 template <
class EvaluationT>
153 using Evaluation = EvaluationT;
158 maxOilSat_ = maxOilSat;
159 regionIdx_ = regionIdx;
169 template <
class OtherCache>
172 regionIdx_ = other.regionIndex();
173 maxOilSat_ = other.maxOilSat();
184 {
return regionIdx_; }
194 { regionIdx_ = val; }
196 const Evaluation& maxOilSat()
const
197 {
return maxOilSat_; }
199 void setMaxOilSat(
const Evaluation& val)
200 { maxOilSat_ = val; }
203 Evaluation maxOilSat_;
214 static void initFromState(
const EclipseState& eclState,
const Schedule& schedule)
216 size_t numRegions = eclState.runspec().tabdims().getNumPVTTables();
219 numActivePhases_ = 0;
220 std::fill_n(&phaseIsActive_[0],
numPhases,
false);
223 if (eclState.runspec().phases().active(Phase::OIL)) {
228 if (eclState.runspec().phases().active(Phase::GAS)) {
233 if (eclState.runspec().phases().active(Phase::WATER)) {
247 assert(numActivePhases_ >= 1 && numActivePhases_ <= 3);
254 gasPvt_ = std::make_shared<GasPvt>();
255 gasPvt_->initFromState(eclState, schedule);
259 oilPvt_ = std::make_shared<OilPvt>();
260 oilPvt_->initFromState(eclState, schedule);
264 waterPvt_ = std::make_shared<WaterPvt>();
265 waterPvt_->initFromState(eclState, schedule);
269 for (
unsigned regionIdx = 0; regionIdx <
numRegions; ++regionIdx) {
283 if (eclState.runspec().co2Storage()) {
284 for (
unsigned regionIdx = 0; regionIdx <
numRegions; ++regionIdx) {
292 const auto& diffCoeffTables = eclState.getTableManager().getDiffusionCoefficientTable();
293 if(!diffCoeffTables.empty()) {
296 diffusionCoefficients_.resize(
numRegions,{0,0,0,0,0,0,0,0,0});
298 for (
unsigned regionIdx = 0; regionIdx <
numRegions; ++regionIdx) {
299 const auto& diffCoeffTable = diffCoeffTables[regionIdx];
300 molarMass_[regionIdx][
oilCompIdx] = diffCoeffTable.oil_mw;
301 molarMass_[regionIdx][
gasCompIdx] = diffCoeffTable.gas_mw;
306 if(diffCoeffTable.gas_in_oil_cross_phase > 0 || diffCoeffTable.oil_in_oil_cross_phase > 0) {
307 throw std::runtime_error(
"Cross phase diffusion is set in the deck, but not implemented in Flow. "
308 "Please default DIFFC item 7 and item 8 or set it to zero.");
326 isInitialized_ =
false;
328 enableDissolvedGas_ =
true;
329 enableVaporizedOil_ =
false;
330 enableVaporizedWater_ =
false;
331 enableDiffusion_ =
false;
342 std::fill_n(&phaseIsActive_[0],
numPhases,
true);
344 resizeArrays_(numPvtRegions);
354 { enableDissolvedGas_ = yesno; }
363 { enableVaporizedOil_ = yesno; }
372 { enableVaporizedWater_ = yesno; }
380 { enableDiffusion_ = yesno; }
387 { gasPvt_ = pvtObj; }
393 { oilPvt_ = pvtObj; }
399 { waterPvt_ = pvtObj; }
413 referenceDensity_[regionIdx][
oilPhaseIdx] = rhoOil;
415 referenceDensity_[regionIdx][
gasPhaseIdx] = rhoGas;
426 for (
unsigned regionIdx = 0; regionIdx <
numRegions; ++ regionIdx) {
448 int activePhaseIdx = 0;
449 for (
unsigned phaseIdx = 0; phaseIdx <
numPhases; ++phaseIdx) {
451 canonicalToActivePhaseIdx_[phaseIdx] = activePhaseIdx;
452 activeToCanonicalPhaseIdx_[activePhaseIdx] = phaseIdx;
456 isInitialized_ =
true;
459 static bool isInitialized()
460 {
return isInitialized_; }
494 throw std::logic_error(
"Phase index " + std::to_string(phaseIdx) +
" is unknown");
513 static constexpr
unsigned oilCompIdx = IndexTraits::oilCompIdx;
517 static constexpr
unsigned gasCompIdx = IndexTraits::gasCompIdx;
520 static unsigned char numActivePhases_;
521 static std::array<bool,numPhases> phaseIsActive_;
526 {
return numActivePhases_; }
532 return phaseIsActive_[phaseIdx];
547 throw std::logic_error(
"Phase index " + std::to_string(phaseIdx) +
" is unknown");
556 throw std::logic_error(
"The water phase does not have any solutes in the black oil model!");
563 throw std::logic_error(
"Phase index " + std::to_string(phaseIdx) +
" is unknown");
579 throw std::logic_error(
"Component index " + std::to_string(compIdx) +
" is unknown");
585 {
return molarMass_[regionIdx][compIdx]; }
613 {
return molarMass_.size(); }
622 {
return enableDissolvedGas_; }
631 {
return enableVaporizedOil_; }
640 {
return enableVaporizedWater_; }
648 {
return enableDiffusion_; }
656 {
return referenceDensity_[regionIdx][phaseIdx]; }
662 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
663 static LhsEval
density(
const FluidState& fluidState,
666 {
return density<FluidState, LhsEval>(fluidState, phaseIdx, paramCache.
regionIndex()); }
669 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
675 return fugacityCoefficient<FluidState, LhsEval>(fluidState,
682 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
686 {
return viscosity<FluidState, LhsEval>(fluidState, phaseIdx, paramCache.
regionIndex()); }
689 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
690 static LhsEval
enthalpy(
const FluidState& fluidState,
693 {
return enthalpy<FluidState, LhsEval>(fluidState, phaseIdx, paramCache.
regionIndex()); }
700 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
701 static LhsEval
density(
const FluidState& fluidState,
708 const LhsEval& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
709 const LhsEval& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
710 const LhsEval& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
716 const LhsEval& Rs = BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
717 const LhsEval& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
725 const LhsEval Rs(0.0);
726 const auto& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
734 const LhsEval& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
735 const LhsEval& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
736 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
745 const LhsEval Rvw(0.0);
746 const LhsEval& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
747 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
755 const LhsEval Rv(0.0);
756 const LhsEval& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
757 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
765 const LhsEval Rv(0.0);
766 const LhsEval Rvw(0.0);
767 const auto& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
774 * waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, saltConcentration);
777 throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
787 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
795 const auto& p = fluidState.pressure(phaseIdx);
796 const auto& T = fluidState.temperature(phaseIdx);
802 const LhsEval& Rs = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState,
oilPhaseIdx, regionIdx);
803 const LhsEval& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
811 const LhsEval Rs(0.0);
812 const LhsEval& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
819 const LhsEval& Rv = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState,
gasPhaseIdx, regionIdx);
820 const LhsEval& Rvw = saturatedVaporizationFactor<FluidState, LhsEval>(fluidState,
gasPhaseIdx, regionIdx);
821 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
831 const LhsEval Rvw(0.0);
832 const LhsEval& Rv = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState,
gasPhaseIdx, regionIdx);
833 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
842 const LhsEval Rv(0.0);
843 const LhsEval& Rvw = saturatedVaporizationFactor<FluidState, LhsEval>(fluidState,
gasPhaseIdx, regionIdx);
844 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
852 const LhsEval Rv(0.0);
853 const LhsEval Rvw(0.0);
854 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
863 *inverseFormationVolumeFactor<FluidState, LhsEval>(fluidState,
waterPhaseIdx, regionIdx);
866 throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
877 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
885 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
886 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
887 const auto& saltConcentration = decay<LhsEval>(fluidState.saltConcentration());
892 const auto& Rs = BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
894 && Rs >= (1.0 - 1e-10)*oilPvt_->saturatedGasDissolutionFactor(regionIdx, scalarValue(T), scalarValue(p)))
896 return oilPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
898 return oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
902 const LhsEval Rs(0.0);
903 return oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
907 const auto& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
908 const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
910 && Rvw >= (1.0 - 1e-10)*gasPvt_->saturatedWaterVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p))
912 && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
914 return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
916 return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
921 const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
923 && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
925 return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
927 const LhsEval Rvw(0.0);
928 return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
933 const auto& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
935 && Rvw >= (1.0 - 1e-10)*gasPvt_->saturatedWaterVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
937 return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
939 const LhsEval Rv(0.0);
940 return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
944 const LhsEval Rv(0.0);
945 const LhsEval Rvw(0.0);
946 return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
949 return waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, saltConcentration);
950 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
961 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
969 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
970 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
971 const auto& saltConcentration = decay<LhsEval>(fluidState.saltConcentration());
974 case oilPhaseIdx:
return oilPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
975 case gasPhaseIdx:
return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
976 case waterPhaseIdx:
return waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, saltConcentration);
977 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
982 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
992 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
993 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
998 const LhsEval phi_oO = 20e3/p;
1001 constexpr
const Scalar phi_gG = 1.0;
1005 const LhsEval phi_wW = 30e3/p;
1023 const auto& R_vSat = gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p);
1027 const auto& R_sSat = oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p);
1030 const auto& x_oOSat = 1.0 - x_oGSat;
1032 const auto& p_o = decay<LhsEval>(fluidState.pressure(
oilPhaseIdx));
1033 const auto& p_g = decay<LhsEval>(fluidState.pressure(
gasPhaseIdx));
1035 return phi_oO*p_o*x_oOSat / (p_g*x_gOSat);
1043 throw std::logic_error(
"Invalid component index "+std::to_string(compIdx));
1059 const auto& R_vSat = gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p);
1062 const auto& x_gGSat = 1.0 - x_gOSat;
1064 const auto& R_sSat = oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p);
1068 const auto& p_o = decay<LhsEval>(fluidState.pressure(
oilPhaseIdx));
1069 const auto& p_g = decay<LhsEval>(fluidState.pressure(
gasPhaseIdx));
1071 return phi_gG*p_g*x_gGSat / (p_o*x_oGSat);
1078 throw std::logic_error(
"Invalid component index "+std::to_string(compIdx));
1093 throw std::logic_error(
"Invalid component index "+std::to_string(compIdx));
1097 throw std::logic_error(
"Invalid phase index "+std::to_string(phaseIdx));
1100 throw std::logic_error(
"Unhandled phase or component index");
1104 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1112 const LhsEval& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1113 const LhsEval& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1114 const LhsEval& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1119 const auto& Rs = BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1121 && Rs >= (1.0 - 1e-10)*oilPvt_->saturatedGasDissolutionFactor(regionIdx, scalarValue(T), scalarValue(p)))
1123 return oilPvt_->saturatedViscosity(regionIdx, T, p);
1125 return oilPvt_->viscosity(regionIdx, T, p, Rs);
1129 const LhsEval Rs(0.0);
1130 return oilPvt_->viscosity(regionIdx, T, p, Rs);
1135 const auto& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1136 const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1138 && Rvw >= (1.0 - 1e-10)*gasPvt_->saturatedWaterVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p))
1140 && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
1142 return gasPvt_->saturatedViscosity(regionIdx, T, p);
1144 return gasPvt_->viscosity(regionIdx, T, p, Rv, Rvw);
1148 const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1150 && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
1152 return gasPvt_->saturatedViscosity(regionIdx, T, p);
1154 const LhsEval Rvw(0.0);
1155 return gasPvt_->viscosity(regionIdx, T, p, Rv, Rvw);
1159 const auto& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1161 && Rvw >= (1.0 - 1e-10)*gasPvt_->saturatedWaterVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
1163 return gasPvt_->saturatedViscosity(regionIdx, T, p);
1165 const LhsEval Rv(0.0);
1166 return gasPvt_->viscosity(regionIdx, T, p, Rv, Rvw);
1170 const LhsEval Rv(0.0);
1171 const LhsEval Rvw(0.0);
1172 return gasPvt_->viscosity(regionIdx, T, p, Rv, Rvw);
1178 return waterPvt_->viscosity(regionIdx, T, p, saltConcentration);
1181 throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1185 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1193 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1194 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1199 oilPvt_->internalEnergy(regionIdx, T, p, BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx))
1200 + p/density<FluidState, LhsEval>(fluidState, phaseIdx, regionIdx);
1204 gasPvt_->internalEnergy(regionIdx, T, p, BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx))
1205 + p/density<FluidState, LhsEval>(fluidState, phaseIdx, regionIdx);
1209 waterPvt_->internalEnergy(regionIdx, T, p, BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx))
1210 + p/density<FluidState, LhsEval>(fluidState, phaseIdx, regionIdx);
1212 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1215 throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1224 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1232 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1233 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1237 case gasPhaseIdx:
return gasPvt_->saturatedWaterVaporizationFactor(regionIdx, T, p);
1239 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1249 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1253 const LhsEval& maxOilSaturation)
1258 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1259 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1260 const auto& So = decay<LhsEval>(fluidState.saturation(
oilPhaseIdx));
1263 case oilPhaseIdx:
return oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p, So, maxOilSaturation);
1264 case gasPhaseIdx:
return gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p, So, maxOilSaturation);
1266 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1278 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1286 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1287 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1290 case oilPhaseIdx:
return oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p);
1291 case gasPhaseIdx:
return gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p);
1293 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1300 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1311 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1328 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1336 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1339 case oilPhaseIdx:
return oilPvt_->saturationPressure(regionIdx, T, BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1340 case gasPhaseIdx:
return gasPvt_->saturationPressure(regionIdx, T, BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1342 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1353 template <
class LhsEval>
1359 return XoG/(1.0 - XoG)*(rho_oRef/rho_gRef);
1366 template <
class LhsEval>
1372 return XgO/(1.0 - XgO)*(rho_gRef/rho_oRef);
1379 template <
class LhsEval>
1385 return XgW/(1.0 - XgW)*(rho_gRef/rho_wRef);
1393 template <
class LhsEval>
1399 const LhsEval& rho_oG = Rs*rho_gRef;
1400 return rho_oG/(rho_oRef + rho_oG);
1407 template <
class LhsEval>
1413 const LhsEval& rho_gO = Rv*rho_oRef;
1414 return rho_gO/(rho_gRef + rho_gO);
1421 template <
class LhsEval>
1427 const LhsEval& rho_gW = Rvw*rho_wRef;
1428 return rho_gW/(rho_gRef + rho_gW);
1434 template <
class LhsEval>
1440 return XoG*MO / (MG*(1 - XoG) + XoG*MO);
1446 template <
class LhsEval>
1452 return xoG*MG / (xoG*(MG - MO) + MO);
1458 template <
class LhsEval>
1464 return XgO*MG / (MO*(1 - XgO) + XgO*MG);
1470 template <
class LhsEval>
1476 return xgO*MO / (xgO*(MO - MG) + MG);
1487 {
return *gasPvt_; }
1497 {
return *oilPvt_; }
1507 {
return *waterPvt_; }
1515 {
return reservoirTemperature_; }
1523 { reservoirTemperature_ = value; }
1525 static short activeToCanonicalPhaseIdx(
unsigned activePhaseIdx) {
1527 return activeToCanonicalPhaseIdx_[activePhaseIdx];
1530 static short canonicalToActivePhaseIdx(
unsigned phaseIdx) {
1533 return canonicalToActivePhaseIdx_[phaseIdx];
1538 {
return diffusionCoefficients_[regionIdx][
numPhases*compIdx + phaseIdx]; }
1542 { diffusionCoefficients_[regionIdx][
numPhases*compIdx + phaseIdx] = coefficient ; }
1547 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
1558 if(!diffusionCoefficients_.empty()) {
1562 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1563 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1569 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1580 static Scalar reservoirTemperature_;
1582 static std::shared_ptr<GasPvt> gasPvt_;
1583 static std::shared_ptr<OilPvt> oilPvt_;
1584 static std::shared_ptr<WaterPvt> waterPvt_;
1586 static bool enableDissolvedGas_;
1587 static bool enableVaporizedOil_;
1588 static bool enableVaporizedWater_;
1589 static bool enableDiffusion_;
1594 static std::vector<std::array<
Scalar, 3> > referenceDensity_;
1595 static std::vector<std::array<
Scalar, 3> > molarMass_;
1596 static std::vector<std::array<
Scalar, 3 * 3> > diffusionCoefficients_;
1598 static std::array<short, numPhases> activeToCanonicalPhaseIdx_;
1599 static std::array<short, numPhases> canonicalToActivePhaseIdx_;
1601 static bool isInitialized_;
1604 template <
class Scalar,
class IndexTraits>
1605 unsigned char BlackOilFluidSystem<Scalar, IndexTraits>::numActivePhases_;
1607 template <
class Scalar,
class IndexTraits>
1608 std::array<bool, BlackOilFluidSystem<Scalar, IndexTraits>::numPhases> BlackOilFluidSystem<Scalar, IndexTraits>::phaseIsActive_;
1610 template <
class Scalar,
class IndexTraits>
1611 std::array<short, BlackOilFluidSystem<Scalar, IndexTraits>::numPhases> BlackOilFluidSystem<Scalar, IndexTraits>::activeToCanonicalPhaseIdx_;
1613 template <
class Scalar,
class IndexTraits>
1614 std::array<short, BlackOilFluidSystem<Scalar, IndexTraits>::numPhases> BlackOilFluidSystem<Scalar, IndexTraits>::canonicalToActivePhaseIdx_;
1616 template <
class Scalar,
class IndexTraits>
1620 template <
class Scalar,
class IndexTraits>
1624 template <
class Scalar,
class IndexTraits>
1626 BlackOilFluidSystem<Scalar, IndexTraits>::reservoirTemperature_;
1628 template <
class Scalar,
class IndexTraits>
1629 bool BlackOilFluidSystem<Scalar, IndexTraits>::enableDissolvedGas_;
1631 template <
class Scalar,
class IndexTraits>
1632 bool BlackOilFluidSystem<Scalar, IndexTraits>::enableVaporizedOil_;
1634 template <
class Scalar,
class IndexTraits>
1635 bool BlackOilFluidSystem<Scalar, IndexTraits>::enableVaporizedWater_;
1637 template <
class Scalar,
class IndexTraits>
1638 bool BlackOilFluidSystem<Scalar, IndexTraits>::enableDiffusion_;
1640 template <
class Scalar,
class IndexTraits>
1641 std::shared_ptr<OilPvtMultiplexer<Scalar> >
1642 BlackOilFluidSystem<Scalar, IndexTraits>::oilPvt_;
1644 template <
class Scalar,
class IndexTraits>
1645 std::shared_ptr<GasPvtMultiplexer<Scalar> >
1646 BlackOilFluidSystem<Scalar, IndexTraits>::gasPvt_;
1648 template <
class Scalar,
class IndexTraits>
1649 std::shared_ptr<WaterPvtMultiplexer<Scalar> >
1650 BlackOilFluidSystem<Scalar, IndexTraits>::waterPvt_;
1652 template <
class Scalar,
class IndexTraits>
1653 std::vector<std::array<Scalar, 3> >
1654 BlackOilFluidSystem<Scalar, IndexTraits>::referenceDensity_;
1656 template <
class Scalar,
class IndexTraits>
1657 std::vector<std::array<Scalar, 3> >
1658 BlackOilFluidSystem<Scalar, IndexTraits>::molarMass_;
1660 template <
class Scalar,
class IndexTraits>
1661 std::vector<std::array<Scalar, 9> >
1662 BlackOilFluidSystem<Scalar, IndexTraits>::diffusionCoefficients_;
1664 template <
class Scalar,
class IndexTraits>
1665 bool 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...
A central place for various physical constants occuring in some equations.
Provides the opm-material specific exception classes.
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:44
Scalar Scalar
The type used for scalar quantities.
Definition: BaseFluidSystem.hpp:49
A fluid system which uses the black-oil model assumptions to calculate termodynamically meaningful qu...
Definition: BlackOilFluidSystem.hpp:141
static constexpr unsigned waterCompIdx
Index of the water component.
Definition: BlackOilFluidSystem.hpp:515
static unsigned numActivePhases()
Returns the number of active fluid phases (i.e., usually three)
Definition: BlackOilFluidSystem.hpp:525
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:1537
static void setReferenceDensities(Scalar rhoOil, Scalar rhoWater, Scalar rhoGas, unsigned regionIdx)
Initialize the values of the reference densities.
Definition: BlackOilFluidSystem.hpp:408
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:1354
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:1394
static constexpr unsigned numComponents
Number of chemical species in the fluid system.
Definition: BlackOilFluidSystem.hpp:510
static void initEnd()
Finish initializing the black oil fluid system.
Definition: BlackOilFluidSystem.hpp:422
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:1422
static constexpr unsigned soluteComponentIndex(unsigned phaseIdx)
returns the index of "secondary" component of a phase (solute)
Definition: BlackOilFluidSystem.hpp:552
static LhsEval saturatedInverseFormationVolumeFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the formation volume factor of a "saturated" fluid phase.
Definition: BlackOilFluidSystem.hpp:962
static LhsEval convertXoGToxoG(const LhsEval &XoG, unsigned regionIdx)
Convert a gas mass fraction in the oil phase the corresponding mole fraction.
Definition: BlackOilFluidSystem.hpp:1435
static LhsEval viscosity(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Calculate the dynamic viscosity of a fluid phase [Pa*s].
Definition: BlackOilFluidSystem.hpp:1105
static Scalar reservoirTemperature(unsigned=0)
Set the temperature of the reservoir.
Definition: BlackOilFluidSystem.hpp:1514
static void setEnableDiffusion(bool yesno)
Specify whether the fluid system should consider diffusion.
Definition: BlackOilFluidSystem.hpp:379
static LhsEval bubblePointPressure(const FluidState &fluidState, unsigned regionIdx)
Returns the bubble point pressure $P_b$ using the current Rs.
Definition: BlackOilFluidSystem.hpp:1301
static LhsEval density(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Calculate the density [kg/m^3] of a fluid phase.
Definition: BlackOilFluidSystem.hpp:701
static LhsEval saturatedDensity(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Compute the density of a saturated fluid phase.
Definition: BlackOilFluidSystem.hpp:788
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:683
static bool isIdealMixture(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal mixture.
Definition: BlackOilFluidSystem.hpp:588
static bool enableVaporizedWater()
Returns whether the fluid system should consider that the water component can dissolve in the gas pha...
Definition: BlackOilFluidSystem.hpp:639
static const char * componentName(unsigned compIdx)
Return the human readable name of a component.
Definition: BlackOilFluidSystem.hpp:568
static const GasPvt & gasPvt()
Return a reference to the low-level object which calculates the gas phase quantities.
Definition: BlackOilFluidSystem.hpp:1486
static constexpr unsigned numPhases
Number of fluid phases in the fluid system.
Definition: BlackOilFluidSystem.hpp:467
static constexpr unsigned solventComponentIndex(unsigned phaseIdx)
returns the index of "primary" component of a phase (solvent)
Definition: BlackOilFluidSystem.hpp:536
static void setWaterPvt(std::shared_ptr< WaterPvt > pvtObj)
Set the pressure-volume-saturation (PVT) relations for the water phase.
Definition: BlackOilFluidSystem.hpp:398
static void setReservoirTemperature(Scalar value)
Return the temperature of the reservoir.
Definition: BlackOilFluidSystem.hpp:1522
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:1250
static size_t numRegions()
Returns the number of PVT regions which are considered.
Definition: BlackOilFluidSystem.hpp:612
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:1186
static constexpr unsigned oilCompIdx
Index of the oil component.
Definition: BlackOilFluidSystem.hpp:513
static bool isLiquid(unsigned phaseIdx)
Return whether a phase is liquid.
Definition: BlackOilFluidSystem.hpp:499
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:690
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:983
static LhsEval convertxoGToXoG(const LhsEval &xoG, unsigned regionIdx)
Convert a gas mole fraction in the oil phase the corresponding mass fraction.
Definition: BlackOilFluidSystem.hpp:1447
static Scalar molarMass(unsigned compIdx, unsigned regionIdx=0)
Return the molar mass of a component in [kg/mol].
Definition: BlackOilFluidSystem.hpp:584
static LhsEval saturatedDissolutionFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the dissolution factor of a saturated fluid phase.
Definition: BlackOilFluidSystem.hpp:1279
static LhsEval dewPointPressure(const FluidState &fluidState, unsigned regionIdx)
Returns the dew point pressure $P_d$ using the current Rv.
Definition: BlackOilFluidSystem.hpp:1312
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:1329
static LhsEval convertXgOToxgO(const LhsEval &XgO, unsigned regionIdx)
Convert a oil mass fraction in the gas phase the corresponding mole fraction.
Definition: BlackOilFluidSystem.hpp:1459
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:1548
static LhsEval saturatedVaporizationFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the water vaporization factor of saturated phase.
Definition: BlackOilFluidSystem.hpp:1225
static bool enableDissolvedGas()
Returns whether the fluid system should consider that the gas component can dissolve in the oil phase...
Definition: BlackOilFluidSystem.hpp:621
static Scalar referenceDensity(unsigned phaseIdx, unsigned regionIdx)
Returns the density of a fluid phase at surface pressure [kg/m^3].
Definition: BlackOilFluidSystem.hpp:655
static Scalar surfaceTemperature
The temperature at the surface.
Definition: BlackOilFluidSystem.hpp:480
static void setOilPvt(std::shared_ptr< OilPvt > pvtObj)
Set the pressure-volume-saturation (PVT) relations for the oil phase.
Definition: BlackOilFluidSystem.hpp:392
static constexpr unsigned waterPhaseIdx
Index of the water phase.
Definition: BlackOilFluidSystem.hpp:470
static LhsEval convertxgOToXgO(const LhsEval &xgO, unsigned regionIdx)
Convert a oil mole fraction in the gas phase the corresponding mass fraction.
Definition: BlackOilFluidSystem.hpp:1471
static LhsEval inverseFormationVolumeFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the formation volume factor of an "undersaturated" fluid phase.
Definition: BlackOilFluidSystem.hpp:878
static constexpr unsigned gasPhaseIdx
Index of the gas phase.
Definition: BlackOilFluidSystem.hpp:474
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:353
static const char * phaseName(unsigned phaseIdx)
Return the human readable name of a fluid phase.
Definition: BlackOilFluidSystem.hpp:483
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:1367
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:670
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:663
static constexpr unsigned gasCompIdx
Index of the gas component.
Definition: BlackOilFluidSystem.hpp:517
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:1408
static bool enableVaporizedOil()
Returns whether the fluid system should consider that the oil component can dissolve in the gas phase...
Definition: BlackOilFluidSystem.hpp:630
static bool isIdealGas(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal gas.
Definition: BlackOilFluidSystem.hpp:600
static unsigned phaseIsActive(unsigned phaseIdx)
Returns whether a fluid phase is active.
Definition: BlackOilFluidSystem.hpp:529
static Scalar surfacePressure
The pressure at the surface.
Definition: BlackOilFluidSystem.hpp:477
static const WaterPvt & waterPvt()
Return a reference to the low-level object which calculates the water phase quantities.
Definition: BlackOilFluidSystem.hpp:1506
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:371
static void initBegin(size_t numPvtRegions)
Begin the initialization of the black oil fluid system.
Definition: BlackOilFluidSystem.hpp:324
static constexpr unsigned oilPhaseIdx
Index of the oil phase.
Definition: BlackOilFluidSystem.hpp:472
static bool enableDiffusion()
Returns whether the fluid system should consider diffusion.
Definition: BlackOilFluidSystem.hpp:647
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:1380
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:362
static void setDiffusionCoefficient(Scalar coefficient, unsigned compIdx, unsigned phaseIdx, unsigned regionIdx=0)
Definition: BlackOilFluidSystem.hpp:1541
static const OilPvt & oilPvt()
Return a reference to the low-level object which calculates the oil phase quantities.
Definition: BlackOilFluidSystem.hpp:1496
static bool isCompressible(unsigned)
Returns true if and only if a fluid phase is assumed to be compressible.
Definition: BlackOilFluidSystem.hpp:596
static void setGasPvt(std::shared_ptr< GasPvt > pvtObj)
Set the pressure-volume-saturation (PVT) relations for the gas phase.
Definition: BlackOilFluidSystem.hpp:386
static Scalar molarMass()
The molar mass in of the component.
Definition: Brine.hpp:80
static Scalar molarMass()
The mass in [kg] of one mole of CO2.
Definition: CO2.hpp:66
A central place for various physical constants occuring in some equations.
Definition: Constants.hpp:41
This class represents the Pressure-Volume-Temperature relations of the gas phase in the black-oil mod...
Definition: GasPvtMultiplexer.hpp:100
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:322
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:96
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:271
This class represents the Pressure-Volume-Temperature relations of the water phase in the black-oil m...
Definition: WaterPvtMultiplexer.hpp:75
The type of the fluid system's parameter cache.
Definition: BlackOilFluidSystem.hpp:152
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:170
void setRegionIndex(unsigned val)
Set the index of the region which should be used to determine the thermodynamic properties.
Definition: BlackOilFluidSystem.hpp:193
unsigned regionIndex() const
Return the index of the region which should be used to determine the thermodynamic properties.
Definition: BlackOilFluidSystem.hpp:183