3.1. Regime I: low-amplitude, locally forced jets

At low, stationary forcing amplitude, we observe the fast development of five prograde jets – in the same direction as the tank's rotation, $u_\phi >0$ (figure 1) – and six retrograde jets. We will refer to this regime as regime I.

To describe this regime, we chose a typical experiment where the pumps power are respectively $P_i = \{ 7, 10, 20, 30, 45, 90\}$ % of their nominal power, corresponding to a forcing amplitude $U_f= 2.4 \times 10^{-3}\ \mathrm {m}\ \textrm {s}^{-1}$ (see appendix B). Figure 4(a) represents the temporal evolution (Hovmöller diagram) of the instantaneous azimuthal mean of the azimuthal component of the velocity $\langle u_\phi \rangle _\phi (\rho ,t)$ – called zonal flow in the following of the paper, whereas mean flow refers to the time-averaged velocity field. The jets develop almost instantaneously (${\sim }10\,t_R$), and reach their saturating amplitude in about 100 $t_R$ (another example is shown in figure 6 for $t < 360\,t_R$). Supplementary movie 2 (azimuthal component of the velocity) and movie 3 (vorticity) illustrate the development of regime I.

Figure 4. Zonal flow and zonal mean potential vorticity in the two experimental regimes represented in figures 3 and 5. (a–c) Regime I. (a) Space–time diagram showing the evolution of the instantaneous zonal flow radial profile $\langle u_\phi \rangle _\phi (\rho ,t)$ during the experiment. Dotted lines: location of the forcing rings. (b) Symbols: instantaneous zonal flow at time $t=2500\ t_R$. Black line: time-averaged zonal flow $\overline {\langle u_\phi \rangle _\phi }(\rho )$ from $t=1500$ to 3000 $t_R$. Dotted lines: location of the forcing rings. (c) Symbols: potential vorticity instantaneous profile at time $t=2500\ t_R$. Grey line: initial potential vorticity profile (hidden behind the symbols). (d–f) Regime II (the same quantities are plotted).

The velocity and vorticity fields obtained after saturation are represented in figures 3(a) and 5(a–c). The retrograde jets are uniform and quasi-axisymmetric, whereas the progade jets are associated with clear non-axisymmetric perturbations. Consistently with the direction of the zonal flow, the anticyclones – negative vorticity $\zeta$ – are located on the outer radius flank of the prograde jets, whereas cyclones are located on their inner radius side. In addition, the prograde jets are thinner than the retrograde ones. These observations highlight the asymmetry between prograde and retrograde jets, generically observed in these type of systems (e.g. Scott Reference Scott2010; Scott & Dritschel Reference Scott and Dritschel2012).

Figure 5. Instantaneous maps in regime I (a–c) and II (d–f ). The black circle marks the tank boundary, and the shaded areas are the shadows where no measurements can be performed. The dashed curves in the top-left quadrant of each subplot represent the forcing rings location. (a–c) Regime I at time $t=24 \ \min = 1800$ $t_R$ and averaged over 1 s. (a) Azimuthal component of the velocity $u_\phi$. (b) Radial component of the velocity $u_\rho$. (c) Vertical component of the vorticity $\zeta$. (d–f ) Regime II at time $t=19\ \min = 1425$ $t_R$ and averaged over 1 s. (d) $u_\phi$. (e) $u_\rho$. (f) $\zeta$. Note the different colour scales for the two regimes.

The saturated zonal flow profile is plotted in figure 4(b) along with its time average, and figure 4(c) shows the zonal mean of the potential vorticity $\langle (\zeta + f)/h \rangle _\phi$. In the absence of dissipation, we expect the material conservation of potential vorticity (PV) (Vallis Reference Vallis2006, § 4.5). In this limit, zonal flows formation can be viewed as a process of mixing of the initial potential vorticity profile $f/h(\rho )$. Dritschel & McIntyre (Reference Dritschel and McIntyre2008) showed that this profile should be turned into a staircase where the prograde jets correspond to steep gradients, and the retrograde jets correspond to weak gradients, i.e. zones of strong mixing. Here, despite the visible segregation of vorticity (figure 5c), the initial vorticity profile is almost not perturbed, showing that zonal jets can exist instantaneously even without this process of potential vorticity mixing. Said differently, this regime is characterized by a local Rossby number $Ro_\zeta = \zeta /f \ll 1$ (see table 1); hence, the initial PV profile is not expected to be strongly modified. The instantaneous root-mean-square velocity, defined as

(3.1)\begin{equation} {u}^\textit{RMS} = \left[ \frac{1}{N} \sum_{i=1}^{N} | {\boldsymbol{u}_i} |^{2} \right]^{1/2}, \end{equation}

where $N$ is the number of PIV velocity vectors, is provided in table 1 along with the global and local Rossby and Reynolds numbers of the flow. Here, $u^{\textit{RMS}}$ is computed from the velocity fields of figure 3. It is considered ‘instantaneous’ in opposition to the same quantity computed after a very long time average $\bar {u}^{\textit{RMS}}$ which will be used later in the paper. Table 1 shows that the flow is barely turbulent in regime I (the local Reynolds number being approximately 100), and highly constrained by rotation given the very small Rossby numbers. Finally, the zonal flow contains $23 \pm 5\,\%$ of the total kinetic energy.

Table 1. Typical instantaneous RMS velocity (3.1), Rossby and Reynolds global and local numbers. Here $R$ is the tank's inner radius $R=0.49$ m, and $f = 2\varOmega = 15.7\ \textrm {rad}\ \textrm {s}^{-1}$. For the local Rossby and Reynolds numbers, we use the distance between two forcing rings as a length scale, $\ell =7.3$ cm. Note that these values correspond to two typical experiments, but may vary in each regime depending on the forcing amplitude.

In this regime, each prograde jet stands right above a forcing ring (see the dashed lines in figures 4 and 5). The only exception is the inner ring (ring C1), which is geometrically constrained due to its small radius and significantly perturbed by the peak at the centre of the bottom plate (figure 2c). Despite this anomaly, the five other forcing rings are clearly associated with a prograde jet. This leads us to hypothesize that the prograde jets are forced locally by prograde momentum convergence towards the forcing radii. The local Reynolds stresses generated by our forcing are then balanced by viscous effects. This mechanism of zonal flow formation is reminiscent of the pioneering experiments of Whitehead (Reference Whitehead1975) and De Verdiere (Reference De Verdiere1979). Whitehead (Reference Whitehead1975) demonstrated that the generation of a train of Rossby waves in a rotating tank with a paraboloidal free surface induces a prograde flow at the radius of the forcing, with two weak retrograde flows on both sides of the forcing. De Verdiere (Reference De Verdiere1979) did the same observation with a forcing consisting in a ring of sink and sources able to be azimuthally translated. Corresponding theoretical studies of Thompson (Reference Thompson1980) and McEwan, Thompson & Plumb (Reference McEwan, Thompson and Plumb1980) then followed and accounted for the mechanism of momentum convergence due to eddy forcing. It is now believed to be the primary mechanism of westerlies formation in the mid-latitude atmosphere (Vallis Reference Vallis2006, chapter 12). This mechanism will be further explored in § 4.

