# Fluctuations in 2D reversibly-damped turbulence

## 1 Introduction

The fundamental laws governing the behaviour of fluids are known, and universally believed to be correctly represented by the Navier-Stokes (NS) equations. However, a clear cut connection between microscopic and macroscopic scales, which definitely justifies such a belief, has still to be made, particularly for systems subjected to nonconservative external fields . If one considers this kind of problems as merely technical, and accepts the NS equations as a valid tool for the description of a fluid’s dynamics, serious difficulties are met nevertheless in the study of the mathematical properties of such equations. For instance, the problem of global existence and regularity of solutions of the NS equations is far from being solved in as general terms as desired. This fact is not of exclusively mathematical interest, because it casts doubts on the effectiveness of the approximation methods devised to extract information from the NS equations. Therefore, what is reasonable and presently feasible is to study simplified problems, and to construct theories connecting such problems with the full NS equations or with the real dynamics of the fluids at hand. In Ref., this is stated as: “It has been recognized that a realistic goal for a statistical theory of turbulence is to determine the equations governing the dynamics of some reduced set of modes which allow calculation of fundamental quantities of the flow”.

The situation is similar, to some extent, to that of molecular dynamics, in nonequilibrium statistical mechanics. There, some progress took place when infinite reservoirs, or driving boundary conditions were replaced by artificial constraints imposed on the bulk dynamics of $`N`$-particle systems, with $`N<\mathrm{}`$ (see e.g. ). This way direct numerical simulations of the particle models become feasible, and dynamical systems theory leads to theoretical predictions  which can be tested in numerical simulations or real experiments. The trade-off of this approach is that the relevant dynamical equations do not seem to be fully justified on physical grounds: no known fundamental force acts on the particles of the system in such a way as to implement the desired constraints. Nonetheless, the results obtained this way are in excellent agreement with experience , and several arguments have been developed to explain why this should be the case. Among such arguments, we deem more convincing those invoking a kind of equivalence of ensembles, which is a well known concept in equilibrium statistical mechanics, but rather new in nonequilibrium statistical mechanics . If verified, the equivalence of ensembles guarantees that different microscopic dynamics result in the same macroscopic behaviour, thus justifying the use of one kind of dynamics or another, depending on which one is more natural to study a given physical problem.<sup>1</sup><sup>1</sup>1Consider, for instance, the success of the idealized dynamics of lattice gas cellular automata in describing quite complex hydrodynamic situations .

Similarly, in Ref., Gallavotti argues that different probability distributions characterizing the steady state properties of a fluid should yield the same values for (some, at least) macroscopic quantities, if such distributions are obtained from different microscopic models of the same fluid. In particular, inspired in part by the mentioned developments in statistical mechanics  and by Ref., Gallavotti conjectures one kind of equivalence between the NS equations and special time reversal invariant equations which he calls GNS, where the G stands for Gaussian constraint.

The paper  introduces a constrained Euler system, which encompasses a portion of the degrees of freedom of the NS equations sufficient to obtain some of the fundamental statistical properties of the fluid. The explicit form of the constrained equations is derived from the observation that in stationary isotropic turbulence, the mean energy in a narrow wave-number shell is nearly constant in time. These equations are similar to those used in nonequilibrium molecular dynamics, for driven particle systems subjected to a “gaussian thermostat” . Therefore, the properties of such particle systems should be observed to some degree in the constrained Euler system of  and, perhaps, in still more general settings. This led Gallavotti  to revive and put under a new light Ruelle’s principle for hydrodynamics . This principle had been previously extended to nonequilibrium statistical mechanincs in , where the first proof of the Gallavotti-Cohen Fluctuation Theorem (GCFT) was obtained, under the assumptions of the Chaotic Hypothesis (CH) quoted below in Section 2.

Tests of the validity of the equivalence conjecture (EC) and of the CH include a Rayleigh-Benard convection experiment by S. Ciliberto and C. Laroche , and numerical simualtions of the GOY shell models by L. Biferale, D. Pierotti and A. Vulpiani . The experiment of  is consistent with the validity of a fluctuation relation  similar to that of the GCFT; while the simualtions of  evidenced a kind of equivalence of different hydrodynamic models: the equivalence of energy cascades. A posteriori, the work by She and Jackson  can also be taken as a verification of the validity of the EC, although apparently not motivated by a general theory. Among other relevant works, notable are Refs., although they do not concern hydrodynamic equations.

In the present paper, we test both the validity of the mentioned fluctuation relation and of the EC for the NS and GNS equations. Because lengthy calculations are needed, we consider two-dimensional systems rather than three-dimensional ones. Our results confirm, and actually indicate possible extensions of both the EC and the fluctuation relation of Ref., consistently with Gallavotti’s predictions for the slope of such (linear) relation. Furthemrore, as explained in Section 5, these results provide an example in which the difference between the GCFT and an identity previously obtained by Evans and Searles  is evident.

## 2 Reversible damping and the equivalence conjecture

The CH has been introduced  in the study of $`N`$-particle systems subjected to nonconservative external forces and to “gaussian” constraints , like dynamical systems of the form

$$\dot{𝐪}_i=𝐩_i/m_i,\dot{𝐩}_i=𝐅_i^{\text{int}}+c_i𝐅^{\text{ext}}\alpha 𝐩_i,\text{for}i=1,\mathrm{},N$$
(1)   
defined on a $`(2dN1)`$-dimensional manifold $`\mathrm{\Omega }IR^{2dN}`$, $`d`$ being the dimension of the physical space. Here, $`(𝐪_i,𝐩_i)`$ is the usual notation for position and momentum of particle $`i`$, $`𝐅_i^{\text{int}}`$ represents the action of the other $`N1`$ particles on particle $`i`$, $`c_i`$ is a charge coupling particle $`i`$ to the external field $`𝐅^{\text{ext}}`$, and $`\alpha 𝐩_i`$ is a dissipation term, which allows the system to reach a stationary state. In particular, the choice

