Differential molecular diffusion (DMD) is a fundamental physical phenomenon that occurs in many fluid flow problems such as turbulent reactive flows. Because DMD is a small-scale event, its modeling is intrinsically challenging, and hence in practical applications, it is more feasible to develop phenomenological models for treating the effect of DMD. In order to develop these phenomenological models, a set of model constraints based on physical observations are needed in order to constrain the model development to yield consistent results with the physical observations. In this work, we adopt an existing power-law Reynolds number scaling of DMD as the model constraints and examine the turbulence modeling requirement of DMD in order to yield the desired scaling. The Reynolds-averaged Navier-Stokes simulations are employed as the modeling framework, and a turbulent mixing layer test case is used as a test case. Perturbation analysis is conducted to examine the model consistency in order to yield the power-law scaling for DMD in the mixing layer test case. It is found that a differential mixing time scale model is needed to yield the power-law scaling, and the commonly used equal mixing time scale model cannot produce the scaling correctly. Numerical simulations of the turbulent mixing problem are also performed to further demonstrate the turbulence modeling requirement for producing the desired power-law scaling of DMD.

## I. INTRODUCTION

In many fluid flow problems involving the transport of multiple components or multiple scalars such as flows in natural water bodies (transport of salinity and heat) and chemically reactive flows (transport of multiple species and energy), molecular diffusion is a fundamental physical process that governs the transport and mixing of different scalars at the molecular scale. The diffusion coefficients of different scalars are usually different, which leads to the phenomenon of differential molecular diffusion (DMD)^{1,2} (or preferential molecular diffusion^{3}). Understanding of DMD and the correct modeling of DMD are thus important for the study of many problems such as the atmospheric/oceanic flows and combustion engines found in transportation vehicles and propulsion devices.

The phenomenon of DMD has significant importance from the standpoint of understanding and accurate modeling of the physics of scalar transport. The transport of the scalar is affected by multi-scale turbulent eddies caused by the cascading of turbulence. The different eddies with different length scales induce spatial and temporal fluctuations of scalars. At the smallest scale, these fluctuations are dissipated by molecular diffusion. A thorough understanding of the different dissipation rates of different scalars due to DMD is crucial in order to completely understand the physics of turbulent scalar transport and mixing. There is plenty of experimental and computational evidence to support the significance of DMD. Drake *et al.*^{4} examined a turbulent hydrogen-air diffusion flame experimentally and numerically. By comparing the measurements and simulation results, they found that DMD significantly affects the transport of major species on the fuel-rich side. In a turbulent piloted jet flame,^{5} the DMD effect is evident in the regions near the jet exit. Barlow, Dunn, and Magnotti^{6} investigated DMD by considering turbulent premixed bluff-body CH_{4}/H_{2} flames. They found that DMD causes a higher local equivalence ratio in the recirculation zone compared to the reactant stream. Besides the experiments, numerous direct numerical simulations (DNSs) have been conducted to demonstrate the importance of DMD. Hilbert and Thévenin^{7} conducted DNS of hydrogen flames with and without DMD and observed a higher temperature prediction by incorporating DMD. Minamoto *et al.*^{8} analyzed the DNS results for a turbulent cross-flow reacting jet flame and found that DMD plays an important role in the preparation of a mixture which favors burning and helps stabilize the flame. The DNS for turbulent sooting flames by Attili *et al.*^{9} demonstrates that DMD leads to a higher total mass of soot precursors which leads to higher total soot mass.

Despite the demonstrated importance of DMD, its treatment in many engineering problems is usually highly simplified because of the difficulty of modeling DMD. The fundamental difficulty arises from the fact that molecular diffusion is a small scale process and it is generally difficult to develop accurate models to incorporate the effect of DMD when the turbulence resolution scale imposed by turbulence models and numerical resolution constraints is much larger than the molecular diffusion length scale. In fact, many models in their original forms typically chose to neglect DMD or completely neglect the effect of molecular diffusion when considering the transport of scalar moments or scalar probability density function (PDF), e.g., the steady flamelet model^{10} and the transported PDF methods.^{11}

