In this contribution, we have performed Large Eddy Simulations (LES) of the dense particle laden turbulent flow
as described by prof. Jennifer Sinclair Curtis (University of Florida, USA) and co-workers. The fluid phase is
modelled with three different LES models, i.e. a Smagorinsky model where van Driest dampening is employed
near the walls, and the more recent WALE and Vreman models. Also, the effect of mesh spacing on the flow
prediction is researched. The particles are modelled by the Discrete Element Model (DEM), in which the particle
collisions are approximated by their local elastic and plastic deformation. The importance of collisions is shown by
comparing the simulations taking particle collisions into account to simulations neglecting this effect. The effect
of the geometry is also researched. Due to the dimensions of the box in which the turbulent pipe flow enters,
it is expected that this dimension has a significant effect on the spreading rate of the turbulent gas-particle jet.
Finally, the particle size is varied, in accordance with the available experimental data. The particles used in the simulations
range from 25 μm to 70 μm as well as a mixture of these. All findings are compared with the experimental data.

General Note:

The International Conference on Multiphase Flow (ICMF) first was held in Tsukuba, Japan in 1991 and the second ICMF took place in Kyoto, Japan in 1995. During this conference, it was decided to establish an International Governing Board which oversees the major aspects of the conference and makes decisions about future conference locations. Due to the great importance of the field, it was furthermore decided to hold the conference every three years successively in Asia including Australia, Europe including Africa, Russia and the Near East and America. Hence, ICMF 1998 was held in Lyon, France, ICMF 2001 in New Orleans, USA, ICMF 2004 in Yokohama, Japan, and ICMF 2007 in Leipzig, Germany. ICMF-2010 is devoted to all aspects of Multiphase Flow. Researchers from all over the world gathered in order to introduce their recent advances in the field and thereby promote the exchange of new ideas, results and techniques. The conference is a key event in Multiphase Flow and supports the advancement of science in this very important field. The major research topics relevant for the conference are as follows: Bio-Fluid Dynamics; Boiling; Bubbly Flows; Cavitation; Colloidal and Suspension Dynamics; Collision, Agglomeration and Breakup; Computational Techniques for Multiphase Flows; Droplet Flows; Environmental and Geophysical Flows; Experimental Methods for Multiphase Flows; Fluidized and Circulating Fluidized Beds; Fluid Structure Interactions; Granular Media; Industrial Applications; Instabilities; Interfacial Flows; Micro and Nano-Scale Multiphase Flows; Microgravity in Two-Phase Flow; Multiphase Flows with Heat and Mass Transfer; Non-Newtonian Multiphase Flows; Particle-Laden Flows; Particle, Bubble and Drop Dynamics; Reactive Multiphase Flows

In this contribution, we have performed Large Eddy Simulations (LES) of the dense particle laden turbulent flow
as described by prof. Jennifer Sinclair Curtis (University of Florida, USA) and co-workers. The fluid phase is
modelled with three different LES models, i.e. a Smagorinsky model where van Driest dampening is employed
near the walls, and the more recent WALE and Vreman models. Also, the effect of mesh spacing on the flow
prediction is researched. The particles are modelled by the Discrete Element Model (DEM), in which the particle
collisions are approximated by their local elastic and plastic deformation. The importance of collisions is shown by
comparing the simulations taking particle collisions into account to simulations neglecting this effect. The effect
of the geometry is also researched. Due to the dimensions of the box in which the turbulent pipe flow enters,
it is expected that this dimension has a significant effect on the spreading rate of the turbulent gas-particle jet.
Finally, the particle size is varied, in accordance with the available experimental data. The particles used in the simu-
lations range from 25 pm to 70 pm as well as a mixture of these. All findings are compared with the experimental data.

Particle-laden turbulent flows can be found in various
industrial and environmental processes. Examples
of such processes are pneumatic transport, fluidised
beds; energy generation of fossil fuels; movement of
soot particles in the atmosphere; the flow of particles
in cyclones and many more. Understanding particle-
particle, particle-fluid and particle-wall interactions is of
outmost importance as a more accurate implementation
of these processes will result in more reliable analysis
tools of such flows. Reliable numerical simulations
will therefore help the optimisation and better design
of industrial processes and also reliable prediction of

environmental processes.

