Ordering of time scales predicts applicability of quasilinear theory in unstable flows

Abstract

We discuss the applicability of quasilinear-type approximations for a turbulent system with a large range of spatial and temporal scales. We consider a paradigm fluid system of rotating convection with vertical and horizontal temperature gradients. In particular, the interaction of rotation with the horizontal temperature gradient drives a ‘thermal wind’ shear flow whose strength is controlled by the horizontal temperature gradient. Varying this parameter therefore systematically alters the ordering of the shearing time scale, the convective time scale and the correlation time scale. We demonstrate that quasilinear-type approximations work well when the shearing time scale or the correlation time scale is sufficiently short. In all cases, the generalised quasilinear approximation systematically outperforms the quasilinear approximation. We discuss the consequences for statistical theories of turbulence interacting with mean gradients.

1. Introduction

Developing a description of nonlinear systems on a wide range of temporal and spatial scales is a significant problem for nonlinear physics. One such canonical example is turbulence (in both fluids and plasmas) where strong nonlinearities lead to interactions on many spatial and temporal scales. In its theoretically simplest manifestation of homogeneous and isotropic turbulence (HIT) the goal is to develop a theory that can describe phenomena as diverse as energy transfer via cascades, development of structures and spatial and temporal intermittency. In this paper, we are concerned with geophysical and astrophysical flows that interact with instabilities, often leading to anisotropic, inhomogeneous states and significant mean flows (Marston & Tobias 2023). At first sight, this problem is less tractable than HIT, because in that case the assumption of homogeneity and isotropy enables the specification of the form of the tensorial representation of transport coefficients. However, the interaction of fluctuating velocities with mean flows may justify approximations that can lead to successful description of the flows.

The approximations to represent anisotropic flows usually involve the prioritisation of certain nonlinear interactions over others (and are formally exact when certain non-dimensional parameters are asymptotically small). These approximations are not only of interest because they give information about which interactions are key in any given physical situation; they may lead to the development of accurate closures for statistical theories describing fluids and plasmas (Farrell & Ioannou 2003; Marston 2011; Marston, Qi & Tobias 2019).

The best-known example of such an approximation is the quasilinear (QL) approximation, which dates back to Malkus (1954). For a historical perspective/description of the QL approximation, see § 2.1 of Marston & Tobias (2023). Malkus proposed that turbulence driven by an instability saturates by modifying the mean state so that it is marginally stable to perturbations (on average). In this framework, the primary interaction is that of the turbulent fluctuations with themselves to modify the mean (back to marginality); the interactions of fluctuations to produce other fluctuations (and a cascade leading to dissipation) are subdominant. Herring provided an early illustration in the context of convection (Herring 1963). However, this idea quite often does not work well, although it can be shown to be formally correct in certain limits (Plumb 1977). In particular, the QL approximation works well when there is a strong separation of time scales, for example in the case of rapid rotation (Scott & Dritschel 2012), strong stratification (Chini, Malecha & Dreeben 2014) or for the cases of first-order smoothing in the mean-field electrodynamics (see e.g. Krause & Raedler 1980; Moffatt & Dormy 2019; Tobias 2021).

Recently, the QL approximation has been extended to consider more sets (triads) of nonlinear, non-local interactions to obtain the generalised quasilinear (GQL) approximation. For details of this approximation, see § 2.1. Briefly, the modes are spatially separated into low and high spectral wavenumbers – below and above a cutoff $\varLambda$ – representing large and small horizontal spatial scales, respectively. Certain triad interactions are prioritised over others. Generalised quasilinear theory has been shown to yield significant improvements over QL for even $\varLambda = 1$ (see e.g. Child et al. 2016; Hernández, Yang & Hwang 2022; Oishi & Baxter 2023).

In any driven, dissipative, turbulent system, it is difficult to predict and control the important time scales in the flow since they emerge as a result of the dynamics; hence, a separation of time scales that leads to effectiveness of QL-type theories (whether QL or GQL) can be difficult to enforce a priori. However, a theory may be developed based on an a posteriori calculation of time scales. In this paper, we propose a system where the ratio of time scales can be modified by smoothly changing an input parameter: the case of convection in the presence of a thermal wind. It is for this system we wish to investigate the utility of QL theories. The paper is organised as follows: in § 2 we introduce QL theories and describe a theory for their applicability. In § 3 we formulate the model and we discuss the results from direct numerical simulation for the time scale ordering in § 4. In § 5 we describe how well the QL and GQL approximations perform – demonstrating that the ordering of the time scales does play a key role in the effectiveness of these approximations. Section 6 presents our conclusions.

