EDP Sciences logo
Submit your paper
Search
Menu
Home All issues Volume 78 (2023) Sci. Tech. Energ. Transition, 78 (2023) 10 Full HTML
Issue
Sci. Tech. Energ. Transition
Volume 78, 2023
Characterization and Modeling of the Subsurface in the Context of Ecological Transition
Article Number 10
Number of page(s) 16
DOI https://doi.org/10.2516/stet/2023004
Published online 12 April 2023

© The Author(s), published by EDP Sciences, 2023

Licence Creative CommonsThis is an Open Access article distributed under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

1 Introduction

In the context of the energy transition, renewable energies like solar and wind are promising solutions to replace fossil fuels [1]. However, such kind of energy production is intermittent. To overcome this issue, these renewable energies can be stored during high-energy production phases, and used when needed [2]. Among all storage technologies, underground storage is well known and provides several advantages such as safety, longevity, and high storage capacity [3]. More specifically, storage in salt caverns is very promising for light gas storage due to very low permeability. The electricity produced from renewables must first be turned into an energy carrier that can be stored, like reactive hydrogen [4], or be used for Compressed Air Energy Storage (CAES) [5].

Salt caverns are developed by injecting water into a geological salt deposit, thus leaching the salt and then withdrawing the brine. Once the desired geometrical volume is obtained, the first gas fill is made by injecting gas and withdrawing the brine. At the end of this so-called “debrining” stage, an amount of brine typically between 1.000 m3 and 10.000 m3 stay at the bottom of the cavern. The cavern is then put into operation. This presence of the residual brine impacts the gas storage operations. First, it can solubilize a part of the injected gas. Second, while it is injected dry, the gas produced contains a higher amount of water, making it necessary to use dehydration surface processes. In the case of compressed air, the withdrawn humid air causes aggressive corrosion conditions for the well and surface piping and equipment. It is thus of primary importance to have a thermodynamic model able to predict phase equilibrium between brine and gas in large conditions of temperature and pressure, and up to saturation in terms of salinity.

Concerning the phase equilibrium model between hydrogen and brines, only a few thermodynamic models have been published. Lopez-Lazaro et al. [6] proposed an approach based on the Soreide and Whitson equation of state [7] calibrated on a few existing experimental and molecular simulation data. Recently, Chabab et al. [8] published new experimental data in temperature and salinity ranges in better accordance with underground storage, and they used both Soreide and Whitson and an electrolyte version of the Cubic Plus Association equation of state [9] to successfully model their data. Li et al. [10] proposed a heterogeneous approach based on the Pitzer activity model and Henry constants for the aqueous phase and a pure-hydrogen equation of state for the vapour phase. This model assumes that the vapour phase behaves ideally, which leads to non-negligible deviations at high pressure, as demonstrated in a previous study [11]. Recently, Zhu et al. [12] also proposed a heterogeneous approach using a Pitzer model for the liquid phase. The fugacity of hydrogen in the vapour phase is calculated with a pure-hydrogen model. A correction is brought to evaluate the water content in this phase. Although accurate on data used for model calibration, these heterogeneous approaches of Li et al. and Zhu et al. are rather devoted to low-pressure conditions and do not allow an accurate calculation of the vapour phase composition. Another heterogeneous approach was recently proposed by Torín-Ollarves and Trusler, based on the Setchenow coefficient and Henry constant correlations [13]. This model fits well with experimental data but remains quite empirical with a large number of parameters. Again, this model does not allow an accurate calculation of the vapour phase composition. Roa Pinto et al. proposed to use an electrolyte version of the PPC-SAFT equation of state and show good agreement with previous models, but also a good capacity in extrapolation when compared to data not used in the calibration database [11]. However, this model exhibits non-negligible deviations in the compressibility factor of pure hydrogen at high pressure, opening the way to improvements to get better accuracy of gas volume produced during storage withdrawal.

Concerning phase equilibrium between air and saturated brines, no specific model has been found in the literature probably due to the lack of experimental data, although molecular simulation is more and more used as a complementary tool to experiments in this field [14]. Nevertheless, we can mention the modeling work recently proposed by Chabab et al. on N2 + brine and O2 + brine systems with the Soreide and Whitson equation of state [15]. This model has been validated on both gas solubility and water content in the vapor phase at a large range of temperatures and salinities. Coupled with their previous parameterization for other gases like H2, CO2, and CH4, this is probably today the most complete thermodynamic model to deal with gas storage in salt caverns.

