Atmospheric Evolution on Low-gravity Waterworlds

We assume that the low-gravity waterworld has a pure water vapor atmosphere and a water reservoir fixed at 40% of the planet's total mass. For a given surface temperature, the atmospheric water content and temperature profile are set by the saturation vapor pressure. The water content of the atmosphere remains orders of magnitude below the world's total water content in even the most extreme cases considered. We assume spherical symmetry; this is required by both the radiative transfer model and the hydrodynamic escape formulation. Our model, and therefore the "habitable zone" we consider, neglects the potential effects of CO₂ cycling (Abbot et al. 2012; Ramirez & Levi 2018), ocean chemistry (Wordsworth & Pierrehumbert 2013; Kite & Ford 2018), and three-dimensional dynamics (Pierrehumbert & Ding 2016; Ding & Pierrehumbert 2018).

It is convenient to parameterize the mass and radius of the waterworld in terms of a single quantity, surface gravity g. This is accomplished using the following scaling relation, accurate for waterworlds in the range of 10⁻² M_⊕ < M < M_⊕ (Sotin et al. 2007):

$$\begin{eqnarray}&&\displaystyle \frac{R}{{R}_{\oplus }}=1.258{\left(\displaystyle \frac{M}{{M}_{\oplus }}\right)}^{0.302}.\end{eqnarray} \tag{ 1 }$$

A schematic of our model is shown in Figure 1. The waterworld gains energy from incoming stellar radiation, and loses it via longwave emission and atmospheric escape. We describe this by the following Equation (2), which is similar to those of Pierrehumbert (2010) and Lehmer et al. (2017). However, we specifically account for radiative surface expansion in both the shortwave and longwave regimes:

$$\begin{eqnarray}&&{\left(\displaystyle \frac{{r}_{\mathrm{SW}}}{{r}_{s}}\right)}^{2}\displaystyle \frac{1}{4}(1-A){F}_{\mathrm{stel}}-({g}_{s}{r}_{s}+L){\rm{\Phi }}={\left(\displaystyle \frac{{r}_{\mathrm{LW}}}{{r}_{s}}\right)}^{2}{F}_{\mathrm{out}}.\end{eqnarray} \tag{ 2 }$$

Here A is the planetary albedo, F_stel and F_out are the incoming shortwave and outgoing longwave fluxes, respectively, g_s and r_s are the surface gravity and radius, respectively, L is the latent heat of vaporization of water, Φ is the escaping mass flux, and r_SW and r_LW are, respectively, the radii of the shortwave absorption and longwave emission surfaces. The critical radius r_c, where the outflow reaches the speed of sound w_c, is far enough above the radiative surfaces that the outflow near the radiative surfaces is negligible. Therefore, radiative transfer and hydrodynamic escape can be treated separately (Appendix). In the low-gravity limit, the energy balance equation is consistent with expressions for volatile loss from comets (Lebofsky 1975; Weissman 1980): since there is no atmosphere, the scaling factors disappear, the outgoing flux reduces to the Stefan-Boltzmann law (${F}_{\mathrm{out}}=\sigma {T}_{s}^{4}$), and the gravity term (${g}_{s}{r}_{s}$) becomes negligible. Conversely, in the high gravity limit, the scaling factors and escape terms are negligible, resulting in a standard runaway greenhouse scenario.

Zoom In Zoom Out Reset image size

Here we obtain r_SW, r_LW, and F_out from line-by-line radiative transfer calculations (Section 2.3), and calculate Φ analytically, assuming hydrodynamic escape (Section 2.4). Then, (2) is solved numerically in instellation-gravity space.

To calculate outgoing longwave radiation and the radii of the radiative surfaces, we use the one-dimensional line-by-line radiative-convective model used previously by Schaefer et al. (2016) and Wordsworth et al. (2017). As in Schaefer et al. (2016), this model uses a truncated version of the HITEMP2010 line list (Rothman et al. 2010), which is accurate below 1000 K. We assume that the atmosphere is in phase equilibrium with the surface reservoir and therefore the temperature profile is given by the saturation vapor pressure curve of water. We limit ourselves to clear-sky calculations, accounting for cloud albedo forcing using the standard assumption (Goldblatt et al. 2013; Goldblatt 2015) of an artificially large surface albedo (0.2). The potential effects of other cloud processes are discussed in Section 4.1.

