1 Introduction

Flame-wall interaction is a common phenomenon in confined combustion chambers such as internal combustion engines and gas turbines (Dreizler and Böhm 2015). The presence of the wall significantly affects the flame structure, and causes the flame to quench, which leads to lower combustion efficiency and increased pollutant formation. In addition, high wall heat fluxes can become a major safety issue.

Head-on quenching (HOQ) and side wall quenching (SWQ) are two generic quenching configurations. Comparatively, looking at global flame characteristics, HOQ results in a smaller quenching distance and higher wall heat flux than SWQ as the flame front impinges vertically on the wall (Boust et al. 2007). This makes HOQ a well-defined limiting case suitable for parametric investigations and consequently this work focuses on the HOQ scenario.

The flow field in HOQ can be either laminar or turbulent. Turbulence adds another level of complexity to HOQ since turbulence-chemistry interaction and turbulence-wall interaction must be considered in addition to flame-wall interaction (Poinsot et al. 1993; Bruneaux et al. 1996; Zhao et al. 2018; Lai and Chakraborty 2016; Lai et al. 2018). Only flame-wall interaction is involved in laminar HOQ, which makes it a suitable configuration to gain deeper insight into flame-wall interaction. Also, results from laminar HOQ can support the development of turbulent near-wall combustion models. Besides, simple configuration such as the laminar HOQ are suitable for complementary experimental and numerical studies. Laminar HOQ has been extensively studied in previous studies (Boust et al. 2007; Hocks et al. 1981; Ezekoye et al. 1992; Popp and Baum 1997; Hasse et al. 2000; Dabireau et al. 2003; Chauvy et al. 2010; Strassacker et al. 2018; Jiang et al. 2019; Bellenoue et al. 2003; Foucher et al. 2003; Sotton et al. 2005; Mann et al. 2014).

Laminar HOQ is always a transient process where the premixed flame propagates through either initially stagnant mixtures (unstrained) or convective flows (strained). In this work, we use the term strain solely to characterize the underlying flow field, which then results in a stretched flame. Until now, both experimental (Bellenoue et al. 2003; Foucher et al. 2003; Sotton et al. 2005) and numerical (Boust et al. 2007; Poinsot et al. 1993; Lai and Chakraborty 2016; Lai et al. 2018; Hocks et al. 1981; Ezekoye et al. 1992; Popp and Baum 1997; Hasse et al. 2000; Dabireau et al. 2003; Chauvy et al. 2010; Jiang et al. 2019) studies related to laminar flame have almost exclusively focused on stagnant and therefore unstrained conditions. Specifically, the configuration of a freely propagating flame impinging on the wall has been studied for different fuel types, wall temperatures, and pressures in numerical simulations (Boust et al. 2007; Poinsot et al. 1993; Lai and Chakraborty 2016; Lai et al. 2018; Hocks et al. 1981; Ezekoye et al. 1992; Popp and Baum 1997; Hasse et al. 2000; Dabireau et al. 2003; Chauvy et al. 2010; Jiang et al. 2019). By contrast, there have been only a very limited number of investigations into laminar strained flow conditions (Bruneaux et al. 1996; Strassacker et al. 2018; Mann et al. 2014), and to the best of our knowledge no numerical studies with detailed chemistry and diffusive transport that conduct a systematic variation of strain rate exist in the literature. Bruneaux et al. (1996) simulated wall quenching of laminar flames in a stagnation line flow using one-step chemistry and simplified diffusion modeling. Mann et al. (2014) conducted the first measurements of the temperature and CO concentration under strained flow conditions. Strassacker et al. (2018) performed corresponding numerical simulations of this configuration using the Reaction-Diffusion-Manifold (REDIM) approach (Bykov and Maas 2007), however the underlying tabulated manifold was generated using gradient information from an unstrained HOQ configuration. The results showed reasonable agreement with the (strained) experimental data.

Previous numerical studies (Ganter et al. 2017; Strassacker et al. 2019) reported that differential diffusion is an important phenomenon during flame-wall interaction. For example, Ganter et al. (2017) reported that when using mixture-averaged diffusion, the simulation results show better agreement with the SWQ experimental data compared to the unity Lewis number assumption. Strassacker et al. (2019) performed an unstrained HOQ simulation, and also confirmed that detailed transport models lead to improved predictions compared to simplified transport models. However, the relevance for varying and especially high-strain conditions is still unknown.

