© The Authors 2024

1 Introduction

Among the >5000 exoplanets discovered thus far, there is a subgroup of ~500 terrestrial-sized worlds (R ≤ 2 R_⊕, M ≤ 10 M_⊕) on extremely close-in orbits of the order of a day or less around their star (Zilinskas et al. 2022). Consequently, they are thought to be tidally locked, with a dayside that experiences permanent illumination, thus exhibiting temperatures that may be in excess of the 1 bar solidi of most silicates (≥1500 K, e.g. Hirschmann 2000). Owing to their extreme characteristics, they are often labelled ultra-short-period (USP) planets, hot-rocky exoplanets (HREs), or lava-ocean planets (LOPs).

Due to the high temperatures, the silicate surface may not only melt, but also partially vaporise (e.g. Léger et al. 2011; Wolf et al. 2023), forming a tenuous vapour shroud commonly called a mineral or silicate atmosphere. Chemical thermodynamic equilibrium between the (molten) surface and the atmosphere is expected, based on the relative convection and crystallisation timescales of fully molten magma oceans (Salvador & Samuel 2023) or their surficial counterparts (Kite et al. 2016); for highly irradiated HREs, these timescales are considered to be long. As a consequence, the atmosphere should carry the imprint of the physicochemical conditions (pressure–temperature-composition) defined at the interface and, therefore (and by extension), information on the mantle composition of the planet, making these worlds promising targets for the study of planetary surfaces, interiors, and evolution.

To this end, several models have been used to explore hypothetical atmospheric compositions and thermal emission of rocky exoplanets through numerical simulations. Fegley & Cameron (1987) introduced the MAGMA code to calculate the partial pressures of stable gas species in high temperature vapours in equilibrium with silicate liquid relevant to the evolution of Mercury. This model was later employed by Schaefer & Fegley (2004), Schaefer & Fegley (2009), and Miguel et al. (2011) for the study of LOPs. According to their findings, alkali metals (Na and K), oxygen (O₂ and O), Fe, and SiO dominate mineral atmospheres, with the specific mixture being a function of irradiation temperature and the composition of the silicate liquid. Ito et al. (2015) refined this model by extending the aforementioned ‘zero-dimensional’ vapour models to a 1D atmosphere, incorporating realistic melt thermodynamics based on the MELTS formalism (Ghiorso & Sack 1995), equilibrium atmospheric chemistry, and full radiative transfer. Their work predicts thermal inversions in mineral atmospheres where temperatures steadily increase with altitude, resulting in strong mid-infrared emission features, especially from SiO. Zilinskas et al. (2022) adopted a similar approach and applied it to known potential LOPs with surface temperatures >1500 K, considering a broader range of molecular species (MgO, SiO₂, etc.). Their study emphasises the significant influence of TiO and MgO interplay on the thermal emission spectrum for compositions akin to the bulk silicate Earth (BSE, Earth’s mantle composition, cf. McDonough & Sun 1995), oceanic and continental crust and komatiite, or Mercury’s mantle bulk composition; these authors utilised the original MAGMA code (Fegley & Cameron 1987; Schaefer & Fegley 2004). Their model was later extended to incorporate the major volatile species such as CO₂, H₂ and H₂O (Zilinskas et al. 2022; Piette et al. 2023).

However, exoplanet compositions may span a broader range than those examined by Zilinskas et al. (2022), which are based on terrestrial rocks. Particularly, some rocky exoplanets might have iron (FeO)-enriched mantles by H₂O-induced oxidation of metal (Elkins-Tanton & Seager 2008) or high oxygen abundances in the planetary source material (Wang et al. 2019). The size of the metallic core is also degenerate with the FeO-content of the mantle when only mass and radius are considered (Dorn et al. 2015). This highlights one key uncertainty: the fraction of oxygen condensing during planet formation (Putirka & Rarick 2019). Small amounts of available oxygen might lead to planets rich in metallic iron (and, thus, large cores), whereas overabundant oxygen could convert more metallic iron into FeO, which would reside in the mantle (Wang et al. 2018, 2022; Spaargaren et al. 2023). This exchange is governed by the chemical availability (chemical potential) of oxygen, assuming the Gibbs free energy of formation of pure O₂ is 0 at 1 bar and T; its standard state: $${\mu }_{{\mathrm{O}}_2}=RT\ln {\mathrm{f}\mathrm{O}}_2,$$(1)

where R is the gas constant, T the temperature, and fO₂ the oxygen fugacity, a measure of the effective partial pressure of O₂, $p{\mathrm{O}}_2={\phi }_{{\mathrm{O}}_2}\cdot f{\mathrm{O}}_2$. Here, φ is the fugacity coefficient. For an ideal gas, φ = 1, such that fugacity and partial pressure are equivalent; at a higher total pressure (P > 1000 bar), they start to diverge markedly. Taking the above results together, one of the key uncertainties in constraining terrestrial exoplanet composition are: (i) the bulk elemental abundances of the major rock forming elements (Fe, Mg, Si) and (ii) the fO₂, as approximated by the FeO/Fe ratio of the planet. The fO₂ at the surface is crucial in shaping the nature of potential atmospheres owing to both the stoichiometry of vaporisation reactions (Visscher & Fegley 2013; Sossi et al. 2019; Jäggi et al. 2021; Wolf et al. 2023) and homogeneous gas phase reactions (Sossi & Fegley 2018; Sossi et al. 2020).

The aforementioned studies either assume a fixed, arbitrary fO₂ (Piette et al. 2023) or rely on codes that try to predict the partial pressure of O₂ in the vapour from the stoichiometry of the vaporisation reactions themselves; these codes include the original MAGMA code (Fegley & Cameron 1987; Schaefer & Fegley 2004; Visscher & Fegley 2013) and the recently published LavAtmos code (van Buchem et al. 2023).

However, the fO₂ calculated according to the stoichiometric method is correct only when (a) all possible reactions that release O₂ in the vapour phase are considered and (b) the set of melt components, along with their thermodynamic properties, is complete and accurate (due to the dependence of fO₂ on melt activities, see Sect. 2.1). The MAGMA and LavAtmos codes are concerned only with Na, K, Ti, Fe, Ca, Al, Mg, and Si self-consistently; yet any other element, for example, chromium, sulfur, carbon or hydrogen, may also be present in natural silicate liquids. In this case, the fO₂ given by the stoichiometric approach would be incorrect. Under thermodynamic equilibrium, the chemical potential, ${\mu }_{{\mathrm{O}}_2}$, must be the same in the atmosphere as in the melt; hence, there is a singular fO₂ that defines the system (Eq. 1). The fO₂ of a magma ocean–atmosphere system is a function of its bulk composition at a given pressure and temperature (e.g. Sossi et al. 2020; Hirschmann 2022). However, the compositions of exoplanets are not known, such that the fO₂ cannot be predicted a priori. Therefore, we consider fO₂ as an independent variable in this work, which resolves both problems.

We aim to study the impact of (hypothetical) exoplanet mantle compositions as well as of oxygen fugacity on the formation, nature and detectability of atmospheres on lava ocean planets. Our new model is developed in a vein similar to that described by Ito et al. (2015), Zilinskas et al. (2022) and Piette et al. (2023), including self-consistent thermodynamics and its effect on atmospheric structure and spectra. This study is structured as follows. Sect. 2 details the methodology, which entails the description of the modified MAGMA code that permits the oxygen fugacity to be controlled as an independent variable, the derivation of the atmospheric speciation, and the atmospheric pressure–temperature structure and emission spectrum. This model is applied to various hypothetical compositions, which we also derive in this section (Sect. 2.6). Section 3 presents the results. Section 4 discusses the implications, and we present our conclusions and summary in Sect. 5.

2 Methods

2.1 Vaporisation of the lava ocean

The first step in modelling a mineral atmosphere is the derivation of the equilibrium pressure and vapour composition at the ocean-atmosphere interface. The initial temperature of the magma ocean is estimated to be that of the dayside temperature of a tidally locked black-body planet, but it is adjusted in subsequent iterations to match the atmospheric base temperature (see Sects. 2.4 and 2.5). To then simulate its vaporisation, we use a modified version of the MAGMA code (Fegley & Cameron 1987; Schaefer & Fegley 2004) that allows us to fix the oxygen fugacity (fO₂) as a free parameter (instead of it determining via stoichiometric vaporisation, as done hitherto). We used the equilibrium constants from Fegley & Cameron (1987), including the one for the reaction SiO₂(l) = Si(g) + 2O, in spite of it having been revised later (Schaefer & Fegley 2004). This is because the Schaefer & Fegley (2004) model underestimates partial pressures of Si-bearing gases in an attempt to improve those for Na and K (which we omit unless otherwise mentioned, due to reasons detailed in Sect. 2.6). We then add the reactions $$\mathrm{F}\mathrm{e}\mathrm{O}(1,\text{ silicate })=\mathrm{F}\mathrm{e}(1,\text{ metal })+\frac{1}{2}{\mathrm{O}}_2,$$(2a) $$2\mathrm{F}\mathrm{e}\mathrm{O}(1,\text{ silicate })+\frac{1}{2}{\mathrm{O}}_2={\mathrm{F}\mathrm{e}}_2{\mathrm{O}}_3(1,\text{ silicate }),$$(2b)

which describe the conversion of the FeO component dissolved in the silicate liquid (l, silicate) to account for the change in silicate melt composition due to effects of fO₂. This leads to (a) the saturation of metallic Fe at low fO₂ as a separate liquid phase (l, metal) and (b) the presence of Fe³⁺ in silicate liquids promoted under oxidising conditions (e.g. Kress & Carmichael 1991; Berry & O’Neill 2021). The restriction of fixed fO₂ marks the chemical network as an open system with respect to O₂, while the number of atoms of non-oxygen elements (Si, Mg, Fe, and so on) is conserved¹. The respective equilibrium constants for the reactions in Eq. (2) were derived from the available reactions in MAGMA, so as to remain self-consistent.

In MAGMA, a series of reactions between assumed melt components (“pseudocompounds”) is balanced to iteratively solve for the activity coefficients of the parent oxide species. The activity of a species in the melt is defined as: $$a_i=x_i{\gamma }_i$$(3)

and relates the mole fraction of component i, x_i in the melt to its activity coefficient γ_i, which is any real positive number. Unfortunately, the implementation limits the range of all γ_i between 0 and 1, whereas, in reality, some activity coefficient (particularly for FeO) might exceed unity (Wolf et al. 2023). We further fix the activity of Fe(l, metal), when stable, at unity, a_Fe = 1 because it constitutes its own phase due to its insolubility in silicate liquids. The fO₂ below which metallic iron can form is computed by exploiting Eq. (2a), re-arranged to solve for the oxygen fugacity: $$f{\mathrm{O}}_2=\frac{K\cdot a_{\mathrm{F}\mathrm{e}\mathrm{O}(1,\text{ silicate })}}{a_{\mathrm{F}\mathrm{e}(1,\text{ metal })}},$$(4)

where K is the equilibrium constant of reaction (2a). MAGMA computes γ_FeO in the silicate; since a_FeO is set by Eq. (4), we can invert Eq. (3) for x_FeO (in the silicate) to solve for the amount (i.e. mole fraction) of metallic Fe (in the metal) with the additional mass balance constraint; $$\mathrm{\Sigma }\mathrm{F}\mathrm{e}={\mathrm{F}\mathrm{e}}^0+\mathrm{F}\mathrm{e}\mathrm{O}+{\mathrm{F}\mathrm{e}}_2{\mathrm{O}}_3,$$(5)

where the ratio xFe₂O₃/xFeO is also fixed at a given P, T, and fO₂ by Eq. (2b). We note that this is incorporated into the iterative scheme of MAGMA in order to account for the decrease in xFeO(l, silicate), and, correspondingly, an increase in the mole fractions (x) of all other oxides. A detailed example of how fO₂ affects the thermochemistry of the system is shown in Appendix A.2. However, some caveats still remain: firstly, in natural systems, other elements can be incorporated into the metal phase, for instance, Si or O, which are not considered here (we discuss the implications in Sect. 4.2). Secondly, we implicitly assume that the metal phase remains well-mixed with the silicate phase and, hence, pFe is buffered at a given pressure and temperature according to the reaction Fe(l, metal) = Fe(g), where aFe(l, metal) = 1. Were the metal to separate, then the silicate melt would record the fO₂ set at the last pressure–temperature of equilibration with the metal. Any subsequent changes in fO₂ would be dictated by Eq. (2b).

At a constant temperature, the formation of a gas species M_μO_ν(g), where M refers to any metal species, and μ and ν are the stoichiometric coefficients, involves two types of reactions. The first is the heterogeneous vaporization of an oxide component M_mO_n(l) from the silicate liquid, occurring between two phases – melt and gas (cf. Sossi et al. 2019): $$M_m{O}_n(l)=mM(\mathrm{g})+\frac{n}{2}{\mathrm{O}}_2(\mathrm{g}).$$(6)

The second type is a homogeneous reaction, i.e. occurring within a single phase (gas): $$M_{\mu }{O}_{\nu }(\mathrm{g})=\mu M(\mathrm{g})+\frac{\nu }{2}{\mathrm{O}}_2(\mathrm{g}).$$(7)

Taken together, the fugacity f M_μO_ν of a gas species M_μO_ν is: $$f{M}_{\mu }{\mathrm{O}}_{\nu }=\frac{K_{\text{het }}^{\alpha }}{K_{\text{hom }}}\cdot a_{M_m{\mathrm{O}}_n(l)}^{\alpha }\cdot f{\mathrm{O}}_2^{\beta },$$(8)

where we define $\alpha :=\frac{\mu }{m}\text{ and }\beta :=\frac{1}{2}(\nu -\mu \frac{n}{m})$; K_het and K_hom are the respective equilibrium constants to Eqs. (6) and (7), adapted from MAGMA, and $a_{M_m{O}_n(l)}$ is the aforementioned activity. Then, fO₂ is kept as a free parameter. The coefficients for evaporation of all gas species considered in this study are given in Table 1.

Finally, we impose the ideal gas law to completely define the fugacity, fG, of any gas G, that is fG = pG, with p being the partial pressure. This is justified by the anticipated small value of the total vapour pressures (≤1000 bar, Ito et al. 2015; Zilinskas et al. 2022; Wolf et al. 2023; van Buchem et al. 2023), for which fugacity and pressure coincide.

The value of fO₂ is itself a function of temperature and pressure, hence, it is useful to normalise absolute fO₂ to a mineral buffer; here, we chose the iron–wüstite (IW) buffer from Hirschmann et al. (2008) since the main redox control occurs via FeO: $$\log f{\mathrm{O}}_2=\mathrm{I}\mathrm{W}(P,T)+\mathrm{\Delta }\mathrm{I}\mathrm{W}.$$(9)

We denote the (logarithmic) deviation from IW (which is Eq. (2a) but for a reaction of pure solid Fe and FeO) as ΔIW. Technically, Eq. (2a) is pressure-dependent, but the anticipated atmospheric pressures of ≤1000 bar have a negligible influence on the volume change (and, therefore, on the pressure sensitivity) of solid-solid buffers.

2.2 Atmospheric speciation

