3.1 Drift velocities from the linear dispersion relation

Although the jet tilting is associated with the tilting of barotropic PV isolines, it is not clear why the jets drift in the presence of a zonal component of PV gradient. In order to address this aspect, we analysed the dispersion relation derived from the two-layer QG model (from (2.1) and (2.2), see appendix B for details) and found frequencies (real parts of the eigenvalues) for eigenmodes matching the jets in the reference solutions. The eigenmodes possess non-zero frequencies, which means that they propagate. The corresponding phase velocities are then compared with the drift velocities estimated from the reference solutions (table 1). The predictions from the dispersion relation are quite close to the fully nonlinear solutions, although the nonlinearities seem to slow the jets down, as the observed drift speeds are slightly smaller than the linear predictions. The difference increases with zonal slope and is almost 30 % in the large-slope simulation. On the basis of these comparisons, we conclude that the linear dispersion relation can be used to predict jet drift speeds, at least for mild zonal slopes. However, the prediction from the linear dispersion relation does not work very well in cases where $O(T_{x})\rightarrow O(\unicode[STIX]{x1D6FD})$ . This comparison is possible because linear dynamics governs a major part of the flow, and the jets correspond to the eigenmodes, which are solutions of the linear problem (Berloff & Kamenkovich Reference Berloff and Kamenkovich2013a ,Reference Berloff and Kamenkovich b ). The main difference is that, unlike zonal eigenmodes over a flat bottom, tilted eigenmodes are not steady but propagate with constant velocities. As we show in the next section, eddies also contribute to the drift by slowing the jets down.

Figure 4.

Sketch of the new coordinate system (new axes are $p$ and $q$ ) that is rotated by an angle $\unicode[STIX]{x1D703}$ (tilt is positive for anticlockwise rotation) and is moving with a speed $c$ (drift speed, positive in $q$ direction). Angle $\unicode[STIX]{x1D703}$ and speed $c$ are directly estimated from the numerical solutions. Green arrows represent the alternating jets; note that the time-mean profile of the jets is a function of $q$ only.

Table 1.

Comparison of drift velocities ( $\text{cm}~\text{s}^{-1}$ ) in the numerical simulations and dispersion relation. Overall, the comparison is good but there is a systematic difference for steeper slopes. $\hat{\boldsymbol{i}}$ and $\hat{\boldsymbol{j}}$ are the zonal and meridional unit vectors, respectively.

3.2 Nonlinear effects and the role of eddies

In the case of purely zonal and stationary jets, the flow can be Reynolds-style decomposed into the time-mean large-scale flow and transient eddies, and the time-mean (or zonally averaged) meridional profiles of nonlinear eddy stress terms can be analysed (e.g. Panetta Reference Panetta1993; Berloff et al. Reference Berloff, Kamenkovich and Pedlosky2009a ). Through a continuous baroclinic instability process, the eddies gain energy from the available potential energy of the system and interact with the jets via exchanges of momentum and PV, and one can study this process by comparing and correlating cross-jet profiles of the eddy forcing and time-mean PV. However, in the case of drifting jets, the analysis is not that straightforward, as the mean state itself is time-dependent. To overcome this problem, we used a new frame of reference in which the jets are stationary and aligned with one of the coordinates. In the new configuration, shown in figure 4, the frame of reference is rotated by the jet tilt angle $\unicode[STIX]{x1D703}$ and moves with a constant speed  $c$ (equal to the jet drift speed). Note that PV conservation holds true irrespective of the frame of reference, as long as there is no additional acceleration. Thus, the governing equations (2.1) and (2.2) can be rewritten in the new frame of reference and interpreted as PV conservation laws (note that $\unicode[STIX]{x1D713}_{i}=\unicode[STIX]{x1D713}_{i}(p,q,t)$ ):

(3.1) $$\begin{eqnarray}\displaystyle & \displaystyle \left[\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}t}-c\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}q}\right]\unicode[STIX]{x1D6F1}_{1}=-\left(\frac{\unicode[STIX]{x2202}\unicode[STIX]{x1D713}_{1}}{\unicode[STIX]{x2202}p}-U_{b}\sin \unicode[STIX]{x1D703}\right)\frac{\unicode[STIX]{x2202}\unicode[STIX]{x1D6F1}_{1}}{\unicode[STIX]{x2202}q}+\left(\frac{\unicode[STIX]{x2202}\unicode[STIX]{x1D713}_{1}}{\unicode[STIX]{x2202}q}-U_{b}\cos \unicode[STIX]{x1D703}\right)\frac{\unicode[STIX]{x2202}\unicode[STIX]{x1D6F1}_{1}}{\unicode[STIX]{x2202}p}+\unicode[STIX]{x1D708}\unicode[STIX]{x1D6FB}^{4}\unicode[STIX]{x1D713}_{1}, & \displaystyle \nonumber\\ \displaystyle & & \displaystyle\end{eqnarray}$$