Approaches using a tabulated manifold, e.g. REDIM- or flamelet-based, are widely-used for numerical simulations. While the flamelet model has been extended to consider strain rate effects in premixed flames (Knudsen et al. 2013; Oijen and Goey 2002; Trisjono et al. 2016; Peters 2000), it has not yet been investigated for near-wall combustion. To consider strain, an additional parameter is introduced. Different parameters including the H radical (Knudsen et al. 2013), CO (Oijen and Goey 2002), H2 (Trisjono et al. 2016) and the dissipation rate (Peters 2000) have been used to represent the strain rate effects. Consequently, we will extend these investigations to near-wall combustion by evaluating the suitability of these parameters to characterize the strain rate effects.

Based on the above literature review, the scope of the current study is to investigate transient HOQ of laminar premixed flames in a counterflow configuration for varying strain rates. The validity of the numerical approach is first evaluated through comparisons with experimental results. Then, strain rate effects on the global quenching quantities (quenching distance, wall heat flux) considering varying wall temperatures are studied. Next, the micro-structure of the quenching flame is analyzed through the local thermo-chemical state in CO/T phase space and near-wall differential diffusion. Finally, guidance is offered for developing a flamelet model for near-wall combustion.

1.1 Governing equations

The strained premixed HOQ problem can be described by the 1D unsteady opposed-flow equations for mass, axial momentum, radial momentum, energy and species (Im 2000), which are based on the conservation equations in cylindrical coordinate, with the assumption that the distributions of axial velocity, scaled radial velocity, temperature and species mass fractions are similar near the centerline, so that they are functions of time and the axial coordinate only: u = u(t,x), V = v/r = V(t,x), T = T(t,x), Yk = Yk(t,x).

Based on the assumption that V = v/r = V(t,x), mass continuity equation can be written as:

$$\frac{\partial \rho }{{\partial t}} + \frac{\partial }{\partial x}\left( {\rho u} \right) + \frac{1}{r}\frac{\partial }{\partial r}(\rho vr){\kern 1pt} {\kern 1pt} {\kern 1pt} { = }{\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} \frac{\partial \rho }{{\partial t}} + \frac{\partial }{\partial x}\left( {\rho u} \right) + 2\rho V = 0$$
(1)

The equation of state:

$$\rho = \frac{PW}{{RT}}$$
(2)

where the total thermodynamic pressure P is the sum of a reference pressure p0 and a varying pressure p, which is the dependent variable.

By differentiating Eq. (2), following relation can be obtained:

$$\frac{\partial \rho }{{\partial t}} = \frac{\rho }{P}\frac{\partial p}{{\partial t}} - \frac{\rho }{T}\frac{\partial T}{{\partial t}} - \rho W\sum\limits_{k} {\frac{1}{{W_{k} }}} \frac{{\partial Y_{k} }}{\partial t}$$
(3)

Substituting Eq. (3) into the Eq. (1), the mass continuity equation can be finally expressed as:

$$\frac{\rho }{P}\frac{\partial p}{{\partial t}} - \frac{\rho }{T}\frac{\partial T}{{\partial t}} - \rho W\sum\limits_{k} {\frac{1}{{W_{k} }}} \frac{{\partial Y_{k} }}{\partial t} + \frac{\partial }{\partial x}(\rho u) + 2\rho V = 0$$
(4)

Similarly, the axial momentum and radial momentum equations can also be simplified using the assumption u = u(t,x) and v/r = V(t,x), which finally read:

$$\rho \frac{\partial u}{{\partial t}} + \rho u\frac{\partial u}{{\partial x}} + \frac{\partial p}{{\partial x}} - 2\mu \frac{\partial V}{{\partial x}} - \frac{4}{3}\frac{\partial }{\partial x}\left( {\mu \frac{\partial u}{{\partial x}}} \right) + \frac{4}{3}\frac{\partial }{\partial x}(\mu V) = 0$$
(5)
$$\begin{gathered} \rho \frac{\partial V}{{\partial t}} + \rho u\frac{\partial V}{{\partial x}} + \rho V^{2} - \frac{\partial }{\partial x}\left( {\mu \frac{\partial V}{{\partial x}}} \right) + \Lambda = 0 \hfill \\ {\kern 1pt} \Lambda \equiv \frac{1}{r}\frac{\partial p}{{\partial r}} = \Lambda (t) \hfill \\ \end{gathered}$$
(6)

