4.1. Numerical set-up

A standard Rogallo (Reference Rogallo1981) code is used to perform a series of experiments on decaying triply periodic homogeneous turbulence in a $(2{\rm \pi} )^3$ box. Equations (3.6) are solved with the dynamic Smagorinsky model described in the previous section. The three velocity components are projected on a Fourier basis, and the nonlinear terms are calculated using a fully dealiased pseudo-spectral method with a 2/3 truncation rule (Canuto et al. Reference Canuto, Hussaini, Quarteroni and Thomas2012). We denote the wavenumber vector by $\boldsymbol {k}$, and its magnitude by $k=|{\boldsymbol {k}}|$. The number of physical points in each direction before dealiasing is $N=128$, so that the highest fully resolved wavenumber is $k_{max}=42$. An explicit third-order Runge–Kutta is used for temporal integration, and the time step is adjusted to a constant Courant–Friedrichs–Lewy number equal to $0.2$. A Gaussian filter whose Fourier expression is

(4.1)\begin{equation} G(k;\varDelta)=\exp(-\varDelta^2k^2/24), \end{equation}

where $\varDelta$ is the filter width (Aoyama et al. Reference Aoyama, Ishihara, Kaneda, Yokokawa, Itakura and Uno2005), is used to evaluate $C_s$ in (3.5). The equations are explicitly filtered at the cutoff filter $\bar {\varDelta }$, and at the test filter $\tilde \varDelta =2\bar {\varDelta }=\varDelta _{g}\sqrt {6}$, where $\varDelta _{g}=2{\rm \pi} /N$ is the grid resolution before dealiasing. Although the explicit cutoff filter is not strictly necessary, it is used for consistency and numerical stability.

The energy spectrum is defined as

(4.2)\begin{equation} E(k;\boldsymbol u)=2{\rm \pi} k^2\langle \hat{u}_i \hat{u}_i^*\rangle_k, \end{equation}

where $\hat {u}_i(\boldsymbol {k})$ are the Fourier coefficients of the $i$-th velocity component, the asterisk is complex conjugation and $\langle \cdot \rangle _k$ denotes averaging over shells of thickness $k\pm 0.5$. The kinetic energy per unit mass is $\mathcal {E}=\langle {u_i u_i}\rangle /{2}=\sum _k E(k,t)$, where $\sum _k$ denotes summation over all wavenumbers.

Meaningful units are required to characterise the simulations. The standard reference for the small scales in Navier–Stokes turbulence are Kolmogorov units, but they are not applicable here because of the absence of molecular viscosity. Instead we derive pseudo-Kolmogorov units using the mean eddy viscosity, $\langle \nu _s\rangle = C_s\bar {\varDelta }^2 \langle |S| \rangle$, and the mean sub-grid energy transfer at the cutoff scale, $\langle \epsilon _s\rangle = C_s\bar {\varDelta }^2 \langle |S|^3 \rangle$. The time and length scale derived from these quantities are $\tau _s=(\langle \nu _s\rangle /\langle \epsilon _s \rangle )^{1/2}$ and $\eta _s=(\langle \nu \rangle ^3_s/\langle \epsilon _s \rangle )^{1/4}$, respectively. We find these units to be consistent with the physics of turbulence, such that the peak of the enstrophy spectrum, $k^2E(k)$, is located at ${\sim }25\eta _s$. In direct numerical simulations, this peak lies at ${\sim }20\eta$.

The larger eddies are characterised by the integral length, velocity and time scales (Batchelor Reference Batchelor1953), respectively defined as

(4.3)\begin{equation} L={\rm \pi} \frac{\displaystyle\sum_k E(k,t)/k}{\displaystyle\sum_k E(k,t)}, \end{equation}

$u'=\sqrt {2\mathcal {E}/3}$ and $T=L/u'$. The Reynolds number based on the integral scale is defined as $Re_L=u' L/\langle \nu _s\rangle$, and the separation of scales in the simulation is represented by the ratio $L/\eta _s$.

In the following we study the energy flux both in Fourier space, averaged over the computational box, and locally in physical space.

The energy flux across the surface of a sphere in Fourier space with wavenumber magnitude $k$ can be expressed as

(4.4)\begin{equation} \varPi(k)=\sum_{q<k} 4{\rm \pi} q^2{Re}\langle \hat{u}^*_i(\boldsymbol{q})\hat{\mathcal{D}}_i(\boldsymbol{q})\rangle_q, \end{equation}

where $\mathcal {D}_i=u_j\partial _j u_i$ and a negative $\varPi$ denotes energy flowing to larger scales.

We use two different definitions of the local interscale energy flux in physical space,

(4.5)\begin{equation} \left.\begin{aligned} &\varSigma(\boldsymbol{x},t; \varDelta) = \tau_{ij}{S}_{ij}, \\ &\varPsi(\boldsymbol{x},t; \varDelta) ={-}u_i\partial_j \tau_{ij}, \end{aligned}\right\}\end{equation}

where the velocity field, the rate-of-strain tensor and $\tau _{ij}$ are filtered with (4.1) at scale $\varDelta$. Both are standard quantities in the analysis of the turbulence cascade (Meneveau & Lund Reference Meneveau and Lund1994; Borue & Orszag Reference Borue and Orszag1998; Aoyama et al. Reference Aoyama, Ishihara, Kaneda, Yokokawa, Itakura and Uno2005; Cardesa et al. Reference Cardesa, Vela-Martín, Dong and Jiménez2015), and are related by $\varPsi =\varSigma - \partial _j (u_i\tau _{ij})$. The second term on the right-hand side of this relation is the divergence of an energy flux in physical space, which has zero mean over the computational box, so that $\langle \varSigma \rangle =\langle \varPsi \rangle$. Positive values of $\varSigma$ and $\varPsi$ indicate that the energy is transferred towards the small scales.

Because the decaying system under study is statistically unsteady, averaging over an ensemble of many realisations is required to extract time-dependent statistics. This is generated using the following procedure. All the flow fields in the ensemble share an initial energy spectrum, derived from a forced statistically stationary simulation, and an initial kinetic energy per unit mass. Each field is prepared by randomising the phases of $\hat {u}_i$, respecting continuity, and integrated for a fixed time $t_{start}$ up to $t_0=0$, where a fully turbulent state is deemed to have been reached, the experiment begins, and statistics start to be compiled. The initial transient, $t_{start}$, common to all the elements in the ensemble, is chosen so that $t_0$ is beyond the time at which dissipation reaches its maximum, and the turbulent structure of the flow is fully developed. Quantities at $t_0$ are denoted by a ‘0’ subindex. Following this procedure, we generate an ensemble of $N_s=2000$ realisations, which are evolved on graphic processing units. The statistics below are compiled over all the members of this ensemble.

It is found that both the small- and large-scale reference quantities defined above vary little across an ensemble prepared in this way, with a standard deviation of the order of 5 % with respect to their mean.

4.2. Turbulence with a reversible model

The basic numerical experiment is conducted as in Carati et al. (Reference Carati, Winckelmans and Jeanmart2001). After preparing the ensemble of turbulent fields at $t=t_0=0$, they are evolved for a fixed time $t_{inv}$, during which some of the initial energy is exported by the model to the unresolved scales. The simulations are then stopped and the sign of the three components of velocity reversed, $\boldsymbol {u}\rightarrow -\boldsymbol {u}$. The new flow fields are used as the initial conditions for the second part of the run from $t_{inv}$ to $2t_{inv}$, during which the flows evolve back to their original state, recovering their original energy and the value of other turbulent quantities.

It is found that $\eta _s$ only changes by 5 % during the decay of the flow, while $\tau _s$ increases by a factor of $1.5$ from $t=0$ to $t_{inv}$. On the other hand, the large-scale quantities $L$ and $u'$ vary substantially as the flow decays. The main parameters of the simulations, and the relation between large- and small-scale quantities, are given in table 1.

Table 1. Main parameters of the simulation averaged over the complete ensemble at time $t_0=0$ and $t_{inv}$. See the text for definitions.

Figure 1(a,b) presents the evolution of the kinetic energy $\mathcal {E}$ and of $C_s$ as a function of time for a representative experiment. We observe a clear symmetry with respect to $t_{inv}$ in both quantities, and we have checked that this symmetry also holds for other turbulent quantities, such as the sub-grid energy fluxes and the skewness of the velocity derivatives. Quantities that are odd with respect to the velocity display odd symmetry in time with respect to $t_{inv}$, and even quantities display even symmetry. In particular, we expect $\boldsymbol {u}(2t_{inv})\simeq -\boldsymbol {u}(0)$.

