28#ifndef EWOMS_RICHARDS_LOCAL_RESIDUAL_HH
29#define EWOMS_RICHARDS_LOCAL_RESIDUAL_HH
31#include <dune/common/fvector.hh>
33#include <opm/material/common/MathToolbox.hpp>
45template <
class TypeTag>
55 enum { contiEqIdx = Indices::contiEqIdx };
56 enum { liquidPhaseIdx = getPropValue<TypeTag, Properties::LiquidPhaseIndex>() };
57 enum { numEq = getPropValue<TypeTag, Properties::NumEq>() };
58 using Toolbox = MathToolbox<Evaluation>;
64 template <
class LhsEval>
66 const ElementContext& elemCtx,
68 unsigned timeIdx)
const
70 const IntensiveQuantities& intQuants = elemCtx.intensiveQuantities(dofIdx, timeIdx);
74 Toolbox::template decay<LhsEval>(intQuants.fluidState().density(liquidPhaseIdx)) *
75 Toolbox::template decay<LhsEval>(intQuants.fluidState().saturation(liquidPhaseIdx)) *
76 Toolbox::template decay<LhsEval>(intQuants.porosity());
83 const ElementContext& elemCtx,
85 unsigned timeIdx)
const
87 const auto& extQuants = elemCtx.extensiveQuantities(scvfIdx, timeIdx);
89 const unsigned focusDofIdx = elemCtx.focusDofIndex();
90 const unsigned upIdx =
static_cast<unsigned>(extQuants.upstreamIndex(liquidPhaseIdx));
92 const IntensiveQuantities& up = elemCtx.intensiveQuantities(upIdx, timeIdx);
96 const Evaluation& rho = up.fluidState().density(liquidPhaseIdx);
97 if (focusDofIdx == upIdx) {
98 flux[contiEqIdx] = extQuants.volumeFlux(liquidPhaseIdx) * rho;
101 flux[contiEqIdx] = extQuants.volumeFlux(liquidPhaseIdx) * Toolbox::value(rho);
109 const ElementContext& elemCtx,
111 unsigned timeIdx)
const
112 { elemCtx.problem().source(source, elemCtx, dofIdx, timeIdx); }
Element-wise calculation of the residual for the Richards model.
Definition: richardslocalresidual.hh:47
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:82
void computeStorage(Dune::FieldVector< LhsEval, numEq > &storage, const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx) const
Evaluate the amount all conservation quantities (e.g. phase mass) within a finite sub-control volume.
Definition: richardslocalresidual.hh:65
void computeSource(RateVector &source, const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx) const
Calculate the source term of the equation.
Definition: richardslocalresidual.hh:108
Declare the properties used by the infrastructure code of the finite volume discretizations.
Defines the common properties required by the porous medium multi-phase models.
Definition: blackoilboundaryratevector.hh:39
typename Properties::Detail::GetPropImpl< TypeTag, Property >::type::type GetPropType
get the type alias defined in the property (equivalent to old macro GET_PROP_TYPE(....
Definition: propertysystem.hh:233
Contains the property declarations for the Richards model.