A Priori Direct Numerical Simulation Modeling of Scalar Dissipation Rate Transport in Head-On Quenching of Turbulent Premixed Flames

ABSTRACT The statistical behavior and modeling of scalar dissipation rate (SDR) transport for head-on quenching of turbulent premixed flames by an inert isothermal wall have been analyzed in the context of Reynolds averaged Navier–Stokes simulations based on three-dimensional simple chemistry direct numerical simulation (DNS) data. It has been found that the density variation, scalar-turbulence interaction, reaction rate gradient, molecular diffusivity gradient, and molecular dissipation terms, i.e., , and , respectively, act as leading order contributors to the SDR transport away from the wall and the turbulent transport and molecular diffusion terms remain negligible in comparison to the other terms. The leading order contributors to the SDR transport have been found to be in a rough equilibrium away from the wall before the quenching is initiated but this equilibrium is not maintained during flame quenching. The predictions of the existing models for the unclosed terms of the SDR transport equation have been assessed with respect to the corresponding quantities extracted from DNS data. No existing models have been found to predict the near-wall behavior of the unclosed terms of the SDR transport equation. The models, which exhibit the most satisfactory performance away from the wall, have been modified to account for near-wall behavior in such a manner that the modified models asymptotically approach the existing model expressions away from the wall.


Introduction
The scalar dissipation rate (SDR) plays a key role in the closure of reaction rate and micro-mixing in premixed turbulent combustion (Bray, 1980). In the context of Reynolds averaged Navier-Stokes (RANS) and large eddy simulations (LES), the SDR is closed using either an algebraic expression in terms of known quantities (Dunstan et al., 2013;Gao et al., 2014Gao et al., , 2015aKolla et al., 2009;Mura and Borghi, 2003;Swaminathan and Bray, 2005) or by solving a transport equation (Chakraborty et al., 2008aGao et al., 2015bGao et al., , 2015cMantel and Borghi, 1994;Mura et al., 2008Mura et al., , 2009. Some of the aforementioned closures (e.g., Dunstan et al., 2013;Gao et al., 2014Gao et al., , 2015aKolla et al., 2009) have already been implemented in RANS and LES of a range of laboratory-scale (Ahmed and Swaminathan, 2014;Butz et al., 2015;Dong et al., 2013;Kolla and Swaminathan, 2010;Langella et al., 2015;Ma et al., 2014;Robin et al., 2010) and industrial (Sadasivuni et al., 2012) configurations, and satisfactory results have been obtained. However, all of the aforementioned modeling of SDR has been carried out for flows away from the wall and the effects of flame quenching by the wall on SDR transport and its modeling are yet to be analyzed in detail. This gap in the open literature has been addressed here by analyzing the SDR transport and its modeling using direct numerical simulations (DNS) data of head-on quenching of turbulent premixed flames by an isothermal inert wall.
To date, limited effort has been directed to the analysis of wall bounded reacting flows using DNS Rutland, 1998, 2002;Bruneaux et al., 1996Bruneaux et al., , 1997Dabireau et al., 2003;Gruber et al., 2010Gruber et al., , 2012Lai and Chakraborty, 2016;Poinsot et al., 1993). Out of the aforementioned DNS analysis, only Bruneaux et al. (1996) and Alshalaan and Rutland (1998) concentrated on the influences of the wall on flame surface density (FSD)-based reaction rate closure of premixed turbulent combustion. Recently, Lai and Chakraborty (2016) modified the SDR-based reaction rate closure by Bray (1980) for the purpose of flame-wall interaction and suggested modifications to an algebraic closure of SDR Kolla et al., 2009) for near-wall effects in the case of head-on quenching of turbulent premixed flames. Algebraic closures of SDR are derived based on the equilibrium of generation and destruction rates of scalar gradients but such an equilibrium may not be maintained under certain conditions (e.g., low Damköhler number combustion, combustion instability, etc.) and it may be necessary to solve a model transport equation for the SDR under those circumstances. The presence of an isothermal wall significantly affects the fluid dynamics and thermal transport Rutland, 1998, 2002;Bruneaux et al., 1996Bruneaux et al., , 1997Dabireau et al., 2003;Gruber et al., 2010Gruber et al., , 2012Lai and Chakraborty, 2016;Poinsot et al., 1993), which gives rise to the inequality between reaction progress variable c and nondimensional temperature T ¼T À T 0 À Á = T ad À T 0 ð Þ (whereT, T 0 , and T ad are the instantaneous dimensional temperature, unburned gas temperature, and adiabatic flame temperature, respectively) even for low Mach number unity Lewis number flames Rutland, 1998, 2002;Bruneaux et al., 1996Bruneaux et al., , 1997Lai and Chakraborty, 2016;Poinsot et al., 1993). Furthermore, turbulence significantly affects the heat flux at the wall Rutland, 1998, 2002;Bruneaux et al., 1996Bruneaux et al., , 1997Dabireau et al., 2003), which is expected to have significant effects on the correlation between reaction rate and scalar gradient, density variation, and scalar gradient alignment with local principal strain rates, which in turn are likely to be manifested in the statistical behavior of the various terms of the SDR transport equation. However, the effects of the wall on the closure of the unclosed terms of the SDR transport equation in the context of RANS have not been analyzed in the existing literature. This deficit is addressed in this article by carrying out three-dimensional (3D) DNS simulations of head-on quenching of statistically planar turbulent premixed flames with different values of Damköhler and Karlovitz numbers (i.e., Da and Ka). The modeling of near-wall physics in the context of RANS is also useful for the purpose of LES as the near-wall modeling in LES is often adopted from RANS, and most LES reduces to RANS in the near-wall region. Furthermore, in the present analysis DNS data has been used to address the inadequacies of the model predictions in the near-wall region, and the modifications to the existing models have been suggested so that they work satisfactorily both near to and away from the wall. Thus, the main objectives of the current analysis are: (1) to identify the near-wall effects on the unclosed terms of the SDR transport equation; and (2) to model the near-wall effects on the unclosed terms of the SDR transport equation.
The rest of the article is organized as follows. The mathematical and numerical background of this work will be presented in the next two sections. Following this, the results will be presented and subsequently discussed. The main findings will be summarized and the main conclusions will be drawn in the final section of this article.

