Combustion and Flames - ČVUT

busyicicleMechanics

Feb 22, 2014 (3 years and 8 months ago)

232 views

Combustion, Flames
and CFD


Rudolf Žitný, Ústav procesní a
zpracovatelské techniky ČVUT FS 2010

HEAT PROCESSES

HP11

CFD analysis of combustion. Transport equations and modeling of chemical reactions.
Non
-
premix combustion and mixture fraction methods. Lagrangian and Eulerian
methods.

HP1
1

Combustion
and Flames

Mark Bradford

Warnatz J.: Combustion. Springer,
Berlin 1996

Homogeneous
reaction in gases

Premixed (
only
one
inlet stream of mixed
fuel and oxidiser)

Non premixed
(separate fuel and
oxidiser inlets)

Laminar
flame

Turbulent
flame

Laminar
flame

Turbulent
flame

Use mixture
fraction method
(PDF)

Use E
BU

(
E
ddy
Br
eak Up models)

Liquid fuels (spray
combustion)

Combustion of
particles (coal)

Lagrangian method
-
trajectories of a
representative set of droplets/particles
in a continuous media

Lagrangian method
-
trajectories of a
representative set of droplets/particles
in a continuous media

HP1
1

Combustion


Premixed combustion

Fuel and oxidiser are mixed first and even then burn

HP1
1

Combustion

Laminar premixed flames are characterised by a narrow flame front (front width depends upon the rate of chemical reaction and

diffusion coefficient). Laminar burning velocity (u
L
) must be greater than the velocity of flow (u
m
) in case of a flat flame front,
otherwise the flame will be extinguished (in fact this is the way how to determine the burning velocity experimentally). Exam
ple
:
flat laminar flame in porous ceramic burner. Nonuniform transversal velocity profile is manifested by a conical flame front,
and

the flame front cone angle depends upon burning velocity u
L
=u
m

sin


. Example: Bunsen burner.

In case of turbulent flow the flame front is wavy and fluctuating (flamelets, the turbulent flame can be viewed as an ensembl
e o
f
premixed laminar flames). Example: Combustion in a cylinder of spark ignition engine, gas turbines in power plants.

Advantage of the premixed combustion: In a gas turbine operating at a fuel
-
lean condition high temperatures and associated
NOx formation are avoided. Soot formation is also suppressed. Disadvantage: Risk of uncontrolled explosion of premixed
reactants.

flamelets

Turbulent premixed flame

u
L

burning velocity
(=u
m
sin

)

flame front

u
m



Bunsen laminar flame

u
m

porous plate

Flat flame

Bunsen burner fully
premixed. Colour
corresponds to radicals
CH (photo Wikipedia)

R.N. Dahms et al. / Combustion and Flame (2011)

HP1
1

Combustion


Nonpremixed combustion

Fuel and oxidiser enter combustion chamber as separated streams. Mixing
and burning take place simultaneously

(nonpremix is called diffusive burning)
.


Turbulent nonpremixed flame

fuel+flue
gas

air+flue
gas

air

stoichiometric
line

Laminar nonpremixed flame

Flame
front

fuel

air

Burned gases

Nonpremixed flames unlike premixed flames do not propagate.

Flame simply cannot propagate into the fuel side because
there is no oxygen, neither to the oxidising region with lack of fuel. Only burned products diffuse to both sides of the flam
e
front.
Flame front position corresponds to stoichiometric ratio fuel/air


there is also the highest temperature.

In turbulent nonpremixed flames the flamelet concept can be used again. Chemistry of nonpremixed combustion is more
complicated because on one side of the flame front is the fuel
-
rich burning accompanied by formation of sooth and on the
air side a fuel lean combustion occurs. Existence of glowing sooth is manifested by yellow luminiscence.

Examples: Industrial burners, combustion engines.

Laminar coflow
(left)
counter flow (right).

HP1
1

Combustion
-
flame

Previous lecture HP10 was devoted to statics (fuel consumption, overall mass and enthalpy balances,…). Questions like:
What is the distance of the flame front disc levitating above a porous ceramic plate? What is the thickness of the flame fron
t?
What is the burning velocity? Which species are generated inside the flame? Concentration profiles? Temperature profiles?
What is the height and width of the flame?... should be addressed to kinetics.



Flame propagation

An idea about processes related to burning velocity and flame propagation can be
obtained from analysis of the simplest 1D cases. There is a similarity with drying:
process is again controlled by diffusion (now diffusion of fuel and products of
burning), and heat fluxes are consumed at a moving interface (flame front, now the
enthalpy changes of chemical reaction substitute the enthalpy of phase changes).

On the other hand the problem is more difficult, because there are much more
transported species and chemical reactions strongly depend upon temperature.
Transport equations are therefore nonlinear and numerical solution is necessary.

HP1
1

Combustion
-
flame