The energy equation can be simplified with the assumption that T = T(t,x), which finally leads to the equation:

$$\rho c_{p} \frac{\partial T}{{\partial t}} + \rho c_{p} u\frac{\partial T}{{\partial x}} - \frac{\partial }{\partial x}\left( {\lambda \frac{\partial T}{{\partial x}}} \right) - \frac{\partial p}{{\partial t}} - u\frac{\partial p}{{\partial x}} + \rho \left( {\sum\limits_{k} {c_{p} Y_{k} V_{k} } } \right)\frac{\partial T}{{\partial x}} + \sum\limits_{k} {\omega_{k} W_{k} h_{k} } = 0$$
(7)

For species equation, assuming Yk = Yk(t,x), it can be written as:

$$\rho \frac{{\partial Y_{k} }}{\partial t} + \rho u\frac{{\partial Y_{k} }}{\partial x} + \frac{\partial }{\partial x}(\rho Y_{k} V_{k} ) - \omega_{k} W_{k} = 0{\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} (k = 1,...,n_{s} )$$
(8)

In Eqs. (1)–(8), t stands for time. x is the axial direction and r is the radial direction. ρ denotes density. T is the temperature. Yk is the mass fraction of the k-th species. W is the mixture-averaged molecular weight while Wk stands for the molecular weight of the k-th species. u represents the axial velocity, while v is the radial velocity, and the scaled radial velocity is denoted as V. Λ is the eigenvalue of the system, which is time-dependent. P stands for the total pressure, which consists of a reference pressure p0 and a varying pressure p. μ is the dynamic viscosity. cp is the specific heat. λ stands for heat conductivity. ωk denotes the molar reaction rate of the k-th species. hk is the enthalpy of formation of the k-th species. Vk is the species’ diffusion velocity, which is considered by both mixture-averaged model and equal diffusivity model (unity Lewis number) in the present study. ns is the total number of species. R stands for the universal gas constant.

The above system of DAEs (differential–algebraic equations) is solved with the in-house solver ULF (Zschutschke et al. 2017). A second-order central differencing scheme is used for spatial discretization. Time integration is implemented with an implicit, variable-order, variable-step method based on the Gear multi-value method (Buzzi Ferraris and Manca 1998) which is suitable for stiff DAE systems. All simulations start from a steady solution as initial condition. GRI-Mech 3.0 (Smith et al. 1999) is used for the combustion chemistry.

2 Numerical implementation

2.1 Computational domain

The configuration studied is schematically shown in Fig. 1. A stoichiometric homogeneous methane/air mixture is injected from a nozzle, and propagates towards the wall. The bulk flow velocity at the nozzle exit is U0. Consistent with the experimental setup (Mann et al. 2014), the premixed flow is ignited downstream by a spark, establishing a strained premixed flame propagating towards the wall.

Fig. 1
figure 1

Schematic of the system configuration for the transient HOQ process

The computational domain extends from the nozzle exit to the wall, which is set to be 10 mm. The coordinate system is defined as shown in Fig. 1 such that x is in axial direction, being zero at the wall surface and positive in the upstream direction, while r is in the radial direction. The domain is initially filled with a stoichiometric homogeneous methane/air mixture at 1 atm. To avoid preheating, the unburnt gas and the wall temperature are set to be equal. The computational domain is discretized by 1200 grid points with local refinement near the wall, resulting a mesh resolution of 5 μm in the near-wall region.

2.2 Boundary conditions

At the inlet, the axial velocity u is set to be constant, which equals to the value at the nozzle exit U0. The scaled radial velocity V is set to zero. The temperature T and species mass fractions Yk are set to be the same of the unburnt gas.

The axial velocity u is set to zero at the wall, while zero-gradient boundary condition is applied for scaled radial velocity V. The wall temperature Tw is set to be constant. Consistently with a previous study (Jiang et al. 2019), the wall is assumed to be inert and consequently a zero-diffusion flux boundary condition ρVkYk = 0 is employed for all species.

