28 #ifndef EWOMS_BLACK_OIL_INTENSIVE_QUANTITIES_HH 29 #define EWOMS_BLACK_OIL_INTENSIVE_QUANTITIES_HH 37 #include <opm/material/fluidstates/CompositionalFluidState.hpp> 38 #include <opm/common/Valgrind.hpp> 40 #include <dune/common/fmatrix.hh> 53 template <
class TypeTag>
56 ,
public GET_PROP_TYPE(TypeTag, FluxModule)::FluxIntensiveQuantities
60 typedef typename GET_PROP_TYPE(TypeTag, DiscIntensiveQuantities) ParentType;
61 typedef typename GET_PROP_TYPE(TypeTag, IntensiveQuantities) Implementation;
63 typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
64 typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
65 typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
66 typedef typename GET_PROP_TYPE(TypeTag, MaterialLaw) MaterialLaw;
67 typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
68 typedef typename GET_PROP_TYPE(TypeTag, PrimaryVariables) PrimaryVariables;
69 typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
70 typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
71 typedef typename GET_PROP_TYPE(TypeTag, FluxModule) FluxModule;
78 enum { waterCompIdx = FluidSystem::waterCompIdx };
79 enum { oilCompIdx = FluidSystem::oilCompIdx };
80 enum { gasCompIdx = FluidSystem::gasCompIdx };
81 enum { waterPhaseIdx = FluidSystem::waterPhaseIdx };
82 enum { oilPhaseIdx = FluidSystem::oilPhaseIdx };
83 enum { gasPhaseIdx = FluidSystem::gasPhaseIdx };
84 enum { dimWorld = GridView::dimensionworld };
85 enum { compositionSwitchIdx = Indices::compositionSwitchIdx };
87 static const bool compositionSwitchEnabled = Indices::gasEnabled;
88 static const bool waterEnabled = Indices::waterEnabled;
91 typedef Opm::MathToolbox<Evaluation> Toolbox;
92 typedef Dune::FieldMatrix<Scalar, dimWorld, dimWorld> DimMatrix;
93 typedef typename FluxModule::FluxIntensiveQuantities FluxIntensiveQuantities;
99 fluidState_.setRs(0.0);
100 fluidState_.setRv(0.0);
110 void update(
const ElementContext& elemCtx,
unsigned dofIdx,
unsigned timeIdx)
112 ParentType::update(elemCtx, dofIdx, timeIdx);
116 const auto& problem = elemCtx.problem();
117 const auto& priVars = elemCtx.primaryVars(dofIdx, timeIdx);
119 unsigned globalSpaceIdx = elemCtx.globalSpaceIndex(dofIdx, timeIdx);
120 unsigned pvtRegionIdx = priVars.pvtRegionIndex();
121 fluidState_.setPvtRegionIndex(pvtRegionIdx);
126 Sw = priVars.makeEvaluation(Indices::waterSaturationIdx, timeIdx);
129 if (compositionSwitchEnabled)
131 if (priVars.primaryVarsMeaning() == PrimaryVariables::Sw_po_Sg)
133 Sg = priVars.makeEvaluation(Indices::compositionSwitchIdx, timeIdx);
134 else if (priVars.primaryVarsMeaning() == PrimaryVariables::Sw_pg_Rv) {
140 Sg -= priVars.makeEvaluation(Indices::solventSaturationIdx, timeIdx);
144 assert(priVars.primaryVarsMeaning() == PrimaryVariables::Sw_po_Rs);
150 Opm::Valgrind::CheckDefined(Sg);
151 Opm::Valgrind::CheckDefined(Sw);
153 Evaluation So = 1.0 - Sw - Sg;
157 So -= priVars.makeEvaluation(Indices::solventSaturationIdx, timeIdx);
159 fluidState_.setSaturation(waterPhaseIdx, Sw);
160 fluidState_.setSaturation(gasPhaseIdx, Sg);
161 fluidState_.setSaturation(oilPhaseIdx, So);
163 asImp_().solventPreSatFuncUpdate_(elemCtx, dofIdx, timeIdx);
166 Evaluation pC[numPhases];
167 const auto& materialParams = problem.materialLawParams(elemCtx, dofIdx, timeIdx);
168 MaterialLaw::capillaryPressures(pC, materialParams, fluidState_);
171 if (priVars.primaryVarsMeaning() == PrimaryVariables::Sw_pg_Rv) {
172 const Evaluation& pg = priVars.makeEvaluation(Indices::pressureSwitchIdx, timeIdx);
173 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
174 fluidState_.setPressure(phaseIdx, pg + (pC[phaseIdx] - pC[gasPhaseIdx]));
178 const Evaluation& po = priVars.makeEvaluation(Indices::pressureSwitchIdx, timeIdx);
179 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
180 fluidState_.setPressure(phaseIdx, po + (pC[phaseIdx] - pC[oilPhaseIdx]));
185 MaterialLaw::relativePermeabilities(mobility_, materialParams, fluidState_);
186 Opm::Valgrind::CheckDefined(mobility_);
189 asImp_().solventPostSatFuncUpdate_(elemCtx, dofIdx, timeIdx);
191 Scalar SoMax = elemCtx.model().maxOilSaturation(globalSpaceIdx);
195 if (priVars.primaryVarsMeaning() == PrimaryVariables::Sw_po_Sg) {
198 if (FluidSystem::enableDissolvedGas()) {
199 const Evaluation& RsSat =
200 FluidSystem::saturatedDissolutionFactor(fluidState_,
204 fluidState_.setRs(RsSat);
207 fluidState_.setRs(0.0);
209 if (FluidSystem::enableVaporizedOil()) {
210 const Evaluation& RvSat =
211 FluidSystem::saturatedDissolutionFactor(fluidState_,
215 fluidState_.setRv(RvSat);
218 fluidState_.setRv(0.0);
220 else if (priVars.primaryVarsMeaning() == PrimaryVariables::Sw_po_Rs) {
223 const auto& Rs = priVars.makeEvaluation(Indices::compositionSwitchIdx, timeIdx);
224 fluidState_.setRs(Rs);
226 if (FluidSystem::enableVaporizedOil()) {
230 FluidSystem::saturatedDissolutionFactor(fluidState_,
235 fluidState_.setRv(RvSat);
238 fluidState_.setRv(0.0);
241 assert(priVars.primaryVarsMeaning() == PrimaryVariables::Sw_pg_Rv);
243 const auto& Rv = priVars.makeEvaluation(Indices::compositionSwitchIdx, timeIdx);
244 fluidState_.setRv(Rv);
246 if (FluidSystem::enableDissolvedGas()) {
250 FluidSystem::saturatedDissolutionFactor(fluidState_,
255 fluidState_.setRs(RsSat);
258 fluidState_.setRs(0.0);
261 typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
262 typename FluidSystem::template ParameterCache<Evaluation> paramCache;
263 paramCache.setRegionIndex(pvtRegionIdx);
264 paramCache.setMaxOilSat(SoMax);
265 paramCache.updateAll(fluidState_);
268 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
269 if (!FluidSystem::phaseIsActive(phaseIdx))
272 const auto& b = FluidSystem::inverseFormationVolumeFactor(fluidState_, phaseIdx, pvtRegionIdx);
273 fluidState_.setInvB(phaseIdx, b);
275 const auto& mu = FluidSystem::viscosity(fluidState_, paramCache, phaseIdx);
276 mobility_[phaseIdx] /= mu;
278 Opm::Valgrind::CheckDefined(mobility_);
282 if (FluidSystem::phaseIsActive(waterPhaseIdx)) {
283 rho = fluidState_.invB(waterPhaseIdx);
284 rho *= FluidSystem::referenceDensity(waterPhaseIdx, pvtRegionIdx);
285 fluidState_.setDensity(waterPhaseIdx, rho);
288 if (FluidSystem::phaseIsActive(gasPhaseIdx)) {
289 rho = fluidState_.invB(gasPhaseIdx);
290 rho *= FluidSystem::referenceDensity(gasPhaseIdx, pvtRegionIdx);
291 if (FluidSystem::enableVaporizedOil()) {
293 fluidState_.invB(gasPhaseIdx) *
295 FluidSystem::referenceDensity(oilPhaseIdx, pvtRegionIdx);
297 fluidState_.setDensity(gasPhaseIdx, rho);
300 if (FluidSystem::phaseIsActive(oilPhaseIdx)) {
301 rho = fluidState_.invB(oilPhaseIdx);
302 rho *= FluidSystem::referenceDensity(oilPhaseIdx, pvtRegionIdx);
303 if (FluidSystem::enableDissolvedGas()) {
305 fluidState_.invB(oilPhaseIdx) *
307 FluidSystem::referenceDensity(gasPhaseIdx, pvtRegionIdx);
309 fluidState_.setDensity(oilPhaseIdx, rho);
313 porosity_ = problem.porosity(elemCtx, dofIdx, timeIdx);
317 Scalar rockCompressibility = problem.rockCompressibility(elemCtx, dofIdx, timeIdx);
318 if (rockCompressibility > 0.0) {
319 Scalar rockRefPressure = problem.rockReferencePressure(elemCtx, dofIdx, timeIdx);
320 Evaluation x = rockCompressibility*(fluidState_.pressure(oilPhaseIdx) - rockRefPressure);
321 porosity_ *= 1.0 + x + 0.5*x*x;
324 asImp_().solventPvtUpdate_(elemCtx, dofIdx, timeIdx);
325 asImp_().polymerPropertiesUpdate_(elemCtx, dofIdx, timeIdx);
330 FluxIntensiveQuantities::update_(elemCtx, dofIdx, timeIdx);
334 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
335 if (!FluidSystem::phaseIsActive(phaseIdx))
338 assert(std::isfinite(Toolbox::value(fluidState_.density(phaseIdx))));
339 assert(std::isfinite(Toolbox::value(fluidState_.saturation(phaseIdx))));
340 assert(std::isfinite(Toolbox::value(fluidState_.temperature(phaseIdx))));
341 assert(std::isfinite(Toolbox::value(fluidState_.pressure(phaseIdx))));
342 assert(std::isfinite(Toolbox::value(fluidState_.invB(phaseIdx))));
344 assert(std::isfinite(Toolbox::value(fluidState_.Rs())));
345 assert(std::isfinite(Toolbox::value(fluidState_.Rv())));
353 {
return fluidState_; }
358 const Evaluation&
mobility(
unsigned phaseIdx)
const 359 {
return mobility_[phaseIdx]; }
365 {
return porosity_; }
382 {
return fluidState_.pvtRegionIndex(); }
390 return fluidState_.viscosity(phaseIdx)*
mobility(phaseIdx);
398 Implementation& asImp_()
399 {
return *
static_cast<Implementation*
>(
this); }
401 FluidState fluidState_;
402 Evaluation porosity_;
403 Evaluation mobility_[numPhases];
Provides the volumetric quantities required for the equations needed by the polymers extension of the...
Definition: blackoilpolymermodules.hh:839
Provides the volumetric quantities required for the equations needed by the solvents extension of the...
Definition: blackoilsolventmodules.hh:879
Definition: baseauxiliarymodule.hh:37
#define GET_PROP_VALUE(TypeTag, PropTagName)
Access the value attribute of a property for a type tag.
Definition: propertysystem.hh:469
Contains the classes required to extend the black-oil model by solvents.
Implements a "taylor-made" fluid state class for the black-oil model.
#define GET_PROP_TYPE(TypeTag, PropTagName)
Access the type attribute of a property for a type tag.
Definition: propertysystem.hh:486
Declares the properties required by the black oil model.
void update(const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx)
Definition: blackoilintensivequantities.hh:110
Evaluation relativePermeability(unsigned phaseIdx) const
Returns the relative permeability of a given phase within the control volume.
Definition: blackoilintensivequantities.hh:387
const Evaluation & mobility(unsigned phaseIdx) const
Returns the effective mobility of a given phase within the control volume.
Definition: blackoilintensivequantities.hh:358
auto pvtRegionIndex() const -> decltype(std::declval< FluidState >().pvtRegionIndex())
Returns the index of the PVT region used to calculate the thermodynamic quantities.
Definition: blackoilintensivequantities.hh:380
Implements a "taylor-made" fluid state class for the black-oil model.
Definition: blackoilfluidstate.hh:49
const FluidState & fluidState() const
Returns the phase state for the control-volume.
Definition: blackoilintensivequantities.hh:352
const Evaluation & porosity() const
Returns the average porosity within the control volume.
Definition: blackoilintensivequantities.hh:364
Contains the classes required to extend the black-oil model by polymer.
Contains the quantities which are are constant within a finite volume in the black-oil model...
Definition: blackoilintensivequantities.hh:54