2. Quasilinear theory

2.1. Generalised quasilinear approximation

The GQL approximation was introduced by Marston, Chini & Tobias (2016) to extend the applicability of QL theories and to enable the derivation of more useful statistical theories. Here, we briefly review the approximation.

We consider the evolution of $\boldsymbol {q}$, a vector of state variables that satisfies a nonlinear evolution equation with quadratic nonlinearities, given by

(2.1)\begin{equation} \partial_t\boldsymbol{q}={\mathcal{L}}[\boldsymbol{q}] +{\mathcal{N}}[\boldsymbol{q},\boldsymbol{q}]. \end{equation}

Here, ${\mathcal {L}}$ is a linear operator and ${\mathcal {N}}$ is the nonlinear operator acting on the state variables. Each physical quantity $q=q_{l}+q_{h}$ can be split into ‘low’ modes (‘$l$’) and ‘high’ modes (‘$h$’), where the separation is achieved via spectral filtering. For a two-dimensional spatial separation in $x$ and $y$ we define the low and high modes via the set of wavevectors given by $\boldsymbol {k}\equiv (k_x,k_y)$

(2.2a)\begin{gather} \textbf{L} \equiv \{{ (k_x,k_y) : |k_x|\leq\varLambda\ \mbox{and}\ |k_y|\leq\varLambda }\}, \end{gather}

(2.2b)\begin{gather}\textbf{H} \equiv \{{(k_x,k_y) : |k_x|>\varLambda\ \mbox{or}\ |k_y|>\varLambda}\}. \end{gather}

The filter is therefore a square in wavenumber space with a size prescribed by the cutoff $\varLambda$ (e.g. see figure 1 showing a choice of $\varLambda =5$). If $\varLambda =0$ then the filter is equivalent to taking a formal mean and GQL reduces to QL. The ‘low’ and ‘high’ variables can then be written as partial Fourier series given by

(2.3a)\begin{gather} q_{\boldsymbol{l}}(x,y,z)=\sum_{\boldsymbol{k}\in\boldsymbol{\mathsf{L}} } {\tilde{q}_{l}}({\boldsymbol{k},z}) \exp\left[{2{\rm \pi} {\rm i}\left({{\frac{k_xx}{L_x}} +{\frac{k_yy}{L_y}} }\right)}\right] \end{gather}

(2.3b)\begin{gather}q_{h}(x,y,z) = \sum_{\boldsymbol{k}\in\boldsymbol{\mathsf{H}} } {\tilde{q}_{h}}({\boldsymbol{k},z}) \exp\left[{2{\rm \pi} {\rm i}\left({{\frac{k_xx}{L_x}} +{\frac{k_yy}{L_y}} }\right)}\right], \end{gather}

where the tilde denotes the coefficients in Fourier space, $L_x$ and $L_y$ are the horizontal grid sizes. In this paper we do not distinguish between the two horizontal directions and we set the coefficient $\varLambda$ to be the same in each direction, although in principle they could differ (for an example of this, see Tobias & Marston 2017).

Figure 1. Schematic of the GQL scheme in horizontal Fourier space, with discrete modes at integer $(k_x,k_y)$. In this example, the set of low modes ($\!\textbf {L}$) is within the orange box with cutoff $\varLambda =5$; and the high modes ($\!\textbf {H}$) include everything outside (cyan). The pink circle represents taking the spectral filter at $K=6$ (for example), during the calculation of kinetic energy transfer functions (5.1)–(5.3).

The corresponding equations of motion for the ‘low’ and ‘high’ modes then take the form

(2.4a)\begin{gather} \partial_t\boldsymbol{q}_{l} ={\mathcal{L}}[\boldsymbol{q}_{l}] +{\mathcal{N}}_{l}[\boldsymbol{q}_{l},\boldsymbol{q}_{l}] +{\mathcal{N}}_{l}[\boldsymbol{q}_{h},\boldsymbol{q}_{h}], \end{gather}

(2.4b)\begin{gather}\partial_t\boldsymbol{q}_{h} ={\mathcal{L}}[\boldsymbol{q}_{h}] +{\mathcal{N}}_{h}[\boldsymbol{q}_{l},\boldsymbol{q}_{h}] +{\mathcal{N}}_{h}[\boldsymbol{q}_{h},\boldsymbol{q}_{l}], \end{gather}

