A Priori Direct Numerical Simulation Modelling of the Curvature Term of the Flame Surface Density Transport Equation for Nonunity Lewis Number Flames in the Context of Large Eddy Simulations

Abstract

A Direct Numerical Simulation (DNS) database of freely propagating statistically planar turbulent premixed flames with Lewis numbers Le ranging from 0.34 to 1.2 has been used to analyse the statistical behaviours of the curvature term of the generalised Flame surface Density (FSD) transport equation, in the context of the Large Eddy Simulation (LES). Lewis number is shown to have significant influences on the statistical behaviours of the resolved and sub-grid parts of the FSD curvature term. It has been found that the existing models for the sub-grid curvature term 𝐶_sg do not capture the qualitative behaviour of this term extracted from the DNS database for flames with Le<<1. The existing models of 𝐶_sg only predict negative values, whereas the sub-grid curvature term is shown to assume positive values within the flame brush for the Le=0.34 and 0.6 flames. Here the sub-grid curvature terms arising from combined reaction and normal diffusion and tangential diffusion components of displacement speed are individually modelled, and the new model of the sub-grid curvature term has been found to capture 𝐶_sg extracted from DNS data satisfactorily for all the different Lewis number flames considered here for a wide range of filter widths.

1. Introduction

Flame Surface Density (FSD) based reaction rate closure is well established in the context of Reynolds Averaged Navier-Stokes (RANS) simulations of turbulent premixed flames [1, 2]. The increased affordability of high performance computing has made Large Eddy Simulation (LES) an alternative simulation tool, where the large-scale physical processes are resolved, but modelling is still required for the subgrid quantities. The FSD-based reaction rate closure has recently been successfully extended for the purpose of LES [3–14]. In LES simulation of premixed combustion, a Favre-filtered reaction progress variable transport equation is solved alongside other filtered conservation equations. The reaction progress variable is defined as 𝑐=(𝑌_𝑅0−𝑌_𝑅)/(𝑌_𝑅0−𝑌_𝑅∞), where 𝑌_𝑅 is the mass fraction of a suitable reactant and the subscripts 0 and ∞ denote the values in the fully unburned and burned gases, respectively. The generalised FSD is defined as Σ_gen=|∇𝑐| [3–14], where the overbar indicates the LES filtering process. The Favre-filtered reaction progress variable transport equation takes the following form: 𝜌𝜕𝑐̃+𝜕𝑡𝜌𝑢_𝑗𝜕𝑐̃𝜕𝑥_𝑗=𝜕𝜕𝑥_𝑗𝜌𝐷𝜕𝑐𝜕𝑥_𝑗+̇𝜕𝑤−𝜕𝑥_𝑗𝜌𝑢_𝑗𝑐−𝑢_𝑗𝑐̃, (1) where 𝑄=𝜌𝑄/𝜌indicates the Favre filtered value of a general variable Q, 𝑢_𝑗 is the velocity component in the jth direction, ρ is the density, D is the molecular diffusivity, and 𝑤̇ is the filtered reaction rate. The first two terms on right hand side of (1) denote the filtered molecular diffusion and reaction rates, respectively, and their combined contribution can be modelled using Σ_gen in the following manner: 𝑤̇+∇⋅𝜌𝐷∇𝑐=𝜌𝑆_𝑑𝑠Σ_gen, (2) where (𝑄)_𝑠=𝑄|∇𝑐|/Σ_gen indicates the surface-weighted filtered value of a general quantity Q and 𝑆_𝑑=𝐷𝑐/𝐷𝑡/|∇𝑐| is the displacement speed, which denotes the speed at which a given c isosurface moves normal to itself with respect to an initially coincident material surface. The generalised FSD Σ_gen is an unclosed quantity and is closed either by using an algebraic expression or by solving a modelled transport equation alongside other conservation equations. The algebraic closure is valid when the generation rate of flame surface area remains in equilibrium with its destruction rate, but this assumption is rendered invalid under unsteady conditions (e.g., combustion instabilities). Under unsteady conditions, it is often advantageous to solve a modelled transport equation ofΣ_gen. The exact transport equation for the generalised FSD Σ_gen is given as [1, 4–7, 9, 10, 12]: 𝜕Σ_gen+𝜕𝑢𝜕𝑡^𝑗Σ_gen𝜕𝑥_𝑗𝜕=−𝑢_𝑖𝑠−𝑢_𝑖Σ_gen𝜕𝑥_𝑖+𝛿_𝑖𝑗−𝑁_𝑖𝑁_𝑗𝜕𝑢_𝑖𝜕𝑥_𝑗𝑠Σ_gen−𝜕𝑆_𝑑𝑁_𝑖𝑠Σ_gen𝜕𝑥_𝑖+𝑆_𝑑𝜕𝑁_𝑖𝜕𝑥_𝑖𝑠Σ_gen, (3) where 𝑁_𝑖=−(𝜕𝑐/𝜕𝑥_𝑖)/|∇𝑐| is the ith component of flame normal vector. The terms on the left hand side of (3) denote transient and mean advection effects, respectively. The first three terms on the right hand side of (3) denote the effects of subgrid convection, flame surface area generation due to fluid-dynamic straining, and flame normal propagation, respectively. The last term of (3) describes the production/destruction of Σ_gen due to flame curvature 𝜅_𝑚=(𝜕𝑁_𝑖/𝜕𝑥_𝑖)/2 and thus referred to as the FSD curvature term [4–7, 9, 10, 12]. It has been found in several previous studies [5–7, 9, 14] that the FSD curvature term remains a leading order contributor to the FSD transport for both unity and nonunity Lewis number turbulent premixed combustion. As the curvature term remains a leading order contributor to the FSD transport, the modelling of (𝑆_𝑑∇⋅𝑁^→)_𝑠Σ_gen is crucial for the transport equation-based FSD closure. The statistical behaviour of (𝑆_𝑑∇⋅𝑁^→)_𝑠Σ_gen is significantly affected by curvature dependence on 𝑆_𝑑 [9, 10, 12]. Earlier a priori Direct Numerical Simulation (DNS) analyses [9, 10, 12] showed that existing models for the subgrid curvature term 𝐶_sg do not adequately capture the qualitative behaviour of this term obtained from DNS data. Moreover, the model parameters for the existing subgrid curvature term 𝐶_sg models are found to be strong functions of the LES filter width Δ [9, 10, 12].

To date, most existing FSD-based models have been proposed for unity Lewis number flames where the differential diffusion of heat and mass has been ignored. The Lewis number is defined as the ratio of thermal diffusivity 𝛼_𝑇 to mass diffusivity 𝐷 (i.e., Le=𝛼_𝑇/𝐷). The effects of Le on the statistical behaviour of the FSD curvature term (𝑆_𝑑∇⋅𝑁^→)_𝑠Σ_gen are yet to be analysed in detail, and this paper aims to bridge this gap in the existing literature. It is worth noting that, in a premixed flame, different species have different values of Lewis number. Thus, specifying a global Lewis number Le characterising the whole combustion process is not straightforward. The Lewis number of the deficient reactant is often considered to be the characteristic Le of the combustion process in question [15, 16]. Moreover, several previous studies [16–29] analysed the effects of differential diffusion of heat and mass by modifying the characteristic Lewis number in isolation, and the same procedure has been adopted here. In the present study, a simplified chemistry-based DNS database of statistically planar turbulent premixed flames with global Lewis numbers ranging from 0.34 to 1.2 has been considered to analyse the statistical behaviour of the FSD curvature term (𝑆_𝑑∇⋅𝑁^→)_𝑠Σ_gen in the context of LES. In this context, the main objectives of this study are as follows: (1)to analyse the statistical behaviours of the subgrid FSD curvature term in the context of LES, for flames with different values of Lewis number; (2)to propose models for different components of the subgrid FSD curvature terms and assess their performances in comparison to the corresponding quantities extracted from DNS data.

