Effects of Mean Inflow Velocity and Droplet Diameter on the Propagation of Turbulent V-Shaped Flames in Droplet-Laden Mixtures

: Three-dimensional carrier phase Direct Numerical Simulations of V-shaped n-heptane spray flames have been performed for different initially mono-sized droplet diameters to investigate the influence of mean flow velocity on the burning rate and flame structure at different axial locations from the flame holder. The fuel is supplied as liquid droplets through the inlet and an overall (i.e., liquid + gaseous) equivalence ratio of unity is retained in the unburned gas. Additionally, turbulent premixed stoichiometric V-shaped n-heptane flames under the same turbulent flow conditions have been simulated to distinguish the differences in combustion behaviour of the pure gaseous phase premixed combustion in comparison to the corresponding behaviour in the presence of liquid n-heptane droplets. It has been found that reacting gaseous mixture burns predominantly under fuel-lean mode and the availability of having fuel-lean mixture increases with increasing mean flow velocity. The extent of flame wrinkling for droplet cases has been found to be greater than the corresponding gaseous premixed flames due to flame-droplet-interaction, which is mani-fested by dimples on the flame surface, and this trend strengthens with increasing droplet diameter. As the residence time of the droplets within the flame decreases with increasing mean inflow velocity, the droplets can survive for larger axial distances before the completion of their evaporation for the cases with higher mean inflow velocity and this leads to greater extents of flame-droplet interaction and droplet-induced flame wrinkling. Mean inflow velocity, droplet diameter and the axial distance affect the flame brush thickness. The flame brush thickens with increasing droplet diameter for the cases with higher mean inflow velocity due to the predominance of fuel-lean gaseous mixture within the flame. However, an opposite behaviour has been observed for the cases with lower mean inflow velocity where the smaller extent of flame wrinkling due to smaller values of integral length scale to flame thickness ratio arising from higher likelihood of fuel-lean combustion for larger droplets dominates over the thickening of the flame front. It has been found that the major part of the heat release arises due to premixed mode of combustion for all cases but the contribution of non-premixed mode of combustion to the total heat release has been found to increase with increasing mean inflow velocity and droplet diameter. The increase in the mean inflow velocity yields an increase in the mean values of consumption and density-weighted displacement speed for the droplet cases but leads to a decrease in turbulent burning velocity. By contrast, an increase in droplet diameter gives rise to decreases in turbulent burning velocity, and the mean values of consumption and density-weighted displacement speeds. Detailed physical explanations have been provided to explain the observed mean inflow velocity and droplet diameter dependences of the flame propagation behaviour.

Vena et al. [8,9] considered turbulent stratified V-flames for iso-octane-air mixtures and indicated that the mixture inhomogeneity plays a key role in the burning rate characteristics and flame propagation rate. These aspects are further strengthened in turbulent droplet-laden combustion due to localised mixture inhomogeneities induced by the evaporation of droplets [10,20]. The mixture inhomogeneity in spray flames is determined by the number density and droplet diameter, which has been demonstrated in several analytical [21][22][23], experimental [24][25][26][27][28][29] and computational  analyses, and this aspect is also valid for turbulent V-flames propagating into droplet-laden mixtures, as demonstrated in a recent analysis by Ozel-Erol et al. [20]. The analysis of Ozel-Erol et al. [20] concentrated on the effects of droplet diameter of the mono-sized droplets supplied ahead of the flame holder for an overall (considering fuel in both gaseous and liquid phases) equivalence ratio of unity (i.e., = 1.0). It was revealed by Ozel-Erol et al. [20] that combustion in spray V-flames with = 1.0 takes place predominantly under fuel-lean premixed mode and the probability of fuel-lean combustion increases with increasing droplet diameter. Moreover, the propensity of fuel-lean combustion decreases with increasing axial distance from the flame holder. It has been found that the influence of the droplet diameter on the density-weighted displacement speed and consumption speed in V-shaped flames is weaker than in the case of spherically expanding flames propagating into droplet-laden mixtures [20].
It is well-known that a variation of the mean inflow velocity changes the flame angle and flame brush in V-shaped flames. Furthermore, a change in mean inflow velocity also alters the ratio of the residence time to the evaporation time of the droplets within Vshaped flames, which is likely to have influences on the mixture composition within the flame. This, in turn, affects the burning rate and flame surface area in turbulent V-shaped flames propagating into droplet-laden mixtures. However, this aspect is yet to be addressed in the existing literature and this paper aims to fill this gap by carrying out carrierphase DNS of turbulent V-shaped flames propagating into mono-sized n-heptane dropletladen mixtures corresponding to an overall equivalence ratio of unity (i.e., = 1.0) for different mean inflow velocities and initial droplet diameters . In this respect, the main objectives of this paper are: (a) to analyse the effects of mean inflow velocity and initial droplet diameters on the flame structure at different axial locations from the flame holder for turbulent V-shaped flames propagating into droplet-laden mixtures with mono-sized droplet diameters.
(b) to explain the effects of mean inflow velocity and initial droplet diameter on consumption speed, density-weighted displacement speed and turbulent burning velocity at different axial locations from the flame holder.
The rest of this paper will be organised in the following manner. The necessary mathematical background is provided in the next section of this paper. This will be followed by a brief discussion on the numerical implementation related to this analysis. The results will be presented and subsequently discussed in Section 4 before the main findings are summarised and conclusions are drawn in the final section of this paper.