The "radiative scaling factors" from (2), ${\left(\tfrac{{r}_{\mathrm{SW}}}{{r}_{s}}\right)}^{2}$ and ${\left(\tfrac{{r}_{\mathrm{LW}}}{{r}_{s}}\right)}^{2}$, are calculated using the method of Goldblatt (2015); we present a summary here. To a good approximation, emission of outgoing longwave radiation and absorption of incoming shortwave radiation occur at the altitude where optical depth $\tau =1$, measured from the top of the atmosphere. These altitudes can be calculated from the output of the radiative transfer model for a given surface gravity g_s and temperature profile using the equation of hydrostatic balance (taking into account the decrease in gravity with altitude, which is nonnegligible for small worlds):

$$\begin{eqnarray}&&\displaystyle \frac{{dp}}{{dz}}=-\rho {g}_{s}{\left(\displaystyle \frac{{r}_{s}}{{r}_{s}+z}\right)}^{2}.\end{eqnarray} \tag{ 3 }$$

Here r_s is the surface radius, z is the altitude, p is pressure, and ρ is mass density. Even for the most rapidly escaping atmospheres we consider, the velocity of the outflow near the radiative surface is small, such that hydrostatic balance can be assumed (Appendix). Once the heights of unity optical depth, as a function of wavenumber (spectral heights) ${z}_{\tau =1}(\nu )$, are obtained, the height of the radiative surface is calculated by a weighted average of ${z}_{\tau =1}(\nu )$ using the relevant irradiance I(ν).

$$\begin{eqnarray}&&{z}_{\mathrm{radiative}\mathrm{surface}}=\displaystyle \frac{\int I(\nu )\times {z}_{\tau =1}(\nu )d\nu }{\int I(\nu )d\nu }.\end{eqnarray} \tag{ 4 }$$

For the longwave emission surface, I(ν) is the outgoing longwave spectral irradiance given by the radiative transfer calculation. For the shortwave emission surface, it is the Planck blackbody spectral irradiance (hereafter, Planck curve) centered around the stellar effective emission temperature. The effect of deviations from blackbody behavior in real stellar spectra is neglected here.

The nature of the water vapor absorption spectrum is such that ${z}_{\tau =1}(\nu )$ shows a decreasing trend with increasing ν. Therefore, this averaging process means that the emission surface is always above the absorption surface. Therefore, the longwave "scaling factor" is larger than its shortwave equivalent; indeed, their separation grows with temperature. This effect is critical to the rest of the study and is illustrated in Figure 2. The spectral heights of τ = 1 and the scaling factors produced by our model are in good agreement with those of Goldblatt (2015), who considered worlds with 0.12 M_⊕ < M < 10 M_⊕.

Zoom In Zoom Out Reset image size

The problem of formulating the hydrodynamic escape from small bodies has been investigated previously. (Kuramoto & Matsui 1994; Pierrehumbert 2010; Catling & Kasting 2017; Lehmer et al. 2017; Zahnle & Catling 2017). Here we simplify the analysis by restricting ourselves to an isothermal escape formulation using the surface temperature T_s. This is because it is simple, easy to work with, and provides an upper bound on the escape rate. We are ultimately interested in water lifetime and habitability, and isothermal escape will provide a conservative estimate for both.

The isothermal sound speed is $\sqrt{\tfrac{{{kT}}_{s}}{m}}\equiv {w}_{c}$, where k is Boltzmann's constant and m is the molecular mass of water. Assuming transonic escape, there is a critical radius r_c at which this speed is first reached. The outflow in spherically symmetric steady state is described by the continuity equation

$$\begin{eqnarray}&&\displaystyle \frac{d}{{dr}}({r}^{2}\rho w)=0\end{eqnarray} \tag{ 5 }$$

and the momentum conservation equation

$$\begin{eqnarray}&&w\displaystyle \frac{{dw}}{{dr}}=-\displaystyle \frac{1}{\rho }\displaystyle \frac{{dp}}{{dr}}-\displaystyle \frac{{GM}}{{r}^{2}},\end{eqnarray} \tag{ 6 }$$

where ρ is density, w is vertical velocity, p is pressure, G is Newton's gravitational constant, and M is the mass of the body. These are combined (e.g., Catling & Kasting 2017; Lehmer et al. 2017) to form the isothermal planetary wind equation

$$\begin{eqnarray}&&({w}^{2}-{w}_{c}^{2})\displaystyle \frac{1}{w}\displaystyle \frac{{dw}}{{dr}}=\displaystyle \frac{2{w}_{c}^{2}}{r}-\displaystyle \frac{{GM}}{{r}^{2}},\end{eqnarray} \tag{ 7 }$$

