The distribution of volatile elements during rocky planet formation

Core segregation and atmosphere formation are two of the major processes that redistribute the volatile elements-hydrogen (H), carbon (C), nitrogen (N), and sulfur (S)-in and around rocky planets during their formation. The volatile elements by definition accumulate in gaseous reservoirs and form atmospheres. However, under conditions of early planet formation, these elements can also behave as siderophiles (i.e., iron-loving) and become concentrated in core-forming metals. Current models of core formation suggest that metal-silicate reactions occurred over a wide pressure, temperature, and compositional space to ultimately impose the chemistries of the cores and silicate portions of rocky planets. Additionally, the solubilities of volatile elements in magmas determine their transfer between the planetary interiors and atmospheres, which has recently come into sharper focus in the context of highly irradiated, potentially molten exoplanets. Recently, there has been a significant push to experimentally investigate the metal-silicate and magma-gas exchange coefficients for volatile elements over a wide range of conditions relevant to rocky planet formation. Qualitatively, results from the metal-silicate partitioning studies suggest that cores of rocky planets could be major reservoirs of the volatile elements though significant amounts will remain in mantles. Results from solubility studies imply that under oxidizing conditions, most H and S are sequestered in the magma ocean, while most N is outgassed to the atmosphere, and C is nearly equally distributed between the atmosphere and the interior. Under reducing conditions, nearly all N dissolves in the magma ocean, the atmosphere becomes the dominant C reservoir, while H becomes more equally distributed between the interior and the atmosphere, and S remains dominantly in the interior.


Introduction
Hydrogen, carbon, nitrogen, and sulfur are abundant volatiles in the atmospheres of terrestrial planets in the Solar System and are expected to similarly dominate the gaseous envelopes of extrasolar rocky planets (Wordsworth and Kreidberg, 2022).Depending on the prevailing thermodynamic and compositional properties, these elements can be present in atmospheres in oxidized (e.g., H 2 O, SO 2 , CO 2 ) or reduced forms (e.g., H 2 , H 2 S, CO) (Gaillard, et al., 2022).A significant proportion of a planet's volatile elements can be locked up in deep reservoirs such as cores and mantles whose bulk compositions are largely determined during early planet formation (Rubie, et al., 2007).Atmospheric outgassing (magma-gas exchange) can occur simultaneously with core-mantle differentiation while planets are hot enough to sustain magma oceans (Elkins-Tanton, 2012;Hirschmann, 2012).Prior to the outgassing of their secondary atmospheres, large rocky planets (>0.2 M E ) that formed early likely accreted primary H and/or He nebular atmospheres that would have interacted with their early magma oceans (Sasaki, 1990;Zahnle, et al., 2010).Although Earth lost its primordial atmosphere (Catling, et al., 2001), studies of exoplanetary systems have revealed that some rocky planets can retain such envelopes throughout their evolution (Misener and Schlichting, 2021).Atmospheres can also be strongly modified by impact-induced shock heating and erosion (Matsui and Abe, 1986;Zahnle, et al., 1988;Davies, et al., 2020;Zahnle, et al., 2020).However, from a geoscientific perspective, the two main volatile redistribution mechanisms are core formation and magma outgassing, which we discuss here with particular emphasis on recent measurements.Volatiles retained within rocky planets' mantles can reside as various mineral or fluid components and participate in geological and biochemical processes (Dasgupta and Hirschmann, 2010;Li, et al., 2013;Mikhail and Sverjensky, 2014;Armstrong, et al., 2015;Hirschmann, 2018).Significant amounts of volatile elements, along with other light elements (e.g., O and Si), are also inferred to be present in the cores of Earth and other rocky planets such a Mars, where they lower the melting temperature and density of iron, leading to partially to fully molten cores capable of sustaining a dynamo (Badro, et al., 2014;Nimmo and Schubert, 2015;Brennan, et al., 2020;Lv and Liu, 2022).
The Earth's core and mantle inherited their elemental abundances during Earth's differentiation.Core-mantle segregation involved chemical reactions between Fe and silicates in a reduced (Fe-alloy saturated) magma ocean (Wade and Wood, 2005;Wood, et al., 2006).Geochemical evidence-the abundance of refractory siderophile elements such as Ni, Co, W and Mo in mantle peridotite-coupled with metal-silicate partitioning experiments have led to inferences about the average conditions of core formation (Ringwood, 1977).The earliest framework for understanding core-mantle differentiation was provided by single-stage models in which metals accumulated at the bottom of the magma ocean before descending into the core (Stevenson, 1981;Stevenson, 1988;Karato and Murthy, 1997;Solomatov, 2007).These models assumed that the entire core and magma ocean equilibrated at a common pressure-temperature-oxygen fugacity and composition (P-T-f O 2 -X) (Jones and Drake, 1986) and led to valuable insights on how redox state and planetary size affect core and mantle chemistry (Righter and Drake, 1996;Wade and Wood, 2005;Corgne, et al., 2009).Early experimental studies of elemental partitioning during core-mantle-differentiation indicated that the mantle's refractory siderophile abundances require metal-silicate partitioning at pressures from 25 to 30 GPa, equivalent to a magma ocean depth of 800 km (Li and Agee, 2001).However, later experiments found that average conditions of 45-55 GPa (~1,500 km depth) and ≥3500 K (Bouhifd and Jephcoat, 2011;Siebert, et al., 2012;Fischer, et al., 2015) could better explain the mantle's Ni and Co abundances.Large impacts may have provided the energy for such deep magma oceans and their associated episodes of core formation (Tonks and Melosh, 1993;Davies, et al., 2020), although rapid pebble accretion (Johansen, et al., 2021) could also have led to high temperatures and mantle melting.
Recent experimental studies done on different groups of elements including highly siderophiles and volatiles have indicated that one set of P-T-f O 2 -X condition cannot account for all the relevant geochemical observables such as the bulk silicate Earth (BSE) contents in highly siderophile and volatile elements.More integrated approaches have been developed which consider the dispersal of impacting materials and degree of chemical equilibration before segregating into a core and a mantle (Rubie, et al., 2003;Rubie, et al., 2007;Deguen, et al., 2014).These multistage approaches can account for some geochemical observables such as the moderately siderophile elemental abundances in the mantles of smaller planetary bodies, which required metal-silicate reactions at lower pressure-temperature conditions (Righter and Drake, 1997;Steenstra, et al., 2016).Though core formation in terrestrial planets is generally thought to occur at pressures and temperatures where metal and silicate phases are largely immiscible, other large-scale processes such as the segregation of a sulfide matte may accompany or follow core segregation (Rubie, et al., 2016;Steenstra, et al., 2020).Additionally, density estimates for rocky exoplanets (Otegi, et al., 2020) and the bulk chemical abundances of disintegrated planetesimals in the atmospheres of polluted white dwarfs (Bonsor, et al., 2020;Bonsor, et al., 2022) indicate that the process of planetary differentiation could be common in planetary systems.
The failure of single-stage homogenous models to find a unique P-T-f O 2 -X condition that can explain all geochemical observables suggests not only that Earth formed in stages and experienced disequilibria, but also that accreting material might have changed over time (Wade, et al., 2012).Indeed, N-body simulations of planetary accretion (Chambers and Wetherill, 1998;O'Brien, et al., 2014) imply that fully formed rocky bodies experience many large impacts during their growth and that radial mixing of inner and outer Solar System materials occurs concomitantly.This heterogeneous framework allows for many episodes of core formation to occur at varying pressures, temperatures, and compositions (Rubie, et al., 2003;Albarede, 2009;Schönbächler, et al., 2010).Crucially, heterogeneous models driven by the mass delivery scenarios of N-body simulations tend to favor the delivery of more volatile-rich materials, ultimately sourced from the outer Solar System, later in the accretion process when the potential for extreme P-T metal-silicate reactions is at its maximum (Halliday, 2013;O'Brien, et al., 2014;Rubie, et al., 2016).However, analyses of small rocky bodies and coupled to astrophysical and geochemical models suggest that volatile elements were present in the inner Solar System well before the later stages of planetary formation (Bar-Nun and Owen, 1998;Busemann, et al., 2006;Marty, 2012;Halliday, 2013;Alexander, 2017;McCubbin and Barnes, 2019;Grewal, et al., 2021c;Deligny, et al., 2021;Grewal, 2022;Grewal and Asimow, 2023).
Another consequence of heterogeneous accretion models is the implication that rocky planets can experience several magma ocean stages throughout their evolution (Elkins-Tanton, 2012).Magmagas interactions during magma ocean stages lead to volatile outgassing at temperature and compositional conditions related, but not equal, to those inferred for core formation (Zahnle, et al., 2010;Hirschmann, 2016;Dasgupta and Grewal, 2019;Lichtenberg, 2022).Volatiles are, by their nature, prone to partitioning into the gas; accordingly, an increasing effort is being invested into developing models to account for the atmospheres of accreting planets as large volatile reservoirs, effectively limiting incorporation into magma oceans or partitioning into cores during planetary accretion.Factors such as oxygen fugacity and magma ocean depth and volatile content strongly affect the resulting atmosphere and residual mantle compositions (Hirschmann, 2016;Grewal, et al., 2021b;Gaillard, et al., 2021;Jackson et al., 2021).Models are now being developed to account for the distribution of volatile elements into atmospheres with varying thermodynamic and compositional properties (Salvador, et al., 2017;Nikolaou, et al., 2019;Lichtenberg, et al., 2021b;Bower, et al., 2022).Although solubility establishes the equilibrium distribution of a volatile element between a primordial atmosphere and its underlying magma ocean, increasing attention is being paid to how kinetic processes, including diffusion and convective stirring within planetary interiors, may affect the outgassing and ingassing of magma oceans (Salvador and Samuel, 2023).
Partition coefficients and solubilities are the fundamental data needed to parameterize core-mantle differentiation and magma ocean-atmosphere exchange models.These parameters have been incorporated into models that are used to understand how volatile elements are fractionated and cycled among major planetary reservoirs-i.e., cores, mantles, and atmospheres-for a range of bodies, including Earth, Mars, Venus, and differentiated planetesimals such as the parent bodies of iron meteorites (Hirschmann, 2012;Wordsworth, 2016;Gaillard, et al., 2021;Jackson et al., 2021;Li, et al., 2021;Grewal, et al., 2022a;Lichtenberg, et al., 2022;Grewal and Asimow, 2023).Although early models were tuned and benchmarked to the Earth, the parameterization made available by recent partitioning and solubility measurements allows the models to be adapted to a wide range of other planetary formation scenarios.Thousands of rocky exoplanets have been discovered by surveys such as the Kepler Mission (Batalha, 2014) and the Transiting Exoplanet Survey Satellite (Kaltenegger, et al., 2019), and atmospheres of super-Earths and sub-Neptunes are now being characterized by the James Webb Space Telescope (JWST) (Mansfield, et al., 2019;Ding and Wordsworth, 2022).In particular, in its first year of observations, JWST has demonstrated its capabilities to constrain the presence or absence of large volatile envelopes around shortperiod rocky exoplanets (Greene, et al., 2023;Ih, et al., 2023;Moran, et al., 2023;Zieba, et al., 2023).After reasonable observational times over the next few observational cycles, JWST will be able to constrain models of the volatile contents and redox states of individual rocky exoplanets (Kempton, et al., 2023;Piette, et al., 2023).It is therefore crucial to continue efforts to model the formation and evolution of rocky planets in these exoplanetary systems (Lichtenberg, et al., 2019;Kite, et al., 2020;Lichtenberg, et al., 2021b;Dorn and Lichtenberg, 2021;Gaillard, et al., 2021;Kite and Schaefer, 2021;Lichtenberg and Krijt, 2021;Schlichting and Young, 2022;Wolf, et al., 2022).Laboratory parameterizations of volatile partitioning and solubilities are increasingly being incorporated into the models which calls for an effort to contextualize the measurements.In this review, we aim to do just that, such that the synergy between laboratory investigations and astronomical characterizations of hot exoplanets can be maximized.