where the ‘$l$’ and ‘$h$’ subscripts on nonlinear operators symbolise low-pass and high-pass filtering. The nonlinear terms in these equations contain the essence of the GQL approximation; certain nonlinear interactions are removed in the evolution equations. In the evolution equations for the low modes the interaction of low modes with high modes to yield low modes $({\mathcal {N}}_{l}[\boldsymbol {q}_{l},\boldsymbol {q}_{h}])$ is removed, whilst in the evolution equation for the high modes two interactions are removed; the interaction of low modes with high modes to yield high modes $({\mathcal {N}}_{h}[\boldsymbol {q}_{l},\boldsymbol {q}_{h}])$ and the interaction of high modes with high modes to yield high modes $({\mathcal {N}}_{h}[\boldsymbol {q}_{h},\boldsymbol {q}_{h}])$. This is the so-called ‘pain-in-the-neck’ or eddy–eddy nonlinear (EENL) term that is ignored in all QL theories.

We note that the nonlinear interactions have been removed in pairs and that GQL is an example of Kraichnan's triad decimation in pairs (Kraichnan 1985). Hence, the equation sets continue to obey conservation laws for all global linear and quadratic invariants (Marston & Tobias 2023).

Two more properties of the system of equations are worthy of note at this point. The first is that, although there is no scope for direct nonlinear interactions of $q_{h}$ with $q_{h}$ to affect $\partial _t q_{h}$ (so-called high/high $\rightarrow$ high interactions), high modes can interact indirectly via the mediation of low modes (so-called eddy scattering off the mean flows). Meanwhile the low modes self-interact fully, equivalent to a model of lower resolution $\varLambda$, but also receive external transfers from high modes. Finally, we emphasise again that, when $\varLambda =0$, the GQL approximation reverts to the QL approximation, where the low modes are formal means and all the spatially fluctuating activity belongs to the high modes. Previous studies have shown that selection of $\varLambda > 0$ may strongly influence the physical realism of QL results, as it can better represent the natural scales of energy exchanges (Child et al. 2016; Marston et al. 2016; Hernández & Hwang 2020).

2.2. Theory for applicability of quasilinear systems: ordering of time scales

Our hypothesis is that the efficacy of QL-type approximations will vary as the ordering of time scales in the turbulence changes. A simple argument is to consider the nonlinear interactions in the momentum equation for the fluctuating turbulent velocity (for QL) or the momentum equation for the ‘high’ modes (for GQL). Here, we give the argument for QL and consider a velocity with mean flow $\bar {{\boldsymbol {u}}}$ and fluctuation ${\boldsymbol {u}}^\prime$. The equation for the fluctuating velocity then reads

(2.5)\begin{equation} \partial_t{{\boldsymbol{u}}^\prime} + {\boldsymbol{u}}^\prime \boldsymbol{\cdot} \boldsymbol{\nabla} \bar{{\boldsymbol{u}}} + \bar{{\boldsymbol{u}}} \boldsymbol{\cdot} \boldsymbol{\nabla} {\boldsymbol{u}}^\prime ={-}[{\boldsymbol{u}}^\prime \boldsymbol{\cdot} \boldsymbol{\nabla} {\boldsymbol{u}}^\prime - \overline{{\boldsymbol{u}}^\prime \boldsymbol{\cdot} \boldsymbol{\nabla} {\boldsymbol{u}}^\prime}] + \nu \nabla^2 {\boldsymbol{u}}^\prime + {\boldsymbol f}^\prime, \end{equation}

where overbar indicates a mean and primed variables are fluctuations. All linear driving and redistribution terms such as buoyancy and Coriolis forces are included in ${\boldsymbol f}^\prime$. Quasilinear theory can be an accurate description of the mean dynamics if the ‘pain-in-the-neck’ term in brackets on the right-hand side of (2.5) is formally smaller than some of the other terms in the equation. What is needed is to provide estimates for the amplitudes of these terms. We note that this is perforce a crude estimation of the importance of nonlinear terms since it is the projection of these terms back onto the mean (or the low modes) that really determines the role of these terms and whether they can be discarded. We shall test whether the ratio of the magnitude of these terms does in fact give a reliable indication of the applicability of QL theory here.