Finally, let us mention that the relaxation dynamics of this regime is consistent with the observation of De Verdiere (Reference De Verdiere1979): when the forcing is stopped, the fluctuating velocities are dissipated more rapidly than the mean flow.

3.2. Regime II: high-amplitude large jets

At high, stationary forcing amplitude, regime I develops as a transient before the system reaches a statistically stationary state with stronger and broader zonal jets, hereafter called regime II.

To describe this regime, we chose a typical experiment where the pumps power are respectively $P_i = \{ 26, 33, 60, 80, 100, 100\}$ % of their nominal power, corresponding to a forcing amplitude $U_f= 4.0 \times 10^{-3}\ \mathrm {m}\ \textrm {s}^{-1}$ (see appendix B). The Hovmöller diagram of this experiment is represented in figure 4(d). The steady jets of the first regime reorganize into three prograde and three retrograde jets. Note that in other experiments, the saturated flow can count four prograde jets instead of three, as will be discussed in § 3.3. This spontaneous transition from regime I to regime II is slow, and the statistically steady state is obtained after a transient of about 800 $t_R$. Furthermore, it involves zonal flows merging events visible in the Hovmöller plots of figures 4(d) and 6. The reorganization of the jets during this transition also shows that they become more independent of the forcing pattern: in the final steady state of regime II, the jets have a typical width which is twice that of the jets in regime I, and their radial position can be shifted compared to the position of the forcing rings. A retrograde flow is even observed above some forcing rings, for instance, above C1 and C3. Thus, in this regime, the system self-organizes at a global scale, and the idea of a direct local forcing is not relevant anymore. Supplementary movie 4 illustrates the development of regime II, and movie 1 shows the particles motion when the system is in steady state.

Figure 6. Experiment illustrating the bistability between regimes I and II: the forcing is initially of low amplitude and stationary, so that we start in regime I ($P_i=\{14,20,46,72,100,100\}\,\%$, $U_f=3.3 \times 10^{-3}\ \mathrm {m}\ \mathrm {s}^{-1}$). At $t=360\ t_R$, a finite perturbation is created by varying the third forcing ring power randomly around its initial value, leading to the transition to regime II ($P_i=\{14,20,46\pm 20,72,100,100\}\,\%$). At $t=1600\ t_R$, the forcing is set back to its initial state but the flow remains in regime II. At time $t=2250\ t_R$, a second finite-amplitude perturbation is performed ($P_i=\{14,20,46,72\pm 20,100,100\}\,\%$) (a). Hovmöller plot: zonal flow profile as a function of time. (b) Total ($\mathcal {K}$), zonal ($\mathcal {K}_z$) and fluctuating ($\mathcal {K}_f$) kinetic energy as a function of time.

The velocity field obtained after saturation is represented in figure 3(b), and the corresponding maps of velocity and vorticity are plotted in figure 5(d–f). The prograde jets are still meandering between cyclones on their right and anticyclones on their left, but these vortices are now large-scale ones. As can be seen in figure 3(a), the vortices forced above the inlets and outlets have a typical diameter of $\sim$3 cm in regime I, whereas in regime II (panel (b)), we observe fewer vortices, with a typical diameter of $\sim$8 cm. The instantaneous RMS velocity (table 1) is about 10 times higher than in the experiment described for regime I. The global and local Rossby numbers are still very small, i.e. the flow is still highly constrained by rotation, but the Reynolds number is multiplied by 10, hence, the flow can now be considered fully turbulent. The fraction of kinetic energy contained in the zonal flow in this experiment reaches $58 \pm 8\,\%$. Figure 4(f) shows that the PV mixing is increased in this second regime and consistently with Dritschel & McIntyre (Reference Dritschel and McIntyre2008), the prograde jets correspond to steepening of the PV profile. But again, the small vorticity of our experiment does not allow an efficient mixing process, though the jets are strong and contains most of the kinetic energy.

3.3. Nature of the transition: a first-order subcritical bifurcation

In this section we investigate the nature of the transition between the two previously described experimental regimes.

Figure 6 shows a Hovmöller diagram representing the evolution of the zonal flow profile $\langle u_\phi \rangle _\phi (\rho ,t)$ during a single experiment as well as the corresponding evolution of the total, zonal and fluctuating kinetic energy defined respectively as

(3.2)\begin{gather} \mathcal{K} = \frac{1}{N} \sum_{i=1}^{N} | \boldsymbol{u}_i | ^{2}, \end{gather}

(3.3)\begin{gather}\mathcal{K}_z = \frac{1}{N} \sum_{i=1}^{N} \langle u_{\phi} \rangle_{\phi,i}^{2}, \end{gather}

(3.4)\begin{gather}\mathcal{K}_f = \mathcal{K}-\mathcal{K}_z, \end{gather}

where $N$ is the number of PIV velocity vectors. The experiment plotted in figure 6 is initialized with a stationary forcing ($U_f=3.3 \times 10^{-3}\ \mathrm {m}\ \mathrm {s}^{-1}$), leading to a steady state in regime I. After 360 $t_R$, this forcing is perturbed at finite amplitude by turning the third ring into a random state. Here, it consists in changing its power every 3 s to random values uniformly distributed in a range centred around ${\pm } 20\,\%$ of its initial power. After such a perturbation, figure 6(a) shows that the system bifurcates towards the second regime through merging events and increasing zonal flow amplitude. Note that without this perturbation, the system would be locked in regime I, as shown by a separate experiment performed with the exact same forcing, at least up to $t=1875\ t_R$. During the transition, the fraction of kinetic energy contained in the zonal flow increases from 21$\pm$7 % to 48$\pm$9 % (figure 6b). This second value is significantly lower than the one mentioned previously for regime II since the forcing of this experiment is weaker. After the transition, the system remains attached to this new steady state even when the forcing is set back to its initial stationary state at time $t=1600\ t_R$. These observations demonstrate that two stable states coexist for this particular forcing and suggest that the transition between the two regimes is of subcritical nature.

To further investigate this bistability, we perform series of experiments where we increase or decrease the forcing step by step. We wait significantly between each step for the system to relax towards a new steady state (typically 20 min, i.e. $1500\ t_R$). We then measure the corresponding mean flow amplitude defined as the RMS velocity computed on a time-averaged velocity field,

(3.5)\begin{equation} \bar{u}^\textit{RMS} = \left[ \frac{1}{N} \sum_{i=1}^{N} | \overline{\boldsymbol{u}_i} |^{2} \right]^{1/2}, \end{equation}