The rest of the paper will be organised as follows. The necessary mathematical background will be provided in the next section. This will be followed by a brief description of the numerical implementation related to the DNS database. Following this, results will be presented and subsequently discussed. The main findings will be summarised, and conclusions will be drawn in the final section.

2. Mathematical Background

The curvature term of the FSD (𝑆_𝑑∇⋅𝑁^→)_𝑠Σ_gen is often decomposed in the following manner [4–7, 9, 10, 12]: 𝑆_𝑑∇⋅𝑁^→_𝑠Σ_gen=𝐶_mean+𝐶_sg, (4) where 𝐶_mean and 𝐶_sg are the resolved and subgrid components of the FSD curvature term, respectively. The resolved curvature term 𝐶_meancan be expressed in three different manners [5, 9, 10, 12]: 𝐶_mean=𝑆_𝑑𝑠𝜕𝑁_𝑖𝑠𝜕𝑥_𝑖Σ_gen𝐶_mean=𝑆_𝑑𝑠𝜕𝑀_𝑖𝜕𝑥_𝑖Σ_gen,𝐶_mean=𝛿_𝑖𝑗−𝑛_𝑖𝑗𝜕𝑆_𝑑𝑠𝑁_𝑖𝑠𝜕𝑥_𝑗Σ_gen, (5a)(5b)(5c) where 𝑀_𝑖=−(𝜕𝑐/𝜕𝑥_𝑖)/|∇𝑐| is the ith component of the resolved flame normal vector. It was demonstrated by Chakraborty and Cant [10, 12] that (5a) provides the best option for the resolved curvature term 𝐶_mean, as it gives rise to the smallest magnitude of 𝐶_sg among all the possibilities shown in (5a)–(5c). Equation (5a) was found to perform the best among the three possibilities shown in (5a), (5b), (5c) for this database. This is advantageous from the perspective of efficient modelling of the FSD curvature term (𝑆_𝑑∇⋅𝑁^→)_𝑠Σ_gen as most of the modelling uncertainty is associated with 𝐶_sg. Moreover, (5a) has also been used for the modelling of (𝑆_𝑑∇⋅𝑁^→)_𝑠Σ_gen in previous LES simulations [5–7, 13]. For the present analysis (5a), (i.e., 𝐶_mean=(𝑆_𝑑)_𝑠[𝜕(𝑁_𝑖)_𝑠/𝜕𝑥_𝑖]Σ_gen) will be considered for the resolved curvature term 𝐶_mean.

It is often useful to decompose the flame displacement speed 𝑆_𝑑=(𝐷𝑐/𝐷𝑡)/|∇𝑐|=[𝑤̇+∇⋅(𝜌𝐷∇𝑐)]/𝜌|∇𝑐| in the following manner for the purpose of modelling the FSD curvature term [9–12, 30, 31]: 𝑆_𝑑=𝑆_𝑟+𝑆_𝑛+𝑆_𝑡,𝑆_𝑟=𝑤̇𝜌||,∇𝑐𝑆_𝑛=^→𝑁⋅∇𝜌𝐷𝑁^→⋅∇𝑐𝜌||,∇𝑐𝑆_𝑡=−𝐷∇⋅𝑁^→=−2𝐷𝜅_𝑚, (6a)(6b)(6c)(6d) where 𝑆_𝑟 and 𝑆_𝑛 are the reaction and normal diffusion components of displacement speed and 𝑆_𝑡 is the tangential diffusion component of displacement speed. The following expression for 𝐶_sg can be obtained using (6a)–(6d) and (5a) (i.e., 𝐶_mean=(𝑆_𝑑)_𝑠[𝜕(𝑁_𝑖)_𝑠/𝜕𝑥_𝑖]Σ_gen): 𝐶_sg=𝐶_sg1+𝐶_sg2=𝑆_𝑑𝜕𝑁_𝑖𝜕𝑥_𝑖𝑠Σ_gen−𝑆_𝑑𝑠𝜕𝑁_𝑖𝑠𝜕𝑥_𝑖Σ_gen, (7) where 𝐶_sg1=𝑆_𝑟+𝑆_𝑛𝜕𝑁_𝑖𝜕𝑥_𝑖𝑠Σ_gen−𝑆_𝑟+𝑆_𝑛𝑠𝜕𝑁_𝑖𝑠𝜕𝑥_𝑖Σ_gen,𝐶_sg2=−𝐷𝜕𝑁_𝑖𝜕𝑥_𝑖²_𝑠Σ_gen−𝐷𝜕𝑁_𝑖𝜕𝑥_𝑖𝑠𝜕𝑁_𝑖𝑠𝜕𝑥_𝑖Σ_gen. (8a)(8b) Equation (8a) indicates that curvature (𝜅_𝑚=∇⋅𝑁^→/2) dependences of (𝑆_𝑟+𝑆_𝑛) and |∇𝑐| significantly influence the statistical behaviour of 𝐶_sg1. Equation (8b) suggests that 𝐶_sg2 is expected to assume negative values throughout the flame brush.

