NIR-driven Moist Upper Atmospheres of Synchronously Rotating Temperate Terrestrial Exoplanets

We use the ROCKE-3D GCM (Way et al. 2017) to obtain 3D atmospheric structures of terrestrial planets with surface water under varying irradiation, focusing on the H₂O profile in the upper atmosphere. ROCKE-3D is a generalization of the ModelE2 GCM (Schmidt et al. 2014), which has been developed for the Earth at NASA's Goddard Institute for Space Studies. The details of the ROCKE-3D GCM are presented in Way et al. (2017).

In the following, we first describe the specific configurations in our experiments to clarify our targets. Then, we also discuss in some detail the convection/cloud parameterizations and the radiation parameterization used in our GCM, as these are of particular importance to our results.

2.1.1. Our Configurations

In our experiments, the radius and the surface gravity of the planet are set to Earth's value.

We only consider aquaplanets, i.e., planets whose surfaces are wholly covered with water. We use a thermodynamic ocean with 50 m depth for the fiducial runs for simplicity. However, in Section 4.1, we also present the results of runs with a fully dynamic ocean, to show that our conclusions are not very sensitive to ocean heat transport.

Our atmosphere has a mean surface pressure of 1 bar, and is composed mainly of N₂, with only 1 ppm CO₂, as in Kopparapu et al. (2016). In Section 4.2, we will briefly discuss the dependence of our results on the mixing ratio of CO₂. We did not include the effects of aerosols and photochemistry.

The model was configured to have a horizontal grid with 4 degree resolution in latitude and 5 degree resolution in longitude, and 40 vertical layers covering up to 0.14 mbar (corresponding to ∼65 km in the case of the Earth).

In our experiments, we mainly changed two parameters: the total incident flux and the spectral type of the star. In the following, we represent the incident flux normalized by the solar constant by S_X. We consider four types of stars: Sun (G2V), HD 22049 (K2V), Kepler-186 (M1V), and GJ 876 (M4V). The spectra of HD 22049 and GJ 876 are obtained from Segura et al. (2003) and Domagal-Goldman et al. (2014), respectively, while the spectrum of Kepler-186 is a modeled spectrum from BT-Settl (Allard et al. 2012) based on the stellar parameters reported in Table S1 of Quintana et al. (2014). The physical parameters of these stars, as well as the references, are summarized in Table 1, and the spectra are displayed in Figure 1.

Zoom In Zoom Out Reset image size

Table 1.  Stellar Properties

Star	Type	${T}_{\mathrm{eff}}$	L [L⊙]	M [M⊙]	R [R⊙]	Reference	Source of Spectrum
Sun	G2V	5772 K	1	1	1	⋯	Lean et al. (2005), SPARC/SOLARISa
HD 22049	K2V	5039 K	0.32	0.82	0.74	Baines & Armstrong (2012)	Segura et al. (2003)
Kepler-186	M1V	3755 K	0.055	0.544	0.523	Torres et al. (2015)	Allard et al. (2012)
GJ 876	M4V	3129 K	0.0122	0.37	0.3761	von Braun et al. (2014)	Domagal-Goldman et al. (2014)

Note.

^a http://solarisheppa.geomar.de/ccmi

Download table as:  ASCIITypeset image

Throughout the paper, we exclusively consider synchronously rotating planets due to their relevance to future transmission observations of close-in systems, as discussed in Section 1. The orbital period (identical to the spin period) is changed according to Kepler's 3rd law, to be consistent with the total incident flux the planet receives and the stellar mass. Thus, the incident flux and the spin period are changed simultaneously for the fiducial set of runs. The orbital/spin periods corresponding to S_X = 1 are shown in Table 2. We also isolate the effect of the spin period in Section 4.3.

Table 2.  Orbital Period, to Which the Rotation Period is Equated, in the Case of S_X = 1

Star	${T}_{\mathrm{eff}}$	Period (days) (SX = 1)
Sun	5772 K	365.3
HD 22049	5039 K	171.6
Kepler-186	3755 K	56.2
GJ 876	3129 K	22.0

