28 #ifndef EWOMS_FLASH_LOCAL_RESIDUAL_HH 29 #define EWOMS_FLASH_LOCAL_RESIDUAL_HH 36 #include <opm/common/Valgrind.hpp> 45 template <
class TypeTag>
48 typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
49 typedef typename GET_PROP_TYPE(TypeTag, EqVector) EqVector;
50 typedef typename GET_PROP_TYPE(TypeTag, RateVector) RateVector;
51 typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
52 typedef typename GET_PROP_TYPE(TypeTag, IntensiveQuantities) IntensiveQuantities;
53 typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
58 enum { conti0EqIdx = Indices::conti0EqIdx };
60 enum { enableDiffusion =
GET_PROP_VALUE(TypeTag, EnableDiffusion) };
66 typedef Opm::MathToolbox<Evaluation> Toolbox;
72 template <
class LhsEval>
74 const ElementContext& elemCtx,
77 unsigned phaseIdx)
const 79 const IntensiveQuantities& intQuants = elemCtx.intensiveQuantities(dofIdx, timeIdx);
80 const auto& fs = intQuants.fluidState();
83 for (
unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
84 unsigned eqIdx = conti0EqIdx + compIdx;
86 Toolbox::template decay<LhsEval>(fs.molarity(phaseIdx, compIdx))
87 * Toolbox::template decay<LhsEval>(fs.saturation(phaseIdx))
88 * Toolbox::template decay<LhsEval>(intQuants.porosity());
91 EnergyModule::addPhaseStorage(storage, elemCtx.intensiveQuantities(dofIdx, timeIdx), phaseIdx);
97 template <
class LhsEval>
99 const ElementContext& elemCtx,
101 unsigned timeIdx)
const 104 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
107 EnergyModule::addSolidHeatStorage(storage, elemCtx.intensiveQuantities(dofIdx, timeIdx));
114 const ElementContext& elemCtx,
116 unsigned timeIdx)
const 120 Opm::Valgrind::CheckDefined(flux);
123 Opm::Valgrind::CheckDefined(flux);
130 const ElementContext& elemCtx,
132 unsigned timeIdx)
const 134 const auto& extQuants = elemCtx.extensiveQuantities(scvfIdx, timeIdx);
136 unsigned focusDofIdx = elemCtx.focusDofIndex();
137 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
139 unsigned upIdx =
static_cast<unsigned>(extQuants.upstreamIndex(phaseIdx));
140 const IntensiveQuantities& up = elemCtx.intensiveQuantities(upIdx, timeIdx);
145 if (upIdx == focusDofIdx) {
147 up.fluidState().molarDensity(phaseIdx)
148 * extQuants.volumeFlux(phaseIdx);
150 for (
unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
151 flux[conti0EqIdx + compIdx] +=
152 tmp*up.fluidState().moleFraction(phaseIdx, compIdx);
157 Toolbox::value(up.fluidState().molarDensity(phaseIdx))
158 * extQuants.volumeFlux(phaseIdx);
160 for (
unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
161 flux[conti0EqIdx + compIdx] +=
162 tmp*Toolbox::value(up.fluidState().moleFraction(phaseIdx, compIdx));
167 EnergyModule::addAdvectiveFlux(flux, elemCtx, scvfIdx, timeIdx);
174 const ElementContext& elemCtx,
176 unsigned timeIdx)
const 178 DiffusionModule::addDiffusiveFlux(flux, elemCtx, scvfIdx, timeIdx);
179 EnergyModule::addDiffusiveFlux(flux, elemCtx, scvfIdx, timeIdx);
186 const ElementContext& elemCtx,
188 unsigned timeIdx)
const 190 Opm::Valgrind::SetUndefined(source);
191 elemCtx.problem().source(source, elemCtx, dofIdx, timeIdx);
192 Opm::Valgrind::CheckDefined(source);
void computeSource(RateVector &source, const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx) const
Calculate the source term of the equation.
Definition: flashlocalresidual.hh:185
void addDiffusiveFlux(RateVector &flux, const ElementContext &elemCtx, unsigned scvfIdx, unsigned timeIdx) const
Adds the diffusive flux at a given flux integration point.
Definition: flashlocalresidual.hh:173
Definition: baseauxiliarymodule.hh:37
Provides the auxiliary methods required for consideration of the energy equation. ...
Definition: energymodule.hh:59
#define GET_PROP_VALUE(TypeTag, PropTagName)
Access the value attribute of a property for a type tag.
Definition: propertysystem.hh:469
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: flashlocalresidual.hh:113
Contains the classes required to consider energy as a conservation quantity in a multi-phase module...
#define GET_PROP_TYPE(TypeTag, PropTagName)
Access the type attribute of a property for a type tag.
Definition: propertysystem.hh:486
Calculates the local residual of the compositional multi-phase model based on flash calculations...
Definition: flashlocalresidual.hh:46
void computeStorage(Dune::FieldVector< LhsEval, numEq > &storage, const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx) const
Evaluate the amount all conservation quantities (e.g.
Definition: flashlocalresidual.hh:98
Declares the properties required by the compositional multi-phase model based on flash calculations...
Provides the auxiliary methods required for consideration of the diffusion equation.
Definition: diffusionmodule.hh:53
void addPhaseStorage(Dune::FieldVector< LhsEval, numEq > &storage, const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx, unsigned phaseIdx) const
Adds the amount all conservation quantities (e.g.
Definition: flashlocalresidual.hh:73
void addAdvectiveFlux(RateVector &flux, const ElementContext &elemCtx, unsigned scvfIdx, unsigned timeIdx) const
Add the advective mass flux at a given flux integration point.
Definition: flashlocalresidual.hh:129
Classes required for molecular diffusion.