2.1 Brief overview of the velocity field

In this section, the main theoretical aspects of Mons et al. (Reference Mons, Cambon and Sagaut2016) are recalled for the sake of clarity and to introduce important notation. The spectral velocity second-order correlation – or spectral Reynolds stress tensor – is defined as

where $\langle \cdot \rangle$ is an ensemble average, $(\cdot )^{\ast }$ the complex conjugate, $\boldsymbol{k}$ and $\boldsymbol{p}$ are wavevectors of modulus the wavenumbers $k$ and $p$ , $\hat{(\cdot )}$ denotes the Fourier transform and $\hat{u} _{i}$ is the spectral fluctuating velocity. Finally, $\hat{\unicode[STIX]{x1D619}}_{ij}$ corresponds to $\unicode[STIX]{x1D619}_{ij}=\langle u_{i}u_{j}\rangle$ in physical space, and its evolution equation in the presence of a mean velocity gradient $A_{ij}=\text{d}U_{i}/\text{d}x_{j}$ is called the Craya equation and reads

where ${\it\nu}$ is the kinematic viscosity, $M_{ij}(\boldsymbol{k})=({\it\delta}_{in}-2{\it\alpha}_{i}{\it\alpha}_{n})A_{nj}$ , with ${\it\alpha}_{i}=k_{i}/k$ , and $T_{ij}^{NL}$ is the total nonlinear transfer, which can be split into a conservative part and a non-conservative one, the latter being linked to slow pressure terms and thus responsible for the RTI mechanism. The kinetic energy and its dissipation rate are defined as

where $E(k,t)$ is the kinetic energy spectrum. In what follows, the $\mathscr{E}$ – $Z$ decomposition is used (Cambon & Rubinstein Reference Cambon and Rubinstein2006) for a flow without helicity,

where $\text{Re}(\cdot )$ denotes the real part, and $\mathscr{E}-E/(4{\rm\pi}k^{2})$ reflects directional anisotropy, i.e. the difference between the energy in a given direction and the spherical average, and is linked to the kinetic energy spectrum through an integral on a spherical surface $S_{k}$ of radius $k$ ,

Then, the complex-valued scalar $Z$ is linked to the polarization anisotropy, i.e. the difference in energy between components, and is related to the spectral Reynolds tensor through

where $N_{i}$ are helical modes (Cambon & Jacquin Reference Cambon and Jacquin1989; Waleffe Reference Waleffe1992) defined as

where $(\boldsymbol{e}^{(1)},\boldsymbol{e}^{(2)},\boldsymbol{e}^{(3)})$ is the Craya–Herring orthonormal right-handed frame of reference associated with a privileged direction $\boldsymbol{n}$ . The first step of the derivation proposed by Mons et al. (Reference Mons, Cambon and Sagaut2016) consists in using EDQNM to close the generalized transfer terms in the evolution equations of $\mathscr{E}$ and $Z$ , equivalent to (2.2).

The second step starts with the following truncated expansions of the spectral Reynolds stress tensor for $\mathscr{E}$ and $Z$ for moderately anisotropic flows,

with $\mathscr{E}_{0}=E(k,t)/(4{\rm\pi}k^{2})$ , and where $\unicode[STIX]{x1D60F}_{ij}^{(dir)}$ and $\unicode[STIX]{x1D60F}_{ij}^{(pol)}$ are deviatoric tensors that describe anisotropy, defined by a spherical integration

Finally, these truncated expansions allow one to derive evolution equations for the spherically averaged descriptors $E$ , $\unicode[STIX]{x1D60F}_{ij}^{(dir)}$ and $\unicode[STIX]{x1D60F}_{ij}^{(pol)}$ (see table 1). After the spherical integration, the $\boldsymbol{k}$ dependence is transformed into a $k$ one, according to the expansions into spherical harmonics (2.8) and (2.9): these expansions allow part of the angular dependence information to be restored.

2.2 Towards the passive scalar and scalar flux Craya equations

In order to derive the spectral evolution equations of the scalar–scalar and velocity–scalar second-order correlations, the procedure described in the previous section for obtaining the Craya equation is used. The scalar field $T$ is decomposed into mean and fluctuating parts $T={\it\Theta}+{\it\theta}$ , with $\langle {\it\theta}\rangle =0$ , and the mean scalar gradient is written

so that ${\it\Theta}={\it\lambda}_{i}x_{i}$ . In physical space, the evolution equation of the scalar fluctuation ${\it\theta}$ in the homogeneous framework is

where $a$ is the scalar diffusivity. The Fourier transform of the previous equation yields

with the convolution product

Then, the spectral scalar–scalar correlation $\mathscr{E}^{T}$ is defined as

which corresponds to the second-order correlation $K_{T}=\langle {\it\theta}^{2}\rangle$ in physical space. Correlation $\mathscr{E}^{T}$ is real and thus satisfies $\mathscr{E}^{T}(\boldsymbol{k})=\mathscr{E}^{T}(-\boldsymbol{k})$ due to Hermitian symmetry. It is linked to the scalar variance spectrum $E_{T}(k,t)$ through a spherical integral

The Craya equation for the three-dimensional scalar variance spectrum is then

where $T^{T,NL}$ is the total nonlinear scalar transfer, which can be written

with $\int (\ldots )\,\text{d}^{3}\boldsymbol{p}$ the integration over the whole spectral domain, and where $S_{i}^{T}(\boldsymbol{k},\boldsymbol{p},t)$ is the scalar three-point third-order spectral correlation,

The scalar flux $F_{i}(\boldsymbol{k},t)$ , or mixed velocity–scalar correlation, introduced in (2.17) is defined as

which corresponds to the second-order correlation $\unicode[STIX]{x1D619}_{i}^{F}=\langle u_{i}{\it\theta}\rangle$ in physical space. The scalar flux Craya equation is then

where $T_{F}^{T,NL}$ is the nonlinear scalar flux transfer, which can be written

with $S_{nm}^{F}(\boldsymbol{k},\boldsymbol{p},t)$ the scalar flux three-point third-order spectral correlation,