Download table as:  ASCIITypeset image

Our fiducial models with a thermodynamic ocean typically equilibrate in ∼30 Earth years in model time. After the steady-state is reached, we average the 3D atmospheric profiles over an integer number of orbits that is approximately equal to 10 Earth years. All of the results presented in this paper are based on the averaged data created this way.

2.1.2. Cloud Scheme

2.1.3. Radiation Scheme

We also simulate transmission spectra based on the outputs of our GCM experiments, setting the modeled planet in a transiting geometry, in order to quantify the absorption features of water vapor. We follow the procedure described in Appendix B of Way et al. (2017), which is briefly summarized below.

First, we obtain the vertical profiles of temperature, water vapor, and cloud particles along the terminator (i.e., dayside–nightside boundary) by linearly interpolating the profiles of the GCM gridpoints adjacent to the terminator. Then, we trace the optical paths of the transmitted ray according to the profile of the refraction index, from the observer toward the stellar disk, in the same way as described in van der Werf (2008), following Misra et al. (2014). After that, the spectral opacity at different points along the path is calculated considering Rayleigh scattering, the absorption by atmospheric molecules, and absorption and scattering by cloud particles, excluding the effects of multiple scattering. The spectral signatures of the planetary atmosphere are discussed in terms of the transit depth ${\rm{\Delta }}F(\lambda )$ where λ is the wavelength, as well as the effective height, ${h}_{\mathrm{eff}}$, which is given as a function of planet radius, ${R}_{{\rm{p}}}$, and stellar radius, ${R}_{\star }$, by

$$\begin{eqnarray}&&{\rm{\Delta }}F(\lambda )=\displaystyle \frac{{({R}_{{\rm{p}}}+{h}_{\mathrm{eff}}(\lambda ))}^{2}}{{R}_{\star }^{2}}.\end{eqnarray} \tag{ 1 }$$

To compute the absorption coefficients of atmospheric molecules, we use line-by-line opacity tables based on the HITRAN2012 database (Rothman et al. 2013), prepared separately from the opacities used in the GCM.

Including the optical effects of clouds is tricky. Strictly speaking, one needs to know the size distribution and the spatial distribution of the cloud particles at the timescale of the planetary transit. In this paper, when calculating transmission spectra, we consider the effect of clouds only in an approximate manner with four assumptions. First, we assume that the radius of the cloud particles ${r}_{\mathrm{cld}}$ is ${r}_{\mathrm{cld}}=30\,{\rm{\mu }}{\rm{m}}$ (regardless of their phase), similar to the typical values of ice cloud particles on the Earth as well as in our GCM simulations. Second, we assume that the scattering cross section is given by the geometric optics limit, i.e., $2\pi {r}_{\mathrm{cld}}^{2}$. Third, we assume that cloud particles are uniformly distributed within each grid cell. Fourthly, we compute the transmitted light based on the optical properties averaged over 10 Earth years, rather than compute it based on the instantaneous optical properties before taking the temporal average. The third and fourth points above imply that we find the effects of clouds to be the largest among the possible realizations given the particle size distribution, because of the nonlinear dependence of transmittance on the optical thickness. Under these assumptions, the absorption coefficient k of the clouds is

$$\begin{eqnarray}&&k=2\pi {r}_{\mathrm{cld}}^{2}\cdot n=3\rho {x}_{\mathrm{cld}}/(2{\rho }_{{{\rm{H}}}_{2}{\rm{O}}}{r}_{\mathrm{cld}}),\end{eqnarray} \tag{ 2 }$$

where n is the number density of cloud particles, ρ is the density of the ambient atmosphere, ${x}_{\mathrm{cld}}$ is the cloud water content provided by the GCM experiments, and ${\rho }_{{{\rm{H}}}_{2}{\rm{O}}}$ is the density of liquid/solid water, which is assumed to be $1\,{\rm{g}}\,{\mathrm{cm}}^{-3}$. The uncertainties in the particle radius, scattering cross section, and the cloud water content could change the absorption coefficient by a factor of a few.