More recently, efforts are being made to incorporate the effect of DMD into models. Kronenburg and Bilger^{12,13} derived the equations that incorporate the DMD effect in the conditional moment closure (CMC) model and proposed a model based on the DNS data to close the additional terms associated with the incorporation of DMD in the CMC model. Their results showed that more accurate NO formation rates are predicted in the near field of a turbulent jet flame. In the transported PDF methods,^{11} an approach to treat DMD was presented by McDermott and Pope,^{14} in which the molecular spatial transport of scalars was modeled by a mean shift (MS) model in the particle composition space. Zhang and Wang^{15} improved the MS model by developing a variance consistent mean shift (VCMS) model to yield fully consistent transport of scalar mean and variance through the particles. In the flamelet models,^{10} the laminar flamelet equations with DMD were first derived by Pitsch and Peters.^{16} Unfortunately, incorporating their laminar flamelet equations with DMD in turbulent flame simulations is not trivial and a straightforward incorporation can lead to a significant over-prediction of the DMD effect.^{17} Wang^{2} argued that this over-prediction is caused by neglecting the dependence of the DMD models on the Reynolds number. Based on this argument, Wang^{2} conducted an analysis of the limiting behaviors of DMD and developed a class of DMD flamelet models to correctly incorporate the turbulence effect on DMD and called them the linear differential diffusion (LDD) model and nonlinear differential diffusion (NDD) model. These recent model developments represent the latest advancement to the existing models toward more accurate modeling of fundamental physics such as DMD in turbulence. The work on the modeling of DMD, however, is still highly inadequate. As mentioned previously, the effect of DMD is non-trivial to model because it is a small-scale behavior and physics-based models are generally difficult to develop. Usually, we rely on phenomenological models to incorporate the effect of DMD in existing modeling frameworks. To develop a phenomenological model, we commonly need to develop a set of model constraints or requirements for the model to be consistent with the physical observations. Currently, consistent model constraints for DMD are generally lacking. In this work, we emphasize on a particular model constraint for the DMD modeling, the Reynolds-number-scaling of DMD, and examine the modeling requirements in order to yield the desired scaling law from many existing studies.

