Previous studies show that at the small scales of stably stratified turbulence, the scale-dependent buoyancy flux reverses sign, corresponding to a conversion of turbulent potential energy (TPE) back into turbulent kinetic energy (TKE). Moreover, the magnitude of the reverse flux becomes stronger with increasing Prandtl number $\textit{Pr}$. Using a filtering analysis we demonstrate analytically how this flux reversal is connected to the mechanism identified in Bragg & de Bruyn Kops (2024 J. Fluid Mech. vol. 991, A10) that is responsible for the surprising observation that the TKE dissipation rate increases while the TPE dissipation rate decreases with increasing $\textit{Pr}$ in stratified turbulence. The mechanism identified by Bragg & de Bruyn Kops, which is connected to the formation of ramp–cliff structures in the density field, is shown to give the scale-local contribution to the buoyancy flux. At the smallest scales this local contribution dominates and explains the flux reversal, while at larger scales a non-local contribution is important. Direct numerical simulations of three-dimensional statistically stationary, strongly stably stratified turbulence confirm the theoretical analysis, and indicate that, while on average the local contribution only dominates the buoyancy flux at the smallest scales, it remains strongly correlated with the buoyancy flux at all scales. The results show that ramp cliffs are not only connected to the reversal of the local buoyancy flux but also the non-local part. At the small scales (approximately below the Ozmidov scale), ramp structures contribute exclusively to reverse buoyancy flux events, whereas cliff structures contribute to both forward and reverse buoyancy flux events.

We examine the dynamic interactions between the large-scale coherent motion and the small-scale turbulence in the passive scalar field of a circular cylinder wake, where the coherent motion exhibits strong periodicity. A combination of four X-wires and four cold wires was used to simultaneously measure the three velocity and temperature fluctuations at nominally the same location. Measurements were taken at $x/d=10$, 20 and 40 in the mean shear plane at Reynolds number 2500, based on the cylinder diameter $d$ and the free-stream velocity. The phase-averaging technique is used to distinguish the large-scale coherent motion from the stochastic motion, enabling the construction of phase-averaged structure functions of the passive scalar in the scale phase plane. The maximum of the coherent scalar $\tilde {\theta }$ closely aligns with the minima of the phase-averaged strain $\langle S \rangle$ and the vortex centre, suggesting that heat is contained within the interior of the vortex. The scale-by-scale distributions of the scalar variance and the streamwise velocity variance exhibit a similar phase dependence associated with the coherent motion. This dependence is perceptible even at the smallest scales. However, as the distance from the cylinder increases, the perceivable scale range decreases and eventually disappears. An expression is formulated to describe the time-averaged second-order structure function of coherent scalar and the time-averaged second-order mixed structure function between the coherent scalar and coherent streamwise velocity at $x/d= 10$ and 20, where the coherent motion is prominent. Furthermore, the scale-by-scale contribution of the coherent scalar variance to the total scalar variance is evaluated. Also, we derive the scale-by-scale scalar variance transport equations that account for the coherent motion in both general and isotropic formulations. Assuming local isotropy, it is found that the equation agrees approximately with the experimental data across all scales at $x/d= 40$. Finally, the differences between the scale-by-scale transport equation for the stochastic scalar variance and that for the stochastic turbulent kinetic energy are discussed.