(3.2) $$\begin{eqnarray}\displaystyle & \displaystyle \left[\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}t}-c\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}q}\right]\unicode[STIX]{x1D6F1}_{2}=-\frac{\unicode[STIX]{x2202}\unicode[STIX]{x1D713}_{2}}{\unicode[STIX]{x2202}p}\frac{\unicode[STIX]{x2202}\unicode[STIX]{x1D6F1}_{2}}{\unicode[STIX]{x2202}q}+\frac{\unicode[STIX]{x2202}\unicode[STIX]{x1D713}_{2}}{\unicode[STIX]{x2202}q}\frac{\unicode[STIX]{x2202}\unicode[STIX]{x1D6F1}_{2}}{\unicode[STIX]{x2202}p}+\unicode[STIX]{x1D708}\unicode[STIX]{x1D6FB}^{4}\unicode[STIX]{x1D713}_{2}-\unicode[STIX]{x1D6FE}\unicode[STIX]{x1D6FB}^{2}\unicode[STIX]{x1D713}_{2}, & \displaystyle\end{eqnarray}$$

where

(3.3) $$\begin{eqnarray}\left.\begin{array}{@{}l@{}}\displaystyle \unicode[STIX]{x1D6F1}_{1}=\overbrace{\unicode[STIX]{x1D6FB}^{2}\unicode[STIX]{x1D713}_{1}+S_{1}(\unicode[STIX]{x1D713}_{2}-\unicode[STIX]{x1D713}_{1})}^{Q_{1}}+\,(\unicode[STIX]{x1D6FD}+S_{1}U_{b})(q\cos \unicode[STIX]{x1D703}+p\sin \unicode[STIX]{x1D703}+ct\cos \unicode[STIX]{x1D703}),\\ \displaystyle \unicode[STIX]{x1D6F1}_{2}=\overbrace{\unicode[STIX]{x1D6FB}^{2}\unicode[STIX]{x1D713}_{2}+S_{2}(\unicode[STIX]{x1D713}_{1}-\unicode[STIX]{x1D713}_{2})}^{Q_{2}}+\,(\unicode[STIX]{x1D6FD}-S_{2}U_{b})(q\cos \unicode[STIX]{x1D703}+p\sin \unicode[STIX]{x1D703}+ct\cos \unicode[STIX]{x1D703})\\ \qquad +\,T_{x}(-q\sin \unicode[STIX]{x1D703}+p\cos \unicode[STIX]{x1D703}-ct\sin \unicode[STIX]{x1D703}).\end{array}\right\}\end{eqnarray}$$

Here, $\unicode[STIX]{x1D6F1}_{i}$ represents the absolute PV (i.e. the sum of PV due to the generated flow $Q_{i}$ and the background PV) in individual layers. Also note that all variables are functions of $p,~q,~t$ and $\unicode[STIX]{x1D735}=(\unicode[STIX]{x2202}/\unicode[STIX]{x2202}p,\unicode[STIX]{x2202}/\unicode[STIX]{x2202}q)$ . The first two terms on the right-hand side of (3.1) and (3.2) contain both linear and nonlinear terms, which are due to the background PV gradients, generated flow and nonlinear eddy effects. Depending on the overall sign of the sum, $Q_{i}$ can be forced or damped, resulting in strengthening or weakening of the jets. The rest of the terms on the right-hand side always smooth the PV gradients by removing energy from the system through viscous dissipation and bottom drag. In order to understand the effects of eddies and linear terms, we derived the time-mean equations from (3.1) and (3.2) by representing the variables as a sum of the time-mean flow ( $\overline{Q}_{i}$ ,  $\overline{\unicode[STIX]{x1D713}}_{i}$ ) and time-dependent eddy field ( $Q_{i}^{\prime }$ ,  $\unicode[STIX]{x1D713}_{i}^{\prime }$ ). This is basically the Reynolds flow decomposition in the moving frame of reference. The equations are then averaged over time (see full derivation in appendix C):

