The turbulent dynamics of Jupiter’s and Saturn’s weather layers: order out of chaos?

Abstract

The weather layers of the gas giant planets, Jupiter and Saturn, comprise the shallow atmospheric layers that are influenced energetically by a combination of incoming solar radiation and localised latent heating of condensates, as well as by upwelling heat from their planetary interiors. They are also the most accessible regions of those planets to direct observations. Recent analyses in Oxford of cloud-tracked winds on Jupiter have demonstrated that kinetic energy is injected into the weather layer at scales comparable to the Rossby radius of deformation and cascades both upscale, mostly into the extra-tropical zonal jets, and downscale to the smallest resolvable scales in Cassini images. The large-scale flow on both Jupiter and Saturn appears to equilibrate towards a state which is close to marginal instability according to Arnol’d’s 2nd stability theorem. This scenario is largely reproduced in a hierarchy of numerical models of giant planet weather layers, including relatively realistic models which seek to predict thermal and dynamical structures using a full set of parameterisations of radiative transfer, interior heat sources and even moist convection. Such models include (amongst others) the Jason GCM, developed in Oxford, which also represents the formation of (energetically passive) clouds of NH₃, NH₄SH and H₂O condensates and the transport of condensable tracers. Recent results show some promise in comparison with observations from the Cassini and Juno missions, but some observed features (such as Jupiter’s Great Red Spot and other compact ovals) are not yet captured spontaneously by most weather layer models. We review recent work in this vein and discuss a number of open questions for future study.

Introduction

The visible atmospheres of Jupiter and Saturn exhibit a bewildering array of features over a huge range of horizontal scales. Both planets are covered in ubiquitous layers of clouds, with the upper layers composed mainly of ammonia ice but with other substances (including water) present at deeper levels (e.g. Irwin 2009). On the largest scales, such clouds are organised into broad bands with meridional widths on the order of 10⁴ km, almost (but not quite) symmetrically disposed about the equator and with boundaries closely aligned with lines of latitude. Such cloud bands have been seen on Jupiter since around 1630 (e.g. see Graney 2010), but also occur in more muted form on Saturn beneath a tenuous layer of upper tropospheric haze. The cloud bands are also now well established to occur in association with an intense and persistent pattern of eastward and westward zonal jet streams, with wind speeds ranging from a few tens of m s⁻¹ up to nearly 400 m s⁻¹ in Saturn’s equatorial jet (e.g. see Ingersoll et al. 2004; Vasavada and Showman 2005; Del Genio et al. 2009).

Within this pattern of zonal jets are found a variety of waves and vortices on scales ranging from a few hundred km up to that of the cloud bands at around 10⁴ km. Such waves and vortices are remarkable both for their diversity and also for their persistence and longevity. The most famous of these features is Jupiter’s Great Red Spot (GRS), a feature which may be traceable back to a large spot discovered in the seventeenth century. But more recent observations have identified many other features that resemble the GRS on both Jupiter and Saturn, though on a somewhat smaller scale. The GRS is a giant anticyclone, but persistent cyclonic vortices have also been observed though tend to have a less compact and circular form (at least for those found at low-mid latitudes) than the corresponding anticyclones. Coherent wave trains are also seen from time to time, though few are as persistent and apparently stable and regular as Saturn’s North Polar Hexagon. Discovered in images from the Voyager spacecraft in the late 1980s (Godfrey 1988), this feature is characterised by a regular wavenumber m = 6 structure around a latitude of 76° N and has been observed most recently by the Cassini spacecraft to persist coherently up to the end of the mission in 2017 and has since continued to be observable from Earth-based telescopes (Hueso et al. 2020).

Most recently, the Juno mission has significantly extended the observational record to include close orbital passes over both poles of Jupiter, where a whole new set of phenomena have been discovered. At high latitudes, Jupiter’s cloud bands and zonal jets are much less evident and the dominant features close to the poles are regular arrays of almost circular, cyclonic vortices surrounding a single cyclone sitting close to the pole itself (Adriani et al. 2018). Infrared images (Adriani et al. 2018) reveal a wealth of small-scale structure associated with such vortices, but most remarkable is the extent to which the eightfold or fivefold cyclone arrays have remained stable and persistent over timescales of at least 1–2 years. Other novel findings from Juno have included the discovery that lightning, most likely associated with intense moist convective storms involving latent heat release from condensation of water vapour, is most common at high latitudes in both the northern and southern hemispheres of the planet. Prior to Juno, the only available information on Jovian lightning had noted a prevalence of such storms in the cyclonic cloud “belts” compared to the anticyclonic “zones” (Little et al. 1999), but could not detect larger scale trends with latitude.

Of particular interest from both the Juno and Cassini missions has been the new insights and constraints provided on the vertical structure of the weather layers of Jupiter and Saturn. Until recently, the “weather layer” has been commonly assumed to include the region of the atmosphere influenced by a combination of solar and thermal radiation, upwelling heat from the deep interior and moist convection, typically spanning a height range down to a pressure level of around 10 bars on both planets. However, both the Juno and Cassini missions have included very close passes in high inclination orbits (in Cassini’s case shortly before the end of the mission when it was deliberately targeted to enter Saturn’s atmosphere) which allowed relatively high order harmonics of each planet’s gravity field to be determined. From a comparison with predictions from simple models of the mass distribution associated with pressure variations in balance with hypothetical patterns of eastward and westward zonal winds, estimates have been obtained of the likely depth scale to which the observed cloud level zonal mean zonal winds could extend into the deep interiors of each planet (Kaspi et al. 2018; Galanti et al. 2019). The values obtained for the depth scale, ~ 3000 km for Jupiter and ~ 9000 km for Saturn, are much deeper than is commonly assumed for the weather layers of either planet (where radiative forcing and moist convection may be prevalent), indicating that the weather layers in these planets may be underlain by a pattern of deeper zonal jets whose vertical:horizontal aspect ratio is close to O(1) (cf Dowling 2019). The presence of these jets, however, is consistent with expectations for the depth of their electrically neutral molecular envelopes (Liu and et al. 2008; Cao and Stephenson 2017), although uncertainties may still be quite large (e.g. see Kong et al. 2018). Juno has also deployed a relatively long wavelength microwave mapping sensor which has, among other results, shown the perturbations in temperature and NH₃ concentration within the GRS to penetrate to depths greater than 300 km (Li et al. 2017). This instrument has also discovered a previously unsuspected, compact and persistent “plume” of NH₃ vapour upwelling at the equator, with relatively low concentrations at higher latitudes (Li et al. 2017).

There is now, therefore, a very considerable range of observations from both Jupiter and Saturn available to guide and constrain models of the structure and dynamics of their weather layers, and to stimulate and challenge their further development. In this article, we review the current state of knowledge and model development for these systems, focusing particularly on questions such as:

• i. What process(es) primarily drive and energise the cloud-level meteorology on Jupiter and Saturn?

• ii. Why does the cloud-level circulation self-organise into patterns of zonal cloud bands and eastward and westward jet streams?

• iii. What factors determine the scales and intensities of these zonal jets?

• iv. What processes sustain the most energetic waves and vortices in the weather layers of these planets?

• v. What are the implications for the transport of heat, momentum and chemical tracers?

  And in particular

