28#ifndef EWOMS_BLACK_OIL_LOCAL_RESIDUAL_HH
29#define EWOMS_BLACK_OIL_LOCAL_RESIDUAL_HH
41#include <opm/material/fluidstates/BlackOilFluidState.hpp>
49template <
class TypeTag>
62 enum { conti0EqIdx = Indices::conti0EqIdx };
67 enum { gasPhaseIdx = FluidSystem::gasPhaseIdx };
68 enum { oilPhaseIdx = FluidSystem::oilPhaseIdx };
69 enum { waterPhaseIdx = FluidSystem::waterPhaseIdx };
71 enum { gasCompIdx = FluidSystem::gasCompIdx };
72 enum { oilCompIdx = FluidSystem::oilCompIdx };
73 enum { waterCompIdx = FluidSystem::waterCompIdx };
74 enum { compositionSwitchIdx = Indices::compositionSwitchIdx };
76 static const bool waterEnabled = Indices::waterEnabled;
77 static const bool gasEnabled = Indices::gasEnabled;
78 static const bool oilEnabled = Indices::oilEnabled;
79 static const bool compositionSwitchEnabled = (compositionSwitchIdx >= 0);
101 template <
class LhsEval>
103 const ElementContext& elemCtx,
108 const IntensiveQuantities& intQuants = elemCtx.intensiveQuantities(dofIdx,
timeIdx);
109 const auto& fs = intQuants.fluidState();
114 if (!FluidSystem::phaseIsActive(
phaseIdx)) {
115 if (Indices::numPhases == 3) {
116 unsigned activeCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::solventComponentIndex(
phaseIdx));
125 unsigned activeCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::solventComponentIndex(
phaseIdx));
134 if (
phaseIdx == oilPhaseIdx && FluidSystem::enableDissolvedGas()) {
135 unsigned activeGasCompIdx = Indices::canonicalToActiveComponentIndex(gasCompIdx);
142 if (
phaseIdx == waterPhaseIdx && FluidSystem::enableDissolvedGasInWater()) {
143 unsigned activeGasCompIdx = Indices::canonicalToActiveComponentIndex(gasCompIdx);
150 if (
phaseIdx == gasPhaseIdx && FluidSystem::enableVaporizedOil()) {
151 unsigned activeOilCompIdx = Indices::canonicalToActiveComponentIndex(oilCompIdx);
158 if (
phaseIdx == gasPhaseIdx && FluidSystem::enableVaporizedWater()) {
169 SolventModule::addStorage(
storage, intQuants);
172 ExtboModule::addStorage(
storage, intQuants);
175 PolymerModule::addStorage(
storage, intQuants);
178 EnergyModule::addStorage(
storage, intQuants);
181 FoamModule::addStorage(
storage, intQuants);
184 BrineModule::addStorage(
storage, intQuants);
187 MICPModule::addStorage(
storage, intQuants);
194 const ElementContext& elemCtx,
205 if (!FluidSystem::phaseIsActive(
phaseIdx))
209 const IntensiveQuantities&
up = elemCtx.intensiveQuantities(
upIdx,
timeIdx);
210 unsigned pvtRegionIdx =
up.pvtRegionIndex();
247 const ElementContext& elemCtx,
252 elemCtx.problem().source(source, elemCtx, dofIdx,
timeIdx);
255 MICPModule::addSource(source, elemCtx, dofIdx,
timeIdx);
266 template <
class UpEval,
class Flu
idState>
269 unsigned pvtRegionIdx,
271 const FluidState&
upFs)
275 unsigned activeCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::solventComponentIndex(
phaseIdx));
277 if (blackoilConserveSurfaceVolume)
284 if (FluidSystem::enableDissolvedGas()) {
285 const auto& Rs = BlackOil::getRs_<FluidSystem, FluidState, UpEval>(
upFs, pvtRegionIdx);
287 unsigned activeGasCompIdx = Indices::canonicalToActiveComponentIndex(gasCompIdx);
288 if (blackoilConserveSurfaceVolume)
293 }
else if (
phaseIdx == waterPhaseIdx) {
295 if (FluidSystem::enableDissolvedGasInWater()) {
296 const auto& Rsw = BlackOil::getRsw_<FluidSystem, FluidState, UpEval>(
upFs, pvtRegionIdx);
298 unsigned activeGasCompIdx = Indices::canonicalToActiveComponentIndex(gasCompIdx);
299 if (blackoilConserveSurfaceVolume)
307 if (FluidSystem::enableVaporizedOil()) {
308 const auto& Rv = BlackOil::getRv_<FluidSystem, FluidState, UpEval>(
upFs, pvtRegionIdx);
310 unsigned activeOilCompIdx = Indices::canonicalToActiveComponentIndex(oilCompIdx);
311 if (blackoilConserveSurfaceVolume)
317 if (FluidSystem::enableVaporizedWater()) {
318 const auto& Rvw = BlackOil::getRvw_<FluidSystem, FluidState, UpEval>(
upFs, pvtRegionIdx);
321 if (blackoilConserveSurfaceVolume)
340 template <
class Scalar>
343 if (blackoilConserveSurfaceVolume)
353 FluidSystem::referenceDensity(waterPhaseIdx, pvtRegionIdx);
357 unsigned activeGasCompIdx = Indices::canonicalToActiveComponentIndex(gasCompIdx);
359 FluidSystem::referenceDensity(gasPhaseIdx, pvtRegionIdx);
363 unsigned activeOilCompIdx = Indices::canonicalToActiveComponentIndex(oilCompIdx);
365 FluidSystem::referenceDensity(oilPhaseIdx, pvtRegionIdx);
Contains the classes required to extend the black-oil model by brine.
Classes required for dynamic convective mixing.
Classes required for molecular diffusion.
Contains the classes required to extend the black-oil model by energy.
Contains the classes required to extend the black-oil model by solvent component.
Contains the classes required to extend the black-oil model to include the effects of foam.
Contains the classes required to extend the black-oil model by MICP.
Contains the classes required to extend the black-oil model by polymer.
Declares the properties required by the black oil model.
Contains the classes required to extend the black-oil model by solvents.
Contains the high level supplements required to extend the black oil model by brine.
Definition blackoilbrinemodules.hh:50
Definition blackoilconvectivemixingmodule.hh:58
Provides the auxiliary methods required for consideration of the diffusion equation.
Definition blackoildiffusionmodule.hh:48
Contains the high level supplements required to extend the black oil model by energy.
Definition blackoilenergymodules.hh:52
Contains the high level supplements required to extend the black oil model.
Definition blackoilextbomodules.hh:59
Contains the high level supplements required to extend the black oil model to include the effects of ...
Definition blackoilfoammodules.hh:60
Calculates the local residual of the black oil model.
Definition blackoillocalresidual.hh:51
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:193
void computeSource(RateVector &source, const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx) const
Calculate the source term of the equation.
Definition blackoillocalresidual.hh:246
static void evalPhaseFluxes_(RateVector &flux, unsigned phaseIdx, unsigned pvtRegionIdx, const ExtensiveQuantities &extQuants, const FluidState &upFs)
Helper function to calculate the flux of mass in terms of conservation quantities via specific fluid ...
Definition blackoillocalresidual.hh:267
static void adaptMassConservationQuantities_(Dune::FieldVector< Scalar, numEq > &container, unsigned pvtRegionIdx)
Helper function to convert the mass-related parts of a Dune::FieldVector that stores conservation qua...
Definition blackoillocalresidual.hh:341
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:102
Contains the high level supplements required to extend the black oil model by MICP.
Definition blackoilmicpmodules.hh:49
Contains the high level supplements required to extend the black oil model by polymer.
Definition blackoilpolymermodules.hh:54
Contains the high level supplements required to extend the black oil model by solvents.
Definition blackoilsolventmodules.hh:58
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242
typename Properties::Detail::GetPropImpl< TypeTag, Property >::type::type GetPropType
get the type alias defined in the property (equivalent to old macro GET_PROP_TYPE(....
Definition propertysystem.hh:235