Integral modelling of turbulent buoyant plumes is crucial for rapid predictions of plume characteristics. While the governing equations are typically derived using self-similarity and a Boussinesq approximation, these assumptions may not hold for plumes originating from finite-area sources with large density ratios. This work evaluates the accuracy of integral-scale models for non-Boussinesq lazy plumes using high-fidelity numerical simulations of turbulent helium plumes. We analyse the plume kinematics by computing vertical fluxes, plume radius and radial profiles, establishing some disparities between common practice and physical accuracy. We identify how the definition of the plume radius changes the perception of the plume structure when the flow is not self-similar and derive a relationship between the flux-based and threshold-based definitions without requiring self-similarity. We then examine the plume dynamics by evaluating the source terms from the governing plume equations. Our results support neglecting diffusive and viscous effects but emphasise the importance of the mean pressure gradient, even in the self-similar regime. Two coefficients need to be modelled: the well-known entrainment coefficient and the lesser-known momentum correction coefficient, which is a correction required for the momentum equation to account for self-similar and slender approximations. The momentum correction coefficient is found to be approximately constant and slightly greater than the assumed value of 1. The standard entrainment coefficient models perform well up to a local Richardson number three times the asymptotic value but overpredict entrainment for larger Richardson numbers. We propose a correction using the known finite limit of entrainment at infinite Richardson number.

The spatial organisation of a passive scalar plume originating from a point source in a turbulent boundary layer is studied to understand its meandering characteristics. We focus shortly downstream of the isokinetic injection ($1.5\leqslant x/\delta \leqslant 3$, $\delta$ being the boundary-layer thickness) where the scalar concentration is highly intermittent, the plume rapidly meanders and breaks up into concentrated scalar pockets due to the action of turbulent structures. Two injection locations were considered: the centre of the logarithmic region and the wake region of the boundary layer. Simultaneous quantitative acetone planar laser-induced fluorescence and particle image velocimetry were performed in a wind tunnel, to measure scalar mixture fraction and velocity fields. Single- and multi-point statistics were compared with established works to validate the diagnostic novelties. Additionally, the spatial characteristics of plume intermittency were quantified using ‘blob’ size, shape, orientation and mean concentration. It was observed that straining, breakup and spatial reorganisation were the primary plume-evolution modes in this region, with little small-scale homogenisation. Further, the dominant role of coherent vortex motions in plume meandering and breakup was evident. Their action is found to be the primary mechanism by which the injected scalar is transported away from the wall in high concentrations (‘large meander events’). Strong spatial correlation was observed in both instantaneous and conditional fields between the high-concentration regions and individual vortex heads. This coherent transport was weaker for wake injection, where the plume only interacts with outer vortex motions. A coherent-structure-based mechanism is suggested to explain these transport mechanisms.

Using pore-resolved direct numerical simulation (DNS), we investigate passive scalar transport at a unit Schmidt number in a turbulent flow over a randomly packed bed of spheres. The scalar is introduced at the domain’s free-slip top boundary and absorbed by the bed, which maintains a constant and uniform scalar value on the sphere surfaces. Eight cases are analysed, which are characterised by friction Reynolds numbers of ${\textit{Re}}_\tau \in [150, 500]$ and permeability Reynolds numbers of ${\textit{Re}}_{{\kern-1pt}K} \in [0.4, 2.8]$, while flow depth-to-sphere-diameter ratios vary within $h/D \in \{ 3, 5, 10 \}$ and the roughness Reynolds numbers lie within $k_s^+ \in [20,200]$. For cases with comparable ${\textit{Re}}_\tau$, the permeable wall behaviour enhances scalar absorption, as indicated by increases in the Sherwood number and the scalar roughness function $\Delta c^+$ over ${\textit{Re}}_{{\kern-1pt}K}$. At progressively higher ${\textit{Re}}_{{\kern-1pt}K}$, the scalar absorption diverges increasingly from the momentum absorption, as its distribution peaks deeper below the crests of the sphere pack and spreads over a wider vertical region. The fixed-value scalar boundary condition emphasises certain similarities between the scalar and velocity fields. Near-interface scalar fluctuations are correlated with streamwise velocity fluctuations, and the turbulent Schmidt number remains close to its value in the free-flow region. Compared with zero-flux scalar boundary conditions, prescribing a uniform scalar value on the sphere surfaces reduces spatial heterogeneity within the pore space, thereby limiting both dispersive transport and the form-induced production of temporal scalar fluctuations.

