# Carbon Storage in Earth's Deep Interior Implied by Carbonate-Silicate-Iron Melt Miscibility

A. Davis, N. Solomatova, R. Caracas, A. Campbell

► To cite this version:

A. Davis, N. Solomatova, R. Caracas, A. Campbell. Carbon Storage in Earth's Deep Interior Implied by Carbonate-Silicate-Iron Melt Miscibility. Geochemistry, Geophysics, Geosystems, 2023, 24 (8), 10.1029/2023GC010896. https://hal-04284726

HAL Id: hal-04284726

https://hal.science/hal-04284726v1

Submitted on 17 Nov 2023

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
# Geochemistry, Geophysics, Geosystems®

## RESEARCH ARTICLE

10.1029/2023GC010896

**Key Points:**

* Carbonate-iron and carbonate-silicate melts have miscibility gaps that close with increasing pressure

* Carbonate-silicate-iron melts may contribute to the existence of ultra-low velocity zones at the core-mantle boundary

* Carbon's affinity for iron indicates that much of Earth's carbon could have been transported to the core during core-mantle differentiation

**Supporting Information:**

Supporting Information may be found in the online version of this article.

**Correspondence to:**

A. H. Davis,
a.h.davis@geo.uio.no

**Citation:**

Davis, A. H., Solomatova, N. V., Caracas, R., & Campbell, A. J. (2023). Carbon storage in Earth's deep interior implied by carbonate-silicate-iron melt miscibility. *Geochemistry, Geophysics, Geosystems*, 24, e2023GC010896. https://doi.org/10.1029/2023GC010896

Received 3 FEB 2023
Accepted 4 AUG 2023

© 2023 The Authors. Geochemistry, Geophysics, Geosystems published by Wiley Periodicals LLC on behalf of American Geophysical Union. This is an open access article under the terms of the Creative Commons Attribution License, which permits use, distribution and reproduction in any medium, provided the original work is properly cited.

<figure><img src="image_2.png" /></figure>

## Carbon Storage in Earth's Deep Interior Implied by Carbonate-Silicate-Iron Melt Miscibility

A. H. Davis¹,², ∗, N. V. Solomatova³, †, R. Caracas⁴,⁵ †, and A. J. Campbell¹

¹Department of the Geophysical Sciences, University of Chicago, Chicago, IL, USA, ²Now at The Centre for Planetary Habitability (PHAB), University of Oslo, Oslo, Norway, ³Ecole Normale Supérieure de Lyon, Laboratoire de Géologie de Lyon, Lyon, France, ⁴PHAB, University of Oslo, Oslo, Norway, ⁵Université de Paris, Institut de Physique du Globe de Paris, CNRS, Paris, France

**Abstract** Carbonate melts have been proposed to exist in the lower mantle, but their interaction with other lower mantle melt compositions is poorly understood. To understand miscibility in the carbonate-silicate-metal melt system, we simulate endmember, binary, and ternary melt mixtures and study how their Gibbs free energies of mixing evolve with pressure. We find that carbonate-metal and carbonate-silicate melts have miscibility gaps that close with increasing pressure, while silicate-metal melts are immiscible at all lower-mantle pressures. Extending this analysis to the core-mantle boundary, we suggest three miscible melt fields near the endmember carbonate, silicate, and iron melt compositions. Analysis of the densities of these miscible melt compositions indicates that some carbonate-rich and some silicate-rich melt compositions are gravitationally stable at the core-mantle boundary and could be candidate compositions to explain ultra-low velocity zones. Additionally, we evaluate the speciation of an example immiscible melt composition at various pressures throughout the mantle and identify reduced carbon species that would be expected to form in the melt. Our analysis reveals that a majority of Earth's carbon could have been transported to the core during core-mantle differentiation and that much of Earth's carbon may be stored in the deep interior today.

**Plain Language Summary** Understanding the storage and cycling of carbon in the Earth's deep interior improves our knowledge of the Earth's formation and evolution throughout geologic time. Carbon-bearing melts are candidate phases for carbon storage in the lower mantle and may react and mix with other melt phases at places such as the core-mantle boundary. The extent of mixing upon reaction is dependent on the thermodynamic properties of the mixture components and determines the possible range of compositions, structures, and densities of multicomponent melt mixtures. To determine possible melt compositions that may arise from mixtures of carbonate, silicate, and iron melts in the lower mantle and their physical and chemical properties, we performed molecular simulations to determine whether a mixture will separate or stay mixed. We find that carbonate melts mix with silicate and iron melts at all lower mantle pressures, and that mixing in the carbonate-silicate-iron melt system increases with pressure. Additionally, we find that certain melt mixtures have densities at the core-mantle boundary that make them candidate compositions to explain ultra-low velocity zones. Finally, we find that carbon has an affinity for iron that leads to the formation of carbide-like structures that may have allowed carbon to become sequestered in the Earth's core during core formation.

## 1. Introduction

Carbon plays a vital role at the Earth's surface in biological and atmospheric processes, but the role of carbon in the deep Earth is less well understood. Carbon and other volatiles in the lower mantle are thought to be remnants of an early magma ocean (Labrosse et al., 2007) or to derive from subducting slabs (Plank & Manning, 2019). However, reports vary on the carbon content of the Earth (Javoy, 1997; McDonough & Sun, 1995), the distribution of carbon between core and mantle (Dasgupta & Walker, 2008; Wood et al., 2013), the phase relations of carbon-bearing phases at depth (Merlini et al., 2012; Oganov et al., 2008), and the reactivity of carbon-bearing phases with the surrounding core and mantle. The existence and role of carbonates in the lower mantle are highly contested. Multiple studies show that carbonates undergo melting (Li et al., 2017), reduce to diamond or iron carbide (Rohrbach & Schmidt, 2011), or decarbonate (Drewitt et al., 2019) before they reach the lower mantle in subducting slabs. However, the stability of carbonate phases depends on a host of thermodynamic variables, including pressure, temperature, and oxygen fugacity. Reports of carbonate melt inclusions in deep Earth
<table><caption>**Table 1**
Melt Compositions Simulated</caption><thead><tr><th>Melt</th><th>Mg</th><th>Si</th><th>C</th><th>O</th><th>Fe</th><th>Total</th></tr></thead><tbody><tr><td>MgCO<sub>3</sub></td><td>24</td><td>0</td><td>24</td><td>72</td><td>0</td><td>120</td></tr><tr><td>MgSiO<sub>3</sub></td><td>24</td><td>24</td><td>0</td><td>72</td><td>0</td><td>120</td></tr><tr><td>Fe</td><td>0</td><td>0</td><td>0</td><td>0</td><td>108</td><td>108</td></tr><tr><td>Mg(C,Si)O<sub>3</sub></td><td>24</td><td>12</td><td>12</td><td>72</td><td>0</td><td>120</td></tr><tr><td>MgCO<sub>3</sub> + Fe</td><td>24</td><td>0</td><td>24</td><td>72</td><td>13</td><td>133</td></tr><tr><td>MgSiO<sub>3</sub> + Fe</td><td>24</td><td>24</td><td>0</td><td>72</td><td>13</td><td>133</td></tr><tr><td>Mg(C,Si)O<sub>3</sub> + Fe</td><td>24</td><td>12</td><td>12</td><td>72</td><td>13</td><td>133</td></tr></tbody></table>

Note. The numbers refer to the number of atoms included in the simulation.

diamonds (Korsakov & Hermann, 2006) as well as petrologic experiments on solid carbonates in lower mantle phase assemblages (Dorfman et al., 2018; Lv et al., 2021) indicate that carbonates and carbonate melts may be stable and present in the lower mantle, plausibly even in the core-mantle boundary region. However, few studies have examined carbonate melt interactions in the lower mantle.