The form of the nonlinear scalar flux transfer (2.22) has been obtained in a recent active scalar study (Burlot et al. Reference Burlot, Gréa, Godeferd, Cambon and Griffond2015). Writing $T_{i}^{F,NL}$ in a manner analogous to the kinetic case (Mons et al. Reference Mons, Cambon and Sagaut2016) is done by introducing

so that the nonlinear scalar transfer reads

where $W_{i}^{F}(\boldsymbol{k},t)=-{\it\alpha}_{i}{\it\alpha}_{m}{\it\tau}_{m}^{F}(\boldsymbol{k},\boldsymbol{p},t)$ represents nonlinear mixed pressure–velocity correlations, from which the RTI of $F_{i}$ is expected, but not prescribed. This RTI mechanism is responsible for the destruction of the scalar flux since it does not exist in isotropic turbulence. The quantity ${\it\tau}_{i}^{F}+{\it\tau}_{i}^{\ast F}$ is a conservative transfer quantity with zero integral over  $\boldsymbol{k}$ .

2.3 Generalized scalar Lin equations and EDQNM closure

The generalized Lin equations for the scalar correlation $\mathscr{E}^{T}$ and the scalar flux $F_{i}$ can be written

Here $T^{T,L}$ and $T_{i}^{F,L}$ are the linear passive scalar and scalar flux transfers, Now, the nonlinear passive scalar and scalar flux transfers $T^{T,NL}$ and $T^{F,NL}$ need to be modelled. The EDQNM approximation was first introduced by Orszag (Reference Orszag1970), Leith (Reference Leith1971) and Orszag (Reference Orszag, Balian and Peube1977) for hydrodynamic turbulence.

Hereafter, an EDQNM closure adapted to homogeneous anisotropic flows is derived similarly to Mons et al. (Reference Mons, Cambon and Sagaut2016). The complete anisotropic form of the second-order spectral tensors is kept in the quasi-normal (QN) formulation, and is therefore expressed in terms of the $\boldsymbol{k}$ -vector descriptors of table 1, whereas the eddy damping (ED), which models the contribution from fourth-order cumulants to the equations that govern triple correlations, is kept in an isotropic relationship. Assuming that the scalar probability distributions are close to normal distributions, the QN approximation is used to model the passive scalar and scalar flux nonlinear transfers. Then, combining isotropic ED, complete QN and Markovianization allows one to write from (2.18) and (2.22)

Here ${\it\theta}_{kpq}^{T}$ and ${\it\theta}_{kpq}^{F}$ are the characteristic times of the third-order scalar and scalar flux correlations, respectively, and ${\it\theta}_{kpq}^{T}$ and ${\it\theta}_{kpq}^{F}$ are obtained by writing the evolution equations of the passive scalar and scalar flux third-order correlations defined in (2.19) and (2.23), respectively. Such an approach for the scalar case has already been used by Bos et al. (Reference Bos, Touil and Bertoglio2005) for HITSG. The eddy-damping terms ${\it\mu}_{2}$ and ${\it\mu}_{3}$ reflect departure from normal laws and, combined with Markovianization, ensure realizability,

with $A_{2}=0$ and $A_{3}=1.3$ . The setting of $A_{2}$ and $A_{3}$ is discussed in Herring et al. (Reference Herring, Schertzer, Lesieur, Newman, Chollet and Larcheveque1982) and Lesieur (Reference Lesieur2008), along with the choice of the eddy-damping terms ${\it\mu}_{2}$ and ${\it\mu}_{3}$ . These two constants are set, based on experimental and theoretical considerations, to recover the Corrsin–Obukhov constant $K_{CO}\simeq 0.8$ . For the kinetic field, the eddy damping used in EDQNM is

as defined by Pouquet et al. (Reference Pouquet, Lesieur, Andre and Basdevant1975) – as an improved version, more general (assessed in two and three dimensions, at various Reynolds numbers) than in Orszag (Reference Orszag1970) – where $A_{1}=0.355$ , so that the Kolmogorov constant is recovered, $K_{0}\simeq 1.3$ . As a result of the QN approximation, $\mathfrak{T}_{i}^{T,QN}$ and $\mathfrak{T}_{ij}^{F,QN}$ can be written

With this closure, and with calculations very similar to the kinetic case, the nonlinear scalar transfer becomes

The single and double primes $^{\prime }$ and $^{\prime \prime }$ refer to wavevectors $\boldsymbol{p}$ and $\boldsymbol{q}$ , respectively; $x$ , $y$ and $z$ are cosines of the angles formed by the triad vectors $(\boldsymbol{k},\boldsymbol{p},\boldsymbol{q})$ ; and $X=Z\text{e}^{2\text{i}{\it\lambda}}$ , $X^{\prime }=Z^{\prime }\text{e}^{2\text{i}{\it\lambda}^{\prime }}$ and $X^{\prime \prime }=Z^{\prime \prime }\text{e}^{2\text{i}{\it\lambda}^{\prime \prime }}$ , where ${\it\lambda}$ , ${\it\lambda}^{\prime }$ and ${\it\lambda}^{\prime \prime }$ are rotation angles of the triad around $\boldsymbol{k}$ , $\boldsymbol{p}$ and $\boldsymbol{q}$ , respectively, and $Z$ has been defined in (2.6). Since the scalar flux exists only in an anisotropic configuration, and is thus a purely anisotropic field, the second part of (2.38) is neglected for moderately anisotropic flows, where all quadratic anisotropic contributions can be neglected. The closure for $F_{i}(\boldsymbol{k})$ is developed with an appropriate decomposition in the following section.

3.1 Spherically averaged descriptors

A decomposition similar to the one of $\mathscr{E}$ in the kinetic case is used for the scalar correlation

with $\mathscr{E}^{(T,dir)}=-15\mathscr{E}_{0}^{T}\unicode[STIX]{x1D60F}_{ij}^{(T)}{\it\alpha}_{i}{\it\alpha}_{j}$ and $\mathscr{E}_{0}^{T}=E_{T}/(4{\rm\pi}k^{2})$ . The following decomposition is chosen for the scalar flux:

