28 #ifndef EWOMS_IMMISCIBLE_BOUNDARY_RATE_VECTOR_HH
29 #define EWOMS_IMMISCIBLE_BOUNDARY_RATE_VECTOR_HH
31 #include <opm/common/Valgrind.hpp>
32 #include <opm/material/constraintsolvers/NcpFlash.hpp>
44 template <
class TypeTag>
47 typedef typename GET_PROP_TYPE(TypeTag, RateVector) ParentType;
48 typedef typename GET_PROP_TYPE(TypeTag, ExtensiveQuantities) ExtensiveQuantities;
49 typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
50 typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
51 typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
52 typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
56 enum { conti0EqIdx = Indices::conti0EqIdx };
59 typedef Opm::MathToolbox<Evaluation> Toolbox;
97 template <
class Context,
class Flu
idState>
98 void setFreeFlow(
const Context& context,
unsigned bfIdx,
unsigned timeIdx,
const FluidState& fluidState)
100 typename FluidSystem::template ParameterCache<typename FluidState::Scalar> paramCache;
101 paramCache.updateAll(fluidState);
103 ExtensiveQuantities extQuants;
104 extQuants.updateBoundary(context, bfIdx, timeIdx, fluidState, paramCache);
105 const auto& insideIntQuants = context.intensiveQuantities(bfIdx, timeIdx);
110 (*this) = Evaluation(0.0);
111 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
113 if (fluidState.pressure(phaseIdx) > insideIntQuants.fluidState().pressure(phaseIdx))
114 density = FluidSystem::density(fluidState, paramCache, phaseIdx);
116 density = insideIntQuants.fluidState().density(phaseIdx);
118 Opm::Valgrind::CheckDefined(density);
119 Opm::Valgrind::CheckDefined(extQuants.volumeFlux(phaseIdx));
123 (*this)[conti0EqIdx + phaseIdx] += extQuants.volumeFlux(phaseIdx)*density;
126 Evaluation specificEnthalpy;
127 Scalar pBoundary = fluidState.pressure(phaseIdx);
128 const Evaluation& pElement = insideIntQuants.fluidState().pressure(phaseIdx);
129 if (pBoundary > pElement)
130 specificEnthalpy = FluidSystem::enthalpy(fluidState, paramCache, phaseIdx);
132 specificEnthalpy = insideIntQuants.fluidState().enthalpy(phaseIdx);
135 Evaluation enthalpyRate = density*extQuants.volumeFlux(phaseIdx)*specificEnthalpy;
136 EnergyModule::addToEnthalpyRate(*
this, enthalpyRate);
140 EnergyModule::addToEnthalpyRate(*
this, EnergyModule::heatConductionRate(extQuants));
143 for (
unsigned i = 0; i < numEq; ++i)
144 Opm::Valgrind::CheckDefined((*
this)[i]);
145 Opm::Valgrind::CheckDefined(*
this);
158 template <
class Context,
class Flu
idState>
162 const FluidState& fluidState)
164 this->
setFreeFlow(context, bfIdx, timeIdx, fluidState);
167 for (
unsigned eqIdx = 0; eqIdx < numEq; ++eqIdx) {
168 Evaluation& val = this->operator[](eqIdx);
169 val = Toolbox::min(0.0, val);
182 template <
class Context,
class Flu
idState>
186 const FluidState& fluidState)
188 this->
setFreeFlow(context, bfIdx, timeIdx, fluidState);
191 for (
unsigned eqIdx = 0; eqIdx < numEq; ++eqIdx) {
192 Evaluation& val = this->operator[](eqIdx);
193 val = Toolbox::max(0.0, val);
201 { (*this) = Evaluation(0.0); }
Contains the quantities which are are constant within a finite volume for the immiscible multi-phase ...
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
#define GET_PROP_TYPE(TypeTag, PropTagName)
Access the type attribute of a property for a type tag.
Definition: propertysystem.hh:486
void setNoFlow()
Specify a no-flow boundary for all conserved quantities.
Definition: immiscibleboundaryratevector.hh:200
Implements a boundary vector for the fully implicit multi-phase model which assumes immiscibility...
Definition: immiscibleboundaryratevector.hh:45
void setOutFlow(const Context &context, unsigned bfIdx, unsigned timeIdx, const FluidState &fluidState)
Specify an outflow boundary.
Definition: immiscibleboundaryratevector.hh:183
void setInFlow(const Context &context, unsigned bfIdx, unsigned timeIdx, const FluidState &fluidState)
Specify an inflow boundary.
Definition: immiscibleboundaryratevector.hh:159
ImmiscibleBoundaryRateVector(const Evaluation &value)
Constructor that assigns all entries to a scalar value.
Definition: immiscibleboundaryratevector.hh:73
void setFreeFlow(const Context &context, unsigned bfIdx, unsigned timeIdx, const FluidState &fluidState)
Specify a free-flow boundary.
Definition: immiscibleboundaryratevector.hh:98