Energy fluxes in quasi-equilibrium flows

TL;DR

This paper investigates the relationship between equilibrium states and turbulent flows in spectrally truncated Euler and Navier-Stokes equations, revealing conditions under which flows resemble equilibrium or turbulence, with implications for energy flux behavior.

Contribution

It provides a theoretical and numerical analysis of how truncated fluid systems transition between equilibrium and turbulent states, highlighting the role of forcing and dissipation scales.

Findings

Small $k_{max}\eta$ flows approach Kraichnan's equilibrium with an effective temperature.

Numerical simulations confirm the model predictions for small $k_{max}\eta$.

A transition from quasi-equilibrium to Kolmogorov turbulence occurs at intermediate $k_{max}\eta$.

Abstract

We examine the relation between the absolute equilibrium state of the spectrally truncated Euler equations (TEE) predicted by Kraichnan (1973) to the forced and dissipated flows of the spectrally truncated Navier-Stokes (TNS) equations. In both of these idealized systems a finite number of Fourier modes is kept contained inside a sphere of radius $k_{max}$ but while the first conserves energy in the second energy is injected by a body-force $\bf{f}$ and dissipated by the viscosity $\nu$. For the TNS system stochastically forced with energy injection rate $\mathcal{I}_\mathcal{E}$ we show, using an asymptotic expansion of the Fokker-Planck equation, that in the limit of small $k_{max}\eta$ (where $\eta=(\nu^3/\mathcal{I}_\mathcal{E})^{1/4}$ the Kolmogorov lengthscale) the flow approaches the absolute equilibrium solution of Kraichnan with such an effective "temperature" so that there is a balance between the energy injection and the energy dissipation rate. Numerical simulations verify the predictions of the model for small values of $k_{max}\eta$. For intermediate values of $k_{max}\eta$ a transition from the quasi-equilibrium "thermal" state to Kolmogorov turbulence is observed. If the forcing is applied at small scales and the dissipation acts only at large scales then the equipartition spectrum appears at all scales for all values of $\nu$. In both cases a finite forward or inverse flux is present even for the cases where the flow is close to the equilibrium state solutions. However, unlike the classical turbulence where an energy cascade develops with a mean energy flux that is large compared to its fluctuations, the quasi-equilibrium state has a mean flux of energy that is subdominant to the large flux fluctuations observed.

Equations

$${\partial }_t\mathbf{u}+{\mathbb{P}}_K[\mathbf{u}\cdot \nabla \mathbf{u}+\nabla p]=0,\nabla \cdot \mathbf{u}=0.$$

$$\mathrm{the}\mathrm{Energy}\mathcal{E}=\frac{1}{2}\int |\mathbf{u}{|}^{2}d{x}^3\mathrm{and}\mathrm{the}\mathrm{Helicity}\mathcal{H}=\frac{1}{2}\int \mathbf{u}\cdot \nabla \times \mathbf{u}d{x}^3.$$

$${{\color[rgb]{0,0,0}E(k)=\frac{1}{2k_{0}}\sum_{k\leq|{\bf k}|<k+k_{0}}|\tilde{\bf u}_{\bf k}|^{2}\quad\mathrm{and}\quad H(k)=\frac{1}{2k_{0}}\sum_{k\leq|{\bf k}|<k+k_{0}}\tilde{\bf u}_{\bf-k}\cdot(i{\bf k}\times\tilde{\bf u}_{\bf k}).}}$$

$$\mathcal{P}(\mathbf{u})={\mathcal{Z}}^{-1}\exp \left[-\alpha \mathcal{E}-\beta \mathcal{H}\right]$$

$$E(k)=\frac{4\pi \alpha k^2}{{\alpha }^2-{\beta }^2{k}^2},H(k)=\frac{4\pi \beta k^4}{{\alpha }^2-{\beta }^2{k}^2}.$$

$${{\color[rgb]{0,0,0}\mathcal{E}=k_{0}\sum_{\bf k}E(k)\quad\mathrm{and}\quad\mathcal{H}=k_{0}\sum_{\bf k}H(k)}}$$

$${\partial }_t\mathbf{u}+{\mathbb{P}}_K[\mathbf{u}\cdot \nabla \mathbf{u}+\nabla P]=\mathbf{f}+\nu \Delta \mathbf{u}.$$

$${\mathbf{h}}_{\mathbf{k}}^s=\frac{\mathbf{k}\times (\hat{\mathbf{e}}\times \mathbf{k})}{\sqrt{2}|k\times (\hat{\mathbf{e}}\times \mathbf{k})|}+is\frac{\hat{\mathbf{e}}\times \mathbf{k}}{\sqrt{2}|\hat{\mathbf{e}}\times k|}$$

$${\partial }_t{\tilde{u}}_{\mathbf{k}}^s={\mathcal{V}}_{\mathbf{k}}^s-\nu k^2{\tilde{u}}_{\mathbf{k}}^s+{\tilde{f}}_{\mathbf{k}}^s$$

$${\mathcal{V}}_{\mathbf{k}}^s=\sum _{\mathbf{p}+\mathbf{q}=\mathbf{k}}\sum _{s_p,s_q}{C}_{\mathbf{k},\mathbf{q},\mathbf{p}}^{s,s_q,s_p}{\tilde{u}}_{\mathbf{p}}^{s_p}{\tilde{u}}_{\mathbf{q}}^{s_q}$$

$$\sum _{s,\mathbf{k}}{\tilde{u}}_{-\mathbf{k}}^s{\mathcal{V}}_{\mathbf{k}}^s=0,\sum _{s,\mathbf{k}}sk{\tilde{u}}_{-\mathbf{k}}^s{\mathcal{V}}_{\mathbf{k}}^s=0$$