Before examining the dependence of the water vapor mixing ratio in the upper atmosphere, we discuss the 3D structures of temperature and water vapor, as we need to know at which altitude the mixing ratio of water vapor becomes relatively uniform and hence can be regarded as representative of the humidity in the upper atmosphere.

Figure 2 shows the structures of temperature and water vapor mixing ratio along the equator for the planets around GJ 876 (M4V) and around the Sun (G2V), respectively, for varying incident flux normalized by the solar constant (S_X). In all panels the substellar point corresponds to $0^\circ$ longitude. In general, the temperature and the water vapor mixing ratio near the surface strongly depend on the distance from the substellar point, but are fairly uniform in the upper part of the atmosphere. Above about 1 mbar, which is located between 35 and 50 km altitude, depending on the atmospheric temperature profiles, the water vapor mixing ratio becomes horizontally uniform within a few tens of percent. Therefore, in the following, we denote the globally averaged water vapor mixing ratio at 1.39 mbar by ${x}_{{{\rm{H}}}_{2}{\rm{O}}}$, and use it as representative of the water vapor abundance in the upper atmosphere.

Zoom In Zoom Out Reset image size

We also observe optically thick clouds around the substellar point and a resulting increase of the planetary albedo as a function of total incident flux (not shown in figures), consistent with Yang et al. (2013, 2014), Kopparapu et al. (2016), and Way et al. (2016). As shown in Figure 2, the maximum surface temperature does not exceed 320 K, even in the highly irradiated cases in our experiments where the water vapor mixing ratio near the model top is as high as ${10}^{-3};$ we will discuss this in the subsequent sections.

In this section, we examine in more detail how the water vapor mixing ratio in the upper atmosphere depends on the incident flux (both the total flux and the spectral distribution).

The left panel of Figure 3 presents ${x}_{{{\rm{H}}}_{2}{\rm{O}}}$ versus S_X. Overall, we find an increase of ${x}_{{{\rm{H}}}_{2}{\rm{O}}}$ as a function of S_X. However, the slope is much more gradual than suggested by 1D models, in particular those of Kasting et al. (1993) and Wordsworth & Pierrehumbert (2013) with a small CO₂ mixing ratio, which typically predict a rapid transition from an Earth-like dry stratosphere to an H₂O-dominated atmosphere. Our result that the water vapor mixing ratio responds more slowly to the incident flux is qualitatively consistent with Yang et al. (2013), who studied synchronously rotating planets using CAM3 with an emphasis on the stability of the climate. Note that we cannot directly compare to their results, as our assumptions do not match perfectly; specifically, Yang et al. (2013) assumed a planet with a larger radius ($R=2{R}_{\oplus }$) and higher gravity ($g=1.4{g}_{\oplus }$), and a 60-day orbital period. The dependence on these other planetary parameters will be investigated in a future paper.

Zoom In Zoom Out Reset image size

Because these are synchronously rotating planets where the surface temperature is stabilized by the optically thick clouds around the substellar point, one might think that this stabilized surface temperature would explain the slow response of water vapor to the increased total incident flux. To see the effect of the surface temperature, we have plotted ${x}_{{{\rm{H}}}_{2}{\rm{O}}}$ as a function of surface temperature around the substellar point, where convection is strongest (Figure 3, middle panel). It is then clearly seen that the water vapor mixing ratio starts to rise at a fairly modest surface temperature and increases even more dramatically as a function of surface temperature than 1D models suggest. This indicates that the increasing surface temperature cannot by itself explain the increased stratospheric water vapor mixing ratios.

Instead, we find that the near-infrared (NIR) component of the incident flux appears to be a good predictor of stratospheric humidity across the full range of stellar spectral types (Figure 3, the right panel). Here, we integrated the star spectra between 0.9 and 3 μm (the interval is somewhat arbitrary) to obtain the NIR component of the incident flux. We will examine the reason for this in the next subsection.