The magnitude of the terms in (2.5) can be calculated by estimating the scales of local spatial and temporal derivatives. These are characterised by a correlation length $\ell _{c}$ (spatial) and time $\tau _{c}$ (temporal) (we shall describe the procedure adopted for calculating these in the results section). It is also useful to define the convective turnover time based on root mean square (rms) velocity fluctuation,

(2.6)\begin{equation} \tau_{o}\approx \frac{\ell_{c}}{u_{rms}}, \end{equation}

and a shear time

(2.7)\begin{equation} \tau_{s} \approx \left(\frac{\mathrm{d}\bar{{\boldsymbol{u}}}}{\mathrm{d}z}\right)^{{-}1}. \end{equation}

The ratio of the EENL term to the shearing term is given by $S = \tau _{s}/\tau _{o}$. If $S$, sometimes termed the Strouhal number, is small then the EENL term is negligible compared with the shear term, and the ‘pain-in-the-neck’ term may safely be relegated. Another interesting ratio of time scales is that which compares the importance the EENL term with the local acceleration $\partial {\boldsymbol {u}}^\prime /\partial t$ for which the time scale is the correlation time. This ratio is measured by the Kubo number $K=\tau _{c}/\tau _{o}$. If $K$ is small then the turbulence rapidly decorrelates and so the EENL term may be discarded. For convective turbulence in the absence of rotation and shear, $K$ is often found to be order unity and so QL theories might not be applicable. However, the presence of rotation can lead to flows and structures that remain correlated on time scales longer than the turnover time of the eddy (sometimes termed coherent structures), whilst the presence of shear can lead to a shortening of the correlation time (as the eddies are sheared apart rapidly before having the opportunity to turn over).

In § 4 we calculate the time scales defined above for a range of parameters and thermal wind strengths, before testing the theory that the ratio of these time scales will provide a measure of how effective a QL (or GQL) theory is. A nice feature of the model we examine is that the time scales can be efficiently re-ordered by varying a single parameter, which is the relative strength of the horizontal temperature gradient.

3. Formulation of the model

The model is set up similarly to that considered by Hathaway & Somerville (1986) and Currie (2014). We consider Boussinesq convection in a local three-dimensional Cartesian domain $(x,y,z)$ rotating with angular velocity $\boldsymbol {\varOmega } = \varOmega \boldsymbol {\hat {\varOmega }}$ at an angle ${\rm \pi} /2-\phi$ to the direction of gravity (${\hat {\boldsymbol {z}}}$), as shown in figure 2. The Boussinesq equations are non-dimensionalised; lengths are non-dimensionalised with the layer depth $d$, whilst time is non-dimensionalised with the thermal diffusion time $\tau = d^2/\kappa$, given thermal diffusivity $\kappa$. The fluid is subjected to both vertical and horizontal temperature gradients, and temperatures are non-dimensionalised with the vertical temperature difference across the layer $\Delta T = \beta d$, where $\beta = |\partial T/\partial z|$ of the conducting state.

Figure 2. Computational domain. The background temperature is illustrated in colour in addition to the geometry of the constant shear thermal wind.

With this non-dimensionalisation we write the non-dimensional basic state temperature as

(3.1)\begin{equation} T_\text{B} = T_0 - z + {T_y} y, \end{equation}

where $(0 \leq x < L_x, 0 \leq y < L_y, 0 \leq z \leq 1)$. The latitudinal temperature gradient drives a thermal wind via the interaction with rotation. The basic state velocity is orientated in the $x$-direction and takes the form of a linear flow with constant shear (Hathaway & Somerville 1986; Currie 2014) that is given by

(3.2)\begin{equation} U_\text{B} ={-}\frac{{{Ra}} {T_y}}{{{Ta}}^{1/2} \sin \phi} \left( z - \frac{1}{2}\right). \end{equation}

Here, the non-dimensional parameters are the Rayleigh number ${{Ra}} \equiv g \alpha d^3 \Delta T/\kappa \nu$, and Taylor number ${{Ta}} \equiv 4 \varOmega ^2 d^4/\nu ^2$. In this paper we hold fixed both the Taylor number ${{Ta}} = 10^5$ and the other non-dimensional parameter, which is the Prandtl number, ${\textit {Pr}} \equiv \nu /\kappa = 1$. The horizontal dimensions are set as $L_x=5$, $L_y=5$. We set the angle $\phi ={\rm \pi} /4$ so the box is placed at midlatitude.

