Probing universal dynamics with topological data analysis in a gluonic plasma

We study nonequilibrium dynamics of SU(2) lattice gauge theory in Minkowski space-time in a classical-statistical regime, where characteristic gluon occupancies are much larger than unity. In this strongly correlated system far from equilibrium, the correlations of energy and topological densities show self-similar behavior related to a turbulent cascade towards higher momentum scales. We employ persistent homology to infer topological features of the gluonic plasma via a hierarchy of simplicial and cubical complexes. All topological observables under investigation are also manifestly gauge invariant and are shown to exhibit self-similar evolution, which relate the spatial and temporal properties of the plasma in terms of universal scaling exponents and functions. The findings may help to understand the early stages of heavy-ion collisions in the limit of high energies, and our methods can also facilitate the topological analysis of other complex systems such as encountered in experiments with ultracold quantum gases.


I. INTRODUCTION
Collective phenomena are ubiquitous during the thermalization of the quark-gluon plasma (QGP) generated in ultrarelativistic collisions of heavy nuclei at the Relativistic Heavy Ion Collider (RHIC) or the Large Hadron Collider (LHC) [1][2][3].Most of the freed quarks and gluons during the first yoctosecond of the collision originate from the small-x regime, giving rise to the notion of gluon saturation and large gluon occupations [4,5].The 'bottom-up' scenario provides a systematic description for the subsequent thermalization of the pre-equilibrium QGP [6].Its early stage is governed by gluon overoccupation at low momenta, such that the weakly coupled system can be accurately mapped onto classicalstatistical Yang-Mills theory [7][8][9][10][11][12], which sets the model for this study.For simplicity we focus on homogeneous initial states in non-expanding Minkowski geometry, while more realistic extensions feature also longitudinal expansion [13].
In this work we compute the nonequilibrium time evolution of an initially over-occupied gluonic plasma.For a real-time SU(2) lattice gauge theory in 3+1 dimensions in the classical-statistical regime, we analyze gaugeinvariant observables using energy and topological densities and their fluctuations, which evade the typical drawbacks of using gauge-fixed field correlation functions.In order to infer non-local geometric information about the plasma's evolution, we use persistent homology [14][15][16].Persistent homology computes topological structures that appear in a hierarchy of simplicial or cubical complexes constructed from the field data and features measures of their dominance (persistence).It behaves stable with respect to perturbations on the data, and can be efficiently calculated [14][15][16].We devise gauge-invariant persistent homology observables and demonstrate that they exhibit self-similar evolution in terms of universal scaling exponents and functions.Identifying universal phenomena in the time evolution of the QGP off equilibrium provides a rich line of understanding the otherwise complicated and theoretically challenging relaxation dynamics of strongly interacting nuclei [3].For the gluonic plasma these include universal time dependences related to nonthermal fixed points [17][18][19], and the universal approach to local thermal equilibrium governed by viscous hydrodynamics [20][21][22].Nonthermal fixed points provide far-from-equilibrium attractor solutions characterized by universal dynamical selfsimilarity.They have first been discussed for relativistic and non-relativistic scalar theories, partly in the context of reheating of the early universe after inflation [23][24][25][26][27].In recent years, dynamical self-similarity has been found experimentally for ultracold Bose gases [28][29][30], which includes evidence for its universality across very different initial conditions.
For the gluonic plasma a dynamically separating hierarchy of momentum scales characterizes transport processes related to nonthermal fixed points.In the infrared, the string tension scale encodes the area dependence of spatial Wilson loops [31,32], going along with condensation far from equilibrium [33].The Debye mass describing the electric screening scale evolves dynamically towards the infrared, too [12,31], and the hard scale indicates energy transport towards the ultraviolet in a universal self-similar process [12,18,31,[34][35][36][37][38][39][40].
Such universal behavior is typically probed using variants of occupation numbers.These are constructed as two-point correlation functions of gluon fields, which are not gauge invariant and contain limited information on the plasma.Therefore, it is instrumental to go beyond this description using gauge invariant observables.Energy densities and their fluctuations are relevant quantities for the hydrodynamical description of the later fluid dynamics of the QGP in heavy ion collisions [41,42].Topological densities are of importance for the structure of the vacuum of quantum chromodynamics (QCD), linked to chiral properties of quarks in the QGP via anomalies and of central relevance for topological effects such as sphaleron transitions [31].For these reasons the study of two-point correlation functions of energy and topological (charge) densities is promising, which involve gauge invariant four-point correlation functions in the gluon fields.
This paper is organized as follows.In Section II we provide details on the lattice setup and discuss energy and topological density correlations, which show clear manifestations of the direct cascade related to energy transport.We argue that dynamical scaling exponents can be understood from occupation numbers and the energymomentum conservation Ward identity.We suggest an explanation for the accurate matching of energy and topological density correlations.In Section III we introduce persistent homology, and discuss the self-similarity of Betti number distributions.We explain how this can heuristically be understood from the behavior of correlation functions in conjunction with the bounded packing of topological features in the constant lattice volume.Finally, in Section IV we conclude and provide an outlook.In the Appendixes we provide details on the lattice setup and occupation numbers, and discuss further theoretical aspects and more of our results for persistent homology.

II. DIRECT CASCADE IN ENERGY AND TOPOLOGICAL DENSITY CORRELATORS
Energy and topological density correlators show dynamical self-similarity indicative of an energy cascade.In this section we first describe the lattice simulations and the examined observables and proceed with an analysis of energy and topological density correlators.

A. Real-time non-Abelian gauge theory on the lattice
We consider SU(N c ) gauge theory for N c = 2 on a cubic lattice with Minkowski metric, spatial extent N 3 s , and temporal and spatial lattice spacings dt and a s , respectively.Spatially periodic boundary conditions are applied.The volume of the spatial lattice Λ s is V = a 3 s N 3 s .We study the theory in temporal-axial gauge, A 0 ≡ 0, and formulate it in terms of su(N c )-valued (chromo-)electric fields E i (t, x) and SU(N c )-valued link variables U i (t, x), i = 1, 2, 3.The latter are related to the gauge fields A i (t, x) via U i (t, x) ≈ exp(iga s A i (t, x)) up to lattice corrections, with g the gauge coupling.
We study the system in the weakly coupled regime, g ≪ 1, and consider highly occupied gluonic initial conditions.Then, at sufficiently early times the full quantum dynamics is accurately described by the classicalstatistical approximation [7][8][9][10][11][12]26].Specifically, we sample over Gaussian initial conditions with variances where averages over transverse polarizations and color degrees are implied.Here the momentum scale Q determines the width of the initial momentum-space distribution.In particular, these correlation functions are related to the distribution function of gluonic quasiparticles as f (t=0, p) ≃ |p|⟨AA⟩(t=0, p) ≃ ⟨EE⟩(t=0, p)/|p|, see Appendix B. Such gluon over-occupation up to a gluon saturation scale is characteristic for a state shortly after the collision of heavy ions [3].More details on how initial conditions are generated can be found in [18,39].The time evolution for the lattice fields E i (t, x) and U i (t, x) results from solving the classical equations of motion.These are the Hamilton equations for the instantaneous classical lattice Hamiltonian =: with elementary spatial plaquette variables and ĵ, k unit lattice vectors in directions j, k.We refer to T 00 (t, x) as the energy density, and emphasize its gaugeinvariance.Specifically, the discretized lattice equations of motion derived from the classical Hamiltonian (2) read [39] U where The anti-Hermitian part of a matrix U is defined as Additionally, Gauss' law has to be satisfied.Once fixed initially, which is accomplished using a projection algorithm onto the constraint surface [56], it is preserved by the time evolution of Equation (5).Finally, observables are computed as the average over sampled initial conditions.We set N s = 512, dt/a s = 0.05 and Qa s = 0.125, and explicitly verified insensitivity of our results to variations of the lattice spacing, specifically upon comparison with results for Qa s = 0.25 with the same lattice volume. 1 Accordingly, no renormalization is applied to any of the computed observables.Results are shown for a single run, which for our large lattice agree with results averaged over multiple runs up to statistical fluctuations.Rescaling the fields A i → gA i and E i → gE i , the equations of motion and initial conditions become independent of the coupling g.This formally corresponds to the classical field limit g → 0 while keeping g 2 f fixed.For more details on the lattice formulation we refer to [39].Based on [32,59,60], we expect similar results as in this work for SU(N c ) theories with N c ≥ 2.
In conjunction with energy densities, we study the dynamics of topological densities q(t, x), which are spacetime integrands of Chern-Simons numbers.In the continuum, we have q ∼ Tr(E • B), with the (chromo-)magnetic field B. On the lattice we define a topological density as where E av (t, x) and B av (t, x) are SU(N c )-valued cloverleaf variants of lattice electric and magnetic fields, averaging contributions of neighboring lattice sites as detailed in Appendix A. We emphasize that q(t, x) defined through Equation ( 8) is gauge-invariant.