This decomposition is consistent with the scalar flux being a solenoidal field ( $k_{i}F_{i}=0$ ). These equations (3.1) and (3.2) are exact in HITSG (Herr et al. Reference Herr, Wang and Collins1996), whereas they are truncations at the first non-trivial order (the second order) of the scalar correlation $\mathscr{E}^{T}$ and scalar flux $F_{i}$ expansions for the cases involving mean velocity gradients (HST and HSTSG). These truncations are consistent with the modelling of $\mathscr{E}$ and $Z$ done in (2.8) and (2.9), respectively.

Hermitian symmetry for the scalar flux, $F_{i}(-\boldsymbol{k})=F_{i}^{\ast }(\boldsymbol{k})$ , is straightforward from the decomposition of a vector into helical modes. Notably, it implies that the vector $\mathscr{E}_{j}^{F}$ is purely real. About the distinction between directional and polarization anisotropies: scalars admit directional anisotropy only. Consequently, $\unicode[STIX]{x1D60F}_{ij}^{(T)}$ , which appears in the expansion of $\mathscr{E}^{T}$ , reflects directional anisotropy, as does $\unicode[STIX]{x1D60F}_{ij}^{(dir)}$ in the decomposition of $\mathscr{E}$ . Then, as $F_{i}$ is a vector, $\mathscr{E}_{j}^{F}$ represents polarization anisotropy. Finally, a solenoidal second-order tensor such as $\hat{\unicode[STIX]{x1D619}}_{ij}$ admits both contributions (Cambon & Rubinstein Reference Cambon and Rubinstein2006). This classification is summarized in table 1.

An additional solenoidal contribution in terms of the Levi-Civita permutation tensor is discussed by O’Gorman & Pullin (Reference O’Gorman and Pullin2005), resulting in the so-called quadrature spectrum. It can be shown that, using a toroidal–poloidal decomposition of a solenoidal vector, the latter contribution is redundant with the projection part in (3.2). Therefore, as the helicity spectrum (Gomez, Politano & Pouquet Reference Gomez, Politano and Pouquet2000) is neglected in the spectral Reynolds stress tensor, this additional solenoidal contribution is consistently discarded.

The scalar anisotropic descriptor $\unicode[STIX]{x1D60F}_{ij}^{(T)}$ and the spherically averaged scalar flux $E_{i}^{F}$ are defined similarly to the kinetic ones $\unicode[STIX]{x1D60F}_{ij}^{(dir)}$ and $\unicode[STIX]{x1D60F}_{ij}^{(pol)}$ using spherical integration,

The scalar anisotropy tensor $\unicode[STIX]{x1D623}_{ij}^{T}$ is then given by

and is analogous in the description of anisotropy to $\unicode[STIX]{x1D623}_{ij}^{(dir)}$ for the velocity field (Mons et al. Reference Mons, Cambon and Sagaut2016). Similar scalar anisotropy descriptors were introduced by Kassinos et al. (Reference Kassinos, Knaepen and Carati2007), with a slightly different convention for structure-based modelling. Finally, the scalar energy $K_{T}$ , or scalar variance, and its dissipation rate are defined as

3.2 Spherical average of the passive scalar and scalar flux fields

The final evolution equations for the spherically averaged scalar descriptors, derived from the previous generalized Lin equations, and briefly called spherically averaged scalar Lin equations, read

Here $S^{T,NL(iso)}$ is the classical nonlinear spherically averaged isotropic scalar transfer,

which is exactly the one derived by Lesieur (Reference Lesieur2008); and $S_{ij}^{T,NL(dir)}$ is the nonlinear spherically averaged directional scalar transfer,

with ${\it\Delta}_{k}$ the domain where $k$ , $p$ and $q$ are the moduli of the triad vectors. The transfer $S^{T,NL(iso)}$ is a conservative transfer, meaning that its integral over $k$ is zero. However, the integral of $S_{ij}^{T,NL(dir)}$ is different from zero (like the directional transfer of the velocity field). From the non-conservative part of $S_{ij}^{T,NL(dir)}$ , an RTI of the passive scalar is expected. Nevertheless, as there is no fluctuating pressure in the passive scalar equation, it cannot return to isotropy on its own, and therefore this mechanism is driven by the dynamics of the velocity field only. For the linear terms, $S^{T,L(iso)}$ is the linear spherically averaged isotropic scalar transfer,

and $S_{ij}^{T,L(dir)}$ is the linear spherically averaged directional scalar transfer,

where $A_{ij}^{+}$ and $A_{ij}^{-}$ respectively denote the symmetric and antisymmetric parts of the velocity mean-gradient matrix, ${\it\lambda}_{i}$ is the mean scalar gradient and $E_{i}^{F}=4{\rm\pi}k^{2}\mathscr{E}_{i}^{F}$ .

Regarding the scalar flux, its nonlinear transfer becomes, at first order in anisotropy,

One can observe that in the framework of moderate anisotropy, only the isotropic part of the velocity field $\mathscr{E}_{0}=E/(4{\rm\pi}k^{2})$ remains. The nonlinear spherically averaged scalar flux transfer is then

One can extract from this nonlinear scalar-flux transfer the mixed pressure–velocity part, which is expected to have an RTI function

And $S_{i}^{F,L}$ is the spherically averaged linear scalar flux transfer,

As a conclusion, the dynamics of a passive scalar field in homogeneous anisotropic turbulence is driven by six spherically averaged compact equations within this anisotropic EDQNM modelling: three for the velocity field first given in Mons et al. (Reference Mons, Cambon and Sagaut2016),

two for the passive scalar (equations (3.7) and (3.8)), and one for the scalar flux (equation (3.9)). The last three are original results of the present work. The different anisotropy descriptors are gathered in table 1.

3.3 Numerical parameters, initial conditions and validity of the model

