Versteeg H., Malalasekra W. An Introduction to Computational Fluid Dynamics: The Finite Volume Method

Подождите немного. Документ загружается.

12.11 REACTION RATE 357

For example, consider the simple reaction

+ O

→ 2H

O (12.57)

For this reaction the total number of species involved is N = 3; the species

involved are M

= H

, M

= O

and M

= H

Considering the left hand side of equation (12.57) we have

Stochiometric coefﬁcient

′

for H

= 2

Stochiometric coefﬁcient

′

for O

= 1

Stochiometric coefﬁcient

′

for H

O = 0, no H

O on the left hand

side of the equation

Considering the right hand side of equation (12.57) we have

Stochiometric coefﬁcient

″

for H

= 0, no H

on the right hand

side of the equation

Stochiometric coefﬁcient

″

for O

= 0, no O

on the right hand

side of the equation

Stochiometric coefﬁcient

″

for H

O = 2

The mass conservation (12.56) can be veriﬁed for this equation as

Left hand side (2 × 2 + 1 × 32) = (2 × 18) Right hand side

The law of mass action (see Kuo, 2005; Turns, 2000), which has been

conﬁrmed by experimental observations, states that the rate of consumption

of a chemical species in a reaction is proportional to the (mathematical) prod-

uct of the concentration of the reacting species, each concentration being

raised to the powers corresponding to the stoichiometric coefﬁcients of the

reactants.

The reaction rate RR (units kmol/m

.s) for the reaction (12.55) is given

by (see Kuo, 2005)

RR ==k

)

v ′

(12.58)

Since species M

appear on both sides of equation (12.55), the net progress

rate of reaction j is given by the production of the species M

minus the

destruction of the species M

==(

″

−

′

)

v ′

(12.59)

Here C denotes the molar concentration of the species in units kmol/m

(or gmol/cm

in CGS units) and the proportionality constant k

is called the

speciﬁc reaction rate constant. It should be noted here that the subscript k

used to identify species is different from this k. The reaction rate constant k

is independent of the concentration and is usually expressed as

= AT

exp − (12.60)

This is known as the Arrhenius law. In this expression A is called the

pre-exponential constant. Parameter

is a temperature exponent and E

the activation energy. These constants are speciﬁc to a given reaction. R

is the universal molar gas constant (= 8.314 kJ/kmol.K) and T is absolute

∏

k=1

∏

k=1

d(C

product

)

ANIN_C12.qxd 29/12/2006 04:44PM Page 357

358 CHAPTER 12 CFD MODELLING OF COMBUSTION

temperature. A note about units is required here. Consider the simple

forward-only reaction

aA + bB → cC + dD (12.61)

The progress rates of reactant consumption and product formation for this

reaction can be evaluated using equation (12.59) as

==−ak

)

(12.62)

==−bk

)

(12.63)

==ck

)

(12.64)

==dk

)

(12.65)

Note that equation (12.61) is a forward-only reaction. When a backward

reaction is also present similar expressions can be written using a backward

reaction rate constant k

. The units of C

and C

are kmol/m

. The units

of the progression rate are kmol/m

.s. Because of the involvement of

exponents in equations (12.59) it can be seen that the units of k

, which is

given by equation (12.60), depend on the stoichiometry of the reaction. For

the reaction (12.61), the sum of the coefﬁcients a and b is called the reaction

order n = a + b. E

has units kJ/kmol. RT has the same units. The exponen-

tial term −E

/RT does not contribute to units. The factor T

has units

depending on the value of

(see Turns, 2000). Therefore for this particular

reaction the constant A has units (m

/kmol)

n−1

/(K

−

s). The units of k

are

/kmol)

n−1

/s. In most cases

is zero, in which case the units required for

are the same as those for A (see e.g. Henson and Salimian, 1984).

For a simple reversible reaction of the form

aA + bB bcg cC + dD (12.66)

the progress rate of consumption of the reactant species A in equation

(12.66) is

==−a[k

)

− k

)

] (12.67)

Similarly for a large scheme of reactions of the form

′

bcg

″

for j = 1, 2,..., m (12.68)

the most general form of the reaction progress rate for species k by reaction

j is given by

= (

″

−

′

) . k

)

′

− k

)

″

(12.69)

∏

k=1

∏

k=1

∑

k=1

∑

k=1

ANIN_C12.qxd 29/12/2006 04:44PM Page 358

12.11 REACTION RATE 359

where k

is the forward reaction rate constant and k

is the backward reaction

rate constant for the reaction j. This kind of compact notation is very useful