$$\alpha =\frac{1}{2K}𝐅^{\text{ext}}\underset{j=1}{\overset{N}{}}c_j𝐩_j,$$
(2)   
where $`K=_{i=1}^N𝐩_i^2/2m_i`$ is the kinetic energy of the system, implies that the internal energy of the system $`E=K+\mathrm{\Phi }`$, $`\mathrm{\Phi }`$ being the internal potential energy, is a constant of the motion. We denote by $`V_t:\mathrm{\Omega }\mathrm{\Omega }`$, $`tIR`$, the flow, so that $`tV_t\gamma `$ represents a solution of Eqs.(1) with initial condition $`\gamma `$. The dynamical system $`(V_t,\mathrm{\Omega })`$ is but an example from a wide class of systems which share remarkable properties . For example, such systems are time reversal invariant, although dissipative. This means that the involution $`i:\{𝐩_j\}_{j=1}^N\{𝐩_j\}_{j=1}^N`$ anticommutes with the time evolution: $`iV_t=V_ti`$ (reversibility); but phase space volumes contract on average (dissipativity), resulting in a multi-fractal stationary state . Such systems have been studied in detail, under the assumption that the following holds :

Chaotic Hypothesis (CH). A chaotic many-particle system or fluid in a stationary state can be regarded, for the purpose of computing macroscopic properties, as a smooth dynamical system with a transitive Axiom-A global attractor. In reversible systems it can be regarded, for the same purpose, as a smooth transitive Anosov system.

This approach led to interesting results for the macroscopic properties of given systems, directly from their microscopic dynamics. Among these are proofs of the positivity of transport coefficients , of the validity of Onsager relations<sup>2</sup><sup>2</sup>2The microscopic models considered in these studies constitute a wide class with some peculiar feature, such as the unusual dissipation terms. Therefore, it is not immediately clear that locutions borrowed from physics –such as “Onsager relations”– can be attributed the same phenomenological meaning they usually have. For instance, a minimal prerequisite for the models to be physically relevant is that their number of particles be large , while this is not assumed in . , and of the fluctuation relation of , whose hydrodynamic version Eq.(26) is given below.

Following Gallavotti’s ideas , we now consider the NS equations for a newtonian incompressible fluid:

$$\dot{𝐮}+(𝐮)𝐮=\frac{1}{\rho }p+𝐠+\nu _1\mathrm{\Delta }𝐮,𝐮=0.$$
(3)   
Here, u is the velocity field, $`\rho `$ is the fluid density, $`p`$ the pressure, g is a constant forcing term and $`\nu _1`$ is the constant viscosity. The curl of this equation gives

$$\dot{\omega }+(𝐮)\omega =(\omega )𝐮+𝐟+\nu _1\mathrm{\Delta }\omega $$
(4)   
in which $`\omega =\times 𝐮`$ is called vorticity and $`𝐟=\times 𝐠`$ represents the forcing term. If we replace $`\nu _1`$ by

$$\beta _Q(𝐮,\omega ,𝐟)=\frac{\left[\omega 𝐟+\omega (\omega )𝐮\right]𝑑𝐱}{\left(\times \omega \right)^2𝑑𝐱},$$
(5)   
we obtain a system whose total enstrophy $`Q=\omega ^2𝑑𝐱`$ is a constant of motion, and which is time reversal invariant in the sense Eqs.(1) are. Similarly to Ref., we refer to Eq.(4) with $`\beta _Q`$ in place of $`\nu _1`$, as to the GNS equations, i.e. the NS equations with a Gaussian constraint. In two spatial dimensions, the vortex-stretching term $`(\omega )𝐮`$ vanishes; hence only a dynamical equation for the third component of the vector $`\omega `$ is necessary. Denoting this component also by $`\omega `$, the vorticity equation reduces to

$$\dot{\omega }=(𝐮)\omega +f+\alpha \mathrm{\Delta }\omega .$$
(6)   
where $`\alpha `$ stands either for $`\nu _1`$ or for $`\beta _Q`$, depending on the case. We impose doubly-periodic boundary conditions on Eq.(6), and rescale the size of the system to $`2\pi `$, making natural a Fourier expansion for $`\omega `$. Such expansion, truncated for numerical implementation, yields the approximation

$$\omega (𝐱)\underset{k_x=N}{\overset{N}{}}\underset{k_y=N}{\overset{N}{}}e^{i𝐤𝐱}\omega _𝐤,𝐤=(k_x,k_y),NIN.$$
(7)   
Substituting Eq.(7) in Eq.(6), we get the equations for the Fourier modes $`\omega _𝐤`$. We write such equations in the more general form:

$$\dot{\omega }_𝐤=r_𝐤+f_𝐤\alpha 𝐤^{2l}\omega _𝐤,k_x,k_y=N,\mathrm{},N,l=0,1,2,\mathrm{}$$
(8)   
where each $`l`$ accounts for the dissipation given by a power of the Laplace operator applied to $`\omega `$. This way, we get Ekman damping for $`l=0`$ ($`\alpha =\nu _0`$), normal viscosity for $`l=1`$ ($`\alpha =\nu _1`$), and a different hyperviscosity $`\nu _l`$ for any $`l>1`$, while only $`l=1`$ was considered in . The quantity $`r_𝐤`$ is the nonlinear interaction term, which reads