Parameterizing experiments
In addition to gravitational segregation, core-mantle differentiation in terrestrial planets can be thought of in terms of the chemical separation of core Fe from mantle oxides and silicates.Thus, laboratory experiments study the chemical transfer of relevant species between metal alloys and silicate/oxide mixtures.The parameterization of metal-silicate partitioning reactions and experimental results has been covered extensively in several recent works (Wood, et al., 2006;Corgne, et al., 2008a;Blanchard, et al., 2017;Huang and Badro, 2018;Chidester, et al., 2022), and only a brief overview is given here.
The equilibrium partitioning of an element between metal and silicate can be thought of as an exchange reaction between molten metal (M) and oxide: where n is the valence of element M. At equilibrium, the partition coefficient D is the ratio of the concentration (here in terms of molar fraction X, but can also be in terms of weight percent, wt.%) of an element in the metal phase to that in the silicate phase: The partition coefficient is related to exchange coefficient of the reaction, K D , by normalizing by the partition coefficient for Fe: That is, K D is the ratio of the partition coefficient of element M to that of iron and is independent of the Fe concentrations during the reaction.
K D is related to the equilibrium constant K of reaction (1), which can be defined in terms of K D and the activity coefficient (γ) of the metal and silicate phases (E.g., Chidester, et al., 2022): Note that the oxygen fugacity of the experiments can be expressed in terms of the Fe-FeO equilibrium: Fe alloy + 1/2 O 2 = FeO silicate melt , i.e., the Iron-Wüstite (IW) redox buffer (E.g., Corgne, et al., 2008a;Campbell, et al., 2009;Chidester, et al., 2022): Where a in this case refers to the activity.
For reaction (1), K is related to the Gibbs free energy of the reaction (E.g., Chapter 11 of Ganguly, 2008): where R is the gas constant, T temperature, P pressure, and ΔHo, ΔSo, and ΔVo are the changes in enthalpy, entropy, and volume of reaction (1), respectively.The natural log of the partition coefficient, log D(M) can then be expressed in terms of the above variables and constants: The constants a, b, c, and d are determined using multivariate least-squares linear regression of the measured Ds, Ps, Ts, and compositions.The parameterization can also include activity terms for the light components in the alloy and silicate melt polymerization.The expression obtained from the regression can then be used to determine partition coefficients over a wide P-T-f O 2 -X space, although care must be taken when regressing because a, b, and c may significantly change with variations in experimental conditions (Walter and Cottrell, 2013).By coupling these functions to core formation models, we can calculate the redistribution of elements between cores and mantles during differentiation.

Core formation models
In the single-stage core formation model, a planet is assumed to differentiate instantaneously and the re-distribution of a chemical species between core and mantle is determined at one set of conditions.Using the formalism of (Rudge, et al., 2010), the mass of species c in the body is conserved as follows: where x b is the bulk composition, x me is the concentration in the mantle embryo, x ce is the concentration in the core of the embryo, and F is the body's core mass fraction, which varies in the Solar System from 0.25 (Mars) to 0.8 (Mercury).The effective partitioning D x in the embryo is defined at one set of P-T-X conditions and assumed constant in single-stage models: In multi-stage core formation models, partition coefficients evolve with the magma ocean depth to account for multiple or continuous differentiation events (Rubie, et al., 2003).The conservation of mass of species x into the accreting planet's core and mantle can integrated numerically from Eqs 9, 10 following from Rudge et al. (2010).Similar formalisms are also developed elsewhere (Rubie, et al., 2011): where x m and x c are the concentrations of species x in the mantle and core respectively of the resultant body.M is the mass of the Earth which changes with time t.D c evolves with P-T-f O 2 -X conditions as the planet accretes.This framework can be adapted to different evolution and accretion paths, including accreting bodies of different sizes and compositions and varying degrees of metalsilicate equilibration between impactors and the magma ocean.
The first term on the right-hand side of Eq. 10 represents material from the embryo's mantle that is added to the Earth's mantle, whereas the second term represents the mantle that reequilibrated with the impactor core.Partial equilibration due to inefficient mixing between large impactors towards the end of accretion (Rubie, et al., 2003;Dahl and Stevenson, 2010) is introduced by the parameter k, which is the fraction of the embryo's core that equilibrates with the planet's mantle.The expressions on the right-hand side of Eq. 11 represent the two paths that the metal can take to the core: (1 − k) Fc ce is the fraction of the embryo's core that is added directly to the planet's core whereas kFD x c m is the fraction that equilibrates with the mantle.Alternatively, partial equilibration can be implemented using a turbulent mixing parameterization derived from fluid dynamic experiments (Deguen, et al., 2014), which is identical in terms of mass conservation to the approach of Rudge et al. (2010), but the single fraction k is replaced by an efficiency factor ξ: where Δ is a metal dilution ratio that is related to the relative densities of the equilibrating metal and silicate (Deguen et al., 2014).

