To better understand the two-way coupling between turbulence and chemistry, the changes in turbulence characteristics through a premixed flame are investigated. Specifically, this study focuses on vorticity, ** ω**, which is characteristic of the smallest length and time scales of turbulence, analyzing its behavior within and across high Karlovitz number (Ka) premixed flames. This is accomplished through a series of direct numerical simulations (DNS) of premixed

*n*-heptane/air flames, modeled with a 35-species finite-rate chemical mechanism, whose conditions span a wide range of unburnt Karlovitz numbers and flame density ratios. The behavior of the terms in the enstrophy,

*ω*

^{2}=

**⋅**

*ω***, transport equation is analyzed, and a scaling is proposed for each term. The resulting normalized enstrophy transport equation involves only a small set of parameters. Specifically, the theoretical analysis and DNS results support that, at high Karlovitz number, enstrophy transport obtains a balance of the viscous dissipation and production/vortex stretching terms. It is shown that, as a result, vorticity scales in the same manner as in homogeneous, isotropic turbulence within and across the flame, namely, scaling with the inverse of the Kolmogorov time scale,**

*ω**τ*

_{η}. As

*τ*

_{η}is a function only of the viscosity and dissipation rate, this work supports the validity of Kolmogorov’s first similarity hypothesis in premixed turbulent flames for sufficiently high

*Ka*numbers. Results are unaffected by the transport model, chemical model, turbulent Reynolds number, and finally the physical configuration.

## I. INTRODUCTION

Turbulent combustion involves the tight two-way coupling between turbulence and chemistry.^{1–4} This includes the individual yet connected processes of turbulence impacting the flame structure^{4–8} and the flame altering the turbulence characteristics;^{9–12} each process is integral to understanding and predicting the behavior of turbulent combustion.^{9} In premixed turbulent combustion, the coupling of the flame and turbulence has been observed to vary with the ratio of the flame time scale (*τ _{F}*) to that of the smallest turbulent eddies (

*τ*

_{η}),

^{4–6,11–14}defined as the Karlovitz number (

*Ka*),

