Numerical investigations of unsteady flame propagation in stepped microtubes

Abstract

Transient numerical simulations with detailed chemistry have been carried out for premixed stoichiometric CH₄–air and H₂–air flames in two-dimensional stepped microtubes for a range of wall heat transfer conditions. Investigations on such configurations are important from the perspective of the design of micro combustion devices, flame arresters and safety in domestic and industrial combustion devices. Similarities in flame propagation characteristics have been brought out through a detailed analysis for both the fuels. Detailed analysis of the propagating flame near the channel step revealed an interesting phenomenon of sudden increase in flame propagation velocities for a certain range of wall heat transfer coefficient, h. A quantitative value of ratio of heat-loss to heat-generation at the contraction has been proposed which helps predict the flame propagation through sudden channel steps.

---

1. Introduction

Recent trends of miniaturization of various small scale devices such as drones, smartphones, micro UAVs, laptops and computers etc. have motivated researchers to develop high energy density alternatives to traditional low energy density electrochemical Li-ion batteries.^1–6 Micro-combustion based systems are one of the promising options due to their inherent advantages of higher heat and mass transfer coefficients, energy densities, and shorter recharge times compared to existing batteries.^1–6 Furthermore, the energy density of hydrocarbon fuels is 20–50 times greater than that of advanced Li-ion batteries.^2–6 Micro-combustors also ensure lesser chemical emissions due to their lower operating temperature range.

A series of studies aimed at understanding flame propagation in microchannels are briefly summarized below. Numerical studies of Kim and Maruta⁷ report the effect of flow velocity profile, wall heat transfer and tube diameters on laminar flame propagation velocities, flame shapes and quenching characteristics. Lee et al.⁸ through their experimental studies have revealed the flame stability limits for methane–air and propane–air mixtures over a wide range of tube diameters and mixture temperatures. Maruta et al.⁹ have reported various stable and unstable flame propagation modes in externally preheated tubes. Kumar et al.³ and Fan et al.⁴ have observed that for a range of operating conditions, rotating flame patterns are formed in radial micro channels. Effect of inlet velocity, mixture equivalent ratio on wall heat losses, flame and wall temperatures were investigated by various researchers.^10,11 Jiang et al.¹² have proposed a new entropy generation rate perspective to analyze the effects of flow velocity, equivalence ratio and heat recuperation in their detailed chemistry premixed flame simulations for H₂–air mixtures. Flame stabilization experiments in diverging micro-channels by Kumar and coworkers^13–15 report the existence of various stable and partially stable flame propagation modes for a range of mixture equivalence ratios and inlet velocities. Various studies on flame propagation in closed tubes of various diameters have been reported in literature.^16–20 The formation of various flame propagation modes, flame evolution and acceleration has been attributed to Darrieus–Landau instability (D–L instability), Taylor instability, vortex motion and hydrodynamic effects of near field flow.^16–18,20

Suddenly diverging channels or stepped micro-combustors used in this study serve as an important configuration for design purpose of flame arresters and, safety in domestic and industrial combustion devices. Various investigations on flame stabilization in these multi-step micro combustors with premixed fuel–air mixtures have been reported in the recent literature.^21–24 They have reported the effect of length of various steps, diameter and material type on flame stability limits. A wire insertion method on a stepped microtube combustor was introduced and numerically investigated by Baigmohammadi et al.²⁵ to improve flame stabilization and combustion process and temperature distribution by activating certain chemical reactions which significantly affect the combustion initiation reaction. Numerical investigations of Gutkowski²⁶ report that the shape of the propagating flame just before a sudden contraction of a stepped tube is expected to govern quenching of the flame at the step. His investigations for a range of tube diameter ratios revealed that dimensions of larger tube determine the flame quenching when the smaller tube diameter is slightly greater than the quenching distance of a given fuel–air mixture. Liu et al.²⁷ introduced annular stepwise diverging tube (ASDT) method to supersede annular diverging tubes (ADT), to measure laminar-burning velocities of premixed flames. The ASDT helped enhance flow uniformity and to formulate flatter flames in majority of the experimental conditions. Liu and Kim²⁸ later introduced the assembled annular stepwise diverging tube (A-ASDT), which provided additional benefits with ease of burner adjustments and control of experimental resolutions.

Although experimental studies are important criteria towards validations of results, however, not all aspects of micro-combustion can be measured and studied through experiments due to difficulties associated with smaller dimensions of systems. Numerical studies allow understanding of the flame propagation characteristics through a more detailed analysis and from a different perspective. Majority of the numerical studies reported in the literature on flame propagation in various channels are restricted to simple single-step reaction models. Furthermore, majority of these studies were conducted with the limiting cases of adiabatic and isothermal walls with little attention to flame propagation through stepped channels. Although these conditions provide insights on general aspects of flame propagation, it is not easy to extrapolate these results to accurately model the realistic and practical thermal boundary conditions existing between an adiabatic and isothermal wall.

The aim of present work is to understand flame propagation in stepped tube configurations with various thermal boundary conditions and to understand transient flame dynamics near the sudden contraction of the step; through the effects of hydrodynamic instability, flame–flow interactions, flame–wall interactions and heat-loss to heat-generation ratio leading to flame quenching. This study also aims at identifying and proposing fuel independent trends and species dependent trends governing flame propagation and quenching in such channels. A suddenly converging geometry is used in this study to investigate the following important aspects related to flame propagation and quenching. (a) A sudden contraction encountered by the flame is a chemically inert and cold surface, resulting in a critical loss of heat from the flame thus hampering the continuous flame propagation through the contraction. (b) The thermal boundary conditions where a flame propagates into the tube with its diameter smaller than the numerically calculated quenching distance of corresponding fuel–air mixture.

