27 #ifndef EWOMS_NCP_LOCAL_RESIDUAL_HH
28 #define EWOMS_NCP_LOCAL_RESIDUAL_HH
42 template <
class TypeTag>
45 typedef typename GET_PROP_TYPE(TypeTag, DiscLocalResidual) ParentType;
46 typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
47 typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
48 typedef typename GET_PROP_TYPE(TypeTag, EqVector) EqVector;
49 typedef typename GET_PROP_TYPE(TypeTag, RateVector) RateVector;
50 typedef typename GET_PROP_TYPE(TypeTag, IntensiveQuantities) IntensiveQuantities;
51 typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
52 typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
57 enum { ncp0EqIdx = Indices::ncp0EqIdx };
58 enum { conti0EqIdx = Indices::conti0EqIdx };
60 enum { enableDiffusion =
GET_PROP_VALUE(TypeTag, EnableDiffusion) };
66 typedef Dune::FieldVector<Evaluation, numEq> EvalEqVector;
67 typedef Dune::BlockVector<EvalEqVector> ElemEvalEqVector;
68 typedef Opm::MathToolbox<Evaluation> Toolbox;
74 template <
class LhsEval>
76 const ElementContext &elemCtx,
81 const IntensiveQuantities &intQuants = elemCtx.intensiveQuantities(dofIdx, timeIdx);
82 const auto &fluidState = intQuants.fluidState();
85 for (
int compIdx = 0; compIdx < numComponents; ++compIdx) {
86 int eqIdx = conti0EqIdx + compIdx;
88 Toolbox::template toLhs<LhsEval>(fluidState.molarity(phaseIdx, compIdx))
89 * Toolbox::template toLhs<LhsEval>(fluidState.saturation(phaseIdx))
90 * Toolbox::template toLhs<LhsEval>(intQuants.porosity());
93 EnergyModule::addPhaseStorage(storage, elemCtx.intensiveQuantities(dofIdx, timeIdx), phaseIdx);
99 template <
class LhsEval>
101 const ElementContext &elemCtx,
106 for (
int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
109 EnergyModule::addSolidHeatStorage(storage, elemCtx.intensiveQuantities(dofIdx, timeIdx));
116 int scvfIdx,
int timeIdx)
const
118 flux = Toolbox::createConstant(0.0);
120 Valgrind::CheckDefined(flux);
123 Valgrind::CheckDefined(flux);
130 int scvfIdx,
int timeIdx)
const
132 const auto &extQuants = elemCtx.extensiveQuantities(scvfIdx, timeIdx);
134 int interiorIdx = extQuants.interiorIndex();
135 for (
int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
138 int upIdx = extQuants.upstreamIndex(phaseIdx);
139 const IntensiveQuantities &up = elemCtx.intensiveQuantities(upIdx, timeIdx);
144 if (upIdx == interiorIdx) {
146 up.fluidState().molarDensity(phaseIdx)
147 * extQuants.volumeFlux(phaseIdx);
149 for (
int compIdx = 0; compIdx < numComponents; ++compIdx) {
150 flux[conti0EqIdx + compIdx] +=
151 tmp*up.fluidState().moleFraction(phaseIdx, compIdx);
156 Toolbox::value(up.fluidState().molarDensity(phaseIdx))
157 * extQuants.volumeFlux(phaseIdx);
159 for (
int compIdx = 0; compIdx < numComponents; ++compIdx) {
160 flux[conti0EqIdx + compIdx] +=
161 tmp*Toolbox::value(up.fluidState().moleFraction(phaseIdx, compIdx));
166 EnergyModule::addAdvectiveFlux(flux, elemCtx, scvfIdx, timeIdx);
173 int scvfIdx,
int timeIdx)
const
175 DiffusionModule::addDiffusiveFlux(flux, elemCtx, scvfIdx, timeIdx);
176 EnergyModule::addDiffusiveFlux(flux, elemCtx, scvfIdx, timeIdx);
186 const ElementContext &elemCtx,
190 Valgrind::SetUndefined(source);
191 elemCtx.problem().source(source, elemCtx, dofIdx, timeIdx);
192 Valgrind::CheckDefined(source);
200 ElemEvalEqVector& storageTerm,
201 const ElementContext& elemCtx,
int timeIdx)
const
204 for (
int dofIdx = 0; dofIdx < elemCtx.numPrimaryDof(0); ++dofIdx) {
205 for (
int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
206 residual[dofIdx][ncp0EqIdx + phaseIdx] =
207 phaseNcp_(elemCtx, dofIdx, timeIdx, phaseIdx);
212 ParentType::evalConstraints_(residual, storageTerm, elemCtx, timeIdx);
223 const auto &fluidState = elemCtx.intensiveQuantities(dofIdx, timeIdx).fluidState();
227 return Toolbox::min(a, b);
234 template <
class Flu
idState>
236 {
return fluidState.saturation(phaseIdx); }
242 template <
class Flu
idState>
247 for (
int i = 0; i < numComponents; ++i)
248 a -= fluidState.moleFraction(phaseIdx, i);
Evaluation phasePresentIneq_(const FluidState &fluidState, int phaseIdx) const
Returns the value of the inequality where a phase is present.
Definition: ncplocalresidual.hh:235
void addAdvectiveFlux(RateVector &flux, const ElementContext &elemCtx, int scvfIdx, int timeIdx) const
Add the advective mass flux at a given flux integration point.
Definition: ncplocalresidual.hh:129
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: ncplocalresidual.hh:115
Classes required for molecular diffusion.
Provides the auxiliary methods required for consideration of the diffusion equation.
Definition: diffusionmodule.hh:46
Evaluation phaseNotPresentIneq_(const FluidState &fluidState, int phaseIdx) const
Returns the value of the inequality where a phase is not present.
Definition: ncplocalresidual.hh:243
void evalConstraints_(ElemEvalEqVector &residual, ElemEvalEqVector &storageTerm, const ElementContext &elemCtx, int timeIdx) const
Set the values of the constraint volumes of the current element.
Definition: ncplocalresidual.hh:199
#define GET_PROP_VALUE(TypeTag, PropTagName)
Access the value attribute of a property for a type tag.
Definition: propertysystem.hh:468
#define GET_PROP_TYPE(TypeTag, PropTagName)
Access the type attribute of a property for a type tag.
Definition: propertysystem.hh:485
Evaluation phaseNcp_(const ElementContext &elemCtx, int dofIdx, int timeIdx, int phaseIdx) const
Returns the value of the NCP-function for a phase.
Definition: ncplocalresidual.hh:218
void computeSource(RateVector &source, const ElementContext &elemCtx, int dofIdx, int timeIdx) const
Calculate the source term of the equation.
Definition: ncplocalresidual.hh:185
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: ncplocalresidual.hh:100
Provides the auxiliary methods required for consideration of the energy equation. ...
Definition: energymodule.hh:54
Definition: baseauxiliarymodule.hh:35
void addDiffusiveFlux(RateVector &flux, const ElementContext &elemCtx, int scvfIdx, int timeIdx) const
Adds the diffusive flux at a given flux integration point.
Definition: ncplocalresidual.hh:172
Details needed to calculate the local residual in the compositional multi-phase NCP-model ...
Definition: ncplocalresidual.hh:43
Declares the properties required for the NCP compositional multi-phase model.
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: ncplocalresidual.hh:75
Contains the classes required to consider energy as a conservation quantity in a multi-phase module...