Characteristics of the turbulent/non-turbulent interface (TNTI) and entrainment in separated and reattaching flows induced by an oscillating fence are investigated using time-resolved particle image velocimetry. Disturbed flows are classified into subcritical, transitional, critical and supercritical cases based on the ratio of the oscillation frequency to the natural vortex shedding frequency. In the recirculation zone, distinct vortices across different cases lead to significant variations in TNTI characteristics. In the subcritical case, the TNTI evolution resembles that in the stationary fence case but with intensified height fluctuations due to the undulation of separated shear layer. For other cases, the mean TNTI height increases with the oscillation frequency, while height fluctuation diminishes. The TNTI thickness varies with nearby vortices, scaling with the Taylor microscale. After the reattachment, TNTI height distributions converge into two groups: subcritical and transitional cases exhibit larger fluctuations and positively skewed probability density functions (PDFs), while critical and supercritical cases show smaller fluctuations and basically symmetric PDFs. The TNTI thickness becomes consistent across various cases, matching the adjacent small-scale vortex size. Besides, the nibbling mechanism of entrainment aligns well with the flow development. The minimum mean entrainment velocity coincides with the strongest prograde vortex while the maximum occurs at $x\approx 1.2x_{{r}}$ (where $x$ denotes the streamwise coordinate and $x_{{r}}$ is the mean reattachment position) in all cases. Engulfment is enhanced near the reattachment location by oscillations in the transitional and critical cases, but is suppressed in the supercritical cases due to the weakness of vortex structures at higher oscillation frequencies.

An imposed constant magnetic field parallel to the interface in the Rayleigh–Taylor framework strongly modifies the dynamics of the flow. The growth rate of the turbulent mixing layer is almost doubled compared with the purely hydrodynamic case, mainly due to a strong reduction of small-scale mixing. Indeed, magnetic tension inhibits the small-scale perturbations from developing, which in turn creates a strong anisotropy with structures elongated in the field direction. Two theoretical predictions for the asymptotic state of the magnetic Rayleigh–Taylor instability (MRTI) are put forward. First, considering the large-scale dynamics, an upper bound for the mixing layer growth rate is obtained. Second, the phenomenology is embedded in a buoyancy–drag equation from which an analytical relation between the growth rate, mixing, anisotropy and induced magnetic fields is derived. Both predictions are successfully assessed with high resolution direct numerical simulations of the Boussinesq–Navier–Stokes equations under the magnetohydrodynamics approximation. These predictions characterize the quasi-self-similar state of the MRTI driven by strong magnetic fields.

This study investigates the influence of free-stream turbulence (FST) and the thrust coefficient ($C_T$) on wind turbine wakes. Wakes generated at $C_T \in \{0.5, 0.7,0.9\}$ are exposed to turbulent inflows with varying FST intensities ($1\,\% \lesssim {\textit{TI}}_{\infty } \lesssim 11\,\%$) and integral length scales ($0.1 \lesssim {\mathcal L}_x/\!D \lesssim 2$, $D$ is the rotor diameter). For high-${\textit{TI}}_{\infty }$ inflows, a flow region in the wake is observed where a mean momentum deficit persists despite the turbulence intensity having already homogenised with that of the free stream, challenging traditional wake definitions. A ‘turning point’ in the mean wake width evolution is identified, beyond which wakes spread at slower rates. Near-field ($x\!/\!D \lesssim 7$) wake growth rate increases with higher ${\textit{TI}}_{\infty }$ and $C_T$, while far-field ($x\!/\!D \gtrsim 15$) wake growth rate decreases with higher ${\textit{TI}}_{\infty }$ – a finding with profound implications for wind turbine wake modelling that also aligns with the entrainment behaviours observed in bluff- and porous-body wakes exposed to FST. Increasing ${\mathcal L}_x$ delays wake recovery onset and reduces the mean wake width, with minimal effect on the spreading rate. Both $C_T$ and FST influence the high- and low-frequency wake dynamics, with varying contributions in the near and far fields. For low-${\textit{TI}}_{\infty }$ and small-${\mathcal L}_x$ inflows, wake meandering is minimal, sensitive to $C_T$ and appears to be triggered by a shear-layer instability. Wake meandering is enhanced for high-${\textit{TI}}_{\infty }$ and large-${\mathcal L}_x$ inflows, with the integral length scale playing a leading role. This emphasises the complex role of FST integral length scale: while increasing ${\mathcal L}_x$ amplifies meandering, it does not necessarily translate to larger mean wake width due to the concurrent suppression of entrainment rate.

