*Article* **Elemental Abundances of Moon Samples Based on Statistical Distributions of Analytical Data**

**Zhiguan Hou 1, Qingjie Gong 1,\*, Ningqiang Liu 1, Biao Jiang 2,\*, Jie Li 1, Yuan Wu 1, Jiaxin Huang <sup>1</sup> and Weixuan Gu <sup>1</sup>**

<sup>1</sup> School of Earth Sciences and Resources, China University of Geosciences, Beijing 100083, China

<sup>2</sup> Institute of Mineral Resources, Chinese Academy of Geological Sciences, Beijing 100037, China

**\*** Correspondence: qjiegong@cugb.edu.cn (Q.G.); jiangbiao334223@163.com (B.J.)

**Abstract:** The successful return of Chang'E-5 (CE5) samples urges the hot topic of the study of the Moon in geochemistry. The elemental data of the analyzed moon samples reported in the literature were collected to determine the elemental abundances in moon samples. Based on 2365 analytical records of moon samples from ten missions of Apollo, Luna, and CE5, elemental abundances of 11 major oxides including Cr2O3, 50 trace elements including Ti, P, Mn, Cr, and 15 rare earth elements (REEs) including Y are derived based on statistical distributions of normal, log-normal, and additive log-ratio transformation, respectively. According to the value of 13.5% CaO content, moon samples are classified into two types, as low-Ca and high-Ca samples, whose elemental abundances are also calculated respectively based on the methods used in the total moon samples. With respect to the mid-ocean ridge basalt (MORB) of the Earth, moon samples (including the Moon, low-Ca, and high-Ca samples) are rich in Cr, REEs, Th, U, Pb, Zr, Hf, Cs, Ba, W, and Be and poor in Na, V, Cu, and Zn in terms of their concentrations, and are enriched in Cr and depleted in Na, K, Rb, P, V, Cu, Zn in spider diagrams. The CE5 sample is a low-Ca type of moon sample and is clearly rich in Ti, Fe, Mn, P, Sc, REEs, Th, U, Nb, Ta, Zr, Hf, Sr, Ba, W, and Be and poor in Mg, Al, Cr, and Ni in terms of their concentrations relative to the moon or the low-Ca samples. If compared with the moon sample, the CE5 sample is also clearly rich in K, REE, and P.

**Keywords:** lunar samples; CE5; low-Ca; high-Ca; spider diagram

## **1. Introduction**

The successful return of Chang'E-5 (CE5) has marked China as the third country to retrieve moon samples after the United States and the former Soviet Union. Studies on the returned samples have been an interesting and hot topic in geochemistry recently [1–11]. Elemental abundance is a basic topic in geochemistry, such as in the Earth [12–16]; therefore, the elemental composition of CE5 samples and moon samples from the Apollo and Luna missions is very attractive to geochemists.

On the elemental compositions of moon samples, Rose et al. [17–21] proposed the compositions from the Apollo missions. Samples in each Apollo mission were divided into the types of soil and rock, and the average elemental concentrations were used to represent their compositions. However, the analyzed samples or sample count from each mission was mostly less than 30, and items of major oxides, trace elements, and rare earth elements (REEs) were limited, such as 11 major oxides (including Cr2O3), 14 trace elements, and 3 REEs (Table 1). In addition, Taylor et al. [22] proposed a set of average elemental compositions in the lunar highland, including 8 major oxides (lack of P2O5 and MnO), 16 trace elements, and 14 REEs based on 4 samples from Apollo 17 and 6 samples from Apollo 16 (Table 1). Korotev [23] proposed the compositions of Apollo 16 soils, including 8 major oxides, 12 trace elements, and 8 REEs based on 8 sampling stations with 2~6 samples in each station (Table 1). Warren and Taylor [24] proposed the compositions of mare basalts

**Citation:** Hou, Z.; Gong, Q.; Liu, N.; Jiang, B.; Li, J.; Wu, Y.; Huang, J.; Gu, W. Elemental Abundances of Moon Samples Based on Statistical Distributions of Analytical Data. *Appl. Sci.* **2023**, *13*, 360. https:// doi.org/10.3390/app13010360

Academic Editors: Dibyendu Sarkar and Andrea L. Rizzo

Received: 30 November 2022 Revised: 22 December 2022 Accepted: 25 December 2022 Published: 27 December 2022

**Copyright:** © 2022 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

and highland regolith on the Moon. The compositions of mare basalts include 8 major oxides, 10 trace elements, and 5 REEs derived from the elemental averages of mare basalts from all the Apollo and Luna missions except Apollo 16 and Luna 20. The sample count of the mare basalts is less than 200 from the Apollo missions and less than 20 from the Luna missions. The compositions of highland regolith, including 9 major oxides (lack of P2O5), 12 trace elements, and 9 REEs, were derived based on the statistical average of soils from Apollo 16, 2 regolith breccias from Apollo 14, and 7 lunar meteorites (Table 1). Although the elemental data from CE5 samples have been reported recently [1–3,5–7], the elemental compositions or abundances of moon samples are incomplete, which were derived based on only a few analytical samples with limited items (Table 1).


**Table 1.** Sample counts and element counts of analyzed moon samples.

In this paper, the analytical data reported in the literature on moon samples from Apollo, Luna, and CE5 missions were compiled firstly. Then the elemental abundances of moon samples were derived based on statistical distributions of the analytical data. Thirdly, the moon samples were classified into two types according to their concentrations of CaO, and the elemental abundances of each type were also calculated. Finally, the elemental compositions of CE5 samples were compared with the newly derived elemental abundances of moon samples.

## **2. Samples and Analytical Methods**

## *2.1. Samples*

The exploration and research of the Moon have been launched since 1957 [25]. In 1969, the Apollo 11 mission realized a manned moon landing and retrieved moon samples, which shifted the research on the Moon from theoretical conception to practical analysis [26–31]. Then in 1970, Luna 16 realized an unmanned moon landing sampling [32–34]. In 2020, Chang'E-5 (CE5) brought back moon samples [35]. So far, humans have retrieved moon samples 10 times, including 6 Apollo missions, 3 Luna missions, and the Chang'E-5 mission (Figure 1). The uppermost few meters of the Moon's crust, from which all the moon samples came, is a layer of loose, highly porous regolith or moon soil [24]. Samples from the Apollo missions contain lumps of rocks, breccias, surface soils, and core soils. The Apollo missions brought back a total of 380.95 kg of samples (Table 2), of which rock and breccia samples were picked up by the astronauts on the surface, and soil samples were mainly extracted by the surface scoop and deep drilling core, with a depth of 2 to 3 m. The Luna missions brought back a total of 0.301 kg of samples (Table 2). The sampling method was mechanical core drilling with a depth of 0.2 to 1.6 m. As for the CE5 mission, it brought back a total of 1.731 kg of lunar samples [1–3,5–7]. The sampling methods were mechanical shovel sampling and drilling sampling, but the drilling samples have not been analyzed or reported until now.

**Figure 1.** The landing sites of Apollo, Luna, and Chang'E-5 missions.


**Table 2.** The counts of analyzed records collected in this paper from each mission with other information.

## *2.2. Analytical Methods*

The Apollo samples were first analyzed by the Lunar Sample Preliminary Examination Team of NASA. It reported a few samples' composition data, including major oxides and some trace elements without REEs for each mission [26–31]. Later, other experts and scholars applied to NASA for samples, of which the analytical composition data have been reported one after another [17–21,36,37]. The samples from the Luna missions were first analyzed and tested by the former Soviet Union [32–34]. Then after exchanging with the United States, the composition data was also analyzed and reported [32–34,38–40]. The composition data of CE5 samples have been analyzed and reported only recently [1–3,5–7].

