2. Characteristic scales in jets

In addition to the inverse cascade of energy, Rhines theorised the existence of a steep $k^{-5}$ power law. This was found in Chekhlov et al. (Reference Chekhlov, Orszag, Sukoriansky, Galperin and Staroselsky1996) to manifest in the vicinity of the meridional wavenumber axes ($k_y$) (Huang, Galperin & Sukoriansky Reference Huang, Galperin and Sukoriansky2001; Sukoriansky et al. Reference Sukoriansky, Galperin and Dikovskaya2002). As the energy along $k_y$ follows different scaling laws to the rest of the spectrum, Sukoriansky et al. (Reference Sukoriansky, Dikovskaya and Galperin2007) divided the total energy spectrum into zonal $E_Z(k)$ and residual $E_R(k)$ (the inverse cascade) components:

(2.1)

where $E(k)$ is the energy density as function of wavenumber $k$, $\epsilon$ is the energy injection rate, and $C_{\beta }$ and $C_K$ are the dimensionless Rhines and Kolmogorov constants, empirically determined to be $C_{\beta } \sim 0.5$ and $C_K \sim 5\unicode{x2013}6$ (Sukoriansky et al. Reference Sukoriansky, Galperin and Dikovskaya2002). Here, $\phi = \arctan (-k_x/k_y)$ is the angle measured anticlockwise from the $k_y$-axis, and $k_x$ is the zonal wavenumber. We chose a threshold $\phi = {\rm \pi}/8$ as the size of the wedge around the $k_y$-axis that contains zonal energy. By equating (2.1a) and (2.1b), Sukoriansky et al. (Reference Sukoriansky, Dikovskaya and Galperin2007) found the transition point in the total spectrum, i.e. the wavenumber at which the energy contained within zonal modes exceeds that of the residual modes:

(2.2)\begin{equation} k_\beta = \left(\frac{C_\beta}{C_K}\right)^{{3}/{10}} \left(\frac{\beta^3}{\epsilon}\right)^{1/5} \approx 0.5\left(\frac{\beta^3}{\epsilon}\right)^{1/5}. \end{equation}

This cross-over wavenumber (Pelinovsky Reference Pelinovsky1978; Maltrud & Vallis Reference Maltrud and Vallis1991; Vallis & Maltrud Reference Vallis and Maltrud1993) marks the point at which the spectrum anisotropises in $\beta$-plane turbulence. Note that this is not to say that for length scales where $k>k_\beta$ the system is isotropic. Instead, we expect energy near $k_x=0$ (the meridional axis) to dominate for length scales exceeding $k_\beta$.

After $k_{\beta }$ is crossed, most of the energy in the system resides in the zonal modes (2.1). By integrating this over wavenumber space, Sukoriansky et al. (Reference Sukoriansky, Dikovskaya and Galperin2007) obtained the total energy of the system, and found that most of the energy resides in the largest occupied scale on the spectrum, namely the Rhines scale given by (1.1). The size of the jets scales according to the largest, most energetic scale, confirming the Rhines (Reference Rhines1975) theory. When there is a sufficiently strong large-scale dissipation such that the energy in the system reaches equilibrium before the flow approaches the domain scale, the Rhines scale settles at the equilibrium wavenumber given by

(2.3)\begin{equation} k_{jet} = k_R = \left(\frac{\beta^2}{8E_{eq}}\right)^{{1}/{4}}, \end{equation}

where $E_{eq}$ is the energy at equilibrium (Danilov & Gurarie Reference Danilov and Gurarie2002; Sukoriansky et al. Reference Sukoriansky, Dikovskaya and Galperin2007). For a system to exhibit jet formation, Sukoriansky et al. (Reference Sukoriansky, Dikovskaya and Galperin2007) argued using a range of numerical experiments that characteristic scales must be well separated as follows:

(2.4)\begin{equation} k_0 \ll k_R \ll k_\beta \ll k_f \ll k_D, \end{equation}

where $k_0$ is the largest allowable scale, $k_f$ is the forcing scale, and energy dissipation occurs at small scales with wavenumber greater than $k_D$. Thus the Rhines scale must be much larger than the cross-over wavenumber to allow significant anisotropy to develop, and must also be much smaller than the domain scale so that domain effects do not interfere. Furthermore, the anisotropic wavenumber must not be close to the forcing scale, so that the inverse cascade may develop. The last inequality arises from classic turbulence arguments (e.g. Kraichnan Reference Kraichnan1967) that the small-scale dissipation must not be significant in the inertial range. Flows obeying (2.4), particularly where $k_\beta /k_R \gtrsim 2$, are described by Galperin et al. (Reference Galperin, Sukoriansky and Dikovskaya2008) to be in the zonostrophic regime. Independently, Scott & Dritschel (Reference Scott and Dritschel2012) found that flows where $k_\beta /k_R \gtrsim 4$ demonstrated sharp eastward and broad westward jet structures. The greater this ratio, the more the jet structures approached the idealised potential vorticity staircase structure (Dritschel & McIntyre Reference Dritschel and McIntyre2008).