Similarly to (Escudié et al. 1999), the strain rate K is defined as the mean axial velocity gradient near the wall, which can be varied by adjusting the inlet velocity. Since the value of K changes during the transient flame propagation, it is calculated based on the initial state of the flow field prior to ignition. It is found that the wall temperature and the diffusion model have only little impact on the exact value of the strain rate. In the present study, four different strain rates of 600 s−1, 1200 s−1, 1800 s−1, and 2400 s−1 are investigated for two distinctively different wall temperatures relevant for confined combustion (Bruneaux et al. 1996) (Tw = 300 K, 600 K), as summarized in Table 1. The wall distance where the axial flow velocity |u| equals the unstretched adiabatic laminar burning velocity sL0 (Tw = 300 K, mixture-averaged model) is also included for reference.

Table 1 Summary of investigated conditions

2.3 Dimensionless parameters

For comparing the current results with the literature and to quantify the effect of strain rate, dimensionless parameters previously used in quenching studies are employed.

The distance between the flame front and the wall (xf) is often expressed in its dimensionless form as a Peclet number (Popp and Baum 1997)

$$Pe = \frac{{x_{f} }}{\delta },$$
(9)

where δ is the thickness of a laminar adiabatic unstretched flame δ = λu/(ρucpsL0), where λu and ρu are the thermal conductivity and density of the unburnt gas, respectively. cp is the specific heat capacity per unit mass. The flame front position corresponds to the location of the maximum heat release rate, which is consistent with previous studies (Dabireau et al. 2003; Palulli et al. 2019).

The instant of quenching (tQ) is defined as the time instant of the maximum wall heat flux \(q_{wall,\max } = (\left. {{{\lambda \partial T} \mathord{\left/ {\vphantom {{\lambda \partial T} {\partial x}}} \right. \kern-\nulldelimiterspace} {\partial x}}} \right|_{wall} )_{\max }\), which is consistent with other numerical studies (Strassacker et al. 2019). The dimensionless form of wall heat flux is \(Q = {{q_{wall} } \mathord{\left/ {\vphantom {{q_{wall} } {[\rho s_{L}^{0} c_{p} (T_{ad} - T_{u} )]}}} \right. \kern-\nulldelimiterspace} {[\rho s_{L}^{0} c_{p} (T_{ad} - T_{u} )]}}\), where Tad is the adiabatic flame temperature and Tu is the temperature of the unburnt gas. The flame-wall distance at quenching instant is defined as quenching distance xf,Q = xf (t = tQ), with the corresponding Peclet number PeQ = xf,Q /δ. Based on the definitions above, the dimensionless time τ is expressed as

$$\tau = \frac{{t - t_{Q} }}{{\delta /s_{L}^{0} }}.$$
(10)

Besides, the wall distance x is normalized as

$$D_{w} = \frac{x}{\delta }.$$
(11)

3 Results and discussion

In this section, the available experimental data is firstly used to validate the numerical solver. Using the validated solver, the effects of strain rates and wall temperature are investigated in a systematic manner next.

3.1 Validation and comparison to experiments

The simulation results are first compared to the experimental data reported in (Mann et al. 2014), which includes the temperature and CO concentration measured through two-photon laser-induced fluorescence (LIF) and nanosecond coherent anti-Stokes Raman spectroscopy of nitrogen (CARS), respectively. Unless specifically indicated, otherwise, all simulations use a mixture-averaged model for the diffusion velocities.