Previous ab initio studies have examined carbonate melts (Koura et al., 1996; Li et al., 2017; Xu et al., 2020), carbon-bearing silicate melts (Bajgain & Mookherjee, 2021; Ghosh & Karki, 2017; Ghosh et al., 2017), carbon and iron-bearing silicate melts (Karki et al., 2020; Solomatova & Caracas, 2021; Solomatova et al., 2019, 2020), and carbon partitioning between silicate and iron melts (Zhang & Yin, 2012). This work and Davis et al. (2022) represent the first ab initio studies with subequal amounts of carbon, silicon, and metal in a melt composition. In this study, we simulate seven total endmember, binary, and ternary melt compositions (Table 1) at pressures between 0 and

200 GPa and at a temperature of 4000 K to estimate melt miscibility, structure, and density. With these results, we determine the degree of mixing in the carbonate-silicate-iron ternary system, the buoyancy of the melts that form, and the chemical species that exist within the melt. We evaluate the viability of carbonate-silicate-metal melt compositions as contributors to ultra-low velocity zones (ULVZs) and the implications for carbon sequestration and distribution throughout the Earth's mantle and core.

## 2. Methods

Ab initio molecular dynamics simulations using the projector-augmented wave method (Kresse & Furthmuller, 1996) of density functional theory were performed with the Vienna ab initio simulation package (Blochl, 1994). We used the generalized gradient approximation in the Perdew-Burke-Ernzerhof form (Perdew et al., 1996) to treat electron exchange and correlation. The kinetic energy cutoffs for the plane-wave expansion of the wavefunctions were set to 600 eV. We used the canonical ensemble (NVT) with a Nosé-Hoover thermostat (Hoover, 1985; Nosé, 1984) with a time step of 1–2 fs for 18–80 ps, depending on the density. The Brillouin zone was sampled at the gamma point. The calculations were spin-polarized at all pressures. A Hubbard $U_{eff}$ ($U-J$) parameter of 4 eV was applied, which enhances the magnetic moment of the Fe atoms and corrects for their volume and coordination environment. The mean-square displacement as a function of time shows a ballistic regime below approximately 1,000 fs, after which the atoms reach a diffusive regime. For the carbonate-silicate-metal melt composition, calculations were run with a minimum of two starting configurations, and the results were averaged. We employ the Universal Molecular Dynamics package for the analysis of the results (Caracas et al., 2021).

We work with seven melts representing endmember (MgCO₃, MgSiO₃, and Fe), binary (Mg(C,Si)O₃, MgCO₃ + Fe, MgSiO₃ + Fe), and ternary (Mg(C,Si)O₃ + Fe) melt compositions (Table 1), with supercells ranging from 108 to 133 atoms. Simulations span a pressure range of 0–200 GPa, and all calculations are performed at 4000 K. Bond distances were determined from the pair distribution functions. The first peak in the pair distribution function marks the radius of the first coordination sphere for the reference atom, and the first minimum translates to the maximum acceptable bond distance for a bonding pair. The fitted minimum values were used in the speciation analysis to determine carbon clusters.

## 3. Melt Miscibilities

The Gibbs free energy of mixing, $\Delta G_{mix}$, determines whether a given solution of melt components will mix or unmix. Negative $\Delta G_{mix}$ values indicate that a mixture is energetically favorable and therefore miscible. Positive $\Delta G_{mix}$ values indicate that a mixture is energetically unfavorable and therefore immiscible. The Gibbs free energy of mixing was estimated using the following equation:

$$ \Delta G_{\text{mix}} = \Delta H_{\text{mix}} - T \Delta S_{\text{mix}} + \int P \Delta V_{\text{mix}} \quad (1) $$
where $\Delta H_{\text{mix}}$ is the enthalpy of mixing, $T$ is the temperature, $\Delta S_{\text{mix}}$ is the entropy of mixing, $P$ is the pressure, and $\Delta V_{\text{mix}}$ is the mixing volume. An example of the contribution of each term to $\Delta G_{\text{mix}}$ on the carbonate-silicate binary is shown in Figure S1 of the Supporting Information S1. We describe the calculation of each term in the equation in the following sections.

### 3.1. Enthalpy of Mixing

To calculate the enthalpy of mixing, we use the equation:

$$ \Delta H_{\text{mix}} = \sum_{i \neq j} \beta_{ij} X_i X_j \qquad (2) $$

where $\beta_{ij}$ represents the binary parameter along the $i-j$ binary, and $X_i$ and $X_j$ represent the mole fractions of the $i$ and $j$ components, respectively. Along each of the three binaries calculated in this study, $X_i$ and $X_j$ are known, but $\beta$ is unknown and is dependent on the chemistry of the mixture. $\Delta H_{\text{mix}}$ values for these compositions are hard to find, making it difficult to select the correct value of $\beta$ along each binary. Solubility numbers, however, are easier to find in the literature, and can be directly compared to numbers extracted from our calculated $\Delta G_{\text{mix}}$ curves (see Section 3.4). Thus, we employ an iterative guess and check method to determine which value of $\beta$ best matches the solubility along a given binary under specific pressure and temperature conditions that can be compared to data from the literature or from our simulations. The selected $\beta$ value is then used to calculate $\Delta H_{\text{mix}}$ (Equation 2), and this $\Delta H_{\text{mix}}$ value is used in all further calculations of $\Delta G_{\text{mix}}$ (Equation 1) under different pressure and temperature conditions. To find appropriate values for $\beta$, we plot $\Delta G_{\text{mix}}$ at 0 GPa and 4000 K along each binary and select $\beta$ values that match the expected degree of mixing (Figure S2 in Supporting Information S1). At 0 GPa, the mixing volume term of the $\Delta G_{\text{mix}}$ equation goes to zero, and the entropy term at 4000 K is the same along all three compositional binaries studied (see Section 3.2). Thus, at 0 GPa and 4000 K, $\Delta G_{\text{mix}}$ across all three binaries is influenced exclusively by the enthalpy term, and by changing the $\beta$ value of $\Delta H_{\text{mix}}$, we directly change the degree of mixing along a given binary. Reasonable degrees of mixing are determined by examining solubility experiments on binary systems in addition to our own simulation results at 0 GPa. A study of orthopyroxene solubility in carbonate melts reports that carbonate melts contain 4 atomic percent silicate at 2 GPa and 1273 K (Kamenetsky & Yaxley, 2015). We select a value of 95 kJ for $\beta$ on the carbonate-silicate join, which leads to silicate solubility of 8 atomic percent in carbonate melts at 4000 K. Silicate-metal melts are immiscible at 0 GPa (Fichtner et al., 2021), and our simulations show groupings of silicon and iron atoms that is suggestive of immiscibility. Thus, we select a $\beta$ parameter of 135 kJ, which leads to limited miscibility (2 atomic percent Fe in silicate melt) at 0 GPa and 4000 K. Experimental reports of carbonate solubility in iron melt are lacking, but our carbonate-metal simulation indicates less miscibility than the carbonate-silicate simulation and more miscibility than the silicate-metal simulation at 0 GPa. We select a value of 115 kJ for the $\beta$ parameter, which leads to 4 atomic percent Fe in the carbonate melt. $\Delta G_{\text{mix}}$ values calculated with the chosen $\beta$ parameters at 0 GPa are plotted in Figure S4 of the Supporting Information S1.

### 3.2. Entropy of Mixing

To calculate the entropy of mixing, we use the ideal entropy of mixing

$$ \Delta S_{\text{mix}} = -R(X_i \ln X_i + X_j \ln X_j + X_k \ln X_k) \qquad (3) $$

where $R$ is the gas constant and $X_i$, $X_j$, and $X_k$ are the mole fractions of the $i$, $j$, and $k$ components, respectively.

### 3.3. Mixing Volumes

