Atmospheric Dynamics on Terrestrial Planets: The Seasonal Response to Changes in Orbital, Rotational, and Radiative Timescales

Thousands of planets outside of the solar system have been detected in the last two decades. These planets vary in their mass, radius, orbital period, eccentricity, solar constant, and more. In this variety of planets, terrestrial planets within their host star habitable zone, are of particular interest due to their potential to harbor life on their surface. Planetary and atmospheric characteristics can change the planet habitability potential (e.g., Spiegel et al. 2009; Kopparapu et al. 2013, 2017; Vladilo et al. 2013; Armstrong et al. 2014), and will have a large impact on the planet's atmospheric dynamics (e.g., Ferreira et al. 2014; Kaspi & Showman 2015; Linsenmeier et al. 2015; Chemke & Kaspi 2017; Penn & Vallis 2018). Studying the atmospheric dynamics dependencies on planetary parameters, will not only improve the understanding of atmospheric physics but may give rise to possible future observables from which we can learn about the planet and its atmosphere characteristics.

This study focuses mainly on seasonal variations of climate. We focus on three timescales that define the climate's seasonal variability: the rotational, orbital, and radiative timescales. This is done through a series of simulations varying the rotation rate (Ω) and three other parameters that are closely related to the radiative forcing, the obliquity (γ), orbital period (ω), and atmospheric surface pressure (p_s). We examine the zonal mean climate, more specifically the zonal mean meridional circulation and surface temperature dependence on the parameters mentioned above. The mean meridional circulation is one of the main circulation features of a terrestrial atmospheres (Vallis 2017), which on Earth, strongly relates to the water cycle and on Titan to the methane cycle (e.g., Mitchell et al. 2006). On Earth, this meridional circulation in the tropics is dominated by the thermally driven Hadley cell up to the subtropics. Beyond that, the meridional circulation is driven by the turbulence forming the Ferrel cell (Vallis 2017). The zonal mean meridional circulation, is described using the zonal mean meridional streamfunction, $\psi =2\pi a\int \mathop{\,v}\limits^{\_}\cos \phi d\sigma$ where a is the planetary radius, $\mathop{\,v}\limits^{\_}$ is the zonal mean meridional wind, ϕ is latitude, and σ = p/p_s is the vertical pressure coordinate normalized by the surface pressure.

The atmospheric circulation dependence on the rotation rate is a result of the Coriolis acceleration importance in the momentum balance. Fast rotation rate results in a multicellular structure of the streamfunction together with multiple jets, while as the rotation rate is slowed down, the cell number decreases and the Hadley cell widens and strengthens (e.g., Walker & Schneider 2006; Chemke & Kaspi 2015a, 2015b; Kaspi & Showman 2015). Slowing down the rotation rate results also in a weaker meridional temperature gradient. This is due to the eddy scale dependence on rotation rate, where eddies become smaller as the rotation rate is increased (Walker & Schneider 2006; Kaspi & Showman 2015), resulting in a less efficient meridional heat transport that increases the meridional temperature gradient.

The other three parameters that this study focuses on are closely related to the radiative forcing. The obliquity (γ), determines the seasonality of the insolation. Nonzero obliquity introduces seasonality to the insolation, and increasing the obliquity toward 90°, shifts poleward the latitude of maximum insolation and increases the insolation meridional gradient (Figure 1). At obliquity 90°, for every degree the planet moves in its orbit the latitude of zero insolation (black contours in Figure 1) also moves a degree, which explains the linearity seen in the zero insolation line in Figure 1(c). Several studies used different types of general circulation models (GCMs) to study the circulation and habitability dependence on the obliquity. There is a general agreement between the models that seasonal variability increases by increasing the obliquity. Though most of these studies focus mainly on the effect of obliquity on the planetary habitability (Williams & Kasting 1997; Spiegel et al. 2009; Armstrong et al. 2014; Ferreira et al. 2014; Linsenmeier et al. 2015; Wang et al. 2016; Nowajewski et al. 2018). Mitchell et al. (2014) studied the seasonality effect on climate, using an idealized parameterization for the seasonality and radiative timescale. In contrast, in this study we examine the climate seasonality response to different physical parameters. This study also provides a systematic study of the seasonal cycle climate dynamical response to changes in the obliquity, and the obliquity effect in a larger parameter space. This study can be viewed as an expansion of Kaspi & Showman's (2015) perpetual equinox study to include the effects of seasonality.

Zoom In Zoom Out Reset image size

Introducing seasonality to the climate system requires taking into account different timescales. The orbital period is a natural timescale of the seasonal cycle, dictating the time over which radiative changes take place. When the orbital period is longer, radiative changes occur over more extended periods giving the atmosphere a longer time to adjust. As radiative changes take place over shorter periods (short orbital periods), a point where the radiative forcing is effectively the annual mean one (Figure 2) is reached. Note that although the orbital period is coupled to the distance from the host star and thus to the stellar flux reaching to the planet, in this study the stellar flux is kept constant while the orbital period varies. This is analogous to studying similar planets that orbit stars with different mass and thus different luminosity. Due to the simplified radiation scheme in our model, effects related to the spectral distribution of the insolation are not taken into account, although it can have some effect on the resulting climate (Goldot et al. 2015; Shields et al. 2016).

Zoom In Zoom Out Reset image size

Another important timescale is the atmospheric radiative timescale, meaning, the time that it takes the atmosphere to adjust to radiative changes. This timescale depends on the atmospheric mass per unit area, defined to be m = p_s/g (Hartmann 2015), where p_s is the surface pressure and g is the surface gravity. The radiative timescale is defined as

$$\begin{eqnarray}&&\tau ={\tau }_{\mathrm{IR}}\displaystyle \frac{{c}_{p}{p}_{s}}{4g{\sigma }_{b}{T}_{e}^{3}},\end{eqnarray} \tag{ 1 }$$

where τ_IR is the infrared optical depth of the atmosphere, c_p is the heat capacity of the air, σ_b is the Stefan–Boltzmann constant, and T_e is the radiative equilibrium temperature. The radiative timescale is linearly proportional to the atmospheric mass (Equation (1)), meaning that as the atmospheric mass is increased, the atmosphere needs a longer time to adjust to the changes in the radiative forcing. The radiative timescale has a strong dependence on the planetary equilibrium temperature T_e, yet because it is not an input parameter of the model and in order to avoid additional complexities as shown in Kaspi & Showman (2015), we focus on varying the radiative timescale by changing only the atmospheric mass.