$$\sum _{s,\mathbf{k}}\frac{\partial }{\partial {\tilde{u}}_{\mathbf{k}}^s}{\mathcal{V}}_{\mathbf{k}}^s=0.$$

$$\frac{\partial }{\partial t}\mathcal{P}+\sum _{s,\mathbf{k}}\frac{\partial }{\partial {\tilde{u}}_{\mathbf{k}}^s}({\mathcal{V}}_{\mathbf{k}}^s\mathcal{P})=\nu \sum _{s,\mathbf{k}}\frac{\partial }{\partial {\tilde{u}}_{\mathbf{k}}^s}(k^2{\tilde{u}}_{\mathbf{k}}^s\mathcal{P})+\sum _{s,\mathbf{k}}{ϵ}_{\mathbf{k}}^s\frac{\partial }{\partial {\tilde{u}}_{\mathbf{k}}^s}\frac{\partial }{\partial {\tilde{u}}_{-\mathbf{k}}^s}\mathcal{P},$$

$$\frac{\partial }{\partial t}⟨\mathcal{E}⟩=-⟨{\mathcal{D}}_{\mathcal{E}}⟩+⟨{\mathcal{I}}_{\mathcal{E}}⟩$$

$$\sum _{s,\mathbf{k}}\frac{\partial }{\partial {\tilde{u}}_{\mathbf{k}}^s}({\mathcal{V}}_{\mathbf{k}}^s{\mathcal{P}}_0)=0.$$

$${\mathcal{P}}_0=f(\mathcal{E})=f\left(\frac{1}{2}\sum _{s,\mathbf{k}}|{\tilde{u}}_{\mathbf{k}}^s{|}^2\right).$$

$$\sum _{s,\mathbf{k}}\frac{\partial }{\partial {\tilde{u}}_{\mathbf{k}}^s}({\mathcal{V}}_{\mathbf{k}}^s{\mathcal{P}}_1)={\nu }^{\prime }\sum _{s,\mathbf{k}}\frac{\partial }{\partial {\tilde{u}}_{\mathbf{k}}^s}(k^2{\tilde{u}}_{\mathbf{k}}^s{\mathcal{P}}_0)+\sum _{s,\mathbf{k}}{ϵ}_{\mathbf{k}}^{\prime s}\frac{\partial }{\partial {\tilde{u}}_{\mathbf{k}}^s}\frac{\partial }{\partial {\tilde{u}}_{-\mathbf{k}}^s}{\mathcal{P}}_0$$

$$\frac{\partial }{\partial {\tilde{u}}_{\mathbf{k}}^s}f(\mathcal{E})=\frac{\partial \mathcal{E}}{\partial {\tilde{u}}_{\mathbf{k}}^s}\frac{\partial f}{\partial \mathcal{E}}=\left(\frac{\partial }{\partial {\tilde{u}}_{\mathbf{k}}^s}\frac{1}{2}\sum _{s_s,\mathbf{q}}{\tilde{u}}_{-\mathbf{q}}^{s_q}{\tilde{u}}_{\mathbf{q}}^{s_q}\right)\frac{\partial f}{\partial \mathcal{E}}={\tilde{u}}_{-\mathbf{k}}^s\frac{\partial f}{\partial \mathcal{E}}$$

$$\sum _{s,\mathbf{k}}\frac{\partial }{\partial {\tilde{u}}_{\mathbf{k}}^s}({\mathcal{V}}_{\mathbf{k}}^s{\mathcal{P}}_1)={\nu }^{\prime }\sum _{s,\mathbf{k}}\left(f+|{\tilde{u}}_{\mathbf{k}}^s{|}^2\frac{\partial f}{\partial \mathcal{E}}\right)k^2+\left(\left[\sum _{s,\mathbf{k}}{ϵ}_{\mathbf{k}}^{\prime s}\right]\frac{\partial f}{\partial \mathcal{E}}+\left[\sum _{\mathbf{k}}{ϵ}_{\mathbf{k}}^{\prime s}|{\tilde{u}}_{\mathbf{k}}^s{|}^2\right]\frac{{\partial }^{2}f}{\partial {\mathcal{E}}^2}\right).$$

$${\nu }^{\prime }\left(\sum _{\mathbf{k}}|\mathbf{k}{|}^2\right)\left(f+\frac{\mathcal{E}}{N}\frac{\partial f}{\partial \mathcal{E}}\right)+\left[\sum _{s,\mathbf{k}}{ϵ}_{\mathbf{k}}^{\prime s}\right]\left(\frac{\partial f}{\partial \mathcal{E}}+\frac{\mathcal{E}}{N}\frac{{\partial }^{2}f}{\partial {\mathcal{E}}^2}\right)=0.$$

$$\frac{\partial }{\partial \mathcal{E}}\left({\nu }^{\prime }{K}^2{\mathcal{E}}^{N}f+{\mathcal{I}}_{\mathcal{E}}^{\prime }{\mathcal{E}}^N\frac{\partial f}{\partial \mathcal{E}}\right)=0$$

$$f(\mathcal{E})={\mathcal{Z}}^{-1}\exp \left(-\frac{\nu K^2}{{\mathcal{I}}_{\mathcal{E}}}\mathcal{E}\right)$$

$$\mathcal{Z}=S_{2N}{2}^{N-1}{\left(\frac{{\mathcal{I}}_{\mathcal{E}}}{\nu K^2}\right)}^N\Gamma (N)$$

$$P(\mathcal{E})=\frac{S_{2N}(2\mathcal{E}{)}^{N-1}}{\mathcal{Z}}\exp \left(-\frac{\nu K^2}{{\mathcal{I}}_{\mathcal{E}}}\mathcal{E}\right).$$

$$⟨\mathcal{E}⟩=\frac{N{\mathcal{I}}_{\mathcal{E}}}{\nu K^2}.$$