Hawkes and Cant [6, 7] modified a version of the Coherent Flamelet Model (CFM) by Candel et al. [2] for the purpose of LES as: 𝐶_sg𝛼=−^𝑁𝛽₁𝑆_𝐿Σ²_gen(1−𝑐), (9) where 𝛼_𝑁=1−(𝑁_𝑘)_𝑠(𝑁_𝑘)_𝑠 is a resolution parameter which vanishes when the flow is fully resolved and 𝛽₁ is a model parameter. Hawkes [5] discussed a possibility of modifying a RANS model proposed by Cant et al. [1] for the purpose of LES as: 𝐶_sg𝐶=−^𝐻𝑆_𝐿Σ²_gen(1−𝑐), (10) where 𝐶_𝐻=𝛼_𝑁𝛽₂(1−(1/3)[1−exp(−10(1−𝑐)𝑘/Σ_gen𝑆_𝐿Δ)], 𝐴=10.0, 𝑢^′_Δ=2𝑘/3 is the subgrid turbulent velocity fluctuation, 𝑘=(𝜌𝑢_𝑖𝑢_𝑖−𝜌𝑢_𝑖𝑢_𝑖)/2𝜌 is the subgrid kinetic energy, and 𝛽₂ is a model parameter. Another model of 𝐶_sg was proposed by Charlette et al. [4]: 𝐶_sg𝛽=−³𝑆_𝐿Σ_gen−|∇𝑐|Σ_gen𝑐(1−𝑐), (11) where 𝛽₃ is a model parameter. The models given by (9)–(11) (henceforth will be referred to as CSGCFM, CSGCPB, and CSGCHAR, resp.) ensure that 𝐶_sg vanishes when the flow is fully resolved (i.e., (𝑁_𝑘)_𝑠(𝑁_𝑘)_𝑠=1.0 and Σ_gen=|∇𝑐|). A priori DNS assessment of the CSGCFM, CSGCPB, and CSGCHAR models and the modelling of 𝐶_sg1 and 𝐶_sg2 will be addressed in Section 4 of this paper.

3. Numerical Implementation

In principle combustion, DNS should account for both three dimensionality of turbulence and detailed chemical mechanism. However, until recently, most combustion DNS studies were carried out either in two dimensions with detailed chemistry or in three dimensions with simplified chemistry due to the limitation of computer storage capacity. Although it is now possible to carry out three-dimensional DNS with detailed chemistry, they remain extremely expensive (e.g., millions of CPU hours and thousands of processors [32]) and the cost of an extensive parametric analysis based on three-dimensional detailed chemistry-based DNS often becomes prohibitive. As the present analysis concentrates on an extensive parametric variation in terms of Lewis number, the chemical mechanism is simplified here by an Arrhenius-type irreversible single-step chemical reaction (i.e., Reactants→Products) following several previous studies [1–12, 14]. It has been found that the strain rate and curvature dependences of 𝑆_𝑑 and |∇𝑐| obtained from three-dimensional simplified chemistry DNS [25–27, 33, 34] are found to be qualitatively similar to the corresponding behaviours obtained from detailed chemistry-based DNS simulations [16, 30, 31, 35]. As the statistical behaviours of the FSD curvature term (𝑆_𝑑∇⋅𝑁^→)_𝑠Σ_gen are strongly dependent on the curvature dependences of 𝑆_𝑑 and |∇𝑐|, the results for this analysis are expected to be valid even for detailed chemistry based simulations at least in a qualitative sense without much loss of generality. Several studies [3–7, 9–12] have concentrated on a priori DNS modelling of FSD based on simplified chemistry in the past and the same approach has been adopted here.

A compressible three-dimensional DNS code SENGA [36] was used for the simulations where the conservation equations of mass, momentum, energy, and species are solved in nondimensional form. A cubic domain of each side equal to 24.1𝛿_th is considered for the present DNS database where 𝛿_th is the thermal flame thickness, which is defined as 𝛿_th=(𝑇_ad−𝑇₀)/Max|∇𝑇|_𝐿, and the subscript 𝐿 refers to quantities in an unstrained planar laminar flame with 𝑇_ad, 𝑇₀, and 𝑇 being the adiabatic flame, unburned gas, and instantaneous gas temperatures, respectively. The computational domain was discretised using a Cartesian grid of 230×230×230 with equal grid spacing in each direction. The grid spacing Δ𝑥 is determined based on the flame resolution, and about 10 grid points are kept within the thermal flame thickness 𝛿_th for all the cases considered here. This grid spacing Δ𝑥 corresponds to 0.73𝜂, where 𝜂 is the Kolmogorov length scale. The boundaries in the mean flame propagation were taken to be partially nonreflecting and were implemented using the Navier-Stokes Characteristic Boundary Conditions (NSCBC) technique [37]. The boundary conditions in the transverse direction were taken to be periodic. The spatial derivatives for the internal grid points were evaluated using a tenth-order central differencing scheme, and the order of differentiation gradually decreases to a one-sided 2nd order scheme at the partially nonreflecting boundaries. The time advancement was carried out using an explicit low storage third-order Runge-Kutta scheme [38].

For the current DNS database, the turbulent velocity field was initialised using a pseudospectral method [39] following the Batchelor-Townsend turbulent kinetic energy spectrum [40]. The flame is initialised using a steady planar unstrained laminar flame solution. The initial values of normalised root mean square (rms) turbulent velocity fluctuation 𝑢′/𝑆_𝐿, integral length scale to thermal flame thickness ratio 𝑙/𝛿_th, heat release parameter 𝜏=(𝑇_ad−𝑇₀)/𝑇₀, Damköhler number Da=𝑙𝑆_𝐿/𝑢′𝛿_th, and Karlovitz number Ka=(𝑢′/𝑆_𝐿)^3/2(𝑙/𝛿_th)^−1/2 are listed in Table 1. According to Peters [41], all the cases considered here can be taken to represent the thin reaction zone regime combustion, as Ka remains greater that unity. Standard values are considered for Prandtl number Pr and the Zel’dovich number 𝛽 (i.e., Pr=0.7 and 𝛽=𝑇_ac(𝑇_ad−𝑇₀)𝑇²_ad=6.0, where 𝑇_ac is the activation temperature).

Table 1 

Initial values of the simulation parameters and non-dimensional numbers relevant to the DNS database.

Under decaying turbulence, DNS simulations should be carried out for a simulation time 𝑡_sim≥max(𝑡_𝑓,𝑡_𝑐) [42], where 𝑡_𝑓=𝑙/𝑢^′ is the initial eddy turn over time and 𝑡_𝑐=𝛿_th/𝑆_𝐿 is the chemical time scale. For this database, the statistics were extracted after about three eddy turn over times (i.e., 3𝑡_𝑓=3𝑙/𝑢^′), which corresponded to one chemical time scale (i.e., 𝑡_𝑐=𝛿_th/𝑆_𝐿). This simulation time remains small but comparable to several studies [3, 24, 28, 43–47] which contributed significantly to the FSD-based modelling in the past. The statistics presented in this paper did not change significantly since halfway through the simulation (i.e., 1.5𝑡_𝑓=1.5𝑙/𝑢^′). The value of 𝑢^′/𝑆_𝐿 in the unburned gas ahead of the flame had decayed by 50% of its initial value when the statistics were extracted. By this time, the normalised integral value 𝑙/𝛿_th had increased to around 1.7 times of its initial value. The values of 𝑢^′/𝑆_𝐿 and 𝑙/𝛿_th at the time statistics were extracted are also representative of the thin reaction zones regime combustion [41]. This DNS database was used extensively earlier for the purpose of RANS modelling [27, 28, 48, 49], and the interested readers are referred to these papers for further details.

The DNS data was explicitly LES filtered using a Gaussian filter kernel in physical space for the purpose of a priori analysis. The filtered quantity 𝑄(𝑥^→,𝑡), is given by 𝑄𝑥^→,𝑡=𝑄𝑥^→−^→𝑟𝐺^→𝑟𝑑^→𝑟, (12) where 𝐺(^→𝑟) is the Gaussian filter kernel, which is defined in the following manner: 𝐺^→𝑟=6𝜋Δ^23/2exp−6^→𝑟⋅^→𝑟Δ². (13)

The filtered quantities of interest were extracted for filter widths Δ ranging from 0.4𝛿_th to 2.4𝛿_th in steps of 0.4𝛿_th. These filter sizes are comparable to the range of Δ used in several previous studies [3, 4, 9–12, 14] for a priori DNS analysis and span a useful range of length scales (i.e., from Δ≈0.4𝛿_th, where the flame is partially resolved, up to 2.4𝛿_th, where the flame becomes fully unresolved and Δ is comparable to the integral length scale).

4. Results and Discussion

The instantaneous isosurfaces of 𝑐 ranging from 0.01 to 0.99 at 𝑡_𝑐=𝛿_th/𝑆_𝐿 are shown in Figure 1, which indicates that the flame wrinkling increases with decreasing Lewis number and this tendency is particularly prevalent for the Le≪1 flames due to thermodiffusive instabilities [17–29]. The unburned reactants diffuse into the reaction zone at a faster rate than the rate at which heat diffuses out in the Le<1 flames. This gives rise to simultaneous presence of high temperature and reactant concentration in the reaction zone for the Le<1 flames, which in turn leads to greater burning rate and flame surface area generation in comparison to the unity Lewis number flame. By contrast, heat diffuses faster than the diffusion rate of reactants into the reaction zone in the case of Le>1, which reduces the burning rate and the rate of flame area generation in comparison to the unity Lewis number flame. The increase in burning rate and flame area generation with decreasing Lewis number can be substantiated by the values of normalised turbulent flame speed 𝑆_𝑇/𝑆_𝐿 and normalised flame surface area 𝐴_𝑇/𝐴_𝐿 which are presented in Table 2. The values of 𝑆_𝑇/𝑆_𝐿 have been evaluated by volume integrating the reaction rate 𝑤̇ using the expression 𝑆_𝑇=(1/𝜌₀𝐴_𝑃)∫_𝑉𝑤̇𝑑𝑉, where 𝐴_𝑃 is the projected area of the flame in the direction of mean flame propagation, while the values of 𝐴_𝑇/𝐴_𝐿 have been evaluated by volume integrating |∇𝑐| (i.e., ∫_𝑉|∇𝑐|𝑑𝑉) under both turbulent and laminar conditions. Table 2 shows that both 𝑆_𝑇/𝑆_𝐿 and 𝐴_𝑇/𝐴_𝐿 increase with decreasing Lewis number, and this effect is particularly prevalent in the flames with Le < 1 due to the presence of thermodiffusive instabilities [17–29]. The increase in flame wrinkling with decreasing Lewis number is also visually evident from the c isosurfaces presented in Figure 1.

Table 2 

The effects of Lewis number on normalised turbulent flame speed 𝑆_𝑇/𝑆_𝐿 and normalised flame surface area 𝐴_𝑇/𝐴_𝐿 after 3.0 initial eddy turn over times.

(a)

(b)

(c)

(d)

(e)

Figure 1 

Instantaneous isosurfaces of 𝑐 ranging from 0.01 to 0.99 at 𝑡=3𝑡𝑓=𝑡𝑐 for cases (a) Le=0.34; (b) Le=0.6; (c) Le=0.8; (d) Le=1.0; (e) Le=1.2. The domain size is 24.1𝛿th×24.1𝛿th×24.1𝛿th.

The variations of (𝑆_𝑑∇⋅𝑁^→)_𝑠Σ_gen, 𝐶_mean, and 𝐶_sg conditionally averaged in bins of 𝑐̃ isosurfaces for cases (a)–(e) are shown in Figure 2 for filter widths Δ=8Δ_𝑚≈0.8𝛿_th and Δ=24Δ_𝑚≈2.4𝛿_th, where Δ_m is the DNS grid size. It is evident from Figure 2 that Le significantly affects the statistical behaviours of the curvature terms. The filter widths Δ=8Δ_𝑚≈0.8𝛿_th, and Δ=24Δ_𝑚≈2.4𝛿_th span a useful range of length scales (i.e., from Δ≈0.8𝛿_th, where the flame is partially resolved, up to 2.4𝛿_th where the flame becomes fully unresolved and Δ is comparable to the integral length scale). In the Le ≪ 1 flames (e.g., cases (a) and (b)), the FSD curvature term (𝑆_𝑑∇⋅𝑁^→)_𝑠Σ_gen behaves as a source term for the major part of the flame brush before assuming negative values towards the burned gas side for Δ=8Δ_𝑚≈0.8𝛿_th. For Δ=24Δ_𝑚≈2.4𝛿_th, the FSD curvature term (𝑆_𝑑∇⋅𝑁^→)_𝑠Σ_gen acts as a source (sink) term towards the unburned (burned) gas side of the flame brush in the Le≪1 flames. In the case of Le≈1.0 flames (i.e., cases (c)–(e)) the curvature term (𝑆_𝑑∇⋅𝑁^→)_𝑠Σ_gen behaves as a sink type term throughout the flame brush for all filter widths. It can be seen from Figure 2 that 𝐶_mean acts as a source (sink) term for cases (a)-(b) ((c)–(e)). The magnitude of 𝐶_mean (𝐶_sg) decreases (increases) with increasing Δ in all cases, and for large filter widths (𝑆_𝑑∇⋅𝑁^→)_𝑠Σ_gen is principally made up of 𝐶_sg. The LES filtering is a convolution process, and the weighted averaging involved in the filtering process leads to a decrease in the magnitude of 𝐶_mean with increasing filter width Δ. The flow becomes increasingly unresolved with increasing filter width Δ, and this is reflected in the rise in 𝐶_sg magnitude with increasing filter width Δ.

(a)

(b)

(c)

(d)

(e)

(f)

(g)

(h)

(i)

(j)

Figure 2 

Variation of (𝑆𝑑∇⋅𝑁→)𝑠Σgen (solid line), 𝐶mean=(𝑆𝑑)𝑠𝜕(𝑁𝑖)s/𝜕𝑥𝑖Σgen (dashed line) and 𝐶sg (asterisk line) conditionally averaged in bins of 𝑐̃ across the flame brush for filter sizes Δ=8Δ𝑚≈0.8𝛿th (top row): (a) Le=0.34 (1st column); (b) Le=0.6 (2nd column); (c) Le=0.8 (3rd column); (d) Le=1.0 (4th column); (e) Le=1.2 (5th column) and for filter size Δ=24Δ𝑚≈2.4𝛿th (bottom row): (f) Le=0.34 (1st column); (g) Le=0.6 (2nd column); (h) Le=0.8 (3rd column); (i) Le=1.0 (4th column); (j) Le=1.2 (5th column). All terms are ensemble averaged on 𝑐̃ isosurfaces in Figure 2 and subsequent cases. All the curvature terms in this and subsequent figures are normalised by 𝑆𝐿/𝛿2th.

The resolved curvature term 𝐶_mean=(𝑆_𝑑)_𝑠𝜕(𝑁_𝑖)_𝑠/𝜕𝑥_𝑖Σ_gen can be seen to capture the behaviour of the curvature term (𝑆_𝑑∇⋅𝑁^→)_𝑠Σ_gen, well at small filter widths (i.e., Δ≤𝛿_th) for flames with Le≈1.0 (i.e., cases (c)–(e)). However, the magnitude of 𝐶_mean decreases with increasing Δ and it does not capture the behaviour of the FSD curvature term (𝑆_𝑑∇⋅𝑁^→)_𝑠Σ_gen for the Le ≪ 1.0 flames (i.e., cases (a) and (b)). The subgrid curvature term, C_sg follows the qualitative behaviour of the FSD curvature term (𝑆_𝑑∇⋅𝑁^→)_𝑠Σ_gen for all filter widths. The subgrid curvature term 𝐶_sg almost entirely makes up the FSD curvature term (𝑆_𝑑∇⋅𝑁^→)_𝑠Σ_gen for Δ≫𝛿_th, and this is especially true for the Le ≪ 1.0 cases (i.e., cases (a) and (b)). It can further be observed from Figure 2 that 𝐶_sg assumes positive values towards the unburned gas side of the flame brush in the Le ≪ 1 flames (e.g., cases (a) and (b)), whereas the existing models for 𝐶_sg allow for only negative values (see (9)–(11)). This suggests that new models for 𝐶_sg are warranted to account for the influences of nonunity Lewis number (i.e., Le≠1.0).

In order to be able to model the subgrid curvature term 𝐶_sg, the decomposition prescribed in (8a)-(8b) has been used here. The variations of ((𝑆_𝑟+𝑆_𝑛)∇⋅𝑁^→)_𝑠Σ_gen, −(𝐷(∇⋅𝑁^→)²)_𝑠Σ_gen=−4(𝐷𝜅²_𝑚)_𝑠Σ_gen, 𝐶_sg1 and 𝐶_sg2 conditionally averaged in bins of 𝑐̃ isosurfaces for cases (a)–(e) are shown in Figure 3 for filter widths Δ=8Δ_𝑚≈0.8𝛿_th and Δ=24Δ_𝑚≈2.4𝛿_th. It is evident from Figure 3 that 𝐶_sg2 remains negative throughout the flame brush for all cases and follows the qualitative behaviour of (−4(𝐷𝜅²_𝑚)_𝑠Σ_gen). A comparison between ((𝑆_𝑟+𝑆_𝑛)∇⋅𝑁^→)_𝑠Σ_gen and −4(𝐷𝜅²_𝑚)_𝑠Σ_gen reveals that −4(𝐷𝜅²_𝑚)_𝑠Σ_gen remains the major contributor to (𝑆_𝑑∇⋅𝑁^→)_𝑠Σ_gen for all the flames at all values of Δ, which is consistent with the expected behaviour in the thin reaction zones regime [41]. The contribution of ((𝑆_𝑟+𝑆_𝑛)∇⋅𝑁^→)_𝑠Σ_gen remains significant for the Le<1 cases (i.e., cases (a), (b) and (c)), but its contribution remains weak in comparison to the magnitude of −4(𝐷𝜅²_𝑚)_𝑠Σ_gen in the Le=1.0 and 1.2 flames (i.e., cases (d) and (e)). Figure 3 demonstrates that 𝐶_sg1 remains close to the magnitude of ((𝑆_𝑟+𝑆_𝑛)∇⋅𝑁^→)_𝑠Σ_gen for all Δ for the Le=1.0 flame (i.e., case (d)), indicating that (𝑆_𝑟+𝑆_𝑛)_𝑠𝜕(𝑁_𝑖)_𝑠/𝜕𝑥_𝑖Σ_gen does not play a major role in capturing the behaviour of ((𝑆_𝑟+𝑆_𝑛)∇⋅𝑁^→)_𝑠Σ_gen. However, there is a significant difference in magnitudes of ((𝑆_𝑟+𝑆_𝑛)∇⋅𝑁^→)_𝑠Σ_gen and C_sg1 for small values of Δ (i.e., Δ<𝛿_th) in the nonunity Lewis number flames (i.e., cases (a)–(c) and (e)), which indicates that (𝑆_𝑟+𝑆_𝑛)_𝑠𝜕(𝑁_𝑖)_𝑠/𝜕𝑥_𝑖Σ_gen plays a key role for small values of filter width in these flames. For large values of filter width (i.e., Δ≫𝛿_th) 𝐶_sg1 remains the major contributor to ((𝑆_𝑟+𝑆_𝑛)∇⋅𝑁^→)_𝑠Σ_gen for all cases considered here, indicating that (𝑆_𝑟+𝑆_𝑛)_𝑠𝜕(𝑁_𝑖)_𝑠/𝜕𝑥_𝑖Σ_gen plays progressively less important role for increasing values of Δ.

(a)

(b)

(c)

(d)

(e)

(f)

(g)

(h)

(i)

(j)

Figure 3 

Variation of ((𝑆𝑟+𝑆𝑛)∇⋅𝑁→)𝑠Σgen (solidline), −4(𝐷𝜅2𝑚)𝑠Σgen (dashed line), 𝐶sg1 (asterisk line), and 𝐶sg2 (squared line) conditionally averaged in bins of 𝑐̃ across the flame brush for Δ=8Δm≈0.8𝛿th (top row): (a) Le=0.34 (1st column); (b) Le=0.6 (2nd column); (c) Le=0.8 (3rd column); (d) Le=1.0 (4th column); (e) Le=1.2 (5th column) and for filter size Δ=24Δ𝑚≈2.4𝛿th (bottom row): (f) Le=0.34 (1st column); (g) Le=0.6 (2nd column); (h) Le=0.8 (3rd column); (i) Le=1.0 (4th column); (j) Le=1.2 (5th column).

Figure 3 shows that there is a significant difference between −4(𝐷𝜅²_𝑚)_𝑠Σ_gen and 𝐶_sg2 for all cases for small values of Δ, and the difference between these quantities decreases with increasing Δ. As most of the contribution of −4(𝐷𝜅²_𝑚)_𝑠Σ_gen remains unresolved for large values of Δ, the subgrid curvature term 𝐶_sg2 remains the major contributor to −4(𝐷𝜅²_𝑚)_𝑠Σ_gen, indicating that (−(𝐷𝜕𝑁_𝑖/𝜕𝑥_𝑖)_𝑠𝜕(𝑁_𝑖)_𝑠/𝜕𝑥_𝑖Σ_gen) plays progressively less important role for increasing values of Δ where the flame is fully unresolved. However, the contribution of (−(𝐷𝜕𝑁_𝑖/𝜕𝑥_𝑖)_𝑠𝜕(𝑁_𝑖)_𝑠/𝜕𝑥_𝑖Σ_gen) remains significant for small values of Δ, where the flame is partially resolved. Figure 3 further shows that the order of magnitudes of both 𝐶_sg1 and 𝐶_sg2 remains comparable and thus accurate modelling of 𝐶_sg1 and 𝐶_sg2 is necessary for precise predictions of 𝐶_sg.

As the range of 𝜅_𝑚 values obtained on a flame surface increases with increasing flame wrinkling, the magnitude of −4(𝐷𝜅²_𝑚)_𝑠 increases with decreasing Le, which in turn leads to increasing magnitude of −4(𝐷𝜅²_𝑚)_𝑠Σ_gen and 𝐶_sg2 (see Figure 3). The positive contribution of 𝐶_sg1 overcomes the negative contribution of 𝐶_sg2 towards the unburned gas side of the flame brush for the Le=0.34 and 0.6 flames (i.e., cases (a) and (b)) and yields a net positive contribution of 𝐶_sg towards the reactant side of the flame brush (see Figure 2).

The statistical behaviours of ((𝑆_𝑟+𝑆_𝑛)∇⋅𝑁^→)_𝑠Σ_gen and 𝐶_sg1 depend on the nature of the curvature 𝜅_𝑚=∇⋅𝑁^→/2 dependences of (𝑆_𝑟+𝑆_𝑛) and |∇𝑐|, and the variation of (𝜅_𝑚)_𝑠 across the flame brush. The correlation coefficients for (𝑆_𝑟+𝑆_𝑛)−𝜅_𝑚 and |∇𝑐|−𝜅_𝑚 for five different 𝑐 isosurfaces across the flame brush for all the cases are shown in Figures 4(a) and 4(b), respectively. For all cases, 𝑆_𝑡=−2𝐷𝜅_𝑚 remains negatively correlated with 𝜅_𝑚 with a correlation coefficient close to (−1.0). However, Figures 4(a) and 4(b) demonstrate that Le significantly affects the curvature 𝜅_𝑚 dependences of (𝑆_𝑟+𝑆_𝑛) and |∇𝑐|. It can be seen from Figures 4(a) and 4(b) that (𝑆_𝑟+𝑆_𝑛) and |∇𝑐| remain positively (negatively) correlated with 𝜅_𝑚 for the Le<1.0 (Le>1.0) flames, whereas both (𝑆_𝑟+𝑆_𝑛) and |∇𝑐| show weak curvature dependences in the unity Lewis number flame. The positive correlations between (𝑆_𝑟+𝑆_𝑛) and 𝜅_𝑚, and between |∇𝑐| and 𝜅_𝑚 strengthen with decreasing Le for the Le<1 flames. The physical explanations for the observed influences of Lewis number on the curvature dependence of (𝑆_𝑟+𝑆_𝑛) and |∇𝑐| have been discussed elsewhere [25–27] and thus will not be discussed in this paper.

(a)

(b)

(c)

(d)

Figure 4 

Correlation coefficients between (a) (𝑆𝑟+𝑆𝑛) and 𝜅𝑚, and between (b) |∇c| and 𝜅𝑚 correlations on 𝑐= 0.1, 0.3, 0.5, 0.7, and 0.9 isosurfaces for cases A–E. Variation of (𝜅𝑚)𝑠×𝛿th with 𝑐̃ across the flame brush for (c) Δ=8Δ𝑚≈0.8𝛿th and (d) Δ=24Δ𝑚≈2.4𝛿th for cases A–E.

The variations of (𝜅_𝑚)_𝑠 conditionally averaged in bins of 𝑐̃ isosurfaces for cases A–E are shown in Figures 4(c) and 4(d) for filter widths Δ=8Δ_𝑚≈0.8𝛿_th and Δ=24Δ_𝑚≈2.4𝛿_th, respectively. It is evident from Figures 4(c) and 4(d) that (𝜅_𝑚)_𝑠 assumes positive (negative) values towards the unburned (burned) gas side of the flame brush. For small values of Δ, the surface-weighted filtered value of curvature (𝜅_𝑚)_𝑠 approaches to 𝜅_𝑚 (i.e., lim_Δ→0(𝜅_𝑚)_𝑠=𝜅_𝑚|∇𝑐|/|∇𝑐|=𝜅_𝑚) and thus the ensemble averaged value of (𝜅_𝑚)_𝑠 remains small for small values of filter width as the ensemble averaged value of 𝜅_𝑚 remains negligible for statistically planar flames. The difference between the ensemble averaged values of (𝜅_𝑚)_𝑠 and 𝜅_𝑚 increases with increasing filter width Δ, as flame wrinkling increasingly takes place at the subgrid level. For the Le=1.0 flame (i.e., case D), the combination of positive (negative) value of (𝜅_𝑚)_𝑠 and weak (𝑆_𝑟+𝑆_𝑛)−𝜅_𝑚 and |∇𝑐|−𝜅_𝑚 correlations gives rise to positive (negative) values of the ensemble averaged values of ((𝑆_𝑟+𝑆_𝑛)∇⋅𝑁^→)_𝑠Σ_gen and 𝐶_sg1 towards the unburned (burned) gas side of the flame brush for all values of Δ. The predominant positive (𝑆_𝑟+𝑆_𝑛)−𝜅_𝑚 and |∇𝑐|−𝜅_𝑚 correlations give rise to positive values of the ensemble averaged values of ((𝑆_𝑟+𝑆_𝑛)∇.𝑁^→)_𝑠Σ_gen and 𝐶_sg1 throughout the flame brush for small values of Δ in the Le=0.34, 0.6, and 0.8 flames. By contrast, negative (𝑆_𝑟+𝑆_𝑛)−𝜅_𝑚 and |∇𝑐|−𝜅_𝑚 correlations (see Figures 4(a) and 4(b)) give rise to negative values of the ensemble averaged values of ((𝑆_𝑟+𝑆_𝑛)∇⋅𝑁^→)_𝑠Σ_gen and 𝐶_sg1 throughout the flame brush for small values of Δ in the Le=1.2 flame. These local dependences are progressively smeared with increasing Δ because of the convolution operation associated with LES filtering process, and this leads to positive (negative) values of ((𝑆_𝑟+𝑆_𝑛)∇⋅𝑁^→)_𝑠Σ_gen and 𝐶_sg1 towards the unburned (burned) gas side of the flame brush for all cases considered here, including the nonunity Lewis number flames where the curvature dependences of (𝑆_𝑟+𝑆_𝑛) and |∇𝑐| are particularly strong.

The dependences of (𝑆_𝑟+𝑆_𝑛)_𝑠 and Σ_gen on 0.5×𝜕(𝑁_𝑖)_𝑠/𝜕𝑥_𝑖Σ_gen are likely to capture some of 𝜅_𝑚 dependences of (𝑆_𝑟+𝑆_𝑛) and |∇𝑐| at small values of filter widths Δ (i.e., Δ<𝛿_th), where the flame is partially resolved. This effect is particularly prevalent in the nonunity Lewis number flames where both (𝑆_𝑟+𝑆_𝑛) and |∇𝑐| are strongly correlated with curvature 𝜅_𝑚 even though the flames are statistically planar in nature. As a result of this, the contribution of (𝑆_𝑟+𝑆_𝑛)_𝑠𝜕(𝑁_𝑖)_𝑠/𝜕𝑥_𝑖Σ_gen remains close to that of ((𝑆_𝑟+𝑆_𝑛)∇⋅𝑁^→)_𝑠Σ_gen for small filter widths (i.e., Δ<𝛿_th) for the non-unity Lewis number flames, which is reflected in the small contribution of 𝐶_sg1 (see Δ=0.8𝛿_th variations in Figures 3(a)–3(c) and Figure 3(e)). The correlation between the resolved quantities (e.g., dependences of (𝑆_𝑟+𝑆_𝑛)_𝑠 and Σ_gen on 0.5×𝜕(𝑁_𝑖)_𝑠/𝜕𝑥_𝑖Σ_gen) weakens with increasing filter width Δ due to smearing of local information. Moreover, physical processes take place increasingly at the subgrid level for Δ≫𝛿_th, and thus (𝑆_𝑟+𝑆_𝑛)_𝑠𝜕(𝑁_𝑖)_𝑠/𝜕𝑥_𝑖Σ_gen does not capture the behaviour of ((𝑆_𝑟+𝑆_𝑛)∇⋅𝑁^→)_𝑠Σ_gen for large filter widths in all cases considered here, including the nonunity Lewis number flames where the curvature dependences of (𝑆_𝑟+𝑆_𝑛) and |∇𝑐| are particularly strong. This leads to 𝐶_sg1≈((𝑆_𝑟+𝑆_𝑛)∇⋅𝑁^→)_𝑠Σ_gen for Δ≫𝛿_th in all cases considered here (see Δ=2.4𝛿_th variations in Figures 3(f)–3(j)). It can be seen from Figure 3 that the positive contribution of 𝐶_sg1 overcomes the negative contribution of 𝐶_sg2 towards the unburned gas side of the flame brush in the Le=0.34 and 0.6 flames, which lead to positive value of 𝐶_sg=𝐶_sg1+𝐶_sg2 towards the unburned gas side for all values of Δ in these cases (see Figure 2). By contrast, negative values of 𝐶_sg2 overcome the positive contributions of 𝐶_sg1 towards the unburned gas side of the flame brush in the Le=0.8, 1.0, and 1.2 flames, which lead to negative values of 𝐶_sg=𝐶_sg1+𝐶_sg2 throughout the flame brush in these cases (see Figure 2).

The subgrid fluctuations of the surface-weighted contributions of (𝑆_𝑟+𝑆_𝑛) and ∇⋅𝑁^→ are scaled here using 𝑆_𝐿 and (Σ_gen−|∇𝑐|), respectively, to propose the following model for 𝐶_sg1: 𝐶_sg1𝛽=−⁴Σ_gen−|∇𝑐|𝑐−𝑐^∗𝑆_𝐿Σ_genexp−𝑎_Σ(1−𝑐)𝑐(1−𝑐)^𝑚, (14) where 𝛽₄,𝑐^∗,𝑎_Σ, and 𝑚 are the model parameters. The function (𝑐−𝑐^∗)/{exp[−𝑎_Σ(1−𝑐))]𝑐(1−𝑐)^𝑚} in (14) is used to capture the correct qualitative behaviour of 𝐶_sg1 throughout the flame brush. In a compressible, LES simulation 𝑐̃ is readily available and 𝑐 needs to be extracted from 𝑐̃.The methodology of extracting 𝑐 from 𝑐̃ in the context of LES was discussed elsewhere [9, 10, 12] and will not be discussed in detail in this paper. The model parameter 𝑐^∗ ensures that the transition from positive to negative value of 𝐶_sg1 takes place at the correct location within the flame brush. The quantity (Σ_gen−|∇𝑐|) vanishes when the flow is fully resolved (i.e., lim_Δ→0(Σ_gen−|∇𝑐|)=lim_Δ→0(|∇𝑐|−|∇𝑐|)=|∇𝑐|−|∇𝑐|=0.0), and thus 𝐶_sg1 becomes exactly equal to zero when the flow is fully resolved (i.e., Δ→0) according to (14). It has been found that 𝑚=1.85 enables (14) to capture the qualitative behaviour of 𝐶_sg1 when the optimum values of 𝑐^∗and 𝑎_Σ are chosen. The optimum values of 𝑐^∗(𝑎_Σ) tend to increase with decreasing (increasing) Δ. The curvature 𝜅_𝑚 dependences of (𝑆_𝑟+𝑆_𝑛) and |∇𝑐| are influenced by Le (see Figures 4(a) and 4(b)), and these local dependences also appear in the resolved scale but their strength diminishes with increasing Δ due to filtering operation. As the resolved and subgrid curvature terms are closely related [9, 10, 12], the qualitative behaviour of 𝐶_sg1 is also affected by the curvature dependences of displacement speed components and scalar gradient at the resolved scale, which leads to the variation of the optimum values of 𝑎_Σ, 𝛽₄, and 𝑐^∗ with Le and Δ. The model parameter 𝛽₄ needs to be deceased for decreasing values of Σ_gen for satisfactory prediction of (14). The prediction of (14) ensemble averaged on 𝑐̃ isosurfaces is compared with the ensemble averaged values of C_sg1 in Figure 5 for all cases considered here for the optimum values of 𝛽₄,𝑐^∗, and 𝑎_Σ for Δ=0.8𝛿_th and Δ=2.4𝛿_th when 𝑚 is taken to be 𝑚=1.85. The optimum values of 𝛽₄,𝑐^∗, and 𝑎_Σ are estimated by calibrating the prediction of (14) with respect to the ensemble averaged values of 𝐶_sg1 obtained from DNS data and the variation of the optimum values of 𝛽₄/Σ_gen, 𝑐^∗, and 𝑎_Σ with Δ/𝛿_th for all cases are shown in Figure 6. The optimum values of 𝛽₄/Σ_gen, 𝑐^∗, and 𝑎_Σ are parameterised here as 𝛽₄Σ_gen=9.81𝛿_th𝑙₁+𝑙₂−𝑙₁1.0+exp[−10.0(Le−1)]^1/2, (15a) where 𝑙₁Δ=1.2^2.79+1.2Δ+𝛿_th^2.79Δ+𝛿_th^2.79;𝑙₂Δ=1.34^0.67+0.53Δ+𝛿_th^0.673.1Δ^0.67+0.1Δ+𝛿_th^0.67, (15b)𝑐^∗=𝑘₁+𝑘₂−𝑘₁1.0+exp−2.0Δ/𝛿_th;𝑎−1.5_Σ=𝑘₄1.0+exp−5.0Δ/𝛿_th,−1.0 (15c) where 𝑘₁=0.75+0.151.0+exp−5.0𝑘₃;𝑘−4.6₂=0.65+0.051.0+exp−9.0𝑘₃,−4.0 (15d)𝑘₄=0.81−0.671.0+exp−5.0𝑘₃;𝑘−4.6₃=Re_Δ^0.83+0.1Δ/𝛿_th^1.73+0.1;Re_Δ=4𝜌₀Δ𝜇₀2𝑘3. (15e)

