## Abstract

The east–west striped pattern of clouds in Jupiter’s weather layer is accompanied by a zonal flow containing 12 eastward-going jet streams alternating in latitude with westward-going jet streams. Based on theory, simulation and observations of the Earth’s oceans and atmosphere, it is conjectured that Jupiter’s weather layer is made of bands of constant potential vorticity (PV), where the interfaces between bands are at the latitudes of the maxima of the eastward-going jet streams. It is speculated that the mixing of PV on Jupiter is analogous to the mixing of salt in the ocean by the Phillips effect, which causes the salt density to form a monotonic ‘staircase’. It is hypothesized that the PV in Jupiter’s weather layer is also a staircase, decreasing from north to south. PV is a function of vorticity, as well as parameters with unknown values, e.g. the vertical stratification and the zonal flow beneath the observable weather layer. Therefore, these hypotheses cannot be tested directly. Using an atmospheric model that contains these unknown parameters, we solved the inverse problem and found values of the unknown parameters (and their uncertainties) that best fit Jovian observations. The unknown parameters influence how the zonal flow interacts with large vortices, e.g. the Great Red Spot (GRS; the largest and longest-lived Jovian vortex, centred at 23°S) and the Oval BA (the second largest vortex, centred at 33°S). Although we found that the PV distribution is approximately piecewise-constant and that the peaks of the eastward-going jet streams are at the latitudes of PV interfaces, there is also a PV interface at 20°S, where there is a westward-going jet stream. We find that the zonal PV is not a monotonic staircase due to the ‘backwards’ interface at 20°S. We show that this backwards interface is necessary to make the GRS nearly round, and that without that interface, the Red Spot would be highly elongated in the east–west direction and probably unstable.

## 1. Introduction

Jupiter’s most prominent, and also most permanent, visual feature is its pattern of horizontally (east–west) striped bands of multi-coloured clouds. This pattern was first observed by Hook [1] nearly 350 years ago. Since the first measurements of the zonal (east–west) winds of Jupiter were extracted from Voyager cloud images [2], it has been known that accompanying the east–west cloud pattern is a system of a dozen alternating jet streams, in which the characteristic peak velocities are approximately 100 m s^{−1}. Although there is approximately 10 per cent temporal variability (on scales from hours to decades) in the maximum amplitude of the winds [3], the latitudes of the locations of their maximum and minimum velocities have remained remarkably constant since the first observations in 1979 [4–8]. Figure 1*a* shows the oscillatory nature of the longitudinally averaged east–west flow in Jupiter’s weather layer (i.e. the layer containing the visible clouds and, therefore, the layer of atmosphere in which the observed zonal velocity resides). The figure shows that the zonal flow is remarkably constant over time.

There is no overarching theory for the formation and maintenance of the winds that can account for all of their features, but in general, there are two leading, but competing, scenarios. According to one scenario, the zonal winds are part of a shallow, nearly two-dimensional system driven by ageostrophic forcing at length scales much smaller than the characteristic latitudinal distances between adjacent jet-stream maxima. In this picture, the large-scale winds are created by an inverse cascade of energy from the small to the large scales with the cascade acting in concert with the anisotropic β-plane effect that makes round vortices into zonal flows [9–17]. In the second scenario of zonal wind formation, the winds are deeply rooted and are coupled to the underlying convective zone, which is the ultimate source of energy of the winds [18–23]. In addition to these general scenarios, there are many other theories: for example, in one theory, the zonal winds are due to entropy maximization [24]; in another theory, the zonal winds are due to tidal forcing of the Jovian atmosphere from the Jovian moons [25–27]; and in another theory, there is a resonant interaction between the tidal forces of the Jovian moons and an atmospheric inertial mode, which drives the winds [28].

In all of the theories listed above, the flow has a well-defined Ertel’s potential vorticity (PV). Dritschel & McIntyre [29] have reviewed the importance of PV mixing with respect to the formation of persistent zonal jets in the Earth’s atmosphere and oceans. They point out the fact that numerical simulations and laboratory experiments both show that mixing often causes the PV to organize into parallel bands with nearly uniform values of PV and that these bands result in multiple zonal jets. In many cases, the PV monotonically increases or decreases in the direction perpendicular to the orientation of the bands and jets so that the PV forms a *staircase*. They argue that there is an analogy between the banding and staircasing of PV and the Phillips effect [30], i.e. the vertical staircasing of the density in a stably stratified fluid due to stirring. It has been suggested by Dritschel & McIntyre and many other authors, especially those that are advocates for the scenario that the Jovian zonal winds are a shallow and nearly geostrophic system, that the zonal flows of Jupiter consist of a series of east–west, staircased bands in which the PV within each band has been homogenized by mixing to a nearly uniform value within each band. The interfaces between the homogenized bands are at the latitudes where the eastward-going jet streams have their maximum velocities.

The earliest attempt to fit the observed zonal velocities of Jupiter to this picture was carried out by Ingersoll *et al.* [31], who used a one-dimensional, local (Cartesian approximation of a spherical shell), quasi-geostrophic approximation of the atmosphere with an infinite Rossby deformation radius so that the quasi-geostrophic PV *q* was
(1.1)
where *ω* is the local vertical (*z*) component of the vorticity, *β* is the longitudinal gradient of the Coriolis force and *y* is the local north–south coordinate. In the *i*th band, if the value of the PV is *q*_{i} and constant, and if the velocity is purely zonal, such that , where is the unit vector in the east (positive *x*) direction, then *v*_{x}=*βy*^{2}/2−*q*_{i}*y*+*c*_{i}, where *c*_{i} is an integration constant. That is, the zonal flow forms a series of linked parabolas. The eastward-going jet streams have their maxima at the interfaces between the bands.

Marcus [32] found a better fit to the Jovian zonal flow over a range of latitudes by using a finite value of the Rossby deformation radius *L*_{r}, so that the quasi-geostrophic PV in the Jovian weather layer was
(1.2)
where *ψ* is the stream function. With a finite value of *L*_{r}, if the PV in the *i*th band has constant value *q*_{i}, and if the velocity is purely zonal, then in the *i*th band, , where *A*_{i} and *y*_{i} are integration constants. As in the case of infinite *L*_{r}, the eastward-going jet streams with finite *L*_{r} have their maxima at the interfaces between the bands.