$$E(k)=\frac{4\pi {\mathcal{I}}_{\mathcal{E}}}{\nu K^2{k}_0^3}{k}^2$$

$$⟨{\mathcal{D}}_{\mathcal{E}}⟩$$

$$k_{tran}=k_{max}{\left(\frac{2C_K}{5}\right)}^{3/11}(k_{max}\eta {)}^{4/11}$$

$$\underset{\nu \to 0}{\lim }\frac{⟨{\mathcal{I}}_{\mathcal{E}}⟩}{⟨2\mathcal{E}{⟩}^{2/3}{k}_0}=0$$

Full text

1 Introduction

Turbulence is a classical example of an out of equilibrium system. In steady state, energy is constantly injected at some scale $\ell_{in}$ while it is dissipated at smaller scales $\ell_{\nu}$ by viscous forces. This process requires a finite flux of energy from the former scale $\ell_{in}$ to the latter $\ell_{\nu}$ that is provided by the well known Kolmogorov-Richardson cascade. Despite the out of equilibrium nature of turbulence there are circumstances where equilibrium dynamics become relevant. This has been claimed to be the case for the scales larger than than the injection scale $\ell>\ell_{in}$. At these scales the energy flux is zero and can possibly be modeled using equilibrium dynamics (Dallas et al., 2015; Cameron et al., 2017; Alexakis & Brachet, 2018). Furthermore, at the smallest scales of the inertial range a so-called ‘bottleneck’ manifests where the power-law slope of the energy spectrum becomes less steep (Donzis & Sreenivasan, 2010; Falkovich, 1994; Martinez et al., 1997; Lohse & Müller-Groeling, 1995). This has been interpreted by Frisch et al. (2008) as an ‘incomplete thermalization’ that becomes asymptotically at equilibrium for hyper-viscous flows when the order of the hyper-viscosity tends to infinity.

Equilibrium dynamics become also relevant in the presence of inverse cascades in finite domains where large scale condensates form (Kraichnan, 1967; Robert & Sommeria, 1991; Naso et al., 2010; Bouchet & Venaille, 2012; Shukla et al., 2016). Finally understanding equilibrium dynamics is important for systems that display a transition of from forward to an inverse cascade (Alexakis & Biferale, 2018; Sahoo et al., 2017; Benavides & Alexakis, 2017; Sozza et al., 2015; Deusebio et al., 2014; Seshasayanan et al., 2014) since the large scale flows in these systems transition from an equilibrium state to an out of equilibrium state. Besides the possible applications, understanding the equilibrium dynamics in turbulence is also a much needed step required before understanding its much harder out-of-equilibrium counterpart. This has led many researchers (Lee, 1952; Hopf, 1952; Kraichnan, 1967, 1973; Orszag, 1977) to investigate the equilibrium state of the truncated Euler equations (TEE) where only a finite number of Fourier modes is kept and are given by:

[TABLE]

Here $\bf u$ is the incompressible velocity field, $p$ is the pressure and $\mathbb{P}_{K}$ is a projection operator that sets to zero all Fourier modes except those that belong to a particular set ‘$K$’ (here chosen to be all wavenumbers inside a sphere centered at the origin with radius $k_{max}$). These equations conserve exactly two quadratic invariants

[TABLE]

In Fourier space these invariants are distributed among the different modes which are quantified by the energy and helicity spherically averaged spectra $E(k),H(k)$ respectively defined as

[TABLE]

Here $\tilde{\bf u}$ is the Fourier transform of $\bf u$ and we have assumed a triple periodic cubic domain of size $2\pi/k_{0}$. The spectra have been divided by the smallest non-zero wavenumber $k_{0}$ so that they have units of energy and helicity density respectively.

At late times this system reaches a statistically steady state whose properties are fully determined by these two invariants. Using Liouville’s theorem and assuming ergodicity Lee (1952) predicted that at absolute equilibrium this system will be such that every state ${\bf u}$ of a given energy $\mathcal{E}$ is equally probable. This is equivalent to the micro-canonical ensemble in statistical physics and it leads to equipartition of energy among all the degrees of freedom (ie among all Fourier amplitudes) and an energy spectrum given by $E(k)\propto k^{2}$. Kraichnan (1973) generalized these results including helicity and assuming a Gaussian equipartition ensemble

[TABLE]

where $\mathcal{P}({\bf u})$ is the probability distribution for the system to be found in the state ${\bf u}$, $\mathcal{E}$ is the energy and $\mathcal{H}$ the helicity given in 2 and $\mathcal{Z}$ a normalization constant. The parameters $\alpha$ and $\beta$ are the equivalent of an inverse temperature and inverse chemical potential in analogy with statistical physics. We note that for the TEE system eq.4 is not exact! This is because eq. 4 allows for fluctuations of the energy $\mathcal{E}$ which are not allowed for the TEE system. However it becomes closer to the true distribution as the number N of fourier modes becomes larger. For no helicity $\beta=0$ this leads to the energy spectrum $E(k)\propto k^{2}$ predicted by Lee. In the presence of helicity one obtains

[TABLE]

The coefficients $\alpha$ and $\beta$ are determined by imposing the conditions

[TABLE]

and $\mathcal{E}$ and $\mathcal{H}$ are the initial energy and helicity respectively. The predictions above have been verified for the truncated Euler system in numerous numerical simulations (Orszag & Patterson, 1977; Cichowlas et al., 2005; Krstulovic et al., 2009; Dallas et al., 2015; Cameron et al., 2017; Alexakis & Brachet, 2018).

