PARTICLE-LADEN FLOW - ERCOFTAC SERIES Phần 6 ppt

Numerical particle tracking studies in a

turbulent round jet

Giordano Lipari, David D. Apsley and Peter K. Stansby

School of Mechanical Aerospace and Civil Engineering - University of Manchester -

Manchester - M60 1QD - UK [email protected]

1 Overview

This paper discusses numerical particle tracking of a 3D cloud of monod￾isperse particles injected within a steady incompressible free round turbulent

jet. With regard to particle-turbulence interaction, the presented modeling

is adequate for dilute suspensions [7], as the carrier and dispersed phase’s

solutions are worked out in two separate steps.

Section 2 describes the solution of the carrier fluid’s Reynolds-averaged

flow. The Reynolds numbers of environmental concern are generally high,

and here the turbulence closure is a traditional k- model ´a la Launder and

Spalding [2] with an ad hoc correction of Pope’s to the  equation to account for

circumferential vortex stretching in a round jet [21]. The resulting mean-flow

and Reynolds-stress fields are discussed in the light of the LDA measurements

by Hussein et al. (1994) with Re ∼ 105 [11].

Section 3 deals with the solution of the dispersed phase. The carrier fluid’s

unresolved turbulence is modeled as a Markovian process. We particularly

refer to the reviews of Wilson, Legg and Thomson (1983) [30] and McInnes

and Bracco (1992) [18]. Clouds of marked fluid particles, rather than traject￾ories, are used for visualizing the dispersing power of fluctuations. As fluc￾tuations in inhomogeneous turbulence are known to entail sizeable spurious

effects, the consistency of the Eulerian and Lagrangian statistics are checked

by comparing the first- and second-order moments of the particle velocity with

the mean flow and Reynolds stresses of the Eulerian solution, as well as the

concentration fields from either solution.

Surprisingly, our tests failed to confirm the full effectiveness of the correc￾tions proposed in either model. The particle spurious mean-velocity vanishes

towards the jet edge, thus abating the unphysical migration towards low￾turbulence regions. However, because of a residual disagreement between the

Lagrangian and Eulerian mean velocities, mass conservation entails concen￾tration profiles that do not follow the anticipated scaling. Possible reasons for

this are discussed in the closing section.

Bernard J. Geurts et al. (eds), Particle Laden Flow: From Geophysical to Kolmogorov Scales, 207–219.

© 2007 Springer. Printed in the Netherlands.

208 Lipari G. Apsley D.D. Stansby P.K.

2 Eulerian modeling of the Reynolds-averaged jet flow

The symbolism for the Reynolds averaging is ui = ui + u

i (i = 1 ... 3).

The Eulerian governing equations are the elliptic Reynolds-averaged mo￾mentum and continuity equations and the transport equations for the turbu￾lence scalars k and . In the Cartesian space and at the steady state, they

read

uj

∂ui

∂xj

= −1

ρ

∂p

∂xi

− ∂Rij

∂xj

+ νt

∂2ui

∂xixj

;

∂uj

∂xj

= 0.

The Reynolds-stress tensor Rij is modeled with the Boussinesq approxima￾tion:

Rij = −νt

∂ui

∂xj

+

∂uj

∂xi



+

2

3

kδij .

The eddy viscosity is based on the scaling νt = Cµk2/, where the fields of

the k and  turbulent scalars are computed with

uj

∂k

∂xj

= −Rij

∂ui

∂xj

−  +

∂

∂xj

ν + νt

σk

 ∂k

∂xj



,

uj

∂

∂xj

= −C1



k

Rij

∂ui

∂xj

− (C2 − C3χ)

2

k +

∂

∂xj

ν + νt

σ

 ∂

∂xj



. (1)

The constants Cµ through C2 and σ take the classic values of Launder and

Spalding [2]. The extra term having C3 in Eq.(1) depends on the vortex￾stretching invariant χ:

χ =

k

2

3 ∂ui

∂xj

− ∂uj

∂xi

 ∂uj

∂xk

− ∂uk

∂xj

 ∂uk

∂xi

− ∂ui

∂xk



,

put forward by Pope (1978) to reconcile the spreading rate of the axial ve￾locity profile [21], for which the uncorrected model would yield 0.11 against

the measured 0.094-0.096. We used a value C3 = 0.5, lower than 0.7, to

match more recent measurements than those originally used by Pope – Fig.

1a. Benefits and limitations of this correction are also discussed in [25].

The equations are solved in dimensionless form by normalization with the

