Disentanglement of inhomogeneous and anisotropic
Rayleigh–Taylor turbulenceTakeshi Matsumoto1

1

Division of Physics and Astronomy, Graduate School of Science, Kyoto University, Japan

Introduction

law or a diﬀerent intermittency. We address this question
with the disentanglement method.

The Kolmogorov’s phenomenology of turbulence, which
is a dimensional analysis with identification of the key
quantity such as the energy dissipation rate, provides us
the intuitive and quantitative understanding of the homogeneous and isotropic turbulence. However, it is very
hard to extend this phenomenology to inhomogeneous andanisotropic turbulent systems. To confront this diﬃculty,
there developed a way to “disentangle” [1] complicatedknots of various inhomogeneous and anisotropic eﬀects
one-by-one. This can be achieved by expanding the statistical quantities like the longitudinal velocity structure
functions with the spherical harmonics Yℓ,m (θ, φ) as⟨{S(p)(x; r)r[ u(x + r, t) − u(x, t) ] ·|r|=∑ ∑∞=504030z

120L(t)10}p ⟩00.5 1 1.5ℓ(p)Sℓ,m (x, r) Yℓ,m (θ, φ),

ρ(z)

(1)

ℓ=0 m=−ℓwhich is also known as the SO(3) decomposition. We
here assume that the system is in a statistically steady(p)state. The coeﬃcient Sℓ,m (x, r) carries information of eachanisotropic sector (ℓ, m) on a specified spatial point x. It
is believed that the coeﬃcient has the scaling behavior(p)ζFig. 1. Left: the horizontally averaged density profiles in
the accelerated frame. The dashed line corresponds to theinitial profile. The domain size is (4π)2 × 16π and theresolution is 1282 × 512. Right: color-coded density fieldshowing typical plumes going up and down (Red: heavy
fluid, Blue: light fluid).

(p)Sℓ,m (x, r) = C(x) × r ℓ,m(2)

if the energy containing scale and the dissipation scale areorder-of-magnitude diﬀerent. This implies that the inhomogeneity does not aﬀect the scaling exponents. Indeed,
for inhomogeneous and anisotropic systems, the isotropic
sector (ℓ = m = 0) is known to follow the Kolmogorov(2)
scaling, e.g., ζ0,0 = 2/3 [1]. With this disentangle method,
one can address many questions such as how does the scal(p)ing exponent ζℓ,m depend on the sector (ℓ, m) ?; is there
universality of these exponent values? and so on [1].
In this paper, as an inhomogeneous and anisotropic turbulent system, we numerically study so-called Rayleigh–Taylor (RT) turbulence, which is a turbulent state initiated by the Rayleigh–Taylor instability between the two
miscible fluids with diﬀerent densities under the eﬀect of
the gravity. If the density diﬀerence of the top and bottom fluids is very small, recent simulation studies of RT
turbulence [2, 3] motivated by its Kolmogorov-type phenomenology [4] have revealed that this particular kind of
RT turbulence is very close to homogeneous and isotropic
turbulence in terms of the scaling laws of the velocity anddensity structure functions (the density behaves the same
way as the passive scalar). What if the density diﬀerence
is not small? We expect that there is a diﬀerent scaling

2Equations of the Rayleigh–Taylor turbulencewith a moderate density diﬀerence

The density diﬀerence of the RT turbulence canbe measured with the Atwood number A = (ρtop −ρbottom )/(ρtop + ρbottom ). In the studies mentioned in theprevious section, the typical value of A was 0.10, wherethe Boussinesq approximation (valid for A ≪ 1) was employed. Here we take a moderate value A = 0.50 and use
one form of the anelastic approximations proposed in [5],thereby neglecting a compressible eﬀect. For this RT turbulence with A = 0.50, there is numerical evidence that its
mixing property is peculiar and diﬀerent from that of homogeneous and isotropic barotropic turbulence [6], which
is why we study this system. The equations of the RT
turbulence with the anelastic approximation are∂t ρ + ∇(ρu)

=0,

∇·u

=

−∇

ρ[∂t u + (u · ∇)u]=−∇(p + ρ0 gz) + Aρ̃g + η∇2 u()η++ ζ ∇(∇ · u).
(5)
3()

D
∇ρ ,ρ

(3)
(4)

〈(δxux)p〉, 〈(δyuy)p〉, 〈(δzuz)p〉

1

r2/3
0.1

r4/3
0.01

p=2
p=40.1110

r10(b)

r2/3

1

r4/3
0.1
p=2, (l,m)=(0,0)
p=4, (l,m)=(0,0)

0.01
0.1

110r
103

r6/5(a)

0.001