Increasing the atmospheric mass expands the outer edge of the habitable zone due to the increase in surface temperatures with the atmospheric mass (e.g., Wordsworth et al. 2010, 2011; Vladilo et al. 2013). Increasing the atmospheric mass also results in flattening of the meridional temperature gradient and lowering of the troposphere height (Goldblatt et al. 2009; Kaspi & Showman 2015; Chemke et al. 2016; Chemke & Kaspi 2017). Several explanations for this increase in the surface temperature have been given. First, increasing the atmospheric mass results in enhancement of the greenhouse effect by pressure broadening of absorption lines overcoming the enhancement of Rayleigh scattering (Goldblatt et al. 2009). Second, by increasing the atmospheric mass, the adiabatic lapse rate is increased, that, in turn, warms the surface and lowers the tropopause (Goldblatt et al. 2009). A third explanation suggested by Chemke & Kaspi (2017) is that increasing the atmospheric mass increases the atmospheric heat capacity that in turn weakens the atmospheric radiative cooling. The decrease in radiative cooling, which is more pronounced in colder latitudes, results in warming of the surface and flattening of the meridional temperature gradient. Chemke & Kaspi (2017) also studied the Hadley circulation response to increasing the atmospheric mass, and found that the Hadley circulation narrows and weakens as the atmospheric mass is increased; this is attributed mainly to the lowering of the tropopause and the flattening of the temperature gradient with the increase in atmospheric mass.

The solar system terrestrial atmospheres of Venus, Earth, Mars, and Titan exhibit significant variability in their different planetary parameters (Table 1), which results in a variety of circulations. For example, Venus has a massive atmosphere, no significant obliquity, and a very slow rotation rate. It has a circulation in its lower atmosphere that is composed of two hemispherically symmetric equator to pole Hadley cells, with no seasonality, while its upper atmosphere is dominated by a day to night side circulation (e.g., Read 2013; Sánchez-Lavega et al. 2017). On Earth, there is a seasonal cycle, in which the Hadley cell transits from a hemispherically symmetric circulation to a strong, wide, winter cross-equatorial cell and a weak and narrow summer cell (e.g., Dima & Wallace 2003).

During Mars' seasonal cycle, its maximum surface temperature shifts from one pole to the other. This strong seasonality is attributed to Mars' thin atmosphere and rocky surface, that also explains Mars' large meridional temperature gradient (McCleese et al. 2010; Read et al. 2015). Titan, on the other hand, although having a long orbital period, has a high atmospheric mass and low surface temperature resulting in a long radiative timescale (Mitchell & Lora 2016). As a result, the surface temperature does not shift significantly off the equator during the seasonal cycle (Lora et al. 2015; Jennings et al. 2016). Also, in contrast to Mars, Titan has an all tropics climate (Mitchell et al. 2006), meaning that its meridional temperature gradient is weak (Jennings et al. 2009). The low meridional temperature gradient on Titan is attributed mainly to its slow rotation rate (Mitchell & Lora 2016; Hörst 2017). The solstice Hadley circulation on both planets is composed of a cross-equatorial cell, with air rising at midlatitudes of the summer hemisphere and descending at midlatitudes of the winter hemisphere (e.g., Lora et al. 2015; Read et al. 2015). Although the Hadley cell extent on both planets is similar, on Mars the extent of the circulation is mainly due to the poleward shift of the maximum surface temperature at solstice, and on Titan, it is mainly due to its slow rotation rate (Guendelman & Kaspi 2018).

Section 2 introduces the model used in this study, focusing on the temperature equation. Section 3 describes how the surface temperature depends on the different parameters. In order to illustrate the effects of the orbital period and atmospheric mass, the case of γ = 54° is studied in detail. Section 4 describes the zonal mean meridional circulation dependence on the study parameters and the dynamical interpretation of this dependence. We then suggest an empirical power-law fit to the model results for the latitudes of the Hadley circulation ascending and descending branches and its strength. Section 5 summarizes and concludes the results with a discussion on the implications of this study results on possible future observables and planetary habitability.

In order to study the climate sensitivity to the parameters mentioned above, we use an idealized, aquaplanet GCM (Frierson et al. 2006), based on the GFDL dynamical core (Anderson et al. 2004), similar to the model described in detail in Kaspi & Showman (2015). The model uses a two-stream gray radiation scheme, where the temperature field is determined by

$$\begin{eqnarray}&&\displaystyle \frac{{DT}}{{Dt}}-\displaystyle \frac{{R}_{d}{T}_{v}w}{{c}_{p}p}={Q}_{r}+{Q}_{c}+{Q}_{b},\end{eqnarray} \tag{ 2 }$$

with the material derivative given by $\tfrac{D}{{Dt}}=\tfrac{\partial }{\partial t}+{\boldsymbol{u}}\cdot {\boldsymbol{\nabla }}$, where ${\boldsymbol{u}}$ = (u, v, w) are the velocities in the longitudinal (λ), meridional (ϕ), and vertical (σ) directions. R_d is the dry gas constant of air, T_v is the virtual temperature, which is the temperature that an air parcel would have in a water vapor free air at a given pressure and density. Q_r, Q_c, and Q_b are the radiative, convective, and boundary layer heating per unit of mass, respectively.

The radiative heating, Q_r is calculated using

$$\begin{eqnarray}&&{Q}_{r}=\displaystyle \frac{g}{{c}_{p}}\displaystyle \frac{\partial }{\partial p}\left(U-D-{R}_{s}\right),\end{eqnarray} \tag{ 3 }$$

where U and D are the upward and downward longwave radiation, respectively, and R_s is the solar shortwave radiation. In order to include seasonality, we calculate R_s using

$$\begin{eqnarray}&&{R}_{s}={S}_{0}\cos \zeta {e}^{-{\tau }_{s}{\sigma }^{2}},\end{eqnarray} \tag{ 4 }$$

where, S₀ = 1360 Wm⁻² is the solar constant, τ_s is the parameter that controls the vertical absorption of the solar radiation, and ζ is the diurnal mean zenith angle, which is given by (Pierrehumbert 2010; Hartmann 2015)

$$\begin{eqnarray}&&\cos \zeta =\displaystyle \frac{h}{\pi }\left(\sin \phi \sin \delta +\displaystyle \frac{1}{h}\cos \phi \cos \delta \sin h\right),\end{eqnarray} \tag{ 5 }$$

where $\cos h=-\tan \phi \tan \delta$ is the daily average hour angle, where for perpetual night h = 0 and for perpetual day h = π. The declination angle is calculated using (Pierrehumbert 2010)

$$\begin{eqnarray}&&\sin \delta =-\sin \kappa \sin \gamma ,\end{eqnarray} \tag{ 6 }$$

where κ is the orbital position relative to the northern hemisphere autumnal equinox angle and γ is the obliquity.

Increasing the obliquity, shifts the maximum insolation poleward and increases the meridional insolation gradient (Figure 1). For short orbital periods, the yearly mean insolation can become important. The latitude of maximum yearly mean insolation flips from the equator to the poles around obliquity 54° (Figure 2(a)). Also, the yearly mean meridional insolation gradient is lower for higher obliquities, reaching a minimum at around γ = 54° (Figure 2(b)). This was shown to be the response of the yearly mean surface temperature (Nowajewski et al. 2018; Kang 2019).

We use a spectral horizontal resolution of T42 (2.8° × 2.8°) and 25 uneven vertical levels for all simulations. The parameters and their ranges are as follows:

• 1.  
  Rotation rate, ${\rm{\Omega }}=\left[\tfrac{1}{16},\tfrac{1}{8},\tfrac{1}{4},\tfrac{1}{2},1\right]$.