nozzle diameter D and jet exit velocity u(0,0) and, exploiting axi-symmetry,

in the polar-cylindrical space (x, r, θ). The origins of both frames are placed

at the jet exit centerline.

The physical domain is the flow’s symmetry half-plane 100 × 20 diameter

long and wide respectively. A pipe protrudes into the domain for 8 diameters.

A structured grid of 200×90 suitably clustered, orthogonal cells is more than

adequate to resolve the expected gradients accurately. A plug-flow profile is

assigned as inflow condition.

The general-purpose in-house research code stream, thoroughly described

in [16], has been used to solve the above equations with a finite-volume

method. Suffice it here to say that the norms of the algebraic-equation re￾siduals could be brought down below the order of 10−13 routinely.

Numerical Particle Tracking in a Round Jet 209

2.1 Results

Fig. 1. Radial profiles in self-similarity variables of: a) ux; b) k; c) Rxx; d) Rrr;

e) Rθθ; f) Rxr. Thin lines: k- results at transects x/D = 55, 74, 83, 92. Bold line:

measured data fit by Hussein et al. [11]. Dashed line: selected profile at x/D = 74

without Pope’s correction. Symbols: measurements from [23] (), [15] (◦), [19] (),[9]

(×), [28] (), [10] (•), [24] () and [8] ().

Centerline values (not shown here) . The normalized axial velocity is expected

to decay as x−1. The inverse quantity u(0.0)u−1

(x.0) increases linearly with a slope

of 0.1564 very close to 0.1538 as measured. The virtual origin at x0 = 1.07D,

less than 4D as measured, implies a shorter zone of flow establishment.

210 Lipari G. Apsley D.D. Stansby P.K.

The turbulent kinetic energy k is to decay as x−2 [2], and the quadratic fit

of the inverse quantity is excellent from x = 10D onwards. Similarly, the rate

of turbulent energy dissipation  should decay as x−4, which is well reproduced

by computation; the fourth-order polynomial fit to the inverse quantity has a

leading-order coefficient of 0.0188 against 0.0208 as measured by Antonia et

al. for Re = 1.5 · 105 [5].

The turbulence timescale T = k/, therefore, increases as x2, e.g. in ac￾cordance with Batchelor’s analysis [6], with values ranging from 5 to 50 time

units. This derived quantity is central to modeling the autocorrelated part in

the fluctuation velocity – Eq. (2).

Radial profiles (Fig. 1). All plots are in self-similarity variables. Bold lines

represent the data fits of the benchmark experiment [11]. Continuous lines

show the computed quantities at selected far-field stations, which do collapse

on a single curve, achieving self-similarity. Dashed lines indicate the k- per￾formance without Pope’s correction.

Symbols are used to report the LDA measurements of high-Re single-phase

jets made available by some authors – Popper et al. (1974) [23], Levy and

Lockwood (1981) [15], Modarres et al. (1984) [19], Fleckhaus et al. (1987) [9],

Tsuji et al. (1988) [28], Hardalupas et al. (1989) [10], Prevost et al. (1996)

[24] and Fan et al. (1997) [8] – prior to studying the two-phase case.

Pope’s correction helps reduce to some extent the discrepancy between

measured and computed flow quantities. A C3-value to match the axial￾velocity spreading rate (the point of ordinate 0.5 in Fig. 1a) worsens the

prediction of the turbulent axial stress Rxx only (Fig. 1c), while those of Rrr,

Rθθ and Rxr improve to match the correct proportion with the scaling vari￾able u2

(x,0) (Fig. 1d-f). The off-axis peaks of Rrr and Rθθ are not supported

by the corresponding measurements though.

Further, the cross-comparison between the experimental data sets reveals

a noticeable disagreement between the benchmark and the two-phase studies

that, except for Fan et al.s, spread less than expected

3 Lagrangian modeling of the particulate cloud

Particles enter the domain at uniformly-distributed random positions on a

pipe cross-section with a chosen input rate N˙ (equal to 100 particles per unit

time as a baseline default). The flow properties at a particles position are

worked out by mapping the Cartesian position (x1, x2, x3) into the compu￾tational grid (x, r) and, then, working out the Reynolds-averaged dependent

variables with a bilinear interpolation. The resulting values are then mapped

back into the Cartesian space with the standard vector/tensor rotation opera￾tions. The local instantaneous fluid velocity ui is then created by summing ui

and u

i as obtained from Sec 2 and 3.1 respectively. The Lagrangian equations

of motion are finally resolved.