Brief Communication

Sensitivity of Calculated Extinction Strain Rate to Molecular

Transport Formulation in Nonpremixed Counterflow Flames

BRADLEY A. WILLIAMS

Chemistry Division, Code 6185 Naval Research Laboratory, Washington, DC 20375–5342, USA

The counterflow flame is a commonly used

geometry for experimental and chemical kinetic

modeling studies of nonpremixed combustion.

Because the flame structure is quasi-one-dimensional and more computationally tractable than

inherently multidimensional flame geometries

[1], detailed chemical kinetic modeling codes

for counterflow flames have come into widespread use.

An important characteristic of counterflow

flames is the extinction strain rate, the maximum velocity gradient which a flame can support and still burn. At strain rates below the

extinction value, there exist three solutions to

the governing equations which determine the

flame structure. Two of the solutions correspond to a stable flame and an essentially

non-reacting cold flow, while the third branch is

an unstable solution which cannot be physically

realized [2]. Various computational approaches

have been used to extend stable solutions onto

the unstable branch, and thus determine the

limits of bistability corresponding to flame extinction and autoignition [2, 3].

The extinction strain rate is important for

modeling turbulent combustion using the laminar flamelet approach [4], and is a figure of

merit for the effectiveness of fire suppressants

in nonpremixed flames [5]. It is desirable for

calculations to accurately predict extinction

strain rates. Substantial variations in predicted

extinction strain rates of nonpremixed methane/

air counterflow flames have been noted between

different chemical kinetic mechanisms [6].

Here, the predicted extinction strain rate is

found to also be sensitive to the treatment of

molecular transport in the computational

model.

The effect of transport on premixed flame

structure and extinction has been recently investigated by Paul and Warnatz [7], and by Ern and

Giovangigli [8, 9]. To our knowledge, the effect

of transport formalism on extinction of nonpremixed counterflow flames has not been previously reported in the literature. Furthermore,

earlier studies reporting computational predictions of extinction strain rates have not always

specified the manner in which the calculation

handled species transport.

In the present study, we employ two computer programs developed for counterflow

flames incorporating the CHEMKIN program

packages developed at Sandia National Laboratories. Initially, the chemical kinetics [10] and

molecular transport [11] routines were applied

to modeling one-dimensional premixed flames

[12], and more recently, to counterflow flames

[13]. The software developed at Sandia, including the counterflow code, OPPDIF, is available

commercially [14]. Independently, Nishioka et

al. adapted the Sandia premixed flame code for

counterflow flames [2]. The code of Nishioka et

al. has the capability to extend the stable flame

solution past the extinction point onto the unstable branch, which the OPPDIF code does

not. Extinction strain rates can be estimated

using OPPDIF by slowly increasing the flow

rates, relative to a previous solution, until the

calculation returns a non-reacting solution [15].

The counterflow flame codes use routines

which calculate molecular diffusivities and viscosities of multicomponent gas mixtures [11].

The two codes employ different versions of the

same routines, both of which permit different

options for treatment of molecular transport.

One choice available to the user is use of

mixture averaged diffusivities (weighted averCorresponding author. E-mail: brad@code6185.nrl.navy.mil age of the binary diffusion coefficients) or of

COMBUSTION AND FLAME 124:330–333 (2001)

0010-2180/01/$–see front matter © 2001 by The Combustion Institute

PII S0010-2180(00)00197-8 Published by Elsevier Science Inc.
multicomponent values obtained by matrix diagonalization. Also, the influence of thermal

diffusion on species transport (Soret effect) may

be included or neglected for the species H and H2.

In the results shown below, we have computed the structure of methane/air flames, on a

domain of 1.5 cm using the code of Nishioka et

al. Plug flow boundary conditions were used

(radial velocity gradients set to zero at the

boundaries). Axial velocities of the fuel and air

flows were chosen to be equal at the two

boundaries, with temperatures at both boundaries set to 300K. Central differencing was used

on the convective term in the species conservation equations.

The reaction mechanism, species thermodynamics functions, and transport parameters of

the Gas Research Institute GRI-Mech 2.11 [16]

were used for all calculations. All nitrogen

species except for N2 were deleted from the

mechanism.

Four series of calculations were performed:

1) multicomponent diffusivity formulation, with

thermal diffusivity considered, 2) multicomponent diffusivity formulation, thermal diffusivity

neglected, 3) mixture-averaged diffusivity formulation, thermal diffusivity considered, 4) mixture-averaged diffusivity formulation, thermal

diffusivity neglected. For the mixture-averaged

calculations, the “correction velocity” formalism was employed, which subtracts from the

calculated diffusion velocities of each species

any net convective velocity produced by the

accumulated errors [11], to ensure species conservation.

Near the extinction condition, the “one-point

temperature-controlling continuation” method

described by Nishioka et al. [2] was used. The

flame temperature at a single point within the

reaction zone was specified as an input condition, and the flow velocities of reactants were

allowed to vary (but were constrained to be

equal to one another) to satisfy the temperature

