Exoplanet secondary atmosphere loss and revival
1 Exoplanet secondary atmosphere loss and revival
Edwin S. Kite and Megan Barnett . Abstract.
The next step on the path toward another Earth is to find atmospheres similar to those of Earth and Venus – high-molecular-weight (secondary) atmospheres – on rocky exoplanets. Many rocky exoplanets are born with thick (> 10 kbar) H -dominated atmospheres but subsequently lose their H ; this process has no known Solar System analog. We study the consequences of early loss of a thick H atmosphere for subsequent occurrence of a high-molecular-weight atmosphere using a simple model of atmosphere evolution (including atmosphere loss to space, magma ocean crystallization, and volcanic outgassing). We also calculate atmosphere survival for rocky worlds that start with no H . Our results imply that most rocky exoplanets orbiting closer to their star than the Habitable Zone that were formed with thick H -dominated atmospheres lack high-molecular-weight atmospheres today. During early magma ocean crystallization, high-molecular-weight species usually do not form long-lived high-molecular-weight atmospheres; instead they are lost to space alongside H . This early volatile depletion also makes it more difficult for later volcanic outgassing to revive the atmosphere. However, atmospheres should persist on worlds that start with abundant volatiles (for example, waterworlds). Our results imply that in order to find high-molecular-weight atmospheres on warm exoplanets orbiting M-stars, we should target worlds that formed H -poor, that have anomalously large radii, or which orbit less active stars. Significance statement.
Earth and Venus have significant atmospheres, but Mercury does not. Thousands of exoplanets are known, but we know almost nothing about the atmospheres of rocky exoplanets. A large fraction of the known rocky exoplanets were formed by a sub-Neptune to super-Earth conversion process during which planets lose most of their H -rich (primary) atmospheres and are reduced in volume by a factor of >2. It remains unclear whether such a gas-rich adolescence increases or decreases the likelihood that super-Earths will subsequently exhibit a H -poor (secondary) atmosphere. Here we show that secondary atmospheres exsolved from the magma ocean are unlikely to be retained by super-Earths, but it is possible for volcanic outgassing to revive super-Earth atmospheres. For M-dwarf planetary systems, super-Earths that have atmospheres close to the star likely were formed with abundant volatiles. Introduction.
The Solar System has three planets – Earth, Venus, and Mars – that have atmospheres derived from H -free solids (secondary atmospheres), and four giant planets whose atmospheres are derived from protoplanetary nebula gas (H -dominated primary 2 atmospheres) (Atreya et al. 1989, Marty et al. 2016). However, this clean separation in process and outcome is apparently unrepresentative of the known exoplanets. The two most common types of known exoplanet, the rocky super-Earth-sized planets (planet radius R pl < 1.6 R ⊕ , where “⊕” is the Earth symbol, planet density ρ pl > 4 g/cc; “super-Earths”) and the gas-shrouded sub-Neptunes ( R pl = 2 - 3 R ⊕ ), are divided by a valley in (planet radius)-(orbital period) space in which planets are less common (Fig. 1) (Fulton et al. 2017). The radius valley can be understood if, and only if, a substantial fraction of planets that are born with thick (>10 kbar) H -dominated primary atmospheres lose those atmospheres and shrink in radius to become rocky super-Earths (Van Eylen et al. 2018). This radius-shrinking process, which carves out the radius valley, may be the way that most super-Earths form. There is no evidence that Earth and Venus underwent this process, and primary atmospheres are not thought to have contributed to the origin of major volatile elements in Earth (Dauphas & Morbidelli 2014). Both Venus and Earth have secondary volatile envelopes (comprised of solid-derived volatiles, including H O and CO ) that are much less massive than the atmospheres of sub-Neptunes. Large surface reservoirs of H O, C species, and N species are essential to life on Earth and to Earth’s habitable climate (Bergin et al. 2015, Catling & Kasting 2017). Does forming with a thick primary atmosphere (sub-Neptune) but ending up as a rocky super-Earth favor the rocky planet ending up with a secondary atmosphere? Do primary atmospheres, in dying, shield high-molecular-weight species from loss during the early era of atmosphere-stripping impacts and intense stellar activity, allowing those constituents to later form a secondary atmosphere? Or does the light (H -dominated) and transient primary atmosphere drag away the higher-molecular weight species? Getting physical insight into the transition from primary to secondary atmospheres is particularly important for rocky exoplanets that are too hot for life. Hot rocky exoplanets are the highest signal/noise rocky targets for upcoming missions such as James Webb Space Telescope (JWST) (Kempton et al. 2018) and so will be the most useful for checking our understanding of this atmospheric transition process. Secondary atmospheres are central to the exoplanet exploration strategy (National Academy of Sciences 2018, Zahnle & Catling 2017). Previous work on the hypothesis that primary atmospheres played a role in forming secondary atmospheres includes that of Eucken in the 1940s (Turekian & Clark 1969), Urey (1951), Cameron and coworkers (e.g. Slattery et al. 1980), Sasaki (1990), and Ozima & Zahnle (1993). Secondary atmosphere formation on exoplanets in the absence of a primary atmosphere has been investigated by (e.g.) Elkins-Tanton & Seager (2008), Dorn et al. (2018), and references therein. Gas survival during planetary volume reduction.
During sub-Neptune to super-Earth conversion, we suppose the planet contains both nebular-derived H and also high- μ species (derived from the planet-forming solid materials) that could form a secondary atmosphere – if retained. Retention of volatiles is controlled by atmospheric loss and atmosphere-interior exchange (Fig. 2). For both sub- 3 Neptunes and super-Earths, silicates (both magma and solid rock) apparently make up most of the planet’s mass (e.g. Dai et al. 2019, Owen and Wu 2017). Atmosphere loss is more difficult for high-molecular-weight volatiles than for H because higher-molecular-weight volatiles such as H O are more easily shielded within the silicate interior and are also more resistant to escape. This can be understood as follows. First, a basic control on atmosphere loss is the ratio, λ , of gravitational binding energy to thermal energy: λ = ( μ / k T ua ) ( G M pl / ( R + z ) ) ≈ 2.3 μ (10 K / T ua )( M pl / 6 M ⊕ )/( ( R + z ) / 2 R ⊕ ) (1) where μ is molecular mass, k is Boltzmann’s constant, T ua is upper-atmosphere temperature, G is the gravitational constant, M pl is planet mass, R is the radius of the silicate planet, and z is atmosphere thickness. When λ ≲ 2, the upper atmosphere flows out to space at a rate potentially limited only by the energy available from upper-atmosphere absorption of light from the star (Watson et al. 1981). If the upper atmosphere absorbs 100 W/m of light from the star, the upper limit on loss is ~1000 bars/Myr. This hydrodynamic outflow ejects the tenuous upper atmosphere at high speed, but the dense lower atmosphere remains close to hydrostatic equilibrium. When λ > 10, the hydrodynamic outflow shuts down. For primary atmospheres, the mean molecular weight μ avg ≈ 2 Da (H ) while for secondary atmospheres, μ avg is ≳10× higher favoring secondary atmosphere retention. Moreover, many secondary-atmosphere constituents (e.g CO ) are much more effective coolants than H , so for secondary atmospheres T ua is lower (e.g. Kulikov et al. 2007, Johnstone et al. 2018). Atmospheric thickness z is also smaller for high- μ avg atmospheres due to their smaller scale height, raising λ . For all atmospheres, proximity to the star increases T ua (lowering λ ) and also increases the upper atmosphere outflow rate when the λ ≲ 2 condition is satisfied. Moreover, closer to the star impacts occur at higher velocities that are more erosive (Kegerreis et al. 2020). Thus we expect massive worlds far from the star to retain atmospheres and low-mass worlds closer to the star to lose them (Zahnle & Catling 2017). This expectation, while physically valid, offers little guidance as to whether or not super-Earths will have atmospheres. To go further, we need to consider the effect of evolving atmospheric composition on loss rate (Eqn. 1) and also track the shielding of volatiles within silicates (magma oceans and solid rock) (Fig. 2). The atmospheres of young sub-Neptunes are underlain by magma oceans (e.g. Vazan et al. 2018; Fig. S1). Equilibrium partitioning of a volatile i between the atmosphere (where it is vulnerable to escape) and the magma ocean (where it is shielded from escape) is given by c i = e i + p i s i m magma + other reservoirs = e i + [ ( e i / A mai ) g mai ( μ avg / μ i ) ] s i m magma + other reservoirs (2) where c i is the total inventory of the volatile, e i is the mass in the atmosphere, p i is the partial pressure of i at the magma-atmosphere interface, A mai (and g mai ) are the area of (and gravitational acceleration at) the magma-atmosphere interface, μ avg is the mean μ of the 4 atmosphere, s i is the solubility coefficient (mass fraction Pa -1 ) of the volatile in the magma, and m magma is the mass of liquid in the magma ocean. Solubility in magma is higher for many secondary-atmosphere constituents (e.g. CH , H S, H O) than for H . When the atmosphere contains both H and a more-soluble, high- μ constituent, the H takes the brunt of atmosphere loss processes that might otherwise remove the high- μ species. Differential solubility also protects the more-soluble species from impact shock (Genda & Abe 2005). Other reservoirs include volatiles stored in the crystal phase. This is an unimportant reservoir early on when the planet has a massive magma ocean, but becomes important later, because these volatiles can be released by volcanic degassing. We neglect volatiles that are stored in the iron core because they do not contribute to the observable atmosphere. Atmospheric loss typically drives magma ocean crystallization (Fig. S1) (except for planets that either orbit uncommonly close to their star, or undergo strong tidal heating, bare-rock planets do not have a global molten-rock layer). The greenhouse effect from the atmosphere keeps the magma ocean liquid for longer, which delays partitioning of the volatile into crystals. This partitioning is described by X xtl = D i X i (3) where X xtl is the concentration of i in the crystal phase, D i is a solid-melt distribution coefficient / partition coefficient, and X i is the concentration of i in the magma. Crystallization enriches the residual melt in volatiles (i.e. D i << 1; e.g. Elkins-Tanton 2008). One last effect can aid retention of high-molecular-weight volatiles. If outflow to space is slow, then the high- μ species can sink back toward the planet (diffusive separation driven by buoyancy) (e.g. Hunten et al. 1987, Hu et al. 2015; Supplementary Information 1g). Eqns. (1)-(3) open two paths from a primary to a secondary atmosphere. A secondary atmosphere can be exsolved as the magma ocean crystallizes. Alternatively, the atmosphere can be revived by volcanic degassing long after the mantle has almost entirely solidified (Fig. 2). We explore these paths below. Model of the transition from primary to secondary atmospheres.
We model the effect of the loss of a thick H -dominated atmosphere on the retention (or loss) of a hypothetical volatile, s . Species s has molecular mass 30 Da, which is within a factor of <1.5 of the molecular mass of all known major secondary-atmosphere constituents. In the uppermost atmosphere, species s is modeled to split into fragments of mass 15 Da for the purpose of calculating whether or not s is effectively entrained by escaping H. We assume that reactions between s and H do not affect the mass inventory of either species; in effect, s is chemically inert. We assume that deviations from thermochemical equilibrium driven by photochemistry in the uppermost atmosphere are reset to thermochemical equilibrium at the high temperature and pressure of the magma-atmosphere interface. While idealized, the model includes many effects that have not 5 previously been incorporated into a model of atmosphere loss, such as a realistic rock melting curve, differential solubility effects, e.t.c. (Supplementary Information). Planet equilibrium temperature T eq (in K), for zero planet albedo, is given by T eq = 278 ( F/F ⊕ ) (4) where F is insolation (W/m ) and the Sun’s insolation at Earth’s orbit, F ⊕ , is 1361 W/m . Upper-atmosphere temperatures T ua > T eq are essential for a Super-Earth atmosphere to flow out to space: this is possible because upper atmospheres efficiently absorb light at wavelength ≲ 100 nm, but do not readily re-emit this light. We adopt s H2 = 2 × 10 -12 Pa -1 (Hirschmann et al. 2012; Supplementary Information 1e), which for initial atmospheric pressure P atm,init = 50 kbar gives a total mass of H (atmosphere plus dissolved-in-magma) of 7 × 10 kg (2% of planet mass). We neglect He, so our model primary atmosphere is slightly more soluble in magma and has a slightly lower molecular weight than in reality. We assume the crystal-melt partition coefficients are D H2 = 0 (for simplicity) and D s = 0.02 (Supplementary Information 1f). Sub-Neptunes have a global shell of magma that freezes during conversion to a Super-Earth (Fig. S1). Planet thermal structure (below the photosphere, which is isothermal by assumption) is as follows (Fig. S2). The temperature at the top of the magma ocean, T mai , is given by the atmosphere adiabat ( T mai / T RCB ) = ( P atm / P RCB ) (γ-1)/ γ (5) where T RCB is the temperature at the radiative-convective boundary within the atmosphere (RCB), P RCB is the pressure at the RCB, and γ is the adiabatic index. We assume T RCB / T eq = 1, that the temperature jump at the magma-atmosphere interface is small, and that the magma ocean is isentropic. (If T RCB / T eq = 1.5 in a 1000 K orbit then the planet cooling timescale is only < 1 Kyr). This basic model absorbs the planet cooling rate and the radiative opacity of the atmosphere into a single parameter, P RCB (Supplementary Information 1c); for more sophisticated models, see e.g.
Marcq et al. (2017). As P atm decreases, T mai cools, and the magma crystallizes. Crystallization starts at great depth and the crystallization front sweeps slowly upward (e.g. Elkins-Tanton 2008). Volatiles enriched near the crystallization front due to the low solubility of volatiles in crystals will be stirred by fast magma currents (speed up to 10 m s -1 , Solomatov 2015) to the near-surface, where they form bubbles that pop and add gas to the atmosphere. Stirring and bubbling can degas the mantle down to at least ~100 GPa depths (e.g. Caracas et al. 2019; Solomotov 2015, and references therein). A smaller portion of s will go into the crystals. This portion is shielded within the rock and available for later volcanic outgassing. We do not include liquid volatiles (e.g. clouds) or fluid-fluid immiscibility. We emphasize results for T eq < 1150 K, cold enough for silicates to condense (Woitke et al. 2018). The model is intended for super-Earths too hot for life. Volcanic outgassing after the magma ocean has completely crystallized is guided by the results of Kite et al. (2009) (Fig. S3) (Supplementary Information 1h). 6 We approximate diffusive separation of H and heavy gases as zero during sub-Neptune to super-Earth conversion. During conversion, escape proceeds too quickly for the s to settle out (Hu et al. 2015) (Supplementary Information 1g). On worlds that are cool enough for life, diffusive separation is important (Supplementary Information 1g). On worlds that are cool enough for life, diffusive separation can allow high- μ species to be retained in the atmosphere while H escapes, including in cases where s is a reducing species that is dissociated into easily-escaping H plus a heavy atom (e.g. H O à τ > 1 Gyr; e.g., Gonnerman & Mukhopadhyay 2009, Saji et al. 2018). Because volcanic degassing of terrestrial planets is so slow/inefficient (e.g. Kaula 1999), we consider magma-ocean exsolution separately from volcanic outgassing. Specifically, we set re-release of s from the solid mantle to zero so long as the exsolved atmosphere is present. With these approximations, for a given insolation L , planet mass M pl , initial dose of s , and P atm at which that dose is applied, the output depends on how much atmosphere has been removed, but not on the speed (or process) of removal. In other words, the equations are time-independent. The small planet evolution sequence.
We first model atmospheric evolution during magma ocean crystallization (Figs. 3-4). We show results for 6 M ⊕ (≈1.6 R ⊕ ), corresponding to the largest (and therefore highest signal/noise) planets that commonly have densities consistent with loss of all H (Rogers 2015). We use a total mass of high- μ volatile ( M s ) = 3 × 10 kg. This corresponds to the near-surface C inventory of Earth and Venus, scaled up to 6 M ⊕ (and it is 2× the mass of Earth’s ocean). We drive the model by decreasing P atm . In Figure 3, as P H2 (blue dashed lines) falls, the magma ocean crystallizes, so that remaining volatiles go into the atmosphere (green), go into the solid mantle (maroon), or escape to space (black). We find that the main controls on the transition on exoplanets from primary to secondary atmospheres are F and the solubility ( s s ) of the high- μ atmosphere constituent. Fig. 3a shows a case with low s s (10 -11 Pa -1 ), for a planet close to its host star ( F/F ⊕ = 240; T eq = 1100 K; typical for Kepler super-Earths). The magma ocean stays fully liquid until the atmospheric mass is reduced by 90%. (Release of dissolved-in-magma H by bubbling is a negative feedback on atmospheric loss.) Crystallization begins at P atm = 5 kbar, and completes at ~200 bars. s is passively entrained (either as atoms or molecules) in the escaping H . We stop the run at P atm = 1 bar, so we do not track the removal of the last bit of the exsolved atmosphere. In this limit of small s s , only 6 × 10 kg × ( s s / 10 -11 Pa -1 ) × ( D s / 0.02 ) is shielded within the solid mantle. Most of the s is lost to space before crystallization begins. The outcome is a bare rock with a volatile-starved solid mantle, incapable of much volcanic outgassing. 7 Fig. 3b shows results for the same s s as in Fig. 3a and a cooler orbit ( F/F ⊕ = 3, T eq = 370 K, intermediate between Mercury and Venus in our Solar System). In the cooler orbit, some crystals are present initially. With no fully-liquid stage, crystallization can start to shield s within crystals as soon as atmospheric loss starts. The s available for later volcanic outgassing in the cooler orbit case is 30× greater. Raising s s to 10 -9 Pa -1 (equivalent to 1 wt% solubility for 100 bars partial pressure) favors secondary atmosphere occurrence (Figs. 3c-3d). More s is dissolved in the magma, and so more s partitions (during crystallization) into the rock, where it is shielded. Because s is much more soluble in the magma than H , very little s is initially in the atmosphere. Therefore, relatively little s is carried away to space during H removal and so more s is available for exsolution during the final crystallization of the magma ocean. This can create a high- μ avg atmosphere, which is easier to retain. Protection by differential solubility is enhanced by the hot orbit ( L/ L ⊕ = 240; Fig. 3c) as opposed to the cool orbit ( L/ L ⊕ = 3; Fig. 3d), because the magma ocean is long-lived for this case and so the volatiles are safely dissolved for longer. As a result, in the hot orbit case, enough high- μ species remain at final magma ocean crystallization to create an atmosphere with high- μ avg (the green line crosses the blue line in Fig. 3c). Such a high- μ avg atmosphere is easier to retain. The small amount of H that remains in the atmosphere at this point can be lost by diffusion-limited escape. Together, F and s s have strong effects on the chance of secondary atmosphere occurrence. Volatiles that are shielded within rock are available for late volcanic outgassing. For high s s the mass of shielded volatiles tends towards the product of the initial inventory of s and the solid-liquid partition coefficient D s (Fig. 4a). This is because almost all of the s is in the magma until the magma ocean has almost completely crystallized. The effect of F/F ⊕ on the mass of shielded volatiles is relatively modest for s s > 3 × 10 -11 Pa -1 . For F/F ⊕ ≤ 25 and low s s < ( μ s /μ H2 ) s H2 , the amount of the high-molecular-weight species that is shielded within rock is proportional to s s . For low s s , the mass of shielded volatiles decreases rapidly for hot orbits. In hot orbits the H wind can carry away s for a long time before the magma ocean cools enough for solidification (and shielding). When μ avg is high, loss to space is less likely (Eqn. 1 and surrounding discussion). Fig. 4b shows how atmospheric μ avg (after crystallization of the magma ocean is complete) depends on F and s s . For s s > 10 -11 Pa -1 , a greater fraction of the s s is stored in the magma than is H . As a result, the atmosphere becomes more s -rich as the magma ocean crystallizes. The enrichment is especially strong for hot orbits, because for hot orbits there is more s available to be exsolved: less s has escaped to space before crystallization completes. For s s < 10 -11 Pa -1 , most of the s is stored in the atmosphere, and so the atmospheric s mixing ratio is near-constant during atmosphere loss. (At F/F ⊕ < 3, diffusive separation favors high- μ avg atmospheres; Supplementary Information 1 g). Cooler orbits favor volcanic outgassing but hotter orbits permit exsolved high-molecular-weight atmospheres. To determine which effect is more important for the chance of seeing a secondary atmosphere on a super-Earth, we turn to a time-dependent model. Where are planets today on the small planet evolution sequence? F XUV plateaus at ~10 -3 × total F for planets around young stars, switching to a power-law decay at < 0.1 Gyr for Solar-mass stars ( F XUV = 3 × 10 -6 F at Earth today) (Supplementary Information 1a). The plateau of high F XUV / F is longer at red dwarf (M) stars (≥ 0.3 Gyr long for ≤ 0.5 M ☉ ) . Therefore, we expect that (for a given T eq ) planets orbiting M-stars will have lost more atmosphere (e.g. Lammer et al. 2007). F XUV drives atmospheric loss (rate d M atm /dt) (Supplementary Information 1b). Nebula-composition atmospheres do not cool efficiently, leading to high T ua that favors hydrodynamic escape (Murray-Clay et al. 2009). For nebular-composition atmospheres, we use: d M atm /d t = ε R pl ( R + z (t) ) F XUV / ( G M pl ) (6) where ε is an efficiency factor. For high- μ atmospheres we adopt the loss fluxes of a CO atmosphere (Tian 2009) as an example of a strong coolant. The model of Tian (2009) and Tian et al. (2009) includes dissociation of CO under high XUV levels, and the escaping material is atomic C and O. Tian (2009) finds low ε for CO atmospheres, and negligible hydrodynamic escape for F XUV < 0.6 W m -2 (= 150× the value on Earth today) (Fig. S4). At intermediate compositions we interpolate using a logistic function (Supplementary Information 1b). Fig. 5 shows atmosphere thickness versus time for a Solar-mass star and P atm,init = 50 kbar. High- μ avg atmospheres can only persist at a narrow range of F . (A similar pattern is seen for low-mass stars; Figure S6). Worlds far from the star receive a low XUV flux and stay as sub-Neptunes. Worlds in hotter orbits lose their atmospheres completely. They may undergo a rapid increase in μ avg (blue to yellow in Fig. 5a), but the resulting high- μ avg atmosphere is almost always short-lived. Why is this? According to the pure-CO model of Tian (2009), d M atm /d t for high molecular weight atmospheres on super-Earths has a threshold at ~150× F XUV / F XUV ,⊕ , below which loss is much slower. But in most cases the primary atmosphere is lost when F XUV is still very high so any exsolved secondary atmosphere is swiftly lost. For low s s , the atmosphere has a composition that is always H -dominated (by number), so exsolved high- μ species are lost in the H wind (Fig. 5b). We conclude that exsolved atmospheres are rare. This conclusion is robust because we use parameters that are favorable for an exsolved atmosphere. For example, we use a high solubility-in-magma for the high- μ species, but our high- μ escape-to-space parameterization is for a species (CO ) whose solubility-in-magma is low. Revival of secondary atmospheres by volcanic outgassing. release at Earth’s mid-ocean ridges (12±2 bars/Gyr; Tucker et al. 2018). The rate of outgassing is adjusted downward to account for loss of volatiles during sub-Neptune-to-super-Earth conversion, assuming 50% of worlds have s s = 10 -9 Pa -1 , and 50% of worlds have s s = 10 -11 Pa -1 . We neglect atmosphere re-uptake to form minerals (e.g., Sleep & Zahnle 2001), because our focus is on worlds that are too hot for aqueous weathering. This omission is conservative relative to our conclusion that volcanically-outgassed atmospheres on hot rocky exoplanets are uncommon. Fig. 6 (gray curve) outlines the region within which, in our model, volcanic outgassing will build up a secondary atmosphere. The line of revival sweeps towards the star over time, because the rate of volcanic degassing falls off more slowly with time than does the star’s XUV flux. Volcanic revival of the atmosphere is difficult for planets around M = 0.3 M ☉ stars, but easier for rocky planets around solar-mass stars. The results are sensitive to changes in s s , and to the choice of XUV flux models (Figs. S9-S11, S17). Fig. S17 also shows results for volcanism at a constant Earth-scaled outgassing rate. Fig. 6 also shows the atmosphere presence/absence line for rocky worlds that start with no H , and with all high-molecular-weight volatiles in the atmosphere (red curves). Such “intrinsically rocky” worlds retain residual secondary atmospheres over a wider range of conditions than do worlds that start as sub-Neptunes. The line of atmosphere loss for these residual atmospheres sweeps further away from the star with time. Discussion.
Our model results are sensitive to the rate of XUV-driven mass loss. The XUV flux of young Solar-mass stars varies between stars of similar age by a factor of 3-10 (e.g. Tu et al. 2015). Stars with low XUV flux are more likely to host planets with atmospheres (Fig. S11). Escape of N or H O is likely faster than the Tian (2009) loss rate estimate (which is for pure-CO atmospheres) that is used in our model (Zahnle et al. 2019, Johnstone et al. 2018). Improved knowledge of escape rate will require more escape rate data and XUV flux data (e.g. Bourrier et al. 2018, Ardila et al. 2018). Currently, XUV flux data as a function of star mass and star age are limited, and upcoming space missions such as SPARCS will gather more data (Ardila et al. 2018, Linsky et al. 2019).
Our model considers only thermal loss. However, solar wind erosion can remove atmospheres that are already thin (Dong et al. 2018). A possible example is Mars over the last ~4 Ga. Planetary dynamos can under some circumstances suppress solar-wind erosion (e.g. Gunnell et al. 2018). Mars had a dynamo prior to 4 Gya, and Earth (but not Venus) has one today. Solubility of gas in magma varies between species. Carbon-bearing volatile species are very insoluble in magma. H O’s solubility in basaltic magma at H O partial pressure 0.03 GPa is ~2 wt%; linearizing, this gives solubility s s = 7 × 10 -10 Pa -1 (e.g. Papale 1997). HCl is even more soluble in magma than H O (Fegley 2020), so a Cl-rich initial composition would have 10 a greater chance of forming an exsolved high- μ atmosphere, although we still predict it would be short-lived. The effects on atmosphere prevalence of 100-fold reduction in solubility are shown in Figs. S9-S11. Better constraints on solubility at T > 2000 K are desirable (e.g. Guillot & Sator 2011). Volcanism on super-Earths should wane over Gyr, according to models. We need more data to test these models. In models, the rate of decrease of volcanism depends on whether or not the planet has plate tectonics, on planet mass, and on mantle composition (e.g. Stevenson 2003, Unterborn et al. 2015, Dorn et al. 2018, Byrne 2019; Supplementary Information). Two effects, of uncertain relative importance, are ignored in our volcanism model. The first is a buffering effect: if the mantle is volatile-rich, then magma is produced more easily (Médard & Grove 2006), but if the mantle is volatile-poor, then the melt rate is reduced (reducing the rate of volatile loss). The second is a potential fine-tuning issue: if volcanic degassing is very rapid, then volatiles will be released from the protective custody of the mantle before atmosphere loss process have lost their bite.
Overall, our choices of M pl , D s , and s s tend to favor the existence of volcanically outgassed atmospheres. Even with these choices atmospheres are usually not stable for planets at T eq > 500 K around M-stars. So, our conclusions are broadly unfavorable for atmospheres on rocky exoplanets at T eq > 500 K around M-stars. However, this conclusion is moderated by the possibility (discussed next) that worlds with abundant H O exist close to the star.
Observational tests.
Atmospheres on rocky exoplanets can now be detected (e.g. Demory et al. 2016, Kreidberg et al. 2019). Theory predicts that retaining an atmosphere should be harder on planets orbiting low-mass stars and the present study extends that prediction to super-Earths that form as sub-Neptunes. A test of this theory would have major implications for habitable zone planets. If this prediction fails, that would suggest that M-star rocky exoplanets formed more volatile-rich than rocky exoplanets orbiting Sun-like stars (Tian & Ida 2015). It is possible that some planets form with more volatiles than can be removed by loss processes (Bitsch et al. 2019). Some models predict formation of hot Super-Earths with 1 wt % - 30 wt% H2O, either by accretion of volatile-rich objects (for example, extrasolar analogs to CI/CM chondrites), or by planet migration (e.g. Raymond et al. 2018, Bitsch et al. 2019). Such volatile-rich worlds are hard to distinguish from bare-rock planets using current data. XUV-driven loss can remove at most a few wt% of an M = 6 M ⊕ planet’s mass for T eq < 1000 K. Our model implies that a planet in the “no atmosphere” zone of Fig. 6 with a JWST detection of H O-dominated atmosphere is more volatile-rich than Venus and Earth. Possible volatile-rich worlds include planets that have radii ~0.2 R ⊕ greater than expected for Earth-composition (Fig S12) (Turbet et al. 2019). Fig. 6 enables the following testable predictions. Since atmospheres close to the star can only persist if the initial H O inventory is high, N /NH should be diluted to very low mixing ratio for such atmospheres. If a super-Earth-sized planet has an atmosphere, then planets at greater semimajor axis in the same system should also have an atmosphere. Starting out 11 as a sub-Neptune is unfavorable for atmosphere persistence, so systems where the planets formed intrinsically rocky should have statistically more atmospheres. Multi-planet systems enable strong tests because uncertainty in time-integrated stellar flux cancels out. Conclusions.
A large fraction of rocky exoplanets on close-in orbits (closer to their star than the Habitable Zone) were born with thick (>10 kbar) H -dominated (primary) atmospheres but have since lost their H . In our model these T eq > 400 K exoplanets almost never transition smoothly to worlds with high-molecular weight atmospheres (Fig. 7). Instead, the high-molecular-weight species are usually carried away to space by the H wind. Volcanic outgassing is an alternative source for a high-molecular-weight atmosphere. Revival of a bare-rock planet by volcanic outgassing gets easier with time, because atmosphere loss slows down rapidly with time, but atmosphere supply by volcanism decays slowly. But volcanic outgassing is also enfeebled by early loss of biocritical volatiles via the H wind. Overall, for a given initial dose of high-molecular-weight species, atmospheres are less likely on hot rocky exoplanets that were born with thick H -dominated atmospheres. Many uncertainties remain, the most important of which is the initial planet volatile content. Within our model, for planets that orbit Solar-mass stars, super-Earth atmospheres are possible at insolations much higher than for planet Mercury. For planets that orbit ~0.3 M ☉ stars, secondary atmospheres at much higher insolation than planet Mercury in our solar system are unlikely unless the planet formed H -poor, or includes a major (~1 wt%) contribution of solids from beyond the water ice line (waterworld). Data availability:
All of the code for this paper, together with instructions to reproduce each of the figures and supplementary figures, can be obtained via the Open Science Framework at https://osf.io/t9h68 or by emailing the lead author.
Acknowledgements:
We thank two reviewers for accurate and useful reviews. We thank B. Fegley, Jr., L. Schaefer, L. Rogers, E. Ford, and J. Bean (discussions). Grants: NASA Exoplanets Research Program (
NNX16AB44G ). Figures.
Figure 1.
The exoplanet abundance histogram (gray band, for orbital periods < 100 days, corrected for detection biases; Fulton & Petigura 2018). Two classes of small exoplanet are seen: volatile-rich sub-Neptunes and rocky super-Earth-sized exoplanets. Sub-Neptune to super-Earth conversion is implied by the data and may be the way that most super-Earths form. 13
Figure 2.
Processes (italics) and reservoirs (upright font) in our model. Atmosphere-interior exchange is central to the transition from primary to secondary atmospheres. Timescales are approximate. 14 (a) s s = 10 -11 Pa -1 , F/F ⊕ = 240 (b) s s = 10 -11 Pa -1 , F/F ⊕ = 3 (c) s s = 10 -9 Pa -1 , F/F ⊕ = 240 (d) s s = 10 -9 Pa -1 , F/F ⊕ = 3 Figure 3.
The small planet evolution sequence, for 6 M ⊕ . Black triangles correspond to total initial inventory of high- μ species. (a) Solubility of high- μ species in magma ( s s ) = 10 -11 Pa -1 , insolation normalized to Earth’s insolation ( F/F ⊕ ) = 240 (planet equilibrium temperature 1095 K). (b) s s = 10 -11 Pa -1 , F/F ⊕ = 3 (planet equilibrium temperature = 360 K). (c) s s = 10 -9 Pa -1 , F/F ⊕ = 240. (d) s s = 10 -9 Pa -1 , F/F ⊕ = 3. (a) (b) Fig. 4. (a)
How the mass of volatiles shielded within the solid mantle depends on F and s s . (b) How atmospheric mean molecular weight (after crystallization of the magma ocean is complete) depends on F and s s . 16 (a) (b) Figure 5.
Time-dependent results for planets orbiting a Solar-mass star. (Fig. S6 shows the results for a 0.3 M ☉ star.) (a) Atmospheric pressure vs. time for s s = 10 -9 Pa -1 , (from top to bottom) F/F ⊕ = {49, 283, 347, 422, 720}, corresponding to 17 planet equilibrium temperature ( T eq ) = {735, 1140, 1200, 1275, 1440} K. (b) As (a), but for s s = 10 -11 Pa -1 . Fig. 6.
Secondary atmosphere presence/absence model output for 6 M ⊕ (higher planet mass favors atmosphere retention). The red lines and gray lines show atmosphere presence/absence contours for two different scenarios. The red lines show atmosphere retention thresholds after 3.0 Gyr for the case where all volatiles are in the atmosphere initially and there is no primary atmosphere; the 16 th and 84 th percentiles are shown, for varying XUV flux (by ±0.4 dex, 1σ; Loyd et al. 2019) relative to the baseline model following the results of Jackson et al. (2012) and Guinan et al. (2016) (Supplementary Information 1a). The red lines move away from the star over time (red arrows). The gray lines show the 16 th and 84 th percentiles for exhibiting an atmosphere after 3.0 Gyr for the case where volcanic outgassing rebuilds the atmosphere from a bare-rock state. The solid gray lines are for stagnant-lid tectonics and the dashed gray lines are for plate tectonics. The lines of atmospheric revival sweep towards the star over time (gray arrows) because the rate of volcanic degassing falls off more slowly with time than does the star’s XUV flux. In each case the atmosphere/no-atmosphere threshold is 1 bar. The black symbols show known planets that may be tested for atmospheres using JWST (Koll et al. 2019). For any individual planet, star-specific XUV-flux estimates, star age, and the planet’s mass, should 18 be combined to make a more accurate estimate than is possible using this overview diagram. Figs. S8 and S16-17 show further details.
Fig. 7.
Graphical summary. The left column corresponds to atmosphere revival (gray lines in Fig. 6) and the right column corresponds to residual atmospheres (red lines in Fig. 6). For worlds in T eq > 400 K orbits that start as sub-Neptunes, formation of a high- μ atmosphere is unlikely, unless the planet starts with abundant high- μ volatiles. References.
Andrault, D., Bolfan-Casanova, N., Nigro, G. L., et al. 2011, Solidus and liquidus profiles of chondritic mantle: Implication for melting of the Earth across its history, Earth Planet. Sci. Lett., 304, 251
Ardila, D., et al. 2018, The Star-Planet Activity Research CubeSat (SPARCS), Proc. AIAA/USU Conference on Small Satellites, Instruments / Science I, SSC18-WKIV-02.
Atreya S., Pollack, J., & Matthews, M. (editors), Origin & Evolution of Planetary and Satellite Atmospheres, U. Arizona Press, 1989
Bergin, E. A., Blake, G. A., Ciesla, F., et al. 2015, Tracing the ingredients for a habitable earth from interstellar space through planet formation, PNAS, 112, 8965.
Bitsch, B.; Raymond, S. N.; Izidoro, A., 2019, Rocky super-Earths or waterworlds: the interplay of planet migration, pebble accretion, and disc evolution, A&A, 624, A109.
Bourrier, V., et al., 2018, Hubble PanCET: an extended upper atmosphere of neutral hydrogen around the warm Neptune GJ 3470b, A&A, 620, A147.
Byrne, P.K., 2019, A comparison of inner Solar System volcanism, Nature Astronomy, doi:10.1038/s41550-019-0944-3.
Caracas, R., et al. 2019, Melt–crystal density crossover in a deep magma ocean, EPSL, 516, 202-211.
Catling, D. C., & Kasting, J. F. 2017, Atmospheric Evolution on Inhabited and Lifeless Worlds (Cambridge: Cambridge Univ. Press)
Dai, F., et al. 2019, Homogeneous Analysis of Hot Earths: Masses, Sizes, and Compositions, ApJ 883, 79
Dauphas, N., Morbidelli, A. 2014. Geochemical and planetary dynamical views on the origin of Earth’s atmosphere and oceans. in Treatise on Geochemistry (2nd Edition), Heinrich D. Holland & Karl K. Turekian, Eds, Elsevier, Oxford, 2014, 1-35.
Demory, B.-O., Gillon, M., Madhusudhan, N., & Queloz, D. 2016, Variability in the super-Earth 55 Cnc e, MNRAS, 455, p.2018.
Dong, C., Jin, M., Lingam, M., et al. 2018, Atmospheric escape from the TRAPPIST-1 planets and implications for habitability, PNAS, 115, 260
Dorn, C., et al., 2018, Outgassing on stagnant-lid super-Earths, A&A, 614, A18
Elkins-Tanton, L., & Seager, S., 2008, Ranges of Atmospheric Mass and Composition of Super-Earth Exoplanets, ApJ, 685, 1237-1246.
Elkins-Tanton 2008, Linked magma ocean solidification and atmospheric growth for Earth and Mars, EPSL 271, 181-191
Elkins-Tanton, L. T. 2011, Formation of early water oceans on rocky planets, Ap&SS, 332, 359
Fegley, B., et al., 2020, Volatile element chemistry during accretion of the Earth, Geochemistry: Chemie der Erde, in press.
Fulton, B. J., Petigura, E. A., Howard, A. W., et al. 2017, The California-Kepler survey. III. A gap in the radius distribution of small planets, AJ, 154, 109
Fulton, B. J., & Petigura, E. A. 2018, The California-Kepler survey. VII. Precise planet radii leveraging Gaia DR2 reveal the Stellar mass dependence of the planet radius gap, AJ, 156, 264
Genda H. & Y. Abe 2005, Enhanced atmospheric loss on protoplanets at the giant impact phase in the presence of oceans, Nature 433, 842
Gonnerman H. & S. Mukhopadhyay, 2009, Preserving noble gases in a convecting mantle, Nature 459(7246):560-3.
Guillot B. & N. Sator 2011, Carbon dioxide in silicate melts: A molecular dynamics simulation study, GCA 75, 1829-1857.
Gupta, A. & H. Schlichting, 2019, Sculpting the valley in the radius distribution of small exoplanets as a by-product of planet formation: the core-powered mass-loss mechanism, Monthly Notices Roy. Astro. Soc. 2019, 487, 24–33.
Hu, R., et al. 2015, Helium atmospheres on warm Neptune-and sub-Neptune-sized exoplanets and applications to GJ 436b, ApJ, 807:8.
Huang, C.X., et al. 2018, Expected Yields of Planet discoveries from the TESS primary and extended missions, arXiv:1807.11129
Hunten et al. 1987, Mass fractionation in hydrodynamic escape, Icarus 69, 532-549.
Johnstone et al. 2018, Upper atmospheres of terrestrial planets: Carbon dioxide cooling and the Earth's thermospheric evolution, A&A 617 A107
Kegerreis, J.A., et al., 2020, Atmospheric Erosion by Giant Impacts onto Terrestrial Planets, arxiv:2002.02977
Kempton, E. M.-R., et al., 2018, A framework for prioritizing the TESS planetary candidates most amenable to atmospheric characterization, Proc. Astron. Soc. Pacific, 130, 993.
Kite, E. S., Manga, M., & Gaidos, E. 2009, Geodynamics and rate of volcanism on massive Earth-like planets, ApJ, 700, 1732
Kreidberg, L., et al. 2019, Absence of a thick atmosphere on the terrestrial exoplanet LHS 3844b, Nature 573, 87–90.
Kulikov, Y.N., et al. 2007, A comparative study of the influence of the active young Sun on the early atmospheres of Earth, Venus, and Mars, Space Science Reviews, 129, 207–243
Lammer, H., et al. 2007, Coronal Mass Ejection (CME) Activity of Low Mass M Stars as An Important Factor for The Habitability of Terrestrial Exoplanets. II. CME-Induced Ion Pick Up of Earth-like Exoplanets in Close-In Habitable Zones, Astrobiology 7, 185-207.
Linsky, J., 2019, Host Stars and their Effects on Exoplanet Atmospheres, Springer.
Marcq, E. et al. 2017, Thermal radiation of magma ocean planets using a 1 ‐ D radiative ‐ convective model of H2O ‐ CO2 atmosphere, JGR-Planets Volume 122, Pages 1539-1553
Marty, B., et al EPSL 2016, Origins of volatile elements (H, C, N, noble gases) on Earth and Mars in light of recent results from the ROSETTA cometary mission, 441, 91-102.
Médard, E., & T. Grove 2006, Early hydrous melting and degassing of the Martian interior JGR-Planets 111, E11003.
Murray-Clay, R., et al. 2009, Atmospheric escape from hot Jupiters, ApJ 693, 23-42.
Owen, J. E., & Wu, Y. 2017, The evaporation valley in the Kepler planets, ApJ, 847, 29
Ozima, M., & K. Zahnle 1993, Mantle degassing and atmospheric evolution: Noble gas view, Geochemical Journal, 27, 185-200.
Papale, P. 1997, Modeling of the solubility of a one-component H2O or CO2 fluid in silicate liquids, Contrib. Mineral. Petrol., 126, 237.
Raymond, S.N.; Boulet, T.; Izidoro, A.; Esteves, L.; Bitsch, B., 2018, Migration-driven diversity of super-Earth compositions, MNRAS Letters, 479, L81-L85
Rogers, L., 2015, Most 1.6 Earth-radius planets are not rocky, ApJ 801:41 (13 pp).
Saji, N.S., et al. 2018, Hadean geodynamics inferred from time-varying 142Nd/144Nd in the early Earth rock record, Geochemical Perspectives Letters, 7, 43-48.
Sasaki, S. 1990, in The Primary Solar-type Atmosphere Surrounding the Accreting Earth: H O-induced High Surface Temperature, ed. H. Newsom & J. H. Jones (Oxford: Oxford Univ. Press), 195
Slattery, W. L.; Decampli, W. M.; Cameron, A. G. W., 1980, Protoplanetary Core Formation by Rain-Out of Minerals, The Moon and the Planets, Volume 23, Issue 4, pp.381-390
Sleep, N., and K. Zahnle, Carbon dioxide cycling and implications for climate on ancient Earth, J. Geophys. Res. Planets 2001, 106, 1373-1399. nd edition, Chapter 9.04, Magma Oceans and Primordial Mantle Differentiation, vol. 9, p. 81-104. Stevenson, D.J., 2003, Styles of mantle convection and their influence on planetary evolution, Comptes Rendus Geoscience, 335, 99-111.
Tian, F. 2009, Thermal escape from super Earth atmospheres in the habitable zones of M stars, ApJL 703:905-909
Tian, F., et al. 2009, Thermal escape of carbon from the early Martian atmosphere, GRL 36, L02205, doi:10.1029/2008GL036513.
Tian, F. 2015, History of water loss and atmospheric O2 buildup on rocky exoplanets near M dwarfs, EPSL, 432, p. 126-132.
Tian, F. & S. Ida 2015, Water contents of Earth-mass planets around M dwarfs, Nature Geoscience 8, 177–180
Tu, L., et al., 2015, The extreme ultraviolet and X-ray Sun in Time: High-energy evolutionary tracks of a solar-like star, A&A, 577, L3
Tucker, J.M., Mukhopadhyay, S., & Gonnerman, H.M., 2018, Reconstructing mantle carbon and noble gas contents from degassed mid-ocean ridge basalts, EPSL, 496, 108-119.
Turbet, M., et al., 2019, Revised mass-radius relationships for water-rich terrestrial planets beyond the runaway greenhouse limit, arXiv:1911.08878.
Turekian, K. K.; Clark, S. P., 1969, Inhomogeneous accumulation of the earth from the primitive solar nebula, Earth and Planetary Science Letters, 6, 346-348.
Unterborn C. T., Johnson J. A., Panero W. R., 2015, Thorium abundances in solar twins and analogs: implications for the habitability of extrasolar planetary systems, ApJ, 806, 139
Urey, H., 1951, The origin and development of the earth and other terrestrial planets, GCA 1, 209-277 van Summeren, J., Conrad, C. P., & Gaidos, E. 2011, Mantle convection, plate tectonics, and volcanism on hot exo-Earths, ApJL, 736, L15
Van Eylen, V., et al. 2018, An asteroseismic view of the radius valley: stripped cores, not born rocky, MNRAS 479, 4786-4795
Vazan, A., et al. 2018, Contribution of the core to the thermal evolution of sub-Neptunes, ApJ 869:163.
Watson, A., et al. 1981, The dynamics of a rapidly escaping atmosphere: applications to the evolution of Earth and Venus, Icarus 48, 150-166
Woitke, P, Helling, Ch., Hunter, G.H., et al., 2018, Equilibrium chemistry down to 100 K. Impact of silicates and phyllosilicates on the carbon to oxygen ratio, Astronomy & Astrophysics, Volume 614, id.A1, 28 pp.
Zahnle, K., & D. Catling 2017, The cosmic shoreline: The evidence that escape determines which planets have atmospheres, and what this may mean for Proxima Centauri b, ApJ 843:122
Zahnle, K., et al. 2019, Strange messenger: A new history of hydrogen on Earth, as told by Xenon, Geochim. Cosmichim. Acta 244, 56-85.
Supplementary Information.
Atmosphere evolution involves many processes (Catling & Kasting 2017). This is (to our knowledge) the first paper on the transition on exoplanets from primary to secondary atmospheres. Therefore our approach is to keep the treatment of each process simple. 22 We consider star ages from 5 Myr, which is a typical time for nebular disk dispersal, up to 8 Gyr. We consider star masses from 0.1-1 M ☉ . We frequently state results for 0.3 M ☉ , because these are the most common host stars for worlds detected by the TESS mission (e.g. Sullivan et al. 2015, Huang et al. 2018). We emphasize planet mass 6 M ⊕ (bare-rock radius ≈ 1.6 R ⊕ ), because these are the largest (and therefore highest signal/noise) worlds that commonly have densities consistent with loss of all H (Rogers 2015). The silicate mass fraction is held fixed at 2/3, with the balance consisting of the Fe/Ni metal core. This silicate mass fraction is based on Solar System data. Our results are only weakly sensitive to reasonable variations in the silicate mass fraction. The total mass of H is small compared to planet mass (Lopez & Fortney 2014). The non-H volatile mass fraction is also assumed to be small, consistent with data (e.g. van Eylen et al. 2018). Candidate drivers for the atmosphere loss that converts sub-Neptunes into super-Earths are photoevaporation, impact erosion, and core luminosity (Inamdar & Schlichting 2016, Zahnle & Catling & 2017, Burger et al. 2018, Ginzburg et al. 2018, Owen & Wu 2017, Biersteker & Schlichting 2019). We focus on photoevaporation in this study. Photoevaporation has been directly observed for sub-Neptune exoplanets (Bourrier et al. 2018), fits almost all of the data (e.g. Van Eylen et al. 2018), and is relatively well-understood (e.g. Murray-Clay et al. 2009). Other mechanisms are important at least for a small number of planetary systems (e.g. Owen & Estrada 2019, Loyd et al. 2019). The rate of photoevaporation is set by the X-ray / Extreme Ultraviolet (XUV) flux, F XUV (e.g. Murray-Clay et al. 2009, but see also Howe et al. 2019). XUV is generated by hot regions near the star’s surface. These regions are heated by magnetic fields whose importance declines as the star sheds angular momentum, through the stellar wind, over time. Thus the star’s XUV luminosity F XUV declines over time. Early in a star’s history, the star’s XUV luminosity L XUV reaches a plateau value (referred to as “saturation”) of between 10 -4 × and 10 -3 × the star’s total, bolometric luminosity ( L) . Estimating L XUV involves combining direct measurements for stars of similar spectral type and age, and interpolation over wavelength regions where star UV is absorbed by the interstellar medium. Data on L XUV is synthesized by Linsky (2019). Lopez & Rice (2018) and Luger & Barnes (2015) propose power-law decays for Sun-like stars. Tu et al. (2015) quantify the variability of X-ray emission for stars of Solar mass. Guinan et al. (2016) compile L XUV data for low-mass stars. We used two approaches to interpolating in L XUV / L as a function of star mass and time. In one approach, we used Figure 5 from Selsis et al. (2007) for the X-ray flux, obtaining the XUV flux from the X-ray flux using the top row of Table 4 in King et al. (2018). Selsis et al. (2007) assume a constant L XUV / L of 10 -3.2 during saturation. Selsis et al. (2007) do not attempt to quantify the variation of X-ray flux between stars of the same age and the same stellar mass. In the other approach, we used Jackson et al.’s (2012) fits for L XUV / L for M ≥ 23 0.5 M ☉ , patching to Guinan et al.’s (2016) fits for M < 0.5 M ☉ . We assumed that Guinan et al.’s (2016) fits applied for M = 0.35 M ☉ and extrapolated to other M-stars. We used standard models for L and star radius (Baraffe et al. 2015) (Fig. S13). The L XUV / L output from the two approaches is shown in Fig. S14. The zones of secondary atmosphere loss and revival from these approaches are shown in Fig. S9-S11. For both approaches, we also calculated how results would vary between stars of different L XUV / L , assuming a 1-standard-error scatter in L XUV / L of ±(0.4-0.5) dex (following Loyd et al. 2019) (Fig. 6, S16). We use an upper wavelength cutoff of 91.2 nm (13.6 eV), corresponding to the lowest-energy photon that can ionize H. However, longer-wavelength light can still contribute to atmospheric escape, especially as μ avg rises. Future work using more sophisticated escape models (e.g. Wang & Dai 2018) might investigate the effect of switching on or off various UV spectral bands at lower energies than 13.6 eV. For pure-H atmospheres, we set d M atm /d t = ε R pl ( R + z ) F XUV / ( G M pl ) (S1) (Eqn. 6 from main text) where the heating efficiency factor, ε = 0.15. ε = 0.15 is a common choice in atmospheric evolution models. Models indicate that ε = 0.1-0.2 (e.g. Shematovich et al. 2014). As the planet’s atmosphere loses mass and cools, it shrinks (Lopez & Fortney 2014). This acts as a negative feedback on loss in part because smaller planets intercept fewer XUV photons. To relate the atmosphere mass M atm to the planet radius ( R + z ), we used the transit-radius tables of Lopez & Fortney (2014). The homopause radius could be bigger than the transit radius by typically 5-25% (Malsky & Rogers 2020). For our pure-high-molecular-mass atmosphere endmember loss rate, we use the loss rate in molecules/sec/cm calculated for a pure-CO atmosphere on a super-Earth by Tian (2009). This is an endmember of a high -μ avg atmosphere that can self-cool efficiently. This loss rate is multiplied by μ avg , and by the area of the surface-atmosphere interface, to get the loss rate in kg/s. Our intent is to bracket the likely down-shift in atmospheric loss as the atmosphere evolves to high -μ avg . We extrapolate beyond the XUV-flux limits shown in Figure 6 of Tian (2009) linearly in loss rate for very large XUV fluxes, and linearly in the log of loss rate for very small XUV fluxes. High -μ avg atmospheres of different composition may shed mass faster (or slower) than the pure-CO case (Johnstone 2020, Johnstone et al. 2018, Zahnle et al. 2019). For atmospheres of intermediate composition, we interpolate between the energy-limited formula (d M atm /dt ∝ F XUV ) to the exponential cutoff in F XUV for high- μ avg atmospheres (d M atm /dt ≈ 0 below a critical F XUV ) found by Tian (2009). Constraints for hydrodynamic loss rates of atmospheres of intermediate composition are limited (Kulikov et al. 2007, Johnston et al. 2018). In the absence of better constraints, we use a logistic curve in log (escape rate) and log (mixing ratio) : 24 X s,atm = ( μ avg – μ p ) / ( μ s – μ p ); (S2) Y = 1 / (1 + exp(- k log X s,atm – log k ) ) (S3) with k = 8 and k = 0.2. We think that these parameter choices understate the X CO2 needed to suppress escape (i.e. favor the survival of high- μ atmospheres), which is conservative in terms of our conclusion that atmosphere survival is unlikely. Next we set Z = log (d M atm /d t ) p - Y log (d M atm /d t ) p - log (d M atm /d t ) s (S4) (d M atm /d t ) mixed = 10 Z (S5) where (d M atm /d t ) p is the energy-limited escape rate with ε = 0.15, (d M atm /d t ) s is from Tian (2009), and (d M atm /d t ) mixed is the loss rate for impure atmospheres. Planet thermal structure in our model is defined by the temperature at the atmosphere’s radiative-convective boundary ( T RCB ) and the temperature at the magma-atmosphere interface, T mai . T RCB and T mai are related by the thickness of the convecting (adiabatic) envelope, which (in pressure units) is ( P atm – P RCB ), and by the adiabatic index γ (held constant at 4/3). We do not consider changes in γ with changes in atmospheric composition. Planets form hot and cool over time. Cooling involves energy transport through the atmosphere by radiation and/or convection. If, at a given optically thick level of the atmosphere, all cooling is accomplished by radiation, then the radiative temperature gradient is d T /dr = – (3 κ ρ / 15 σT ) ( L int /4 π r ) (S6) (e.g. Rogers et al. 2011), where κ is the local Rosseland-mean opacity, ρ is the local atmospheric density, T is local temperature, σ is the Stefan constant, and ( L int /4 π r ) is the internal luminosity corresponding to the cooling of the planet (where r is the distance from the center of the planet). If the radiative d T /d r from Eqn. (S6) exceeds the convective adiabatic temperature gradient, then convection takes over and d T /d r is reduced to a value slightly greater than the adiabatic lapse rate. Thus, as the planet cools, and L int decreases, the radiative zone expands ( P RCB moves to greater depths) (e.g. Bodenheimer et al. 2018, Vazan et al. 2018). We set P RCB = 10 bars and neglect P RCB change with time, which means our magma oceans are slightly longer-lived than they would be in reality. On the other hand, we use the 1 Gyr T eq value throughout the calculation, when for cool stars T eq will be higher earlier in the planet’s life. This simplification tends to shorten model magma ocean lifetime relative to reality. We assume full redistribution of energy absorbed on the dayside to the entire (4π) area of the planet. Our model planets are also isothermal with altitude above the radiative- 25 convective boundary. In effect, we neglect the temperature difference between the planet photosphere (level at which the optical depth is approximately unity) and the radiative-convective boundary. At the top of the atmosphere we assume zero albedo. This is reasonable for cloud-free atmospheres of planets orbiting cool stars, although for Solar-mass stars the albedo will be ~0.2, considering only CO and H O (g) opacity (Pluriel et al. 2019). Our focus is interior to the habitable zone and for thick atmospheres. For such worlds atmospheric retention would result in magma oceans that last Gyr (Hamano et al. 2013, Vazan et al. 2018, Fig. S1). Thus, interior to the habitable zone, atmosphere loss is the rate-limiter for magma ocean crystallization (Hamano et al. 2013, Hamano et al. 2015, Lebrun et al. 2013). In principle our model could lead to an equilibrium situation where a few-hundred-bar high- µ atmosphere sustains a magma ocean indefinitely. However such outcomes in our model require fine-tuning. During sub-Neptune-to-super-Earth conversion, planet cooling is driven by atmospheric removal. As the overburden pressure is removed, the remaining atmosphere cools adiabatically. This causes a temperature gradient between the magma and the atmosphere, and heat flow from the magma into the atmosphere soon warms the atmosphere because the atmosphere heat capacity is small compared to the magma ocean heat capacity. This raises T RCB , and because T RCB / T eq > 1 cools the planet quickly, this cools the planet. All of these processes occur continuously in our simplified model: we assume the magma temperature tracks the loss of atmosphere. In effect we assume that convection swiftly reduces super-adiabatic temperature gradients and that the boundary layer temperature contrast at the top of the liquid magma (viscosity < 10 -1 Pa s) is small. More detailed models of sub-Neptune thermal evolution include Rogers et al. (2011), Howe & Burrows (2015), Vazan et al. (2018), and Bodenheimer et al. (2018). For cooler orbits in our model, some crystals are present initially. This incomplete-magma-ocean initial condition is appropriate if the planet is assembled by impacts of objects that are Mars-sized or smaller (e.g. planetesimals, pebbles, or planetary embryos), and the planet cools down to the adiabat between impacts.
Magma ocean crystallization shields dissolved volatiles from loss. The magma ocean crystallizes (losing mass to crystals) as the planet cools. The magma ocean mass at a given T mai depends on conditions deep within the silicate layer (higher T and much higher P than at the magma-atmosphere interface). The silicate solidus ( T vs. P for 0% melt fraction) and liquidus ( T vs. P for 100% melt fraction) are curved in T-P space. As a result, for a steady decline in T mai , the magma mass will decrease first quickly and then slowly. To track this effect, we used the solidus, the liquidus, and the adiabat of Andrault et al. (2011) (Fig. S15). We assume that the melt fraction increases linearly with increasing temperature between the solidus and the liquidus, and neglect chemical fractionation during crystallization. We also neglect the reduction in the solidus T due to volatile enrichment as crystallization proceeds. More detailed models of magma ocean crystallization include 26 Nikolau et al. (2019) and Katyal et al. (2019). Our model predicts whole mantle melting for T mai ≥ 3000 K with negligible melt below T mai ~ 1750 K. Melting curves differ between experimenters and differ depending on assumed mantle composition (e.g. Andrault et al. 2017, Miyazaki & Korenaga 2019). To relate pressure to depth, we use the pressure-density curve from Dziewonski & Anderson (1981). To get the radius at the magma-atmosphere interface we assume ( R pl / R ⊕ ) = ( M pl / M ⊕ ) (Valencia et al. 2006). We integrate downward from the magma-atmosphere interface. Starting at P atm and T mai , we follow the adiabat until we reach the solidus. The melt mass is the mass of melt (if any) between the magma-atmosphere interface and the depth of the magma solidus. If this melt mass exceeds 2/3 of planet mass, the melt mass is set equal to 2/3 of planet mass. The effect of the weight of the atmosphere on the location of the solidus is included in our calculations, but turns out to be unimportant. We assume full re-equilibration between magma and atmosphere as the atmosphere is removed (efficient degassing). Degassing is efficient if the cooling magma ocean is fully convective (Ikoma et al. 2018), which can be understood as follows. Whether or not a fully liquid magma ocean undergoing whole-magma-ocean convection will degas as the H envelope is removed and the liquid cools depends on (a) the number of times each magma parcel cycles through the upper thermal boundary layer during cooling, and (b) the ratio of the degassed boundary layer thickness to the thermal boundary layer thickness. Suppose the degassed boundary layer thickness is equal to the maximum depth at which bubbles can form (this assumes swift bubble ascent, and that supersaturation is unimportant). For sub-Neptune to super-Earth conversion, this is (up to 10s kbar)/( ρ magma g ) = up to 50 km. The thermal boundary layer thickness for a fully liquid magma ocean will be ≪ O dissolves in magma in proportion to the square root of partial pressure at low pH O, becoming more linear at high pH O (e.g. Stolper 1982, Matsui & Abe 1986). We adopt s H2 = 2 × 10 -12 Pa -1 . The basis for this is as follows. Hirschmann et al. (2012) report results from laboratory experiments on basaltic melts. Extrapolating to melts of peridotitic composition, they state that “at 1 GPa in the presence of pure H , the molecular H concentration [in the melt] is 0.19 wt%.” Fegley et al. (2020, their Table 5) compile solubilities of various volatiles in molten silicates. We ignore non-ideal solubility (due to non-ideal fugacity of the gas in the atmosphere), because it is not very important at the relatively low P atm we consider (Kite et al. 2019). We also ignore joint-solubility effects because relevant data are not available. Few measurements have been made of the temperature dependence of volatile solubility in magma at T > 2000 K (e.g. Fegley & Schaefer 2014, Guillot & Sator 2011) and we neglect T dependence of volatile solubility. 27 Solid-melt distribution coefficients, D , for water partitioning into nominally anhydrous mantle minerals vary from 10 -4 to 10 -1 (Table S1 in Elkins-Tanton 2008), but for the most common mantle minerals are 0.02 or less. Here, “nominally anhydrous” means minerals that do not have H in their chemical formula. We use D = 0.02, which is at the high end of experimentally observed range (corresponding to water portioning into pyroxene). We neglect saturation limits. C distribution coefficients are much lower, from 2 × 10 -4 to 7 × 10 -3 according to Elkins-Tanton (2008) and with low (1-5 ppmw) saturation limits. Key initial conditions are that the quantity of the high molecular weight volatile initially in the crystal phase is zero ( s xtl = 0), and e s = C s /(1 + ( μ p / μ s )( g mai / A mai ) m magma s s ) (S7) which assumes that the atmosphere is dominantly μ p at the zeroth timestep. We restrict ourselves to F/F ⊕ < 10 . For F/F ⊕ > 10 , silicate vapor-pressure equilibrium atmospheres develop (e.g. Kite et al. 2016). Our model loops through the following steps: (Step 1) Lose some atmosphere. à (Step 2) Post-escape outgassing. à (Step 3) Crystallize (sequestering some volatiles, but increasing the concentration of volatiles in the magma). à (Step 4) Post-crystallization outgassing à Loop back to Step 1. Step 1. Lose atmosphere. A down-step in atmospheric pressure is prescribed. No fractionation is permitted between H and the high- μ species during this escape (see section 1g). The H and the high- μ species are lost in proportion to their bulk abundance in the atmosphere. Step 2. Post-escape outgassing . The magma and atmosphere now re-equilibrate through outgassing of both H and s until the partial pressures of H and of s are in equilibrium with the concentrations of H and s dissolved in the magma. The partial pressure of H depends on the mass of H in the atmosphere and on the mean molecular weight of the atmosphere, which is affected by the amount of s in the atmosphere. Similarly the partial pressure of s depends on the mass of s in the atmosphere and on the mean molecular weight of the atmosphere, which is affected by the amount of p in the atmosphere. From Eqn. (2) in the main text, we obtain the coupled equations C s = e s + e s (( e s + e p )/( e s / μ s + e p / μ p )) k / μ s (S8) C p = e p + e p (( e s + e p )/( e s / μ s + e p / μ p )) k / μ p (S9) where C s (kg) is the total on-planet mass of s , C p (kg) is the total on-planet mass of H ( p ), e s is the kg of s in the atmosphere, e p is the kg of H in the atmosphere, the second term on the right-hand-side of (S8) is the mass (kg) of s dissolved in the magma, the second term on 28 the right-hand-side of (S9) is the mass (kg) of H dissolved in the magma, μ s is the molecular weight of s, μ p is the molecular weight of H , and k = ( g mai / A mai ) m magma s s (S10) k = ( g mai / A mai ) m magma s p (S11) where g mai and A mai are (respectively) the gravity at, and area of, the magma-atmosphere interface for ( R mai / R ⊕) = ( M pl / M ⊕ ) , m magma is calculated as described in section 1d above, s s is the solubility of s in the magma, and s p is the solubility of H in the magma. The simultaneous equations (S10)-(S11) have two solutions, and the physical solution is e p = (C p k + C s k + (k (C p k μ s - k (C p k μ p + 2 C p k μ p μ s + C p μ s + 2 C p C s k k μ p μ s - 2 C p C s k μ p + 4 C p C s k μ p μ s + 4 C p C s k μ p μ s - 2 C p C s k μ s + 2 C p C s μ p μ s + C s k μ s + 2 C s k μ p μ s + C s μ p ) - C s k μ p + 2 C s k μ s + C p k μ p + C s k k μ s ))/(2 (k μ p - k μ s + k μ p - k k μ s )) + (k k (C p k μ s - k (C p k μ p + 2 C p k μ p μ s + C p μ s + 2 C p C s k k μ p μ s - 2 C p C s k μ p + 4 C p C s k μ p μ s + 4 C p C s k μ p μ s - 2 C p C s k μ s + 2 C p C s μ p μ s + C s k μ s + 2 C s k μ p μ s + C s μ p ) - C s k μ p + 2 C s k μ s + C p k μ p + C s k k μ s ))/(2 (k μ p - k μ s + k μ p - k k μ s )))/(k + k k ) (S12) e s = -(C p k μ s - k (C p k μ p + 2 C p k μ p μ s + C p μ s + 2 C p C s k k μ p μ p - 2 C p C s k μ p + 4 C p C s k μ p μ s + 4 C p C s k μ p μ s - 2 C p C s k μ s + 2 C p C s μ p μ s + C s k μ s + 2 C s k μ p μ s + C s μ p ) - C s k μ p + 2 C s k μ s + C p k μ p + C s k k μ s )/(2 (k μ p - k μ s + k μ p - k k μ s )) (S13) We recalculate the other variables by inserting the new e s and e p into Eqns. (S8) and (S9). The coupling of Eqn. (S8) and Eqn. (S9) via μ avg has two important effects that disfavor secondary atmospheres on exoplanets. During the loss of the primary atmosphere, the μ avg effect increases the amount of the high- μ constituent in the atmosphere by a factor of 15 (= μ s /μ p ), increasing the amount of the high- μ constituent that is lost to space. Moreover, as μ avg rises from 2 to 30, because the partial pressure of the high- μ constituent scales with the number ratio of s in the atmosphere, f s = e s / ( e p + e s )( μ avg /μ p ), the number of moles of s in the atmosphere for a given saturation vapor pressure of s decreases by a factor of 15 ( = μ s /μ p ). This negative feedback reduces the number of worlds that smoothly transition from a primary to a secondary atmosphere. Step 3. Crystallize. The magma mass is forced to a lower value, consistent with the lower T mai associated with the atmospheric pressure that was reduced in Step 1. The freshly-produced crystals gain an s concentration (in mass fraction) of D s × X s , where X s is the mass fraction of the volatile in the magma. The crystals are not remixed (by assumption) on the timescale of magma-ocean crystallization. s within crystals is relatively safe from release into the unsafe surface environment during sub-Neptune-to-super-Earth conversion, 29 because solid-state creep speeds within the solid mantle are 10 -9 – 10 -8 m s -1 , versus up to 10 m s -1 in the magma. As a result, the concentration of the crystals is independent of the previous partitioning of the volatiles into the solid rock. Partitioning of H into the crystals is not considered. Step 4. Post-crystallization outgassing. Because X s < 1, the magma ocean is enriched in s after the crystallization step. Therefore we re-equilibrate the magma and the atmosphere taking account of the now-higher concentration of s ( X s ) in the magma. The processes involved in Steps 1-4 all occur on timescales short compared to the magma ocean lifetime. To map the rocky planet evolution sequence onto time we combine the environmental drivers from (1a), the loss rates from (1b), and the small planet evolution sequence from (1f). T RCB is fixed in the crystallization calculation, but T eq evolves with star age. We deal with this by setting T RCB throughout the calculation equal to T eq at a fixed age of 1 Gyr. Thus, planets orbiting cool stars in our model have shorter-lived magma oceans than in reality. Diffusive separation of H from high- µ species is important at lower XUV flux (e.g. Hu et al. 2015, Wordsworth et al. 2018, Saito & Kuramoto 2018, Odert et al. 2018, Malsky & Rogers 2020). However, diffusive separation is washed out at the high XUV fluxes required if photoevaporation is to cause sub-Neptune-to-super-Earth conversion (Tian 2015, Catling & Kasting 2017). To demonstrate this, we first define the (atomic) H flux, F , on the assumption that the overwhelming majority of the upper atmosphere is comprised of H: F = ε F XUV ( R + z ) / (4 G M pl m ) (S14) in units of molecules s -1 m -2 . Here, ε = 0.15 is the efficiency of conversion of XUV energy into atmospheric escape, R is solid-planet radius, z is atmosphere thickness (which we set to 0.4 R to correspond to a small sub-Neptune), G is the gravitational constant, M pl is planet mass, and m is the mass (kg) of the hydrogen atom. From this we obtain the crossover mass m c (Tian 2015) m c = m + k T F /( b g X ) (S15) where T is the temperature of the heterosphere, and b = 4.8 × 10 T m -1 s -1 is the bindary diffusion coefficient, corresponding to neutral H and neutral O, used by Tian (2015). Including the effect of ionization of H would decrease the effective value of b (Hu et al. 2015) and so strengthen our conclusion that diffusive separation is not important during sub-Neptune to super-Earth conversion for Kepler ’s super-Earths. We set T = 10 K, the same value used by Hu et al. (2015) and Malsky & Rogers (2020). 30 The flux, F , of the high- µ species is zero when m c < m . When m c > m the flux F is given by F = F ( X ( m c – m )) / ( X ( m c – m )) (S16) Substituting in values for 6 M ⊕ and atomic mass 15 Da we find F XUV = 0.7 W/m at crossover, 1.3 W/m for 50% reduction in escaping flux due to weight differences, and 6.1 W/m for escaping flux of high-molecular weight species at 90% of the flux if there was no fractionation at all. All three thresholds are shown as dashed vertical lines on Fig. S5 and Fig. S6. In the context of sub-Neptune to super-Earth conversion, most of the atmosphere mass is removed at higher F XUV . Most exoplanets either stay as sub-Neptunes or lose their atmosphere too quickly for the XUV flux to drop to levels that permit >10% differences in loss rate due to fractionation. This can be understood as follows. If the upper atmosphere absorbs 1 W/m of light from the star, the upper limit on loss rate for a 6 M ⊕ world is ~10 bars/Myr, falling to 1.5 bars/Myr at ε = 0.15. Since sub-Neptune-to-super-Earth conversion takes ~300 Myr, and ~5 × 10 bars of H need to be removed (including H that is initially dissolved in the magma), F XUV must be ~100 W/m . Using the saturation threshold of Selsis et al. (2007), ( F XUV / F = 10 -3.2 ) this corresponds to 160 kW/m , F/F ⊕ = 120, for distance from the star = 0.1 AU for a Sun-like star. This is the typical distance of a Kepler super-Earth from its host star (the abundance of 1.6 R ⊕ ~ 6 M ⊕ super-Earths drops off rapidly at larger separations (Fulton & Petigura 2018). At such high XUV fluxes, fractionation causes negligible differences in the loss rate relative to the no-fractionation case. Therefore, our approximation of no fractionation during sub-Neptune to super-Earth conversion is valid. However, fractionation by diffusion is important can affect the composition of late-stage volcanically outgassed atmospheres. For example, the volcanic H outgassing flux is calculated to be less than the diffusively-limited H escape rate unless H stays at a very low level (e.g. Catling et al. 2001, Batalha et al. 2016, Zahnle et al. 2019). This effect maintains volcanically outgassed atmospheres at high μ av , even if outgassing of H by serpentinization (Sleep et al. 2004) and/or volcanic outgassing is large. Fractionation by diffusion is also much more important for Habitable Zone super-Earth-sized planets. For F bol = 1361 W/m (the Solar flux at 1 AU today), even at saturation F XUV is 0.9 W/m , which is only marginally capable of lifting atomic mass = 15 Da species out of the atmosphere. This favors the occurrence of secondary atmospheres on Habitable Zone exoplanets. Finally, exposure to moderate F XUV for Gyr may convert a H -dominated sub-Neptune atmosphere into a He-enriched sub-Neptune, provided the planet does not start with too much H (Hu et al. 2015, Malsky & Rogers 2020). No exoplanets have yet been confirmed to be volcanically active. We use a simple model of the geodynamics and rate of volcanism on super-Earth-sized planets (Kite et al. 2009). The model includes parameterized mantle convection, thermal evolution, and volcanism for both plate tectonics and stagnant lid modes. The model is tuned to reproduce Earth’s 31 present-day rate of volcanism in planet tectonics mode. Kite et al. (2009) consider three different melting models. In the present study, we use the results from the widely-used melting model of Katz et al. (2003). Many more sophisticated models of rocky exoplanet thermal evolution, geodynamics, and volcanic activity exist (e.g. Dorn et al. 2018, and references therein). The procedure for estimating the rate of volcanism starts from the results of Kite et al. (2009) (Fig. S3). This is an Earth-tuned parameterized mantle convection and rate of volcanism model that predicts the rate of volcanism versus time for either stagnant lid mode and plate tectonics mode. The rate of CO release at Earth’s mid-ocean ridges (12±2 bars/Gyr; Tucker et al. 2018) is multiplied by the computed rate of volcanism, in units of Earth’s present-day rate, to obtain the predicted rate of CO release on the rocky exoplanet. This is adjusted downward in proportion to any early loss of high-molecular-weight volatiles to space. For example, if a model exoplanet has lost 70% of its high-molecular-weight volatiles to space, then the rate of volcanic outgassing is reduced to 30% of the rate predicted by comparing the volumetric volcanic flux from the Earth-tuned code to that of the Earth. Fig. 6 compares the broad zone of worlds with atmosphere for planets that form without thick H atmospheres to the much more restricted zone of volcanically revived atmospheres for worlds that form with thick H atmospheres. This calculation assumes s s = 10 -9 Pa -1 (for which ~90% of high- μ volatiles are loss to space during the loss of the H -dominated atmosphere). Atmospheric ingassing is neglected because worlds well inside the habitable zone are too hot for weathering reactions to sequester atmophiles in the crustal rocks. Volatile release by volcanic outgassing is neglected, by construction, on worlds that start with all their volatiles in the atmosphere. Forming with a secondary atmosphere is advantageous for subsequently exhibiting a secondary atmosphere relative to the total-loss-and-subsequent-revival process. The stagnant-lid and plate-tectonics predictions diverge greatly after ~4 Gyr, when volcanism shuts down on the stagnant-lid model planet (Fig. S8). The plate-tectonics and stagnant lid predictions otherwise appear similar in this Earth-scaled model (Kite et al. 2009). However, this overall similarity of plate and stagnant-lid predictions is specific to the Earth-tuned model of Kite et al. (2009): at least one other model predicts much stronger suppression of volcanism on stagnant-lid super-Earth-sized rocky planets (e.g. Dorn et al. 2018). Tidal locking, by itself, has little effect on volcanism (e.g. van Summeren et al. 2011). So long as we do not mechanistically understand volcanism versus time for Solar System worlds (Byrne et al. 2019), why Earth’s mantle took so long to mix, nor why Earth’s mantle is not yet completely outgassed, we think that it is appropriate to use basic models to predict rate of volcanism versus time for rocky exoplanets. Robust constraints (for non-tidally-heated planets) are decline of mantle temperature by 50-200K per Gyr, with a concomitant decline in the potential for melting (Stevenson 2003). We assume the rate of degassing is proportional to the rate of magmatism. This is a simplification. Even low-degree partial melts can effectively extract almost all the volatiles from the full mass of silicate mantle that sweeps through the partial melt zone. This is 32 because volatiles partition strongly into the melt during partial melting (i.e., D i ≪ outgassing rate: a factor drawn randomly from a lognormal distribution centered on 0 and with a standard deviation of 2 bars/Gyr is added to the nominal value of 12 bars/Gyr. (3) Solubility of the high-molecular-weight species in the magma: alternated between values of 10 -11 Pa -1 and 10 -9 Pa -1 . We used 1 bar atmospheric pressure as the threshold for atmosphere presence/absence. The model results are more sensitive to changes in the XUV flux than they are to changes in the rate of volcanic outgassing. This is due to the exponential cutoff in the XUV-flux-driven loss rate of CO (Tian 2009). The inventory of CO even on today’s Earth is uncertain. For example, the review of Lee et al. (2019) reports a range of estimates for the C stored in the non-sedimentary rocks of Earth’s continental crust, from 4.2 × 10 Gton C to 2.6 × 10 Gton C. This increases the 33 likelihood that some hot rocky exoplanets have C inventories large enough to retain an observable atmosphere.
No exoplanets have yet been confirmed to contain >1 wt% H O. Modeling papers tracking the formation, migration, and (for habitable-zone worlds) climate evolution of worlds with >1 wt% H O include Kite & Ford (2018), Bitsch et al. (2019), and references therein.
Another way to increase volatile mass is to ballast the hydrogen with oxygen obtained from the reaction FeO magma + H à Fe + H O (g) ; i.e., endogenic water which does not require planet migration (Kite et al. 2020). Code availability. Everything used to make this paper can be obtained for unrestricted further use by emailing the lead author. Supplementary References.
Andrault, D., Bolfan-Casanova, N., Bouhifd, M. A., et al. 2017, Toward a coherent model for the melting behavior of the deep Earth's mantle, Phys. Earth & Planet Interiors, 265, 67
Baraffe, I., et al. 2015, New evolutionary models for pre-main sequence and main sequence low-mass stars down to the hydrogen-burning limit, Astronomy & Astrophysics, Volume 577, id.A42.
Biersteker, J. B., & Schlichting, H. E. 2019, Atmospheric mass-loss due to giant impacts: the importance of the thermal component for hydrogen–helium envelopes, MNRAS, 485, 4454
Bodenheimer, P., Stevenson, D. J., Lissauer, J. J., & D’Angelo, G. 2018, New Formation Models for the Kepler-36 System, ApJ, 868, 138
Catling, D.C., K. J. Zahnle, and C. P. McKay, Biogenic methane, hydrogen escape, and the irreversible oxidation of early Earth, Science, 293, 839-843, 2001.
Dziewonski, A. M., & Anderson, D. L. 1981, Preliminary reference Earth model, Phys. Earth Planet. Int., 25, 297
Fegley, B., & L. Schaefer 2014 Chemistry of the Earth's Earliest Atmosphere. Treatise on Geochemistry, 2nd edition, chapter 6.3.
Ginzburg, S., et al., Super-Earths, in Formation, Evolution, and Dynamics of Young Solar Systems, Astrophys. & Space Sci. Library, v.445, edited by M. Pessah & O. Gressel, Springer.
Guinan et al. 2016, Living with a red dwarf: rotation and X-ray and ultraviolet properties of the halo population Kapteyn's star, ApJ, 821, Issue 2, article id. 81.
Hamano, K. et al. 2015, Lifetime and spectral evolution of a magma ocean with a steam atmosphere: its detectability by future direct imagingApJ, 806, Issue 2, article id. 216.
Howe, A. R., & Burrows, A. 2015, Evolutionary models of super-Earths and mini-Neptunes incorporating cooling and mass loss, ApJ, 808, 150
Howe, A., F. Adams and M. Meyer 2019, Survival of Primordial Planetary Atmospheres: Photodissociation Driven Mass Loss, arXiv:1912.08820
Huang, C. X., Burt, J., Vanderburg, A., et al. 2018, TESS Discovery of a Transiting Super-Earth in the pi Mensae System, ApJL, 868, L39
Ikoma, M., Elkins-Tanton, L., Hamano, K., Suckale, J. 2018. Water partitioning in planetary embryos and protoplanets with magma oceans, Space Sci. Rev. 214, 76.
Inamdar, N. K., Schlichting, H. E. 2016. Stealing the gas: giant impacts and the large diversity in exoplanet densities, ApJ 817, L13.
Johnstone, C. et al. 2019, Extreme hydrodynamic losses of Earth-like atmospheres in the habitable zones of very active stars, Astronomy & Astrophysics, Volume 624, id.L10.
Johnstone, C. 2020, Hydrodynamic Escape of Water Vapor Atmospheres near Very Active Stars, Astrophysical Journal, Volume 890, Issue 1, id.79.
Katyal, N., et al. 2019, Evolution and Spectral Response of a Steam Atmosphere for Early Earth with a Coupled Climate–Interior Model, ApJ 875:31
Katz, R. F., Spiegelman, M., & Langmuir, C. H. 2003, A new parameterization of hydrous mantle melting, Geochem. Geophys. Geosyst. 4, 1073
King, G.W., et al., 2018, The XUV environments of exoplanets from Jupiter-size to super-Earth MNRAS 478, 1193-1208.
Kite, E. S., Fegley, B., Jr., Schaefer, L., & Gaidos, E. 2016, Atmosphere-interior exchange on hot, rocky exoplanets, ApJ, 828, 80.
Kite, E.S., & Ford, E.B., 2018, Habitability of exoplanet waterworlds, ApJ 864, 75.
Kite, E. S., Fegley, B., Schaefer, L., & Ford, E. B., 2019, Superabundance of exoplanet sub-Neptunes explained by fugacity crisis, ApJL, 87, L33.
Kite, E. S., Fegley, B., Schaefer, L., & Ford, E. B., 2020, Atmosphere origins on exoplanet sub-Neptunes, ApJ, 891, 111.
Kulikov et al. 2007, A comparative study of the influence of the active young Sun on the early atmospheres of Earth, Venus, and Mars, Space Science Reviews, 129, 207-243.
Lebrun, T., et al. 2013, Thermal evolution of an early magma ocean in interaction with the atmosphere, JGR-Planets 118, 1155-176.
Lee, C., Jiang, H., Dasgupta, R., & Torres, M. (2019). A Framework for Understanding Whole-Earth Carbon Cycling. In B. Orcutt, I. Daniel, & R. Dasgupta (Eds.), Deep Carbon: Past to Present (pp. 313-357). Cambridge: Cambridge University Press.
Lopez, E., & Fortney, J., 2014, Understanding the mass-radius relation for sub-Neptunes: Radius as a proxy for composition, ApJ 721:1.
Lopez, E. D., & Rice, K. 2018, How formation time-scales affect the period dependence of the transition between rocky super-Earths and gaseous sub-Neptunes and implications for η ⊕ , MNRAS, 479, 5303-5311. Loyd, R.O.P., et al., 2019, Current Population Statistics Do Not Favor Photoevaporation over Core-powered Mass Loss as the Dominant Cause of the Exoplanet Radius Gap, arXiv:1912.12305.
Luger, R. & R. Barnes 2015, Extreme water loss and abiotic O2 buildup on planets throughout the habitable zones of M dwarfs, Astrobiology, 15(2), 119-143.
Malsky, I., & L. Rogers 2020, Coupled Thermal and Compositional Evolution of Photo Evaporating Planet Envelopes, arXiv:2002.06466.
Mansfield, M. et al. 2019, Identifying Atmospheres on Rocky Exoplanets through Inferred High Albedo, ApJ , 886:141.
Matsui, T. & Y. Abe, 1986, Evolution of an impact-induced atmosphere and magma ocean on the accreting Earth, Nature 319, 303–305
Miyazaki, Y., & Korenaga, J. 2019, On the timescale of magma ocean solidification and its chemical consequences: 2. Compositional differentiation under crystal accumulation and matrix compaction, J. Geophys. Res. (Solid Earth), 124, 3399
Odert, P., et al., 2018, Escape and fractionation of volatiles and noble gases from Mars-sized planetary embryos and growing protoplanets, Icarus 307, 327-346.
Olson, P. L., & Sharp, Z. D. 2019, Nebular atmosphere to magma ocean: A model for volatile capture during Earth accretion, Phys. Earth Planet. Int., 294, 106294.
Owen, J., & B. C. Estrada 2019, Testing exoplanet evaporation with multitransiting systems, arXiv:1912.01609.
Owen, J. E., & Wu, Y. 2017, The evaporation valley in the Kepler planets, ApJ, 847, 29.
Otegi, J.F., et al. 2019, Revisited mass-radius relations for exoplanets below 120 M ⊕ , A&A, 643, A43. Pluriel, W., et al. 2019, Modeling the albedo of Earth-like magma ocean planets with H2O-CO2 atmospheres, Icarus 317, 583-590.
Rogers, L. A., Bodenheimer, P., Lissauer, J. J., Seager, S. 2011, Formation and structure of low-density exo-Neptunes, ApJ, 738, 59.
Saito, H., & K. Kuramoto, 2018, Formation of a hybrid-type proto-atmosphere on Mars accreting in the solar nebula, MNRAS, 475, 1274-1287.
Selsis, F., et al., 2007, Habitable planets around the star Gliese 581?, 476, 1373-1387.
Shematovich, V. I., Ionov, D. E., & Lammer, H. 2014, Heating efficiency in hydrogen-dominated upper atmospheres, A&A, 571, A94.
Sleep, N., et al., 2004, H2-rich fluids from serpentinization: geochemical and biotic implications, PNAS 101, 12818-12823.
Solomon, S.C. 1978, On volcanism and thermal tectonics on one ‐ plate planets, GRL, 5, 461-464. Stolper, E. 1982, The speciation of water in silicate melts, GeoCoA, 46, 2609.
Tu, L., et al., 2015, The extreme ultraviolet and X-ray Sun in Time: High-energy evolutionary tracks of a solar-like star, A&A 577, L3.
Valencia, D., O’Connell, R. J., & Sasselov, D. 2006, Internal structure of massive terrestrial planets, Icarus, 181, 545.
Van Eylen, V., Agentoft, C., Lundkvist, M. S., et al. 2018, An asteroseismic view of the radius valley: stripped cores, not born rocky, MNRAS, 479, 4786.
Vazan, A., Ormel, C. W., Noack, L., Dominik, C. 2018. Contribution of the core to the thermal evolution of sub-Neptunes, ApJ, 869, 163.
Wang, L., & F. Dai, 2018, Evaporation of Low-mass Planet Atmospheres: Multidimensional Hydrodynamics with Consistent Thermochemistry, ApJ, 860:175.
Wordsworth, R.D., et al. 2018, Redox evolution via gravitational differentiation on low-mass planets: Implications for abiotic oxygen, water loss, and habitability, AJ 155:195.
Supplementary Figures. Fig. S1.
How magma ocean mass increases with atmospheric thickness. Output from a toy model of sub-Neptune thermal structure (Kite et al. 2020). Dashed lines correspond to magma ocean mass, labeled in Earth-masses of magma, for a volatile-free planet mass of 6 M ⊕ . Colored lines correspond to temperatures at the magma-atmosphere interface of 1500 K (maroon), 2000 K (red), 3000 K (orange), and 4000 K (yellow). Magma ocean masses in excess of 2 Earth masses are not plotted because this corresponds to a magma pressure range that is not well explored by experiments. Fig. S2.
Planet thermal structure (reproduced from Kite et al., 2020). Temperature versus pressure plot, showing adiabats within the atmosphere and the magma (black line). T rheo (white line) corresponds to the temperature of the rheological transition (~40% melt fraction) for rock (dashed white line at levels where no rock is present). RCB = radiative-convective boundary. T mai = temperature at the magma-atmosphere interface. For real planets the equilibrium temperature ( T eq ), the effective temperature ( T eff ), and the temperature at the RCB ( T RCB ), are related by T eq < T eff ≲ T RCB . In this paper we make the approximation that these three temperatures are equal (Supplementary Information 1c). 37
Fig. S3. (Reproduced from Kite et al. 2009). Rate of volcanism per unit mass on massive Earth-like planets undergoing stagnant lid convection, normalized to calculated rate on a plate-tectonic Earth, for Katz et al. (2003) melting model. Dark gray shaded regions correspond to mantle temperatures associated with very intense volcanism, too high for a reliable crustal thickness calculation. Contours are at 0, 0.5, 1, 2, 5, and 10 times Earth’s present-day rate of volcanism (~4 × 10 −19 s −1 , equivalent to 24 km yr −1 , in plate-tectonics mode). For example, a 6 M ⊕ planet on the “2” contour erupts 6 × 2 × 24 ≈ 300 km yr −1 . Fig. S4.
Atmospheric loss rate for a pure CO atmosphere, combining results for super-Earths (Tian 2009) and for Mars (Tian et al. 2009, scaled to 1 AU). The vertical dashed line corresponds to 150× Earth’s present-day XUV flux. 38 Fig. S5.
The same evolutionary tracks as in Fig. 5, shown in F XUV – μ avg coordinates. F/F ⊕ increases from lowest track ( F/F ⊕ = 49, T eq = 735 K) to uppermost track ( F/F ⊕ = 720, T eq = 1440 K). The dashed lines highlight (from top to bottom) the threshold of F XUV below which all of the high-molecular-mass species is retained by the planet; the F XUV corresponding to a 50% reduction in the no-fractionation loss rate of the high-molecular-mass species; and the F XUV corresponding to escape of the high-molecular-mass species at 90% of the rate at which no loss would occur (Supplementary Information 1g). Right of the rightmost dashed line, fractionation protects the constituents of the secondary atmosphere, and left of the leftmost dashed line, fractionation is much less important. Calculations are done assuming an atomic wind of H entraining a 1% mixing ratio of atoms of mass 15 Da. 39 (a) (b)
Fig. S6. (a)
As Fig. 5 (top panel), but for a planet orbiting a star of 0.3 Solar masses. Time-dependent results. Atmospheric pressure vs. time for planets for (from top to bottom)
F/F ⊕ = {49, 283, 347, 422, 720}, corresponding to planet equilibrium temperature ( T eq ) = {735, 1140, 1200, 1275, 1440} K. (b) As Fig. S5 but for a planet orbiting a star of 0.3 Solar masses. The same evolutionary tracks as in the top panel, shown in F XUV – μ avg coordinates. F/F ⊕ increases from lowest track ( F/F ⊕ = 49, T eq = 735 K) to uppermost track ( F/F ⊕ = 720, T eq = 1440 K). The dashed lines highlight (from top to bottom) the threshold of F XUV below which all of the high-molecular-mass species is retained by the planet; the F XUV corresponding to a 50% reduction in the no-fractionation loss rate of the high-molecular-mass species; and the F XUV corresponding to escape of the high-molecular-mass species at 90% of the rate at which no loss would occur (Supplementary Information 1g). Right of the rightmost dashed line, fractionation protects the constituents of the secondary atmosphere, and left of the leftmost dashed line, fractionation is much less important. Fractionation calculations are done assuming an atomic wind of H entraining a 1% mixing ratio of atoms of mass 15 Da. a) b) c) d) Fig. S7.
Details of planet evolution for the tracks shown in Fig. 5a. Solar-mass star, 6 M ⊕ . Line thickness corresponds to insolation, with the thickest lines corresponding to the greatest insolation. Results are shown for F/F ⊕ = 49 (blue), F/F ⊕ = 283 (red), F/F ⊕ = 347 (yellow), F/F ⊕ = 422 (purple), and F/F ⊕ = 720 (green), corresponding to T eq = {735, 1140, 1200, 1275, 1440} K, respectively. (a) Planet radius versus time. (b)
Planet equilibrium temperature versus time. (c)
Planet mass loss rate versus time (solid lines), and the mass loss rate that the planet would have if the atmosphere was composed entirely of CO based on the calculations of Tian (2009) (dashed lines). The planets corresponding to the yellow, purple, and green lines evolve to a high -μ atmosphere, but only in the yellow case is this atmosphere long-lived. (d) The fraction of high -μ shielded within solid rock. This is zero for the two worlds that remain as sub-Neptunes (because no solid rock forms), and very similar for the three worlds that transition to super-Earths (full magma crystallization). 41 a) b) c) d) Fig. S8.
Atmospheric thickness versus time for the Selsis et al. (2007) XUV model. (a-b) volcanic revival, s s = 10 -9 Pa -1 , Kite et al. (2009) volcanism model. There is a minor trend in plate-tectonics mode for planets around Solar-mass stars at T eq < T eq . This is due to the smaller initial volume of the ( s -protecting) magma ocean at low T eq . (c-d) Start-with-all-volatiles-in-the-atmosphere, residual atmosphere case. 42
Fig. S9.
Deterministic calculation for atmosphere presence/absence after 3.0 Gyr. Deterministic calculation with a fixed XUV flux for each star mass, combining the fits of Jackson et al. (2012) and Guinan et al. (2016). The red dash-dot line corresponds to the line of vanished atmospheres for planets that have all volatiles in the atmosphere initially. This line moves away from the star over time. The thick gray solid and dashed lines correspond to the line of atmospheric revival by volcanic outgassing for planets that lose all atmosphere during transition from a sub-Neptune to a super-Earth. These lines of atmospheric revival sweep towards the star over time because the rate of volcanic degassing falls off more slowly with time than does the star’s XUV flux. This chart assumes initial volatile supply is independent of F and similar to Earth and Venus. Increasing volatile supply will move thresholds to higher F . The thin gray lines reduce the solubility of the high-molecular-weight volatile in the magma from 10 -9 Pa -1 to 10 -11 Pa -1 . Selected exoplanets overplotted. Note that because most planets are too small to have Earth-like composition and 6 M ⊕ , this plot is optimistic for atmospheric survival (higher planet mass favors atmosphere retention). 43 Fig. S10.
Deterministic calculation for atmosphere presence/absence after 3.0 Gyr. Deterministic calculation with a fixed XUV flux for each star mass, following Selsis et al. (2007). The red dash-dot line corresponds to the line of vanished atmospheres for planets that have all volatiles in the atmosphere initially. This line moves away from the star over time. The thick gray solid and dashed lines correspond to the line of atmospheric revival by volcanic outgassing for planets that lose all atmosphere during transition from a sub-Neptune to a super-Earth. These lines of atmospheric revival sweep towards the star over time because the rate of volcanic degassing falls off more slowly with time than does the star’s XUV flux. This chart assumes initial volatile supply is independent of F and similar to Earth and Venus. Increasing volatile supply will move thresholds to higher F . The thin gray lines reduce the solubility of the high-molecular-weight volatile in the magma from 10 -9 Pa -1 to 10 -11 Pa -1 . Selected exoplanets overplotted. Note that because most planets are too small to have Earth-like composition and 6 M ⊕ , this plot is optimistic for atmospheric survival (higher planet mass favors atmosphere retention). 44 Fig. S11.
Small planet atmosphere presence/absence diagram after 3.0 Gyr. Reduction in XUV flux by a factor of 10, or, equivalently, decreasing the efficiency of XUV-driven loss by a factor of 10, relative to the baseline based on the estimates of Selsis et al. 2007. This increases the range of exoplanets for which a volcanically supported atmosphere is possible. The red dash-dot line corresponds to the line of vanished atmospheres for planets that have all volatiles in the atmosphere initially. This line moves away from the star over time. The dark gray solid and dashed lines correspond to the line of atmospheric revival by volcanic outgassing for planets that lose all atmosphere during transition from a sub-Neptune to a super-Earth. The lines of atmospheric revival sweep towards the star over time because the rate of volcanic degassing falls off more slowly with time than does the star’s XUV flux. This chart assumes initial volatile supply is independent of F and similar to Earth and Venus. Increasing volatile supply will move thresholds to higher F . The light gray lines reduce the solubility of the high-molecular-weight volatile in the magma from 10 -9 Pa -1 to 10 -11 Pa -1 . Note that because most planets are too small to have Earth-like composition and 6 M ⊕ , this plot is optimistic for atmospheric survival (higher planet mass favors atmosphere retention). 45 Fig. S12.
Planet densities from Otegi et al. (2020). The blue line corresponds to density for Earth-like composition, assuming scaling of mass with planet radius as ( R / R ⊕ ) = ( M / M ⊕ ) . The red line corresponds to 80% of that density. Many super-Earth sized planets plot below the red line (for example HD 3167 b), indicating either a volatile envelope (atmosphere) or alternatively a very low (Mg+Fe)/Si ratio relative to that of the Earth. Super-Earths with the same or lower density as Earth (e.g. HD 3167 b) are excellent candidates for having elevated-molecular-weight secondary atmospheres. Planets plotting in the low-density, large radius clump in the bottom right are sub-Neptunes, retaining thick H -dominated atmospheres. Fig. S13.
Stellar luminosity, L , versus time from the model of Baraffe et al. (2015) for star masses ranging from 0.1 – 1.0 M ☉ , and normalized to the Sun’s luminosity at the present-day L ☉ . Contours drawn at intervals of 0.1 M ☉ . The blue line highlights 0.1 M ☉ . The red line highlights 0.3 M ☉ . The black line highlights 1.0 M ☉ . 46 (a) (b) Fig. S14. F XUV at an insolation corresponding to a circular orbit 0.1 AU from today’s Sun, relative to the F XUV at Earth’s orbit today. (a) for the model of Selsis et al. (2007), and (b) combining the fits of Jackson et al. (2012) and Guinan et al. (2016). The dashed lines highlight (from top to bottom) the threshold of F XUV below which all of the high-molecular-mass species is retained by the planet; the F XUV corresponding to a 50% reduction in the no-fractionation loss rate of the high-molecular-mass species; and the F XUV corresponding to escape of the high-molecular-mass species at 90% of the rate at which no loss would occur. Above the top dashed line, fractionation protects the constituents of the secondary atmosphere, and below the lowest dashed line, fractionation is much less important. Calculations are done assuming an atomic wind of H entraining a 1% mixing ratio of atoms of mass 15 Da, and planet mass 6 M ⊕ . 47 Fig. S15.
Andrault et al. (2011) melt curves. Colored lines are magma adiabats. Blue lines are melt fraction curves, where the lowest one is the solidus. We find that the effect of atmospheric overburden pressure on the solidus temperature is relatively small. X axis corresponds to depth within a hypothetical Earth-mass planet (the curves correspond to lines in P-T space, so they map to smaller depth on higher-gravity worlds). 48
Fig. S16.
As Fig. 6, substituting in an XUV-flux-vs.-time parameterization following Selsis et al. (2007). Secondary atmosphere presence/absence model output for 6 M ⊕ (higher planet mass favors atmosphere retention). The dashed red lines show the lines of atmosphere retention after 3.0 Gyr for the case where all volatiles are in the atmosphere initially and there is no primary atmosphere; the 16 th and 84 th percentiles are shown, for varying XUV flux (by ±0.4 dex, 1σ; Loyd et al. 2019) relative to the baseline model following the results of Jackson et al. (2012) and Guinan et al. (2016) (Supplementary Information 1a). These lines move away from the star over time (red arrows). The gray lines show the 16 th percentile for exhibiting an atmosphere after 3.0 Gyr for the case where volcanic outgassing rebuilds the atmosphere from a bare-rock state (the 84 th percentile is at T eq < 400 K). The solid gray lines are for stagnant-lid tectonics and the dashed gray lines are for plate tectonics. The lines of atmospheric revival sweep towards the star over time (gray arrows) because the rate of volcanic degassing falls off more slowly with time than does the star’s XUV flux. In each case the atmosphere/no-atmosphere threshold is 1 bar. The black symbols show known planets that may be tested for atmospheres using JWST (Koll et al. 2019). For any individual planet, star-specific XUV-flux estimates, star age, and the planet’s mass, should be combined to make a more accurate estimate than is possible using this overview diagram. (a) (b) Fig. S17.
Secondary atmosphere presence/absence diagrams for 6 M ⊕ (higher planet mass favors atmosphere retention). The dash-dot lines correspond to the line of atmosphere vanishing for planets that have all volatiles in the atmosphere initially. The solid and dashed lines correspond to the line of atmospheric revival by volcanic outgassing for planets that lose all atmosphere during transition from a sub-Neptune to a super-Earth. The line of revival sweeps towards the star over time because the rate of volcanic degassing falls off more slowly with time than does the star’s XUV flux. The locations of these curves change depending on model assumptions. Increasing volatile supply will move thresholds to higher L . Allowing volatile supply to increase with decreasing T eq , which is realistic, will steepen gradients with L . The gray line is for a constant Earth-scaled outgassing rate (72 bars/Gyr). (b). As (a), but showing changes over time. The gray lines are for 72 bars/Gyr outgassing at different star masses; from top to bottom, {0.3, 0.6, 1.0} M ☉ ..