Mixing volumes were calculated by taking weighted averages of molar volumes of individual melt components. To calculate molar volumes, we first fit either second- or third-order Birch-Murnaghan equations of state for our simulated melt compositions (Figure 1a). We also included the pyrolite and pyrolite + 8CO compositions from Solomatova et al. (2019) for comparison. The equations of state fit parameters are reported in Table S8 of the Supporting Information S1. The fit parameters reveal that melts with a carbonate component are highly compressible, which is in agreement with previous studies of carbon-bearing melts (Ghosh & Karki, 2017; Ghosh et al., 2007; Sakamaki et al., 2011). There is a significant covariance between $K_0'$, $K_0$, and $V_0$ values for
<figure><img src="image_3.png" /><figcaption>**Figure 1.** Density, molar volume, and mixing volume data for all simulated melt compositions are plotted as a function of pressure. (a) Melt densities are fit to second- or third-order Birch-Murnaghan equations of state. (b) Molar volumes. Due to the non-stoichiometric nature of the ternary melt composition, we report volumes for all melt compositions per mole of atoms. (c) Mixing volumes for all binary and ternary melt compositions over the pressure range of the Earth's lower mantle.</figcaption></figure>

all melts, which is common in finite strain equations of state. One continuous equation of state was fit across multiple structural transitions, which stem from gradual coordination changes in the melt. However, the majority of coordination changes occur between 0 and 20 GPa, below the pressure regime of interest. In the pressure regime of the lower mantle, the fits closely match the data. Using the Birch-Murnaghan equation of state fits, we calculate molar volumes of the melts at pressures from 0 to 200 GPa (Figure 1b). Due to the non-stoichiometric nature of the melt mixtures, we calculate volumes per mole of atoms instead of per formula unit, allowing the molar volumes to be directly compared. Iron and iron-bearing melts have the largest molar volumes, while pyrolite melts have the smallest. The densities and molar volumes of each melt composition are reported in Tables S1–S7 of the Supporting Information S1.

To calculate the mixing volumes, we compared the molar volume of our simulated mixture with the weighted average of the molar volumes of the mixture components. The magnitude of the mixing volume indicates the nonideality of a melt mixture and the degree of interaction between melt components. As the pressure derivative of $\Delta G_{\text{mix}}$, the mixing volume is the tendency of the mixture to become more or less energetically favorable with changing pressure. Thus, the sign of the mixing volume is suggestive of miscibility in multicomponent mixtures. This is especially true at high pressures, where the mixing volume term dominates the contribution to $\Delta G_{\text{mix}}$. For instance, a composition that is 50% MgSiO$_3$ and 50% Fe has a $\Delta G_{\text{mix}}$ value of 58 kJ at 136 GPa and 4000 K (see Figure 3). Of the 58 kJ, 47 kJ is from the mixing volume component, accounting for 82% of the contribution to $\Delta G_{\text{mix}}$. Mixtures with positive mixing volumes become larger upon mixing and become less stable with increasing pressure, enforcing immiscibility. The mixing volumes for the four multicomponent melts in this study were plotted as a function of pressure in Figure 1c. Over the range of the lower mantle, Mg(C,Si)O$_3$ and MgCO$_3$ + Fe have negative mixing volumes. MgSiO$_3$ + Fe and Mg(C,Si)O$_3$ + Fe have positive mixing volumes.

To expand our analysis to any composition in the MgCO$_3$–MgSiO$_3$–Fe ternary system, the mixing volume data are fit to the following power series multicomponent mixing model (Ganguly, 2001; Wohl, 1946, 1953):

$$ \Delta V_{\text{mix}} = \sum_{i \neq j} X_i X_j (W_{ij}^G X_{ji} + W_{ji}^G X_{ij}) + \sum_{i \neq j, \neq k} X_i X_j X_k C_{ijk} \quad (4) $$

where $\Delta V_{\text{mix}}$ is the mixing volume, the $W_{ij}^G$'s are the binary interaction parameters, and $C_{ijk}$ is the ternary interaction term. $X_i$, $X_j$, and $X_k$ are the mole fractions of the i, j, and k components, and $X_{ij}$ and $X_{ji}$ are the projected mole fractions of the i and j components in the binary join $i - j$. $X_{ij}$ is given analytically by $\frac{1}{2}(1 + X_i - X_j)$. Similar to our calculation of molar volumes, in this analysis, we set the mole fractions by counting the number of atoms of each component rather than the number of formula units. Thus, each iron atom in the model is compared to an average atom of either MgCO$_3$ or MgSiO$_3$. For example, in our simulated ternary melt composition (12MgSiO$_3$ + 12MgCO$_3$ + 13Fe), the mole fractions are $X_{\text{Fe}} = 0.1$, $X_{\text{MgCO}_3} = 0.45$, and $X_{\text{MgSiO}_3} = 0.45$, which are derived from the number of atoms of each component. Given the limited number of simulations, in this model, we assume that mixing along the binaries is symmetric, and thus $W_{ij}^G = W_{ji}^G$.

Figures 2a–2c display mixing volumes along the binary joins MgCO$_3$–MgSiO$_3$, MgCO$_3$–Fe, and MgSiO$_3$–Fe at 24, 50, and 136 GPa,
<figure><img src="image_4.png" /><figcaption>**Figure 2.** Mixing volumes for all binary and ternary melts are plotted as a function of melt composition. (a–c) Mixing volumes for binary solutions of (a) MgCO₃ and MgSiO₃; (b) MgSiO₃ and Fe; and (c) MgCO₃ and Fe at 24, 50, and 136 GPa. Simulated compositions are marked by data points, and lines are fits to the multicomponent mixing model. Both the MgCO₃ and MgSiO₃ binary and the MgCO₃ and Fe binary have negative mixing volumes at lower mantle pressures, suggesting miscibility. The MgSiO₃ and Fe binary has positive mixing volumes across all lower mantle pressures, suggesting immiscibility. (d–f) Mixing volumes for ternary solutions of MgCO₃, MgSiO₃, and Fe at (d) 24 GPa; (e) 50 GPa; and (f) 136 GPa. Positive mixing volumes are shaded red and negative mixing volumes are shaded blue. The simulation compositions used to fit the mixing model are marked by the data points, and the ternary composition always lies within the immiscible region of the plot.</figcaption></figure>

all at the temperature of 4000 K, calculated using the fit parameters reported in Table S9 of the Supporting Information S1. Along the carbonate-silicate join, mixing volumes are negative at all lower mantle conditions and become more negative with increasing pressure, suggesting continuous solubility. Along the silicate-metal join, mixing volumes are positive across all lower mantle pressure conditions, and decrease with pressure. At pressures beyond those of the Earth's mantle, we would expect silicate and metal melt to become miscible. Along the carbonate-metal join, mixing volumes are always negative but become less negative with increasing pressure, suggesting potential immiscibility beyond the core-mantle boundary pressure. Additionally, the magnitude of the mixing volumes represents the degree of interaction between the melts. Generally, silicate-metal melts have the most interaction, followed by carbonate-metal melts, and carbonate-silicate melts. Carbonate-silicate melt interaction terms are very small, even at their most negative point at 136 GPa, indicating that this mixture is close to ideal.

Figures 2d–2f show the calculated mixing volumes for ternary compositions. Note that the mixing model was fit using four data points (three binary compositions and one ternary composition), which are labeled in Figures 2d–2f. At 24 GPa, melts with greater than 50% carbonate have negative mixing volumes. For melts less than 50% carbonate, mixing volumes are more negative with increasing iron percentage and more positive with increasing silicate percentage. With increasing pressure, the negative mixing volume regime shrinks and the positive mixing volume regime grows to cover more iron- and carbonate-rich parts of the ternary plot. By 136 GPa, only compositions that are greater than 70% carbonate and compositions close to the carbonate-metal and carbonate-silicate binaries have negative mixing volumes. Additionally, mixing volume magnitudes decrease with increasing pressure, indicating that these melts tend to become more ideal with increasing pressure. This conclusion is supported by the trends in the binary and ternary interaction parameters (Figure S3 and Table S9 in
<figure><img src="image_5.png" /><figcaption>**Figure 3.** $\Delta G_{\text{mix}}$ for all binary and ternary melts are plotted as a function of melt composition. (a–c) $\Delta G_{\text{mix}}$ for binary solutions of (a) MgCO$_3$ and MgSiO$_3$; (b) MgSiO$_3$ and Fe; and (c) MgCO$_3$ and Fe at 24, 50, and 136 GPa. MgCO$_3$ and MgSiO$_3$ have limited miscibility near the endmembers that increases with pressure, and MgSiO$_3$ and Fe are immiscible at all pressures studied. MgCO$_3$ and Fe have limited miscibility at 24 GPa, and the miscibility gap closes with increasing pressure. (d–f) $\Delta G_{\text{mix}}$ for ternary solutions of MgCO$_3$, MgSiO$_3$, and Fe at (d) 24 GPa; (e) 50 GPa; and (f) 136 GPa. Positive values of $\Delta G_{\text{mix}}$ are shaded red and negative values are shaded blue. Solid gold lines outline miscible melt compositions, two coexisting melt regions, and three coexisting melt regions, labeled with 1, 2, or 3, respectively, and are estimated from the $\Delta G_{\text{mix}}$ values on the binaries. Dashed gold lines are tie lines in the two melt regions. Data points represent compositions simulated in this study.</figcaption></figure>

Supporting Information S1). With increasing pressure, the interaction parameters trend toward 0, indicating that the interaction between melt components becomes increasingly less important with depth.