The main analytical methods used to determine the composition data of moon samples are X-ray fluorescence spectrometry (XRF), inductively coupled plasma mass spectrometry (ICP-MS), instrumental neutron activation analysis (INNA), and electron microprobe analysis (EPMA), which can analyze many items such as major oxides, some trace elements, and/or REEs simultaneously [17,36–38]. In addition, other analytical methods such as isotope dilution analyses (IDA), emission spectrography (ES), and radiochemical neutron activation (RCNA) were adopted to analyze the composition of Rb, Sr, Zr, Re, Au, etc. [41–43]. The analytical quality, such as the detection limits, precisions, accuracies, relative errors, and relative standard deviations, were illustrated in the original literature. In general, the relative errors were less than 5% for major oxides and 10% for trace elements, including REEs.

All the analytical data collected in this paper are from moon samples brought back from the Apollo, Luna, and CE5 missions (excluding meteorites found on the Earth) and have been reported in the literature. The analyzed samples are soils, breccias, and rocks, and the counts of analytical samples (or records) in each mission or landing site are listed in Table 2. There are 2365 records of analytical data of moon samples used in this paper in total.

## **3. Statistical Methods and Results**

#### *3.1. Elemental Abundances in Moon Samples*

In order to calculate the elemental abundances in moon samples, a total of 2365 records of analytical data were used with equal weight, ignoring the different missions or landing sites, sample types, and analytical methods from the literature.

According to the rule of "The contents of major oxides commonly obey normal distribution", the average of the analytical data excluding outliers was calculated as the abundance for each major oxide here. Firstly, the average (Avg) and standard deviation (Std) of each oxide's data were calculated for all the samples (the count of samples is labeled as n0). Then, the boundary values of Avg ± 3Std were used to delete the outliers repeatedly until no outlier data were found. Finally, the average of each oxide's data without outliers (the count of samples is labeled as n) was calculated and viewed as its abundance. The abundances of 11 major oxides (including Cr2O3) along with counts of samples (n0 and n) are listed in Table 3, and statistical histograms of major oxides are illustrated in Figure 2. The sum of the 11 major oxides (or Total in Table 2) is 99.10%, which is close to the closure value of 100%.

**Figure 2.** Histograms of major oxides with their counts of samples (n0 is the count of total analytical data and n is the count of analytical data without outliers). The histogram of K2O contents is drawn without the three highest values of 5.8, 3.2, and 3.1, and the histogram of P2O5 contents is without the two highest values of 1.38 and 1.1.


**Table 3.** Elemental abundances in moon samples with their counts of samples (n0 and n).

Note: The units of major oxides and trace elements (including REEs) are % and μg/g, respectively, except Ag, Au, Cd, Re, Ru, Rh, Pd, Os, Ir, Pt, which are in ng/g.

According to the rule of "The contents of trace elements commonly obey log-normal distributions", the geometric average of analytical data excluding outliers was calculated as the abundance for each trace element. The elemental abundances of 50 trace elements (including Ti, P, Mn, and Cr), along with their counts of samples (n0 and n), are also listed in Table 3.

REE pattern is a useful tool for traceability or provenance in geochemistry [44,45], and the key signature is the variation trend of the pattern, which is dependent on REE concentrations. Therefore, the covariation of REE abundances needs to be considered. Here, the additive log-ratio (alr) transformation method [46,47] was adopted to calculate the REE (including Y) abundances, and Yb was selected as the denominator to calculate the ratios for other REEs. Yb was selected as the denominator of the alr transformation method because it not only has the largest count of analyzed samples (n0 = 1680 and n = 1652) but also obeys the log-normal distribution well relative to the other REEs. According to the geometric average of Yb without outliers, 5.84 μg/g was set as its abundance of moon samples and was used to calculate the abundances of the other REEs (including Y). The abundances of 15 REEs (including Y) with their counts of samples (n0 and n) are also listed in Table 3.

## *3.2. Elemental Abundances in Moon Samples with Low-Ca and High-Ca*

Although elemental abundances in moon samples were derived based on statistical distributions such as normal, log-normal, and alr-normal distributions, some items clearly deviate from their ideal distributions, such as CaO, FeO, Al2O3, TiO2, P2O5, Cr2O3, etc., as shown in Figure 2. It is worth mentioning that the distribution of CaO contents is near a bimodal distribution (Figure 2). In order to derive more meaningful elemental abundances, moon samples were classified into two types, as low-Ca and high-Ca samples, based on the CaO content boundary of 13.5% used in this paper.

In the total 2365 analyzed moon samples, there were only 1979 records with valid CaO contents. According to the content boundary of 13.5% CaO, there were 1468 records classified as low-Ca samples and 511 as high-Ca samples. Therefore, moon samples with low Ca were about three-quarters of the total moon samples collected in this paper, and samples with high Ca were about one-quarter.

With respect to the low-Ca moon samples and high-Ca moon samples, we used the same statistical methods as adopted for the total moon samples to calculate the elemental abundances of each type separately. The abundances of each type are also listed in Table 3, along with their counts of samples (n0 and n).

## **4. Discussion**

## *4.1. Geochemical Signatures of Elemental Abundances*

Based on the elemental abundances of the moon, low-Ca, and high-Ca samples in Table 2, we compared their geochemical signatures with the elemental abundances of carbonaceous chondrite CI [12], primitive mantle [13], bulk oceanic crust [14], mid-ocean ridge basalt (including those of Atlantic, India, and Pacific [14]), continental crust (including total continental crust, lower continental crust, middle continental crust, and upper continental crust [15]), and rocks of China (including acidic rock, intermediate rock, and basic rock of China [16]) and found moon samples are more close to the mid-ocean ridge basalt (MORB) of the Earth. Here, only the comparison results with the MORB were illustrated to derive the geochemical signatures of moon samples.

## 4.1.1. Major Oxides

According to the illustration method of spider diagrams, the 11 major oxides were first sequenced descending on their abundances of moon samples. Then, moon samples were normalized based on the MORB of the Earth. Finally, the spider diagrams of major oxides of moon samples were derived and illustrated in Figure 3.

**Figure 3.** Spider diagrams of major oxides of moon samples normalized based on MORB. Abundances of the MORB are from White and Klein [14]. Elements in the horizontal axis are the abbreviations of their major oxides listed in Table 2.

With respect to the MORB, the moon sample (labeled as Moon in Figure 3) is clearly rich in Cr2O3 and TiO2 and poor in K2O in terms of their concentrations. The moon sample with low-Ca (labeled as Low-Ca in Figure 3) is also clearly rich in Cr2O3 and TiO2 and poor in K2O, like the moon sample. While the moon sample with high-Ca (labeled as High-Ca in Figure 3) is clearly rich in Cr2O3 and Al2O3 and poor in Na2O, TiO2, MnO, and FeO in terms of their concentrations. In the three abundances of moon samples, the moon sample with low-Ca is closer to the MORB, except for the clear signature of higher concentrations of Cr2O3 and TiO2 and lower concentrations of K2O.

With respect to the MORB, the moon samples (including Moon, low-Ca, and high-Ca in Figure 3) are clearly enriched in Cr and depleted in Na in the diagrams. Here, the terms of rich and poor are used for concentrations (comparison between/among samples), and the terms of enriched and depleted are used for diagrams (comparison among elements in the same sample).

## 4.1.2. REEs

The REE patterns of the three abundances of moon samples (including Moon, low-Ca, and high-Ca) are illustrated in Figure 4.