at the specified point. The specified temperature in each continuation run was typically

lowered by some 5 K from the value obtained at

the same grid point in the previous solution.

Figure 1 shows the peak flame temperature as

a function of strain rate for a nonpremixed

methane/air counterflow flame. The four series

of calculations shown in the plot are identical

except in the treatment of molecular transport.

The calculations were performed on meshes

having between 160 and 260 grid points, which

were periodically refined as the strain rate was

varied. For all calculations, local strain rates

were determined by numerical differentiation of

the calculated axial velocity with respect to

position. We estimate the uncertainty of the

strain rate determination due to the pointwise

differentiation and the finite mesh at 6 5 s-1.

Use of different transport treatments produces significant differences in the predicted

extinction strain rates. The lowest and the highest of the values differ by almost 100 s-1, ranging

from approximately 450 s-1 to 550 s-1. For

comparison, most experimental measurements

cluster around 400 s-1 [17]. The calculations

using the multicomponent transport formulation predict a lower extinction strain rate than

do those using the mixture-averaged approximation. For both types of calculations, inclusion

of the Soret effect for H and H2 increases the

predicted extinction strain rate. Additionally,

we have compared methane/air extinction strain

rates estimated using the OPPDIF code for

multicomponent and mixture-averaged transport formulations (neglecting the Soret effect),

and find a similar sensitivity of extinction strain

rate to transport treatment.

Figure 2 shows plots of temperature and H

atom mole fraction as a function of position for

Fig. 1. Relationship between maximum flame temperature

and peak air-side strain rate using the counterflow flame

code of Ref. [2] with different options specified in the

transport formulation. The extinction strain rates are given

by the maximum values obtained as the peak temperature

decreases.

CALCULATED EXTINCTION STRAIN RATE SENSITIVITY 331
each of the four transport treatments at a strain

of 351 s-1. The temperature profile has a similar

peak location and shape for all the transport

treatments, as does the H atom profile. The

peak H atom concentration varies substantially

between the calculations with a difference of

about 20% between the highest value (mixture

averaged transport with thermal diffusivity) and

the lowest value (multicomponent transport, no

thermal diffusivity). The other two cases have

peak H atom concentrations midway between

these two values.

The present results demonstrate that use of

different transport formulations has a substantial effect (on the order of 20%) on the predicted extinction strain rate of a nonpremixed

methane/air flame. This implies that accurate

treatment of molecular transport is a prerequisite to obtaining credible computational predictions of extinction strain rates of nonpremixed

flames. The extension of these findings to other

fuel/oxidizer systems is beyond the scope of the

present work, but comparison of hydrogen and

methane fuels in premixed flames [8, 9] shows

that the former is far more sensitive to transport

treatment.

The effect of transport also appears to be, in

general, much more significant in a nonpremixed than in a premixed flame. A series of

calculations (using the PREMIX code) of the

burning velocity of a mixture of 15.8% O2, 7.9%

CH4, and 76.3% N2 yielded values which varied

by less than 3% (between 22.5 and 23.1 cm/s)

between the four treatments of transport considered here. Similarly, Ern and Giovangigli [9]

found differences of only about 2% between the

calculated extinction strain rates of premixed

methane/air counterflow flames using mixtureaveraged, multicomponent, and multicomponent-Soret transport formulations. The only

nonpremixed flame modeled in Ref. [9] was an

axisymmetric methane-air co-flow flame, in

which the Soret effect was reported to have only

a minor influence on flame structure. This result

may suggest an influence of strain rate (typically

higher in a counterflow than in a co-flow flame)

on the relative importance of transport.

Using the mathematical formulation employed by the CHEMKIN codes, there is a

significant penalty in computational time in

using the multicomponent transport formulation. Both Ern and Giovangigli [8, 9] and Blasenbrey et al. [18], however, have proposed

alternate formulations of multicomponent

transport to reduce the computational expense.

The appropriate treatment of the Soret effect

is a more problematic issue. The transport

routines used for the CHEMKIN codes can

include the Soret effect only for H and H2, and

neglect the Dufour effect (influence of species

diffusion on thermal transport) entirely. The

approach used by Ern and Giovangigli is to

include both the Soret and Dufour effects for all

species. In Ref. [9] Ern and Giovangigli comment that under certain conditions, the effect of

thermal diffusion considered for all species can

be opposite the effect when considered only for

H and H2. On the other hand, Paul and Warnatz [7] include the Soret effect only for H and

H2 on the grounds that the uncertainties are too

Fig. 2. Profiles of temperature and atomic hydrogen mole

fraction in the vicinity of the reaction zone using the

different molecular transport treatments. All four calculations used identical reactant flows; the maximum air-side

strain rate for all calculations was 351 s-1.

332 B. A. WILLIAMS
large to justify its inclusion for other species. If

the Soret effect does need to be included for all

species, the transport formulations available in

the CHEMKIN routines will have to be augmented in order to obtain credible predictions