where the ideal gas law has also been used. At the critical radius w = w_c and (7) yields

$$\begin{eqnarray}&&{r}_{c}=\displaystyle \frac{{GM}}{2{w}_{c}^{2}}.\end{eqnarray} \tag{ 8 }$$

Typically, r_c ≳ 10r_s. Following Lehmer et al. (2017), we integrate (7) from r_s to r_c and neglect the small kinetic energy term near the surface. This yields the mass flux, Φ:

$$\begin{eqnarray}&&{\rm{\Phi }}={\rho }_{s}{w}_{s}={\rho }_{s}{w}_{c}{\left(\displaystyle \frac{{r}_{c}}{{r}_{s}}\right)}^{2}\exp \left(-\displaystyle \frac{1}{2}+\displaystyle \frac{{GM}}{{w}_{c}^{2}}\left(\displaystyle \frac{1}{{r}_{c}}-\displaystyle \frac{1}{{r}_{s}}\right)\right)\end{eqnarray} \tag{ 9 }$$

Φ is a function of T_s via w_c. The rate of mass loss is calculated as

$$\begin{eqnarray}&&\dot{M}=4\pi {r}_{s}^{2}{\rm{\Phi }}\end{eqnarray} \tag{ 10 }$$

and the water lifetime τ_water is given by

$$\begin{eqnarray}&&{\tau }_{\mathrm{water}}=\displaystyle \frac{{M}_{\mathrm{water}}}{\dot{M}}.\end{eqnarray} \tag{ 11 }$$

Terrestrial-mass worlds can experience the phenomenon known as the "runaway greenhouse" (Komabayasi 1967; Ingersoll 1969; Kasting 1988; Nakajima et al. 1992; Goldblatt & Watson 2012): at high enough surface temperatures, total emission converges toward a limit. If the incoming stellar flux exceeds this threshold, there will be a runaway surface temperature increase, boiling off the oceans—making the world uninhabitable. Therefore, the runaway greenhouse is believed to set the ultimate inner edge of the habitable zone. (Kasting et al. 1993; Goldblatt & Watson 2012; Goldblatt et al. 2013; Kopparapu et al. 2013)

Based on (2) and the preceding analysis, it is clear that low-gravity waterworlds should behave differently. Although the flux per unit area from the effective emission surface will still reach a limit as temperature increases, the emission surface itself can keep gaining area relative to the absorption surface (Figure 2). If the ratio between the longwave and shortwave scaling factors does not reach a limit before escape becomes important, a runaway greenhouse limit will not be reached.

We investigate this quantitatively by studying the steady-state instellation F_stel as a function of surface temperature. Equation (2) can be rewritten in the form

$$\begin{eqnarray}&&{F}_{\mathrm{stel}}=\displaystyle \frac{{\left(\tfrac{{r}_{\mathrm{LW}}}{{r}_{s}}\right)}^{2}{F}_{\mathrm{out}}+({g}_{s}{r}_{s}+L){\rm{\Phi }}}{{\left(\tfrac{{r}_{\mathrm{SW}}}{{r}_{s}}\right)}^{2}\tfrac{1}{4}(1-A)}.\end{eqnarray} \tag{ 12 }$$

F_stel(T) is plotted in Figure 3, for a given g = 2.3 ms⁻², corresponding to M ≃ 0.08 M_⊕. This value is illustrative because the escape rate increases slowly enough with T that the radiative surface expansion effect can be clearly observed. Water lifetimes are also added: they are a nonlinear function of temperature. For comparison, we also plot the results of two other radiation assumptions: the same line-by-line calculation with radiative surface expansion ignored, and blackbody emission from the surface ($\sigma {T}_{s}^{4}$) until the runaway greenhouse limit is reached. The latter approach was used by Lehmer et al. (2017), to study Europa- and Ganymede-like worlds. Hydrodynamic escape is assumed to be governed by (9) in all cases. This is responsible for the rapid increase in F_stel(T) as higher temperatures are approached. We note that such states with high mass flux are "quasi-steady": the rapid expulsion of water into space is not sustainable indefinitely.

Zoom In Zoom Out Reset image size

