My Project
Loading...
Searching...
No Matches
H2OAirFluidSystem.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_H2O_AIR_SYSTEM_HPP
28#define OPM_H2O_AIR_SYSTEM_HPP
29
30#include "BaseFluidSystem.hpp"
32
38
40
41#include <cassert>
42
43namespace Opm {
44
54template <class Scalar,
55 //class H2Otype = SimpleH2O<Scalar>,
56 class H2Otype = TabulatedComponent<Scalar, H2O<Scalar> >>
58 : public BaseFluidSystem<Scalar, H2OAirFluidSystem<Scalar, H2Otype> >
59{
61 typedef BaseFluidSystem <Scalar, ThisType> Base;
62 typedef ::Opm::IdealGas<Scalar> IdealGas;
63
64public:
65 template <class Evaluation>
66 struct ParameterCache : public NullParameterCache<Evaluation>
67 {};
68
70 typedef H2Otype H2O;
72 typedef ::Opm::Air<Scalar> Air;
73
75 static const int numPhases = 2;
76
78 static const int liquidPhaseIdx = 0;
80 static const int gasPhaseIdx = 1;
81
83 static const char* phaseName(unsigned phaseIdx)
84 {
85 switch (phaseIdx) {
86 case liquidPhaseIdx: return "liquid";
87 case gasPhaseIdx: return "gas";
88 };
89 throw std::logic_error("Invalid phase index "+std::to_string(phaseIdx));
90 }
91
93 static bool isLiquid(unsigned phaseIdx)
94 {
95 //assert(0 <= phaseIdx && phaseIdx < numPhases);
96 return phaseIdx != gasPhaseIdx;
97 }
98
100 static bool isCompressible(unsigned phaseIdx)
101 {
102 //assert(0 <= phaseIdx && phaseIdx < numPhases);
103 return (phaseIdx == gasPhaseIdx)
104 // ideal gases are always compressible
105 ? true
106 :
107 // the water component decides for the liquid phase...
109 }
110
112 static bool isIdealGas(unsigned phaseIdx)
113 {
114 return
115 (phaseIdx == gasPhaseIdx)
117 : false;
118 }
119
121 static bool isIdealMixture(unsigned /*phaseIdx*/)
122 {
123 //assert(0 <= phaseIdx && phaseIdx < numPhases);
124 // we assume Henry's and Rault's laws for the water phase and
125 // and no interaction between gas molecules of different
126 // components, so all phases are ideal mixtures!
127 return true;
128 }
129
130 /****************************************
131 * Component related static parameters
132 ****************************************/
133
135 static const int numComponents = 2;
136
138 static const int H2OIdx = 0;
140 static const int AirIdx = 1;
141
143 static const char* componentName(unsigned compIdx)
144 {
145 switch (compIdx)
146 {
147 case H2OIdx: return H2O::name();
148 case AirIdx: return Air::name();
149 };
150 throw std::logic_error("Invalid component index "+std::to_string(compIdx));
151 }
152
154 static Scalar molarMass(unsigned compIdx)
155 {
156 return
157 (compIdx == H2OIdx)
159 : (compIdx == AirIdx)
161 : 1e30;
162 }
163
164
170 static Scalar criticalTemperature(unsigned compIdx)
171 {
172 return
173 (compIdx == H2OIdx)
175 : (compIdx == AirIdx)
177 : 1e30;
178 }
179
185 static Scalar criticalPressure(unsigned compIdx)
186 {
187 return
188 (compIdx == H2OIdx)
190 : (compIdx == AirIdx)
192 : 1e30;
193 }
194
200 static Scalar acentricFactor(unsigned compIdx)
201 {
202 return
203 (compIdx == H2OIdx)
205 : (compIdx == AirIdx)
207 : 1e30;
208 }
209
210 /****************************************
211 * thermodynamic relations
212 ****************************************/
213
220 static void init()
221 {
222 if (H2O::isTabulated)
223 init(/*tempMin=*/273.15,
224 /*tempMax=*/623.15,
225 /*numTemp=*/50,
226 /*pMin=*/-10,
227 /*pMax=*/20e6,
228 /*numP=*/50);
229 }
230
242 static void init(Scalar tempMin, Scalar tempMax, unsigned nTemp,
243 Scalar pressMin, Scalar pressMax, unsigned nPress)
244 {
245 if (H2O::isTabulated) {
246 H2O::init(tempMin, tempMax, nTemp,
247 pressMin, pressMax, nPress);
248 }
249 }
250
252 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
253 static LhsEval density(const FluidState& fluidState,
254 const ParameterCache<ParamCacheEval>& /*paramCache*/,
255 unsigned phaseIdx)
256 {
257 assert(phaseIdx < numPhases);
258
259 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
260 LhsEval p;
261 if (isCompressible(phaseIdx))
262 p = decay<LhsEval>(fluidState.pressure(phaseIdx));
263 else {
264 // random value which will hopefully cause things to blow
265 // up if it is used in a calculation!
266 p = - 1e30;
267 Valgrind::SetUndefined(p);
268 }
269
270
271 LhsEval sumMoleFrac = 0;
272 for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx)
273 sumMoleFrac += decay<LhsEval>(fluidState.moleFraction(phaseIdx, compIdx));
274
275 if (phaseIdx == liquidPhaseIdx)
276 {
277 // assume ideal mixture: Molecules of one component don't discriminate
278 // between their own kind and molecules of the other component.
279 const LhsEval& clH2O = H2O::liquidDensity(T, p)/H2O::molarMass();
280
281 const auto& xlH2O = decay<LhsEval>(fluidState.moleFraction(liquidPhaseIdx, H2OIdx));
282 const auto& xlAir = decay<LhsEval>(fluidState.moleFraction(liquidPhaseIdx, AirIdx));
283
284 return clH2O*(H2O::molarMass()*xlH2O + Air::molarMass()*xlAir)/sumMoleFrac;
285 }
286 else if (phaseIdx == gasPhaseIdx)
287 {
288 LhsEval partialPressureH2O =
289 decay<LhsEval>(fluidState.moleFraction(gasPhaseIdx, H2OIdx))
290 *decay<LhsEval>(fluidState.pressure(gasPhaseIdx));
291
292 LhsEval partialPressureAir =
293 decay<LhsEval>(fluidState.moleFraction(gasPhaseIdx, AirIdx))
294 *decay<LhsEval>(fluidState.pressure(gasPhaseIdx));
295
296 return H2O::gasDensity(T, partialPressureH2O) + Air::gasDensity(T, partialPressureAir);
297 }
298 throw std::logic_error("Invalid phase index "+std::to_string(phaseIdx));
299 }
300
302 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
303 static LhsEval viscosity(const FluidState& fluidState,
304 const ParameterCache<ParamCacheEval>& /*paramCache*/,
305 unsigned phaseIdx)
306 {
307 assert(phaseIdx < numPhases);
308
309 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
310 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
311
312 if (phaseIdx == liquidPhaseIdx)
313 {
314 // assume pure water for the liquid phase
315 // TODO: viscosity of mixture
316 // couldn't find a way to solve the mixture problem
317 return H2O::liquidViscosity(T, p);
318 }
319 else if (phaseIdx == gasPhaseIdx)
320 {
321 /* Wilke method. See:
322 *
323 * See: R. Reid, et al.: The Properties of Gases and Liquids,
324 * 4th edition, McGraw-Hill, 1987, 407-410 or
325 * 5th edition, McGraw-Hill, 2000, p. 9.21/22
326 */
327 LhsEval muResult = 0;
328 const LhsEval mu[numComponents] = {
331 };
332
333 for (unsigned i = 0; i < numComponents; ++i) {
334 LhsEval divisor = 0;
335 for (unsigned j = 0; j < numComponents; ++j) {
336 LhsEval phiIJ =
337 1 +
338 sqrt(mu[i]/mu[j]) * // 1 + (mu[i]/mu[j]^1/2
339 std::pow(molarMass(j)/molarMass(i), 1./4.0); // (M[i]/M[j])^1/4
340
341 phiIJ *= phiIJ;
342 phiIJ /= std::sqrt(8*(1 + molarMass(i)/molarMass(j)));
343 divisor += decay<LhsEval>(fluidState.moleFraction(phaseIdx, j))*phiIJ;
344 }
345 const auto& xAlphaI = decay<LhsEval>(fluidState.moleFraction(phaseIdx, i));
346 muResult += xAlphaI*mu[i]/divisor;
347 }
348 return muResult;
349 }
350 throw std::logic_error("Invalid phase index "+std::to_string(phaseIdx));
351 }
352
354 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
355 static LhsEval fugacityCoefficient(const FluidState& fluidState,
356 const ParameterCache<ParamCacheEval>& /*paramCache*/,
357 unsigned phaseIdx,
358 unsigned compIdx)
359 {
360 assert(phaseIdx < numPhases);
361 assert(compIdx < numComponents);
362
363 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
364 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
365
366 if (phaseIdx == liquidPhaseIdx) {
367 if (compIdx == H2OIdx)
368 return H2O::vaporPressure(T)/p;
369 return BinaryCoeff::H2O_Air::henry(T)/p;
370 }
371
372 // for the gas phase, assume an ideal gas when it comes to
373 // fugacity (-> fugacity == partial pressure)
374 return 1.0;
375 }
376
378 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
379 static LhsEval binaryDiffusionCoefficient(const FluidState& fluidState,
380 const ParameterCache<ParamCacheEval>& /*paramCache*/,
381 unsigned phaseIdx,
382 unsigned /*compIdx*/)
383 {
384 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
385 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
386
387 if (phaseIdx == liquidPhaseIdx)
389
390 assert(phaseIdx == gasPhaseIdx);
392 }
393
395 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
396 static LhsEval enthalpy(const FluidState& fluidState,
397 const ParameterCache<ParamCacheEval>& /*paramCache*/,
398 unsigned phaseIdx)
399 {
400 const auto& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
401 const auto& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
402 Valgrind::CheckDefined(T);
403 Valgrind::CheckDefined(p);
404
405 if (phaseIdx == liquidPhaseIdx)
406 {
407 // TODO: correct way to deal with the solutes???
408 return H2O::liquidEnthalpy(T, p);
409 }
410
411 else if (phaseIdx == gasPhaseIdx)
412 {
413 LhsEval result = 0.0;
414 result +=
415 H2O::gasEnthalpy(T, p) *
416 decay<LhsEval>(fluidState.massFraction(gasPhaseIdx, H2OIdx));
417
418 result +=
419 Air::gasEnthalpy(T, p) *
420 decay<LhsEval>(fluidState.massFraction(gasPhaseIdx, AirIdx));
421 return result;
422 }
423 throw std::logic_error("Invalid phase index "+std::to_string(phaseIdx));
424 }
425
427 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
428 static LhsEval thermalConductivity(const FluidState& fluidState,
429 const ParameterCache<ParamCacheEval>& /*paramCache*/,
430 unsigned phaseIdx)
431 {
432 assert(phaseIdx < numPhases);
433
434 const LhsEval& temperature =
435 decay<LhsEval>(fluidState.temperature(phaseIdx));
436 const LhsEval& pressure =
437 decay<LhsEval>(fluidState.pressure(phaseIdx));
438
439 if (phaseIdx == liquidPhaseIdx)
440 return H2O::liquidThermalConductivity(temperature, pressure);
441 else { // gas phase
442 const LhsEval& lambdaDryAir = Air::gasThermalConductivity(temperature, pressure);
443
444 const LhsEval& xAir =
445 decay<LhsEval>(fluidState.moleFraction(phaseIdx, AirIdx));
446 const LhsEval& xH2O =
447 decay<LhsEval>(fluidState.moleFraction(phaseIdx, H2OIdx));
448 LhsEval lambdaAir = xAir*lambdaDryAir;
449
450 // Assuming Raoult's, Daltons law and ideal gas
451 // in order to obtain the partial density of water in the air phase
452 LhsEval partialPressure = pressure*xH2O;
453
454 LhsEval lambdaH2O =
455 xH2O*H2O::gasThermalConductivity(temperature, partialPressure);
456 return lambdaAir + lambdaH2O;
457 }
458 }
459};
460
461} // namespace Opm
462
463#endif
A simple class implementing the fluid properties of air.
The base class for all fluid systems.
Material properties of pure water .
Binary coefficients for water and nitrogen.
Relations valid for an ideal gas.
A parameter cache which does nothing.
A generic class which tabulates all thermodynamic properties of a given component.
Some templates to wrap the valgrind client request macros.
static Evaluation gasThermalConductivity(const Evaluation &, const Evaluation &)
Specific heat conductivity of steam .
Definition Air.hpp:217
static Scalar criticalPressure()
Returns the critical pressure of .
Definition Air.hpp:92
static Evaluation gasDensity(const Evaluation &temperature, const Evaluation &pressure)
The density of at a given pressure and temperature [kg/m^3].
Definition Air.hpp:102
static Scalar molarMass()
The molar mass in of .
Definition Air.hpp:80
static const char * name()
A human readable name for the .
Definition Air.hpp:60
static bool gasIsIdeal()
Returns true iff the gas phase is assumed to be ideal.
Definition Air.hpp:72
static Evaluation gasEnthalpy(const Evaluation &temperature, const Evaluation &)
Specific enthalpy of liquid water with 273.15 K as basis.
Definition Air.hpp:180
static Scalar criticalTemperature()
Returns the critical temperature of .
Definition Air.hpp:86
static Evaluation gasViscosity(const Evaluation &temperature, const Evaluation &)
The dynamic viscosity of at a given pressure and temperature.
Definition Air.hpp:137
The base class for all fluid systems.
Definition BaseFluidSystem.hpp:46
ScalarT Scalar
The type used for scalar quantities.
Definition BaseFluidSystem.hpp:51
static Evaluation henry(const Evaluation &temperature)
Henry coefficent for air in liquid water.
Definition H2O_Air.hpp:55
static Evaluation liquidDiffCoeff(const Evaluation &temperature, const Evaluation &)
Lacking better data on water-air diffusion in liquids, we use at the moment the diffusion coefficient...
Definition H2O_Air.hpp:101
static Evaluation gasDiffCoeff(const Evaluation &temperature, const Evaluation &pressure)
Binary diffusion coefficent for molecular water and air.
Definition H2O_Air.hpp:70
static Scalar acentricFactor()
Returns the acentric factor of the component.
Definition Component.hpp:110
static void init(Scalar, Scalar, unsigned, Scalar, Scalar, unsigned)
A default routine for initialization, not needed for components and must not be called.
Definition Component.hpp:60
A fluid system with a liquid and a gaseous phase and water and air as components.
Definition H2OAirFluidSystem.hpp:59
H2Otype H2O
The type of the water component used for this fluid system.
Definition H2OAirFluidSystem.hpp:70
static LhsEval density(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx)
Calculate the density [kg/m^3] of a fluid phase.
Definition H2OAirFluidSystem.hpp:253
static const int AirIdx
The index of the air component.
Definition H2OAirFluidSystem.hpp:140
static const int numComponents
Number of chemical species in the fluid system.
Definition H2OAirFluidSystem.hpp:135
static Scalar molarMass(unsigned compIdx)
Return the molar mass of a component in [kg/mol].
Definition H2OAirFluidSystem.hpp:154
static const int gasPhaseIdx
The index of the gas phase.
Definition H2OAirFluidSystem.hpp:80
static bool isIdealGas(unsigned phaseIdx)
Returns true if and only if a fluid phase is assumed to be an ideal gas.
Definition H2OAirFluidSystem.hpp:112
static const char * componentName(unsigned compIdx)
Return the human readable name of a component.
Definition H2OAirFluidSystem.hpp:143
static const int liquidPhaseIdx
The index of the liquid phase.
Definition H2OAirFluidSystem.hpp:78
static LhsEval binaryDiffusionCoefficient(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx, unsigned)
Calculate the binary molecular diffusion coefficient for a component in a fluid phase [mol^2 * s / (k...
Definition H2OAirFluidSystem.hpp:379
static const int H2OIdx
The index of the water component.
Definition H2OAirFluidSystem.hpp:138
static bool isCompressible(unsigned phaseIdx)
Returns true if and only if a fluid phase is assumed to be compressible.
Definition H2OAirFluidSystem.hpp:100
static const int numPhases
Number of fluid phases in the fluid system.
Definition H2OAirFluidSystem.hpp:75
static const char * phaseName(unsigned phaseIdx)
Return the human readable name of a fluid phase.
Definition H2OAirFluidSystem.hpp:83
static LhsEval viscosity(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx)
Calculate the dynamic viscosity of a fluid phase [Pa*s].
Definition H2OAirFluidSystem.hpp:303
static bool isIdealMixture(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal mixture.
Definition H2OAirFluidSystem.hpp:121
static void init(Scalar tempMin, Scalar tempMax, unsigned nTemp, Scalar pressMin, Scalar pressMax, unsigned nPress)
Initialize the fluid system's static parameters using problem specific temperature and pressure range...
Definition H2OAirFluidSystem.hpp:242
::Opm::Air< Scalar > Air
The type of the air component used for this fluid system.
Definition H2OAirFluidSystem.hpp:72
static Scalar criticalTemperature(unsigned compIdx)
Critical temperature of a component [K].
Definition H2OAirFluidSystem.hpp:170
static bool isLiquid(unsigned phaseIdx)
Return whether a phase is liquid.
Definition H2OAirFluidSystem.hpp:93
static LhsEval enthalpy(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx)
Given a phase's composition, temperature, pressure and density, calculate its specific enthalpy [J/kg...
Definition H2OAirFluidSystem.hpp:396
static Scalar acentricFactor(unsigned compIdx)
The acentric factor of a component [].
Definition H2OAirFluidSystem.hpp:200
static void init()
Initialize the fluid system's static parameters.
Definition H2OAirFluidSystem.hpp:220
static Scalar criticalPressure(unsigned compIdx)
Critical pressure of a component [Pa].
Definition H2OAirFluidSystem.hpp:185
static LhsEval thermalConductivity(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx)
Thermal conductivity of a fluid phase [W/(m K)].
Definition H2OAirFluidSystem.hpp:428
static LhsEval fugacityCoefficient(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx, unsigned compIdx)
Calculate the fugacity coefficient [Pa] of an individual component in a fluid phase.
Definition H2OAirFluidSystem.hpp:355
static Evaluation liquidDensity(const Evaluation &temperature, const Evaluation &pressure, bool extrapolate=false)
The density of pure water in at a given pressure and temperature.
Definition H2O.hpp:687
static const Scalar criticalTemperature()
Returns the critical temperature of water.
Definition H2O.hpp:95
static Evaluation gasDensity(const Evaluation &temperature, const Evaluation &pressure)
The density of steam in at a given pressure and temperature.
Definition H2O.hpp:562
static Evaluation vaporPressure(Evaluation temperature)
The vapor pressure in of pure water at a given temperature.
Definition H2O.hpp:141
static Evaluation gasViscosity(const Evaluation &temperature, const Evaluation &pressure)
The dynamic viscosity of steam.
Definition H2O.hpp:790
static Evaluation liquidEnthalpy(const Evaluation &temperature, const Evaluation &pressure)
Specific enthalpy of liquid water .
Definition H2O.hpp:237
static Evaluation liquidThermalConductivity(const Evaluation &temperature, const Evaluation &pressure)
Thermal conductivity of water (IAPWS) .
Definition H2O.hpp:844
static const Scalar acentricFactor()
The acentric factor of water.
Definition H2O.hpp:89
static bool gasIsIdeal()
Returns true iff the gas phase is assumed to be ideal.
Definition H2O.hpp:627
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
static bool liquidIsCompressible()
Returns true iff the liquid phase is assumed to be compressible.
Definition H2O.hpp:546
static Evaluation liquidViscosity(const Evaluation &temperature, const Evaluation &pressure, bool extrapolate=false)
The dynamic viscosity of pure water.
Definition H2O.hpp:815
static Evaluation gasThermalConductivity(const Evaluation &temperature, const Evaluation &pressure)
Thermal conductivity of water (IAPWS) .
Definition H2O.hpp:864
static Evaluation gasEnthalpy(const Evaluation &temperature, const Evaluation &pressure)
Specific enthalpy of water steam .
Definition H2O.hpp:186
static const char * name()
A human readable name for the water.
Definition H2O.hpp:77
Relations valid for an ideal gas.
Definition IdealGas.hpp:38
A parameter cache which does nothing.
Definition NullParameterCache.hpp:40
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30
Definition H2OAirFluidSystem.hpp:67