• vi. What aspects of the observable flows can be accounted for solely by the dynamics of a shallow weather layer?

In this article we consider each of these questions, with a particular focus on the energetics of the circulation and the extent to which many of the observations can be captured in a numerical model that is limited to the upper troposphere and stratosphere of a Jupiter-like gas giant planet. The “Energetics of eddies and zonal jets” section considers first the strong zonal mean circulation and current ideas for why it takes the form observed, including some of the implications for the apparent stability or instability of the zonal jets. The “Global circulation models of gas giant weather layers” section reviews the recent development of relatively complex, 3D time-dependent numerical models of jovian weather layers, with examples of some recent successes but also some of their limitations. We conclude in the “Discussion and conclusion” section with a brief outlook for future research in this area.

Energetics of eddies and zonal jets

The maintenance of the zonal jets on Jupiter by the action of atmospheric waves and eddies has been established since the time of the Voyager encounters with these planets in the early 1980s (Beebe et al. 1980; Ingersoll et al. 1981). Velocity measurements from manually tracking cloud motions in Voyager ISS images allowed the estimation of zonally averaged horizontal eddy momentum fluxes $\left( {\overline{u'v'} } \right)$ and the shear structure of the zonal mean zonal flow $\left( {d\bar{u}/dy} \right)$, from which the rate of transfer of kinetic energy into the mean zonal flow could be computed. The measured rate of transfer (typically 1–3.5 × 10⁻⁴ W kg⁻¹) was initially greeted with some scepticism (e.g. Sromovsky et al. 1982), on the grounds that the manual method of velocity tracking used by Ingersoll et al. (1981) could inadvertently introduce biases in the calculation. This was important, not least because were such a large apparent transfer rate to extend significantly from Jupiter’s cloud tops into the deeper atmosphere (e.g. down to pressure levels of ~ 3 bars) it would imply that more than 10% of the incoming thermal energy from solar and internal heating was being converted into zonal kinetic energy.

More recent work has employed more advanced, automated image correlation techniques to measure velocity fields for both Jupiter (Salyk et al. 2006; Asay-Davis et al. 2009; Choi and Showman 2011; Galperin et al. 2014) and Saturn (Del Genio and Barbara 2012) using ISS images from the Cassini spacecraft. These have enabled estimates of the eddy-zonal flow kinetic energy transfers at the cloud tops which should be relatively free of some of the observational biases of the manual tracking work. Nevertheless, the results have broadly confirmed the earlier estimates for Jupiter, with an overall conversion rate of around 10⁻⁴ W kg⁻¹ at latitudes up to 50° N and S. For Saturn, estimates seem closer to around half of this value (Del Genio and Barbara 2012). The question remains, however, as to how deep this conversion might extend into the abyssal layers on both planets, and hence what is the overall thermodynamic efficiency of the conversion of thermal to zonally averaged mechanical energy.

Resolution of this question (and the corresponding question as to the rate of production of KE from thermally generated potential energy), however, would require global measurements of velocity and temperature structure that extend much deeper into both atmospheres than has been feasible so far. However, it seems clear that at least the extra-tropical zonal jets on Jupiter and Saturn owe their origin largely to the nonlinear rectification of momentum transfer by eddies. But this raises other related questions such as:

• Do eddies of all scales and locations act in this way to drive zonal flows?

• What is the primary source of the eddies found at the cloud tops?

• What factors set the dominant scales of the zonal jets?

• More generally, how do turbulent cascades of energy and enstrophy operate in gas giant atmospheres?

KE spectra

An important first step in addressing the first of these questions is to determine the kinetic energy spectrum of motions in the atmospheres of Jupiter and Saturn. A number of attempts to estimate spectra for Jupiter’s atmosphere have been made in the past (e.g. Mitchell 1982; Choi and Showman 2011), though most of them have focused mainly on estimating these only in the zonal direction. These have provided hints that Jupiter’s atmosphere might exhibit evidence for KE spectra resembling the classical Kolmogorov k^−5/3 spectrum over a certain range of scales, prompting comparisons with expectations from classical theories for 3D or 2D turbulence. However, it is clear from even a superficial inspection of images of Jupiter that the motions in the atmosphere are highly anisotropic with quite different structures in the zonal and meridional directions. A proper analysis should therefore include a consideration of both the zonal and meridional directions, for which the natural decomposition uses the set of spherical harmonics (Choi and Showman 2011; Galperin et al. 2014).

Figure 1 shows an example of the KE spectrum of Jupiter’s atmosphere, derived from automated tracking measurements of cloud motions in Cassini images during the fly-by in December 2000 by Galperin et al. (2014) between latitudes ± 50^o, in which the KE was projected onto spherical harmonics (see Appendix A of Galperin et al. 2014 for detailed analysis of resolution and measurement uncertainties). This shows KE decomposed into the zonal mean component (for zonal wavenumber m =0, shown with a thin line) and the other (m ≠0) components (thick line). Clear differences are seen between the zonal and residual spectra, with the zonal spectrum roughly following an n⁻⁵ trend for n >20 and somewhat flatter (nearer n^−5/3) at larger scales. The residual spectrum, however, is quite different, being almost flat for n <50 and then steepening towards a form close to the Kolmogorov slope n^−5/3.

Fig. 1

Kinetic energy spectrum of Jupiter’s cloud layer (in units of Ω²R²), obtained by Galperin et al. (2014) from Cassini images. Spectra are projected onto spherical harmonics (with total wavenumber index n) and decomposed into the zonal mean flow E_Z(n) (1, thin line) and residual KE E_R(n) (2, thick line), such that total energy E_T(n) = E_Z(n) + E_R(n). The analysis used the CB2 filter images with a native resolution of 0.05° pixel⁻¹ and a correlation box size of 0.95°, corresponding to a spectral resolution limit of around 1° or n ≈ 200

This spectrum shows few features in common with that of either an idealised 3D or 2D turbulent fluid. For the latter case, one would expect to see an isotropic spectrum with two segments (e.g. Davidson 2015), one with the Kolmogorov form n^−5/3 for n < n_I (where n_I represents the scale at which KE is injected, forcing the flow), where the spectrum is dominated by an upscale (spectrally local) cascade of kinetic energy, and the other with a slope of n⁻³ for n > n_I, dominated by the downscale transfer of enstrophy (relative or potential). In Jupiter’s case, the spectrum is clearly anisotropic, at least at large scales, but for which a Kolmogorov-like segment is only found at small scales in the isotropic component of the residual (eddy) spectrum E_R (and arguably at large scales in the zonal spectrum E_Z). At large scales, the KE is heavily dominated by the zonal mean jet circulation, with more than 90% of the total kinetic energy residing in the zonal jets, at least at the cloud tops (Galperin et al. 2014).

The steep n⁻⁵ character of the zonal spectrum led Sukoriansky et al. (2002) to suggest the relevance of the zonostrophic paradigm for geostrophic turbulence, forced at a relatively small scale and dissipated across a range of scales. This represents a generalisation of the classical Kolmogorov–Batchelor–Kraichnan theory to take account of the effects of spherical planetary curvature and the consequent poleward gradients of planetary vorticity and dominance of dispersive Rossby waves at large scales (cf. Rhines 1975). The consequences for the expected KE spectrum (at least for a predominantly barotropic flow) lead to an anisotropic spectrum with zonal and residual components of the spectrum of the form:

$$ E_{\text{R}} = C_{\text{K}} \varepsilon^{2/3} \left( {\frac{n}{R}} \right)^{ - 5/3} , $$

(1)

$$ E_{\text{Z}} = C_{\text{Z}} \left( {\frac{\varOmega }{R}} \right)^{2} \left( {\frac{n}{R}} \right)^{ - 5} \quad {\text{for}}\left( {n < n_{\beta } } \right), $$

(2)

where C_K and C_Z are dimensionless constants (O(1) in each case), R is the planetary radius and $\varepsilon$ is the upscale kinetic energy transfer rate. The total wavenumber n_β is the so-called anisotropy wavenumber, defined in a spherical domain by

$$ \frac{{n_{\beta } }}{R} = \left[ {\frac{{\varOmega^{3} }}{{R^{3} \varepsilon }}} \right]^{1/5} , $$

(3)

such that the spectrum is anisotropic for n < n_β and the two spectra cross over at n = n_β. This assumes that the turbulent flow is forced at an even smaller scale n_f for which n_f> n_β.

The zonal spectrum E_Z for both Jupiter and Saturn appear to fit the zonostrophic form in Eq. (2) above in both slope and amplitude, given a value for $C_{\text{Z}} \approx 0.5$ (e.g. Sukoriansky et al. 2002; Galperin et al. 2014). For E_R on Jupiter, taking C_K to be approximated by the Kolmogorov–Kraichnan constant $C_{\text{K}} \approx 6$, the observed spectrum for n > n_β ~ 57 suggests an estimate for $\left| \varepsilon \right|$ of around 10⁻⁵ m² s⁻³ (Galperin et al. 2014). However, this method of estimating $\varepsilon$ does not distinguish between upscale or downscale transfers, leaving open the interpretation of the spectrum as zonostrophic in the sense defined by Sukoriansky et al. (2002) and subsequent work.

Spectral transfers and structure functions

To determine the direction as well as the magnitude of the transfer of energy or enstrophy (squared relative vorticity ½ ζ²) between different scales within the spectrum, it is necessary to consider higher order covariances of velocities than the quadratic product for the energy spectrum itself. This is commonly analysed within the spectral framework described e.g. by Frisch (1995) or Burgess et al. (2013), for which

$$ \frac{{\partial {\mathcal{E}}_{n} }}{\partial t} + \varPi_{n} = - 2\upsilon \varOmega_{n} + {\mathcal{F}}_{n} , $$

(4)

where $\varepsilon_{n}$ represents the cumulative energy contained in scales larger than wavenumber n, ν is the kinematic viscosity and Ω_n is the cumulative enstrophy over wavenumbers ≤ n, so the first term on the right hand side of (4) represents viscous dissipation, and ${\mathcal{F}}_{n}$ represents the forcing of kinetic energy at wavenumbers ≤ n. The term $\varPi_{n}$ represents the flux of energy between wavenumbers through wavenumber n, due to turbulent interactions between different scales. The sign of $\varPi_{n}$ is such that $\varPi_{n} > 0$ implies a direct energy flux towards smaller scales and $\varPi_{n} < 0$ implies an inverse transfer towards larger scales. $\varPi_{n}$ is computed (for the rotational, non-divergent component of the flow) from the nonlinear tendency for enstrophy (i.e. the squared vorticity, $1/2 \zeta_{n}^{m*} \zeta_{n}^{m}$, where $\zeta_{n}^{m}$ represents the (m, n) component of the vertical component of relative vorticity), as represented for flows decomposed into spherical harmonics, $Y_{n}^{m} = P_{n} \left( {\cos \theta } \right)e^{im\varphi }$ (where m is the zonal wavenumber index and n the total wavenumber index), defined by

$$ J_{n} = - \frac{1}{4}\mathop \sum \limits_{{m = - n}}^{n} \zeta _{n}^{{m*}} \left( {\user2{u}.{\mathbf{\nabla }}\zeta } \right)_{n}^{m} + \zeta _{n}^{m} \left( {\user2{u}.{\mathbf{\nabla }}\zeta } \right)_{n}^{{m*}} , $$

(5)

where ()^* denotes the complex conjugate. This represents the rate of change of enstrophy with total wavenumber n due to quadratic nonlinear interactions with other wavenumber components with zonal index m. The corresponding spectral tendency for kinetic energy, $I_{n} ,$ can be obtained from the relationship

$$ I_{n} = \frac{{a^{2} }}{{n\left( {n + 1} \right)}}J_{n} , $$

(6)

and the (isotropic) spectral fluxes of energy ($\varPi_{n}$) and enstrophy ($H_{n}$) are then obtained from

$$ H_{n + 1} = - \mathop \sum \limits_{\ell = 0}^{n} J_{\ell } ; $$

(7a)

$$ \varPi_{n + 1} = - \mathop \sum \limits_{\ell = 0}^{n} I_{\ell } . $$

(7b)

Similar tendencies and fluxes for the divergent, irrotational components of the flow can also be obtained via the tendency of the squared horizontal divergence, also decomposed into spherical harmonics.

Profiles of the spectral fluxes of KE and enstrophy in Jupiter’s cloud tops were obtained recently by Young and Read (2017) from the same two-dimensional, horizontal wind vector fields derived from tracking cloud features by Galperin et al. (2014) as mentioned above. Some of the results are illustrated in Fig. 2. Figure 2a shows the profile of the spectral enstrophy flux, indicating a positive (forward) cascade at almost all scales down to the resolution of the measurements (around 0.5° in latitude and longitude) with no evidence of an inertial range (i.e. no plateau in H_n).

Fig. 2

Profiles of spectral fluxes of enstrophy (squared vorticity) and kinetic energy with total horizontal wavenumber and wavelength at Jupiter’s cloud tops, as derived from Cassini ISS images of Jupiter during the 2000 fly-by by Young and Read (2017). a Enstrophy flux, b total (rotational) kinetic energy flux and c kinetic energy flux due only to eddy–eddy interactions. The direction of energy or enstrophy transfer is indicated by the arrows. G14g and G14s refer, respectively, to the global mosaic and single image cloud tracking analyses of Galperin et al. (2014) while C11 refers to the global mosaic analysis of the same images by Choi and Showman (2011)

The total kinetic energy transfer is shown in Fig. 2b and indicates a much more complicated picture. In particular, it shows a negative (upscale) flux between wavelengths of 2 × 10³–2 × 10⁴ km, with a hint of an inertial range plateau around 10⁴ km. The flux of KE changes sign, however, at two scales, (i) at a wavelength ~ 2–3 × 10⁴ km where energy converges, and (ii) around a wavelength of 2 × 10³ km, from which energy diverges, forming a direct forward cascade at smaller scales. Figure 2c shows the same computation including only eddy–eddy interactions (i.e. with the zonal mean components removed). This shows a qualitatively similar pattern, but with an amplitude 5–10 times smaller. This essentially shows that the strongest interactions are direct transfers between eddies and zonal flows, though the spectrally local, eddy–eddy cascade is of a similar form. The overall picture is consistent with the generation of KE as eddies of a scale close to 2000 km, and an accumulation of KE at scales around 20,000 km, close to the scales of the large-scale zonal bands.