B. Correlations of energy and topological densities
Connected two-point correlation functions of the energy density T 00 (t, x) are constructed via with the volume-averaged energy density which remains constant in time in the simulations as required by total energy conservation.We Fouriertransform to lattice momentum space in x-direction, C T 00 (t, px ) = Ns−1 ∆x=0 C(t, ∆x=(∆x, 0, 0)) e −i px∆x , (11) with px ∈ 2π/(a s N s ){−N s /2, . . ., N s /2 − 1} lattice momentum, corresponding to physical momentum p x = (2/a s ) sin(p x a s /2).Finally, we set In continuum this correlator corresponds to where T 00 (t, p) = x T 00 (t, x) exp(−ipx) with the notation x = V d 3 x and pi = dp i /(2π).The topological density two-point correlation function ⟨(q(t, p x )) 2 ⟩ c is defined analogously.
In Figure 1 we show results for energy and topological density correlators, compared across times.We find that they start off flat in momentum space and rapidly decay above the momentum scale Q.This indicates an equidistribution of low momentum fluctuations of energy and topological densities in momentum space.After a short early-time dynamics (not displayed), the correlators remain constant in shape for times above Qt = 100 and solely decrease in overall numbers while shifting to larger momenta.This is indicative of a nonthermal fixed point.We test for the characteristic dynamical self-similar scaling according to the scaling ansatz for suitably chosen scaling exponents α T 00 , β and a reference time t ′ in the time interval when Equation (14) describes the dynamics accurately.In the right column of Figure 1 we display correlations rescaled according to the scaling ansatz (14), where we suggestively set β = −1/7 and α T 00 = −1/7.We find accurate matching of the rescaled correlations for the entire considered time range from Qt = 100 to Qt = 1800.The same observations apply to topological density correlations.In fact, both types of correlations agree with each other nearly exactly.We argue for this at the end of this subsection, and primarily discuss T 00 correlations in the following, noting that similar arguments apply for topological density correlations.
One can heuristically motivate the scaling ansatz (14) by following a kinetic quasi-particle picture.The Fouriertransformed local energy density fluctuations for sample i can be approximated on an individual sample level as Here ω p ≃ |p| is the gluon dispersion and f i (t, p) denotes a distribution function of gluon occupation numbers computed for a single classical-statistical sample.It has previously been demonstrated that the gluon occupancy ⟨f i (t, p)⟩ shows dynamical self-similarity for initial conditions similar to ours, scaling as [12,18,[37][38][39]] even for single samples as for our simulations, see Appendix B. The scaling is reminiscent of a direct energy cascade to higher momenta since the momenta of hard modes grow as (t/t ′ ) 1/7 Q while the energy density remains constant during this process, ⟨ p T 00 i (t, p)⟩ = const.Due to Equation (15) one expects that the exponent β = −1/7 governs the dynamics at hard momenta also for ⟨(T 00 (t, p x )) 2 ⟩ c .
The scaling exponent α T 00 can be understood from the energy-momentum conservation Ward identity.We work in a continuum space-time with spatial volume V in a general spatial dimension d.In [61] via metric variations the identity has been derived for correlations of the stress-energy tensor T µν (x), connected in the T µν (x i ).Integrating x 1 over V , we find using Gauss' theorem and periodic boundary conditions, such that no spatial boundary terms occur, Expectation values in our simulations exhibit approxi-mate homogeneity and isotropy in space, such that where ∆x = x − y and ∆y = y − x.We note that Analogous to Equation ( 14) we assume a scaling ansatz in position space for T 00 -correlations, Together with classical-statistical simulations corresponding to symmetric operator ordering, denoted with the substitution ∆x ′ = (t/t ′ ) −β ∆x and exploiting the empty boundary of the lattice torus due to periodic boundary conditions.This is solved by α ′ = −dβ.We obtain the scaling behavior of ⟨(T 00 (t, p x ))2 ⟩ c via momentum integration over directions 2, . . ., d and a Fourier transformation: such that α T 00 = β, which is consistent with Figure 1.It remains to discuss the similarity of energy density T 00 ∼ Tr(E 2 + B 2 ) and topological density q ∼ Tr(E • B) correlations visible in Figure 1 and giving rise to identical dynamical scaling behavior.While their one-point functions are different, T 00 = (0.106 ± 8.4 the error indicates fluctuations in time), their variances agree.This suggests that the colour and spatial directions of E(t, x) and B(t, x) are statistically independent and, by spatial homogeneity, identically distributed (i.i.d.).This is consistent with the observed ⟨Tr(E 2 + B 2 )⟩ > 0, ⟨Tr(E • B)⟩ = 0, and identical higher-order connected correlation functions.Phrased differently, the space-time and colour components of the field strength tensor F µν (t, x) are i.i.d. up to antisymmetry of the indices.Similar observations have been reported for two-point correlation functions of electric and magnetic fields [37].Note that this indicates suppressed contributions from non-trivial topological excitations such as sphalerons [31].
The agreement of energy and topological densities provides further evidence for the universality of the nonthermal fixed point across different observables O with respect to its scaling exponent β.The amplitude exponent α O is linked to β by constraints such as energymomentum conservation.In the next section we provide a complementary study of this universality by analyzing non-local geometric observables instead of correlation functions.

