Document [original]

Jiahan Wang, Jörn Sesterhenn, Wolf-Christian Müller

Coherent structure detection and the inverse

cascade mechanism in two-dimensional

Navier–Stokes turbulence

Open Access via institutional repository of Technische Universität Berlin

Document type

Journal article | Accepted version

(i. e. final author-created version that incorporates referee comments and is the version accepted for

publication; also known as: Author’s Accepted Manuscript (AAM), Final Draft, Postprint)

This version is available at

https://doi.org/10.14279/depositonce-18247

Citation details

Wang, J., Sesterhenn, J., & Müller, W.-C. (2023). Coherent structure detection and the inverse cascade

mechanism in two-dimensional Navier–Stokes turbulence. In Journal of Fluid Mechanics (Vol. 963).

Cambridge University Press (CUP). https://doi.org/10.1017/jfm.2023.313.

This article has been published in a revised form in Journal of Fluid Mechanics

https://doi.org/10.1017/jfm.2023.313. This version is free to view and download for private research and study

only. Not for re-distribution or re-use. © The Author(s), 2023.

cbed This work is licensed under a Creative Commons Attribution-NonCommercial- NoDerivatives 4.0

International license: https://creativecommons.org/licenses/by-nc-nd/4.0/

Coherent structure detection and the inverse

cascade mechanism in two-dimensional

Navier-Stokes turbulence

Jiahan Wang1, Jörn Sesterhenn2and Wolf-Christian Müller1†

1Zentrum für Astronomie und Astrophysik, ER 3-2, Technische Universität Berlin,

Hardenbergstr. 36a, 10623 Berlin, Germany

2Technische Mechanik und Strömungsmechanik, Universität Bayreuth, Universitätsstrasse 30,

95440 Bayreuth, Germany

(Received xx; revised xx; accepted xx)

Coherent structures in two-dimensional Navier-Stokes turbulence are ubiquitously

observed in nature, experiments and numerical simulations. The present study

conducts a comparison between several structure detection schemes based on the

Okubo-Weiss criterion, the vorticity magnitude, and Lagrangian coherent struc-

tures (LCSs), focusing on the inverse cascade in two-dimensional hydrodynamic

turbulence. A recently introduced vortex scaling phenomenology [B. H. Burgess,

R. K. Scott, J. Fluid Mech., 811:742–756, 2017] allows the quantification of the

respective thresholds required by these methods based on physical properties of

the flow. The resulting improved comparability allows to identify characteristic

relative differences in the detection sensitivity between the employed structure

†wolf-c[email protected]

detection techniques. With respect to the inverse cascade of energy, coherent

structures contribute, as expected, substantially less to the cross-scale flux than

the residual incoherent parts of the flow although the energetically dominant

coherent structures lead to an important large-scale deformation of the energy

spectrum. This cascade inactivity can be understood by an increased misalignment

of strain-rate and subgrid stress tensors within coherent structures. At the same

time, the structures exhibit strong and localised nonlinear cross-scale interactions

that appear to stabilize them. We quantify and interpret the resulting shape

preservation of coherent structures in terms of a multi-scale gradient approach [G.

L. Eyink, J. Fluid Mech., 549:191–214, 2006] as the depletion of strain rotation

and vorticity gradient stretching while the dynamics of the residual fluctuations

are consistent with the vortex thinning picture.

1. Introduction

The formation of structures in turbulent flows is ubiquitously observed in nature

for example in atmospheric flows or oceans. In two-dimensional hydrodynamic

turbulence, the structure formation process is associated with the inverse cascade

of kinetic energy, tranferring this quantity from smaller to ever larger spatial

scales. This well-known phenomenon has already been predicted in the seminal

papers of Kraichnan (1967), Leith (1968) and Batchelor (1969) (KLB) and ob-

served in numerous simulations (e.g. Lilly 1971; Frisch & Sulem 1984; Maltrud

& Vallis 1993; Boffetta & Musacchio 2010) and experiments (e.g. Rutgers 1998;

Paret & Tabeling 1998; Chen et al. 2006). We have chosen this physical system for

the present work since it allows to study naturally emergent coherence typically

appearing in the form of structurally rather simple vortices or combinations of

those. These vortical structures are embedded in a statistically isotropic and tur-

bulent two-dimensional flow which is conveniently accessible to direct numerical

simulation (DNS) and measurement.

An intuitive and commonly accepted defining characteristic of coherence is

persistence for a finite time, which leaves room for more detailed specification. A

mathematically unique definition would not only be beneficial for fluid mechanics

research but also for related disciplines such as astrophysics. There, the problem

of the non-universality regarding vortex identification has been pointed out by

Canivete Cuissa & Steiner (2020) and Yadav et al. (2021) with respect to studies

of the solar atmosphere.

A number of methods for the detection of coherent structures exist that are

often built on different specifications of coherence. In fact, most of the comparative

studies in the literature focus on a specific class of coherence specification, e.g.

Jeong & Hussain (1995) studied various vortex criteria, Hadjighasem et al. (2017)

compared different techniques for Lagrangian coherent structure (LCS) identifi-

cation and Taira et al. (2017) discussed numerous mode decomposition methods.

A comparison and meaningful evaluation of different detection strategies and of

their respective coherence specification requires fiducial physical properties of the

considered flow which can be related to the detected structures.

In the present work, we consider three coherence detection schemes – two

vorticity based, one of Lagrangian type – that use the Okubo-Weiss criterion

(Okubo 1970; Weiss 1991), the vorticity magnitude, and the finite-time Lyapunov

exponent (FTLE) field determining LCSs (Haller & Yuan 2000; Haller 2015).

We compare the detection results by investigating the physical properties of

the coherent and the residual (non-coherent) structures in the two-dimensional

turbulent system.

The work of Ouellette (2012) follows a similar approach, which has led to several

experimental studies (see Liao & Ouellette 2013; Kelley et al. 2013). The low

Reynolds numbers attained in these two-dimensional experiments (Re = 185 (Liao

& Ouellette 2013) and Re = 220 (Kelley et al. 2013)) however do not allow

for, e.g., an adequate investigation of cross-scale turbulent interactions in the

framework of the KLB similarity ansatz. This motivates the use of DNSs in the

present investigation.

More specifically, we consider spectral nonlinear cross-scale fluxes of energy, as

well as their scale-filtered correspondents in configuration space. We also compare

with previous results based on a multi-scale gradient ansatz and employ a re-

cently proposed vortex scaling phenomenology for two-dimensional Navier-Stokes

turbulence (Burgess & Scott 2017, 2018) which applies dimensional arguments

to impose physically motivated constraints onto coherent vortices. The results

obtained by these theoretical approaches serve as physical reference points for

the comparison of structure detection techniques and the interpretation of their

results in relation to the inverse turbulent cascade of energy.

This paper is structured as follows: Section 2 presents the decomposition of

the flow into coherent and residual contributions. Section 3 briefly introduces

coherent structure specifications. Section 4 describes the applied diagnostics and

theoretical concepts. Section 5 presents numerical methods and the parameters

used for simulations and analysis. The main results are presented in section 6. A

conclusion is given in section 7.

2. Physical model and flow decomposition

We consider Navier-Stokes turbulence on a two-dimensional 2π-periodic square of

size Agoverned by the differential equations,

∂tu+ (u·∇)u=−∇p+ν∇2u,(2.1)

∇·u= 0,(2.2)

where u= (ux, uy),pand νare the velocity, pressure and kinematic viscosity,

respectively. The kinetic energy per unit mass E= (1/2A)RAu2dA and the

enstrophy Ω= (1/2A)RAω2dA defined with the vorticity ω=∂xuy−∂yux, are

inviscid invariants in a two-dimensional configuration. The kinetic energy exhibits

an inverse cascade, transferring energy from small to large length scales, contrary

to the enstrophy, which exhibits a direct cascade.

We employ a decomposition of the total vorticity field into a coherent part, ωc,

and a residual/incoherent contribution, ωr, (cf. Ohkitani 1991) to carry out the

analysis of different schemes for coherence detection (cf. section 3 below):

ω=ωc+ωr,

ωc(x) = 









ω(x), (x)⩾thr,

0, (x)< thr,

(2.3)

ωr(x) = 









0, (x)⩾thr,

ω(x), (x)< thr.

(2.4)

Based on the particular specification of coherence, a physical characteristic of the

flow, (x), serves as an indicator of this property, turning the detection into a

thresholding procedure with a fixed threshold, thr. In order to improve compara-

bility of different detection schemes, it is important to gauge their thresholds with

respect to a physical property of the flow (see section 3.4 below). Technical details

of the decomposition are pointed out in appendix A.1. The coherent velocity field,

uc= (uc,x, uc,y), and the residual velocity ur= (ur,x, ur,y), are approximated by

inverting ωc/r =∇ × uc/r in Fourier space, a procedure symbolically represented

by the operator ∇×∇−2, which similarly has been employed in several related

works (see Benzi et al. 1986, 1988; Borue 1994; Okamoto et al. 2007; Yoshimatsu

et al. 2009; Vallgren 2011; Burgess & Scott 2018):

uc/r =−∇ × (∇−2ωc/r),

u=uc+ur.(2.5)

This enables straightforward access to various decomposed turbulent fields and

related quantities such as the Fourier spectrum of kinetic energy per unit mass. It

is defined as E(k) = Pk|ˆ

u(k)|2/2with the wavevector k= (kx, ky)and for the