Figure 1. Evolution of (a) mean kinetic energy, normalised with the initial energy $\mathcal {E}_0$, and (b) $C_s\bar \varDelta ^2$ as a function of time (arbitrary units). The dashed line in (a,b) is $t=t_{inv}$, and the dotted lines in (a) are $t=t_{stats}$, used in § 6. (c) Energy spectrum $E(k)$. ——, $t=0$ and $t=2t_{inv}$; – – –, $t=t_{inv}$. The diagonal line is $E(k)\propto k^{-5/3}$. (d) Spectrum of the error between the velocity fields at $t=2t_{inv}$ and $t=0$, as a function of wavenumber. The dashed line denotes $k^4$.

Figure 1(c) shows the energy spectrum $E(k)$ at times $t=0$, $t_{inv}$ and $2t_{inv}$. Although there are no observable differences between the initial and final spectra, we quantify this difference in scale using the spectrum of the error between the velocity field at $t=0$ and minus the velocity field at $2t_{inv}$,

(4.6)\begin{equation} \rho(k)=\frac{2E(k;\boldsymbol{u}(0)+\boldsymbol{u}(2t_{inv}))}{E(k;\boldsymbol{u}(0)) + E(k;\boldsymbol{u}(2t_{inv}))}, \end{equation}

which is presented in figure 1(d). The spectrum of the error is of the order of $0.01$ in the cutoff wavenumber and decreases as $k^4$ with the wavenumber, confirming that the flow fields in the forward and backward evolution are similar, except for the opposite sign and minor differences in the small scales. This suggests that the energy cascade is microscopically reversible in the inertial scales. Sustained reverse cascades are possible in the system, even if only direct ones are observed in practice, supporting the conclusion that the one-directional turbulence cascade is an entropic (probabilistic) effect, unrelated to the presence of an energy sink at the small scales.

This conclusion is in agreement with the empirical evidence that the kinetic energy dissipation in turbulent flows is independent of the Reynolds number (Sreenivasan Reference Sreenivasan1984). This phenomenon, known as the dissipative anomaly, exposes the surrogate nature of kinetic energy dissipation, which is controlled by large-scale dynamics through the energy cascade process (Taylor Reference Taylor1935; Kolmogorov Reference Kolmogorov1941). In agreement with the dissipative anomaly, we will present evidence that the energy dissipation is a consequence, rather than a cause, of the energy cascade.

4.3. The reverse cascade without the model

The validity of the conclusions in the previous section depends on the ability of the model system to represent turbulence dynamics. In particular, since the object of our study is the turbulence cascade rather than the SGS model, it is necessary to assess whether the reversibility properties of the system, and the presence of a sustained backwards energy cascade, stem from the SGS model, or whether intrinsic turbulent mechanisms are involved. In this subsection we show that the model injects energy at the smallest resolved scales, but that the energy travels back to the large scales due to inertial mechanisms. In the next subsection we further demonstrate that the inverse cascade can exist for some time even for irreversible formulations of the sub-grid model.

If the SGS model contributed substantially to the presence of a sustained inverse cascade, its removal during the backward leg of the simulation would have an immediate effect on the system, resulting in the instantaneous breakdown of the reverse cascade. This is shown not to be the case by an experiment in which the model is removed at $t_{inv}$ by letting $C_s=0$. The system then evolves according to the Euler equations

(4.7)\begin{equation} \left. \begin{aligned} &\partial_t u_i+u_j\partial_j u_i={-}\partial_i p, \\ &\partial_i u_i=0, \end{aligned}\right\}\end{equation}

which are conservative and include only inertial forces. The resulting evolution of the inertial energy flux across Fourier scales, $\varPi (k,t)$, defined in (4.4), is displayed in figure 2(a) as a function of wavenumber and time.