with *τ*_{η} being the Kolmogorov time scale *τ*_{η} ≡ (*ν*/*ϵ*)^{1/2}, where *ν* is the kinematic viscosity and *ϵ* is the dissipation rate. *τ _{F}* is evaluated as

*τ*=

_{F}*l*/

_{F}*S*with

_{L}*S*the laminar flame speed and

_{L}*l*the laminar flame thickness, defined here as

_{F}*l*= (

_{F}*T*−

_{b}*T*)/|∇

_{u}*T*|

_{max}. It is common to evaluate the Karlovitz number using

*τ*

_{η}in the unburnt flow; this quantity is referred to here as the unburnt Karlovitz number (

*Ka*). Previous studies have also shown the importance of the local Karlovitz number (as opposed to

_{u}*Ka*) in determining the internal flame structure.

_{u}^{15}As the Karlovitz number is controlled by the smallest turbulent scales, to understand the effects on the flame by turbulence, it is important to describe the behavior of these turbulent scales within the flame. More fundamentally, the significance of

*τ*

_{η}within the Karlovitz number relies on the assumption that the smallest turbulent scales depend on

*ν*and

*ϵ*alone (Kolmogorov’s first similarity hypothesis), which has not been tested within premixed flames. This study focuses on the behavior of the smallest turbulent scales within high

*Ka*premixed flames.

To study the effects of the flame on the incoming turbulence, previous studies have considered various quantities. Many focused on quantities related to the turbulent kinetic energy (TKE) transport equation, often for modeling purposes.^{14,16–19} Though studying the TKE provides valuable insight into the integral scales, unfortunately, it does not describe the behavior of the smallest scales. Other studies focused on the rate of strain tensor, *S* = 1/2(∇** u** + ∇

*u*^{T}), for its relevance to the transport equation of the reaction progress variable scalar dissipation rate

^{20}and its relation to flame stretch.

^{20–23}The rate of strain tensor is often studied for the alignment of its principle axes with various quantities such as the flame normal and vorticity. Used in this way, the rate of strain tensor provides detailed information on the dynamics of the flame turbulence interaction, but does not directly provide information about the smallest turbulence scales. Study of the smallest turbulent scales can be accomplished through the energy spectrum;

^{24}however, evaluating the energy spectrum along a line in physical space, which crosses different regions (for instance burnt and unburnt) of a curved, turbulent flame, does not precisely measure the evolution of the energy spectrum as a function of progress through the flame. In order to study the behavior of the smallest turbulent scales, vorticity and the terms in its transport equation may be investigated. Previous studies focused on the vorticity vector,

**, for the importance of the vortex-stretching mechanism in the energy cascade and its appearance in the scalar gradient transport equation.**

*ω*^{11–13,25}Vorticity is also studied due to its relation to the smallest turbulent scales, as, like dissipation, it has been shown to scale with the smallest turbulent eddies in homogeneous, isotropic turbulence.

^{26–28}Studies on vorticity often consider enstrophy,

*ω*

^{2}=

**⋅**

*ω***, and its transport equation, where each term is associated with a specific physical process: vortex stretching/production, dilatation, baroclinic torque, and viscous dissipation. For these reasons, enstrophy provides insight into the behavior of the smallest turbulent scales along with the processes affecting their behavior.**

*ω*Several relevant conclusions have been made in previous work with respect to the evolution of vorticity.^{11,12,25,29} Chakraborty^{11} analyzed direct numerical simulations (DNS) using one-step chemistry with $ Ka u \u2217 $ up to 13 for the alignment of vorticity with the principle axes of the local strain rate tensor, where $ Ka u \u2217 $ is the unburnt Karlovitz number defined as $ Ka u \u2217 = ( l F / l ) 1 / 2 ( u o / S L ) 3 / 2 $ for which *l* and *u _{o}* are the integral length and velocity scales, respectively. (The Karlovitz number is often reported using $ Ka u \u2217 $, which is obtained from

*Ka*by assuming

_{u}*ν*=

*S*.) They also found the vortex-stretching term is on average positive even with different flame density ratios and Lewis numbers. Lipatnikov

_{L}l_{F}*et al.*

^{25}considered DNS at low unburnt Karlovitz numbers ($ Ka u \u2217 =.2$–.3). They found that dilatation and baroclinic torque are important in the transport of enstrophy and observed that large density ratios resulted in an average production of vorticity in the flame. The only relevant study at high Karlovitz number is attributed to Hamlington

*et al.*

^{12}who performed several simulations of high

*Ka*($ Ka u \u2217 =3$–125) premixed H

_{u}_{2}-air flames varying the turbulence intensity. These simulations relied on numerical viscosity using an implicit large eddy simulation (ILES) framework

^{30}and employed one-step chemistry with unity Lewis number assumption. In their work, the vorticity magnitude was observed to be reduced by heat release for low turbulence intensities, while at high turbulence intensities, the flame weakly affected the vorticity magnitude. Recently, Poludnenko

^{29}discussed the magnitude of terms in the vorticity equation for moderately high values of the Karlovitz number (

*Ka*= 7–30) from simulations of H

_{u}_{2}-air premixed turbulent flames also relying on numerical viscosity in an ILES framework and using one-step chemistry with unity Lewis number assumption. It was observed that for the higher value of

*Ka*(

_{u}*Ka*= 30), vorticity production was largely isotropic and had a similar magnitude through the flame as in the reactants, while at a lower Karlovitz number (

_{u}*Ka*= 7), anisotropic vorticity production and baroclinic torque were observed to play a dominant role in the transport of vorticity. As viscosity is critical to the behavior of the smallest turbulent scales, it is relevant to consider if these observations of enstrophy in simulations relying on numerical viscosity are impacted by temperature dependent molecular viscosity, which increases through a flame.

_{u}Despite these previous contributions, it still remains unclear how the flame affects enstrophy at high *Ka _{u}*, in what manner do terms in the enstrophy transport equation vary with parameters such as

*Ka*and the flame density ratio, and under what conditions Kolmogorov’s first similarity hypothesis is valid within premixed turbulent flames. Additionally, many of the previous studies are numerical and employ different models for the chemistry and species transport. As noted by Chakraborty,

_{u}^{11}it is unclear if the models have different effects on the observed results. For example, one-step chemical models assume a single chemical pathway, whereas many engine-relevant applications use complex hydrocarbon fuels with many chemical pathways. At high

*Ka*, when the turbulence is expected to disrupt the complex internal structure of the flame, the need to retain this structure in order to capture accurately the effects of the flame on the turbulence has not been studied. This knowledge is critical as it would support selecting models with minimal yet sufficient detail in order to numerically investigate the effects of the flame on turbulence accurately and efficiently. Finally, as viscosity is critical to the behavior of the smallest turbulent scales, it is unclear if these scales are accurately represented by an ILES framework which relies on numerical viscosity.

_{u}Considering the above discussion, the primary goal of the present work is to understand the changes in the smallest turbulent scales within the turbulent premixed flame brush through study of the mean enstrophy, focusing on high *Ka _{u}* flames. This information will be used to assess Kolmogorov’s first similarity hypothesis, i.e., whether the small turbulent scales depend on

*ϵ*and

*ν*alone. At high Karlovitz number, the time and length scales of the smallest turbulence scales are both smaller than those of the flame. As a result, the behavior and important mechanisms are expected to be different than for low

*Ka*flames, in line with the different observed behaviors based on

_{u}*Ka*.

_{u}^{11,12,25}The second goal of this study is to assess model effects on the observed behavior, focusing on chemical and transport models. These results have implications on large eddy simulations (LES), which rely on modeling the behavior of the small turbulent scales.

To accomplish these goals, a series of DNS are performed varying the Karlovitz number and flame density ratio in order to investigate the manner in which enstrophy and the terms in its transport equation vary through premixed flames. The simulations are of a slightly lean *n*-heptane/air flame modeled with a finite-rate chemical model and constant non-unity Lewis numbers. Several additional DNS are performed varying the chemical model and species Lewis numbers to investigate the effects of modeling assumptions on the previous results.

In Sec. II, the governing equations, physical conditions, and numerical setup are described. In Section III, the results are introduced and in Section IV analysis of the enstrophy transport equation is performed. This is followed by Section V which explores the effects of different modeling assumptions on the behavior of enstrophy. Section VI provides an extension of these results to higher Reynolds numbers and comparison with the recent high Karlovitz number slot Bunsen flames of Sankaran *et al.*^{31} Finally, the results are then discussed in light of ILES and LES modeling.

## II. CONFIGURATION OF DIRECT NUMERICAL SIMULATIONS

### A. Physical configuration

The present study considers statistically stationary, statistically planar premixed turbulent *n*-heptane/air flames at a slightly lean equivalence ratio (*ϕ* = 0.9) and atmospheric pressure. The three-dimensional domain has an inflow and outflow at the left and right *x* boundaries, respectively, and periodic boundary conditions in the *y* and *z* directions (Fig. 1). The height and width of the channel are equal and denoted as *L*, while the length, *L _{x}*, is equal to 11

*L*. Based on previous studies,

^{33,34}the turbulence integral length scale is expected to be proportional to the domain width, specifically

*l*= 0.19

_{o}*L*, which is used for the values in Tables I and II and further discussed in Section III A. Tables I and II also provide the length scale $L$, defined as

which yields similar values as *l _{o}*. Here,

*κ*is the wavenumber and

*E*(

*κ*) is the two-dimensional three component velocity spectrum calculated in the unburnt gas using a single

*y*-

*z*plane and averaged over time. A separate DNS is performed of relatively weak, homogeneous, isotropic, triply periodic box turbulence and is used to generate the inflow condition. The mean inflow velocity is constant for each case and set to a value which approximates the turbulent flame speed, allowing for an arbitrary long run-time. This configuration lacks any mean shear so that the effects of the flame on the turbulence may be specifically studied.

. | A . | B . | C . | C′ . | C^{∗}
. | D . |
---|---|---|---|---|---|---|

T (K) _{u} | 298 | 298 | 800 | 500 | 298 | 800 |

ρ/_{u}ρ _{b} | 7.8 | 7.8 | 3.3 | 4.9 | 7.8 | 3.3 |

L (mm) _{x} | 25.6 | 25.6 | 16.8 | 18.7 | 25.6 | 16.8 |

L (mm) | 2.33 | 2.33 | 1.53 | 1.70 | 2.33 | 1.53 |

Grid | 11 × 128^{3} | 11 × 128^{3} | 11 × 128^{3} | 11 × 146^{3} | 11 × 240^{3} | 11 × 220^{3} |

Re _{t} | 83 | 190 | 170 | 290 | 390 | 380 |

Ka _{u} | 70 | 220 | 200 | 650 | 640 | 750 |

u′/S _{L} | 9 | 18 | 19 | 38 | 37 | 45 |

l/_{o}l _{F} | 1.1 | 1.1 | 1.2 | 1.0 | 1.1 | 1.2 |

$L/ l F $ | 1.3 | 1.2 | 1.3 | 1.0 | 1.1 | 1.2 |

η (_{u}μm) | 16 | 9 | 7 | 4.6 | 5.1 | 3.5 |

S (m/s) _{L} | 0.36 | 0.36 | 2.3 | 0.86 | 0.36 | 2.3 |

l (mm) _{F} | 0.39 | 0.39 | 0.25 | 0.32 | 0.39 | 0.25 |

. | A . | B . | C . | C′ . | C^{∗}
. | D . |
---|---|---|---|---|---|---|

T (K) _{u} | 298 | 298 | 800 | 500 | 298 | 800 |

ρ/_{u}ρ _{b} | 7.8 | 7.8 | 3.3 | 4.9 | 7.8 | 3.3 |

L (mm) _{x} | 25.6 | 25.6 | 16.8 | 18.7 | 25.6 | 16.8 |

L (mm) | 2.33 | 2.33 | 1.53 | 1.70 | 2.33 | 1.53 |

Grid | 11 × 128^{3} | 11 × 128^{3} | 11 × 128^{3} | 11 × 146^{3} | 11 × 240^{3} | 11 × 220^{3} |

Re _{t} | 83 | 190 | 170 | 290 | 390 | 380 |

Ka _{u} | 70 | 220 | 200 | 650 | 640 | 750 |

u′/S _{L} | 9 | 18 | 19 | 38 | 37 | 45 |

l/_{o}l _{F} | 1.1 | 1.1 | 1.2 | 1.0 | 1.1 | 1.2 |

$L/ l F $ | 1.3 | 1.2 | 1.3 | 1.0 | 1.1 | 1.2 |

η (_{u}μm) | 16 | 9 | 7 | 4.6 | 5.1 | 3.5 |

S (m/s) _{L} | 0.36 | 0.36 | 2.3 | 0.86 | 0.36 | 2.3 |

l (mm) _{F} | 0.39 | 0.39 | 0.25 | 0.32 | 0.39 | 0.25 |

. | B_{1}
. | B_{Tab,1}
. | B_{OS,1}
. | $ B Tab , 1 4 $ . |
---|---|---|---|---|

T (K) _{u} | 298 | 298 | 298 | 298 |

ρ/_{u}ρ _{b} | 7.8 | 7.8 | 7.3 | 7.8 |

L (mm) _{x} | 25.6 | 25.6 | 25.6 | 60.6 |

L (mm) | 2.33 | 2.33 | 2.33 | 9.32 |

Grid | 11 × 128^{3} | 11 × 128^{3} | 11 × 128^{3} | 2574 × 512^{2} |

Re _{t} | 190 | 190 | 190 | 1150 |

Ka _{u} | 280 | 280 | 250 | 280 |

u′/S _{L} | 21 | 21 | 22 | 33 |

l/_{o}l _{F} | 1.0 | 1.0 | 1.2 | 4 |

$L/ l F $ | 1.0 | 1.0 | 1.2 | 3.7 |

η (_{u}μm) | 9 | 9 | 9 | 9 |

S (m/s) _{L} | 0.29 | 0.29 | 0.27 | 0.29 |

l (m) _{F} | 0.43 | 0.43 | 0.36 | 0.43 |

. | B_{1}
. | B_{Tab,1}
. | B_{OS,1}
. | $ B Tab , 1 4 $ . |
---|---|---|---|---|

T (K) _{u} | 298 | 298 | 298 | 298 |

ρ/_{u}ρ _{b} | 7.8 | 7.8 | 7.3 | 7.8 |

L (mm) _{x} | 25.6 | 25.6 | 25.6 | 60.6 |

L (mm) | 2.33 | 2.33 | 2.33 | 9.32 |

Grid | 11 × 128^{3} | 11 × 128^{3} | 11 × 128^{3} | 2574 × 512^{2} |

Re _{t} | 190 | 190 | 190 | 1150 |

Ka _{u} | 280 | 280 | 250 | 280 |

u′/S _{L} | 21 | 21 | 22 | 33 |

l/_{o}l _{F} | 1.0 | 1.0 | 1.2 | 4 |

$L/ l F $ | 1.0 | 1.0 | 1.2 | 3.7 |

η (_{u}μm) | 9 | 9 | 9 | 9 |

S (m/s) _{L} | 0.29 | 0.29 | 0.27 | 0.29 |

l (m) _{F} | 0.43 | 0.43 | 0.36 | 0.43 |

The turbulence and temperature in the reactants are varied between simulations to investigate the effects of both the unburnt Karlovitz number and flame density ratio. All necessary information about the different simulations is provided in Tables I and II, where *Re _{t}* =

*u*′

*l*/

_{o}*ν*,

*u*′, and

*η*are the turbulent Reynolds number, rms velocity fluctuation, and Kolmogorov length scale, respectively, all calculated in the unburnt gas. The cases in this study are based on the previous work of Savard

_{u}*et al.*

^{32}and Lapointe

*et al.*

^{15}Specifically, cases B and B

_{1}are based on the simulation performed by Savard

*et al.*

^{32}while cases A, C, C′, and D are based on the simulations performed by Lapointe

*et al.*

^{15}These simulations are repeated in this study with a slightly modified turbulence forcing method, discussed in Section II C. Cases C

^{∗}, B

_{Tab,1}, $ B Tab , 1 4 $, and B

_{OS,1}are new to this study. Cases B

_{1}, B

_{Tab,1}, B

_{OS,1}, and $ B Tab , 1 4 $ are performed to test the effects of the transport models, chemical models, and Reynolds number and use the same method described here unless specifically stated otherwise in Sections V or VI A where they are discussed.

Between the cases studied, the unburnt Karlovitz number varies by an order of magnitude (*Ka _{u}* = 70–750) and the unburnt temperature spans practically relevant conditions (

*T*= 298–800 K). Cases A, B, and C

_{u}^{∗}have the same density ratio and are used to test the effects of

*Ka*independently, while the pairs B, C and C

_{u}^{∗}, C′ have the same

*Ka*and are used to test the density ratio independently. While they have the same value of

_{u}*Ka*, cases C′ and C

_{u}^{∗}have different values of

*Re*. Figure 2 shows

_{t}*Ka*and the density ratio for each case as well as their location on the Peters’ regime diagram. These conditions span the transition from the thin to broken/distributed reaction zone regimes.

_{u}To ensure a statistical steady-state, each case is run initially for at least 13 eddy turnover times (*τ _{o}* =

*k*/

*ϵ*, where

*k*is the TKE) to remove any initial transient effects. After this period, data are collected for over 25

*τ*in order to provide sufficient statistical samples. During the simulation, data are collected at a constant rate of approximately 0.5

_{o}*τ*, in order for each data file to represent an independent statistical sample. Further specifications of the simulation conditions are listed in Tables I and II.

_{o}### B. Governing equations

In this study, we solve the low-Mach number reacting flow equations,^{35,36} which include equations for the conservation of mass, conservation of momentum, species transport, and temperature,

In these equations, *ρ* is the density, ** u** is the velocity vector,

*P*is the hydrodynamic pressure, and

**is an applied forcing term. The viscous stress tensor is defined as**

*f* where *μ* is the mixture dynamic viscosity, and *I* is the identity tensor. In the species equations, *Y _{i}* is the mass fraction of species

*i*, $ \omega \u0307 i $ is the species chemical source term, and

*j*_{i}is the species diffusion mass flux vector defined as

where $ u c =\u2212\u2211 D i Y i X i \u2207 X i $ is the correction velocity, *X _{i}* is the species mole fraction, and

*D*is the species molecular diffusivity. In the temperature equation,

_{i}*T*is the temperature,

*α*is the mixture thermal diffusivity, $ \omega \u0307 T $ is the heat source term defined as $ \omega \u0307 T = \u2212 1 c p \u2211 h i \omega \u0307 i $, where

*h*(

_{i}*T*) is the species enthalpy,

*c*

_{p,i}is the species heat capacity, and

*c*is the mixture heat capacity. These equations are combined with the ideal gas law as the equation of state,

_{p}*ρ*=

*P*/

_{o}W*RT*with 1/

*W*= ∑

*Y*/

_{i}*W*. Here, R is the universal gas constant,

_{i}*W*is the species molecular weight, and

_{i}*P*is the thermodynamic pressure.

_{o}The *n*-heptane/air chemistry is modeled with a reduced finite-rate chemical model containing 35 species and 217 reactions (forward and backward counted separately).^{37} An additional reduction of the mechanism from Ref. 37 was performed removing aromatic species, justified by the slightly lean conditions. Constant non-unity Lewis numbers are employed, determined as the Lewis numbers of each species in the burnt region in a one-dimensional unstretched premixed flame simulation using full transport. The species Lewis numbers used in the present simulations are the same as those listed in the work of Savard and Blanquart.^{38} This chemical and transport model has been tested against experimental data and numerical results using full transport (mixture-averaged formulation). Good agreement was found in species mass fractions and species chemical source terms through the flame as well as the laminar flame speed across a range of equivalence ratios.^{38}

The governing equations are solved in this study using the low-Mach, variable density, reacting flow solver NGA.^{39} The chemical source term time integration is performed using a recently developed iterative, semi-implicit method which allows numerical time steps limited here only by the convective Courant-Friedrichs-Lewy (CFL) number, while remaining second order accurate in time and free of lagging errors.^{40} The integration uses an approximation of the diagonal of the chemical Jacobian as the preconditioner, which is calculated at negligible computational cost. The overall scheme is second-order accurate in space and time while discretely conserving kinetic energy.^{39} Scalar transport is performed with the Bounded QUICK scheme, BQUICK.^{41,42} The numerical resolution is designed to resolve all relevant physical length scales of the turbulence and flame, given by considering the criteria *κ _{max}η* > 1.5

^{43}and a minimum of 20 grid points per

*l*.

_{F}^{40}Additional numerical parameters are also provided in Tables I and II.

### C. Turbulence forcing

A variety of configurations have been used in previous DNS studies of premixed turbulent combustion. In one approach, mean shear is present in the flow which generates turbulence through the flame. For example, Sankaran *et al.*^{31,44} studied slot Bunsen flames, which develop downstream of the burner exit, while being statistically stationary. Dunstan *et al.*^{45} investigated V-flames, again, with spatially developing statistics while being statistically stationary. Finaly, Kolla *et al.*^{24} simulated rectangular slot-jet premixed flames which are statistically homogeneous in the plane of the flame, while developing in time. For the purpose of the current study, these configurations present two drawbacks. First, they require directing computational resources towards large domains to contain the mean shear, which limits the highest attainable Karlovitz number for given computational resources. Second, these flows develop in space or time which increases the complexity of isolating the effects of the flame on turbulence. To alleviate these difficulties, several studies have considered another class of configurations, where a statistically planar, statistically steady flame interacts with decaying homogeneous isotropic turbulence. This was used, for example, in the work of Lipatnikov *et al.*^{25} and Chakraborty *et al.*^{14} at low to moderate *Ka _{u}* ($ Ka u \u2217 =\u2009$0.2–0.3 and $ Ka u \u2217 =\u2009$0.5–13, respectively). However, at high Karlovitz numbers, the relative fast decay of turbulence compared to the flame transit time inhibits sustaining high Karlovitz numbers through the flame.

^{32}

Due to the above considerations, most direct numerical simulations of high Karlovitz number premixed flames considered a statistically planar flame and used turbulence forcing to prevent the decay of TKE. This configuration is similar to that of Lipatnikov *et al.*^{25} and Chakraborty *et al.*,^{14} but the decay of turbulence is prevented. Various researchers have used this configuration to investigate both the dynamics of the flame^{15,32,38,46–48} and that of the turbulence.^{12,13,29} For instance, Aspden *et al.*^{4} investigated distributed burning in lean hydrogen flames, Poludnenko and Oran^{3} studied the mechanisms impacting the turbulent flame speed, and Hamlington *et al.*^{13} studied the intermittency of enstrophy.

The present study takes a similar approach while selecting linear forcing, over spectral forcing. It is relevant here to summarize the analysis of Lundgren^{49} on the origins and benefits of linear forcing. By applying a Reynolds decomposition of velocity to the momentum equation, namely, ** u** =

**′ + 〈**

*u***〉, where 〈**

*u***〉 denotes the ensemble average, the transport equation for the fluctuating velocity,**

*u***′, may be derived. This equation contains the term ∇〈**

*u***〉 ⋅**

*u***′, which drives the fluctuating velocity as a result of mean shear. In the TKE transport equation, this term becomes 〈**

*u***′ ⋅ ∇〈**

*u***〉 ⋅**

*u***′〉 and represents the energy transfer from mean shear into turbulent kinetic energy. It is important to highlight that in practical turbulent reacting flows, the mean shear term is non-zero throughout the flame. In other words, energy production occurs in the preheat zone as well as on the burnt side and throughout the flame. Following this analysis, linear forcing in its basic form appends the term,**

*u**A*

**′, to the momentum equation. Momentum is thus injected across all wavenumbers and in proportion to the local velocity, just as in shear flows, with production primarily occurring at the large scales.**

*u*^{33}

Even though some energy (albeit very small) is injected at the small scales, linear forcing has been found to maintain the correct dynamics of these small turbulent scales.^{33} As shown by Carroll and Blanquart,^{50} all forcing techniques (spectral or linear) produce the exact same small scales behavior as characterized by the second order structure function, namely,

Equation (9) was derived analytically by considering a Taylor expansion of the forced Karman-Howarth equation and did not presume any form of the forcing term. This theoretical behavior was also confirmed by numerical simulations, as shown in the same paper. Rosales and Meneveau^{33} also concluded that linear forcing leaves the small scale behavior essentially unaltered.

Appropriate large scale characteristics are also observed. It is important to note that linear forcing does not impose isotropy nor enforces any arbitrary energy spectrum (*κ*^{−5/3} or otherwise) and allows these characteristics to evolve naturally. As a result, the decay rate of the energy spectrum, i.e., the *n* in *E*(*κ*) ∝ *κ*^{−n}, was found to be smaller than 5/3 and to depend on the Reynolds number, a result consistent with a large body of experimental studies.^{50,51} This is shown in Fig. 3 which presents the two-dimensional three component velocity spectra in the unburnt gas for each of the present cases, calculated in a single *y*–*z* plane (at *x* = 1.5*L*) and averaged over time. Superimposed are the energy spectrum decay rates measured experimentally in decaying grid turbulence^{51} for the high Reynolds number limit, *κ*^{−5/3}, and for a comparable *Re*_{λ} of case C^{∗}, *κ*^{−1.35}. Agreement with the latter demonstrates that the current simulations present the correct scaling of the energy spectrum and capture accurately the energy cascade at the large scales.

In the current simulations, linear forcing, modified from the work of Carroll and Blanquart^{34} by subtracting the Favre average velocity, is employed by appending the following term to the momentum equation:

The definition of ** f** differs from that of Savard

*et al.*

^{32}and Lapointe

*et al.*

^{15}by only using the instantaneous local density, not the planar average. This simplifies the TKE and enstrophy budget analysis, but does not alter the results. Here,

*A*is the forcing parameter and is equal to

*A*= 1/(2

*τ*),

_{o}*k*is the target TKE given by the expression $ k o = ( 27 / 2 ) l o 2 A 2 $, and

_{o}*k*(

*x*,

*t*) is the instantaneous Favre

*y*-

*z*planar averaged TKE, defined as $k ( x , t ) = u \u2033 \u22c5 u \u2033 \u02dc /2$. The planar Favre average for an arbitrary field

*ψ*is defined as $ \psi \u0303 = ( \rho \psi ) \xaf / \rho \xaf $ with $ \psi \u2033 =\psi \u2212 \psi \u0303 $ where $ \psi \xaf $ here is the standard planar average,

with $ \psi \u2032 =\psi \u2212 \psi \xaf $. In order to avoid negative velocities at the inflow and outflow, forcing is not applied near these boundaries but within the range 0.5*L* to 8*L* (Fig. 1).

## III. FORCING PERFORMANCE AND ENSTROPHY BUDGET

This section presents an analysis of the effects of the forcing in the simulations. Next, the method of conditional averaging is introduced and global properties of the turbulent flames are presented. This is followed by an overview of the terms in the enstrophy transport equation and qualitative observations of the enstrophy transport budget.

### A. Forcing analysis

The implemented method of turbulence forcing intends to maintain a constant TKE throughout domain (where the forcing is applied). Figure 4(a) displays the planar averaged TKE for case B, which, as expected, is nearly constant in the region of forcing. The TKE shows a slight dip near the location of the flame and then returns to the imposed value behind the flame. This decrease, though small, is in line with the experimental results of Cheng *et al.*^{52} which contain a similar trend.

Next we consider the planar averaged dissipation rate, defined as $ \u03f5 \xaf = \tau \u2032 : S \u2033 \xaf / \rho \xaf $.^{53} Figure 4(b) shows that this quantity is constant in the region of forcing. This may be explained through the TKE transport equation for the case of statistical stationarity and homogeneity in the *y* and *z* directions,

where *u* is the velocity in the *x* direction and *e*_{x} is the unit vector in the *x* direction. The first four terms in this equation, namely, the left hand side (LHS), dissipation, forcing, and dilatation, respectively, are plotted in Fig. 4(c) along with the residual, which represents the cumulative magnitude of all the remaining terms. Within the region where the forcing is active, dissipation and forcing are the dominant terms, even through the flame. Therefore, the TKE transport equation reduces to the balance: $ \u03f5 \xaf =2A k o $, which is indeed the value obtained by $ \u03f5 \xaf $ in Fig. 4(b).

Previous studies in 3D periodic homogeneous, isotropic turbulence found that linear forcing results in an integral length scale, $l= ( 2 k / 3 ) 3 / 2 / \u03f5 \xaf $, which is proportional to the domain size,^{34} namely, about 0.19*L*. The integral length scale for the present configuration is fairly constant through the domain and acquires a value of approximately 0.16*L* (Fig. 4(d)). It is important to note that the effective integral length scale, *l*, differs only slightly from the *a priori* integral length scale, *l _{o}*, used in setting up the simulations (see Section II A). The difference is the result of a slightly smaller TKE in the present work than observed for HIT.

### B. Conditional averaging

This study reports quantities primarily through conditional averages in order to present their behavior *through* the flame. Turbulence quantities are expected to vary through the flame based upon the local thermodynamic properties of the fluid (such as density and viscosity). In a curved and instantaneously transient flame, these quantities correlate better with a flame progress variable, denoted as *C* (Fig. 5(a)), compared to the spatial coordinate *x*. For this reason, averages are conditioned on *C* and denoted as

for a given field *ψ*. As each fluid property (especially density) collapses to a single curve as function of *C*, Reynolds and Favre averages are virtually identical in *C* space. The form of the progress variable chosen is *C* = *Y*_{H2O} + *Y*_{H2} + *Y*_{CO2} + *Y _{CO}*, as it tracks the flame evolution through the preheat and reaction zones. Additionally, its maximum value was found to exhibit little dependence on the unburnt temperature.

^{54}The progress variable range is standardized by considering $ C \u02c6 =C/ C max $ so that 0 represents the reactants and 1 represents the products. Conditional averages are performed excluding the domain outside the region of turbulence forcing due to the presence of very weak turbulence.

By use of this conditional averaging, we define the local Kolmogorov time, length, and velocity scales as

the local Karlovitz number as *Ka*(*C*) = *τ _{F}*/

*τ*

_{η}(

*C*), and the dissipation rate conditioned on

*C*as 〈

*ϵ*|

*C*〉 = 〈

*τ*′ :

*S*″|

*C*〉/〈

*ρ*|

*C*〉. Additionally, we define a quantity involving the flame density change,

*γ*(

*C*) = Δ

*ρ*/〈

*ρ*|

*C*〉, where Δ

*ρ*=

*ρ*−

_{u}*ρ*, which will be used in the subsequent analysis. Further references to the Karlovitz number and the Kolmogorov scales are written as

_{b}*τ*

_{η},

*η*,

*u*

_{η}, and

*Ka*, where the dependence on

*C*is implicit.

It is relevant to note that conditioning on *C* rather than *x* can highlight aspects of the turbulence transformation. This is demonstrated through comparing the dissipation rate conditioned on *x* and *C* for two cases, B and D. As mentioned earlier, *ϵ* conditioned on *x* is constant in the region of forcing for both cases (Fig. 4(b)). However, *ϵ* is able to vary as a function of *C* within each plane, which is observed for case B, and D to a much lesser extent (Fig. 5(b)). As the Karlovitz number increases, the variation in 〈*ϵ*|*C*〉 deceases, as shown by case D.

### C. Global flame properties

Two of the most important quantities to characterize a turbulent flame are the turbulent flame speed, *S _{T}*, and brush thickness,

*l*, which are reported here. The time-dependent turbulent flame speed,

_{T}*S*, and brush thickness,

_{T}*l*, are presented in Fig. 6. The turbulent flame speed is defined as

_{T} where $ \omega \u0307 F $ is the fuel source term and the subscript *u* indicates quantities evaluated in the unburnt gas. The flame brush is defined as the volume of fluid for which the normalized progress variable is between $ C \u02c6 =0.05$ and $ C \u02c6 =0.8$ divided by the cross-sectional area,

This definition applied to the laminar flames under consideration produces nearly the same flame thickness as the thermal width used above, *l _{F}* = (

*T*−

_{b}*T*)/|∇

_{u}*T*|

_{max}. As shown in Fig. 6, the turbulent flame speed and flame brush thickness change in time, but are statistically constant. This supports that the present simulations have achieved a statistically steady-state. All simulations have been performed for at least 25

*τ*, but several were run for a longer period of time. Table III summarizes the time averaged turbulent flame speeds and brush thicknesses. As the Karlovitz number increases, these quantities increase as well. Lapointe

_{o}*et al.*

^{15}further investigated the behavior of these quantities (for a fixed ratio of integral length scale to flame thickness), specifically showing that they vary with the reaction zone Karlovitz number. While employing a different definition of the turbulent flame speed and flame brush thickness, Sankaran

*et al.*

^{31}also observed an increase in these quantities with the Karlovitz number.

### D. Vorticity overview

The enstrophy, *ω*^{2} = ** ω** ⋅

**, transport equation is derived from the momentum equation as**

*ω* where *D*/*Dt* is the material or total derivative. Each term on the right hand side is associated with a specific physical processes: vortex stretching, dilatation, baroclinic torque, viscous dissipation, and forcing, respectively. Vortex stretching, viscous dissipation, and forcing are active in constant density flows, while dilatation and baroclinic torque arise here only due to the presence of the flame. The change of fluid properties (such as density and viscosity) within a premixed flame alters the enstrophy of the incoming turbulence through these five terms.

To begin the discussion of enstrophy, a budget analysis of its transport equation is performed for case B (Fig. 7) to demonstrate behavior common to all cases tested. The other cases show similar trends, but baroclinic torque and dilatation are observed to have even smaller magnitudes relative to the other terms (with the exception of case A for which these two terms are slightly larger). In Fig. 7, the terms are conditionally averaged on *x* specifically to illustrate the flow configuration. In order to minimize numerical errors, the calculation of each term is performed using compact finite difference stencils and minimal spatial interpolations, accounting for the spatial staggering of variables. The flow enters the domain on the left side of the figure with very weak turbulence; vortex stretching, viscous dissipation, and forcing then increase due to the onset of turbulence forcing (*x*/*l _{F}* ≃ 3). These terms maintain fairly constant values in the unburnt region, and subsequently decrease through the flame to smaller values in the burnt region. The flame location may be approximated by the peak density gradient which occurs near