We solve for perturbations to the basic state, which satisfy

(3.3)\begin{gather} \frac{\mathrm{D} {\boldsymbol{u}}}{\mathrm{D} t} +{{Ta}} {\textit{Pr}} \boldsymbol{\hat{\varOmega}} \times {\boldsymbol{u}} ={-}{\textit{Pr}} \boldsymbol{\nabla} p + {{Ra}} {\textit{Pr}} \theta {\hat{\boldsymbol{z}}} + {\textit{Pr}} \nabla^2 {\boldsymbol{u}}, \end{gather}

(3.4)\begin{gather}\boldsymbol{\nabla}\boldsymbol{\cdot} {\boldsymbol{u}} = 0, \end{gather}

(3.5)\begin{gather}\frac{\mathrm{D}}{\mathrm{D} t} (T_\text{B} + \theta) = \nabla^2 \theta, \end{gather}

where the total convective derivative is given by $\mathrm {D}/\mathrm {D}t \equiv (\partial /\partial t + (U_\text {B} {\hat {\boldsymbol {x}}} + {\boldsymbol {u}}) \boldsymbol {\cdot } \boldsymbol {\nabla })$. The boundary conditions at $z=0,1$ are impenetrable, stress free (on the velocity perturbations to the thermal wind – so that the total stress is fixed on the boundaries) and isothermal. The system is periodic in the horizontal directions.

4. Direct numerical solutions and calculation of time scale ordering

We initially perform fully nonlinear three-dimensional direct numerical simulations (DNS) of (3.3)–(3.5) using the fully open source Dedalus pseudo-spectral code (Burns et al. 2020). We use Fourier bases in the horizontal ($x$ and $y$) directions, and a Chebyshev representation in $z$. We perform simulations for three different choices of latitudinal temperature gradient (${T_y} = 0$, $-0.5$ and $-2$) leading to different thermal wind strengths. We also consider DNS for five different choices of Rayleigh number to calculate the time scale ordering. Two of these Rayleigh number choices (${{Ra}} = 4\times 10^4$ and $2\times 10^5$) are further used to compare the fully nonlinear results with QL and GQL calculations, leading to six different dynamical regimes for that comparison. Resolutions are $(n_x,n_y,n_z) = (64,64,64)$ for ${ {Ra}} = 4\times 10^4$; $(64,64,128)$ for ${{Ra}} = 2\times 10^5$, ${T_y}=0$; and $(128,128,64)$ for ${{Ra}} = 2\times 10^5$, ${T_y}=-2$.

Figure 3 shows snapshots of the out-of-plane ($y$-component of the) vorticity in the $(x,z)$ plane for four fiducial parameter choices. These show the typical structure of the eddies in the convection simulations. For all the parameters considered, the solution settles down into a statistically steady state after an initial transient; for each case the temporal behaviour of the statistically steady state is chaotic. The convection takes the form of cellular structures and, unsurprisingly, as the vertical thermal driving (as measured by ${{Ra}}$) is increased the level of turbulence increases and the structures become stronger, more time dependent and stay coherent for less time (as will be shown later). Moreover, as the latitudinal thermal gradient is increased in magnitude, the thermal wind strengthens (as shown by (3.1)) and the solutions become more sheared, with a strongest horizontal vorticity.

Figure 3. Results of DNS. (a,b) The $y$ component of vorticity for ${T_y} = 0$: ${{Ra}} = 4\times 10^4$ (a), ${{Ra}} = 2\times 10^5$ (b); (c,d) $y$ component of vorticity for ${T_y} = -2$: ${{Ra}} = 4\times 10^4$ (c), ${{Ra}} = 2\times 10^5$ (d).

Figure 4(a–c) shows how the relevant time scales of the flow, $\tau _{s}$, $\tau _{o}$ (given by (2.6)–(2.7)) and $\tau _{c}$, change as a function of ${{Ra}}$ and ${T_y}$. The correlation length and correlation time are calculated using the procedure of Käpylä et al. (2006), which is briefly summarised here. The calculation of the correlation time uses timeseries of velocity components for each spatial location. We deduct the local temporal mean to obtain the fluctuations, form the temporal autocorrelation function of the fluctuating velocities at times $t$ and $t+\tau$ and normalise by the temporal variance. The local correlation time at point $(x,y,z)$ is the first root where this autocorrelation drops from unity to half its value (similar to calculations by Käpylä et al. 2006). The global correlation time $\tau _{c}$ is then given by the globally representative peak of the distribution of local correlation times.