Far above the lava ocean, chemical, and thermal equilibrium with the underlying magma ocean is lost, and, as such, the temperature and pressure in the atmosphere differ from those at the interface. Here, we assume that all atmospheric layers have the same bulk elemental composition as the vapour in equilibrium with the melt, namely, the atmosphere is assumed to be well-mixed, even though this might not be the case in the presence of expected thermal inversions (cf. Ito et al. 2015). Further, the condition of thermochemical equilibrium in each atmospheric layer is invoked. Following these restrictions, we can compute the gas speciation with the FastChem code (Stock et al. 2018, 2022), which produces look-up tables for the thermodynamic properties of gas species, spanning a range of 10⁻⁸−l0³ bar and 1500–6000 K (sufficient to model mineral atmospheres). The reactions involved in FastChem are displayed in Table B.1. Although reaction kinetics are ignored, the high temperatures (≥2000 K) are likely sufficient to establish equilibrium. For the purposes of this study, we neglected the effect of photochemistry.

2.3 Opacities

2.3.1 Atoms

The atomic species allowed in the modelled atmospheres are Si, Mg, Fe, O, and Ti (in a few cases, also Na & K). Their opacities were obtained from The Data & Analysis Center for Exoplanets (DACE)², which also provided the tabulated opacities as a function of wavelength. We utilised the line lists from Kurucz (2017) for all atomic species. The spectra on DACE were generated with the opacity calculator HELIOS-K (Grimm & Heng 2015; Grimm et al. 2021), using a resolution of Δν = 0.01 cm⁻¹. The lines were treated as Voigt-profiles where the line wings are cut at a distance of 1000 000 Δν cm⁻¹ from the line center, effectively amounting to no line cut. This procedure almost certainly overestimates the effect of the line wing on the opacity, but the strong resonant lines of many metal gases do not merit the assumption of line cuts (e.g. Mg at 285 nm). Ideally, they would need to be modelled with other line shape theories, as done for Na in Zilinskas et al. (2022). We used atomic opacities only for a pressure of 10⁻⁸ bar, noting that their opacities on DACE are only available at 10⁻⁸ bar.

2.3.2 Molecules

We included the molecules SiO (Yurchenko et al. 2021), SiO₂ (Owens et al. 2020), MgO (Li et al. 2019), O₂ (Gordon et al. 2017), and TiO (McKemmish et al. 2019). Their opacities were obtained from the Exoclimes simulation platform³ which hosts opacities prior to their upload to DACE. The line lists for most molecules stem from ExoMol⁴ (Tennyson et al. 2016), except for O₂, which originates from HITRAN⁵ (Rothman et al. 2005, 2009; Gordon et al. 2017, 2022). Contrary to the atoms, the Voigt profiles have a line wing cut at 100 cm⁻¹ from the line core (smaller than for the atomic species, as they lack broad resonance lines), but are similarly sampled at a resolution of 0.01 cm⁻¹. Most molecular species here have line list data available over the range 1500–4500 K, except for SiO₂, for which no data beyond 3000 K exist.

Further, we use opacities sampled over a pressure range 10⁻⁸−10³ bar. Pressure induces an increase in the width of a spectral line due to collisions with the surrounding atmosphere, an effect known as pressure broadening which becomes more relevant for atmospheres of higher pressure. When approximating line shapes with Voigt profiles, the pressure effect enters through an increase in the half-width Γ_L of the line (cf. Grimm et al. 2021). ExoMol models the half-width as (Tennyson et al. 2016): $${\mathrm{\Gamma }}_L\propto {\mathrm{\Gamma }}_{ref}{\left(\frac{T_{ref}}{T}\right)}^n\left(\frac{P}{P_{ref}}\right),$$(10)

where Γ_ref and n are the broadening parameters and depend on the nature of the ambient gas. Different lines of the same species might have different coefficients. The broadening coefficients that are available are usually valid only for interactions with H₂, He, or air as background gases. However, for the metal and metal oxide species relevant to this study, no broadening coefficients (regardless of atmospheric background) are available. We therefore make do with the existing default parameters from ExoMol, Γ_ref = 0.07 cm⁻¹/bar and n = 0.5, which are applied to all lines of all molecules indiscriminately.

Here, O₂ is treated differently because its line-lists were obtained from the HITRAN database. In this case, the line broadening is assumed to stem from collisions between molecules of the same species, and the half-width takes the following form: $${\mathrm{\Gamma }}_L\propto {\mathrm{\Gamma }}_{ref}{\left(\frac{T_{ref}}{T}\right)}^n\left(\frac{{\mathrm{\Gamma }}_{self}{P}_{self}}{P_{ref}}\right).$$(11)

Here, the broadening parameter Γ_self denotes the self-broadening, the reference pressures and temperatures are P_ref and T_ref, and the temperature exponent is n.

An example for the opacities of major gas species in mineral atmospheres, sampled at 2500 K and 0.1 bar (broadening included) is shown in Fig. 1. We acknowledge that our treatment of line broadening is a simplification necessitated by a lack of data. However, we also consider a case with and without broadening and highlight the differences in the computed spectra (cf. Sect. 4.7 and Appendix F).

2.4 Radiative transfer

The atmospheric structure is computed self-consistently with the radiative transfer code HELIOS (Malik et al. 2017)⁶. This code treats the atmosphere on a log₁₀-equidistant grid confined between the pressure at the top-of-the-atmosphere (TOA), fixed at 10⁻⁸ bar, coinciding with the lower limit on the available opacities, see Sect. 2.3, and the bottom-of-the-atmosphere (BOA), which equals the vapour pressure atop the lava ocean. HELIOS uses the tabulated atmospheric speciation and the specified opacities to construct each layer’s transmissivity 𝒯. The mixing of the opacity of each species within the 0.01 cm⁻¹ wavelength-wide spectral bin is calculated via the random overlap method (Lacis & Oinas 1991) unless the maximum opacity of the species to be added is 1% of the total opacity; in that case, the correlated-k method is used (Malik et al. 2017). Effectively, this assumes that the transmissivity is a product of the individually opacities weighted by the mixing ratio of the respective gas species (Amundsen et al. 2017). Then, 𝒯 is used to solve the radiative transfer problem by iteratively adjusting the temperature profile to the radiation flux (cf. Malik et al. 2017), starting from an isothermal atmosphere as the initial guess.

