Line data Source code
1 : //* This file is part of the MOOSE framework
2 : //* https://mooseframework.inl.gov
3 : //*
4 : //* All rights reserved, see COPYRIGHT for full restrictions
5 : //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 : //*
7 : //* Licensed under LGPL 2.1, please see LICENSE for details
8 : //* https://www.gnu.org/licenses/lgpl-2.1.html
9 :
10 : #pragma once
11 :
12 : #include "PorousFlowFluidStateMultiComponentBase.h"
13 :
14 : class BrineFluidProperties;
15 : class SinglePhaseFluidProperties;
16 : class Water97FluidProperties;
17 :
18 : /**
19 : * Specialized class for brine and CO2 including calculation of mutual
20 : * solubility of the two fluids using a high-accuracy fugacity-based formulation.
21 : *
22 : * For temperatures 12C <= T <= 99C, the formulation is based on
23 : * Spycher, Pruess and Ennis-King, CO2-H2O mixtures in the geological
24 : * sequestration of CO2. I. Assessment and calculation of mutual solubilities
25 : * from 12 to 100C and up to 600 bar, Geochimica et Cosmochimica Acta, 67, 3015-3031 (2003),
26 : * and Spycher and Pruess, CO2-H2O mixtures in the geological sequestration of CO2. II.
27 : * Partitioning in chloride brine at 12-100C and up to 600 bar, Geochimica et
28 : * Cosmochimica Acta, 69, 3309-3320 (2005).
29 : *
30 : * For temperatures 109C <= T <= 300C, the formulation is based on
31 : * Spycher and Pruess, A Phase-Partitioning Model for CO2-Brine Mixtures at Elevated
32 : * Temperatures and Pressures: Application to CO2-Enhanced Geothermal Systems,
33 : * Transport in Porous Media, 82, 173-196 (2010)
34 : *
35 : * As the two formulations do not coincide at temperatures near 100C, a cubic
36 : * polynomial is used in the intermediate temperature range 99C < T < 109C to
37 : * provide a smooth transition from the two formulations in this region.
38 : *
39 : * Notation convention
40 : * Throughout this class, both mole fractions and mass fractions will be used.
41 : * The following notation will be used:
42 : * yk: mole fraction of component k in the gas phase
43 : * xk: mole fraction of component k in the liquid phase
44 : * Yk: mass fraction of component k in the gas phase
45 : * Xk: mass fraction of component k in the liquid phase
46 : */
47 : class PorousFlowBrineCO2 : public PorousFlowFluidStateMultiComponentBase
48 : {
49 : public:
50 : static InputParameters validParams();
51 :
52 : PorousFlowBrineCO2(const InputParameters & parameters);
53 :
54 : virtual std::string fluidStateName() const override;
55 :
56 : void thermophysicalProperties(Real pressure,
57 : Real temperature,
58 : Real Xnacl,
59 : Real Z,
60 : unsigned int qp,
61 : std::vector<FluidStateProperties> & fsp) const override;
62 :
63 : void thermophysicalProperties(const ADReal & pressure,
64 : const ADReal & temperature,
65 : const ADReal & Xnacl,
66 : const ADReal & Z,
67 : unsigned int qp,
68 : std::vector<FluidStateProperties> & fsp) const override;
69 :
70 : /**
71 : * Mole fractions of CO2 in brine and water vapor in CO2 at equilibrium
72 : *
73 : * In the low temperature regime (T <= 99C), the mole fractions are calculated directly,
74 : * while in the elevated temperature regime (T >= 109C), they are calculated iteratively.
75 : * In the intermediate regime (99C < T < 109C), a cubic polynomial is used to smoothly
76 : * connect the low and elevated temperature regimes.
77 : *
78 : * @param pressure phase pressure (Pa)
79 : * @param temperature phase temperature (K)
80 : * @param Xnacl NaCl mass fraction (kg/kg)
81 : * @param[out] xco2 mole fraction of CO2 in liquid
82 : * @param[out] yh2o mole fraction of H2O in gas
83 : */
84 : void equilibriumMoleFractions(const ADReal & pressure,
85 : const ADReal & temperature,
86 : const ADReal & Xnacl,
87 : ADReal & xco2,
88 : ADReal & yh2o) const;
89 :
90 : /**
91 : * Mass fractions of CO2 in brine and water vapor in CO2 at equilibrium
92 : *
93 : * @param pressure phase pressure (Pa)
94 : * @param temperature phase temperature (K)
95 : * @param Xnacl NaCl mass fraction (kg/kg)
96 : * @param[out] Xco2 mass fraction of CO2 in liquid (kg/kg)
97 : * @param[out] Yh2o mass fraction of H2O in gas (kg/kg)
98 : */
99 : void equilibriumMassFractions(const ADReal & pressure,
100 : const ADReal & temperature,
101 : const ADReal & Xnacl,
102 : ADReal & Xco2,
103 : ADReal & Yh2o) const;
104 :
105 : /**
106 : * Mass fractions of CO2 and H2O in both phases, as well as derivatives wrt
107 : * PorousFlow variables. Values depend on the phase state (liquid, gas or two phase)
108 : *
109 : * @param pressure phase pressure (Pa)
110 : * @param temperature phase temperature (K)
111 : * @param Xnacl NaCl mass fraction (kg/kg)
112 : * @param Z total mass fraction of CO2 component
113 : * @param[out] PhaseStateEnum current phase state
114 : * @param[out] FluidStateProperties data structure
115 : */
116 : void massFractions(const ADReal & pressure,
117 : const ADReal & temperature,
118 : const ADReal & Xnacl,
119 : const ADReal & Z,
120 : FluidStatePhaseEnum & phase_state,
121 : std::vector<FluidStateProperties> & fsp) const;
122 :
123 : /**
124 : * Thermophysical properties of the gaseous state
125 : *
126 : * @param pressure gas pressure (Pa)
127 : * @param temperature temperature (K)
128 : * @param[out] FluidStateProperties data structure
129 : */
130 : void gasProperties(const ADReal & pressure,
131 : const ADReal & temperature,
132 : std::vector<FluidStateProperties> & fsp) const;
133 :
134 : /**
135 : * Thermophysical properties of the liquid state
136 : *
137 : * @param pressure liquid pressure (Pa)
138 : * @param temperature temperature (K)
139 : * @param Xnacl NaCl mass fraction (kg/kg)
140 : * @param[out] FluidStateProperties data structure
141 : */
142 : void liquidProperties(const ADReal & pressure,
143 : const ADReal & temperature,
144 : const ADReal & Xnacl,
145 : std::vector<FluidStateProperties> & fsp) const;
146 :
147 : /**
148 : * Gas saturation in the two-phase region
149 : *
150 : * @param pressure gas pressure (Pa)
151 : * @param temperature phase temperature (K)
152 : * @param Xnacl NaCl mass fraction (kg/kg)
153 : * @param Z total mass fraction of CO2 component
154 : * @param FluidStateProperties data structure
155 : * @return gas saturation (-)
156 : */
157 : ADReal saturation(const ADReal & pressure,
158 : const ADReal & temperature,
159 : const ADReal & Xnacl,
160 : const ADReal & Z,
161 : std::vector<FluidStateProperties> & fsp) const;
162 :
163 : /**
164 : * Gas and liquid properties in the two-phase region
165 : *
166 : * @param pressure gas pressure (Pa)
167 : * @param temperature phase temperature (K)
168 : * @param Xnacl NaCl mass fraction (kg/kg)
169 : * @param Z total mass fraction of NCG component
170 : * @param qp quadpoint for capillary presssure
171 : * @param[out] FluidStateProperties data structure
172 : */
173 : void twoPhaseProperties(const ADReal & pressure,
174 : const ADReal & temperature,
175 : const ADReal & Xnacl,
176 : const ADReal & Z,
177 : unsigned int qp,
178 : std::vector<FluidStateProperties> & fsp) const;
179 :
180 : /**
181 : * Fugacity coefficients for H2O and CO2 for T <= 100C
182 : * Eq. (B7) from Spycher, Pruess and Ennis-King (2003)
183 : *
184 : * @param pressure gas pressure (Pa)
185 : * @param temperature temperature (K)
186 : * @param co2_density CO2 density (kg/m^3)
187 : * @param[out] fco2 fugacity coefficient for CO2
188 : * @param[out] fh2o fugacity coefficient for H2O
189 : */
190 : void fugacityCoefficientsLowTemp(const ADReal & pressure,
191 : const ADReal & temperature,
192 : const ADReal & co2_density,
193 : ADReal & fco2,
194 : ADReal & fh2o) const;
195 :
196 : ///@{
197 : /**
198 : * Fugacity coefficients for H2O and CO2 at elevated temperatures (100C < T <= 300C).
199 : * Eq. (A8) from Spycher & Pruess (2010)
200 : *
201 : * @param pressure gas pressure (Pa)
202 : * @param temperature temperature (K)
203 : * @param co2_density CO2 density (kg/m^3)
204 : * @param xco2 mole fraction of CO2 in liquid phase (-)
205 : * @param yh2o mole fraction of H2O in gas phase (-)
206 : * @param[out] fco2 fugacity coefficient for CO2
207 : * @param[out] fh2o fugacity coefficient for H2O
208 : */
209 : void fugacityCoefficientsHighTemp(const ADReal & pressure,
210 : const ADReal & temperature,
211 : const ADReal & co2_density,
212 : const ADReal & xco2,
213 : const ADReal & yh2o,
214 : ADReal & fco2,
215 : ADReal & fh2o) const;
216 :
217 : ADReal fugacityCoefficientH2OHighTemp(const ADReal & pressure,
218 : const ADReal & temperature,
219 : const ADReal & co2_density,
220 : const ADReal & xco2,
221 : const ADReal & yh2o) const;
222 :
223 : ADReal fugacityCoefficientCO2HighTemp(const ADReal & pressure,
224 : const ADReal & temperature,
225 : const ADReal & co2_density,
226 : const ADReal & xco2,
227 : const ADReal & yh2o) const;
228 : ///@}
229 :
230 : /**
231 : * Activity coefficient of H2O
232 : * Eq. (12) from Spycher & Pruess (2010)
233 : *
234 : * @param temperature temperature (K)
235 : * @param xco2 mole fraction of CO2 in liquid phase (-)
236 : * @return activity coefficient
237 : */
238 : ADReal activityCoefficientH2O(const ADReal & temperature, const ADReal & xco2) const;
239 :
240 : /**
241 : * Activity coefficient of CO2
242 : * Eq. (13) from Spycher & Pruess (2010)
243 : *
244 : * @param temperature temperature (K)
245 : * @param xco2 mole fraction of CO2 in liquid phase (-)
246 : * @return activity coefficient
247 : */
248 : ADReal activityCoefficientCO2(const ADReal & temperature, const ADReal & xco2) const;
249 :
250 : /**
251 : * Activity coefficient for CO2 in brine. From Duan and Sun, An improved model calculating
252 : * CO2 solubility in pure water and aqueous NaCl solutions from 257 to 533 K and from 0 to
253 : * 2000 bar, Chem. Geol., 193, 257-271 (2003)
254 : *
255 : * Note: this is not a 'true' activity coefficient, and is instead related to the molality
256 : * of CO2 in water and brine. Nevertheless, Spycher and Pruess (2005) refer to it as an
257 : * activity coefficient, so this notation is followed here.
258 : *
259 : * @param pressure phase pressure (Pa)
260 : * @param temperature phase temperature (K)
261 : * @param Xnacl salt mass fraction (kg/kg)
262 : * @return activity coefficient for CO2 in brine
263 : */
264 : ADReal activityCoefficient(const ADReal & pressure,
265 : const ADReal & temperature,
266 : const ADReal & Xnacl) const;
267 :
268 : /**
269 : * Activity coefficient for CO2 in brine used in the elevated temperature formulation.
270 : * Eq. (18) from Spycher and Pruess (2010).
271 : *
272 : * Note: unlike activityCoefficient(), no pressure dependence is included in
273 : * this formulation
274 : *
275 : * @param temperature phase temperature (K)
276 : * @param Xnacl salt mass fraction (kg/kg)
277 : * @return activity coefficient for CO2 in brine
278 : */
279 : ADReal activityCoefficientHighTemp(const ADReal & temperature, const ADReal & Xnacl) const;
280 :
281 : /**
282 : * Equilibrium constant for H2O
283 : * For temperatures 12C <= T <= 99C, uses Spycher, Pruess and Ennis-King (2003)
284 : * For temperatures 109C <= T <= 300C, uses Spycher and Pruess (2010)
285 : * For temperatures 99C < T < 109C, the value is calculated by smoothly interpolating
286 : * the two formulations
287 : *
288 : * @param temperature temperature (K)
289 : * @return equilibrium constant for H2O
290 : */
291 : ADReal equilibriumConstantH2O(const ADReal & temperature) const;
292 :
293 : /**
294 : * Equilibrium constant for CO2
295 : * For temperatures 12C <= T <= 99C, uses Spycher, Pruess and Ennis-King (2003)
296 : * For temperatures 109C <= T <= 300C, uses Spycher and Pruess (2010)
297 : * For temperatures 99C < T < 109C, the value is calculated by smoothly interpolating
298 : * the two formulations
299 : *
300 : * @param temperature temperature (K)
301 : * @return equilibrium constant for CO2
302 : */
303 : ADReal equilibriumConstantCO2(const ADReal & temperature) const;
304 :
305 : /**
306 : * Partial density of dissolved CO2
307 : * From Garcia, Density of aqueous solutions of CO2, LBNL-49023 (2001)
308 : *
309 : * @param temperature fluid temperature (K)
310 : * @return partial molar density (kg/m^3)
311 : */
312 : ADReal partialDensityCO2(const ADReal & temperature) const;
313 :
314 : virtual Real totalMassFraction(
315 : Real pressure, Real temperature, Real Xnacl, Real saturation, unsigned int qp) const override;
316 :
317 : /**
318 : * The index of the salt component
319 : * @return salt component number
320 : */
321 1 : unsigned int saltComponentIndex() const { return _salt_component; };
322 :
323 : /**
324 : * Henry's constant of dissolution of gas phase CO2 in brine. From
325 : * Battistelli et al, A fluid property module for the TOUGH2 simulator for saline brines
326 : * with non-condensible gas, Proc. Eighteenth Workshop on Geothermal Reservoir Engineering (1993)
327 : *
328 : * @param temperature fluid temperature (K)
329 : * @param Xnacl NaCl mass fraction (kg/kg)
330 : * @return Henry's constant (Pa)
331 : */
332 : ADReal henryConstant(const ADReal & temperature, const ADReal & Xnacl) const;
333 :
334 : /**
335 : * Enthalpy of dissolution of gas phase CO2 in brine calculated using Henry's constant
336 : * From Himmelblau, Partial molal heats and entropies of solution for gases dissolved
337 : * in water from the freezing to the near critical point, J. Phys. Chem. 63 (1959).
338 : * Correction due to salinity from Battistelli et al, A fluid property module for the
339 : * TOUGH2 simulator for saline brines with non-condensible gas, Proc. Eighteenth Workshop
340 : * on Geothermal Reservoir Engineering (1993).
341 : *
342 : * @param temperature fluid temperature (K)
343 : * @param Xnacl NaCl mass fraction (kg/kg)
344 : * @return enthalpy of dissolution (J/kg)
345 : */
346 : ADReal enthalpyOfDissolutionGas(const ADReal & temperature, const ADReal & Xnacl) const;
347 :
348 : /**
349 : * Enthalpy of dissolution of CO2 in brine calculated using linear fit to model of
350 : * Duan and Sun, An improved model calculating CO2 solubility in pure water and aqueous NaCl
351 : * solutions from 273 to 533 K and from 0 to 2000 bar, Chemical geology, 193, 257--271 (2003).
352 : *
353 : * In the region of interest, the more complicated model given in Eq. (8) of Duan and Sun
354 : * is well approximated by a simple linear fit (R^2 = 0.9997).
355 : *
356 : * Note: as the effect of salt mass fraction is small, it is not included in this function.
357 : *
358 : * @param temperature fluid temperature (K)
359 : * @return enthalpy of dissolution (J/kg)
360 : */
361 : ADReal enthalpyOfDissolution(const ADReal & temperature) const;
362 :
363 : /**
364 : * Mole fractions of CO2 in brine and water vapor in CO2 at equilibrium in the low
365 : * temperature regime (T <= 99C).
366 : *
367 : * @param pressure phase pressure (Pa)
368 : * @param temperature phase temperature (K)
369 : * @param Xnacl NaCl mass fraction (kg/kg)
370 : * @param[out] xco2 mole fraction of CO2 in liquid
371 : * @param[out] yh2o mass fraction of mole in gas
372 : */
373 : void equilibriumMoleFractionsLowTemp(const ADReal & pressure,
374 : const ADReal & temperature,
375 : const ADReal & Xnacl,
376 : ADReal & xco2,
377 : ADReal & yh2o) const;
378 :
379 : /**
380 : * Function to solve for yh2o and xco2 iteratively in the elevated temperature regime (T > 100C)
381 : *
382 : * @param pressure gas pressure (Pa)
383 : * @param temperature fluid temperature (K)
384 : * @param Xnacl NaCl mass fraction (kg/kg)
385 : * @param co2_density CO2 density (kg/m^3)
386 : * @param[out] xco2 mole fraction of CO2 in liquid phase (-)
387 : * @param[out] yh2o mole fraction of H2O in gas phase (-)
388 : */
389 : void solveEquilibriumMoleFractionHighTemp(Real pressure,
390 : Real temperature,
391 : Real Xnacl,
392 : Real co2_density,
393 : Real & xco2,
394 : Real & yh2o) const;
395 :
396 : protected:
397 : /// Check the input variables
398 : virtual void checkVariables(Real pressure, Real temperature) const;
399 :
400 : /**
401 : * Cubic function to smoothly interpolate between the low temperature and elevated
402 : * temperature models for 99C < T < 109C
403 : *
404 : * @param temperature temperature (K)
405 : * @param f0 function value at T = 372K (99C)
406 : * @param df0 derivative of function at T = 372K (99C)
407 : * @param f1 function value at T = 382K (109C)
408 : * @param df1 derivative of function at T = 382K (109C)
409 : * @param[out] value value at the given temperature
410 : * @param[out] deriv derivative at the given temperature
411 : */
412 : void smoothCubicInterpolation(
413 : Real temperature, Real f0, Real df0, Real f1, Real df1, Real & value, Real & deriv) const;
414 :
415 : ///@{
416 : /**
417 : * The function A (Eq. (11) of Spycher, Pruess and Ennis-King (2003) for T <= 100C,
418 : * and Eqs. (10) and (17) of Spycher and Pruess (2010) for T > 100C)
419 : *
420 : * @param pressure gas pressure (Pa)
421 : * @param temperature fluid temperature (K)
422 : * @param Xnacl NaCl mass fraction (kg/kg)
423 : * @param co2_density CO2 density (kg/m^3)
424 : * @param xco2 mole fraction of CO2 in liquid phase (-)
425 : * @param yh2o mole fraction of H2O in gas phase (-)
426 : * @param[out] A the function A
427 : * @param[out] B the function B
428 : */
429 : void funcABHighTemp(Real pressure,
430 : Real temperature,
431 : Real Xnacl,
432 : Real co2_density,
433 : Real xco2,
434 : Real yh2o,
435 : Real & A,
436 : Real & B) const;
437 :
438 : void funcABHighTemp(Real pressure,
439 : Real temperature,
440 : Real Xnacl,
441 : Real co2_density,
442 : Real xco2,
443 : Real yh2o,
444 : Real & A,
445 : Real & dA_dp,
446 : Real & dA_dT,
447 : Real & B,
448 : Real & dB_dp,
449 : Real & dB_dT,
450 : Real & dB_dX) const;
451 :
452 : void funcABLowTemp(const ADReal & pressure,
453 : const ADReal & temperature,
454 : const ADReal & co2_density,
455 : ADReal & A,
456 : ADReal & B) const;
457 : ///@}
458 :
459 : /// Salt component index
460 : const unsigned int _salt_component;
461 : /// Fluid properties UserObject for water
462 : const BrineFluidProperties & _brine_fp;
463 : /// Fluid properties UserObject for the CO2
464 : const SinglePhaseFluidProperties & _co2_fp;
465 : /// Molar mass of water (kg/mol)
466 : const Real _Mh2o;
467 : /// Inverse of molar mass of H2O (mol/kg)
468 : const Real _invMh2o;
469 : /// Molar mass of CO2 (kg/mol)
470 : const Real _Mco2;
471 : /// Molar mass of NaCL
472 : const Real _Mnacl;
473 : /// Molar gas constant in bar cm^3 /(K mol)
474 : const Real _Rbar;
475 : /// Temperature below which the Spycher, Pruess & Ennis-King (2003) model is used (K)
476 : const Real _Tlower;
477 : /// Temperature above which the Spycher & Pruess (2010) model is used (K)
478 : const Real _Tupper;
479 : /// Minimum Z - below this value all CO2 will be dissolved. This reduces the
480 : /// computational burden when small values of Z are present
481 : const Real _Zmin;
482 : /// Henry's coefficeients for CO2
483 : const std::vector<Real> _co2_henry;
484 : };
|
