28 #ifndef EWOMS_PVS_INTENSIVE_QUANTITIES_HH
29 #define EWOMS_PVS_INTENSIVE_QUANTITIES_HH
36 #include <opm/material/constraintsolvers/ComputeFromReferencePhase.hpp>
37 #include <opm/material/constraintsolvers/MiscibleMultiPhaseComposition.hpp>
38 #include <opm/material/fluidstates/CompositionalFluidState.hpp>
39 #include <opm/common/Valgrind.hpp>
41 #include <dune/common/fvector.hh>
42 #include <dune/common/fmatrix.hh>
55 template <
class TypeTag>
60 ,
public GET_PROP_TYPE(TypeTag, FluxModule)::FluxIntensiveQuantities
62 typedef typename GET_PROP_TYPE(TypeTag, DiscIntensiveQuantities) ParentType;
64 typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
65 typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
66 typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
67 typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
68 typedef typename GET_PROP_TYPE(TypeTag, MaterialLaw) MaterialLaw;
69 typedef typename GET_PROP_TYPE(TypeTag, MaterialLawParams) MaterialLawParams;
70 typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
71 typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
72 typedef typename GET_PROP_TYPE(TypeTag, FluxModule) FluxModule;
74 enum { switch0Idx = Indices::switch0Idx };
75 enum { pressure0Idx = Indices::pressure0Idx };
78 enum { enableDiffusion =
GET_PROP_VALUE(TypeTag, EnableDiffusion) };
80 enum { dimWorld = GridView::dimensionworld };
82 typedef Opm::MathToolbox<Evaluation> Toolbox;
83 typedef Opm::MiscibleMultiPhaseComposition<Scalar, FluidSystem, Evaluation> MiscibleMultiPhaseComposition;
84 typedef Opm::ComputeFromReferencePhase<Scalar, FluidSystem, Evaluation> ComputeFromReferencePhase;
86 typedef Dune::FieldVector<Scalar, numPhases> PhaseVector;
87 typedef Dune::FieldVector<Evaluation, numPhases> EvalPhaseVector;
88 typedef Dune::FieldMatrix<Scalar, dimWorld, dimWorld> DimMatrix;
90 typedef typename FluxModule::FluxIntensiveQuantities FluxIntensiveQuantities;
96 typedef Opm::CompositionalFluidState<Evaluation, FluidSystem>
FluidState;
108 void update(
const ElementContext& elemCtx,
unsigned dofIdx,
unsigned timeIdx)
110 ParentType::update(elemCtx, dofIdx, timeIdx);
111 EnergyIntensiveQuantities::updateTemperatures_(fluidState_, elemCtx, dofIdx, timeIdx);
113 const auto& priVars = elemCtx.primaryVars(dofIdx, timeIdx);
114 const auto& problem = elemCtx.problem();
119 Evaluation sumSat = 0.0;
120 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
121 fluidState_.setSaturation(phaseIdx, priVars.explicitSaturationValue(phaseIdx, timeIdx));
122 Opm::Valgrind::CheckDefined(fluidState_.saturation(phaseIdx));
123 sumSat += fluidState_.saturation(phaseIdx);
125 Opm::Valgrind::CheckDefined(priVars.implicitSaturationIdx());
126 Opm::Valgrind::CheckDefined(sumSat);
127 fluidState_.setSaturation(priVars.implicitSaturationIdx(), 1.0 - sumSat);
134 const MaterialLawParams& materialParams =
135 problem.materialLawParams(elemCtx, dofIdx, timeIdx);
137 MaterialLaw::capillaryPressures(pC, materialParams, fluidState_);
140 const Evaluation& p0 = priVars.makeEvaluation(pressure0Idx, timeIdx);
141 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
142 fluidState_.setPressure(phaseIdx, p0 + (pC[phaseIdx] - pC[0]));
148 typename FluidSystem::template ParameterCache<Evaluation> paramCache;
149 unsigned lowestPresentPhaseIdx = priVars.lowestPresentPhaseIdx();
150 unsigned numNonPresentPhases = 0;
151 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
152 if (!priVars.phaseIsPresent(phaseIdx))
153 ++numNonPresentPhases;
157 if (numNonPresentPhases == numPhases - 1) {
160 Evaluation sumx = 0.0;
161 for (
unsigned compIdx = 1; compIdx < numComponents; ++compIdx) {
162 const Evaluation& x = priVars.makeEvaluation(switch0Idx + compIdx - 1, timeIdx);
163 fluidState_.setMoleFraction(lowestPresentPhaseIdx, compIdx, x);
168 fluidState_.setMoleFraction(lowestPresentPhaseIdx, 0, 1 - sumx);
173 ComputeFromReferencePhase::solve(fluidState_, paramCache,
174 lowestPresentPhaseIdx,
180 unsigned numAuxConstraints = numComponents + numNonPresentPhases - numPhases;
181 Opm::MMPCAuxConstraint<Evaluation> auxConstraints[numComponents];
184 unsigned switchIdx = 0;
185 for (; switchIdx < numPhases - 1; ++switchIdx) {
186 unsigned compIdx = switchIdx + 1;
187 unsigned switchPhaseIdx = switchIdx;
188 if (switchIdx >= lowestPresentPhaseIdx)
191 if (!priVars.phaseIsPresent(switchPhaseIdx)) {
192 auxConstraints[auxIdx].set(lowestPresentPhaseIdx, compIdx,
193 priVars.makeEvaluation(switch0Idx + switchIdx, timeIdx));
198 for (; auxIdx < numAuxConstraints; ++auxIdx, ++switchIdx) {
199 unsigned compIdx = numPhases - numNonPresentPhases + auxIdx;
200 auxConstraints[auxIdx].set(lowestPresentPhaseIdx, compIdx,
201 priVars.makeEvaluation(switch0Idx + switchIdx, timeIdx));
207 MiscibleMultiPhaseComposition::solve(fluidState_, paramCache,
208 priVars.phasePresence(),
218 Scalar myNan = std::numeric_limits<Scalar>::quiet_NaN();
219 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
220 fluidState_.setEnthalpy(phaseIdx, myNan);
229 MaterialLaw::relativePermeabilities(relativePermeability_,
230 materialParams, fluidState_);
231 Opm::Valgrind::CheckDefined(relativePermeability_);
234 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
235 mobility_[phaseIdx] =
236 relativePermeability_[phaseIdx] /
fluidState().viscosity(phaseIdx);
239 porosity_ = problem.porosity(elemCtx, dofIdx, timeIdx);
240 Opm::Valgrind::CheckDefined(porosity_);
243 intrinsicPerm_ = problem.intrinsicPermeability(elemCtx, dofIdx, timeIdx);
246 FluxIntensiveQuantities::update_(elemCtx, dofIdx, timeIdx);
249 EnergyIntensiveQuantities::update_(fluidState_, paramCache, elemCtx, dofIdx, timeIdx);
252 DiffusionIntensiveQuantities::update_(fluidState_, paramCache, elemCtx, dofIdx, timeIdx);
254 fluidState_.checkDefined();
261 {
return fluidState_; }
267 {
return intrinsicPerm_; }
273 {
return relativePermeability_[phaseIdx]; }
278 const Evaluation&
mobility(
unsigned phaseIdx)
const
279 {
return mobility_[phaseIdx]; }
285 {
return porosity_; }
289 Evaluation porosity_;
290 DimMatrix intrinsicPerm_;
291 Evaluation relativePermeability_[numPhases];
292 Evaluation mobility_[numPhases];
Provides the volumetric quantities required for the calculation of molecular diffusive fluxes...
Definition: diffusionmodule.hh:146
const Evaluation & relativePermeability(unsigned phaseIdx) const
Returns the relative permeability of a given phase within the control volume.
Definition: pvsintensivequantities.hh:272
#define GET_PROP_VALUE(TypeTag, PropTagName)
Access the value attribute of a property for a type tag.
Definition: propertysystem.hh:469
const Evaluation & porosity() const
Returns the average porosity within the control volume.
Definition: pvsintensivequantities.hh:284
void update(const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx)
Definition: pvsintensivequantities.hh:108
Contains the classes required to consider energy as a conservation quantity in a multi-phase module...
Opm::CompositionalFluidState< Evaluation, FluidSystem > FluidState
The type of the object returned by the fluidState() method.
Definition: pvsintensivequantities.hh:96
const DimMatrix & intrinsicPermeability() const
Returns the intrinsic permeability tensor a degree of freedom.
Definition: pvsintensivequantities.hh:266
#define GET_PROP_TYPE(TypeTag, PropTagName)
Access the type attribute of a property for a type tag.
Definition: propertysystem.hh:486
const FluidState & fluidState() const
Returns the phase state for the control-volume.
Definition: pvsintensivequantities.hh:260
const Evaluation & mobility(unsigned phaseIdx) const
Returns the effective mobility of a given phase within the control volume.
Definition: pvsintensivequantities.hh:278
Provides the volumetric quantities required for the energy equation.
Definition: energymodule.hh:539
Contains the quantities which are are constant within a finite volume in the compositional multi-phas...
Definition: pvsintensivequantities.hh:56
Declares the properties required for the compositional multi-phase primary variable switching model...
Classes required for molecular diffusion.