(a)

(b)

(c)

(d)

(e)

(f)

(g)

(h)

(i)

(j)

Figure 5 

Variations of 𝐶sg1 (solid line) and 𝐶sg2 (dashed line) conditionally averaged in bins of 𝑐̃ across the flame brush along with the predictions of (14) (vertical line) and (16) (crosses line) for Δ=8Δ𝑚≈0.8𝛿th (top row): (a) Le=0.34 (1st column); (b) Le=0.6 (2nd column); (c) Le=0.8 (3rd column); (d) Le=1.0 (4th column); (e) Le=1.2 (5th column) and for filter size Δ=24Δ𝑚≈2.4𝛿th (bottom row): (f) Le=0.34 (1st column); (g) Le=0.6 (2nd column); (h) Le=0.8 (3rd column); (i) Le=1.0 (4th column); (j) Le=1.2 (5th column).

(a)

(b)

(c)

(d)

(e)

Figure 6 

Variations of the model parameters 𝛽1 (◯), 𝛽2(□), 𝛽3(▵), 𝛽4/Σgen (down-pointing triangle), 𝛽5(+), 𝑎Σ(*), and 𝑐∗(×) with Δ for: (a) Le=0.34; (b) Le=0.6; (c) Le=0.8; (d) Le=1.0; (e) Le=1.2.