In Figure 4, we show the vertical profiles of temperature, water vapor mixing ratio, short-wave heating rate, and upward velocity, averaged over the region within a distance $0.2\,{R}_{\oplus }$ from the substellar point. The variability of the profiles in the individual gridboxes that contribute to the average is negligible compared to the changes due to varying S_X. Rising motion exists at most altitudes in the substellar region. In the troposphere, where moist convection occurs, large-scale vertical motion is relatively insensitive to S_X. Above the convective layers, though, large-scale upward motion increases strongly with incident flux, concurrently with the increase in water vapor mixing ratio and heating rate in the upper atmosphere.

Zoom In Zoom Out Reset image size

Here, in the upper atmosphere we see the interplay among the radiative heating, water abundance, and the vertical air transport: when the incident flux increases, the radiative heating due to the absorption of short-wave radiation by water vapor and potential cloud particles, mainly in the NIR regime, also increases. The air then reacts by rising and adiabatically cooling. The induced upward motion of the atmosphere transports water vapor from the moist lower atmosphere to the upper atmosphere, lowering the vertical gradient of water vapor and increasing the moisture at high altitude, further increasing the short-wave heating.

To be more specific, the evolution of the profile satisfies the conservation equation for the potential temperature:

$$\begin{eqnarray}&&\displaystyle \frac{T}{\theta }\left(\displaystyle \frac{\partial \theta }{\partial t}+({\boldsymbol{v}}\cdot {\rm{\nabla }})\,\theta \right)=\displaystyle \frac{{Q}_{r}+{Q}_{l}}{\rho {c}_{p}},\end{eqnarray} \tag{ 3 }$$

or

$$\begin{eqnarray}&&\displaystyle \frac{T}{\theta }\displaystyle \frac{\partial \theta }{\partial t}=-\displaystyle \frac{T}{\theta }\left(u\displaystyle \frac{\partial \theta }{\partial x}+v\displaystyle \frac{\partial \theta }{\partial y}\right)-\displaystyle \frac{T}{\theta }\left(\omega \displaystyle \frac{\partial \theta }{\partial p}\right)+\displaystyle \frac{{Q}_{r}}{\rho {c}_{p}}+\displaystyle \frac{{Q}_{l}}{\rho {c}_{p}},\end{eqnarray} \tag{ 4 }$$

where p is the pressure, Q_r and Q_l are the radiative and latent heating rates, c_p is the specific heat of the atmosphere at constant pressure, ρ is atmospheric density, and $\{u,v,\omega \}$ are the conventional velocity components toward the east ($x-$axis), north ($y-$axis), and upward ($z-$axis; expressed here in pressure coordinates as $\omega \equiv {dp}/{dt}$), respectively. Given that the horizontal advection (the first term on the right side) has a negligible contribution since the horizontal gradients are small (Figure 2), and that the latent heating (the last term on the right side) is also negligible above the altitude to which convection penetrates, the increase of the radiative heating in the upper atmosphere must be balanced by the adiabatic cooling by vertical advection. The contribution of each term of the right side of Equation (4) in the steady-state is shown in Figure 5, for the substellar regions of low and high-irradiation cases with GJ 876 and the Sun. It is clearly seen that the net radiative heating in high-irradiation cases is compensated by the cooling by the vertical advection.

Zoom In Zoom Out Reset image size

Figure 6 puts this upward motion in the substellar regions into a global context by presenting the velocity fields in the altitude–longitude plane along the equator for the low and high-irradiation cases for GJ 876 and the Sun. Here, the vertical components of the vectors are multiplied by 1000 to make them more visible, because large-scale vertical motions are about three orders of magnitude weaker than horizontal motions. The development of the large-scale circulation above the troposphere is clearly seen.

Zoom In Zoom Out Reset image size

Altogether, Figures 4–6 suggest that the stratospheric dynamical response to increased stellar radiation is similar for stars of different spectral type, but for a warmer star the onset of strong upward motion occurs at a much larger incident flux because less of the instellation is in the NIR where the radiation is more efficiently absorbed by water vapor and clouds. This is the reason for the behavior seen in the right panel of Figure 3.