• 2.  
  Obliquity, $\gamma =\left[10^\circ ,20^\circ ,30^\circ ,40^\circ ,50^\circ ,60^\circ ,70^\circ ,80^\circ ,90^\circ \right]$.
• 3.  
  Atmospheric mass, ${p}_{s}=\left[\tfrac{1}{2},1,2,5,10,30\right]$.
• 4.  
  Orbital period, $\omega =\left[\tfrac{1}{8},\tfrac{1}{4},\tfrac{1}{2},1,2,4\right]$.

Note that in this study, the units for Ω, p_s, and ω are normalized to the corresponding Earth-like values, where for convenience we take ω = 1 to be 360 days instead of 365.

Also note that in this study, given this vast parameter range, we simulate six different subspaces. We select the subspaces by varying two parameters at a time, holding the other two with Earth-like values, except for the obliquity, that is held at 30°. The range chosen for the different parameters is discussed in Section 5.

Comparing the model with an Earth-like configuration and Earth observational reanalysis shows a general similarity between the model results and observations, with some differences. The model exhibits stronger seasonality, a stronger meridional streamfunction, and a colder climate compared with observations (Figure 3). Taking into account that the model is highly idealized, where, for example, clouds, land, or ice effects on the climate are neglected, the resulting climate with Earth-like parameters is satisfying. In addition, as the aim of this study is to examine the climate sensitivity to the planetary parameters to leading order, a perfect reconstruction of Earth's climate is not necessary. Note that in all relevant figures in this study, for convenience, the time coordinate is represented by the year fraction.

Zoom In Zoom Out Reset image size

Visualizing the climate dependence on a broad set of parameters becomes complex when the seasonal cycle is taken into account. A good form to visualize it in are figures such as Figure 4. This type of figure demonstrates how the climate, and its seasonality, change in response to changes in the parameters in a continuous form. The reader should approach this figure as one plot following the axis that represents the varied parameters. For orientation, we highlight in blue the simulation with Earth-like parameters.

Zoom In Zoom Out Reset image size

Increasing the obliquity shifts the maximum surface temperature toward the summer hemisphere pole, increases the meridional temperature gradient, and results in a sharper transition between summer and winter, where at high obliquities the transition seasons are less pronounced (Figure 4). This temperature response to changes in the obliquity is a result of the insolation dependence on the obliquity (Figure 1).

In long orbital periods, the surface temperature follows more closely the solar insolation pattern (black contours in Figure 4). As a result, the climate in long orbital periods exhibits strong seasonality, with maximum surface temperatures reaching to the summer pole (bottom row in Figure 4). Short orbital periods (top row in Figure 4), exhibit weak to no seasonality, and the surface temperature takes a similar shape to the yearly mean insolation. The surface temperature dependence on the orbital period is a result of the orbital period being the timescale in which radiative changes occur. This means that long orbital periods allow radiative changes to take place over longer periods, giving the atmosphere more time to adjust to these changes, resulting in a pronounced seasonality. Alternatively, when the orbital period is short, the atmosphere does not have enough time to adjust to the seasonal radiative changes, resulting in the solar forcing being effectively the yearly mean forcing.

Faster rotation rates result in a stronger meridional temperature gradient (Figure 5). This is explained by the eddy scale decreasing with faster rotation rate, which results in a less efficient meridional heat transport, that in turn increases the meridional temperature gradient (Walker & Schneider 2006; Kaspi & Showman 2015).

Zoom In Zoom Out Reset image size

By increasing the atmospheric mass the surface temperature increases and the meridional temperature gradient decreases (Figure 5), in agreement with previous studies (e.g., Goldblatt et al. 2009; Chemke et al. 2016; Chemke & Kaspi 2017). Using the temperature equation (Equation (2)), Chemke & Kaspi (2017) showed that the main component in the temperature equation that is strongly affected by the atmospheric mass is the radiative term (Equation (3)), which can be rewritten as

$$\begin{eqnarray}&&{Q}_{r}=\displaystyle \frac{1}{{c}_{p}}\displaystyle \frac{\partial }{\partial m}\left(U-D-{R}_{s}\right).\end{eqnarray} \tag{ 7 }$$

Writing Q_r in this form highlights its dependence on the atmospheric mass (m = p/g), or the atmospheric heat capacity, c_pdm. Equation (7) shows that an increase in the atmospheric heat capacity lowers the radiative cooling effect, which in turn results in a warmer surface climate and a flatter meridional temperature gradient (Chemke & Kaspi 2017).

Focusing on the atmospheric mass dependence (Figure 5) shows that the atmospheric mass has only a small effect on the surface temperature seasonality. More specifically, increasing the atmospheric mass has only a small effect on the latitude of maximum surface temperature and on the temperature temporal variations. The dominant effect of increasing the atmospheric mass, is, as mentioned, a warmer surface climate and flatter meridional temperature gradient. This is puzzling, as the radiative timescale of the atmosphere strongly depends on the atmospheric mass (Equation (1)).

In order to understand why the atmospheric mass has little effect on seasonality, it is interesting to take a closer look at the γ ≈ 54° case. This is a special case, as it is a turning point in the yearly mean forcing from maximum radiation at the equator to maximum radiation at the poles. Also, it is a minimum point in the yearly mean insolation gradient (Figure 2). As a result of this unique obliquity value, in a short orbital period, the expected surface temperature will have a relatively small meridional gradient as the effective forcing in a short orbital period resembles the yearly mean forcing (Figure 4).

A priori assuming all other parameters stay constant, decreasing the orbital period or increasing the atmospheric mass should have a similar effect on the seasonality as both relate to the adjustment timescale of the atmosphere. For this reason, comparing the effects of increasing the atmospheric mass to decreasing the orbital period, can help to explain why increasing the atmospheric mass has a small effect on seasonality. A quantitative comparison shows that increasing the atmospheric mass mainly warms the surface, with significant warming of the cold latitudes that flatten the meridional temperature gradient (Figures 6(b), 7(a)). In contrast, a shorter orbital period tones down the seasonality, by symmetrically decreasing the temperature in warm latitudes and increasing them in cold latitudes (Figures 6(c), 7(b)). This asymmetric effect of atmospheric mass on the surface temperature can be explained by the radiative cooling dependence on the atmospheric mass (Equation (7)). As mentioned, increasing the atmospheric mass decreases the radiative cooling, which is more effective at cold latitudes. As a result, there is more warming of cold latitudes than cooling at warm latitudes as the atmospheric mass is increased. This decrease in radiative cooling with the atmospheric mass has a more pronounced effect that diminishes the seasonal effect of the atmospheric mass.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The coupled effect of a short orbital period together with high atmospheric mass results in a climate with close to a latitudinally uniform surface temperature throughout the seasonal cycle for a planet with γ = 54°. More specifically, the difference between the yearly maximum and minimum surface temperature is ∼6° K (Figure 6(d)). This suggests that not all high obliquity planets will experience a strong seasonal cycle, with fast and strong transitions between warm and cold temperatures, and the possibility that a high obliquity planet will have an equable climate that is habitable exists. It is important to note that diurnal cycle effects are neglected, and although they might become important for short orbital periods, they are out of the scope of this study. Salameh et al. (2018) showed that diurnal effects on the surface temperature become important when the ratio between the orbital period (ω) and the rotation period (P_Ω) is ω/P_Ω ≈ 10. In this case, ω/P_Ω = 45, so neglecting diurnal effects is justified.