In this work, we propose a thermodynamic model able to predict both phase equilibrium: H2 + brine and Air + brine. Air is considered a mixture of 80% mol. of nitrogen and 20% mol. of oxygen. Thus, binary systems N2 + brine and O2 + brine must be studied first, before evaluating the model in prediction on the multi-component mixture air + brine. Brine is assumed to be made of water and NaCl since this salt is the most abundant in salt caverns. The model proposed is a revision of the initial model of Roa Pinto et al. [11]. It is based on the e-PPC-SAFT equation of state: due to its strong physical background, this model has a good capability of extrapolation, which is an important requirement in the present case of scarce experimental data. This model will be parameterized to reproduce not only gas solubility in brine but also the composition of the vapor phase. The capability of the model to predict accurately the water content of the vapor phase is a key benefit compared to existing modeling works based on a heterogeneous approach [10, 12, 13]. Furthermore, the accuracy of a SAFT-based thermodynamic model, which includes specific terms to deal with association and electrolyte interactions, is expected to be at least the same as a cubic-based model such as the Soreide and Whitson equation of state, and to have better predictive behavior. Another key benefit of the model proposed in this work is to be based on an ion-parameterization approach, and not a salt-parameterization approach like with the Soreide and Whitson approach. Consequently, such a model allows calculating specifically the liquid fugacity of each ion individually, which can be of prior importance if this model has to be further used in a reactive equilibria scheme where ions are involved, for instance in typical fluid-rock interaction calculations.

In the first part of this article, the e-PPC-SAFT framework will be recalled, with a specific focus on the model parameters. In the second part, a review of available experimental data is proposed. The systems investigated cover pure components (hydrogen, nitrogen, oxygen), brine, and mixtures of light gases with pure water and brine. Then, in the third part, the model parameterization is carried out for both pure components and mixtures. In the fourth part, the model is validated on its capability to predict data and properties not used in the parameterization database, and on its capability to compute multicomponent mixtures such as air only from binary considerations. Finally, an application of this model will be presented to calculate gas moisture during typical Compressed Air Energy Storage injection-withdrawal cycles.

2 The e-PPC-SAFT model

The e-PPC-SAFT equation of state [16] is used to model the different forces acting in an electrolytic multiphasic system. It gives the reduced residual Helmholtz energy A res as a sum of different contributions, each contribution representing a specific interacting force between particles:(1)where R is the ideal gas constant and T the temperature.

The first four contributions (hard-sphere, chain, dispersion, and association) derive directly from the original PC-SAFT version of Gross and Sadowski [17]. The hard sphere term describes the repulsive interactions between hard spheres. Its expression is the following:(2)with(3)where, ρ is the molar density and total volume of the fluid, NA the Avogadro’s number (6.02214 × 1023 mol−1), is the hard sphere diameter, ε i is the dispersive energy parameter, λ is the softness parameter which is equal to 0.12 (except for water in the model used in this work: λ i  = 0.203) and x i is the molar fraction of component i.

The chain term is used to model the chain formation between segments. This contribution is given by:(4)where m i is the chain length parameter and is the radial distribution function given by Boublik [18] and Mansoori et al. [19]:(5)and(6)

The dispersive contribution describes the short-range Van der Waals attraction between segments. The corresponding expression is:(7)with(8) (9)where is the reduced radial distance around a segment and Z hc is the compressibility factor of hard chain. Mixing parameters and are given by:(10) (11)where ε ij and are evaluated using the Lorentz–Berthelot mixing rules:(12)where k ij is an empirical interaction binary parameter.

The association term describes the strong short-range attractions due to association forces between association sites. It derives from the Wertheim theory [20] and is given by [21]:(13)with(14)and(15)

is the fraction of association sites non-bonded to site A of molecule i. M i is the number of association sites in molecule i, and the association force between site A of molecule i and site B of molecule j. It is characterized by two binary parameters: the association energy () and is the association volume (), determined from the CR-1 mixing rule:(16)where w ij and u ij are empirical interaction binary parameters.

From this original version of the PC-SAFT equation of state, Nguyen-Huynh et al. [22] added a specific polar term from the Jog and Chapman theory [23] to include interactions between dipoles and quadrupoles. The expression of this term is the following:(17) (18) (19)where μ i is the dipole moment of the molecule and x pi its polar fraction. and are correlations whose expressions can be found in [23]. The quadripolar moment Q i of the molecule and its quadrupolar fraction can also be taken into account, following equations given in [23].

The addition of a Non-Additive Hard-Sphere (NAHS) term to the PPC-SAFT equation of state was proposed by Trinh et al. to improve prediction with hydrogen-containing mixtures [24] by correcting the additivity assumption of the hard spheres. The expression of this term is:(20)where is the radial distribution function of square-well potential [25] and d ij the cross-diameter calculated by:(21)where l ij is an empirical binary interaction parameter.

The last two terms of equation (1) are specific to electrolytic systems. For long-range electrolytic interactions, the use of MSA (Mean Spherical Approximation) is used as suggested by many authors [11, 16, 26–32], although an approach based on Debye–Hückel theory is expected to give a similar performance [33]. Its expression is given by:(22)where ε r is the dielectric constant of the solution, ε o the dielectric constant in the vacuum (8.85418 × 10−12 C mol−1), n i the number of moles of ion i, e is the electron charge (1.60218 × 10−19 C), Z i the ion charge, the solvated ion diameter and Γ the shielding parameter calculated as follows:(23)

