**Simulation of H2-Air Non-Premixed Flame Using Combustion Simulation Technique to Reduce Chemical Mechanisms**

Kazui Fukumoto and Yoshifumi Ogami *Ritsumeikan University, 1-1-1 Noji-Higashi, Kusatsu, Shiga Japan* 

#### **1. Introduction**

356 Fluid Dynamics, Computational Modeling and Applications

Nebot, E.; Casanueva, J. F.; Casanueva, T. & Sales, D. (2007). *Model for fouling deposition on* 

Ocłon, P. & Łopata, S. (2011). *Modelling of the flow distribution inside the collectors of the high* 

Quan, Z., Chen, Y. & Chongfang, M. A. (1997). *Experimental Study of Fouling on Heat Transfer* 

Reddy, J. N. (1993). *Introduction to Finite Element Method*, McGraw Hill, ISBN10: 0070513554,

Rocha, L. A. O.; Saboya, F. E. M. & Vargas, J.V.C. (1997). *A comparative study of elliptical and* 

Saboya, F. E. M. & Sparrow, E. M. (1974). *Local and Average Transfer Coefficients for One-Row* 

Saboya, S. M. & Saboya, F. E. M. (1981). *Transfer Coefficients for Plate Fin and Elliptical Tube* 

Smith, I. M & Griffiths, D. V. (2004). *Programming the Finite Element Method*, Willey, ISBN-10

Taler, D. (2009). *Dynamics of Tube Heat Exchangers*, Monograph 193, AGH University of Science

Taler, J. & Duda, P. (2006). *Solving Direct and Inverse Heat Conduction Problems*, Springer,

White F. (2003). *Fluid Mechanics*, McGraw Hill, ISBN - 10 0071215662, ISBN - 13 978-

Zienkiewicz, O. C. & Taylor, R. L. (2000). *The Finite Element Method 5th Ed*., Elsevier, Vol.1-3,

– 400, Cracow University of Technology press, , Cracow

Journal of Heat and Fluid Flow, pp. 247–252, Vol. 18, Issue 2

and Technology Press, ISSN 0867-6631, Cracow, Poland (in Polish)

Engineering, pp. 535 – 540, Vol. 16, No. 4

8817 , pp. 410-419, vol. 31

ISBN- 13: 978-0070513556

Engineering, pp. 153-162

ISBN 978-3-540-33470-5

ISBN 978-0-7506-5055-7

047084969X, ISBN-13 978- 0470849699

272, Vol. 96

0071215664

*power plant steam condensers cooled with seawater: Effect of water velocity and tube material*, International Journal of Heat and Mass Transfer, pp. 3351–3358, Vol. 50 Nishiyama, H.; Ota, T. & Matsumo, T. (1988). *Heat Transfer and Flow around an Elliptic* 

*Cylinders in Tandem Arrangement*, JSME International Journal Series II, ISSN 0914-

*performance heat exchanger*, Technical Transactions Mechanics 4-M/2011/B, pp. 393

*Surface During Forced Convective Heat Transfer*, Chinese Journal of Chemical

*circular sections in one- and two-row tubes and plate fin heat exchangers*, International

*Plate Fin and Tube Heat Exchanger Configurations*, Journal of Heat Transfer, pp. 265-

*Heat Exchangers*, Proceedings of the VI Brazilian Congress of Mechanical

Combustion simulation that uses computational fluid dynamics (CFD) has been widely adopted as the design tool for combustion equipment. Because flow inside such equipment is generally turbulent, turbulence and combustion models are needed to simulate combustion; many combustion simulations have been performed to verify a system's internal state, such as velocity, pressure, mole fractions of chemical species, and temperature.

Combustion simulation of a confined impinging jet reactor has been performed by the large eddy simulation (LES) model serving as the turbulence model and the presumed probability density function (PDF) serving as the combustion model (Daniele, 2009). The analysis shows that a confined impinging jet reactor is indeed an interesting device because of its high mixing efficiency and absence of stagnant and recirculation zones.

Under the condition of moderate or intense low-oxygen dilution (MILD), the effect of H2 on H2-CH4 turbulent non-premixed flames was investigated with the improved standard *k* model as the turbulence model and with the eddy dissipation concept (EDC) model (Amir et al., 2010). Simulation results show that H2 addition to CH4 leads to improved mixing, increase in turbulent kinetic energy decay along the flame axis, increase in flame entrainment, higher reaction intensities, and increase in mixture ignitability and rate of heat release.

Although combustion simulation was considered to be an efficient designing tool, considerable computational time was needed to calculate the chemical reaction. Combustion models that detail chemical mechanisms require reaction calculations involving *n*dimensional ordinary differential equations (ODEs) that are solved according to the number of chemical species involved. Therefore, reducing computation time for the combustion simulation is a significant problem. If computation time could be easily reduced according to the required prediction accuracy, we would be able to obtain the results more quickly. For example, O is a significant species whose mass fraction is necessary to compute the amount of NO present; therefore, the accuracy of the mass fraction of O cannot be neglected. To determine this mass fraction with sufficient accuracy, it is necessary to build a reduced mechanism including O. Generally, a quasi-steady state or partial equilibrium is assumed when the reduced mechanism is built (Warrants et al., 2006).

The chemical equilibrium method (Nagai et al., 2002) does not use reaction equations; instead, the equilibrium composition of a chemical system is determined by minimizing the

Simulation of H2-Air Non-Premixed Flame

TNF web site (Barrow, 2003b), and *T* was computed as

at constant pressure of chemical species *J* .

**2.1 Combustion simulation technique** 

(including the overall reaction).

where

species *J* , *YJ*

is the turbulent viscosity,

Using Combustion Simulation Technique to Reduce Chemical Mechanisms 359

, *<sup>J</sup> i T*

represents the density, *ui* represents the gas velocity of *i* components,

enthalpy, rad *q* is the source term of the radiation, *<sup>i</sup> x* is the coordinate in the *i* direction,

*<sup>h</sup>* is the effective Prandtl number,

and *T* is the temperature. The symbols "~" and "–" placed over the variables indicate density weighted and time average, respectively. rad *q* was modeled as described on the

> <sup>0</sup> , *J J <sup>J</sup> J*

*Cp*

*h Yh*

where 0 *T* 298.15 . *<sup>J</sup>*<sup>0</sup> *h* represents the standard enthalpy of the formation of chemical

The proposed technique is explained using the example of a H2-air turbulent non-premixed flame. Figure 1 illustrates the scheme of the proposed technique. In this technique, chemical species are divided, roughly, into two groups: group a—the species categorized by the chemical equilibrium, and group b—the species categorized by the chemical kinetics

**a Species obtained by Look-up table of** 

chemical equations

H2 O H OH O2 H2O HO2 H2O2 NO HNO N2O

Partial species was improved by

**b Species modified by chemical equations**

In the proposed technique, the mass fraction equation of chemical species *J* was used as

*j j*

*x x* 

� ( , ) ( ) *<sup>j</sup> <sup>J</sup> jJ <sup>j</sup> J J*

*uY j uY w*

 

(5)

**chemical equilibrium**

H2 H O OH

Fig. 1. Overview of proposed technique

N NO2 N2

*x* 

0

is the mass fraction of chemical species *<sup>J</sup>* , and *<sup>C</sup> <sup>J</sup> <sup>p</sup>* is the mean specific heat

(3)

*h* is the

is the heat conductivity,

(4)

*t*

*q*

*T T*

Gibbs free energy, which is subject to the conservation of the chemical elements involved in the combustion. The composition is obtained by solving a system of simultaneous nonlinear equations with five unknowns (one for each of the elemental species involved in that study, namely, C, H, O, and N, plus one). In another study, a partial non-premixed flame of CH4 was modeled using the vortex method in combination with the chemical equilibrium method (Ogami et al., 2010). Although this method did not use any chemical equations, the calculated temperature and main products were in good agreement with experimental results. In some cases, the data obtained from the simulation concerning the minor products did not agree with the experimental results.

A combustion simulation technique based on a combination of the chemical equilibrium method and the EDC model was developed in a previous study and validated by simulation of a H2-air turbulent non-premixed flame (Fukumoto et al., 2010) and a CO-H2-air turbulent non-premixed flame (Fukumoto et al., 2009). An advantage of the technique is the ease with which a reduced chemical mechanism can be built according to the accuracy requirement for the chemical species. Thus, the technique can predict intermediate species with high accuracy when a reduced mechanism for the minor species is built. The key to this proposed technique is the incorporation of the concept of equilibrium into the chemical kinetics, which decreases the dimensionality of the chemical reaction calculation and reduces the computational time for a combustion simulation.

This chapter presents a technique for combustion simulation based on a chemical equilibrium method and an EDC model (Magnussen, 1989, 2005). The advantage of the proposed technique is that it makes it easy for us to construct a reduced mechanism according to the accuracy requirements of the chemical species. The proposed technique is validated through a simulation of a H2-air non-premixed flame (Barrow et al., 2003a). The results computed using the proposed technique are then compared with experimental data and computational data obtained by the EDC model.

#### **2. Method**

Ansys FLUENT version 13.0 (Ansys, Inc., 2011) was used in this study. The temperature, enthalpy, and reaction rate of the chemical species were calculated by a user-defined function (UDF) in FLUENT. The UDF is a function that can be programmed to load dynamically with the FLUENT solver in order to enhance the standard features of FLUENT. Moreover, the momentum equation, turbulence model, and mass fraction equation were computed using Ansys FLUENT.

To solve the momentum equation, velocity and pressure were coupled by the SIMPLE method; the basic Reynolds stress model with a function of the enhanced wall treatment in Ansys FLUENT (Ansys, Inc., 2010) was chosen as the turbulence model. The temperature was obtained by computation of the following energy equation:

$$\frac{\partial(\overline{\rho}\overline{u}\_i\overline{h})}{\partial\mathbf{x}\_i} = -\frac{\partial}{\partial\mathbf{x}\_i}\left(\overline{q}\_f + \overline{\rho}\overline{\overline{u}\_i'}h'\right) + \overline{q}\_{\text{rad}'} \tag{1}$$

$$
\rho \overline{\mu}\_i^{\overline{\text{h}}} \dot{\mu}^{\cdot} = -\frac{\mu\_t}{\sigma\_{\text{h}}} \frac{\partial \tilde{\text{h}}}{\partial \mathbf{x}\_i} \,' \tag{2}
$$

$$
\overline{\eta}\_{\parallel} = -\mathcal{X} \frac{\partial \overline{T}}{\partial \mathbf{x}\_i} \,' \tag{3}
$$

where represents the density, *ui* represents the gas velocity of *i* components, *h* is the enthalpy, rad *q* is the source term of the radiation, *<sup>i</sup> x* is the coordinate in the *i* direction, *t* is the turbulent viscosity, *<sup>h</sup>* is the effective Prandtl number, is the heat conductivity, and *T* is the temperature. The symbols "~" and "–" placed over the variables indicate density weighted and time average, respectively. rad *q* was modeled as described on the TNF web site (Barrow, 2003b), and *T* was computed as

$$
\overline{T} = T\_0 + \frac{\tilde{h} - \sum\_{f} \tilde{Y}\_f h\_{f0}}{C p\_f},
\tag{4}
$$

where 0 *T* 298.15 . *<sup>J</sup>*<sup>0</sup> *h* represents the standard enthalpy of the formation of chemical species *J* , *YJ* is the mass fraction of chemical species *<sup>J</sup>* , and *<sup>C</sup> <sup>J</sup> <sup>p</sup>* is the mean specific heat at constant pressure of chemical species *J* .

#### **2.1 Combustion simulation technique**

358 Fluid Dynamics, Computational Modeling and Applications

Gibbs free energy, which is subject to the conservation of the chemical elements involved in the combustion. The composition is obtained by solving a system of simultaneous nonlinear equations with five unknowns (one for each of the elemental species involved in that study, namely, C, H, O, and N, plus one). In another study, a partial non-premixed flame of CH4 was modeled using the vortex method in combination with the chemical equilibrium method (Ogami et al., 2010). Although this method did not use any chemical equations, the calculated temperature and main products were in good agreement with experimental results. In some cases, the data obtained from the simulation concerning the minor products

A combustion simulation technique based on a combination of the chemical equilibrium method and the EDC model was developed in a previous study and validated by simulation of a H2-air turbulent non-premixed flame (Fukumoto et al., 2010) and a CO-H2-air turbulent non-premixed flame (Fukumoto et al., 2009). An advantage of the technique is the ease with which a reduced chemical mechanism can be built according to the accuracy requirement for the chemical species. Thus, the technique can predict intermediate species with high accuracy when a reduced mechanism for the minor species is built. The key to this proposed technique is the incorporation of the concept of equilibrium into the chemical kinetics, which decreases the dimensionality of the chemical reaction calculation and reduces the

This chapter presents a technique for combustion simulation based on a chemical equilibrium method and an EDC model (Magnussen, 1989, 2005). The advantage of the proposed technique is that it makes it easy for us to construct a reduced mechanism according to the accuracy requirements of the chemical species. The proposed technique is validated through a simulation of a H2-air non-premixed flame (Barrow et al., 2003a). The results computed using the proposed technique are then compared with experimental data

Ansys FLUENT version 13.0 (Ansys, Inc., 2011) was used in this study. The temperature, enthalpy, and reaction rate of the chemical species were calculated by a user-defined function (UDF) in FLUENT. The UDF is a function that can be programmed to load dynamically with the FLUENT solver in order to enhance the standard features of FLUENT. Moreover, the momentum equation, turbulence model, and mass fraction equation were

To solve the momentum equation, velocity and pressure were coupled by the SIMPLE method; the basic Reynolds stress model with a function of the enhanced wall treatment in Ansys FLUENT (Ansys, Inc., 2010) was chosen as the turbulence model. The temperature was

> ( , ) ( ) *<sup>i</sup> J i*

> > � , *<sup>t</sup> <sup>i</sup>*

*<sup>h</sup> u h*

*u h q uh q*

*i i*

*x x*

�

*h i*

*x*

rad

(1)

(2)

did not agree with the experimental results.

computational time for a combustion simulation.

and computational data obtained by the EDC model.

obtained by computation of the following energy equation:

**2. Method** 

computed using Ansys FLUENT.

The proposed technique is explained using the example of a H2-air turbulent non-premixed flame. Figure 1 illustrates the scheme of the proposed technique. In this technique, chemical species are divided, roughly, into two groups: group a—the species categorized by the chemical equilibrium, and group b—the species categorized by the chemical kinetics (including the overall reaction).

In the proposed technique, the mass fraction equation of chemical species *J* was used as

$$\frac{\partial}{\partial \mathbf{x}\_{j}} (\overline{\rho} \overline{u}\_{j} \tilde{Y}\_{l}) = -\frac{\partial}{\partial \mathbf{x}\_{j}} \left( \overline{j}\_{jl} + \overline{\rho} \overline{u}\_{j}^{\Box} Y\_{l} \right) + \overline{w}\_{l'l'} \tag{5}$$

Simulation of H2-Air Non-Premixed Flame

where *J* indicates a*J* .

where 2.1377 *C*

where *C* 0.408 

structures; b

beq

When

where b *wJ*

b *wJ*

b *YJ*

.

turbulent kinetic energy.

.

b

and temperature in each computational cell:

*J*

*w*

was modeled as

\* 

Using Combustion Simulation Technique to Reduce Chemical Mechanisms 361

<sup>3</sup> 3 0 (1.0 ) . *YY Y <sup>J</sup> J J*

 

> 1 4

1 2 *C* , 

b b beq b b

*YY Y Y*

*J J J J*

*J J JJ*

0 \* 0

b beq b b b

\* 0 \* 0

( ) , ( 0.0),

*Y Y Y Y Y*

, ) 0.0) ( ) max( , (

*J*

represents the viscosity, the turbulent dissipation rate, and *k* the

<sup>2</sup> *<sup>C</sup>* , *<sup>k</sup>* 

The proposed technique was needed to compute *wJ* in Equation (5), and different expressions were used to compute , *wJ* depending on the group that the species belonged to. The reaction rate of chemical species in group b can be obtained with the EDC model. The source term in Equation (5) for these modified chemical species was computed as

(12)

(13)

(15)

*YJ* , the mass fraction

, the chemical species b*J* in the fine

(16)

(14)

\*

<sup>0</sup> *YJ* , the chemical species b*J* in the surrounding gas; and beq

Since the lookup table included the mass fraction of the chemical species in group b, the production resulting from using the values take from the table must be taken into account.

 was obtained by calculating the detailed chemical mechanism. The EDC model assumes that combustion occurs in a constant pressure reactor. The following equation was integrated under the initial conditions obtained from the current mass fractions, density,

> b b \* , *J J dY w dt*

 can be obtained by computing detailed chemical mechanisms based on the Arrhenius equation. Equation (16) is a stiff ODE, and therefore, the robust computational code CVODE

Although the reaction rates of fuels can be computed by Equation (15), the flame sometimes blows out. That is why we used the following equations to estimate the reaction rates of the

represents the reaction rate of chemical species b*J* in the fine structures, and

0 *YJ* , these equations were similar to what was seen in the EDC model.

*Y*

 

was given as

2

 

 

 

where b *wJ* represents the reaction rate in group b; b *YJ*

of the chemical species b*J* obtained from the lookup table.

(Hindmarsh et al., 2009) was used to acquire the steady solution.

2

*J*

$$
\overline{\dot{j}}\_{\parallel l} = -\overline{\rho}D\_{l} \frac{\partial \overline{\Psi}\_{l}}{\partial \mathbf{x}\_{j}} \,' \, \tag{6}
$$

$$
\overline{\rho}\overline{\mu}\_{\dot{I}}^{\dot{I}}Y\_{\dot{I}} = -\frac{\mu\_{t}}{\sigma\_{m}}\frac{\partial \tilde{Y}\_{\dot{I}}}{\partial \mathbf{x}\_{i}}\tag{7}
$$

where *wJ* represents the reaction rate of chemical species *J* ; *DJ* , the effective diffusion coefficient of chemical species *J* ; and *<sup>m</sup>* , the effective Schmidt number. The eddy dissipation concept (Magnussen, 1989; Gran et al., 1996; Magnussen, 2005) was used to model *wJ* , and *DJ* was estimated using Smooke's model. To obtain *wJ* in Equation (5), calculations of the chemical reaction and chemical equilibrium were performed; the mass fraction equation can be computed by assignment of a value to *wJ* .

To perform the chemical equilibrium calculation, we needed the mass fraction of elemental species *I* and enthalpy in group b. We therefore added the mass fraction equation of elemental species *J* in group a:

$$\frac{\partial(\overline{\rho}\overline{u}\_{j}\overline{Z}\_{l})}{\partial\mathbf{x}\_{j}} = \frac{\partial}{\partial\mathbf{x}\_{j}}\left(\frac{\mu\_{l}}{\sigma\_{m}}\frac{\partial\overline{Z}\_{l}}{\partial\mathbf{x}\_{i}}\right) + \overline{w}\_{l,r} \tag{8}$$

where *ZI* represents the mass fraction of elemental species *<sup>I</sup>* . *wI* is the source term that exchanges the mass in group b for that in group a:

$$
\overline{w}\_{l} = -\mathcal{M}\_{l} \left[ \sum\_{f\_{\mathbf{b}}} \frac{a\_{\frac{\eta\_{\mathbf{b}}}{\eta\_{\mathbf{b}}}} \overline{w}\_{f\_{\mathbf{b}}}}{\mathcal{M}\_{f\_{\mathbf{b}}}} \right], \tag{9}
$$

where *MI* represents the atomic weight of species *I* , b*J* the chemical species *J* in group b, <sup>b</sup> *IJ a* the number of elemental species *I* included in chemical species b*J* , b *wJ* the reaction rate of chemical species b*J* , and b *MJ* the mole weight of chemical species b*J* . The enthalpy in group a was approximated as

$$
\tilde{h}\_{\mathbf{a}} = \tilde{h} \sum\_{\mathbf{l}} \tilde{Z}\_{\mathbf{l}\prime} \tag{10}
$$

where a *h* represents the enthalpy in group a. The source term in Equation (5) for the equilibrium species a*J* was modeled as

$$\overline{\boldsymbol{w}}\_{\boldsymbol{I}\_{\boldsymbol{u}}} = \frac{\overline{\rho} (\boldsymbol{\xi}^{\star})^{2}}{\boldsymbol{\pi}^{\star}} (\boldsymbol{Y}\_{\boldsymbol{I}\_{\text{aug}}} - \boldsymbol{Y}\_{\boldsymbol{I}\_{\text{u}}}^{0})\_{\prime} \tag{11}$$

where a*J* is the chemical species *J* in group a, a *wJ* is the reaction rate of chemical species a*<sup>J</sup>* , \* is the fine scale, is the time scale of the reaction, and aeq *YJ* is the mass fraction of chemical species a*J* obtained from the lookup table, which was prepared before the combustion simulation. The lookup table stored the results of the chemical equilibrium calculation when *ZI* and a *<sup>h</sup>* were varied, and these values were also used to extract the equilibrium values from the lookup table.

a <sup>0</sup> *YJ* in Equation (11) was estimated as

$$
\tilde{Y}\_{\!\!\!\!\/} = \xi^3 Y\_{\!\!\!\/\!\/}^\* + (\mathbf{1}.0 - \xi^3) Y\_{\!\!\!\/\!\/}^0. \tag{12}
$$

where *J* indicates a*J* .

360 Fluid Dynamics, Computational Modeling and Applications

*jJ J*

*u Y*

fraction equation can be computed by assignment of a value to *wJ* .

exchanges the mass in group b for that in group a:

The enthalpy in group a was approximated as

equilibrium values from the lookup table.

<sup>0</sup> *YJ* in Equation (11) was estimated as

and a *<sup>h</sup>*

represents the enthalpy in group a.

coefficient of chemical species *J* ; and

elemental species *J* in group a:

where *ZI*

where a *h*

a

is the fine scale,

calculation when *ZI*

*<sup>Y</sup> j D*

� , *<sup>J</sup> <sup>t</sup> j J*

where *wJ* represents the reaction rate of chemical species *J* ; *DJ* , the effective diffusion

dissipation concept (Magnussen, 1989; Gran et al., 1996; Magnussen, 2005) was used to model *wJ* , and *DJ* was estimated using Smooke's model. To obtain *wJ* in Equation (5), calculations of the chemical reaction and chemical equilibrium were performed; the mass

To perform the chemical equilibrium calculation, we needed the mass fraction of elemental species *I* and enthalpy in group b. We therefore added the mass fraction equation of

( ) ( ) , *j I <sup>t</sup> <sup>I</sup>*

represents the mass fraction of elemental species *<sup>I</sup>* . *wI* is the source term that

*J J*

 

where *MI* represents the atomic weight of species *I* , b*J* the chemical species *J* in group b, <sup>b</sup> *IJ a* the number of elemental species *I* included in chemical species b*J* , b *wJ* the reaction

> <sup>a</sup> , *<sup>I</sup> I*

a aeq a 2 <sup>0</sup> ( ) ( ), *w YY J JJ*

where a*J* is the chemical species *J* in group a, a *wJ* is the reaction rate of chemical species a*<sup>J</sup>* , \*

chemical species a*J* obtained from the lookup table, which was prepared before the combustion simulation. The lookup table stored the results of the chemical equilibrium

b b b b

*a w*

*M*

, *IJ <sup>J</sup>*

*j jmi u Z Z*

*xxx*

*I I*

*w M*

rate of chemical species b*J* , and b *MJ* the mole weight of chemical species b*J* .

The source term in Equation (5) for the equilibrium species a*J* was modeled as

 

 

*I*

(8)

(9)

*h hZ* (10)

(11)

is the time scale of the reaction, and aeq *YJ* is the mass fraction of

were varied, and these values were also used to extract the

*w*

� , *J*

(6)

*<sup>m</sup>* , the effective Schmidt number. The eddy

(7)

*j*

*x*

*m i Y*

*x*

\* was modeled as

$$\xi^\* = \mathbf{C}\_{\xi} \left( \frac{\mu \epsilon}{\overline{\rho} k^2} \right)^{\mathbf{\varUpsilon}},\tag{13}$$

where 2.1377 *C* . represents the viscosity, the turbulent dissipation rate, and *k* the turbulent kinetic energy. was given as

$$
\tau^\* = \mathbf{C}\_{\tau} \left( \frac{\mu}{\overline{\rho}\epsilon} \right)^{\frac{1}{2}},
\tag{14}
$$

where *C* 0.408 .

The proposed technique was needed to compute *wJ* in Equation (5), and different expressions were used to compute , *wJ* depending on the group that the species belonged to. The reaction rate of chemical species in group b can be obtained with the EDC model. The source term in Equation (5) for these modified chemical species was computed as

$$
\overline{w}\_{I\_{\rm b}} = \begin{cases}
\overline{\rho} (\overline{\xi}^\*)^2 \left[ \max (Y\_{I\_{\rm b}}^\* - Y\_{I\_{\rm b}}^0, Y\_{I\_{\rm b\rm w}}) \right], & (Y\_{I\_{\rm b}}^\* - Y\_{I\_{\rm b}}^0 \ge 0.0) \\
\overline{\rho} (\overline{\xi}^\*)^2 \left[ Y\_{I\_{\rm b}}^\* + Y\_{I\_{\rm b\rm w}} - Y\_{I\_{\rm b}}^0 \right], & (Y\_{I\_{\rm b}}^\* - Y\_{I\_{\rm b}}^0 < 0.0),
\end{cases} \tag{15}
$$

where b *wJ* represents the reaction rate in group b; b *YJ* , the chemical species b*J* in the fine structures; b <sup>0</sup> *YJ* , the chemical species b*J* in the surrounding gas; and beq *YJ* , the mass fraction of the chemical species b*J* obtained from the lookup table.

Since the lookup table included the mass fraction of the chemical species in group b, the production resulting from using the values take from the table must be taken into account. When beq 0 *YJ* , these equations were similar to what was seen in the EDC model.

b *YJ* was obtained by calculating the detailed chemical mechanism. The EDC model assumes that combustion occurs in a constant pressure reactor. The following equation was integrated under the initial conditions obtained from the current mass fractions, density, and temperature in each computational cell:

$$\frac{d\boldsymbol{Y}\_{f\_b}^\*}{dt} = \frac{\boldsymbol{w}\_{f\_b}^\*}{\rho}\tag{16}$$

where b *wJ* represents the reaction rate of chemical species b*J* in the fine structures, and b *wJ* can be obtained by computing detailed chemical mechanisms based on the Arrhenius equation. Equation (16) is a stiff ODE, and therefore, the robust computational code CVODE (Hindmarsh et al., 2009) was used to acquire the steady solution.

Although the reaction rates of fuels can be computed by Equation (15), the flame sometimes blows out. That is why we used the following equations to estimate the reaction rates of the

Simulation of H2-Air Non-Premixed Flame

were the concentrations after combustion.

extended Zeldovich mechanism as

the reaction rate of the NO was expressed as

alternately, the PDF was related to two variables:

were taken into account using the proposed technique.

and *<sup>a</sup> h*

lookup table.

chosen as

1 and 

**2.4 Numerical conditions** 

Using Combustion Simulation Technique to Reduce Chemical Mechanisms 363

where *<sup>J</sup> h* represents the enthalpy of chemical species *J* , and *<sup>r</sup> h* represents the enthalpy of reactants (before combustion). When we found *T* using Equation (22), the obtained results

FLUENT version 13.0 has a function for computing the amount of the thermal NOX using the

 N OH H NO. (25) Because we used the basic function in FLUENT, we describe only the outline of the method to estimate the amount of NO here (Warnatz et al., 2006; Ansys, Inc, 2010). The reaction rate constants used in FLUENT were selected based on the evaluation of Hanson and Salimian (Hanson et al., 1984). To calculate the NO concentration in turbulence combustion, the timeaveraged formation rate of NO must be computed at each computational cell. Methods of modeling the mean turbulent reaction rate of NO can be based on the presumed PDF model (Peng et al., 2009). FLUENT adopted the PDF model to predict the mass fraction of NO, and

NO NO 1 1 *w wP d* ( )

NO NO 1 2 1 2 *w w P dd* ( , 

represent the PDFs of the variables. In this study, the mass fraction of O2 and *T* were

mechanism includes intermediate species such as O, N, H, and OH. The mass fractions of these species can be obtained by using the proposed technique. Thus, only the reaction rate of NO was computed by the presumed PDF model, while those of the intermediate species

A H2-air non-premixed flame was simulated to verify the accuracy of the proposed technique. The numerical conditions were based on Barrow's experimental data (Barrow et al., 2003a). Figure 2 shows the sketch of the computational domain used in this study. H2 *<sup>Y</sup>*

where *w*NO represents the instantaneous reaction rate of NO, and 1 *P*( )

 

<sup>2</sup> , respectively. As shown in Equations (23)–(25), the extended Zeldovich

*J*

The equilibrium value *YJeq* was obtained from reading the lookup table using *ZH*

**2.3 NO prediction using the extended Zeldovich mechanism** 

, *<sup>J</sup> J r*

. Then, a simple multi-linear interpolation was used to extract information from the

*hX h* (22)

O N NO N, <sup>2</sup> (23)

N O O NO, <sup>2</sup> (24)

, (26)

) , (27)

 and 1 2 *P*(,) 

, *ZO*

 , *ZN* ,

fuels in this study. If we assumed that the reaction rates of fuels were determined with an overall reaction such as 1(kg)fuel (k *st* g)oxidant (1 )(k *st* g)products , \* *<sup>Y</sup>*fu can be modeled based on the expression of the overall reaction rate as

$$Y\_{\rm fu}^{\*} = \tilde{Y}\_{\rm fu} - \Delta Y\_{\rm fu} \tag{17}$$

$$
\Delta Y\_{\rm fu} = A\_{\rm fu} \min(\tilde{Y}\_{\rm fu}, \frac{\tilde{Y}\_{\rm ox}}{\rm st}). \tag{18}
$$

The subscripts fu and ox represent the fuels and oxidants, respectively. *st* represents the theoretical weight of the oxidants, and *A*fu is the mixing parameter of fuels set to 1.0 (equal to the value in the original EDC model).

#### **2.2 Chemical equilibrium method and reading lookup table**

When preparing the lookup table, we developed a different program from the UDF based on the chemical equilibrium method, by minimization of the Gibbs free energy. This method assumes that the equilibrium composition of a chemical system can be determined by minimizing the Gibbs free energy. Under the minimum conditions of Gibbs free energy, the following simultaneous nonlinear equations are obtained:

$$\left[X\_0 \sum\_{I} a\_{I\parallel} \exp\left[-\sum\_{K} \frac{a\_{K\parallel} \beta\_I}{RT} - \left(\frac{G\_I}{RT} + \ln p\right)\right] \right] = d\_{I\prime} \tag{19}$$

$$\sum\_{I} \exp\left[-\sum\_{K} \frac{a\_{Kj}\beta\_{I}}{RT} - \left(\frac{G\_{I}}{RT} + \ln p\right)\right] = 1,\tag{20}$$

where *X*0 represents the molar concentration, *IJ a* the number of elemental species *I* included in chemical species *J* , *KJ a* the number of elemental species *K* included in chemical species *J* , *<sup>I</sup>* the negative chemical potential only when the system is in chemical equilibrium, *R* the gas constant, *GJ* the Gibbs free energy of chemical species *J* , *p* the pressure, and *<sup>I</sup> d* the concentration of elemental species *I* .

When *<sup>I</sup> d* and *T* were assigned as constraints, *<sup>I</sup>* and *X*0 were obtained by solving the simultaneous nonlinear Equations (19) and (20). In this study, the computational code KINSOL (Collier et al., 2009) was used to solve these equations, and the concentration of each species at given *T* was obtainable by

$$X\_{I} = X\_{0} \exp - \left[ \sum\_{I} \frac{a\_{I\parallel} \beta\_{I}}{RT} - \left( \frac{G\_{I}}{RT} + \ln p \right) \right] \tag{21}$$

where *XJ* represents the molar concentration of chemical species *J* , and thermodynamic properties such as *GJ* and *<sup>J</sup> h* were computed as functions of the temperature (Mcbride et al., 2002). According to this procedure, the concentration of the chemical species *J* at given *T* was acquired. However, the additional constraint that the enthalpy was constant before and after combustion was required so that the chemical equilibrium composition in combustion could be found. Thus, *T* must be found under the constraint

fuels in this study. If we assumed that the reaction rates of fuels were determined with an overall reaction such as 1(kg)fuel (k *st* g)oxidant (1 )(k *st* g)products , \* *<sup>Y</sup>*fu can be

ox fu fu fu min( , ). *<sup>Y</sup> YA Y*

The subscripts fu and ox represent the fuels and oxidants, respectively. *st* represents the theoretical weight of the oxidants, and *A*fu is the mixing parameter of fuels set to 1.0 (equal

When preparing the lookup table, we developed a different program from the UDF based on the chemical equilibrium method, by minimization of the Gibbs free energy. This method assumes that the equilibrium composition of a chemical system can be determined by minimizing the Gibbs free energy. Under the minimum conditions of Gibbs free energy, the

> <sup>0</sup> exp ln , *KJ I J IJ I*

> > e , xp ln 1 *KJ I J*

*<sup>p</sup> RT RT* 

*X a p d RT RT*

*a G*

where *X*0 represents the molar concentration, *IJ a* the number of elemental species *I* included in chemical species *J* , *KJ a* the number of elemental species *K* included in

equilibrium, *R* the gas constant, *GJ* the Gibbs free energy of chemical species *J* , *p* the

simultaneous nonlinear Equations (19) and (20). In this study, the computational code KINSOL (Collier et al., 2009) was used to solve these equations, and the concentration of each

*I*

combustion could be found. Thus, *T* must be found under the constraint

*X X <sup>p</sup> RT RT* 

where *XJ* represents the molar concentration of chemical species *J* , and thermodynamic properties such as *GJ* and *<sup>J</sup> h* were computed as functions of the temperature (Mcbride et al., 2002). According to this procedure, the concentration of the chemical species *J* at given *T* was acquired. However, the additional constraint that the enthalpy was constant before and after combustion was required so that the chemical equilibrium composition in

0exp ln , *IJ I J*

*a G*

(21)

  *G*

(19)

(20)

*<sup>I</sup>* and *X*0 were obtained by solving the

*<sup>I</sup>* the negative chemical potential only when the system is in chemical

*st*

\* *Y Y* fu fu fu *Y* (17)

(18)

modeled based on the expression of the overall reaction rate as

**2.2 Chemical equilibrium method and reading lookup table** 

following simultaneous nonlinear equations are obtained:

pressure, and *<sup>I</sup> d* the concentration of elemental species *I* .

*J*

*J K*

*J K*

to the value in the original EDC model).

chemical species *J* ,

species at given *T* was obtainable by

When *<sup>I</sup> d* and *T* were assigned as constraints,

$$\sum\_{l} h\_{l} X\_{l} = h\_{r},\tag{22}$$

where *<sup>J</sup> h* represents the enthalpy of chemical species *J* , and *<sup>r</sup> h* represents the enthalpy of reactants (before combustion). When we found *T* using Equation (22), the obtained results were the concentrations after combustion.

The equilibrium value *YJeq* was obtained from reading the lookup table using *ZH* , *ZO* , *ZN* , and *<sup>a</sup> h* . Then, a simple multi-linear interpolation was used to extract information from the lookup table.

#### **2.3 NO prediction using the extended Zeldovich mechanism**

FLUENT version 13.0 has a function for computing the amount of the thermal NOX using the extended Zeldovich mechanism as

$$\text{O} + \text{N}\_2 \Leftrightarrow \text{NO} + \text{N}\_\prime \tag{23}$$

$$\text{N} + \text{O}\_2 \Leftrightarrow \text{O} + \text{NO}\_2 \tag{24}$$

$$\text{N} + \text{OH} \Leftrightarrow \text{H} + \text{NO} \tag{25}$$

Because we used the basic function in FLUENT, we describe only the outline of the method to estimate the amount of NO here (Warnatz et al., 2006; Ansys, Inc, 2010). The reaction rate constants used in FLUENT were selected based on the evaluation of Hanson and Salimian (Hanson et al., 1984). To calculate the NO concentration in turbulence combustion, the timeaveraged formation rate of NO must be computed at each computational cell. Methods of modeling the mean turbulent reaction rate of NO can be based on the presumed PDF model (Peng et al., 2009). FLUENT adopted the PDF model to predict the mass fraction of NO, and the reaction rate of the NO was expressed as

$$\overline{w}\_{\rm NO} = \int w\_{\rm NO} P(\phi\_1) d\phi\_1 \,\tag{26}$$

alternately, the PDF was related to two variables:

NO NO 1 2 1 2 *w w P dd* ( , ) , (27)

where *w*NO represents the instantaneous reaction rate of NO, and 1 *P*( ) and 1 2 *P*(,) represent the PDFs of the variables. In this study, the mass fraction of O2 and *T* were chosen as 1 and <sup>2</sup> , respectively. As shown in Equations (23)–(25), the extended Zeldovich mechanism includes intermediate species such as O, N, H, and OH. The mass fractions of these species can be obtained by using the proposed technique. Thus, only the reaction rate of NO was computed by the presumed PDF model, while those of the intermediate species were taken into account using the proposed technique.

#### **2.4 Numerical conditions**

A H2-air non-premixed flame was simulated to verify the accuracy of the proposed technique. The numerical conditions were based on Barrow's experimental data (Barrow et al., 2003a). Figure 2 shows the sketch of the computational domain used in this study. H2 *<sup>Y</sup>*

Simulation of H2-Air Non-Premixed Flame

0.7,

 *<sup>h</sup>* and 0.7 

based on the EDC model.

0.06

0.6 / *x L* were slightly larger than the experimental data.

0.12

0.18

Mixture fraction

Fig. 3. Axial mixture fraction

0.24

0.3

Using Combustion Simulation Technique to Reduce Chemical Mechanisms 365

Barrow's detailed mechanism (Barrow et al., 1999); and H2 *A* 100, the computational data obtained by the proposed technique with H2 *A* 100 . Although H2 *A* 100 was seemingly large, choosing this value means that H2 is mixed with O2 and is immediately consumed. Fig. 3 shows the axial profile of the mixture fraction. Because the parameters 2 *C* 1.80,

fraction obtained by the experiment, the computational results were close to the experimental data at 0.2 / *x L* . However, as shown in Fig. 3, those at *x L* / 0.2 were slightly smaller than the experimental data. The mixture fraction was related to mass diffusion, which in turn depended on turbulence models. It can be seen from this result that mass diffusion was not similar to the experimental data at / 0.2 *x L* (near the nozzle), and the turbulence model would have to be modified to obtain the closer result. However, the objective in this section is to confirm the validity of the proposed technique. To that end, we performed the simulation using the EDC model for reference. The EDC model has been widely applied to the simulation of combustion equipment (Stefanidis et al., 2006; Minotti et al., 2010), and the results obtained using the EDC model have been considered to be reliable for the engineered approach. The prediction of the proposed technique must approach that of the EDC model, because the expression of the reaction rate in the proposed technique was

<sup>0</sup> 0.2 0.4 0.6 0.8 <sup>1</sup> <sup>0</sup>

*x/L*

The axial temperature *T* and mole fraction of H2O are shown in Fig. 4. The mixture fractions of the computational results are in good agreement with each other, although those at

Fig. 5 shows the axial mole fraction of H2 and O2. Both the computational results of H2 and O2 were lower than those of the experimental data. It was found that the simulations slightly overestimated the consumption of H2 and O2 as compared with those of the experiment. The mole fraction of H2 at 0.3 / 0.6 *x L* obtained by the proposed technique was slightly larger than that of the EDC model. The reaction rate of H2 was determined by the overall reaction, choosing H2 *A* 100 in the case of the proposed technique; the rate was computed

EXP

EDC

Proposed (AH2=100)

*<sup>h</sup>* of the Reynolds stress model were chosen to match the axial mixture

in the fuel was 1.0; O2 *<sup>Y</sup>* , N2 *<sup>Y</sup>* , and H O2 *<sup>Y</sup>* in air were 0.231, 0.761, and 0.00720, respectively. Since O2 *<sup>Y</sup>* , N2 *<sup>Y</sup>* , and H O2 *<sup>Y</sup>* belong to group b, � *Z*<sup>H</sup> 0.008 , � *Z*<sup>O</sup> 0.23772 , and � *Z*<sup>N</sup> 0.76148 were assigned in air. The fuel was injected from a nozzle with a radius of 3.75 mm. The velocity of the fuel and ambient air were 296 and 1.0 m/s, respectively. The temperature of the fuel and ambient air was set to 300 K. In all cases, a second-order upwind difference was applied to the governing equations of the momentum, mass fraction of chemical species *J* and elemental species *I* , and the energy equation. The two-dimensional axisymmetric solver in FLUENT was chosen using the axis indicated by the centerline in Fig. 2. The Reynolds stress model was chosen as the turbulence model; 2 *C* 1.80 , 0.7 *<sup>h</sup>* , and *<sup>m</sup>* 0.7 of the model parameters were used to adjust the axial profile of the mixture fraction of the experimental data. In this study, the detailed chemical mechanism of H2 suggested by Barrow (Barrow et al., 1999) was used to compute the reaction rates of the modified species in the proposed technique, and those of the chemical species taken into account by the EDC model.

Fig. 2. Sketch of computational domain

#### **3. Results**

#### **3.1 Prediction of mixture fraction, temperature, and mole fraction**

First, to validate the prediction accuracy of the proposed technique, we compared the obtained results with those of the EDC model and the experiment. The lookup table was prepared before the combustion simulation, and it included the following species: H2, O2, H2O, O, H, OH, HO2, H2O2, NO, HNO, N2O, N, NO2, and N2. In Figs. 3–18, *x* represents the axial distance from the exit of the nozzle; *r*, the radial distance from the centerline of the nozzle; *L* , the approximate visible flame length, which is equal to 675 mm in this flame; *D* , the diameter of the nozzle, which is equal to 3.75 mm in this flame; *T*, the temperature; EXP, the experimental data; EDC, the computational data predicted by the EDC model with

in the fuel was 1.0; O2 *<sup>Y</sup>* , N2 *<sup>Y</sup>* , and H O2 *<sup>Y</sup>* in air were 0.231, 0.761, and 0.00720, respectively.

*Z*<sup>N</sup> 0.76148 were assigned in air. The fuel was injected from a nozzle with a radius of 3.75 mm. The velocity of the fuel and ambient air were 296 and 1.0 m/s, respectively. The temperature of the fuel and ambient air was set to 300 K. In all cases, a second-order upwind difference was applied to the governing equations of the momentum, mass fraction of chemical species *J* and elemental species *I* , and the energy equation. The two-dimensional axisymmetric solver in FLUENT was chosen using the axis indicated by the centerline in Fig. 2. The Reynolds stress model was chosen as the turbulence model; 2 *C* 1.80 , 0.7

 *<sup>m</sup>* 0.7 of the model parameters were used to adjust the axial profile of the mixture fraction of the experimental data. In this study, the detailed chemical mechanism of H2 suggested by Barrow (Barrow et al., 1999) was used to compute the reaction rates of the modified species in the proposed technique, and those of the chemical species taken into

*Z*<sup>H</sup> 0.008 , �

1350

Centerline

Atmosphere

*Z*<sup>O</sup> 0.23772 , and

*r*

Outlet

*x*

*<sup>h</sup>* , and

Since O2 *<sup>Y</sup>* , N2 *<sup>Y</sup>* , and H O2 *<sup>Y</sup>* belong to group b, �

 **Air**

160

1.875

**3.1 Prediction of mixture fraction, temperature, and mole fraction** 

First, to validate the prediction accuracy of the proposed technique, we compared the obtained results with those of the EDC model and the experiment. The lookup table was prepared before the combustion simulation, and it included the following species: H2, O2, H2O, O, H, OH, HO2, H2O2, NO, HNO, N2O, N, NO2, and N2. In Figs. 3–18, *x* represents the axial distance from the exit of the nozzle; *r*, the radial distance from the centerline of the nozzle; *L* , the approximate visible flame length, which is equal to 675 mm in this flame; *D* , the diameter of the nozzle, which is equal to 3.75 mm in this flame; *T*, the temperature; EXP, the experimental data; EDC, the computational data predicted by the EDC model with

**Air** 250

**1e-5 m/s**

150

�

account by the EDC model.

**O2:0.231 N2:0.761**

**Fuel**

**Fuel**

**H2:1.0**

**3. Results** 

**296 m/s**

**Mass fraction**

Fig. 2. Sketch of computational domain

**H2O:0.00720**

**1.0 m/s**

Wall

**Air**

Barrow's detailed mechanism (Barrow et al., 1999); and H2 *A* 100, the computational data obtained by the proposed technique with H2 *A* 100 . Although H2 *A* 100 was seemingly large, choosing this value means that H2 is mixed with O2 and is immediately consumed.

Fig. 3 shows the axial profile of the mixture fraction. Because the parameters 2 *C* 1.80, 0.7, *<sup>h</sup>* and 0.7 *<sup>h</sup>* of the Reynolds stress model were chosen to match the axial mixture fraction obtained by the experiment, the computational results were close to the experimental data at 0.2 / *x L* . However, as shown in Fig. 3, those at *x L* / 0.2 were slightly smaller than the experimental data. The mixture fraction was related to mass diffusion, which in turn depended on turbulence models. It can be seen from this result that mass diffusion was not similar to the experimental data at / 0.2 *x L* (near the nozzle), and the turbulence model would have to be modified to obtain the closer result. However, the objective in this section is to confirm the validity of the proposed technique. To that end, we performed the simulation using the EDC model for reference. The EDC model has been widely applied to the simulation of combustion equipment (Stefanidis et al., 2006; Minotti et al., 2010), and the results obtained using the EDC model have been considered to be reliable for the engineered approach. The prediction of the proposed technique must approach that of the EDC model, because the expression of the reaction rate in the proposed technique was based on the EDC model.

Fig. 3. Axial mixture fraction

The axial temperature *T* and mole fraction of H2O are shown in Fig. 4. The mixture fractions of the computational results are in good agreement with each other, although those at 0.6 / *x L* were slightly larger than the experimental data.

Fig. 5 shows the axial mole fraction of H2 and O2. Both the computational results of H2 and O2 were lower than those of the experimental data. It was found that the simulations slightly overestimated the consumption of H2 and O2 as compared with those of the experiment. The mole fraction of H2 at 0.3 / 0.6 *x L* obtained by the proposed technique was slightly larger than that of the EDC model. The reaction rate of H2 was determined by the overall reaction, choosing H2 *A* 100 in the case of the proposed technique; the rate was computed

Simulation of H2-Air Non-Premixed Flame

0.002

Fig. 6. Radial temperature and mole fraction of OH

0.004

0.006

Mole fraction

model or the proposed technique.

mainly from turbulence model accuracy.

0.1

0.2

Mixture fraction

Fig. 7. Radial mixture fraction at *x*/*L* = 1/8

0.3

0.008

0.01

Using Combustion Simulation Technique to Reduce Chemical Mechanisms 367

and we found that the mole fraction of OH obtained using the proposed technique was

0 0.2 0.4 0.6 0.8 1

*x/L*

Fig. 7 shows the radial mixture fraction at / 1 /8 *x L* . Although the results obtained by the EDC model and the proposed technique were in good agreement, they were not close to the experimental data. As mentioned above, the difference of the mixture fraction between the simulations and experiment was related to the turbulence model, and not to either the EDC

As we can observe from Fig. 8, the peaks of the temperature and the mole fraction of H2O predicted with the use of the EDC and proposed models were at approximately / 4 *r D* , while that of the experimental data was at approximately / 2.9 *r D* . The deference of coordinates of the peaks between the simulation data and the experimental data derived

<sup>0</sup> <sup>2</sup> <sup>4</sup> <sup>6</sup> <sup>8</sup> <sup>0</sup>

*r/D*

EXP

EDC

Proposed (AH2=100)

EXP

EDC

Proposed (AH2=100)

OH

qualitatively similar to that obtained by the EDC model and experimental data.

from the Barrow's detailed mechanism in the case of the EDC model. The overall reaction with choosing H2 *A* 100 in the proposed technique caused the difference of H2 at 0.3 / 0.6 *x L* .

Fig. 4. Axial temperature *T* and mole fraction of H2O

Fig. 5. Radial temperature and mole fraction of H2 and O2

The mole fraction of OH is plotted in Fig. 6. Generally, reaction equations including OH are needed to obtain the mole fraction of OH, but it was not required in the proposed technique. (However, the reaction equation of OH can be taken into account in the proposed technique using Equations (15)–(16) to improve the prediction accuracy.)

As shown in Fig. 6, the result of OH was comparatively small as compared with major species such as H2, O2, and H2O; the amount of the radical species, such as H, O, and OH, were generally very subtle. For this reason, it was important to estimate the tendency of yield and consumption of the mole fraction of OH. In Fig. 6, the coordinates of the peaks of the mole fraction of OH obtained by the simulations were similar to the experimental data,

from the Barrow's detailed mechanism in the case of the EDC model. The overall reaction with choosing H2 *A* 100 in the proposed technique caused the difference of H2 at

Proposed (AH2=100)

<sup>0</sup> 0.2 0.4 0.6 0.8 <sup>1</sup> <sup>0</sup>

*x/L*

EXP

EDC

<sup>0</sup> 0.2 0.4 0.6 0.8 <sup>1</sup> <sup>0</sup>

*x/L*

The mole fraction of OH is plotted in Fig. 6. Generally, reaction equations including OH are needed to obtain the mole fraction of OH, but it was not required in the proposed technique. (However, the reaction equation of OH can be taken into account in the proposed technique

As shown in Fig. 6, the result of OH was comparatively small as compared with major species such as H2, O2, and H2O; the amount of the radical species, such as H, O, and OH, were generally very subtle. For this reason, it was important to estimate the tendency of yield and consumption of the mole fraction of OH. In Fig. 6, the coordinates of the peaks of the mole fraction of OH obtained by the simulations were similar to the experimental data,

H2O

Proposed (AH2=100)

*T*

600

O2

1200

1800

*T* [K]

2400

3000

0.3 / 0.6 *x L* .

0.2

0.2

0.4

0.6

Mole fraction

0.8

1

Fig. 4. Axial temperature *T* and mole fraction of H2O

H2

Fig. 5. Radial temperature and mole fraction of H2 and O2

using Equations (15)–(16) to improve the prediction accuracy.)

0.4

0.6

Mole fraction

0.8

1

EXP

EDC

and we found that the mole fraction of OH obtained using the proposed technique was qualitatively similar to that obtained by the EDC model and experimental data.

Fig. 6. Radial temperature and mole fraction of OH

Fig. 7 shows the radial mixture fraction at / 1 /8 *x L* . Although the results obtained by the EDC model and the proposed technique were in good agreement, they were not close to the experimental data. As mentioned above, the difference of the mixture fraction between the simulations and experiment was related to the turbulence model, and not to either the EDC model or the proposed technique.

As we can observe from Fig. 8, the peaks of the temperature and the mole fraction of H2O predicted with the use of the EDC and proposed models were at approximately / 4 *r D* , while that of the experimental data was at approximately / 2.9 *r D* . The deference of coordinates of the peaks between the simulation data and the experimental data derived mainly from turbulence model accuracy.

Fig. 7. Radial mixture fraction at *x*/*L* = 1/8

Simulation of H2-Air Non-Premixed Flame

mole fraction of H2O are in good agreement.

0.002

Fig. 10. Radial mole fraction of OH at x/L = 1/8

0.01

0.02

0.03

Mixture fraction

Fig. 11. Radial mixture fraction at *x*/*L* = 1/2

0.04

0.05

0.06

0.004

Mole fraction

0.006

0.008

Using Combustion Simulation Technique to Reduce Chemical Mechanisms 369

The mixture fraction and temperature as well as mole fraction of H2O are shown in Figs. 11 and 12, respectively. Both the computational data of the mixture fraction, temperature, and

0 2 4 6 8

*r/D*

EXP

EDC

Proposed (AH2=100)

<sup>0</sup> <sup>4</sup> <sup>8</sup> <sup>12</sup> <sup>16</sup> <sup>20</sup> <sup>0</sup>

*r/D*

Fig. 13 shows the radial mole fraction of H2 and O2 at *x* /L 1/2 . Although the mole fractions of O2 obtained by the proposed model and EDC model were almost similar, the result of H2 computed using the proposed technique was not similar to that of the EDC model. This was because the proposed technique used H2 *A* 100 to estimate the reaction rate of H2, while the reaction rate was computed using the detailed chemical mechanism in the EDC model. However, the difference of the results between the proposed technique and

EXP

EDC

OH

Proposed (AH2=100)

Fig. 8. Radial temperature and mole fraction of H2O at *x*/*L* = 1/8

Fig. 9 shows the radial mole fraction of H2 and O2. The mole fraction of O2 obtained by the proposed technique was lower than that predicted by the EDC model. In the proposed technique, the mole fraction of O2 was obtained from the lookup table, which was prepared using the chemical equilibrium method. For this reason, consumption of O2 was sometimes overestimated at points where the reaction was in the state of strong nonequilibrium.

Fig. 9. Radial mole fraction of H2 and O2 at *x*/*L* = 1/8

As shown in Fig. 10, the coordinate values of the peak of the mole fraction of OH obtained by the simulations and experimental data are located at approximately / 4.1 *r D* and *r D*/ 2.7 , respectively. These values were slightly different from what was seen in the experiment. The coordinate of the peak mole fraction of OH was close to that of the mole fraction of H2O and *T*. It can be seen from this tendency that OH was the indicative species of extinction and ignition phenomena.

EXP

EDC

*T*

Proposed (AH2=100)

<sup>0</sup> <sup>2</sup> <sup>4</sup> <sup>6</sup> <sup>8</sup> <sup>0</sup>

*r/D*

Fig. 9 shows the radial mole fraction of H2 and O2. The mole fraction of O2 obtained by the proposed technique was lower than that predicted by the EDC model. In the proposed technique, the mole fraction of O2 was obtained from the lookup table, which was prepared using the chemical equilibrium method. For this reason, consumption of O2 was sometimes overestimated at points where the reaction was in the state of strong

EXP

EDC

<sup>0</sup> <sup>2</sup> <sup>4</sup> <sup>6</sup> <sup>8</sup> <sup>0</sup>

*r/D*

As shown in Fig. 10, the coordinate values of the peak of the mole fraction of OH obtained by the simulations and experimental data are located at approximately / 4.1 *r D* and *r D*/ 2.7 , respectively. These values were slightly different from what was seen in the experiment. The coordinate of the peak mole fraction of OH was close to that of the mole fraction of H2O and *T*. It can be seen from this tendency that OH was the indicative species

O2

Proposed (AH2=100)

H2O

Fig. 8. Radial temperature and mole fraction of H2O at *x*/*L* = 1/8

H2

500

1000

1500

*T* [K]

2000

2500

3000

0.2

0.2

of extinction and ignition phenomena.

Fig. 9. Radial mole fraction of H2 and O2 at *x*/*L* = 1/8

0.4

0.6

Mole fraction

0.8

1

0.4

0.6

Mole fraction

nonequilibrium.

0.8

1

The mixture fraction and temperature as well as mole fraction of H2O are shown in Figs. 11 and 12, respectively. Both the computational data of the mixture fraction, temperature, and mole fraction of H2O are in good agreement.

Fig. 10. Radial mole fraction of OH at x/L = 1/8

Fig. 11. Radial mixture fraction at *x*/*L* = 1/2

Fig. 13 shows the radial mole fraction of H2 and O2 at *x* /L 1/2 . Although the mole fractions of O2 obtained by the proposed model and EDC model were almost similar, the result of H2 computed using the proposed technique was not similar to that of the EDC model. This was because the proposed technique used H2 *A* 100 to estimate the reaction rate of H2, while the reaction rate was computed using the detailed chemical mechanism in the EDC model. However, the difference of the results between the proposed technique and

Simulation of H2-Air Non-Premixed Flame

0.002

Fig. 14. Radial mole fraction of OH at *x*/*L* = 1/2

0.01

0.2

0.4

Mole fraction

H2O

Fig. 16. Radial temperature and mole fraction of H2O at x/*L* = ¾

*T*

0.6

0.8

0.02

0.03

Mixture fraction

Fig. 15. Radial mixture fraction at *x*/*L*= 3/4

0.04

0.05

0.004

0.006

Mole fraction

0.008

0.01

Using Combustion Simulation Technique to Reduce Chemical Mechanisms 371

OH

0 4 8 12 16 20

*r/D*

EXP

EDC

Proposed (AH2=100)

<sup>0</sup> <sup>3</sup> <sup>6</sup> <sup>9</sup> <sup>12</sup> <sup>15</sup> <sup>0</sup>

*r/D*

EXP

EDC

Proposed (AH2=100)

<sup>0</sup> <sup>3</sup> <sup>6</sup> <sup>9</sup> <sup>12</sup> <sup>15</sup> <sup>0</sup>

*r/D*

600

1200

1800

*T* [K]

2400

3000

EXP

EDC

Proposed (AH2=100)

EDC model were small as compared to the maximum mole fraction of H2. Basically, the maximum mole fraction of H2 was equal to 1.0 taken at / 0 *x L* ; the difference at *r D*/ 0.0 between both computational results was about 0.03.

Fig. 12. Radial temperature and mole fraction of H2O at *x*/*L* = 1/2

As shown in Fig. 14, the peaks of the mole fractions of OH obtained by the simulations and that of the experiment were located at approximately / 7.6 *r D* and / 6.4 *r D* , respectively. The coordinates of the peaks obtained by the simulations were similar to each other, but they were slightly different from that of the experiment.

As shown in Figs. 15, 16, and 17, computational results for EDC and proposed models were in good agreement, but were either slightly greater or less than the experimental data, which showed a similar tendency to that of the other area, such as the radical direction at *x L* / 1/8 and / 1 / 2. *x L*

Fig. 13. Radial mole fraction of H2 and O2 at *x*/*L* = 1/2

Fig. 14. Radial mole fraction of OH at *x*/*L* = 1/2

EDC model were small as compared to the maximum mole fraction of H2. Basically, the maximum mole fraction of H2 was equal to 1.0 taken at / 0 *x L* ; the difference at

EXP

EDC

Proposed (AH2=100)

<sup>0</sup> <sup>4</sup> <sup>8</sup> <sup>12</sup> <sup>16</sup> <sup>20</sup> <sup>0</sup>

*r/D*

As shown in Fig. 14, the peaks of the mole fractions of OH obtained by the simulations and that of the experiment were located at approximately / 7.6 *r D* and / 6.4 *r D* , respectively. The coordinates of the peaks obtained by the simulations were similar to each

As shown in Figs. 15, 16, and 17, computational results for EDC and proposed models were in good agreement, but were either slightly greater or less than the experimental data, which showed a similar tendency to that of the other area, such as the radical direction at

<sup>0</sup> <sup>4</sup> <sup>8</sup> <sup>12</sup> <sup>16</sup> <sup>20</sup> <sup>0</sup>

*r/D*

O2

Proposed (AH2=100)

EXP

EDC

500

1000

1500

*T* [K]

2000

2500

3000

*r D*/ 0.0 between both computational results was about 0.03.

*T*

Fig. 12. Radial temperature and mole fraction of H2O at *x*/*L* = 1/2

other, but they were slightly different from that of the experiment.

0.1

0.2

0.3

Mole fraction

*x L* / 1/8 and / 1 / 2. *x L*

0.05

Fig. 13. Radial mole fraction of H2 and O2 at *x*/*L* = 1/2

0.1

0.15

Mole fraction

0.2

H2

0.25

0.3

H2O

0.4

0.5

0.6

Fig. 15. Radial mixture fraction at *x*/*L*= 3/4

Fig. 16. Radial temperature and mole fraction of H2O at x/*L* = ¾

Simulation of H2-Air Non-Premixed Flame

OH modified using Equations (15)–(16).

yield and consumption of NO.

0.0004

0.0008

0.0012

Mole fraction

0.0016

0.002

**3.2 Prediction of NO** 

Using Combustion Simulation Technique to Reduce Chemical Mechanisms 373

We show the application of the proposed technique to the simulation of NO in this section. In Figs. 19–30, "Proposed" indicates computational data obtained by the proposed technique with H2 *A* 100 ; "Proposed (H, O, OH)" indicates the computational data obtained by the proposed technique with H2 *A* 100 and reaction rates of H, O, and OH modified using Equations (15)–(16). Other symbols in Figs. 19–30 are as they were in the previous section. Two lookup tables were prepared before the combustion simulation. One table includes the species H2, O2, H2O, O, H, OH, HO2, H2O2, NO, HNO, N2O, N, NO2, and N2; it was used for the case of the proposed technique with H2 *A* 100 . The other table includes the species H2, O2, H2O, HO2, H2O2, NO, HNO, N2O, N, NO2, and N2; it was prepared for the case of the proposed technique H2 *A* 100 and reaction rates of H, O, and

As can be seen in Fig. 19, the peaks of the mole fraction of NO obtained by the simulation and that of the experiment were located at approximately / 0.82 *x L* . Those of the EDC model and proposed technique with modified reaction rates of H, O, and OH were bigger than the experimental data, while the peak obtained using the proposed technique without modified reaction rates was similar to that of the experiment. On the other hand, as shown in Fig. 20, the mole fraction of NO obtained by the proposed technique without the modified reaction rates was not similar to the experimental data. Thus, we did not determine which cases were most accurate in these simulations. As can be observed from Figs. 19 and 20, the amount of NO was very small as compared to the major species such as H2, O2, and H2O. Therefore, the most important aim of the NO prediction was to simulate the tendency of

0 0.4 0.8 1.2 1.6 2

*x/L*

Fig. 19. Axial mole fraction of NO predicted with the extended Zeldovich mechanism

at roughly / 6 *r D* , and the mole fraction of NO decreased at 6 / 15 *r D* .

In Figs. 20 and 21, the coordinates of the peaks obtained by both the simulations were in good agreement. Moreover, as shown in Fig. 22, the peaks of the mole fractions were located

 EXP Proposed

EDC

Proposed (H, O, OH)

Fig. 17. Radial mole fraction of H2 and O2 at *x*/*L* = 3/4

Fig. 18. Radial mole fraction of OH at *x*/*L* = 3/4

The radial mole fractions of OH at / 3 / 4 *x L* are depicted in Fig. 18. The computational results were overestimated as compared with the experimental data, and the mole fraction of OH at 9 / *r D* obtained using the proposed technique was similar to that obtained using the EDC model.

In this portion of the experiment, we found that the results of the mixture fraction, major species, and temperature obtained by the proposed technique were roughly similar to those predicted by the EDC model. The mole fraction of OH was not quantitatively similar to that of the EDC model; however, the coordinate values of the peaks of mole fraction of OH predicted by the proposed technique and EDC model were in good agreement with one another. Thus, the tendency of yield and consumption of OH was predicted with the proposed technique, which did not use the chemical equation of OH.

The computations show slight differences in the experimental data, even when the EDC model was applied to the combustion simulation. Moreover, we found that it is important to accurately predict the non-premixed flame using not only the combustion model but also the turbulence model. For this reason, further improvement of the combustion and turbulence models is a significant task for combustion simulation.

#### **3.2 Prediction of NO**

372 Fluid Dynamics, Computational Modeling and Applications

<sup>0</sup> <sup>3</sup> <sup>6</sup> <sup>9</sup> <sup>12</sup> <sup>15</sup> <sup>0</sup>

*r/D*

EXP

EDC

0 3 6 9 12 15

*r/D*

The radial mole fractions of OH at / 3 / 4 *x L* are depicted in Fig. 18. The computational results were overestimated as compared with the experimental data, and the mole fraction of OH at 9 / *r D* obtained using the proposed technique was similar to that obtained using

In this portion of the experiment, we found that the results of the mixture fraction, major species, and temperature obtained by the proposed technique were roughly similar to those predicted by the EDC model. The mole fraction of OH was not quantitatively similar to that of the EDC model; however, the coordinate values of the peaks of mole fraction of OH predicted by the proposed technique and EDC model were in good agreement with one another. Thus, the tendency of yield and consumption of OH was predicted with the

The computations show slight differences in the experimental data, even when the EDC model was applied to the combustion simulation. Moreover, we found that it is important to accurately predict the non-premixed flame using not only the combustion model but also the turbulence model. For this reason, further improvement of the combustion and

proposed technique, which did not use the chemical equation of OH.

turbulence models is a significant task for combustion simulation.

O2

EXP

EDC

H2

Proposed (AH2=100)

Proposed (AH2=100)

0.06

Fig. 17. Radial mole fraction of H2 and O2 at *x*/*L* = 3/4

0.002

Fig. 18. Radial mole fraction of OH at *x*/*L* = 3/4

the EDC model.

0.004

Mole fraction

0.006

0.008

0.12

0.18

Mole Fraction

0.24

0.3

We show the application of the proposed technique to the simulation of NO in this section. In Figs. 19–30, "Proposed" indicates computational data obtained by the proposed technique with H2 *A* 100 ; "Proposed (H, O, OH)" indicates the computational data obtained by the proposed technique with H2 *A* 100 and reaction rates of H, O, and OH modified using Equations (15)–(16). Other symbols in Figs. 19–30 are as they were in the previous section. Two lookup tables were prepared before the combustion simulation. One table includes the species H2, O2, H2O, O, H, OH, HO2, H2O2, NO, HNO, N2O, N, NO2, and N2; it was used for the case of the proposed technique with H2 *A* 100 . The other table includes the species H2, O2, H2O, HO2, H2O2, NO, HNO, N2O, N, NO2, and N2; it was prepared for the case of the proposed technique H2 *A* 100 and reaction rates of H, O, and OH modified using Equations (15)–(16).

As can be seen in Fig. 19, the peaks of the mole fraction of NO obtained by the simulation and that of the experiment were located at approximately / 0.82 *x L* . Those of the EDC model and proposed technique with modified reaction rates of H, O, and OH were bigger than the experimental data, while the peak obtained using the proposed technique without modified reaction rates was similar to that of the experiment. On the other hand, as shown in Fig. 20, the mole fraction of NO obtained by the proposed technique without the modified reaction rates was not similar to the experimental data. Thus, we did not determine which cases were most accurate in these simulations. As can be observed from Figs. 19 and 20, the amount of NO was very small as compared to the major species such as H2, O2, and H2O. Therefore, the most important aim of the NO prediction was to simulate the tendency of yield and consumption of NO.

Fig. 19. Axial mole fraction of NO predicted with the extended Zeldovich mechanism

In Figs. 20 and 21, the coordinates of the peaks obtained by both the simulations were in good agreement. Moreover, as shown in Fig. 22, the peaks of the mole fractions were located at roughly / 6 *r D* , and the mole fraction of NO decreased at 6 / 15 *r D* .

Simulation of H2-Air Non-Premixed Flame

more reactive because of the increasing *T* of the EDC model*.* 

0.0003

was also larger than that computed by the EDC model.

0.0003

0.0006

0.0009

Mole fraction

Fig. 23. Axial mole fraction of O

0.0012

0.0015

computing the reaction rates of H, O, and OH can be omitted in this flame.

Fig. 22. Radial mole fraction of NO at *x*/*L* = 3/4

0.0006

0.0009

Mole fraction

0.0012

0.0015

Using Combustion Simulation Technique to Reduce Chemical Mechanisms 375

computed by the proposed technique with modified reaction rates of H, O, and OH, as can be seen in Fig. 24. Thus, the first reaction of the extended Zeldovich mechanism became

0 3 6 9 12 15

*r/D*

0 0.4 0.8 1.2 1.6 2

*x/L*

Proposed

EDC

O

Proposed (H, O, OH)

On the other hand, as can be seen in Figs. 25 and 26, either the mole fraction of O or *T* obtained using the proposed technique with modified reaction rates of H, O, and OH were larger than those of the EDC model. Therefore, the radial mole fraction of NO at / 1 /8 *x L*

We found that the amount of NO was determined by the amount of O and *T*, and those of prediction accuracy were so important of obtaining the accurate mole fraction of NO. In the proposed technique, the reaction rates of the intermediate species such as H, O, and OH were computed, and the mole fractions can be obtained taking into account the chemical equations. Generally, when the reaction rates were modified using the chemical equations, the prediction accuracy of the chemical species was also increased (Fukumoto et al., 2010). However, even without the modified reaction rates, the proposed technique well simulated the mole fraction of O, and the improvement of the modified reaction rates was relatively small. Therefore, the

 EXP Proposed

EDC

Proposed (H, O, OH)

Fig. 20. Radial mole fraction of NO at *x*/*L*=1/8 predicted with the extended Zeldovich mechanism

It can be seen from Figs. 19–22 that the amounts of the NO prediction were not quantitatively correct; however, the peaks of the mole fractions obtained by the simulations were similar. Thus, we found that the tendency of yield and consumption of NO was simulated by all the cases in this part of the experiment.

The prediction accuracy of O is significant in the computation of the mass fraction of NO because the first reaction of the extended Zeldovich mechanism, expressed as O + N2 N + NO, is the rate-limiting step (Warnatz et al., 2006). Thus, the amount of NO increases with an increase in O. As shown in Fig. 23, the axial mole fraction of O obtained using the proposed technique with modified reaction rates of H, O, and OH was larger than that of the EDC model.

Fig. 21. Radial mole fraction of NO at *x*/*L* = 1/2

However, in Fig. 19, the axial mole fraction of NO obtained by the proposed technique with modified reaction rates of H, O, and OH was smaller than that of the EDC model. This was because axial *T* computed by the EDC model was overestimated as compared to that

0.00003

simulated by all the cases in this part of the experiment.

0.0002

Fig. 21. Radial mole fraction of NO at *x*/*L* = 1/2

0.0004

0.0006

Mole fraction

0.0008

0.001

0.00006

0.00009

Mole fraction

mechanism

the EDC model.

0.00012

0.00015

0 2 4 6 8

*r/D*

0 4 8 12 16 20

*r/D*

However, in Fig. 19, the axial mole fraction of NO obtained by the proposed technique with modified reaction rates of H, O, and OH was smaller than that of the EDC model. This was because axial *T* computed by the EDC model was overestimated as compared to that

 EXP Proposed

EDC

NO

Proposed (H, O, OH)

Fig. 20. Radial mole fraction of NO at *x*/*L*=1/8 predicted with the extended Zeldovich

It can be seen from Figs. 19–22 that the amounts of the NO prediction were not quantitatively correct; however, the peaks of the mole fractions obtained by the simulations were similar. Thus, we found that the tendency of yield and consumption of NO was

The prediction accuracy of O is significant in the computation of the mass fraction of NO because the first reaction of the extended Zeldovich mechanism, expressed as O + N2 N + NO, is the rate-limiting step (Warnatz et al., 2006). Thus, the amount of NO increases with an increase in O. As shown in Fig. 23, the axial mole fraction of O obtained using the proposed technique with modified reaction rates of H, O, and OH was larger than that of

 EXP Proposed

EDC

Proposed (H, O, OH)

NO

computed by the proposed technique with modified reaction rates of H, O, and OH, as can be seen in Fig. 24. Thus, the first reaction of the extended Zeldovich mechanism became more reactive because of the increasing *T* of the EDC model*.* 

Fig. 22. Radial mole fraction of NO at *x*/*L* = 3/4

On the other hand, as can be seen in Figs. 25 and 26, either the mole fraction of O or *T* obtained using the proposed technique with modified reaction rates of H, O, and OH were larger than those of the EDC model. Therefore, the radial mole fraction of NO at / 1 /8 *x L* was also larger than that computed by the EDC model.

We found that the amount of NO was determined by the amount of O and *T*, and those of prediction accuracy were so important of obtaining the accurate mole fraction of NO. In the proposed technique, the reaction rates of the intermediate species such as H, O, and OH were computed, and the mole fractions can be obtained taking into account the chemical equations. Generally, when the reaction rates were modified using the chemical equations, the prediction accuracy of the chemical species was also increased (Fukumoto et al., 2010). However, even without the modified reaction rates, the proposed technique well simulated the mole fraction of O, and the improvement of the modified reaction rates was relatively small. Therefore, the computing the reaction rates of H, O, and OH can be omitted in this flame.

Fig. 23. Axial mole fraction of O

Simulation of H2-Air Non-Premixed Flame

0.0005

2100

0.0005

0.001

Mole fraction

Fig. 29. Radial mole fraction of O at *x*/*L*=3/4

0.0015

0.002

2200

*T* [K]

2300

2400

2500

0.001

Mole fraction

Fig. 27. Radial mole fraction of O at *x*/*L*=1/2

0.0015

0.002

Using Combustion Simulation Technique to Reduce Chemical Mechanisms 377

0 4 8 12 16 20

*r/D*

Proposed

EDC

Proposed (H, O, OH)

<sup>0</sup> <sup>4</sup> <sup>8</sup> <sup>12</sup> <sup>16</sup> <sup>20</sup> <sup>2000</sup>

Fig. 28. Radial mole fraction of *T* at *x*/*L*=1/2 ranging from 2000 K to 2500 K

*r/D*

Proposed

EDC

O

Proposed (H, O, OH)

0 3 6 9 12 15

*r/D*

T

Proposed

EDC

O

Proposed (H, O, OH)

Fig. 24. Axial *T* ranging from 2000 K to 2800 K

Fig. 25. Radial mole fraction of O at *x*/*L* = 1/8

Fig. 26. Radial mole fraction of *T* at *x*/*L* = 1/8 ranging from 1800 K to 2400 K

Fig. 27. Radial mole fraction of O at *x*/*L*=1/2

<sup>0</sup> 0.4 0.8 1.2 1.6 <sup>2</sup> <sup>2000</sup>

*x/L*

Proposed

O

EDC

Proposed (H, O, OH)

0 2 4 6 8

*r/D*

T

Proposed

EDC

Proposed (H, O, OH)

<sup>0</sup> <sup>2</sup> <sup>4</sup> <sup>6</sup> <sup>8</sup> <sup>1800</sup>

Fig. 26. Radial mole fraction of *T* at *x*/*L* = 1/8 ranging from 1800 K to 2400 K

*r/D*

T

Proposed

EDC

Proposed (H, O, OH)

2200

Fig. 24. Axial *T* ranging from 2000 K to 2800 K

0.0006

1900

2000

2100

2200

2300

2400

Fig. 25. Radial mole fraction of O at *x*/*L* = 1/8

0.0012

0.0018

Mole fraction

0.0024

0.003

2400

2600

2800

Fig. 28. Radial mole fraction of *T* at *x*/*L*=1/2 ranging from 2000 K to 2500 K

Fig. 29. Radial mole fraction of O at *x*/*L*=3/4

Simulation of H2-Air Non-Premixed Flame

0

We draw the following conclusions:

http://www.ansys.com/

equation of OH.

**5. References** 

Fig. 31. Comparison of computational time

species, so that the computational time can be minimized.

fraction, and *T* was similar to that of the EDC model.

technique with the extended Zeldovich mechanism.

2010), pp. 11324–11331, ISSN: 0360-3199

2

4

Time [s]

**4. Conclusion** 

6

8

10

Using Combustion Simulation Technique to Reduce Chemical Mechanisms 379

68.5 %

Proposed Proposed (H, O, OH) EDC

In this chapter, we presented a combustion simulation technique based on a chemical equilibrium method combined with the EDC model. The proposed technique makes it easy to construct a reduced mechanism according to the accuracy requirements of the chemical

The prediction accuracy of the proposed technique for the mixture fraction, major mole

 The coordinates of the peak OH mole fractions predicted by the proposed technique and EDC model were in good agreement. Thus, the tendency of yield and consumption of OH can be predicted using the proposed technique, which did not use the chemical

The tendency of yield and consumption of NO was simulated by the proposed

 A lower number modified species led to lower computational time. However, the prediction accuracy also relied on the number of modified species included. Modified

Amir, M. & Sadegh, T. (2010). Effect of Hydrogen on Hydrogen-Methane Turbulent Non-

Barlow, R S. & Smith, N S A. (1999). Nitric Oxide Formation in Dilute Hydrogen Jet Flames:

Barlow, R S. (2003b). Radiation Models, In: *International Workshop on Measurement and* 

premixed Flame under MILD Condition. *Hydrogen Energy*, Vol. 35, No. 20, (October

Isolation of the Effects of Radiation and Turbulence-Chemistry Submodels. *Combustion and Flame*, Vol. 117, No. 1, (April 1999), pp. 4–31, ISSN: 0010–2180 Barlow, R S. (2003a). Sandia H2/He Flame Data-Release 2.0, In: *International Workshop on* 

*Measurement and Computation of Turbulent Non-premixed Flames*, June 2011, Available

*Computation of Turbulent Non-premixed Flames*, June 2011, Available from: <

species should be included according to the required prediction accuracy.

Ansys, Inc. (2011). Ansys FLUENT, (2010), Available from: http://www.ansys.com/ Ansys, Inc. (2010). *ANSYS FLUENT Theory Guide*, Ansys, Inc. Retrieved from

from: < http://www.sandia.gov/TNF/simplejet.html>

http://www.sandia.gov/TNF/radiation.html>

99.5 %

Fig. 30. Axial *T* at *x*/*L*=3/4 ranging from 2000 K to 2500 K

#### **3.3 Influence of number of modified species and computational time**

Due to the time required to compute reaction rates of these species, the computational time is dependent on the number of the chemical species in group a. Our proposed technique enables the shortening of these computational times. The purpose of this section is to discuss the prediction accuracy for the minor species and computational time required by the proposed technique. The numerical condition of the EDC model conformed to that presented in the previous sections, and the other two cases of the proposed technique were performed in this part of the experiment. In the first case of the proposed technique, the reaction rate of H2 was computed by the estimated by Equations (17)–(18). In the second case, those of H, O, and OH were modified, and the reaction rate of H2 was estimated.

Fig. 31 shows a comparison of the computational time in the case of the EDC model and two cases of the proposed technique. Each computational time per iteration was measured after the solution was converged and was measured during the 100 iterations of performing the reaction calculation and reading the lookup table. The CPU used in these simulations was the Intel Core i7-980X Processor Extreme Edition with 12M Cache, and 3.33 GHz. GFORTRAN 4.4.4 was chosen as a compiler, and CVODE 2.6 (Hindmarsh et al., 2009) was adopted to solve the ODEs of the reaction calculation. As can be seen in Fig. 31, the computational times for the EDC model—the two proposed techniques with modified H, O, and OH (second case) and with no modified species (first case)—were 8.99 s, 2.83 s, and 0.04 s, respectively. We found that the lower the number of included modified species, the fewer were the demands on lower computational time. On the other hand, the prediction accuracy also depended on the number of modified species included. Therefore, when the proposed technique is used, the number of modified species should be determined according to the required prediction accuracy. In the proposed technique, the computational time required increased with an increase in the number of modified species. However, the proposed technique has an advantage: We only chose the modified species and then the reduced mechanism was automatically built from the detailed mechanism.

Fig. 31. Comparison of computational time

## **4. Conclusion**

378 Fluid Dynamics, Computational Modeling and Applications

Proposed

EDC

Proposed (H, O, OH)

<sup>0</sup> <sup>3</sup> <sup>6</sup> <sup>9</sup> <sup>12</sup> <sup>15</sup> <sup>2000</sup>

Due to the time required to compute reaction rates of these species, the computational time is dependent on the number of the chemical species in group a. Our proposed technique enables the shortening of these computational times. The purpose of this section is to discuss the prediction accuracy for the minor species and computational time required by the proposed technique. The numerical condition of the EDC model conformed to that presented in the previous sections, and the other two cases of the proposed technique were performed in this part of the experiment. In the first case of the proposed technique, the reaction rate of H2 was computed by the estimated by Equations (17)–(18). In the second case, those of H, O, and OH were modified, and the reaction rate

Fig. 31 shows a comparison of the computational time in the case of the EDC model and two cases of the proposed technique. Each computational time per iteration was measured after the solution was converged and was measured during the 100 iterations of performing the reaction calculation and reading the lookup table. The CPU used in these simulations was the Intel Core i7-980X Processor Extreme Edition with 12M Cache, and 3.33 GHz. GFORTRAN 4.4.4 was chosen as a compiler, and CVODE 2.6 (Hindmarsh et al., 2009) was adopted to solve the ODEs of the reaction calculation. As can be seen in Fig. 31, the computational times for the EDC model—the two proposed techniques with modified H, O, and OH (second case) and with no modified species (first case)—were 8.99 s, 2.83 s, and 0.04 s, respectively. We found that the lower the number of included modified species, the fewer were the demands on lower computational time. On the other hand, the prediction accuracy also depended on the number of modified species included. Therefore, when the proposed technique is used, the number of modified species should be determined according to the required prediction accuracy. In the proposed technique, the computational time required increased with an increase in the number of modified species. However, the proposed technique has an advantage: We only chose the modified species and then the reduced mechanism was automatically built from the detailed

*r/D*

2100

Fig. 30. Axial *T* at *x*/*L*=3/4 ranging from 2000 K to 2500 K

**3.3 Influence of number of modified species and computational time** 

2200

*T* [K]

of H2 was estimated.

mechanism.

2300

2400

T

2500

In this chapter, we presented a combustion simulation technique based on a chemical equilibrium method combined with the EDC model. The proposed technique makes it easy to construct a reduced mechanism according to the accuracy requirements of the chemical species, so that the computational time can be minimized.

We draw the following conclusions:


#### **5. References**


**16** 

*Italy* 

Claudio Bruno

**Nuclear Propulsion** 

*DIMA, University of Rome (La Sapienza), Roma* 

Nuclear propulsion (NP) concepts go back to the very end of WW II. Scientists informed about the effects of the US atomic bomb thought of exploiting its energy release for applications like commercial electric power generation, but also rockets and space flight [Shepherd and Cleaver, 1948, 1949; Bussard and DeLauer, 1958]. However, space flight was still considered science fiction, and the military had to deal with more concrete things, like the Cold War. Thus, besides power generation, second stages of ICBM, submarine propulsion, long range and long duration airplanes and missiles became the focus of nuclear

It was the second-stage and airplane application that drove R&D in nuclear propulsion. With the advent of reliable ICBM (the Atlas missile) and lighter fission and thermonuclear warheads, a nuclear-powered second stage became no longer necessary. Airplane applications were found impractical: the Convair NB-36 required such a heavy lead shield for the crew that testing and operation were much restricted. Nuclear-powered missiles were easier to design, e.g., project PLUTO, but still far more complicated compared to conventional. The Soviets investigated airplanes and rockets powered by nuclear power as

The history of NP can be found in [Czysz and Bruno, 2009, Chapter 7; Lawrence, 2008; Lawrence et al, 1995; Gunn and Ehresman, 2003; Dewar, 2004] and will not be reported here. Basic technology is also discussed in the references above, in particular reactor design is in [Lawrence et al, 1995]. In the US investing in NP technology by NASA and LASL stopped in 1972, but was kept alive by the USAF, and to some extent also at Los Alamos Science Laboratories (LASL, now LANL) who focused on unconventional concepts. In the former Soviet Union, now Russia, research was slowed down, but never completely stopped,

The present work focuses rather synthetically on NP based on fission. (Fusion has been demonstrated, but the net power gain budget is still negative, and power generation is still far in the future, so NP based on this technology will not be discussed.) Note that only reluctantly space agencies or international study groups like ISECG [ISECG] are beginning to include *fission-based NP* in their exploration scenarios. The cost of R&D and experimenting with NP technology is in fact much higher than with conventional rocket propulsion, and until now MIR first, and then ISS, have absorbed most of the human exploration R&D moneys. After completing the ISS, space agencies are now discussing inhouse and among themselves what should be the next goal of *human* space exploration. This

because it was thought indispensable for future space missions.

**1. Introduction** 

energy applications.

well, and discarded them too.

Collier, A M. & Hindmarsh, A C. (2009). KINSOL version 2.6.0, In: *SUNDIALS web site*, 30 June 2011, Available from

https://computation.llnl.gov/casc/sundials/main.html/


https://computation.llnl.gov/casc/sundials/main.html/