Further boundary conditions are the stellar spectrum at TOA (which can either be a black-body or a (user-provided) real or simulated spectrum. Either spectrum is then scaled to match an arbitrary imposed irradiation temperature, T_irr (Sect. 2.5) and the internal temperature which dictates the heat flux, F, at the BOA. The internal temperature, T_int, is given by: $$T_{int}={\left(\frac{F_{out}-F_{in}}{\sigma }\right)}^{\frac{1}{4}},$$(12)

where σ is the Stefan–Boltzmann constant, and reflects the energy balance between that released by the planet to space (F_out) and that being shed from the interior (F_in). We leave this parameter at its default setting of 30 K. We neglect any convective adjustment to the temperature-pressure profile, since previous studies have found the atmosphere to be characterised by a thermal inversion and is thus strongly stratified (Ito et al. 2015; Zilinskas et al. 2022).

The converged solution yields the pressure–temperature profile of the atmosphere. Various useful byproducts are generated, such as the spectrum and atmospheric transmission function. However, the BOA-temperature changes during the run, hence the melt temperature would no longer correspond to the initial guess, namely, the substellar temperature, and is iterated accordingly. This requires subsequent calls to both the vaporisation routine (Sect. 2.1), FastChem and HELIOS, until the temperature of the melt has converged; at this point, the equilibrium partial pressures of the gases species stabilise and the atmosphere becomes stationary. In our case, we consider ΔT ≤ 35 K acceptable, given the wide range of temperatures considered. Upon convergence, the surface pressure may have changed by an order of magnitude from its initial guess. However, the irradiation temperature of the planet is maintained throughout a run.

Fig. 1 Opacities of major gas species sampled at $\mathrm{\Delta }\nu =0.01{\mathrm{c}\mathrm{m}}^{-1},\nu =\frac{c}{\lambda }$, at 2500 K and 0.1 bar (a typical value for the photosphere) as a function of wavelength. Note that the spectral lines of atomic O are omitted since are so weak that they would not be visible in the plot.

2.5 Temperature

The energy balance of the planet is set by the stellar radiation input, which predominately depends on the orbital separation d. HELIOS uses the stellar spectrum and d to infer the radiation flux at TOA. However, describing the problem via the planet’s irradiation temperature, T_irr, instead is more intuitive, and generalizes the result for planets that orbit larger or smaller stars of sufficiently similar spectrum. We find the orbital separation, d, which matches the desired T_irr by: $$d=\sqrt{\mathfrak{f}\cdot (1-A_B)}{R}_{\star }{\left(\frac{T_{\star }}{T_{\mathrm{i}\mathrm{r}\mathrm{r}}}\right)}^2,$$(13)

where A_B is the bond albedo (assumed to be 0, namely, the molten surface is a perfect black body, Essack et al. 2020) and T_* and R_* the effective temperature and radius of the host star (note: here we also assume the star to be a black body, which is an acceptable assumption to derive d). Also, $\mathfrak{f}$ is the dilution factor, a technical necessity to account for any redistribution of energy due to geometrical or dynamical processes. A tidally locked planet with no heat redistribution requires $\mathfrak{f}=\frac{2}{3}$, imposed by the spherical geometry of the planet; such an approximation is warranted for LOPs based on their expected thin atmospheres with little heat redistribution (Nguyen et al. 2020); however, it relies on the parallel beam approximation for the incident stellar light, which cannot be wholly descriptive for planets as close-in as LOPs (e.g. Léger et al. 2011; Carter 2018; Nguyen et al. 2020). Luckily, the choice of $\mathfrak{f}$ is not of relevance for our study of 1D-atmospheric profiles parameterised by T_irr, unless otherwise mentioned. We refer the interested reader to Hansen (2008) for a more in-depth discussion of $\mathfrak{f}$.

2.6 Exoplanet compositions

The compositions of hypothetical exoplanetary mantles were estimated from stellar compositions from the Hypatia catalogue (Hinkel & Unterborn 2018). To do so, their elemental abundances (in dex) were converted to elemental fractions. The bulk abundances of the major elements Si, Fe, Al, Ca, Mg, and O were then scaled to the observed Sun-to-Earth ratios (e.g. Wang et al. 2019) to mimic fractionation attending planetary accretion. Oxygen is then stoichiometrically distributed among the elements in the order Ca-Al-Mg-Si-Fe, based on the order of their respective metal-metal oxide equilibrium constants at 1 bar. The amount of O might be insufficient to oxidise all Fe and Si, in which case they will form the metallic core. However, this model is merely an approximation since the exact partitioning of O, Si and Fe between the mantle and core depends on the pressure and temperature of equilibration (Rubie et al. 2011); for exoplanets, this property is unknown.

This procedure results in 5296 planets, which is too many to model individually. Therefore, we used a Gaussian mixture model (GMM) to extract the representative exoplanetary mantle and core compositions (see Appendix C). To the oxides listed above, we added TiO₂ due to its strong impact on the atmospheric structure and spectrum (Zilinskas et al. 2022). Given that Ti is a refractory element with a similar condensation temperature to AI (Lodders & Fegley 1998), we adopted a constant TiO₂/Al₂O₃ ratio as is observed in BSE: ~0.0578 by molar abundance or ~0.0451 by mass (McDonough & Sun 1995). The constructed planetary compositions are listed in Table 2 and their position in the exoplanetary compositional space is indicated in Fig. 2.

We neglected Na₂O and K₂O in this analysis, as they are not refractory species; hence there is no orthodox way to estimate their abundance in exoplanet compositions. Thus, they were treated separately (see Sect. 4.1).

Fig. 2 Spread in SiO2/MgO and FeO/MgO ratios in the mantles of prospective exoplanet compositions according to the Hypatia database. See Sect. 2.6 for modelling details. Contours denote the distribution withing the planets that contain FeO in their mantles, and define the levels where the distribution contains 16%, 50%, and 84% of probability mass (from inner to outer contour). The white dots represent hypothetical exoplanet compositions selected by a GMM. The mantle ratios SiO2/MgO and FeO/MgO are shown for the solar system planets Earth (McDonough & Sun 1995), Mars (Khan et al. 2022) and Mercury (Nittler et al. 2018). Venus is assumed to be similar to Earth.

3 Results

We explored the possible parameter space of hot, rocky exoplanet (HRE) atmospheres on a grid with three axes representing the composition, temperature, and oxygen fugacity. For the composition, we used the planetary endmembers listed in Table 2. We assumed they represent fully molten magma oceans. Then, fO₂ covered the range ΔIW ∈ (−6, −4, −2, 0, 2, 4, 6), by some 12 orders of magnitude. This range is motivated by the trend observed among Solar System objects; the present-day upper mantle of the Earth has ΔIW+3.5 (Frost & McCammon 2008), Mars is roughly at ΔIW, whereas Earth and Mercury are inferred to have undergone core formation at ΔIW-2 and ΔIW-5, respectively (Righter et al. 2006). The irradiation temperatures were set to 2000, 2500, 3000, and 3500 K, which covers the temperature range of HREs (cf. Zilinskas et al. 2022). In total, this amounts to a grid of 5 × 7 × 4 = 140 HRE atmospheres. We further assumed a 1 M_⊕ and 1 R_⊕ planet, orbiting a l R_⊙ star. The stellar spectrum is that of the Sun (Gueymard 2004), but beyond 2.9 μm, the spectrum is assumed to be a black body of 5770 K temperature in order to remove erratic noise. The extrapolation was performed with the startool program, which is part of the HELIOS package (Malik et al. 2017). Most stars with potential lava ocean planets are “Sun-like” FGK stars (Zilinskas et al. 2022).

3.1 Vapour chemistry

3.1.1 Effect of temperature

Temperature dictates the total amount of outgassing; the pressure of any gas species increases exponentially with temperature. Because all evaporation reactions have similar enthalpies (curves are approximately parallel in Fig. 3, left column, see also Sossi & Fegley 2018), relative changes in the stabilities of gas species are minor, but not negligible. As shown by Schaefer & Fegley (2009); Miguel et al. (2011); Schaefer et al. (2012), we find that low-temperature atmospheres (~2000 K) are rich in Fe and some Mg and SiO. FeO may dominate over MgO (and under oxidising conditions over SiO as well) when the melt is comparatively cold and sufficiently iron enriched (a few wt%). At higher temperatures (≥2500 K), atmospheres become more SiO-rich at the expense of Fe, irrespective of composition (provided iron is present, archetypes CORL, ARS and TERRA), even when the FeO content of the mantle is strongly enhanced (~35 wt%, as in CORL).

3.1.2 Effect of composition

Unlike temperature, the effect of composition on the out-gassed atmospheres is, to a first-order, linear, where p_i ∝ a_i (cf. Eq. (8)). Hence, as the mole fractions of the major oxide components in exoplanetary mantle compositions typically vary by a factor of 2–3 (Table 2), the partial pressures of the gas species containing a given element vary by the same factor. This neglects any potential change in the activity coefficient (γ; Eq. (3)) with varying bulk melt composition, but this is a secondary effect in all cases (see Wolf et al. 2023, for more discussion). Following this logic, since most of the planetary archetypes have similar Mg/Si ratios, the equilibrium partial pressures (and hence the total abundances) of their respective gas species (Si, SiO, SiO₂, Mg, and MgO) are nearly identical. The only exceptions are compositions that are extremely deficient in SiO₂, such as XTREM. In this case, pSiO is second only to pMg in the vapour for the temperature range considered (≤3500 K). The abundance of iron, as FeO, on the other hand differs among our compositions, and thus the partial pressures of Fe-bearing species vary between high-FeO (CORL) and low-FeO archetypes (TERRA) by a factor of roughly 2–3, as expected. FeO-absent planetary archetypes (HERM & XTREM) have no Fe or FeO in their atmosphere.

Fig. 3 Vapour composition directly above the magma ocean as a function of temperature (left) and oxygen fugacity (right) for three representative mantle compositions in our sample: the composition for a coreless planet, CORL (top), the Earth-analogue TERRA (center) and the extremely silicon-depleted XTREM (bottom); see Table 2 for the respective compositions.

3.1.3 Effect of oxygen fugacity

Unlike the effect of temperature or composition, fO₂ can significantly change the nature and abundance of the outgassed atmosphere, by virtue of the fact that it varies by orders of magnitude. The most prominent effects are listed below:

• The partial pressures of SiO, Mg, and Fe scale as log $p_i\propto -\frac{1}{2}\log f{\mathrm{O}}_2$. Hence, their partial pressures are low for highly oxidised atmospheres (ΔIW > 3.5), but predominate under reducing to highly reducing conditions (ΔIW < 0). Their relative proportions, however, are nearly independent of fO₂, because their vaporisation reactions are defined by the same stoichiometry (cf. Eq. (8), Table 1);

• The partial pressures of MgO(g), SiO₂(g), and FeO(g) – with the latter having limitations that are discussed below – are independent of fO₂. Hence, they are less abundant than their oxygen-poor counterparts; Mg(g), Si(g), SiO(g), Fe(g) for highly reducing to oxidising conditions (ΔIW < 3.5, Fig. 3) but become the dominant form of the respective element in the gas phase at highly oxidising conditions (ΔIW > 3.5). Their partial pressures relative to other gas species (i.e. O₂ and O), however, remain relatively low at these conditions;

• The partial pressure of FeO(g) declines as FeO(l) is reduced to its metallic form in the liquid (i.e. Fe(l)) under highly reducing conditions (≤ΔIW-2) following reaction (2a). Since aFe = 1 when it precipitates (see Sect. 2.1), the partial pressure is fixed by pFe = K · aFe = K (with K the equilibrium constant of the vapourisation reaction Fe(l, metal) = Fe(g)) for all conditions more reducing than ~ΔIW ≤ −2 (dependent on temperature and composition). Thus, the partial pressure of Fe(g) becomes independent of fO₂ (Fig. 3), and highly reduced atmospheres are dominated by SiO and Mg, even if the melt composition is iron-rich. On the other hand, a minor decrease in pFe and pFeO is observed in highly oxidising atmospheres (≥ΔIW+3.5), where FeO in the melt is oxidised to Fe₂O₃ (Eq. (2b));

• When the melt is oxidising (ΔIW ≥ 2), the atmospheres are dominated by O and O₂, and log P ~ log pO₂ ∝ log fO₂. Other species become increasingly diluted.

The generalities in the trends of p vs. ΔIW detailed above are valid for all temperatures tested, but their absolute values (as well as the position of cross-overs between the partial pressures of species) depend on temperature (Fig. 3). The reason lies in the fact that reaction (8) depends on absolute fO₂, while the partial pressures displayed in Fig. 3 are plotted relative to the IW buffer, the fO₂ of which increases with temperature.

3.2 Surface pressure

Here, we report the pressure at the base of the atmosphere after it has achieved radiative equilibrium (Fig. 4). We find it to be a strong function of irradiation temperature and oxygen fugacity, but not composition. A 500 °C increase in temperature roughly increases the atmospheric pressure by a factor ~20, irrespective of composition and fO₂. Oxygen fugacity can result in a ~5-log unit change in total pressure for a given temperature and composition. We note a particular trend: $${\log }_{10}{P}_{\text{surface }}\propto \{\begin{array}{ll}-\frac{1}{2}{\log }_{10}{f}_{{\mathrm{O}}_2}& \mathrm{\Delta }\mathrm{I}\mathrm{W}\le 1,\\ {\log }_{10}{f}_{{\mathrm{O}}_2}& \mathrm{\Delta }\mathrm{I}\mathrm{W}\ge 1,\end{array}$$(14)

which is a direct result of the systematics described in Sect. 3.1.3. In particular, this comes from the fact that atmospheres below ΔIW-1 are mainly comprised of metal-bearing gas species (Fe, Mg, SiO), whose partial pressures are proportional to fO₂ ^−1/2 (Eq. (6)); whereas above ΔIW+1, O₂ predominates and hence the total pressure ~fO₂ (see also Fig. 3). However, oxidised atmospheres do not seem to follow this trend at an irradiation temperature of 2500 K (purple lines); their pressures ramp up suddenly for ΔIW≥ +2, with XTREM transitioning the earliest. This is a result of a greenhouse effect induced by MgO (see Sect. 3.3), which affects all oxidised atmospheres (ΔIW ≥ + 2) and, due to the hotter T_BOA of the atmosphere, enhances outgassing.

The minimum in atmospheric pressure is typically attained in the ΔIW~0–2 range; pressures are ~10^−4.5, 10⁻³, 10⁻² and 10⁻¹ bar for 2000, 2500, 3000, and 3500 K, respectively. Similar pressures at these temperatures are found by Miguel et al. (2011); Zilinskas et al. (2022); van Buchem et al. (2023). They do not explicitly report fO₂, so in order to compare, we estimate their fO₂ to be ΔIW+3.36, ΔIW+2.53 and ΔIW+1.93 at 2000 K, 2500 K and 3000 K, respectively, supporting the validity of our vaporisation model⁷.

The highest pressures are achieved by ultra-hot (3500 K) and highly oxidised (Δ IW>3.5) atmospheres, reaching up to ~3 GPa at Δ IW+6; more moderate temperatures of 2500 K lead to surface pressures of ~500 − 1000 bar in case of ΔIW+6. We note that because all gases are assumed ideal, only qualitative conclusions should be drawn from models at ΔIW+6. Highly reducing atmospheres (ΔIW ≤ −2) reach pressures of up to 145 bar for 3500 K irradiation, but are more tenuous at lower temperatures, around 0.1–10 bar. This indicates that mineral atmospheres are not necessarily as tenuous as previously thought (i.e. ~3 · 10⁻⁴, 1.4 · 10⁻², 2 · 10⁻¹ and 2 bar at 2000, 2500, 3000, and 3500 K, respectively, Visscher & Fegley 2013); however, extreme temperatures and oxygen fugacities are required to achieve pressures beyond 1 bar.

Fig. 4 Pressure above a lava ocean as a function of composition, irradiation temperature, and oxygen fugacity (expressed relative to the iron-wüstite buffer, IW), compared to atmospheric pressures of solar-system objects (grey, dashed). Source for pressures: Catling & Kasting (2017).

3.3 Atmospheric structure

The 140 simulated atmospheres can be broadly grouped by temperature and fO₂ regime based on their speciation, optical thickness, and pressure–temperature structure. The composition was found to exert only minor influences on most characteristics, so it is ignored in the classification. The resulting classes are labelled “cold and reducing,” “cold and oxidising,” “hot and reducing,” and “hot and oxidising,” based on their P–T-profiles (shown in Fig. 5) and their speciation & physical characteristics (e.g. Fig. 6). From each class, we selected a representative end-member in fugacity and temperature, which are labelled A–D in Fig. 5; the chosen melt composition is TERRA. The atmospheric characteristics of these representative endmembers are displayed in Fig. 6 and will be discussed in the following. We note, however, that this is only a broad classification, mostly based on our discrete sampling in oxygen fugacity, temperature, and composition space; in reality, all these cases transition into each other continuously.

In order to discuss the radiation physics of these atmospheres, we introduce the concept of optical depth or (as outlined here) its change along a single path of length, z, taken by a lightray as $$\mathrm{\Delta }\tau =\kappa \rho z,$$(15)

where κ is the opacity in units of area per mass, e.g. cm² g⁻¹ (see Sect. 2.3), ρ is the mass density, and z is the thickness of the absorbing layer. The transmissivity, t, is then defined by the following function (the Beer-Lambert-law): $$t=\frac{I_{\mathrm{o}\mathrm{u}\mathrm{t}}}{I_{\mathrm{i}\mathrm{n}}}=e^{-\mathrm{\Delta }\tau },$$(16)

where I_in/out are the incoming and outgoing intensity of the attenuated ray. These simple considerations apply to every layer of the atmosphere. We make use of τ in Fig. 6 to demonstrate how the radiation interacts with the atmospheric chemistry. Additionally, we will define the integrated transmissivity: $$t_k={\displaystyle \prod _{i=0}^k}{t}_i,$$(17)

where t_k is the integrated transmissivity of layer k and t_i is the transmissivity of atmospheric layer i, with i = 0 being the top of the atmosphere. This describes which fraction of a radiation beam coming from space has already been absorbed by the time it reaches layer k (or, by invoking time symmetry, how strongly the outgoing radiation from layer k has been attenuated). Setting t_k = 50% provides us with a tentative definition of the photosphere, displayed by the white ribbons in Fig. 6A1–D1.

We further define the (volumetric) mixing ratio of a species i in a gas as (Catling & Kasting 2017): $$v_i=\frac{N_i}{N},$$(18)

where N_i is the number of moles of species i in a gas, and N is the total mole number of the gas.

3.3.1 Cold and reducing

Outgassed atmospheres of ‘cold’ (~2000 K) lava planets of intermediate to reducing redox state (below ΔIW+2) are predominantly composed of SiO, Si, MgO, and Fe (if available) and monoatomic O (Fig. 6A3); the latter being as abundant as Si(g) throughout the atmosphere as a result of dissociation of SiO, the main carrier of oxygen. Free electrons from thermal ionisation become abundant towards the upper atmosphere, making it a potential conductor, albeit the atmosphere remains electrically neutral due to the remnant cations, which are predominately singly ionised (cf. Fig. 6A3 and Table B.1).

The relative stability of SiO throughout the atmospheric column has a strong effect on the pressure–temperature profile: since it is an effective absorber in the UV, the atmosphere becomes almost fully opaque at short wavelengths (≤0.3 μm; Fig. 6A1). The absorbed heat has to be re-emitted in the infrared, which is mostly transparent except for the strong bands of SiO at 4.5 and 9 μm. This requires the upper layers to heat up in order to achieve radiative balance by emission through the only available lines, and as a result, the atmospheric profile develops a thermal inversion (Malik et al. 2019; Gandhi & Madhusudhan 2019). In this class, the inversion extends throughout the entire atmospheric column (Fig. 6A2). Similar inversions were also found by Ito et al. (2015) and Zilinskas et al. (2022).

Atmospheric profiles of this class are largely independent of melt composition. This can be explained by SiO being the most important (and abundant) absorber or emitter in all reducing atmospheres (see Sect. 3.1), with Si, Mg, Fe, and O gases being consigned to play minor roles. It is worth mentioning that the emission of SiO in the infrared blocks important spectral bands (located at 9 μm) for the characterisation of silicate surfaces (Hu et al. 2012; Fortin et al. 2022). Otherwise, the atmosphere remains transparent throughout the infrared (see Fig. 6A1). The relative transparency of these atmospheres mean that the melt temperature closely approximates the irradiation temperature.

Fig. 5 Temperature–pressure profiles of all simulated atmospheres. Each subfigure contains all spectra (light grey), with the ones corresponding to the Tirr/ΔIW of the given row and column being highlighted. Color corresponds to composition (Table 2). The black dotted, vertical lines denote the irradiation temperature of the planet. Points A–D refer to a selected set of atmospheres that are shown in more detail in Fig. 6.

3.3.2 Hot and reducing

The atmospheres of hot (T_irr ≥ 2500 K) and reducing (ΔIW ≤ 0) planets are similar to their cooler counterparts (abundant Mg, Fe, Si and O), but they all have uniformly higher partial pressures, including SiO(g), due to the higher temperatures (Miguel et al. 2011; Wolf et al. 2023) (Fig. 6B3). However, high temperatures also cause SiO to dissociate more readily into weaker UV absorbers Si and O, which makes the upper layers of the atmosphere (p ≤ 10⁻⁵ bar) slightly more transparent than in cooler, equally reducing planets (i.e. the photosphere is situated at ~10⁻⁶ bar in the hot case vs. ~10⁻⁷ bar in the cold case, Fig. 6A1 vs. B1). Consequently, the capacity to absorb UV radiation is reduced, resulting in the development of(nearly) isothermal parts in the upper atmospheres; the isotherms in XTREM and HERM are more pronounced due to their even lower SiO content. Similar to the cold and reduced case, free electrons and cations are abundant in this layer.

Closer to the surface, at pressures higher than ~10⁻⁴ bar, SiO becomes the predominant species, again leading to a strong thermal inversion and hence to a temperature decrease towards the surface. Contrary to the cold & reducing case, this inversion is confined to intermediate atmospheric layers, and it cools the bottom layers of the atmospheres below the irradiation temperature of the planet (cf. Fig. 5).

Underneath the inversion, at pressures ≥0.1 bar, the higher pressures and cooler temperatures stabilise molecules such as MgO, FeO, and TiO. In particular, MgO is a strong, grey absorber in the infrared, and in combination with SiO it traps light in essentially all wavelength ranges, leading to optically thick lower layers (see also Zilinskas et al. 2022). Diffusive transport of radiation is initiated, as evidenced by the (approximately) isothermal P–T-profiles in the lower atmosphere (Fig. 5 upper right and Fig. 6B2), which becomes more pronounced in hotter, more reducing planets.

Fig. 6 Atmospheric structure of representative combinations of ΔIW and Tirr, selected from Figs. 5A–D. The left column displays the optical depth in every atmospheric layer (color gradient), as well as an indication of the photospere (white). The latter is defined where the atmosphere has absorbed 16% (upper dotted), 50% (solid) and 84% (lower dotted) of all incoming light. The center displays the atmospheric pressure–temperature profile, with the irradiation temperature Tirr highlighted as vertical dotted line. On the right, the (volume) mixing ratios of major chemical species are shown.

3.3.3 Cold and oxidising

Cold (~2000 K) and oxidising (≥ΔIW+4) lava ocean planets produce atmospheres mainly composed of monatomic O. Diatomic O₂ is also present, but dissociates in the layers above the lava ocean as temperatures increase (Fig. 6C3). Other species are scarce (with mixing ratios ≤10⁻⁴), though chief among them are Fe, SiO and Mg, whereas Si is less abundant compared to the reducing class of atmospheres (due to the homogeneous reaction $\mathrm{S}\mathrm{i}(\mathrm{g})+\frac{1}{2}{\mathrm{O}}_2(\mathrm{g})\to \mathrm{S}\mathrm{i}\mathrm{O}(\mathrm{g})$ being favoured). Additionally, FeO reaches relatively high mixing ratios in the lower atmospheres (with the potential to exceed SiO), but dissociates in the upper layers into Fe(g) and O(g). Free electrons and ions increase in abundance in the upper atmosphere, but the ionisation degree remains lower than in a reducing atmosphere.

Since SiO is no longer present in high mixing ratios (≤ 10⁻⁴), its role as a UV absorber is diminished, and only becomes relevant at relatively high pressures (above 10⁻⁶ bar) in the lower atmosphere. Furthermore, its number density remains low enough such that the infrared bands (weaker than the UV bands) only cause a minor attenuation of the in- and outgoing radiation and thus remain negligible (Fig. 6C1). Thermal inversions are still present, but weaken as the melt-vapour system moves towards ΔIW+6, a result of lower p_SiO. Due to the comparatively low temperature of the planet (~2000 K), the atmosphere never becomes thick enough for SiO, MgO, and TiO to act as efficient absorbers, rendering it relatively transparent throughout the entire infrared wavelength range, implying that the surface is observable from space. The transition from relatively transparent to fully opaque at the 9 μm feature – important for surface characterisation – occurs between ΔIW+2 and ΔIW+4.

Interestingly, this atmospheric class is sensitive to composition (see Fig. 5, lower-left). In particular, the upper atmospheres of the iron-free archetypes HERM and XTREM remain cooler than atmospheres with iron vapour (but are still in excess of the irradiation temperature). There is also a correlation between the iron content of an atmosphere and the temperature of the upper layers under highly oxidising conditions (>ΔIW+3.5). As SiO(g) has significantly lower mixing ratios than in the reducing cases, the relative opacity of Fe(g) becomes more important in shaping the upper atmosphere, thus imparting the difference in melt FeO content to the P–T-profile of the atmosphere. However, the majority of differences between HERM & XTREM versus CORL-ARS-TERRA lies within the lower SiO₂(l) abundance of the former group, which leads to lower mixing ratios of SiO(g) and therefore weaker inversions (Fig. 5).

3.3.4 Hot and oxidising

Hot and oxidising atmospheres extend to significantly higher pressures than their cooler counterparts (in most cases >1 bar; cf. Fig. 4), but are similar in chemistry (dominated by O₂ and O, with low pSiO), transparency (the photosphere is located at comparatively high pressures, 10⁻² bar) and the weakening thermal inversions with increasing fO₂. The atmospheric inversions almost disappear at ΔIW+6 and T_irr ≥ 3000 K, becoming near isothermal. This is due to the dilution of SiO by O and O₂, as well as suppressed pSiO in absolute terms.

Further down the atmospheric column, the same effect leaves MgO(g) as the most efficient absorber in the atmosphere; its opacities span over nearly all wavelengths except for the shortwave bands (≤400 nm) and are intrinsically higher than for SiO (see Fig. 1); hence, since MgO(g) is present in similar mixing ratios as is SiO at lower altitudes, it causes the atmosphere to become an efficient infrared emitter, leading to a flattened photosphere (Fig. 6D1). Below the photosphere, the atmosphere is optically thick and radiation must be transferred diffusively.

Another striking feature that is not seen in other classes is a distinct cold trap, expressed in all but the hottest and most oxidising cases. It occurs in the mid-atmosphere just below the photosphere at 0.001–1 bar (Fig. 6D2), and arises from a translucent optical window spanning the range 0.3–0.5 μm (Fig. 6D1). The window enables a substantial portion of stellar radiation to penetrate deeper atmospheric layers, heating them from below. Since the atmosphere is optically thick in IR at the pressure level of the window, it lacks the ability to re-emit to space, resulting in radiation entrapment that ultimately causes a greenhouse effect. The window is a consequence of a low SiO mixing ratio; in its wavelength range, there is a particular drop in the opacity continuum (cf. Fig. 1) unless either SiO or TiO appear in large mixing rations, both of which have low partial pressures in oxidising atmospheres.

Similar to their cooler counterparts, hot and oxidising atmospheres show some response to the iron content in their upper atmospheres, for the same reason. It is only in ultra-hot and oxidising cases (T_irr ≥ 3000, ΔIW+6) that this effect nearly vanishes, together with the thermal inversion due to the absolute dominance of oxygen.

3.4 Spectra

Due to their proximity to the star and expected tidal locking, atmospheres of LOPs are most readily studied during secondary occultations. This involves capturing a time series of spectra during the brief moments when the planet transits behind its host star. In the intervals just before and after the eclipse, the planet’s dayside is revealed, resulting in a combined flux of both the planet (F_p) and the star (F_s). The observed brightness difference is a function of the flux ratio between the planet and planet+star combination, which can be expressed as (e.g. Seager 2010): $$\frac{F_p}{F_s}=\frac{{\epsilon }_p{R}_p{}^2}{{\epsilon }_s{R}_s{}^2}$$(19)

where R_p and R_s are the radii of the planet and star, respectively, and ϵ_p and ϵ_s denote their respective spectral exitance (i.e. the energy flux emitted per wavelength and surface area, given here in erg s⁻² cm⁻³). This unique scenario facilitates the separation of the planet’s spectral information from that of its host star. The observation of LOPs during the secondary eclipse is favourable due to the fact that: (i) their atmospheres are expected to have higher molecular weights (M) than Hot Jupiters or Sub-Neptunes, which renders transit spectroscopy more challenging because of their low scale height, H, $$H=\frac{RT}{\mathrm{M}\mathrm{g}},$$(20)

where R is the gas constant and g the acceleration due to gravity; furthermore, (ii) the anticipated condensation of vapour toward the terminator, leaving little to no atmosphere along the planetary limb for study (e.g. Nguyen et al. 2020). Throughout this work, we interchangeably refer to the depth of the secondary eclipse as the secondary eclipse depth, occultation depth, or planet-to-star flux ratio. However, using the secondary occultation depth means that observables unrelated to the spectral exitance of the planet, such as the stellar emission and the planetary radius, are convoluted with the result (see Eq. (19)). However, the spectral existance of the modelled atmospheres is expected to vary negligibly with the characteristics of the planet, and the occultation depth of a larger planet with identical T_irr, fO₂, composition, and host-star type to the simulations shown here (1 M_⊕, 1 R_⊕) is expected to scale-up by $R_p^2$, where R_p is in R_⊕.

3.4.1 General features

Secondary occultation depths in Fig. 7 reveal that mineral atmospheres show features over a vast range of spectral domains. It should be borne in mind that, although the flux ratio increases as a function of wavelength, the absolute flux declines beyond the black body emission peak, making features beyond ~20 μm increasingly challenging to detect. Many spectral features, largely derived from molecular species (MgO, SiO, and SiO₂) are located in the mid-infrared (MIR, loosely defined here as 3–25 μm). Also, TiO has some strong bands that appear similar in shape and strength to the SiO features (cf. Fig. 1), however, owing to its low atmospheric mixing ratios they never become relevant. All these species should appear in atmospheres derived from any melt composition, as both Si and particularly Mg are lithophile. Moreover, this region is also accessible to modern or planned missions such as the JWST or Ariel.

Further, we find a plethora of lines in the ultraviolet (UV), visible (VIS), and near-infrared (NIR) up to ~1 μm. A considerable contribution comes from neutral atoms, which are particularly abundant at high altitudes in mineral atmospheres (see Fig. 6). For the most part, these lines are produced by Mg and Fe (300–520 nm), with some contribution from MgO (480 nm). The UV is dominated by emission from SiO, where it has its strongest feature (see Fig. 1); but also the strong doublet line of Mg (280 nm) is located there, as are some lines of Fe. In the visible and very near infrared (500 nm-1 μm), a mixture of MgO and TiO lines are present. The latter has stronger opacities per unit mass (Fig. 1), but MgO has higher mixing ratios and therefore contributes more than does TiO. Some weaker lines of Fe might also contribute in the NIR, but their opacities are lower than those of the aforementioned molecular species and therefore do not generate observable features at the resolution used to generate spectra by our model.

3.4.2 Effect of composition

Few substantial differences in the flux ratio of the MIR spectra can be ascribed to the composition of the underlying melt, regardless of fO₂ or temperature (Fig. 7). Compositions TERRA, ARS and CORL are nearly identical, even though the iron content varies significantly between them. HERM shows some slight reduction in the strength of the SiO features, yet follows the other compositions closely despite its considerably lower SiO₂/MgO ratio (~0.863 instead of ~1.16, the value for TERRA, CORL and ARS). Only the extremely reduced planetary archetype XTREM is distinct from its peers, as it has weaker SiO emission at both 4.5 and 9 μm owing to its low SiO₂ content. This observation can be reconciled with the predictions from the vaporisation model, Eq. (6), namely that pSiO should rely on the activity (and therefore concentration) of SiO₂ in the melt.

In the UVIS, the effect of composition under constant fO₂ (ΔIW = −2) and temperature (T_irr = 2500 K) on the emission spectrum are shown in Fig. 8; other T_irr would induce similar effects with stronger lines (compare Fig. 7), while the effect of fO₂ is discussed in Sect. 3.4.3. One of the major compositional controls is the emission from 330 to 460 nm by Fe(g), which differs markedly between the two groups of planetary mantles that contain FeO (CORL, ARS, TERRA) and those that do not (HERM, XTREM). We also noticed a consistent (albeit subtle) decrease in the Fe line intensity with increasing FeO(l). We note that CORL has highest pFe abundance and TERRA the lowest.

The other species reacting to melt composition is TiO via its broadband emission from 0.5–1 μm, which varies by a factor of ~3 from XTREM to TERRA. However, while the variation within the archetypes CORL, ARS, TERRA, and HERM corresponds to the abundance of TiO₂ in the melt (Table 2), XTREM shows an unusually low TiO(g) emission despite its higher abundance of TiO₂. This suppression of emission is probably not realistic⁸, but serves to illustrate the effect of changing TiO₂(l). The increased flux of XTREM in the MgO band beyond 1 μm (Fig. 8) is a result of increased heating of the lower atmosphere due to lower TiO(g).

3.4.3 Effect of oxygen fugacity

Oxygen fugacity has a considerable impact on the MIR flux ratio for the planets on our grid (Fig. 7). Most obviously, the intensities of the SiO features at 4.5 μm and 9 μm are more strongly expressed at reducing conditions (ΔIW ≤ 0), whereas they vanish under oxidising conditions (ΔIW ≥ 2) to the point where atmospheres at ΔIW ≥ 4 mimic a black body. Factors controlling the intensity of SiO emission features include:

• increased opacity through increased mixing ratios of SiO(g), where pSiO increases as a function of fO₂^−0.5 (see Sects. 2.1 and 3.1);

• increased opacity through pressure broadening: the total pressures are higher for reducing- (ΔIW ≤ 0) with respect to intermediate (0 ≤ ΔIW ≤ 2) atmospheres since the partial pressures of the monatomic gases of other major elements (Fe, Mg) also depend on fO₂ ^−0.5 (see Sect. 3.2);

• thermal inversions, which are more prevalent in reducing-and intermediate atmospheres, but weaker in oxidising atmospheres (see Sect. 3.3). It is this property that produces SiO emission features in the first place.

The corollary is that atmospheres produced in equilibrium with a more reduced silicate liquid produce stronger SiO(g) emission, while atmospheres derived from oxidised magmas (>ΔIW+2) have barely distinguishable SiO emission features relative to the black body baseline (Fig. 7).

The second most striking feature in the MIR is the difference in MgO background emission (yellow fields, Fig. 7) between oxidising (blue) and reducing (yellow) cases. Because pMgO is independent of fO₂, it is diluted by SiO(g), Mg(g) and Fe(g) (if present) in reducing atmospheres. Therefore, it can achieve high mixing ratios only at higher pressures, which pushes the photosphere in the MgO-bands to lower layers cooled by TiO(g) emission, thus decreasing emission in MgO bands (see Appendix D). This results in enhanced contrasts between the photosphere (dominated by the MgO background) and the SiO features at 4.5 and 9 μm, from essentially nil at ΔIW+6 and +4 to ~10 ppm (Earth-sized planet) at ΔIW-6. Larger and hotter planets should express stronger absolute contrasts; a 2 R_⊕, 2500 K planet might reach ~40 ppm contrast, a similarly sized 3000 K planet could achieve ~60 ppm (cf. Eq. (19)). The linearity in the spectral contrast over this wavelength range holds promise for the deduction of fO₂ in such atmospheres; a more detailed discussion can be found in Sect. 4.2.

Previously, it has been shown that, in addition to the SiO(g) emission feature, SiO₂(g) can be detected from an absorption feature at 7 μm in MIR-observations of LOPs with JWST (Zilinskas et al. 2022) and can carry information on the atmospheric fO₂ (Wolf et al. 2023). We find the SiO₂ 7 μm absorption feature to typically exist in oxidising atmospheres (ΔIW +2 to +4) of 2500–3000 K, which corresponds to the fO₂ expected in Zilinskas et al. (2022). In more reducing systems it is expressed as an emission feature that is roughly independent of fO₂ (as expected from Eq. (8)) and could thus serve as a baseline to compare the SiO(g) feature to, in order to derive fO₂. However, it is of invariably weaker intensity than the adjacent MgO and SiO features, making it comparatively difficult to quantify.

The UVIS/NIR spectrum is also affected by fO₂ (Fig. 9). However, even though the planet strongly radiates via these lines, the secondary eclipse depth is low (≤2 ppm for an Earth-sized planet) due to the brightness of the star in this wavelength range, though it increases with temperature (Fig. 7). As highlighted in Fig. 9, the iron lines between 330 to 460 nm respond to changes in fO₂. This could make them potential tools to trace fO₂, however, FeO contents can vary greatly owing to its extraction into planetary cores and iron line intensity in this range should be coupled with other fO₂-sensitive spectral feature(s). Moreover, the intensity of iron lines is not a linear function of logfO₂, owing to the saturation in Fe metal at low fO₂ (see Sect. 3.1), meaning maximum emission occurs at ΔIW-2 (shifting towards ΔIW+0 for higher temperatures, cf. Fig. 7). The partial pressure of Fe(g) decreases with fO₂−^0.5, such that oxidising atmospheres (≥ΔIW+2) have no remarkable iron lines and instead an emission ‘bump’ is observed, where hotter, lower layers radiate through the optical window formed by the paucity of SiO(g) (see Sect. 3.3).

Another notable difference is caused by the presence of TiO(g). Towards the visible wavelength range, we find that its emission (~0.45–1 μm) is suppressed for highly oxidising atmospheres (ΔIW ≥ 4), due to its low mixing ratios, a result of lower partial pressures (see Eq. (8)) and its dilution by O₂ in the atmosphere. Reducing conditions, on the other hand, produce strong TiO features (except for XTREM; see Sect. 3.4.2); however the variation in intensity as a function of fO₂ is not linear.

Lastly, the Mg and MgO peaks at 495 and 520 nm, respectively, grow in intensity as conditions become more reducing. Their ratio changes: in reducing atmospheres, the Mg(g) feature is more pronounced, whereas MgO(g) shows the higher peak under oxidising conditions. However, we note that the resolution of our opacities (0.01 cm⁻¹) is not sufficient to fully capture the shape of such narrow, highly peaked lines.

Fig. 7 Spectra of all 140 simulations of this study. Each column contains all spectra of equal irradiation temperature, and each row constant composition. In each box, all spectra of constant irradiation temperature are plotted (grey), but only the ones corresponding to the composition indicated in the row are coloured. The colour-coding is made according to their oxygen fugacity (yellow most reducing, brown intermediate, light blue most oxidising). Species responsible for important spectral features are indicated by the pointers, important wavebands by coloured patches.

Fig. 8 Effect of composition on the spectrum of a mineral atmosphere at constant fO2 relative to the IW buffer (ΔIW-2) and irradiation temperature (2500 K) for an Earth-sized planet. Shown are the stellar spectrum (top), the spectrum of planetary emission (center, coloured) and the secondary occultation depth Fp/F*, the ratio between the two (bottom, colours). Flux ratios for outgassed atmospheres may overlap, as evident for the CORL and ARS compositions. Only the ultra-SiO2-poor composition XTREM is slightly distinct by showing less emission in the TiO band. Both iron-free compositions (HERM and XTREM) do not express iron lines.

Fig. 9 Effect of fO2 on the spectrum of a mineral atmosphere at constant melt composition (TERRA) and irradiation temperature (2500 K) for an Earth-sized planet. Shown are the stellar spectrum (top), the spectrum of planetary emission (center, coloured) and the secondary occultation depth Fp/F*, the ratio between the two (bottom, colours). The emission from metallic species (Fe, Mg, MgO, TiO) in the visible distorts the overall shape of the spectral energy distribution (SED) compared to a pure black body. This effect is pronounced in reducing atmospheres, whereas oxidised ones follow a black body more closely except for a handful of strong lines.

4 Discussion

4.1 Alkali metals

We deliberately studied alkali-free systems, even though these elements have been theorised to occur in abundance in lava planet atmospheres (Fegley & Cameron 1987; Schaefer & Fegley 2004; Léger et al. 2011; Ito et al. 2015; Zilinskas et al. 2022). However, their volatility during planet-forming processes (Larimer 1967; Sossi et al. 2019; Gellissen et al. 2019) means that their abundances in the planet cannot be estimated reliably from those in the star, in contrast to more refractory elements such as Si, Mg, Fe, Al, Ca, and Ti. If Na₂O and K₂O were present in the melt, the resulting atmospheres would become rich in their respective metal-bearing gases: Na, K, NaO, and KO (Fegley & Cameron 1987; Schaefer & Fegley 2004; Zilinskas et al. 2022; Wolf et al. 2023; van Buchem et al. 2023). These gases would have high mixing ratios, particularly for low temperature (~ 1500 K) planets with tenuous atmospheres, and would absorb in the visible at their respective doublet lines around 589 and 769 nm.

In order to quantitatively determine the effect of alkalis on the compositions, P–T profiles and emission spectra of mineral atmospheres, the BSE composition from McDonough & Sun (1995), nearly identical to TERRA, Table 2, but containing 0.36 wt% Na₂O and 0.029 wt% K₂O, were used and modelled at ΔIW-4, IW and ΔIW+4 (Fig. 10). As a consequence, the planet exhibits strong but narrow emission lines, particularly of Na, in the visible & NIR (Fig. 10, left). This causes a slight cooling effect due to the increased emissivity of the upper atmosphere (see also Zilinskas et al. 2022); its magnitude depends on the temperature and fO₂ of the planet. It is most pronounced for cold planets (T_irr ~ 2000 K, all redox states) and the ΔIW = 0 case at 2500 K. Monatomic Na becomes the dominant gas in such atmospheres (Fig. A.1 because $p\mathrm{N}\mathrm{a}\propto f{\mathrm{O}}_2^{-0.25}$; however, it is overwhelmed by SiO at lower fO₂, and, at higher fO₂, by O₂). This also explains their higher surface pressures relative to the alkali-free cases (Fig. 10, right panel), consistent with the findings of Fegley & Cameron (1987); Miguel et al. (2011); Ito et al. (2015); Zilinskas et al. (2022). Hence, alkalis contribute to an increase in the total pressures of all atmospheres below T_irr ~ 2500 K. For planets with T_irr ≥ 2500 K, the overall impact on MIR spectra remains limited, and hence does not jeopardise the conclusions drawn earlier.

The volatility of sodium and potassium in the context of lava ocean evaporation means that a larger mass fraction of their total budget resides in the atmosphere relative to the other rock-forming elements considered here. Thus, continuous removal of atmospheric mass through condensation or escape to space would deplete Na and K in the residual magma more rapidly than for other elements, leaving behind the more refractory species (cf. Schaefer & Fegley 2009; Kite et al. 2016; Nguyen et al. 2020; Erkaev et al. 2023). Together with the weaker dependence of p(Na,K) on fO₂^{−0 25} compared to other species like SiO(g) or Fe(g), this indicates that the spectral features of alkali metals are (when taken on their own) poor indicators of fO₂.

Fig. 10 Influence of alkali metals (Na, K) on structure and spectrum of mineral atmospheres of an Earth-sized planet. On the left hand side, the difference in the planet-to-star flux ratio between an atmosphere with and without alkali metals is shown. The composition is TERRA (similar to bulk silicate Earth, BSE) with 0.36 wt% Na2O and 0.029 wt% K2O added (McDonough & Sun 1995).

4.2 Degeneracy between melt SiO₂/MgO ratio and fO₂

In Sect. 3.4.2, it was pointed out that both increasing aSiO₂(l) and decreasing fO₂ will result in more intense emission from SiO(g) features. To quantify the relative importance of these two variables, we ran some simple models of MgO-SiO₂ liquids at various compositions along the binary join, together with their respective atmospheric gas species, Mg, MgO, Si, SiO, SiO₂, O, and O₂. The planet was again assumed to be Earth-sized, irradiated at 2500 K, and fO₂ was sampled at ΔIW−4, IW, and ΔIW +4. The resulting MIR-spectra are displayed in Fig. 11, where they are compared to the native spectrum of TERRA (blue) at the given conditions. As expected, the strength of the SiO feature correlates with the SiO₂ content of the melt; the SiO(g) emission features at 9 μm at any T or fO₂ are near-absent (<5 ppm) in SiO₂-poor compositions (SiO₂/MgO < 1), and become more intense for atmospheres above more silica-rich melts. MgO shows no clear features in its spectral bands, contrary to what was seen in Sect. 3.4.3 and Fig. 7. Since it acts as a near-grey absorber, its “feature” in the MgO-band (as observed in the TERRA spectrum) is due to TiO(g), which acts to cool the lower atmosphere where MgO is spectrally active (see Fig. 11, right column; the P–T-profile of TERRA in blue compared to the TiO-free binary cases), leading to less thermal emission (Zilinskas et al. 2022).

To quantify the intensity of the SiO(g) feature, we introduce the contrast ratio Φ by comparing the intensity of its strongest feature at 9 μm to that of an adjacent MgO band at 6 μm (purple and grey fields in Fig. 11). We define the contrast ratio as: $$\mathrm{\Phi }=\frac{{\int }_{{\lambda }_{\mathrm{S}\mathrm{i}\mathrm{O},0}}^{{\lambda }_{\mathrm{S}\mathrm{i}\mathrm{O},1}}{F}_p\mathrm{d}\lambda }{{\int }_{{\lambda }_{\mathrm{M}\mathrm{g}\mathrm{O},0}}^{{\lambda }_{\mathrm{M}\mathrm{g}\mathrm{O},1}}{F}_p\mathrm{d}\lambda }\cdot \frac{{\int }_{{\lambda }_{\mathrm{M}\mathrm{g}\mathrm{O},0}}^{{\lambda }_{\mathrm{M}\mathrm{g}\mathrm{O},1}}{F}_{\star }\mathrm{d}\lambda }{{\int }_{{\lambda }_{\mathrm{S}\mathrm{i}\mathrm{O},0}}^{{\lambda }_{\mathrm{S}\mathrm{i}\mathrm{O},1}}{F}_{\star }\mathrm{d}\lambda },$$(21)

where λ_SiO,0 denotes the lower end of the SiO-band, and λ_SiO,1 the upper (likewise for MgO). Here, we use λ_MgO,0 = 5 μm and λ_MgO,1 = 6.7 μm, λ_SiO,0 = 7.5 μm, and λ_SiO,1 = 10 μm.

As shown in Fig. 12, the contrast ratio Φ at constant T_irr is sensitive to both fO₂ and composition. We recovered the earlier discussed suppression of the SiO(g) feature in oxidising atmospheres of a fixed composition, but also find that they can still be forced into emission by increasing the SiO(g) content, which occurs when the lava has a SiO₂/MgO ratio significantly higher than the bulk mantle (as indicated by growing Φ, see Fig. 12). In contrast, the reducing endmember (ΔIW−4) shows a turnover for high SiO₂/MgO ratios (≥3), which is the result of a growing inversion due to SiO(g), cooling the photosphere (cf. Fig. 11).

At fixed fO₂ (i.e. along the colored lines in Fig. 12), Φ varies by a maximum amount of ~0.20 over the entire range of SiO₂/MgO studied here (cases ΔIW−4 and ΔIW0, while ΔIW+4 shows a significantly smaller ΔΦ), amounting to ΔF_p/F_* ~ 5.4 ppm for a1 R_⊕ planet (21.53 ppm for 2 R_⊕). However, over the range of plausible bulk compositions – estimated by the 16th and 84th percentile in the hypothetical exoplanet SiO₂/MgO ratio – we observe a limited contrast ratio change of maximally ~0.07 (achieved by ΔIW−4), amounting to ΔF_p/F_* ≤ 1.9 ppm for a 1 R_⊕ planet or ≤ 5.88 ppm for 2 R_⊕. This means that its more challenging to separate between the majority of expected exoplanet bulk compositions, but it is possible to identify outliers (provided fO₂ is known), such as the SiO₂ depleted XTREM or (hypothetical) SiO₂-enriched surfaces derived from partial melting of the bulk mantle (Zilinskas et al. 2022); for instance, SiO₂/MgO~5.7–7 for mid-ocean ridge basalts (Klein 2003) and SiO₂/MgO~16.6 for continental crust (Wedepohl 1995).

A change in fO₂ at fixed SiO₂/MgO=0.98 (the 50% percentile) from ΔIW−4 to ΔIW+4 induces variation of ΔΦ ~ 0.21 (~5.65 ppm) for an Earth-sized planet or ΔF_p/F_* ~ 22.6 ppm for R_p = 2 R_⊕. For higher SiO₂/MgO ratios, we find equal or larger ΔΦ and ΔF_p/F_*, while for lower SiO₂/MgO ratios, the Δ values become smaller. Therefore, the effect of fO₂ on Φ exceeds that of composition by a factor of at least 2 for the plausible range of mantle compositions.

Taken together, this highlights fO₂ as a strong diagnostic tool, yet it is degenerate with composition in terms of their effect on Φ, which is only exacerbated once compositions with SiO₂/MgO outside the main bulge of the distribution (Fig. 12, top) are considered. In order to break the degeneracy, spectral features of other species, such as SiO₂(g) should be investigated. Equation (8) shows that pSiO₂ is independent of fO₂, but relates to the activity (and, hence, concentration) of SiO₂ in the melt. In Fig. 11, we see that its features at 7.5 and 34 μm are indeed sensitive to the SiO₂/MgO ratio; unfortunately, the variation with melt composition within the 7.5 μm feature is comparatively weak, while the (more reliable) 34 μm feature is beyond the reach of modern telescopes, including JWST.

The above considerations were made for a binary SiO₂-MgO liquid, while real melts also incorporate other oxides; therefore, we compare the contrast ratios of the planetary archetypes to the simple binary liquid (markers in Fig. 12). We find great agreement for the intermediate to oxidising atmospheres, but highly reducing cases (ΔIW−4) show larger contrast ratios for multi-component silicate liquids than in the binary. This is related to the aforementioned effect of TiO(g), which weakens MgO emission, therefore enhancing Φ (see Fig. 11, Sect. 3.4.3, and Appendix D for more detail). Additionally, TiO(g) seems to increase the degeneracy between melt SiO₂ and fO₂, as shown by the steep rise of Φ for the archetypes (markers in Fig. 12) compared to the pure MgO-SiO₂ liquid. However, a multicomponent melt releases additional redox sensitive species into the atmosphere (e.g. Fe(g) and TiO(g) that could help to constrain fO₂; see Sect. 3.4.3) and, by lifting the degeneracy, the SiO₂ content of the melt as well.

Fig. 11 Planet-to-star flux ratio spectra of atmospheres above a binary MgO-SiO2 melt, with varying mass ratios of MgO/SiO2 (shown in inset). We tested the conditions Tirr = 2500 K for an Earth-sized planet at a) Δ IW−4, b) IW and c) Δ IW+4.