Figure 5 shows that (14) satisfactorily predicts 𝐶_sg1 when 𝑚 is taken to be 𝑚=1.85, and the optimum values of 𝛽₄,𝑐^∗, and 𝑎_Σ are used. According to the parameterisation given by (15a)–(15e), 𝛽₄ increases with decreasing Le, as the effects of chemical reaction strengthen with decreasing values of Lewis number (see Table 2). Moreover, 𝛽₄/Σ_gen, 𝑐^∗, and 𝑎_Σ approach to asymptotic values for large values of Δ and turbulent Reynolds number based on LES filter width Re_Δ=4𝜌₀2𝑘/3Δ/𝜇₀, where 𝜌₀ and 𝜇₀ are the unburned gas density and viscosity, respectively.

Here, the contribution of (𝐷𝜅²_𝑚)_𝑠−(𝐷𝜕𝑁_𝑖/𝜕𝑥_𝑖)_𝑠𝜕(𝑁_𝑖)_𝑠/𝜕𝑥_𝑖 is scaled with (Ξ_Δ−1)^𝑛𝑆_𝐿Σ_gen (i.e., (𝐷𝜅²_𝑚)_𝑠−(𝐷𝜕𝑁_𝑖/𝜕𝑥_𝑖)_𝑠𝜕(𝑁_𝑖)_𝑠/𝜕𝑥_𝑖∼(Ξ_Δ−1)^𝑛𝑆_𝐿Σ_gen) where the sub-grid fluctuations of 𝐷 are taken to scale with 𝑆_𝐿/Σ_gen. The above relations are utilised here to propose a model for 𝐶_sg2 in the following manner: 𝐶_sg2𝛽=−⁵𝑆_𝐿Ξ_Δ−1^𝑛Σ²_gen𝑐(1−𝑐), (16) where Ξ_Δ=Σ_gen/|∇c| is the wrinkling factor [8, 11, 43, 50, 51], 𝛽₅ and 𝑛 are the model parameters, and 𝑐(1−𝑐) is used to capture the correct qualitative behaviour of 𝐶_sg2. The subgrid curvature term 𝐶_sg2 vanishes when the flow is fully resolved according to (16), (i.e., lim_Δ→0Ξ=lim_Δ→0Σ_gen/|∇𝑐|=lim_Δ→0|∇𝑐|/|∇𝑐|=|∇𝑐|/|∇𝑐|=1.0). It has been found that (16) satisfactorily captures the behaviour of 𝐶_sg2 throughout the flame brush for 𝑛=1.0 in all cases when a suitable value of 𝛽₅ is used. The variation of the global mean optimum values of 𝛽₅ with Δ/𝛿_th is shown in Figure 6 for all cases considered here. The optimum values of 𝛽₅ have been parameterised here in the following manner: 𝛽₅=𝑚(Le)Re_ΔRe_Δ×𝑟+1.0₁+𝑟₂−𝑟₁1.0+exp−5.0Re_Δ−𝑟₃, (17a) where 𝑟₁𝑟=1.6^41.23+6.247.17𝑟₄^1.23;𝑟+0.26₂𝑟=1.88^42.27+5.928.47𝑟₄^2.27+0.47𝑟₃=35.0erfexp5.3r₄;𝑟−1.0₄=ΔΔ+𝛿_th𝑚(Le)=𝑟₅+1.0−𝑟₅1.0+exp[−10.0(Le−1.0)]^1/4;𝑟₅𝑟=0.46^45.22+4.538.0𝑟₄^5.22.+2.96 (17b)(17c)(17d)