*x*/

*l*= 22 in Fig. 7. Dilatation and baroclinic torque peak within the flame, but they have a smaller magnitude than the other three terms (for all the present cases). Finally, each term approaches zero as the forcing subsides (

_{F}*x*/

*l*≃ 48) prior to the outflow.

_{F}Next, the transformation of vorticity is qualitatively observed by plotting the vorticity magnitude through the flame. Figure 8 displays instantaneous 2D contours of cases A, B, C, C^{∗}, and D. The change in vorticity is dramatic as the magnitude is greatly suppressed through the flame. Preliminary observation suggests that the vorticity is reduced to a lesser extent in cases C and D, which has a similar *Ka _{u}* as B and C

^{∗}, respectively, but a higher unburnt temperature. The vorticity magnitude is significantly altered by the flame for all values of

*Ka*tested, in contrast to the observations of Hamlington

_{u}*et al.*

^{12}who found the vorticity magnitude varied little through the flame at similar

*Ka*values. Possible reasons for these differences will become more apparent in Sec. IV.

_{u}## IV. ANALYSIS AND RESULTS

In this section, scaling estimates for each term of the enstrophy equation, namely, vortex stretching, dilatation, baroclinic torque, viscous dissipation, and forcing, are derived to explain their variation through the flame and across reactant conditions. Predictions from this analysis are tested using results from the present DNS.