where $\overline {\cdot }$ denotes the time average over the whole duration of the record once the flow has reached the statistically steady state. Typically, the time average is performed over 200 to 1000 $t_R$ depending on the duration of the record. Figure 7(a) represents the mean flow amplitude (3.5) as a function of the forcing amplitude $U_f$ as defined in appendix B. Typical maps of the time-averaged zonal velocity are represented in panel (b). For low values of the forcing amplitude, regime I is observed with the six prograde jets structure and $\bar {u}^{\textit{RMS}} \sim 2.5 \times 10^{-3} \ \mathrm {m}\ \mathrm {s}^{-1}$. As the forcing amplitude increases, the jets structure does not change but their amplitude increases smoothly. When the forcing further increases, a sharp transition occurs around $U_f \approx 3.32\times 10^{-3}\ \mathrm {m}\ \mathrm {s}^{-1}$ corresponding to a bifurcation from regime I to regime II: both the jets size and amplitude increase abruptly ($\bar {u}^{\textit{RMS}} \sim 7.5\times 10^{-3}\ \mathrm {m}\ \mathrm {s}^{-1}$). Once in regime II, the amplitude of the jets continues to increase with the forcing amplitude. When the forcing amplitude is gradually decreased, the bifurcation from regime II to regime I is again abrupt, but obtained at a lower forcing $U_f \approx 3.11\times 10^{-3}\ \mathrm {m}\ \mathrm {s}^{-1}$. These hysteresis experiments confirm that the two regimes coexist in a given forcing range $U_f \in [3.11,3.32]\times 10^{-3}\ \mathrm {m}\ \mathrm {s}^{-1}$. The particular forcing of the experiment represented in figure 6 ($U_f=3.3\times 10^{-3} \ \mathrm {m}\ \mathrm {s}^{-1}$) belongs to the bistable range in which the first regime is metastable. In § 4 we propose a model to explain this hysteresis phenomenon.

Figure 7. (a) Experimental hysteresis loop. The time-averaged RMS velocity of the flow $\bar {u}^{\textit{RMS}}$ (3.5) is plotted as a function of the forcing amplitude $U_f$ (see appendix B). The top and right axes correspond to associated non-dimensional quantities. On the right axis, we use the Reynolds number based on the mean flow RMS velocity, $\overline {Re}=\bar {u}^{\textit{RMS}}\ell /\nu$ ($\ell =7.3\ \textrm {cm}$ is the distance between two forcing rings). For the top axis, we use the typical velocity expected at the transition $U_t$ (see § 4.4 and (4.34)). For each curve, the forcing is either increased (reddish) or decreased (bluish) step by step. The different colours are different experiments. The shaded area is the bistable zone. The grey points in regime II correspond to experiments initialized with the exact same forcing and for which saturation leads to three possible jets configurations represented in panel (b). The last point of the yellow curve is in configuration 2, but may not have reached its stationary state. (b) Time-averaged zonal velocity maps in regime I and II. In regime II multiple steady states can be obtained for the same forcing. The configurations 1, 2, 3 correspond to the points labelled accordingly in panel (a).

Finally, we note a significant variability in the mean flow amplitude in regime II. The grey points in figure 7(a), located at $U_f=4\times 10^{-3}\ \mathrm {m}\ \mathrm {s}^{-1}$, correspond to nine experiments where we apply the exact same forcing ($P_i=\{26,33,60,80,100,100\}\%$), starting from solid body rotation. Despite the similarity of the forcing, the flow may evolve towards different statistically steady states where the mean flow amplitude and scale are roughly the same, but the position of the jets differ. Based on 15 realizations, we have identified three different steady states with permutations between the location of the prograde jets, as represented in figure 7(b). The last point of the yellow curve in figure 7(a) is in configuration 2. It is located below the others probably because it had not reached its steady state when the measurements were performed (500 $t_R$ after the forcing change in contrast to $1500\ t_R$ for the grey points). Note that we have never observed spontaneous transitions between these three states, even during day-long experiments (38 000 $t_R$). The origin of this multi-stability is beyond the scope of the present study, but will be the focus of future investigations. In particular, its link with the multi-stability observed recently in numerical simulations of stochastically forced barotropic turbulence (Bouchet, Nardine & Tangarife Reference Bouchet, Nardine and Tangarife2019a; Bouchet, Rolland & Simonnet Reference Bouchet, Rolland and Simonnet2019b) should be addressed.

It is of interest to compute the transition rates between the two regimes. On figure 8 we plot the evolution of the total kinetic energy for transitions from regime I to regime II and vice versa in order to compute the corresponding time scales. Transitions in the direction $\textrm {II}\rightarrow \textrm {I}$ (decreasing power) are accompanied by an exponential decay of the total kinetic energy

(3.6)\begin{equation} \mathcal{K}=\mathcal{K}_0 + (\mathcal{K}_\infty-\mathcal{K}_0) \,\textrm{e}^{-t/\tau}, \end{equation}

where $\mathcal {K}_0$ is the kinetic energy in the initial steady state, and $\mathcal {K}_\infty$ the kinetic energy reached in the final steady state after the transition. We plot in figure 8(b) the time evolution of the normalized kinetic energy for three transitions with different initial and/or final steady states. Despite these differences, it is clear that the three transitions have the same characteristic time $\tau _{\textrm {II}\rightarrow \textrm {I}} \approx 150\ t_R$. The picture is different for transitions in the direction $\textrm {I}\rightarrow \textrm {II}$ (figure 8c). The kinetic energy increases in a non-trivial way before saturating in a new steady state. We compare this evolution for three transitions starting from the same initial state in the first regime (red dot I in figure 8(a)), but evolving towards three different steady states in regime II ($\textrm {II}_{1,2,3}$ in figure 8(a)). Note that $\textrm {I}\rightarrow \textrm {II}_1$ corresponds to a finite-amplitude perturbation of a steady state in regime I inside of the bistable range. We plot in figure 8(c) the kinetic energy normalized the same way as for $\textrm {II}\rightarrow \textrm {I}$ transitions. This time, the curves do not collapse. We observe that the closer (in terms of forcing amplitude) the second steady state, the longer the transition. The transition $\textrm {I}\rightarrow \textrm {II}_3$ is, for instance, about three times faster than the transition $\textrm {I}\rightarrow \textrm {II}_1$.

Figure 8. Evolution of the total kinetic energy $\mathcal {K}$ during transitions from regime $\textrm {I} \rightarrow \textrm {II}$ (red curves) and $\textrm {II} \rightarrow \textrm {I}$ (blue curves). The normalized total kinetic energy is plotted as a function of time in units of rotation period. (a) Qualitative location of the transitions on the hysteresis loop. (b) Transitions from $\textrm {II}\rightarrow \textrm {I}$. The time is initialized at the moment when the forcing was changed from a super-critical one to a subcritical one. The normalized kinetic energy decays on a time scale $\tau$ obtained from an exponential fit (lines). (c) Transitions from $\textrm {I}\rightarrow \textrm {II}$, starting from the same initial state.