**Figure 4.** REE patterns of moon samples normalized to the MORB. Abundances of the MORB are from White and Klein [14].

With respect to the MORB, the REE patterns of the Moon and low-Ca samples (Figure 4) are near flat, except for the clear negative Eu anomaly. While the pattern of high-Ca (Figure 4) is tilted right for light REEs and nearly flat for heavy REEs. Therefore, the Moon and low-Ca samples are closer to the MORB than the high-Ca samples in the REE pattern, except for the clear negative Eu anomaly. With respect to the absolute concentrations of REEs, the three abundances of moon samples (including Moon, low-Ca, and high-Ca) are all higher or richer than those of the MORB. The descending sequence of total REE concentrations is (CE5 discussed in the following), low-Ca, Moon, high-Ca, and the MORB (Figure 4).

#### 4.1.3. Trace Elements

Trace elements are illustrated using the spider diagram suggested by Sun and Mc-Donough [48], with 31 elements, including K, P, and Ti. In order to avoid the repetition of REEs, the Ce and Eu were deleted from the spider diagram in which Ce and Eu are following the La and Sm, respectively. Therefore, a total of 29 elements (including K, P, and Ti) were used to draw the spider diagrams of the moon samples, which are illustrated in Figure 5.

**Figure 5.** Spider diagrams of moon samples normalized to the MORB. Abundances of the MORB are from White and Klein [14].

With respect to the MORB in terms of their concentrations, the moon samples (including the Moon, low-Ca, and high-Ca) are all clearly rich in Cs, Ba, W, Th, U, and Pb (Figure 5).

With respect to the MORB, spider diagrams of moon samples (including the Moon, low-Ca, and high-Ca) are nearly flat (variations are limited to one order of magnitude), except for the clear depletion in Rb, K, and P (Figure 5). Furthermore, the diagram of high-Ca is also clearly depleted in Ti (Figure 5).

## 4.1.4. Other Trace Elements

Except for the aforementioned major oxides, REEs, and trace elements, there are only eight remaining elements of Sc, V, Co, Ni, Cu, Zn, Be, and B, which are reported abundances both of moon samples and the MORB. Here we supplement Ti, Cr, Mn, and Fe to the eight elements to form a series of the first transition elements plus Be and B. Therefore, 12 elements were used to form the spider diagrams of the moon samples which are illustrated in Figure 6.

**Figure 6.** Spider diagrams of moon samples on first transition elements plus Be and B normalized to the MORB. Abundances of the MORB are from White and Klein [14].

With respect to the MORB in terms of their concentrations, the moon samples (including the Moon, low-Ca, and high-Ca) are all clearly rich in Cr, Ni, and Be and poor in V, Cu, and Zn (Figure 6).

With respect to the MORB, the spider diagrams of moon samples (including the Moon, low-Ca, and high-Ca) are nearly flat (variations are limited to one order of magnitude), except for the clear enrichment in Cr and depletion in V, Cu, and Zn (Figure 6). In addition, the diagram for the high-Ca moon sample is also clearly enriched in Ni and B (Figure 6).

In summary, the geochemical signatures of major oxides, REEs, and trace elements in moon samples (including the Moon, low-Ca, and high-Ca) are rich in Cr, REEs, Th, U, Pb, Zr, Hf, Cs, Ba, W, and Be and poor in Na, V, Cu, and Zn in terms of their concentrations, and are enriched in Cr and depleted in Na, K, Rb, P, V, Cu, and Zn in the spider diagrams, relative to the MORB. Among the total moon samples, low-Ca samples, and high-Ca samples, the low-Ca samples are closer to the total moon samples in concentration or patterns (or diagrams), which is not only illustrated in Figure 3 to Figure 6 but also consistent with the counts of analytical data with equal weighting used in this study. Except the 46 elements discussed in this paper, the other 26 elemental abundances in moon samples are not discussed here because of the data lack on the MORB.

## *4.2. Chang'E-5 Samples*

## 4.2.1. Elemental Concentrations

Here we compile the analyzed elemental data of CE5 samples reported recently. The averages are used to represent the elemental concentrations of the CE5 samples which include 33 analyzed records, including two repetitions reported by different authors. If the elemental concentrations were only reported by Zong et al. [7], the suggested concentrations by Zong et al. [7] are used in this paper.

The calculated elemental concentrations of the CE5 samples are also listed in Table 3, including 11 major oxides, 15 REEs (including Y), and 28 trace elements (including Ti, P, Mn, and Cr).

## 4.2.2. Geochemical Signatures

The content of CaO of CE5 is 11.36%, which is lower than the boundary value of 13.5% for low-Ca and high-Ca moon samples. Therefore, the CE5 sample is the low-Ca type of moon sample. In terms of Ti content, moon samples can be divided into three types: high Ti (TiO2 ≥ 6%), low Ti (1% ≤ TiO2 < 6%), and very low Ti (TiO2 < 1%) [11]. The content of TiO2 of CE5 is 5.43%, which indicates that the CE5 sample is the low Ti type of moon sample.

According to the illustrations of moon samples, the geochemical signatures of the CE sample are also illustrated in Figure 3 to Figure 6.

With respect to the MORB, the CE5 sample is clearly rich in Cr2O3, TiO2, and FeO and poor in Na2O (Figure 3) in major oxides in terms of their concentrations, and is also enriched in Cr2O3, TiO2, and FeO and depleted in Na2O in the diagrams. The REE pattern of CE5 is tilted right slowly, except for the clear negative Eu anomaly, and its REE concentrations are clearly higher than those of the MORB (Figure 4). In the spider diagrams (Figures 5 and 6), the CE5 sample is clearly rich in Cs, Ba, W, Th, U, Nb, Ta, Pb, Zr, Hf, Ti, Li, Sr, Cr, and Be and poor in V, Cu, and Zn in terms of their concentrations, and is enriched in Cr and depleted in Rb, K, P, V, Cu, and Zn in the spider diagrams. Among the total moon sample, low-Ca sample, and high-Ca sample, the CE5 sample is closer to the low-Ca sample in concentrations or patterns (or diagrams), which are illustrated in Figure 3 to Figure 6. This is consistent with the low-Ca type of moon samples discriminated on CaO content, as mentioned previously.

With respect to the total moon sample, the CE5 sample is clearly rich in Ti, Fe, Mn, P, Sc, K, REEs, Th, U, Nb, Ta, Zr, Hf, Sr, Ba, W, and Be and poor in Mg, Al, Cr, and Ni in terms of their concentrations. From this view, the CE5 sample is the KREEP type of moon sample [49,50] because it is rich in K, REEs, and P relative to the moon sample. However, the contents of K2O in the CE5 sample and the low-Ca sample are 0.175% and 0.171%, respectively, which are almost the same within the relative error. From this view, the CE5 sample is the non-KREEP type of moon sample [1] because of the non-enrichment of K relative to the low-Ca moon sample.

Except for the aforementioned studies on the CE5, more and more articles are being published [51–54] on hot and interesting topics, which are very helpful in promoting the research of the Moon. The elemental abundances of the CE5 sample in this paper will be improved on with more studies in the near future.

## **5. Conclusions**


**Author Contributions:** Conceptualization, Z.H., Q.G., B.J., J.L. and Y.W.; Methodology, Q.G. and N.L.; Formal analysis, Z.H.; Data curation, Z.H., N.L., B.J., J.L., Y.W., J.H. and W.G.; Writing—original draft, Z.H.; Writing—review & editing, Q.G. and N.L. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research received no external funding.