The time evolution of the kinetic spectra $(E,E\unicode[STIX]{x1D60F}_{ij}^{(dir)},E\unicode[STIX]{x1D60F}_{ij}^{(pol)})$ and the passive scalar and scalar flux spectra $(E_{T},E_{T}\unicode[STIX]{x1D60F}_{ij}^{(T)},E_{i}^{F})$ are obtained by solving $20$ coupled integro-differential equations using a third-order Runge–Kutta scheme with implicit treatment of diffusion terms. The wavenumber space is discretized using a logarithmic mesh $k_{i+1}=10^{1/f}k_{i}$ for $1\leqslant i\leqslant n$ , where $n$ is the number of modes and $f=17$ is the number of discrete points per decade. This mesh spans from $k_{min}$ to $k_{max}=10k_{{\it\eta}}$ , where $k_{{\it\eta}}=({\it\epsilon}/{\it\nu}^{3})^{1/4}$ is the Kolmogorov wavenumber. Such a resolution at high Reynolds numbers allows a fine description of the viscous range, with $k_{min}$ higher than the ones usually used in DNS. The time step is monitored by defining a constant Courant–Friedrichs–Lewy (CFL) number and is obtained by considering the characteristic time scales of scalar and kinetic dynamics. The values of the physical and numerical parameters are given for simulations at large and small Reynolds numbers in table 2. If not mentioned otherwise, the initial conditions are isotropic ( $F_{i}(k,t=0)=0$ ) with the kinetic energy spectrum $E(k,t=0)$ proposed by Pope (Reference Pope2000) and Meyers & Meneveau (Reference Meyers and Meneveau2008),

where $f_{L}$ and $f_{{\it\eta}}$ are shape functions for large and small scales, respectively,

and ${\it\sigma}$ is the infrared slope. Simulations have shown that initial moderately anisotropic conditions only slightly impact the early dynamics, without altering the asymptotic states. Thus, classical isotropic initial conditions are used such as $E_{T}(k,0)=E(k,0)$ , avoiding the introduction of supplementary numerical parameters.

In what follows, the Prandtl number is equal to one, meaning that $a={\it\nu}$ , and integrated quantities correspond to quantities obtained by integration over the whole wavenumber space, such as the kinetic energy.

It has been mentioned several times that the present model is valid for moderate anisotropy as quadratic contributions in terms of anisotropy descriptors have been neglected. In this view, criteria need to be given to quantify the upper bounds of the anisotropy that the model can handle at the levels of the velocity and scalar fields. Such criteria, or realizability conditions, can be derived directly from (2.8) and (3.1), to ensure that the kinetic and scalar energy densities $\mathscr{E}$ and $\mathscr{E}^{T}$ are positive. This yields

where $\mathscr{L}_{i}$ and $\mathscr{L}_{i}^{T}$ are eigenvalues of $\unicode[STIX]{x1D60F}_{ij}^{(dir)}$ and $\unicode[STIX]{x1D60F}_{ij}^{(T)}$ , respectively. In the following simulations, these realizability conditions are a posteriori systematically verified as satisfied. Moreover, it is worth noting that the violation of (3.23) immediately yields negative values for the spectra $E(k,t)$ and $E_{T}(k,t)$ . Hence, moderate anisotropy could be roughly defined as the case for which the intensity of the mean-field gradients is such that the realizability conditions are satisfied.

4.1 Homogeneous isotropic turbulence with a mean scalar gradient (HITSG)

In this section, the emphasis is put on the specific case where the scalar flux appears due to a uniform mean gradient ${\it\lambda}_{j}=(0,0,-{\it\Lambda})$ with ${\it\Lambda}>0$ , whereas the velocity field is purely isotropic and decreases with time. This configuration has been widely investigated (Bos et al. Reference Bos, Touil and Bertoglio2005; O’Gorman & Pullin Reference O’Gorman and Pullin2005; Bos Reference Bos2014): turbulent eddies bring the hot fluid to the cold parts of the flow (and the opposite), thus creating a turbulent flux.

4.2 Homogeneous shear turbulence (HST)

In this section, qualitative comparisons and general results are given for the passive scalar field submitted to a shear without any scalar gradient. This particular framework has not received much attention in experiments or DNS, unlike the case with both velocity and scalar gradients that will be addressed in the next section. Firstly, the case of a mean shear imposed at short times only and then released is investigated, and finally the emphasis is put on the passive scalar field submitted to a sustained shear. For a shear-driven flow, the mean velocity gradient matrix is written $A_{13}=\text{d}U_{1}/\text{d}x_{3}=-S$ with $S>0$ . The dimensionless time $St$ refers to the quantity of anisotropy injected in the flow. In both cases, the scalar spectrum displays a $k^{-5/3}$ inertial–convective range as revealed in figure 5.

Figure 5. (a) Scalar spectrum $E_{T}(k,t)$ in a sustained shear flow, with $St=50$ ( $Re_{{\it\lambda}}=10^{4}$ ) and ${\it\sigma}=2$ , along with the scalar integral and Kolmogorov wavenumbers $k_{T}$ and $k_{{\it\eta}}$ , as a function of the wavenumber $k$ . (b) Scalar decay exponents ${\it\alpha}_{T}$ in high-Reynolds-number regime, with $St=1$ for ${\it\sigma}=2$ and ${\it\sigma}=4$ , as a function of the Reynolds number $Re_{{\it\lambda}}$ . Symbols represent the theoretical values expected from (4.17): ○ is for ${\it\sigma}={\it\sigma}_{T}=2$ , and ▫ is for ${\it\sigma}={\it\sigma}_{T}=4$ .

4.3 Homogeneous shear turbulence with a mean scalar gradient (HSTSG)

This final section focuses on homogeneous shear turbulence with a mean scalar gradient (HSTSG). The emphasis is put on the impact of both mean velocity and scalar gradients on the scalar flux and the passive scalar.

4.4 Conclusion on the comparisons

The present anisotropic EDQNM modelling for the passive scalar and its flux has been compared with several experiments and DNS. The agreement between asymptotic values of integrated quantities, and even sometimes in the transitory state, validate the modelling introduced in §§ 2 and 3. Nevertheless, it has to be mentioned that the present model seems to slightly under-estimate the turbulent Prandtl number.

5.1 Decay and growth laws in HITSG

The results presented here are an extension of an analysis previously applied to the passive scalar in HIT (Briard et al. Reference Briard, Gomez, Sagaut and Memari2015). The theoretical decay exponents of kinetic and scalar integrated quantities are recalled in table 3.

