28 #ifndef EWOMS_RICHARDS_LOCAL_RESIDUAL_HH 29 #define EWOMS_RICHARDS_LOCAL_RESIDUAL_HH 41 template <
class TypeTag>
44 typedef typename GET_PROP_TYPE(TypeTag, EqVector) EqVector;
45 typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
46 typedef typename GET_PROP_TYPE(TypeTag, RateVector) RateVector;
47 typedef typename GET_PROP_TYPE(TypeTag, IntensiveQuantities) IntensiveQuantities;
48 typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
49 typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
51 enum { contiEqIdx = Indices::contiEqIdx };
52 enum { liquidPhaseIdx =
GET_PROP_VALUE(TypeTag, LiquidPhaseIndex) };
54 typedef Opm::MathToolbox<Evaluation> Toolbox;
60 template <
class LhsEval>
62 const ElementContext& elemCtx,
64 unsigned timeIdx)
const 66 const IntensiveQuantities& intQuants = elemCtx.intensiveQuantities(dofIdx, timeIdx);
70 Toolbox::template decay<LhsEval>(intQuants.fluidState().density(liquidPhaseIdx))
71 *Toolbox::template decay<LhsEval>(intQuants.fluidState().saturation(liquidPhaseIdx))
72 *Toolbox::template decay<LhsEval>(intQuants.porosity());
79 const ElementContext& elemCtx,
81 unsigned timeIdx)
const 83 const auto& extQuants = elemCtx.extensiveQuantities(scvfIdx, timeIdx);
85 unsigned focusDofIdx = elemCtx.focusDofIndex();
86 unsigned upIdx =
static_cast<unsigned>(extQuants.upstreamIndex(liquidPhaseIdx));
88 const IntensiveQuantities& up = elemCtx.intensiveQuantities(upIdx, timeIdx);
92 const Evaluation& rho = up.fluidState().density(liquidPhaseIdx);
93 if (focusDofIdx == upIdx)
94 flux[contiEqIdx] = extQuants.volumeFlux(liquidPhaseIdx)*rho;
96 flux[contiEqIdx] = extQuants.volumeFlux(liquidPhaseIdx)*Toolbox::value(rho);
103 const ElementContext& elemCtx,
105 unsigned timeIdx)
const 106 { elemCtx.problem().source(source, elemCtx, dofIdx, timeIdx); }
Intensive quantities required by the Richards model.
Definition: baseauxiliarymodule.hh:37
void computeSource(RateVector &source, const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx) const
Calculate the source term of the equation.
Definition: richardslocalresidual.hh:102
#define GET_PROP_VALUE(TypeTag, PropTagName)
Access the value attribute of a property for a type tag.
Definition: propertysystem.hh:469
#define GET_PROP_TYPE(TypeTag, PropTagName)
Access the type attribute of a property for a type tag.
Definition: propertysystem.hh:486
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: richardslocalresidual.hh:78
void computeStorage(Dune::FieldVector< LhsEval, numEq > &storage, const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx) const
Evaluate the amount all conservation quantities (e.g.
Definition: richardslocalresidual.hh:61
Element-wise calculation of the residual for the Richards model.
Definition: richardslocalresidual.hh:42
Calculates and stores the data which is required to calculate the flux of fluid over a face of a fini...