Stoichiometric H₂–air and CH₄–air are used for the present study because of their broad usage, and a large difference between their laminar burning velocity and quenching distance. Although methane is a popular fuel for commercial interests, hydrogen has recently gained popularity due to its CO₂ free and soot-free combustion. Enriching hydrocarbon fuels with H₂ gas is also considered beneficial for many practical applications. Furthermore, its entire detailed chemistry can be modeled with 21 elementary reactions, which makes numerical investigations with hydrogen computationally less expensive. H₂–air also has a small quenching distance, allowing the development of smallest micro-combustors.

The paper is organized in three major sections. Details of the computational domain and various models, such as domain details, reaction modeling, boundary conditions and solution approach are presented in Section 2. Results of flame propagation, extinction are discussed in Section 3. The outcomes of the present work are summarized in last section of this paper.

2. Computational details

2.1 Computational domain

The details of the computational domain used for the present studies are shown in Fig. 1. It consists of a smaller tube with diameter, d and a length of 30 mm joined with a larger tube of diameter D and length 30 mm. Entire domain is considered with an axis of symmetry about the centerline of the tube. The tube diameters are kept different for both the fuels as it is intended to keep the smaller tube diameter close to flame quenching distance (d_q) of respective fuels. The quenching distances of CH₄–air and H₂–air are determined computationally for isothermal case with T_w = 300 K. They were found to be 3.0 mm and 0.74 mm, respectively using the detailed reaction mechanisms for CH₄ and H₂ fuels.^6,23 For present numerical work, a smaller tube diameter, d is chosen as 0.7 mm for H₂–air studies, a value slightly smaller than the quenching distance of H₂–air mixtures. This is further scaled directly with diameter ratio −d_q(CH₄)/d_q(H₂) to obtain d for CH₄–air studies as 2.84 mm. D is chosen to fix contraction ratio (D/d) = 2 for both the fuels.

Fig. 1 Schematic of computational domain and grid (a) details of computational domain of suddenly converging channel (b) adapted domain near the flame front in CH₄/air studies (c) zoomed view of adapted region.

For H₂–air studies, a uniform grid of square cells with dimensions equal to 15 μm was chosen. An average thermal flame thickness of 0.3 mm for H₂–air case ensured with at least 20 points within the reaction zone to sufficiently resolve all the features of the propagating flame front. For CH₄–air flames, a dynamically adaptive grid system was chosen. Here, the initial grid consisting of uniform square cells with dimensions equal to 0.1 mm was chosen. The criterion for grid refinement and coarsening was normalized temperature gradient in the computational domain being greater or less than 0.3. The grid refining process was carried out in three levels with the regions of highest gradient being refined to the finest square cells of 12.5 μm. Since the flame (region of highest temperature gradient) is unsteady and propagating throughout the domain, the dynamic gradient adaptation, which adapts the grid after every 30 time steps (Δt = 10⁻⁵ s) ensures that region of interest is refined just before the flame enters it and again coarsened when the flame leaves it. All these details of the adapted grid and computational domain are shown in Fig. 1a–c.

2.2 Governing equations

The present study is carried out for stoichiometric CH₄–air and H₂–air mixtures with laminar flow in a suddenly contracting 2-D axisymmetric channel. Several assumptions were made to simplify the present analysis. Work due to pressure and viscous forces (viscous heating) has been neglected. Since, it was intended to study the effect of thermal boundary conditions on the flame shape, gravity effects were neglected. The species transport model assumes only volumetric reactions with no effect of wall reactions (chemically inert walls), since wall reactions play important role on flame extinction through radical quenching only at very low pressures²⁹ and high wall temperatures (>1000 K).³⁰ The overall effect of thermal radiation is also neglected as it results in a slight reduction in the maximum temperature in the domain³¹ along with a slight shift in the position of the stabilized flame to a downstream location.¹⁴ The effect of radiation is neglected as it is not expected to significantly affect the flame propagation in present studies. The governing equations are summarized below.



where v⃑ = (u₁, u₂) or (u, v) is the velocity vector with x and y components.

In the above equations ρ, p, μ, h, k, M_j and Y_j represent fluid density, pressure, viscosity, enthalpy, thermal conductivity, molecular weight of j^th species and mass fraction of j^th species, respectively. is the diffusive mass flux of j^th species due to species and temperature gradients. It is calculated from Maxwell–Stefan equations of multi-component diffusion using local mass diffusivity and thermal diffusion coefficients from kinetic theory relations. R_j is the volumetric reaction rate of j^th species, calculated using finite rate chemistry. The thermodynamic properties such as viscosity, specific heats, and thermal conductivity of individual species are assumed to vary with temperature and calculated using kinetic theory, NASA piecewise polynomials and kinetic theory relations, respectively. The fluid specific heat is calculated as mass-fraction weighted average of specific heats of all species and the fluid thermal conductivity, viscosity are calculated using ideal gas mixing laws. The flow is assumed to be incompressible as the flow velocities are very small (Ma < 0.01).

2.3 Chemical reaction modeling

The H₂–air studies were carried out using 9 species and 21-step detailed mechanism developed by Connaire et al.³² A smaller domain and relatively less number of reactions reduce the requirements of the computational resources. For CH₄–air mixtures, due to large number of species and hundreds of elementary reactions, a reduced form of the established GRI 3.0 (ref. 33) detailed mechanism has been used in this study. The reduced form of GRI 3.0 predicts burning velocity at 1 atm and 300 K with an error less than 4% compared to full GRI 3.0 mechanism. The detailed chemistry is efficiently modeled using 27 species, and 159 reactions. NO chemistry is completely ignored because the emphasis of the present studies is on understanding the flame propagation in small tubes with a step in the channel.

2.4 Boundary and initial conditions