The predictions of (16) ensemble averaged on 𝑐̃ isosurfaces are compared with ensemble averaged values of 𝐶_sg2 in Figure 5 for all cases at Δ=0.8𝛿_th and Δ=2.4𝛿_th, which show that (16) satisfactorily predicts the statistical behaviour of 𝐶_sg2 when n is taken to be 𝑛=1.0 and the optimum value of 𝛽₅ is used. According to (17a)–(17d), 𝛽₅ approaches to asymptotic values for large values of Δ and Re_Δ=4𝜌₀2𝑘/3Δ/𝜇₀.

Equations (13) and (15a)–(15e) can be combined to propose a model for 𝐶_sg in the following manner: 𝐶_sg=−𝛽₄Σ_gen−|∇𝑐|𝑐−𝑐^∗𝑆_𝐿Σ_genexp−𝑎_Σ(1−𝑐)𝑐(1−𝑐)^𝑚−𝛽₅𝑆_𝐿Ξ_Δ−1^𝑛Σ²_gen𝑐(1−𝑐). (18)

The above model will henceforth be referred to CSGNEW model in this paper. Equation (18) allows for a positive contribution of 𝐶_sg through the contribution of −𝛽₄(Σ_gen−|∇𝑐|)(𝑐−𝑐^∗)𝑆_𝐿Σ_gen/{exp[−𝑎_Σ(1−𝑐)]𝑐(1−𝑐)^𝑚}, which is absent in the CSGCAND, CSGCANT, and CSGCHAR models. The predictions of the CSGCAND, CSGCANT, CSGCHAR, and CSGNEW models for Δ=0.8𝛿_th and Δ=2.4𝛿_th are compared with 𝐶_sg obtained from DNS in Figure 7 for the optimum values of 𝛽₁, 𝛽₂, 𝛽₃, and 𝛽₅. The optimum values of 𝛽₁, 𝛽₂, and 𝛽₃ are estimated by calibrating the models based on the ensemble averaged value of 𝐶_sg obtained from DNS data. The variations of the optimum values of 𝛽₁, 𝛽₂, and 𝛽₃ with Δ for all cases are also shown in Figure 6. It is evident from Figure 6 that 𝛽₁, 𝛽₂, 𝛽₃, and 𝛽₅ remain greater than unity for all cases. This is found to be consistent with the realisability analysis by Hawkes and Cant [52]. Figure 6 further demonstrates that the optimum values of 𝛽₁, 𝛽₂, and 𝛽₃ change appreciably with increasing Δ, which is consistent with earlier findings [9, 10, 12]. Moreover, optimum values of 𝛽₁, 𝛽₂, and 𝛽₃ for a given Δ are affected by Le (see Figure 6). It is worth noting that parameterisation of the optimum values of 𝛽₁, 𝛽₂ and 𝛽₃ also yields complex relations similar to (15a)–(15e) and (17a)–(17d). However, such parameterisation is not presented here because the CSGCAND, CSGCANT, and CSGCHAR models do not capture the qualitative behaviour of 𝐶_sg for the Le=0.34 and 0.6 flames.

