1. Introduction
When turbulence occurs in environmental flows, it is often affected by both mean-shear and stable stratification (Vallis Reference Vallis2006; Ferrari & Wunsch Reference Ferrari and Wunsch2009; Wyngaard Reference Wyngaard2010; Zorzetto, Bragg & Katul Reference Zorzetto, Bragg and Katul2018; Ayet et al. Reference Ayet, Katul, Bragg and Redelsperger2020), leading to sheared, stably stratified turbulence (SSST). Mean-shear (vertical gradient of horizontal flow) produces turbulence, with hairpin vortices, streaks and strong fluctuations in all three directions (Lee, Kim & Moin Reference Lee, Kim and Moin1990; Pope Reference Pope2000; Davidson Reference Davidson2004). Stable stratification suppresses fluctuations in the vertical direction, and if sufficiently strong generates a quasi-two-dimensional flow behaviour (Riley & Lelong Reference Riley and Lelong2000; Riley & Lindborg Reference Riley and Lindborg2012). Stable stratification also provides a restoring buoyancy force that enables the propagation of internal waves and the formation of quasi-horizontal ‘pancake’ vortical structures in the flow (Davidson Reference Davidson2004). Since mean-shear and stable stratification have competing effects, and produce different flow structures, the dynamics of SSST are rich and complex. Moreover, the flow is often ‘patchy’, with turbulent and non-turbulent regions interspersed, depending on the local competition between shear and buoyancy. Indeed, when the local shear in the flow is strong enough for the local Richardson number to be sufficiently small, the local flow may undergo Kelvin–Helmholtz instabilities that can evolve due to nonlinearity into turbulent motion (Riley & de Bruyn Kops Reference Riley and de Bruyn Kops2003).
Understanding and modelling SSST is an active area of research with many open questions. One such issue concerns understanding the mixing efficiencies in SSST and their parametric dependence (Peltier & Caulfield Reference Peltier and Caulfield2003; Portwood, de Bruyn Kops & Caulfield Reference Portwood, de Bruyn Kops and Caulfield2019), which are vital for predicting mixing in oceans (Jayne Reference Jayne2009; Gregg et al. Reference Gregg, D'Asaro, Riley and Kunze2018). Another vital area is to understand the properties of SSST across its range of dynamical scales, and the physical mechanisms that govern the fluxes of turbulent kinetic and potential energy between scales. Not only is this important for a basic understanding of the flow, it is also of crucial importance for developing large-eddy simulation (LES) models for SSST. A number of studies have considered the effect of stable stratification on the multiscale properties of turbulence using Fourier analysis and considering the average behaviour of the flow in terms of the energy spectrum, the mean buoyancy and mean interscale energy transfer terms (Riley & de Bruyn Kops Reference Riley and de Bruyn Kops2003; Lindborg Reference Lindborg2006; Almalkie & de Bruyn Kops Reference Almalkie and de Bruyn Kops2012). The study of Riley & de Bruyn Kops (Reference Riley and de Bruyn Kops2003) showed that the horizonal energy spectrum exhibits a 
$k_h^{-5/3}$ scaling for wavenumbers smaller than the Ozmidov wavenumber 
$k_O$ (corresponding to the wavenumber at which inertial and buoyancy forces are of the same order), where 
$k_h$ is the horizontal wavenumber, and the results indicated a downscale energy transfer of kinetic energy in the flow. This motivated (Lindborg Reference Lindborg2006), who confirmed that strongly stratified turbulent flows exhibit a downscale cascade of turbulent kinetic and potential energy on average, and developed phenomenological predictions similar in spirit to Kolmogorov's 1941 theory (Kolmogorov Reference Kolmogorov1941). The observation of a downscale energy cascade was contrary to predictions that had been made by Gage (Reference Gage1979) and Lilly (Reference Lilly1983) based on the assumption that in the limit of strong stratification, the flow should behave as two-dimensional turbulence. The basic reason why the predictions of Gage (Reference Gage1979) and Lilly (Reference Lilly1983) failed is that, as shown in Billant & Chomaz (Reference Billant and Chomaz2001), as the strength of the stratification increases and the vertical velocity of the flow is suppressed, the vertical length scale of the flow also reduces in such a way that the terms in the dynamical equations associated with vertical motion always remain 
$O(1)$, hence the flow never becomes two-dimensional.
Many questions remain, however, regarding the multiscale properties of SSST. For example, how do fluctuations of the SSST flow about its mean-field state behave? What are the mechanisms of the turbulent kinetic and potential energy transfers among scales, and to what extent are these transfers coupled to each other, and to fluctuations in the local buoyancy? How do the sub-grid scale terms in SSST contribute to the flow dynamics across scales? To address these and other questions, in this paper we will use a filtering-based approach (Liu, Meneveau & Katz Reference Liu, Meneveau and Katz1994; Eyink Reference Eyink1995a,Reference Eyinkb; Meneveau & Katz Reference Meneveau and Katz2000) wherein the velocity and density fields are considered for different filtering lengths 
$\ell$. The filtering-based approach has been applied previously to the study of energy and potential enstrophy cascades in the context of rotating, stratified flows (Aluie & Kurien Reference Aluie and Kurien2011). We apply the filtering approach to explore the turbulent kinetic and potential energies at different scales of SSST, and the various processes that determine their behaviour. Using a filtering approach allows us to consider spatially local couplings between the processes that control the flow energetics at different scales, and to consider fluctuations in the flow as well as the mean-field behaviour. Moreover, such an analysis of the multiscale properties of SSST can provide insights for developing LES models of SSST, since LES models are often developed in the physical space, rather than Fourier space, context. LES modelling of stratified turbulence is an underdeveloped area (Khani & Waite Reference Khani and Waite2015), and new insights into the dynamics of stratified turbulence across scales could aid in the development of appropriate sub-grid models.
The outline of the paper is as follows. In § 2, we introduce the equations governing the turbulent kinetic and potential energy at different scales, consider the behaviour of the mean-field state of the flow at different scales, and discuss the physical mechanisms governing the energy transfer among scales in SSST. In § 3, we describe the direct numerical simulations (DNS) dataset used in this study, and in § 4, we present and discuss results from the DNS for the quantities introduced in § 2 that govern the turbulent kinetic and potential energy at different scales. Finally, in § 5, we draw conclusions from the study and discuss areas for future investigation.
2. Theoretical considerations
2.1. Governing equations for scale-dependent energy fields
We consider the case where the mean velocity gradient 
$\gamma$ and the mean density gradient 
$\zeta$ (note that 
$\zeta <0$) are constant in space and time, such that the total fluid velocity vector can be written as 
$\boldsymbol {U}=z\gamma \boldsymbol {e}_x+\boldsymbol {u}$, where 
$\boldsymbol {u}$ is the fluctuating component of the velocity, and the total density can be written as 
$\rho =\rho _r +z\zeta +\rho '$, where 
$\rho _r$ is a constant reference density, and 
$\rho '$ is the fluctuating component of the density. While 
$\boldsymbol {e}_x$ is the unit vector in the direction of the mean velocity (
$\langle \boldsymbol {U}\rangle =z\gamma \boldsymbol {e}_x$), 
$\boldsymbol {e}_z$ is the unit vector pointing in the direction opposite to the gravitational acceleration 
$\boldsymbol {g}$ (i.e. the vertical direction), and 
$\boldsymbol {e}_y\equiv \boldsymbol {e}_z\times \boldsymbol {e}_x$.
The filtering operator to be used is defined for an arbitrary field 
$\boldsymbol {a}(\boldsymbol {x},t)$ as 
$\tilde {\boldsymbol {a}}(\boldsymbol {x},t)\equiv \int _{\mathbb {R}^3}\boldsymbol {a}(\boldsymbol {x}-\boldsymbol {y},t)\,\mathcal {G}_\ell (\boldsymbol {y})\,{\rm d}\boldsymbol {y}$, where 
$\mathcal {G}_\ell$ is an isotropic filtering kernel with length scale 
$\ell$. In order for the energy fields introduced below to be strictly non-negative, which is a physical requirement, a necessary and sufficient condition is that the filtering kernel 
$\mathcal {G}_\ell$ be non-negative (Vreman, Geurts & Kuerten Reference Vreman, Geurts and Kuerten1994). We will use the isotropic Gaussian filter 
$\mathcal {G}_\ell (\boldsymbol {y})\equiv (2{\rm \pi} \ell ^2)^{-3/2}\exp (-\|\boldsymbol {y}\|^2/2\ell ^2)$, which satisfies this condition. Further considerations regarding the implications of using different filtering kernels can be found in Sadek & Aluie (Reference Sadek and Aluie2018).
Assuming that 
$|\rho '|/\rho \ll 1$, the governing equation for 
$\boldsymbol {u}$ is then the Boussinesq–Navier– Stokes equation coupled with an advection–diffusion equation for 
$\rho '$ (Vallis Reference Vallis2006). The filtered versions of these equations are
$$\begin{gather} D_t\tilde{\boldsymbol{u}}={-}(1/\rho_r)\,\boldsymbol{\nabla}\tilde{p}- \boldsymbol{\nabla}\boldsymbol{\cdot} \boldsymbol{\tau}+2\nu\, \boldsymbol{\nabla} \boldsymbol{\cdot}\tilde{\boldsymbol{\mathsf{s}}}- N\tilde{\phi} \boldsymbol{e}_z+\tilde{\boldsymbol{F}}, \end{gather}$$
$$\begin{gather}D_t\tilde{\phi}={-}\boldsymbol{\nabla} \boldsymbol{\cdot} \boldsymbol{\varSigma}+\kappa\,\nabla^2\tilde{\phi} +N \tilde{u}_z+\tilde{f}, \end{gather}$$
where 
${D_t\equiv \partial _t+\tilde {\boldsymbol {u}}\boldsymbol {\cdot } \boldsymbol {\nabla }}$, 
$\tilde {\phi }\equiv g\tilde {\rho '}/N\rho _r$ is the scaled density field (with dimensions of a velocity), 
$g\equiv \|\boldsymbol {g}\|$, 
$\nu$ and 
$\kappa$ are the kinematic viscosity and thermal diffusivity, respectively, 
$N\equiv \sqrt {-g\zeta /\rho _r}$ is the buoyancy frequency, 
$\boldsymbol {\tau }\equiv \widetilde {\boldsymbol {uu}}-\tilde {\boldsymbol {u}}\tilde {\boldsymbol {u}}$ is the sub-grid stress tensor, 
$\tilde {\boldsymbol{\mathsf{s}}}\equiv (\boldsymbol {\nabla }\tilde {\boldsymbol {u}}+[\boldsymbol {\nabla }\tilde {\boldsymbol {u}}]^{\rm T})/2$ is the filtered strain rate, and 
$\boldsymbol {\varSigma }\equiv \tilde {\boldsymbol {u}\phi }-\tilde {\boldsymbol {u}}\tilde {\phi }$ is the sub-grid vector for the scaled density field. The terms 
$\boldsymbol {F}$ and 
$f$, whose filtered forms appear in the equations above, are the forcing terms associated with mean velocity gradient applied to the flow, and for the SSST under consideration are given by
$$\begin{gather} \boldsymbol{F}={-}\gamma(\boldsymbol{e}_z\boldsymbol{\cdot} \boldsymbol{x}) (\boldsymbol{e}_x \boldsymbol{\cdot} \boldsymbol{\nabla})\boldsymbol{u} -\gamma\boldsymbol{e}_x(\boldsymbol{e}_z\boldsymbol{\cdot} \boldsymbol{u}), \end{gather}$$
$$\begin{gather}f={-}\gamma(\boldsymbol{e}_z\boldsymbol{\cdot} \boldsymbol{x})(\boldsymbol{e}_x \boldsymbol{\cdot} \boldsymbol{\nabla})\phi. \end{gather}$$
Following Germano (Reference Germano1992), the turbulent kinetic energy (TKE) in a given region of size 
$\ell$ is 
$\widetilde {\boldsymbol {u}\boldsymbol {\cdot } \boldsymbol {u}}/2$, and this may be split up into the TKE at scales 
$\geq \ell$ (‘large scales’) denoted by 
$E_K\equiv \|\tilde {\boldsymbol {u}}\|^2/2$, and the TKE at scales 
$< \ell$ (‘small scales’) denoted by 
$e_K\equiv (\widetilde {\|\boldsymbol {u}\|^2}-\|\tilde {\boldsymbol {u}}\|^2)/2=\mathrm {tr}[\boldsymbol {\tau }]/2$. The equations governing 
$E_K$ and 
$e_K$ may be derived from (2.1) and are given by
$$\begin{gather} D_t E_K={-}(1/\rho_r)\,\boldsymbol{\nabla}\boldsymbol{\cdot}(\tilde{\boldsymbol{u}}\tilde{p}){-}\boldsymbol{\nabla} \boldsymbol{\cdot}(\tilde{\boldsymbol{u}}\boldsymbol{\cdot} \boldsymbol{\tau})\,{-}\,\varPi_K+2\nu\,\boldsymbol{\nabla} \boldsymbol{\cdot}(\tilde{\boldsymbol{u}}\boldsymbol{\cdot} \tilde{\boldsymbol{\mathsf{s}}})\,{-}\,2\nu\,\|\tilde{\boldsymbol{\mathsf{s}}}\|^2\,{-}\,N\tilde{{u}}_z\tilde{\phi}+\tilde{\boldsymbol{F}}\boldsymbol{\cdot} \tilde{\boldsymbol{u}}, \end{gather}$$
$$\begin{gather}D_t e_K=\boldsymbol{\nabla} \boldsymbol{\cdot}\boldsymbol{\mathsf{T}}_\boldsymbol{\mathsf{K}} +\mathcal{B}+ \varPi_K-\varepsilon_K+\mathcal{F}_K. \end{gather}$$
In these equations, the scale-to-scale TKE flux is defined as 
${\varPi _K\equiv -\boldsymbol {\tau :}\tilde {\boldsymbol{\mathsf{s}}}}$, such that 
$\varPi _K>0$ corresponds to a transfer of TKE from the large to the small scales. The small-scale buoyancy term is 
$\mathcal {B}\equiv -N(\widetilde {u_z\phi }-\tilde {{u}}_z\tilde {\phi })$, the small-scale TKE dissipation rate is 
$\varepsilon _K\equiv 2\nu (\widetilde {\|\boldsymbol{\mathsf{s}}\|^2}-\|\tilde {\boldsymbol{\mathsf{s}}}\|^2)$, the small-scale forcing is 
$\mathcal {F}_K\equiv \widetilde {\boldsymbol {F} \boldsymbol {\cdot } \boldsymbol {u}}-\tilde {\boldsymbol {F}}\boldsymbol {\cdot }\tilde {\boldsymbol {u}}$, and the small-scale TKE transport term involves
\begin{equation} \boldsymbol{\mathsf{T}}_\boldsymbol{\mathsf{K}}\equiv{-}(1/\rho_r)(\widetilde{\boldsymbol{u}p} -\tilde{\boldsymbol{u}}\tilde{p})-(1/2)(\widetilde{\boldsymbol{u} \|\boldsymbol{u}\|^2}-\tilde{\boldsymbol{u}} \widetilde{\|\boldsymbol{u}\|^2})+(\tilde{\boldsymbol{u}}\boldsymbol{\cdot} \boldsymbol{\tau})+2\nu(\widetilde{\boldsymbol{u}\boldsymbol{\cdot} \boldsymbol{\mathsf{s}}} - \tilde{\boldsymbol{u}}\boldsymbol{\cdot}\tilde{\boldsymbol{\mathsf{s}}}). \end{equation}
The turbulent potential energy (TPE) in a given region of size 
$\ell$ is 
$\widetilde {\phi \phi }/2$, and this may be split up into the amount contained in the large scales denoted by 
$E_P\equiv \tilde {\phi }\tilde {\phi }/2$, and the amount contained in the small scales denoted by 
$e_P\equiv (\widetilde {\phi \phi }-\tilde {\phi }\tilde {\phi })/2$. The equations governing 
$E_P$ and 
$e_P$ may be derived from (2.2) and are given by
$$\begin{gather} D_t E_P={-}\boldsymbol{\nabla} \boldsymbol{\cdot}(\tilde{\phi} \boldsymbol{\varSigma})-\varPi_P +\tfrac{1}{2}\kappa\,\nabla^2\tilde{\phi}\tilde{\phi}-\kappa\,\| \boldsymbol{\nabla}\tilde{\phi}\|^2 +N \tilde{\phi}\tilde{u}_z+\tilde{\phi}\tilde{f}, \end{gather}$$
$$\begin{gather}D_t e_P=\boldsymbol{\nabla} \boldsymbol{\cdot}\boldsymbol{\mathsf{T}}_\boldsymbol{\mathsf{P}}-\varepsilon_P-\mathcal{B} +\varPi_P+\mathcal{F}_P. \end{gather}$$
In these equations, the scale-to-scale TPE flux is defined as 
$\varPi _P\equiv -\boldsymbol {\varSigma } \boldsymbol {\cdot }\boldsymbol {\nabla }\tilde {\phi }$, such that 
$\varPi _P>0$ corresponds to a transfer of TPE from the large to the small scales. The small-scale TPE dissipation rate is 
$\varepsilon _P\equiv \kappa (\widetilde {\|\boldsymbol {\nabla }\phi \|^2}-\|\boldsymbol {\nabla }\tilde {\phi }\|^2)$, the small-scale scalar forcing is 
$\mathcal {F}_P\equiv \tilde {\phi f}-\tilde {\phi }\tilde {f}$, and the small-scale scalar transport involves
\begin{equation} \boldsymbol{\mathsf{T}}_\boldsymbol{\mathsf{P}}\equiv (\kappa/2)\,\boldsymbol{\nabla}(\widetilde{\phi\phi}- \tilde{\phi}\tilde{\phi})+\tilde{\phi}\boldsymbol{\varSigma}. \end{equation}
The equations for 
$E_K$ and 
$E_P$ shown above ((2.5) and (2.8)) were also derived in the previous study of Aluie & Kurien (Reference Aluie and Kurien2011).
In the following, attention will be given to the behaviour of the terms that contribute to 
$D_t e_K$ and 
$D_t e_P$. Note that while 
$e_K$ and 
$e_P$ are referred to as the small-scale TKE and TPE, respectively, they correspond to the TKE and TPE contained in all scales 
$<\ell$, so that when 
$\ell$ exceeds the integral length scale of the flow, 
$e_K$ and 
$e_P$ actually contain the contributions from the largest dynamical scales in the flow. Furthermore, while we choose to focus on the small-scale fields, the questions that we are interested in exploring can also be addressed from a different vantage point by considering the dynamics of the large-scale fields.
A schematic is shown in figure 1 to illustrate the various energy pathways in the flow according to the equations presented above.
Figure 1. Schematic to illustrate the various energy pathways in the flow. The large-scale TKE field gains energy from the forcing 
$\tilde {\boldsymbol {F}}\boldsymbol {\cdot }\tilde {\boldsymbol {u}}$, and loses TKE irreversibly through the large-scale TKE dissipation-rate term 
$\varepsilon _{KE}$. The large-scale TPE field gains energy from the forcing 
$\tilde {\phi }\tilde {f}$, and loses TPE irreversibly through the large-scale TPE dissipation-rate term 
$\varepsilon _{PE}$. The small-scale TKE field gains energy from the forcing term 
$\mathcal {F}_K$, can lose/gain energy reversibly to/from the large-scale TKE field 
$E_K$ through the interscale TKE flux term 
$\varPi _K$, can lose/gain energy reversibly to/from the small-scale TPE field 
$e_P$ through the buoyancy term 
$\mathcal {B}$, and loses TKE irreversibly through the TKE dissipation-rate term 
$\varepsilon _K$. The transport term 
$\boldsymbol {\nabla } \boldsymbol {\cdot }\boldsymbol{\mathsf{T}}_\boldsymbol{\mathsf{K}}$ moves 
$e_K$ conservatively around in space and so is neither a source nor a sink for 
$e_K$. The small-scale TPE field gains energy from the forcing term 
$\mathcal {F}_P$, can lose/gain energy reversibly to/from the large-scale TPE field 
$E_P$ through the interscale TPE flux term 
$\varPi _P$, can lose/gain energy reversibly to/from the small-scale TKE field 
$e_K$ through the buoyancy term 
$\mathcal {B}$, and loses TPE irreversibly through the TPE dissipation-rate term 
$\varepsilon _P$. The transport term 
$\boldsymbol {\nabla } \boldsymbol {\cdot } \boldsymbol{\mathsf{T}}_\boldsymbol{\mathsf{P}}$ moves 
$e_P$ conservatively around in space and so is neither a source nor a sink for 
$e_P$.
2.2. Length scales in SSST
As discussed in Portwood et al. (Reference Portwood, de Bruyn Kops and Caulfield2019), for the velocity field there are four important length scales for the SSST flow under consideration. The first is the large-eddy length scale, which may be characterized by 
$L=\mathcal {E}_K^{3/2}/\langle \epsilon _K\rangle$, where the angle brackets 
$\langle \ \rangle$ denote an ensemble average and 
$\mathcal {E}_K\equiv \lim _{\ell \to \infty }\langle e_K\rangle$ is the total (i.e. involving contributions from all scales) mean TKE in the flow, and 
$\langle \epsilon _K\rangle \equiv \lim _{\ell \to \infty }\langle \varepsilon _K\rangle$ is the total mean TKE dissipation rate (note that here we are assuming an unbounded flow that is suitable for the SSST under consideration). The second length scale is the Kolmogorov length scale 
$\eta =(\nu ^3/\langle \epsilon _K\rangle )^{1/4}$, which characterizes the scale at which viscous and inertial scales in the flow are of the same order. The third is the Ozmidov length scale 
$\ell _O=(\langle \epsilon _K\rangle /N^3)^{1/2}$, which characterizes the scale at which buoyancy and inertial forces are of the same order. The fourth is the Corrsin length scale 
$\ell _C=(\langle \epsilon _K\rangle /\gamma ^3)^{1/2}$, which characterizes the scale at which mean-shear and inertial forces are of the same order. Different non-dimensional flow parameters may be related to these length scales (assuming 
$\nu =\kappa$):
$$\begin{gather} Ri\equiv\frac{N^2}{\gamma^2}=(\ell_C/\ell_O)^{4/3}, \end{gather}$$
$$\begin{gather}Fr\equiv\frac{\langle\epsilon_K\rangle}{N\mathcal{E}_K}= (\ell_O/L)^{2/3}, \end{gather}$$
$$\begin{gather}Re_b\equiv\frac{\langle\epsilon_K\rangle}{\nu N^2}= (\ell_O/\eta)^{4/3}, \end{gather}$$
$$\begin{gather}Re_s\equiv\frac{\langle\epsilon_K\rangle}{\nu \gamma^2}= (\ell_C/\eta)^{4/3}. \end{gather}$$
In order for the flow to become turbulent, the Richardson number must satisfy 
$Ri< O(1)$ (and therefore 
$\ell _C< O(\ell _O)$), and in order for stratification to have an impact on the flow, the Froude number must satisfy 
$Fr< O(1)$. In order for there to exist a range of scales 
$\ell \ll \ell _O$ that are not affected by buoyancy, the buoyancy Reynolds number must satisfy 
$Re_b\gg 1$. Correspondingly, in order for there to exist a range of scales 
$\ell \ll \ell _C$ that are not affected by the mean-shear, the shear Reynolds number must satisfy 
$Re_s\gg 1$. These Reynolds numbers are, however, related, since 
$Re_s=Re_b\,Ri$.
For the scalar field, the large scale is 
$L_P\equiv \mathcal {E}_P^{3/2}\langle \epsilon _K\rangle ^{1/2}/\langle \epsilon _P\rangle ^{3/2}$, where 
$\mathcal {E}_P\equiv \lim _{\ell \to \infty }\langle e_P\rangle$ is the total mean TPE in the flow, and 
$\langle \epsilon _P\rangle \equiv \lim _{\ell \to \infty }\langle \varepsilon _P\rangle$ is the total mean TPE dissipation rate. The small scale is the Batchelor scale 
$\eta _B\equiv (\nu \kappa ^2/\langle \epsilon _K\rangle )^{1/4}$. In this paper, we focus on flows with 
$\nu =\kappa$ such that 
$\eta _B=\eta$.
2.3. Mean-field behaviour
We now turn to consider the contributions to the equations governing the mean energy fields 
$\langle e_K\rangle$ and 
$\langle e_P\rangle$. For SSST, the mean transport equations for these quantities reduce to
$$\begin{gather} 0=\langle\mathcal{B}\rangle-\langle\varepsilon_K\rangle+\langle\varPi_K \rangle+\langle\mathcal{F}_K\rangle, \end{gather}$$
$$\begin{gather}0={-}\langle\mathcal{B}\rangle-\langle\varepsilon_P\rangle+\langle \varPi_P\rangle. \end{gather}$$
In the following discussion, we will first consider the case where the forcing is confined to the large scales (such as in Lindborg Reference Lindborg2006), a flow that we refer to here as forced stably stratified turbulence (FSST), for which we introduce the scale 
$\ell _F$ as the scale below which the forcing plays a sub-leading role in in the flow (note that 
$\ell _F\leq O(L)$). This behaviour will then be compared to that of SSST to understand the differences.
Focusing first on the TKE, in the limit 
$\ell /L \to \infty$ , 
$\langle \varPi _K\rangle \to 0$, 
$\langle \varepsilon _K\rangle \to \langle \epsilon _K\rangle$, 
$\langle \mathcal {F}_K\rangle \to \langle \mathcal {F}_K\rangle ^{\infty }$, 
$\langle \mathcal {B}\rangle \to \langle \mathcal {B}\rangle ^\infty$ , so that we have
\begin{equation} \langle\mathcal{F}_K\rangle^\infty\sim{-}\langle\mathcal{B} \rangle^\infty+\langle\epsilon_K\rangle, \end{equation}
reflecting a balance between the total injection of TKE by the forcing and the total energy lost due to viscous dissipation and conversion to TPE. The behaviour of (2.15) and (2.16) as the scale 
$\ell$ decreases depends on the dynamical scales of the system that were introduced in § 2.2.
For FSST, 
$\langle \varPi _K\rangle \sim -\langle \mathcal {B}\rangle +\langle \epsilon _K\rangle$ when 
$\ell _F\gg \ell \gg \eta$, and for 
$\ell _O/\ell _F\to \infty$ (neutrally buoyant), this would correspond to a TKE cascade 
$\langle \varPi _K\rangle \sim \langle \epsilon _K\rangle$. However, for a stably stratified flow where 
$\langle \mathcal {B}\rangle ^{\infty }<0$, since buoyancy effects reduce with decreasing scale, 
$\langle \varPi _K\rangle \sim -\langle \mathcal {B}\rangle +\langle \epsilon _K\rangle$ implies that the TKE flux 
$\langle \varPi _K\rangle$ will actually reduce as 
$\ell$ decreases, until it approaches a constant value in the regime 
$\ell _O>\ell \gg \eta$, as discussed in Riley & Lindborg (Reference Riley and Lindborg2012) and Kumar, Chatterjee & Verma (Reference Kumar, Chatterjee and Verma2014). As a result of this, at scales where buoyancy is active (i.e. 
$\ell \geq \ell _O$), the TKE flux in stratified turbulence cannot be in the form of a cascade (which would require a constant energy flux). This is simply a reflection of the fact that at these scales, TKE is being lost as it is passed down to smaller scales due to conversion of TKE to TPE.
For SSST, the behaviour of 
$\langle \varPi _K\rangle$ is quite different. In this case, 
$\langle \mathcal {F}_K\rangle$ operates down to the scale 
$\ell _C$, so because 
$\ell _C<\ell _O$, the behaviour 
$\langle \varPi _K\rangle \sim -\langle \mathcal {B}\rangle +\langle \epsilon _K\rangle$ never emerges, and the regime that emerges instead for 
$Ri\ll 1$ is 
$\langle \varPi _K\rangle \sim -\langle \mathcal {F}_K\rangle +\langle \epsilon _K\rangle$ for 
$\ell \gg \eta$. In this regime, the TKE flux is again not in the form of a cascade since 
$\langle \mathcal {F}_K\rangle$ depends on 
$\ell$. Indeed, based on the behaviour of the co-spectrum (Katul et al. Reference Katul, Porporato, Manes and Meneveau2013), we might expect the following behaviour to emerge at high 
$Re$ in SSST:
\begin{equation} \langle\varPi_K\rangle\sim{-}\gamma^2\langle\epsilon_K \rangle^{1/3}\ell^{4/3}+\langle\epsilon_K\rangle, \quad\text{for } L\gg \ell\gg\eta. \end{equation}
Therefore, in this range 
$\langle \varPi _K\rangle$ increases as 
$\ell$ decreases, until it asymptotes to the constant flux cascade regime 
$\langle \varPi _K\rangle \sim \langle \epsilon _K\rangle$ for 
$\ell _C\gg \ell \gg \eta$. This is the opposite behaviour to that discussed in Kumar et al. (Reference Kumar, Chatterjee and Verma2014) for FSST, where the forcing is confined to scales 
$\gg \ell _O$ and where 
$\langle \varPi _K\rangle$ reduces as 
$\ell$ decreases towards 
$\ell _O$.
In the above discussion, it has been assumed that 
$\langle \varepsilon _K\rangle \sim \langle \epsilon _K\rangle$ for 
$\ell \gg \eta$. While this is true for isotropic turbulence, it will not apply in general in strongly stratified flows with 
$Fr\ll 1$. This is because when 
$Fr\ll 1$, the flow structures are highly anisotropic, so while the horizontal length scale may be large enough for viscous effects to be negligible, the vertical length scale of the structure may be small enough for viscous effects to be important (Riley & Lindborg Reference Riley and Lindborg2012). In such a situation, even if 
$\langle \varPi _K\rangle \sim \langle \varepsilon _K\rangle$ emerges in some range of scales, this would not correspond to an inertial cascade since 
$\langle \varepsilon _K\rangle$ would depend on 
$\ell$.
Concerning the TPE field, for flows with 
$\nu =\kappa$, we have for 
$\ell \gg \eta$ that
\begin{equation} \langle\varPi_P\rangle\sim\langle\mathcal{B}\rangle+\langle\epsilon_P\rangle. \end{equation}
Hence in the regime 
$\ell >\ell _O$, as 
$\ell$ decreases, 
$\langle \varPi _P\rangle$ increases, while for 
$\ell _O\gg \ell \gg \eta$, the constant flux cascade regime 
$\langle \varPi _P\rangle \sim \langle \epsilon _P\rangle$ emerges. Therefore, there cannot be a TPE cascade in the strict sense at scales where density is an active scalar, i.e. 
$\ell >\ell _O$, since the scale dependency of 
$\langle \mathcal {B}\rangle$ at these scales leads to a non-constant flux of TPE at these scales. This observation seems to be in conflict with the well-known Bolgiano–Obukhov (BO) scaling (Bolgiano Reference Bolgiano1959; Obukhov Reference Obukhov1959) that has been proposed for stratified turbulent flows, in which 
$\langle \varPi _K\rangle$ is supposed to decay with reducing 
$\ell$ due to buoyancy as 
$\langle \varPi _K\rangle \propto \ell ^{4/5}$ while 
$\langle \varPi _P\rangle$ is constant, 
$\langle \varPi _P\rangle \sim \langle \epsilon _P\rangle$. This behaviour, however, seems problematic. For example, in FSST, for 
$\ell _F\gg \ell \gg \eta$,
$$\begin{gather} \langle\varPi_K\rangle\sim{-}\langle\mathcal{B}\rangle+\langle \epsilon_K\rangle, \end{gather}$$
$$\begin{gather}\langle\varPi_P\rangle\sim\langle\mathcal{B}\rangle+ \langle\epsilon_P\rangle. \end{gather}$$
It would seem that BO scaling could emerge only if 
$\langle \epsilon _P\rangle \gg \langle \epsilon _K\rangle =O(|\langle \mathcal {B}\rangle | )$, as noted in Alam, Guha & Verma (Reference Alam, Guha and Verma2019), so that we could have 
$\langle \varPi _K\rangle \sim -\langle \mathcal {B}\rangle +\langle \epsilon _K\rangle$ in which 
$\langle \varPi _K\rangle$ reduces as 
$\ell$ reduces while 
$\langle \varPi _P\rangle \sim \langle \epsilon _P\rangle$. However, DNS of FSST show 
$\langle \epsilon _P\rangle \lesssim \langle \epsilon _K\rangle$ over a range of 
$Fr$ and 
$Re_b$ (Lindborg Reference Lindborg2006; de Bruyn Kops Reference de Bruyn Kops2015) that includes the ‘moderately stratified’ regime 
$Re\gg 1$ and 
$Fr\approx 1$, for which it has been argued that BO should apply (Alam et al. Reference Alam, Guha and Verma2019). In SSST, BO scaling cannot emerge because there are by definition no scales at which the dynamics are given by a balance between inertia and buoyancy forces, since the mean-shear production is larger than the buoyancy at all scales at which it is active, reflected in that SSST requires 
$\ell _C<\ell _O$.
2.4. Mechanisms governing the TKE and TPE fluxes
Since stable stratification acts to two-dimensionalize the flow, it had been conjectured that such flows might feature an inverse energy cascade (Gage Reference Gage1979; Lilly Reference Lilly1983). However, it has now been demonstrated that this does not occur, but that there is a forward/downscale cascade of both kinetic and potential energy (Riley & de Bruyn Kops Reference Riley and de Bruyn Kops2003; Lindborg Reference Lindborg2005, Reference Lindborg2006; Waite & Bartello Reference Waite and Bartello2006; Brethouwer et al. Reference Brethouwer, Billant, Lindborg and Chomaz2007; Aluie & Kurien Reference Aluie and Kurien2011). The TKE cascade in homogeneous, isotropic turbulence (HIT) is also downscale on average in three dimensions; however, the underlying mechanisms driving the energy transfers in HIT and SSST might be quite different. In the context of HIT, it has long been thought that the key mechanism driving the TKE cascade is that of vortex stretching (VS) (Taylor Reference Taylor1922, Reference Taylor1938; Tennekes & Lumley Reference Tennekes and Lumley1972; Davidson Reference Davidson2004; Doan et al. Reference Doan, Swaminathan, Davidson and Tanahashi2018). However, recent studies have demonstrated quantitatively that while VS plays an important role, the largest contribution to the TKE cascade comes from the dynamical process of the self-amplification of the strain-rate field (Carbone & Bragg Reference Carbone and Bragg2020; Johnson Reference Johnson2020, Reference Johnson2021). An important question is how this understanding applies to SSST, where effects such as internal waves and mean-shear can play a role (whether directly or indirectly) in how TKE is transferred between scales, as well as the question of the mechanism driving the TPE transfer.
In Johnson (Reference Johnson2020, Reference Johnson2021), a powerful, exact relationship was derived for 
$\varPi _K$ that assumes only that the filtering kernel 
$\mathcal {G}_\ell$ used in constructing 
$\tilde {\boldsymbol {u}}$ is Gaussian. The result is
$$\begin{gather} \varPi_K=\varPi_K^{F,SSA}+\varPi_K^{F,VS}+\varPi_K^{SG,SSA}+\varPi_K^{SG,VS}+\varPi_K^{SG,C}, \end{gather}$$
$$\begin{gather}\varPi_K^{F,SSA}={-}\ell^2\tilde{\boldsymbol{\mathsf{s}}}\boldsymbol{:}(\tilde{\boldsymbol{\mathsf{s}}}\boldsymbol{\cdot}\tilde{\boldsymbol{\mathsf{s}}}), \end{gather}$$
$$\begin{gather}\varPi_K^{F,VS}=(1/4)\ell^2 \tilde{\boldsymbol{\mathsf{s}}} \boldsymbol{:}(\tilde{\boldsymbol{\omega}} \tilde{\boldsymbol{\omega}}), \end{gather}$$
where superscript 
$F$ denotes that the quantity depends on only the filtered fields, while superscript 
$SG$ denotes that the quantity depends on the sub-grid fields as well as the filtered fields. Explicit integral formulae for the 
$SG$ terms can be found in Johnson (Reference Johnson2020, Reference Johnson2021); we do not quote them here as they will not be considered in detail in the present paper. Note that in Johnson (Reference Johnson2020, Reference Johnson2021), the 
$F$ contributions are referred to as the ‘local’ contributions, while the 
$SG$ contributions are referred to as ‘non-local’ contributions. While the use of this terminology is appropriate in studies of turbulence (L'vov & Falkovich Reference L'vov and Falkovich1992; Eyink Reference Eyink2005; Aluie & Eyink Reference Aluie and Eyink2009), we prefer to avoid it, because typically in turbulence, a local energy flux is taken to imply an energy flux dominated by the interaction between scales of similar size, which is not necessarily a property possessed by 
$\varPi _K^{F,SSA}$ and 
$\varPi _K^{F,VS}$, which in principle may involve the interaction of any scales in the flow of size 
$\geq \ell$ (Aluie & Eyink Reference Aluie and Eyink2010).
The term 
$\varPi _K^{F,SSA}$ describes the contribution to the TKE flux arising from the process of strain self-amplification (SSA) whereby the filtered strain-rate field 
$\tilde {\boldsymbol{\mathsf{s}}}$ interacts with itself due to nonlinearity to either amplify (if 
$\tilde {\boldsymbol{\mathsf{s}}}\boldsymbol {:}(\tilde {\boldsymbol{\mathsf{s}}} \boldsymbol {\cdot }\tilde {\boldsymbol{\mathsf{s}}})<0$) or else suppress (if 
$\tilde {\boldsymbol{\mathsf{s}}}\boldsymbol {:}(\tilde {\boldsymbol{\mathsf{s}}} \boldsymbol {\cdot }\tilde {\boldsymbol{\mathsf{s}}})>0$) the strain-rate magnitude 
$\|\tilde {\boldsymbol{\mathsf{s}}}\|$. The term 
$\varPi _K^{F,VS}$ describes the contribution to the TKE flux arising from the process of VS whereby the filtered strain-rate field 
$\tilde {\boldsymbol{\mathsf{s}}}$ either amplifies (if 
$\tilde {\boldsymbol{\mathsf{s}}} \boldsymbol {:}(\tilde {\boldsymbol {\omega }}\tilde {\boldsymbol {\omega }})>0$) or else suppresses (if 
$\tilde {\boldsymbol{\mathsf{s}}} \boldsymbol {:}(\tilde {\boldsymbol {\omega }}\tilde {\boldsymbol {\omega }})<0$) the enstrophy 
$\|\tilde {\boldsymbol {\omega }}\|^2$ through the process of vortex stretching (or compression).
The average of the contribution from the filtered field 
$\varPi _K^{F}=\varPi _K^{F,SSA}+\varPi _K^{F,VS}$ in (2.22) is similar to the expression for the TKE flux in a two-point Kármán–Howarth equation derived using filtering and an asymptotic expansion (Carbone & Bragg Reference Carbone and Bragg2020). Moreover, due to the relation of Betchov (Betchov Reference Betchov1956; Eyink Reference Eyink2006), 
$\langle \varPi _K^{F,SSA}\rangle =3\langle \varPi _K^{F,VS}\rangle$ for incompressible homogeneous turbulence, so that the contribution of SSA to 
$\langle \varPi _K^{F}\rangle$ is three times larger than that from VS. The sub-grid contributions 
$\varPi _K^{SG,SSA}$ and 
$\varPi _K^{SG,VS}$ are similar to 
$\varPi _K^{F,SSA}$ and 
$\varPi _K^{F,VS}$ except that in the 
$SG$ contributions, the filtered strain rate 
$\tilde {\boldsymbol{\mathsf{s}}}$ acts to amply the strain rate at sub-grid scales (for 
$\varPi _K^{SG,SSA}$) and vorticity at sub-grid scales (for 
$\varPi _K^{SG,VS}$). DNS data in Johnson (Reference Johnson2020, Reference Johnson2021) showed that 
$\langle \varPi _K^{SG,SSA}\rangle \approx \langle \varPi _K^{SG,VS}\rangle$ in the inertial range of isotropic turbulence, so that overall, the SSA mechanism contributes more to the average energy cascade 
$\langle \varPi _K\rangle$ than VS does, contrary to the traditional explanation according to which VS is seen as the key mechanism driving the energy cascade (Taylor Reference Taylor1938; Tennekes & Lumley Reference Tennekes and Lumley1972; Davidson Reference Davidson2004; Doan et al. Reference Doan, Swaminathan, Davidson and Tanahashi2018).
The result in (2.22) also applies to SSST because it makes no assumption about the flow dynamics. As such, in SSST, SSA and VS will still be the key dynamical processes governing the TKE flux. Moreover, the relation of Betchov (Reference Betchov1956) still applies (since it assumes only incompressibility and homogeneity; it does not assume isotropy). Therefore, at least with respect to the filtered field contribution 
$\langle \varPi _K^{F}\rangle$, SSA still contributes three times as much as VS to the total TKE flux. The relative contribution of 
$\langle \varPi _K^{SG,SSA}\rangle$ and 
$\langle \varPi _K^{SG,VS}\rangle$ to 
$\langle \varPi _K^{SG}\rangle$, however, may differ from that in HIT. Moreover, the actual behaviour of the SSA and VS processes in SSST may differ appreciably from that in isotropic turbulence, since, for example, in SSST, internal waves can contribute to the behaviour of 
$\tilde {\boldsymbol{\mathsf{s}}}$ and 
$\tilde {\boldsymbol {\omega }}$, and anisotropy in the flow can modify the alignments between 
$\tilde {\boldsymbol{\mathsf{s}}}$ and 
$\tilde {\boldsymbol {\omega }}$ (see Sujovolsky & Mininni Reference Sujovolsky and Mininni2020) which affects the VS process. Some of these more involved questions will be the subject of a forthcoming work. A key point to be explored here is the relative contributions of 
$\varPi _K^{F}$ and 
$\varPi _K^{SG}$ to 
$\varPi _K$, and the correlations between these terms, which is important to understand for LES modelling of SSST.
Following the procedure outlined in Johnson (Reference Johnson2020, Reference Johnson2021), a result analogous to (2.22) can be derived for 
$\varPi _P$:
$$\begin{gather} \varPi_P=\varPi_P^{F}+\varPi_P^{SG,S}+\varPi_P^{SG,V}, \end{gather}$$
$$\begin{gather}\varPi_P^{F}={-}\ell^2\tilde{\boldsymbol{\mathsf{s}}}\boldsymbol{:}(\boldsymbol{\nabla}\tilde{\phi}\, \boldsymbol{\nabla}\tilde{\phi}), \end{gather}$$
$$\begin{gather}\varPi_P^{SG,S}={-}\boldsymbol{\nabla}\tilde{\phi}\boldsymbol{\cdot} \int_0^{\ell^2} \tau_\beta(\boldsymbol{\nabla} \tilde{\phi}^{\sqrt{\alpha}}, \tilde{\boldsymbol{\mathsf{s}}}^{\sqrt{\alpha}})\,{\rm d}\alpha, \end{gather}$$
$$\begin{gather}\varPi_P^{SG,V}={-}\boldsymbol{\nabla}\tilde{\phi}\boldsymbol{\cdot} \int_0^{\ell^2} \tau_\beta(\boldsymbol{\nabla}\tilde{\phi}^{\sqrt{\alpha}}, \tilde{\boldsymbol{r}}^{\sqrt{\alpha}})\,{\rm d}\alpha, \end{gather}$$
where 
$\widetilde {(\cdot )}^{\sqrt {\alpha }}$ denotes filtering at scale 
$\sqrt {\alpha }$ (rather than scale 
$\ell$), and for arbitrary first-order 
$\boldsymbol {p}$ and second-order 
$\boldsymbol {q}$ tensor fields, 
$\tau _\beta$ is defined as 
$\tau _\beta (\boldsymbol {p},\boldsymbol {q})\equiv \widetilde {\boldsymbol {p} \boldsymbol {\cdot }\boldsymbol {q}}^{\beta }-\tilde {\boldsymbol {p}}^{\beta }\boldsymbol {\cdot }\tilde {\boldsymbol {q}}^{\beta }$, with 
$\beta \equiv \sqrt {\ell ^2-\alpha }$, and 
${\boldsymbol{\mathsf{r}}}\equiv (\boldsymbol {\nabla }{\boldsymbol {u}}-[\boldsymbol {\nabla }{\boldsymbol {u}}]^{\rm T})/2$ is the rotation-rate tensor. The rotational motion in the flow makes no explicit contribution to the filtered flux 
$\varPi _P^{F}$ because 
$\tilde {\boldsymbol{\mathsf{r}}}\boldsymbol {:}(\boldsymbol {\nabla }\tilde {\phi }\, \boldsymbol {\nabla }\tilde {\phi })=0$, although it does affect 
$\varPi _P^{F}$ implicitly since rotation in the fluid affects the alignment between 
$\tilde {\boldsymbol{\mathsf{s}}}$ and 
$\boldsymbol {\nabla }\tilde {\phi }$.
The contribution from the filtered field 
$\varPi _P^{F}$ describes the flux of TPE associated with the amplification (if 
$\tilde {\boldsymbol{\mathsf{s}}}\boldsymbol {:}(\boldsymbol {\nabla }\tilde {\phi }\, \boldsymbol {\nabla }\tilde {\phi })<0$) or suppression (if 
$\tilde {\boldsymbol{\mathsf{s}}}\boldsymbol {:}(\boldsymbol {\nabla }\tilde {\phi }\, \boldsymbol {\nabla }\tilde {\phi })>0$) of 
$\|\boldsymbol {\nabla }\tilde {\phi }\|$ due to the filtered strain rate 
$\tilde {\boldsymbol{\mathsf{s}}}$. The contribution 
$\varPi _P^{F}$ is similar in form to 
$\varPi _K^{F,VS}$, in that both depend upon the strain-rate field 
$\tilde {\boldsymbol{\mathsf{s}}}$ acting to amplify or suppress another dynamical field, rather than SSA in which the strain-rate field amplifies or suppresses itself. However, the negative sign appearing in 
$\varPi _P^{F}$ but not in 
$\varPi _K^{F,VS}$ leads to a significant difference in how the strain rate contributes to these energy fluxes. In particular, if we define 
$\tilde {\lambda }_i$ and 
$\tilde {\boldsymbol {v}}_i$ as the 
$i$th ordered eigenvalues and eigenvectors of 
$\tilde {\boldsymbol{\mathsf{s}}}$, then we may write the expressions as (see Carbone & Bragg Reference Carbone and Bragg2020; Johnson Reference Johnson2021; Ballouz & Ouellette Reference Ballouz and Ouellette2018) 
$\langle \varPi _K^{F,VS}\rangle = (1/4)\ell ^2\sum _{i=1}^3\langle \tilde {\lambda }_i\,\|\tilde {\boldsymbol {\omega }}\|^2 \cos ^2\theta _{\omega,i}\rangle$ and 
$\langle \varPi _P^{F}\rangle = -\ell ^2\sum _{i=1}^3\langle \tilde {\lambda }_i\,\|\boldsymbol {\nabla }\tilde {\phi }\|^2 \cos ^2\theta _{\boldsymbol {\nabla }\phi,i}\rangle$, where 
$\theta _{\omega,i}$ denotes the angle between 
$\tilde {\boldsymbol {\omega }}$ and 
$\boldsymbol {v}_i$, and similarly for 
$\theta _{\boldsymbol {\nabla }\phi,i}$. In view of this, if 
$\langle \varPi _K^{F,VS}\rangle >0$ so that VS contributes to the average downscale TKE flux, then 
$\langle \tilde {\lambda }_i\,\|\tilde {\boldsymbol {\omega }}\|^2 \cos ^2\theta _{\omega,i}\rangle$ must be dominated by contributions from the extensional eigendirections of 
$\tilde {\boldsymbol{\mathsf{s}}}$. In isotropic turbulence and for 
$\ell /\eta \to 0$, 
$\tilde {\boldsymbol {\omega }}$ aligns preferentially with 
$\boldsymbol {v}_2$, but the contribution to 
$\sum _{i=1}^3\langle \tilde {\lambda }_i\,\|\tilde {\boldsymbol {\omega }}\|^2 \cos ^2\theta _{\omega,i}\rangle$ from 
$\tilde {\lambda }_1$ is larger than that from 
$\tilde {\lambda }_2$ since 
$\tilde {\lambda }_1$ tends to be much larger than the positive values of 
$\tilde {\lambda }_2$ (Tsinober Reference Tsinober2001; Buaria, Bodenschatz & Pumir Reference Buaria, Bodenschatz and Pumir2020). This is all the more true for 
$\ell /\eta$ in the inertial range, where the contribution from 
$\tilde {\lambda }_1$ to 
$\sum _{i=1}^3\langle \tilde {\lambda }_i\,\|\tilde {\boldsymbol {\omega }}\|^2 \cos ^2\theta _{\omega,i}\rangle$ is much larger than that from 
$\tilde {\lambda }_2$ (Carbone & Bragg Reference Carbone and Bragg2020). In an analogous way, if 
$\langle \varPi _P^{F}\rangle = -\ell ^2\sum _{i=1}^3\langle \tilde {\lambda }_i\|\boldsymbol {\nabla }\tilde {\phi }\|^2 \cos ^2\theta _{\boldsymbol {\nabla }\phi,i}\rangle >0$, then 
$\sum _{i=1}^3\langle \tilde {\lambda }_i\,\|\boldsymbol {\nabla }\tilde {\phi }\|^2 \cos ^2\theta _{\boldsymbol {\nabla }\phi,i}\rangle$ must be dominated by contributions from the compressive eigendirections of 
$\tilde {\boldsymbol{\mathsf{s}}}$. Hence, while the VS contribution to the downscale TKE flux is generated by the amplification of 
$\tilde {\boldsymbol {\omega }}$ due to the extensional straining motions in the flow, the downscale TPE flux is generated by the amplification of 
$\boldsymbol {\nabla }\tilde {\phi }$ due to compressional straining motions in the flow. Similar behaviour was also suggested previously based on an LES model for scalar gradients in turbulence (Leonard Reference Leonard1997; Higgins, Parlange & Meneveau Reference Higgins, Parlange and Meneveau2004).
The sub-grid TPE fluxes associated with the strain-rate field 
$\varPi _P^{SG,S}$ and rotation-rate field 
$\varPi _P^{SG,V}$ do not have simple interpretations (unlike 
$\varPi _K^{SG,SSA}$ and 
$\varPi _K^{SG,VS}$), but are related to the amplification of the scalar gradients at multiple scales. It is interesting, however, that the rotation-rate (and therefore vorticity) field makes an explicit contribution to the TPE flux through the sub-grid term 
$\varPi _P^{SG,V}$, even though it makes no explicit contribution to the filtered flux 
$\varPi _P^{F}$.
The discussion above highlights that the physical mechanisms governing the TKE and TPE fluxes in stratified turbulence are quite different, so that while both are positive on average (e.g. Lindborg Reference Lindborg2006), it may not be reasonable to model them in similar ways. We will explore this further by considering the statistical correlation between 
$\varPi _K$ and 
$\varPi _P$, which will provide insights into the extent to which they might be modelled using similar sub-grid closures in the context of LES.
3. Direct numerical simulations
In the next section, data from DNS will be analysed by computing various terms in (2.6) and (2.9) to understand the processes and mechanisms controlling the behaviour of the TKE and TPE across scales in SSST. The DNS data used are from the data set presented in Portwood et al. (Reference Portwood, de Bruyn Kops and Caulfield2019) and Portwood, de Bruyn Kops & Caulfield (Reference Portwood, de Bruyn Kops and Caulfield2022), which we summarize here. In the DNS, the unfiltered versions of (2.1) and (2.2) are solved with constant mean velocity gradient 
$\gamma$ and mean density gradient 
$\zeta$ using the Fourier pseudospectral scheme described in de Bruyn Kops (Reference de Bruyn Kops2015) and Almalkie & de Bruyn Kops (Reference Almalkie and de Bruyn Kops2012), but with the mean-shear term handled using an integrating factor (Brucker et al. Reference Brucker, Isaza, Vaithianathan and Collins2007; Chung & Matheou Reference Chung and Matheou2012; Sekimoto, Dong & Jiménez Reference Sekimoto, Dong and Jiménez2016).
In order to generate a statistically stationary flow, a method similar to that of Taylor et al. (Reference Taylor, Deusebio, Caulfield and Kerswell2016) is adopted wherein the Richardson number of the flow is adjusted via 
$g$ using a mass–spring damper control system with a target value for the total TKE 
$(1/2)\|\boldsymbol {u}\|^2$. With this method, the Richardson and Froude numbers are emergent in the flow, rather than imposed. In the absence of stratification, the TKE in a homogeneous turbulent shear flow grows with time (while the flow scales remain smaller than the domain size). With stratification, mechanisms such as spontaneous shear instabilities that form in the flow due to the layering provide for a downscale TKE flux. A statistically stationary, homogeneous, sheared, stably stratified flow therefore exists at a special equilibrium point where mechanisms generating a downscale TKE flux prevent the growth of TKE that would occur in the flow in the absence of stratification. Put another way, coupling of the fluctuating momentum to the density fields provides for additional energy dissipation mechanisms, which allows the sheared flow to attain a steady state.
The flow is well-resolved, with isotropic grid spacing and maximum wavenumber 
$k_{max}$ satisfying 
$k_{max}\eta \approx 2$. In order to allow for the development of large anisotropic flow scales, a large domain of size 
$\mathcal {L}_x/\mathcal {L}_y=2$, 
$\mathcal {L}_x/\mathcal {L}_z=4$, 
$\mathcal {L}_x/L=40$ was used, where 
$\mathcal {L}_x,\mathcal {L}_y,\mathcal {L}_z$ are the dimensions of the domain in the 
$x,y,z$ directions. In the 
$x$ direction, 3072 grid points were used. Table 1 summarizes the parameters in the DNS.
Table 1. Table of parameters in DNS.
Figure 2 shows snapshots of the spanwise, streamwise and vertical velocity components normalized by 
$\sqrt {2\mathcal {E}_K/3}$, as well as the fluctuating density field 
$\rho '$ normalized by 
$\sqrt {\langle \rho '\rho '\rangle }$, for a two-dimensional plane in the streamwise and vertical directions. The snapshots illustrate clearly both the strong anisotropy in the flow and the inclination of the flow structures relative to the horizontal due to the mean velocity shear in the flow. That preferential orientation is, however, absent apparently for the vertical component corresponding to the direction in which the buoyancy force acts. Figure 3 shows the same quantities but for the filtered velocity and density fields, which also reveal strong anisotropy and preferential orientation of the flow structures. The difference between the fields visualized in figures 2 and 3 corresponds to the contribution from the small-scale fields.
Figure 2. Snapshots of velocity and density fields from the DNS showing a plane in the streamwise and vertical directions. Normalized quantities are: (a) spanwise component 
$u_y/\sqrt {2\mathcal {E}_K/3}$, (b) streamwise component 
$u_x/\sqrt {2\mathcal {E}_K/3}$, (c) vertical component 
$u_z/\sqrt {2\mathcal {E}_K/3}$, (d) fluctuating density 
$\rho '/\sqrt {\langle \rho '\rho '\rangle }$. Values go from (red, blue), which correspond to 
$(-3,3)$ for the normalized quantities, centred at 
$\text {white}=0$. Velocities are normalized using the same scale 
$\sqrt {2\mathcal {E}_K/3}$ to highlight the anisotropic partitioning of kinetic energy in the flow.
Figure 3. Snapshots of the filtered (with 
$\ell =60\eta$) velocity and density fields from the DNS showing a plane in the streamwise and vertical directions. Normalized quantities are filtered: (a) spanwise component 
$\tilde {u}_y/\sqrt {2{E}_K/3}$, (b) streamwise component 
$\tilde {u}_x/\sqrt {2{E}_K/3}$, (c) vertical component 
$\tilde {u}_z/\sqrt {2{E}_K/3}$, (d) fluctuating density 
$\widetilde {\rho '}/\sqrt {\langle \widetilde {\rho '}\widetilde {\rho '}\rangle }$. Values go from (red, blue), which correspond to 
$(-3,3)$ for the normalized quantities, centred at 
$\text {white}=0$. Velocities are normalized using the same scale 
$\sqrt {2{E}_K/3}$ to highlight the anisotropic partitioning of kinetic energy in the flow.
4. Results and discussion
4.1. Mean-field behaviour
We begin by considering the behaviour of 
$\langle e_K\rangle$, 
$\langle e_P\rangle$ and the diagonal components of 
$\langle \boldsymbol {\tau }\rangle /2$ (which correspond to the TKE associated with different components of 
$\boldsymbol {u}$) as a function of filter length 
$\ell$ in order to understand how energy is partitioned in the flow at different scales. The results in figure 4(a) show that at all scales, 
$\langle e_K\rangle >\langle e_P\rangle$, with 
$\langle e_K\rangle /\langle e_P\rangle =O(10)$ at the large scales of the flow. The flow is therefore far from a state of equipartition of large-scale energy among the TKE and TPE fields, unlike the behaviour that is thought to emerge for strongly stratified flows with 
$Fr\ll 1$ (Billant & Chomaz Reference Billant and Chomaz2001). Below the Corrsin scale 
$\ell _C=(\langle \epsilon _K\rangle /\gamma ^3)^{1/2}$, the difference between 
$\langle e_K\rangle$ and 
$\langle e_P\rangle$ reduces, but remains significant even for 
$\ell /\eta \ll 1$. When considering the components 
$\langle \boldsymbol {\tau }\rangle /2$, it is seen that although the total TKE is much larger than the TPE, the TKE associated with particular components of the velocity field is comparable to the TPE. In particular, 
$\langle \tau _{zz}\rangle /2\approx 0.8\langle e_P \rangle$ and 
$\langle \tau _{yy}\rangle /2\approx 1.3\langle e_P\rangle$ at the large scales. As such, even though the energy contained in the TPE field is small compared to that in the total TKE field, it is of the same order as that contained in the spanwise (
$y$) and vertical (
$z$) directions of the flow, and therefore plays an important energetic role in the system. In terms of scaling, 
$\langle e_K\rangle$ and 
$\langle e_P\rangle$ both show the expected behaviour 
$\langle e_K\rangle \propto \ell ^2$ and 
$\langle e_P\rangle \propto \ell ^2$ at 
$\ell /\eta \leq O(1)$, where viscous effects lead to smooth velocity and density fields. Outside of this, a well-defined scaling regime does not emerge, due to the limited Reynolds number of the flow (as well as the fact that 
$\ell _C/\eta$ is too small for an inertially dominated regime to emerge).
Figure 4. (a) Plot of mean small-scale TKE 
$\langle e_K\rangle$, TPE 
$\langle e_P\rangle$, and diagonal components of 
$\langle \boldsymbol {\tau }\rangle /2$, normalized by total energy 
$\langle E_T\rangle \equiv \lim _{\ell /\eta \to \infty }[\langle e_K\rangle +\langle e_P\rangle ]$, as a function of filter scale 
$\ell$. The thick dotted line indicates scaling 
$\propto \ell ^2$. (b) Plot of mean small-scale TKE 
$\langle \varepsilon _K\rangle$ and TPE 
$\langle \varepsilon _P\rangle$ dissipation rates, normalized by the total turbulent energy dissipation rate 
$\langle \epsilon _T\rangle \equiv \lim _{\ell /\eta \to \infty }[\langle \varepsilon _K\rangle +\langle \varepsilon _P\rangle ]$. The thin vertical dotted lines from right to left are 
$L/\eta,\ell _O/\eta,\ell _C/\eta =126.3,43.9,10.4$, respectively.
In figure 4(b) we consider the mean small-scale turbulent kinetic 
$\langle \varepsilon _K\rangle$ and potential 
$\langle \varepsilon _P\rangle$ energy dissipation rates. For 
$\ell /\eta \to \infty$, these satisfy 
$\langle \varepsilon _K\rangle \to \langle \epsilon _K\rangle$ and 
$\langle \varepsilon _P\rangle \to \langle \epsilon _P\rangle$, and the results show that these are in the ratio 
$\langle \epsilon _K\rangle /\langle \epsilon _P\rangle \approx 5$. Both 
$\langle \varepsilon _K\rangle$ and 
$\langle \varepsilon _P\rangle$ are approximately independent of 
$\ell$ only down to 
$\ell /\eta = O(10)$, consistent with the usual observation that the Kolmogorov scale underestimates the scale at which viscous forces become important (Pope Reference Pope2000).
We now turn to consider the contributions to the equations governing 
$\langle e_K\rangle$ and 
$\langle e_P\rangle$. The results in figure 5 show that the flow exhibits the 
$\ell /L\to \infty$ asymptotic behaviour 
$\langle \mathcal {F}_K\rangle ^\infty \sim -\langle \mathcal {B}\rangle ^\infty +\langle \epsilon _K\rangle$ for 
$\ell /\eta \gtrsim O(100)$, which corresponds to 
$\ell /L \gtrsim O(1)$. As 
$\ell /\eta$ is decreased below 
$O(100)$, 
$\langle \varPi _K\rangle$ begins to increase and becomes positive, indicating a downscale flux of TKE, as has been observed previously for stratified turbulent flows (Riley & de Bruyn Kops Reference Riley and de Bruyn Kops2003; Lindborg Reference Lindborg2005, Reference Lindborg2006; Waite & Bartello Reference Waite and Bartello2006; Brethouwer et al. Reference Brethouwer, Billant, Lindborg and Chomaz2007). However, it does not become significant until 
$\ell /\eta \lesssim O(10)$. Moreover, when it reaches its peak value, it does not exhibit the behaviour 
$\langle \varPi _K\rangle /\langle \epsilon _K\rangle \sim 1$ (the maximum value obtained is 
$\langle \varPi _K\rangle /\langle \epsilon _K\rangle \approx 0.46$) that would be expected for non-stratified isotropic turbulence. This is because 
$\langle \varPi _K\rangle$ cannot grow significantly until 
$\ell$ is small enough such that 
$\langle \mathcal {F}_K\rangle \ll \langle \mathcal {F}_K\rangle ^\infty$, which occurs only for 
$\ell \leq O(\ell _C)$. Once this regime is obtained, since 
$\ell _C<\ell _O$, the TKE balance becomes 
$\langle \varPi _K\rangle \sim \langle \varepsilon _K\rangle$, a behaviour that can be observed in figure 5 to be approached at small scales. However, in the present flow, this behaviour does not give rise to an inertial TKE cascade since at the scales where 
$\langle \varPi _K\rangle \sim \langle \varepsilon _K\rangle$, 
$\langle \varepsilon _K\rangle$ is a decreasing function of 
$\ell$. To observe an inertial TKE cascade with 
$\langle \varPi _K\rangle \sim \langle \epsilon _K\rangle$ would require considering a flow possessing a range of scales 
$\ell _O>\ell _C\gg \ell \gg \eta$. The results in figure 5 show that the contribution to the TKE flux coming from the filtered field dynamics, i.e. 
$\langle \varPi _K^F\rangle$, makes a significant contribution to 
$\langle \varPi _K\rangle$ at scales where 
$\langle \varPi _K\rangle$ plays a significant role in the small-scale TKE budget equation, and also that the contribution involving the sub-grid fields 
$\langle \varPi _K^{SG}\rangle$ is also significant, just as for isotropic turbulence (Johnson Reference Johnson2020, Reference Johnson2021). At sufficiently small 
$\ell$, 
$\langle \varPi _K\rangle \approx \langle \varPi _K^F\rangle$, consistent with the exact limiting behaviour 
$\lim _{\ell /\eta \to 0}\varPi _K^F\to \varPi _K$.
Figure 5. Plot of terms in the average small-scale TKE budget equation (2.15). Thin dotted lines from right to left are 
$L/\eta,\ell _O/\eta,\ell _C/\eta =126.3,43.9,10.4$, respectively.
The mean buoyancy term 
$\langle \mathcal {B}\rangle$ is negative at all scales, indicating a mean transfer of TKE to TPE, though the magnitude of 
$\langle \mathcal {B}\rangle$ is sub-leading compared to the other terms shown in figure 5. Nevertheless, buoyancy plays a key role in the flow because its magnitude is significant compared to the processes governing the vertical motion of the flow. As expected, the magnitude of 
$\langle \mathcal {B}\rangle$ begins to reduce significantly below the Ozmidov scale 
$\ell _O$.
Concerning the TPE behaviour, the results in figure 6 exhibit the 
$\ell /L\to \infty$ asymptotic behaviour 
$\langle \mathcal {B}\rangle ^\infty \sim -\langle \epsilon _P\rangle$ for 
$\ell /\eta \gtrsim O(100)$, which corresponds to 
$\ell /L \gtrsim O(1)$. The behaviour of 
$\langle \varPi _P\rangle$ is very similar to that of 
$\langle \varPi _K\rangle$ in that 
$\langle \varPi _P\rangle$ does not become significant until 
$\ell /\eta \lesssim O(10)$, and does not exhibit the behaviour 
$\langle \varPi _P\rangle /\langle \epsilon _P\rangle \sim 1$ (the maximum value obtained is 
$\langle \varPi _P\rangle /\langle \epsilon _P\rangle \approx 0.61$) that would be expected for a passive scalar advected in isotropic turbulence. Whereas the growth of 
$\langle \varPi _K\rangle$ is constrained in SSST to scales 
$\ell \leq O(\ell _C)$, as discussed previously, the growth of 
$\langle \varPi _P\rangle$ is constrained in SSST to scales 
$\ell \leq O(\ell _O)$. In our DNS, these two ranges are similar, so the ranges of scales over which 
$\langle \varPi _K\rangle$ and 
$\langle \varPi _P\rangle$ are active are similar. However, in an SSST flow with 
$Ri\lll 1$ and 
$Re\ggg 1$, this could lead to an interesting scenario where 
$\langle \varPi _P\rangle /\langle \epsilon _P\rangle \sim 1$ while 
$\langle \varPi _K\rangle /\langle \epsilon _K\rangle \approx 0$ for 
$\ell _O\gg \ell \gg \ell _C$. The results in figure 6 also show that like the TKE results, the contribution to the TPE flux coming from the filtered field dynamics, i.e. 
$\langle \varPi _P^F\rangle$, makes a significant contribution to 
$\langle \varPi _P\rangle$ at scales where 
$\langle \varPi _P\rangle$ plays a significant role in the small-scale TPE budget equation, but that the contribution involving the sub-grid fields 
$\langle \varPi _P^{SG}\rangle$ is also significant. At sufficiently small 
$\ell$, 
$\langle \varPi _P\rangle \approx \langle \varPi _P^F\rangle$, consistent with the exact limiting behaviour 
$\lim _{\ell /\eta \to 0}\varPi _P^F\to \varPi _P$.
Figure 6. Plot of terms in the average small-scale TPE budget equation (2.16). Thin dotted lines from right to left are 
$L/\eta,\ell _O/\eta,\ell _C/\eta =126.3,43.9,10.4$, respectively.
4.2. Fluctuations about the mean-field
In order to understand the energetics of the flow beyond its mean-field behaviour, we will now consider the probability density functions (PDFs) of various quantities. We begin in figures 7(a,b) with the PDFs of the normalized small-scale energies 
$e_K/\langle e_K\rangle$ and 
$e_P/\langle e_P\rangle$ for different filtering lengths 
$\ell /\eta$, with the largest filter scale considered around half of the integral length scale 
$L$. At the small scales, the PDFs are highly non-Gaussian, such that the TKE and TPE fields exhibit frequent large fluctuations about their mean-field behaviour. As the filter length increases, the heavy tails reduce, which is associated with the dominance of the linear mean-shear and buoyancy forces at the large scales of the flow, and the decreasing strength of the nonlinear term. The results also show that the TPE field exhibits stronger fluctuations from the mean-field behaviour than the TKE field (although the difference is not that large), with the difference most notable at intermediate scales. This is consistent with the fact that scalar fields are known to be more intermittent in turbulent flows because they lack the pressure gradient in their dynamics that regulates large fluctuations in the velocity gradient field.
Figure 7. Plots of PDFs of (a) 
$e_K/\langle e_K\rangle$, (b) 
$e_P/\langle e_P\rangle$ and (d) 
$\mathcal {B}/\langle \mathcal {B}\rangle$ for different filter lengths 
$\ell$. In (c), the solid lines correspond to the PDFs of 
$\varPi _K/\langle \varPi _K\rangle$, and the dashed lines correspond to the PDFs of 
$\varPi _P/\langle \varPi _P\rangle$. Different colours/symbols correspond to different 
$\ell /\eta$ as indicated by the legend in (a).
In figure 7(c), we show the PDFs of 
$\varPi _K/\langle \varPi _K\rangle$ and 
$\varPi _P/\langle \varPi _P\rangle$ for different filtering lengths 
$\ell /\eta$. At larger scales, the PDFs have a large variance, showing that at these scales, the fluxes of TKE and TPE can exceed their mean values significantly. The PDFs are also close to being symmetric, so that the probabilities of upscale and downscale fluxes of TKE are similar, and the same for TPE. As one moves to smaller scales where 
$\langle \varPi _K\rangle$ and 
$\langle \varPi _P\rangle$ both increase, the probability of large fluctuations about the mean-field behaviour reduces. However, there is still a significant probability of observing regions of the flow where the local values of 
$\varPi _K$ and 
$\varPi _P$ exceed significantly their mean values in magnitude. As 
$\ell$ is reduced, the mode of the PDF remains close to zero, while the mean and skewness both become positive, with the probability of downscale fluxes significantly exceeding that of upscale fluxes. The PDFs for TKE and TPE are similar in shape, and at scales where 
$\langle \varPi _K\rangle$ and 
$\langle \varPi _P\rangle$ play an important role in the TKE and TPE scalewise budgets, the fluctuations of 
$\varPi _K$ about 
$\langle \varPi _K\rangle$ are larger than those of 
$\varPi _P$ about 
$\langle \varPi _P\rangle$.
In figure 7(d), we show the PDFs of 
$\mathcal {B}/\langle \mathcal {B}\rangle$ for different filtering lengths 
$\ell /\eta$. At larger scales where 
$-\langle \mathcal {B}\rangle$ is the dominant source term in the mean TPE budget, the fluctuations about 
$\langle \mathcal {B}\rangle$ are not that strong, and remarkably, the probability of observing 
$\mathcal {B}/\langle \mathcal {B}\rangle <0$, which would correspond to conversion of TPE to TKE, is zero. For smaller scales, the probability of observing 
$\mathcal {B}/\langle \mathcal {B}\rangle <0$ becomes finite, and increases with decreasing 
$\ell$. The significant fluctuations of 
$\mathcal {B}$ about its mean value indicate that the effects of buoyancy could be felt at scales considerably smaller than the mean-field description suggests.
The significant increase in the width of the PDFs observed in figures 7(c,d) for 
$\varPi _K/\langle \varPi _K\rangle$ and 
$\varPi _P/\langle \varPi _P\rangle$ as 
$\ell /\eta$ increases, and for 
$\mathcal {B}/\langle \mathcal {B}\rangle$ as 
$\ell /\eta$ decreases, is due mainly to the fact that the standard deviations of the variables are much larger than their mean values in these ranges of 
$\ell /\eta$. This is illustrated in figure 8 where we show the standard deviations of the variables 
$\varPi _K, \varPi _P, -\mathcal {B}$ normalized by their means, as functions of scale 
$\ell /\eta$.
Figure 8. Plots of 
$\sqrt {\langle \xi -\langle \xi \rangle \rangle ^2}/\langle \xi \rangle$ for 
$\xi =\varPi _K$, 
$\xi =\varPi _P$ and 
$\xi =-\mathcal {B}$ as functions of filter scale 
$\ell /\eta$ to illustrate how the standard deviation of a variable compares to its mean value at different scales.
In figure 9, we show the joint PDFs of 
$\varPi _P/\langle \varPi _P\rangle$ and 
$\varPi _K/\langle \varPi _K\rangle$ for different filtering lengths 
$\ell /\eta$. The results reveal that there is a very weak positive correlation between 
$\varPi _P$ and 
$\varPi _K$ across all scales, and the correlation is largest at 
$\ell /\eta =16$. Indeed, the probability of observing a downscale TKE flux with an upscale TPE flux, and vice versa, is quite significant. This relatively weak correlation can be understood in view of the discussion in § 2.4 regarding the mechanisms governing the TKE and TPE fluxes. It is common to model TKE and TPE fluxes using similar closures but with a coefficient that accounts for their differences only through the Prandtl number effects. However, the weak correlation between 
$\varPi _K$ and 
$\varPi _P$ revealed in figure 9 shows that such an approach is not suitable, and that the closures need to reflect to some extent the different physical processes dominating these fluxes, which leads to the relatively weak correlation between 
$\varPi _K$ and 
$\varPi _P$. In the context of LES, the contributions to 
$\varPi _K$ and 
$\varPi _P$ from the filtered field dynamics, namely 
$\varPi _K^F\equiv \varPi _K^{F,SSA}+\varPi _K^{F,VS}$ and 
$\varPi _P^F$ (see § 2.4), are closed since they depend on only the filtered field quantities, not the sub-grid fields. It is therefore of particular interest to understand how these terms contribute to the total TKE and TPE fluxes in the flow. To explore this, in figure 10 we consider the joint PDF of 
$\varPi _K$ and 
$\varPi _K^F$ for different filtering lengths 
$\ell /\eta$. The filtered fluxes satisfy 
$\lim _{\ell /\eta \to 0}\varPi _K^F\to \varPi _K$, and the results for 
$\ell /\eta =0.25$ are approaching this regime, showing a very small probability of events deviating from the state 
$\varPi _K=\varPi _K^F$. As 
$\ell /\eta$ is increased, there remains a significant positive correlation between the variables 
$\varPi _K$ and 
$\varPi _K^F$, indicating that 
$\varPi _K^F$ makes an important contribution to the total flux 
$\varPi _K$, consistent with the mean-field results in figure 5. However, there is a significant probability of events where either 
$\varPi _K^F$ and 
$\varPi _K$ have significantly different magnitudes, or else even have opposite signs. In terms of the physical mechanisms discussed in § 2.4, this means that at a given scale 
$\ell$, if the strain and vorticity fields are being amplified at that scale by nonlinearity such that 
$\varPi _K^F>0$, this may nevertheless not lead to a downscale flux of TKE associated with 
$\varPi _K>0$ since at scales smaller than 
$\ell$, the strain and vorticity fields may be experiencing suppression (rather than amplification) due to nonlinearity yielding 
$\varPi _K^{SG}<0$, and if this is strong enough, then 
$\varPi _K=\varPi _K^F+\varPi _K^{SG}<0$. The spread of the PDF about the line 
$\varPi _K=\varPi _K^F$ also implies that the relative contribution of 
$\varPi _K^F$ to 
$\varPi _K$ is similar both during events where 
$\varPi _K\sim \langle \varPi _K\rangle$ and in large fluctuations where 
$|\varPi _K|\ll \langle \varPi _K\rangle$ or 
$|\varPi _K|\gg \langle \varPi _K\rangle$. For 
$\ell /\eta =60$, the results show that 
$\varPi _K$ and 
$\varPi _K^F$ are almost uncorrelated, and at this scale, 
$\langle \varPi _K\rangle /\langle \varPi _K^F\rangle \approx 43$, so that the filtered field makes a small, uncorrelated contribution to 
$\varPi _K$, which is itself very small at 
$\ell /\eta =60$. The results for the joint PDF of 
$\varPi _P$ and 
$\varPi _P^F$ are shown in figure 11. These are very similar to the TKE flux results, with 
$\varPi _P^F\sim \varPi _P$ observed at the smallest scales, and a significant correlation between the two at intermediate scales. The correlation between 
$\varPi _P^F$ and 
$\varPi _P$ is, however, slightly stronger than that between 
$\varPi _K$ and 
$\varPi _K^F$.
Figure 9. Contour plot of the logarithm of the joint PDF of 
$\varPi _P/\langle \varPi _P\rangle$ and 
$\varPi _K/\langle \varPi _K\rangle$ for (a) 
$\ell /\eta =0.25$, (b) 
$\ell /\eta =6$, (c) 
$\ell /\eta =16$, (d) 
$\ell /\eta =60$. Colours correspond to the logarithm of the PDF, and the correlation coefficient is shown at the top of each plot.
Figure 10. Contour plot of the logarithm of the joint PDF of 
$\varPi _K/\langle \varPi _K\rangle$ and 
$\varPi _K^F/\langle \varPi _K^F\rangle$ for (a) 
$\ell /\eta =0.25$, (b) 
$\ell /\eta =6$, (c) 
$\ell /\eta =16$, (d) 
$\ell /\eta =60$. Colours correspond to the logarithm of the PDF, and the correlation coefficient is shown at the top of each plot.
Figure 11. Contour plot of the logarithm of the joint PDF of 
$\varPi _P/\langle \varPi _P\rangle$ and 
$\varPi _P^F/\langle \varPi _P^F\rangle$ for (a) 
$\ell /\eta =0.25$, (b) 
$\ell /\eta =6$, (c) 
$\ell /\eta =16$, (d) 
$\ell /\eta =60$. Colours correspond to the logarithm of the PDF, and the correlation coefficient is shown at the top of each plot.
In figure 12, we consider the joint PDF of 
$\mathcal {B}$ and 
$\varPi _K$ for different filter lengths 
$\ell$. The results show that as 
$\ell$ is decreased, the PDF reorients from being extended along the 
$\varPi _K$ axis, to being extended along the 
$\mathcal {B}$ axis. Furthermore, the peak of the PDF (indicating the mode) is along 
$\varPi _K/\langle \varPi _K \rangle \approx 0$ for small scales and 
$\mathcal{B}/\langle \mathcal{B} \rangle >0$ for larger scales, which is associated with the different levels of skewness of these random variables at these scales. This corresponds to the transition from larger scales where buoyancy plays a more dominant role in the TKE energetics, to smaller scales where the nonlinear energy flux plays a dominant role. As also observed in figure 7(d), at larger scales the probability of observing 
$\mathcal {B}/\langle \mathcal {B}\rangle <0$ is very low, suggesting that convective motion is very rare at these scales. As 
$\ell$ decreases, however, the probability of observing 
$\mathcal {B}/\langle \mathcal {B}\rangle <0$ increases significantly. At 
$\ell /\eta =16$, 
$-\langle \mathcal {B}\rangle \approx \langle \varPi _K\rangle$, and the results in figure 12 show that at this scale, convective motion 
$\mathcal {B}/\langle \mathcal {B}\rangle <0$ can occur, but the PDF is strongly skewed towards stably stratified regions that have 
$\mathcal {B}/\langle \mathcal {B}\rangle >0$. Moreover, for this scale the PDF is stretched significantly along the 
$\varPi _K$ axis, showing that fluctuations of 
$\varPi _K$ are considerably stronger than those of 
$\mathcal {B}$ about their respective mean values, which are approximately equal in magnitude at this scale. The results for the joint PDF of 
$\mathcal {B}$ and 
$\varPi _P$ are similar to those for the joint PDF of 
$\mathcal {B}$ and 
$\varPi _K$, therefore for brevity we do not show them.
Figure 12. Contour plot of the logarithm of the joint PDF of 
$\mathcal {B}/\langle \mathcal {B}\rangle$ and 
$\varPi _K/\langle \varPi _K\rangle$ for (a) 
$\ell /\eta =0.25$, (b) 
$\ell /\eta =6$, (c) 
$\ell /\eta =16$, (d) 
$\ell /\eta =60$. Colours correspond to the logarithm of the PDF, and the correlation coefficient is shown at the top of each plot.
Finally, given the importance of the velocity gradient dynamics for the TKE flux, in figure 13 we plot the results for the joint PDF of the filtered velocity gradient invariants 
$Q\equiv -{\boldsymbol {\nabla }\tilde {\boldsymbol {u}}\boldsymbol {:}\boldsymbol {\nabla }\tilde {\boldsymbol {u}}}/2$ and 
$R\equiv -{(\boldsymbol {\nabla }\tilde {\boldsymbol {u}}\boldsymbol {\cdot }\boldsymbol {\nabla }\tilde {\boldsymbol {u}}) \boldsymbol {:}\boldsymbol {\nabla }\tilde {\boldsymbol {u}}}/3$ at different filter scales 
$\ell /\eta$. In Danish & Meneveau (Reference Danish and Meneveau2018), this joint PDF was considered for isotropic turbulence, and they found that the classic sheared-drop shape of the isoprobability lines is preserved as 
$\ell /\eta$ is varied over the range considered, while other statistical characterizations of the velocity gradients also reveal behaviour that is qualitatively similar as the filter scale is increased (Tom, Carbone & Bragg Reference Tom, Carbone and Bragg2020). For 
$\ell /\eta =0.25$, the shape of the PDF in figure 13 is identical to that observed in isotropic turbulence (Meneveau Reference Meneveau2011), with the contours extending down the right Viellefosse tail, and a strong preference for quadrants 
$Q>0$, 
$R<0$ and 
$Q<0$, 
$R>0$. However, the results in figure 13 show that at sufficiently large scales, the shape of the PDF changes from the sheared-drop shape and becomes more symmetric. Such behaviour has also been observed for sufficiently large scales in homogeneous isotropic turbulence (Chertkov, Pumir & Shraiman Reference Chertkov, Pumir and Shraiman1999; Naso & Pumir Reference Naso and Pumir2005) and rotating turbulence (Naso & Godeferd Reference Naso and Godeferd2012), as well as compressible turbulence (Wang et al. Reference Wang, Wan, Chen, Xie, Zheng, Wang and Chen2020). Consequently, the change in the shape of the PDF contours shown in figure 13 as 
$\ell /\eta$ is increased could be caused by the fact that the largest filter lengths considered are approaching large/integral scales where nonlinearity weakens, or could perhaps be due to the dominance of linear forces (e.g. buoyancy and mean-shear) at these scales.
Figure 13. Contour plot of the logarithm of the joint PDF of 
$Q$ and 
$R$ for (a) 
$\ell /\eta =0.25$, (b) 
$\ell /\eta =6$, (c) 
$\ell /\eta =16$, (d) 
$\ell /\eta =60$. Colours correspond to the logarithm of the PDF.
5. Conclusions
In this paper, we have analysed the scale-dependent TKE and TPE in sheared, stably stratified turbulence (SSST) using a filtering approach, where the flow has constant mean velocity gradient and mean density gradient. Equations for the scale-dependent TKE and TPE are to explore the competing effects in the flow as well as the physical mechanisms governing the TKE and TPE fluxes between scales. Various quantities in these equations were then evaluated using data from DNS of SSST, with attention given to both the mean-field behaviour of the flow, as well as fluctuations about this mean-field state.
In terms of the mean-field properties, while the mean TKE 
$\langle e_K\rangle$ is larger than the mean TPE 
$\langle e_P\rangle$ by an order of magnitude at the large scales, the difference between them reduces with reducing scale 
$\ell$. Nevertheless, the TPE plays an important energetic role in the system, and is larger than the vertical component of TKE at all scales in the flow. The mean small-scale TKE dissipation rate 
$\langle \varepsilon _K\rangle$ is also significantly larger than the potential dissipation rate 
$\langle \varepsilon _P\rangle$. The mean TKE and TPE fluxes between scales, 
$\langle \varPi _K\rangle$ and 
$\langle \varPi _P\rangle$, respectively, do not reveal a cascade regime in the flow due to the impact of the mean-shear down to relatively small scales in the flow. The contributions to these fluxes from the filtered fields are also shown to be significant at scales where the energy fluxes play significant roles in the energy budget equations; however, the contributions from the sub-grid fields are also important, just as has been observed for isotropic turbulence.
To understand the flow energetics beyond the mean-fields, PDFs of various quantities (normalized by their mean values) have been studied and discussed. The PDFs of small-scale TKE and TPE are highly non-Gaussian at the smallest scales, indicating large fluctuations of the TKE and TPE about their mean-field values. The TPE shows stronger fluctuations from its mean-field value than the TKE, which is consistent with the known result that scalar fields are more intermittent than velocity fields in turbulence, since scalar fields lack a pressure gradient term in their equation that can act to suppress large fluctuations. As the filter scale is increased, these PDFs approach a Gaussian distribution at scales where the linear forces in the flow (mean-shear and buoyancy) dominate the dynamics. At larger scales in the flow, buoyancy seems to convert TKE to TPE always, with the PDF of the buoyancy term 
$\mathcal {B}$ showing almost zero probability of locally convective events where 
$\mathcal {B}>0$. The probability of locally convective regions increases as the filter scale is decreased, however, the probability remains small at scales downward of the Ozmidov scale. The TKE and TPE fluxes between scales are on average positive (downscale flux), with their instantaneous values correlated positively, and increasingly so with decreasing filter scale. However, the correlation is quite weak, and this is because the physical mechanisms governing the TKE and TPE fluxes are quite different, as discussed in our theoretical analysis in § 2.4. Indeed, the joint PDFs of the TKE and TPE fluxes reveal a significant probability of events where the TKE and TPE fluxes may have very different magnitudes, and even opposite signs. Finally, the PDF of the principal invariants of the filtered velocity gradients (the so-called 
$Q,R$ invariants; Tsinober Reference Tsinober2001) reveals the classical sheared-drop shape for the smallest filter scales. However, the PDF becomes increasingly symmetric about 
$R=0$ as the filter scale increases. This change in shape may be caused by the fact that the largest filter lengths considered are approaching large/integral scales where nonlinearity weakens, or perhaps also due to the dominance of the linear forces associated with buoyancy and mean-shear at these scales.
This contrasts with the behaviour observed in isotropic turbulence (Danish & Meneveau Reference Danish and Meneveau2018) and is due to the increasingly dominant effects of the linear mean-shear and buoyancy effects over the nonlinear inertial effects at these scales. This symmetry that emerges means that vortex stretching and compression become equally probable, as well as strain self-amplification and self-suppression.
One important direction for future work is to compute and explore the sub-grid scale contributions to the TKE and TPE interscale fluxes. Exploring these terms, probing in more detail the physical processes that govern them, and analysing their relation to the filtered scale flow would provide the insights required for building LES models of SSST. Work on this area will be the subject of a forthcoming article. It would also be desirable to consider a DNS of SSST flow with sufficiently large 
$Re_b$ to observe a constant flux energy cascade at scales smaller than the Ozmidov and Corrsin scales, something that could not be observed in our DNS. This would enable a better assessment of the extent to which the flow at these scales is similar to the idealized case of isotropic turbulence, and the extent to which large deviations about the mean-field behaviour allow the effects of stratification to be felt even in this range. A flow with a large scale separation between the integral length and the Ozmidov scale would also allow for a better understanding of the way TKE and TPE are exchanged in this regime, and the mechanisms governing the TKE and TPE interscale fluxes in this regime, where buoyancy plays a key role. Exploring such flows is, however, very challenging computationally, and it may be some time before a flow in such a regime can be simulated with DNS.
Acknowledgements
Specifically, the Comet and Expanse clusters were used under allocation CTS170009. The DNS database was generated with funding from the US Office of Naval Research via grant N00014-15-1-2248. High performance computing resources for the DNS were provided through the US Department of Defense High Performance Computing Modernization Program by the Army Engineer Research and Development Center and the Army Research Laboratory under Frontier Project FP-CFD-FY14-007.
The contributions by G.D.P. were performed under the auspices of the US Department of Energy by Lawrence Livermore National Laboratory under contract DE-AC52-07NA27344.
Funding
A.D.B. and X.Z. acknowledge support through a National Science Foundation (NSF) CAREER award no. 2042346. This work used the Extreme Science and Engineering Discovery Environment (XSEDE), which is supported by NSF grant no. ACI-1548562 (Towns et al. Reference Towns2014).
Declaration of interests
The authors report no conflict of interest.
Open access