The use of a spherical harmonic decomposition for this analysis assumes that this overall picture for transfers of KE and enstrophy between scales applies across the whole planet, and neglects any local variations in the structure of the energy budget. Young and Read (2017) did, however, also verify the picture derived from spectral fluxes using 3rd order structure functions, which derive the sense and magnitude of energy transfers between scales more directly. The longitudinal structure function, $F_{L3} \left( r \right)$, is defined by

$$ F_{L3} \left( r \right) = \delta u_{L}^{3} , $$

(8)

where $\delta u_{L} \left( r \right) = u_{L} \left( {\varvec{x} + \left| \varvec{r} \right|} \right) - u_{L} \left( \varvec{x} \right)$, angle brackets denote an average over all x and $u_{L}$ is the velocity component lying along the path joining x and x + r. Positive values of $F_{L3} \left( r \right)$ indicate a positive correlation between separating particle paths and kinetic energy, and hence a tendency to transfer KE to larger scales (e.g. Frisch 1995; Davidson 2015). For a classical Kolmogorov–Batchelor–Kraichnan (KBK) turbulent spectrum with a k^−5/3 slope, the third-order structure function is expected to scale as $F_{L3} \left( r \right)\sim \varPi r$ (e.g. Davidson 2015). Young and Read’s (2017) results show a change in sign of $F_{L3} \left( r \right)$ at values of r around 2000–3000 km for both local and global velocity maps towards a nearly linear segment with r, confirming the suggested reversal of the KE cascade towards a direct, forward cascade at scales smaller than 2000–3000 km.

The nature of the forward cascade at scales < 2000 km is not well characterised. Such a direct cascade (with a spectrum roughly consistent with the KBK slope of − 5/3; see “KE spectra” section above) is not consistent with the classical paradigm of 2D or QG turbulence (e.g. Charney 1971; Salmon 1980; Davidson 2015), although a somewhat similar dual cascade pattern is also observed in the Earth’s oceans (e.g. Scott and Wang 2005). Internal gravity waves may play an important role if the atmosphere is stably stratified in this region, but observations are insufficient at present to characterise this segment of the cascade with any confidence.

The significance of this injection scale is important with respect to the dominant mechanism that generates the eddies that energise the turbulent cascades. At around 2000–3000 km, this is highly suggestive of the first baroclinic Rossby deformation radius $L_{d}$. This has been estimated in several previous studies of Jupiter (e.g. Achterberg and Ingersoll 1989; Ingersoll and Kanamori 1995; Ingersoll et al. 2004; Vasavada and Showman 2005; Read et al. 2006) to lie around $L_{d} \sim 2000$ km (to within a factor ~ 2), consistent with the scale of the KE source found by Young and Read (2017).

Based on Galileo and Cassini observations of Jupiter, Ingersoll et al. (2000) and Gierasch et al. (2000) have suggested the potential significance of moist convection as an energetically important process for generating eddies in Jupiter’s atmosphere, estimating that up to 50% of Jupiter’s internal heat might be transported upwards through moist convective thunderstorms in the troposphere. However, the typical scale of Jovian thunderstorms is quite a lot smaller than $L_{d}$, even when aggregated into storm clusters, suggesting that this may not be the principal source of KE in the large-scale turbulent KE cascades measured at Jupiter’s cloud tops. More likely would be the possibility that a form of baroclinic instability is responsible for the generation of eddies on scales comparable with $L_{d}$ (e.g. Vallis 2017), in a similar fashion to the Earth’s oceans (e.g. Scott and Wang 2005), although the style of instability may well be rather different from that found more familiarly in terrestrial atmospheres. On Earth, for example, baroclinic instability tends to be strongly enhanced through horizontal thermal gradients at the lower boundary. Such a boundary does not exist on the gas giant planets, however, at least within the upper weather layers. But from the viewpoint of the potential vorticity structure, the tropopause may act as a weak interface, allowing the development of baroclinically unstable modes that are enhanced in amplitude near the top of the tropopause through thermal gradients along the tropopause itself (e.g. Conrath et al. 1981; Kaspi and Flierl 2007; Read et al. 2020). Many uncertainties remain, however, because the vertical structures of the atmospheres of Jupiter and Saturn are relatively unexplored in observations.

Potential vorticity structure and jet stability

The vorticity structure of the zonal jets on Jupiter and Saturn is of particular interest in relation to their likely stability to breaking up into meanders and/or vortex-like eddies. As we discuss below, it may also hold clues as to why the zonal spectrum, E_z(n), adopts the steep form ~ n⁻⁵. A long-standing conundrum has been the observation that the latitudinal profile of $\bar{u}\left( y \right)$ (where y is the northward coordinate) appears to violate the well known Rayleigh–Kuo necessary (though not sufficient) condition for instability, namely, that $\beta - \partial^{2} \bar{u}/\partial y^{2} = \partial \zeta_{a} /\partial y$ changes sign across the domain, where β is the northward planetary vorticity gradient due to spherical planetary curvature and $\zeta_{a}$ is the absolute vorticity. This is clearly seen in Fig. 3b, which shows the latitudinal profiles of northward vorticity gradient (solid lines) and $\beta = 2\varOmega \cos \theta /R$ (dashed line). The fact that the solid line crosses the dashed one demonstrates that the Rayleigh–Kuo criterion is violated, principally in association with westward zonal jets.

Fig. 3

Observed structure in a zonal velocity, b northward gradient of relative vorticity and c quasi-geostrophic potential vorticity of Jupiter’s cloud-level zonal jets, as derived by Read et al. (2006) from Cassini cloud tracking measurements (Porco et al. 2003) and infrared thermal retrievals from the Cassini Composite Infrared Spectrometer (CIRS)

This tendency for many zonal jets on Jupiter and Saturn was noted in analyses of the early Voyager observations (e.g. Ingersoll et al. 1981), though was explained at the time as indicating that aspects of the vertical structure of the flow (which could not be measured directly) might ameliorate the change of sign of $\partial \zeta_{a} /\partial y$ (see also Scott and Dunkerton 2017). A more complete stability measure is related to the potential vorticity structure, based on quasi-geostrophic (QG) theory, for which the instability criterion equivalent to the Rayleigh–Kuo condition is the Charney–Stern–Pedlosky (CSP) criterion (e.g. Vallis 2017), namely, that $\partial Q/\partial y$ should change sign within the domain, where

$$ Q = \varvec{k}.\nabla \times \varvec{u}_{g} + f_{0} \frac{\partial }{\partial z}\left[ {\frac{1}{{N^{2} \left( z \right)}}\frac{{\partial \psi_{g} }}{\partial z}} \right] + \beta y $$

(9)

is the QG potential vorticity, $\varvec{u}_{g}$ is the geostrophic velocity, $\psi_{g}$ the geostrophic stream function and z is the log-pressure pseudo-height coordinate if Q is computed on pressure surfaces. The corresponding profile of Q for Jupiter was derived from a combination of cloud tracked winds and retrieved temperatures from Cassini observations by Read et al. (2006) and is shown in Fig. 3(c). Despite the addition now of realistic (albeit low vertical resolution) thermal information to include the second “stretching” term on the RHS of Eq. (9), it is evident that $\partial Q/\partial y$ also changes sign with latitude, suggesting that even the CSP instability criterion is also satisfied. An important caveat for both these criteria is that they are necessary conditions, but not sufficient to guarantee instability (e.g. see Vallis 2017). However, the fact remains that the zonal velocity on both gas giant planets robustly satisfies both criteria in repeated measurements.