There have been some substantial efforts to support a power-law Reynolds-number-scaling of the effect of DMD in multi-component turbulent flows. One of the simplest ways to quantify the effect of DMD is perhaps the *z*_{αβ} parameter,^{1,2,18} e.g., *z*_{HC} = *ξ*_{H} − *ξ*_{C}, the difference of the mixture fractions *ξ*_{H} and *ξ*_{C} in terms of the mass fractions of elements hydrogen H and carbon C, respectively. Bilger and Dibble^{1} first suggested that the mean of *z*_{αβ}, $z\u0303\alpha \beta $ , where the tilde symbol denotes Favre averaging, and the root mean square (RMS) of *z*_{αβ}, *z*_{αβ,RMS} both follow a power-law scaling $Ret\u22121$ in turbulence. The Reynolds number *Re*_{t} used for the DMD scaling study is usually defined as a local turbulent Reynolds number, e.g., *Re*_{t} = $vl$/*ν*, where $v$ is with the integral velocity scale $v=2k/3$, *l* is the integral length scale $l=2k3/3/\epsilon $, *k* is the turbulent kinetic energy, $\epsilon $ is the turbulent kinetic energy dissipation rate, and *ν* is the kinematic viscosity.^{2,19} A second power-law scaling of $z\alpha \beta ,RMS\u223cRet\u22120.25$ was reported by Kerstein, Cremer, and McMurtry,^{20} Nilsen and Kosály,^{21} and Ulitsky, Vaithianathan, and Collins,^{22} in non-reacting flows, based on a more rigorous theoretical analysis. Han *et al.*^{18} and Han and Wang^{19} conducted a rigorous analysis of the scaling of DMD in a series of DNS flames (the Sandia CO/H_{2} temporally evolving jet flames).^{23} They concluded that there is evidence to support the existing scaling of

when a local turbulent Reynolds number is used for the scaling analysis but there is also some evident deviation from these scaling results in these DNS flames.^{19} It is hypothesized that these scaling results can be observed only statistically and the probability to observe them is the highest, as suggested by Han and Wang.^{19} The theoretical scaling results in Eqs. (1) and (2) provide a first-order approximation to the scaling in real problems, which is useful for developing reduced-order models for DMD. The deviation from the theoretical scaling which can be derived from idealized homogeneous isotropic turbulence^{20–22} is arguably caused by many reasons such as insufficiently high Reynolds numbers and deviation from idealized turbulence.^{19} Based on the existing work, there appears to be little evidence to support the scaling of $z\alpha \beta ,RMS\u223cRet\u22121$ by Bilger and Dibble,^{1} and there seems a unique Reynolds-number-scaling for both non-reacting and reacting turbulence in terms of $z\u0303\alpha \beta $ and *z*_{αβ,RMS}. This power-law scaling can be suitably used as the model constraints for the development of phenomenological DMD models. So far, none of the existing DMD models have been thoroughly examined in terms of its consistency with the power-law scaling in Eqs. (1) and (2). More importantly, the requirements for turbulence modeling in order to yield the desired power-law scaling are also not well developed. In this work, we aim to examine the modeling requirements for the turbulent transport modeling of the scalar mean and variance in the Reynolds-averaged Navier-Stokes (RANS) modeling context in order to yield consistency with the DMD power-law scaling. The work is expected to be important to enable turbulent transport models consistent with physical observations as well as to provide a theoretical basis for the future development of fully consistent DMD models.

The rest of the paper is organized as follows. Section II summarizes the models for turbulence and scalar transport used in this work. Section III introduces a turbulent mixing layer problem for the examination of DMD modeling. Section IV examines the consistency of the DMD model for yielding the desired power-law Reynolds-number-scaling by using a perturbation analysis. Numerical simulations of the turbulent mixing layer while incorporating DMD are conducted in Sec. V to further examine the model consistency requirements. The conclusions are drawn in Sec. VI.

## II. CONSISTENT TURBULENT MODELING OF SCALAR TRANSPORT

We aim to model the effect of DMD by using the RANS simulations. A simple turbulence model, the *k* − $\epsilon $ model,^{24,25} is employed in this study for describing the turbulent transport. The governing equations for the flow, turbulence, and scalars are generally written as

where the overline “–” denotes the Reynolds averaging, “″” denotes fluctuations, *ũ*_{i} is the Favre mean velocity, $\rho \xaf$ is the mean density, $P\xaf$ is the mean pressure, *k* is the turbulent kinetic energy, $\epsilon $ is the turbulent kinetic energy dissipation rate, $Sij=1/2(\u2202ui\u0303/\u2202xj+\u2202uj\u0303/\u2202xi)$ is the mean strain rate tensor, *δ*_{ij} is the Kronecker delta, $Y\u0303\alpha $ is the mean mass fraction for the *α*th scalar, $Y\alpha \u2032\u2032Y\beta \u2032\u2032\u0303$ is the covariance between the *α*th scalar and *β*th scalar, *ν* is the mixture molecular viscosity, *D*_{α} is the molecular diffusivity of the *α*th scalar, *ν*_{t} is the turbulent eddy viscosity, *σ*_{k} and *σ*_{$\epsilon $} are the turbulent Prandtl number for *k* and $\epsilon $, respectively, *D*_{t} is the turbulent eddy diffusivity, and $\chi \u0303\alpha \beta $ is the cross scalar dissipation rate defined as

A summation rule is applied to the repeated indices in all the equations discussed in this work except the scalar indices *α* and *β*. When *α* = *β*, the covariance Eq. (8) reduces to the equation for the scalar variance $Y\alpha \u2032\u20322\u0303$. The number of scalars considered in the above equations is arbitrary. No reaction source term is considered in the scalar equations in order to have a non-reacting system for a focused study on the examination and modeling of DMD. It is expected that the modeling requirements for yielding the desired DMD power-law scaling in the non-reacting system are applicable to reacting problems as well since the same power-law Reynolds-number-scaling has been reported in both non-reacting^{20,21} and reacting systems.^{19}

To close the governing equations, *ν*_{t} and *D*_{t} are modeled as *ν*_{t} = *C*_{μ}*k*^{2}/$\epsilon $ and *D*_{t} = *ν*_{t}/*Sc*_{t} by using the *k* − $\epsilon $ model,^{24,25} where *C*_{μ} = 0.09 and *Sc*_{t} is the turbulent Schmidt number. We write cross scalar dissipation rate as

where *τ*_{αβ} is a mixing time scale for the covariance $Y\alpha \u2032\u2032Y\beta \u2032\u2032\u0303$. Conventionally, the mixing time scale is modeled by assuming a constant ratio between the mixing time scale and the turbulence integral time scale,^{24}

where *C*_{ϕ} is a constant usually taken to be 2.0.^{24} In this model, the mixing time scales for the different scalars are equal, and we name it the equal mixing time scale (EMTS) model. The EMTS model for the mixing time scales leads to an important question which is whether it is able to yield the desired power-law Reynolds-number-scaling in Eqs. (1) and (2). For comparative studies, we employ another model for the mixing time scale as follows:

in which *Sc*_{αβ} is an nominal molecular Schmidt number for the scalars *α* and *β*, *Sc*_{αβ} = 2/(1/*Sc*_{α} + 1/*Sc*_{β}),^{22} which implies that a nominal molecular diffusivity for scalars *α* and *β* is *D*_{αβ} = (*D*_{α} + *D*_{β})/2. This mixing time scale model is adapted from Kerstein, Cremer, and McMurtry,^{20} Nilsen and Kosály,^{21} and Ulitsky, Vaithianathan, and Collins.^{22} The model constant *C* in the model is specified to be *C* = 2.0^{22} and *f* is a function of the molecular Schmidt number *Sc*_{αβ},^{26}

This model equation is valid for all values of *Sc*_{αβ}. It is noted that the value of *Sc*_{αβ} does not affect the Reynolds-number dependence ($\u223cRet\u22121/2$) in Eq. (12) and hence it does not affect the Reynolds-number scaling of DMD in Eqs. (1) and (2). The model equation (13) has different behaviors when *Sc*_{αβ} ≤ 1 and *Sc*_{αβ} > 1. This behavior difference is caused by the difference of the scalar dissipation at the small scale.^{26} When *Sc*_{αβ} < 1, the scalar variance dissipates at the Oboukov-Corrsin scale which is larger than the Kolmogorov scale and the scalar variance spectrum in the inertial subrange follows the *κ*^{−5/3} scaling.^{27,28} When *Sc*_{αβ} > 1, the scalar variances dissipate at the Batchelor scale which is smaller than the Kolmogorov scale, and the scalar variance spectrum in the viscous subrange follows the *κ*^{−1} scaling.^{29} The model in Eq. (12) allows difference of the mixing time scale for the different scalars. We name this model as the differential mixing time scale (DMTS) model.

## III. TURBULENT MIXING LAYER PROBLEM

A statistically one-dimensional transient mixing layer (as shown in Fig. 1) with constant density^{2} is employed as the test case to examine the models. The turbulence in the mixing layer is forced^{2} so that the turbulent kinetic energy *k* remains constant, *k* = 1 m^{2}/s^{2}. The kinematic viscosity *ν* is set to be 1.6 × 10^{−5} m^{2}/s so that the turbulence is fully determined by the Reynolds number *Re*_{t}. Other turbulence parameters can be obtained from the work of Wang,^{2} i.e., the turbulent kinetic energy dissipation rate $\epsilon $ = 2*k*^{2}/3*νRe*_{t}, the turbulence length scale $l=3/2k\nu Ret$, the turbulence time scale *τ* = *k*/$\epsilon $ = 3*νRe*_{t}/(2*k*), and the turbulent viscosity $\nu t=C\mu k2/\epsilon =32C\mu \nu Ret$. For scalars in the mixing layer test case, the molecular diffusivity *D*_{α} and the turbulent diffusivity *D*_{t} are determined by the molecular Schmidt number *Sc*_{α} and the turbulent Schmidt number *Sc*_{t}, respectively, i.e., *D*_{α} = *ν*/*Sc*_{α} and *D*_{t} = *ν*_{t}/*Sc*_{t}. Different scalars have different *Sc*_{α} but with the same *Sc*_{t} = 0.7. The transient mixing layer is confined in a domain *x* ∈ [−*L*, *L*].

To simplify the analysis, Eqs. (7) and (8) are first non-dimensionalized. With *y* = *x*/*l* and *s* = *t*/*τ*, we recast Eqs. (7) and (8) into the following for the one-dimensional mixing layer test case which is unsteady and statistically one-dimensional:

where *ξ*_{α} is the mixture fraction defined based on the *α*th scalar as $\xi \alpha =Y\alpha \u2212Y\alpha ,\u2212L/(Y\alpha ,L\u2212Y\alpha ,\u2212L)$, *Y*_{α,±L} are the boundary values of *Y*_{α}, $\xi \alpha \u2032\u2032\xi \beta \u2032\u2032\u0303$ is the mixture fraction covariance, and $\chi ^\alpha \beta =\chi \u0303\alpha \beta \tau $ is the dimensionless cross scalar dissipation rate. The equations for mass, momentum, and turbulence have been eliminated in the homogeneous turbulence mixing layer problem. The DMD effect in the mixing layer test case can be quantified by using the *z* parameter as $z\u0303\alpha \beta =\xi \u0303\alpha \u2212\xi \u0303\beta $ and $z\alpha \beta ,RMS=\xi \alpha \u2032\u20322\u0303+\xi \beta \u2032\u20322\u0303\u22122\xi \alpha \u2032\u2032\xi \beta \u2032\u2032\u03031/2$.

## IV. PERTURBATION ANALYSIS OF TURBULENT MIXING LAYER

Perturbation analysis is performed on the model equations in Sec. III for the mixing layer problem to examine the scaling of the effect of DMD with respect to the Reynolds number. In Sec. IV A, perturbation analysis is first performed to examine the scaling of $z\u0303\alpha \beta $, and then in Secs. IV B and IV C, perturbation analysis is performed to examine the scaling of *z*_{αβ,RMS} with the EMTS model and DMTS model discussed in Sec. II, respectively.

### A. Perturbation analysis for $z\u0303\alpha \beta $

To perform perturbation analysis for $z\u0303\alpha \beta $, we rewrite Eq. (14) as

where *ϑ*_{α} = 2*Sc*_{t}/(3*C*_{μ}*Re*_{t}*Sc*_{α}). For a high *Re*_{t}, *ϑ*_{α} is much less than one, and hence the solution to Eq. (16) can be viewed as a small perturbation to the solution of Eq. (16) with *ϑ*_{α} = 0. We denote the solution to Eq. (16) with *ϑ*_{α} = 0 as $\xi \u0303\alpha y,s=G\alpha ,0(y,s)$. An approximate solution to Eq. (16) can then be written as

where *p* is a key parameter to be determined for the scaling analysis and *G*_{α,m} are the coefficients of the Taylor series that are independent of *Re*_{t}. It is noted that *p* must be positive because $lim\mathit{\vartheta}\alpha \u21920\mathit{\vartheta}\alpha mp=0$ for *m* > 0 so that $\xi \u0303\alpha y,s=G\alpha ,0(y,s)$ when *ϑ*_{α} = 0.

By equating the coefficients of the $\mathit{\vartheta}\alpha mp$ terms on both sides, we can readily obtain the governing equations for $G\alpha ,my,s$ which can be used to determine the scaling parameter *p*. It can be readily found that 1/*p* must be a positive integer to admit solutions to Eq. (18). Otherwise, from Eq. (18), we can get the following two contradicting equations, by equating the coefficients of $\mathit{\vartheta}\alpha mp$ terms with *m* = 0 and *m* = 1 + 1/*p*:

When 1/*p* is a positive integer, i.e., *p* = 1/*N* where *N* is a positive integer, by collecting the first *N* + 1 terms from Eq. (18), we can obtain Eq. (19) and the following equations:

Considering that the initial and boundary conditions for the mean mixture fraction $\xi \u0303\alpha (y,s)$ are in general independent of *Re*_{t} and hence are independent of *ϑ*_{α}, we can obtain the initial and boundary conditions for *G*_{α,m} in Eqs. (21) and (22) as

where ±*L* are the boundary locations of the mixing layer test case in Sec. III. With these initial and boundary conditions, a trivial solution to Eq. (21) can be obtained as

The results in Eq. (24) can also be generalized to the conditions of $iN<m<i+1N$ for any positive integer *i*. Consequently, the solution of $\xi \u0303\alpha $ in Eq. (17) is reduced to

Similarly, the solution of mixture fraction $\xi \u0303\beta $ (*β* ≠ *α*) can be obtained as

The DMD parameter $z\u0303\alpha \beta \u2009$ is then obtained as

Also, it is noted that $G\alpha ,my,s\u2261G\beta ,my,s$ (for any *m*) since they have the same governing equations in Eqs. (21) and (22) and boundary conditions in Eq. (23). Thus $z\u0303\alpha \beta $ can be estimated as, by using the definition for *ϑ*_{α} and *ϑ*_{β} and by keeping the leading order terms in Eq. (27),

Since $G\alpha ,Ny,s$ is independent of *Re*_{t}, $z\u0303\alpha \beta $ is found to be inversely proportional to *Re*_{t}, i.e., $z\u0303\alpha \beta \u223cRet\u22121$. This result shows that the turbulence model outlined in Sec. II is capable of yielding the desired power-law Reynolds-number-scaling for $z\u0303\alpha \beta $ in Eq. (1). It is worthwhile to mention that we obtain the scaling result in Eq. (28) through a perturbation analysis. Alternatively, we can derive the analytical solution to Eq. (16) directly and seek its asymptote at the limit of *Re*_{t} → *∞*, which has already been done by Wang.^{2} The perturbation analysis presented here is a more general approach that can be used for both the scalar mean and covariance, the latter of which does not permit analytical solutions. In Secs. IV B and IV C, we continue to perform the perturbation analysis for *z*_{αβ,RMS} with the EMTS model in Eq. (11) and the DMTS model in Eq. (12).

### B. Perturbation analysis for *z*_{αβ,RMS} with equal mixing time scale model

To perform perturbation analysis for *z*_{αβ,RMS}, we rewrite Eq. (15) as

where *ϑ*_{αβ} = *Sc*_{t}(*Sc*_{α} + *Sc*_{β})/(3*C*_{μ}*Re*_{t}*Sc*_{α}*Sc*_{β}). For a high *Re*_{t}, *ϑ*_{αβ} is much less than one, and hence the solution to Eq. (29) can be viewed as a small perturbation to the solution of Eq. (29) with *ϑ*_{αβ} = 0. We denote the solution to Eq. (29) with *ϑ*_{αβ} = 0 as $\xi \alpha \u2032\u2032\xi \beta \u2032\u2032\u0303y,s=F\alpha \beta ,0(y,s)$. An approximate solution to Eq. (29) can then be written as

where *q* is a parameter to be determined for the scaling analysis of *z*_{αβ,RMS} and *F*_{αβ,m} are the coefficients of the Taylor series. Substituting Eqs. (10), (25), (26), and (30) into Eq. (29), we obtain

in which *τ*_{αβ} can be modeled with either the EMTS model in Eq. (11) or the DMTS model in Eq. (12) in Sec. II. Following the same procedure as shown in Sec. IV A, we can find that 1/*q* must be a positive integer, i.e., *q* = 1/*N* where *N* is a positive integer.

We first consider the EMTS model for *τ*_{αβ} in Eq. (11) for the perturbation analysis of *z*_{αβ,RMS}. Equating the coefficients of the $\mathit{\vartheta}\alpha \beta mq$ terms for the first 2*N* + 1 terms in Eq. (31), we have

Considering that the initial and boundary conditions for the covariance $\xi \alpha \u2032\u2032\xi \beta \u2032\u2032\u0303y,s$ are in general independent of *Re*_{t} and hence are independent of *ϑ*_{αβ}, we can obtain the initial and boundary conditions for *F*_{α,m} in Eqs. (33)–(35) as

With these initial and boundary conditions, a trivial solution to Eq. (33) can be obtained as

This solution can be generalized to $iN<m<i+1N$ for any non-negative integer *i* by considering more terms in the above analysis. Consequently, the solution of $\xi \alpha \u2032\u2032\xi \beta \u2032\u2032\u0303$ in Eq. (30) reduces to

By setting *α* = *β* in Eq. (38), we can obtain the solutions for the variance $\xi \alpha \u2032\u20322\u0303$ and $\xi \beta \u2032\u20322\u0303$ as

It is noted that the equations for *F*_{αβ,0}(*y*, *s*), *F*_{αα,0}(*y*, *s*), and *F*_{ββ,0}(*y*, *s*) are exactly the same, and without losing generality, we can assume that they have the same initial and boundary conditions for the mixing layer test case in Sec. III. As a result, they have identical solutions

We introduce a new quantity *H*_{αβ}(*y*, *s*) as

By combining the equations for *F*_{αβ}, *F*_{αα}, and *F*_{ββ} in Eq. (34), we can obtain the equation for *H*_{αβ}(*y*, *s*) as

From Eq. (36), we can obtain the initial and boundary conditions for *H*_{αβ} as

Consequently, the DMD parameter *z*_{αβ,rms} can be approximated as, by combining the results in Eqs. (38)–(40) and keeping the leading order terms,

in which the aforementioned solution *H*_{αβ} = 0 has been used. This result shows that $z\alpha \beta ,RMS2$ scales as $z\alpha \beta ,RMS2\u2009\u223c\u2009Ret\u22122$ and $z\alpha \beta ,RMS\u2009\u223c\u2009Ret\u22121$. This scaling is obtained by using the EMTS model in Eq. (11). Obviously, this scaling is inconsistent with Eq. (2) obtained from previous theoretical studies,^{19–21} which indicates that the EMTS model in Eq. (11) is incorrectly formulated for modeling DMD.

### C. Perturbation analysis for *z*_{αβ,RMS} with differential mixing time scale model

We next perform the perturbation analysis for *z*_{αβ,RMS} with the DMTS model in Eq. (12). When the DMTS model is used, following a similar procedure in Sec. IV B, we can readily find that if 1/*q* is not equal to *N*, where *N* is an even positive integer, no solutions are permitted. Thus, 1/*q* must be an even positive integer in Eq. (30).

By equating the coefficients of $\mathit{\vartheta}\alpha \beta mq$ for the first *N* terms in Eq. (31), we have Eq. (32) and the following equations:

in which $C\alpha \beta =C\varphi C[3C\mu Sc\alpha Sc\beta /(Sct(Sc\alpha +Sc\beta ))]12$.

By using the same arguments for obtaining the solution in Eq. (37) in Sec. IV B, we can obtain a trivial solution to Eq. (46) as

This solution can be generalized to $iN/2<m<i+1N/2$ for any non-negative integer *i* by considering more terms in the above analysis. Consequently, the solution of $\xi \alpha \u2032\u2032\xi \beta \u2032\u2032\u0303$ in Eq. (30) reduces to

By letting *α* = *β* in Eq. (49), we can obtain the solutions for the variance $\xi \alpha \u2032\u20322\u0303$ and $\xi \beta \u2032\u20322\u0303$ as

The DMD parameter *z*_{αβ,RMS} can then be obtained as, by combining the results in Eqs. (49)–(51),

This shows that $z\alpha \beta ,RMS2$ scales as $z\alpha \beta ,RMS2\u223cRet\u22120.5$ and $z\alpha \beta ,RMS\u223cRet\u22120.25$, which is consistent with the theoretical result^{20–22} in Eq. (2). This consistent scaling is obtained by using the DMTS model in Eq. (12), indicating a need to properly construct turbulence models for yielding the desired power-law Reynolds-number-scaling of DMD. The DMTS in Eq. (12) is such a model that meets this need, while the conventional EMTS model in Eq. (11) does not.

## V. NUMERICAL SIMULATIONS OF TURBULENT MIXING LAYER

Perturbation analysis has been performed in Sec. IV to examine the scaling of $z\u0303\alpha \beta $ and *z*_{αβ,RMS} with respect to *Re*_{t} with both the EMTS and DMTS models. In this section, RANS simulations of the turbulent mixing layer test case are performed to further examine the scaling results.

In the RANS simulations, a second-order central difference scheme is used to discretize Eqs. (7) and (8) in space and a second-order predictor-corrector scheme is used for the temporal discretization. The computational domain for the mixing layer is divided into 256 non-uniform grids (with the dimensionless grid size Δ*y*/*l* = 0.0625) based on a convergence study to ensure that the numerical error is smaller than *O*(10^{−6}). The time step size Δ*t* is specified by ensuring the Fourier number *Fo* = max_{α}(*D*_{T} + *D*_{α})Δ*t*/Δ*y*^{2} < 0.25. The initial and boundary conditions are specified according to the work of Wang.^{2} Five scalars are considered in the numerical simulations, with the Lewis number *Le*_{α} = [0.1, 0.5, 1.0, 2.0, 100]. The Prandtl number *Pr* is set to be *Pr* = 1 and hence *Sc*_{α} = *PrLe*_{α} = [0.1, 0.5, 1.0, 2.0, 100]. More details about the numerical simulations can also be found in the work of Wang.^{2}