Figure 4. Time scales of the DNS. (a–c) Overturning, correlation and shear time scales $\tau _{o}$, $\tau _{c}$, $\tau _{s}$ as a function of ${{Ra}}$ for ${T_y} = 0, -0.5, -2$ (left to right).

Spatial correlation scales in three principal directions are obtained more simply from the spatially integrated autocorrelation function of velocity fluctuations and again calculating where this drops by a factor of half. The correlation length is then given by the Pythagorean mean of these directed correlation lengths, i.e. $\ell _{c}=( x_{c}^2+y_{c}^2+z_{c}^2)^{1/2}$.

Robustly measuring all three characteristic time scales from DNS ideally requires each scale be resolved by a large number of temporal snapshots, and also requires a simulation prolonged enough to show many cycles.

In the case of no thermal wind (regular rotating convection), shown in figure 4(a), the extremely weak and intermittent shear leads to a shear time $\tau _{s}$ that is long compared with both the turnover time $\tau _{o}=\ell _{c}/u_{rms}$ and correlation time $\tau _{c}$. As the Rayleigh number increases, the stronger flows lead to a decrease of the shear time. For non-sheared convection, the correlation time of the convective structures is roughly of the same order as the turnover time, as is expected; cells stay coherent for the order of a turnover time. Again, as ${{Ra}}$ in increased both the correlation time and turnover time of the stronger convection decrease, remaining much smaller than the shear time scale. Hence, $S$ is large and $K\sim {{O}}(1)$ here; QL theory is not expected to be successful for any of these cases.

For the case with strong thermal wind ${T_y}=-2$, shown in figure 4(c) the situation has changed dramatically. In this case the shear is strong and the shear time scale is short. Moreover, interestingly for this case the correlation time of the convective eddies can now be shorter than their turnover time at high ${{Ra}}$ – the eddies are ripped apart by the strong shear before they can turn over once. Hence, we are in the regime where $\tau _{s} \ll \tau _{o}$ and $\tau _{c} < \tau _{o}$. Hence, for this case, the Strouhal number $S \ll 1$ and the Kubo number is small ($K<1$). Our prediction, therefore, is that QL-type theories will work better for this case with a strong thermal wind.

For the case with ${T_y}=-0.5$, all time scales are roughly of the same order for all values of the Rayleigh number; both the Strouhal number and the Kubo number are therefore of order unity and there is no compelling reason why QL or GQL should provide an accurate description of the statistics of the flow.

5. Effectiveness of quasilinear theories and nonlinear transfers

In order to assess the effectiveness of the QL approximations, we apply GQL approximations to calculations at two different values of the Rayleigh number for three different choices of the thermal wind parameter. We consider GQL at different wavenumber threshold cutoffs $\varLambda$; i.e. $k_x,k_y \le \varLambda$ are considered low modes whilst all others are considered high modes, as discussed in § 2.1. We present results for $\varLambda =0$, a choice that corresponds to QL evolution in both horizontal directions. Higher cutoffs set at $\varLambda =1,5,10$ yield GQL evolution in the horizontal.

In figure 5 we calculate the total kinetic energy as a function of ${{Ra}}$, ${T_y}$ and cutoff threshold $\varLambda$, and compare with results from the full DNS. The figure shows the average kinetic energy in the saturated state as a fraction of that obtained for DNS, for the various degrees of QL and GQL approximation (i.e. as $\varLambda$ varies). The first thing to note is that all QL approximations (for all the parameters) overestimate the kinetic energy in convection. This is unsurprising and is easily explained; the lack of the EENL term usually leads to underestimated dissipation (and hence overestimated energy) owing to the missing local cascades; it is these cascades that are important for carrying the energy that is input into the system at large or intermediate scales down to the dissipative scale where it can be decay. The absence of the local cascade makes this process significantly less efficient. For straight rotating convection (no thermal wind; ${T_y}=0$), all the GQL approximations perform poorly until $\varLambda =10$. As $\varLambda$ increases, the approximation does improve – and the kinetic energy more closely matches that found for the DNS solutions. This demonstrates that GQL performs better than QL ($\varLambda =0$). This agrees with previous studies that demonstrate how GQL can lead to improved performance over QL by including energy transfer between high modes via non-local (in spectral space) interactions (Child et al. 2016; Tobias & Marston 2017). This eddy-scattering effect is important, as it allows energy transfer between high modes via scattering off the low modes.

