[

# [

[
###### Abstract

We investigate the statistical properties, based on numerical simulations and analytical calculations, of a recently proposed stochastic model for the velocity field of an incompressible, homogeneous, isotropic and fully developed turbulent flow. A key step in the construction of this model is the introduction of some aspects of the vorticity stretching mechanism that governs the dynamics of fluid particles along their trajectories. An additional further phenomenological step aimed at including the long range correlated nature of turbulence makes this model dependent on a single free parameter that can be estimated from experimental measurements. We confirm the realism of the model regarding the geometry of the velocity gradient tensor, the power-law behaviour of the moments of velocity increments (i.e. the structure functions) including the intermittent corrections and the existence of energy transfer across scales. We quantify the dependence of these basic properties of turbulent flows on the free parameter and derive analytically the spectrum of exponents of the structure functions in a simplified non dissipative case. A perturbative expansion in power of shows that energy transfer, at leading order, indeed take place, justifying the dissipative nature of this random field.

A dissipative random velocity field for fully developed fluid turbulence]A dissipative random velocity field for fully developed fluid turbulence ]Rodrigo M. Pereira, Christophe Garban, Laurent Chevillard

## 1 Introduction

Fluid turbulence is an archetypal phenomenon belonging to out-of-equilibrium and non-linear classical physics. Starting probably with the work of Reynolds, the complex and multiscale nature of turbulent velocity fluctuations is usually apprehended in a statistical way. In this spirit, Kolmogorov (1941) proposed in his seminal article a dimensional based argument explaining the spatial two-point correlation structure of velocity fluctuations, i.e. the 2/3-law (see for instance classical textbooks Batchelor 1953; Tennekes & Lumley 1972; Frisch 1995; Pope 2000), as it was observed in early experimental measurements of laboratory flows. Furthermore, he derived rigorously from the Navier-Stokes equations, using the stationary solution of the von Kármán and Howarth equation, the behaviour at infinite Reynolds number and vanishing scale of the third order moment of velocity increments, known as the 4/5-law (Frisch 1995), reminiscent of a non vanishing mean energy transfer across scales. This gives solid ground to the following phenomenology of three dimensional homogeneous, isotropic and incompressible turbulence: energy is injected at large scale , the so-called integral length scale (typically the mesh size of a grid generated turbulence in a wind tunnel, or the typical scale of propellers, etc.), that is then transferred to smaller scales via a direct cascading process, until it is dissipated by viscosity.

In this context, a more complete picture could be given while proposing a stochastic representation of a velocity field able to reproduce in probability law the formerly described spatial structure of turbulence. In other words, we ask whether it is possible to build a random vector field , incompressible, statistically homogeneous and isotropic, seen as a statistically stationary solution of the Navier-Stokes equations, that reproduces in particular the observed 2/3 and 4/5-law.

The very first idea would be to consider a Gaussian approximation. This was first considered by Kolmogorov himself, and the respective process belongs to the more general class of fractional Brownian motions (Mandelbrot & Van Ness 1968). Unfortunately, such a Gaussian model fails to reproduce the observed and derived mean energy transfer encoded in the third-order moment of the velocity increments, as previously mentioned. However, an underlying Gaussian velocity field is an appealing starting point, and we will see in the following how to modify it in order to obtain a more realistic picture that includes energy transfer.

Hereafter, we consider homogeneous, isotropic and incompressible velocity fields in three dimensional space, i.e. and , with a regularizing scale that plays the role, in a schematic way, of the Kolmogorov dissipative length scale. In a Gaussian framework, let us call such a field. To fully determine this Gaussian velocity field, we have to prescribe its velocity components covariance, that will in particular take into account the self-similar law of Kolmogorov (i.e. the -law). Homogeneity, isotropy, incompressibility and the self-similar property can be modelled as a stochastic integral (Robert & Vargas 2008) in the following way

 ug,ϵ(x)=−∫R3φL(x−z)x−z|x−z|5/2−Hϵ∧W(dz) , (1.0)

where is a zero-average vector Gaussian white noise whose components are independent with variance equal to the infinitesimal volume , and is a large-scale cut-off of characteristic extension (i.e. the integral length scale) ensuring a finite variance of this random velocity field. It is chosen radially symmetric to ensure isotropy, that is, for any vector , . The singular kernel is regularized over the scale such that is proportional to when (see section 3 and Robert & Vargas 2008 for further details). This gives differentiability to the vector field for any . The vector product entering equation (1) recalls the structure of the Biot-Savart law and ensures incompressibility (i.e. a divergence-free vector field). Using standard rescaling techniques (as done in Robert & Vargas 2008, using similar techniques from Mandelbrot & Van Ness 1968), it can be shown that the limiting process , which corresponds in a turbulent context to the limit of vanishing viscosity, is a finite variance random vector field when the so-called Hurst (or Hölder) exponent is strictly positive, i.e. . Let us introduce the velocity increments in order to make a connection with turbulence phenomenology and comment on the free parameter entering the definition of the field from equation (1). As is usually done in the turbulence literature (see for instance classical textbooks such as Batchelor 1953; Frisch 1995), it is convenient to define the longitudinal and transverse velocity increments. Note and , the projections of the vector onto the direction of and, correspondingly, onto any perpendicular direction. The longitudinal and transverse velocity increments are given by

 δ∥ℓuϵ(x)=uϵ∥(x+ℓ/2)−uϵ∥(x−ℓ/2) and δ⊥ℓuϵ(x)=uϵ⊥(x+ℓ/2)−uϵ⊥(x−ℓ/2) . (1.0)

