CO2 transports in the Earth’s interior play a crucial role in understanding the deep carbon cycle and the global climate changes. Currently, CO2 transports inside of the Earth under extreme condition of pressure and temperature have not been understood well. In this study, the molecular dynamics (MD) calculations were performed to study CO2 transports under different CO2 pressures in slit-like magnesite pores with different pore sizes at 350~2500 K and 3~50 GPa are presented. Diffusion of CO2 in magnesite was improved as the temperature increases but showed the different features as a function of pressure. The diffusion coefficients of CO2 in magnesite were found in the range of . Magnesite with the pore size of 20~25 Å corresponds to the highest transports. Anisotropic diffusion of CO2 in magnesite may help to understand the inhomogeneous distribution of carbon in the upper mantle. The time of CO2 diffusion from the mantle to Earth surface was estimated to be around several tens of Ma and has an important effect on deep carbon cycle. The simulation of CO2 transports based on the Earth condition provides new insights to revealing the deep carbon cycle in the Earth’s interiors.

1. Introduction

Prediction of the changes in the Earth’s climate due to the increase of carbon dioxide (CO2) concentrations needs to deeply understand the chemical and physical processes that control the carbon cycle [1]. Generally, the carbon cycle involves the deep and shallow subcycles of carbon. So far, intensive attention has been paid to the shallow carbon cycle that is mainly related to the atmosphere, oceans, and shallow crustal. However, there is relatively little known about the deep carbon cycle. Due to the Earth’s deep interior that contains more than 90% of the total Earth’s carbon [2, 3], the carbon transport and storage in the deep Earth play a key role in determining the global carbon cycle [1, 49]. Therefore, it is very necessary to understand the carbon cycle in the Earth’s deep interior.

Carbon in the mantle exists in the different oxidized forms (e.g., CO, CO2, and CO32−) and the different states (solid, molten carbonate, fluid, etc.) [1013]. Magnesite (MgCO3) is the most common carbonate phase in the subducted slabs [14]. High-pressure experiments suggest that magnesite with a rhombohedral structure can be stable up to 115 GPa at 2000~3000 K which may be the major carbon host in the mid-lower mantle [1417].

CO2 is an important carbon-bearing form in the Earth’s interior. About 40 megatons of CO2 emits from the mantle to the surface every year by diffuse outgassing [18]. Therefore, CO2 plays a fundamental role in determining the global carbon cycle. In situ Raman spectroscopy in laser-heated diamond anvil cells revealed that CO2 molecular was stable up to 34 GPa at ~1770 K [6, 19]. Therefore, the transport properties of CO2 in the deep Earth under high pressure and high temperature are a vital link to understand the deep carbon cycle.

Except for high-pressure and high-temperature (high P-T) experiments, molecular dynamics (MD) simulations are powerful tools to investigate the molecular-scale dynamic properties of materials ([2026]; Hu et al. 018). There are some achievements on researching diffusion of CO2 in different materials by MD simulation [27]. For example, the CO2 diffusion coefficients of crystalline polystyrene at 1 atm. are (at 333.15 K), (at 343.15 K), and (at 353.15 K) [28]. Diffusion coefficients of supercritical CO2 and its mixtures monotonously decreased as its increasing molar density at 317.5 K [29]. The diffusion coefficients of H2, CH4, CO, O2, and CO2 in water near the critical point (600~670 K, 250 atm) were investigated by MD simulation [30], and CO2 can greatly influence the diffusion of mobile species in clay minerals [31]. However, most of those works focused on model structures, which are different from the potential minerals in the mantle. In addition, previous works were mainly performed under the static conditions without concerning high temperature and high pressure. Accordingly, the transport properties of CO2 under high temperature-pressure in the Earth’s interior conditions have not yet been identified. In this work, the diffusion of CO2 in magnesite was investigated under the high temperature-pressure conditions mimicking the extreme environment in the mantle.

2. Methods

2.1. Crystalline Models