A uniform temperature profile and axial velocity (u = 1 cm s⁻¹ in CH₄–air; u = 30 cm s⁻¹ in H₂–air) is assigned at the channel inlet, where the premixed fuel–air mixture (T_in = 300 K) enters in stoichiometric proportions. A small flow velocity at the mixture inlet helps in improved solution convergence. The flame dynamics remains unaffected due to this small mixture velocity at the inlet, as the respective flame propagation velocities are relatively much higher. Ambient pressure is assigned at the channel outlet with Neumann boundary conditions. Various flame propagation studies in micro-tubes by Pizza et al.³⁴ and Fan et al.³⁵ with H₂/air and CH₄/air mixtures report the formation of axisymmetric flames at low mixture velocities. Therefore, the axisymmetric condition was applied at the channel axis (center). Walls are considered to be chemically inert and negligible thickness.^12,29–31 No slip boundary condition is applied at the channel walls as the flow Knudsen number is very small and varies between 0.3 × 10⁻⁴ and 1.5 × 10⁻⁴. Thermal boundary conditions at the walls are varied between the limiting cases of adiabatic (h = 0) and isothermal (T_w = constant), along with different convective heat transfer coefficients, h to model the intermediate wall heat transfer conditions. Ignition is provided to the fuel–air mixture through a 0.3 mm thick, 2400 K temperature flat ignition patch at the exit plane of the channel to ensure uniformity in the obtained results.^7,8,36,37

2.5 Numerical solution approach

Simulations were carried out using a general purpose CFD code, Fluent 15.0,³⁸ which allows import of detailed chemistry of fuel–air mixtures in CHEMKIN format. Pressure–velocity coupling algorithm, PISO³⁹ was used to carry out the transient computations for both H₂–air and CH₄–air cases. The chemistry integration was performed using ISAT integration method.⁴⁰ This method speeds up the numerical simulations significantly by tabulating the accessed composition space region “on the fly” (in situ) with error control (adaptive tabulation) with sufficiently small error tolerance of 10⁻⁴. Detailed simulations were carried out to understand the effect of various spatial and temporal discretization schemes and the results are shown in Fig. 2. It is clear from Fig. 2 that the heat release rate and OH species profile are predicted well by various schemes. Therefore, second order accurate schemes are used for spatial and temporal discretization in the present work. The solution is considered converged when the scaled residuals of mass, momentum, species and energy equations drop to 10⁻⁵ and there was no appreciable change in the residuals further. A time step (Δt, s) of 10⁻⁵ and 10⁻⁶ is taken for CH₄–air and H₂–air mixtures, respectively. An average CPU time of 90 hours was required for each CH₄–air and 50 hours for H₂–air flames on a 16 core Xeon server.

Fig. 2 Comparison of (a) heat release rate (b) OH species mole fraction at a position 4d_q from the channel exit plane for various grids and spatial and temporal discretization schemes. A_T1S2 – first order in time and second order in space and 2 level adaptive grid, B_T2S2 – second order in time and second order in space and 3 level adaptive grid to 0.0125 mm level, C_T2S1 – second order in time and first order in space without adaptive grid (0.025 mm) D_T2S2 – second order in time and second order in space and without adaptive grid (0.0125 mm) (c) predicted laminar burning velocity for different adaptive and non-adaptive grids (d) flame structure of a laminar flame obtained from present computations.

2.6 Grid independence studies

Grid independence for H₂–air studies was carried out by comparing the species profiles and heat release rates (HRR) of the chosen grid with the finest square grids of size 15 μm.⁴¹ Detailed grid independence studies for H₂–air mixtures in similar meso-scale channels have been reported by authors earlier.⁴¹ Grid independence studies for CH₄–air were carried out by comparing the heat release rates, OH species profiles and laminar-burning velocities (S^o_L) for various grids and the results are shown in Fig. 2. It can be noted from Fig. 2a and b that the heat release rate and OH species profile predicted with grids of size 25 μm and 12.5 μm completely overlap with each other. There are sufficient grid points within the propagating flame indicating that the reaction zone is sufficiently resolved. The laminar burning velocity is obtained by noting the flame propagation velocity (S_L) for adiabatic walls with free slip velocity condition in sufficiently small tubes. The flame physics with these boundary conditions corresponds to a freely propagating adiabatic flame. The location of the mean temperature contour in the domain (average of maximum and minimum temperature in the domain) is defined as the flame position and its rate of change as flame propagation speed, S_L. Other criteria based on maximum CH₃ concentration and maximum HRR, also yield same results. The flame propagation speed, S^o_L obtained using the finest grid and the adaptive grid chosen for present work is equal to 36 cm s⁻¹, which is equal to laminar burning velocity (S^o_L) of stoichiometric CH₄–air mixtures reported in recent studies. For instance, 35.3 cm s⁻¹ has been reported as the laminar burning velocity for stoichiometric CH₄–air mixtures in the studies of Akram and Kumar,¹⁵ Akram et al.⁴² and Liu and Kim.²⁷ Flame structure depicting the variation of major species for adiabatic methane flames in sufficiently small tubes is shown in Fig. 2d. The flame thickness obtained from present studies compares well with that reported by Kim and Maruta⁷ using single-step mechanism of CH₄–air mixtures.

3. Flame propagation characteristics near contraction

3.1 Variations in flame shape angle