In equation (1), the Born term is also classically used to model short-range ion solvation forces [11, 16, 27]. It is defined by:(24)where is the Born diameter for each ion. This diameter is determined from experimental Gibbs energies of solvation [34], through the equation resulting from equation (25):(25)

Both MSA and Born terms involve the dielectric constant of the aqueous solution. In this work, it is calculated using the model of the dielectric constant of pure water proposed by Schreckenberg et al. [27] to which the salt correction of Pottel is added [35]:(26)with(27)where N T is the total number of moles, V is the volume of the system, d v and d T are constants with a value for the water of 0.3777 dm3/mol and 1403 K respectively, and(28)

Finally, Table 1 provides a summary of what parameters are used for each term of the e-PPC-SAFT model.

Table 1

e-PPC-SAFT pure-component and binary parameters.

Knowing the residual Helmholtz energy from equation (1), the fugacity coefficient of each component i in a given phase p can be calculated from:(29)where Vp is the volume of the phase p, is the number of mole of component i, P is the pressure and n is the total mole number in phase p.

For all component i, phase equilibrium is calculated by equalizing fugacities of i in liquid and vapor phases, under mass balance constraint:(30)where x and y stand for the liquid and vapor molar fraction, respectively, and and the vectors of compositions of these phases.

For light gases showing a very low solubility in the aqueous phase, it is common to express solubility in terms of the Henry constant. It is defined as:(31)where fi is the fugacity of i in the liquid or vapor phase. From an equation of state, it can be calculated by:(32)where is the saturation pressure of the solvent, and the superscript stands for infinite dilution of i in the solvent. It can be noticed that Henry’s constants computed with equation (32) are for real mixtures and not salt-free mixtures.

The computation of the Henry constant is carried out in two steps. First, a “bubble point” flash of the real solvent (i.e. with mole fraction including ions) is carried out. Then, the fugacity coefficient at infinite dilution for the gas is calculated in this real solvent. Other equilibria computed in this work (more specifically to plot pressure-composition and temperature-composition diagrams) are performed using a diphasic TP flash, always using the mole fractions of the real mixture, including ions. This flash algorithm equals the fugacity of all species (including ions) between liquid and vapor phases with respect to mass balance. It is worth noticing that the Born term used in the SAFT model makes the liquid fugacity coefficient of ions extremely low, resulting thus in a vapor composition very close to zero.

3 Review of experimental data

3.1 Pure components (H2, N2, O2)

In the context of underground gas storage, it is of primary importance for economic evaluations to calculate with very good accuracy the volume of the gas produced. Consequently, the parameterization of the light gases involved in this study (hydrogen, nitrogen, and oxygen) is carried out to reproduce experimental compressibility factors Z in temperature and pressure ranges of industrial interest:(33)

The compressibility factors of all these pure components are collected from the NIST database [36] for temperatures ranging from 273 to 473 K and pressures ranging from 1 to 300 bar.

3.2 Brine

The brine is considered to be an aqueous NaCl solution since halite is the most abundant salt found in salt caverns. NaCl is assumed fully dissociated in Na+ and Cl− ions. One of the most important properties to be correctly modeled is the Mean Ionic Activity Coefficient (MIAC) of NaCl salt since it describes the deviation from the ideality of the aqueous solution. MIAC is defined from a single ion activity coefficient by:(34)

The e-PPC-SAFT equation of state computes fugacity coefficients φi, and not directly activity coefficients. The relationship between both for component i is:(35)where the exponent * is related to the reference state (infinite dilution for ions).

The MIAC experimental data covering a temperature range from 273 K to 473 K and salinities up to saturation are taken from [37]. Note that experimental ion activity coefficients are very often given in the molality scale, while the equation of state computes MIAC from equation (34) in the molar fraction scale. A conversion can be done using the molar fraction of the solvent:(36)where the exponents (m) and (x) stand for molality and molar fraction scale, respectively.

The composition of the gas produced during storage cycles is also of primary importance to designing surface facilities. More specifically, the water content of the vapor phase must be accurately known. Thus, special attention to our modeling is also paid to correctly reproduce the vapor pressure of the brine. Experimental data for this property on the same temperature and salinity ranges are taken from [38, 39]. Finally, a set of experimental brine densities [40] are also used to parameterize the model.

3.3 Light gas + pure water systems

Concerning hydrogen, Lopez-Lazaro et al. published recently a review of experimental hydrogen solubility in pure water data [6]: around 300 data from 273 to 623 K were collected and reconciliation in terms of Henry constants. These data were fitted following the empirical correlation suggested by Trinh et al. [41]:(37)where Tr is the reduced temperature of pure water (Tr = T/TC where TC is the critical temperature of pure water: 647 K). is the vapor pressure of the solvent, taken from the DIPPR database [42]. a, b, c are three empirical parameters given in Table 2.

Table 2

