Paper No 7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30-June 4, 2010
Effects of drift velocity on the interaction timescales between heavy particles and SGS
eddies in the Langevin equation for LES
Guodong Jin', Guo-Wei He' and Jian Zhang'
LNM, Institute of Mechanics, Chinese Academy of Sciences, Beijing 100190, China
gdjin@lnm.imech.ac.cn, hgw@ lnm.imech.ac.cn, /lhii-n!Jiiji7lll iniiKchll Ci
Keywords: Keywords: large-eddy simulation, Langevin equation, particle-eddy interaction, crossing trajectory effect
Abstract
The interaction timescale between the heavy particle and the subgrid scale (SGS) eddy is an important parameter for the
closure of the Langevin equation used to recover the SGS fluid velocity seen by heavy particles in Large-Eddy Simulation
(LES). In addition to the filter width and particle inertia, the interaction timescale is also strongly affected by the averaged
relative velocity or the drift velocity between the heavy particle and SGS eddies caused by the gravitational force or other body
forces. In this paper, we study the effects of the drift velocity on the particle- SGS eddy interaction timescales using an
isotropic turbulent field (Rek=102) obtained by Direct Numerical Simulation (DNS). We derive the analytical formulae of the
timescales of the velocity components parallel and normal to the drift velocity in the SGS turbulent field and validate them by
directly calculating the timescales in a DNS field. It is observed that the proposed formulae reasonably fit the numerical results,
and that with the increase of the drift velocity, the timescale of the velocity component parallel to the drift velocity 6TLp1
decreases due to the crossing trajectory effect, the continuity effect further decreases the timescale of the velocity component
normal to the drift velocity 6TLp2. In the limit of very high drift velocity, 6TLp2=0.5 6TLpl.
Introduction
Large-eddy simulation (LES) is widely used to simulate the
particle-laden turbulent flows with the Eulerian-Lagrangian
approach. In LES, only large-scale eddies of the carrier fluid
flow are explicitly resolved and the effects of unresolved,
small or subgrid scale (SGS) eddies on the large-scale ones
are modelled using a SGS model. The SGS turbulent flow
field itself is not available. The trajectories of particles are
tracked by integrating the equations of individual particles.
The effects of SGS turbulent velocity field on particle
motion have been studied by Wang and Squires (1996),
Armenio et al. (1999), Yamamoto et al. (2001), Shotorban
and Mashayek (2006), Fede and Simonin (2006), Berrouk et
al. (2007), Bini and Jones (2008), and Pozorski and Apte
(2009), amongst others. In order to include the effects of
SGS eddies on inertial particle motions, the stochastic
differential equation (SDE), that is, the Langevin stochastic
equation is introduced to model the SGS fluid velocity seen
by inertial particles (Fede et al., 2006; Shotorban and
Mashayek, 2006; Berrouk et al., 2007; Bini and Jones, 2008;
Pozorski and Apte, 2009). However, the accuracy of such a
Langevin equation model depends primarily on the
prescription of the SGS fluid velocity autocorrelation time
seen by an inertial particle or the inertial particle-SGS eddy
interaction timescale (denoted by 6TLp). From the theoretical
point of view, 6TLp differs significantly from the Lagrangian
fluid velocity correlation time (Reeks, 1977; Wang and
Stock, 1993), and this carries the essential nonlinearity in
the statistical modeling of particle motion. In previous
studies, 6TLp is modeled either by the fluid SGS Lagrangian
timescale (Fede et al., 2006; Shotorban and Mashayek,
2006; Pozorski and Apte, 2009; Bini and Jones, 2008) or by
a simple extension of the timescale obtained from the full
flow field (Berrouk et al., 2007).
We have studied the subtle and non-monotonic dependence
of 6TLp on the filter width and particle Stokes number using
a flow field obtained from Direct Numerical Simulation
(DNS) and proposed an empirical closure model for
6TLp(Jin et al., 2010). However, the gravitational settling
effect was neglected. Under the effect of gravity, a heavy
particle moves relatively to the surrounding fluid with a
drift velocity, crossing the trajectories of fluid elements and
interacting with different small-scale eddies. As a result, we
expect that the timescale of the SGS velocity seen by the
heavy particle decreases as the settling velocity increases.
This is known as the crossing trajectory effect for a heavy
particle in turbulent flows. In the context of LES, the
interaction is restricted to SGS turbulent eddies.
Furthermore, a related effect is that the timescale in the
vertical direction is larger than that in the horizontal
direction due to the continuity effect in an incompressible
turbulent flow. These effects related to the gravity have been
studied by Yudine (1959), Csanady (1963), Reeks (1977),
Wells and Stock (1983), and Wang and Stock (1993),
amongst others, for the full turbulent flow field. The
gravitational effect is important when the particle settling
velocity is comparable to the rms velocity fluctuation of the
SGS flow field. Using the mean fluid-particle relative
Paper No
velocity as the characteristic velocity, and the longitudinal
and transverse subgrid Eulerian length scales as the
characteristic length scales respectively, Fede, Simonin and
Villedieu (2007) derive a model to describe the crossing
trajectory effect on the timescale of subgrid turbulence
along the trajectory of inertial particles.
Following Wang & Stock's idea (1993) in full-scale
turbulent field, we derive the analytical formulae of the
timescales of the velocity components parallel and normal
to the drift velocity in SGS turbulent field and validate them
by directly calculating the timescales in a DNS field.
Nomenclature
kof Cutoff wavenumber
u Fluid velocity in wavenumber space (ms-1)
u Fluid velocity in physical space (ms-1)
u Filtered fluid velocity in physical space (ms-1)
u SGS fluid velocity in physical space (ms-1)
u, I Modelled velocity using the Langevin equation
St Stokes number, St= r / ST
wo Gravitational velocity
T Particle relaxation time
5Lf The longitudinal lengthscale
6f(r) The normalized SGS longitudinal fluid spatial
velocity correlation
g) The SGS transverse fluid velocity spatial
correlation function
5T1, Particle-SGS eddy interaction timescale parallel
to the gravitational direction
STLp2 Particle-SGS eddy interaction timescale vertical
to the gravitational direction
ST, Particle-SGS eddy interaction timescale under
zero gravity
Governing Equations
The one-way coupling assumption is applied, under which
the particle motion is controlled by turbulence and the
effects of particle on fluid motion is negligible.
The forced isotropic turbulent flow field in a periodic box of
side 227 is simulated using the Direct Numerical Simulation
(DNS) method. The Navier-Stokes equations with a
large-scale random forcing term is solved using a
pseudo-spectral method in the Fourier space,
+vk2 >(k,t)= P(k)F(uxo)+f(k,), (1)
where f(k,t) is the Fourier coefficients or the fluid
velocity in wavenumber space, u and o are full-scale
fluid velocity and vorticity in physical space, F denotes a
Fourier transformation, the projection tensor Pm = m -
kkm/k2 (j, m=1, 2,3) projects F(uxo) onto the plane
normal to k and eliminates the pressure term in
Navier-Stokes equations. The flow is driven and maintained
7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30-June 4, 2010
by a random artificial force f(k,t) which is non-zero only
at low wavenumbers kl <, in Fourier space. The
details for the DNS can be found in Wang & Maxey(1993).
In the filtered-DNS (FDNS), the filtered velocity U is
calculated by truncating the high-wavenumber Fourier
modes
k*f
u(x,t)= i(k,t)el'k- (2)
Ikl=ko
where kc is the cutoff wavenumber in FDNS, k0 = 1 is
the lowest wavenumber. The subgrid scale velocity is then
u (x, t) = u(x,t) (x,t)
The DNS of turbulent flow field with the grid resolution
2563 in a cube with each side 2n is performed. The main
parameters of the flow field are listed in Table 1.
Table 1 Flow parameters in the DNS 2563
Reynolds number Rek 102.05
r.m.s. fluid velocity u 19.34
Dissipation rate a 3554.4
Minimum length scale q 0.0135
Minimum timescale z 0.0037
Minimum velocity scale vK 3.6272
Eulerian integral timescale TE 0.050
Lagrangian integral timescale TL 0.037
Viscosity v 0.0488
Using different cutoff wavenumbers, k we can obtain
different SGS flow fields. Figure 1 shows the energy
spectrum of the simulated flows in DNS 2563 and the
filtered DNS at different cutoff locations. Figure 2 shows
the decomposition of a vorticity field with the cutoff wave
number kr =21.
102 I
102
101
S10 0
10-
10-2
10-2
Figure 1: The energy spectrum of the simulated flows in
2563 DNS. The vertical dashed lines represent different
cutoff locations in the filtered-DNS.
Paper No
500
400
300
200
100
50
a
0 1 2 3 4 5 6
(a) full scalevorticity contour
(a) full scale vorticity contour K +
P500
400
300
200
100
50
Io 4
o
U 1 2 3 4 5 6
xvorticity contour
(b) filtered vorticity contour 5(7 +
SGS
500
400
300
200
100
7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30-June 4, 2010
The dispersed phase is composed of Np spherical heavy
particles with a diameter smaller than the Kolmogorov
length scale, dp = 0.5q. In general, there are several different
forces acting on a heavy particle suspended in a nonuniform
and unsteady turbulent flow field. When the ratio of particle
density to fluid density is much larger than one, the equation
of motion for a heavy particle can be approximated as
dx (t)
= v,(t, (4)
dt
(5)
where Xp(t) and Vp(t) are particle position and velocity at
time t, wo is the particle Stokes terminal velocity under
gravity in otherwise quiescent fluid, wo = gzp, g is
gravitational acceleration, Tp is the particle Stokes response
time, representing the timescale for particle to respond to
the variation of surrounding fluid, u(Xp(t), t) is the fluid
velocity seen by the particle which is obtained from DNS
flow field by a three-dimensional six-point Lagrangian
interpolation scheme, f(Rep) is the nonlinear drag correction
coefficient,
f(Re,)=1+ 0.15Re687, (6)
which is determined by the instantaneous value of the
particle Reynolds number
Re =lu-vld/v. (7)
Similar to the previous studies (Fede et al., 2006), the full
fluid velocity seen by an inertial particle is modeled using
an extended, stochastic Langevin equation as
du:{ .. d] x iT' - U )dt+ _k_, T )
3 'ST )
where the superscript + in Eq. (8) denotes the modeled full
scale velocity, of which one part is from the resolved
velocity in LES, here, from the filtered DNS flow field in
this study, and the other part is from the Langevin equation.
The modeled full scale velocity is different from the real
solution of the Navier-Stokes equations in DNS. 4 is a
Gaussian random variable of zero mean and unit variance,
ksGS, p is the fluid SGS kinetic energy seen by an inertial
particle, which may not be equal to the fluid SGS kinetic
energy averaged over the whole space, kSGs. 6TLp,, is the
timescale between heavy particles and SGS eddies along the
direction of velocity component u, 6TLp,, is obtained by
integrating the correlation, 6RLp,,, of the SGS fluid velocity
component, u,, seen by an inertial particle. Here, 6RLp,,, is
calculated as
ui (Xp (to), to)u (Xp((to + r), to + r)
u'2 ((to,t)
6
40
0 1 2 3 4 5 6
I^TM.:
-- -- - - - - -
x
(c) SGS vorticity contour a? +627
Figure 2: Vorticity contours of the full scale flow field,
filtered flow field and SGS flow field at z = 7i in DNS
2563, where kc = 21 and ox = 5 +) my = y +)y2.
It is noted that i denotes the velocity component instead of
summation.
Gravitational effects on the interaction timescale
between inertial particle and SGS velocity field
In general situation, both the inertia and the drift velocity
are not zero for inertial particles suspended in turbulent flow
fields, the two effects can be included in 6RLP,,(r) by
extending the hypothesis of Csanady(1963), which can be
used to transfer 6RLp,,(r) from SGS fluid velocity temporal
correlation to SGS fluid velocity spatial correction when the
58R,P (Zr)
dvp(t (_ (u(xp (t), t) Vp (t))f + wo
7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30-June 4, 2010
drift velocity increases from 0 to oo. The hypothesis of
Csanady (1963) states that 6RLp,1(T) along the drift velocity
is a constant on the ellipses
2 2
+ = const (10)
(T,)2 (6L,)2
where (w) is the mean settling velocity of inertial particles
in turbulent flow field, 8T,, is the particle-SGS eddy
interaction timescale under the condition of zero drift
velocity, 8Ls is the longitudinal lengthscale of the SGS
velocity field. The dependence of 6TLp on the filter width
and particle Stokes number can be expressed as (Jin et al,
2003)
T,
- (0.444-0.77kcs)exp -L In ]2 +I
-0.2
0
i qko-0.284
f
------- rko=0.216
-.-.-.-.-.-.. rlkf-0.135
------------ exp(-r/6Lf)
. .
, ~ ~~ I
4 8
r/6Lf
12 16
-(1- f)exp ,
exp 5.15)
where f = 8T, / 8TE, 6T, and 8TE are the Lagrangian and
Eulerian timescales of the SGS flow field. kcf can be
approximately expressed as kcs = /A and A is the
lengthscale of the filter width in physical space. Eq. (11)
captures all the main characteristics of the dependence of
6TLp on St and kef. Eq. (10) makes a smooth transition of
6RLp,I(T) from SGS fluid velocity temporal correlation to
SGS fluid velocity spatial correlation when (w) changes
from 0 to oo. The longitudinal lengthscale 8Lf is defined
as
L f- 8f (r)dr, (12)
, f( \
< ux (xo,Yo,zo)ux(xo + r, yo,o) >
S
The normalized SGS longitudinal fluid spatial velocity
correlation at different cut-off locations is plotted in
Figure 3. As shown in Figure 3, 6f(r) from the DNS
result has negative loops and oscillations in the tail. As a
first order approximation, we use an exponential relation to
approach 6f(r) in the SGS flow field,
8f(r) = exp (14)
For isotropic turbulence, as the continuity equation is
satisfied kfi = 0 for every wavenumber, k. The SGS
transverse fluid velocity spatial correlation function g(r) can
be obtained using the continuity relation,
r d(6f (r))
g(r)= 6 f(r) + (15)
2 dr
From Eq.(14) and Eq.(15), we can obtain
8g(r)= 1-- exp(---) (16)
S26Lf ) 6Lf
The SGS Lagrangian and Eulerian fluid velocity time
correlation functions are shown in Figure 4 and can be
approximately expressed as exponential functions,
Figure 3: Normalized SGS longitudinal fluid spatial
velocity correlation at different cutoff locations. The
abscissa has been normalized by the longitudinal integral
scale. The dash-dot-dotted line denotes the exponential
approximation.
1
0.8
- 0.6
0.4
- 0.2
0 0.01 0.02 0.03 0.04
Figure 4: SGS Eulerian and Lagrangian fluid velocity
temporal correlation at qkcf = 0.135 Reverse to full scale
flow field case, the SGS Eulerian temporal correlation
decorrelates faster than its Lagrangian counterpart in SGS
because the small scale eddies are convected by
energy-containing large eddies in turbulence.
6R, () =exp(- ) (17)
8T,
8R,(T) =exp(- )
aTE
With the expressions for SGS fluid velocity statistics,
R, (r) 8RE (r), 6f(r), 8g(r), we now consider the
SGS fluid velocity correlation seen by an inertial particle
with time delay r The mean spatial displacement parallel
Paper No
)
)1l\
7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30-June 4, 2010
to the gravity direction in the time interval r for an inertial
particle is r =(w) This displacement is due to both
gravity and turbulence preferential sweeping as discussed by
Wang & Maxey(1993). This relative displacement causes a
loss of SGS fluid velocity spatial correlation seen by the
inertial particle and should be reflected in 6Rp, (r). As
stated by Csanady (1963), in the case of
8Lf /(w) <
=Sf((w)r) parallel to the gravity direction and
8R p,2 (() = 8R Lp () = 6g(r)= Sg((w) r) perpendicular to
gravity direction. If replacing r/SL, and r with
r J/ST, +(w) SL, and (wr in Eq. (14) and Eq. (16)
respectively, we can obtain
SRL = exp(-r 1/ST,2 +(w) /L) (19)
5RL = 1- 26 exp(-rl/ST +(w) L) (20)
2Lf )
Integrating Eq. (19) and Eq.(20), we can get inertial
particle-SGS eddy interaction timescales along the direction
of the drift velocity STLP and normal to the direction of
the drift velocity STLp2 as
1
5T = I =
1 (w)V
ST2, 6L
J_ ^
STL2
(w) 1
( v) 1 (22)
28L1 1 (xv)
ST2 8+L
Validation of the interaction timescales between
inertial particle and SGS velocity field using DNS
data
In this section, we will validate the proposed relations of Eq.
(21) and Eq. (22) using the DNS data. The SGS Lagrangian
fluid velocity correlation along an inertial particle can be
calculated using Eq. (9). Figure 5 plots the comparison
between the analytical results and numerical results from the
DNS flow field.
It is shown that the proposed formulae well fit the numerical
results, and with the increase of the drift velocity, inertial
particle-SGS eddy interaction timescale in the direction of
body force STp1 decreases due to the crossing trajectory
effect (solid line and squares), the continuity effect further
decreases the inertial particle-SGS eddy interaction
timescale normal to the direction of body force (dashed line
0.006'
- 0.004
&
o
Go
0.002 -
~TY----o
: ^o^^
0 50 100 150 200
Figure 5: Variation of longitudinal and transverse
particle-SGS eddy interaction time scale STL,, and 6T,2.
and circles). In the limiting of very high drift velocity,
(w) o STp2 -> 1/26Tl .
Conclusions
This paper studies the effects of the drift velocity on the
interaction timescales between heavy particles and SGS
eddies in the condition of gravity.
The proposed formulae well fit the numerical results, and
with the increase of the drift velocity, inertial particle-SGS
eddy interaction timescale in the direction of body force
STL, decreases due to the crossing trajectory effect, the
continuity effect further decreases the inertial particle-SGS
eddy interaction timescale normal to the direction of body
force. In the limiting of very high drift velocity, (w) oo,
8Tp, > 1/2 TL .
The proposed formulae will be used to study LES of
turbulent dispersion of inertial particles under gravity for
next study.
Acknowledgements
This work was supported by CAS (KJCX-SW-L08), 973
Program of China (2007CB814800), NSFC (10732090,
10702074), the LNM initial funding for young investigators
and SRF for ROCS, SEM. G-D. Jin would like to thank
Professor L.-P Wang at the Department Mechanical
Engineering, University of Delaware for his help in the
course of this study.
References
Armenio, V, Piomelli, U., Fiorotto, V Effect of the subgrid
scales on particle motion. Phys. Fluids, 11, 3030-3042
(1999)
Paper No
S Analytical 6TLp1, Eq.(21)
S----- Analytical 6TLp2, Eq.(22)
SD Numerical 6TLpl
O Numerical STLp2
I-
Paper No 7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30-June 4, 2010
Berrouk, A.S., Laurence, D., Riley, J. J., Stock, D.E. Yudine, M. I. Physical considerations on heavy-particle
Stochastic modelling of inertial particle dispersion by diffusion. Adv. Geophys., 6, 185-191(1959)
subgrid motion for LES of high Reynolds number pipe flow.
J. Turbul., 8 (50), 1-20 (2007)
Bini, M., Jones, W. P Large-eddy simulation of
particle-laden turbulent flows. J. Fluid Mech., 614, 207-252
(2008)
Csanady, GT. Turbulent diffusion of heavy particles in the
atmosphere. J. Atmos. Sci., 20, 201-208 (1963)
Fede, P, Simonin, O. Numerical study of the subgrid fluid
turbulence effects on the statistics of heavy colliding
particles. Phys. Fluids, 18, 045103 (2006)
Fede P, Simonin O., Villedieu P, Crossing trajectory effect
on the subgrid turbulence seen by solid heavy particles, 6th
Int. Conf. on Multiphase Flow, ICMF2007, Leipzig,
Germany, July 9-13, 2007
Fede, P, Simonin, O., Villedieu, P, Squires, K.D. Stochastic
modeling of turbulent subgrid fluid velocity along inertial
particle trajectories. In: Proceedings of the 2006 CTR
Summer Program
Jin, G. D., He, G-W., Wang L.-P and Zhang J. Subgrid scale
fluid velocity timescales seen by inertial particles in
large-eddy simulation of particle-laden turbulence. Int. J.
Multiphase Flow, 36, 432-437(2006)
Pozorski, J., Apte, S.V Filtered particle tracking in isotropic
turbulence and stochastic modeling of subgrid-scale
dispersion. Int. J. Multiphase Flow, 35, 118-128(2009)
Reeks, M.W.. On the dispersion of small particles suspended
in an isotropic turbulent fluid. J. Fluid Mech., 83, 529-546
(1977)
Shotorban, B., Mashayek, F. A stochastic model for particle
motion in largeeddy simulation. J. Turbul., 7 (18), 1-13
(2006)
Wang L.-P, Maxey M.R. Settling velocity and concentration
distribution of heavy particles in a forced isotropic and
homogeneous turbulence. J. Fluid Mech., 256, 27-68 (1993)
Wang, Q., Squires, K.D. Large eddy simulation of
particle-laden turbulent channel flow. Phys. Fluids, 8,
1207-1223 (1996)
Wang, L.-P, Stock, A.D. Dispersion of heavy particles by
turbulent motion. J. Atmos. Sci., 50 (13), 1897-1913 (1993)
Wells, M.R., Stock, D.E. The effects of crossing trajectories
on the dispersion of particles in a turbulent flow. J. Fluid
Mech., 136, 31-62(1983)
Yamamoto, Y, Potthoff, M., Tanaka, T, Kajishima, T., Tsuji,
Y. Large eddy simulation of turbulent gas-particle flow in a
vertical channel: effect of considering inter-particle
collisions. J. Fluid Mech., 442, 303-334(2001)
__
__ |