We also note that the contribution of radiative heating by cloud particles above convective layers is not negligible. In order to see the radiative effect of clouds above convective layers, we created models that are comparable to the fiducial models, but where the influence of cloud particles above the 237 mbar level (∼11 km for the runs shown below) has been turned off in the radiation calculations (i.e., cloud particles in the upper atmosphere are transparent to radiation). Figure 7 displays the comparison between such "transparent high clouds" models (dashed black, dashed gray) and fiducial models (red, pink), for S_X = 1.0 and 1.2 with the incident spectrum of GJ 876. While these high clouds have a minimal effect on the surface temperature and atmospheric profiles in the convection layers, it is seen that radiative heating above the convection layer is substantially suppressed without the radiative effect of clouds, leading to a colder, and drier upper atmosphere. Note that we still see heating attributed to water vapor in the "transparent high clouds" models; for example, the water vapor mixing ratio and short-wave heating at 20–30 km altitude increase considerably from S_X = 1.0 to S_X = 1.2, even in the absence of cloud radiative effects, but the sensitivity to S_X is weaker at higher altitudes. We also ran equivalent simulations for other types of stars and found that the overall effects of the "transparent high clouds" assumption are very similar.

Zoom In Zoom Out Reset image size

In summary, atmospheric vertical motion and radiative heating by water vapor and cloud particles interact with each other through positive feedbacks, and thus a modest increase of the incident radiative flux, S_X, can result in an order-of-magnitude or more increase of the water vapor mixing ratio in the upper atmosphere, while the surface temperature increases only modestly due to the stabilization by clouds. This picture is verified by a good correlation between the water vapor mixing ratio at 1.39 mbar and the NIR portion of the incident flux as shown in the right panel of Figure 3, insensitive to the stellar spectral types. Because the NIR flux is an external parameter that can be found from direct observations of the host star, it would be a useful guide to predict the high-altitude water vapor mixing ratio for a planet with surface liquid water.

In this subsection, we show the dependence of the transmission spectra on the incident flux, taking the modeled planets around one of our M-type stars, GJ 876, as examples. Figure 8 shows the transmission spectra for the S_X = 0.6, 1.0 and 1.2 cases, where the water vapor mixing ratio at 1.39 mbar is roughly 10⁻⁶, 10⁻⁴ and $3\times {10}^{-3}$, respectively. The corresponding vertical profiles in the substellar region are shown in Figure 4, and note that the water vapor mixing ratio is even higher at medium altitudes (10–30 km). Although refraction due to the planetary atmosphere is included, its effect on the transmission spectra of planets around late-type stars like GJ 876 is minor (Bétrémieux & Kaltenegger 2014; Misra et al. 2014).

Zoom In Zoom Out Reset image size

Without the opacity due to clouds (colored dashed lines in Figure 8), the effective altitude of the water vapor absorption features increases from ∼15 km to ∼40 km when S_X increases from 0.6 to 1.2. This level of effective altitude of water signatures is comparable to the strongest features (CO₂) in the modeled visible-to-NIR transmission spectrum of the Earth (solid black line), and substantially larger than the H₂O signatures of the Earth. Note that the feature around 2.7 μm is the overlap between H₂O and CO₂, and in our simulation cases with 1 ppm of CO₂ the feature predominantly comes from H₂O, while in the case of the Earth it mostly comes from CO₂. Thus, for the highly irradiated cases shown here, H₂O signatures will be a natural target when trying to observe the strongest atmospheric signatures of planets similar to the Earth.

However, we also find clouds in the upper atmosphere around the terminator, and the inclusion of clouds substantially raises the baseline of the transmission spectra. (Note that the high clouds discussed in the previous sections are those around the substellar point, while what matter here for the transmission spectra are those at the terminator.) These clouds at the terminator presumably form from water vapor that is transported from the substellar point and cools radiatively or adiabatically. Although these tenuous clouds at high altitude would be optically thin vertically, the grazing geometry of transmission observations causes them to significantly affect the transmission spectra due to the long optical path length. After including the effect of clouds, however, we still find the overall increase of water vapor signatures as the incident flux increases. Compared to the standard Earth model, or the low-irradiation case (water vapor mixing ratio $\sim {10}^{-6}$) with the same assumptions for the cloud optical properties, the signal of a high-irradiation case could be larger by a factor of a few.

