The decay of turbulence in rotating flows

0.91 MB·English

by  T. Teitelbaum

The decay of turbulence in rotating flows
Tomas Teitelbaum1 and Pablo D. Mininni1,2
1 Departamento de F´ısica, Facultad de Ciencias Exactas y Naturales,
Universidad de Buenos Aires, IFIBA and CONICET,
Ciudad Universitaria, 1428 Buenos Aires, Argentina.
1
1 2 NCAR, P.O. Box 3000, Boulder, Colorado 80307-3000, U.S.A.
0 (Dated: January 25, 2011)
2
Wepresentaparametricspacestudyofthedecayofturbulenceinrotatingflowscombiningdirect

n numerical simulations, large eddy simulations, and phenomenological theory. Several cases are
a considered: (1)theeffectofvaryingthecharacteristic scaleoftheinitialconditionswhencompared
J with the size of the box, to mimic “bounded” and “unbounded” flows; (2) the effect of helicity

4 (correlation between the velocity and vorticity); (3) the effect of Rossby and Reynolds numbers;
2 and(4)theeffectofanisotropyintheinitialconditions. InitialconditionsincludetheTaylor-Green
vortex, the Arn’old-Beltrami-Childress flow, and random flows with large-scale energy spectrum

] proportionaltok4. Thedecaylawsobtainedinthesimulationsfortheenergy,helicity,andenstrophy
n
in each case can be explained with phenomenological arguments that consider separate decays
y
for two-dimensional and three-dimensional modes, and that take into account the role of helicity
d
and rotation in slowing down the energy decay. The time evolution of the energy spectrum and
-
u developmentofanisotropies inthesimulationsarealsodiscussed. Finally,theeffectofrotation and
l helicity in theskewness and kurtosis of the flow is considered.
f
.
s
c
i I. INTRODUCTION dimensional(3D)modes,andevolveundertheir ownau-
s
tonomous dynamics. Moreover, in that limit the av-
y
eraged equation for the slow modes splits into a 2D
h Naturepresentsseveralexamplesofrotatingflows. Ro-
p Navier-Stokes equation for the vertically-averaged hor-
tationinfluenceslarge-scalemotionsintheEarth’satmo-
[ izontal velocity and a passive scalar equation for the
sphereandoceans,aswellasconvectiveregionsofthesun
vertically-averaged vertical velocity. Although simula-
2 and stars. Rotation is also important in many industrial
tions of forced rotating flows [6] and of ideal truncated
v flows, such as turbo machinery, rotor-craft,and rotating
5 channels. In a rotating system, the Coriolis force, lin- rotatingflows[15,16]usingperiodicboundaryconditions
9 show good agreement with these predictions for small
ear in the velocity, modifies the flow nonlinear dynamics
5 enough Rossby numbers, for long times the decoupling
when strong enough. In its presence, the Navier-Stokes
1 ofslowand fastmodes seems to break down. Also, some
equation becomes a multi-scale problem with a “slow”
.
9 authors[17]arguethatinunboundeddomainsnodecou-
time scale τ ∼ L/U associated with the eddies at a
0 L pling is achievable even for Ro=0.
characteristicscaleL(U isacharacteristicvelocity),and
0
1 a “fast” time scale τΩ ∼ 1/Ω ∼ τLRo associated with Ofparticularimportanceformanygeophysicalandas-
: inertial waves. The dimensionless Rossby number Ro is trophysical problems is how turbulence decays in time.
v
the ratio of advection to Coriolis forces, and measures The problem is also important for laboratory experi-
i
X the influence of rotationupon the nonlinear dynamics of ments,asitprovides,e.g.,onewaytomeasurechangesin
the system (decreasing as rotation becomes dominant). the energy dissipation rate associated with the presence
r
a Resonant wave theory provides a framework to study ofrotation. Evenintheabsenceofrotation,thedecayof
the effect of rapid rotation in turbulence [1–5]. The sep- isotropic turbulence proves difficult to tackle because of
arationbetweenthefastandslowtime scalesresultsina thedifferentdecaylawsobtaineddependingonboundary
selection of the resonant triadic interactions as the ones andinitialconditions. As anexample,for boundedflows
responsible for the energy transfer among scales. As a (i.e., flows for which the initial characteristic length is
result, energy transfer and dissipation is substantially closetothesizeofthevessel)the energydecaysas∼t−2
decreased in the presence of strong rotation [3]. The (see, e.g., [13, 18, 19]). For unbounded flows (i.e., flows
resonant condition is also responsible for the transfer of in an infinite domain, or in practice, flows for which the
energy towards two-dimensional (2D) slow modes, driv- initialcharacteristiclength is muchsmallerthanthe size
ing the flow to a quasi-2D state [3, 4] (this result is of- of the vessel) a ∼ t−10/7 [20, 21] or ∼ t−6/5 [13, 22, 23]
ten referred to in the literature as the “dynamic Taylor- decay law is observed depending on whether the initial
Proudman theorem”, see e.g., [6]). The development of energyspectrumatlargescalesbehavesas∼k4 or∼k2,
anisotropy and reduction of the energy transfer and dis- respectively.
sipation rates has been verified in numerical simulations In the presence of rotation, the decay of turbulence
[7–11] and experiments [12, 13]. becomes substantially richer, with decay rates depend-
Similar arguments (see, e.g., [14]) indicate that in ing not only on whether turbulence is bounded or un-
the limit of fast rotation (small Rossby number) the bounded and on the initial spectrum at large scales, but
slow 2D modes decouple from the remaining fast three- also, e.g., on the strength of background rotation (see
2
[24]). A detailed experimental study of this dependence The structure of the paper is as follows. In Sec. II we
canbe foundinRefs. [12,13],wherethe authorsstudied introducethe equations,describetheDNSandLES,and
the energy decay of grid-generated turbulence in a ro- give details of the initial conditions and different param-
tating tank using particle image velocimetry, and found eters used. Section III presents phenomenological argu-
different decay laws depending upon the rate of rotation ments to obtain decay laws in turbulent flows with and
andthesaturation(ornot)ofthecharacteristicsizeofthe without rotation. The phenomenological predictions are
largest eddies. For large Rossby number they reported then compared with the numerical results in Sec. IV,
a decay ∼ tα with exponent α ≈ −1.1 for non rotating which presents the energy, helicity, and enstrophy decay
turbulence (consistent with the value of −6/5 predicted in all runs. The spectral evolution and development of
in [22, 23] for ∼ k2 unbounded turbulence), which later anisotropy in the flows is discussed at the end of that
turned to α ≈ −2 after the largest eddies grew to the section. The effect of initial anisotropy in the decay is
experiment size. For small Rossby number this energy consideredinSec.V. Astatisticalanalysisincludingevo-
decay rate became smaller saturating at ≈ −1. Simi- lution of skewness and kurtosis is presented in Sec. VI.
lar results were reported in simulations in [10], were a Finally, Sec. VII gives our conclusions.
decrease from ≈ −10/7 to −5/7 was observed as rota-
tion wasincreased. These results are consistentwith the
reduction of the energy transfer discussed above. They II. EQUATIONS AND MODELS
also reported a steeper energy spectrum together with
positive vorticity skewness for the rotating flows, and A. Equations
anisotropic growth of integral scales (see also [7, 25]).
The decay rate of rotating turbulence also seems to
The evolution of an incompressible fluid in a rotating
depend on the helicity content of the flow. Helicity (the
frame is described by the Navier-Stokes equation with
alignment between velocity and vorticity) is an ideal in-
the Coriolis force,
variantof the equations of motion (a quantity conserved
in the inviscid limit) with intriguing properties. In the ∂ u+ω×u+2Ω×u=−∇P+ν∇2u, (1)
t
absence of rotation, the presence of helicity does not
modify theenergydecayrate[19,26,27]northedissipa- and the incompressibility condition,
tion rate [28]. However, in the presence of rotation [27]
reported a further decrease of the energy transfer when ∇·u=0, (2)
both rotation and helicity were present, and [19] showed
that the decay rate of bounded rotating flows changes where u is the velocity field, ω =∇×u is the vorticity,
drastically depending on whether helicity is present or the centrifugal term is included in the total pressure per
not. unit of mass P, and ν is the kinematic viscosity (uni-
In this paper we conduct a detailed study of param- form density is assumed). The rotation axis is in the z
eter space of rotating helical flows, taking into account direction, so Ω=Ωzˆ, where Ω is the rotation frequency.
(1) the effect of varying the characteristic scale of the As mentioned in the introduction, these equations are
initial conditions when compared with the size of the solved numerically using two different methods: DNS,
box, (2) the effect of helicity, (3) the effect of Rossby andLES using a dynamicalsubgrid-scalespectralmodel
and Reynolds numbers, and (4) the effect of anisotropy of rotating turbulence that also takes into account the
in the initial conditions. The numerical study uses a helicity cascade. All simulations were performed in a
two pronged approachcombining direct numerical simu- threedimensionalperiodicboxoflength2π,usingdiffer-
lations(DNS) andlargeeddy simulations(LES).Several ent spatial resolutions ranging from 963 grid points for
initialconditionsareconsidered,althoughwhenthechar- thelowestresolutionLESrunsupto5123 forthehighest
acteristic initial scale of the flow is smaller than the size resolution DNS.
ofthe domain, wefocus only onthe casewith largescale
initial energy spectrum ∼ k4. The different decay laws
obtained(whichinsomecasescoincidewithpreviousex- B. Models
perimentalornumericalobservations,whileinothersare
new) are explained using phenomenological arguments, In a DNS, all spatial and time scales (up to the dis-
and we classify the results depending on the relevant ef- sipation scale) are explicitly resolved. The simulations
fects on each case. wereperformedusingaparallelizedpseudo-spectralcode
After studying the decay laws, we study the evolution [29, 30] with the two-third rule for dealiasing. As a re-
of anisotropy, of skewness and kurtosis, and the forma- sult, the maximum wave number resolved in the DNS is
tion of columnar structures in the flow. We consider k =N/3whereN isthelinearresolution;toproperly
max
how helicity affects the evolution of skewness and kur- resolve the dissipative scales the condition k /k < 1
η max
tosis, andassociatepeaks observedin the time evolution must be satisfied during all simulations, where k is the
η
of these quantities with the dynamics of the columnar dissipation wave number. In practice, this condition is
structures. more stringent if reliable data about velocity gradients
3
TABLE I: Parameters used in the simulations: kinematic viscosity ν, rotation rate Ω, Reynolds number Re, Rossby number
Ro, micro-Rossby number Roω, initial relative helicity h, relative helicity at the time of maximum dissipation h∗, and time
of maximum dissipation t∗. The values of Re, Ro, and Roω are always given at t∗. The last column succinctly describes the
initial energy spectrum E(k): the power law followed by the spectrum, the range of scales where this power law is satisfied,
and the flow (TG for Taylor-Green, ABC for Arn’old-Beltrami-Childress, and RND for random).
Run ν Ω Re Ro Roω h h∗ t∗ Initial E(k)
D256-1 1.5×10−3 0 450 ∞ ∞ 0 9×10−10 1.26 k−4 (4-14) TG
D256-2 1.5×10−3 4 550 0.12 1.28 0 −1×10−8 1.06 k−4 (4-14) TG
D256H-1 1.5×10−3 0 600 ∞ ∞ 0.95 0.34 2.28 k−4 (4-14) ABC
D256H-2 1.5×10−3 4 830 0.08 0.80 0.95 0.65 2.25 k−4 (4-14) ABC
D512-1 7×10−4 4 1100 0.12 1.82 0 7×10−9 0.88 k−4 (4-14) TG
D512-2 8.5×10−4 0 420 ∞ ∞ 8×10−5 8×10−4 0.60 k4 (1-14) RND
D512-3 8.5×10−4 10 450 0.10 0.95 4×10−3 4×10−3 0.70 k4 (1-14) RND
D512H-1 7×10−4 4 1750 0.08 1.15 0.95 0.44 1.70 k−4 (4-14) ABC
D512H-2 8×10−4 0 440 ∞ ∞ 0.90 0.38 0.94 k4 (1-14) RND
D512H-3 8×10−4 10 530 0.07 0.70 0.99 0.5 1.50 k4 (1-14) ABC
L96-1 8.5×10−4 0 550 ∞ ∞ 0.03 0.02 0.30 k4 (1-14) RND
L96-2 8.5×10−4 2 540 0.42 2.90 −0.03 −0.02 0.30 k4 (1-14) RND
L96-3 8.5×10−4 4 540 0.21 1.45 −0.03 −0.02 0.30 k4 (1-14) RND
L96-4 8.5×10−4 6 550 0.14 0.95 −0.03 −0.02 0.30 k4 (1-14) RND
L96-5 8.5×10−4 8 550 0.11 0.73 −0.03 −0.02 0.30 k4 (1-14) RND
L96-6 8.5×10−4 10 530 0.08 0.65 0.03 0.02 0.30 k4 (1-14) RND
L96H-1 8×10−4 0 500 ∞ ∞ 0.90 0.51 0.70 k4 (1-14) RND
L96H-2 8.5×10−4 10 540 0.08 0.63 0.90 0.70 0.70 k4 (1-14) RND
L96H-3 8.5×10−4 10 490 0.08 0.60 0.99 0.80 1.15 k4 (1-14) ABC
L192-1 2×10−4 0 1200 ∞ ∞ −7×10−3 −6×10−3 0.10 k4 (1-30) RND
L192-2 2×10−4 10 1100 0.22 1.65 −7×10−3 −6×10−3 0.13 k4 (1-30) RND
L192H-1 2×10−4 0 950 ∞ ∞ 0.90 0.60 0.30 k4 (1-30) RND
L192H-2 2×10−4 10 1000 0.20 1.60 0.94 0.71 0.38 k4 (1-30) ABC
L192HA-1 2×10−4 10 1200 0.16 1.40 0.90 0.56 0.50 k4 (1-25) RND
L192HA-2 2×10−4 10 1300 0.14 1.35 0.90 0.59 0.46 k4 (1-25) RND
L192HA-3 2×10−4 10 1300 0.15 1.35 0.90 0.58 0.45 k4 (1-25) RND
andhigh-orderstatisticsoftheflowareneeded(see,e.g., scales of scales smaller than a cut-off scale is modeled
[31],wheretheyindicatek /k <0.5forsuchstudies). withsimplifiedequations. Tothisendweusethespectral
η max
modelderivedin[32]forisotropichelicalandnon-helical
Thedissipationwavenumberasafunctionoftimewas
turbulence, and its extension to the rotating case [33].
computedforallsimulationsintwodifferentways: asthe
The model is based on the eddy damped quasi-normal
Kolmogorov dissipation wave number for isotropic and
homogeneous turbulence k = (ǫ/ν3)1/4 = (hω2i/ν2)1/4 Markovian (EDQNM) closure to compute eddy viscos-
η
(whereǫistheenergydissipationrateandhω2ithemean ityandeddynoise,andassumesunresolvedscales(scales
smallerthanthecut-off)areisotropic. Botheddyviscos-
square vorticity), and as the wavenumber where the en-
ity and noise are computed considering the contribution
strophy spectrum peaks. The Kolmogorov dissipation
from the energy and the helicity spectra (see, e.g., [34]
wave number was found to be always larger than the
for another subgrid model that takes into account the
wavenumber where dissipationpeaks, andin the follow-
effect of helicity). The model adapts dynamically to the
ing we therefore only consider the Kolmogorov scale as
inertial indices of the resolved energy and helicity spec-
it gives a more stringent condition on the resolution.
tra, and as a result it is well suited to study rotating
In all DNS discussed below, the ratio k /k was
η max turbulence for whichthe scalinglawsare notwellknown
≤ 0.7 at the time of maximum dissipation (t ≈ 1 to
andmaydependontheRossbynumber. Foravalidation
t ≈ 3 depending on the simulation), 0.2 to 0.5 at t ≈ 10
of the LES againstDNS results, the reader is referredto
(when the self-similar decay starts in most runs), and
[32, 33].
monotonously decreases to values between 0.05 to 0.2 at
t≈100. Thespatialresolutionsusedwere2563 and5123
The subgrid model starts by applying a spectral fil-
grid points, and an explicit second-order Runge-Kutta
ter to the equations; this operation consists in truncat-
method was used to evolve in time, with a Courant- ing all velocity components at wave vectors k such that
Friedrichs-Levy (CFL) number smaller than one. |k| = k > k , where k is the cut-off wave number. One
c c
In the LES approach, only the large scales are explic- then models the transfer between the large (resolved)
itly resolved,while the statistical impact onthe resolved scales and the small (subgrid unresolved) scales of the
4
flowbyaddingeddyviscosityandeddynoisetotheequa- wherethebracketsdenotespatialaverage. Wealsodefine
tions for the resolved scales. These are obtained solving the relative helicity as
the EDQNMequationsfor estimatedenergyandhelicity
H
spectra in the subgrid range. To this end, an interme- h= , (6)
h|u||ω|i
diate range is defined, lying between k′ and k (in most
c c
cases kc′ =kc/3), where the energy spectrum is assumed which is bounded between −1 and 1 and can be inter-
to present a power-law behavior possibly followed by an preted as the mean cosine of the angle between the ve-
exponential decrease. As an example, for the energy the locity and the vorticity.
following expression is used: To control the net amount of relative helicity in the
initial conditions we consider three different flows: a su-
E(k,t)=E0k−αEe−δEk, kc′ ≤k<kc . (3) perposition of Taylor-Green(TG) vortices [35],
The coefficients αE, δE and E0 are computed at each uTG = Usin(k0x)cos(k0y)cos(k0z)xˆ−
time step doing a mean square fit of the resolved energy
spectrum. The spectrum is extrapolated to the unre- Ucos(k0x)sin(k0y)cos(k0z)yˆ, (7)
solved scales using these coefficients, and the EDQNM a superposition of Arn’old-Beltrami-Childress (ABC)
equations are solved. Then one solves the Navier-Stokes flows [36],
equation (1) with an extra term on the r.h.s. which in
spectral space takes the form uABC = [Bcos(k0y)+Asin(k0z)]xˆ+
[Ccos(k0z)+Asin(k0x)]yˆ+
−ν(k|k ,t)k2u(k,t), (4)
c [Acos(k0x)+Bsin(k0y)]zˆ, (8)
where u(k,t) is the Fourier transform of the velocity
(withA=0.9,B =1.1,andC =1),andasuperposition
field u(x,t), −k2 is the Laplacian in Fourier space, and
ofFouriermodeswithrandomphases(RND)inwhichwe
ν(k|k ,t)isaneddyviscosityproportionaltotheratioof
c usethealgorithmdescribedin[37]tocontroltherelative
the so-calledabsorptionterms in EDQNM to the energy
helicity. Ineachcase,theflowsaresuperposedinarange
(andhelicity)spectrum. Eddynoiseisaddedinasimilar
ofwavenumbersasdescribedinTableI,andwithglobal
manner (for more details, see [32, 33]).
amplitudesforeachwavenumbertogivethedesiredslope
LES runs using this model have a resolution of either
in the initial energy spectrum.
963or1923gridpoints. Apseudo-spectralmethodisalso
The TG vortex is non-helical (h = 0) and has no en-
used, but without dealiasing, resulting in the maximum
ergy in the k = 0 modes, whose amplification in the
z
wave number k = N/2. As in the DNS, an explicit
max rotatingcases(see below)canthus be attributedonly to
second-order Runge-Kutta method is used to evolve in
abidimensionalizationprocess. TheTGvortexwasorig-
time.
inally motivated as an initial condition which leads to
Parameters for all sets of runs are listed in Table I.
rapid development of small spatial scales, and also mim-
DNS runs are labeled with a D, followed by the linear
ics the von K´arm´an flow between two counter-rotating
resolution,aletter“H”iftherunhashelicity,aletter“A”
disks used in several experiments [38]. The ABC flow is
if the initial energy spectrum is anisotropic,and the run
an eigenfunction of the curl operator and as a result has
number. LES runs start with an L, followed by numbers
maximum helicity (h = ±1 depending on the sign of k0,
and letters using the same convention as in the DNS.
when only one value of k0 is excited), whereas the RND
flow allows us to tune the amount of initial relative he-
licity between −1 and 1 as well as the initial anisotropy.
C. Initial conditions and definitions
Whengeneratingthe flows,twodifferentinitialenergy
spectrawereconsidered. Tostudyinitialconditionswith
As mentioned in the introduction, we are interested characteristic length close to the size of the computa-
in the decay laws obtained in the system depending on tional domain, a spectrum E(k) ∼ k−4 for k ∈ [4,14]
properties of the initial conditions and the amount of (followed by exponential decay) was imposed. To study
rotation. Inparticular,wewillvarytheinitialamountof initial conditions with length smaller than the domain
helicity, the initial energy containing scale (with respect size, we imposed a spectrum E(k) ∼ k4 for k ∈ [1,14],
to the largest available scale in the box), the shape of [1,25], or [1,30](also followed by exponential decay). In
the energyspectrum, andthe strengthofturbulence and thelattercase,thecharacteristiclengthcangrowintime
of rotation as controlled by the Reynolds and Rossby as the spectrum peaks around k = 14, 25, or 30. This
numbers respectively. allows us to mimic (at least for a finite time before the
HelicityisanidealinvariantoftheNavier-Stokesequa- characteristic length reaches the domain size) the decay
tion which measures the alignment betweenvelocity and of unbounded flows. The characteristic length will be
vorticity. If zero, the initial conditions are mirror- associated in the following with the flow integral scale,
symmetric, so it also measures the departure from a which is defined as
mirror-symmetric state. We define the net helicity as
k−1E(k)
L=2π k , (9)
H =hu·ωi, (5) E(k)
P k
P
5
where E(k) is the isotropic energy spectrum. spectrum that depends only on k and k , which relates
k ⊥
Simulations in Table I are alsocharacterizedby differ- to the reduced energy spectra as follows:
ent Reynolds and Rossby numbers. The Reynolds and
Rossby numbers in the Table are defined as E(kk)= e(kk,k⊥), (15)
Xk⊥
UL
Re= (10) and
ν
E(k )= e(k ,k ). (16)
and ⊥ k ⊥
Xkk
U
Ro= , (11)
2ΩL
respectively. Of importance is also the micro-scale III. TIME EVOLUTION - PHENOMENOLOGY
Rossby number (see e.g., [7])
Wepresentnowphenomenologicalargumentsthatwill
ω
Roω = , (12) become handy to understand the different decay rates
2Ω
thatareobservedinoursimulationsaswellasinprevious
which can be interpreted as the ratio of the convective studies. Some of the arguments are well known, while
to the Coriolis acceleration at the Taylor scale. The others are new, and we quote previous derivations when
Rossby number Ro must be small enough for rotation needed.
to affect the turbulence, while the micro-Rossbynumber
Roω mustbe largerthanonefor scramblingeffects ofin-
A. Non-rotating flows
ertial wavesnot to completely damp the nonlinear term,
which would lead to pure exponential viscous energy de-
cay [5]. In all runs in Table I, Ro and Roω are one order 1. Bounded
of magnitude apart at the time of maximum enstrophy
t∗, and this interval is roughly preservedthroughout the From the energy balance equation
simulations.
dE
Here and in the following, the isotropic energy spec- ∼ǫ (17)
dt
trumisdefinedbyaveraginginFourierspaceoverspher-
ical shells where ǫ is the energy dissipation rate, Kolmogorov phe-
nomenology leads to
1
E(k,t)= 2 u∗(k,t)·u(k,t), (13) dE E3/2
k≤|kX|<k+1 dt ∼ L , (18)
whereu(k,t)istheFouriertransformofthevelocityfield, where E =E(t) ∼ kE(k) and L is an energy containing
and the asterisk denotes complex conjugate. Other two length scale. For bounded flows where L ∼ L0 (L0 is
spectra can also be used to characterize anisotropy. thesizeofthebox),Eq.(18)becomesdE/dt∼E3/2/L0,
Onthe onehand,the so-called“reduced”energyspec- resulting in the self-similar decay [18, 39, 40]
tra E(k ) and E(k ) are defined averaging in Fourier
⊥ k E(t)∼t−2. (19)
spaceovercylindersandplanesrespectively. Morespecif-
ically, the reduced energy spectra as a function of wave
numbers k with k = (k ,k ,0), and k with k =
⊥ ⊥ x y k k 2. Unbounded
(0,0,k ), are defined by computing the sum above over
z
all modes in the cylindrical shells k ≤ |k | < k +1
⊥ ⊥ ⊥
andoverplanesk ≤|k |<k +1respectively(isotropic In unbounded flows, a similarity solution of Eq. (18)
k k k
requires some knowledge of the behavior of the energy
and reduced spectra for the helicity are defined in the
containing scale L, which is in turn related to the evo-
sameway). Fromthereducedspectra,perpendicularand
lution of E(k) for low wave numbers. In the case of an
parallel integral scales can be defined; e.g., for the per-
initial large scale spectrum ∼k4, the quasi-invariance of
pendicular direction,
Loitsyanski’s integral I (see [21, 41]) leads, on dimen-
k−1E(k ) sionalgrounds,toI ∼L5E,andreplacinginEq.(18)we
L⊥ =2π k⊥ ⊥ ⊥ . (14) get Kolmogorov’sresult [20]
E(k )
P k⊥ ⊥
E(t)∼t−10/7. (20)
On the other hand, morPe information of the spectral
anisotropy can be obtained studying the axisymmetric A different decay law is obtained if an initial ∼k2 spec-
energy spectrum e(k ,k ) (see, e.g., [3, 5]). Assuming trum is assumed for low wave numbers [22, 23]. In the
k ⊥
theflowisaxisymmetric,thethree-dimensionalspectrum following,we will consider only the bounded or the ∼k4
canbe integratedaroundthe axisofrotationto obtaina unbounded cases.
6
FIG. 2: Energy decay for different values of Ω from 0 to 10
for unbounded, non-helical runs L96-1, L96-2, L96-3, L96-4,
L96-5, and L96-6. The decay becomes slower with increasing
rotation rate. We also show t−10/7 and t−5/7 laws as refer-
ences.
B. Rotating flows: isotropic arguments
1. Bounded
In the case of solid-body rotationwithout net helicity,
a spectra E(k) ∼ ǫ1/2Ω1/2k−2 is often assumed at small
scales (i.e., wave numbers larger than the integral wave
number) [9, 42–45]. Replacing in the balance equation,
this spectrum leads to
dE E2
∼ . (21)
dt L2Ω
For bounded flows L∼L and we get [13, 44, 46]
o
E(t)∼t−1. (22)
In helical rotating flows the small-scale energy spec-
trum takes a different form. The direct transfer is dom-
inated by the helicity cascade. In this case we can write
the helicity flux as δ ∼h /(Ωτ2) where h is the helicity
l l l
at the scale l, and τ the eddy turnover time [11]. Con-
l
stancy of δ leads to small scale spectra E(k)∼ k−n and
FIG. 1: (a) Energy decay for non-rotating unbounded runs.
H(k)∼kn−4,wheren=5/2obtainsforthecaseofmax-
Non-helical runs D512-2 (solid), L96-1 (dashed), L192-1
(dash-dotted),andhelicalrunsD512H-2(solid,thick),L96H- imum helicity [11]. Further use of dimensional analysis
1 (dashed, thick), and L192H-1 (dash-dotted, thick) are leads to E(k) ∼ ǫ1/4Ω5/4k−5/2 for the energy spectrum
shown. A −10/7 slope is shown as a reference. The inset in terms of the energy dissipation rate, and replacing in
showstheenergydecayfornon-rotatingboundedrunsD256- the balance equation we get
1 (solid), and D256H-1 (solid, thick). (b) Enstrophy decay
for the same runs, with a −17/7 slope shown as a reference. dE E4
Theinsetshowstheenstrophydecayintheboundedruns. (c) ∼ . (23)
dt L6Ω5
Helicity decay in theunboundedhelical runsof (a);theinset
showsthehelicitydecayfortheboundedhelicalrunD256H-1.
For L∼L then [19]
o
E(t)∼t−1/3. (24)
7
2. Unbounded
For non-helical unbounded flows with E(k) ∼ k4 at
small wave numbers we can again make use of the con-
stancy of I in Eq. (21), leading to [25]
E(t)∼t−5/7. (25)
For helical flows, assuming I remains constant in
Eq. (23), we obtain [47]
E(t)∼t−5/21. (26)
C. Rotating flows: anisotropic arguments
ThedecaylawsobtainedforrotatingflowsinEqs.(22),
(24), and (25), have been reported in experiments or in
simulations [12, 13, 19, 24]. However,the analysis above
is based on the isotropic energy spectrum and on the
quasi-invarianceoftheisotropicLoitsyanskiintegral. For
an anisotropic flow, other quantities are expected to be
quasi-invariants during the decay instead [41, 48, 49].
Rotating flows tend to become quasi-2D, and the as-
sumptionofanaxisymmetricenergyspectrumseemsnat-
ural considering the symmetries of the problem. If there
is no dependence on wavenumbers on the paralleldirec-
tion, the energy spectrum for small values of k can be
⊥
expanded as (see, e.g., [50])
E(k⊥)≈Lk⊥−1+Kk⊥+I2Dk⊥3 +··· (27)
We will be interested in the following coefficients:
K = hu·u′irdr, (28)
Z
and
I2D = hu·u′ir3dr, (29)
Z
where hu·u′i is the two-point correlation function for
spatial increments r perpendicular to the rotation axis.
If the correlation function decays fast enough for large
FIG. 3: (a) Energy decay for rotating unbounded runs (Ω= values of r [48], these quantities can be expected to be
10). Non-helical runs D512-3 (solid), L96-6 (dashed), L192-2 quasi-invariants during the decay, respectively for initial
(dash-dotted),andhelicalrunsD512H-3(solid,thick),L96H- large-scale energy spectra ∼ k and ∼ k3, in the same
⊥ ⊥
3 (dashed, thick), L96H-2 (dash-tripe-dotted, thick), and way I is quasi-conserved during the decay of isotropic
L192H-2 (dash-dotted, thick) are shown. At late times, the flowswithaninitiallarge-scale∼k4 energyspectrum. A
non-helicalrunsdecayslightlyfasterthant−5/7,whilethehe- detailedproofofthe conservationofK for rotatingflows
licalrunsareclosetoa−1/3decay. Theinsetshowsbounded
canbefoundin[49];itisadirectconsequenceofthecon-
non-helicalrunsD256-2(solid)andD512-1(dashed),andhe-
servationoflinearmomentuminthe directionparallelto
licalrunsD256H-2(solid,thick)andD512H-1(dashed,thick).
(b) Enstrophy decay for the same runs, with a −12/7 slope the rotation axis. It is worth pointing out that these
quantities were also shown to be conserved in other sys-
shown as a reference. The inset shows the enstrophy decay
in the bounded runs. (c) Helicity decay in the unbounded tems: proofs of the conservationof K and I2D for quasi-
helical runs; boundedhelical runsare shown in theinset. geostrophicflows can be found in [50]. In practice, these
quantitiesareonlyapproximatelyconservedinnumerical
8
FIG. 5: Evolution of I/I(0) (thick lines) and of I2D/I2D(0)
(thinlines)forrunsD512-3(solid)andL96-6(dashed). While
in both runs I2D mantains an approximately constant value,
I growths monotonically and during the self-similar energy
decay increases by one order of magnitude.
simulations, see e.g., the approximate constancy of I2D
and K reported for rotating flows in [47].
AspervirtueofthedecaytheRossbynumberdecreases
withtime, wewillfurtherassumeforourphenomenolog-
icalanalysisthat2Dand3Dmodesareonlyweaklycou-
pled,andwriteequationsfortheenergyinthe2Dmodes,
E2D. Inthenon-helicalcase,ifK remainsapproximately
constant with K ∼ E2DL2⊥L0k (where L0k is the size of
the box in the direction parallelto Ω), then Eq. (21) for
the 2D modes becomes
dE2D ∼ E23DL0k, (30)
dt KΩ
which leads to a decay
E2D(t)∼t−1/2. (31)
Alternatively,constancyofI2D ∼E2DL4⊥L0k inEq.(21)
for the 2D modes leads to
5/2 1/2
dE2D E2DL0k
∼ , (32)
dt I1/2Ω
2D
and
E2D(t)∼t−2/3. (33)
The same arguments can be extended to the helical
FIG. 4: Energy decay for E3D (solid) and E2D (dashed) for
rotatingcaseusingEq.(23). IfconstancyofKisassumed
runswith rotation: (a)Unboundednon-helicalD512-3(thin)
andL96-6 (thick);E3D ∼t−10/7 andE2D ∼t−2/3 decaysare we get
idnidtiiocnasteLd9.6H(b-)2;UEn2bDouinsdceldosehetloicat−l1w/2it.h(rca)nUdonmbouinnidtieadl hcoenli-- dE2D ∼ E27DL30k, (34)
cal with ABC initial conditions D512H-3 (thin) and L96H-3 dt Ω5K3
(thick);E2D isclosetot−1/3. (d)Boundedhelical withABC and
initial conditions D512H-1.
E(t)∼t−1/6. (35)
9
Finally, constancy of I2D leads to Asaresult,the“bounded”runsarehereonlybrieflycon-
sidered to study the time evolution of global quantities
dE2D E21D1/2L30/k2 (energy, enstrophy, and helicity), and to compare with
∼ , (36) thepredictionobtainedinthecorrespondingcasesinthe
dt I3/2Ω5
2D phenomenological analysis.
and
E(t)∼t−2/9. (37) A. Non-rotating flows
These decay laws will be important to analyze the evo- Numerical results for non-rotating, bounded and un-
lution of the energy in the simulations discussed in the bounded flows are shown in Fig. 1. In the unbounded
next section. case (runs with an initial energy spectrum ∼ k4 peak-
ing at k = 14 in the DNS and 963 LES, and peaking
at k = 30 in the 1923 LES), the runs show a decay for
D. Enstrophy decay the energy close to ∼ t−10/7 independently of the pres-
enceofhelicityornot(notetherunsalsospanarangeof
Reynolds numbers from Re ≈ 420 to 1200). The decay
From any of the previous energy decay laws, one can
alsocomputelawsfortheenstrophydecayΩ(t)= ω2 /2 is consistentwith the predictiongivenbyEq.(20) foran
initial ∼k4 energy spectrum.
usingtheisotropicenergybalanceequationandreplacing
ǫ=νΩ(t), which results in Ω(t)=ν−1dE/dt. Fro(cid:10)m (cid:11)this The enstrophy decay is also consistent with this law,
as expressed by Eq. (38), decaying close to ∼ t−17/7 in
equation, for every solution for which the energy decays
as E(t)∼tα, the enstrophy decay results all cases. In the absence of rotation, helicity only delays
the onset of the self-similar decay by retarding the time
Ω(t)∼tα−1. (38) when the maximum of enstrophy takes place, as already
reported in [27] and [53]. This is more clearly seen in
Althoughrotatingflowsareanisotropic,the enstrophy the DNS; see, e.g., the time of the peak of enstrophy for
ispredominantlyasmall-scalemagnitudeandwewillsee runs D512-3 and D512H-3 in Fig. 1(b). Finally, in the
that this isotropic argument gives good agreement with helicalruns,helicityseemstodecayastheenstrophy,just
thenumericalresultsforrotatingandnon-rotatingflows. slightly slower than the ∼t−17/7 law.
Since helicity is related with the energy and the enstro- Similarresultsareobservedforboundedflows,i.e., for
phy only througha Schwartzinequality,no simple decay initial conditions with the initial energy containing scale
laws can be derived in its case using these phenomeno- close to the size of the box (runs with a ∼k−4 spectrum
logical arguments. from k = 4 to 14, peaking at k = 4). In this case, all
runs are consistent with a ∼ t−2 decay for the energy
(seethe insets ofFig.1)inagreementwithEq.(19),and
IV. TIME EVOLUTION - NUMERICAL a decay for the enstrophy close to ∼ t−3 in agreement
RESULTS with Eq. (38). In the helical runs, helicity decays again
slightlyslowerthantheenstrophy,butclosetothe∼t−3
power law.
We present here the results for the energy, enstrophy,
and helicity decay obtainedin the numericalsimulations
listed in Table I, classifying them as rotating or non-
B. Rotating flows
rotating, bounded or unbounded (in the sense that the
initial integral scale is smaller than the size of the box),
and helical or non-helical. 1. Global quantities
Concerning the terminology of “bounded” and “un-
bounded” used to describe the numerical simulations, it As rotation is increased, the simulations show a shal-
is important to note that confinement effects in a ro- lower power law in the energy decay. As an illustration,
tating flow go beyond a saturation of the integral scale Fig. 2 shows the energy decay rate in simulations of un-
when it grows to the box size. Confinement also selects bounded non-helical flows with increasing rotation rate
a discrete set of inertial waves which are normal modes Ω. As reported in previous numerical simulations [24]
of the domain, and boundaries can introduce dissipation and experiments [12, 13], as Ω increases the decay slows
throughEkmanlayers. Thelattereffectisnotpresentin down until reaching a saturated decay for Ro≈0.1. We
our numerical simulations with periodic boundary con- will focus in the following in simulations with Rossby
ditions. Finally, it was shown in [51] (see also [52]) that number small enough to observe this saturated decay,
thesmallnumberofFouriermodesavailableintheshells althoughnotsosmallthattherotationquenchesallnon-
with wave number k ≈ 1 gives rise to poor represen- linear interactions giving only exponential decay. A de-
tation of isotropy and of the integral scale in runs for tailed study of the transition between the non-rotating
which the integral scale approaches 1/5 of the box size. and rotating cases can be found in [12].
10
FIG. 7: Axisymmetric energy spectrum e(k ,k )/sinθ for
k ⊥
different times for run L96-6 (non-helical, Ω=10, initial en-
ergy spectrum ∼k4 peaking at k=14).
well as the effect of anisotropy in the initial conditions
which is specially relevant for this particular case. The
enstrophy and helicity show a decay close to ∼ t−12/7.
FIG. 6: (a) Evolution of the isotropic energy spectrum E(k) Note that in the presence of rotation, helicity not only
forL96-6(non-helical,Ω=10, initial∼k4 spectrumpeaking slows down the occurrence of the peak of enstrophy as
atk=14)fromt=5tot=100withtimeincrements∆t=5.
already reported in [27], but it also changes the energy
Inset: reduced perpendicular energy spectrum E(k ) for the
⊥ decayafterthispeak. Theenstrophydecayisnotaffected
same times. (b) Evolution of the isotropic energy spectrum
forL96H-2(helical, Ω=10, initial∼k4 spectrumpeakingat by the presence of helicity.
Overall, the case of constrained runs shows a similar
k = 14) at the same times, with the reduced perpendicular
energy spectrum in theinset. scenario, with a significant slow down of the decay rates
inthe presence ofrotation,and withan extraslow down
of energy decay in the presence of helicity (see the in-
sets in Fig. 3). Rotating non-helical flows are close to
Figure 3 shows the energy, enstrophy, and helicity de-
E(t) ∼ t−1, Ω(t) ∼ t−2, and H(t) ∼ t−2, while heli-
cay in simulations of rotating flows with and without
cal flows in this case display a shallower decay in the
helicity, in the unbounded andbounded cases (the latter
energy consistent with E ∼ t−1/3 as predicted by the
in the insets). The energy decay in the unbounded non-
phenomenological arguments that take into account the
helical runs (thin lines in Fig. 3) is slightly steeper than
what Eq. (25) predicts (E ∼t−5/7). A better agreement effectofhelicityintheenergyspectrumofrotatingturbu-
isobservedfortheenstrophy,whichisclosertoa∼t−12/7 lence. As inthe unboundedcase,the presenceofhelicity
does not affect the decay rate of enstrophy.
law. As will be shown next, the agreement between the
phenomenologicalargumentsforthe energyandthe sim-
ulations is improved if the decay of 2D and 3D modes is
considered separately. 2. Anisotropic global quantities
Alternatively, the unbounded helical runs (thick lines
inFig. 3) showfor the energya ∼t−1/3 decayor steeper Althoughtheimpactofrotationandhelicityintheen-
(although shallower than ∼t−5/7). Runs with ABC ini- ergydecayis clear,thepredictionsgivenbytheisotropic
tial conditions tend to develop a clearer power law de- phenomenologicalargumentsinSec.IIIBdonotcoincide
cay and to be closer to a ∼ t−1/3 decay than runs with in all cases with the results from the simulations. This
random helical modes. Again, these differences can be can be ascribed to the fact that these arguments assume
explained considering the decay of 2D and 3D modes, as an isotropic energy scaling, while rotation breaks down