In the TEE however there is no exchange of energy with external sources or sinks and it is thus harder to make contact with the more realistic systems mentioned in the beginning of the introduction such as the scales larger than the forcing scale in a turbulent flow. It was shown recently in Alexakis & Brachet (2018) that, in some cases, although the large scales are close to an equilibrium state there is still exchange of energy with the smaller turbulent and forcing scales generating energy fluxes (from the forced scales to the large scales and from the large scales to the turbulent scales). It appears thus that, even in the presence of sources and sinks, equilibrium dynamics can still be relevant. In this work we examine further this possibility by looking at the truncated Navier-Stokes (TNS) equations where there is constant energy injection and dissipation like in the regular Navier-Stokes Equation but the system is limited to a finite number of Fourier modes as in TEE.

We show analytically in the next section that for weak energy injection and weak viscosity so the $k_{max}\eta\ll(k_{0}/k_{max})^{11/4}$ (where $\eta$ the kolmogorov lengthscale) the system indeed reaches a quasi-equilibrium state whose probability distribution $\mathcal{P}({\bf u})$ can be calculated. We verify and extend these results using direct numerical simulations in section 3. Our conclusions are presented in the last section.

2 Asymptotic expansion

We consider the truncated Navier-Stokes (TNS) equations

[TABLE]

where ${\bf u}$ is the incompressible velocity field, $P$ is the pressure, $\nu$ is the viscosity and $\bf f$ is a forcing function. The domain is a $2\pi$ periodic cube so that the smallest non-zero wavenumber is $k_{0}=1$. The projection operator $\mathbb{P}_{K}$ sets to zero all Fourier modes with wavenumbers outside the sphere of radius $k_{max}$. In total there are $N\simeq\frac{4\pi}{3}(k_{max}/k_{0})^{3}$ Fourier wavenumbers inside this sphere. In order to proceed it helps to write the truncated Navier Stokes equation in Fourier space using the Craya-Lesieur-Herring decomposition (Craya (1958); Lesieur (1972); Herring (1974)) where every Fourier mode is written as the sum of two modes one with positive helicity and one with negative helicity $\tilde{\bf u}_{\bf k}=\tilde{u}_{\bf k}^{+}{\bf h}^{+}_{\bf k}+\tilde{u}_{\bf k}^{-}{\bf h}_{\bf k}^{-}$. The two vectors $\bf h_{k}^{\pm}$ are given by:

[TABLE]

where $\bf\hat{e}$ is an arbitrary unit vector. The sign index $s=\pm 1$ indicates the sign of the helicity of ${\bf h}^{s}_{\bf k}$. The basis vectors ${\bf h}^{s}_{\bf k}$ are eigenfunctions of the curl operator in Fourier space such that $i{\bf k\times h}^{s}_{\bf k}=s|{\bf k}|{\bf h}^{s}_{\bf k}$. They satisfy ${\bf h}^{s}_{\bf k}\cdot{\bf h}^{s}_{\bf k}=0$ and $({\bf h}^{s}_{\bf k})^{*}\cdot{\bf h}^{s}_{\bf k}=1$, where the complex conjugate of ${\bf h}^{s}_{\bf k}$ is given by $({\bf h}^{s}_{\bf k})^{*}={\bf h}^{-s}_{\bf k}={\bf h}^{s}_{-\bf k}$. They form a complete base for incompressible zero-mean vector fields. This decomposition has been extensively used and discussed in the literature (Cambon & Jacquin, 1989; Waleffe, 1992; Chen et al., 2003; Biferale et al., 2012; Moffatt, 2014; Alexakis, 2017; Sahoo et al., 2017). Note that since every $\bf\tilde{u}_{\bf k}$ is described by two complex amplitudes $\tilde{u}^{\pm}_{\bf k}$ that satisfy $(\tilde{u}^{\pm}_{\bf k})^{*}=\tilde{u}^{\pm}_{-{\bf k}}$ there are in total $2N$ independent degrees of freedom. The truncated Navier-Stokes using the helical decomposition can then be written as

[TABLE]

where the nonlinear term $\mathcal{V}^{s}_{\bf k}$ is written as the convolution

[TABLE]

and the tensor $C^{s_{k},s_{q},s_{p}}_{{\bf k},{\bf q},{\bf p}}$ is given by $C^{s_{k},s_{q},s_{p}}_{{\bf k},{\bf q},{\bf p}}=\frac{1}{2}(s_{q}q-s_{p}p)({\bf h}_{-{\bf k}}^{s_{\bf k}}\cdot{\bf h}_{\bf q}^{s_{\bf q}}\times{\bf h}_{\bf p}^{s_{\bf p}}).$ The nonlinearity $\mathcal{V}^{s}_{\bf k}$ satisfies the following relations

[TABLE]

that correspond to the energy and helicity conservation respectively and

[TABLE]

This last relation indicates that phase space volume is conserved by the nonlinearity (ie it satisfies a Liouville condition). We will assume that the forcing is written as $\tilde{f}^{s}_{\bf k}={\epsilon}^{s}_{\bf k}\xi^{s}_{\bf k}$ where $\xi^{s}_{\bf k}$ are random complex amplitudes that are statistically independent, normally distributed and delta-correlated in time such that $\langle\xi^{s}_{\bf k}(t)\xi^{s^{\prime}}_{-{\bf q}}(t^{\prime})\rangle=2\delta_{s,s^{\prime}}\,\delta_{{\bf q},{\bf k}}\,\delta(t-t^{\prime})$. With this choice each forcing mode injects energy to the system on average at rate $\epsilon^{s}_{\bf k}$. Then the Fokker-Plank equation for the probability density $\mathcal{P}({\bf u})$ in the 2N-dimensional space of all complex amplitudes $\tilde{u}^{\pm}_{\bf k}$ is given by

[TABLE]

where the sum is over all 2N modes $\tilde{u}^{s}_{\bf k}$. Multiplying by $\mathcal{E}=\frac{1}{2}\sum_{s,{\bf q}}|\tilde{u}^{s}_{\bf q}|^{2}$, integrating over the phase-space volume and using integration by parts we obtain the energy balance equation