III. DIRECT CASCADE IN PERSISTENT HOMOLOGY
Persistent homology provides a quantitative means to extract topological structures from point clouds of data, and is based on the construction of a hierarchy of combinatorial objects like simplicial or cubical complexes.We will start this section with an intuitive introduction to the utilized alpha complexes and their persistent homology in Section III A, followed by a discussion of dynamical self-similarity in persistent homology in Section III B. We then study the dynamics of topological structures in point clouds computed from energy and topological densities in Section III C. We find clear indications for selfsimilar scaling associated to the energy cascade, which again reveals universality across observables.Moreover, we confirm a geometric relation between appearing exponents in persistent homology.In Section III D we discuss persistence ratios of topological structures, i.e., spatially scale-invariant measures of their dominance, which form a set of dynamical invariants beyond self-similarity in the gluonic plasma.
Similarly to the findings for alpha complexes reported in this section, in Appendix F we discuss so-called cubical complexes, which give access to the persistent homology of level sets without spatial metric information entering the analysis, and are particularly suitable for pixelized data.In this approach the application of coarse-graining is inevitable to exclude lattice artefacts.Again, we find evidence for self-similarity associated to the energy cascade.Hence, the universal dynamics also encompasses the persistent homology of cubical complexes.Since nonthermal fixed points are a scale-dependent phenomenon by means of self-similar transport across length scales, we focus on the results for alpha complexes here, which seem to reveal self-similar scaling for a wider time interval than cubical complexes.

