26 #ifndef EWOMS_BLACK_OIL_PRIMARY_VARIABLES_HH
27 #define EWOMS_BLACK_OIL_PRIMARY_VARIABLES_HH
33 #include <dune/common/fvector.hh>
35 #include <opm/material/constraintsolvers/NcpFlash.hpp>
36 #include <opm/material/fluidstates/CompositionalFluidState.hpp>
44 template <
class TypeTag>
49 typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
50 typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
51 typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
52 typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
53 typedef typename GET_PROP_TYPE(TypeTag, MaterialLaw) MaterialLaw;
54 typedef typename GET_PROP_TYPE(TypeTag, MaterialLawParams) MaterialLawParams;
60 enum { gasPressureIdx = Indices::gasPressureIdx };
61 enum { waterSaturationIdx = Indices::waterSaturationIdx };
62 enum { switchIdx = Indices::switchIdx };
66 enum { gasPhaseIdx = FluidSystem::gasPhaseIdx };
67 enum { waterPhaseIdx = FluidSystem::waterPhaseIdx };
68 enum { oilPhaseIdx = FluidSystem::oilPhaseIdx };
72 enum { gasCompIdx = FluidSystem::gasCompIdx };
73 enum { waterCompIdx = FluidSystem::waterCompIdx };
74 enum { oilCompIdx = FluidSystem::oilCompIdx };
76 typedef typename Opm::MathToolbox<Evaluation> Toolbox;
77 typedef Dune::FieldVector<Scalar, numComponents> ComponentVector;
79 static_assert(numPhases == 3,
"The black-oil model assumes three phases!");
80 static_assert(numComponents == 3,
"The black-oil model assumes three components!");
83 enum SwitchingVarMeaning {
92 Valgrind::SetUndefined(*
this);
94 temperature_ = FluidSystem::surfaceTemperature;
103 Valgrind::SetUndefined(switchingVarMeaning_);
105 temperature_ = FluidSystem::surfaceTemperature;
113 , temperature_(FluidSystem::surfaceTemperature)
114 , switchingVarMeaning_(value.switchingVarMeaning_)
115 , pvtRegionIdx_(value.pvtRegionIdx_)
118 void setPvtRegionIndex(
int value)
119 { pvtRegionIdx_ = value; }
121 unsigned pvtRegionIndex()
const
122 {
return pvtRegionIdx_; }
124 SwitchingVarMeaning switchingVarMeaning()
const
125 {
return switchingVarMeaning_; }
127 void setSwitchingVarMeaning(SwitchingVarMeaning newMeaning)
128 { switchingVarMeaning_ = newMeaning; }
133 template <
class Flu
idState>
134 void assignMassConservative(
const FluidState &fluidState,
135 const MaterialLawParams &matParams,
136 bool isInEquilibrium =
false)
138 typedef typename std::remove_reference<typename FluidState::Scalar>::type ConstEvaluation;
139 typedef typename std::remove_const<ConstEvaluation>::type FsEvaluation;
140 typedef typename Opm::MathToolbox<FsEvaluation> FsToolbox;
144 for (
int phaseIdx = 1; phaseIdx < numPhases; ++phaseIdx) {
145 Valgrind::CheckDefined(fluidState.temperature(0));
146 Valgrind::CheckDefined(fluidState.temperature(phaseIdx));
148 assert(fluidState.temperature(0) == fluidState.temperature(phaseIdx));
154 if (isInEquilibrium) {
161 typename FluidSystem::ParameterCache paramCache;
162 paramCache.setRegionIndex(pvtRegionIdx_);
165 typedef Opm::NcpFlash<Scalar, FluidSystem> NcpFlash;
166 typedef Opm::CompositionalFluidState<Scalar, FluidSystem> FlashFluidState;
167 FlashFluidState fsFlash;
168 fsFlash.setTemperature(FsToolbox::value(fluidState.temperature(0)));
169 for (
int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
170 fsFlash.setPressure(phaseIdx, FsToolbox::value(fluidState.pressure(phaseIdx)));
171 fsFlash.setSaturation(phaseIdx, FsToolbox::value(fluidState.saturation(phaseIdx)));
172 for (
int compIdx = 0; compIdx < numComponents; ++compIdx)
173 fsFlash.setMoleFraction(phaseIdx, compIdx, FsToolbox::value(fluidState.moleFraction(phaseIdx, compIdx)));
176 for (
int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
177 Scalar rho = FluidSystem::template density<FlashFluidState, Scalar>(fsFlash, paramCache, phaseIdx);
178 fsFlash.setDensity(phaseIdx, rho);
182 ComponentVector globalMolarities(0.0);
183 for (
int compIdx = 0; compIdx < numComponents; ++compIdx) {
184 for (
int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
185 globalMolarities[compIdx] +=
186 fsFlash.saturation(phaseIdx) * fsFlash.molarity(phaseIdx, compIdx);
194 NcpFlash::template solve<MaterialLaw>(fsFlash, paramCache, matParams, globalMolarities);
203 template <
class Flu
idState>
206 typedef typename std::remove_reference<typename FluidState::Scalar>::type ConstEvaluation;
207 typedef typename std::remove_const<ConstEvaluation>::type FsEvaluation;
208 typedef typename Opm::MathToolbox<FsEvaluation> FsToolbox;
210 temperature_ = FsToolbox::value(fluidState.temperature(gasPhaseIdx));
214 (*this)[gasPressureIdx] = FsToolbox::value(fluidState.pressure(gasPhaseIdx));
215 (*this)[waterSaturationIdx] = FsToolbox::value(fluidState.saturation(waterPhaseIdx));
217 bool gasPresent = (fluidState.saturation(gasPhaseIdx) > 0.0);
218 bool oilPresent = (fluidState.saturation(oilPhaseIdx) > 0.0);
220 if ((gasPresent && oilPresent) || (!gasPresent && !oilPresent))
222 switchingVarMeaning_ = GasSaturation;
223 else if (oilPresent) {
225 if (FluidSystem::enableDissolvedGas())
226 switchingVarMeaning_ = GasMoleFractionInOil;
228 switchingVarMeaning_ = GasSaturation;
233 if (FluidSystem::enableVaporizedOil())
234 switchingVarMeaning_ = OilMoleFractionInGas;
236 switchingVarMeaning_ = GasSaturation;
241 if (switchingVarMeaning() == GasSaturation)
242 (*this)[switchIdx] = FsToolbox::value(fluidState.saturation(gasPhaseIdx));
243 else if (switchingVarMeaning() == GasMoleFractionInOil)
244 (*
this)[switchIdx] = FsToolbox::value(fluidState.moleFraction(oilPhaseIdx, gasCompIdx));
246 assert(switchingVarMeaning() == OilMoleFractionInGas);
247 (*this)[switchIdx] = FsToolbox::value(fluidState.moleFraction(gasPhaseIdx, oilCompIdx));
257 bool adaptSwitchingVariable()
259 Scalar pg = (*this)[Indices::gasPressureIdx];
261 if (switchingVarMeaning() == GasSaturation) {
263 Scalar Sw = (*this)[Indices::waterSaturationIdx];
264 Scalar Sg = (*this)[Indices::switchIdx];
265 Scalar So = 1 - Sw - Sg;
267 if (Sg < 0.0 && So > 0.0 && FluidSystem::enableDissolvedGas()) {
270 Scalar xoGsat = FluidSystem::saturatedOilGasMoleFraction(temperature_,
273 setSwitchingVarMeaning(GasMoleFractionInOil);
274 (*this)[Indices::switchIdx] = xoGsat;
278 if (So < 0.0 && Sg > 0.0 && FluidSystem::enableVaporizedOil()) {
282 Scalar xgOsat = FluidSystem::saturatedGasOilMoleFraction(temperature_,
285 setSwitchingVarMeaning(OilMoleFractionInGas);
286 (*this)[Indices::switchIdx] = xgOsat;
291 else if (switchingVarMeaning() == OilMoleFractionInGas) {
293 Scalar Sw = (*this)[Indices::waterSaturationIdx];
294 Scalar xgO = (*this)[Indices::switchIdx];
295 Scalar xgOsat = FluidSystem::saturatedGasOilMoleFraction(temperature_,
300 setSwitchingVarMeaning(GasSaturation);
301 (*this)[Indices::switchIdx] = 1 - Sw;
308 assert(switchingVarMeaning() == GasMoleFractionInOil);
311 Scalar xoG = (*this)[Indices::switchIdx];
312 Scalar xoGsat = FluidSystem::saturatedOilGasMoleFraction(temperature_,
317 setSwitchingVarMeaning(GasSaturation);
318 (*this)[Indices::switchIdx] = 0.0;
331 ParentType::operator=(other);
332 temperature_ = other.temperature_;
333 switchingVarMeaning_ = other.switchingVarMeaning_;
334 pvtRegionIdx_ = other.pvtRegionIdx_;
340 for (
int i = 0; i < numEq; ++i)
348 SwitchingVarMeaning switchingVarMeaning_;
349 unsigned char pvtRegionIdx_;
Represents the primary variables used by the a model.
#define GET_PROP_VALUE(TypeTag, PropTagName)
Access the value attribute of a property for a type tag.
Definition: propertysystem.hh:468
Represents the primary variables used by the black-oil model.
Definition: blackoilprimaryvariables.hh:45
Represents the primary variables used by the a model.
Definition: fvbaseprimaryvariables.hh:41
Declares the properties required by the black oil model.
Definition: baseauxiliarymodule.hh:35
void assignNaive(const FluidState &fluidState)
Assign the primary variables "somehow" from a fluid state.
Definition: fvbaseprimaryvariables.hh:98