7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010
A Lagrangian approach for quantifying the segregation of inertial particles in
incompressible turbulent flows
E. Meneguz* and M. W. Reeks*
Mechanical and Systems Engineering, University of Newcastle, Newcastle upon Tyne, NE1 7RU, UK
Elena.Meneguz@ncl.ac.uk and Mike.Reeks@ncl.ac.uk
Keywords: turbulent flows, inertial particles, segregation
Abstract
This paper is concerned with the quantification of preferential concentration of heavy particles in incompressible
turbulent flow. We exploit the so called Full Lagrangian Approach (Osiptsov (2000)) based on the definition of a
unit deformation tensor and we show how this can be used to evaluate the statistical properties of the compressibility
of the particle velocity field and of the particle concentration. The main objective of this work is the comparison
of results in a DNS of statistically stationary homogeneous isotropic turbulence with those of a three dimensional
\ Imlici ik flow composed of 200 random Fourier modes. We observe qualitatively the same trend in both flow fields: in
particular there exists a critical value of the Stokes number (the dimensionless particle relaxation time) below which
the compressibility of the particle velocity field remains negative in the course of time, indicating that the process of
segregation goes on indefinitely (until the particles collide). In addition, we calculate the probability density function
of the log J(t) I and we show that its evolution in time is close to a Gaussian distribution.
Introduction
There have been numerous studies devoted to the seg
regation of particles/droplets in turbulent flows, given
the relevance of this phenomenon in many environmen
tal and industrial activities (warm rain initiation (Shaw
(2003)) and formation and growth of PM10 particu
lates in the atmosphere are only two of many exam
ples). Early experiments and simulations (e.g. Crowe
et al. (1993)) have shown that not only the particles seg
regate reaching a maximum when the particle response
time is approximately equal to the timescale of the tur
bulent structure (i.e. the Stokes number ~ 1), but also
the suspended particles seem to cluster into regions of
high strain rate and expelled from regions of high vor
ticity. Maxey and his coworkers (e.g. Maxey (1987))
showed that the gravitational settling of particles in ho
mogeneous turbulence is enhanced due to preferential
sweeping in the direction of gravity. Since then there
have been many studies to understand and quantify this
segregation process, e.g. the seminal studies by Sun
daram and Collins (1997) and Wang et al. (1998) to
quantify the influence of segregation on twoparticle dis
persion and the process of particle agglomeration.
In more recent years it seems that there has been a ten
dency, even though from different viewpoints, in mod
eling the turbulent flow as the sum of two contributions
that we could attempt here to generalize as follows: an
underlying velocity field that accounts for all particle
particle and fluidparticle two point spatial correlations
and a spatially uncorrelated component that makes par
ticles with inertia to cross their trajectory and yet having
different velocities one to another (F6vrier et al. (2005),
Reeks (2004), Falkovick and Pumir (2004) and Wilkin
son et al. (2007)).
A different mechanism is the one proposed by Chen et
al. (2006) and Coleman and Vassilicos (2009), who have
demonstrated that there exists a strong correlation be
tween the segregation of inertial particles and the pres
ence of zero acceleration points in DNS of turbulent
flows. In addition, it has been hypothesized that in
Kinematic Simulations (KS) where some features like
sweeping of small scales by the large ones are not easily
reproduced heavy particles anticluster with zero veloc
ity points (Chen et al. (2006)). In addition, the motion of
particles in terms has been studied in terms of dynami
cal systems theory as in Bec (2005) and Wilkinson et
al. (2007) who obtained an analytical expression for the
Lyapunov exponents associated with the motion of iner
tial particles in physical space.
In the present work we exploit the Full Lagrangian Ap
proach (FLA) (see Osiptsov (2000), Reeks (2004) and
Healy and Young (2005)) that consists of calculating the
size of an infinitesimally small volume occupied by a
group of particles, along the trajectory of one single par
ticle. This immediately yields the concentration of parti
cles along the trajectory, since the inverse of the volume
occupied by a fixed number of particles corresponds to
the particle concentration by definition. It is worth not
ing that in the previous studies mentioned, the FLA was
mainly used to calculate the particle concentration at
some point and for a given time, while we recently ex
tended the methodology further (IJzermans et al. (2010)
and IJzermans et al. 2's i'i.. Because it was the first
time this has been done, at first we needed to test the va
lidity of our results in rather simple flow fields, as we did
in the the work above cited. We then went on to validate
our results in a more complex three dimensional turbu
lent velocity field where no modelling is applied but the
Navier Stokes equation are solved directly, this being the
main object of the present work.
The paper is structered as follows: in the next section
we describe the methodology used, including details of
the characteristics of the flow fields used in the simu
lations, the governing equations of motion for particles,
the deformation tensor, its evolution equations and re
lated statistics. The third section is devoted to the dis
cussion of findings and open questions currently under
investigation; we summarize our work and prospects for
future studies in the fourth section.
Method
Details of the simulations are described in the following
subsection. We then discuss the governing equation
of motion for the inertial particles and we illustrate in
details the Full Lagrangian Method Osiptsov (2000).
Fluid velocity field
The threedimensional incompressible Navier Stokes
equations are solved using a pseudospectral method
performed on a triplyperiodic cubic domain of length
L = 27r. To remove the aliasing error so introduced, the
so called "3/2" rule is implemented, according to Peyret
(2002). In order to keep the statistically stationary
turbulence, forcing is applied in a similar way to Goto
(2008), that is by fixing the magnitude of the Fourier
components of the vorticity corresponding to the lowest
wavenumbers. Based on a 1283 grid, we achieve a
Re = 65. Other relevant parameters of the simulation
are contained in Table 1. Finally, for what concerns the
discretization in time we used a fourth order Lawson
scheme which treats the linear terms of the differential
equation with exact integration (integrating factor, see
Canuto et al. (1991)) and the non linear ones with an
explicit method (fourth order RungeKutta scheme).
For the details of the kinematic simulations, we refer to
7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010
Table 1: Numerical and flow parameters of DNS.
r (Kolmogorov length scale) 0.032840 m
L (integral length scale) 1.2 m
t, (eddy turnover time) 1.04 s
Ax (grid space) 0.049 m
v cinematicc viscosity) 0.00923 rm2/s
T, (Kolmogorov time scale) 0.116847 s 1
kmax,, 2.1
Uzermans et al. (2010). Herewith we only state that the
incompressible velocity field is obtained by imposing
a Kraichnan energy spectrum as in Kraichnan (2000)
and Spelt and Biesheuvel (1997) which provides a good
description of turbulence at low Reynolds number. As
well as in our DNS simulations, also in this case we
impose periodic boundary conditions by using a three
dimensional periodic box.
Lagrangian description of the particle phase
We consider a pointparticle approach with the dispersed
phase described as a continuum. Under the assumption
that the density of the particle is much higher than the
density of the carrier flow, that is pp/p > 1, and neglect
ing Brownian motion and any body force, the equation
of motion of heavy particles is reduced to:
dx dv 1
v, (u v), (1)
dt dt Tp
where Tp is the particle relaxation time defined as Tp 
2. ,', /9vp for a spherical particle of radius a,.
The Full Lagrangian Approach (see Osiptsov (2000),
Reeks (2004) and Healy and Young (2005)) is based
upon a second order deformation tensor Jj
d9xz(xo, t)/dxoj where by definition Jij: J det(Jij).
xo represents the position of the centre of the infinites
imally small volume surrounding each particle at some
initial time t 0 Since we are interested in the time
dependent behaviour of such parameter, we then calcu
late the first and second derivatives for each component
of Jij. Rearranging, we obtain the following evolution
equations:
dt Jdti Tp Ji
(2)
We choose as initial conditions Jij(0) = ij and
ij (0) dui (xo, 0)/Oxj. That being said, the choice of
the initial conditions does not affect the long term statis
tics. We note that by definition it holds:
IJW (t) 1(t), (3)
where c(t) is the particle number concentration normal
ized so that at time t 0 it is equal to J 1 1 From
7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010
Eq. 3 it is trivial to deduce that J'a = c a for any
a E R. We now define two different type of averages,
a particle ensemble average denoted by the brackets (K),
and a spatial average defined over a certain domain Q2,
denoted by the overbar It holds the following relation:
( 1) A ((x)c(x)dx F, (4)
By recursively substituting 4 = c and recalling Eq. 3,
one can calculate any spatially averaged moment of the
particle number concentration provided that the number
of particle trajectories is sufficiently large inside a par
ticular volume. The latter condition is easily verified in
the cases here examinated where we are dealing with a
periodic domain. This yields to the following equation:
C (Il1 "), V aER. (5)
The Full Lagrangian Approach can also be used to deter
mine the compressibility of the particle phase. Starting
again from Eq. 3, we can write:
d d
d In J d In n, (6)
dt dt
Because dn/dt
nV v, Eq. 6 can be rewritten as:
dt
hInl J V v. (7)
Results and discussion
The result of our simulations carried out in the kinematic
model are discussed in details in IJzermans et al. (2010).
The primary objective of this work is the validation of
such outcomes in a more complex flow field in which the
dynamical NavierStokes equations are solved explic
itly; this means that some characteristic features of real
flows, such as the sweeping of the small scales by the
large ones, are in principle taken into account. For this
reason we focus on the three dimensional DNS of ho
mogeneous isotropic turbulence described in the method
section. We let the flow run for a few eddy turnover
times until it reaches a steady state; the corresponding
fluid velocity on the grid nodes is then saved and used
as initial conditions for when we inject the particles. At
time t = 0, we release randomly a group of 100.000
heavy particles with a uniform distribution inside a pe
riodic cubic domain of length L = 27r. Trajectories are
then computed in time using a fourth order RungeKutta
numerical scheme; the velocity of the fluid at the particle
position is obtained with a 6th order polynomial interpo
lation routine.
To start with, we calculate the spatial average of the mo
ments of the particle number concentration ca as a func
tion of time. For the sake of brevity, we illustrate only
ca=3
a=2
 =0
10 15 20 25
Figure 1: Spatial average of the moments of the particle
number concentration in the 3D DNS, St=0.4.
the case of St = 0.4 (Fig. 1), where this latter is defined
as the ratio between the particle relaxation time and the
Kolmogorov time scale of the flow, i.e. St = Tp/T. We
compare it with for the result obtained in our KS model
in Fig. 2 where the Stokes number is defined as the ra
tio between the particle response time and the typical
timescale of the flow. As we previously observed (IJz
ermans et al. (2009)), the curves exhibit an exponential
trend, and the order of magnitude is proportional to the
order of the moment. What is really important here, is
to observe the intermittency of the lines which might be
explained by the occurance of "singularities", i.e. events
for which JI = 0 that happen when the particle number
concentration becomes large enough to be considered in
finite. We have previously attributed this phenomenon to
the so called Random Uncorrelated Motion Uzermans et
al. (2010)), although the contribution that this can make
is very limited for St < 1. We are currently investigat
ing the distribution of singularities and its evolution in
space and time.
We then calculate the compressibility of the particle
phase averaged over all the particles versus time as Fig.
3 displays. Among the St numbers considered, the peak
of the segregation seems to happen for St 0.7, this
corresponding to the minimum of d < Inl J > /dt. We
observe that the curves are highly intermittent.
In Fig. 4 we also plot the long time limit of
the compressibility, that is lirmtot1 lnl J(t)l
limnt, d(ln J(t) )/dt limrnt, Vv as a function of
three different particle response times which correspond
to St 0.7, St 1 and St 4. For comparison,
in Fig. 5 we plot the same quantity obtained in our KS
7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010
1020
1015
0 1010
105
S0
. 2 8..... ...... .. ... ........  2
0 2 4 6 8 10 12 14 16 0
t
Figure 2: Spatial average of the moments of the particle
number concentration in the 3D kinematic simulation,
St=0.4.
simulations. Clearly there is a good qualitatively agree
ment between the two. Quantitatively, the magnitude of
the segregation seems to be higher in the case of DNS,
and this enhancement is probably due to a broader even
though still limited distribution of scales. On the other
hand, how do we justify the positive values of the com
pressibility? We do not have an answer to such question
yet. The calculation of an average quantity only pro
vides us with a picture of the net effect between com
pression/expansion and what effect prevails in the course
of time, but it is likely that having a positive compress
ibility on average does not mean that the segregation is
zero everywhere. The same concept applies to negative
compression. Given the nature of Jij(t), every time a
particle experience a positive compressibilty the elemen
tal volume that sorrounds that particle is dilated. How
ever, it is questionable whether Jij (t) is able to measure
a sort of "local diffusion", where by this term we intend
the displacement with respect to a fixed point. It would
be somehow intuitive to argue that the vorticity plays an
important role. Nevertheless, we recently constructed a
simple flow field model of homogeneous isotropic tur
bulence made of pairs of counter rotating vortices (IJz
ermans et al. (2010)) in which the linear velocity field
contains no vorticity, and yet we were able to observe
positive compressibility for approximately St > 1.5. In
addition, Chen et al. (2006) have found that in their two
dimensional flow field simulated kinematically, inertial
particles cluster regardless of the amount of vorticity and
strain in the flow, the ratio among those two quantities
remaining 1 contrary to the traditional view of particles
expelled from regions of high vorticity and segregating
SSt=O 7
St=1
3 St=O 1
# St=10
+*u *
.. A
4~~
A
10 20 30 40 50 60
t/eta
eta
Figure 3: Particle averaged compressibility as a function
of time for four different Stokes number in the 3D DNS.
into regions of strain. We are currently in the process of
investigating what causes the compressibility to become
positive.
Finally, we calculated the probability density function
(PDF) for the occurrence of loglJ(t)\ (see Fig. 6 for
St 1). We have found out that the more time in
creases the more it approaches a Gaussian with negative
mean as illustrated by the standardized central moments
of the distribution for the three different timescale con
sidered. We can then confirm what was argued by Reeks
(2004), i.e. that a Gaussian convection diffusion process
can accurately describe the dispersion of log J(t) .
Conclusions
In this paper we have studied the dispersion of inertial
particles in a three dimensional DNS of homogeneous
isotropic turbulent flow, and we have benchmarked our
results with a three dimensional model of turbulent flows
kinematically simulated. With the use of the Full La
grangian Approach, we have been able to show that there
is not a qualitative difference in the statistics of the par
ticle concentration and of the compressibility.
We felt the need to use DNS calculations to introduce a
few effects not included in our kinematic simulations: a
broader range of scales (we now have (L/rI 126 in com
parison with the "nearly single scale" KS model) and
other additional features such as vortex stretching and
interaction between the modes. As we expected, a wider
distribution of scales has not significantly altered the re
sults for J(t) essentially because this is a local param
eter whose value only depends on the gradients of the
I U
7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010
St=1
St=O 7
E
S05
1 5
2
0 10 20 30 40 50 60 70 80
t/t
Figure 4: Longtime approximation for the particle av
eraged compressibility as a function of time for three
different Stokes number in the 3D DNS.
velocity field. In the context of our study this is not re
garded as a limitation as many interparticle effects such
as collisions and twoparticle separation are local effects
as well. Therefore, this authorizes us to carry out further
investigations in kinematic models that have the advan
tage of not being very heavy computationally (in con
trast to other techniques such as DNS).
We are currently in the process of gaining a better un
derstanding of the process of segregation especially re
garding the nature of two effects: positive values for the
compressibility above a certain St number and the fre
quency of singularities and their distribution in time and
space.
Acknowledgements
The present work has been carried out under EPSRC
GB's financial support. The authors would like to ac
knowledge Stuart Coleman and Christos Vassilicos, De
partment of Aeronautics, Imperial College London for
providing the basic version of the DNS code described
in the paper.
References
Osiptsov A.N., Lagrangian modelling of dust admixture
in gas flows, Astrophysics and Space Science, Vol. 274,
pp. 377386, 2000
Shaw R.A., Particleturbulence interactions in atmo
spheric clouds, Ann. Rev. Fluid Mech., Vol. 35, pp. 183
227, 2003
05
04 j
03
o 02
A 01
S,
5 10 15
time*
Figure 5: Longtime approximation for the particle av
eraged compressibility as a function of time for three
different Stokes number in the 3D KS.
Crowe C.T., Chung J.N. and Troutt T.R., Chapter 18. In
Particulate TwoPhase Flow (ed. M.C.Roco), Vol. 626,
pp. 11, Heinemann, Oxford, 1993
Maxey M.R., The gravitational settling of aerosol parti
cles in homogeneous turbulenceand random flow fields,
J. Fluid Mech., Vol. 174, pp. 441465, 1987
Sundaram S. and Collins L.R.,Collision statistics in an
isotropic particleladen turbulent suspension. part 1 di
rect numerical simulations, J. Fluid Mech, Vol. 335,
pp.75109
Wang L.P., Wexler S. and Zhou Y., Statistical mechani
cal descriptions of turbulent coagulation, Physics of Flu
ids, Vol. 10(10), pp. 26472651
F6vrier P, Simonin O. and Squires K.D., Partitioning of
particle velocities in gassolid turbulent flows into a con
tinuous field and a spatially uncorrelated random dis
tribution; theoretical formalism and numerical study, J.
Fluid Mech., Vol. 553, pp. 146, 2005
Reeks M.W., Simulation of particle diffusion, segrega
tion, and intermittency in turbulent flows, in Proc. of IU
TAM Symposium on Computational Modelling of Dis
perse Multiphase Flow, (ed. S. Balachandar A. Pros
peretti), pp. 2130, Elsevier, 2004
Wilkinson M., Mehelig B., Ostlund S. and Duncan K.P.,
Unmixing in random flows, Physics of Fluids, Vol.19,
113303,2007
. time=25, mean=2 97, skewst= 15, kurtt=3 10
0 time=10, mean=1 19, skewst=0 28, kurtst3 40
time=l, mean=0 18, skewst2 22, kurtt=1327
Lo
C 04r
0
30 25 20 15 10
IoglJI
7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010
Denissenko P., Falkovich G. and Lukaschuk S., How
waves affect teh distribution of particles that float on a
liquid surface, Phys. Rev. Letters, Vol. 97, pp. 244501,
2006
Falkovich G. and Pumir A., Sling effect in collisions of
water droplets in turbulent clouds, American Metereo
logical Society, pp. 44974505, 2007
Goto S., A physical mechanism of the energy cascade in
homogeneous isotropic turbulence, J. Fluid Mech., Vol.
605, pp. 355366, 2008
Canuto C., Hussaini M.Y., Quarteroni A. and Zang T.,
Spectral methods in Fluid Dynamics, Springer, 1991

5 0 5 10 Peyret R., Spectral methods for incompressible visocus
flow, Springer, 2002
Figure 6: PDF of occurrence of loglJI for St=l.
Chen L., Goto S. Vassilicos J.C., Turbulent clustering of
stagnation points and inertial particles, J. Fluid Mech.,
Vol. 553, pp. 143/154, 2006
Bec J., Multifractal concentrations of inertial particles in
smooth random flow, J. Fluid Mech., Vol. 528, pp. 255
277, 2005
Healy D.P and Young J.B., Full lagrangian methods for
calculating particle concentration fields in dilute gas
particle flows, Proc. Roy. Soc. London A: Mathematical,
Physical and Engineering Sciences, Vol. 461(2059), pp.
21972225, 2005
Kraichnan R.H., Diffusion by a random velocity field,
Physics of Fluids, Vol. 13, pp. 2231, 2000
Spelt PD. amd Biesheuvel A., On the motion of gas
bubbles in homogeneous isotropic turbulence, J. Fluid
Mech, Vol. 336, pp. 221224, 1997
Coleman S.W. and Vassilicos J.C., A unified sweepstick
mechanism to explain particle clustering in two and
threedimensional homogeneous, isotropic turbulence,
Physics of Fluids, Vol. 21, pp. 113301113310, 2009
IJzermans R.H.A., Meneguz E. and Reeks M.W., Segre
gation of particles in incompressible random flows: sin
gularities, intermittency and random uncorrelated mo
tion (RUM), J. Fluid Mech., in press
IJzermans R.H.A., Reeks M.W., Meneguz E., Picciotto
M., Soldati A., Measuring segregation of inertial parti
cles in turbulence by a full Lagrangian approach, Physi
cal Review E, Vol. 80, pp. 015302, 2009
nJ