A partial differential equation governing the evolution of the joint probability distribution of multicomponent flow observations, drawn randomly from one or more control volumes, is derived and applied to examples involving irreversible mixing. Unlike local probability density methods, this work adopts an integral perspective by regarding a control volume as a sample space with an associated probability distribution. A natural and general definition for the boundary of such control volumes comes from the magnitude of the gradient of the sample space distribution, which can accommodate Eulerian or Lagrangian frames of reference as particular cases. The formulation exposes contributions made by uncertain or stochastic boundary fluxes and internal cross-gradient mixing in the equation governing the observables’ joint probability distribution. Advection and diffusion over a control volume’s boundary result in source and drift terms, respectively, whereas internal mixing, in general, corresponds to the sign-indefinite diffusion of probability density. Several typical circumstances for which the corresponding diffusion coefficient is negative semidefinite are identified and discussed in detail. The framework is a natural setting for examining available potential energy, the incorporation of uncertainty into bulk models, and establishing a link with the Feynman–Kac formula and Kolmogorov equations that are used to analyse stochastic processes.

We explore the fundamental flow structure of temporally evolving inclined gravity currents with direct numerical simulations. A velocity maximum naturally divides the current into inner and outer shear layers, which are weakly coupled by momentum and buoyancy exchanges on time scales that are much longer than the typical time scale characterising either layer. The outer layer evolves to a self-similar state and can be described by theory developed for a current on a free-slip slope (Van Reeuwijk et al. 2019, J. Fluid Mech., vol. 873, pp. 786–815) when expressed in terms of outer-layer properties. The inner layer evolves to a quasi-steady state and is essentially unstratified for shallow slopes, with flow statistics that are virtually indistinguishable from fully developed open channel flow. We present the classic buoyancy–drag force balance proposed by Ellison & Turner (1959, J. Fluid Mech., vol. 6, pp. 423–448) for each layer, and find that buoyancy forces in the outer layer balance entrainment drag, while buoyancy forces in the inner layer balance wall friction drag. Using scaling laws within each layer and a matching condition at the velocity maximum, the entire flow system can be solved as a function of the slope angle, in good agreement with the simulation data. We further derive an entrainment law from the solution, which exhibits relatively high accuracy across a wide range of Richardson numbers, and provides new insights into the long runout of oceanographic gravity currents on mild slopes.

Mixing and heat transfer rates are typically enhanced in high-pressure transcritical turbulent flow regimes. This is largely due to the rapid variation of thermophysical properties near the pseudo-boiling region, which can significantly amplify velocity fluctuations and promote flow destabilisation. The stability conditions are influenced by the presence of baroclinic torque, primarily driven by steep, localised density gradients across the pseudo-boiling line; an effect intensified by differentially heated wall boundaries. As a result, enstrophy levels increase compared with equivalent low-pressure systems, and flow dynamics diverge from those of classical wall-bounded turbulence. In this study the dynamic equilibrium of these instabilities is systematically analysed using linear stability theory. It is shown that under isothermal wall transcritical conditions, the nonlinear thermodynamics near the pseudo-boiling region favour destabilisation more readily than in subcritical or supercritical states; though this typically requires high-Mach-number regimes. The destabilisation is further intensified in non-isothermal wall configurations, even at low Brinkman and significantly low Mach numbers. In particular, the sensitivity of neutral curves to Brinkman number variations, along with the modal and non-modal perturbation profiles of hydrodynamic and thermodynamic modes, offer preliminary insight into the conditions driving early destabilisation. Notably, a non-isothermal set-up (where walls are held at different temperatures) is found to be a necessary condition for triggering destabilisation in low-Mach, low-Reynolds-number regimes. For the same Brinkman number, such configurations accelerate destabilisation and enhance algebraic growth compared with isothermal wall cases. As a consequence, high-pressure transcritical flows exhibit increased kinetic energy budgets, driven by elevated production rates and reduced viscous dissipation.