Of course, this picture of exactly homogenized bands of PV with sharp interfaces between them is a great simplification of the actual Jovian flow. Even using the quasi-geostrophic approximation, the above analysis ignores the effects on the PV in the observable weather layer due to an underlying layer of denser fluid. Accounting for that underlying layer by using the one-and-a-half layer quasi-geostrophic model makes the PV in the upper, observable layer equal to [33]
(1.3)
(1.4)
where
(1.5)
and where *g* is the reduced gravity in the weather layer, *h*(*x*,*y*,*t*) is the free-surface height of the upper layer (which acts like a stream function because *ψ*≡*gh*/*f*_{0}), *h*_{b}(*y*) is the height of the rigid bottom layer of fluid, **v** is the weather layer’s horizontal velocity in system III,^{1} is the local vertical unit vector, *f*(*y*) is the Coriolis parameter and *f*_{0} is the local value of *f*(*y*). Note that ∇^{2}*ψ*=*ω* and *ψ*_{b}(*y*)≡*gh*_{b}/*f*_{0} is the stream function of the (unobservable) east–west flow in the underlying layer, which is assumed to be a function of latitude but not longitude or time.

Support for the shallow-layer theory with small-scale forcing as the source of the Jovian wind comes from numerical simulations (cf. Scott & Polvani [14]). Almost all numerical simulations of the quasi-geostrophic equations of motion on a β-plane with a large-scale Rayleigh (or Ekman) dissipation, and a sufficiently large-amplitude, small-scale, homogeneous and isotropic forcing show the formation of zonal flows. Figure 2 shows the PV and zonal velocity of a typical simulation with finite Rossby deformation radius at late times when the flow has become statistically steady. The PV has formed a nearly monotonic staircase as a function of latitude *y* and is nearly independent of longitude *x*. The widths (latitudinal distance between adjacent peaks of the eastward-going jet streams) of the zonal flows in the numerical simulations are approximately equal to 2*πL*_{r}, in fair agreement with Jupiter’s zonal widths (figure 1; assuming that the atmosphere’s vertical stratification and pressure scale-height are uniform over the planet at the weather layer, the deformation radius is proportional to , where *θ* is the latitude). In these initial-value numerical simulations, the mechanism that creates the monotonic steps of PV is the same mixing mechanism that Phillips hypothesized creates monotonic steps of salinity in the upper layer of the ocean [30]. Initially, the PV, *q* as defined in equation (1.2), in the numerical simulations is independent of the east–west coordinate *x*, linear in *y* and equal to *βy* because the initial stream function and vorticity of the flow is zero. As the small-scale forcing builds up the flow, the resulting velocity starts to mix the PV. At random (and unequally spaced) latitudes *y*_{i}, the mixing causes the PV to depart from its initial value and have a near uniform value of *q*=*βy*_{i} in a band centred on latitude *y*_{i}. Thus, the PV forms *bands* of approximately homogenized PV where the boundaries of the bands, by definition, are the latitudes of the near-discontinuities in the PV and are the locations of the peaks of the eastward-going zonal jets. The PV of the late-time, statistically steady flow is approximately *staircased*, decreasing from north to south. The fact that the PV profiles at different values of *x* (longitude) in figure 2 are nearly coincident indicates that the PV is nearly independent of *x* (as is the stream function), and, therefore, the flow is approximately zonal (east–west). During the mixing, the circulation of the initial *q*(*x*,*y*,*t*) is locally conserved, such that , where *L*_{x} is the computational domain size in the *x*-direction and Δ_{i} is the latitudinal width of the *i*th band. Here, *q*(*x*,*y*,*t*=0)=*βy*, so *q*_{i}=(*βy*_{i}). As time goes on, the widths of the mixed bands become broader. Occasionally, neighbouring bands merge together, and homogenize the PV throughout the larger merged band such that the circulation of the PV is locally conserved. Where the homogenized bands of PV abut each other, there are sharp interfaces of PV. The end result is a series of steps of PV *q*_{i} that are monotonic in *y*. Note that the final piecewise-constant PV profile is monotonic in *y* in these simulations *only because the initial distribution of PV was also monotonic*. Our numerical experiments show that in many cases, when the zonally averaged initial value of *q*(*x*,*y*,*t*) is not monotonic in *y*, the final, steady-state values of *q*_{i} are not monotonic in *y*_{i}.

All current formation scenarios have major difficulties. For example, although the longitudinal spacing between adjacent peaks of the eastward-going jet-stream maxima in figure 2 are in good agreement with Jupiter’s zonal flows, the zonal velocities are an order of magnitude smaller than the observed Jovian velocities. The purpose of this paper is not to discuss the failures or successes of any specific model of zonal flow formation on Jupiter, but rather to answer more general questions about the current and observable distribution of PV in Jupiter’s atmosphere. From the observations: (i) do Jupiter’s zonal flows consist of bands of homogenized PV, (ii) if homogenized bands exist, are the interfaces of PV located at the eastward-going jet streams, and (iii) if bands of homogenized PV exist, is the PV in a monotonic staircase, decreasing from north to south?

The remainder of the paper is organized as follows. In §2, we review the solutions of an *inverse problem* where, by fitting the parameters of a reduced-parameter model to observations of the velocities in Jupiter’s weather layer, we can deduce the PV distribution of the zonal flow. In §3, we show how those results answer the three questions posed in the previous paragraph. In §4, we review zone–vortex interactions and show why those interactions reveal the zonal PV profile. We also show qualitatively why the observed shape of the Great Red Spot (GRS) necessitates that the zonal PV does not form a staircase. Our conclusions are in §5.

## 2. Inverse calculations

A commonly used approximation for the layer of atmosphere containing Jupiter’s observable clouds, zonal flow and vortices is the one-and-a-half layer quasi-geostrophic model, in which the PV and velocity are given by equations (1.4) and (1.5) and for which the PV is advectively conserved,
(2.1)
The only independent dynamical variable in equations (1.4)–(2.1) is the free-surface height *h*; *q* and **v** are slaved to *h* through equations (1.4)–(1.5). The influence of an anomalous region of PV, such as a Jovian vortex, falls off exponentially fast far from the vortex with a characteristic decay distance of *L*_{r}. Therefore, far from the vortices, we would expect the far-field velocity to be zonal, independent of *x* and *t* and a function of *y*. It is this far-field zonal velocity , shown in figure 1, in which we are interested in determining whether its PV is in nearly uniform bands with sharp interfaces and whether those bands form monotonic steps. Previous researchers [32,34,35] parametrized the influence of deep layers on the Jovian weather layer with the height *h*_{b}(*y*) of the rigid bottom topography. However, it is more useful (but equivalent) for us to use the PV of the far-field zonal flow . From equation (1.4), there is a one-to-one relation between *h*_{b}(*y*) and ,
(2.2)
where is the free-surface height of the far-field zonal flow, defined by
(2.3)
and where the integral of over the computational domain is arbitrarily defined to be zero.^{2} Because can be measured, *β* and *f*_{0} are known for the domain of interest, the above equation shows that *h*_{b}(*y*) can be parametrized in terms of and *L*_{r}. We cannot use an initial-value solver to solve equations (1.4)–(2.1) for *h* (and thus for **v** and *q*), with parametrized by (via equations (2.2) and (2.3)) because the values of *L*_{r} and (the objective of this study) are unknown. A standard, but tedious, method of dealing with this difficulty would be to solve the inverse problem: (i) guess a set of values for *L*_{r} and , (ii) solve the initial-value equations (1.4)–(2.1) with these parameters until the velocity converges to a statistically steady state, (iii) if the resulting velocity looks like the Jovian zonal velocity, then stop, but if is does not, then (iv) vary the guessed parameter values of *L*_{r} and , and go back to step (ii). When the resulting velocity looks like the actual Jovian zonal velocity, the parameter values of *L*_{r} and used in the initial-value equations are the ‘best-fit’ parameter values, and we could use the best-fit value to determine if is a staircase, as hypothesized. Unfortunately, this method cannot work in the case when the Jovian flow is *exactly* zonal, i.e. with a *q* and **v** that are functions of *y* only, and with **v** in the *x*-direction. In this case, **v**⋅∇*q*≡0. Thus, the steady-state form of equation (2.1) is trivially satisfied for **v** equal to the observed zonal velocity and for *any* values of *L*_{r} and , and therefore, all profiles of are best fits because zonal flows are not affected by (or ).

However, as shown below, (and ) have a strong effect on *non-zonal* velocities and in particular on the Jovian vortices. We can, therefore, exploit an earlier study [36] that was used for finding best fits for Jupiter’s GRS and Oval BA (Jupiter’s second largest vortex) to determine between the latitudes of approximately 9° S and approximately 39° S. We used several sets of observations to find best-fit parameter values, but they all resulted in similar results, so we restrict our discussion here to a review of our results based on observations we took on April 2006 with the Hubble Space Telescope (HST). Our study looked for steady or for uniformly translating (in *x*) solutions of equations (1.4)–(2.1). Most studies in which the velocity fields of the GRS or Oval BA were determined were carried out manually. That is, clouds were assumed to be passive tracers and their displacements over time were used to determine velocities. In laboratory flows, manual methods for extracting velocities were replaced by automated methods such as particle image velocimetry or coherent image velocimetry (CIV) [37,38] to extract velocities automatically using correlation methods. Those automated methods and similar correlation methods [39] fail to work with Jovian images taken more than 42 min apart because the Jovian velocities distort the clouds so much that the automated methods are unable to find correlations between the images. The advantage of longer time intervals between the images is that the extracted velocities have smaller uncertainties, and the advantage of automated methods over manual methods is that the velocity fields obtained of the GRS and Oval BA contain millions of vectors, while the manually extracted field contains hundreds or at most a few thousand. We recently developed a refinement of CIV, called the advected corrected coherent image velocimetry algorithm (ACCIV) that allowed us to use images separated by 10 h to extract automatically the velocities of the GRS and the surrounding zonal flows at 23° S. With ACCIV and CIV, we also extracted the velocities of the Oval BA and its surrounding zonal flows at 34° S [3]. The characteristic velocities of these vortices and zones are approximately 100 m s^{−1}, and the uncertainty of the measurements is approximately ±5 m s^{−1}.

Rather than solving the *inverse problem* with an initial-value code, as outlined above, we use a *steady-state finder*, which is a code that quickly computes the steady or uniformly translating solutions to equations (1.4)–(2.1). This is a much more efficient method and is described in Shetty *et al.* [35,36]. Rather than using generic values in searching for the best-fit parameter values, we use the genetic algorithm as described in Zohdi [40] and Shetty & Marcus [36] to find the values of the unknown parameters that minimize the *L*_{1} norm between the velocities computed with the numerical model and the observed velocities (see §3). Rather than search over all possible parametrizations of the Jovian atmosphere between approximately 9° S and approximately 39° S, we separately solve for the flow near the latitudes of the GRS and of the Oval BA.

To model the flow near the GRS (and the Oval BA), we use the reduced-parameter model that is shown schematically in figure 3. In this reduced-parameter model we *assume* that far from the vortices, the flow forms bands of uniform PV (three bands are shown in figure 3). We exploit the fact that the GRS has been shown by several groups of authors to be a compact region of anomalous PV (cf. the references cited in the review article [32]). Here, we parametrize each vortex by assuming that it consists of a set of nested patches of anomalous PV (two nested regions are shown in figure 3). In our reduced model of the GRS, the PV of the vortex is compact and consists of two nested patches of uniform PV with values and . Different values of PV are indicated by different shades of grey in figure 3. The outer boundaries of these patches have infinitesimal widths, and the diameters of the east–west principal axes of the inner and outer patches are (*D*_{x})_{2} and (*D*_{x})_{1}, respectively. In our models, the GRS and Oval BA are each embedded in a zonal flow that consists of three uniform patches of PV. Far from the vortex, the zonal flow in figure 3 becomes parallel to the east–west axis, and the latitudes of the two interfaces between these three patches asymptote to *y*_{E} and *y*_{W}. At these two interfaces, the ‘jumps’ or changes in the values of the PV are Δ*Q*_{E} and Δ*Q*_{W}, respectively. These interfaces have widths of *δ*_{E} and *δ*_{W}. The GRS or Oval BA, represented by the compact, nested patches of PV, deform the contours of the interfaces between the patches of zonal PV so that they bend around the vortex. In our model, the locus of the four contours between the five patches of uniform PV are all computed numerically and determined such that the flow (which is uniquely determined by the PV) is a steady solution to the quasi-geostrophic equations of motion in a frame translating eastward with velocity *U* (where *U* is uniquely determined by the equations of motion).

In figure 3, the regions of anomalous PV representing the vortices are embedded within the Jovian zonal flows, but their PV is not simply superposed on the zonal bands with *straight* interfaces. A set of east–west zonal winds that consists of a set of bands of uniform PV with straight edges could exist as an equilibrium solution to equations (1.4)–(2.1). However, the *superposition* of those bands of PV with the anomalous PV of a GRS could never be a steady (or uniformly translating) solution of equations (1.4)–(2.1) because the velocity produced by the anomalous PV of a GRS (an extension of the way in which the Biot–Savart law produces circulation around a vortex) distorts the PV of the zonal bands, so that the interfaces between the PV zonal bands would become curved. In a steady (or uniformly translating) flow, the interfaces between regions of uniform PV must be coincident with streamlines. It is, therefore, impossible for an interface of the zonal PV to intersect the outermost streamline of the GRS, so it is necessary that the interfaces between the PV of the zonal flow bend around the GRS, as in the schematic in figure 3. Thus, the *shapes* of each patch, interface and streamline are *not* free parameters, but are computed by the steady- (or uniformly translating) state finding code.

The fast-solver that computes (uniformly translating) steady-state solutions of the quasi-geostrophic equations (1.4)–(2.1) requires slightly different information to initialize than a standard initial-value solver. Instead of requiring an initial velocity field, the steady-state finding code requires that we specify: (i) (which is known by observations), (ii) the PV along the north–south line at the extreme right-hand side of the computational domain, which is effectively , the function for which we seek a best fit, and (iii) the values of the PV along an east–west line that passes through the vortex (which are unknown parameters). The PV along these two lines remains fixed as the steady-state finder is iterated. For the GRS, we choose that east–west line to be the latitude 23° S, and for the Oval BA, we choose it to be latitude 33° S (based on observations, these latitudes are the principal east–west axes of these two vortices).

In this study, the GRS or Oval BA are each represented by at most two nested patches of PV because we found that having more than two did not improve the match between the best-fit solution and the observed ACCIV-extracted velocities.

To model , we assumed that it was nearly piecewise-constant, so we carried out calculations in which consisted of *N*_{z} bands of uniform PV, and in which the interfaces between each band were smoothed over some finite distance (a separate distance for each interface). However, we found that for both the GRS and Oval BA, having just one interface south of the vortex centre and one north of the vortex centre was sufficient and that more interfaces did not improve the fit of the best-fit solution to the observations. Thus, the zonal flows embedding both the GRS and the Oval BA were modelled with three bands of PV, as shown in figure 3. Near the latitude of the interface of zonal PV south (or north) of a vortex centre, there is an eastward-going (or westward-going) jet stream. We model that interface as a jump in PV of strength Δ*Q*_{E} (or Δ*Q*_{W}), located asymptotically far from the vortex at latitude *y*_{E} (or *y*_{W}), and with a latitudinal width *δ*_{E} (or *δ*_{W}), such that for the zonal flow for either the GRS or the Oval BA is discretized as
(2.4)

Thus, our model for the flow near the GRS (or for the Oval BA) consists of 11 parameters: *L*_{r}, (*D*_{x})_{1}, (*D*_{x})_{2}, *q*^{VOR}_{1}, *q*^{VOR}_{2}, *y*_{E}, *y*_{W}, Δ*Q*_{W}, Δ*Q*_{E}, *δ*_{E} and *δ*_{W}. Once the 11 parameter values are specified, the fast-solver finds a uniformly translating equilibrium solution to equations (1.4)–(2.1) with those values. The fast-solver computes the shapes of the interfaces along with the value of *U*, the rate at which the vortex drifts eastward with respect to system III.

## 3. The best-fit solutions

The results of the search for the best-fit parameter values as determined by Shetty & Marcus [36] are shown in table 1. The primary measure for determining that a set of parameters was the best fit was that that the set minimized the *L*_{1} norm of the difference between the velocity field computed with those parameter values and the observed velocity vectors as extracted from the images by the ACCIV method. For both the GRS and the Oval BA, the differences were compared at approximately 1×10^{6} grid points that covered the vortices and the surrounding zonal flows. The uncertainties listed in table 1 for each of the 11 best-fit parameters for the GRS and Oval BA were computed with a boot-strap Monte Carlo method [36,41], which is a standard method for finding uncertainties when the distribution function of the uncertainty at each measurement point is not known.

Quantitatively, the numerical solutions of the GRS and Oval BA evaluated with the best-fit parameter values in table 1 were determined to be good by showing that the *L*_{1} norm difference between the velocity fields of the numerical solutions and the velocities extracted from the observations is small [36]. Qualitatively, those numerical solutions can be seen to be good by examining figures 4 and 5, which show comparisons between the velocities from the best-fit models and the observed velocities along the two principal axes of the GRS and Oval BA. Despite the fact that the model has only 11 free parameters, the best-fit models are good fits to the observations due to the facts that the models are true equilibrium solutions to the equations of motion and that tens of thousands of sets of parameter values were examined to determine the set of values that best fit the observations.

The results that are relevant with respect to the question of the staircasing of the PV for Jupiter are the values of *y*_{E}, *y*_{W}, Δ*Q*_{W}, Δ*Q*_{E}, *δ*_{E} and *δ*_{W} for the GRS and Oval BA. The (up to an arbitrary gauge constant) or far-field zonal PV as a function of latitude between 9° S and 39° S is shown in figure 6. The answer to the question ‘do Jupiter’s zonal flows consist of bands of homogenized PV’ is *partially yes*. If the zonal flow had not formed regions with homogenized PV, and instead the PV was a smoothly varying function of latitude, then the Jovian PV profile near the GRS (or Oval BA) would be poorly approximated by a piecewise-constant function with only three ‘pieces’ with two interfaces between them. The approximation of a continuous PV profile by a piecewise-constant function would improve if the number of discontinuities *N*_{z}, or interfaces, increased (and a near perfect fit as ). However, our numerical experiments found that piecewise-constant PV profiles with *N*_{z}>2 did not improve the fit for the zonal flow near the GRS (and with *N*_{z}>1 did not improve the fit for the zonal flow near the Oval BA), thereby demonstrating that the PV profile in each of our models was well-approximated by regions of uniform PV separated by sharp interfaces. The evidence that the Jovian PV interfaces are sharp is demonstrated by the fact that the widths *δ*_{E} and *δ*_{W} of the interfaces of the best-fit models of both the GRS and Oval BA are much smaller than the latitudinal widths of the bands of homogenized PV (figure 6).

The answer to the question ‘do Jupiter’s zonal flows consist of bands of homogenized PV’ is only a partial, rather than a full, *yes* is due to our definition of a ‘band’. Based on numerical simulations such as those in figure 2, we had defined a band to be the region between two adjacent eastward-going jet streams. Figure 6 and the fact that Δ*Q*_{W}≠0 for the GRS show that with this definition of band, the PV in a band is not uniform. By now re-defining band to mean the region between an eastward-going jet stream and a neighbouring westward-going jet stream, we can state that we have found that the PV within a band is homogenized over the latitudes that we have examined.

The answer to the question ‘if homogenized bands exist, are the interfaces of the PV located at the eastward-going jet streams’ is *no*. Figure 6 shows that *y*_{E} and *y*_{W} for the GRS are at the latitudes of its ambient eastward-going and westward-going jet streams, respectively, and that *y*_{E} for the Oval BA is at the latitude of the eastward-going jet stream to its south.

The answer to the question ‘if bands of homogenized PV exist, is the PV in a monotonic staircase, decreasing from north to south’ is *no* as indicated by the negative and large-amplitude value of Δ*Q*_{W} for the GRS and as shown in figure 6.

It is the fact that Δ*Q*_{W}<0, derived from our solution of the inverse problem, which shows that Jupiter’s zonal PV does not form a monotonic staircase. Because the workings of the inverse problem are not immediately transparent, we show qualitatively in §4 that Δ*Q*_{W}<0 is necessary for the GRS not to be highly elongated in the east–west direction.

## 4. The aspect ratio of the Great Red Spot and non-monotonic zonal potential vorticity

Here, we show how affects the GRS’s *aspect ratio*, which is defined to be north–south minor diameter of *q*^{VOR}_{1}, the GRS’s outer PV anomaly (the white region in figure 3), divided by its east–west major diameter, (*D*_{x})_{1}. We show that the aspect ratio is a strong function of and that if did not have a large *backwards jump of PV* (defined to be a PV interface where the sign of the PV change makes the staircase non-monotonic) at approximately the latitude of the westward-going jet stream just north of the centre of the GRS, the GRS would be highly elongated in its east–west direction. Thus, even without the inverse problem, we can argue that Jupiter’s zonal PV does not form a monotonic staircase.

### (a) Zone–vortex interactions

Using the one-and-a-half layer quasi-geostrophic approximation, Van Buskirk [42] and Shetty *et al.* [35] showed that the aspect ratio of a vortex embedded in a shearing zonal flow is primarily a function of the ratio of the magnitude of the PV of the vortex to the shear of the local zonal flow. To a first approximation, the aspect ratio of a vortex that is a steady or uniformly translating solution of the one-and-a-half layer quasi-geostrophic equations behaves as it does for a Moore–Saffman vortex, i.e. a vortex of uniform vorticity embedded in a zonal flow with uniform shear and with 1/*L*_{r}≡*h*_{b}≡*ψ*_{b}≡0. With no shear, a Moore–Saffman vortex is round, and as the shear increases, the vortex is stretched in the direction of the zonal flow [43]. For a vortex in which the PV consists of two nested regions, such as the vortex shown schematically in figure 3, when [(*D*_{x})_{1}−(*D*_{x})_{2})]>2*L*_{r}, the aspect ratio is primarily a function of *q*^{VOR}_{1} divided by the local zonal shear. The PV in the vortex’s core, due to *q*^{VOR}_{2}, and the other parameters of the model do not strongly affect the aspect ratio because the influence of a delta-function patch of PV falls off exponentially at distances greater than *L*_{r} from the delta function [32]. It would, therefore, appear that the aspect ratio of the GRS is a function primarily of *q*^{VOR}_{1} and the shear of the observed zonal flow , and that the aspect ratio should be independent of . However, due to zone–vortex interactions, the previous sentence is misleading.

To understand the interactions between vortices and zonal flows, it is first necessary to understand that the total PV of the flow is not a *linear* superposition of the far-field zonal PV and the PV anomaly of a vortex. For example, the total PV *q* of our best-fit GRS model (which is a uniformly translating solution of equations (1.4)–(2.1)) is shown in figure 7*a*(i), and it is *not* the sum of the best-fit PV of the far-field zonal flow (figure 7*a*(ii)) as defined in equation (2.2) and the best-fit PV anomaly *q*^{GRS} of the GRS (figure 7*a*(iv)), which we define to be the PV due to *q*^{VOR}_{2} and *q*^{VOR}_{1}. To understand the nonlinearity mathematically, define the *interaction PV* as
(4.1)
The interaction PV is found by numerically computing a steady (or uniformly translating) solution of equations (1.4)–(2.1), thereby finding the parameters , *q*^{GRS}(*x*,*y*) and the total PV *q*(*x*,*y*). Then, *q*^{INT}(*x*,*y*) is computed using the definition in equation (4.1). In general, equilibrium is only possible with a non-zero *q*^{INT}(*x*,*y*) (figure 7*a*(iii)). To understand the PV interaction from a physical point of view, consider the locations in latitude *y* of the PV interfaces with jumps Δ*Q*_{E} and Δ*Q*_{W} as functions of longitude *x*. These locations for the flow around the GRS are the two streamlines shown as thick, white curves in figure 8. The velocity created by the PV anomaly of the GRS causes those streamlines to bend around the GRS, rather than intersect it. (The outer boundary of the GRS’s PV anomaly is also a streamline, and streamlines cannot intersect, except at stagnation points where the velocity must be zero.) Figure 7 shows that *q*^{INT}(*x*,*y*) (figure 7*a*(iii)) is created by the bending of the PV jumps of the zonal flow around the GRS. That is, *q*^{INT}(*x*,*y*) is what needs to be added to the straight streamlines of the far-field zonal flow (figure 7*a*(ii)) so that the zonal flow bends around the GRS.

Figure 7 shows that *q*^{INT} consists of two compact regions of PV anomaly just north and south of the GRS. By definition, *q*^{INT} asymptotes to zero far from the GRS. For the GRS, we found from the solution to the inverse problem that Δ*Q*_{E}>0 and Δ*Q*_{W}<0, so the two patches of *q*^{INT} are both anti-cyclonic (counter-clockwise) like the GRS. The velocities created by the total PV and its three components are shown in figure 7(*b*). The velocity created by embeds the GRS in anti-cyclonic shear; however, the velocity produced by the two anti-cyclones that make up *q*^{INT} embeds the GRS in a cyclonic (clockwise) shear, such that the total shear that the GRS feels from is nearly zero. If the PV of the far-field zonal flow decreased monotonically from the North to South Pole, then Δ*Q*_{W} would be non-negative. A positive Δ*Q*_{W}, together with a positive Δ*Q*_{E}, would make *q*^{INT} into a dipole with an anti-cyclonic patch south of the GRS and a cyclonic patch north of it. This dipole would not significantly reduce the shear in which the GRS is embedded, but rather produce a westward-going jet near the centre of the GRS, which would make the GRS drift rapidly to the west.

With no interface in the zonal PV (i.e. if the zonal PV of Jupiter were everywhere constant), *q*^{INT} would be zero, the GRS would only feel the shear , and (using the values of *q*^{VOR}_{2}, *L*_{r} and Δ*Q*_{E} from table 1—but see the next subsection) the GRS would be much more elongated. Thus, it is the interfaces of the zonal PV and, in particular, the non-monotonicity of , that prevents the GRS from being highly elongated by the local zonal shear .

### (b) ‘Traits’ of the Great Red Spot

In a previous study of reduced-parameter models [35], we showed that the values of *q*^{VOR}_{1}, *L*_{r}, *y*_{E} and Δ*Q*_{E} can be approximately determined without formally solving the inverse problem. We determined the *dependence* of the features or traits of the GRS velocity field on the model’s parameters. For example, we found that the maximum north–south amplitude of the velocity on the east–west principal axis of the GRS (as shown in figure 4*a*) is strongly dependent on the product and is insensitive to all of the other parameters of the model. In hindsight, this lack of dependence on the other parameters is not surprising. First, there is no obvious mechanism for the east–west zonal flow to affect the north–south velocity along the principal axis (because the *q*^{INT} PV is farther than *L*_{r} from the principal axis). Secondly, the longitudes of the two locations of the extrema of the north–south velocity along the GRS’s principal axis are approximately located at distances ±(*D*_{x})_{1}/2 from the centre of the GRS. Thus, because [(*D*_{x})_{1}−(*D*_{x})_{2}]/2>*L*_{r}, *q*^{VOR}_{2} and (*D*_{x})_{2} do not strongly affect the north–south velocity at these locations.

We also found that the north–south velocity along the latitude of the principal east–west axis decayed approximately exponentially in *x* outside the region of anomalous PV (i.e. for |*x*|>(*D*_{x})_{1}, where the centre of the GRS is defined to be at *x*=0). This observation was consistent with the fact that the velocity produced by a delta function of PV in the one-and-a-half layer quasi-geostrophic approximation (equations (1.4)–(2.1)) falls off exponentially, with an e-folding length *L*_{r}, with distance from the delta function [32]. Fitting the north–south velocity along the principal axis for |*x*|>(*D*_{x})_{1}, we found a value of *L*_{r} that was consistent with the value in table 1. Thus, both *q*^{VOR}_{1} and *L*_{r} can be found without relying on solutions of the inverse problem.

From figure 7, it is apparent that the shear on the GRS from could be decreased by either increasing the value of Δ*Q*_{E} or by making Δ*Q*_{W} negative and increasing its absolute value. Thus, it would appear plausible that the effective shear of the zonal flow on the GRS could be reduced so that the GRS did not become elongated and maintained its observed aspect ratio by sufficiently increasing Δ*Q*_{E}. However, we found (and discuss in §4*c*) that by making Δ*Q*_{E} large, the model’s location of the eastward-going jet stream south of the centre of the GRS differs significantly from its observed location.

### (c) Ingredients of the model of the Great Red Spot

One way to illustrate the effect of each of the parameters of the reduced-parameter model of the GRS on the flow, and which illustrates the need for Δ*Q*_{W} to be negative, is to examine the computed best-fit flow when we start with a reduced-parameter model of the GRS with very few parameters, and then see how the fit to observations improves as we add more parameters to the model. In particular, we can examine how the fit changes when we include or exclude the parameter Δ*Q*_{W} from the model. We begin with a model of the GRS that has only one PV interface (i.e. with and Δ*Q*_{W}≡Δ*Q*_{E}≡0). In this model, *L*_{r}, *q*^{VOR}_{1} and (*D*_{x})_{1} are the only parameters of the model. The uniformly translating solution to equations (1.4)–(2.1) using this very constrained model that best fits the observations is shown in figure 9. This model was computed by minimizing the *L*_{1} norm of the difference between the computed and observed velocities. Figure 9*b* shows that this highly constrained model accurately captures the north–south velocities along the latitude of the principal east–west axis exterior to the GRS (i.e. for |*x*|>(*D*_{x})_{1}), which is consistent with our arguments that the east–west zonal flow does not strongly affect the north–south velocity along the latitude that passes through the centre of the GRS. That figure is also consistent with our arguments that the value of *q*^{VOR}_{2} does not strongly affect the magnitude of the peak north–south velocities along the principal east–west axis and that the locations of the north–south peak velocities are at *x*=±(*D*_{x})_{1}/2. However, because , the constrained best-fit model in figure 9 does a poor job at capturing the quiet core of the GRS. Also, figure 9*c* shows that the fit between the computed and observed east–west velocities along the longitude of the north–south, minor principal axis is poor. Most relevant to this discussion is that figure 9*a* shows that the aspect ratio of the modelled GRS is 2.2, rather than 1.6, which shows that the GRS is significantly elongated in its east–west direction.

If we now loosen our constraints and compute the best-fit model with two PV interfaces for the PV anomaly of the GRS (i.e. *q*^{VOR}_{2} and (*D*_{x})_{2} are also parameters of the model), but still have no PV interfaces in the zonal flow (i.e. with Δ*Q*_{W}≡Δ*Q*_{E}≡0), then the best-fit model with these constraints (figure 10) captures the observed north–south flow along the entire latitude line that passes through the centre of the GRS (figure 10*b*). However, the model’s east–west velocities still fail to fit the observed velocities (figure 10*c*), and the GRS is still elongated in its east–west direction with an aspect ratio of 2.4 (figure 10*a*). The east–west elongation of the GRS model in figure 9 is less than the elongation of the GRS model in figure 10 because the PV in the core of the former GRS model is greater than that in the core of the latter.

As figure 11 shows, to capture the observed north–south velocities along the latitude of the major axis, and the observed east–west velocities along the portion of the minor axis longitude to the *south* of the centre of the GRS, it is necessary to have a reduced-parameter model with two PV interfaces for the GRS anomaly (as in figure 10) and also an interface in the zonal flow *south* of the centre of the GRS (i.e. with Δ*Q*_{E}, *y*_{E} and *δ*_{E} as model parameters in addition to all of the parameters used in the model in figure 10). The best-fit model in this constrained model puts the location of the peak of the eastward-going jet stream at its observed latitude. The maximum velocity of the observed eastward-going jet stream is less than that of the model’s maximum velocity (figure 11*c*). However, as discussed by Asay-Davis *et al.* [3], the way in which velocities are extracted from the satellite observations tends to ‘round’ the velocity at the peaks significantly, so it is possible that the model’s maximum velocity of the eastward-going jet stream is more accurate than that of the observations. Figure 11*a* shows that the aspect ratio of the GRS in this model is 2.1, which is less than the aspect ratios of the models in figures 9 and 10 due to the positive value of Δ*Q*_{E}. However, this model’s aspect ratio is still high with respect to its observed value of 1.6.

Using the constrained model whose best-fit solution is shown in figure 11 (i.e. a model without an interface in the zonal PV north of the GRS’s centre), we can examine vortices that are uniformly translating solutions of equations (1.4)–(2.1) but are not the best-fit solution (i.e. that do not minimize the *L*_{1} norm difference between their computed and observed velocities). For this family of vortices, if we increase the value of Δ*Q*_{E} sufficiently (and allow the values of the other parameters of the model to vary), we can find vortices that correctly reproduce the observed aspect ratio of the GRS, in addition to the correct north–south velocities along the latitude line that passes through the centre of the GRS. However, in these models, the computed location of the maximum velocity of the eastward-going jet stream differs significantly from the observed latitude. Only by including an interface of the zonal PV *north* of the centre of the GRS (i.e. using the full 11-parameter model whose best-fit solution is shown in figures 4, 6 and 8) can we correctly reproduce the observed east–west velocities, the north–south velocities *and* the observed aspect ratio of the GRS. For that model, Δ*Q*_{W} is negative, as shown in table 1 and figure 6.

## 5. Conclusions

### (a) Results

Exploiting the manner in which the PV of Jupiter’s east–west zonal flows interact with the planet’s two largest vortices, the GRS and the Oval BA, we have determined Jupiter’s zonal PV between latitudes 9° S and 39° S. We found, as speculated by other authors and as shown by many numerical simulations, that the Jovian PV consists of a series of east–west bands in which the PV is nearly uniform. Also, as speculated, we found that the eastward-going jet streams at 25.6° S and 35.3° S are located at the latitudes of sharp interfaces between the bands of constant PV. The *jump* or change in PV at these interfaces have the same sign and that sign is consistent with the PV forming a monotonic staircase from the North Pole, where the PV is large and positive, to the South Pole, where it is large and negative. However, contrary to most previous hypotheses and simulations, we have shown that there is also a sharp interface of zonal PV at the latitude of the westward-going jet stream at 19.9° S. The jump in PV at 19.9° S is *backwards*, i.e. it has the opposite sign as the other jumps in PV at the zonal interfaces. Thus, between 9° S and 39° S, the Jovian zonal PV is not a monotonic staircase.

### (b) Reasons for a non-monotonic staircase of potential vorticity

One unlikely explanation of why our analysis did not find a monotonic staircase of PV at 19.9° S is that the interface at 19.9° S in figure 6 is spurious. Near the Equator, the quasi-geostrophic approximation needed for equations (1.4)–(2.1) breaks down, and it could be that our analysis near the Equator is therefore incorrect. However, 19.9° S is not so close to the Equator that such a breakdown would be expected. Moreover, there must be some physical mechanism that keeps the GRS from being elongated by the local zonal shear, and other than a backwards jump in zonal PV, no other mechanism has been proposed. Therefore, we believe that the non-monotonicity in the Jovian zonal PV is real, rather than a symptom of bad modelling.

Another possible explanation for the non-monotonic PV interface of at 19.9° S is that Jupiter’s flows are constrained by the Taylor–Proudman theorem so that the zonal winds are ‘columns’ aligned with the rotation axis. Busse [19] argued that the zonal flows are Taylor columns confined between an outer sphere at radius *R* (the weather layer) and an inner sphere of radius *R*_{B} at the boundary of Jupiter’s solid core. Schneider & Liu [23] argued that *R*_{B} is much larger and is at the radius where the atmosphere becomes electrically conducting and the magnetic field retards the flow. With this Taylor-column model of the zonal winds, *q* and are the components of the potential vorticity and far-field potential vorticity in the direction of the rotation axis (rather than the spherical radial components). At most latitudes, equations (1.4)–(2.1) are still valid, and is still well-approximated by equation (2.2). However, now *β*≡0 in equation (2.2) because the axis of rotation is everywhere parallel to the direction of the PV. Furthermore, *h*_{b}(*y*) is primarily due to the spherical geometry that constrains the Taylor columns between *R* and *R*_{B}. A latitude of special interest is *θ*_{T} where the vertically aligned tangent cylinder of the sphere of radius *R*_{B} intersects *R*, i.e. *θ*_{T} is the dividing latitude in the weather layer: for |*θ*|>*θ*_{T}, the bottom boundary of a Taylor column of fluid is at radius *R*_{B}; for |*θ*|<*θ*_{T}, the Taylor column lies outside the radius *R*_{B}, and columns extend uninterrupted from the northern hemisphere at radius *R* to the southern hemisphere at radius *R*. Equation (2.2) is only approximately valid near *θ*_{T} because of the discontinuity in *h*_{b}(*y*) at *θ*_{T} where the value of *h*_{b}(*y*) doubles. Therefore, if the zonal flows are Taylor columns, our result in figure 6 is only approximately valid near *θ*_{T}. None the less, it is worthwhile considering whether *θ*_{T} is equal to 19.9°, and whether the discontinuity in *h*_{b}(*y*) at *θ*_{T} causes a barrier to PV mixing there, resulting in a PV interface. A value of *θ*_{T}=19.9° implies that , which is close to the value of 0.96 used by Schneider & Liu [23]. If there were no meridional PV mixing on Jupiter, and if Jupiter rotated as a solid body, then the abrupt change in *h*_{b} at *θ*_{T} would create a jump in the value of at *θ*_{T}, and it is reasonable to ask whether some remnant of that jump currently exists at 19.9° S. The answer is probably *no*. The change in PV at *θ*_{T} of a solid-body rotating Jupiter is not backwards (i.e. the change in PV from the Pole to *θ*_{T} because of the topographic PV has the same sign as the jump in PV at *θ*_{T}). Therefore, the jump in PV at *θ*_{T} for a solid-body rotating Jupiter has the *opposite sign* of the jump in PV that we found at 19.9° S and that was shown in figure 6.

We believe that the most plausible explanation for the non-monotonicity of the Jovian zonal PV is its ‘initial condition’. Both the theoretical argument (Phillips effect) for, and the numerical simulations of (cf. Scott & Polvani [14]), staircasing require that the initial distribution of PV be approximately monotonic. If the ‘initial’ Jupiter had solid-body rotation and if its relative vorticity in a spherical weather layer were dominated by its radial component, then the initial PV would be monotonic in latitude decreasing from north to south.^{3} However, for the Jovian atmosphere, what do we mean by its initial condition? The complex processes by which the final atmosphere formed and spun up to its present rotation rate are not known. If, for example, the Jovian atmosphere formed with a strong Equatorial jet stream, then the initial condition of the zonal PV might not have been monotonic in latitude. Our numerical simulations of the evolution of zonal PV (similar to the evolution of the flow shown in figure 2), show that the late-time zonal PV is a function of the initial conditions. For example, the number of bands of nearly uniform PV in the late-time flow changes as a function of the initial conditions [12], and the late-time zonal PV can be non-monotonic for some non-monotonic initial conditions. Therefore, future theoretical and computational studies of Jupiter’s formation and ‘spin-up’, or of the formation of its zonal flows, especially those that use global circulation models, should be required to reproduce and explain the current non-monotonicity of the zonal PV near its Equator.

## Footnotes

One contribution of 9 to a Theme Issue ‘Dynamical barriers and interfaces in turbulent flows’.

↵1 System III is the rotating coordinate system based on a 9 h 55 m 29.711 s day, in which many researchers report the atmospheric velocities of Jupiter; it is based on the rotation rate of Jupiter’s magnetosphere.

↵2 A ‘gauge’ constant can be added to ,

*h*_{b}and/or*h*, representing a shift in origin or the vertical coordinate. Such a shift has no effect on the velocity or any other observable quantity because*h*acts as a stream function and only its spatial derivative enters into the dynamics. A change in gauge will change the PV by a constant value. We have arbitrarily chosen the gauge for*h*_{b}such that the PV along the southern boundary of our computational domain is zero.↵3 In the model in which the zonal flows are Taylor columns, the Ertel PV of a solid-body rotating Jupiter is not monotonic: the magnitude of the PV gradually decreases from the North Pole to (

*θ*_{T})^{+}; decreases by a factor of 2 from (*θ*_{T})^{+}to (*θ*_{T})^{−}; increases back to the value it had at (*θ*_{T})^{+}at a latitude of ; and becomes infinite at the Equator. However, 10° is so close to the Equator, the stratification due to gravity is nearly perpendicular to the planet’s rotation axis, and due to the stratification, it is unlikely that the Taylor-column hypothesis is valid there.

- This journal is © 2011 The Royal Society