My Project
Loading...
Searching...
No Matches
Spe5FluidSystem.hpp
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*/
27#ifndef OPM_SPE5_FLUID_SYSTEM_HPP
28#define OPM_SPE5_FLUID_SYSTEM_HPP
29
30#include "BaseFluidSystem.hpp"
32
35
37
38#include <string_view>
39
40namespace Opm {
41
56template <class Scalar>
58 : public BaseFluidSystem<Scalar, Spe5FluidSystem<Scalar> >
59{
61
62 typedef typename ::Opm::PengRobinsonMixture<Scalar, ThisType> PengRobinsonMixture;
63 typedef typename ::Opm::PengRobinson<Scalar> PengRobinson;
64
65 static const Scalar R;
66
67public:
69 template <class Evaluation>
70 struct ParameterCache : public Spe5ParameterCache<Evaluation, ThisType>
71 {};
72
73 /****************************************
74 * Fluid phase parameters
75 ****************************************/
76
78 static const int numPhases = 3;
79
81 static const int gasPhaseIdx = 0;
83 static const int waterPhaseIdx = 1;
85 static const int oilPhaseIdx = 2;
86
88 typedef ::Opm::H2O<Scalar> H2O;
89
91 static std::string_view phaseName(unsigned phaseIdx)
92 {
93 static const std::string_view name[] = {
94 "gas",
95 "water",
96 "oil",
97 };
98
99 assert(phaseIdx < numPhases);
100 return name[phaseIdx];
101 }
102
104 static bool isLiquid(unsigned phaseIdx)
105 {
106 //assert(0 <= phaseIdx && phaseIdx < numPhases);
107 return phaseIdx != gasPhaseIdx;
108 }
109
115 static bool isCompressible(unsigned /*phaseIdx*/)
116 {
117 //assert(0 <= phaseIdx && phaseIdx < numPhases);
118 return true;
119 }
120
122 static bool isIdealGas(unsigned /*phaseIdx*/)
123 {
124 //assert(0 <= phaseIdx && phaseIdx < numPhases);
125 return false; // gas is not ideal here!
126 }
127
129 static bool isIdealMixture(unsigned phaseIdx)
130 {
131 // always use the reference oil for the fugacity coefficents,
132 // so they cannot be dependent on composition and they the
133 // phases thus always an ideal mixture
134 return phaseIdx == waterPhaseIdx;
135 }
136
137 /****************************************
138 * Component related parameters
139 ****************************************/
140
142 static const int numComponents = 7;
143
144 static const int H2OIdx = 0;
145 static const int C1Idx = 1;
146 static const int C3Idx = 2;
147 static const int C6Idx = 3;
148 static const int C10Idx = 4;
149 static const int C15Idx = 5;
150 static const int C20Idx = 6;
151
153 static std::string_view componentName(unsigned compIdx)
154 {
155 static const std::string_view name[] = {
156 H2O::name(),
157 "C1",
158 "C3",
159 "C6",
160 "C10",
161 "C15",
162 "C20"
163 };
164
165 assert(0 <= compIdx && compIdx < numComponents);
166 return name[compIdx];
167 }
168
170 static Scalar molarMass(unsigned compIdx)
171 {
172 return
173 (compIdx == H2OIdx)
175 : (compIdx == C1Idx)
176 ? 16.04e-3
177 : (compIdx == C3Idx)
178 ? 44.10e-3
179 : (compIdx == C6Idx)
180 ? 86.18e-3
181 : (compIdx == C10Idx)
182 ? 142.29e-3
183 : (compIdx == C15Idx)
184 ? 206.00e-3
185 : (compIdx == C20Idx)
186 ? 282.00e-3
187 : 1e30;
188 }
189
193 static Scalar criticalTemperature(unsigned compIdx)
194 {
195 return
196 (compIdx == H2OIdx)
198 : (compIdx == C1Idx)
199 ? 343.0*5/9
200 : (compIdx == C3Idx)
201 ? 665.7*5/9
202 : (compIdx == C6Idx)
203 ? 913.4*5/9
204 : (compIdx == C10Idx)
205 ? 1111.8*5/9
206 : (compIdx == C15Idx)
207 ? 1270.0*5/9
208 : (compIdx == C20Idx)
209 ? 1380.0*5/9
210 : 1e30;
211 }
212
216 static Scalar criticalPressure(unsigned compIdx)
217 {
218 return
219 (compIdx == H2OIdx)
221 : (compIdx == C1Idx)
222 ? 667.8 * 6894.7573
223 : (compIdx == C3Idx)
224 ? 616.3 * 6894.7573
225 : (compIdx == C6Idx)
226 ? 436.9 * 6894.7573
227 : (compIdx == C10Idx)
228 ? 304.0 * 6894.7573
229 : (compIdx == C15Idx)
230 ? 200.0 * 6894.7573
231 : (compIdx == C20Idx)
232 ? 162.0 * 6894.7573
233 : 1e30;
234 }
235
239 static Scalar criticalMolarVolume(unsigned compIdx)
240 {
241 return
242 (compIdx == H2OIdx)
244 : (compIdx == C1Idx)
246 : (compIdx == C3Idx)
248 : (compIdx == C6Idx)
250 : (compIdx == C10Idx)
252 : (compIdx == C15Idx)
254 : (compIdx == C20Idx)
256 : 1e30;
257 }
258
262 static Scalar acentricFactor(unsigned compIdx)
263 {
264 return
265 (compIdx == H2OIdx)
267 : (compIdx == C1Idx)
268 ? 0.0130
269 : (compIdx == C3Idx)
270 ? 0.1524
271 : (compIdx == C6Idx)
272 ? 0.3007
273 : (compIdx == C10Idx)
274 ? 0.4885
275 : (compIdx == C15Idx)
276 ? 0.6500
277 : (compIdx == C20Idx)
278 ? 0.8500
279 : 1e30;
280 }
281
287 static Scalar interactionCoefficient(unsigned comp1Idx, unsigned comp2Idx)
288 {
289 unsigned i = std::min(comp1Idx, comp2Idx);
290 unsigned j = std::max(comp1Idx, comp2Idx);
291 if (i == C1Idx && (j == C15Idx || j == C20Idx))
292 return 0.05;
293 else if (i == C3Idx && (j == C15Idx || j == C20Idx))
294 return 0.005;
295 return 0;
296 }
297
298 /****************************************
299 * Methods which compute stuff
300 ****************************************/
301
310 static void init(Scalar minT = 273.15,
311 Scalar maxT = 373.15,
312 Scalar minP = 1e4,
313 Scalar maxP = 100e6)
314 {
315 PengRobinsonParamsMixture<Scalar, ThisType, gasPhaseIdx, /*useSpe5=*/true> prParams;
316
317 // find envelopes of the 'a' and 'b' parameters for the range
318 // minT <= T <= maxT and minP <= p <= maxP. For
319 // this we take advantage of the fact that 'a' and 'b' for
320 // mixtures is just a convex combination of the attractive and
321 // repulsive parameters of the pure components
322
323 Scalar minA = 1e30, maxA = -1e30;
324 Scalar minB = 1e30, maxB = -1e30;
325
326 prParams.updatePure(minT, minP);
327 for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
328 minA = std::min(prParams.pureParams(compIdx).a(), minA);
329 maxA = std::max(prParams.pureParams(compIdx).a(), maxA);
330 minB = std::min(prParams.pureParams(compIdx).b(), minB);
331 maxB = std::max(prParams.pureParams(compIdx).b(), maxB);
332 };
333
334 prParams.updatePure(maxT, minP);
335 for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
336 minA = std::min(prParams.pureParams(compIdx).a(), minA);
337 maxA = std::max(prParams.pureParams(compIdx).a(), maxA);
338 minB = std::min(prParams.pureParams(compIdx).b(), minB);
339 maxB = std::max(prParams.pureParams(compIdx).b(), maxB);
340 };
341
342 prParams.updatePure(minT, maxP);
343 for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
344 minA = std::min(prParams.pureParams(compIdx).a(), minA);
345 maxA = std::max(prParams.pureParams(compIdx).a(), maxA);
346 minB = std::min(prParams.pureParams(compIdx).b(), minB);
347 maxB = std::max(prParams.pureParams(compIdx).b(), maxB);
348 };
349
350 prParams.updatePure(maxT, maxP);
351 for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
352 minA = std::min(prParams.pureParams(compIdx).a(), minA);
353 maxA = std::max(prParams.pureParams(compIdx).a(), maxA);
354 minB = std::min(prParams.pureParams(compIdx).b(), minB);
355 maxB = std::max(prParams.pureParams(compIdx).b(), maxB);
356 };
357
358 PengRobinson::init(/*aMin=*/minA, /*aMax=*/maxA, /*na=*/100,
359 /*bMin=*/minB, /*bMax=*/maxB, /*nb=*/200);
360 }
361
363 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
364 static LhsEval density(const FluidState& fluidState,
365 const ParameterCache<ParamCacheEval>& paramCache,
366 unsigned phaseIdx)
367 {
368 assert(phaseIdx < numPhases);
369
370 return fluidState.averageMolarMass(phaseIdx)/paramCache.molarVolume(phaseIdx);
371 }
372
374 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
375 static LhsEval viscosity(const FluidState& /*fluidState*/,
376 const ParameterCache<ParamCacheEval>& /*paramCache*/,
377 unsigned phaseIdx)
378 {
379 assert(phaseIdx <= numPhases);
380
381 if (phaseIdx == gasPhaseIdx) {
382 // given by SPE-5 in table on page 64. we use a constant
383 // viscosity, though...
384 return 0.0170e-2 * 0.1;
385 }
386 else if (phaseIdx == waterPhaseIdx)
387 // given by SPE-5: 0.7 centi-Poise = 0.0007 Pa s
388 return 0.7e-2 * 0.1;
389 else {
390 assert(phaseIdx == oilPhaseIdx);
391 // given by SPE-5 in table on page 64. we use a constant
392 // viscosity, though...
393 return 0.208e-2 * 0.1;
394 }
395 }
396
398 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
399 static LhsEval fugacityCoefficient(const FluidState& fluidState,
400 const ParameterCache<ParamCacheEval>& paramCache,
401 unsigned phaseIdx,
402 unsigned compIdx)
403 {
404 assert(phaseIdx <= numPhases);
405 assert(compIdx <= numComponents);
406
407 if (phaseIdx == oilPhaseIdx || phaseIdx == gasPhaseIdx)
409 paramCache,
410 phaseIdx,
411 compIdx);
412 else {
413 assert(phaseIdx == waterPhaseIdx);
414 return
415 henryCoeffWater_(compIdx, fluidState.temperature(waterPhaseIdx))
416 / fluidState.pressure(waterPhaseIdx);
417 }
418 }
419
420protected:
421 template <class LhsEval>
422 static LhsEval henryCoeffWater_(unsigned compIdx, const LhsEval& temperature)
423 {
424 // use henry's law for the solutes and the vapor pressure for
425 // the solvent.
426 switch (compIdx) {
427 case H2OIdx: return H2O::vaporPressure(temperature);
428
429 // the values of the Henry constant for the solutes have
430 // are faked so far...
431 case C1Idx: return 5.57601e+09;
432 case C3Idx: return 1e10;
433 case C6Idx: return 1e10;
434 case C10Idx: return 1e10;
435 case C15Idx: return 1e10;
436 case C20Idx: return 1e10;
437 default: throw std::logic_error("Unknown component index "+std::to_string(compIdx));
438 }
439 }
440};
441
442template <class Scalar>
443const Scalar Spe5FluidSystem<Scalar>::R = Constants<Scalar>::R;
444
445} // namespace Opm
446
447#endif
The base class for all fluid systems.
A central place for various physical constants occuring in some equations.
Implements the Peng-Robinson equation of state for a mixture.
Specifies the parameter cache used by the SPE-5 fluid system.
Class implementing cubic splines.
The base class for all fluid systems.
Definition BaseFluidSystem.hpp:47
Scalar Scalar
The type used for scalar quantities.
Definition BaseFluidSystem.hpp:52
static const Scalar R
The ideal gas constant [J/(mol K)].
Definition Constants.hpp:45
static const Scalar criticalTemperature()
Returns the critical temperature of water.
Definition H2O.hpp:95
static std::string_view name()
A human readable name for the water.
Definition H2O.hpp:77
static Evaluation vaporPressure(Evaluation temperature)
The vapor pressure in of pure water at a given temperature.
Definition H2O.hpp:141
static const Scalar criticalMolarVolume()
Returns the molar volume of water at the critical point.
Definition H2O.hpp:113
static const Scalar acentricFactor()
The acentric factor of water.
Definition H2O.hpp:89
static const Scalar criticalPressure()
Returns the critical pressure of water.
Definition H2O.hpp:101
static const Scalar molarMass()
The molar mass in of water.
Definition H2O.hpp:83
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
The mixing rule for the oil and the gas phases of the SPE5 problem.
Definition PengRobinsonParamsMixture.hpp:60
void updatePure(const FluidState &fluidState)
Update Peng-Robinson parameters for the pure components.
Definition PengRobinsonParamsMixture.hpp:82
const PureParams & pureParams(unsigned compIdx) const
Return the Peng-Robinson parameters of a pure substance,.
Definition PengRobinsonParamsMixture.hpp:205
Scalar a() const
Returns the attractive parameter 'a' of the Peng-Robinson fluid.
Definition PengRobinsonParams.hpp:50
Scalar b() const
Returns the repulsive parameter 'b' of the Peng-Robinson fluid.
Definition PengRobinsonParams.hpp:57
Implements the Peng-Robinson equation of state for liquids and gases.
Definition PengRobinson.hpp:59
The fluid system for the oil, gas and water phases of the SPE5 problem.
Definition Spe5FluidSystem.hpp:59
static const int oilPhaseIdx
Index of the oil phase.
Definition Spe5FluidSystem.hpp:85
static bool isIdealMixture(unsigned phaseIdx)
Returns true if and only if a fluid phase is assumed to be an ideal mixture.
Definition Spe5FluidSystem.hpp:129
static const int H2OIdx
Index of the water component.
Definition Spe5FluidSystem.hpp:144
static const int numComponents
Number of chemical species in the fluid system.
Definition Spe5FluidSystem.hpp:142
static const int C15Idx
Index of the C15 component.
Definition Spe5FluidSystem.hpp:149
static const int numPhases
Number of fluid phases in the fluid system.
Definition Spe5FluidSystem.hpp:78
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 Spe5FluidSystem.hpp:399
static Scalar criticalPressure(unsigned compIdx)
Critical pressure of a component [Pa].
Definition Spe5FluidSystem.hpp:216
static LhsEval density(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &paramCache, unsigned phaseIdx)
Calculate the density [kg/m^3] of a fluid phase.
Definition Spe5FluidSystem.hpp:364
static Scalar criticalTemperature(unsigned compIdx)
Critical temperature of a component [K].
Definition Spe5FluidSystem.hpp:193
static Scalar criticalMolarVolume(unsigned compIdx)
Molar volume of a component at the critical point [m^3/mol].
Definition Spe5FluidSystem.hpp:239
static Scalar acentricFactor(unsigned compIdx)
The acentric factor of a component [].
Definition Spe5FluidSystem.hpp:262
static LhsEval viscosity(const FluidState &, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx)
Calculate the dynamic viscosity of a fluid phase [Pa*s].
Definition Spe5FluidSystem.hpp:375
static void init(Scalar minT=273.15, Scalar maxT=373.15, Scalar minP=1e4, Scalar maxP=100e6)
Initialize the fluid system's static parameters.
Definition Spe5FluidSystem.hpp:310
static const int C1Idx
Index of the C1 component.
Definition Spe5FluidSystem.hpp:145
static bool isCompressible(unsigned)
Returns true if and only if a fluid phase is assumed to be compressible.
Definition Spe5FluidSystem.hpp:115
static const int C20Idx
Index of the C20 component.
Definition Spe5FluidSystem.hpp:150
static const int C6Idx
Index of the C6 component.
Definition Spe5FluidSystem.hpp:147
static const int C10Idx
Index of the C10 component.
Definition Spe5FluidSystem.hpp:148
static Scalar interactionCoefficient(unsigned comp1Idx, unsigned comp2Idx)
Returns the interaction coefficient for two components.
Definition Spe5FluidSystem.hpp:287
static const int waterPhaseIdx
Index of the water phase.
Definition Spe5FluidSystem.hpp:83
static Scalar molarMass(unsigned compIdx)
Return the molar mass of a component in [kg/mol].
Definition Spe5FluidSystem.hpp:170
::Opm::H2O< Scalar > H2O
The component for pure water to be used.
Definition Spe5FluidSystem.hpp:88
static const int C3Idx
Index of the C3 component.
Definition Spe5FluidSystem.hpp:146
static bool isLiquid(unsigned phaseIdx)
Return whether a phase is liquid.
Definition Spe5FluidSystem.hpp:104
static const int gasPhaseIdx
Index of the gas phase.
Definition Spe5FluidSystem.hpp:81
static std::string_view componentName(unsigned compIdx)
Return the human readable name of a component.
Definition Spe5FluidSystem.hpp:153
static std::string_view phaseName(unsigned phaseIdx)
Return the human readable name of a fluid phase.
Definition Spe5FluidSystem.hpp:91
static bool isIdealGas(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal gas.
Definition Spe5FluidSystem.hpp:122
Specifies the parameter cache used by the SPE-5 fluid system.
Definition Spe5ParameterCache.hpp:47
Scalar molarVolume(unsigned phaseIdx) const
Returns the molar volume of a phase [m^3/mol].
Definition Spe5ParameterCache.hpp:202
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30
The type of the fluid system's parameter cache.
Definition Spe5FluidSystem.hpp:71