Exoplanet Volatile Carbon Content as a Natural Pathway for Haze Formation

The compositions of bodies in the solar system point to an array of chemical environments that were present during the formation of planets and their building blocks. Among these, a fundamental change in chemistry occurred in the solar nebula, the protoplanetary disk that circled our own Sun, at so-called "snow lines." These chemical transitions mark the region outside of which the pressures and temperatures are such that a given molecular species would exist as a solid and inside of which that solid sublimates to the vapor. An important consequence of these locations is that corresponding species would be abundant in solids that form outside a transition point, while these components would be scarce in solids located inside.

Of particular importance to planet formation are the locations of the water and CO snow lines as these have traditionally been considered to be the primary carriers of oxygen and carbon in protoplanetary disks (Öberg et al. 2011). Importantly, the CO snow line is located many tens of astronomical units (au) from the star (Qi et al. 2013), only where temperatures are <30 K, low enough for CO ice to form. This would seemingly make many planets, including the known rocky exoplanets, most of which are found closer to their host star than the Earth, carbon poor.

This picture assumes that all carbon is locked up in CO (or similarly volatile species such as CH₄ or CO₂). However, recent work has shown that a significant amount of carbon in the interstellar medium (ISM), up to 60% of cosmic carbon (Mishra & Li 2015), is carried by refractory organics (Bergin et al. 2015; Gail & Trieloff 2017). These organics, hereafter called "soot," are predominantly macromolecular in nature and comprised of hydrocarbons/organic species (Alexander et al. 2013) and are the product of disequilibrium reactions in the interstellar medium and/or outer protoplanetary disk. These soots will remain solid at temperatures up to ∼500 K (Li et al. 2021) and have the important property that when they are heated above their destruction temperature, they decompose into simpler, more volatile species. That is, their vaporization is irreversible. This leads to the concept of the "soot line" in a protoplanetary disk, a location close to the star, outside of which refractory carbon would be available for incorporation into solid planetary materials, but inside of which it is absent (Kress et al. 2010; Li et al. 2021). A unique property of the soot line is that any carbon contained in vapor that mixes outward remains in the gas and does not freeze out again as expected around traditional snow lines (Ros & Johansen 2013).

Planets that form outside of the soot line can thus be carbon rich, leading to highly reducing conditions during their early evolution, particularly if they formed interior to the water snow line (and thus did not have access to another major solid hydrogen/oxygen carrier). Figure 1 indicates that such planets very well may exist, showing the temperature profiles of millimeter-sized pebbles extrapolated from ALMA measurements of protoplanetary disks of ages associated with potential incipient planet formation; the corresponding locations of the soot and snow lines are also shown. The area in between these two locations marks the portion of the disk where planetary materials would be relatively rich in carbon but chemically reduced because of the preservation of refractory organics but loss of water to the gas. A histogram of known Earth-like planets and super-Earths and their semimajor axes is also shown, with many being found in this important region. If these planets formed predominately from solids from this region, then they would form from carbon-rich/water-poor material. It has been suggested that many of these systems formed at larger distances and migrated inwards at earlier stages (Ida & Lin 2008; Coleman & Nelson 2014; Izidoro et al. 2017). This calls into question the correspondence shown in Figure 1. However, other models argue for "in situ" formation (e.g., Lee et al. 2014; Batygin & Morbidelli 2023), and as we discuss below, at earlier stages the nebular gas is warmer with a more distant soot line. Observational tests of these competing ideas are thus strongly desired.

Zoom In Zoom Out Reset image size

Given the above, one may expect the Earth to have formed with a high abundance of carbon, yet it is severely carbon depleted (Bergin et al. 2015). However, there are ways that carbon can be lost from solids during planet formation. Li et al. (2021) demonstrated that temperatures early in disk evolution can be much higher than shown in Figure 1 provided the mass accretion rate from the disk to the star is sufficiently high. This could push the soot line out to beyond 1 au during the first 1 Myr; its low carbon inventory could reflect that much of Earth's primary materials were assembled during this phase (Li et al. 2021). Alternatively, Hirschmann et al. (2021) argued that heating of planetesimals from radioactive decay (²⁶Al in our solar system) can reach high enough temperatures to destroy the organics and drive off volatiles. If Earth's building blocks formed later than suggested by Li et al. (2021), the low carbon content of the Earth could then be a result of the differentiation and thermal metamorphism of its progenitors, which we readily see in the meteorite record (see also Grewal 2022).

