28#ifndef EWOMS_BLACK_OIL_RATE_VECTOR_HH
29#define EWOMS_BLACK_OIL_RATE_VECTOR_HH
31#include <dune/common/fvector.hh>
33#include <opm/material/common/MathToolbox.hpp>
34#include <opm/material/common/Valgrind.hpp>
35#include <opm/material/constraintsolvers/NcpFlash.hpp>
58template <
class TypeTag>
60 :
public Dune::FieldVector<GetPropType<TypeTag, Properties::Evaluation>,
61 getPropValue<TypeTag, Properties::NumEq>()>
73 enum { numEq = getPropValue<TypeTag, Properties::NumEq>() };
74 enum { numComponents = getPropValue<TypeTag, Properties::NumComponents>() };
75 enum { conti0EqIdx = Indices::conti0EqIdx };
76 enum { enableEnergy = getPropValue<TypeTag, Properties::EnableEnergy>() };
77 enum { enableSolvent = getPropValue<TypeTag, Properties::EnableSolvent>() };
78 enum { enablePolymer = getPropValue<TypeTag, Properties::EnablePolymer>() };
79 enum { enablePolymerMolarWeight = getPropValue<TypeTag, Properties::EnablePolymerMW>() };
80 enum { enableFoam = getPropValue<TypeTag, Properties::EnableFoam>() };
81 enum { enableBrine = getPropValue<TypeTag, Properties::EnableBrine>() };
82 enum { enableBioeffects = getPropValue<TypeTag, Properties::EnableBioeffects>() };
83 using Toolbox = MathToolbox<Evaluation>;
84 using ParentType = Dune::FieldVector<Evaluation, numEq>;
88 { Valgrind::SetUndefined(*
this); }
102 void setMassRate(
const ParentType& value,
unsigned pvtRegionIdx = 0)
104 ParentType::operator=(value);
107 if constexpr (getPropValue<TypeTag, Properties::BlackoilConserveSurfaceVolume>()) {
108 if (FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) {
109 (*this)[FluidSystem::canonicalToActiveCompIdx(FluidSystem::gasCompIdx)] /=
110 FluidSystem::referenceDensity(FluidSystem::gasPhaseIdx, pvtRegionIdx);
112 if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)) {
113 (*this)[FluidSystem::canonicalToActiveCompIdx(FluidSystem::oilCompIdx)] /=
114 FluidSystem::referenceDensity(FluidSystem::oilPhaseIdx, pvtRegionIdx);
116 if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) {
117 (*this)[FluidSystem::canonicalToActiveCompIdx(FluidSystem::waterCompIdx)] /=
118 FluidSystem::referenceDensity(FluidSystem::waterPhaseIdx, pvtRegionIdx);
120 if constexpr (enableSolvent) {
122 (*this)[Indices::contiSolventEqIdx] /=
123 solventPvt.referenceDensity(pvtRegionIdx);
137 ParentType::operator=(value);
140 for (
unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
141 (*this)[conti0EqIdx + compIdx] *= FluidSystem::molarMass(compIdx, pvtRegionIdx);
145 (*this)[Indices::contiSolventEqIdx] *= solventPvt.molarMass(pvtRegionIdx);
147 if constexpr (enablePolymer) {
148 if constexpr (enablePolymerMolarWeight) {
149 throw std::logic_error(
"Set molar rate with polymer weight tracking not implemented");
155 if constexpr (enableFoam) {
156 throw std::logic_error(
"setMolarRate() not implemented for foam");
159 if constexpr (enableBrine) {
160 throw std::logic_error(
"setMolarRate() not implemented for salt water");
163 if constexpr (enableBioeffects) {
164 throw std::logic_error(
"setMolarRate() not implemented for bioeffects (biofilm/MICP)");
168 if constexpr (getPropValue<TypeTag, Properties::BlackoilConserveSurfaceVolume>()) {
169 if (FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) {
170 (*this)[FluidSystem::canonicalToActiveCompIdx(FluidSystem::gasCompIdx)] /=
171 FluidSystem::referenceDensity(FluidSystem::gasPhaseIdx, pvtRegionIdx);
173 if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)) {
174 (*this)[FluidSystem::canonicalToActiveCompIdx(FluidSystem::oilCompIdx)] /=
175 FluidSystem::referenceDensity(FluidSystem::oilPhaseIdx, pvtRegionIdx);
177 if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) {
178 (*this)[FluidSystem::canonicalToActiveCompIdx(FluidSystem::waterCompIdx)] /=
179 FluidSystem::referenceDensity(FluidSystem::waterPhaseIdx, pvtRegionIdx);
181 if constexpr (enableSolvent) {
182 (*this)[Indices::contiSolventEqIdx] /=
183 solventPvt.referenceDensity(pvtRegionIdx);
191 template <
class Flu
idState,
class RhsEval>
194 const RhsEval& volume)
196 for (
unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
197 (*this)[conti0EqIdx + compIdx] =
198 fluidState.density(phaseIdx) *
199 fluidState.massFraction(phaseIdx, compIdx) *
207 template <
class RhsEval>
210 for (
unsigned i = 0; i < this->size(); ++i) {
Contains the classes required to extend the black-oil model by bioeffects.
Contains the classes required to extend the black-oil model by brine.
Contains the classes required to extend the black-oil model to include the effects of foam.
Contains the classes required to extend the black-oil model by polymer.
Contains the classes required to extend the black-oil model by solvents.
Contains the high level supplements required to extend the black oil model by brine.
Definition: blackoilbrinemodules.hh:56
Contains the high level supplements required to extend the black oil model to include the effects of ...
Definition: blackoilfoammodules.hh:58
Contains the high level supplements required to extend the black oil model by polymer.
Definition: blackoilpolymermodules.hh:64
Scalar molarMass() const
Definition: blackoilpolymermodules.hh:553
Implements a vector representing mass, molar or volumetric rates for the black oil model.
Definition: blackoilratevector.hh:62
BlackOilRateVector(Scalar value)
Definition: blackoilratevector.hh:93
void setMolarRate(const ParentType &value, unsigned pvtRegionIdx=0)
Set a molar rate of the conservation quantities.
Definition: blackoilratevector.hh:134
BlackOilRateVector & operator=(const RhsEval &value)
Assignment operator from a scalar or a function evaluation.
Definition: blackoilratevector.hh:208
void setVolumetricRate(const FluidState &fluidState, unsigned phaseIdx, const RhsEval &volume)
Set a volumetric rate of a phase.
Definition: blackoilratevector.hh:192
void setMassRate(const ParentType &value, unsigned pvtRegionIdx=0)
Set a mass rate of the conservation quantities.
Definition: blackoilratevector.hh:102
BlackOilRateVector()
Definition: blackoilratevector.hh:87
Contains the high level supplements required to extend the black oil model by solvents.
Definition: blackoilsolventmodules.hh:68
static const SolventPvt & solventPvt()
Definition: blackoilsolventmodules.hh:363
Declare the properties used by the infrastructure code of the finite volume discretizations.
Defines the common properties required by the porous medium multi-phase models.
Definition: blackoilbioeffectsmodules.hh:43
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:233