The cospectrum $\mathscr{F}$ is continuously destroyed by the decay of the velocity field, and is consequently also decaying. Theoretical decay exponents can be derived similarly to what Comte-Bellot & Corrsin did for the passive scalar. Two assumptions, based on physical arguments for high Reynolds numbers, are needed. Firstly, figure 2 revealed that the cospectrum infrared dynamics is driven only by the kinetic one. Then, we assume that the main contribution to the decay mechanism comes from the inertial range and consequently that only the $k^{-7/3}$ inertial range should be taken into account to determine the decay exponents of $\mathscr{F}$ . Secondly, in the case of Batchelor turbulence, when the permanence of the large-eddy hypothesis is broken, backscatter parameters $p$ and $p_{T}$ are introduced for the velocity and scalar fields, respectively, in HIT (Briard et al. Reference Briard, Gomez, Sagaut and Memari2015). Since $\mathscr{F}$ is the spectral counterpart of the velocity–scalar cross-correlation, its backscatter parameter $p_{\mathscr{F}}$ should contain both effects. Hence, one defines $p_{\mathscr{F}}=(p+p_{T})/2$ for ${\it\sigma}=4$ (so that $p_{\mathscr{F}}=0$ for ${\it\sigma}\leqslant 3$ ). This gives $p_{\mathscr{F}}=0.4075$ in the case $Pr=1$ for ${\it\sigma}=4$ . These reasonable assumptions yield

Injecting in this equation the decay exponents of ${\it\epsilon}$ and $k_{L}$ recalled in table 3, and using $p_{\mathscr{F}}$ , one obtains the theoretical decay exponent of the cospectrum for high Reynolds numbers,

In the low-Reynolds-number regime, when the cospectrum survives, the production term $P_{\mathscr{F}}$ drives the dynamics of the velocity–scalar correlation. The RTI term ${\it\Pi}_{\mathscr{F}}$ and the cospectrum dissipation rate ${\it\epsilon}_{\mathscr{F}}$ can be considered as negligible. Indeed, considering the form of the cospectrum dissipation (4.4), this quantity reaches its highest values in the inertial range, which tends to disappear for low Reynolds numbers. Therefore, the evolution equation (4.5) of the cospectrum energy becomes $(\text{d}K_{\mathscr{F}}/\text{d}t)\sim P_{\mathscr{F}}=2{\it\Lambda}K/3$ . This immediately yields

which does not depend on the scalar gradient ${\it\Lambda}$ .

These theoretical decay exponents for the cospectrum for high and low Reynolds numbers are assessed numerically in figure 11. The agreement is excellent, even for the particular case of Batchelor turbulence. Such laws could not have been observed before since it requires high Reynolds numbers in decaying turbulence, whereas most of the DNS of HITSG consider a stationary velocity field at moderate Reynolds numbers. These theoretical results, assessed numerically, give further insights into the prediction of high-Reynolds-number decay in HITSG. From (5.2) and (5.3), it follows that $K_{\mathscr{F}}$ does not decay for ${\it\sigma}=1$ . In this case, the Reynolds number $Re_{{\it\lambda}}$ remains constant and so the dynamics of the inertial range, on which ${\it\alpha}_{\mathscr{F}}$ is based, remains unchanged.

The dissipation rate ${\it\epsilon}_{\mathscr{F}}$ is now investigated. Probably because the scalar flux is a purely anisotropic quantity, which implies that ${\it\epsilon}_{\mathscr{F}}$ is not conserved in the inertial range, unlike ${\it\epsilon}$ and ${\it\epsilon}_{T}$ , one cannot express its time dependence as a power law. Nevertheless, as the inertial range disappears for low Reynolds numbers, obtaining the decay exponent $n_{{\it\epsilon}_{\mathscr{F}}}$ of ${\it\epsilon}_{\mathscr{F}}$ is possible in this regime, considering the evolution equation (4.5), so that

The agreement between this theoretical result and numerical simulations is displayed in figure 11 as well.

Figure 11. Decay exponents of the cospectrum energy ${\it\alpha}_{\mathscr{F}}$ (a) and of the cospectrum dissipation $n_{{\it\epsilon}_{\mathscr{F}}}$ (b) for different infrared slopes ${\it\sigma}$ , as a function of the Reynolds number $Re_{{\it\lambda}}$ . Symbols represent theoretical predictions: ▫ for high Reynolds numbers, given by (5.2); ○ for low Reynolds numbers, given by (5.3) and (5.4).

The emphasis is now put on the effect of the mean scalar gradient ${\it\Lambda}$ on the passive scalar itself. The scalar spectrum still displays a $k^{-5/3}$ inertial–convective range despite the mean gradient, as revealed in figure 1, and this was obtained experimentally as well (Mydlarski & Warhaft Reference Mydlarski and Warhaft1998). From the evolution equation of the scalar variance (4.19), it follows that the dynamics of $K_{T}$ is driven by the production term in the presence of a mean scalar gradient ${\it\Lambda}$ . Using the previous results regarding the decay exponents of $K_{\mathscr{F}}$ , the theoretical exponent of $K_{T}$ in the presence of a mean scalar gradient ${\it\alpha}_{T}^{{\it\Lambda}}$ can be found. In the high-Reynolds-number regime, one has

and for low Reynolds numbers,