### 3.4. Gibbs Free Energy of Mixing

With the equations and approximations describing $\Delta H_{\text{mix}}$, $\Delta S_{\text{mix}}$, and $\Delta V_{\text{mix}}$, we determine how $\Delta G_{\text{mix}}$ evolves along binary and ternary joins. $\Delta G_{\text{mix}}$ is plotted along the carbonate-silicate, carbonate-metal, and silicate-metal binaries in Figures 3a–3c. MgCO$_3$ and MgSiO$_3$ melts demonstrate limited miscibility at all lower mantle pressures, with a large miscibility gap. Note that although $\Delta G_{\text{mix}}$ is negative for many mixture compositions plotted along the MgCO$_3$ – MgSiO$_3$ binary, not all mixture compositions are miscible (see Figure 4 for the translation of the Gibbs free energy curve to a binary phase diagram). In the immiscible region, two melt compositions coexist, and these compositions are determined by the common tangent of the $\Delta G_{\text{mix}}$ curves. These tangents are quasi-horizontal and have support points that are very close to the minima of the free energy. The shared tangents of the curves are plotted in Figure 4 and are at ~9 and ~91 mol percent MgCO$_3$ at all pressures examined for the carbonate-silicate binary. $\Delta G_{\text{mix}}$ decreases with pressure, suggesting eventual closing of the miscibility gap at higher pressures than those reached by the Earth's mantle. MgSiO$_3$ and Fe are immiscible at all lower mantle pressures. The metallic character of the pure Fe melt makes it incompatible with the insulating character of the molten silicate melts, and any iron that is dissolved in the silicate is always incorporated as an ionic phase, FeO or Fe$_2$O$_3$. Similar to the MgSiO$_3$ and MgCO$_3$ binary, the MgCO$_3$ and Fe binary also has a miscibility gap that begins to close with increasing pressure. At 24 GPa, the two coexisting melt compositions are at 6 and 94 mol percent MgCO$_3$, but by 136 GPa, the two coexisting melt compositions are at 38 and 62 MgCO$_3$ mole percentage (Figure 4).
<figure><img src="image_6.png" /><figcaption>**Figure 4.** Phase diagram of binary melt mixtures. The central region inside each of the pairs of lines is the miscibility gap, where it is more energetically favorable to form two melt immiscible melt compositions. The miscibility gap closes with increasing pressure for carbonate-metal mixtures and widens for silicate-metal mixtures. Although the carbonate-silicate miscibility gap remains mostly constant at lower-mantle conditions, a decreasing $\Delta G_{\text{mix}}$ suggests the eventual closing of the miscibility gap at higher pressures.</figcaption></figure>

The ternary diagrams (Figures 3d–3f) show a range of miscibilities that expand with increasing pressure. As in the case of the binary phase diagrams, the miscibility fields for the ternary mixture compositions are set by the shared tangent in the $\Delta G_{\text{mix}}$ curves, and not solely by where $\Delta G_{\text{mix}}$ is calculated to be negative. These miscible melt fields, outlined by solid gold lines, are extended from $\Delta G_{\text{mix}}$ values along each of the binaries. There are three miscibility fields, and each is located near an endmember composition. As the carbonate-metal and the carbonate-silicate miscibility gaps close with pressure, the miscibility fields grow to accommodate more mixing. In between the miscibility fields are two-melt regions, and compositions that fall in these regions will exsolve two immiscible melts. Dashed gold lines represent example tie lines in these regions. The central triangle is the three-melt region, and compositions that fall in this region will exsolve three immiscible melt compositions. With increasing pressure, the two-melt regions grow and the three-melt region shrinks, indicating the overall increase of miscibility in this system at high pressure. It is important to note that the chosen value for $\Delta H_{\text{mix}}$ affects the miscibilities of the melt mixtures. We estimate values for $\Delta H_{\text{mix}}$ based on chemical speciation and literature experiments, but without additional constraints, there is some ambiguity in the selected value. The blue regions of the plot indicate compositions with negative $\Delta G_{\text{mix}}$ values, and thus, show a possible range of miscible compositions that are available under smaller $\Delta H_{\text{mix}}$ values. Nonetheless, the results based on these reasonable estimates of $\Delta H_{\text{mix}}$ illustrate the plausibility of reduced immiscibility with increasing pressure in the carbonate-silicate-metal system, such that an Fe-rich carbonate melt and a carbon-rich Fe melt would be expected to segregate from other phases at the base of the mantle.

## 4. Carbon-Bearing Clusters

As evidenced by Figures 3d–3f, many carbonate-silicate-iron melt compositions are immiscible at lower mantle conditions, even at the high-pressure conditions of the core-mantle boundary. Although ab initio molecular dynamics cannot model phase separation because of size effects, the clustering of species observed in our simulations is suggestive of the process of melt separation. In this section, we identify the species that segregate in an example ternary melt mixture through speciation analysis and determine individual cluster densities to understand how elements distribute through the lower mantle.

The example melt composition is the ternary melt composition simulated for this study (10% Fe, 45% $\text{MgCO}_3$, 45% $\text{MgSiO}_3$). A more complete speciation analysis of this melt composition is reported by Davis et al. (2022). In
<figure><img src="image_7.png" /><figcaption>**Figure 5.** Examples of types of carbon-bearing clusters and their abundances. (a–c) Examples of isolated carbon-iron-oxygen clusters classified as carbonate, carbide, and polymer from 74 GPa to 4000 K. (d) Abundances of carbonate, carbide, and polymer clusters with increasing pressure.</figcaption></figure>

general, large increases in C-Fe and C-C bonding with pressure at the expense of C-O bonding are observed. No evidence of Fe-Si bonding is found at any pressure. These tendencies are in good agreement with the predicted miscibilities for binary solutions (Figures 3a–3c). In the melt, large carbon-carbon clusters and, at higher pressures, carbon-iron clusters form and have limited interaction with the silicate melt network, indicating the types of carbon-bearing melt species we might expect to segregate from a silicate melt. However, the extent of carbon-iron interaction is difficult to quantify. Carbon-iron clusters often consist of iron atoms surrounding a polymerized carbon core, which could be classified as either a diamond seed nucleus or an iron carbide cluster. Previous simulations (Davis et al., 2022; Karki et al., 2020; Solomatova et al., 2019) show that we would expect carbon to bond to O, Fe, C, and Si. For each carbon atom, we would expect 12% (13/108) of the bonds to be C-Fe bonds as there are 13 iron atoms from available 108 coordinating anions for carbon (12 silicon, 72 oxygen, 11 out of the 12 carbon, and 13 iron). Additionally, we anticipate that 67% (72/108) of carbon bonds are to oxygen. Starting from these estimates based on the statistical sampling, we classify the carbon-based clusters in the melt. Bond abundances greater than the abundances predicted from statistical sampling indicate that there is a chemical preference for the bonding element. In our cluster analysis, we classify carbon clusters with C-O bond abundances greater than 67% as carbonates and carbon clusters with C-Fe bond abundances greater than 12% as carbides. Clusters with less C-O and C-Fe abundances than expected from statistical sampling are classified as carbon polymers. The abundances of the different types of clusters at 1, 74, and 148 GPa are plotted in Figure 5d. At all three pressures, clusters of each type are formed, but the relative abundances of the cluster types evolve with pressure. Carbonates are the most abundant cluster type at 1 GPa and account for 49% of the total clusters, but that number drops to 31% at 74 GPa and 30% at 148 GPa. Carbides almost match the number of carbonate clusters
<figure><img src="image_8.png" /><figcaption>**Figure 6.** Carbon clusters identified from three separate snapshots at (a) 74 GPa and (b) 148 GPa and plotted by the fraction of bonds that are C-O and the fraction of bonds that are C-Fe. Points circled with a dotted line correspond to carbonates, with a solid line correspond to carbide, and with a dashed line correspond to polymers. Data points are colored by the density of the cluster. Clusters denser than MgSiO₃ (blue) will sink and clusters less dense than MgSiO₃ (red) will float.</figcaption></figure>

