27 #ifndef OPM_BLACK_OIL_FLUID_SYSTEM_HPP
28 #define OPM_BLACK_OIL_FLUID_SYSTEM_HPP
45 #include <opm/parser/eclipse/EclipseState/EclipseState.hpp>
46 #include <opm/parser/eclipse/EclipseState/Tables/FlatTable.hpp>
47 #include <opm/parser/eclipse/EclipseState/Tables/TableManager.hpp>
60 template <class FluidSystem, class FluidState, class LhsEval>
61 LhsEval getRs_(typename std::enable_if<!HasMember_Rs<FluidState>::value, const FluidState&>::type fluidState,
65 decay<LhsEval>(fluidState.massFraction(FluidSystem::oilPhaseIdx, FluidSystem::gasCompIdx));
66 return FluidSystem::convertXoGToRs(XoG, regionIdx);
69 template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
70 auto getRs_(
typename std::enable_if<HasMember_Rs<FluidState>::value,
const FluidState&>::type fluidState,
72 -> decltype(decay<LhsEval>(fluidState.Rs()))
73 {
return decay<LhsEval>(fluidState.Rs()); }
75 template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
76 LhsEval getRv_(
typename std::enable_if<!HasMember_Rv<FluidState>::value,
const FluidState&>::type fluidState,
80 decay<LhsEval>(fluidState.massFraction(FluidSystem::gasPhaseIdx, FluidSystem::oilCompIdx));
81 return FluidSystem::convertXgOToRv(XgO, regionIdx);
84 template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
85 auto getRv_(
typename std::enable_if<HasMember_Rv<FluidState>::value,
const FluidState&>::type fluidState,
87 -> decltype(decay<LhsEval>(fluidState.Rv()))
88 {
return decay<LhsEval>(fluidState.Rv()); }
90 template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
91 LhsEval getSaltConcentration_(
typename std::enable_if<!HasMember_saltConcentration<FluidState>::value,
92 const FluidState&>::type,
96 template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
97 auto getSaltConcentration_(
typename std::enable_if<HasMember_saltConcentration<FluidState>::value,
const FluidState&>::type fluidState,
99 -> decltype(decay<LhsEval>(fluidState.saltConcentration()))
100 {
return decay<LhsEval>(fluidState.saltConcentration()); }
110 template <
class Scalar,
class IndexTraits = BlackOilDefaultIndexTraits>
121 template <
class EvaluationT>
124 typedef EvaluationT Evaluation;
129 maxOilSat_ = maxOilSat;
130 regionIdx_ = regionIdx;
140 template <
class OtherCache>
143 regionIdx_ = other.regionIndex();
144 maxOilSat_ = other.maxOilSat();
155 {
return regionIdx_; }
165 { regionIdx_ = val; }
167 const Evaluation& maxOilSat()
const
168 {
return maxOilSat_; }
170 void setMaxOilSat(
const Evaluation& val)
171 { maxOilSat_ = val; }
174 Evaluation maxOilSat_;
185 static void initFromState(
const EclipseState& eclState,
const Schedule& schedule)
187 size_t numRegions = eclState.runspec().tabdims().getNumPVTTables();
190 numActivePhases_ = 0;
191 std::fill_n(&phaseIsActive_[0],
numPhases,
false);
194 if (eclState.runspec().phases().active(Phase::OIL)) {
199 if (eclState.runspec().phases().active(Phase::GAS)) {
204 if (eclState.runspec().phases().active(Phase::WATER)) {
218 assert(numActivePhases_ >= 1 && numActivePhases_ <= 3);
224 gasPvt_ = std::make_shared<GasPvt>();
225 gasPvt_->initFromState(eclState, schedule);
229 oilPvt_ = std::make_shared<OilPvt>();
230 oilPvt_->initFromState(eclState, schedule);
234 waterPvt_ = std::make_shared<WaterPvt>();
235 waterPvt_->initFromState(eclState, schedule);
239 for (
unsigned regionIdx = 0; regionIdx <
numRegions; ++regionIdx) {
253 if (eclState.runspec().co2Storage()) {
254 for (
unsigned regionIdx = 0; regionIdx <
numRegions; ++regionIdx) {
262 const auto& diffCoeffTables = eclState.getTableManager().getDiffusionCoefficientTable();
263 if(!diffCoeffTables.empty()) {
266 diffusionCoefficients_.resize(
numRegions,{0,0,0,0,0,0,0,0,0});
268 for (
unsigned regionIdx = 0; regionIdx <
numRegions; ++regionIdx) {
269 const auto& diffCoeffTable = diffCoeffTables[regionIdx];
270 molarMass_[regionIdx][
oilCompIdx] = diffCoeffTable.oil_mw;
271 molarMass_[regionIdx][
gasCompIdx] = diffCoeffTable.gas_mw;
276 if(diffCoeffTable.gas_in_oil_cross_phase > 0 || diffCoeffTable.oil_in_oil_cross_phase > 0) {
277 throw std::runtime_error(
"Cross phase diffusion is set in the deck, but not implemented in Flow. "
278 "Please default DIFFC item 7 and item 8 or set it to zero.");
296 isInitialized_ =
false;
298 enableDissolvedGas_ =
true;
299 enableVaporizedOil_ =
false;
300 enableDiffusion_ =
false;
311 std::fill_n(&phaseIsActive_[0],
numPhases,
true);
313 resizeArrays_(numPvtRegions);
323 { enableDissolvedGas_ = yesno; }
332 { enableVaporizedOil_ = yesno; }
340 { enableDiffusion_ = yesno; }
347 { gasPvt_ = pvtObj; }
353 { oilPvt_ = pvtObj; }
359 { waterPvt_ = pvtObj; }
373 referenceDensity_[regionIdx][
oilPhaseIdx] = rhoOil;
375 referenceDensity_[regionIdx][
gasPhaseIdx] = rhoGas;
386 for (
unsigned regionIdx = 0; regionIdx <
numRegions; ++ regionIdx) {
408 int activePhaseIdx = 0;
409 for (
unsigned phaseIdx = 0; phaseIdx <
numPhases; ++phaseIdx) {
411 canonicalToActivePhaseIdx_[phaseIdx] = activePhaseIdx;
412 activeToCanonicalPhaseIdx_[activePhaseIdx] = phaseIdx;
416 isInitialized_ =
true;
419 static bool isInitialized()
420 {
return isInitialized_; }
454 throw std::logic_error(
"Phase index " + std::to_string(phaseIdx) +
" is unknown");
480 static unsigned char numActivePhases_;
481 static std::array<bool,numPhases> phaseIsActive_;
486 {
return numActivePhases_; }
492 return phaseIsActive_[phaseIdx];
507 throw std::logic_error(
"Phase index " + std::to_string(phaseIdx) +
" is unknown");
516 throw std::logic_error(
"The water phase does not have any solutes in the black oil model!");
523 throw std::logic_error(
"Phase index " + std::to_string(phaseIdx) +
" is unknown");
539 throw std::logic_error(
"Component index " + std::to_string(compIdx) +
" is unknown");
545 {
return molarMass_[regionIdx][compIdx]; }
573 {
return molarMass_.size(); }
582 {
return enableDissolvedGas_; }
591 {
return enableVaporizedOil_; }
599 {
return enableDiffusion_; }
607 {
return referenceDensity_[regionIdx][phaseIdx]; }
613 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
614 static LhsEval
density(
const FluidState& fluidState,
617 {
return density<FluidState, LhsEval>(fluidState, phaseIdx, paramCache.
regionIndex()); }
620 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
626 return fugacityCoefficient<FluidState, LhsEval>(fluidState,
633 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
637 {
return viscosity<FluidState, LhsEval>(fluidState, phaseIdx, paramCache.
regionIndex()); }
640 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
641 static LhsEval
enthalpy(
const FluidState& fluidState,
644 {
return enthalpy<FluidState, LhsEval>(fluidState, phaseIdx, paramCache.
regionIndex()); }
651 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
652 static LhsEval
density(
const FluidState& fluidState,
659 const LhsEval& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
660 const LhsEval& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
661 const LhsEval& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
667 const LhsEval& Rs = BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
668 const LhsEval& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
676 const LhsEval Rs(0.0);
677 const auto& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
685 const LhsEval& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
686 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv);
694 const LhsEval Rv(0.0);
695 const auto& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv);
702 * waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, saltConcentration);
705 throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
715 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
723 const auto& p = fluidState.pressure(phaseIdx);
724 const auto& T = fluidState.temperature(phaseIdx);
730 const LhsEval& Rs = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState,
oilPhaseIdx, regionIdx);
731 const LhsEval& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
739 const LhsEval Rs(0.0);
740 const LhsEval& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
747 const LhsEval& Rv = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState,
gasPhaseIdx, regionIdx);
748 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv);
756 const LhsEval Rv(0.0);
757 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv);
766 *inverseFormationVolumeFactor<FluidState, LhsEval>(fluidState,
waterPhaseIdx, regionIdx);
769 throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
780 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
788 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
789 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
790 const auto& saltConcentration = decay<LhsEval>(fluidState.saltConcentration());
795 const auto& Rs = BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
797 && Rs >= (1.0 - 1e-10)*oilPvt_->saturatedGasDissolutionFactor(regionIdx, scalarValue(T), scalarValue(p)))
799 return oilPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
801 return oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
805 const LhsEval Rs(0.0);
806 return oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
810 const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
812 && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
814 return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
816 return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv);
820 const LhsEval Rv(0.0);
821 return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv);
824 return waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, saltConcentration);
825 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
836 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
844 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
845 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
846 const auto& saltConcentration = decay<LhsEval>(fluidState.saltConcentration());
849 case oilPhaseIdx:
return oilPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
850 case gasPhaseIdx:
return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
851 case waterPhaseIdx:
return waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, saltConcentration);
852 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
857 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
867 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
868 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
873 const LhsEval phi_oO = 20e3/p;
876 const Scalar phi_gG = 1.0;
880 const LhsEval phi_wW = 30e3/p;
898 const auto& R_vSat = gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p);
902 const auto& R_sSat = oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p);
905 const auto& x_oOSat = 1.0 - x_oGSat;
907 const auto& p_o = decay<LhsEval>(fluidState.pressure(
oilPhaseIdx));
908 const auto& p_g = decay<LhsEval>(fluidState.pressure(
gasPhaseIdx));
910 return phi_oO*p_o*x_oOSat / (p_g*x_gOSat);
918 throw std::logic_error(
"Invalid component index "+std::to_string(compIdx));
934 const auto& R_vSat = gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p);
937 const auto& x_gGSat = 1.0 - x_gOSat;
939 const auto& R_sSat = oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p);
943 const auto& p_o = decay<LhsEval>(fluidState.pressure(
oilPhaseIdx));
944 const auto& p_g = decay<LhsEval>(fluidState.pressure(
gasPhaseIdx));
946 return phi_gG*p_g*x_gGSat / (p_o*x_oGSat);
953 throw std::logic_error(
"Invalid component index "+std::to_string(compIdx));
968 throw std::logic_error(
"Invalid component index "+std::to_string(compIdx));
972 throw std::logic_error(
"Invalid phase index "+std::to_string(phaseIdx));
975 throw std::logic_error(
"Unhandled phase or component index");
979 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
987 const LhsEval& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
988 const LhsEval& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
989 const LhsEval& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
994 const auto& Rs = BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
996 && Rs >= (1.0 - 1e-10)*oilPvt_->saturatedGasDissolutionFactor(regionIdx, scalarValue(T), scalarValue(p)))
998 return oilPvt_->saturatedViscosity(regionIdx, T, p);
1000 return oilPvt_->viscosity(regionIdx, T, p, Rs);
1004 const LhsEval Rs(0.0);
1005 return oilPvt_->viscosity(regionIdx, T, p, Rs);
1010 const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
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);
1020 const LhsEval Rv(0.0);
1021 return gasPvt_->viscosity(regionIdx, T, p, Rv);
1027 return waterPvt_->viscosity(regionIdx, T, p, saltConcentration);
1030 throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1034 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1042 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1043 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1048 oilPvt_->internalEnergy(regionIdx, T, p, BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx))
1049 + p/density<FluidState, LhsEval>(fluidState, phaseIdx, regionIdx);
1053 gasPvt_->internalEnergy(regionIdx, T, p, BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx))
1054 + p/density<FluidState, LhsEval>(fluidState, phaseIdx, regionIdx);
1058 waterPvt_->internalEnergy(regionIdx, T, p)
1059 + p/density<FluidState, LhsEval>(fluidState, phaseIdx, regionIdx);
1061 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1064 throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1073 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1077 const LhsEval& maxOilSaturation)
1082 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1083 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1084 const auto& So = decay<LhsEval>(fluidState.saturation(
oilPhaseIdx));
1087 case oilPhaseIdx:
return oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p, So, maxOilSaturation);
1088 case gasPhaseIdx:
return gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p, So, maxOilSaturation);
1090 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1102 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1110 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1111 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1114 case oilPhaseIdx:
return oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p);
1115 case gasPhaseIdx:
return gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p);
1117 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1124 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1135 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1152 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1160 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1163 case oilPhaseIdx:
return oilPvt_->saturationPressure(regionIdx, T, BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1164 case gasPhaseIdx:
return gasPvt_->saturationPressure(regionIdx, T, BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1166 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1177 template <
class LhsEval>
1183 return XoG/(1.0 - XoG)*(rho_oRef/rho_gRef);
1190 template <
class LhsEval>
1196 return XgO/(1.0 - XgO)*(rho_gRef/rho_oRef);
1203 template <
class LhsEval>
1209 const LhsEval& rho_oG = Rs*rho_gRef;
1210 return rho_oG/(rho_oRef + rho_oG);
1217 template <
class LhsEval>
1223 const LhsEval& rho_gO = Rv*rho_oRef;
1224 return rho_gO/(rho_gRef + rho_gO);
1230 template <
class LhsEval>
1236 return XoG*MO / (MG*(1 - XoG) + XoG*MO);
1242 template <
class LhsEval>
1248 return xoG*MG / (xoG*(MG - MO) + MO);
1254 template <
class LhsEval>
1260 return XgO*MG / (MO*(1 - XgO) + XgO*MG);
1266 template <
class LhsEval>
1272 return xgO*MO / (xgO*(MO - MG) + MG);
1283 {
return *gasPvt_; }
1293 {
return *oilPvt_; }
1303 {
return *waterPvt_; }
1311 {
return reservoirTemperature_; }
1319 { reservoirTemperature_ = value; }
1321 static short activeToCanonicalPhaseIdx(
unsigned activePhaseIdx) {
1323 return activeToCanonicalPhaseIdx_[activePhaseIdx];
1326 static short canonicalToActivePhaseIdx(
unsigned phaseIdx) {
1329 return canonicalToActivePhaseIdx_[phaseIdx];
1334 {
return diffusionCoefficients_[regionIdx][
numPhases*compIdx + phaseIdx]; }
1338 { diffusionCoefficients_[regionIdx][
numPhases*compIdx + phaseIdx] = coefficient ; }
1343 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
1354 if(!diffusionCoefficients_.empty()) {
1358 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1359 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1365 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1376 static Scalar reservoirTemperature_;
1378 static std::shared_ptr<GasPvt> gasPvt_;
1379 static std::shared_ptr<OilPvt> oilPvt_;
1380 static std::shared_ptr<WaterPvt> waterPvt_;
1382 static bool enableDissolvedGas_;
1383 static bool enableVaporizedOil_;
1384 static bool enableDiffusion_;
1389 static std::vector<std::array<
Scalar, 3> > referenceDensity_;
1390 static std::vector<std::array<
Scalar, 3> > molarMass_;
1391 static std::vector<std::array<
Scalar, 3 * 3> > diffusionCoefficients_;
1393 static std::array<short, numPhases> activeToCanonicalPhaseIdx_;
1394 static std::array<short, numPhases> canonicalToActivePhaseIdx_;
1396 static bool isInitialized_;
1399 template <
class Scalar,
class IndexTraits>
1400 unsigned char BlackOilFluidSystem<Scalar, IndexTraits>::numActivePhases_;
1402 template <
class Scalar,
class IndexTraits>
1403 std::array<bool, BlackOilFluidSystem<Scalar, IndexTraits>::numPhases> BlackOilFluidSystem<Scalar, IndexTraits>::phaseIsActive_;
1405 template <
class Scalar,
class IndexTraits>
1406 std::array<short, BlackOilFluidSystem<Scalar, IndexTraits>::numPhases> BlackOilFluidSystem<Scalar, IndexTraits>::activeToCanonicalPhaseIdx_;
1408 template <
class Scalar,
class IndexTraits>
1409 std::array<short, BlackOilFluidSystem<Scalar, IndexTraits>::numPhases> BlackOilFluidSystem<Scalar, IndexTraits>::canonicalToActivePhaseIdx_;
1411 template <
class Scalar,
class IndexTraits>
1415 template <
class Scalar,
class IndexTraits>
1419 template <
class Scalar,
class IndexTraits>
1421 BlackOilFluidSystem<Scalar, IndexTraits>::reservoirTemperature_;
1423 template <
class Scalar,
class IndexTraits>
1424 bool BlackOilFluidSystem<Scalar, IndexTraits>::enableDissolvedGas_;
1426 template <
class Scalar,
class IndexTraits>
1427 bool BlackOilFluidSystem<Scalar, IndexTraits>::enableVaporizedOil_;
1429 template <
class Scalar,
class IndexTraits>
1430 bool BlackOilFluidSystem<Scalar, IndexTraits>::enableDiffusion_;
1432 template <
class Scalar,
class IndexTraits>
1433 std::shared_ptr<OilPvtMultiplexer<Scalar> >
1434 BlackOilFluidSystem<Scalar, IndexTraits>::oilPvt_;
1436 template <
class Scalar,
class IndexTraits>
1437 std::shared_ptr<GasPvtMultiplexer<Scalar> >
1438 BlackOilFluidSystem<Scalar, IndexTraits>::gasPvt_;
1440 template <
class Scalar,
class IndexTraits>
1441 std::shared_ptr<WaterPvtMultiplexer<Scalar> >
1442 BlackOilFluidSystem<Scalar, IndexTraits>::waterPvt_;
1444 template <
class Scalar,
class IndexTraits>
1445 std::vector<std::array<Scalar, 3> >
1446 BlackOilFluidSystem<Scalar, IndexTraits>::referenceDensity_;
1448 template <
class Scalar,
class IndexTraits>
1449 std::vector<std::array<Scalar, 3> >
1450 BlackOilFluidSystem<Scalar, IndexTraits>::molarMass_;
1452 template <
class Scalar,
class IndexTraits>
1453 std::vector<std::array<Scalar, 9> >
1454 BlackOilFluidSystem<Scalar, IndexTraits>::diffusionCoefficients_;
1456 template <
class Scalar,
class IndexTraits>
1457 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:112
static unsigned numActivePhases()
Returns the number of active fluid phases (i.e., usually three)
Definition: BlackOilFluidSystem.hpp:485
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:1333
static void setReferenceDensities(Scalar rhoOil, Scalar rhoWater, Scalar rhoGas, unsigned regionIdx)
Initialize the values of the reference densities.
Definition: BlackOilFluidSystem.hpp:368
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:1178
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:1204
static void initEnd()
Finish initializing the black oil fluid system.
Definition: BlackOilFluidSystem.hpp:382
static constexpr unsigned soluteComponentIndex(unsigned phaseIdx)
returns the index of "secondary" component of a phase (solute)
Definition: BlackOilFluidSystem.hpp:512
static LhsEval saturatedInverseFormationVolumeFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the formation volume factor of a "saturated" fluid phase.
Definition: BlackOilFluidSystem.hpp:837
static LhsEval convertXoGToxoG(const LhsEval &XoG, unsigned regionIdx)
Convert a gas mass fraction in the oil phase the corresponding mole fraction.
Definition: BlackOilFluidSystem.hpp:1231
static LhsEval viscosity(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Calculate the dynamic viscosity of a fluid phase [Pa*s].
Definition: BlackOilFluidSystem.hpp:980
static const unsigned oilPhaseIdx
Index of the oil phase.
Definition: BlackOilFluidSystem.hpp:432
static Scalar reservoirTemperature(unsigned=0)
Set the temperature of the reservoir.
Definition: BlackOilFluidSystem.hpp:1310
static void setEnableDiffusion(bool yesno)
Specify whether the fluid system should consider diffusion.
Definition: BlackOilFluidSystem.hpp:339
static const unsigned waterPhaseIdx
Index of the water phase.
Definition: BlackOilFluidSystem.hpp:430
static LhsEval bubblePointPressure(const FluidState &fluidState, unsigned regionIdx)
Returns the bubble point pressure $P_b$ using the current Rs.
Definition: BlackOilFluidSystem.hpp:1125
static LhsEval density(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Calculate the density [kg/m^3] of a fluid phase.
Definition: BlackOilFluidSystem.hpp:652
static LhsEval saturatedDensity(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Compute the density of a saturated fluid phase.
Definition: BlackOilFluidSystem.hpp:716
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:634
static bool isIdealMixture(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal mixture.
Definition: BlackOilFluidSystem.hpp:548
static const unsigned numComponents
Number of chemical species in the fluid system.
Definition: BlackOilFluidSystem.hpp:470
static const char * componentName(unsigned compIdx)
Return the human readable name of a component.
Definition: BlackOilFluidSystem.hpp:528
static const GasPvt & gasPvt()
Return a reference to the low-level object which calculates the gas phase quantities.
Definition: BlackOilFluidSystem.hpp:1282
static constexpr unsigned solventComponentIndex(unsigned phaseIdx)
returns the index of "primary" component of a phase (solvent)
Definition: BlackOilFluidSystem.hpp:496
static void setWaterPvt(std::shared_ptr< WaterPvt > pvtObj)
Set the pressure-volume-saturation (PVT) relations for the water phase.
Definition: BlackOilFluidSystem.hpp:358
static void setReservoirTemperature(Scalar value)
Return the temperature of the reservoir.
Definition: BlackOilFluidSystem.hpp:1318
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:1074
static size_t numRegions()
Returns the number of PVT regions which are considered.
Definition: BlackOilFluidSystem.hpp:572
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:1035
static bool isLiquid(unsigned phaseIdx)
Return whether a phase is liquid.
Definition: BlackOilFluidSystem.hpp:459
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:641
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:858
static LhsEval convertxoGToXoG(const LhsEval &xoG, unsigned regionIdx)
Convert a gas mole fraction in the oil phase the corresponding mass fraction.
Definition: BlackOilFluidSystem.hpp:1243
static Scalar molarMass(unsigned compIdx, unsigned regionIdx=0)
Return the molar mass of a component in [kg/mol].
Definition: BlackOilFluidSystem.hpp:544
static LhsEval saturatedDissolutionFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the dissolution factor of a saturated fluid phase.
Definition: BlackOilFluidSystem.hpp:1103
static LhsEval dewPointPressure(const FluidState &fluidState, unsigned regionIdx)
Returns the dew point pressure $P_d$ using the current Rv.
Definition: BlackOilFluidSystem.hpp:1136
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:1153
static LhsEval convertXgOToxgO(const LhsEval &XgO, unsigned regionIdx)
Convert a oil mass fraction in the gas phase the corresponding mole fraction.
Definition: BlackOilFluidSystem.hpp:1255
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:1344
static bool enableDissolvedGas()
Returns whether the fluid system should consider that the gas component can dissolve in the oil phase...
Definition: BlackOilFluidSystem.hpp:581
static Scalar referenceDensity(unsigned phaseIdx, unsigned regionIdx)
Returns the density of a fluid phase at surface pressure [kg/m^3].
Definition: BlackOilFluidSystem.hpp:606
static Scalar surfaceTemperature
The temperature at the surface.
Definition: BlackOilFluidSystem.hpp:440
static void setOilPvt(std::shared_ptr< OilPvt > pvtObj)
Set the pressure-volume-saturation (PVT) relations for the oil phase.
Definition: BlackOilFluidSystem.hpp:352
static LhsEval convertxgOToXgO(const LhsEval &xgO, unsigned regionIdx)
Convert a oil mole fraction in the gas phase the corresponding mass fraction.
Definition: BlackOilFluidSystem.hpp:1267
static LhsEval inverseFormationVolumeFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the formation volume factor of an "undersaturated" fluid phase.
Definition: BlackOilFluidSystem.hpp:781
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:322
static const char * phaseName(unsigned phaseIdx)
Return the human readable name of a fluid phase.
Definition: BlackOilFluidSystem.hpp:443
static const unsigned numPhases
Number of fluid phases in the fluid system.
Definition: BlackOilFluidSystem.hpp:427
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:1191
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:621
static const unsigned gasPhaseIdx
Index of the gas phase.
Definition: BlackOilFluidSystem.hpp:434
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:614
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:1218
static bool enableVaporizedOil()
Returns whether the fluid system should consider that the oil component can dissolve in the gas phase...
Definition: BlackOilFluidSystem.hpp:590
static bool isIdealGas(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal gas.
Definition: BlackOilFluidSystem.hpp:560
static unsigned phaseIsActive(unsigned phaseIdx)
Returns whether a fluid phase is active.
Definition: BlackOilFluidSystem.hpp:489
static Scalar surfacePressure
The pressure at the surface.
Definition: BlackOilFluidSystem.hpp:437
static const WaterPvt & waterPvt()
Return a reference to the low-level object which calculates the water phase quantities.
Definition: BlackOilFluidSystem.hpp:1302
static void initBegin(size_t numPvtRegions)
Begin the initialization of the black oil fluid system.
Definition: BlackOilFluidSystem.hpp:294
static bool enableDiffusion()
Returns whether the fluid system should consider diffusion.
Definition: BlackOilFluidSystem.hpp:598
static const unsigned gasCompIdx
Index of the gas component.
Definition: BlackOilFluidSystem.hpp:477
static const unsigned oilCompIdx
Index of the oil component.
Definition: BlackOilFluidSystem.hpp:473
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:331
static const unsigned waterCompIdx
Index of the water component.
Definition: BlackOilFluidSystem.hpp:475
static void setDiffusionCoefficient(Scalar coefficient, unsigned compIdx, unsigned phaseIdx, unsigned regionIdx=0)
Definition: BlackOilFluidSystem.hpp:1337
static const OilPvt & oilPvt()
Return a reference to the low-level object which calculates the oil phase quantities.
Definition: BlackOilFluidSystem.hpp:1292
static bool isCompressible(unsigned)
Returns true if and only if a fluid phase is assumed to be compressible.
Definition: BlackOilFluidSystem.hpp:556
static void setGasPvt(std::shared_ptr< GasPvt > pvtObj)
Set the pressure-volume-saturation (PVT) relations for the gas phase.
Definition: BlackOilFluidSystem.hpp:346
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:86
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:276
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:123
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:141
void setRegionIndex(unsigned val)
Set the index of the region which should be used to determine the thermodynamic properties.
Definition: BlackOilFluidSystem.hpp:164
unsigned regionIndex() const
Return the index of the region which should be used to determine the thermodynamic properties.
Definition: BlackOilFluidSystem.hpp:154