Note that the exponents do not depend explicitly on the scalar gradient ${\it\Lambda}$ . The agreement between these theoretical expressions of ${\it\alpha}_{T}^{{\it\Lambda}}$ and numerical simulations is presented in figure 12 in both high- and low-Reynolds-number regimes. For high Reynolds numbers, the scalar variance grows in time whatever ${\it\sigma}$ is, whereas for low Reynolds numbers, it decays only for ${\it\sigma}=4$ . This can be explained using some physical arguments. The theoretical prediction (5.6) of ${\it\alpha}_{T}^{{\it\Lambda}}$ is based on the fact that the dynamics of $K_{T}$ is driven by the production term $2{\it\Lambda}K_{\mathscr{F}}$ , and $K_{\mathscr{F}}$ is greater for smaller ${\it\sigma}$ . Consequently, for high infrared slope such as ${\it\sigma}=4$ , the spectrum $E(k,t)$ , and thus $\mathscr{F}(k,t)$ , has less energy in large scales, resulting in a weak production term for the scalar variance that does not balance dissipation. This is consistent with HITSG experimental results at moderate $Re_{{\it\lambda}}$ , where the scalar variance can grow or decay depending on the intensity of ${\it\Lambda}$ : with a weak ${\it\Lambda}$ , $K_{T}$ still decays (nevertheless more slowly than in HIT). Therefore, a strong mean scalar gradient ${\it\Lambda}$ corresponds to a small infrared slope ${\it\sigma}$ , or equivalently corresponds to energetic large scales.

In the particular case of Saffman turbulence ( ${\it\sigma}=2$ ), the value ${\it\alpha}_{T}^{{\it\Lambda}}=4/5$ was already found by Chasnov (Reference Chasnov1995). There, the decay and growth laws of passive and active scalar fields, with and without buoyancy or mean gradient, are studied. Power laws for the active scalar field were assessed by large-eddy simulations (LES). Hence, the present EDQNM simulations validate the power laws for the passive scalar field, with an explicit dependence on the initial large-scale conditions ${\it\sigma}$ . Consequently, this analysis can be considered as an extension of Chasnov (Reference Chasnov1995) regarding the growth of a passive scalar in HITSG.

Figure 12. Growth and decay exponents of the scalar energy ${\it\alpha}_{T}^{{\it\Lambda}}$ (a) and of the scalar dissipation rate $n_{{\it\epsilon}_{T}}^{{\it\Lambda}}$ (b) for different infrared slopes ${\it\sigma}$ , as a function of the Reynolds number $Re_{{\it\lambda}}$ . Symbols represent theoretical predictions: ▫ for high Reynolds numbers, given by (5.5); ○ for low Reynolds numbers, given by (5.6).

Then, from the scalar variance evolution equation (4.19), the scalar dissipation should evolve in time similarly as the cospectrum energy,

in both high- and low-Reynolds-number regimes. This is assessed in figure 12 as well. All the new results regarding the growth or decay predictions for $K_{\mathscr{F}}$ , $K_{T}$ , ${\it\epsilon}_{\mathscr{F}}$ and ${\it\epsilon}_{T}$ are gathered in table 4.

As a conclusion, the dynamics of a scalar field forced by a mean gradient is fully dominated by the decaying isotropic kinetic field. Indeed, the initial scalar large-scale condition ${\it\sigma}_{T}$ never appears in the theoretical decay and growth exponents. The velocity field completely leads the dynamics of the flow. This is consistent with a similar result obtained by De Marinis et al. (Reference De Marinis, Chibbaro, Meldi and Sagaut2013) in HIT, where the passive scalar field experiences a Joule heat production. Therefore, it can be concluded that, with a production mechanism, the velocity field completely dominates the passive scalar dynamics.

With these new algebraic exponents for $K_{T}$ and $K_{\mathscr{F}}$ , one can show that for high or low Reynolds numbers, i.e. when both ${\it\alpha}_{\mathscr{F}}$ and ${\it\alpha}_{T}^{{\it\Lambda}}$ are recovered, the correlation ${\it\rho}_{w{\it\theta}}$ previously investigated is independent of time. Indeed, its algebraic exponent is then ${\it\alpha}_{\mathscr{F}}-({\it\alpha}+{\it\alpha}_{T}^{{\it\Lambda}})/2=0$ . This means that for moderate Reynolds numbers ( $Re_{{\it\lambda}}\sim 100$ ) not high enough to ensure that ${\it\alpha}_{\mathscr{F}}$ and ${\it\alpha}_{T}^{{\it\Lambda}}$ are verified, ${\it\rho}_{w{\it\theta}}$ varies in time. This could be the reason for the scattering of the obtained values for ${\it\rho}_{w{\it\theta}}$ in DNS and experiments, pointed out in § 4.

5.2 Return to isotropy in HITSG

The small-scale RTI is briefly investigated for the HITSG configuration, at the level of the scalar second-order moments. Since only the third component of the scalar gradient is non-zero, this is an axisymmetric configuration, meaning that the scalar anisotropy indicators satisfy $2\unicode[STIX]{x1D623}_{11}^{T}=2\unicode[STIX]{x1D623}_{22}^{T}=-\unicode[STIX]{x1D623}_{33}^{T}$ and $2\unicode[STIX]{x1D60F}_{11}^{(T)}=2\unicode[STIX]{x1D60F}_{22}^{(T)}=-\unicode[STIX]{x1D60F}_{33}^{(T)}$ . In figure 13, the $\unicode[STIX]{x1D623}_{ij}^{T}$ are presented: they become constant in both Saffman and Batchelor turbulence, which is qualitatively the same behaviour as the $\unicode[STIX]{x1D623}_{ij}$ in a sustained shear flow. The spectral anisotropy tensors $\unicode[STIX]{x1D60F}_{ii}^{(T)}$ reveal that small scales of the scalar second-order moments completely return to isotropy for high Reynolds numbers $Re_{{\it\lambda}}\sim 10^{3}$ . Moreover, it has been pointed out in experiments and DNS (Pumir Reference Pumir1994; Danaila et al. Reference Danaila, Le Gal, Anselmet, Plaza and Pinton1999) that, at the level of the scalar third-order moments, some anisotropy remains in the small scales. This is not incompatible with small-scale isotropic second-order moments, as shown recently by Bos (Reference Bos2014) in HITSG.

Figure 13. (a) Diagonal components of $\unicode[STIX]{x1D623}_{ij}^{T}$ , defined in (3.5), for both Saffman (black) and Batchelor (grey) turbulence, as a function of time $t$ . (b) Diagonal components of $\unicode[STIX]{x1D60F}_{ij}^{(T)}$ , defined in (3.3), at $Re_{{\it\lambda}}=10^{3}$ for ${\it\sigma}=2$ , along with the Kolmogorov and scalar integral wavenumbers $k_{{\it\eta}}$ and $k_{T}$ , as a function of the wavenumber $k$ .