These differences highlight an asymmetry in the transition mechanism depending on its direction. We propose the following interpretation. In the case $\textrm {II}\rightarrow \textrm {I}$ the transition resembles a classical relaxation following a linear dissipation process. If the Reynolds stresses sustaining the strong jets abruptly decrease when the forcing decreases, the linear friction dominates the zonal flow evolution equation $\partial _t \langle u_\phi \rangle _\phi = -\alpha \langle u_\phi \rangle _\phi$ (see (4.1) and (A10)), where $\alpha$ is a linear friction, and the zonal flow decreases exponentially. We then expect the time scale of this transition to be of the order of the Ekman friction time scale $\tau _E=\alpha ^{-1}=\varOmega ^{-1}E^{-1/2}$. In our case, $\tau _E\approx$ 206 rotation periods, which is consistent with $\tau _{\textrm {II}\rightarrow \textrm {I}}\approx 150\ t_R$ determined previously. On the contrary, in the case $\textrm {I}\rightarrow \textrm {II}$, we expect a nonlinear mechanism leading to a non-trivial increase of the zonal flow amplitude. Contrary to the linear friction, this mechanism depends on the forcing amplitude – as may be intuited from the theoretical model developed in § 4. The higher the forcing, the faster we expect the transition to occur.

4.1. Linear model for the Reynolds stresses

In this section we determine the Reynolds stresses divergence $\mathcal {R}$ and its sensitivity to the zonal flow. To do so, we compute the linear response to a stationary forcing on the $\beta$- plane, in the presence of a background zonal flow $U$. Besides, we adopt a local approach by assuming a length scale separation between the wavelength of the forcing and the spatial variations of the zonal flow. We also assume homogeneity by considering an infinite fluid domain in both directions. This approach allows us to

1. (i) forget about the geometrical effects inherent to the cylindrical geometry and work in equivalent 2-D local cartesian coordinates $(x,y)$ (see figure 1);

2. (ii) assume that the background flow $U$ is constant in $(x,y)$, which is only true locally, inside of a single jet.

For the basis $(\boldsymbol {e}_x,\boldsymbol {e}_y,\boldsymbol {e}_z)$ to be direct, with $\boldsymbol {e}_z$ downward and $\boldsymbol {e}_x$ zonal, in the same direction as $\boldsymbol {e}_\phi$, $\boldsymbol {e}_y$ has to be oriented towards the axis of rotation (figure 1). We denote $\boldsymbol {u}=(u,v)_{ \boldsymbol {e}_x, \boldsymbol {e}_y}$ the 2-D cartesian velocity components, and $\zeta = \partial _x v - \partial _y u$ the associated vorticity. The $\beta$-plane barotropic vorticity equation (4.1) reduces to

(4.10)\begin{equation} \frac{ \partial \zeta}{ \partial t} + u \frac{ \partial \zeta}{ \partial x} + v \frac{ \partial \zeta}{ \partial y} + \beta v + \alpha \zeta = \nu \nabla^{2} \zeta + q(x,y), \end{equation}

with $\alpha$ defined by (4.3) and

(4.11)\begin{equation} \beta = - \frac{f}{h} \frac{\mathrm{d}h}{\mathrm{d}\kern 0.05em y}. \end{equation}

Note that in this cartesian framework $\beta$ is now positive $(\mathrm {d}_yh <0)$. We have added an arbitrary stationary forcing $q(x,y)$ representing a vorticity source. We linearize this equation around a uniform background zonal flow $\boldsymbol {U}=U\,\boldsymbol {e}_x$ by setting $\boldsymbol {u}=\boldsymbol {U}+\boldsymbol {u'}$ with $| \boldsymbol {u'} | \ll | \boldsymbol {U} |$, and keeping only first-order terms:

(4.12)\begin{equation} \frac{ \partial \zeta'}{ \partial t} + U \frac{ \partial \zeta'}{ \partial x} + \beta v' + \alpha \zeta' = \nu \nabla^{2} \zeta' + q(x,y). \end{equation}

We drop the primes in the following and define the streamfunction $\psi$ ($u=- \partial _y\psi$, $v= \partial _x\psi$ and $\zeta = \nabla ^{2} \psi$) such that

(4.13)\begin{equation} \displaystyle \frac{ \partial}{ \partial t} \nabla^{2} \psi + U \frac{ \partial}{ \partial x}\nabla^{2}\psi + \beta \frac{ \partial \psi}{ \partial x} + \alpha \nabla^{2} \psi - \nu \nabla^{2} \nabla^{2} \psi = q(x,y). \end{equation}

We perform a spatial Fourier transform of this equation in both $x$ and $y$ leading to

(4.14)\begin{equation} \frac{ \partial \hat{\psi}}{ \partial t} + \left[ \mathrm{i}\omega(k,l) + \omega_E(k,l) \right] \hat{\psi} = - \frac{\hat{q}(k,l)}{k^{2}+l^{2}}, \end{equation}

where $\hat {\psi }$ and $\hat {q}$ are the Fourier coefficients associated with $\psi$ and $q$, and $\boldsymbol {k}=(k,l)_{\boldsymbol {e}_x,\boldsymbol {e}_y}$ the wave vector. We denote $\omega$ the Rossby waves angular frequency, Doppler shifted by the advection by the zonal flow

(4.15)\begin{equation} \omega = kU - \frac{k \beta}{k^{2}+l^{2}}, \end{equation}

(Vallis Reference Vallis2006, p. 230) and $\omega _E$ the damping rate due to the viscous dissipation in the bulk and the bottom friction

(4.16)\begin{equation} \omega_E = \alpha + \nu(k^{2}+l^{2}). \end{equation}

Note that there is no gravity effects (or deformation radius) in the Rossby waves dispersion relation because we make the rigid lid approximation (appendix A). This is justified in the present work since the short wave limit is valid for the forced waves. The solution to (4.14) with the initial condition $\hat {\psi }(k,l,t=0)=0$ is

(4.17)\begin{equation} \hat{\psi}(k,l,t)=\frac{-\hat{q}(k,l)}{(k^{2}+l^{2})(\mathrm{i}\omega(k,l)+\omega_E)} \left[ 1 - \exp({-(\mathrm{i}\omega+\omega_E)t}) \right]. \end{equation}

The inverse Fourier transform $\mathcal {F}^{-1}$ of $\psi$ can be computed numerically to retrieve the physical streamfunction $\psi$ at a time $t$ for a given forcing. Similarly, the vorticity and velocity components can be computed using

