28 #ifndef EWOMS_PVS_BOUNDARY_RATE_VECTOR_HH 29 #define EWOMS_PVS_BOUNDARY_RATE_VECTOR_HH 34 #include <opm/common/Valgrind.hpp> 45 template <
class TypeTag>
48 typedef typename GET_PROP_TYPE(TypeTag, RateVector) ParentType;
49 typedef typename GET_PROP_TYPE(TypeTag, ExtensiveQuantities) ExtensiveQuantities;
50 typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
51 typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
52 typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
53 typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
58 enum { conti0EqIdx = Indices::conti0EqIdx };
62 typedef Opm::MathToolbox<Evaluation> Toolbox;
87 template <
class Context,
class Flu
idState>
88 void setFreeFlow(
const Context& context,
unsigned bfIdx,
unsigned timeIdx,
const FluidState& fluidState)
90 typename FluidSystem::template ParameterCache<typename FluidState::Scalar> paramCache;
91 paramCache.updateAll(fluidState);
93 ExtensiveQuantities extQuants;
94 extQuants.updateBoundary(context, bfIdx, timeIdx, fluidState, paramCache);
95 const auto& insideIntQuants = context.intensiveQuantities(bfIdx, timeIdx);
100 (*this) = Evaluation(0.0);
101 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
102 Evaluation meanMBoundary = 0;
103 for (
unsigned compIdx = 0; compIdx < numComponents; ++compIdx)
105 fluidState.moleFraction(phaseIdx, compIdx)*FluidSystem::molarMass(compIdx);
108 if (fluidState.pressure(phaseIdx) > insideIntQuants.fluidState().pressure(phaseIdx))
109 density = FluidSystem::density(fluidState, paramCache, phaseIdx);
111 density = insideIntQuants.fluidState().density(phaseIdx);
113 for (
unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
115 if (fluidState.pressure(phaseIdx) > insideIntQuants.fluidState().pressure(phaseIdx))
116 molarity = fluidState.moleFraction(phaseIdx, compIdx)*density/meanMBoundary;
118 molarity = insideIntQuants.fluidState().molarity(phaseIdx, compIdx);
122 (*this)[conti0EqIdx + compIdx] += extQuants.volumeFlux(phaseIdx)*molarity;
126 Evaluation specificEnthalpy;
127 if (fluidState.pressure(phaseIdx) > insideIntQuants.fluidState().pressure(phaseIdx))
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 (
unsigned i = 0; i < numEq; ++i)
142 Opm::Valgrind::CheckDefined((*
this)[i]);
149 template <
class Context,
class Flu
idState>
153 const FluidState& fluidState)
155 this->
setFreeFlow(context, bfIdx, timeIdx, fluidState);
158 for (
unsigned eqIdx = 0; eqIdx < numEq; ++eqIdx) {
159 Evaluation& val = this->operator[](eqIdx);
160 val = Toolbox::min(0.0, val);
167 template <
class Context,
class Flu
idState>
171 const FluidState& fluidState)
173 this->
setFreeFlow(context, bfIdx, timeIdx, fluidState);
176 for (
unsigned eqIdx = 0; eqIdx < numEq; ++eqIdx) {
177 Evaluation& val = this->operator[](eqIdx);
178 val = Toolbox::max(0.0, val);
186 { (*this) = Evaluation(0.0); }
PvsBoundaryRateVector(const Evaluation &value)
Definition: pvsboundaryratevector.hh:73
Definition: baseauxiliarymodule.hh:37
void setNoFlow()
Specify a no-flow boundary for all conserved quantities.
Definition: pvsboundaryratevector.hh:185
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
void setOutFlow(const Context &context, unsigned bfIdx, unsigned timeIdx, const FluidState &fluidState)
Specify an outflow boundary.
Definition: pvsboundaryratevector.hh:168
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
void setFreeFlow(const Context &context, unsigned bfIdx, unsigned timeIdx, const FluidState &fluidState)
Specify a free-flow boundary.
Definition: pvsboundaryratevector.hh:88
void setInFlow(const Context &context, unsigned bfIdx, unsigned timeIdx, const FluidState &fluidState)
Specify an inflow boundary.
Definition: pvsboundaryratevector.hh:150
Declares the properties required for the compositional multi-phase primary variable switching model...
Implements a rate vector on the boundary for the fully implicit compositional multi-phase primary var...
Definition: pvsboundaryratevector.hh:46