27 #ifndef OPM_OIL_PVT_MULTIPLEXER_HPP
28 #define OPM_OIL_PVT_MULTIPLEXER_HPP
37 #include <opm/input/eclipse/EclipseState/EclipseState.hpp>
38 #include <opm/input/eclipse/EclipseState/Runspec.hpp>
42 #define OPM_OIL_PVT_MULTIPLEXER_CALL(codeToCall) \
43 switch (approach_) { \
44 case OilPvtApproach::ConstantCompressibilityOilPvt: { \
45 auto& pvtImpl = getRealPvt<OilPvtApproach::ConstantCompressibilityOilPvt>(); \
49 case OilPvtApproach::DeadOilPvt: { \
50 auto& pvtImpl = getRealPvt<OilPvtApproach::DeadOilPvt>(); \
54 case OilPvtApproach::LiveOilPvt: { \
55 auto& pvtImpl = getRealPvt<OilPvtApproach::LiveOilPvt>(); \
59 case OilPvtApproach::ThermalOilPvt: { \
60 auto& pvtImpl = getRealPvt<OilPvtApproach::ThermalOilPvt>(); \
64 case OilPvtApproach::BrineCo2Pvt: { \
65 auto& pvtImpl = getRealPvt<OilPvtApproach::BrineCo2Pvt>(); \
69 case OilPvtApproach::NoOilPvt: \
70 throw std::logic_error("Not implemented: Oil PVT of this deck!"); \
73 enum class OilPvtApproach {
77 ConstantCompressibilityOilPvt,
94 template <
class Scalar,
bool enableThermal = true>
100 approach_ = OilPvtApproach::NoOilPvt;
101 realOilPvt_ =
nullptr;
106 , realOilPvt_(realOilPvt)
117 case OilPvtApproach::LiveOilPvt: {
118 delete &getRealPvt<OilPvtApproach::LiveOilPvt>();
121 case OilPvtApproach::DeadOilPvt: {
122 delete &getRealPvt<OilPvtApproach::DeadOilPvt>();
125 case OilPvtApproach::ConstantCompressibilityOilPvt: {
126 delete &getRealPvt<OilPvtApproach::ConstantCompressibilityOilPvt>();
129 case OilPvtApproach::ThermalOilPvt: {
130 delete &getRealPvt<OilPvtApproach::ThermalOilPvt>();
133 case OilPvtApproach::BrineCo2Pvt: {
134 delete &getRealPvt<OilPvtApproach::BrineCo2Pvt>();
138 case OilPvtApproach::NoOilPvt:
149 void initFromState(
const EclipseState& eclState,
const Schedule& schedule)
151 if (!eclState.runspec().phases().active(Phase::OIL))
155 if (eclState.runspec().co2Storage())
156 setApproach(OilPvtApproach::BrineCo2Pvt);
157 else if (enableThermal && eclState.getSimulationConfig().isThermal())
158 setApproach(OilPvtApproach::ThermalOilPvt);
159 else if (!eclState.getTableManager().getPvcdoTable().empty())
160 setApproach(OilPvtApproach::ConstantCompressibilityOilPvt);
161 else if (eclState.getTableManager().hasTables(
"PVDO"))
162 setApproach(OilPvtApproach::DeadOilPvt);
163 else if (!eclState.getTableManager().getPvtoTables().empty())
164 setApproach(OilPvtApproach::LiveOilPvt);
166 OPM_OIL_PVT_MULTIPLEXER_CALL(pvtImpl.initFromState(eclState, schedule));
172 { OPM_OIL_PVT_MULTIPLEXER_CALL(pvtImpl.initEnd()); }
178 { OPM_OIL_PVT_MULTIPLEXER_CALL(
return pvtImpl.numRegions());
return 1; }
184 { OPM_OIL_PVT_MULTIPLEXER_CALL(
return pvtImpl.oilReferenceDensity(regionIdx));
return 700.; }
189 template <
class Evaluation>
191 const Evaluation& temperature,
192 const Evaluation& pressure,
193 const Evaluation& Rs)
const
194 { OPM_OIL_PVT_MULTIPLEXER_CALL(
return pvtImpl.internalEnergy(regionIdx, temperature, pressure, Rs));
return 0; }
199 template <
class Evaluation>
201 const Evaluation& temperature,
202 const Evaluation& pressure,
203 const Evaluation& Rs)
const
204 { OPM_OIL_PVT_MULTIPLEXER_CALL(
return pvtImpl.viscosity(regionIdx, temperature, pressure, Rs));
return 0; }
209 template <
class Evaluation>
211 const Evaluation& temperature,
212 const Evaluation& pressure)
const
213 { OPM_OIL_PVT_MULTIPLEXER_CALL(
return pvtImpl.saturatedViscosity(regionIdx, temperature, pressure));
return 0; }
218 template <
class Evaluation>
220 const Evaluation& temperature,
221 const Evaluation& pressure,
222 const Evaluation& Rs)
const
223 { OPM_OIL_PVT_MULTIPLEXER_CALL(
return pvtImpl.inverseFormationVolumeFactor(regionIdx, temperature, pressure, Rs));
return 0; }
228 template <
class Evaluation>
230 const Evaluation& temperature,
231 const Evaluation& pressure)
const
232 { OPM_OIL_PVT_MULTIPLEXER_CALL(
return pvtImpl.saturatedInverseFormationVolumeFactor(regionIdx, temperature, pressure));
return 0; }
237 template <
class Evaluation>
239 const Evaluation& temperature,
240 const Evaluation& pressure)
const
241 { OPM_OIL_PVT_MULTIPLEXER_CALL(
return pvtImpl.saturatedGasDissolutionFactor(regionIdx, temperature, pressure));
return 0; }
246 template <
class Evaluation>
248 const Evaluation& temperature,
249 const Evaluation& pressure,
250 const Evaluation& oilSaturation,
251 const Evaluation& maxOilSaturation)
const
252 { OPM_OIL_PVT_MULTIPLEXER_CALL(
return pvtImpl.saturatedGasDissolutionFactor(regionIdx, temperature, pressure, oilSaturation, maxOilSaturation));
return 0; }
261 template <
class Evaluation>
263 const Evaluation& temperature,
264 const Evaluation& Rs)
const
265 { OPM_OIL_PVT_MULTIPLEXER_CALL(
return pvtImpl.saturationPressure(regionIdx, temperature, Rs));
return 0; }
270 template <
class Evaluation>
272 const Evaluation& pressure,
273 unsigned compIdx)
const
275 OPM_OIL_PVT_MULTIPLEXER_CALL(
return pvtImpl.diffusionCoefficient(temperature, pressure, compIdx));
return 0;
278 void setApproach(OilPvtApproach appr)
281 case OilPvtApproach::LiveOilPvt:
285 case OilPvtApproach::DeadOilPvt:
289 case OilPvtApproach::ConstantCompressibilityOilPvt:
293 case OilPvtApproach::ThermalOilPvt:
297 case OilPvtApproach::BrineCo2Pvt:
301 case OilPvtApproach::NoOilPvt:
302 throw std::logic_error(
"Not implemented: Oil PVT of this deck!");
314 {
return approach_; }
317 template <OilPvtApproach approachV>
318 typename std::enable_if<approachV == OilPvtApproach::LiveOilPvt, LiveOilPvt<Scalar> >::type& getRealPvt()
324 template <OilPvtApproach approachV>
325 typename std::enable_if<approachV == OilPvtApproach::LiveOilPvt, const LiveOilPvt<Scalar> >::type& getRealPvt()
const
328 return *
static_cast<LiveOilPvt<Scalar>*
>(realOilPvt_);
331 template <OilPvtApproach approachV>
332 typename std::enable_if<approachV == OilPvtApproach::DeadOilPvt, DeadOilPvt<Scalar> >::type& getRealPvt()
335 return *
static_cast<DeadOilPvt<Scalar>*
>(realOilPvt_);
338 template <OilPvtApproach approachV>
339 typename std::enable_if<approachV == OilPvtApproach::DeadOilPvt, const DeadOilPvt<Scalar> >::type& getRealPvt()
const
342 return *
static_cast<DeadOilPvt<Scalar>*
>(realOilPvt_);
345 template <OilPvtApproach approachV>
346 typename std::enable_if<approachV == OilPvtApproach::ConstantCompressibilityOilPvt, ConstantCompressibilityOilPvt<Scalar> >::type& getRealPvt()
349 return *
static_cast<ConstantCompressibilityOilPvt<Scalar>*
>(realOilPvt_);
352 template <OilPvtApproach approachV>
353 typename std::enable_if<approachV == OilPvtApproach::ConstantCompressibilityOilPvt, const ConstantCompressibilityOilPvt<Scalar> >::type& getRealPvt()
const
356 return *
static_cast<ConstantCompressibilityOilPvt<Scalar>*
>(realOilPvt_);
359 template <OilPvtApproach approachV>
360 typename std::enable_if<approachV == OilPvtApproach::ThermalOilPvt, OilPvtThermal<Scalar> >::type& getRealPvt()
363 return *
static_cast<OilPvtThermal<Scalar>*
>(realOilPvt_);
366 template <OilPvtApproach approachV>
367 typename std::enable_if<approachV == OilPvtApproach::ThermalOilPvt, const OilPvtThermal<Scalar> >::type& getRealPvt()
const
370 return *
static_cast<const OilPvtThermal<Scalar>*
>(realOilPvt_);
373 template <OilPvtApproach approachV>
374 typename std::enable_if<approachV == OilPvtApproach::BrineCo2Pvt, BrineCo2Pvt<Scalar> >::type& getRealPvt()
377 return *
static_cast<BrineCo2Pvt<Scalar>*
>(realOilPvt_);
380 template <OilPvtApproach approachV>
381 typename std::enable_if<approachV == OilPvtApproach::BrineCo2Pvt, const BrineCo2Pvt<Scalar> >::type& getRealPvt()
const
384 return *
static_cast<const BrineCo2Pvt<Scalar>*
>(realOilPvt_);
387 const void* realOilPvt()
const {
return realOilPvt_; }
389 bool operator==(
const OilPvtMultiplexer<Scalar,enableThermal>& data)
const
391 if (this->
approach() != data.approach())
395 case OilPvtApproach::ConstantCompressibilityOilPvt:
396 return *
static_cast<const ConstantCompressibilityOilPvt<Scalar>*
>(realOilPvt_) ==
397 *
static_cast<const ConstantCompressibilityOilPvt<Scalar>*
>(data.realOilPvt_);
398 case OilPvtApproach::DeadOilPvt:
399 return *
static_cast<const DeadOilPvt<Scalar>*
>(realOilPvt_) ==
400 *
static_cast<const DeadOilPvt<Scalar>*
>(data.realOilPvt_);
401 case OilPvtApproach::LiveOilPvt:
402 return *
static_cast<const LiveOilPvt<Scalar>*
>(realOilPvt_) ==
403 *
static_cast<const LiveOilPvt<Scalar>*
>(data.realOilPvt_);
404 case OilPvtApproach::ThermalOilPvt:
405 return *
static_cast<const OilPvtThermal<Scalar>*
>(realOilPvt_) ==
406 *
static_cast<const OilPvtThermal<Scalar>*
>(data.realOilPvt_);
407 case OilPvtApproach::BrineCo2Pvt:
408 return *
static_cast<const BrineCo2Pvt<Scalar>*
>(realOilPvt_) ==
409 *
static_cast<const BrineCo2Pvt<Scalar>*
>(data.realOilPvt_);
415 OilPvtMultiplexer<Scalar,enableThermal>& operator=(
const OilPvtMultiplexer<Scalar,enableThermal>& data)
417 approach_ = data.approach_;
419 case OilPvtApproach::ConstantCompressibilityOilPvt:
420 realOilPvt_ =
new ConstantCompressibilityOilPvt<Scalar>(*
static_cast<const ConstantCompressibilityOilPvt<Scalar>*
>(data.realOilPvt_));
422 case OilPvtApproach::DeadOilPvt:
423 realOilPvt_ =
new DeadOilPvt<Scalar>(*
static_cast<const DeadOilPvt<Scalar>*
>(data.realOilPvt_));
425 case OilPvtApproach::LiveOilPvt:
426 realOilPvt_ =
new LiveOilPvt<Scalar>(*
static_cast<const LiveOilPvt<Scalar>*
>(data.realOilPvt_));
428 case OilPvtApproach::ThermalOilPvt:
429 realOilPvt_ =
new OilPvtThermal<Scalar>(*
static_cast<const OilPvtThermal<Scalar>*
>(data.realOilPvt_));
431 case OilPvtApproach::BrineCo2Pvt:
432 realOilPvt_ =
new BrineCo2Pvt<Scalar>(*
static_cast<const BrineCo2Pvt<Scalar>*
>(data.realOilPvt_));
442 OilPvtApproach approach_;
446 #undef OPM_OIL_PVT_MULTIPLEXER_CALL
This class represents the Pressure-Volume-Temperature relations of the liquid phase for a CO2-Brine s...
This class represents the Pressure-Volume-Temperature relations of the oil phase without dissolved ga...
This class represents the Pressure-Volume-Temperature relations of the oil phase without dissolved ga...
This class represents the Pressure-Volume-Temperature relations of the oil phas with dissolved gas.
This class implements temperature dependence of the PVT properties of oil.
This class represents the Pressure-Volume-Temperature relations of the liquid phase for a CO2-Brine s...
Definition: BrineCo2Pvt.hpp:57
This class represents the Pressure-Volume-Temperature relations of the oil phase without dissolved ga...
Definition: ConstantCompressibilityOilPvt.hpp:42
This class represents the Pressure-Volume-Temperature relations of the oil phase without dissolved ga...
Definition: DeadOilPvt.hpp:45
This class represents the Pressure-Volume-Temperature relations of the oil phas with dissolved gas.
Definition: LiveOilPvt.hpp:52
This class represents the Pressure-Volume-Temperature relations of the oil phase in the black-oil mod...
Definition: OilPvtMultiplexer.hpp:96
unsigned numRegions() const
Return the number of PVT regions which are considered by this PVT-object.
Definition: OilPvtMultiplexer.hpp:177
Evaluation internalEnergy(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure, const Evaluation &Rs) const
Returns the specific enthalpy [J/kg] oil given a set of parameters.
Definition: OilPvtMultiplexer.hpp:190
Evaluation inverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure, const Evaluation &Rs) const
Returns the formation volume factor [-] of the fluid phase.
Definition: OilPvtMultiplexer.hpp:219
Evaluation saturatedViscosity(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: OilPvtMultiplexer.hpp:210
OilPvtApproach approach() const
Returns the concrete approach for calculating the PVT relations.
Definition: OilPvtMultiplexer.hpp:313
Evaluation saturatedInverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure) const
Returns the formation volume factor [-] of the fluid phase.
Definition: OilPvtMultiplexer.hpp:229
Evaluation saturatedGasDissolutionFactor(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure) const
Returns the gas dissolution factor [m^3/m^3] of saturated oil.
Definition: OilPvtMultiplexer.hpp:238
Evaluation viscosity(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure, const Evaluation &Rs) const
Returns the dynamic viscosity [Pa s] of the fluid phase given a set of parameters.
Definition: OilPvtMultiplexer.hpp:200
Evaluation saturationPressure(unsigned regionIdx, const Evaluation &temperature, const Evaluation &Rs) const
Returns the saturation pressure [Pa] of oil given the mass fraction of the gas component in the oil p...
Definition: OilPvtMultiplexer.hpp:262
Evaluation saturatedGasDissolutionFactor(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure, const Evaluation &oilSaturation, const Evaluation &maxOilSaturation) const
Returns the gas dissolution factor [m^3/m^3] of saturated oil.
Definition: OilPvtMultiplexer.hpp:247
const Scalar oilReferenceDensity(unsigned regionIdx) const
Return the reference density which are considered by this PVT-object.
Definition: OilPvtMultiplexer.hpp:183
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 implements temperature dependence of the PVT properties of oil.
Definition: OilPvtThermal.hpp:50