Parameters of the correlation fitting the experimental Henry constant in pure water (Eq. (37)).

A similar review is proposed in this work for the solubility of nitrogen and oxygen in pure water. For both nitrogen and oxygen, experimental data are collected and converted in terms of the Henry constant. In the literature, the solubility data can be directly reported as Henry constant, or also in terms of Kuenen coefficient Si, Ostwald coefficient Li, Bunsen coefficient αi, and molar fraction as a function of temperature and pressure (“TPx” data). For this last kind of data, we retain here only data at low pressure, since conversion into Henry constant requires to be close to infinite dilution. The reference of experimental data collected is given in Table 3. The conversion rules in Henry constant recommended by IUPAC [43] are:(38)where Tθ and Pθ are standard temperature (273.15 K) and pressure (1.013 bar), is the liquid molar volume of pure water at T, taken from the DIPPR databank, and MS the molar mass of pure water. From this molar fraction and equation (32), Henry constants are then estimated by:(39)where it is assumed that component i is pure in the vapor phase (composition of water is neglected), and is the fugacity coefficient of pure component i at temperature T and pressure P calculated with e-PPC-SAFT model and parameters of Table 4. To check the assumption that water in the vapor phase is negligible, the pure water vapor pressure is calculated for each experimental data and compared to the experimental total pressure. It is found that on average, the pure water vapor pressure represents around 3% of the total pressure, which is less of the experimental uncertainty. This suggests that this assumption is valid.

Table 3

Experimental data used for model parameterization and model validation.

Table 4

Pure component parameters in the e-PPC-SAFT model. Parameters in bold are those adjusted in this work. For O2, N2, and H2, parameters are adjusted on compressibility factors. For ions, parameters are adjusted on MIAC, brine vapor pressures, and brine densities.

The reconciliation experimental data are plotted in Figure 1. An average experimental uncertainty is estimated to be ±10%. These data are fitted using empirical equation (38) (bold line in Fig. 1). The parameters of this correlation are given in Table 2. Note that for O2, it can be observed that a few experimental data published are not consistent with most of the data. These inconsistent data have not been considered to establish the correlation.

thumbnail Fig. 1

Experimental Henry constants of N2 and O2 in pure water. Blue diamonds: from direct Henry constant data. Red squares: from Ostwald coefficient data. Green triangles: from Kuenen coefficient data. Purple crosses: from Bunsen coefficient data. Blue stars: from PTx data. The solid line is a fit of all data with correlation 38. The dotted lines are guides indicated ±10% from the correlation.

3.4 Light gas + brine systems

Experimental hydrogen solubility data in saline water are relatively scarce in literature and cover only a limited range of temperature (up to 333 K) and salinities (2 mol/kg) [44–47]. Lopez-Lazaro et al. reviewed these data and proposed additional ones at higher temperatures (up to 500 K) from molecular simulations [6]. This set of solubility data converted into the Henry constant is used to parameterize our model. Very recently, Chabab et al. [8] and Torín-Ollarves and Trusler [13] published new solubility data for temperatures up to 423 K and salinities up to 5 mol/kg. These new data are used in this work as a validation dataset for the model.

In contrast to hydrogen, the experimental data of oxygen solubility in salted water cover a wider range of temperature (up to 473 K) and salinities (5.8 mol/kg) [48–51]. As indicated in Table 3, the recent experimental oxygen solubility data of Chabab et al. [52], not used in the parameterization database, will be used to validate the model.

Concerning nitrogen, experimental solubility data in saline water are very scarce in the literature. The published Henry constants (or data which can be converted in Henry constant following Eq. (39)) concern a very restricted temperature range (up to 298 K) [53, 54], and cannot be used alone to parameterize our model in an underground gas storage context. Thus, in contrast to other systems, the binary interaction parameter between N2 and ions is directly fitted to match “PTx” experimental data at elevated pressure, covering a wider range of temperature (303 to 398 K) and salinities (up to 5.5 mol/kg) [55, 56].

4 Model parameterization

4.1 Pure components (H2, N2, O2)

For hydrogen, nitrogen, and oxygen, the experimental compressibility factors are used to adjust the dispersion energy and segment diameter of each of these light gas, while segment numbers are taken from previous studies [92, 93]. Table 4 gives the optimized parameters, and Figure 2 compares calculated values and reference data.

thumbnail Fig. 2

Compressibility factor of H2, N2, and O2. Symbols: reference date from NIST database [36]. Lines: this model.

The average deviations between model results and the reference data are equal to 0.1%, 0.2%, and 0.1% for hydrogen, nitrogen, and oxygen, respectively. These deviations are less than the uncertainties reported by the NIST for the reference data (0.2%, 0.5%, and 1.2%, respectively). The parameters fitted are consistent with other PC-SAFT parameter sets available in the literature [94–97].

4.2 Brine