at 1 GPa, at 48% of the total, and are the most abundant cluster type at 74 and 148 GPa, at 55% and 64% of the total, respectively. Polymers are always the least abundant cluster type. They increase in abundance from 3% to 14% from 1 to 74 GPa, and then decrease in abundance to 5% at 148 GPa. In Davis et al. (2022), we noticed that the majority of changes in C-O and C-C bond abundances occur in the first 25 GPa. Thus, we expect diamond formation to peak around 25 GPa, and this expectation is reflected in the relative increase in polymer formation between 1 and 74 GPa. Similarly, we expect carbonate cluster abundance to decrease rapidly in the first 25 GPa before plateauing, and this result is also observed. Finally, the large and linearly increasing number of carbide clusters matches the speciation results in both Davis et al. (2022) and Solomatova et al. (2019), which report linear increases in C-Fe bond abundances with increasing pressure.

The composition and the volume of the carbon clusters determine their relative density within the mantle. Using the Bader charge analysis algorithm (Henkelman et al., 2006; Sanville et al., 2007; Tang et al., 2009; Yu & Trinkle, 2011), we calculate the volumes of individual atoms within carbon clusters to determine cluster densities. The densities of example carbon clusters isolated at 74 and 148 GPa are plotted in Figure 6. The selected clusters are grouped according to their classification as a carbonate, carbide, or polymer. We directly compare the density of the cluster to the calculated density of MgSiO₃ melt under the same conditions. Carbide clusters are much denser than MgSiO₃ melt, and with enough time and aggregation, we expect these clusters to segregate from the multicomponent melt and sink to the core. Similarly, polymers are slightly denser than MgSiO₃ melt. Carbonate clusters are lighter than MgSiO₃ melt, and we expect these clusters to be buoyant within the mantle.

**5. Implications**

From the miscibility analysis, we found three miscible melt compositional fields: carbonate-rich, silicate-rich, and iron-rich melts (Figures 3d–3f). We consider the densities of these melt compositions to determine their buoyancy in the lower mantle and to evaluate their implications for carbon distribution and sequestration in the lower mantle and core. Densities of liquids calculated from molecular dynamics methods have been shown to systematically deviate from experimental values depending on the approximation used for the exchange correlation functional (Zhang et al., 2013; Zhao et al., 2014). However, relative comparisons of density between calculated melts are useful, provided the same approximations are made. As an example, density differences from MgSiO₃ melt at the core-mantle boundary are plotted in Figure 7. Here, miscible melt compositions could be formed from a deep Earth carbonatite melt interacting with iron melt at the core-mantle boundary. Compositions are shaded in red, white, and blue to represent buoyant, neutrally buoyant, and dense compositions, respectively, as compared to the density of MgSiO₃ melt, which is used as a proxy for the lower mantle composition. Three groups of miscible compositions emerge. Compositions outlined in red are iron-rich and denser than MgSiO₃. We expect these compositions to sink into the outer core, dragging carbon and silicon out of the mantle and enriching the outer core with light elements over time. Compositions outlined in blue are carbonate-rich and buoyant. These compositions are anticipated to rise through the mantle and return
<figure><img src="image_9.png" /><figcaption>**Figure 7.** Densities of ternary melt compositions at 136 GPa compared to the density of MgSiO₃ melt. Melts denser than MgSiO₃ and less dense than pure Fe are colored blue, and would sink in the mantle. Melts less dense than MgSiO₃ are colored red, and would float. Miscible melt compositions that would sink, float, and be neutrally buoyant are outlined with blue, red, and green lines, respectively.</figcaption></figure>

carbon to shallower depths. Finally, the green compositions are neutrally buoyant and thus gravitationally stable at the core-mantle boundary. These compositions could serve as possible contributors to ULVZs. ULVZs have multiple proposed explanations, including FeSi formed through core-mantle reactions (Mergner et al., 2021), hydrous phases such as (Al,Fe)OOH (Thompson et al., 2021), Fe-rich post perovskite (Garnero & McNamara, 2008), Fe-rich (Mg,Fe) O (Solomatova et al., 2016; Wicks et al., 2010), and patches of partial melt (Williams & Garnero, 1996). Partial melt is a likely explanation for ULVZs due to the 3:1 ratio of S-to-P wave velocity reduction (Garnero & McNamara, 2008; Williams & Garnero, 1996). Therefore, buoyantly neutral melt compositions, such as the carbonate-silicate-metal melt compositions calculated in this study, could serve as one possible explanation and contribute to the ULVZs.

Within the immiscible melt compositions, carbon, carbon-iron, and carbon-oxygen clusters form (Figure 5). Given enough time and aggregation, we expect the carbon and carbon-iron clusters to exsolve from the melt, as has been previously suggested (Dasgupta & Hirschmann, 2010; Karki et al., 2020; Mysen et al., 2011; Stagno et al., 2013). In our example ternary melt composition, the majority of the clusters formed at 148 GPa are carbide (64%), and the propensity for carbon to bond with iron indicates the high siderophily of carbon under these thermodynamic conditions. Carbide clusters are denser than the surrounding mantle (Figure 6). Thus, a significant amount of Earth's carbon contained in the lower mantle may bond with iron and sink to the core, matching previous ab initio predictions of carbon's fate under reduced conditions in the lower mantle (Karki et al., 2020; Rohrbach & Schmidt, 2011). This not only prevents carbon from being recycled back to the Earth's surface but also changes the evolution of the core composition. An increasingly carbon-rich core composition would evolve to have density, sound velocities, and electrical and thermal conductivity more similar to the carbon-rich alloys Fe₃C and Fe₇C₃ (Fiquet et al., 2009; Ghosh & Karki, 2017; Wood et al., 2013). Moreover, given the chemical preference of carbon to be bonded to iron, we propose that during core formation, iron droplets that segregate from the magma ocean and fall downwards would constitute strong attraction basins for carbon. In this way, the magma ocean would be leached of its carbon. After the Moon-forming impact, metal and silicate melt would be well-mixed, and siderophile elements, such as carbon, would be segregated with iron into the core, supporting the idea that carbon is a candidate element to explain the density deficit in the core (Prescher et al., 2015; Solomatova et al., 2019).

