Evolution of a Stratified Turbulent Cloud under Rotation

1. Introduction

Localized turbulence is common in both atmospheric and oceanic flows, which are strongly influenced by the Earth’s rotation and stable density stratification [1,2,3,4,5]. Understanding the evolution of these turbulent patches under the combined effects of rotation and stratification is crucial for gaining deeper insights into various geophysical phenomena. A key mechanism behind the generation of such turbulence is the local breakdown of internal waves [6,7]. In particular, the resulting turbulence often exhibits horizontal scales that are significantly larger than their vertical counterparts. Empirical observations from deep waters and the equatorial pycnocline also reveal pronounced differences in horizontal and vertical scales of turbulent patches [2,3]. Furthermore, mixing induced by internal wave breaking tends to weaken local stratification. Such reduced stratification is also observed in turbulent regions within the ocean mixed layer [8]. With these characteristics in mind, our present study is concerned with a horizontally extended turbulent cloud within a vertically stably stratified, rotating fluid system, where the effects of rotation dominate over the stratification. Our focus is on understanding the formation of flow structures and analysing the associated energy dispersion patterns.

While numerous studies have highlighted the dominance of waves in spatially uniform and stationary scenarios of rotating/stratified turbulence [9,10,11], the dynamics of localized turbulent patches has also attracted considerable attention. In particular, under rotation, these turbulent patches often give rise to columnar vortices parallel to the axis of rotation, associated with inertial waves [12,13,14,15]. For example, in their oscillating-grid turbulence experiments, Dickinson and Long [13] found that, for a Rossby number close to unity, the boundary between perturbed and quiescent regions moves at speeds related to certain wave velocities. Similarly, Davidson et al. [14] showed experimentally that the leading edges of columnar vortices under rotation align with inertial wave group velocities. Using direct numerical simulations (DNS), Ranjan and Davidson [15] observed these columnar structures, attributed their formation to inertial waves, and emphasized their role in significant energy dissipation from localized turbulence.

Analogous to turbulent patches under rotation, the presence of internal wave radiation is clearly evident in localized turbulence within stratified fluids [16,17,18,19], such as wakes caused by object motion [20,21]. The flow structures in these scenarios are predominantly wave-dominated. In addition, these waves contribute significantly to energy reorganization, highlighting their fundamental role in the atmosphere and oceans [22]. Anisotropic structures and wave propagation are observed in both rotating and stratified turbulent flows. These phenomena highlight the complex interplay between rotation and stratification in geophysical flows [23,24,25]. While many studies have investigated the turbulent patch behavior in rotating and stratified fluids, the emphasis has often been on aspects such as lateral intrusions, Thorpe scales and the vertical expansion of the patch [26,27,28,29]. However, the waves emitted by these disturbances often receive less attention. One notable exception [30] used DNS to study a buoyant cloud subjected to vertical rotation and horizontal gravity, reminiscent of equatorial conditions. It was confirmed that the formation of columnar structures was due to inertial waves arising from the buoyant cloud.

Building on previous work, we explore an under-explored area: the evolution of a horizontal turbulent cloud under vertical rotation and stratification. This has significant relevance in geophysical contexts such as the high latitudes of the atmosphere and oceans [31] and the deep Arctic waters [32,33]. In particular, our study examines scenarios where rotation effects outweigh stratification, in contrast to many typical geophysical flows where stratification predominates. In this context, we address three key questions: First, how does stratification affect the resulting flow structures, as opposed to purely rotating scenarios? Second, are these structures influenced by inertial-gravity waves, and how does different stratification affect their formation? Third, what fraction of the cloud’s energy is transferred to waves at different stratification levels?

2. The evolution of a single stratified eddy under rotation

In this section, we study the evolution of a single stratified eddy under rotation. This is motivated by the fact that a turbulent cloud can be considered as a sea of randomly oriented vortex blobs with different scales. First, we introduce internal-inertial gravity waves under the Boussinesq approximation. Second, we solve analytically the initial value problem of a compact vortex blob in a rapidly rotating stratified environment. Finally, a numerical simulation is performed to validate the analytical solution. It is shown that as the eddy evolves, part of the energy is retained in the eddy by the potential vorticity (PV) mode, while the other part is carried away from the eddy by inertial-gravity waves.