(4.18)\begin{equation} \left. \begin{gathered} \zeta(x,y,t) = \mathcal{F}^{-1}\left( -(k^{2}+l^{2})\hat{\psi}(k,l,t) \right), \\ u(x,y,t) = \mathcal{F}^{-1}\left( -\mathrm{i}l\hat{\psi}(k,l,t) \right), \\ v(x,y,t) = \mathcal{F}^{-1}\left( \mathrm{i}k\hat{\psi}(k,l,t) \right). \end{gathered}\right\} \end{equation}

The Reynolds stresses term $\mathcal {R}(y,t;U)$ is then easily computed as

(4.19)\begin{equation} \mathcal{R}(y,t;U) = -\frac{ \partial}{ \partial y} \langle uv \rangle_x = \langle v\zeta \rangle_x. \end{equation}

The sign difference in the second equality compared to expression (4.8) in cylindrical coordinates comes from $\boldsymbol {e}_y$ pointing inward whereas $\boldsymbol {e}_\rho$ is pointing outward.

4.2. Comparison of the linear model with experimental results

To confirm that the reduced QG approximation is an appropriate model for the experiment, we compare the linear solution with the very beginning of experiments where only one forcing ring is turned on. Note that a good agreement is expected since in our experiments, $Ro_\zeta = \zeta /f \ll 1$ (table 1), which is the main assumption of the QG approximation. To carry this comparison, we first set all the model parameters to the experimental ones, that is,

(4.20)\begin{equation} \left. \begin{gathered} \displaystyle \beta = -\frac{f}{h} \frac{\mathrm{d}h}{\mathrm{d} y} \sim 50\ \mathrm{m}^{-1}\ {\mathrm{s}^{-1}}, \\ \alpha = \frac{1}{2} f E^{1/2} \sim 5.6\times 10^{-3}\ \mathrm{s}^{-1}, \\ U = 0\ \mathrm{m}\ {\mathrm{s}}^{-1}. \end{gathered}\right\} \end{equation}

When the pump is activated, the fluid is at rest in the rotating frame (solid body rotation). In addition, we design the forcing term to mimic the experimental forcing on the chosen ring: for the third ring, it corresponds to a line of 18 vortices (9 cyclones and 9 anticyclones) regularly spaced onto a perimeter of $2{\rm \pi} R_3 \approx 0.688$ m. We chose to represent each vortex by a Gaussian one, leading to

(4.21)\begin{equation} q(x,y)= q_m \sum_{i=1}^{18} (-1)^{i} \exp\left(-\frac{(x-x_i)^{2} + (y-y_i)^{2}}{r_v^{2}}\right), \end{equation}

with $q_m$ the forcing amplitude, $r_v$ the radius of the vortices and $(x_i,y_i)$ the centre location of each vortex. The vortices’ radius is set to $2/5$ times the spacing between two vortices based on the experimental measurements. To estimate the forcing amplitude $q_m$ that we should use to better represent the experimental regime, we measure the vorticity linear growth rate above the forcing injection points and adjust $q_m$ so that the growth rate obtained with the linear model is comparable to the experimental one. This method leads us to use $q_m = 0.5\ \mathrm {s}^{-2}$.

Figure 9 shows a comparison between an experiment and the linear model 2 s after the third forcing ring was turned on at its maximum power. The experimental flow and the linear solution are qualitatively and quantitatively very close. They both exhibit a westward stretching of each vortex, westward meaning in the retrograde direction compared to the background rotation (decreasing $\phi$ in the experiment, decreasing $x$ in the cartesian model). The dispersive nature of the Rossby waves emitted by the vortex are responsible for this chevron pattern pointing eastward, as explained by Firing & Beardsley (Reference Firing and Beardsley1976), Flierl (Reference Flierl1977), Chan & Williams (Reference Chan and Williams1987) in the case of isolated vortices. In particular, the long Rossby waves which propagate westward faster than short waves are responsible for the westward stretching. This response can also be understood as the deformation of each vortex into a so-called $\beta$-plume. Stommel (Reference Stommel1982) first described $\beta$-plumes when trying to understand the circulation induced by rising water from hydrothermal vents in the Pacific, and his theory was further developed by Davey & Killworth (Reference Davey and Killworth1989) who considered the evolution of buoyancy sources on a $\beta$-plane. In both cases, the convergence or divergence of fluid around the perturbation generates cyclonic or anticyclonic motions which are subsequently elongated westward due to the emission of Rossby waves (see the review and experiments by Afanasyev & Ivanov (Reference Afanasyev and Ivanov2019)). As a consequence, an east–west asymmetry in the radial component develops, which is clear in figure 9(b,f): the advection of the background potential vorticity leads to a weakening of the flow on the west (left) side of each vortex, and a strengthening on their east (right) side, for both cyclones and anticyclones. More interesting for us is the fact that this asymmetry leads to a prograde momentum convergence towards the region of generation of the vortices as demonstrated by figure 9(d,h). Figure 10 compares the Reynolds stresses divergence profile in the experiment and in the linear solution. There is indeed an eastward acceleration of the zonal flow in the forcing region, located between two westward acceleration regions. This mechanism thus explains the experimental regime I, i.e. the formation of prograde jets flanked by two retrograde jets above each forcing ring. Vallis (Reference Vallis2006, chapter 12) provides an overview of this mechanism in the framework of barotropically forced surface westerlies in the atmosphere. Our experiment then shows that this mechanism is robust since the generated prograde jets persist at later times, even in the nonlinearly saturated regime.

Figure 9. Comparison between the experimental flow and the linear QG model. (a–d) Experimental flow measured 2 s after the third forcing was turned on. The experimental data, originally obtained in a cylindrical geometry, is remapped for the sake of comparison with our cartesian model. Here $R_3$ is the radius of the third forcing ring. Only eight vortices are visualized, but the third ring contains 18 vortices. (e–h) Solution of the linear model (4.17) at time $t=$ 2 s. The model parameters ($\alpha$, $\beta$, $q_m$, $k_f$ and $U$) are estimated from the experimental parameters. (a,e) Zonal velocity perturbation. (b,f) Radial velocity perturbation (note that $u_\rho$ should be compared with $-v$). (c,g) Vorticity. (d,h) Reynolds stresses divergence (zonal flow acceleration).

Figure 10. Comparison of the Reynolds stresses divergence between the experiment and the linear model solution. (a) Experimental Reynolds stresses divergence (4.8). (b) Reynolds stresses divergence computed from the linear model (4.19). Both the experiment and the linear model show the generation of a prograde jet at the forcing location, flanked with two retrograde jets.

4.3. Resonance of the Reynolds stresses and associated feedback

The previous section shows that our experiment can be successfully described by a simple 2-D QG model incorporating only the $\beta$-effect and the bottom friction, at least in the linear regime. We can thus use this model to investigate the feedback that the zonal flow can have on the Reynolds stresses divergence $\mathcal {R}$ and study the possibility of bistability. This is the goal of the present section.

