7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30June 4, 2010
Prediction of Tjunction mixing via Lagrangian tracking of tracer particles
A. Dehbi* and F. de Crecyt
Laboratory for ThermalHydraulics, Paul Scherrer Institut, Villigen 5232, Switzerland
t Laboratoire de Mecanique des Fluides et d'Acoustique, Ecole Centrale de Lyon, 69130 Ecully, France
abdel.dehbi@psi.ch and florence.decrecy@eclyon.fr
Keywords: Lagrangian simulations, particleladen fluid, CFD, mixing
Abstract
Turbulent fluid mixing downstream of a Tjunction is predicted by a hybrid EulerLagrange approach whereby tracer particle
trajectories are computed and mixing of the streams deduced from the relative concentration of particles originating from the
two branches. Both the RANS Reynolds Stress Model (RSM) and timeaveraged Detached Eddy Simulations (DES) are used to
obtain the mean flow field. The fluctuating fluid velocities seen by the particles are computed from the normalized Langevin
Continuous Random Walk (CRW) stochastic model. Simulation results are compared to experimental data on mixing of two
isothermal streams consisting of tap and deionized water, respectively. Initial simulations with RANSReynolds Stress Model
(RSM) yield strong underprediction of the mixing. Closer look at the results shows that the Reynolds stresses, which are
required inputs to the CRW, are poorly predicted with RSM in comparison with previous Large Eddy Simulations. Detached
Eddy Simulations (DES) were subsequently performed to provide the flow field, and the DESCRW model predictions compare
well with the mixing data.
TJunctions are commonly found in a variety of applications
which involve mixing of fluids of different properties. (e.g.
in chemical or nuclear plants). Safety issues are often
encountered in Tjunction flows, for example, thermal
stripping, i.e. piping surface damage, is a well known
phenomenon that is linked to large temperature fluctuations
which are generated when two fluid streams of significantly
different temperatures mix with one another in a Tjunction.
These fluctuations can cause cyclic thermal stresses and
eventually initiate cracking at the wall (Gelineau et al., 2002,
Chapuliot et al., 2005).
Owing to the complex turbulent nature of the flow, CFD
investigations of Tjunction mixing have focused in recent
years on the use of Large Eddy Simulations (LES).
Examples of recent LES studies of Tjunctions can be found
in Hu & Kazimi (2006), Lee et. al (2I ** ) or Niceno et al.
(lI= II ). Well resolved LES is principle captures most of the
turbulent structures which govern turbulent mixing, and
allows the determination of mixing dynamics parameters
such as oscillation frequencies and amplitudes.
In practice, mixing is obtained by simply solving for an
additional scalar equation (e.g. temperature) in the
framework of RANS or LES. In this study however, we use
Lagrangian tracer particle tracking instead. This study
therefore amounts to a further validation of the Langevin
Continuous Random Walk (CRW) model for particle
transport in turbulent fields (Dehbi, 2008). Both the
complicated flow in a Tjunction and the use of inertialess
particles render this a challenging test for the CRW model.
The carrier fluid mean field will be computed using both
RANSRSM and Detached Eddy Simulations (DES).
2. Geometry and boundary conditions of the
simulated experiments
The investigations reported in this paper aim at simulating
isothermal Tjunction mixing tests conducted at the Paul
Scherrer Institut (Zboray et. al, 2007, Walker et al., (2009)).
The tests yielded a detailed database that is ideally suited for
CFD code validation. The experimental setup is shown in
Figure 1. The piping is made from Plexiglass with an inner
diameter of 51 mm. The main and branch pipes have length
of 1000 mm and 500 mm, respectively. Tap water flows in
main branch pipe, whereas deionized water flows in the
branch pipe. The two water streams have very different
electrical conductivities, which allows one to measure
mixing with a specially designed wire mesh sensor (WMS).
The WMS consists of 16 parallel transmitter wires and 16
parallel receiving wires, both having a pitch of 3 mm. The
transmitting and receiving wires cross perpendicularly as
shown in Figure 2. Hence data is collected in 16 times 16
points of a given cross section. The WMS can be displaced
in any plane downstream of the Tjunction to collect mixing
data as a function of distance. The data processed from the
WMS signals allows the determination the instantaneous,
mean and standard deviation of the "mixing scalar", which
is proportional to the concentration of the main branch fluid.
It has value of 1 for the high conductivity main branch fluid
and 0 for the low conductivity for the side branch fluid.
1. Introduction
A drift correction term (s 1)
CD ra coefficient ()
d particle diameter (m)
F drgprunit mass (N kg )
g gravity acceleration (m s )
k turbulent kinetic energy (m's )
Re Renods number ()
U velocity (m/s)
u fluctuating velocity (m/s)
u* friction veloct (m/s)
t time (s)
x length scale (m), along mixing pipe
v lnth scale (m), normal to ie, symmetry lane
z lnth scale (m) normal to pp
Greek letters
p constants
a kinetic energy dissipation rate (m's
rl forcing moment ()
p viscosity (Pa s)
v kinematic viscosity (m's )
o rms of velocity, m s
p density (kg m
z time scale (s)
z, wall shear stress (Pa)
Gaussian random numbers ()
Superscripts
+ dimensionless
Sub scis
f fluid
dimension
L Lagrangian
p pricle
1streamwise direction
2 wall normal direction
3 spamvise direction
Figure 2: Schematic of the wire mesh sensor (WMS)
The principles of the WMS are explained in the paper by
Prasser et al. (1998). The test matrix consists of a dozen tests
with various combinations of flow rates. Only three, namely
tests T l, T2, and T3 (see Table 1) are simulated in this
investigation. These cover equal flow rates at high and low
Reynolds numbers, as well as a case where the main flow is
2.5 times larger than the branch flow.
Main branch Side branch Mixing pipe
flow flow Reynolds
Test number
kg/s kg/s
Tl 0.5 0.5 44000
T2 0.1 0.1 8800
T3 0.5 0.2 30700
Table 1: Experimental conditions of the simulated tests
7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30June 4, 2010
1:y
Ilinil Ibuncll
anon 111111'~ ourb 11
in, unun b
Figure 1: Schematic of the Tjunction test facility
top viw
2. Stochastic Lagrangian tracking of tracer
particles
Mixing in this investigation is determined via Lagrangian
tracking of tracer particles which are injected at the inlet
faces of the two branches and subsequently tracked. The
inertia of the particles has to be small enough in order for
them to faithfully follow the fluid flow and big enough to get
a reasonable integration time step. Since the carrier flow is
quite turbulent, the mean flow velocity alone will not give
satisfactory predictions of the mixing, and therefore
turbulent fluid fluctuations have to be modeled.
In an Euler/Lagrange approach, the tracking requires the
Nomenclature
A c
7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30June 4, 2010
equation in boundary layers, Durbin (1984) and Thomson
(1984), proposed that the Langevin equation be normalized
to account for strongly inhomogeneous turbulence.
Following Iliopoulos et al. :* IIi), the normalized Langevin
equation along the ith coordinate is written as:
specification of the instantaneous fluid velocity at particle
location. The latter can be decomposed as follows:
U= L+u
In a RANS framework, U is considered the timeaveraged
fluid velocity, whereas in an LES/DES approach, U can be
viewed as the timeaveraged value of the spatially filtered
velocity.
The determination of the particle paths necessitates also the
specification of the fluctuating velocity part u, and this is
described by the stochastic Langevin model outlined later in
this section,
We consider a rigid pointwise small inertia particle is
assumed to be entrained in a turbulent flow at isothermal
conditions. The only forces acting on the particle are taken
to be drag and gravity. Brownian diffusion is ignored as
particles in this study have diameters greater than 1 pLm. The
lift force is also neglected as the particles are of the same
unit density as the carrier fluid. The vector force balance on
a spherical particle reduces then to:
dU 18pu Re
= FD pp,df2 D pU)+ (2)
where U is the fluid velocity, U, the particle velocity, p, the
particle density, d, the particle geometric diameter, 4L the
fluid molecular viscosity, g the gravity acceleration vector,
and Re, the particle Reynolds number defined as:
dnur
at
In the above, u, is the fluid fluctuating velocity component,
o, the rms of velocity zL ai Lagrangian time scale, dll,
a succession of uncorrelated random forcing terms, and A,
the mean drift correction term which ensures the wellmixed
criterion (Thomson, 1987). For tracer particles with
vanishing inertia (Stokes number 0), A, can be written as
follows (Bocksell & Loth, 2006):
ttuu
a("
A, (6)
axJ
where repeated summation is performed on the j index only.
Dehbi (2008) has shown that the normalized Langevin
equations can be cast as follows for the bodyfitted
streamwise, normal, and spanwise directions of the
boundary layer (v' > 100):
u, dt 12
< ) +  d5z
d(a
u,
d( 1)
u, dt 2Z ao,
( ) +  g + dt
d~u3 _ dt df3(9
where the d ,'s are taken to be an uncorrelated succession of
Gaussian random numbers with zero mean and variance dt.
2.2 The Langevin equation in the bulk
In the bulk region, for which y 2100, turbulence can be
expected to be both anisotropic and inhomogeneous given
the complex Tjunction Hlow. No simplifications of the drift
coefficient can be made, since one is away from the
boundary layer, and hence the Langevin equations for the
bulk region can then be expressed in the computational
domain as follows:
d(u' ) u( 2 +, d( +A 4 dt (10)
o, o, Fr L.
v being the fluid kinematic viscosity. The drag coefficient iS
computed in the ANSYSFluent code (2008) from the
following equation:
CD 1
Re, Re:
where the p's are constants which apply to spherical
particles tor wide ranges of Re,.
The stochastic Langevin equations defining the fluctuating
velocity field along a particle track are presented briefly.
The domain is subdivided in two regions: the boundary layer
region with strongly anisotropic turbulence, and a bulk
region with generally anisotropic and inhomogeneous
turbulence. The Langevin equations will take different
forms depending on the location of the particle.
2.1 The Langevin equation in boundary layers
To improve the predictive capabilities of the Langevin
u, dt /2
( ) +, d + A, dt
o, z, V
u,
d~uz) uz 2 dg+AzA, dt
o, o, z, y
The drift terms A, as used as defined by equation (6) without
Iu, tt +drl, +Ardt
+ x dt
further approximation.
2.3 Specification of Eulerian rms and time scales
The flow in the boundary layer is modeled in the RANS
framework, and the prevailing Eulerian statistics are
assumed to be reasonably well approximated by those given
by Direct Numerical Simulation (DNS) investigations of
fully developed channel flows. This approximation yielded
reasonable estimations of the particle deposition rates in a
variety of geometries (Dehbi, 2008). The channel flow DNS
data by Marchioli et al. (13) were curvefitted to give the
rms of velocity by ratios of polynomials of order 3 to 5 and
match the data with a correlation coefficient better than
0.99:
+_ (13)
U bk7+k
y' is the wall distance in dimensionless units defined as:
.F y
(14)
y is the particle distance to the nearest wall, u* the friction
velocity derived from the wall shear stress Tw and wall fluid
density pc as follows:
u "
9/ P
(15)
The crossterm uzu2 1, in equation 7 was curvefitted in
similar fashion.
For the estimation of the fluid Lagrangian time scale 2L'
Bocksell and Loth (2006) have performed DNS calculations
in the boundary layer and showed that the Lagrangian time
scales in all directions are nearly equal and quite well
approximated by the fits obtained by Kallio and Reeks
(1989) and given by:
7,= 10 y < 5
z: = 7.122 + 0.5731 y
0.00129 y+2 5
7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30June 4, 2010
2.4 Coupling the CRW model to the ANSYSFluent code
The ANSYSFluent CFD code (ANSYS, 2008) provides the
mean flow parameters as well as a module to integrate the
particle equations of motion. The Langevin model described
earlier was implemented in ANSYSFluent as a User
Defined Function (UDF) subroutine which supplies the
trajectory calculation module with the fluctuating fluid
velocity seen by a particle at each time step. Details of the
coupling between the ANSYSFluent CFD code and the
stochastic Langevin model are provided by Dehbi (2008).
3. Results: RANSRSM simulations
3.1 Flow field
Prior to computing the tracer particle trajectories, the flow
field must be specified. In a first approach, the Reynolds
Averaged Navier Stokes (RANS) is attempted, because of
its modest computational requirements. The Reynolds
Stress Model (RSM) available in the CFD package
ANSYSFluent 12 (ANSYS, 2008) is used as it is
considered the most physically sound model in the RANS
approach, especially for flows with anisotropy and
curvature. In the RSM formulation, the Boussinesq
assumption is not used, but rather an equation is employed
for each of the Reynolds stresses, hence allowing the
computation of possible turbulence anisotropies. The
simulated domain consists of a Tjunction section with 10
diameter lengths in each of the branch.
The corresponding RANSRSM grid is made of fully
hexahedral cells, with grid clustering in the wall region to
ensure a y+ of about 12 for the walladjacent cells, and
1520 cells inside the boundary layer (up a y+ of 200), in
accordance with the Best Practice Guidelines for boundary
layer resolved simulations (ERCOFTAC, 2000). Grid
refinement was higher in the mixing region and coarser
away from it. The ANSYS preprocessor Gambit was used
to build three hexahedral meshes: a coarse, a medium and a
fine mesh. These have 0.34, 1.1 and 3.7 million cells,
respectively.
Test Tl is first simulated in the RANS approach, and it is
characterized by equal main and branch flows of 0.5 kg/s.
The Reynolds number is 22 000 for inlet and branch pipes,
44 000 for outlet pipe. A fully developed velocity profile is
imposed on both inlets, and second order discretization is
used. The steadystate solution is driven to tight
convergence whereby all scaled residuals are smaller than
104, and area weighted average wall shear stress and wall
y+ reach iterationindependent plateaus. For mixing
purposes, the ydirection which is the horizontal direction
normal to the main mixing pipe is most important, because it
has the steepest mixing scalar gradients. Thus results are
displayed along that direction throughout this study.
The streamwise xvelocity profiles along the yaxis at axial
distance x = 7. 1 cm (L/D=1.42) are shown in Figure 3 for all
where :
The Lagrangian time scale in the bulk is calculated from the
total turbulent kinetic energy k and dissipation rate E as
follows:
2 k
, E
A value of 14 for Co provides good agreement with time
scales computed DNS investigations (Mito and Hanratty,
2002).
7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30June 4, 2010
three meshes. All three simulations capture the recirculation
zone (negative velocity), and it can be concluded that grid
independence is achieved. The RANS simulations are
subsequently conducted with the medium grid having 1.1
million cells.
0.025
0.02
0.015
0 g1
0.005 
0.005
0.01
0.015 

0 0.5 1 1.5 2
x Velocity [m/s] x=151mm
Figure 5: x velocity at x=15.1 cm. Comparison of
RANSRSM and LES
0.5 1
Velocity [m/s]
Figure 3: x velocity vs. y at x=7.1 cm. Effect of grid
resolution
The accuracy of the mean flow results is checked against
previous LES computations by Niceno et al.( 2009). Hence,
results for the mean xvelocity along the yaxis are shown in
Figures 4 and 5 for two axial positions, namely x = 5.1 cm
(L/D=1) and x = 15.1 cm (L/D=3). The RSM and LES
profile are in reasonable agreement, and hence the RSM
flow field was subsequently used along with the CRW
model to deduce the mixing via Lagrangian tracking of
tracer particles.
Particles having diameter 20 pLm and density 1000 kg/m3 are
used as tracers. The corresponding particle dimensionless
relaxation time is:
_' r, u~ P~dZ u (19)
2, is the particle relaxation time and 4L the fluid molecular
viscosity. Taking the maximum mean velocity to be 1.5 m/s,
the maximum friction velocity is of the order of 0. 1 m/s.
Hence, the particle dimensionless relaxation time 2,' is at
most 0.2, which means the used particles act as tracers.
Particles that hit wall (less than 1% in all simulations) are
removed from the domain and considered lost.
0.025 _
0.02
0.015 
0.01 
0 a0
0 
0.005 
0.01 
0.015 
0.02 
100000 particles uniformly distributed across the two inlet
faces are injected and subsequently tracked. This sample
was found to be large enough not to affect the results.
Particles are injected for a total of 10 flowthrough times T,
patil t rah heotltT corresponding to the average time it takes for a fluid
paticl toreac th outet.At a long enough time to
establish stationary conditions, typically t=5T, one freezes
the particles tracks and particle concentration is obtained by
simply counting the number of particles in a small volume at
any desired cross section.
S RSM
SLES To get information on mixing in a given plane, a cylinder of
1.5 volume equaling the pipe cross sectional area times a width
of 0.02 cm is defined at different axial locations
(measurement locations) and discretized into 16 by 16
identical cubes to replicate the measurement distribution of
the WMS. The mixing profiles are then obtained by simply
parisn of dividing the number of particles originating from the main
branch by the total number of particles (branch+main).
Contour profiles of the mixing scalar as a function of the
axial distance x (0.051m, 0. 191 m and 0.311 m) are shown
in Figure 6 for both the simulation and the experimental data.
Plots of the mixing scalar along the central yaxis line are
also shown.
0 0.5 1
x Velocity [m/s] x=51mm
Figure 4: x velocity
RANSRSM and LES
at x
5.1 cm. Com
3.2 Lagrangian particle tracking
Experimental data, x=0 051m
0028
0 01
6
4
001
002 001 0 001 002
y [cm]
Experimental data, x=0 191m
002
8
0 01
0~6
0 01 4
002
002 001 0 001 002
y [cm]
Simulation, x=0 051m
00208
0 01
06
04
001
002 001 0 001 002
y [cm]
Simulation, x=0 191m
002
08
0 01
0 06
00
002 001 0 001 002
y [cm]
Comparison along y axis, x=0 051m
experiment
1 computation
002 001 0 0 01 0 02
y[m]
Comparison along y axis, x=0 191m
experiment
1 computation
0 8
002 001 0 0 01 0 02
y[m]
Comparison along y axis, x=0 311m
1 experiment
computation
002 001 0 0 01 0 02
y[m]
Experimental data, x=0 311m
002
8
0 01
0 01 6
2
002
002 001 0 001 002
y [cm]
Simulation, x=0 311m
002
08
0 01
000
002 001 0 001 002
y [cm]
7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30June 4, 2010
Figure 6: Mixing scalar for Test Tl at different axial locations:
Left: experimental data
Middle: RANSRSM simulation
Right: comparison of data and simulation along yaxis
7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30June 4, 2010
4. Detached Eddy Simulations (DES)
4.1 Fluid flow
More advanced models need to be used to provide not only
accurate mean flow predictions, but also accurate Reynolds
stresses. One obvious approach is the Large Eddy
Simulation (LES) method which aims at directly resolving
the large anisotropic structures while modeling the smaller
isotropic turbulence dissipating eddies. However, properly
resolving the nearwall region at moderate Re numbers still
demands large CPU resources, especially if many
simulations are to be performed, as is the case here. An
alternative way is the Detached Eddy Simulation (DES)
(Spalart et al, 1997) approach which has become an
attractive alternative to LES for a wide range of applications.
DES is a hybrid RANSLES method which uses unsteady
RANS models for the nearwall regions and LES treatment
in the bulk region. The motivation of such an approach is
that RANS models in the boundary layer region are accurate
enough and require only moderate grid resolution to yield
acceptable results. In the bulk region, LES is used to
properly resolve large coherent and generally anisotropic
structures which drive particle motion. In the bulk LES
region, the remaining unresolved fraction of turbulent
kinetic energy is modeled using so called subgrid scale
(SGS) models. The model employed in this investigation is
the Realizable ke.
A DES computation was conducted to simulate Test 1 as in
the RANS study discussed earlier. A hexahedral grid with
1.5 million meshes is constructed with appropriate grading
to resolve the boundary layer. The grid has adequate
resolution to capture the most of the turbulent kinetic energy
in the bulk LES region. The time step is 0.0005 s, which
provides a Courant number that is less than 1 throughout the
domain, hence ensuring timeaccurate results. After
stationary conditions were reached, all required Eulerian
statistics were gathered for 3 flowthrough times. Figures 9
and 10 show the symmetry plane instantaneous and mean
xvelocity, respectively. As in previous investigations, a
recirculation region forms immediately downstream of the
Tjunction, due to momentum exchange between the two
steams. By the time the flow reaches about x/D= 2, the
recirculation zone disappears and the turbulent mixing
becomes dominant.
Figure 11 shows comparison of the DES results and the LES
(Niceno, 2009) for the rms of wallnormal velocity. It can
thus be concluded that the DES provides reasonably close
results to LES for the Reynolds stresses that control
turbulent diffusion of tracerlike particle.
Subsequently, two additional experimental tests were
simulated, i.e. Test 2 and Test 3. The same grid was used in
these simulations that have a smaller Reynolds numbers,
thus ensuring that the turbulent scales are adequately
resolved.
In the experimental data, the scalar profile in not perfectly
symmetric with respect to the yaxis. This is due to a small
temperature difference between the streams, which causes
one stream to slightly tumn over the other as the flow
progresses into the mixing pipe. On the other hand, the
simulation predicts roughly symmetric profiles for the scalar,
as itshould.
It is observed that the predicted and experimental mixing are
comparable up to x=0.051 m (L/D=1). This is not surprising,
as mixing very near the Tjunction is dictated to a great
extend by inertia rather than turbulence. Beyond this point,
the simulation becomes continuously less accurate, as
turbulent mixing becomes dominant.
To see why turbulent mixing is not well predicted by the
RSM model despite reasonably accurate mean field, we go
back to the Langevin model equations and the important
drift correction factor A (equation 6). This factor ensures
that tracer particles remain fully mixed if inj ected uniformly,
i.e. that they don't undergo the socalled "spurious
drift" which will cause them to artificially accumulate in
low turbulence regions. If one considers, to a first order, that
the diagonal Reynolds stresses dominate, the factor A for the
ydirection reduces to 802 By Where (32 is the rms of
velocity in the wallnormal ydirection. We thereafter
compare the LES and RSM rms of the yvelocity along the
yaxis at x= 0.051 m and x = 0. 151 m in Figures 7 and 8,
respectively, and observe that they are quite different. Since
turbulent dispersion of fluid particles is governed by the
gradients of the turbulent kinetic energy (socalled
turbophoresis), it is natural to expect that the turbulent
mixing is poorly predicted by the RSM if one takes the LES
results as a reference.
Closer look at Figures 7 and 8 shows that for the LES
profiles, the wallnormal velocity gradient is essentially
negative for y > 0, i.e. for the region where the flow is
originating from the main branch. This means that the tracer
particles will on average experience a push towards the
center of the pipe, hence promoting mixing. Conversely, the
rms gradient is by and large positive for y < 0.01, i.e. for
the region where the flow is originating from the side branch.
Hence tracers in that region will experience a push towards
the center of the pipe, enhancing mixing as well. In the RSM
simulations, the magnitudes of the gradients are smaller, and
the sign of the gradient goes sometimes in the
"wrong direction", i.e. toward enhancing demixing, clearly
and unphysical result.
In conclusion, this preliminary investigation highlights the
fact that accurate predictions of the mixing require not only
an accurate mean flow field but also a good estimation of
the Reynolds stresses, which is not achieved with the
RANSRSM approach. Hence a better model of turbulence
is required.
7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30June 4, 2010
04
0.2
00
0.
Z
Figure 7: rms of y velocity at x
RANS and LES
5.1 cm. Comparison of
7 0.000
> 0.005
0.010
0.015
0.020
0.025i,
0.00 0.05 0.10 0.15 0.20 0.25
v, [m/s]
Figure 11: rms of y velocity at x=15. 1 cm. Comparison of
DES vs. LES profiles
Figure 8: rms of y velocity at x=15.1 cm. Comparison of
RANS and LES
17
I
02
0.4
I o
Z L
Figure 9: Instantaneous xvelocity at symmetry plane
4.2 Results: Lagrangian particle tracking
Lagrangian particle tracking is conducted in similar fashion
as in the RANS simulation, with the mean field and mean
Reynolds stresses based on the timeaveraged data of the
DES computation. Again, particles having unit density and
diameter 20 pLm are used as tracers. For Test 2 with equal
carrier fluid flow rates (0.1 kg/s), 100000 particles
uniformly distributed across the two inlet faces are injected
and subsequently tracked. For Test 3, 250000 particles for
the main branch and 100000 particles for the side branch are
injected, keeping therefore the proportion between the fluid
flow rates at the branches inlets. Apart from these
adjustments, the procedure to compute the particle
COncentrations distribution at different downstream sections
is identical to what was described in the RANSCRW
procedure earlier.
Figure 12 shows the comparison between the DESCRW
and the data for Test (0.5 kg/s, 0.5 kg/s). We note a
significant improvement in the predictions compared with
the RANS simulations. This is largely due more accurate
computations of the Reynolds stresses. It is nonetheless
Figure 10: Timeaveraged xvelocity at symmetry plane
7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30June 4, 2010
clear that the DESCRW computations still slightly
underpredict the mixing rate, a feature that was also seen in
the LES investigation (Niceno, 2009). Still, the predictions
are largely within the rms band of the experimental data.
Figures 12 and 13 show contour plots as well as profiles of
the mixing scalar along the yaxis for tests T2 and T3. Again,
the results show that the DESCRW model predicts well the
experimental mixing data both qualitatively and
quantitatively. In particular, simulations are able to capture
the rapid mixing when the main and side branches have
comparable Hlow rates (Test Tl and T2), and the slow
mixing when the branches have significantly different flow
rates.
9. Conclusions
Turbulent fluid mixing downstream of a Tjunction is
predicted by a hybrid EulerLagrange approach whereby
tracer particle trajectories are computed and mixing of the
streams deduced from the local concentration of particles
originating from the two branches. Both the RANSRSM
and timeaveraged DES are used to obtain the mean flow
field. The fluctuating fluid velocities seen by the particles
are computed from the normalized Langevin CRW
stochastic model which is well suited to tackle flow
inhomogeneities.
Simulation results are compared to experimental data on
mixing of two isothermal streams consisting of tap and
deionized water, respectively. Three cases are simulated,
covering equal Hlow rates at high and low Reynolds numbers,
as well as a case in which the main Hlow is 2.5 times larger
than the branch flow.
Initial simulations with RANSRSM yield strong
underprediction of the mixing. Closer look at the results
shows that the Reynolds stresses, which are one of the inputs
to the CRW, are poorly predicted with RSM in comparison
with previous LES computations. DES simulations were
then performed to provide the Eulerian field, and the model
predictions compare well with the whole range of data.
It can therefore be concluded that the Langevin CRW for
particle transport reproduces realistic fluid velocity
fluctuations in complex flows whenever the underlying Hluid
mean Hlow and Reynolds stresses are computed accurately
by the fluid turbulence model.
Simulation, x=0 051m
002
8
0 01
6
4
001
2
002
002 001 0 001 002
y [cm]
Simulation, x=0 191m
002
8
0 01
6
002 001 0 001 002
y [cm]
Comparison along y axis, x=0 051m
experiment
1 computation
02~y~
002 001 0 0 01 0 02
y[m]
Comparison along y axis, x=0 191m
experiment
computation
002 001 0 0 01 0 02
y[m]
Experimental data, x=0 051m
002
08
0 01
06
04
001
02
00
002 001 0 001 002
y [cm]
Experimental data, x=0 191m
002
08
0 01
06
00
002 001 0 001 002
y [cm]
Experimental data, x=0 311m
002
08
0 01
06
0 04
00102
0 02
00
002 001 0 001 002
y [cm]
Simulation, x=0 311m
002
8
0 01
6
0 4
001 2
0 02
002 001 0 001 002
y [cm]
Comparison along y axis, x=0 311m
experiment
computation
0T 8
002 001 0 0 01 0 02
y[m]
Figure 12: Mixing scalar for Test Tl (0.50.5 kg/s) at different axial locations:
Left: experimental data
Middle: DES simulation
Right: comparison of data and simulation along yaxis
7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30June 4, 2010
Experimental data, x=0 051m
00208
0 01
06
04
0 01
02
00
002 001 0 001 002
y [cm]
Experimental data, x=0 191m
002
08
0 01
04
00102
002 001 0 001 002
y [cm]
Experimental data, x=0 311m
002
08
0 02
000
002 001 0 001 002
y [cm]
Simulation, x=0 051m
00208
0 01
06
04
0 01
02
00
002 001 0 001 002
y [cm]
Simulation, x=0 191m
002
08
0 01
04
00102
002 001 0 001 002
y [cm]
Simulation, x=0 311m
0 02
8
0 02
002 001 0 001 002
y [cm]
Comparison along y axis, x=0 051m
illexperiment
1 computation
002 001 0 0 01 0 02
y[m]
Comparison along y axis, x=0 191m
1 experiment
1 computation
002 001 0 0 01 0 02
y[m]
7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30June 4, 2010
Figure 13: Mixing scalar for Test T2 (0.1
Left: experimental data
Middle: DES simulation
0. 1 kg/s) at different axial locations:
Right: comparison of data and simulation along yaxis
Experimental data, x=0 051m
002
08
0 01
06
04
001
02
002
002 001 0 001 002
y [cm]
Simulation, x=0 051m
002
8
0 01
6
4
001
2
002
002 001 0 001 002
y [cm]
Simulation, x=0 191m
0028
0 01
6
4
001
2
002
002 001 0 001 002
y [cm]
Simulation, x=0 311m
002
8
0 01
0~6
4
0 01
2
002
002 001 0 001 002
y [cm]
Comparison along y axis, x=0 051m
1 experiment
computation
002 001 0 0 01 0 02
y[m]
Comparison along y axis, x=0 191m
experiment
1 computation
002 001 0 0 01 0 02
y[m]
Comparison along y axis, x=0 311m
experiment
computation
1
002 001 0 0 01 0 02
y[m]
Experimental data, x=0 191m
00208
0 01
06
04
001
02
00
002 001 0 001 002
y [cm]
Experimental data, x=0 311m
002
08
0 01
0 06
04
0 01
02
00
002 001 0 001 002
y [cm]
7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30June 4, 2010
Figure 13: Mixing scalar for Test T3 (0.50.2 kg/s) at different axial locations:
Left: experimental data
Middle: DES simulation
Right: comparison of data and simulation along yaxis
7th International Conference on Multiphase Flow
ICMF 2010, Tampa, FL USA, May 30June 4, 2010
Durbin, P.A., Comments on the papers by Wilson et al.
(1981) and Legg and Raupach (1982). Boundary Layer
Meteorology 29, 409411 (1984).
Thomson, D.J., Random walk modeling of dispersion in
inhomogeneous turbulence. Quarterly Journal of the Royal
Meteorological Society, 110,11071120 (1984)
Iliopoulos I., Mito Y., Hanratty T.J., A stochastic model for
solid particle dispersion in a nonhomogeneous turbulent field.
International Joumnal of Multiphase Flow 29, 375394
Thompson, D.J., Criteria for the selection of stochastic
models of particle trajectories in turbulent flows. J. Fluid
Mech. 180, 529559 (1987)
Bocksell T.L., Loth E., Stochastic modeling of particle
diffusion in a turbulent boundary layer. International Joumnal
of Multiphase Flow, 32, 12341253 (2006)
Marchioli, C., Picciotto, M., Soldati, A., Influence of
gravity and lift on particle velocity statistics and transfer
rates in turbulent vertical channel flow, Intemnational Joumnal
of Multiphase Flow 33, 227251 (2007)
Kallio, G.A., Reeks, M.W., A numerical simulation of
particle deposition in turbulent boundary layers.
International Joumnal of Multiphase Flow 3, 433446 (1989)
Mito Y. Hanratty T.J., Use of a modified Langevin equation
to describe turbulent dispersion of fluid particles in a channel
flow. Flow, Turbulence and Combustion, 68, 126 !* II I)
ERCOFTAC. European Research Community on Flow
Turbulence and Combustion (ERCOFTAC) Best Practice
Guidelines, Version 1 (2000)
Spalart, P.R., Jou, W.H., Strelets, M., Allmaras, S.R.,
Comments on the feasibility of LES for wings, and on a
hybrid RANS/LES approach. Advances in LES/DNS, First
AFOSR International Conference on DNS/LES. (1997)
Acknowledgements
The authors are grateful to Dr. R. Zboray for making
available his experimental data base.
References
Gelineau, O., Sperandio, M., Simoneau, J.P., Hamy, J.M.,
Roubin, P.H.L., Thermomechanical, and Thermalhydraulic
Analyses of a Tjunction using Experimental Data, in
IAEATECDOC1318, 3162 ('1 'II)
Chapuliot, S., Gourdin, C., Payen, T., Magnaud, J.P.,
Monavon, A., 200 5. Hydrothermalmechanic al analy si s of
thermal fatigue in a mixing tee. Nucl. Eng. Des. 235,
575596 (2005)
Hu, L.M., Kazimi M. S., LES benchmark study of high cycle
temperature fluctuations caused by thermal striping in a
mixing tee, Intemnational Joumnal of Heat and Fluid Flow, 27,
5464 (2006)
Lee, J.I., Hu, L.W., Saha, P., Kazimi, M. S.,
Numerical analysis of thermal striping induced high cycle
thermal fatigue in a mixing tee, Nuclear Engineering and
Design, 239, 833839 (2009)
Niceno B., Zboray R. and Prasser HM., Experimental and
Numerical Mixing Studies Relevant for Prediction of
Thermal Loads in Tjunctions, Proceedings 6th Intemnational
Symposium on Turbulence, Heat and Mass Transfer, Rome,
Italy (2009)
Dehbi, A., Turbulent particle dispersion in arbitrary
wallbounded geometries: a coupled CFDLangevin
equation based approach, Intemnational Journal of Multiphase
Flow, 34, 819828 (2008)
Zboray, R., Manera, A., Niceno, B., Prasser, H.M.,
Investigations on Mixing Phenomena in SinglePhase Flows
in a Tjunction Geometry,NURETH12, Pittsburgh
Pennsylvania, USA, September 30October 4, paper 71,
(2007)
Walker, C., Simiano, M., Zboray, R. Prasser H.M.,
Investigations on mixing phenomena in singlephase flow in
a Tjunction geometry, Nuclear Engineering and Design,
239,116126 (2009)
Prasser, H.M., Bottger, A., Zschau, J., A new
electrodemesh tomograph for gasliquid flow. Flow Meas.
Instrument. 9, 111119 (1998)
ANSYS Fluent 12 Users Guide, 2008.
