28 #ifndef OPM_BLACK_OIL_FLUID_STATE_HH
29 #define OPM_BLACK_OIL_FLUID_STATE_HH
42 template <class FluidState>
43 unsigned getPvtRegionIndex_(typename std::enable_if<HasMember_pvtRegionIndex<FluidState>::value,
44 const FluidState&>::type fluidState)
45 {
return fluidState.pvtRegionIndex(); }
47 template <
class Flu
idState>
48 unsigned getPvtRegionIndex_(
typename std::enable_if<!HasMember_pvtRegionIndex<FluidState>::value,
49 const FluidState&>::type)
54 template <class FluidSystem, class FluidState, class LhsEval>
55 auto getInvB_(typename std::enable_if<HasMember_invB<FluidState>::value,
56 const FluidState&>::type fluidState,
59 -> decltype(decay<LhsEval>(fluidState.invB(phaseIdx)))
60 {
return decay<LhsEval>(fluidState.invB(phaseIdx)); }
62 template <
class Flu
idSystem,
class Flu
idState,
class LhsEval>
63 LhsEval getInvB_(
typename std::enable_if<!HasMember_invB<FluidState>::value,
64 const FluidState&>::type fluidState,
66 unsigned pvtRegionIdx)
68 const auto& rho = fluidState.density(phaseIdx);
69 const auto& Xsolvent =
70 fluidState.massFraction(phaseIdx, FluidSystem::solventComponentIndex(phaseIdx));
74 *decay<LhsEval>(Xsolvent)
75 /FluidSystem::referenceDensity(phaseIdx, pvtRegionIdx);
80 template <class FluidState>
81 auto getSaltConcentration_(typename std::enable_if<HasMember_saltConcentration<FluidState>::value,
82 const FluidState&>::type fluidState)
83 {
return fluidState.saltConcentration(); }
85 template <
class Flu
idState>
86 auto getSaltConcentration_(
typename std::enable_if<!HasMember_saltConcentration<FluidState>::value,
87 const FluidState&>::type)
97 template <
class ScalarT,
99 bool enableTemperature =
false,
100 bool enableEnergy =
false,
101 bool enableDissolution =
true,
102 bool enableBrine =
false,
103 unsigned numStoragePhases = FluidSystem::numPhases>
106 enum { waterPhaseIdx = FluidSystem::waterPhaseIdx };
107 enum { gasPhaseIdx = FluidSystem::gasPhaseIdx };
108 enum { oilPhaseIdx = FluidSystem::oilPhaseIdx };
110 enum { waterCompIdx = FluidSystem::waterCompIdx };
111 enum { gasCompIdx = FluidSystem::gasCompIdx };
112 enum { oilCompIdx = FluidSystem::oilCompIdx };
115 typedef ScalarT Scalar;
116 enum { numPhases = FluidSystem::numPhases };
117 enum { numComponents = FluidSystem::numComponents };
130 Valgrind::CheckDefined(pvtRegionIdx_);
132 for (
unsigned storagePhaseIdx = 0; storagePhaseIdx < numStoragePhases; ++ storagePhaseIdx) {
133 Valgrind::CheckDefined(saturation_[storagePhaseIdx]);
134 Valgrind::CheckDefined(pressure_[storagePhaseIdx]);
135 Valgrind::CheckDefined(density_[storagePhaseIdx]);
136 Valgrind::CheckDefined(invB_[storagePhaseIdx]);
139 Valgrind::CheckDefined((*enthalpy_)[storagePhaseIdx]);
142 if (enableDissolution) {
143 Valgrind::CheckDefined(*Rs_);
144 Valgrind::CheckDefined(*Rv_);
148 Valgrind::CheckDefined(*saltConcentration_);
151 if (enableTemperature || enableEnergy)
152 Valgrind::CheckDefined(*temperature_);
160 template <
class Flu
idState>
163 if (enableTemperature || enableEnergy)
166 unsigned pvtRegionIdx = getPvtRegionIndex_<FluidState>(fs);
169 if (enableDissolution) {
170 setRs(BlackOil::getRs_<FluidSystem, FluidState, Scalar>(fs, pvtRegionIdx));
171 setRv(BlackOil::getRv_<FluidSystem, FluidState, Scalar>(fs, pvtRegionIdx));
175 setSaltConcentration(BlackOil::getSaltConcentration_<FluidSystem, FluidState, Scalar>(fs, pvtRegionIdx));
177 for (
unsigned storagePhaseIdx = 0; storagePhaseIdx < numStoragePhases; ++storagePhaseIdx) {
178 unsigned phaseIdx = storageToCanonicalPhaseIndex_(storagePhaseIdx);
186 setInvB(phaseIdx, getInvB_<FluidSystem, FluidState, Scalar>(fs, phaseIdx, pvtRegionIdx));
197 { pvtRegionIdx_ =
static_cast<unsigned short>(newPvtRegionIdx); }
203 { pressure_[canonicalToStoragePhaseIndex_(phaseIdx)] = p; }
209 { saturation_[canonicalToStoragePhaseIndex_(phaseIdx)] = S; }
214 void setPc(
unsigned phaseIdx,
const Scalar&
pc)
215 { pc_[canonicalToStoragePhaseIndex_(phaseIdx)] =
pc; }
222 totalSaturation_ = value;
233 assert(enableTemperature || enableEnergy);
235 (*temperature_) = value;
246 assert(enableTemperature || enableEnergy);
248 (*enthalpy_)[canonicalToStoragePhaseIndex_(phaseIdx)] = value;
254 void setInvB(
unsigned phaseIdx,
const Scalar& b)
255 { invB_[canonicalToStoragePhaseIndex_(phaseIdx)] = b; }
261 { density_[canonicalToStoragePhaseIndex_(phaseIdx)] = rho; }
283 { *saltConcentration_ = newSaltConcentration; }
289 {
return pressure_[canonicalToStoragePhaseIndex_(phaseIdx)]; }
295 {
return saturation_[canonicalToStoragePhaseIndex_(phaseIdx)]; }
300 const Scalar&
pc(
unsigned phaseIdx)
const
301 {
return pc_[canonicalToStoragePhaseIndex_(phaseIdx)]; }
308 return totalSaturation_;
316 if (!enableTemperature && !enableEnergy) {
317 static Scalar tmp(FluidSystem::reservoirTemperature(pvtRegionIdx_));
321 return *temperature_;
330 const Scalar&
invB(
unsigned phaseIdx)
const
331 {
return invB_[canonicalToStoragePhaseIndex_(phaseIdx)]; }
340 const Scalar&
Rs()
const
342 if (!enableDissolution) {
343 static Scalar
null = 0.0;
357 const Scalar&
Rv()
const
359 if (!enableDissolution) {
360 static Scalar
null = 0.0;
373 static Scalar
null = 0.0;
377 return *saltConcentration_;
389 {
return pvtRegionIdx_; }
395 {
return density_[canonicalToStoragePhaseIndex_(phaseIdx)]; }
404 {
return (*enthalpy_)[canonicalToStoragePhaseIndex_(phaseIdx)]; }
413 {
return (*enthalpy_)[canonicalToStoragePhaseIndex_(phaseIdx)] -
pressure(phaseIdx)/
density(phaseIdx); }
424 const auto& rho =
density(phaseIdx);
426 if (phaseIdx == waterPhaseIdx)
427 return rho/FluidSystem::molarMass(waterCompIdx, pvtRegionIdx_);
430 rho*(
moleFraction(phaseIdx, gasCompIdx)/FluidSystem::molarMass(gasCompIdx, pvtRegionIdx_)
431 +
moleFraction(phaseIdx, oilCompIdx)/FluidSystem::molarMass(oilCompIdx, pvtRegionIdx_));
447 {
return FluidSystem::viscosity(*
this, phaseIdx, pvtRegionIdx_); }
456 if (compIdx == waterCompIdx)
461 if (compIdx == waterCompIdx)
463 else if (compIdx == oilCompIdx)
464 return 1.0 - FluidSystem::convertRsToXoG(
Rs(), pvtRegionIdx_);
466 assert(compIdx == gasCompIdx);
467 return FluidSystem::convertRsToXoG(
Rs(), pvtRegionIdx_);
472 if (compIdx == waterCompIdx)
474 else if (compIdx == oilCompIdx)
475 return FluidSystem::convertRvToXgO(
Rv(), pvtRegionIdx_);
477 assert(compIdx == gasCompIdx);
478 return 1.0 - FluidSystem::convertRvToXgO(
Rv(), pvtRegionIdx_);
483 throw std::logic_error(
"Invalid phase or component index!");
493 if (compIdx == waterCompIdx)
498 if (compIdx == waterCompIdx)
500 else if (compIdx == oilCompIdx)
501 return 1.0 - FluidSystem::convertXoGToxoG(FluidSystem::convertRsToXoG(
Rs(), pvtRegionIdx_),
504 assert(compIdx == gasCompIdx);
505 return FluidSystem::convertXoGToxoG(FluidSystem::convertRsToXoG(
Rs(), pvtRegionIdx_),
511 if (compIdx == waterCompIdx)
513 else if (compIdx == oilCompIdx)
514 return FluidSystem::convertXgOToxgO(FluidSystem::convertRvToXgO(
Rv(), pvtRegionIdx_),
517 assert(compIdx == gasCompIdx);
518 return 1.0 - FluidSystem::convertXgOToxgO(FluidSystem::convertRvToXgO(
Rv(), pvtRegionIdx_),
524 throw std::logic_error(
"Invalid phase or component index!");
530 Scalar
molarity(
unsigned phaseIdx,
unsigned compIdx)
const
539 for (
unsigned compIdx = 0; compIdx < numComponents; ++ compIdx)
540 result += FluidSystem::molarMass(compIdx, pvtRegionIdx_)*
moleFraction(phaseIdx, compIdx);
548 {
return FluidSystem::fugacityCoefficient(*
this, phaseIdx, compIdx, pvtRegionIdx_); }
553 Scalar
fugacity(
unsigned phaseIdx,
unsigned compIdx)
const
562 static unsigned storageToCanonicalPhaseIndex_(
unsigned storagePhaseIdx)
564 if (numStoragePhases == 3)
565 return storagePhaseIdx;
567 return FluidSystem::activeToCanonicalPhaseIdx(storagePhaseIdx);
570 static unsigned canonicalToStoragePhaseIndex_(
unsigned canonicalPhaseIdx)
572 if (numStoragePhases == 3)
573 return canonicalPhaseIdx;
575 return FluidSystem::canonicalToActivePhaseIdx(canonicalPhaseIdx);
578 ConditionalStorage<enableTemperature || enableEnergy, Scalar> temperature_;
579 ConditionalStorage<enableEnergy, std::array<Scalar, numStoragePhases> > enthalpy_;
580 Scalar totalSaturation_;
581 std::array<Scalar, numStoragePhases> pressure_;
582 std::array<Scalar, numStoragePhases> pc_;
583 std::array<Scalar, numStoragePhases> saturation_;
584 std::array<Scalar, numStoragePhases> invB_;
585 std::array<Scalar, numStoragePhases> density_;
586 ConditionalStorage<enableDissolution,Scalar> Rs_;
587 ConditionalStorage<enableDissolution, Scalar> Rv_;
588 ConditionalStorage<enableBrine, Scalar> saltConcentration_;
589 unsigned short pvtRegionIdx_;
A fluid system which uses the black-oil model assumptions to calculate termodynamically meaningful qu...
A simple class which only stores a given member attribute if a boolean condition is true.
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
Provides the OPM_UNUSED macro.
Some templates to wrap the valgrind client request macros.
Implements a "tailor-made" fluid state class for the black-oil model.
Definition: BlackOilFluidState.hpp:105
void setDensity(unsigned phaseIdx, const Scalar &rho)
\ brief Set the density of a fluid phase
Definition: BlackOilFluidState.hpp:260
void setInvB(unsigned phaseIdx, const Scalar &b)
\ brief Set the inverse formation volume factor of a fluid phase
Definition: BlackOilFluidState.hpp:254
const Scalar & temperature(unsigned) const
Return the temperature [K].
Definition: BlackOilFluidState.hpp:314
void checkDefined() const
Make sure that all attributes are defined.
Definition: BlackOilFluidState.hpp:127
const Scalar & Rs() const
Return the gas dissulition factor of oil [m^3/m^3].
Definition: BlackOilFluidState.hpp:340
Scalar fugacityCoefficient(unsigned phaseIdx, unsigned compIdx) const
Return the fugacity coefficient of a component in a fluid phase [-].
Definition: BlackOilFluidState.hpp:547
Scalar molarDensity(unsigned phaseIdx) const
Return the molar density of a fluid phase [mol/m^3].
Definition: BlackOilFluidState.hpp:422
Scalar molarVolume(unsigned phaseIdx) const
Return the molar volume of a fluid phase [m^3/mol].
Definition: BlackOilFluidState.hpp:440
const Scalar & totalSaturation() const
Return the total saturation needed for sequential.
Definition: BlackOilFluidState.hpp:306
void assign(const FluidState &fs)
Retrieve all parameters from an arbitrary fluid state.
Definition: BlackOilFluidState.hpp:161
Scalar internalEnergy(unsigned phaseIdx) const
Return the specific internal energy [J/kg] of a given fluid phase.
Definition: BlackOilFluidState.hpp:412
void setPc(unsigned phaseIdx, const Scalar &pc)
Set the capillary pressure of a fluid phase [-].
Definition: BlackOilFluidState.hpp:214
void setRs(const Scalar &newRs)
Set the gas dissolution factor [m^3/m^3] of the oil phase.
Definition: BlackOilFluidState.hpp:268
void setEnthalpy(unsigned phaseIdx, const Scalar &value)
Set the specific enthalpy [J/kg] of a given fluid phase.
Definition: BlackOilFluidState.hpp:244
Scalar moleFraction(unsigned phaseIdx, unsigned compIdx) const
Return the mole fraction of a component in a fluid phase [-].
Definition: BlackOilFluidState.hpp:489
void setRv(const Scalar &newRv)
Set the oil vaporization factor [m^3/m^3] of the gas phase.
Definition: BlackOilFluidState.hpp:276
const Scalar & Rv() const
Return the oil vaporization factor of gas [m^3/m^3].
Definition: BlackOilFluidState.hpp:357
const Scalar & saturation(unsigned phaseIdx) const
Return the saturation of a fluid phase [-].
Definition: BlackOilFluidState.hpp:294
Scalar fugacity(unsigned phaseIdx, unsigned compIdx) const
Return the fugacity of a component in a fluid phase [Pa].
Definition: BlackOilFluidState.hpp:553
void setSaltConcentration(const Scalar &newSaltConcentration)
Set the salt concentration.
Definition: BlackOilFluidState.hpp:282
unsigned short pvtRegionIndex() const
Return the PVT region where the current fluid state is assumed to be part of.
Definition: BlackOilFluidState.hpp:388
void setPressure(unsigned phaseIdx, const Scalar &p)
Set the pressure of a fluid phase [-].
Definition: BlackOilFluidState.hpp:202
void setSaturation(unsigned phaseIdx, const Scalar &S)
Set the saturation of a fluid phase [-].
Definition: BlackOilFluidState.hpp:208
const Scalar & enthalpy(unsigned phaseIdx) const
Return the specific enthalpy [J/kg] of a given fluid phase.
Definition: BlackOilFluidState.hpp:403
const Scalar & pressure(unsigned phaseIdx) const
Return the pressure of a fluid phase [Pa].
Definition: BlackOilFluidState.hpp:288
Scalar viscosity(unsigned phaseIdx) const
Return the dynamic viscosity of a fluid phase [Pa s].
Definition: BlackOilFluidState.hpp:446
const Scalar & pc(unsigned phaseIdx) const
Return the capillary pressure of a fluid phase [-].
Definition: BlackOilFluidState.hpp:300
void setPvtRegionIndex(unsigned newPvtRegionIdx)
Set the index of the fluid region.
Definition: BlackOilFluidState.hpp:196
Scalar molarity(unsigned phaseIdx, unsigned compIdx) const
Return the partial molar density of a component in a fluid phase [mol / m^3].
Definition: BlackOilFluidState.hpp:530
const Scalar & invB(unsigned phaseIdx) const
Return the inverse formation volume factor of a fluid phase [-].
Definition: BlackOilFluidState.hpp:330
Scalar density(unsigned phaseIdx) const
Return the density [kg/m^3] of a given fluid phase.
Definition: BlackOilFluidState.hpp:394
Scalar massFraction(unsigned phaseIdx, unsigned compIdx) const
Return the mass fraction of a component in a fluid phase [-].
Definition: BlackOilFluidState.hpp:452
Scalar averageMolarMass(unsigned phaseIdx) const
Return the partial molar density of a fluid phase [kg / mol].
Definition: BlackOilFluidState.hpp:536
void setTemperature(const Scalar &value)
Set the temperature [K].
Definition: BlackOilFluidState.hpp:231
void setTotalSaturation(const Scalar &value)
Set the total saturation used for sequential methods.
Definition: BlackOilFluidState.hpp:220
const Scalar & saltConcentration() const
Return the concentration of salt in water.
Definition: BlackOilFluidState.hpp:370