It can be shown, for , that the limiting Gaussian random field , as we will recall in this article, is scale-invariant in the sense that moments of velocity increments, the so-called structure functions, behave as power-laws, i.e. for ,

 E[(δ∥ℓug)2q]∼ℓ→0Cg,∥2qℓ2qH and E[(δ⊥ℓug)2q]∼ℓ→0Cg,⊥2qℓ2qH,

with and two strictly positive constants that are universal in the sense that they depend on the order and only on the value of the cut-off function at the origin (and not on its entire shape). Based on dimensional arguments, Kolmogorov phenomenology predicts , that is , justifying the terminology of the 2/3-law. We gather the proofs in appendix B. Roughly speaking, the scale-invariance property comes from the singular power-law shape of the kernel entering the definition of the Gaussian velocity field in equation 1. We can see that, as far as Kolmogorov’s 2/3-law is concerned, we can give a clear meaning (equation 1) to a stochastic representation of a turbulent homogeneous, isotropic and incompressible velocity field. Unfortunately, such a representation is too naive to reproduce the 4/5-law, that requires a non vanishing (strictly negative) third order moment of velocity increments, whereas the Gaussian random field (1) is such that

 E[(δ∥ℓug)3]=E[(δ⊥ℓug)3]=0.

As noted, the Gaussian velocity field (1) fails to give a realistic picture of K41 phenomenology since the Gaussian structure leads to vanishing odd-order correlators and thus, to a vanishing mean energy transfer across scales. Furthermore, higher-order even correlators are themselves poorly predicted. This intrinsic non-Gaussianity of the small scales was observed in early experimental measurements (Batchelor 1953; Monin & Yaglom 1971). This was taken into account by Kolmogorov (1962) and Obukhov (1962) (hereafter referred to as KO62) while refining the K41 theory and setting a very peculiar statistical structure of the dissipation field. This is the so-called intermittency, i.e. multifractal phenomenon (see Frisch 1995 for a review on the subject).

Various approaches were developed in the past to provide a stochastic representation of such a dissipation field, starting probably with the discrete cascades models initiated by the Russian school (Monin & Yaglom 1971, see also Benzi et al. 1993; Arneodo et al. 1998). We will prefer here to follow a continuous version of these discrete cascade models that ensures homogeneity, known as the limit-lognormal model of Mandelbrot (1972), rigorously studied in the framework of Gaussian multiplicative chaos by Kahane (1985) (see Rhodes & Vargas 2014 for recent developments on this matter). The aim is to model the dissipation field as a lognormal process with a long range correlation structure of the fluctuations, as observed in experiments (Monin & Yaglom 1971; Gagne & Hopfinger 1979; Antonia et al. 1981). Gaussian multiplicative chaos consists in defining such a scalar lognormal process as the exponential of a Gaussian field with logarithmic covariance, i.e. , with being the integral length scale. It is then possible to give a clear meaning to the scalar field (Kahane 1985), where is a dimensionless free parameter of the theory. In particular, in three dimensional space, it is then possible to show that the local average of this scalar field over a ball of size is a well-posed random field whose moments are scale invariant in the sense that behaves as , with a non-linear (quadratic) function of the order . As we said, such a construction is defined up to a dimensionless free parameter known as the intermittency coefficient that can be precisely estimated on experimental signals (see for instance Chevillard et al. 2012 and references therein).

Until now, the statistical properties of turbulence that we have mentioned concern mainly the fluctuations of the longitudinal velocity profile which is accessible with traditional experimental techniques, hot-wire anemometry in particular, and do not characterize the vector nature of the velocity field. For example, at this stage, nothing is said about the peculiar correlation structure of the components of the velocity gradients tensor . We are thus asking if, furthermore, it is possible to build up a differentiable velocity field (at a finite ) consistent with two important properties of the velocity gradient tensor that are (i) the teardrop shape of the joint density of the invariants and (Tsinober 2001; Wallace 2009; Meneveau 2011) and (ii) the preferential alignment of the vorticity vector with the eigenframe of the rate-of-strain matrix.

Going beyond the Gaussian approximation (1) is a difficult matter since the mathematical theory of non-Gaussian processes is far more sophisticated. In this direction, some recent attempts by Çağlar (2007) and Hedevang & Schmiegel (2014) are interesting but it is not clear whether these vector fields exhibit energy transfer. Let us also mention the iterative procedure of Rosales & Meneveau (2008) that gives a realistic picture but which is not explicit, making analytical results out of reach at the present time. In a one-dimensional context, several models have been proposed in the literature in order to apply the discrete cascade models to reproduce synthetically the observed fluctuations of longitudinal velocity profiles, including a model for energy transfer (Juneja et al. 1994) with additional parameters and propositions to extend to spatio-temporal (Biferale et al. 1998) and Lévy-based (Schmiegel et al. 2004) stochastic representations. In a different spirit, Nawroth & Peinke (2004) propose to reconstruct velocity time series, starting from a time series at a given (small) scale and assuming a Markov property in scale. As far as we know, Robert & Vargas (2008) are the first to have proposed a compressible velocity field with non-symmetrical probability laws. To generalize their approach to incompressible velocity fields, they propose to modify the Gaussian field (1) in order to include energy transfer and intermittency effects. This was done while disturbing the vector white noise field by the scalar multifractal measure given by the multiplicative chaos. Unfortunately, for symmetry reasons, Robert & Vargas (2008) show that this incompressible intermittent velocity field has a vanishing mean energy transfer. It is tempting to think that the present picture is too heuristic to represent the complex local structure of turbulence. A further step in this direction was proposed by Chevillard et al. (2010) in which the Euler equations, and more precisely the vorticity stretching mechanism, is used in order to motivate the exponentiation of a homogeneous field of isotropic symmetric trace-free Gaussian matrices, that eventually lead to energy transfer.