Magnesite has a space group of R3c. Its structure is composed of the alternative layers of Mg and CO3 groups. Among them, one C atom and three neighbor O atoms lie in the same plane that is parallel to the - plane. Under the ambient condition, the cell parameters of magnesite are as follows: , , , [32]. A supercell of magnesite was established. A vacuum layer was added between two CO3 planes to build slit-like pores. The pore size (H, Å) was determined by the height of a vacuum layer. The two-layer slit-like pore structures of magnesite with electrically neutral were built and simulated (Figure 1).

2.2. Equilibrium Adsorption Configurations: The Grand Canonical Monte Carlo Method (GCMC) Simulation

The GCMC employed the μVT ensemble. The chemical potential , the volume (), and the temperature () were kept constant, while the number of molecules was allowed to fluctuate [33]. GCMC simulations were carried out by sorption code as implemented in Materials Studio 2017. The structures of CO2 in magnesite pores with energy minimized adsorption were obtained under the different conditions. During the GCMC simulation, CO2 molecules within the framework were randomly rotated and translated. Only were the positions and orientations of CO2 molecules sampled. Each configuration was treated as a rigid body. The COMPASS force field was used to describe the hypersurface of the potential energy on which the atomic nuclei moved [34]. The Ewald method was adopted to calculate the Coulomb force interaction. The atom interaction-based method was used to calculate van der Waals force with the Lennard-Jones potential. The cut-off distance was 0.15 Å. The number of load step in each calculated run is , in which the balance step number is and the process step number is . This GCMC method was successfully used to calculate the absorption of CO2 in magnesite in our previous works [35], and its validity also was proved by different works [31, 3640].

2.3. Molecular Dynamics (MD) Simulations

MD simulations were carried out by the Forcite code as implemented in Materials Studio 2017. The equilibrium configurations from GCMC calculations were further used in MD simulations. The MD simulations were performed with a time step of 1 fs (up to 1000 ps) using an NPT ensemble at 350~2500 K and 3~50 GPa. Interlayer molecular configurations were sampled every 1000 steps for further analysis. The temperature was kept constant using a velocity scale thermostat algorithm, and the convergence criterion is 1 K. The Ewald summation method was employed with an Ewald accuracy of  kcal/mol. Like the GCMC simulations, the COMPASS force field was used for describing the potential energy hypersurface [34].

2.4. Calculations for Diffusivities

The mean square displacement (MSD) is a mode of the displacement of particles followed over time. It can be defined aswhere is the position of the atom at time , and the angular bracket denotes that the quantity is an ensemble average property. In an equilibrium ensemble, the average squared displacement is independent of the time , which can be averaged out, leaving the mean square displacement over an interval :

If there are equivalent particles, the MSD can be further averaged as

If the particle is diffusing, the MSD becomes linear as a function of time (diffusive regime). The self-diffusion coefficient () can be described by the slope according to Einstein equation:where is the number of diffusion molecules, is the time, and is the displacement vector of the th molecule at time . The angular bracket denotes that the quantity is an ensemble average property.

3. Results

3.1. The Diffusivities at High Pressure and High Temperature

Diffusions of CO2 in magnesite pore were calculated with the pore size of 25 Å and the CO2 pressure of 50 MPa at high pressures (3~50 GPa) and high temperatures (350~2500 K); the results are shown in Figure 2. Diffusion of CO2 in magnesite increased as temperature increases but decreased as pressure increases. Specifically, the ( in 10-12 m2/s) increased from 1.67 at 350 K to 4.45 at 2500 K (3 GPa), while it decreased from 4.45 at 3 GPa to 3.29 at 50 GPa (2500 K). However, the change trend of the diffusion with temperature is different under different pressures. At low pressure (<9 GPa), the diffusion coefficient consistently increased with increasing the temperature at the whole range explored in this work. However, as the pressure increased, the diffusion coefficients showed the different trends when the temperature is higher or lower than 900 K.

For the temperature higher than 900 K, the values showed a good linear relation with the , which corresponds well with the Arrhenius law:where is the preexponential factor, is the gas constant, and is the activation energy for diffusion. Fitting parameters at 900~2500 K and 3~50 GPa are shown in Figure 3 and Table 1. The decreased from 4.929 kJ/mol at 3 GPa to 4.501 kJ/mol at 50 GPa. The increased from 25.2 kJ/mol at 3 GPa to 60.7 kJ/mol at 50 GPa, indicating that diffusion of CO2 in magnesite becomes more difficult under high pressure.