(3.4) $$\begin{eqnarray}\displaystyle & \displaystyle \frac{\unicode[STIX]{x2202}\overline{Q}_{1}}{\unicode[STIX]{x2202}t}=\underbrace{-\unicode[STIX]{x1D735}\boldsymbol{\cdot }\overline{(\boldsymbol{u}_{1}^{\prime }Q_{1}^{\prime })}}_{\text{eddy forcing}}+\underbrace{\overbrace{(c+U_{b}\sin \unicode[STIX]{x1D703})\frac{\unicode[STIX]{x2202}\overline{Q}_{1}}{\unicode[STIX]{x2202}q}}^{L_{1}^{c}}\overbrace{-\,(\unicode[STIX]{x1D6FD}+S_{1}U_{b})\sin \unicode[STIX]{x1D703}\,\,\overline{u}_{1}}^{L_{1}^{m}}}_{\text{linear stress terms}}+\underbrace{\overbrace{\unicode[STIX]{x1D708}\unicode[STIX]{x1D6FB}^{2}\overline{\unicode[STIX]{x1D701}}_{1}}^{D_{1}^{\unicode[STIX]{x1D708}}}}_{\text{dissipation}},\qquad & \displaystyle\end{eqnarray}$$

(3.5) $$\begin{eqnarray}\displaystyle & \displaystyle \frac{\unicode[STIX]{x2202}\overline{Q}_{2}}{\unicode[STIX]{x2202}t}=\overbrace{-\unicode[STIX]{x1D735}\boldsymbol{\cdot }\overline{(\boldsymbol{u}_{2}^{\prime }Q_{2}^{\prime })}}^{\text{eddy forcing}}+\overbrace{\underbrace{c\frac{\unicode[STIX]{x2202}\overline{Q}_{2}}{\unicode[STIX]{x2202}q}}_{L_{2}^{c}}\underbrace{-\,[(\unicode[STIX]{x1D6FD}-S_{2}U_{b})\sin \unicode[STIX]{x1D703}+T_{x}\cos \unicode[STIX]{x1D703}]\,\overline{u}_{2}}_{L_{2}^{m}}}^{\text{linear stress terms}}+\overbrace{\underbrace{\unicode[STIX]{x1D708}\unicode[STIX]{x1D6FB}^{2}\overline{\unicode[STIX]{x1D701}}_{2}}_{D_{2}^{\unicode[STIX]{x1D708}}}\underbrace{-\,\unicode[STIX]{x1D6FE}\overline{\unicode[STIX]{x1D701}}_{2}}_{D_{2}^{\unicode[STIX]{x1D6FE}}}}^{\text{dissipation}},\qquad & \displaystyle\end{eqnarray}$$

where $\boldsymbol{u}_{i}$ is the new velocity ( $\boldsymbol{u}_{i}=(u_{i},v_{i})=(-\unicode[STIX]{x2202}\unicode[STIX]{x1D713}_{i}/\unicode[STIX]{x2202}q,\unicode[STIX]{x2202}\unicode[STIX]{x1D713}_{i}/\unicode[STIX]{x2202}p)$ ) and $\unicode[STIX]{x1D701}_{i}$ is the new relative vorticity ( $\unicode[STIX]{x1D701}_{i}=\unicode[STIX]{x1D6FB}^{2}\unicode[STIX]{x1D713}_{i}$ ) in each layer. Here, we retain the time derivatives in the time-mean PV (in the rest of the paper, ‘mean PV’ refers to the mean PV of the developed flow only i.e.  $\overline{Q}_{i}$ ), but the terms vanish in the steady state and are kept here only for clarity. As seen in (3.4) and (3.5), the nonlinear terms (referred to hereafter as ‘eddy forcing’ or ‘nonlinear stress term’) and the linear terms (‘linear stress terms’ $L_{i}^{c}+L_{i}^{m}$ , and ‘dissipation terms’ $D_{i}^{\unicode[STIX]{x1D708}}+D_{i}^{\unicode[STIX]{x1D6FE}}$ ) balance each other in the stationary equilibrated state.

Figure 5.

Cross-jet profiles for the medium-slope simulation averaged in the moving frame in the interval 10 000–20 000 days. (a,b) Eddy forcing ( $NL_{i}$ , smoothed with 10-point moving average) and linear stress terms plus dissipation terms ( $L_{i}$ ). (c,d) Radiation stress term, inertial stress term and dissipation terms. (e,f) Reynolds stress term and form stress term (smoothed with 10-point moving average). Top and bottom rows are for the top and bottom layers, respectively. In order to normalise the plots, the stress terms (units are in $\text{s}^{-2}$ ) are multiplied by factors of $5\times 10^{12}$ and $2\times 10^{13}$ in the top and bottom layers, respectively. Dashed blue curves in (a,b,e,f) represent the time-mean PV profiles ( $\overline{Q}_{i}$ , normalised), which are multiplied by factors of $2.2\times 10^{5}$ and $8.8\times 10^{5}$ in the top and bottom layers, respectively. The eddies act against the jets, whereas the jets are forced by the linear stress terms.

