My Project
Loading...
Searching...
No Matches
Brine_CO2.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_CO2_HPP
29#define OPM_BINARY_COEFF_BRINE_CO2_HPP
30
33
34namespace Opm {
35namespace BinaryCoeff {
36
41template<class Scalar, class H2O, class CO2, bool verbose = true>
42class Brine_CO2 {
43 typedef ::Opm::IdealGas<Scalar> IdealGas;
44 static const int liquidPhaseIdx = 0; // index of the liquid phase
45 static const int gasPhaseIdx = 1; // index of the gas phase
46
47public:
55 template <class Evaluation>
56 static Evaluation gasDiffCoeff(const Evaluation& temperature, const Evaluation& pressure, bool extrapolate = false)
57 {
58 //Diffusion coefficient of water in the CO2 phase
59 Scalar k = 1.3806504e-23; // Boltzmann constant
60 Scalar c = 4; // slip parameter, can vary between 4 (slip condition) and 6 (stick condition)
61 Scalar R_h = 1.72e-10; // hydrodynamic radius of the solute
62 const Evaluation& mu = CO2::gasViscosity(temperature, pressure, extrapolate); // CO2 viscosity
63 return k / (c * M_PI * R_h) * (temperature / mu);
64 }
65
72 template <class Evaluation>
73 static Evaluation liquidDiffCoeff(const Evaluation& /*temperature*/, const Evaluation& /*pressure*/)
74 {
75 //Diffusion coefficient of CO2 in the brine phase
76 return 2e-9;
77 }
78
96 template <class Evaluation>
97 static void calculateMoleFractions(const Evaluation& temperature,
98 const Evaluation& pg,
99 const Evaluation& salinity,
100 const int knownPhaseIdx,
101 Evaluation& xlCO2,
102 Evaluation& ygH2O,
103 bool extrapolate = false)
104 {
105 OPM_TIMEFUNCTION_LOCAL();
106 Evaluation A = computeA_(temperature, pg, extrapolate);
107
108 /* salinity: conversion from mass fraction to mol fraction */
109 Evaluation x_NaCl = salinityToMolFrac_(salinity);
110
111 // if both phases are present the mole fractions in each phase can be calculate
112 // with the mutual solubility function
113 if (knownPhaseIdx < 0) {
114 Evaluation molalityNaCl = moleFracToMolality_(x_NaCl); // molality of NaCl //CHANGED
115 Evaluation m0_CO2 = molalityCO2inPureWater_(temperature, pg, extrapolate); // molality of CO2 in pure water
116 Evaluation gammaStar = activityCoefficient_(temperature, pg, molalityNaCl);// activity coefficient of CO2 in brine
117 Evaluation m_CO2 = m0_CO2 / gammaStar; // molality of CO2 in brine
118 xlCO2 = m_CO2 / (molalityNaCl + 55.508 + m_CO2); // mole fraction of CO2 in brine
119 ygH2O = A * (1 - xlCO2 - x_NaCl); // mole fraction of water in the gas phase
120 }
121
122 // if only liquid phase is present the mole fraction of CO2 in brine is given and
123 // and the virtual equilibrium mole fraction of water in the non-existing gas phase can be estimated
124 // with the mutual solubility function
125 if (knownPhaseIdx == liquidPhaseIdx)
126 ygH2O = A * (1 - xlCO2 - x_NaCl);
127
128 // if only gas phase is present the mole fraction of water in the gas phase is given and
129 // and the virtual equilibrium mole fraction of CO2 in the non-existing liquid phase can be estimated
130 // with the mutual solubility function
131 if (knownPhaseIdx == gasPhaseIdx)
132 //y_H2o = fluidstate.
133 xlCO2 = 1 - x_NaCl - ygH2O / A;
134 }
135
139 template <class Evaluation>
140 static Evaluation henry(const Evaluation& temperature, bool extrapolate = false)
141 { return fugacityCoefficientCO2(temperature, /*pressure=*/1e5, extrapolate)*1e5; }
142
151 template <class Evaluation>
152 static Evaluation fugacityCoefficientCO2(const Evaluation& temperature, const Evaluation& pg, bool extrapolate = false)
153 {
154 OPM_TIMEFUNCTION_LOCAL();
155 Valgrind::CheckDefined(temperature);
156 Valgrind::CheckDefined(pg);
157
158 Evaluation V = 1 / (CO2::gasDensity(temperature, pg, extrapolate) / CO2::molarMass()) * 1.e6; // molar volume in cm^3/mol
159 Evaluation pg_bar = pg / 1.e5; // gas phase pressure in bar
160 Evaluation a_CO2 = (7.54e7 - 4.13e4 * temperature); // mixture parameter of Redlich-Kwong equation
161 Scalar b_CO2 = 27.8; // mixture parameter of Redlich-Kwong equation
162 Scalar R = IdealGas::R * 10.; // ideal gas constant with unit bar cm^3 /(K mol)
163 Evaluation lnPhiCO2;
164
165 lnPhiCO2 = log(V / (V - b_CO2));
166 lnPhiCO2 += b_CO2 / (V - b_CO2);
167 lnPhiCO2 -= 2 * a_CO2 / (R * pow(temperature, 1.5) * b_CO2) * log((V + b_CO2) / V);
168 lnPhiCO2 +=
169 a_CO2 * b_CO2
170 / (R
171 * pow(temperature, 1.5)
172 * b_CO2
173 * b_CO2)
174 * (log((V + b_CO2) / V)
175 - b_CO2 / (V + b_CO2));
176 lnPhiCO2 -= log(pg_bar * V / (R * temperature));
177
178 return exp(lnPhiCO2); // fugacity coefficient of CO2
179 }
180
189 template <class Evaluation>
190 static Evaluation fugacityCoefficientH2O(const Evaluation& temperature, const Evaluation& pg, bool extrapolate = false)
191 {
192 OPM_TIMEFUNCTION_LOCAL();
193 const Evaluation& V = 1 / (CO2::gasDensity(temperature, pg, extrapolate) / CO2::molarMass()) * 1.e6; // molar volume in cm^3/mol
194 const Evaluation& pg_bar = pg / 1.e5; // gas phase pressure in bar
195 const Evaluation& a_CO2 = (7.54e7 - 4.13e4 * temperature);// mixture parameter of Redlich-Kwong equation
196 Scalar a_CO2_H2O = 7.89e7;// mixture parameter of Redlich-Kwong equation
197 Scalar b_CO2 = 27.8;// mixture parameter of Redlich-Kwong equation
198 Scalar b_H2O = 18.18;// mixture parameter of Redlich-Kwong equation
199 Scalar R = IdealGas::R * 10.; // ideal gas constant with unit bar cm^3 /(K mol)
200 Evaluation lnPhiH2O;
201
202 lnPhiH2O =
203 log(V/(V - b_CO2))
204 + b_H2O/(V - b_CO2) - 2*a_CO2_H2O
205 / (R*pow(temperature, 1.5)*b_CO2)*log((V + b_CO2)/V)
206 + a_CO2*b_H2O/(R*pow(temperature, 1.5)*b_CO2*b_CO2)
207 *(log((V + b_CO2)/V) - b_CO2/(V + b_CO2))
208 - log(pg_bar*V/(R*temperature));
209 return exp(lnPhiH2O); // fugacity coefficient of H2O
210 }
211
212private:
218 template <class Evaluation>
219 static Evaluation salinityToMolFrac_(const Evaluation& salinity) {
220 OPM_TIMEFUNCTION_LOCAL();
221 const Scalar Mw = H2O::molarMass(); /* molecular weight of water [kg/mol] */
222 const Scalar Ms = 58.44e-3; /* molecular weight of NaCl [kg/mol] */
223
224 const Evaluation X_NaCl = salinity;
225 /* salinity: conversion from mass fraction to mol fraction */
226 const Evaluation x_NaCl = -Mw * X_NaCl / ((Ms - Mw) * X_NaCl - Ms);
227 return x_NaCl;
228 }
229
235 template <class Evaluation>
236 static Evaluation moleFracToMolality_(const Evaluation& x_NaCl)
237 {
238 // conversion from mol fraction to molality (dissolved CO2 neglected)
239 return 55.508 * x_NaCl / (1 - x_NaCl);
240 }
241
249 template <class Evaluation>
250 static Evaluation molalityCO2inPureWater_(const Evaluation& temperature, const Evaluation& pg, bool extrapolate = false)
251 {
252 OPM_TIMEFUNCTION_LOCAL();
253 const Evaluation& A = computeA_(temperature, pg, extrapolate); // according to Spycher, Pruess and Ennis-King (2003)
254 const Evaluation& B = computeB_(temperature, pg, extrapolate); // according to Spycher, Pruess and Ennis-King (2003)
255 const Evaluation& yH2OinGas = (1 - B) / (1. / A - B); // equilibrium mol fraction of H2O in the gas phase
256 const Evaluation& xCO2inWater = B * (1 - yH2OinGas); // equilibrium mol fraction of CO2 in the water phase
257 return (xCO2inWater * 55.508) / (1 - xCO2inWater); // CO2 molality
258 }
259
269 template <class Evaluation>
270 static Evaluation activityCoefficient_(const Evaluation& temperature,
271 const Evaluation& pg,
272 const Evaluation& molalityNaCl)
273 {
274 OPM_TIMEFUNCTION_LOCAL();
275 const Evaluation& lambda = computeLambda_(temperature, pg); // lambda_{CO2-Na+}
276 const Evaluation& xi = computeXi_(temperature, pg); // Xi_{CO2-Na+-Cl-}
277 const Evaluation& lnGammaStar =
278 2*molalityNaCl*lambda + xi*molalityNaCl*molalityNaCl;
279 return exp(lnGammaStar);
280 }
281
290 template <class Evaluation>
291 static Evaluation computeA_(const Evaluation& temperature, const Evaluation& pg, bool extrapolate = false)
292 {
293 OPM_TIMEFUNCTION_LOCAL();
294 const Evaluation& deltaP = pg / 1e5 - 1; // pressure range [bar] from p0 = 1bar to pg[bar]
295 Scalar v_av_H2O = 18.1; // average partial molar volume of H2O [cm^3/mol]
296 Scalar R = IdealGas::R * 10;
297 const Evaluation& k0_H2O = equilibriumConstantH2O_(temperature); // equilibrium constant for H2O at 1 bar
298 const Evaluation& phi_H2O = fugacityCoefficientH2O(temperature, pg, extrapolate); // fugacity coefficient of H2O for the water-CO2 system
299 const Evaluation& pg_bar = pg / 1.e5;
300 return k0_H2O/(phi_H2O*pg_bar)*exp(deltaP*v_av_H2O/(R*temperature));
301 }
302
311 template <class Evaluation>
312 static Evaluation computeB_(const Evaluation& temperature, const Evaluation& pg, bool extrapolate = false)
313 {
314 OPM_TIMEFUNCTION_LOCAL();
315 const Evaluation& deltaP = pg / 1e5 - 1; // pressure range [bar] from p0 = 1bar to pg[bar]
316 const Scalar v_av_CO2 = 32.6; // average partial molar volume of CO2 [cm^3/mol]
317 const Scalar R = IdealGas::R * 10;
318 const Evaluation& k0_CO2 = equilibriumConstantCO2_(temperature); // equilibrium constant for CO2 at 1 bar
319 const Evaluation& phi_CO2 = fugacityCoefficientCO2(temperature, pg, extrapolate); // fugacity coefficient of CO2 for the water-CO2 system
320 const Evaluation& pg_bar = pg / 1.e5;
321 return phi_CO2*pg_bar/(55.508*k0_CO2)*exp(-(deltaP*v_av_CO2)/(R*temperature));
322 }
323
331 template <class Evaluation>
332 static Evaluation computeLambda_(const Evaluation& temperature, const Evaluation& pg)
333 {
334 OPM_TIMEFUNCTION_LOCAL();
335 static const Scalar c[6] =
336 { -0.411370585, 6.07632013E-4, 97.5347708, -0.0237622469, 0.0170656236, 1.41335834E-5 };
337
338 Evaluation pg_bar = pg / 1.0E5; /* conversion from Pa to bar */
339 return
340 c[0]
341 + c[1]*temperature
342 + c[2]/temperature
343 + c[3]*pg_bar/temperature
344 + c[4]*pg_bar/(630.0 - temperature)
345 + c[5]*temperature*log(pg_bar);
346 }
347
355 template <class Evaluation>
356 static Evaluation computeXi_(const Evaluation& temperature, const Evaluation& pg)
357 {
358 OPM_TIMEFUNCTION_LOCAL();
359 static const Scalar c[4] =
360 { 3.36389723E-4, -1.98298980E-5, 2.12220830E-3, -5.24873303E-3 };
361
362 Evaluation pg_bar = pg / 1.0E5; /* conversion from Pa to bar */
363 return c[0] + c[1]*temperature + c[2]*pg_bar/temperature + c[3]*pg_bar/(630.0 - temperature);
364 }
365
372 template <class Evaluation>
373 static Evaluation equilibriumConstantCO2_(const Evaluation& temperature)
374 {
375 OPM_TIMEFUNCTION_LOCAL();
376 Evaluation temperatureCelcius = temperature - 273.15;
377 static const Scalar c[3] = { 1.189, 1.304e-2, -5.446e-5 };
378 Evaluation logk0_CO2 = c[0] + temperatureCelcius*(c[1] + temperatureCelcius*c[2]);
379 Evaluation k0_CO2 = pow(10.0, logk0_CO2);
380 return k0_CO2;
381 }
382
389 template <class Evaluation>
390 static Evaluation equilibriumConstantH2O_(const Evaluation& temperature)
391 {
392 OPM_TIMEFUNCTION_LOCAL();
393 Evaluation temperatureCelcius = temperature - 273.15;
394 static const Scalar c[4] = { -2.209, 3.097e-2, -1.098e-4, 2.048e-7 };
395 Evaluation logk0_H2O =
396 c[0] + temperatureCelcius*(c[1] + temperatureCelcius*(c[2] + temperatureCelcius*c[3]));
397 return pow(10.0, logk0_H2O);
398 }
399
400};
401
402} // namespace BinaryCoeff
403} // namespace Opm
404
405#endif
Relations valid for an ideal gas.
Some templates to wrap the valgrind client request macros.
Binary coefficients for brine and CO2.
Definition Brine_CO2.hpp:42
static Evaluation fugacityCoefficientH2O(const Evaluation &temperature, const Evaluation &pg, bool extrapolate=false)
Returns the fugacity coefficient of the H2O component in a water-CO2 mixture.
Definition Brine_CO2.hpp:190
static Evaluation gasDiffCoeff(const Evaluation &temperature, const Evaluation &pressure, bool extrapolate=false)
Binary diffusion coefficent [m^2/s] of water in the CO2 phase.
Definition Brine_CO2.hpp:56
static Evaluation fugacityCoefficientCO2(const Evaluation &temperature, const Evaluation &pg, bool extrapolate=false)
Returns the fugacity coefficient of the CO2 component in a water-CO2 mixture.
Definition Brine_CO2.hpp:152
static void calculateMoleFractions(const Evaluation &temperature, const Evaluation &pg, const Evaluation &salinity, const int knownPhaseIdx, Evaluation &xlCO2, Evaluation &ygH2O, bool extrapolate=false)
Returns the mol (!) fraction of CO2 in the liquid phase and the mol_ (!) fraction of H2O in the gas p...
Definition Brine_CO2.hpp:97
static Evaluation liquidDiffCoeff(const Evaluation &, const Evaluation &)
Binary diffusion coefficent [m^2/s] of CO2 in the brine phase.
Definition Brine_CO2.hpp:73
static Evaluation henry(const Evaluation &temperature, bool extrapolate=false)
Henry coefficent for CO2 in brine.
Definition Brine_CO2.hpp:140
static Evaluation gasViscosity(Evaluation temperature, const Evaluation &pressure, bool extrapolate=false)
The dynamic viscosity [Pa s] of CO2.
Definition CO2.hpp:207
static Scalar molarMass()
The mass in [kg] of one mole of CO2.
Definition CO2.hpp:70
static Evaluation gasDensity(const Evaluation &temperature, const Evaluation &pressure, bool extrapolate=false)
The density of CO2 at a given pressure and temperature [kg/m^3].
Definition CO2.hpp:193
static const Scalar molarMass()
The molar mass in of water.
Definition H2O.hpp:83
Relations valid for an ideal gas.
Definition IdealGas.hpp:38
static const Scalar R
The ideal gas constant .
Definition IdealGas.hpp:41
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30