Fig. 3 shows the normalized heat release rate contours of the flame shape at various axial locations as flame propagates through the tube for a wall heat transfer condition of h = 37 W m⁻² K⁻¹. It is clear from this figure that the flame shape angle remains constant beyond 3d_q distance from the exit plane as flame shape is fully developed. Therefore, the flame shape in terms of flame shape angle has been used in the present work to discuss various characteristics related to flame propagation in the channels. All the results presented in this paper are obtained upstream of 3d_q location from the exit plane to ensure that it represents the behavior of a fully developed propagating flame. Fig. 4 shows the normalized heat release rate contours for CH₄–air and H₂–air flames with different wall heat transfer conditions. There are six contour levels for each of the CH₄–air flame front shown in Fig. 4a–e, each representing a multiple of 0.17 of maximum heat release rate. The innermost contour corresponds to a maximum normalized value of 1.0. Five contour levels are plotted for H₂–air flames as shown in Fig. 4(i)–(v), each representing a multiple of 0.2 times of the maximum heat release rate. A notable characteristic is the flame shape angle, θ, defined with respect to the flame shape as shown in Fig. 4 for various positively and negatively stretched flames observed for a range of wall heat transfer conditions. This characteristic will help in quantifying the change in flame shape much better as compared to other parameters. The angle made by the flame tip along the maximum heat release rate contour is defined as flame shape angle. A positive or negative flame shape angle is defined with respect to the vertical line as shown in Fig. 4b and d. As the wall heat transfer coefficient is varied for both CH₄–air and H₂–air mixtures, the flame shape changes from a positive flame shape to a negative flame shape. The transition in flame shape occurs for CH₄–air flame at h ∼ 43 W m⁻² K⁻¹ and for H₂–air flame, the transition has been observed to occur in h = 450–500 W m⁻² K⁻¹ range. This flame transition is expected to significantly alter the flame propagation characteristics as discussed in the later sections.

Fig. 3 Normalized heat release rate contours depicting the flame shape at different axial locations for h = 37 W m⁻² K⁻¹ condition. HRR contours are normalized to 4.05 × 10⁹ W m⁻³.

Fig. 4 Normalized heat release rate contours depicting flame shape angles at different thermal boundary conditions (a)–(e) CH₄–air flames (a) h = 0 W m⁻² K⁻¹, (b) h = 25 W m⁻² K⁻¹, (c) h = 42 W m⁻² K⁻¹, flames (d) h = 43 W m⁻² K⁻¹, (e) h = 50 W m⁻² K⁻¹; (i)–(v) H₂–air flames (i) h = 0 W m⁻² K⁻¹, (ii) h = 100 W m⁻² K⁻¹, (iii) h = 450 W m⁻² K⁻¹, (iv) h = 500 W m⁻² K⁻¹ (v) h = 800 W m⁻² K⁻¹ (contours normalized to maximum heat release rates – (a) 4.66 × 10⁹ W m⁻³ (b) 4.33 × 10⁹ W m⁻³ (c) 4.05 × 10⁹ W m⁻³ (d) 3.9 × 10⁹ W m⁻³ (e) 4.01 × 10⁹ W m⁻³ (i) 2.25 × 10¹⁰ W m⁻³ (ii) 2.23 × 10¹⁰ W m⁻³ (iii) 1.86 × 10¹⁰ W m⁻³ (iv) 1.79 × 10¹⁰ W m⁻³ (v) 1.893 × 10¹⁰ W m⁻³; contours plotted for x = −4d_q in CH₄/flames and x = 18 mm in H₂/air flames).

Many researchers have argued that flame stretch and hence the flame shape angle can be closely correlated to Darrieus–Landau instability.^7,16–18 The shape of the stabilized flame varies with tube diameter even when same boundary conditions such as velocity profile and wall heat transfer conditions are maintained. These transitions in flame shape are difficult to observe during experimental studies due to thermal stabilization effects.^34,43 However, these transitions have been observed to appear in the numerical studies of Kim and Maruta⁷ and these instabilities lead to flame stretching, even for adiabatic walls with free stream boundary conditions.

Fig. 5 and 6, show the variation of the flame shape angle across the channel step in CH₄–air flames for two different thermal wall boundary conditions of h = 0, 25 W m⁻² K⁻¹. The normalized heat release rate contours are plotted in the figures with their maximum values summarized below the figure. The flame angle decreases as the flame enters the step of the channel. The change in the flame shape angle is significant for h = 25 W m⁻² K⁻¹ case as shown in Fig. 6. A similar trend can also be observed for H₂–air flames as shown in Fig. 7 for adiabatic wall conditions. The flame shape angle varies from a value of 36° in the larger tube to ∼21° in the smaller tube. As the flame enters the smaller tube, the maximum temperature in the domain decreases because of increased heat loss due to increased surface area-to-volume ratio. This decrease in peak flame temperature reduces the density jump across the flame and hence the effect of D–L instability decreases. The hydrodynamic instability being an integral part of causing a change in the flame shape reduces, resulting in a further decrease in the flame shape angle for various fuel–air mixtures.

Fig. 5 Flame stretch angle variation before and after contraction for CH₄/air flames in adiabatic tube (a) after entering the contraction (b) before entering the contraction (contours normalized to maximum heat release rates: (a) 4.82 × 10⁹ W m⁻³ (b) 4.66 × 10⁹ W m⁻³).