In addition to carbide clusters, we observed the formation of carbon polymers in our simulated ternary melt composition, which could be precursors for diamonds. Our simulations reveal a possible mechanism for diamond formation, where carbon polymers exsolve from a silicate melt. Previously, this formation mechanism was observed in
oxygen-deficient carbon-bearing silicate melts (Ghosh et al., 2017), and the addition of iron in our simulations may actually increase carbon polymerization (Belonoshko et al., 2015). Our previous speciation analysis of this melt composition (Davis et al., 2022) indicates that carbon-carbon bond formation reaches a peak around 25 GPa, beyond which it plateaus. C-Fe bonding, however, increases linearly and with increasing depth, the percentage of polymers decreases as carbide clusters are preferentially formed. Thus, our analysis indicates a diamond formation zone around 25 GPa. This depth in the Earth matches reports of diamonds with a deep Earth origin, which are returned from either the transition zone or the top of the lower mantle (Smith et al., 2016; Stachel et al., 2005). Ultradeep diamonds are thought to crystallize from a carbon-rich melt composition (Haggerty, 1999), which matches the formation mechanism outlined by our speciation analysis. Additionally, this formation mechanism is thought to play a role in eclogitic (E-type) diamonds (Gao et al., 2017; Haggerty, 1999). E-type diamonds crystallize within an eclogitic source rock, and have been reported to form at depths ranging from 150 to >660 km (Kirkley et al., 1991; Meyer, 1985; Moore & Gurney, 1985; Shirey & Shigley, 2013). They tend to form at higher pressure and temperature conditions than peridotitic (P-type) diamonds (Deines et al., 1993, 2000; Shirey & Shigley, 2013). E-type diamonds tend to have a wide range of δ¹³C values and are a good match for the carbon isotopic values of diamonds returned from the transition zone (Shirey & Shigley, 2013). Thus, the crystallization of a melt composition similar to our studied carbonate-silicate-metal melt, which forms abundant carbon polymers around the depth of the transition zone, may be a formation mechanism for ultradeep E-type diamonds. Our cluster analysis (Figure 5) also indicates that polymers are formed even at the core-mantle boundary, and if these polymers aggregate to form diamonds, these diamonds may be brought to the surface by deep mantle plumes. Diamonds with a lowermost mantle origin may be identified through compositional analysis of fluid inclusions. Diamonds containing silicate-poor metal-rich carbonate melt compositions that fall into the miscible melt regions indicated in Figure 3f would indicate a core-mantle boundary origin and would provide evidence for carbonate-silicate-metal melt reactions in the lowermost mantle.

Finally, we examine the carbon distribution at pressure and temperature conditions of the core-mantle boundary to provide some insight into possible carbon distributions between core and mantle phases. We examine an equimolar composition (i.e., 1/3 Fe, 1/3 MgCO₃, and 1/3 MgSiO₃), which falls into the three-melt region at the center of the ternary plot (Figure 3f). Of the three melts that exsolve from this composition, 31% is a carbonate-rich melt, 37% is an iron-rich melt, and 32% is a silicate-rich melt, where the melt compositions that exsolve are given by the corners of the miscible melt fields in Figure 3f. From mass balance calculations, we determine that for this case, 90% of the carbon is contained in outer-core compositions (carbonate-rich and iron-rich melts) and 10% is contained in a potential ULVZ composition (silicate-rich melt). In fact, carbon is distributed to varying degrees between outer-core and ULVZ-type compositions for the majority of compositions in this ternary system, and buoyant carbonate-rich melt phases only form when there is less than ~5% Fe in the system. At the core-mantle boundary where iron melt is abundant, we anticipate that carbon is preferentially stored in the lower-mantle and core phases, indicating that the ultimate fate of Earth's carbon may be storage in the deep interior.

## 6. Conclusions

Carbonates are important compounds in the crust and upper mantle and may play a role in the lower mantle as well. Carbonate melts in the deep Earth may react with silicates and metals, especially at the core-mantle boundary where these phases are abundant. The chemical and physical properties of the melts that form from these reactions have important consequences for the distribution and storage of carbon in the deep Earth. Ab initio molecular dynamics simulations of carbonate-silicate-iron melt compositions allow for the examination of melt miscibilities, densities, and speciation. We find that carbonate-silicate and carbonate-iron melts have miscibility gaps that close with increasing pressure, and that carbonate-iron melts have the highest affinity for mixing. Silicate-iron melts are immiscible at all lower mantle pressures. By expanding our analysis to the ternary carbonate-silicate-iron system, we find that three miscible melt fields exist near each of the endmember compositions (Fe-rich, carbonate-rich, and silicate-rich melts). Iron-rich melts are dense and sink into the core, providing a mechanism to enrich the outer core with light elements, such as carbon, oxygen, and silicon. Silicate-rich melts are neutrally buoyant and sit at the core-mantle boundary, providing one possible explanation for the existence of ULVZs. Carbonate-rich melts, depending on their iron content, may sink into the core, remain at the core-mantle boundary, or rise through the mantle. Thus, depending on the composition that forms through the reaction of carbonate, silicate, and iron phases, carbon may be stored in the deep Earth in core- or
ULVZ-type compositions or may return to shallower depths. The majority of melt compositions have densities that classify them as core- or ULVZ-type compositions, indicating that the fate of carbon may be to be stored in the Earth's deep interior. Finally, the speciation of carbonate-silicate-iron melts indicates that carbon polymers, iron carbides, and carbonate clusters are formed in the melt, and the relative proportions of these clusters at various pressures indicate carbon's changing affinity for the other elements. Iron carbides, which are favorably formed at higher pressures, indicate carbon's increasingly siderophile nature with depth. Carbon polymers, when aggregated, could form diamonds and are abundant at transition zone pressures, indicating a propensity for diamond formation in and around the transition zone. The distribution of carbon throughout the Earth's interior is a complicated topic, affected by many thermodynamic variables, including pressure, temperature, composition, and oxygen fugacity. More experimental and computational studies of carbonate melts and their interactions with other phases at lower-mantle and especially core-mantle boundary conditions will help elucidate the role of carbon in the Earth's deep interior.

## Conflict of Interest

The authors declare no conflicts of interest relevant to this study.

## Data Availability Statement

All simulation output presented in this manuscript is publicly available at (Davis et al., 2023).

## References

* Bajgain, S. K., & Mookherjee, M. (2021). Carbon bearing aluminosilicate melt at high pressure. *Geochimica et Cosmochimica Acta*, 312, 106–123. https://doi.org/10.1016/j.gca.2021.07.039

* Belonoshko, A. B., Lukinov, T., Rosengren, A., Bryk, T., & Litasov, K. D. (2015). Synthesis of heavy hydrocarbons at the core-mantle boundary. *Scientific Reports*, 5(1), 18382. https://doi.org/10.1038/srep18382

* Blochl, P. E. (1994). Projector augmented-wave method. *Physical Review B*, 50(24), 17953–17979. https://doi.org/10.1103/PhysRevB.50.17953

* Caracas, R., Kobsch, A., Solomatova, N. V., Li, Z., Soubiran, F., & Hernandez, J.-A. (2021). Analyzing melts and fluids from ab initio molecular dynamics simulations with the UMD package. *Journal of Visualized Experiments*, (175). https://doi.org/10.3791/61534-v

* Dasgupta, R., & Hirschmann, M. M. (2010). The deep carbon cycle and melting in Earth's interior. *Earth and Planetary Science Letters*, 298(1–2), 1–13. https://doi.org/10.1016/j.epsl.2010.06.039

* Dasgupta, R., & Walker, D. (2008). Carbon solubility in core melts in a shallow magma ocean environment and distribution of carbon between the Earth's core and the mantle. *Geochimica et Cosmochimica Acta*, 72(18), 4627–4641. https://doi.org/10.1016/j.gca.2008.06.023

* Davis, A. H., Solomatova, N. V., Campbell, A. J., & Caracas, R. (2022). The speciation and coordination of a deep Earth carbonate-silicate-metal melt. *Journal of Geophysical Research: Solid Earth*, 127(3), e2021JB023314. https://doi.org/10.1029/2021JB023314

* Davis, A. H., Solomatova, N. V., Caracas, R., & Campbell, A. J. (2023). Output files for carbon storage in the Earth's deep interior implied by carbonate-silicate-metal melt miscibility [Dataset]. Zenodo. https://doi.org/10.5281/zenodo.7093216

* Deines, P., Harris, J. W., & Gurney, J. J. (1993). Depth-related carbon isotope and nitrogen concentration variability in the mantle below the Orapa kimberlite, Botswana, Africa. *Geochimica et Cosmochimica Acta*, 57(12), 2781–2796. https://doi.org/10.1016/0016-7037(93)90390-I