Mathematical background
For the present analysis a single-step Arrhenius-type irreversible chemical reaction is considered for the purpose of computational economy as 3D DNS simulation with detailed chemistry remains prohibitively expensive for a detailed parametric analysis (Chen et al., 2009). The scalar field is represented by a reaction progress variable c, which can be defined in terms of a suitable product mass fraction Y P in the following manner: where the values in the unburned and burned gases are denoted by 0 and 1, respectively. According to Eq. (1), c increases monotonically from zero in the unburned gas to unity in the fully burned products. Bray (1980) derived the following closure for the mean reaction rate of reaction progress variable _ ω in terms of SDRε c ¼ ρDÑc 00 Á Ñc 00 =" ρ in the following manner for high Damköhler number (i.e., Da ) 1) flames: where ρ is the gas density; D is the progress variable diffusivity; f b c ð Þ is the burning mode probability density function (pdf); " Q,Q ¼ ρQ=" ρ, and Q 00 ¼ Q ÀQ are the Reynolds averaged, Favre averaged, and Favre fluctuation values of a general quantity Q; and the subscript L is used to refer to unstrained laminar flame quantities. The value of c m remains between 0.7 and 0.9 for typical hydrocarbon-air flames (Bray, 1980). Although Eq.
(2) was originally proposed for high Damköhler number (i.e., Da ) 1), it was subsequently shown by Chakraborty and Cant (2011) that Eq. (2) remains valid for low Damköhler number (i.e., Da <1) combustion as long as the flamelet assumption holds true. It is evident from Eq. (2) that e ε c needs to be modeled in order to close the mean reaction rate _ ω. Using transport equation of c, i.e., ρ Dc=Dt ð Þ¼ _ ω þ Ñ Á ðρDÑc), it is possible to derive a transport equation of the SDRε c ¼ ρDÑc 00 Á Ñc 00 =" ρ, which takes the following form (Borghi, 1990;Gao et al., 2014;2015b, 2015cMantel and Borghi, 1994;Mura and Borghi, 2003;Swaminathan and Bray, 2005): The terms on the left-hand side of Eq. (3a) represent the effects of unsteadiness and mean advection, respectively. The terms D 1 and T 1 denote the molecular diffusion and turbulent transport of e ε c , respectively. The term T 2 is the density variation term due to heat release, whereas the turbulence-scalar interaction term T 3 arises from the alignment of Ñc with local principal strain rates. The term T 4 signifies the chemical reaction contribution to the SDR transport and ÀD 2 ð Þ is the molecular dissipation term. The f D ð Þ term arises due to diffusivity gradients. The terms T 1 , T 2 , T 3 , T 4 , ÀD 2 ð Þ, and f D ð Þ are unclosed and need modeling for transport equation-based closure. The statistical behavior and the modeling of these terms will be discussed in detail in the "Results and Discussion" section of this article.