[TABLE]

where $\mathcal{D}_{\mathcal{E}}=\nu\sum_{s,{\bf k}}k^{2}\langle|\tilde{u}^{s}_{\bf k}|^{2}\rangle$ is the energy dissipation, $\mathcal{I}_{\mathcal{E}}=\sum_{s,{\bf k}}\epsilon^{s}_{\bf k}$ is the injection rate and the brackets stand for the average $\left\langle f\right\rangle\equiv\int f\mathcal{P}dU$ where $dU$ stands for the phase space volume element $dU=\prod_{s,{\bf k}}d\tilde{u}^{s}_{\bf k}$.

We are interested in the limit that the energy injection and dissipation rate are a small perturbation to the thermalized fluctuations. We thus set $\nu_{k}=\delta\nu_{k}^{\prime}$, $\epsilon^{s}_{\bf k}=\delta\epsilon^{\prime s}_{\bf k}$ where $\delta\ll 1$ is a small parameter. We then expand $\mathcal{P}$ in power series of $\delta$ as $\mathcal{P}=\mathcal{P}_{0}+\delta\mathcal{P}_{1}+\dots$ . We are going to also consider the long time limit so we can neglect the time derivative. To zeroth order we then have

[TABLE]

The equation above implies that $\mathcal{P}_{0}$ is constant along the trajectories in the phase space followed by solutions of the truncated Euler equations. These trajectories are expected to be chaotic for large $N$ and since this is such a high dimensional space we can also conjecture that these trajectories are space filling (ie ergodic) in the subspace constrained by the invariants of the system. In other words we assume that the trajectory will pass arbitrarily close to any point that has the same energy and helicity as the initial conditions. In this case $\mathcal{P}_{0}$ is determined by the energy and helicity of the system $\mathcal{P}({\bf u})=f(\mathcal{E},\mathcal{H})$. For the present work however we are going to neglect the second invariant the helicity and assume dependence only on the energy. We then write the solution of eq. (15) as

[TABLE]

To next order we then get

[TABLE]

Substituting eq. (16) and using the chain rule for the derivatives

[TABLE]

we obtain for the function $f(\mathcal{E})$

[TABLE]

To obtain a closed equation for $f(\mathcal{E})$ we average eq. (19) over the volume $dU_{\mathcal{E}}$ of all points in phase space of energy between $\mathcal{E}$ and $\mathcal{E}+d\mathcal{E}$. This consists of a spherical shell in the 2N-dimensional phase space of radius $\sqrt{2\mathcal{E}}$. Averaging over this volume leads the sum in the left hand side to drop out because it is a divergence and the trajectories determined by $\mathcal{V}^{s}_{\bf k}$ stay inside the shell. The volume integrals of terms independent of $\tilde{u}^{s}_{\bf k}$ are proportional to the shell volume $dU_{\mathcal{E}}=S_{2N}(2\mathcal{E})^{N-1}d\mathcal{E}$ where $S_{2N}$ is the surface of an unit radius 2N-dimensional sphere. Terms proportional to $|\tilde{u}^{s}_{\bf k}|^{2}$ result due to symmetry $\int_{dU_{\mathcal{E}}}|\tilde{u}^{s}_{\bf k}|^{2}dU=(2N)^{-1}\int_{dU_{\mathcal{E}}}\sum_{s,{\bf k}}|\tilde{u}^{s}_{\bf k}|^{2}dU=(2N)^{-1}S_{2N}(2\mathcal{E})^{N}$. This leads to

[TABLE]

If we set $K^{2}=\sum_{s,{\bf k}}|{\bf k}|^{2}\simeq 8\pi k_{max}^{5}/(5k_{0}^{3})$ and $\mathcal{I}_{\mathcal{E}}^{\prime}=\left[\sum_{s,{\bf k}}\epsilon^{\prime s}_{{\bf k}}\right]$ then by multiplying by $\mathcal{E}^{N-1}$ the equation simplifies to

[TABLE]

that has the bounded solution:

[TABLE]

where $\mathcal{Z}$ is a normalization constant that imposes $\int\mathcal{P}({\bf u})d{\bf u}=1$ and is given by:

[TABLE]

with $\Gamma$ being the Gamma function. We have thus recovered the Kraichnan distribution of eq. (4) with $\beta=0$ and inverse temperature given by $\alpha=\nu K^{2}/\mathcal{I}_{\mathcal{E}}$. Note that $\alpha$ depends only on the ratio of $\nu^{\prime}$ and $\mathcal{I}^{\prime}_{\mathcal{E}}$ and thus is independent of $\delta$ and we have thus dropped the primes.

It worth restating that $\mathcal{P}({\bf u})=f\left(\frac{1}{2}\sum_{s,{\bf k}}|\tilde{u}^{s}_{\bf k}|^{2}\right)$ given in eq.22 expresses the probability that the system finds itself in the particular state ${\bf u}$ with energy $\mathcal{E}=\frac{1}{2}\sum_{s,{\bf k}}|\tilde{u}^{s}_{\bf k}|^{2}$. If we would like to find the probability $P(\mathcal{E})$ of finding the system in any state ${\bf u}$ of energy $\mathcal{E}$ we need to average over all states ${\bf u}$ that have energy $\mathcal{E}$. This leads to the chi-distribution for the energy

[TABLE]

For large $N$ the distribution $P(\mathcal{E})$ in 24 is highly peaked at the mean energy

[TABLE]