For the brine (mixture of H2O, Na+, and Cl−), special attention should be paid to accurately reproduce MIAC close to NaCl saturation, since it corresponds to conditions encountered in salt caverns. It requires tuning not only pure ion parameters but also binary parameters between cation and anion since at high ionic concentration, ion pairing can be observed [98]. For pure ion parameters, we consider the hard sphere diameter equal to the Pauling diameter, and a segment number equal to 1.

In the e-PPC-SAFT framework, ions are assumed to interact with other particles only through association forces, and not through dispersion forces. Consequently, the energy dispersion parameter of ions is equal to 0. Finally, there are three parameters per ion to be adjusted: association energy, association volume, and MSA diameter. To better mimic ion-pairing formation at high salt concentration, the cross-binary parameters on association energy wij and association volume uij (Eq. (16)) between Na+ and Cl− are also adjusted. A total of seven parameters is then used for brine parameterization. The e-PPC-SAFT parameters of pure water are directly taken from a previous study [16]. Finally, the optimized parameters are reported in Table 4 (for pure component parameters), and on Table 5 (for binary parameters). Figure 3 compares experimental MIAC [37], vapor pressures [38, 39], and densities [40] with calculated values.

thumbnail Fig. 3

NaCl MIAC, vapor pressure, and density of NaCl aqueous solutions. Symbols: experimental data [37–40]. Lines: this model.
thumbnail Fig. 4

Henry constant of H2 (diamonds), N2 (triangles), and O2 (circles) in pure water. Symbols: experimental data (Eq. (38)). Lines: this model.
Table 5

Summary of binary parameters adjusted in this work for the e-PPC-SAFT model.

A very good agreement is found in the whole range of temperatures and salinities. The average deviations between experimental and calculated MIAC, vapor pressures, and densities are equal to 2.0%, 1.3%, and 1.4%, respectively.

4.3 Light gas + pure water systems

Correlation 37 is used as reference data to adjust interaction binary parameters of the e-PPC-SAFT model between H2 and H2O, N2 and H2O, and O2 and H2O. This optimization is done for temperatures ranging from 273 to 473 K in order to focus on conditions of interest for gas underground storage. As mentioned in Table 1, two binary parameters can be tuned: the kij and lij parameters, involved in the cross energy dispersion and the cross diameter of the NAHS term, respectively. A previous study based on molecular simulations of the H2 + H2O system showed that the lij parameter has a stronger influence on Henry’s constant calculations [6]. Consequently, only this binary parameter is adjusted, and the kij parameter is set to 0. A preliminary optimization using a constant value of the lij parameter led to important deviations between experimental and calculated values. Furthermore, the model was not able to reproduce the maximum in the Henry constant curves observed experimentally. Thus, it was decided to introduce a temperature dependency for this binary parameter:(4)

The adjusted , and parameters are given in Table 5. It is highly recommended to use this correlation only in the subcritical domain (temperature less than 647 K) to avoid the too-large value of lij. Figure 3 presents the experimental data calculated data. This optimization leads to a good agreement, with an average deviation between 1.8%, 1.1%, and 1.3% for systems H2 + H2O, N2 + H2O, and O2 + H2O, respectively.

4.4 Light gas + brine systems

In the e-PPC-SAFT model, ions are assumed to interact with solvent and neutral solutes only through repulsive and associative forces (no dispersion interaction). Thus, the lij interaction parameter of the NAHS is the only binary parameter available to tune the model. This binary parameter between Na+ and X (X = H2, O2 or N2) and between Cl− and X are taken equal since it is not possible to distinguish these two interactions independently in the brine. For H2 + salted water and O2 + salted water systems, the optimization of this parameter shows that a constant value (given in Tab. 5) is sufficient to reach a good accuracy as shown in Figure 5 and Figure 6 with an average deviation between experimental and calculated Henry constant of 3% and 8%, respectively. Note that for oxygen, experiments are not reported for fixed salinity values. Thus, the graphs in Figure 6 show experimental data gathered in specific salinity ranges.

thumbnail Fig. 5

Henry constant of hydrogen in salted water. Lines: this model. Symbols: experimental data from the parameterization database (Table 3). Blue diamonds: 0.5 mol/kg. Red squares: 1 mol/kg. Green triangles: 2 mol/kg.
thumbnail Fig. 6

Henry constant of oxygen in salted water. Lines: this model for salinity of 0.9 mol/kg (blue), 3 mol/kg (green), and 5.3 mol/kg (orange). Symbols: experimental data. Blue diamonds: from 0.6 to 1 mol/kg. Green triangles: from 2.8 to 4 mol/kg. Orange circles: from 5.1 to 5.8 mol/kg.

For the N2 + saline water system, a temperature dependency is found necessary to reach an acceptable average deviation (9%). The optimized parameters are reported in Table 5, and Figure 7 shows a comparison between experimental and calculated data.

thumbnail Fig. 7

N2 solubility in salted water at 303 K (left) and 398 K (right). Lines: this model. Symbols: experimental data. Blue triangles: 1 mol/kg. Green circles: 2.7 mol/kg. Orange diamonds: 5.5 mol/kg. Purple squares: 4.2 mol/kg.