An alternative interpretation was offered by Dowling (1993), Stamp and Dowling (1993) and Dowling (1993), who noted that a stronger (nonlinear) stability criterion was provided in the so-called second stability theorem of Arnol’d (Arnol’d 1966)—hereafter “Arnol’d II”. Although still only a necessary condition, it does allow for the possibility of flows exhibiting a change of sign of $\partial Q/\partial y$ while remaining marginally stable to non-zonal perturbations. For this situation, the condition for stability reduces to

$$ \left( {\frac{{\bar{u} - \alpha }}{\partial Q/\partial y}} \right) \ge L_{d}^{2} , $$

(10)

where α is a constant representing a shift in the frame of reference and $\partial Q/\partial y$ can take either sign. Marginal stability corresponds to the equality condition in Eq. (10).

Based on some earlier diagnostic calculations of the vorticity structure of Jupiter’s winds (assuming an analogy with rotating shallow water theory), Dowling (1993), Stamp and Dowling (1993) and Dowling (1993) suggested that both Jupiter’s and Saturn’s weather layers were close to marginally stable conditions, despite violating the Rayleigh–Kuo condition for stability. The essence of this approach relates to the ability of Rossby waves propagating within a local region with a given value of $\bar{u}\;{\text{and}}\; \partial Q/\partial y$ to phase lock with a parallel wave train with a different value of $\bar{u}\;{\text{and}}\;\partial Q/\partial y$. Marginal instability implies just one value of α that would allow the gravest (deepest) Rossby wave modes to phase-lock, hence implying a unique reference frame in which instability (through interaction between pairs of Rossby wave trains) is possible.

A unified interpretation of these two stability theorems has recently been developed by Dowling (2014, 2020) through the introduction of a dimensionless parameter, Ma, analogous to the Mach number in supersonic flows. Ma is defined with respect to the phase speed of the gravest Rossby waves that can propagate within a shear flow and takes the general form

$$ {\text{Ma}} = - \kappa_{1}^{2} \frac{{d\psi_{G} }}{dQ}, $$

(11)

where $\kappa_{1}^{2}$ is the first baroclinic eigenvalue connecting the (eddy component of the) QG potential vorticity with the stream function and tends to $L_{d}^{ - 2}$ for small Burger number (= [L_d/R]², where R is the planetary radius). For finite L_d in the long wave limit, the phase speed of Rossby waves $c \to - \left( {\overline{\partial Q/\partial y} } \right)L_{d}^{2}$, so that “Ma” takes the form of a (signed) Mach-type number, $- \left( {u - \alpha } \right)/c$, since $\left( {u - \alpha } \right)/\left[ {\partial Q/\partial y} \right] \to - d\psi_{G} /dQ$. Thus, Eq. (10) becomes equivalent to requiring Ma > 1 for stability while the CSP criterion implies stability for Ma  < 0. Instability is therefore expected under either criteria for Ma⁻¹ > + 1 with marginal stability/instability corresponding to Ma = + 1 (see Dowling (2019, 2020) for more details).

This unique characteristic of a marginally unstable flow was later exploited by Read et al. (2009a, b) to determine the implied reference frame for marginal instability on both Jupiter and Saturn, based on local computations of the zonal mean Q and $\bar{u}$. The results showed that the marginally unstable reference frame on Jupiter was indistinguishable from the System III frame, determined from the rotation of its magnetic field. For Saturn, however, the results indicated a significantly different frame of reference that rotated several minutes faster than the System III frame determined from Voyager observations (period = 10 h 39 min 24 ± 7 s; Desch and Kaiser 1981). With a value of 10 h 34 min 13 ± 20 s, it implies a shift in apparent zonal velocity on Saturn’s equator of around 85 m s⁻¹, significantly altering the shape of its zonal velocity profile towards a more symmetrical distribution of $\bar{u}$ between eastward and westward jets, more like Jupiter.

This result was controversial at the time, because Cassini measurements of Saturn’s magnetic field obtained a rather different rotation period to Voyager of 10 ^h 47 ^m 6 ^s. Moreover, the magnetic rotation period measured by Cassini continued to vary throughout the Cassini mission (e.g. Ye et al. 2016), strongly indicating that neither these measurements nor Voyager’s were accurately sensing Saturn’s “true” internal rotation period. Coincidentally, however, another method of estimating Saturn’s internal rotation rate was also published in 2007 (Anderson and Schubert 2007), using measurements of its oblateness and gravity field, which also inferred an even faster rotation period of 10 h 32 min 35 ± 13 s. These results have received some further support from more recent measurements of Saturn’s gravity field and oblateness (Helled et al. 2015) and measurements of perturbations in Saturn’s rings (Mankovitch et al. 2019), which, respectively, suggest internal rotation periods of 10 h 32 min 45 ± 46 s and 10 h 33 min 38 s, though with the latter having fairly large uncertainties ~ ±2 min. Nevertheless, it seems clear that, assuming Saturn has a well defined internal rotation period, it must be significantly faster than the System III frame determined by either Voyager or Cassini radio observations. Moreover, the influence of this internal rotation seems to be felt in the dynamical structure of the flow, even in the weather layers of both Jupiter and Saturn.

One further implication of the gas giant weather layers being close to a marginally unstable state relates to the n⁻⁵ form of the zonal KE spectrum. Planetary rotation allows for the existence of inflection points in the latitudinal profile of $\bar{u}$ while remaining marginally stable under the Rayleigh-Kuo criterion. This may be consistent with a zonal KE spectrum with an index of at least − 4 (e.g. Huang et al. 2001), although this makes a number of simplifying assumptions (notably of linearity and a uniform value of β). Huang et al. (2001) argue that the n⁻⁵ form of the zonal KE spectrum could represent a saturated flow, in which nonlinearities limit the amplitude of both $\bar{u}$ and its vorticity. The observations, therefore, that both Jupiter and Saturn exhibit an n⁻⁵ spectrum in their zonal KE of the form in Eq. (2) appears to be consistent with the notion that their zonal mean flow structure is defined by a marginally unstable, saturated state.

Global circulation models of gas giant weather layers

One of the most stringent tests of our understanding of a complex system such as Jupiter’s or Saturn’s weather layer is the extent to which we can reproduce observed structures and behaviour, and even predict features not yet observed, in a model based on fundamental principles and assumptions, i.e. without overt empirical “tuning”. In this regard, the task of developing such models for Jupiter and Saturn is still at a fairly primitive stage, although a rich hierarchy of dynamical models of varying sophistication and complexity have appeared in the past ~ 40 years.

Idealised models