A given steady-state instellation will allow for a habitable (long-lived liquid water) state if the surface temperature is larger than 273 K and the water inventory is long-lived enough. In Figure 3 and the rest of this study we choose 10⁹ yr as a cutoff; this choice is further discussed in Section 3.2. In our model, where radiative surface expansion is considered, F_stel(T) never reaches a limit. Therefore, there is no runaway greenhouse. This leads to a relatively wide habitable zone. In contrast, neglecting radiative surface expansion but still using the line-by-line model would predict an extremely thin habitable zone, and assuming blackbody emission until the runaway greenhouse would predict no habitable zone at all.

We now expand the picture of Section 3.1 and Figure 3 into gravity space, allowing the instellation to be set by the orbital distance. We choose a water lifetime of 10⁹ yr as an inner habitability boundary. While the principal requirement for the water lifetime habitability boundary is that it feasibly allows the evolution of life, this particular choice of timescale also ensures that such worlds have a reasonable chance of being habitable at the time of observation—even when the system's age is uncertain. With this definition, the "longer-lived states" in our model consistently have surface temperatures below proposed upper bounds for survival of life (Bains et al. 2015), suggesting that temperature is not additionally prohibitive.

As can be observed from Figure 2, stellar emission temperature plays a role in determining the size of the shortwave scaling factor. We thus consider two limiting cases: an M-star (modeled on AD Leonis) and a G-star (modeled on the Sun). The resulting steady-state temperatures and the 10⁹ yr water lifetime contour are shown in Figure 4. We label five qualitatively different climate states: long-lived snowball, long-lived temperate, short-lived snowball, short-lived temperate, and short-lived hot. The mass of Ganymede and a 10⁷ yr water lifetime contour are plotted additionally for reference.

Zoom In Zoom Out Reset image size

The "habitable zone" is defined by the area enclosed by the 273 K temperature contour and the 10⁹ yr water lifetime contour. The conservative maximum width in the low-gravity regime is ∼0.08 au for G-stars and ∼0.008 au for M-stars. Relaxing the 10⁹ yr habitability boundary would extend the habitable zone both spatially and to lower gravities. Worlds with appropriate geochemical cycling could of course have a wider habitable zone, because the outer edge would move outwards.

The width of the habitable zone shown in Figure 4 decreases with gravity, ultimately disappearing. This marks the transition from "planet-like" (habitable at some stellar flux) to "comet-like" (never habitable regardless of stellar flux). For worlds with a mass–radius scaling accurately described by (1), our model conservatively gives this boundary at g ≃ 1.48 ms⁻², or M/M_⊕ ≃ 0.0268. The boundary is essentially independent of stellar spectral type, because the differences in the shortwave scaling factor are negligible for waterworld surface temperatures near 273 K.

Although our model predicts a habitable climate state, such a state is not necessarily easy to access. Figure 4 exhibits a range of possible climate states: short-lived hot, temperate or snowball, and long-lived temperate or snowball. The ice-albedo feedback can hamper transitions from snowball states to temperate states: this has already been demonstrated for terrestrial-mass worlds (Yang et al. 2017). We can incorporate the ice-albedo feedback into our model using the simple albedo step function