In order to compute the eddy forcing and linear stress terms, we first linearly interpolated the streamfunction field on a rectangular grid rotated by the jet tilt angle $\unicode[STIX]{x1D703}$ and aligned the streamfunction snapshots in the moving frame of reference at every time step. We then computed the stress terms using the relations (3.4) and (3.5) and averaged them over time and along the jets. Cross-jet profiles of these linear and nonlinear stress terms for the medium-slope solution are shown in figure 5. As expected from the equations, the linear and nonlinear stress terms balance each other (figure 5 a,b). (There are small deviations from zero due to the data interpolation. An interpolation scheme with higher-order accuracy could be used, but we found linear interpolation to be sufficient for our purposes. In fact, we also computed the stress terms in the original $xy$ coordinate system and interpolated the stress field, but the resulting differences were found to be insignificant.) Relative contributions of the individual linear terms are also analysed (see figure 5 c,d): the $L^{c}$ and $L^{m}$ terms (hereafter referred to as ‘radiation stress term’ and ‘inertial stress term’, respectively) in the top layer are almost four times larger in magnitude than the eddy forcing, whereas all these terms are of comparable strength in the bottom layer. In the top layer, the radiation and inertial stress terms largely balance each other, and the difference between them is compensated by the eddy forcing. In both layers, the dissipation terms are significantly small and can be safely neglected. In terms of kinetic energy, the top layer dominates over the bottom layer and is two to three times more energetic (not shown). Thus, this supports the hypothesis that linear modes control a major part of the flow dynamics and hence the drift.

Table 2.

Correlation coefficients between cross-jet profiles of the time-mean linear and nonlinear stress terms, and the time-mean PV ( $\overline{Q}_{i}$ ) in both layers for the medium-slope simulation. The linear stress terms support the jets, whereas the eddy forcing opposes them.

It is important to note that the eddy forcing tends to smooth the mean PV profile in both layers and counteract the linear stress terms (figure 5 a,b). We computed correlation coefficients between the cross-jet profiles of the time-mean PV and stress terms. The linear stress terms are positively correlated (correlation coefficients ${\sim}0.6$ ) with the time-mean PV, whereas the eddy forcing terms show negative correlation ( ${\sim}-0.3$ ) with the time-mean PV (table 2). This is in contrast to previous studies on oceanic jets (Berloff et al. Reference Berloff, Kamenkovich and Pedlosky2009a ,Reference Berloff, Kamenkovich and Pedlosky b ; Chen et al. Reference Chen, Kamenkovich and Berloff2015) and jets in planetary atmospheres (Ingersoll et al. Reference Ingersoll, Beebe, Mitchell, Garneau, Yagi and Müller1981; Panetta Reference Panetta1993; Ingersoll et al. Reference Ingersoll, Gierasch, Banfield, Vasavada and Team2000), where the jets are essentially eddy-driven. To clarify this matter, we decomposed the eddy forcing into the Reynolds stress term ( $Rs$ ) and form stress term ( $Fs$ ):

(3.6) $$\begin{eqnarray}\displaystyle & \displaystyle Rs_{i}=-\unicode[STIX]{x1D735}\boldsymbol{\cdot }\overline{(\boldsymbol{u}_{i}^{\prime }\unicode[STIX]{x1D701}_{i}^{\prime })}, & \displaystyle\end{eqnarray}$$

(3.7) $$\begin{eqnarray}\displaystyle & \displaystyle Fs_{i}=-\unicode[STIX]{x1D716}_{i}S_{i}\unicode[STIX]{x1D735}\boldsymbol{\cdot }\overline{(\boldsymbol{u}_{i}^{\prime }(\unicode[STIX]{x1D713}_{2}^{\prime }-\unicode[STIX]{x1D713}_{1}^{\prime }))}, & \displaystyle\end{eqnarray}$$

