26 #ifndef EWOMS_FLASH_LOCAL_RESIDUAL_HH
27 #define EWOMS_FLASH_LOCAL_RESIDUAL_HH
41 template <
class TypeTag>
44 typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
45 typedef typename GET_PROP_TYPE(TypeTag, EqVector) EqVector;
46 typedef typename GET_PROP_TYPE(TypeTag, RateVector) RateVector;
47 typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
48 typedef typename GET_PROP_TYPE(TypeTag, IntensiveQuantities) IntensiveQuantities;
49 typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
54 enum { conti0EqIdx = Indices::conti0EqIdx };
56 enum { enableDiffusion =
GET_PROP_VALUE(TypeTag, EnableDiffusion) };
62 typedef Opm::MathToolbox<Evaluation> Toolbox;
68 template <
class LhsEval>
70 const ElementContext &elemCtx,
75 const IntensiveQuantities &intQuants = elemCtx.intensiveQuantities(dofIdx, timeIdx);
76 const auto &fs = intQuants.fluidState();
79 for (
int compIdx = 0; compIdx < numComponents; ++compIdx) {
80 int eqIdx = conti0EqIdx + compIdx;
82 Toolbox::template toLhs<LhsEval>(fs.molarity(phaseIdx, compIdx))
83 * Toolbox::template toLhs<LhsEval>(fs.saturation(phaseIdx))
84 * Toolbox::template toLhs<LhsEval>(intQuants.porosity());
87 EnergyModule::addPhaseStorage(storage, elemCtx.intensiveQuantities(dofIdx, timeIdx), phaseIdx);
93 template <
class LhsEval>
95 const ElementContext &elemCtx,
100 for (
int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
103 EnergyModule::addSolidHeatStorage(storage, elemCtx.intensiveQuantities(dofIdx, timeIdx));
110 const ElementContext &elemCtx,
116 Valgrind::CheckDefined(flux);
119 Valgrind::CheckDefined(flux);
126 const ElementContext &elemCtx,
130 const auto &extQuants = elemCtx.extensiveQuantities(scvfIdx, timeIdx);
132 int interiorIdx = extQuants.interiorIndex();
133 for (
int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
135 int upIdx = extQuants.upstreamIndex(phaseIdx);
136 const IntensiveQuantities &up = elemCtx.intensiveQuantities(upIdx, timeIdx);
141 if (upIdx == interiorIdx) {
143 up.fluidState().molarDensity(phaseIdx)
144 * extQuants.volumeFlux(phaseIdx);
146 for (
int compIdx = 0; compIdx < numComponents; ++compIdx) {
147 flux[conti0EqIdx + compIdx] +=
148 tmp*up.fluidState().moleFraction(phaseIdx, compIdx);
153 Toolbox::value(up.fluidState().molarDensity(phaseIdx))
154 * extQuants.volumeFlux(phaseIdx);
156 for (
int compIdx = 0; compIdx < numComponents; ++compIdx) {
157 flux[conti0EqIdx + compIdx] +=
158 tmp*Toolbox::value(up.fluidState().moleFraction(phaseIdx, compIdx));
163 EnergyModule::addAdvectiveFlux(flux, elemCtx, scvfIdx, timeIdx);
170 const ElementContext &elemCtx,
174 DiffusionModule::addDiffusiveFlux(flux, elemCtx, scvfIdx, timeIdx);
175 EnergyModule::addDiffusiveFlux(flux, elemCtx, scvfIdx, timeIdx);
182 const ElementContext &elemCtx,
186 Valgrind::SetUndefined(source);
187 elemCtx.problem().source(source, elemCtx, dofIdx, timeIdx);
188 Valgrind::CheckDefined(source);
Classes required for molecular diffusion.
Provides the auxiliary methods required for consideration of the diffusion equation.
Definition: diffusionmodule.hh:46
void addAdvectiveFlux(RateVector &flux, const ElementContext &elemCtx, int scvfIdx, int timeIdx) const
Add the advective mass flux at a given flux integration point.
Definition: flashlocalresidual.hh:125
void computeSource(RateVector &source, const ElementContext &elemCtx, int dofIdx, int timeIdx) const
Calculate the source term of the equation.
Definition: flashlocalresidual.hh:181
#define GET_PROP_VALUE(TypeTag, PropTagName)
Access the value attribute of a property for a type tag.
Definition: propertysystem.hh:468
void computeStorage(Dune::FieldVector< LhsEval, numEq > &storage, const ElementContext &elemCtx, int dofIdx, int timeIdx) const
Evaluate the amount all conservation quantities (e.g. phase mass) within a finite sub-control volume...
Definition: flashlocalresidual.hh:94
Calculates the local residual of the compositional multi-phase model based on flash calculations...
Definition: flashlocalresidual.hh:42
#define GET_PROP_TYPE(TypeTag, PropTagName)
Access the type attribute of a property for a type tag.
Definition: propertysystem.hh:485
void computeFlux(RateVector &flux, const ElementContext &elemCtx, int scvfIdx, int timeIdx) const
Evaluates the total mass flux of all conservation quantities over a face of a sub-control volume...
Definition: flashlocalresidual.hh:109
void addDiffusiveFlux(RateVector &flux, const ElementContext &elemCtx, int scvfIdx, int timeIdx) const
Adds the diffusive flux at a given flux integration point.
Definition: flashlocalresidual.hh:169
Provides the auxiliary methods required for consideration of the energy equation. ...
Definition: energymodule.hh:54
Definition: baseauxiliarymodule.hh:35
void addPhaseStorage(Dune::FieldVector< LhsEval, numEq > &storage, const ElementContext &elemCtx, int dofIdx, int timeIdx, int phaseIdx) const
Adds the amount all conservation quantities (e.g. phase mass) within a single fluid phase...
Definition: flashlocalresidual.hh:69
Declares the properties required by the compositional multi-phase model based on flash calculations...
Contains the classes required to consider energy as a conservation quantity in a multi-phase module...