3. Geometric eddy ellipse formulation

The Reynolds decomposed shallow-water momentum equation is given by

(3.1)\begin{equation} \frac{\partial \bar{\boldsymbol{u}}}{\partial t} + \overline{f\boldsymbol{k}\times\boldsymbol{u}} ={-}\boldsymbol{\nabla}\boldsymbol{\cdot}\boldsymbol{\mathsf{T}} - g\,\boldsymbol{\nabla}\bar{\eta} + \bar{\boldsymbol{F}} + \bar{\boldsymbol{D}}, \end{equation}

where $\boldsymbol {u}\equiv (u,v)$ is the horizontal velocity, with zonal and meridional velocity components denoted $u$ and $v$, respectively, $f = f_0+\beta y$ is the Coriolis frequency under the $\beta$-plane approximation, $\boldsymbol {k}$ is the vertical unit vector, $g$ is the vertical component of the acceleration due to gravity, $\eta$ is the free surface of the fluid, $\boldsymbol {\nabla }$ is the horizontal gradient operator, and $\boldsymbol {F}$ and $\boldsymbol {D}$ are the horizontal forcing and dissipation mechanisms, respectively. We choose here to define the ‘eddying’ portion of the flow as fluctuations with respect to the zonal mean since jets formed in geostrophic turbulence are statistically zonal. The zonal mean is denoted with an overline, and primes indicate fluctuations. This is a crude assumption because in this fully developed turbulent system, there are many interacting scales of motion, so it is difficult to separate exactly what constitutes transient eddies and the long-lived zonal structure. This is something that we will explore. The advective terms vanish upon zonal averaging. We find that eddies influence the mean flow through the divergence of the eddy velocity correlation tensor:

(3.2)\begin{equation} \boldsymbol{\mathsf{T}} = \boldsymbol{\mathsf{T}}(y) = \begin{pmatrix} \overline{u'u'} & \overline{u'v'}\\ \overline{u'v'} & \overline{v'v'} \end{pmatrix}. \end{equation}

Equation (3.2) may be separated into isotropic (trace-only) and anisotropic (trace-free) parts:

(3.3)\begin{equation} \boldsymbol{\mathsf{T}} = K\boldsymbol{\mathsf{I}} + \boldsymbol{\mathsf{E}}, \end{equation}

where $\boldsymbol{\mathsf{I}}$ is the identity matrix, and the zonally averaged kinetic energy density associated with the eddying portion of the flow field is given by

(3.4)\begin{equation} K = K(y) = \tfrac{1}{2}(\overline{u'^2} + \overline{v'^2}), \end{equation}

and we have introduced the eddy momentum stress tensor

(3.5)\begin{equation} \boldsymbol{\mathsf{E}} = \begin{pmatrix} M & N\\ N & -M \end{pmatrix}. \end{equation}

The normal stress difference is given by

(3.6)\begin{equation} M = M(y) = \tfrac{1}{2}(\overline{u'^2} - \overline{v'^2}), \end{equation}

and the shear stress is given by

(3.7)\begin{equation} N = N(y) = \overline{u'v'}. \end{equation}

In a quasi-geostrophic, incompressible fluid system, $\bar {v} = 0$. Thus the dynamical evolution of (3.1) is governed by the divergence of the sheer stress:

(3.8)\begin{equation} \frac{\partial \bar{u}}{\partial t} ={-}\frac{\partial N}{\partial y} \end{equation}

in the absence of forcing and dissipation. The Reynolds stresses (3.6) and (3.7) describe the single-point correlations between eddy velocities $u'$ and $v'$. In the geometric eddy ellipse formulation, if there are strong correlations, then the mean pattern traced in $u'\unicode{x2013}v'$ space will be elliptical, and if no correlations exist, then the trace will be circular. We can recover these ellipse properties from the covariance tensor given by (3.2) by performing a principal axes decomposition (Preisendorfer Reference Preisendorfer1988), diagonalising (3.5) under a rotation through angle $\theta$, following Waterman & Lilly (Reference Waterman and Lilly2015). This is equivalent to rotating the ellipse such that its semi-major axis aligns with the zonal direction. This allows us to rewrite (3.5) as

(3.9)\begin{equation} \boldsymbol{\mathsf{E}} = L\begin{pmatrix} \cos{2\theta} & \sin{2\theta}\\ \sin{2\theta} & -\cos{2\theta} \end{pmatrix}, \end{equation}

where we have defined

(3.10)\begin{equation} L = L(y) = \sqrt{M^2+N^2}. \end{equation}