Numerical implementation
Three-dimensional DNS simulations of head-on quenching of premixed turbulent flames have been carried out using a compressible DNS code, SENGA (Jenkins and Cant, 1999), where the conservation equations of mass, momentum, energy, and reaction progress variable are solved in nondimensional form in a coupled manner. The simulations have been carried out for a domain of size 70:6δ Z Â 35:2δ Z Â 35:2δ Z (where δ Z ¼ α T 0 =S L is the Zel'dovich flame thickness with α T0 and S L being the thermal diffusivity of unburned gas and unstrained laminar burning velocity, respectively) with the largest side aligned with the mean direction of flame propagation (i.e., -ve x 1 -direction here). The simulation domain is discretized using a uniform Cartesian grid of size 512 Â 256 Â 256, ensuring 10 grid points across the thermal flame thickness δ th ¼ T ad À T 0 ð Þ =MaxjÑTj L . Most existing analyses of flame-wall interaction have been carried out for a constant wall temperature boundary condition Rutland, 1998, 2002;Bruneaux et al., 1996Bruneaux et al., , 1997Dabireau et al., 2003;Gruber et al., 2010Gruber et al., , 2012Poinsot et al., 1993). The same approach has been adopted here. Cooling of the walls is necessary in most combustors because the burned gas temperature is often higher than the melting point of the combustor material, and thus the constant wall temperature boundary condition has a practical relevance. A no-slip isothermal inert wall with temperature T w ¼ T 0 is taken to be placed on the left-hand side boundary in the x 1 -direction and the boundary opposite to the wall is taken to be partially nonreflecting. The mass flux normal to the wall is considered to be zero, which can be written in terms of reaction progress variable as @c=@x 1 ð Þ x 1 ¼0 ¼ 0 at the wall. The boundaries in x 2 and x 3 directions are taken to be periodic. The spatial discretization for the internal grid points has been carried out by a 10th-order central difference scheme but the order of differentiation gradually drops to a one-sided 2nd-order scheme at the nonperiodic boundaries (Jenkins and Cant, 1999). The temporal advancement has been carried out by an explicit low storage 3rd-order Runge-Kutta scheme (Wray, 1990). The reactive field is initialized by a steady unstrained planar laminar premixed flame solution. The flame is initially placed at a location where the influence either from or to the wall is marginal (in these cases T ¼ 0:9 isosurface for the unstrained planar laminar flame solution is kept at a distance of 20δ Z away from the wall), so that flame gets enough time to evolve in the presence of turbulence before interacting with the wall. For the purpose of initializing the velocity field, an initially homogeneous isotropic field of turbulent velocity fluctuations is generated using a pseudo-spectral method (Rogallo, 1981) following the Batchelor-Townsend spectrum (Batchelor and Townsend, 1948), but the velocity components at the wall u 1 , u 2 , and u 3 are specified to be zero to ensure a no-slip condition. This velocity field is allowed to evolve for an initial eddy turn-over time (i.e., t e ¼ l=u 0 ) before it interacts with the flame.
The initial values of normalized rms turbulent velocity fluctuation u 0 =S L , the ratio of turbulent integral length scale to thermal flame thickness l=δ th for the turbulent velocity field away from the wall are listed in Table 1 along with the corresponding values of Damköhler number Da ¼ lS L =δ th u 0 , Karlovitz number Ka ¼ ðu 0 =S L Þ 3=2 ðl=δ th Þ À1=2 , and turbulent Reynolds number Re t ¼ ρ 0 u 0 l=μ 0 , where ρ 0 and μ 0 are the unburned gas density and viscosity, respectively. It can be seen from Table 1 that the cases A, C, and E (B, C, and D) have the same values of Da (Ka) and Ka (Da) is modified to bring about the changes in Re t . The global Lewis number is taken to be unity for all cases considered here (i.e., Le ¼ 1:0). Standard values are chosen for Prandtl number Pr and ratio of specific heats γ (i.e., Pr ¼ 0:7 and γ ¼ 1:4). For the present analysis, both the heat release parameter τ ¼ T ad À T 0 ð Þ =T 0 and Zel'dovich number β ¼ E ac T ad À T 0 ð Þ =RT 2 ad are taken to be 6.0 (i.e., τ ¼ 6:0 and β ¼ 6:0). The gaseous mixture is assumed to follow the ideal gas laws and the flame Mach number Ma ¼ S L = ffiffiffiffiffiffiffiffiffiffi γRT 0 p is taken to be 0.014. The simulations for turbulent cases have been carried out for a time when the maximum, mean, and minimum values of wall heat flux assume identical values following the flame quenching. The simulation time remains different for different cases but the simulations for all cases were continued for t ! 12δ Z =S L , where 12δ Z =S L corresponds to 21, 30, 21, 15, and 21 initial eddy turnover times (i.e., t e ¼ l=u 0 ) for cases A-E, respectively. The nondimensional grid spacing next to the wall y þ ¼ u τ Δx=ν remains smaller than unity for all turbulent cases (the maximum value of y þ has been found to be 0.93 during the course of the simulation), where u τ ¼ ffiffiffiffiffiffiffiffiffi ffi τ w =ρ p , τ w , and ν are the friction velocity, mean wall shear stress, and kinematic viscosity, respectively. For y þ ¼ u τ Δx=ν % 0:93, the minimum normalized wall normal distance u τ x 1 =ν of T ¼ 0:9 isosurface has been found to be about 15:0 for the quenching flames considered here. This implies that the flame does not necessarily quench within the viscous sublayer.
The DNS data has been ensemble averaged on the transverse plane (i.e., x 2 and x 3 direction) at a given x 1 location for the purpose of extracting Reynolds/Favre averaged quantities. The Reynolds/Favre averaged quantities have been assessed for statistical convergence by comparing the values obtained based on full sample size with the corresponding values evaluated using the distinct half of the domain in the transverse direction, and the agreement has been found to be satisfactory. In the next section, the results based on full sample size will be presented for the sake of conciseness.

Flame wall interaction
The distributions of reaction progress variable c, nondimensional temperature T, and nondimensional reaction rate _ ω Â δ Z =ρ 0 S L in central x 1 À x 3 plane are shown in Figure 1. For unity Lewis number flame, c and T are identical to each other under adiabatic low Mach number conditions. This can be observed when the flame is away from the wall prior to the flame quenching. However, the adiabatic condition is not maintained in the case of isothermal boundary condition, and the equality between c and T does not hold during flame quenching. It can be seen from Figure 1 that the flame wrinkling increases t Da À1=2 when the flame is away from the wall.
This can be confirmed from the values of normalized flame surface area A T =A L (where flame area A is evaluated using the volume integral ð V Ñc j jdV and the subscripts T and L are used to refer to turbulent and laminar conditions) presented in Table 2 at different time instants. Table 2 indicates that an increase in u 0 =S L promotes a high extent of flame area generation at early times (i.e., t 4δ Z =S L ) when the flame is away from the wall. However, a high extent of flame wrinkling for large values of u 0 =S L enables the flame elements to reach close to the wall, which in turn leads to early initiation of flame quenching. This can be confirmed from Table 2, which shows that smaller values of A T =A L have been obtained for the cases with higher initial u 0 =S L at later times (e.g., t ! 8δ Z =S L ) due to earlier initiation of flame quenching. The aforementioned behavior can further be confirmed by comparing the distributions of c, T, and _ ω shown in Figure 1. It was shown previously (Lai and Chakraborty, 2016) that the flame quenches once T ¼ 0:9 isosurface reaches a certain nondimensional distance x 1 =δ Z normal to the wall and this nondimensional distance is often referred to as the minimum Peclet number Pe min (Lai and Chakraborty, 2016;Poinsot et al., 1993). For the unity Lewis number flames considered here the minimum Peclet number Pe min ð Þ L for laminar premixed flame head-on quenching is found to be 2.83 (Lai and Chakraborty, 2016), which is consistent with several previous experimental analyses (Huang et al., 1986;Jarosinsky, 1986;Vosen et al., 1984). It was demonstrated by Lai and Chakraborty (2016) that the minimum value of wall Peclet number Pe min for head-on quenching of turbulent Le ¼ 1:0 flames assumes approximately the same values as in the case of the planar laminar flame quenching. Lai and Chakraborty (2016) showed that _ ω drops significantly and eventually disappears for x 1 =δ z Pe min ð Þ L , which can be confirmed from the distributions of _ ω shown here. Interested readers are referred to Lai and Chakraborty (2016) for further discussion on temporal evolutions of wall heat flux and wall Peclet number in the current configuration for different values of Lewis number.
The present article will only concentrate on the closures of the unclosed terms of the SDR transport equation in the case of head-on quenching. Algebraic closure of SDR for head-on quenching has been addressed elsewhere (Lai and Chakraborty, 2016) and will not be repeated here. However, for the sake of completeness, the variation of normalized SDRε þ c (i.e.,ε c Â δ Z =S L ) with normalized wall normal distance x 1 =δ Z at different time instants is shown in Figure 2 for all cases considered here. It can be seen from Figure 2 that the distribution ofε þ c broadens with increasing u 0 =S L due to an increase in flame wrinkling. A comparison between Figures 1 and 2 reveals that the regions of high reaction Table 2. List of normalized flame surface area A T =A L at different stages of flame quenching for all cases considered here. progress variable gradient shifts towards the wall and thus the value of x 1 =δ Z , at which the peak value of normalized SDR is obtained, decreases as the flame propagates towards the wall. However, it can be discerned from the reaction progress variable c and _ ω distributions in Figure 1 that the gradient of c disappears upon flame quenching and thus the magnitude of normalized SDRε þ c drops significantly once the quenching starts andε þ c eventually vanishes (see t ¼ 8δ Z =S L and 10δ Z =S L ) after flame quenching. The near-wall behavior of SDRε c can be understood by examining the terms in the SDR transport equation, i.e., Eq. (3a). The variations of normalized values of the terms on the right-hand side of the SDR transport equation, Eq. (3a), i.e., D 1 , T 1 , T 2 , T 3 , T 4 , ÀD 2 ð Þ, 10δ Z =S L for cases A-E. The x 1 -axis is shown in log scale for the inset. Blue and black vertical lines indicate and f D ð Þ, with normalized wall normal distance x 1 =δ Z at different time instants are shown in Figure 3. The following observations can be made from Figure 3 regarding the SDRε c transport in the near-wall region: • The terms T 2 , T 3 , T 4 , ÀD 2 ð Þ, and f D ð Þ remain the leading order contributors to the SDRε c transport for all cases. However, the cases with high values of u 0 =S L exhibit high magnitudes of these terms.
• The chemical reaction rate term T 4 acts as a dominant source away from the wall for all cases considered here. However, the magnitude of T 4 drops significantly and eventually vanishes following flame quenching (e.g., t ! 8δ Z =S L Þ. • The dissipation term ÀD 2 ð Þ acts as the dominant sink for all cases but its magnitude also decreases in the near-wall region as a result of flame quenching. • The density variation term T 2 acts as the second dominant source term away from the wall for all cases. As a result of large heat loss through the wall, the reaction rate vanishes in the near-wall region, which leads to a weak/zero contribution of T 4 . However, T 2 assumes non-negligible values in the near-wall region because of significant variations of density in this region due to temperature change even in the absence of chemical reaction rate. • The scalar-turbulence interaction term T 3 acts predominantly as a sink term for cases A-D, but this term locally exhibits positive values in case E at early times. However, it assumes negative values as time progresses before vanishing altogether. • The diffusivity variation term f D ð Þ acts as a sink in the near-wall region but this term assumes small positive values away from the wall. This term especially plays the role of a dominant sink in the near-wall region. • It is evident from Figure 3 that the magnitudes of the molecular diffusion term D 1 and the turbulent transport term T 1 remain negligible in comparison to the magnitudes of T 2 , T 3 , T 4 , ÀD 2 ð Þ, and f D ð Þ. Swaminathan and Bray (2005) proposed the following scaling estimates for the unclosed terms of SDR transport equation where the velocity fluctuations and length scales associated with scalar fluctuations are scaled with respect to S L and δ Z , respectively, whereas the gradients of Favre/Reynolds averaged quantities are scaled with respect to integral length scale l: An alternative order of magnitude analysis was proposed by Mantel and Borghi (1994) involving the rms turbulent velocity fluctuation u 0 and Taylor micro-scale λ to scale the velocity fluctuations and the gradients of the fluctuating quantities, respectively, according to Tennekes and Lumley (1972): Furthermore, the following scaling estimate for ÀD 2 ð Þ in nonreacting flows can be obtained if the second derivatives of fluctuating quantities are scaled using the Kolmogorov length scale, η: Equation (4) indicates that T 2 , T 3 , T 4 , ÀD 2 ð Þ, and f D ð Þ are indeed expected to be leading order contributors to the SDR transport, and the magnitudes of D 1 and T 1 are expected to be negligible in comparison to T 2 , T 3 , T 4 , ÀD 2 ð Þ, and f D ð Þ for all cases irrespective of the value of Da and Re t . This can be substantiated from the behaviors of T 2 , T 3 , T 4 , ÀD 2 ð Þ, and f D ð Þ away from the wall for all cases considered here. It can be seen from Figure 3 that for all cases T 2 , T 3 , T 4 , and ÀD 2 ð Þ indeed scale with ρ 0 S 2 L =δ 2 Z away from the wall but the magnitudes of all of the unclosed terms of the SDR transport equation decrease once the quenching starts. The magnitude of f D ð Þ remains smaller than the terms T 2 ; T 3 ; T 4 , and ÀD 2 ð Þ, but it cannot be ignored. However, the magnitudes of all of the terms decay with time once the flame starts to interact with the wall and eventually vanishes after flame quenching. A comparison of the magnitudes of T 2 ; T 3 ; T 4 , ÀD 2 ð Þ, and f D ð Þ in Figure 3 reveals that the terms T 3 ; T 4 , ÀD 2 ð Þ, and f D ð Þ in the near-wall region assume negligible values in comparison to T 2 , which acts to generate e ε c in the near-wall region even where _ ω vanishes. Furthermore, it can be seen from Figure 3 that the terms (T 2 þ T 3 þ T 4 þ f D ð ÞÞ and ÀD 2 ð Þ remain roughly in equilibrium when the flame is away from the wall, but such an equilibrium is not maintained once the flame starts to quench. The modeling of the unclosed terms T 1 , T 2 , T 3 , T 4 , ÀD 2 ð Þ, and f D ð Þ will be discussed next in this article. Most models for the unclosed terms of the SDR transport equation have been proposed based on a-priori analysis of DNS data in canonical configurations, where there is no mean shear (Chakraborty et al., 2008aSwaminathan, 2010, 2013;Mantel and Borghi, 1994;Mura et al., 2008Mura et al., , 2009). The physical mechanisms that affect the statistical behavior of the various terms of the SDR transport remain mostly unchanged in the presence of the wall because SDR statistics are principally governed by small-scale physics, which are largely independent of mean-scale forcing. For example, dilatation rate is expected to play a key role in turbulent premixed combustion, even for small values of Mach number, irrespective of the presence of the wall. The statistical behavior of the scalar-turbulence interaction term T 3 is governed by the alignment of Ñc with the local principal strain rates, irrespective of the presence of the wall. Thus, many of the underlying physical mechanisms that affect the SDR transport remain unchanged in the presence of the wall. Therefore, the models that were originally proposed for premixed flames without walls account for some of the physical mechanisms, which are still valid in the presence of the wall. Furthermore, it has recently been found that the models for the unclosed terms of the SDR transport equation, which were originally proposed based on DNS data for turbulent premixed flames without any mean shear, also provide satisfactory predictions for a configuration (i.e., rod-stabilized V-flame), where mean shear is present (Gao et al., 2015c). Thus, the models that were originally proposed based on a-priori DNS analysis of premixed turbulent flames in canonical configuration without walls and mean shear have been considered as a starting point of this analysis because it is easier to modify the existing models to account for near-wall behavior, rather than switching to completely different models in the near-wall region.
Modeling of turbulent transport term T 1 It can be seen from Eqs. (4a) and (4b) that the magnitude of T 12 is expected to be negligible in comparison to that of T 11 for high values of Re t . Hence, T 1 in practical high Re t turbulent flows can be approximated as: In order to model turbulent transport T 1 , it is essential to model the turbulent flux ρu 00 j ε c , which is often modeled for passive scalar mixing using a gradient flux model for ρu 00 j ε c : where μ t ¼ " ρC μk 2 =ε is the eddy viscosity, C μ ¼ 0:09 model constant, and σ ε is the turbulent Schmidt number with turbulent kinetic energyk ¼ ρu 00 i u 00 i = 2" ρ ð Þ and its dis- However, it is well known that turbulent fluxes of the quantities related to turbulent scalar gradient (e.g., FSD AE ¼ Ñc j j and SDRε c ) exhibit counter-gradient (gradient) when turbulent scalar flux ρu 00 i c 00 shows counter-gradient (gradient) type behavior (Chakraborty and Cant, 2009;Swaminathan, 2010, 2013;Veynante et al., 1997). One obtains counter-gradient transport when the velocity jump due to flame normal acceleration dominates over turbulent velocity fluctuations and vice versa (Chakraborty and Cant, 2009;Veynante et al., 1997).  proposed a model (referred to as the T1CS model here) for ρu 00 i ε c in terms of ρu 00 i c 00 , which is capable of predicting both gradient and counter-gradient type transport: where the model parameters are given by λ c ¼ 2 and Φ ¼ 0:5 . Figure 4 shows that the T1CS model underpredicts the turbulent flux ρu 00 1 ε c in a region where x 1 =δ Z Pe min ð Þ L . As the flame propagates towards the quenching zone, the SDR flux ρu 00 1 ε c exhibits mainly gradient-type transport. The T1CS model assumes the transition of ρu 00 j ε c =ρu 00 j c 00 from a positive to a negative value atc % 0:5 by using the factor Φ Àc ð Þ. This transition is no longer valid in the near-wall region because of predominantly gradient transport and is also due to the absence of the effects of flame normal acceleration as a result of flame quenching. The T1CS model has been revised here in the following manner:  parameterized in such a manner that Eq. (5d) approaches Eq. (5c) away from the wall (i.e., x 1 =δ Z ) Pe min ð Þ L ). It can be seen from Figure 4 that the revised model captures both qualitative and quantitative behaviors of ρu 00 1 ε c both away from and close to the wall.

Modeling of density variation term T 2
The fluid density ρ in low Mach number combustion is given by (Bray et al., 1985): The nondimensional temperature T can be equated to the reaction progress variable c (e.g., T ¼ c) for globally adiabatic unity Lewis number flames in absence of the wall. Thus, under the aforementioned condition, ρ and " ρ can be expressed in terms of c andc as: Using Eq. (6b), T 2 can be expressed in the following manner (Chakraborty et al., 2008aSwaminathan and Bray, 2005): According to above equation, an alternative expression for T 2 can be obtained in the following manner for unity Lewis number flames: @c 00 @x j @u k @x k @c @x j |fflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflffl} According to the scaling argument by Swaiminathan and Bray (2005), the terms T 21 , T 22 , and T 23 can be scaled in the following manner: where the velocity fluctuation, gradients of fluctuation, and its mean quantities are scaled using S L , δ Z , and l, respectively. Equation (6e) suggests that T 22 and T 23 become negligible in comparison to T 21 for high values ofRe t . Consequently, T 2 can be taken to scale with ρ 0 S 2 L =δ 2 Z (Swaminathan and Bray, 2005). For unity Lewis number flame, the dilatation rate Ñ Áũ is scaled as Ñ Áũ , τS L =δ th (Chakraborty and Swaminathan, 2007a, 2007bChakraborty et al. 2008aSwaminathan and Bray, 2005), and SDR scales with c , S L =δ th (Chakraborty et al., 2008aSwaminathan, 2010, 2013;Gao et al., 2015c). The aforementioned scaling arguments have been taken into account in the model proposed by Chakraborty et al. (2008a) for T 2 , which takes the following form (T2CS): where C T 2 ¼ B T 2 = 1 þ Ka L ð Þ 0:5 is the model parameter with B T 2 ¼ 2:0 and Ka L % S L ð Þ À3=2ε δ th ð Þ 1=2 is the local Karlovitz number (Chakraborty et al., 2008aSwaminathan, 2010, 2013;Gao et al., 2015c). The Karlovitz number Ka L dependence of C T 2 accounts for weakening of T 2 magnitude for large values of Ka L (Chakraborty et al., 2008aSwaminathan, 2010, 2013;Gao et al., 2015c) due to diminished heat release effects as the broken reaction zones regime is approached (Peters, 2000).
The variations of T 2 with normalized wall normal distance x 1 =δ Z at different time instants are shown in Figure 5 for all cases considered here. Figure 5 shows that T 2 acts as a source term and assumes higher magnitudes before quenching than after quenching because most of the density variation occurs due to the chemical heat release. The magnitude of T 2 also diminishes as reaction rate sharply reduces close to the wall, but nonzero value of T 2 has been observed during quenching because of the density variation driven by the temperature change across the thermal boundary layer on the isothermal wall. It can be seen from Figure 5 that the T2CS model satisfactorily predicts T 2 extracted from DNS data before flame quenching (i.e., when the flame is away from the wall). However, the T2CS model gives rise to significant overprediction of T 2 when flame-wall interaction takes place. Here, the T2CS model has been modified in the following manner to account for flame-wall interaction: þ 1 account for the effects of the wall. For unity, Lewis number flamesc ÀT À Á vanish away from the wall butc ÀT À Á assumes nonzero values only in the near-wall region. This type ofc ÀT À Á dependence was used by Bruneaux et al. (1997) in the context of FSD modeling and the same approach has been adopted here. The exponent A 3 rises in the near-wall region, which acts to mimic the reduction of T 2 magnitude as a result of weakening of heat release effects arising from flame quenching. It is worth noting that A 2 and A 3 approach 1.0 and 0.5, respectively, away from the wall and Eq. (6g) becomes identical to the T2CS model, i.e., Eq. (6f), for x 1 =δ Z ) Pe min ð Þ L . It can be seen from Figure 5 that Eq. (6f) reduces the overprediction of T 2 , and its predictions are in better agreement with DNS data than the T2CS model in the near-wall region when the flame interacts with the wall. Furthermore, the prediction of Eq. (6g) becomes identical to the T2CS model away from the wall before the flame quenching is initiated.

Modeling of the turbulent scalar interaction term T 3
The variations of T 31 , T 32 , and T 33 with x 1 =δ Z at different time instants for cases A-E are shown in Figure 6. It can be seen that T 31 and T 33 assume predominantly negative values when the flame is away from the wall. The contribution of T 32 remains a dominant contribution to T 3 . However, as flame interacts with the wall, the magnitude of T 32 drops significantly. At the quenching stage, the magnitudes of T 31 and T 33 gradually become comparable to that of T 32 in the near-wall region. The local Damköhler number Da L 1 kS L =εδ th drops close to the wall due to the combination of the decay of turbulent kinetic energyk and a sharp increase of dissipation of turbulent kinetic energy in the near-wall region. This drop in Da L leads to an enhancement in magnitudes of T 31 and T 33 according to Eq. (4a). At the final stage of quenching, T 33 assumes positive values because of predominantly negative values of @ũ 1 =@x 1 as a result of the reversal of the direction of the flow (after quenching flow is directed towards the wall in contrast to the flow away from the wall before quenching). The statistical behavior of T 3 can also be explained by using the scalar-turbulence interaction contribution Λ: Figure 6. Variation of T 3 Â δ 2 Z =ρ 0 S 2 L with x 1 =δ Z along with its components T31, T32, and T33 at t = 4δ Z =S L ; 6δ Z =S L ; 8δ Z =S L ; and 10δ Z =S L for cases A-E. The x 1 -axis is shown in log scale for the inset.
Λ ¼ À 2ρD @c @x i @u i @x j @c @x j ¼ À2ρ e α cos 2 θ α þ e β cos 2 θ β þ e γ cos 2 θ γ À Á N c where e α , e β , and e γ are the most extensive, intermediate, and most compressive principal strain rates and θ α , θ β , and θ γ are the angles between Ñc and the eigenvectors associated with e α , e β , and e γ , respectively. According to Swaminathan and Bray (2005), the following scaling relation can be obtained: It can be deduced from Eq. (7b) that the quantity À 2ρD @c=@x j À Á @ũ i =@x j À Á @c=@x j À Á remains negligible in comparison to the contributions from T 31 , T 32 , and T 33 . It can be seen from Eq. (7a) that a preferential alignment of Ñc with e α (e γ ) leads to negative (positive) contributions of Λ and T 3 . Several previous analyses Swaminathan, 2007a, 2007b;Swaminathan and Grout, 2006) indicated that Ñc aligns preferentially with the most extensive principal strain rate e α , when the strain rate a chem induced by flame normal acceleration dominates over the effects of turbulent straining a turb . By contrast, a preferential alignment of Ñc with e γ occurs when turbulent straining a turb overwhelms the influences of strain rate a chem arising from flame normal acceleration. Scaling a chem and a turb by τS L =δ th and u 0 =l, alternatively u 0 =λ, respectively yields ða chem =a turb Þ , τDa, alternatively ða chem =a turb Þ , τDa=Re 1 2 t , τ=Ka Chakraborty and Swaminathan, 2007a, 2007b, which suggests that a chem is likely overwhelm a turb for large values of Da and/or τ. For the cases considered here, a chem dominates over a turb in spite of Da<1 due to the large value of τ. Thus, Ñc predominantly aligns with e α for all cases, which predominantly gives rise to negative values of T 3 , except at the final stage of flame quenching when T 3 assumes positive values close to the wall due to positive values of T 33 due to flow direction reversal. From the aforementioned discussion it is clear that the models of T 3 components need to account for Ñc alignment characteristics with local principal strain rates.
For the present analysis four existing models of T 31 have been considered for model comparison, which are listed as: T31MB (Mantel and Borghi, 1994): T31M1 (Mura et al., 2009): T31M2 (Mura et al., 2009): whereñ f ¼ Ñc= Ñc j j is a local flamelet normal vector.
T31CS : j c 00 @c @x j À C C τ Á S L " ρε c <ñ f Áx j > @c @x jc where Da Ã L ¼ S L ρ 0k = δ th " ρε ð Þ is the local density-weighted Damköhler number, C 1 ¼ 0:5 and C 2 ¼ 1:3Ka 2 L = 1 þ Ka L ð Þ 2 are the model parameters, and the model parameter C C is expressed as: It was demonstrated by  that the T31CS model captures both the qualitative and quantitative behaviors of T 31 for a range of values of Da, Ka, and Re t in the absence of the wall. The variations of T 31 with normalized wall normal distance x 1 =δ Z are shown in Figure 7 along with the predictions of the T31MB, T31M1, T31M2, and T31CS models for different time instants for all cases considered here. All of these models underpredict the magnitude of T 31 and the extent of this underprediction is particularly severe for the T31MB and T3M1 models. The agreement of the T3M2 and T3CS model predictions with DNS data remain better than the other models when the flame is away from the wall (i.e., t 4δ Z =S L ) and also when the quenching starts. Nonetheless, the T31M2 and T31CS models do not adequately predict T 31 extracted from DNS data in the near-wall region. The T31CS model starts to underpredict the DNS data once the quenching is initiated. In order to address this deficiency, the following modification to the T31CS model has been proposed here: j c 00 @c @x j À C C τ Á S L " ρε c <ñ f Áx j > @c @x jc are the model parameters, which account for the wall effects. The model parameter A 4 acts to reduce the overprediction of T 31 magnitude once the flame starts to interact with the wall (see t ! 6δ Z =S L in Figure 7). The model parameter A 5 becomes active in the near-wall region wherecÞT, and is responsible for damping the magnitude of T 31 close to the near-wall region. The term À C C τ Á S L " ρε c <ñ f Áx j > @c=@x j À Á c 1:5 is necessary to accurately predict T 31 away from the wall. However, it has a strong dependence onc 1:5 , andc changes rapidly in the near-wall region due to the interaction of the flame with the wall. The involvement of A 5 in Eq. (8f) reduces thisc 1:5 dependence close to the wall. It can be seen from Figure 7 that the model given by Eq. (8f) provides better performance than the other available models and thus is recommended for T 31 modeling. Mantel and Borghi (1994) proposed the following model for T 32 (T32MB): where A e ¼ 0:9 is a model parameter. Mura et al. (2008) proposed the following models for T 32 based on a-priori analysis of DNS data for high Da flames: Figure 7. Variation of T 31 Â δ 2 Z =ρ 0 S 2 L with x 1 =δ Z along with the predictions by the T31MB, T31CS, T31M1, and T31M2 models and Eq. (8f) (new model) at t = 4δ Z =S L ; 6δ Z =S L ; 8δ Z =S L ; and 10δ Z =S L for cases A-E. The x 1 -axis is shown in log scale for the inset.
where A M1 ¼ 1:0, C A ¼ 0:6, and C B ¼ 1:6 are the model constants and Da L ¼ S Lk =δ thε is the local Damköhler number. Equations (9b) and (9c) will henceforth be referred to as the T32M1 and T32M2 models in this article.  proposed a model for T 32 , which includes nonunity Lewis number effects (T32CS): Þε c accounts for the destruction of scalar gradient as a result of preferential alignment between Ñc and e α [see Eq. (7a)], and the local Karlovitz number Ka L dependence of C Ã 4 accounts for weakening of Ñc alignment with e α with increasing Karlovitz number due to diminishing influence of flame normal acceleration for high Karlovitz number combustion. The involvement of Le ÀP in the second term on the right-hand side of Eq. (9d) allows for strengthening of flame normal acceleration with decreasing Le . The involvement of 1 Àc ð Þ Φ Le ð Þ ensures that the qualitative behavior of T 32 variation withc is adequately captured Swaminathan, 2010, 2013). The predictions of the T32MB, T32M1, T32M2, and T32CS models are compared to T 32 extracted from DNS data in Figure 8. It can be seen from the variations of T 32 with normalized wall normal distance x 1 =δ Z in Figure 8 that the T32MB model fails to predict the negative values of T 32 for all cases. The performances of the T32M1 and T32M2 models remain comparable but their predictions remain an order of magnitude smaller than the corresponding quantity extracted from DNS data. The agreement between the T32CS model prediction with DNS data is better than the other models before quenching is initiated (i.e., when the flame remains away from the wall by the wall) in spite of overpredictions of the magnitude of the negative values of T 32 for the flames considered here. However, the quantity ðε=kÞ assumes large values close to the wall (due to augmentation ofεand decay ofk in the vicinity of the wall), which leads to severe overprediction of T 32 by all of the models in the near-wall region. In order to capture the near-wall behavior of T 32 the T32CS model has been modified in the following manner: Here, the modified parameters are A 6 ¼ 0:  the wall. Furthermore, the modifications suggested in Eq. (9e) reduce the extent of the overprediction of the magnitude of the negative values of T 32 in comparison to the T32CS model. Mantel and Borghi (1994) proposed the following model for T 33 : Chakraborty and Swaminathan (2007b) proposed an alternative model for T 33 : , and ϕ ¼ τc 1 Àc ð Þ= 1 þ τc ð Þ are the model parameters. Mura et al. (2009) proposed the following models for T 33 : The models given by Eqs. (10c)-(10e) will henceforth be referred to as the T33M1, T33M2, and T33M3, respectively. The predictions of the T33MB, T33CS, T33M1, T33M2, and T33M3 models are compared with DNS data in Figure 9. It can be seen from Figure 9 that the T33MB model overpredicts the magnitude of the negative contribution of T 33 for cases A and B when the flame is away from the wall. However, it performs satisfactorily away from the wall in cases C-E before flame quenching. However, the model T33MB underpredicts the magnitude of T 33 in the near-wall region once the flame quenching is initiated. It can be seen from Figure 9 that the T33M1 model satisfactorily predicts T 33 for cases A and B when the flame is away from the wall as well as at the quenching stage. The models T33M2 and T33M3 overpredict the magnitude of the negative value of T 33 when the flame is away from the wall (e.g., t 4δ Z =S L ), nonetheless, the T33M2 model underpredicts whereas the T33M3 model significantly overpredicts the magnitude of T 33 during the final stage of quenching.
It is worth noting that the T33M1 model is consistent with the scaling arguments given by Eqs. (4a) and (4b). By contrast, the T33M2 model is consistent with the scaling given by Eq. (4b) (i.e., T 33 , O ρ 0 u 0 2 =l 2 À Á ), but the scaling arguments by Swaminathan and Bray (2005) Moreover, the T33M3 model scales as ρ 0 S 2 L =δ 2 Z Â U ref =u 0 Â Da=Re t and ρ 0 u 0 2 =l 2 Â Da=Re t according to the scaling arguments by Swaminathan and Bray (2005) and Mantel and Borghi (1994), respectively, which are different from the scalings of T 33 given by Eqs. (4a) Swaminathan and Bray (2005) and thus the T33M3 model underpredicts the magnitude of T 33 away from the wall for the thin reaction zones regime flames (i.e., Ka > 1) considered here. It can be seen from Figure 9 that the performance of the T33CS model remains comparable to that of T33M1 for high turbulent Reynolds number cases (i.e., cases C-E) when the flame is away from the wall. However, the T33CS model offers a more accurate prediction Figure 9. Variation of T 33 Â δ 2 Z =ρ 0 S 2 L with x 1 =δ Z along with the predictions by the T33MB, T33CS, T33M1, T33M2, and T33M3 models and Eq. (10f) (new model) at t = 4δ Z =S L ; 6δ Z =S L ; 8δ Z =S L ; and 10δ Z =S L for cases A-E. The x 1 -axis is shown in log scale for the inset. than the T3M1 model for cases A and B before the flame interacts with the wall. However, the T33CS model underpredicts the magnitude of T 33 at the final stage of quenching. This inadequacy is addressed here by the following modification: Þþ2 are the model parameters. The parameter A 8 and the involvement of expc ÀT À Á in C Ã 5 increase the magnitude of the model prediction in the near-wall region where the magnitude of T 33 is underpredicted by the T33CS model. The model parameter A 8 asymptotically approaches unity and C Ã 5 approaches C 5 away from the wall wherec ¼T. It can be seen from Figure 9 that the new model given by Eq. (10f) predicts the quantitative behavior of T 33 more satisfactorily than the other available models both away from and near to the wall.
Modeling of the combined reaction rate, dissipation, and diffusivity gradient terms The transport equation of N c ¼ DÑc Á Ñc can be rearranged in the following manner (Chakraborty et al., 2008a;Gao et al., 2015c): =ρ Ñc j j is the local flame displacement speed (Chakraborty et al., 2008a) andñ ¼ ÀÑc= Ñc j j is the local flame normal vector. Consequently, the combined contribution of the terms D 1 , T 4 , ÀD 2 ð Þ, and f D ð Þ can be written as (Chakraborty et al., 2008a: wherem ¼ ÀÑc= Ñc is the resolved flame normal, and S r ¼ _ ω= ρ Ñc j j ð Þ and S n 1 N Á Ñ ρDÑ Á Ñc À Á = ρ Ñc j j ð Þ are the reaction and normal diffusion components of the displacement speed, respectively (Echekki and Chen, 1999;Peters et al., 1998). It can be seen with Eq. (11b) that the net contribution of D 1 þ T 4 À D 2 þ f D ð Þ ½ specifies the effect due to flame normal propagation and flame curvature. Followed by previous analyses (Chakraborty et al., 2008aSwaminathan, 2010, 2013;Mantel and Borghi, 1994), it might be more convenient to model the net contribution of D 1 þ T 4 À D 2 þ f D ð Þ ½ rather than its individual components. Mantel and Borghi (1994) proposed the following model for D 1 þ T 4 À D 2 ½ : where β 1 ¼ 4:2 and C ¼ 0:1 are the model parameters. Since, D 1 is a close term, it is more convenient to model only T 4 À D 2 þ f D ð Þ ½ . Chakraborty et al. (2008a) proposed the following model for T 4 À D 2 ½ : where β 2 ¼ 6:7 is a model parameter. It is worth noting that f D ð Þ was ignored by Mantel and Borghi (1994) and Chakraborty et al. (2008a) and thus the models given by Eqs. (11c) and (11d) will not be discussed further in this article. Recently Gao et al. (2015c) extended the model given by Eq. (11d) in the following manner (i.e., T4D2CS model): where β V ¼ 6:7 is a model parameter. The predictions of the T4D2CS model are compared with DNS data in Figure 10, which shows the variation of T 4 À D 2 þ f D ð Þ ½ with normalized wall normal distance x 1 =δ Z at different time instances. The net contribution of remains negative when the flame is away from the wall (i.e., x 1 =δ Z ) Pe min ð Þ L ). However, a weakly positive contribution of T 4 À D 2 þ f D ð Þ ½ is observed in the region given by Pe min ð Þ L x 1 =δ Z 2 Pe min ð Þ L . The magnitude of negative contribution of T 4 À D 2 þ f D ð Þ ½ increases significantly in the region given by x 1 =δ Z Pe min ð Þ L during flame quenching. This increase in the sink contribution of T 4 À D 2 þ f D ð Þ ½ arises due to significant f D ð Þ contribution in the near-wall region (see Figure 4). It can be seen from Figure 10 that the T4D2CS model is able to capture T 4 À D 2 þ f D ð Þ ½ obtained from the DNS data satisfactorily when the flame is away from the wall (i.e., x 1 =δ Z ) Pe min ð Þ L ). However, the T4D2CS model does not adequately capture the qualitative and quantitative behaviors of T 4 À D 2 þ f D ð Þ ½ in the near-wall region. Here, the T4D2CS model has been modified in the following manner in order to improve the near-wall predictions: The corresponding model parameters are given as: Θ ¼ 0:5 1 À erf x 1 δ Z À 1:9 Pe min ð Þ The model parameter A 9 has been introduced in order to capture the augmentation of negative contribution of T 4 À D 2 þ f D ð Þ ½ in the near-wall region, whereas A 10 is responsible for changing the sign of the model prediction and Θ makes sure this change in sign becomes active at x 1 =δ z % 1:9 Pe min ð Þ L . However, Θ vanishes away from the wall (i.e., Figure 10. Variation of T 4 À D 2 þ f D ð Þ ½ Â δ 2 Z =ρ 0 S 2 L with x 1 =δ Z along with the predictions by the T4D2CS model and Eq. (11f) (new model) at t = 4δ Z =S L ; 6δ Z =S L ; 8δ Z =S L ; and 10δ Z =S L for cases A-E. The x 1 -axis is shown in log scale for the inset.
x 1 =δ Z ) Pe min ð Þ L ) and thus Eq. (11f) becomes identical to the T4D2CS model, i.e., Eq. (11e), which can be substantiated from Figure 10, where the predictions of Eq. (11f) are also shown. It can be seen from Figure 10 that the model given by Eq. (11f) adequately predicts the augmentation of the magnitude of negative contribution of T 4 À D 2 þ f D ð Þ ½ in the near-wall region and also captures a slightly positive value of T 4 À D 2 þ f D ð Þ ½ in the region given by Pe min ð Þ L x 1 =δ Z 2 Pe min ð Þ L . Thus, Eq. (11f) can be used for modeling T 4 À D 2 þ f D ð Þ ½ both close to and away from the wall.