**Institutional Review Board Statement:** Not applicable.

**Informed Consent Statement:** Not applicable.

**Data Availability Statement:** Not applicable.

**Acknowledgments:** We thank the authors who reported the geochemical data of moon samples used in this manuscript. We greatly appreciate the comments from the anonymous reviewers for their valuable suggestions to improve the quality of this manuscript.

**Conflicts of Interest:** The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

## **References**


**Disclaimer/Publisher's Note:** The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

**Jia Zhang 1,\*, Shide Mao <sup>2</sup> and Zeming Shi <sup>1</sup>**

<sup>2</sup> School of Earth Sciences and Resources, China University of Geosciences, Beijing 100083, China

**\*** Correspondence: jiazhang202211@163.com

**Abstract:** An equation of state (EOS) of CH4-N2 fluid mixtures in terms of Helmholtz free energy has been developed by using four mixing parameters, which can reproduce the pressure-volumetemperature-composition (*PVTx*) and vapor-liquid equilibrium (VLE) properties of CH4-N2 fluid mixtures. The average absolute deviation of all the *PVTx* data available up to 673.15 K and 1380 bar from this EOS is 0.38%. Combining this EOS of CH4-N2 fluid mixtures and the EOS of CH4-CO2 and CO2-N2 fluid mixtures in our previous work, an EOS of CO2-CH4-N2 fluid mixtures has been developed, which is named ZMS EOS. The ZMS EOS can calculate all thermodynamic properties of ternary CO2-CH4-N2 fluid mixtures and the average absolute deviation of the *PVTx* data from the ZMS EOS is 0.40% for the CO2-CH4-N2 system. The ZMS EOS can be applied to calculate excess enthalpies of CO2-CH4-N2 fluid mixtures, predict the solubility of CO2-CH4-N2 fluid mixtures in brine and water, and quantitatively estimate the impact of the impurities (CH4 and N2) on the CO2 storage capacity in the CO2 capture and storage (CCS) processes. The ZMS EOS can also be applied to calculate the isochores of CO2-CH4-N2 system in the studies of fluid inclusions. All Fortran computer codes and Origin drawing projects in this paper can be obtained freely from the corresponding author.

**Keywords:** CO2-CH4-N2 system; equation of state; CCS; excess enthalpies; fluid inclusion

## **1. Introduction**

The CO2, CH4 and N2 are important natural fluids, and non-aqueous CO2–CH4–N2 mixtures are often reported in the studies of fluid inclusions [1–3]. In recent decades, the amount of greenhouse gas CO2 in the atmosphere has gradually increased because of the development of industry. Numerous studies have shown that CO2 capture and storage (CCS) is an effective method to reduce the amount of CO2 in the atmosphere [4–6]. However, CO2 in industrial exhaust gases is usually impure, mixing with other components, such as CH4 and N2 [7,8]. The impurities can affect the design of the CCS processes [9–11]. Therefore, predicting thermodynamic properties of the CO2-CH4-N2 fluid mixtures, especially the pressure-volume-temperature-composition (*PVTx*) and vapor-liquid equilibrium (VLE) properties at different temperatures and pressures, is of great significance for the related CCS and fluid inclusion studies [12–15].

Although experimental thermodynamic data are the most reliable in the corresponding applications, they are costly and time-consuming to obtain. On the other hand, the predictive EOS is a better choice for us. In the past century, the cubic EOSs (e.g., Peng-Robinson (PR) EOS [16], Soave–Redlich–Kwong (SRK) EOS [17]) originated from the van der Waals EOS [18] are highly preferred for the CCS researchers because of their simplicity. Based on the standard PR EOS, some more accurate EOSs have been presented for the CCS mixtures, such as a model that integrates the classical PR EOS with advanced mixing rules, called the "Peng-Robinson + residual part of excess Helmholtz energy model" [19], and the Enhanced-Predictive-PR78 (E-PPR78) EOS [20]. However, the cubic EOSs cannot well reproduce the *PVTx* properties of fluids under high pressure-temperature conditions, as can be seen in

**Citation:** Zhang, J.; Mao, S.; Shi, Z. A Helmholtz Free Energy Equation of State of CO2-CH4-N2 Fluid Mixtures (ZMS EOS) and Its Applications. *Appl. Sci.* **2023**, *13*, 3659. https:// doi.org/10.3390/app13063659

Academic Editor: Andrea L. Rizzo

Received: 20 January 2023 Revised: 6 March 2023 Accepted: 7 March 2023 Published: 13 March 2023

**Copyright:** © 2023 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

Section 2. Based on Heyen EOS [21], Heyen et al. [22] simulated the phase equilibria of the CH4-CO2 system below 50 ◦C and 100 bar; Darimont and Heyen [23] predicted the phase equilibria of the CO2-N2 system below 20 ◦C and 90 bar. Because cubic EOSs cannot well reproduce the *PVTx* properties, Thiery et al. [24,25] and Thiery and Dubessy [26] modeled the phase equilibria of the CO2-CH4-N2 system including liquid, vapor and solid by the cubic EOSs and the *PVTx* properties by the Lee–Kesler correlation [27]. However, these models have the disadvantage of using two EOSs that may cause incoherency of various fluid parameters [28].

In recent years, the Statistical Associating Fluid Theory (SAFT) EOS [29] from Wertheim's thermodynamic perturbation theory (TPT) has become popular for CCS fluid mixtures. Because it is based on statistical mechanics and has an explicit physical meaning, it has more reliable extrapolation and more accurate predictive capability than traditional empirical models. Many successful modifications of the SAFT-based EOSs have been presented, such as the EOS of perturbed-chain SAFT (PC-SAFT) [30,31], variable range SAFT (SAFT-VR) [32], and Lennard-Jones SAFT (LJ-SAFT) [33]. However, the SAFT EOSs cannot well reproduce experimental VLE data in the near-critical region of the CO2-N2 fluid mixture, as can be seen in the previous work [34]. Up to now, the EOS in terms of Helmholtz free energy is also commonly used to calculate the thermodynamic properties of CCS fluid mixtures [35–39], because the EOS in the form of Helmholtz free energy has high precision, wide temperature-pressure range, and can calculate all thermodynamic properties of the fluids by the derivation from the Helmholtz free energy EOS, such as enthalpy, heat capacity, entropy, etc. The derivation from the Helmholtz free energy EOS is not complex, as can be seen later.

In our previous work [34,40], the EOS in terms of dimensionless Helmholtz free energy for the CH4-CO2 and CO2-N2 fluid mixtures have been established by using four mixing parameters, which have been applied to the studies of the CH4-CO2 and CO2-N2 fluid inclusions. In this work, firstly, the EOS in terms of dimensionless Helmholtz free energy of the CH4-N2 fluid mixture has been developed by the same approach. The ZMS EOS is obtained by combining the binary interaction parameters of CH4-CO2, CO2-N2 and CH4-N2 systems. Experimental data available for the binary CH4-N2 and ternary CO2-CH4-N2 fluid mixtures are used to verify the accuracy of the ZMS EOS. Then the ZMS EOS is applied to calculate the excess enthalpies, the solubility of CO2-CH4-N2 gas mixtures in brine and water, the impact of impurities (CH4 and N2) on the CO2 storage capacity, and the isochores of the CO2-CH4-N2 fluid inclusions.

## **2. The ZMS EOS**

The EOS of the CO2-CH4-N2 mixtures is in terms of dimensionless Helmholtz free energy *α*, which is represented by

$$
\mathfrak{a} = \mathfrak{a}^0 + \mathfrak{a}^r \tag{1}
$$

