Geophysical turbulence and the duality of the energy flow across scales
GGeophysical turbulence and the duality of the energy flow across scales
A. Pouquet ∗ , and R. Marino Computational and Information Systems Laboratory, NCAR, Boulder CO 80307, USA. Department of Applied Mathematics, University of Colorado at Boulder, Boulder, CO 80309, USA.
The ocean and the atmosphere, and hence the climate, are governed at large scale by interactionsbetween pressure gradient, Coriolis and buoyancy forces. This leads to a quasi-geostrophic balancein which, in a two-dimensional-like fashion, the energy injected by solar radiation, winds or tidesgoes to large scales in what is known as an inverse cascade. Yet, except for Ekman friction, energydissipation and turbulent mixing occur at small scale implying the formation of such scales associatedwith breaking of geostrophic dynamics through wave-eddy interactions [1, 2] or frontogenesis [3, 4],in opposition to the inverse cascade. Can it be both at the same time? We exemplify here thisdual behavior of energy with the help of three-dimensional direct numerical simulations of rotatingstratified Boussinesq turbulence. We show that efficient small-scale mixing and large-scale coherencedevelop simultaneously in such geophysical and astrophysical flows, both with constant flux asrequired by theoretical arguments, thereby clearly resolving the aforementioned contradiction.
PACS numbers: 47.32.Ef, 47.55.Hd, 47.27.-i, 47.27.ek
Geostrophic balance, in which nonlinearities are ne-glected, leads to simplified quasi-bi-dimensional behaviorwith energy flowing to large scales, and reduced small-scale dissipation, contrary to observations [5]: verticalmixing can decrease water density, contributing to the(upward) closing of the ocean global circulation [1]. It isidentified with breaking of internal gravity waves [6], andit can potentially control the amplitude of the mesoscales.Such flows are neither three-dimensional (3D) nortwo-dimensional (2D), since at small scales, 3D eddiesmay prevail. Considering the system dimensionality D S proved essential when examining critical phenomenawhich simplify in higher dimensions, due to more modeinteractions as D S grows. Fluid turbulence is vastly dif-ferent in two or three dimensions, because of the strongconstraint imposed by the new 2D invariants (such as theintegrated powers of vorticity). This leads to energy flow-ing towards the largest scales, ending up in a condensate[7]; it can take the form of features such as jets, observedin the atmosphere of planets, or in the oceans as stri-ations [8]. Thus, geophysical turbulence is anisotropic,quasi-2D at large scale and quasi-3D at small scale [9].However, traditional three-dimensional homogeneousisotropic turbulence (HIT) is known to break structures(meso-scale eddies, clouds) into progressively smaller en-tities which will be dissipated at small scale, enhancingmixing of tracers such as pollutants [10] or biota [11].Whereas the fate of energy in 3D is modeled through anenhanced viscosity ν turb >
0, the 2D evolution leadingto large-scale structures can be related to a destabilizingtransport coefficient, e.g. ν turb ≤
0. Since the directionof the cascade is known to affect the amount of energyavailable to irreversible processes of dissipation and mix-ing, it is thus an essential parameter in the overall energybudget of the atmosphere and ocean [12].A transition from 2D to 3D in turbulence has beeninvestigated in various contexts. For example, is therea critical dimension for which ν turb changes sign, indica-tive of a change of behavior in the overall flow dynamics? Using two-layer quasi-geostrophic (QG) models with bot-tom friction, it was shown recently that when adding, ina somewhat ad-hoc fashion, a horizontal eddy-viscositymimicking coupling to smaller scales and thereby pre-sumably changing locally the sign of ν turb , both a directand inverse energy cascades were obtained [13].More formally, starting from two-point turbulence clo-sure, space dimensionality appears through incompress-ibility. The critical dimension that separates 2D from3D behavior can be computed and is found to be ≈ . R Π = | (cid:15) I /(cid:15) D | ;it is found to be a smooth monotonic function of D S , ina fashion similar to critical phenomena, thus providinga path between 2D and 3D behavior. In order to modelthe anisotropy of geophysical flows, one can alternativelyintroduce an anisotropic scale contraction/dilation. Thisallows to break the geostrophy constraint by consideringexplicitly the production of horizontal vorticity by hori-zontal or vertical eddies; it leads to a fractal dimensionof turbulence, close to 2.55 for stratified flows [17].Furthermore, an inverse energy cascade can also occurin 3D-HIT. On the one hand, when restricting nonlinearinteractions in 3D to those between helical waves of thesame polarization, energy is found to flow to large scale,with helicity (velocity-vorticity correlations) populatingthe small scales [18]. In reality, cross-polarization inter-actions dominate, but the tendency for strong inversetransfer is clearly displayed in this restricted model.On the other hand, taking a purely 2D input of energyand a fluid with a variable aspect ratio A r , energy againhas an increased tendency to flow to large scales as A r becomes small, with a transition at A r ≈ / A r is de-fined as the ratio of the vertical resolution to the forcing a r X i v : . [ phy s i c s . f l u - dyn ] S e p scale) [19]. A clear dual energy cascade obtains, with R Π a decreasing function of A r . Also, inverse transfer inthick layers (with now A r ≈ .
78) is observed experimen-tally, the suppression of vertical motions being attributedto interactions with vertical shear for eddies whose time-scale is larger than the characteristic shear time [20].These are idealized physical systems, modeling com-plex fluids under rather restrictive conditions. However,the link between large scales and small scales (or non-local interactions between Fourier modes) is embodiedin coherent structures such as chlorophyl filaments [21],water vapor, ozone, temperature or salinity tracer fronts,and in magnetohydrodynamics, current sheets, plasmoidsand Alfv´en vortices [22]. These structures have one di-mension comparable to the integral scale of the flow orlarger, and one close to the dissipative scale. One elementaltering the way such structures arise and evolve is theideal invariants, and in particular whether or not theyinvolve gradients. Finally, if one expects the symmetriesof the primitive equations to recover at small scale, us-ing a statistical argument based on the large number ofmodes, this recovery may be impeded by the presenceof large-scale shear [23]. For example, direct couplingbetween large scales (at which the inertio-gravity wavesreside) and small scales (at which turbulence resides) wasdemonstrated in [24], providing a progressive destructionof shear layers together with propagation, over the layerdepth, of efficient mixing induced by the turbulence.Stratified turbulence is not 2D in the traditional sense:it has strong vertical shearing [9, 25–29], allowing for theefficient creation of small scales, as well as of large scalesin the presence of rotation [30]. What is perhaps notwell recognized is that the 3D Boussinesq equations, in-cluding rotation and stratification as in the atmosphereand oceans, can produce both large scale and small scaleenergy excitation, both with constant flux. Numerousnumerical studies suffer from a lack of resolving both thelarge and the small eddies: because of the inherent costof such computations, a divide-and-conquer approach hasbeen successfully followed, analyzing either the direct orthe inverse cascade, but not convincingly both. Fluxes ofenergy to large scales and to small scales become compa-rable for strong rotation [31], as well as in the presenceof stratification [32]. However, in all these studies, thesmallness of the forcing wavenumber ( ≈ Methods:
Oceanic turbulence is studied in the idealizedcontext of the incompressible stably stratified rotatingBoussinesq primitive equations, with u the velocity and θ the density fluctuations in units of velocity. Solid-bodyrotation of strength Ω (with f = 2Ω ) is imposed in thevertical ( z ) direction with unit vector ˆ z , as well as anti-aligned gravity g ; isotropic three-dimensional forcing F Run
Re F r Ro
N/f R B R Π α
10a 5000 0.020 0.08 4 2.0 5.77 -3.9910b 5000 0.045 0.18 4 10.1 2.70 -2.9310c 5000 0.060 0.24 4 18.0 1.36 -2.3410d 4000 0.040 0.08 2 6.4 9.04 -3.9910e 5000 0.090 0.18 2 40.5 1.62 -2.1215a 8000 0.100 0.20 2 80.0 1.08 -1.87TABLE I. List of the runs done on cubic grids of n p points,with 10 & 15 standing for n p = 1024 & 1536 respectively.All runs use a random force in the wavenumber band k F ∈ [10 , Re, F r and Ro are the Reynolds, Froude and Rossbynumbers, with N/f = Ro/F r and R B = ReF r the buoyancyReynolds number. R Π = (cid:15) I /(cid:15) D is the ratio of the direct tothe inverse flux of energy in the vicinity of k F (1 < k < (cid:15) I , 11 < k <
20 for (cid:15) D ); it is computed on spectra averagedover 10 turn-over times τ NL = L F /U , in the range 12 22. Finally, α is the best fit for the small-scale kineticenergy spectral index; note the significant decrease of α withincreasing Re and R B . All large-scale indices, computed for k < k F , are close to 5 / is included; ∇ · u = 0 ensures incompressibility: ∂ t u − ν ∆ u + N θ ˆ z + F + ∇ p − f u × ˆ z = − u · ∇ u , (1) ∂ t θ − κ ∆ θ − N w = − u · ∇ θ , (2) w being the vertical velocity, p the pressure, ν the viscos-ity, and κ = ν the thermal diffusivity. The square Brunt-V¨ais¨al¨a frequency is given by N = − ( g/θ )( d ¯ θ /dz ),where d ¯ θ/dz is the imposed background stratification,assumed to be linear and constant. In the ideal case( ν = 0 , F = 0), the total (kinetic plus potential) energy E T = (cid:10) | u | + θ (cid:11) = E V + E P is conserved and thepoint-wise potential vorticity P V = − N ω z + f ∂ z θ + ω ·∇ θ is a material invariant. No modeling of small-scale dy-namics is included.The numerical code, GHOST (Geophysical High OrderSuite for Turbulence), uses a pseudo-spectral method andis tri-periodic, with n p grid points; it is parallelized witha hybrid MPI/Open-MP method and scales linearly up to98,000 processors for grid of up to 6144 points [33]. Forc-ing is introduced in the momentum equation as a randomfield centered in the wavenumber band k F ∈ [10 , L = 2 π ,corresponding to a minimum wavenumber k min =1; thesmallest resolved scale is 2 π/k max = 6 π/n p . Initial con-ditions are zero for the density and random for u .Three dimensionless parameters characterize the flow:the Reynolds number Re = U L F /ν , the Rossby num-ber Ro = U / [ L F f ] and the Froude number, F r = U / [ L F N ]; U is the rms velocity, L F = 2 π/k F is theforcing scale; finally, (cid:15) V ≡ dE V /dt = − (cid:104) u · F (cid:105) is thekinetic energy injection rate. Note that in order to re-solve the Ozmidov scale, at which the eddy turn-overtime and 1 /N become equal and isotropisation recovers,one can show that R B ≥ R B = ReF r is the a 200 400 600 800 1000 1200 1400200400600800100012001400 −0.6−0.4−0.200.20.40.6 b 200 400 600 800 1000 1200 1400200400600800100012001400 −0.6−0.4−0.200.20.40.6 FIG. 1. Horizontal (xy, (a) ) and vertical (xz, (b) ) two-dimensional cuts of the vertical velocity for run 15a at thelatest time, with R B ≈ 80 and a small-scale spectrum slightlysteeper than a Kolmogorov law. The axes are labeled in termsof grid spacing, and the forcing scale corresponds to roughly145 in these units. Observe the large-scale structures, with asize of up to a third of the overall flow (or more in the fila-ments), arising from the inverse cascade, together with super-imposed intense small scale eddies (e.g., at x ≈ , y ≈ buoyancy Reynolds number. Runs are performed with2 < R B ≤ 120 (see Table ). Whether the Ozmidov scaleis properly resolved or not may well alter the efficiencyof mixing, and the properties of stratified turbulence, asadvocated in [27] and as also observed here.The right-hand sides of equations (1, 2) are used to de-rive the evolution of the total (kinetic + potential) energydensity. Taking its Fourier transform (denoted by ˆ . , with (cid:63) denoting complex conjugate) gives access to the spec-tral transfer which, upon integration over wavenumber, yields the total isotropic energy flux Π T = Π V + Π P :Π V ( k ) = (cid:90) kk min T V ( q ) dq , T V ( q ) = − (cid:88) C q ˆ u (cid:63) q · (cid:92) ( u · ∇ u ) q with C q the shell q ≤ | q | < q + 1. An expression for Π P can be written in a similar fashion. Note that in theseBoussinesq runs, the eventual change of sign of energyfluxes at a “zero-crossing” wavenumber is given by k F since the forcing is added at that scale.Results: Fig.1 shows full 2D cuts of vertical velocityin the vertical and horizontal for run 15a; the forcingis roughly 1/10th of the box and one clearly observesboth intense small-scale features where dissipation oc-curs, and organized patches significantly larger than theforcing scale, indicative of the dual flux of energy.Results concerning scale-to-scale distribution inFourier space are displayed in Fig.2 for runs with N/f =2 and 4, with the fluxes Π T ( k ) (right) being averaged for10 turn-over times after the peak of dissipation t p ≈ . k < k F ), with a negative flux, and with an ap-proximate k − / scaling [30], as expected from classicaltheory of two-dimensional (2D) turbulence [7, 34]. Thisinverse cascade to large scales in 2D was demonstratedusing e.g. two-point closures of turbulence [35], or morerecently high-resolution numerical simulations [36].These runs also have a clear direct energy cascade ( k >k F ), with a constant positive flux. Spectral indices α aredefined through E V ( k ) ∼ k α where the fit is performedin the inertial range of wavenumber, k F < k < k diss with k diss ≈ k max marking the onset of the dissipation range.These exponents (see Table ), vary between ≈ . 99 and ≈ . 87, the steeper the lower Re and R B . The shallowerspectrum is close to a Kolmogorov solution α Kol = 5 / R B (see [31]for the rotating case).The inset in Fig. 2 gives the temporal variation of E V (solid lines) and (scaled) dissipation D V = 2 ν (cid:10) | ω | (cid:11) (dashed lines). The steady energy increase, after an ini-tial transient, is typical of inverse cascades; The variationof the ratio of inverse to direct flux with the buoyancyReynolds number is indicative of the increased effective-ness of turbulence as R B grows. One can expect thisratio to decrease as N/f increases since no inverse cas-cade occurs in the purely stratified case [30].Such direct cascades of energy in rotating stratifiedturbulence have been analyzed using theoretical closuremodels of turbulence [37]. Dual cascades were also foundwhen examining AVISO altimeter data for the Kuroshiocurrent [13], with values of R Π approaching those ofoceanic data for the largest imposed turbulent (horizon-tal) viscosity. Whereas these authors conclude to someambiguity in the interpretation of their results due tothe necessary filtering of the data, our DNS of the Boussi-nesq equations unambiguously show that dual energy cas-cades are realistic outcomes in a geophysical setting. The − − − k E v ( k ) Ev , 5 × (cid:105) (cid:157) | (cid:116) | (cid:156) k − a Re=4,000 Fr=0.04 Ro=0.08Re=5,000 Fr=0.09 Ro=0.18Re=8,000 Fr=0.1 Ro=0.2 (cid:190) k − (cid:190) k − (cid:190) k − (cid:111) NL − − − k (cid:87) T / (cid:161) V b Re=4,000 Fr=0.04 Ro=0.08 Re=5,000 Fr=0.09 Ro=0.18 Re=8,000 Fr=0.1 Ro=0.2 Re=5,000 Fr=0.02 Ro=0.08 Re=5,000 Fr=0.045 Ro=0.18 Re=5,000 Fr=0.06 Ro=0.24 FIG. 2. (Color Online.) (a): Kinetic energy spectra for run 10d (red), 10e (blue), 15a (black), all with N/f = 2 and increasing R B . The straight lines with different power laws are given as indications. In the bottom inset are shown the temporal evolutionof the kinetic energy for the same runs (solid lines), together with their (scaled) dissipation (dashed lines) 5 × ν (cid:10) | ω | (cid:11) , with ω = ∇ × u the vorticity. The spectra, not averaged in time, are shown at t/τ NL ∼ 22, whereas the peak of dissipation occursfor all the runs around t/τ NL ∼ . 3, time after which the energy starts to grow. (b): Total (kinetic plus potential) energyfluxes normalized by energy input (cid:15) V = (cid:104) u · F (cid:105) for the same runs, as well as run 10a (magenta-dash), 10b (green-dash) and 10c(cyan-dash) for which N/f = 4. higher values of R Π found in our runs likely reflect thefact that buoyancy is not dominant in our DNS, with N/f ≤ 4. However, we note that the abyssal southernocean at mid latitudes has N/f as low as 4 or 5 andshows considerable mixing [1, 38].Conclusion and discussion: We have shown in this pa-per that a dual (direct and inverse) constant flux en-ergy cascade is present in rotating stratified turbulence,thereby resolving the paradox noted by some authors(see, eg., [4, 13]) and thus adding credence to having bothgeostrophic balance and anomalous transport in geophys-ical turbulence. The computations clearly point out thepossibility of the co-existence in the ocean and the atmo-sphere of idealized large-scale dynamics dominated byquasi-geostrophic motions, together with the productionof small scales, essential to mixing [38].More computations and data analysis are required tocategorize in a quantitative way the mixing efficiencyone can expect in such flows. For example, the varia-tion of R Π with the relevant dimensionless parameters,is an open problem which requires huge numeral as wellas observational resources. Sub-grid scale modeling ofsmall-scale dynamics may be introduced to study thisphenomenon in a parametric fashion (see e.g. [39] forrotating flows). However, there are some indications ofa dual flux, using quasi-geostrophy [13], or in more com-plex settings using a numerical oceanic model applied tothe California coastal current [40]. This somewhat para- doxical behavior of the energy directivity can be under-stood if one recalls that triadic energetic exchanges canbe either positive or negative, and it is a delicate balancebetween the two that determines the overall sign of theflux, as also found for helical flows [18].Physical descriptions beyond the Boussinesq equationscan be used in modeling geophysical turbulence. For ex-ample, one can consider the evaporatively-driven (as op-posed to radiatively driven) configurations of stratocu-mulus clouds, in which case the buoyancy term is alteredby the existence of a threshold (in saturation mixturefraction), leading to a nonlinear equation of state. Simi-lar phenomena may occur in the oceans, for which thereis a complex set of state relations between temperature,density and salinity which may lead to distorted isopyc-nal surfaces. However, using the Boussinesq framework,it is clear that, beyond the energy cascades with small-scale or (exclusive) large-scale constant fluxes, other –mixed– solutions are found that explain how the oceanicand atmospheric systems are in quasi-geostrophic bal-ance at large scale and yet have a sufficient productionof small scales leading to enhanced mixing. We are thankful to B. Galperin and C. Herbert forfruitful discussions. This work is supported by NSFthrough grant CMG/1025183, and NCAR. Computa-tions were performed at NCAR (ASD), NSF/XSEDETGPHY-100029 & 110044, and INCITE/ DOE DE-AC05-00OR22725. [1] J. R. Ledwell, E. T. Montgomery, K. L. Polzin, L. St-Laurent, R. W. Schmitt, and J. M. Toole, Nature ,179 (2000).[2] J. Vanneste, Ann. Rev. Fluid Mech. , 147 (2013).[3] B. Hoskins and F. Bretherton, J. Atmos. Sci. , 11(1972).[4] M. Molemaker, J. McWilliams, and X. Capet, J. FluidMech. , 35 (2010).[5] G. Ivey, K. Winters, and J. Koseff, Ann. Rev. FluidMech. , 169 (2008).[6] M. Nikurashin, G. K. Vallis, and A. Adcroft, NatureGeosci. , 48 (2013).[7] R. Kraichnan and D. Montgomery, Rep. Prog. Phys. ,547 (1980).[8] B. Galperin, H. Nakano, H.-P. Huang, and S. Sukorian-sky, Geophys. Res. Lett. , L13303 (2004).[9] P. Sagaut and C. Cambon,Homogeneous Turbulence Dynamics (Cambridge Uni-versity Press, Cambridge, 2008).[10] B. I. Shraiman and E. Siggia, Nature , 639 (2000).[11] P. Klein and G. Lapeyre, Ann. Rev. Mar. Sci. , 351(2009).[12] R. Ferrari and C. Wunsch, Ann. Rev. Fluid Mech. ,253 (2009).[13] B. Arbic, K. Polzin, R. Scott, J. Richman, and J. Shriver,J. Phys. Oceano. , 283 (2013).[14] U. Frisch, M. Lesieur, and P. Sulem, Phys. Rev. Lett. , 895 (1976).[15] J. D. Fournier and U. Frisch, Phys. Rev. A , 747(1978).[16] T. Bell and M. Nelkin, Phys. Fluids , 345 (1977).[17] S. Lovejoy and D. Schertzer, Multifractal Cascades andthe Emergence of Atmospheric Dynamics, CambridgeUniversity Press, Cambridge (2012).[18] L. Biferale, S. Musacchio, and F. Toschi, Phys. Rev.Lett. , 164501 (2012).[19] A. Celani, S. Musacchio, and D. Vincenzi, Phys. Rev.Lett. , 184506 (2010).[20] H. Xia, D. Byrne, G. Falkovich, and M. Shats, NaturePhys. , 321 (2011). [21] A. Davis and X.-H. Yan, Geophys. Res. Lett. , L17304(2004).[22] D. Sundkvist, V. Krasnoselskikh, P. Shukla, A. Vaivads,M. Andr´e, S. Buchert, and H. R`eme, Nature , 825(2005).[23] A. Pumir and B. I. Shraiman, Phys. Rev. Lett. , 3114(1995).[24] D. C. Fritts, L. Wang, and J. Werne, Geophys. Res. Lett. , 396 (2009).[25] P. Billant and J.-M. Chomaz, Phys. Fluids , 1645(2001).[26] E. Lindborg, J. Fluid Mech. , 207 (2006).[27] G. Brethouwer, P. Billant, E. Lindborg, and J.-M.Chomaz, J. Fluid Mech. , 343 (2007).[28] M. Waite and P. Smolarkiewicz, J. Fluid Mech. , 239(2008).[29] P. Billant, A. Deloncle, J.-M. Chomaz, and P. Otheguy,J. Fluid Mech. , 396 (2010).[30] R. Marino, P. Mininni, D. Rosenberg, and A. Pouquet,EuroPhys. Lett. , 44006 (2013).[31] P. Mininni, D. Rosenberg, and A. Pouquet, J. FluidMech. , 263 (2012).[32] H. Aluie and S. Kurien, Eur. Phys. Lett. , 44006(2011).[33] P. Mininni, D. Rosenberg, R. Reddy, and A. Pouquet,Parallel Computing , 316 (2011).[34] G. Boffetta and R. Ecke, Ann. Rev. Fluid Mech. , 427(2012).[35] A. Pouquet, M. Lesieur, J. C. Andr´e, and C. Basdevant,J. Fluid Mech. , 305 (1975).[36] G. Boffetta, J. Fluid Mech. , 253 (2007).[37] B. Galperin and S. Sukoriansky, Ocean Dynamics ,1319 (2010).[38] K. J. Heywood, A. N. Garabato, and D. Stevens, Nature , 1011 (2002).[39] A. Pouquet, J. Baerenzung, P. Mininni, D. Rosenberg,and S. Thalabard, J. Phys., Conf. Series, Europ. Turb.Conf. Proc. ETC13, K. Bajer Ed. , 042015 (2011).[40] X. Capet, J. McWilliams, M. Molemaker, andA. Shchepetkin, J. Phys. Ocean.38