There are various frameworks in which gas-solid
flows can be predicted in the Eulerian-Lagrangian
framework, e.g. Large Eddy Simulation (LES), Direct
Numerical Simulation (DNS) and the Reynolds Aver-
aged Navier-Stokes (RANS) method. RANS methods
do not resolve the small scales, which are important
for turbulent flows, but only consider the averaged
quantities of the transport equations. The smaller scales
are modelled by appropriate closure models. DNS
methods, on the other hand, offer much higher accuracy
in solving the Navier-Stokes equations at the expense
of huge computational time; currently, DNS can only
solve flows of low Reynolds (Re) numbers, which are
outside of most the engineering and industrial interests.
Although the computational effort for LES is still
very high, it is considerably lower than for DNS and
therefore has become very fashionable for analysis of
flows in research and is an emerging tool in industrial
problems.

LES solves the Navier-Stokes equations up to a
particular length-scale due to the application of a filter.
Length-scales smaller than the cut-off width (A) of the
filter are modelled. This cut-off width is an indication
of the size of eddies that are retained in the computa-
tions and the size of eddies that are rejected Versteeg
and Malalasekera (2007). Due to the filtering of the
Navier-Stokes equations, models are required to resolve
the sub-grid-scale (SGS) stresses to account for the
effect of unresolved scales on the convective momentum
transport. In this paper, three different models are used
and compared with one another. These models are
the well known Smagorinsky model where van Driest
damping is employed near the wall; the Wall-Adapting
Local-Eddy (WALE) viscosity model proposed by
Nicoud et al Nicoud and Durcos (1999); and finally
the Vreman model proposed by Vreman Vreman 21 11 14).

The particle interactions are modelled by the Dis-
crete Element Model (DEM), first applied by Cundall
and Strack Cundall and Strack (1979). The particle-
particle and particle-wall collisions are approximated
by the elastic and plastic deformation occurring during
collision. Such deformation is mathematically described
by a spring-dashpot-slider model, or soft-sphere model,
as proposed by Tsuji et al. Tsuji et al. (1991).

Due to the large number of particles, the compu-
tational time and memory space required for these
simulations is very large. This is because computations
for each particle take place for every timestep. On
the other hand, DEM has a very high accuracy as no

empirical data, besides from material properties, are
required to predict the particle-particle and particle-wall
interactions. The model is also fully deterministic -
making an accurate comparison of the effect of initial
conditions possible. To describe the fluid-particle inter-
actions, the Wen and Yu Wen and Yu (1966) drag force
has been employed in this work. Although accounting
for the volume fraction effects on the drag force is
probably not of large importance in the test-cases
simulated in the current research work, it might have an
effect on the particle-clustering. The forces resulting
from the fluid-particle interactions, the particle-particle
interactions and the particle-wall interactions are added
and used to track the individual particles employing the
Verlet integration scheme.

Theory

Fluid-phase modelling

The equations for the fluid-phase are similar to the
single-phase, except for including a local volume frac-
tion and including source terms accounting for the pres-
ence of particles. The governing equations are the mass
continuity and conservation of momentum, which are
given by equations (1) & (2) respectively,

where af is the local volume fraction of the fluid phase;
Sf the source term linear in the velocity field; T3 the
additional source terms; and 7)3 the stress tensor of the
fluid, given by

' 9 vi a vi (
Tf 8 'I xi V

2f ) a (3)
3 8x"V

where pf is the shear and Af represent the bulk viscosity
of the fluid.

Subgrid Scale Models

Large Eddy Simulation (LES) is used to solve the
mass-continuity and momentum equations, models are
therefore needed to determine the subgrid-scale (SGS)

stresses. Three models are used and compared with ex-
perimental data. These are the Smagorinsky model with
van Driest damping near the wall; the WALE model
Nicoud and Durcos (1999); and finally the Vreman
model Vreman (21' '4). It is worthwhile mentioning that
the models are compared with one another by consid-
ering two different geometries, by varying the length of
the inlet pipe. Predictions have been made for an inlet
pipe of length L = 20.0 and of L = 90.0. In addition,
different mesh spacing is researched in order to examine
the effect on the flow prediction.

Van Driest Wall damping
It is well-known that the Smagorisnky model is not
suitable to account for the effect of walls, as the strong
velocity gradient perpendicular to the wall is present
due to the no-slip boundary condition at the wall and
not a gradient driving the turbulence. Hence, near the
wall the so-called van Driest dampening is employed
which is based upon a simple mixing length model.