where *α*<sup>0</sup> and *α<sup>r</sup>* are the ideal-gas part and the residual part of dimensionless Helmholtz free energy. *α*<sup>0</sup> and *α<sup>r</sup>* are defined by

$$\boldsymbol{\alpha}^{0} = \sum\_{i=1}^{n} \boldsymbol{\pi}\_{i} \left[ \boldsymbol{\alpha}\_{i}^{0} + \ln \boldsymbol{\pi}\_{i} \right] \tag{2}$$

$$\alpha^r = \sum\_{i=1}^n x\_i a\_i^r(\delta, \pi) + \alpha^E(\delta, \pi, \pi) \tag{3}$$

where *n* is the number of components in the mixture, *xi* is the mole fraction of component *i* in the mixture, the superscripts "0", "r" and "E" denote the ideal-gas part, the residual part, and the excess of dimensionless Helmholtz free energy, respectively. The subscript *i* denotes the component *i*. α<sup>r</sup> *<sup>i</sup>* can be calculated by the EOSs of pure CO2, CH4 and N2 fluids [41–43], which are in terms of dimensionless Helmholtz free energy and are recommended as the standard EOSs of pure CO2, CH4 and N2 fluids by the National Institute of Standards and Technology (NIST).

*δ* and *τ* of Equation (3) are the reduced mixture density (*δ* = *ρ*/*ρc*) and the inverse reduced mixture temperature (*τ* = *Tc*/*T*), where *ρ* and *T* are the density and temperature of the mixture, and *ρ<sup>c</sup>* and *Tc* are the pseudo-critical density and temperature of the mixture. *ρ<sup>c</sup>* and *Tc* are defined by

$$\rho\_c = \left[\sum\_{i=1}^n \frac{\mathbb{x}\_i}{\rho\_{ci}} + \sum\_{i=1}^{n-1} \sum\_{j=i+1}^n \mathbb{x}\_i \mathbb{x}\_j \mathbb{\zeta}\_{ij}\right]^{-1} \tag{4}$$

$$T\_c = \sum\_{i=1}^{n} \mathbf{x}\_i T\_{ci} + \sum\_{i=1}^{n-1} \sum\_{j=i+1}^{n} \mathbf{x}\_i^{\mathcal{S}\_{ij}} \mathbf{x}\_j \varsigma\_{ij} \tag{5}$$

where *xj* is the mole fraction of component *j* in the mixture, *ρci* and *Tci* are the critical density and temperature of component *i*, respectively. The critical temperatures and densities of the three components considered in this study are listed in Table 1. *ζij*, *βij* and *ςij* are binary parameters associated with components *i* and *j*.

**Table 1.** Critical parameters of pure fluids.


*α<sup>E</sup>* of Equation (3) takes the general form developed by Lemmon and Jacobsen (LJ-1999 EOS) [35].

$$\alpha^{E}(\delta,\pi,\boldsymbol{x}) = \sum\_{i=1}^{n-1} \sum\_{j=i+1}^{n} \boldsymbol{x}\_{i}\boldsymbol{x}\_{j}\boldsymbol{F}\_{i\bar{j}} \sum\_{k=1}^{10} N\_{k} \delta^{d\_{k}} \boldsymbol{x}^{t\_{k}} \tag{6}$$

where *Fij* is a binary parameter associated with components *i* and *j*, *N*k, *d*<sup>k</sup> and *t*<sup>k</sup> are the general parameters independent of fluids, which can be obtained from Lemmon and Jacobsen [35]. There are four binary parameters (*ζij*, *βij*, *ςij* and *Fij*) in the above equations. Binary interaction parameters of CH4-CO2 and CO2-N2 have been determined in previous work [34,40]. In this work, the binary interaction parameters of CH4-N2 systems are obtained by the nonlinear regression of selected experimental data. Since experimental *PVTx* data available are much more than experimental VLE data, only part of new and accurate experimental *PVTx* data are selected to fit the parameters to make the *PVTx* and VLE data keep similar weights in the fitting. In this work, the Levenberg–Marquardt algorithm of the nonlinear least squares method is used for fitting, which is an efficient and widely used mathematical optimization technique. The fitting condition is to minimize the weighted sum of squares of the calculation errors of the equation on the selected experimental *PVTx* and VLE data. Compared to Newton's method (e.g., The Newton– Raphson method used in REFPROP [44]), it combines the advantages of two numerical minimization algorithms: the gradient descent method and the Gauss–Newton method, and it has good convergence and robustness. Regressed parameters of the CH4-N2 system are listed in Table 2, which also includes the parameters of CH4-CO2 and CO2-N2 systems from previous studies.

**Table 2.** Binary interaction parameters of the CO2-CH4-N2 system.


Note: Subscripts i and j refer to the first component and the second component, respectively.

## *2.1. The Binary CH4-N2 Mixture*

The calculated deviations of the *PVTx* and VLE properties from each experimental data set [45–66] by this EOS for the CH4-N2 fluid mixture are listed in Table 3. Detail calculation method for the *PVTx* and VLE properties can be found in Mao et al. [40]. The deviations of the *PVTx* data are the density deviations, and the temperature and pressure of the experimental *PVTx* data are up to 673.15 K and 1380 bar, respectively. The deviations of VLE property are the deviations of the vapor-liquid phase equilibrium composition, and the experimental VLE composition covers the full range. The average absolute deviation of all *PVTx* data from this EOS is 0.38% for the CH4-N2 fluid mixture.



Note: Deviations of the PVTx data are the density deviations, and the deviations of VLE data are the equilibrium composition deviations. AAD is the abbreviation of average absolute deviation.

Figure 1 shows the deviations between this EOS and the experimental density data of the CH4-N2 mixture [49,52,54,64] at different pressures. It can be seen in Figure 1 that most of the density deviations are within ±1%, with or close to experimental accuracy. The comparisons between experimental density and the calculated densities by this EOS are shown in Figure 2, where the calculated densities from the PR and SRK EOSs [67] are also included. The errors of three EOSs are also plotted in Figure 2. Figure 2 shows the PR and SRK EOSs cannot well reproduce the *PVTx* properties, especially under high-pressure conditions. In contrast, the calculated densities from this EOS are in good agreement with experiment data [68].

Figure 3 compares the vapor-liquid phase equilibrium curves calculated from this EOS with experimental data [55,57,59,66] at three temperatures of 170 K, 115 K and 110 K. It can be seen from Figure 3 that the calculated vapor-liquid phase equilibrium curves from this EOS are consistent with experimental data, indicating that this EOS can accurately calculate the vapor-liquid phase equilibrium of CH4-N2 mixture.

**Figure 1.** Density deviations for the CH4-N2 system at different pressures. Experimental data(triangles [64], stars [52], rounds [54], squares [49]): *P* is the pressure, *ρ*cal is the calculated density, and *ρ*exp is the experimental density.

**Figure 2.** Experimental and calculated densities for the 0.8 mole CH4 + 0.2 mole N2 mixture at different temperatures. The experimental density data (rounds) [68], the PR (green dash lines) and SRK (blue dash lines) EOSs [67], this work (red lines): *P* is the pressure, *ρ* is the density, and *T* is temperature.