where $\unicode[STIX]{x1D716}_{1}=-\unicode[STIX]{x1D716}_{2}=1$ , and $i=1$ , 2 for the top and bottom layers, respectively. Note that the Reynolds stress term and form stress term are simply the convergence of eddy fluxes of relative vorticity and buoyancy, respectively. As shown in figure 5(e,f), the form stress term tends to smooth the PV gradients (correlation coefficients between the form stress term and mean PV profiles are ${\sim}-0.9$ in both layers); hence, the form stress term opposes the jets. On the other hand, the Reynolds stress term shows a positive (no) correlation with the mean PV in the top (bottom) layer, which means that only the upper-layer Reynolds stress term forces the jets – see also the discussion in Thompson & Young (Reference Thompson and Young2007). The former observation is comforting in terms of theoretical aspects, because the form stress term is expected to be against the mean PV in a baroclinically unstable flow system forced with an eastward background flow (Berloff et al. Reference Berloff, Kamenkovich and Pedlosky2009a ). The Reynolds stress term and form stress term are comparable in magnitude in the top layer; however, in the bottom layer, the magnitude of the Reynolds stress term is quite small in comparison to the form stress term, and the effect of the form stress term is more prominent. In the absence of topography, the Reynolds stress term, which is responsible for maintaining the jets in the eastward background shears, dominates over the form stress term and, overall, the eddy forcing is positively correlated with the mean PV. However, here, in the presence of zonal bottom slope, the form stress term overcomes the Reynolds stress term and the overall impact of the eddies is the opposite. This is also in agreement with the linear predictions, which say that any zonal variation in bottom topography enhances the baroclinic instability (Hart Reference Hart1975; Chen & Kamenkovich Reference Chen and Kamenkovich2013). In a study of the North Atlantic jets, Kamenkovich, Berloff & Pedlosky (Reference Kamenkovich, Berloff and Pedlosky2009) also found both kinds of eddy feedbacks (see also Melnichenko et al. Reference Melnichenko, Maximenko, Schneider and Sasaki2010; Barthel et al. Reference Barthel, Hogg, Waterman and Keating2017; Youngs et al. Reference Youngs, Thompson, Lazar and Richards2017). It is important to note that, if the eddy field is computed by high-pass time filtering in a stationary frame of reference, then the resulting eddy field would also include the drifting jets. Thus, the eddy forcing would include terms interpreted as linear stress terms in this study. This may complicate comparison between our work and some earlier studies.

Figure 6.

Cross-jet profiles for the medium-slope simulation averaged in the moving frame in the interval 10 000–20 000 days. (a,b) Eddy PV flux $\overline{v_{i}^{\prime }Q_{i}^{\prime }}$ , eddy relative vorticity flux $\overline{v_{i}^{\prime }\unicode[STIX]{x1D701}_{i}^{\prime }}$ and eddy buoyancy flux $\unicode[STIX]{x1D716}_{i}S_{i}\,\overline{v_{i}^{\prime }(\unicode[STIX]{x1D713}_{2}^{\prime }-\unicode[STIX]{x1D713}_{1}^{\prime })}$ across the jets (smoothed with 10-point moving average). (c) Heat diffusivity. The thin green curves in (a,b) represent the mean flow (normalised) in individual layers, and actual magnitudes of the velocities are shown in (c).

We also analysed the cross-jet profiles of the eddy PV fluxes $\overline{v_{i}^{\prime }Q_{i}^{\prime }}$ , relative vorticity fluxes $\overline{v_{i}^{\prime }\unicode[STIX]{x1D701}_{i}^{\prime }}$ and buoyancy fluxes $\unicode[STIX]{x1D716}_{i}S_{i}\,\overline{v_{i}^{\prime }(\unicode[STIX]{x1D713}_{2}^{\prime }-\unicode[STIX]{x1D713}_{1}^{\prime })}$ across the jets (figure 6 a,b). The eddy PV fluxes have most of the contribution from the eddy buoyancy fluxes (also referred to as ‘eddy heat fluxes’), which are uniformly negative (positive) in the top (bottom) layer. In order to understand the role of the eddy buoyancy fluxes, we computed the eddy heat diffusivity  $\unicode[STIX]{x1D707}_{h}$ ,

(3.8) $$\begin{eqnarray}\unicode[STIX]{x1D707}_{h}=\frac{\unicode[STIX]{x1D716}_{i}S_{i}\,\overline{v_{i}^{\prime }(\unicode[STIX]{x1D713}_{2}^{\prime }-\unicode[STIX]{x1D713}_{1}^{\prime })}}{-\unicode[STIX]{x2202}_{q}\overline{b}_{i}}=\frac{\unicode[STIX]{x1D716}_{i}S_{i}\,\overline{v_{i}^{\prime }(\unicode[STIX]{x1D713}_{2}^{\prime }-\unicode[STIX]{x1D713}_{1}^{\prime })}}{-\unicode[STIX]{x2202}_{q}(\unicode[STIX]{x1D716}_{i}S_{i}(\overline{\unicode[STIX]{x1D713}}_{2}-\overline{\unicode[STIX]{x1D713}}_{1})+\unicode[STIX]{x1D716}_{i}S_{i}U_{b}q\cos \unicode[STIX]{x1D703})}=\frac{\overline{v_{i}^{\prime }(\unicode[STIX]{x1D713}_{2}^{\prime }-\unicode[STIX]{x1D713}_{1}^{\prime })}}{\overline{u}_{2}-\overline{u}_{1}-U_{b}\cos \unicode[STIX]{x1D703}},\end{eqnarray}$$