Figure 14. Evolution of integrated quantities in HSTSG, as a function of the accumulated anisotropy $St$ , for ${\it\sigma}=2$ and $Re_{{\it\lambda}}(St=60)=4\times 10^{4}$ . (a) Scalar anisotropy indicators $\unicode[STIX]{x1D623}_{ij}^{T}$ defined in (3.5), and scalar shear rapidity $S_{R}^{T}={\it\epsilon}_{T}/(K_{T}S)$ . (b) Ratios ${\it\epsilon}_{\mathscr{F}}/(K_{\mathscr{F}}S)$ and ${\it\epsilon}_{\mathscr{F}}^{S}/(K_{\mathscr{F}}^{S}S)$ .

5.3 Growth laws of the energies in HSTSG

The HSTSG configuration is now considered at very high Reynolds numbers. As mentioned before, the kinetic energy eventually grows exponentially in shear flows at the rate ${\it\gamma}=0.34$ , defined in (4.18), within the present anisotropic EDQNM closure in a completely homogeneous framework (Briard et al. Reference Briard, Gomez, Mons and Sagaut2016; Mons et al. Reference Mons, Cambon and Sagaut2016). In this section, the emphasis is put on the growth of the scalar variance and its interactions with the scalar flux. Some additional results about the passive scalar and the scalar flux are presented, which may be of interest for one-point modelling, such as identifying negligible quantities, and their corresponding mechanisms, at high shear rates.

Figure 15. Evolution of integrated quantities in HSTSG, as a function of the accumulated anisotropy $St$ , for ${\it\sigma}=2$ and $Re_{{\it\lambda}}(St=60)=4\times 10^{4}$ . (a) Exponential growth of the kinetic, scalar, cospectrum and streamwise flux energies. The cospectrum and streamwise flux dissipation rates ${\it\epsilon}_{\mathscr{F}}$ and ${\it\epsilon}_{\mathscr{F}}^{S}$ are displayed in grey. (b) The exponential growth rates ${\it\gamma}_{\mathscr{F}}$ , ${\it\gamma}_{\mathscr{F}}^{S}$ and ${\it\gamma}_{T}$ .

Firstly, the scalar anisotropy tensors $\unicode[STIX]{x1D623}_{ij}^{T}$ are presented in figure 14. As in the HST case without scalar gradient, the scalar indicators reach constant values for large $St$ , and the ratio ${\it\Lambda}/S$ impacts only the short-time dynamics without modifying the asymptotic state. There is a noteworthy similarity with the behaviour of $\unicode[STIX]{x1D623}_{ij}=\unicode[STIX]{x1D619}_{ij}/(2K)-{\it\delta}_{ij}/3$ in shear flows. An interesting feature is that the ratios ${\it\epsilon}_{\mathscr{F}}/(K_{\mathscr{F}}S)$ and ${\it\epsilon}_{\mathscr{F}}^{S}/(K_{\mathscr{F}}^{S}S)$ tend to zero for large $St$ , whereas their kinetic and scalar counterparts do not. This means that the linear effects of shear are preponderant over nonlinear exchanges at the scalar flux level: this is in agreement with figure 8, where the scalar flux transfers of energy are gathered at large scales, dominated by linear mechanisms. The evolution of the scalar, cospectrum and streamwise flux energies are given by

where $K_{T}$ , $K_{\mathscr{F}}$ and $K_{\mathscr{F}}^{S}$ , respectively, have been defined in (3.6), (4.4) and (4.24). Here is what happens simultaneously: The cross-correlation $\unicode[STIX]{x1D619}_{13}$ produces $K_{\mathscr{F}}^{S}$ through the mean scalar gradient ${\it\Lambda}$ . Then $\unicode[STIX]{x1D619}_{13}$ brings energy to the transverse component $\unicode[STIX]{x1D619}_{33}$ thanks to nonlinear redistribution, which causes $K_{\mathscr{F}}$ to grow as well through ${\it\Lambda}$ . Finally, $K_{\mathscr{F}}$ provokes the growth of $K_{T}$ . The main result, displayed in figure 15, is that the scalar variance $K_{T}$ , which was exponentially decreasing in HST, now grows exponentially in HSTSG. Its scalar growth rate is ${\it\gamma}_{T}={\it\gamma}=0.34$ . This is qualitatively in agreement with the experimental work of Ferchihi & Tavoularis (Reference Ferchihi and Tavoularis2002), where the scalar variance and kinetic energy grow exponentially at the same rate in the presence of both mean velocity and scalar gradients. Moreover, both the cospectrum energy $K_{\mathscr{F}}$ and the streamwise flux energy $K_{\mathscr{F}}^{S}$ grow exponentially with the rate ${\it\gamma}$ as well. This is consistent with the constant scalar flux correlations ${\it\rho}_{u_{i}{\it\theta}}$ obtained in § 4 for sufficiently high $St$ . From (5.10), it is possible to determine the analytical expressions of the cospectrum, streamwise flux and scalar exponential growth rates ${\it\gamma}_{\mathscr{F}}$ , ${\it\gamma}_{\mathscr{F}}^{S}$ and ${\it\gamma}_{T}$ . Using the fact that ${\it\epsilon}_{\mathscr{F}}/(K_{\mathscr{F}}S)\rightarrow 0$ and ${\it\epsilon}_{\mathscr{F}}^{S}/(K_{\mathscr{F}}^{S}S)\rightarrow 0$ , one has The agreement between the asymptotic values of these quantities and the value $0.34$ expected is presented in figure 15. Moreover, simulations show that ${\it\gamma}_{\mathscr{F}}$ , ${\it\gamma}_{\mathscr{F}}^{S}$ and ${\it\gamma}_{T}$ do not depend on large-scale initial conditions ${\it\sigma}$ . All the results regarding the exponential rates of the energies and their dissipation rates in HST and HSTSG are gathered in table 5.

5.4 Discussion on the streamwise flux spectrum $\mathscr{F}_{S}(k,t)$ scaling

