28#ifndef EWOMS_PVS_BOUNDARY_RATE_VECTOR_HH
29#define EWOMS_PVS_BOUNDARY_RATE_VECTOR_HH
34#include <opm/material/common/Valgrind.hpp>
45template <
class TypeTag>
55 enum { numEq = getPropValue<TypeTag, Properties::NumEq>() };
56 enum { numPhases = getPropValue<TypeTag, Properties::NumPhases>() };
57 enum { numComponents = getPropValue<TypeTag, Properties::NumComponents>() };
58 enum { conti0EqIdx = Indices::conti0EqIdx };
59 enum { enableEnergy = getPropValue<TypeTag, Properties::EnableEnergy>() };
62 using Toolbox = Opm::MathToolbox<Evaluation>;
87 template <
class Context,
class Flu
idState>
88 void setFreeFlow(
const Context& context,
unsigned bfIdx,
unsigned timeIdx,
const FluidState& fluidState)
90 ExtensiveQuantities extQuants;
91 extQuants.updateBoundary(context, bfIdx, timeIdx, fluidState);
92 const auto& insideIntQuants = context.intensiveQuantities(bfIdx, timeIdx);
93 unsigned focusDofIdx = context.focusDofIndex();
94 unsigned interiorDofIdx = context.interiorScvIndex(bfIdx, timeIdx);
99 (*this) = Evaluation(0.0);
100 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
102 if (fluidState.pressure(phaseIdx) > insideIntQuants.fluidState().pressure(phaseIdx)) {
103 if (focusDofIdx == interiorDofIdx)
104 density = fluidState.density(phaseIdx);
106 density = Opm::getValue(fluidState.density(phaseIdx));
108 else if (focusDofIdx == interiorDofIdx)
109 density = insideIntQuants.fluidState().density(phaseIdx);
111 density = Opm::getValue(insideIntQuants.fluidState().density(phaseIdx));
113 for (
unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
115 if (fluidState.pressure(phaseIdx) > insideIntQuants.fluidState().pressure(phaseIdx)) {
116 if (focusDofIdx == interiorDofIdx)
117 molarity = fluidState.molarity(phaseIdx, compIdx);
119 molarity = Opm::getValue(fluidState.molarity(phaseIdx, compIdx));
121 else if (focusDofIdx == interiorDofIdx)
122 molarity = insideIntQuants.fluidState().molarity(phaseIdx, compIdx);
124 molarity = Opm::getValue(insideIntQuants.fluidState().molarity(phaseIdx, compIdx));
128 (*this)[conti0EqIdx + compIdx] += extQuants.volumeFlux(phaseIdx)*molarity;
132 Evaluation specificEnthalpy;
133 if (fluidState.pressure(phaseIdx) > insideIntQuants.fluidState().pressure(phaseIdx)) {
134 if (focusDofIdx == interiorDofIdx)
135 specificEnthalpy = fluidState.enthalpy(phaseIdx);
137 specificEnthalpy = Opm::getValue(fluidState.enthalpy(phaseIdx));
139 else if (focusDofIdx == interiorDofIdx)
140 specificEnthalpy = insideIntQuants.fluidState().enthalpy(phaseIdx);
142 specificEnthalpy = Opm::getValue(insideIntQuants.fluidState().enthalpy(phaseIdx));
144 Evaluation enthalpyRate = density*extQuants.volumeFlux(phaseIdx)*specificEnthalpy;
145 EnergyModule::addToEnthalpyRate(*
this, enthalpyRate);
151 EnergyModule::addToEnthalpyRate(*
this, EnergyModule::thermalConductionRate(extQuants));
154 for (
unsigned i = 0; i < numEq; ++i)
155 Opm::Valgrind::CheckDefined((*
this)[i]);
162 template <
class Context,
class Flu
idState>
166 const FluidState& fluidState)
168 this->
setFreeFlow(context, bfIdx, timeIdx, fluidState);
171 for (
unsigned eqIdx = 0; eqIdx < numEq; ++eqIdx) {
172 Evaluation& val = this->operator[](eqIdx);
173 val = Toolbox::min(0.0, val);
180 template <
class Context,
class Flu
idState>
184 const FluidState& fluidState)
186 this->
setFreeFlow(context, bfIdx, timeIdx, fluidState);
189 for (
unsigned eqIdx = 0; eqIdx < numEq; ++eqIdx) {
190 Evaluation& val = this->operator[](eqIdx);
191 val = Toolbox::max(0.0, val);
199 { (*this) = Evaluation(0.0); }
Provides the auxiliary methods required for consideration of the energy equation.
Definition: energymodule.hh:50
Implements a rate vector on the boundary for the fully implicit compositional multi-phase primary var...
Definition: pvsboundaryratevector.hh:47
PvsBoundaryRateVector & operator=(const PvsBoundaryRateVector &value)=default
PvsBoundaryRateVector(const PvsBoundaryRateVector &value)=default
void setOutFlow(const Context &context, unsigned bfIdx, unsigned timeIdx, const FluidState &fluidState)
Specify an outflow boundary.
Definition: pvsboundaryratevector.hh:181
void setInFlow(const Context &context, unsigned bfIdx, unsigned timeIdx, const FluidState &fluidState)
Specify an inflow boundary.
Definition: pvsboundaryratevector.hh:163
PvsBoundaryRateVector()
Definition: pvsboundaryratevector.hh:65
PvsBoundaryRateVector(const Evaluation &value)
Definition: pvsboundaryratevector.hh:73
void setNoFlow()
Specify a no-flow boundary for all conserved quantities.
Definition: pvsboundaryratevector.hh:198
void setFreeFlow(const Context &context, unsigned bfIdx, unsigned timeIdx, const FluidState &fluidState)
Specify a free-flow boundary.
Definition: pvsboundaryratevector.hh:88
Contains the classes required to consider energy as a conservation quantity in a multi-phase module.
Definition: blackoilboundaryratevector.hh:37
typename Properties::Detail::GetPropImpl< TypeTag, Property >::type::type GetPropType
get the type alias defined in the property (equivalent to old macro GET_PROP_TYPE(....
Definition: propertysystem.hh:235
Declares the properties required for the compositional multi-phase primary variable switching model.