where $\overline{b}_{i}$ represents the time-mean buoyancy in the layer. Note that the notation is the same as in (3.6) and (3.7). Although the heat diffusivity relations are different in the top and bottom layers, they are mathematically equivalent and would result in the same cross-jet profile of the heat diffusivity in both layers. As shown in figure 6(c), the heat diffusivity is positive everywhere. This means that the eddies advect buoyancy anomaly in the direction of decreasing $\overline{b}$ , resulting in the smoothing of mean PV gradients, or, in other words, the eddy buoyancy fluxes are down-gradient. The minima in the eddy heat diffusivity profile coincide with the westward jets, whereas the eddy heat diffusivity is strongest just south of the eastward jets. Since the largest eddy flux covers the eastward jet core, the dynamics is profoundly different from the barotropic ‘PV staircases’ characterised by partial eddy transport barriers located on the eastward jet cores (e.g. Dritschel & McIntyre Reference Dritschel and McIntyre2008). Also, the magnitude of the eddy heat diffusivity increases with zonal slope. On the other hand, the eddy relative vorticity fluxes force the jets (figure 6 a,b), as they are strongly correlated with the mean flow (see expression (15.21) in Vallis Reference Vallis2017). Since eddy heat diffusivity and PV diffusivity are not the same, eddy PV diffusivity was also computed in each layer and the PV diffusivity magnitudes were found to be largely positive. Overall, the eddy buoyancy fluxes dominate over the eddy relative vorticity fluxes; hence, the eddy PV fluxes are down-gradient. The roles of these fluxes are the same, as in the case of stationary zonal jets (Panetta Reference Panetta1993; Thompson & Young Reference Thompson and Young2007), and only the relative magnitudes of the relative vorticity and buoyancy fluxes change because of the zonal asymmetry.

Given this, there must exist a critical value of zonal slope where the effects of the Reynolds stress term and form stress term exactly counter each other, and eddies start to oppose the jets beyond that critical slope. We computed correlation coefficients between the eddy forcing and mean PV profiles in simulations with different zonal slope magnitudes (table 3). Over a flat bottom, jets are supported by the eddy forcing and this is clearly captured by a strong positive correlation in the top layer for $T_{x}=0$ . As the zonal slope increases, the eddy forcing correlation coefficient in the top layer changes sign from positive to negative. On the other hand, correlation coefficients in the bottom layer are negative in all cases. This is mostly because only the upper-layer Reynolds stress term is responsible for maintaining the jets (Thompson & Young Reference Thompson and Young2007), whereas the lower-layer Reynolds stress term is small in magnitude and does not have any significant effect (see figure 5 f). It should be noted that the critical slope value computed here is not unique, and it would differ in different parameter regimes. Moreover, periodic boundary conditions may also affect the value of critical slope.

Looking at both the eddy forcing and eddy PV fluxes, we conclude that the vertical structure plays a significant part in deciding the role of eddies. Oceanic eddies are complex, as they can either force or damp a large-scale circulation. The primary role of the eddies is to stabilise the system and bound the alternating jets. In the case of a flat bottom, eddies force the jets and compensate for the energy removed through dissipation. The bounding mechanism is due to over-localisation of the eddies and the resulting loss of the eddy forcing efficiency (Berloff et al. Reference Berloff, Karabasov, Farrar and Kamenkovich2011). In contrast, in the presence of zonal slope, the jets are no longer parallel to the imposed background flow and can be directly forced by the vertical shear, as seen in the inertial stress term $L^{m}$ . Here, the eddies act against the jets to balance the linear forcing, so that the jets can be stabilised and bounded.

Table 3.

Correlation coefficients between the cross-jet profiles of the time-mean nonlinear stress terms, and the time-mean PV ( $\overline{Q}_{i}$ ) in both layers in different simulations (all parameters were kept the same, as given in § 2.1, except  $T_{x}$ ). As the zonal slope is increased, eddies change their behaviour from being jet-supporting to jet-opposing.

3.3 Jet drift mechanism

