26 #ifndef EWOMS_NCP_INTENSIVE_QUANTITIES_HH
27 #define EWOMS_NCP_INTENSIVE_QUANTITIES_HH
34 #include <opm/material/constraintsolvers/NcpFlash.hpp>
35 #include <opm/material/fluidstates/CompositionalFluidState.hpp>
36 #include <opm/material/constraintsolvers/CompositionFromFugacities.hpp>
38 #include <dune/common/fvector.hh>
39 #include <dune/common/fmatrix.hh>
49 template <
class TypeTag>
54 ,
public GET_PROP_TYPE(TypeTag, FluxModule)::FluxIntensiveQuantities
56 typedef typename GET_PROP_TYPE(TypeTag, DiscIntensiveQuantities) ParentType;
58 typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
59 typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
60 typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
61 typedef typename GET_PROP_TYPE(TypeTag, MaterialLaw) MaterialLaw;
62 typedef typename GET_PROP_TYPE(TypeTag, MaterialLawParams) MaterialLawParams;
63 typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
64 typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
66 typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
67 typedef typename GET_PROP_TYPE(TypeTag, FluxModule) FluxModule;
72 enum { enableDiffusion =
GET_PROP_VALUE(TypeTag, EnableDiffusion) };
73 enum { fugacity0Idx = Indices::fugacity0Idx };
74 enum { saturation0Idx = Indices::saturation0Idx };
75 enum { pressure0Idx = Indices::pressure0Idx };
76 enum { dimWorld = GridView::dimensionworld };
78 typedef Opm::CompositionFromFugacities<Scalar, FluidSystem, Evaluation>
79 CompositionFromFugacitiesSolver;
80 typedef Opm::CompositionalFluidState<Evaluation, FluidSystem,
81 enableEnergy> FluidState;
82 typedef Dune::FieldVector<Evaluation, numComponents> ComponentVector;
83 typedef Dune::FieldMatrix<Scalar, dimWorld, dimWorld> DimMatrix;
86 typedef typename FluxModule::FluxIntensiveQuantities FluxIntensiveQuantities;
95 void update(
const ElementContext &elemCtx,
99 ParentType::update(elemCtx,
102 ParentType::checkDefined();
104 typename FluidSystem::ParameterCache paramCache;
105 const auto &priVars = elemCtx.primaryVars(dofIdx, timeIdx);
108 Evaluation sumSat = 0;
109 for (
int phaseIdx = 0; phaseIdx < numPhases - 1; ++phaseIdx) {
110 const Evaluation& val = priVars.makeEvaluation(saturation0Idx + phaseIdx, timeIdx);
111 fluidState_.setSaturation(phaseIdx, val);
114 fluidState_.setSaturation(numPhases - 1, 1.0 - sumSat);
115 Valgrind::CheckDefined(sumSat);
118 EnergyIntensiveQuantities::updateTemperatures_(fluidState_, elemCtx, dofIdx, timeIdx);
121 const auto &problem = elemCtx.problem();
122 const MaterialLawParams &materialParams =
123 problem.materialLawParams(elemCtx, dofIdx, timeIdx);
125 Evaluation capPress[numPhases];
126 MaterialLaw::capillaryPressures(capPress, materialParams, fluidState_);
128 const Evaluation& pressure0 = priVars.makeEvaluation(pressure0Idx, timeIdx);
129 for (
int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
130 fluidState_.setPressure(phaseIdx, pressure0 + (capPress[phaseIdx] - capPress[0]));
134 for (
int compIdx = 0; compIdx < numComponents; ++compIdx)
135 fug[compIdx] = priVars.makeEvaluation(fugacity0Idx + compIdx, timeIdx);
138 const auto *hint = elemCtx.thermodynamicHint(dofIdx, timeIdx);
139 for (
int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
142 for (
int compIdx = 0; compIdx < numComponents; ++compIdx) {
144 const Evaluation& moleFracIJ = hint->fluidState().moleFraction(phaseIdx, compIdx);
147 fluidState_.setMoleFraction(phaseIdx, compIdx, moleFracIJ);
151 CompositionFromFugacitiesSolver::guessInitial(fluidState_,
157 CompositionFromFugacitiesSolver::solve(fluidState_, paramCache,
162 porosity_ = problem.porosity(elemCtx, dofIdx, timeIdx);
163 Valgrind::CheckDefined(porosity_);
166 MaterialLaw::relativePermeabilities(relativePermeability_,
167 materialParams, fluidState_);
170 for (
int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
172 const Evaluation& mu = FluidSystem::viscosity(fluidState_, paramCache, phaseIdx);
173 fluidState_.setViscosity(phaseIdx, mu);
175 mobility_[phaseIdx] = relativePermeability_[phaseIdx]/mu;
179 intrinsicPerm_ = problem.intrinsicPermeability(elemCtx, dofIdx, timeIdx);
182 FluxIntensiveQuantities::update_(elemCtx, dofIdx, timeIdx);
185 EnergyIntensiveQuantities::update_(fluidState_, paramCache, elemCtx, dofIdx, timeIdx);
188 DiffusionIntensiveQuantities::update_(fluidState_, paramCache, elemCtx, dofIdx, timeIdx);
197 {
return fluidState_; }
203 {
return intrinsicPerm_; }
209 {
return relativePermeability_[phaseIdx]; }
215 {
return mobility_[phaseIdx]; }
221 {
return porosity_; }
228 #if !defined NDEBUG && HAVE_VALGRIND
229 ParentType::checkDefined();
231 Valgrind::CheckDefined(porosity_);
232 Valgrind::CheckDefined(relativePermeability_);
234 fluidState_.checkDefined();
239 DimMatrix intrinsicPerm_;
240 FluidState fluidState_;
241 Evaluation porosity_;
242 Evaluation relativePermeability_[numPhases];
243 Evaluation mobility_[numPhases];
Provides the volumetric quantities required for the energy equation.
Definition: energymodule.hh:530
Classes required for molecular diffusion.
const Evaluation & mobility(int phaseIdx) const
ImmiscibleIntensiveQuantities::mobility.
Definition: ncpintensivequantities.hh:214
#define GET_PROP_VALUE(TypeTag, PropTagName)
Access the value attribute of a property for a type tag.
Definition: propertysystem.hh:468
NcpIntensiveQuantities()
Definition: ncpintensivequantities.hh:89
#define GET_PROP_TYPE(TypeTag, PropTagName)
Access the type attribute of a property for a type tag.
Definition: propertysystem.hh:485
void checkDefined() const
IntensiveQuantities::checkDefined.
Definition: ncpintensivequantities.hh:226
const Evaluation & porosity() const
ImmiscibleIntensiveQuantities::porosity.
Definition: ncpintensivequantities.hh:220
Contains the quantities which are are constant within a finite volume in the compositional multi-phas...
Definition: ncpintensivequantities.hh:50
const FluidState & fluidState() const
ImmiscibleIntensiveQuantities::fluidState.
Definition: ncpintensivequantities.hh:196
Definition: baseauxiliarymodule.hh:35
void update(const ElementContext &elemCtx, int dofIdx, int timeIdx)
IntensiveQuantities::update.
Definition: ncpintensivequantities.hh:95
const DimMatrix & intrinsicPermeability() const
ImmiscibleIntensiveQuantities::intrinsicPermeability.
Definition: ncpintensivequantities.hh:202
Provides the volumetric quantities required for the calculation of molecular diffusive fluxes...
Definition: diffusionmodule.hh:138
Declares the properties required for the NCP compositional multi-phase model.
const Evaluation & relativePermeability(int phaseIdx) const
ImmiscibleIntensiveQuantities::relativePermeability.
Definition: ncpintensivequantities.hh:208
Contains the classes required to consider energy as a conservation quantity in a multi-phase module...