28 #ifndef EWOMS_NCP_LOCAL_RESIDUAL_HH 29 #define EWOMS_NCP_LOCAL_RESIDUAL_HH 36 #include <opm/common/Valgrind.hpp> 45 template <
class TypeTag>
48 typedef typename GET_PROP_TYPE(TypeTag, DiscLocalResidual) ParentType;
49 typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
50 typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
51 typedef typename GET_PROP_TYPE(TypeTag, EqVector) EqVector;
52 typedef typename GET_PROP_TYPE(TypeTag, RateVector) RateVector;
53 typedef typename GET_PROP_TYPE(TypeTag, IntensiveQuantities) IntensiveQuantities;
54 typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
55 typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
60 enum { ncp0EqIdx = Indices::ncp0EqIdx };
61 enum { conti0EqIdx = Indices::conti0EqIdx };
63 enum { enableDiffusion =
GET_PROP_VALUE(TypeTag, EnableDiffusion) };
69 typedef Dune::FieldVector<Evaluation, numEq> EvalEqVector;
70 typedef Dune::BlockVector<EvalEqVector> ElemEvalEqVector;
71 typedef Opm::MathToolbox<Evaluation> Toolbox;
77 template <
class LhsEval>
79 const ElementContext& elemCtx,
82 unsigned phaseIdx)
const 84 const IntensiveQuantities& intQuants = elemCtx.intensiveQuantities(dofIdx, timeIdx);
85 const auto& fluidState = intQuants.fluidState();
88 for (
unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
89 unsigned eqIdx = conti0EqIdx + compIdx;
91 Toolbox::template decay<LhsEval>(fluidState.molarity(phaseIdx, compIdx))
92 * Toolbox::template decay<LhsEval>(fluidState.saturation(phaseIdx))
93 * Toolbox::template decay<LhsEval>(intQuants.porosity());
96 EnergyModule::addPhaseStorage(storage, elemCtx.intensiveQuantities(dofIdx, timeIdx), phaseIdx);
102 template <
class LhsEval>
104 const ElementContext& elemCtx,
106 unsigned timeIdx)
const 109 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
112 EnergyModule::addSolidHeatStorage(storage, elemCtx.intensiveQuantities(dofIdx, timeIdx));
119 const ElementContext& elemCtx,
121 unsigned timeIdx)
const 125 Opm::Valgrind::CheckDefined(flux);
128 Opm::Valgrind::CheckDefined(flux);
135 const ElementContext& elemCtx,
137 unsigned timeIdx)
const 139 const auto& extQuants = elemCtx.extensiveQuantities(scvfIdx, timeIdx);
141 unsigned focusDofIdx = elemCtx.focusDofIndex();
142 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
145 unsigned upIdx =
static_cast<unsigned>(extQuants.upstreamIndex(phaseIdx));
146 const IntensiveQuantities& up = elemCtx.intensiveQuantities(upIdx, timeIdx);
151 if (upIdx == focusDofIdx) {
153 up.fluidState().molarDensity(phaseIdx)
154 * extQuants.volumeFlux(phaseIdx);
156 for (
unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
157 flux[conti0EqIdx + compIdx] +=
158 tmp*up.fluidState().moleFraction(phaseIdx, compIdx);
163 Toolbox::value(up.fluidState().molarDensity(phaseIdx))
164 * extQuants.volumeFlux(phaseIdx);
166 for (
unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
167 flux[conti0EqIdx + compIdx] +=
168 tmp*Toolbox::value(up.fluidState().moleFraction(phaseIdx, compIdx));
173 EnergyModule::addAdvectiveFlux(flux, elemCtx, scvfIdx, timeIdx);
180 const ElementContext& elemCtx,
182 unsigned timeIdx)
const 184 DiffusionModule::addDiffusiveFlux(flux, elemCtx, scvfIdx, timeIdx);
185 EnergyModule::addDiffusiveFlux(flux, elemCtx, scvfIdx, timeIdx);
195 const ElementContext& elemCtx,
197 unsigned timeIdx)
const 199 Opm::Valgrind::SetUndefined(source);
200 elemCtx.problem().source(source, elemCtx, dofIdx, timeIdx);
201 Opm::Valgrind::CheckDefined(source);
204 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
205 source[ncp0EqIdx + phaseIdx] =
206 phaseNcp(elemCtx, dofIdx, timeIdx, phaseIdx);
213 template <
class LhsEval = Evaluation>
217 unsigned phaseIdx)
const 219 const auto& fluidState = elemCtx.intensiveQuantities(dofIdx, timeIdx).fluidState();
220 typedef typename std::remove_const<typename std::remove_reference<decltype(fluidState)>::type>::type FluidState;
222 typedef Opm::MathToolbox<LhsEval> LhsToolbox;
224 const LhsEval& a = phaseNotPresentIneq_<FluidState, LhsEval>(fluidState, phaseIdx);
225 const LhsEval& b = phasePresentIneq_<FluidState, LhsEval>(fluidState, phaseIdx);
226 return LhsToolbox::min(a, b);
234 template <
class Flu
idState,
class LhsEval>
235 LhsEval phasePresentIneq_(
const FluidState& fluidState,
unsigned phaseIdx)
const 237 typedef Opm::MathToolbox<typename FluidState::Scalar> FsToolbox;
239 return FsToolbox::template decay<LhsEval>(fluidState.saturation(phaseIdx));
246 template <
class Flu
idState,
class LhsEval>
247 LhsEval phaseNotPresentIneq_(
const FluidState& fluidState,
unsigned phaseIdx)
const 249 typedef Opm::MathToolbox<typename FluidState::Scalar> FsToolbox;
253 for (
unsigned i = 0; i < numComponents; ++i)
254 a -= FsToolbox::template decay<LhsEval>(fluidState.moleFraction(phaseIdx, i));
void computeStorage(Dune::FieldVector< LhsEval, numEq > &storage, const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx) const
Evaluate the amount all conservation quantities (e.g.
Definition: ncplocalresidual.hh:103
Details needed to calculate the local residual in the compositional multi-phase NCP-model ...
Definition: ncplocalresidual.hh:46
void addDiffusiveFlux(RateVector &flux, const ElementContext &elemCtx, unsigned scvfIdx, unsigned timeIdx) const
Adds the diffusive flux at a given flux integration point.
Definition: ncplocalresidual.hh:179
Definition: baseauxiliarymodule.hh:37
Provides the auxiliary methods required for consideration of the energy equation. ...
Definition: energymodule.hh:59
#define GET_PROP_VALUE(TypeTag, PropTagName)
Access the value attribute of a property for a type tag.
Definition: propertysystem.hh:469
Contains the classes required to consider energy as a conservation quantity in a multi-phase module...
#define GET_PROP_TYPE(TypeTag, PropTagName)
Access the type attribute of a property for a type tag.
Definition: propertysystem.hh:486
Declares the properties required for the NCP compositional multi-phase model.
Provides the auxiliary methods required for consideration of the diffusion equation.
Definition: diffusionmodule.hh:53
LhsEval phaseNcp(const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx, unsigned phaseIdx) const
Returns the value of the NCP-function for a phase.
Definition: ncplocalresidual.hh:214
void computeSource(RateVector &source, const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx) const
Calculate the source term of the equation.
Definition: ncplocalresidual.hh:194
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: ncplocalresidual.hh:118
void addAdvectiveFlux(RateVector &flux, const ElementContext &elemCtx, unsigned scvfIdx, unsigned timeIdx) const
Add the advective mass flux at a given flux integration point.
Definition: ncplocalresidual.hh:134
void addPhaseStorage(Dune::FieldVector< LhsEval, numEq > &storage, const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx, unsigned phaseIdx) const
Adds the amount all conservation quantities (e.g.
Definition: ncplocalresidual.hh:78
Classes required for molecular diffusion.