The scaling of the streamwise flux $\mathscr{F}_{S}$ is now addressed. This is of particular interest since $\langle u_{1}{\it\theta}\rangle$ only exists when both mean velocity and scalar gradients are active, such as in atmospheric flows. Wyngaard & Coté (Reference Wyngaard and Coté1972) first proposed

by assuming that $\mathscr{F}_{S}$ depends on ${\it\epsilon}$ , $k$ , $S$ and ${\it\Lambda}$ (this expression can also be found starting from $\mathscr{F}_{S}\sim {\it\Lambda}{\it\epsilon}^{1/3}k^{-7/3}$ and replacing ${\it\epsilon}^{1/3}$ by its expression as a function of the shear scale ${\it\epsilon}^{1/3}\sim k^{-2/3}S$ ). This is not satisfactory, as there is no dependence on ${\it\epsilon}$ . In a recent paper, Knaus & Pantano (Reference Knaus and Pantano2009) studied reactive and non-reactive scalar flux spectra with DNS and found that a $k^{-7/3}$ range was a satisfactory behaviour. Such a scaling

assumes that ${\it\epsilon}_{\mathscr{F}}^{S}$ is conserved, which is false, since $\mathscr{F}_{S}$ is purely anisotropic and does not exist in HIT. With the present modelling, it is revealed in figure 16 that there is a good agreement with the $k^{-23/9}$ scaling predicted by Bos & Bertoglio (Reference Bos and Bertoglio2007),

Numerically, the constant is found to be $C_{\mathscr{F}}^{S}\simeq 1.5$ . This scaling relies on the assumption that one needs to take into account anisotropy for the form of the dissipation rate ${\it\epsilon}$ . Starting with the scaling of the cospectrum, ${\it\epsilon}^{1/3}$ is replaced by ${\it\epsilon}_{ij}^{1/3}(k)$ , so that

Then, one has ${\it\epsilon}_{ij}(k)E(k)=3{\it\phi}_{ij}(k){\it\epsilon}$ , where ${\it\phi}_{ij}$ is the spherical average of the spectral Reynolds tensor. The $k^{-23/9}$ law is recovered using the scaling proposed by Lumley (Reference Lumley1967), ${\it\phi}_{13}\sim S{\it\epsilon}^{1/3}k^{-7/3}$ .

Figure 16. (a) Kinetic, cospectrum and streamwise flux spectra $E$ , $\mathscr{F}$ and $\mathscr{F}_{S}$ , as a function of the normalized wavenumber $k/k_{L}$ . (b) Associated spectral scaling: the horizontal dashed lines correspond to $-5/3$ , $-7/3$ and $-23/9$ , as a function of the normalized wavenumber $\log (k/k_{L})$ . Both at $St=60$ , where $Re_{{\it\lambda}}=4\times 10^{4}$ , and ${\it\sigma}=2$ .

5.5 Return to isotropy in HSTSG

The small-scale RTI is finally addressed for the HSTSG configuration, at the level of the scalar second-order moments. In figure 17, the scalar anisotropy tensors $\unicode[STIX]{x1D60F}_{ij}^{(T)}$ are presented at large $St$ at two different Reynolds numbers: local isotropy is recovered for the diagonal components of $\unicode[STIX]{x1D60F}_{ij}^{(T)}$ , which are zero at small scales, whereas the extra-diagonal component $\unicode[STIX]{x1D60F}_{13}^{(T)}$ shows that some residual anisotropy still remains at small scales. This confirms that, in the presence of shear, and in agreement with most of the DNS and experiments, some anisotropy persists at the scalar small scales, even at the second-order moments level. Furthermore, the Reynolds number is found to have a non-negligible impact on the small-scale anisotropy, as previously noticed by Garg & Warhaft (Reference Garg and Warhaft1998) for the velocity field: indeed, the reduction of small-scale anisotropy, especially for $\unicode[STIX]{x1D60F}_{13}^{(T)}$ , is noteworthy from $Re_{{\it\lambda}}=2400$ (already higher than the $Re_{{\it\lambda}}$ reached in DNS) to $Re_{{\it\lambda}}=2\times 10^{4}$ . The shear wavenumber $k_{S}=\sqrt{S^{3}/{\it\epsilon}}$ (Corrsin Reference Corrsin1958) is displayed as well: for wavenumbers $k>k_{S}$ , nonlinear effects are dominant, consistent with the RTI of small scales, whereas for $k<k_{S}$ , linear effects are stronger and carry most of the anisotropy. A last remark is that the scalar gradient seems to smooth the scalar large-scale anisotropy: indeed, for shear-driven flows without mean scalar gradient, anisotropy is gathered around the scalar peak of energy $k_{T}$ and is weaker in the infrared range (see figure 6). Whereas for HSTSG, anisotropy progressively increases from moderate to large scales.

Therefore, the conclusion about the RTI of passive scalar small scales is twofold. Firstly, and consistently with what is observed for shear flows, an increasing Reynolds number is found to reduce small-scale anisotropy for the scalar second-order moments. Nevertheless, some anisotropy remains even at the highest Reynolds numbers reached here ( ${\sim}10^{4}$ ). Secondly, simulations indicate that, at high Reynolds numbers, scalar small scales are still slightly anisotropic in shear-driven flows (with or without an additional mean scalar gradient), whereas they completely return to isotropy when there is only a mean scalar gradient.

Figure 17. Scalar anisotropy descriptors $\unicode[STIX]{x1D60F}_{ij}^{(T)}$ , defined in (3.3), along with the Kolmogorov, shear and scalar integral wavenumbers $k_{{\it\eta}}$ , $k_{S}=\sqrt{S^{3}/{\it\epsilon}}$ and $k_{T}$ , for ${\it\sigma}=2$ , at two different Reynolds numbers, as a function of the normalized wavenumber $k/k_{T}$ for: (a)  $Re_{{\it\lambda}}(St=40)=2400$ ; and (b)  $Re_{{\it\lambda}}(St=50)=2\times 10^{4}$ . The zoom represents the small scales before $k_{{\it\eta}}$ where anisotropy persists.