### A. Vortex stretching

In dimensional form, the vortex stretching term can vary by more than an order of magnitude across the flame and varies by several orders of magnitude between the different cases (Fig. 9(a)). Scaling of this term requires estimates for vorticity and the rate of strain tensor.

Scaling the rate of strain tensor is accomplished through the viscous dissipation rate. Using the definitions of the dissipation rate and the Kolmogorov scales, we may write

by assuming *μ* = *μ*(*C*) and neglecting spatial gradients of the mean flow. By this, the rate of strain tensor scales with the Kolmogorov time scale, a quantity related to the turbulence characteristics, as well as the velocity divergence, a quantity related to the flame characteristics. The velocity divergence may be rewritten through the continuity equation as −(*Dρ*/*Dt*)/*ρ* and estimated using the density jump across the flame and the laminar flame time scale,

The magnitudes of these two components of the rate of strain tensor are then compared through the ratio,

As this ratio gets larger, the component related to the turbulence increases in magnitude compared to the component related to the flame. Therefore, in the present configuration of high Karlovitz number (*Ka _{u}* > 70), we estimate the magnitude of

*S*with 1/

*τ*

_{η}. This result is consistent with homogeneous, isotropic turbulence where $ S \u2032 : S \u2032 \xaf $ also scales as $ ( 1 / \tau \eta 2 ) $.

^{43}As a spatial gradient of velocity, like the rate of strain tensor, vorticity is estimated as 1/

*τ*

_{η}and enstrophy with $1/ \tau \eta 2 $. Previous experimental and numerical work supports a correlation of vorticity with the Kolmogorov time scale under conditions of homogeneous, isotropic turbulence.

^{26–28}The above analysis results in the scaling,

The same scaling was obtained in the case of homogeneous isotropic turbulence by Tennekes and Lumley.^{55}

Normalized according to the above expression, the vortex stretching terms for each case collapse to a fairly constant value close to 0.15, which is the same value obtained in DNS of homogeneous, isotropic turbulence (Fig. 9(b)). Near $ C \u02c6 =0.95$, the values for cases A and B decrease below 0.15; and at this point, the local values of *Ka*/*γ* are 2.1 and 6.5, respectively, which are the lowest values among all the cases. As the validity of Eq. (22) depends on Eq. (21) being large, this observed decrease may be due to the low values of *Ka*/*γ* and a breakdown of the proposed scaling. Otherwise, the proposed normalization accurately represents the changes in this term through the flame and across all runs. The small magnitude of the constant obtained for this term (0.15) may be due to the preferential alignment of vorticity with the second eigenvector of the rate of strain tensor,^{12} whose eigenvalue is known to be small.^{56,57} The success of the normalization suggests that within the flame the vortex stretching term behaves similar to homogeneous, isotropic turbulence in the limit of high Karlovitz number.

### B. Dilatation

Dilatation in the present configuration is due only to the effects of the flame. Therefore, it will tend toward zero away from the flame. This is shown in Fig. 10(a), where dilatation decreases going towards the largest and smallest values of the progress variable and the peak value varies by nearly six orders of magnitude between cases. Scaling the dilatation term requires estimates of enstrophy and the divergence of velocity, which are analyzed in the previous subsection on vortex stretching. This leads to the following scaling,

Normalized in this manner, results from all cases obtain a similar trend and magnitude (Fig. 10(b)). Considering the wide range of conditions tested, the curves in Fig. 10(b) are sufficiently similar to support that Eq. (23) captures the scaling of dilatation. Additionally, the peak value obtained is near unity, supporting that the normalization captures not only the scaling of this term but the magnitude as well.

### C. Baroclinic torque

Baroclinic torque, like dilatation, is only present due to the density variations within the flame. Likewise, it tends toward zero in the reactants and products. This is shown in Fig. 11(a), which also demonstrates that the peak value varies by orders of magnitude between cases. Though analytically equivalent, baroclinic torque is calculated as ** ω** ⋅ ∇ × (∇

*P*/

*ρ*) rather than