Here, $L$ is the excess energy in the direction of the major axis of the ellipse compared to the minor axis, and $\theta$ is given by

(3.11)\begin{equation} \theta = \theta(y) = \frac{1}{2}\arctan{\left(\frac{N}{M}\right)}, \end{equation}

which is the tilt of the ellipse with respect to the zonal direction. The geometric eddy ellipse is not a description of the shape of individual eddies themselves, but a description of their net fluxes. Through examining the distribution of eddy ellipses and observing their patterns, we can determine the direction and magnitude of these fluxes (Waterman & Hoskins Reference Waterman and Hoskins2013). Flows may demonstrate significant anisotropy between the zonal and meridional components of their flow velocities, such that $M\neq 0$, but this may not necessarily give rise to a zonal momentum tendency. This kind of behaviour corresponds to neutral ellipse tilts, where $N=0$. Conversely, when $M=0$, flow velocities do not demonstrate any significant anisotropy between the zonal and meridional components, but there may be a zonal momentum tendency where $N\neq 0$. An eddy ellipse tilted at $\theta = \pm ({{\rm \pi} }/{4})$ is achieved in this instance.

By examining the eddy ellipse shapes, we can tell which eddy distributions give rise to momentum fluxes into the mean flow, and which do not. An example of an eddy shape that demonstrates a zonal momentum flux is the banana-shaped eddy depicted by Wardle & Marshall (Reference Wardle and Marshall2000) (see their figure 5), where the eddy shape is such that momentum is fluxed eastwards. Examining the distribution of the geometric ellipses visualises these flux directions more intuitively than examining the velocities. Circular eddy ellipses are achieved only when $M=0$ and $N=0$, and there is no significant anisotropy, either through correlations between $u'$ and $v'$ or where $\overline {u'u'}\neq \overline {v'v'}$. In this scenario, only the isotropic component of the tensor given by (3.3) remains. Examples of other eddy ellipse distributions that give rise to zonal momentum fluxes can be found in figures 3 and 4 of Tamarin et al. (Reference Tamarin, Maddison, Heifetz and Marshall2016). In this paper, we will evolve our flow using the barotropic vorticity equation given by

(3.12)\begin{equation} \frac{\partial\zeta}{\partial t} ={-} u\,\frac{\partial \zeta}{\partial x} - v\, \frac{\partial \zeta}{\partial y} - \beta v + F - r\zeta - \mathcal{D}\,\nabla^4\zeta, \end{equation}

where $r$ and $\mathcal {D}$ are the Rayleigh friction and small-scale biharmonic dissipation parameters, respectively. Here, the zonal and meridional velocity components are given by

(3.13)\begin{equation} (u,v) \equiv \left(-\frac{\partial\psi}{\partial y},\frac{\partial\psi}{\partial x}\right), \end{equation}

written here in terms of the stream function $\psi$. The relative vorticity is given by

(3.14)\begin{equation} \zeta \equiv \left(\frac{\partial v}{\partial x} - \frac{\partial u}{\partial y}\right) \equiv \nabla^2\psi, \end{equation}

and we introduce rotation into our system through the $\beta$-plane. Rayleigh friction is chosen so that the flow reaches statistical equilibrium before the energy approaches the largest scales, such that characteristic scales are well separated according to (2.4).

The Reynolds-averaged inviscid unforced form of the barotropic vorticity equation is

(3.15)\begin{equation} \frac{\partial \bar{\zeta}}{\partial t} ={-}\boldsymbol{\nabla}\boldsymbol{\cdot}(\overline{\boldsymbol{u}'\zeta'}), \end{equation}

where eddies influence the mean flow in the barotropic vorticity equation through the divergence of an eddy vorticity flux $\overline {\boldsymbol {u}'\zeta '}$. From this, we can relate (3.2) to the zonal mean eddy vorticity flux in (3.15) through the Taylor–Bretherton identity (Taylor Reference Taylor1915; Bretherton Reference Bretherton1966; Plumb Reference Plumb1986):

(3.16)\begin{equation} \overline{\boldsymbol{u}'\zeta '} = \left(-\frac{\partial M}{\partial y}, -\frac{\partial N}{\partial y}\right). \end{equation}

Thus only the components of the eddy momentum stress tensor given by (3.5) appear explicitly in the eddy forcing of the mean barotropic vorticity equation.

4. Model

We will solve (3.12) numerically on a zonally periodic, meridionally bounded channel model with zonal and meridional dimensions $\pm L_x/2$ and $\pm L_y/2$, respectively, where $L_x = 6000\,\mathrm {km}$ and $L_y = 3000\,\mathrm {km}$. A laterally bounded model is chosen so that jets do not migrate meridionally. We must ensure that the appropriate boundary conditions are applied to maintain consistency between (3.12) and the shallow-water equations from which it is derived. We determine $\psi$ by inverting (3.14), and obtain $u,v$ using (3.13). There must be no flow through the meridional walls so $v=0$ such that $\psi '=0$. This imposes no conditions on the zonal mean flow on the boundary, so we are free to choose a mixture of boundary conditions: no-slip on the zonal mean modes, and free-slip boundary conditions on the non-mean modes (e.g. Esler Reference Esler2008; Williams et al. Reference Williams, Haine, Read, Lewis and Yamazaki2009) given by