2.1. Internal-inertial gravity waves

The Boussinesq set of equations for a linearly stratified fluid under system rotation can be written as:

$$\nabla ·\mathbf{u}=0,$$

(1)

$$\frac{\partial \mathbf{u}}{\partial t}+\mathbf{u}·\nabla \mathbf{u}=-\frac{1}{{\rho }_0}\nabla p+2\mathbf{u}\times \Omega -N\varphi {\mathbf{e}}_z+\nu {\nabla }^2\mathbf{u},$$

(2)

$$\frac{\partial \varphi }{\partial t}+\mathbf{u}·\nabla \varphi =N{u}_z+\kappa {\nabla }^2\varphi ,$$

(3)

where $\mathbf{u}$ is the velocity vector, $u_z$ is the z-component of the velocity, p is the modified pressure incorporating a centrifugal term, $\Omega =\Omega {\mathbf{e}}_z$ represents the rotation vector, N is the Brunt-Väisälä frequency, and $\nu$ and $\kappa$ are the kinematic viscosity and the diffusion coefficient, respectively. We define $\varphi =\left(g/{\rho }_{0}N\right){\rho }^{\prime }$, which has the dimension of the velocity, and ${\rho }^{\prime }$ is the density perturbation from the ambient density.

Inertial-gravity waves can be obtained after a linearization of the Boussinesq equations. Neglecting the molecular diffusion and the second-order terms of $\mathbf{u}$, $\mathbf{w}$ and $\varphi$ one obtains:

$$\frac{\partial \omega }{\partial t}+\nabla \times (2\Omega \times \mathbf{u})=-N\nabla \varphi \times {\mathbf{e}}_z,$$

(4)

$$\frac{\partial \varphi }{\partial t}=N{u}_z,$$

(5)

from which one finds [see [34], p. 56]

$$\frac{{\partial }^2}{\partial t^2}{\nabla }^2{u}_z+4{\Omega }^2\frac{{\partial }^2{u}_z}{\partial z^2}+N^2\left(\frac{{\partial }^2{u}_z}{\partial x^2}+\frac{{\partial }^2{u}_z}{\partial y^2}\right)=0.$$

(6)

This equation admits inertial-gravity waves with the dispersion relation given as

$$\varpi =\pm \sqrt{N^2{k}_h^2+4{\Omega }^2{k}_z^2}/k,$$

(7)

where $\varpi$ is the frequency, while $k_h=\sqrt{k_x^2+k_y^2}$ and $k=\sqrt{k_x^2+k_y^2+k_z^2}$ are the horizontal wavenumber and the total wavenumber, respectively. The group velocity is ${\mathbf{c}}_g=\nabla \varpi$, whose z-component is

$$c_{g,z}=\frac{\partial \varpi }{\partial k_z}=\pm \frac{\left(4{\Omega }^2-N^2\right)k_h^2{k}_z}{\sqrt{N^2{k}_h^2+4{\Omega }^2{k}_z^2}{k}^3}.$$

(8)

2.2. Analytical study of a single eddy

Following Davidson et al. [14], we consider the evolution of an eddy in a stratified fluid under rotation, assuming that the process is axisymmetric with respect to the eddy axis. In this context, in a cylindrical coordinate system $(r,\theta ,z)$, an axisymmetric velocity field $\mathbf{u}$ can be decomposed into azimuthal and poloidal components

$$\mathbf{u}=\left(\Gamma /r\right){\mathbf{e}}_{\theta }+\nabla \times \left[\left(\psi /r\right){\mathbf{e}}_{\theta }\right],$$

(9)

where $\Gamma =u_{\theta }r$ is the angular momentum and $\psi$ is the Stokes streamfunction. We substitute Equation (9) into the linearized inviscid vorticity equation (4) and the linearized non-diffusive Equation (5) to obtain

$$\frac{\partial \Gamma }{\partial t}=2\Omega \frac{\partial \psi }{\partial z},\frac{\partial }{\partial t}\left(r{\omega }_{\theta }\right)=2\Omega \frac{\partial \Gamma }{\partial z}+Nr\frac{\partial \varphi }{\partial r},\frac{\partial \varphi }{\partial t}=\frac{N}{r}\frac{\partial \psi }{\partial r},$$

