| Line | Branch | Exec | Source |
|---|---|---|---|
| 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 | // SPDX-FileCopyrightInfo: Copyright © DuMux Project contributors, see AUTHORS.md in root folder | ||
| 5 | // SPDX-License-Identifier: GPL-3.0-or-later | ||
| 6 | // | ||
| 7 | /*! | ||
| 8 | * \file | ||
| 9 | * \ingroup Binarycoefficients | ||
| 10 | * \brief Binary coefficients for CO2 and brine. | ||
| 11 | */ | ||
| 12 | #ifndef DUMUX_BINARY_COEFF_BRINE_CO2_HH | ||
| 13 | #define DUMUX_BINARY_COEFF_BRINE_CO2_HH | ||
| 14 | |||
| 15 | #include <dune/common/math.hh> | ||
| 16 | |||
| 17 | #include <type_traits> | ||
| 18 | #include <dune/common/debugstream.hh> | ||
| 19 | #include <dumux/common/parameters.hh> | ||
| 20 | #include <dumux/material/components/brine.hh> | ||
| 21 | #include <dumux/material/components/h2o.hh> | ||
| 22 | #include <dumux/material/components/co2.hh> | ||
| 23 | #include <dumux/material/idealgas.hh> | ||
| 24 | |||
| 25 | namespace Dumux::BinaryCoeff { | ||
| 26 | |||
| 27 | /*! | ||
| 28 | * \ingroup Binarycoefficients | ||
| 29 | * \brief Binary coefficients for brine and CO2. | ||
| 30 | */ | ||
| 31 | template<class Scalar, class CO2, bool verbose = true> | ||
| 32 | class Brine_CO2 { | ||
| 33 | using H2O = Components::H2O<Scalar>; | ||
| 34 | using IdealGas = Dumux::IdealGas<Scalar>; | ||
| 35 | static constexpr int lPhaseIdx = 0; // index of the liquid phase | ||
| 36 | static constexpr int gPhaseIdx = 1; // index of the gas phase | ||
| 37 | |||
| 38 | public: | ||
| 39 | /*! | ||
| 40 | * \brief Binary diffusion coefficient \f$\mathrm{[m^2/s]}\f$ of water in the CO2 phase. | ||
| 41 | * | ||
| 42 | * According to B. Xu et al. (2002) \cite xu2003 <BR> | ||
| 43 | * \param temperature the temperature \f$\mathrm{[K]}\f$ | ||
| 44 | * \param pressure the phase pressure \f$\mathrm{[Pa]}\f$ | ||
| 45 | */ | ||
| 46 | 74102554 | static Scalar gasDiffCoeff(Scalar temperature, Scalar pressure) | |
| 47 | { | ||
| 48 |
7/14✓ Branch 0 taken 17 times.
✓ Branch 1 taken 74102537 times.
✓ Branch 3 taken 17 times.
✗ Branch 4 not taken.
✓ Branch 6 taken 17 times.
✗ Branch 7 not taken.
✓ Branch 9 taken 17 times.
✗ Branch 10 not taken.
✓ Branch 12 taken 17 times.
✗ Branch 13 not taken.
✓ Branch 15 taken 17 times.
✗ Branch 16 not taken.
✗ Branch 17 not taken.
✗ Branch 18 not taken.
|
74102605 | static const bool hasGasDiffCoeff = hasParam("BinaryCoefficients.GasDiffCoeff"); |
| 49 |
1/2✓ Branch 0 taken 74102554 times.
✗ Branch 1 not taken.
|
74102554 | if (!hasGasDiffCoeff) //in case one might set that user-specific as e.g. in dumux-lecture/mm/convectivemixing |
| 50 | { | ||
| 51 | //Diffusion coefficient of water in the CO2 phase | ||
| 52 | 74102554 | constexpr Scalar PI = 3.141593; | |
| 53 | 74102554 | constexpr Scalar k = 1.3806504e-23; // Boltzmann constant | |
| 54 | 74102554 | constexpr Scalar c = 4; // slip parameter, can vary between 4 (slip condition) and 6 (stick condition) | |
| 55 | 74102554 | constexpr Scalar R_h = 1.72e-10; // hydrodynamic radius of the solute | |
| 56 | 74102554 | const Scalar mu = CO2::gasViscosity(temperature, pressure); // CO2 viscosity | |
| 57 | 74102554 | return k / (c * PI * R_h) * (temperature / mu); | |
| 58 | } | ||
| 59 | else | ||
| 60 | { | ||
| 61 | ✗ | static const Scalar D = getParam<Scalar>("BinaryCoefficients.GasDiffCoeff"); | |
| 62 | ✗ | return D; | |
| 63 | } | ||
| 64 | } | ||
| 65 | |||
| 66 | /*! | ||
| 67 | * \brief Binary diffusion coefficient \f$\mathrm{[m^2/s]}\f$ of CO2 in the brine phase. | ||
| 68 | * | ||
| 69 | * \param temperature the temperature \f$\mathrm{[K]}\f$ | ||
| 70 | * \param pressure the phase pressure \f$\mathrm{[Pa]}\f$ | ||
| 71 | */ | ||
| 72 | 74102554 | static Scalar liquidDiffCoeff(Scalar temperature, Scalar pressure) | |
| 73 | { | ||
| 74 | //Diffusion coefficient of CO2 in the brine phase | ||
| 75 |
7/14✓ Branch 0 taken 17 times.
✓ Branch 1 taken 74102537 times.
✓ Branch 3 taken 17 times.
✗ Branch 4 not taken.
✓ Branch 6 taken 17 times.
✗ Branch 7 not taken.
✓ Branch 9 taken 17 times.
✗ Branch 10 not taken.
✓ Branch 12 taken 17 times.
✗ Branch 13 not taken.
✓ Branch 15 taken 17 times.
✗ Branch 16 not taken.
✗ Branch 17 not taken.
✗ Branch 18 not taken.
|
74102605 | static const bool hasLiquidDiffCoeff = hasParam("BinaryCoefficients.LiquidDiffCoeff"); |
| 76 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 74102554 times.
|
74102554 | if (!hasLiquidDiffCoeff) //in case one might set that user-specific as e.g. in dumux-lecture/mm/convectivemixing |
| 77 | return 2e-9; | ||
| 78 | else | ||
| 79 | { | ||
| 80 | ✗ | const Scalar D = getParam<Scalar>("BinaryCoefficients.LiquidDiffCoeff"); | |
| 81 | ✗ | return D; | |
| 82 | } | ||
| 83 | } | ||
| 84 | |||
| 85 | /*! | ||
| 86 | * \brief Returns the _mol_ (!) fraction of CO2 in the liquid | ||
| 87 | * phase and the _mol_ (!) fraction of H2O in the gas phase | ||
| 88 | * for a given temperature, pressure, CO2 density and brine | ||
| 89 | * salinity. | ||
| 90 | * | ||
| 91 | * Implemented according to Spycher and Pruess (2005) \cite spycher2005 <BR> | ||
| 92 | * applying the activity coefficient expression of Duan and Sun (2003) \cite duan2003 <BR> | ||
| 93 | * and the correlations for pure water given in Spycher, Pruess and Ennis-King (2003) \cite spycher2003 <BR> | ||
| 94 | * | ||
| 95 | * \param temperature the temperature \f$\mathrm{[K]}\f$ | ||
| 96 | * \param pg the gas phase pressure \f$\mathrm{[Pa]}\f$ | ||
| 97 | * \param salinity the salinity \f$\mathrm{[kg \ NaCl / kg \ solution]}\f$ | ||
| 98 | * \param knownPhaseIdx indicates which phases are present | ||
| 99 | * \param xlCO2 mole fraction of CO2 in brine \f$\mathrm{[mol/mol]}\f$ | ||
| 100 | * \param ygH2O mole fraction of water in the gas phase \f$\mathrm{[mol/mol]}\f$ | ||
| 101 | */ | ||
| 102 | 104176086 | static void calculateMoleFractions(const Scalar temperature, | |
| 103 | const Scalar pg, | ||
| 104 | const Scalar salinity, | ||
| 105 | const int knownPhaseIdx, | ||
| 106 | Scalar &xlCO2, | ||
| 107 | Scalar &ygH2O) | ||
| 108 | { | ||
| 109 | |||
| 110 | 104176086 | const Scalar A = computeA_(temperature, pg); | |
| 111 | |||
| 112 | /* salinity: conversion from mass fraction to mol fraction */ | ||
| 113 | 208352172 | const Scalar x_NaCl = salinityToMoleFrac_(salinity); | |
| 114 | |||
| 115 | // if both phases are present the mole fractions in each phase can be calculate | ||
| 116 | // with the mutual solubility function | ||
| 117 | 104176086 | if (knownPhaseIdx < 0) | |
| 118 | { | ||
| 119 | 208352172 | const Scalar molalityNaCl = molFracToMolality_(x_NaCl); // molality of NaCl //CHANGED | |
| 120 | 208352172 | const Scalar m0_CO2 = molalityCO2inPureWater_(temperature, pg); // molality of CO2 in pure water | |
| 121 | 104176086 | const Scalar gammaStar = activityCoefficient_(temperature, pg, molalityNaCl);// activity coefficient of CO2 in brine | |
| 122 | 104176086 | const Scalar m_CO2 = m0_CO2 / gammaStar; // molality of CO2 in brine | |
| 123 | 104176086 | xlCO2 = m_CO2 / (molalityNaCl + 55.508 + m_CO2); // mole fraction of CO2 in brine | |
| 124 | 104176086 | ygH2O = A * (1 - xlCO2 - x_NaCl); // mole fraction of water in the gas phase | |
| 125 | } | ||
| 126 | |||
| 127 | // if only liquid phase is present the mole fraction of CO2 in brine is given and | ||
| 128 | // and the virtual equilibrium mole fraction of water in the non-existing gas phase can be estimated | ||
| 129 | // with the mutual solubility function | ||
| 130 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 104176086 times.
|
104176086 | if (knownPhaseIdx == lPhaseIdx) |
| 131 | ✗ | ygH2O = A * (1 - xlCO2 - x_NaCl); | |
| 132 | |||
| 133 | // if only gas phase is present the mole fraction of water in the gas phase is given and | ||
| 134 | // and the virtual equilibrium mole fraction of CO2 in the non-existing liquid phase can be estimated | ||
| 135 | // with the mutual solubility function | ||
| 136 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 104176086 times.
|
104176086 | if (knownPhaseIdx == gPhaseIdx) |
| 137 | ✗ | xlCO2 = 1 - x_NaCl - ygH2O / A; | |
| 138 | 104176086 | } | |
| 139 | |||
| 140 | /*! | ||
| 141 | * \brief Returns the fugacity coefficient of the CO2 component in a water-CO2 mixture | ||
| 142 | * (given in Spycher, Pruess and Ennis-King (2003) \cite spycher2003 ) | ||
| 143 | * | ||
| 144 | * \param T the temperature \f$\mathrm{[K]}\f$ | ||
| 145 | * \param pg the gas phase pressure \f$\mathrm{[Pa]}\f$ | ||
| 146 | */ | ||
| 147 | 104176086 | static Scalar fugacityCoefficientCO2(Scalar T, Scalar pg) | |
| 148 | { | ||
| 149 | 127538034 | const Scalar V = 1 / (CO2::gasDensity(T, pg) / CO2::molarMass()) * 1.e6; // molar volume in cm^3/mol | |
| 150 | 104176086 | const Scalar pg_bar = pg / 1.e5; // gas phase pressure in bar | |
| 151 | 104176086 | const Scalar a_CO2 = (7.54e7 - 4.13e4 * T); // mixture parameter of Redlich-Kwong equation | |
| 152 | 104176086 | constexpr Scalar b_CO2 = 27.8; // mixture parameter of Redlich-Kwong equation | |
| 153 | 104176086 | constexpr Scalar R = IdealGas::R * 10.; // ideal gas constant with unit bar cm^3 /(K mol) | |
| 154 | |||
| 155 | using std::log; | ||
| 156 | using std::exp; | ||
| 157 | using std::pow; | ||
| 158 | 416704344 | const Scalar lnPhiCO2 = log(V / (V - b_CO2)) + b_CO2 / (V - b_CO2) - 2 * a_CO2 / (R | |
| 159 | 416704344 | * pow(T, 1.5) * b_CO2) * log((V + b_CO2) / V) + a_CO2 * b_CO2 | |
| 160 | 208352172 | / (R * pow(T, 1.5) * b_CO2 * b_CO2) * (log((V + b_CO2) / V) | |
| 161 | 208352172 | - b_CO2 / (V + b_CO2)) - log(pg_bar * V / (R * T)); | |
| 162 | |||
| 163 | 104176086 | return exp(lnPhiCO2); // fugacity coefficient of CO2 | |
| 164 | } | ||
| 165 | |||
| 166 | /*! | ||
| 167 | * \brief Returns the fugacity coefficient of the H2O component in a water-CO2 mixture | ||
| 168 | * (given in Spycher, Pruess and Ennis-King (2003) \cite spycher2003 ) | ||
| 169 | * | ||
| 170 | * \param T the temperature \f$\mathrm{[K]}\f$ | ||
| 171 | * \param pg the gas phase pressure \f$\mathrm{[Pa]}\f$ | ||
| 172 | */ | ||
| 173 | 208352172 | static Scalar fugacityCoefficientH2O(Scalar T, Scalar pg) | |
| 174 | { | ||
| 175 | 255076068 | const Scalar V = 1 / (CO2::gasDensity(T, pg) / CO2::molarMass()) * 1.e6; // molar volume in cm^3/mol | |
| 176 | 208352172 | const Scalar pg_bar = pg / 1.e5; // gas phase pressure in bar | |
| 177 | 208352172 | const Scalar a_CO2 = (7.54e7 - 4.13e4 * T);// mixture parameter of Redlich-Kwong equation | |
| 178 | 208352172 | constexpr Scalar a_CO2_H2O = 7.89e7;// mixture parameter of Redlich-Kwong equation | |
| 179 | 208352172 | constexpr Scalar b_CO2 = 27.8;// mixture parameter of Redlich-Kwong equation | |
| 180 | 208352172 | constexpr Scalar b_H2O = 18.18;// mixture parameter of Redlich-Kwong equation | |
| 181 | 208352172 | constexpr Scalar R = IdealGas::R * 10.; // ideal gas constant with unit bar cm^3 /(K mol) | |
| 182 | |||
| 183 | using std::log; | ||
| 184 | using std::pow; | ||
| 185 | using std::exp; | ||
| 186 | 625056516 | const Scalar lnPhiH2O = log(V / (V - b_CO2)) + b_H2O / (V - b_CO2) - 2 * a_CO2_H2O | |
| 187 | 416704344 | / (R * pow(T, 1.5) * b_CO2) * log((V + b_CO2) / V) + a_CO2 | |
| 188 | 416704344 | * b_H2O / (R * pow(T, 1.5) * b_CO2 * b_CO2) * (log((V + b_CO2) | |
| 189 | 416704344 | / V) - b_CO2 / (V + b_CO2)) - log(pg_bar * V / (R * T)); | |
| 190 | |||
| 191 | 208352172 | return exp(lnPhiH2O); // fugacity coefficient of H2O | |
| 192 | } | ||
| 193 | |||
| 194 | private: | ||
| 195 | /*! | ||
| 196 | * \brief Returns the molality of NaCl \f$\mathrm{[mol \ NaCl / kg \ water]}\f$ for a given mole fraction | ||
| 197 | * \param salinity the salinity \f$\mathrm{[kg \ NaCl / kg \ solution]}\f$ | ||
| 198 | */ | ||
| 199 | static Scalar salinityToMoleFrac_(Scalar salinity) | ||
| 200 | { | ||
| 201 | 104176086 | constexpr Scalar Mw = H2O::molarMass(); // molecular weight of water [kg/mol] | |
| 202 | 104176086 | constexpr Scalar Ms = 58.8e-3; // molecular weight of NaCl [kg/mol] | |
| 203 | |||
| 204 | 104176086 | const Scalar X_NaCl = salinity; | |
| 205 | // salinity: conversion from mass fraction to mol fraction | ||
| 206 |
1/2✓ Branch 0 taken 104176086 times.
✗ Branch 1 not taken.
|
104176086 | const Scalar x_NaCl = -Mw * X_NaCl / ((Ms - Mw) * X_NaCl - Ms); |
| 207 | return x_NaCl; | ||
| 208 | } | ||
| 209 | |||
| 210 | /*! | ||
| 211 | * \brief Returns the molality of NaCl \f$\mathrm{(mol \ NaCl / kg \ water)}\f$ | ||
| 212 | * for a given mole fraction \f$\mathrm{(mol \ NaCl / mol\ solution)}\f$ | ||
| 213 | * | ||
| 214 | * \param x_NaCl mole fraction of NaCL in brine \f$\mathrm{[mol/mol]}\f$ | ||
| 215 | */ | ||
| 216 | static Scalar molFracToMolality_(Scalar x_NaCl) | ||
| 217 | { | ||
| 218 | // conversion from mol fraction to molality (dissolved CO2 neglected) | ||
| 219 | 104176086 | return 55.508 * x_NaCl / (1 - x_NaCl); | |
| 220 | } | ||
| 221 | |||
| 222 | /*! | ||
| 223 | * \brief Returns the equilibrium molality of CO2 \f$\mathrm{(mol \ CO2 / kg \ water)}\f$ for a | ||
| 224 | * CO2-water mixture at a given pressure and temperature | ||
| 225 | * | ||
| 226 | * \param temperature the temperature \f$\mathrm{[K]}\f$ | ||
| 227 | * \param pg the gas phase pressure \f$\mathrm{[Pa]}\f$ | ||
| 228 | */ | ||
| 229 | 104176086 | static Scalar molalityCO2inPureWater_(Scalar temperature, Scalar pg) | |
| 230 | { | ||
| 231 | 104176086 | const Scalar A = computeA_(temperature, pg); // according to Spycher, Pruess and Ennis-King (2003) | |
| 232 | 104176086 | const Scalar B = computeB_(temperature, pg); // according to Spycher, Pruess and Ennis-King (2003) | |
| 233 | 104176086 | const Scalar yH2OinGas = (1 - B) / (1. / A - B); // equilibrium mol fraction of H2O in the gas phase | |
| 234 | 104176086 | const Scalar xCO2inWater = B * (1 - yH2OinGas); // equilibrium mol fraction of CO2 in the water phase | |
| 235 | 104176086 | return (xCO2inWater * 55.508) / (1 - xCO2inWater); // CO2 molality | |
| 236 | } | ||
| 237 | |||
| 238 | /*! | ||
| 239 | * \brief Returns the activity coefficient of CO2 in brine for a | ||
| 240 | * molal description. According to Duan and Sun (2003) \cite duan2003 <BR> | ||
| 241 | * given in Spycher and Pruess (2005) \cite spycher2005 <BR> | ||
| 242 | * | ||
| 243 | * \param temperature the temperature \f$\mathrm{[K]}\f$ | ||
| 244 | * \param pg the gas phase pressure \f$\mathrm{[Pa]}\f$ | ||
| 245 | * \param molalityNaCl molality of NaCl \f$\mathrm{(mol \ NaCl / kg \ water)}\f$ | ||
| 246 | */ | ||
| 247 | 104176086 | static Scalar activityCoefficient_(Scalar temperature, Scalar pg, Scalar molalityNaCl) | |
| 248 | { | ||
| 249 | 104176086 | const Scalar lambda = computeLambda_(temperature, pg); // lambda_{CO2-Na+} | |
| 250 | 208352172 | const Scalar xi = computeXi_(temperature, pg); // Xi_{CO2-Na+-Cl-} | |
| 251 | 208352172 | const Scalar lnGammaStar = 2 * lambda * molalityNaCl + xi * molalityNaCl | |
| 252 | 104176086 | * molalityNaCl; | |
| 253 | using std::exp; | ||
| 254 | 104176086 | return exp(lnGammaStar); // molal activity coefficient of CO2 in brine | |
| 255 | } | ||
| 256 | |||
| 257 | /*! | ||
| 258 | * \brief Returns the parameter A for the calculation of | ||
| 259 | * them mutual solubility in the water-CO2 system. | ||
| 260 | * Given in Spycher, Pruess and Ennis-King (2003) \cite spycher2003 <BR> | ||
| 261 | * | ||
| 262 | * \param T the temperature \f$\mathrm{[K]}\f$ | ||
| 263 | * \param pg the gas phase pressure \f$\mathrm{[Pa]}\f$ | ||
| 264 | */ | ||
| 265 | 208352172 | static Scalar computeA_(Scalar T, Scalar pg) | |
| 266 | { | ||
| 267 | 208352172 | const Scalar deltaP = pg / 1e5 - 1; // pressure range [bar] from p0 = 1bar to pg[bar] | |
| 268 | 208352172 | const Scalar v_av_H2O = 18.1; // average partial molar volume of H2O [cm^3/mol] | |
| 269 | 208352172 | constexpr Scalar R = IdealGas::R * 10; | |
| 270 | 208352172 | const Scalar k0_H2O = equilibriumConstantH2O_(T); // equilibrium constant for H2O at 1 bar | |
| 271 | 208352172 | const Scalar phi_H2O = fugacityCoefficientH2O(T, pg); // fugacity coefficient of H2O for the water-CO2 system | |
| 272 | 208352172 | const Scalar pg_bar = pg / 1.e5; | |
| 273 | using std::exp; | ||
| 274 | 208352172 | return k0_H2O / (phi_H2O * pg_bar) * exp(deltaP * v_av_H2O / (R * T)); | |
| 275 | } | ||
| 276 | |||
| 277 | /*! | ||
| 278 | * \brief Returns the parameter B for the calculation of | ||
| 279 | * the mutual solubility in the water-CO2 system. | ||
| 280 | * Given in Spycher, Pruess and Ennis-King (2003) \cite spycher2003 <BR> | ||
| 281 | * | ||
| 282 | * \param T the temperature \f$\mathrm{[K]}\f$ | ||
| 283 | * \param pg the gas phase pressure \f$\mathrm{[Pa]}\f$ | ||
| 284 | */ | ||
| 285 | 104176086 | static Scalar computeB_(Scalar T, Scalar pg) | |
| 286 | { | ||
| 287 | 104176086 | const Scalar deltaP = pg / 1e5 - 1; // pressure range [bar] from p0 = 1bar to pg[bar] | |
| 288 | 104176086 | constexpr Scalar v_av_CO2 = 32.6; // average partial molar volume of CO2 [cm^3/mol] | |
| 289 | 104176086 | constexpr Scalar R = IdealGas::R * 10; | |
| 290 | 208352172 | const Scalar k0_CO2 = equilibriumConstantCO2_(T); // equilibrium constant for CO2 at 1 bar | |
| 291 | 208352172 | const Scalar phi_CO2 = fugacityCoefficientCO2(T, pg); // fugacity coefficient of CO2 for the water-CO2 system | |
| 292 | 104176086 | const Scalar pg_bar = pg / 1.e5; | |
| 293 | using std::exp; | ||
| 294 | 312528258 | return phi_CO2 * pg_bar / (55.508 * k0_CO2) * exp(-(deltaP | |
| 295 | 312528258 | * v_av_CO2) / (R * T)); | |
| 296 | } | ||
| 297 | |||
| 298 | /*! | ||
| 299 | * \brief Returns the parameter lambda, which is needed for the | ||
| 300 | * calculation of the CO2 activity coefficient in the brine-CO2 system. | ||
| 301 | * Given in Spycher and Pruess (2005) \cite spycher2005 <BR> | ||
| 302 | * \param T the temperature \f$\mathrm{[K]}\f$ | ||
| 303 | * \param pg the gas phase pressure \f$\mathrm{[Pa]}\f$ | ||
| 304 | */ | ||
| 305 | 104176086 | static Scalar computeLambda_(Scalar T, Scalar pg) | |
| 306 | { | ||
| 307 | constexpr Scalar c[6] = { -0.411370585, 6.07632013E-4, 97.5347708, | ||
| 308 | -0.0237622469, 0.0170656236, 1.41335834E-5 }; | ||
| 309 | |||
| 310 | using std::log; | ||
| 311 | 104176086 | const Scalar pg_bar = pg / 1.0E5; /* conversion from Pa to bar */ | |
| 312 | 208352172 | return c[0] + c[1] * T + c[2] / T + c[3] * pg_bar / T + c[4] * pg_bar | |
| 313 | 208352172 | / (630.0 - T) + c[5] * T * log(pg_bar); | |
| 314 | } | ||
| 315 | |||
| 316 | /*! | ||
| 317 | * \brief Returns the parameter xi, which is needed for the | ||
| 318 | * calculation of the CO2 activity coefficient in the brine-CO2 system. | ||
| 319 | * Given in Spycher and Pruess (2005) \cite spycer2005 <BR> | ||
| 320 | * \param T the temperature \f$\mathrm{[K]}\f$ | ||
| 321 | * \param pg the gas phase pressure \f$\mathrm{[Pa]}\f$ | ||
| 322 | */ | ||
| 323 | static Scalar computeXi_(Scalar T, Scalar pg) | ||
| 324 | { | ||
| 325 | constexpr Scalar c[4] = { 3.36389723E-4, -1.98298980E-5, | ||
| 326 | 2.12220830E-3, -5.24873303E-3 }; | ||
| 327 | |||
| 328 | 104176086 | Scalar pg_bar = pg / 1.0E5; /* conversion from Pa to bar */ | |
| 329 | 104176086 | return c[0] + c[1] * T + c[2] * pg_bar / T + c[3] * pg_bar / (630.0 - T); | |
| 330 | } | ||
| 331 | |||
| 332 | /*! | ||
| 333 | * \brief Returns the equilibrium constant for CO2, which is needed for the | ||
| 334 | * calculation of the mutual solubility in the water-CO2 system | ||
| 335 | * Given in Spycher, Pruess and Ennis-King (2003) \cite spycher2003 <BR> | ||
| 336 | * \param T the temperature \f$\mathrm{[K]}\f$ | ||
| 337 | */ | ||
| 338 | static Scalar equilibriumConstantCO2_(Scalar T) | ||
| 339 | { | ||
| 340 | 104176086 | const Scalar TinC = T - 273.15; //temperature in °C | |
| 341 | constexpr Scalar c[3] = { 1.189, 1.304e-2, -5.446e-5 }; | ||
| 342 | 104176086 | const Scalar logk0_CO2 = c[0] + c[1] * TinC + c[2] * TinC * TinC; | |
| 343 | using std::pow; | ||
| 344 | 208352172 | return pow(10, logk0_CO2); | |
| 345 | } | ||
| 346 | |||
| 347 | /*! | ||
| 348 | * \brief Returns the equilibrium constant for H2O, which is needed for the | ||
| 349 | * calculation of the mutual solubility in the water-CO2 system | ||
| 350 | * Given in Spycher, Pruess and Ennis-King (2003) \cite spycher2003 <BR> | ||
| 351 | * \param T the temperature \f$\mathrm{[K]}\f$ | ||
| 352 | */ | ||
| 353 | 208352172 | static Scalar equilibriumConstantH2O_(Scalar T) | |
| 354 | { | ||
| 355 | 208352172 | const Scalar TinC = T - 273.15; //temperature in °C | |
| 356 | constexpr Scalar c[4] = { -2.209, 3.097e-2, -1.098e-4, 2.048e-7 }; | ||
| 357 | 416704344 | const Scalar logk0_H2O = c[0] + c[1] * TinC + c[2] * TinC * TinC + c[3] | |
| 358 | 208352172 | * TinC * TinC * TinC; | |
| 359 | using std::pow; | ||
| 360 | 416704344 | return pow(10, logk0_H2O); | |
| 361 | } | ||
| 362 | }; | ||
| 363 | |||
| 364 | /*! | ||
| 365 | * \brief Old version of binary coefficients for CO2 and brine. | ||
| 366 | * Calculates molfraction of CO2 in brine according to Duan and Sun 2003 | ||
| 367 | * molfraction of H2O has been assumed to be a constant value | ||
| 368 | * For use with the actual brine_co2_system this class still needs to be adapted | ||
| 369 | */ | ||
| 370 | template<class Scalar, class CO2Impl, bool verbose = true> | ||
| 371 | class Brine_CO2_Old | ||
| 372 | { | ||
| 373 | using H2O = Components::H2O<Scalar>; | ||
| 374 | using Brine = Components::Brine<Scalar,H2O>; | ||
| 375 | using CO2 = CO2Impl; | ||
| 376 | using IdealGas = Dumux::IdealGas<Scalar>; | ||
| 377 | |||
| 378 | public: | ||
| 379 | /*! | ||
| 380 | * \brief Returns the _mole_ (!) fraction of CO2 in the liquid | ||
| 381 | * phase at a given temperature, pressure and density of | ||
| 382 | * CO2. | ||
| 383 | * \param temperature the temperature \f$\mathrm{[K]}\f$ | ||
| 384 | * \param pg the gas phase pressure \f$\mathrm{[Pa]}\f$ | ||
| 385 | * \param rhoCO2 density of CO2 | ||
| 386 | */ | ||
| 387 | static Scalar moleFracCO2InBrine(Scalar temperature, Scalar pg, Scalar rhoCO2) | ||
| 388 | { | ||
| 389 | // regularisations: | ||
| 390 | if (pg > 2.5e8) { | ||
| 391 | pg = 2.5e8; | ||
| 392 | } | ||
| 393 | if (pg < 2.e5) { | ||
| 394 | pg = 2.e5; | ||
| 395 | } | ||
| 396 | if (temperature < 275.) { | ||
| 397 | temperature = 275; | ||
| 398 | } | ||
| 399 | if (temperature > 600.) { | ||
| 400 | temperature = 600; | ||
| 401 | } | ||
| 402 | |||
| 403 | const Scalar Mw = H2O::molarMass(); /* molecular weight of water [kg/mol] */ | ||
| 404 | const Scalar Ms = 58.8e-3; /* molecular weight of NaCl [kg/mol] */ | ||
| 405 | |||
| 406 | const Scalar X_NaCl = Brine::salinity(); | ||
| 407 | /* salinity: conversion from mass fraction to mole fraction */ | ||
| 408 | const Scalar x_NaCl = -Mw * X_NaCl / ((Ms - Mw) * X_NaCl - Ms); | ||
| 409 | |||
| 410 | // salinity: conversion from mole fraction to molality | ||
| 411 | const Scalar mol_NaCl = -55.56 * x_NaCl / (x_NaCl - 1); | ||
| 412 | |||
| 413 | const Scalar A = computeA_(temperature, pg); /* mu_{CO2}^{l(0)}/RT */ | ||
| 414 | const Scalar B = computeB_(temperature, pg); /* lambda_{CO2-Na+} */ | ||
| 415 | const Scalar C = computeC_(temperature, pg); /* Xi_{CO2-Na+-Cl-} */ | ||
| 416 | const Scalar pgCO2 = partialPressureCO2_(temperature, pg); | ||
| 417 | const Scalar phiCO2 = fugacityCoeffCO2_(temperature, pgCO2, rhoCO2); | ||
| 418 | |||
| 419 | using std::log; | ||
| 420 | using Dune::power; | ||
| 421 | const Scalar exponent = A - log(phiCO2) + 2*B*mol_NaCl + C*power(mol_NaCl,2); | ||
| 422 | |||
| 423 | using std::exp; | ||
| 424 | const Scalar mol_CO2w = pgCO2 / (1e5 * exp(exponent)); /* paper: equation (6) */ | ||
| 425 | |||
| 426 | const Scalar x_CO2w = mol_CO2w / (mol_CO2w + 55.56); /* conversion: molality to mole fraction */ | ||
| 427 | //Scalar X_CO2w = x_CO2w*MCO2/(x_CO2w*MCO2 + (1-x_CO2w)*Mw); /* conversion: mole fraction to mass fraction */ | ||
| 428 | |||
| 429 | return x_CO2w; | ||
| 430 | } | ||
| 431 | |||
| 432 | private: | ||
| 433 | /*! | ||
| 434 | * \brief computation of \f$\mathrm{[mu_{CO2}^{l(0)}/RT]}\f$ | ||
| 435 | * \param T the temperature \f$\mathrm{[K]}\f$ | ||
| 436 | * \param pg the gas phase pressure \f$\mathrm{[Pa]}\f$ | ||
| 437 | */ | ||
| 438 | static Scalar computeA_(Scalar T, Scalar pg) | ||
| 439 | { | ||
| 440 | static const Scalar c[10] = { | ||
| 441 | 28.9447706, | ||
| 442 | -0.0354581768, | ||
| 443 | -4770.67077, | ||
| 444 | 1.02782768E-5, | ||
| 445 | 33.8126098, | ||
| 446 | 9.04037140E-3, | ||
| 447 | -1.14934031E-3, | ||
| 448 | -0.307405726, | ||
| 449 | -0.0907301486, | ||
| 450 | 9.32713393E-4, | ||
| 451 | }; | ||
| 452 | |||
| 453 | const Scalar pg_bar = pg / 1.0E5; /* conversion from Pa to bar */ | ||
| 454 | const Scalar Tr = 630.0 - T; | ||
| 455 | |||
| 456 | using std::log; | ||
| 457 | return | ||
| 458 | c[0] + | ||
| 459 | c[1]*T + | ||
| 460 | c[2]/T + | ||
| 461 | c[3]*T*T + | ||
| 462 | c[4]/Tr + | ||
| 463 | c[5]*pg_bar + | ||
| 464 | c[6]*pg_bar*log(T) + | ||
| 465 | c[7]*pg_bar/T + | ||
| 466 | c[8]*pg_bar/Tr + | ||
| 467 | c[9]*pg_bar*pg_bar/(Tr*Tr); | ||
| 468 | } | ||
| 469 | |||
| 470 | /*! | ||
| 471 | * \brief computation of B | ||
| 472 | * | ||
| 473 | * \param T the temperature \f$\mathrm{[K]}\f$ | ||
| 474 | * \param pg the gas phase pressure \f$\mathrm{[Pa]}\f$ | ||
| 475 | */ | ||
| 476 | static Scalar computeB_(Scalar T, Scalar pg) | ||
| 477 | { | ||
| 478 | const Scalar c1 = -0.411370585; | ||
| 479 | const Scalar c2 = 6.07632013E-4; | ||
| 480 | const Scalar c3 = 97.5347708; | ||
| 481 | const Scalar c8 = -0.0237622469; | ||
| 482 | const Scalar c9 = 0.0170656236; | ||
| 483 | const Scalar c11 = 1.41335834E-5; | ||
| 484 | |||
| 485 | const Scalar pg_bar = pg / 1.0E5; /* conversion from Pa to bar */ | ||
| 486 | |||
| 487 | using std::log; | ||
| 488 | return | ||
| 489 | c1 + | ||
| 490 | c2*T + | ||
| 491 | c3/T + | ||
| 492 | c8*pg_bar/T + | ||
| 493 | c9*pg_bar/(630.0-T) + | ||
| 494 | c11*T*log(pg_bar); | ||
| 495 | } | ||
| 496 | |||
| 497 | /*! | ||
| 498 | * \brief computation of C | ||
| 499 | * | ||
| 500 | * \param T the temperature \f$\mathrm{[K]}\f$ | ||
| 501 | * \param pg the gas phase pressure \f$\mathrm{[Pa]}\f$ | ||
| 502 | */ | ||
| 503 | static Scalar computeC_(Scalar T, Scalar pg) | ||
| 504 | { | ||
| 505 | const Scalar c1 = 3.36389723E-4; | ||
| 506 | const Scalar c2 = -1.98298980E-5; | ||
| 507 | const Scalar c8 = 2.12220830E-3; | ||
| 508 | const Scalar c9 = -5.24873303E-3; | ||
| 509 | |||
| 510 | const Scalar pg_bar = pg / 1.0E5; /* conversion from Pa to bar */ | ||
| 511 | |||
| 512 | return | ||
| 513 | c1 + | ||
| 514 | c2*T + | ||
| 515 | c8*pg_bar/T + | ||
| 516 | c9*pg_bar/(630.0-T); | ||
| 517 | } | ||
| 518 | |||
| 519 | /*! | ||
| 520 | * \brief computation of partial pressure CO2 | ||
| 521 | * | ||
| 522 | * We assume that the partial pressure of brine is its vapor pressure. | ||
| 523 | * \warning: Strictly this is assumption is invalid for CO2 because the | ||
| 524 | * mole fraction of CO2 in brine can be considerable | ||
| 525 | * | ||
| 526 | * \param temperature the temperature \f$\mathrm{[K]}\f$ | ||
| 527 | * \param pg the gas phase pressure \f$\mathrm{[Pa]}\f$ | ||
| 528 | */ | ||
| 529 | static Scalar partialPressureCO2_(Scalar temperature, Scalar pg) | ||
| 530 | { | ||
| 531 | return pg - Brine::vaporPressure(temperature); | ||
| 532 | } | ||
| 533 | |||
| 534 | /*! | ||
| 535 | * \brief The fugacity coefficient of CO2 for a CO2-H2O mixture. | ||
| 536 | * | ||
| 537 | * \param temperature the temperature \f$\mathrm{[K]}\f$ | ||
| 538 | * \param pg the gas phase pressure \f$\mathrm{[Pa]}\f$ | ||
| 539 | * \param rhoCO2 the density of CO2 for the critical volume \f$\mathrm{[kg/m^3]}\f$ | ||
| 540 | */ | ||
| 541 | static Scalar fugacityCoeffCO2_(Scalar temperature, | ||
| 542 | Scalar pg, | ||
| 543 | Scalar rhoCO2) | ||
| 544 | { | ||
| 545 | static const Scalar a[15] = { | ||
| 546 | 8.99288497E-2, | ||
| 547 | -4.94783127E-1, | ||
| 548 | 4.77922245E-2, | ||
| 549 | 1.03808883E-2, | ||
| 550 | -2.82516861E-2, | ||
| 551 | 9.49887563E-2, | ||
| 552 | 5.20600880E-4, | ||
| 553 | -2.93540971E-4, | ||
| 554 | -1.77265112E-3, | ||
| 555 | -2.51101973E-5, | ||
| 556 | 8.93353441E-5, | ||
| 557 | 7.88998563E-5, | ||
| 558 | -1.66727022E-2, | ||
| 559 | 1.3980, | ||
| 560 | 2.96000000E-2 | ||
| 561 | }; | ||
| 562 | |||
| 563 | // reduced temperature | ||
| 564 | const Scalar Tr = temperature / CO2::criticalTemperature(); | ||
| 565 | // reduced pressure | ||
| 566 | const Scalar pr = pg / CO2::criticalPressure(); | ||
| 567 | |||
| 568 | // reduced molar volume. ATTENTION: Vc is _NOT_ the critical | ||
| 569 | // molar volume of CO2. See the reference! | ||
| 570 | const Scalar Vc = IdealGas::R*CO2::criticalTemperature()/CO2::criticalPressure(); | ||
| 571 | const Scalar Vr = | ||
| 572 | // molar volume of CO2 at (temperature, pg) | ||
| 573 | CO2::molarMass() / rhoCO2 | ||
| 574 | * | ||
| 575 | // "pseudo-critical" molar volume | ||
| 576 | 1.0 / Vc; | ||
| 577 | |||
| 578 | // the Z coefficient | ||
| 579 | const Scalar Z = pr * Vr / Tr; | ||
| 580 | |||
| 581 | const Scalar A = a[0] + a[1] / (Tr * Tr) + a[2] / (Tr * Tr * Tr); | ||
| 582 | const Scalar B = a[3] + a[4] / (Tr * Tr) + a[5] / (Tr * Tr * Tr); | ||
| 583 | const Scalar C = a[6] + a[7] / (Tr * Tr) + a[8] / (Tr * Tr * Tr); | ||
| 584 | const Scalar D = a[9] + a[10] / (Tr * Tr) + a[11] / (Tr * Tr * Tr); | ||
| 585 | |||
| 586 | using std::log; | ||
| 587 | using std::exp; | ||
| 588 | const Scalar lnphiCO2 = | ||
| 589 | Z - 1 - | ||
| 590 | log(Z) + | ||
| 591 | A/Vr + | ||
| 592 | B/(2*Vr*Vr) + | ||
| 593 | C/(4*Vr*Vr*Vr*Vr) + | ||
| 594 | D/(5*Vr*Vr*Vr*Vr*Vr) | ||
| 595 | + | ||
| 596 | a[12]/(2*Tr*Tr*Tr*a[14])* | ||
| 597 | ( | ||
| 598 | a[13] + 1 - | ||
| 599 | ( a[13] + 1 + | ||
| 600 | a[14]/(Vr*Vr) | ||
| 601 | )*exp(-a[14]/(Vr*Vr))); | ||
| 602 | |||
| 603 | return exp(lnphiCO2); | ||
| 604 | } | ||
| 605 | |||
| 606 | }; | ||
| 607 | |||
| 608 | } // end namespace Dumux::BinaryCoeff | ||
| 609 | |||
| 610 | #endif | ||
| 611 |