Let us recall how to include some aspects of the vorticity stretching mechanism in the present picture. The Euler equations reads, in the vorticity formulation,

 DωDt=∂ω∂t+(u⋅∇)ω=\mathsfbiSω ,

where is the velocity field, solution of the Euler equation and given by the Biot-Savart law, i.e.

 u(x)=−14π∫x−z|x−z|3∧ω(z)dz,

and the deformation field is defined as the symmetric part of the velocity gradient tensor, namely , where stands for matrix transpose. In incompressible flows, the deformation field is fully determined by the vorticity field and the explicit form reads (Constantin 1994; Majda & Bertozzi 2002)

 \mathsfbiS(x)=38πP.V.∫[(x−z)⊗[(x−z)∧ω(z)]|x−z|5+[(x−z)∧ω(z)]⊗(x−z)|x−z|5]dz , (1.0)

where the integral is understood as a Cauchy Principal Value (P.V.) and is the tensor product, i.e. . The first underlying idea of Chevillard et al. (2010) is to study the implication of a linearization of the previous formulation of the Euler equation on the velocity field generated by the stretching of an initial Gaussian vorticity field (with a K41 structure) by the initial deformation field. This first motivates the use of the exponentiation of a Gaussian random field of symmetric matrices, although it was not expected from this short time study of the Euler equations to reproduce the peculiar intermittent nature of the velocity field. This structure was introduced by hand using the Gaussian multiplicative chaos that is naturally obtained while modifying the integration kernel of the deformation field (equation 1.0). This heuristic procedure, motivated by the short time dynamics of the Euler equations, leads to the following proposition of a velocity field representing a realistic local structure of turbulence:

 uϵ(x)=−∫R3φL(x−z)x−z|x−z|5/2−Hϵ∧eγ\mathsfbiXϵ(z)W(dz) , (1.0)

where is an isotropic trace-free symmetric random matrix, whose structure recalls the one of the deformation field (1.0), given explicitly by a tensor Wiener integral that we will specify later. The non dimensional constant governs the level of intermittency. Let us finally remark that a crucial step of this construction, as dictated by the short-time dynamics of the Euler equations, is the intrinsic dependence of this statistically isotropic matrix on the vector white noise . We can see, given a Hurst exponent that we will take to be to be consistent with K41 phenomenology, that the proposed stochastic model (1.0) does depend on a single free parameter that can be determined empirically. Therefore, if this vector field is to provide a realistic picture of the local structure of turbulence, this unique free parameter should govern at the same time both the intermittency phenomenon and the physics of the energy transfer, which is, as far as we know, a new type of relationship between these phenomena. We will indeed derive from a perturbative approach (section 6) that the third order moment of velocity increments is proportional to the scale, with a multiplicative factor that is itself proportional to this free parameter .

The purpose of this article is to go beyond the results obtained by Chevillard et al. (2010) in which the field (1.0) has been proposed for the first time and studied mostly numerically for a single value of the intermittency coefficient representing in a satisfactory manner the statistical properties of turbulence. As we will see in the following quick description of the various sections of the article, the proposed new material include (i) an extensive numerical study of the statistical properties of the velocity field at the smallest resolutions we were able to reach, for several values of the free parameter , (ii) an analytical derivation of the spectrum of exponents of the structure functions in the asymptotic limit of vanishing resolutions of a simplified ersatz of the field named and (iii) a perturbative approach for small able to capture some aspects of the energy transfer taking place while reconsidering the field (equation 1.0).

In section 2, we set our notation and define the field of random matrices .

In section 3, we describe the numerical procedure in order to obtain realizations of the velocity field (1.0). In short, is simulated in a periodic box of size . We rely then on the discrete Fourier transform to perform the convolutions. The matrix exponential is evaluated at each point of space using a Padé approximant with scaling and squaring. The Fast Fourier Transform (FFT) algorithm is used in its fully parallel form. We study then the numerical properties of the velocity field based on realizations up to collocation points.

In section 4, we use these numerical simulations to compute the joint density of the invariants and at various intermittency coefficients and discuss their comparison with what is obtained in laboratory and numerical flows. Similarly, we show the preferential alignment of vorticity with the intermediate eigendirection of the eigenframe of the deformation, and quantify its dependence on .

Section 5 is devoted to a joint numerical and analytical study of the intermittency phenomenon observed in the velocity field (1.0). We will indeed observe that this field is intermittent (in a sense that we will make precise in the devoted section), and its level of intermittency is given in terms of the coefficient . A rigorous derivation of the behaviour of the structure functions of the velocity field is mathematically very demanding, and even obtaining the variance of the components is a difficult task. The reason is related to the strong correlation between the exponentiated Gaussian field of matrices and the underlying vector white noise . To obtain analytical results, we study an ersatz, which has the same structure as the proposed field (1.0) but assuming the independence of the matrix and vector fields. We will call this case the independent case and note the respective velocity field . We will show in section 5 that indeed, when properly renormalized, the velocity field converges towards a finite-variance process when and we will compute its respective structure functions, obtaining

with and two strictly positive constants. Note that, asymptotically, higher order longitudinal and transverse structure functions share similar scaling behaviour. Note also that we do not obtain perfect power-laws since an additional logarithmic factor appears in the asymptotic behaviour. This factor is related to the matrix nature of the chaos and was already observed in Chevillard et al. (2013). This former scale dependence of structure functions is based on an exact calculation for , and has been extended to higher orders based on a conjecture proposed in Chevillard et al. (2013). Thus, this field allows us to understand the intermittent corrections to the scaling behaviour with respect to the Gaussian case . The independence assumption leads on the other side to vanishing third- and more generally odd-orders structure functions, namely

 E(δ∥ℓuind)3=0,