For simplicity, we now forget about the specific geometry of the experimental forcing, and use a generic forcing consisting of a doubly periodic array of vortices with a wavelength comparable to the experimental one, i.e.

(4.22)\begin{equation} q(x,y)=q_m \cos(k_f x) \cos(l_f y), \end{equation}

with $q_m=0.5\ \mathrm {s}^{-2}$ the forcing amplitude and $k_f=l_f$ the forcing wavenumber. In the following, we work with $k_f=63\ \mathrm {rad}\ \mathrm {m}^{-1}$ corresponding to a typical forcing wavelength of 10 cm. We show in appendix C that the forcing scale has only a small influence on the physical mechanism presented. The linear response obtained with this forcing is the superposition of the response to each forcing line, with a prograde acceleration above each forcing horizontal line, and retrograde accelerations in between.

To study the feedback of the zonal flow, we solve for the stationary linear response to this forcing ((4.17) with $t \rightarrow +\infty$) for various amplitudes of the background zonal flow $U$. We recall the local nature of our analysis: the background zonal flow $U$ should be seen as the – uniform – zonal flow at the core of a jet. For each solution, we extract the streamfunction amplitude $\vert \hat {\psi } \vert$ and the Reynolds stresses divergence at $y=0$, and represent them in figure 11 as a function of the zonal flow $U$. Given the amplitude of the streamfunction

(4.23)\begin{equation} \vert \hat{\psi} \vert = \frac{q_m}{k_f^{2} + l_f^{2}} \left( {\omega(k_f,l_f)^{2}+\omega_E(k_f,l_f)^{2}} \right)^{-1/2}, \end{equation}

we expect a resonance of the linear response when $\omega =0$, or, in other words, when the directly forced Rossby waves are stationary, i.e.

(4.24)\begin{equation} U = \frac{\beta}{k_f^{2} + l_f^{2}} = -c, \end{equation}

where $c$ is the directly forced Rossby waves phase speed in the absence of zonal flow (Vallis Reference Vallis2006, p. 542). The resonant amplification of the response amplitude is visible in figure 11(a). Then, the amplitude of the Reynolds stresses is expected to be proportional to the squared amplitude of the streamfunction (see (4.18) and (4.19)), leading to

(4.25)\begin{equation} \vert \mathcal{R}\vert (U) \propto \vert \hat{\psi} \vert^{2} \propto \frac{1}{\omega^{2} + \omega_E^{2}} \propto \frac{1}{\left(1+\dfrac{U}{c}\right)^{2} + \gamma^{2} } , \end{equation}

where $\gamma$ is a non-dimensional parameter characterizing the Rossby waves damping

(4.26)\begin{equation} \gamma^{2} = \left( \frac{{\omega_E}}{k_f c}\right)^{2}. \end{equation}

For the problem to be analytically tractable, we chose to model the resonant curve $\mathcal {R}(U)$ plotted in figure 11(b) with a parametrized Lorentzian. We thus forget about the spatial structure of the momentum flux convergence, and set

(4.27)\begin{equation} \mathcal{R}(U) = \mathcal{R}_m ~\frac{1}{\gamma^{2} + \left(1+\dfrac{U}{c}\right)^{2}}. \end{equation}

Doing so, we focus on the amplitude of the response rather than on the details due to our particular choice of forcing. The important physical effects of the zonal flow, the $\beta$-effect and the friction, are contained in the Lorentzian and influence the position and flatness of the resonant peak. Figure 11(b) shows that the amplitude of the Reynolds stresses is indeed largely dominated by this resonant amplification. Hence, we do not loose any important feature by modelling $\mathcal {R}$ with (4.27) which has the advantage of making the problem analytically tractable.

Figure 11. (a) Streamfunction amplitude as a function of the background zonal flow $U$. Each symbol is a solution of the linear model with a different background zonal flow. The black line is the corresponding analytical amplitude in the case of a doubly periodic forcing (4.23). (b) Reynolds stresses at $y=0$ as a function of the background zonal flow $U$. Each symbol is a solution of the linear model with a different background zonal flow. The black line represents the best fit of the Lorentzian given by (4.27) where the amplitude $\mathcal {R}_m$ is a free parameter. In both panels, the vertical dashed black line shows the location of $U=-c$, where $c$ is the forced Rossby wave phase speed (4.24).

It is now clear that the momentum flux can lead to abrupt transitions: on the left side of the resonant peak, any increase of a zonal flow $U$ leads to an increased prograde momentum convergence, and an increased acceleration of this zonal flow. We have thus identified a potential positive feedback mechanism, provided that it is not cancelled by the negative feedback of the viscous dissipation. Note that since the Rossby waves propagate in the westward direction, this resonance can only occur in an eastward jet ($U>0$).

4.4. Stationary solutions and linear stability

The linear QG model showed the resonant amplification of the wave-induced Reynolds stresses when the zonal flow is such that the directly forced Rossby waves are stationary. Our goal is now to verify whether this feedback of the zonal flow can explain the transition and bistability observed in our experiment.

Consistently with our local approach, we consider a minimal model where the zonal flow $U(t)$ is assumed to be only time dependent, with no spatial modulation. Such a uniform zonal flow is sustained by the Reynolds stresses divergence $\mathcal {R}$ and dissipated by the linear friction due to the Ekman pumping, i.e.

(4.28)\begin{equation} \frac{ \partial U}{ \partial t} = \mathcal{R}(U) - \alpha U, \end{equation}

with $\mathcal {R}(U)$ the Lorentzian given by (4.27). The stationary solutions of this equation are the roots of the third-order polynomial

(4.29)\begin{equation} P(U)=U^{3} + 2c U^{2} + (1+\gamma^{2})c^{2} U - \frac{\mathcal{R}_m c^{2}}{\alpha}. \end{equation}

Depending on the sign of the discriminant of $P$, one or three stationary solutions can exist, as represented in figure 12. We denote $U_1$, $U_2$ and $U_3$ those three solutions such that $U_1<U_2<U_3$. For three stationary solutions to exist, i.e. bistability to be possible, a necessary but not sufficient condition is

(4.30)\begin{equation} \gamma^{2} < \tfrac{1}{3}. \end{equation}

When this condition is satisfied, the sufficient condition for three solutions to exist is that

(4.31)\begin{equation} \mathcal{R}_m \in [\mathcal{R}_{1}, \mathcal{R}_{2}] = -\tfrac{2}{27}\alpha c \big[9\gamma^{2}+1-\sqrt{\varGamma}, 9\gamma^{2}+1+\sqrt{\varGamma}\big], \end{equation}