and numerical solution is necessary
. There are few exceptions: Zeldovic and Frank
-
Kamenetskii (1938) described the flame propagation of a steady premixed planar
flame by a transport equation of fuel (mass fraction

F
(z)) and by transport
equation for enthalpy (temperature T(z)) assuming first
-
order chemical reaction
(fuel

product) and similarity of transport equations (assuming that the diffusion
coefficient D is approximately the same as the temperature diffusion coefficient
a
).
Analysis of these equations reveals the fact that a steady solution (stationary
position of the flame front) exists only as soon as the flow velocity u (=u
F

flame
velocity) is

u
F

T
-
temperature


F

-

fuel

z

10 20 [mm]

)
exp(
1
RT
E
A
D
u
reaction
reaction
diffusion
F





Result: The flame velocity increases when the diffusion
coefficient or the chemical reaction rate increases.

E is activation energy and
A preexponential factor
characterising rate of the
chemical reaction
(Arrhenius term).

HP1
1

Combustion
-
flame

0
)
(
1
2
2





r
r
D
r
r
fuel
dif

R



r

Some information about the flame front can be obtained from the solution of
diffusion equation. For example the flame front position around the burning fuel
particle in air can be estimated from the radial mass fraction profile
s

obtained by
solution of

the

diffusion equation


fueĺ
=
1


ox
=
1

r
R
r
r
R
R
r
ox
fuel











1

)
(
)
(








R
D
m
R
R
D
m
dif
ox
dif
fuel



)
(
…and now: position of front (distance

) is determined by stoichiometry of
chemical reaction (1 kg of fuel consumes
s

kg of oxidiser), therefore
and assuming the same diffusion coefficient in fuel and in air

fuel ox
sm m

/
R s


all fuel is
consumed in the
thin flame front

HP1
1

Combustion
-
flame


Flame size

The size (length) of a nonpremixed burning jet can be estimated from the velocity
of jet and from the diffusion time necessary for mixing in the transversal direction.
Burke and Schumann (1928), Ind.Eng.Chem 20:998


h

r
-
mixing
length

Velocity of circular jets decreases with the distance from nozzle (u
~1/x) however
this approximation assumes a constant velocity at axis and also cylindrical form
of jet, having radius r of fuel nozzle. Mixing length (r) is associated with the
mixing time by theory of penetration depth

t
D
r
dif


and the height of flame h is estimated as the product of the mixing time and axial
velocity

dif
fuel
dif
fuel
fuel
D
V
D
r
r
V
t
r
V
ut
h
2
2
2
2











This is qualitatively correct conclusion stating that the flame length is proportional
to flowrate and indirectly proportional to diffusion coefficient

(Warnatz 1996 documents the diffusion coefficient effect by comparing the height of flame at hydrogen combustion that is
2.5 shorter than the flame height at carbon monoxide combustion)

HP1
1

Combustion
-
flame
-
explosion

Previous examples concern steady combustion processes (time independent
solutions). Ignition and explosions are always time dependent processes.

Explosion
is
characterized by
a
positive feedback, when a temperature increase
releases more heat than could be removed. The simplest demonstration is a
perfectly mixed chamber filled by reaction mixture

(Semenov 1928)


)
(
)
exp(
w
a
reaction
p
T
T
S
RT
E
A
h
dt
dT
c








Q
reaction
(T)

Q
removed
(T)

There exists a critical temperature T
crit

when Q
reaction
(T
crit
)=Q
removed
(T
crit
). For
T
>T
crit

begins the so called run
-
away and temperature increases without any
control (this is explosion).

Remark: the term
deflagration

is used for combustion and flame velocity determined by diffusion (this velocity is rather
small, of the order 1 m/s). The term
detonation

concerns explosions, driven be pressure gradients (propagation
velocities are very high of the order 1000 m/s, because speed of sound is high in burnt gases).

HP1
1

Chemical reactions

HP1
1

Rate of chemical reactions

When discussing combustion statics a detailed knowledge of actual chemical
reactions is not needed


released heat can be evaluated only from the
composition of fuel and burned gas (this composition can be obtained for
example by analyzers of flue gas). This is because according to the Hess’s
law the enthalpy change of a reaction (enthalpy of products minus enthalpy of
reactants) is independent of a specific reaction mechanism:

0 0 0
298,298,298
enthalpy released by combustion
of moles of reactants
( ) ( ) ( )
R
r P f P R f R
P R
h h h

 
    
 
where are stoichiometric coefficients of reactants (R) and
products (P), is the standard enthalpy of formation
corresponding to the formation of compound from elements in their standard
state (at temperature 298 K and atmospheric pressure).


What is substantial: It is sufficient to evaluate only the overall chemical
reaction and it is not necessary to sum up reaction enthalpies of many actual
elementary reactions.




R
P


0
0
,
(
~
)
,

h
f
compound
0
298
Combustion of

Methane

Example: Gas burner operating at atmospheric pressure consumes mass
flowrate of CH
4

0.582 kg/s. Calculate mass flowrate of flue gas. Assuming
temperature of fuel and air 298K and temperature of flue gases 500 K
calculate heating power.

4 2 2
4 2 2 2
0 0 0
,298,298,298
1, 2, 1, 2
( ) 75, ( ) 393, ( ) 242
2 2
4 2 2 2
CH O CO H O
f CH f CO f H O
h h h
CH O CO H O
   
     
        
 

4
0
0 0
298
298,298 4
( )
802
( ) ( ) 75 1 393 2 242 802 / 50
0.016
R
R i f i R
i
CH
h
MJ
h h kJ mol CH h
M kg


               

Flue gas production

1kg CH
4

needs 4 kg O
2
. Mass fraction of oxygen in air is 0.233.

2 2
4 2 2 4 2 4
2 2
0.767
(1 ) (5 4 ) 0.582(5 4 ) 10.57/
0.233
N N
fg CH O N CH O CH
O O
m m m m m m m kg s
 
 
          
Heat removed

4
0
( ) 50 0.582 0.01057 (500 298) 27
R CH fg p fg
Q h m m c T T MW
        
c
p

0.001 MJ/kgK

Air

CH
4

fg
-
flue gas

Q

HP1
1

Remark: methane combustion is described by more
than 200 elementary chemical reactions. Because we
need only overall heat it is sufficient to use only one
overall chemical reaction.

HP1
1

Rate of chemical reactions

Released heat calculated from the equation


corresponds to a certain number of reacting moles of reactants [J/mole].

So that we could calculate the
power

of heat released in a unit volume [W/m
3
] it is
necessary to know number of moles of reactants consumed by each chemical
reaction in unit volume and in unit time, therefore it is necessary to know the time
rate of concentration

0 0 0
298,298,298
( ) ( ) ( )
r P f P R f R
P R
h h h
 
    
 
dt
d
M
dt
d
M
dt
A
d
dt
dc
A
A
A
A
A


1
1
]
[



where c
A

[A] is the molar concentration of specie A [kmol/m
3
], M
A

is the molecular
mass [kg/kmol],

A

is the mass concentration [kg/m
3
] and

A

is the mass fraction of
the specie A in a mixture.

HP1
1

Rate of chemical reactions

If we need to
calculate

the composition inside a flame and the composition of
flue gas or if it is necessary to evaluate temperature distribution inside a
combustion chamber, the rate of chemical reactions must be known.

Reaction rate of a general reaction A+B+…


D+E+… can be expressed as

k

...
]
[
]
[
]
[
b
a
B
A
k
dt
A
d


where [A] is molar concentration of the specie A [kmol/m
3
], exponents a,b,…
are reaction orders with respect to the species A,B,… (only reactants, not
products!) and the sum of all exponents a+b+…=n is the overall reaction
order. K is the rate coefficient of reaction that depends upon temperature and
pressure.

Example:

Combustion of hydrogen with chlorine H
2
+Cl
2

2HCl

Reaction is of the first order with respect to hydrogen and of order 0.5 with
respect to chlorine (overall reaction order n=1.5).

5
.
0
2
2
]
][
[
]
[
Cl
H
k
dt
HCl
d

HP1
1

Rate of chemical reactions

Previous example demonstrates that the reaction orders need not be integers, they
are usually rational, sometimes even negative, numbers. Even worse: the rate
equation of overall chemical reaction has a limited validity, reaction orders and
reaction coefficients depend upon temperature and pressure.

The reason for this complicated behavior is in complexity of chemical reactions. Overall reactions should be decomposed to
several elementary reactions, and in the previous case instead of H
2
+Cl
2

2HCl the following set of four elementary reactions
should be considered
.


1
2
3
4
2
2
2
2
2
2 Cl


k
k
k
k
Cl Cl
Cl
Cl H HCl H
H Cl HCl Cl


  
  
Elementary reactions can be either unimolecular (for example the reaction 1, decomposition of Cl
2

molecule to radicals Cl) or
bimolecular (reactions 2,3,4). Unimolecular elementary reactions are always of the first order (reaction exponent
a
=1) and
bimolecular are of the second order (exponents
a=b
=1). Actual rate of overall chemical reaction is evaluated by simultaneous
solution of the four differential equations with the rate coefficients k
1
,k
2
,k
3
,k
4
. Solution of these differential equations can be
simplified by assuming that the radicals Cl and H are at chemical equilibrium and the time derivatives of their concentration

ca
n
be neglected. Then the differential equations describing the time changes of radicals reduce to algebraic equations and the
concentration of [H] and [Cl] can be eliminated from the set of elementary reactions. Rate equations are therefore expressed
only in terms of [Cl
2
] and [H
2
], giving previously presented result

5
.
0
2
2
]
][
[
]
[
Cl
H
k
dt
HCl
d

Remark: Not all bimolecular reactions are elementary reactions (H
2
+Cl
2

2HCl is an example). Elementary bimolecular
reactions are characterized by presence of highly reactive radicals as reactants.

HP1
1

Elementary reactions

Rate equations of uni
-

and bi
-
molecular elementary reactions are characterized
by reaction orders that are identical with the stoichiometric coefficients 0,1, or 2

2






A
A
B
A
B
A
A
A
kc
dt
dc
D
A
A
c
kc
dt
dc
dt
dc
D
B
A
kc
dt
dc
D
A












Kinetics of a general system of R elementary reactions (r=1,2,…,R) of
S species (s=1,2,…,S) can be described as







S
s
s
R
r
ri
r
i
rs
c
k
dt
dc
1
)
,
0
max(
1


i=1,2,…,S

Explanation will be shown in the following
example (H
2
+Cl
2

2HCl)


stoichiometric coefficient of
the specie s in the reaction
r (negative for reactant,
positive for product)

HP1
1

Elementary reactions

There is 5 species S
1
(Cl
2
), S
2
(H
2
), S
3
(Cl), S
4
(H), S
5
(HCl) and four reactions

1
2
3
4
1 3
3 1
3 2 5 4
4 1 5 3
2
2
k
k
k
k
S S
S S
S S S S
S S S S


  
  
max(0,)
1
1
1
1
1 0 0 0 0 0 0 2 0 0 1 0 0 1 0
1 1 2 3 4 5 2 1 2 3 4 5 4 1 2 3 4 5
2
1 1 2 3 4 1 4
( 1)( ) (1)( ) ( 1)( )
rs
S
R
r r s
r
s
dc
k c
dt
k c c c c c k c c c c c k c c c c c
k c k c k c c





 
     
   


1
1
1
0
1
1
1
1
1
0
0
0
2
0
1
0
0
2
0
1
45
44
43
42
41
35
34
33
32
31
25
24
23
22
21
15
14
13
12
11














































Let us apply rate equations for the concentration of chlorine Cl
2

(i=1)

matrix of stoichiometric
coefficients

I hope that it is comprehensible… all other rate equations are similar.

HP1
1

Temperature dependence

Reaction
constants k depend upon temperature according to Arrhenius law

a
E
RT
k Ae


Activation energy of
chem.reaction
J/mol

Preexpontial factor
(the same unit as k)

Relative amount of molecules having kinetic energy greater than E
a

is given by
Maxwell distribution
of energies

/
E RT
e

The greater is temperature
the more molecules have
energy greater than E

Collision theory: only the molecules having energy greater than E
a

react at mutual
collision.

Reactions having high activation energy proceed at high temperatures and vice
versa (thermal NOx reaction with E
a
~300 kJ/mol starts at about 1500 K, while
prompt NOx reaction with
E
a
~100 kJ/mol is important at relatively low temperature
about 1000 K)

HP1
1

Simplified reaction mechanism

Reality is such that even the simple looking reaction like hydrogen combustion
2H
2
+O
2

2H
2
O requires about 40 elementary reactions for a satisfactory
reaction mechanism. Combustion of hydrocarbons is even more complicated
and is represented by hundreds and thousands of elementary reactions
(
oxidation of the simplest hydrocarbon
CH
4

is described by
277 reactions among
49 chemical species). Simultaneous solution of thousands differential equations
is too expensive for simulation of practical systems with 3 dimensional
inhomogeneous turbulent flows. There exist several techniques how to reduce
the complicated system, generally based upon exclusion of reactions that are
very fast and are in partial equilibrium (QSSA Quasi
-
Steady
-
State
-
Assumption,
ILDM Intrinsic Low Dimensional Method, see Peters 2000). Methane oxidation
can be in this way reduced to global mechanism consisting in only four or even
two steps

2
2
2
2
4
2
1
2
2
3
CO
O
CO
O
H
CO
O
CH





HP1
1

Simplified reaction mechanisms are applied also for description of NOx formation

(important for design of burners from the point of view of emissions). While the
conce
n
trations of CO
2
, H
2
O can be predicted by assuming equilibrium, NO
reactions are slow and far from equilibrium at short residence times in burning
zone (some ms, time scales of reactions are thousand times longer).


Thermal NOx (Zeldovic mechanism)





Prompt NOx (Fenimore)




Fuel NOx

Simplified reaction mechanism

H
NO
OH
N
O
NO
O
N
N
NO
N
O
k
k
k















3
2
1
2
2
Source of NO from atmospheric nitrogen formed by
the following
reactions
:

The first reaction (attack of oxygen radical to N
2
) has very high activation
energy (319 kJ/mol) and proceeds only at very high temperatures.

Prompt NO is formed also from atmospheric N
2

but at the flame front from
radicals CH (fuel rich flame)

NO
N
HCN
N
CH





...
2
Fuel
-
bound nitrogen is converted into NO mainly in coal
combustion (there is approximately 1% chemically bound
N in coal)

O
N
NO
N
H
NO
OH
N






2
CFD
-

TURBULENT FLAMES

Habisreuter P. et al.
Mathematical modeling of turbulent swirling flames

in Cox M.et al: Combustion Residues: Current, Novel and Renewable
Applications. J.Wiley 2008


CFD
-

THERMOCHEMISTRY

Kubota, N. (2007) Combustion Wave Propagation, in Propellants and
Explosives:
Thermochemical Aspects of Combustion
, Wiley Weinheim

HP1
1

CFD mod
e
l
s


CFD
-

TURBULENT COMBUSTION

Peters N.: Turbulent combustion. Cambridge Univ.Press, 2000


Alexander

HP1
1

Flow modelling


Navier Stokes

( ) 0
u
t



 

( )
m
Du
p u s
Dt
 
    
Combustion modelling is always based upon
the
solution of Navier Stokes
equations for vector of velocity field (even if it is combustion of liquid or solid
fuels). There are 4 unknowns: 3 components of velocity u
x
,u
y
,u
z

and pressure.
Therefore 4 partial dif
f
erential equations must be solved

simultaneously

Continuity equation
, stating
that the accumulation of mass is the
divergence of mass flux

Momentum balance
, stating
that the mass x acceleration equals
the sum of forces acting upon a fluid
element (pressure forces + viscous
forces + volumetric forces)

Solution of these equations (velocity field) describes convective transport of
fuel, oxidiser as well
as
products of burning inside a combustion chamber.
These equations are quite general, describe incompressible/compressible
fluids, laminar as well as turbulent flow regime, but…

HP1
1

Flow modelling


Navier Stokes


laminar as well as turbulent flow
s
, but

in the turbulent regime it would be
necessary to model all details of chaotic motion of turbulent eddies. Large and
slowly rotating energetic eddies decay to smaller and smaller ones down to tiny
Kolmogorov eddies so small that viscous effects prevail and dissipate all kinetic
energy to heat. Extremely fine mesh is necessary for resolution of Kolmogorov
eddies and number of nodes N
nodes

rapidly increases with the Reynolds number.
Also the the range of time scales is very broad due to high frequency of turbulent
fluctuations of small eddies (10 kHz).

9/4
3/4
_
Re
Re
nodes
time steps
N
N


Very small time step of Direct Numerical Simulation and number of nodes 10
9

make standard numerical methods of finite volumes or elements prohibitively
expensive for the solution of engineering problems (including combustion).

HP1
1

The most frequently used technique for modeling
of
turbulent flow
s

consists in
averaging of fluctuating quantities (velocities, density, pressure, temperature and
concentrations). Filtration of time fluctuations (substituting actual noisy time
course

at a point x,y,z

by average) is the principle of
RANS

(Reynolds Averaging
Navier Stokes) family of numerical methods while
the
spatial fluctuations filtering
is a characteristic of
LES

(Large Eddy Simulation) class.

It is assumed that filtering the time fluctuations smooth out the spatial fluctuations
and vice versa (because turbulent fluctuations are correlated).

/2
/2
( ) (,)/
t
t
u x u x d
 


 

Δ

t

u

'
u
u
u


Turbulent flow RANS
/LES

Mean
value

Temporary fluctuation

HP1
1

Applying
time
average operator to each term of continuity equation results in
equation which is exactly the same as the original equation just only instead of
actual velocities the mean (filtered) velocities are substituted.
T
ime averaging
applied to the momentum equation
also
results to
a similar

equation for the mean
velocities but a new term, Reynolds turbulent stresses

t
, appear
s

on the right
hand side

Turbulent flow RANS

( )
t m
Du
p u s
Dt
  
     
( )
ef m
Du
p u s
Dt
 
    
Boussinesque hypothesis suggests that the Reynolds turbulent stresses are
again proportional to the gradient of mean velocity as in the case of laminar
flow but with a different proportionality coefficient
-

turbulent viscosity

Resulting equation can be solved by the same numerical techniques as in the
laminar flows, but using variable viscosity

ef
(k,

)
.

HP1
1

Turbulent flow RANS k
-


2
t
k
C

 


Effective viscosity can be expressed as the sum of actual (molecular) viscosity


and
turbulent viscosity

t

representing momentum transport by turbulent eddies. This
contribution (

t
) depends upon the magnitude of velocity fluctuation (or square of this
magnitude, expressed by
the kinetic energy of turbulence
k

[m
2
/s
2
]) and also upon
the rate of how the kinetic energy of turbulent eddies is converted to smaller eddies
and in the final instance to heat. This second characteristic is the dissipation of
energy of turbulence



2
/s
3
]
(the unit of


follows from

the unit of k, because


is the change of k per second).

Turbulent viscosity can be expressed by the model

where C


is a constant (usually 0.09). There is no mystery in the correlation: it can
be derived easily by dimensional analysis (unit of viscosity is Pa.s) assuming
power law formula and comparing exponents at kg, m, s

t
k
  
 

2 2 2 2 3
2 3 3 2 2 3
[.] [ ][ ][ ]
Ns kg kg m m kg m
Pa s
m ms m s s s
      
    
 

 
   
 
 
1
1 2 2 3
1 2 3

 
 

   
 
1
1
2




 

HP1
1

Turbulent flow RANS k
-


1
2
2
1
( ) (( ) ) ( 2:)
t
t
u C E E C
t k
f f
 


 
   




      

( ) (( ) ) 2:
t
t
k
k
ku k E E
t


 




      

2 2
':'''
ij ij
e e e e
 

 
 
2 2 2
1 1
''(''')
2 2
i i
k u u u v w
   

Transport equation for the kinetic energy of
turbulence
k

[m
2
/s
2
]


Transport equation for dissipation of energy
of turbulence



2
/s
3
]

Production term (generator of
turbulence)

Dissipation
term

'
'
1
'( )
2
j
i
ij
i j
u
u
e
x x


 
 
1
( ( ) )
2
T
E u u
   
It is assumed that k and


are scalar properties, that can be calculated from

HP1
1

Turbulent flow RANS

The choice k
-


and the two transport equations for k and


is not by far the only way
how to characterize a state of turbulence.

Two equation models calculate turbulent viscosity from the pairs of turbulent characteristics
k
-

, or
k
-

,
or

k
-
l
(Rotta 1986)

k (kinetic energy of turbulent fluctuations) [m
2
/s
2
]



(dissipation of kinetic energy) [m
2
/s
3
]


(specific dissipation energy) [1/s]

l (length scale)


[m]

2
t
t
k
k
C

 

 



Wilcox

(1998)

k
-
omega model
, Kolmogorov (1942)

Jones,Launder (1972),
Launder Spalding

(1974)

k
-
epsilon model

These pairs of turbulent characteristics can describe not only the turbulent viscosity
but also transport coefficients (diffusivity, thermal conductivity) and the rate of
micromixing R
m

(kg of component mixed in unit volume per second) that controls
the rate of combustion, e.g.






C
R
k
C
R
m
m



as follows from
dimensional analysis (the
same way as turbulent
viscosity)

HP1
1

Turbulent flow RANS

All the presented CFD
RANS models (and many
others), are available in
commercial programs like
Fluent.

Knowing velocity field and
scalar characteristics of
turbulence it is possible to
analyzed transport of
individual species,
chemical reactions, and
temperatures…


i



mass fraction of specie i in mixture
[kg of i]/[kg of mixture]





i

mass concentration of specie [kg of i]/[m
3
]


Mass balance of species

(
for each specie one transport equation)

( ) ( ) ( )
i i i i
i
u
t
S
  

   

Rate of production
of specie i [kg/m
3
s]

Production of species is controlled by


Diffusion of reactants (micromixing)


t
diffusion

(diffusion time constant)


Chemistry (rate equation for perfectly mixed reactants)


t
reaction

(reaction constant)

Damkohler number

diffusion
reaction
t
Da
t

exp( )
a
i i
i i
E
S A
RT
S C
k
 


 

Da<<1

Da>>1

Reaction controlled by kinetics (Arrhenius)

Turbulent diffusion controlled combustion

Because only micromixed
reactants can react

HP1
1

Combustion


mass balances

HP1
1

Combustion


mass balances

max(0,)
1
1
rs
S
R
i
i i r ri s
r
s
d
S M k c
dt






 


The source term S
i

[kg/(m
3
s)] is the sum of contributions of all chemical reactions
where the specie i exists either as a reactant or a product (when stoichiometric
coefficient

ri

of specie i in reaction r is nonzero)


S
-
is number of species (i=1,2,…,S) and R is number of chemical reactions.

Enthalpy balance is written for mixture of all species (result
-
temperature field)

( ) ( ) ( )
h
h h
S
u h
t
 

   

Sum of all reaction enthalpies of all reactions

,
enthalpy of formation
/
h i f i i
i
S S h M
 

It holds only for reaction without
phase changes h
~ c
p
T

Energy transport must be solved together with the fluid flow equations (usually
using turbulent models, k
-

, RSM,…). Special attention must be paid to
radiative energy transport (not discussed here, see e.g.
P1
-
model
,
DTRM
-
discrete transfer radiation,…). For modeling of chemistry and transport of
species there exist many different methods but only one
-

mixture fraction
method for the non
-
premixed combustion will be discussed in more details.

HP1
1

Combustion
-

enthalpy

3 3
,
kg/kmol
enthalpy of
W/m kg/(m s)
formation J/kmol
/
h i f i i
i
S S h M
 

HP1
1

Combustion
-

enthalpy

The production term S
h

[W/m
3
] is the sum of reaction enthalpies over all exo
-

or
endothermic reactions.

It can be calculated from a simplified reaction scheme (from reduced number of
overall reactions, reaction enthalpies and reaction progresses of these reactions)
but also from the previously calculated S
i

values (kg of specie i generated in unit
volume and unit time due to all chemical reactions) and corresponding enthalpies
of formation

HP1
1

Mixture Fraction method

Bingham

Let us consider non
-
premixed combusion, and fast chemical reactions
(paraphrased as “
What is mixed is burned

(or at equilibrium)”)

m
fuel

m
oxidiser

Flue gases

Mass balance of fuel



Mass balance of oxidant

( ) ( ) ( )
fuel fuel fuel fuel
u S
t
  

   

( ) ( ) ( )
ox ox ox ox
u S
t
  

   

Calculation of fuel and oxidiser
consumption is the most
difficult part. Mixture fraction
method is the way, how to
avoid it

3
of produced fuel
[ ]
kg
s m

Mass fraction of
oxidiser (e.g.air)

Mass fraction of fuel
(e.g.methane)

HP1
1

Mixture Fraction method


Stoichiometry

1 kg of fuel +
s

kg of oxidiser


(1+s) kg of product

a
nd subtracting previous equations

( ) ( ) ( )
( ) ( ) ( )
fuel fuel fuel fuel
ox ox ox ox
s s s s
u S
t
u S
t
  
  

    


    

( ) ( ) ( )
fuel ox
s
u S S
t
 

      

Introducing new variable

fuel ox
s
 
  
This term is ZERO
due to stoichimetry

HP1
1

Mixture Fraction method

-
S
fuel

is kg of fuel consumed in unit
volume per one second (at a given place
x,y,z and time t).

Mixture fraction
f

is defined as linear function of


normalized in such a way that
f=0 at oxidising stream and f=1 in the fuel stream

Resulting

transport equation for the mixture fraction
f

is
without any source term

( ) ( ) ( )
f fu f
t
 

  

,0
0
1 0,1,0
fuel ox ox
fuel ox
s
f
s
  
 
 

 
  
Mixture fraction is
a
property that is CONSERVED, only dispersed and transported
by convection.
f

can be interpreted as a concentration of a key element (for
example carbon). And because it was assumed that „what is mixed is burned“ the
information about the carbon concentration at a place x,y,z bears information
about all other participating species
, see next slide…



ox

is the mass fraction of
oxidiser at an arbitrary point
x,y,z, while

ox,0

at inlet (at the
stream 0)

HP1
1

Mixture Fraction method

Knowing
f

we can calculate mass fraction of fuel and oxidiser at any place x,y,z

,0
,1,0
ox
st
fuel ox
oichio
f
s

 


For example the mass fraction of fuel is calculated as


ox,1
(fuel rich region, oxidiser is consumed =
0)

1
1
stoichio
stoichio
stoichio
fuel fuel
f
f
f
f
f

 

  

(fuel lean region)

0 0
stoichio
fuel
f
f

  
The concept can be generalized assuming that chemical reactions are at
equlibrium

f



i


mass fraction of species is calculated from equilibrium constants





(evaluated from Gibbs energies)

At the point x,y,z where f=f
stoichio

are all reactants consumed
(therefore

ox
=

fuel
=0)

HP1
1

Mixture Fraction method

Equilibrium depends upon concentration of
the
key component (upon
f
) and
temperature. Mixture fraction
f

undergoes turbulent fluctuations and these
fluctuations are characterized by probability density function p(f). Mean value of
mass fraction, for example
the
mass fraction of fuel is to
be
calculated from this
distribution

1
0
( ) ( )
fuel fuel
f p f df
 


Mass fraction corresponding to an
arbitrary value of mixture fraction is
calculated from equlibrium constant


Probability density function, defined in
terms of
mean and variance of f

0 fmean 1

p

Frequently used


distribution

1 1
1
1 1
0
(1 )
( )
(1 )
p q
p q
f f
p f
f f df
 
 




Variance of f is calculated from another transport equation

2 2 2 2 2
(') (') (') ( )'
f g t d
f f u f C f C f
t k

   

      

HP1
1

Mixture Fraction method

Final remark: In the case, that

fuel
is a linear function of f, the mean value of
mass fraction


fuel

can be evaluated directly from the mean value of f (and it is
not necessary to identify probability density function p(f), that is to solve the
transport equation for variation of f). Unfortunately the relationship


fuel
(f) is
usually highly nonlinear.

1
0
( ) ( ) ( )
fuel fuel fuel
f p f df f
  
 

HP1
1

Mixture Fraction method

m
fuel

1
| | ( )
2
D D
du
m F
dt
F c A u v u v


  
Lagrangian method
: trajectories, heating and evaporation of droplets
injected from a nozzle are calculated.

Sum of all forces acting to liquid
droplet moving in continuous
fluid (fluid velocity
v

is calculated
by solution of NS equations)

Relative velocity
(fluid
-
particle)

Drag force

Drag coefficient c
D
depends upon Reynolds number

0.687
5
24 3
(1 Re) Re 5 Oseen
Re 16
24
(1 0.15Re ) Re 800 Schiller Nauman
Re
0.4 1000< Re 3.1
0 Newton
D
D
D
c
c
c
  
  
 
1 10
4

10
5

Re

c
D

Newton’s region
c
D
=0.44

Effect of cloud
(

c

volume
fraction of dispersed phase
-
gas)

3.7
0
/
D D c
c c


HP1
1

Combustion


liquid fuel

Evaporation of fuel droplet

Diffusion from droplet surface to gas:

2
2
05 0.33
( ) ( )
6
2 0.6Re
p g dif fg s
dm d D
D Sh D m m
dt dt
Sh Sc

  
  
 
Mass fraction of
fuel at surface

Sherwood
number

Schmidt number =

/D
dif

Ranz Marshall correlation for mass transport

HP1
1

Combustion


liquid fuel

dif
D
D
Sh


CFD examples

Vít Kermes, Petr Belohradský, Jaroslav Oral, Petr Stehlík:
Testing of gas and liquid fuel burners for power and process
industries. Energy, Volume 33, Issue 10, October 2008,
Pages 1551
-
1561



Lempicka

HP1
1

CFD


gas burner

HP1
1

CFD


gas burner

Zeldovic

N
O
x


t
hermal NOx

HP1
1

CFD


gas burner

HP1
1

CFD


gas burner


Low computational cost CFD analysis of thermoacoustic oscillations

Applied Thermal Engineering
,
Volume 30, Issues 6
-
7
,
May 2010
,
Pages 544
-
552

Andrea Toffolo, Massimo Masi, Andrea Lazzaretto


Timing of the pressure fluctuation due to the acoustic mode at 36 Hz and of the heat released by
the fuel injected through the main radial holes (if the heat is released

in the gray zones, the
necessary condition stated by Rayleigh criterion is satisfied and thermoacoustic oscillations at this
frequency are likely to grow).

HP1
1

CFD
thermoacoustic burner

HP
1
1

EXAM

HP
1
1

Combustion

chemistry and CFD

What is important
(at least for exam)

HP
1
1

0 0 0
298,298,298
0
products (with positive enthalpy of forma
tion
stoichiometric coefs. ) at standard press
ure
and temperature 298K
( ) ( ) ( )
R
P
r P f P R f R
P R
h h h


 

    
 
Heat released by reaction of reactants with stoichiometric coefficients

R

max(0,)
1
1
1
rs
S
R
i i
r ri s
r
s
i
dc d
k c
dt M dt






 


stoichiometric coefficient of
the specie s in the reaction
r (negative for reactant,
positive for product)

Kinetics of R
-

elementary chemical reactions (c
i

is molar concentration of specie i).
This rate equation holds only for elementary reactions (mono and bimolecular reactions) and not for overall reactions.

H
NO
OH
N
O
NO
O
N
N
NO
N
O
k
k
k















3
2
1
2
2
Example: complicated reaction mechanism of NO formation from atmospheric nitrogen can be simplified to only 3
elementary bimolecular reactions (Zeldovic mechanism of NOx formation)

OH
N
O
N
N
O
NO
c
c
k
c
c
k
c
c
k
dt
dc
3
2
1
2
2



reaction constant depends upon
temperature (Arrhenius term
~exp(
-
E
r

/RT), where E
r

is
activation energy

of reaction r)

What is important
(at least for exam)

HP
1
1

reaction heat corresponding
to all chemical reactions
( ) ( )
h
S
hu h

   
production of specie i
convective transport diffusion by turbule
nt eddies
( ) ( )
i i i i
u S
 
     
Concentration of species c
i

(or mass fractions

i
)

in the combustion chamber
follows from the mass balance (transport equation for specie i)

max(0,)
1
1
rs
S
R
i
i i r ri s
r
s
d
S M k c
dt






 


The source term S
i

corresponds to all chemical reactions where the specie i
exists either as a reactant or a product (when stoichiometric coefficient

ri

is
nonzero)


Temperature field follows from enthalpy balance where all chemical equations
(and their reaction enthalpies) should be included

What is important
(at least for exam)

HP
1
1

CFD (Computer Fluid Dynamics) solution is based upon numerical solution of
previous transport equations for concentrations (or

i
) and temperatures (or
enthalpies h) together with the Navier Stokes equations and continuity equation
for velocity and pressure. In case of turbulent flow it is necessary to solve also
additional transport equation for turbulent characteristics (for example k and

).

This complicated system of partial differential equations can be sometimes
simplified: In the case of non
-
premixed combustion (mixing of fuel and air
proceeds not in front of but inside the combustion chamber) and assuming
practically instantaneous chemical reactions (concept “what is mixed is burned”),
the transport equations for all species can be substituted by only one transport
equation for “mixture fraction f” (mass fraction of a key element, usually carbon).

convection of f dispersion of f
( ) ( )
fu f

  
There is no source term in the transport equation for the mixture fraction f
because elements on contrary to species are conserved. Concentrations of
individual species

i
(x,y,z) can be derived from the calculated f(x,y,z) and from
stoichiometry (or assuming immediately established chemical equilibrium).