Hydrogen
Hydrogen is the most cosmochemically abundant of the volatile elements and is found in significant concentrations in a wide range of meteorites (Alexander, 2019a;Alexander, 2019b;McCubbin and Barnes, 2019;Lewis, et al., 2022;Peterson, et al., 2023).The D/H ratios of various terrestrial, meteoritic, and planetary reservoirs have been crucial in constraining the origins of hydrogen on Earth, which measurements indicate could be due to contributions from different reservoirs accreted at different times throughout the Solar System's history (Albarede, et al., 2013;Halliday, 2013).Studies of hydrogen metal-silicate partitioning behavior (and solubility, discussed in section 3.1.1)have recently become important in the context of sub-Neptune exoplanets, which are thought to have retained thick primodrial H and/or He envelopes (Fulton, et al., 2017;Lichtenberg, et al., 2021a;Dorn and Lichtenberg, 2021;Rogers and Owen, 2021;Schlichting and Young, 2022).There is also isotopic evidence that the Earth and Mars could have had an early H-rich atmosphere (Tian, et al., 2005;Sharp, 2017;Pahlevan, et al., 2022).Hydrogen could also have been ingassed into Earth's magma ocean if it was in contact with such an atmosphere, and may have ultimately partitioned into the core (Wu, et al., 2018;Young, et al., 2023).Likewise, the cores of H-rich exoplanets could also have inherited nebular H, which would have influenced their total water abundances (Kimura and Ikoma, 2022).
Though hypothesized to be a low-density core component, its tendency to diffuse and exsolve makes measuring the metal-silicate partitioning of hydrogen difficult, and data remain sparse (See compilation in Figure 1); it readily forms FeH x at pressures above 5 GPa, but decomposes at ambient conditions (Badding, et al., 1991;Okuchi, 1997).Okuchi (1997) indirectly (i.e., by estimating the volume of gas bubbles) determined D H at 7.5 GPa and 1,200 °C-1,500 °C in large-volume press experiments to be ~1 on average and proposed that H could be a primary light element in Earth's core.More recently, large-volume press experiments have yielded new estimates for D H up to 20 GPa, as constrained by a combination of analytical techniques including electron probe microanalysis (EPMA), nuclear microprobe resonance (NMR), and electron recoil detection analysis (ERDA) inferred from the latter experiments (Clesi, et al., 2018;Malavergne, et al., 2019).
Resultant D H values from these two studies are lower than unity (lithophile behavior), but thermodynamic modeling also suggested that higher pressure conditions could make H more siderophile.It is important to note that these studies measured sample compositions after quenching, and it is possible that hydrogen diffused out of their quench products, leading to large biases on D H . Carbon saturation in the metals of the latter experiments could also have contributed to lowering D H .
Measurements in laser-heated DAC experiments indirectly inferred H content from the unit cell expansion with in situ X-ray diffraction (Tagawa, et al., 2021) and found H to be siderophile, with D H > 29 at 30-60 GPa and 3,100-4600 K, in agreement with theoretical calculations (Yuan and Steinle-Neumann, 2020).When coupled with a core formation model, those values imply that up to 0.6 wt.% H could have been incorporated into the Earth's core and >0.15 wt.% into the cores of planets more massive than 0.1M E if water was accreted during the main stage of accretion.

Carbon
The presence of carbon in Earth's and other planetary cores has been inferred from its affinity for iron at ambient conditions and its depletion in the BSE relative to carbonaceous chondrites (CIs) (Hirschmann, 2016).The results of metal-silicate partitioning experiments at pressures up to 15 GPa in large-volume presses (Dasgupta, et al., 2013;Chi, et al., 2014;Malavergne, et al., 2019) indicate that although carbon is a highly siderophile (D C up to 10 5 ) TABLE 1 Recommended values of volatile element partition coefficients (D x ) over a range of P-T-f O 2 -X conditions relevant to different magma ocean depth/ planetary size and redox conditions.Functional forms from the literature are used to calculate these partition coefficients, see Supplementary Appendix SA (H (Malavergne, et al., 2019;Tagawa, et al., 2021), N (Grewal, et al., 2021b), C (Malavergne, et al., 2019;Fischer, et al., 2020;Blanchard, et al., 2022), S (Boujibar, et al., 2014;Suer, et al., 2017).Note that these recommended values based on specific set of conditions and different values will be obtained if conditions/variables in the functions are changed.Ds are rounded off to the nearest significant decimal.Uncertainties on the Ds are dependent on the uncertainties in the fits to the functional forms and can be as large as ± 20%. at lower pressures conditions (~1 GPa), it will partition significantly less into Fe as pressure increases (e.g., D C ≈ 10 at 15 GPa) (See compilation in Figure 2).Models based on these results, suggest that the cores of smaller planetesimals could retain wt.% levels of carbon, whereas those of larger bodies would be less enriched (<1 wt.% C) while their mantles would become progressively more C-rich during accretion (Kuwahara, et al., 2021).Though D C can increase as fO 2 decreases from ΔIW −1 to −3 (Malavergne, et al., 2019), under highly reducing conditions, Si will partition strongly into iron, limiting C dissolution in cores of very reduced bodies.This chemical trend has led to the suggestion that the Ureilite parent body or Mercury could have become carbon-saturated, possibly leading to a graphite-rich crust (Vander Kaaden and McCubbin, 2015;Keppler and Golabek, 2019;Steenstra and van Westrenen, 2020).The presence of sulfur and nitrogen in core materials can also affect the partitioning of carbon into the Fe (Grewal, et al., 2019b;Jackson et al., 2021), and immiscible C-and S-rich layers could have been contemporaneous in the cores of some planetesimals (Corgne, et al., 2008b).

Planetesimals
Experimental results at higher P-T conditions in laser-heated DACs show that carbon has a lowered affinity for iron at the conditions of a deep terrestrial magma ocean (Fischer, et al., 2020;Blanchard, et al., 2022).When coupled with models of Earth's core formation and accretion, these partition coefficients imply that less than 1 wt.% C could be present in the cores of Earth and similar-sized bodies.Complementary mantle compositions could be hundreds of ppm, depending on the timing of carbon accretion.

Nitrogen
Nitrogen is the most abundant component of Earth's atmosphere and its presence in meteoritic material indicates that Hydrogen partition coefficients (log D H ) plotted as a function of oxygen fugacity, parameterized here in terms of the iron-wüstite (ΔIW) buffer (A) and experimental temperature (B).Data sources: (Okuchi, 1997;Clesi, et al., 2018;Malavergne, et al., 2019;Tagawa, et al., 2021).
it was present in early planetary building blocks (Marty, 2012;Grewal, et al., 2021b;Grewal, 2022;Grewal and Asimow, 2023).Its depletion in Earth's mantle relative to C (Marty, 2012) suggests that N could have been incorporated into the Earth's core.Furthermore, it has been noted from early solubility experiments that N is increasingly incorporated into iron with increasing pressures of up to 7 GPa (Roskosz, et al., 2013;Speelmanns, et al., 2018).More recent metal-silicate partitioning experiments in large-volume presses (Dalou, et al., 2017;Grewal, et al., 2019a;Dalou, et al., 2019) found N to be siderophile (10 < D N < 31) at 1-7 GPa and up to 1800 °C (See Figure 3); its partitioning behavior being sensitive to P, T, f O 2 , and the presence of light elements such as S in Fe (Grewal, et al., 2019a).The sensitivity of N partitioning to oxygen fugacity suggests that it dissolves as nitride (N 3− ) in silicate melts and as neutral N in core metals (Dalou et al., 2019;Grewal et al., 2020), and this change in oxidation state makes N partitioning relatively sensitive to oxygen fugacity compared to the other volatile elements.
Nitrogen is poorly soluble in magmas at oxygen fugacities around IW (Libourel et al., 2003).Thus, magmatically active bodies can lose a large fraction of their N to their atmosphere.Combined differentiation, outgassing, and accretion models show that protoplanets that differentiated early could have been depleted in N (Grewal, et al., 2021b).However, if differentiation occurred later at the embryonic stage, nitrogen reservoirs could have been maintained within cores while the mantles would have remained depleted, a scenario which could explain the N signature of the terrestrial mantle (Grewal, et al., 2021c).Limited laser-heated DAC data demonstrate that N remains siderophile up to 26 GPa and ~3500 K (Jackson et al., 2021), further suggesting that it could be sequestered into larger planetary cores.This behavior is qualitatively different from C and S, which both become less siderophile at increased pressures and temperatures.This change in siderophile behavior among the volatile elements has significant implications for the S-C-N reservoirs on planetary embryos and larger planets (Jackson et al., 2021).The measurements compiled in Figure 3 highlight that D N increases with increasing f O 2 .Obtaining measurements at higher P-T conditions will be valuable to further investigations of N incorporation into the cores of large planets.