**Figure 3.** Vapor-liquid phase equilibria of the CH4-N2 system at different temperatures, (**a**) *P*-*x*N2/*y*N2 figure at 170 K; (**b**) *P*-*x*CH4/*y*CH4 figure for 170 K; (**c**) *P*-*x*N2/*y*N2 figure for 115 K; (**d**) *P*-*x*CH4/ *y*CH4 figure for 115 K; (**e**) *P*-*x*N2/*y*N2 figure for 110 K; (**f**) *P*-*x*CH4/*y*CH4 figure for 110 K. Experimental data (rounds [55], triangles [57], squares [59]), this work (solid lines): *x*N2 and *y*N2 are the mole fractions of N2 in the liquid and vapor phases, respectively *x*CH4 and *y*CH4 are the mole fractions of CH4 in the liquid and vapor phases, respectively.

## *2.2. The Ternary CO2-CH4-N2 Mixture*

Combining binary interaction parameters of the CH4-CO2, CO2-N2 and CH4-N2 systems, the EOS can also predict the thermodynamic properties of the CO2-CH4-N2 mixture. Here, the experimental *PVTx* and vapor-liquid equilibrium data are used to verify the accuracy of this EOS for the ternary mixture.

The calculated *PVTx* deviations for the ternary CO2-CH4-N2 mixture from each experimental data set [1,64,69–71] are given in Table 4, where the temperature and pressure of the *PVTx* data are up to 573.15 K and 1000 bar, and the composition almost covers the full range. The average absolute deviation of all *PVTx* data from this EOS is 0.40% for the CO2-CH4-N2 fluid mixture. Figure 4 shows the comparison between experimental data [71] and calculated densities for the CO2-CH4-N2 mixture at different temperatures. It can be seen from Figure 4 that calculated densities from this EOS are in good agreement with experiment data, indicating that this EOS has a good predictive ability for the CO2-CH4-N2 mixture.

**Figure 4.** Comparisons between experimental and calculated densities for the CO2-CH4-N2 mixture at different temperatures. (**a**) *X*CO2 = 0.0998, *X*CH4 = 0.6625, *X*N2 = 0.2477; (**b**) *X*CO2 = 0.9949, *X*CH4 = 0.02, *X*N2 = 0.0301; (**c**)*X*CO2 = 0.1504, *X*CH4 = 0.1004, *X*N2 = 0.7492; (**d**) *X*CO2 = 0.3346, *X*CH4 = 0.3319, *X*N2 = 0.3335. Experimental data (squares [71]), this work (solid lines): *X*CO2 is the mole fraction of CO2 in the mixture, *X*CH4 is the mole fraction of CH4 in the mixture, *X*N2 is the mole fraction of N2 in the mixture.


**Table 4.** Calculated *PVTx* deviations for the ternary CO2-CH4-N2 mixture.

Note: Deviations of the *PVTx* data are the density deviations. *N* is the Number of data points, AAD is the abbreviation of average absolute deviation.

The vapor-liquid phase equilibrium properties of the ternary CO2-CH4-N2 system are shown in Figure 5, where the curves are calculated from this EOS, and the experimental data are from the literature of [72–74]. As can be seen from Figure 5, all the experimental data points in the non-critical region agree well with this EOS, but in the near-critical region (Figure 5c), the calculated values deviate largely from experimental data [72–74].

**Figure 5.** Vapor-liquid phase equilibria for the CO2-CH4-N2 system at different temperatures and pressures, (**a**) *T* = 220 K, *P* = 60.8 bar; (**b**) *T* = 230 K, *P* = 62.05 bar; (**c**) *T* = 230 K, *P* = 86.19 bar; (**d**) *T* = 270 K, *P* = 45.6 bar. Experimental data (squares [73], rounds [72], triangles [74]), this work (solid lines).

#### **3. The Applications of the ZMS EOS**

## *3.1. Calculating Excess Enthalpies*

The excess enthalpy (*HE*) is an important thermodynamic property of the mixture for the mixing and separation processes, which is defined by

$$H^E = \left\{ H\_m - \sum\_i x\_i H\_i \right\}\_{P,T} \tag{7}$$

where *Hm* is the enthalpy of the mixture, *Hi* is the enthalpy of pure component *i*, and *xi* is the mole fraction of component *i* in the mixture. According to the thermodynamic relationship between Helmholtz free energy and enthalpy, *Hm* and *Hi* are given by

$$H\_{\rm m} = RT\left(1 + \pi\left(\alpha\_{\tau}^{0} + \alpha\_{\tau}^{r}\right) + \delta\alpha\_{\delta}^{r}\right) \tag{8}$$

$$H\_{\dot{l}} = RT\left(1 + \tau\_{\dot{l}}\left(\alpha\_{\tau \dot{i}}^{0} + \alpha\_{\tau \dot{i}}^{r}\right) + \delta\_{\dot{l}}\alpha\_{\delta \dot{i}}^{r}\right) \tag{9}$$

where *α*<sup>0</sup> *<sup>τ</sup>* = *∂α*<sup>0</sup> *∂τ δ* , *α<sup>r</sup> <sup>τ</sup>* = *∂α<sup>r</sup> ∂τ δ* , *α<sup>r</sup> <sup>δ</sup>* = *∂α<sup>r</sup> ∂δ τ* , *α*<sup>0</sup> *<sup>τ</sup><sup>i</sup>* = *∂α*<sup>0</sup> *i ∂τ<sup>i</sup> δi* , *α<sup>r</sup> <sup>τ</sup><sup>i</sup>* = *∂α<sup>r</sup> i ∂τ<sup>i</sup> δi* and *αr <sup>δ</sup><sup>i</sup>* = *∂α<sup>r</sup> i ∂δ<sup>i</sup> τi* . *δ<sup>i</sup>* and *τ<sup>i</sup>* are the reduced density and inverse reduced temperature of pure component i, which are defined by

$$
\delta\_i = \rho\_i / \rho\_{ci} \tag{10}
$$

$$
\pi\_{\bar{l}} = T\_{c\bar{l}} / T \tag{11}
$$

where *ρ<sup>i</sup>* is the density of component *i*. Based on Equations (8)–(11), *Hm* can be calculated by this EOS of the CO2-CH4-N2 fluid mixtures and *Hi* can be calculated by the abovementioned equations of pure CO2, CH4 and N2 fluids.

The calculated excess enthalpy curves of the CH4-CO2 and CO2-N2 mixtures are, respectively, shown in Figures 6 and 7, where the experimental data [75,76] are included for comparison. The following EOSs are also included for comparison: the standard Peng– Robinson (PR) EOS, either optimized (optimal *kij*) or not (*kij* = 0); the PR EOS with the residual part of the Wilson excess Helmholtz energy model (PR + EOS/*aE*,*Wilson res* ) [19]. It can be seen from Figures 6 and 7 that this EOS is more accurate than the standard PR EOS, either optimized (optimal *kij*) or not (*kij* = 0) at most cases. In general, the PR + EOS/*aE*,*Wilson res* and this EOS shows a similar predictive capability for the excess enthalpies of the CH4-CO2 mixture. However, this EOS is slightly more accurate than the PR + EOS/*aE*,*Wilson res* for the CO2-N2 mixture at low pressures.

Figure 8 compares the excess enthalpy curves of the CH4-N2 mixture calculated from the ZMS EOS with experimental data [77], and enthalpies calculated from the Enhanced-Predictive-PR78 (E-PPR78) EOS [20] are also added for comparison. Figure 8a shows the ZMSEOS are significantly more accurate than the (E-PPR78) EOS at high and middle pressures. Figure 8b shows the enthalpies of mixing calculated from the E-PPR78 EOS are slightly better than those of the ZMS EOS at low and middle pressures, but the enthalpies of mixing calculated from the ZMS EOS are better than those from the E-PPR78 EOS at *P* = 101.33 bar.

**Figure 6.** Comparisons between experimental and calculated mixing enthalpies for the CH4-CO2 mixture, (**a**) *T* = 283.15 K, *P* = 10.13bar; (**b**) *T* = 283.15 K, *P* = 30.40 bar; (**c**) *T* = 313.15 K, *P* = 10.13 bar; (**d**) *T* = 313.15 K, *P* = 30.40 bar. Experimental data (squares [75]). PR*(kij* = 0) (red dash curves); PR (optimal *k*ij) (green dash curves); PR + EOS/*aE*,*Wilson res* (blue dash curves); This work (black solid curves): *h<sup>E</sup>* is the mixing enthalpy.

**Figure 7.** *Cont*.

**Figure 7.** Comparisons between experimental and calculated mixing enthalpies for the CO2-N2 mixture, (**a**) *T* = 313.15 K, *P* = 10.13bar; (**b**) *T* = 313.15 K, *P* = 30.40 bar; (**c**) *T* = 313.15 K, *P* = 81.06 bar; (**d**) *T* = 313.15 K, *P* = 121.59 bar. Experimental data (black squares [76]). PR (*k*ij = 0) (red dash curves); PR (optimal *<sup>k</sup>*ij) (green dash curves); PR + EOS/*aE*,*Wilson res* (blue dash curves); This work (black solid curves).

**Figure 8.** Comparisons between experimental and calculated mixing enthalpies for the CH4-N2 mixture at 195.15 K and 253.15 K. (**a**) *T* = 195.15 K; (**b**) *T* = 253.15 K. Experimental data (squares [77]), this work (solid curves); the Enhanced-Predictive-PR78 (E-PPR78) EOS (dash curves).

## *3.2. Calculating the Solubility of CO2-CH4-N2 Mixtures in Aqueous Electrolyte Solution*

The solubility of the CO2-CH4-N2 gas mixtures in the electrolyte solutions can provide the quantitative assessment for the storage of CO2 in deep saline aquifers. According to the thermodynamic principle, when the CO2-CH4-N2 gas mixtures reach dissolution equilibrium in the electrolyte solution, the chemical potential of component *i* in the liquid phase (*μ<sup>L</sup> <sup>i</sup>* ) and its chemical potential in the vapor phase (*μ<sup>V</sup> <sup>i</sup>* ) are equal. The chemical potential can be written in terms of fugacity in the vapor phase and activity in the liquid phase

$$
\mu\_i^L = \mu\_i^{L(0)} + RT\ln m\_i \gamma\_i \tag{12}
$$

$$
\mu\_i^V = \mu\_i^{V(0)} + RT \ln \eta\_i y\_i P \tag{13}
$$

where *i* is the gas component in vapor mixtures, and *P* is the pressure in bar. *yi* is the molar fraction of *i* in the vapor phase, *ϕ<sup>i</sup>* is the fugacity coefficient of *i* in the vapor phase and *mi* is the solubility of *i* in the liquid phase in mol/kg. *γ<sup>i</sup>* is the activity coefficient of component *i* in the liquid phase. *μL*(0) *<sup>i</sup>* and *<sup>μ</sup>V*(0) *<sup>i</sup>* are the standard chemical potential of *i* in the liquid and vapor phase, respectively.

At the dissolution phase equilibrium, *μ<sup>V</sup> <sup>i</sup>* = *<sup>μ</sup><sup>L</sup> <sup>i</sup>* . From Equations (12) and (13), we can obtain

$$
\ln m\_i = \ln y\_i P + \ln \rho\_i - \frac{\mu\_i^{L(0)} - \mu\_i^{V(0)}}{RT} - \ln \gamma\_i \tag{14}
$$

Since the water content in the vapor phase at equilibrium is generally small, it has a negligible effect on the fugacity coefficient. When calculating the fugacity coefficient, the water content in the vapor phase can be ignored, the vapor phase can be approximated as the CO2-CH4-N2 system and the mole fraction of gas component *i* in the vapor phase is expressed as *yi*. Consequently, *ϕ<sup>i</sup>* can be calculated from the ZMS EOS developed in this work.

Pitzer activity coefficient model [78] is chosen to calculate *γ<sup>i</sup>*

$$
\ln \gamma\_i = \sum\_{\mathfrak{c}} 2\lambda\_{i-\mathfrak{c}} m\_{\mathfrak{c}} + \sum\_{\mathfrak{a}} 2\lambda\_{i-\mathfrak{a}} m\_{\mathfrak{a}} + \sum\_{\mathfrak{c}} \sum\_{\mathfrak{a}} \zeta\_{i-\mathfrak{c}-\mathfrak{a}} m\_{\mathfrak{c}} m\_{\mathfrak{a}} \tag{15}
$$

where λ and *ζ* are second-order and third-order interaction parameters, respectively; c and *a* refer to cation and anion, respectively. Substituting Equation (15) into Equation (14) yields

$$\begin{split} lm\boldsymbol{m}\_{i} &= \quad lm\boldsymbol{y}\_{i}\boldsymbol{P} + lm\boldsymbol{\varrho}\_{i} - \frac{\mu\_{i}^{L(0)} - \mu\_{i}^{V(0)}}{R\boldsymbol{I}} - \sum\_{c} 2\lambda\_{i-c}m\_{c} \\ &- \sum\_{a} 2\lambda\_{i-a}m\_{a} - \sum\_{c} \sum\_{a} \zeta\_{i-c-a}m\_{c}m\_{a}a \end{split} \tag{16}$$

As can be seen from Equation (17), *mi* is related to the difference between *<sup>μ</sup>L*(0) *<sup>i</sup>* and *μV*(0) *<sup>i</sup>* , and is not related to the specific value of *<sup>μ</sup>L*(0) *<sup>i</sup>* or *<sup>μ</sup>V*(0) *<sup>i</sup>* . Therefore, to simplify the model, *μV*(0) *<sup>i</sup>* is assumed to be zero. Equation (16) can be simplified as

$$\begin{split} l\ln m\_{i} &= -\ln y\_{i}P + \ln \phi\_{i} - \frac{\mu\_{i}^{L(0)}}{RT} - \sum\_{c} 2\lambda\_{i-c}m\_{c} \\ &- \sum\_{a} 2\lambda\_{i-a}m\_{a} - \sum\_{c} \sum\_{a} \zeta\_{i-c-a}m\_{c}m\_{a} \end{split} \tag{17}$$

where *<sup>μ</sup>L*(0) *i RT* , *λ* and *ζ* are the function of temperature and pressure, and the parameters can be determined by the solubility model of pure gases (CO2, CH4, N2). Since the 1980s, many solubility models for CO2, CH4 and N2 gases in pure water and brine have been developed [79–91], most of which do not exceed 473 K and 1000 bar. Among them, the solubility models of CH4, CO2 and N2 established by Mao and co-workers [83,86,91] are applicable to a wider range of temperature, pressure, and salinity with higher accuracy. They have been chosen as the solubility models of pure CH4, CO2 and N2 gases in this work. Table 5 lists the applicable temperature, pressure, and salinity ranges of these solubility models for pure gases.

**Table 5.** Ranges of pure gas solubility models.


To verify the accuracy of the predictions, the experimental data for the solubility of the CO2-CH4-N2 gas mixtures are compared with the calculated results. Figure 9 compares the experimental data on the solubility of the CH4-CO2 mixture in pure water [92,93] for the temperature range of 324.5–375.5 K and the pressure range of 100–750 bar, and it can be seen that the agreement is very good.

The solubility experimental data of the CO2-N2 gas mixture in pure water [94] and the calculated results are compared in Figure 10, where the temperature range is 308.15–318.15 K and the pressure range is 80–160 bar. As can be seen from Figure 10 the predictive results are in good agreement with the experimental data. The average absolute deviations between the calculated N2 and CO2 solubility of this model and the experimental data are 6.04% and 1.52%, respectively. Figure 11 compares the solubility

**Figure 9.** Comparisons between experimental and calculated solubilities for the CH4-CO2 mixture in water. (**a**) *mCH*<sup>4</sup> -*yCH*<sup>4</sup> figure at 344.15 K; (**b**) *mCO*<sup>2</sup> - *yCH*<sup>4</sup> figure at 344.15 K;(**c**) *mCH*<sup>4</sup> -*yCH*<sup>4</sup> figure at 373.5 K; (**d**) *mCO*<sup>2</sup> - *yCH*<sup>4</sup> figure at 375.5 K; (**e**) *mCH*<sup>4</sup> -*yCH*<sup>4</sup> figure at 324.5 K; (**f**) *mCO*<sup>2</sup> - *yCH*<sup>4</sup> figure at 324.5 K. Experimental data (rounds [92], squares [93]), this work (solid curves): *mCO*<sup>2</sup> is the solubility of CO2, *mCH*<sup>4</sup> is the solubility of CH4.

**Figure 10.** Comparisons between experimental and calculated solubilities for the CO2-N2 mixture in water. (**a**) *mN*<sup>2</sup> - *yN*<sup>2</sup> figure at 308.15 K; (**b**) *m*CO2 - *yN*<sup>2</sup> figure at 308.15 K; (**c**) *mN*<sup>2</sup> - *yN*<sup>2</sup> figure at 318.15 K; (**d**) *m*CO2 - *yN*<sup>2</sup> figure at 318.15 K. Experimental data (squares [94]), this work (solid curves): *mN*<sup>2</sup> is the solubility of N2.

**Figure 11.** Comparisons between experimental and calculated solubilities for the CO2-N2 mixture in saline water. (**a**) *mN*<sup>2</sup> - *yN*<sup>2</sup> figure; (**b**) *m*CO2 - *yN*<sup>2</sup> figrure. Experimental data (black squares [94]), this work (solid curves).

## *3.3. The Impact of Impurities (CH4 and N2) on the CO2 Storage Capacity*

CH4 and N2 are non-condensable impurities, which reduce the CO2 storage capacity in geological formations. Based on the EOS of the CO2-CH4-N2 fluid mixtures, the impact of CH4 and N2 on CO2 storage capacity can be calculated quantitatively by the normalized storage capacity proposed by Wang et al. [9]. For a certain reservoir with a certain volume, the normalized storage capacity of the CO2-CH4-N2 gas mixtures is shown as

$$\frac{A}{A\_{\rm CO\_2}} = \frac{\rho}{\rho\_{\rm CO\_2}(1 + \frac{M\_{\rm CH\_4}}{M\_{\rm CO\_2}} + \frac{M\_{N\_2}}{M\_{\rm CO\_2}})} \tag{18}$$

where *A* and *A*CO2 are the mass of CO2 in the mixtures and the mass of pure CO2 under the same volume, respectively. *ρ*CO2 and *ρ* are the density of the pure CO2 and CO2-CH4-N2 gas mixtures, respectively; *M*CH4 , *MN*<sup>2</sup> and *M*CO2 are the mass of CH4, N2, and CO2 in the mixture, respectively. *A*/*A*CO2 is the ratio of the mass of CO2 per unit volume in the mixture to that in the pure state. The value of *A*/*A*CO2 can be viewed as the normalized storage capacity for CO2 (i.e., the storage capacity for structural trapping of CO2). For pure CO2, the normalized storage capacity (*A*/*A*CO2 ) is 1. For the CO2-CH4-N2 gas mixtures of a given composition, *M*CH4 , *MN*<sup>2</sup> and *M*CO2 are also known. *ρ*CO2 can be calculated by the above-mentioned equation of pure CO2 fluid and *ρ* can be calculated by the ZMS EOS.