(10)

where ${\nabla }_{*}^2\psi =\left(r\partial /\partial r\right)\left(r^{-1}\partial \psi /\partial r\right)+{\partial }^2\psi /\partial z^2=-r{\omega }_{\theta }$, ${\nabla }_{*}^2$ being the Stokes operator. Combining the equations in (10) yields

$$\frac{\partial }{\partial t}\left(\frac{{\partial }^2}{\partial t^2}{\nabla }_{*}^2\Gamma +\left[{\left(2\Omega \right)}^2\frac{{\partial }^2}{\partial z^2}+N^2\left(r\frac{\partial }{\partial r}\right)\left(\frac{\partial }{r\partial r}\right)\right]\Gamma \right)=0.$$

(11)

Equation (11) can be readily solved by using the Hankel-cosine transform

$${\widehat{u}}_{\theta }=\frac{1}{2{\pi }^2}{\int }_0^{\infty }{\int }_0^{\infty }r{u}_{\theta }{J}_1\left(k_{r}r\right)cos\left(k_{z}z\right)\mathrm{d}r\mathrm{d}z,$$

(12)

where $J_1$ is the Bessel function of the first kind of order 1. Here, $k_r$ and $k_z$ represent the wavenumbers in the r-direction and the z-direction, respectively, with $0\le k_r,k_z<\infty$. Given the initial conditions ${\widehat{u}}_{\theta }={\widehat{u}}_{\theta }^{\left(0\right)}$, $\psi =0$ and $\varphi =0$, we can derive that

$${\widehat{u}}_{\theta }=\left[{\widehat{u}}_{\theta }^{\left(0\right)}/\left(N^2{k}_r^2+{\left(2\Omega \right)}^2{k}_z^2\right)\right]\left[N^2{k}_r^2+{\left(2\Omega \right)}^2{k}_z^{2}cos\left(\sqrt{N^2{k}_r^2+{\left(2\Omega \right)}^2{k}_z^2}t/k\right)\right],$$

(13)

from which we obtain

$$\begin{array}{cc}\hfill u_{\theta }=& \underset{\mathrm{PV}\mathrm{mode}}{\underbrace{4\pi {\int }_0^{\infty }{\int }_0^{\infty }\left(N^2{k}_r^3{\widehat{u}}_{\theta }^{\left(0\right)}/{\alpha }^2\right)J_1\left(k_{r}r\right)cos\left(k_{z}z\right)\mathrm{d}{k}_r\mathrm{d}{k}_z}}\hfill \\ \hfill & +\underset{\mathrm{inertial}-\mathrm{gravity}\mathrm{waves}}{\underbrace{2\pi {\int }_0^{\infty }{\int }_0^{\infty }\left[{\left(2\Omega \right)}^2{k}_r{k}_z^2{\widehat{u}}_{\theta }^{\left(0\right)}/{\alpha }^2\right]J_1\left(k_{r}r\right)\left[cos\left(k_{z}z-\alpha t/k\right)+cos\left(k_{z}z+\alpha t/k\right)\right]\mathrm{d}{k}_r\mathrm{d}{k}_z}},\hfill \end{array}$$

(14)

where we introduce $\alpha =\sqrt{N^2{k}_r^2+{\left(2\Omega \right)}^2{k}_z^2}$ for conciseness. In the absence of stratification (i.e. when $N=0$), Equations (13) and (14) are consistent with the results presented in Davidson et al. [14], which consider only the rotation of the system. However, the inclusion of stratification introduces a new time-independent term and changes the frequencies of the previously identified wave-like components. Note that the phase velocities of the wave-like components, given by $\pm \alpha /k$, closely match those of the inertial-gravity waves according to equation (7), provided we replace $k_r$ with $k_h$. This alignment strongly suggests that these wave-like components are indeed manifestations of the inertial-gravity waves. In a periodic domain with spatially uniform rotating stratified flows, Smith and Waleffe [35] identified the linear eigenmodes as two inertial-gravity waves and a PV mode, suggesting that the time-independent term in (14) is consistent with the zero-frequency PV mode. However, our solution shows how these PV mode and inertial-gravity waves manifest in a localized, non-uniform turbulent background - an aspect not addressed in previous studies. Later in Section 2.3, we will explore the associated flow structures and elucidate their role in energy dispersion.