Fig. 12 Φ as function of composition for a binary SiO2-MgO melt (lines, colored by fO2), assuming Tirr = 2500 K, 1 M⊕ and 1 R⊕, compared to Φ for a atmosphere-free, pure black body planet (black dashed line). The Φ of the multicomponent-melts of the archetypes in Table 2 are overlaid (markers). The distribution of SiO2/MgO of hypothetical exoplanet mantles is indicated as histogram (top), with the 16, 50 and 80 percentile shown as grey dotted lines.

4.3 Link between melt chemistry and redox state

At a given temperature and pressure, composition, and fO₂ are linked through the ratio of the activities of two oxidation states of the same element (e.g. aFeO/aFe₂O₃, aFe/aFeO or aSi/aSiO₂, and so on). In our model, the redox-sensitive species are Fe(l), FeO(l), and Fe₂O₃(l).

Fe₂O₃(l), like FeO(l), is a component in the silicate melt and becomes abundant (>10%) only above Δ IW+6 at 2500 K, and thus has minimal influence on the atmospheric P–T-profile or spectrum (see Appendix A). However, at reducing conditions, Eq. (2a) predicts that oxidised iron, FeO(l), should no longer be stable in large quantities; instead, it forms metallic iron Fe(l) + 0.5O₂(g). We treated Fe(l) as a distinct, pure phase where a_Fe = 1, which fixes the partial pressure of Fe(g) over the lava (Sects. 2.1 and 3.1). In reality, however, the metal phase may not be pure. As oxygen fugacity decreases, iron melt may also incorporate silicon (or other elements in low quantities) and, at high temperatures, oxygen (e.g. O’Neill et al. 1998). The incorporation of Si into metal is predicted, based on 1 bar thermodynamic data, to occur at highly reducing conditions (below ΔIW-6), but SiO₂ reduction is favoured with increasing T (O’Neill et al. 1998; Gessmann et al. 2001; Ricolleau et al. 2011). The reduction of SiO₂ into Si would lower the SiO₂ concentration of the silicate melt; once a fraction of SiO₂ is reduced to its metallic form: Si(l, metal) + O₂ ↔ SiO₂(l, silicate) (with activity aSi in the metal); thus, it must hold that $p\mathrm{S}\mathrm{i}\propto a\mathrm{S}\mathrm{i},p\mathrm{S}\mathrm{i}\mathrm{O}\propto a\mathrm{S}\mathrm{i}\cdot f{\mathrm{O}}_2^{1/2}$ and pSiO₂ ∝ aSi ·fO₂. If we were to treat the formation of Si-metal in a manner analogous to Fe-metal (i.e. it forms its own, separate metallic phase), then aSi= 1 and we would find constant vapour pressures of Si(g), while pSiO and pSiO₂ would drop for increasingly reducing conditions. This would leave the vapour increasingly dominated by Mg(g) as conditions become more reducing. In reality, Si and Fe-metal would form an alloy, and the activity of Si will be lower than 1 (Lacaze & Sundman 1991) and non-trivially interlinked with aSiO₂. Simultaneously, the activity of Fe in the metal phase would likely fall as more foreign species are dissolved into the metal. This would reduce the partial pressure of iron in the vapour, leaving atmospheres sourced from highly reducing melts increasingly iron-poor.

