My Project
Loading...
Searching...
No Matches
Brine_H2.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*/
28#ifndef OPM_BINARY_COEFF_BRINE_H2_HPP
29#define OPM_BINARY_COEFF_BRINE_H2_HPP
30
32
33namespace Opm {
34namespace BinaryCoeff {
35
40template<class Scalar, class H2O, class H2, bool verbose = true>
41class Brine_H2 {
42 static const int liquidPhaseIdx = 0; // index of the liquid phase
43 static const int gasPhaseIdx = 1; // index of the gas phase
44
45public:
56 template <class Evaluation>
57 static Evaluation calculateMoleFractions(const Evaluation& temperature,
58 const Evaluation& pg,
59 const Evaluation& salinity,
60 bool extrapolate = false)
61 {
62 // Convert salinity to molality from mass fraction
63 Evaluation X_NaCl = salinityToMolality_(salinity);
64
65 // All intermediate calculations
66 Evaluation lnYH2 = moleFractionGasH2_(temperature, pg);
67 Evaluation lnPg = log(pg / 1e6); // Pa --> MPa before ln
68 Evaluation lnPhiH2 = fugacityCoefficientH2(temperature, pg, extrapolate);
69 Evaluation lnKh = henrysConstant_(temperature);
70 Evaluation PF = computePoyntingFactor_(temperature, pg);
71 Evaluation lnGammaH2 = activityCoefficient_(temperature, X_NaCl);
72
73 // Eq. (4) to get mole fraction of H2 in brine
74 Evaluation xH2 = exp(lnYH2 + lnPg + lnPhiH2 - lnKh - PF - lnGammaH2);
75 return xH2;
76 }
77
84 template <class Evaluation>
85 static Evaluation computePoyntingFactor_(const Evaluation& temperature, const Evaluation& pg)
86 {
87 // PF is approximated as a polynomial expansion in terms of temperature and pressure with the following
88 // parameters (Table 4)
89 static const Scalar a[4] = {6.156755, -2.502396e-2, 4.140593e-5, -1.322988e-3};
90
91 // Eq. (16)
92 Evaluation pg_mpa = pg / 1.0e6; // convert from Pa to MPa
93 Evaluation PF = a[0]*pg_mpa/temperature + a[1]*pg_mpa + a[2]*temperature*pg_mpa + a[3]*pg_mpa*pg_mpa/temperature;
94 return PF;
95 }
96
105 template <class Evaluation>
106 static Evaluation activityCoefficient_(const Evaluation& temperature, const Evaluation& salinity)
107 {
108 // Linear approximation in temperature with following parameters (Table 5)
109 static const Scalar a[2] = {0.64485, 0.00142};
110
111 // Eq. (17)
112 Evaluation lnGamma = (a[0] - a[1]*temperature)*salinity;
113 return lnGamma;
114 }
115
121 template <class Evaluation>
122 static Evaluation henrysConstant_(const Evaluation& temperature)
123 {
124 // Polynomic approximation in temperature with following parameters (Table 2)
125 static const Scalar a[5] = {2.68721e-5, -0.05121, 33.55196, -3411.0432, -31258.74683};
126
127 // Eq. (13)
128 Evaluation lnKh = a[0]*temperature*temperature + a[1]*temperature + a[2] + a[3]/temperature
129 + a[4]/(temperature*temperature);
130 return lnKh;
131 }
132
140 template <class Evaluation>
141 static Evaluation moleFractionGasH2_(const Evaluation& temperature, const Evaluation& pg)
142 {
143 // Need saturaturated vapor pressure of pure water
144 Evaluation pw_sat = H2O::vaporPressure(temperature);
145
146 // Eq. (12)
147 Evaluation lnyH2 = log(1 - (pw_sat / pg));
148 return lnyH2;
149 }
150
159 template <class Evaluation>
160 static Evaluation fugacityCoefficientH2(const Evaluation& temperature,
161 const Evaluation& pg,
162 bool extrapolate = false)
163 {
164 // Convert pressure to reduced density and temperature to reduced temperature
165 Evaluation rho_red = H2::reducedMolarDensity(temperature, pg, extrapolate);
166 Evaluation T_red = H2::criticalTemperature() / temperature;
167
168 // Residual Helmholtz energy, Eq. (7) in Li et al. (2018)
169 Evaluation resHelm = H2::residualPartHelmholtz(T_red, rho_red);
170
171 // Derivative of residual Helmholtz energy wrt to reduced density, Eq. (73) in Span et al. (2018)
172 Evaluation dResHelm_dRedRho = H2::derivResHelmholtzWrtRedRho(T_red, rho_red);
173
174 // Fugacity coefficient, Eq. (8) in Li et al. (2018)
175 Evaluation lnPhiH2 = resHelm + rho_red * dResHelm_dRedRho - log(rho_red * dResHelm_dRedRho + 1);
176
177 return lnPhiH2;
178 }
179
185 template <class Evaluation>
186 static Evaluation gasDiffCoeff(const Evaluation& temperature, const Evaluation& pressure)
187 {
188 // atomic diffusion volumes
189 const Scalar SigmaNu[2] = { 13.1 /* H2O */, 7.07 /* H2 */ };
190 // molar masses [g/mol]
191 const Scalar M[2] = { H2O::molarMass()*Scalar(1e3), H2::molarMass()*Scalar(1e3) };
192
193 return fullerMethod(M, SigmaNu, temperature, pressure);
194 }
195
196private:
201 template <class Evaluation>
202 static Evaluation salinityToMolality_(const Evaluation& salinity) {
203 // Molar mass NaCl
204 const Scalar M_NaCl = 58.44e-3;
205
206 // Convert
207 const Evaluation X_NaCl = salinity / ((1 - salinity) * M_NaCl);
208 return X_NaCl;
209
210 }
211}; // end class Brine_H2
212
213} // end namespace BinaryCoeff
214} // end namespace Opm
215#endif
Estimate binary diffusion coefficents in gases according to the method by Fuller.
Evaluation fullerMethod(const Scalar *M, const Scalar *SigmaNu, const Evaluation &temperature, const Evaluation &pressure)
Estimate binary diffusion coefficents in gases according to the method by Fuller.
Definition FullerMethod.hpp:56
Binary coefficients for brine and H2.
Definition Brine_H2.hpp:41
static Evaluation calculateMoleFractions(const Evaluation &temperature, const Evaluation &pg, const Evaluation &salinity, bool extrapolate=false)
Returns the mol (!) fraction of H2 in the liquid phase for a given temperature, pressure,...
Definition Brine_H2.hpp:57
static Evaluation activityCoefficient_(const Evaluation &temperature, const Evaluation &salinity)
Returns the activity coefficient of H2 in brine which is needed in calculation of H2 solubility in Li...
Definition Brine_H2.hpp:106
static Evaluation gasDiffCoeff(const Evaluation &temperature, const Evaluation &pressure)
Binary diffusion coefficent [m^2/s] for molecular water and H2 as an approximation for brine-H2 diffu...
Definition Brine_H2.hpp:186
static Evaluation henrysConstant_(const Evaluation &temperature)
Returns Henry's constant of H2 in brine which is needed in calculation of H2 solubility in Li et al (...
Definition Brine_H2.hpp:122
static Evaluation moleFractionGasH2_(const Evaluation &temperature, const Evaluation &pg)
Returns mole fraction of H2 in gasous phase which is needed in calculation of H2 solubility in Li et ...
Definition Brine_H2.hpp:141
static Evaluation computePoyntingFactor_(const Evaluation &temperature, const Evaluation &pg)
Returns the Poynting Factor (PF) which is needed in calculation of H2 solubility in Li et al (2018).
Definition Brine_H2.hpp:85
static Evaluation fugacityCoefficientH2(const Evaluation &temperature, const Evaluation &pg, bool extrapolate=false)
Calculate fugacity coefficient for H2 which is needed in calculation of H2 solubility in Li et al (20...
Definition Brine_H2.hpp:160
static Evaluation vaporPressure(Evaluation temperature)
The vapor pressure in of pure water at a given temperature.
Definition H2O.hpp:141
static const Scalar molarMass()
The molar mass in of water.
Definition H2O.hpp:83
static Evaluation residualPartHelmholtz(const Evaluation &T_red, const Evaluation &rho_red)
The residual part of Helmholtz energy.
Definition H2.hpp:464
static Scalar criticalTemperature()
Returns the critical temperature of molecular hydrogen.
Definition H2.hpp:83
static Evaluation reducedMolarDensity(const Evaluation &temperature, const Evaluation &pg, bool extrapolate=false)
Calculate reduced density (rho/rho_crit) from pressure and temperature.
Definition H2.hpp:380
static Evaluation derivResHelmholtzWrtRedRho(const Evaluation &T_red, const Evaluation &rho_red)
Derivative of the residual part of Helmholtz energy wrt.
Definition H2.hpp:498
static constexpr Scalar molarMass()
The molar mass in of molecular hydrogen.
Definition H2.hpp:77
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30