5 Model validation

For each system, a set of experimental data is excluded from the model parameterization database and is used to validate our model. The model is validated not only on its capability to correctly predict data not used in the parameterization database, but also to extrapolate to mixtures. This last validation is of primary importance to accurately evaluates phase equilibrium between brine and air for CAES units. Air is considered in this work as a mixture of N2 (80%mol.) and O2 (20%mol.). Table 3 provides a synthetic view of the references of experimental data used for both model parameterization and model validation.

For H2 + H2O, N2 + H2O, and O2 + H2O systems, the model predicts very satisfactory solubility data not used in the parameterization database, with an average deviation of 5%. The comparisons between experimental and predicted bubble curves are reported in the Supplementary Material (Fig. S1). More interestingly, the model predicts very well the water content in the vapor phase, as illustrated in Figure 8, with an average deviation of 3.5%. Finally, the extension of the model to multicomponent systems is also validated by its capacity to accurately predict air solubility in water in a large range of temperatures and pressure, as illustrated in Figure 9, with an average deviation of 7%.

thumbnail Fig. 8

Dew curves of H2 + H2O (left) and N2 + H2O (right). Lines: predictions with this model. Symbols: experimental data not used in the parameterization database [57]. Red squares: 366 K, green triangles: 422 K; blue diamonds: 310 K.
thumbnail Fig. 9

Predictions of Air + H2O bubble curves at 1 bar (left) and 98 bar (right). Lines: this model. Symbols: experiments not used in parameterization database.

As for gas + pure water systems, the model proposed for gas + salted water systems is validated on its capability to predict data not included in the parameterization database as indicated in Table 3. For the H2 + H2O + NaCl system, recent hydrogen solubility data are used for validation [8]. As illustrated in Figure 10, our model predicts with good accuracy these data in the whole temperature range (323 to 373 K) and salinity range (1–5 mol/kg). The average deviation between experiments and prediction is about 5%, while the experimental uncertainty estimated by the authors is very close (4%).

thumbnail Fig. 10

Predictions of hydrogen solubility in salted water at 323 K (left) and 373 K (right). Lines: this model. Symbols: experimental data not used in the parameterization database [8]. Orange diamonds: 1 mol/kg. Green triangles: 3 mol/kg. Blue circles: 5 mol/kg.

For O2 + H2O + NaCl system, the recent experimental solubility data of Chabab et al. [52], not used in the parameterization database, are used to validate the model. As illustrated in Figure 11, the model predictions are found in good agreement with these experiments, with an average deviation of 7%. It can be nevertheless noticed that predictions are less accurate at the highest salinity (deviation of 15% at 4 mol/kg).

thumbnail Fig. 11

Predictions of oxygen solubility in salted water at 323 K (left) and 373 K (right). Lines: this model. Symbols: experimental data not used in the parameterization database [52]. Purple squares: 0.5 mol/kg. Orange diamonds: 1 mol/kg. Gray triangles: 4 mol/kg.

For N2 + H2O + NaCl system, due to the lack of experimental data, all the existing data have been used for the model parameterization. Finally, for Air + H2O + NaCl, no data have been found in the literature to validate model predictions. It could be interesting for further studies to acquire new experimental data on such a system and compare them to our predictions.

6 Application to air storage in salt caverns

The e-PPC-SAFT equation of state parameterized in this work is used to evaluate some thermodynamic properties in the case of Compressed Air Energy Storage (CAES) in the salt cavern. Air is assumed to be a mixture of 80% N2 + 20% O2. The cavern considered in this study has a volume of 315.000 m3, with a residual volume of the brine of 20.000 m3. This brine is made of saturated halite (water + NaCl brine, with a salt concentration of 6 mol/kg). The following cycles are considered for this CAES: a cycle lasts for one day, constituted by 12 h of withdrawing during which both temperature and pressure decrease to the minimum storage conditions, and 12 h of injection during which temperature and pressure reach the maximum storage conditions. The pressure and temperature ranges have been obtained by using the GUSTv2 software [99, 100], modeling the air (20% O2 and 80% N2) thermodynamic behavior in the cavern and the heat transfer with the salt mass: the fluctuation of pressure and temperature obtained are 141–166 bar and 35–50 °C, respectively. The mass transfer between the air and the residual brine is not considered in this model. Injection and withdrawing flows are assumed constant for each cycle. An illustration of temperature and pressure variations during an 8 days operation is provided in Supplementary Material (Fig. S2).