$$r_𝐤=\underset{𝐩}{}\underset{𝐪}{}\frac{𝐪^{}𝐩}{𝐪^2}\omega _𝐪\omega _𝐩\delta _{𝐤,𝐩+𝐪};𝐪^{}=(q_y,q_x).$$
(9)   
Equations (8) determine the dynamics of the $`(2N+1)^2`$ complex modes $`\omega _𝐤`$, in a phase space $`\mathrm{\Omega }^{(2N+1)^2}`$ whose dimension is only $`2N(N+1)`$, because of the reality condition $`\omega _𝐤=\omega _𝐤^{}`$ and of the absence of the $`𝐤=(0,0)`$ mode.

Differently from , we consider all the quadratic global quantities, such as the energy $`E=\omega \mathrm{\Delta }^1\omega 𝑑𝐱`$, the enstrophy $`Q=\omega ^2𝑑𝐱`$, the palinstrophy $`P=(\omega )^2𝑑𝐱`$, the “hyperpalinstrophy” $`H=(\mathrm{\Delta }\omega )^2𝑑𝐱`$, etc, instead of considering $`E`$ and $`Q`$ only. Using the spectral notation, we can write any one of them as

$$Q_m=\underset{𝐤}{}𝐤^{2m}\omega _𝐤^{}\omega _𝐤,$$
(10)   
where $`Q_1=E`$, $`Q_0=Q`$, $`Q_1=P`$, $`Q_2=H`$ etc. These quantities satisfy the evolution equations

$$\dot{Q}_m=2F_m2\alpha Q_{l+m}+2R_{m,0},m=1,0,1,2,\mathrm{}$$
(11)   
where we have defined for convenience of notation

$`F_m`$ $`=`$ $`{\displaystyle \underset{𝐤}{}}𝐤^{2m}\omega _𝐤^{}f_𝐤,`$   
$`R_{m,l}`$ $`=`$ $`{\displaystyle \underset{𝐤}{}}{\displaystyle \underset{𝐩}{}}{\displaystyle \underset{𝐪}{}}\omega _𝐤\omega _𝐩\omega _𝐪\delta _{𝐤+𝐩+𝐪,0}(𝐤^{2l}+𝐩^{2l}+𝐪^{2l}){\displaystyle \frac{𝐪^{}𝐩}{𝐪^2}}𝐤^{2m}`$ (12)   
and $`l`$ depends on the kind of viscosity under consideration. The term $`R_{m,0}`$ vanishes for $`m=1`$ and $`m=0`$, because of the symmetries of the summands: this corresponds to the conservation of energy and enstrophy in inviscid two-dimensional flow. Now, if we take

$$\alpha =\beta _{l,m}=\frac{F_m+R_{m,0}}{Q_{l+m}},$$
(13)   
instead of the constant viscosity coefficient $`\nu _l`$ in Eq.(8), $`Q_m`$ becomes a constant of motion. In this case, we refer to Eqs.(8) as to the “cut-off GNS equations”. <sup>3</sup><sup>3</sup>3Equations (8), with $`\alpha =\beta _{l,m}`$, for an arbitrary choice of $`l,m`$, cannot be interpreted in terms of Gauss’ principle, while that interpretation was allowed in the cases of Ref.. An important quantity is $`\sigma _\alpha `$, the negative of the divergence in phase space of the right hand side of Eqs.(8),

$$\sigma _\alpha =\underset{𝐤}{}\frac{\delta \dot{\omega }_𝐤}{\delta \omega _𝐤}=\alpha \underset{𝐤}{}𝐤^{2l}+\underset{𝐤}{}𝐤^{2l}\omega _𝐤\frac{\delta \alpha }{\delta \omega _𝐤},$$
(14)   
which may be identified  with the entropy production rate of the system. Then, for $`\alpha =\nu _l`$ we get

$$\sigma _l^{NS}=\nu _l\underset{k_x=N}{\overset{N}{}}\underset{k_y=N}{\overset{N}{}}𝐤^{2l}0^{2l}$$
(15)   
while for $`\alpha =\beta _{l,m}`$ (the GNS case) we have

$$\sigma _{l,m}^{GNS}=\beta _{l,m}\left(\frac{\sigma _\nu }{\nu _l}2\frac{Q_{2l+m}}{Q_{l+m}}\right)+\frac{F_{m+l}+R_{m,l}}{Q_{l+m}}.$$
(16)   
The sums in Eq.(15) evaluate to

$`(2N+1)^21`$ $`\mathrm{for}`$ $`l=0(\mathrm{Ekman}\mathrm{damping}),`$   
$`{\displaystyle \frac{2}{3}}(2N+1)^2(N+1)N`$ $`\mathrm{for}`$ $`l=1(\mathrm{normal}\mathrm{viscosity}),`$ (17)   
$`{\displaystyle \frac{2}{45}}(14N^2+14N3)(2N+1)^2(N+1)N`$ $`\mathrm{for}`$ $`l=2(\mathrm{hyperviscosity}),\mathrm{etc}.`$   
Clearly, $`\sigma _l^{NS}`$ is a constant, while $`\sigma _{l,m}^{GNS}`$ depends on the vorticity field, hence on time. Nevertheless, $`\beta _{l,m}`$ could fluctuate around an average value $`\beta _{l,m}`$ in such a rapid fashion as to be unresolved on macroscopic scales, making the macroscopic behaviour of the GNS model identical to one with given constant viscosity. This is our interpretation of the EC of  which for $`l=1,m=0`$ states that:

Equivalence Conjecture (EC). The stationary probability distributions of the NS equations and of the GNS equations are equivalent in the limit of large Reynolds number, provided the enstrophy $`Q`$ and the viscosity $`\nu `$ are so related that $`\sigma _l^{NS}`$ and the average of $`\sigma _{l,m}^{GNS}`$ are equal.