Figure 5. Results of QL/GQL simulations. Time-averaged kinetic energy (KE) in saturation normalised by DNS kinetic energy as a function of GQL cutoff $\varLambda$ for ${T_y}= -2, -0.5, 0$.

The fact that GQL does not produce reasonable results until $\varLambda =10$ shows that the rotating convection system is not very amenable to QL-type approximations – even if eddy scattering off the mean flow is included – and that the EENL plays an important role. This is in agreement with our arguments in the previous section, based on the ratio of time scales; there is not much scope for neglecting the EENL term a priori in favour of either of the advective terms in the equation for the fluctuating velocity.

As $|{T_y}|$ is increased, the QL approximations begin to perform better. In all cases, however, the QL approximation ($\varLambda =0$) significantly overestimates the kinetic energies. However, when the shear time is small compared with the overturning time, GQL (even with $\varLambda =1$, the least elaborate version of the theory) provides a good approximation to the full system – at least in representing global quantities such as the kinetic energy.

Figure 6 demonstrates how well the various levels of approximation reproduce the mean flow profiles (where the average is taken horizontally). Again, for straight convection, the means are poorly represented by the QL theories. However, as the separation of time scales becomes more pronounced and the EENL term becomes subdominant, GQL, even with $\varLambda =1$, provides an accurate representation of the saturated mean flow structures. It is precisely in these regimes that QL theories become effective, as predicted by the theory in § 2.1.

Figure 6. Results of QL/GQL simulations. Vertical ($z$) profiles of time and $(x,y)$-averaged $u$ velocity perturbations as a function of GQL cutoff $\varLambda$ for ${T_y}= -2, -0.5, 0$.

The QL approximation and its generalisations are examples of approximation via constrained triad decimation in pairs (Kraichnan 1985). The removal of certain interactions prevents certain routes for energy transfer. The effectiveness of the approximation presumably depends on how important these transfers were in the fully nonlinear DNS and how much damage removing them does to the dynamics and statistics of the system.

To reveal the preferred paths for transfer under the various degrees of approximation, we calculate energy transfer functions for the kinetic energy using methods outlined in Verma (2004), Alexakis, Mininni & Pouquet (2005, 2007), Lesur & Longaretti (2011), Favier, Silvers & Proctor (2014) and Currie & Tobias (2020). Briefly, we calculate the transfer in horizontal wavenumber at a fixed $z$ location, so that each physical variable such as $q=q(x,y,z)$ has a corresponding Fourier transform $\tilde {q}=\tilde {q}(k_x,k_y,z)$. We denote a circular shell in Fourier space (e.g. pink in figure 1) where wavenumbers are near absolute value $K$

(5.1)\begin{equation} \textbf{K}\equiv \{{ (k_x,k_y) : K<\sqrt{k_x^2+k_y^2}\leq K+1 }\}.\end{equation}

This leads to a $K$-filtered representation

(5.2)\begin{equation} q_{K}(x,y,z)\equiv \sum_{\boldsymbol{k}\in\boldsymbol{\mathsf{K}}} \tilde{q}(k_x,k_y,z) \exp\left[{ 2{\rm \pi} {\rm i}\left({ {\frac{k_x x}{L_x}}+{\frac{k_y y}{L_y}} }\right) }\right], \end{equation}

which forms an orthogonal basis. Given the $K$-filtered vector velocity $\boldsymbol {u}_{K}\equiv (u_{K},v_{K},w_{K})$, the $K$-component of kinetic energy is $E_{K}={\frac 12}\boldsymbol {u}_{K}\boldsymbol {\cdot }\boldsymbol {u}_{K}$. The equation for $\partial E_{K}/\partial t$ leads to the energy ‘transfer function’ from modes at scale $Q$ to modes at a different scale $K$

(5.3)\begin{equation} \mathcal{T}(Q,K)\equiv{-}\iiint{\mathrm{d}\kern0.7pt x\,\mathrm{d}y\,\mathrm{d}z\,\boldsymbol{u}_{K} \boldsymbol{\cdot} (\boldsymbol{u}\boldsymbol{\cdot}\boldsymbol{\nabla})\boldsymbol{u}_{Q}}. \end{equation}