As $N$ tends to infinity $P(\mathcal{E})$ becomes asymptotically a delta function centered at $\langle\mathcal{E}\rangle$. The mean energy of any mode $\tilde{u}_{\bf k}^{s}$ is given by $\frac{1}{2}\langle|\tilde{u}_{\bf k}^{s}|^{2}\rangle=\langle\mathcal{E}\rangle/2N$ Averaging over spherical shells then leads to the thermal equipartition spectrum

[TABLE]

Finally using eq. (22) one can calculate the energy dissipation

[TABLE]

and verify that the energy balance relation in eq. (14) is satisfied. The results indicate therefore that for small viscosity the truncated system will converge to the absolute equilibrium solutions of such “temperature” $1/\alpha$ so that the viscous dissipation balances the energy injection rate!

3 Numerical simulations

In this section we test the results of the previous section and extend our investigation beyond the asymptotic limit using direct numerical simulations of the TNS system of eq. (7). The simulations were performed using the ghost code (Mininni et al., 2011) that is a pseudospectral code with 2/3 de-alliasing and a second order Runge-Kutta. For all runs the energy injection rate was fixed to unity and the integration times were sufficiently long so that steady states were reached. The forcing used is random and white in time as the one discussed in the previous section. It is limited to a spherical shell of wavenumbers satisfying $(k_{{}_{F}}\leq|{\bf k}|\leq k_{{}_{F}}^{\prime})$.

Three cases were examined. In the first case small resolution runs were performed on a cubic domain with $N_{G}=32$ grid points in each direction. These runs due to their small size allow to some extend a direct investigation of probability distribution function $\mathcal{P}({\bf u})$. In the second case the simulations were performed on a larger grid $N_{G}=256$ numerical grid that, after de-aliasing, leads to $k_{max}=85$ and was forced at large scales $(k_{{}_{F}}=1,k_{{}_{F}}^{\prime}=2)$. These simulations demonstrate the transition from a forward cascade to a quasi-equilibrium state predicted in the section before.

In the third case the simulations were designed to demonstrate the presence of an inverse flux in the thermalized state. A smaller grid was used $N_{G}=128$ with $k_{max}=42$. The energy injection was at large wavenumbers $(k_{{}_{F}}=31,k_{{}_{F}}^{\prime}=35)$ while we replaced $\nu\nabla^{2}{\bf u}$ by the modified viscous term $\nu\nabla^{2}\mathbb{P}_{Q_{{}_{D}}}[{\bf u}]$ that acts only a particular spherical set of small wavenumbers $Q_{{}_{D}}$ satisfying $(1\leq|{\bf k}|\leq 4)$. With this setup energy is forced to be transported inversely from the forced wavenumbers to the dissipation wave numbers.

3.1 Small resolution runs

For these small resolution runs the energy injection rate was fixed to unity and the viscosity was set to $\nu=10^{-4}$. The maximum wavenumber for the grid $N_{G}=32$ that was used was $k_{max}=10$. Despite the small value of $k_{max}$ the total number of Fourier modes in this system is still quite high $N=5040$. It is thus still impossible to verify the predictions for $\mathcal{P}({\bf u})$ in full detain for such a high dimensional space. Nonetheless in the left panel of figure 1 we plot the probability distribution function for the energy $P(\mathcal{E})$ based on the direct numerical simulations plotted with solid brown line and based on the analytic predictions of 24. The two curves overlap indicating that both the mean value and fluctuations around it are correctly captured. Note that for this large value of $N$ the distribution is highly peaked however there is still a finite variance of $\mathcal{E}$ indicating that $\mathcal{E}$ is a fluctuating quantity at difference with the TEE where $\mathcal{E}$ is fixed by the initial conditions.

To further explore the validity of the result 22 on the right panel of figure 1 we plot the pdfs of the real and imaginary part of the modes $\tilde{u}^{s}_{\bf k}$ for the three wavenumbers ${\bf k}=(0,2,0),(0,4,0),(0,8,0)$, and for $s=\pm 1$. The results of the previous section predict that these amplitudes follow a Gaussian distribution with the same variance. There are 16 curves in total in the left panel of 1 that have the same variance and perfectly overlap with the Gaussian shown by the dashed line in agreement the prediction of eq. 22. Further more we calculated the elements of the co-variance matrix for these modes $\Sigma{i,j}=\langle X_{i}X_{j}\rangle$ (where $X_{i}$ stand for the mode amplitudes $\Re(\tilde{u}^{s}_{\bf k})$ and $\Im(\tilde{u}^{s}_{\bf k})$). The off-diagonal elements $i\neq j$ are two orders of magnitude smaller than the diagonal elements $i=j$. This indicates that the Fourier amplitudes $\tilde{u}^{s}_{\bf k}$ are independent variables with Gaussian distribution further verifying the predictions of the previous section.

3.2 Quasi-equilibrium state and forward flux

In the previous case although it was possible to verify some of the predictions of the previous section on $\mathcal{P}({\bf u})$. However the limited range of wavenumbers did not allow to test the predictions on the energy spectrum. To that end we used a series of simulations on a larger grid ($N_{G}=256$) with $\mathcal{I}_{\mathcal{E}}=1$ varying the viscous coefficient $\nu$ from $\nu=2\cdot 10^{-1}$ to $\nu=10^{-10}$. We must note here that the time to reach saturation from zero initial conditions is proportional to $T\propto\mathcal{E}/\mathcal{I}_{\mathcal{E}}$ and can be very large for small values of $\nu$. For this reason the runs with small $\nu$ started with random initial conditions with energy close to the one predicted. The same runs were repeated with slightly smaller or large energy to make sure all runs converged to the same point.

The resulting energy spectra are shown in the left panel of figure 2 for ten values of $\nu$. Dark colors indicate large values of $\nu$ while bright values indicate small values of $\nu$.