In Ref. the large Reynolds number plays the role of the thermodynamic limit, in which limit $`\beta _{l,m}\sigma _l^{NS}/\nu _l`$ in Eq.(16) is expected to dominate over the other terms, making the equality $`\sigma _l^{NS}=\sigma _{l,m}^{GNS}`$ imply $`\beta _{l,m}=\nu _l`$, for $`l=1,m=0`$. We further conjecture that the same should hold for any $`l`$ and $`m`$.

How can these reasonings be justified? In the first place the solutions to the NS equations are expected to take a form, in the $`\nu _l0`$ limit, which is independent of the viscosity and of the boundary conditions, and therefore universal. This is because the flow is expected to develop on characteristic length scales which are much smaller than the system size but much larger than the small scale at which the dissipation is effective. This limit is then natural in hydrodynamics, as it is desirable for the EC to have a universal rather than particular character. Moreover, in standard computational practice, the limit $`\nu _l0`$ is connected to the $`N\mathrm{}`$ limit: a numerical simulation of turbulence is considered reliable when the dissipative scale $`1/k_d`$ is well resolved ($`N>Lk_d`$, $`L`$ being a length scale of the system). The value $`k_d`$ is estimated by applying Kolmogorov’s theory, which assumes a power law scaling in a range of $`k`$ called inertial; this scale becomes smaller as $`\nu _l0`$. However, the GNS case is more garbled, as it is not known how the scaling is affected by the truncation and by the reversible dissipation. Arguments based on power counts shed little light on this, at least for the easier cases in which $`R_{m,l}=0`$: for increasing $`N`$, $`\sigma _\nu ^{NS}`$ behaves as

$$\sigma _l^{NS}=\nu _l\underset{k_x=N}{\overset{N}{}}\underset{k_y=N}{\overset{N}{}}𝐤^{2l}2\nu _l\pi _1^Nk^{2l+1}𝑑k\nu _l\frac{N^{2l+2}1}{2l+2}.$$
(18)   
In order to derive scaling relations in $`N`$ for the quantities involved in the definition of $`\sigma _{l,m}^{GNS}`$, we can assume that the shell energy $`E(k)`$, defined as

$$E(k)=\underset{[|𝐤|]=k}{}\frac{|\omega _𝐤|^2}{𝐤^2},$$
(19)   
where $`[]`$ denotes the nearest integer, scales in $`k`$ as $`k^\lambda `$. Such a power law provides at least a standard term of comparison. The reference theory for $`2D`$ turbulent flows, due to Batchelor, Kraichnan and Leith , predicts $`\lambda =3`$. This slope would appear only in the inertial range, and a steeper decay of $`E(k)`$ is expected in the dissipative range $`k>k_d`$. It can then be inferred that, for $`N`$ large, but smaller than $`Lk_d`$, $`Q_n`$ scales in $`N`$ as

$$Q_n=\underset{k_x=N}{\overset{N}{}}\underset{k_y=N}{\overset{N}{}}𝐤^{2n}|\omega _𝐤|^2\underset{k}{}k^{2n+2}E(k)2\pi _1^Nk^{2n+\lambda +2}𝑑k\frac{N^{2n+\lambda +3}1}{2n+\lambda +3}.$$
(20)   
This implies that $`Q_n`$ converges for large $`N`$ if $`\lambda <2n+3`$, while it would diverge in the opposite case. This divergence might indeed not be observed, as the limit $`N\mathrm{}`$ has to be taken before $`\nu _l0`$, i.e. a dissipative range is always included. As for $`F_m`$ in Eq.(13), we argue that for a given forcing concentrated on a few (or just one) $`𝐤`$ vectors, the term $`F_m`$ can depend on $`N`$ only if the magnitudes of $`|\omega _𝐤|`$ or the phase between $`\omega _𝐤`$ and $`f_𝐤`$ change when adding more ultraviolet modes. If we assume that the (large scale) vorticty field is not significantly altered by the (much smaller) additional large $`N`$ modes, then also $`F_m`$ should not change significantly while increasing $`N`$. Eq.(13) therefore says that

$$\beta _{l,m}\frac{1}{Q_{l+m}}\frac{1}{|N^{2m+2l+\lambda +3}1|},m=1,0.$$
(21)   
Our results reported in Section 4 show that $`\beta _{l,m}\nu _l`$, for all $`N`$, so we should have $`\lambda <2m2l3`$ (e.g. $`\lambda <3`$ for constrained energy and $`\lambda <5`$ for constrained enstrophy and normal viscosity), provided that $`E(k)k^\lambda `$ and within the inertial range. With this constraint for $`\lambda `$, we have