Figure 2. (a) Energy flux $\varPi (k)$ as a function of the time after removing the SGS model, $t/T_0$: $\square$, solid black line, $0$; $\bullet$, solid red line, $0.1$; $\blacktriangleright$, solid blue line, $0.15$; $\blacktriangle$, solid magenta line, $0.25$. Energy flux normalised with $\epsilon _0={u'_0}^3/L_0$. (b) Reversal time $t_{rev}/T_0$ at which $\varPi =0$ as a function of wavenumber $kL_0$. Initial condition before the inversion of the velocities and the removal of the model obtained with: ——, red, reversible dynamic SGS model ($N=128$); – – –, irreversible spectral SGS model ($N=512$). The solid black line corresponds to $t_{rev}/T_0=2{(kL_0)}^{-2/3}$.

At $t_{inv}$, energy fluxes are negative at all scales and energy flows towards the large scales. Shortly after, at $t-t_{inv}=0.1T_0$, the inverse cascade begins to break down at the small scales, while it continues to flow backwards across most wavenumbers. This continues to be true at $t=0.15T_0$, where the direct cascade at the smallest scales coexists with a reverse cascade at the larger scales. The wavenumber separating the two regimes moves progressively towards larger scales, and the whole cascade becomes direct after $t\simeq 0.25T_0$.

Figure 2(b) shows the dependence on the wavenumber of the time $t_{rev}(k)$ at which $\varPi (k,t_{rev})=0$. It is significant that, for wavenumbers that can be considered inertial, $kL_0\gg 1$, it follows that $t_{rev}/T_0\approx 2(k L_0)^{-2/3}$, which is consistent with the Kolmogorov (Reference Kolmogorov1941) self-similar structure of the cascade. If we assume a Kolmogorov spectrum, $E(k)=C_K \epsilon _s^{2/3} k^{-5/3}$, where $C_K\approx 1.5$ is the Kolmogorov constant (Pope Reference Pope2001) and $\epsilon _s$ is the energy flux through wavenumber $k$, which is equal to the energy dissipation, the energy above a given wavenumber is

(4.8)\begin{equation} \mathcal{E}(k)=\int^{\infty}_{k} E(q)\,\mathrm{d}q \approx \tfrac{3}{2}C_K \epsilon_s^{2/3} k^{{-}2/3}. \end{equation}

When the SGS model is removed, an inverse cascade can only be sustained by drawing energy from this reservoir. Independently of whether the cascade can maintain its integrity from the informational point of view, the time it takes for a flux $-\epsilon _s$ to drain (4.8) can be expressed as

(4.9)\begin{equation} \frac{t_{max}}{T_0}\approx \frac{3}{2} C_K \left( \frac{u'^3}{L_0\epsilon_s} \right)^{1/3} (k L_0)^{{-}2/3} \approx 2 (k L_0)^{{-}2/3}, \end{equation}

where we have used values from table 2. As shown in figure 2(b), this approximation is only 50 % larger than the observed reversal time for the initial conditions prepared with the reversible SGS model.

Table 2. Auto-correlation coefficient among scales of $\varSigma$ and $\varPsi$, evaluated for: ‘All’, the complete field; ‘$+$’, conditioned to positive energy transfer events at $5\varDelta _g$; ‘$-$’, conditioned to negative energy transfer events at $5\varDelta _g$.

The simplest conclusion is that the model just acts as a source to provide the small scales with energy as they are being depleted by the flux towards larger eddies. If this source is missing, the predominant forward cascade reappears, but the inertial mechanisms are able to maintain the reverse cascade process as long as energy is available.

4.4. The effect of irreversible models

We have repeated the same experiment for an irreversible spectral SGS model (Métais & Lesieur Reference Métais and Lesieur1992) with a higher resolution $(N=512)$. Results of $t_{rev}(k)$ for this experiment are shown by the dashed line in figure 2(b). The behaviour is similar to the experiment in the previous section, but the wider separation of scales in this simulation allows the preservation of an inverse cascade for times of the order of the integral time. This experiment confirms that the microscopic reversibility of the inertial scales is independent of the type of dissipation, and holds as long as the information of the direct cascade process is not destroyed by the LES model. This experiment suggests that microscopic reversibility should also hold in the inertial range of fully developed Navier–Stokes turbulence. Note that the good agreement of the reversal time with the energy estimate in (4.9) implies that $t_{rev}$ is very close to the maximum possible for a given energy.

5.1. The geometry of phase space

In the first place we confirm, by perturbing the inverse cascade, that in the neighbourhood of each inverse evolution there exist a dense distribution of phase-space trajectories that also display inverse dynamics. Each inverse trajectory is ‘unstable’, in the sense that it is eventually destroyed when perturbed, but inverse dynamics is easily found considerably far from the original phase-space trajectories, even for distances of the order of the size of the accessible phase space, suggesting that inverse and direct trajectories lie in separated regions of phase space, rather than mostly being intertwined in the same neighbourhood. Hereafter we refer to the phase space of the Fourier coefficients, and consider, unless stated, experiments with the reversible SGS model.

In this experiment we generate a backwards initial condition by integrating the equations of motion until $t_{inv}$, changing the sign of the velocities $\boldsymbol {u}^\iota =-\boldsymbol {u}(t_{inv})$, and introducing a perturbation, $\delta \boldsymbol {u}$. The perturbed flow field, $\boldsymbol {u}^p=\boldsymbol {u}^\iota +\delta \boldsymbol {u}$, is evolved and compared with the unperturbed trajectory, $\boldsymbol {u}^\iota (t)$. We choose the initial perturbation field in Fourier space as

(5.1)\begin{equation} \delta \hat{u}_i({\boldsymbol{k}})= \hat{u}^\iota_i({\boldsymbol{k}})\cdot \mu \exp(\textrm{i} \phi), \end{equation}

where $\phi$ is a random angle with uniform distribution between $0$ and $2{\rm \pi}$ and $\mu$ is a real parameter that sets the initial energy of the perturbation. The same $\phi$ is used for the three velocity components, so that incompressibility also holds for $\boldsymbol {u}^p$. Figure 3(a,b) shows the evolution of the mean kinetic energy of the perturbed field, $\mathcal{E}_p=\langle u^p\rangle^2/2$, and of the perturbation field, $\delta \mathcal {E}=\langle \delta \boldsymbol {u}^2\rangle /{2}$, for $\mu =0.1$, $0.01$ and $0.001$. The energy of the perturbed fields increases for a considerable time in the three cases, demonstrating the presence of inverse energy transfer and negative mean eddy viscosity, $C_s<0$. Even under the strongest initial perturbation, $\mu =0.1$, the inverse cascade is sustained for approximately $0.5T_0$. Note that $\delta \mathcal {E}$ grows exponentially at a constant rate for some time under the linearised dynamics of the system, indicating that perturbations tend to the most unstable linear perturbation (Benettin et al. Reference Benettin, Galgani, Giorgilli and Strelcyn1980), which is independent of the form of the initial perturbation.

Figure 3. (a) Evolution of the energy of the perturbed fields $\mathcal {E}_p$ and (b) of $\delta \mathcal {E}=1/2\langle \delta {\boldsymbol {u}}^2\rangle$, normalised with $\mathcal {E}_0$, for: ——, $\mu =0.1$; – – –, red, $\mu =0.01$; $-\cdot -$, blue, $\mu =0.001$; $\cdots \cdots \cdot$, energy of the unperturbed trajectory, $\mathcal {E}(t)/\mathcal {E}_0$. (c) Perturbed trajectories for $\mu =0.1$ and $\mu =0.01$, represented in the energy-dissipation space. Symbols as in (a). The solid circle represents the beginning of the inverse trajectory at $t_{inv}$, the empty circles represent the end of the perturbed trajectories at $2t_{inv}$, and the solid square represents the end of the direct trajectory at time $t_{inv}$. The dotted lines and arrows represent the direct and inverse trajectories. (d) Evolution of the probability density function of $\epsilon _s$ as a function of time starting from random initial conditions at $t_{start}$. From left to right: $(t-t_{start})/T_0=: 0.0$; $0.008$; $0.017$; $0.033$; $0.067$; $1.0$. $\epsilon _s$ normalised with the ensemble average at $T_0$.

As shown in figure 3(b), by the time the inverse cascade is destroyed and the energy of the perturbed fields has evolved to a maximum and $C_s=0$, the energy of the perturbation is comparable to the total energy in all cases, $\delta \mathcal {E}\sim \mathcal {E}$, which indicates that the states that separate inverse and direct dynamics are considerably far in phase space from the unperturbed trajectories, at least when the initial perturbations are small. If we estimate distances as the square root of the energy ($L_2$ norm of the velocity field), the maximum of $\mathcal {E}_p$ is located approximately at $\delta \mathcal {E}^{1/2}\simeq 0.4\mathcal {E}_0^{1/2}$ from the unperturbed phase-space trajectory. In figure 3(c) we represent the evolution of the perturbed trajectories in the dissipation-energy space. This plot intuitively shows that a typical direct trajectory is located far apart in phase space from inverse trajectories, and that perturbed inverse trajectories must first cross the set of states with zero dissipation, $C_s=0$, before developing a direct cascade.

In figure 3(d) we show the evolution of the probability density function of the sub-grid energy transfer, $\epsilon _s$, in the direct evolution of the ensemble, starting at $t_{start}$, when the initial condition is a random field with a prescribed energy spectrum. The initial probability distribution of $\epsilon _s$ is symmetric and has zero mean, indicating that the probabilities of direct and inverse evolutions are similar. As the system evolves in time and develops the turbulence structure, the probability of finding a negative value of $\epsilon _s$ decreases drastically. At $t-t_{start}=0.008T_0$, we do not observe any negative values of $\epsilon _s$ in the ensemble. The standard deviation of $\epsilon _s$ is small compared with the mean at $t-t_{start}=T_0$, which indicates the negligible probability of observing negative values of $\epsilon _s$ after a time of the order of an eddy-turnover time.

Although decaying turbulent flows have no attractor in the strict sense, i.e. a set of phase-space points to which the system is attracted and where it remains during its evolution, we use the term here to denote the set of phase-space trajectories with a fully turbulent structure and a direct energy cascade. As we have shown, this set of states is in fact ‘attracting’, but the flow does not remain indefinitely turbulent due to its decaying nature. Let us note that it is possible to construct a reversible turbulent system with a proper turbulent attractor by including a reversible forcing in the large scales, for instance, a linear forcing that keeps the total energy of the system constant in time.

Inverse evolutions sustained in time are almost only accessible from forward-evolved flow fields with changed sign. These trajectories lie in the antiattractor, which is constructed by applying the transformation $\boldsymbol {u}\rightarrow -\boldsymbol {u}$ to the turbulent attractor. However, we have shown in this analysis that inverse trajectories exist in a larger set of states outside the antiattractor, which can be escaped by perturbing the reversed flow fields. The destruction of the inverse cascade under perturbations reflects the negligible probability of inverse evolutions, and indicates that they can only be accessed temporarily. Moreover, the small standard deviation of the sub-grid energy transfer with respect to the mean in the turbulent attractor reflects the negligible probability that the system spontaneously escapes the attractor and develops and inverse energy cascade.

The reversible turbulent system under study provides access to trajectories outside the turbulence attractor, which contains almost exclusively direct cascades. By studying these inverse trajectories, we identify the physical mechanisms that lead to the prevalence of direct over inverse cascades.

5.2. The structure of local energy fluxes in physical space

Up to now we have only dealt with the mean transfer of energy over the complete domain. In this section we characterise the energy cascade as a local process in physical space. We study here two markers of local energy transfer in physical space, $\varSigma (\boldsymbol {x},t;\varDelta )$ and $\varPsi (\boldsymbol {x},t;\varDelta )$, previously defined in (4.5), calculated at filter scales $\varDelta =5\varDelta _{g}$ and $\varDelta =10\varDelta _{g}$. Figure 4(a,b) show the probability distribution of these quantities in the direct evolution. We observe wide tails in the probability distribution of the two quantities, and skewness towards positive events, which is more pronounced for $\varSigma$ than for $\varPsi$. Energy fluxes change sign under $\boldsymbol {u}\rightarrow -\boldsymbol {u}$, and, as shown in figure 4, their odd-order moments are in general non-zero and have a definite sign in turbulent flows. This sign denotes statistical irreversibility, i.e, the privileged temporal direction of the system in an out-of-equilibrium evolution. The spatial average, which marks the direction of the cascade, is an important example, $\langle \varSigma \rangle =\langle \varPsi \rangle >0$.

Figure 4. (a,b) Probability density function of local markers of the energy cascade without volume averaging: (a) $\varSigma$; (b) $\varPsi$. Symbols correspond to: $\blacktriangle$, dashed red line, $\varSigma (5\varDelta _{g})$; $\blacksquare$, dashed blue line, $\varSigma (10\varDelta _{g})$; $\blacktriangleright$, solid magenta line, $\varPsi (5\varDelta _{g})$; $\bullet$, solid black line, $\varPsi (10\varDelta _{g})$. The dashed–dotted lines represent $-\langle \psi \rangle /\sigma _V$ for quantities calculated at $5\varDelta _{g}$, and the open symbols represent the probability density function of $\langle \varSigma (5\varDelta _{g})\rangle _V$ and $\langle \varPsi (5\varDelta _{g})\rangle _V$ for $V=(16\varDelta _{g})^3$.

To characterise the non-local spatial properties of the local energy fluxes, we define the correlation coefficient between two scalar fields, $\psi$ and $\zeta$, as

(5.2)\begin{equation} \mathcal{S}_{\psi,\zeta}(\Delta \boldsymbol{x})=\frac{\langle \psi'(\boldsymbol{x} + \Delta \boldsymbol{x}) \zeta'(\boldsymbol{x})\rangle}{\sqrt{\langle \psi'^2\rangle \langle\zeta'^2\rangle}}, \end{equation}

where we subtract the spatial average from quantities marked with primes. Due to isotropy, $\mathcal {S}$ only depends on the spatial distance, $\Delta x=|\Delta \boldsymbol {x}|$. We also define the auto-correlation coefficient $\mathcal {S}_{\psi ,\psi }(\Delta x)$ and the auto-correlation length as

(5.3)\begin{equation} \ell_\psi=\int^\infty_0 \mathcal{S}_{\psi,\psi}(\xi){\,\textrm{{d}}} \xi, \end{equation}

which measures the typical length over which $\psi$ decorrelates, and is related to the typical size of events in $\psi$. The correlation lengths of $\varSigma$ and $\varPsi$ are proportional to the filter size $\varDelta$, $\ell _\varSigma =0.47 \varDelta$ and $\ell _\varPsi =0.51 \varDelta$.

We use $\langle \varSigma \rangle _V$ and $\langle \varPsi \rangle _V$ to study the spatial structure of energy fluxes, where $\langle \cdot \rangle _V$ represents the volume-averaging operation over a sphere of volume $V$. The standard deviation of the probability distribution of the volume-averaged fields,

(5.4)\begin{equation} \sigma_V=\sqrt{\langle\left(\langle\psi\rangle_V-\langle\psi\rangle\right)^2\rangle}, \end{equation}

is shown in figure 5(a) as a function of the averaging volume, where $\sigma _0=(\langle \psi ^2\rangle - \langle \psi \rangle ^2)^{1/2}$ is the standard deviation of the test field without volume averaging. When the averaging volume is sufficiently large, $V/\ell ^3\gtrsim 10^3$, the standard deviation becomes inversely proportional to the square root of $\mathcal {N}=V/\ell ^3$ in all cases, where $\mathcal {N}$ is a measure of the number of independent energy transfer events within the averaging volume. These results suggest statistical independence of the events within the averaging volumes and some degree of space locality in the energy cascade.

Figure 5. (a) Plot of $\sigma _V$ as a function of $V/\ell ^3$, where the dashed line represents $\sigma _V\sim \mathcal {N}^{-1/2}$. Symbols correspond to: $\blacktriangle$, dashed red line, $\varSigma (5\varDelta _{g})$; $\blacksquare$, dashed blue line, $\varSigma (10\varDelta _{g})$; $\blacktriangleright$, solid magenta line, $\varPsi (5\varDelta _{g})$; $\bullet$, dashed black line, $\varPsi (10\varDelta _{g})$. (b) Correlation coefficient between $\varSigma$ and $\varPsi$ as a function of the averaging volume: $\blacksquare$, solid dash red, $\mathcal {S}_{\langle \varSigma \rangle _V,\langle \varPsi \rangle _V}(0)$ at $5\varDelta _{g}$; $\bullet$, solid blue line, $\mathcal {S}_{\langle \varSigma \rangle _V,\langle \varPsi \rangle _V}(0)$ at $10\varDelta _{g}$. Error bars are calculated from the average standard deviation of $\log \mathcal {P}$ when partitioning the complete dataset in four subsets. (c) Interscale auto-correlation coefficient as a function of the averaging volume: $\blacktriangle$, solid black line, $\mathcal {R}_{\langle \varSigma \rangle _V}(5\varDelta _{g},10\varDelta _{g})$; $\blacktriangleright$, solid magenta line, $\mathcal {R}_{\langle \varPsi \rangle _V}(5\varDelta _{g},10\varDelta _{g})$; $\blacktriangle$, dashed black line, $\mathcal {R}^-_{\langle \varSigma \rangle _V}(5\varDelta _{g},10\varDelta _{g})$; $\blacktriangleright$, solid dash magenta, $\mathcal {R}^-_{\langle \varPsi \rangle _V}(5\varDelta _{g},10\varDelta _{g})$. (d) Asymmetry function, $\mathcal {P}$, as a function of $\langle \psi \rangle ^2/\sigma ^2_V$. Symbols as in (a). The dashed line is the exact solution of $\mathcal {P}$ for a Gaussian, which approaches $\log \mathcal {P}=1/2{\langle \psi \rangle ^2}/{\sigma _V^2}$ for large ${\langle \psi \rangle ^2}/{\sigma _V^2}$. (e, f) Probability density function of local markers of the energy cascade averaged at scale $V^{1/3}=1.2L_0=32\varDelta _g$: (e) $\varSigma$; ( f) $\varPsi$. The dashed–dotted lines represent $-\langle \psi \rangle /\sigma _V$ for quantities calculated at $5\varDelta _{g}$. Symbols as in (a).

Despite this locality, $\varSigma$ and $\varPsi$ are very different pointwise, which might suggest that local energy fluxes are not uniquely defined, and that only global averages are robust with respect to the particular definition of fluxes, $\langle \varSigma \rangle =\langle \varPsi \rangle$. In fact, the correlation between $\varSigma$ and $\varPsi$ for $\Delta x=0$ is low, $\mathcal {S}_{\varSigma ,\varPsi }(0)\sim 0.05$, at both filter scales, $5\varDelta _{g}$ and $10\varDelta _{g}$. However, when these quantities are averaged, the correlation coefficient increases rapidly with the averaging volume. In figure 5(b) we show the dependence of $\mathcal {S}_{\langle \varSigma \rangle _V, \langle \varPsi \rangle _V}(0)$ on $V$. For averaging volumes $V\sim (4\ell )^3$, the correlation increases to approximately $0.7$ for both filter widths, indicating that $\varSigma$ and $\varPsi$ are similar when averaged over volumes of the order of their cubed correlation length.

5.3. The structure of local energy fluxes in scale space

We extend this analysis to scale space by considering the correlation of the energy fluxes at different scales, defined as

(5.5)\begin{equation} \mathcal{R}_\psi(\varDelta_1,\varDelta_2)=\frac{\langle \psi'(\boldsymbol{x};\varDelta_1)\psi'(\boldsymbol{x};\varDelta_2)\rangle}{\sqrt{\langle (\psi'(\varDelta_1))^2\rangle \langle(\psi'(\varDelta_2))^2\rangle}}, \end{equation}

and calculate the same quantity conditioning the averages to positive or negative energy transfer events at scale $\varDelta _1$. We denote the conditional correlations by $\mathcal {R}^+$ when conditioning to $\psi (\varDelta _1)>0$, and $\mathcal {R}^-$ when conditioning to $\psi (\varDelta _1)<0$. These interscale auto-correlation coefficients are presented in table 2. We find correlation values of approximately $0.7$ for scale increments of $2$ in both quantities. Table 2 also includes results for $\varDelta _2=20\varDelta _g$. It shows that all fluxes substantially decorrelate for $\varDelta _2/\varDelta _1=4$ (Cardesa et al. Reference Cardesa, Vela-Martín and Jiménez2017), and that this is especially true for the backscatter.

Figure 5(c) shows the dependence of $\mathcal {R}_{\langle \varSigma \rangle _V}(5\varDelta _{g}, 10\varDelta _{g})$ and $\mathcal {R}_{\langle \varPsi \rangle _V}(5\varDelta _{g},10\varDelta _{g})$ on the averaging volume $V$. The interscale auto-correlation increases with the averaging volume when we consider all direct events (not shown). On the other hand, the correlations conditioned to backscatter increase up to volumes of the order of $V^{1/3}\approx 2\ell$, and then decrease. Beyond $V^{1/3}\approx 4\ell$, the number of averaged inverse energy transfer events is not enough to compute reliable correlations.

5.4. Physical-space estimates of the probability of inverse cascades

We have shown in § 5.1 that the probability of spontaneously observing an inverse cascade over extended regions of a turbulent flow is negligible, making it impractical to quantify the probability of such evolutions by direct observation. However, as shown in § 5.2, we frequently observe local inverse energy transfer events over restricted regions of physical space. We attempt to estimate the probability of observing an inverse cascade using the integral asymmetry function,

(5.6)\begin{equation} \mathcal{P}=\frac{P(\langle\psi\rangle_{V}>0)}{P(\langle\psi\rangle_{V}<0)}, \end{equation}

which compares the probability of observing a direct to an inverse cascade in a volume $V$. This approach follows the methodology of the fluctuation relations (Evans & Searles Reference Evans and Searles2002) and its local versions (Ayton, Evans & Searles Reference Ayton, Evans and Searles2001; Michel & Searles Reference Michel and Searles2013).

We evaluate the integral asymmetry function of $\varSigma$ and $\varPsi$ for different averaging volumes, and display it as a function of $\langle \psi \rangle ^2/\sigma _V^2$ in figure 5(d). The statistical independence of the energy transfer events reported in § 5.2 suggests that, when the averaging volume is large enough, the integral asymmetry function should behave as that of a Gaussian distribution with non-zero mean, which for large $\langle \psi \rangle ^2/\sigma _V^2$ is

(5.7)\begin{equation} \log \mathcal{P}\simeq2\langle \psi \rangle^2/\sigma_V^2. \end{equation}

If we assume that the number of independent events in a volume $V$ is $\mathcal {N}\sim V/\ell ^3$, we can conclude that $\log P \sim \mathcal {N} \sim V/\ell ^3$, and the probability of direct over inverse cascades increases exponentially with the number of independent energy transfer events considered. The exponential dependence in (5.7) is well satisfied in figure 5(d), but the prefactor is not, and $\log \mathcal {P}>2\langle \psi \rangle ^2/\sigma _V^2$. Since the prefactor is a property of the Gaussian distributions that arise from the application of the central limit theorem, this failure can be traced to the strong dependence of $\mathcal {P}$ on the negative tails of the statistical distributions of $\varSigma$ and $\varPsi$, which are not Gaussian. The fast growth of $\log \mathcal {P}$ in figure 5(d) signals that the probability of inverse-transfer events decreases with volume averaging faster than what could be expected for independent events, in agreement with the evidence in figure 5(c) that inverse energy transfer is shallow in scale space. In figure 4(a,b) we show the probability density function of $\langle \varSigma (5\varDelta _g) \rangle _V$ and $\langle \varSigma (5\varDelta _g) \rangle _V$ for $V=(16\varDelta _g)^3$. The negative tails decrease considerably with the averaging, and do not collapse with those of $\varSigma (10\varDelta _g)$, and $\varPsi (10\varDelta _g)$, while the positive tails are less affected. This asymmetry in the tails persists for large averaging volumes, as shown in figure 5(e, f), where we display the probability distribution of $\langle \varSigma \rangle _V$ and $\langle \varPsi \rangle _V$ for a volume of the order of the integral scale of the flow, $V^{1/3}=1.2L_0=32\varDelta _g$. Neither of the two probability distributions are Gaussian, specially in its negative tails.

These results point to essential differences between inverse and direct energy transfer events, which seem to emerge from different dynamical processes. While direct energy fluxes are correlated with fluxes at larger scales, inverse energy fluxes are not, suggesting that local backscatter, measured as $\varSigma <0$ or $\varPsi <0$, does not necessarily imply a net inverse cascade of energy. This is corroborated by the strong depletion of inverse energy transfer events with volume averaging, which suggest that backscatter might be a consequence of space-local conservative fluxes, which cancel out when averaged, rather than of interscale interactions. Let us note that $\varSigma$ and $\varPsi$ differ in fact by a spatial flux, and that the intensity and probability of backscatter events are a specially noticeable difference between the two definitions. These space-local fluxes seem to cause a higher probability of backscatter in $\varPsi$ than in $\varSigma$, suggesting that they are also, to some extent, the cause of backscatter in $\varSigma$. In view of these results, it is unlikely that the probability of inverse cascades can be quantified from the probability of local backscatter. Even if we could do so by analysing a massive database, and by devising a local definition of the energy transfer free of spatial fluxes, our results suggest that the energy cascade is sufficiently out of equilibrium for these predictions to have little practical importance.

6.1. Dynamics of the invariants of the velocity gradient tensor

Let $\mathcal {A}_{ij}=\partial _{j} u_i$ be the velocity gradient tensor of an incompressible flow. The second and third invariants of $\mathcal {A}_{ij}$ are

(6.1)\begin{gather} Q ={-}\tfrac{1}{2}\mathcal{A}_{ij}\mathcal{A}_{ji}=\tfrac{1}{4}\omega_i\omega_i-\tfrac{1}{2}S_{ij}S_{ij}, \end{gather}

(6.2)\begin{gather}R ={-}\tfrac{1}{3}\mathcal{A}_{ij}\mathcal{A}_{jk}\mathcal{A}_{ki}={-}\tfrac{1}{3}S_{ij}S_{jk}S_{ki} - \tfrac{1}{4}\omega_iS_{ij}\omega_j, \end{gather}

where $\omega _i$ is the $i$-th component of the vorticity vector. Considering $R$ and the discriminant $D=27/4R^2+Q^3$, the flow can be classified in four different topological types: positive $D$ corresponds to rotating topologies and negative $D$ to saddle-node topologies; negative $R$ accounts for topologies with a stretching principal direction and positive $R$ for topologies with a compressing principal direction. These invariants are also related to the vorticity vector and the rate-of-strain tensor, such that $Q$ indicates the balance between enstrophy and strain, denoted by $|S|^2=2S_{ij}S_{ij}$ and $|\omega |^2=\omega _i\omega _i$, and $R$ represents the balance between vortex stretching and strain self-amplification, $\omega _iS_{ij}\omega _j$ and $S_{ij}S_{jk}S_{ki}$. These terms appear in the evolution equations of $|S|^2$ and $|\omega |^2$,

(6.3)\begin{gather} \tfrac{1}{4}D_t |S|^2 ={-}S_{ij}S_{jk}S_{ki}-\tfrac{1}{4}\omega_iS_{ij}\omega_j-S_{ij}\partial_{ij}p, \end{gather}

(6.4)\begin{gather}\tfrac{1}{2}D_t |\omega|^2= \omega_iS_{ij}\omega_j, \end{gather}

where $D_t=\partial _t+u_j\partial _j$ is the substantial derivative along a Lagrangian trajectory. Note that for simplicity we have considered the evolution of an inviscid flow for simplicity.

The statistical distribution of $Q$ and $R$ in turbulent flows has a typical teardrop shape, which is shown in figure 6(e). For ease of reference, we have divided the $Q$–$R$ plane into four quadrants, where $q_1$ corresponds to $Q>0$ and $R>0$, $q_2$ to $Q>0$ and $R<0$, $q_3$ to $Q<0$ and $R<0$, and $q_4$ to $Q<0$ and $R>0$. The teardrop shape is characterised by a lobe in $q_2$, where the enstrophy is dominant over the strain and the vortex stretching over the strain self-amplification, and a tail in $q_4$, which is known as the Vieillefosse (Reference Vieillefosse1984) tail, and represents dominant strain self-amplification over vortex stretching in strain-dominated regions. The high absolute values of $Q$ in the Vieillefosse tail and in the upper semiplane indicate the spatial segregation of $|S|$ and $|\omega |$. Low absolute values of $Q$ do not in general imply low values of $|S|$ or $|\omega |$, but rather that $|S|\sim |\omega |$. We actually find that $|S|\sim \langle |S| \rangle$ in regions where $Q\sim 0$ and $R\sim 0$. For a detailed interpretation of the $Q$–$R$ plane, we refer the reader to Tsinober (Reference Tsinober2000).

Figure 6. Line arrow, solid line, conditional mean trajectories in the $Q$–$R$ plane for the inverse (a–d) and direct (e–h) evolutions due to different contributions. (i–l) Modulus of the probability transport velocities $|\boldsymbol {\varPhi }|$ in the $Q$–$R$ plane in the direct evolution due to the different contributions. (a,e,i) All contributions $\boldsymbol {\varPhi }_T$, (b, f, j) SGS model $\boldsymbol {\varPhi }_M$, (c,g,k) restricted Euler $\boldsymbol {\varPhi }_E$ and (d,h,l) non-local component of the pressure Hessian $\boldsymbol {\varPhi }_P$. Contours of the probability density function of $Q$ and $R$ contain 0.9 and 0.96 of the total data. All quantities are normalised using $t_Q={\langle Q^2 \rangle }^{1/4}$.

The dynamics of the velocity gradients can be statistically represented by the probability transport velocities in the $Q$–$R$ plane, which are obtained from the average Lagrangian evolution of $Q$ and $R$ (Ooi et al. Reference Ooi, Martín, Soria and Chong1999). Taking spatial derivatives in the evolution equations of the velocity field described by (3.6) and considering (6.1), (6.2), the Lagrangian evolution of the invariants reads as

(6.5)\begin{gather} {D_tQ}={-}3R+\mathcal{A}_{ij}\mathcal{H}_{ji}, \end{gather}

(6.6)\begin{gather}{D_tR}=\tfrac{2}{3}Q^2-\mathcal{A}_{ij}\mathcal{A}_{jk}\mathcal{H}_{ki}, \end{gather}

where $\mathcal {H}_{ij}=\mathcal {H}^P_{ij}+\mathcal {H}^M_{ij}$ comprises the contribution of the non-local component of the pressure Hessian and the SGS model, $\mathcal {H}_{ij}^P={\partial _{ij} p -((\partial _{kk}p)\delta _{ij})/{3}}$ and $\mathcal {H}_{ij}^M=\partial _j{M_i}$. Here, $M_i$ is the $i$-th component of the SGS model in the momentum equation and $\delta _{ij}$ is the Kronecker's delta. The total probability transport velocity is defined as

(6.7)\begin{equation} \boldsymbol{\varPhi}_T=\{\langle {D_tQ}\rangle_C t_Q^3, \langle {D_tR}\rangle_C t_Q^4\}, \end{equation}

where $\langle \cdot \rangle _C$ is the probability conditioned to $Q$ and $R$, and $t_Q=1/\langle Q^2 \rangle ^{1/4}$ is a characteristic time extracted from the standard deviation of $Q$. Let us note that the probability density flux in the $Q$–$R$ plane is $P(Q,R)\boldsymbol {\varPhi }_T$, where $P(Q,R)$ is the joint probability density of $Q$ and $R$. The probability transport velocities are integrated to yield conditional mean trajectories (CMTs) in the $Q$–$R$ plane (Ooi et al. Reference Ooi, Martín, Soria and Chong1999), which are shown also in figure 6(e).

To analyse the dynamics of the invariants, we decompose the CMTs into the contribution of the different terms in the evolution equation of $Q$ and $R$,

(6.8)\begin{gather} \boldsymbol{\varPhi}_E=\{{-}3R, 2/3Q^2 \}, \end{gather}

(6.9)\begin{gather}\boldsymbol{\varPhi}_P=\{\langle \mathcal{A}_{ij}\mathcal{H}_{ji}^P \rangle_{C}, \langle -\mathcal{A}_{ij}\mathcal{A}_{jk}\mathcal{H}_{ki}^P\rangle_{C}\}, \end{gather}

(6.10)\begin{gather}\boldsymbol{\varPhi}_M=\{\langle \mathcal{A}_{ij}\mathcal{H}_{ji}^M \rangle_{C}, \langle -\mathcal{A}_{ij}\mathcal{A}_{jk}\mathcal{H}_{ki}^M \rangle_{C}\}, \end{gather}

where $\boldsymbol {\varPhi }_T=\boldsymbol {\varPhi }_E+\boldsymbol {\varPhi }_P+\boldsymbol {\varPhi }_M$, and the normalisation with $t_Q$ has been dropped for simplicity. This decomposition separates the strictly space-local dynamics of restricted Euler (RE) (Cantwell Reference Cantwell1992), $\boldsymbol {\varPhi }_E$, which depends exclusively on $Q$ and $R$, and includes advection and the local action of the pressure Hessian, from the non-local action of the pressure Hessian and the SGS model, which are included in $\boldsymbol {\varPhi }_P$ and $\boldsymbol {\varPhi }_M$, respectively. This separation into local and non-local dynamics is convenient to justify the prevalence of direct over inverse cascades. In order to compare the importance of the different terms, we consider the norm of the probability transport velocities,

(6.11)\begin{equation} |\boldsymbol{\varPhi}|=\sqrt{(t_Q^3\langle {D_tQ}\rangle_C)^2 + (t_Q^4\langle{D_tR}\rangle_C)^2}. \end{equation}

Statistics of $Q$ and $R$ and their CMTs have been compiled for two different times in the ensemble of realisations, $t_{stats}=(1 \pm 0.12)t_{inv}$, which correspond to the direct and inverse evolutions and are marked in figure 1(a). If we assume that the decay of the system is self-similar, conclusions drawn from this analysis should be independent of $t_{stats}$. The total probability transport velocity, $\boldsymbol {\varPhi }_T$, has been calculated directly from the expression of the substantial derivative, computing separately the temporal and convective derivatives. Appropriate numerical methods have been used in all the computations (Lozano-Durán, Holzner & Jiménez Reference Lozano-Durán, Holzner and Jiménez2015). Results show that the number of realisations in our database is sufficient to converge the statistics.

The different contributions to the CMTs and the norm of the probability transport velocities in the inverse and direct evolutions are shown in figures 6(a–h) and 6(i–l). In the direct evolution, the CMTs develop an average rotating cycle characteristic of turbulent flows, which exposes causality between the different configurations of the flow in an average sense. The CMTs rotate clockwise: from dominant vortex stretching to dominant vortex compression in the upper semiplane, to the Vieillefosse tail in $q_4$, and again to dominant vortex stretching in enstrophy-dominated regions.

The CMTs move from negative to positive $R$, and finally to the Vieillefosse tail, due to the effect of the RE dynamics. This trend would cause the appearance of infinite gradients in a finite time (Vieillefosse Reference Vieillefosse1984), but the non-local effect of the pressure Hessian and the SGS model, or viscosity in the case of direct numerical simulations, prevent it by bringing the CMTs back to $q_3$ and restarting the cycle. For closed steady states, CMTs describe a closed cycle (Lozano-Durán et al. Reference Lozano-Durán, Holzner and Jiménez2016), while in our decaying flow they spiral inwards due to the action the SGS model, which contracts the probability distribution, resembling the action of the viscosity in direct numerical simulations (Meneveau Reference Meneveau2011). The CMTs of the non-local component of the pressure Hessian evolve from positive to negative values of $R$, following the same behaviour observed by Chevillard et al. (Reference Chevillard, Meneveau, Biferale and Toschi2008) and Meneveau (Reference Meneveau2011). The norm of the probability transport velocities reveals that the RE dynamics is mostly counteracted by the non-local component of the pressure Hessian in regions where $Q$ and $R$ are large. These observations are in agreement with Lüthi, Holzner & Tsinober (Reference Lüthi, Holzner and Tsinober2009), and evidence a secondary role of the model in the dynamics of intense gradients.

In the inverse evolution we identify substantial changes. As explained in § 4, time reversal changes only the sign of quantities which are odd with the velocity: $Q$ remains unaltered while $R$ changes sign, leading to an inverse teardrop shape. The Vieillefosse tail lies now in $q_3$, forming an antitail, and a higher probability of vortex compression appears in enstrophy-dominated regions. This transformation also affects the CMTs, and is most relevant in the effect of the SGS model, which now expands the probability distribution and leads to an average outward spiralling. The behaviour of $\boldsymbol {\varPhi }_E$ and $\boldsymbol {\varPhi }_P$ is similar to the direct evolution in the upper semiplane, but undergoes fundamental changes in the lower semiplane. In the direct evolution the non-local component of the pressure Hessian counteracts the RE dynamics, preventing the formation of intense gradients in the Vieillefosse tail. Conversely, in the inverse evolution, $\boldsymbol {\varPhi }_P$ favours the growth of intense gradients in the antitail, while the RE opposes it.

This analysis compares the structure of the ‘attractor’ with that of the ‘antiattractor’, and can be simply derived by considering the effect that the transformation $\boldsymbol {u}\rightarrow - \boldsymbol {u}$ has on the dynamics of the velocity gradients. However, we have identified in § 5.1 inverse trajectories outside the antiattractor. We will show that, for these trajectories, this analysis yields non-trivial results.

6.2. Asymmetry in the $Q$–$R$ space

The statistical irreversibility of the cascade is clearly manifest in the statistics of the invariants through the asymmetry of the probability distribution with respect to $R=0$. To characterise statistical irreversibility in the $Q$–$R$ plane, we use an asymmetry function

(6.12)\begin{equation} \varpi(Q,R)=\log \frac{P^+(Q,R)}{P^-(Q,R)}, \end{equation}

where $P^+$ and $P^-$ denote the probability density of $Q$ and $R$ in the direct and inverse evolutions, respectively. Figure 7(a) shows the distribution of $\varpi (Q,R)$, while 7(b) shows the absolute values of the asymmetry function averaged along the $R$ axis,

(6.13)\begin{equation} \langle\varpi(Q)\rangle_R=\frac{\displaystyle\int|\varpi(Q,R)|{\,\textrm{{d}}}{R}}{\displaystyle\int {\,\textrm{{d}}} R}. \end{equation}

Most of the temporal asymmetry of the velocity gradients is related to $Q<0$, where the strain is dominant over the vorticity. In enstrophy-dominated regions vortex stretching is on average dominant in the direct evolution, but the probability of vortex compression is not negligible. Both processes also occur in the inverse cascade, in which vortex compression is more probable but coexists with vortex stretching. On the other hand, the structure of strain-dominated regions depends strongly on the direction of the system in time. Intense velocity gradients in $Q<0$ are either in $q_3$ or $q_4$, depending on the direction of the cascade, but not simultaneously in both quadrants. The dynamics of strain-dominated regions are more fundamentally related to the direction of the system in time than the dynamics of enstrophy-dominated regions.

Figure 7. (a) Colour plot of the asymmetry function $\varpi (Q,R)$. Isocontours contain $0.7$ and $0.85$ of the total data. Dashed lines represent the contours of the asymmetry function when $R\rightarrow -R$. (b) Asymmetry function averaged along the $R$ axis, $\langle \varpi (Q)\rangle _R$, for $R$ in the interval $(-2,2)$. Quantities are normalised using $t_Q={\langle Q^2 \rangle }^{1/4}$.

6.3. Energy fluxes in the $Q$–$R$ plane: inverse evolutions outside the ‘antiattractor’

We now analyse the relation between the local structure of the flow and the energy cascade by conditioning the statistics of the local energy fluxes to the invariants of the filtered velocity gradients. First, we study the experiments on the inverse cascade without the model, which are presented in § 4.3. These experiments provide non-trivial inverse evolutions outside the ‘antiattractor’, for which this analysis identifies the structures that support the inverse energy cascade.

We calculate the invariants of the filtered velocity gradients and their CMTs at scale $\check {\varDelta }=5\varDelta _g$. We use a Gaussian filter (4.1), and calculate the CMTs using the substantial derivative of the filtered field, $\check {D}_t=\partial _t+\check {u}_i\partial _i$. We analyse these quantities in the inverse evolutions when the model is removed, and connect them to the evolution of the local energy fluxes in physical space, $\varSigma (\boldsymbol {x},t)$, and in scale space, $\varPi (k,t)$. The analysis of $\varPsi$ yields qualitatively similar results.

Figure 8(a–h) show the probability distribution of the invariants of the unfiltered (a–d), and filtered (e–h) velocity gradients, and the their CMTs, in the inverse evolution without the model. Figure 8(i–l) show the averaged $\varSigma (5\varDelta _g)$ conditioned to the invariants of the filtered velocity gradients.

Figure 8. Conditional mean trajectories in the $Q$–$R$ plane of the inverse evolution without the model (see § 4.3) for (a–d) the unfiltered case and (e–h) filtered at $\check {\varDelta }=5\varDelta _g$, and (i–l) conditional average of $\varSigma$ at scale $\check {\varDelta }=5\varDelta _g$, at different times: (a,e) $t/T_0=0$; (b, f) $t/T_0=0.1$; (c,g) $t/T_0=0.15$; (d,h) $t/T_0=0.25$. ——, solid magenta line, contours of the probability density function of $Q$ and $R$ containing $0.8$ and $0.9$ of the data. Here $Q$ and $R$ are normalised with $t_Q=\langle Q(t_{inv})^2\rangle ^{1/4}$. Conditional-averaged energy fluxes in (i–l) are normalised with the absolute value of $\langle \varSigma \rangle$ at each time.

Initially, the statistics of the filtered and unfiltered velocity gradients are similar, and resemble those of the inverse evolutions. The average $\varSigma$ conditioned to $Q$ and $R$ is predominantly negative and most intense in $q_3$ and $q_1$, suggesting a relevant role of the antitail and of vortex compression in the inverse cascade.

At time $0.1T_0$ after the inversion, the antitail of the unfiltered gradients has contracted, while the antitail of the filtered gradients remains unaltered. We do not yet see changes in the upper semiplane, where the probability distribution remains similar to the initial state in both cases. According to $\varPi (k,t)$ in figure 2(a), at this time the direct cascade has not yet regenerated at the small scales. At the filter scale, $\varPi (k,t)$ and $\langle \varSigma \rangle$ are still negative, and have a similar magnitude to the initial state.

At $0.15T_0$, the unfiltered velocity gradients start developing dominant vortex stretching, represented by a prominent lobe in $q_2$, and a regular Vieillefosse tail. The antitail, although reduced, is still observable in the filtered gradients and $\langle \varSigma \rangle <0$. At this stage, as shown in figure 2(a), an inverse cascade at the filter scale, and a direct cascade in the small scales coexist. Finally, at $0.25T_0$, when the filtered velocity gradients develop structures in the Vieillefosse tail, the direct cascade regenerates at the filter scale and $\langle \varSigma \rangle >0$.

These results indicate a strong connection between the direction of the cascade and the orientation of the Vieillefosse tail. First, the conditional averages of $\varSigma$ in the $Q$–$R$ plane show that, at $0.25T_0$, when the energy cascade starts to regenerate at the filter scale, direct local energy fluxes are most intense in the $q_4$ quadrant. Second, the differences in the invariants of the filtered velocity gradients between $0.1T_0$ and $0.15T_0$ are only significant in the Vieillefosse tail, but the inverse energy fluxes are reduced by a half, $\langle \varSigma (0.15T_0)\rangle \simeq 0.5\langle \varSigma (0.1T_0)\rangle$.

The symmetry in the statistics of enstrophy-dominated regions reported in § 6.2 also holds during the transition from inverse to direct phase-space trajectories. At the filter scale, the probability distribution or the CMTs in the upper semiplane do not substantially change during the regeneration of the energy cascade.

6.4. Energy fluxes in the $Q$–$R$ plane: direct evolutions in the turbulent ‘attractor’

We perform the same analysis for the direct evolutions with the SGS model. Figure 9(a) shows the average $\varSigma$ conditioned to $Q$ and $R$. To discriminate between direct and inverse energy transfer, we condition the probability distribution of $Q$ and $R$ to the sign of the local energy fluxes, and present the results in figure 9(b,c). We quantify in table 3 the relevance of each quadrant to the energy cascade by calculating the average of $\varSigma$ conditioned to $Q$ and $R$ belonging to a quadrant, $\langle \varSigma \rangle _{q_i}=\langle \varSigma |(Q,R)\in q_i\rangle$, and the total contribution of each quadrant to the energy fluxes, $\{\varSigma \}_{q_i}=\langle \varSigma \rangle _{q_i}P_{q_i}$, where $P_{q_i}$ is the probability that $Q$ and $R$ belong to $q_i$, so that $\langle \varSigma \rangle =\sum _{i=1,4}\{\varSigma \}_{q_i}$. We have also calculated these quantities conditioning the averages to the direction of the energy fluxes, being

(6.14)\begin{equation} \langle \varSigma \rangle^+_{q_i}=\langle \varSigma |\varSigma>0,(Q,R)\in q_i\rangle, \end{equation}

the average energy transfer conditioned to the quadrant $q_i$ and to $\varSigma >0$. Accordingly, $\{\varSigma \}^+_{q_i}=\langle \varSigma \rangle ^+_{q_i}P^+_{q_i}$ is the total contribution of the quadrant $q_i$ to the positive energy fluxes. Here $P^+_{q_i}$ is the probability that $Q$ and $R$ belongs to $q_i$ and that $\varSigma >0$. Statistics conditioned to negative fluxes are marked with a minus superscript.

Figure 9. (a) Conditional average of $\varSigma$ at scale $\check {\varDelta }=5\varDelta _g$ in the $Q$–$R$ plane of the filtered velocity gradients for the direct evolution with the model. (b) Average $\varSigma$ over $Q$ and $R$ conditioned to $\varSigma >0$. (c) Absolute value of the average $\varSigma$ over $Q$ and $R$ conditioned to $\varSigma <0$. Values of $\varSigma$ normalised with $\langle \varSigma \rangle$ at $t_{stats}$. ——, solid magenta line, contours of the probability density function of $Q$ and $R$ containing $0.85$ and $0.9$ of the data. In (c), $-\cdot -$ marks isocontours of $Q$ and $R$ conditioned to negative $\varSigma$ for the inverse evolutions. The dashed line corresponds to $D=0$. Here $Q$ and $R$ are normalised with $t_Q={\langle Q^2 \rangle }^{1/4}$.

Table 3. Conditional averages of energy transfer and contributions to the total energy fluxes of each quadrant in the $Q$–$R$ plane during the direct evolution. In the first table we have considered all data, in the second only direct energy transfer events, and in the last only backscatter events. Quantities are normalised with the average $\varSigma$ in the first case, and with the averages conditioned to positive and negative fluxes in the other two cases.

The most intense direct energy transfer events are related to structures in the Vieillefosse tail, and $q_4$ has the highest average energy transfer, $\langle \varSigma \rangle _{q_4}\sim 1.5\langle \varSigma \rangle$, and contributes with approximately $50\,\%$ of the total energy transfer. Also structures in $q_2$ are relevant to the cascade, with a contribution of around $30\,\%$ to the total fluxes. Conversely, the contribution of $q_1$ to direct energy transfer is negligible. When we condition these statistics to direct energy fluxes, we obtain qualitatively similar results. The statistics of $Q$ and $R$ conditioned to positive $\varSigma$ yield almost a complete teardrop, except for a lower probability of events in $q_1$.

When we condition to negative energy transfer events, the average inverse energy transfer is roughly similar in the four quadrants, but the $q_1$ quadrant contributes the most to the total, with a $56\,\%$. This is a consequence of inverse energy transfer events being mostly located in $q_1$. As shown in figure 9(c), the statistics of $Q$–$R$ conditioned to inverse energy transfer resemble an inverse teardrop without the antitail. In the same figure we have plotted the probability density function of $Q$ and $R$ conditioned to inverse energy transfer in the inverse evolutions. By time symmetry, the antitail in $q_3$ contributes the most to inverse energy transfer in the inverse cascade, but this tail is absent in the direct evolution.

6.5. An entropic argument for the prevalence of direct energy fluxes

We have presented in the previous section compelling evidence that topologies in the antitail are strongly connected to the inverse cascade, but that these topologies have negligible probability in the direct cascade. Here we provide arguments that explain the low probability of these topologies, and, consequently, of inverse cascades.

Consider the evolution of $R$ in regions where the strain is dominant over the enstrophy, $|S|^2\sim \langle |S|^2\rangle \gg |\omega |^2$. Neglecting the contribution of the enstrophy, we obtain

(6.15)\begin{equation} R\approx{-}\tfrac{1}{3}S_{ij}S_{jk}S_{ki}={-}\alpha_1\alpha_2\alpha_3, \end{equation}

where $\alpha _1\leqslant \alpha _2\leqslant \alpha _3$ are the eigenvalues of the rate-of-strain tensor. This is a good approximation to the Vieillefosse tail and the antitail, where $|S|^2\sim 100|\omega |^2$. Due to compressibility $\alpha _1 + \alpha _2 + \alpha _3=0$, and the sign of $R$ is determined by the sign of the intermediate eigenvalue. Topologies in the Vieillefosse tail are characterised by $\alpha _2>0$, while in the antitail $\alpha _2<0$. We describe the evolution of $R$ in terms of the evolution of $\alpha _2$, whose equation reads as

(6.16)\begin{equation} {D_t\alpha_2}={-}\alpha_2^2-\left(\frac{\varLambda}{3}+\varGamma_2\right), \end{equation}

where $\varLambda =\partial _{kk}p\approx -(|S|^2)/{2}=-\alpha ^2_1-\alpha ^2_2-\alpha ^2_3$ is the trace of the pressure Hessian, which is local, and $\varGamma _2$ is the contribution of $\mathcal {H}^P_{ij}=\partial _{ij} p-(\partial _{kk}p\delta _{ij})/{3}$ to the evolution of $\alpha _2$, which is non-local. We have discarded quadratic terms in the vorticity and the action of the SGS model, which we show above to have a negligible contribution to the dynamics of intense gradients (Nomura & Post Reference Nomura and Post1998). Disregarding the action of the non-local component of the pressure Hessian, $\varGamma _2$, the equation reads as

(6.17)\begin{equation} D_t\alpha_2=\tfrac{1}{3}(\alpha_1^2+\alpha_3^2-2\alpha_2^2)>0, \end{equation}

where we have used the fact that, from their definition, $\alpha _1^2>\alpha _2^2$ and $\alpha _3^2>\alpha _2^2$. The evolution of any initial value of $\alpha _2$ results in the growth of the intermediate eigenvalue, leading to positive values of $R$. This is the effect of the RE dynamics, which depletes the antitail in the direct cascade. Left to itself, it generates intense strain in the Vieillefosse tail, and is independent on the non-local structure of the velocity gradients. From (6.17) we conclude that the only inertial mechanism capable of preventing the rate-of-strain tensor from developing a positive intermediate eigenvalue is the non-local action of the pressure Hessian through $\varGamma _2$. This analysis is consistent with the results presented in § 6.1, where we show that this term sustains the antitail in the inverse cascade.

The non-local pressure Hessian depends on the complete flow field. In the absence of boundaries it can be expressed as a singular integral, which depends on all the points of the domain (Speziale, Sarkar & Gatski Reference Speziale, Sarkar and Gatski1991; Ohkitani & Kishiba Reference Ohkitani and Kishiba1995). As a consequence, the non-local pressure Hessian is in general decorrelated from the local dynamics of the velocity gradients (She et al. Reference She, Jackson, Orszag, Hunt, Phillips and Williams1991). Figure 10(a) shows the joint probability distribution of the two components of $\boldsymbol {\varPhi }_P=\{\mathcal {A}_{ij}\mathcal {H}^P_{ji},-\mathcal {A}_{ij}\mathcal {A}_{jk}\mathcal {H}^P_{ki}\}$ in the centre of the $Q$–$R$ plane, and in the Vieillefosse tail during the inverse and direct evolutions. See figure 10(b) for reference. In regions where gradients are not strong, the non-local component of the pressure Hessian appears decorrelated from the local dynamics, with a large scatter of the data with respect to the mean. In the Vieillefosse tail the non-local component of the pressure Hessian counteracts the RE dynamics and the action of the non-local component of the pressure Hessian over $Q$ and $R$ appears more correlated. When time is reversed, these correlations are necessary to sustain the antitail. As shown in § 6.3, the contraction of the antitail is the first identifiable process in the transition from the inverse to the direct cascade, which indicates that these correlations are quickly destroyed under perturbations.

Figure 10. (a) Joint probability density function of the components of the probability transport velocities due to the non-local part of the pressure Hessian, $\{\mathcal {A}_{ij}\mathcal {H}^P_{ji},-\mathcal {A}_{ij}\mathcal {A}_{jk}\mathcal {H}^P_{ki}\}$, conditioned to different areas in the $Q$–$R$ plane: $-\cdot -$, $|Q|<1$ and $|R|<1$, dashed region in (b); – – –, dashed blue line, $Q<-1$ in the direct evolution, shaded region in (b); ——, solid red line, $Q<-1$ in the inverse evolution, shaded region in (b). Quantities normalised with $t_Q$.

While the generation of topologies in the Vieillefosse tail is a direct consequence of the local self-amplification of the rate-of-strain tensor, the generation of topologies in the antitail requires a global configuration of the complete flow field that counteracts this amplification through the non-local action of the pressure Hessian. This scenario appears intuitively unlikely, suggesting that only specially organised flows are able to produce such action (Popper Reference Popper1956). This is the case of the initial conditions used for the inverse evolutions, which conserve all the information of the direct cascade process. In particular, they retain the correlations between the rate-of-strain tensor and the non-local contribution of the pressure Hessian in the Vieillefosse tail, thus inhibiting the breakdown of the antitail in the inverse cascade.