missing all the physics of energy transfer and showing that the intrinsic correlation between the matrix and vector fields in the velocity field (equation 1.0) is crucial to reproduce non vanishing third order moment. Nonetheless, we show numerically that, to fourth order, the ersatz and the full velocity field (1.0) share similar intermittent properties. Furthermore, an analytical study that takes into account finite-scale corrections is performed in order to interpret with high precision the numerical results. From the behaviour of the velocity increment flatness that we will define later on (equation 5.0), taking into account non trivial finite-scale corrections, we are led to propose the very particular value for turbulent applications in order to be consistent with experimental measurements and numerical simulations.

Section 6 is devoted to the physics of energy transfer. As we said, a rigorous study of the statistical properties of the proposed velocity field (1.0) is a difficult task. In order to discuss the important physics of the energy transfer, as required by the -law (Frisch 1995), we will rely on a perturbative analysis of this field, at a finite resolution , using the intermittency parameter as the small parameter, which is indeed the case as far as turbulence is concerned. We show that such a perturbative expansion of the longitudinal third order velocity structure function is given by

 E(δ∥ℓuϵ)3=γDϵ(ℓ)ℓ3H+oϵ(γ),

where stands for a term that depends on but depends on a higher power of than 1 (typically this term is of order by symmetry). We are then able to show that the dominating term linked to converges when towards a non trivial function which is such that

 limϵ→0Dϵ(ℓ)=D(ℓ)⟶ℓ→0D(0)=D,

with a constant. Numerical simulations show indeed such a linear behaviour, with , of the third order structure function with both the intermittency coefficient and the scale when is small and when we use the Hurst of K41, namely . This shows, up to first order in , that the proposed velocity field (1.0) exhibits energy transfer according to Kolmogorov phenomenology.

We gather in section 7 our conclusion and perspectives.

## 2 Notations and basic properties of the velocity field

In what follows, will denote the Kronecker delta and the Levi-Civita symbol. We adopt Einstein’s convention of sum over repeated indices, unless explicitly stated, and we note that .

The full vector field (1.0) reads, with index notation,

 uϵi(x)=∫ϕϵik(x−z)(eγ\mathsfbiXϵ(z))klWl(dz), (2.0)

where the kernel encodes the structure of the underlying Gaussian velocity field (1) and is given by

 ϕϵik(x)=−ϵijkφL(x)xj|x|52−Hϵ,

and the following matrix field built from the very same vector white noise that enters the construction of the underlying structure:

 \mathsfbiXϵ(x)=√1532π∫|x−y|⩽Lx−y|x−y|7/2ϵ⊗[(x−y)∧W(dy)]+[(x−y)∧W(dy)]⊗x−y|x−y|7/2ϵ, (2.0)

which is inspired by the tensor structure of the rate-of-strain matrix (equation 1.0) that stretches the vorticity vector along its path. We will motivate the use of the multiplicative factor when we give the variance and covariance of the elements of the matrix . At this stage, remark that this matrix is Gaussian since it is defined through a linear operation on a Gaussian measure (equation 2.0). It is indeed symmetric, and it is easy to check that it is trace free, according to

 tr(\mathsfbiXϵ)=√1532π∫|x−y|⩽L2x−y|x−y|7/2ϵ⋅[(x−y)∧W(dy)]=0.

The free parameter entering the velocity field (equation 2.0) plays the same role as the free parameter used in the field given by equation 12 of Chevillard et al. (2010), and their relation is .

The Gaussian white vector field , for , follows the following rules of calculation. For any suitable deterministic function , such that it is integrable along its diagonal, we have

 E∫f(x,y)Wi(dx)=∫f(x,y)E[Wi(dx)]=0,

and

### 2.1 Covariance structure of the field of isotropic matrices

Let us first show that the matrix field (2.0) is indeed homogeneous and isotropic. Homogeneity of the field follows from the convolution with the homogeneous white measure . Consider now a rotation matrix such that , where the identity matrix. Then, it can be shown that for any rotation matrix , we have . The equality in law stands for equality in probability. This shows that, in that sense, the matrix field is statistically isotropic.

As an important further characterization of the homogeneous field of matrices , we want to obtain its covariance structure, component by component. Let us first remark that all of the elements are of zero mean, which follows from the definition of the field as a convolution with a zero-mean white noise . We gather all the proofs of the following results in annex A.

We have seen that the field of matrices is statistically isotropic. Recall that each element is a Gaussian random variable, and is a symmetric matrix. Thus, the covariance structure of its elements is given by the general framework developed in Chevillard et al. (2013). We recall several key properties of this random matrix.

The first property of the elements of is the divergence of their variance with the regularizing parameter . Henceforth, we focus only on the element of the matrix. See annex A for a general discussion on the statistical behaviour of the other elements. Defining the variance of this element as , then it is easy to obtain its asymptotic behaviour when as

 σ2ϵ=E[(Xϵ11)2]∼ϵ→0lnLϵ . (2.0)

Thus, the variance of the elements of diverge logarithmically with . This situation is classically encountered in the context of multiplicative chaos (see a review on this topic by Rhodes & Vargas 2014). Similarly, the covariance of the element can be computed and we find, taking first the limit and then looking for an equivalent at small distances,

 σ2|x−y|=limϵ→0E[Xϵ11(x)Xϵ11(y)]∼|x−y|→0lnL|x−y| . (2.0)

In other words, the Gaussian random variable converges when towards a random Gaussian distribution whose covariance behaves logarithmically at small distances. We remark that the factor entering the definition of (2.0) ensures a unit-factor in front of the logarithmic behaviours seen in equations (2.0) and (2.0). The very purpose of the theory of multiplicative chaos (Rhodes & Vargas 2014) is to give a meaning to the exponential of such a field.