The analysis in the previous sections predicts that the linear dynamics controls the jet drift. In order to understand this physically, the linear stress terms in expressions (3.4) and (3.5) need to be examined more carefully. Here, the inertial stress term, which is the forcing term for the mean PV, is proportional to the developed mean flow. By construction, the maxima of the mean along-jet velocity and mean PV profiles do not align and are significantly offset. Similarly, the radiation stress term, which partially compensates against the inertial stress term, is also significantly shifted and, in fact, exactly orthogonal to the mean PV profile (i.e. zero correlation), as it is proportional to the first derivative of $\overline{Q}_{i}$ (see table 4 for the maximum/minimum correlation coefficients and corresponding lead/lag phases). In contrast, the eddy forcing maxima are almost aligned with the mean PV maxima, and there is a negligible improvement in correlations with the time-lagged profiles (table 5). The maximum increase in local PV occurs at the locations of the along-jet velocity maxima (figure 7). To account for PV conservation along with this local PV increase, the jets drift southwards. In other words, the off-core linear forcing pushes the jets in the meridional direction. The same is true for northward drifting jets in the case of a negative zonal slope, where the relative shifts between the mean PV and mean flow are simply reversed.

Table 4.

Correlation coefficients and maximum/minimum lag-correlation coefficients between the time-mean cross-jet profiles of the linear stress terms and the time-mean PV ( $\overline{Q}_{i}$ ) in both layers for the medium-slope simulation. Here, the lag period of $2\unicode[STIX]{x03C0}$ is equivalent to the meridional width of one pair of jets.

Figure 7.

Sketch of cross-jet profiles of the time-mean PV (blue thick) and mean flow (red-brown thin) in both layers in numerical simulations with positive zonal slopes in topography (relative magnitudes are not to scale). The maxima of the mean flow and mean PV profiles do not align, and there is a large shift. Red dashed arrows indicate the tendency of PV to increase/decrease at the velocity maxima due to the mean flow, and blue solid arrows indicate the direction of the drift. Vertical dotted lines represent the zero lines.

Table 5.

Correlation coefficients and maximum/minimum lag-correlation coefficients between the time-mean cross-jet profiles of Reynolds stress term and form stress term, and the time-mean PV ( $\overline{Q}_{i}$ ) in both layers for the medium-slope simulation. Here, the lag period of $2\unicode[STIX]{x03C0}$ is equivalent to the meridional width of one pair of jets.

The drift mechanism also agrees with the hypothesis by Boland et al. (Reference Boland, Thompson, Shuckburgh and Haynes2012). They suggest that the jet drift is a result of PV conservation, and the along-jet PV advection is compensated by a meridional drift. From (3.4) and (3.5), it is clear that the inertial stress term is basically advection of a PV component by the mean flow, and the jets drift to compensate for this advection (see expression of the radiation stress term). Melnichenko et al. (Reference Melnichenko, Maximenko, Schneider and Sasaki2010) studied the role of eddy fluxes in the formation of striations in the Pacific. They also showed that linear advection terms appear in the time-averaged PV budget due to the non-zero angle between the background flow and tilted striations. Here, too, the striations tend to drift across the mean PV contours, which may be due to the presence of the linear advection terms. It is natural to ask why eddies do not stop this drift, as the eddy forcing maxima are also not aligned with the mean PV maxima in a way that would support a northward drift and counteract the southward propagation of the jets. Indeed, the eddies do tend to slow the jets down (see table 1); however, for the topographic slopes considered, the eddy effect is small and cannot overcome the drift induced by the mean flow. On the other hand, the eddy effect may be significant for steeper slopes, which was seen in Thompson & Richards (Reference Thompson and Richards2011) and Chen et al. (Reference Chen, Kamenkovich and Berloff2015).

It is important to note here that a sloping topography is not the only way to generate drifting jets, which are also observed in systems forced with a meridionally non-uniform stress term (Williams Reference Williams2003; Chan, Plumb & Cerovecki Reference Chan, Plumb and Cerovecki2007; Chemke & Kaspi Reference Chemke and Kaspi2015b ). In these systems, the drift cannot be explained by the off-core linear forcing mechanism described above, which is only valid for tilted, drifting jets over topography, and is thought to be caused by meridionally asymmetric eddy momentum fluxes or a residual flow in the meridional direction. Srinivasan (Reference Srinivasan2013) suggests that the jet drift in general is caused by the breaking of reflectional (or mirror) symmetry in a system. Most of the studies that use idealised models on the $\unicode[STIX]{x1D6FD}$ -plane are forced with either a uniform zonal flow or an isotropic small-scale random forcing. These choices preserve the mirror symmetry and, thus, generate steady zonal jets. The mirror symmetry can be broken, even in a barotropic model, e.g. by using an anisotropic forcing term, and in response drifting jets can be formed (see chap. 3 in Srinivasan Reference Srinivasan2013). In our study, the mirror symmetry is broken by the zonal slope in topography, resulting in drifting jets. However, there are differences between these two scenarios. The drift mechanisms are different in these two cases, and also, the drifting jets with anisotropic forcing are purely zonal, whereas the jets formed over zonally sloping topography are tilted. The concept of the mirror symmetry can potentially predict the drifting behaviour in various systems, but it is not clear, whether it can be also quantitative.