A. Introducing alpha complexes and persistent homology
Let X ⊂ R 3 be a point cloud, i.e., a finite set of distinct points. 2 Alpha complexes of X are computationally useful geometric simplicial complexes describing X.Their topology can be extracted algorithmically, which is facilitated by their combinatorial properties.Here we focus on an intuitive description, and in Appendix C we deliver more mathematical details related to persistent homology.For a thorough introduction to topological data analysis we refer to the literature [14][15][16].
Simplices of different dimensions are illustrated in Figure 2(a).We refer to a point as a 0-simplex, a line as a 1-simplex, a triangle as a 2-simplex, and a filled tetrahedron as a 3-simplex.A simplicial complex is a collection of simplices closed under taking boundaries.We construct alpha complexes of X. Trivially, we can connect two points in X by a straight line, i.e., a 1-simplex.Similarly, through any three points in X we can draw a unique circle and through any four points in X we can draw a unique sphere.We call a circle or sphere empty if all points of X lie on or outside of it.We can associate to any two, three and four points a 1-, 2-or 3-simplex connecting them, and a radius parameter.The latter can be defined as half of the length of the line between the two points, as the radius of the unique circle through the three points, or as the radius of the unique sphere through the four points, respectively.We finally define the alpha complex α r (X) of X with radius r to consist of empty simplices of arbitrary dimension ≤ 3 with radius ≤ r.The construction works analogously in an arbitrary number of dimensions.Alpha complexes are simplicial complexes.
In Figure 2(b) we give two-dimensional examples of alpha complexes of increasing radii (indicated by gray disks around the points) for a cloud consisting of points sampled randomly from a circle with Gaussian noise added to their positions.We note that simplices which enter at larger radii have visually larger area.For a finite point cloud there are finitely many different alpha complexes α ri (X) with r i < r j for i < j, i, j = 1, . . ., K for an integer K.They form a filtration: The final, 'full' alpha complex Del(X) is also known as the Delaunay complex or Delaunay triangulation.That alpha complexes of increasing radii are nested can be seen in Figure 2(b).
As we sweep through the filtration of alpha complexes α r (X), topology changes appear.We observe in Figure 2(b) that for the lowest radius the alpha complex consists of many distinct connected components.Upon increasing the radius the connected components successively merge with each other, such that the large loop gets born in the alpha complex.The large loop persists for a comparably large range of radii, until it dies when filled entirely by triangles.Its topology does not 'see' its thickening upon increasing the radius.We describe the topology of α r (X) for a specific radius r by homology, which can be explicitly computed.Homology can describe holes of different dimensions in complexes.A dimension-0 hole is a connected component, a dimension-1 hole is an unfilled, planar-like loop, and a dimension-2 hole is an empty void in the complex.
Persistent homology describes the changing topology of the entire filtration of alpha complexes at once.It associates to a dimension-k hole a birth radius r b and a death radius r d > r b , such that the hole is present in the alpha complexes for all radii r ∈ [r b , r d ), with r b minimal and r d maximal.We can describe the persistent homology of the filtration of alpha complexes {α r (X)} r by the collection of all birth-death pairs (r b , r d ), which we call the persistence diagram, denoted Dgm k (X).In Figure 2(c) we display the dimension-1 persistence diagram of the noisy, approximately circular point cloud.We note that many features appear near the diagonal r b = r d , i.e., have persistence ratio r d /r b near 1, represented by the tiny, short-lived unfilled loops which appear at smaller radii in Figure 2 Persistent homology has a number of useful theoretical properties.It is stable, i.e., perturbations of the input point cloud result in only slight changes of corresponding persistence diagrams with respect to suitable metrics [62,63].Furthermore, most observables computed with persistent homology behave well with regard to statistical limits, and also large volume asymptotics exist, including notions of ergodicity [57,58].The efficient computation of persistent homology is facilitated by a variety of suitable libraries for different programming languages, see [14] for an overview.We employ the versatile GUDHI library [64], and compute periodic alpha complexes to take the spatially periodic boundary conditions of our lattice into account.