* Deines, P., Viljoen, F., & Harris, J. W. (2000). Implications of the carbon isotope and mineral inclusions record for the formation of diamonds in the mantle underlying a mobile belt: Ventia, South Africa. *Geochimica et Cosmochimica Acta*, 65(5), 813–838. https://doi.org/10.1016/S0016-7037(00)00569-X

* Dorfman, S. M., Badro, J., Nabiei, F., Prakapenka, V. B., Cantoni, M., & Gillet, P. (2018). Carbonate stability in the reduced lower mantle. *Earth and Planetary Science Letters*, 489, 84–91. https://doi.org/10.1016/j.epsl.2018.02.035

* Drewitt, J. W. E., Walter, M. J., Zhang, H., McMahon, S. C., Edwards, D., Heinen, B. J., et al. (2019). The fate of carbonate in oceanic crust subducted into Earth's lower mantle. *Earth and Planetary Science Letters*, 511, 213–222. https://doi.org/10.1016/j.epsl.2019.01.041

* Fichtner, C. E., Schmidt, M. W., Liebske, C., Bouvier, A.-S., & Baumgartner, L. P. (2021). Carbon partitioning between metal and silicate melts during Earth accretion. *Earth and Planetary Science Letters*, 554, 116659. https://doi.org/10.1016/j.epsl.2020.116659

* Fiquet, G., Badro, J., Gregoryanz, E., Fei, Y., & Occelli, F. (2009). Sound velocity in iron carbide (Fe₂C) at high pressure: Implications for the carbon content of the Earth's inner core. *Physics of the Earth and Planetary Interiors*, 172(1–2), 125–129. https://doi.org/10.1016/j.pepi.2008.05.016

* Ganguly, J. (2001). In C. A. Geiger (Ed.), *Solid solutions in silicate and oxide systems* (pp. 37–70). Etvs University Press. (European Mineralogical Union Notes in Mineralogy).

* Gao, J., Niu, J., Qin, S., & Wu, X. (2017). Ultradeep diamond originate from deep subducted sedimentary carbonates. *Science China Earth Sciences*, 60(2), 207–217. https://doi.org/10.1007/s11430-016-5151-4

* Garnero, E. J., & McNamara, A. K. (2008). Structure and dynamics of Earth's lower mantle. *Science*, 320(5876), 626–628. https://doi.org/10.1126/science.1148028

* Ghosh, D. B., Bajgain, S. K., Mookherjee, M., & Karki, B. B. (2017). Carbon-bearing silicate melt at deep mantle conditions. *Scientific Reports*, 7(1), 848. https://doi.org/10.1038/s41598-017-00918-x

* Ghosh, D. B., & Karki, B. B. (2017). Transport properties of carbonated silicate melt at high pressure. *Science Advances*, 3(12). https://doi.org/10.1126/sciadv.1701840

* Ghosh, S., Ohtani, E., Litasov, K., Suzuki, A., & Sakamaki, T. (2007). Stability of carbonated magmas at the base of the Earth's upper mantle. *Geophysical Research Letters*, 34(22), L22312. https://doi.org/10.1029/2007GL031349

* Haggerty, S. E. (1999). A diamond trilogy: Superplumes, supercontinents, and supernovae. *Science*, 285(5429), 851–860. https://doi.org/10.1126/science.285.5429.851

## Acknowledgments

We thank the anonymous reviewers for their constructive criticism and discussion of the paper. This research was supported with a Graduate Research Fellowship (DGE 1746045) to AHD and with a National Science Foundation GROW Grant to AHD, by the European Research Council (ERC) under the European Union's Horizon 2020 research and innovation program (Grant Agreement No. 681818-IMPACT to RC), and by the National Science Foundation with Grant EAR—1651017 to AJC. We also acknowledge access to the PSMN supercomputer at ENS Lyon to perform all calculations.
Henkelman, G., Arnaldsson, A., & Jónsson, H. (2006). A fast and robust algorithm for Bader decomposition of charge density. *Computational Materials Science*, **36**(3), 354–360. https://doi.org/10.1016/j.commatsi.2005.04.010

Hoover, W. G. (1985). Canonical dynamics: Equilibrium phase-space distributions. *Physical Review A: General Physics*, **31**(3), 1695–1697.
https://doi.org/10.1103/PhysRevA.31.1695

Javoy, M. (1997). The major volatile elements of the Earth: Their origin, behavior, and fate. *Geophysical Research Letters*, **24**(2), 177–180.
https://doi.org/10.1029/96GL03931

Kamenetsky, V. S., & Yaxley, G. M. (2015). Carbonate-silicate liquid immiscibility in the mantle propels kimberlite magma ascent. *Geochimica et Cosmochimica Acta*, **158**, 48–56. https://doi.org/10.1016/j.gca.2015.03.004

Karki, B. B., Ghosh, D. B., & Banjara, D. (2020). Mixed incorporation of carbon and hydrogen in silicate melts under varying pressure and redox conditions. *Earth and Planetary Science Letters*, **549**, 116520. https://doi.org/10.1016/j.epsl.2020.116520

Kirkley, M. B., Gurney, J. J., & Levinson, A. A. (1991). Age, origin, and emplace of diamonds: Scientific advances in the last decade. *Gems and Gemology*, **27**(1), 2–25. https://doi.org/10.5741/gems.27.1.2

Korsakov, A. V., & Hermann, J. (2006). Silicate and carbonate melt inclusions associated with diamonds in deeply subducted carbonate rocks. *Earth and Planetary Science Letters*, **241**(1–2), 104–118. https://doi.org/10.1016/j.epsl.2005/10.037

Koura, N., Kohara, S., Takeuchi, K., Takahashi, S., Curtiss, L. A., Grimsditch, M., & Saboungi, M. (1996). Alkali carbonates: Raman spectroscopy, ab initio calculations, and structure. *Journal of Molecular Structure*, **382**(3), 163–169. https://doi.org/10.1016/0022-2860(96)09314-3

Kresse, G., & Furthmuller, J. (1996). Efficient iterative schemes for ab initio total-energy calculations using a plane-wave basis set. *Physical Review B: Condensed Matter*, **54**(16), 11169–11186. https://doi.org/10.1103/physrevb.54.11169

Labrosse, S., Hernlund, J. W., & Coltice, N. (2007). A crystallizing dense magma ocean at the base of the Earth's mantle. *Nature*, **450**(7171), 866–869. https://doi.org/10.1038/nature06355

Li, Z. Y., Li, J., Lange, R., Liu, J. C., & Mintzer, B. (2017). Determination of calcium carbonate and sodium carbonate melting curves up to Earth's transition zone pressures with implications for the deep carbon cycle. *Earth and Planetary Science Letters*, **457**, 395–402. https://doi.org/10.1016/j.epsl.2016.10.027

Lv, M., Dorfman, S. M., Badro, J., Borensztajn, S., Greenberg, E., & Prakapenka, V. B. (2021). Reversal of carbonate-silicate cation exchange in cold slabs in Earth's lower mantle. *Nature Communications*, **12**(1), 1712. https://doi.org/10.1038/s41467-021-21761-9

McDonough, W. F., & Sun, S. S. (1995). The composition of the Earth. *Chemical Geology*, **120**(3–4), 223–253. https://doi.org/10.1016/0009-2541(94)00140-4

Mergner, V., Kupenko, I., Spiekermann, G., Petitgirard, S., Libon, L., Chariton, S., et al. (2021). Sound velocities in FeSi at lower mantle conditions and the origin of ultra-low velocity zones. *Geophysical Research Letters*, **48**(14), e2020GL092257. https://doi.org/10.1029/2020GL092257

Merlini, M., Crichton, W. A., Hanfland, M., Gemmi, M., Müller, H., Kupenko, I., & Dubrovinsky, L. (2012). Structures of dolomite at ultrahigh pressure and their influence on the deep carbon cycle. *Proceedings of the National Academy of Sciences of the United States of America*, **109**(34), 13509–13514. https://doi.org/10.1073/pnas.1201336109