For each temperature/pressure couple of the cycle, the e-PPC-SAFT model is used to evaluate the compressibility factor Z of the vapor phase. The evolution of this factor over the cycles is illustrated in Figure S3 of the Supplementary Material. This factor enables computing the total quantity of gas mixture in this phase. The e-PPC-SAFT is also used to compute phase equilibrium for each point of the cycle in order to determine more specifically the water content of the vapor phase. The knowledge of this gas water content during CAES operation is of primary importance, since water may increase corrosion issues in surface facilities and promote hydrate formation in pipes. Hence, Figure 12 illustrates the evolution of the water content in the gas produced during the cycles, as well as the corresponding mass of gaseous water present in the cavern. In order to highlight the effect of salinity, two calculations are carried out: in a hypothetical salt-free system, and a saturated halite brine, representative of salt caverns currently operated. The quantity of air solubilized in the brine is also plotted in Figure 13.

thumbnail Fig. 12

Prediction of the water content of the vapor phase (top), and its corresponding mass according to the dimensions of the cavern (bottom). Solid line: in a salt cavern, with residual brine at its bottom. Dashed line: in a hypothetical cavern with residual (salt-free) water.
thumbnail Fig. 13

Prediction of the quantity of air solubilized in the residual brine. Solid line: in a salt cavern. Dashed line: in a hypothetical cavern with residual (salt-free) water.

These figures show that in a salt cavern, with residual saturated brine at its bottom, the quantity of air solubilized is lower than what would be found in a hypothetical cavern with residual pure water, as a result of the “salting-out” effect. In presence of salt, water molecules mainly solvate ions and are less available to solubilize gas, leading also to lower volatility. As a result of the mass balance, the water content in the vapor phase is 30% less in salt cavern than in pure water. Note that the calculated quantities result only from thermodynamic considerations and neglect kinetic aspects, and more specifically the time needed to reach the thermodynamic equilibrium. Temperature gradients in the cavern are also neglected, although due to vaporization heat being absorbed in the brine-gas interface only, at the bottom of the cavern, while condensation heat release is distributed on the whole cavern volume, some conditions can create a cold gas layer just above the brine, preventing natural convection and may increase the time to reach this equilibrium in the whole cavern [101]. These assumptions necessarily cause an overestimation of the water content for short injection/withdrawing cycles. Nevertheless, this approach provides a first conservative estimation of gas moisture, useful to design surface dehydration units [102].

7 Conclusion

In this work, a parameterization of the e-PPC-SAFT equation of state has been carried out to predict phase equilibrium in salt caverns for hydrogen or compressed air storage. The amount of gas stored in the cavern being a key economical factor, the pure component parameters of the light gases (H2, N2, O2) have been fitted to accurately match reference compressibility factors. Deviations are in the same range of experimental uncertainties. The brine properties have also to be carefully modeled, to accurately predict gas solubilities, phase density, and water volatility. Brine is assumed to be a mixture of water, Na+, and Cl−, and model parameters have been fitted to match experimental mean ionic activity coefficients, brine vapor pressure, and brine densities. The average deviation between experimental and calculated data are 2.0, 1.3, and 1.4%, respectively, for temperatures ranging from 273 to 473 K, and salinities up to saturation (6 mol/kg). The solubility of the light gases is first investigated in pure water. For N2 and O2, a new Henry constant correlation is proposed based on a literature review of existing data, while for H2 an existing recent correlation is used. These correlations are used to fit the binary interaction parameters of the e-PPC-SAFT model between H2 + H2O, N2 + H2O, and O2 + H2O. Calculated Henry constants in pure water reproduce reference experimental correlation with good accuracy (average deviation less than 2%).

A similar approach is carried out for systems involving saline water. A literature review of the Henry constant of light gases in saline water is carried out, and binary interaction parameters of the e-PPC-SAFT model between H2 + ion, N2 + ion, and O2 + ion have been fitted to match experimental data. Average deviations between 3% and 9% are obtained, in the same range of experimental uncertainties.

This parameterized model is validated on its capability to accurately predict data or properties not used in the parameterization database, but also on its capability of predicting mixture properties. Thus, it has been used to predict Air + pure water phase equilibrium with good accuracy (deviation less than 7% on air molar fraction solubilized), the air being considered as a mixture of 80%mol. N2 + 20%mol. O2. The model has also been used to predict air + saline water phase equilibrium to evaluate the water content of the vapor phase for an application of CAES in a salt cavern. It can be noticed that there are no experimental data available for air + saline water phase equilibrium. Thus, it could be useful in further studies to acquire new experimental data and to compare the predictions generated by this new model. It is also worth noticing that the case study investigated in this work only deals with a thermodynamic approach. Kinetic and diffusive effects are neglected, but this approach can nevertheless be used as a first conservative evaluation to design surface dehydration units.

Supplementary material

thumbnail Figure S1.

Bubble curves of systems H2 + H2O (top left), N2 + H2O (top right) and O2 + H2O (bottom). Lines: predictions by the model. Symbols: experimental data not used in the model parameterization. Red squares: 366 K. Green triangles: 422 K. Blue diamonds: 373 K. Purple triangles: 436 K.
thumbnail Figure S2.

Temperature and pressure variations during 8 cycles of CAES unit. Solid line: temperature. Dashed lines: pressure.
thumbnail Figure S3.