We can evaluate the detectability of these features for idealized cases where the only source of observational uncertainty is photon noise. For an Earth-sized planet 10 pc away, assuming a 6.5 m telescope (having the James Webb Space Telescope in mind) with 40% throughput and a wavelength resolution ${ \mathcal R }=20$, absorption features with 20 km differential effective altitude in this wavelength range (e.g., the feature around 2.7 μm for the cases with clouds) would become detectable (5σ) after ∼100 hr of integration with an M4 star (assuming blackbody radiation with ${R}_{\star }=0.4\,{R}_{\odot }$, T = 3500 K, similar to GJ 876), or after ∼15 hr with an M8 star (with ${R}_{\star }=0.1\,{R}_{\odot }$, T = 2500 K, similar to TRAPPIST-1 planets), respectively.

So far, we have used a thermodynamic ocean model with zero ocean heat transport and 50 m depth. Indeed, this simplified treatment has been regularly used in previous 3D atmospheric modeling studies of exoplanets. A more realistic approach is to include the dynamics of the ocean. Therefore, we check whether our results are sensitive to the prescription of the ocean by comparing our fiducial results and those obtained using a dynamic ocean module.

For this purpose, we ran the models for GJ 876 assuming a 900 m deep dynamic ocean. The models with the dynamic ocean typically equilibrated after a hundred to several hundred Earth years. In the same way as the fiducial models, we averaged the atmospheric profiles over approximately 10 Earth years after the model reached a steady-state.

Figure 9 shows ${x}_{{{\rm{H}}}_{2}{\rm{O}}}$ for GJ 876 with varying S_X both for our fiducial runs and those using the dynamic ocean module. There are minor differences between the two, with the maximum difference being less than an order of magnitude. The models with a dynamic ocean tend to have somewhat less ${x}_{{{\rm{H}}}_{2}{\rm{O}}}$ than the models with zero ocean heat transport. This appears to be related to the fact that ocean heat transport from the dayside to the nightside lowers the substellar surface temperature and weakens convection there, and thus there is less transport of water vapor to the upper atmosphere.

Zoom In Zoom Out Reset image size

Overall, however, our fiducial models with a thermodynamic ocean give a reasonable order-of-magnitude estimate of water vapor mixing ratio in the upper atmosphere.

Our fiducial models assumed 1 ppm CO₂, which was an arbitrary choice. Although a full investigation of the dependence on the CO₂ mixing ratio is beyond the scope of this paper, here we briefly mention how much the results are changed if we assume an Earth-like CO₂ mixing ratio. The results of the models with 300 ppm of CO₂ are plotted in Figure 9. As shown, this level of variation of CO₂ has a very small effect on the water vapor profile. Larger amounts of CO₂ up to that of a CO₂-dominated atmosphere could have a more substantial effect on the 3D profiles (Wordsworth & Pierrehumbert 2013).

In our fiducial set of experiments, we changed the planetary spin period together with the stellar flux for a given star, because we consider synchronously rotating planets and the orbital period (which is equal to the spin period) must change to be consistent with the total incident flux and the host star mass according to Kepler's 3rd law. Thus, strictly speaking, the change in the water vapor mixing ratio presented in Figure 3 includes the effect of both the varying incident flux and the varying spin period. Meanwhile, it is well-known that the spin period has first-order effects on the atmospheric circulation and hence can potentially affect ${x}_{{{\rm{H}}}_{2}{\rm{O}}}$.

In order to isolate the effects of rotation, we performed experiments where we varied the spin period while keeping the incident flux fixed at S_X = 1 with the GJ 876 spectrum, but allowing the planet to remain in synchronous rotation. The effect of varying spin period on ${x}_{{{\rm{H}}}_{2}{\rm{O}}}$ is presented in Figure 10. It is considerably smaller than the dependence on S_X seen in Figure 3, i.e., the increase in stratospheric water vapor as the planet is moved closer to its star is mostly an intrinsic sensitivity to the stellar flux rather than an effect of the changing planet rotation as the orbital period changes.