The normalized storage capacity of the CO2-CH4-N2 fluid mixture calculated by the ZMS EOS is plotted in Figure 12. Figure 12a shows the normalized storage capacity of a given composition at different temperatures. With the increase in pressure, the normalized storage capacity decreases at first and then increases, and there is a minimum point. When the temperature changes, the pressure corresponding to the lowest point of the normalized storage also changes. Figure 12b shows the normalized storage capacity of different compositions at the same temperature, which indicates that the content of impurities is much larger and the normalized storage capacity is much smaller.

**Figure 12.** Normalized CO2 storage capacity calculated by this EOS of CO2-CH4-N2 mixture. (**a**) the Normalized capacity-*P* figure at different temperatures for the same composition; (**b**) the Normalized capacity-*P* figure at differents compositions for the same temperatures.

## *3.4. Isochores of the CO2-CH4-N2 fluid Inclusions*

In the studies of fluid inclusions, the isochores (pressure-temperature relation at constant density and composition) are frequently used to estimate the trapping temperatures and pressures.

Based on the EOS of the CO2-CH4-N2 fluid mixtures, isochores of the CO2-CH4-N2 inclusions can be calculated by the following equation:

$$P(\delta, \tau, \mathbf{x}) = \rho RT[1 + \delta a\_{\delta}^{r}] \tag{19}$$