The different diffusion properties of CO2 in magnesite under different pressure and temperature (9 GPa and 900 K) might be related to the phase transition of CO2 and stability of magnesite under different pressure and temperature. The CO2 might transform from fluid phase to Pa3 molecular structure, Pbcn, P42/mnm, and I42d with increasing pressure [4144]; however, the stability of CO2 under high temperature is not yet well understood. For example, at ambient temperature, Pa3 CO2 undergoes a pressure-induced transformation to the Cmca symmetry at about 9-12 GPa; when pressure is larger than 12 GPa, there are two molecular phases: P42/mnm phase for  K and Pbcn phase for  K [43, 44]. At the same time, the presence of CO2 might affect the stability of magnesitite [45, 46].

3.2. Diffusion Anisotropy

Figure 4 shows an anisotropic feature of CO2 diffusion in magnesite. For the magnesite structure, it has the same -axis and -axis; therefore, the diffusion in the - plane of magnesite is the same. Accordingly, we only calculated diffusion coefficient along -axis and -axis direction under high T-P. The diffusion coefficients () of CO2 along the -axis and -axis varied differently with temperature and pressure (Figure 4). At temperature of 350~500 K, the ( along -axis direction) was larger than the ( along -axis direction) at a low pressure (<9 GPa) and then the became smaller than the with increasing pressure. When temperature was higher than 700 K, the was larger than the , but the differences became very small with increasing pressure. Anisotropy of the CO2 diffusion in magnesite decreased with increasing pressure and temperature. At 25 GPa, anisotropy increased from 3% at 2500 K to 60% at 350 K. For comparison, at 3 GPa, the anisotropy increased from 3.5% at 2500 K to 280% at 350 K. As the pressure was higher than 25 GPa, the difference between and was very small. The results indicated diffusion of CO2 in magnesite is anisotropic at a low pressure and is isotropic under high pressure.

Anisotropy of CO2 diffusion in magnesite is important to understand the anisotropic properties of deep carbon cycle. A study on olivine-hosted melt inclusions from the Mid-Atlantic Ridge indicated that carbon content in the upper mantle is highly heterogeneous, varying by almost two orders of magnitude globally [47]. Volatiles in the mantle can be hosted either in independent phases or as dissolved components in minerals. The content of magnesite in the upper mantle is 1000~8000 ppm [17, 48]. Under the upper mantle ( and ), diffusion of CO2 in magnesite shows an obvious anisotropy. Therefore, this diffusion difference of CO2 probably leads to the heterogeneous transport and storage of carbon in the mantle.

3.3. The Diffusion under Different CO2 Pressures

The diffusion of CO2 in magnesite is closely related to its content in the mantle. Therefore, values of CO2 in magnesite under different CO2 pressures were calculated in order to understand the effect of CO2 content on its diffusion behavior. Figure 5 shows the values of CO2 in magnesite with the CO2 pressures of 10, 30, 50, 70, and 100 MPa under the 900~2500 K and 3~50 GPa. The diffusibility increased as the CO2 pressure increased. The averages of values increased by 1.5%-4.8% as the CO2 pressure increased from 10 MPa to 100 MPa under the different pressures. As detected by Raman spectroscopy, the diffusion coefficient of CO2 in aqueous solutions decreases significantly at high CO2 concentrations [49]. Our results here show different diffusion behavior from the diffusivity of CO2 in aqueous solutions under geologic carbon sequestration conditions. In addition, the effects of CO2 pressure on its diffusion show a slight decrease with increasing the temperature. Compared with the change caused by temperature or pressure (orders of magnitude), the effect of CO2 pressure on its diffusion is limited.

3.4. The Diffusion under Different Magnesite Pore Sizes