$$\frac{Q_{2l+m}}{Q_{l+m}}=\{\begin{array}{cc}𝒪(N^{2l})\mathrm{if}0>\lambda +2m+2l+3>2l\hfill & \\ 𝒪(N^0)\mathrm{if}\lambda +2m+2l+3<2l\hfill & \end{array}$$
(22)   
The power count in eq. (16) then shows that for $`l0`$ we have

$$\sigma _{l,m}^{GNS}=\frac{\beta _{l,m}}{\nu _l}𝒪(N^{2l+2}),$$
(23)   
so that for large $`N`$ (depending on $`\nu _l`$) $`\sigma _{l,m}^{GNS}`$ is proportional to $`\beta _{l,m}`$, in the leading order. Indeed this is correct for any $`N`$, if $`l=0`$, and numerical inspection shows that, in all our simulations (made with $`l=1`$ and $`m=1,0,1,2`$), the ratio $`\nu _l\sigma _{l,m}^{GNS}/(\sigma _l^{NS}\beta _{l,m})`$ differs from one by no more than a few parts per thousand, already at truncations with $`N`$ as low as 3. We can therefore fix values of $`\nu _l`$ and number of modes $`N`$, and check the equality $`\beta _{l,m}=\nu _l`$ needed for the macroscopic equivalence of the NS and GNS systems, from which the equality $`\sigma _l^{NS}=\sigma _{l,m}^{GNS}`$ follows.

Whether we take the thermodynamic limit or not, there is a formal similarity between the cut-off GNS equations and Eqs.(1), which may then enjoy similar properties. Of these, the validity of the GCFT of  is perhaps the most striking one, and in  Gallavotti has argued that a similar relation (introduced in ) should hold for the fluctuations of $`\sigma _{l,m}^{GNS}`$, defined by Eq.(16) with $`l=1`$ and $`m=0`$. Because of Eq.(23), a similar fluctuation law should also hold for the fluctuations of $`\beta _{l,m}`$ (possibly for any $`l,m`$), which has a more direct physical interpretation than $`\sigma _{l,m}^{GNS}`$, since $`\beta _{l,m}`$ represents the viscosity in the GNS equations.

To test these ideas in our simulations of a $`2D`$ fluid, we need the following definitions. Let $`tV_t\omega `$ be the time evolution for an initial vorticity field $`\omega `$, take an overall simulation time $`T`$ which must be adequately longer than the characteristic time of the fluctuations of $`\beta _{l,m}`$, and assume that relaxation to a statistically stationary state takes place in a time $`t_0T`$. Then, introduce the time average

$$\beta _{l,m}=\frac{1}{T}_0^T\beta _{l,m}(V_t\omega )𝑑t.$$
(24)   
We omit to indicate the dependence of $`\beta _{l,m}`$ on $`T`$ and on $`\omega `$, assuming that $`T`$ is large enough, and the dynamics ergodic enough, that larger simulation times and different randomly chosen initial fields would produce negligible differences in the end result. Let us subdivide the time interval $`[t_0,T]`$, in a number of subintervals of length $`\tau `$, and consider the quantities

$$\overline{\beta }_{l,m}^\tau (i)=\frac{1}{\beta _{l,m}\tau }_{t_0+(i1)\tau }^{t_0+i\tau }\beta _{l,m}(V_t\omega )𝑑t;i=1,\mathrm{},\frac{Tt_0}{\tau }.$$
(25)   
These values, arranged in a histogram, allow us to construct the probability distribution $`\pi _{l,m}^\tau `$ of the fluctuations of $`\beta _{l,m}`$ integrated over a time $`\tau `$. Then, a test of the validity of the fluctuation relation Eq.(4.1) of  amounts in our setting to verify that the values

$$C_\beta ^{l,m}(\tau ,p)=\frac{1}{\tau \beta _{l,m}p}\mathrm{log}\frac{\pi _{l,m}^\tau (p)}{\pi _{l,m}^\tau (p)},$$
(26)   
obtained for different values of $`\tau `$ and $`p`$ (in the support of $`\pi _{l,m}^\tau `$), do converge for growing $`\tau `$ to a definite value $`c_\beta ^{l,m}`$ independent of $`p`$. The fluctuation relation of  is then a statement on the form of the measured p.d.f. $`\pi _{l,m}^\tau `$ for large $`\tau `$. We found that $`C_\beta ^{l,m}(\tau ,p)`$ can be fitted by a $`p`$-independent function of $`\tau `$, $`C_\beta ^{l,m}(\tau )`$ say. Therefore, the odd part of $`\mathrm{log}(\pi _{l,m}^\tau (p))`$ appears to be linear in $`p`$ at all $`\tau `$. Then, to assess the validity of the fluctuation relation, we verify whether the slope $`C_\beta ^{l,m}(\tau )`$ of the linear function of $`p`$

$$\frac{1}{\tau \beta _{l,m}}\left[\mathrm{log}\left(\pi _{l,m}^\tau (p)\right)\mathrm{log}\left(\pi _{l,m}^\tau (p)\right)\right]=C_\beta ^{l,m}(\tau )p$$
(27)   
does converge to a given $`c_\beta ^{l,m}`$, as $`\tau `$ is increased. This linear form is shared by many kinds of distributions, including distributions with exponential tails and Gaussian distributions, giving the erroneous impression that a trivial connection exists between the fluctuation relations of  and the central limit theorem (see point 4 of Section 5).

## 3 Description of the numerical method

We performed a number of numerical experiments, integrating in time several vorticity fields and using a standard protocol for the NS equations . We used a two-dimensional $`2/3`$ dealiased pseudospectral code on the periodical square $`[0,2\pi ]\times [0,2\pi ]`$. The resolution of such simulations is traditionally expressed by the number of modes $`M\times M`$ over which the convolution term is computed; the number $`M`$ is often a power of 2 for FFT convenience. With the dealiasing procedure, the vorticity field is reconstructed in physical space using only the modes with $`M/3k_x,k_yM/3`$. The formulas given above have therefore to be read with $`N=M/3`$: for example, $`M=32`$ implies a system of 221 complex coupled modes. Time advancement in our algorithm is accomplished by a fourth order Runge–Kutta integrator. This NS code is easily modified to integrate the GNS equations; even the convolution term $`R_{m,l}`$ in (12) can be computed via a multiplication in physical space, as is done according to the pseudospectral method for the nonlinear term $`n_k`$ of eq. (9). However, in the GNS case the analytical integration of the linear viscous term, often done to allow longer integration steps for the NS equations, is not possible. Moreover, the mere substitution of the constant $`\nu _l`$ by $`\beta _{l,m}`$, to be computed at each substep of the Runge–Kutta algorithm, does not guarantee the exact conservation of the constrained global quantity $`Q_m`$. This is not a real problem, because the relative error $`ϵ=\left|Q_m(t_0)Q_m(t_0+T)\right|/Q_m(t_0)`$ made in the conservation of $`Q_m`$ can be controlled by an adequate reduction of the integration time steps. In Table 1, the values of $`ϵ`$ for each run are explicitly given, confirming that our accuracy in the conservation of $`Q_m`$ is satisfactory.

To study the properties of statistically stationary, nonequilibrium states both for the NS and for the GNS systems, we need a criterion to judge if a stationary state has been reached in a simulation. We found that, for the NS system, the approach to a steady state can be assessed by looking at the time series of the quantity $`Q_m`$ which is going to be constrained in the GNS system. Fluctuations of such quantities in the NS system have amplitudes and correlation times which vary from case to case, but they take relatively short times to settle around an average value, if $`\nu _l`$ is not too small and the forcing is not too large.<sup>4</sup><sup>4</sup>4This and the fact that sufficient numerical accuracy requires very small integration time steps if $`\nu _l`$ is small and $`f`$ large, determines the range of parameters which we can investigate. As an empirical criterion, we then assume that the NS system has reached a stationary state when $`Q_m`$ fluctuates around a definite value. Figures 1 and 2 show how this is typically achieved.

To perform our NS simulations, we have first to choose the values of the viscosity coefficient $`\nu _l`$, of the force $`f_𝐤`$, and of the number of active modes, i.e. of $`N`$. The correspondig GNS system is obtained keeping the same $`f_𝐤`$ and $`N`$, discarding $`\nu _l`$, and choosing the global quantity to be fixed which, in our simulations, was $`E`$ or $`Q`$ or $`P`$ or $`H`$. Once this is done, we let the NS system evolve towards a steady state using a “buildup” approach: starting from a small random initialization of the $`\omega _𝐤`$ modes, we numerically integrate the field for a sufficiently long time, $`t_0`$ say. Simulations are considered sufficiently long in this context, when they are long with respect to the correlation times of the fluctuating series and to the time required to stabilize the average of the considered $`Q_m`$. A snapshot of the vorticity field at an arbitrary time $`t>t_0`$, is then used as initial datum for the GNS runs.

In order to perform longer runs, we minimized the number of modes of the vorticity, i.e. the resolution. Some numerical experiments led us to conclude that, given a force $`f_𝐤`$, an adequate resolution should include the spectral modes up to at least $`2𝐤`$. Thus, for instance, we performed many of our runs with forcings on $`𝐤=(\pm 3,\pm 4)`$ and $`𝐤=(\pm 4,\pm 3)`$ at a resolution of $`M=32`$, which, taking into account the dealiasing, implies a truncation of the Fourier expansions at $`N=10`$. Such resolution allows us to carry on simulations for up to some million timesteps, i.e. to total simulation times $`T`$ of the order of $`10^5`$ time units, whereas fluctuations show characteristic times of the order of ten units. Runs at higher resolution converge to a steady state with approximately the same values of the chosen $`Q_m`$, in roughly the same time and, in general, with similar fluctuation amplitudes. The values of the Reynolds numbers (computed as $`\frac{E}{\nu _1\sqrt{Q}}`$) are not particularly relevant in our context; however, they are observed to be of the order of $`10^2`$ up to $`10^4`$, in the various runs. The complete list of runs performed, inclusive of simulation parameters and results, is presented in Table 1.

## 4 Tests: equivalence conjecture and fluctuation formula

Our GNS systems, with our choice of parameters, take a reasonable computation time to reach a stationary state, i.e. a state in which the fluctuations of $`\beta _{l,m}`$ occur about its time average $`\beta _{l,m}`$. This quantity is found to coincide (apart from small errors, and from the cases with insufficient statistics) with the value of the viscosity $`\nu _l`$ of the NS equations, if the quantities $`E,Q,P`$ and $`H`$ are fixed to the value measured in the final NS snapshot. This fact is quite remarkable if we consider that the only information we pass from the NS to the GNS system is the NS field at an arbitrary time $`t`$ after the stationary state has been reached, and that two consecutive snapshots may yield rather different values of the chosen $`Q_m`$. In doing so, we assumed that any phase space point of the steady state trajectory carries enough information for the GNS run, and our results validate this approach.

The fact that good agreement is obtained when imposing the conservation of a quantity like $`E`$ or $`Q`$ is not entirely surprising, and might be connected with the notion of “rugged invariants”, i.e. of global quantities which are approximately conserved in spontaneous turbulent decay. While the constancy of such quantities has been imposed for structural reasons here, that was earlier the key assumption of a number of theories of turbulence. Without attempting to make a review, we remark that the stationarity of $`E`$, together with that or not of $`Q`$, were the basis of various models of either stationarily forced or freely decaying two-dimensional turbulence. Among these we quote: a) the equilibrium statisical mechanics of Fourier modes, with fixed $`E`$ and $`Q`$, due to Kraichnan ; b) the minimum enstrophy theory, in which the decay is ruled by a minimization process of $`Q`$ at fixed $`E`$ ; c) the maximal entropy theory, in which the decay involves again a lowering of $`Q`$ at constant $`E`$, but the decrease in $`Q`$ is explained probabilistically as loss of information due to mixing . In our cases, however, the hierarchy of conserved quantities implied in these theories seems to be absent: the equivalence between the stationary state of the GNS and of the NS equations is found for all constrained $`Q_m`$. On the other hand, the amplitude of the fluctuations of $`\beta _{l,m}`$ increases with $`m`$, as can be seen comparing Figs.1 and 2.