For the large values of $\nu$ the simulations are ‘well resolved’ in the sense that one can observe clearly the dissipation range where there is an exponential decrease of the energy spectrum. For $\nu=0.01$ and $\nu=0.001$ (second and third dark line from the bottom) one can also see the formation of an inertial range that displays a negative power-law close to the Kolmogorov prediction $E(k)=C_{K}\mathcal{I}_{\mathcal{E}}^{2/3}k^{-5/3}$ (where $C_{K}\simeq 1.6$ is the Kolmogorov constant Sreenivasan (1995); Donzis & Sreenivasan (2010)). As the value of $\nu$ is decreased a bottleneck at large wave numbers appears and energy starts to pile up at the smallest scales of the system. As the value of $\nu$ is decreased further this bottleneck appears to take the form of a positive power-law close to the thermal equilibrium prediction $E(k)\propto\mathcal{A}\,k^{2}$. This thermal spectrum occupies more wavenumbers as the viscosity is decreased until all wavenumbers follow this scaling. At the smallest value of $\nu$ the result is compared with the asymptotic result obtained in the previous section. The proportionality coefficient $\mathcal{A}$ can be estimated from eq. 26 to be $\mathcal{A}=5\mathcal{I}_{\mathcal{E}}/(2\nu k_{max}^{5})$ that guaranties that the energy balance condition $2\nu\int_{0}^{k_{max}}k^{2}E(k)dk=\mathcal{I}_{\mathcal{E}}$ is satisfied. Matching the thermal with the Kolmogorov spectrum we obtain that the transition occurs at the wavenumber

[TABLE]

where $\eta=(\nu^{3}/\mathcal{I}_{\mathcal{E}})^{1/4}$ is the Kolmogorov length-scale. The right panel of fig. 2 shows a comparison of this estimate with $k_{tran}$ measured from the spectra as the wavenumber at which $E(k)$ obtains its minimum. The scaling agrees very well with the results from the simulations.

In summary for values of $k_{max}\eta$ larger than $k_{0}/k_{max}$ (ie $k_{0}\eta\gg 1$) the flow is laminar displaying an exponential spectrum. For $1\ll k_{max}\eta\ll k_{0}/k_{max}$ there is the formation of an inertial range where the Kolmogorov spectrum dominates followed by the dissipation range. If $(k_{0}/k_{max})^{11/4}\ll k_{max}\eta\ll 1$ (so that $k_{tran}\gg k_{0}$) then the dissipative range no longer exist and there is coexistence of the Kolmogorov spectrum followed by a thermalized spectrum at large wavenumbers. Finally for $k_{max}\eta\ll(k_{0}/k_{max})^{11/4}$ the system is in the quasi-equilibrium thermal state predicted in the previous section. These result are summarized in the sketch in fig. 3.

The transition from a Kolmogorov spectrum to a thermal one resembles a lot the time evolution of the TEE studied in Cichowlas et al. (2005). For the TEE at early times a $k^{-5/3}$ energy spectrum develops as energy is transferred to larger and larger wavenumbers. When the maximum wavenumber $k_{max}$ is reached, the thermalized energy spectrum starts to develop displaying at intermediate times both spectral slopes $k^{-5/3}$ and $k^{2}$. The difference with the present runs is that in the Euler case the transition occurs as time is increased while in the present case we only consider the steady state, and vary the value of $k_{max}\eta$.

Similarities can also be found with the recent work on a time-reversible version of the Navier-Stokes (Gallavotti, 1996). Using shell models Biferale et al. (2018) and three dimensional simulations Shukla et al. (2019). of the time-reversible Navier-Stokes these authors found similar transitions from a Kolmogorov to a thermal quasi-equilibrium with the formation of both spectra depending on the parameter regime. The results were interpreted in terms of a phase transition, a possibility that could further explored for the present work as well.

We note that the time averaged flux from large scales to small scales is constant in the inertial range. What varies as we change the value of $k_{max}\eta$ is the amplitude of the fluctuations of the flux around this mean value of the flux. This is displayed in figure 4 where the mean flux is shown with dark line together with the instantaneous fluxes at different times with bright colors for four different values of $\nu$. For $\nu$ in the Kolmogorov turbulence regime the fluctuations of the flux are concentrated around the mean value without large deviations from it. As the value of the viscosity is decreased the fluctuations at large wavenumbers are increased until finally the fluctuations of the mean flux are orders of magnitude larger than the averaged flux for all wave numbers. The pdfs of the energy flux for different wavenumbers and values of viscosity are shown in fig. 5. Note that for a flow in equilibrium all third order quantities like the flux have zero value. In this case the mean value of the flux remains fixed but the variance of the fluctuations as $\nu$ is decreased (right panel of fig. 5) or $k$ is increased (left panel of fig. 5) increase. Thus compared to the variance of the fluctuations the mean flux becomes negligible and thus the flow is in quasi-equilibrium.

3.3 Quasi-equilibrium state and inverse flux

In most three dimensional turbulent flows the scales larger than the forcing scale are close to an equilibrium state (Alexakis & Brachet, 2018). In many instances however there is a change of dynamics at large enough scales and the flow is constrained to two dimensional dynamics (like for example turbulence in thin layers or in rotating flows) where energy tends to cascade inversely (Alexakis & Biferale, 2018). There is thus a transition from a thermal state with zero energy flux, to a state that has a finite inverse flux of energy. Close to such transition (unless the transition is discontinuous) the system has to be close to the equilibrium state with an inverse energy flux. With this motivation we examine the case where the forcing is located at small scales while the dissipation is limited in the large scales so that there is an net inverse transfer of energy.

We have thus performed simulations on a $N_{G}=128$ numerical grid that leads to a $k_{max}=48$ forcing at wavenumbers in the range $(k_{{}_{F}}=31,k_{{}_{F}}^{\prime}=45)$ with unit energy injection rate. At difference with the previous simulations the dissipation limited only to small wavenumbers that satisfy $1\leq|{\bf k}|\leq k^{\prime}=4$. The energy is thus injected at small scales and dissipated at large scales forcing that way an inverse transfer of energy.

