Statistics of vorticity alignment with local strain rates in turbulent premixed flames

The instantaneous alignment of the vorticity vector with local principal strain rates is analysed for statistically planar turbulent premixed flames with different values of heat release parameter and global Lewis number spanning different regimes of combustion. It has been shown that the vorticity vector predominantly aligns with the intermediate principal strain rate in turbulent premixed flames, irrespective of the regime of combustion, heat release parameter and Lewis number. However, the relative alignment of vorticity with the most extensive and compressive principal strain rates changes based on the underlying combustion conditions. Detailed physical explanations are provided for the observed behaviours of vorticity alignment with local principal strain rates. It has been shown that heat release due to combustion significantly affects the alignment of vorticitywith local principal strain rates. However, the mean contribution of the vortex-stretching term in the transport equation of enstrophy remains positive for all cases considered here, irrespective of the nature of the vorticity alignment. © 2014 Elsevier Masson SAS. All rights reserved.


Introduction
The alignment of the vorticity vector with local principal strain rates is of fundamental importance for the understanding and modelling of turbulent fluid motion, as the alignment statistics directly affect the nature of the vortex-stretching mechanism [1]. It has been demonstrated in several previous studies that the vorticity vector instantaneously aligns with the intermediate eigenvector of strain rate tensor for non-reacting turbulence [2][3][4][5][6][7][8][9][10][11][12]. However, relatively limited attention was given to the analysis of alignment of vorticity with local strain rates in the case of turbulent reacting flows [13][14][15]. In many applications (e.g. Spark Ignition (SI) engines and industrial gas turbines), the fuel and oxidiser are homogeneously mixed prior to the combustion process (i.e. premixed combustion). Thus, the understanding of vorticity alignment with local principal strain rates is of fundamental interest for the development of high-fidelity models which can, in turn, contribute * Tel.: +44 0191 208 3570; fax: +44 0191 208 8600.
E-mail address: nilanjan.chakraborty@newcastle.ac.uk. to the design of new generation energy-efficient and environmentfriendly combustion devices. The analysis of Nomura and Elghobashi [13], Boratov et al. [14] and Jaberi et al. [15] concentrated on vorticity alignment with local principal strain rates for non-premixed flames where fuel and oxidiser are completely separated from each other prior to the combustion process. Recently, Hamlington et al. [16] analysed vorticity statistics in premixed combustion based on numerical solutions of reactive systems. The analysis by Nomura and Elghobashi [13] demonstrated that the vorticity vector aligns with the intermediate principal strain rate in non-premixed flames similar to non-reacting turbulent flows but vorticity in non-premixed flames shows appreciable probabilities of local alignment with the most extensive principal strain rate. The analysis by Boratov et al. [14] on non-premixed flame DNS data reveals that the extent of vorticity alignment with the most extensive principal strain rate increases in the regions where the magnitude of strain rate dominates over the vorticity magnitude. By contrast, vorticity shows preferential alignment with the intermediate principal strain rate in the regions where the vorticity magnitude dominates over the strain rate magnitude. The analysis by Jaberi et al. [15] further demonstrated that the alignment of  [15] was carried out for constant volume homogeneous turbulence but their findings were found to be qualitatively similar to the results by Nomura and Elghobashi [13] for non-premixed combustion in the presence of inhomogeneous turbulence. Moreover, Jaberi et al. [15] showed that vorticity remains mostly perpendic-ular to the most compressive principal strain rate in both reactive and non-reactive regions of non-premixed turbulent combustion. The computational analysis by Hamlington et al. [16] reveals that vorticity alignment with local principal strain rates in turbulent premixed flames is qualitatively similar to the previous findings in the context of non-premixed combustion (i.e. preferential alignment with the intermediate principal strain rate; negligible alignment with the most compressive principal strain rate and an increased alignment with the most extensive principal strain rate in the heat releasing zone). In general it has been found that the alignment of vorticity with the most extensive principal strain rate increases in the heat releasing zone in turbulent reacting flows [13][14][15][16]. Moreover, Hamlington et al. [16] demonstrated that the alignment of vorticity vector with local principal strain rates is expected to show a behaviour similar to a passive non-reacting isothermal turbulent flow when the root-mean-square (rms) velocity fluctuation u ′ remains much greater than the unstrained laminar burning velocity S L . However, the statistical behaviour of vorticity alignment with local principal strain rate changes significantly in comparison to passive turbulent conditions when u ′ remains either small than or comparable to S L . Detailed physical explanations for the u ′ /S L dependence of vorticity alignment have also been provided by Hamlington et al. [16] who indicated that the alignment of vorticity with flame normal vector plays a key role and vorticity tends to align with local flame normal vector in the reaction zone for small values of u ′ /S L , which in turn increases alignment of vorticity with the most extensive principal strain rate in the thin reaction zones regime flames with unity Lewis number Le (i.e. the ratio of thermal diffusivity to mass diffusivity). The relative sizes of flame thickness, and the Kolmogorov length scale, significantly affect the nature of combustion process in turbulent premixed combustion [17,18]. Energetic turbulent eddies cannot penetrate into the flame structure for the regime of combustion where Kolmogorov length scale remains greater than the flame thickness. This regime is commonly referred to as the corrugated flamelets regime [17]. By contrast, energetic eddies can penetrate into the preheat zone and can perturb the reaction-diffusion balance when the flame thickness remains greater than the Kolmogorov length scale but the reaction zone thickness needs to be smaller than the Kolmogorov length scale to prevent flame extinction. This regime of combustion is referred to as the thin reaction zones regime of combustion [17]. The regime of combustion is often characterised by a non-dimensional number known as the Karlovitz number Ka, which scales as Ka ∼ δ 2 th /η 2 [17]. Thus Ka < 1(Ka > 1) is associated with the corrugated flamelets (thin reaction zones) regime of combustion. Alternatively Ka can be taken to scale as Ka ∼ v 2 η /S 2 L , which suggests that the velocity induced by thermal expansion due to heat release is likely to play an important role in the corrugated flamelets regime (where Ka < 1). By contrast, turbulent velocity fluctuations are likely to mask the effects of heat release in the thin reaction zones regime (where Ka > 1). Thus it can be expected that the regime of combustion is likely to have an appreciable effect on the vorticity statistics in turbulent premixed flames. The analysis of Hamlington et al. [16] concentrated on the regimes of premixed turbulent combustion with unity Lewis number (i.e. Le = 1.0) where δ th remains greater than η (i.e. Ka > 1) but differences in vorticity alignment statistics between Ka < 1 and Ka > 1 combustion are yet to be analysed in detail. Moreover, the analysis by Hamlington et al. [16] was carried out based on implicit Large Eddy Simulations (LES) where the numerical scheme was responsible for the dissipation of sub-grid turbulent kinetic energy. In contrast to LES, which inherently allows for approximations, Direct Numerical Simulations (DNS) resolve all the relevant length and time scales of turbulence without any recourse to turbulence modelling. However, most of the relevant length scales are resolved in LES and as a result the implications of modelling inaccuracies are much less severe than in Reynolds Averaged Navier-Stokes (RANS) simulations. Thus, the vorticity alignment with local principal strain rates reported by Hamlington et al. [16] based on LES is found to be consistent with the present DNS findings for the thin reaction zones regime flames with unity Lewis number (see Section 3).
A recent DNS analysis [19] indicates that the characteristic Lewis number Le of turbulent premixed combustion has significant influences on the alignment statistics of reactive scalar gradients and thus Le is expected to have appreciable effects on vorticity alignment statistics. The effects of the ratio of flame thickness to the Kolmogorov length scale (characterised by Karlovitz number where the repeated indices in Eq. (1) and subsequent equations implicitly indicate a summation process according to the standard tensor notation. The first term on the right hand side of Eq. (1) is the vortex-stretching term. The second term on the right hand side originates from the non-zero dilatation rate (i.e. ∇ • ⃗ u = ∂u k /∂x k ) due to thermal expansion, whereas the third term denotes the molecular diffusion of vorticity. The fourth term on the right hand side depicts vorticity generation/destruction by baroclinic effects. The fifth and sixth terms originate due to density change caused by the chemical reaction, and the last term f 1 (µ) involves the term involving spatial gradients of viscosity.
Multiplying Eq. (1) by ω i yields the transport equation of Ω = ω i ω i /2, which is given by: where the penultimate term on right hand side (i.e. −µ(∂ω i /∂x j ) (∂ω i /∂x j )) is responsible for molecular dissipation of Ω, and f 2 (µ) accumulates the terms involving spatial derivatives of viscosity.
The statistical behaviour of the first term on the right hand side (i.e. ρω i ∂u i /∂x j ω j ) originates due to vortex stretching, and the statistical behaviour of this term depends on the alignment of ω i with strain rate tensor e ij = 0.5(∂u i /∂x j + ∂u j /∂x i ), which can be verified from the following expression: It is evident from Eq. (3) that the magnitude and sign of ρω i e ij ω j depend on the statistical behaviours of |cos α| , |cos β| and |cos γ |, and on the relative magnitudes of e α , e β and e γ . The alignment of vorticity with local principal strain rates will be discussed in detail in Section 3 of this paper. The enstrophy Ω is closely related to the dissipation rate of turbulent kinetic energy [20], which plays a key role in the closure of turbulent Reynolds stresses [20], and serves as an input parameter to the turbulent combustion models [21,22]. Thus it is important to analyse the vorticity alignment statistics (e.g. vortex-stretching term in Eq. (3)) for the modelling dissipation rate of turbulent kinetic energy.

Mathematical background & numerical implementation
In the present analysis, the alignment statistics of vorticity have been studied based on three-dimensional DNS data of statistically planar freely propagating turbulent premixed flames. A wellknown DNS database [23] has been used for analysing vorticity alignment in the corrugated flamelets regime. The simulations representing the thin reaction zones regime were carried out using a compressible DNS code called SENGA [24]. Ideally, both the three-dimensionality of turbulence, and the detailed structure of the flame should be accounted for in combustion DNS studies. However, accounting for detailed chemical mechanism along with three-dimensionality of turbulent flow field remains extremely computationally expensive even by contemporary standards. Therefore, the chemical aspect of combustion in this analysis is accounted for by a single-step Arrhenius-type reaction (i.e. Reactants → Products). In the context of simplified chemistry, a reaction progress variable c uniquely represents the species field. The reaction progress variable c is defined in terms of the mass in such a manner that c rises monotonically from 0 to 1 from unburned gases to fully burned gases. The reaction rate of c takes the following form, under the assumptions of single-step chemistry: In the present study, the specific heats of all species are taken to be equal and independent of temperature, and Fick's law has been used for mass diffusion. The Soret and Dufour effects are neglected. Moreover, thermo-physical properties such as, viscosity µ, thermal conductivity k and the density-weighted mass diffusivity ρD are taken to be equal for all species, and independent of temperature, following several previous studies [24][25][26][27]. The initial values of normalised root-mean-square (rms) turbulent velocity fluctuation u ′ /S L , integral length scale to flame thickness ratio l/δ th , turbulent Reynolds number Re t = ρ 0 u ′ l/µ 0 , Damköhler number   Table 1, along with the values of heat release parameter τ = (T ad − T 0 )/T 0 and Lewis number Le. The thermal flame thick- where the subscript L refers to the unstrained laminar flame quantities. Unstrained planar laminar flame simulation is performed on different grids to ensure grid-independent solution and the maximum value of was evaluated based on these laminar flame simulations. Table 1 suggests that case A represents the corrugated flamelets regime combustion (i.e. Ka < 1), whereas the other cases represent the thin reaction zones regime (i.e. Ka > 1) according to the regime diagram by Peters [17]. Among the cases considered here, the differential diffusion of heat and mass has been ignored in cases A, B and F, which is reflected in the unity value of Lewis number. The rate of mass diffusion is taken to be faster (slower) than the thermal diffusion rate in cases C-E (case G), which is reflected in the Le < 1 (Le > 1) value in these cases. In a real combustion process, different species have different values of Lewis number but often the Lewis number of the deficient reactant is taken to be the characteristic Lewis number of the underlying combustion process [28]. In several previous studies [29][30][31][32][33][34][35][36], the characteristic Lewis number was modified independently of other parameters in order to analyse the effects of differential diffusion of heat and mass in isolation. The same approach has been adopted in this analysis.
It is evident from Table 1 that the initial values of turbulent Reynolds number Re t are different but comparable for cases A and B-G. The difference between Re t for cases A and B-G (i.e. Re t = 56.7 in case A and Re t = 47.0 in cases B-G) is not large enough to induce any major Reynolds number dependence. Cases B and F are identical to each other except the value of heat release parameter τ . If a comparison between the results for cases B and F reveals that a 1.5 fold increase τ (i.e. τ = 3.0 in case B and τ = 4.5 in case F) does not significantly modify the qualitative nature of vorticity dynamics in the thin reaction zones regime (which will be substantiated later in this paper), it can be inferred that the differences in vorticity alignment statistics between cases A and B originate due to differences in the combustion regimes because the difference in τ between cases A and B (τ = 2.3 in case A and τ = 3.0 in case B) is not expected to play a major role.
There have been several previous analyses where either cases A, B and F or the simulations parameters similar to cases A, B and F have been used to demonstrate and analyse the fundamental differences in flame propagation [18], turbulent kinetic energy [37] and scalar dissipation rate [22] and Flame Surface Density [38,39] closures, scalar gradient alignment [40] and conditional velocity statistics [41] between the corrugated flamelets and thin reaction zones regimes and the same approach has been followed here.
For case A, inflow and outflow boundaries are specified in the direction of mean flame propagation and transverse directions are considered to be periodic. The inlet turbulent velocity field is specified by scanning a plane though a frozen turbulent velocity fluctuation field. The outflow boundary is specified using the Navier-Stokes Characteristic Boundary Conditions (NSCBC) outlined by Poinsot and Lele [42]. The boundaries in the direction of the mean flame propagation (i.e. negative x 1 -direction in the present configuration) for cases B-G are taken to be partially non-reflecting in nature, following NSCBC [42] technique, and transverse boundaries are considered to be periodic. According to NSCBC technique, the inviscid part of the Navier-Stokes equations is specified using a characteristic analysis based on a locally one- 3 and A 5 = p + ρu 1 a, and the wave velocities λ i are given by λ 1 = u 1 − a, λ 2 = λ 3 = λ 4 = u 1 and λ 5 = u 1 + a. The wave amplitude variations L i directed into the domain need to be specified at the partially non-reflecting boundary and the wave amplitude variations for the outgoing waves are estimated from internal solutions. For the incoming waves the following expressions for L i : where σ i s are the relevant relaxation parameters [42]. The values of ∂τ 1j /∂x 1 (for j = 1, 2 and 3), ∂q T 1 /∂x 1 and ∂q C 1 /∂x 1 are considered to be zero at the partially non-reflecting boundaries [42], where q T 1 = −k(∂T /∂x 1 ) and q C 1 = −ρD∂c/∂x 1 are the heat and mass fluxes respectively.
Interested readers are referred to Ref. [42] for a more detailed discussion on boundary conditions. In DNS, the turbulent flow field is simulated without any physical approximation, and the grid spacing remains smaller than the Kolmogorov length scale. The grid size for DNS is chosen in such a manner that the grid spacing satisfies ∆ ≤ min(η, δ th /10). For all cases considered here the grid spacing is determined by the resolution of the flame, and at least about 10 grid points are kept within δ th . As turbulence decays rapidly with the increase in kinematic viscosity in the burned gases (even in the present case due to the decrease in gas density in the burned gas), the flame-turbulence interaction is mainly characterised by the turbulent processes in the unburned gas side [25,26]. Thus, the findings based on constant thermophysical properties are likely to be qualitatively valid even when temperature dependent properties are considered. This matter has been discussed earlier in detail by Poinsot et al. [25] and Louch and Bray [26]. Standard values are taken for Zel'dovich number β Z , Prandtl number Pr and the ratio of specific heats γ G (i.e. β Z = 6.0, The turbulent velocity fluctuation is initialised by a homogeneous isotropic incompressible flow field which is generated using a pseudo-spectral method [43] following the Batchelor-Townsend spectrum [44]. The flame is initialised using unstrained planar laminar flame simulations. In case A all the convective terms are time advanced in an explicit manner using a low-storage 3rd order Runge-Kutta scheme [45]. By contrast, all the diffusive/viscous terms in case A are time-advanced implicitly using the Crank-Nicolson scheme. For cases B-G all the terms are timeadvanced in an explicit manner using a low-storage 3rd order Runge-Kutta scheme [45]. Similar procedure has been followed in several previous analyses [18,19,24,[28][29][30][31][32][33][34][35][36][37][38][39][40][41][46][47][48][49][50][51]. In all cases flame-turbulence interaction takes place under decaying turbulence. Under decaying turbulence, simulations should be carried out for t sim , where t f = l/u ′ is the initial eddy turn over time and t c = δ th /S L is the chemical time scale. The simulation in case A was run for about 4 initial eddy turn over times (t sim ∼ 4t f ∼ 4l/u ′ ), whereas cases B-G were run for a time equivalent to 3.34t f . The aforementioned simulation times remain either greater than (case A) or equal to (cases B-G) one chemical time scale t c = δ th /S L , and are comparable to several previous studies [18,19,24,[28][29][30][31][32][33][34][35][36][37][38][39][40][41][46][47][48][49][50][51]. As the typical values of eddy turn over time t f is much greater in case A than in cases B-G, the numerical treatment for case A is different in comparison to cases B-G for the purpose of computational economy.
The turbulent kinetic energy and its dissipation rate in the unburned gas ahead of the flame were not varying significantly with time when statistics were extracted. The values of u ′ /S L in the unburned reactants ahead of the flame at the time when statistics were extracted decreased by about 52% and 50% of its initial value in cases A and B-G respectively. The values of l/δ th have increased from their initial values by a factor of about 1.10 and 1.7 for cases A and cases B-G respectively, but there are still enough turbulent eddies in each side of the computational domain. All the statistics presented in Section 3 correspond to t = 4.0l/u ′ (t = 3.34l/u ′ ) for case A (cases B-G) but the qualitative nature of the statistics was found to have remained unchanged since t = 2.0 l/u ′ for all cases.

Distribution of vorticity
The distributions of vorticity magnitude mid-plane are shown in Fig. 1  (2.14, 1.44, 1.15 and 0.95) times that in case F respectively [36]. This is found to be consistent several previous analyses on the effects of non-unity Lewis number [19,[29][30][31][32][33][34][35][36]52,53]. The enhanced (reduced) burning rate in the Le < 1 (Le > 1) flames give rise strengthening (weakening) of dilatation rate and velocity jump across the flame in comparison to the corresponding unity Lewis number flame, and the effects of flame generated turbulence strengthen (weaken) with decreasing (increasing) Le. These effects were discussed elsewhere [19,35,36] in detail by this author (see Fig. 5 of Ref. [19] for the pdfs of ∂u i /∂x i for cases C, E and G. The values of (∂u i /∂x i ) and velocity jump across the flame for cases C-G are shown in Fig. 3b, e of Ref. [35]. The augmentation of flame generated turbulence with decreasing Le can be seen from Fig. 2 of Ref. [36]) and interested readers are referred to these papers for further information. It has also been demonstrated by the present author that the alignment of ∇c with local principal strain rates is also significantly affected by the enhanced heat release effects at small values of Le and interested readers are referred to Ref. [19] for more information. √ ω i ω i × δ th /S L in the central x 1 − x 2 plane when the statistics were extracted for (a-g) cases A-G (see Table 1 for attributes of cases A-G). The white lines indicate reaction progress variable contours from = 0.1 to 0.9 from left to right in steps of 0.1. Jaberi et al. [15] (see Fig. 17 in Ref. [15]) and thus are not shown here for the sake of brevity. It can be seen from Fig. 1a that there is a large separation between the vorticity magnitudes at the inlet boundary and within the flame in case A. In case A the turbulent velocity fluctuation values are much smaller than other cases whereas integral length scale is much larger than other cases (see Table 1). The dilatation rate ∂u i /∂x i can be scaled as (∂u i /∂x i ) ∼ τ S L /δ th ∼ τ Da(u ′ /λ)/Re 1/2 t [18,19,40] and thus dilatation rate dominates over turbulent straining (which can be taken to scale with u ′ /λ [20]) in case A due to large values of Da. This is reflected in the preferential alignment of ∇c with e α , which was demonstrated and discussed elsewhere [40]. The strong dilatation rate field suppresses the vorticity magnitude (see the second term of Eq. (2)) as the flow approaches towards the flame. However, there are localised spots where vorticity magnitude increases across the flame due to large velocity gradients induced by strong velocity jump across the flame (see Ref. [36]) in case A. It is also worth noting that Ka ∼ δ 2 th /η 2 < 1 for case A, which indicates that this flame is not expected to encounter energetic turbulent eddies, which can be substantiated from Fig. 1. By contrast,   The flames with Le ≈ 1.0 within the thin reaction zones regime (i.e. cases B, E-G) do not impart much influence on the vorticity field in the upstream of the flame. As the initial turbulent velocity field is the same in cases B, and E-G, and the thermo-chemistry is similar as well, they exhibit almost identical reaction progress variable and vorticity fields. However, the extent of wrinkling and burning increases with decreasing Le in the Le ≪ 1.0 flames (e.g. cases C and D) due to thermo-diffusive instabilities [29][30][31][32][33][34][35][36]52,53]. The strong velocity jump across the flame and high dilatation rate magnitude due to high rate of heat release alter the vorticity field significantly in cases C and D in comparison to cases B, E-G. The flame-generated turbulence locally augments vorticity magnitude across the flame in the Le ≪ 1 flames (e.g. cases C and D). The effects of turbulent combustion regime (i.e. Ka) and Lewis number Le on flame-generated turbulence (e.g. turbulent kinetic energy generation) have been analysed elsewhere [36,37], and will not be repeated here for the sake of conciseness.