It is important to note, however, that accretion rates through disks vary by orders of magnitude (e.g., Hartmann et al. 2016), and thus, in many protoplanetary systems the soot line would be very close to the star and well inside the location where planets are found. Even if not, carbon-rich organic pebbles can be replenished by inward drift from regions that never saw the inside of the soot line, a process that potentially did not occur in the solar system as a result of Jupiter's formation (Kruijer et al. 2020). Further, not all planetary systems are expected to form with as high an abundance of ²⁶Al that our solar system had (if any at all; Ciesla et al. 2015; Lichtenberg et al. 2019), which would limit the heating that planetesimals experienced prior to their accretion into planets. Thus, it is possible, perhaps even likely, that a significant fraction of the Earth-sized planets and super-Earths were assembled from rocky materials with large carbon inventories.

This is the thesis that we explore in this Letter. Here we focus on worlds that are silicate rich but have greater volatile inventories than seen in the Earth with 0.1–1.0 wt% refractory carbon present in their mantles. We will show that one implication of this composition is that atmospheric hazes, which appear to be present in numerous systems (e.g., Kreidberg et al. 2014; Crossfield & Kreidberg 2017; Gao et al. 2020; Dymont et al. 2022), would be a natural outcome. In Section 2 we provide the baseline model of the mantle composition. In Section 3 we explore the atmospheric composition of these planets and apply a photochemical model of haze formation. In Section 4 we present basic predictions of this model and discuss the implications.