(4.1)\begin{equation} \left.\begin{gathered} \frac{\partial\bar{\zeta}}{\partial y} = \frac{\partial^3\bar{\zeta}}{\partial y^3} = 0,\\ \zeta' = \frac{\partial^2 \zeta'}{\partial y^2} = 0. \end{gathered}\right\} \end{equation}

Solving for $\psi$ will also produce $\bar {\psi } = \mathrm {constant}$ on each of the lateral boundaries. These constants are, in general, different. The first is determined by setting $\bar {u}=0$ on each of the lateral boundaries, and the second is found through mass conservation. Setting $\bar {u}=0$, together with the boundary conditions (4.1), trivially satisfies consistency conditions between the barotropic vorticity equation and the shallow-water equations through conserving circulation (McWilliams Reference McWilliams1977; Graef & Müller Reference Graef and Müller1996).

We discretise (3.12) on a channel model with $513\times 257$ grid points such that the zonal and meridional grid spacing is given by $\Delta x = \Delta y = 11.7\,\mathrm {km}$, and the time step is $\Delta t = 112.5\,\mathrm {s}$. The flow is spun up from rest and forced using a stochastic isotropic forcing function (e.g. Chekhlov et al. Reference Chekhlov, Orszag, Sukoriansky, Galperin and Staroselsky1996) with a constant enstrophy injection rate $\xi = 8\times 10^{-18}\,\mathrm {s}^{-3}$ at the forcing scale given by $k_f = 64$, such that the constant energy injection rate is $\epsilon = 1.78\times 10^{-8}\,\mathrm {m}^2\,\mathrm {s}^{-3}$. The scale-selective biharmonic diffusion term is used to remove energy at small scales with a coefficient $\mathcal {D} = 5\times 10^{10}\,\mathrm {m}^4\,\mathrm {s}^{-1}$, and the flow is damped through Rayleigh friction with a coefficient $r = 1\times 10^{-8}\,\mathrm {s}^{-1}$. We choose $\beta = 4\beta _0$, where $\beta _0 = 2.29\times 10^{-11}\,\mathrm {m}^{-1}\,\mathrm {s}^{-1}$ is the terrestrial equatorial value. Though the $\beta$-value that we are using is much greater than the terrestrial value, the parameters chosen possess dynamic similarity to terrestrial problems. Note that all wavenumbers quoted have been scaled by $2{\rm \pi} /L_y$.

5. A zonostrophic jet

The parameters for this barotropic channel model have been chosen in such a way that the flow is in the zonostrophic regime and exhibits multiple zonal jets, obeying the scaling relationships presented in Sukoriansky et al. (Reference Sukoriansky, Dikovskaya and Galperin2007). In this section, we will discuss these characteristic scales in detail, and show how they relate to the structure and stability of the jets.

We observe the evolution of the flow as it is forced isotropically and spun up from rest, eventually approaching its equilibrium jet configuration (figure 1). Early in the flow evolution, zonal structures emerge that merge together as time evolves. What we are essentially observing is the evolution of the Rhines scale as the system increases in energy and $U$ increases. Using (2.3), we calculate the final destination of the Rhines scale, the mode number $k_{eq} \approx 4$. This agrees well with the observed number of peak to peak jets, in which there are $5$ sharp eastward jets interleaved by $4$ broad westward jets. We also note a subtle interplay between the energy and jet stability. As noted in Vallis & Maltrud (Reference Vallis and Maltrud1993), sharp eastward and broad westward jet profiles correspond to barotropically stable structures. Rearranging a jet configuration then requires breaking the stability, and we see a fast change in energy as this occurs. So as the flow progresses, there are a few metastable jet configurations that the flow assumes, before approaching its final configuration.

Figure 1. (a) Hovmöller diagram of $\bar {u}$ normalised by the instantaneous maximum zonal velocity. Superimposed is the evolution of the domain average energy as a function of time (red). (b) The zonal time average velocity profile, $[\bar {u}]$ (thick black) and the time average divergence of momentum flux, $-r^{-1}\,\partial _y N$ (thin orange). The vertical black line in (a) marks the beginning of the $20\,000\,\mathrm {h}$ period over which the time averages in (b) are taken.

The tendency in the zonal mean flow is governed by a divergence of the shear stress (3.8). Upon introducing Rayleigh friction, the flow will equilibrate as energy is removed from the largest scales. We introduce this in (5.2), by setting