**⋅ (∇**

*ω**ρ*× ∇

*P*)/

*ρ*

^{2}to reduce errors in the numerical evaluation. Scaling of baroclinic torque requires estimates for the gradients of pressure and density. In this estimate, only the magnitude and not alignment of these two vectors is considered.

Analogous to the previous analysis for the rate of strain tensor, we consider pressure decomposed into turbulence and flame related quantities and note that pressure gradients scale as *ρu*^{2}/*l*. Using quantities related to the flame produces the scaling $\rho S L 2 / l F $. Using turbulent quantities, this scaling is largest at the smallest turbulent length scales, providing the scaling $\rho u \eta 2 /\eta $. The ratio of these scalings is given by

where *ν*(*C*) = *α*(*C*)^{2}*S _{L}l_{F}*. Given the present definition of

*l*,

_{F}*α*(

*C*) commonly varies from about 0.1 to 1. Therefore, in the present case of high

*Ka*, the pressure gradients scale with the turbulence quantities, $\rho u \eta 2 /\eta $.

Next, we estimate the gradient of density. Analogous to the temporal derivative of density, we scale spatial density gradients with the density change across the flame, Δ*ρ*, divided by the laminar flame thickness so that $ \u2202 \rho \u2202 x \u2243 \Delta \rho l F $. The following scaling for baroclinic torque is obtained:

When normalized according to the above expression, the variation in the peak value between the six cases reduces from five orders of magnitude to a factor of 2 (Fig. 11). The peak value is small and varies from 0.1 at the lower Karlovitz numbers to a constant value of about 0.04 for sufficiently high Karlovitz numbers. This transition is made more clear in Fig. 12(a), which shows for each case the peak value of the normalized baroclinic torque versus the local Karlovitz number. To possibly reduce this variability with *Ka*, two alternative scalings are proposed: multiplying the proposed scaling with $1/ Ka u $ and replacing *l _{F}* with

*l*(Figs. 12(b) and 12(c), respectively). Both reduce the overall variability; however, the first is somewhat arbitrary and the second requires an

_{T}*a priori*expression for

*l*based on

_{T}*l*and

_{F}*Ka*for further scaling analysis, of which no adequate expression is known. The original scaling is used subsequently for its theoretical basis, suitability for scaling analysis, and consistent behavior at high Karlovitz number, which is the focus of the present study. As the present simulations only include pressure gradients from hydrodynamic pressure fluctuations, the overpressures discussed by Poludnenko

_{u}^{29}and observed to play a role in the production of baroclinic torque are not included here. However, their effect on vorticity transport was observed to be lower at higher Karlovitz numbers,

^{29}suggesting their relevance to the transport of enstrophy diminishes as

*Ka*increases. Their effect at high

*Ka*on vorticity and anisotropy is the subject of future work.

To further investigate the mechanisms controlling the magnitude of baroclinic torque, the normalized gradients of density and pressure are calculated (Fig. 13). Though baroclinic torque decreases with *Ka _{u}*, these quantities do not. The gradient of pressure increases with

*Ka*, due to finite

_{u}*Re*effects, and the density gradient primarily varies with the flame density ratio. Additionally, as both have a magnitude close to unity, the small peak magnitude of the normalized barolinic torque is not explained by the behavior of these quantities. Next, the alignment of the density and pressure gradient vectors is calculated. Figure 13(c) shows for case B that these vectors are either preferentially parallel or anti-parallel. Either cases corresponds to a cross product of zero, so that the preferential alignment modulates the magnitude of baroclinic torque. This alignment of the pressure and density gradients is not accounted for in the proposed scaling (Eq. (25)), hence the small magnitude observed in Fig. 11(b).

### D. Viscous dissipation

Viscous dissipation, like vortex stretching, is present in constant density flow. In dimensional form, it can vary by nearly two orders of magnitude across the flame and it varies by many orders of magnitude between the cases tested (Fig. 14(a)). The primary component of viscous dissipation, ∇ × ((∇ ⋅ *τ*)/*ρ*), scales as *νu*/*l*^{3}. Given a Kolmogorov turbulent cascade, this is again largest for the smallest scales of turbulence resulting in the scaling,