respective length scale `∼k−1. Fourier-transformed quantities are denoted by ‘ˆ’

and the sum over all wavevectors located on a wavenumber shell, k⩽|k|< k+1, is

indicated as Pk. According to the KLB phenomenology, the spectrum possesses

scaling properties for the inverse kinetic energy and direct enstrophy cascade

ranges, which are

E(k)∼2/3

Ik−5/3, k kf,(2.6)

E(k)∼η2/3

Ik−3, k kf,(2.7)

with kfthe forcing wavenumber at which energy and enstrophy are injected

with an injection rate Ior ηI=k2

fI, respectively. Here, we are interested in

a decomposition of the kinetic energy spectrum as

E(k) = Ec(k) + Er(k) + Ecr(k),(2.8)

where Ec(k) = Pk|ˆ

uc(k)|2/2and Er(k) = Pk|ˆ

ur(k)|2/2are associated to

spectral contributions from purely coherent and residual regions, respectively, and

Ecr(k) = PkRe[ˆ

u∗

c(k)·ˆ

ur(k)] the spectrum resulting from mixed-contributions

of coherent and residual parts, with ‘ ∗’ denoting the complex conjugate.

In the following, we introduce the identification schemes considered for (x)in

equations (2.3) and (2.4) and the choice of the corresponding threshold values

thr.

3. Coherence specifications

In general, schemes for the detection of coherent structures can be grouped into

several categories including threshold methods, modal decomposition methods

such as proper orthogonal decomposition (POD) (Holmes et al. 2012), dynamic

mode decomposition (DMD) (Rowley et al. 2009; Schmid 2010) or spectral proper

orthonal decomposition (SPOD) (Towne et al. 2018), and wavelet methods (see

Okamoto et al. 2007; Yoshimatsu et al. 2009; Farge & Schneider 2015). In this

work, we focus on threshold methods, in particular the approaches based on the

Okubo-Weiss (OW) criterion (Okubo 1970; Weiss 1991), the vorticity magnitude

(VM), and the LCSs (Haller & Yuan 2000; Haller 2015). These schemes are

straightforwardly employed using equations (2.3) and (2.4). The thresholds for

the VM and for the LCS based structure detection are chosen with the help of

vortex scaling (see section 3.4 below).

3.1. Okubo-Weiss (OW)/Q-criterion

A frequently applied quantity for structure identification is the Eulerian velocity

gradient tensor ∇u, which is often investigated in decomposed form ∇u=S+W,

with S= (∇u+ (∇u)T)/2the symmetric strain-rate tensor and W= (∇u−

(∇u)T)/2the skew-symmetric spin tensor. The usage of invariants of ∇u, e.g. the

eigenvalues or the trace, and of its tensor decomposition have led to numerous

identification schemes, see, e.g. (Hunt et al. 1988; Chong et al. 1990; Jeong &

Hussain 1995; Hua & Klein 1998; Zhou et al. 1999; Chakraborty et al. 2005).

However, all of these methods face the problem of objectivity (Haller 2005; Haller

et al. 2016), i.e. they lack invariance under certain transformations of the frame of

reference which combine rotation and translation. Thus for the sake of simplicity,

we restrict ourselves to the well-known Q-criterion (Hunt et al. 1988), whose two-

dimensional equivalent resembles the OW criterion. It is defined as

Q=1

2|W|2− |S|2=1

21

2ω2− |S|2,(3.1)

where for Q > 0vortex dominated/elliptical regions and for Q < 0strain

dominated/hyperbolic regions are detected. Thus, (x) = Q(x)and thr = 0

are set in equations (2.3) and (2.4).

3.2. Vorticity magnitude (VM)

Coherent structures in two-dimensional flows are often most clearly visible in

the spatial distribution of the vorticity. Thus, an intuitive approach is to set

(x) = |ω(x)|in equations (2.3) and (2.4). Furthermore, the vorticity is closely

connected to the Lagrangian-averaged vorticity deviation (LAVD) method (Haller

et al. 2016), which is an objective detection criterion.

3.3. Lagrangian coherent structures (LCSs)

LCSs take the evolution of the flow field into account by determining the pair-

dispersion characteristic of passively advected Lagrangian tracers (Haller & Yuan

2000; Haller 2015). Thus, they reveal structures in the flow, which are neither

captured by the vorticity ωnor variants of the velocity gradient tensor ∇u. To

this end, the flowmap Ft

t0(x0) = x(t;t0,x0)is considered, with x0= (x0, y0)the

initial position at time t0. The detection of LCSs can be realised by determining

the finite-time Lyapunov exponent (FTLE) field which is given by

Λt

t0(x0) = 1

Teddy

log qλC

2(x0),(3.2)

with λC

2the largest eigenvalue of the Cauchy-Green strain tensor Ct

t0(x0) =

∇Ft

t0(x0)T∇Ft

t0(x0). The FTLE is interpreted as a local measure of stretching

and can be calculated forward and backward in time. Thus, the values are set

to (x) = Λt0+Teddy

t0(x)for the forward-in-time and (x) = Λt0−Teddy

t0(x)for the

backward-in-time case. Please note, that for the FTLE case the roles of ωcand

ωrare switched in equations (2.3) and (2.4), meaning that small FTLE values

correspond to coherent regions, contrary to the VM |ω(x)|. This is because large

FTLE values isolate coherent regions as illustrated in figures 4 (c) and (d). To

our knowledge no condition exists for the flowmap integration time. Hence, we

suggest setting it to the large-eddy turnover time Teddy according to section 5,

which is typically the longest characteristic correlation time scale of the system.

Further numerical details for the FTLE calculation are discussed in appendix A.2.

Although more refined LCS approaches exist, for our purposes the FTLE

yields sufficient insight into the flow physics, as high-valued FTLE regions, which

are visually perceived as sharp ridges in the flow, are supposed to materially

separate dynamically distinct domains with different transport characteristics.

For example, these domains mark areas of zero cross-scale energy fluxes in low

Reynolds number systems (cf. Kelley et al. 2013). Furthermore, forward-in-time

FTLE (f-FTLE) ridges are associated with repelling LCSs and backward-in-

time FTLE (b-FTLE) ridges to attracting LCSs, indicating stable and unstable

manifolds in the flow in the sense of dynamical systems theory.

3.4. Determining the threshold: vortex scaling

Two of the three detection schemes considered here include free threshold param-

eters which complicate a meaningful comparison of the detection methods and the

physical interpretation of the detection results. In order to achieve comparability

between the three coherence specifications, the VM and LCS schemes are gauged

by making use of the above-mentioned vortex scaling phenomenology. This

model, which we briefly summarize here for completeness, provides a physically

motivated diagnostic signature which we use as a reference for the highly non-

trivial threshold choice of thr in equations (2.3) and (2.4). The phenomenology

characterizes coherent structures by their vortex area Ain configuration space

instead of the classical wavenumber dependence in Fourier space. Therefore, a

time-dependent vortex number density distribution n(A, t)is defined, which yields

the number of coherent vortices per unit area for a certain vortex area Aat time

t. The model is based on the first three moments of nω2

vwith the vortex intensity

ω2

v. They are the vortex energy Ev, vortex enstrophy Zvand vortex number

Nv, respectively. All three quantities are assumed to be approximately conserved

during the spatial growth of an ‘average’ vortex of area A. The number density

is anticipated to follow a power-law tαiA−riwith exponents αiand ridetermined

via the conservation of Ev,Zvand Nv. The range of areas is divided into a thermal

bath regime Af⩽A<A−, an intermediate scaling regime A−<A<A+and a

front of the vortex population A+< A ⩽Amax, respectively, where A−and A+

are transitional areas, Afthe forcing-scale area, and Amax the maximum vortex

area.

In this model, the thermal bath is associated with the equilibration of the flow

with the continuous forcing, which injects energy at a constant rate generating

small-scale vorticity. This leads to an A-independent flux of Evin A-space. The

intermediate scaling regime consists of a self-similar distribution of vortex sizes. It

is assumed that the enstrophy lost through filament shedding during merger and

aggregation processes is replaced by the enstrophy injection such that the vortex

enstrophy Zvis also approximately conserved. In the front regime, vortices are

expected to be large and distant from each other, such that merging events rarely

occur. Thus, approximately conserving the vortex number Nv. Based on these

conservation assumptions, the scaling laws of the number density for varying area

regimes are derived as (see Burgess & Scott 2017, 2018)

n(A, t)∼











A−3, Af⩽A < A−,

t−1A−1, A−< A < A+,

t5A−6, A+< A ⩽Amax.

(3.3)

We take the best achievable agreement with the three regime subdivision (3.3)

as a reference to gauge the threshold values in (2.3) and (2.4). Please note that this

qualitative level of agreement mainly relies on the assumption that the emergence

and the evolution of coherence are asymptotically self-similar for sufficiently large

scale-separation between the regions of the forcing and the large scales of the

system under consideration. This can only be fulfilled up to a rather modest

approximate level in turbulence DNS. In the present work, the scaling exponents

are considered relative to each other. Thus, their absolute numerical values are

not of principal importance to the investigation. They are nevertheless mentioned

above for completeness.

4. Diagnostic methods for the inverse cascade

The inverse cascade of kinetic energy corresponds to a cross-scale energy flux of

which we distinguish coherent and residual contributions from three perspectives:

(i)spectrally in Fourier space, (ii)scale-filtered in configuration space which com-

bines the aspects of spatial scale and position, and (iii)via a multi-scale gradient

(MSG) approach (Eyink 2006b) which adds scale locality and the differentiation

between involved physical processes.

4.1. Spectral flux

The temporal evolution of the energy spectrum is straightforwardly obtained from

the Navier-Stokes equations (2.1) and (2.2) as:

∂tE(k) + T(k) = D(k) + F(k),(4.1)

with the nonlinear transfer term T(k) = PkRe[ˆ

u∗(k)·\

(u·∇u)(k)]. Kinetic

energy is provided to the flow by a forcing term +fuon the right-hand side

of the Navier-Stokes momentum balance (2.1). Thus, the energy source term is

determined as F(k) = PkRe[ˆ

u∗(k)·ˆ

fu(k)], which is equivalent to the energy

injection rate Iwhen summed over all Fourier wavenumbers. In order to allow for

a statistically stationary state, kinetic energy accumulating at the largest length

scales of the flow due to the inverse cascade has to be continuously extracted

from the system. For this purpose a large scale damping term −dαuis added to

the right-hand side of equation (2.1). The energy sink D(k) = Dν(k) + Dα(k)is

split into two dissipative contributions, where Dν(k) = −2νk2E(k)is the viscous

dissipation active on small length scales and Dα(k) = −2dαE(k)introduces

friction active on large length scales. These terms are equivalent to the energy

dissipation rate on viscous scales νand on large scales α, respectively, when

summed over all wavenumbers. Details on the numerical implementation of +fu

and −dαuare given in the text around equation (5.1) in section 5 .

The spectral cross-scale energy flux Z(k)is obtained by summing the transfer

term over consecutive shells

Z(k) =

k0=0

T(k0) =

k0=0 Xk0Re hˆ

u∗(k)·\

(u·∇u)(k)i,(4.2)

and corresponds to the flux of energy from scales smaller than kto scales larger

than k. The influence of the coherent and residual contributions with regard to

the cascade mechanism is measured by the decomposition

Z(k) = X

α,β,γ∈{c,r}

Zα,β,γ(k),(4.3)

Zα,β,γ(k) =

k0=0

Tα,β,γ(k0) =

k0=0 Xk0Re hˆ

u∗

α(k)·\

(uβ·∇uγ)(k)i,(4.4)

which results in eight independent flux contributions. We investigate the homoge-

neous fluxes originating from purely coherent Zc,c,c(k)and residual components

Zr,r,r(k), and the mixed flux arising through coherent-residual interactions as

Zcr(k) = Z(k)−Zc,c,c(k)−Zr,r,r(k).

4.2. Spatial flux distribution

A complementary formulation of the cross-scale energy flux which captures its

local structure in configuration space and which enables a detailed analysis

regarding its spatial distribution is obtained by a scale-filter approach, cf., e.g.,

(Ouellette 2012). For the i-th component of the velocity vector the filter operation

at a length scale `is given by

ui(x) = ZG`(r)ui(x+r)d2r,(4.5)

where we choose G`as a smooth, non-negative, spatially well localised filter

kernel with unit integral. Because we are interested in the spatial distribution

of the cross-scale flux, the locality aspect of the filter is crucial for the accurate

localisation of flux contributions in configuration space. Hence, in this work we

employ a Gaussian filter: ˆ

G`(k) = exp(−k2`2/24) to achieve sufficient filter

locality in Fourier space as well as in configuration space. The temporal evolution

of the filtered kinetic energy, E=|u|2/2, is given by (see Pope 2000)

∂tE+∂iqi=−ν−Z, (4.6)

with ∂ithe partial derivative of the i-th component, where we use the Einstein

summation convention. Additionally, qi=uiE+uj(pδij +τij −2νSij)contains

the nonlinear spatial transport and the viscous dissipation of the filtered large-

scale kinetic energy with δij the Kronecker delta function, ν= 2νSijSij the

viscous dissipation from the filtered velocity field, and Z=−Sij τij the spatial

cross-scale flux term representing the exchange of kinetic energy between the

known filtered fields and the fluctuations which have been depleted by the filter

operation in the filtered numerical system. The flux term is an inner product

of the (filtered) strain-rate tensor Sij = (∂jui+∂iuj)/2and the subgrid stress

tensor τij =uiuj−uiuj, that expresses the stresses exerted by the depleted

fluctuations. Please note that we are referring to the deviatoric (trace-free) stress

term ˚

τ=τ−(1/2)tr(τ)I, with tr the trace operator and Ithe unit matrix. For

the remainder we will write τinstead of ˚

τ. The production term Zis equally

understood as a spatial cross-scale flux term. Note that choosing a sharp spectral

filter instead of a smooth Gaussian filter will lead to the equality between the

spatial average of the production term and the spectral flux in equation (4.2) as