(5.1)\begin{equation} \bar{\boldsymbol{D}} ={-}r\bar{\boldsymbol{u}}. \end{equation}

Since $\bar {v}=0$, the dynamics is governed solely by the zonal component of the shallow-water equation,

(5.2)\begin{equation} \frac{\partial\bar{u}}{\partial t} ={-} \frac{\partial N}{\partial y} -r\bar{u}. \end{equation}

When the flow reaches a steady state,

(5.3)\begin{equation} \bar{u} ={-} \frac{1}{r}\,\frac{\partial N}{\partial y}. \end{equation}

So the jets will be maintained by the divergence of the shear stresses. Of course, in our case we also have a biharmonic diffusion term, so this relationship is approximate. We examine this relationship in figure 1(b). The time average (denoted by square brackets) is taken over the last $20\,000\,\textrm {h}$ of the flow evolution. Over this period, the divergence in the shear stress broadly matches the zonal mean structure with some small variation. Shorter time intervals produce less agreement, with the flow demonstrating short-time-scale variations that produce a tendency in the mean jet profile such that the relationship given by (5.3) no longer holds. The fact that the divergence of the shear stress produces steady jet patterns after a long-time mean suggests that a scale separation may exist between scales producing the jet patterns and the scales producing temporal variations.

We now demonstrate that these jets are characteristic of those in the literature and lie within the zonostrophic regime. The energy spectrum $E(k_x,k_y)$ demonstrates a cascade of energy towards the largest scales (figure 2a), and shows the characteristic ‘dumbbell’ shape that was first described in Maltrud & Vallis (Reference Maltrud and Vallis1991). This is associated with the anisotropisation of the flow due to the presence of Rossby waves. A large amount of energy is seen to have built up along the $k_y$-axis, which can be understood by examining the angularly averaged spectra (figure 2b).

Figure 2. (a) The time average energy spectrum $\log E(k_x,k_y)$. (b) The time and angular averaged spectra: total $E(k)$ (solid black), zonal $E_Z(k)$ (solid orange), and residual $E_R(k)$ (turquoise). The vertical dashed lines from right to left are the forcing scale $k_f$ (plum red), cross-over wavenumber $k_\beta$ (blue), and Rhines scale $k_R$ (green). The jet profile spectrum is marked with a thin yellow line. The slopes of the theoretical residual spectrum ($\propto k^{5/3}$, dashed red) and zonal spectrum ($\propto k^{-5}$, dot-dashed red) are also given. The time average is taken over a $20\,000\,\mathrm {h}$ equilibrium period.

The total spectrum (figure 2b), initially traces the Kolmogorov slope given by $k^{-5/3}$, associated with isotropic dynamics. We calculate the cross-over wavenumber given by (2.2): $k_\beta \approx 16$. We find that as the cascade approaches $k_\beta$, $E(k)$ steepens and possesses several sharp spikes of energy that trace the $k^{-5}$ slope, confirming this as the approximate point at which the spectrum anisotropises. This corresponds to the energy build up along the $k_y$-axis in $E(k_x,k_y)$ (figure 2a).

In splitting the spectrum into the zonal and residual spectra, we find that at smaller scales, most of the energy is concentrated in the residual spectrum, i.e. the region of $E(k_x,k_y)$ that lies outside of the $k_y$-axis, and traces the $k^{-5/3}$ slope. The zonal spectrum also initially traces the $k^{-5/3}$ slope, which is thought to be indicative of the absence of Rossby waves at these scales. The zonal spectrum steepens and surpasses the residual spectrum at $k \approx k_\beta$. At the largest scales, the energy contained in the zonal modes dwarfs that of the residual such that the amplitude of the zonal spectrum is approximately equal to that of the total spectrum (figure 2a). This possesses sharp, discrete spikes of energy that are harmonics of the largest energy-containing scales, i.e. mode numbers $4$ and $5$, corresponding to the number of westward and eastward jets, respectively (Danilov & Gurarie Reference Danilov and Gurarie2004).

Examining figure 2(b), we see that all characteristic scales of interest (i.e. $k_R\approx 4$, $k_\beta \approx 16$ and $k_f=64$) are well separated according to (2.4). However, we note that the biharmonic diffusion causes a strong suppression of the enstrophy cascade throughout the inertial range, so the last inequality in (2.4) is not obeyed by this system.

We also examine the jet profile spectrum, which is the spectrum corresponding to the zonal mean jet structure (figure 1b). Here, we see that the spectrum associated with the jet profile possesses the discrete harmonics observed in the zonal and total energy spectrum (figure 2b). The asymmetric sharp eastward and broad westward jets are required for the jet profile to be barotropically stable (Vallis & Maltrud Reference Vallis and Maltrud1993); we argue that these discrete harmonics are then a necessary consequence of barotropic stability.

