# Point-vortex statistical mechanics applied to turbulence without vortex stretching

Tong Wu, Tomos David, Wouter J T Bos

► To cite this version:

Tong Wu, Tomos David, Wouter J T Bos. Point-vortex statistical mechanics applied to turbulence without vortex stretching. Journal of Statistical Mechanics: Theory and Experiment, 2023, 2023, pp.113203. 10.1088/1742-5468/ad063a. hal-04310023

HAL Id: hal-04310023

https://hal.science/hal-04310023v1

Submitted on 27 Nov 2023

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
Point-vortex statistical mechanics applied to turbulence without vortex stretching

Tong Wu,¹ Tomos David,² and Wouter J. T. Bos¹

1) LMFA-CNRS, Ecole Centrale de Lyon, Université de Lyon, Ecully,
France

2) Univ Lyon, ENS de Lyon, Univ Claude Bernard, CNRS,
Laboratoire de Physique (UMR CNRS 5672), F-69342 Lyon,
France

(Dated: 14 September 2023)

**Abstract:** In turbulent systems with inverse cascades, energy will pile up at large scales if no large-scale sink is present. We observe that in forced-dissipative three-dimensional turbulence from which vortex-stretching is removed, such condensation is observed, associated with an inverse cascade of helicity. The large-scale structure of this condensate is characterized by a hyperbolic sine relation between vorticity and velocity, analogous to the sinh relation between vorticity and stream function observed in freely decaying 2D turbulence in periodic domains. We generalize a 2D point-vortex statistical mechanics approach to our 3D system. It is shown that the predictions of this approach are in agreement with observations of both the forced-dissipative system, after appropriate averaging, and of the freely decaying system.

**Keywords:** Vortex-stretching, condensation, statistical mechanics, turbulence
# I. INTRODUCTION

In two-dimensional (2D) turbulence at high Reynolds numbers, nonlinear interactions yield a transfer of energy towards large scales¹. If energy is continuously injected, and in the absence of a large-scale dissipation mechanism, energy eventually accumulates in low wave numbers, leading to a condensate¹. This condensation phenomenon has been observed numerically and, in 2D periodic square domains, the physical manifestation of this condensation is the generation of a large-scale counter-rotating vortex pair²,³. Such large-scale coherent structures in turbulent flows are at the heart of a number of geophysical processes⁴⁻⁶.

Large-scale structuring in 2D and quasi-2D systems has been extensively studied using statistical mechanics⁷⁻⁹. Two different statistical mechanics approaches are the Lee-Kraichnan equilibrium statistical mechanics¹,¹⁰ and the Robert-Sommeria-Miller (RSM) approach¹¹⁻¹³, respectively. Kraichnan predicted the shape of energy spectra associated with the absolute equilibrium state of a truncated Fourier system taking into account the two inviscid invariants surviving Galerkin truncation: energy and enstrophy¹. However, a weakness of Kraichnan’s theory is that it cannot predict structuring and does not take into account the higher-order moments of the vorticity. For instance, an observation that is not captured by such an approach is the following: In freely evolving 2D turbulence in square periodic domains, at long times, vorticity $\omega$ and stream function $\Psi$ are locally related by a hyperbolic sine function, corresponding to a steady solution of the 2D Euler-equations. More precisely, observations in numerical simulations are consistent with a functional relation

$$ \omega = -c^{-1} \sinh(b\Psi), \qquad (1) $$

where c and b are constants ($c > 0, b < 0$)¹⁴. This relation was predicted by applying statistical mechanics to a point-vortex system. Indeed, building upon the ideas of Onsager¹⁵, Joyce and Montgomery obtained the above quantitative prediction¹⁶.

Navier-Stokes turbulence is a dissipative system, whereas the point-vortex dynamics constitute an ideal finite-dimensional Hamiltonian system. The agreement of theory with observations of freely evolving Navier-Stokes turbulence is, therefore, not a trivial fact. However, point-vortex descriptions have allowed important progress in the understanding of two-dimensional turbulence⁸,¹⁷,¹⁸.
It is even less straightforward that these approaches should work when a turbulent system is continuously stirred by a body force. In such forced-dissipative systems, inviscid invariants are continuously modified so that the application of equilibrium statistical mechanics is not well founded in principle. Nevertheless, on average such approaches might work. Indeed, recently it was observed that maximum entropy theory could be applied to forced-dissipative geostrophic turbulence<sup>19,20</sup>. Another indication that statistical mechanics might be applied to averaged turbulent systems is the success of the application of statistical mechanics of axisymmetric turbulence<sup>21-23</sup> to observations of experiments of turbulent flows, which are only axisymmetric on average<sup>24</sup>.

In the present work, we continue this line of research and show how point-vortex statistical mechanics can be applied to forced-dissipative and freely decaying three-dimensional (3D) turbulence scenarios, however, without vortex stretching. We have recently established a framework for incompressible turbulence without vortex stretching<sup>25,26</sup>. In cases where no forcing or dissipation is present, the governing dynamics can be described by:

$$ \frac{\partial \omega}{\partial t} + \mathbf{u} \cdot \nabla \omega = -\nabla P_{\omega}, \qquad (2) $$