### 2.2 Homogeneity and isotropy

Let us now show that the velocity field (2.0), defined with the field of matrices (2.0), is indeed homogeneous and isotropic. Again, homogeneity of the vector field follows from the convolution with the homogeneous field . Consider again a rotation matrix . Then, for any rotation matrix , we have . Thus, the velocity field is statistically isotropic. As a consequence, the velocity field is of zero-mean, i.e., for any

 Euϵ=0.

## 3 Numerical procedure

As we will see in the following, a rigorous derivation of the statistical properties of the velocity field is a difficult matter. When analytical results are not possible, we will rely on numerical simulations. To do so, we perform a numerical approximation of in the periodic domain . Define as the number of collocation points in one direction. We will typically present results for , using a fully parallelized algorithm of the Fast Fourier Transform (Frigo & Johnson 2005). The elementary volume is given by . Standard algorithms allow to generate independent realizations of a zero-mean Gaussian variable of variance in order to define the vector white noise . The elements of the matrix exponential entering the construction are calculated with the Expokit tool (Sidje 1998) using the irreducible rational Padé approximant. The remaining convolutions are performed in the Fourier space.

We choose as an isotropic cut-off function the following , compactly supported function

 φL(x)=e−|x|2L2−|x|21|x|⩽L,

and we will consider the particular value in order to obtain a couple of integral scales in our simulations. The precise shape of this function is not important, besides its characteristic length scale , and only the large scale statistical quantities such as the variance depend on it. We will see that at small scales, explored as an example by velocity increments, only the value at the origin matters.

As a regularization mechanism, we use the following regularized norm

 |x|2ϵ=|x|2+ϵ2. (3.0)

This regularization procedure makes the continuous field (2.0) differentiable, and divergence free in particular. In the discrete approximation, we cannot choose arbitrarily small, since it is bounded from below by the finiteness of the smallest accessible scale . As mentioned, plays the role, in a schematic way, of the Kolmogorov scale, and therefore should depend on the Reynolds number. As far as the Gaussian field is concerned (equation 1), we can relate to the kinematic viscosity and the Hurst exponent in such a way that the average dissipation per unit of mass remains finite and strictly positive when (Chevillard 2015). In the following, we will work numerically at a finite viscosity, and we will be interested theoretically in the asymptotic limit , which corresponds to the infinite Reynolds number limit. We thus have to take greater than . When , then the numerical field is smooth and gradients are well approximated. In particular, in standard deviation, the divergence of the field is much smaller than the gradient of one of its components. When , the numerical field is rough, and gradients are poorly approximated. In other words, the divergence of the field can be of the order of the gradient of one of its components (in standard deviation). We are also interested in simulations where the inertial range is wide, i.e. we would like to maximize the ratio . In section 4, since we will focus on velocity gradients, we will use . In the following sections, we will use . Once again, the precise regularization procedure is not important as long as is of order at the origin, and equal to at a distance from the origin. We can rigorously show that this is the case for the matrix multiplicative chaos (Chevillard et al. 2013).

## 4 Statistical structure of the velocity gradient tensor

As an important characterization of the velocity gradient tensor , we study its two non vanishing invariants. For instance, the second invariant is given by

 Q=−12tr(\mathsfbiA2)=14|ω|2−12tr(\mathsfbiS2) (4.0)

where is the vorticity vector and the symmetric part of , and can be interpreted as the competition between enstrophy and dissipation (per unit viscosity). Then, positive represents rotation-dominated regions and negative dissipation-dominated regions. Analogously, the third invariant is given by

 R=−13tr(\mathsfbiA3)=−14ωiSijωj−13tr(\mathsfbiS3) (4.0)

representing competition between enstrophy production and dissipation production. See Tsinober (2001); Wallace (2009); Meneveau (2011) for a discussion on this topic. We simulate the vector field for four different values of intermittency coefficients , with , and (see discussion in section 3) and represent the numerical estimation of the joint density of the invariants and in figure 1.

As is well known, a Gaussian velocity field corresponding to (1), or equivalently the velocity field (1.0) with , predicts a joint density of the invariants symmetrical with respect to the line. This is what we obtain in figure 1(a). In figures 1(b,c,d), we study the effect of increasing . We can see that the bigger the value of , the more elongated is the joint density along the right tail of the zero discriminant (or Vieillefosse) line, where (Vieillefosse 1982). We will see that increasing corresponds to increasing the level of intermittency. As justified in section 5, we choose the very particular value for turbulence applications whose corresponding joint density of and is displayed in figure 1(c).

Another striking property of turbulence is the preferential alignment of vorticity with the strain eigendirection associated to the intermediate eigenvalue. We refer again to Tsinober (2001); Wallace (2009); Meneveau (2011) for further discussions. We represent in figure 2 the probability density of the cosine of the angle between vorticity and the eigenvectors of the strain. Figure 2(b) indicates the preferential alignment of vorticity with the correct eigenvector, as observed already in Chevillard et al. (2010) for a single value of . Here, we can see that this alignment is governed by the intermittency coefficient : no preferential alignment is observed when , as expected from a Gaussian velocity field, and this preferential alignment increases with increasing . We observe also in figure 2(a) that the density of the preferential orthogonality of vorticity with the eigendirection associated to the smallest (negative) eigenvalue is barely sensitive to , except in the Gaussian case . As for the angle between vorticity and the eigendirection associated to the biggest (positive) eigenvalue (figure 2(c)), as observed in real flows, the density is almost flat, with a slight dependence on the parameter , showing no preferential orientation.