We parameterize the meridional surface temperature difference by

$$\begin{eqnarray}&&{{\rm{\Delta }}}_{H}=\displaystyle \frac{\max ({T}_{s})-\min ({T}_{s})}{\mathrm{mean}({T}_{s})},\end{eqnarray} \tag{ 8 }$$

where T_s is the surface temperature in the southern hemisphere summer. Δ_H is a measure for the meridional surface temperature gradient and represents the value of the temperature difference between the warm and cold latitudes in the extreme season. Generally, Δ_H increases as the rotation rate, obliquity, and orbital period are increased, or when the atmospheric mass is decreased (Figure 8). However, a clear exception is found for short orbital periods (ω < 1/4), were a minimum in Δ_H is found around obliquity ∼50° (Figure 8(c)). This minimum corresponds to the minimum in the yearly mean insolation meridional difference around obliquity ∼54° (Figure 2(b)), indicating that in short orbital periods, the climate becomes similar to the yearly mean climate.

Zoom In Zoom Out Reset image size

4.1. Simulation Results

The Hadley circulation is a thermally driven component of the zonal mean meridional circulation, with air rising at warm latitudes and descending at colder ones. As air rises, it condenses, creating a zone highly populated with clouds and intense precipitation. On Earth, it is called the intertropical convergence zone (ITCZ). The width, namely the latitude of ascending and descending branches and the seasonality of the Hadley circulation varies between the solar system planets. Venus and Titan are extreme examples, where the mean meridional circulation is composed mainly of the Hadley circulation, which reaches high latitudes (e.g., Roe 2012; Sánchez-Lavega et al. 2017). Numerous theories were suggested for predicting the positions of the ascending and descending branches, and a good summary of the different theories can be found in Faulk et al. (2017).

The Hadley circulation is closely related to the water cycle on Earth, and the methane cycle on Titan (Mitchell et al. 2006; Schneider et al. 2012). As mentioned, the ascending branch of the circulation correlates to the region of intense precipitation and is a region highly populated with clouds. In contrast, the descending branch region is clear of clouds and generally associated with low precipitation and desert areas. The cloud and desert zones might be good candidates for future observables through their very different albedos. Determining the detection methods for the circulation effects is out of the scope of this study.

The seasonality of the meridional circulation and its strength increase with increasing obliquity (Figure 9). During solstice of strong seasonal cases, the zonal mean meridional streamfunction is composed of a strong cross-equatorial winter Hadley cell, with no summer cell and weak transition seasons. This Hadley cell response follows the surface temperature response (Figures 4 and 9). However, for Earth-like rotation rate cases, even in cases where the maximum surface temperature is at the pole (fourth row in Figure 4), the ascending branch does not reach the pole (bottom row of Figure 9). This result is in agreement with Faulk et al. (2017), who showed that even in an eternal solstice case (for Earth-like forcing), where the maximum surface temperature is at the pole, the ITCZ stays at low latitudes, and does not shift to the position of maximum surface temperature. However, by slowing down the rotation rate, a wider and stronger circulation emerges, extending from the summer hemisphere pole to around latitude ∼60° in the winter hemisphere (Figure 9). This result indicates that the rotation rate is a limiting factor for the width of the circulation, a result that is consistent with the axisymmetric theory (Guendelman & Kaspi 2018; Hill et al. 2019; Singh 2019).

Zoom In Zoom Out Reset image size

The circulation narrows and weakens in response to an increase in atmospheric mass, and an opposite response occurs when the orbital period is increased (Figure 10). Similar to the response of the circulation to the obliquity changes, the circulation response to the orbital period and atmospheric mass changes also correlates with the surface temperature response (Figures 4, 5, and 10). More specifically, increase in seasonality and temperature gradient (where both increases by decreasing the atmospheric mass or increasing the obliquity and orbital period) results in a wider and stronger Hadley circulation (Figures 4, 5, and 8–10).

Zoom In Zoom Out Reset image size

4.2. Theoretical Arguments

Several theories for the ITCZ position and the Hadley cell ascending branch latitude have been suggested. Among them, Neelin & Held (1987) theorized that air will ascend close to the latitude of maximum surface temperature or around the maximum of the low level moist static energy (MSE),

$$\begin{eqnarray}&&m={L}_{e}q+{gz}+{c}_{p}T,\end{eqnarray} \tag{ 9 }$$

where L_e is the latent heat of vaporization, q is the specific humidity, and z is the geopotential height. Kang et al. (2008) suggested that the position of the maximum precipitation correlates with the energy flux equator (EFE), which is the latitude where the vertically integrated MSE flux ($\overline{{mv}}$) is zero.

The axisymmetric theory (Held & Hou 1980; Lindzen & Hou 1988; Caballero et al. 2008; Guendelman & Kaspi 2018) neglects any eddy contribution to the momentum balance. Assuming angular momentum conservation at the top branch of the circulation, the angular momentum conserving wind is

$$\begin{eqnarray}&&{u}_{m}={\rm{\Omega }}a\displaystyle \frac{{\cos }^{2}{\phi }_{1}-{\cos }^{2}\phi }{\cos \phi },\end{eqnarray} \tag{ 10 }$$

where a is the planetary radius and ϕ₁ is the latitude of the ascending branch (see a more detailed derivation in the Appendix). Assuming also that gradient-wind balance holds, meaning

$$\begin{eqnarray}&&{fu}+\displaystyle \frac{{u}^{2}\tan \phi }{a}=-\displaystyle \frac{1}{a}\displaystyle \frac{\partial {\rm{\Phi }}}{\partial \phi },\end{eqnarray} \tag{ 11 }$$

with Φ = p/ρ₀, where p is pressure and ρ₀ is a reference density used in the Boussinesq approximation (Lindzen & Hou 1988). Using Equations (10) and (11) we can derive the vertical averaged potential temperature field in gradient balance with angular momentum conserving wind (Lindzen & Hou 1988; Guendelman & Kaspi 2018)

$$\begin{eqnarray}&&\displaystyle \frac{\overline{\theta }(\phi )-\overline{\theta }({\phi }_{1})}{{\theta }_{0}}=-\displaystyle \frac{{\delta }_{H}}{{R}_{t}}\displaystyle \frac{{\left({\sin }^{2}\phi -{\sin }^{2}{\phi }_{1}\right)}^{2}}{{\cos }^{2}\phi },\end{eqnarray} \tag{ 12 }$$

where $\overline{\theta }$ , the potential temperature, is the temperature that an air parcel would have if it were brought adiabatically to some reference pressure (Hartmann 2015),

