27#ifndef OPM_BLACK_OIL_FLUID_SYSTEM_HPP
28#define OPM_BLACK_OIL_FLUID_SYSTEM_HPP
36#include <opm/common/TimingMacros.hpp>
60OPM_GENERATE_HAS_MEMBER(Rv, )
61OPM_GENERATE_HAS_MEMBER(Rvw, )
62OPM_GENERATE_HAS_MEMBER(Rsw, )
63OPM_GENERATE_HAS_MEMBER(saltConcentration, )
64OPM_GENERATE_HAS_MEMBER(saltSaturation, )
66template <class FluidSystem, class FluidState, class LhsEval>
67LhsEval getRs_(typename std::enable_if<!HasMember_Rs<FluidState>::value, const FluidState&>::type fluidState,
71 decay<LhsEval>(fluidState.massFraction(FluidSystem::oilPhaseIdx, FluidSystem::gasCompIdx));
72 return FluidSystem::convertXoGToRs(XoG, regionIdx);
75template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
76auto getRs_(
typename std::enable_if<HasMember_Rs<FluidState>::value,
const FluidState&>::type fluidState,
78 ->
decltype(decay<LhsEval>(fluidState.Rs()))
79{
return decay<LhsEval>(fluidState.Rs()); }
81template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
82LhsEval getRv_(
typename std::enable_if<!HasMember_Rv<FluidState>::value,
const FluidState&>::type fluidState,
86 decay<LhsEval>(fluidState.massFraction(FluidSystem::gasPhaseIdx, FluidSystem::oilCompIdx));
87 return FluidSystem::convertXgOToRv(XgO, regionIdx);
90template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
91auto getRv_(
typename std::enable_if<HasMember_Rv<FluidState>::value,
const FluidState&>::type fluidState,
93 ->
decltype(decay<LhsEval>(fluidState.Rv()))
94{
return decay<LhsEval>(fluidState.Rv()); }
96template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
97LhsEval getRvw_(
typename std::enable_if<!HasMember_Rvw<FluidState>::value,
const FluidState&>::type fluidState,
101 decay<LhsEval>(fluidState.massFraction(FluidSystem::gasPhaseIdx, FluidSystem::waterCompIdx));
102 return FluidSystem::convertXgWToRvw(XgW, regionIdx);
105template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
106auto getRvw_(
typename std::enable_if<HasMember_Rvw<FluidState>::value,
const FluidState&>::type fluidState,
108 ->
decltype(decay<LhsEval>(fluidState.Rvw()))
109{
return decay<LhsEval>(fluidState.Rvw()); }
111template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
112LhsEval getRsw_(
typename std::enable_if<!HasMember_Rsw<FluidState>::value,
const FluidState&>::type fluidState,
116 decay<LhsEval>(fluidState.massFraction(FluidSystem::waterPhaseIdx, FluidSystem::gasCompIdx));
117 return FluidSystem::convertXwGToRsw(XwG, regionIdx);
120template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
121auto getRsw_(
typename std::enable_if<HasMember_Rsw<FluidState>::value,
const FluidState&>::type fluidState,
123 ->
decltype(decay<LhsEval>(fluidState.Rsw()))
124{
return decay<LhsEval>(fluidState.Rsw()); }
126template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
127LhsEval getSaltConcentration_(
typename std::enable_if<!HasMember_saltConcentration<FluidState>::value,
128 const FluidState&>::type,
132template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
133auto getSaltConcentration_(
typename std::enable_if<HasMember_saltConcentration<FluidState>::value,
const FluidState&>::type fluidState,
135 ->
decltype(decay<LhsEval>(fluidState.saltConcentration()))
136{
return decay<LhsEval>(fluidState.saltConcentration()); }
138template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
139LhsEval getSaltSaturation_(
typename std::enable_if<!HasMember_saltSaturation<FluidState>::value,
140 const FluidState&>::type,
144template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
145auto getSaltSaturation_(
typename std::enable_if<HasMember_saltSaturation<FluidState>::value,
const FluidState&>::type fluidState,
147 ->
decltype(decay<LhsEval>(fluidState.saltSaturation()))
148{
return decay<LhsEval>(fluidState.saltSaturation()); }
158template <
class Scalar,
class IndexTraits = BlackOilDefaultIndexTraits>
169 template <
class EvaluationT>
172 using Evaluation = EvaluationT;
177 maxOilSat_ = maxOilSat;
178 regionIdx_ = regionIdx;
188 template <
class OtherCache>
191 regionIdx_ = other.regionIndex();
192 maxOilSat_ = other.maxOilSat();
203 {
return regionIdx_; }
213 { regionIdx_ = val; }
215 const Evaluation& maxOilSat()
const
216 {
return maxOilSat_; }
218 void setMaxOilSat(
const Evaluation& val)
219 { maxOilSat_ = val; }
222 Evaluation maxOilSat_;
233 static void initFromState(
const EclipseState& eclState,
const Schedule& schedule);
253 { enableDissolvedGas_ = yesno; }
262 { enableVaporizedOil_ = yesno; }
271 { enableVaporizedWater_ = yesno; }
280 { enableDissolvedGasInWater_ = yesno; }
287 { enableDiffusion_ = yesno; }
294 { gasPvt_ = pvtObj; }
300 { oilPvt_ = pvtObj; }
306 { waterPvt_ = pvtObj; }
325 static bool isInitialized()
326 {
return isInitialized_; }
366 static constexpr unsigned oilCompIdx = IndexTraits::oilCompIdx;
370 static constexpr unsigned gasCompIdx = IndexTraits::gasCompIdx;
373 static unsigned char numActivePhases_;
374 static std::array<bool,numPhases> phaseIsActive_;
379 {
return numActivePhases_; }
385 return phaseIsActive_[phaseIdx];
399 {
return molarMass_[regionIdx][compIdx]; }
427 {
return molarMass_.size(); }
436 {
return enableDissolvedGas_; }
446 {
return enableDissolvedGasInWater_; }
455 {
return enableVaporizedOil_; }
464 {
return enableVaporizedWater_; }
472 {
return enableDiffusion_; }
480 {
return referenceDensity_[regionIdx][phaseIdx]; }
486 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
487 static LhsEval
density(
const FluidState& fluidState,
490 {
return density<FluidState, LhsEval>(fluidState, phaseIdx, paramCache.
regionIndex()); }
493 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
499 return fugacityCoefficient<FluidState, LhsEval>(fluidState,
506 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
510 {
return viscosity<FluidState, LhsEval>(fluidState, phaseIdx, paramCache.
regionIndex()); }
513 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
514 static LhsEval
enthalpy(
const FluidState& fluidState,
517 {
return enthalpy<FluidState, LhsEval>(fluidState, phaseIdx, paramCache.
regionIndex()); }
524 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
525 static LhsEval
density(
const FluidState& fluidState,
532 const LhsEval& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
533 const LhsEval& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
534 const LhsEval& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
540 const LhsEval& Rs = BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
541 const LhsEval& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
549 const LhsEval Rs(0.0);
550 const auto& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
558 const LhsEval& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
559 const LhsEval& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
560 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
569 const LhsEval Rvw(0.0);
570 const LhsEval& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
571 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
579 const LhsEval Rv(0.0);
580 const LhsEval& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
581 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
589 const LhsEval Rv(0.0);
590 const LhsEval Rvw(0.0);
591 const auto& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
598 const LhsEval& Rsw = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState,
waterPhaseIdx, regionIdx);
599 const LhsEval& bw = waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
604 const LhsEval Rsw(0.0);
607 * waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
610 throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
622 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
630 const auto& p = fluidState.pressure(phaseIdx);
631 const auto& T = fluidState.temperature(phaseIdx);
637 const LhsEval& Rs = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState,
oilPhaseIdx, regionIdx);
638 const LhsEval& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
646 const LhsEval Rs(0.0);
647 const LhsEval& bo = oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
654 const LhsEval& Rv = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState,
gasPhaseIdx, regionIdx);
655 const LhsEval& Rvw = saturatedVaporizationFactor<FluidState, LhsEval>(fluidState,
gasPhaseIdx, regionIdx);
656 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
666 const LhsEval Rvw(0.0);
667 const LhsEval& Rv = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState,
gasPhaseIdx, regionIdx);
668 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
677 const LhsEval Rv(0.0);
678 const LhsEval& Rvw = saturatedVaporizationFactor<FluidState, LhsEval>(fluidState,
gasPhaseIdx, regionIdx);
679 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
687 const LhsEval Rv(0.0);
688 const LhsEval Rvw(0.0);
689 const LhsEval& bg = gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
699 const auto& saltConcentration = decay<LhsEval>(fluidState.saltConcentration());
700 const LhsEval& Rsw = saturatedDissolutionFactor<FluidState, LhsEval>(fluidState,
waterPhaseIdx, regionIdx);
701 const LhsEval& bw = waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
708 *inverseFormationVolumeFactor<FluidState, LhsEval>(fluidState,
waterPhaseIdx, regionIdx);
712 throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
723 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
732 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
733 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
738 const auto& Rs = BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
740 && Rs >= (1.0 - 1e-10)*oilPvt_->saturatedGasDissolutionFactor(regionIdx, scalarValue(T), scalarValue(p)))
742 return oilPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
744 return oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
748 const LhsEval Rs(0.0);
749 return oilPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rs);
753 const auto& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
754 const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
756 && Rvw >= (1.0 - 1e-10)*gasPvt_->saturatedWaterVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p))
758 && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
760 return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
762 return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
767 const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
769 && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
771 return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
773 const LhsEval Rvw(0.0);
774 return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
779 const auto& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
781 && Rvw >= (1.0 - 1e-10)*gasPvt_->saturatedWaterVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
783 return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
785 const LhsEval Rv(0.0);
786 return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
790 const LhsEval Rv(0.0);
791 const LhsEval Rvw(0.0);
792 return gasPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rv, Rvw);
796 const auto& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
798 const auto& Rsw = BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
800 && Rsw >= (1.0 - 1e-10)*waterPvt_->saturatedGasDissolutionFactor(regionIdx, scalarValue(T), scalarValue(p), scalarValue(saltConcentration)))
802 return waterPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p, saltConcentration);
804 return waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
807 const LhsEval Rsw(0.0);
808 return waterPvt_->inverseFormationVolumeFactor(regionIdx, T, p, Rsw, saltConcentration);
810 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
823 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
832 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
833 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
834 const auto& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
837 case oilPhaseIdx:
return oilPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
838 case gasPhaseIdx:
return gasPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p);
839 case waterPhaseIdx:
return waterPvt_->saturatedInverseFormationVolumeFactor(regionIdx, T, p, saltConcentration);
840 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
845 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
855 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
856 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
861 const LhsEval phi_oO = 20e3/p;
864 constexpr const Scalar phi_gG = 1.0;
868 const LhsEval phi_wW = 30e3/p;
886 const auto& R_vSat = gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p);
890 const auto& R_sSat = oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p);
893 const auto& x_oOSat = 1.0 - x_oGSat;
895 const auto& p_o = decay<LhsEval>(fluidState.pressure(
oilPhaseIdx));
896 const auto& p_g = decay<LhsEval>(fluidState.pressure(
gasPhaseIdx));
898 return phi_oO*p_o*x_oOSat / (p_g*x_gOSat);
906 throw std::logic_error(
"Invalid component index "+std::to_string(compIdx));
922 const auto& R_vSat = gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p);
925 const auto& x_gGSat = 1.0 - x_gOSat;
927 const auto& R_sSat = oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p);
931 const auto& p_o = decay<LhsEval>(fluidState.pressure(
oilPhaseIdx));
932 const auto& p_g = decay<LhsEval>(fluidState.pressure(
gasPhaseIdx));
934 return phi_gG*p_g*x_gGSat / (p_o*x_oGSat);
941 throw std::logic_error(
"Invalid component index "+std::to_string(compIdx));
956 throw std::logic_error(
"Invalid component index "+std::to_string(compIdx));
960 throw std::logic_error(
"Invalid phase index "+std::to_string(phaseIdx));
963 throw std::logic_error(
"Unhandled phase or component index");
967 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
976 const LhsEval& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
977 const LhsEval& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
982 const auto& Rs = BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
984 && Rs >= (1.0 - 1e-10)*oilPvt_->saturatedGasDissolutionFactor(regionIdx, scalarValue(T), scalarValue(p)))
986 return oilPvt_->saturatedViscosity(regionIdx, T, p);
988 return oilPvt_->viscosity(regionIdx, T, p, Rs);
992 const LhsEval Rs(0.0);
993 return oilPvt_->viscosity(regionIdx, T, p, Rs);
998 const auto& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
999 const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1001 && Rvw >= (1.0 - 1e-10)*gasPvt_->saturatedWaterVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p))
1003 && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
1005 return gasPvt_->saturatedViscosity(regionIdx, T, p);
1007 return gasPvt_->viscosity(regionIdx, T, p, Rv, Rvw);
1011 const auto& Rv = BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1013 && Rv >= (1.0 - 1e-10)*gasPvt_->saturatedOilVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
1015 return gasPvt_->saturatedViscosity(regionIdx, T, p);
1017 const LhsEval Rvw(0.0);
1018 return gasPvt_->viscosity(regionIdx, T, p, Rv, Rvw);
1022 const auto& Rvw = BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1024 && Rvw >= (1.0 - 1e-10)*gasPvt_->saturatedWaterVaporizationFactor(regionIdx, scalarValue(T), scalarValue(p)))
1026 return gasPvt_->saturatedViscosity(regionIdx, T, p);
1028 const LhsEval Rv(0.0);
1029 return gasPvt_->viscosity(regionIdx, T, p, Rv, Rvw);
1033 const LhsEval Rv(0.0);
1034 const LhsEval Rvw(0.0);
1035 return gasPvt_->viscosity(regionIdx, T, p, Rv, Rvw);
1040 const LhsEval& saltConcentration = BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1042 const auto& Rsw = BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx);
1044 && Rsw >= (1.0 - 1e-10)*waterPvt_->saturatedGasDissolutionFactor(regionIdx, scalarValue(T), scalarValue(p), scalarValue(saltConcentration)))
1046 return waterPvt_->saturatedViscosity(regionIdx, T, p, saltConcentration);
1048 return waterPvt_->viscosity(regionIdx, T, p, Rsw, saltConcentration);
1051 const LhsEval Rsw(0.0);
1052 return waterPvt_->viscosity(regionIdx, T, p, Rsw, saltConcentration);
1056 throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1060 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1068 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1069 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1074 oilPvt_->internalEnergy(regionIdx, T, p, BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx))
1075 + p/density<FluidState, LhsEval>(fluidState, phaseIdx, regionIdx);
1079 gasPvt_->internalEnergy(regionIdx, T, p,
1080 BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1081 BlackOil::template getRvw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx))
1082 + p/density<FluidState, LhsEval>(fluidState, phaseIdx, regionIdx);
1086 waterPvt_->internalEnergy(regionIdx, T, p,
1087 BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1088 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx))
1089 + p/density<FluidState, LhsEval>(fluidState, phaseIdx, regionIdx);
1091 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1094 throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1103 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1111 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1112 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1113 const auto& saltConcentration = decay<LhsEval>(fluidState.saltConcentration());
1117 case gasPhaseIdx:
return gasPvt_->saturatedWaterVaporizationFactor(regionIdx, T, p, saltConcentration);
1119 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1129 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1133 const LhsEval& maxOilSaturation)
1139 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1140 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1141 const auto& So = decay<LhsEval>(fluidState.saturation(
oilPhaseIdx));
1144 case oilPhaseIdx:
return oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p, So, maxOilSaturation);
1145 case gasPhaseIdx:
return gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p, So, maxOilSaturation);
1146 case waterPhaseIdx:
return waterPvt_->saturatedGasDissolutionFactor(regionIdx, T, p,
1147 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1148 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1160 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1169 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1170 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1173 case oilPhaseIdx:
return oilPvt_->saturatedGasDissolutionFactor(regionIdx, T, p);
1174 case gasPhaseIdx:
return gasPvt_->saturatedOilVaporizationFactor(regionIdx, T, p);
1175 case waterPhaseIdx:
return waterPvt_->saturatedGasDissolutionFactor(regionIdx, T, p,
1176 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1177 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1184 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1195 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1212 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar>
1220 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1223 case oilPhaseIdx:
return oilPvt_->saturationPressure(regionIdx, T, BlackOil::template getRs_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1224 case gasPhaseIdx:
return gasPvt_->saturationPressure(regionIdx, T, BlackOil::template getRv_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1225 case waterPhaseIdx:
return waterPvt_->saturationPressure(regionIdx, T,
1226 BlackOil::template getRsw_<ThisType, FluidState, LhsEval>(fluidState, regionIdx),
1227 BlackOil::template getSaltConcentration_<ThisType, FluidState, LhsEval>(fluidState, regionIdx));
1228 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1239 template <
class LhsEval>
1245 return XoG/(1.0 - XoG)*(rho_oRef/rho_gRef);
1252 template <
class LhsEval>
1258 return XwG/(1.0 - XwG)*(rho_wRef/rho_gRef);
1265 template <
class LhsEval>
1271 return XgO/(1.0 - XgO)*(rho_gRef/rho_oRef);
1278 template <
class LhsEval>
1284 return XgW/(1.0 - XgW)*(rho_gRef/rho_wRef);
1292 template <
class LhsEval>
1298 const LhsEval& rho_oG = Rs*rho_gRef;
1299 return rho_oG/(rho_oRef + rho_oG);
1306 template <
class LhsEval>
1312 const LhsEval& rho_wG = Rsw*rho_gRef;
1313 return rho_wG/(rho_wRef + rho_wG);
1320 template <
class LhsEval>
1326 const LhsEval& rho_gO = Rv*rho_oRef;
1327 return rho_gO/(rho_gRef + rho_gO);
1334 template <
class LhsEval>
1340 const LhsEval& rho_gW = Rvw*rho_wRef;
1341 return rho_gW/(rho_gRef + rho_gW);
1347 template <
class LhsEval>
1353 return XgW*MG / (MW*(1 - XgW) + XgW*MG);
1359 template <
class LhsEval>
1365 return XwG*MW / (MG*(1 - XwG) + XwG*MW);
1371 template <
class LhsEval>
1377 return XoG*MO / (MG*(1 - XoG) + XoG*MO);
1383 template <
class LhsEval>
1389 return xoG*MG / (xoG*(MG - MO) + MO);
1395 template <
class LhsEval>
1401 return XgO*MG / (MO*(1 - XgO) + XgO*MG);
1407 template <
class LhsEval>
1413 return xgO*MO / (xgO*(MO - MG) + MG);
1424 {
return *gasPvt_; }
1434 {
return *oilPvt_; }
1444 {
return *waterPvt_; }
1452 {
return reservoirTemperature_; }
1460 { reservoirTemperature_ = value; }
1462 static short activeToCanonicalPhaseIdx(
unsigned activePhaseIdx);
1464 static short canonicalToActivePhaseIdx(
unsigned phaseIdx);
1468 {
return diffusionCoefficients_[regionIdx][
numPhases*compIdx + phaseIdx]; }
1472 { diffusionCoefficients_[regionIdx][
numPhases*compIdx + phaseIdx] = coefficient ; }
1477 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
1488 if(!diffusionCoefficients_.empty()) {
1492 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
1493 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
1499 default:
throw std::logic_error(
"Unhandled phase index "+std::to_string(phaseIdx));
1504 static void resizeArrays_(std::size_t
numRegions);
1506 static Scalar reservoirTemperature_;
1508 static std::shared_ptr<GasPvt> gasPvt_;
1509 static std::shared_ptr<OilPvt> oilPvt_;
1510 static std::shared_ptr<WaterPvt> waterPvt_;
1512 static bool enableDissolvedGas_;
1513 static bool enableDissolvedGasInWater_;
1514 static bool enableVaporizedOil_;
1515 static bool enableVaporizedWater_;
1516 static bool enableDiffusion_;
1521 static std::vector<std::array<
Scalar, 3> > referenceDensity_;
1522 static std::vector<std::array<
Scalar, 3> > molarMass_;
1523 static std::vector<std::array<
Scalar, 3 * 3> > diffusionCoefficients_;
1525 static std::array<short, numPhases> activeToCanonicalPhaseIdx_;
1526 static std::array<short, numPhases> canonicalToActivePhaseIdx_;
1528 static bool isInitialized_;
1531template <
class Scalar,
class IndexTraits>
1532unsigned char BlackOilFluidSystem<Scalar, IndexTraits>::numActivePhases_;
1534template <
class Scalar,
class IndexTraits>
1535std::array<bool, BlackOilFluidSystem<Scalar, IndexTraits>::numPhases> BlackOilFluidSystem<Scalar, IndexTraits>::phaseIsActive_;
1537template <
class Scalar,
class IndexTraits>
1538std::array<short, BlackOilFluidSystem<Scalar, IndexTraits>::numPhases> BlackOilFluidSystem<Scalar, IndexTraits>::activeToCanonicalPhaseIdx_;
1540template <
class Scalar,
class IndexTraits>
1541std::array<short, BlackOilFluidSystem<Scalar, IndexTraits>::numPhases> BlackOilFluidSystem<Scalar, IndexTraits>::canonicalToActivePhaseIdx_;
1543template <
class Scalar,
class IndexTraits>
1547template <
class Scalar,
class IndexTraits>
1551template <
class Scalar,
class IndexTraits>
1553BlackOilFluidSystem<Scalar, IndexTraits>::reservoirTemperature_;
1555template <
class Scalar,
class IndexTraits>
1556bool BlackOilFluidSystem<Scalar, IndexTraits>::enableDissolvedGas_;
1558template <
class Scalar,
class IndexTraits>
1559bool BlackOilFluidSystem<Scalar, IndexTraits>::enableDissolvedGasInWater_;
1562template <
class Scalar,
class IndexTraits>
1563bool BlackOilFluidSystem<Scalar, IndexTraits>::enableVaporizedOil_;
1565template <
class Scalar,
class IndexTraits>
1566bool BlackOilFluidSystem<Scalar, IndexTraits>::enableVaporizedWater_;
1568template <
class Scalar,
class IndexTraits>
1569bool BlackOilFluidSystem<Scalar, IndexTraits>::enableDiffusion_;
1571template <
class Scalar,
class IndexTraits>
1572std::shared_ptr<OilPvtMultiplexer<Scalar> >
1573BlackOilFluidSystem<Scalar, IndexTraits>::oilPvt_;
1575template <
class Scalar,
class IndexTraits>
1576std::shared_ptr<GasPvtMultiplexer<Scalar> >
1577BlackOilFluidSystem<Scalar, IndexTraits>::gasPvt_;
1579template <
class Scalar,
class IndexTraits>
1580std::shared_ptr<WaterPvtMultiplexer<Scalar> >
1581BlackOilFluidSystem<Scalar, IndexTraits>::waterPvt_;
1583template <
class Scalar,
class IndexTraits>
1584std::vector<std::array<Scalar, 3> >
1585BlackOilFluidSystem<Scalar, IndexTraits>::referenceDensity_;
1587template <
class Scalar,
class IndexTraits>
1588std::vector<std::array<Scalar, 3> >
1589BlackOilFluidSystem<Scalar, IndexTraits>::molarMass_;
1591template <
class Scalar,
class IndexTraits>
1592std::vector<std::array<Scalar, 9> >
1593BlackOilFluidSystem<Scalar, IndexTraits>::diffusionCoefficients_;
1595template <
class Scalar,
class IndexTraits>
1596bool BlackOilFluidSystem<Scalar, IndexTraits>::isInitialized_ =
false;
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
Some templates to wrap the valgrind client request macros.
The base class for all fluid systems.
Definition: BaseFluidSystem.hpp:46
Scalar Scalar
The type used for scalar quantities.
Definition: BaseFluidSystem.hpp:51
A fluid system which uses the black-oil model assumptions to calculate termodynamically meaningful qu...
Definition: BlackOilFluidSystem.hpp:160
static constexpr unsigned waterCompIdx
Index of the water component.
Definition: BlackOilFluidSystem.hpp:368
static unsigned numActivePhases()
Returns the number of active fluid phases (i.e., usually three)
Definition: BlackOilFluidSystem.hpp:378
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:1467
static void setReferenceDensities(Scalar rhoOil, Scalar rhoWater, Scalar rhoGas, unsigned regionIdx)
Initialize the values of the reference densities.
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:1240
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:1293
static constexpr unsigned numComponents
Number of chemical species in the fluid system.
Definition: BlackOilFluidSystem.hpp:363
static void initEnd()
Finish initializing the black oil fluid system.
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:1335
static LhsEval convertXwGToxwG(const LhsEval &XwG, unsigned regionIdx)
Convert a gas mass fraction in the water phase the corresponding mole fraction.
Definition: BlackOilFluidSystem.hpp:1360
static LhsEval saturatedInverseFormationVolumeFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the formation volume factor of a "saturated" fluid phase.
Definition: BlackOilFluidSystem.hpp:824
static void initBegin(std::size_t numPvtRegions)
Begin the initialization of the black oil fluid system.
static LhsEval convertXoGToxoG(const LhsEval &XoG, unsigned regionIdx)
Convert a gas mass fraction in the oil phase the corresponding mole fraction.
Definition: BlackOilFluidSystem.hpp:1372
static LhsEval viscosity(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Calculate the dynamic viscosity of a fluid phase [Pa*s].
Definition: BlackOilFluidSystem.hpp:968
static Scalar reservoirTemperature(unsigned=0)
Set the temperature of the reservoir.
Definition: BlackOilFluidSystem.hpp:1451
static unsigned soluteComponentIndex(unsigned phaseIdx)
returns the index of "secondary" component of a phase (solute)
static LhsEval convertRswToXwG(const LhsEval &Rsw, unsigned regionIdx)
Convert a gas dissolution factor to the the corresponding mass fraction of the gas component in the w...
Definition: BlackOilFluidSystem.hpp:1307
static void setEnableDiffusion(bool yesno)
Specify whether the fluid system should consider diffusion.
Definition: BlackOilFluidSystem.hpp:286
static LhsEval bubblePointPressure(const FluidState &fluidState, unsigned regionIdx)
Returns the bubble point pressure $P_b$ using the current Rs.
Definition: BlackOilFluidSystem.hpp:1185
static LhsEval density(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Calculate the density [kg/m^3] of a fluid phase.
Definition: BlackOilFluidSystem.hpp:525
static LhsEval saturatedDensity(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Compute the density of a saturated fluid phase.
Definition: BlackOilFluidSystem.hpp:623
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:507
static bool isIdealMixture(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal mixture.
Definition: BlackOilFluidSystem.hpp:402
static bool enableVaporizedWater()
Returns whether the fluid system should consider that the water component can dissolve in the gas pha...
Definition: BlackOilFluidSystem.hpp:463
static unsigned solventComponentIndex(unsigned phaseIdx)
returns the index of "primary" component of a phase (solvent)
static constexpr unsigned numPhases
Number of fluid phases in the fluid system.
Definition: BlackOilFluidSystem.hpp:333
static void setWaterPvt(std::shared_ptr< WaterPvt > pvtObj)
Set the pressure-volume-saturation (PVT) relations for the water phase.
Definition: BlackOilFluidSystem.hpp:305
static void setReservoirTemperature(Scalar value)
Return the temperature of the reservoir.
Definition: BlackOilFluidSystem.hpp:1459
static void setEnableDissolvedGasInWater(bool yesno)
Specify whether the fluid system should consider that the gas component can dissolve in the water pha...
Definition: BlackOilFluidSystem.hpp:279
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:1130
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:1061
static constexpr unsigned oilCompIdx
Index of the oil component.
Definition: BlackOilFluidSystem.hpp:366
static const char * componentName(unsigned compIdx)
Return the human readable name of a component.
static bool isLiquid(unsigned phaseIdx)
Return whether a phase is liquid.
Definition: BlackOilFluidSystem.hpp:352
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:514
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:846
static LhsEval convertxoGToXoG(const LhsEval &xoG, unsigned regionIdx)
Convert a gas mole fraction in the oil phase the corresponding mass fraction.
Definition: BlackOilFluidSystem.hpp:1384
static Scalar molarMass(unsigned compIdx, unsigned regionIdx=0)
Return the molar mass of a component in [kg/mol].
Definition: BlackOilFluidSystem.hpp:398
static LhsEval saturatedDissolutionFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the dissolution factor of a saturated fluid phase.
Definition: BlackOilFluidSystem.hpp:1161
static LhsEval dewPointPressure(const FluidState &fluidState, unsigned regionIdx)
Returns the dew point pressure $P_d$ using the current Rv.
Definition: BlackOilFluidSystem.hpp:1196
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:1213
static LhsEval convertXwGToRsw(const LhsEval &XwG, unsigned regionIdx)
Convert the mass fraction of the gas component in the water phase to the corresponding gas dissolutio...
Definition: BlackOilFluidSystem.hpp:1253
static LhsEval convertXgOToxgO(const LhsEval &XgO, unsigned regionIdx)
Convert a oil mass fraction in the gas phase the corresponding mole fraction.
Definition: BlackOilFluidSystem.hpp:1396
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:1478
static LhsEval saturatedVaporizationFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the water vaporization factor of saturated phase.
Definition: BlackOilFluidSystem.hpp:1104
static bool enableDissolvedGas()
Returns whether the fluid system should consider that the gas component can dissolve in the oil phase...
Definition: BlackOilFluidSystem.hpp:435
static Scalar referenceDensity(unsigned phaseIdx, unsigned regionIdx)
Returns the density of a fluid phase at surface pressure [kg/m^3].
Definition: BlackOilFluidSystem.hpp:479
static Scalar surfaceTemperature
The temperature at the surface.
Definition: BlackOilFluidSystem.hpp:346
static void setOilPvt(std::shared_ptr< OilPvt > pvtObj)
Set the pressure-volume-saturation (PVT) relations for the oil phase.
Definition: BlackOilFluidSystem.hpp:299
static constexpr unsigned waterPhaseIdx
Index of the water phase.
Definition: BlackOilFluidSystem.hpp:336
static LhsEval convertxgOToXgO(const LhsEval &xgO, unsigned regionIdx)
Convert a oil mole fraction in the gas phase the corresponding mass fraction.
Definition: BlackOilFluidSystem.hpp:1408
static LhsEval inverseFormationVolumeFactor(const FluidState &fluidState, unsigned phaseIdx, unsigned regionIdx)
Returns the formation volume factor of an "undersaturated" fluid phase.
Definition: BlackOilFluidSystem.hpp:724
static LhsEval convertXgWToxgW(const LhsEval &XgW, unsigned regionIdx)
Convert a water mass fraction in the gas phase the corresponding mole fraction.
Definition: BlackOilFluidSystem.hpp:1348
static constexpr unsigned gasPhaseIdx
Index of the gas phase.
Definition: BlackOilFluidSystem.hpp:340
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:252
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:1266
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:494
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:487
static constexpr unsigned gasCompIdx
Index of the gas component.
Definition: BlackOilFluidSystem.hpp:370
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:1321
static const WaterPvt & waterPvt()
Return a reference to the low-level object which calculates the water phase quantities.
Definition: BlackOilFluidSystem.hpp:1443
static bool enableVaporizedOil()
Returns whether the fluid system should consider that the oil component can dissolve in the gas phase...
Definition: BlackOilFluidSystem.hpp:454
static bool isIdealGas(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal gas.
Definition: BlackOilFluidSystem.hpp:414
static unsigned phaseIsActive(unsigned phaseIdx)
Returns whether a fluid phase is active.
Definition: BlackOilFluidSystem.hpp:382
static const char * phaseName(unsigned phaseIdx)
Return the human readable name of a fluid phase.
static Scalar surfacePressure
The pressure at the surface.
Definition: BlackOilFluidSystem.hpp:343
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:270
static constexpr unsigned oilPhaseIdx
Index of the oil phase.
Definition: BlackOilFluidSystem.hpp:338
static const OilPvt & oilPvt()
Return a reference to the low-level object which calculates the oil phase quantities.
Definition: BlackOilFluidSystem.hpp:1433
static bool enableDissolvedGasInWater()
Returns whether the fluid system should consider that the gas component can dissolve in the water pha...
Definition: BlackOilFluidSystem.hpp:445
static bool enableDiffusion()
Returns whether the fluid system should consider diffusion.
Definition: BlackOilFluidSystem.hpp:471
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:1279
static std::size_t numRegions()
Returns the number of PVT regions which are considered.
Definition: BlackOilFluidSystem.hpp:426
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:261
static void setDiffusionCoefficient(Scalar coefficient, unsigned compIdx, unsigned phaseIdx, unsigned regionIdx=0)
Definition: BlackOilFluidSystem.hpp:1471
static const GasPvt & gasPvt()
Return a reference to the low-level object which calculates the gas phase quantities.
Definition: BlackOilFluidSystem.hpp:1423
static bool isCompressible(unsigned)
Returns true if and only if a fluid phase is assumed to be compressible.
Definition: BlackOilFluidSystem.hpp:410
static void setGasPvt(std::shared_ptr< GasPvt > pvtObj)
Set the pressure-volume-saturation (PVT) relations for the gas phase.
Definition: BlackOilFluidSystem.hpp:293
This class represents the Pressure-Volume-Temperature relations of the gas phase in the black-oil mod...
Definition: GasPvtMultiplexer.hpp:102
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:316
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:97
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:254
This class represents the Pressure-Volume-Temperature relations of the water phase in the black-oil m...
Definition: WaterPvtMultiplexer.hpp:82
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:245
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition: Exceptions.hpp:30
The type of the fluid system's parameter cache.
Definition: BlackOilFluidSystem.hpp:171
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:189
void setRegionIndex(unsigned val)
Set the index of the region which should be used to determine the thermodynamic properties.
Definition: BlackOilFluidSystem.hpp:212
unsigned regionIndex() const
Return the index of the region which should be used to determine the thermodynamic properties.
Definition: BlackOilFluidSystem.hpp:202