Zoom In Zoom Out Reset image size

Nonetheless, the spin effect is non-negligible; we find that ${x}_{{{\rm{H}}}_{2}{\rm{O}}}$ becomes smaller as the spin period becomes very short, i.e., increasing planet spin partly offsets the effect of increasing stellar flux. Our interpretation of this tendency is as follows. At slow rotation (long spin period) the day–night circulation dominates the atmosphere. As the spin period becomes shorter, however, the atmosphere develops high-altitude zonal winds, as demonstrated in previous papers (Merlis & Schneider 2010; Edson et al. 2011; Kopparapu et al. 2016). The super-rotating zonal wind seen in the equatorial region in the shorter period runs (not shown) efficiently transports water vapor from the dayside to the nightside, limiting the cloudiness in the substellar region, and resulting in a higher surface temperature (Kopparapu et al. 2016). Consistently, in the runs with a shorter spin period, the convection at the substellar point is deeper, and the temperature at the top of the convection where it moistens the atmosphere is colder, and thus the convection is less able to moisten the atmosphere there. Moreover, the water vapor transported upward by convection is more efficiently transported to the nightside by the zonal wind, as noted above; the water vapor mixing ratio on the nightside increases by about an order of magnitude when the rotation period decreases from 365 Earth days to 4 Earth days. These processes are likely to lead to the smaller upper atmosphere water vapor mixing ratio observed in our shorter period runs.

We have shown that the radiative heating caused by absorption of the incident NIR radiation by water vapor is one of the key processes governing the upper atmosphere humidity, and the water vapor continuum absorption is partly responsible for this. There is significant uncertainty in the strength of both the self- and foreign-broadened water vapor continua at NIR wavelengths (Campargue et al. 2016; Shine et al. 2016). Our fiducial models used the CAVIAR water vapor continuum (Ptashnik et al. 2011), which from recent laboratory measurements appear to overestimate particularly the self-broadened water vapor continuum (Shine et al. 2016). Another popular continuum model, MT_CKD (Mlawer et al. 2012), or its predecessor CKD (Clough et al. 1989), appears to be more in line with, or a bit below, the recent measurements of the self-broadened continuum. Measurements of the foreign-broadened continuum are sparse, but more consistent with CAVIAR as MT_CKD (and CKD) appears to underestimate the strength of this continuum (Shine et al. 2016).

In order to evaluate the impact of the water vapor continuum model on our results, we have also run simulations with the CKD 2.4 water vapor continuum. As CAVIAR and CKD 2.4 approximately bracket the available measurements, this should give a good indication of the sensitivity of our results to the adopted continuum model.

Figure 11 compares ${x}_{{{\rm{H}}}_{2}{\rm{O}}}$ with the CAVIAR model (our fiducial model) and those with the CKD 2.4 model, indicating that the effect of the water vapor continuum is very small. We also compare in Figure 12 the vertical profiles of temperature and water vapor in the substellar region for the most highly irradiated (i.e., most moist) runs with the GJ 876 and Solar spectra. The resulting profiles overall closely match.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

When the atmospheric structure is fixed, the difference in short-wave fluxes due to the different water vapor continua is observed predominantly in the lower atmosphere (<10 km), as expected since the strength of continuum absorption is ∝P², where P is the pressure. Using the P–T and specific humidity profiles at the substellar point from our fiducial run with the spectrum of GJ 876 and S_X = 1.2, we find that the net short-wave clear-sky flux at the surface is larger with the CKD continuum compared to the CAVIAR continuum by about 7%, with the net surface clear-sky short-wave flux with the CAVIAR model being 838 W m⁻². However, letting the system evolve, we find that such a difference is adjusted by small changes in temperature and humidity profiles and cloud patterns, leaving the basic atmospheric structure similar.

Based on these results, we are confident that the correlation between the humidity in the upper atmosphere and the amount of NIR incident radiation shown in Figure 3 is not sensitive to the choice of water vapor continuum model.