with $\varGamma = (1-3\gamma ^{2})^{3}$. Physically, the first condition (4.30) means that bistability can never exist if the Rossby waves are too strongly damped. The second condition shows that even when the friction is not too high, three stationary solutions exist only for a given range of the forcing amplitude $\mathcal {R}_m$. As represented in figure 12, if the forcing is too high, only the super-resonant solution $U_3$ exists. Conversely, if the forcing is too weak, only the low-amplitude solution $U_1$ can exist.

Figure 12. Visualization of the stationary solutions ($U_1,U_2,U_3$) of the zonal flow evolution equation (4.28). Except the forcing amplitude, all the parameters are fixed and equal to the experimental parameters ($\alpha = 5.6\times 10^{-3}\ \mathrm {s}^{-1}$, $k_f=l_f= 63\ \mathrm {rad}\ \mathrm {m}^{-1}$, $c=6.3\times 10^{-3}\ \mathrm {m}\ \mathrm {s}^{-1}$, $\gamma ^{2}=1.99\times 10^{-4}$). (a) Illustrative case of a small forcing amplitude $\mathcal {R}_m < \mathcal {R}_{1}$ and a small friction or a too high friction. Note that $U_1$ is very small, but not zero. (b) Illustrative case of a forcing in the range $\mathcal {R}_m \in [\mathcal {R}_{1}, \mathcal {R}_{2}]$ and a small friction. (c) Illustrative case of a high forcing amplitude $\mathcal {R}_m > \mathcal {R}_{2}$ and a small friction. (d) Amplitude of the three stationary solutions as a function of the forcing amplitude $\mathcal {R}_m$. The black dots represent the unstable solution whereas the white ones are stable. The shaded area is the bistable range, bounded by $\mathcal {R}_1$ and $\mathcal {R}_2$ (4.31).

To investigate the linear stability of the stationary solutions, we go back to the zonal flow evolution equation (4.28). We linearize the nonlinear operator $\mathcal {R}$ around the stationary state $U_s$ and denote $U'(t)=U(t)-U_s$ the perturbed zonal flow, to obtain the perturbations evolution equation