Another interesting fact is that in the NS equations we can fix the value of $`f`$ and of $`\nu _l`$ independently of each other, while the amplitude of the fluctuations of $`\beta _{l,m}`$ grows as the ratio $`f/\nu _l`$ increases (see Fig.3). This is different from the case of particle systems, where the dissipative term is proportional to the forcing term. Also, differently from the case of particle systems, the quantity $`\dot{\omega }_k`$ in the cut-off GNS equations has $`k^{2l}\beta _{l,m}\omega _k`$ instead of just $`\beta _{l,m}\omega _k`$. It is striking, then, that the equality between $`\beta _{l,m}`$ and $`\nu _l`$ is verified irrespective of the value of $`f`$, as seen varying $`f`$ over several orders of magnitude.

To check the validity of the fluctuation relation Eq.(27), we considered a number of long runs, which we cut in segments of given length $`\tau `$, for various values of $`\tau `$, for different realizations of the forcing term, and for different values of the viscosity $`\nu _1`$. Moreover, we considered four different constraints in the cut-off GNS equations: $`E`$, $`Q`$, $`P`$ or $`H=`$constant. Sample histograms of the distribution of the averages of $`\overline{\beta }_{1,1}^\tau `$ over the times $`\tau `$ are shown in Fig.4 for constant $`E`$ and a few different values of $`\tau `$. Clearly, the distributions become more and more peaked around the mean $`\beta _{1,1}`$ as $`\tau `$ is increased (large fluctuations of $`\beta _{1,1}`$ are wiped out by the average). Such distributions deviate however from gaussian, and the non-gaussian character is not mildened by the averaging; the actual computation shows that the kurtosis $`(\overline{\beta }_{l,m}^\tau 1)^4/(\overline{\beta }_{l,m}^\tau 1)^2^2`$ is significantly different from 3 (the value for the gaussian distribution), and that this difference seems to increase both with $`\tau `$ and $`f`$ (see Fig.4). This is not exactly the case for particle systems (at moderate forcings), whose probability distributions can be interpolated quite well by Gaussians, although they should not be Gaussian distributions . Similarly to previous works  for particle systems, we do not separate in time the adjacent segments to decorrelate them. The reason is that, in all the cases we considered, the only effect of the decorrelation was the worsening of data statistics. Our test consisted in computing the left hand side of Eq.(27) for many values of $`p`$, in verifying the consistency of the results with a linear law of slope $`C_\beta ^{l,m}(\tau )`$, and in extrapolating $`C_\beta ^{l,m}(\tau )`$ to large $`\tau `$. Given the observed proportionality between $`\beta _{l,m}`$ and $`\sigma _{l,m}^{GNS}`$, Eq.(27) holds for the distributions of $`\sigma _{l,m}^{GNS}`$, with the slope $`C_\beta ^{l,m}(\tau )`$ replaced by $`C_\sigma ^{l,m}(\tau )=C_\beta ^{l,m}(\tau )\nu _l/\sigma _l^{NS}`$.