6. Eddy tilts

In the previous section, we have seen how the parameters chosen have led to jets forming in the zonostrophic regime, with characteristic scales adequately well separated according to (2.4) and an equilibrium jet profile consisting of sharp eastward and broad westward jets representing a barotropically stable structure. We will now turn our attention to examining eddy–mean-flow interactions using the geometric eddy ellipse formulation, attempting to identify the scales responsible for fluxing momentum into the jet structure, balancing the mean flow.

Tamarin et al. (Reference Tamarin, Maddison, Heifetz and Marshall2016) examined a barotropically unstable zonal jet profile that successively strengthened or weakened as instabilities on the jet flanks decayed and grew, respectively. Here, the jet structures that form as a result of spectral anisotropisation by the $\beta$ effect have barotropically stable jet profiles that satisfy the Rayleigh–Kuo stability criterion. We expect the eddy ellipse patterns corresponding to the jets to persist in time in order to maintain the jet structure. It is peculiar that we observe the relationship given by (5.3) only after a relatively long time averaging, when the jet structures themselves are persistent and stable. To understand this, we examine the Hovmöller plots of the momentum fluxes $N, M$ and the ellipse quantity $\theta$ in figure 3 over a $20\,000\,\textrm {h}$ period in which the energy is in a statistical equilibrium. Though the energy has reached equilibrium, we see that there are significant time-dependent processes occurring in the shear stress $N$, where momentum is fluxed north and south in an alternating pattern (figure 3a). The tendency in $\bar {u}$ is governed by (5.2), so these time-dependent processes result in short-term fluctuations in jet intensity. We see, however, that over a time average, the shear stress possesses a pattern that is consistent with the jet structure (figure 3b). In particular, we find that eastward jets are flanked by a negative flux of eastward momentum to their north and a positive flux of eastward momentum to their south. The converse is true for westward jets. This leads to a characteristic pattern in the stress that is half a jet width offset from the jet structure.

Figure 3. Hovmöller diagrams of: (a) the shear stress $N(y)$ and (b) its normalised time average $[N(y)]/10^{-3}$ (blue solid); (c) the normal stress difference $M(y)$ and (d) its normalised time average $[M(y)]/10^{-2}$ (blue solid); and (e) the corresponding ellipse tilt angle $\theta (y)$ and (f) its time average $[\theta ]$ normalised by ${\rm \pi} /2$ (blue solid). The normalised zonal-time average zonal velocity $[\bar {u}]/U_{max}$ is also plotted for comparison (black dashed), where $U_{max}$ is the maximum value of $[\bar {u}]$. The time averages are taken over a $20\,000\,\mathrm {h}$ equilibrium period.

The difference in normal stresses $M$ is in general an order of magnitude larger than $N$ (figure 3c). The temporal evolution of $M$ shows a zonal structure that is time-dependent and pulsing regularly. The pulses are always positive and coincide with regions where the $N$ changes sign in time. These pulses are located at some of the jet cores and at domain walls, but do not appear to correlate with any strengthening pattern in the zonal mean structure itself. We speculate that the intermittent signal of pulses in $M$ can be attributed to the action of the Rayleigh friction, which becomes important as the jet strengthens. In both the time evolution and the time average, $M$ is always non-negative, such that $\overline {u'^2}> \overline {v'^2}$. Also, although there is a zonal jet pattern in the normal stress field, this structure does not match the zonal mean jet pattern in $\bar {u}$ (figure 3d).

The resulting Hovmöller plot of the eddy tilt angle $\theta$, calculated from $M$ and $N$ using (3.11), does not look particularly instructive (figure 3e). It demonstrates a time alternating pattern of net northward and net southward momentum fluxing which possesses a meridional structure that is strongest at the jet cores, though it is not obvious how these fluxes interact with the jet. However, the time average shows a very clear eddy–mean-flow interaction, consistent with the jet pattern, where eddy ellipses are tilted towards the jet structure in an alternating pattern on the flanks of the jets (figure 3f). We see that there are time-dependent momentum fluxes that dominate the flow as discovered in Huang & Robinson (Reference Huang and Robinson1998), which mask a more subtle momentum flux that supports the jet structure. The flow converges to this structure only approximately over a 20 000 h period. Our task is therefore to find exactly what scales are responsible for supporting the jet.

6.1. Zonal filter

In the previous section, we calculated a number of characteristic scales. We found that scales where $k > k_\beta$ are dominated by energy in the residual spectrum (2.1b), and scales where $k < k_\beta$ are dominated by energy accumulation along the $k_y$-axis corresponding to the development of the zonal spectrum (2.1a). The latter we associate with the formation of jets. We also know that in the steady state, the number of jets is equal to the Rhines scale $k_R$. Consequently, we might expect that eddies supporting the jet structures will reside between the Rhines scale $k_R$ and the anisotropy scale $k_\beta$. Since the flow favours zonal modes in the zonostrophic regime, we will develop a filter that separates the velocity correlations that contribute to the eddy stresses $M$ and $N$ by zonal scale.

