26 #ifndef EWOMS_RICHARDS_INTENSIVE_QUANTITIES_HH
27 #define EWOMS_RICHARDS_INTENSIVE_QUANTITIES_HH
31 #include <opm/material/fluidstates/ImmiscibleFluidState.hpp>
33 #include <dune/common/fvector.hh>
34 #include <dune/common/fmatrix.hh>
44 template <
class TypeTag>
47 ,
public GET_PROP_TYPE(TypeTag, FluxModule)::FluxIntensiveQuantities
49 typedef typename GET_PROP_TYPE(TypeTag, DiscIntensiveQuantities) ParentType;
50 typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
51 typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
52 typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
53 typedef typename GET_PROP_TYPE(TypeTag, MaterialLaw) MaterialLaw;
54 typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
55 typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
56 typedef typename GET_PROP_TYPE(TypeTag, FluxModule) FluxModule;
58 typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
59 enum { pressureWIdx = Indices::pressureWIdx };
60 enum { numPhases = FluidSystem::numPhases };
61 enum { liquidPhaseIdx =
GET_PROP_VALUE(TypeTag, LiquidPhaseIndex) };
63 enum { dimWorld = GridView::dimensionworld };
65 typedef typename FluxModule::FluxIntensiveQuantities FluxIntensiveQuantities;
66 typedef Dune::FieldMatrix<Scalar, dimWorld, dimWorld> DimMatrix;
67 typedef Dune::FieldVector<Scalar, numPhases> ScalarPhaseVector;
68 typedef Dune::FieldVector<Evaluation, numPhases> PhaseVector;
69 typedef Opm::MathToolbox<Evaluation> Toolbox;
73 typedef Opm::ImmiscibleFluidState<Evaluation, FluidSystem>
FluidState;
78 void update(
const ElementContext &elemCtx,
int dofIdx,
int timeIdx)
80 ParentType::update(elemCtx, dofIdx, timeIdx);
82 const auto& T = elemCtx.problem().temperature(elemCtx, dofIdx, timeIdx);
83 fluidState_.setTemperature(T);
86 typedef typename GET_PROP_TYPE(TypeTag, MaterialLaw) MaterialLaw;
87 const auto &problem = elemCtx.problem();
88 const typename MaterialLaw::Params &materialParams =
89 problem.materialLawParams(elemCtx, dofIdx, timeIdx);
90 const auto &priVars = elemCtx.primaryVars(dofIdx, timeIdx);
97 fluidState_.setSaturation(liquidPhaseIdx, 1.0);
98 fluidState_.setSaturation(gasPhaseIdx, 0.0);
100 MaterialLaw::capillaryPressures(pC, materialParams, fluidState_);
105 const Evaluation& pW = priVars.makeEvaluation(pressureWIdx, timeIdx);
107 Toolbox::max(elemCtx.problem().referencePressure(elemCtx, dofIdx, 0),
108 pW + (pC[gasPhaseIdx] - pC[liquidPhaseIdx]));
113 fluidState_.setPressure(liquidPhaseIdx, pW);
114 fluidState_.setPressure(gasPhaseIdx, pN);
117 MaterialLaw::saturations(sat, materialParams, fluidState_);
118 fluidState_.setSaturation(liquidPhaseIdx, sat[liquidPhaseIdx]);
119 fluidState_.setSaturation(gasPhaseIdx, sat[gasPhaseIdx]);
121 typename FluidSystem::ParameterCache paramCache;
122 paramCache.updateAll(fluidState_);
125 const Evaluation& mu = FluidSystem::viscosity(fluidState_, paramCache, liquidPhaseIdx);
126 fluidState_.setViscosity(liquidPhaseIdx, mu);
127 fluidState_.setViscosity(gasPhaseIdx, 1e-20);
130 const Evaluation& rho = FluidSystem::density(fluidState_, paramCache, liquidPhaseIdx);
131 fluidState_.setDensity(liquidPhaseIdx, rho);
132 fluidState_.setDensity(gasPhaseIdx, 1e-20);
135 MaterialLaw::relativePermeabilities(relativePermeability_, materialParams, fluidState_);
138 for (
int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
139 mobility_[phaseIdx] = relativePermeability_[phaseIdx]/fluidState_.viscosity(phaseIdx);
142 porosity_ = problem.porosity(elemCtx, dofIdx, timeIdx);
145 intrinsicPerm_ = problem.intrinsicPermeability(elemCtx, dofIdx, timeIdx);
148 FluxIntensiveQuantities::update_(elemCtx, dofIdx, timeIdx);
155 {
return fluidState_; }
161 {
return porosity_; }
167 {
return intrinsicPerm_; }
173 {
return relativePermeability_[phaseIdx]; }
179 {
return mobility_[phaseIdx]; }
182 FluidState fluidState_;
183 DimMatrix intrinsicPerm_;
184 Evaluation relativePermeability_[numPhases];
185 Evaluation mobility_[numPhases];
186 Evaluation porosity_;
void update(const ElementContext &elemCtx, int dofIdx, int timeIdx)
Definition: richardsintensivequantities.hh:78
Opm::ImmiscibleFluidState< Evaluation, FluidSystem > FluidState
The type returned by the fluidState() method.
Definition: richardsintensivequantities.hh:73
#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
const Evaluation & porosity() const
Returns the average porosity within the control volume.
Definition: richardsintensivequantities.hh:160
const Evaluation & mobility(int phaseIdx) const
Returns the effective mobility of a given phase within the control volume.
Definition: richardsintensivequantities.hh:178
const DimMatrix & intrinsicPermeability() const
Returns the intrinsic permeability tensor a degree of freedom.
Definition: richardsintensivequantities.hh:166
Intensive quantities required by the Richards model.
Definition: richardsintensivequantities.hh:45
Definition: baseauxiliarymodule.hh:35
const Evaluation & relativePermeability(int phaseIdx) const
Returns the relative permeability of a given phase within the control volume.
Definition: richardsintensivequantities.hh:172
const FluidState & fluidState() const
Returns the phase state for the control-volume.
Definition: richardsintensivequantities.hh:154
Contains the property declarations for the Richards model.