$$\begin{eqnarray}&&{R}_{t}=\displaystyle \frac{8\pi \rho {GH}{\delta }_{H}}{3{{\rm{\Omega }}}^{2}a},\end{eqnarray} \tag{ 13 }$$

is the thermal Rossby number (Guendelman & Kaspi 2018), where ρ is the planet's mean density, G is the universal gravitational constant, and δ_H is the meridional fractional change of the radiative equilibrium temperature that is given by

$$\begin{eqnarray}&&\displaystyle \frac{{\overline{\theta }}_{e}}{{\theta }_{0}}=1+\displaystyle \frac{{\delta }_{H}}{3}(1-3{\left(\sin \phi -\sin {\phi }_{0}\right)}^{2}),\end{eqnarray} \tag{ 14 }$$

where ϕ₀ is the latitude of maximum θ_e. In order to find the predicted ascending and descending branches of the circulation, we need to assume that the circulation is closed and temperature is continuous at the edge of the cells, e.g.,

$$\begin{eqnarray}&&{\int }_{{\phi }_{1}}^{{\phi }_{j}}(\overline{\theta }-{\overline{\theta }}_{e})\cos \phi d\phi =0,\end{eqnarray} \tag{ 15 }$$

$$\begin{eqnarray}&&\overline{\theta }({\phi }_{j})={\overline{\theta }}_{e}({\phi }_{j}),\end{eqnarray} \tag{ 16 }$$

where j = w, s represents the position for the winter and summer descending branch. If we assume that the circulation is composed of only one pole-to-pole cell we can use

$$\begin{eqnarray}&&{\int }_{-\pi /2}^{\pi /2}(\overline{\theta }-{\overline{\theta }}_{e})\cos \phi d\phi =0,\end{eqnarray} \tag{ 17 }$$

$$\begin{eqnarray}&&\overline{\theta }(\pi /2)={\overline{\theta }}_{e}(\pi /2),\end{eqnarray} \tag{ 18 }$$

to give

$$\begin{eqnarray}&&{R}_{t}=\displaystyle \frac{1}{3\sin {\phi }_{0}-1}.\end{eqnarray} \tag{ 19 }$$

We can translate Equation (19) to conditions on R_t and ϕ₀ for having a pole-to-pole circulation. First, because R_t > 0 we need ${\phi }_{0}\gt \arcsin (1/3)$, second, ϕ₀ = 90° represents the minimal value of R_t to have a pole-to-pole circulation, meaning that in order to have a pole-to-pole circulation we need R_t ≥ 0.5. Solving numerically for different values of R_t and ϕ₀ we see that this condition holds (Figure 11).

Zoom In Zoom Out Reset image size

This limit on R_t suggests that planets with small R_t (usually fast rotating planets, such as Earth and Mars), will not be able to reach a state with a pole-to-pole Hadley circulation. In contrast, planets with high R_t such as Titan, need a relatively small ϕ₀ in order to reach a pole-to-pole circulation (Guendelman & Kaspi 2018).

One advantage of the axisymmetric theory, compared to other theories, is that the axisymmetric theory gives a prediction for both edges of the circulation, ascending and descending, using one set of arguments. However, as mentioned, this theory neglects any eddy contribution, which may play an important role (Walker & Schneider 2005, 2006). Other theories that do take into account eddy contribution have been suggested (e.g., Held 2000; Korty & Schneider 2008; Kang et al. 2008; Levine & Schneider 2015).

4.3. Matching Theory to Simulations

An alignment of the streamfunction and angular momentum, suggests that the streamfunction conserves its zonal mean angular momentum and the circulation is less influenced by eddies and is more axisymmetric. As a result, cases where the meridional streamfunction follows the zonal mean angular momentum contours, the axisymmetric theory becomes more appropriate. Indeed, in some cases, mainly cases with strong seasonality, an alignment between the streamfunction and zonal mean angular momentum is noticed (Figure 12).

Zoom In Zoom Out Reset image size

The streamfunction becomes aligned with the angular momentum contours both by slowing down the planetary rotation rate and by increasing seasonality, namely, having a strong cross-equatorial circulation (Figure 12). The seasonal alignment is related to a previously studied regime transition from an eddy mediated equinox cell to an axisymmetric solstice cell (Bordoni & Schneider 2008, 2010; Merlis et al. 2013; Geen et al. 2018). Both alignments are a result of the weakening of eddy momentum flux convergence, where the seasonal alignment due to shielding by the upper level easterlies and the rotation rate alignment is due to a global weakening of the eddy momentum flux convergence as the rotation rate is slowed down (Faulk et al. 2017).

These alignments indicate that axisymmetric arguments (e.g., Lindzen & Hou 1988) can be used to explain the Hadley cell dependence on the different parameters, mainly the circulation dependence on the thermal Rossby number. Note that although this theory is very idealized, it was found to give a good approximation for Earth Hadley cell width (Held & Hou 1980), for the position of the Hadley cell ascending branch on Titan and Mars (Guendelman & Kaspi 2018), and can give a physical interpretation to changes in the width and strength of the circulation as a response to changes in planetary parameters (Kaspi & Showman 2015; Chemke & Kaspi 2017). As mentioned above, according to the axisymmetric theory, the circulation depends on the thermal Rossby number (R_t) and the latitude of maximum warming (ϕ₀). The theory predicts that the circulation will become wider and stronger as the troposphere height (H), the meridional radiative equilibrium temperature gradient (δ_H) and ϕ₀ are increased or by slowing down the rotation rate (e.g., Guendelman & Kaspi 2018).

Different from the axisymmetric theory, where the radiative equilibrium temperature, more specifically, δ_H, ϕ₀, and the troposphere height H are predetermined (Lindzen & Hou 1988), in our model these parameters are not predetermined, so a parameterization using the model results for these parameters is needed. The surface meridional temperature difference Δ_H (Equation (8)) is used to parameterize the meridional radiative equilibrium temperature gradient (δ_H). For the latitude of maximum radiative equilibrium temperature, we choose the latitude of maximum surface temperature. For the troposphere height parameterization, the height of the circulation is used (Walker & Schneider 2006), being the height where the circulation is 10% of its maximum. This parameterization for the tropopause height seems to be a more appropriate one than the WMO lapse rate definition (Reichler et al. 2003), as the focus of this study is on the circulation itself. Note that the choice of parameterizations, although with reason, is not unique, and other physically consistent parameterizations can be used.

Comparing between the model results for the ascending branch (Figure 13, shading), taken to be the latitude where the streamfunction gets to 5% of its maximum value (Walker & Schneider 2006; Faulk et al. 2017), with the calculated prediction from the axisymmetric theory (dashed contours in Figures 13 and 15(b)), shows a general agreement with some numerical misfit. Although the correlation between the axisymmetric prediction and the model results is not perfect, its existence, together with the general agreement between the contour shapes (Figure 13) suggests that the use of axisymmetric arguments are appropriate to leading order. Yet there is a definite numerical misfit, which can be a result of neglecting the eddy contribution. Another possible explanation for this misfit is a problem in the parameterization used to calculate the theoretical prediction, such as the use of the surface temperature to parameterize δ_H. The surface temperature is the final result of the dynamical and radiative processes, so using it to parameterize δ_H that in the theory context represents only the radiative forcing is problematic, so a separation between the two contributions might be in need of a more appropriate parameterization.