$$\begin{eqnarray}A({T}_{s})=\left\{\begin{array}{ll}{\alpha }_{I} & {T}_{s}\lt 273\ {\rm{K}}\\ {\alpha }_{L} & {T}_{s}\geqslant 273\ {\rm{K}}\end{array}\right.,\end{eqnarray} \tag{ 13 }$$

where α_I is the albedo of the icy (snowball) state, and α_L is the albedo when there is surface liquid water. Hysteresis plots for different choices of α_I are shown in Figure 5. We observe that a snowball state experiencing a stellar flux-driven deglaciation generally bypasses the long-lived state entirely, except for very low α_I values. Although the mechanism setting the inner edge of the habitable zone is different, the conclusion of the habitable state likely being bypassed upon stellar flux-driven deglaciation is the same as that of Yang et al. (2017).

Zoom In Zoom Out Reset image size

Beyond those considered in this study, other factors play a role in determining the habitability of low-gravity waterworlds. For example, we have restricted ourselves to pure water vapor atmospheres; in reality, a variety of species could be present, depending on interior and atmospheric chemistry. We did not consider the potential effects of clouds except via an increased planetary albedo. Cloud greenhouse forcing could reduce outgoing radiation as it does for terrestrial-mass worlds, and clouds in an expanded low-gravity atmosphere could increase absorbed flux by scattering more stellar radiation toward the surface. Depending on whether cloud effects contribute more to cooling or warming, the size of the habitable zone presented in this study could either increase or decrease. The net effect should be carefully assessed in the future using three-dimensional models (e.g., Pierrehumbert & Ding 2016). Finally, as in our solar system, low-gravity waterworlds may often be found as exomoons around giant planets. Then, tidal forces are important for habitability (Heller & Barnes 2013), becoming especially deleterious around M-stars (Zollinger et al. 2017).

Even if a habitable state exists, it is not necessarily "accessible," as we showed in Section 3.4. Worlds with a warm history have access to the entire habitable zone derived here; the proto-atmosphere phenomenon discussed by Kuramoto & Matsui (1994) could make a warm start more likely. However, for worlds with a cold history, stellar-flux-driven deglaciation into a long-lived state only occurs for artificially low ice-albedo values. Somewhat more favorable conditions might be found around M-stars, where the ice-albedo feedback is weakened (Joshi & Haberle 2012; Shields et al. 2013). However, such systems experience another problem: the dramatic pre-main-sequence evolution of the star may strip potentially habitable worlds of their water early on (Ramirez & Kaltenegger 2014; Luger & Barnes 2015; Tian & Ida 2015). Further research is required to integrate all of these considerations.

In this study we have used the conventional but limited definition of habitability as the ability to maintain surface liquid water on sufficient timescales. While our model categorizes worlds of Europa's or Ganymede's size as "comet-like" (unable to maintain surface liquid water for gigayear timescales regardless of stellar flux), that does not necessarily make them uninhabitable. For example, there is a large body of work on how life could survive on Europa within a subsurface ocean, supported by energy from a geologically active interior and favorable chemistry (Reynolds et al. 1983; Gaidos et al. 1999; McCollom 1999; Chyba 2000). Detection of such life remotely would be extremely challenging. Thus the question of life on "comet-like" worlds ultimately remains open, if unanswerable at present.

Observation of low-gravity waterworlds outside of our solar system is possible in principle. The smallest exoplanet discovered to date is the rocky Kepler 37-b, with a radius of 0.303 R_⊕ (Barclay et al. 2013). This is already smaller than all of the low-gravity waterworlds modeled in this study. Regarding exomoons, a previous search using the Kepler Space Telescope (Kipping et al. 2015) was sensitive to masses as low as 1.7 M_Ganymede ≃ 0.043 M_⊕, which includes most of the mass range for habitable low-gravity waterworlds shown in Figure 4. Unfortunately, due to the nature of transit photometry, low-gravity worlds can currently only be found in orbits with relatively small semimajor axes; this is likely to preclude long-lived water. Correctly classifying low-gravity waterworlds presents another issue: under current observational uncertainties in this size regime (e.g., Jontof-Hutter et al. 2015), mass–radius determinations alone are likely to be insufficient (Sotin et al. 2007).

The phenomena described in this study will have observable effects. Radiative surface expansion increases transit depths; planetary density will be underestimated by a factor of ${r}_{\mathrm{SW}}^{3}/{r}^{3}$, a decrease of up to 15% for temperate worlds. Clouds or escaping dust will amplify this: the latter already provides a means to observe the much smaller exocomets (Kiefer et al. 2014; Rappaport et al. 2017), and has also been proposed to explain anomalously low-density sub-Neptunes (Wang & Dai 2019). Escaping H atoms could produce a detectable Lyα signal; for exomoons, these atoms would likely accumulate in a torus around the host planet (Lehmer et al. 2017). Transmission spectroscopy has already been successfully used to study Earth-mass planets (e.g., de Wit et al. 2016); for low-gravity waterworlds, the highly expanded atmosphere will amplify the variations in such a spectrum. Coupling the model developed here with a forward model for observables would be an interesting follow-up to this study, and would aid future detection efforts.

Further work could consider more complicated models of hydrodynamic escape. Since the isothermal escape approximation represents an upper bound, such work is likely to widen the habitable zone presented in Figure 4.

The two key phenomena governing evolution on low-gravity waterworlds, radiative surface expansion, and hydrodynamic escape, are treated separately in our model. This is possible while the outflow near the radiative surface is negligible, which is the case for the parameter space regime we consider. While our model is sufficient to estimate constraints on habitability, the scenario of nonnegligible outflow near the radiative surface remains of physical interest, and should be explored further.

Finally, the planet–comet transition calculated in this study represents a hard lower gravity boundary for habitability in the conventional sense. An upper gravity boundary is likely given by the onset of runaway gas accretion, which distinguishes super-Earths from gas giants (e.g., Ginzburg et al. 2016). Such boundaries, when used in tandem, could be a useful tool when estimating a total number of potentially habitable worlds.