The earliest numerical circulation models made minimal assumptions concerning the vertical structure of Jupiter’s atmosphere, treating it essentially as a single, shallow fluid layer on the surface of a sphere (Williams 1978, 1979, 1985; Williams and Yamagata 1984; Williams and Wilson 1988), based on the rotating shallow water equations. Later work led to some account being taken of the likely presence of a deep abyssal layer (assumed to be relatively dynamically passive, though potentially with deep zonal jets) underlying a shallow (but dynamically active) weather layer, forming effectively a 1½ layer (or “1¾ layer”; e.g. Flierl et al. 2019, Dowling 2020) model (e.g. Ingersoll and Cuong 1981; Dowling and Ingersoll 1988, 1989). This class of model was able to capture certain features and phenomena that were consistent with some observations, such as the formation of patterns of parallel jet streams and the short-term behaviour of some waves and oval vortices, though they were strongly idealised and lacked a representation of a number of realistic physical processes to sustain or dissipate circulation features or their vertical structure. Nevertheless, this class of model has continued to be employed in a number of more recent theoretical studies intended to investigate particular hypotheses for the origin and stability of various features (e.g. Li et al. 2006; Showman 2007; O'Neill et al. 2015; Thomson and McIntyre 2016, Thomson 2020).

Intermediate complexity models

Because of their computational expense, as well as a lack of detailed observations of the vertical structure of gas giant atmospheres, researchers have been relatively slow to apply full scale 3D numerical models, based on the meteorological Primitive Equations, to this region of gas giant atmospheres. The earliest attempts to apply such 3D models were still quite heavily idealised (Williams 1997, 2002, 2003) though did seek to impose plausible patterns of diabatic forcing and vertical structure to a shallow weather layer. Results did confirm the role of baroclinically unstable eddies in generating patterns of multiple, upper surface-enhanced zonal jets which tended to drift in latitude when heating was confined to upper levels.

Around the same time, Dowling et al. (1998) published the first results from their Explicit Planetary Isentropic Coordinate (EPIC) model, formulated (as the acronym would suggest) using isentropic coordinates in the vertical direction. This was initially applied to a number of idealised problems, mostly in regional rather than global domains, to study the dynamics of artificially initialised vortices resembling Neptune’s Great Dark Spot (LeBeau and Dowling 1998) though has more recently been applied to other features such as Jupiter’s equatorial hot spots (Showman and Dowling 2000) and White Ovals (Morales-Juberías et al. 2003), Saturn’s ribbon waves and polar hexagon (Morales-Juberías et al. 2011) and other giant planet vortices (e.g. Brueshaber et al. 2019).

Other weather layer models have appeared which use a well established dynamical core to solve the dynamical equations but also have used idealised thermal relaxation or mechanical forcing to drive circulations, either in a limited area domain (e.g. Yamazaki et al. 2004; Zuchowski et al. 2009a, based on a dynamical core from the UK Met Office) or globally (e.g. Lian and Showman 2008; O'Gorman and Schneider 2008; Medvedev et al. 2013, Thomson and Vallis 2019). These models were able to capture certain features relevant to gas giant meteorology, such as the interaction of waves and zonal flows, though were more diagnostic than predictive. Lian and Showman (2008), for example, showed that forcing localised within a shallow weather layer could nonetheless lead to a much deeper response e.g. in the pattern of zonal winds. Lian and Showman (2010) also introduced a very simple parameterisation to represent the effects of moist convective storms as a random pattern of mass sources and sinks, suggesting the feasibility of driving the zonally banded circulation through nonlinear interactions between moist convection and zonal jets in a rapidly rotating, stratified spherical atmosphere.

Thomson and Vallis (2019) describe versions of the Isca model, based on the GFDL FMS spectral dynamical core, set up in either highly simplified 1½ level or intermediate-complexity configurations intended to emulate Jupiter-like conditions. Like the EPIC model, Isca can be configured in an “n + ¾ layer” form with prescribed deep zonal jets in an abyssal layer at the base of the model domain, that have allowed comparison of simplified 1¾ level model results (e.g. Thomson and McIntyre 2016; Thomson 2020) with simulations with more realistic vertical structure (Thomson and Vallis 2019). But few attempts have been made in any of these studies so far accurately to match the conditions and forcing processes in the troposphere or stratosphere of Jupiter or Saturn.

More realistic weather layer models

Most recently, several groups have attempted to develop global circulation models that include a more realistic set of parameterisations that represent forcing and dissipation processes that (arguably!) more closely resemble those found in the upper tropospheres and stratospheres of Jupiter and/or Saturn, and even in some cases including the ice giant planets (Uranus and Neptune).

The EPIC model, referenced above in the previous section, is a notable example which has also been extended recently to use hybrid isentropic-pressure coordinates (Dowling et al. 2006) and more realistic parameterisations e.g. for cloud microphysics (Palotai 2008). However, most recent studies using EPIC to date have not implemented realistic radiative forcing in the upper troposphere or stratosphere.

Lian and Showman (2010) used the MITgcm dynamical core, forced by a combination of thermal relaxation towards a temperature structure inspired by the observed structures of the gas and ice giant atmospheres, trending towards a neutrally stable deep interior (down to a depth ~ 100 bar), with latent heat release associated with the condensation of a water vapour tracer. The results demonstrated the spontaneous formation of multiple, extra-tropical, eddy-driven zonal jets, together with reasonably strong equatorial jets that reproduced the observed directions for conditions corresponding to either Jupiter, Saturn or Neptune. The equilibrated flow also exhibited a pattern of zonal jets that violated the Charney–Stern stability criterion but in most cases, remained close to a marginally unstable Ma⁻¹ ≈ 1 state with respect to Arnol’d’s second stability theorem, much as discussed in Sect. "Potential vorticity structure and jet stability" above.

In a similar vein, Schneider and Liu (2009) used the pseudo-spectral dynamical core from the Flexible Modeling System (FMS) of Princeton’s GFDL to obtain simulations of all of the gas and ice giant weather layers (down to a pressure of ~ 3 bar). In contrast to Lian and Showman (2010), however, they used a semi-gray radiative transfer parameterisation with insolation of the observed magnitude (nominally for all four gas and ice giants but for a planet with a fixed Earth-like obliquity) and realistic upward heat fluxes from the interior but no latent heating or moist convection. A pattern of surface drag was also used to emulate magnetohydrodynamic drag at depth to close the energy budget. This model was also able to capture patterns of extra-tropical zonal jets on reasonably realistic spatial scales and amplitudes for each planet, together with equatorial jets of the appropriate sign.

They also computed the exchange of KE between different scales, using the same spectral flux diagnostic as discussed above in the “Spectral transfers and structure functions” section. This is illustrated in Fig. 4a, which is taken from Schneider and Liu (2009) and clearly shows a negative (upscale) KE flux between wavenumber indices 10 < n < 100 turning into a positive (downscale) flux for n ≥120, much as found in the observations of Young and Read (2017). This calculation also shows most of the KE spectral flux residing in the eddy-zonal flow interactions, with only around 20% of the flux in the eddy–eddy components and no obvious evidence of a plateau in spectral flux indicative of an inertial range. Subsequent work (Liu and Schneider 2010, 2015) investigated mechanisms for zonal jet formation in their simulations, highlighting the roles of baroclinic eddies in driving the extra-tropical zonal jets and the relative magnitudes of the interior and solar thermal fluxes in each case in determining the predominant direction (though not necessarily the magnitude) of the equatorial jet, and investigating the scaling of jet strength with model parameters.

Fig. 4

Spherical harmonic spectral fluxes of KE in model simulations of a Jupiter-like weather layer, obtained by a Schneider and Liu (2009) and b in the Jason model (Young et al. 2019a, 2020). Fluxes in a are averaged over 3 bar in pressure whereas in b this is evaluated (like the observations) at a pressure of 1 bar

Most recently, Spiga et al. (2020) and Cabanes et al. (2020) have developed a new model, based on a new dynamical core on an icosahedral grid (DYNAMICO), that includes a more realistic multi-band (correlated-k) radiative transfer parameterisation which computes heating rates specific to the composition of Saturn’s atmosphere (Guerlet et al. 2014). The model also takes account of a realistic heat flux from the deep interior, together with parameterisations for turbulent dissipation and dry convective adjustment, though without moist convection or latent heating. Like the model of Schneider and Liu (2009), the domain only extends to 3 bar pressure with a similar pattern of “MHD drag”. Nevertheless, the model was able to capture many features of Saturn’s observed thermal structure and meteorology, including its extra-tropical multiple zonal jets and polar vortices. The equatorial jet obtained, however, was not particularly realistic in its magnitude though was at least in the observed direction.

The Jason weather layer model

Young et al. (2019a, b) have recently developed a new GCM for Jupiter’s weather layer, designated Jason, which attempts to include a reasonably full range of moderately realistic parameterisations as well as high spatial resolution in its dynamical core. Like the model of Lian and Showman (2008), the dynamical core is based on MITgcm (in its latitude–longitude grid configuration), but follows Schneider and Liu (2009) in utilising a semi-gray radiative transfer scheme that is calibrated to reproduce approximately the observed vertical distribution of solar and infrared radiative fluxes. A uniform internal heat flux is input at the bottom of the domain (at ~ 18 bar) and the model includes a uniform, weak linear drag at the bottom boundary to emulate “MHD drag” (including in the tropics). The model can also include the effects of latent heat release by condensation of a water vapour tracer and parameterisations of both dry and moist convection (using the heat engine parameterisation of Zuchowski et al. 2009b). As in the earlier work of Zuchowski et al. (2009c), tracers can be advected for water, NH₃ and NH₄SH and their condensable cloud aerosols, allowing the model to develop its own layered cloud structures which can then be transported by the dynamics and form features that can resemble observed cloud structures. An example of the modelled distribution of water clouds is illustrated in Fig. 5, which shows a snapshot of water ice column density across the planet. Clouds clearly organise into zonal bands with sharp edges, rendering visible small-scale travelling waves propagating along zonal jets.

Fig. 5

Visualisation of a snapshot of the field of water cloud ice column density from a simulation of the weather layer circulation of Jupiter at a horizontal resolution of 0.7° × 0.7° in latitude × longitude using the Jason model (Young et al. 2019b). This simulation does not include moist convection

Like the other realistic models mentioned above, Jason exhibits the spontaneous formation of a pattern of zonal jets, whose scale and magnitude depends somewhat on which forcing processes and parameterisations are active, as well as on model resolution (e.g. see Fig. 6). Most runs so far have been carried out, as in the other models in the “More realistic weather layer models” section, using a resolution of approximately 0.7° × 0.7° in latitude and longitude. This resolution, or higher, is necessary when simulating Jupiter’s meteorology in order to resolve adequately the first baroclinic Rossby radius on Jupiter, L_d, which on Jupiter is around 1000–3000 km (Achterberg and Ingersoll 1989; Vasavada and Showman 2005; Read et al. 2006) and represents the dominant energy-containing scales for baroclinic instabilities on a planet like Jupiter. To date, however, the Jason model has not imposed a pattern of deep zonal jets at the bottom of the model domain, although the pattern of eddy-driven zonal jets does seem to extend to the (fixed pressure) base of the model where it spontaneously develops a corresponding geostrophically balanced topography in geopotential height. This captures some effects of deep zonal jets on the vorticity dynamics of the resolved flow in the model weather layer, though more passively than in some other models (such as Dowling and Ingersoll 1989; Thomson and Vallis 2019).

Fig. 6

Zonal mean zonal velocity cross-section, averaged over the final 100 days of a run using the Jason Jupiter model with semi-gray radiative transfer and moist convection active (assuming $4 \times$ solar composition of H₂O). White contours indicate zonal velocities in excess of 30 m s⁻¹ and the tropopause (where $\partial T/\partial p$ changes sign) is shown as a grey contour

As evident in spectral KE flux calculations from the Jason simulations (e.g. see Fig. 4b and Young et al. 2020), zonal jets arise through predominantly upscale nonlinear transfers of KE from eddies into zonal flows. In the Jason simulations, this takes place between horizontal wavelengths of 20,000–30,000 km (comparable to the latitudinal spacing of the zonal jets and the so-called Rhines scale π(2U/β)^1/2; Hide 1969; Rhines 1975) and ~ 2000–3000 km (comparable to L_d). At scales smaller than 2000 km, the KE transfer is direct, towards small scales, much as found in the observations and in the simulations of Schneider and Liu (2009), consistent with a source of KE at scales around 3000 km.

The origin of this KE is indicated from preliminary calculations of the spectral conversion rate of potential energy (PE) to KE, $C_{n}^{E}$, which represents the spectrally decomposed covariance of vertical velocity and density (Young et al. 2020). Preliminary results suggest a significant conversion of APE to KE at scales around 3000 km. If confirmed, this would be consistent with the interpretation that baroclinic eddies on scales comparable with L_d are dominant, not only in generating large-scale structures and zonal jets through upscale nonlinear interactions, but also in energising the KE cascade by conversion from potential energy at small scales in a scenario reminiscent of the QG baroclinic turbulence paradigm of Salmon (1980). Unlike this classical paradigm, however, the inverse KE cascade to large scales is also accompanied by a forward cascade of KE to even smaller scales. Such a dual KE cascade is also found, for example, in the Earth’s ocean circulation (e.g. Scott and Wang 2005). The nature of this forward cascade on scales smaller than the deformation radius is still not well understood, either for Jupiter or the Earth’s oceans, but seems likely to involve three-dimensional processes (such as internal gravity waves) related to the stratification. Such processes should be explored in more detail in future modelling work and observations.

Although baroclinic instability may dominate the dynamics, at least at extra-tropical latitudes, moist convection is also active in runs that include the moist convection parameterisation. The distribution of moist convective events, as captured by the parameterisation, is far from uniform, however. Figure 7 shows an example from a Jason simulation in which the Zuchowski et al. (2009b) parameterisation is active (Young and Read 2020), which shows (a) the net heating/cooling rates due to moist convection as a function of latitude and height, and (b) the vertical distribution of moist convective adjustment during the run. This clearly demonstrates that the most frequent and intense moist convection takes place at high latitudes, with weak heating and cooling close to the water condensation level, around 3–4 bars in this case, and stronger warming at the detrainment levels around 200–600 hPa. This distribution corresponds remarkably well to the observed distribution of lightning, as found recently by Brown et al. (2018) from Juno observations. This tendency to favour moist convection at high latitudes may reflect the stabilising effect of absorbed solar radiation at low latitudes, which tends to suppress deep convection down to levels of a few bars (though not beneath this). One aspect of the observations that is not well captured in the Jason simulations, however, is the significant asymmetry between the northern and southern hemispheres on Jupiter, for which the north seems to have been much more active recently than the south. Recent Juno results have also noted a strong hemispheric asymmetry between the northern and southern hemispheres (Moore et al. 2018), which might suggest that the circulation at deep levels exerts an influence on the weather layer that should be allowed for, perhaps through a “¾-layer” extension, in future work.

Fig. 7

Distribution of moist convection in a Jason simulation of Jupiter’s weather layer at a resolution of 0.7° × 0.7° in latitude × longitude (Young and Read 2020), showing (a) the net rate of change of humidity as a function of latitude and pressure, and (b) the fraction of convecting grid points as a function of pressure, showing the confinement of moist convection between ~ 4 bars and 200 hPa

The modelled transport and distribution of ammonia vapour is another quantity that is represented in the Jason model (Young et al. 2019b), and can be compared with recent Juno MicroWave Radiometer (MWR) observations (Li et al. 2017) and recent ground-based radio observations (de Pater et al. 2019). Figure 8 shows the distribution in latitude and height (pressure) of the NH₃ vapour mass mixing ratio 100–200 days after initialisation in a latitudinally uniform state. This simulation clearly captures a strong, ascending plume of ammonia from the dense reservoir in thermochemical equilibrium at depth, capped at around 500 hPa where NH₃ vapour tends to condense into clouds. A similar pattern of strongly ascending NH₃ vapour was found by Li et al. (2017) in the Juno MWR observations, where a narrow (~ 5° wide) equatorial ammonia plume was found to rise up from a reservoir at p >20 bar, capped at around 700 hPa. In the Jason simulations, however, the narrow NH₃ plume does not persist indefinitely, but tends to spread in latitude as the circulation evolves, equilibrating with a much more diffuse equatorial maximum in NH₃ MMR than observed. This almost certainly indicates that, although the model exhibits some realistic tendencies in the organisation of the circulation, there are processes missing (and as yet poorly understood) that either keep the equatorial NH₃ confined to a narrow ascending plume or deplete the upper tropospheric NH₃ concentration at mid-latitudes (e.g. Guillot 2020). It will be of interest, however, to see if future observations (e.g. from further ground-based radio observations; cf de Pater et al. 2019) continue to confirm the persistence of Jupiter’s compact equatorial NH₃ plume.

Fig. 8

A map in latitude and height (pressure) of the mass mixing ratio of NH₃ vapour obtained 100–200 days after tracers were injected into a Jason model simulation of Jupiter’s weather layer, without H₂O latent heating and moist convection (Young et al. 2019b)

Discussion and conclusions

The results highlighted in this review paint an increasingly clear picture of how the circulation of Jupiter’s and Saturn’s weather layers work, although many questions still remain unresolved.

Global circulation in weather layers

Thermodynamic forcing in these weather layers seems to be dominated mainly by a combination of large-scale differential solar heating down to pressures of a few bars and upwelling heat from the planetary interior, balanced by radiative cooling to space. The dynamical response, in the form of baroclinic instabilities and active but complex and anisotropic turbulent cascades of energy and enstrophy, leads to the generation of a pattern of eddy-driven zonal jets with forcing that is probably strongest in the upper troposphere. The response on large horizontal scales, however, seems to penetrate to much deeper levels, as evident in the gravity harmonic results from Juno and Cassini (Kaspi et al. 2018; Galanti et al. 2019).

Somewhat remarkably, our most recent models suggest that moist convection, although energetically potentially significant (Gierasch et al. 2000), does not seem to play a direct role in driving large-scale motions, at least at low-mid latitudes. Young and Read (2020) find that parameterised moist convection acts to intensify the circulations produced mainly by solar and interior heating towards more realistic intensities but does not seem to result in any new phenomena, at least on global scales (although it may be more important on local or regional scales). This may, of course, be due to the limitations of the parameterisation schemes and limited spatial resolution in the simulations. But the disparity of scales between the typical sizes of individual moist convective storms and the inferred scale of energy injection in the observed turbulent KE cascade (at least for Jupiter; Young and Read 2017, but also likely for Saturn; see Cabanes et al. 2020) would seem to favour baroclinic eddies over convective storms as the principal energising features for large scale flows. However, it will be important in future to establish from observations whether a similar situation prevails on Saturn and the other giant planets.

The main response to such forcing in the weather layer is in the form of the dominant pattern of zonal jets, at least in the extra-tropics, the strength of which may be limited by a tendency for the jet flows to adjust to a state of marginal barotropic instability. Arnol’d’s 2nd stability criterion would seem to provide a significant constraint (Stamp and Dowling 1993) that has proved effective in pointing towards a particular reference frame in which marginal instability is satisfied for both Jupiter and Saturn (Read et al. 2009a, b). It will be of great interest in future work to explore whether this marginally unstable configuration also occurs on the ice giant planets, Uranus and Neptune, although observations of a sufficient quality to test this hypothesis may need to await a new mission to these remote planets.

Atmospheric models

The most recent generation of numerical weather layer models, using physically based parameterisations constrained by observations, would seem to offer some significant promise in capturing many, but by no means all, observed features of these atmospheres. The demands on model resources, however, are challenging, requiring quite high horizontal (< 1° × 1° in latitude × longitude) and vertical resolution, deep domains (at least to 10–20 bar) and long integration times (> 200,000 simulated rotation periods) in order to reach thermal and dynamical equilibrium with realistic heating and cooling rates. Recent results, however, show such approaches to be rewarded by the reproduction of reasonably realistic patterns of extra-tropical zonal jets, cloud bands and certain aspects of meridional tracer transport.

Virtually all weather layer models published to date, however, seem to fail to reproduce a number of major features in observations, at least spontaneously. The large oval vortices on Jupiter and Saturn, such as the Great Red Spot and White Ovals on Jupiter, do not develop spontaneously in most of these model simulations without a finite-amplitude initial perturbation, although a few models have proven able to develop an equilibrated GRS-like vortex from small amplitude perturbations of a different shape and scale; e.g. Dowling and Ingersoll (1989); see Fig. 17 and associated text). This raises some intriguing questions as to the physical origin of these features and might suggest either some missing physical processes in many models or inadequate spatial resolution. Deep-rooted hot spots might offer one possible explanation, given the appearance of isolated vortical structures in some deep convection models (e.g. see Heimpel et al. 2015). The formation of stable arrays of polar vortices or symmetric wave patterns, such as Saturn’s North Polar Hexagon, have also yet to be demonstrated in a realistically maintained circulation (though see Morales-Juberías et al. (2015) for a promising approach). The equatorial jets in all four giant planets pose other significant challenges to such models, for which the question remains unresolved as to whether these jets require a deep circulation (e.g. Yano et al. 2005; Kaspi et al. 2009; Heimpel et al. 2005) or can be sustained solely within a shallow weather layer (e.g. Yamazaki et al. 2005; Liu and Schneider 2011).

More generally the connection between shallow weather layers and deep-seated circulations and convection remains poorly understood. Deep models are currently unable to capture the complex structures of the weather layer while shallow weather layer models do not represent a physically realistic lower boundary. The development of coupled deep convection and weather layer models, similar to the coupled atmosphere–ocean models for the Earth’s climate, may ultimately be the only way to address this adequately, although this approach also presents major technical challenges.