Vorticity alignment with local principal strain rates
In turbulent premixed combustion the reacting scalar field is often characterised in terms of c and thus it is useful to present the vorticity alignment statistics for different values of c across the flame. In the context of single step chemistry the reaction zone is confined to the region given by 0.7 ≤ c ≤ 0.9 [18,32]. Thus the effects of dilatation rate and heat release are particularly strong in the region given by 0.7 ≤ c ≤ 0.9 [19]. As the alignment of the vorticity vector with local strain rates can be characterised by the magnitude of the cosines of the angles α, β and γ (see Eq. (3)

Difference in vorticity alignment behaviour between the different regimes of combustion for unity Lewis number flames
The vorticity vector aligns predominantly with the intermediate principal strain rate e β towards the unburned gas side of the flame (e.g. c = 0.1), irrespective of the heat release parameter τ , and the regime of combustion. This is consistent with previous findings in the context of non-reacting [2][3][4][5][6][7][8][9][10][11][12] and reacting [13][14][15][16] turbulent flows. However, the alignment of vorticity vector with e α and e γ changes within the flame and the vorticity alignment with local principal strain rates differs from one case to the other. The vorticity vector shows a significant extent of alignment with the most compressive principal strain rate e γ for the major part of the flame, except the unburned gas side in the corrugated flamelets regime flame (i.e. case A) considered here, and this tendency increases towards the heat releasing zone of the flame (i.e. 0.7 ≤ c ≤ 0.9). The vorticity vector does not align with the most extensive principal strain rate e α throughout the flame in case A. This preferential alignment of vorticity in case A with e γ is contrary to the findings in the context of both reacting [13][14][15][16] and nonreacting flows [2][3][4][5][6][7][8][9][10][11][12]. It is worth noting that Hamlington et al. [16] demonstrated that the alignment of vorticity with local principal strain rates exhibits deviation from passive turbulent flow conditions (i.e. predominant vorticity alignment with e β ) for small values of u ′ /S L , which is qualitatively consistent with the findings of Hamlington et al. [16]. However, Hamlington et al. [16] did not report preferential alignment of vorticity with e γ as in case A, but the values of Da and Ka for case A are significantly different from the values considered in Ref. [16].
The vorticity vector preferentially aligns with e β throughout the flame in cases B and F, which represent the thin reaction zones regime combustion but in case F vorticity shows a relatively higher level of alignment with the most extensive principal strain rate in the heat releasing zone (i.e. 0.7 ≤ c ≤ 0.9) than in case B. The analysis of Hamlington et al. [16] demonstrated that the alignment of vorticity with e α decreases with weakening of heat release effects for unity Lewis number thin reaction zones regime flames., which is responsible for smaller extent of vorticity alignment with e α in case B than in case F due to smaller value of τ in case B.
Jimenez [4] showed that vorticity tubes with high vorticity magnitude induce velocity field in such a manner that the most extensive and compressive principal strain rates occur on the planes orthogonal to the plane in which these vortex tubes are placed. Thus, the high vorticity regions are preferentially aligned with the intermediate principal strain rate e β in non-reacting turbulent flows. This behaviour has been observed in all the flames in both unburned and burned gas sides of the flame for all cases. Moreover, in case B, where combustion is taking place in the thin reaction zones regime (i.e. v η ≫ S L ) for a small value of τ , the presence of the flame does not alter the vorticity alignment, in comparison to the familiar behaviour obtained in the context of non-reacting turbulent flows (i.e. preferential alignment of vorticity with e β ).
It can be seen from Figs. 2-4 that the variation of τ between cases B and F does not have any major influence on the qualitative behaviour of the alignment of vorticity with local principal strain rates in the thin reaction zones regime. By contrast, the vorticity alignments with local principal strain rates for cases A and B are markedly qualitatively different although the ratio of τ between cases B and A is comparable to the ratio of τ between cases F and B. As the turbulent Reynolds number Re t for cases A, B and F are comparable to each other (Re t in case A is 56.7, whereas this value is 47 in cases B and F), it is unlikely that this small variation of Re t is responsible for the differences in vorticity alignment between cases A and B, F, and the observed differences in vorticity alignment statistics between these cases arise due to different combustion regimes (i.e. case A represents a typical corrugated flamelets regime flame, whereas cases B and F represent typical thin reaction zones regime flames).

Difference in vorticity alignment behaviour in response to global Lewis number Le
The alignment of vorticity with principal strain rates for the decreases from the unburned gas side to the heat releasing zone of the flame. The extent of vorticity alignment with e α was shown to increase in the heat releasing zone in the previous analyses on turbulent non-premixed [13][14][15] and premixed [16] combustion. Similar behaviour has been observed for the unity Lewis number flames within the thin reaction zones regime (i.e. cases B and F).
It was shown in Ref. [19] that ∇c aligns with e α in the reaction zone in cases E-G and a greater level of vorticity alignment with e α is consistent with the findings of Hamlington et al. [16], which demonstrated that vorticity tends to show significant alignment with e α in the regions of intense heat release where vorticity aligns with ∇c in the thin reaction zones regime flames with unity Lewis number. However, the alignment of vorticity with e α (e γ ) decreases (increases) in the heat releasing zone for the high Damköhler number corrugated flamelets regime flame (i.e. case A) and in the thin reaction zones regime flames with Le ≪ 1 (i.e. cases C and D).
The alignment of vorticity with e γ was not observed in the previous analyses [2][3][4][5][6][7][8][9][10][11][12][13][14][15][16]. The analysis of Hamlington et al. [16] ignored the effects of baroclinic torque on vorticity alignment behaviour under the conditions of intense heat release (see Eq. (25) in Ref. [16]) and it will be demonstrated later in the paper that baroclinic torque might play a key role in the vorticity alignment statistics in Le ≪ 1 flames (e.g. cases C and D) and also in the flame representing the corrugated flamelets regime combustion (i.e. case A). A comparison between Figs. 2-4 reveals that the extent of vorticity alignment with e α decreases with decreasing Le and this trend is particularly prevalent for flames with Le ≪ 1 (e.g. cases C and D). By contrast, the extent of vorticity alignment with e β and e γ increases with decreasing Le in the reaction zone 0.7 ≤ c ≤ 0.9 (see cases C-G in Figs. 3 and 4) but no major difference is observed in the preheat zone (i.e. c ≤ 0.5) in response to Le. The physical explanations for the aforementioned vorticity alignment behaviour will be provided next in this paper.

Physical explanations for the observed vorticity alignment
Figs. 2-4 suggest that the regime of combustion and Lewis number have significant influences on the alignment of the vorticity vector with local principal strain rates. The most extensive principal strain rate e α in turbulent premixed flames can be scaled with the strain rate induced by chemical reaction a chem . The strain rate induced by flame normal acceleration due to thermal expansion is referred to as the strain rate induced by chemical reaction a chem in this analysis. In an unstrained laminar flame the velocity jump across the flame (i.e. difference between velocities in products and reactants) is given by τ S L . As this velocity difference takes place over the flame thickness δ th , the strain rate induced by chemical heat release can be scaled as a chem ∼ τ S L /δ th [19,40]. It can be shown for low Mach number unity Lewis number flames ∂u i /∂x i can be expressed as (∂u i /∂x i ) = τ (ρS d ) |∇c| /ρ 0 [18,32] where S d = (Dc/Dt)/ |∇c| is the displacement speed (i.e. the speed with which the flame surface moves normal to itself with respect to an initially coincident material surface). Scaling ρS d |∇c| by ρ 0 S L /δ th enables one to scale (∂u i /∂x i ) in the following manner (∂u i /∂x i ) ∼ τ S L /δ th . It has been shown in a previous analysis [19] that a similar scaling can also be applied to turbulent nonunity Lewis number flames in the following manner a chem ∼ ϕ 1 (Le)τ S L /δ th and (∂u i /∂x i ) ∼ ϕ 2 (Le)τ S L /δ th where the functions ϕ 1 (Le) and ϕ 2 (Le) increase with decreasing Le to account for strengthening of heat release effects with decreasing Le [19].
The strain rate induced by turbulence a turb can be scaled with tangential strain rate a turb ∼ a T = (δ ij − N i N j )∂u i /∂x j . In turbulent premixed flames the chemical strain rate can be scaled using the dilatation rate ∂u i /∂x i as: a chem ∼ ∂u i /∂x i . Thus the ratio between a chem to a turb can be scaled as: a chem /a turb ∼ (∂u i /∂x i )/a T . The variations of ⟨∂u i /∂x i ⟩/⟨a T ⟩ with c across the flame are shown in Fig. 5 for all cases considered here. It is evident from Fig. 5 that ⟨∂u i /∂x i ⟩/⟨a T ⟩ assumes values much greater than unity in cases A, C and D for the major part of the flame. In cases B, E-G the quantity ⟨∂u i /∂x i ⟩/⟨a T ⟩ assumes much smaller values than the values obtained in cases A, C and D. A comparison between cases B and F reveals an increase in τ yields a greater value of ⟨∂u i /∂x i ⟩/⟨a T ⟩. This suggests that the differences between ⟨∂u i /∂x i ⟩/⟨a T ⟩ values between cases A and B arise principally due to the differences in the combustion regimes. Fig. 5 suggests that ⟨∂u i /∂x i ⟩/⟨a T ⟩ assumes higher values in the corrugated flamelets regime flame than in the thin reaction zones regime flames with comparable values of Re t (see cases A, B and F). As Re t scales as: Re t ∼ Da 2 Ka 2 in the Le = 1.0 flames [17], it can be inferred by comparing the values of ⟨∂u i /∂x i ⟩/⟨a T ⟩ in cases A, B and F that ⟨∂u i /∂x i ⟩/⟨a T ⟩ strengthens with increasing Da (as Re t remains comparable for cases A, B and F but Ka < 1 in case A and Ka > 1 in cases B and F). The strengthening of ⟨∂u i /∂x i ⟩/⟨a T ⟩ in the corrugated flamelets regime is consistent with earlier findings of Chakraborty and Swaminathan [40] in the context of alignment statistics of ∇c. Fig. 5 further indicates that the peak value of ⟨∂u i /∂x i ⟩/⟨a T ⟩ decreases with increasing Le, which is found to be consistent with previous findings of Chakraborty et al. [19]. A comparison between Figs. 2-5 reveals that the cases with large values of ⟨∂u i /∂x i ⟩/⟨a T ⟩ (e.g. cases A, C and D) exhibit mostly perpendicular alignment of vorticity with e α . However, a one-toone correspondence between the increase in ⟨∂u i /∂x i ⟩/⟨a T ⟩ and decreasing collinear alignment of vorticity with e α on a given c isosurface may not always be achieved because Fig. 5 shows the ratio of mean quantities, whereas the pdfs in Fig. 2 show the instantaneous behaviour. In non-reacting turbulent flows vorticity does not show any preferential (either collinear or perpendicular) alignment with e α . Moreover, a comparison between Figs. 2 and 5 indicates that the vorticity alignment with e α shows similar qualitative behaviour as that of non-reacting turbulent flows for The effects of dilatation rate on vorticity alignment can further be elucidated by the variation of the mean values of |cos α|, |cos β| and |cos γ | conditional on dilatation rate (∂u i /∂x i ) within the progress variable range given by 0.1 ≤ c ≤ 0.9, which are shown in Fig. 6a-g for cases A-G respectively. Fig. 6a-g show that vorticity alignment with e α decreases in the regions of high magnitude of (∂u i /∂x i ) although this trend is relatively weak for the Le = 1.2 flame where the rate of burning is the weakest [19,[29][30][31][32][33][34][35][36]52,53]. By contrast, the alignment of vorticity with e γ increases at the regions of high dilatation rate. The extent of vorticity alignment with e β also decreases with increasing (∂u i /∂x i ) for Le = 0.8 and 1.0 flames (i.e. cases A, B, E and F). However, no systematic trend has been observed between |cos β| and (∂u i /∂x i ) for Le = 0.34, 0.6 and 1.2 flames (i.e. cases C, D and G).
The mean values of   cos α p   ,   cos β p   and   cos γ p   conditional on dilatation rate (∂u i /∂x i ) within the progress variable range given by 0.1 ≤ c ≤ 0.9 are shown in Fig. 7a-g for cases A-G respectively. It is evident from Fig. 7a-g that ∇p remains imperfectly aligned with both e β and e γ , and the extent of ∇p alignment with e α increases for large values of (∂u i /∂x i ). The density gradient ∇ρ is closely related to ∇c (i.e. |∇ρ| ∼ τ ρ 2 |∇c| /ρ 0 and an equality holds for low Mach number unity Lewis number flames). Previous analyses [19,40,49,54] indicated that ∇c aligns preferentially e α in high Da unity Lewis number and low Da (i.e. Da < 1) subunity Lewis number Le ≪ 1 flames (e.g. cases A, C and D) where the effects of (∂u i /∂x i ) are particularly strong, so ∇ρ also predominantly aligns with e α in these cases. This suggests that the cross product between ∇p and ∇ρ is expected to have significant contributions in the directions of e β and e γ in high Da unity Lewis number and low Da (i.e. Da < 1) sub-unity Lewis number Le ≪ 1 flames (e.g. cases A, C and D), whereas a weak contribution is expected in the direction of e α . Thus, the effects of the baroclinic torque (i.e. ε ijk (∂ρ/∂x j )(∂p/∂x k )(1/ρ)) are likely to be weak in the direction of e α in cases A, C and D. A comparison between Figs. 2-4 and 6 also reveals that vorticity aligns significantly with e β and e γ in the regions of flame where the vorticity alignment with e α is weak in high Da unity Lewis number and low Da (i.e. Da < 1) subunity Lewis number Le ≪ 1 flames (e.g. cases A, C and D).
It is worth noting that the intensity of heat release increases with decreasing Lewis number which in turn increases the magnitude of dilatation rate with decreasing Le (see Fig. 5 in Ref. [19] and Fig. 3b in Ref. [36]). However, dilatation rate is one of the major manifestations of thermal expansion and flame normal acceleration also strengthens with decreasing Le. This strengthening of flame normal acceleration and dilatation rate also acts to aid the flame-generated turbulence in premixed flames and this effect is particularly strong in Le ≪ 1 flames (see Refs. [35,36]). Thus, the qualitative behaviours of the variation of |cos α|, |cos β| and |cos γ | (   cos α p   ,   cos β p   and   cos γ p   ) conditional on dilatation rate (∂u i /∂x i ) for cases C-G are not expected to be identical because dilatation rate also affects the flame-generated turbulence and turbulence intensity, which in turn have significant influences of the alignment behaviours of vorticity and pressure gradient with local principal strain rates. Fig. 1 of Ref. [6] showed imperfect vorticity alignment between ∇p and e β (and also with e γ to some extent) for non-reacting turbulence and the present results are in qualitative agreement with that finding but the extent of alignment of ∇p with the eigenvectors corresponding to e α , e β and e γ in turbulent premixed flames are significantly different in comparison to non-reacting turbulent flows. The analysis of Jaberi et al. [15] demonstrated that ∇ρ and ∇p show predominant collinear alignment in reacting flows without heat release but ∇ρ and ∇p tend to show high probability of mutually perpendicular alignment in the presence of heat release. As the extents of ∇c and ∇p alignments with e α increase for large values of (∂u i /∂x i ) (see Fig. 7 and Ref. [54]), both ∇ρ ∼ (−τ ρ 2 ∇c/ρ 0 ) and ∇p show greater extent of collinear alignment in the premixed flames analysed here than in the nonpremixed flames with heat release analysed by Jaberi et al. [15]. Lee et al. [55] reported that the extent of vorticity alignment with all the principal strain rates (i.e. e α , e β and e γ ) decreases with increasing (∂u i /∂x i ) for non-reacting compressible flows, which is in contrast to the present observations. In non-reacting compressible flows ∇ρ does not have any preferential alignment with e α unlike turbulent premixed flames. Thus the baroclinic torque does not remain weak in a preferred direction at high magnitudes of (∂u i /∂x i ) in non-reacting compressible flows. It is worth noting that dilatation rate effects on vorticity alignment has been observed here for small values of Mach number where the dilatation rate effects principally originate due to the thermal expansion as a result of highly exothermic chemical reaction. Moreover, dilatation rate is predominantly positive in turbulent premixed flames due to thermal expansion in contrast to nonreacting turbulent compressible flows (see Fig. 6a-g).
The pdfs of e β /e α for all cases are shown in Fig. 8 whereas the pdfs of the ratio of the most compressive principal strain rate to the most extensive principal strain rate (i.e. e γ /e α ) for all cases are shown in Fig. 9.   cos α p   ,   cos β p   and   cos γ p   ) conditional on normalised dilatation rate ∂u i /∂x i × δ th /S L in the region corresponding to 0.1 ≤ c ≤ 0.9 for cases (a-g) A-G.
finding e α ≫   e β   and e α ≫   e γ   are significant for high Da unity Lewis number and low Da (i.e. Da < 1) sub-unity Lewis number Le ≪ 1 flames (i.e. cases A, C and D), which indicates a significant probability of finding (∂u i /∂x i ) ≈ e α (compare between Figs. 10 and 11) in these cases. The probability of (∂u i /∂x i ) ≈ e α decreases with increasing (decreasing) Le(Da). This can further be substantiated from Fig. 12 where the pdfs of (∂u i /∂x i )/e α are shown. Fig. 12 shows that the probability of finding (∂u i /∂x i ) ≈ e α is relatively smaller in low Damköhler number (i.e. Da < 1) flames with Le ≈ 1.0 (i.e. cases B, E-G) than in high Da unity Lewis number   e β   /e α ) on the reaction progress variable c = 0.1, 0.3, 0.5, 0.7 and 0.9 isosurfaces, for cases (a-g) A-G.
(i.e. ε ijk (∂ρ/∂x j )(∂p/∂x k )(1/ρ)) effects in Eq. (1) are likely to be weak in the direction of e α . This suggests that there is no major vorticity generating mechanism in the direction of e α in cases A, C and D and thus vorticity does not align with e α for the major portion of the flame in these cases. Therefore, vorticity alignment with e α remains unstable in cases A, C and D in the region of the flame where the probability of finding (∂u i /∂x i ) ≈ e α is significant. Fig. 12 shows that the probability of finding (∂u i /∂x i ) ≈ e α is overwhelming throughout the flame in high Da unity Lewis number combustion (e.g. case A) and thus vorticity does not align with e α a b c d e f g Fig. 11. Pdfs of the ratio of the magnitude of the most compressive principal strain rate to the most extensive principal strain rate (i.e.   e γ   /e α ) on the reaction progress variable c = 0.1, 0.3, 0.5, 0.7 and 0.9 isosurfaces, for cases (a-g) A-G. The baroclinic torque plays a major contributor to the vorticity generation in reacting flows [15,16]. As the contribution of ρ | ⃗ ω| (e α − ∇ • ⃗ u) cos α ≈ ρ | ⃗ ω| (e α − e α ) cos α ≈ 0 and the baroclinic torque (i.e. ε ijk (∂ρ/∂x j )(∂p/∂x k )(1/ρ)) is weak in the direction of e α , only viscous diffusion of vorticity remains major contributor in this direction in cases A, C and D. This further indicates that vorticity alignment is unstable in the direction of e α in cases A, C and D. As both ∇p and ∇ρ assume large values in Le ≪ 1 flames (e.g. cases C and D) [19,35,36], the effects of baroclinic torque in the directions of e β and e γ strengthen with Table 2 Ratio of the magnitudes of the mean values of e β cos 2 β to e γ cos 2 γ (i.e.   ⟨e β cos 2 β⟩   /   ⟨e γ cos 2 γ ⟩   ) conditional on c = 0.1, 0.3, 0.5, 0.7 and 0.9 isosurfaces for cases A-G.
Case t /Da < 1) will behave in a similar manner to case A, as Re t ≫ 1 in most practical engineering applications. However, it is theoretically possible to obtain Da ≫ 1 even in the thin reaction zones regime combustion (i.e. Ka ∼ Re where the vorticity alignment is expected to qualitatively similar to that of case A. The qualitative nature of vorticity alignment in low Damköhler number (i.e. Da < 1) thin reaction zones regime combustion depends significantly on the global Lewis number Le, as demonstrated by cases C-G considered here. The effects of Lewis number on heat release rate and dilatation rate ∂u i /∂x i statistics are dependent on relative balances of thermal and mass diffusion and thus are mostly unaffected by the regime of combustion. As a result, the extent of alignment of vorticity with e α (e β and e γ ) is expected to decrease (increase) with decreasing Le in the corrugated flamelets regime as well.
Figs. 2-4 clearly demonstrate that the vorticity alignment behaviour changes with c for all the cases. In cases A, C and D the dilatation rate effects in the heat releasing region of the flame (i.e. 0.7 ≤ c ≤ 0.9) are strong enough to alter the background turbulent flow field. However, the effects of dilatation rate are not strong for both unburned and burned gas side of the flame [19] and thus the vorticity alignment behaviour shows similarity to the well-known non-reacting turbulent flow alignment behaviour (i.e. preferential alignment with e β ) on both sides of the flame. In the thin reaction zones regime flames with Le ≈ 1.0 (e.g. cases B, E-G) v η ≫ S L and thus these flames imparts relatively less influence on the background turbulent fluid motion in comparison to Da > 1 unity Lewis number combustion (e.g. case A). The effects of heat release are strong enough to influence the background fluid motion in Le ≪ 1 flames even in the thin reaction zones regime.
The pdfs of e β /e α in Fig. 8 reveal that e β /e α predominantly assumes positive values throughout the flame, indicating predominantly positive values of e β , as e α is the most extensive (i.e. positive) principal strain rate. Eq. (3) suggests that the alignment of vorticity with e α (e γ ) tends to yield positive (negative) values of Λ. The ratio of the mean values of the magnitudes of e β cos 2 β to e γ cos 2 γ (i.e.   ⟨e β cos 2 β⟩   /   ⟨e γ cos 2 γ ⟩   ) on five different c isosurfaces across the flame are presented in Table 2, which demonstrates that the magnitude of the mean contribution of e β cos 2 β supersedes the magnitude of mean contribution of e γ cos 2 γ for the major portion of the flame ensuring a positive mean contribution of Λ throughout the flame irrespective of the nature of vorticity alignment with local principal strain rates. This indicates that the mean contribution of the vortex-stretching term Λ in the transport equation of Ω remains positive, irrespective of the nature of vorticity alignment with local principal strain rates.

Conclusions
A three-dimensional simplified chemistry DNS database of statistically planar turbulent premixed flames, with a range of different values of Karlovitz number, heat release parameter and Lewis number, has been used to analyse the alignment of the vorticity vector with local principal strain rates. It has been found that vorticity aligns with the intermediate principal strain rate in all cases considered here, which is consistent with previous analyses for both non-reacting and reacting flows [2][3][4][5][6][7][8][9][10][11][12][13][14][15][16].
Moreover, vorticity vector also shows non-negligible alignment with the most extensive principal strain rate and this alignment increases in the heat releasing zone for high values of Karlovitz number (i.e. Ka > 1) for flames with Le ≈ 1.0, which is consistent with previous analyses on turbulent reacting flows [13][14][15][16]. However, it has been found that the vorticity alignment with the most extensive principal strain rate is unstable for flames where the probability of finding (∂u i /∂x i ) ≈ e α is significant, which can be realised for small values of Le and also for large (small) values of Da (Ka) in case of unity Lewis number flames. The absence of vorticity alignment with e α , in low Le and small Karlovitz number unity Lewis number flames, and alignment of vorticity with e γ are in contrast to previous findings [13][14][15][16]. The intermediate principal strain rate predominantly assumes positive values, which makes the mean contribution of the vortex-stretching term in the transport equation of Ω = ω i ω i /2 to be positive for all cases considered here. The parameters such as the regime of combustion (which is characterised by Da and Ka), Lewis number Le and heat release parameter τ have major influences on the strength of the dilatation rate field in comparison to turbulent straining. The effects of dilatation rate on vorticity alignment statistics are significant in turbulent premixed flames even for small values of Mach number. The alignment of vorticity with the most extensive and intermediate principal strain rates decrease with increasing dilatation rate within the flame, which was observed in compressible non-reacting turbulent flows at much higher Mach number [55]. By contrast, the extent of alignment of vorticity with the most compressive principal strain rate increases at high dilatation rate locations in premixed turbulent flames unlike non-reacting compressible flows [55]. The dilatation rate effects strengthen significantly in comparison to the effects of turbulent straining with increasing Da and τ , and with decreasing Le. It has been found that there is a large probability of finding (∂u i /∂x i ) ≈ e α for large (small) values of Da (Le), and under this situation vorticity ceases to align with e α . One needs strong chemical reaction to have conditions given by e α ≫   e β   and e α ≫   e γ   so the vorticity alignment observed in high Da and Le ≪ 1 flames are unlikely to be obtained for compressible non-reacting turbulent flows. Moreover, (∂u k /∂x k ) is identically zero for incompressible flows so the observed change in vorticity alignment with e α cannot be seen in turbulent incompressible flows.
Hamlington et al. [16] analysed the effects of turbulence intensity on vorticity alignment statistics for unity Lewis number flames in the thin reaction zones regime. Given the similarity between the present findings and the results reported in Hamlington et al. [16], it can be expected that major qualitative differences in terms of the effects of turbulence intensity on vorticity alignment statistics are unlikely to be observed in unity Lewis number thin reaction zones regime flames. However, the vorticity alignment statistics in the corrugated flamelets regime and non-unity Lewis number flames show marked difference in comparison to the results reported in Hamlington et al. [16]. Thus, it necessitates further analyses of the effects of turbulence intensity of vorticity alignment with local principal strain rates in flames with Le ̸ = 1 and for flames representing the corrugated flamelets regime (i.e. Ka < 1). These analyses are beyond the scope of present study but will form the basis for future investigations. Moreover, it should also be noted that the effects of detailed chemistry and transport are not considered in the present analysis, and therefore future research in these directions is necessary for a more comprehensive analysis of vorticity alignment characteristics.