Figure 2 shows the CO mole fraction XCO over the temperature T for different positions, such plots in phase space were originally presented in (Mann et al. 2014) and are also used here for comparison. Also indicated are the two branches for CO formation and oxidation, which are characteristic for premixed flames, more details can be found in (Kosaka et al. 2018). Experimental results are shown as scatter data and a corresponding mean; the variance from the fluctuations is shown as an error bar. At different locations, the trends in the variation of CO against T are similar. At the beginning, the CO level increases along with the temperature (CO formation branch). After the peak value is reached, CO decreases due to oxidation. The CO level drops with initially increasing and later decreasing temperature depending on the wall distance. However, non-negligible differences between the simulations and experiments can be observed for positions near the wall, which is consistent with the findings in (Strassacker et al. 2018). This can be seen mostly for the oxidation branch, e.g. Fig. 2b, c, but also for the formation branch in the very near-wall position, as in Fig. 2a. These differences can be partially attributed to experimental uncertainties for wall distance determination. The wall in the experiments is convex which is part of a sphere with a radius of Rwall = 300 mm. The maximum value of the uncertainty is estimated to be 100 μm (Dreizler A 2019). Thus, we additionally show simulation results shifted by 100 μm, which show very good agreement with experiments for all wall distances. It is evident that in the very near-wall region, the shift has a significant effect on the CO-T variation, see e.g. at x = 0.1 mm. However, such effects become almost insignificant for locations far away from the wall, e.g. for x = 1.7 mm, because the wall influence on the flame structure decreases with increasing distance from the wall. Further away from the wall, the wall influence can be neglected and both CO formation and oxidation branches correspond to a freely-propagating flame. For these regions, discrepancies between numerical results and experimental data, especially in the formation branch, may be attributed to the chemistry mechanism and the influence of initial ignition which is generated by placing an artificial time-decaying energy source near the inlet, with the latter effect becoming negligible in the near-wall region.

Fig. 2
figure 2

CO mole fraction over temperature for different positions until 0.5 ms after quenching (numerical results are obtained using mixture-averaged diffusion model). The error bar is obtained according to the experimental fluctuations (Kosaka et al. 2018)

Figure 3a shows the evolution of the temperature as a function of the dimensionless time τ and normalized wall distance Dw. It can be clearly seen that as the flame propagates towards the wall, the near-wall temperature increases. After quenching, the temperature first decreases in the near-wall region, and then spreads towards regions away from the wall. The corresponding evolution of the CO mole fraction is shown in Fig. 3b. It can be seen that CO does not vanish during or even after quenching, and it remains high within the quenching distance, which is consistent with the experimental data as reported above.

Fig. 3
figure 3

Dw-τ plots of a temperature and b CO mole fraction for experimental setup. The horizontal line represents the normalized quenching distance. The vertical line represents the quenching instant

In summary, the comparison between the numerical and experimental results show favorable agreement. Furthermore, our numerical results for this reference case are comparable to previously reported works, and differences to experimental results can at least partially be attributed to uncertainties in the determination of the wall position. In the following section, the influences of strain on the quenching process both in terms of global quantities and local thermochemical structures are investigated.

3.2 Global quenching parameters and local thermo-chemical state

As indicated above the quenching distance is the flame-wall distance at time instant of quenching. It is a key parameter quantifying the thickness of the unburnt layer. Dimensionless quenching distances (PeQ) for different strain rates and different wall temperatures, see Table 1, are plotted in Fig. 4a, where the subscript Q denotes the instant of quenching. In addition to the mixture-averaged diffusion model, results obtained with the unity Lewis number assumption are shown for comparison. Overall, PeQ decreases with an increasing strain rate, which is a direct consequence of the underlying flow pushing the flame towards the wall. However, for a higher wall temperature, these strain rate effects on PeQ become less pronounced. For example, for Tw = 300 K, PeQ decreases from 2.66 to 1.96 for a strain variation of 600 s−1 to 2400 s−1 using a mixture-averaged diffusion model. By contrast, for Tw = 600 K and the same diffusion model, the maximum difference in PeQ is only 0.01. Comparing PeQ for different wall temperatures for the same strain rate, it can be seen that PeQ decreases as the wall temperature increases. This specific trend is consistent with previous findings in the literature for unstrained HOQ (Norden 1965). The gap becomes narrower as the strain rate increases. Finally, for both Tw = 300 K and Tw = 600 K, the unity Lewis number assumption leads to an underestimation of PeQ compared to the mixture-averaged model, which are mainly due to the differences in absolute values of quenching distance xf and laminar burning velocity s 0L . This confirms that a unity Lewis number is not suitable for laminar HOQ of methane-air flames.

Fig. 4
figure 4

Normalized global quantities: a quenching distances, b quenching wall heat flux

