Prévia do material em texto
Chemical Engineering Science 59 (2004) 3477–3490
www.elsevier.com/locate/ces
PDF modelling of turbulent non-premixed combustion
with detailed chemistry
Haifeng Wang∗, Yiliang Chen
Department of Thermal Science and Energy Engineering, University of Science and Technology of China, Heifei,
Anhui 230027, People’s Republic of China
Abstract
Although the signi/cant advantage for the probability density function (PDF) methods of the exact treatment of chemical reactions in
turbulent combustion problems, a detailed chemistry mechanism (e.g., the GRI mechanism) has not been implemented in the practical
calculations by now due to the prohibitive computation of PDF methods. In this work, a detailed mechanism (GRI-Mech 3.0, consisting
of 53 species and 325 elemental reactions) is /rstly incorporated into the PDF calculation of a turbulent non-premixed jet 6ame (Sandia
Flame D). The 6ow is formulated in the boundary layer form. The joint composition PDF closure level is applied and a multiple-time-scale
(MTS) k–� turbulence model is combined for the closure of turbulent transport terms. The molecular mixing process is modelled by
the Euclidean minimum spanning tree (EMST) mixing model. The solutions are obtained by using the space marching algorithm for
turbulence equations and node-based Monte Carlo particle method for PDF evolution equation. The chemical reaction source terms are
integrated directly. Extensive comparisons between the predictions and the measurements are made, which involve radial pro/les of mean
and rms (root mean square), conditional mean, scatter plots of scalars and conditional PDF distribution etc. The 6ame structures are well
represented by the present calculation, including intermediate species (e.g. CO and H2) mass fractions, pollutant NO emission and local
extinction.
? 2004 Elsevier Ltd. All rights reserved.
Keywords: Turbulent non-premixed combustion; PDF modelling; Detailed chemistry; Multiple-time-scale k–� model; EMST mixing model
1. Introduction
Turbulent non-premixed combustion takes place widely
in chemical engineering and industrial applications. A
key issue arising from combustion problems is the strong
interactions between turbulent 6uctuations and chemical
reactions, which make the failure of conventional moment
closure to calculate turbulent combustion problems. During
the past two decades, several closure methods have been ad-
vanced to treat such interactions, e.g. the probability density
function (PDF) methods (Pope, 1985, 1994; Dopazo, 1994),
steady and unsteady 6amelet models (Peters, 1984; Pitsch
et al., 1998), conditional moment closure (CMC)
(Klimenko, 1990; Bilger, 1993) and one dimensional tur-
bulence (ODT) model (Hewson and Kerstein, 2002) etc.
In the aspect of experimental investigations, several
kinds of burners have been designed and measured in TNF
∗ Corresponding author. Tel.: +86-551-3601650;
fax:+86-551-3606459.
E-mail address: wanghf@mail.ustc.edu.cn (H. Wang).
workshop (1996), e.g. piloted jet 6ames, bluF-body stabi-
lized jet 6ames and swirling jet 6ames, to investigate the
in6uence of turbulence on 6ame structures. These 6ames
have become standard test cases for the combustion mod-
els and have been simulated extensively (see proceedings
of TNF workshop, 1996). Piloted jet 6ames are geometri-
cally simplest, while it still re6ects strong interactions be-
tween turbulence and /nite rate chemistry. Therefore, they
are preferred for the validation of combustion models. One
of the piloted jet 6ames (Sandia 6ame D, Barlow and Frank
(1998), experimental data sets are available from the web-
site of TNF workshop, 1996) is mainly concerned with in
the present work.
The piloted jet 6ame (6ame D) has been extensively sim-
ulated, including 6amelet simulations (Pitsch et al., 1998;
Coelho and Peters, 2001; Pitsch, 2002), CMC simulations
(Roomina and Bilger, 2001) and PDF simulations (Lindstedt
et al., 2000; Tang et al., 2000; Xiao et al., 2000; Xu and
Pope, 2000) etc (more simulations can also be found in the
proceedings of TNF workshop, 1996). Finite rate chemistry
can all be handled by above three models, while detailed
0009-2509/$ - see front matter ? 2004 Elsevier Ltd. All rights reserved.
doi:10.1016/j.ces.2004.05.015
mailto:wanghf@mail.ustc.edu.cn
3478 H. Wang, Y. Chen / Chemical Engineering Science 59 (2004) 3477–3490
mechanisms are only easily incorporated into the 6amelet
models and CMC. For PDF methods, it is not tractable to
perform arbitrary detailed mechanism because of its heavy
computational cost. Reducing computational time is a main
research task of PDF methods. Usually reduced mecha-
nisms, e.g. for methane oxidation, reduced four-step mecha-
nism (Peters and Kee, 1987; Bilger et al., 1990), are applied
in the PDF calculations (Chen et al., 1989; Jones and Kakhi,
1997, 1998). Saxena and Pope (1999) incorporate a skele-
tal mechanism in the PDF methods by using the in situ
adaptive tabulation (ISAT) algorithm (Pope, 1997) which
eKciently decreases the computational time. In Sung et al.
(1998), it is recognized that, with the increasing compu-
tational capacity, conventional reduced mechanisms (four-
or /ve-step) are not necessary, and eForts to develop aug-
mented reduced mechanisms (ARMs) are worthwhile to
be made. An ARM for methane oxidation consisting of 16
species and 12-step reactions is generated by Sung et al.
(1998) from GRI-Mech 1.2 (GRI Website, 1995) and
ARMs containing nitrogen chemistry are then developed
based on GRI-Mech 3.0 (GRI Website, 1995) by Sung et al.
(2001). Xu and Pope (2000) and Tang et al. (2000) incor-
porate ARM in their PDF calculations of 6ame D by using
ISAT algorithm (Pope, 1997; Saxena and Pope, 1999). The
6ame structures are well predicted by their calculations, in-
cluding the pro/le of CO in the fuel-rich region. Usually the
CO concentration in the fuel-rich region is under-predicted
by conventional reduced mechanisms (Chen et al., 1989)
and skeletal mechanism (Saxena and Pope, 1999), which
shows the superiority of ARMs. Similar improvement is
also reported by Lindstedt et al. (2000), where a compre-
hensive reduced mechanism containing 16 independent, 4
dependent and 28 steady-state species is adopted.
Although the progress of reduced mechanisms and im-
proved predictions of PDF methods with developed reduced
mechanisms, a detailed mechanism has not been imple-
mented in the practical PDF calculations and the examina-
tion of detailed mechanisms performance in the PDF simula-
tions is desired. Chen (Proceedings of third TNF workshop,
1996) simulated 6ame D using joint scalar PDFmethod with
a detailed mechanism (GRI-Mech 1.2) on a parallel cluster,
while poor results were obtained, which may be attributed
to the small number of particles (50 particles/cell) used for
the Monte Carlo simulation (Pope, 1981). In the present
work, we concentrate on the implementation of a detailed
mechanism (GRI-Mech 3.0) in the PDF modelling of 6ame
D within acceptable computational time.
In the following sections, the experimental conditions for
6ame D are outlined /rstly. The modelling methods are
described subsequently, including turbulence model, scalar
joint PDF evolution equation and small-scale mixing model
etc. Then the numerical solution methods and speci/cations
of boundary conditions are discussed. The numerical re-
sults are presented and compared with experimental data
in the next section, and conclusions are drawn in the /nal
section.
D
F
ue
l j
et
P
ilo
te
d
fla
m
e
A
ir
co
flo
w
DP
Fig. 1. Burner geometry of Sandia 6ame D.
2. Sandia Flame D
The piloted methane/air turbulent non-premixed jet 6ame
(Sandia 6ame D of Barlow and Frank, 1998) is chosen as
numerical test case. The geometry of the burner is illustrated
in Fig. 1. The jet 6ows consist of three parts, the fuel jet, the
piloted 6ame and the air co-6ow. The inner diameter of the
fuel jet nozzle is D=7:2×10−3 m and the outer diameter of
the annular piloted 6ame is DP=18:4×10−3 m. The fuel is
a mixture of air and methanewith the ratio 3:1 by volume.
The temperature of the fuel is 294 K, and the bulk velocity
of the fuel jet equals 49:6 m=s (±2 m=s). The annular piloted
6ame is a lean mixture (� = 0:77) of C2H2, H2, air, CO2
and N2 with the same nominal enthalpy and equilibrium
composition as methane/air at the same equivalence ratio.
The bulk velocity of the pilot is 11:4 m=s (±0:5 m=s). The
air co-6ow temperature is 291 K, and the velocity is 0:9 m=s.
The velocity /elds of 6ame D are measured by Schneider
et al. (2003) and the scalar /elds are measured by Barlow
and Frank (1998). The measured scalars include tempera-
ture, mixture fraction, CH4, O2, CO2, H2O, CO, H2, OH and
NO. The mixture fraction � is de/ned as
�=
0:5(YH − YH;2)=WH + 2(YC − YC;2)=WC
0:5(YH − YH;2)=WH + 2(YC − YC;2)=WC : (1)
The stoichiometric mixture fraction for 6ameD is �st=0:351.
3. Modelling methods
Due to the high expense of PDF methods in dealing with
detailed chemical reaction mechanism, modelling methods
should be well designed to accommodate the computational
H. Wang, Y. Chen / Chemical Engineering Science 59 (2004) 3477–3490 3479
cost. In this work, the boundary layer approximation is
adopted to reduce the 6ow to parabolic type 6ow. The
joint scalar closure level PDF method is employed to re-
duce the dimension of the joint PDF and to facilitate the
usage of node-based Monte-Carlo algorithm (Pope, 1981)
which makes the solution easier. For the model closure, a
multiple-time-scale (MTS) k–� turbulence model of Kim
and Chen (1989) is introduced to calculate the turbulence
/eld, and the Euclidean minimum spanning trees (EMST)
mixing model developed by Subramaniam and Pope (1998)
is employed to model the small-scale mixing term appearing
in the PDF evolution equation. Chemical reactions are de-
scribed by the newest version of GRImechanism, GRI-Mech
3.0 (GRI Website, 1995), which consists of 53 species and
325 elemental reactions and contains nitrogen chemistry.
3.1. Parabolized Navier–Stokes equations
Parabolic type 6ow can be described by the parabolized
Navier–Stokes (PNS) equations (Rubin and Tannehill,
1992) which omit the turbulent transport along the stream
direction. Since the PNS equations can be solved by the
space marching algorithm which only needs to scan the
computational domain once to obtain the solutions and
does not need iteration, the computational time will be
greatly saved in contrast to general iteration method. Hence
the PNS equations are preferred in the present models. In
the cylindrical coordinate system, the PNS equations are
written as
@( M�ũ)
@x
+
1
r
@(r M�ṽ)
@r
= 0; (2)
@( M�ũũ)
@x
+
1
r
@(r M�ũṽ)
@r
=
1
r
@
@r
(
r�t
@ũ
@r
)
+ (�∞ − M�)g: (3)
The transport equation for radial momentum is not given
here since it is not necessary for the numerical solution (see
Section 4).
The turbulent eddy-viscosity �t can be modelled by stan-
dard models, such as k–� model and Reynolds stress model
(RSM), while model coeKcients should be adjusted to ob-
tain correct jet spreading rate, e.g. the c�2 is adjusted from
1.92 to 1.8 (Lindstedt et al., 2000). Here, an alternative tur-
bulence model, multiple-time-scale (MTS) k–� model (Kim
and Chen, 1989), is introduced to the present modelling of
turbulent axisymmetric jet non-premixed 6ame which does
not need coeKcient adjustment.
3.2. MTS k–� turbulence model
In the conventional k–� model, only the production and
dissipation of turbulent kinetic energy are accounted for,
whereas in the MTS k–� model of Kim and Chen (1989),
the in6uence of cascade of turbulent kinetic energy is also
involved except for the production and dissipation. The con-
cept of the MTS k–� model is that the turbulent kinetic en-
ergy spectrum is partitioned into two regions, production
range and dissipation range. The total turbulent kinetic en-
ergy k contains two parts, the energy of large eddies in
production range kp and the energy of /ne-scale eddies in
dissipation range kt (k = kp + kt). The large eddies turbu-
lent energy is generated by the mean 6ow instability and
cascades to /ner eddies. The energy transfer rate �p is in-
troduced to describe the eddies cascade. The /ne-scale ed-
dies turbulent energy is dissipated into thermal energy by
the viscous forces, which are characterized by the turbulent
kinetic energy dissipation rate �t (or �). The modelled trans-
port equations for kp, kt , �p and �t are as follows (Kim and
Chen, 1989).
@( M�ũkp)
@x
+
1
r
@(r M�ṽkp)
@r
=
1
r
@
@r
(
r
�t
�kp
@kp
@r
)
+ M�(P − �p); (4)
@( M�ũ�p)
@x
+
1
r
@(r M�ṽ�p)
@r
=
1
r
@
@r
(
r
�t
��p
@�p
@r
)
+
M�
kp
(cp1P2 + cp2P�p − cp3�2p); (5)
@( M�ũkt)
@x
+
1
r
@ (r M�ṽkt)
@r
=
1
r
@
@r
(
r
�t
�kt
@kt
@r
)
+ M�
(
�p − �t
)
; (6)
@( M�ũ�t)
@x
+
1
r
@(r M�ṽ�t)
@r
=
1
r
@
@r
(
r
�t
��t
@�t
@r
)
+
M�
kt
(ct1�2 + ct2�p�t − ct3�2t ); (7)
where P = �t= M�(@ũ=@r)2. The model constants in Eqs. (5)–
(7) are �kp=�kt=0:75, ��p=��t=1:15, cp1=0:21, cp2=1:24,
cp3 = 1:84, ct1 = 0:29, ct2 = 1:28, ct3 = 1:66.
The turbulent eddy viscosity is de/ned as
�t = M�c�fk2=�p; (8)
where c�f = 0:09. Based on the MTS turbulence model, a
time scale for turbulent kinetic energy dissipation can be
de/ned from large-eddy quantities as
�= kp=�p: (9)
3480 H. Wang, Y. Chen / Chemical Engineering Science 59 (2004) 3477–3490
3.3. Joint scalar PDF transport equation
Suppose there are �-dimension scalar vector �(x; t) =
{’i(x; t); i = 1; : : : ; �} in the 6ow /eld, and at position x
and at time t, the scalar joint PDF is f( ; x; t), where =
{ i; i = 1; : : : ; �} is the sample space vector corresponding
to �(x; t). Using the standard technology of Pope (1985),
the evolution equation for f( ; x; t) can be derived as
@
@x
( M�ũf̃) +
1
r
@
@r
(r M�ṽf̃)
=− 1
r
@
@r
(r〈v′′|�= 〉 M�f̃)− @
@ i
( M�Sif̃)
+
@
@
{
M�
〈
1
�
@(rJ ir)
r@r
∣∣∣∣�= 〉 f̃} : (10)
The /rst term on the right-hand side of Eq. (10) denotes
the transport of f( ; x; t) along the radial direction in the
physical space that is induced by turbulent 6uctuations. The
gradient transport assumption (Pope, 1985; Lindstedt et al.,
2000) is often used to model this term as
− 1
r
@
@r
(〈v′′|�= 〉 M�f̃) = 1
r
@
@r
(
r
�t
�f
@f̃
@r
)
; (11)
where �f = 0:9.
The second term on the right-hand side of Eq. (10) de-
notes the change of f( ; x; t) in sample space due to chem-
ical reaction. Obviously, the chemical reaction term is in
closed form, which enables the exact treatment of arbitrarily
complicated chemical reaction mechanisms. In this work,
the chemical reaction is described by GRI-Mech 3.0 for
methane oxidation and nitrogen oxides formation.
The third term on the right-hand side of Eq. (10) is the
conditional mean of molecular transport, namely small-scale
mixing term. It is not closed and has become one of the
main eForts for the study of PDF method. The often used
mixing model are IEM (interaction by exchange with the
mean) (IEM)model (Pope, 1985) Curl’s model and its mod-
i/cations (Pope, 1985; Lindstedt et al., 2000) etc. Recently,
Subramaniam and Pope (1998) developed a new mixing
model based on EMST which will be used in the present
simulation.
The radiative heat losses are not taken into account in
Eq. (10) due to their relative small impact on the tempera-
ture /eld of gaseous 6ames. Meanwhile the existing radia-
tion models, such as the optically thin limit model (Barlow
et al., 2001), are too simple to describe the interactions be-
tween turbulence and radiation (Li and Modest, 2002) and
the implementation of such radiation model will spend too
much extra CPU time. Therefore, for simplicity, the radia-
tion losses are ignored in the present 6ame temporarily.
The PDF evolution equation (10) is usually solved by
Monte Carlo particle method (Pope, 1981, 1985). The PDF
f( ; x; t) is represented by an ensemble of N particles in
terms of �∗j (j=1; : : : ; N ). Each particle evolves in physical
and sample space and the Favre average of scalars can be
estimated from the ensemble average of the particlesas
�̃(x; t) =
1
N
N∑
j=1
�∗j (x; t): (12)
3.4. EMST mixing model
One common drawback of popular mixing models, such
as the above-mentioned IEM model and Curl’s models, is
that they do not ful/ll the localness principle of mixing mod-
els (Subramaniam and Pope, 1998). With them, the particles
will mix homogeneously with each other in each grid cell
like a homogeneous volume reaction system. When they are
applied to diFusion combustion problem, the fuel and the
oxidizer might penetrate the reaction zone, which results in
the non-physical mixing of cold fuel and cold oxidizer. How-
ever, in fact, only those particles in the physical space neigh-
borhood can in6uence the local small-scale mixing behav-
ior, non-adjacent particle should not aFect each other. Since
the scalar /eld is continuous in physical space, the principle
of localness in physical space is equivalent to the localness
in scalar space. In order to remedy the drawback of conven-
tional mixing models, Subramaniam and Pope (1998) devel-
oped an EMST based mixing model, where the localness of
multiple scalars is de/ned through the Euclidean minimum
spanning tree constructed in scalar space, and the princi-
ple of localness is ful/lled through interactions of neighbor
particles in scalar space. The good performance of EMST
model has been shown by Xu and Pope (2000) in their PDF
calculations of turbulent non-premixed 6ames.
A complete description of the EMST model can be found
in Subramaniam and Pope (1998). Here, we brie6y describe
the implementation of the model. At any time, given the
ensemble of N particles in the grid cell, a subset of NT
particles is chosen from the particle ensemble for mixing,
according to an age property associated with each particle.
An EMST is constructed on this subset of NT particles. This
tree containsNT particles (nodes) andNT−1 edges, and each
particle is connected with at least one neighbor particle. The
mixing is simulated through the evolution of each particle
�∗j (j = 1; : : : ; NT ) as
wj
d�∗j
dt
=−(
NT−1∑
)=1
B)
×{(�∗j − �∗n)),jm) + (�∗j − �∗m)),jn)}; (13)
where the )th edge of the tree connects the particle pair
(m); n)). The speci/cation of the model constants ( and
B) is described in Subramaniam and Pope (1998). In the
determination of (, a mixing time-scale �’ is needed which
is usually related to the turbulence time-scale � (Eq. (9))
as
�’ = �=c’; (14)
H. Wang, Y. Chen / Chemical Engineering Science 59 (2004) 3477–3490 3481
where c’ is the mechanical to scalar time-scale ratio. Exper-
iments (Beguier et al., 1978; Panchapakesan and Lumely,
1993; Dai et al., 1995; Wang and Tong, 2002; Anderson
and Bremhorst, 2002) and DNS (Eswaran and Pope, 1988;
Johansson and Wikstrom, 1999) show that c’ is nearly a
constant. The standard value 2.0 is usually used for c’
(Spalding, 1971). Sometimes it performs better with small
tuning of the standard value, e.g. with 1.5 by Xu and Pope
(2000) and with 2.3 by Lindstedt et al. (2000). However,
actually, c’ is 6ow dependent, to account for the c’ vari-
ation in the 6ow /eld, Hidouri et al. (2003) have ever
tried a variable c’ method in the modelling of turbulent jet
6ame, while no substantial improvement is shown. There-
fore the constant method is still adopted in the present
calculation and the standard value c’ = 2:0 is chosen for
simplicity.
4. Numerical method
The space marching algorithm is employed to solve the
problem. The axial coordinate x is viewed as time-like
coordinate. The governing equations for mean velocity
(Eq. (3)) and turbulence (Eqs. (4)–(7)) are discretized
by FV (/nite-volume) method, namely upwind diFerence
for convection terms, central diFerence for diFusion terms
and implicit scheme for axial convection. The obtained
algebraic equations are solved by using iteration method.
The radial momentum equation is not solved and the radial
velocity is deduced from the continuum equation (6). The
computational domain along the radial direction is divided
into 50 cells, most of which are laid out around the cen-
terline. The step size along axial direction is /xed to be
Qx=D=0:05 within the whole computational domain except
the beginning, where Qx=D=0:005 is used temporarily.
Total about 2000 steps are marched from nozzle to down-
stream x=D∼ 80. Since the implicit scheme is used for the
turbulence equations discretization, there is no step size
restriction for turbulence /elds calculation. The time-like
step size can be speci/ed as Qt =Qx=ũ.
The PDF evolution equation (10) is solved by the
node-based Monte Carlo particle method (Pope, 1981)
which /xes the particle ensemble on grid cells. The deriva-
tives with respect to physical space coordinate in Eq. (10)
are /rstly discretized by FV method with explicit scheme
for axial convection. Through a time-splitting scheme,
the evolution of the PDF can be separated into three pro-
cesses, convection and diFusion process C, small-scale
mixing process M , and chemical reaction process R. Since
the explicit scheme is used for Eq. (10), the step size
Qx for three processes is restricted. To match the step
size speci/ed in the turbulence /elds, some sub-steps are
involved inside these processes according to CFL condi-
tions, e.g. nR sub-steps for R, nM sub-steps for M , and nC
sub-steps for C. Therefore the discrete form of Eq. (10) is
written as
f̃|x+Qx =
nR∏
i=1
(I +QtiR) ·
nM∏
j=1
(I +QtjM)
×
nC∏
k=1
(I +QtkC) · f̃|x; (15)
where
∑nR
i=1 Qti=
∑nM
j=1 Qtj=
∑nC
k=1 Qtk =Qt. In each grid
cell, 400 particles are allocated, and total 400× 50=20000
particles are located in the stream-cross section.
The boundary conditions need to be speci/ed. Most
of them are determined according to experimental con-
ditions (Schneider et al., 2003; Barlow and Frank, 1998;
Xu and Pope, 2000). The turbulent kinetic energy k at
the inlet boundary is estimated from measured Reynolds
normal stress. The turbulent kinetic energy of large eddies
in production range kp is assumed to be 80% of k. The
determination of the inlet condition for �p and �t (or �)
needs special attention since experimental data for them
is hardly available, and Merci et al. (2002) found that the
computational results are very sensitive to inlet condition
of �. We assume that �p is equal to �t at the inlet boundary
and estimate their radial pro/les through mixing length lm
as �p= �t= c
3=4
�f · k3=2=lm, where lm is given following Merci
et al. (2002) as:
lm
Dh
=
{
1− exp
[
−2× 106 ·
(
y
Dh
)3]}
×
[
1
15
−
(
1
2
− y
Dh
)4]
: (16)
The thermodynamic, chemical kinetic parameters and trans-
port properties are provided by CHEMKIN-II package
(Kee et al., 1989). The particle evolution equations due
to chemical reaction are solved by subroutine LSODE
(Radhakrishnan and Hindmarsh, 1993). EMST mixing
model is implemented using Fortran package from Ren
et al. (2002), where the numerical weight wj (j=1; : : : ; N )
of particles is set to be unity.
The calculation begins at x=D = 0. At /rst, The particle
evolutions of convection and diFusion, mixing, and reac-
tion process are solved in turn to obtain the reactive scalar
/elds. The mean density is calculated from the obtained tem-
perature and composition concentrations through ideal gas
state relation, and is substituted into the turbulence /elds to
compute mean velocity and turbulent eddy viscosity. When
the velocity and scalar /elds are all obtained, the compu-
tational plane will be advanced by Qx along axial direc-
tion and new computations continue. The cycle will not stop
until it reaches downstream x=D ∼ 80. The required CPU
time for present calculation is about a week on a PC (AMD
Athlon(tm) XP 2400+, 2:0 GHz, 512 MB memory).
3482 H. Wang, Y. Chen / Chemical Engineering Science 59 (2004) 3477–3490
5. Results and discussion
In this section, the numerical predictions of 6ame D are
presented and compared to the measurements of Barlow
and Frank (1998) and Schneider et al. (2003) extensively.
The explored results involve radial pro/les of mean and
rms (root mean square), conditionalmean, scatter plots of
scalars and conditional PDF distribution, etc. Due to space
limitation, only the results at three diFerent axial locations,
viz. x=D = 7:5, 15 and 30, are presented.
The predicted and measured radial pro/les of mean ax-
ial velocity, mean mixture fraction and rms mixture frac-
tion (�′′ ≡ �̃′′2
1=2
) are compared in Fig. 2. The predictions
of mean axial velocity and mixture fraction are in good
agreement with experimental data. The jet spreading rate is
well predicted. The general over-prediction of the spread-
ing rate of round jet by conventional k–� is remedied by
the adopted MTS k–� model, which shows the superior-
ity of the MTS model. The overall agreement between the
predicted and measured rms mixture fraction is also good,
while at x=D=7:5, the rms mixture fraction is over-predicted
by 40% and at x=D = 30, it is under-predicted by 30%.
The evolution of scalar variance is mainly determined by
small-scale mixing which is not closed in the PDF evolu-
tion equation. Therefore, the large discrepancy between the
predicted and measured rms mixture fraction is put down
to the mixing model. The modelling of mixing is the most
diKcult task for PDF methods, and by now no existing mix-
ing model can ful/l all the performance criteria of mixing
models (Subramaniam and Pope, 1998). Moreover, the only
single scalar time-scale for all scalars and the unwarranted
r/
D
0
1
2
x/D=7.5
r/
D
0
2
4
x/D=15
u (m/s)
r/
D
0 50
0
3
6 x/D=30
~
r/
D
0
1
2
x/D=7.5
r/
D
0
2
4
x/D=15
ξ
r/
D
0 0.5 1
0
3
6 x/D=30
~
r/
D
0
1
2
x/D=7.5
r/
D
0
2
4
x/D=15
ξ
r/
D
0 0.2
0
3
6 x/D=30
"
Fig. 2. Predicted (solid lines) and measured (symbols) radial pro/les of
mean axial velocity, mean mixture fraction and rms mixture fraction.
x/D=7.5 ;r=6mm
1
x/D=7.5 ;r=9mm
C
P
D
F
0 0.25 0.5 0.75
0.25 0.5 0.75
0.25 0.5 0.75 0.25 0.5 0.75 0.25 0.5 0.75
0.25 0.5 0.75 0.25 0.5 0.75
0.25 0.5 0.75 0.25 0.5 0.75
0
5
10
x/D=7.5 ;r=3mm
x/D=15 ;r=6mm
1
x/D=15 ;r=12mm
C
P
D
F
0
0
5
10
x/D=15 ;r=2mm
x/D=30 ;r=9mm
C
P
D
F
0
0
5
10
x/D=30 ;r=3mm
1
x/D=30 ;r=21mm
ξ
ξ
ξ ξ ξ
ξ ξ
ξ ξ
Fig. 3. Predicted (solid lines) and measured (dashed lines) mixture fraction
PDFs conditional on particular physical positions.
constant mechanical to scalar time-scale ratio c’ also cause
unpredictable error. Based on these situations, the observed
discrepancy between predicted and measured rms mixture
fraction is not strange and the error is acceptable. The pre-
diction of the radial position corresponding to the peak rms
mixture fraction is same as the measurement, which also
indicates the accurate prediction of the jet spreading rate.
The evolution trend of the predicted rms mixture fraction is
similar to the prediction of Xu and Pope (2000) in which
the same mixing model, viz. EMST model, was used. Fig. 3
compares the predicted and measured mixture fraction PDFs
at diFerent axial and radial locations. At the same axial dis-
tance, from inner to outer, the peak of mixture fraction PDF
moves from low mixture fraction to higher one. The predic-
tions reasonably represent this trends and the overall agree-
ment between the prediction and measurement is good. The
6ow and mixing /eld of 6ame D is well predicted by present
models, which provides the prerequisite to the accurate pre-
diction of the reactive scalar /elds.
The radial pro/les of mean temperature and mass frac-
tions of reactive species in 6ame D are shown in Figs. 4–6.
The predicted mean temperature in Fig. 4 and mass fractions
of major species , such as CH4, O2 in Fig. 4 and CO2, H2O
in Fig. 5, are in good agreement with the experimental data
of Barlow and Frank (1998), including the magnitude and
the radial positions corresponding to peak values. The good
prediction of the temperature indicates that, for the present
6ame, the impact of radiative heat losses is minor and the
omission of that is acceptable. The predicted mass fractions
of radicle OH in Fig. 5 and intermediate species CO, H2
H. Wang, Y. Chen / Chemical Engineering Science 59 (2004) 3477–3490 3483
r/
D
0
1
2
x/D=7.5
0
2
4
x/D=15
T(K)
0 1000
0
3
6 x/D=30
~
r/
D
0
1
2
x/D=7.5
0
2
4
x/D=15
YCH4
0 0.1
0
3
6 x/D=30
~
r/
D
0
1
2
x/D=7.5
r/
D
r/
D
r/
D
r/
D
r/
D
0
2
4
x/D=15
YO2
r/
D
0 0.1 0.2
0
3
6 x/D=30
~
Fig. 4. Predicted (solid lines) and measured (symbols) radial pro/les of
mean temperature and mean mass fractions of CH4 and O2.
r/
D
r/
D
r/
D
r/
D
r/
D
r/
D
r/
D
r/
D
r/
D
0
1
2
x/D=7.5
0
2
4
x/D=15
YCO2
0 0.1
0
3
6 x/D=30
~
0
1
2
x/D=7.5
0
2
4
x/D=15
YH2O
0 0.1
0
3
6 x/D=30
~
0
1
2
x/D=7.5
0
2
4
x/D=15
YOH
0 0.002
0
3
6 x/D=30
~
Fig. 5. Predicted (solid lines) and measured (symbols) radial pro/les of
mean mass fractions of CO2, H2O and OH.
in Fig. 6 are also in reasonable agreement with measure-
ments. Whereas the pollutant NO emission is over-predicted
by 30% in Fig. 6. The omission of radiative heat losses may
lead to the over-prediction of NO emission, while, for the
present 6ame, this in6uence may be ignored due to the good
prediction of temperature. Barlow et al. (2001) found that
the NO formation in laminar opposed-6ow partially pre-
mixed methane/air 6ames is also over-predicted by using
r/
D
r/
D
r/
D
r/
D
r/
D
r/
D
r/
D
r/
D
r/
D
0
1
2
x/D=7.5
0
2
4
x/D=15
YCO
0 0.05
0
3
6 x/D=30
~
0
1
2
x/D=7.5
0
2
4
x/D=15
YH2
0 0.002
0
3
6 x/D=30
~
0
1
2
x/D=7.5
0
2
4
x/D=15
YNO
0 5E-05
0
3
6 x/D=30
~
Fig. 6. Predicted (solid lines) and measured (symbols) radial pro/les of
mean mass fractions of CO, H2 and NO.
0
0.3
0.6
0.9
0
0.3
0.6
0.9
0
0.3
0.6
0.9
0
0.3
0.6
0.9
0
0.3
0.6
0.9
0
0.3
0.6
0.9
0
0.3
0.6
0.9
0
0.3
0.6
0.9
0
0.3
0.6
0.9
x/D=7.5
x/D=15
T(K)
0 1000 2000
x/D=30
x/D=7.5
x/D=15
YCH4
0 0.1
x/D=30
x/D=7.5
x/D=15
YO2
0 0.2
x/D=30
ξ ξ ξ
ξ
ξ
ξ
ξ
ξ
ξ
Fig. 7. Predicted (solid lines) and measured (symbols) conditional mean
temperature and mass fractions of CH4 and O2 versus mixture fraction.
GRI-Mech 3.0 and including the radiative heat losses, which
implies that the over-prediction of NO formation is a com-
mon drawback of GRI-Mech 3.0. Excluding this drawback,
we may say that the NO emission from 6ame D is well pre-
dicted by present models.
To further explore the turbulence–chemistry interactions
directly, the predicted conditional mean temperature and
species mass fractions versus mixture fraction are compared
with the experimental data in Figs. 7–9. The agreements be-
tween the predictions and measurements are very good for
3484 H. Wang, Y. Chen / Chemical Engineering Science 59 (2004) 3477–3490
0
0.3
0.6
0.9
0
0.3
0.6
0.9
0
0.3
0.6
0.9
0
0.3
0.6
0.9
0
0.3
0.6
0.9
0
0.3
0.6
0.9
0
0.3
0.6
0.9
0
0.3
0.6
0.9
0
0.3
0.6
0.9
x/D=7.5
x/D=15
YCO2
0 0.1
x/D=30
x/D=7.5
x/D=15
YH2O
0 0.1
x/D=30
x/D=7.5
x/D=15
YOH
0 0.003
x/D=30
ξ
ξ
ξ ξ
ξ
ξ ξ
ξ
ξ
Fig. 8. Predicted (solid lines) and measured (symbols) conditional mean
mass fractions of CO2, H2O and OH versus mixture fraction.
0
0.3
0.6
0.9
0
0.3
0.6
0.9
0
0.3
0.6
0.9
0
0.3
0.6
0.9
0
0.3
0.6
0.9
0
0.3
0.6
0.9
0
0.3
0.6
0.9
0
0.3
0.6
0.9
0
0.3
0.6
0.9
x/D=7.5
x/D=15
YCO
0 0.05
x/D=30
x/D=7.5
x/D=15
YH2
0 0.003
x/D=30
x/D=7.5
x/D=15
YNO
0 0.0001
x/D=30
ξ
ξ
ξ ξ
ξ
ξ ξ
ξ
ξ
Fig. 9. Predicted (solid lines) and measured (symbols) conditional mean
mass fractions of CO, H2 and NO versus mixture fraction.
conditional mean temperature in Fig. 7 and major species
(CH4, O2 in Fig. 7 and CO2, H2O in Fig. 8) mass fractions.
The predicted shape of radicle OH conditional mean mass
fraction is somewhat narrower than the experimental data,
and the peak value of OH at x=D=7:5 is a bit over-predicted,
while the overall agreements between predictions and mea-
surements are still satisfying. The intermediate species CO,
H2 in Fig.9 are also well predicted. It has been reported
that the intermediate species in the fuel-rich region (�¿�st)
tend to be over-predicted by both steady 6amelet model and
/rst-order CMC (Barlow et al., 2001; Coelho and Peters,
2001; Roomina and Bilger, 2001). The reaction progress
in the region �¿�st seems to be over-estimated by steady
Fig. 10. Scatter plot of measured (left) and predicted (right) temperature
versus mixture fraction. Lines: stretched laminar 6amelet (a= 100 s−1).
6amelet model and /rst-order CMC, while it is well repre-
sented by PDF model. The conditional mean mass fraction
of NO is also over-predicted like its radial pro/le shown in
Fig. 6. With regard to the mentioned drawback of GRI-Mech
3.0, the NO emission is also well represented by present
models. The fairly good agreements between the predicted
and measured conditional means in mixture fraction space
demonstrate the capacity of PDF model to represent the
turbulence–chemistry interactions. In the physical space,
the agreements between the predictions and measurements
(shown in Figs. 4–6) is similar to those in mixture fraction
space. This bene/ts from the good predictions of 6ow and
mixing /eld except for the capacity of PDF model.
In the node-based Monte Carlo algorithm, the scalar joint
PDF f̃( ; x; r) is represented by an ensemble of particles
located at position (x; r). Each particle can be viewed as a
realization of the turbulent reactive 6ow at (x; r). Therefore,
the particle ensemble is comparable to the instantaneous
measurements of 6ame D. It is informative to compare the
scatter plot of the particle ensemble with the instantaneous
experimental data of scalars, although the comparisons have
limitations (Xu and Pope, 2000).
Fig. 10 shows the scatter plot of measured and predicted
temperature versus mixture fraction at three axial locations
x=D= 7:5, 15 and 30. For the sake of reference, a stretched
laminar 6amelet pro/le under a moderate strain rate (a =
100 s−1) is also plotted in the /gure. The 6amelet pro/le
is obtained by solving the steady 6amelet equations (Peters,
1984). The overall agreement between the experimental data
H. Wang, Y. Chen / Chemical Engineering Science 59 (2004) 3477–3490 3485
Fig. 11. Scatter plot of measured (left) and predicted (right) CH4
mass fraction versus mixture fraction. Lines: stretched laminar 6amelet
(a = 100 s−1).
and the prediction is comparatively good. Most of the points
distribute around the 6amelet line, while some points de-
parture far away from the 6amelet line, e.g. the temperature
near �st drops to 1000 K. The suppressed temperature in the
reaction zone indicate the 6ame local extinction (Barlow and
Frank, 1998). This is the re6ection of the strong turbulence
–chemistry interactions. The amount of local extinction can
be estimated from the number of points that departure away
from the 6amelet line. According to the experimental data
(left in Fig. 10), the local extinction shows low probability
at x=D=7:5 and gets more at x=D=15, while the amount of
local extinction decreases at x=D= 30, which indicates that
some extinguished samples return to the 6amelet state again,
namely re-ignition process. The local extinction process is
well represented by present models. While the amount of
local extinction seems to be a bit over-predicted at x=D=30
and the re-ignition process is postponed to downstream. In
the center of the scatter plot of predicted temperature, There
is a void hole surrounded by the extinguished samples and
by samples near the 6amelet line, where almost no sam-
ples locate. The experimental data does not show this phe-
nomenon. This may come from the “stranding” problem of
EMST mixing model (Subramaniam and Pope, 1998).
Fig. 11 shows the scatter plot of measured and predicted
CH4 mass fraction versus mixture fraction. In the fuel-lean
region (�¡�st), the CH4 mass fraction of 6amelet is almost
zero, while some measured instantaneous CH4 mass fraction
has non-zero value. The oxidizer of 6ame D is air, so the
Fig. 12. Scatter plot of measured (left) and predicted (right) CO
mass fraction versus mixture fraction. Lines: stretched laminar 6amelet
(a = 100 s−1).
fuel existing in the fuel-lean region is penetrated from the
fuel side through the reaction zone. These samples must be
extinguished locally and are transported through molecular
diFusion to the oxidizer side (If not, they will be consumed
within the reaction zone). The prediction of scatter plot of
CH4 mass fraction is in good agreement with the experiment.
The number of predicted non-zero samples in the fuel-lean
region are more than the experiment at x=D = 30, which is
consistent with the over-predicted amount of local extinction
at the same axial location.
The scatter plot of measured and predicted CO mass frac-
tion versus mixture fraction is illustrated in Fig. 12. The
measured scatter plot of CO mass fraction shows much more
dispersion than that of temperature (see Fig. 10). Most sam-
ples departure from the 6amelet line. The prediction repre-
sents this property fairly well. However the prediction shows
conglomeration of some samples at some places. This is not
re6ected by the experimental data and may result from the
“stranding” problem of EMST mixing model.
Fig. 13 shows the scatter plot of measured and predicted
NO mass fraction versus mixture fraction. The strain rate
of referenced stretched laminar 6amelet here is diFerent
from the former, viz. a = 800 s−1 instead of 100 s−1, be-
cause the NO pro/le of 6amelet under a = 100−1 is much
higher than the NO emission in 6ame D. The NO emission
is over-predicted, and the reason has been discussed in the
conditional mean of NO mass fraction in Fig. 9. There is
also a void hole in the center of the predicted scatter plot
of NO just like the scatter plot of temperature in Fig. 10.
3486 H. Wang, Y. Chen / Chemical Engineering Science 59 (2004) 3477–3490
Fig. 13. Scatter plot of measured (left) and predicted (right) NO
mass fraction versus mixture fraction. Lines: stretched laminar 6amelet
(a = 800 s−1).
The experiment of NO does not exhibit such structure and
the “stranding” phenomenon of EMST mixing model is also
thought to be responsible for this.
Other species such as O2, CO2, H2O, H2 and OH scatter
plot are also well represented by present models (/gures not
shown).
Fig. 14 shows the CPDFs of temperature and CO2, H2O
mass fractions at three diFerent axial locations x=D=7:5, 15
and 30. The mixture fraction interval used for the estimation
of CPDFs is also shown in the /gure for diFerent scalars, e.g.
(0:3¡�¡ 0:4) for temperature. From x=D= 7:5 to 30, the
shape of measured temperature CPDF become thinner and
the peak value shifts toward higher temperature. This trends
are well predicted by present models. However, the predicted
CPDF shapes are more thinner than measurements, which
results in a higher peak value of CPDF. The distributions
of CPDFs of CO2 and H2O in Fig. 14 are similar to that of
temperature.
The measured and predicted CPDFs of CO, H2 and OH
mass fractions are illustrated in Fig. 15. The moving trend
of the peak value of the measured CPDFs of CO and H2
is same as that of temperature in Fig. 14, and it is well
represented by the predictions. The predicted shapes of the
CPDFs of CO and H2 are also thinner than measurements. At
x=D=30, the peak value of CO and H2 CPDFs is very high,
which is in accordance with the conglomeration of samples
in the scatter plot of CO in Fig. 12. At this axial location,
a bimodal shape of the CO and H2 CPDFs is observed,
which is not re6ected by the measurements. The bimodal
0.05 0.1 0.15
PDF (0.3<ξ<0.4)
YCO2 YCO2
0.05 0.1 0.15
0
50
100
EXP (0.3<ξ<0.4)
0.05 0.1 0.15
PDF (0.35<ξ<0.45)
YH2O YH2O
0.05 0.1 0.15
0
50
100
EXP (0.35<ξ<0.45)
T(K)
1500 2000
PDF (0.3<ξ<0.4)
T(K)
C
P
D
F
C
P
D
F
C
P
D
F
1000 1500 2000
0
0.01
x/D=7.5
x/D=15
x/D=30
EXP (0.3<ξ<0.4)
Fig. 14. Conditional PDFs of measured (left) and predicted (right) tem-
perature and CO2, H2O mass fractions.shape of CPDFs corresponds to the void hole structure in
the center of the scatter plots (see Fig. 12). The peak value
of measured OH CPDF moves oppositely to the lower side
when the axial location moves downstream, which is well
represented by the prediction. The CPDF of NO is shown in
Fig. 16. Since the NO emission is over-predicted by
GRI-Mech 3.0, the predicted peak value locates at higher
concentration. The predicted CPDF shapes of NO are bi-
modal, which is in accordance with the void hole structure
in the center of scatter plot of NO in Fig. 13. The moving
trend of peak value is also well re6ected by the prediction,
namely to higher concentration with the increase of the
axial distance from the nozzle.
The 6ame structures of 6ame D are well predicted by
present models, namely MTS k–� model, joint scalar PDF
model, EMST mixing model and GRI-Mech 3.0 summarily.
Compared to the literature, our predictions of 6ame D are
very similar to those of Xu and Pope (2000). In their work,
the ingredients of used models consist of a joint velocity–
composition–turbulence frequency PDF model, EMST mix-
ing model and ARMmechanism of Sung et al. (1998). They
report good numerical predictions and make extensive com-
parisons between their predictions and experimental data.
According to the comparisons with their results, the present
predictions validate the reduced mechanism ARM used by
H. Wang, Y. Chen / Chemical Engineering Science 59 (2004) 3477–3490 3487
YCO
0.05 0.1
PDF (0.43<ξ<0.53)
YH2
0.0025 0.005
PDF (0.48<ξ<0.58)
YH2
C
P
D
F
0.0025 0.005
0
500
1000
1500
EXP (0.48<ξ<0.58)
YOH
0.0025 0.005
PDF (0.28<ξ<0.36)
YOH
C
P
D
F
0 0.0025 0.005
0
500
1000
EXP (0.28<ξ<0.36)
YCO
C
P
D
F
0.05 0.1
0
50
100
x/D=7.5
x/D=15
x/D=30
EXP (0.43<ξ<0.53)
Fig. 15. Conditional PDFs of measured (left) and predicted (right) CO,
H2 and OH mass fractions.
YNO
5E-05 0.0001
x/d=7.5
x/d=15
x/d=30
PDF (0.33< ξ<0.41)
YNO
C
P
D
F
0 5E-05 0.0001
0
50000 EXP (0.33< ξ<0.41)
Fig. 16. Conditional PDFs of measured (left) and predicted (right) NO
mass fraction.
Xu and Pope (2000) through present predictions using de-
tailed mechanism GRI-Mech 3.0, and the capacity of the
adopted scalar joint PDF model with MTS k–�model is also
demonstrated by the more sophisticated self-contained joint
velocity-composition-turbulence frequency PDFmodel used
by Xu and Pope (2000). The same mixing model (EMST
model) is adopted by the two predictions. Due to model
limitation, some phenomena, such as the void hole in the
center of the scatter plot and conglomeration of samples,
are brought on in Xu and Pope (2000) and here, which is
not re6ected by the experimental data. Further research on
the small-scale process is needed and sophisticated mixing
model is desired to represent the turbulence–chemistry in-
teractions subtly. More comparisons of predictions with the
literature is restricted by the limited report of numerical re-
sults, e.g. Lindstedt et al. (2000), Tang et al. (2000), Xiao
et al. (2000) etc.
Another comparison to the literature is about the compu-
tational time. It should be noted that the comparison has lim-
itations since too much ingredients are involved, e.g. hard-
ware, software, adopted modelling methods and numerical
parameters etc, although the comparison is informative. The
novelty of present work is the incorporation of the detailed
reaction mechanism (GRI-Mech 3.0) in the PDF calculation
of 6ame D, which is the most hugest mechanism ever imple-
mented in the practical calculations of turbulent combustion
problems by using PDF methods. The required CPU time for
the present calculation is about a week on a PC. In Xu and
Pope (2000), the reported CPU time is about 3 weeks un-
der their modelling strategies. The elliptic solution method
is adopted in their work, so iterations are needed, while in
the present predictions, lower PDF closure level is adopted
and the space marching algorithm is employed, which scans
the computational domain only once to obtain the solutions.
Therefore, the present computational cost is much less than
Xu and Pope (2000). Our solution methods are very simi-
lar to those used by Lindstedt et al. (2000), which consist
of Reynolds stress model (RSM), joint scalar PDF model,
modi/ed Curl’s model, reduced mechanisms (consisting of
16 independent, 4 dependent and 28 steady-state scalars),
and distributed particle Monte Carlo algorithm. The required
CPU time for the modelling strategies of Lindstedt et al.
(2000) is of the order of 3 weeks, that is also more than the
time required for the present calculation with detailed chem-
istry. In Lindstedt et al. (2000), about 10,000 axial steps are
used to cover the 6ame length of ∼ 80 jet diameters. How-
ever, the present axial step size for turbulence /elds calcula-
tion and scalar /elds calculation is diFerent. Fixed step size
is used for turbulence /elds, viz. Qx=D=0:05, and adaptive
sub-steps are used within the /xed step for the scalar /elds
(see Section 4). Total about 2000 steps are used for the tur-
bulence /eld. The number of sub-steps of diFerent evolution
process, e.g. convection and diFusion, mixing, and reaction,
and of diFerent grid cells along the radial direction is dif-
ferent. Therefore, the adaptive sub-steps will certainly re-
duce the computational cost compared to the uniform steps
(e.g. Lindstedt et al., 2000). Although the comparison about
computational time is not very strict, it shows the impor-
tance of the matching of diFerent closure elements and the
feasible implementation of detailed chemistry in the PDF
calculation.
In addition, only Sandia Flame D has been successfully
computed by the present model. Flame E and F, where Xu
and Pope (2000) and Lindstedt et al. (2000) claim some suc-
cess with PDF method, has not yet been calculated with the
detailed chemistry (GRI-Mech 3.0) and seem a bit diKcult
to be conducted due to their large amount of local extinc-
tion which is resulted from the even intenser interactions
between turbulence and chemistry in these 6ame. The
3488 H. Wang, Y. Chen / Chemical Engineering Science 59 (2004) 3477–3490
prediction of Flame F is very sensitive to the mechanical
to scalar time-scale ratio c’ in Xu and Pope (2000). It
seems that c’, or broadly the mixing, is very crucial to the
modelling of the 6ames with intense interactions between
turbulence and chemistry, especially Flame F. More eForts
are desired on this issue, such as the sensitivity of PDF
predictions to mixing models (Ren and Pope, 2004) and
c’, the determination of c’ and about variable c’ in the
6ow /eld etc. Moreover the performance of detailed mech-
anism GRI-Mech 3.0 in these 6ames will also need to be
examined in the future.
6. Conclusions
In this paper, the detailed reaction mechanism, GRI-Mech
3.0, is incorporated into the PDF modelling of turbulent
non-premixed jet 6ame (Sandia 6ame D) for the /rst time.
The adopted models contain theMTS k–�model, joint scalar
PDF method and EMST small-scale mixing model. The 6ow
is reduced to a parabolic type problem and is described
by parabolized Navier–Stokes (PNS) equations. The space
marching algorithm is adopted to solve this problem, which
reduces the computational cost greatly. The large /xed step
size for turbulence /elds calculation is adopted and adaptive
sub-steps are used for scalar /elds calculation. The required
CPU time for the calculation is about a week on a modern
PC, which is acceptable for present computational capacity.
The numerical results are presented in detail and are com-
pared to the experimental data extensively. The predictions
are in fairly good agreement with the measurements, includ-
ing the radial pro/les of mean and rms, conditional mean,
scatter plots and conditional PDFs. The shown and com-
pared quantities involve axial mean velocity, mean and rms
mixture fraction, temperature and mass fractions of species
CH4; O2; CO2; H2O, CO, H2, OH and NO. The NO emis-
sion is over-predicted by GRI-Mech 3.0. Due to the limita-tion of EMST model, some phenomena are resulted in, such
as the void hole in the center of the scatter plots and con-
glomeration of samples, which is not re6ected by measure-
ments. Sophisticated mixing models are desired to improve
the capacity of PDF model further.
Notation
a strain rate
B) EMST model constant
c�f constant coeKcient in Eq. (8)
cpj, ctj model constants in Eqs. (5) and (7),
respectively (j = 1–3)
c’ mechanical to scalar time-scale ratio
D diameter of fuel jet nozzle
DP diameter of annular piloted 6ame
C convection and diFusion process
Dh hydraulic diameter
f( ; x; t) scalar joint probability density func-
tion (PDF)
g acceleration of gravity
J ir molecular transport 6ux of ith scalar
along radial direction
kp turbulent kinetic energy of large ed-
dies in production range
kt turbulent kinetic energy of
/ne-scale eddies in dissipation
range
k turbulent kinetic energy, =kp + kt
lm mixing length
M small-scale mixing process
nC;M ;R number of sub-steps for process C,
M and R, respectively
N number of particles
NT number of particles in the EMST
P production rate of turbulent kinetic
energy
R chemical reaction process
Si reaction rate of ith scalar
t time
u axial velocity
v radial velocity
w numerical weight of particle
WX molecular or atomic weight of
species or element X
x coordinate vector
x, r cylindrical coordinate
y normal distance from the nearest
solid boundary
YX mass fraction of species or element
X
〈·〉 mathematical expectations
Greek letters
� dimension of scalar vector
( model parameter controlling the
variance decay rate
, Kronecker delta
Qt time step size
Qx step size along axial direction
�p energy transfer rate from production
range to dissipation range
�, �t dissipation rate of turbulent kinetic
energy
�t turbulent eddy-viscosity
� mixture fraction
� density
�kp, �kt , ��p, ��t , �f turbulent Prandtl number for kp, kt ,
�p, �t and f, respectively
� time scale of turbulent kinetic en-
ergy dissipation
�’ mixing time-scale
H. Wang, Y. Chen / Chemical Engineering Science 59 (2004) 3477–3490 3489
� equivalence ratio
’i ith scalar
i sample space variable of ith scalar
� scalar vector (={’i; i = 1; : : : ; �})
�∗ scalar vector of particle ensemble
sample space scalar vector
(={ i; i = 1; : : : ; �})
Subscripts
1 fuel jet 6ow
2 oxidizer 6ow
st stoichiometric condition
∞ in/nity
Superscripts
− conventional average
∼ Favre average
′′ 6uctuation from Favre average; rms
(root mean square)
Acknowledgements
This work was supported by the Special Funds for
Major State Basic Research Projects (G1999022207), PR
China and the National Natural Science Foundation of
China (50206021). We are grateful to Dr. R.S. Barlow and
Dr. A. Dreizler to provide us the experimental data for
scalar /elds and velocity /elds of 6ame D, respectively.
References
Anderson, S.M., Bremhorst, K., 2002. Investigation of the 6ow /eld of a
highly heated jet of air. International Journal of Heat and Fluid Flow
23, 205–219.
Barlow, R.S., Frank, J.H., 1998. EFects of turbulence on species mass
fractions in methane/air jet 6ames. In: 27th Symposium (International)
on Combustion. The Combustion Institute, Pittsburgh, PA,
pp. 1087–1095.
Barlow, R.S., Karpetis, A.N., Frank, J.H., Chen, J.-Y., 2001. Scalar
pro/les and NO formation in laminar opposed-6ow partially premixed
methane/air 6ames. Combustion and Flame 127, 2102–2118.
Beguier, C., Dekeyser, I., Launder, B.E., 1978. Ratio of scalar and velocity
dissipation time scales in shear 6ow turbulence. Physics of Fluids 21,
307–310.
Bilger, R.W., 1993. Conditional moment closure for turbulent reacting
6ow. Physics of Fluids A 5, 436–444.
Bilger, R.W., Starner, S.H., Kee, R.J., 1990. On reduced mechanisms
for methane-air combustion in nonpremixed 6ames. Combustion and
Flame 80, 135–149.
Chen, J.-Y., Kollmann, W., Dibble, R.W., 1989. PDF modeling of
turbulent nonpremixed methane jet 6ames. Combustion Science and
Technology 64, 315–346.
Coelho, P.J., Peters, N., 2001. Unsteady modelling of a piloted methane/air
jet 6ame based on the Eulerian particle 6amelet model. Combustion
and Flame 124, 444–465.
Dai, Z., Tseng, L.K., Faeth, G.M., 1995. Velocity/mixture fraction
statistics of round, self-preserving, buoyant turbulent plumes. Journal
of Heat Transfer-Transactions of the ASME 117, 918–926.
Dopazo, C., 1994. Recent developments in PDF methods. In: Libby, P.A.,
Williams, F.A. (Eds.), Turbulent Reacting Flows. Academic, New
York (Chapter 7).
Eswaran, V., Pope, S.B., 1988. Direct numerical simulations of the
turbulent mixing of a passive scalar. Physics of Fluids 31, 506–520.
GRI-Mech website, 1995. http://www.me.berkeley.edu/gri mech/.
Hewson, J.C., Kerstein, A.R., 2002. Local extinction and reignition in
nonpremixed turbulent CO=H2=N2 jet 6ames. Combustion Science and
Technology 174, 35–66.
Hidouri, A., Gazzah, M.H., Ticha, H.B., et al., 2003. Dynamic and scalar
turbulent 6uctuation in a diFusion 6ame of an-axisymmetric methane
jet into air. Computational Mechanics 31, 253–261.
International Workshop on Measurement and Computation of
Turbulent Nonpremixed Flames, 1996. http://www.ca.sandia.gov/TNF/
abstract.html.
Johansson, A.V., Wikstrom, P.M., 1999. DNS and modelling of passive
scalar transport in turbulent channel 6ow with a focus on scalar
dissipation rate modelling. Flow, Turbulence and Combustion 63,
223–245.
Jones, W.P., Kakhi, M., 1997. Application of the transported PDF
approach to hydrocarbon-air turbulent jet diFusion 6ames. Combustion
Science and Technology 129, 393–430.
Jones, W.P., Kakhi, M., 1998. PDF modeling of /nite-rate chemistry
eFects in turbulent nonpremixed jet 6ames. Combustion and Flame
115, 210–229.
Kee, R.J., Rupley, F.M., Miller, J.A., 1989. CHEMKIN-II: a fortran
chemical kinetics package for the analysis of gas-phase chemical
kinetic. Sandia Report, SAND89-8009.
Kim, S.-W., Chen, C.-P., 1989. A multiple-time-scale turbulence model
based on variable partitioning of the turbulent kinetic energy spectrum.
Numerical Heat Transfer, Part B 16, 193–211.
Klimenko, A.Y., 1990. Multicomponent diFusion of various admixtures
in turbulent 6ow. Fluid Dynamics 25, 327–334.
Li, G.N., Modest, M.F., 2002. Application of composition PDF methods
in the investigation of turbulence–radiation interactions. Journal of
Quantitative Spectroscopy and Radiative Transfer 73, 461–472.
Lindstedt, R.P., Louloudi, S.A., Vaos, E.M., 2000. In: Proceedings of
the Combustion Institute. The Combustion Institute, Pittsburgh, PA,
Vol. 28, pp. 149–156.
Merci, B., Dick, E., Vierendeels, J., Langhe, C.D., 2002. Determination
of � at inlet boundaries. International Journal of Numerical Methods
for Heat and Fluid Flow 12, 65–80.
Panchapakesan, N.R., Lumely, J.L., 1993. Turbulence measurements in
axisymmetric jets of air and helium. Part 2. Helium jet. Journal of
Fluid Mechanics 246, 225–247.
Peters, N., 1984. Laminar diFusion 6amelet models in non-premixed
turbulent combustion. Progress in Energy and Combustion Science 10,
319–339.
Peters, N., Kee, R.J., 1987. The computation of stretched laminar
methane-air diFusion 6ames using a reduced four-step mechanism.
Combustion and Flame 68, 17–29.
Pitsch, H., 2002. Improved pollutant predictions in large-eddy simulations
of turbulent non-premixed combustion by considering scalar dissipation
rate 6uctuations. In: Proceedings of the Combustion Institute. The
Combustion Institute, Pittsburgh, PA, Vol. 29, pp. 1971–1978.
Pitsch, H., Chen, M., Peters, N., 1998. Unsteady 6amelet modeling
of turbulent hydrogen-air diFusion 6ames. In: 27th Symposium
(International) on Combustion. The Combustion Institute, Pittsburgh,
PA, pp. 1057–1064.
Pope, S.B., 1981. A Monte Carlo method for the PDF equations of
turbulent reactive 6ow. Combustion Science and Technology 25,
159–174.
Pope, S.B., 1985. PDF methods for turbulent reactive 6ows. Progress in
Energy and Combustion Science 11, 119–192.
Pope, S.B., 1994. Lagrangian PDF methods for turbulent 6ows. Annual
Review of Fluid Mechanics 26, 23–63.
http://www.me.berkeley.edu/gri_mech/http://www.ca.sandia.gov/TNF/abstract.html
http://www.ca.sandia.gov/TNF/abstract.html
3490 H. Wang, Y. Chen / Chemical Engineering Science 59 (2004) 3477–3490
Pope, S.B., 1997. Computationally eKcient implementation of combustion
chemistry using in situ adaptive tabulation. Combustion Theory and
Modelling 1, 41–63.
Radhakrishnan, K., Hindmarsh, A.C., 1993. Description and use of
LSODE, the Livermore solver for ordinary diFerential equations.
Livermore National Laboratory Report, UCRL-1D-113855.
Ren, Z., Subramaniam, S., Pope, S.B., 2002. Implementation of the EMST
mixing model. http://mae.cornell.edu/∼laniu/emst.
Ren, Z.Y., Pope, S.B., 2004. An investigation of the performance of
turbulent mixing models. Combustion and Flame 136, 208–216.
Roomina, M.R., Bilger, R.W., 2001. Conditional moment closure (CMC)
predictions of a turbulent methane-air jet 6ame. Combustion and Flame
125, 1176–1195.
Rubin, S.G., Tannehill, J.C., 1992. Parabolized/reduced Navier–Stokes
computational techniques. Annual Review of Fluid Mechanics 24,
117–144.
Saxena, V., Pope, S.B., 1999. PDF simulations of turbulent combustion
incorporating detailed chemistry. Combustion and Flame 117,
340–350.
Schneider, C.H., Dreizler, A., Janicka, J., Hassel, E.P., 2003. Flow /eld
measurements of stable and locally extinguishing hydrocarbon-fuelled
jet 6ames. Combustion and Flame 135, 185–190.
Spalding, D.B., 1971. Concentration 6uctuations in a round turbulent free
jet. Chemical Engineering Science 26, 95–107.
Subramaniam, S., Pope, S.B., 1998. A mixing model for turbulent reactive
6ows based on Euclidean minimum spanning trees. Combustion and
Flame 115, 487–514.
Sung, C.J., Law, C.K., Chen, J.-Y., 1998. An augmented reduced
mechanism for methane oxidation with comprehensive global
parametric validation. In: 27th Symposium (International) on
Combustion. The Combustion Institute, Pittsburgh, PA, pp. 295–304.
Sung, C.J., Law, C.K., Chen, J.-Y., 2001. Augmented reduced mechanism
for NO emission in methane oxidation. Combustion and Flame 125,
906–919.
Tang, Q., Xu, J., Pope, S.B., 2000. Probability density function
calculations of local extinction and NO formation in piloted-jet
turbulent methane/air 6ame. In: Proceedings of the Combustion
Institute. The Combustion Institute, Pittsburgh, PA, Vol. 28,
pp. 133–139.
Wang, D., Tong, C., 2002. Conditionally /ltered scalar dissipation,
scalar diFusion, and velocity in a turbulent jet. Physics of Fluids 14,
2170–2185.
Xiao, K., Schmidt, D., Maas, U. 2000. A study of local 6ame structures in
piloted jet diFusion 6ame using a probability density function method.
In: Proceedings of the Combustion Institute. The Combustion Institute,
Pittsburgh, PA, Vol. 28, pp. 157–165.
Xu, J., Pope, S.B., 2000. PDF calculations of turbulent nonpremixed
6ames with local extinction. Combustion and Flame 123, 281–307.
http://mae.cornell.edu/~laniu/emst
PDF modelling of turbulent non-premixed combustionwith detailed chemistry
Introduction
Sandia Flame D
Modelling methods
Parabolized Navier--Stokes equations
MTS k--epsilon turbulence model
Joint scalar PDF transport equation
EMST mixing model
Numerical method
Results and discussion
Conclusions
Acknowledgements
References