hZ(x)i=Z(k= 2π/`), if the wavenumber kis chosen according to the filtering

length scale `.

The strain-rate and stress tensor are further decomposed into coherent, residual

and mixed contributions

Z(x) = X

α,β,γ∈{c,r}

Zα,β,γ(x) = X

α,β,γ∈{c,r}

Sα:τβ,γ,(4.7)

with Sα= (∇uα+(∇uα)T)/2and τβ,γ =uβuγ−uβuγ. We propose the following

three-part decomposition

Z(x) = −S:τ=−Sc:τc,c

| {z }

Zc(x)

−Sr:τr,r

| {z }

Zr(x)

−Sc:τr,r −Sr:τc,c −Sc: (τc,r +τr,c)−Sr: (τc,r +τr,c)

| {z }

Zcr(x)

(4.8)

where Zc(x)consists of purely coherent and Zr(x)of purely residual contribu-

tions, and Zcr(x)is the flux contribution originating from mixed interactions.

Because the spatial cross-scale flux consists of an inner product of two tensors,

the analysis of angle alignments between tensor eigenframes is possible. Thus, a

polar decomposition leads to the following expression for the total and decom-

posed fluxes (see Eyink 2006b; Fang & Ouellette 2016)

Z(x) = −2σ(x)λ(x) cos(2δθ(x)),(4.9)

Zc/r(x) = −2σc/r(x)λc/r(x) cos(2δθc/r(x)),(4.10)

respectively. The positive eigenvalues of the strain-rate and subgrid stress tensors

are σand λ, respectively, and the angle between their corresponding eigenvectors

is δθ as illustrated in figure 1 (a). The same definitions are used for the eigenvalues

and angles of coherent and residual parts, which are indicated by the indices

cand r, respectively. The cosine of the rotation angle between strain-rate and

stress tensors, cos(2δθ), can be understood as an efficiency of the cross-scale

energy transfer (Fang & Ouellette 2016). Therefore, a detailed analysis of angle

distributions from coherent δθcand residual parts δθris conducted in section 6.1.

The mixed cross-scale flux Zcr in equation (4.8) is a very complex object due to

the heterogeneous subgrid stress tensors, τc,r and τr,c, which are not symmetric

and thus not straightforward to interpret. Only the sum of τc,r +τr,c yields a

symmetric stress quantity. Thus, the mixed cross-scale flux contribution consists

of a sum of four different physical contributions: (i)Exertion of residual stress on

coherent strain-rate, (ii)exertion of coherent stress on residual strain-rate, (iii)

exertion of mixed stress on coherent strain-rate, and (iv)exertion of mixed stress

on residual strain-rate. For conciseness of this paper, we abstain from analysing

all the single contributions of this mixed flux regarding their rotation angles, and

focus on the sum of all four contributions altogether.

4.3. Multi-scale gradient (MSG) flux expansion

As a final extension of the flux analysis, the locality between strain-rate tensors

on varying scales is analysed according to the second-order MSG approach (Eyink

2006a,b). For that, a second filtering operation is defined as

u(b)

i(x) = ZG`b(r)ui(x+r)d2r,(4.11)

where G`bfilters out contributions from all scales smaller than `b=λ−b`, with a

geometric factor λ > 1. This leads to the band-pass filtered velocity

u[b]

i=









u(b)

i−u(b−1)

i, b ⩾1,

ui, b = 0,

(4.12)

representing contributions from a band of length scales between `band `b−1. The

filtering operation leads to the multi-scale property of the MSG expanded cross-

scale flux approach. The multi-gradient nature comes from a Taylor expansion of

(a) (b)

(c)

Scale ux e ciency Vortex thinning

Figure 1: Overview of different angles used in this work. (a): Angle δθ between the strain-rate

Sand subgrid stress tensor τ. (b): Angle δα[b]between the large-scale strain-rate tensor S(0)

and the band-pass filtered strain-rate tensor S[b]. Remade from (Eyink 2006b). (c): Angle δβ[b]

between the contractile direction of the large-scale strain-rate tensor S(0) and the band-pass

filtered vorticity gradient vector ∇ω[b].

the velocity increments δu(r;x) = u(x+r)−u(x)with separation vector r. The

technical details for the derivation of the second-order MSG flux are outlined in

appendix A.3 and yield (see Eyink 2006b):

ZMSG

∗=−S(0) :τMSG

∗

=ZMSG

∗(S(0) :˜

S[b],S(0) :S[b],(∇ω[b])TS(0)(∇ω[b]))

b=0

(Z[b]

SR +Z[b]

DSR +Z[b]

DSM +Z[b]

V GS)−Z(nb)

F SF .(4.13)

The parameter b∈N0denotes the level of scale locality of the respective MSG

flux contributions Z[b]

SR,Z[b]

DSR,Z[b]

DSM and Z[b]

V GS, meaning that for low b-values

contributions from strongly scale local interactions are measured, whereas contri-

butions of non-local interactions are obtained for larger values. The total number

of filter bands is denoted as nb. The inner products between tensors, as well as

matrix vector products are expressible in polar coordinates as

S(0) :˜

S[b]=−σ(0)σ[b]sin(2δα[b]),(4.14)

S(0) :S[b]=σ(0)σ[b]cos(2δα[b]),(4.15)

(∇ω[b])TS(0)(∇ω[b]) = −σ(0)|∇ω[b]|2cos(2δβ[b]),(4.16)

where σ(0) and σ[b]are the positive eigenvalues of the strain-rate tensors S(0)

and S[b], respectively, with α(0) and α[b]the angles between their corresponding

eigenvectors to a fixed orthogonal frame of reference, and ˜

S[b]the skew-strain-

rate matrix rotated counterclockwise by π/4to the original strain matrix S[b].

According to figure 1 (b), δα[b]=α[b]−α(0) is the rotation angle between the

large-scale tensor S(0) and the subfilter-scale tensors S[b]. Figure 1 (c) shows δβ[b],

which is the angle between the vorticity gradient vector ∇ω[b]and the eigenvector

of S(0) corresponding to the negative eigenvalue. The latter is equivalent to its

contractile direction.

The second-order MSG flux can be subdivided into four flux channels, in which

the investigation of the angles δα[b]and δβ[b]directly illuminates the proposed

vortex thinning picture (Eyink 2006b; Xiao et al. 2009):

•The strain rotation (SR) Z[b]

SR is equivalent to the first-order MSG expansion

and relates to the following physical picture: A small-scale vortex ω[b]embedded

in a large-scale strain-rate field S(0), as illustrated in figure 1 (b), is stretched

along the positive and compressed along the negative eigendirection of the strain.

This leads to an elliptical shape inducing a shear layer and thus a small-scale

strain rotated with δα[b]=±π/4towards the large-scale strain, depending on the

sign of the vorticity.

•The differential strain rotation (DSR) Z[b]

DSR contains a Newtonian stress-

strain relation of the form τ[b]=−ν[b]

TS[b], with negative eddy-viscositiy ν[b]

According to figure 1 (b), the elliptically-shaped vortex still possesses the same

area, but the circumference increases leading to a loss of energy, due to Kelvin’s

theorem of the conservation of circulation Γ=Hu·ds. As a consequence, the

small-scale stress τ[b]exerts negative work on the large-scale strain S(0) because

of its parallel alignment to that strain. This results in an energy transfer towards

larger length scales.

•The vorticity gradient stretching (VGS) Z[b]

V GS is a measure of the elongation

of vortex lines. The angle δβ[b]between the vorticity gradient vector ∇ω[b]and

the contractile direction of the large-scale strain S(0) measures the alignment of

the stretching direction to the vorticity isolines. Thus, a higher tendency for this

alignment increases the rate of stretching parallel to the isolines, as depicted in

figure 1 (c), leading to a thinning of the vortex.

•The differential strain magnification (DSM) Z[b]

DSM contains, similar to the

SR term, a skew-Newtonian stress-strain relation with skew-eddy-viscositiy γ[b]

It measures the logarithmic rate of strain increase, when moving in the direction

of increasing vorticity. According to Xiao et al. (2009), this term is generally

expected to be smaller, as we can confirm in the results of section 6.2.

Although the vortex thinning picture is not necessarily associated with single

coherent vortices, but rather with the whole vorticity ensemble itself, we intend to

measure the influence of coherent regions and their residual backgrounds regarding

this mechanism. This is achieved by the following three-part decomposition:

ZMSG

∗=−S(0)

c:τMSG

∗,c,c

| {z }

ZMSG

∗,c

−S(0)

r:τMSG

∗,r,r

| {z }

ZMSG

∗,r

−S(0)

c:τMSG

∗,r,r −S(0)

r:τMSG

∗,c,c −S(0)

c: (τMSG

∗,c,r +τMSG

∗,r,c )−S(0)

r: (τMSG

∗,c,r +τMSG

∗,r,c )

| {z }

ZMSG

∗,cr

(4.17)

with ZMSG

∗,c and ZMSG

∗,r the purely coherent and residual second-order MSG flux

expansions, respectively, and ZMSG

∗,cr the flux contribution originating from the

mixed interactions. For the reasons similar to those already given in section 4.2,

the same form of heterogeneous stresses, τMSG

∗,c,r and τMSG

∗,r,c , appears in the mixed

MSG expanded flux. To limit the scope of this paper, we refrain from an in-depth

analysis of the vortex thinning angles δα[b]and δβ[b]for the mixed MSG flux

contribution. However, the decomposition of the MSG expanded flux into purely

coherent and residual parts implies a decomposition of the different flux channels

as well

ZMSG

∗,c/r =

b=0

(Z[b]

SR,c/r +Z[b]

DSR,c/r +Z[b]

DSM,c/r +Z[b]

V GS,c/r)−Z(nb)

F SF,c/r.(4.18)

This leads to the analysis of different angles between strain-rate tensors and

vorticity gradient vectors for varying scale localities, set by b, originating from

coherent and residual components

S(0)

c:˜

S[b]

c=−σ(0)

cσ[b]

csin(2δα[b]

c),(4.19)

S(0)

r:˜

S[b]

r=−σ(0)

rσ[b]

rsin(2δα[b]

r),(4.20)

S(0)

c:S[b]

c=σ(0)

cσ[b]

ccos(2δα[b]

c),(4.21)

S(0)

r:S[b]

r=σ(0)

rσ[b]

rcos(2δα[b]

r),(4.22)

(∇ω[b]

c)TS(0)

c(∇ω[b]

c) = −σ(0)

c∇ω[b]

c

2cos(2δβ[b]

c),(4.23)

(∇ω[b]

r)TS(0)

r(∇ω[b]

r) = −σ(0)

r∇ω[b]

r

2cos(2δβ[b]

r).(4.24)

The variables are interpreted in the same fashion as above for the total field but

now with respect to the coherent (index c) and residual (index r) contributions.

We present an analysis of the thinning effects in section 6.2, for which Z[b]

SR,c/r,

Z[b]

DSR,c/r,Z[b]

V GS,c/r and their corresponding angles δα[b]

c/r,δβ[b]

c/r are the relevant

quantities measuring the thinning tendencies of purely coherent and residual

parts, respectively.

5. Numerical methods and parameters

Equation (2.1) and (2.2) are solved in Fourier space using the equivalent and

numerically more favourable vorticity representation. The differential equation

includes a small-scale forcing term, ˆ

fω, and a large-scale damping function, −ˆ

dωˆω,

yielding

∂tˆω+\

[(u·∇)ω] = νk2ˆω+ˆ

fω−ˆ

dωˆω. (5.1)

It is solved by a pseudospectral approach, with a second-order trapezoidal

Leapfrog time integration scheme, and a 2/3-dealiasing method (Canuto

et al. 1988). The forcing components of the velocity field ˆ

fu,x and ˆ

fu,y

are drawn from Gaussian normal distributions and they are afterwards

projected onto the solenoidal components ˆ

fu,j = (δij −kikj/k2)ˆ

fu,i to satisfy

the incompressibility condition. The forcing term is then constructed as

fω(k) = i kxˆ

fu,y(k)−kyˆ

fu,x(k)and applied at a wavenumber of kf= 200.

In order to avoid the accumulation of energy at large scales due to the inverse

cascade, a large-scale linear damping term with a Gaussian damping factor

dω(k) = αωexp(−(k−k0,ω)2/(2σ2

ω)) is employed. The parameters for the large-

scale friction factor αω, the center of the Gaussian damping profile k0,ω and

its variance σ2

ωare given in table 2 in appendix B. We solve the system at a

resolution of 40962in the square periodic domain 2π×2π.

The large-eddy turnover time is estimated as Teddy =Lint/urms, with Lint =

Rk−1E(k)dk/ RE(k)dk the integral length scale and urms the root-mean square

velocity, which is also used for the integration time of the passive tracers in the

FTLE/LCS calculation in equation (3.2). Our results are taken after reaching a

statistically stationary state, as confirmed in figure 2 (a). They are averaged over

100 snapshots equidistantly distributed over roughly 20Teddy.

A discussion of the chosen values of the energy injection rate Iand the general

system parameters for the present numerical setup can be found in appendix B,

which is related to the characteristics of structure formation, the kinetic energy

spectrum and the cross-scale kinetic energy flux. There, we conclude that run1

(table 2 in appendix B) is the best choice for the purpose of our present study. The

spatial vorticity distribution in figure 2 (b) exhibits a clearly developed population

of visually distinguishible vortices or coherent structures. The kinetic energy

spectrum in figure 2 (c) deviates from the theoretically expected k−5/3scaling

due to finite-size effects discussed in appendix B, but we deem it to be more

adequate for the subsequent analysis due to a more clearly discernible structuring

of the flow. The cross-scale kinetic energy and enstrophy fluxes shown in figure

2 (d) possess sufficiently extended ranges of inverse and direct spectral transfer.

This facilitates the cross-scale flux decompositions in sections 6.1 and 6.2.

5.1. Structure detection

As already mentioned in the previous sections 3.2 and 3.3, the threshold choice for

the VM, the f-FTLE and the b-FTLE criterion to sample coherent regions from

the vorticity distribution is not straightforward. Therefore, the threshold is chosen

such that the three-subregime structure of the coherent vortex number density,

n(A), in equation (3.3) is realized most clearly. The system studied by Burgess &

Scott (2017, 2018) assumed stationarity by imposing an integral length scale far

below the largest length scales of the system domain. On the contrary, our system

is in a statistically stationary state with constant kinetic energy Eand enstrophy

Ωaccording to figure 2 (a). Therefore, we analyse the scaling sensitivity of the

number density n(A)only in dependence of the coherent area Awithout the time

t, as presented in figures 3 (b), (e), (h) and (k).

Our comparison of coherence specifications and detection methods begins with

the VM criterion. This approach is technically closest to the detection method

chosen in Burgess & Scott (2017, 2018) and thus expected and observed to yield

the best agreement with the findings published there and with the corresponding

asymptotic scaling laws of the vortex number density, n(A). We vary the VM

threshold introduced in section 3.2 with three different values thr = 1.1ωrms,

0.9ωrms,0.7ωrms, for which the total areas occupied by the coherent regions in

figure 3 (a) amount to 3.4%,6.1%,13.5%, respectively. Filamentary structures

occur with increasing area occupation extending the thermal bath regime in figure

3 (b). The number density approximately exhibits the phenomenologically pro-

0 100 200 300

t/T

eddy

E, Ω(×10

−3

)

(a)

Ω

/max(

)

(b)

110

−8

−6

−4

−2

E(k)

−5/3

−3

(c)

110

−3

−1

E(k)k

5/3

−2/3

110

-1

-0.8

-0.6

-0.4

-0.2

0.2

(k)/

(d)

0.2

0.4

0.6

0.8

Ω

(k)/

-1

-0.75

-0.5

-0.25

0.25

0.5

0.75

Figure 2: Flow observables of run1 (see table 2 in appendix B). (a): Time evolution of the

kinetic energy Eand enstrophy Ω. The enstrophy is divided by 103for better visualisation.

(b): Vorticity ωfrom a 10242region, which is 6.25% of the total physical domain. The colorbar

is normalised to ω/ max |ω|for better visualisation of emerging structures. The upper circle

highlights a vortex pair and the lower circle marks a single large-scale vortex. Physical

quantities for both of these regions are further analysed in figure 4. (c): Kinetic energy

spectrum E(k)with inset showing the compensated spectrum E(k)k5/3−2/3

I, where the black

dashed line indicates a value of CE= 6.69 predicted by the test-field model (TFM) closure of

Kraichnan (1971). (d): Normalised cross-scale kinetic energy flux ZE(k)/I, where the inset

shows the normalised enstrophy flux ZΩ(k)/ηI.

posed asymptotic A−3scaling for the highest area occupation. The approximate

power-law deteriorates as the threshold is raised and the detected set comprises of

fewer structures. As the VM criterion explicitly detects regions of high vorticity,

the most energetic regions of the flow are included in the coherent part, which

is evident from the decomposed kinetic energy spectra in figure 3 (c). Most of

the energy is concentrated in Ec(k), whereas the energy in the residual spectrum

Er(k)decreases for larger length scales with scalings which become flatter than

k−5/3for increasing area occupations. The clearly discernible intermediate region

exhibits an increased fluctuation level for lower thresholds since larger and thus

fewer structures tend to be detected. This is consistent with larger deviations of

the measured scaling exponent from the predicted value than in the case of the

thermal bath as seen in table 1.

In contrast to the VM method, the OW criterion only permits one possible

threshold thr = 0. The detected coherent regions in figure 3 (d) amount to an

area occupation of 5.3%. The thermal bath range obtained via the OW method

drops steeper with increasing area Athan for the VM technique (cf. table 1). This

indicates a preference of the OW criterion for intense vortices at the cost of lesser

vortical structures in the interval Af⩽A<A−. The intermediate range of the

number density n(A)in figure 3 (e) has the best fulfillment of the A−1scaling

compared to the other methods according to table 1. This is not surprising, since

the OW criterion favours circularly shaped structures, which are also detected by

the modified vorticity threshold criterion used by Burgess & Scott (2017, 2018).

The residual energy spectrum Er(k)in figure 3 (f) possesses a scaling closer to

k−5/3in the inertial range from k≈10 −200 versus the total spectrum E(k).

In comparison, the coherent spectrum Ec(k)has a much shorter k−5/3scaling

range from k≈10 −30, which becomes shallower for increasing wavenumbers.

This indicates that the lacking energetic self-similarity of largest-scale coherent

structures detected by the OW criterion pollutes the scaling of the theoretically

expected KLB spectrum of total energy (see appendix B). The residual energy

whose scaling is not suffering from this specific finite-size effect of the numerical

simulation exhibits good agreement with KLB expectations, see also, e.g., (Borue

1994; Scott 2007; Vallgren 2011; Burgess & Scott 2018).

With regard to the sensitivity of the LCS detection, presented in section 3.3,

we set the thresholds of the f-FTLE sampling scheme to thr = 0.87Λt0+Teddy

t0,rms ,

0.905Λt0+Teddy

t0,rms ,0.94Λt0+Teddy

t0,rms , which amount to coherent area occupations of 3.1%,

6.0%,11.2%, respectively. The thresholds for the b-FTLE are set to thr =

0.88Λt0−Teddy

t0,rms ,0.915Λt0−Teddy

t0,rms ,0.94Λt0−Teddy

t0,rms , resulting in coherent area occupations

of 3.3%,6.8%,10.9%, respectively. The crinkly-shaped structures occuring in the

domain, as shown in figures 3 (g) and (j), lead to the flattening of the thermal bath

regime and a simultaneous steepening of the intermediate range with increasing

area occupation, as depicted in figures 3 (h) and (k). Note that the different

structural shape also results in a more extended thermal bath regime and a

diminishing intermediate range with increasing area occupation in contrast to

the OW and VM criterion. Comparing f-FTLE and b-FTLE fields, similar larger-

sized structures are detected by both fields, whereas smaller-sized structures are

detected at differing positions. This is not surprising, since forward- and backward-

in-time LCSs, are associated with different fluid dynamics, i.e. repelling and

attracting manifolds in a dynamical systems sense, respectively (Haller & Yuan

2000; Haller 2015). The residual energy spectra Er(k)in figures 3 (i) and (l) are

closer to a k−5/3scaling throughout the entire inertial range compared to the

total spectrum E(k), indicating a pollution of the KLB spectrum by the coherent

structures similar to the results of the OW criterion.

The large-scale damping required for the achieving stationarity impacts the

apparition of large-scale structures as discussed in appendix B. This leads to the

anomalous and partially polluted power-laws of the front regime in figures 3 (b),

(e), (h) and (k) compared to the phenomenologically expected exponent of A−6

(Burgess & Scott 2017, 2018), where a large-scale damping mechanism has not

been employed. Another reason for the deviation is the resolution of 81922in

the conducted DNS of Burgess & Scott (2017, 2018), which is higher compared

to our present study. Nevertheless, the overall trend of the three-part number

density n(A)is satisfied in all of these definitions. Moreover, independent of the

definition, coherent structures tend to distort the theoretically predicted KLB

scaling as illustrated in figures 3 (c), (f), (i) and (l).

In summary and in comparison to the VM detection technique the OW and the

LCS approaches tend to detect similar large-scale coherence features of the flow,

in spite of their distinct underlying coherence specifications. The different char-

acterizations of coherence lead to a reduced sensitivity for small-scale turbulent

structures in the case of OW, while the physically most involved LCS method

tends to detect a surplus of small-scale structures compared to the most simple

VM specification. We proceed by setting the respective detection parameters

such that the detection signatures as shown in figures 3 (b), (e), (h) and (k)

become most similar to each other, i.e. an area occupation of 6.1% for VM, of

6.0% for f-FTLE, and of 6.8% for b-FTLE. The resulting scalings in table 1

reflect a rough overall agreement with the three expected scaling regimes of the

number density n(A). Small adjustments of these thresholds do not result into

(a)

2.6 3.0 3.4 3.8 4.2

log

(A)

-5

-4

-3

-2

-1

log

(n(A)/N)

−3

−1

−6

(b)

area

occupation: 3.4%

area occupation: 6.1%

area occupation: 13.5%

(k)k

5/3

(c)

(d)

2.6 3.0 3.4 3.8

log

(A)

-5

-4

-3

-2

-1

−3

−1

−6

(e)

area

occupation: 5.3%

(f)

(g)

f-FTLE

2.6 3.0 3.4 3.8

log

(A)

-5

-4

-3

-2

-1

−3

−1

−6

(h)

area

occupation: 3.1%

area occupation: 6.0%

area occupation: 11.2%

(i)

(j)

b-FTLE

2.6 3.0 3.4 3.8

log

(A)

-5

-4

-3

-2

-1

−3

−1

−6

(k)

area

occupation: 3.3%

area occupation: 6.8%

area occupation: 10.9%

(l)

Figure 3: Coherent structures, number density and kinetic energy spectra obtained by the VM

(a)-(c), OW criterion (d)-(f), f-FTLE (g)-(i) and b-FTLE (j)-(l), respectively. The coherent

areas are color-coded matching the normalised number densities n(A)/N and the decomposed

energy spectra, with Nthe total number of detected structures. Color-code: (i)lowest area

occupation: green regions corresponding to the green lines of the number densitiy and

decomposed energy spectra plots, (ii)intermediate area occupation: green+blue regions

corresponding to the blue lines and (iii)largest area occupation: green+blue+red regions

corresponding to the red lines. The OW criterion has only one permitted threshold, hence only

a single area occupation is shown. All the energy spectra are compensated by k−5/3, with the

black solid line indicating the total energy spectrum E(k), and the colored solid and dashed

lines denoting the coherent Ec(k)and residual Er(k)spectra, respectively. The black dotted

horizontal line visualises deviations from the k−5/3scaling.

qualitatively significant differences regarding the inverse cascade analysis, which

is conducted in the next section.

Method log10(Af) log10(A−) log10(A+) log10(Amax) (i) (ii) (iii)

VM 2.6 3.44 3.94 4.14 −2.92 ±0.09 −1.38 ±0.11 −5.15 ±0.86

OW 2.57 3.42 3.72 3.95 −3.41 ±0.07 −1.01 ±0.13 −4.63 ±0.66

f-FTLE 2.57 3.53 3.66 3.73 −2.63 ±0.01 0.11 ±0.31 −24.45 ±2.24

b-FTLE 2.56 3.45 3.58 3.68 −2.8±0.02 0.08 ±0.31 −23.79 ±1.68

Table 1: Logarithmic power-law scalings of the three-part vortex number densities n(A)

corresponding to structures detected by the VM (thr = 0.9ωrms, area occupation: 6.1%), OW

criterion, f-FTLE (thr = 0.905Λt0+Teddy

t0,rms , area occupation: 6.0%) and b-FTLE

(thr = 0.915Λt0−Teddy

t0,rms , area occupation: 6.8%), respectively. The scaling ranges are denoted

as: (i)thermal bath Af⩽A < A−,(ii)intermediate A−< A < A+, and (iii)front

A+< A ⩽Amax. Scaling exponents based on linear regression are given together with the

respective standard deviation. The transitional areas A−and A+are estimated from figures 3

(b), (e), (h) and (k). For simplicity, the forcing-scale area Afand the maximum vortex area

Amax are chosen as the first and last datapoints of the respective abscissae of n(A).

6. Relation of coherent structures to the inverse cascade process

Coherent features of the vorticity field ω(shown in figure 2 (b)) contain most

of its kinetic energy, as inferred from the spatial distribution from the absolute

velocity |u|in figure 4 (a), with the largest fraction located in the vicinity of the

vortex core. For the vortex pair in figure 2 (b), the energy increases towards their

separatrix where the maximum value is reached.

The f-FTLE and b-FTLE fields, Λt0+Teddy

t0and Λt0−Teddy

t0, are illustrated in figure

4 (c) and (d), respectively, where high FTLE values are potential candidates

for LCSs. The LCS method considered here detects coherent vortices by the

characteristic patterns of two-particle dispersion dynamics perpendicular to the

identified material lines that are shown in the figure. From this perspective,

the approach senses the imprint of a coherent structure on the surrounding

flow rather than detecting specific differential (OW) or amplitude markers (VM)

associated with coherence. LCS can thus be regarded as complementary to the

two other methods considered here. The LCS based on the f-FTLE field displays

more pronounced small-scale fluctuations perpendicular to the respective material

line as compared to the b-FTLE field. This reflects the different repelling and

attracting dynamics expected along forward- and backward-in-time LCSs, while

the difference between the Lagrangian scheme used for the f-FTLE and the semi-

Lagrangian approach employed for the b-FTLE (see appendix A.2) can play a

role here as well. There is a strong visual correlation between ridges in both

FTLE fields with the boundaries of vortices observed in the vorticity field. This

corroborates the above choice of sampling the vorticity distribution with the f-

FTLE and b-FTLE according to equations (2.3) and (2.4).

6.1. Cross-scale flux efficiency

When aiming at establishing a link between the nonlinear cross-scale flux and

the detected structures in configuration space, the temporal derivative of energy,

dE/dt, see figure 4 (b), is only of limited value as it does not yield clear localized

signatures correlated with coherent vortices. In contrast, the spatial distribution

of the cross-scale energy flux Zexhibits intense quadrupolar structures in coherent

regions reflecting their high level of symmetry. This is shown in figure 4 (e) and

has been observed by Xiao et al. (2009) and Liao & Ouellette (2013) as well.

Thus, coherent regions lead to large local contributions to the cross-scale flux,

but do not necessarily generate significant contributions to the spatially averaged

net inverse flux hZi.

|/|

rms

(a) |

(

)

(

)

rms

(b) dE/dt

+ T

eddy

(d) Λ

− T

eddy

Z/Z

rms

(e) Z

δθ/(π/2)

(f) δθ

-4

-2

-8

-4

0.25

0.5

0.75

Figure 4: Various physical quantities of the same region as shown in figure 2 (b). (a): Absolute

velocity |u|. (b): Energy rate dE/dt. (c): f-FTLE Λt0+Teddy

t0. (d): b-FTLE Λt0−Teddy

t0. (e):

Cross-scale kinetic energy flux Zobtained from a smooth Gaussian filter with a filtering

wavenumber of k= 60. (f): Angle between strain-rate and subgrid stress tensor δθ obtained

from the same filter as in (e).

The cross-scale energy flux is decomposed into its coherent, residual, and mixed

contributions, for the spectral cross-scale fluxes ZE

c,c,c(k),ZE

r,r,r(k),ZE

cr(k)(equa-

tion (4.4)) and the spatially averaged cross-scale flux distributions in configuration

space hZc(x)i,hZr(x)i,hZcr(x)i(equation (4.8)), respectively. Both cross-scale

flux measurements are shown in figure 5, because the ZE(k)representation is

usually favored for quantifying the existence of turbulent cascade activity, whereas

the Z(x)representation is commonly used to reveal the spatial structure of the

flux. We observe that the fluxes of all coherence definitions are qualitatively

similar, although the OW criterion in figure 5 (b), and the f-/b-FTLE in figures

5 (c) and (d), respectively, possess higher quantitative similarity in contrast

to the VM scheme. For the latter, the coherent flux ZE

c,c,c(k)is slightly higher

for all scales as shown in figure 5 (a) which shows the limits of the applied

simple gauge criterion. This is because the VM favourably extracts higher-valued

vorticity regions compared to the other detection schemes, which are attributed

to larger-sized structures according to the coherent spectrum Ec(k)in figure 3

(c). Therefore, the highest coherent flux contributions are rather in the lower

wavenumber range, with a decreasing contribution towards higher wavenumbers.

The structures detected by the f-FTLEs and b-FTLEs have vanishing coherent

flux contributions, which are close to zero for the entire inertial range. This is in

agreement with both forward- and backward-in-time LCSs, having the tendency

to collectively inhibit the energy transfer among scales (Kelley et al. 2013). The

circularly-shaped structures identified by the OW criterion also exhibit the same

behavior of inhibiting the cross-scale flux. These observations contribute to an

alternative definition of coherence in a turbulence context in the sense that energy

within these structures tends to remain rather closely at their given length scales

without cascading across scales.

Generally, the merging of coherent vortices has been one of the most appeal-

ing physical mechanisms for the inverse cascade for quite some time. However,

independent of the definition and as shown above, the coherent part of the flux

c,c,c(k)has an overall low negative contribution throughout the inertial range.

Therefore, merging effects of coherent structures have a weak influence to the

overall inverse cascade, which also has been pointed out by Xiao et al. (2009).

The residual flux ZE

r,r,r(k)remains negative throughout the inertial range as

well, peaking at small-scales, close to the forcing wavenumbers, and decaying

roughly logarithmically towards larger-scales. This is due to the contribution of

detected smaller-sized structures, according to the residual spectra Er(k)already

shown in figures 3 (c), (f), (i) and (l) above. We propose in section 6.2 that the

negative contribution of the residual flux is attributed to a stronger impact of the

thinning mechanism on smaller scales.

Furthermore, a substantial amount of the net negative flux on each length scale

originates from the mixed coherent-residual interactions ZE

cr(k), which roughly

stays at a constant level throughout the inverse flux region. We abstain from

ascribing this flux contribution to more specific physical dynamics as this flux

results from the complex interplay between coherent/residual/mixed stresses and

coherent/residual strain-rates as already mentioned in section 4.2.

Because the Fourier cross-scale flux is determined via spatial integration, an

investigation in configuration space is reasonable for gaining further insight. The

coherent part Zc(x)reconstructed from ωcmostly consists of highly ordered

quadrupolar structures and is exemplarily illustrated for the VM criterion in figure

6 (a). Compared to that, the ωr-reconstructed residual part Zr(x)in figure 6 (b)

has rather complex and unordered spatial features. Similar spatial characteristics

for coherent and residual parts are obtained for the other coherence-detection

methods. Also, qualitatively similar findings are gained for filtering lengths below

the damping-dominated and above the forcing-dominated length scales, hence we

employ a filtering wavenumber of k= 60 for the remaining analysis.

Another observation, in the context of the present work, is that coherent

structures, independent of their detection method, have much higher local cross-

scale flux amplitudes compared to their residual counterparts. Using the VM

criterion as an example, the ratio between the maximum absolute values of the ωc-

reconstructed and the ωr-reconstructed spatial cross-scale fluxes is max(|Zc|)

max(|Zr|)≈44.

Similar ratios are obtained for the other coherence definitions. However, high

amplitudes are not necessarily responsible for a high net inverse flux. As a general

observation, the probability density function (PDF) of the total cross-scale flux

Z, used for reference for the comparison to the various coherence schemes shown

in figures 7 (a)-(d), centralises around zero with a slightly negative skewness of

˜µ3=−0.93. Furthermore, the PDF has a kurtosis of ˜µ4= 336 corresponding

to a very flat distribution, which indicates that rare high negative amplitude

events contribute to the total net negative flux. The PDFs of the coherent

parts Zcof all coherence definitions exhibit larger tails to both positive and

negative values compared to the total field, with kurtosis values of 404,430,

503 and 547 for the VM, OW criterion, f-FTLE, and b-FTLE, respectively. This

suggests that coherent structures are responsible for the high spatial fluctuations

of the total cross-scale flux distribution, although these strong fluctuations do

not sustain a cascade as already seen before with the overall low Fourier cross-

-1

-0.8

-0.6

-0.4

-0.2

0.2

⟨Z(

)⟩/ε

(k)/ε

(a) VM (b) OW

2π/ℓ, k

-1

-0.8

-0.6

-0.4

-0.2

0.2

⟨Z(

)⟩/ε

, Z

(k)/ε

⟨Z(

)⟩

(k)

⟨Z

(

)⟩

c, c, c

(k)

⟨Z

(

)⟩

r, r, r

(k)

⟨Z

(

)⟩

(k)

2π/ℓ, k

(d) b-FTLE

Figure 5: Comparison of cross-scale kinetic energy fluxes of varying coherence definitions,

obtained from the spectral formulation ZE(k)(dashed lines) and the spatially averaged

formulation hZ(x)i(solid lines) obtained by Gaussian filtering for varying filter length scales `.

The grey vertical line indicates the wavenumber k= 60, at which the subsequent analysis of

the spatial cross-scale flux distribution is conducted.

scale flux contributions in figures 5 (a)-(d). In contrast, the PDFs of the residual

contributions Zrpossess lower tails to both positive and negative values. The

negative skewness values of −0.76,−1.4,−2.01, and −0.96 for the VM, OW

criterion, f-FTLE, and b-FTLE, respectively, lead to a measurable net inverse

(a) Z

Z/Zrms

(b) Z

-8

-4

Figure 6: Reconstructed spatial cross-scale kinetic energy flux distributions of the (a) purely

coherent Zc(x)and (b) purely residual Zr(x)parts obtained by the VM criterion for a

filtering wavenumber of k= 60, with each normalised by its root-mean square value for

visualisation purposes. The absolute values of Zcare much higher compared to Zrwith a ratio

max(|Zc|)

max(|Zr|)≈44.

cascade of the residual part. However, due to the overall flat nature of the residual

PDFs, the residual cascade is also driven by rare high negative amplitude events.

The angle distributions between stress and strain-rate tensors as a flux transfer

efficiency measure, motivated by equations (4.9) and (4.10), is shown in figure

4 (f), where the largest part of the marked coherent regions display angles close

to δθ =π/4. The corresponding misalignment of strain-rate and subfilter stress

tensors results in a lower nonlinear flux efficiency. A clearer picture of the overall

cascade direction tendencies is obtained from the PDFs of δθ,δθcand δθrshown

in figures 7 (e)-(h). The PDF of the total field δθ clearly shifts towards values

smaller than π/4, with skewness and kurtosis values of ˜µ3= 0.06 and ˜µ4= 3.79,

respectively. This indicates that the strain-rate and stress tensors, Sand τ, tend

to positively align, which lead to the overall net negative flux. However, the shape

of all the PDFs of the coherent parts δθcare mostly symmetric and centered at

a value of π/4with skewness values 0.06,0.02,0.0, and 0.03 for the VM, OW

criterion, f-FTLE and b-FTLE, respectively. Thus, the coherent strain-rate and

stress tensors, Scand τc,c, have the tendency to fully misalign, resulting in a lower

cross-scale flux efficiency, despite their high local flux amplitudes. On the contrary,

the PDFs of the residual parts δθrskew even more to the right compared to the

total field δθ, with skewness values of 0.14,0.13,0.12, and 0.11 for the VM, OW

criterion, f-FTLE and b-FTLE, respectively. Additionally, for the OW criterion in

figure 7 (f), for the f-FTLE in figure 7 (g), and for the b-FTLE in figure 7 (h), the

left tails of the residual parts δθrare even higher and the right tails even lower

compared to their corresponding coherent counterparts δθc. These tendencies lead

to a higher stress-strain tensor alignment behavior and thus a higher cross-scale

flux efficiency of the residual compared to the coherent part.

We conclude that although strong nonlinear cross-scale flux interactions are

found within the quadrupolar structures of figure 6 (a), these high local flux

amplitudes are not responsible for an actual cascade process. They are rather

responsible for sustaining the structures’ coherence in a turbulent environment

by keeping the associated energy close to specific length scales.

6.2. Coherent shape preservation and residual thinning mechanism

A possible reason for the small coherent Fourier cross-scale flux contributions is

the shape preservation of coherent structures in combination with the enhanced

flux efficiency of the residual part due to increased vortex thinning. To verify this

hypothesis, the MSG expanded cross-scale energy flux in equation (4.13) and its

decomposed form in equations (4.17) and (4.18) are investigated in more detail

with a filtering length of `=π/15, a geometric factor of λ= 2 and a total number

−6

−4

−2

(

)

, P

(

)

, P

(

)

(a) VM

−2

6 ⋅ 10

−2

(

δθ

)

, P

(

δθ

)

, P

(

δθ

)

(e) VM

δθ

−6

−4

−2

(

)

, P

(

)

, P

(

)

(b) OW

−2

6 ⋅ 10

−2

(

δθ

)

, P

(

δθ

)

, P

(

δθ

)

(f) OW

−6

−4

−2

(

)

, P

(

)

, P

(

)

-FTLE

−2

6 ⋅ 10

−2

(

δθ

)

, P

(

δθ

)

, P

(

δθ

)

(g) f

-FTLE

-45 -30 -15 0 15 30 45

−6

−4

−2

(

)

, P

(

)

, P

(

)

(d) b

-FTLE

0 0.25 0.5

δθ

−2

6 ⋅ 10

−2

(

δθ

)

, P

(

δθ

)

, P

(

δθ

)

(h) b

-FTLE

Figure 7: All quantities are filtered at a wavenumber of k= 60. (a)-(d): PDFs of spatial

cross-scale flux distributions Z,Zcand Zrwith each normalised by its standard deviation

ZSD,ZSD

cand ZSD

r, respectively. (e)-(h): PDFs of rotation angle distributions between

strain-rate and subgrid stress tensors δθ,δθcand δθr.

of filter bands of nb= 5. For a detailed analysis the flux fractions Q[b]on different

subfilter-scales are divided into

Q[b]=hZ[b]

SRi+hZ[b]

DSRi+hZ[b]

DSM i+hZ[b]

V GSi

hZMSG

∗i=Q[b]

SR +Q[b]

DSR +Q[b]

DSM +Q[b]

V GS,

(6.1)

where Q[b]

SR,Q[b]

DSR,Q[b]

DSM and Q[b]

V GS are the flux fractions resulting from strain

rotation (SR), differential strain rotation (DSR), differential strain magnification

(DSM) and vorticity gradient stretching (VGS), respectively. All the flux channels

are further decomposed into their coherent (c), residual (r) and mixed (cr)

contributions, e.g. the fraction of the SR is decomposed as Q[b]

SR =Q[b]

SR,c +Q[b]

SR,r +

Q[b]

SR,cr. All these contributions are illustrated in figure 8.

As already discussed by Xiao et al. (2009), the MSG flux ZMSG

∗exhibits the

locality behavior predicted by the test-field model (TFM) closure of Kraichnan

(1971). Thus, it is not surprising in figure 8 that the majority of the contributions

do not result from strongly local interactions b= 0, but rather from non-

local interactions b∈[1,3], with a decreasing influence from highly non-local

interactions b⩾4, with bbeing a measure of the scale locality of nonlinear

interactions. The first-order MSG flux only consists of the SR term Z[b]

SR, which

has no contribution to the strongly local interactions (Q[0]

SR = 0) (Xiao et al.

2009; Eyink 2006b). Thus, second-order contributions are necessary to sufficiently

capture the cross-scale flux. The differential strain magnification Z[b]

DSM has the

lowest overall contribution and is also non-existent for strongly local interactions

(Q[0]

DSM = 0) (Xiao et al. 2009; Eyink 2006b). Hence, its effect originating from

coherent and residual parts are not further analysed with respect to strain-rate

alignment properties.

In figure 8 the relative contributions of the coherent, residual and mixed parts

to each flux channel are nearly independent of the scale locality band the specific

coherence definition. The simple VM method represents a notable exception

generally yielding increased coherent flux contributions at the cost of the residual

fluxes. As already mentioned above, this is the consequence of the simple gauge

criterion chosen in section 5.1 that is not capable of eliminating all differences

between the detection methods. Test simulations (not shown) indicate that the

VM method is more robust than the LCS techniques to variations of the threshold

value with regard to the scaling of the vortex number density, but exhibits a

stronger and monotonic impact of the threshold on the relative importance of the

different MSG flux fractions than the LCS methods. The gauging procedure thus

leads to an effective relative threshold decrease for the VM method compared

to the other more complex coherence specifications, illustrating the difficulty

of neutralizing all differences between the coherence specifications. The flux

fractions mostly reflect the decomposed Fourier cross-scale flux contributions

already presented in figure 5 above. For example, the f-FTLE criterion has small

coherent cross-scale flux contributions ZE

c,c,c(k)throughout the inertial range in

figure 5 (c), and the coherent fractions Q[b]

SR,c,Q[b]

DSR,c,Q[b]

DSM,c, and Q[b]

V GS,c of the

MSG flux in figure 8 are small for all of the flux channels as well, except for the VM

method as mentioned above. We associate the lack of substantial SR and VGS of

the coherent structures with their shape preservation characteristic. In contrast,

the residual fractions Q[b]

SR,r,Q[b]

DSR,r,Q[b]

DSM,r, and Q[b]

V GS,r have a much higher

contribution in all these MSG flux channels, which is also reflected by the residual

cross-scale flux ZE

r,r,r(k). Additionally, the coherent cross-scale flux of the VM

criterion is decreasing, while the residual flux is increasing for higher wavenumbers

in the inertial range in figure 5 (a). This is also reflected by the MSG flux fractions

in figure 8. There, the coherent fractions of all MSG flux channels decrease, while

the residual fractions increase for lower scale locality. These observations are a first

sign that vortex thinning plays a minor role in coherent structures and is more

active in the residual counterpart instead. Lastly, similar to the decomposed cross-

scale fluxes in figures 5 (a)-(d), the mixed interactions Q[b]

SR,cr,Q[b]

DSR,cr,Q[b]

DSM,cr,

and Q[b]

V GS,cr contribute to almost half of the fractions in all MSG flux channels in

figure 8. However, this does not imply that coherent-residual interactions lead to

an enhanced thinning mechanism, since the MSG cross-scale flux of the mixed

contribution ZMSG

∗,cr in equation (4.17) consists of four different heterogeneous

strain-rate and stress components. We leave the analysis of these heterogeneous

interactions for future investigations.

According to the numerical studies of Xiao et al. (2009) and as described in

section 4.3, the vortex thinning mechanism is geometrically quantified by the

rotation angles δα[b]and δβ[b], which are illustrated in figures 1 (b) and (c),

respectively. Conditioning the angle δα[b]onto the band-pass filtered vorticities

ω[b]allows us to obtain statistics of the rotation angle behavior between the large

scale strain-rate S(0) and the band-pass filtered strain-rates S[b], and thus reveals

the physical mechanism of the MSG flux contribution originated from the SR

term Z[b]

SR in equation (4.13). A conditioning of the angle δα[b]onto the band-

pass filtered eddy-viscosities ν[b]

Tallows the analysis of the alignment behavior

between the large scale strain-rate S(0) and the band-pass filtered stress τ[b]=

−ν[b]

TS[b], which comes from a Newtonian stress-strain relation as already described

in section 4.3. From this we are able to infer the physical mechanism of the

DSR term Z[b]

DSR in equation (4.13), which quantifies the direction of the stress

work exertion on the strain-rate field from small to large length scales. Finally,

analysing the δβ[b]angle reveals the influence of the VGS mechanism described

01234 5

0.05

0.1

0.15

0.2

0.25

[b]

f-FTLE

b-FTLE

[

]

[

]

[

]

[

]

DSR

[

]

DSR

[

]

DSR

[

]

DSM

[

]

DSM

[

]

DSM

[

]

VGS

[

]

VGS

[

]

VGS

Figure 8: Second-order MSG flux fractions for the length scale bands in the range of b∈[0,5]

originating from SR Q[b]

SR (green), DSR Q[b]

DSR (violet), DSM Q[b]

DSM (blue), and VGS Q[b]

V GS

(orange) contributions, decoded into their corresponding coherent (subscript c), residual

(subscript r), and mixed (subscript cr) parts. The varying coherence criterions are compared

for every bin the following order: VM, OW criterion, f-FTLE and b-FTLE.

by Z[b]

V GS in equation (4.13). We extend the angle analysis by investigating the

angles δα[b]

c/r and δβ[b]

c/r of the purely coherent and residual contributions based on

the equations (4.19) - (4.24). According to figures 9 (b), (d) and (f) all coherence

definitions provide similar results regarding the above mentioned angle statistics.

Thus, the following observations and conclusions are made for coherent structures

in general, independent of the concrete definition.

The closer the angles between the strain-rate tensors δα[b]are to values of ±π/4

(their sign reflecting the sign of ω[b]), the higher the contribution of the SR term

Z[b]

SR to the inverse cascade. This rotation angle tendency is exemplarily illustrated

in figure 9 (a) for the negative vorticity condition −ω[b], for which the peak of the

PDF P(δα[b]|ω[b]<0) gradually shifts towards +π/4with decreasing scale locality,

implied by the increasing b-values (cf. Xiao et al. 2009). This gradual peak shift

quantifies SR and the accompanied transformation of circularly shaped structures

into shear layers. Thus, investigating the presence of SR is one possibility to show

the existence of the thinning process for the total field. Therefore, if we measure

the shift of the PDF peak, we are able to determine the SR effect separately for the

purely coherent and residual parts. We achieve this by evaluating the conditional

expectation values of the angles δα[b],δα[b]

cand δα[b]

ras

Ehδα[b]|(ω[b]<0)i,Ehδα[b]

c|(ω[b]

c<0)i,Ehδα[b]

r|(ω[b]

r<0)i.(6.2)

Figure 9 (b) shows that the conditional expectation values for the coher-

ent structures, independent of the concrete definition, are approximately

Eδα[b]

c|(ω[b]

c<0)≈0. This means that there is only a minor gradual shift

of the PDF peak for coherent structures with decreasing scale locality (increasing

b-values). This leads to the conclusion that SR and therefore thinning effects are

barely present within structures of the coherent part. Hence, coherent structures

are generally not turned into shear layers and tend to preserve their shape. On the

contrary, the residual part exhibits an increased expected value Eδα[b]

r|(ω[b]

r<0)

with decreasing scale locality. Therefore, the background residual field is prone

to the development of shear layers and thus has a higher thinning tendency

compared to the coherent structures in the system.

The DSR term Z[b]

DSR has an increasing contribution to the overall inverse

cascade, if the angle is δα[b]≷π/4for positive or negative eddy-viscosities ±ν[b]

respectively. Therefore, the band-pass filtered stress tensors τ[b]should exhibit

parallel or anti-parallel alignments to the large-scale strain-rate tensor S(0) based

on the eddy-viscosity’s sign in order to maximise the Z[b]

DSR contribution. The

alignment tendency reveals the energy transfer between scales during the vortex

thinning procedure. Thus, negative work exertion of the small-scale stress τ[b]on

the large-scale strain-rate field S(0) is understood as a general energy transfer

from small to large length scales enabling an overall inverse energy cascade,

whereas positive work exertion allows for a direct energy cascade. The PDF of the

angles exemplarily conditioned onto the positive eddy-viscosity P(δα[b]|ν[b]

T>0)

for the total field is presented in figure 9 (c) and shows that the angle becomes

δα[b]> π/4for b⩾3(cf. Xiao et al. 2009). This reveals the negative work

exertion of the small-scale stress τ[b]on the large-scale strain-rate tensor S(0),

as depicted in figure 1 (b), and leads to the conclusion that the vortex thinning

mechanism facilitates the energy transfer from small to large length scales for the

entire field. The conditional expectation values of the angles δα[b],δα[b]

cand δα[b]

Ehδα[b]|(ν[b]

T>0)i,Ehδα[b]

c|(ν[b]

T,c >0)i,Ehδα[b]

r|(ν[b]

T,r >0)i,(6.3)

are determined to measure the shift of the maxima towards δα[b]> π/4for

increasing b-values in figure 9 (d) and are also used to determine the alignment

tendencies of coherent and residual parts. The absolute values δα[b]instead

of the original angle δα[b]are considered, because the PDFs in figure 9 (c) are

symmetric. As a result, figure 9 (d) shows that the negative work exertion is

present within the residual parts for b⩾3as the angles become δα[b]

r>

π/4. Although Ehδα[b]

c|(ν[b]

T,c >0)iincreases for larger b, the angles for the

coherent structures remain at δα[b]

c< π/4independent of the scale locality. This

means that the coherent part even has the tendency that the coherent band-pass

filtered stress tensors τ[b]

cexert positive instead of negative work on the coherent

large-scale strain-rate tensor S(0)

c. Thus, energy within coherent structures is not

transferred from the small-scale stress to the large-scale strain-rate due to the

weak contributions of the thinning mechanism.

Lastly, the VGS term Z[b]

V GS is dependent on the angles δβ[b]between the

contractile direction of the large-scale strain-rate S(0) and the vorticity gradi-

ents ∇ω[b]. These angles approach zero δβ[b]→0for decreasing scale locality

(increasing b), as presented in figure 9 (e) (cf. Xiao et al. 2009), which shows the

contribution Z[b]

V GS for the total field. This quantifies the presence of the vorticity

isoline deformation along the stretching direction of the large-scale strain-rate

field, as depicted in figure 1 (b) for the total field. The alignment angles δβ[b]

c/r

between the purely coherent/residual vorticity gradients ∇ω[b]

c/r and the purely

coherent/residual contractile direction of the strain-rate tensor S(0)

c/r are used

to measure the effect of VGS in the coherent and residual parts of the flow

respectively. The following expected values

Ehδβ[b]i,Ehδβ[b]

ci,Ehδβ[b]

ri,(6.4)

for the total, coherent and residual fields are illustrated in figure 9 (f). The

expected value for the residual angle Eδβ[b]

rdecreases for large scale separations

(increasing b), indicating the presence of the VGS in the residual part and thus

is another indicator of the thinning mechanism. For the coherent part Eδβ[b]

cis

close to π/4for all values of b. This implies the physical picture that the coherent

large-scale strain-rate tensor S(0)

cis not distorting the vorticity isolines, which

ultimately leads to a shape preservation of coherent structures.

In conclusion, the small cascade efficiency of coherent structures is caused by

their shape preserving nature determined by the depletion of Z[b]

SR,c and Z[b]

V GS,c

terms independent of the scale separations bin combination with the positive work

exertion of the small-scale stresses on the large scale strain-rate, as captured by

the Z[b]

DSR,c term. In contrast, the higher flux efficiency in the residual parts of

the flow is caused by the enhanced thinning mechanism quantified by the Z[b]

SR,r

and Z[b]

V GS,r contributions and the negative work exertion captured in the Z[b]

DSR,r

term.

7. Conclusions

The present investigation deals with the nonlinear dynamics of coherent and resid-

ual structures in pseudospectral DNS of two-dimensional Navier-Stokes turbulence

forced at small scales. This involves the application of three commonly applied

coherence detection schemes based on the OW criterion, on the VM, and on LCSs

to identify and to isolate the respective flow components. Among these threshold-

based detection methods, the VM technique and the LCS approach are gauged by

using statistical properties of detected structures to improve comparability of the

detection results. Using this setup, (i)the performance of the employed detection

methods is discussed in relation to each other, (ii)the coherent and residual

contributions to the cross-scale energy flux of the inverse turbulent cascade are

analyzed in spectral Fourier as well as in configuration space to study their role

for the characteristics and for the dynamics of the respective parts of the flow.

We find that, i), under application of the chosen gauge criterion and in compar-