However, we acknowledge that our model cannot capture the formation of a multi-component liquid alloy as of yet. Further studies on the thermodynamics of silicate melts have to be performed, and the activities of more species should be calibrated.

Should metal formation occur in a magma ocean, as envisaged herein, then its high density relative to silicate liquid promotes its segregation and sinking to form a metallic core. Thus, the system is no longer chemically closed, and the iron-loving elements are extracted from the magma ocean. Hence, the compositions HERM and XTREM (Table 2) should be descriptive of planetary melts under highly reducing conditions (ΔIW ≤ −2).

4.4 Implications for geochemical characterisation of lava ocean planets

Thus far, the mantle composition and oxygen fugacity have been treated as independent variables, though they may be linked as established above. Since the planetary archetypes in Table 2 were constructed with a fixed total amount of oxygen, we can place a crude estimate on fO₂ set during core-formation of these hypothetical endmembers (see Appendix C for more details). The resulting fugacities are shown in Table 3 and were used to generate atmospheric structure and spectra of said composition-fO₂ combinations for T_irr = 2500 K and 1 R_⊕. Results are shown in Fig. 13.

We find the intensity of the SiO(g) feature to increase in the order XTREM-CORL-ARS-TERRA-HERM, which is mostly dictated by the descending fO₂. XTREM is the only exception; it expresses the weakest SiO(g) feature despite its highly reducing nature (ΔIW-5.82). This is a consequence of the previously established effects of SiO₂/MgO-ratio and fO₂: the extremely low SiO₂ in XTREM archetype leads to a low pSiO. HERM also has its SiO(g) feature slightly diminished (cf. Fig. 7), but due to its more reducing nature it expresses the strongest feature overall. TERRA only shows a slightly weaker feature; this echoes the degeneracy seen between low SiO₂ and oxidising conditions (Sect. 4.2 and Fig. 12). However, the question of whether compositions with SiO₂/MgO as low as XTREM are fully molten is uncertain due to the prevalence of MgO, which has a high melting temperature of 3098 K (Dubrovinsky & Saxena 1997). If it is excluded, we find that, despite the degeneracy, the strength of the SiO feature relates to fO₂ and rightfully identifies the most reducing compositions HERM and TERRA as such, while the more oxidised ARS and CORL show correspondingly weaker features.

