28 #ifndef EWOMS_BLACK_OIL_LOCAL_RESIDUAL_HH 29 #define EWOMS_BLACK_OIL_LOCAL_RESIDUAL_HH 42 template <
class TypeTag>
45 typedef typename GET_PROP_TYPE(TypeTag, IntensiveQuantities) IntensiveQuantities;
46 typedef typename GET_PROP_TYPE(TypeTag, ExtensiveQuantities) ExtensiveQuantities;
47 typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
48 typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
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, EqVector) EqVector;
52 typedef typename GET_PROP_TYPE(TypeTag, RateVector) RateVector;
53 typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
55 enum { conti0EqIdx = Indices::conti0EqIdx };
60 enum { gasPhaseIdx = FluidSystem::gasPhaseIdx };
61 enum { oilPhaseIdx = FluidSystem::oilPhaseIdx };
62 enum { waterPhaseIdx = FluidSystem::waterPhaseIdx };
64 enum { gasCompIdx = FluidSystem::gasCompIdx };
65 enum { oilCompIdx = FluidSystem::oilCompIdx };
66 enum { waterCompIdx = FluidSystem::waterCompIdx };
67 enum { compositionSwitchIdx = Indices::compositionSwitchIdx };
69 static const bool compositionSwitchEnabled = Indices::gasEnabled;
70 static const bool waterEnabled = Indices::waterEnabled;
72 static constexpr
bool blackoilConserveSurfaceVolume =
GET_PROP_VALUE(TypeTag, BlackoilConserveSurfaceVolume);
74 typedef Opm::MathToolbox<Evaluation> Toolbox;
83 template <
class LhsEval>
85 const ElementContext& elemCtx,
87 unsigned timeIdx)
const 90 const IntensiveQuantities& intQuants = elemCtx.intensiveQuantities(dofIdx, timeIdx);
91 const auto& fs = intQuants.fluidState();
95 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
96 if (!FluidSystem::phaseIsActive(phaseIdx))
99 unsigned activeCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::solventComponentIndex(phaseIdx));
100 LhsEval surfaceVolume =
101 Toolbox::template decay<LhsEval>(fs.saturation(phaseIdx))
102 * Toolbox::template decay<LhsEval>(fs.invB(phaseIdx))
103 * Toolbox::template decay<LhsEval>(intQuants.porosity());
105 storage[conti0EqIdx + activeCompIdx] += surfaceVolume;
108 if (phaseIdx == oilPhaseIdx && FluidSystem::enableDissolvedGas()) {
109 unsigned activeGasCompIdx = Indices::canonicalToActiveComponentIndex(gasCompIdx);
110 storage[conti0EqIdx + activeGasCompIdx] +=
111 Toolbox::template decay<LhsEval>(intQuants.fluidState().Rs())
116 if (phaseIdx == gasPhaseIdx && FluidSystem::enableVaporizedOil()) {
117 unsigned activeOilCompIdx = Indices::canonicalToActiveComponentIndex(oilCompIdx);
118 storage[conti0EqIdx + activeOilCompIdx] +=
119 Toolbox::template decay<LhsEval>(intQuants.fluidState().Rv())
125 if (!blackoilConserveSurfaceVolume) {
126 unsigned pvtRegionIdx = intQuants.pvtRegionIndex();
128 unsigned activeWaterCompIdx = Indices::canonicalToActiveComponentIndex(waterCompIdx);
129 storage[conti0EqIdx + activeWaterCompIdx] *=
130 FluidSystem::referenceDensity(waterPhaseIdx, pvtRegionIdx);
132 if (compositionSwitchEnabled) {
133 unsigned activeGasCompIdx = Indices::canonicalToActiveComponentIndex(gasCompIdx);
134 storage[conti0EqIdx + activeGasCompIdx] *=
135 FluidSystem::referenceDensity(gasPhaseIdx, pvtRegionIdx);
137 unsigned activeOilCompIdx = Indices::canonicalToActiveComponentIndex(oilCompIdx);
138 storage[conti0EqIdx + activeOilCompIdx] *=
139 FluidSystem::referenceDensity(oilPhaseIdx, pvtRegionIdx);
142 SolventModule::addStorage(storage, intQuants);
145 PolymerModule::addStorage(storage, intQuants);
153 const ElementContext& elemCtx,
155 unsigned timeIdx)
const 157 assert(timeIdx == 0);
161 const ExtensiveQuantities& extQuants = elemCtx.extensiveQuantities(scvfIdx, timeIdx);
162 unsigned focusDofIdx = elemCtx.focusDofIndex();
163 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
164 if (!FluidSystem::phaseIsActive(phaseIdx))
167 unsigned upIdx =
static_cast<unsigned>(extQuants.upstreamIndex(phaseIdx));
168 const IntensiveQuantities& up = elemCtx.intensiveQuantities(upIdx, timeIdx);
169 if (upIdx == focusDofIdx)
170 evalPhaseFluxes_<Evaluation>(flux, phaseIdx, extQuants, up);
172 evalPhaseFluxes_<Scalar>(flux, phaseIdx, extQuants, up);
175 if (!blackoilConserveSurfaceVolume) {
176 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
177 unsigned activeCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::solventComponentIndex(phaseIdx));
178 unsigned upIdx =
static_cast<unsigned>(extQuants.upstreamIndex(phaseIdx));
179 const IntensiveQuantities& up = elemCtx.intensiveQuantities(upIdx, timeIdx);
180 unsigned pvtRegionIdx = up.pvtRegionIndex();
181 Scalar refDensity = FluidSystem::referenceDensity(phaseIdx, pvtRegionIdx);
182 flux[conti0EqIdx + activeCompIdx] *= refDensity;
187 SolventModule::computeFlux(flux, elemCtx, scvfIdx, timeIdx);
190 PolymerModule::computeFlux(flux, elemCtx, scvfIdx, timeIdx);
197 const ElementContext& elemCtx,
199 unsigned timeIdx)
const 202 elemCtx.problem().source(source, elemCtx, dofIdx, timeIdx);
206 template <
class UpEval>
207 void evalPhaseFluxes_(RateVector& flux,
209 const ExtensiveQuantities& extQuants,
210 const IntensiveQuantities& up)
const 212 unsigned compIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::solventComponentIndex(phaseIdx));
213 const auto& fs = up.fluidState();
215 Evaluation surfaceVolumeFlux =
216 Toolbox::template decay<UpEval>(fs.invB(phaseIdx))
217 * extQuants.volumeFlux(phaseIdx);
219 flux[conti0EqIdx + compIdx] +=
223 if (phaseIdx == oilPhaseIdx && FluidSystem::enableDissolvedGas()) {
224 flux[conti0EqIdx + Indices::canonicalToActiveComponentIndex(gasCompIdx)] +=
225 Toolbox::template decay<UpEval>(fs.Rs())
231 if (phaseIdx == gasPhaseIdx && FluidSystem::enableVaporizedOil()) {
232 flux[conti0EqIdx + Indices::canonicalToActiveComponentIndex(oilCompIdx)] +=
233 Toolbox::template decay<UpEval>(fs.Rv())
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.
#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 computeSource(RateVector &source, const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx) const
Calculate the source term of the equation.
Definition: blackoillocalresidual.hh:196
Contains the high level supplements required to extend the black oil model by solvents.
Definition: blackoilprimaryvariables.hh:47
Calculates the local residual of the black oil model.
Definition: blackoillocalresidual.hh:43
void computeFlux(RateVector &flux, const ElementContext &elemCtx, unsigned scvfIdx, unsigned timeIdx) const
Evaluates the total mass flux of all conservation quantities over a face of a sub-control volume...
Definition: blackoillocalresidual.hh:152
Contains the high level supplements required to extend the black oil model by polymer.
Definition: blackoilpolymermodules.hh:64
void computeStorage(Dune::FieldVector< LhsEval, numEq > &storage, const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx) const
Evaluate the amount all conservation quantities (e.g.
Definition: blackoillocalresidual.hh:84
Contains the classes required to extend the black-oil model by polymer.