We explore the consequences for the atmospheres of planets that form in this critical region where the planets are chiefly comprised of refractories and organics. The bulk of meteoritic organic material is insoluble in typical solvents and is thought to be macromolecular in form with a typical composition of C₁₀₀H_75–79O_11–17N_3–4S_1–3 (that is, normalized relative to 100 carbon atoms) (Alexander et al. 2017; Glavin et al. 2018). We assume this material is representative of the soot composition. In the interstellar medium, up to 60% of cosmic carbon (Mishra & Li 2015) is carried by refractory organics, and the bulk refractory organic carbon composition in cometary material is comparable to that of ISM material (Bergin et al. 2015; Gail & Trieloff 2017). As such, the refractory organic carbon content of planet-building materials forming beyond the soot line is expected to be high and likely comparable to that in comets. In this case, Comet 67P, which is similar to Comet Halley, had refractory organic carbon content that is ∼6× that of CI chondrites (Bardyn et al. 2017), which have 2–4 wt% carbon (Pearson et al. 2006). Models of dust emission in protoplanetary disk systems also commonly assume carbonaceous dust comprised of refractory organics is present in abundances consistent with the interstellar medium (Pollack et al. 1994; D'Alessio et al. 2001; Birnstiel et al. 2018). Further chemical evidence of soot destruction at the soot line may have been recently detected in the JWST spectra of a young protoplanetary disk (Tabone et al. 2023).

We therefore explore outcomes where the composition of the starting material is 0.1 and 1.0% wt% soot, with the majority of the mass in these planets being the more refractory silicates and metals. A final state with 0.1%–1% soot by mass is conservative in that it assumes substantial volatile loss during formative stages (Hirschmann et al. 2021) from the assumed initial values of 12–24 wt%. While rich in reduced carbon, these silicate-dominated rocky bodies are not the extreme "carbon planets" discussed by other researchers (e.g., Madhusudhan et al. 2012; Unterborn et al. 2014).

We first assume that the planet accreted without water as this represents an extreme case, but we also explore solutions where some water is provided in the form of hydrous silicates. Even in the most water-poor cases, water in the mantle is generated via the oxygen in silicates and the hydrogen from organics/atmosphere; this water can be released to the atmosphere (this is discussed in the Appendix A). We perform several calculations to predict the end-state compositions of nascent planetary atmospheres (see Appendix). We adopt rocky masses of 0.3 M_⊕, 1.0 M_⊕, and 3 M_⊕. Planets are assumed to accumulate a nebular atmosphere consisting of hydrogen that increases with their mass, following Stökl et al. (2015), which can be significant when the planet mass exceeds that of Earth. For all objects, we first determine the geochemical equilibrium between the outgassing from a molten mantle and the overlying atmosphere, using the outgassing model of Gaillard et al. (2022). Based on the assumed ratio of refractories to organics we calculate the oxygen fugacity and the compositions of coexisting atmosphere and molten volatile-bearing silicate. In each case, the resulting mantle contains significant amounts of carbon. The destruction of the soot results in the formation of C species dissolved in the molten silicate and in some cases, graphite, along with outgassed C-O-H vapor.

The results from these equilibrium calculations, which represent the base of the overlying atmosphere, are given in Table 1. For the lower-mass planets (0.3 M_⊕ and 1.0 M_⊕) their initial atmospheres are dominated by H₂ and CO, but significant (a few %+) amounts of CH₄ are typically present. The presence of a massive hydrogen envelope in the 3 M_⊕ planet changes the underlying equilibrium such that the atmospheres are H₂ and CH₄ dominated. Essentially, an appreciable fraction of C released from the mantle is processed into methane in these atmospheres. The stability of CH₄ rather than CO, even at high temperatures, is a product of the high pressure of these thick atmospheres, as the reaction CO + 3 H₂ $\longrightarrow$ CH₄ + H₂O has a negative volume change. Surprisingly, the inclusion of water, across all modeled planetary masses, does not alter this composition as the high atmospheric carbon content continues to favor some hydrocarbon production. These results demonstrate that abundant hydrocarbons are outgassed to the base atmosphere across a range of planet masses and even in the case of high (compared to the Earth) water content. This is notable as the ingredients for haze production have been linked to the presence of CH₄ and other simple hydrocarbons, which are processed by UV photons in the upper atmosphere (e.g., Miller-Ricci Kempton et al. 2012; Morley et al. 2013; Kawashima & Ikoma 2018; Lavvas et al. 2019). We note that other pathways of haze production have been suggested for CO and CO₂ dominated atmospheres (He et al. 2020), alongside the potential importance of trace species, e.g., sulfur (Zahnle et al. 2016; Gao et al. 2017). In this case we concentrate on CH₄, which we show can be quite abundant, even in cases for which temperatures would nominally convert CH₄ into CO or CO₂.

Table 1. Model Assumptions and Base Atmosphere Properties/Composition

% Soot a	% H2O a	Mp	Mp,soot	M ${}_{p,{{\rm{H}}}_{2}}$	P	log10(f${}_{{{\rm{O}}}_{2}}$)	N2	H2O	H2	CO2	CO	CH4	Element Fractions
(by mass)		(M⊕)			(MPa)		(Mixing Ratios)						H	O	C
1.0	0.0	0.30	0.003	8.7 × 10−8	84.08	−10.47	0.00	0.04	0.23	0.03	0.66	0.03	0.314	0.351	0.335
1.0	0.0	1.00	0.01	1.5 × 10−4	98.47	−10.55	0.00	0.05	0.33	0.02	0.52	0.07	0.459	0.271	0.270
1.0	0.0	3.00	0.03	1.4 × 10−2	2970.00	−14.80	0.00	0.00	0.59	0.00	0.00	0.41	0.874	0.000	0.126
0.1	0.0	0.30	0.003	8.7 × 10−8	71.28	−10.19	0.00	0.00	0.01	0.05	0.94	0.00	0.007	0.508	0.484
0.1	0.0	1.00	0.01	1.5 × 10−4	70.52	−10.26	0.01	0.01	0.07	0.04	0.87	0.00	0.082	0.471	0.446
0.1	0.0	3.00	0.03	1.4 × 10−2	866.10	−14.80	0.00	0.00	0.95	0.00	0.00	0.05	0.977	0.000	0.023
0.1	1.0	0.30	0.003	8.7 × 10−8	109.34	−10.36	0.00	0.08	0.50	0.01	0.29	0.10	0.668	0.165	0.167
0.1	1.0	1.00	0.01	1.5 × 10−4	188.16	−10.43	0.00	0.06	0.42	0.01	0.28	0.21	0.679	0.135	0.186
0.1	1.0	3.00	0.03	1.4 × 10−2	1120.00	−14.84	0.00	0.00	0.96	0.00	0.00	0.04	0.983	0.000	0.017

^a Percent mass added to M_p .

Download table as:  ASCIITypeset image

3.1. Baseline Atmospheric Composition

The results from Table 1 refer to the composition at the base of atmosphere, but tying this to haze production requires quantitative models of chemistry and irradiation of the upper atmosphere. To determine whether we expect high CH₄ abundances to persist to altitude, we first investigate the atmospheric composition in thermochemical equilibrium as a function of atmospheric temperature and pressure (methods are described in the Appendix). The results from these calculations are provided in Figure 2. Here we include the atmospheric temperature as a proxy for the planet's orbital distance. The left-hand panel provides the atmospheric composition as a function of planet mass and atmospheric temperature for planets with 0.1% soot by mass. For lower-mass planets, a greater diversification of the carbon content is seen, resulting in CO- and CO₂-dominated atmospheres. Despite this, at lower temperatures, some transitioning to CH₄ is observed, yielding significant CH₄/CO₂ ratios. This is of interest as the simultaneous detection of CH₄ and CO₂ has been suggested as a biosignature (Krissansen-Totton et al. 2018); for these planets it may instead be a natural outcome of formation. However, for planets with ≤1 M_⊕, detailed calculations of mantle evolution are needed to understand the impact of higher soot content in the upper mantle (compared to the Earth) and the resulting effect on atmospheric composition and evolution.

Zoom In Zoom Out Reset image size

Our solutions in Table 1 refer to young planets, but there can be significant mantle and atmospheric composition evolution over billions of years. Most directly this would be the loss of the primary H₂ atmosphere and development of a secondary atmosphere. Thus, these models are not necessarily predictive for the composition of observed, and evolved, exoplanets. Among our modeled planets, the 3 M_⊕ is dominated by its hydrogen envelope and is expected to experience the least atmospheric evolution. We therefore focus on these planets as representative of an evolved outcome of sub-Neptune atmosphere formation. These planets are more abundant in our collection of known exoplanets and are more representative of those atmospheres likely to be characterized by JWST in the next few years. The right-hand panel of Figure 2 shows the equilibrium atmospheric composition for a 3 M_⊕ planet with variable soot and water content at a pressure of 1 mbar (approximately the pressure probed by transmission spectroscopy measurements). Here, we find that the atmospheres are hydrogen/methane dominated and in some instances with significant concentrations of acetylene (C₂H₂) and ethylene (C₂H₄). These results show that methane and other hydrocarbons can persist at high abundance, even in high-temperature and low-pressure conditions, due to the effectively elevated C/O ratio provided by the soot-rich mantle.

3.2. Implementation of Haze Model

Models of haze formation have been developed for exoplanetary atmospheres based upon irradiation of methane and other carriers. We apply one such model including chemical kinetics, photochemistry, and haze formation to our 3 M_⊕ planet with 0.1% soot and no water. We stress that this hydrocarbon-based haze model is for illustrative purposes. Our calculations show that these atmospheres will be methane rich. But nitrogen and sulfur are carried alongside carbon within soot (Alexander et al. 2012). Thus, other chemical solutions for hazes are possible. Regardless, methane will be present in these systems in abundance.

The baseline haze model is discussed in Appendix B. We specifically model the planet at an equilibrium temperature of 600 K placed in orbit around an M-dwarf host star to align its properties with sub-Neptune exoplanet targets that will be observed with JWST during its first year of operations. The resulting chemical abundance profiles are presented in Figure 3, which demonstrate that these atmospheres readily produce hazes via hydrocarbon polymerization channels.

Zoom In Zoom Out Reset image size

4.1. Current Exoplanet Landscape

4.2. Predictions for JWST and Mantle Water Content

The predicted transmission spectrum of a volatile rich world from 0.3 to 30 μm is shown in Figure 3. The baseline spectrum (no hazes) exhibits strong methane features. In the presence of hazes, the spectrum is significantly muted, in line with existing observations of sub-Neptune atmospheres (Knutson et al. 2014; Kreidberg et al. 2014; Libby-Roberts et al. 2020). These models are for a 3 M_⊕ planet, and existing spectroscopic data on atmospheres currently do not detect discernible atmospheric features toward these lower-mass planets (de Wit et al. 2018; Diamond-Lowe et al. 2020; Libby-Roberts et al. 2022). The expectation of our model results is that similar features would be anticipated for slightly more massive planets (e.g., 5 M_⊕).

We find that even adding comparable amounts of soot and water does not dilute the presence of a rich, methane-dominated atmosphere (Figure 2). To determine the robustness of this model, we vary the water content within our baseline model for a 3 M_⊕, 600 K equilibrium temperature planet, with 0.1% soot and no water. The water content in the mantle of this planet is 0.23 wt% with an atmospheric water mixing ratio of 8 ×10⁻⁴ and a total water column depth of 1.5 ×10²⁶ molecules/cm². The atmospheric mixing ratio of water provides the floor value for varying the water content (i.e., 1×O). This might occur if some material is supplied to the young silicate-/soot-rich planet from beyond the water ice line. To model this we raise the mixing ratio of water by integer units noted in Figure 2 but maintain oxygen mass balance in the calculation by lowering other carriers uniformly. We otherwise maintain the temperature–pressure (T–P) profile, eddy diffusion coefficient, and stellar spectrum as for the modeling described in Section 3.2. This model is not self-consistent as we are not modeling all the geochemical steps described in Section 3.1 and do not account for the potential effect of this additional water vapor on the physical/chemical profile. Rather, we vary the atmospheric water content to capture the impact of additional oxygen in the atmosphere on haze production and approximately simulate scenarios in which additional water is provided from beyond the snow line.

The impact of decreasing haze production on transmission spectra is shown in Figure 4. Our atmospheric modeling suggests that substantial H₂O incorporation during planet formation that translates to additional atmospheric oxygen exceeding a fortyfold increase, such as what would be expected for planets born beyond the snow line, is needed to throttle haze production in this model of a soot-rich planet. Overall, we find that these results indicate high atmospheric methane concentration and ample haze production, and transmission spectra dominated by haze and methane features are robust even when considerable excess oxygen is added to the atmosphere from water incorporation at the time of formation.

Zoom In Zoom Out Reset image size

4.3. Volatile, Rich Worlds

This research comes from an interdisciplinary collaboration initiated by the Integrated NSF Support Promoting Interdisciplinary Research and Education Program through grant AST1344133. Additional funding has been provided by National Aeronautics and Space Administration grants 80NSSC19K0959 (to M.M.H.), XRP NNX16AB48G (to G.A.B.), XRP 80NSSC20K0259 (to E.A.B. and F.J.C.), and 80GSFC21M0002 (to S.T.B); and National Science Foundation Grant AAG 2009095 (to E.M.-R.K. and supporting D.J.T.). S.T.B. also received support from the GSFC Sellers Exoplanet Environments Collaboration (SEEC), which is funded in part by the NASA Planetary Science Division's Internal Scientist Funding Model. This project is supported, in part, by funding from Two Sigma Investments, LP to EAB. Any opinions, findings, and conclusions or recommendations expressed in this material are those of the authors and do not necessarily reflects the views of Two Sigma Investments, LP.

Facility: JWST. -

Software: HELIOS (Malik et al. 2017, 2019), Exo-Transmit (Kempton et al. 2017; Teal et al. 2022).

Initial calculations begin with different fractions of C–H–O soot (C₁₀₀H₇₇O₁₅) and silicate, assuming that the silicate initially contains 16 wt.% FeO, a comparatively oxidized assumption similar to the mantle of Mars.

We consider a planet with Earth-like proportions of core (33% by mass) and silicate (67%), and for the condensed portion of the planet, the mass–radius relationship is given by the parameterization a + b M_⊕+c ln (M_⊕), where a = 0.9868, b = 0.0231, and c = 0.2599 are empirical coefficients taken from mass–radius relationships of Earth-like planets from Santerne et al. (2018) and Zeng et al. (2019). This relationship allows for the calculation of the gravitational acceleration at the surface for each planetary mass.

For equilibration between molten silicate and overlying atmosphere, we employ the thermodynamic outgassing model presented in Gaillard et al. (2022). ⁹ This model calculates partial pressures of outgassed species (H₂, H₂O, CO₂, CO, CH₄) based on assumed total mantle mass, volatile content, temperature (1773 K for our calculations), and oxygen fugacity, but the critical values, passed to the atmospheric calculations described below, are the total elemental masses of outgassed elements (principally C–H–O).

In the calculation of Gaillard et al. (2022), the calculated masses of silicate outgassed volatiles do not conserve oxygen mass balance. This is because volatile inputs are taken only as oxidized species (H₂O, CO₂) but output as both oxidized and reduced C–H–O species. This required several adjustments. First, input of reduced volatiles, such as hydrogen gas or "soot" required adjustments according to reactions such as H₂ + FeO = Fe + H₂O.

Second, because the Gaillard et al. (2022) calculator assumes that the silicate is effectively an infinite reservoir, leaving oxygen fugacity unchanged even as oxidized input species are converted to a combination of oxidized and reduced species, we took an iterative approach. After each iterative step, we recalculated the concentration of FeO in the silicate by enforcing O mass balance and then calculated ${{\rm{f}}}_{{O}_{2}}$ based on an empirical curve derived from Frost et al. (2008): ΔIW = 0.8763 $\cdot \mathrm{ln}$(FeO) − 3.80 (IW = iron-wüstite buffer), where FeO is in units of weight percent. Oxygen fugacities were bounded at a minimum of ΔIW = −6, below which melt FeO is effectively zero, making mass balance ineffective. Convergence was accepted when the resulting ${{\rm{f}}}_{{O}_{2}}$ (oxygen fugacity) differed by less than 0.03 log units from the previous iteration, and generally four to five iterative steps were required.

The temperature for the principal calculations was selected as 1773 K because these are conditions close to the experimental constraints on volatile solubilities in silicate liquid employed by the calculator. Greater temperatures are expected at the surfaces of magma oceans, particularly for larger planets, and this will affect volatile speciation. For example, at higher temperatures, CO and H₂ are favored relative to CH₄ and H₂O. This effect is not so consequential for the purposes of the present calculation as the information from the outgassing calculation that is passed to the atmospheric calculations described below is that of elemental abundances rather than gaseous species. To illustrate the temperature sensitivity of our calculations, we include one calculation at 2773 K for the case of 0.1% soot and 1 wt% H₂O. As shown in Table A1, the resulting elemental abundances of the initial outgassed atmosphere are not consequentially different for the same bulk composition at 1773 K.

Table A1. Model Assumptions and Base Atmosphere Properties/Composition with Variable Temperatures

% Soot a	% H2O a	Mp	Mp,soot	M ${}_{p,{{\rm{H}}}_{2}}$	P	T	log10(f${}_{{{\rm{O}}}_{2}}$)	Element Fractions
(by Mass)		(M⊕)			(MPa)	(K)		H	O	C
0.1	1.0	0.3	0.003	8.7 × 10−8	109.3	1773	−10.36	0.668	0.165	0.167
0.1	1.0	0.3	0.003	8.7 × 10−8	105.9	2773	−4.85	0.666	0.167	0.167
0.1	1.0	1.0	0.01	1.5 × 10−4	188.2	1773	−10.43	0.679	0.135	0.186
0.1	1.0	1.0	0.01	1.5 × 10−4	186.0	2773	−4.85	0.674	0.137	0.189
0.1	1.0	3.00	0.03	1.4 × 10−2	1120.0	1773	−14.84	0.983	0.000	0.017
0.1	1.0	3.0	0.03	1.4 × 10−2	1120.0	2773	−9.07	0.983	0.000	0.017

^a Percent mass added to M_p.

Download table as:  ASCIITypeset image

A rigorous treatment of volatile reservoirs developed during planetary differentiation would include sequestration of C and H in the core as both elements are siderophile (Hirose et al. 2019; Li et al. 2021). Developing a model that includes partitioning of C and H into the core would require addressing the processes of accretion and core segregation that precede the initial conditions of our model, which begins with a full-formed planet from which the core has already segregated. However, the absence of this reservoir does not have a strong effect on the modeling we present because the amount of soot in each calculation is a free variable. For a fixed supply of volatile materials (i.e., soot, ice) to a planet of a given mass, the capture of volatiles by the core would reduce their masses residing in the mantle and atmosphere, thereby ultimately diminishing the total atmospheric pressure. Core segregation does not change the speciation of carbon in the mantle, which occurs as accessory phases including diamond, graphite, iron carbide, and C–H–O fluid (e.g., Frost et al. 2008). Therefore, for the purposes of investigating methane outgassing and haze production, a model that incorporates core segregation and a greater amount of accreted soot would give essentially the same results as a model with no core segregation and a smaller amount of soot.

To model the composition of the portion of the atmosphere that would be observable via spectroscopic techniques with JWST, we apply two different methodologies. The first is a chemical equilibrium calculation of atmospheric abundances as a function of atmospheric pressure. The second is a chemical kinetics calculation that includes photolysis reactions and photochemical production of important hydrocarbon haze precursors and thus haze.

For the chemical equilibrium calculation, we start from the surface composition at the atmosphere–mantle boundary (Table 1), calculated as described in Appendix A. We then derive the underlying elemental abundances (i.e., H₂O → 2 H + 1 O) of H, C, N, O, S, and Ar. From these abundances, we re-derive thermochemical equilibrium as a function of temperature and pressure using the Gibbs free energy minimization techniques described in Mbarek & Kempton (2016). We perform our calculations over a pressure range of 1 μbar–100 bar and for temperatures from 300 to 1200 K for a set of 69 molecules made up of H, C, N, O, and S (and Ar). An example of the resulting chemical equilibrium abundances versus pressure is shown in Figure 5. We additionally show the chemical equilibrium abundances at a pressure of 1 mbar (approximately the pressure level probed in transmission spectroscopy) for planets across our entire model domain in Figure 2 in the main text.

Zoom In Zoom Out Reset image size

For the chemical kinetics modeling, we first must generate realistic T–P profiles for the atmospheres in question. (This step is unnecessary for the chemical equilibrium modeling, described above, because in that case the chemical composition depends uniquely on the local temperature and pressure of the gas rather than the full vertical T–P profile.) We use the open-source HELIOS ¹⁰ code (Malik et al. 2017, 2019) to calculate temperature–pressure profiles in radiative convective equilibrium. We generate T–P profiles for the 3 M_⊕ planet, which for reasons already discussed in Section 3.1, is the scenario for which we believe our modeled atmospheres are most representative of the evolved planets that will typically be observed with JWST. We model planets with equilibrium temperatures of 600, 900, and 1200 K (we focus on the 600 K model in the main text), set by selecting the planet's orbital semimajor axis assuming zero albedo and fully efficient day–night heat redistribution. The pressure at the bottom of the atmosphere is set to 10³ bar. The host star properties and spectrum are selected to match the M-dwarf star GJ 876 (T_eff = 3300 K, R_⋆ = 0.367R_⊙) as representative of a typical system that would be observed with JWST.

The resulting HELIOS T–P profiles are then passed into a chemical kinetics code to calculate atmospheric abundances of gas-phase species and hydrocarbon haze as a function of altitude. As the impact of photodissociation is particularly pronounced at low pressures beyond the pressure cutoffs commonly used in radiative transfer models (here: 10⁻⁷ bar), we extrapolate the HELIOS T–P profiles as isothermal to 10⁻⁹ bar. We use the version of the Atmos photochemistry code described in Harman et al. (2022), with the addition of carbon-bearing species and chemical reactions up to C–4 (C₃H₂, C₃H₃, C₃H₄, C₄H₂, C₄H₃, C₄H₅) and nitrogen-bearing species and reactions (N₂, N, NH, NH₂, NH₃, N₂H, N₂H₂, N₂H₃, CN, NCO, HCN, HNO, HNCO, NO, H₂CN, HC₃N, C₂H₃CN, CH₂NH, CH₂NH₂, CH₃NH₂, CH₂CN, CH₃CN) sourced from Tsai et al. (2021). We additionally account for the formation of organic haze using the fractal haze model from Wolf & Toon (2010) and Arney et al. (2016, 2017) adapted for a H₂-dominated atmosphere (Parmentier et al. 2013). Haze formation is primarily initiated by CH₄ photolysis, which catalyzes the formation of complex organic molecules in the atmosphere. Our chemical network cannot capture the full complexity of reactions occurring among all of these high-order hydrocarbon molecules. We instead follow a common practice of selecting lower-order haze "precursor" species from our chemical network that are formed high up in the atmosphere. For the current work we select polyacetylene (C_2nH₂; e.g., Allen et al. 1980; Wilson & Atreya 2003; Lavvas et al. 2008) and allene (CH₂CCH₂) polymerization (Pavlov et al. 2001) pathways, both proceeding through reactions with the ethynyl radical C₂H, and a nitrogen-bearing copolymer pathway based on cyanoacetylene HC₃N (Lavvas et al. 2008; Krasnopolsky 2009) for haze production:

$$\begin{eqnarray}&&{{\rm{C}}}_{4}{{\rm{H}}}_{2}+{{\rm{C}}}_{2}{\rm{H}}\to \mathrm{Polymer}+{\rm{H}}\end{eqnarray} \tag{ B1 }$$

$$\begin{eqnarray}&&{\mathrm{CH}}_{2}{\mathrm{CCH}}_{2}+{{\rm{C}}}_{2}{\rm{H}}\to \mathrm{Polymer}+{\rm{H}}\end{eqnarray} \tag{ B2 }$$

$$\begin{eqnarray}&&{\mathrm{HC}}_{3}{\rm{N}}+{{\rm{C}}}_{4}{{\rm{H}}}_{3}\to \mathrm{Co}\mbox{--}\mathrm{Polymer}.\end{eqnarray} \tag{ B3 }$$

We assume a 100% conversion efficiency into haze. Once hazes form in the photochemistry model, they scatter and absorb incoming UV photons, which ultimately self-regulates the formation of additional haze. Aerosol particles form as Mie scatterers that grow and coagulate into fractal aggregate particles composed of monomers of a fixed size of 50 nm. Haze optical properties for spherical and fractal aggregate particles were calculated with the mean field approximation model described in Rannou et al. (1999) and Botet et al. (1997) assuming Titan tholin complex refractive indices from Khare et al. (1984). The irradiating host star is again selected to be GJ 876, using its UV spectrum from the MUSCLES catalog (France 2016; France et al. 2016). ¹¹ We assume a uniform Eddy diffusion coefficient of K_zz = 6 × 10⁸ cm² s⁻¹, similar in range as previous studies (Kawashima & Ikoma 2018; Tsai et al. 2021; Harman et al. 2022). While the choice of K_zz influences particle coagulation and atmospheric mixing, we forgo a detailed discussion and note that all atmospheric models we generated produced significant amounts of haze for the complete range of K_zz values we tested (5 × 10⁸–5 × 10¹⁰ cm² s⁻¹). The atomic composition determined above in the chemical equilibrium modeling was scaled to preserve the relative abundance ratios while introducing a solar metallicity abundance of He, which Atmos uses as a (required) nonreactive filler gas. The planet's gravity at the 10³ bar level and radius were set to 1481.86 cm s⁻² and 1.41 R_⊕, respectively.

Finally, we model the transmission spectra of the resulting atmospheres. For this we use the Exo-Transmit code (Kempton et al. 2017), as modified in Teal et al. (2022), to generate transmission spectra from the vertical abundance profiles output by the chemical kinetics code. Haze opacities are included in this version of Exo-Transmit, which depend on the haze particle radius. We use an identical set of hydrocarbon haze optical properties for all haze particles in the atmosphere, regardless of which of the three precursor formation pathways generated the haze. In Figure 3 of the main text, we show versions of the transmission spectra with the haze opacity included and removed, emphasizing the impact of hazes on muting/obscuring spectral features.