B. Self-similarity in persistent homology
How can dynamical self-similarity manifest itself in persistent homology?Let X(t) ⊂ R d be a family of time-dependent point clouds.In [43] the dimension-k persistence pair distribution has been introduced as In general, its expectation value ⟨P k ⟩(t, r b , r d ) exists and is not a sum of Dirac δ-functions anymore [58,65].We note that ⟨P k ⟩(t, r b , r d ) has support only for r d > r b .We say that the expected persistence pair distribution ⟨P k ⟩(t, r b , r d ) scales self-similarly in time if (26) for suitable scaling exponents η 1 , η 2 and arbitrary reference time t ′ in the self-similar regime.This scaling ansatz is similar to the position space scaling ansatz for the T 00 -correlator given in Equation (21).
From ⟨P k ⟩ we can compute many persistent homology descriptors.For instance, the expected number of dimension-k holes at radius r, known as the Betti number, can be computed as With Equation ( 26) Betti numbers scale as ⟨P k ⟩ also yields the total number of homology classes that appear in the filtration at arbitrary radius as which scales as The exponent η 2 thus encodes the power-law decline of the number of persistent homology classes for given birth and death radii r b , r d .Similarly, we can compute from ⟨P k ⟩ the average death radius of dimension-k holes as The exponent η 1 describes the dynamical power-law blow-up of length scales associated to persistent homology.For later use we note that the normalized distribution of persistence ratios π = r d /r b can be computed from ⟨P k ⟩ as For self-similar time evolutions all Π k (t, π) constitute invariants of motion.It has been conjectured in [66] that the distributions Π k (t, π) are in fact universal across different processes to generate i.i.d.point clouds, thus invariant beyond self-similar time evolutions.We will confirm this in Section III D.
To study the self-similar behavior of integrated quantities such as the total number of homology classes n k (t) in lattice simulations, persistent homology classes of sizes near the lattice spacing a s and near the lattice size a s N s need to be explicitly excluded from the integrations.The employed procedure for persistence ratio distributions is described in Section III D.

C. Dynamics of homology classes in Betti numbers
Alpha complexes and their persistent homology require point clouds as input.Motivated by the findings of Section II B, we examine sublevel sets of energy and topological densities, constructed for T 00 i (t, x) and q i (t, x), and computed from an individual classical-statistical sample i as For each such point cloud at time t we compute the persistent homology of its filtration of alpha complexes.
In the classical-statistical approximation we approximate expectation values of persistent homology descriptors as ensemble averages of observables computed for individual samples [43].
In Figure 3 we show results for Betti number distributions of energy and topological density sublevel sets.We set ν T 00 = −0.081/Q 4 and ν q = −0.13/Q 4 , such that the point clouds comprise 20, 000 − 50, 000 points which reflect the dynamics of local minima.For comparison, T 00 and q have been given at the end of Section II B. We notice that dimension-0 Betti numbers monotonously decrease with increasing radii, which originates from connected components merging and thus decreasing in number.Dimension-1 and -2 Betti numbers show distinct peaks.We observe peaks appearing at larger radii the larger the dimension of the holes is.Dimension-0 holes are present at lowest radii, required to merge with each other to form dimension-1 holes, which in turn need to die (getting filled with triangles) to form dimension-2 holes, and thus giving rise to the hierarchy.Focusing on the inset figures for T 00 i and q i sublevel sets, we notice a continuous motion of holes to lower radii with time, across all dimensions.The numbers of holes increase, while the shapes of Betti number distributions qualitatively remain invariant.This suggests that topological structures the alpha complexes decrease in size in the course of time.
The main figures are rescaled according to Equation (28).We propose the exponents η 1 = −1/7, η 2 = −5/7.For these we observe that the rescaled distributions coincide to good accuracy, which is indicative of the validity of the self-similar scaling ansatz for persistent homology in Equation (26).The relation η 2 = 5η 1 is an example of the general packing relation first observed in [43] and rigorously proven in [58], and originates from the geometric bound that the constant volume puts on the number of holes.If holes decrease in size, more of them fit on the lattice.
We explicitly checked for the independence of the scaling behavior and qualitative shapes of Betti number distributions from the choice of ν T 00 and ν q in nearby regimes, see Appendix D. We verified an approximate insensitivity from infrared and ultraviolet lattice cutoffs, such that the non-local features that persistent homology captures are well-resolved and sufficiently dense.
We conclude that the persistent homology of sublevel sets of energy and topological densities reveals selfsimilar scaling related to the energy cascade with the same spatial scaling exponent η 1 = β as for correlation functions.Yet, while energy densities are bounded by positivity, T 00 i ≥ 0, this is not the case for topological densities.This generally results in different functional shapes through the involved point cloud construction.This is to be contrasted with correlation functions.Still, with regard to scaling exponents, we find that the universality of the self-similar dynamics related to the energy cascade encompasses persistent homology observables.

D. Persistence ratio distributions
In Figure 4 we display distributions of persistence ratios, Π k (t, π), for the energy and topological density sublevel sets X T 00 ,ν T 00 ,i (t) and X q,νq,i (t).The distributions Π k (t, π) are particularly sensitive to lattice artefacts.We explicitly removed the two most prominent lattice artefacts from the computed distributions by excluding in dimension 1 holes with π = √ 2 and π = 4/3, and in dimension 2 holes with π = 3/2 and π = 9/8.Square roots of quotients of small integers are unique to specific point configurations on the lattice originating from vertices of one up to few elementary lattice voxels.We checked that further lattice artifacts in this sense do not contribute significantly to the distributions shown in Figure 4.
The distributions Π k (t, π) decrease fast and have support roughly up to π ≈ 2 in dimension 1 and up to π ≈ 1.6 in dimension 2. Crucially, the distributions remain constant in time up to statistical fluctuations.Although showing more fluctuations at smaller Π k (t, π)values, this also applies to the initial time (indicated in gray).Thus, the time translation invariance of Π k (t, π) is not due to the specifics of a self-similar time evolution.The persistence ratio distributions Π k (t, π) rather constitute approximate invariants of motion for all π.
As already mentioned in Section III B, it has been conjectured that persistence ratio distributions are universal across i.i.d.point cloud generation processes [66].We regard the constancy in time of Π k (t, π) in Figure 4 as a manifestation of this.An even stronger universality conjecture for a distribution of transformed π-values has been provided in the same work [66].This applies to our data, too, as we detail in Appendix E.

