28#ifndef EWOMS_BLACK_OIL_BOUNDARY_RATE_VECTOR_HH
29#define EWOMS_BLACK_OIL_BOUNDARY_RATE_VECTOR_HH
31#include <opm/material/common/Valgrind.hpp>
32#include <opm/material/constraintsolvers/NcpFlash.hpp>
44template <
class TypeTag>
62 enum { conti0EqIdx = Indices::conti0EqIdx };
63 enum { contiEnergyEqIdx = Indices::contiEnergyEqIdx };
93 template <
class Context,
class Flu
idState>
97 const FluidState& fluidState)
110 if (!FluidSystem::phaseIsActive(
phaseIdx)) {
127 using RhsEval =
typename std::conditional<std::is_same<typename FluidState::Scalar, Evaluation>::value,
128 Evaluation, Scalar>::type;
137 for (
unsigned i = 0; i < tmp.size(); ++i)
138 (*
this)[i] += tmp[i];
141 if constexpr (enableEnergy) {
146 density = fluidState.density(
phaseIdx);
168 if constexpr (enableSolvent) {
169 (*this)[Indices::contiSolventEqIdx] =
extQuants.solventVolumeFlux();
170 if (blackoilConserveSurfaceVolume)
171 (*this)[Indices::contiSolventEqIdx] *=
insideIntQuants.solventInverseFormationVolumeFactor();
173 (*
this)[Indices::contiSolventEqIdx] *=
insideIntQuants.solventDensity();
177 if constexpr (enablePolymer) {
178 (*this)[Indices::contiPolymerEqIdx] =
extQuants.volumeFlux(FluidSystem::waterPhaseIdx) *
insideIntQuants.polymerConcentration();
181 if constexpr (enableMICP) {
182 (*this)[Indices::contiMicrobialEqIdx] =
extQuants.volumeFlux(FluidSystem::waterPhaseIdx) *
insideIntQuants.microbialConcentration();
183 (*this)[Indices::contiOxygenEqIdx] =
extQuants.volumeFlux(FluidSystem::waterPhaseIdx) *
insideIntQuants.oxygenConcentration();
184 (*this)[Indices::contiUreaEqIdx] =
extQuants.volumeFlux(FluidSystem::waterPhaseIdx) *
insideIntQuants.ureaConcentration();
188 LocalResidual::adaptMassConservationQuantities_(*
this,
insideIntQuants.pvtRegionIndex());
191 if constexpr (enableEnergy)
195 for (
unsigned i = 0; i < numEq; ++i) {
196 Valgrind::CheckDefined((*
this)[i]);
198 Valgrind::CheckDefined(*
this);
205 template <
class Context,
class Flu
idState>
209 const FluidState& fluidState)
216 Scalar&
val = this->operator[](
eqIdx);
217 val = std::min<Scalar>(0.0,
val);
224 template <
class Context,
class Flu
idState>
228 const FluidState& fluidState)
235 Scalar&
val = this->operator[](
eqIdx);
236 val = std::max( Scalar(0),
val);
244 { (*this) = Scalar(0); }
253 template <
class Context,
class Flu
idState>
263 if constexpr (enableEnergy) {
267 (*this)[contiEnergyEqIdx] +=
extQuants.energyFlux();
270 for (
unsigned i = 0; i < numEq; ++i)
271 Valgrind::CheckDefined((*
this)[i]);
272 Valgrind::CheckDefined(*
this);
Contains the classes required to extend the black-oil model by energy.
Contains the quantities which are are constant within a finite volume in the black-oil model.
Implements a boundary vector for the fully implicit black-oil model.
Definition blackoilboundaryratevector.hh:46
BlackOilBoundaryRateVector()
Default constructor.
Definition blackoilboundaryratevector.hh:75
BlackOilBoundaryRateVector(Scalar value)
Definition blackoilboundaryratevector.hh:81
void setOutFlow(const Context &context, unsigned bfIdx, unsigned timeIdx, const FluidState &fluidState)
Specify an outflow boundary.
Definition blackoilboundaryratevector.hh:225
BlackOilBoundaryRateVector(const BlackOilBoundaryRateVector &value)=default
Copy constructor.
void setFreeFlow(const Context &context, unsigned bfIdx, unsigned timeIdx, const FluidState &fluidState)
Specify a free-flow boundary.
Definition blackoilboundaryratevector.hh:94
void setNoFlow()
Specify a no-flow boundary for all conserved quantities.
Definition blackoilboundaryratevector.hh:243
void setInFlow(const Context &context, unsigned bfIdx, unsigned timeIdx, const FluidState &fluidState)
Specify an inflow boundary.
Definition blackoilboundaryratevector.hh:206
void setThermalFlow(const Context &context, unsigned bfIdx, unsigned timeIdx, const FluidState &boundaryFluidState)
an energy flux that corresponds to the thermal conduction from
Definition blackoilboundaryratevector.hh:254
Contains the high level supplements required to extend the black oil model by energy.
Definition blackoilenergymodules.hh:52
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242
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