The spectra for four different values of $\nu$ are shown in the left panel of fig. 6. The forcing and dissipation shells are indicated in the graph. In this case a spectrum $E(k)=\mathcal{A}k^{2}$ forms for all values of $\nu$ with small changes at the dissipation wave numbers for large values of $\nu$. For small values o $\nu$ the amplitude of the spectrum $\mathcal{A}$ is such so that the energy balance is satisfied that leads to $\mathcal{A}\simeq 5\mathcal{I}_{\mathcal{E}}/(2\nu k^{\prime 5})$ where $k^{\prime}$ is the maximum wavenumber where the dissipation acts. The dashed line shows this prediction for the smallest value of $\nu$. The flux fluctuations shown in the right panel of the same figure (fig. 6) are always dominant. We note however that their averaged value leads to a constant negative flux in the $k^{2}$ inertial range where no forcing or dissipation are present. It is worth noting that in these quasi-equilibrium states the direction of the transfer of energy is not determined by the nonlinear term as in Kolmogorov turbulence but only by the location of the energy source and sink. If regular viscosity is used, although a thermalized state is still reached as predicted by the previous section, there is no net inverse flux of energy. This is demonstrated in fig. 7 where the energy spectra and the energy flux is shown from simulations with regular viscosity forced at small scales.

4 Conclusions

In this work we examined the spectrally truncated Navier-Stokes equations flows that are close to equilibrium. We showed analytically that in the limit of small viscosity the statistically steady state of these flows converge to the Kraichnan (1973) solutions. We note that our prediction is not just for the spectrum in eq. (26) but for the full probability distribution $\mathcal{P}({\bf u})$ for the system to find itself in a state ${\bf u}$ given in eq. (22).

The derivation was based on two major assumptions. First we assumed ergodicity for the solutions of the TEE. This assumption appears in most calculations of classical statistical physics that although it can only be proved in very few systems it appears to be a plausible one for many systems with large numbers of degrees of freedom. For the TEE it appears at least to be in agreement with the results of numerical simulations. The second assumption we made was to neglect the effect of the second invariant: the helicity. This assumption was made in order to simplify the (rather involved) calculation. Had we kept the effect of helicity then the zeroth order solution $\mathcal{P}_{0}({\bf u})$ would have been reduced to a function of two variables $f(\mathcal{E},\mathcal{H})$ and we would have ended up with an elliptic partial differential equation to solve for $f$. However the presence of helicity would break the spherical symmetry in phase space that allowed us to calculate the involved integrals. The calculation is still feasible but much more lengthier and we leave it for future work.

The numerical investigation verified our analytical results and demonstrated that the Fourier amplitudes $\tilde{u}_{\bf k}^{s}$ indeed become independent Gaussian variables, and that the energy distribution and energy spectrum approach that of the analytic predictions as the asymptotic limit is reached. Furthermore, the numerical investigation also shed light on how the TNS system transitions from the classical Kolmogorov turbulence state to the thermalized solutions of Kraichnan (1973) as the $k_{max}\eta$ is varied and led to precise predictions on when the transition takes place. The transitions are compactly summarized in fig. 3 where the three different states “laminar”, “turbulent” and “thermal” as well as the transitions from one state to the other are clearly marked. Finally it was also demonstrated that these quasi equilibrium states are present when the dissipation is localized in the small wave numbers and the forcing at large wavenumbers forming an inverse flux of energy.

It is worth noting that a fixed amplitude flux (positive or negative) was always present in our simulations and that it was determined by the injection rate. However, as the quasi-equilibrium state is approached, the amplitude of the velocity fluctuations increases. These velocity fluctuations then lead to fluctuations of the energy flux to have variance that is much larger then the mean value making this mean flux a sub-dominant. Furthermore, the mean energy in the quasi-equilibrium state scales like $\langle\mathcal{E}\rangle\propto 1/\nu$ (see eq. 25) and therefore in the small viscosity limit the energy dissipation normalized by $u_{rms}^{3}k_{0}$

[TABLE]

becomes zero. Thus, opposed to Kolmogorov turbulence, there is no finite (normalized) energy dissipation rate in the zero viscosity limit. We have thus to distinguish the processes of energy transfer in these quasi-equilibrium states (sometimes referred as “warm cascades”) from the energy cascade that is met classical turbulent flows. These states are dominated by fluctuations and the direction of the energy transfer is not determined by the properties of the non-linearity but only the location of the sources and sinks of energy and in general is a non-local process (see fig. 14 in Alexakis & Brachet (2018)).

There are many directions in which the present results can be pursued further. First of all including the effect of helicity is crucial to have a complete description of the system. Moreover carrying out the calculation at the next order so that statistics of the fluxes can also be calculated would be equally desirable. Finally, extending these results to two-dimensional flows, where the equilibrium states can take the form of large scale condensates, is another possible direction. Such calculations, although considerably longer than the ones presented here, should still be feasible and we hope to address them in our future work.

Acknowledgements.

This work was granted access to the HPC resources of MesoPSL financed by the Region Ile de France and the project Equip@Meso (reference ANR-10-EQPX-29-01) of the programme Investissements d’Avenir supervised by the Agence Nationale pour la Recherche and the HPC resources of GENCI-TGCC & GENCI-CINES (Project No. A0050506421) where the present numerical simulations have been performed. This work has also been supported by the Agence nationale de la recherche (ANR DYSTURB project No. ANR-17-CE30-0004). This work was also supported by the research Grant No. 6104-1 from Indo-French Centre for the Promotion of Advanced Research (IFCPAR/CEFIPRA).