in handling large systems of kinetic schemes containing many equations

(more than 100 or so) and is used in common kinetic solution algorithms

such as CHEMKIN (Kee et al., 1996).

In complex chemical schemes involving m equations the total rate of

production rate L

of a certain species k is the sum of individual rates of each

equation producing the species k, i.e.

= L

(12.70)

where L

is the total reaction rate for the species k and L

is the reaction

rate for species k for the speciﬁc reaction j.

These reaction rates are used in the transport equations for each of the

species in the solution of combustion problems. When reaction rate expres-

sions are written for each species participating in the mechanism this yields

a system of ﬁrst-order differential equations describing the progress of a

chemical system with time. Given the initial conditions (initial concentration

of the reacting species and temperature), additional equations for conservation

of mass, momentum and energy can be written, and the system of differential

equations can be numerically integrated to obtain a solution. In the solution

of such systems, when one or more variables change rapidly in comparison

with other variables, the resulting system is said to be stiff. Chemical kinetic

solution packages use special routines to deal with such stiff systems (see Kee

et al., 1996; Radhakrishnan and Pratt, 1988).

The reaction rates are used as source terms (see later) in transport equa-

tions for species mass fractions, which have the units kg/m

.s. To convert to

required units we multiply L

, which is in kmol/m

.s, by the molecular

weight (MW). After converting to mass rate units the usual reaction rate

source term used in species transport equations may be written as

= (MW )

(12.71)

where L

is obtained from (12.70) and (12.69) using details of kinetic data

(values of constants A,

, E

etc.) from Baulch et al. (1994), Turns (2000),

Smith et al. (2003), Kuo (2005), Gardiner (1984) and many other publica-

tions cited later. In some of the early combustion literature tabulated kinetic

data have been tabulated in CGS units. Care must be given in using these

datasets, and appropriate SI unit conversion is required when incorporating

reaction rates.

A further point to note is the link between forward and backward reaction

rate constants. When the forward reaction rate constant (k

) is known, the

backward rate constant (k

) for the jth reaction can be calculated through the

equilibrium constant introduced in section 12.8:

(12.72)

kkC

KRTp

fj bj M

kj kj

( ) ( /)

()

∑

∑∑

−

∏

νν

″′

″

′

′″

products

reactants

∑

j =1

ANIN_C12.qxd 29/12/2006 04:44PM Page 359

360 CHAPTER 12 CFD MODELLING OF COMBUSTION

Figure 12.1 Schematic

representation of the reaction

rate expression for a single-

step reaction

For example, one of the reactions known as the Zel’dovich mechanism for

thermal NO formation (discussed later) is

NO + O → N + O

The rate coefﬁcient (from Turns, 2000) for the forward reaction is

= 3.80 × 10

1.0

exp(−20820/T) cm

/gmol.s

At 2300 K, the value is k

= 1.029 × 10

/gmol.s. For this reaction

(∑

′

−∑

″

) = 0. The K

-value for the above reaction at 2300 K is

= 1.94 × 10

−4

. Using equation (12.72) the backward reaction rate

coefﬁcient is obtained as

= k

==5.278 × 10

/(gmol s)

Further details of chemical kinetics, the fundamental development of reac-

tion rate expressions, reaction kinetic data and how they are used in packages

like CHEMKIN are available in Turns (2000) and Kuo (2005).

To illustrate the behaviour of the reaction rate expression for a simple

reaction of the form (12.61) the rate of consumption of species A can be

written as

= A′[Y

]

exp(−E/R

T) (12.73)

where A′=(MW )

a AT

It should be noted that constant A′ contains density and molecular weight

terms, which arise from the conversion from molar concentrations (12.69) to

mass fractions (as used in CFD computations). According to this expression,

it can be seen that the rate increases exponentially with temperature. As

the reaction progresses, the mass fractions of reactants decrease as they get

consumed. Therefore the reaction rate decreases as species A and B are con-

sumed. The typical behaviour of the rate expression is shown in Figure 12.1.

(MW )

1.024 × 10

1.94 × 10

− 4

ANIN_C12.qxd 29/12/2006 04:44PM Page 360

12.13 REDUCED MECHANISMS 361

The reaction rate increases almost exponentially with temperature (T) at

ﬁrst because of the dominant role of the term exp(−E/RT). As T increases

towards the ﬁnal temperature T

burnt

, the reaction rate drops in spite of the

exponential term because [Y

]

decreases rapidly as reactants A and B

species are consumed in the reaction. This demonstrates that the availability