The calculated isochores of the CO2-CH4-N2 inclusions at two different compositions are plotted in Figure 13, from which it can be seen that the isochores of the CO2-CH4-N2 inclusions are a bit curved.

**Figure 13.** Calculated isochores of the CO2-CH4-N2 mixture by this EOS, (**a**) *X*CO2 = 0.8, *X*CH4 = 0.1, *X*N2 = 0.1; (**b**) *X*CO2 = 0.7, *X*CH4 = 0.15, *X*N2 = 0.15: *V*<sup>m</sup> is the molar volume.

This section may be divided into subheadings. It should provide a concise and precise description of the experimental results, their interpretation, as well as the experimental conclusions that can be drawn.

## **4. Conclusions**

A fundamental EOS for the Helmholtz free energy of the CH4-N2 mixture has been developed by using four binary interaction parameters. Comparisons with experimental *PVTx* and VLE data available showed that the EOS can satisfactorily reproduce the experimental volumetric and vapor-liquid phase equilibria data of binary CH4-N2 mixtures up to 673.15 K and 1380 bar, with or close to experimental accuracy. Combining this EOS of the CH4-N2 fluid mixtures and the EOS of the CH4-CO2 and CO2-N2 fluid mixtures developed in our previous work, an EOS of ternary CO2-CH4-N2 fluid mixtures has been presented, which is named ZMS EOS. The ZMS EOS can be applied to calculate excess enthalpies, predict the solubility of the CO2-CH4-N2 gas mixtures in water and brine, estimate the impurities CH4 and N2 on the CO2 storage capacity, and calculate the isochores of the CO2-CH4-N2 fluid mixtures.

**Author Contributions:** J.Z.: Conceptualization, Software, Investigation, Writing—original. S.M.: Supervision, Methodology, Writing—review, and editing. Z.S.: Investigation, Software, Validation. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research was funded by the National Natural Science Foundation of China (Grant No. 42073056) and Everest Scientific Research Program of Chengdu University of Technology (80000- 2021ZF11419).

**Institutional Review Board Statement:** Not applicable.

**Informed Consent Statement:** Not applicable.

**Data Availability Statement:** Not applicable.

**Conflicts of Interest:** The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

## **References**


**Disclaimer/Publisher's Note:** The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.