Our results show that the distributions of the values $`\overline{\beta }_{l,m}^\tau `$ are in all cases consistent with the linear law Eq.(27). The fit to a straight line, performed according to standard least-squares minimization, is always remarkably good, as in the cases shown in Fig.5. The error bars are derived from the errors on the histograms, which, according to standard statistical analysis, are estimated as the square root of the bin count. The value of the linear regression coefficient vs. $`\tau `$ is plotted in Fig.6 for several cases. The slope of these lines decreases with $`\tau `$, and converges to a finite value at large $`\tau `$ when $`E`$, $`Q`$, $`P`$ or $`H`$ are fixed. Therefore, the fluctuation relation is valid in these cases, although presently we cannot compare our results with Gallavotti’s theoretical prediction  which links $`c_\sigma ^{l,m}=lim_\tau \mathrm{}C_\beta ^{l,m}(\tau )\nu _l/\sigma _l^{NS}`$ to the Lyapunov spectrum of the system at hand. However, consistently with these predictions, we obtain $`0<c_\sigma ^{l,m}<1`$. We conclude this section noting that the dependence of $`C_\beta ^{l,m}(\tau )`$ on $`\tau `$ (cf. Fig.6) is different from that observed in the work , while it is similar to that observed in the particle systems of Ref..

## 5 Conclusions

1. We investigated the equivalence between the NS and GNS equations. Such equivalence is verified for a variety of cases, and actually appears more robust than expected in . Indeed, $`\beta _{1,m}`$ approaches $`\nu _1`$ at least when $`E`$ or $`Q`$ or $`P`$ or $`H`$ are fixed in the cut-off GNS equations (cf. Table 1). Moreover, neither the thermodynamic limit nor a large number of modes is needed to obtain this result. This leads us to conclude that the EC can be extended to the case of constrained $`Q_m`$, for any $`m1`$. The equivalence of the NS and GNS equations can be equally expressed in terms of $`\beta _{l,m}`$ and $`\nu _l`$, or of $`\sigma _{l,m}^{GNS}`$ and $`\sigma _l^{NS}`$, due to the proportionality between $`\sigma _{l,m}^{GNS}`$ and $`\beta _{l,m}`$ expressed by Eq.(23). We preferred to use $`\beta _{l,m}`$, which is directly connected with a quantity of ordinary hydrodynamic interest: the viscosity of the fluid $`\nu _l`$. Moreover, the values of $`c_\beta ^{l,m}`$ appear to cluster around a single value, independent of the resolution of the run and of the parent viscosity, as can be inferred from Table 1 using $`c_\sigma ^{l,m}=c_\beta \nu _l/\sigma _\nu ^{NS}`$. This is not the case for $`c_\sigma ^{l,m}`$.