Figure 2 shows the profiles of the mean mixture fraction $\xi \u0303\alpha $, variance $\xi \alpha \u2032\u20322\u0303$, covariance $\xi \alpha \u2032\u2032\xi \beta \u2032\u2032\u0303$, $z\u0303\alpha \beta $, and *z*_{αβ,RMS} predicted by the EMTS model and the DMTS model with the different values of Lewis number in the mixing layer with *Re*_{t} = 100. In the results, the scalar index *β* is *β* = 3 with *Le*_{3} = 1 and *α* varies from 1 to 5. From the plot, it is observed that the time scale models (EMTS or DMTS) do not affect the results for the mean mixture fraction $\xi \u0303\alpha $ and the mean DMD $z\u0303\alpha \beta $, expectedly. For the variance and covariance, evident difference can be observed between the different model results with the difference being most obvious for large deviations of Lewis number from *Le* = 1 (e.g., *Le* = 0.1). Consequently, the predictions *z*_{αβ,RMS} by using the different mixing time scale models (EMTS or DMTS) are widely different from each other, as shown in Fig. 2.

To examine the prediction of the Reynolds-number-scaling of DMD (in terms of $z\u0303\alpha \beta $ and *z*_{αβ,RMS}), RANS simulations are performed with a wide range of different *Re*_{t} ∈ [10^{2}, 10^{6}]. The peak values of the predicted $z\u0303\alpha \beta $ and *z*_{αβ,RMS} by using EMTS or DMTS are plotted against *Re*_{t} in Fig. 3. From the figure, we can see that with both the EMTS model and the DMTS model, the scaling of $z\u0303\alpha \beta \u223cRet\u22121$ is correctly captured, while for *z*_{αβ,RMS}, the EMTS model predicts $z\alpha \beta ,RMS\u223cRet\u22121$ and DMTS model predicts $z\alpha \beta ,RMS\u223cRet\u22120.25$. These results are consistent with the perturbation analysis results in Sec. IV. The results further confirm the necessity of developing proper mixing time scale models for yielding the desired power-law Reynolds-number-scaling of DMD and the DMTS model in Eq. (12) is one of these models that is capable of yielding such a scaling consistently.

## VI. CONCLUSIONS

A fundamental physical phenomenon, DMD, is studied in this work and the focus is placed on model constraints for the development of phenomenological models for DMD. An existing power-law Reynolds-number-scaling of DMD is adopted as a model constraint for the model assessment. A turbulent mixing layer test case is considered. Perturbation analysis is carefully conducted to examine the consistency of different models for yielding the desired power-law Reynolds-number-scaling of DMD. It is found that the EMTS model is unable to yield the desired scaling, while the DMTS model is found to be able to yield such a scaling consistently. This finding is important to guide the future model development and simulations to be fully consistent with theoretical results and physical observations. Numerical simulations of the mixing layer test case are also conducted to further demonstrate the model requirements in order to yield the desired power-law Reynolds-number-scaling.

## ACKNOWLEDGMENTS

Acknowledgment is made to the Donors of the American Chemical Society Petroleum Research Fund (Grant No. 53781-DNI9) for support of this research. This paper is based upon work supported by the National Science Foundation under Grant No. CBET-1336075. This research was supported in part through computational resources provided by Information Technology at Purdue University, West Lafayette, Indiana.

## REFERENCES

_{4}/H

_{2}flames

_{2}kinetics