of reactants (concentrations) and temperature both play major roles in deter-

mining the reaction rate.

There are many types of reaction: reactions of various orders (ﬁrst-order,

second-order etc.), consecutive reactions, competitive reactions, opposing

reactions, chain reactions, chain branching reactions etc. For oxidation of

various fuels detailed reaction mechanisms and appropriate rate constants

for those reactions are available in the literature. For example, see Gardiner

(1984), Drake and Blint (1988), Dryer (1991), Smooke (1991), Peters (1993),

Turns (2000), Seshadri and Williams (1994), Warnatz et al. (2001), the

GRI 3.0 mechanism (Smith et al., 2003 and references therein), and also

the San Diego mechanism (2003, http://maemail.ucsd.edu/combustion/

cermech/).

Detailed mechanisms for methanol and the CO/H

system are

available in Dryer (1991). This reference also gives one-step, two-step and

four-step mechanisms for a variety of hydrocarbon fuels. Peters (1993) is

a valuable source book which has contributions from many other authors

and contains details of multi-step and reduced mechanisms for a range

of fuels. A detailed mechanism for methane combustion, NO formation,

single-step and multi-step reaction mechanisms for common hydrocarbons

are also available in Turns (2000). A mechanism containing 46 steps for

methane combustion has been reported in Smooke et al. (1986). A very well-

known and widely used mechanism for methane combustion is GRI 2.11

(Bowman et al., 1996). The most recent version of this mechanism is GRI 3.0

(Smith et al., 2003), which includes a detailed mechanism for NO formation.

However, it has been observed in some test cases that the previous version

GRI 2.11 appears to give better NO predictions (at the time of writing this

text) than the latest version GRI 3.0 (see Kim and Huh, 2002). Further

details of detailed mechanisms can be found in Seshadri and Williams

(1994), Warnatz et al. (2001), and also at the Lawrence Livermore Laboratory