The non-dimensional wall heat flux is another key quantity in flame-wall interaction. Following the definition above and consistent with previous works, QQ denotes the normalized wall heat flux at the quenching instant. Figure 4b shows the normalized quenching wall heat flux for different strain rates and different wall temperatures. Overall, a higher strain rate results in a larger wall heat flux. This is closely related to the variation in the PeQ number. It is interesting to note that the variation in QQ is more distinct for smaller strain rates. At higher wall temperatures, the normalized wall heat flux almost stays constant at a value of 0.66 for varying strain rates, which is consistent with the unstrained value reported in (Popp and Baum 1997). The difference between the results calculated using unity Lewis number assumption and mixture-averaged diffusion becomes significant for lower wall temperatures and higher strain rates, which is mainly caused by the differences in laminar burning velocity s 0L , the adiabatic flame temperature Tad, as well as the near-wall temperature gradient. This finding again confirms that the unity Lewis number assumption is unsuitable here and therefore it will not be considered further in the remainder of the study.

After the discussion of the global quantities, we now turn our attention to the local thermochemical state in vicinity of the wall and how it is affected by strain. Figure 5 presents CO-T variations under different strain rates for different positions (Dw = 1, 2, 3), and two wall temperatures of 300 K and 600 K are shown. The positions are chosen to be within the quenching region, see the discussion on the Peclet number above. The results are recorded until τ = 0.5. It is seen that CO profiles vary strongly at the different positions and for the two wall temperatures. Closer to the wall, both CO formation and oxidation happen under lower temperature conditions.

Fig. 5
figure 5

CO-T plot for different positions until 0.5 dimensionless time after quenching

From these results it is evident that strain rate effects on the thermo-chemical state are especially obvious for Tw = 300 K. At a higher strain rate, both the formation and the oxidation branch move to a higher temperature. Close to the wall, the strain rate has a larger influence on the formation branch than on the oxidation branch. Further away from the wall, the strain rate effect on CO-T decreases and remains only important for the oxidation branch. This can be confirmed by looking specifically at Dw = 3, the formation branches for different strain rates almost overlap.

Summarizing the analyses of quenching distances, wall heat flux and CO-T variation, the strain rate has a significant effect at lower wall temperatures. Contrarily, the effect almost diminishes at higher wall temperatures and this applies both for the global quantities and the local microstructure. Therefore, in the subsequent discussion, we focus on cases with Tw = 300 K to investigate the impact of strain on differential diffusion.

3.3 Strain rate effects on differential diffusion

As a first step, we quantify the impact of differential diffusion on the quenching process for the case Tw = 300 K and K = 600 s−1. Then, we analyze the influence of varying strain rates.

For this analysis, we use the differential diffusion parameter ZHC to quantify the degree of differential diffusion

$$Z_{{{\text{HC}}}} = F_{{\text{H}}} - F_{{\text{C}}} ,$$
(12)

where FH = (YH − YH,1)/(YH,2 − YH,1) and FC = (YC − YC,1)/(YC,2 − YC,1) (Barlow et al. 2000). Y stands for the mass fraction of the respective element (H or C), and subscripts 1 and 2 refer to pure fuel (here methane) and air (i.e. YC,2 = YH,2 = 0), respectively. Larger absolute values of ZHC indicate stronger differential diffusion. Note that the main focus of the present study is the variation of ZHC, rather than the ZHC value itself, which might strongly differ for other fuels.

Figure 6 shows the ZHC profile as a function of the normalized wall distance for different times. Symbols on the lower abscissa denote the flame front position at the respective time (indicated by the color). At different times, the minimum (largest absolute) value of ZHC always occurs ahead of the position of maximum heat release towards the preheat zone, i.e. between the reaction zone and the wall. As the flame approaches the wall, the minimum of ZHC moves towards the wall and differential diffusion within the quenching distance becomes more pronounced. The smallest value of ZHC is found just prior to quenching, and it exceeds the value of the propagating flame by a factor of up to 2 in terms of absolute numbers. This is consistent with the findings above concerning the change of the local thermochemical state during quenching.

Fig. 6
figure 6

Differential diffusion parameter profile at different times for Tw = 300 K, K = 600 s−1

For the subsequent analysis of differential diffusion during HOQ, the balance equation for differential diffusion parameter ZHC is considered. Using Eq. (12), the total derivative can be expressed as

$$\rho \frac{{DZ_{{HC}} }}{{Dt}} = \frac{1}{{Y_{{H,2}} - Y_{{H,1}} }}\rho \frac{{DY_{H} }}{{Dt}} - \frac{1}{{Y_{{C,2}} - Y_{{C,1}} }}\rho \frac{{DY_{C} }}{{Dt}}.$$
(13)