Overall, as far as velocity gradient statistics are concerned, the velocity field predicts non trivial facts of fluid turbulence. In this picture, at least for the range of values studied, the dependence on this parameter is weak. We will study in the following section the influence of the parameter on the scaling of structure functions, where it will play a key role.

## 5 Numerical and theoretical study of the intermittent properties

As we have seen in section 4, at a finite (or finite Reynolds number), the velocity field (1.0) predicts realistic velocity gradient statistics. In particular, at any , we reproduce the teardrop shape of the -plane (figure 1) and the preferential alignment of vorticity (figure 2). The precise value of representing realistic turbulent statistics has not been selected yet. This is the purpose of this section.

We focus now on the intermittency phenomenon, and explore the influence of the parameter on the anomalous scaling of the structure functions (Frisch 1995). To do so, we will discuss this phenomenon based on the flatness of velocity increments (defined in equation 1), namely

 F∥,⊥ϵ(ℓ)=E(δ∥,⊥ℓuϵ)4[E(δ∥,⊥ℓuϵ)2]2, (5.0)

in both the longitudinal (i.e. ) and transverse cases (). We perform first simulations of the field for , and several . We choose in order to maximize the extend of the inertial range, and represent the results of the estimation of the flatness in figure 3.

### 5.1 Numerical estimations

We display in figure 3(a) (respectively figure 3c) the flatness of the longitudinal (transverse) velocity increments as a function of the scale for selected values of the parameter , including the Gaussian case . We indeed notice that in the Gaussian case the flatnesses and do not depend on scale and equal 3. Then, as increases, for any fixed scale, the flatness increases. We observe a power law with the scale within a limited range, i.e.

 F∥,⊥ϵ(ℓ)∼ϵ≪ℓ

the regularization over the scale polluting a large part of the accessible scales. We represent in figure 4 the obtained value for as a function of . We indeed see that the scaling exponents decrease when increases, showing the augmentation on the level of intermittency Frisch (1995). We observe also, in this context, that the level of intermittency of transverse velocity increments is higher than the one observed on longitudinal velocity increments. In laboratory and numerical flows (see for instance Chevillard et al. 2012), we find a universal behaviour (independent of the flow geometry and the Reynolds number) for the longitudinal case, with . Thus, turbulence statistics seem to be well reproduced for a very particular small value of the parameter , which was already found in Chevillard et al. (2010).

The present velocity field (1.0) is an example of random process that exhibits a higher level of intermittency for the transverse case than for the longitudinal case. This is indeed a surprising effect, also observed in real flows (see discussions in Dhruva et al. 1997; Chen et al. 1997; Grauer et al. 2012). In our model, only an analytical study could give a clear answer to this observed discrepancy, in particular in the asymptotic limit . Unfortunately, the underlying mathematical structure of this field is subtle, the strong correlation between the field and the white measure is difficult to handle. In the following, we will study both numerically and theoretically an ersatz of the velocity field (1.0) that follows the same rules of construction, except that the field of matrices is built independently of the underlying white measure . This ersatz is amenable to exact derivations of its statistical properties.

### 5.2 The hypothesis of independence

Consider now the following velocity field

 uind,ϵi(x)=1cϵ∫ϕϵik(x−z)(eγ\mathsfbiXϵ(z))klWl(dz), (5.0)

where again

 ϕϵik(x)=−ϵijkφL(x)xj|x|52−Hϵ,

and the following matrix field:

 \mathsfbiXϵ(x)=√1532π∫|x−y|⩽Lx−y|x−y|7/2ϵ⊗[(x−y)∧W′(dy)]+[(x−y)∧W′(dy)]⊗x−y|x−y|7/2ϵ (5.0)

where now the vector noise is independent of the vector noise of the underlying Gaussian structure (5.0), namely for any , and any components and , we have . As we show in appendix C, the field (5.0) needs to be renormalized in order to converge, when , towards a finite-variance process. This deterministic normalization constant (which diverges when ) is given by

 c2ϵ=13E[tr e2γ\mathsfbiXϵ]. (5.0)

This is a standard way to renormalize the multiplicative chaos and this situation is well understood as far as multiplicative chaos is concerned (Chevillard et al. 2013; Rhodes & Vargas 2014). We can see that the difference between the full vector field (1.0) and the one just mentioned (5.0) is that the matrix free field and the vector white noise are independent. This is a strong simplification but we can get simple exact results.

#### 5.2.1 Numerical simulations

We perform similar simulations of the velocity field (5.0) as we did for the (dependent) velocity field (2.0) and compare the estimation of velocity increment flatnesses. We display our results in figures 3(b,d) and 4(b).

Overall, at any value of the parameter , both velocity fields and share qualitatively similar levels of intermittency. This gives confidence in explaining the intermittent properties of (2.0) using results from the intermittent nature of the velocity field ersatz (5.0), which is amenable to analytical derivation. The main difference, as far as intermittency is concerned, comes from the behaviour of the flatness of transverse velocity increments. We observe indeed that in the independent case, the observed power-laws of the flatnesses and (equation 5.0) have the same scaling exponent (5.0), a property that we show in the following sections. Also, as we will see, the velocity field ersatz does not exhibit energy transfer. In other words, is a good ersatz to study the intermittent nature of the velocity field , but fails at describing the physics of energy transfer. This reveals also the great importance of building up the matrix field from the very same white measure in order to predict energy transfer. We will come back to this point in section 6.

#### 5.2.2 Mean and covariance

Since and are independent, we easily get (see annex C)

 E[uind,ϵi(x)]=0.

Thus, the vector field (5.0) has zero average, as expected from an isotropic vector field. As for the covariance, we find

 E[uind,ϵi(0)uind,ϵp(h)]=(ϕϵik⋆ϕϵpk)(h),