(a)

(b)

(c)

(d)

(e)

(f)

(g)

(h)

(i)

(j)

Figure 7 

Variation of 𝐶𝑠𝑔 (solid line) conditionally averaged in bins of 𝑐̃ across the flame brush along with the predictions of CSGCAND (asterisk line), CSGCANT (dashed line), CSGCHAR (squared line) and CSGNEW (triangle line) for filter sizes Δ=8Δ𝑚≈0.8𝛿th (top row): (a) Le=0.34 (1st column); (b) Le=0.6 (2nd column); (c) Le=0.8 (3rd column); (d) Le=1.0 (4th column); (e) Le=1.2 (5th column) and for filter size Δ=24Δ𝑚≈2.4𝛿th (bottom row): (f) Le=0.34 (1st column); (g) Le=0.6 (2nd column); (h) Le=0.8 (3rd column); (i) Le=1.0 (4th column); (j) Le=1.2 (5th column).

It can further be seen from Figure 7 that the CSGCHAR model tends to overpredict the negative values of 𝐶_sg towards the unburned gas side in cases C–E (𝐿𝑒 = 0.8, 1.0 and 1.2 flames), and this behaviour becomes more prominent with increasing filter size. It is clear from Figure 7 that for Δ=24Δ_𝑚=2.4𝛿_th, the CSGCHAR model predicts the maximum magnitude of 𝐶_sg near the middle of the flame whereas the actual maximum magnitude of 𝐶_sg is attained slightly towards the burned gas side. The CSGCAND and CSGCANT models give comparable performance for optimum values of 𝛽₁ and 𝛽₂ in cases C–E. However, the CSGCAND and CSGCANT models do not satisfactorily capture the qualitative behaviour of 𝐶_sg and underpredict (overpredict) the magnitude of 𝐶_sg towards the burned gas side (middle) of the flame brush in the Le = 0.8, 1.0 and 1.2 flames. Figure 7 demonstrates that the CSGNEW model captures the qualitative behaviour of 𝐶_sg in a better manner than the CSGCAND and CSGCANT models and the quantitative agreement between 𝐶_sg and the CSGNEW model remains better than the CSGCAND, CSGCANT, and CSGCHAR models in all cases for all values of Δ when optimum values of 𝛽₄, 𝛽₅, 𝑎_Σ, and 𝑐^∗ are used.