IV. CONCLUSIONS
In this work we have studied energy and topological densities in real-time SU(2) lattice gauge theory simulations in the classical-statistical regime with over-occupied initial conditions.We found that energy and topological density correlation functions reveal self-similarity related to a direct energy cascade.Corresponding scaling expo-nents can be understood from a kinetic quasi-particle picture (β) and the energy-momentum conservation Ward identity (α(β)).
The approximate agreement of energy and topological density correlations hints at independently and identically distributed colour and space directions of E(t, x) and B(t, x).Differences between energy and topological density correlations could have been due to topological excitations such as sphalerons.While topological excitations typically manifest themselves in coarsening dynamics due to their mutual annihilation, in our study both energy and topological density correlations instead showed a dynamical refinement of associated structures.This suggests that topological excitations are barely of relevance for the dynamics in gluonic plasmas at the length scales probed by the energy and topological densities for the examined over-occupied initial conditions.
In order to look for self-similar dynamics beyond local correlation functions, we constructed non-local observables based on persistent homology.Specifically, we examined the filtration of alpha complexes of energy and topological density sublevel sets, and the filtration of cubical complexes of energy and topological densities (in Appendix F).Even for topological observables like Betti number distributions we observed self-similarity.The associated scaling exponents can be understood from the correlation functions (η 1 = β) and the packing relation (η 2 (η 1 )), which follows from the bounded number of topological features for a given non-expanding lattice volume.We emphasize that all observables investigated are gauge invariant, and the universality related to the direct energy cascade encompasses local and non-local observables.In particular, we have demonstrated in Appendix F that for cubical complex filtrations of sub-or superlevel sets where the spatial metric information does not enter, the packing relation applies as soon as there is any bound on the filtration.This goes beyond previous findings for alpha complexes [43], and is of relevance to any persistent homology study of critical phenomena in or out of equilibrium.
More generally, persistent homology can provide a sensitive probe for non-local structures in lattice field data, as our results and persistent homology studies of the partly intricate phase structures of diverse physical models with non-local order parameters have been demonstrating [47][48][49][50][51][52][53].Non-local excitations including topological defects can be important for a wide range of dynamical phenomena ranging from heavy-ion collisions to ultracold quantum gases [33,67].Beyond that, higherorder correlation functions are numerically hard to access but can contain important information in strongly correlated systems.Persistent homology can thus provide valuable observables in statistical contexts sensitive to higher-order and non-local correlations, in line with the results of the present study.
While we concentrated on gluonic plasma dynamics, understanding far-from-equilibrium universality classes is of relevance for a wide range of applications in complex many-body systems such as experiments with ul-tracold quantum gases.Our study opens a pathway for their characterization including topological and non-local quantifiers.as [39] f EE (t, p) with color decomposition E j (t, p) = E a j (t, p)T a of the Fourier-transformed electric field, T a generators of su(N c ), p ≡ |p|, and the transversal polarization projector P ij T (p) defined as in [39].For simplicity we ignore effects from hard thermal loop self-energies, in particular screening masses, since we focus on hard momentum modes while they would alter correlators primarily at soft momenta.
In Figure 5 we display the occupation numbers for both definitions in Equation (B1) across the entire time range investigated.Without dynamical rescaling, we observe an approximate power-law decline ∼ p −4/3 in momentum space at lower momenta, which for f Ė Ė appears overlaid by an approximate ∼ p −2 behavior of the correlator.Above Q, a steep decline towards zero appears for all times.We note that for both definitions occupations shift towards larger momenta in the course of time.For dynamical rescaling we suggestively set exponents β = −1/7 and α = −4/7 = 4β, consistent with Equation (16).In particular, the same β describes occupation numbers and all dynamical observables considered in the main text.
Let C be a simplicial complex such as the alpha complexes α r (X) constructed in the main text.The same constructions apply for cubical complexes.We consider chain complexes and homology groups with coefficients in Z 2 , such that the k-th chain complex We can study the topology of C using the homology groups H k (C), which capture similar topological information compared to homotopy groups of C, but are in general not the same.Technically, elements of H k (C), called homology classes, are equivalence classes of kcycles modulo higher-dimensional boundary contributions.We may intuitively think of homology classes as independent holes.Their number is the Z 2 -dimension of H k (C), called k-th Betti number, (C2)