Slit-like magnesite pores of 7.5, 10, 12.5, 15, 20, 25, and 30 Å were constructed. Figure 6 shows diffusions of CO2 in these models. It is clear that firstly increases and then decreases with increasing pore sizes. Under the same simulated condition of temperature and pressure, the diffusion has the highest ability in magnesite with a pore size of 20~25 Å. The differences of the diffusion caused by the change of pore size become smaller with increasing pressure.

4. Discussion

Plenty of CO2 is emitted from the deep Earth to the surface by diffusion [18, 47] and made significant influence on the Earth’s carbon cycle. Therefore, it is critical to understand the diffusion process of CO2 in the Earth’s interior. The diffusion coefficient of CO2 in magnesite is at 350~2500 K and 3~50 GPa. CO2 diffusion in crystalline polystyrene was at 1 atm and 333.15~353.15 K [28]. The diffusion coefficients of CO2 in silicate melt (Jadeite+albite+quartz) are at 1323~1673 K and 500 MPa [6]. Under the similar condition of temperature and pressure, our calculated diffusion coefficient of CO2 is smaller than that in crystalline polystyrene, however is larger than that of CO2 in silicate melt. The results indicate the faster diffusion ability of CO2 in magnesite crystal than in silicate melt.

Another critical issue on deep carbon cycle is the distribution of the diffusion ability of CO2 along the Earth’s geothermal profile. From Figure 7, we can estimate diffusion distance of CO2 in magnesite is 0.017~0.933 m each year. Accordingly, CO2 diffusion from the mantle to the surface takes around dozens of Ma. Carbon hosted in the mantle ultimately affected the surface carbon budgets over billion-year time scales [50, 51]. Similarly, isotopic results indicated three time-scale carbon cycles in the entire upper mantle: 1: >120 km, 5~10 Ma; 2: 120~410 km, 48~53 Ma; and 3: 410~660 km, >60 Ma [5256]. The diffusion of CO32- would be three to four orders of magnitude slower than that of molecular CO2 [6]. Therefore, CO2 diffusion plays the important role in determining the deep carbon cycles. For the first time, we present the CO2 diffusion behavior in real minerals (magnesite) under the Earth’s interior conditions. The results here provide new data to understand the deep carbon cycle.

The above results indicate that temperature, pressure, and pore size of magnesite have important effect on CO2 diffusion. The excess adsorption amounts of CO2 in the pores of magnesite decrease with increasing temperature but increase with increasing pressure and pore size [35]. So interior circumstance changes from deep to shallow of the Earth might drive the transport of CO2 in its interior, such as the changes of stress and porosity caused by the movement of fault or tectonic slab. However, the grain or mineral boundaries are always the important channel for fluid migration in the Earth’s interior, which also needs to be paid more attention.

5. Conclusion

For understanding the Earth’s deep carbon cycle, diffusion of CO2 in magnesite pore with different CO2 pressure and pore sizes was calculated under 350~2500 K and 3~50 GPa. The diffusion coefficient of CO2 in magnesite increases with increasing temperature but decreases with increasing pressure. The increases with increasing pressure, indicating the diffusion of CO2 in magnesite becomes more difficult under high pressure. CO2 has the strongest diffusion in the magnesite when the pore size is 20~25 Å, and the diffusion slightly increases with increasing CO2 pressure. Diffusion of CO2 in magnesite shows obvious anisotropy and might be responsible for the inhomogeneous distribution of carbon in the upper mantle. The diffusion coefficient of CO2 in magnesite is in the range of under the pressure and temperature investigated here. Our calculated results indicate that diffusion of CO2 is faster in magnesite than in silicate melt. The results of CO2 diffusion in magnesite under high pressure and high temperature provide new insights to understand the carbon cycle in the Earth’s interiors.

Data Availability

All data included in this study are available upon request by contact with the corresponding author.

Conflicts of Interest

The authors declare that they have no conflicts of interest.


This work was supported by the Foundation of the Key Laboratory of Earthquake Forecasting and United Laboratory of High-Pressure Physics and Earthquake Science, the Institute of Earthquake Forecasting, CEA (Grant No. 2019IEF0101-1, 2019IEF0502, 2017KLEP03, and 2019HPPES06), and the National Natural Science Foundation of China (Grant No. 41762008).