28 #ifndef EWOMS_RICHARDS_BOUNDARY_RATE_VECTOR_HH 29 #define EWOMS_RICHARDS_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;
54 enum { contiEqIdx = Indices::contiEqIdx };
55 enum { liquidPhaseIdx =
GET_PROP_VALUE(TypeTag, LiquidPhaseIndex) };
57 typedef Opm::MathToolbox<Evaluation> Toolbox;
81 template <
class Context,
class Flu
idState>
82 void setFreeFlow(
const Context& context,
unsigned bfIdx,
unsigned timeIdx,
const FluidState& fluidState)
84 typename FluidSystem::template ParameterCache<typename FluidState::Scalar> paramCache;
85 paramCache.updateAll(fluidState);
87 ExtensiveQuantities extQuants;
88 extQuants.updateBoundary(context, bfIdx, timeIdx, fluidState, paramCache);
89 const auto& insideIntQuants = context.intensiveQuantities(bfIdx, timeIdx);
94 (*this) = Evaluation(0.0);
96 unsigned phaseIdx = liquidPhaseIdx;
98 if (fluidState.pressure(phaseIdx) > insideIntQuants.fluidState().pressure(phaseIdx))
99 density = FluidSystem::density(fluidState, paramCache, phaseIdx);
101 density = insideIntQuants.fluidState().density(phaseIdx);
105 (*this)[contiEqIdx] += extQuants.volumeFlux(phaseIdx) * density;
108 for (
unsigned i = 0; i < numEq; ++i) {
109 Opm::Valgrind::CheckDefined((*
this)[i]);
111 Opm::Valgrind::CheckDefined(*
this);
118 template <
class Context,
class Flu
idState>
122 const FluidState& fluidState)
124 this->
setFreeFlow(context, bfIdx, timeIdx, fluidState);
127 for (
unsigned eqIdx = 0; eqIdx < numEq; ++eqIdx) {
128 Evaluation& val = this->operator[](eqIdx);
129 val = Toolbox::min(0.0, val);
136 template <
class Context,
class Flu
idState>
140 const FluidState& fluidState)
142 this->
setFreeFlow(context, bfIdx, timeIdx, fluidState);
145 for (
unsigned eqIdx = 0; eqIdx < numEq; ++eqIdx) {
146 Evaluation& val = this->operator[](eqIdx);
147 val = Toolbox::max(0.0, val);
155 { (*this) = Evaluation(0.0); }
Intensive quantities required by the Richards model.
Definition: baseauxiliarymodule.hh:37
RichardsBoundaryRateVector(const Evaluation &value)
Definition: richardsboundaryratevector.hh:67
#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 setOutFlow(const Context &context, unsigned bfIdx, unsigned timeIdx, const FluidState &fluidState)
Specify an outflow boundary.
Definition: richardsboundaryratevector.hh:137
void setInFlow(const Context &context, unsigned bfIdx, unsigned timeIdx, const FluidState &fluidState)
Specify an inflow boundary.
Definition: richardsboundaryratevector.hh:119
void setFreeFlow(const Context &context, unsigned bfIdx, unsigned timeIdx, const FluidState &fluidState)
Specify a free-flow boundary.
Definition: richardsboundaryratevector.hh:82
Implements a boundary vector for the fully implicit Richards model.
Definition: richardsboundaryratevector.hh:44
void setNoFlow()
Specify a no-flow boundary for all conserved quantities.
Definition: richardsboundaryratevector.hh:154