(4.32)\begin{equation} \frac{ \partial U'}{ \partial t} = \left.\frac{\mathrm{d} \mathcal{R}}{\mathrm{d}U} \right|_{U_s} U' - \alpha U'. \end{equation}

We seek $U'$ under the form $U'=U'_0\,\textrm {e}^{\sigma t}$, where $\sigma$ is the growth rate of the perturbation. Substituting into (4.32) leads to

(4.33)\begin{align} \sigma = \left.\frac{\mathrm{d} \mathcal{R}}{\mathrm{d}U} \right|_{U_s} - \alpha. \end{align}

The stationary solution $U_s$ is unstable if and only if the growth rate $\sigma >0$. This condition is always verified for the second stationary solution $U_2$ (the sub-resonant one), whereas $U_1$ and $U_3$ are stable stationary solutions.

We plot in figure 12(d) the amplitude of the stationary solutions obtained for varying forcing amplitudes $\mathcal {R}_m$, all the other parameters being fixed and equal to the experimental parameters ($\alpha = 5.6\times 10^{-3}\ \mathrm {s}^{-1}$, $k_f=l_f= 63\ \mathrm {rad}\ \mathrm {m}^{-1}$, $c=6.3\times 10^{-3}\ \mathrm {m}\ \mathrm {s}^{-1}$, $\gamma ^{2}=1.99\times 10^{-4}$). The three stationary solutions branches are visible, and coexist only in a given range of forcing amplitude bounded by $\mathcal {R}_1$ (purple line) and $\mathcal {R}_2$ (blue line) given by (4.31). This figure also demonstrates the bistability and possibility of an hysteresis: for an experiment with increasing forcing, the transition $U_1 \rightarrow U_3$ happens for $\mathcal {R}_m=\mathcal {R}_2$. At this forcing, there is a saddle-node bifurcation ($S_2$) through which $U_1$ loses its stability. But if the forcing amplitude is then decreased, the observed solution will remain $U_3$ until the forcing reaches $\mathcal {R}_1 < \mathcal {R}_2$ where there is another saddle-node bifurcation ($S_1$), and we go back to the lower branch solution $U_1$. Our model predicts that close to the transition, the amplitude of the zonal flow on the lower branch is of $U\approx 2\ \textrm {mm}\ \mathrm {s}^{-1}$, and $U\approx -c=6.3\ \mathrm {mm}\ \mathrm {s}^{-1}$ for the upper branch. Finally, this model allows us to define a typical velocity expected at the transition, to compare with the forcing RMS velocity $U_f$ used to characterize the experimental hysteresis, on figure 7. From (4.31), the Reynolds stresses at the transition are typically of $\mathcal {R}_{m,t} \sim \alpha \vert c\vert (1+9\gamma ^{2})$. A typical transition velocity can be obtained supposing that the forcing is balanced by friction $U_t \sim \mathcal {R}_{m,t}/\alpha$, leading to

(4.34)\begin{align} U_t \sim \vert c \vert (1+9\gamma^{2}). \end{align}

If we use $U_f^{*}=U_f/U_t$ as a non-dimensional forcing, then the transition should always occur at $U_f^{*}$ of order unity. The dependence on the $\beta$-effect, the friction and the forcing scale are incorporated in $c$ and $\gamma$. Note that the width of the bistable zone will however vary depending on $\gamma$.

4.5. Comparison of the experimental to theoretical transition

In this section we report additional experimental observations that support the mechanism of Rossby waves resonance.

First, the amplitudes of the zonal flow expected in the two steady states $U_1$ and $U_3$ are very close to the experimental ones, where $U_1$ represents regime I and $U_3$ regime II, as can be seen in figure 7. Another way to see it is by comparing the zonal flow amplitude with the phase speed of the – non-advected – directly forced Rossby waves: if it is higher than $-c$, the system is in the super-resonant steady state, whereas if it is lower than $-c$, the system is sub-resonant. In the experimental set-up, with $k_f \approx 63\ \mathrm {rad} \ \mathrm {m}^{-1}$ and $\beta \approx 50\ \mathrm {m}^{-1}\ \mathrm {s}^{-1}$, the phase speed of the directly forced Rossby waves is $c \approx -6.3\ \mathrm {mm}\ \mathrm {s}^{-1}$, without advection. First, figure 13 demonstrates that the forcing indeed excites Rossby waves. The radial component of the velocity at a given radius exhibits patterns that move

Figure 13. Rossby waves excited by the forcing. The third forcing ring ($R_3=0.214$ m) is turned on at its maximum power from time $t=0$ to 75 $t_R$. (a) Time evolution of the radial component of the velocity at a fixed radius $\rho =0.17$ m. (b) Same data band-passed filtered between 0.06 and 0.12 Hz. The slope of the dashed line is the non-Doppler-shifted phase speed of the Rossby waves excited by the forcing (4.24) $c \approx -6.3\ \mathrm {mm}\ \textrm {s}^{-1}$. The slope of the full line is the Doppler-shifted phase speed $\langle u_\phi \rangle _\phi + c$. It increases in absolute value since the zonal flow increases when the forcing is turned on, and at the chosen radius the flow is retrograde, as represented on the top panel.

in the retrograde direction (i.e. decreasing $\phi$), with the same wavelength as the forcing and at the Doppler-shifted speed $\langle u_\phi \rangle _\phi + c$. Then, given this typical phase speed $c$, figure 7 shows that in regime I the zonal flow is sub-resonant ($\bar {u}^{\textit{RMS}} < -c$), and super-resonant in regime II ($\bar {u}^{\textit{RMS}} > -c$), which is consistent with our model. Furthermore, the model predicts that when the forcing is decreased and $U_3$ (regime II) loses its stability ($S_2$ in figure 12(d)), the zonal flow is quasi-resonant ($\langle u_\phi \rangle _\phi \approx -c$ ) which is again compatible with the measured velocity at the transition $\textrm {II}\rightarrow \textrm {I}$ in figure 7. Note that if this equilibrium close to resonance can exist ($U_3$), our analysis only gives an explanation for its origin, ultimately, the nonlinearities are responsible for locking the system into a near-resonant state. Note finally that the bistability range and zonal flow amplitudes are slightly varying with the range of possible forcing wavenumbers for our experiment (figure 16(b) in appendix C). For a closer match between the predicted zonal flow amplitudes and the measured ones, we would choose $k_f = 57\ \mathrm {rad}\ {\mathrm {m}^{-1}}$, which is reasonable given the uncertainty in the relevant forcing scale in our set-up.

Second, we have not yet discussed the inherent asymmetry between the prograde and the retrograde jets in our model. Since the Rossby waves propagate in the retrograde direction, their only way to become stationary and resonate with the fixed forcing is within a prograde jet. This implies that the transition only occurs in the prograde jets which is again consistent with our experimental observations. Indeed, the prograde jets govern the dynamics, by increasing in amplitude and merging, whereas the retrograde flow seems to adjust passively to the transition (e.g. figure 6). If we suppose that the forcing imparts no net angular momentum to the fluid, then the integral of the angular momentum per unit mass, $\rho (u_\phi + \varOmega \rho )$, over the domain should be constant. If an eastward acceleration is produced by the transition then this convergence of prograde angular momentum should be balanced by negative angular momentum elsewhere (see the Reynolds stresses profile in figure 10). The retrograde flows in our experiment seem to arise as such, which is consistent with the fact that the retrograde flow is smooth whereas the prograde one strongly interacts with vortices (figure 3).

Third, our model implies that the faster the Rossby waves, the stronger the forcing will have to be to reach the transition (increasing $U_t$), and conversely for slower Rossby waves. We performed experiments at 80 RPM and 60 RPM (not shown) instead of the 75 RPM for which the bottom plate was designed. In these experiments, $\beta$ is no longer uniform, but slightly varying with radius. It is higher at 80 RPM due to the increased curvature of the free surface ($\beta _{80}\in [57,73]\ \mathrm {m}^{-1}\ {\mathrm {s}^{-1}}$ with a mean at $65.5\ \mathrm {m}^{-1}\ \mathrm {s}^{-1}$) and weaker at 60 RPM ($\beta _{60}\in [20,30]\ \mathrm {m}^{-1}\ \mathrm {s}^{-1}$ with a mean at $22.8 \ \mathrm {m}^{-1}\ \mathrm {s}^{-1}$). As a consequence, the forced Rossby waves are respectively faster and slower ($|c|=$ 8.3 and $2.9\ \mathrm {mm}\ \mathrm {s}^{-1}$). Consistently with our model, we observed that for a similar forcing ($U_f=4.0\ \mathrm {mm}\ \mathrm {s}^{-1}$), regime II is obtained at 75 RPM whereas regime I is observed at 80 RPM. Conversely, for a forcing at which regime I is observed at 75 RPM ($U_f=2.9\ \mathrm {mm}\ \mathrm {s}^{-1}$), regime II is obtained at 60 RPM. The transition thus occurs at larger forcing amplitudes for an increased $\beta$-effect. For a more quantitative view of the sensitivity to the $\beta$-effect, we refer the reader to appendix C.

Finally, we wish to discuss the local aspect of our model relatively to the experiment. Our model only explains the local feedback mechanism inside of a given prograde jet. How the global system responds is a different question. Based on the Rhines scale $L_R \sim (u^{\textit{RMS}}/\beta )^{1/2}$ (Rhines Reference Rhines1975), where we recall that $u^{\textit{RMS}}$ takes into account all the components of the flow (3.1), we do expect an increase of the jets width during the transition, since the RMS eddy velocity increases. However, our model does not explain why the jets merge during the transition. In addition, this model does not rule out the possibility of the coexistence of regime I and regime II flows side by side. For instance, the most external forcing ring in our experiment (above C6) has a greater pressure loss because of the number of hoses (38) and it is possible that on this ring, the forcing is never super-resonant. This suggests that regional stable equilibria, where both regimes can be locally sustained in distinct regions of space, may exist. Finally, we observe that the merger events during the transition are associated with a radial shift of the jets leading to an uncorrelation between the jets position and the forcing rings. This further suggest that a local approach will not be sufficient to explain the saturation in regime II. Instead, it may be relevant to adopt a global approach based, for instance, on the turbulent properties and energy transfer of anisotropic turbulence on a $\beta$-plane. This type of approach has led to the development of the theory of zonostrophic turbulence (Galperin et al. Reference Galperin, Sukoriansky, Young, Chemke, Kaspi, Read and Dikovskaya2019, and references therein). Due to its fast rotation and owing to the Taylor–Proudman theorem, the flow is quasi two dimensional and may bear an inverse turbulent energy cascade. Because of the $\beta$-effect and associated Rossby waves, the energy transfer becomes anisotropic and redirected towards zonal currents. Ultimately, the large-scale drag halts the expansion of the inverse cascade (Sukoriansky et al. Reference Sukoriansky, Galperin and Dikovskaya2002, Reference Sukoriansky, Dikovskaya and Galperin2007). Determining whether such a theory is valid to explain the nonlinear saturation in regime II is beyond the scope of the present work, but will be investigated in a separate study. Note also that the spectral analysis which can be found in Cabanes et al. (Reference Cabanes, Aurnou, Favier and Le Bars2017), Cabanes, Favier & Le Bars (Reference Cabanes, Favier and Le Bars2018) for the previous version of the experiment and corresponding DNS supports the idea that the second regime is close to zonostrophic turbulence. To conclude, we wish to underline that Cabanes et al. (Reference Cabanes, Aurnou, Favier and Le Bars2017) probably only observed regime II in their experimental set-up because their forcing was of larger amplitude than in the present study, thus, probably always super-resonant.