blackoilratevector.hh
Go to the documentation of this file.
1// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2// vi: set et ts=4 sw=4 sts=4:
3/*
4 This file is part of the Open Porous Media project (OPM).
5
6 OPM is free software: you can redistribute it and/or modify
7 it under the terms of the GNU General Public License as published by
8 the Free Software Foundation, either version 2 of the License, or
9 (at your option) any later version.
10
11 OPM is distributed in the hope that it will be useful,
12 but WITHOUT ANY WARRANTY; without even the implied warranty of
13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 GNU General Public License for more details.
15
16 You should have received a copy of the GNU General Public License
17 along with OPM. If not, see <http://www.gnu.org/licenses/>.
18
19 Consult the COPYING file in the top-level source directory of this
20 module for the precise wording of the license and the list of
21 copyright holders.
22*/
28#ifndef EWOMS_BLACK_OIL_RATE_VECTOR_HH
29#define EWOMS_BLACK_OIL_RATE_VECTOR_HH
30
31#include <dune/common/fvector.hh>
32
33#include <opm/material/common/MathToolbox.hpp>
34#include <opm/material/common/Valgrind.hpp>
35#include <opm/material/constraintsolvers/NcpFlash.hpp>
36
44
45#include <stdexcept>
46
47namespace Opm {
48
58template <class TypeTag>
60 : public Dune::FieldVector<GetPropType<TypeTag, Properties::Evaluation>,
61 getPropValue<TypeTag, Properties::NumEq>()>
62{
67
72
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>;
85
86public:
87 BlackOilRateVector() : ParentType()
88 { Valgrind::SetUndefined(*this); }
89
93 BlackOilRateVector(Scalar value) : ParentType(Toolbox::createConstant(value))
94 {}
95
102 void setMassRate(const ParentType& value, unsigned pvtRegionIdx = 0)
103 {
104 ParentType::operator=(value);
105
106 // convert to "surface volume" if requested
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);
111 }
112 if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)) {
113 (*this)[FluidSystem::canonicalToActiveCompIdx(FluidSystem::oilCompIdx)] /=
114 FluidSystem::referenceDensity(FluidSystem::oilPhaseIdx, pvtRegionIdx);
115 }
116 if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) {
117 (*this)[FluidSystem::canonicalToActiveCompIdx(FluidSystem::waterCompIdx)] /=
118 FluidSystem::referenceDensity(FluidSystem::waterPhaseIdx, pvtRegionIdx);
119 }
120 if constexpr (enableSolvent) {
121 const auto& solventPvt = SolventModule::solventPvt();
122 (*this)[Indices::contiSolventEqIdx] /=
123 solventPvt.referenceDensity(pvtRegionIdx);
124 }
125 }
126 }
127
134 void setMolarRate(const ParentType& value, unsigned pvtRegionIdx = 0)
135 {
136 // first, assign molar rates
137 ParentType::operator=(value);
138
139 // then, convert them to mass rates
140 for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
141 (*this)[conti0EqIdx + compIdx] *= FluidSystem::molarMass(compIdx, pvtRegionIdx);
142 }
143
144 const auto& solventPvt = SolventModule::solventPvt();
145 (*this)[Indices::contiSolventEqIdx] *= solventPvt.molarMass(pvtRegionIdx);
146
147 if constexpr (enablePolymer) {
148 if constexpr (enablePolymerMolarWeight) {
149 throw std::logic_error("Set molar rate with polymer weight tracking not implemented");
150 }
151
152 (*this)[Indices::contiPolymerEqIdx] *= PolymerModule::molarMass(pvtRegionIdx);
153 }
154
155 if constexpr (enableFoam) {
156 throw std::logic_error("setMolarRate() not implemented for foam");
157 }
158
159 if constexpr (enableBrine) {
160 throw std::logic_error("setMolarRate() not implemented for salt water");
161 }
162
163 if constexpr (enableBioeffects) {
164 throw std::logic_error("setMolarRate() not implemented for bioeffects (biofilm/MICP)");
165 }
166
167 // convert to "surface volume" if requested
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);
172 }
173 if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)) {
174 (*this)[FluidSystem::canonicalToActiveCompIdx(FluidSystem::oilCompIdx)] /=
175 FluidSystem::referenceDensity(FluidSystem::oilPhaseIdx, pvtRegionIdx);
176 }
177 if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) {
178 (*this)[FluidSystem::canonicalToActiveCompIdx(FluidSystem::waterCompIdx)] /=
179 FluidSystem::referenceDensity(FluidSystem::waterPhaseIdx, pvtRegionIdx);
180 }
181 if constexpr (enableSolvent) {
182 (*this)[Indices::contiSolventEqIdx] /=
183 solventPvt.referenceDensity(pvtRegionIdx);
184 }
185 }
186 }
187
191 template <class FluidState, class RhsEval>
192 void setVolumetricRate(const FluidState& fluidState,
193 unsigned phaseIdx,
194 const RhsEval& volume)
195 {
196 for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
197 (*this)[conti0EqIdx + compIdx] =
198 fluidState.density(phaseIdx) *
199 fluidState.massFraction(phaseIdx, compIdx) *
200 volume;
201 }
202 }
203
207 template <class RhsEval>
208 BlackOilRateVector& operator=(const RhsEval& value)
209 {
210 for (unsigned i = 0; i < this->size(); ++i) {
211 (*this)[i] = value;
212 }
213 return *this;
214 }
215};
216
217} // namespace Opm
218
219#endif
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