My Project
Loading...
Searching...
No Matches
GenericOilGasFluidSystem.hpp
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 Copyright 2023 SINTEF Digital, Mathematics and Cybernetics.
5
6 This file is part of the Open Porous Media project (OPM).
7
8 OPM is free software: you can redistribute it and/or modify
9 it under the terms of the GNU General Public License as published by
10 the Free Software Foundation, either version 2 of the License, or
11 (at your option) any later version.
12
13 OPM is distributed in the hope that it will be useful,
14 but WITHOUT ANY WARRANTY; without even the implied warranty of
15 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16 GNU General Public License for more details.
17
18 You should have received a copy of the GNU General Public License
19 along with OPM. If not, see <http://www.gnu.org/licenses/>.
20
21 Consult the COPYING file in the top-level source directory of this
22 module for the precise wording of the license and the list of
23 copyright holders.
24*/
25
26#ifndef OPM_GENERICOILGASFLUIDSYSTEM_HPP
27#define OPM_GENERICOILGASFLUIDSYSTEM_HPP
28
29#include <opm/common/OpmLog/OpmLog.hpp>
30
32#include <opm/material/fluidsystems/PTFlashParameterCache.hpp> // TODO: this is something else need to check
34
35#include <cassert>
36#include <string>
37#include <string_view>
38
39#include <fmt/format.h>
40
41namespace Opm {
42
43
52 template<class Scalar, int NumComp>
53 class GenericOilGasFluidSystem : public BaseFluidSystem<Scalar, GenericOilGasFluidSystem<Scalar, NumComp> > {
54 public:
55 // TODO: I do not think these should be constant in fluidsystem, will try to make it non-constant later
56 static const int numPhases=2;
57 static const int numComponents = NumComp;
58 static const int numMisciblePhases=2;
59 // \Note: not totally sure when we should distinguish numMiscibleComponents and numComponents.
60 // Possibly when with a dummy phase like water?
61 static const int numMiscibleComponents = NumComp;
62 // TODO: phase location should be more general
63 static constexpr int oilPhaseIdx = 0;
64 static constexpr int gasPhaseIdx = 1;
65
66 template <class ValueType>
70
72 std::string name;
73 Scalar molar_mass;
74 Scalar critic_temp;
75 Scalar critic_pres;
76 Scalar critic_vol;
77 Scalar acentric_factor;
78
79 ComponentParam(const std::string_view name_, const Scalar molar_mass_, const Scalar critic_temp_,
80 const Scalar critic_pres_, const Scalar critic_vol_, const Scalar acentric_factor_)
81 : name(name_),
82 molar_mass(molar_mass_),
83 critic_temp(critic_temp_),
84 critic_pres(critic_pres_),
85 critic_vol(critic_vol_),
86 acentric_factor(acentric_factor_)
87 {}
88 };
89
90 template<typename Param>
91 static void addComponent(const Param& param)
92 {
93 // Check if the current size is less than the maximum allowed components.
94 if (component_param_.size() < numComponents) {
95 component_param_.push_back(param);
96 } else {
97 // Adding another component would exceed the limit.
98 const std::string msg = fmt::format("The fluid system has reached its maximum capacity of {} components,"
99 "the component '{}' will not be added.", NumComp, param.name);
100 OpmLog::note(msg);
101 // Optionally, throw an exception?
102 }
103 }
104
105 static void init()
106 {
107 component_param_.reserve(numComponents);
108 }
109
115 static Scalar acentricFactor(unsigned compIdx)
116 {
117 assert(isConsistent());
118 assert(compIdx < numComponents);
119
120 return component_param_[compIdx].acentric_factor;
121 }
127 static Scalar criticalTemperature(unsigned compIdx)
128 {
129 assert(isConsistent());
130 assert(compIdx < numComponents);
131
132 return component_param_[compIdx].critic_temp;
133 }
139 static Scalar criticalPressure(unsigned compIdx)
140 {
141 assert(isConsistent());
142 assert(compIdx < numComponents);
143
144 return component_param_[compIdx].critic_pres;
145 }
151 static Scalar criticalVolume(unsigned compIdx)
152 {
153 assert(isConsistent());
154 assert(compIdx < numComponents);
155
156 return component_param_[compIdx].critic_vol;
157 }
158
160 static Scalar molarMass(unsigned compIdx)
161 {
162 assert(isConsistent());
163 assert(compIdx < numComponents);
164
165 return component_param_[compIdx].molar_mass;
166 }
167
172 static Scalar interactionCoefficient(unsigned /*comp1Idx*/, unsigned /*comp2Idx*/)
173 {
174 assert(isConsistent());
175 // TODO: some data structure is needed to support this
176 return 0.0;
177 }
178
180 static std::string_view phaseName(unsigned phaseIdx)
181 {
182 static const std::string_view name[] = {"o", // oleic phase
183 "g"}; // gas phase
184
185 assert(phaseIdx < numPhases);
186 return name[phaseIdx];
187 }
188
190 static std::string_view componentName(unsigned compIdx)
191 {
192 assert(isConsistent());
193 assert(compIdx < numComponents);
194
195 return component_param_[compIdx].name;
196 }
197
201 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
202 static LhsEval density(const FluidState& fluidState,
203 const ParameterCache<ParamCacheEval>& paramCache,
204 unsigned phaseIdx)
205 {
206 assert(isConsistent());
207 assert(phaseIdx < numPhases);
208
209 if (phaseIdx == oilPhaseIdx || phaseIdx == gasPhaseIdx) {
210 return decay<LhsEval>(fluidState.averageMolarMass(phaseIdx) / paramCache.molarVolume(phaseIdx));
211 }
212
213 return {};
214 }
215
217 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
218 static LhsEval viscosity(const FluidState& fluidState,
219 const ParameterCache<ParamCacheEval>& paramCache,
220 unsigned phaseIdx)
221 {
222 assert(isConsistent());
223 assert(phaseIdx < numPhases);
224 // Use LBC method to calculate viscosity
225 return decay<LhsEval>(ViscosityModel::LBC(fluidState, paramCache, phaseIdx));
226 }
227
229 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
230 static LhsEval fugacityCoefficient(const FluidState& fluidState,
231 const ParameterCache<ParamCacheEval>& paramCache,
232 unsigned phaseIdx,
233 unsigned compIdx)
234 {
235 assert(isConsistent());
236 assert(phaseIdx < numPhases);
237 assert(compIdx < numComponents);
238
239 return decay<LhsEval>(PengRobinsonMixture::computeFugacityCoefficient(fluidState, paramCache, phaseIdx, compIdx));
240 }
241
243 static bool isCompressible([[maybe_unused]] unsigned phaseIdx)
244 {
245 assert(phaseIdx < numPhases);
246
247 return true;
248 }
249
251 static bool isIdealMixture([[maybe_unused]] unsigned phaseIdx)
252 {
253 assert(phaseIdx < numPhases);
254
255 return false;
256 }
257
259 static bool isLiquid(unsigned phaseIdx)
260 {
261 assert(phaseIdx < numPhases);
262
263 return (phaseIdx == 0);
264 }
265
267 static bool isIdealGas(unsigned phaseIdx)
268 {
269 assert(phaseIdx < numPhases);
270
271 return (phaseIdx == 1);
272 }
273 private:
274 static bool isConsistent() {
275 return component_param_.size() == NumComp;
276 }
277
278 static std::vector<ComponentParam> component_param_;
279 };
280
281 template <class Scalar, int NumComp>
282 std::vector<typename GenericOilGasFluidSystem<Scalar, NumComp>::ComponentParam>
283 GenericOilGasFluidSystem<Scalar, NumComp>::component_param_;
284}
285#endif // OPM_GENERICOILGASFLUIDSYSTEM_HPP
The base class for all fluid systems.
Specifies the parameter cache used by the SPE-5 fluid system.
The base class for all fluid systems.
Definition BaseFluidSystem.hpp:47
Scalar Scalar
The type used for scalar quantities.
Definition BaseFluidSystem.hpp:52
A two phase system that can contain NumComp components.
Definition GenericOilGasFluidSystem.hpp:53
static bool isIdealMixture(unsigned phaseIdx)
Returns true if and only if a fluid phase is assumed to be an ideal mixture.
Definition GenericOilGasFluidSystem.hpp:251
static LhsEval fugacityCoefficient(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &paramCache, unsigned phaseIdx, unsigned compIdx)
Calculate the fugacity coefficient [Pa] of an individual component in a fluid phase.
Definition GenericOilGasFluidSystem.hpp:230
static Scalar interactionCoefficient(unsigned, unsigned)
Returns the interaction coefficient for two components.
Definition GenericOilGasFluidSystem.hpp:172
static Scalar criticalTemperature(unsigned compIdx)
Critical temperature of a component [K].
Definition GenericOilGasFluidSystem.hpp:127
static bool isCompressible(unsigned phaseIdx)
Returns true if and only if a fluid phase is assumed to be compressible.
Definition GenericOilGasFluidSystem.hpp:243
static Scalar criticalVolume(unsigned compIdx)
Critical volume of a component [m3].
Definition GenericOilGasFluidSystem.hpp:151
static Scalar criticalPressure(unsigned compIdx)
Critical pressure of a component [Pa].
Definition GenericOilGasFluidSystem.hpp:139
static bool isIdealGas(unsigned phaseIdx)
Returns true if and only if a fluid phase is assumed to be an ideal gas.
Definition GenericOilGasFluidSystem.hpp:267
static LhsEval density(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &paramCache, unsigned phaseIdx)
Calculate the density [kg/m^3] of a fluid phase.
Definition GenericOilGasFluidSystem.hpp:202
static bool isLiquid(unsigned phaseIdx)
Return whether a phase is liquid.
Definition GenericOilGasFluidSystem.hpp:259
static std::string_view phaseName(unsigned phaseIdx)
Return the human readable name of a fluid phase.
Definition GenericOilGasFluidSystem.hpp:180
static Scalar acentricFactor(unsigned compIdx)
The acentric factor of a component [].
Definition GenericOilGasFluidSystem.hpp:115
static Scalar molarMass(unsigned compIdx)
Return the molar mass of a component in [kg/mol].
Definition GenericOilGasFluidSystem.hpp:160
static LhsEval viscosity(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &paramCache, unsigned phaseIdx)
Calculate the dynamic viscosity of a fluid phase [Pa*s].
Definition GenericOilGasFluidSystem.hpp:218
static std::string_view componentName(unsigned compIdx)
Return the human readable name of a component.
Definition GenericOilGasFluidSystem.hpp:190
Specifies the parameter cache used by the SPE-5 fluid system.
Definition PTFlashParameterCache.hpp:48
Scalar molarVolume(unsigned phaseIdx) const
Returns the molar volume of a phase [m^3/mol].
Definition PTFlashParameterCache.hpp:199
Implements the Peng-Robinson equation of state for a mixture.
Definition PengRobinsonMixture.hpp:41
static LhsEval computeFugacityCoefficient(const FluidState &fs, const Params &params, unsigned phaseIdx, unsigned compIdx)
Returns the fugacity coefficient of an individual component in the phase.
Definition PengRobinsonMixture.hpp:74
Definition LBC.hpp:40
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30
Definition GenericOilGasFluidSystem.hpp:71