3.4 Interpreting the linear forcing

In the presence of a zonally sloping topography, the jets are maintained by the inertial stress term, which is proportional to the mean flow, and the eddies act against the jets. It is very tempting to assume here that the forcing is due to the mean flow itself; however, that is not feasible according to the laws of physics. The jets cannot force themselves and simultaneously lose energy to the eddies. Here, the jets are directly forced by the imposed vertical shear in the presence of the zonal slope, and the intensity of the linear forcing depends on the mean flow itself. The zonal bottom slope acts like a coupling between the tilted jets and background flow, and the jets gain energy directly from the background flow. This balance between the jets and eddies only works in forced dissipative turbulent systems. If the background flow and dissipation terms were removed, the jets would dissipate and lose all their energy to eddies, thus leading to a very different flow dynamics. However, a detailed analysis of that is beyond the scope of this paper and is being pursued as a separate piece of work.

3.5 Energy equilibration and the long-time flow adjustment

In addition to tilting and drifting of the jets, another difference made by the zonal component of PV gradient is that there is a significant increase in the equilibrated energy of the system (figure 2 e). Arbic & Flierl (Reference Arbic and Flierl2004) showed that the eddy energy of the equilibrium state greatly increases when the vertical shear is oriented not in the east–west direction. The increase is proportional to the angle between the direction of the vertical shear and PV isolines (see also Smith Reference Smith2007). In this work, we tilt the PV isolines by introducing the zonal slope in bottom topography, which basically introduces a non-zero angle between the background PV and vertical shear, and thus increases the energy of the equilibrated state. This can also be understood by analysing the eddy field. In the presence of the zonal slope, eddies grow through baroclinic instability (by extracting energy from the imposed background shear) and also gain energy from the mean flow itself, i.e. the jets in this case. In this situation, both the jets and eddies tend to be more energetic.

Figure 8.

Flow properties in the transient state of the medium-slope solution. The time evolutions of: (a)  $\text{barotropic velocity}=(1/(H_{1}+H_{2}))(H_{1}u_{1}+H_{2}u_{2})$ ; (b)  $\text{baroclinic velocity}=u_{1}-u_{2}$ (cross-jet profiles, averaged along the jets, are aligned in the moving frame of reference by assuming a constant drift velocity; note that legend units are in days); and (c) correlation coefficients between the cross-jet profiles of eddy forcing and PV, and the fraction of the total energy in the barotropic component. Here, we used the full flow field to compute the eddy forcing, because it was not possible to compute a mean flow in the transient state.

Moreover, it takes a long time for the system to attain equilibrium. The continuous increase of the total energy (roughly in the interval 1000–4000 days, figure 2 e) is definitely not during the time period when the jets are being formed (§ 2.1), but it has to do with some secondary effects. To investigate this, we first look at the time evolution of the barotropic and baroclinic velocity components in the transient period (figure 8 a,b). Both the barotropic and baroclinic velocity fields increase at similar rates, and, in fact, the ratio of energy in the barotropic component to the total energy (barotropic plus baroclinic,  $E_{BT}+E_{BC}$ ) remains almost constant during this period (blue curve in figure 8 c).

One possible explanation behind this slow energy increase could be the direct forcing by the background flow once the jets are formed. The inertial stress term depends on the mean flow; hence, the linear forcing becomes prominent only after the formation of the jets. Note that the jet formation occurs due to the upscale transfer of energy by mesoscale eddies that receive energy from the background flow. However, at a later stage, the eddies start acting against the linear stress terms and oppose the jets. So, there must be a time interval when the eddies change their behaviour from being jet-supporting to jet-opposing. To confirm this hypothesis, we computed correlation coefficients between the cross-jet eddy forcing and PV profiles during the transient period (figure 8 c). The time evolution of the correlation coefficient shows a decreasing trend in the top layer, where the correlation coefficient changes from a small positive value to negative values with time. On the other hand, the correlation remains negative in the bottom layer. This is probably because the eddy forcing in the bottom layer is largely due to the form stress term, thus producing a negative correlation at all times (also read the discussion of table 3). The results are in accord with our hypothesis, at least in the more active top layer that offers an explanation of the long time scale in the transient period.