where y is the length from the wall, in the case of pipe
flow given by,
y (R -r) (5)
and pT represents the eddy viscosity, which is

t c2 dU
fiT i t dy

where l~ is the mixing length. The mixing lei
scale is approximated as

/mix = ly
Near the wall, however, ,ix is not well described
Equation (7). This is beacuse the mixing length me
'seems incapable of reproducing the y+3 variation r
the wall' Grifoll and Giralt (2000). The mixing length
therefore, modified by the van Driest damping funct
This function dampens the effect of viscosity near
wall. Hence equation (7) is modified to equation (8)

mix ts-[1- c /A o]

where y+ is the dimensionless distance defined as,

and A+ is a constant, typically A+ 26. Although
the model seems very crude, it has proven to be quite
successful in many types of wall bounded turbulent flows.

Smagorinsky model
The Smagorinsky model assumes that the local SGS
stresses are proportional to the local rate of strain of the
resolved flow Versteeg and Malalasekera (2007). Thus,
the stresses are given as

Tij = -2psGsSij + 3T16Sij (10)

where Sj = + %j psGs is estimated using
equation (11)

PSGS = Pf (CsGSA)2 2 (11)

where Cscs is the Smagorinsky constant, A is the LES
filter width.

The Smagorinsky van Driest constant The
Smagorinsky model, however, does not accurately rep-
resent the eddy viscosity near the wall. As shown by
Weickert et al Weickert et al. (2009), the CSGS is modi-
fied as,