mechanisms link (http://www-cms.llnl.gov/combustion/combustion2.html).

It should be noted that rate constants of various mechanisms have been

derived on the basis of computational estimations and comparisons with

measurable combustion properties such as ﬂame speeds, temperatures and

species mass fractions obtained in controlled experiments. For large systems

of chemical equations dedicated computer programs are required to solve

chemical kinetic problems. CHEMKIN (Kee et al., 1996), for example, is a

widely used software package for such problems, and many commercial CFD

codes allow the incorporation of its information relating to chemical reactions.

The computational cost of chemical kinetics evaluations and the associated

species transport equations is substantial, and grows rapidly as the reaction

mechanisms become more elaborate and detailed. Consequently, efforts have

Detailed

mechanisms

12.12

Reduced

mechanisms

12.13

ANIN_C12.qxd 29/12/2006 04:44PM Page 361

been made to develop more practical reaction schemes involving fewer reac-

tions to represent combustion of basic fuels. These simpliﬁed schemes are

called reduced mechanisms and consist of a few key equations to predict

major and important minor species (see e.g. Dryer, 1991; Seshadri and

Williams, 1994). For example, the detailed mechanism for H

oxidation given

in Conaire et al. (2004) consists of 19 reactions. Massias et al. (1999) give a

ﬁve-step reduced mechanism for H

combustion including NO formation

(Table 12.1).

362 CHAPTER 12 CFD MODELLING OF COMBUSTION

Table 12.1

bcg 2H

+ O bcg H

3 4H + O

bcg 2H

+ O

bcg 2NO

+ O bcg N

The detailed chemical–kinetic mechanism for methane combustion

shown in Seshadri and Williams (1994) consists of 39 elementary reactions

involving 17 chemical species including N

. A ﬁve-step reduced mechanism

by Hewson and Bollig (1996) for CH

combustion is shown in Table 12.2.

The same reference gives a 52-step, 13-species reaction detailed mechanism

for nitrogen chemistry (species NO, NH

, HCN etc.) and a reduced six-step

reduced mechanism for nitrogen chemistry.

Table 12.2

1CH

+ 2H + H

O bcg 4H

+ CO

O + CO bcg H

+ CO

32Hbcg H

+ 3H

bcg 2H

O + 2H

+ 2CO bcg C

+ O

Rate constants and further data for the reduced mechanisms shown

in Tables 12.1 and 12.2 can be obtained from the above cited references.

Extensive studies have been made of most practically important combustion

reactions, resulting in a vast body of publications with proposals for detailed

and reduced mechanisms. Some detailed mechanisms involve many equa-

tions – the above-mentioned GRI 3.0 mechanism for natural gas combus-

tion including NO

chemistry, for example, consists of 325 reactions and

53 species. The accuracy of reduced mechanisms has been compared with

detailed mechanisms for many practically important combustion reactions

(see e.g. Barlow et al., 2001; Massias et al., 1999). The choice of a reduced

mechanism is always a trade-off between computational cost and accuracy

required. The user should carefully select a suitable mechanism depend-

ing on the particular application and available resources. When results are

interpreted it should also be noted that some reduced mechanisms have

been tuned for speciﬁc temperature ranges and speciﬁc types of combustion

processes.

ANIN_C12.qxd 29/12/2006 04:44PM Page 362

12.14 GOVERNING EQUATIONS FOR COMBUSTING FLOWS 363

So far we have discussed thermodynamics, chemical equilibrium and chem-

ical kinetics without any reference to ﬂow conditions. In many combustion

situations ﬂuid ﬂow is an integral part of the combustion process. In non-

premixed combustion situations such as furnaces, fuel and air streams are

mixed by ﬂuid ﬂow and turbulence, and the resulting combustion tempera-

tures, species concentrations and distribution of species are very much con-

trolled by ﬂuid ﬂow. Similarly, in premixed combustion, e.g. a spark ignition

IC engine, the geometry, ﬂuid ﬂow and turbulence created by the induction

and compression processes prior to ignition play an important role in the

combustion characteristics of the engine. We now present the transport

equations governing gaseous fuel combustion. We quote all equations in

compact sufﬁx notation, introduced ﬁrst in Chapter 3.

The governing equations of continuity and momentum developed in

Chapter 2 can be used in unchanged form.

Continuity

The continuity equation is

+ (

) = 0 (12.74)

It should be noted that density in combusting ﬂows is a variable, and depends

on pressure, temperature and species concentration.

Momentum equations

As in other ﬂows the velocity ﬁeld is governed by momentum equations

(

) + (

) =− + +F

(12.75)

where

is the viscous stress tensor and F

is the body force (which includes

gravity):

+−

(12.76)

The ﬁrst transport equation that is speciﬁc to reacting ﬂows is the equation

for the conservation of mass of a species k. It can be easily developed from

the general transport equation (2.39) by setting

= Y

. In sufﬁx notation the

statement is as follows.

Transport equations for species (k)

(

) + (

) =

+ M

(12.77)

Rate of Net rate of Net rate of Net rate of

change of

decrease of mass

increase of mass

mass of of species k due of species k due of species k due

species k to convection to diffusion to sources

∂

∂τ

∂

∂ρ

∂

Governing

equations for

combusting flows

12.14

ANIN_C12.qxd 29/12/2006 04:44PM Page 363

364 CHAPTER 12 CFD MODELLING OF COMBUSTION

In equation (12.77) D

is the species diffusion coefﬁcient (units m

/s). The

volumetric rate of generation (or destruction) of a species due to chemical

reactions appears as the source (or sink) term M

in each of their transport

equations. For kinetically controlled combustion this term takes the form of

(12.71) or (12.73). A wide variety of different models are available for phys-

ically controlled combustion. The most important ones will be discussed in

the remaining sections of this chapter.

It is common practice to assume a single diffusion coefﬁcient for all

species. This simpliﬁes equation (12.77) to

(

) + (

) =

D + M

(12.78)

Whilst the single D assumption is not always accurate (and may in fact be

quite inaccurate) it is very attractive, since it enables far-reaching simpliﬁca-

tion of combustion calculations.

Energy equation

In combusting ﬂows the temperature depends on the thermodynamic state

and the composition of the mixture. Some combustion models do not require

a transport equation for enthalpy: for example, in the laminar ﬂamelet model

the temperature is obtained from the laminar ﬂamelet library curves. Other

combustion submodels, however, require the solution of the transport equa-

tion for enthalpy. Chemical energy is released as heat during combustion,

and the resulting enthalpy is obtained by solving its transport equation:

(

h) + (

h) =+

− h

++S

rad

(12.79)

Net rate of Net rate of

Rate of

Net rate of increase of increase of

increase of increase of

change of

decrease of enthalpy due enthalpy due to

enthalpy due enthalpy due

enthalpy

enthalpy due

to diffusion

mass diffusion

to pressure

to radiative

to convection along gradients along gradients of

work heat transfer

of enthalpy species concentration

In equation (12.79) the source term of the transport equation for enthalpy

rad

is the radiation loss or gain. Viscous energy dissipation is normally

assumed to be negligible in low Mach number combusting ﬂows. Here h is

the mixture enthalpy per unit mass and h

is the speciﬁc enthalpy of species

k, and the summation is carried out over all N species considered in the

chosen reaction mechanism. The mixture Prandtl number is

, and Sc

the species Schmidt number, Sc

≡

The Prandtl number (

) is deﬁned as

Rate of momentum transport

Rate of energy transport

∂

∑

k =1

∂

ANIN_C12.qxd 29/12/2006 04:44PM Page 364

12.14 GOVERNING EQUATIONS FOR COMBUSTING FLOWS 365

The Lewis number is deﬁned as

≡=

The Schmidt number is deﬁned as

= Le

If a single diffusion coefﬁcient is used, i.e. D

= D for k = 1, 2,..., N, as in

equation (12.78), equation (12.79) can be written as

(

h) + (

h) =+− 1 h

++S

rad

(12.80)

This can be rewritten as

(

h) + (

h) =+− 1 h

++S

rad

(12.81)

Therefore for the case where the Lewis number is unity the enthalpy equa-

tion simpliﬁes to

(

h) + (

h) =++S

rad

(12.82)

For low-speed ﬂows

∂

t can be neglected, so under the assumptions of

low-speed ﬂow, single diffusion coefﬁcient and unity Lewis number, the

enthalpy equation has exactly the same form as the general transport equa-

tion (2.39). Furthermore, if the radiation source term S

rad

is also small, the

enthalpy is a conserved or passive scalar.

It should be noted that the thermal conductivity of individual species is a

function of temperature and can be evaluated using a polynomial expression

in the form k = a + bT + cT

+ dT

, where coefﬁcients are available in Reid

et al. (1987), for example. The binary diffusion coefﬁcient D

between two

species depends on the concentration of components, temperature and pres-

sure. Key equations which can be used to calculate D

and then a single

diffusion coefﬁcient D

(diffusion coefﬁcient of a single species to the rest

of mixture) are also available in the above reference. D

also depends on

individual species concentrations, temperature and pressure. The reader is

also referred to Paul and Warnatz (1998), who describe transport models

for species properties. In combustion modelling the required diffusion co-

efﬁcients and other transport properties are calculated (in detailed chemistry

∂

∑

k =1

∂

∑

k=1

∂

Rate of momentum transport

Rate of mass transport

Rate of energy transport

Rate of mass transport

ANIN_C12.qxd 29/12/2006 04:44PM Page 365

models) using packages such as CHEMKIN (Kee et al., 1996). Just to quote

some rough values, thermal conductivities of H

, H

O and CO

at 1000 K are

0.4 W/mK, 0.097 W/mK and 0.068 W/mK respectively. Approximate

values of binary diffusion coefﬁcient for mixtures with N

at 1000 K reported

in Paul and Warnatz (1998) are: for H–N

, 10 cm

/s; for OH–N

, 2 cm

/s;

and for O

–N

, 1.5 cm

/s. Lewis numbers for species in an H

/air mixture

in a one-dimensional laminar ﬂame calculation evaluated at the maximum

species mass fraction location are shown in Table 12.3. Apart from species H

and H

, other species have Lewis numbers close to 1.0. These values are

presented for information only, and the actual values should be evaluated

case by case depending on local conditions. The reader is also referred to

Clarke (2002) for further details on the calculation of Lewis numbers for

various fuel/air mixtures.

366 CHAPTER 12 CFD MODELLING OF COMBUSTION

Table 12.3 Lewis numbers for species evaluated in H

/air combustion

HOHOH

0.920 1.052 0.200 0.726 0.700 0.175 0.995 1.005

Other relationships

The total of mass fractions of fuel, oxidant and inert species is equal to 1, so

= 1 (12.83)

where k represents species.

The temperature can be calculated from the enthalpy by means of

T = (12.84)

where ∆h

is the enthalpy of combustion.

The average value of the speciﬁc heat O

is deﬁned as follows:

= c

dT (12.85)

Here

≡ Y

p,k

and c

p,k

is the speciﬁc heat of species k. Polynomial ﬁts for temperature-

dependent c

p,k

for various species are available in various texts: see for

example Cengel and Boles (2002).

The local density of the mixture is dependent on pressure, reactant and

product concentrations and on the mixture temperature. Its value can be

calculated from the following equation of state:

all species

∑



(T − T

ref

)

fuel species

h −

∑

∆h

all species

∑

ANIN_C12.qxd 29/12/2006 04:44PM Page 366