28 #ifndef EWOMS_BLACK_OIL_BOUNDARY_RATE_VECTOR_HH 29 #define EWOMS_BLACK_OIL_BOUNDARY_RATE_VECTOR_HH 31 #include <opm/common/Valgrind.hpp> 32 #include <opm/material/constraintsolvers/NcpFlash.hpp> 43 template <
class TypeTag>
46 typedef typename GET_PROP_TYPE(TypeTag, RateVector) ParentType;
47 typedef typename GET_PROP_TYPE(TypeTag, ExtensiveQuantities) ExtensiveQuantities;
48 typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
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, Indices) Indices;
56 enum { conti0EqIdx = Indices::conti0EqIdx };
80 template <
class Context,
class Flu
idState>
84 const FluidState& fluidState)
86 typename FluidSystem::ParameterCache paramCache;
87 paramCache.updateAll(fluidState);
89 ExtensiveQuantities extQuants;
90 extQuants.updateBoundary(context, bfIdx, timeIdx, fluidState, paramCache);
91 const auto& insideIntQuants = context.intensiveQuantities(bfIdx, timeIdx);
97 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
98 Scalar meanMBoundary = 0;
99 for (
unsigned compIdx = 0; compIdx < numComponents; ++compIdx)
100 meanMBoundary += fluidState.moleFraction(phaseIdx, compIdx)
101 * FluidSystem::molarMass(compIdx);
104 if (fluidState.pressure(phaseIdx) > insideIntQuants.fluidState().pressure(phaseIdx))
105 density = FluidSystem::density(fluidState, paramCache, phaseIdx);
107 density = insideIntQuants.fluidState().density(phaseIdx);
109 for (
unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
111 if (fluidState.pressure(phaseIdx) > insideIntQuants.fluidState().pressure(phaseIdx))
113 fluidState.moleFraction(phaseIdx, compIdx)
117 molarity = insideIntQuants.fluidState().molarity(phaseIdx, compIdx);
121 (*this)[conti0EqIdx + compIdx] += extQuants.volumeFlux(phaseIdx) * molarity;
126 for (
unsigned i = 0; i < numEq; ++i) {
127 Opm::Valgrind::CheckDefined((*
this)[i]);
129 Opm::Valgrind::CheckDefined(*
this);
136 template <
class Context,
class Flu
idState>
140 const FluidState& fluidState)
142 this->
setFreeFlow(context, bfIdx, timeIdx, fluidState);
146 for (
unsigned eqIdx = 0; eqIdx < numEq; ++eqIdx) {
147 Scalar& val = this->operator[](eqIdx);
148 val = std::min<Scalar>(0.0, val);
155 template <
class Context,
class Flu
idState>
159 const FluidState& fluidState)
161 this->
setFreeFlow(context, bfIdx, timeIdx, fluidState);
165 for (
unsigned eqIdx = 0; eqIdx < numEq; ++eqIdx) {
166 Scalar& val = this->operator[](eqIdx);
167 val = std::max( Scalar(0), val);
175 { (*this) = Scalar(0); }
Definition: baseauxiliarymodule.hh:37
Contains the quantities which are are constant within a finite volume in the black-oil model...
#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
BlackOilBoundaryRateVector()
Default constructor.
Definition: blackoilboundaryratevector.hh:62
void setNoFlow()
Specify a no-flow boundary for all conserved quantities.
Definition: blackoilboundaryratevector.hh:174
Implements a boundary vector for the fully implicit black-oil model.
Definition: blackoilboundaryratevector.hh:44
void setFreeFlow(const Context &context, unsigned bfIdx, unsigned timeIdx, const FluidState &fluidState)
Specify a free-flow boundary.
Definition: blackoilboundaryratevector.hh:81
BlackOilBoundaryRateVector(Scalar value)
Definition: blackoilboundaryratevector.hh:68
void setOutFlow(const Context &context, unsigned bfIdx, unsigned timeIdx, const FluidState &fluidState)
Specify an outflow boundary.
Definition: blackoilboundaryratevector.hh:156
void setInFlow(const Context &context, unsigned bfIdx, unsigned timeIdx, const FluidState &fluidState)
Specify an inflow boundary.
Definition: blackoilboundaryratevector.hh:137