As with several quantities discussed prior, in the case of high *Ka*, viscous dissipation scales with the turbulent quantities, specifically,

In the high Reynolds number limit of homogeneous, isotropic turbulence, Tennekes and Lumley^{55} propose a similar scaling.

Normalized according to the above expression, viscous dissipation obtains a fairly constant value through the flame and across conditions (Fig. 14(b)). As *Ka _{u}* increases between cases A, B, and C

^{∗}, the variation in the normalized quantity decreases. This suggests that the normalization better characterizes the behavior of this term as the Karlovitz number increases, which aligns with the chosen normalization representing the high

*Ka*limit. Additionally, the constant tends toward the value obtained for homogeneous, isotropic turbulence as

*Ka*increases.

_{u}### E. Forcing term

In dimensional form, the forcing term varies by several orders of magnitude between the cases tested (Fig. 15(a)). Assuming the TKE is everywhere equal to the imposed value, *k _{o}*, the forcing term is identically

*Aω*

^{2}. As stated in Section II C,

*A*is equal to

*A*= (2

*τ*)

_{o}^{−1}and is imposed as a parameter of the simulation. Therefore, the forcing term scales as

Normalized in this manner, the forcing term obtains a nearly constant value through the flame and across conditions (Fig. 15(b)). It is observed that the normalized forcing term has less variation as the Karlovitz number increases.

### F. Normalized enstrophy transport equation

The above scaling estimates are used to propose a normalization of the entire enstrophy transport equation, which is then given by

where $ t \u02c6 =t/ \tau F $, $ \omega \u02c6 2 =\u3008 \omega 2 |C\u3009 \tau \eta 2 $, and the Damköhler number is *Da* = *τ _{o}*/

*τ*. Again,

_{F}*Ka*is the local Karlovitz number,

*Ka*(

*C*). As written above, $ T \u02c6 1 $ is vortex stretching, $ T \u02c6 2 $ is dilatation, $ T \u02c6 3 $ is baroclinic torque, $ T \u02c6 4 $ is viscous dissipation, and $ T \u02c6 5 $ is the forcing with each component normalized according to the above discussion. The preceding analysis in this section shows that the normalized terms obtain nearly constant values (or a constant peak value) through the flame and across conditions with $ T \u02c6 1 \u22430.15$, $0< T \u02c6 2 \u22720.7$, $0< T \u02c6 3 \u22720.1$, $\u22120.3< T \u02c6 4 <\u22120.15$, and $ T \u02c6 5 \u22431$, which are shown to be valid at high

*Ka*. Equation (29) is not intended to be precisely solved for enstrophy transport, but provides the scaling and approximate magnitude of each term. It also represents the local balance of vortex stretching, viscous dissipation, and forcing in the present high

_{u}*Ka*cases, as shown in Fig. 7.

_{u}The scaling in Eq. (29) suggests that as the Karlovitz number increases, vortex stretching and dissipation increase in magnitude relative to baroclinic torque, dilatation, and forcing. This is consistent with the observed relative decrease of baroclinic torque and dilatation as *Ka* increased in the results of Hamlington *et al.*^{12} It is also consistent with the behavior observed at low *Ka _{u}* where baroclinic torque and dilation contribute significantly to the transport of vorticity.

^{25}Going further, using Eq. (29), we may predict the relative magnitude of vortex stretching, dilatation, baroclinic torque, and viscous dissipation as a function of the Karlovitz number (Fig. 16). It is important to stress that these predictions are only accurate when the proposed scaling is valid, i.e., at high Karlovitz numbers. Thus, dashed lines in Fig. 16 are used to reflect the uncertainty in extrapolating below the lowest value of

*Ka*tested in this work. Lines in Fig. 16 are calculated using

*α*and

*γ*from the results of case B evaluated at $ C \u02c6 =0.5$, which is near the location of peak dilatation and baroclinic torque. Symbols represent the simulation results with each term again evaluated at $ C \u02c6 =0.5$. The numerical and theoretical values show good agreement, and this agreement improves as the Karlovitz number increases. Using these results, vortex stretching and viscous dissipation are both predicted to be larger than dilatation and baroclinic torque when

*Ka*≳ 20 for T

_{u}= 298. Sufficiently above this value of

*Ka*, it is predicted that vortex stretching and dissipation dominate the transport of enstrophy.

### G. Summary

In the limit of high *Ka*, enstrophy transport results in a local balance between production and dissipation. This is also the case for homogeneous, isotropic turbulence.^{55,58} This confirms theoretically the hypothesis that, in high *Ka* premixed flames, the mean enstrophy behaves locally in a similar manner to constant density, homogeneous, isotropic turbulence. This implies that enstrophy should scale with the Kolmogorov time scale, which is confirmed in Fig. 17, with a proportionality constant close to unity. From homogeneous, isotropic turbulence, it can be shown that

as well as

so that a normalized value of unity is expected.^{43} Variations in the normalized enstrophy decrease as the Karlovitz number increases, which is expected as the normalizations are chosen for the limit of high *Ka*. Therefore, for sufficiently high *Ka _{u}*, the mean enstrophy has the same value as in homogeneous, isotropic turbulence given the local

*ϵ*and

*ν*. This demonstrates that it is the change in kinematic viscosity, as opposed to the effect of the density change through dilatation and baroclinic torque, which drives the enstrophy transformation through the flame.

These results shed light on the validity of Kolmogorov’s first similarity hypothesis within premixed flames, since enstrophy is characteristic of the smallest turbulent scales. Though the fluid properties and turbulence characteristics vary widely across the flame and between the present cases, enstrophy is found to vary only as a function of *τ*_{η}, or equivalently *ν* and *ϵ* alone. This supports that Kolmogorov’s first similarity hypothesis is valid for sufficiently high Karlovitz number premixed flames. In contrast, this result suggests that Kolmogorov’s first similarity hypothesis may not be valid for low *Ka _{u}* premixed flames, but remains to be confirmed. While aspects of these results may not be surprising, high

*Ka*behavior is found to be fundamentally different than low

_{u}*Ka*behavior. For example, Kolla

_{u}*et al.*

^{24}performing DNS of low Karlovitz number flames, show a collapse of the viscous scales in the TKE spectrum when normalizing using the flame thickness. At low

*Ka*, the impact of the flame is more than to simply increase the viscosity, which is the primary effect in high

*Ka*flames, as demonstrated here.

## V. MODEL TESTING

The above conclusions are obtained using finite-rate chemistry and non-unity Lewis numbers. Turbulent combustion research often employs other chemical and transport models. In this section, the effects of chemical and transport models on the conclusions of this work are tested. Case B is chosen for these tests as the Karlovitz number is high enough such that substantial flame broadening and chemical source term variation occurs, but low enough such that differential diffusion effects are not eliminated^{6,32,38} providing a rigorous test case for both transport and chemical models.

### A. Effects of unity-Lewis number assumption

To test the effects of transport models, case B is repeated setting all Lewis numbers to unity, referred to as B_{1}. These two cases are compared primarily to investigate changes in the mean enstrophy, as this study focuses on the transport of this quantity.

Figures 18(a) and 18(b) show that varying the transport model has negligible effects on the dimensional and normalized mean enstrophy. The agreement between cases B and B_{1} can be understood by considering the vortex stretching and viscous dissipation terms. Figure 18(c) shows that vortex stretching is nearly identical between cases B and B_{1}; similar agreement is found for viscous dissipation. As vortex stretching and viscous dissipation dominate the transport of enstrophy at high Karlovitz number, the mean enstrophy evolves similarly in the two cases. Good agreement is also observed in the remaining terms, such as dilatation (Fig. 18(d)). These results suggest that, at high Karlovitz number, the transport model has negligible effects on the smallest turbulent scales. This is in contrast to the importance of the transport model in representing the effects of the turbulence on the flame. For example, Savard and Blanquart^{38} found that with non-unity *Le* transport (compared with unity *Le* transport) the fuel chemical consumption rate in case B exhibited more local extinctions and its mean value was reduced by 40%. This variation in the fuel source term was found to significantly reduce the turbulent flame speed (*S _{T}*/

*S*) compared to the unity