Persistent homology groups
Let {C r } r≥0 be a filtration of complexes, such as the alpha complex filtration considered in the main text.We compute all their individual homology groups {H k (C r )} r .In addition, the filtration contains inclusion maps C r → C s for all r ≤ s.These induce maps on homology groups, The ι r,s k map a homology class in H k (C r ) either to one in H k (C s ), if it is still present for C s , or to zero, if corresponding (potentially deformed) cycles appear as boundaries in H k (C s ).Further, non-trivial cokernels can appear for ι r,s k : new homology classes may appear in C s , which are not present in C r .Then, s can be chosen such that for sufficiently small ϵ > 0, We call the collection {(H k (C r ), ι r,s k )} r≤s a persistence module.It is tame, if Equation (C4) holds only for finitely many distinct s-values.
By the structure theorem of persistent homology [71,72], any tame persistence module is isomorphic to its persistence diagram, i.e., the collection of all its birthdeath pairs (r b , r d ), r b < r d ∈ R ∪ {∞}.The same birthdeath pair may appear multiple times.q i (t, x) before computing the persistent homology of their superlevel set cubical complexes.Thus, higher lattice momentum scales than π/(8a s ) do not enter the persistent homology analysis.Concerning the infrared, we have verified an approximate extensive scaling of Betti numbers and other persistent homology quantifiers with respect to the lattice volume.

Cubical complexes
A cubical complex is a collection of cubes of different dimensions, which is closed under taking boundaries, similarly to simplicial complexes.We explicitly construct the filtration of cubical complexes for sub-and superlevel sets of lattice functions at a given time t.Cubical complexes are well-suited to describe sub-and superlevel sets of lattice functions, or pixelized data more generally.Let C be the full cubical complex of the lattice, consisting of one 3-cube x + [−1/2, 1/2] 3 for each spatial lattice point x ∈ Λ s .Such a 3-cube comes with all its faces, edges and vertices, as required to have C closed under taking boundaries.On our lattice, a 3-cube is a cube of side length a s , a 2-cube is a square of side length a s , a 1-cube is an edge of length a s , and a 0-cube is a point.In the literature, such cubes are called elementary [73].Analogously, any 1-cube is contained in the boundaries of four 2-cubes, and any 0-cube is contained in the boundaries of six 1-cubes.We inductively apply Equation (F1) to construct T t,i for lower-dimensional cubes from higherdimensional ones, until T t,i is defined on all C.This construction is called the lower star filtration.
We define cubical complexes corresponding to lattice sublevel sets of T 00 i (t, x) at time t as These are closed under taking boundaries, thus indeed cubical complexes.We define superlevel set cubical complexes as reflecting the structure of lattice superlevel sets.For topological densities the construction is the same; simply replace T 00 i (t, x) by q i (t, x).We are interested in the persistent homology of the filtrations of complexes {C T 00 i (t, ν)} ν and {D T 00 i (t, ν)} ν , as ν is swept through.Indeed, we notice that they define filtrations, e.g.C T 00 i (t, ν) ⊆ C T 00 i (t, µ) whenever ν ≤ µ.The previous constructions of persistent homology apply also in this setting.The persistent homology of cubical complex filtrations can be efficiently calculated [73].We again use GUDHI [64] to compute persistent homology of the cubical complex filtration, and use periodic cubical complexes to take spatially periodic boundary conditions into account.
In Figure 9(a) we illustrate the meaning of cubical complex homology classes.In Figure 9(b) we sketch the meaning of birth b and death d of homology classes as the filtration of superlevel sets is swept through.We note that the birth of dimension-0 homology classes happens at function values of local maxima.Their persistences quantify their dominances, dying when merging at saddle points with tails of other local maxima.Dimension-1 homology classes can be thought of as describing potentially deformed volcano-like features in function landscapes, and analogously for higher-dimensional homology classes.

Betti number distributions
In Figure 10 we show results for Betti number distributions of all dimensions for both the energy and the topological density superlevel set filtrations using cubical complexes.Insets show numbers without, and main plots with dynamical rescaling.Note that energy densities are shown with the average subtracted.We notice from the unrescaled figures that topological features of both energy and topological densities give rise to a hierarchy of peaks between dimensions.For superlevel sets of a function whose values scatter approximately symmetrically around zero, connected components appear first at positive filtration parameters.Multiple such connected components need to merge to form a dimension-1 feature, which is born around zero.Typically at negative filtration parameters, multiple dimension-1 features die to give birth to a dimension-2 feature.
In the course of time, throughout dimensions and energy and topological densities topological features shift towards zero filtration parameters suggesting a homogenization of structures.Further, their number increases, indicative of their dynamical refinement.In the main plots we suggestively rescaled the figures with η 1 = −1/7 and η 2 = −3/7 = 3η 1 .We notice approximate matching of the rescaled curves.Notice that aside of the initial condition we here show data from Qt = 300 onwards, while in the main text we show data starting with Qt = 100.Betti number distributions do not allow sufficiently consistent rescaling before Qt = 300 (not displayed).
The exponent η 1 = −1/7 is consistent with the earlier findings of η 1 for alpha complexes and β for correlation functions.However, here it describes the dynamical refinement of structures in the space of function values, not with respect to spatial sizes of structures.How can the exponents still agree?We recall from Section II B that the correlators ⟨(T 00 (t, p x )) 2 ⟩ c and ⟨(q(t, p x )) 2 ⟩ c scale overall as ∼ t β for spatial scaling exponent β = −1/7.If such correlators contribute predominantly to local fluctuations in T 00 i (t, x) and q i (t, x), we can heuristically find that topological features scale in the space of function values as ∼ t β , too.We take our observation of η 1 = −1/7 for the superlevel set filtration  as an indication for this.
The number of topological features which appear in T 00 i superlevel sets is bounded by total energy conservation in a similar way that the number of topological features for alpha complexes is bounded by the constant system volume.This effectively provides a one-dimensional constraint for the persistent homology of the superlevel set filtration.Then, for cubical complexes of T 00 i superlevel sets the packing relation yields η 2 = 3η 1 [43,58], consistent with the data.
The apparent similarity of energy and topological density cubical complexes in persistent homology can again be understood from independent spatial and colour directions of E(t, x) and B(t, x) as proposed in Section II B. Yet, differences are visible in particular for the left-sided tail of dimension-2 Betti numbers, where the Betti number distribution has larger support for topological rather than energy densities.

