27 #ifndef OPM_WATER_PVT_MULTIPLEXER_HPP 28 #define OPM_WATER_PVT_MULTIPLEXER_HPP 33 #define OPM_WATER_PVT_MULTIPLEXER_CALL(codeToCall) \ 34 switch (approach_) { \ 35 case ConstantCompressibilityWaterPvt: { \ 36 auto& pvtImpl = getRealPvt<ConstantCompressibilityWaterPvt>(); \ 40 case ThermalWaterPvt: { \ 41 auto& pvtImpl = getRealPvt<ThermalWaterPvt>(); \ 46 OPM_THROW(std::logic_error, "Not implemented: Water PVT of this deck!"); \ 54 template <
class Scalar,
bool enableThermal = true>
60 enum WaterPvtApproach {
68 approach_ = NoWaterPvt;
75 delete &getRealPvt<ConstantCompressibilityWaterPvt>();
78 case ThermalWaterPvt: {
79 delete &getRealPvt<ThermalWaterPvt>();
93 void initFromDeck(
const Deck& deck,
const EclipseState& eclState)
95 bool enableWater = deck.hasKeyword(
"WATER");
100 && (deck.hasKeyword(
"WATDENT")
101 || deck.hasKeyword(
"VISCREF")))
102 setApproach(ThermalWaterPvt);
103 else if (deck.hasKeyword(
"PVTW"))
106 OPM_WATER_PVT_MULTIPLEXER_CALL(pvtImpl.initFromDeck(deck, eclState));
108 #endif // HAVE_OPM_PARSER 111 { OPM_WATER_PVT_MULTIPLEXER_CALL(pvtImpl.initEnd()); }
117 { OPM_WATER_PVT_MULTIPLEXER_CALL(
return pvtImpl.numRegions());
return 1; }
122 template <
class Evaluation>
124 const Evaluation& temperature,
125 const Evaluation& pressure)
const 126 { OPM_WATER_PVT_MULTIPLEXER_CALL(
return pvtImpl.viscosity(regionIdx, temperature, pressure));
return 0; }
131 template <
class Evaluation>
133 const Evaluation& temperature,
134 const Evaluation& pressure)
const 135 { OPM_WATER_PVT_MULTIPLEXER_CALL(
return pvtImpl.inverseFormationVolumeFactor(regionIdx, temperature, pressure));
return 0; }
137 void setApproach(WaterPvtApproach appr)
144 case ThermalWaterPvt:
149 OPM_THROW(std::logic_error,
"Not implemented: Water PVT of this deck!");
161 {
return approach_; }
164 template <WaterPvtApproach approachV>
165 typename std::enable_if<approachV == ConstantCompressibilityWaterPvt, Opm::ConstantCompressibilityWaterPvt<Scalar> >::type& getRealPvt()
171 template <WaterPvtApproach approachV>
172 typename std::enable_if<approachV == ConstantCompressibilityWaterPvt, const Opm::ConstantCompressibilityWaterPvt<Scalar> >::type& getRealPvt()
const 178 template <WaterPvtApproach approachV>
179 typename std::enable_if<approachV == ThermalWaterPvt, Opm::WaterPvtThermal<Scalar> >::type& getRealPvt()
185 template <WaterPvtApproach approachV>
186 typename std::enable_if<approachV == ThermalWaterPvt, const Opm::WaterPvtThermal<Scalar> >::type& getRealPvt()
const 193 WaterPvtApproach approach_;
197 #undef OPM_WATER_PVT_MULTIPLEXER_CALL This class implements temperature dependence of the PVT properties of water.
Definition: WaterPvtThermal.hpp:55
Definition: Air_Mesitylene.hpp:33
This class represents the Pressure-Volume-Temperature relations of the gas phase without vaporized oi...
This class represents the Pressure-Volume-Temperature relations of the gas phase without vaporized oi...
Definition: ConstantCompressibilityWaterPvt.hpp:48
Evaluation viscosity(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure) const
Returns the dynamic viscosity [Pa s] of the fluid phase given a set of parameters.
Definition: WaterPvtMultiplexer.hpp:123
This class represents the Pressure-Volume-Temperature relations of the water phase in the black-oil m...
Definition: WaterPvtMultiplexer.hpp:55
WaterPvtApproach approach() const
Returns the concrete approach for calculating the PVT relations.
Definition: WaterPvtMultiplexer.hpp:160
This class implements temperature dependence of the PVT properties of water.
unsigned numRegions() const
Return the number of PVT regions which are considered by this PVT-object.
Definition: WaterPvtMultiplexer.hpp:116
Evaluation inverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure) const
Returns the formation volume factor [-] of the fluid phase.
Definition: WaterPvtMultiplexer.hpp:132