1 Introduction
The quintessential paradigm for investigating the fundamentals of rotating thermally forced flows is provided by rotating Rayleigh–Bénard convection (RRBC) in a horizontal layer rotating about a vertical axis with constant angular velocity
$\unicode[STIX]{x1D6FA}$
, i.e. convection in a layer of Boussinesq fluid confined between flat horizontal rigidly rotating upper and lower boundaries maintaining a destabilizing temperature jump
$\unicode[STIX]{x0394}T>0$
. Of particular relevance to the dynamics of stellar and planetary interiors, planetary atmospheres and terrestrial oceans is the regime of geostrophic turbulence where fluid motions are sufficiently constrained by rotation to enforce pointwise balance between the pressure gradient and the Coriolis force, otherwise known as geostrophic balance. This balance is characteristic of rapidly rotating systems for which the convective Rossby number is small,
$$\begin{eqnarray}Ro\equiv \frac{\sqrt{g\unicode[STIX]{x1D6FC}\unicode[STIX]{x0394}T/H}}{2\unicode[STIX]{x1D6FA}}=\sqrt{\frac{Ra}{Pr}}E\ll 1.\end{eqnarray}$$
This number denotes the ratio of the rotation time scale to the free-fall or free-rise time scale for a parcel of fluid with temperature difference
$\unicode[STIX]{x0394}T$
relative to the ambient fluid. Here,
$g$
denotes the acceleration due to gravity,
$\unicode[STIX]{x1D6FC}$
is the thermal expansion coefficient and
$H$
is the layer depth. The second equality rewrites this definition in terms of quantities familiar from studies of Rayleigh–Bénard convection: the Rayleigh number
$Ra$
, the Ekman number
$E$
and the Prandtl number
$Pr$
, given by
$$\begin{eqnarray}Ra=\frac{g\unicode[STIX]{x1D6FC}\unicode[STIX]{x0394}TH^{3}}{\unicode[STIX]{x1D708}\unicode[STIX]{x1D705}},\quad E=\frac{\unicode[STIX]{x1D708}}{2\unicode[STIX]{x1D6FA}H^{2}},\quad Pr=\frac{\unicode[STIX]{x1D708}}{\unicode[STIX]{x1D705}}.\end{eqnarray}$$
These measure respectively the strength of the thermal forcing, the importance of viscous diffusion relative to rotation and the thermometric properties of the fluid through its kinematic viscosity
$\unicode[STIX]{x1D708}$
and thermal diffusivity
$\unicode[STIX]{x1D705}$
. Since
$Ra$
must be large to generate turbulence, the presence of geostrophic turbulence requires that
$E$
be correspondingly small. This regime can be accessed by exploring the simultaneous limits
$Ra\rightarrow \infty$
,
$E\rightarrow 0$
such that
$Ro\ll 1$
. Unfortunately, this regime is inaccessible to both direct numerical simulations (DNS) of the Navier–Stokes equations (NSEs) and laboratory investigations. To date, the lowest achievable Ekman numbers are in the neighbourhood of
$E=O(10^{-7})$
, whereas an adequate exploration of geostrophic turbulence requires Ekman numbers that are much lower.
An alternative and fruitful approach that has recently been advanced (Julien & Knobloch Reference Julien and Knobloch2007) employs an asymptotic reformulation of the NSEs for incompressible thermal convection valid in the limit
$Ro\downarrow 0$
to derive a reduced system of partial differential equations (PDEs) called the non-hydrostatic quasi-geostrophic equations (NH-QGEs). This reduced system filters out fast inertial waves and thin Ekman boundary layers, and is therefore amenable to extensive numerical explorations. These have been validated qualitatively by DNS studies at moderately low
$Ro$
and
$E$
(Favier, Silvers & Proctor Reference Favier, Silvers and Proctor2014; Guervilly, Hughes & Jones Reference Guervilly, Hughes and Jones2014; Stellmach et al.
Reference Stellmach, Lischper, Julien, Vasil, Cheng, Ribeiro, King and Aurnou2014; Plumley et al.
Reference Plumley, Julien, Marti and Stellmach2016), and have enabled a comprehensive mapping of the
$Ra$
–
$Pr$
space (Julien et al.
Reference Julien, Rubio, Grooms and Knobloch2012; Rubio et al.
Reference Rubio, Julien, Knobloch and Weiss2014). Figure 1 illustrates volume renderings of the vorticity and streamfunction fields in the regime of geostrophic turbulence in a domain of unit aspect ratio in the horizontal. A remarkable feature of this state is the presence of a strong inverse energy cascade resulting in a box-scale condensate in the form of a vortex dipole. The feature appears to be barotropic (i.e. depth-independent), as demonstrated by, and most noticeable in, the barotropic vorticity and streamfunction fields (b,c). The barotropic dynamics satisfies the two-dimensional (2D) barotropic vorticity equation forced baroclinically by the underlying depth-dependent geostrophic turbulence and damped by viscosity (Rubio et al.
Reference Rubio, Julien, Knobloch and Weiss2014). The energetics of this process can be viewed as a two-way barotropic–baroclinic interaction: the barotropic dynamics is directly forced by and extracts energy from the convective (baroclinic) dynamics. This interaction can be highly efficient in that the barotropic vortex is capable of growing to large amplitudes with little impact on the underlying geostrophic turbulence, as measured by the small adjustment in convective (baroclinic) kinetic energy when the pathway to exciting the barotropic manifold is switched on at
$t=0$
(Rubio et al.
Reference Rubio, Julien, Knobloch and Weiss2014). At the same time, the baroclinic fluctuations are aligned by the barotropic flow, leading to a self-sustaining process.
Volume render of geostrophic turbulence at
$RaE^{4/3}=90,Pr=1$
. Top view of (a) the total vertical vorticity
$\unicode[STIX]{x1D701}$
, (b) the barotropic vorticity
$\langle \unicode[STIX]{x1D701}\rangle$
and (c) the barotropic streamfunction
$\langle \unicode[STIX]{x1D6F9}\rangle$
.
The generation of large-scale structure in turbulent flows is primarily investigated in 2D, and focuses on the 2D NSEs with damping (provided by Rayleigh friction) and dissipation (provided by viscosity), and driven by externally imposed noise, usually taken to be white (Smith & Yakhot Reference Smith and Yakhot1994; Bouchet & Simonnet Reference Bouchet and Simonnet2009; Frishman, Laurie & Falkovich Reference Frishman, Laurie and Falkovich2017). Numerical study of this system in a periodic domain with a square aspect ratio also realizes condensation into a box-scale vortex structure. When the aspect ratio becomes elongated, turbulent jets oriented parallel to the short side form instead of a box-scale vortex. As shown recently (Frishman et al. Reference Frishman, Laurie and Falkovich2017), these jets may be populated by large numbers of prominent vortices embedded in an anisotropic turbulent background state. This type of condensation process has also been examined using ideas from equilibrium statistical mechanics (Bouchet & Simonnet Reference Bouchet and Simonnet2009; Bouchet & Venaille Reference Bouchet and Venaille2012), which predicts a transition from a box-scale vortex to a jet state as the aspect ratio increases and the domain becomes rectangular. Both the simulations and theory find that jets are already present when the elongation is of order 10 %. However, the statistical approach describes only box-scale structures and so cannot examine the finer details of the turbulent jets it predicts. In addition, the 2D system is a driven dissipative system, and any predictions from equilibrium statistical mechanics have to be treated with caution despite the similarities between the predictions and the numerical simulations.
In the present paper, we also identify a transition between a box-scale vortex dipole and jets, and find that jets first appear when the elongation is of order 10 %. However, our system is quite different from the 2D damped noise-driven NSEs studied by Smith & Yakhot (Reference Smith and Yakhot1994), Bouchet & Simonnet (Reference Bouchet and Simonnet2009) and Frishman et al. (Reference Frishman, Laurie and Falkovich2017) in that our equations are fully three-dimensional (3D) and the fluctuations driving the condensation process have to be determined self-consistently with the vortices or jets they produce. Thus, the noise process is both anisotropic and non-white, and the physics behind the condensation process necessarily differs. Our conclusion, elaborated further below, is that the condensation process is highly robust, with respect to both the physics behind the fluctuations and the substantially different nature of the governing equations themselves.
2 The non-hydrostatic reduced equations
A complete derivation and discussion of the NH-QGEs is presented in Sprague et al. (Reference Sprague, Julien, Knobloch and Werne2006), Julien & Knobloch (Reference Julien and Knobloch2007) and Julien et al. (Reference Julien, Rubio, Grooms and Knobloch2012). The equations are obtained as the leading-order reduction of the incompressible NSEs based on a multiscale asymptotic expansion in
$Ro=E^{1/3}\equiv \unicode[STIX]{x1D716}\ll 1$
employing a small scale
$L=\unicode[STIX]{x1D716}H\ll 1$
as well as the large vertical scale
$H$
. For the case of stress-free upper and lower boundaries, the leading-order velocity field
$\boldsymbol{u}\equiv (\boldsymbol{u}_{\bot },W)$
is in geostrophic balance, i.e.
$\widehat{\boldsymbol{z}}\times \boldsymbol{u}_{\bot }=-\unicode[STIX]{x1D735}_{\bot }p$
. It follows that the horizontal velocity field is non-divergent, with
$\boldsymbol{u}_{\bot }\equiv (u,v,0)=(-\unicode[STIX]{x2202}_{y}\unicode[STIX]{x1D6F9},\unicode[STIX]{x2202}_{x}\unicode[STIX]{x1D6F9},0)$
, where the pressure
$p\equiv \unicode[STIX]{x1D6F9}$
is the geostrophic streamfunction. The reduced NH-QGEs governing the motion of the fluid are
$$\begin{eqnarray}\displaystyle & \displaystyle D_{t}^{\bot }\unicode[STIX]{x1D701}-\unicode[STIX]{x2202}_{Z}W=\unicode[STIX]{x1D6FB}_{\bot }^{2}\unicode[STIX]{x1D701}, & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle D_{t}^{\bot }W+\unicode[STIX]{x2202}_{Z}\unicode[STIX]{x1D6F9}=\frac{\widetilde{Ra}}{Pr}\unicode[STIX]{x1D6E9}^{\prime }+\unicode[STIX]{x1D6FB}_{\bot }^{2}W, & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle D_{t}^{\bot }\unicode[STIX]{x1D6E9}^{\prime }+W\unicode[STIX]{x2202}_{Z}\overline{\unicode[STIX]{x1D6E9}}=\frac{1}{Pr}\unicode[STIX]{x1D6FB}_{\bot }^{2}\unicode[STIX]{x1D6E9}^{\prime }, & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x2202}_{\unicode[STIX]{x1D70F}}\overline{\unicode[STIX]{x1D6E9}}+\unicode[STIX]{x2202}_{Z}(\overline{W\unicode[STIX]{x1D6E9}^{\prime }})=\frac{1}{Pr}\unicode[STIX]{x2202}_{ZZ}\overline{\unicode[STIX]{x1D6E9}}, & \displaystyle\end{eqnarray}$$
capturing respectively the evolution of vertical vorticity
$\unicode[STIX]{x1D701}=-\unicode[STIX]{x1D6FB}_{\bot }^{2}\unicode[STIX]{x1D6F9}$
, vertical velocity
$W$
and temperature
$\unicode[STIX]{x1D6E9}=\overline{\unicode[STIX]{x1D6E9}}+\unicode[STIX]{x1D716}\unicode[STIX]{x1D6E9}^{\prime }$
at the reduced Rayleigh number
$\widetilde{Ra}=Ra\unicode[STIX]{x1D716}^{4}$
for a given Prandtl number
$Pr$
. The temperature is decomposed into a mean (horizontally averaged) component
$\overline{\unicode[STIX]{x1D6E9}}$
evolving on the slow time scale
$\unicode[STIX]{x1D70F}=\unicode[STIX]{x1D716}^{2}t$
and a small fluctuating component
$\unicode[STIX]{x1D6E9}^{\prime }$
. Here,
$D_{t}^{\bot }\equiv \unicode[STIX]{x2202}_{t}+\boldsymbol{u}_{\bot }\boldsymbol{\cdot }\unicode[STIX]{x1D735}_{\bot }$
denotes the horizontal material derivative. The system is solved with impenetrable stress-free fixed-temperature boundary conditions,
$$\begin{eqnarray}\overline{\unicode[STIX]{x1D6E9}}=1,W=\unicode[STIX]{x1D6E9}^{\prime }=0\text{ at }Z=0,\quad \text{and}\quad \overline{\unicode[STIX]{x1D6E9}}=0,W=\unicode[STIX]{x1D6E9}^{\prime }=0\text{ at }Z=1.\end{eqnarray}$$
The NH-QGEs are discretized in the horizontal and vertical spatial directions using a sparse Fourier–Chebyshev spectral decomposition (Julien & Watson Reference Julien and Watson2009). They are then time evolved using a third-order semi-implicit explicit Runge–Kutta scheme.
3 Results and discussion
In the following, we present results for RRBC obtained from a series of simulations of the NH-QGE system (2.1)–(2.5) performed within the geostrophic turbulence regime at
$\widetilde{Ra}=90$
,
$Pr=1$
. The horizontal aspect ratio 1 :
$\unicode[STIX]{x1D6E4}$
is varied from 1 : 1 to 1 : 6. All cases exhibit similar efficiency in heat transport as measured by the Nusselt number, namely
$Nu=36.84\pm 1.97$
.
Top view of a suite of simulations using the NH-QGEs for RRBC with impenetrable stress-free fixed-temperature boundaries. All simulations were performed at
$\widetilde{Ra}=90$
,
$Pr=1$
. (a) The total vorticity
$\unicode[STIX]{x1D701}$
, (b) the barotropic vorticity
$\langle \unicode[STIX]{x1D701}\rangle$
and (c) the barotropic streamfunction
$\langle \unicode[STIX]{x1D6F9}\rangle$
. In units of a fixed
$x$
-horizontal box scale
$L_{x}=10L_{c}$
(where
$L_{c}=4.82E^{1/3}H$
is the critical wavelength for convection), the domain aspect ratio
$(1:\unicode[STIX]{x1D6E4})L_{x}$
increases from 1 : 1 (upper row) to 1 : 1.1 (second row) and incrementally from 1 : 2 to 1 : 6 (remaining rows). It should be noted that
$x$
is in the vertical direction with
$y$
plotted horizontally. Banded flow of alternating cyclonicity when
$\unicode[STIX]{x1D6E4}=2$
, and a doubling, tripling and quadrupling of the band structure at
$\unicode[STIX]{x1D6E4}=3$
,
$5$
and
$6$
respectively is evident. The second row (
$\unicode[STIX]{x1D6E4}=1.1$
) shows that the transition to banded flow occurs upon even a small departure from isotropy. The resolution in the
$x,y,Z$
directions is
$144\times 144\unicode[STIX]{x1D6E4}\times 192$
.
3.1 Visualizations
Figure 2 depicts the top views of volume renderings of the total vertical vorticity
$\unicode[STIX]{x1D701}$
(a), the barotropic vorticity
$\langle \unicode[STIX]{x1D701}\rangle$
(b) and the barotropic streamfunction
$\langle \unicode[STIX]{x1D6F9}\rangle$
(c). Here,
$\langle \cdots \rangle$
indicates an average in the vertical. The barotropic vortex dipole present at aspect ratio 1 : 1 (top row) is replaced by a state of approximately parallel (i.e. banded) flow consisting of an alternating sequence of cyclonic and anticyclonic vortical bands. This transition first occurs at an aspect ratio of approximately 1 : 1.1 (second row) and is most prominently revealed in the rendering of the barotropic streamfunction;
$\langle \unicode[STIX]{x1D701}\rangle$
exhibits greater spatial complexity due to higher spatial derivatives. We see that strong vortical eddies persist within a band of given cyclonicity and that small-scale geostrophic turbulence is globally advected and organized by the banded large-scale flow. The latter resembles the small-scale filamentary structures in figure 1(a). From a more global perspective, the aspect ratio
$\unicode[STIX]{x1D6E4}$
provides a selection mechanism for the number of alternating bands or jets. It can be seen that an increase from
$\unicode[STIX]{x1D6E4}=2$
to
$\unicode[STIX]{x1D6E4}=3$
results in a doubling of the number of alternating bands from one to two. The latter persists for
$\unicode[STIX]{x1D6E4}=4$
before losing stability to a state of three alternating bands at an aspect ratio of
$\unicode[STIX]{x1D6E4}=5$
and then four at
$\unicode[STIX]{x1D6E4}=6$
.
Inspection of the horizontal velocity fields offers another viewpoint for interpreting the large-scale barotropic structure (figure 3). In the following, we refer to the velocity components
$u=-\unicode[STIX]{x2202}_{y}\unicode[STIX]{x1D6F9}$
and
$v=\unicode[STIX]{x2202}_{x}\unicode[STIX]{x1D6F9}$
as parallel and transverse respectively. The parallel velocity field (a,c,e,g) clearly shows a state resembling Kolmogorov flow, i.e. alternating bands of unidirectional but oppositely directed turbulent flow (a,e). Averaging along the
$x$
direction reveals a jet structure that is sawtooth in shape in the
$y$
or transverse direction (dashed lines, c,g), with small root mean square (RMS) fluctuations (solid lines, c,g). This unidirectional jet structure is found at all sufficiently non-square aspect ratios (figure 2). Observations of
$v$
(b,d,f,h) indicate coexistence with a weak meandering transverse jet (b,f). Averaging in the
$y$
direction (dashed lines, d,h) shows that the transverse mean flow is substantially weaker than the parallel flow (by a factor of approximately 10). The RMS fluctuations about this profile are large, indicating that the transverse mean flow is highly oscillatory (solid lines, d,h). Figure 3(i,j) illustrates the persistence of the parallel and transverse flows with time through Hovmöller diagrams. Here, we see that after an initial transient, stable parallel jets are formed while the antisymmetric transverse flow exhibits random switching.
(a–h) Parallel and transverse barotropic velocity fields
$\langle u\rangle$
(a,c,e,g) and
$\langle v\rangle$
(b,d,f,h). Aspect ratios of 1 : 2 and 1 : 5 are depicted in (a,b,c,d) and (e,f,g,h) respectively. Here,
$\overline{f}^{d}$
indicates averaging of
$f$
in the
$d$
direction. The average and RMS fluctuations of
$\langle u\rangle$
in the parallel direction, i.e. averaging in the
$x$
direction, reveal a Kolmogorov-like velocity profile, and the fluctuations of
$\langle v\rangle$
in the transverse direction reveal the presence of a meandering jet (red dashed lines). The amplitude of the RMS fluctuations about these profiles is shown by solid black lines. (i,j) Hovmöller diagrams of (i)
$\overline{\langle u\rangle }^{x}$
and (j)
$\overline{\langle v\rangle }^{y}$
for an aspect ratio of 1 : 5 at time intervals
$\unicode[STIX]{x0394}t=0.5$
.
(a) Term by term decomposition of the barotropic vorticity equation (3.3) showing unbounded growth in the barotropic mode for
$\unicode[STIX]{x1D6E4}=6$
(black) and
$\unicode[STIX]{x1D6E4}=1$
(grey). (b) One-dimensional barotropic (solid) and baroclinic (dashed) kinetic energy spectra as a function of
$k/k_{Lx}$
for
$\unicode[STIX]{x1D6E4}=6$
, with the corresponding results for
$\unicode[STIX]{x1D6E4}=1$
in grey. Here,
$k_{Lx}=2\unicode[STIX]{x03C0}/L_{x}$
denotes the box wavenumber in the
$x$
direction. Peak energy scales
$k/k_{Lx}=4/6$
for
$\unicode[STIX]{x1D6E4}=6$
and
$k/k_{Lx}=1$
for
$\unicode[STIX]{x1D6E4}=1$
characterize the jet scale.
3.2 Energetics
In the absence of dissipation, and like the incompressible NSEs, all QG systems conserve the volume-averaged energy
${\mathcal{E}}$
and the pointwise potential vorticity
${\mathcal{P}}{\mathcal{V}}$
. For the NH-QGEs (2.1)–(2.4), these are given by
$$\begin{eqnarray}\displaystyle & \displaystyle {\mathcal{E}}=\left\langle \left(\frac{1}{2}\overline{|\unicode[STIX]{x1D735}_{\bot }\unicode[STIX]{x1D6F9}|^{2}+W^{2}}\right)+\frac{\widetilde{Ra}}{Pr}\left(\overline{\unicode[STIX]{x1D6E9}^{\prime 2}}+\unicode[STIX]{x1D716}^{-2}\left(\overline{\unicode[STIX]{x1D6E9}}-z\right)^{2}\right)\right\rangle , & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle {\mathcal{P}}{\mathcal{V}}=\unicode[STIX]{x1D701}-J\left[W,\frac{\unicode[STIX]{x1D6E9}^{\prime }}{\unicode[STIX]{x2202}_{Z}\overline{\unicode[STIX]{x1D6E9}}}\right]+\unicode[STIX]{x2202}_{Z}\left(\frac{\unicode[STIX]{x1D6E9}^{\prime }}{\unicode[STIX]{x2202}_{Z}\overline{\unicode[STIX]{x1D6E9}}}\right), & \displaystyle\end{eqnarray}$$
where
$J[f,g]=\unicode[STIX]{x2202}_{x}f\unicode[STIX]{x2202}_{y}g-\unicode[STIX]{x2202}_{y}f\unicode[STIX]{x2202}_{x}g$
. As a consequence of strong vertical motion, with
$|W|\sim |\boldsymbol{u}_{\bot }|$
,
${\mathcal{E}}$
and
${\mathcal{P}}{\mathcal{V}}$
are not solely functionals of the geostrophic streamfunction
$\unicode[STIX]{x1D6F9}$
: the second and third terms in (3.2) represent ageostrophic baroclinic contributions to
${\mathcal{P}}{\mathcal{V}}$
. Sole functional dependence on
$\unicode[STIX]{x1D6F9}$
may be recovered from the NH-QGE system only in the limit of strong stratification where
$\unicode[STIX]{x2202}_{Z}\overline{\unicode[STIX]{x1D6E9}}\rightarrow \unicode[STIX]{x1D716}^{-1}$
,
$\unicode[STIX]{x2202}_{Z}\rightarrow \unicode[STIX]{x1D716}^{-1}$
and
$W\rightarrow \unicode[STIX]{x1D716}$
. Here, hydrostatic balance in (2.2) implies
$\unicode[STIX]{x1D6E9}^{\prime }\rightarrow (Pr/\widetilde{Ra})\unicode[STIX]{x2202}_{Z}\unicode[STIX]{x1D6F9}$
and the classical hydrostatic quasi-geostrophic equation is recovered. In this case, the volume-averaged potential enstrophy
$\langle \overline{{\mathcal{P}}{\mathcal{V}}^{2}}\rangle$
becomes the second conserved quantity required to guarantee a dual cascade.
As already noted, however, the dynamics within the barotropic subspace provides an alternative pathway for an inverse cascade. Depth-averaging of the vertical vorticity equation (2.1) gives the barotropic vorticity equation (BVE),
$$\begin{eqnarray}\unicode[STIX]{x2202}_{t}\langle \unicode[STIX]{x1D701}\rangle +J[\langle \unicode[STIX]{x1D6F9}\rangle ,\langle \unicode[STIX]{x1D701}\rangle ]=-\langle J[\unicode[STIX]{x1D6F9}^{\prime },\unicode[STIX]{x1D701}^{\prime }]\rangle +\unicode[STIX]{x1D6FB}_{\bot }^{2}\langle \unicode[STIX]{x1D701}\rangle ,\end{eqnarray}$$
indicating that the material growth of
$\langle \unicode[STIX]{x1D701}\rangle$
depends on the net balance between the two terms on the right-hand side, i.e. between the baroclinic forcing and the barotropic viscous dissipation. In the absence of forcing and damping, the conserved quantities are the volume-averaged barotropic energy and enstrophy,
$$\begin{eqnarray}{\mathcal{E}}_{bt}={\textstyle \frac{1}{2}}\langle \overline{|\unicode[STIX]{x1D735}_{\bot }\langle \unicode[STIX]{x1D6F9}\rangle |^{2}}\rangle ,\quad {\mathcal{Z}}_{bt}=\overline{\langle \unicode[STIX]{x1D701}\rangle ^{2}}=\langle \overline{(\unicode[STIX]{x1D6FB}_{\bot }^{2}\langle \unicode[STIX]{x1D6F9}\rangle )^{2}}\rangle .\end{eqnarray}$$
Both are sole functionals of
$\langle \unicode[STIX]{x1D6F9}\rangle$
, suggesting a dual cascade. For all aspect ratios, we find that the evolution of the volume-averaged kinetic energy is similar to the isotropic case. Figure 4(a) illustrates the decomposition of the right-hand side of the BVE (3.3) as a function of time for
$\unicode[STIX]{x1D6E4}=1$
and
$\unicode[STIX]{x1D6E4}=6$
. The figure clearly demonstrates that in both cases the convective forcing
$\langle J[\unicode[STIX]{x1D6F9}^{\prime },\unicode[STIX]{x1D701}^{\prime }]\rangle$
and damping
$\unicode[STIX]{x1D6FB}_{\bot }^{2}\langle \unicode[STIX]{x1D701}\rangle$
are saturated but unbalanced, resulting in unbounded growth of the large-scale barotropic mode. However, despite this similarity, notable distinctions exist in the energetics. Figure 4(b) shows the corresponding kinetic energy spectra as a function of the renormalized horizontal wavenumber
$k/k_{Lx}$
, where
$k=\sqrt{k_{x}^{2}+k_{y}^{2}}$
and
$k_{Lx}=2\unicode[STIX]{x03C0}/L_{x}$
are associated with the box dimension
$L_{x}$
. Here,
$$\begin{eqnarray}E(k)=\int _{0}^{2\unicode[STIX]{x03C0}}E(\boldsymbol{k})k\,\text{d}\unicode[STIX]{x1D719}_{k},\end{eqnarray}$$
where in polar representation
$k=|\boldsymbol{k}|$
and
$\tan \unicode[STIX]{x1D719}_{k}=k_{y}/k_{x}$
. The spectra have also been decomposed into barotropic (bt) and baroclinic (bc) components. For both
$\unicode[STIX]{x1D6E4}=1$
and
$\unicode[STIX]{x1D6E4}=6$
(and indeed all intermediate cases), the barotropic signature (solid curves) exhibits a steep power law with
$E_{bt}(k)\sim k^{-3}$
, while the baroclinic signature (dashed curves) gains dominance at higher wavenumber and exhibits a shallower instantaneous power
$E_{bc}(k)\sim k^{-5/3}$
(Rubio et al.
Reference Rubio, Julien, Knobloch and Weiss2014). The total kinetic energy spectrum
$E_{bt}+E_{bc}$
exhibits a steep to shallow transition in the power law exponents, a result reminiscent of the Nastrom–Gage spectrum observed in atmospheric and oceanic measurements (Nastrom & Gage Reference Nastrom and Gage1985). In the present RRBC case, however, the
$k^{-3}$
barotropic spectrum is a consequence of the large-scale condensate (Smith & Waleffe Reference Smith and Waleffe1999). For
$\unicode[STIX]{x1D6E4}=1$
(grey curve), we observe that the most energetic barotropic scale is
$k/k_{Lx}=1$
, indicating that the large-scale condensate (the vortex dipole) has reached the box scale, the largest scale possible. For anisotropic aspect ratios, the most energetic scale is that associated with the unidirectional jet occurring at a scale intermediate to the box dimensions, i.e.
$k_{Ly}<k<k_{Lx}$
. For
$\unicode[STIX]{x1D6E4}=6$
, where four jets are observed, this occurs at
$k/k_{Lx}=4/6$
(see the solid line in figure 4
b), with a steep decline in power from its peak to the largest box scale
$2\unicode[STIX]{x03C0}/k_{Ly}$
. These larger scales are associated with the weaker meandering transverse jet (see figure 3
b,f). For
$2\leqslant \unicode[STIX]{x1D6E4}\leqslant 6$
, we observe maximal power in the barotropic energy spectra at
$1/2\leqslant k/k_{Lx}\leqslant 4/6$
, or, equivalently,
$3/2\leqslant L/L_{x}\leqslant 2$
.
The flow of energy associated with the generation of large-scale dynamics in the barotropic subspace can be determined from the nonlinear advection term in (3.3). Accordingly, we detail how power is transferred to the horizontal wavenumber
$\boldsymbol{k}$
through triadic interactions involving wavenumbers
$\boldsymbol{p}$
and
$\boldsymbol{q}$
such that
$\boldsymbol{p}+\boldsymbol{q}+\boldsymbol{k}=\mathbf{0}$
. We define the transfer functions
$$\begin{eqnarray}\displaystyle & \displaystyle T_{\boldsymbol{k}\boldsymbol{p}\boldsymbol{q}}=b_{\boldsymbol{p}\boldsymbol{q}}\text{Re}[\langle \hat{\unicode[STIX]{x1D6F9}}_{\boldsymbol{k}}\rangle \langle \hat{\unicode[STIX]{x1D6F9}}_{\boldsymbol{p}}\rangle \langle \hat{\unicode[STIX]{x1D6F9}}_{\boldsymbol{q}}\rangle ]\unicode[STIX]{x1D6FF}_{\boldsymbol{p}+\boldsymbol{q}+\boldsymbol{k},\mathbf{0}}, & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle b_{\boldsymbol{p}\boldsymbol{q}}=b_{\boldsymbol{q}\boldsymbol{p}}={\textstyle \frac{1}{2}}(p^{2}-q^{2})(p_{x}q_{y}-p_{y}q_{x}), & \displaystyle\end{eqnarray}$$
where
$\text{Re}$
denotes the real part and
$\unicode[STIX]{x1D6FF}_{\boldsymbol{p}+\boldsymbol{q}+\boldsymbol{k},\mathbf{0}}$
is the Kronecker delta function. Anisotropy is handled by replicating a
$1:\unicode[STIX]{x1D6E4}$
barotropic field
$\unicode[STIX]{x1D6E4}$
times to form a periodic square of size
$(\unicode[STIX]{x1D6E4}L_{x})^{2}$
. Owing to this periodic extension, the power map of a 2D Fourier transform of the barotropic field contains the non-integer wavenumber array
$(k_{x},k_{y})=(2\unicode[STIX]{x03C0}/L_{x})(i/\unicode[STIX]{x1D6E4},j/\unicode[STIX]{x1D6E4})$
, with
$i,j=0,1,\ldots ,N_{x}\unicode[STIX]{x1D6E4}$
. The array is sparse due to zero row entries corresponding to wavenumbers
$k_{x}$
with non-integer values
$i/\unicode[STIX]{x1D6E4}$
that do not fulfil the periodicity of the
$(0,L_{x})$
domain. To avoid the impact of these zero entries when averaging over annular or spherical shells in wavenumber space, we invoke coarse-graining by averaging over blocks of
$\unicode[STIX]{x1D6E4}^{2}$
entries: the block associated with integer wavelengths
$(2\unicode[STIX]{x03C0}/L_{x})(p,q)$
is indexed by
$$\begin{eqnarray}\displaystyle i_{x}=p-\frac{i-1}{\unicode[STIX]{x1D6E4}},\quad i=1,\unicode[STIX]{x1D6E4};\quad i_{y}=q-\frac{j-1}{\unicode[STIX]{x1D6E4}},\quad j=1,\unicode[STIX]{x1D6E4}. & & \displaystyle\end{eqnarray}$$
The result of coarse-graining is a power map array that contains integer wavenumbers
$(k_{x},k_{y})=(2\unicode[STIX]{x03C0}/L_{x})(p,q)$
, with
$p,q=0,1,\ldots ,N_{x}$
. The coarse-grained transfer map
$$\begin{eqnarray}\displaystyle T_{kp}=\int k\,\text{d}\unicode[STIX]{x1D719}_{k}\int p\,\text{d}\unicode[STIX]{x1D719}_{p}\mathop{\sum }_{\boldsymbol{q}}T_{\boldsymbol{k}\boldsymbol{p}\boldsymbol{q}} & & \displaystyle\end{eqnarray}$$
details the transfer of energy from wavenumber
$p$
to
$k$
. Figure 5 illustrates three cases:
$\unicode[STIX]{x1D6E4}=1$
and the coarse-grained transfer maps for
$\unicode[STIX]{x1D6E4}=3$
and
$6$
. The power signature in the super- and sub-off-diagonal lines in all barotropic self-interaction
$T_{kp}$
maps indicates the existence of a forward or direct cascade, i.e. direct transfer of spectral power from low to high
$k$
at constant wavenumber
$p$
. The non-local inverse cascade occurs for
$p\gg k$
and corresponds to the direct transfer of power from the high
$p$
wavenumbers to
$k\approx 1$
. Likewise, when
$p\ll k$
, energy is extracted from high
$k$
wavenumbers and transferred to
$p\approx 1$
. It should be recalled that for
$\unicode[STIX]{x1D6E4}>1$
, the wavenumber
$(p,k)=(1,1)$
associated with the large-scale structure is associated with subharmonic wavenumbers identified in (3.8) with
$p,q=1$
.
Spectral transfer maps of barotropic self-interaction for aspect ratios of 1 : 1 (a), 1 : 3 (b) and 1 : 6 (c), showing how energy is transferred from wavenumbers
$p$
to wavenumbers
$k$
. The results in (b,c) have been coarse-grained to manage the non-integer wavenumbers arising from the anisotropy of the domain.
3.3 Potential vorticity
From the full NH-QGE system, the barotropic and baroclinic potential vorticities are given by
$$\begin{eqnarray}{\mathcal{P}}{\mathcal{V}}_{bt}=\left\langle \unicode[STIX]{x1D701}\right\rangle -\left\langle J\left[W,\frac{\unicode[STIX]{x1D6E9}^{\prime }}{\unicode[STIX]{x2202}_{Z}\overline{\unicode[STIX]{x1D6E9}}}\right]\right\rangle ,\quad {\mathcal{P}}{\mathcal{V}}_{bc}=\unicode[STIX]{x1D701}^{\prime }-J\left[W,\frac{\unicode[STIX]{x1D6E9}^{\prime }}{\unicode[STIX]{x2202}_{Z}\overline{\unicode[STIX]{x1D6E9}}}\right]^{\prime }+\unicode[STIX]{x2202}_{Z}\left(\frac{\unicode[STIX]{x1D6E9}^{\prime }}{\unicode[STIX]{x2202}_{Z}\overline{\unicode[STIX]{x1D6E9}}}\right).\end{eqnarray}$$
Details of the characteristics of these two quantities are illustrated in figure 6, with results for aspect ratios of 1 : 2 and 1 : 5. A snapshot of
${\mathcal{P}}{\mathcal{V}}_{bt}$
as a function of
$y$
indicates that it is dominated by the barotropic vorticity
$\langle \unicode[STIX]{x1D701}\rangle$
. This signal oscillates in
$y$
, reflecting the spontaneous generation of a Kolmogorov-like flow. Consistent with the inverse energy cascade, we find that the magnitude of this signal is unbounded in time. In contrast, the baroclinic contribution to
${\mathcal{P}}{\mathcal{V}}_{bt}$
remains bounded and fluctuates randomly about zero without coherence. This result suggests that the barotropic dynamics is essentially linear. In contrast,
${\mathcal{P}}{\mathcal{V}}_{bc}$
and its component terms, which are of roughly equal magnitude, vary rapidly about zero mean and saturate in time (not shown). The finding that the dissipation and the forcing are both weak suggests that the barotropic manifold is amenable to the application of equilibrium statistical mechanics (Bouchet & Venaille Reference Bouchet and Venaille2012).
The barotropic potential vorticity
${\mathcal{P}}{\mathcal{V}}_{bt}$
and its decomposition for aspect ratios of 1 : 2 (a) and 1 : 5 (b). The
$\langle \unicode[STIX]{x1D701}\rangle$
term dominates, suggesting the linear relation
${\mathcal{P}}{\mathcal{V}}_{bt}\approx \langle \unicode[STIX]{x1D701}\rangle$
.
4 Conclusion
In this paper, we have confirmed, following earlier work (Julien et al.
Reference Julien, Rubio, Grooms and Knobloch2012; Rubio et al.
Reference Rubio, Julien, Knobloch and Weiss2014), that geostrophic turbulence is unstable to the formation of large-scale vortices, and investigated the properties of this state when the doubly periodic domain we use in the horizontal changes from square to rectangular. We have shown that with increasing domain anisotropy, the large-scale vortex structure is replaced by a shear flow of Kolmogorov type parallel to the shorter side and superposed on the turbulent state. The flow has jet-like structure with a well-defined characteristic scale comparable to the short box scale that forms spontaneously. It is neither externally imposed nor the result of imposed external stirring (Smith & Yakhot Reference Smith and Yakhot1994; Bouchet & Simonnet Reference Bouchet and Simonnet2009; Frishman et al.
Reference Frishman, Laurie and Falkovich2017) – here, the flow is maintained by 3D turbulent fluctuations that are determined self-consistently. The jets all have a well-defined mean separation that depends on the aspect ratio but undergo meander that may be intrinsic or driven by fluctuations. Evidently, jet formation does not require the presence of a
$\unicode[STIX]{x1D6FD}$
term.
These predictions resemble qualitatively the predictions from both numerical simulations of the damped 2D NSEs driven by imposed white noise (Smith & Yakhot Reference Smith and Yakhot1994; Bouchet & Simonnet Reference Bouchet and Simonnet2009; Frishman et al.
Reference Frishman, Laurie and Falkovich2017) and those from equilibrium statistical mechanics for 2D flows (Bouchet & Venaille Reference Bouchet and Venaille2012). Since the noise is imposed in the simulations, it cannot respond to the box-scale structures that evolve. This is also the case in the equilibrium statistical description, where noise must be assumed to be present to drive the system to equilibrium but the equilibrium reached is independent of the noise. The equilibrium statistical mechanics approach is based on maximizing the entropy, defined in terms of the vertical vorticity, subject to constraints derived from the inviscid equations of motion. The variational problem leads to a monotonically increasing relation between the most probable values of the vorticity
$\bar{\unicode[STIX]{x1D701}}$
and the streamfunction
$\bar{\unicode[STIX]{x1D6F9}}$
. We assume that these values are those realized by the flow, and hence drop the overbars in the discussion that follows. Whether the predictions of this type of theory are relevant to a 3D forced dissipative flow such as geostrophic turbulence remains a question, however. Figure 7 provides evidence that a
$\unicode[STIX]{x1D701}(\unicode[STIX]{x1D6F9})$
relation in fact exists for both
$\unicode[STIX]{x1D6E4}=1$
and
$\unicode[STIX]{x1D6E4}\neq 1$
, and is qualitatively similar to that found in 2D hydrodynamics. Indeed, if we suppose, following Bouchet & Simonnet (Reference Bouchet and Simonnet2009) and Bouchet & Venaille (Reference Bouchet and Venaille2012), a relation of the form
$\unicode[STIX]{x1D701}=a_{2}\unicode[STIX]{x1D6F9}+a_{4}\unicode[STIX]{x1D6F9}^{3}+\cdots \!,$
we conclude from figure 7(a) that in both cases
$a_{2}<0$
,
$a_{4}<0$
. The non-zero value of
$a_{4}$
is important in determining the nature of the transition from the box-scale dipole present when
$\unicode[STIX]{x1D6E4}=1$
to the jet-like flow parallel to the shorter side present when
$\unicode[STIX]{x1D6E4}\gtrsim 1$
, a transition whose presence here accords with both simulations of the damped 2D NSEs driven by noise (Bouchet & Simonnet Reference Bouchet and Simonnet2009; Frishman et al.
Reference Frishman, Laurie and Falkovich2017) and the prediction from equilibrium statistical mechanics (Bouchet & Simonnet Reference Bouchet and Simonnet2009; Bouchet & Venaille Reference Bouchet and Venaille2012). The latter predicts the presence of a phase transition from the dipole state to a parallel shear flow already at small values of the elongation
$\unicode[STIX]{x1D6E4}-1$
(Bouchet & Simonnet Reference Bouchet and Simonnet2009; Bouchet & Venaille Reference Bouchet and Venaille2012), a result qualitatively similar to both the 2D simulations and our self-consistent 3D simulations in figure 2 (second row:
$\unicode[STIX]{x1D6E4}=1.1$
).
Scatter plots of
$\langle \unicode[STIX]{x1D701}\rangle$
versus
$\langle \unicode[STIX]{x1D6F9}\rangle$
for aspect ratios of 1 : 1 (a) and 1 : 5 (b), suggesting a functional relation (black curve) obtained by pointwise averaging. The plots indicate that
$a_{2}a_{4}>0$
.
While the correspondence between our numerical results and equilibrium statistical mechanics theory is encouraging, it is also evident that our system is not an equilibrium system – it is a forced dissipative system. In this case, as already discussed, we expect to see an inverse energy cascade towards large scales and a direct entropy cascade to small scales. However, if both the forcing and the dissipation are appropriately weak on the large scales of interest, which we believe to be the case, the resulting non-equilibrium states that are observed are nevertheless expected to be close to the equilibrium states identified in the equilibrium statistical mechanics approach, a prediction corroborated in DNS of a stochastically driven vorticity equation in 2D (Smith & Yakhot Reference Smith and Yakhot1994; Bouchet & Simonnet Reference Bouchet and Simonnet2009; Frishman et al. Reference Frishman, Laurie and Falkovich2017). In particular, the transition from a dipole flow to a parallel shear flow with increasing elongation, predicted by the statistical mechanics approach, persists into the non-equilibrium regime (see figure 23 of Bouchet & Venaille Reference Bouchet and Venaille2012).
The above theory (Bouchet & Venaille Reference Bouchet and Venaille2012) has been developed for 2D flows. Our system is, in contrast, fully 3D, at least on small scales. However, the large-scale barotropic mode studied here obeys 2D dynamics driven by stochastic baroclinic forcing, and therefore shares many of the properties of these flows discussed by Bouchet & Venaille (Reference Bouchet and Venaille2012) and Frishman et al. (Reference Frishman, Laurie and Falkovich2017). In particular, we expect that our flows should also undergo a transition with increasing anisotropy from a dipole flow to a parallel shear flow, a prediction confirmed in our simulations. We emphasize that these are performed on an asymptotically reduced model valid in the limit
$Ro\rightarrow 0$
(Julien et al.
Reference Julien, Rubio, Grooms and Knobloch2012; Rubio et al.
Reference Rubio, Julien, Knobloch and Weiss2014). Given that the predictions of this model have been confirmed in subsequent simulations of the primitive equations at
$E=10^{-7}$
(Stellmach et al.
Reference Stellmach, Lischper, Julien, Vasil, Cheng, Ribeiro, King and Aurnou2014),
$E=10^{-6}$
(Favier et al.
Reference Favier, Silvers and Proctor2014) and
$E=5\times 10^{-6}$
(Guervilly et al.
Reference Guervilly, Hughes and Jones2014), we expect that our findings will also apply to the NSEs at sufficiently low Ekman numbers. That this is in fact the case is demonstrated by Guervilly & Hughes (Reference Guervilly and Hughes2017), who show that when
$E=10^{-5}$
, jets also appear when the elongation exceeds 10 %. The characteristic horizontal scales of the multi-jet states for larger elongations are also comparable. However, there are some differences too, in that Guervilly and Hughes also observe quite long-lived cyclonic vortices within their jets, whereas our system exhibits both cyclonic and anticyclonic vortices in equal numbers. This is a consequence of the reflection symmetries
$x\rightarrow -x$
,
$y\rightarrow -y$
of the NH-QGE system that are present at leading order in the limit
$Ro\rightarrow 0$
(
$E\rightarrow 0$
).
The above discussion indicates that the formation of large-scale structures, be they box-scale vortices or jets, is a very robust phenomenon, independent of the details of the fluctuations driving the system, independent of the specific system studied and even independent of its dimensionality, provided only that the flow is strongly anisotropic. Thus, large-scale vortices are also present in 3D non-rotating systems, provided that the fluid layer is sufficiently thin, thereby forcing the turbulent flow to be anisotropic on large scales (Xia et al. Reference Xia, Punzmann, Falkovich and Shats2008; Xia, Shats & Falkovich Reference Xia, Shats and Falkovich2009). We expect that in the presence of doubly periodic boundary conditions, this system will also undergo a transition with increasing domain anisotropy from the large-scale vortex state observed in square domains to a parallel shear flow, although this is of course difficult to confirm in experiments carried out in bounded domains. In thicker layers, the large-scale vortex reduces vertical motion, rendering the system susceptible to upscale energy cascade (Xia et al. Reference Xia, Byrne, Falkovich and Shats2011) and reinforcing the vortex. In this case, too, we expect a transition to large-scale shear flow with increased anisotropy much as occurs in convectively driven turbulence (Goluskin et al. Reference Goluskin, Johnston, Flierl and Spiegel2014). Both conjectures are amenable to confirmation by DNS.
Acknowledgements
This work was supported in part by the National Science Foundation under grants DMS-1317666, EAR-1620649 (K.J.) and DMS-1317596 (E.K.), and the NASA Earth and Space Science Fellowship Program (M.P.). The authors are grateful to C. Guervilly and D. W. Hughes for sharing their results prior to publication.