Persistence distributions
In Figure 11 we display the distribution of absolute persistences, d − b, for energy and topological density superlevel sets.Again, insets show figures without rescaling, which we discuss first.We notice that dimension-0 distributions decline approximately exponentially with persistence values.Dimension-1 distributions have more restricted support only up to T 00 − T 00 , q ≃ 0.05 Q 4 .Persistence distributions of dimension-2 features of en-ergy densities have more restricted support compared to topological densities, whose dimension-2 persistence distributions moreover show approximate exponential behavior.This indicates that energy density persistences have an upper bound due to the positivity of T 00 i , which is not the case for q i .There, the approximate symmetry among dimension-0 and dimension-2 features hints at a symmetric distribution of topological density values around zero, which in particular features local lumps, similarly to the deductions in [53].
Let us now consider absolute persistences defined as integrals of the averaged dimension-k persistence pair distributions over the death parameter, (F4) The scaling ansatz Equation ( 26) then leads to the selfsimilar scaling P k (t, p) = (t/t ′ ) η1−η2 P k (t ′ , (t/t ′ ) −η1 p) .(F5) Accordingly, we have rescaled the curves in the main plots in Figure 11 using the same exponents η 1 = −1/7 and η 2 = −3/7 as for the Betti numbers shown in Figure 10.We observe approximate agreement among the rescaled curves mostly up to statistical fluctuations, which is indicative of dynamical self-similarity.We suggest that the deviations for dimension-2 persistences of energy density sublevel sets are due to the positivity bound.

Figure 1 .
Figure 1.Connected two-point correlation functions of energy density T 00 (t, px) (top) and topological density q(t, px) (bottom).Variants without (left) and with dynamical rescaling (right) are shown.Gray: Initial condition.

Figure 2 .
Figure 2. (a): Simplices of different dimensions; the 3-simplex is to be regarded as filled.(b): Alpha complexes of increasing radii.The point cloud in (b) consists of points sampled randomly from a circle with Gaussian noise added on their positions.(c): Persistence diagram of dimension-1 holes of the point cloud.Dashed line: the diagonal r b = r d .
(b).One prominent feature appears at large r d , highlighted in Figure2(c).Having a high persistence ratio r d /r b ≃ 5, this feature is represented by the large circular structure present in the point cloud.In summary, persistent homology provides a quantitative means to extract topological structures from point cloud data, and incorporates persistence measures of their dominance.

Figure 4 .
Figure 4. Normalized persistence ratio distributions of the alpha complex filtration for T 00 (t, x) (top) and q(t, x) sublevel sets (bottom), with ν T 00 = −0.081/Q 4 and νq = −0.13/Q 4 .Distributions at the initial time are displayed in gray.
C k (C) of C consists of formal sums of chains of ksimplices with Z 2 -coefficients.The boundary operator ∂ k : C k (C) → C k−1 (C) is defined to map a chain of k-simplices to its boundary, which is a (k − 1)-chain.Since boundaries of chain boundaries are empty, on has ∂ k−1 • ∂ k = 0. We define the cycle group Z k (C) := ker(∂ k ) consisting of all closed k-chains, i.e., k-chains without boundary.We further define the boundary group B k (C) := im(∂ k+1 ) consisting of all k-chains which are boundaries of (k + 1)-chains.We find B k (C) ⊆ Z k (C) as subgroups, such that we can define their quotient groups, H k (C) := Z k (C)/B k (C) , (C1) called homology groups.

Figure 9 .
Figure 9. (a): Homology classes of different dimensions in cubical complexes, the enclosed volume is indicated in red.(b): Persistent homology classes of dimensions 0 and 1 arising from superlevel sets of a function with 2-dimensional domain.Figure adapted from [53].

Figure 10 .
Figure 10.Betti number distributions of the T 00 (t, x) (top) and the q(t, x) superlevel set filtrations (bottom) using cubical complexes.The initial time is displayed in gray.The scaling exponents are set to η1 = −1/7 and η2 = 3η1 = −3/7.The insets show the original distributions without rescaling; the inset axes have the same ranges as the axes in the main plots.

Figure 11 .
Figure 11.Persistence distributions of the T 00 (t, x) (top) and the q(t, x) superlevel set filtrations (bottom) using cubical complexes.The initial time is displayed in gray.The scaling exponents are set to η1 = −1/7 and η2 = 3η1 = −3/7.The insets show the distributions without rescaling; the inset axes have the same ranges as the main plots.