where $\mathbf{u}$ represents the velocity with $\nabla \cdot \mathbf{u} = 0$, $\omega = \nabla \times \mathbf{u}$ denotes the vorticity, and $P_{\omega}$ is a pressure term ensuring the vorticity to remain solenoidal. Indeed, Eq. (2) is obtained by removing the vortex-stretching term, $\omega \cdot \nabla \mathbf{u}$, from the curl of the 3D Euler equation. In previous work<sup>27</sup>, we illustrated that the forced turbulent system without vortex stretching exhibits an inverse helicity cascade (and also a forward enstrophy cascade) bearing similarities to the behavior observed in 2D turbulence. However, different from 2D turbulence, the no-vortex-stretching system does not conserve kinetic energy. This inverse transfer led us to consider the following questions: In the absence of dissipation at large scales, will this system manifest a condensation structure? Can we predict the properties of such a structure using point-vortex theory?

In this paper, we show that, in forced-dissipative turbulence without vortex stretching, condensation indeed occurs at the scale of the box. We will also investigate the system when it is freely decaying from this condensation state. Moreover, we attempt to predict the characteristics of large-scale structures in both the forced-dissipative and freely decaying cases using the principle of point-vortex statistical mechanics.

The rest of this paper is organized as follows. In Sec. II, we predict properties of
large-scale structures in 3D turbulence without vortex stretching by using a generalized
point-vortex model. In Sec. III, we present the numerical setup. Then, in Sec. IV we report
on the assessment of the theoretical results. Finally, Sec. V presents the conclusions.

II. ANALYTICAL CONSIDERATIONS

At large but finite Reynolds numbers, a 2D Navier-Stokes fluid will relax to a quasi-static maximum entropy state where most of the energy is concentrated at large scales. When this state is attained, a sinh relation between vorticity and stream function emerges (Eq. (1)). This was observed numerically<sup>14,28</sup> and explained analytically<sup>11,29</sup> by a point-vortex model. The original idea to use such a point-vortex model in this context came from Onsager. He considered a Hamiltonian system of *N* parallel point-vortices representing an incompressible and inviscid 2D fluid<sup>15</sup>. Then Joyce and Montgomery used this idea to study the 2D electrostatic guiding center plasma<sup>16,30</sup>. The system of interacting point-vortices is equivalent to the system which describes the interaction of long, uniformly charged rods aligned with a uniform magnetic field. Joyce and Montgomery considered *N* rods of charge +*e* and *N* more of charge −*e* and predicted the emergence of the hyperbolic sine relation.

We investigate whether there is a similar functional relation in our turbulent system
without vortex stretching. We generalize the point-vortex model to a 3D system. Thereto
we propose that a 3D vorticity field can be considered as a combination of 6*N* point-vortices,
the axes of which are each parallel to one of the three coordinate axes. We suppose that in
each direction there are *N* point-vortices of positive vorticity +1, and *N* ones of negative
vorticity −1. Without loss of generality, we use non-dimensional vorticities +1 and −1 here.
Similar to Ref.<sup>16</sup>, we imagine the total volume *V* = *L*<sup>3</sup> subdivided into small cells of volume
Δ ≪ *V*. The cells are large enough, however, to contain many particles. We call *N*<sub>*i*,*x*</sub><sup>+</sup>, *N*<sub>*i*,*y*</sub><sup>+</sup>,
*N*<sub>*i*,*z*</sub><sup>+</sup> and *N*<sub>*i*,*x*</sub><sup>−</sup>, *N*<sub>*i*,*y*</sub><sup>−</sup>, *N*<sub>*i*,*z*</sub><sup>−</sup> the number of positive and negative point-vortices in each direction
inside cell *i*. We can consider that one small cell constitutes a basic unit, which means
the resultant quantities (such as velocity, vorticity, etc.) of all point vortices in a cell can
be considered as the quantity at the center point of the cell. Furthermore, it might seem
artificial that all point vortices are aligned with the coordinate axes. However, since all
cells contain a large number of vortices, we can without loss of generality assume that three
vortices, in the three directions, represent the three components of one point vortex. This
would reduce the total number of vortices in the description by a factor of three, which will not make any difference in the asymptotic limit, where all $N_{i,a}^{\pm}$ ($a$ indicates $x$, $y$, or $z$) tend to infinity. Then the vorticity field can be written as