We consider numerically a Lagrangian view of turbulent mixing in two-layer stably stratified parallel shear flow. By varying the ratio of shear layer depth to density interface thickness, these flows are prone to either a primary Kelvin–Helmholtz instability (KHI) or to a primary Holmboe wave instability (HWI). These instabilities are conventionally thought to mix qualitatively differently; by vortical ‘overturning’ of the density interface induced by KHI, or by turbulent ‘scouring’ on the edges of the density interface induced by HWI. By tracking Lagrangian particles in direct numerical simulations, so that the fluid buoyancy sampled along particle paths provides a particular Lagrangian measure of mixing, we investigate the validity of this overturning/scouring classification. The timing of mixing events experienced by particles inside and outside the interface is qualitatively different in simulations exhibiting KHI and HWI. The root mean square (r.m.s.) buoyancy for particles that start with the same buoyancy is actually larger for HWI-associated flows than for KHI-associated flows for the same bulk Richardson number $Ri_b$, implying heterogeneous mixing along particle paths for HWI. The number of particles starting close to the mid-plane of the interface which experience a change in sign in the local fluid buoyancy (and hence end up on the opposite side of the mid-plane after mixing) is compared for KHI and HWI in flows with various $Ri_b$. Perhaps surprisingly, for HWI with a large $Ri_b$, more than half of the particles that start near the mid-plane end up on the opposite side of the mid-plane.

We investigate the energetics of mixing induced by a continuously supplied dense current (density $\rho _0$) propagating beneath a lighter ambient fluid (density $\rho _a$) along a horizontal rigid boundary within a rectangular domain. The flow fields are computed using direct numerical simulations (DNS) performed with the Nek5000 spectral element solver. Mixing is quantified through the temporal evolution of the background potential energy, which exhibits a linear increase over time. This linear trend enables the definition of a dimensionless mixing parameter $\gamma$, representing the rate of background potential energy growth. The value of $\gamma$ depends on the initial density contrast for a fixed volumetric discharge at the source, characterised by the dimensionless source Froude number. The results reveal a non-monotonic dependence of $\gamma$ on the source Froude number, highlighting a complex interaction between flow forcing and mixing efficiency. We find that, under the assumption of uniform mixing along the current’s length, a fraction $\gamma /2$ of the total supplied energy is invested in mixing along a horizontal distance equal to the height of the inlet.

Space–time correlations of velocity and high-Schmidt-number ($Sc \approx 2000$) passive scalar fields are investigated in turbulent pipe flow using particle image velocimetry and planar laser-induced fluorescence, respectively. Both the velocity and scalar fields exhibit characteristic elliptical patterns in their respective space–time correlations. The elliptic approximation model, originally developed for the velocity field, is applied to estimate convection and sweeping velocities for both fields. In both fields, the convection velocity decreases, while the sweeping velocity increases, along the pipe radius. The convection velocity ratio between the scalar and velocity fields shows that high-Schmidt-number scalar fluctuations are advected faster than the velocity fluctuations. Similarly, the sweeping velocity of the scalar fluctuations is found to be larger than that of the velocity fluctuations. Furthermore, the high-Schmidt-number scalar is found to decorrelate more rapidly than the corresponding velocity, with the scalar Taylor microscale distinctly smaller than the velocity Taylor microscale.