ison with the VM scheme the OW method exhibits a bias towards largest-scale

and most energetic structures. In contrast, the LCS scheme shows an increased

susceptibility for small-scale coherence as compared to the VM method. Both

tendencies can largely be neutralized by adjusting the free parameters of the VM

and the LCS methods to yield the same scale-dependency of the vortex number

density as the OW specification.

-0.5 -0.25 0 0.25 0.5

δα

[b]

/π

−2

2 ⋅ 10

−2

P(δα

[b]

|ω

[b]

< 0)

(

b=1

b=2

b=3

b=4

b=5

-0.5 -0.25 0 0.25 0.5

[b]

/π

−2

2 ⋅ 10

−2

P(δα

[b]

|ν

[b]

> 0)

(c)

0 0.2

5 0.5 0.75 1

δβ

[

/(π/

−2

8 ⋅ 10

−3

P(δ

[

)

(e)

12345

0.02

0.04

0.06



[

(δα

[b]

/π)|(ω

[b]

< 0)

]

(b)

tot

coherent

,VM

residual,VM

coherent,OW

residual,OW

coherent,f

-FTLE

residual,f-FTLE

coherent,b-FTLE

residual,b-FTLE

12345

0.15

0.2

0.25



[

(|δα

[b]

|/π)|(ν

[b]

> 0)

]

(

)

12345

0.44

0.46

0.48

0.5



[

/(π/2)

]

(f)

Figure 9: (a): PDFs of strain-rate tensor angles δα[b]conditioned onto the negative vorticities

ω[b]. (b): Expected values of δα[b],δα[b]

cand δα[b]

rconditioned onto the negative vorticities ω[b],

ω[b]

cand ω[b]

r, respectively, for varying coherence definitions. (c): PDFs of strain-rate tensor

angles δα[b]conditioned onto the positive eddy-viscosities ν[b]

T. (d): Expected values of 



δα[b]



,



δα[b]



and 



δα[b]



conditioned onto the positive eddy-viscosities ν[b]

T,ν[b]

T,c and ν[b]

T,r, respectively,

for varying coherence definitions. (e): PDFs of angles between the large-scale strain-rate tensor

and vorticity gradient vector δβ[b]. (f): Expected values of δβ[b],δβ[b]

cand δβ[b]

rfor varying

coherence definitions.

With respect to the role of detected coherent structures for turbulence dynam-

ics, (ii), we find that they are responsible for a pollution of the phenomenologically

expected spectral scaling of the kinetic energy spectrum E(k)∼k−5/3at largest

spatial scales. This finding is supported by the largely unaffected k−5/3-scaling

observed in the energy spectrum of the residual (incoherent) fraction of the turbu-

lent fluctuations. The observation suggests the possible use of coherence detection

and decomposition in DNS of homogeneous turbulence for the reduction of the

large-scale condensation of inversely cascading quantities for physical systems

that feature inverse cascades, e.g. two-dimensional Navier-Stokes turbulence or

magnetohydrodynamic turbulence.

The application of a spatial scale-filter approach for the analysis of the nonlinear

dynamics of coherent and residual parts of the turbulent flow indicates a high

nonlinear activity within coherent structures. The finding shows that coherent

structures in two-dimensional Navier-Stokes turbulence are in general dynamically

sustained while the spatial structure of the dynamics yields a shape-preserving

depletion of the nonlinear cross-scale flux with regard to the entire structure. This

is in agreement with the observed coherent Fourier cross-scale energy fluxes and

the low flux efficiency due to the high misalignment tendencies of coherent strain-

rate and subgrid stress tensors. The shape preservation of coherent structures in

this case is verified by employing the MSG expansion of the coherent spatial cross-

scale energy fluxes that exhibits a clear depletion of the deformation processes

that are scale-flux generating. These findings suggest to employ the depletion of

the MSG contributions of SR and VGS as markers for structural coherence in

two-dimensional turbulent flows.

The inverse cascade is instead driven by a combination of (i)interactions

entirely among residual fluctuations and of (ii)nonlinear interactions between

coherent structures and residual fluctuations. The former contribution is strongest

at small scales while the latter dominates at large scales. This suggests that two

different physical processes are responsible for the respective energy fluxes. For

the first contribution the dominant physical process has recently been introduced

as vortex thinning. This is in line with the enhanced alignment properties of the

residual strain-rate and subgrid stress tensors, yielding a high flux efficiency. The

second contribution is the dominant flux contribution and stays on a relatively

high and roughly constant level thoughout the inverse flux region. We abstain from

ascribing this flux contribution to more specific physical dynamics as multiple

factors may be determining its characteristics due to the heterogeneous character

of the interacting strain-rate and and stress tensor fields. Further work is presently

being pursued along these directions.

Acknowledgements

The authors thank B. Beck, R. Mäusle, J. Reiss, and J.-M. Teissier for fruitful

discussions. This work was supported by the German Research Foundation (DFG)

within the Research Training Group GRK 2433. Computing resources from the

Max Planck Computing and Data Facility (MPCDF) are also acknowledged.

Declaration of Interests. The authors report no conflict of interest.

Appendix A. Technical details

A.1. Vorticity decomposition

After thresholding according to one of the investigated coherence criteria in the

fashion of equations (2.3) and (2.4), the nonzero scalar values of the vorticity field

are clustered. All adjacent neighboring nonzero pixels (maximum of four neighbors

for each pixel) are grouped into separately countable and connected clusters.

Then, clusters whose number of pixels are below the forcing area of Af=π(`f/2)2,

with `f= 2π/kfthe forcing length scale, are excluded and discarded from the

coherent field. After that, a 5×5smooth Gaussian filter is applied to avoid

regularity issues caused by sharp boundaries. This leads to the coherent vorticity

field ωsmooth

c. The residual field is obtained by subtraction,

ωr=ω−ωsmooth

c.(A 1)

A.2. Efficient FTLE calculation

From a technical viewpoint, integrating larger amounts of Lagrangian tracers

to sufficiently resolve a turbulent flow in space and time requires substantial

computational resources. Thus, we utilize an efficient numerical technique, as

suggested by Finn & Apte (2013), to simultaneously calculate the flowmap Ft

forward and backward in time by employing a Lagrangian and a semi-Lagrangian

scheme, respectively. The forward-in-time/Lagrangian scheme advects passive

tracers for an integration time Taccording to the underlying velocity field,

which generate the forward-in-time flowmap, Ft0+T

t0. The backward-in-time/semi-

Lagrangian scheme, introduced by Leung (2011), is based on the solution of level-

set equations and constructs the backward-in-time flowmap, Ft0−T

t0, by tracking

coordinates on an Eulerian grid.

For an increased temporal resolution of the FTLE fields, a flowmap composition

method, proposed by Brunton & Rowley (2010), is applied as

Ft0+T

t0=IFt0+Nh

t0+(N−1)h◦. . . ◦ IFt0+2h

t0+h◦Ft0+h

t0,(A 2)

where ◦is the composition operator, Nthe number of substeps, h= 100 the

flowmap substep and Ithe interpolation operator. A second-order Runge-Kutta

integration scheme is applied for particle integration and a cubic interpolation

scheme, suggested by Staniforth & Côté (1991), is used for the interpolation

operator and the mapping of particles between grid points.

A.3. MSG flux derivation

The multi-gradient nature of the approach arises from a Taylor expansion

δu(b,m)(r;x) =

p=1

p!(r·∇)pu(b)(x)(A 3)

of the filtered velocity increments δu(b)(r;x) = u(b)(x+r)−u(b)(x)with sep-

aration vector r, which is subject to the filtering operation given by equation

(4.11). The filtered subgrid stress tensor τ(b)can then be entirely expressed by the

Taylor expanded velocity increments instead of the original increments δu(b)(r).

This yields (cf. Eyink 2006a)

τ(b,m)=ZGl(r)δu(b,m)(r)δu(b,m)(r)d2r

−ZGl(r)δu(b,m)(r)d2rZGl(r)δu(b,m)(r)d2r,(A 4)

where τ(b,m)is a multi-scale and multi-gradient expression for the stress.

It can be shown that the MSG expanded stress converges as limm→∞ τ(b,m)=

τ(b)in the L1-norm. Nevertheless, for increasing scale indices b, which corre-

sponds to adding finer-scale structures, a growing amount of space gradients of

higher-order, m, is required to approximate τ(b,m)≈τ(b). Therefore, a coherent-

subregions approximation (CSA) approach is suggested by Eyink (2006a) enabling

the approximation of the MSG stress by low-order gradients m. As a result, the

CSA corrected MSG stress is obtained, which is more accurately approximated as

τ(b,m)

∗≈τ(b,m)with fewer gradients m. According to Eyink (2006b) and Xiao et al.

(2009), the best approximation with regard to the original cross-scale flux term

is reached for expansions up to the second-order in gradients m= 2. For the sake

of brevity we define τMSG

∗=τ(b,2)

∗and ZMSG

∗=Z(b,2)

∗, which are the MSG-CSA

stress and cross-scale flux, respectively, expanded to second order, giving

ZMSG

∗=−S(0) :τMSG

∗

b=0 "−C[b]

2`2

bω[b](S(0) :˜

S[b])

| {z }

Z[b]

+C[b]

8`4

ν[b]

z }| {

∇ω[b]·∇α[b])(S(0) :S[b])

| {z }

Z[b]

DSR

−C[b]

16 `4

γ[b]

z }| {

∇ω[b]·∇ln σ[b])(S(0) :˜

S[b])

| {z }

Z[b]

DSM

+C[b]

32 `4

b(∇ω[b])TS(0)(∇ω[b])

| {z }

Z[b]

V GS

−(∇ψ(nb)

∗)TS(0)(∇ψ(nb)

∗)

| {z }

Z(nb)

F SF

,(A 5)

with the fluctuation stream function ψ(nb)

∗(cf. Eyink 2006b) and the coefficients

C[b]

p(cf. appendix C in Eyink 2006a). The flux contribution from the fluctua-

tion stream function (FSF) Z(nb)

F SF is similarly interpreted as a vorticity gradient

stretching but considered much smaller in magnitude due to cancellations from

spatial averaging. Additionally, it possesses a positive mean as shown by Xiao et al.

(2009), thus a further analysis with regard to the inverse cascade mechanism is

neglected for this term.

Appendix B. Sensitivity of energy injection rate, KLB theory and

finite-size effects

The present coherent structure analysis is based on DNS of the two-dimensional

Navier-Stokes equations. In this section, we briefly discuss the choice of system

parameters, with regard to the dynamics of structure formation and turbulence

statistics.

In contrast to other studies (see Borue 1994; Babiano et al. 1987; Maltrud

& Vallis 1993; Danilov & Gurarie 2001; Vallgren 2011; Burgess & Scott 2018),

the present work does not employ a hyperviscous dissipation term to avoid the

accompanying unphysical distortion of the spatial structure of the flow. Although

this shortens the scaling-range of the enstrophy cascade, the resulting spectral

extent still suffices for our purposes, as shown by the spectra and enstrophy fluxes

in figures 2 (c) and (d), and figure 10. We have found no further implications for

other diagnostics relevant in the context of the present investigation.

The formation of structures is dependent on the energy injection rate Iand

therefore we vary its rate for fixed viscosity νand large-scale friction values

(k0,ω,σ2

ω,αω), which also affects the ratio of energy dissipated at largest scales

with rate αand at viscous scales with rate ν. Due to the unavoidably limited

spectral bandwidth of the numerical simulations, it is not possible to fulfill both

requirements of the KLB picture at the same time, namely a ratio of α/I= 1,

such that all the injected energy is dissipated at large-scales, in combination

with an exact power-law ∼k−5/3for the energy scaling-range. Even DNS at

significantly higher numerical resolution (see e.g. Boffetta & Musacchio 2010)

exhibit smaller deviations from the asymptotic scaling exponent. We observe, that

simultaneously only one of the two characteristics is approximately realisable with

sufficient accuracy. Therefore, three simulation runs are performed, as listed in

table 2.

Next to numerical resolution and the large-scale damping required for quasi-

Run Iα/Iν/IRe[105]Teddy

run1 13.05 0.67 0.33 7.87 0.1

run2 3.5 0.6 0.4 4.36 0.2

run3 0.88 0.47 0.53 2.34 0.4

Table 2: Simulation parameters of the two-dimensional hydrodynamic turbulence system with

a resolution of 40962, forcing wavenumber kf= 200, large-scale friction factor αω= 500, center

of Gaussian damping profile k0,ω = 0.1, variance of the Gaussian σ2

ω= 1, viscosity ν= 5 ·10−5

and timestep dt = 5 ·10−6.

110

−8

−6

−4

−2

E(k)

−5/3

−3

run2

run3

110

−3

−1

E(k)k

5/3

−2/3

110

-1

-0.8

-0.6

-0.4

-0.2

0.2

(k)/

0.2

0.4

0.6

0.8

Ω

(k)/

Figure 10: Kinetic energy spectra E(k)of run2 and 3 with varying kinetic energy injection

rates (left) with inset showing the compensated spectra E(k)k5/3−2/3

I, where the black dashed

line indicates a value of CE= 6.69 predicted by the TFM closure of Kraichnan (1971).

Normalised cross-scale kinetic energy ZE(k)/I(top right) and enstrophy fluxes ZΩ(k)/ηI

(bottom right) from run2 and 3, respectively.

stationarity of the flow, the large-scale driving, i.e. the energy injection rate I,

exerts an important influence on the system.

According to figure 2 (c), run1 with the highest energy injection rate exhibits

the strongest deviation from the k−5/3scaling, but the best developed normalised

cross-scale energy ZE(k)/Iand enstrophy fluxes ZΩ(k)/ηI(figure 2 (d)) close

to values of −1and 1in the energy and enstrophy inertial ranges, respectively.

Run3 with the lowest energy injection rate in figure 10 is the closest to fulfill

the power-law but displays weaker nonlinear fluxes. Although claims exist that

relate the large-scale steepening of the energy spectrum to the large-scale damping

leading to the formation of large-scale coherent vortices at largest scales (see Borue

1994; Danilov & Gurarie 2001), other studies show that vortex formation already

occurs on smaller scales and is not entirely an artifact due to hypofriction effects

(see Babiano et al. 1987; Vallgren 2011; Burgess & Scott 2017). Therefore, it

appears to be plausible that the formation of coherent structures in configuration

space is an inherent property of two-dimensional turbulence not captured by

the wavenumber-based KLB phenomenology. This structure formation property

has the tendency to pollute the scaling of the energy spectrum, which is further

discussed in section 5.

With decreasing energy injection rate the vorticity field has less distinct coher-

ent features as shown in figure 2 (b) and figures 11 (a) and (b), where single vortex

structures become less intense and gradually dissolve into the residual background,

similarly observed by Burgess & Scott (2017). At the same time the probability

density function (PDF) of the vorticity becomes flatter with increasing injection

rates according to figure 11 (c), where the increasing values at the tails of the

PDFs are associated to the large vorticity values in the visible vortex cores. For

the present analysis run1 is used, due to the most visible presence of coherent

structures and the stronger spectral cross-scale flux in the inverse cascade regime

compared to the other DNS. However, the remaining DNS, run2 and run3, also

lead to qualitatively similar results. This suggests that the similarity scaling is not

the determining factor for the inverse cascade dynamics of coherent structures.

The large-scale damping certainly has a strong and intended effect on the large-

scale energetics that naturally extends over a limited spectral range into the

smaller-scale dynamics. The high level of isotropy of the damping mechanism,

however, is an effort to avoid additional severe nuisances such as violent and

random artificial straining or other unwanted anisotropic processes. Although

comparison with similar works without applied large-scale damping suggests that

the damping does not lead to even more severe artefacts with regard to coherence

dynamics as already inflicted by the finite-size of the system. Such unwanted side-

effects in particular with regard to the higher-order MSG results can of course not

be fully ruled out. Furthermore, the impact of numerical resolution, as observed

in test simulations with a monotonically decreasing number of collocation points

down to 2562reveals that the three-regime signature of the vortex number density

already becomes less discernible at a resolution of 20482. The qualitative results

obtained via the MSG expansion, in contrast, remain unchanged and observable

for all tested resolutions.

REFERENCES

Babiano, A., Basdevant, C., Legras, B. & Sadourny, R. 1987 Vorticity and passive-scalar

dynamics in two-dimensional turbulence. Journal of Fluid Mechanics 183, 379–397.

Batchelor, G. K. 1969 Computation of the energy spectrum in homogeneous two-dimensional

turbulence. The Physics of Fluids 12, 233–239.

Benzi, R., Paladin, G., Patarnello, S., Santangelo, P. & Vulpiani, A. 1986

Intermittency and coherent structures in two-dimensional turbulence. Journal of Physics

A: Mathematical and General 19, 3771–3784.

(a) run2

/max(

)

(b) run3

-10 -5 0 5 10

ω/ω

rms

−6

−5

−4

−3

−2

−1

P(ω/ω

rms

)

(c)

run1: 

:7.

run2: 

:5.51

run3: 

:4.16

-1

-0.75

-0.5

-0.25

0.25

0.5

0.75

Figure 11: (a), (b): Vorticities of run2 and 3 from a 10242region, which is 6.25% of the total

physical domain. The colorbar is normalised to ω/ max |ω|for a better visual comparison of

structures from varying injection rates. (c): The kurtosis of the vorticity PDFs from runs1-3

are denoted as ˜µ4.

Benzi, R., Patarnello, S. & Santangelo, P. 1988 Self-similar coherent structures in two-

dimensional decaying turbulence. Journal of Physics A: Mathematical and General 21,

1221–1237.

Boffetta, G. & Musacchio, S. 2010 Evidence for the double cascade scenario in two-

dimensional turbulence. Physical Review E 82, 016307.

Borue, V. 1994 Inverse energy cascade in stationary two-dimensional homogeneous turbulence.

Physical Review Letters 72 (10), 1475–1478.

Brunton, S. L. & Rowley, C. W. 2010 Fast computation of finite-time Lyapunov exponent

fields for unsteady flows. Chaos 20, 017503.

Burgess, B. H. & Scott, R. K. 2017 Scaling theory for vortices in the two-dimensional

inverse energy cascade. Journal of Fluid Mechanics 811, 742–756.

Burgess, B. H. & Scott, R. K. 2018 Robustness of vortex populations in the two-dimensional

inverse energy cascade. Journal of Fluid Mechanics 850, 844–874.

Canivete Cuissa, J. R. & Steiner, O. 2020 Vortices evolution in the solar atmosphere: A

dynamical equation for the swirling strength. Astronomy and Astrophysics 639, A118.

Canuto, C., Hussaini, M. Y., Quarteroni, A. & Zang, T. A. 1988 Spectral methods in

fluid dynamics. Springer-Verlag Berlin Heidelberg.

Chakraborty, P., Balachandar, S. & Adrian, R. J. 2005 On the relationships between

local vortex identification schemes. Journal of Fluid Mechanics 535, 189–214.

Chen, S., Ecke, R. E., Eyink, G. L., Rivera, M., Wan, M. & Xiao, Z. 2006 Physical

mechanism of the two-dimensional inverse energy cascade. Physical Review Letters 96,

084502.

Chong, M. S., Perry, A. E. & Cantwell, B. J. 1990 A general classification of three-

dimensional flow fields. Physics of Fluids A: Fluid Dynamics 2(5), 765–777.

Danilov, S. & Gurarie, D. 2001 Forced two-dimensional turbulence in spectral and physical

space. Physical Review E 63 (6), 061208.

Eyink, G. L. 2006aMulti-scale gradient expansion of the turbulent stress tensor. Journal of

Fluid Mechanics 549, 159–190.

Eyink, G. L. 2006bA turbulent constitutive law for the two-dimensional inverse energy cascade.

Journal of Fluid Mechanics 549, 191–214.

Fang, L. & Ouellette, N. T. 2016 Advection and the efficiency of spectral energy transfer

in two-dimensional turbulence. Physical Review Letters 117 (10), 104501.

Farge, M. & Schneider, K. 2015 Wavelet transforms and their applications to MHD and

plasma turbulence: a review. Journal of Plasma Physics 81.

Finn, J. & Apte, S. V. 2013 Integrated computation of finite-time Lyapunov exponent fields

during direct numerical simulation of unsteady flows. Chaos 23, 013145.

Frisch, U & Sulem, P. L. 1984 Numerical simulation of the inverse cascade in two-dimensional

turbulence. The Physics of Fluids 27, 1921–1923.

Hadjighasem, A., Farazmand, A., Blazevski, D., Froyland, G. & Haller, G. 2017 A

critical comparison of Lagrangian methods for coherent structure detection. Chaos 27,

053104.

Haller, G. 2005 An objective definition of a vortex. Journal of Fluid Mechanics 525, 1–26.

Haller, G. 2015 Lagrangian coherent structures. Annual Review of Fluid Mechanics 47, 137–

162.

Haller, G., Hadjighasem, A., Farazmand, M. & Huhn, F. 2016 Defining coherent vortices

objectively from the vorticity. Journal of Fluid Mechanics 795, 136–173.

Haller, G. & Yuan, G. 2000 Lagrangian coherent structures and mixing in two-dimensional

turbulence. Physica D: Nonlinear Phenomena 147, 352–370.

Holmes, P., Lumley, J. L., Berkooz, G. & Rowley, C. W. 2012 Turbulence, coherent

structures, dynamical systems, and symmetry, 2nd edn. Cambridge University Press.

Hua & Klein 1998 An exact criterion for the stirring properties of nearly two-dimensional

turbulence. Physica D 113, 98–110.

Hunt, J. C. R., Wray, A. A. & Moin, P. 1988 Eddies, streams and convergence zones in

turbulent flows. Proceedings of the Summer Program 1988 pp. 193–208.

Jeong, J. & Hussain, F. 1995 On the identification of a vortex. Journal of Fluid Mechanics

285, 69–94.

Kelley, D. H., Allshouse, M. R. & Ouellette, N. T. 2013 Lagrangian coherent structures

separate dynamically distinct regions in fluid flows. Physical Review E 88, 013017.

Kraichnan, R. H. 1967 Inertial ranges in two-dimensional turbulence. The Physics of Fluids

10, 1417–1423.

Kraichnan, R. H. 1971 Inertial-range transfer in two- and three-dimensional turbulence.

Journal of Fluid Mechanics 47, 525–535.

Leith, C. E. 1968 Diffusion approximation for two-dimensional turbulence. The Physics of

Fluids 11, 671–673.

Leung, S. 2011 An Eulerian approach for computing the finite time Lyapunov exponent. Journal

of Computational Physics 230, 3500–3524.

Liao, Y. & Ouellette, N. T. 2013 Spatial structure of spectral transport in two-dimensional

flow. Journal of Fluid Mechanics 725, 281–298.

Lilly, D. K. 1971 Numerical simulation of developing and decaying two-dimensional turbulence.

Journal of Fluid Mechanics 45, 395–415.

Maltrud, M. E. & Vallis, G. K. 1993 Energy and enstrophy transfer in numerical simulations

of two-dimensional turbulence. Physics of Fluids A: Fluid Dynamics 5(7), 1760–1775.

Ohkitani, K. 1991 Wave number space dynamics of enstrophy cascade in a forced two-

dimensional turbulence. Physics of Fluids A: Fluid Dynamics 3, 1598–1611.

Okamoto, N., Yoshimatsu, K., Schneider, K., Farge, M. & Kaneda, Y. 2007

Coherent vortices in high resolution direct numerical simulation of homogeneous isotropic

turbulence: a wavelet viewpoint. Physics of Fluids 19, 115109.

Okubo, A. 1970 Horizontal dispersion of floatable particles in the vicinity of velocity

singularities such as convergences. Deep Sea Research and Oceanographic Abstracts 17 (3),

445–454.

Ouellette, N. 2012 On the dynamical role of coherent structures in turbulence. Comptes

Rendus Physique 13, 866–877.

Paret, J & Tabeling, P. 1998 Intermittency in the two-dimensional inverse cascade of energy:

experimental observations. Physics of Fluids 10, 3126–3136.

Pope, S. B. 2000 Turbulent flows. Cambridge University Press.

Rowley, C. W., Mezić, I., Bagheri, S., Schlatter, P. & Henningson, D. S. 2009

Spectral analysis of nonlinear flows. Journal of Fluid Mechanics 641, 115–127.

Rutgers, M. A. 1998 Forced 2d turbulence: experimental evidence of simultaneous inverse

energy and forward enstrophy cascades. Physical Review Letters 81, 2244–2247.

Schmid, P. J. 2010 Dynamic mode decomposition of numerical and experimental data. Journal

of Fluid Mechanics 656, 5–28.

Scott, R. K. 2007 Nonrobustness of the two-dimensional turbulent inverse cascade. Physical

Review E 75, 046301.

Staniforth, A. & Côté, J. 1991 Semi-Lagrangian integration schemes and their application

to environmental flows. Monthly Review Weather 119 (9), 2206–2223.

Taira, K., Brunton, S. L., Dawson, S. T. M., Rowley, C. W., Colonius, T., McKeon,

B. J., Schmidt, O. T., Gordeyev, S., Theofilis, V. & Ukeiley, L. S. 2017 Modal

analysis of fluid flows: an overview. AIAA Journal 55, 4013–4041.

Towne, A., Schmidt, O. T. & Colonius, T. 2018 Spectral proper orthogonal decomposition

and its relationship to dynamic mode decomposition and resolvent analysis. Journal of

Fluid Mechanics 847, 821–867.

Vallgren, A. 2011 Infrared Reynolds number dependency of the two-dimensional inverse

energy cascade. Journal of Fluid Mechanics 667, 463–473.

Weiss, J. 1991 The dynamics of enstrophy transfer in two-dimensional hydrodynamics. Physica

D48, 273–294.

Xiao, Z., Wan, M., Chen, S. & Eyink, G. L. 2009 Physical mechanism of the inverse

energy cascade of two-dimensional turbulence: a numerical investigation. Journal of Fluid

Mechanics 619, 1–44.

Yadav, N., Cameron, R. H. & Solanki, S. K. 2021 Vortex flow properties in simulations

of solar plage region: evidence for their role in chromospheric heating. Astronomy and

Astrophysics 645, A3.

Yoshimatsu, K., Kondo, Y., Schneider, K., Okamoto, N., Hagiwara, H. & Farge, M.

2009 Wavelet-based coherent vorticity sheet and current sheet extraction from three-

dimensional homogeneous magnetohydrodynamic turbulence. Physics of Plasmas 16,

082306.

Zhou, J., Adrian, R. J., Balachandar, S. & Kendall, T. M. 1999 Mechanisms for

generating coherent packets of hairpin vortices in channel flow. Journal of Fluid Mechanics

387, 353–396.