$$ \left\{ \begin{aligned} \omega_x(\mathbf{r}) &= \sum_i (N_{i,x}^+ - N_{i,x}^-) \delta(\mathbf{r} - \mathbf{r}_i), \\ \omega_y(\mathbf{r}) &= \sum_i (N_{i,y}^+ - N_{i,y}^-) \delta(\mathbf{r} - \mathbf{r}_i), \\ \omega_z(\mathbf{r}) &= \sum_i (N_{i,z}^+ - N_{i,z}^-) \delta(\mathbf{r} - \mathbf{r}_i), \end{aligned} \right. \qquad (3) $$

where $\mathbf{r}_i$ indicates the position of cell $i$. The associated velocity field can be obtained using Biot-Savart's law

$$ \mathbf{u}(\mathbf{x}) = \frac{1}{4\pi} \int_{\mathbb{R}^3} \frac{\omega(\mathbf{y}) \times (\mathbf{x} - \mathbf{y})}{|\mathbf{x} - \mathbf{y}|^3} d\mathbf{y}. \qquad (4) $$

Substituting Eq. (3) into Eq. (4), yields expressions for the velocity components in each direction

$$ \left\{ \begin{aligned} u_x(\mathbf{r}) &= \frac{1}{4\pi} \sum_i \frac{(N_{i,y}^+ - N_{i,y}^-)(r_z - r_{i,z}) - (N_{i,z}^+ - N_{i,z}^-)(r_y - r_{i,y})}{|\mathbf{r} - \mathbf{r}_i|^3}, \\ u_y(\mathbf{r}) &= \frac{1}{4\pi} \sum_i \frac{(N_{i,z}^+ - N_{i,z}^-)(r_x - r_{i,x}) - (N_{i,x}^+ - N_{i,x}^-)(r_z - r_{i,z})}{|\mathbf{r} - \mathbf{r}_i|^3}, \\ u_z(\mathbf{r}) &= \frac{1}{4\pi} \sum_i \frac{(N_{i,x}^+ - N_{i,x}^-)(r_y - r_{i,y}) - (N_{i,y}^+ - N_{i,y}^-)(r_x - r_{i,x})}{|\mathbf{r} - \mathbf{r}_i|^3} \end{aligned} \right. \qquad (5) $$

with $r_{i,x}$, $r_{i,y}$, $r_{i,z}$ the $x$, $y$, and $z$-component of $\mathbf{r}_i$ and $r_x$, $r_y$, $r_z$ the $x$, $y$, and $z$-component of $\mathbf{r}$.

This point-vortex model allows us to define an entropy using Boltzmann's formula. The entropy of this 3D point-vortex system can be written as

$$ S = \ln \mathcal{W} \qquad (6) $$

with

$$ \mathcal{W} = \frac{(6N)!}{\prod_i \{N_{i,x}^+! N_{i,x}^-! N_{i,y}^+! N_{i,y}^-! N_{i,z}^+! N_{i,z}^-!\}}. \qquad (7) $$

The quantity $\mathcal{W}$ indicates the number of microstates associated with the macrostate where there are $N_{i,x}^+$, $N_{i,y}^+$, $N_{i,z}^+$ positive and $N_{i,x}^-$, $N_{i,y}^-$, $N_{i,z}^-$ negative point-vortices inside cell $i$. We assume $N_{i,x}^+$, $N_{i,y}^+$, $N_{i,z}^+$, $N_{i,x}^-$, $N_{i,y}^-$, $N_{i,z}^-$ to be large enough for Stirling's formula to be valid. We then obtain

$$ S \approx 6N\ln(6N) - 6N + \sum_i (-N_{i,x}^+ \ln N_{i,x}^+ + N_{i,x}^+ - N_{i,x}^- \ln N_{i,x}^- + N_{i,x}^- - N_{i,y}^+ \ln N_{i,y}^+ + N_{i,y}^+ \\ - N_{i,y}^- \ln N_{i,y}^- + N_{i,y}^- - N_{i,z}^+ \ln N_{i,z}^+ + N_{i,z}^+ - N_{i,z}^- \ln N_{i,z}^- + N_{i,z}^-) $$
(8)

Similar to a 2D system, we assume the final state of 3D turbulence without vortex stretching to be close to the maximum entropy state. We will therefore attempt to determine this state by maximizing the entropy under the constraints imposed by the dynamical equation. Indeed, our system (2) conserves enstrophy and helicity²⁶. To proceed, we need therefore to consider formulas of these two invariants expressed as a function of the variables governing the 3D point-vortex system. The total helicity of the system is given by

$$ H = \frac{1}{2} \int \boldsymbol{u} \cdot \boldsymbol{\omega} d\boldsymbol{r} = \frac{1}{8\pi} \sum_{i,j} \frac{1}{|\boldsymbol{r}_j - \boldsymbol{r}_i|^3} \left( \begin{aligned} & (N_{j,x}^+ - N_{j,x}^-) (N_{i,y}^+ - N_{i,y}^-)(r_{j,z} - r_{i,z}) - (N_{i,z}^+ - N_{i,z}^-)(r_{j,y} - r_{i,y}) + \\ & (N_{j,y}^+ - N_{j,y}^-) (N_{i,z}^+ - N_{i,z}^-)(r_{j,x} - r_{i,x}) - (N_{i,x}^+ - N_{i,x}^-)(r_{j,z} - r_{i,z}) + \\ & (N_{j,z}^+ - N_{j,z}^-) (N_{i,x}^+ - N_{i,x}^-)(r_{j,y} - r_{i,y}) - (N_{i,y}^+ - N_{i,y}^-)(r_{j,x} - r_{i,x}) \end{aligned} \right). \quad (9) $$

And enstrophy of this 3D point-vortex system can be written as

$$ W = \sum_{i} (1)^2 N_{i,x}^{+} + \sum_{i} (-1)^2 N_{i,x}^{-} + \sum_{i} (1)^2 N_{i,y}^{+} + \sum_{i} (-1)^2 N_{i,y}^{-} + \sum_{i} (1)^2 N_{i,z}^{+} + \sum_{i} (-1)^2 N_{i,z}^{-}. \quad (10) $$

The maximum entropy state is now obtained by solving the variational problem,

$$ \delta S - \gamma \delta H - \alpha \delta W = 0 \quad (11) $$

with $\gamma$, $\alpha$ Lagrange multipliers. Similar to what was assumed for 2D systems, we suppose the numbers of positive or negative point-vortices to be independent of each other. Thus, the partial derivative of helicity with respect to $N_{j,x}^+$ is

$$ \frac{\partial H}{\partial N_{j,x}^{+}} = \frac{1}{8\pi} \sum_{i} \frac{(N_{i,y}^{+} - N_{i,y}^{-})(r_{j,z} - r_{i,z}) - (N_{i,z}^{+} - N_{i,z}^{-})(r_{j,y} - r_{i,y})}{|\boldsymbol{r}_j - \boldsymbol{r}_i|^3} = \frac{u_x(\boldsymbol{r}_j)}{2}. \quad (12) $$

Other derivatives have similar forms. Then, substituting Eq. (8), (9) and (10) into Eq. (11) and deriving the integrand with respect to $N_{i,x}^+$ and $N_{i,x}^-$ yields the expressions

$$ \begin{aligned} \ln N_{i,x}^{+} + \alpha + \frac{\gamma}{2} u_x(\boldsymbol{r}_i) &= 0, \\ \ln N_{i,x}^{-} + \alpha - \frac{\gamma}{2} u_x(\boldsymbol{r}_i) &= 0. \end{aligned} \quad (13) $$

Results in the other two directions are similar. Hence, the number of point-vortices in cell $i$ with vorticity vectors along the $x$-direction are

$$ \begin{aligned} N_{i,x}^{+} &= \exp\left(-\left(\alpha + \frac{\gamma}{2} u_x(\boldsymbol{r}_i)\right)\right), \\ N_{i,x}^{-} &= \exp\left(-\left(\alpha - \frac{\gamma}{2} u_x(\boldsymbol{r}_i)\right)\right), \end{aligned} \quad (14) $$
at the maximum entropy state. Using Eq. (3), we find

$$ 
\begin{aligned}
\omega_x &= \exp(-\alpha)(\exp(-\frac{\gamma}{2}u_x) - \exp(\frac{\gamma}{2}u_x)) \\
&= 2\exp(-\alpha)\sinh(-\frac{\gamma}{2}u_x).
\end{aligned} 
\qquad (15) 
$$

After applying similar derivations in the other two directions, the vorticity components in the y- and z-directions are

$$ 
\begin{aligned}
\omega_y &= 2 \exp(-\alpha) \sinh(-\frac{\gamma}{2}u_y), \\
\omega_z &= 2 \exp(-\alpha) \sinh(-\frac{\gamma}{2}u_z).
\end{aligned} 
\qquad (16) 
$$

This sinh relation between vorticity and velocity, Eqs. (15) and (16) are the principal theoretical results of the present investigation.

### III. NUMERICAL METHOD AND SETUP

Direct numerical simulations (DNSs) are performed using a standard pseudo-spectral solver with a third-order Adams-Bashforth time integration scheme³¹, which was modified to remove vortex stretching²⁶. Our computational domain is a cubic periodic box of size $L = 2\pi$. Aliasing errors are removed using the 2/3 rule. DNSs are executed on grids of size 256³.

The Biot-Savart operator is applied in Fourier space to the vorticity equation in order to obtain an equation for the velocity field. Then the equation for Fourier coefficients of velocity $\hat{u}$, which is solved in our code, is

$$ \frac{\partial \hat{u}}{\partial t} + \frac{i}{k^2} \mathbf{k} \times \mathcal{F}[(\mathbf{u} \cdot \nabla)\mathbf{\omega}] = \hat{f} - \nu k^{2\eta} \hat{u} \qquad (17) $$

with $\mathbf{k}$ the wave vector, $k = |\mathbf{k}|$ the wave number, $\mathcal{F}[\bullet]$ the Fourier transform, $\hat{f}$ a forcing term, $\nu$ the (hyper)viscosity, and $\eta$ the hyperviscosity parameter. The physical space equivalent of (17) (Eq. (2)) contains a pressure term which allows to impose incompressibility. Equivalently, in Fourier-space, the incompressibility $(\nabla \cdot \mathbf{u} = 0)$ is contained in $\frac{i}{k^2} \mathbf{k} \times \mathcal{F}[\bullet]$ which projects the nonlinearity of the evolution-equation onto a plane perpendicular to the wave-vector.

We consider two cases. The first one is a forced-dissipative system. For the second case, we eliminate the force and allow the system to decay freely. Since we are concentrating on the condensation state without focusing too much on the dissipation range, we reduce
the size of the latter by using hyperviscosity and we use therefore the value $\eta = 4$ (normal viscous dissipation corresponds to $\eta = 1$). The parameter $\nu$ represents the hyperviscosity rate which is set to $10^{-13}$ in our simulations. We choose an injection mechanism which keeps the energy constant in a narrow wavenumber range$^{32,33}$ $29.5 < k < 31.5$ around the forcing wavenumber $k_f = 30.5$ at a level $E(k) = 10^{-3}$. A convenient property of this type of forcing, with respect to the present investigation, is that it does not only inject energy in the system, but also helicity$^{27}$.

The initial energy spectrum is chosen as $E(k) = 5 \times 10^{-7}$. The initial value of the energy is small enough to allow a clear observation of the energy condensation process. Complex phases of $\hat{u}$ are set randomly at the initial moment.

# IV. NUMERICAL RESULTS

We first show the energy condensation process in the forced system in Sec. IV A. Subsequently, in Sec. IV B, we eliminate the force and allow the system to decay freely from the condensation state. In both cases, we verify our predictions from the point-vortex model by illustrating the emergence of a hyperbolic sine relation between velocity and vorticity.

## A. Self-organization in the forced system

We start by forcing the system from a random initial condition, dissipating mainly the small scales. In Fig. 1, we show the time evolution of the helicity flux normalized by its dissipation rate $\epsilon_H$, defined as $\epsilon_H = \int 2\nu k^{2\eta} H(k) dk$, during the condensate phase from $t = 5$ to $t = 160$. In turbulence governed by the reduced Navier-Stokes equation Eq. (17), the flux of helicity is defined as $\Pi_H(k) = - \int_{\Sigma_k} \text{Re} [\hat{N}(k) \cdot \hat{u}^*(k)] dk$ with $\hat{N}(k) = \mathcal{F}[-(u \cdot \nabla)\omega]$ the non-linear term in the vorticity equation. Re[·] represents the real part of the quantity in brackets. Asterisk (*) denotes the complex conjugate. $\Sigma_k$ is the spherical domain in Fourier space consisting of all wave vectors with $||k|| \le k$. The normalized flux is negative at scales smaller than the forcing wavenumber $k_f = 30.5$, implying the helicity cascades from small towards large scales during this phase.

In Fig. 2 (a) and (b), we show the time evolution of energy and helicity spectra, respectively. A peak in the energy and helicity spectra appears at large scales after $t = 60$. Indeed,
<figure><img src="image_2.png" /><figcaption>Fig. 1. Flux of helicity normalized by the dissipation rate of helicity during the force-dissipative phase at $t = 5, 10, 20, 30, 40, 50, 60, 70, 80, 100, 130, 160$. The coloring of the curves evolves in time from dark to light.</figcaption></figure>

<figure><img src="image_3.png" /><figcaption>Fig. 2. Time evolution of wavenumber spectra of the kinetic energy (a) and helicity (b) during the force-dissipative phase. The time-instants are the same as indicated in Fig. 1.</figcaption></figure>

since helicity is conserved by the nonlinearity and there is only very weak damping present at scales $k < k_f$, the physical process leading to the build-up at the large scales of the system is associated with the helicity which is transferred to large scales. After arriving at the largest scales of the system, no dissipation mechanism is able to absorb the helicity injected by the forcing. The helicity and its associated energy piles then up, leading to large spectral
<figure><img src="image_4.png" /><figcaption>Fig. 3. Time evolution of integral quantities in three directions during the forced phase. $\langle\bullet\rangle$ indicates the volume average. (a) Energy. (b) Helicity. (c) Enstrophy.</figcaption></figure>

peaks at the small wavenumbers. The spectra of kinetic energy and helicity for intermediate wavenumbers during the inverse cascade of helicity are theoretically proportional to $k^{-7/3}$ and $k^{-4/3}$, as expected for the present system$^{27}$. In unmodified 3D turbulence, the helicity is transferred in the opposite direction (to small scales). However, dimensional arguments lead in that case also to the possibility of a $-7/3$ scaling$^{34}$, but this is only observed if the forcing scheme is tuned such that the flow becomes extremely helical$^{35}$. Due to the limitation of computing resources, the inertial range is not sufficiently large to clearly display these theoretical scalings.
<figure><img src="image_5.png" /><figcaption>Fig. 4. Visualizations of condensate-like structure at $t = 160$. (a) Energy. (b) Helicity.</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>Fig. 5. Scatter plot of y-components velocities and vorticities at the condensation state at $t = 160$.</figcaption></figure>

In Fig. 3, we show the time evolution of kinetic energy, helicity, and enstrophy in the three directions, respectively. After $t \approx 60$, kinetic energy, helicity, and enstrophy in the y direction become larger than in the other two directions. At this time the condensation process starts. The associated condensate-like structure is visualized in Fig. 4. This structure is constituted of highly anisotropic columnar structures in the y direction. A parallel can here be drawn with 2D turbulence where the condensate-like structure is a pair of counter-rotating vortices.

After reaching this condensate state, we assess the relation between velocities and vor-
<figure><img src="image_7.png" /><figcaption>Fig. 6. Velocity components in three directions at the condensation state at $t = 160$. (a) $u_x$. (b) $u_y$. (c) $u_z$.</figcaption></figure>

tivities as predicted by the point-vortex model in Sec. II. Fig. 5 shows the scatter plot of $u_y$ against $\omega_y$ for the condensation state at $t = 160$. A tendency of a sinh-relation can be observed, but a considerable noise hinders us from identifying a clear functional relation-ship. In order to reduce the noise, instead of investigating the relation between velocity and
vorticity components directly, we focus on their averaged values

$$ \begin{aligned} \langle u_y \rangle_y (x, z) &= \frac{1}{M_y} \sum_y u_y(x, y, z), \\ \langle \omega_y \rangle_y (x, z) &= \frac{1}{M_y} \sum_y \omega_y(x, y, z), \end{aligned} \qquad (18) $$

where $u_y$ and $\omega_y$ indicate velocity and vorticity in the y-direction, $\langle \cdot \rangle_y$ indicates a space average in the y-direction and $M_y$ is the number of grid points in the y-direction ($M_y = 256$ in our simulations). Similar definitions can be defined in the other directions. As sinh is not a linear function, theoretically, we can not derive $\langle \omega_y \rangle_y = 2 \exp(-\alpha) \sinh(-\frac{\gamma}{2} \langle u_y \rangle_y)$ from Eq. (15) and (16). But if $u_y$ is nearly constant along the y-direction, we can get an approximation $\langle \sinh(-\frac{\gamma}{2} u_y) \rangle_y \approx \sinh(-\frac{\gamma}{2} \langle u_y \rangle_y)$. An instantaneous visualization of three velocity components for the condensation state at $t = 160$ is shown in Fig. 6. We observe that $u_x$, $u_y$ and $u_z$ depend only weakly on the x, y and z directions respectively. Thus, the averaged values are appropriate for investigating the hyperbolic sine relation. Furthermore, the structure of $u_y$ is similar to a dipole, while structures of $u_x$ and $u_z$ are closer to a unidirectional flow. Note that dipoles and unidirectional flows are two condensation states of 2D turbulence caused by bifurcations of the stochastic Navier-Stokes equations$^6$.

The values of the two Lagrange multipliers $\alpha$ and $\gamma$ can be determined a priori. Firstly, we introduce the averaged quantities associated with helicity and entrophy, as

$$ \begin{aligned} H_{ave} &= \frac{1}{2} \sum_{y,z} \langle u_x \rangle_x (y,z) \langle w_x \rangle_x (y,z) + \frac{1}{2} \sum_{x,z} \langle u_y \rangle_y (x,z) \langle w_y \rangle_y (x,z) + \\ & \frac{1}{2} \sum_{x,y} \langle u_z \rangle_z (x,y) \langle w_z \rangle_z (x,y), \end{aligned} \qquad (19) $$

$$ W_{ave} = \frac{1}{2} \sum_{y,z} \langle w_x \rangle_x^2 (y,z) + \frac{1}{2} \sum_{x,z} \langle w_y \rangle_y^2 (x,z) + \frac{1}{2} \sum_{x,y} \langle w_z \rangle_z^2 (x,y) $$

respectively. We substitute the relation $\langle \omega_m \rangle_m = 2 \exp(-\alpha) \sinh(-\frac{\gamma}{2} \langle u_m \rangle_m)$ into Eq. (19) where $m$ indicates $x$, $y$ and $z$, so that $H_{ave}$ and $W_{ave}$ can be written as functions of averaged velocities $\langle u_m \rangle_m$. At each moment, exact values of $\langle u_m \rangle_m$, $H_{ave}$, and $W_{ave}$ can be obtained from the simulation data. Then the two unknowns $\alpha$ and $\gamma$ can be calculated from these two equations of $H_{ave}$ and $W_{ave}$. At $t = 160$, the analytically predicted curve is $\langle \omega_m \rangle_m = 1.5207 \sinh(0.3779 \langle u_m \rangle_m)$ in our simulation.

In Fig. 7, we show scatter plots of averaged components of velocities and vorticities in three directions for the condensation state at $t = 160$. A hyperbolic sine relation is
<figure><img src="image_8.png" /><figcaption>(a)</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>(b)</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>(c)</figcaption></figure>

Fig. 7. Scatter plots of averaged components of velocities and vorticities in three directions at the condensation state at $t = 160$. The yellow lines present the analytically predicted curve $\langle\omega_m\rangle_m = 1.5207 \sinh(0.3779 \langle u_m\rangle_m)$. (a) X-components. (b) Y-components. (c) Z-components.

observed between y-components $\langle u_y\rangle_y$ and $\langle \omega_y\rangle_y$. The yellow lines in Fig. 7 indicate the analytically predicted curve, which overlaps the data points in the direction where most energy is contained. This collapse verifies our analytical prediction of the sinh relation in the y-direction. However, for other directions, we observe a linear instead of a sinh relation between velocities and vorticities, as shown in Fig. 7(a) and (c). **The Taylor expansion of sinhUNITARY OPERATOR** for small $x$ yields $cx$ with $c$ a constant number, which indicates that the form of sinh($cx$) resembles a straight line when $x$ varies in a small range around 0. We insist that
these functional relations are not direct fits of the sinh relation to the data but obtained
evaluating the averaged quantities $\langle u_m \rangle_m$, $H_{ave}$ and $W_{ave}$ only.

Moreover, we note that the $y$ direction is not necessarily the direction which contains most energy. Further details can be found in Appendix A, where we present results from additional simulations.

## B. Freely decaying state

At the condensation state at $t = 160$, the external force is removed from the system, and we indicate in this section the moment when the removal of force commences as time $t = 0$.

In Fig. 8, we show the time evolution of kinetic energy, helicity, and enstrophy in the three directions, respectively. During the free-decay phase, most of kinetic energy, helicity, and enstrophy remain concentrated in the $y$ direction. As shown in Fig. 8(a), energy displays an increasing trend during the dissipation phase, especially in the $y$ direction, which is not violating any conservation laws because energy is not conserved in turbulence without vortex stretching$^{26}$. From Fig. 8(b), we can see that the value of helicity exhibits very little change over time since most of the helicity cascades towards small wavenumbers and is only weakly dissipated by the hyperviscosity. Fig. 8(c) shows the decrease of enstrophy.

In Fig. 9(a) and (b), we show the time evolution of energy and helicity spectra, respectively. We observe that the energy and helicity at small scales decay during the dissipation. The insets in Fig. 9 provide a detailed view around the first modes, revealing that energy and helicity at large scales increase for a duration $t \le 16$ after the external force is removed. The large-scale structures persist after eliminating the force, as shown in Fig. 10.

In Fig. 11(a), we show scatter plots of $y$-components of velocities and vorticities, along with the analytically predicted black curves. Here, no-averaged values $u_p$, $H$, and $W$ are used to calculate Lagrange multipliers $\alpha$ and $\gamma$. During the dissipation phase, noise at small scales decreases, and we observe a clear hyperbolic sine relation between $u_y$ and $\omega_y$. Thus, we can infer that the noise originates from the small-scale motions, and the sinh relation is characteristic of the large-scale coherent structures. Our analytical prediction of the sinh relation in the $y$ direction is confirmed by the overlapping of the predicted curve and data points. For other directions, similar to the forced phase, we observe a linear instead of a sinh relation between velocities and vorticities, as shown in Fig. 11(b).
<figure><img src="image_11.png" /><figcaption>Fig. 8. Time evolution of integral quantities contained in three directions during the freely decaying phase. (a) Energy. (b) Helicity. (c) Enstrophy.</figcaption></figure>

V. CONCLUSIONS

In this paper, we have illustrated a large-scale condensation of energy in 3D turbulence without vortex stretching. As we presented in previous work, helicity cascades from small to large scales²⁷. During the cascade, energy is carried by helicity towards the smallest wavenumbers. In the absence of friction at large scales, energy accumulates and forms a large-scale structure. And this large-scale structure persists when the force is removed.

This condensate of 3D turbulence without vortex stretching behaves like that of 2D turbulence. We showed analytically a generalized 3D point-vortex model to predict a hy-
<figure><img src="image_12.png" /><figcaption>Fig. 9. Time evolution of wavenumber spectra of kinetic energy (a) and helicity (b) during the freely decaying phase at $t = 0, 1, 2, 4, 8, 16$. The insets in (a) and (b) provide a detailed view around the first modes.</figcaption></figure>

perbolic sine relation between vorticity and velocity in the direction which contains most energy. In a forced system, after the condensate occurs, the predicted hyperbolic sine relation is observed between averaged vorticity and velocity, which suggests that the statistical mechanics approach can be applied to averaged systems, as also observed in experiments$^{24}$. At the freely decaying phase, the relation between vorticity and velocity very clearly relaxes to the predicted hyperbolic sine function, similar to the sinh relation between vorticity and stream function in freely decaying 2D turbulence.

ACKNOWLEDGEMENT

We acknowledge the China Scholarship Council (CSC) for the financial support. All simulations were carried out using the facilities of the PMCS2I (École Centrale de Lyon). For the purpose of Open Access, a CC-BY public copyright licence has been applied by the authors to the present document and will be applied to all subsequent versions up to the Author Accepted Manuscript arising from this submission.
<figure><img src="image_13.png" /><figcaption>Fig. 10. Visualization of energy during the freely decaying phase. (a) $t = 1$. (b) $t = 2$. (c) $t = 4$. (d) $t = 8$. (e) $t = 16$.</figcaption></figure>
<figure><img src="image_14.png" /><figcaption>Fig. 11. Scatter plots of the components of velocities and vorticities during the freely decaying phase. (a) $u_y$ versus $\omega_y$ at time instants $t = 0, 1, 2, 4, 8, 16$ from dark to light. (b) x- and z-components of velocities and vorticities at $t = 16$. Black lines in (a) and (b) are the analytically predicted curves at $t = 16$, given by $\omega_y = 1.511 \sinh(0.3718u_y)$.</figcaption></figure>

## Appendix A: Supplemental simulations

To assess the robustness of the results presented above, we also executed eight supplemental DNS-runs using the same initial energy spectrum but different initial complex phases of $\hat{u}$. We find that the $y$ direction is not necessarily the direction which contains most energy. Furthermore, in one of the eight simulations, energy in two of the three directions is comparable and larger than that in the third direction. In this case the sinh relation is found in the two directions containing larger energy, presented in Fig. 12(b-d).

In all our nine simulations we have thus observed the same type of symmetry breaking where two components of the kinetic energy take the same value, the other being either larger or smaller. No observations are reported where all three components tend to the same value.

## REFERENCES

¹R. H. Kraichnan. Inertial ranges in two-dimensional turbulence. *Phys. Fluids*, 10(7):1417–1423, 1967.
<figure><img src="image_15.png" /><figcaption>(a)</figcaption></figure>

<figure><img src="image_16.png" /><figcaption>(b)</figcaption></figure>

<figure><img src="image_17.png" /><figcaption>(c)</figcaption></figure>

<figure><img src="image_18.png" /><figcaption>(d)</figcaption></figure>

Fig. 12. Simulation results of the case where energy in two of the three directions are larger than that in the third direction. (a) Mean square values of velocity components in three directions. (b-d) Scatter plots of averaged components of velocities and vorticities in three directions at $t = 126$. (b) Relation between $\langle u_x \rangle_x$ and $\langle \omega_x \rangle_x$. (c) Relation between $\langle u_y \rangle_y$ and $\langle \omega_y \rangle_y$. (d) Relation between $\langle u_z \rangle_z$ and $\langle \omega_z \rangle_z$. Yellow lines in (b), (c) and (d) are the analytically predicted curves $\langle \omega_m \rangle_m = 1.5391 \sinh(0.3479 \langle u_m \rangle_m)$.

<footnote n="2">²L. M. Smith and V. Yakhot. Finite-size effects in forced two-dimensional turbulence. *J. Fluid Mech.*, 274:115, 1994.</footnote>

<footnote n="3">³M. Chertkov, C. Connaughton, I. Kolokolov, and V. Lebedev. Dynamics of energy condenstation in two-dimensional turbulence. *Phys. Rev. Lett.*, 99:084501, 2007.</footnote>
⁴A. van Kan and A. Alexakis. Condensates in thin-layer turbulence. *J. Fluid Mech.*, 864:490–518, 2019.

⁵H. Xia, D. Byrne, G. Falkovich, and M. Shats. Upscale energy transfer in thick turbulent fluid layers. *Nat. Phys.*, 7(4):321–324, 2011.

⁶F. Bouchet and E. Simonnet. Random changes of flow topology in two-dimensional and geophysical turbulence. *Phys. Rev. Lett.*, 102(9):094504, 2009.

⁷F. Bouchet and A. Venaille. Statistical mechanics of two-dimensional and geophysical flows. *Phys. Rep.*, 515:227, 2012.

⁸J. G. Esler and T. L. Ashbee. Universal statistics of point vortex turbulence. *J. Fluid Mech.*, 779:275–308, 2015.

⁹G. Gauthier, M. T. Reeves, X. Yu, A. S. Bradley, M. A. Baker, T. A. Bell, H. Rubinsztein-Dunlop, M. J. Davis, and T. W. Neely. Giant vortex clusters in a two-dimensional quantum fluid. *Science*, 364(6447):1264–1267, 2019.

¹⁰T. D. Lee. On some statistical properties of hydrodynamical and magnetohydrodynamical fields. *Q. Appl. Math*, 10:69, 1952.

¹¹R. Robert and J. Sommeria. Relaxation towards a statistical equilibrium state in two-dimensional perfect fluid dynamics. *Phys. Rev. Lett.*, 69(19):2776, 1992.

¹²J. Miller. Statistical mechanics of Euler equations in two dimensions. *Phys. Rev. Lett.*, 65:2137–2140, 1990.

¹³G. L. Eyink and K. R Sreenivasan. Onsager and the theory of hydrodynamic turbulence. *Rev. Mod. Phys.*, 78(1):87, 2006.

¹⁴D. Montgomery, W. H. Matthaeus, W. T. Stribling, D. Martinez, and S. Oughton. Relaxation in two dimensions and the “sinh-Poisson” equation. *Phys. Fluids A*, 4:3, 1992.

¹⁵L. Onsager. Statistical hydrodynamics. *Il Nuovo Cimento*, 6:279, 1949.

¹⁶G. Joyce and D. Montgomery. Negative temperature states for the two-dimensional guiding center plasma. *J. Plasma Phys.*, 10:107, 1973.

¹⁷Y. B. Pointin and T. S. Lundgren. Statistical mechanics of two-dimensional vortices in a bounded container. *Phys. Fluids*, 19:1459, 1976.

¹⁸J. B. Taylor, M. Borchardt, and P. Helander. Interacting vortices and spin-up in two-dimensional turbulence. *Phys. Rev. Lett.*, 102:124505, 2009.

¹⁹T. W. David, L. Zanna, and D. P. Marshall. Eddy-mixing entropy and its maximization in forced-dissipative geostrophic turbulence. *J. Stat. Mech. Theory Exp.*, 2018(7):073206,
2018.

²⁰M. S. Singh and M. E. O’Neill. The climate system and the second law of thermodynamics. *Rev. Mod. Phys.*, 94(1):015001, 2022.

²¹N. Leprovost, B. Dubrulle, and P. H. Chavanis. Dynamics and thermodynamics of axisymmetric flows: Theory. *Phys. Rev. E*, 73:046308, 2006.

²²S. Thalabard, B. Dubrulle, and F. Bouchet. Statistical mechanics of the 3d axisymmetric Euler equations in a Taylor–Couette geometry. *J. Stat. Mech.*, 2014(1):P01005, 2014.

²³B. Qu, W. J. T. Bos, and A. Naso. Direct numerical simulation of axisymmetric turbulence. *Phys. Rev. Fluids*, 2:094608, 2017.

²⁴R. Monchaux, F. Ravelet, B. Dubrulle, A. Chiffaudel, and F. Daviaud. Properties of steady states in turbulent axisymmetric flows. *Phys. Rev. Lett.*, 96:124502, 2006.

²⁵W. J. T. Bos. Three-dimensional turbulence without vortex stretching. *J. Fluid Mech.*, 915:A121, 2021.

²⁶T. Wu and W. J. T. Bos. Statistical mechanics of the Euler-equations without vortex stretching. *J. Fluid Mech.*, 929, 2021.

²⁷T. Wu and W. J. T. Bos. Cascades of enstrophy and helicity in turbulence without vortex stretching. *Phys. Rev. Fluids*, 7:094601, 2022.

²⁸D. Montgomery, X. Shan, and W. H. Matthaeus. Navier–Stokes relaxation to Sinh–Poisson states at finite Reynolds numbers. *Phys. Fluids A: Fluid Dynamics*, 5(9):2207–2216, 1993.

²⁹R. Robert and J. Sommeria. Statistical equilibrium states for two-dimensional flows. *J. Fluid Mech.*, 229:291, 1991.

³⁰D. Montgomery and G. Joyce. Statistical mechanics of “negative temperature” states. *Phys. Fluids*, 17:1139, 1974.

³¹A. Delache, C. Cambon, and F. Godeferd. Scale by scale anisotropy in freely decaying rotating turbulence. *Phys. Fluids*, 26(2):025104, 2014.

³²S. Chen, G. D. Doolen, R. H. Kraichnan, and Z. She. On statistical correlations between velocity increments and locally averaged dissipation in homogeneous turbulence. *Phys. Fluids A: Fluid Dynamics*, 5(2):458–463, 1993.

³³L. P. Wang, S. Chen, J. G. Brasseur, and J. C. Wyngaard. Examination of hypotheses in the Kolmogorov refined turbulence theory through high-resolution simulations. Part 1. Velocity field. *J. Fluid Mech.*, 309:113–156, 1996.

³⁴A. Brissaud, U. Frisch, J. Léorat, M. Lesieur, and A. Mazure. Helicity cascades in fully
developed isotropic turbulence. *Phys. Fluids*, 16(8):1366-1367, 1973.

<footnote n="35"><sup>35</sup>M. Kessar, F. Plunian, R. Stepanov, and G. Balarac. Non-Kolmogorov cascade of helicity-driven turbulence. *Phys. Rev. E*, 92(3):031004, 2015.</footnote>