Twenty-Sixth Symposium (International) on Combustion/The Combustion Institute, 1996/pp. 915–921
TWO-DIMENSIONAL MATHEMATICAL MODELING OF LAMINAR, PREMIXED, METHANE-AIR COMBUSTION ON AN EXPERIMENTAL SLOT BURNER D. BRADLEY, P. H. GASKELL, K. C. KWAN and M. J. SCOTT Department of Mechanical Engineering University of Leeds Leeds LS2 9JT, UK
A two-dimensional model of the premixed laminar burning of methane-air has been applied successfully to combustion on a slot burner. Computing times could be excessive and a reduced reaction scheme of Peters was found to be satisfactory for nonrich combustion. CARS temperature measurements on such a burner were in good agreement with predictions, as were the flame heights. Heat loss by conduction from the hot gases to the burner surface are shown, both theoretically and experimentally, to be appreciable at low flow rates. Although conditions for flashback could be computed, it was not possible to compute those for liftoff. Nevertheless, the conditions that give rise to it can be deduced. The flame curvature and aerodynamic strain contributions to the flame stretch rate provide insights into the liftoff phenomena and these are quantified. Their relationship to local burning velocities enables Markstein lengths also to be quantified for each of the two contributions. As the flow rate increases, negative curvature stretch becomes dominant at the flame tip, whereas at the base of the flame the aerodynamic strain is dominant and negative. The values of Markstein length are such as to induce stability at the tip where the local burning velocity increases with flow rate, and instability at the base where the local burning velocity decreases with flow rate. These effects are reinforced by the diffusion of H2 toward the tip of the unreacted gases and by heat loss to the burner surface. Consequently, flame liftoff originates at the base.
Introduction The slot burner, with a steady, premixed, laminar flame is of practical importance in small boilers and space heaters. Key design parameters are the conditions for flashback and liftoff, and the relationship between the flow velocity and flame height. Other performance characteristics include the heat transfer rate from the flame to the burner surface, the conditions that lead to flickering and instability at the tip, and the interactions between flames on adjacent slots. An explanation of these rests upon a fundamental understanding of combustion, including the influence of the flow field on the flame and the effects of flame stretch. Slots usually are long and, consequently, mathematical models can be two-dimensional [1]. The application of such a steady-state model to lean and stoichiometric laminar burning of premixed methane-air is described, together with the associated experiments. A number of chemical kinetic schemes has been investigated and a compromise has been reached between kinetic accuracy and acceptable computing time. Predictions of flame heights are in good agreement with those experimentally measured, as are the predictions of gas temperatures when compared with those measured by the CARS
technique. Of particular interest are the conditions that lead to flame liftoff, and the role of flame stretch. Both the curvature stretch and the aerodynamic strain contributions to the flame stretch rate are computed over the entire leading edge of the flame surface. It would appear that aerodynamic straining at the base of the flame is the main cause of liftoff, while the tip of the flame, which experiences negative curvature stretch, is comparatively stable.
The Mathematical Model Governing Equations Horizontal and vertical coordinates are x and y, while the velocities are u and m, respectively. The global mass conservation equation is
915
](uq) ](mq) ` 40 ]x ]y
(1)
The density is q and the mass fraction of the species i is Yi, so that the species mass conservation equation is
916
LAMINAR PREMIXED FLAMES
in which Cp is the specific heat at constant pressure, k the thermal conductivity, T, the temperature, hi, the specific enthalpy, and Wi the molecular mass of species i. The solution domain and the boundary conditions for adiabatic flow are shown in Fig. 1. Numerical rounding errors in Yi were reduced by obtaining YN2 by difference N
YN2 4 1 1
o i 1 4
Yi
(6)
i?N2
The ideal gas law was assumed and, for the values of transport properties, the simplifications of the classical Chapman-Enskog procedure described by Dixon-Lewis [2] were adopted. Diffusion velocities were expressed by the three terms in the series proposed by Oran and Boris [3]. In order to discuss flame stability, it is necessary to derive an expression for the flame stretch rate. A number of expressions have been proposed for this [4] and all have their context of convenience. From the general expression derived by Matalon [5], it can be shown that the flame stretch rate is given by 1 dA ]u ]u ]m 4 cos2 a 1 sin a cos a ` A dt ]x ]y ]x
1
Fig. 1. Solution domain and boundary conditions.
` sin2 a ]((u ` ui)qYi) ]((m ` mi)qYi) ` 4w ˙i ]x ]y
(2)
in which ui and mi are diffusion velocities and w˙i is the chemical source term. The two momentum conservation equations are ](quu) ](qum) ]p ` 41 ]x ]y ]x `
] ]u ] ]u g ` g ]x ]x ]y ]y
1 2
1 2
(3)
](qum) ](qmm) ]p ] ]m ` 41 ` g ]x ]y ]y ]x ]x
1 2
`
] ]m g ` qg ]y ]y
1 2
(4)
in which p is the pressure, g the viscosity, and g the gravitational constant. The energy conservation equation is Cp
](uqT) ](mqT) ] ]T ] ` Cp 1 k 1 ]x ]y ]x ]x ]y
1 2
]T
]p
]p
]
1k ]y2 4 u ]x ` m ]y 1 ]x ]
N
N
N
1q o Y u h 2 1 ]y 1q o Y m h 2 1 o h W w˙ i41
i i i
i41
i i i
i41
i
i
i
(5)
2
]a ]m ]a ` un sin a 1 cos a ]y ]y ]x
1
2
(7)
in which A is the area of a flame surface element, t is the time, un the laminar burning velocity influenced by flame stretch, and a is the angle shown in Fig. 1. The first three terms in Equation 7 comprise the aerodynamic straining, and the last two the curvature stretch. The distinction between the two contributions to the total flame stretch rate is important in that they are each associated with a different Markstein length [4]. Chemical Kinetics and Numerical Solution Four chemical kinetic schemes were investigated. The first was the “full” reaction mechanism of Dixon-Lewis [2], comprising 53 elementary reactions with 18 chemical species. This proved to be expensive in computer time, and a “reduced” mechanism was sought. This led to a second scheme proposed by Peters [6] that retains the 18 major reactions in the “short” scheme of Miller et al. [7]. In this scheme, more detailed computations [8] had led Peters to achieve simplification through steady-state approximations for OH, O, HO2, CH3, CH2O, and CHO, and various partial equilibria. It was found that for lean mixtures up to stoichiometric, there was good agreement between profiles of species and temperatures computed by this and by the DixonLewis scheme. As a result, the Peters scheme was adopted for most of the computations because of its
TWO-DIMENSIONAL SLOT BURNER LAMINAR FLAMES
combination of sufficient accuracy and shorter computing time. A third scheme of Peters and Williams [9] employed all 31 elementary reactions of Ref. 7, with some assumptions of partial equilibrium. This scheme, however, yielded burning velocities that were excessively high, a condition that, no doubt, could have been remedied with different numerical values for some of the rate constants. The final scheme of Mauss and Peters [10] employed 40 reactions in the C1 mechanism. In this, there are no partial equilibrium, but some steady-state, assumptions. Numerical convergence with this final scheme was very slow, the computing times prohibitive, and it was abandoned. The present study does not involve rich mixtures, and C1 mechanisms are sufficiently adequate. Numerical solutions were obtained with control volume discretization and central differencing for diffusive fluxes. For the convective terms, the boundednesspreserving CCCT scheme was employed [11]. The standard, line-by-line iterative method took the form of the penta-diagonal matrix algorithm (PDMA) for matrices that included five coefficients. This is necessary for the higher-order CCCT scheme, but a tridiagonal matrix can be used for lower-order schemes. The SIMPLE algorithm [12] was used to solve for the velocity and pressure fields. Grid dependency of the solutions was checked and it was found that 80 grid points were adequate in the vertical direction for a flame height of about 5.5 mm, with 40 points in the horizontal direction for a slot width of 3 mm. For stability, it was necessary to extend the computational domain in the y direction to a height of at least four times the flame height, although the grid spacing could be increased outside the reaction domain. For the computation of flame stretch rates, a finer grid was necessary and a grid size 140 2 140 was employed. Full computational details are given in Ref. 13. Solutions were obtained for an initial pressure of 1 atm, with inlet gas and burner surface temperatures of 290 K.
Experimental Measurements The burner had five slots, each 50 mm in length and 3 mm in width. The measurement volume was midway along the length of the central slot. The two adjacent slots created symmetrical conditions and the outer slots prevented air entrainment. The center planes of the slots were pitched 3 mm apart. The burner was made of five 6-mm-thick brass plates, with a 3-mm-thick end plate. Air entrainment was eliminated by shields of flowing nitrogen. The premixture entered at the base of the burner and flowed through 1-mm-diameter holes drilled in a 12-mmthick ceramic block, to distribute the flow uniformly. The flow developed further aerodynamically as it
917
flowed along a square-sectioned steel tube 195 mm in length into the head of the slotted burner. The integrated CARS laser system was supplied by Spectron Laser System, Ltd. This comprised a Qswitched, frequency-doubled Nd:YAG laser with an energy of 350 mJ, and a 35-mJ broadband dye laser. Data collection was via a CCD camera. The laser system and receiving optics remained stationary while the burner was traversed perpendicular to the beams. A Dell System 433 PC was used in conjunction with the Harwell CARP-PC and QUICK-PC software to generate and match nitrogen spectra. A single measurement point in the flame yielded 100– 150 individual spectra, and the standard deviation of the measurements was about 80 K. The flame geometry was ideal for CARS measurements in that the length of the measuring volume, which could be approximately 3 mm in the direction of the laser beams [14], could be aligned parallel to the slot and along an isotherm. Results and Discussion Temperatures and Flame Heights There was good agreement between measured and predicted flame heights and temperatures. The latter were well predicted through the reaction zone. Some temperatures, for an equivalence ratio f of 0.84, are presented in Fig. 2 for two different cold mean flow velocities, vm, in the slot. Temperatures are plotted against vertical height above the burner surface for three different horizontal distances from the axis of the slot. Measured values are shown by the symbols. Those of burned gas at the burner tip tended to be lower than those predicted, suggesting some underestimation of the heat loss to the top of the burner, which was computed as 5% of the total energy at the higher flow rates and 30% at the lowest flow rates. This was confirmed by thermocouple temperature measurements at the top surface. These were 334 K and 343 K at the high and low flow rates, respectively, above the assumed boundary condition of 290 K. Heat loss cannot completely explain the lower measured temperature at the flame tip in Fig. 2a. Another relevant factor is the difficulty of aligning the measuring volume along an isotherm in this region of high gradients. Any presence of nonrepresentative colder gas, because of the nonlinearity of the CARS technique, results in a disproportionately lower measured temperature. The conductive heat flux to the horizontal burner surface is given by (k]T/]y)y40, and computed values of it across the surface for the two flow rates of the stoichiometric mixture are shown in Figs. 3b and 4b. Shown in Figs. 3a and 4a, along with the velocity vectors, are the associated contours of volumetric heat release rate and, beneath these, two isotherms.
918
LAMINAR PREMIXED FLAMES
Fig. 2. Comparison of measured (points) and computed (line) temperature profiles, with z 4 0.84 and (a) mm 4 2.0 m s11; (b) mm 4 0.6 m s11.
Note the differences in height scales and that the lower velocity flame is located closer to the burner surface, with an increase in the heat flux to it. Flame height is defined as the axial distance above the top of the burner surface to the point of maximum heat release rate. Over the range of the computations, it was found to vary linearly with the mean cold flow velocity. For a 2-mm slot width, the gradients of flame height with respect to the mean velocity are 5.0 2 1013, 4.6 2 1013, and 3.6 2 1013 s11, for values of f of 0.75, 0.84, and 1.0, respectively. For the 3-mm slot width, these values are increased by 50%. The lifted height is defined as the minimum distance between the top burner surface and the base of the flame where the volumetric heat release rate is 2 GW m13. Like flame height, it increases linearly with cold flow velocity, but at a lower rate than the former. The computations show flashback to occur at mean velocities of 0.30, 0.34, and 0.40 m s11, for values of f of 0.75, 0.84, and 1.0, respectively. Attempts to predict liftoff were less successful due to numerical instabilities (the origin of which is probably physical) at the associated higher flow rates. The maximum mean cold flow burner velocities to give converged solutions are 2.0, 2.1, and 2.75 m s11, again in ascending order of f. The corresponding
Fig. 3. Methane-air, f 4 1.0, mm 4 0.4 m s11. (a) Twodimensional fields of heat release (GW m13), cold isotherm (300 K), and velocity; (b) corresponding stretch rates at the 300-K isotherm, and heat flux to burner surface.
experimental liftoff velocities are 3.2, 4.7, and 5.6 m s11. Flame Structure At the lower flow rates, photographs of the flame on the slot burner reveal a continuous flame sheet, with fairly equal luminosity at the tips and bases of the flame. Similarly, the computed contours of volumetric heat release rate are continuous, although there is some reduction in this at the bases. The velocity vectors in Fig. 3a show diverging flow toward the tip, and for these conditions the flame stretch, measured on the 300-K isotherm, predominantly arises from aerodynamic strain at the tip, where it is positive. Figure 3b shows it to decrease toward the base where it becomes negative, a consequence of the low-velocity, converging flow above
TWO-DIMENSIONAL SLOT BURNER LAMINAR FLAMES
919
the angle a. At the base, the two contributions to flame stretch are more comparable, but, as will be shown, the aerodynamic strain term is more dominant and negative. The computations show there to be little diffusive migration of CH4 away from the center of the inverted “V” of unburned gas to the reaction zones on the flanks. Any diffusive stratification of the mixture arises from diffusion of H2 from the flanks, resulting in a significant concentration of H2 ahead of the reaction zone at the tip. These aspects are illustrated in Fig. 5 where mole fraction contours are presented for both CH4 and H2, for two different sets of conditions. Any H2 enhancement at the base is far less significant than it is at the tip. Markstein Lengths The effect of flame stretch rate on the burning velocity, un, is expressed by [4] ul 1 un 4 Lsas ` Lcacb
Fig. 4. Methane-air, f 4 1.0, mm 4 2.75 m s11. (a) Twodimensional fields of heat release (GW m13), cold isotherm (300 K), and velocity; (b) corresponding stretch rates at the 300-K isotherm, and heat flux to burner surface.
the burner surface between slots. The burning velocities decrease in harmony with the strain rate. As the flow rate increases, both the intensity of the luminosity and the computed heat release rate diminish at the base, much more than at the tip. At the highest flow rates, the flame stretch rate is predominantly due to flame curvature at the tip where it is negative and large, and to aerodynamic strain at the base where it remains negative. At the tip, the burning velocity increases strikingly as the curvature stretch becomes more negative, whereas the burning velocity decreases toward the base as the aerodynamic strain rate also decreases. Figure 4 typifies these conditions, although it is difficult to determine the stretch rate terms accurately in equation 7 due to the sensitivity of the terms involving gradients to
(8)
in which ul is the stretch-free laminar burning velocity, L a Markstein length, a an associated component of the flame stretch rate; subscript s indicates aerodynamic strain and subscript c indicates curvature stretch. Values of the two Markstein lengths have been evaluated separately for different equivalence ratios by the application of equation 8 to a region of the flame where one of the terms on the right, in turn, is dominant over the other, both at the tip and along the flanks. The equation was applied to the front of the flame along the 300-K isotherm. Values of un with a dominant stretch rate component were plotted against that stretch and extrapolated to zero stretch to give the value of ul. The slope of the straight line plot gives the value of the Markstein length, and mean values are given in Table 1. Also tabulated are the rather more exactingly computed and more accurate values for symmetric, one-dimensional, spherical flames, taken from Ref. 4. Although the agreement between the values is not good, Markstein lengths are notoriously difficult to quantify precisely. What is more important for present purposes is that Ls is negative and Lc is positive. Another aspect is that attribution of numerical values of Markstein lengths to a given mixture suggests they are physicochemical parameters. However, diffusive stratification which is quite different for the slot burner from spherical flames and which clearly depends upon the flame geometry, can alter the mixture composition ahead of the reaction zone. Flame Stability and Liftoff Values of un peak at the flame tip, slowly decrease down the flank, then decrease more sharply at the base of the flame. For example, with f 4 1.0, mm 4 1.0 m s11, and for a 2-mm slot, the value of un is
920
LAMINAR PREMIXED FLAMES
Fig. 5. Mole fractions of (a) CH4 and (b) H2 for f 4 1.0, mm 4 1.8 m 11, and f 4 0.84, mm 4 1.4 m s11. Slot width 4 3 mm. TABLE 1 Computed values of Markstein lengths (mm)
f Ls (present) (Ref. 4) Lc (present) (Ref. 4)
0.75
0.84
1.0
10.82 10.28 0.03 0.10
10.65 10.23 0.01 0.10
10.25 10.19 0.01 0.10
1.25 m s11 close to the tip, 0.35 m s11 at the center of the flanks, and 0.1 m s11 close to the base. These changes are entirely commensurate with a positive value of Lc at the tip, which with the negative curvature stretch increases un there, and a large negative value of Ls at the base which with the slightly negative aerodynamic strain decreases un. Two other influences that reinforce these trends are the diffusion of H2 into unreacted gases at the tip and the conductive loss to the burner surface at the base, which reduces the burning velocity there more than
would be expected from the values of Ls and the aerodynamic strain. As the flow rate increases toward liftoff, the dominant negative curvature stretch at the tip, because of the positive value of Lc, leads to a further increase in burning velocity, with un $ ul. In contrast, at the base of the flame, as the flow rate increases so does the negative aerodynamic strain, and the burning velocity is further reduced. The computations also show the volumetric heat release rate to decrease significantly at the base. The flattened base of the flame disappears, a sideways and outward flow develops between the burner and the base of the flame, that bypasses the flame. Ultimately, the burning velocity will become insufficient for the flame to propagate against the gas flow. Thus, liftoff would originate at the base of the flame, notwithstanding the stable characteristics of the tip. This is in line with the observations of Law et al. [15] who were unable to observe tip opening with lean CH4 flames because the flame lifted off before such an opening could occur. However, this behavior of methane-air flames cannot be
TWO-DIMENSIONAL SLOT BURNER LAMINAR FLAMES
generalized to other fuels, particularly those of high molecular mass and other equivalence ratios. A positive value of Ls would alter the liftoff characteristics.
Conclusions 1. A two-dimensional laminar combustion model with a short reaction scheme from Peters has been implemented successfully for lean methaneair mixtures. 2. CARS temperature measurements on a slot burner are in fairly good agreement with model predictions. 3. The heat loss to the burner surface increases as the gas flow decreases and ranges from 5 to 30% of the energy in the fuel. 4. Flame height and lift increase linearly with flow rate. 5. Burning velocities and flame stretch rates have been evaluated at the front of the flame along the 300-K isotherm. Aerodynamic strain and flame curvature Markstein lengths have been obtained from these data. 6. Increasing negative stretch curvature, associated with an increasing flow rate, increases un locally and stabilizes the tip of the flame. Simultaneously, increasingly negative aerodynamic strain at the base of the flame decreases un locally and eventually leads to liftoff. 7. These effects are reinforced by the diffusion of H2 toward the tip and heat loss at the base. Acknowledgments The authors thank British Gas plc, particularly Drs. M. Davies, S. Hasko, S. C. Taylor, and D. B. Smith, and EPSRC for their support.
921
REFERENCES 1. Somers, L. M. T. and De Goey, L. P. H., Combust. Sci. Technol. 108:121–132 (1995). 2. Dixon-Lewis, G., in Combustion Chemistry, (W. C. Gardiner, Ed.), Springer-Verlag, New York, 1984. 3. Oran, E. S. and Boris, J. P., Prog. Energy Combust. Sci. 7:1 (1981). 4. Bradley, D., Gaskell, P. H., and Gu, X. J., Combust. Flame 104:176–198 (1996). 5. Matalon, M., Combust. Sci. Technol. 31:169 (1983). 6. Peters, N., in Numerical Simulation of Combustion Phenomena, Lecture Notes in Physics, 241 (R. Glowinski, B. Larrouturou, and R. Temam, Eds.), SpringerVerlag, Berlin, 1985, pp. 90–109. 7. Miller, J. A., Kee, R. J., Smooke, M. D., and Grcar, J. F., Western States Section of the Combustion Institute, Spring Meeting, 1984, Paper WSS/CI 84-10. 8. Warnatz, J., Eighteenth Symposium (International) on Combustion, The Combustion Institute, Pittsburgh, 1981, pp. 369–384. 9. Peters, N. and Williams, F. A., Combust. Flame 68:185–207 (1987). 10. Mauss, F. and Peters, N., in Reduced Kinetic Mechanisms for Applications in Combustion Systems, Lecture Notes in Physics, m 15, (N. Peters and B. Rogg, Eds.), Springer-Verlag, Berlin, 1993, pp. 58–75. 11. Gaskell, P. H. and Lau, A. K. C., Int. J. Num. Methods Fluids 8:617 (1988). 12. Patankar, S. V. and Spalding, D. B., Int. J. Heat Mass Transfer 15:1787 (1972). 13. Kwan, K. C., “Mathematical Modelling of Premixed, Laminar Methane-Air Flames,” Ph.D. Thesis, University of Leeds, 1994. 14. Bradley, D., Lawes, M., Scott, M. J., Sheppard, C. G. W., Greenhalgh D., and Porter, F. M., TwentyFourth Symposium (International) on Combustion, The Combustion Institute, Pittsburgh, 1992, pp. 527– 535. 15. Law, C. K., Ishizuka, S., and Cho, P., Combust. Sci. Technol. 28:89–96 (1982).