where is the correlation product defined in equation (B 0). We can see that the covariance structure of the vector field (5.0) is the same as the one obtained from the underlying Gaussian field (1). In particular, the field (5.0) converges in a sense when , and it has same variance and covariance as the underlying Gaussian field (see annex B for the properties of the covariance of the underlying Gaussian field). We will note the corresponding limiting process as . Similarly, the Gaussian field and the field share the same second order structure functions (both longitudinal and transverse).

#### 5.2.3 Fourth order structure function and Flatnesses

We define the velocity increment as

 δℓuind,ϵi=uind,ϵi(ℓ/2)−uind,ϵi(−ℓ/2)=1cϵ∫Φϵ,ℓik(z)(eγ\mathsfbiXϵ(z))klWl(dz), (5.0)

where we have defined the even function

 Φϵ,ℓik(x)=ϕϵik(x+ℓ/2)−ϕϵik(x−ℓ/2).

Of special interest is the fourth order structure function (longitudinal and transverse cases), and the respective flatnesses, i.e. and (5.0). We get (no summation over repeated index implied),

 E(δℓ uϵ,indi)4= 3c4ϵ∫Φϵ,ℓik1(z2)Φϵ,ℓik2(z2)Φϵ,ℓik3(z4)Φϵ,ℓik4(z4)E[(e2γ\mathsfbiXϵ(z2))k1k2(e2γ\mathsfbiXϵ(z4))k3k4]dz2dz4, (5.0)

entering therefore the covariance of the matrix multiplicative chaos , which is analytically derived in Chevillard et al. (2013). In annex C, we show that the fourth-order structure function converges when if we choose , and behaves asymptotically in the limit of vanishing scale as

where the multiplicative constant is derived in annex C. This shows that longitudinal and transverse fourth order structure functions have the same scaling behaviours. More precisely, we obtain

 E(δ∥ℓuind)4∼ℓ→0Cind,∥4ℓ4H−4γ2ln1ℓ and E(δ⊥ℓuind)4∼ℓ→0Cind,⊥4ℓ4H−4γ2ln1ℓ

where the different constants and are also given in annex C. This shows that the velocity field (5.0), built assuming independence of and , is intermittent, and the respective flatnesses (5.0) behave as power-laws times a logarithmic correction with the scale (see annex C):

 F∥(ℓ)∼ℓ→0Cind,∥4(Cind,∥2)2ℓ−4γ2ln1ℓ and F⊥(ℓ)∼ℓ→0Cind,⊥4(Cind,⊥2)2ℓ−4γ2ln1ℓ. (5.0)

Thus, according to this asymptotic prediction (5.0), the power-law exponents of the flatness are the same and related to the parameter , namely .

As we can observe in figure 4, this asymptotic prediction performs poorly against our numerical results. We will see in section 5.2.5 that this quantitative discrepancy can be explained while taking into account finite scale corrections, as it is necessarily done while fitting power-laws of flatnesses estimated in numerical simulations.

#### 5.2.4 Heuristics for higher order structure functions

It is easy to see, and shown in annex C, that all odd order structure functions vanish under the hypothesis of independence: the assumption of independence prevents us from studying the physics of energy transfer. Let us study now the even higher-order structure functions and let us consider the -order moment of velocity increments , for , given in equation (C 0). Doing so, we are left considering the -correlator of the matrix exponential of :

 E[n∏q=1(e2γ\mathsfbiXϵ(z2q))k2q−1k2q].

For calculations are tedious, but, based on a conjecture of Chevillard et al. (2013) and an appropriate range of orders , we can write that the structure functions behave as

 limℓ→0lnE(δ∥ℓuind,∥,⊥)2qlnℓ=ζind2q,

with a similar spectrum of exponents for both longitudinal and transverse structure functions given by a quadratic function of the order , namely

 ζindq=qH−q(q−2)2γ2,

showing that indeed, the parameter fully determines the intermittent properties of the velocity field.

#### 5.2.5 Finite size corrections and interpretation of numerics

As far as flatnesses are concerned, we have found in our numerical simulations, for which results are displayed in figures 3(b,d), a power-law exponent , defined in (5.0), much smaller (in magnitude) than our prediction (5.0). This is shown in figure 4(b). In this section we propose to explain this surprising fact while taking into account finite-size corrections based on our predictions before looking at the asymptotic limit .

We have seen while deriving the flatnesses of velocity increments (section 5.2.3 and annex C) that the covariance of the matrix chaos enters the expression of the fourth-order structure function, as shown in (5.2.3). The mathematical theory developed in Chevillard et al. (2013) allows us not only to derive its asymptotic structure in the double limits and then as we have already seen, but also in the limit at a finite scale . Recall that in the double limits, we have found a remaining logarithmic correction to the power-laws of the flatnesses (see equation 5.0). The purpose of this section is to explore the behaviour of the matrix chaos covariance when the scale is finite, after taking the limit , that eventually leads to intermittency with logarithmic corrections at vanishing scales .

Clearly, the dependence on the scale of the flatnesses is linked to this matrix chaos covariance. Under the hypothesis of isotropy, we can show that the matrix chaos covariance depends on only two scalar functions and defined as

We give the expressions of and in (C 0) and (C 0), that are intermediate steps before taking the limit . As we show in annex C, the quantities and are responsible for the intermittent correction of the fourth order structure functions, the remaining kernels entering the expressions of the full structure functions (5.2.3) participate mainly to the scaling . When renormalizing by the square of the second-order structure function, defining thus the respective flatnesses, we focus only on the intermittent corrections. Thus, we will approximate the flatness exponent by the logarithmic derivative of the contribution associated to and , to write

 ∂lnF∥,⊥(ℓ)∂lnℓ=β∥,⊥(γ)≈∂ln[3f(ℓ)+6g(ℓ)]∂lnℓ. (5.0)