Zoom In Zoom Out Reset image size

Figure 13 summarizes the effect of the different study parameters on the Hadley cell ascending branch latitude. Increasing the orbital period and the obliquity or decreasing the rotation rate and atmospheric mass, results in the poleward shift of the ascending branch. These trends of the ascending branch position with the atmospheric mass, orbital period, and obliquity, correlate well with the trends of H and Δ_H. More specifically, increasing the obliquity and orbital period or decreasing the atmospheric mass results in an increase in H and Δ_H (Figures 8 and 14). The increase in H and Δ_H correlates with the poleward shift of the ascending branch, in agreement with axisymmetric arguments (Held & Hou 1980). However, increasing the rotation rate increases both the temperature gradient and tropopause height (Figures 8 and 14), and in contrast to the axisymmetric arguments, the circulation narrows. Nonetheless, this is not a violation of the axisymmetric theory; first, because the rotation rate is itself an important parameter in this theory, where by increasing the rotation rate the circulation becomes narrower and this effect dominants over the effect of the other parameters (Guendelman & Kaspi 2018). Second, the increase in the temperature gradient with faster rotation rate is associated with the decrease of the eddy scale with the rotation rate (e.g., Kaspi & Showman 2015), which is a dynamical eddy effect; thus it does not contradict axisymmetric arguments.

Zoom In Zoom Out Reset image size

In agreement with Faulk et al. (2017) this study suggests that the only form in which the ascending branch can reach the summer pole is by slowing down the rotation rate (Figure 13). This strong dependence of the circulation on the rotation rate comes mainly from the dominance of the Coriolis force in the momentum balance. A similar strong dependence on the rotation rate also arises in the axisymmetric theory, where it predicts that the width of the circulation is a function of thermal Rossby number that depends strongly on the rotation rate (Guendelman & Kaspi 2018).

Comparing the model results with the EFE (e.g., Kang et al. 2008; Bischoff & Schneider 2014; Schneider et al. 2014), taken here to be the the latitude where the vertically integrated MSE flux ($\overline{{mv}}$) reaches 5% of its maximum, although somewhat dispersed, it shows a good correlation with the position of the ascending branch (Figure 15(b)). These results are in agreement with Wei & Bordoni (2018) who showed that there is a lag between the EFE and the ITCZ seasonal cycle; however, when taking seasonal time average, this lag disappears.

Zoom In Zoom Out Reset image size

Another predictor is the latitude of low level maximum MSE (e.g., Neelin & Held 1987), ϕ_MSE, where it correlates well with the ascending branch as long as ϕ_MSE stays around low latitudes. However, there are points where ϕ_MSE is at the pole, but the ascending branch stays at low latitudes (Figure 15(c)). These points correspond to cases with relative fast rotation rates (close to Earth-like rotation rate), where the maximum surface temperature is at the pole, but ϕ₁ does not reach the pole, as the fast rotation rate acts to limit the extension of the circulation (Faulk et al. 2017; Guendelman & Kaspi 2018).

4.4. Power-law Fit for the Hadley Circulation Response

In order to get a relation between the circulation response to the studied parameters, we use a power-law fit of the form $\propto {p}_{s}^{a}{\gamma }^{b}{{\rm{\Omega }}}^{c}{\omega }^{d}$ for the latitude of the ascending branch, descending branch, and the circulation strength. The power-law coefficients are determined by minimizing the sum of squared residuals, e.g., the sum of squared differences between the model function and the fitted power-law function. For the fit process, all saturated points were left out.

The fit for the ascending branch, Φ₁ is

$$\begin{eqnarray}&&{{\rm{\Phi }}}_{1}\propto {p}_{s}^{-0.13}{\gamma }^{0.71}{{\rm{\Omega }}}^{-0.62}{\omega }^{0.28},\end{eqnarray} \tag{ 20 }$$

with R² = 0.87.

Although possible correlations between the different parameters are ignored, this power-law fit does reasonably well as a first-order approximation for both the ascending and descending branches (Figures 13 and 16). The fit for the descending branch is

$$\begin{eqnarray}&&{{\rm{\Phi }}}_{w}\propto {p}_{s}^{-0.13}{\gamma }^{0.14}{{\rm{\Omega }}}^{-0.44}{\omega }^{0.13},\end{eqnarray} \tag{ 21 }$$

with R² = 0.93.

Zoom In Zoom Out Reset image size

Another important aspect is the meridional circulation strength (Figure 17), for which the best fit is

$$\begin{eqnarray}&&{{\rm{\Psi }}}_{M}\propto {p}_{s}^{-0.43}{\gamma }^{0.27}{{\rm{\Omega }}}^{-0.39}{\omega }^{0.30},\end{eqnarray} \tag{ 22 }$$

with R² = 0.87.

Zoom In Zoom Out Reset image size

The first step in studying this empirical fit is to examine the power sign of each parameter. The signs show that the circulation becomes stronger and wider as we increase the obliquity and orbital period (positive sign), and also by decreasing the atmospheric mass and rotation rate (negative sign), in agreement with the results shown in this study.

The widening and strengthening of the circulation with the obliquity and orbital period correlates with the increase in the seasonality and meridional temperature gradient that results by increasing these parameters (Figures 4 and 8). Although there is also some contribution from the increase in the tropopause height with these parameters (Figure 14), the tropopause height increase is not as substantial as the increase in the temperature gradient.

Comparing the strength of the circulation and latitude of ascending and descending branch powers of p_s shows a clear difference, where the atmospheric mass strongly affects the circulation strength but has little effect on the width of the circulation. The reduction in the meridional temperature difference (Figure 8) and of the troposphere height (Figure 14) together with the small to no influence on the maximum surface temperature position (Figure 5) can explain the difference between the dependence of the circulation strength and width on the atmospheric mass. Our results are consistent with Chemke & Kaspi (2017), which found that in the perpetual equinox case, ${{\rm{\Phi }}}_{w}\propto {p}_{s}^{-0.16}$ and ${{\rm{\Psi }}}_{m}\propto {p}_{s}^{0.41}$.

Previous studies of the Hadley cell dependence on the rotation rate suggested some empirical and theoretical power laws. Caballero et al. (2008), using a small angle approximation showed that the axisymmetric theory predicts that Φ₁ ∝ Ω^−2/3, which is close to the result presented here for the ascending branch (Equation (20)). However, Caballero et al. (2008) gave the same scaling for Φ_w which is different from the result here (Equation (21)), indicating, that the axisymmetric arguments cannot help to fully explain the descending branch position in our simulations and there is a need for other processes. This conclusion is strengthened by the fact that the descending branch in our simulations does not go further than latitude ∼60°, which is not a feature of the axisymmetric theory. Faulk et al. (2017) studied the ITCZ position in the seasonal case, where the insolation was kept Earth-like, and only the rotation rate was changed and found that the position of the ITCZ is proportional to Ω^−0.63, in agreement with these results. Other studies gave different theoretical and empirical fits; however, those studies focused mainly on the perpetual equinox case (e.g., Held & Hou 1980; Walker & Schneider 2006).