Mathematical Background
In the interest of an extensive parametric analysis, a modified single-step Arrheniustype chemical reaction [51] is considered for n-heptane-air combustion. This modified single-step chemical reaction provides a realistic gaseous phase equivalence ratio dependence of the unstrained laminar burning velocity ( ) and adiabatic flame temperature ( ) in hydrocarbon-air flames, and therefore has been chosen to keep the computational cost within reasonable limits. In this chemical mechanism both the activation energy and the heat of combustion are taken to be functions of gaseous phase equivalence ratio and interested readers are referred to Refs. [20,52] where it was shown that the variations of the normalised unstrained laminar burning velocity / and adiabatic flame temperature with for the present thermo-chemistry agree well with both previously reported experimental data and detailed chemistry simulations. It is worth noting that several previous DNS analyses [30,31,[33][34][35][36][37]40,[43][44][45][46][47][48][49][50] used simple chemistry for simulating combustion of droplet-laden mixtures and the same approach has been adopted in this analysis. For the purpose of simplicity, the Lewis numbers of all species are assumed to be unity and all species in the gaseous phase are taken to be perfect gases following several previous studies. The ratio of specific heats ( = / = 1.4, where and are the gaseous specific heats at constant pressure and constant volume, respectively) and Prandtl number ( = / = 0.7 where is the dynamic viscosity and is the thermal conductivity of the gaseous phase) are taken to have standard values. A Lagrangian approach has been adopted to transport the position ⃗ , velocity ⃗ , diameter and temperature of the individual droplets in the following manner [30,31,[33][34][35][36][37]40,[43][44][45][46][47]: where is the gaseous phase temperature, is the latent heat of vaporization and , and stand for the relaxation timescales for droplet velocity, diameter and temperature, respectively, and they are given by [30,31,[33][34][35][36][37]40,[43][44][45][46][47] Here, is the value of fuel mass fraction at the droplet surface, which needs the knowledge of the partial pressure of the fuel vapour at the droplet surface . This is achieved by using the Clausius-Clapeyron equation which leads to the following expressions [30,31,[33][34][35][36][37]40,[43][44][45][46][47]: In Equation (4), is the gas constant, is the boiling point of the fuel at the reference pressure and is assumed to be , and and are the molecular weights of gaseous mixture and fuel, respectively. The knowledge of is necessary to calculate the Spalding number . The Lagrangian tracking of liquid fuel droplets are coupled with the Eulerian treatment of the carrier gaseous mixture with additional source terms in the gaseous transport equations in the following manner [20,36,37,40,[43][44][45][46][47]: In the context of Equation (5), = {1, , , , } and = {1, , , , } for the conservation equations of mass, momentum, energy, and mass fractions, respectively and = / is the generalised diffusion coefficient (where is the kinematic viscosity, and is an appropriate Schmidt number corresponding to ) associated with = {1, , , } whereas is given by = for = , and is the component of fluid velocity. The term ̇ in Equation (5) is the reaction rate contribution, ̇ is an appropriate source term and ̇ refers to the term associated with the coupling between Lagrangian particulate and Eulerian carrier phases. A tri-linear interpolation has been employed for the purpose of the evaluation of ̇ from the droplet location ⃗ to the eight surrounding nodes [20,36,37,[40][41][42][43][44][45][46][47].
It is useful to define a reaction progress variable , based on oxygen mass fraction, and mixture fraction, = ( − / + / ) ( + / ) ⁄ , in the following manner in order to identify the reaction zone [20,36,37,[40][41][42][43][44][45][46][47]: Here, = 0.233 represents the oxygen mass fraction in air, whereas = 1.0 is the fuel mass fraction in the pure fuel stream. For n-heptane, C H as the fuel in the analysis, = 3.52 denotes the stoichiometric mass ratio of oxygen to fuel and = 0.0621 and = 0.0621 represent the corresponding stoichiometric fuel mass fraction and stoichiometric mixture fraction, respectively. The transport equation of the reaction progress variable takes the following form [20,43,[45][46][47]: where represents the molecular diffusivity, ̇ is the reaction rate of reaction progress variable, ̇ , is the source/sink term arising due to droplet evaporation and ̇ represents the cross-scalar dissipation term arising due to mixture inhomogeneity. The reaction rate of the reaction progress variable ̇ is given by [20,43,[45][46][47]: In Equation (8), ̇ is the reaction rate of oxygen and ̇ , in Equation (7) is calculated as [20,43,[45][46][47]: Here, ̇ = ( ̇ − ̇ / )/( + ⁄ ) is the source/sink term in the mixture fraction transport equation with ̇ = Γ (1 − ) and ̇ = −Γ being the droplet source/sink terms in the mass fraction transport equations for fuel and oxygen, respectively and Γ is the source term in the mass conservation equation due to evaporation [20,43,[45][46][47]. The term associated with mixture inhomogeneity ̇ is expressed as [20,43,[45][46][47]: Equation (7) can alternatively be written in the kinematic form for a given isosurface in the following manner [20,43,[45][46][47]: Here stands for the displacement speed, which is given by [20,43,[45][46][47]53,54]: Equation (12) suggests that the density-weighted displacement speed can be defined as * = / where is the unburned gas density. Another alternative flame speed, known as the consumption speed , is often used in the combustion literature, which is defined as [20,47]: where is the elemental distance in the local flame normal direction.
In SENGA+, the governing equations of mass, momentum, energy and species of the gaseous phase are solved in non-dimensional form. The spatial discretisation for the internal grid points is carried out using a 10th order central difference scheme, but the order of differentiation gradually decreases to a one-sided 2nd order scheme at the non-periodic boundaries. An explicit low-storage 3rd order Runge-Kutta scheme [60] is used for explicit time-advancement. The boundaries in the stream-wise direction ( -direction) are taken to be inflow and outflow and are specified using to the Navier-Stokes Characteristic Boundary Conditions (NSCBC) technique [61]. The remaining transverse directions (i.e., and ) are considered to be periodic and thus no explicit boundary conditions are needed for and directions. The simulation domain is taken to be a cube of = = = 30 with = ( − )/ max| ∇ | being the thermal flame thickness of the unstrained stoichiometric laminar premixed flame. A uniform Cartesian grid of 384 is used to discretise the aforementioned computational domain, which ensures about 10 grid points within . A commercial software COSILAB (Rotexo-Softpredict-Cosilab) [62] is utilised to obtain steady-state planar spray flame solutions as previously done by Neophytou and Mastorakos [63], and these laminar flame solutions are utilised to initialise density, species mass fractions, temperature, flame normal velocity and droplet properties for three different initial values of droplet diameter (i.e., ⁄ = 0.04, 0.05 and 0.06) for an overall equivalence ratio of unity (i.e., = 1.0). All the V-flame cases taken to have a moderate level of inlet turbulence intensity (i.e., / is the integral length scale based on turbulent kinetic energy and its dissipation rate ) for the flames considered here remains moderate, the value of remains comparable to that used in experiments by Vena et al. [8,9] and previous V-premixed flame DNS analyses by Dunstan et al. [15,16]. A small value of , ⁄ is chosen for this study, as it allows for the analysis of flame-droplet interaction, and this aspect is eclipsed by flame-turbulence interaction for large values of , ⁄ [26,29,45]. A frozen field of homogeneous isotropic turbulence is scanned by a plane in order to specify turbulent velocity fluctuations at the inlet plane following Taylor's hypothesis [64], as done in several previous analyses [15,16,20,65,66]. The flame holder is located at approximately 10 from the inlet plane. Three different axial distances (i.e., position A, B and C) from the flame holder (i.e., 9 , 13 and 17 ) are considered for the statistics presented in this paper, which are taken after 2.0 where = / is the flow through time. A Gaussian weighting function ( ) = [− /2 ] following Dunstan et al. [5,16] is used to impose the reactant mass fractions (i.e., and ) and the non-dimensional temperature of gaseous phase within the flame holder volume, where and are adjustable constants to obtain the desired flame holder radius (i.e., = 0.3 ).
The initial droplet number density is estimated to be 1.23 ≤ ( ) ⁄ ≤ 1.85 in the unburned gas for the cases considered here and the liquid volume fraction is well below 0.01%. Individual droplets are tracked using a Lagrangian approach in the context of point-source assumption which is extensively used for DNS of diluted spray flames [30][31][32][33][34]36,37,[40][41][42][43][44][45][46][47]. The point source formulation is valid when the droplet diameter is smaller than both the Kolmogorov length scale and the grid size used in DNS. For the inlet turbulence intensity (i.e., , ⁄ = 2.0) considered here, the ratio of initial droplet diameter to the Kolmogorov length scale is ⁄ = 0.09, 0.11, 0.13 for the cases of ⁄ = 0.04, 0.05, 0.06, respectively. This ratio of the initial droplet diameter to grid spacing remains comparable to several previous analyses [31,34,36,[39][40][41][42]. Furthermore, the diameter of the droplets interacting with the flame in the reaction zone reduces significantly (i.e., at least 50% smaller than the initial size) due to the relatively high volatility of n-heptane. The sub-grid evaporation model assumption is not expected to have an important impact on the analysis of flame-droplet interaction, as flame-droplet interaction occurs at larger scales than the initial size of the droplets. Therefore, it can be expected that the qualitative nature of the findings is not sensitive to the evaporation modelling used in this analysis. Additionally, validity of the point source assumption for droplets with sizes smaller than the Kolmogorov length scale is supported by a recent comparative analysis by Haruki et al. [67] with respect to the evaporation characteristics obtained from fully resolved multiphase simulations. Carrier phase DNS with point source and fully Eulerian phase-DNS were considered by de Chaisemartin et al. [68] to compare their performance in terms of polydispersity and droplet crossing. The results by de Chaisemartin et al. [68] showed a good agreement between gaseous fuel mass fraction fields obtained from these two approaches. They [68] also reported that the simulations with the full Eulerian approach were 10 times more expensive than those with the point source approach and a large-scale parallelisation capability for 3D simulations is needed for the Eulerian phase DNS for sprays. Moreover, combustion simulations involving fully resolved dispersed phase are currently in a primitive stage and need further validation to be used on a routine basis. Thus, a point source assumption for the dispersed phase (i.e., droplets) has been considered in this analysis.
The Stokes number for this problem can be estimated based on turbulent flow time- ⁄ is the particle time scale and /√ is the turbulent time scale) and it has been found that remains smaller than 0.13 for the cases considered here. An alternative Stokes number can be defined based on the chemical time scale / , (where is the thermal diffusivity in the unburned gas) as The group number helps to categorise combustion process of droplet-laden mixtures as individually burning droplets ( ≪ 1.0) and external sheath combustion ( ≫ 1.0) [69]. According to Chiu and Liu [69] the group number is expressed as where and are the Lewis and Schmidt numbers, respectively, is the number of droplets in a specified volume. For all droplet cases, is estimated to much greater than unity which is indicative of the external sheath combustion [30,70]. A recent analysis by Zhao et al. [70] indicated that external sheath combustion is the most likely mode of burning under the conditions analysed in this paper.

Reacting Mixture Composition
A comparison of gaseous mixture composition for droplet cases with , ⁄ = 5.0 and 10 can be made for different initial droplet diameters (i.e., ⁄ = 0.04, 0.05, 0.06) from the distributions of gaseous phase equivalence ratio in the midplane presented in Figure 1. The droplets residing on the plane are illustrated by grey dots. Inlet boundary on the left side of the domain supplies the fuel in the form of liquid droplets with a temperature of = 300K. The mass flux of liquid fuel at the inlet plane is determined based on the values of and . Evaporation of the liquid droplets provides gaseous fuel to maintain the V-shaped flame. Figure 1 shows that gaseous equivalence ratio in the unburned gas region can reach at the stoichiometric value (i.e., = 1.0) for small droplet cases with initial ⁄ = 0.04 due to their high evaporation rates, however, mostly remains smaller than unity (i.e.,  , and consequently in these cases less fuel vapour becomes available in the gaseous phase than in the cases with smaller mean flow velocity. Furthermore, large droplets, which reach at the flame and escape through the flame without complete evaporation, eventually evaporate completely in the high temperature zone in the burned gas (i.e., ≈ 1.0). The resulting evaporated fuel vapour diffuses from the high temperature burned gas region into the flame and mixes with excess air diffusing from the unburned gas side to give rise to local pockets of stoichiometric mixture (i.e., = 1.0) where diffusion mode of combustion is obtained. The predominance of fuellean combustion in the gaseous phase for higher mean inflow velocity cases can be veri- , the gaseous phase combustion nominally takes places in the thin reaction zones regime [71]. The predominance of low Damköhler number combustion can be substantiated from the PDFs of reaction progress variable for different values of Favre-averaged reaction progress variable . Figure 3 shows the PDFs of reaction progress variable for the Favre-averaged reaction progress variable values ̃= 0.1, 0.5 and 0.9 at locations A, B and C, whereas the contours of ̃= 0.1, 0.5 and 0.9 for the cases considered here are shown in Figure 4. In this paper, the Favre averaging is conducted by averaging in time (i.e., over one flow-through time = / ) and space (i.e., spanwise z-direction). The PDFs of are supposed to show a bimodal distribution with peaks at = 0.0 and = 1.0 for high values of Damköhler number [72]. However, Figure 3 reveals that the PDFs of do not exhibit bimodal distribution for the most cases considered in this study. The departure of PDF of from bimodal distribution with impulses at = 0.0 and = 1.0 shows the tendency of low Damköhler number combustion for the cases considered in this study. Figure 3 further demonstrates that the PDFs of for the cases considered here show mostly mono-modal distributions and this tendency strengthens with increasing mean flow velocity.

Flame Brush Thickness and Flame Angle
The small values of Damköhler number and high values of Karlovitz number are indicative of the thickened flame regime of combustion [71] and local occurrences of flame thickening can be discerned from the contours of in Figure 1. This along with the flame wrinkling determines the flame brush thickness in the V-flame cases considered in this analysis. The flame brush thickness ( . . ) can be defined as the distance between ̃= 0.1 and 0.9 isosurfaces for a given location, as shown in Figure 4. The values of flame brush thickness, ( . . ) / for V-flame cases under droplet-laden mixtures are presented and compared to the corresponding stoichiometric turbulent premixed gaseous flame cases in Table 1. It can be seen from Table 1 that the flame brush thickness ( . . ) / increases in the downstream and decreases with increasing droplet diameter for the cases with low mean flow velocity. On the contrary, flame brush becomes wider for large droplets for the cases with , ⁄ = 10.0. According to turbulent diffusion velocity, flame brush thickness for a V-shaped flame is correlated with the distance from flame holder and flow parameters (i.e., bulk mean velocity , root-mean-square ′ and integral length scale ) [11,64,73]. Kheirkhah and Gulder [11] demonstrated the effects of equivalence ratio on flame brush thickness under different flow conditions (i.e., / ) for V-shaped premixed flames. The flame brush thickness can alternatively be estimated using the expressions: , ̃= 1/max |∇ | and , ̅ = 1/ max |∇ ̅ |. Reported values of , ̃ and , ̅ in Table 1 [74]. This suggests that ( ) (resp. / ( ) ) is expected to increase (resp. decrease) with increasing droplet diameter due to greater likelihood of fuel-lean combustion (where The changes in flame brush thickness in the flow direction lead to the variation of flame angle for V-shaped flames. Furthermore, slightly curved nature of ̃ isosurfaces yields to significant differences in flame angle in the mean flow direction from location A to C. The contours of Favre averaged reaction progress variable presented in Figure 4 also demonstrate that the variations of mean flow velocity and droplet diameter give rise to significant differences in the angles made by the tangent on the ̃= 0.1, 0.5 and 0.9 contours with the mean flow direction (i.e., x-direction). The angles ̃ . , ̃ . and ̃ . for ̃= 0.1, 0.5 and 0.9 isosurfaces are provided in Table 2 for locations A, B and C where these angles are calculated using the slopes of the tangent on the respective ̃ isosurfaces at a given axial location. It can be seen from Figure 4 that flame angles decrease with increasing mean inflow velocity, as expected, and the flame angles ̃ . , ̃ . and ̃ . for gaseous premixed flame are greater than those for droplet cases for ̃ . in the case of small droplet diameters. This behaviour originates due to differences in the evolution of the flame brush thickness in the flow direction (see Table 1). The flame angles at a given plane remain comparable for different droplet diameters for both mean inlet velocities considered here.

Reaction Zone Structure
The variation of flame angle in the mean flow direction for different droplet sizes is consistent with the evolution of the distribution of reacting mixture composition at various ̃ isosurfaces. Figure 5 shows the PDFs of gaseous equivalence ratio on ̃= 0.1, 0.5 and 0.9 isosurfaces at different axial locations for both mean inflow velocities. Figure 5 indicates that the variation of the mean flow velocity causes significant differences in the PDFs of at ̃= 0.1. Evaporated fuel in droplet cases with higher mean flow velocity burns leaner mixture at ̃= 0.1 than that with lower mean flow velocity and the differences between the PDF profiles for different mean velocities decrease in the mean flow direction from A to C for ̃= 0.1. Additionally, distribution for ̃= 0.5 becomes comparable for a given droplet diameter with different mean flow velocities at all locations considered here. On the contrary, PDFs of for higher mean flow velocity cases attain a peak value at ≈ 1.0 towards the burned gas side of the flame brush (e.g., ̃= 0.9) which implies that spray flames with vector. In relation to these definitions, the flame normal vector points towards the reactants and a positive (negative) curvature value indicates an element of the flame surface, which is convex (concave) to the reactants. It is evident from Figure 6 that flame surfaces for the gaseous premixed flames remain smooth irrespective of the mean inflow velocity and exhibit wrinkles induced by the turbulent fluid motion. Furthermore, small droplets with , ⁄ = 5.0 rarely interact with the flame because they mostly evaporate before reaching the flame surface. However, a decrease in advection time in comparison to the evaporation timescale with increasing , ⁄ enhances the number of droplets interacting with the flame farther in the mean flow direction. This leads to greater propensity to obtain dimples on the flame surface for cases with higher mean inflow velocity. Moreover, flame wrinkling as a result of the dimples arising from flame-droplet interaction increases with increasing droplet diameter due to the slow evaporation rate of large droplets for a given mean inflow velocity. The larger droplets survive longer in the mean flow direction because of their slower evaporation rate and thus the dimples become denser for larger droplet diameters and higher mean inflow velocities. This behaviour is consistent with previous experimental [29] and numerical [45][46][47] findings for spherically expanding flames propagating into droplet-laden mixtures.  The aforementioned droplet-induced flame surface deformation can be quantified with the help of PDFs of flame surface curvature . Figure 7 shows the PDFs of local normalised flame curvature × on the = 0.8 isosurface at different planes (i.e., A, B and C) for all droplet cases and the corresponding gaseous premixed flame cases considered here. The reaction rate of reaction progress variable reaches its maximum value close to = 0.8 isosurface for the present thermo-chemistry [20], therefore thisisosurface is assumed as the flame surface for the rest of this study. Although the PDFs of local normalised flame curvature × are shown only for = 0.8 isosurface in Figure 7, they show qualitatively similar behaviour for other values of and thus are not presented for the sake of brevity. Figure 7 reveals that the PDFs of assume peak value at zero curvature at all locations for both gaseous premixed and droplet-laden flames. For droplet cases, the PDF profiles of × broaden both for negative and positive values in comparison to the corresponding gaseous premixed flames due to flame-droplet interactions. The PDFs of × for droplet cases widen with increasing mean inflow velocity, whereas an opposite trend is observed for gaseous premixed flames because the flame becomes increasingly narrow for increasing mean inflow velocity. It can be seen from Figures 1 and 4, and Table 2 that the flame angle shrinks with increasing , ⁄ and thus the gaseous premixed flames shows greater range of both positive and negative curvature values in the case of smaller mean velocity because of greater likelihood of the flame interaction with a range of different turbulent flow conditions. It has already been explained earlier that the greater availability of droplet within the flame brush increases the possibility flame-droplet interaction for higher mean inlet velocities and this yields wider distributions of × for all the droplet cases at all locations (i.e., A, B and C) with a higher possibility of finding long negative and positive tails. Differences in PDF of between two mean flow velocities for a given droplet diameter increase in the streamwise direction from A to C, since some of the droplets complete their evaporation before reaching the flame, especially farther in the streamwise direction for lower mean flow velocity cases. Moreover, the widening of × of droplet flame cases in comparison to the corresponding gaseous premixed flames increases with increasing droplet diameter because these droplets survive longer within the flame due to their slower evaporation rates and thus induce greater extents of flame-droplet interaction.
The wrinkling of −isosurfaces in these flames, in turn, have implications on the burning rate statistics in premixed flames. The mode of combustion can be characterised by the Flame index FI = (∇ • ∇ )/(|∇ ||∇ |) [75]. The FI takes a positive value for a premixed mode of combustion while FI is negative for non-premixed combustion mode. Figure 8 shows the percentages of heat release rate arising from premixed and non-premixed modes of combustion for all droplet cases considered here. It is evident from Figure  8 that premixed mode of combustion is the dominant contributor to the total heat release for all droplet cases at all locations considered in this analysis. However, the percentage of heat release associated with non-premixed combustion reaches non-negligible levels for large droplet cases particularly at higher mean flow velocity. The slow evaporation of large droplets gives rise to high extent of mixture inhomogeneity in the gaseous phase and this increases the likelihood of non-premixed combustion for large droplet cases. In the case of higher mean flow velocity, larger number of droplets pass through the flame and eventually evaporate in the burned gas region (i.e., ≈ 1.0) generating locally fuel rich regions. Unburned fuel vapour from these fuel-rich pockets in the burned gas region eventually diffuses into the flame and mixes with the air diffusing from the predominantly fuel-lean unburned gas side to give rise to local diffusion mode of burning. Therefore, the extent of heat release due to non-premixed mode of combustion increases with increasing , ⁄ . As non-premixed combustion occurs where stoichiometric mixture is obtained, the mild peaks at ≈ 1.0 in the PDFs of in Figure 2 are predominantly associated with non-premixed mode of combustion in these droplet cases.

Statistical Analysis of Flame Speeds
With the information that heat release predominantly takes place in premixed mode, it is worthwhile to analyse the flame propagation statistics in terms of density-weighted displacement speed * and consumption speed . Therefore, it is instructive to analyse the mean behaviour of the terms, which drive the statistical behaviour of densityweighted displacement speed * . Figure 9 shows the mean values of the terms ̇ , ̇ , , ̇ and ∇ • ( ∇ ) conditional on for locations A, B and C for the cases with different mean inflow velocities. Figure 9 indicates that the mean molecular diffusion term ∇ • ( ∇ ) assumes positive values towards the preheat zone (i.e., < 0.5) but becomes negative in the reaction zone (i.e., 0.5 ≤ ≤ 0.9) in contrast with the positive mean values of ̇ for all cases considered here. Furthermore, since stratified mixtures with < 1.0 and > 1.0 burn slower in the droplet cases than in the stoichiometric (i.e., = 1.0) premixed flame, the mean values of ̇ for droplets cases is expected to be smaller than the gaseous premixed case. Overwhelmingly leaner combustion trend in the case of larger droplets reduces to the mean values of ̇ with increasing and the differences in the mean values of ̇ between the gaseous premixed and droplet cases are greater for , ⁄ = 10.0 than those of the cases with , ⁄ = 5.0 due to the predominance of fuel-lean combustion for higher inflow velocities (see Figures 2 and 5). Although, the mean values of terms associated with droplet evaporation ̇ , and mixture inhomogeneity ̇ are negligible compared with the magnitudes of the mean values of ̇ and ∇ • ( ∇ ) throughout the flame in the cases with lower mean flow velocity, these mean values (i.e., ̇ , and ̇ ) slightly increase for higher mean flow velocity cases. Figure 9 further shows that the mean value of |∇ | remains comparable with the mean value of , |∇ | for the gaseous premixed V-flames with different mean flow velocities. However, the mean values of |∇ | = ∇ • ( ∇ ) + ̇ + ̇ , + ̇ have been found to be smaller than those of , |∇ | for all droplet cases considered in this analysis. It can also be observed from Figure 9 that the statistical behaviours of |∇ | and , |∇ | are not significantly influenced by the variation of the mean flow velocity for V-shaped flames. As the inverse of the maximum value of |∇ | can be taken to be a measure of the flamelet thickness (i.e., ~1/ |∇ | ) [43,45,76], a decreasing trend of , |∇ | with increasing droplet diameter and mean inflow velocity for droplet cases suggests that the flame thickens with increasing droplet diameter and mean inflow velocity due to increased propensity for fuel-lean combustion. As the premixed combustion takes place in the flamelet regime, the background turbulent motion does not affect the flame structure and as a result the mean values of ̇ , ∇ • ( ∇ ) and , |∇ | are not affected by the mean inflow velocity.
Although the statistical behaviour of density-weighted displacement speed * is strongly linked with the evolution of the aforementioned terms (i.e., ∇ • ( ∇ ), ̇ , ̇ , , ̇ ), consumption speed is dependent only on ̇ and the flame thickness. Figure 10 shows the normalised mean consumption speed / , and density-weighted displacement speed, * / , calculated on = 0.8 isosurface at locations A, B and C for both mean inflow velocities. A decreasing trend in / , is observed for increasing droplet diameter for both inflow velocities. In addition, / , for the premixed gaseous flame cases is found to be greater than the droplet cases for both mean inflow velocities, which is consistent with the behaviour of the mean ̇ presented in Figure 9. However, the mean value of / , for droplet cases increases with increasing mean inflow velocity. It is important to note from Figure 9 that the mean value of ̇ and |∇ | conditional upon do not change with increasing , ⁄ irrespective of the location for premixed flames and therefore / , remains mostly unaffected by the location and mean inflow velocity for gaseous premixed flame cases. It has already been demonstrated in the context of the discussion of Figure 9 that the flame thickens and the mean value of ̇ decreases with increasing , ⁄ and droplet diameter, and this reduction in mean reaction rate magnitude with increasing droplet diameter is responsible for the decreasing trend of / , with increasing . Although the mean value of ̇ conditional upon for droplet cases decreases with increasing   shows the same qualitative trend and thus is not shown here) and the normalised ratio of volume-integrated product formation rate to the flame surface area = ∫ ̇ /[ , ∫ |∇ | ] are reported in Table 3 for V-shaped flames with , ⁄ = 5.0 and 10.0. Table 3 reveals that increases with increasing droplet diameter for the cases with higher mean inflow velocity while values remain comparable for droplet cases with lower mean flow velocity. It can further be observed that for lower mean flow velocity case, of the gaseous premixed case is larger than those of the droplet cases with initial ⁄ = 0.04 and 0.05, however, this value remains smaller than that for the case with initial ⁄ = 0.06 in the case of , ⁄ = 10.0. This behaviour originates due to greater extent of droplet-induced flame wrinkling for large droplet cases. The relative magnitudes of flame area of a given ̃ isosurfaces for droplet and gaseous premixed flame cases follow the same qualitative trend as that of irrespective of the mean inflow velocity. The increase in the mean flow direction yields smaller values of Ω for all droplet cases due to higher probability of finding < 1 in the reacting mixture composition in these cases. Furthermore, ∫ ̇ changes in proportion to that of flame surface area for the gaseous premixed flames according to Damköhler's first hypothesis [77] and therefore remains close to unity in these cases. A marginal change is observed in = ∫ ̇ / [ , ∫ |∇ | ] values when the mean flow velocity increases. However, < 1 in the droplet cases is indicative of fuel-lean combustion because a laminar stoichiometric premixed flame is expected to have = 1. The increased propensity of fuel-lean combustion leads to a slight decrease in with increasing , ⁄ . Table 3. Time averaged normalised values of = ∫ |∇ | , , ̃ , Ω = ∫ | ̇ | / the cases with higher mean flow velocity due to the predominance of fuel-lean gaseous mixture within the flame. However, an opposite behaviour has been observed for the cases with lower mean inflow velocity where the smaller extent of flame wrinkling dominates over the effects of the thickening of the flame front. Liquid droplets interact with the flame and lead to dimples on the flame surface. Under higher mean flow velocity conditions, the PDFs of local curvature values broaden as a result of the significantly increased droplet-induced flame wrinkling. It has been demonstrated that the major part of total heat release arises from the premixed mode of combustion. However, the contribution of non-premixed combustion to the total heat release becomes stronger for larger droplets with higher mean flow velocity and this trend is relatively strong at the locations closer to the flame holder due to the greater availability of droplets interacting with the flame. The mean values of the terms which govern the statistical characteristic of displacement speed exhibit comparable magnitudes for different mean flow velocities irrespective of the distance from flame holder. The mean values of consumption speed diminish with increasing droplet diameters irrespective of the mean inflow velocity. Furthermore, under small mean inflow velocity conditions, the mean values of density weighted displacement speed have been found to be comparable for different droplet diameters at all locations. This trend changes under higher mean inflow conditions and the mean values of density weighted displacement speed increase for larger droplet cases at the locations close to the flame holder. The increased probability of fuel-lean combustion leads to a decrease in turbulent burning velocity with increasing mean inflow velocity.
The present analysis has been conducted using a single step chemistry and thus the investigation carried out in this paper needs to be extended for detailed chemistry DNS for further validation of the current findings although several previous analyses demonstrated that the flame speed statistics and flame structure obtained from simple chemistry simulations [36,65,78,79] remain at least qualitatively consistent with detailed chemistry computations [41,53,54,80] for premixed turbulent flames. Moreover, the physical mechanisms, which are responsible for the variations observed based on two mean inlet velocities in this analysis are unlikely to change for a different value of , ⁄ , and thus, the qualitative trends in response to the changes in , ⁄ can be inferred from the results presented here. However, further analyses for different mean inlet velocities will be necessary for deeper physical understanding.