Vertically integrated transfer functions should ideally be antisymmetric about the line $Q=K$. In practice, the transfers at low $K$ and low $Q$ are noisy.

Figure 7 shows some examples of these $\mathcal {T}(Q,K)$ maps (temporally averaged); here the colour scales are proportional to absolute values, but coloured to indicate the sign (blue for $\mathcal {T}<0$; red for $\mathcal {T}>0$). Figure 7 shows these transfer functions for two different choices of ${T_y}=0, -2$ (${{Ra}}$ is fixed at the higher value of ${{Ra}} = 2 \times 10^5$) and two different choices of $\varLambda =0, 1$ and for DNS. In all cases the QL version cuts off most of the channels for interaction and so no cascades are available for energy transfer. Utilising GQL provides channels for non-local interactions and energy is scattered among modes.

Figure 7. Spectral transfer function $\mathcal {T}(k, p)$ between two wavenumbers $(k, p)$. Panels (a–c) show ${T_y} = 0$; (d–f) show ${T_y} = -2$.

The transfer function figure clearly shows that when $|{T_y}|$ is small, almost all transfers occurring in triad interactions between fluctuating modes and the EENL term are important. This is most obvious from the presence of energy transfer just above and below the diagonal. The local cascade is therefore dominant for turbulent flows with no mean flows. Note that in the case ${T_y}=0$, QL saturates via interaction of turbulent fluctuations with the mean flow (figure 7a), even though this is not the dominant mechanism in the DNS (c); the system has to saturate by altering the means, rather than utilising a cascade to increase dissipation.

However, when ${T_y}=-2$ nearly all of the energy is transferred via the mean flow in DNS; energy is transferred between turbulent modes in the cascade because of the interaction with the mean. It is this interaction that is more readily captured via QL descriptions; with GQL with $\varLambda =1$ giving a nearly perfect representation of the transfers for this case. A similar improvement in the accuracy of QL was found for the case of Langmuir turbulence which also breaks horizontal isotropy (Skitka, Marston & Fox-Kemper 2020).

6. Conclusions

In this paper, we investigated how the effectiveness of the QL and GQL approximations in describing turbulent behaviour depends on the ordering of time scales in the flow. We derived and tested the hypothesis that the effectiveness of such approximations depends on the ordering of time scales in the turbulence. If the advective time of the turbulence is long compared with either a shearing time or a correlation time, then the QL-type approximations will perform better.

We test these ideas by examining a rotating convective system in the presence of a thermal wind. This system is carefully chosen as the ordering of the time scales is naturally selected by the amplitude of the latitudinal temperature gradient (for other parameters fixed); the strength of the thermal wind can be systematically varied so that the time scales are ordered differently for different parameter sets. Using full DNS, as well as QL and GQL simulations, we confirm our hypothesis; as the thermal wind is increased and the shear time scale becomes small compared with the overturning and correlation time scales then the accuracy of QL theories becomes better. Moreover, we demonstrate that GQL systematically improves on QL dynamical descriptions, owing to its ability to scatter energy between high modes via an interaction with low modes.

What we have shown here is that, by constructing a fluid system with a tuneable ordering of time scales, nonlinear interactions on a wide range of spatial and temporal scales can be understood and modelled via restricted modal equation sets. Although we have shown this for a turbulent fluid dynamic system, the conclusion applies to a wide variety of quadratically nonlinear systems. As the restricted modal equation set consists of a reduction of interactions via triad decimation in pairs, it conserves all global quadratic invariants in the limit of no driving and dissipation. We believe therefore that our approach is very general and that progress can be made on similar nonlinear problems both in fluid dynamics and plasmas with complicated spatial and temporal interactions.

Finally, we comment that QL (and GQL models) are useful as they enable the derivation of statistical models for turbulent flows. The absence of the EENL term enables the exact closure of the system at second order for statistical theories. Such (inhomogeneous and anisotropic) models for QL systems (sometimes known as second-order cumulant expansion (CE2) or stochastic structural stability theories (Farrell & Ioannou 2003; Marston 2011; Marston et al. 2019)) are effective in a number of important geophysical and astrophysical situations. Excitingly, the statistical theory for GQL (termed GCE2) has recently been introduced (Nivarti, Marston & Tobias 2023) and shown to improve over CE2 in a wide range of circumstances.