_{L}*Le*flame. While the transport model may affect quantities central to flame propagation, the resulting impact on the mean enstrophy is found here to be negligible. It is important to note that this analysis is in the absence of thermo-diffusive instabilities, and these results may not extend to flames, such as lean hydrogen/air, where these instabilities occur.

### B. Effect of chemical models

To investigate the implications of chemical models, case B_{1} is repeated using two alternative chemical mechanisms: one-step (B_{OS,1}) and tabulated chemistry (B_{Tab,1}). In order to focus on the effects of the chemical model, unity Lewis numbers are used in all three simulations in order to eliminate effects due to differences in the transport models.

Using one-step chemistry, only the reactants and products, *n*-C_{7}H_{16}, O_{2}, H_{2}O, and CO_{2}, are transported. The chemical source terms are determined through a single irreversible, one-step reaction. The same equation is solved as in the work of Lapointe *et al.*,^{59} but with a different fuel. The rate constant, *A*, and activation temperature, *T _{a}*, are chosen so that the 1D unstretched laminar flame speed and thickness closely match that of the finite-rate mechanism (

*A*= 6.0 × 10

^{9}m

^{3}/kg s and

*T*= 15 000 K). As this model excludes intermediate species, the progress variable is defined as the sum of the chemical product mass fractions,

_{a}*C*=

*Y*

_{H2O}+

*Y*

_{CO2}. One-step chemistry requires less computation resources compared to finite-rate chemistry as fewer species must be transported.

Using flamelet generated manifolds (FGM), or tabulated chemistry, the flame is modeled through a single progress variable, *C*, with which fluid properties and the chemical source term are tabulated.^{60–63} Tabulation is performed using the 1D unstretched flame solution of the finite-rate chemical model. Once again, the progress variable is defined as the sum of H_{2}O, H_{2}, CO_{2}, and CO mass fractions and is governed by the transport equation

The use of tabulated chemistry relaxes the computational cost compared to finite-rate chemistry as a single scalar must be transported. Time integration of the chemical source terms is performed explicitly for tabulated and one-step chemistry.

The results obtained with each model are compared primarily through the mean enstrophy as, again, the transport of enstrophy is the focus of this study. Figures 19(a) and 19(b) show that the chemical models induce little differences in the dimensional and normalized mean enstrophy. This is again explained by the agreement of vortex stretching (Fig. 19(c)) and dissipation among the chemical models, which control the transport of enstrophy in high *Ka* premixed flames. As one-step has a slightly different definition of *C* and value of *C _{max}*, a small shift in the curves in $ C \u02c6 $ space and differences near $ C \u02c6 =1$ are expected. A greater effect is observed for the dilatation (Fig. 19(d)) and baroclinic torque terms with one-step chemistry, though this has negligible effects in the transport of enstrophy as their magnitude is small relative to vortex stretching and viscous dissipation. The agreement observed in the mean enstrophy further emphasizes that the smallest turbulent scales evolve through the flame with the local fluid properties, independent of the tested transport or chemical models.

## VI. RESULTS EXTENSION

In this section, the above results are compared with a higher Reynolds number simulation and the recent slot Bunsen flames of Sankaran *et al.*,^{31} demonstrating the applicability of the present conclusions to these alternative configurations.

### A. Higher Reynolds numbers

While this study considers a wide range of Karlovitz numbers, the integral length scale in each DNS discussed thus far is limited to approximately the flame thickness. It is important to consider whether these results are applicable to larger values of *l _{o}*/

*l*, or equivalently, larger Reynolds numbers. This is investigated, first, theoretically through the vorticity spectrum, Ω(

_{F}*κ*) and, second, practically through an additional DNS with a larger integral length scale.

In the present low *Re* DNS, all the enstrophy is contained in scales smaller than the flame thickness. For a larger *Re*, the portion of enstrophy contained in scales smaller than the flame thickness may be determined by considering the model spectrum of Pope for the high *Re* limit,^{43}

and evaluating the vorticity spectrum as Ω(*κ*) = 2*κ*^{2}*E*(*κ*). Using this model spectrum, the wavenumber *κ*^{∗} corresponding to *κ*^{∗} = 2*π*/*l _{F}* for which 80% of the enstrophy is contained in scales smaller than the flame thickness may be found through the integral

Considering an infinite *Re*, at least 80% of the enstrophy is contained in scales smaller than the flame if *l _{F}*/

*η*> 40. This is approximately the condition of case B. Said otherwise, if case B had an arbitrarily large

*Re*but the same

*Ka*, then only 20% of the enstrophy would be in scales larger than the flame thickness. In summary, any high Reynolds number turbulent flame with the same Karlovitz number of the present cases would have nearly the same fraction of vorticity contained in scales smaller than the flame. The present results should thus remain valid for larger integral to flame length scale ratios.

_{u}An additional DNS, labeled case $ B Tab , 1 4 $, is performed to verify if the results are indeed independent of the Reynolds number. This simulation has the exact same conditions as case B_{Tab,1}, but a higher Reynolds number, *Re _{t}* = 1150. The increase in

*Re*is accomplished by increasing both

_{t}*L*as well as

*u*′, in order to maintain the same

*Ka*number, resulting in a four-fold increase in the integral length scale (

_{u}*l*/

_{o}*l*= 4 and

_{F}*u*′/

*S*= 33). Figure 20 presents instantaneous density contour plots from B

_{L}_{Tab,1}and $ B Tab , 1 4 $. While still maintaining sufficient length for the turbulent flame brush,

*L*is set to 6.5

_{x}*L*with the forcing applied between

*L*= 0.5

_{x}*L*and

*L*= 4.5

_{x}*L*. To further minimize the computational cost, the grid is stretched in the

*x*-direction near the x-outflow where the forcing is not applied. Utilizing the results of Section V, which found little effect of the transport and chemical models on the vorticity transformation, this much larger simulation employs tabulated chemistry with unity Lewis number transport. This, again, reduces the computational cost of the simulation (compared to finite-rate chemistry), while still retaining the physics necessary to this study. Due to the larger domain size, each data file represents significantly more independent statistical samples than case B

_{Tab,1}. Therefore, after an initial transient, it is only necessary to collect data over 6

*τ*.

_{o}Cases B_{Tab,1} and $ B Tab , 1 4 $ are compared, as before, by considering the mean enstrophy, vortex stretching, and dilatation terms. Figure 21 shows that increasing the turbulent Reynolds number introduces only small changes to these quantities. Particularly, the mean enstrophy exhibits a similar variation through the flame. This agreement, along with the preceding theoretical analysis, supports that the conclusions of this study are independent of the turbulent Reynolds number (or equivalently the integral to flame length scale ratio). Combined with the results of Section V, this reinforces the primary conclusion of this study. To summarize, first, the unity Lewis number case, B_{1}, has significantly less local extinction than case B.^{15,38} Yet, there is virtually no effect observed on the behavior of the smallest turbulent scales (Fig. 18). Second, B_{Tab,1} and B_{OS,1} have entirely different chemical models and responses to straining and wrinkling. Even more, in the case of tabulated chemistry, all chemical source terms are only functions of the progress variable, and hence remain the same regardless of flow straining and flame wrinkling. Despite these difference, no effects were observed on enstrophy (Fig. 19). Third, case $ B Tab , 1 4 $ has a 4 times larger integral length scale, intentionally introducing straining and wrinkling on a much larger scale (as shown in Fig. 20). Yet again, virtually no effect is observed on the smallest turbulent scales (Fig. 21). These results support the conclusion of the manuscript that, given a sufficiently high Karlovitz number, it is not the large flame dynamics which determines the smallest turbulent scales, but only the local kinematic viscosity and dissipation rate. This is essentially a restatement of Kolmogorov’s first similarity hypothesis.

### B. Slot Bunsen flame

The simulations in the present study consider a statistically planar flame with zero mean shear and forced turbulence. While providing a configuration favorable to fundamental study, other simulations and most practical devices operate under different conditions. To investigate the applicability of the current results and conclusions to alternative configurations, the present results are compared to the recent high Karlovitz number slot Bunsen flames of Sankaran *et al.*^{31} The reader is referred to Ref. 31 for complete details about the DNS. Only a brief overview is provided here.

Sankaran *et al.*^{31} simulated preheated (*T _{u}* = 800 K) methane/air premixed slot Bunsen flames at atmospheric pressure. The slot, containing unburnt reactants, extends in the spanwise,

*z*, direction and is bounded in the transverse,

*y*, direction by a heated, laminar coflow of burnt products. The difference in inlet velocities between the slot and laminar coflow introduces a strong mean shear. The flame and flow spatially develop in the streamwise,

*x*, direction. The present comparison considers case C, which has the highest unburnt Karlovitz number and turbulent Reynolds number among their cases. To remove any confusion with our case C, case C from Sankaran

*et al.*

^{31}is referred to as S-C in the following. The chemistry is modeled using a reduced finite rate chemical mechanism with 13 transported species and 4 additional species assumed in quasi-steady state. Species transport is performed with constant non-unity Lewis numbers.

Downstream of the slot, the Karlovitz number decreases along the centerline. This gives rise to different values of *Ka* through the flame within different *y*–*z* planes, as shown in Fig. 22. As we use a different definition of *Ka _{u}* than used by Sankaran