However, the composition–fugacity combinations shown here may not be fully reflective of the actual exoplanet population. In the Solar System, Mercury has a higher SiO₂ abundance than predicted from solar abundances despite its low fO₂ of ~IW-5.4 (Namur et al. 2016), a combination that would result in enhanced SiO(g) features. Mercury-like compositions do, however, still occur in the stellar compositional range with the star-to-planet composition conversion employed (see Fig. 2). Additionally, the fugacities in Table 3 are expected for the core–mantle interface; the mantle, provided it contains FeO, may become more oxidised towards the surface due to evolution of the redox gradient throughout the mantle (Armstrong et al. 2019; Deng et al. 2020). This effect is exacerbated for coreless planets, which should be of oxidising nature (Elkins-Tanton & Seager 2008). FeO-free (reducing) mantles are less well studied, but the redox conditions might be set by species such as sulfur or carbon, as could have been the case on Mercury (Namur et al. 2016). In this case, the partial pressures of O₂ and mineral gases can still be predicted properly by using fO₂ as independent variable, even though the atmospheric composition and spectra would correspondingly change (cf. Jäggi et al. 2021).

Fig. 13 Secondary eclipse depth and atmospheric profiles produced by the planet/melt compositions from Table 2, evaluated with the respective fO2 (Table 3).

4.5 Determination of the oxygen fugacity of exoplanet targets with MIRI

As shown by Zilinskas et al. (2022) and Piette et al. (2023), the SiO₂(g) and SiO(g) features at 7 and 9 μm can be observed with JWSTs mid-infrared instrument (MIRI). To quantitatively examine the potential to constrain oxygen fugacity from such an observation, we constructed synthetic mineral atmospheres from the evaporation of a bulk silicate Earth-like mantle composition (TERRA), together with the orbital and physical parameters of K2-141 b, 55 Cancri e, TOI-431, and HD 231885 b; all potential lava planets (although recent JWST observations hint at a CO-CO₂ atmosphere on 55 Cnc e, Hu et al. 2024). Stellar spectra were taken from Zilinskas et al. (2022)⁹ and the physical and orbital parameters can be found in Table 4. We tested fO₂ between −6 and 6 log units relative to IW, assuming no heat redistribution, namely, with a dilution factor of $\mathfrak{f}=\frac{2}{3}$ (Hansen 2008). As close-in planets are most likely tidally locked, this is a reasonable assumption.

The resulting spectra are shown in Fig. 14A and analysed in Fig. 14B in terms of the flux contrast between the MgO and SiO bands, Φ (Eq. (21)). This provides merely a simplistic approach; actual observations should be analysed by atmospheric retrievals instead (e.g. Piette et al. 2023), but we expect them to mainly rely on the contrast between the MgO and SiO feature to constrain the redox state.

The simulated spectra in Fig. 14A show that the wavelength-dependent secondary eclipse depth varies according to the planetary and stellar characteristics, predominantly the relative size between planet and star, as well as their respective temperatures (cf. Eq. (19)). In particular, the model for planet K2-141 b shows large flux ratios in the MIR (up to ~250 ppm), despite its moderate size of 1.51 R_⊕ as it orbits a small and cooler K-type star (roughly 4000 ≤ T_* ≤ 5000, see Table 4). TOI-431 b, orbiting a similar star, is the smallest and coldest planet in this sample, yet it still produces an appreciable secondary occultation depth for identical reasons as K2-141 b. Finally, 55 Cnc e produces a strong signal due to its large size, while the similarly sized HD 213885 b (the hottest planet in this sample) orbits a star slightly larger and hotter than the Sun, thus exhibiting a weaker signal overall.

To distinguish between various redox states, however, the relative intensity of SiO-MgO emission features (the contrast ratio, Φ) is shown to be the most precise. The value of Φ shows a near-monotonic decrease as a function of fO₂ at constant temperature and composition, for all planets (Fig. 14B). Although K2-141 shows the deepest secondary occulation depth, its maximum Φ is significantly lower compared to the other planets considered (Fig. 14B1); TOI-431 b also exhibits a reduced contrast ratio (Fig. 14B2), yet it is not as extreme. This is attributable to the weaker thermal inversion expected in mineral atmospheres of planets orbiting K-stars (Zilinskas et al. 2022, and cf. Fig. E.1, this study). These weakened inversions are the result of less short-wave radiation being available to drive heating of the upper and intermediate atmosphere, where the SiO feature emits. The planets around the ‘Solar-type’ stars 55 Cnc and HD 213885 have stronger atmospheric thermal inversions and therefore produce a more marked contrast ratio between the SiO and the MgO features (Figs. 14B3 and B4).

To determine the extent to which such contrast ratios are observable, we constructed mock observations with JWST’s MIRI instrument in LRS mode using pandexo (Batalha et al. 2017), assuming the observations are stacked after one, three, and five occultations. The precision with which fO₂ can be inverted depends not only on the flux ratio, but also the absolute flux, for which the observation of the brightest target (55-Cnc-e) results in the smallest uncertainties per occultation (Fig. 14B3). In this (best-case scenario), five occultations would permit distinction between ΔIW cases to a precision of ~±l log unit, whereas for all other planets (B1, B2 and B4), only highly reducing (≤ΔIW−2) or highly oxidising (≥ΔIW+4) conditions could be ruled out.

The uncertainty on the feature strength strongly relates to the emission-spectroscopy-metric (ESM, Kempton et al. 2018): $$\mathrm{E}\mathrm{S}\mathrm{M}=4.29\cdot {10}^6\cdot \frac{B_{7.5}\left(T_{\text{day }}\right)}{B_{7.5}\left(T_{\star }\right)}\cdot {\left(\frac{R_p}{R_{\star }}\right)}^2\cdot {10}^{-m_K/5}$$(22)

where B(T) is the Planck function, T_* is the temperature of the host star, T_day is the dayside temperature of the planet (i.e. assuming $\mathfrak{f}=\frac{2}{3}$), R_p and R_* are the radii of planet and star, respectively, and m_K is the K-band magnitude. It captures the quality of a target for observation, in this case with the JWST. As seen in Table 4, 55 Cnc e has by far the highest ESM, and thus requires fewer occultations to reach a given signal-to-noise ratio. All other targets have smaller ESMs, and thus larger uncertainties in Fig. 14. The controlling factor is the distance of the target to the observer; the ESM scales proportional to ${10}^{-m_K}/5$ (Kempton et al. 2018), introducing a variance of a factor of ~6 between 55 Cnc e and the other planets. The other parameters in Eq. (22) are not as significant as they are similar between all planets, highlighting the importance of observing nearby targets.

Fig. 14 Emission spectra simulated for various existing lava planets, K2-141-b (green), 55-Cnc-e (blue), TOI-431-b (purple) and HD 213885-b (orange), assuming constant composition (TERRA, Table 2) and their corresponding stellar spectra (Zilinskas et al. 2022). Two diagnostic features for fO2 determination are highlighted: MgO at 6 μm and SiO at 9 μm. The contrast between the fluxes of the respective bins are shown on the right, for each planet, respectively, once as the contrast expected from a perfect observation (solid coloured lines) and the respective standard deviations (1σ) for various number of observed transits with MIRI (black, simulated with pandexo Batalha et al. 2017).

4.6 Beyond JWST

Within the atmospheres explored, we find that in addition to the intensity of the 9 μm SiO feature relative to the MgO background at ~6 μm, the Mg–MgO pair (as well as the lines of TiO and Fe in the visible) are promising with respect to the determination of fO₂. These lines correspond to a high spectral emittance of the planet (cf. Fig. 9) but a low flux ratio. However, the use of high-resolution cross-correlation spectroscopy (HRCCS) might prove fruitful in detecting element-oxide pairs in the atmosphere as well as species not included in this study that (potentially) have lines in the visible and are abundant in the melt: Cr–CrO–CrO₂, Mn–MnO, Fe–FeO, Ni–NiO, Al–AlO, or Ca–CaO, for example. Some of these species have already been observed with HRCCS to occur in hot Jupiters (Hoeijmakers et al. 2019). The relevant wavelength range of 0.3–1 micron is within the optical window of Earth’s atmosphere and is therefore less contaminated by telluric lines (e.g. Smette et al. 2015). To study the redox state of the outgassing melt, it is imperative to measure the abundances of the gases; while many studies in HRCCS aim solely at detection of species, they also allow for the retrieval of atmospheric abundances (Brogi & Line 2019).

4.7 Assumptions and limitations

To complete this study, a plethora of assumptions and simplifications had to be made. They are listed here in subjective order from most to least severe:

1. We treated the atmosphere as a 1D gas column with a single irradiation temperature. However, the hemispherical asymmetry implied by tidal locking could cause a continuous temperature decrease away from the substellar point, such that there is no single atmospheric temperature (Zieba et al. 2022). Furthermore, we cannot account for temperature variations caused by day-to-night winds (Castan & Menou 2011; Nguyen et al. 2020, 2022). However, the impact of geometry or atmospheric circulation on emission spectra is beyond the scope of this study;

2. We lack the opacity of one of the main atmospheric species, FeO(g), therefore neither the models shown here nor elsewhere (Ito et al. 2015; Zilinskas et al. 2022; Piette et al. 2023) are able to account for its potential effect on atmospheric structure and spectra. Additionally, it may reveal another channel through which to derive fO₂ by comparing the intensities of the spectral features generated by Fe(g) and FeO(g) or to link the spectra to mantle FeO content;

3. Pressure broadening coefficients for the species considered in this study are not well understood. When available, broadening coefficients are only valid for a background atmosphere composed of hydrogen, helium, or air (Tennyson et al. 2016). However, ignoring this effect might induce larger errors. We discuss this in more detail in Appendix F, where we find that the intensities of the SiO features are reduced by up to 5.36 ppm at the SiO-9μm feature (roughly 50% relative, compare with Fig. 7). Further, the atomic lines of all elemental species (Fe, Mg, Si, Ti, and O) are pressure-insensitive in our model. For the highly abundant Fe(g), this implies that its lines are likely saturated and will not respond linearly to changes in atmospheric iron abundance, hindering the retrievability of the planets iron content. Systematic studies of line broadening in different background atmospheres should be undertaken;

4. Melts at high oxygen fugacities (ΔIW ≥ 4) produce thick atmospheres (>1000 bar) that are oxygen-rich and, in extremely hot cases (3500 K irradiation temperature), they may even reach the GPa (10⁴ bar) range. At such pressures, our assumption of the ideal gas law is no longer correct, requiring an equation of state to properly calculate fugacities (e.g. Belonoshko & Saxena 1991);

5. More accurate models for the thermodynamics of liquid solutions are required. MAGMA might under-predict the activities of certain melt species and, hence, their partial pressure above the lava. On the other hand, VapoRock (Wolf et al. 2023) applies a more realistic thermodynamic model, but was not calibrated on compositions foreign to terrestrial rocks;

6. In our models (Sect. 2.2), we assume a well-mixed atmosphere (in terms of elemental abundances) despite the lack of thermal convection. However, diffusion could homogenise the atmosphere. Diffusion rates in a binary gas are well approximated by the Chapman–Enskog equation (Chapman & Cowling 1990), which shows that the diffusivity of a dilute species in an atmosphere is proportional to 1/P (P = total pressure). Diffusion is therefore faster at lower pressures and, hence, more efficient at higher altitudes. The spectroscopically active regions are typically below 0.1 bar, so diffusion should mix the relevant atmospheric layers relatively well. On the other hand, the lower atmosphere is frequently cooler than the upper atmosphere (Fig. 5), impeding mixing, either by advection or by diffusion, owing to higher P;

7. We assume that the melt can be considered as a well-mixed, effectively infinite reservoir. This implies that its chemical composition and oxygen fugacity does not change during vaporisation, a condition that might be violated by finite (or unmixed) melts. This may arise once the amount of vapour can no longer be considered infinitesimally small compared to the amount of liquid. Depending on T_irr and fO₂, the minimal required depth of the lava ocean can be crudely estimated as ranging from nil to ~1000 km (see Appendix G). The depths required will be correspondingly higher if the masses of other atmospheric gases (e.g. H, C) are present (Charnoz et al. 2023);

8. We assume no condensation. Contrary to more reducing atmospheres, oxidised ones above ΔIW+2 develop significant cold traps in lower layers (cf. Fig. 5, which could induce condensation and formation of cloud layers. Cold traps may lead to cloud formation if the P–T curve for a given condensation reaction intersects that of the atmosphere (Mahapatra et al. 2017; Herbort et al. 2020) and the atmospheric P–T profile might change in response. The net effect would be to trap easily condensable species in the lower atmosphere, leaving only the most volatile (e.g. K or Na) to diffuse into the upper atmosphere. Furthermore, hot lower layers with overlying cold strata are prone to vertical convection, which we do not consider in this study.

5 Summary and conclusions

We have devised a coupled framework for silicate melt vaporisation (a modified version of MAGMA that allows oxygen fugacity to be set as an independent variable), gas speciation (FastChem), and radiative transfer (HELIOS) to self-consistently compute the atmospheric structure and emission spectrum of lava ocean planets (LOPs). This model was used to test how (hypothetical) exoplanet temperatures, compositions, and oxygen fugacities influence atmospheric structure and resultant UV-, optical, and infrared spectra. The sensitivity of spectral features in the forward models were then compared to mock observations of four LOPs, in order to identify spectral features that, through inversion, may be used to place constraints on terrestrial exoplanet geochemistry. Our key findings are the following:

• The lava oceans redox state substantially shapes atmospheric chemistry. Under reducing conditions (below ~ΔIW), the atmosphere is dominated by metal-bearing gases such as SiO, Mg, and Fe, while above ΔIW+2, O and O₂-rich atmospheres predominate. The partial pressure of Fe(g) remains constant under strongly reducing conditions (≤ΔIW−2) due to the formation of metallic Fe. Gaseous TiO, a strong absorber, is mostly found in reducing atmospheres. Redoxneutral molecules such as MgO(g), SiO₂(g), and FeO(g) have constant partial pressures for all fO₂ for a given temperature and composition, and occur in comparatively low abundance;

• The total atmospheric pressure scales approximately as p ∝ fO₂ ^−0.5 for ΔIW ≤ 1, and p ∝ fO₂ for ΔIW ≥ 1. This implies that very reducing (ΔIW ≤ −4) and very oxidising (ΔIW ≥ 4) systems produce thicker atmospheres, sometimes in excess of 1 bar, while a well-defined pressure minimum occurs around IW to ΔIW+2, depending on temperature. Total atmospheric pressures generated at a given T and fO₂ depend only weakly on composition, unless the planet is cold (≤2500 K) and alkali metals are present;

• Atmospheric P–T profiles are nearly independent of melt composition, but are sensitive to fO₂. Reducing atmospheres show strong thermal inversions, resulting in hot upper atmospheres and increasing SiO emission thereby hindering cloud formation. Oxidised atmospheres have less pronounced inversions, to the point where (in some extreme cases) they have become nearly isothermal, e.g. ΔIW ≥ 6 and T = 3500 K). The latter can also form cold traps in their photospheres that are potentially sufficient to induce cloud formation;