As we have seen also, a key quantity that enters the expression of and (equations C 0 and C 0) is the covariance of the diagonal elements of , i.e. (equation C 0), and we will write it as

 σ2ℓ=ln(Lℓ)+α, (5.0)

where is a constant independent of the scale (it is more generally a bounded function of the scale, but we will neglect this functional dependence). Obviously, the constant is negligible in front of the logarithm (5.0) when . It is not the case when is finite. We have estimated this constant in our numerical simulation and we find (data not shown). Using equation (5.0), we evaluate the logarithmic derivative at the scale , using the form of the covariance given in (5.0) with thus . We display the result of this fit, as a function of , in figure 4(b). The comparison with numerical data is fairly good at low values of and deteriorates at higher values. Several remarks can be made at this stage to justify the level of adequacy of our fit with numerical data. First, the model used to fit our data has been obtained in the limit of vanishing resolutions , whereas it can remain some finite- corrections when looking at a numerical simulation. Secondly, the parameter rigorously should be seen as a bounded function of the scales . We took it as a constant for the sake of simplicity. Thirdly, relation (5.0) is only an approximation, and there could be additional finite scale corrections related to the underlying Gaussian velocity field. Recall indeed, as shown in appendix B, that the Gaussian velocity field exhibits exact power laws in the double limits and , whereas we focus here on finite scale corrections. On the theoretical side, these corrections are difficult to obtain and will depend on the precise shape of the large-scale cut-off function entering the definition of the velocity field (5.0). Finally, let us note that the influence of the free parameter entering equation (5.0) is only quantitatively substantial for higher values of . We explain in this way the surprising fact that the estimated values of based on our simulations differ from the asymptotic prediction by taking into account finite-size corrections.

Under the independence assumption, we are thus able to quantify finite-size corrections to the scalings. Going back to the full vector field (1.0), we observe in figure 4(a) that, likewise, longitudinal and transverse structure functions seem to be affected by finite size corrections. We furthermore observe a slight difference between the longitudinal and transverse cases: according to our numerical simulations, it seems that transverse intermittency corrections are bigger than the longitudinal ones. The underlying strong correlation between the chaos and the white noise prevents us from deriving analytically the asymptotic regime and thus, we cannot conclude at this stage whether this slight difference will remain at vanishing resolution and vanishing scale.

## 6 Energy transfer: skewness phenomenon

Let us now turn back to the full (i.e. correlated) vector field (1.0) that we recall here for convenience,

 uϵi(x)=∫ϕϵik(x−z)(eγ\mathsfbiXϵ(z))klWl(dz), (6.0)

where

 ϕϵik(x)=−ϵijkφL(x)xj|x|52−Hϵ,

and the following field (2.0) is built from the very same vector white noise that enters the definition of the velocity field (6.0):

 \mathsfbiXϵ(x)=√1532π∫|x−y|⩽Lx−y|x−y|7/2ϵ⊗[(x−y)∧W(dy)]+[(x−y)∧W(dy)]⊗x−y|x−y|7/2ϵ.

The same vector noise enters both the velocity field equation (6.0) and the matrix field implying peculiar correlations, fully given by the correlators , that read, for ,

 Γpϵ(x−y) =E[\mathsfbiXϵ(x)Wp(dy)dy] =√1532π⎡⎣x−y|x−y|7/2ϵ⊗[(x−y)∧ep]+[(x−y)∧ep]⊗x−y|x−y|7/2ϵ⎤⎦ , (6.0)

where we have defined the unit vector with . Let us also remark that at a given finite ,

 Γpϵ(0)=0.

The components of the matrices are noted .

Obtaining the exact statistical properties of the velocity field (6.0) is a difficult task. The very peculiar correlation between and , fully encoded in the correlator (6), and the non-commutative nature of the field of matrices make the calculation out of reach at the present time. As an example, we do not know today how to perform such a calculation even for the variance of the field . Instead, in order to interpret the numerical evidences of energy transfer observed in Chevillard et al. (2010), we propose to do a simpler calculation, namely a perturbative expansion in power of that we hope will capture several key ingredients of the physics of energy transfer. Making such an expansion prevents the analysis of intermittent corrections since, in nature, the intermittency phenomenon cannot be treated, as far as we know, with such an expansion.

As in the independent case (5.0), we expect for the full vector field (6.0) a normalization constant such that is of finite variance. Recall that the elements of are Gaussian random variables whose variance diverges logarithmically with (cf. equation 2.0). Indeed, for the independent case, we have shown that the velocity field has to be normalized by a constant that itself diverges with . It has been shown in Chevillard et al. (2013) that the multiplicative chaos has to be renormalized in order to define a proper random variable (see also Rhodes & Vargas 2014 for a review on this topic). As far as (6.0) is concerned, since we do not know how to get the variance, we cannot make such a normalization constant explicit, but is expected to be of the order of (as it has been proved for the chaos in Chevillard et al. 2013) or (as we have shown for the independent case). In the following perturbative expansion, we expect then a contribution of order (and all the following even powers of ) from this constant. Since we do not know it explicitly, we will limit ourselves to a first order expansion in , which will not have any contribution from this possible unknown renormalizing constant. We will see then that such a first order expansion exhibit energy transfer.

### 6.1 First order expansion of the covariance

At a given finite , the covariance of the vector field is given by

 E[uϵi(0)uϵp(h)]=∫ϕϵik(−z1)ϕϵpq(h−z2)E[(eγ\mathsfbiXϵ(z1))kl(eγ\mathsfbiXϵ(z2))qrWl(dz1)Wr(dz2)].

Expanding the matrix exponentials up to first order gives

 E[(e