2.3. Numerical validation of the analytical results

To assess the applicability of the analytical results in Section 2.2 and to develop a better understanding of the evolution of a single eddy in a rotating stratified fluid, a simple initial velocity field of the Gaussian-eddy form is chosen:

$$\mathbf{u}=\Lambda rexp\left(-\frac{r^2+z^2}{{\delta }^2}\right){\mathbf{e}}_{\theta },$$

(15)

for which Davidson et al. [14] have given an analytical solution when only the rotation is considered. In this equation, $\Lambda$ is the characteristic angular rotation rate and $\delta$ is the characteristic size. The corresponding Rossby and Froude numbers are defined by $Ro=\Lambda /\left(2\Omega \right)$ and $Fr=\Lambda /N$, respectively. To validate our analytical solutions, we chose parameters $Ro=0.02$ and $Fr=0.16$, corresponding to a vortex characterized by $\Lambda =0.422$ and $\delta =0.125$. The analytical solution is then obtained by substituting (15) into (14). We then performed a DNS to study the evolution of this Gaussian eddy. Equations (1)-(3) are solved using a parallelized pseudo-spectral code in a ${512}^3$ periodic box of size $l_{box}=2\pi$, where the fourth-order Runge-Kutta time-stepping scheme is employed. The linear terms caused by rotation and stratification, together with the viscous and diffusive terms, are integrated exactly using an integrating factor technique. A combination of phase-shifting and truncation is used to de-alias the nonlinear terms [see e.g. [36].

Figure 1(a) and (b) show contours of $u_{\theta }^2$ (normalized by its maximum value) in the $x-z$ plane passing through the axis of the Gaussian eddy at $t/t_f=13.5$. These are results from the analytical solution and the DNS respectively. Here $t_f=1/\left(2\Omega \right)$ denotes the rotation time scale. The DNS result, which takes into account the full Navier-Stokes equations, is still in close agreement with the analytical counterpart. A prominent region of $u_{\theta }^2$ near the origin stands out, which is absent in the unstratified scenario [Figure 3 [30]. This region is indicative of the PV mode, while the other intense regions represent inertial-gravity waves.

To illustrate the evolution of the flow structures and their corresponding energy dispersion, Figure 2 shows iso-surfaces of the vertical vorticity, ${\omega }_z=0.1\Lambda$, color-coded with $\varphi$, at times $t/t_f=13.5$ and $16.9$. As the system begins to evolve, a central vortex appears at the origin, flanked by two vortices, one above and one below (Figure 2(a)). Over time, the central vortex undergoes minor changes, while the flanking vortices become smaller and their centers move further away from the origin (Figure 2(b)). This behavior is in stark contrast to the purely rotating case [30], where only two columnar structures emerge, progressively moving away from each other as they elongate. The shifts in the flow structures can be attributed to the interplay between the Coriolis and buoyancy forces. In Figure 2(a), all three cyclonic vortices ($w_z>0$) are influenced by the Coriolis force, causing them to stretch horizontally. For the central vortex, $\varphi$ has negative values at the top and positive ones at the base. This distribution results in a vertical stretch from the buoyancy force, which counteracts the Coriolis force due to the incompressibility of the fluid. As a result, this vortex remains relatively stable, maintaining a nearly consistent energy distribution over time. For each flanking vortex, $\varphi$ retains its sign, showing larger absolute values for smaller magnitudes of $\left|z\right|$. This implies that the buoyancy force is acting to move these vortices away from the origin and to compress them vertically. The interaction of the Coriolis and buoyancy forces causes these vortices to recede and decrease in size, resulting in energy dispersion and potential-kinetic energy exchange. It is worth noting that a similar phenomenon to the central vortex has been documented in homogeneous rotating stratified turbulence, where it correlates with a pronounced exchange between kinetic and potential energy [37].

3. A stratified turbulent cloud under rotation

In this section, we move from the study of a single eddy to the nonlinear evolution of a turbulent cloud. Using DNS, we study the flow structures and observe characteristics consistent with inertial-gravity waves. Finally, we calculate the energy radiated from the turbulent cloud by these waves to understand their role in energy transfer for localized turbulent patches.

3.1. DNS of a turbulent cloud

Following the procedure described by Davidson [38], we performed DNS on fully developed homogeneous isotropic turbulence to obtain the initial velocity field for the turbulent cloud study. We then selected a velocity field corresponding to a specific moment characterized by an integral length scale of $l_0\approx 0.04l_{box}$ and a Reynolds number $R{e}_0=u_0{l}_0/\nu \approx 169$, where $u_0$ is the r.m.s. velocity. To numerically generate a horizontal turbulent cloud, the resulting velocity field was spatially filtered with

$$f\left(z\right)=exp\left(-{10}^a{\left|z-\pi \right|}^m\right),$$

(16)

as shown in Figure 3(a). Here we set $a=4$ and $m=9$, resulting in a cloud thickness of $l_c=0.13l_{box}$, approximately three times the integral length scale $l_0$. The buoyancy field was initialized with $\varphi =0$ throughout the computational domain, which ensures a uniformly linear density profile, indicating no initial mixing within the cloud. Figure 3(b) shows the initial turbulent cloud as indicated by the velocity module, and its spectrum in the x-y plane, averaged over its central part in the z-direction, is shown in Figure 3(c).

To investigate the effect of different stratification on the evolution of a turbulent cloud under rotation, we conducted DNS with a constant Rossby number, $Ro=0.11$, and five different Froude numbers: $Fr=\infty$ (corresponding to $N=0$, a purely rotating scenario), $0.88$, $0.44$, $0.22$ and $0.11$. The simulation parameters are given in Table 1, where $Ro=u_0/\left(2\Omega l_0\right)$ and $Fr=u_0/\left(N{l}_0\right)$. With the exception of the purely rotating case R0.11, the simulations are denoted by their Froude numbers, e.g. the case F0.88 corresponds to $Fr=0.88$.

3.2. Evolution of flow structures

To show the evolution of the turbulent cloud structures, Figure 4 shows iso-surfaces of $u_z$ for all simulations at different times. For the case R0.11 under pure rotation, Figure 4(a) shows the emergence of vertical columnar structures from the turbulent cloud, which extend into the quiescent region, in agreement with previous experimental results [14] and numerical study [15]. For the weakly stratified case F0.88, Figure 4(b) shows subtle differences in the flow structures at $t/t_f=8.4$ and more pronounced variations at $t/t_f=12.7$ compared to the case R0.11.

For the moderately stratified cases F0.44 and F0.22, Figure 4(c) and (d) show that the columnar structures tend to deviate from a purely vertical orientation as they extend into the quiescent region. This tilt becomes more pronounced with increasing stratification. A comparison of Figure 4(a)-(e) shows that, at equivalent time intervals, the vertical growth of these structures decreases with decreasing $Fr$. This suggests that strong density stratification limits the vertical expansion of the cloud, consistent with the observation in experiments [27].

In the case F0.11 (Figure 4(e)), the rotation frequency and the Brunt-Väisälä frequency are identical. As a result, no structures emerge from the turbulent cloud, suggesting minimal energy emission from the cloud. This behavior mirrors the case without both rotation and stratification, as documented in [15]. A plausible interpretation is that the stratification strength is sufficient to completely inhibit the vertical expansion of the turbulent cloud in the presence of system rotation. It is noteworthy that Equation (7) for this scenario yields $\varpi =2\Omega =N$, implying uniform frequencies of inertial-gravity waves in all directions. Thus, in the absence of both rotation and stratification, the system exhibits isotropy.

3.3. Are the flow structures formed by inertial-gravity waves?

Columnar structures from a turbulent cloud under rotation have been identified as inertial waves [15], while the large-scale pancake structures from a stratified turbulent cloud are recognized as internal gravity waves [22]. This section examines whether the flow structures discussed in Section 3.2 from the stratified turbulent cloud under rotation can be attributed to inertial-gravity waves.

In all simulations, the growth of flow structures in different directions leads to the expansion of turbulent clouds in the vertical direction, which can be quantified by tracking the boundary of the $|u_z|$ iso-surfaces with time. The reason we choose $u_z$ to indicate structure is that it satisfies the wave equation as described by (6). In particular, given a threshold $|u_z^b|$, the upper and lower boundaries of the cloud can be obtained from the positions where $|u_z|=|u_z^b|$. We can therefore calculate the vertical extent of the cloud, $h(x,y,t)$, for each point on the $x-y$ plane and time t. The mean wave-cloud thickness $\overline{h}\left(t\right)$ is defined as the average of $h(x,y,t)$ over x and y. Figure 5 shows the mean wave-cloud thickness for all the simulations where we choose $|u_z^b|=0.35u_0$. The results remain robust when the threshold varies around this value, and this threshold is used consistently throughout the analysis in this paper. For the purely rotating case R0.11, the mean wave-cloud thickness shows a linear growth with time, resulting from the inertial waves with $k_z=0$ traveling at the group velocity $c_g=2\Omega ⁄k_h$[15]. For the case F0.88 with weak stratification, the mean wave-cloud thickness is close to that in the case R0.11 at $t\le 0.2$, indicating that the stratification does not affect the cloud before $t=0.2$. This is because that the linear time scale corresponding to the rotation, $t_f=1/\left(2\Omega \right)$, is smaller than that corresponding to the stratification, $t_N=1/N$. When $t>0.2$, the mean wave-cloud thickness shows a smaller linear growth rate. For the cases F0.44 and F0.22 with stronger stratification, the behavior of the mean wave-cloud thicknesses can be interpreted as the superposition of a linear growth and an oscillation with a fixed frequency over time. In addition, the frequency is larger in the case F0.22. The mean wave-cloud thickness in the case F0.11 is almost constant with time, consistent with the observation in Figure 4(e).

To investigate whether the flow structures are closely related to the inertial-gravity waves, we compare the results from the DNS with those predicted by the linear theory. Since the initial turbulent cloud can be viewed as a distribution of random vortex blobs of different shapes and sizes, we consider the linear evolution of a Gaussian vortex satisfying Equation (15) under the same Rossby and Froude numbers. We set $\Lambda =2.31$, $\delta =0.125$, and utilize $|u_z^b|=0.002\Lambda \delta$ to determine the vertical extent of the eddy, $h_e$. This extent is defined as the largest vertical distance between the positions where $|u_z|=|u_z^b|$. Figure 5(b) shows that the results for $h_e$ align qualitatively with those of $\overline{h}$ in Figure 5(a), suggesting that the expansion of the turbulent cloud may be prominently influenced by the linear mechanism.

We now show quantitatively that the flow structures are closely related to inertial-gravity waves. For a vortex blob in the initial turbulent cloud, suppose it is characterized by a length scale of $1/k$. Then its rescaled group velocity in the vertical direction is

$${\tilde{c}}_{g,z}=k{c}_{g,z}=\pm \frac{2\Omega (1-{(N/\left(2\Omega \right))}^2)(k_z/k_h)}{\sqrt{{(N/\left(2\Omega \right))}^2+{(k_z/k_h)}^2}(1+{(k_z/k_h)}^2)}.$$

(17)

Given the ratio of the Brunt-Väisälä frequency to the rotation frequency $N/\left(2\Omega \right)$, it is determined that when

$${\left(\frac{k_z}{k_h}\right)}^2=\frac{1}{4}\left[-{(N/\left(2\Omega \right))}^2+\sqrt{{(N/\left(2\Omega \right))}^4+8{(N/\left(2\Omega \right))}^2}\right],$$

(18)

${\tilde{c}}_{g,z}$ reaches the maximum ${\tilde{c}}_{g,z}^{max}$, which is the vertical growth rate of the turbulent cloud predicted by the linear inertial-gravity waves. Furthermore, the corresponding value of $k_z/k_h={(k_z/k_h)}_0$ characterizes the tilt angle of the flow structure. The vertical growth rate of the turbulent cloud in the DNS can be determined from the growth rate of the mean wave-cloud thickness, $\overline{h}\left(t\right)$, as shown in Figure 5(a). The rate is calculated using a linear fitting method. In Figure 6(a) we present two vertical growth rates: ${\alpha }_{DNS}$ derived from DNS and ${\alpha }_{LT}$ derived from linear theory. To facilitate comparison, both rates are normalized by their values at $N/2\Omega =0$. Figure 6(b) shows the angle between the flow structure and the vertical direction, $\theta =arctan\left[{(k_z/k_h)}_0\right]$, predicted by linear theory for different $N/\left(2\Omega \right)$. The results are consistent with Figure 4, where the angle becomes larger as $N/\left(2\Omega \right)$ increases. For quantitative comparison, we also plot cones in the third columns of Figure 4(b), (c) and (d), where the angle between the cone surface and the vertical axis is equal to $\theta$. It is notable that there is a strong agreement between the DNS results and the predictions of linear theory. However, for the case F0.22, both $\alpha$ and $\theta$ are slightly underestimated by the linear theory. This discrepancy could be due to the interactions between waves and turbulence, especially given the small vertical expansion of the cloud in this case. It is important to note that in scenarios where both rotation and stratification are present, the frequency of the inertial-gravity waves (7) tied to ${\tilde{c}}_{g,z}^{max}$ is not zero. This explains the oscillatory component observed in $\overline{h}$ for the cases F0.44 and F0.22. When only rotation or stratification is present, the maximum group velocity always coincides with inertial (or gravity) waves of zero frequency. Under these circumstances, tools such as two-dimensional energy spectrum analysis can be used effectively, given the consistent phase of $u_z$[15] or $u_x$ [22] at the edge of the turbulent cloud.

3.4. Wave-dominated and turbulence-dominated regions

In Section 3.3 it is highlighted that the flow structures emerging from the turbulent cloud show a strong association with inertial-gravity waves. In stratified turbulence, PV has been used to identify turbulent/non-turbulent interfaces [22,39]. Our study, however, uses Lagrangian particle tracking to measure the vertical extent of turbulent advection [15]. Lagrangian particles act as passive tracers, with their position $\mathbf{r}$ at time t given by $\mathrm{d}\mathbf{r}/\mathrm{d}t=\mathbf{u}(\mathbf{r},t)$. In this relationship, $\mathbf{u}$ represents the velocity of the fluid. To determine $\mathbf{u}$ at the particle positions, we use a 6th-order Lagrange interpolation. Furthermore, the positions of these particles are integrated using a 2nd-order Adams-Bashforth scheme.

Tracer particles with a number $N_p=65536$ are initialized randomly from a uniform distribution throughout the initial turbulent cloud, the vertical range of which is denoted as $[z_0,z_1]$. Let $N_p^{out}$ be the number of particles outside of the range of the initial turbulent cloud. The r.m.s. distance of these particles to the initial turbulent cloud can be calculated as

$$d_{rms}=\sqrt{\frac{1}{N_p^{out}}\sum _{p=1}^{N_p^{out}}{d}_p^2},$$

(19)

where

$$d_p=\left\{\begin{array}{cc}{z}_p-z_1,\hfill & z_p>z_1\hfill \\ z_0-z_p,\hfill & z_p<z_0\hfill \end{array}\right.$$

(20)

represents the distance of the particle p to the initial turbulent cloud. We consider $[z_0-d_{rms},z_1+d_{rms}]$ as the vertical range of the turbulence-dominated region and the other domain as the wave-dominated region. Figure 7 shows the ratio of particles outside the initial turbulent cloud, $N_p^{out}/N_p$, and $d_{rms}$ normalized with the initial cloud thickness $l_c$. For the case R0.11, both $N_p^{out}/N_p$ and $d_{rms}/l_c$ increase with time, while for the cases F0.88 and F0.44 both quantities first increase and then decrease after reaching the maximum values, probably due to the stable stratification. However, in the cases F0.22 and F0.11, where the stratification is strong, both quantities immediately reach stationary states around small values.

Table 2 presents the maximum values of $N_p^{out}/N_p$ and $d_{rms}/l_c$ for all cases, where $max(d_{rms}/l_c)$ can also be interpreted as the maximum vertical expansion of the initial turbulent cloud. It is observed that both quantities become smaller with increasing stratification, indicating the limitation of vertical turbulent advection by stratification.

Once the regions dominated by waves and turbulence are identified, we can determine the energy of the turbulent cloud carried away by the inertial-gravity waves, denoted as $E^{wave}$. This energy is equal to the energy within the wave-dominated region. Figure 8(a) shows that for cases where $Ro<Fr$, the ratio $E^{wave}/E$ exhibits a monotonic increase with time. In particular, for a given time point, a decrease in $Fr$ results in smaller values of $E^{wave}/E$. This suggests that under stronger stratification, inertial-gravity waves contain a smaller fraction of the energy. In contrast, for the scenario where $Ro=Fr$, the ratio $E^{wave}/E$ oscillates with a fixed pattern. Its values fluctuate between nearly zero and a peak of about $15\%$, indicating that the inertial-gravity waves only sporadically extract a small fraction of energy from the turbulent cloud. Thus, while inertial waves account for a significant fraction of the energy in a system dominated by rotation alone, this energy contribution diminishes with the introduction of stratification. Figure 8(b) shows the ratio of potential energy to total energy in the wave-dominated region, $E_P^{wave}/E^{wave}$. As stratification increases, potential energy becomes increasingly dominant and varies more frequently with time. This highlights the increasing importance of the kinetic-potential energy exchange and suggests that the dominant inertial-gravity waves exhibit greater phase velocities under increased stratification. It is worth noting that a similar kinetic-potential energy exchange phenomenon has been reported in homogeneous rotating stratified turbulence [40].

4. Discussion and Conclusions

In this study, we investigate the evolution of a stratified turbulent cloud influenced by rotation. By analyzing the behavior of a single eddy under both rotation and stratification, we observe clear deviations from the purely rotating scenario. In particular, within the stratified environment, we detect a spatial manifestation of a known PV mode, characterized as a stable vortex at the origin. Furthermore, the system transitions from exhibiting inertial waves to producing inertial-gravitational waves, whose propagation relations depend on both rotation and stratification.

Numerical simulations have been performed to investigate the nonlinear evolution of the turbulent cloud influenced by different degrees of stratification under rotation. In cases where $Ro<Fr$, columnar structures emerge spontaneously from the turbulent cloud. These structures are vertically oriented in the case of pure rotation, but deviate from the vertical as stratification is introduced. As $Fr$ decreases, the tilt angle increases and the rate of cloud expansion decreases. Conversely, in the $Ro=Fr$ system, the turbulent cloud remains devoid of distinct structures, reminiscent of scenarios without rotation and stratification. The flow structures are observed to have characteristics consistent with inertial-gravity waves. This observation is supported by the agreement between DNS and linear theory, which give consistent predictions for the vertical growth rates of the turbulent cloud and the tilt angles of these structures. Using Lagrangian particle tracking, we distinguish between wave-dominated and turbulence-dominated regions. It is found that the inertial gravity waves emitted by the cloud can transport a substantial fraction of the energy to the neighbouring quiescent fluid, which decreases with increasing stratification. This suggests the central role of inertial-gravity waves in energy transfer processes within atmospheres and oceans, particularly around regions of localized turbulence.

The present study elucidates the effects of stratification on localized turbulence under system rotation. First, compared to the purely rotating case, the presence of stratification introduces the zero-frequency PV mode, whose energy distribution barely changes with time. While this mode is not new, our finding is its distinct spatial manifestation within a non-uniformly turbulent background. Second, as the stratification intensifies, the dominant inertial-gravity waves exhibit increased tilting and reduced vertical group velocity. In combination, these phenomena limit the energy emitted by the initial turbulent cloud. It should be noted that this study focuses on the parameter regime where $Ro<Fr$. While the $Ro>Fr$ regime is thought to dominate in the Earth’s atmosphere and oceans [6], the $Ro<Fr$ regime is relevant for several geophysical flows, such as in the deep sea where the influence of stratification is weak but rotation effects are prominent [41,42,43,44,45]. Furthermore, the observed flow structures can be conceptualized as Fourier modes with phase coherence. This provides an alternative perspective: the evolution of these turbulent clouds could be seen as the evolving phase coherence of Fourier modes influenced by both rotation and stratification, which is an avenue for future research.