• The primary spectral features in the mid-infrared include SiO peaks at 4.5 and 9 μm and the MgO “grey” background around 6 μm. Both fO₂ and the SiO₂ content of the melt influence their relative intensity, positively impacting SiO(g) emission. However, within plausible mantle SiO₂ ranges, fO₂-induced variation surpasses composition-induced changes. TiO(g) further decreases emission in the MgO band under reducing conditions, intensifying the contrast between SiO(g) and MgO(g) features, a trend amplified with increasing melt TiO₂ content;

• A subordinate factor influencing the intensity of SiO(g) mission is pressure-broadening. Hence, precise broadening theories, presently lacking, are essential for correlating atmospheric spectra with geochemical features of the underlying planet;

• The oxygen fugacity can be determined via observations made in the MIR. The diagnostic feature sensitive to fO₂ is the ratio between the line strength of SiO at 9μm and the emission in the MgO-band ~6μm, which could constrain fO₂ ~± 1 log unit with JWST observations, depending on the brightness of the target and the number of occultations. Systems that appear most promising for this endeavour are hot (≥2500 K) super-Earths orbiting nearby Sun-like stars (i.e. those with a high emission spectroscopy metric);

• Mineral atmospheres express strong lines in the visible due to evaporating neutral metal atoms, such as Fe, Mg, and Si, as well as some molecules such as TiO and MgO. They also carry the imprint of the melts fO₂ and could be potentially identified by ground based telescopes through high-resolution cross-correlation spectroscopy (HRCCS).

The MAGMA code can be obtained from Bruce Fegley Jr. upon reasonable request. The other codes used in this study are distributed under permissive software licences and can be obtained from the respective sources (see below).

Data availability

The model pipeline can be found on github, at https://github.com/ExPlanetology/phaethon. The atmospheric structure, chemistry, and spectra data can be found on zenodo at https://doi.org/10.5281/zenodo.13837367.

Acknowledgements

We are grateful to Laura Schaefer for a rigorous and constructive review that resulted in substantial improvements to the work. We appreciate the editorial handling and input of the editor, Emmanuel Lellouch. We thank Bruce Fegley Jr. for providing the MAGMA code, Lukas Carmichael for providing the exoplanet composition data and Matteo Brogi for our discussions on the possibility of HRCCS for small exoplanets. This work was supported by the Swiss National Science Foundation (SNSF) through an Eccellenza Professorship (203668) and the Swiss State Secretariat for Education, Research and Innovation (SERI) under contract No. MB22.00033, a SERI-funded ERC Starting grant “2ATMO” to P.A.S. Parts of this work have been carried out within the framework of the National Centre of Competence in Research (NCCR) PlanetS supported by the SNSF under grant 51NF40_205606. This publication makes use of The Data & Analysis Center for Exoplanets (DACE), which is a facility based at the University of Geneva (CH) dedicated to extrasolar planets data visualisation, exchange and analysis. DACE is a platform of the Swiss NCCR PlanetS, federating Swiss expertise in Exoplanet research. The DACE platform is available at https://dace.unige.ch. This work made use of the following codes: HELIOS (Malik et al. 2017, 2019), FastChem (Stock et al. 2018), Pandexo (Batalha et al. 2017), MAGMA (Fegley & Cameron 1987; Schaefer & Fegley 2004), VapoRock (Wolf et al. 2023), numpy (Harris et al. 2020), scipy (Virtanen et al. 2020), pandas (Wes McKinney 2010), matplotlib (Hunter 2007), seaborn (Waskom 2021) and scikit-learn (Pedregosa et al. 2011).

Appendix A Elaboration on outgassing

Here, we treat the outgassing from a bulk silicate Earth (BSE) composition (McDonough & Sun 1995), spanning the range of ΔIW-6 to ΔIW+8 to highlight (a) the effect of alkalis on the atmospheric pressure and composition and (b) the evolution of the concentration and activity of FeO and Fe₂O₃ as a function of redox state.

A.1 Alkali metals

As can be seen in the upper panel of Fig. A.1, Na(g) is prevalent under reducing and highly reducing conditions (ΔIW≤ 0), but reacts less sensitively to fO₂ than SiO(g) (cf. Eq. 8 and Table 1). Hence, it becomes the dominant species in atmospheres of reducing to intermediate redox states (ΔIW-2 to +2), significantly increasing the total atmospheric pressure compared to the alkali free case. This effect diminishes as temperature increases, as SiO(g) is more temperature sensitive (Wolf et al. 2023, and Sec. 3.1.1 this study). The increase in total pressure can be seen in the atmospheric profiles in Fig. 10. Under oxidising and highly oxidising conditions (ΔIW≥ + 2), Na(g) and NaO(g) are the more relevant metal gas compared to species of Si, Mg and Fe. As was seen in Fig. 10, Na(g) maintains some spectral features in the visible+NIR even under oxidising conditions, which is contrary to e.g. Fe(g) and Mg(g) (cf. Fig. 9). Potassium K(g) and KO(g) do not play significant roles in terms of abundances or spectra for BSE-like abundances (cf. Fig. 10).

A.2 Iron speciation and melt thermodynamics

As described in Sect. 2.1, the speciation of iron changes with fO₂; both within the silicate melt phase as Fe³⁺ and Fe²⁺, and as Fe⁰, which forms a new phase because it is insoluble in silicate melts. Consider a chemical system with fO₂ as an independent variable (Fig. A.1). In such a system, the total moles of Fe, Mg, Al, Ca, Na, K, Ti and Si are constant as a function of fO₂, while the number of moles of O is permitted to vary in accordance with fO₂ to satisfy Eqs. 2a and 2b. When fO₂ is sufficiently low that the system is saturated in Fe⁰ (i.e. xFe > 0, as given by Eq. 2a), the concentration (mole fraction) of FeO in the melt must decrease to satisfy the constraint that ∑Fe = constant (Eq. 5, see Fig. A.1, lower panel). In so doing, the mole fractions of other oxides in the silicate melt increase by a factor (1 − xFe_new)/(1 − xFe_old). This is observable in the small increases in xCaO and xAl₂O₃ in Fig. A.1 for a BSE-like composition. This effect would be exacerbated for compositions initially richer in FeO, like ARS or CORL. We emphasise that, although their abundances in the silicate melt phase increase, the moles of all the elements in the bulk system (except O) are conserved (Fig. A.1, lower center); FeO is simply replaced by Fe as fO₂ decreases, thereby maintaining the same number of moles of non-oxygen elements in the system.

In terms of activities, aFeO is low under reducing conditions due to the aforementioned formation of Fe-metal (see Sect. 2.1), but increases as conditions become more oxidising. The activity of FeO stabilises once Fe-metal destabilizes. This transition point is temperature and composition dependent, and for the BSE composition at 2500 K located at ~ΔIW-3. This value is too low; for these conditions, the transition should be closer to ΔIW-2 (Frost et al. 2008), likely a consequence of the crude thermodynamic model applied by MAGMA. aFeO shows a second “roll-over” at very high fO₂, at which more and more FeO is converted into Fe₂O₃ (indicated by the increase in the Fe₂O₃ activity, and decrease in a_FeO). This suppresses the outgassing rate of iron, as shown in Fig. A.1, where FeO(g) shows a roll-over at very oxidising states (≥ΔIW+6). This would yield very oxidising atmopsheres increasingly iron poor. However, it remains to be investigated if planets with such oxidising mantles can form, and if they could maintain such massive O₂-atmospheres in excess of 1 GPa.

Fig. A.1 Vapour pressures and thermodynamics of a silicate melt of BSE composition. Upper: Similar to Fig. 3, but for a melt of BSE composition (including Na2O). The dominance of Na(g) at intermediate redox states (ΔIW-3 to ΔIW+2) is evident, which is consistent with the findings of other studies (Schaefer & Fegley 2004; Miguel et al. 2011; Wolf et al. 2023; van Buchem et al. 2023). Upper center: Concentration of oxides in the silicate melt. FeO concentration deplets under reducing conditions, as Fe is formed. This is evident in Lower center: here, the bulk system composition is shown, highlighting the conversion of FeO to Fe. Lower: Activity of Fe and Fe2O3 in the melt, compared to the far more dominant SiO2 and MgO.

Appendix B FastChem reactions

The reactions allowed in FastChem are listed in Table B.1; this includes all species that can form in the mineral atmospheres of this study. The “basic” reactions cover the chemistry available in all atmospheres, while the “extended” selection is only used in specific use cases, i.e. Sec. 4.1. Species that also have an associated opacity and outgass significantly are highlighted in bold. Al, AlO, Ca and CaO for example would also be available opacity sources, but outgass in such low mixing ratios that we do not consider them. FeO however, despite its prevalence in mineral atmospheres, lacks line lists and therefore cannot be included.

Appendix C Determination of hypothetical exoplanet compositions

To derive the hypothetical exoplanet compositions we utilise (Table 2), we used a Gaussian mixture model (GMM), which can identify distinctive archetypes within the population, as depicted in Fig. 2. GMMs characterize a given population of d-dimensional datapoints by approximating it with a distribution composed of k d-dimensional Gaussians, with d being the number of parameters. For the problem at hand, it is best to transform the composition of mantle (wt% of SiO₂, MgO, CaO, Al₂O₃ and FeO) and core (wt% of Fe and Si) into a logratio w.r.t to MgO (mantle) or Fe (core) in order to remove the restriction of having to sum to 100%, while the core-mass-fraction CMF was also log-transformed in order to avoid the boundaries at 0 and 1. A challenge emerges from the oversimplified condensation model, which produces a trifurcation in the resultant planetary classifications: while ‘Earth-like’ planets (e.g. ARS, TERRA) exhibiting a well-defined parameter set, core-less planets are characterized by an absence of a core and its composition (e.g. CORL in our dataset), whereas planets formed under reducing conditions (HERM & XTREM) lack mantle FeO. This implies that d varies between the subclasses, and we therefore have to treat each one individually before collecting the final sample of exoplanet archetypes.

While the choice of k is arbitrary, one can ascertain an “optimal” number by varying it and evaluating the Bayesian Information Criterion (BIC) (Schwarz 1978): $$\mathrm{B}\mathrm{I}\mathrm{C}=-2\log L\widehat{(k)}+\log N\cdot d$$(C.1)

Here, N represents the number of samples, and $\hat{L}$ denotes the likelihood as determined by the respective algorithm in use (in this case, the GaussianMixtureModel from scikit-learn, Pedregosa et al. 2011). Typically, one would anticipate the BIC to be minimized for the number of Gaussian components that best fits the data without overfitting. However, the distribution of hypothetical exoplanets is not distinctly Gaussian, rendering the BIC less decisive in identifying an optimal k. Instead, we primarily examine the gradient of the BIC and identify points where it flattens, signifying a potential optimal number of Gaussian components.

We use k = 2 for coreless and reduced, FeO-free planets and k = 6 for regular terrestrial planets with an iron core and mantle FeO. The determination of these node counts was informed by the earlier analysis of the BIC, although we use less components to fit the reduced sub-population than predicted by consideration of the BIC (5 components). In total, we obtain 10 Gaussian modes, each characterized by a multidimensional mean representing a descriptive exoplanet composition.

From the resulting set of 10 potential candidates, we filter out compositions that either have a low weight (contributing minimally to the overall population) or are too similar to another composition in the set (e.g. marginal difference in their mantle chemistry). This process yields a final sample of five planets that exhibit maximum diversity. These compositions, detailed in Table 2, are named according to their characteristics: CORL represents the population of coreless planets; TERRA serves as an analog for Earth; ARS, richer in mantle iron but with a smaller core than Earth, can be likened to a “Super-Exo-Mars”; HERM, fully depleted in FeO yet with a relatively high SiO₂/MgO ratio, bears semblance to Mercury albeit with a smaller core; and XTREM, unlike any observed in the Solar System, is extremely reduced, leading to Si entering the core rather than the mantle, thereby enriching the mantle with MgO.

Additionally, with the core composition (Table 2), we are able to estimate the oxygen fugacity at core-mantle equilibrium (Table 3). To do so requires assumptions to be made as to the equilibrium conditions at the core-mantle interface. Here we assume core formation takes place at a fictive p-T; 2000 K and 1 bar, at which the thermodynamic properties of the pure end-members (Fe, Si; Chase 1998) and of the mixtures (Fe-Si alloys; Lacaze & Sundman 1991) are well known, facilitating calculation of fO₂. Depending on the Fe-Si-O element balance, fO₂ is determined accordingly:

1. The ratio of FeO(silicate) to Fe(metal) sets the fugacity at the core-mantle boundary: $$\mathrm{\Delta }\mathrm{I}\mathrm{W}=2\log \left(\frac{{\gamma }_{\mathrm{F}\mathrm{e}\mathrm{O}}\cdot x_{\mathrm{F}\mathrm{e}\mathrm{O}}}{{\gamma }_{Fe}\cdot x_{\mathrm{F}\mathrm{e}}}\right),$$(C.2)

   following the law of mass action. We set γ_FeO = 1.3 (O’Neill & Eggins 2002), γ_Fe = 0.9 which accounts for the dilution of pure Fe with Ni, which is reduced at higher fO₂ than is Fe.

2. The subsample of the planets that have insufficient oxygen to oxidise all of the Si present in the bulk planet need to be treated separately, since some Si partitions into the metallic core. Hence, the oxygen fugacity is set to the Si-SiO₂ (dubbed here the ‘silicon-quartz’ buffer, abbrev. SiQ): $$\mathrm{\Delta }\mathrm{S}\mathrm{i}\mathrm{Q}=\log \left(\frac{{\gamma }_{{\mathrm{S}\mathrm{i}\mathrm{O}}_2}\cdot x_{{\mathrm{S}\mathrm{i}\mathrm{O}}_2}}{{\gamma }_{Si}\cdot x_{Si}}\right).$$(C.3)

   We set ${\gamma }_{{\mathrm{S}\mathrm{i}\mathrm{O}}_2}=0.8$ (Wolf et al. 2023), γ_Si = 0.25 for XTREM and γ_Si = 0.01 for HERM (Lacaze & Sundman 1991). The fO₂ value given by Eq. C.3 can be expressed relative to the IW-buffer by adding a correction factor (ΔIW at 2000 K, 1 bar): $$\mathrm{\Delta }\mathrm{I}\mathrm{W}=\mathrm{\Delta }\mathrm{S}\mathrm{i}\mathrm{Q}-\mathrm{7.327.}$$(C.4)