Sulfur
Sulfur's affinity for iron at ambient conditions and its presence in iron meteorites imply that it could be a component of protoplanets and planetary cores (Dreibus and Wanke, 1985;Chabot, 2004;Gounelle and Zolensky, 2014).Combined with its cosmochemical abundance, these factors have been used to argue that S is one of the light alloying components in Earth's and other planetary cores.Its incorporation into Fe has been studied extensively in an attempt to explain the density deficit in Earths's core (i.e., relative to that of pure iron) and the properties of other planetary cores (Fei, et al., 1995;Morard, et al., 2013;Boujibar, et al., 2020;Brennan, et al., 2020).S is also known to affect the behaviors of other chemical species in various geological settings, particularly chalcophile elements (Jana and Walker, 1997;Mahan, et al., 2018) and C (Li, et al., 2016b;Tsuno, et al., 2018).Experiments in largevolume presses up to 30 GPa (Li and Agee, 2001;Boujibar, et al., 2014); showed the increasing siderophile tendency of S with increasing pressure, lending support to its incorporation into the cores of planetesimals and planetary embryos.For example, Mars, which is the size of a planetary embryo and more volatile rich than Earth, could contain up to ~20 wt.% S in its core (Brennan, et al., 2020) while it has been speculated that Mercury's core could be surrounded by an FeS layer (Namur, et al., 2016).Experiments in laser-heated DACs ~29-100 GPa, up to 5300 K (Suer, et al., 2017;Mahan, et al., 2018;Jackson et al., 2018;Chidester, et al., 2022); found D S values that are an order of magnitude lower on average than large-volume press results (ranging from less than 10-100), suggesting that sulfur's siderophile tendency does increase with pressure, but is also strongly suppressed at high temperatures (see compilation in Figure 4).When the results of laser-heated DAC experiments are included in fitting the functional form for log D S , the entropy term is an order of magnitude larger than if the fit included only the large-volume pressure cell results.This is generally
the case with laser-heated DAC studies of other elements.Redox and the chemical speciation of S have also been shown to play an important role in its partitioning behavior (Mavrogenes and O'Neill, 1999;Chidester, et al., 2022).
When the functional form for D S is incorporated into models of Earth's core formation (Boujibar, et al., 2014;Suer, et al., 2017), less than 2 wt.% S can partition into the core, in agreement with cosmochemical (McDonough and Sun, 1995) and geophysical constraints (Badro, et al., 2014).However, depending on the initial bulk sulfur content and/or the timing of sulfur accretion and the efficiency of core-mantle equilibration, more sulfur could be incorporated into planetary cores.Predicted sulfur compositions of the Earth's mantle could range from a few hundred to thousands of ppm also depending on the timing of sulfur accretion (Rubie, et al., 2016;Suer, et al., 2017), overlapping with geochemical observations suggesting that the BSE contains ~200 ppm S (McDonough and Sun, 1995;Wang and Becker, 2013).It is possible that magma ocean crystallization leads to an increase in S concentrations (Rubie, et al., 2016).Above the sulfur capacity at sulfide saturation (SCSS) limit, precipitation and segregation of a sulfide-rich matte could occur on a planet-wide scale in both small and large bodies (O'Neill and Mavrogenes, 2002;Steenstra, et al., 2020), leading to a late pulse of core formation if the sulfide phase can mobilize to the core-mantle boundary.

Experimental limitations
In the last decade, advances in experimental and analytical techniques have enabled high-quality measurements of metalsilicate partitioning coefficients and the solubilities of volatile elements.Nonetheless, the interpretation of some of these measurements remain controversial.
Although not specific to volatile-element partitioning, ambiguity surrounds the interpretation of quench textures in both large-volume press and DAC experimental samples.Dendritic textures and overgrowths of quench rims in run products of large-volume press experiments lead to difficulty in determining stable liquid compositions.Likewise, the origin of metallic inclusions in the silicate portions of DAC run products has been controversial.Figure 5 shows examples of a large volume and laser-heated DAC experiment post quench.Large metal inclusions are usually excluded from analyses of silicate melt compositions, but smaller nanoscale inclusions (nanonuggets) are difficult to avoid, particularly with EPMA which integrates over several to tens of cubic microns.To assess this issue, NanoSIMS has been used to obtain highly resolved measurements for S and C (Suer, et al., 2017;Fischer, et al., 2020;Blanchard, et al., 2022).TEM measurements have also been performed to help resolve this controversy in recent studies (Fischer, et al., 2015;Suer, et al., 2017;Suer, et al., 2021), which suggest that nanoparticles form upon rapid quench but were dissolved during equilibrium melting.However, if nano-metal inclusions are contaminants, they could lower measured partition coefficients for the volatiles and other elements, significantly impacting the results of models that utilize these measurements.
De-volatilization and the loss of volatile species during quenching can lead to large uncertainties on experimental measurements.The degree to which a high pressure-high temperature assemblage is preserved depends on the quench rate of the experimental apparatus: example quench rates are 175 °C/s in piston-cylinder experiments, ~800 °C/s in multi-anvil experiments, and thousands of degrees per microsecond in DAC experiments.Magnesium-rich silicates (analogs of mafic terrestrial magma oceans) are difficult to quench to a homogenous glass, and some works have thus used basaltic and andesitic silicate compositions (Suer, et al., 2017).Quenching is particularly important for hydrogen, which is known to escape from silicate melts and metals upon recovery at ambient conditions (Okuchi, 1997).Therefore, it has been particularly difficult to obtain reliable partitioning measurements for hydrogen.Recently, the combination of EPMA, NMR and ERDA have been used successfully to measure partition coefficients for H and C (Clesi, et al., 2018;Malavergne, et al., 2019).
The large concentrations of volatile elements used in the iron alloys of some experiments can also affect partition coefficients.Such concentrations are typically much higher than those observed in nature and can be in violation of Henry's law.Thus, the application of such partition coefficients in models may introduce additional, and difficult to quantify, uncertainties (Kuwahara, et al., 2019;Grewal, et al., 2021a;Grewal, et al., 2022b;Shi, et al., 2022).
Statistical regressions of measurements are used to obtain constants for the functional description of partitioning behavior.Differences have been noted between partitioning behaviors in large-volume press and DAC experiments.These have been attributed to changes in chemical behaviors in the different P-T regimes probed by these two types of apparatuses.This effect is perhaps best documented for S and can be seen in the different magnitudes of the regression constants obtained across different experimental techniques (Boujibar, et al., 2014;Suer, et al., 2017).Thus, a single functional form might not be sufficient to fully describe datasets spanning large P-T ranges.In addition, linear functional relationships might not fully describe the convolution of the IW redox buffer with the behaviors of other elements involved, and further modeling efforts are needed to deconvolve these effects.For these reasons, it is not recommended that the functional forms be extrapolated beyond the ranges of the measurements on which they are based.

Volatile element solubilities and magma-atmosphere interactions
Volatile elements distribute themselves between metal, silicate, and gas during planetary formation (Chao, et al., 2021;Lichtenberg, et al., 2022).Metal reacts with magma before its segregation to the core, whereas gas and magma continuously exchange volatiles at the atmosphere-magma ocean interface.Thus, it is expected that planetary bodies continuously redistribute their volatile element budgets as they grow and undergo new periods of core formation and mantle melting.In the previous section, our focus was on the sequestration of volatile elements within planetary cores via metal segregation in magma oceans.In this section, however, our focus is to understand how the presence of an atmosphere can affect the distribution of volatile elements during core formation, because any gaseous molecule present in an atmosphere subtracts from the budget available to metal-silicate partitioning.

Gas solubilities in magma oceans
Solubility is perhaps the most central parameter governing the chemical interaction of a primordial atmosphere and its underlying magma ocean.More specifically, solubility quantifies the relationship between the fugacity of a gas species and its corresponding concentration in a condensed phase.This can be expressed as: where X i is the concentration of element i dissolved in the condensed phase, f α i is the fugacity of element i in the gas (or fluid) raised to an exponent that relates to the stoichiometry of the dissolution reaction(s), and S i is the solubility of element i in the condensed phase.Solubility depends on the temperature and pressure of the reaction as well as the compositions of the condensed and gaseous phases.In our example, the condensed phase is the magma.Fugacity can be further deconvolved as: where γ i and P i are the fugacity coefficient and the partial pressure of element i, respectively.The high temperatures and moderate pressures in atmospheres suggest that the gases present will only have moderate deviations from ideality, and fugacity coefficients are often assumed to be unity, such that fugacity and partial pressure are equal.Eqs 13, 14 reveal the overarching role of atmospheric pressure on atmosphere-magma interactions, as it is partial pressures that drive the dissolution of volatile species into the magma.
The stability of a gas species depends on the prevailing P-T-X conditions within the gas phase.For example, the stability of CO with respect to CO 2 can be expressed as: Le Chatelier analysis indicates that more oxidized conditions favor CO 2 stability compared to CO.In isolation, this reaction also predicts that CO 2 will be favored under lower temperature and higher pressure given that the reactants have higher entropy and volume, but all major gas species should be considered together when making specific predictions regarding how species concentration is affected by changing P-T conditions for any bulk gas composition.
The example above stresses the importance of oxygen fugacity in determining the stability of gas species.We note, however, that it is not precisely known how oxygen fugacity varies throughout a growing planet.It is known that, at depth, magma oceans are saturated in liquid Fe alloy and are therefore strongly reduced environments (well below the IW oxygen buffer, as inferred by the low residual FeO concentrations in the bulk silicate portions of planets), but isochemical decompression (or compression) of a magma can drive homogeneous reactions that cause autoreduction or auto-oxidation (Zhang, et al., 2018;Armstrong, et al., 2019).Moreover, magma ocean crystallization should lead to fO 2 variations given the different partitioning of ferric and ferrous Fe between minerals and melts.This implies that the oxygen fugacity associated with metal-silicate reactions need not be the same as that associated with gas-magma reactions, particularly because the magma ocean may become increasingly chemically isolated form the atmosphere as solidification progresses (Bower et al., 2022).However, our goal is not to review this dynamic and the potential resulting redox stratification within planetary magma oceans (Hirschmann, 2012).Rather, we introduce this concept to stress that a wide range of oxygen fugacities are likely to impact volatile element solubilities associated with planetary formation.We correspondingly organize our review of gas solubilities below by element, and within each element by oxygen fugacity.We also emphasize that our goal is not to provide an unabridged review of volatile element solubilities in magmas, but we simply seek to highlight the basic controls on volatile element solubilities and to provide reference solubility values that may serve to accelerate the development of magma ocean-atmosphere interaction models.

Hydrogen
When reacted with oxygen, hydrogen gas forms water as (similar to Reaction 15): Relative to H 2 , gaseous water molecules have a wide stability field under geological conditions, as demonstrated by the ratio of f H 2 O/f H 2 plotted as a function of ΔIW at 2273 K (Figure 7A).The crossover of f H 2 O/f H 2 ratio occurs near IW and is nearly independent of temperature.We note that we use thermodynamic data from the NIST WebBook to calculate equilibrium constants of reactions between gases.
Water vapor can dissolve into magma as either OH or H 2 O (molecular) following the reactions: Reactions 17 and 18 predict that low f H 2 O values should favor the dissolution of OH; however, with increasing OH concentration, relatively more H 2 O will dissolve as molecular H 2 O. Indeed, water solubility has been experimentally shown to scale with f H 2 O 0.5 at lower concentrations, but with f H 2 O at higher concentration.The transition from OH-to H 2 O-dominated water solubility typically occurs at the wt.% level (Stolper, 1982), which requires f H 2 O > 100 bars, whereas terrestrial planets are typically estimated to contain between 100 and 1,000 ppm H 2 O integrated over their silicate reservoirs, atmospheres, and oceans (Peslier, 2010;Halliday, 2013).This implies that only modest H 2 O pressures could have been present in primordial atmospheres (~1-10 bars on average), unless large amounts of H were lost from the atmospheres (Catling, et al., 2001).Moreover, the products in Reaction 18 are favored at higher temperatures (Nowak and Frontiers in Earth Science frontiersin.org10 Behrens, 1995), and, given that magma oceans exist at higher temperatures than modern magmatic systems, it is expected that OH is a major component of water solubility in growing planets, even under relatively reducing conditions (Figure 6).
Hydroxyl forms chemical bonds when dissolved in silicate melts, lowering melt viscosity and implying that water reacts with bridging oxygens in Reaction 18 to depolymerize the melt network.Magma ocean liquids are depolymerized, meaning that they have few bridging oxygens to support OH dissolution.Studies of water solubility in depolymerized liquids highlight that water also bonds with non-bridging oxygens, preferentially forming complexes with higher field strength cations (Mysen and Virgo, 1986;Xue and Kanzaki, 2004).It is therefore expected that water should dominantly dissolve as OH at the low fH 2 O conditions relevant to magma ocean liquids.A recent work suggests that peridotitic magmas contain relatively more OH than H 2 O compared to more polymerized magmas at equal total H 2 O contents, such that H 2 O solubility in peridotite liquids will remain dominated by OH up to 5 wt.% H 2 O (Bondar, et al., 2023).
Experimental determinations of water solubility in magma oceans have historically been precluded by difficulties associated with working in ultramafic systems, although recent advances in rapid-quench multi-anvil and laser levitation techniques are highly promising (Bondar, et al., 2023;Sossi, et al., 2023).In the absence of well-established data, we rely on models that allow for the P-T-X dependencies of water solubility determined within the realm of quenchable liquids to be extrapolated to the P-T-X realm of magma oceans.Several models have been developed (Newman and Lowenstern, 2002;Papale, et al., 2006;Iacovino, et al., 2021), and we use MAGMASAT (Ghiorso and Gualda, 2015) for our example because it is implemented in the popular MELTs algorithm and considers both OH and H 2 O dissolution.
with PH 2 O in bars.The fact that MAGMASAT predicts that H2O solubility scales close to PH 2 O 0.5 is consistent with OH dominating solubility up to 0.6 wt.% H2O dissolved in peridotitic magma.Such predictions require significant extrapolation because most water solubility determinations are limited to basaltic (or move evolved) compositions and typical magmatic temperatures (<1673 K) and should therefore be treated with caution.Nevertheless, it is clear that, relative to other gases (see following subsections), water vapor remains highly soluble in ultramafic magma, even under magma ocean conditions.Support for MAGMASAT model prediction comes from recent results from a laser levitation apparatus study (Sossi, et al., 2023).This apparatus enables quenching of ultramafic glasses reacted with controlled f H 2 O at temperatures that are directly applicable to magma oceans.Current data collected in systems with f H 2 O/ f H 2 > 1 yield a similar solubility relationship to that predicted above: assuming a molar absorption coefficient (ε 3550 ) of 5.1 m 2 /mol.Under more reducing conditions, water vapor is progressively destabilized to produce H 2 .Spectral evidence demonstrates that H 2 gas dissolves as molecular H 2 in magmas (Luth, et al., 1987).As a neutral, non-polar molecule, H 2 is thought to dissolve into the ionic porosity of the melt structure, following the well-documented behavior of noble gases (Carroll and Stolper, 1993).The radius of H 2 is similar to that of helium, and indeed, experiments demonstrate that H 2 solubility is also similar to that of helium, at least for basaltic melts (Hirschmann, 2012).The solubility of helium, or hydrogen, has not be experimentally investigated in peridotitic liquids, but their solubilities can be estimated by extrapolating the relationship of ionic porosity and solubility for helium, assuming that hydrogen solubility continues to track with helium towards more depolymerized compositions on a molar basis.To do this, we apply the algorithm developed for helium solubility in silicate melt by (Iacono-Marziano, et al., 2010) to estimate H 2 solubility at 1 bar and 2273 K. Taking a melt composition equal to the BSE yields the relationship: Experiments demonstrate that temperature has a modest effect on the solubilities of noble gases (Jambon, et al., 1986) although data are limited to relatively cool conditions compared to those associated with magma oceans (<1873 K).It has been suggested that higher temperatures should favor the dissolution of H 2 into magmas (Chachan and Stevenson, 2018;Kite, et al., 2020), and this effect could be significant given the wide temperature range over which magma oceans likely exist.The model of (Iacono-Marziano, et al., 2010) predicts an order of magnitude increase in helium solubility upon heating from 1773 to 3272 K, also suggesting the importance of high temperatures in promoting H 2 dissolution into magmas.
We now apply the solubility laws for H to calculate its equilibrium distribution between a magma ocean and a primordial atmosphere.Solving the equilibrium distribution for any volatile element is accomplished by solving a system of equations that include the solubility laws, statements of mass balance, and the equilibrium constants for gas phase equilibria to define the relative stabilities of gas species as a function of magma ocean-atmosphere interface P-T-X conditions.Mass balance statements use total mass as a constraint, but the mass balance constraint could also be the pressure of a gas species in the atmosphere or the concentration of the relevant volatile element in the magma ocean.
To connect solubility laws to mass distribution, melt concentration is scaled to mass by multiplying concentration by the mass of the assumed magma ocean.Partial pressure (fugacity) is scaled to atmospheric mass as: where P x is the partial pressure of gas x in the atmosphere, A is the surface area of the planet, g is the gravitational constant, and r is the ratio of mass of the volatile element of interest in the gas molecule to the full mass of the gas molecule (e.g., 1 for H 2 and 2/18 for H 2 O).It is clear that the distribution of a volatile element depends on the size of the planet (A), its gravitational field (g), and the depth (or mass) of the magma ocean, so these also must be specified.We provide as an example the system of equations needed to solve the distribution of H.We apply a similar approach for other volatile elements in their respective sections.
where K is the equilibrium constant for Reaction 24, M H is the total mass of H in the atmosphere-magma ocean system, and M MO is the mass of the magma ocean.The equilibrium constant is defined by the P-T conditions selected and the associated Gibbs energy change for the reaction (calculated from data tabulated on the NIST WebBook and the ideal gas law).The factors of 100 are required in the leftmost terms on the righthand side to convert from wt. % to wt. fraction.The system is solved for the partial pressures of H 2 and H 2 O (PH 2 and PH 2 O).This assumes ideal gas behavior; high temperature-moderate pressure systems (hundreds of bars) only have small deviations from ideality.For example, the fugacity coefficient of steam is 0.91 at 1073 K and 500 bars (Helgeson and Kirkham, 1974) and that of CO 2 is 1.2 at 1500 K and 500 bars (Mel'nik, 1972), taking the standard state to be 1 bar.Higher temperatures serve only to force fugacity coefficients closer to unity.CO 2 -H 2 O mixtures do have excess energies of mixing, but again, this effect is generally small at the high temperatures and moderate pressures of volatile-dominated atmospheres overlying magma oceans.The solution may require iteration between the solved pressure of the atmosphere and that assumed to calculate K until there is agreement.We note that this approach can also be extended to predict core chemistry by adding a term that includes a partition coefficient (with implied P-T-X conditions for metal-silicate chemical exchange) along with the available masses of silicate and metal to be reacted to define the mass of the volatile element sequestered in the core for the specified P-T-X conditions.
Figure 8 plots the distributions of H species between atmospheres and magma oceans as a function of f O 2 (relative to the IW oxygen buffer at 1 bar) and the total abundance of H.The temperature dependance of H 2 and H 2 O stability is largely canceled by referencing f O 2 to the IW buffer.For simplicity, we assume a total atmospheric pressure of 100 bars and a temperature of 2773 K.An assumed constant total pressure is supported by the relatively constant and high pressure that C species provide given the abundance of C in the BSE, as will be demonstrated in the next subsection.
The crossover for equal partial pressures of H 2 and H 2 O occurs near ΔIW−1 at 100 bars and 2273 K (Figure 8A).Water (steam) is much more soluble than H 2 , and this shifts the crossover point for equal H 2 O and H 2 concentrations in magma oceans to extremely reducing conditions (<ΔIW−6), such that H 2 O may dominate H dissolution in magma oceans over a wide range of accretion scenarios despite the higher partial pressures of H 2 in atmospheres (Figure 8B).In all scenarios, the magma ocean contains more H by mass than the atmosphere (Figure 8C).H preferentially partitions into the magma ocean under oxidizing conditions, but as f O 2 drops and f H 2 rises, the distribution becomes nearly equal, indicating the increased volatility of H in reduced systems.The small masses of H in primordial atmospheres indicate that they have only a small capacity to limit the incorporation of H into cores and determine the H abundances of the bulk silicate reservoirs (i.e., materials later derived from the magma ocean-atmosphere system, that includes mantle, crust, atmosphere, and oceans).

Carbon
Under high f O 2 , CO 2 is expected to be the dominant C-bearing gas species (Reaction 15).Spectroscopic results have identified that C is present as carbonate groups in magmas when reacted with CO 2 gas (Mysen and Boettcher, 1975;Fine and Stolper, 1986), and this observation suggests that the following set of reactions control CO 2 solubility in magmas: From these reactions, C solubility should linearly scale with CO 2 fugacity (α = 1 in Eq. 13), at least up to moderate pressures, and indeed this has been experimentally documented (Dixon, et al., 1995).
Carbon solubility increases with increasing alkali and alkali Earth metal contents in magmas, which is consistent with the stability of carbonate mineral-like complexes in magmas (Lesne, et al., 2011;Duncan, et al., 2017).Higher temperatures favor the dissolution of CO 2 in the melt as carbonate (Konschak and Keppler, 2014).Relatively large amounts of CO 2 (molecular) can also dissolve into more evolved (higher SiO 2 ) melt compositions, following Reaction 25.Higher temperatures also favor reactants in Reaction 25 (Konschak and Keppler, 2014), and it is therefore possible that molecular CO 2 contributes to CO 2 solubility within magma oceans.
Sufficient experimental data have been collected to permit the development of models for CO 2 solubility in magmas (Dixon and Stolper, 1995;Papale, et al., 2006;Iacono-Marziano, et al., 2012).We again apply MAGMASAT to predict CO 2 solubility in peridotitic liquids over magma ocean P-T conditions (Gualda and Ghiorso, 2015) (Figure 9).Solubility scales with PCO 2, as suggested by Reactions 25 and 26 MAGMASAT predicts that higher temperatures are predicted to decrease solubility, but the effect is relatively minor, with solubility decreasing by a factor of 2 between 1,500 and 3000 K. Experimental work on mafic systems also supports a relatively small effect of temperature; although there is evidence that CO 2 solubility increases with increasing temperature, this is opposite to the predictions of MAGMASAT.Least squares fitting of MAGMASAT predictions of CO 2 solubility at 2273 K for 1, 10, 50, and 100 bars yields the following relationship for a BSE magma composition: Under more reducing conditions, f CO 2 will decrease while f CO increases for any given temperature and pressure.With decreasing f CO 2 , the carbonate and CO 2 contents of a magma will correspondingly decrease, while C species associated with the dissolution of CO will increase.
The exact nature of C species dissolved in magmas in equilibrium with CO remains relatively uncertain.Higher than expected C solubilities were observed in reduced, lunar analog melts rich in Fe, leading to the suggestion that Fe-carbonyl groups could be a major component of CO solubility (Wetzel, et al., 2013).Other spectroscopic analyses of melts, with or without Fe, reacted with variable f CO have not found additional evidence that Fe-carbonyl is a significant component of C solubility.Rather, there is building evidence that single CO groups are the dominant species (Armstrong, et al., 2015;Yoshioka, et al., 2019).These groups may be molecular CO or CO species complexed to other components in the melt.Without a more detailed understanding for how CO interacts with magma, extrapolation of CO solubility determinations to magma ocean conditions will introduce some uncertainty.Nonetheless, it is clear that C becomes

FIGURE 8
Distributions of H between atmosphere and magma ocean as functions of fO 2 (ΔIW).Low-H (H equal to a 1 bar atmosphere of H 2 or H 2 O) and high-H scenarios (H equal to a 10 bar atmosphere of H 2 or H 2 O) are considered.Masses of H in the low-H and high-H scenarios are broadly similar to those inferred for the bulk silicate reservoirs of terrestrial bodies within the solar system.All scenarios assume an atmosphere of 100 bars total pressure, a temperature of 2773 K, an Earth-sized planet, and a full mantle magma ocean.(A) Predictions of partial pressures of H 2 and H 2 O gases in equilibrium with a magma ocean.PH 2 O is nearly constant across the range of fO 2 considered, whereas PH 2 increases under more reducing conditions.(B) Predictions of H species concentrations for magma oceans in equilibrium with partial pressures of gases in (A).(C) Mass ratios of H in the atmosphere versus H dissolved in magma oceans.Oxidizing conditions force H 2 O stability and nearly all H to dissolve into magma oceans.Reducing conditions enable a more equal distribution of H between atmosphere and magma ocean.less soluble in magma as f CO 2 falls and f CO rises.Yoshioka et al. (2019) developed the following equation to predict the solubility of C (wt. %) in MORB magmas for a given f CO: The unit for f CO is bar.Note that the solubility of CO is expressed as reduced carbon, not CO, so statements of mass balance that incorporate Eq. 28 must take this into account.
The distribution of volatile carbon species between atmospheres and magma oceans are not predicted to be strongly dependent on the amount of carbon present in the combined system as the solubilities of all relevant species scale with their respective fugacities, with exponents close to unity.
Further decreases in f O 2 will progressively destabilize CO, possibly resulting in graphite/diamond/carbide saturation if overall C contents are sufficiently high.Decreases in f O 2 also potentially help to stabilize methane or other CH-species if a source of H is also present (Figure 6): and CH 4 can then dissolve into the magma following: Once dissolved in a magma, methane may undergo further reactions to produce other hydrocarbon species.Spectroscopic observations identify CH 4 , or possibly methyl groups, dissolved in reduced, H-bearing silicate melts, supporting the importance of Reaction 29 (Mysen and Yamashita, 2010).Methane presumably dissolves into the ionic porosity of magma, following the example of H 2 above, but how melt-reactive hydrocarbon species interact with the silicate network has not been well explored, and their importance to magma oceans is therefore not as well established as it is for other species containing C or H.
Figure 10 plots the distribution of C species between atmospheres and magma oceans as a function of f O 2 and the total abundance of C.Only solubility related to CO 2 and CO is considered.We again assume a temperature of 2273 K, but total pressure is now determined by the combined pressure of CO and CO 2 .Atmospheric pressure therefore decreases as CO is converted to CO 2 with increasing f O 2 .We explore a high-C scenario (mass equivalent of 500 bars of C species in the atmosphere-magma ocean system) and a low-C scenario (mass equivalent 50 bars of C species in the atmosphere-magma ocean system).For context, present-day Venus has an atmosphere that is nearly 100 bars of CO 2 .The crossover for partial pressures of CO and CO 2 is near IW, but higher pressures stabilize CO 2 gas (Eq.28) and shift the crossover to more reducing conditions (Figure 10A).The crossover for the melt concentration of C related to CO and CO 2 dissolution occurs below IW because CO 2 species are more soluble than those related to CO (Figure 10B).Magma oceans are the larger C reservoir when CO 2 dominates solubility, whereas the atmosphere is the larger reservoir when CO dominates solubility (Figure 10C).The ability of the atmosphere to dominate the C budget therefore appears restricted to conditions more reducing than IW.
The discussion above for carbon assumes that the magma ocean remains undersaturated with respect to solid forms of C. High f CO (and f CO 2 ) and reducing conditions promote graphite or diamond saturation following the general reaction: Once stabilized, a buffering phase such as diamond or graphite is a new volatile reservoir that must be accounted for.The physical stability of the buffering phase must also be considered.For example, graphite is expected to be buoyant in magma ocean liquids and could therefore rise to the surface (Keppler and Golabek, 2019), whereas volatile elements dissolved in alloys or stabilized in intermetallic compounds (carbide, nitrides, sulfides) may be dense and could therefore sink within a magma ocean to ultimately join the core.

Nitrogen
N 2 is the dominant N-bearing species in the gas phase above magma oceans (Boulliung, et al., 2020) independent of fO 2 (Dalou, et al., 2022;Sossi, et al., 2023).However, depending on fO 2 , the speciation of N in silicate melts simulating magma ocean conditions falls into two domains.Under oxidized conditions, N dissolves as molecular N 2 with a valence state of 0 and its solubility is primarily dependent on pressure following Henry's Law (Libourel, et al., 2003;Dalou, et al., 2019;Grewal, et al., 2020;Bernadou, et al., 2021).At reduced conditions, N dissolves as ions with a −3 valence state (Libourel et al., 2003;Dalou et al., 2019;Mosenfelder et al., 2019;Boulliung et al., 2020;Grewal et al., 2020;Bernadou et al., 2021).Depending on fH 2 , these ions could either be anhydrous nitrides (N 3− ) or hydrous amines (NH 2− ) and ammonia (NH 3 ) molecules (Mosenfelder, et al., 2019;Grewal, et al., 2020).N 2 and NH 3 molecules physically dissolve into the ionic porosity of silicate melts whereas ionic nitrides and amides chemically dissolve into the silicate melt network by displacing O 2− from the silicate melt network (Libourel et al., 2003).Whereas N 2 , NH 3 , and NH 2− species have been observed by Raman and FTIR (fourier transform infrared) spectroscopy in quenched silicate glasses (Dalou, et al., 2019;Mosenfelder, et al., 2019;Grewal, et al., 2020) (Figure 6), the presence of anhydrous N 3− has only been observed by Raman by Dalou et al. (2022) at very low fO 2 .N solubility as N 2 can be represented as: Even though hydrous N-H species are widely observed in experimental silicate melts simulating magma ocean conditions, the dissolution of N solely as N 3− can adequately fit the observed N abundances in reduced silicate melts in both hydrous and anhydrous conditions, and can be represented as: The solubility of N as molecular N 2 under oxidized conditions (Eq.31) scales with fN 2 and is independent of fO 2 , whereas the solubility of N as N 3− under reduced conditions (Eq.32), in addition to scaling with (fN 2 ) with an exponent between 0.5 and 0.75.Although the exact fO 2 for the transition of N speciation depends upon several thermodynamic parameters like T, P, fluid composition, and melt composition, available evidence suggests that this transition likely occurs near the IW redox buffer for silicate melts applicable to shallow magma ocean conditions (Boulliung et al., 2020;Bernadou et al., 2021).
The solubility of N in analog magma ocean silicate melts has been extensively studied by experimental studies over the last two decades for a wide range of fO 2 , P, and T conditions.Libourel et al. (2003) calibrated the two-species model (dissolution of N as N 2 and N 3− in oxidized and reduced conditions, respectively) with experimental data at 1 bar to determine N solubility in basaltic silicate melts as a function of P N and fO 2 : The first term on the righthand side accounts for N 2 solubility, and the second term accounts for N 3− solubility.Bernadou, et al. (2021) showed that the formalism for N solubility in the silicate melt represented by Eq. 33 remains almost unaltered from 1 to 3,000 bars in the C-H-O-N system.
Figure 11 plots the distribution of N as a function of fO 2 between atmospheres and magma oceans.We explore a high-N scenario (mass equivalent of 10 bars of N species in the atmosphere-magma ocean system) and a low-N scenario (mass equivalent of 1 bar of N species in the atmosphere-magma ocean system).Atmospheres start to lose significant N while magma oceans start to gain significant N below IW (Figure 11A and 9 (center)).Figure 11B shows that there are essentially two regimes for N: i) above IW, nearly all N resides in atmospheres, as N 2 is relatively insoluble in melts; and ii) below IW −3, nearly all N resides in magma oceans, demonstrating the highly soluble nature of nitride complexes.Between these two regimes is a relatively narrow transition zone, although the fO 2 values within this transitional zones overlap with those inferred for magma oceans from studies of core formation.Atmospheres are the largest N reservoir when N 2 controls solubility to ~IW −2, but magma oceans quickly dominate the N budget as fO 2 drops.

Sulfur
It has long been known (Fincham and Richardson, 1954) that at fO 2 below the quartz-fayalite-magnetite buffer (<IW +3), S dissolves as S 2− by displacing O 2− from the anion sublattice.Therefore, at oxygen fugacities relevant to magma oceans (generally <IW), S solubility in the silicate melts can be represented as: The equilibrium of this equation can be described as: ln K eq ln aS 2− + ln fO 0.5 2 − ln fS 0.5 2 − ln aO 2− (36) The concentration of O 2− anions in the silicate melt greatly exceeds those of other anions, including S 2− .O 2− concentration or activity is assumed to be constant and Eq.36 can be modified as (assuming the activity coefficient of S  Frontiers in Earth Science frontiersin.orgcalibrated with experimental data for FeO-poor magmas (Namur, et al., 2016), their utility for atmosphere-magma ocean equilibration models is limited due the lack of a fS 2 term.
The empirical equation has FeO and fO 2 terms whose effects are counteract each other in Fe metal-saturated systems.Because magma ocean silicate melts are in equilibrium with metal during core-mantle differentiation, their FeO contents are thought to be lower than those of terrestrial basaltic magmas.Therefore, below IW, the effect of the fO 2 term dominates such that the S solubility in the silicate melt increases at increasingly reduced conditions.
Sulfur solubility is typically expressed as a function of fS 2 (Eq.37), but S 2 is not necessarily the dominant S-bearing gas species over the P-T-X conditions of magma ocean-atmosphere exchange.Indeed, S is unique in its tendency to form gas molecules with other volatile elements.Hydrogen sulfide is produced by S 2 when exposed to elevated fH 2 , and COS is produced by S 2 when exposed to elevated fCO.This dynamic makes predicting the atmospheremagma ocean distribution of S correspondingly more complicated than for H, C, or N.
Figure 12 plots the stability of S-bearing species for an atmosphere with 100 bars of total pressure.The gas phase is 10 mol% H, 0.03 mol% S, and the balance is C species.Temperature is fixed at 2273 K, fO 2 ranges between ΔIW −4 and ΔIW +2.Hydrogen sulfide is the most abundant S-bearing gas, followed by COS, under reducing conditions, whereas SO 2 is most abundant under oxidizing conditions.S 2 remains a minority species over the entire range of fO 2 .
Figure 12B plots the mass of S in magma oceans for our nominal scenario (10 mol% H, 0.03 mol% S, and the balance is C species), a high-S scenario (10 mol% H, 0.3 mol% S, and the balance C is species), and a high-H scenario (30 mol% H, 0.03 mol% S, and the balance C species).The mass of S in the magma ocean is calculated using Eq.37.All atmospheres are 100 bars total, and equilibrium is calculated for 2273 K.The horizontal solid lines bracket estimated S masses for the BSE in Figure 12B.Sulfur is highly soluble in magmas, and this manifests as the much larger masses of S dissolved in magma oceans compared to the equilibrium masses of S present in atmospheres (Figure 12A).
Eq. 37 assumes that HS − species do not contribute to S solubility, as the compositional dependence of S solubility does not include a term related to fH 2 .With this assumption, the stabilization of H 2 S gas with increasing fH 2 under reducing conditions acts to make S more volatile.An example of this is provided in Figure 12 as H scenario.There is debate as to whether HS − species are significant for S solubility (Baker and Moretti, 2011), so Eq.37 may underestimate S solubility in H 2 -rich systems.Nevertheless, it is clear that S is relatively soluble in melts and that primordial atmospheres contain only a small amount of S.

Beyond an equilibrium model
The analysis above reviewed solubility laws and mass balance statements to calculate the equilibrium model for the distribution of H, C, N, and S between magma oceans and atmospheres.But this framework is not complete and much work is needed to make new experimental measurements of solubility under the P-T-X conditions of magma oceans, including high temperatures, high pressures, ultramafic melt compositions, and highly reducing conditions.Nevertheless, applying an equilibrium model with a single P-T-X condition for chemical exchange within a system as large and dynamic as a combined magma ocean and atmosphere must fall short because P-T-X conditions vary with time and heliocentric location as planets grow.
As magma oceans cool, the interface temperature with the atmosphere will change (Lichtenberg, et al., 2021a;Lichtenberg, et al., 2022).The ability of a magma ocean and atmosphere to remain in equilibrium relies on mass exchange across their

FIGURE 11
Distributions of N between atmosphere and magma ocean as functions of fO 2 (ΔIW).Low-N (N equal to a 1 bar atmosphere of N 2 ) and high-N scenarios (N equal to a 10 bar atmosphere of N 2 ) are considered.Masses of N in the low-N and high-N scenarios are broadly similar to those inferred for the bulk silicate reservoirs of terrestrial bodies within the solar system.All scenarios assume a temperature of 2773 K, an Earth-sized planet, and a full mantle magma ocean.(A) Predictions of partial pressures of N 2 gas in equilibrium with a magma ocean.(B) Predictions of N species concentrations for magma oceans in equilibrium with partial pressures of N 2 gas in (C).(C) Mass ratios of N in the atmosphere versus N dissolved in magma oceans.Under oxidizing conditions, N nearly all N is present in the atmosphere, whereas reducing conditions force nearly all N to dissolve into magma oceans.
interface, but also the nature and associated kinetics of this exchange remain areas of study (Salvador, et al., 2017).The diffusion of volatiles from the interface may be important here, in addition to the convective patterns of the underlying magma ocean and its associated timescales and length scales for mixing (Solomatov, 2007).With continued cooling, magma oceans crystallize, leading to a new era of evolution in which volatiles will partition between minerals and melts, depending on the elemental exchange between the magma and its overlying atmosphere, and hence the formation of a protocrust (Bower, et al., 2022).Some minerals may be major hosts for volatiles, but many likely exclude volatiles from their structures.Crystallization should therefore lead to net outgassing.Outgassing may be facilitated by bubble nucleation, bubble growth, and the diffusion of dissolved volatiles in the magma towards bubbles (Ikoma, et al., 2018).Again, kinetics may be important in determining the chemical distribution of volatiles within the magma ocean-atmosphere system.Crystallization also offers the possibility of trapping buoyant magmas in the cumulate piles of magma oceans.Trapped magma provides a physical, rather than a chemical, mechanism for retaining volatiles in magma oceans (Hier-Majumder and Hirschmann, 2017) and is highly relevant for inhomogeneous crystallization paths in Earth-or super-Earth-sized exoplanets (Moore and Cowan, 2020).The amount of trapped liquid is also an expression of kinetic processes and the competition of melt migration timescales against crystallization timescales.Finally, recent work suggests that large impacts may result in transient, highly energetic periods in which the distinction between atmosphere and surface may vanish (Caracas and Stewart, 2023).Upon cooling, these atmospheres would shift from being silicate-rich to rich in moderately volatile elements.Ultimately, all of these processes must be accounted for to understand how much nature deviates from the equilibrium framework.In light of associated degeneracy, information from exoplanetary systems may provide additional constraints on the redox diversity and crystallization pathways of rocky planets within and beyond the Solar System.

Summary and outlook
The cores of rocky planets could host significant amounts of volatile elements.However, a major finding from the recent metal-silicate partitioning literature is that high pressures and temperatures can have opposing effects on the partition coefficients, and on average can lower the siderophile tendencies of some elements including volatiles S and C.This could result in decreased incorporation of these elements into Fe cores during late-stage core formation in larger rocky planets.Though partial equilibration of large impactors during the giant impact phase of planetary growth could result in the transfer of unequilibrated volatile-rich materials to cores.More measurements at high P-T are needed to ascertain trends for N and H.Both N and S appear to show increased affinity for iron while C partitions less into the metal at more oxidizing conditions.These change in siderophile behaviors of the volatile elements has significant implications for the volatile reservoirs as planets grow.When coupled with accretion models, the metal-silicate partitioning measurements support heterogeneous accretion scenarios whereby volatiles were delivered towards the later stages of core formation.This is in contrast with the increasing evidence for the volatile-rich nature of inner Solar System planetesimals which could have strongly influenced the volatile budgets of terrestrial planets.Though most of the studies reviewed here were applied to the formation of the Earth and other rocky Solar System bodies, it has been inferred that core segregation likely occur during the formation of rocky planets in extrasolar planetary systems and in planets larger than the Earth such as super-earths and sub-neptunes (E.g., Otegi, et al., 2020).Calculations have indicated efficient mixing between

FIGURE 12
Distributions of sulfur between atmosphere and magma ocean as functions of fO 2 (ΔIW).(A) An example calculation of S gas species partial pressures for a 100-bar atmosphere with a nominal composition of 10 mol% H, 0.03 mol% S, and the balance is C species at 2273 K. (B) Predicted concentrations of S (dissolved in magma oceans for a nominal atmosphere, a high-S atmosphere (10 mol% H, 0.3 mol% S, and the balance is C), and a high-H atmosphere (30 mol% H, 0.03 mol% S, and the balance is C species).The horizontal lines bracket the estimated range of S concentrations in the bulk silicate Earth (McDonough and Sun, 1995) (C) Mass ratios of S in the atmosphere versus S dissolved in magma oceans.All scenarios predict that most S is present within magma oceans.
iron and rocky materials (Wahl and Militzer, 2015) and, similarly, volatiles and rocky materials (Dorn and Lichtenberg, 2021;Kovačević, et al., 2022;Vazan, et al., 2022) would occur in such planets, and it would be beneficial to obtain measurements at higher P-T and a wider range of compositions than the current dataset to better understand metal-rock and volatile-rock differentiation mechanisms in super-Earths.
In the magma ocean stage of growth, significant amounts of volatiles can dissolve into the silicate parts of planets, in which solubilities and redox state determine the availability of volatiles to outgassing.The recent solubility studies indicate that under oxidizing conditions, most H and S are sequestered in the magma ocean, most N is outgassed to the atmosphere, and C is nearly equally distributed between the atmosphere and the interior.Under reducing conditions, nearly all N dissolves in the magma ocean, the atmosphere becomes the dominant C reservoir, H becomes more equally distributed between the interior and the atmosphere, and S remains dominantly in the interior.Coupling of measurements with microphysical to large-scale models is a direction for future work that could provide further clarity on some of these multi-scale processes.Finally, the understanding gained from these studies can guide the interpretation of new measurements of exoplanetary atmospheres, which can in turn provide additional constraints to assess the bulk compositions, formation histories, and the abilities of rocky planets to develop conditions suitable for the origin and sustainability of life as we know it.
Recent first principles/ab initio molecular dynamics and density functional theory (DFT) simulations have been used to gain further insights on how volatile elements are distributed among cores, magmas, and atmospheres.For example (Solomatova and Caracas, 2019;Davis, et al., 2022), found C coordination in silicate melts to be largely pressure dependent, implying that terrestrial magma oceans could have contained oxidized carbon polymers.The latter study also found strong clustering between C and Fe, which increases as a function of pressure, suggesting C could be present in cores.Molecular dynamics studies of Fe-N compounds lead to the conclusion that Earth's core likely contains very low amounts of nitrogen (Bajgain, et al., 2019), though this does not rule out higher N contents in the cores of other planets.DFT studies have been used to argue for the preferential incorporation of H in the core (Li, et al., 2020;Yuan and Steinle-Neumann, 2020).If validated by meaurements, ab initio methods could provide a way to estimate partition coefficients and solubilities at a wider range of conditions than are achievable with current experiments.A Bayesian inference approach might also provide a way to extrapolate these quantities to conditions beyond the measurements (e.g., Gaffney, et al., 2022).

FIGURE 5
FIGURE 5Scanning electron microscope images of quench products from metal-silicate partitioning experiments on a large-volume press (A)(Dalou, et al.,  2017)  and laser-heated diamond anvil cell (B).Bright regions indicate material with higher atomic number such as Fe-Ni alloys.The quenched silicate regions of both samples contain small inclusions.In the case of diamond anvil cell experiments, inclusions can be as small as the nanometer scale.

FIGURE 6
FIGURE 6Schematic contours of H, C, and N speciation determined by Raman spectroscopy in basaltic glasses synthesized in piston-cylinder experiments at 1-2 GPa and 1,400 °C-1,600 °C, shown as a function of fO 2 (relative to the IW buffer) and H content (as equivalent H 2 O, wt.%).Raman and FTIR measurements are from Hirschmann (2012),Armstrong et al. (2015),Dalou et al. (2019Dalou et al. ( , 2022)), and Grewal et al.(2020).Blank areas represent conditions at which literature data are currently lacking.This figure demonstrates that H, C, and N on their own may be present as multiple species at a given fO 2 or water content, and that they combine to form different molecules in natural and synthetic glasses.fO 2 of Earth and Mars formation are indicated.Sulfur speciation is not shown because the effect of H on S speciation has not yet been determined at fO 2 relevant to magma ocean conditions.

Figure 7
Figure 7 plots predicted H 2 O solubilities from MAGMASAT for a peridotitic liquid (BSE of (McDonough and Sun, 1995) as a function of PH 2 O (assumed ideal) and temperature.Solubility scales close to PH 2 O 0.5 and is essentially independent of temperature, although recent experimental data suggest that water solubility decreases with temperature.Fitting all predictions together yields the following equation: H 2 O wt.% 0.0575•PH 2 O 0.5158(19)

FIGURE 7
FIGURE 7MAGMASAT predictions of H 2 O solubility in an ultramafic melt at different temperatures.The curve is fit to the 2273 K data and is given by Eq. 19.

FIGURE 9
FIGURE 9MAGMASAT predictions of CO 2 solubility in an ultramafic melt as a function of PCO 2 at different temperatures.The line is fit to the 2273 K data and is given by Eq. 27.

FIGURE 10
FIGURE 10Distributions of C between atmosphere and magma ocean as functions of fO 2 (ΔIW).Low-C (C equal to a 50 bar atmosphere of CO or CO 2 ) and high-C scenarios (C equal to a 500 bar atmosphere of CO or CO 2 ) are considered.Masses of C in the low-C and high-C scenarios are broadly similar to those inferred for the bulk silicate reservoirs of terrestrial bodies within the solar system.All scenarios assume a temperature of 2773 K, an Earth-sized planet, and a full mantle magma ocean.Atmospheric pressure varies with the total pressure of C-bearing gas molecules.(A) Predictions of partial pressures of CO and CO 2 gases in equilibrium with a magma ocean.PCO 2 dominates above IW, whereas PCO dominates below IW. (B) Predictions of C species concentrations for magma oceans in equilibrium with partial pressures of gases in (A).(C) Mass ratios of C in the atmosphere versus C dissolved in magma oceans.Under oxidizing conditions, C is nearly equally distributed, but reducing conditions make the atmosphere the dominant C reservoir.