The balance equations for elemental mass fractions are

$$\rho \frac{{DY_{H} }}{{Dt}} = - \nabla \cdot \left( {\sum\nolimits_{{k = 1}}^{{n_{s} }} {\frac{{a_{{kH}} W_{H} }}{{W_{k} }}} \rho Y_{k} V_{k} } \right),$$
(14)
$$\rho \frac{{DY_{C} }}{{Dt}} = - \nabla \cdot \left( {\sum\nolimits_{{k = 1}}^{{n_{s} }} {\frac{{a_{{kC}} W_{C} }}{{W_{k} }}} \rho Y_{k} V_{k} } \right),$$
(15)

where W stands for the molecular weight. Vk is the diffusion velocity of species k. akH is the number of atoms of element H in a molecule of species k, akC is the number of atoms of element C in a molecule of species k.

For a unity Lewis number assumption, the diffusion term for ZHC can be expressed as \(- \nabla \cdot (\rho D\nabla Z_{{HC}} )\), where D is the diffusion coefficient of ZHC and all species. Combining Eqs. (13)–(15), the governing equation for ZHC can be rearranged as

$$\begin{gathered} \rho \frac{{DZ_{{HC}} }}{{Dt}} - \nabla \cdot (\rho D\nabla Z_{{HC}} ) \hfill \\ \quad = \frac{{ - \nabla \cdot \left( {\sum\nolimits_{{k = 1}}^{{n_{s} }} {\frac{{a_{{kH}} W_{H} }}{{W_{k} }}} \rho V_{k} Y_{k} } \right)}}{{Y_{{H,2}} - Y_{{H,1}} }} - \frac{{ - \nabla \cdot \left( {\sum\nolimits_{{k = 1}}^{{ns}} {\frac{{a_{{kC}} W_{C} }}{{W_{k} }}} \rho V_{k} Y_{k} } \right)}}{{Y_{{C,2}} - Y_{{C,1}} }} - \nabla \cdot (\rho D\nabla Z_{{HC}} ) \hfill \\ \end{gathered}$$
(16)

The right-hand side of Eq. (16) represents the deviation from a scalar with unity Lewis number and therefore it is a direct measure for the importance of differential diffusion.

Figure 7 presents the evolution of the source term at different times up to quenching for K = 600 s−1 (a) and K = 1800 s−1 (b), respectively. For reference, the ZHC profiles are also plotted. Prior to quenching, the shape of the source term remains similar for all times with negative values just ahead of the flame front and positive values around the reaction zone, which is the expected profile for propagating flames without wall influence. As long as the flame is far away from the wall, the source term near the wall is zero, see e.g. τ = − 5 for both strain rates. As the flame approaches the wall, the source term in the near wall region begins to decrease until the instant of quenching, when the positions of minimum source term and minimum ZHC (see also Fig. 6) are both located at the wall. Compared to the propagating flames, both the source term and ZHC decrease to smaller values (larger absolute values).

Fig. 7
figure 7

Source term (dotted line) and ZHC (solid line) profiles at different times for a K = 600 s−1, and b K = 1800 s−1

Both maximum and minimum ZHC source terms are amplified by higher strain rate for the propagating flames, resulting in a small difference in minimum and maximum ZHC value. However, when the flame interacts with wall, the minimum ZHC value for K = 600 s−1 at τ = 0 is significantly smaller compared to K = 1800 s−1 meaning that strain dampens differential diffusion near the wall during quenching. In contrast, the source term decreases with increasing strain rate, which reflects the larger magnitude of species’ gradients near the wall.

This is further analyzed in Fig. 8, where the minimum ZHC value (ZHCmin) is plotted as a function of time, as is the location Dw|ZHC = ZHCmin. When the flame is far away from the wall, the minimum ZHC value is almost constant during propagation, confirming the rather small influence of strain as discussed above. However, when the flame enters the near-wall region (see discussion on Peclet number above), a minimum of ZHCmin occurs right before quenching with a slight delay (τ values closer to zero) for higher strain rates. Further, the higher the strain rate, the later starts the departure from the propagating flame and this is also reflected in the trajectory of Dw |ZHC = ZHCmin. The location of the minimum ZHC value occurs later in the quenching region for higher strain rate. In fact, this explains why ZHCmin does increase with strain. Although the source term is amplified by strain, the smaller residence time in the quenching layer leads to dampening of differential diffusion.