Fig. 6 Flame stretch angle variation before and after contraction for CH₄/air flames in case of h = 25 W m⁻² K⁻¹ (a) after entering the contraction (b) before entering the contraction (contours normalized to maximum heat release rates: (a) 3.7 × 10⁹ W m⁻³ (b) 4.33 × 10⁹ W m⁻³.

Fig. 7 Flame stretch angle variation before and after contraction for H₂/air flames in case of adiabatic walls (a) after entering the contraction (b) before entering the contraction (d_q – numerically calculated quenching distance) (contours normalized to maximum heat release rates: (a) 2.12 × 10¹⁰ W m⁻³ (b) 2.25 × 10¹⁰ W m⁻³).

3.2 Variations in flame propagation velocity with wall conditions

Fig. 8 shows the variation of flame propagation velocity for different wall heat transfer conditions. It is interesting to note that for the cases of h = 0, 25 W m⁻² K⁻¹, the flame propagates steadily and suddenly accelerates as it is about to reach the channel step (<2d_q). As the flame enters the smaller tube, it settles to a flame propagation velocity smaller than initial flame propagation velocity. It is further interesting to note that for the case of h = 42 and 43 W m⁻² K⁻¹ case, the flame acceleration is completely absent and flame simply propagates into the smaller tube with a resulting flame propagation velocity smaller than initial flame propagation velocity. For h = 50 W m⁻² K⁻¹, the flame tries to propagate into the smaller tube, however leads to a complete flame quenching.

Fig. 8 Variation of S_L along the axis at different thermal boundary conditions for CH₄/air flames (a) varying h conditions (b) varying isothermal wall temperature conditions.

Fig. 9a–e, shows the HRR profiles upstream of contraction point for various wall heat transfer conditions. It is clear from this figure that the propagating flame undergoes a transition from a positive flame shape angle to a negative shape angle as the wall heat transfer coefficient increases. The extent of acceleration of the propagating flame shown in Fig. 8 depends on the flame shape just downstream of channel step. It is clear from Fig. 9c that a slight transition or change in flame shape towards negative occurs at h = 42 W m⁻² K⁻¹ condition. This transition has negligible effect on overall flame propagation speed. For the cases of h = 43–50 W m⁻² K⁻¹, the transition has significant variation in flame shape angle and flame propagation speed. The sudden acceleration in the flame speed as it approaches the channel step disappears for conditions of h ≥ 42 W m⁻² K⁻¹. Due to this, flame appears to be quenching as it enters the channel step. Similar transition in the shape of H₂–air flames has been observed to occur for h ≥ 450 W m⁻² K⁻¹.

Fig. 9 HRR contours depicting variation of flame shapes at different wall heat transfer conditions of CH₄–air flames (a)–(e) flame position at 4d_q downstream of contraction (i)–(v) flame position at 1d_q downstream of contraction (contours normalized to maximum heat release rates – (a) 4.67 × 10⁹ W m⁻³ (b) 4.33 × 10⁹ W m⁻³ (c) 4.05 × 10⁹ W m⁻³ (d) 3.86 × 10⁹ W m⁻³ (e) 4.01 × 10⁹ W m⁻³ (i) 4.37 × 10⁹ W m⁻³ (ii) 4.35 × 10⁹ W m⁻³ (iii) 4.05 × 10⁹ W m⁻³ (iv) 4.05 × 10⁹ W m⁻³ (v) 4.05 × 10⁹ W m⁻³).

Fig. 8b shows the variation of the flame propagation speeds for the case of isothermal wall conditions with T_w = 300, 350 and 450 K. It is interesting to note that for isothermal walls, the acceleration in flame speed is relatively gradual near the channel step. The step of the channel acts as additional wall area, preheating the mixture (T_mix = 300 K) in the region near to the contraction. The flame acceleration for T_w = 350 K and 450 K case is due to the fact that as the flame approaches the channel step, the mixture is preheated in the vicinity of contraction due to larger wall area. Due to higher mixture temperature, T_u, the preheated mixture results in higher flame propagation speed S_L.

Fig. 10 and 11 show the HRR profiles and flow streamlines for two cases of flame propagating in a larger channel and flame entering the channel step. The flow streamlines enter the stretched flame at larger angles (β = 30°) near the channel step as compared to the case, when the flame is away from the channel step (β = 60°) and the flow lines are relatively horizontal. In other words, the flux of fuel–air mixture entering the stretched flame increases suddenly as the flame reaches the channel step. The flame also gets stretched further (indicated by θ) to accommodate the increased flux of fuel–air mixture. Hence, θ increases significantly near the channel step. The combined effect of higher flame stretch (higher flame shape angle) and higher mixture mass-flux results in higher flame acceleration near the contraction. This effect for CH₄–air flames can be easily visualized in Fig. 10 and 11 for adiabatic (h = 0) and h = 25 W m⁻² K⁻¹ wall heat transfer conditions. The extent of flame acceleration decreases for higher h values. For instance, a decrease in θ values can be observed along with a smaller change in angle β as flame nears the channel step. Hence, a reduced supply of fuel–air mixture entering normally to the stretched flame is responsible for a reduction in flame acceleration due to a decreases in θ as shown in Fig. 10 and 11.

Fig. 10 HRR contours and streamlines for CH₄/air flames with adiabatic walls, depicting the change in axial mass flux and stretch angle as flame near contraction (a) larger stretch with increased axial mass flux near contraction (b) smaller stretch and lesser axial mass flux at 4d_q downstream and away from contraction (HRR contours normalized to maximum value of (a) 4.37 × 10⁹ W m⁻³ (b) 4.66 × 10⁹ W m⁻³).

Fig. 11 HRR contours and streamlines at h = 25 W m⁻² K⁻¹ for CH₄/air flames, depicting the change in axial mass flux and stretch angle as flame near contraction (a) larger stretch more mass flux near contraction (b) smaller stretch and lesser mass flux at 4d_q downstream and away from contraction (HRR contours normalized to maximum value of (a) 4.38 × 10⁹ W m⁻³ (b) 4.33 × 10⁹ W m⁻³).

Fig. 12 and 13 shows the HRR contours for the wall heat transfer conditions of h = 50 W m⁻² K⁻¹ and T_w = 300 K. It is clear from these figures that the flame at −2d_q location has a negative flame shape angle. As the flame nears the contraction, the flame continues to retain the shape. The flow streamlines show that the axial mass flux of the fuel–air mixture entering the flame is relatively small as compared to the case shown in Fig. 10 and 11. Due to this, the flame does not accelerate at it approaches the channel step. This leads to increased heat-loss from the propagating flame to the channel walls and results in flame quenching for both the conditions.

Fig. 12 HRR contours and stream lines at h = 50 W m⁻² K⁻¹ for CH₄/air flames, depicting the change in axial flux and stretch angle as flame near contraction (a) smaller stretch and axial mass flux near contraction (b) larger stretch and lesser axial mass flux at −2d_q contraction (HRR contours normalized to maximum value of (a) 4.33 × 10⁹ (b) 4.01 × 10⁹ W m⁻³).

Fig. 13 HRR contours and stream lines at T_w = 300 K for CH₄/air flames, depicting the change in axial mass flux and stretch angle as flame near contraction (a) smaller stretch and axial mass flux near contraction (b) larger stretch and lesser axial mass flux at −2d_q downstream and away from contraction (HRR contours normalized to maximum value of (a) 4.33 × 10⁹ W m⁻³ (b) 4.33 × 10⁹ W m⁻³).

3.3 Variations in heat-loss to heat-generation ratio

It was earlier proposed by Gutkowski²⁶ that quenching of the flame is governed by the flame shape near the contraction, with the tulip shape flames being more exposed to losses and hence leading to flame extinction in the smaller channels. The heat-loss from the flame through the walls, q_l depends on thermal boundary condition and the flame shape near the contraction. The interaction of the flame with chemically inert cold walls of channel step also plays an important role in flame propagation. The heat generated within the flame, q_s due to combustion is a property of fuel–air mixture, volume of the reaction zone, mixture temperature and mixture equivalence ratio. Therefore, the variation of heat-loss to heat-generation ratio, q_l/q_s near the contraction is expected to play an important role in governing the flame quenching. It is further discussed in this study that a particular upper limit value of q_l/q_s near the contraction may help predict the possibility of continuation of flame propagation in the smaller channel.

Fig. 14 shows the variation of the heat-loss to heat-generation ratio for both CH₄–air and H₂–air flames at different wall heat transfer conditions. It is interesting to note from Fig. 14a that for all conditions with h ≤ 43 W m⁻² K⁻¹ there is a sudden increase in the q_l/q_s due to increased contact of the propagating flame with the walls of the smaller channel. For all conditions of h ≤ 43 W m⁻² K⁻¹, the q_l/q_s ratio is always less than 3.53 for CH₄–air flames and it results in flame propagation in smaller channels. For h = 50 W m⁻² K⁻¹, the q_l/q_s > 3.53 and it leads to flame quenching. An isothermal wall condition of T_w = 350 K results in higher q_l/q_s ratio, leading to flame extinction in the smaller channel.

Fig. 14 Variation of heat-loss to heat-generation ratio near the contraction for different thermal boundary conditions (a) CH₄/air flames (b) H₂/air flames.

It can be seen in Fig. 14b that a q_l/q_s ratio less than 3.75 for H₂–air flames results in flame propagation through the channel step and can be considered as an upper limit for flame propagation through sudden channel step. Therefore, this value can be considered as an upper limit for H₂/air flames, which corresponds to h = 1000 W m⁻² K⁻¹, where propagating flame continues for about 5 mm before quenching. It is interesting to note that the value of q_l/q_s ratio for CH₄–air flames is 3.53 at h = 42 W m⁻² K⁻¹ and ∼3.75 for H₂–air flames at h = 1000 W m⁻² K⁻¹. More computations were carried out to understand the effect of varying D/d ratio on q_l/q_s ratio and flame propagation into smaller channels. A smaller tube diameter, d slightly greater than quenching distance, d_q with T_w = 300 K is considered. A suitably smaller D/d ratio = 1.14 was chosen and it was observed that the flame propagated through the contraction successfully with q_l/q_s ratio being much smaller, ∼1.28 for H₂–air and ∼1.38 for CH₄–air mixtures. These values compare well with the work of Gutkowski²⁶ where a value of 1.5 for C₃H₈–air flames, with a D/d ratio = 1.16 has been reported. Further, the heat losses are reported to be maximum for a stoichiometric mixture;¹⁰ therefore a q_l/q_s ratio of 3.75 and 3.5 may be considered as upper limits for H₂–air and CH₄–air flames, respectively for walls with no-slip boundary conditions.

Another notable observation for CH₄–air flames is that contrary to the expectation, at h = 42 W m⁻² K⁻¹, a slightly higher q_l/q_s ratio of 3.53 was observed as compared to 3.48 at h = 43 W m⁻² K⁻¹. This may be explained by recalling that h = 42 W m⁻² K⁻¹, had a positive flame shape angle and is more exposed to heat-losses than h = 43 W m⁻² K⁻¹, which transitions to a negative flame shape angle by the time it reaches the contraction. The positive flame shape angles or nearly tulip shaped flames being more exposed to quenching than mushroom shaped flames are in good agreement with the numerical results of Gutkowski.²⁶

3.4 Characteristic follower species in CH₄–air flames

To obtain the profiles of the characteristic follower species at a particular flame position, the maximum value of mass fraction of a particular species (Y_i) is obtained from the axial distribution along centerline for CH₄–air flames with h = 43 W m⁻² K⁻¹ thermal wall boundary condition. The entire distribution is normalized with their respective maximum values and plotted with respect to the non-dimensionalized flame position as shown in Fig. 15 and 16. It is interesting to note that various species of CH₄–air flames can be grouped into three different categories depending on the characteristic trends they follow – (a) the S_L follower species, i.e. the species which follow the axial variation of S_L (Fig. 15a), (b) the flame temperature follower species, i.e. the species which follow the trend of flame temperature (Fig. 15b), and (c) the heat release rate (HRR) follower species, i.e. the species which follow the trend of axial variation of maximum heat release rate (Fig. 16). It is noteworthy that while the trends of flame temperature and HRR are followed without any lag; in S_L follower species, the species follow S_L variation with a time lag of about 2 ms. The following results were computed for the case of h = 43 W m⁻² K⁻¹, because at this condition, the variation in flame shape angle and flame propagation speed is significant and trends are easily distinguishable. The trends of these species were also found to be true for other extreme wall heat transfer conditions of adiabatic and isothermal walls as well. Since, these trend follower species are also present in other hydrocarbon fuels, they are also expected to follow the trends for other hydrocarbon fuels such as propane and butane. However, while comparing the results of CH₄ and H₂–air flames, the common follower species being O and OH, the expected trend of these species for HRR was compared and shown in Fig. 16b. Other species in H₂–air flame were not observed to follow the temperature and flame speed profiles. Therefore, it can be concluded that the proposed species can be used to predict the global flame propagation properties such as flame propagation speed, flame temperature and heat release rate variation as the flame propagates in complex channel configurations.

Fig. 15 Characteristics follower species for h = 43 W m⁻² K⁻¹ condition in CH₄/air flames (a) S_L follower species (b) flame temperature follower species.

Fig. 16 HRR follower species (a) for wall h = 43 W m⁻² K⁻¹ in CH₄/air flames (b) for h = 0 W m⁻² K⁻¹ in H₂/air flames.

4. Conclusions

Unsteady numerical simulations with detailed chemistry have been carried out in stepped microtube combustors to understand the effect of various wall heat transfer conditions on flame propagation and flame extinction. A wide range of thermal boundary conditions between adiabatic, isothermal and heated isothermal walls have been studies using detailed chemistry for stoichiometric CH₄–air and H₂–air mixtures. It has been observed that the flame shape angle decreases significantly when the flame passes through the contraction. The density jump ratio across the flames is an important parameter governing the magnitude of hydrodynamic instability or D–L instability. The average temperature in the domain decreases due to increased wall heat losses resulting in reduced density ratio across the flames. A decrease in the magnitude of this instability results in a reduced flame shape angle of the propagating flame. A sudden acceleration of highly stretched flames near the channel step is due to sudden increase in mass flux entering the flame near the channel step. The ratio of heat-loss to heat-generation is a significant parameter governing the flame propagation in stepped channels. Quantified values of this ratio have been proposed for both CH₄ and H₂ fuels, which will help predict the continuous flame propagation or extinction in such channels. Three groups of follower species for CH₄–air flames, each following different characteristics of flame propagation velocity, flame temperature and maximum heat release rate, have been identified. For H₂–air flames, the identified follower species were observed to follow maximum heat release rate matched with that of CH₄–air flames.

Nomenclature

Adiabatic walls	Ad
Smaller tube diameter (mm)	d
Larger tube diameter (mm)	D
Quenching distance or its numerical value (mm) (mm)	dq
Convective heat transfer coefficient (W m^−2 K^−1)	h
Heat release rate (W m^−3)	HRR
In situ adaptive tabulation	ISAT
Flame propagating speed (m s^−1)	SL
Burnt gas temperature (K)	Tb
Unburnt gas temperature (K)	Tu
Wall temperature of isothermal walls (K)	Tw
Species mass fraction	Yi
Time step (s)	τ

Acknowledgements

This work was partially supported by the DST, Govt of India under joint Indo-Russian international research support wide grant no. INT/RFBR/P-179. The work was partially supported by project No. 13-06-0223 of FEFU and RFBR project No. 14-08-92695.

References

1. C. Fernandez-Pello, Micropower generation using combustion: issues and approaches, Proc. Combust. Inst., 2002, 29, 883–889.
2. Y. Ju and K. Maruta, Microscale combustion: Technology development and fundamental research, Prog. Energy Combust. Sci., 2011, 37(6), 669–715.
3. S. Kumar, K. Maruta and S. Minaev, Experimental investigations on the combustion behavior of methane–air mixtures in a microscale radial combustor configuration, J. Micromech. Microeng., 2007, 17, 900–908.
4. A. Fan, S. Minaev, S. Kumar, W. Liu and K. Maruta, Experimental study on flame pattern formation and combustion completeness in a radial microchannel, J. Micromech. Microeng., 2007, 17, 2398–2406.
5. A. Veeraragavan and C. P. Cadou, Flame speed predictions in planar micro/mesoscale combustors with conjugate heat transfer, Combust. Flame, 2011, 158(11), 2178–2187.
6. R. van Noorden, The rechargeable revolution: A better battery, Nature, 2014, 507(7490), 26–28.
7. N. Kim and K. Maruta, A numerical study on propagation of premixed flames in small tubes, Combust. Flame, 2006, 146, 283–301.
8. S. T. Lee and C. H. Tsai, Numerical investigation of steady laminar flame propagation in a circular tube, Combust. Flame, 1994, 99, 484–490.
9. K. Maruta, T. Kataoka, N. Kim, S. Minaev and R. Fursenko, Characteristics of combustion in a narrow channel with a temperature gradient, Proc. Combust. Inst., 2005, 30(2), 2429–2436.
10. L. Junwei and Z. Beijing, Experimental investigation on heat loss and combustion in methane/oxygen micro-tube combustor, Appl. Therm. Eng., 2008, 28, 707–716.
11. W. Yang, C. Deng, J. Zhou, J. Liu, Y. Wang and K. Cen, Experimental and numerical investigations of hydrogen–air premixed combustion in a converging diverging micro tube, Int. J. Hydrogen Energy, 2014, 39, 3469–3476.
12. D. Jiang, W. Yang and K. J. Chua, Entropy generation analysis of H₂/air premixed flame in micro-combustors with heat recuperation, Chem. Eng. Sci., 2013, 98, 265–272.
13. B. Khandelwal and S. Kumar, Experimental investigations on flame stabilization behavior in a diverging micro channel with premixed methane–air mixtures, Appl. Therm. Eng., 2010, 30, 2718–2723.
14. S. Kumar, Numerical studies on flame stabilization behavior of premixed methane–air mixtures in diverging mesoscale channels, Combust. Sci. Technol., 2011, 183, 779–801.
15. M. Akram and S. Kumar, Experimental studies on dynamics of methane–air premixed flame in meso-scale diverging channels, Combust. Flame, 2011, 158, 915–924.
16. V. Bychkov, V. Y. Akkerman, G. Fru, A. Petchenko and L. E. Eriksson, Flame acceleration in the early stages of burning in tubes, Combust. Flame, 2007, 150(4), 263–276.
17. D. M. Valiev, V. Y. Akkerman, M. Kuznetsov, L. E. Eriksson, C. K. Law and V. Bychkov, Influence of gas compression on flame acceleration in the early stage of burning in tubes, Combust. Flame, 2013, 160(1), 97–111.
18. C. Clanet and G. Searby, On the “tulip flame” phenomenon, Combust. Flame, 1996, 105(1), 225–238.
19. H. Xiao, R. W. Houim and E. S. Oran, Formation and evolution of distorted tulip flames, Combust. Flame, 2015, 162(11), 4084–4101.
20. M. Matalon and P. Metzener, The propagation of premixed flames in closed tubes, J. Fluid Mech., 1997, 336, 331–350.
21. U. W. Taywade, A. A. Deshpande and S. Kumar, Thermal performance of a microcombustor with heat recirculation, Fuel Process. Technol., 2013, 109, 179–188.
22. A. A. Deshpande and S. Kumar, On the formation of spinning flames and combustion completeness for premixed fuel–air mixtures in stepped tube combustors, Appl. Therm. Eng., 2013, 51, 91–101.
23. B. Khandelwal, A. A. Deshpande and S. Kumar, Experimental studies on flame stabilization in a three step rearward facing configuration based micro channel combustor, Appl. Therm. Eng., 2013, 58, 363–368.
24. B. Khandelwal, G. P. S. Sahota and S. Kumar, Investigations into the flame stability limits in a backward step micro scale combustor with premixed methane–air mixtures, J. Micromech. Microeng., 2010, 20, 095030.
25. M. Baigmohammadi, S. S. Sadeghi, S. Tabejamaat and J. Zarvandi, Numerical study of the effects of wire insertion on CH₄ (methane)/air pre-mixed flame in a micro combustor, Energy, 2013, 54, 271–284.
26. A. Gutkowski, Numerical analysis of effect of ignition methods on flame behavior during passing through as sudden contraction near quenching conditions, Appl. Therm. Eng., 2013, 54, 202–211.
27. Z. Liu, M. J. Lee and N. I. Kim, Direct prediction of laminar burning velocity using an adapted annular stepwise diverging tube, Proceedings of Combustion Institute, 2013, 34, 755–762.
28. Z. Liu and N. I. Kim, An assembled annular stepwise diverging tube for the measurement of laminar burning velocity and quenching distance, Combust. Flame, 2014, 161, 1499–1506.
29. Y. Kizaki, H. Nakamura, T. Tezuka, S. Hasegawa and K. Maruta, Effect of radical quenching on CH₄/air flames in a micro flow reactor with a controlled temperature profile, Proc. Combust. Inst., 2015, 35, 3389–3396.
30. C. M. Miesse, R. I. Masel, C. D. Jensen, M. A. Shannon and M. Short, Submillimeter-scale combustion, AIChE J., 2004, 50, 3206–3214.
31. D. G. Norton and D. G. Vlachos, Combustion characteristics and flame stability at the microscale: a CFD study of premixed methane/air mixtures, Chem. Eng. Sci., 2003, 58, 4871–4882.
32. M. O. Connaire, H. J. Curran, J. M. Simmie, W. J. Pitz and C. K. Westbrook, A comprehensive modeling study of hydrogen oxidation, Int. J. Chem. Kinet., 2004, 36, 603–622.
33. G. P. Smith, D. M. Golden, M. Frenklach, N. W. Moriarty, B. Eiteneer, M. Goldenberg, C. T. Bowman, R. K. Hanson, S. Song, W. C. Gardiner Jr, V. V. Lissianski, and Z. Qinhttp://www.me.berkeley.edu/gri_mech/, 1999.
34. G. Pizza, C. E. Frouzakis, J. Mantzaras, A. G. Tomboulides and K. Boulochos, Dynamics of premixed hydrogen/air flames in microchannels, Combust. Flame, 2008, 152, 433–450.
35. A. Fan, K. Maruta, H. Nakamura and W. Liu, Experimental investigation of flame pattern transitions in a heated radial micro-channel, Appl. Therm. Eng., 2012, 47, 111–118.
36. A. Fan, J. Wan, Y. Liu, B. Pi, H. Yao and W. Liu, Effect of bluff body shape on the blow-off limit of hydrogen/air flame in a planar micro-combustor, Appl. Therm. Eng., 2014, 62(1), 13–19.
37. K. T. Kim, D. H. Lee and S. Kwon, Effects of thermal and chemical surface-flame interaction on flame quenching, Combust. Flame, 2006, 146, 19–28.
38. ANSYS® FLUENT, Release 15, Theory Guide. Chapter 7 – Species Transport and Finite Rate Chemistry, ANSYS, Inc, 2013.
39. R. I. Issa, B. Ahmadi-Befrui, K. R. Beshay and A. D. Gosman, Solution of the implicitly discretised reacting flow equations by operator-splitting, J. Comput. Phys., 1991, 93, 388–410.
40. S. B. Pope, Computationally efficient implementation of combustion chemistry using in situ adaptive tabulation, Combust. Theory Modell., 1997, 1, 41–63.
41. A. Nair, V. R. Kishore and S. Kumar, Dynamics of premixed hydrogen–air flames in microchannels with a wall temperature gradient, Combust. Sci. Technol., 2015, 187(10), 1620–1637.
42. M. Akram, P. Saxena and S. Kumar, Laminar Burning velocity of methane–air mixtures at elevated temperatures, Energy & Fuels, 2013, 27, 3460–3466.
43. V. V. Bychkov and M. A. Libermann, Dynamics of stability of premixed flames, Phys. Rep., 2000, 325, 115–237.

---