3. Some planets in our sample might be coreless. This occurs if there is sufficient oxygen to oxidise all Fe, as is the case for CORL, then the FeO/FeO_1.5 ratio dictates fO₂ (e.g. O’Neill et al. 2018). The simplified condensation model does not account for the formation of Fe₂O₃. As such, we use Eq. (C.2) to provide a minimum estimate for the fO₂ at which the silicate liquid would be in equilibrium with pure Fe. For this planet archetype, it yields ΔIW ≥ −0.89.

Appendix D Effect of TiO on spectrum

TiO(g) is a strong shortwave absorber and thus might modify atmospheric spectra significantly, hence posing as a potential probe for melt composition. To test its effect, we used the composition TERRA and variend the TiO₂ abundance from 0.1 to 5 wt%. Results are displayed in Fig. D.1. As can be seen, increasing levels of TiO₂ produce stronger emission in the visible and near-infrared, but reduce emission in the mid-infrared globally. TiO predominantly absorbs/emits in the range 0.45–1 micron, coinciding with maxima in stellar emissions for G-type stars, thus it absorbs stellar radiation before it reaches lower layers. As a consequence, the contrast between the MgO and SiO bands is enhanced as the lower atmosphere – where MgO resides – is cooled. Reducing systems (ΔIW−4) show lower emission in MgO(g) and SiO(g) at 4.5 μm, but SiO(g) at 9 μm is not as strongly affected, possibly a combination of the 9 μm feature being more opaque (cf. Fig. 1) and SiO(g) being more abundant. This enhances the contrast Φ (Eq. 21) with increasing melt TiO₂ content for this redox regime. Intermediate conditions (ΔIW 0) do not adjust the shape of their spectral features to the TiO₂ abundance, but lower overall emissivity. Oxidising conditions (here ΔIW+4) suppress outgassing of TiO(g), therefore they are hardly affected.

Fig. D.1 Spectra of planets with enhanced TiO2 melt content. The underlying melt composition is assumed to be TERRA, and the TiO2 melt abundance is enhanced according to the main text. The planet under study is Earth-sized and irradiated at 2500 K.

However, TiO₂ is generally not expected to occur in large abundances in planetary mantles due to its overall low availability in the protoplanetary nebula (Lodders & Fegley 1998); it can be, however, accumulated at the surface via fractional melting or crystallisation due to its nature as an incompatible element. The lunar mare basalts for example can contain several wt% of T1O₂ (up to ~13 wt%, Giguere et al. 2000). On Earth, crustal material (formed through fractional melting of the mantle) is enriched in SiO₂ and TiO₂ while depleted in MgO, which would generate stronger SiO features (Sec. 3.4). As already noted by Zilinskas et al. (2022), such a combination could be used to distinguish evaporating bulk mantle from crustal material.

Appendix E Stellar spectra influence on pressure–temperature profile

The major cause behind the prominent inversions found earlier is the absorption of shortwave radiation (VIS and NIR), and reemission in MIR. As G-type stars – which we assumed as the fiducial case in our model – have their emission maximum in this wavelength regime, they readily cause thermal inversions (Zilinskas et al. 2022). However, some of the planets tested in Sect. 4.5, namely K2-141 b and TOI-431 b, orbit K-type stars. The spectral maxima of these stars lies deeper in the infrared, thus limiting the UVIS input to their planets atmosphere. As a consequence, the thermal inversion of mineral atmospheres around cooler stars are less pronounced (see Fig. E.1), and the SiO(g) feature is weakened (Zilinskas et al. 2022). This implies that retrieval of fO₂, as outlined in Sect. 4.5, becomes more challenging for planets around cooler stars.

Fig. E.1 Pressure-temperature profiles of atmospheres in Fig. 6.

Appendix F Effect of pressure broadening

Even though the broadening coefficients used in this study are estimated values, neglecting broadening altogether would result in even less accurate results. To quantify the extent of the broadening effect as a function of fO₂, we repeat an analysis for a 1 R_⊕, 1 M_⊕ planet with T_irr = 2500 K and the TERRA composition for the fugacity range in ΔIW from −6 to +6, shown in Fig. F.1. The opacities were constructed as outlined in Sect. 2.3.2, but restricted at a single pressure value of p = 10⁻⁸ bar. Temperature (Doppler) broadening is allowed. The unbroadenend spectra (thick lines) show significantly smaller SiO(g) features, but remain similar to the broadened lines (dashed), indicated by the difference between the respective cases. Pressure-temperature profiles change only noticeably for the highly oxidised cases ΔIW +4 and +6.

Fig. F.1 Spectra of a bulk silicate Earth-like (TERRA) planet of 1 R⊕ (as per Fig. 7) and 2500 K when no pressure broadening is considered (bold lines) vs. the broadened spectra (dashed) seen in the main paper (see Fig. 7). The unbroadenend lines are sampled at a constant pressure of 10−8 bar. The SiO features are significantly less pronounced when pressure broadening is ignored.

Appendix G Considerng whether the lava ocean can buffer the atmosphere

We assumed that the melt can be considered as an infinite reservoir compared to the vapour. To investigate this assumption, we compare the column density for an (assumed isothermal) atmosphere: $$N_{{\mathrm{O}}_2,atmo}=\frac{p_{{\mathrm{O}}_2}}{m_{{\mathrm{O}}_2}\cdot g},$$(G.1)

where p is the pressure (Pa), m the mass of the average atmospheric particle (kg) and g the gravitational acceleration (m/s²), with the column density of a melt-layer: $$N_{{\mathrm{O}}_2,melt}=\frac{\rho h}{m_{melt}}\cdot x_{{\mathrm{O}}_2},$$(G.2)

where ρ is the density of the melt (kg/m³), h is the depth of the melt column (m), and m is its average molecular weight (kg), and $x_{{\mathrm{O}}_2}$ being the molar fraction of O₂ in the melt. Approximate values for each of these parameters are: m_atmo = 32, m_melt based on melt composition but typically ~51 g/mol, ρ ~ 2600 kg/m³ and g ~ 10 m/s², $p_{{\mathrm{O}}_2}$ is based on the evaporation model at given temperature and fO₂ conditions. $x_{{\mathrm{O}}_2}$ = 0.25, based on its rough stoichiometry in the melt (which is by mole number half composed of O). The melt composition is assumed to be BSE (bulk silicate Earth, similar to the TERRA composition in our study); m_melt is therefore ~51.9 · u, where u is the atomic mass unit.

By requiring that $N_{{\mathrm{O}}_2,melt}\ge N_{{\mathrm{O}}_2,atmo}$, we can estimate the depth of the melt column that can buffer the atmosphere, labelled here as the minimal necessary ocean depth (MNOD). MNOD is obtained by finding the value h for which $N_{{\mathrm{O}}_2,melt}=N_{{\mathrm{O}}_2,atmo}$; therefore, the actual ocean depth has to exceed MNOD by at least an order of magnitude. Results are displayed in G.1.

While this estimate is incredibly simplistic and qualitative at its best, it still indicates that a.) cooler melts should be more capable of buffering the atmosphere, given the low column density of the vapour and b.) oxidising lava oceans must be of greater depth than reducing ones in order to buffer the vapour fO₂, based on the large column density of the massive O+O₂ atmospheres.

Fig. G.1 Minimal necessary ocean depth (MNOD) to buffer the atmosphere, in meters.

References

All Tables

All Figures

	Fig. 1 Opacities of major gas species sampled at $\mathrm{\Delta }\nu =0.01{\mathrm{c}\mathrm{m}}^{-1},\nu =\frac{c}{\lambda }$, at 2500 K and 0.1 bar (a typical value for the photosphere) as a function of wavelength. Note that the spectral lines of atomic O are omitted since are so weak that they would not be visible in the plot.
In the text

Fig. 2 Spread in SiO2/MgO and FeO/MgO ratios in the mantles of prospective exoplanet compositions according to the Hypatia database. See Sect. 2.6 for modelling details. Contours denote the distribution withing the planets that contain FeO in their mantles, and define the levels where the distribution contains 16%, 50%, and 84% of probability mass (from inner to outer contour). The white dots represent hypothetical exoplanet compositions selected by a GMM. The mantle ratios SiO2/MgO and FeO/MgO are shown for the solar system planets Earth (McDonough & Sun 1995), Mars (Khan et al. 2022) and Mercury (Nittler et al. 2018). Venus is assumed to be similar to Earth.

In the text

	Fig. 3 Vapour composition directly above the magma ocean as a function of temperature (left) and oxygen fugacity (right) for three representative mantle compositions in our sample: the composition for a coreless planet, CORL (top), the Earth-analogue TERRA (center) and the extremely silicon-depleted XTREM (bottom); see Table 2 for the respective compositions.
In the text

	Fig. 4 Pressure above a lava ocean as a function of composition, irradiation temperature, and oxygen fugacity (expressed relative to the iron-wüstite buffer, IW), compared to atmospheric pressures of solar-system objects (grey, dashed). Source for pressures: Catling & Kasting (2017).
In the text

Fig. 5 Temperature–pressure profiles of all simulated atmospheres. Each subfigure contains all spectra (light grey), with the ones corresponding to the Tirr/ΔIW of the given row and column being highlighted. Color corresponds to composition (Table 2). The black dotted, vertical lines denote the irradiation temperature of the planet. Points A–D refer to a selected set of atmospheres that are shown in more detail in Fig. 6.

In the text

Fig. 6 Atmospheric structure of representative combinations of ΔIW and Tirr, selected from Figs. 5A–D. The left column displays the optical depth in every atmospheric layer (color gradient), as well as an indication of the photospere (white). The latter is defined where the atmosphere has absorbed 16% (upper dotted), 50% (solid) and 84% (lower dotted) of all incoming light. The center displays the atmospheric pressure–temperature profile, with the irradiation temperature Tirr highlighted as vertical dotted line. On the right, the (volume) mixing ratios of major chemical species are shown.

In the text

Fig. 7 Spectra of all 140 simulations of this study. Each column contains all spectra of equal irradiation temperature, and each row constant composition. In each box, all spectra of constant irradiation temperature are plotted (grey), but only the ones corresponding to the composition indicated in the row are coloured. The colour-coding is made according to their oxygen fugacity (yellow most reducing, brown intermediate, light blue most oxidising). Species responsible for important spectral features are indicated by the pointers, important wavebands by coloured patches.

In the text

Fig. 8 Effect of composition on the spectrum of a mineral atmosphere at constant fO2 relative to the IW buffer (ΔIW-2) and irradiation temperature (2500 K) for an Earth-sized planet. Shown are the stellar spectrum (top), the spectrum of planetary emission (center, coloured) and the secondary occultation depth Fp/F*, the ratio between the two (bottom, colours). Flux ratios for outgassed atmospheres may overlap, as evident for the CORL and ARS compositions. Only the ultra-SiO2-poor composition XTREM is slightly distinct by showing less emission in the TiO band. Both iron-free compositions (HERM and XTREM) do not express iron lines.

In the text

Fig. 9 Effect of fO2 on the spectrum of a mineral atmosphere at constant melt composition (TERRA) and irradiation temperature (2500 K) for an Earth-sized planet. Shown are the stellar spectrum (top), the spectrum of planetary emission (center, coloured) and the secondary occultation depth Fp/F*, the ratio between the two (bottom, colours). The emission from metallic species (Fe, Mg, MgO, TiO) in the visible distorts the overall shape of the spectral energy distribution (SED) compared to a pure black body. This effect is pronounced in reducing atmospheres, whereas oxidised ones follow a black body more closely except for a handful of strong lines.

In the text

	Fig. 10 Influence of alkali metals (Na, K) on structure and spectrum of mineral atmospheres of an Earth-sized planet. On the left hand side, the difference in the planet-to-star flux ratio between an atmosphere with and without alkali metals is shown. The composition is TERRA (similar to bulk silicate Earth, BSE) with 0.36 wt% Na2O and 0.029 wt% K2O added (McDonough & Sun 1995).
In the text

	Fig. 11 Planet-to-star flux ratio spectra of atmospheres above a binary MgO-SiO2 melt, with varying mass ratios of MgO/SiO2 (shown in inset). We tested the conditions Tirr = 2500 K for an Earth-sized planet at a) Δ IW−4, b) IW and c) Δ IW+4.
In the text

Fig. 12 Φ as function of composition for a binary SiO2-MgO melt (lines, colored by fO2), assuming Tirr = 2500 K, 1 M⊕ and 1 R⊕, compared to Φ for a atmosphere-free, pure black body planet (black dashed line). The Φ of the multicomponent-melts of the archetypes in Table 2 are overlaid (markers). The distribution of SiO2/MgO of hypothetical exoplanet mantles is indicated as histogram (top), with the 16, 50 and 80 percentile shown as grey dotted lines.

In the text

	Fig. 13 Secondary eclipse depth and atmospheric profiles produced by the planet/melt compositions from Table 2, evaluated with the respective fO2 (Table 3).
In the text

Fig. 14 Emission spectra simulated for various existing lava planets, K2-141-b (green), 55-Cnc-e (blue), TOI-431-b (purple) and HD 213885-b (orange), assuming constant composition (TERRA, Table 2) and their corresponding stellar spectra (Zilinskas et al. 2022). Two diagnostic features for fO2 determination are highlighted: MgO at 6 μm and SiO at 9 μm. The contrast between the fluxes of the respective bins are shown on the right, for each planet, respectively, once as the contrast expected from a perfect observation (solid coloured lines) and the respective standard deviations (1σ) for various number of observed transits with MIRI (black, simulated with pandexo Batalha et al. 2017).

In the text

Fig. A.1 Vapour pressures and thermodynamics of a silicate melt of BSE composition. Upper: Similar to Fig. 3, but for a melt of BSE composition (including Na2O). The dominance of Na(g) at intermediate redox states (ΔIW-3 to ΔIW+2) is evident, which is consistent with the findings of other studies (Schaefer & Fegley 2004; Miguel et al. 2011; Wolf et al. 2023; van Buchem et al. 2023). Upper center: Concentration of oxides in the silicate melt. FeO concentration deplets under reducing conditions, as Fe is formed. This is evident in Lower center: here, the bulk system composition is shown, highlighting the conversion of FeO to Fe. Lower: Activity of Fe and Fe2O3 in the melt, compared to the far more dominant SiO2 and MgO.

In the text

	Fig. D.1 Spectra of planets with enhanced TiO2 melt content. The underlying melt composition is assumed to be TERRA, and the TiO2 melt abundance is enhanced according to the main text. The planet under study is Earth-sized and irradiated at 2500 K.
In the text

	Fig. E.1 Pressure-temperature profiles of atmospheres in Fig. 6.
In the text

	Fig. F.1 Spectra of a bulk silicate Earth-like (TERRA) planet of 1 R⊕ (as per Fig. 7) and 2500 K when no pressure broadening is considered (bold lines) vs. the broadened spectra (dashed) seen in the main paper (see Fig. 7). The unbroadenend lines are sampled at a constant pressure of 10−8 bar. The SiO features are significantly less pronounced when pressure broadening is ignored.
In the text

	Fig. G.1 Minimal necessary ocean depth (MNOD) to buffer the atmosphere, in meters.
In the text