28 #ifndef EWOMS_BLACK_OIL_FLUID_STATE_HH
29 #define EWOMS_BLACK_OIL_FLUID_STATE_HH
33 #include <opm/common/Valgrind.hpp>
34 #include <opm/common/Unused.hpp>
35 #include <opm/common/ErrorMacros.hpp>
36 #include <opm/common/Exceptions.hpp>
48 template <
class TypeTag>
51 typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
52 typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
54 enum { waterPhaseIdx = FluidSystem::waterPhaseIdx };
55 enum { gasPhaseIdx = FluidSystem::gasPhaseIdx };
56 enum { oilPhaseIdx = FluidSystem::oilPhaseIdx };
58 enum { waterCompIdx = FluidSystem::waterCompIdx };
59 enum { gasCompIdx = FluidSystem::gasCompIdx };
60 enum { oilCompIdx = FluidSystem::oilCompIdx };
63 typedef Evaluation Scalar;
64 enum { numPhases = FluidSystem::numPhases };
65 enum { numComponents = FluidSystem::numComponents };
78 Opm::Valgrind::CheckDefined(pvtRegionIdx_);
80 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
81 Opm::Valgrind::CheckDefined(saturation_[phaseIdx]);
82 Opm::Valgrind::CheckDefined(pressure_[phaseIdx]);
83 Opm::Valgrind::CheckDefined(invB_[phaseIdx]);
86 Opm::Valgrind::CheckDefined(Rs_);
87 Opm::Valgrind::CheckDefined(Rv_);
95 template <
class Flu
idState>
96 void assign(
const FluidState& fs OPM_UNUSED)
101 void setPvtRegionIndex(
unsigned newPvtRegionIdx)
102 { pvtRegionIdx_ =
static_cast<unsigned short>(newPvtRegionIdx); }
104 void setPressure(
unsigned phaseIdx,
const Evaluation& p)
105 { pressure_[phaseIdx] = p; }
107 void setSaturation(
unsigned phaseIdx,
const Evaluation& S)
108 { saturation_[phaseIdx] = S; }
110 void setInvB(
unsigned phaseIdx,
const Evaluation& b)
111 { invB_[phaseIdx] = b; }
113 void setDensity(
unsigned phaseIdx,
const Evaluation& rho)
114 { density_[phaseIdx] = rho; }
116 void setRs(
const Evaluation& newRs)
119 void setRv(
const Evaluation& newRv)
122 const Evaluation& pressure(
unsigned phaseIdx)
const
123 {
return pressure_[phaseIdx]; }
125 const Evaluation& saturation(
unsigned phaseIdx)
const
126 {
return saturation_[phaseIdx]; }
128 const Evaluation& temperature(
unsigned phaseIdx OPM_UNUSED)
const
129 {
return temperature_; }
131 const Evaluation& invB(
unsigned phaseIdx)
const
132 {
return invB_[phaseIdx]; }
134 const Evaluation& Rs()
const
137 const Evaluation& Rv()
const
140 unsigned short pvtRegionIndex()
const
141 {
return pvtRegionIdx_; }
143 bool phaseIsPresent(
unsigned phaseIdx)
const
144 {
return saturation_[phaseIdx] > 0.0; }
149 Evaluation density(
unsigned phaseIdx)
const
150 {
return density_[phaseIdx]; }
152 Evaluation molarDensity(
unsigned phaseIdx)
const
154 const auto& rho = density(phaseIdx);
156 if (phaseIdx == waterPhaseIdx)
157 return rho/FluidSystem::molarMass(waterCompIdx, pvtRegionIdx_);
160 rho*(moleFraction(phaseIdx, gasCompIdx)/FluidSystem::molarMass(gasCompIdx, pvtRegionIdx_)
161 + moleFraction(phaseIdx, oilCompIdx)/FluidSystem::molarMass(oilCompIdx, pvtRegionIdx_));
165 Evaluation molarVolume(
unsigned phaseIdx)
const
166 {
return 1.0/molarDensity(phaseIdx); }
168 Evaluation viscosity(
unsigned phaseIdx)
const
169 {
return FluidSystem::viscosity(*
this, phaseIdx, pvtRegionIdx_); }
171 Evaluation enthalpy(
unsigned phaseIdx OPM_UNUSED)
const
173 OPM_THROW(Opm::NotImplemented,
174 "The black-oil model does not support energy conservation yet.");
177 Evaluation internalEnergy(
unsigned phaseIdx OPM_UNUSED)
const
179 OPM_THROW(Opm::NotImplemented,
180 "The black-oil model does not support energy conservation yet.");
183 Evaluation massFraction(
unsigned phaseIdx,
unsigned compIdx)
const
187 if (compIdx == waterCompIdx)
192 if (compIdx == waterCompIdx)
194 else if (compIdx == oilCompIdx)
195 return 1.0 - FluidSystem::convertRsToXoG(Rs_, pvtRegionIdx_);
197 assert(compIdx == gasCompIdx);
198 return FluidSystem::convertRsToXoG(Rs_, pvtRegionIdx_);
203 if (compIdx == waterCompIdx)
205 else if (compIdx == oilCompIdx)
206 return FluidSystem::convertRvToXgO(Rv_, pvtRegionIdx_);
208 assert(compIdx == gasCompIdx);
209 return 1.0 - FluidSystem::convertRvToXgO(Rv_, pvtRegionIdx_);
214 OPM_THROW(std::logic_error,
215 "Invalid phase or component index!");
218 Evaluation moleFraction(
unsigned phaseIdx,
unsigned compIdx)
const
222 if (compIdx == waterCompIdx)
227 if (compIdx == waterCompIdx)
229 else if (compIdx == oilCompIdx)
230 return 1.0 - FluidSystem::convertXoGToxoG(FluidSystem::convertRsToXoG(Rs_, pvtRegionIdx_),
233 assert(compIdx == gasCompIdx);
234 return FluidSystem::convertXoGToxoG(FluidSystem::convertRsToXoG(Rs_, pvtRegionIdx_),
240 if (compIdx == waterCompIdx)
242 else if (compIdx == oilCompIdx)
243 return FluidSystem::convertXgOToxgO(FluidSystem::convertRvToXgO(Rv_, pvtRegionIdx_),
246 assert(compIdx == gasCompIdx);
247 return 1.0 - FluidSystem::convertXgOToxgO(FluidSystem::convertRvToXgO(Rv_, pvtRegionIdx_),
253 OPM_THROW(std::logic_error,
254 "Invalid phase or component index!");
257 Evaluation molarity(
unsigned phaseIdx,
unsigned compIdx)
const
258 {
return moleFraction(phaseIdx, compIdx)*molarDensity(phaseIdx); }
260 Evaluation averageMolarMass(
unsigned phaseIdx)
const
262 Evaluation result(0.0);
263 for (
unsigned compIdx = 0; compIdx < numComponents; ++ compIdx)
264 result += FluidSystem::molarMass(compIdx, pvtRegionIdx_)*moleFraction(phaseIdx, compIdx);
268 Evaluation fugacityCoefficient(
unsigned phaseIdx,
unsigned compIdx)
const
269 {
return FluidSystem::fugacityCoefficient(*
this, phaseIdx, compIdx, pvtRegionIdx_); }
271 Evaluation fugacity(
unsigned phaseIdx,
unsigned compIdx)
const
274 fugacityCoefficient(phaseIdx, compIdx)
275 *moleFraction(phaseIdx, compIdx)
280 static const Evaluation temperature_;
281 std::array<Evaluation, numPhases> pressure_;
282 std::array<Evaluation, numPhases> saturation_;
283 std::array<Evaluation, numPhases> invB_;
284 std::array<Evaluation, numPhases> density_;
287 unsigned short pvtRegionIdx_;
290 template <
class TypeTag>
291 const typename BlackOilFluidState<TypeTag>::Evaluation BlackOilFluidState<TypeTag>::temperature_ =
292 BlackOilFluidState<TypeTag>::FluidSystem::surfaceTemperature;
void checkDefined() const
Make sure that all attributes are defined.
Definition: blackoilfluidstate.hh:75
void assign(const FluidState &fs OPM_UNUSED)
Retrieve all parameters from an arbitrary fluid state.
Definition: blackoilfluidstate.hh:96
Declares the properties required by the black oil model.
Implements a "taylor-made" fluid state class for the black-oil model.
Definition: blackoilfluidstate.hh:49