of extinction strain rates of nonpremixed

flames.

Since the extinction strain rate of nonpremixed flames appears to be a property particularly sensitive to molecular transport, care is

required in its modeling. Both the multicomponent formulation and the Soret effect significantly influence the predicted value, and the

computational methodology must be specified

to allow meaningful comparisons between modeling studies. The need for detailed consideration of the accuracy of the transport coefficients and the computational methodology, as

done in Refs. [7–9], is underscored by the

present results.

The author thanks M. Nishioka for sharing the

counterflow flame code described in Ref. [2]. This

work was supported by the US Department of

Defense’s Next Generation Fire Suppression

Technology Program funded by the DoD Strategic

Environmental Research and Development Program.

REFERENCES

1. Kee, R. J., Miller, J. A., Evans, G. H., and DixonLewis, G., Proc. Combust. Inst. 22:1479–1494 (1988).

2. Nishioka, M., Law, C. K., and Takeno, T., Combust.

Flame104:328–342 (1996).

3. Giovangigli, V. and Smooke, M. D., Combust. Sci.

Technol.53:23–49 (1987).

4. Candel, S., Veynante, D., Lacas, F., Darabiha, N., and

Rolon, C., Combust. Sci. Technol.98:245–264 (1994).

5. Papas, P., Fleming, J. W., and Sheinson, R. S., Proc.

Combust. Inst.26:1405–1412 (1996).

6. Tanoff, M. A., Dobbins, R. R., Smooke, M. D.,

Burgess, Jr., D. R., Zachariah, M. R., Tsang, W., and

Westmoreland, P.R., Halon Options Technical Working Conference, Albuquerque, NM, May 6–8, 1997,

pp. 116–127.

7. Paul, P. H. and Warnatz, J., Proc. Combust. Inst.27:

495–504 (1998).

8. Ern, A. and Giovangigli, V., Combust. Theor. Model.2:

349–372 (1998).

9. Ern, A. and Giovangigli, V., Combust. Sci. Technol.149:157–181 (1999).

10. Kee, R. J., Rupley, F. M., and Miller, J. A., Chemkin II:

A Fortran Chemical Kinetics Package for the Analysis of

Gas-Phase Chemical Kinetics, Sandia Report SAND89-

8009, Sandia National Laboratories, Livermore, CA,

1989.

11. Kee, R. J., Dixon–Lewis, G., Warnatz, J., Coltrin,

M. E., Miller, J. A., and Moffat, H. K., A Fortran

Computer Code Package for the Evaluation of Gas–

Phase, Multicomponent Transport Properties,Report

SAND86-8246B, Sandia National Laboratories, Livermore, CA, 1986.

12. Kee, R. J., Grcar, J. F., Smooke, M. D., Miller, J. A.,

and Meeks, E., A Program for Modeling Steady, Laminar, One-Dimensional Premixed Flames,Report

SAND85-8240, Sandia National Laboratories, Livermore, CA, 1985.

13. Lutz, A. E., Kee., R. J., Grcar, J. F., and Rupley, F. M.,

OPPDIF: A Fortran Program of Computing Opposedflow Diffusion Flames, Report SAND96-8243, Sandia

National Laboratories, Livermore, CA, 1996.

14. Kee, R. J., Rupley, F. M., Miller, J. A., Coltrin, M. E.,

Grcar, J. F., Meeks, E., Moffat, H. K., Lutz, A. E.,

Dixon-Lewis, G., Smooke, M. D., Warnatz, J., Evans,

G. H., Larson, R. S., Mitchell, R. E., Petzold, L. R.,

Reynolds, W. C., Caracotsios, M., Stewart, W. E., and

Glarborg, P., CHEMKIN Collection, Release 3.5, Reaction Design, Inc., San Diego, CA, 1999.

15. Pitts, W. M., Yang, J. C., Huber, M. L., and Blevins, L.

G., Characterization and Identification of Super-Effective Thermal Fire Extinguishing Agents—First Annual

Report, Report NISTIR 6414, National Institute of

Standards and Technology, Gaithersburg, MD, 1999.

16. Bowman, C. T., Hanson, R. K., Gardiner, W. C.,

Lissianski, V., Frenklach, M., Goldenberg, M., and

Smith, G. P., GRI-Mech 2.11 — An Optimized Detailed

Chemical Reaction Mechanism for Methane Combustion and NO Formation and Reburning, GRI Report

GRI-97/0020, 1997, Gas Research Institute, Chicago.

http://www.ME.Berkeley.edu/gri2mech/.

17. Zegers, E. J. P., Williams, B. A., Fisher, E. M.,

Fleming, J. W., and Sheinson, R. S., Combust.

Flame121:471–487 (2000).

18. Blasenbrey, T., Schmidt, D., and Maas, U., Proc.

Combust. Inst. 27:505–511 (1998).

Received 1 March 2000; revised 19 July 2000; accepted

8 August 2000

CALCULATED EXTINCTION STRAIN RATE SENSITIVITY 333