The ability to detect and characterize terrestrial exoplanets has improved in the last years, with more planets detected and future missions planned to detect and characterize terrestrial planets and their atmospheres. The planets span over a wide range of orbital configurations, planetary parameters, and atmospheric characteristics. Studying the circulation dependence on different planetary parameters is important in order to understand the climate on these planets and can also give rise to possible future observables.

This study focuses on the climate dependence on four parameters: obliquity, orbital period, and atmospheric mass, which strongly relate to the radiative forcing, together with the rotation rate, that has significant dynamical importance. In order to study the climate dependence on these four parameters, we use an idealized GCM, where the simplifications in this model allow us to separate the different dynamical effects of each parameter without introducing other complexities such as land, ice, or cloud feedback.

The focus is on the zonal mean circulation, specifically the surface temperature and the Hadley circulation. Understanding the surface temperature response to changes in different planetary parameters is important in order to assess the planet's habitability potential. For example, a highly seasonal climate with large and fast transitions between cold to warm temperature reduces the habitability potential of a planet, even if the planet's mean surface temperature is in the habitable range, and surface water can be present. Also, it is essential to understand the water cycle response to changes in the planetary parameters; therefore, we study the Hadley circulation response. The Hadley circulation will also dictate ground features of the planet, for example, on Earth the distinction between the wet tropical regions, located around the ascending branch of the Hadley cell and the deserts located at the descending branch of the Hadley cell. Also on Titan methane lakes around the pole and ice deserts around the equator and midlatitudes can be related to the Hadley circulation (e.g., Hayes et al. 2008; Aharonson et al. 2009).

Nonzero obliquity introduces seasonality to the solar insolation, and as the obliquity increases toward 90°, the seasonality strengthens and the insolation meridional gradient increases (Figure 1). When introducing seasonality to the climate system, different timescales that relate to the radiative forcing become relevant. The natural seasonal timescale is the orbital period, which is the time over which radiative changes take place. Short orbital period means that the radiative changes take place over a short amount of time and for a short enough orbital period, the resulting surface temperature is closer to the yearly mean climate. Long orbital periods result in a longer time for radiative changes to take place. As a result, the surface temperature in long orbital periods will have a similar structure to the insolation (Figure 4).

The atmospheric radiative adjustment timescale, the time that the atmosphere needs to adjust to changes in the radiative forcing, depends on the atmospheric mass (Equation (1)). However, increasing the atmospheric mass also decreases the radiative cooling, which is a dominant effect on the temperature, and it masks over the seasonal effect of atmospheric mass, and as a result, the seasonality is weakly affected by changes in the atmospheric mass (Figures 5 and 7). That being said, we can still see some effect of the atmospheric mass on the seasonality, more specifically, the equatorward shift of the maximum surface temperature for high enough atmospheric mass (Figure 5). The γ ≈ 54° case is an interesting one, mainly because it exhibits a minimum in the yearly mean insolation (Figure 2). As a result of this minimum, in a short orbital period and a high atmospheric mass planet, the climate becomes equable, with no abrupt seasonal transitions (Figure 6(d)), this can increase the planets' habitability potential.

The response of the circulation to changes in obliquity, orbital period, and atmospheric mass follow the response of the temperature. Meaning that by increasing the obliquity and orbital period, the circulation becomes wider and stronger (Figures 9 and 10) due to the increase in seasonality and the meridional temperature gradient (Figure 8). Decreasing the atmospheric mass also results in a wider and stronger circulation (Figure 10), attributed mainly to the increase in the meridional temperature gradient and the increase in the troposphere height (Chemke & Kaspi 2017). The effect of the atmospheric mass on the circulation strength is stronger than the effect on its width (Equations (20)–(22)). This also correlates to the temperature response, as mentioned, the atmospheric mass does not contribute significantly to seasonality changes, but does change the meridional temperature gradient (Figure 8) and tropopause height (Figure 14).

It is important to note that even in a very seasonal case, where the maximum surface temperature is at the summer pole, the circulation does not extend to the pole if the rotation rate is Earth-like (Figure 9), and only happens when the rotation rate is reduced (Guendelman & Kaspi 2018). Slowing down the rotation rate also flattens the temperature gradient due to the higher heat transport efficiency of larger eddies. In addition, lowering the rotation rate widens the region where the weak temperature gradient limit holds, i.e., where the Coriolis term is less dominant, and gravity waves are acting to smooth the temperature gradient, resulting in wider tropics (Sobel et al. 2001; Raymond & Zeng 2005). The widening of the circulation by lowering the rotation rate can be attributed to the decrease in the eddy momentum flux convergence with the rotation rate (Faulk et al. 2017).

Summarizing the circulation dependence on the different parameters, we suggest an empirical power-law fit for the ascending and descending branches and the circulation strength (Equations (20)–(22), respectively). Although this fit is for a specific model in a specific configuration, it gives a first-order approximation for the dependence of the circulation on these four parameters.

Although theoretically the range of the orbital period, rotation rate, and atmospheric mass can be even larger, the range that this study covers represents the different seasonal regimes of the diurnal mean climate. The chosen values for the orbital period covers climates ranging from close to an annual mean climate for the short orbital periods, to a strong seasonal climate where the surface temperature follows the insolation closely (Figure 4). The atmospheric mass range shows both the strong radiative cooling and the relatively weak radiative timescale effects (Figure 5). For the rotation rate, we mainly cover the range for Earth-like and slower rotation rate, which covers the circulation response from a narrow, hemispherically symmetric Hadley circulation to a wide cross-equatorial circulation (Figure 9). Faster rotation rates are discarded as it enters into a different dynamical regime (Chemke & Kaspi 2015a, 2015b; Kaspi & Showman 2015).

In future observations, with some spatial resolution or other ways of inferring cloud cover, regions highly populated with clouds can possibly be associated with the ascending branch of the Hadley cell, and a sector that is cloud depleted can be related to the descending branch of the Hadley cell. Detecting such temporal changes in these cloud regions' latitude can give an indication about the planetary orbital and atmospheric characteristics.

We thank Rei Chemke for his help in the model configuration. We also appreciate the helpful comments from the anonymous reviewer. We acknowledge support from the Israeli Science foundation (grant 1819/16) and the Helen Kimmel Center for Planetary Science at the Weizmann Institute of Science.

Following the derivation in Lindzen & Hou (1988) closely, we start with the zonal mean angular momentum per unit mass, defined to be

$$\begin{eqnarray}&&M=({\rm{\Omega }}{a}^{2}\cos \phi +{ua})\cos \phi ,\end{eqnarray} \tag{ 23 }$$