S(p)
l,m(x,r)Here ρ, u is the density and the velocity; D is the diﬀusion
coeﬃcient associated with the miscibility of the two fluids;
ρ0 is the reference density ρ0 = (ρtop + ρbottom )/2; p isthe pressure; g = (0, 0, −g) is the gravitational acceleration; η and ζ are viscosity coeﬃcients. They are solved in a
triply periodic domain with the standard spectral method.
Mainly the result of the domain size (4π)2 × 32π with resolution 2562 × 2048 is reported here. The other parametervalues are η/ρ0 = D = 10−3 , ζ = 0.0. We use the 2ndorder Adams–Bashforth method for the temporal integration.
The initial density is given by a step function of thevertical coordinate z with small two-dimensional perturbations, whose profile is shown as the dashed line in Fig.1.
Notice that the RT turbulence is unsteady. The simulations are stopped when the turbulent-zone length (denotedas L(t) in Fig.1) reaches 60 % of the vertical domain sizein order to suppress eﬀect of the boundary condition. Repeating this process with several diﬀerent initial density
perturbations, we obtain an ensemble of the RT turbulenceover which the average is calculated. In the calculation of
the structure functions, the position x and increment r
are both taken to be inside the turbulent zone, which iswithin L(t) shown in Fig.1.
Disentanglement of the structure functions

(c)

r2

|S(4)
l,m(x,r)|1p
In practice, the coeﬃcients Sℓ,m(x, t) in the SO(3) decomposition Eq.1, are numerically obtained as follows.
Firstly, we calculate the structure function on the Cartesian grid points of r. Secondly, we interpolate the structure function data onto the Gauss-Legendre points for asuitably chosen set of r = |r|’s. Lastly, we calculate the coeﬃcients from the data on the Gauss-Legendre points withspherepack, a software library of the spherical harmonics.
We begin with Fig.2 (a) showing that the structure
functions of the longitudinal velocity components without using the disentanglement method. Here, in addition to the ensemble average, the spatial average with respect to the position x inside the turbulent zone is taken.
No scaling behavior is seen, which is mainly due to thelow Reynolds number. However the isotropic components(p)S0,0 (x, r) shown in Fig.2 (b) exhibit better scaling qualities, whose exponents are consistent to the Kolmogorov
scaling law. This demonstrates power of the disentanglement method. Here x is taken to be the middle of the(4)turbulent zone. The anisotropic sector S2,m of the fourthorder structure function is plotted in Fig.2 (c). Again thescaling r2 is observed, which coincides with the scaling of
anisotropic barotropic turbulence. The same plots of otheranisotropic sectors indicate that the scaling law of the longitudinal velocity of the RT turbulence is identical to thatof the barotropic turbulence on the contrary to our expectation that a diﬀerent scaling law might exist for the RT
turbulence with not-small Atwood number. The inhomogeneity eﬀect is studied by changing the position x in theturbulent zone, which shows that the scaling exponents are
the same for diﬀerent x as expected.
For the structure functions of the density, it appears
that they obey the same isotropic and anisotropic scaling
law of the passive scalar advected by the barotropic turbulence, which is again contrary to our expectation. However
for the density at large scales there is an indication of a different scaling law from the Obukhov-Corrsin scaling law,which will be discussed at the conference.

0.1
0.01
(l,m)=(2,0)(l,m)=(2,1)(l,m)=(2,2)

0.001
0.0001
0.1

110r

Fig. 2. (a) Second and fourth-order longitudinal velocity structure functions of the RT turbulence at late time
without the disentanglement. The increments r are alongwith the x, y and z axes (the longer curves correspondto the z-axis case). The power laws r2/3 and r4/3 correspond to the Kolmogorov law and r6/5 corresponds to theBolgiano–Obukhov law of the stratified turbulence, which
is another scaling candidate for the RT turbulence. (b)Their isotropic sectors extracted with the disentanglement
method. (c) The anisotropic sectors ℓ = 2 for the 4th-orderstructure function.
Acknowledgement
A part of the simulation was performed with SX8 in
the Yukawa Institute of Kyoto Univerisity. This work wassupported by the Grant-in-Aid for Young Scientists (B)No. 2174290 from the JSPS.
References[1] Biferale and Procaccia, Phys. Rep., 414, 43–164, 2005.
[2] Boﬀetta et al., Phys. Rev. E, 79, 065301(R), 2009.
[3] Matsumoto, Phys. Rev. E, 79, 055301(R), 2009.
[4] Chertkov, Phys. Rev. Lett., 91, 115001, 2003.
[5] Joseph, Eur. J. .Mech. B, 9, 565–596, 1990.
[6] Livescu and Ristorcelli, J. Fluid Mech., 591, 43–71, 2007.