Non-universal transitions to two-dimensional turbulence

The occurrence of system-scale coherent structures, so-called condensates, is a well-known phenomenon in two-dimensional turbulence and is a consequence of the inverse energy cascade – the energy transfer from small to large scales that is a characteristic of two-dimensional turbulence. Here, the onset of the inverse energy cascade and the ensuing condensate formation are investigated as a function of the magnitude of the force and for different types of forcing. Random forces with constant mean energy input lead to a supercritical transition, while forcing through a small-scale linear instability results in a subcritical transition with bistability and hysteresis. That is, the transition to two-dimensional turbulence is non-universal. For the supercritical case we quantify the effect of large-scale friction on the value of the critical exponent and the location of the critical point.


Introduction
Two-dimensional (2-D) and quasi-2-D flows occur at the macro-and mesoscale in a variety of physical systems. Examples include plasma flow in the solar tachocline (Spiegel & Zahn 1992), Earth's atmosphere near the tropopause (Nastrom, Gage & Jasperson 1984;Gage & Nastrom 1986;Falkovich 1992), stratified layers in the oceans (Vallis 2006), laboratory experiments using electrolyte layers (Sommeria 1986) and soap films (Vorobieff, Rivera & Ecke 1999) and, more recently, also dense bacterial suspensions, where the collective motion of microswimmers induces patterns of mesoscale vortices (Dombrowski et al. 2004;Dunkel et al. 2013;Gachelin et al. 2014). A characteristic feature of 2-D turbulence is the occurrence of an inverse energy cascade (Kraichnan 1967), whereby kinetic energy is transferred from small to large scales. Kraichnan's theoretical prediction has been verified in experimental studies of thin fluid layers (Sommeria 1986;Paret & Tabeling 1997, 1998 and through numerical simulations (Lilly 1969;Frisch & Sulem 1984;Verron & Sommeria 1987), to name only the first few. A detailed overview on 2-D turbulence can be found in the review article by Boffetta & Ecke (2014). In confined systems, 892 A18-2 M. Linkmann, M. Hohmann and B. Eckhardt this self-organisation can result in the formation of large-scale coherent structures (Kraichnan 1967;Hossain, Matthaeus & Montgomery 1983;Sommeria 1986;Smith & Yakhot 1993), so-called condensates (Smith & Yakhot 1993). These may emerge in different forms depending on the geometry and boundary conditions, e.g. as vortex monopoles in the case of experimental conditions, i.e. wall-bounded flows (Sommeria 1986;Paret & Tabeling 1998;Molenaar, Clercx & van Heijst 2004;van Heijst, Clercx & Molenaar 2006), if the Rayleigh damping, that is the friction between the moving fluid and the bottom plate of the confining experimental apparatus, is small, or as vortex dipoles or jets (Smith & Yakhot 1993;Bouchet & Simonnet 2009;Frishman, Laurie & Falkovich 2017) in the case of periodic boundary conditions. The inverse energy cascade in 2-D turbulence is connected with an additional inviscid conservation law, that of enstrophy. Physically, it is a consequence of the anisotropic stretching of small-scale vortices by large-scale strain followed by an alignment of the small-scale velocity field around the stretched vortex with the large-scale strain field, thereby reinforcing the latter (Chen et al. 2006). However, inverse cascades and thus condensates are not specific to 2-D phenomena. They occur whenever fluctuations in one spatial coordinate are suppressed, as is the case in thin fluid layers (Sommeria & Verron 1984;Sommeria 1986;Paret & Tabeling 1997, 1998Shats, Xia & Punzmann 2005;Xia, Shats & Falkovich 2009;Celani, Musacchio & Vincenzi 2010;Xia et al. 2011;Musacchio & Boffetta 2017) or, for instance, in the presence of rapid rotation (Deusebio et al. 2014;Rubio et al. 2014;Gallet 2015), stratification (Sozza et al. 2015) or both (Marino et al. 2013), and in the presence of a strong uniform magnetic field (Gallet & Doering 2015) for weakly conducting flows. Another, fully three-dimensional (3-D) mechanism that leads to inverse energy transfer is breaking of mirror symmetry (Waleffe 1993;Biferale, Musacchio & Toschi 2012). In magnetohydrodynamic turbulence, the latter can result in the formation of magnetic condensates through large-scale dynamo action or the inverse cascade of magnetic helicity (Frisch et al. 1975;Pouquet, Frisch & Léorat 1976). Finally, spectral condensation also occurs in toroidal confined plasmas, to the effect that an analogy between 2-D turbulence and toroidal plasma turbulence exists, at least at the level of theoretical models (Horton & Hasegawa 1994) and in phenomenological terms (Shats et al. 2005).
Inverse energy transfer can thus occur in different physical systems, and one could imagine that the onset thereof may depend on the details of the system, such as the dimensionality or the presence of a magnetic field, for instance. Smooth, supercritical and subcritical transitions between non-equilibrium statistically steady states have indeed been observed in this context. In 3-D rotating domains for example, the nature of the transition between forward and inverse energy transfer with respect to the rotation rate depends on the mechanism by which the condensate saturates (Seshasayanan & Alexakis 2018). In the case of weak or vanishing friction with side or bottom walls, the two saturation scenarios are: (i) saturation by viscous effects as in two dimensions, where the condensate becomes sufficiently energetic for the upscale flux to be balanced by viscous dissipation (Chan, Mitra & Brandenburg 2012), or (ii) saturation by local cancellation of the rotation rate by the counter-rotating vortex that forms part of the condensate (Alexakis 2015). In case (i) the transition is supercritical (Seshasayanan & Alexakis 2018), and in (ii) it is subcritical (Alexakis 2015;Yokoyama & Takaoka 2017;Seshasayanan & Alexakis 2018), showing bistability and hysteresis (Yokoyama & Takaoka 2017). Similar results have been obtained if the magnitude of the forcing is used as a control parameter at a fixed value of the rotation rate (Yokoyama & Takaoka 2017), with random and static forcing both resulting in a subcritical transition. The latter was interpreted as evidence in support of universality. Hysteretic transitions and bistable scenarios also occur in thin layers as a function of the layer thickness (van Kan & Alexakis 2019). Subcriticality in the transition to condensate formation in rapidly rotating Rayleigh-Bénard convection has been connected with non-local energy transfer from the driven scales into the condensate due to persistent phase correlations (Favier, Guervilly & Knobloch 2019).
In summary, transitions in cascade directions from direct to inverse and vice versa have received considerable attention in recent years, Alexakis & Biferale (2018) provide a comprehensive overview thereof. Further to this, certain aspects of condensate dynamics, such as circulation reversals, that occur for weak Rayleigh damping (Sommeria 1986;Molenaar et al. 2004), and vortex breakdown due to viscous boundary layers (Molenaar et al. 2004), have been investigated experimentally and numerically. Spontaneous transitions between condensates and disordered states have been observed in toroidal plasma turbulence (Shats et al. 2005). However, transitions to purely 2-D turbulence have only been studied in the context of wave turbulence described in terms of the Gross-Pitaevsky equation (Vladimirova, Derevyanko & Falkovich 2012), and in active matter. In the former, the transition depends on the details of the small-scale driving, i.e. it is non-universal. In the latter, spatio-temporal chaos and classical 2-D turbulence with a condensate are connected by a subcritical transition (Linkmann et al. 2019(Linkmann et al. , 2020. Here, we extend this work and focus on the transition to two-dimensional turbulence as a function of the intensity and the type of driving, and in the presence of large-scale friction. Conceptually, the 2-D geometry differs substantially from thin layers or rapidly rotating 3-D domains, as the energy transfer is now purely inverse while in the latter two cases 2-D and 3-D dynamics, with the corresponding cascade directions, co-exist. That is, the transition investigated here does not occur between two non-equilibrium statistically steady states with different multiscale dynamics. Instead, in two dimensions one state has multiscale dynamics and the other is a spatio-temporally chaotic state concentrated at the driven scales. In that state the nonlinear interscale transfer is too weak to excite motion at scales outside the driven range of scales. Hence the transition in two dimensions is towards and away from multiscale dynamics, not between different types of such. By means of direct numerical simulations we show that the nature of the transition depends on the type of driving: it is supercritical for random forcing and subcritical if the driving is given by a small-scale linear instability. In the former case we also explore the effect of large-scale friction on the location of the critical point and the value of the critical exponent.

Numerical details
We consider the 2-D Navier-Stokes equations for incompressible flow in a square domain V embedded in the xy-plane with periodic boundary conditions. In this case, the Navier-Stokes equations can be written in vorticity form where u = (u x (x, y), u y (x, y), 0) is the velocity field per unit mass, ω the non-vanishing component of its vorticity ∇ × u = (0, 0, ω), ν the kinematic viscosity, α 0 the Rayleigh damping coefficient and f a solenoidal body force. The subscripts x, y and z denote the respective components of a 3-D vector field. We carry out direct numerical simulations of (2.1) on V = [0, 2π] 2 using the standard pseudospectral method (Orszag 1969) for spatial discretisation in conjunction with full dealiasing by truncation following the 2/3rds rule (Orszag 1971). The initial data consist of random, Gaussian distributed vorticity fields. Owing to the focus on condensate formation and its dependence on large-scale friction, the friction coefficient α was small or set to zero in some of the simulations. In the latter case the condensate saturates on a viscous time scale (Chan et al. 2012;Linkmann et al. 2020) with the consequence that the simulations need to be evolved for a long time in order to obtain statistically stationary states. Similarly long transients occur for low values of α. As such, it was necessary to compromise on resolution, and the simulations were run using 256 2 -512 2 grid points.
In order to study the transitions, we conduct a parameter study for a stochastic, Gaussian distributed and δ-in-time correlated force f s that is applied at scales corresponding to a wavenumber interval [k min , k max ], and compare the results with those obtained with a forcing that is linear in the velocity field (Linkmann et al. 2019(Linkmann et al. , 2020 where γ k is a spherically symmetric Galerkin projector Here, k = |k|, while· denotes the Fourier transform and ν IN > 0 an amplification factor, such that the driving occurs through a linear instability in the wavenumber interval Our previous studies (Linkmann et al. 2019(Linkmann et al. , 2020, where linear forcing through a piecewise constant function as in (2.3) was introduced, were carried out in that context. In order to mimic the functional form of previously proposed single-equation models (Wensink et al. 2012;Słomka & Dunkel 2015), which feature a hyperviscous term, an additional dissipation term ν 2 ω had been used at small scales, i.e. at k > k min , resulting in the following equation where ν 2 > ν for k > k min and zero otherwise. For both f s and f l , statistically stationary states are eventually reached, where the spatio-temporally averaged energy dissipation, ε, balances the spatio-temporally averaged energy input, ε IN , with · V,t = · V t denoting the combined spatial and temporal average. For Gaussiandistributed and δ(t)-correlated random forcing, ε IN is known a priori (Novikov 1965) ε ] with k min = 33 and k max = 40 for all simulations, α the large-scale friction parameter, Re = UL/ν the Reynolds number with respect to the root-mean-square velocity U, the integral length scale L = 2/U 2 ∞ 0 dkE(k)/k and the kinematic viscosity ν. The latter was set to ν = 0.0005 for all simulations. The Reynolds number at the driven scales is Re f = U f L f /ν, with U f = ( kmax kmin dkE(k)) 1/2 denoting the velocity at the driven scales and L f = 2π/(k min + k max ). The large-eddy turnover time is denoted by τ = L/U.
That is, the energy input is a control parameter rather than an observable in simulations using f s . Details of the simulations are summarised in table 1.
For the linear forcing, the energy input is is the energy spectrum. Equation (2.1) with f = f l and the aforementioned enhanced small-scale damping has been solved numerically by Linkmann et al. (2019Linkmann et al. ( , 2020 in the context of transitions to large-scale pattern formation in dense suspensions of active matter. Here, we compare our simulations listed in table 1 against the data of Linkmann et al. (2019Linkmann et al. ( , 2020, as summarised in table 2. All simulations are evolved for several thousand large-eddy turnover times τ = L/U, where U is the root-meansquare velocity and L = 2/U 2 ∞ 0 dkE(k)/k the integral length scale, with E(k) = E(k, t) t .

Random forcing
Before reporting on the results from the parameter study varying F, we briefly discuss dynamical and statistical properties of the simulated flows using three example cases with F = 0.11, F = 0.14 and F = 0.23. In order to facilitate the direct comparison with our previously published results using linear driving, the following presentation and discussion of the example cases is structured similarly to those discussed in Linkmann et al. (2019Linkmann et al. ( , 2020. Time series of the kinetic energy E(t) = |u| 2 V /2 and visualisations of vorticity field samples taken during statistically steady evolution corresponding to the three example  (Linkmann et al. 2019(Linkmann et al. , 2020. The number of grid points in each coordinate is denoted by N, the viscosity by ν, and ν IN , k min , k max and ν 2 are the parameters in (2.2)-(2.4). The Reynolds number Re = UL/ν is based on ν, the root-mean-square velocity U and the integral length scale L = 2/U 2 ∞ 0 dkE(k)/k, with ν = 1.1 × 10 −3 for N = 256 and ν = 1.7 × 10 −5 for N = 1024. The Reynolds number at the driven scales is Re f = U f L f /ν, with U f = ( kmax kmin dkE(k)) 1/2 denoting the velocity at the driven scales and L f = 2π/(k min + k max ). The large-scale friction parameter α = 0 and (ν + ν 2 )/ν = 10 for all simulations. Averages in the statistically stationary state are computed from at least 1800 snapshots separated by one large-eddy turnover time τ = L/U. cases are shown in figure 1. For F = 0.23 a condensate consisting of two counterrotating vortices has formed. The remaining cases do not show large-scale structure formation. The time evolution of E(t) and the representative vorticity fields of the three example cases are qualitatively similar to those obtained with linear forcing discussed in our previous work. Quantitative differences between the mean energy levels reported here and in Linkmann et al. (2020) are due to the choice of parameter values. Please note that the example cases only serve to provide a qualitative overview of the data, which are discussed quantitatively in § 4. Figure 2 presents energy spectra (a) and normalised fluxes (b) for F = 0.11, F = 0.14 and F = 0.23. A scaling range characterised by a scaling exponent of the energy spectrum close to the Kolmogorov value −5/3 and a nearly wavenumber-independent flux only forms at the largest value of F. For smaller F the flux tends to zero rapidly for k < k min , hence dissipation cannot be negligible in this wavenumber range. In all cases the maximum and minimum values of the normalised flux do not add up to unity, which indicates that a substantial amount of energy is dissipated directly in the driving range. Interestingly, for intermediate values of F, the scaling exponent of the energy spectrum is still close but slightly larger than −5/3. For the smallest value of F the energy spectrum scales linearly with k for k < 7, indicative of energy equipartition among Fourier modes in this wavenumber range. A similar transition in the energy spectra in statistically stationary 2-D turbulence occurs if the condensate is avoided through a strong drag term (Tsang & Young 2009), in the sense that the extent of the −5/3 scaling range decreases with increasing large-scale friction and a power law with positive exponent appears at low wavenumbers. However, as a drag term alters the scale-by-scale energy balance, it breaks the zero-flux equilibrium condition that underlies linear scaling in two dimensions, the low-wavenumber scaling in the presence of drag is expected to differ from the absolute equilibrium scaling observed here for α = 0. Indeed, in the former case the spectra are much steeper (Tsang & Young 2009).
Condensates are known to affect inertial-range physics in terms of the properties of the third-order structure function (Xia et al. 2008) and the scaling of the energy spectrum in the inertial range of scales (Chertkov et al. 2007). The spectral slopes observed here for the random forcing case and in the presence of a condensate are  Deviations of the spectral exponent from the Kolmogorov value occur in a variety of turbulent systems. For a modified version of the Kuramoto-Sivashinski equation that allowed systematic deviations from inertial transfer, Bratanov et al. (2013) showed by semi-analytical and numerical means that non-universal power laws arise in spectral intervals where the ratio of linear and nonlinear time scales is wavenumber independent. As strong condensates result in a significant contribution of linear terms to the dynamics, a similar analysis could potentially lead to further insights on the non-universal scaling exponents in 2-D turbulence.

Non-universal transitions
The transition to 2-D turbulence as a function of the energy input has so far only been investigated for a single-equation model of active matter (Linkmann et al. 2019(Linkmann et al. , 2020. Here, we now study the transition for a different kind of forcing and in the presence of large-scale friction, as condensates also occur in the presence of drag (Sommeria 1986;Paret & Tabeling 1998;Danilov & Gurarie 2001;Molenaar et al. 2004;van Heijst et al. 2006;Tsang & Young 2009). Figure 3(a) presents the shell-averaged amplitude of the lowest Fourier modes, i.e. the square root of the average energy at the largest scale, A 1 = √ E(k) k=1 , as a function of F from the parameter study for the random, Gaussian-distributed and δ(t)-correlated forcing. Three main observations can be made from the figure 3(a). First, there is a clear transition point, below which A 1 0 and above which A 1 grows with increasing F, indicating the formation of a condensate and thus the onset of sustained inverse cascade, i.e. 2-D turbulence. Second, the data appear to be continuous at the critical point F c with a possibly discontinuous first derivative. The critical point is approached from above by a power law where F c = F c (α) and γ = γ (α) depend on the value of the large-scale friction coefficient. For α = 0 the functional form A 1 (F) corresponds to the upper branch of the normal form of a supercritical pitchfork bifurcation, that is γ = 1/2. Third, for F F c the amplitude A 1 grows linearly with F in all cases. Equivalently, E(k) k=1 grows linearly with ε IN , which is expected for a sizeable condensate as most of the dissipation should then take place at the largest scales where δk is the grid spacing in Fourier space. The dependence of F c and γ on the large-scale friction coefficient α is further quantified in figure 4. As can be seen in the figure, the approach to the critical point described by the exponent γ , is strongly and nonlinearly dependent on the level of large-scale friction, while the location of the critical point varies little. A least-squares fit of A 1 against F places the critical point at F c = 0.135 for α = 0. For α = 0.0005 we have F c = 0.136 and γ = 0.72, α = 0.001 results in F c = 0.138 and γ = 0.78 and α = 0.005 corresponds to F c = 0.146 and γ = 1. According to the discussion in the previous paragraph, γ = 1 is an asymptotic value in the sense that higher exponents are not expected.
The type of transition is very different if the driving occurs through a small-scale linear instability. Figure 3(b) presents the results of the parameter study carried out by Linkmann et al. (2019Linkmann et al. ( , 2020 as a function of ν IN for the linear forcing specified in (2.2). In contrast to the randomly forced case, the transition is now subcritical (Linkmann et al. 2019(Linkmann et al. , 2020 as evidenced by a discontinuity in the data and the clearly visible hysteresis loop. The latter is discussed by Linkmann et al. (2020) in further detail, with figure 3 showing the same result as figure 5 of the aforementioned publication, using A 1 instead of E(k) |k=1 . As the hysteresis loop is small, one may  expect that the transitions happen at comparable values of a forcing-scale Reynolds number Re f = U f L f /ν, where L f is a length scale that corresponds to the middle wavenumber in the driven range, L f = 2π/(k min + k max ), and U f is the root-mean-square velocity in that range of scales This is indeed the case, the transition occurs at Re f ≈ 20 in the subcritical case (Linkmann et al. 2019) and at Re f ≈ 10 in the supercritical case studied here. Non-universality in the transition to condensate formation also occurs in the Gross-Pitaevsky model of wave turbulence: for small-scale driving by a local-in-scale linear instability, a series of symmetry-breaking sharp transitions occur as function of increasing wave action (Vladimirova et al. 2012). The first statistical symmetry that breaks with increasing condensate growth is isotropy, followed by twofold, threefold and fourfold symmetries. Interestingly, such symmetry breaking does not occur if the driving is realised through a small-scale random process. That is, final turbulent states with different statistical symmetries are obtained for different types of small-scale energy input. In weak turbulence, differences between turbulent states originating from the type of driving arise in the context of information theory through different entropy extraction rates (Falkovich & Shavit 2019). As a non-equilibrium steady state, turbulence is a driven-dissipative system, that is, its phase-space measure will be non-uniform, detailed balance is broken and the information-theoretic entropy of that measure (here, the differential entropy of the phase-space measure with respect to the Lebesgue measure) becomes time dependent. Only driving and dissipation contribute to this time dependence, and Falkovich & Shavit (2019) showed that the rate of change in (differential) entropy depends explicitly on the phase-space measure for random forcing while it is independent thereof if the driving is given by a local-in-scale linear instability. In other words, the information context of the system depends on the type of driving, which may be expected, especially concerning comparisons with random forcing. Here, we did not observe any difference in statistical symmetry between randomly forced 2-D turbulence and 2-D turbulence generated by a small-scale instability, however, this may well be because the condensates attained here and in our previous work are moderate in amplitude.

Conclusions
We here study the formation of the condensate and thereby the transition to 2-D turbulence as a function of the type and amplitude of the forcing. Direct numerical simulations show that the condensate does not appear gradually but in a phase transition. For prescribed energy dissipation the transition is second order, and both the critical point and the critical exponent depend on the value of the large-scale friction coefficient. In this context, we point out that ε does not depend on α for the random, δ(t)-correlated forcing used here, as is the case for time-independent forcing such as Kolmogorov flow (Tsang & Young 2009). However, a series of test simulations using time-independent forcing led to similar results (S. Musacchio & G. Boffetta, private communication). When the forcing is due to a small-scale instability as inspired by continuum models of active matter, the transition is first order (Linkmann et al. 2019(Linkmann et al. , 2020. The phase transitions separate two markedly different types of 2-D dynamics: in 2-D turbulence, energy input is predominantly balanced by large-scale dissipation either in the condensate or through Rayleigh friction, and intermediate scales follow an inertial cascade; in spatio-temporally chaotic states where no condensation occurs, dissipation is spread over the intermediate scales and the properties of the energy transfer are noticeably different and non-universal. In summary, the transition to 2-D turbulence is non-universal in the sense that (i) the type of transition depends on the type of forcing, and (ii) the details of the transition for a given type of forcing depend on other system parameters such as large-scale friction. The presence of these non-universalities naturally motivates questions concerning their origin. Results from rapidly rotating Rayleigh-Bénard convection (Favier et al. 2019) suggest that the hysteretic transition in the linearly forced case may be related to persistent phase correlations between the driven scales and the condensate. Random forcing precludes such a scenario. Further questions concern the theoretical predictions on the dependence of the critical exponent γ on the level of large-scale friction. The value γ = 1 is plausible for strong linear damping by the same argument that predicted a linear dependence of the energy in the condensate on the energy input. Finally, it remains to be seen if symmetry-breaking transitions between condensates as in the Gross-Pitaevsky equation (Vladimirova et al. 2012) also occur for the 2-D Navier-Stokes equations.