26 #ifndef EWOMS_IMMISCIBLE_BOUNDARY_RATE_VECTOR_HH
27 #define EWOMS_IMMISCIBLE_BOUNDARY_RATE_VECTOR_HH
29 #include <opm/material/common/Valgrind.hpp>
30 #include <opm/material/constraintsolvers/NcpFlash.hpp>
42 template <
class TypeTag>
45 typedef typename GET_PROP_TYPE(TypeTag, RateVector) ParentType;
46 typedef typename GET_PROP_TYPE(TypeTag, ExtensiveQuantities) ExtensiveQuantities;
47 typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
48 typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
49 typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
50 typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
54 enum { conti0EqIdx = Indices::conti0EqIdx };
57 typedef Opm::MathToolbox<Evaluation> Toolbox;
95 template <
class Context,
class Flu
idState>
96 void setFreeFlow(
const Context &context,
int bfIdx,
int timeIdx,
const FluidState &fluidState)
98 typename FluidSystem::ParameterCache paramCache;
99 paramCache.updateAll(fluidState);
101 ExtensiveQuantities extQuants;
102 extQuants.updateBoundary(context, bfIdx, timeIdx, fluidState, paramCache);
103 const auto &insideIntQuants = context.intensiveQuantities(bfIdx, timeIdx);
108 (*this) = Toolbox::createConstant(0.0);
109 for (
int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
111 if (fluidState.pressure(phaseIdx) > insideIntQuants.fluidState().pressure(phaseIdx))
112 density = FluidSystem::density(fluidState, paramCache, phaseIdx);
114 density = insideIntQuants.fluidState().density(phaseIdx);
116 Valgrind::CheckDefined(density);
117 Valgrind::CheckDefined(extQuants.volumeFlux(phaseIdx));
121 (*this)[conti0EqIdx + phaseIdx] += extQuants.volumeFlux(phaseIdx)*density;
124 Evaluation specificEnthalpy;
125 Scalar pBoundary = fluidState.pressure(phaseIdx);
126 const Evaluation& pElement = insideIntQuants.fluidState().pressure(phaseIdx);
127 if (pBoundary > pElement)
128 specificEnthalpy = FluidSystem::enthalpy(fluidState, paramCache, phaseIdx);
130 specificEnthalpy = insideIntQuants.fluidState().enthalpy(phaseIdx);
133 Evaluation enthalpyRate = density*extQuants.volumeFlux(phaseIdx)*specificEnthalpy;
134 EnergyModule::addToEnthalpyRate(*
this, enthalpyRate);
138 EnergyModule::addToEnthalpyRate(*
this, EnergyModule::heatConductionRate(extQuants));
141 for (
int i = 0; i < numEq; ++i)
142 Valgrind::CheckDefined((*
this)[i]);
143 Valgrind::CheckDefined(*
this);
156 template <
class Context,
class Flu
idState>
157 void setInFlow(
const Context &context,
int bfIdx,
int timeIdx,
158 const FluidState &fluidState)
160 this->
setFreeFlow(context, bfIdx, timeIdx, fluidState);
163 for (
int eqIdx = 0; eqIdx < numEq; ++eqIdx) {
164 Evaluation& val = this->operator[](eqIdx);
165 val = Toolbox::min(0.0, val);
178 template <
class Context,
class Flu
idState>
179 void setOutFlow(
const Context &context,
int bfIdx,
int timeIdx,
180 const FluidState &fluidState)
182 this->
setFreeFlow(context, bfIdx, timeIdx, fluidState);
185 for (
int eqIdx = 0; eqIdx < numEq; ++eqIdx) {
186 Evaluation &val = this->operator[](eqIdx);
187 val = Toolbox::max(0.0, val);
195 { (*this) = Toolbox::createConstant(0.0); }
void setNoFlow()
Specify a no-flow boundary for all conserved quantities.
Definition: immiscibleboundaryratevector.hh:194
#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
ImmiscibleBoundaryRateVector(const Evaluation &value)
Constructor that assigns all entries to a scalar value.
Definition: immiscibleboundaryratevector.hh:71
Provides the auxiliary methods required for consideration of the energy equation. ...
Definition: energymodule.hh:54
void setFreeFlow(const Context &context, int bfIdx, int timeIdx, const FluidState &fluidState)
Specify a free-flow boundary.
Definition: immiscibleboundaryratevector.hh:96
ImmiscibleBoundaryRateVector(const ImmiscibleBoundaryRateVector &value)
Copy constructor.
Definition: immiscibleboundaryratevector.hh:80
Definition: baseauxiliarymodule.hh:35
ImmiscibleBoundaryRateVector()
Definition: immiscibleboundaryratevector.hh:61
Implements a boundary vector for the fully implicit multi-phase model which assumes immiscibility...
Definition: immiscibleboundaryratevector.hh:43
void setOutFlow(const Context &context, int bfIdx, int timeIdx, const FluidState &fluidState)
Specify an outflow boundary.
Definition: immiscibleboundaryratevector.hh:179
void setInFlow(const Context &context, int bfIdx, int timeIdx, const FluidState &fluidState)
Specify an inflow boundary.
Definition: immiscibleboundaryratevector.hh:157
Contains the quantities which are are constant within a finite volume for the immiscible multi-phase ...