Meyer, H. O. A. (1985). Genesis of diamond: A mantle saga. *American Mineralogist*, **70**, 344–355.

Moore, R. O., & Gurney, J. J. (1985). Pyroxene solid solution in garnets included in diamond. *Nature*, **318**(6046), 553–555. https://doi.org/10.1038/318553a0

Mysen, B. O., Kumamoto, K., Cody, G. D., & Fogel, M. L. (2011). Solubility and solution mechanisms of C-O-H volatiles in silicate melt with variable redox conditions and melt composition at upper mantle temperatures and pressures. *Geochimica et Cosmochimica Acta*, **75**(20), 6183–6199. https://doi.org/10.1016/j.gca.2011.07.035

Nosé, S. (1984). A unified formulation of the constant temperature molecular dynamics methods. *The Journal of Chemical Physics*, **81**(1), 511–519. https://doi.org/10.1063/1.447334

Oganov, A. R., Ono, S., Ma, Y., Glass, C. W., & Garcia, A. (2008). Novel high-pressure structures of $MgCO_3$, $CaCO_3$, and $CO_2$ and their role in Earth's lower mantle. *Earth and Planetary Science Letters*, **273**(1–2), 38–47. https://doi.org/10.1016/j.epsl.2008.06.005

Perdew, J. P., Burke, K., & Ernzerhof, M. (1996). Generalized gradient approximation made simple. *Physical Review Letters*, **77**(18), 3865–3868. https://doi.org/10.1103/PhysRevLett.77.3865

Plank, T., & Manning, C. E. (2019). Subducting carbon. *Nature*, **574**(7778), 343–352. https://doi.org/10.1038/s41586-019-1643-z

Prescher, C., Dubrovinsky, L., Bykova, E., Kupenko, I., Glazyrin, K., Kantor, A., et al. (2015). High Poisson's ratio of Earth's inner core explained by carbon alloying. *Nature Geoscience*, **8**(3), 220–223. https://doi.org/10.1038/ngeo2370

Rohrbach, A., & Schmidt, M. W. (2011). Redox freezing and melting in the Earth's deep mantle resulting from carbon-iron redox coupling. *Nature*, **472**(7342), 209–212. https://doi.org/10.1038/nature09899

Sakamaki, T., Ohtani, E., Urakawa, S., Terasaki, H., & Katayama, Y. (2011). Density of carbonated peridotite magma at high pressure using an X-ray absorption method. *American Mineralogist*, **96**(4), 553–557. https://doi.org/10.2138/am.2011.3577

Sanville, E., Kenny, S. D., Smith, R., & Henkelman, G. (2007). Improved grid-based algorithm for Bader charge allocation. *Journal of Computational Chemistry*, **28**(5), 899–908. https://doi.org/10.1002/jcc.20575

Shirey, S. B., & Shigley, J. E. (2013). Recent advances in understanding the geology of diamonds. *Gems and Gemology*, **49**(4), 188–222. https://doi.org/10.5741/GEMS.49.4.188

Smith, E. M., Shirey, S. B., Nestola, F., Bullock, E. S., Wang, J., Richardson, S. H., & Wang, W. (2016). Large gem diamonds from metallic liquid in Earth's deep mantle. *Science*, **354**(6318), 1403–1405. https://doi.org/10.1126/science.aal1303

Solomatova, N., Caracas, R., & Cohen, R. (2020). Carbon speciation and solubility in silicate melts. In C. E. Manning, J.-F. Lin, & W. L. Mao (Eds.), *Carbon in Earth's interior*. https://doi.org/10.1002/9781119508229.ch16

Solomatova, N. V., & Caracas, R. (2021). Buoyancy and structure of volatile-rich silicate melts. *Journal of Geophysical Research: Solid Earth*, **126**(2), e2020JB021045. https://doi.org/10.1029/2020JB021045

Solomatova, N. V., Caracas, R., & Manning, C. E. (2019). Carbon sequestration during core formation implied by complex carbon polymerization. *Nature Communications*, **10**(1), 789. https://doi.org/10.1038/s41467-019-08742-9

Solomatova, N. V., Jackson, J. M., Sturhahn, W., Wicks, J. K., Zhao, J., Toellner, T. S., et al. (2016). Equation of state and spin crossover of (Mg,Fe)O at high pressure, with implications for explaining topographic relief at the core-mantle boundary. *American Mineralogist*, **101**(5), 1084–1093. https://doi.org/10.2138/am-2016-5510

Stachel, T., Brey, G. P., & Harris, J. W. (2005). Inclusions in sublithospheric diamonds: Glimpses of deep Earth. *Elements*, **1**(2), 73–78. https://doi.org/10.2113/gselements.1.2.73

Stagno, V., Ojwang, D. O., McCammon, C. A., & Frost, D. J. (2013). The oxidation state of the mantle and the extraction of carbon from Earth's interior. *Nature*, **493**(7430), 84–88. https://doi.org/10.1038/nature11679

Tang, W., Sanville, E., & Henkelman, G. (2009). A grid-based Bader analysis algorithm without lattice bias. *Journal of Physics: Condensed Matter*, **21**(8), 084204. https://doi.org/10.1088/0953-8984/21/8/084204
Thompson, E. C., Campbell, A. J., & Tsuchiya, J. (2021). Elastic properties of the pyrite-type FeOOH-AlOOH system from first-principles calculations. *Geochemistry, Geophysics, Geosystems*, 22(5), e2021GC009703. https://doi.org/10.1029/2021GC009703

Wicks, J. K., Jackson, J. M., & Sturhahn, W. (2010). Very low sound velocities in iron-rich (Mg,Fe)O: Implications for the core-mantle boundary region. *Geophysical Research Letters*, 37(15), L15304. https://doi.org/10.1029/2010GL043689

Williams, Q., & Garnero, E. J. (1996). Seismic evidence for partial melt at the base of Earth's mantle. *Science*, 273(5281), 1528–1530. https://doi.org/10.1126/science.273.5281.1528

Wohl, K. (1946). Thermodynamic evaluation of binary and ternary liquid systems. *Transactions of the American Institute of Chemical Engineers*, 42, 215–249.

Wohl, K. (1953). Thermodynamic evaluation of binary and ternary liquid systems. *Chemical Engineering Progress*, 49, 218–221.

Wood, B. J., Li, J., & Shahar, A. (2013). Carbon in the core: Its influence on the properties of core and mantle. *Reviews in Mineralogy and Geochemistry*, 75(1), 231–250. https://doi.org/10.2138/rmg.2013/75.8

Xu, M., Jing, Z., Bajgain, S. K., Mookherjee, M., Van Orman, J. A., Yu, T., & Wang, Y. (2020). High-pressure elastic properties of dolomite melt supporting carbonate-induced melting in deep upper mantle. *Proceedings of the National Academy of Sciences of the United States of America*, 117(31), 18285–18291. https://doi.org/10.1073/pnas.2004347117

Yu, M., & Trinkle, D. R. (2011). Accurate and efficient algorithm for Bader charge integration. *The Journal of Chemical Physics*, 134(6), 064111. https://doi.org/10.1063/1.3553716

Zhang, Y., & Yin, Q. Z. (2012). Carbon and other light element contents in the Earth's core based on first-principles molecular dynamics. *Proceedings of the National Academy of Sciences of the United States of America*, 109(48), 19579–19583. https://doi.org/10.1073/pnas.1203826109

Zhang, Z., Stixrude, L., & Brodholt, J. (2013). Elastic properties of MgSiO₃-perovskite under lower mantle conditions and the composition of the deep Earth. *Earth and Planetary Science Letters*, 379, 1–12. https://doi.org/10.1016/j.epsl.2013.07.034

Zhao, G., Mu, H. F., Tan, X. M., Wang, D. H., & Yang, C. L. (2014). Structural and dynamical properties of MgSiO₃ melt over the pressure range 200-500 GPa: Ab initio molecular dynamics. *Journal of Non-Crystalline Solids*, 385, 169–174. https://doi.org/10.1016/j.jnoncrysol.2013.11.024