A Lagrangian description of bubble swarms has largely eluded both experimental and numerical efforts. Now, in a tour de force of deep-learning-enabled optical tracking measurements, Huang et al. (2025 J. Fluid. Mech. 1014, R1) have managed to follow the three-dimensional trajectories of $10^5$ deforming and overlapping bubbles within a swarm, perhaps for long enough to witness their approach to the diffusive limit. Their results reveal that bubble swarms exhibit a dispersion law strikingly reminiscent of classical Taylor dispersion in isotropic turbulence, but with an earlier, undulatory transition from the ballistic-to-diffusive regime. Huang et al. (2025 J. Fluid Mech. 1014, R1), have helped close the loop on our understanding of Lagrangian bubble dispersion – from self-stirring swarms to bubbles in isotropic turbulence.

In typical nature and engineering scenarios, such as supernova explosion and inertial confinement fusion, mixing flows induced by hydrodynamic interfacial instabilities are essentially compressible. Despite their significance, accurate predictive tools for these compressible flows remain scarce. For engineering applications, the Reynolds-averaged Navier–Stokes (RANS) simulation stands out as the most practical approach due to its outstanding computational efficiency. However, existing RANS studies focus primarily on cases where the compressible effect plays an insignificant role in mixing development, with quite limited attention given to scenarios with significant compressibility influence. Moreover, most of the existing RANS mixing models demonstrate significantly inaccurate predictions for the latter. This study develops a novel compressible RANS mixing model by incorporating physical compressibility corrections into the $K$–$L$–$\gamma$ mixing transition model recently proposed by Xie et al. (J. Fluid Mech. 1002, 2025, A31). Specifically, taking the density-stratified Rayleigh–Taylor mixing flows as representative compressible cases, we first analyse the limitations of the existing model for compressible flows, based on high-fidelity data and local instability criteria. Subsequently, the equation of state for a perfect gas is employed to derive comprehensive compressibility corrections. The crucial turbulent composition and heat fluxes are integrated into the closure of the key turbulent mass flux term of the turbulent kinetic energy equation. These corrections enable the model to accurately depict compressible mixing flows. Systematic validations confirm the efficacy of the proposed modelling scheme. This study offers a promising strategy for modelling compressible mixing flows, paving the way for more accurate predictions in complex scenarios.

Understanding the interplay between buoyancy and fluid motions within stably stratified shear layers is crucial for unravelling the contribution of flow structures to turbulent mixing. In this study, we examine statistically the local relationship between stratification and fluid deformation rate in wave and turbulent regimes, using experimental datasets obtained from a stratified inclined duct (SID) containing fluids of different densities that form an exchange flow. We introduce rotational and shear components of varying strength within the vorticity and a family of coherent gradient Richardson numbers ($Ri_C$), ratios related to the buoyancy frequency and the strength of either the rotational or shearing motion. Conditional statistical analysis reveals that both shear and stratification intensity affect the probability distribution of the $Ri_C$, with extreme events occurring more frequently in areas of weak stratification. In the wave regime, we identify the persistence of fast-spin vortices within the strongly stratified density interface. However, scouring of the density interface is primarily driven by shearing motions, with baroclinic torque making a notable contribution to enstrophy transport. In the turbulent regime, rigid-body rotations occur at significantly shorter time scales than that associated with the local buoyancy frequency, making them more disruptive to stratification than shear. Additionally, correlation analysis reveals that irrotational strain distorts stable stratification similarly to shearing motions, but is weaker than both shearing and rotational motions and less likely to have a time scale longer than that related to the buoyancy frequency. Moreover, we observed that the interplay between rotational and shearing motions intensifies as stratification increases. Finally, a comparison of length scales along the shear layers highlights the $Ri_C$ as a valuable measure of the relative sizes of different motions compared with the Ozmidov scale and shows that stratification can influence sub-Ozmidov scales through baroclinic torque. This study highlights the critical impact of the type, strength and location of fluid deformations on localised mixing, providing new insights into the role of rotational motions in shear-driven stratified flows.