5. Conclusions

The LES modelling of the curvature term (𝑆_𝑑∇⋅𝑁^→)_𝑠Σ_gen of the generalised FSD Σ_gen transport equation has been addressed here using a simplified chemistry-based DNS database of freely propagating statistically planar turbulent premixed flames with Lewis number Le ranging from 0.34 to 1.2. The statistical behaviours of the subgrid curvature term 𝐶_sg for a range of different values of Δ have been analysed in terms of its contributions 𝐶_sg1 and 𝐶_sg2 arising from the combined reaction and normal diffusion component and tangential diffusion components of displacement speed (i.e., (𝑆_𝑟+𝑆_𝑛) and 𝑆_𝑡=−2𝐷𝜅_𝑚), respectively. The Lewis number is shown to have significant influences on the statistical behaviours of the resolved and subgrid components of the FSD curvature term. Detailed physical explanations have been provided for the observed filter size and Lewis number dependences of the different components of (𝑆_𝑑∇⋅𝑁^→)_𝑠Σ_gen. Models have been identified for individual components of the subgrid curvature term (i.e., 𝐶_sg1 and 𝐶_sg2), and the performances of these models have been compared to the corresponding quantities extracted from DNS data. It has been found that the new models for 𝐶_sg1 and 𝐶_sg2 satisfactorily capture the statistical behaviours of the corresponding terms extracted from DNS data. It has been found that the existing models for the subgrid curvature term 𝐶_sg do not satisfactorily capture the qualitative behaviour of the corresponding quantity extracted from DNS data for all the flames considered here. This problem is particularly prevalent for flames with small values of Lewis number (e.g., Le=0.34 and 0.6) where 𝐶_sg locally assumes positive values, whereas the existing models can only predict negative values of 𝐶_sg. The performance of the newly proposed model for 𝐶_sg has been found to be better than the existing models, and it has been shown to capture positive contributions of 𝐶_sg for the Le≪1 flames. The present analysis has been carried out using a DNS database with moderate value of Re_𝑡 in the absence of the effects of detailed chemistry and transport. As simplified chemistry-based DNS qualitatively captures the curvature 𝜅_𝑚=∇⋅𝑁^→/2 and strain rate dependences of 𝑆_𝑑 and |∇c| obtained from detailed chemistry based simulations, it can be expected that the statistical behaviours of the curvature term (𝑆_𝑑∇⋅𝑁^→)_𝑠Σ_gen presented in this paper will be valid at least in a qualitative sense in the context of detailed chemistry. However, the quantitative values of the model parameters (i.e., 𝛽₄, 𝛽₅, 𝑎_Σ, and 𝑐^∗) may need to be altered in the presence of detailed chemistry. Thus, three-dimensional DNS data with detailed chemistry and experimental data at higher values of Re_𝑡 will be necessary for more comprehensive modelling of (𝑆_𝑑∇⋅𝑁^→)_𝑠Σ_gen and 𝐶_sg in the context of LES. Moreover, the newly proposed models need to be implemented in LES simulations for the purpose of a posteriori assessments. However, this is kept beyond the scope of this paper. Several previous studies [3–7, 9–12, 43–49] concentrated purely on the model development based on a priori analysis of DNS data and the same approach has been adopted here. Implementation of the newly developed models in LES simulations will form the basis of future investigations.

Acknowledgment

The financial support by EPSRC, UK, is gratefully acknowledged.