Nicoud et al Nicoud and Durcos (1999) proposed
igth the wall-adapting local eddy-viscosity (WALE) model.
With this model, the non-physical behaviour of the eddy
viscosity as predicted by the Smagorinsky model near
(7) the wall is avoided. This is because the WALE model
Sby is of higher order compared to the Smagorinsky model,
)del which is 0(1). In this way, the eddy viscosity at the wall
ear is zero and hence, the only viscosity contribution at the
h is, wall is the molecular viscosity. The model for the eddy
ion. viscosity is given in Nicoud and Durcos (1999) and We-
the ickert et al. (2009),
Id d j
(SG) s pf (C'A)2 (S (14)
(8) (S (+14

C, is related to the Smagorinsky constant via,

C. = (, G ) (S -ad)3)
(SiJSi (Sdd
w-ithl -- 4j Y-- and./.-- -

with^d 1 +^)1 ._9 a .

Vreman model

Vreman Vreman (2 '14) proposed an eddy viscosity
model used in LES, which is also of higher order accu-
racy and produces the correct limit of the eddy viscosity
near the walls. This model assumes that the dissipation
is relatively small in transitional and near-wall regions.
The eddy viscosity is calculated as

LES utilises low-pass filters so that only the large scale
eddies are resolved, whereas the smaller scale eddies are
removed. In the aforementioned models, a filter width
is used in order to define the cut-off lengthscale. In this
paper the volume of the cells used in the meshing is di-
rectly related to the filter width (A), i.e. A .
As for successful LES calculations, the computational
cell aspect ratio has to remain near unity, this is a good
approximation.

Particle modelling

A discrete element model (DEM) proposed by Cun-
dall and Strack Cundall and Strack (1979) is used to
model the particles. The individual trajectories of the
particles are determined in the Langrangian framework,
where particle collisions are modelled via the soft-
sphere model proposed by Tsuji et al Tsuji et al. (1991),
which accounts for some non-reversible deformation.

Discrete Element Models (DEM) The interactions of
particles with other particles and walls are of dynamic
nature. This is because the particle movements are
essentially defined by the partcle-particle interactions;
particle-wall interactions; particle-fluid interactions
and/or body forces Kuang et al. (2008). The trajectories
of individual particles are considered (i.e. described
in a Langrangian framework) and Newton's 2"d law is
solved for each individual particle, accounting for the
fluid-particle, particle-particle, and particle-wall inter-
actions, and approximating the integral with the Verlet
algorithm. The flowchart shown in figure 1 illustrates
the DEM algorithm. The flowchart describes the logic
of how forces between individual particles with other
particles and/or walls of the domain are determined. In
addition, it illustrates how forces between particles and
the fluid are computed.

where mp is the mass of the particle; Ip is the momen-
tum of rotational inertia; Tp the torque of the particle;
wp is rotational velocity; vp is the translational velocity;
and 3 is the drag function as proposed by Wen and Yu
Wen and Yu (1966), where the reciprocal of the Eulerian
fluid-particle timescale is given by

and CD represents the coefficient of drag for an individ-
ual particle and af represents the fluid volume fraction.
The right hand side terms of equation (17) are outlined
in table 1.

Implementation of particle collisions The particle-
particle and particle-wall interactions as taken into
account in this work are assumed to be local; i.e. no
long-range forces are included. For establishing the
nature of the interaction, each particle pair could be
interrogated. However, this would lead to a scaling of
the computational effort with N2, N being the number
of particles. Instead, a particle-mesh algorithm is
adopted, in which each of the particles is assigned a cell
in the particle mesh based upon its location. Using this
particle-mesh, each particle is only tested for interaction
against particles located in the same or directly neigh-
bouring particle mesh cells. Although there is some
additional computational effort and a slightly more
complex algorithm required for this approach, it leads
to a scaling of N logN with the number of particles,
making it a lot more favourable for large numbers of
particles.

The particle collisions are modelled by the soft-
sphere model as described by Cundall and Strack
Cundall and Strack (1979). In brief, this model uses
a spring-dashpot-slider arrangement to describe the
particle behaviour before, during and after a collision.
The damping coefficient used, described by Tsuji et
al Tsuji (1994) & Tsuji et al. (1991) is, however, used
since it is related to the coefficient of restitution.

The normal and tangential contact forces are given
by the sum of forces due to the spring and dashpot.
From Hertzian contact theory the normal and tangential
contact forces are,

Fij = (-kn6 ,i. G.n)n (20)

Ftij = -kt6t G (21)
where G is the velocity vector of particle i relative to
particle j, Get is the slit velocity vector at the con-

tact point. The subscripts n and t represent the normal
and tangential components respectively and 6 is the dis-
placement, or overlap, of particles i and j during colli-
sion. Three parameters are required by the soft-sphere
model; stiffness (k), damping coefficient (rI) and the co-
efficient of friction (p). p is a well known empirical
quantity, stiffness and the damping coefficient must be
estimated using equations

/n = a VMkn

*t = aV'Mkt64

where M = n and m is the mass of the particle.
Note that for wall collisions m oc, hence M = n.
The relationship between a and the coefficient of restitu-
tion is well defined by Tsuji et al Tsuji et al. (1991). The
spring constants are, based upon elastic deformation,

k 4 1 4(1 1
3 E

- )
Ey

Sir rj

S2 aT 2-a, 1 -V r+rj) 2' I
kt 8( + (rir) 6n (25)
Sj / \ irj /

where r is the radius of the particles a the Poisson's ra-
tio; E the Young's modulus, G the shear modulus (given
by G = 2(+E ) and 6, the magnitude of the deforma-
tion in the normal direction. Note as above, when a par-
ticle reaches a wall, r, oc, hence',+' -1.

Description of test cases

The objective of the simulations performed in this work
is to investigate the behaviour of a turbulent two-phase
flow containing borosilicate glass beads with diameters
of 25 pm and 70 pm. Simulations with three different
test cases are pursued, which are summarised in table 2.

Geometric data The main test-case involves a
particle-laden gas jet issuing into a test chamber, as is
described in the experiment ?. The dimensions of the ge-
ometry as described in the experimental work are sum-
marised in table 3. Refer to figure 2 for the lay-out of
the experiment.
For the simulations, two geometries have been used
for the simulation cases, a pipe of length = 20 and a
pipe of length = 90. Both pipes exit in a chamber of
width and depth of = 6.5 and length of = 13. As
the chamber in which the jet exits is significantly smaller
than the chamber used in the experimental work, full slip
boundary conditions and a zero pressure gradient have
been used as boundary conditions on the chamber walls,
to minimize the influence of their presence.

Copper pipe, L/D = 90

Test Chamber

S0.300 m

0.470 m

0.470 m

Figure 2: Test chamber dimensions

Properties of the particles The physical properties of
the borosilicate glass beads as used in the simulations are
summarised in table 4.
'coefficient of restitution for inelastic particle-particle collisions
2coefficient of restitution for inelastic particle-wall collisions

Properties of the fluid The fluid used in the simula-
tions is air at standard conditions and its properties are
summarised in table ??.

Computational mesh

For each simulated geometry, two computational meshes
were created. All meshes are suitable for resolving the
wall boundary layer, and have 3 mesh points within the
y+ 5 layer. The coarsest geometry for the pipe with
length = 20 has 55,000 computational cells, and the
refined mesh for the same geometry has 260,000 cells.
The coarse mesh of the geometry with the pipe length of
S 90 has 200,000 mesh cells and the fine mesh has
around 600,000 mesh cells. Because resolving all the
individual particle trajectories takes around 85% of the
computational time, the refinement of the mesh does not
have a significant impact on the computational time.

Discretisation

The simulations performed in this paper use a finite vol-
ume scheme, a second order backward Euler time dis-
cretisation and a second order accurate central scheme
to approximate the advective terms in the continuity and
momentum equations. The solving procedure is made
parallel with the MPI libraries. During the simulations,
the CFL number is kept constant at 0.4, leading to a
slight variation in time-step over time.

Boundary and Initial Conditions

Experimental data for the single-phase flow are reported
at the pipe exit. This reported velocity profile is used to
determine the exact mass flow and the centerline veloc-
ity to be specified at the inlet of the simulations. Due to
mass conservation, the mass flow rate at the pipe exit is
constant in the pipe. Using equations (26) & (27), two
simultaneous equations are set up to determine the von
Karman (K) and B constants.
/R-6,, R
T= pf 27rrVl(r) dr + pf 27rrV (r) dr
(26)

U, U, ln(RU ) + B (27)
K Vf

Using the calculated values of K and B and equations
(28) and (30), the inlet velocity profile is produced.

VI(r) In((R
K

r) )+ B,
Vf

valid for 0 < Irl < (R 6,), where 6, is the linear-
sublayer,

2
V (r) = (R- r) (30)
vf
valid for (R 6,) < Ir < R. The frictionvelocity (u,)
is also required to determine the inlet velocity profile.
This is calculated using the following relationship,

U7 = (31)
P 1

where the shear wall stress is estimated using the friction
factor for a smooth pipe at Re = 8,400.

1 2 (3
Tw fpfV2 (32)
8

Note that the friction factor is estimated using Prandtl's
friction law for a fully developed turbulent flow Pope
(2000) for a smooth pipe. The values for the above quan-
tities are given in table 6.
The matched velocity profile is set as an inflow
condition at the pipe inlet. In addition, in a number of
simulations a synthesized turbulent component is added
to this flow. The synthesized turbulence is generated by
sampling isotropic fluctuations, using the Fourier modes

of the fully developed spectrum. The number of Fourier
modes which are considered (N) was 600, the turbulent
kinetic energy (qm) is 0.5, the boundary layer thickness
(Lt) is 0.0014 m and the timescale correlation of the
turbulence was estimated at T 2.48 x 10 4 s.

In the chamber a laminar flow of constant velocity
0.06 ms 1 through the top of the box is defined. In the
experiment of Hadinoto et al Hadinoto et al. (2005),
the top of the box is left open, allowing in order to
allow for air entrainment from the surrounding flow.
Unfortunately, the exact characteristics of this air flow
are not known, making it difficult to predict the correct
spreading-rate of the jet as it leaves the pipe.
The walls of the pipe are enforced with a no-slip
condition for the velocity and pressure is extrapolated
onto the walls from the fluid. The pressure boundary
conditions were interpolated from downstream at the
inlet of the pipe and the inlet of the chamber walls. It is
found that interpolation from downstream improves the
results over the more common zero-gradient boundary
conditions for pressure.
As the chamber in which the jet exits in the simulations
is significantly smaller than the chamber used in the
experimental work, full slip boundary conditions and
a zero pressure gradient have been used as boundary
conditions on the chamber walls, to minimize the
influence of their presence.

Particle Injection The particles are injected at the
pipe inlet (noted as L/D = 0) with a solids loading
very near of m = 1.0. The injection timing is set so as
to compensate for the difference in the solids loading at
each injection. This will produce an overall solids load-
ing very near to 1.0 spread out through the geometry.
Assuming a uniform placement of particles at the inlet,
the location of each particle at the pipe inlet is known. A
random component is added to the spatial coordinates of
the particle. Also, the particles are given a random com-

ponent on all velocity directions. Some form of initial
dispersion, caused by velocity fluctuations in the radial
and azimuthal direction, is crucial to obtain an accurate
representation of the flow behaviour of the particles.

70 pm particles A particle size distribution (PSD) is
used for the 70 pm particles (the particle diameter range
was between 62 pm to 88 pm). The number based distri-
bution used was uniform (i.e. the number based proba-
bility to find a particle is the same). Since the PSD must
be volume based, the cumulative density function (cdf)
of the number based distribution is raised to 1/3 in or-
der to obtain the volume based cdf. Note that the solids
loading for the 70 pm particles is 1.0.

25 pm and 70 pm particles The procedure de-
scribed for the 70 pm particles is repeated for the bi-
nary mixture. A PSD is used with a cdf of 1.'. 25
pm and 5% 70 pm. This is because a solids loading
of 0.5 is required for each fraction. Therefore the ratio
of the number of particles for each fraction is given by,
n7 = 0.05.
25 (dp70

Results and Discussion

Up to now, simulations for test-cases 2 and 3 have been
run, although not all test-cases have sufficient statistics
to draw firm conclusions. The simulations are run on 8
cores in parallel.

Comparison of mesh refinement

Two different meshes were used to simulate the small
pipe geometry; the coarse mesh having approximately
55,000 mesh cells and the fine mesh having approx-
imately 260,000 mesh cells, meaning a refinement of
about 1.75 in all directions. Figure 3 shows the effect of
mesh refinement on the predicted centerline velocity at
the pipe exit. It can be seen that the mesh refinement has
very little effect. A similar conclusion can be drawn for
the results presented in Figure 4, showing the magnitude
of the fluctuations at the pipe exit, and Figure 5 showing
the centerline velocity downstream from the pipe exit. It
can be seen that especially the fine mesh needs to be run
for longer, as the statistics are not yet steady.

Comparison of SGS LES models

Three different LES models for closing the SGS stresses
are employed. The results are shown in Figures 6,7, 8 &
9. In Figure 16 all models seem to overpredict the mass
flow rate at the pipe exit. However, a good correlation
between the experimental results is achieved especially

Figure 3: The simulation results of the mean axial ve-
locity at the pipe exit (S 0) from the fine and coarse
mesh, compared with the single phase and multiphase
data.

Figure 4: The simulation results of the root mean square
of the velocity fluctuations at the pipe exit ( = 0)
from the fine and coarse mesh, compared with the single
phase and multiphase data.

the Vreman model. This is better illustrated in Figure 8
where the Vreman model is in good agreement with the
experimental results (particularly near the pipe exit).
The axial and radial rms velocities, however, are not pre-
dicted well by the models. This is because there are not
enough statistical data and also this comparison is for
the short pipe.

Effect of pipe length

Fluid Statistics There is a strong effect of the length
of the pipe on the fluid velocity statistics. Figure 10
shows the mean axial flow velocity directly at the pipe
exit for the short pipe ( 20) and the long pipe
( 90) compared to experimental data. The data of

Figure 6: The simulation results of the average center-
line velocity at the pipe exit ( =' 0) as a result from the
different SGS closure models compared to experimental
data.

the long pipe is a very good match to the experimen-
tal data. For both cases, it seems that the mass flow is
slightly overestimated at the inlet of the pipe. Figure 11
shows the root mean square of the velocity fluctuations
in the axial direction resulting from the simulation with
the short pipe compared, the simulation with the long
pipe and the comparison with experimental data. Both
predictions are reasonable and there is no clear differ-
ence. Figure 12 shows the root mean square of the veloc-
ity fluctuations in the radial direction resulting from the
simulation with the short pipe compared, the simulation
with the long pipe and the comparison with experimen-
tal data. Both predictions underestimate the magnitude
of the fluctuations.

Figure 7: The simulation results of the root mean square
of the velocity velocity fluctuations at the pipe exit exit
( = 0) as a result from the different SGS closure mod-
els compared to experimental data.

Figure 8: The simulation results of the average center-
line velocity after the pipe exit (from = 0) as a result
from the different SGS closure models compared to ex-
perimental data.

Particle Statistics There is also a strong effect of
the pipe length on the particle velocity statistics; the re-
sults possibly show that the longest pipe, with length
S= 90 is not long enough to achieve fully developed
flow for the particles. Figure 13 shows the mean parti-
cle velocity in the axial direction for the long pipe, the
short pipe, and the experimental data. There is a strong
improvement for the longer pipe, but the match is not
completely satisfactorily.
Figures 14 and 15 show the root mean square of the par-
ticle velocity in the axial and radial directions, respec-
tively. Although the trends in both figures seem to be
correctly predicted, there seems to be an offset in the
predictions compared to the experimental data.

Figure 10: The simulation results of the average fluid
axial velocity component at the pipe exit ( = 0) from
the short pipe ( 20) and the long pipe ( 90)
compared to experimental data.

Bimodal particle mixture

A case employing the Smagorinsky SGS model includ-
ing van Driest Dampening was performed including a
binary particle mixture, with particles of dp 25pm
and dp 70pm with a 50:50 relative volumetric pro-
portion with a total mass loading of 1.0.

Fluid velocity statistics Figures 16 and 17 show
the mean fluid velocity in the axial direction and the
root mean square of the fluid velocity fluctuations in the
axial direction, respectively. Both are compared to the
experimental findings. The results show a fairly good
prediction.

The prediction of the root mean square of the ve-
locity fluctuations, compared to the experiments, as

Figure 11: The simulation results of the root mean
square of the fluid fluctuating axial velocity component
at the pipe exit ( = 0) from the short pipe ( = 20)
and the long pipe ( = 90) compared to experimental
data.

shown in figure 18 is not as good, probably due to a lack
of statistical values.

Particle statistics As the experiments have been
performed with phase doppler anemometry (PDA), the
velocity statistics of the small and the large particles in
the binary mixture have been determined separately. Al-
though very similar, some changes have been found in
the statistics of the small particles compared to the large
ones experimentally.

Conclusions

This paper presents Large Eddy Simulation (LES)
modelling work and comparing this work to the ex-
perimental data as presented by prof. J. Curtis and
co-workers. The simulations have been run without any
empirical parameters and material properties have been
chosen to accurately reflect the properties of the fluid
and the particles.

Generally, the simulations are too short to draw
strong conclusions. However, a few conclusions can
be drawn from the results presented so far. The first
conclusion is that the simulation results seem to be
fairly independent of mesh resolution. Changing the
mesh from 55,000 cells to 260,000 cells hardly affects
the results. Moreover, the results from the WALE and
Vreman LES models accurately capture the presence of
the wall, without the neccesity of applying van Driest
dampening. The application of such wall dampening
is computationally expensive in parallel simulation runs.

Figure 12: The simulation results of the root mean
square of the fluid fluctuating radial velocity component
at the pipe exit ( = 0) from the short pipe ( = 20)
and the long pipe ( 90) compared to experimental
data.

An important aspect of the simulations, is the size
of the geometry. The simulation with the shortest
geometry, with a pipe of 1 20, is too short for
the fluid and particles to be fully developed. There
are also concerns that the chamber in which the pipe
exits is too small to reflect an accurate spreading of the
jet. The uncertainty of the flow in the chamber in the
experimental work makes it difficult to match with the
simulations.

Future work

The most important goal for the future is to perform
longer simulations, so better statistics can be obtained
for the largest geometry. Once this is achieved, the ef-
fect of varying the following properties in the simulation
will be performed:

The size of the chamber. The size of the cham-
ber will affect the spreading rate of the jet. If this
will also affect the dispersion of particles will be
researched.

The van Driest dampening parameters at the wall.
This may have an important effect on the center line
velocity.

The magnitude of the SGS constants. All the SGS
models contain a constant, and although in theory
the magnitude of the constant should not effect the
flow characteristics this is likely not to be true for
this case, which has a fairly low Re number.

The properties of the particles. The fluid and par-
ticle velocities are highly affected by the collisions

Figure 13: The simulation results of the average axial
particle velocity component at the pipe exit ( = 0)
from the short pipe (1 = 20 and the long pipe ( = 90)
compared to experimental data.

Figure 14: The simulation results of the root mean
square of the particle fluctuating axial velocity com-
ponent at the pipe exit ( = 0) from the short pipe
( = 20 and the long pipe ( 90) compared to ex-
perimental data.

Figure 18: The simulation results of the root mean
square of the fluid fluctuating radial velocity component
at the pipe exit (j = 0) of the flow containing the binary
mixture of particles compared to experimental data.