28 #ifndef EWOMS_RICHARDS_INTENSIVE_QUANTITIES_HH 29 #define EWOMS_RICHARDS_INTENSIVE_QUANTITIES_HH 33 #include <opm/material/fluidstates/ImmiscibleFluidState.hpp> 35 #include <dune/common/fvector.hh> 36 #include <dune/common/fmatrix.hh> 46 template <
class TypeTag>
49 ,
public GET_PROP_TYPE(TypeTag, FluxModule)::FluxIntensiveQuantities
51 typedef typename GET_PROP_TYPE(TypeTag, DiscIntensiveQuantities) ParentType;
52 typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
53 typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
54 typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
55 typedef typename GET_PROP_TYPE(TypeTag, MaterialLaw) MaterialLaw;
56 typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
57 typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
58 typedef typename GET_PROP_TYPE(TypeTag, FluxModule) FluxModule;
60 typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
61 enum { pressureWIdx = Indices::pressureWIdx };
62 enum { numPhases = FluidSystem::numPhases };
63 enum { liquidPhaseIdx =
GET_PROP_VALUE(TypeTag, LiquidPhaseIndex) };
65 enum { dimWorld = GridView::dimensionworld };
67 typedef typename FluxModule::FluxIntensiveQuantities FluxIntensiveQuantities;
68 typedef Dune::FieldMatrix<Scalar, dimWorld, dimWorld> DimMatrix;
69 typedef Dune::FieldVector<Scalar, numPhases> ScalarPhaseVector;
70 typedef Dune::FieldVector<Evaluation, numPhases> PhaseVector;
71 typedef Opm::MathToolbox<Evaluation> Toolbox;
75 typedef Opm::ImmiscibleFluidState<Evaluation, FluidSystem>
FluidState;
87 void update(
const ElementContext& elemCtx,
unsigned dofIdx,
unsigned timeIdx)
89 ParentType::update(elemCtx, dofIdx, timeIdx);
91 const auto& T = elemCtx.problem().temperature(elemCtx, dofIdx, timeIdx);
92 fluidState_.setTemperature(T);
95 typedef typename GET_PROP_TYPE(TypeTag, MaterialLaw) MaterialLaw;
96 const auto& problem = elemCtx.problem();
97 const typename MaterialLaw::Params& materialParams =
98 problem.materialLawParams(elemCtx, dofIdx, timeIdx);
99 const auto& priVars = elemCtx.primaryVars(dofIdx, timeIdx);
106 fluidState_.setSaturation(liquidPhaseIdx, 1.0);
107 fluidState_.setSaturation(gasPhaseIdx, 0.0);
108 ScalarPhaseVector pC;
109 MaterialLaw::capillaryPressures(pC, materialParams, fluidState_);
114 const Evaluation& pW = priVars.makeEvaluation(pressureWIdx, timeIdx);
116 Toolbox::max(elemCtx.problem().referencePressure(elemCtx, dofIdx, 0),
117 pW + (pC[gasPhaseIdx] - pC[liquidPhaseIdx]));
122 fluidState_.setPressure(liquidPhaseIdx, pW);
123 fluidState_.setPressure(gasPhaseIdx, pN);
126 MaterialLaw::saturations(sat, materialParams, fluidState_);
127 fluidState_.setSaturation(liquidPhaseIdx, sat[liquidPhaseIdx]);
128 fluidState_.setSaturation(gasPhaseIdx, sat[gasPhaseIdx]);
130 typename FluidSystem::template ParameterCache<Evaluation> paramCache;
131 paramCache.updateAll(fluidState_);
134 const Evaluation& mu = FluidSystem::viscosity(fluidState_, paramCache, liquidPhaseIdx);
135 fluidState_.setViscosity(liquidPhaseIdx, mu);
136 fluidState_.setViscosity(gasPhaseIdx, 1e-20);
139 const Evaluation& rho = FluidSystem::density(fluidState_, paramCache, liquidPhaseIdx);
140 fluidState_.setDensity(liquidPhaseIdx, rho);
141 fluidState_.setDensity(gasPhaseIdx, 1e-20);
144 MaterialLaw::relativePermeabilities(relativePermeability_, materialParams, fluidState_);
147 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
148 mobility_[phaseIdx] = relativePermeability_[phaseIdx]/fluidState_.viscosity(phaseIdx);
151 porosity_ = problem.porosity(elemCtx, dofIdx, timeIdx);
154 intrinsicPerm_ = problem.intrinsicPermeability(elemCtx, dofIdx, timeIdx);
157 FluxIntensiveQuantities::update_(elemCtx, dofIdx, timeIdx);
164 {
return fluidState_; }
170 {
return porosity_; }
176 {
return intrinsicPerm_; }
182 {
return relativePermeability_[phaseIdx]; }
187 const Evaluation&
mobility(
unsigned phaseIdx)
const 188 {
return mobility_[phaseIdx]; }
192 DimMatrix intrinsicPerm_;
193 Evaluation relativePermeability_[numPhases];
194 Evaluation mobility_[numPhases];
195 Evaluation porosity_;
Definition: baseauxiliarymodule.hh:37
#define GET_PROP_VALUE(TypeTag, PropTagName)
Access the value attribute of a property for a type tag.
Definition: propertysystem.hh:469
const Evaluation & relativePermeability(unsigned phaseIdx) const
Returns the relative permeability of a given phase within the control volume.
Definition: richardsintensivequantities.hh:181
const DimMatrix & intrinsicPermeability() const
Returns the intrinsic permeability tensor a degree of freedom.
Definition: richardsintensivequantities.hh:175
#define GET_PROP_TYPE(TypeTag, PropTagName)
Access the type attribute of a property for a type tag.
Definition: propertysystem.hh:486
const Evaluation & porosity() const
Returns the average porosity within the control volume.
Definition: richardsintensivequantities.hh:169
const Evaluation & mobility(unsigned phaseIdx) const
Returns the effective mobility of a given phase within the control volume.
Definition: richardsintensivequantities.hh:187
Intensive quantities required by the Richards model.
Definition: richardsintensivequantities.hh:47
Contains the property declarations for the Richards model.
Opm::ImmiscibleFluidState< Evaluation, FluidSystem > FluidState
The type returned by the fluidState() method.
Definition: richardsintensivequantities.hh:75
const FluidState & fluidState() const
Returns the phase state for the control-volume.
Definition: richardsintensivequantities.hh:163
void update(const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx)
Definition: richardsintensivequantities.hh:87