26 #ifndef EWOMS_BLACK_OIL_INTENSIVE_QUANTITIES_HH
27 #define EWOMS_BLACK_OIL_INTENSIVE_QUANTITIES_HH
31 #include <opm/material/fluidstates/CompositionalFluidState.hpp>
32 #include <dune/common/fmatrix.hh>
42 template <
class TypeTag>
45 ,
public GET_PROP_TYPE(TypeTag, FluxModule)::FluxIntensiveQuantities
47 typedef typename GET_PROP_TYPE(TypeTag, DiscIntensiveQuantities) ParentType;
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, FluidSystem) FluidSystem;
52 typedef typename GET_PROP_TYPE(TypeTag, MaterialLaw) MaterialLaw;
53 typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
54 typedef typename GET_PROP_TYPE(TypeTag, PrimaryVariables) PrimaryVariables;
55 typedef typename GET_PROP_TYPE(TypeTag, BlackOilFluidState) FluidState;
56 typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
57 typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
58 typedef typename GET_PROP_TYPE(TypeTag, FluxModule) FluxModule;
62 enum { waterCompIdx = FluidSystem::waterCompIdx };
63 enum { oilCompIdx = FluidSystem::oilCompIdx };
64 enum { gasCompIdx = FluidSystem::gasCompIdx };
65 enum { waterPhaseIdx = FluidSystem::waterPhaseIdx };
66 enum { oilPhaseIdx = FluidSystem::oilPhaseIdx };
67 enum { gasPhaseIdx = FluidSystem::gasPhaseIdx };
68 enum { dimWorld = GridView::dimensionworld };
70 typedef Opm::MathToolbox<Evaluation> Toolbox;
71 typedef Dune::FieldMatrix<Scalar, dimWorld, dimWorld> DimMatrix;
72 typedef typename FluxModule::FluxIntensiveQuantities FluxIntensiveQuantities;
78 void update(
const ElementContext &elemCtx,
int dofIdx,
int timeIdx)
80 ParentType::update(elemCtx, dofIdx, timeIdx);
82 fluidState_.setTemperature(elemCtx.problem().temperature(elemCtx, dofIdx, timeIdx));
84 const auto& problem = elemCtx.problem();
85 const auto& priVars = elemCtx.primaryVars(dofIdx, timeIdx);
87 int pvtRegionIdx = priVars.pvtRegionIndex();
90 Evaluation Sw = priVars.makeEvaluation(Indices::waterSaturationIdx, timeIdx);
93 if (priVars.switchingVarMeaning() == PrimaryVariables::GasSaturation)
94 Sg = priVars.makeEvaluation(Indices::switchIdx, timeIdx);
95 else if (priVars.switchingVarMeaning() == PrimaryVariables::OilMoleFractionInGas)
98 assert(priVars.switchingVarMeaning() == PrimaryVariables::GasMoleFractionInOil);
102 fluidState_.setSaturation(waterPhaseIdx, Sw);
103 fluidState_.setSaturation(gasPhaseIdx, Sg);
104 fluidState_.setSaturation(oilPhaseIdx, 1 - Sw - Sg);
107 Evaluation pg = priVars.makeEvaluation(Indices::gasPressureIdx, timeIdx);
110 typedef typename GET_PROP_TYPE(TypeTag, MaterialLaw) MaterialLaw;
111 Evaluation pC[numPhases];
112 const auto &materialParams = problem.materialLawParams(elemCtx, dofIdx, timeIdx);
113 MaterialLaw::capillaryPressures(pC, materialParams, fluidState_);
114 for (
int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
115 fluidState_.setPressure(phaseIdx, pg + (pC[phaseIdx] - pC[gasPhaseIdx]));
116 if (fluidState_.pressure(phaseIdx) < 1e5) {
117 OPM_THROW(Opm::NumericalProblem,
118 "All pressures must be at least 1 bar.");
123 const auto& T = fluidState_.temperature(oilPhaseIdx);
126 for (
int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
127 for (
int compIdx = 0; compIdx < numComponents; ++compIdx)
128 fluidState_.setMoleFraction(phaseIdx, compIdx, 0.0);
131 fluidState_.setMoleFraction(oilPhaseIdx, oilCompIdx, 1.0);
132 fluidState_.setMoleFraction(gasPhaseIdx, gasCompIdx, 1.0);
133 fluidState_.setMoleFraction(waterPhaseIdx, waterCompIdx, 1.0);
137 if (priVars.switchingVarMeaning() == PrimaryVariables::GasSaturation) {
140 Evaluation xoG = 0.0;
141 if (FluidSystem::enableDissolvedGas())
142 xoG = FluidSystem::saturatedOilGasMoleFraction(T, pg, pvtRegionIdx);
143 fluidState_.setMoleFraction(oilPhaseIdx, gasCompIdx, xoG);
144 fluidState_.setMoleFraction(oilPhaseIdx, oilCompIdx, 1 - xoG);
146 Evaluation xgO = 0.0;
147 if (FluidSystem::enableVaporizedOil())
148 xgO = FluidSystem::saturatedGasOilMoleFraction(T, pg, pvtRegionIdx);
149 fluidState_.setMoleFraction(gasPhaseIdx, gasCompIdx, 1 - xgO);
150 fluidState_.setMoleFraction(gasPhaseIdx, oilCompIdx, xgO);
152 else if (priVars.switchingVarMeaning() == PrimaryVariables::GasMoleFractionInOil) {
155 const auto& xoG = priVars.makeEvaluation(Indices::switchIdx, timeIdx);
157 fluidState_.setMoleFraction(oilPhaseIdx, gasCompIdx, xoG);
158 fluidState_.setMoleFraction(oilPhaseIdx, oilCompIdx, 1 - xoG);
160 Evaluation xgO = 0.0;
161 if (FluidSystem::enableVaporizedOil())
162 xgO = FluidSystem::saturatedGasOilMoleFraction(T, pg, pvtRegionIdx);
163 fluidState_.setMoleFraction(gasPhaseIdx, gasCompIdx, 1 - xgO);
164 fluidState_.setMoleFraction(gasPhaseIdx, oilCompIdx, xgO);
167 assert(priVars.switchingVarMeaning() == PrimaryVariables::OilMoleFractionInGas);
171 const auto& xgO = priVars.makeEvaluation(Indices::switchIdx, timeIdx);
173 fluidState_.setMoleFraction(gasPhaseIdx, gasCompIdx, 1 - xgO);
174 fluidState_.setMoleFraction(gasPhaseIdx, oilCompIdx, xgO);
176 Evaluation xoG = 0.0;
177 if (FluidSystem::enableDissolvedGas())
178 xoG = FluidSystem::saturatedOilGasMoleFraction(T, pg, pvtRegionIdx);
179 fluidState_.setMoleFraction(oilPhaseIdx, gasCompIdx, xoG);
180 fluidState_.setMoleFraction(oilPhaseIdx, oilCompIdx, 1 - xoG);
184 MaterialLaw::relativePermeabilities(relativePermeability_, materialParams, fluidState_);
185 Valgrind::CheckDefined(relativePermeability_);
187 typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
188 typename FluidSystem::ParameterCache paramCache;
189 paramCache.setRegionIndex(pvtRegionIdx);
190 paramCache.updateAll(fluidState_);
193 for (
int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
194 const auto& rho = FluidSystem::density(fluidState_, paramCache, phaseIdx);
195 fluidState_.setDensity(phaseIdx, rho);
197 const auto& mu = FluidSystem::viscosity(fluidState_, paramCache, phaseIdx);
198 fluidState_.setViscosity(phaseIdx, mu);
200 mobility_[phaseIdx] = relativePermeability_[phaseIdx] / mu;
204 porosity_ = problem.porosity(elemCtx, dofIdx, timeIdx);
208 Scalar rockCompressibility = problem.rockCompressibility(elemCtx, dofIdx, timeIdx);
209 if (rockCompressibility > 0.0) {
210 Scalar rockRefPressure = problem.rockReferencePressure(elemCtx, dofIdx, timeIdx);
211 Evaluation x = rockCompressibility*(pg - rockRefPressure);
212 porosity_ *= 1.0 + x + 0.5*x*x;
216 intrinsicPerm_ = problem.intrinsicPermeability(elemCtx, dofIdx, timeIdx);
220 FluxIntensiveQuantities::update_(elemCtx, dofIdx, timeIdx);
224 for (
int phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
225 assert(std::isfinite(Toolbox::value(fluidState_.density(phaseIdx))));
226 assert(std::isfinite(Toolbox::value(fluidState_.saturation(phaseIdx))));
227 assert(std::isfinite(Toolbox::value(fluidState_.temperature(phaseIdx))));
228 assert(std::isfinite(Toolbox::value(fluidState_.pressure(phaseIdx))));
229 assert(std::isfinite(Toolbox::value(fluidState_.viscosity(phaseIdx))));
230 assert(std::isfinite(Toolbox::value(relativePermeability_[phaseIdx])));
231 for (
int compIdx = 0; compIdx < numComponents; ++ compIdx)
232 assert(std::isfinite(Toolbox::value(fluidState_.moleFraction(phaseIdx, compIdx))));
234 assert(std::isfinite(intrinsicPerm_.frobenius_norm()));
235 assert(std::isfinite(Toolbox::value(porosity_)));
243 {
return fluidState_; }
249 {
return intrinsicPerm_; }
255 {
return relativePermeability_[phaseIdx]; }
261 {
return mobility_[phaseIdx]; }
267 {
return porosity_; }
270 FluidState fluidState_;
271 Evaluation porosity_;
272 DimMatrix intrinsicPerm_;
273 Evaluation relativePermeability_[numPhases];
274 Evaluation mobility_[numPhases];
const DimMatrix & intrinsicPermeability() const
Returns the intrinsic permeability tensor a degree of freedom.
Definition: blackoilintensivequantities.hh:248
#define GET_PROP_VALUE(TypeTag, PropTagName)
Access the value attribute of a property for a type tag.
Definition: propertysystem.hh:468
const Evaluation & porosity() const
Returns the average porosity within the control volume.
Definition: blackoilintensivequantities.hh:266
#define GET_PROP_TYPE(TypeTag, PropTagName)
Access the type attribute of a property for a type tag.
Definition: propertysystem.hh:485
void update(const ElementContext &elemCtx, int dofIdx, int timeIdx)
Definition: blackoilintensivequantities.hh:78
Declares the properties required by the black oil model.
const FluidState & fluidState() const
Returns the phase state for the control-volume.
Definition: blackoilintensivequantities.hh:242
Contains the quantities which are are constant within a finite volume in the black-oil model...
Definition: blackoilintensivequantities.hh:43
Definition: baseauxiliarymodule.hh:35
const Evaluation & relativePermeability(int phaseIdx) const
Returns the relative permeability of a given phase within the control volume.
Definition: blackoilintensivequantities.hh:254
const Evaluation & mobility(int phaseIdx) const
Returns the effective mobility of a given phase within the control volume.
Definition: blackoilintensivequantities.hh:260