*et al.*,

^{31}the exact numerical values are expected to be slightly different. These results are compared with the present cases A and B, as they have comparable Karlovitz numbers. In case S-C, the abrupt decrease in

*Ka*near $ C \u02c6 =1$ is likely due to the laminar coflow bounding the turbulent flame. Additionally, the S-C

^{3/4}curve does not extend to $ C \u02c6 =0$ because there are no more pure unburnt reactants. The simulation results may be parameterized by two variables, one to represent the spatial progress in the

*x*direction and another for the progress through the flame. Investigating different planes in

*x*provides a similar function as varying

*Ka*between cases in the present configuration. Therefore, comparison is performed by analyzing the data at three

_{u}*y*–

*z*planes, corresponding to 1/4, 1/2, and 3/4 of the domain length, denoted as S-C

^{1/4}, S-C

^{1/2}, and S-C

^{3/4}, respectively. These are the same locations where Sankaran

*et al.*

^{31}primarily report results. Quantities are then conditionally averaged on the progress variable

*C*in these individual planes to characterize the variation through the flame. The progress variable is again defined as the sum of H

_{2}O, H

_{2}, CO

_{2}, and CO mass fractions. All quantities are computed based on the fluctuating velocity after the mean flow has been subtracted.

Once again, the normalized enstrophy is considered to test the dependence of the smallest turbulence scales on *ϵ* and *ν*. Figure 23 shows good agreement between these two very different configurations. Most notably, the normalized enstrophy is nearly constant with values slightly below unity within the flame (Fig. 23(a)), as observed in the present forced cases. The normalized vortex stretching term is also nearly constant through the flame (Fig. 23(b)), again supporting the proposed scaling. Similar to *Ka* in Fig. 22, low turbulence levels close to the products (near $ C \u02c6 =0.9$) due to the laminar coflow may be responsible for differences in these quantities, particularly for the plane closest to the slot inlet, S-C^{1/4}. Additionally, the curves corresponding to the S-C^{3/4} plane do not extend to $ C \u02c6 =0$ as there are no pockets of unburnt reactants left.

Despite some differences in the curves shown in Fig. 23, the proposed scaling remains valid in the turbulent Bunsen flame. This is made clearer in Fig. 24 by comparing the dimensional enstrophy and vortex stretching to the expected scaling through the flame ($0.05< C \u02c6 <0.9$). In summary, the present simulations and the results from the slot Bunsen flame of Sankaran *et al.*^{31} follow the proposed scalings, and no effects of the physical configuration could be identified.

## VII. DISCUSSION

The implications of the results of this study on LES modeling and ILES are discussed in this section.

### A. LES

While the development of subgrid scale (SGS) models is beyond the scope of this present work, the above results provide/suggest directions for modeling high Karlovitz number premixed flames using LES. This work found that within high *Ka* premixed flames, enstrophy scales in the same manner as in constant density, homogeneous, isotropic turbulence. Once again, enstrophy is characteristic of the smallest turbulent scales and these scales are smaller than the flame thickness at high *Ka*. Therefore, this suggests the scales smaller than the flame thickness may be modeled in the same manner as in homogeneous isotropic turbulence within the framework of large eddy simulations. This corresponds to a LES filter width smaller than the flame thickness: *l _{F}* > Δ >

*η*. While it is unclear if this result holds for filter widths larger than the flame thickness, in Section VI A, we show that the majority of enstrophy is contained in scales smaller than the flame thickness. This suggests the possibility that models based on vorticity

^{64,65}need not be altered for use in high

*Ka*premixed flames. Testing and validation is required to confirm these hypotheses.

### B. ILES

The results of this study have implications on numerical methods, such as ILES,^{30} which rely on numerical viscosity for all or part of the energy dissipation in simulations of premixed turbulent combustion. This discussion builds upon the current work and the recent work of Lapointe *et al.*,^{15} both of which shed light on the importance of temperature dependent viscosity in premixed turbulent combustion.

Recall that one conclusion of the present work is that, at high *Ka*, the smallest turbulent scales throughout the flame are controlled entirely by the local dissipation rate and the local molecular viscosity; this is Kolmogorov’s first similarity hypothesis. Given the importance of viscosity, two simulations performed with different viscosity models would produce correspondingly different transformation of the enstrophy field. More precisely, if the numerical viscosity does not vary through the flame like the physical molecular viscosity, the smallest turbulence scales, and enstrophy, will transform through the flame in an unphysical manner. For instance, in the work of Hamlington *et al.*^{12} which used an ILES framework (with zero molecular viscosity), it was observed that the flame had little effect on the mean vorticity at high turbulence intensities. This behavior is in contrast to the present work, which demonstrated large changes through the flame, and the behavior is likely due to the absence of temperature dependent molecular viscosity which varies through the flame.

A further implication of these results relates to the response of the flame chemical reaction zone to the incoming turbulent flow. Lapointe *et al.*^{15} found that the reaction zone behavior in high *Ka _{u}* premixed flames is controlled by the local reaction zone Karlovitz number,

*Ka*

_{δ}= (

*δ*/

*η*)

^{2}, where

*δ*is the reaction zone thickness, rather than the unburnt Karlovitz number. As

*η*increases significantly through the flame via its dependence on

*ϵ*and

*ν*, any discrepancies in these quantities due to modeling simplification (such as used in ILES) would result in a different value of

*Ka*

_{δ}. In other words, the treatment of molecular viscosity through the flame directly impacts the behavior of the chemical reaction zone. Then, it is doubtful that the distinction between thin reaction zones, broken reaction zones, and distributed reaction zones can be investigated without proper physical molecular viscosity. This discussion highlights certain limitations of ILES given recent results which demonstrate the importance of temperature dependent viscosity in premixed turbulent combustion.

## VIII. CONCLUSION

The primary goal of the present work was to understand the changes in the smallest turbulent scales within a high Karlovitz number turbulent premixed flame brush via the study of the mean enstrophy. For this purpose, a series of high Karlovitz number direct numerical simulations were performed spanning a wide range of Karlovitz numbers and flame density ratios. These simulations of *n*-heptane/air premixed flames employed a reduced finite-rate chemical mechanism and non-unity Lewis number transport. Using the DNS results, a theoretical scaling analysis was performed to estimate the magnitude of each term in the enstrophy transport equation. As a result, a normalized enstrophy transport equation was proposed for high Karlovitz number. Several conclusions are as follows.

The proposed normalized enstrophy transport equation involves a small set of parameters so that the relative magnitude of vortex stretching, dilatation, baroclinic torque, and viscous dissipation may be predicted as a function of the Karlovitz number and flame conditions.

In the limit of high *Ka*, vortex stretching and viscous dissipation dominate the behavior of enstrophy. A balance of these two terms is also observed in homogeneous, isotropic turbulence. As a consequence, given a sufficiently high *Ka*, the mean enstrophy scales in the same manner as homogeneous, isotropic turbulence, namely, with the inverse of the Kolmogorov time scale squared, $1/ \tau \eta 2 $. Therefore, for sufficiently high *Ka _{u}*, the mean enstrophy obtains the same value as in homogeneous, isotropic turbulence given the local

*ϵ*and

*ν*. These conclusions were found to be independent of the Reynolds number and the physical configuration (i.e., forced or unforced).

As *τ*_{η} is only a function of *ϵ* and *ν*, this conclusion supports the validity of Kolmogorov’s first similarity hypothesis in high *Ka* premixed flames. In contrast, this suggests that in moderate to low *Ka* flames, this hypothesis may not be valid and different characteristics of the turbulence through the flame is expected. Further work is required to test the validity of the remaining Kolmogorov hypotheses in high Karlovitz number premixed flames.

Several additional DNS were performed altering the transport and chemical models to test their effects on the transformation of enstrophy. It was found that unity Lewis number transport provided sufficient detail of the flame structure to capture the relevant effects of the flame on turbulence. Additionally, simulations using one-step and tabulated chemistry captured the behavior of the mean enstrophy. This suggests that future research studying the effects of the flame on turbulence does not need to include details of finite-rate chemistry and differential diffusion in numerical simulations to accurately capture the behavior of the smallest turbulent scales. The present study does however underscore the importance of using physical molecular viscosity as opposed to relying on numerical viscosity.

## Acknowledgments

This research was supported by the Department of Defense [Air Force Office of Scientific Research] under Award No. (FA9550-12-1-0144). Views, opinions of, and endorsements by the author(s) do not reflect those of the US Army or the Department of Defense. The authors also gratefully acknowledge the ARCS Los Angeles Chapter and Fonds de Recherche du Quebec—Nature et Techonologies for their funding and support. This work also greatly benefited from a collaboration with Dr. Ramanan Sankaran, Dr. Evatt R. Hawkes, Dr. Chun Sang Yoo, Dr. Jacqueline H. Chen, along with Dr. Hemanth Kolla who shared with the authors the DNS data of Ref. 31.

## REFERENCES

_{7}H

_{16}premixed turbulent flame

_{7}H

_{16}premixed turbulent flames