Evolution of the vapour phase compressibility factor during the cycles, calculated from the e-PPC-SAFT model.

References

All Tables

Table 1

e-PPC-SAFT pure-component and binary parameters.

In the text
Table 2

Parameters of the correlation fitting the experimental Henry constant in pure water (Eq. (37)).

In the text
Table 3

Experimental data used for model parameterization and model validation.

In the text
Table 4

Pure component parameters in the e-PPC-SAFT model. Parameters in bold are those adjusted in this work. For O2, N2, and H2, parameters are adjusted on compressibility factors. For ions, parameters are adjusted on MIAC, brine vapor pressures, and brine densities.

In the text
Table 5

Summary of binary parameters adjusted in this work for the e-PPC-SAFT model.

In the text

All Figures

thumbnail Fig. 1

Experimental Henry constants of N2 and O2 in pure water. Blue diamonds: from direct Henry constant data. Red squares: from Ostwald coefficient data. Green triangles: from Kuenen coefficient data. Purple crosses: from Bunsen coefficient data. Blue stars: from PTx data. The solid line is a fit of all data with correlation 38. The dotted lines are guides indicated ±10% from the correlation.
In the text
thumbnail Fig. 2

Compressibility factor of H2, N2, and O2. Symbols: reference date from NIST database [36]. Lines: this model.
In the text
thumbnail Fig. 3

NaCl MIAC, vapor pressure, and density of NaCl aqueous solutions. Symbols: experimental data [37–40]. Lines: this model.
In the text
thumbnail Fig. 4

Henry constant of H2 (diamonds), N2 (triangles), and O2 (circles) in pure water. Symbols: experimental data (Eq. (38)). Lines: this model.
In the text
thumbnail Fig. 5

Henry constant of hydrogen in salted water. Lines: this model. Symbols: experimental data from the parameterization database (Table 3). Blue diamonds: 0.5 mol/kg. Red squares: 1 mol/kg. Green triangles: 2 mol/kg.
In the text
thumbnail Fig. 6

Henry constant of oxygen in salted water. Lines: this model for salinity of 0.9 mol/kg (blue), 3 mol/kg (green), and 5.3 mol/kg (orange). Symbols: experimental data. Blue diamonds: from 0.6 to 1 mol/kg. Green triangles: from 2.8 to 4 mol/kg. Orange circles: from 5.1 to 5.8 mol/kg.
In the text
thumbnail Fig. 7

N2 solubility in salted water at 303 K (left) and 398 K (right). Lines: this model. Symbols: experimental data. Blue triangles: 1 mol/kg. Green circles: 2.7 mol/kg. Orange diamonds: 5.5 mol/kg. Purple squares: 4.2 mol/kg.
In the text
thumbnail Fig. 8

Dew curves of H2 + H2O (left) and N2 + H2O (right). Lines: predictions with this model. Symbols: experimental data not used in the parameterization database [57]. Red squares: 366 K, green triangles: 422 K; blue diamonds: 310 K.
In the text
thumbnail Fig. 9

Predictions of Air + H2O bubble curves at 1 bar (left) and 98 bar (right). Lines: this model. Symbols: experiments not used in parameterization database.
In the text
thumbnail Fig. 10

Predictions of hydrogen solubility in salted water at 323 K (left) and 373 K (right). Lines: this model. Symbols: experimental data not used in the parameterization database [8]. Orange diamonds: 1 mol/kg. Green triangles: 3 mol/kg. Blue circles: 5 mol/kg.
In the text
thumbnail Fig. 11

Predictions of oxygen solubility in salted water at 323 K (left) and 373 K (right). Lines: this model. Symbols: experimental data not used in the parameterization database [52]. Purple squares: 0.5 mol/kg. Orange diamonds: 1 mol/kg. Gray triangles: 4 mol/kg.
In the text
thumbnail Fig. 12

Prediction of the water content of the vapor phase (top), and its corresponding mass according to the dimensions of the cavern (bottom). Solid line: in a salt cavern, with residual brine at its bottom. Dashed line: in a hypothetical cavern with residual (salt-free) water.
In the text
thumbnail Fig. 13

Prediction of the quantity of air solubilized in the residual brine. Solid line: in a salt cavern. Dashed line: in a hypothetical cavern with residual (salt-free) water.
In the text
thumbnail Figure S1.

Bubble curves of systems H2 + H2O (top left), N2 + H2O (top right) and O2 + H2O (bottom). Lines: predictions by the model. Symbols: experimental data not used in the model parameterization. Red squares: 366 K. Green triangles: 422 K. Blue diamonds: 373 K. Purple triangles: 436 K.
In the text
thumbnail Figure S2.

Temperature and pressure variations during 8 cycles of CAES unit. Solid line: temperature. Dashed lines: pressure.
In the text
thumbnail Figure S3.

Evolution of the vapour phase compressibility factor during the cycles, calculated from the e-PPC-SAFT model.
In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.