where Ω is the planetary rotation rate, a is the planetary radius, ϕ is latitude, and u is the zonal mean zonal wind. The angular momentum of air that starts at rest from latitude ϕ₁, the Hadley circulation ascending branch, is

$$\begin{eqnarray}&&M={\rm{\Omega }}{a}^{2}{\cos }^{2}{\phi }_{1}.\end{eqnarray} \tag{ 24 }$$

Assuming that the upper branch of the Hadley circulation conserves its angular momentum, we can write an expression for the angular momentum conserving wind,

$$\begin{eqnarray}&&{u}_{m}={\rm{\Omega }}a\displaystyle \frac{{\cos }^{2}{\phi }_{1}-{\cos }^{2}\phi }{\cos \phi }.\end{eqnarray} \tag{ 25 }$$

Assuming also that the zonal flow is in gradient-wind balance

$$\begin{eqnarray}&&{fu}+\displaystyle \frac{{u}^{2}\tan \phi }{a}=-\displaystyle \frac{1}{a}\displaystyle \frac{\partial {\rm{\Phi }}}{\partial \phi },\end{eqnarray} \tag{ 26 }$$

where Φ = p/ρ₀, we can evaluate this balance (Equation (26)) at heights H and 0 and subtract the two to give

$$\begin{eqnarray}&&{fu}(H)+\displaystyle \frac{u{\left(H\right)}^{2}\tan \phi }{a}=-\displaystyle \frac{1}{{\rho }_{0}a}\displaystyle \frac{\partial }{\partial \phi }\left[{\rm{\Phi }}(H)-{\rm{\Phi }}(0)\right],\end{eqnarray} \tag{ 27 }$$

where we assume that $| u(0)| \ll | u(H)|$. Assuming hydrostatic balance

$$\begin{eqnarray}&&\displaystyle \frac{\partial {\rm{\Phi }}}{\partial z}=\displaystyle \frac{g}{{\theta }_{0}}\theta ,\end{eqnarray} \tag{ 28 }$$

where θ is the potential temperature and θ₀ is some reference potential temperature, and integrating the hydrostatic balance (Equation (28)) with respect to z gives

$$\begin{eqnarray}&&\left[{\rm{\Phi }}(H)-{\rm{\Phi }}(0)\right]=\displaystyle \frac{g}{{\theta }_{0}}\overline{\theta },\end{eqnarray} \tag{ 29 }$$

where $\overline{\theta }={\int }_{0}^{H}\theta {dz}/H$. Substituting u(H) = u_m and Equation (29) into Equation (27) gives

$$\begin{eqnarray}&&\displaystyle \frac{1}{{\theta }_{0}}\displaystyle \frac{\partial \overline{\theta }}{\partial \phi }=\displaystyle \frac{{{\rm{\Omega }}}^{2}{a}^{2}}{{gH}}\left(\displaystyle \frac{\tan \phi }{{\cos }^{2}\phi }{\cos }^{4}{\phi }_{1}-\sin \phi \cos \phi \right).\end{eqnarray} \tag{ 30 }$$

Integrating Equation (30) with respect to ϕ yields

$$\begin{eqnarray}&&\displaystyle \frac{\overline{\theta }(\phi )-\overline{\theta }({\phi }_{1})}{{\theta }_{0}}=-\displaystyle \frac{{{\rm{\Omega }}}^{2}{a}^{2}}{2{gH}}\displaystyle \frac{{\left({\sin }^{2}\phi -{\sin }^{2}{\phi }_{1}\right)}^{2}}{{\cos }^{2}\phi },\end{eqnarray} \tag{ 31 }$$

which is similar to Equation (12). In order to get the exact expression as in Equation (12), one needs to use g = 4πGρa/3 where ρ is the planet's mean density and G is the universal gravitational constant and the definition for the thermal Rossby number

$$\begin{eqnarray}&&{R}_{t}=\displaystyle \frac{8\pi \rho {GH}{\delta }_{H}}{3{{\rm{\Omega }}}^{2}a},\end{eqnarray} \tag{ 32 }$$

where δ_h is meridional fractional change of the radiative equilibrium (Equation (33)). Substituting R_t in Equation (31) gives Equation (12) exactly.

To close the set of equations, we parameterize the thermal forcing using Newtonian relaxation to an equilibrium potential temperature of the form

$$\begin{eqnarray}&&\displaystyle \frac{{\overline{\theta }}_{e}}{{\theta }_{0}}=1+\displaystyle \frac{{\delta }_{H}}{3}(1-3{\left(\sin \phi -\sin {\phi }_{0}\right)}^{2}),\end{eqnarray} \tag{ 33 }$$

with ϕ₀ being the latitude of maximum ${\overline{\theta }}_{e}$. In order to find the predicted latitude of ascending and descending branches of the circulation we need to assume that the circulation is thermally closed and that the temperature is continuous at the edge of the cells, e.g.,

$$\begin{eqnarray}&&{\int }_{{\phi }_{1}}^{{\phi }_{j}}(\overline{\theta }-{\overline{\theta }}_{e})\cos \phi d\phi =0,\end{eqnarray} \tag{ 34 }$$

$$\begin{eqnarray}&&\overline{\theta }({\phi }_{j})={\overline{\theta }}_{e}({\phi }_{j}),\end{eqnarray} \tag{ 35 }$$

where j = w, s represents the position for the winter and summer descending branch respectively.

Of particular interest is the pole-to-pole circulation case. In this scenario Equations (34) and (35) translate to

$$\begin{eqnarray}&&{\int }_{-\pi /2}^{\pi /2}(\overline{\theta }-{\overline{\theta }}_{e})\cos \phi d\phi =0,\end{eqnarray} \tag{ 36 }$$

$$\begin{eqnarray}&&\overline{\theta }(\pi /2)={\overline{\theta }}_{e}(\pi /2),\end{eqnarray} \tag{ 37 }$$

writing these equations explicitly, with ϕ₁ = π/2 gives

$$\begin{eqnarray}&&\begin{array}{l}{\displaystyle \int }_{-\pi /2}^{\pi /2}\left[1+\displaystyle \frac{{\delta }_{H}}{3}(1-3{\left(\sin \phi -\sin {\phi }_{0}\right)}^{2})\right]\cos \phi d\phi \\ \ \ =\ {\displaystyle \int }_{-\pi /2}^{\pi /2}\left[\displaystyle \frac{\theta (\pi /2)}{{\theta }_{0}}-\displaystyle \frac{{\delta }_{H}}{{R}_{t}}\displaystyle \frac{{\left({\sin }^{2}\phi -1\right)}^{2}}{{\cos }^{2}\phi }\right]\cos \phi d\phi ,\end{array}\end{eqnarray} \tag{ 38 }$$

$$\begin{eqnarray}&&\displaystyle \frac{\overline{\theta }(\pi /2)}{{\theta }_{0}}=1+\displaystyle \frac{{\delta }_{H}}{3}(1-3{\left(1-\sin {\phi }_{0}\right)}^{2}).\end{eqnarray} \tag{ 39 }$$

After some simple algebraic manipulations of these equations, one can obtain the expression given in Equation (19).