We can exploit the homogeneity in $x$ to find the contributions of eddy velocities of different zonal scales to the eddy stress quantities $M$ and $N$. Since our flow is periodic in $x$, we take a Fourier transform of the eddy velocities in the $x$-direction, and substitute these into the eddy velocity correlation tensor (3.2). Thus

(6.1) \begin{equation} \boldsymbol{\mathsf{T}}(y) = \frac{2}{L_x}\begin{pmatrix} \displaystyle\int_1^{k_D} |\widetilde{u'}|^2 \,{\rm d} k_x & \displaystyle\int_1^{k_D} \widetilde{u'}\widetilde{v'}^\ast \,{\rm d} k_x\\ \displaystyle\int_1^{k_D} \widetilde{v'}\widetilde{u'}^\ast \,{\rm d} k_x & \displaystyle\int_1^{k_D} |\widetilde{v'}|^2 \,{\rm d} k_x \end{pmatrix}, \end{equation}

where the tilde represents the Fourier transform in the zonal direction, the asterisk is the complex conjugate, and we take the limit of integration from $k_x=1$ as the eddying velocities denoted by a prime have the zonal mean subtracted. From this, we can rewrite $M$ as

(6.2a)$$\begin{gather} M = M_l + M_h, \end{gather}$$

(6.2b)$$\begin{gather}\text{where}\quad M_l =\frac{1}{L_x}\int_1^{k_c-1} \left(\widetilde{|u'|}^2 - \widetilde{|v'|}^2\right){\rm d} k_x \end{gather}$$

(6.2c)$$\begin{gather}\text{and}\quad M_h = \frac{1}{L_x}\int_{k_c}^{k_D} \left(\widetilde{|u'|}^2 - \widetilde{|v'|}^2\right){\rm d} k_x, \end{gather}$$

and $N$ as

(6.3a)$$\begin{gather} N = N_l + N_h, \end{gather}$$

(6.3b)$$\begin{gather}\text{where}\quad N_l = \frac{1}{L_x}\int_0^{k_c-1} \widetilde{u'}\widetilde{v'}^\ast \,{\rm d} k_x \end{gather}$$

(6.3c)$$\begin{gather}\text{and}\quad N_h = \frac{1}{L_x}\int_{k_c}^{k_D} \widetilde{u'}\widetilde{v'}^\ast \,{\rm d} k_x, \end{gather}$$

where $k_c$ is some threshold wavenumber, and the subscripts $l$ and $h$ represent low-pass and high-pass filters, respectively.

We examine the scales associated with the eddy–mean-flow relationship given by (5.3), separating $N$ into low-pass ($N_l$) and high-pass ($N_h$) filtered signals as in (6.3a), for different cut-off wavenumbers $k_c$ (figure 4). The low-pass filter on $N$ is only able to capture the jet pattern in its divergence for large values of $k_c$, with the structure beginning to emerge at relatively high wavenumbers $k_c \geq 32$ (figure 4a). Below this, the signal has an amplitude comparable to the jet scale, but the pattern does not match the jet profile. In contrast, the jet pattern is evident for nearly all $k_c$, in the high-pass filter, with a faint, low-amplitude jet pattern still present when wavenumbers $k\geq k_f$ are retained (figure 4b). As $k_c$ varies, the transition between the noisy low wavenumbers and the appearance of the jet pattern in the high wavenumbers is relatively smooth. Though we cannot pinpoint a sharp transition scale, what is clear is that the divergence of $N$ does not capture the jet pattern for scales where $k_c \leq k_\beta$. This is counterintuitive because this is the scale associated with spectral anisotropy; when $k>k_\beta$, it is thought that isotropic dynamics dominate. Choosing $k_c=k_\beta$ does appear to be somewhat representative of the transition between the two regimes of eddy–mean-flow interactions, so we will proceed using this choice.

Figure 4. The time-average zonal velocity profile $[\bar {u}]$ over the $20\,000\,\mathrm {h}$ equilibrium period (black dashed) compared to (a) low-pass filtered divergence of the shear stress, $-r^{-1}\partial _yN_l$, and (b) high-pass filtered divergence of the shear stress, $-r^{-1}\partial _yN_h$. Here, $N_l$ is the low-pass filter on $N$, retaining zonal wavenumbers from $k_x=1$ up to cut-off wavenumber $k_c = 4,8,16,32,64,128,256$ (thin light grey to thick dark grey); $N_h$ is the corresponding high-pass filter from the same cut-off wavenumber $k^c_x$ to $k_x=256$ (thick dark grey to thin light grey) such that $N = N_l+N_h$ for each $k_c$. The special case $k_c = k_\beta = 16$ is shown as the thick blue line.