Fig. 8
figure 8

Minimum ZHC-value (solid line) within Dw = 0–20 and its location Dw |ZHC = ZHCmin (dashed line) plotted over time

3.4 Guidance to near-wall flamelet modeling

In strained premixed flamelet modelling, an additional parameter is introduced to characterize strain rate effects (Knudsen et al. 2013; Oijen and Goey 2002; Trisjono et al. 2016; Peters 2000). It should be sensitive to strain rate, and it is crucial that there is no overlap between different profiles, so that a unique value can be found for a table look-up. This section mainly focuses on finding a suitable strain parameter for near-wall flamelet modeling.

In Fig. 9, three different scalars previously used for strain parameterization, YH, YCO and \(Y_{{{\text{H}}_{2} }}\), are plotted against progress variable. Here, the progress variable (PV) is defined as PV =  \(Y_{{{\text{CO}}_{2} }}\) + \(Y_{{{\text{H}}_{2} {\text{O}}}}\)  + YCO + \(Y_{{H_{2} }}\) (Knudsen et al. 2013). Two time instants are chosen: τ = − 1 represents the time before quenching, while the other time τ = 0 corresponds to the quenching instant. As shown in Fig. 8a, before quenching, both YH and \(Y_{{H_{2} }}\) show a strong sensitivity to the strain rate, while YCO varies only little. It is hard to distinguish YCO between different strain rates, so YCO is not a suitable choice for the time before quenching. At the instant of quenching, YH does not show a monotonic variation with the strain rate as before. Different profiles overlap, so a lack of uniqueness will be encountered for tabulation. At the quenching instant, \(Y_{{H_{2} }}\)is a suitable strain rate parameter. Although the variation trend is in opposition to the time before quenching, the monotonicity over the strain rate is retained. Based on the above analysis for the cases studied, H2 is suitable to parameterize strain rate for the near-wall flamelet modelling.

Fig. 9
figure 9

Different species plotted against progress variable for different strain rates prior and at quenching

4 Conclusions

The present numerical study investigated the transient head-on quenching (HOQ) of stretched premixed laminar methane-air flames using detailed chemistry and transport. While previous studies in the literature focused on unstrained conditions, this work is the first study of laminar HOQ systematically varying the strain rate for different wall temperatures. After validating the numerical approach against experimental data, this work specifically focused on four aspects.

  1. 1.

    The impact of strain rate on the global quantities: The normalized quenching distance expressed as Peclet number PeQ decreases with increasing strain rate, while the normalized wall heat flux increases. This can be directly attributed to the underlying flow field supporting flame propagation towards the wall. Modelling species’ diffusion with a unity Lewis number is not sufficient and this applies for low and high strain rates. The impact of strain is generally higher for low wall temperatures while it is negligible at 600 K.

  2. 2.

    The impact of strain on the near-wall microstructure: the local flame structure is analyzed with CO-T profiles in phase space at various wall distances. Strain influences both the CO formation and oxidation branches. Close to the wall, strain affects mostly CO formation while further away, the oxidation branch is shifted to higher temperatures with increasing strain. The impact of strain is generally more pronounced for lower wall temperatures and it is negligible at 600 K, which is consistent with the findings for the global quantities above.

  3. 3.

    The impact of strain on near-wall differential diffusion: the effect of non-equal diffusivities is quantified by the parameter ZHC. Differential diffusion is substantially enhanced in the quenching layer compared to propagating flames. While strain has little influence on the minimum/maximum ZHC values in the propagating flame, higher strain rates do suppress the increase of differential diffusion in the near-wall quenching layer. This aspect is further quantified by studying the source term in the transport equation for ZHC. While the source term is amplified by strain, the main effect is the reduction of residence time in the near-wall region for higher strain rates.

  4. 4.

    Parameterization of strain rate for manifold-based tabulation approaches: It was shown in previous works that the effect of strain on the structure of a stretched premixed flame can be parameterized by a reactive scalar. This finding is extended here for quenching flames. Results show that the H2 is well-suited for parameterizing strain rate for propagation and quenching, while other scalars previously used for propagating flames are not applicable during quenching. While the parameterization should only be used for the conditions tested, our findings can be considered a good starting point when looking at other conditions.