Conclusions
The SDRε c transport and its modeling in the context of RANS have been analyzed for headon quenching of turbulent premixed flame by an inert isothermal wall based on 3D simple chemistry DNS data. It has been found that an increase in u 0 =S L leads to an increase in the magnitudes of the unclosed terms of the SDR transport equation. For all cases the terms arising from density variation, scalar-turbulence interaction, reaction rate gradient, molecular diffusivity gradient, and molecular dissipation, i.e., T 2 ; T 3 ; T 4 ; f D ð Þ, and ÀD 2 ð Þ, remain the leading order contributors to the SDRε c transport away from the wall, and the turbulent transport and molecular diffusion terms remain negligible in comparison to the aforementioned leading order terms. A rough equilibrium between (T 2 þ T 3 þ T 4 þ f D ð ÞÞ and ÀD 2 ð Þhas been observed away from the wall before the quenching is initiated, but this equilibrium is not maintained at the vicinity of the wall during flame quenching. The performances of previously proposed models for T 1 ; T 2 ; T 31 ; T 32 ; T 33 , and T 4 À D 2 þ f D ð Þ ð Þ have been assessed with respect to the corresponding quantities extracted from DNS data. It has been found that the aforementioned models do not adequately predict the near-wall behavior of the unclosed terms of the SDR transport equation. The models, which exhibit the most promising performance away from the wall, have been modified to account for the near-wall behavior in such a manner that they asymptotically approach the existing model expressions away from the wall. Although the functional form of the modeling parameters have been proposed so that they follow the asymptotic behavior in terms of Da L ; Ka L , and x 1 =δ Z , it is likely that they will need to be modified when datasets with a larger range of Re t with detailed chemistry will be explored. It is worth noting that several previous DNS analyses on flame-wall interaction used a simple chemical mechanism Rutland, 1998, 2002;Bruneaux et al., 1996Bruneaux et al., , 1997Poinsot et al., 1993) and the same approach has been followed here. Moreover, all existing SDR transport closures have been proposed based on simple chemistry DNS data Chakraborty and Swaminathan, 2007a, 2007bKolla et al., 2009;Mura et al., 2008Mura et al., , 2009Swaminathan and Bray, 2005). Since statistical behaviors of Ñc j j have been adequately captured by single-step chemistry [see Chakraborty and Klein (2008) and Chakraborty et al. (2008b for scalar gradient statistics based on both simple and detailed chemistry DNS data], it can be expected that the findings will at least be qualitatively valid in the context of detailed chemistry and transport. However, a different wall boundary condition (e.g., adiabatic wall boundary condition) may lead to the modification of some of the wall corrections proposed here, but this analysis is beyond the scope of this article. Moreover, Lewis number may have some influence on the modeling of SDR transport but the qualitative nature of the present findings is unlikely to be modified . This necessitates comprehensive experimental and DNS-based validations at high values of Re t for accurate estimation of the model parameters. Furthermore, the proposed models need to be implemented in actual RANS simulations to assess their predictive capabilities. Some of the aforementioned issues will form the foundation of future analyses.