In figure 5 we examine the Hovmöller diagrams for the low-pass filtered stresses $M_l$ and $N_l$, and the associated tilt angle; these appear almost identical to the diagrams for the unfiltered quantities in figure 3, where they are dominated by noisy, time-dependent processes.

Figure 5. Hovmöller diagrams of the low-pass filtered, large zonal scale: (a) shear stress $N_l(y)$ and (b) its normalised time average $[N_l(y)]/10^{-3}$ (blue solid); (c) normal stress difference $M_l(y)$ and (d) its normalised time average $[M_l(y)]/10^{-2}$ (blue solid); and (e) the corresponding ellipse tilt angle $\theta _l(y)$ and (f) its time average $[\theta _l]$ normalised by ${\rm \pi} /2$ (blue solid). The normalised zonal-time average zonal velocity $[\bar {u}]/U_{max}$ is also plotted for comparison (black dashed), where $U_{max}$ is the maximum value of $[\bar {u}]$. The time averages are taken over a $20\,000\,\mathrm {h}$ equilibrium period.

The Hovmöller diagram of $N_l$ demonstrates the same alternating pattern of momentum flux observed in the shear stress. Unlike the unfiltered case, the time mean over $N_l$ (figure 5b) does not average to reveal a correlation with the jet structure. Thus we can confirm that the scales retained in the low-pass filter are not responsible for the time-mean structure in $N$ (figure 3a). The zonal pulses observed in the Hovmöller diagrams of $M$ (figure 3c) are present in $M_l$ (figure 5c). We see in the time mean of $M_l$ (figure 5d) that these scales provide the zonal structure observed in the time mean of $M$ (figure 3d). As we have noted before, these pulses are not entirely consistent with the jet pattern, though they do broadly resemble some jet structure. They do, however, seem to correlate with the alternating streaking pattern in $N_l$ (figure 5a). The Hovmöller diagram of $\theta _l$ (figure 5e) demonstrates the same alternating north and south flux direction found in $\theta$ (figure 3e), with no discernible tilting pattern in the time mean (figure 5f).

Momentum flux quantities filtered at high pass, where wavenumbers $k \geq k_c = k_\beta$ are retained (figure 6a), demonstrate the momentum flux structure that we have been searching for.

Figure 6. Hovmöller diagrams of the high-pass filtered, small and intermediate zonal scale: (a) shear stress $N_h(y)$ and (b) its normalised time average $[N_h(y)]/10^{-3}$ (blue solid); (c) normal stress difference $M_h(y)$ and (d) its normalised time average $[M_h(y)]/10^{-2}$ (blue solid); and (e) the corresponding ellipse tilt angle $\theta _h(y)$ and (f) its time average $[\theta _h]$ normalised by ${\rm \pi} /2$ (blue solid). The normalised zonal-time average zonal velocity, $[\bar {u}]/U_{max}$ is also plotted for comparison (black dashed), where $U_{max}$ is the maximum value of $[\bar {u}]$. The time averages are taken over a $20\,000\,\mathrm {h}$ equilibrium period.

From the Hovmöller diagram of $N_h$, we see that high-pass filtering reveals an alternating band structure that does not vary in time, with the exception of some noise (figure 6a). This zonal structure is consistent with the zonal jets. The strengths of these bands are much smaller than the magnitudes of the streaks in $N_l$ (figure 5a), which is why they are not apparent in the unfiltered signal $N$. We can see that these bands are a half-jet-width offset from the jet profile, providing the zonal structure observed in the time mean of $N$ with a magnitude that supports the jets (figure 6b).

The Hovmöller diagram of $M_h$ demonstrates small, negative, persistent band structures that align with the jet structure (figure 6c). The sign of $M_h$ is negative everywhere inside the domain, and is positive only at the domain walls. Also, $M_h$ has its highest amplitude, and is most negative within the jet cores. The pattern in the time mean correlates with the jet structure, though the jet structure is purely zonal. This suggests that the eddies at the scales supporting the jet have strong meridional velocity fluxes (figure 6d).

The high-pass tilt angle $\theta _h$ demonstrates a persistent banded pattern consistent with the jet (figure 6e). The pattern in the time mean is consistent with the jet structure (figure 6f) and also reveals some more subtle features. Jets must demonstrate a sharp eastward and broad westward jet structure for barotropic stability (Vallis & Maltrud Reference Vallis and Maltrud1993). A sharp eastward jet requires a sharp change in tilt angle, changing from strongly negative on its northern flank to strongly positive to the south of the jet maximum. A broad westward jet requires a less sharp change in the tilt angle, gradually changing from positive to negative from the north to the south. This results in the sawtooth shape observed in the time-mean structure of the tilt angle.