{'fraction_non_alphanumeric': 0.03506642389911659, 'fraction_numerical': 0.011632611774226179, 'mean_word_length': 4.643958135109419, 'pattern_counts': {'":': 0, '<': 5, '<?xml version=': 0, '>': 1, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 3, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Absorption of waves by large scale winds in stratified turbulence

The atmosphere is a nonlinear stratified fluid in which internal gravity waves are present. These waves interact with the flow, resulting in wave turbulence that displays important differences with the turbulence observed in isotropic and homogeneous flows. We study numerically the role of these waves and their interaction with the large scale flow, consisting of vertically sheared horizontal winds. We calculate their space and time resolved energy spectrum (a four-dimensional spectrum), and show that most of the energy is concentrated along a dispersion relation that is Doppler shifted by the horizontal winds. We also observe that when uniform winds are let to develop in each horizontal layer of the flow, waves whose phase velocity is equal to the horizontal wind speed have negligible energy. This indicates a nonlocal transfer of their energy to the mean flow. Both phenomena, the Doppler shift and the absorption of waves traveling with the wind speed, are not accounted for in current theories of stratified wave turbulence.


I. INTRODUCTION
Turbulent flows are often pictured as highly disorganized, with energy being transferred from large scale motions to small scale eddies. But in some cases the opposite can happen. In 1967 Kraichnan [1] developed the theory of the inverse energy cascade, in which nonlinearities transfer energy towards larger structures in a selforganization process. The theory had a huge impact in oceanography and meteorology, and in the basic physical understanding of turbulence. However, an inverse cascade is not the only possible mechanism by which energy can be transferred from small to large scales, and in some flows other mechanisms can result in the generation of large scale flows.
An important example is given by stratified flows. Stratification plays a key role in the dynamics of the oceans and the atmosphere. As an example, internal gravity waves arising from it are responsible for transfer of momentum between different regions of the atmosphere [2], and for transport in the oceans [3,4]. Much effort has been put in characterizing internal gravity waves and their relation with geophysical flows [5][6][7].
Several observational studies suggest there is a coupling between the mean (or large-scale) flow and internal gravity waves [8][9][10]. As a result, various models for the interaction between a wave field and a mean background flow have been put forward [11][12][13]. These models, which were mainly developed in the context of atmospheric sciences and oceanography, predict that the dispersion relation of the waves is Doppler shifted by the background flow. Furthermore, they also describe an instability by which waves whose horizontal phase velocity matches that of the horizontal wind speed on a given layer of the fluid, are absorbed by the background flow. This is the so-called critical layer (CL) instability. A fraction of the energy in the waves is then transferred towards the mean flow, and the rest is dissipated. However, most CL models are either linear or perturbative, have the large scale shear externally imposed, lack nonlinear coupling between the waves themselves as well as with small-scale eddies, and cannot fully describe the nonlinear dynamics of a turbulent geophysical flow. Important observational studies of the atmosphere [14] and the ocean [15], and laboratory experiments [16,17] focusing on the properties of single waves passing through, show Doppler shift and indications of CL instability, although they lack the presence of a fully developed turbulent superposition of waves. In the same spirit, numerical studies of Doppler shift and CL instability [12,13,18,19] that solve nonlinear atmospheric models, focus on single wave packets interacting with a background flow which is externally imposed, and do not take into account the turbulent nature of geophysical flows nor the nonlinear interaction between the internal gravity waves themselves and with the eddies present in the flow (see however [18], where the full nonlinear equations are solved, although a single wavepacket and background flow are imposed). This can be understood as identifying waves in a smooth flow can be done by observing the time evolution of the system, while extraction of the waves in a disordered background requires space and time information.
On the other hand, fluid dynamics studies, both theoretically and numerically (see, e.g., [20][21][22][23]), focused mainly on the fully nonlinear aspect of stratified turbulence, on how energy is distributed among scales, and on the development of anisotropy and of flat "pancake" structures in the flow. These studies remark that the interaction between waves and eddies is of fundamental importance in stratified turbulence, but have neglected the effects of Doppler shifts or CL instability. An important topic which is subject to an active debate is the existence of an inverse cascade in purely stratified flows. As already mentioned, an inverse cascade is a self-similar process by which energy is transferred nonlinearly from small scales to larger scales, so energy gets accumulated in structures with larger correlation length than that of the injection mechanism. Inverse cascades have been ex-tensively studied in two dimensional turbulence [1], and have also been observed in rotating turbulence (including rotating stratified turbulence) [24][25][26]. But results for purely stratified turbulence are unclear. It is known that stratified flows can generate large-scale vertically sheared horizontal winds (VSHW) [27], although the mechanism involved is not entirely understood. While evidence of energy flow toward large scales has been reported [28,29], several authors have argued that an inverse cascade is not possible [21,30] based on statistical mechanical arguments.
Recently, weak turbulence theory [31,32] has also been used to study stratified turbulence, although only considering the coupling between the waves and not with the eddies. While this theory can give information on the formation of the nonlinear energy cascade, it cannot characterize interactions of the wave field with a mean flow. Some extensions that consider nonlocal interactions and the development of zonal flows have been proposed to address this, and are specially relevant in the context of quasi-geostrophic turbulence [33], and plasma turbulence in tokamaks [34,35] to explain the development of large-scale flows.
In this article, we identify a mechanism working within a stably stratified turbulent flow that couples the wave field with the large-scale VSHW, thus gaining a better understanding of the role of waves and the generation of large structures in turbulent flows. We obtain for the first time the four-dimensional space and time resolved energy spectra, allowing us to uniquely identify the wave components of the total flow and to study their dynamics. We show that waves are Doppler shifted by the VSHW, and that while a direct nonlinear energy cascade is present, energy can also be transferred nonlocally from the smallscale waves to the large-scale flow though CL instabilities. Finally, we also show that the energy spread created by the Doppler shift is not uniform, a fact that results from the turbulent nature of the dynamically evolving large-scale flow.

II. STABLY STRATIFIED FLOWS
In this section we present the equations we solve numerically that describe a stably stratified flow, introduce characteristic time scales, and give a brief introduction to the linear theory of gravity waves in the presence of background shear.

A. The Boussinesq equations
We consider the Boussinesq equations describing an incompressible stratified flow, where u is the velocity field (with ∇ · u = 0), θ is the potential temperature fluctuations, p is the pressure normalized by the mean fluid density, ν is the kinematic viscosity, κ is the thermal diffusivity, N is the Brunt-Väisälä frequency (associated with the vertical background stratification), and F is an external mechanical forcing.

B. Characteristic timescales
The important characteristic timescales that come into play in this system are the wave period τ ω , and the nonlinear turnover time τ N L . The wave period is simply given by τ ω ∝ 1/ω 0 . The nonlinear turnover time is the time taken by an eddie at a certain scale to transfer its energy through nonlinear interactions to smaller scale eddies locally in wavenumber space. From dimensional analysis, this timescale can be estimated as where ∝ 1/k is a lengthscale in the inertial range, u is the characteristic velocity of eddies at scale , and E(k) is the energy spectrum. If the energy spectrum in stratified turbulence follows a power law E(k) ∝ k −3 (see, e.g., [20]), then it is easy to show that the turnover time is independent of the scale, and becomes τ N L ∝ L/u rms at all scales.

C. Gravity waves in a background flow
In the presence of a background horizontal flow U(z), the frequency measured by an observer is shifted by the Doppler effect, resulting in a frequency Considering for simplicity the unidirectional case with U = U (z)x, and substituting u = U + u in Eqs. (1) and (2) with u = w(z)e ik ⊥ (x−ct)ẑ a planar internal gravity wave propagating inx with phase velocity c = ω/k, we obtain after linearization an equation for the vertical velocity amplitude [37] This equation has a singularity when c = U (z), which is associated with the existence of a CL in which the Reynolds stress tensor has a discontinuity [37]. As a result, a wave traveling through the fluid an approaching the layer with c = U (z) develops faster and faster fluctuations in its vertical velocity amplitude, becomes unstable, and finally degenerates into turbulence, transferring a fraction of its energy to the background flow, while the rest is dissipated. These latter processes cannot be properly described by the linear theory. In practice, the instability only occurs if the Richardson number, Ri = N 2 /(∂ z U ) 2 , is greater than 1/4 [37]. Going beyond the linear regime, the detailed numerical study in [18] considers the energy balance when a wave packet crosses a CL. The authors were able to quantify that more than a third of the energy in the waves is nonlocally transferred to the mean background flow, while the rest of the energy in the waves is transferred towards the small scale turbulence where it is finally dissipated.

III. NUMERICAL SIMULATIONS
Direct numerical simulations of Eqs. (1) and (2) were performed in a cubic periodic box using the parallel code GHOST [38][39][40] that uses a pseudospectral method to estimate spatial derivatives, and evolves the system in time using a second-order Runge-Kutta method. Two forcing functions were used: Taylor-Green forcing with F = f 0 (sin x cos y cos zx − cos x sin y cosẑ), and a randomly generated isotropic three-dimensional forcing acting at k = 1. As we will consider frequency spectra, in both cases the forcing was kept constant in time to prevent exciting spurious timescales in the system, and to prevent disrupting the development of VSHW.
A Taylor-Green flow excites two-dimensional motions and has been used in previous studies of stably stratified turbulence to mimic atmospheric motions [41], while isotropic forcing injects energy directly into vertical motions and thus excites stronger gravity waves. Some of the effects we show are disrupted when VSHW cannot fully develop, so in the next section we focus first on the isotropically forced case with strong stratification. Then, in Sec. IV C we present results comparing both forcings and explain in more detail the impact of Taylor-Green forcing in the development of the CL instability. In that section, we also analyze the effects of varying the level of stratification.
The simulations can be characterized with dimensionless numbers. All simulations have Reynolds number Re = Lu rms /ν ≈ 9700 (where L is the forcing scale and u rms is the r.m.s. velocity), and Prandtl number Pr = ν/κ = 1. For each forcing function, we performed three simulations with varying stratification N , resulting in Froude numbers Fr = u rms /(N L) ≈ 0.04, 0.02, and 0.01 (with decreasing Fr indicating stronger stratification), and buoyancy Reynolds number Re b = ReFr 2 ≈ 15, 3.9 and 0.97. Each simulation was started from the fluid at rest, so neither the large scale flow nor the small scale waves were previously imposed, and the simulations were continued for 12 large-scale turnover times to let the system reach a turbulent steady state. Then, each simulation was continued for other 12 turnover times saving all fields with high temporal cadence, resulting in 1200 outputs of each field that sample the Brunt-Väisälä frequency 40 times per period when Fr ≈ 0.01. This, previously unmatched temporal resolution, allows us to resolve four-dimensional spectra simultaneously in space and in time, although at a computational cost that allows us to reach only moderate spatial resolutions with 512 3 spatial grid points. Nonetheless, as we show, the spatial resolution we use is appropriate, and recent studies also indicate that properties of fully developed turbulence can be identified even at moderate resolutions [42].

A. Spatial spectral analysis
Traditional characterization of turbulent flows is done using one-dimensional spectra. In Fig. 1 we present the parallel spectrum E(k ) for the three simulations with isotropic forcing. Only as a reference, we also show in Fig. 1 a power law ∝ k −3 . The reader can find detailed spatial characterizations of stratified turbulence, which go beyond the aim of this work and at higher spatial resolution, in recent studies (see, e.g., [21,22]).
As a better way to characterize energy distribution among scales in the presence of anisotropy, Fig. 2 shows the spatial axisymmetric energy spectrum e(k ⊥ , k ), normalized by sin θ with θ = arctan(k ⊥ /k ) to obtain circular isocountours in the case of an isotropic flow. In Fig. 2, two spectra are presented, for the simulations with isotropic forcing and Fr ≈ 0.02 and ≈ 0.01. As stratification is increased, energy distribution becomes more anisotropic, with energy being preferentially trans- Fr ≈ 0.01. Note the spectral anisotropy resulting from stratification. The solid line indicates the modes with wave period equal to the eddy turnover time, τω = τNL. In (b) a ridge is formed close to this curve, indicating energy is transferred towards modes with lower k ⊥ but the transfer is halted when τω ≈ τNL. Modes with τω < τNL (i.e., modes below the solid line) are often called wave modes, as these modes have the wave period as the fastest timescale. However, a proper characterization of waves requires space and time resolved spectra.
ferred towards modes with smaller k ⊥ (and, as a result, larger wave period) [20][21][22]27]. However, for Fr ≈ 0.01 energy accumulates near the modes with wave period (τ ω ∝ 1/ω 0 ) equal to the nonlinear turnover time (τ NL ∝ L/u rms ), forming a ridge. As the energy transfer mechanism is often given by the shortest timescale [43], modes below the curve τ ω = τ NL (those with wave period shorter than the turnover time) are associated with waves [21,22,27]. Modes above the curve τ ω = τ NL (and in particular, modes with k ⊥ = 0) are often called vortical modes, as for these modes ω 0 ≈ 0. The fraction of the energy in these modes decreases with decreasing Fr, in good agreement with the observed accumulation near τ ω = τ NL for large Fr. The slow down of the transfer as the energy reaches the ridge is compatible with critical balance arguments [31], and is also of great importance in weak wave turbulence theories [32] that require the energy of the system to remain in weakly interacting waves. Note however that this distinction between waves and vortical modes in Fig. 2, based on the characteristic time of each mode, is only approximate. A strict discrimination between waves and eddies requires the four-dimensional spectrum.

B. Spatio-temporal analysis
Precise identification of the waves, and of their role in the dynamics, requires both space and time information. Figure 3 shows different cuts of the frequency and wavenumber spectrum E θ (k x , k y , k z , ω) for k x = 0, and for either k z = 0 or k z = 10 in the simulation with the strongest stratification and isotropic forcing. As internal gravity waves couple vertical motions with temperature fluctuations, we consider the spectrum of potential energy E θ to isolate the waves more easily. In Fig. 3 there is no significant accumulation of energy in the modes that satisfy the dispersion relation given by Eq. (3), i.e., in modes that could be directly associated with internal gravity waves. Instead, most of the energy (≈ 80% for k y ∈ [10, 80]) falls inside a wider region defined by the Doppler shifted Eq. (5) with −0.4 ≤ U y ≤ 0.4. This suggests the waves are propagating through layers with non-zero mean horizontal velocity (i.e., VSHW). Indeed, a probability density function (PDF) of the Cartesian components of the horizontally averaged horizontal velocity U obtained directly from the same simulation shows two peaks at ≈ ±0.4 (see the inset in Fig. 4). Also, the PDF of the point-wise horizontal velocity in the entire box is approximately Gaussian with variance close to 0.4.
It is important to remind that we force the fluid from rest at the large scales; this generates the winds, and by nonlinear interactions energy cascades towards the smaller scales. Then, the cycle completes itself when the small scale waves interact with the large scale winds; we are not imposing flows at different scales, or imposing a background flow, as sometimes done to study interaction of waves with shear flows. Instead, we just let turbulent interactions develop by themselves. In similar simulations of turbulent stratified flows, waves have been directly observed before by measuring the frequency spectrum E(ω), or by measuring E(ω) for a few Fourier modes [23]. This study identified peaks at the frequency of the internal gravity waves, but from the frequency spectrum a relation between frequency ω and wavenumber k cannot be obtained without assuming the system is dominated by the waves (and then the relation of these quantities is given by the dispersion relation), or domi- Fraction of the energy F that is contained within the two Doppler shifted branches as a function of the wavenumber. In the inertial range, ≈ 80% of the energy corresponds to Doppler shifted waves. nated by vortical motions (and then the relation is given by sweeping, see, e.g., [43]). Figure 3 gives direct evidence that most of the energy is in waves. And the waves satisfy Eq. horizontal wind. However, there is more power in modes close to the curve ω 0 + U y k y (with U y = 0.4) than in modes close to ω 0 − U y k y . This energy distribution cannot be explained by a preferential direction in U (see the PDFs in the inset in Fig. 4), by viscous effects (which introduce damping but no frequency shift), or by nonlinear corrections to the dispersion relation, as internal gravity waves are exact nonlinear solutions of the Boussinesq equations [36]. The preferential concentration is instead compatible with a CL instability. In Fig. 3 we show as a reference the area (shaded in transparent gray or orange) where ω < U y k y with U y = 0.4. As in this simulation −0.4 U y 0.4 (once again, see the PDFs in the inset in Fig. 4), for all modes in that area there is some layer with U y such that the layer can be a CL. When the waves approach these layers, c ≈ U , and from Eq. (6) their vertical velocity starts fluctuating more rapidly, undergoing an instability and transferring a fraction of their energy to the flow [13,18] (i.e., the waves are absorbed). Indeed, there is almost no energy in the wave modes in this area. This mechanism can only act if Ri > 1/4. We verified the point-wise value of Ri in the simulations is larger than 1/4. For the run with Fr ≈ 0.04, the minimum of Ri is ≈ 14. Although this value is rather large for atmospheric flows, it has been reported in other simulations of stratified turbulence, and indicate that CL instability can be the reason for the nonlocal formation of large scale structures in stratified turbulence observed in [27][28][29].
Internal gravity waves in a stratified fluid couple the temperature with the vertical component of the velocity. We now show that the four-dimensional power spectrum of the vertical velocity displays the same features as the spectrum of the potential energy shown above. Figure  5 shows a cut of the power spectrum of the vertical velocity E z (k x = 0, k y , k z = 0, ω) in the simulation with stronger stratification and isotropic forcing. The same features found in the four-dimensional spectrum of the potential energy can be found in this figure, including the Doppler shift and the defect of energy in the modes that have frequency compatible with CL instability. In fact, the spectra are practically indistinguishable. As in the case of the temperature, a large fraction of the energy in vertical motions is also associated with wave motions, (≈ 80% of the energy between 10 < k y < 80 is inside the fan corresponding to Doppler shifted waves, see the inset in Fig. 5). When the spectrum of horizontal velocity is considered instead, none of these signatures can be found (not shown).
As a final and independent verification of the excess of waves traveling with the flow observed in Fig. 3, we resort to a statistical analysis. As most of the energy in the simulation with Fr ≈ 0.01 is in the waves, we can assume the vertical velocity is a superposition of traveling waves ∝ e i(k·x−ωt) . We can then compute U · ∇u z /(|U||∇u z |), which gives an estimation of how the wave vector (i.e., the propagation direction for a pure planar wave) is aligned with the horizontal flow. Figure 4 shows the PDF of this alignment (averaged over every layer in the fluid) minus the PDF of the sine of a uniformly distributed angle θ ∈ [0, 2π) (i.e., of randomly aligned fields with no privileged direction), for the simulation with stronger stratification and isotropic forcing. The result indicates a preference towards an alignment of horizontal gradients of the vertical velocity with the mean flow, as there is an excess for +1 when compared with −1, and a deficit (compared with uniformly distributed angles) for the case in which the two fields are perpendicular. Interestingly, a preference of stratified flows towards developing one sign of velocity gradients (resulting in non-Gaussian PDFs) was reported before in [44].
An interesting fact, that sheds more light on how these mechanisms work, is that these effects are not observed in purely rotating flows [26,43]. While these flows can also generate strong horizontal velocity fields, waves tend to travel in the vertical direction (as compared to the stratified case, where they tend to travel horizontally), and no Doppler shift develops as the waves are perpendicular to the large-scale flow.

C. Dependence on level of stratification and forcing
The results discussed so far correspond to the simulation with randomly generated isotropic forcing and Fr ≈ 0.01. A similar distribution of energy in frequency and wavenumber is observed in the simulations with same  Taylor-Green forcing, which prevents the formation of a mean wind in each horizontal layer, the CL instability is not observed, although Doppler spreading of the waves still takes place. These results are also presented in this section.
As a comparison with the case shown in Fig. 3, Fig. 6(a) shows the space and time resolved spectrum E θ (k x = 0, k y , k z = 0, ω) for a simulation still forced with randomly generated forcing, but with Fr ≈ 0.02. This spectrum can be directly compared with Fig. 3(a). Decreasing the stratification results in a spectrum that bears great resemblance with the spectrum for Fr ≈ 0.01, with non-uniform Doppler spreading of the waves and the excess of waves with ω ≈ ω 0 + U y k y . The defect of energy in the modes compatible with CL instability is also visible in this simulation, although the damping of these modes is slightly weaker than in the case with stronger stratification. This is compatible with the smaller fraction of energy in the waves as stratification is decreased; see the inset of Fig. 6(a). Also, this is compatible with the fact that as stratification is decreased, the power in the large-scale motions (small wavenumbers) also decreases  Fig. 5. The insets show in each case the fraction of the energy F that is contained within the two Doppler shifted branches. In (b), the inset on the right also shows a schematic representation of the Taylor-Green forcing. In (a), note the decrease in the level of stratification of the system still preserves the main features of the spectrum. In (b), they Taylor-Green forcing injects most of the energy in vortical modes (ω ≈ 0), the energy in the waves is therefore smaller, and the flow geometry prevents the CL instability from developing.
(see Fig. 1). These trends are further confirmed by the simulation with weaker stratification and Fr ≈ 0.04.
When Taylor-Green forcing is used instead, the broadening of the dispersion relation by Doppler shift can still be observed; see Fig. 6(b) for Fr ≈ 0.02, specially for wavenumbers 0 < k y 10. However, this effect is dimmed by a large concentration of energy in modes with ω ≈ 0. This is to be expected as the Taylor-Green forc-ing, given by F = f 0 (sin x cos y cos zx − cos x sin y cosẑ), consists of two counter rotating vortices in the horizontal velocity, and excites no vertical motions directly (see a schematic diagram in the inset of Fig. 6(b)). As a result, the fraction of energy in wave modes is much smaller than in the case with isotropic forcing (see also a direct measurement of this in the inset of Fig. 6(b)). Moreover, the imposition of the Taylor-Green vortices in each layer by the external forcing seems to disrupt the development of a non-zero mean wind in each horizontal layer, weakening also the development of the CL instability. Indeed, the damping of the energy in modes consistent with CL instability in the four-dimensional spectrum of potential energy in Fig. 6(b) is almost inexistent.
To ease with the understanding of the behavior caused by the two different forcings, and to illustrate the presence of VSHW in the flows, in Fig. 7 we present vertical slices (at constant x) of the potential temperature and horizontal velocity fields for the two simulations. Although the simulations have the same Fr and Re, the flows are significantly different. On the one hand, there are smaller scale structures and overturning in both fields for the Taylor-Green forced simulation. On the other hand, the simulation with isotropic forcing shows a more dominant large-scale horizontal flow, and smaller scale (although smoother) fluctuations in the temperature. This is consistent with the fact that when forcing isotropically there is more energy in the waves, as seen in Fig. 6. Another way to quantify how much energy is in the VSHW is to compare the ratio of the energy in the modes with ω = 0 (i.e., in the modes with k ⊥ = 0) to the total energy. For Fr = 0.02, the Taylor-Green simulation has ≈ 0.009% of its total energy in those modes, while the isotropically forced one has ≈ 0.03%, indicating that the isotropically forced simulation has stronger VSHW. Note however that both simulations have a broad energy spectrum, indicating that turbulence in both flows may be of different nature, with the case with more energy in wave motions being smoother and possibly closer to a wave turbulence regime.
The comparison between the two different forcing functions further indicates that stratified turbulent flows may display different behavior depending on whether the mechanism used to excite the turbulence allows or prevents the development of large-scale vertically sheared horizontal winds. Interestingly, the case studied here that injects more energy directly into the waves, is also the case in which horizontal winds and CL instability more clearly develop, two features that are not considered in wave turbulence theories.

V. CONCLUSIONS
In a turbulent flow, waves and individual instability events cannot be easily identified. As a result, previous numerical and observational studies of Doppler shift and CL instability focused on analyzing single wave pack- ets travelling through a background flow. Our analysis, based on computation of a four-dimensional spectrum with high temporal and spatial resolution, allowed us to study these phenomena in turbulent flows, and to identify direct evidence of their occurrence.
With these tools we showed that Doppler shift and CL instability occur naturally in a stratified disordered flow, as a result of the interaction of the waves with the horizontal winds. This indicates the CL instability can be one of the mechanism behind the formation of large scales structures in stratified flows, often observed in simulations but whose origin is unclear [27][28][29]. Moreover, although Doppler shift is observed in all forcing functions considered, development of the CL instability requires the external forcing not to disrupt the development of mean horizontal winds. Theories of stratified wave turbulence should take these effects into account. The mechanism, and the tools presented here, can be also relevant in quasi-geostrophic turbulence [33] and plasma turbulence [34,35], where zonal flows are also known to develop.