2. The fluctuation relation Eq.(27) is seen to hold when $`E,Q,P,H`$ are constrained. In these cases, the left hand side of Eq.(27) is perfectly fitted by straight lines at all values of $`\tau `$, (see e.g. Fig.5) with slopes converging to definite values $`c_\sigma ^{l,m}`$ in the limit of large $`\tau `$ (cf. Fig.6 and Table 1). This at one time validates the CH in the framework of hydrodynamics, and yields a new result related to the NS equations. Moreover, the values $`c_\sigma ^{l,m}`$ are consistent with Gallavotti’s predictions (cf. Eq.(4.1) of ), which imply $`0<c_\sigma ^{l,m}1`$.

3. Our results concern statistically stationary states of the dynamics which, from a fluid–dynamicist point of view, are achieved at very late times. The resulting flows are substantially dominated, in absence of a proper large-scale dissipation mechanism, by vorticity structures of the size of the computational domain. These large scale structures are seen to be unsteady, fact which is by itself responsible for the fluctuations of the system. In such conditions the mechanism of transfer of energy toward small $`k`$ has already pushed almost all the avaliable energy to the largest accessible scales; the formation of a self similar direct enstrophy cascade or of an inverse energy cascade are prevented. If one is concerned with universal inertial range properties of the turbulence, such kind of states would be inappropriate because they are dependent on the finite size boundary conditions, hence not universal. Our point of view is indeed different, and we do not consider the poor resolution or the lack of a universal turbulent cascade as shortcomings. We are interested in the statistical behavior of a given dynamical system, and the computational limitations on $`N`$ do not hinder our investigation. In this spirit, we treated the viscosity coefficient as a mere parameter of the simulation; we did not match resolution and viscosity in order to guarantee the resolution of a dissipative range. Again, the point is that we investigate special properties of the vorticity modes which, because of their dynamical nature, should not depend on the level of resolution.

4. Some confusion is sometime made on the nature of the GCFT of Ref., because the p.d.f. obtained for the fluctuations of the entropy production rate in particle systems are indistinguishable to the eye from Gaussian distributions. This led some to believe that the statement of a linear law for the odd part of the p.d.f. of the fluctuations in particle and hydrodynamic systems is a quite generic result and just a manifestation of the validity of the central limit theorem. This belief, however, is not well founded. In the first place, the relation investigated here and the GCFT concern large and not small deviations, hence Gaussian distributions should not be expected. For instance, these distributions have finite support for $`m=1,0`$, because the values of $`\beta _{l,m}`$ are bounded, being a ratio between linear and quadratic functions of the vorticity field. Furthermore, differently from , the kurtosis of our distributions are manifestly different from the Gaussian value $`3`$ (cf. Fig.4) at finite $`\tau `$’s, and we haven’t observed any convergence to $`3`$ for growing $`\tau `$. Even in the case of , in which a Gaussian interpolation of the data is possible, a connection with the central limit theorem is far from obvious.

5. Two recent papers,  and , considered the verification of the EC and of the fluctuation relation in fluid dynamics, suggested by Gallavotti . The present work strengthens these earlier findings, along with new results. In , one experimental fluctuating time series, representative of a heat flux, is analyzed. The p.d.f. of this series is evidently non-gaussian, but its odd part can be fitted by a linear relation, as in our cases. The asymptotic value of the slope of the fitting lines is not given. In contrast, we were able to explore several cases and to check how the asymptotic $`c_\sigma ^{l,m}`$ are reached.

The work of , instead, is mainly concerned with the relation between ordinary shell models and reversibly-damped counterparts. The authors find differences in the statistics produced by the two versions of the shell models, as well as limitations on the values which the parameters can take for one kind of equivalence to hold. Such equivalence is found in the correct reproduction of the energy cascade and of the multiscaling of the structure functions. No analysis of the GCFT is attempted there, at variance with our paper. It is however interesting to compare our and their approach, as it seems that the dynamical equivalence between “normal” and “reversibly constrained” dynamical systems representing fluid turbulence is not always as complete as in the NS-GNS cases.

6. One final remark concerns the debate on the connection between the GCFT and the so called Evans-Searles identity (ESI) . The ESI is a relation which concerns time reversible dynamical systems, and the Liouville mesaure $`\mu _L`$ on the phase space $`\mathrm{\Omega }`$ of such systems. Hence it can also be applied to our truncated GNS equations. In particular, let $`E_p\mathrm{\Omega }`$ be the subset of initial conditions of trajectories along which the phase space contraction after an evolution of any length $`T`$ is $`e^{p\sigma _{l,m}^{GNS}T}`$. Then, the ESI implies that 

$$\mathrm{log}\mu _L(E_p)\mathrm{log}\mu _L(E_p)=p\sigma _{l,m}^{GNS}T.$$
(28)   
This equation is formally similar to Eq.(27) re-written in terms of $`\sigma _{l,m}^{GNS}`$, for $`T=\tau `$ and $`C_\sigma ^{l,m}(\tau )=1`$. Now, the fluctuation relation of the GCFT also has a slope of $`1`$, in the restricted domain of dynamical systems with dense attractors not considered here. This led some to believe that the ESI and the GCFT describe the same quantities, namely the fluctuations of the phase space contraction. However, our results explicitly show that the ESI cannot describe these fluctuations, because we have $`c_\sigma ^{l,m}<1`$. The ESI, instead, concerns the relative probability of independent “trajectory histories” emanating from the Liouville distribution.

## Acknowledgements

We owe special gratitude to G. Gallavotti for suggesting this problem, and for illuminating discussions throughout the time we devoted to this work. We also thank G. Boffetta, F. Bonetto A. Celani and A. Vulpiani for stimulating and constructive criticism. We thank the Institute for Cosmogeophysics of the CNR, Torino, for access to its computer facilities. L.R. acknowledges support from GNFM-CNR (Italy) and through the EC contract ERBCHRXCT940460.