The turbulent transport of momentum, energy and passive scalar is investigated in the flow around a rectangular cylinder of aspect ratio 5 : 1 – a geometry representative of separating and reattaching flows from sharp-edged bodies. The study is based on direct numerical simulation (DNS) conducted at Reynolds numbers up to ${\textit{Re}} = 14\,000$, based on the cylinder thickness, with Schmidt number fixed at ${\textit{Sc}} = 0.71$. At this Reynolds number, the flow exhibits features of asymptotic high-${\textit{Re}}$ behaviour. Budgets of mean momentum, Reynolds stresses, mean scalar and scalar fluxes provide a detailed view of the underlying transport mechanisms. The mean momentum balance elucidates the role of turbulence in entraining free stream fluid, promoting shear-layer reattachment, sustaining backflow in the recirculation region and regulating wake dynamics through large-scale vortex shedding. The leading-edge shear layer is the main site of turbulence production, with energy injected into streamwise fluctuations and redistributed to cross-flow components by pressure–strain interactions. As ${\textit{Re}}$ increases, vertical fluctuations increasingly return energy to the mean upward flow, stabilising the separation bubble height. Turbulent transport dominates scalar redistribution. Scalar fluxes are primarily generated by interactions between Reynolds stresses and scalar gradient, and modulated by pressure-scalar gradient effects. An a priori evaluation of eddy-viscosity and diffusivity models quantifies the misalignment between modelled and DNS-resolved stress and flux tensors, as well as the inhomogeneity of eddy transport coefficients. This analysis deepens the understanding of transport phenomena in bluff-body flows approaching the asymptotic regime, and underpins the validation and improvement of turbulence models for separating and reattaching flows.

We analyse direct numerical simulations of homogeneous, forced, stably stratified turbulence to study how the pressure–strain and pressure scrambling terms are modified as stability is increased from near neutral to strongly stratified conditions. We decompose the pressure into nonlinear and buoyancy components, and find that the buoyancy part of the pressure–strain correlation changes sign to promote large-scale anisotropy at strong stability, unlike the nonlinear component, which always promotes large-scale isotropy. The buoyancy component of the pressure scrambling term is positive semidefinite and increases monotonically with stability. As its magnitude becomes greater than the nonlinear component (which is negative), the overall scrambling term generates buoyancy flux at very strong stability. We apply quadrant analysis (in the pressure-gradient space) to these correlations to study how contributions from the four quadrants change with stability. Furthermore, we derive exact relationships for the volume-averaged buoyancy components of these correlations which reveal (i) the buoyancy component of the pressure–strain correlation involves a weighted sum of the vertical buoyancy flux cospectrum, so counter-gradient buoyancy fluxes contribute to enhanced anisotropy by transferring vertical kinetic energy into horizontal kinetic energy; (ii) the buoyancy component of the pressure scrambling term involves a weighted sum of the potential energy spectrum; (iii) the weighting factor accentuates contributions from layered motions, which are a prominent feature of strongly stratified flows. These expressions apply generally to all homogeneous stratified flows independent of forcing and across all stability conditions, explaining why these effects have been observed for both forced and sheared stably stratified turbulence simulations.

Layer formation can occur within stratified fluids, often associated with the effect of ‘double diffusion’ where the fluid buoyancy depends on two components with differing molecular diffusivities (e.g. heat and salt in seawater). However, since layering also occurs in one-component stratified fluids, the generation mechanism for layers is often unclear. In this paper, we present a framework that unifies multiple-layer generation mechanisms across both one- and two-component stratified fluids. We demonstrate how these mechanisms can be assessed using simulations of double-diffusive intrusions. Our simulations illustrate the importance of the negative turbulent diffusivity for buoyancy in contributing to layer formation.