7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010
Optimal Eulerian model for the simulation of dynamics and coalescence of
alumina particles in solid propellant combustion
Frangois Doisneau*t, Fr6d6rique LaurentNegret, Angelo Murrone*,
Joel Dupays* and Marc Massott
Department of Fundamental and Applied Energetics, ONERA, 91120 Palaiseau, FRANCE
SLaboratoire EM2C, Ecole Centrale Paris and CNRS, 92295 ChatenayMalabry, FRANCE
francois.doisneau@onera.fr, frederique.laurent@em2c.ecp.fr, angelo.murrone@onera.fr,
joel.dupays@onera.fr and marc.massot@em2c.ecp.fr
Keywords: Dilute Polydisperse Spray; Coalescence; Solid Propulsion; Alumina; Eulerian MultiFluid Model
Abstract
The accurate simulation of polydisperse sprays encountering coalescence in unsteady gaseous flows is a crucial
issue for solid rocket booster optimization. Indeed, the internal flow of the engine depends strongly on the alumina
droplet size distribution, which spreads up with coalescence. Yet solving for unsteady twophase flows with a
high dimensional phase space is a challenge for both modelling and scientific computing. The usual Lagrangian
approaches lead to a very high computational cost or to a low resolution level and they induce coupling difficulties to
the Eulerian gaseous phase description. A wide range of Eulerian models have been recently developed to describe
the disperse liquid phase at a lower cost and with an easier coupling to the carrier gaseous phase. Among these
models, the MultiFluid model allows the detailed description of polydispersity and size/velocity correlations by
separately solving fluids of sizesorted droplets, the socalled sections. On the one hand, the existing first order
description of the size distribution in each section provides simple and fast resolution for coalescence. On the other
hand, a second order method allows to reduce the number of sections required to capture accurately coalescence and
to use elaborate droplet collision modelling, yet at the cost of heavier computation algorithms. This paper seeks to
conclude on computational time and precision of both methods in order to choose the most efficient configuration
for multidimensional unsteady rocket chamber simulations. Its objective is threefold: first, to validate the second
order method by comparing simulations to reference solutions and dedicated experimental measurements conducted
at ONERA, second to study the efficiency and robustness of both methods on coalescing sizeconditioned dynamics
problems, third, to draw some firm conclusions about the necessity to use first order or second order methods in order
to capture the physics of solid propulsion configurations.
Introduction
Twophase flows constituted of a gaseous phase and a
disperse condensed phase play a key role in many in
dustrial and scientific applications : fuel atomization
and evaporation by Diesel injectors; fluidized beds; gas
bubbles in oil or boiling water pipeworks; dynamics of
planet formation in solar nebulae; etc. In all these ap
plications the disperse phase is composed of particles of
various sizes that can eventually coalesce or aggregate,
breakup, evaporate and have their own inertia and size
conditioned dynamics. So the importance of polydisper
sion is obvious for a full modelling of these phenomena.
In rocket boosters, aluminum powder is frequently
used as solid propellant additive to increase specific im
pulse. Unlike the other ingredients, aluminum parti
cles can burn in a sigmlik.iii portion of the chamber
and produce a condensed liquid disperse phase of alu
mina. This disperse phase encounters drag forces, coa
lescence and heat exchanges (Dupays et al. 2000). In the
nozzle, the droplets accelerate suddenly and cool down
with the gas, becoming solid particles which breakup
in such velocity gradients (Dukowicz 1980; O'Rourke
1981; Amsden et al. 1989). Thus, the disperse phase
strongly interacts with the gaseous flow field during its
way throughout the engine. It contributes to the booster
performance loss, particularly via a decrease in nozzle
efficiency; droplets are the source of slag material that
may remain in the engine during firing, causing insula
tion erosion in high concentration zones; some of the
droplets, mainly the ones with high inertia which end up
in the eventual aftdome region around the submerged
nozzle, induce sloshing motion of this molten liquid
slag and can lead to control problems and possible ve
hicle instability. In such harsh conditions of pressure,
temperature and velocity, solid propulsion experiments
consume high technology materials and offer poor mea
surement output, especially on the disperse phase. The
abundance of physical phenomena involved makes the
models difficult to scale (T6th 2008). Regarding the
prohibitive cost of experiments, numerical simulation is
the only available tool for optimizing rocket engines.
Until now, complete 3D computations were achieved
at the cost of drastic physics simplifications, allowing
little selfreliance towards validation experiments. Yet
time has come for comprehensive simulations, including
advanced gas/droplet/structure coupled models, to give
predictive answers.
Focusing on the dynamics of the alumina cloud in the
combustion chamber is sufficient to evaluate specific
impulse loss in the nozzle and slag material accumula
tion. Therefore the study takes place in the combustion
chamber, where the principal physical processes that
must be accounted for are : transport in real space,
acceleration of droplets due to drag conditioned by
size and coalescence, leading to polydispersity. These
processes are remarkably sensitive to size distribution.
We will therefore choose a model accurate as regards
the size distribution and work with nonevaporating
sprays throughout the paper, keeping in mind the broad
application fields related to this study. By spray, we
denote a disperse liquid phase constituted of droplets
carried by a gaseous phase. We consider the specific
case of dilute sprays i.e. where the liquid volume
fraction is much smaller than one.
The retained approach called "mesoscopic" or some
times "kinetic" in reference to the kinetic theory of
gases describes the droplets as a cloud of point parti
cles for which the exchanges of mass, momentum and
heat are described using a statistical point of view, with
eventual correlations : a finite set of global properties
such as size of spherical droplets, velocity of the cen
ter of mass, temperature are modelled so that the total
phase space is usually highdimensional. More details
about the droplets, such as angular momentum, can be
predicted by increasing the size of the phase space : it
is established that refined droplets models can be used
as long as they do not include history terms. Williams
(1958, 1985) proposed a transport equation based on
kinetic theory that has proven to be useful for treating
polydisperse, dilute and moderately dense liquid sprays.
7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010
Such an equation describes the evolution of the number
density function or NDF of the spray due to the drag
force of the gaseous phase and the dropletdroplet in
teraction of coalescence (Kuentzmann 1973; Hylkema
1999; Laurent and Massot 2001; Laurent et al. 2004).
There are several strategies in order to solve the liq
uid phase and the major challenge in numerical simula
tions is to account for the strong coupling between all
the involved processes. A first choice is to approximate
the NDF by a sample of discrete numerical parcels of
particles of various sizes through a LagrangianMonte
Carlo approach (Dukowicz 1980; O'Rourke 1981; Ams
den et al. 1989; Hylkema 1999; Rilger et al. 2000). It is
called Direct Simulation MonteCarlo method (DSMC)
by Bird (1994) and is generally considered to be the
most accurate for solving Williams equation; it is spe
cially suited for direct numerical simulations (DNS)
since it does not introduce any numerical diffusion, the
particle trajectories being exactly resolved. Its main
drawback is the coupling of a Eulerian description for
the gaseous phase to a Lagrangian description of the dis
perse phase, thus offering limited possibilities of vec
torization/parallelization and implicitation. Besides, it
brings another difficulty associated with the repartition
of the mass, momentum and heat source terms at the
droplet location onto the Eulerian grid for the gas de
scription. Moreover for unsteady computations of poly
disperse sprays, a large number of parcels in each cell
of the computational domain is generally needed, thus
yielding large memory requirement and CPU cost.
As an alternative, the Eulerian MultiFluid model,
furthered by Laurent and Massot (2001); Laurent et al.
(2004) from the ideas of Greenberg et al. (1993), relies
on the derivation of a semikinetic modelling from the
Williams equation using a moment method for velocity,
but keeping the continuous size distribution function
(Laurent and Massot 2001). This distribution function
is then discretized using a "finite volume approach" in
size phase space that yields conservation equations for
mass, momentum (and eventually other properties such
as number, energy) of droplets in fixed size intervals
called sections, each of them constituting a different
"fluid". Please note that integrating on a continuous size
variable in each section is a key aspect: while other Eu
lerian approaches often consider discrete droplet sizes
gathered into .L,,,s, which cannot account for the
new droplet sizes created by coalescence (except scarce
examples such as in Vasenin et al. (1995)), continuous
size approaches such as the sectional method hereafter
described are the only Eulerian methods handling
coalescence naturally and rigorously. After integration
on the sections, the resulting conservation equations
are similar to those of the pressureless gas dynamics
(Bouchut (1994); Zel'dovich (1970)) and lead to singu
lar behaviors such as delta shocks and vacuum zones.
Wellsuited numerical methods are thus required, some
specific schemes being presented in De ( I.iiciii.illiiii
(2009). The model finally requires closure relations for
the phenomena accounted for by the Williams equation.
We refer to Abramzon and Sirignano (1989); Laurent
and Massot (2001); Laurent et al. (2004) for detailed
droplet models for which the MultiFluid model can
be easily extended. These extensions do not have any
impact on the conclusions of the present study.
The Eulerian MultiFluid model has proven its capa
bility to simulate the sizeconditioned dynamics of poly
disperse sprays including coalescence with a first order
resolution method of the size distribution in each sec
tion (Laurent and Massot 2001) even if the number of
sections with the original method has to be large for ac
curacy purposes. In the context of evaporating sprays,
several second order methods have been developed. A
second order method using exponential size distributions
in the sections is discussed in Laurent (2006). These
two methods are in fact particular cases of a general
method requiring n moments and introduced in Massot
et al. (2010). The first order method can easily solve
coalescence (Laurent et al. 2004) while the second or
der method has been extended to coalescence at the cost
of some more elaborate algebra in Dufour (2005). It can
yet easily include advanced models such as collision and
coalescence efficiencies developed in Hylkema (1999);
Hylkema and Villedieu (1998).
For the purpose of this paper, the first and second
order methods have been implemented in a research
code solving dilute sprays in a pseudo 2D nozzle with
oneway coupling to the gaseous phase. The physical
features have yet been simplified : the drag force is
modelled by a Stokes law; the temperature, compo
sition, density and viscosity of the gaseous phase are
assumed to be constant and uniform; as a consequence
unsteady heating of the droplets is negligible so that
their temperature and density are also constant. This
configuration is used to compare the two methods
solving the dynamics and size evolution of a lognormal
distributed spray. Moreover it fully validates the second
order method by providing detailed comparisons to
lagrangian and MultiFluid reference solutions. The
research code also allows to simulate the complete
dynamics of an experiment on coalescence (D'Herbigny
and Villedieu 2001) which had been initially designed
to validate collision efficiency models with an averaged
analytical formula. This work is thus a continuation
of the experiment exploitation and provides precious
validation for such twophase flow models, which too
often lack experimental back up. Moreover this test
case confirms the robustness of the method towards
7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010
sharp or monodisperse distributions. These validations
globally show the compliance of the second order
MultiFluid method to the features required by rocket
booster simulations i.e. accuracy on polydispersity and
dynamics, advanced coalescence models, robustness
and fast computation. In addition, the second order
MultiFluid model had previously been featured in
CEDRE, a multiphysics 3D industrial code developed
at ONERA, the French aerospace lab. CEDRE pro
vides fully coupled aerothermochemical resolution
for energetic problems. It includes advanced models
and can therefore compute drag force, evaporation and
coalescence terms depending on all the gas parameters.
The paper thus presents rocket booster simulations
performed with an optimized version of the second
order MultiFluid to achieve the proof of the efficiency
of such models in a complex physical background.
The paper is organized as follows. Section one is ded
icated to the derivation of the MultiFluid method as a
theoretical framework : the origin and assumptions of
the corresponding coalescence models are detailed, as
well as the resulting coalescence formulas for first and
second order methods. In section two, we evaluate the
accuracy/speed compromise provided by the two Multi
Fluid methods implemented in the research code on a
simple pseudo2D validation test case used in Laurent
et al. (2004). In section three, we validate the second or
der MultiFluid method and its collision and coalescence
efficiency models on a 1D coalescence efficiency exper
iment conducted at ONERA and detailed in D'Herbigny
and Villedieu (2001). In section four, and to illustrate the
ability of the second order model to simulate in a reason
able time a complete solid propulsion test case, we give
some numerical results in a modelled chamber and noz
zle using the CEDRE code before concluding the study.
1 Mesoscopic Eulerian spray modelling : two
MultiFluid methods
In this section, we introduce the framework of our
study : the kinetic description of the disperse phase and
the derivation of two Eulerian resolution methods. The
separation between the two methods appears when pre
suming the size distributions on the sections. The origin
and assumptions of the corresponding coalescence mod
els are widely detailed.
1.1 A kinetic description : the Williams
equation
Number density function When focusing on poly
dispersity, the size parameter p of droplets is of capi
tal importance but its natural expression depends on the
7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010
phenomena: volume v is relevant towards conservation
of matter, surface S towards evaporation and radius r
towards impact factor for instance. Since we assume
spherical droplets, equivalence relation v = 4ir3
1 S3
6, S 2 allows the size to be expressed in this paper with
the most comfortable notations.
Let us define the NDF function f9 of the spray, where
f (t, x, u, y)dxdydu denotes the averaged number of
droplets (in a statistical sense), at time t, in a volume
of size dx around a space location x, with a velocity
in a duneighborhood of u and with a size in a do
neighborhood of p. As for the size parameter conver
sions, we shall keep in mind that f (r)dr fs(S)dS
f" (v)dv and use the implicit notation f.
Williams equation The evolution of the spray is de
scribed by the Williams transport equation (Williams
1958). For a non evaporating case, it reads :
Jf + u Oxf + 0 (F f) = (1)
where F is the drag force per unit mass and F is the
coalescence source term.
As an illustration, the Stokes law models drag force
F due to the velocity difference with the gaseous phase
and reads :
F(t, x, u, S)
U(t,x) u
TTp(S)
r (S)
plS
18Pi
187 p9
where U is the gas velocity, pg represents its dynamic
viscosity and pi is the droplet (supposedly liquid) mate
rial density.
Coalescence operator The kinetic modelling for
the collision operator leading to coalescence is taken
from (Hylkema and Villedieu 1998). We then assume :
[Cl] We only take into account binary collisions (small
volume fraction of the liquid phase)
[C2] The mean collision time is very small compared to
the intercollision time
[C3] During coalescence, mass and momentum are pre
served.
Thus, r Q+ Q where Q+ and Q respectively
correspond to the quadratic integral operators associated
with creation and destruction of droplets due to coales
cence. These quadratic operators read (Hylkema and
Villedieu 1998; Hylkema 1999):
Q + i f (t, X, U (, *, U), V (, V*))x
2 J1u*,v*n[o,v]
f(t, x, u*, v*) (u u u* , v*)Jdv*du*
Q j(, x, u, v)f(t, x, u*, w) x
T(u u*, v, v*)dv*du*
(2)
where v(v, v*) vv* and u = vu** arethepre
collisional parameters, J is the Jacobian of the transform
(v, u) > (v u ) : J = (v /v)d with d the dimension
of the velocity phase space and B(u u*, v, v*) is the
collision/coalescence probability kernel which reads :
T(uu*, V, J*) = &(uu* ,v,*) iuu* B(v, *)
In this kernel, B(v, v*) = (r + r*)2 is the impact pa
rameter and = con.coal accounts for the collision
efficiency (col and the coalescence efficiency (coal. The
modelling of these efficiencies is detailed in the follow
ing two paragraphs.
Collision efficiency Collision efficiency Ecol is a
probability factor modelling the correlation of droplet
velocities immediately before a collision. It is equivalent
to consider that droplets may dodge each other due to
the gas flow surrounding them. Collision efficiency laws
thus require knowledge of local gas parameters such as
density p, or viscosity pg.
In the case of unbalanced droplet sizes, D'Herbigny
and Villedieu (2001) select two collision efficiency mod
els : the LangmuirBlodgett model, and the Beard
Grover model (Hylkema 1999; Achim 1999). These
laws strongly depend on the bigger droplet Reynolds
number and on a dimensionless number k which read,
when taking rl the smaller radius and r2 the bigger ra
dius :
Re 2P2 U U2 k pl u2 ul
sg 9'Pgr2
The number k is the ratio of the smaller droplet relax
ation time T, = 2pj to its residence time in the bigger
droplet influence zone T7 = .r
U2lU
In the case of low Reynolds numbers, Langmuir and
Blodgett (Langmuir 1948) numerically get the following
expression :
Er (k)
E (k)
0
31n(2k ) 2
+ 2(k 1.214))
if k < 1.214
otherwise
whereas for high Reynolds, they get:
SE (k)
E2 (k)
0 if k < 0.0833
(k+o5)2 otherwise
For intermediate cases they assume the following inter
polation :
FLB El(k) (Re/60))E (k)
Col(k, Re) 1 + Re/60 + 1 + Re/60
Beard and Grover (1974) suggest to increase the ac
curacy of formula (3) which results from a simple inter
polation between two limits. For this purpose, they use a
numerical solution of the incompressible NavierStokes
equations in order to determine the gaseous flow sur
rounding the bigger droplet depending on the Reynolds
number. They can therefore evaluate precisely the forces
on the smaller droplet and compute its trajectory. For
Re E [0,400] and ri < r2, it yields:
4
colj(k, Re) = [arctan(max(H(k, Re),0))] (4)
where
H = 0.1465 + 1.302Z 0.607Z2 + 0.293Z3
Z ln(k/ko)
ko= exp(0.1007 0.358 n(Re) + 0.0261[ln(Re)]2)
Typical collision efficiency values with these two laws
are in Fig. 1.
S05
04
S03
02
01
0 05 1 15 2 25 3
Velocity difference u2 ui (m/s)
Figure 1: Two collision laws selected for booster appli
cations with r2 150pm (red : rl = 2/m,
green: rl = 3pm; blue : rl 4pm)
Coalescence efficiency Coalescence efficiency
(coal is the probability for two droplets to merge after
a collision. Since the droplets can bounce on each other
or separate by reflexion or stretching if the remaining ki
netic energy of the new droplet is too high, coalescence
efficiency depends a priori on the velocity parameters,
the droplet material viscosity and surface tension, etc.
Most simple laws such as the BrazierSmith model
(Achim 1999; BrazierSmith et al. 1972) are commonly
used. Especially, the CEDRE code provides this model
7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010
but it is not used in the following simulations. There
fore, the BrazierSmith model is not detailed any further.
The formalism and the associated assumptions needed
to derive the Eulerian MultiFluid models are introduced
in Laurent and Massot (2001). We shall now recall the
two main steps which are the semikinetic derivation and
the sectional integration in order to precisely introduce
the coalescence terms.
1.2 Semikinetic model
Velocity coherence assumption In a first step we
reduce the size of the phase space, considering only the
moments of order zero and one in the velocity variable
at a given time, a given position and for a given droplet
size : n = fdu and u f ufdu/n which depends
on (t, x, S). In order to close the system, the following
assumptions are introduced:
[H1] For a given droplet size, at a given point (t, x),
there is only one characteristic averaged velocity
u(t, x, S).
[H2] The dispersion around the averaged velocity
u(t, x, S) is zero in each direction, whatever the
point (t, x, S).
It is equivalent to presume the following NDF :
f(t, x, u, S) = n(t, x, S)6(u u(t, x, S))
The monokinetic ]i p' ,lle',i [H1] is equivalent to reduc
ing the velocity distribution support to a one dimensional
submanifold parameterized by droplet size. It is fairly
correct when the Stokes number is small enough for the
droplets to follow the flow (Massot 2007). Yet it pre
vents dropletcrossing, which can occur in solid propul
sion during ejection from fast vortices or high speed
parietal injection.
The zerodispersion h] p, llwc'is [H2] is justified when
neither turbulence nor brownian phenomena interact
with the disperse phase.
These ]h.lp'lllw, thus induce strong limitations for
rocket booster simulation and will be discussed in the
conclusion of this paper.
Equations This step leads to a system of conser
vation equations called the semikinetic model; it is a
monokinetic model, saying that at a fixed position and
time and for a fixed size, there is only one possible ve
locity. The semikinetic system reads :
3tn + x, (nu) = cQ Q,
9t(nu) + 9x (nu (F u) = nF + Q'
Semikinetic coalescence operator In the semi
kinetic framework, the coalescence operator yields the
evolution rate of the zeroth and first order moments of
the velocity phase space. These two terms read, when
omitting the (t, x) dependency :
Q+ = n(o (,, *))n(V*)B3(o I*) do*
n 2 J v* [0,]
Qn n(v) n(o*)B3(, o*)I do*
,*e[o,+oo[
J2 v* [O ,+oo
Q.n(O) n(v*)B(vo, v*)Iu d*[
*C[ o,+ oo[
(6)
where I, I,, I and 1+ are the partial collisional inte
grals. To preserve the monokinetic and zerodispersion
assumptions, these integrals are computed with average
velocities and therefore read :
I+ =lu(v*) u( ) (I(V*) u( ), v, v*)
I =u(o) U() &(u() U(o), o, o)
Su(<) (w) Vu(o*) (u() u( ), o, o*)
fU(ou) U(W)1 CO(W) U(O), V, ow)
4U U(V) U(V) U(*) C(IUM U(V*), V, V*)
1.3 MultiFluid model
In a second step we choose a discretization 0 = So <
S1 < ... < SN for the droplet size phase space and
we average the system of conservation laws over each
fixed size interval [Sk 1, Sk [, called section. The set of
droplets in one section can be seen as a "fluid" for which
conservation equation are written, the sections exchang
ing mass and momentum. In order to close the system,
the following assumptions are introduced :
[H3] In one section, the characteristic averaged velocity
does not depend on the size of the droplets.
[H4] The form of n as a function of S is supposed to
be independent of t and x in a given section, thus
decoupling the evolution of the mass concentration
of droplets in a section from the repartition in terms
of sizes.
These assumptions are equivalent to presuming the NDF
in velocity and size inside each section k :
SCr [u(t, X ix, S)
VS [ 1, S [ n(t, x, S)
Uk(t,x)
Kk(t, x, S)
7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010
Since we have chosen a constant value for the velocity
distribution in each section, MultiFluid methods based
on assumption [H3] feature a first order velocity conver
gence with the number of sections. It induces a conve
nient integration simplification but it can generate tra
jectory modes. For coarse size discretizations indeed,
it assumes the same velocity for droplets of remarkably
different sizes in a given section.
As for [H4], it allows to reduce the size distribution
information of each section at (t, x) to a set of moments
of S, the number of which depends on the choice of the
(K k) k set of functions. The difference between the two
methods we are about to introduce is yet the choice of
their form.
In the socalled first order method, a constant IK
function in each section yields a first order granulome
try convergence with the number of sections. Whereas in
the second order method, a twocoefficient 2Kk function
yields a second order granulometry convergence. The
choice of refining the size distribution description is of
course related to the need of a fine resolution of polydis
persity, as told in introduction.
1.4 First order method
First order size assumption The first order Multi
Fluid method assumes, instead of [H4], the following
notation in each section k :
n(t, x, S) = k (t, x) k(S)
where mk, is the mass concentration of droplets in the
kth section, in such a way that
isk
*Sk1
k(S) Pi S/dS = 1.
6 v
Such an approach only focuses on one moment of the
distribution in the size variable : the moment in terms of
mass is chosen because it is conserved by coalescence.
Please note that the distribution on the last section is
a decreasing exponential with a fixed coefficient. This
choice allows the final section to treat the bigger droplets
but requires not to have a ,igiilik.iIl part of the mass.
This is a major limitation compared to the second order
method.
Equations The conservation equations for the kth
section result from the integration of the mass moment
of the semikinetic system (5) in each section k and
reads :
dtmuk+9x(mkUk) u=k
9t (m kUk )+ x (mk Uk 0gUk)
 M C
k
mkFk+Cmu
_1Cmu
(7)
(7)
vi
";,
vj+: \ 
0 vi Vi+01 v*
Figure 2: Coalescence terms are integrated on the
"(i, j, k) doii.nmiii
where Fe is the "classical" averaged drag force per unit
mass on a section (Greenberg et al. 1993; Laurent and
Massot 2001):
SU a 1 s 1 l (S)S3/2
F T T Sk 61/ I (S)
The total system thus counts twice as much equations as
the number of sections.
First order coalescence terms In each section
equation, the creation coalescence terms result from a
double integration : on the whole colliding partner size
space at the kinetic level (eq. 2) and on the concerned
section at the MultiFluid level. Yet the second depen
dency will not coincide with the section after mapping
the natural variables of the two precursor colliding part
ners. Thus splitting the integration domain thanks to
the section continuity yields elementary integrals Qijk,
triply indexed with the two precursor section numbers
(i, j) and the destination section number (k). These do
mains are illustrated in Fig. 2.
In the particular case of the first order method, con
sidering [H2] and assuming co,, = coal 1 (Laurent
and Massot 2001), the coalescence integrals Qijk take
the following form after factorizing the mass moments
Qi*JI J ( JI 2}+ r < 7,. ,,ddr*dr+
Q dv*<+v *[vV2+l4 3
f k= 1 K(r*)wlI (r )w(r*+ r )2 pir dr*dr
S +v*+v'C [vk,vk+l[
As for the disappearance terms, they can also be com
puted as sums of the elementary creation integrals and
must be so to ensure the conservation of matter and en
ergy. After some algebra, the coalescence terms 1C',
7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010
Ct 1C+ and 'C" read:
1(k) N
crn Zm Zm luj (1Q^I,'Q^I,)
i=1 j=1
N N
jl iI
1(k) N
1CrU+ i= i j *.j I (U dQ+ + IQk)
ik1 : 1
N
1Cku = meUfk E mj u,
j1
N
i=l
(8)
The integrands of the Qsjk integrals depend only on size
parameters r,r* and r' thanks to [H4]. This allows the
Qijk integrals to be precalculated as soon as the section
limits and the (IKsk) are given, i.e. once and for all at
the beginning of a simulation.
1.5 Second order method
Second order size assumption The second order
MultiFluid model is based on a twocoefficient expo
nential approximation of the size distribution in each
section. This means that [H4] reads, for S E [S k1, Sk [:
,((t, x, S) a (x, t)exp(bke(x, t)S)
where (a (x, t), bk (x, t)) I ensures
S kl 2kI(t, x, S)dS = nk(t, x)
Sk 2(t, S) S3/2dS = m(t, x)
S _k1 6
The choice of an exponential function ensures the pos
itivity of the distribution function whatever scheme is
chosen. It also aims at reducing the number of sections
and is well suited for evaporation, which requires mass
flux information at the section boundary. On the other
hand, coalescence becomes more difficult to compute if
one wants to respect realizability conditions on the sec
tions i.e. the fact that (ma, nk) couples are conditioned
by the section boundaries since they drive information
on average droplet volume.
Equations The conservation equations for the kth
section now read:
_2Cn
'Ck
C
Ot(mkUk)+9x. (mu uku )= m C+)U)C "+ 2C"T
(9)
The total system thus counts three times as much equa
tions as the number of sections.
dtnl+d,(ncul) 2Cn+
tnfk+9x(ncukc)=+C.
09mfe+89x (mUk) k
Second order coalescence terms Since the time
and space dependency of the sizedistribution functions
Kk(t, x, S) is no longer factorizable as was mk in the
first order method, the Qijk integrals must be computed
at each time step in each cell on its "(i, j, k) domain".
So it is numerically interesting to integrate the colli
sion/coalescence efficiency at the same time. That is
why this method is "less inappropriate" for efficiency
models. Let us take the notation :
i (t, x, r*r u u*) = (t, x, S*)l (t, x, S )x
(r* + r )2((r*, ru, u_ U* ) i Uj
The coalescence integrals have a different homogeneity
than in the first order method. They now include the
number, mass or momentum information and read :
Q ijt J 'x(t, r, u,
iQ Jj*+V[Vk,Vk+1
2Q* [4 x r[
J Jv*+v'>[vk,vk+1[
11 ij (t, x, r*, r ,
J Jv*+v'[ vk,vk+ [
u*l)dr*dr'
u*) ,. ddr*dr
U*) w4Tplr dr*dr
3
The coalescence source terms 2+', k2 2C(', ,
2C" +and 2C are still written as direct sums of the
"(i, j, ) iinici.iI, but therefore read:
I(k) N N N
k 11 2. ZZk, ki kji
i lj 1 j i= 1
I(k) N N N
2CrE EZQ 2Q ZZ EE1(Q,+ t Q<;)
i=i"=1 j =i 1
i KH
2cmu *Zu 2, 2cm
kc 1 ZZ (Uj2QIti~uUk6,2CTUk k k
j=l kKL
(10)
1.6 Numerical strategy for coalescence
Towards computation, the main difference with the first
order method is the time and space dependency of the in
tegrals ^2 >,, 2Q1 and Qij. The second order method
requires indeed a costly numerical integration for each
time step and space cell. We have chosen to per
form these integration with a 2D 5 point NewtonCotes
quadrature which requires i3 (t, x, r*, r u u*) to
be computed 52 x N,2 x Ne,1 times per time step. The
accuracy of this method on such exponential functions
has been tested but will not be discussed any further.
Yet performing this integration with an advanced col
lision or coalescence efficiency model is now natural and
hardly more costly. That's why these efficiencies have
only been included in the second order method. Laurent
et al. (2004) introduces a general form for coalescence
7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010
efficiency formulas to be suitable for a precalculated
first order method but no collision efficiency models can
be easily implemented since they require local gas pa
rameters. That is why the first order method here ex
cludes such collision efficiency models.
2 Dynamic Study : the Nozzle Test Case
As a main criterion, we want to compare the two Multi
Fluid methods on their ability to describe the dynam
ics of a coalescing cloud. We therefore need a well
suited test case, inducing coalescence as well as size
conditioned dynamics and difficult enough to highlight
the limitations of the methods.
For that purpose, we have chosen the configuration
which was used in Laurent and Massot (2001) to vali
date the first order method compared to a reference La
grangian solver : a 2D axisymmetrical conical decel
erating nozzle, designed in such a way that it admits,
for oneway coupling spray dynamics a selfsimilar so
lution. We though inject a lognormal distribution instead
of a bimodal one.
2.1 Definition of configuration
The chosen configuration is stationary 2D axisymmetri
cal in space and 1D in droplet size. It is described in
detail in Laurent et al. (2004). Hence, only its essential
characteristics are given here.
For the problem to be onedimensional in space, con
ditions for straight trajectories are used and are compat
ible with the assumption of an incompressible gas flow.
This leads to the following expression for the gaseous
axial velocity U, and the reduced radial velocity U,/r,
for z _> z :
U U(z)
S1(Zo) U, U,
.21(zo)
r Z Z
where zo > 0 is the coordinate of the nozzle entrance
and the axial velocity U(zo) at the entrance is fixed.
The trajectories of the droplets are also assumed straight
since their injection velocity is colinear to the one of the
gas. This assumption is only valid when no coalescence
occurs. However, even in the case of coalescence, it is
valid in the neighborhood of the centerline.
Because of the deceleration of the gas flow in the con
ical nozzle, droplets are slowing down too, however at a
rate depending on their size and inertia. This will induce
coalescence. The deceleration at the entrance of the noz
zle is taken as a(zo) = 2U(zo)/zo; it is chosen large
enough so that the velocity difference developed by the
various sizes of droplets is important. Laurent and Mas
sot (2001) chose rather large values, as well as strong de
celeration, leading to extreme cases: U(zo) 5 m.s1,
7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010
zo = 5 cm : these values generate a very strong cou
pling between coalescence and droplet dynamics. These
severe conditions make the test cases under considera
tion very efficient tools for the numerical evaluation of
the two Eulerian models. Finally please note that no
efficiency law is used in this section so that the sec
ond order method remains comparable to the first order
method, which lacks this sort of model. We thus assume
Ecoi coal 1 in this section.
2.2 Droplet initial distribution
The lognormal distribution is useful to characterize alu
mina particles in solid propergol rocket boosters. We
therefore choose a lognormal distribution on the surface
variable, without loss of generality. The droplets are
constituted of liquid alumina, their initial velocity is the
one of the gas, their initial temperature, fixed at the equi
librium temperature 3600 K (corresponding to an infinite
conductivity model), does not change along the trajecto
ries.
For the sake of the first order method accuracy, we
do not want to transfer too much mass in the last sec
tion. The initial injected mass concentration is therefore
taken as mo 1.0609 kg.m 3. For the same reason, the
lognormal parameters are set to SLN 1600pm/2 and
8LN = 1.5 which corresponds to a sharp distribution
centered on a radius of 11.3pm. We shall yet keep in
mind that the second order method can accurately treat
cases where the final section hosts a ,igniliik.i portion
of the mass whereas the first order method would need to
increase the number of sections. The nozzle test case has
thus been designed to meet optimal first order specifica
tions i.e. where as much sections as possible get filled
but the final section remains almost empty.
We now choose the tested numbers of sections : a 5
section test case illustrates what happens with a coarse
discretization. A 13 section and a 25 section test cases
prove the convergence of the method. They are com
pared to a 53 first order run which we use as a reference
solution. The first order indeed has been fully validated
compared to a lagrangian test case in Laurent and Mas
sot (2001).
2.3 Distribution spread up results
The processing of polydispersion gives a first indica
tion on the method's accuracy. We therefore compare
the mass concentration distributions at the nozzle output
to the 53 section reference test case in Fig. 3. With 5
sections, the first order strongly overestimates the size
growth while the second order underestimates it. With
13 sections, the trend is the same but the error is smaller.
Finally for 25 sections, we consider both methods are
1 5 sections
002
a o
.o 008
.0
00 006
0
" 004
U 002
0
 008
O
I 13 sections 
I 25 sections I
0 10 20 30 40 50 60
Radius (pm)
Figure 3: Mass concentration distribution at the noz
zle's end computed with the first and second
order methods (empty and filled symbols) and
reference (black)
roughly converged. The growth overestimation of the
first order method for poor size discretization brings
about major consequences for the spray dynamics.
2.4 Spray dynamics results
Let us now appreciate the consequences of polydisper
sion treatment on dynamics. Fig. 4 shows the evolution
of mass and number concentrations along the center
line of the nozzle computed with both Eulerian methods
(5,13,25 sections) compared to the reference solutions
(53 section first order method and lagrangian reference
test case). For the first order case, number concentra
tions are computed considering section average droplet
volumes given by integration of the (1k)k functions :
data are therefore redundant with mass concentration
data. We emphasize the fact that the convergence of the
7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010
cl
o
.q
o 02
o
C
C
E,,
0
Z
1 05
C
5Y 1
0
c
0
1 05
:
C
 09
C
0 85 L
005
01 015
Position (m)
Figure 4: Number and mass concentration along nozzle
for 5, 13, 25 sections and reference solution
(black)
Groups (im) 5 sections 13 sections 25 sections
Gi=[0, 12.5] 1 1 to 3 1 to 6
G2[12.5, 25] 2 4 to 6 7 to 12
G3[25, 37.5] 3 7 to 9 13 to 18
G4=[37.5, 50] 4 10 to 12 19 to 24
Gs (> 50p/m) 5 13 25
Table 1: Composition of the five section groupings
two methods is achieved for 25 section simulations as
we can see in Fig. 5, which is a zoom of Fig. 4.
To compare precisely the effect of polydispersion on
dynamics, let us now consider five size intervals : they
correspond to the sections in the 5 section case and are
composed of section groupings in the other cases as il
lustrated in table 1. The evolution of the mass concen
tration of these groupings along the nozzle is given in
Fig. 6. It is there obvious that the 5 section first order
error on coalescence is severe, especially in the fifth and
last grouping G5 where no mass should be found conse
quently to the problem parametrization. Finally, the sec
ond order method reveals to have acceptable error with
as few as 5 sections.
016 018 02 022 024
Position (m)
Figure 5: Zoom on the convergence of mass concen
tration at the nozzle's end (black : rough la
grangian validation)
5 sections
13 sections
25 sections
53 sections
I 1st order MF 2nd order MF
Is 10s
1.5 s 50 s
2s 180 s
5 s n.a.
02 025 Table 2: Computational time for 10, 000 iterations on a
2.66GHz Intel Core 2 Duo CPU
2.5 Computational time
To conclude this study, table 2 recaps the duration of the
different runs. It is obvious that the first order method
is much faster since it is highly optimized thanks to the
precalculation. Even when considering the 5 section
second order run which is fairly as accurate as the 13
section first order it is 6 time slower.
Yet two restrictions must be added about the use of the
first order with a thinner size space discretization. First,
it is not possible to meet such a high computation speed
with specific collision/coalescence laws. Second for typ
ical solid propulsion simulations, each section requires
the resolution of an unsteady 3D pressureless Eulertype
system so that solving for coalescence is no longer the
limiting step. Thus increasing the number of section be
comes a costly operation, restraining the profit of the
first order method.
3 Coalescence Experimentation : the
D'Herbigny Test Case
In this section, we provide further validation for the sec
ond order MultiFluid method, by comparing it to a 1D
experiment on collision/coalescence efficiency models.
Moreover this configuration provides tough testing for
Eulerian models because it deals with bimodal size dis
7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010
Empty First order
SFilled Second order
0 / A 5 sections
S 13 sections
S25 sections
005 / Reference
Figure 6: Evolution of mass concentration in the four
finite section groupings and in the last section
for 5, 13 and 25 sections
tributions. Indeed, modal distributions (i.e. monodis
perse or discrete distributions) are rather suited for the
Lagrangian point of view while they introduce stiff size
distributions in Eulerian modelling.
3.1 The D'Herbigny collision efficiency
experiment
In the D'Herbigny experiment, coalescence is studied
through the growth of a bigger droplet falling through
a fog of smaller droplets. Details about the experimental
device and conditions can be found in D'Herbigny and
Villedieu (2001) and in Fig. 7.
In this experiment, the average injection radius of the
bigger droplets is r2 150pm thanks to piezoelectric
injector which is fairly accurate. The smaller droplets,
composing the fog, have a radius ri e [2/m, 4pm].
The two types of droplets have an approximately con
stant velocity difference of 3m/s along the tunnel. One
can therefore compute the values of the k parameter used
in the LangmuirBlodgett and BeardGrover collision ef
ficiency models. They are given in table 3.
Moreover, the collisional Reynolds number Re is
around 60 which is comparable to the values in the Ari
ane 5 P230 booster. This range of values (moderate Re
and low k) is supposed to be favorable to the Beard
Grover law.
Figure 7: D'Herbigny experimental device
rl 2pm 2, 5pm 3pm 3, 5pm 4pm
k 0.78 1.22 1.75 2.39 3.12
Table 3: Collision parameter k values in the
D'Herbigny experiment
3.2 Results
Droplet radius growth study D'Herbigny and
Villedieu (2001) provide an analytic formula to esti
mate the Sauter mean radius growth Ar2 of the bigger
droplets when assuming a small growth and a constant
velocity difference :
Ar2 = Q&collCv
where C, is the volume fraction of the fog, coll the av
erage collision efficiency and a results from the integra
tion of the impact parameter. As a first step, we thus
perform the simulation with a constant velocity of 3m/s
to fit the analytical result. Experimental, analytical and
simulation results are given in Fig. 8.
The first conclusion is that collision efficiency laws
have a strong impact on radius growth, here undermining
it up to 50%.
The second conclusion is a positive validation of
both laws implemented in the second order Multi
Fluid Method. According to D'Herbigny and Villedieu
(2001), this constant velocity study also concludes
that LangmuirBlodgett law is most suited for rocket
chamberlike conditions.
Method numerical stability The third conclusion
of these simulations is that the method allows robust and
accurate simulation of stiff distributions. Indeed the sim
ulations are converged for a rough number of 100 itera
tions and provide a selfsimilar solution which means
low diffusion in the size phase space.
Monodisperse injector
0
0
0
Transparent
o  tunnel (5m)
o
Granulometry
measurement
 Fog machine (r, Ci, ,
Tranquilizing Grid
a
7 i~8a
30
n LangmuirBlodgett
25 O BeardGrover
20
. 15
0 10 20 30 40
Fog concentration C, (ppm)
50 60
Figure 8: Radius growth Ar2 after 5m through a C,
concentration fog (empty red : Analytical;
filled green : Second order simulation with
constant velocity; black : Experimental).
3.3 Discussion on velocity
As a second step, let us revise the constant velocity as
sumption. On the one hand, droplets do reach a terminal
velocity, i.e. the drag force compounds the weight but
this limit velocity increases with the droplet size since
the drag force increases like the droplet surface while
the weight increases like the droplet volume. On the
other hand, collisions with a static fog induce a momen
tum dilution i.e. a significant slow down. We can see the
evolution of the effective velocity as calculated by the
code when gravity, drag force and coalescence momen
tum transfer are enabled in Fig. 9. The droplet velocity
tends to be slightly smaller than the terminal velocity
because of momentum dilution.
2 25
2
15
a
0 1 2 3 4 5
Position (m)
Figure 9: Velocity along the tunnel with gravity, drag
force and coalescence momentum transfer
(blue triangles) and terminal velocity for the
corresponding Sauter Mean Radius (black)
Yet the consequences of velocity on coalescence are
not obvious. First let us neglect efficiency laws. The co
alescence rate per unit time linearly increases with the
7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010
velocity difference (eq. 2) but linearly decreases with
the time spent in the tunnel so that it vanishes in the fi
nal growth formula (eq. 11). In fact radius growth for a
fixed length travel only depends on the number of colli
sions i.e. fog concentration and travel length, assuming
a fairly constant cross section. However when consider
ing efficiency laws, the velocity difference plays a role
since it controls the way smaller droplets can dodge the
bigger ones. For instance eq. 11 depends linearly on an
averaged collision efficiency factor strongly conditioned
by velocity difference as shown in Fig. 1.
The new simulation results with full dynamics is
given in Fig. 10. It shows that the Beard Grover
efficiency correlation is not so bad while the aver
age velocity might have been slightly overestimated in
D'Herbigny and Villedieu (2001). As a final restric
tion, please note that no coalescence efficiency law has
been used (coal 1) which suggests not to conclude on
which collision efficiency law to use.
30
o LangmuirBlodgett
25 BeardGrover
S 10 20 30 40
Fog concentration C, (ppm)
50 60
Figure 10: Radius growth Ar2 with velocity resolu
tion (filled blue : Second order simulation;
black: Experimental)
3.4 Conclusion on the model and methods
The previous results confirm that Eulerian MultiFluid
methods can be used to simulate accurately the size
distribution evolution of a coalescing spray and its size
conditioned dynamics. The first order method provides,
with a reasonable number of sections (between 10 and
20), extremely fast results since it is highly optimized.
The second order MultiFluid method provides good re
sults with a very coarse sizespace discretization (as few
as 5 sections). If both methods can include coalescence
efficiency models, only the second order method easily
provides collision efficiency models. The correspond
ing laws can significantly reduce the particle growth rate
which is experimentally observed so that these efficiency
models are essential to capture the physics of coales
cence.
t
Y sL~~
=I~C1
7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010
4 Rocket booster simulation
In this section, we conclude the study by achieving a
2D complex simulation of a rocket booster chamber and
nozzle which experiences the meeting of a parietal in
jection of gas and small alumina droplets resulting from
the combustion of propergol and an axial twophase flow
carrying bigger droplets supposedly coming from pri
orly injected, burnt and coalesced particles. The upper
part of the chamber generating these bigger droplets is
not solved.
4.1 The CEDRE code
The CEDRE code has been developed at ONERA for
several years to provide a multiphysics comprehensive
tool for energetic problems in aeronautics. The main
feature is the combination of a non structured spatial
solver with models as varied as multiphase flows, multi
species chemistry, thermal conduction, radiation, wall
film models, etc. The coupling can be one or two
way. As for the disperse phase, CEDRE includes a La
grangian and two Eulerian methods i.e. a multiclass and
the second order MultiFluid method studied here. Only
the Lagrangian and the MultiFluid methods provide co
alescence models.
4.2 Configuration
The 2D rocket booster simulation is performed on a sim
plified configuration, yet featuring the main difficulties
of solid propulsion typical flows i.e. parietal injection
and a supersonic nozzle. It takes place on a deformed
structured 1500 cell mesh (Fig. 11, top). When seeking
a twophase stationary flow solution, one usually takes
a converged gas flow field as an initial state. With a to
tal rate of flow of 10 kg.s 1.m 2 from the wall (pro
pellant combustion) and the head end (upstream flow),
this gas flow involves extremely high velocity gradients
in the nozzle (Fig. 11, middle), which will induce size
conditioned droplet dynamics as studied in section two.
The simulation strictly speaking starts when injecting
the disperse phase : a monomodal 5om radius wall in
jection represents the droplets resulting from recently
burnt aluminum particles, directly expelled from the pro
pellant grain, and another roughly monomodal distri
bution around 20/m is injected on the axis to model
the previously coalesced droplets (Fabignon et al. 2003).
For the purpose of our simulation, the injected volume
fractions (more or less 5.10 5) approximate correctly
typical rocket booster conditions and preserve the dilute
spray p 1 ,liL'wNi
wall injection (5p/m)
I I I I I I I I I I
S2000
1181
0697
0412
0243
0143
0085
0050
S2000
1181
0697
0412
0243
0143
0085
0050
Figure 11: Top : Deformedstructured 1500 cell mesh (Ar
rows : injection zones); Middle : Gaseous Mach
number without particle; Bottom : with particles
4.3 Results
For the disperse phase, we choose a 5 section discretiza
tion (0, 12.5, 25, 37.5 and 50 pm) as in section two. A
10 2 s time interval is required to allow the first droplets
to reach the end of the nozzle. With a 10 6 s time step
(10, 000 iterations), we perform a 1 h single processor
AIX platform computation. The stationary volume frac
tions for the five sections are displayed in Fig. 12.
First, the disperse phase has an impact on the gas flow
since we perform a twoway coupling simulation. Spe
cific impulse loss can indeed be observed with the nozzle
Mach number decrease in Fig. 11, bottom. Second, we
do note bigger particle creation in Fig. 12 : coalescence
occurs as soon as the two injected types of droplets meet
since the third section gets filled with ignilik.lI vol
ume fractions ( 10 5) in the chamber. Third, bigger
droplets are created in the nozzle eventhough accelera
tion there induces a strong volume fraction dilution. In
pursuance of section two, section 5 is left almost empty
with a volume fraction around 10 10 (the colormap is
different from the other section ones); yet a different
discretization could be used to activate this section thus
increasing size accuracy.
4.4 Remaining Numerical difficulties in solid
propulsion
The final simulation of this paper illustrates the ability
of MultiFluid models to provide in a reasonable time
fairly accurate information on polydisperse coalescing
sprays in complex geometry gas flows. The above re
sults however suggest several possible improvements,
among which reducing computational time and consid
ering more physical phenomena.
Typical advanced rocket booster simulations are un
steady and performed on more complex meshes : they
can include from 10, 000 cells to millions of cells for
I Soctinn 1 [0 19 lii m I
Figure 12: Dispersed phase volume fractions
3D meshes. Highly optimized code and multiprocessor
computations (with several spatial domains for instance)
are therefore required. In our constant time step sim
ulation, the CFL condition is driven by the flow in the
nozzle but the gap between the chamber and the nozzle
stay times makes most of the early time steps too short
and therefore uselessly numerous. An adaptative time
step method considering both phases including operator
splitting could prevent from illsuited time steps.
Physically, high velocity gradients in the nozzle in
duce droplet secondary breakup which strongly under
mines droplet growth. Since it is out of the scope of this
paper, this phenomenon has been neglected but its mod
elling reaches two requirements : specific impulse com
putation is subject to the precise knowledge of droplet
size after the nozzle throat and experimental data on alu
mina size distribution mainly consist of ejected material
which has encountered complete breakup process. De
veloppement of fragmentation models for the Eulerian
MultiFluid method has been achieved in Dufour (2005).
Finally, the main drawback of such methods comes
from the monokinetic ]ihp ,liLsis : it forbids droplet
crossings which should occur on the centerline of our
simulation and induces instead a droplet accumulation
artifact due to the momentum averaging on the sym
metry axis. On more advanced solid propulsion un
steady configurations, dedicated to hydrodynamic insta
bility studies for instance, highinertia droplets ejected
from vortices need a multivelocity treatment which ad
vanced high order moment methods can provide as de
veloped in Kah et al. (2010).
1 0x10
S37x10
1 4x10'
5 2x10
1 9x10
7 2x10(
2 7x10
10x10O
1 0x10
S37x10'
1 4x10
5 2x10
1 9x10
7 2x10
1 4x10
52x10(
7 19x10
72x10'
27x10
10x lO
10x10
37x10(
14x10(
522x10(
1 9x10
7 2x10
27x10
10xlO0
10x10
73x10
52x10
38x10
27x10
19x10
14x10
10xlO
SSection 2: [12.5, 3]i .
Section 3: [25, 37.5]/um
7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010
5 Conclusions
In this paper, we provide a comprehensive validation of
a Eulerian, second order in size, method for solving two
phase polydisperse coalescing flows. The conclusion of
this study is that the existing method is able to simulate
accurately the dynamics of such flows with less sections
than the first order method as demonstrated in section
two yet in a slower way. It is however the only Eule
rian MultiFluid method including validated advanced
collision efficiency models, which are crucial for rocket
booster studies. Finally, the method is practically effi
cient in rocket booster simulation contexts; moreover re
ducing the number of required sections can be profitable
to the computational time of such simulations. The re
maining limitations underlined in the final section re
quire though further developments but the method can
bear such improvements.
6 Acknowledgments
The present research was done thanks to a Ph. D. Grant
from DGA, Ministry of Defence (M. S. Amiet, Techni
cal Monitor). The authors also thank G. Dufour and T.
Fontfreyde for their preliminary contribution to the sec
ond order method in the CEDRE code.
References
B. Abramzon and W.A. Sirignano. Droplet vaporization
model for spray combustion calculations. Int. J. Heat
Mass Transfer, 32:16051618, 1989.
P. Achim. Simulation de collisions, coalescence et rup
ture de gouttes par une approche lagrangienne: appli
cation aux moteurs a propergol solide. PhD thesis, Uni
versit6 de Rouen, 1999.
A. A. Amsden, P. J. O'Rourke, and T. D. Butler. Kiva ii,
a computer program for chemically reactive flows with
sprays. Technical Report LA11560MS, Los Alamos
National Laboratory, 1989.
K. V. Beard and S. N. Grover. Numerical collision ef
ficiencies for small raindrops colliding with micron size
particles. J. of atmospheric sciences, 31:543550, 1974.
G. A. Bird. Molecular gas dynamics and the direct sim
ulation of gas flows. Oxford Science Publications, 42,
1994.
F Bouchut. On zero pressure gas dynamics. In Advances
in kinetic theory and computing, pages 171190. World
Sci. Publishing, River Edge, NJ, 1994.
I Q rtinn , f,,m I
PR. BrazierSmith, S.G. Jennings, and J. Latham. The
interaction falling water drops : coalescence. Proceed
ings of the Royal Society, 326:393408, 1972.
S. De Chaisemartin. Moddles euleriens et simulation
numerique de la dispersion turbulente de brouillards qui
s'evaporent. PhD thesis, Ecole Centrale de Paris, 03
2009. http://tel.archivesouvertes.fr/tel00443982/en/.
F. X. D'Herbigny and P. Villedieu. Etude exp6rimentale
et num6rique pour la validation d'un module de coales
cence. Technical Report RF1/05166 DMAE, ONERA,
2001.
G. Dufour. Modelisation multifluide eulerienne pour les
ecoulements diphasiques a inclusions dispersees. PhD
thesis, Universit6 Paul Sabatier Toulouse III, 2005.
J. K. Dukowicz. A particlefluid numerical model for
liquid sprays. J. Comput. Phys., 35(2):229253, 1980.
J. Dupays, Y. Fabignon, P. Villedieu, G. Lavergne, and
J. L. Estivalezes. Some aspects of twophase flows
in solidpropellant rocket motors. In Solid Propellant
Chemistry, Combustion, and Motor Interior Ballistics,
volume 185 of Progress in Astronautics and Aeronau
tics. AIAA, 2000.
Y. Fabignon, O. Orlandi, J.F. Trubert, D. Lambert, and
J. Dupays. Combustion of aluminum particles in solid
rocket motors. AIAA Paper 20034807, July 2023 2003.
In 39th AIAA/ASME/SAE/ASEE Joint Propulsion Con
ference and Exhibit, Huntsville, Tx.
J. B. Greenberg, I. Silverman, and Y. Tambour. On the
origin of spray sectional conservation equations. Com
bustion and Flame, 93, 1993.
J. J. Hylkema. Modelisation cinetique et simulation
numerique d'un brouillard dense de gouttelettes. Appli
cation auxpropulseurs a poudre. PhD thesis, Ecole Nat.
Sup6rieure de l'A6ronautique et de l'Espace, 1999.
J. J. Hylkema and P. Villedieu. A random particle
method to simulate coalescence phenomena in dense liq
uid sprays. In Lecture Notes in Physics, volume 515,
pages 488493, Arcachon, France, 1998. Proc. 16th Int.
Conf. on Num. Meth. in Fluid Dyn.
D Kah, F Laurent, L Freret, S De Chaisemartin, R
Fox, J Reveillon, and M Massot. Eulerian quadrature
based moment models for dilute polydisperse evap
orating sprays. Flow Turbulence and Combustion,
pages 126, 04 2010. http://hal.archivesouvertes.fr/hal
00449866/en/.
P. Kuentzmann. Aerothermochimie des suspensions.
GautierVillars Editeur, 1973.
7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010
I. Langmuir. The production of rain by a chain reaction
in cumulous clouds at temperatures above freezing. J.
Meteor, 5:17192, 1948.
F. Laurent. Numerical analysis of Eulerian multifluid
models in the context of kinetic formulations for dilute
evaporating sprays. M2AN Math. Model. Numer. Anal.,
40(3):431468, 2006.
F. Laurent and M. Massot. Multifluid modeling of lam
inar polydispersed spray flames: origin, assumptions
and comparison of the sectional and sampling methods.
Comb. Theory and Modelling, 5:537572, 2001.
F. Laurent, M. Massot, and P. Villedieu. Eulerian multi
fluid modeling for the numerical simulation of coales
cence in polydisperse dense liquid sprays. J. Comp.
Phys., 194:505543, 2004.
M. Massot. Eulerian Multifluid modelsforpolydisperse
evaporating sprays, chapter III of "MultiPhase React
ing Flows: Modelling and Simulation", pages 79123.
CISM Courses and Lectures No. 492. Springer, 2007.
Editors D.L. Marchisio and R. O. Fox,.
M. Massot, F. Laurent, D. Kah, and S. De Chaise
martin. A robust moment method for evaluation of the
disappearance rate of evaporating sprays. to appear
in SIAM J. of App. Math., 2010. http://hal.archives
ouvertes.fr/hal00332423/en/.
P. O'Rourke. Collective drop ifh. i on vaporizing liquid
sprays. PhD thesis, Los Alamos National Laboratory
87545, University of Princeton, 1981.
M. Rilger, S. Hohmann, M. Sommerfeld, and
G. Kohnen. Euler/Lagrange calculations of turbulent
sprays: the effect of droplet collisions and coalescence.
Atomization and Sprays, 10(1):4781, 2000.
B. T6th. Twophase flow investigation in a coldgas
solid rocket motor model through the study of the slag
accumulation process. PhD thesis, Universit6 Libre de
Bruxelles, 2008.
I. M. Vasenin, R. K. Narimanov, A. A. Glazunov, N. E.
Kuvshinov, and V. A. Ivanov. Twophase flows in the
nozzles of solid rocket motors. J. Propulsion and Power,
11(4):583592, 1995.
F A. Williams. Spray combustion and atomization.
Phys. Fluids, 1:541545, 1958.
F A. Williams. Combustion Theory (Combustion Sci
ence and Engineering Series). ed F A Williams (Read
ing, MA: AddisonWesley), 1985.
Y. B. Zel'dovich. Gravitational instability : an approxi
mate theory for large density perturbations. Astronomy
and Astrophysics, 5:8489, 1970.
