26 #ifndef EWOMS_RICHARDS_LOCAL_RESIDUAL_HH
27 #define EWOMS_RICHARDS_LOCAL_RESIDUAL_HH
39 template <
class TypeTag>
42 typedef typename GET_PROP_TYPE(TypeTag, EqVector) EqVector;
43 typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
44 typedef typename GET_PROP_TYPE(TypeTag, RateVector) RateVector;
45 typedef typename GET_PROP_TYPE(TypeTag, IntensiveQuantities) IntensiveQuantities;
46 typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
47 typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
49 enum { contiEqIdx = Indices::contiEqIdx };
50 enum { liquidPhaseIdx =
GET_PROP_VALUE(TypeTag, LiquidPhaseIndex) };
52 typedef Opm::MathToolbox<Evaluation> Toolbox;
58 template <
class LhsEval>
60 const ElementContext &elemCtx,
64 const IntensiveQuantities &intQuants = elemCtx.intensiveQuantities(dofIdx, timeIdx);
68 Toolbox::template toLhs<LhsEval>(intQuants.fluidState().density(liquidPhaseIdx))
69 *Toolbox::template toLhs<LhsEval>(intQuants.fluidState().saturation(liquidPhaseIdx))
70 *Toolbox::template toLhs<LhsEval>(intQuants.porosity());
76 void computeFlux(RateVector &flux,
const ElementContext &elemCtx,
77 int scvfIdx,
int timeIdx)
const
79 const auto &extQuants = elemCtx.extensiveQuantities(scvfIdx, timeIdx);
81 int interiorIdx = extQuants.interiorIndex();
82 int upIdx = extQuants.upstreamIndex(liquidPhaseIdx);
84 const IntensiveQuantities &up = elemCtx.intensiveQuantities(upIdx, timeIdx);
88 const Evaluation& rho = up.fluidState().density(liquidPhaseIdx);
89 if (interiorIdx == upIdx)
90 flux[contiEqIdx] = extQuants.volumeFlux(liquidPhaseIdx)*rho;
92 flux[contiEqIdx] = extQuants.volumeFlux(liquidPhaseIdx)*Toolbox::value(rho);
99 const ElementContext &elemCtx,
102 { elemCtx.problem().source(source, elemCtx, dofIdx, timeIdx); }
Element-wise calculation of the residual for the Richards model.
Definition: richardslocalresidual.hh:40
#define GET_PROP_VALUE(TypeTag, PropTagName)
Access the value attribute of a property for a type tag.
Definition: propertysystem.hh:468
void computeSource(RateVector &source, const ElementContext &elemCtx, int dofIdx, int timeIdx) const
Calculate the source term of the equation.
Definition: richardslocalresidual.hh:98
#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: richardslocalresidual.hh:76
Intensive quantities required by the Richards model.
Calculates and stores the data which is required to calculate the flux of fluid over a face of a fini...
Definition: baseauxiliarymodule.hh:35
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: richardslocalresidual.hh:59