Displacement autocorrelation functions for strong anomalous diffusion: A scaling form, universal behavior, and corrections to scaling

Displacement autocorrelation functions for strong anomalous diffusion: A scaling form, universal behavior, and corrections to scaling / Vollmer, JÃ1⁄4rgen; Rondoni, Lamberto; Tayyab, Muhammad; Giberti, Claudio; Mej('(i))aMonasterio, Carlos. In: PHYSICAL REVIEW RESEARCH. ISSN 2643-1564. STAMPA. 3:013067(2021), pp. 1-14. Original Displacement autocorrelation functions for strong anomalous diffusion: A scaling form, universal behavior, and corrections to scaling

The statistical properties of transport are given by the probability density of the displacement P( x, t ). However, in most situations this density is unknown and transport is studied in terms of the asymptotic behavior of its moments in time where · denotes the average over an ensemble of trajectories, p ∈ R is the order of the moment, and the function γ (p) is called the spectrum of the moments of the displacement. The exponent η = γ (2) characterizes the mean-square displacement. Transport is called subdiffusive for 0 < η < 1, diffusive for η = 1, superdiffusive for 1 < η < 2, and ballistic when η = 2.
Transport is termed scale invariant when the probability density of displacements x admits a scaling P( x, t ) = t −ν F ( x/t ν ) with ν constant, i.e., when all moments of the displacement are characterized by a single scale t ν . The spectrum of the moments of the displacement will then be a linear function of p: γ (p) = ν p.
When the spectrum γ (p) is nonlinear, transport is called strong anomalous diffusion [18]. This phenomenon has been observed in a variety of simple stochastic systems [19], Hamiltonian systems with mixed phase space [18], polygonal billiard channels [20,21], billiards with infinite horizon [22][23][24], one-dimensional maps [25,26] and intermittent maps [18,27], running sand piles [28], stochastic models of inhomogeneous media [29,30], the diffusion in laser-cooled atoms [31], in experiments on the mobility of particles inside living cancer cells [32], particles passively advected by dynamical membranes [33], and in the bulk-mediated diffusion on lipid bilayers [34]. A widely investigated case is the one in which different scalings of the bulk and the tail of the probability distribution give rise to a piecewise linear form of the scaling exponent γ (p), Strong anomalous diffusion is believed to be generic [35][36][37][38] for dynamics with fat-tailed waiting-time distributions. In recent years its dynamical basis has been analyzed through generalizations of the central limit theorem and non-normalizable densities [17,35,38,39]. At the same time, numerous studies have addressed the relation between the properties of deterministic dynamics and transport. Roughly speaking, chaos is commonly associated with fast decay of correlations and normal diffusion, while nonchaotic dynamics often leads to anomalous transport and slow decay of correlations [1,3,12,21,26]. Stochastic processes, on the other hand, are often considered to resemble chaotic dynamics, but they can give rise to normal as well as anomalous diffusion, depending on their correlations decay rate [40,41]. Thus, numerous questions remain open [3,12,15,29,[42][43][44][45]. In particular, the asymptotic behavior of correlation functions is not understood in general, although it is relevant, e.g., to distinguish transport processes that arise from different microscopic mechanisms but have the same moments [15,36,[46][47][48].
In Ref. [26] a deterministic map named slicer map (SM) was introduced to shed light on these issues. It was shown that for an appropriate matching of parameters its moments of the displacement, of order 2 and higher, scale in time like those of the Lévy-Lorentz gas (LLg) [29], a transport model [26,49] where strong anomalous diffusion emerges as a consequence of scattering in an environment with quenched disorder. In Ref. [49] it was found that the matching of exponents also entails the equality of the large-time scaling of the displacement autocorrelation function, making two entirely different systems hardly distinguishable on the level of the statistics of displacements.
Here, we introduce the fly-and-die (FnD) model as a minimal exactly solvable deterministic model where strong anomalous diffusion emerges due to a light front, i.e., ballistic trajectories that did not undergo transitions up to any finite time t. Unlike the SM the FnD is continuous in time and space. Otherwise, it closely mimics the asymptotic transport properties of the SM (and hence of the LLg). We present analytic expressions for its spectrum of the moments γ (p) and its displacement autocorrelation function φ(t, t + h). Our main result is a scaling form expressing the autocorrelation function as a function of the reduced time h/t, This result is interesting for two reasons: (i) There are only very few analytic results about position-position autocorrelation function for strong anomalous diffusion. To our knowledge they have only been calculated for Lévy walks [50,51] and the slicer map [49], and those studies did not introduce the scaling form, Eq. (4). (ii) We propose that Eq. (4) is universally valid. It applies for every system showing strong anomalous diffusion where the mean-square displacement is governed by the light front. This claim is underpinned in the second part of the paper by showing that the correspondence between the transport properties of the FnD, the SM, and the LLg extends to a much broader class of systems. To this end we compare the analytical dependencies of the FnD with the scaling function for Lévy walks that is derived here based on analytical results of Ref. [51], and with numerical results for systems where analytical treatments are beyond reach. For all dynamics the scaling form, Eq. (4), provides a very good data collapse. However, careful inspection reveals that there are small systematic corrections to scaling with system-specific features, or related to nonasymptotic effects. This paper is organized as follows: In Sec. II, we introduce the FnD dynamics and obtain analytical expressions for the spectrum of moments and the autocorrelation function. In Sec. III, we revisit the slicer map (Sec. III A), Lévy walks (Sec. III B), the Lévy Lorentz gas (Sec. III C), the Lorentz gas (Sec. III D), and a family of polygonal billiard channels (Sec. III E), comparing their spectrum of moments, γ (p), with that of the FnD dynamics. The scope of this discussion is strictly confined to material required to address the correlation function in a self-contained manner. We provide links to the literature for relevant further work. In Sec. IV, the comparison concerns the scaling of the autocorrelation function φ(t 1 , t 2 ), which is not treated elsewhere. Section V highlights predictions of the FnD model that apply universally for all models from the perspective of scaling theory and the Buckingham-Pi theorem, and it addresses model-specific features. Section VI summarizes our findings.

II. THE FLY-AND-DIE MODEL
In this section we introduce the fly-and-die (FnD) model.

A. Dynamics
The FnD dynamics addresses the motion of a particle on the semi-infinite line [0, ∞). Different trajectories are labeled by initial conditions x 0 ∈ [0, 1]. Starting at x 0 the particle initially moves along the positive x axis with unit velocity (it flies). At a time t c (x 0 ) it stops and never resumes motion (it dies). Hence, at time t a particle is found in position Anomalous transport emerges when the distribution of the flight times t c (x 0 ) has a power-law tail. Here, we consider the case with ξ and b positive constants, and initial conditions x 0 uniformly distributed in the interval [0, 1]. The minimum flight time is t M = b 1/ξ . In simple words, the FnD dynamics evolves an ensemble of initial conditions in the unit interval. Each point is associated to a trajectory that flies independently and without turning back for a time that is determined by the initial condition x 0 via the function t c (x 0 ). The distribution of the flight times t c determines the transport properties of the FnD dynamics.

B. Moments of displacement
For the FnD dynamics the spectrum of the moments of the displacement is obtained straightforwardly: Initial conditions, x 0 , that perform a flight longer than t, lie in the interval [0, (t M /t ) ξ ). Therefore, adopting the uniform distribution of initial conditions in the unit interval, the probability to fly for a time larger than t is given by In the following we only consider times t t M . For p = ξ the pth moment of the displacement x(t ) takes the form For p = ξ an analogous calculation yields This expression can also be found as the p → ξ limit of Eq. (7a), by observing that At large times, t b 1/ξ , the spectrum of the moments takes the form: The FnD exhibits piecewise-linear scaling in the statistics of the displacement with an exponent γ (p) that is zero for p < ξ and linear, γ (p) = p − ξ , for p > ξ. The mean-square displacement, | x| 2 , approaches a constant value (localization) for ξ > 2, it grows logarithmically for ξ = 2, and it grows as a power law with a mean-square displacement scaling exponent η = 2 − ξ for 0 < ξ < 2. Then, ξ ∈ (0, 1) corresponds to superdiffusive behavior, ξ = 1 to normal diffusion, and ξ ∈ (1, 2) to subdiffusion: The FnD gives rise to all anomalous transport regimes.

C. Displacement autocorrelation function
We now turn our attention to the autocorrelation function Without loss of generality, we set t 2 > t 1 . Accordingly, we split the integration range in three intervals 0 < x 0 < P(> t 2 ): The trajectory is still flying at time t 2 , hence x(t 1 ) = t 1 and x(t 2 ) = t 2 .
Performing calculations analogous to the derivation of Eq. (7a), one finds for η > 0 that (8) It is convenient to normalize this expression by (9) Subtracting one and denoting the time difference by h = t 2 − t 1 , one obtains In the large-time limit t M /t 1 → 0 for a fixed value of h/t 1 , i.e., for a finite ratio t 2 /t 1 , this entails a scaling form Here, again the η = 1 case is found as the η → 1 limit of the η = 1 case. Equation (11a) connects the large-time behavior of the correlation function to the dependence of the meansquare displacement. It has the following asymptotic scaling for small and large values of h/t 1 , As anticipated in Eq. (4) the scaling form, Eq. (11), represents the time dependence of the displacement autocorrelation function φ(t 1 , t 2 ) in terms of a single dimensionless time ratio 013067-3

D. A remark on the equivalence of statistics
Equation (7c) entails that the moments of the FnD dynamics show a piecewise-linear behavior, as in the case of strong anomalous diffusion.
For p > ξ the moments of the FnD dynamics grow in time with an exponent γ (p) = p − ξ . A linear behavior with slope one for γ (p) is typical when the contributions to the higher cumulants are asymptotically dominated by ballistic flights or light fronts [18,37,38]. For Lévy processes this was explicitly discussed in Ref. [36]. Recently, it has been suggested that the scaling holds in general due to a "big jump principle" [37,38]. The mean-square displacement of the FnD dynamics follows this behavior when ξ < 2. In the next section we will explicitly verify that the matching of the value of the mean-square displacement will then also fix all other moments in the largep branch of the piecewise-linear behavior of γ (p).
For the FnD dynamics, the moments p < ξ are constant in time. This is a peculiarity of the dynamics that is nongeneric: The FnD dynamics does not mimic the lower moments of displacement for other dynamics. This is a strong point of the FnD model. It allows us to clearly pin down features of dynamics that derive from rare ballistic flights. The FnD prediction, Eq. (11), for the displacement correlation function φ(t 1 , t 2 ) is of particular interest. It suggests a scaling form, Eq. (4), that turns out to be astonishingly robust in vastly different settings.

III. MODELS AND MOMENTS
In this section, we investigate how the spectrum of the moments of the displacement in dynamics with light fronts can be mapped to the spectrum of the FnD. First, we address the relation of the FnD to the slicer map (SM). Then, we consider a Lévy-walk model that is the most widely studied model featuring strong anomalous diffusion. Subsequently, we discuss a continuous-time stochastic process, the Lévy-Lorentz gas (LLg), a chaotic billiard with infinite horizon, the periodic Lorentz gas (LG), and some polygonal billiard channels (PBC). For each system we also briefly revisit key features of its dynamics.

A. The slicer map (SM)
The SM is a one-parameter deterministic exactly solvable model S α : [0, 1] × Z → [0, 1] × Z defined by [26,49]: For all integers m the so-called "slicer" m (α) is defined by For 1/2 < x < 1 each iteration of the map increases the values of m by one, until x > m (α). Subsequently, the trajectory enters a stable period-two cycle, oscillating back and forth between the two neighboring sites m and m − 1. Similarly, for 0 < x < 1/2 each iteration of the map decreases the values of m by one, until x < − m (α), and then the trajectory enters a stable period-two cycle. The SM was inspired by the dynamics of polygonal billiards, which have vanishing Lyapunov exponents, because their trajectories separate substantially only at a countable set of points [11,20,21]. Analogously, the distance between two points x 1 and x 2 of the SM jumps discontinuously when they reach a cell m where m (α) ∈ [x 1 , x 2 ]. The corners of polygons act as slicers of the bundle of initial conditions (cf. the discussion of polygonal billiards in Sec. III E).
The moments of the displacement of an ensemble where particles are initially distributed uniformly in the interval [0,1] were obtained in Ref. [26], . In the large-time limit, n 2 1/α , this expression has the asymptotic behavior The SM exhibits strong anomalous diffusion with a piecewiselinear spectrum γ (p) and threshold order p c = α. The SM exhibits all transport regimes upon varying the parameter α. It features subdiffusion for 2 α > 1, diffusion for α = 1, and superdiffusion for 1 > α > 0. Comparing Eq. (14) with Eq. (7c), we see that the spectrum of the moments of the displacement of the FnD dynamics and the SM coincide when taking ξ = α and b = t ξ M = 2. These expressions for ξ and b can also be found by matching the mean-square displacement, i.e., the expressions for p = 2, and observing that the two branches of γ (p) have slope zero and one, respectively.

B. Lévy walks (LW)
Lévy walks refer to dynamics where an agent walks with constant speed v along straight line segments of length r that are sampled from a probability distribution with a powerlaw tail [16,19,36]. In the simplest case the dynamics is one dimensional. At the end of each segment there is a random decision taken to go left or right, and the length r is sampled from the probability distribution where the exponent β LW > 0 characterizes the power-law decay of the distribution, and r 0 > 0 is the minimum distance between scatterers.
The moments of displacement of an ensemble of particles that start a LW at x = 0 were obtained in Refs. [19,39], For β LW < 1 the LW exhibits anomalous diffusion with ballistic scaling, η = 2. For β LW > 1 it features strong anomalous diffusion with a piecewise-linear spectrum γ (p) and threshold order p c = β LW . For β LW > 2 there is subdiffusion with exponent η = 2/β LW , and for 1 < β LW < 2 superdiffusion with exponent 3 − β LW . Comparing Eq. (16) with Eq. (7c), we see that the spectrum of the moments of the displacement of the FnD dynamics and the LW coincide for large moments p, when one matches parameters as ξ = −1 + β LW . In the superdiffusive regime this follows from matching the mean-square displacement and acknowledging that superdiffusion emerges in the models due to rare, very long ballistic trajectories, i.e., that γ (p) has slope one.

C. The Lévy-Lorentz gas (LLg)
The Lévy-Lorentz gas (LLg) was introduced in Ref. [44] (and studied further in, e.g., Refs. [29,45]) as a onedimensional model of anomalous transport in semiconductor devices. In the LLg model a particle is randomly scattered with probability 1/2 at randomly fixed positions on the line. Between two consecutive collisions the particle moves at constant velocity ±v. The distances r between neighboring scatterers are sampled from a Lévy distribution with probability density i.e., the same probability distribution also adopted for the LW, cf. Eq. (15). The spectrum of the moments of the displacement was derived under the assumption that higher moments are dominated by the light front, an assumption also denoted as single long jump principle. For nonequilibrium initial conditions, i.e., an ensemble of trajectories starting all at the same scatterer position, it entails [29] (18) The LLg shows strong anomalous diffusion with mean-square displacement r 2 (t ) ∼ t η growing as a power law in time with exponent It exhibits superdiffusive transport for 0 < β LLg < 3/2 and diffusive transport for β LLg 3/2. Unlike the FnD dynamics and the SM, it has no subdiffusive transport regime. On the other hand, the LLg dynamics is much more complex as compared to the former models. Only its moments of displacement, Eq. (18), have been analytically expressed. In Ref. [26] it was established that the exponent η characterizing the mean-square displacement of the LLg, Eq. (19), agrees with the relation η = 2 − α for the SM dynamics if For β LLg < 3/2 the mean-square displacement, p = 2, lies in one of the large-p branches of γ (p), such that the matching of the mean-square displacement entails an agreement of the spectrum of the moments of the displacement for p > β LLg when 0 < β LLg < 1 and for p > 2β LLg − 1 when 1 < β LLg < 3/2. On the other hand, for β LLg > 3/2 the matching only works for the mean-square displacement because We conclude that the matching of the mean-square displacement extends to the large-p branch of the spectrum γ (p) of the LLg by taking ξ = α. The matching of the parameter α of the SM with β LLg of the LLg involves three branches, Eq. (20), while there is a single linear relation required for the matching of the SM and the FnD dynamics. The authors of [26] took this nontrivial matching as evidence that the equivalence is not complete.

D. The Lorentz gas (LG)
Billiards with infinite horizon are a class of deterministic dynamics where the length of free flights is not bounded. The resulting light fronts induce rare ballistic excursions of the displacement and hence strong anomalous diffusion [22]. A distinguished example is the Lorentz gas (LG), which is a paradigmatic model of deterministic diffusion [52,53].
The periodic LG consists of an array of circular scatterers of radius R periodically placed on a triangular lattice of the plane. The separation between nearest-neighbor scatterers is set to = 8 3 cos π 6 so that the horizon is infinite for R < 1 (see inset of Fig. 1). The dynamics consists of free flights of point particles of unit velocity that are started in a region close to the origin and undergo specular collisions with the scatterers.
With infinite horizon the variance of trajectory lengths diverges because the probability of a trajectory segment with no collision and length between l and l + dl decays as l −3 . In a seminal work Bleher showed that for every periodic configuration of scatterers with infinite horizon the displacement scaled by √ t ln(t ) converges in distribution to Gaussian statistics [54]. This weak superdiffusion is hard to observe numerically as the time scales at which it sets in are physically unobservable. Since Bleher's result the existence of this weak superdiffusion in infinite horizon billiard tables has generated a great deal of research [22,[55][56][57][58]. Recently in Ref. [59] a refined central limit theorem in which the √ t ln(t ) scaling to the displacement is substituted by a rescaled Lambert function has shown to appropriately describe the mean-square displacement, at least for some geometries of the infinite horizon periodic LG (see also Ref. [60]).
At infinite horizon (R < 1), it was observed that the LG exhibits strong anomalous diffusion [22], with spectrum When the horizon is finite (R 1), transport is diffusive, and the moments of displacement scale in time with exponent γ (p) = p/2 for any order p. For R < 1, the mean-square displacement is at the crossover of the two branches of the spectrum of the moments of the displacement: For moments p < 2 the spectrum shows normal diffusive transport, while it is dominated by the ballistic excursion for p > 2. The latter branch agrees with the FnD dynamics with ξ = 1. This is shown in Fig. 1, where results are reported for numerical simulations of 10 6 trajectories lasting a total time of 10 6 , with randomly chosen initial positions r(0) and velocities v(0). Solid squares show the spectrum γ (p) for the infinite-horizon case with R = 0.1 and open squares for the finite-horizon case with R = 1.1.

E. Polygonal billiard channels (PBC)
The scatterers of the Lorentz gas are dispersive such that its dynamics is chaotic. Drastically different are the polygonal billiard channels (PBC) [20,21,23,61]. Their dynamics is that of point particles undergoing specular reflections with upper and lower straight walls, periodically arranged in one direction. The geometry of the elementary cell is specified by four parameters as shown in the inset of Fig. 2, where we follow the notation that has been used in Ref. [21]. When y b + y t < H the channel has an infinite horizon, meaning that a trajectory can move for arbitrarily large distances without colliding with a wall. We consider in this case an ensemble of particles that start at a random position in a random direction in a unit cell of the periodic channel. Similarly to the action of the SM, a bundle of nearby trajectories only separates when "sliced" by a corner. The resulting subpopulations of the bundle will follow qualitatively different paths after they hit separate line segments of the wall. Analogously to the LG case, the probability density of trajectory segments of length l scales as ∼l −3 . For all PBCs there are families of trajectories which travel arbitrarily large distances in a given direction without reversing their motion [20,21,23], effectively acting as light fronts.
In other aspects the dynamics of the PBC differ substantially from the LG. Firstly, in the PBC the separation of trajectories produced by the channel corners is not exponentially fast. Secondly, the ballistic excursions appear not only as the result of an infinite horizon but also because of the poor mixing due to nondefocusing collisions [20,21].
In Fig. 2 we show the spectrum γ (p) for three different PBC with x = 1, y t = 0.77, y b = 0.45 and different widths H (previously considered in Ref. [20]). The spectrum γ (p) is described by with p c = 2, 3, and 4 for H = 1.27, 1.17, and 1.07, respectively.
For H = 1.27 the PBC has an infinite horizon, and we observe the same spectrum of the moments of the displacement as for the LG. Quite interesting is the behavior of PBCs with finite horizon, H = 1.17 (green open circles) and H = 1.07 (red open squares). All trajectories experience collisions with the walls within a maximum finite time. However, the bundles of trajectories that mimic light fronts have a weaker impact on the spectrum of the moments of the displacement [21]. In particular they do not dominate the mean-square displacement. Moments dominated by the light front arise for p > p c > 2. For H = 1.17 the ballistic scaling sets in at p c = 3 and for H = 1.07 at p c = 4. Therefore, matching the mean-square displacement with that of the FnD dynamics does not provide information on the spectrum of the moments of the displacement for any other value of p. On the other hand, for ξ = p c /2 the FnD dynamics will describe the spectrum for p > p c .

F. Summary
The large-p branch of the spectrum of the moments of the displacement of strong anomalous diffusion can be matched with the FnD dynamics by an appropriate choice of its parameter ξ . When the mean-square displacement falls into this branch, i.e., when 2 p c , the matching is obtained by equating the exponent η of the mean-square displacement of the considered model with 2 − ξ for the FnD dynamics and observing that the large p branch of γ (p) has slope one for systems with light fronts.

IV. DISPLACEMENT AUTOCORRELATIONS
In this section, we discuss the scaling form, Eq. (4), of the displacement autocorrelation functions. We verify that the scaling holds for the analytical expressions known for the SM [49] and the LW [51]. Subsequently, we analyze the autocorrelation functions for the LLg, the LG, and the PBC, where the correlation functions are not analytically known. A data collapse based on Eq. (4) works in all cases. For small h/t the systems even follow the functional form of the function C(h/t ) of the FnD, up to system-specific corrections for the billiard systems. For large h/t the functions C(h/t ) for the different systems may differ. In all cases the comparison to the expression Eq. (11) for the FnD provides valuable insights into the scaling of the autocorrelations.

A. The slicer map (SM)
The autocorrelation function of the SM, φ(m, n) = (x m − x 0 ) (x n − x 0 ) , was obtained in Ref. [49]. In the limit of large times m and n one has For m = n, the expression reduces to the mean-square displacement: The autocorrelation function of the SM, Eq. (23), exhibits the same scaling as Eq. (11) after identifying η = 2 − α, as established for the spectrum of moments in Eq. (14). The constant terms are subdominant for η > 0 and n > m 1.  [51] (details are given below). There is a perfect data collapse for all values of β LW . For h/t 1 10 −2 there even is quantitative agreement between the LW and the FnD prediction, Eq. (11a).
For β LW = 0.1 we followed 2 × 10 6 trajectories for a time 10 9 , and for β LW = 0.3 we followed 5 × 10 5 trajectories to the same time 10 9 . This provides data for N 16 that allows us to explore the correlation function in the range of reduced times, 10 −8 h/t 1 10 5 , i.e., for 13 orders of magnitude.
For all values β LW < 1 the LW shows ballistic transport with exponent η = 2 such that Eq. (11a) reduces to the linear function C(h/t 1 ) = h/t 1 . For h/t 1 1 the data follow the predicted linear scaling. However, for larger h/t 1 it is better described by power laws with exponents 0.9 and 0.7, respectively (not shown). In this range the data characterize the decay of correlations.
Froemberg and Barkai [51] showed that in the ballistic regime the correlation function involves a sum of three incomplete Bessel functions of t 1 /t 2 with parameters that depend on η and prefactors proportional to t 1 t 2 , t 2 2 , and t 2 1 , respectively. This entails the scaling form, Eq. (4), since φ(t 1 , t 1 ) ∼ t 2 1 . Moreover, in their Eq. (15) they provide the asymptotic scaling of φ(t 1 , t 1 + h) for small and large h/t 1 . In combination with the mean-square displacement their result provides for η = 2 1 .

(25a)
For h/t 1 1 this result predicts the observed crossover to power laws with exponent 1 − β LW = 0.9 and 0.7, respectively. This exponent depends on β LW and clearly FnD tells nothing about this crossover. However, for h/t 1 1 the functional dependences agree. This agreement is remarkable because for ballistic transport the LW shows single-scale ballistic scaling. There is not even strong anomalous diffusion.
For β LW = 1.2 we followed 1.2 × 10 5 trajectories to a time 10 8 , providing data covering 11 orders of magnitude, 10 −7 h/t 1 10 4 . Equation (16) provides η(β LW = 1.2) = 1.8 and accordingly Eq. (11b) predicts a crossover from a power law with exponent 1 for h/t 1 1 to an exponent 0.8 for h/t 1 1. The upper left inset of Fig. 3 shows the ratio of the data and the FnD prediction. For small h/t 1 the data approach the FnD asymptotics with deviations of only a few percent that lie within the scatter of our data. For large h/t 1 the data lie below the FnD dependence. In the main panel one sees that for h/t 1 10 the function is constant to within our numerical resolution.
For β LW = 1.5 the numerical simulations are still more challenging: We only followed 9 × 10 4 trajectories until time 10 7 , providing data for 10 −6 h/t 1  According to Ref. [51] the correlation function of a LW with length and velocity scales r 0 = v = 1 takes the following form in the superdiffusive regime After substituting η = 3 − β LW this provides the following scaling form for the LW correlation function For small h/t 1 this is identical to the FnD expression, Eq. (11a), and for large h/t 1 it accounts for the observed saturation of the scaling function. We conclude that an important qualitative prediction of the FnD model holds for all data sets: The displacement autocorrelation function of the LW admits a data collapse in the form of Eq. (4). Moreover, the FnD provides an accurate, parameter-free description for small values of the reduced time h/t 1 . On the other hand, for h/t 1 1 the functional dependence of the reduced correlation functions for FnD and the LW differ. Indeed, for a strongly ergodic system the correlation φ(t 1 , t 2 ) = x(t 1 ) x(t 2 ) should factorize in the limit of t 1 , t 2 , h/t 1 → ∞ and decay to zero for a symmetric dynamics, x(t 1 )x(t 2 ) → x(t 1 ) x(t 2 ) = 0. As a consequence, φ(t 1 , t 2 )/φ(t 1 , t 1 ) − 1 → −1 for large h/t 1 . However, neither the FnD nor the LW enjoy such ergodic features [16].

C. The Lévy-Lorentz gas (LLg)
For the LLg comprehensive numerical simulations for the autocorrelation function, Eq. (3), have been reported in Ref. [49]. For the present study we performed additional simulations where we evaluated the correlation functions for the times t 1 = 10 m/3 and t 2 = t 1 + 10 n/3 with 9 m, n N. Figure 4 shows the scaling representation, Eq. (4), for three different parameter values of the LLg, β LLg = 0.1 (upper curve, N = 24), β LLg = 0.3 (middle curve, N = 24), and β LLg = 0.6 (lower curve, N = 21), respectively. The value of the time t 1 and h are indicated by the same choice of rainbow colors and symbols as in Fig. 3. Again the data for different β LLg are shifted by successive factors of ten.
As suggested by Eq. (4) all correlations collapse in a master curve when φ(t 1 , t 2 )/φ(t 1 , t 1 ) − 1 is plotted as function of h/t 1 . The FnD expression for the asymptotic scaling, Eq. (11a), is provided by solid gray lines.
For β LLg = 0.1 we performed simulations where we followed 10 3  For increasing values of β LLg the trajectories perform more collisions per unit time. Hence, it becomes harder to reach large times in the simulations, but at the same time one needs fewer realizations of disorder to achieve faithful averages. For β LLg = 0.3 we followed 10 3 trajectories in each of 1150 realizations of the scatterers up to a time 10 9 . This allows us to follow the correlation function in the range of reduced times, 10 −8 h/t 1 10 5 . Estimates of the statistical errors are smaller than the points for most data. The mean-square displacement grows with η = 1.93 such that the two branches have slopes 1 and η − 1 = 0.93, respectively. The change of slope is then becoming noticeable for our numerical data. Also in this case the data lie on the prediction over the full parameter range of h/t 1 that is 13 orders of magnitude wide.
For β LLg = 0.6 we followed 10 3 trajectories in each of 260 realizations of the scatterers up to a time of 10 8 . This allows us to explore the correlation function in the range of reduced times, 10 −7 h/t 1 10 4 , i.e., for 11 orders of magnitude. The mean-square displacement grows with η = 1.775 such that the two branches have clearly different slopes, 1 and 0.775, respectively. For h/t 1 1 the data follow exactly the predicted linear scaling. However, for larger h/t 1 they tend to lie slightly below the prediction.
In the inset to the upper left we show the ratio of the data for β LLg = 0.1 and the FnD prediction. There is a perfect match for h/t 1 10 −2 , and the data deviate systematically towards smaller values for larger h/t 1 . However, the deviation is smaller when both t 1 and h are large. The same is observed for β LLg = 0.3 with a maximum deviation of 50%. The inset to the lower right shows the corresponding plot for β LLg = 0.6. Also in this case there is a very good match for h/t 1 10 −2 , and the data deviate systematically towards smaller values for larger h/t 1 . The maximum deviation reaches a factor of five for h/t 1 10 2 , and we see the same trend, that it becomes smaller again for still larger values of h/t 1 .
In accordance with our expectation for a correction to scaling the deviations differ for different times t 1 at a fixed ratio h/t 1 , as indicated by the different colors of the symbols. In contrast, according to Eq. (10) the corrections to scaling in the FnD dynamics approach one from above when t 1 is increased.
For all data sets an important qualitative prediction of the FnD model holds: The displacement autocorrelation function of the LLg obeys the scaling form Eq. (4). Moreover, the FnD expression, Eq. (11), provides a rather accurate, parameter-free description of the two-variable correlation function φ(t 1 , t 2 ) in terms of the time ratio h/t 1 . There are deviations in particular at intermediate times. They constitute a nonuniversal correction to Eq. (11a), and they behave qualitatively differently in the LLg and the FnD dynamics.

D. The Lorentz gas (LG)
We have numerically computed the autocorrelation function of the Lorentz gas φ(t 1 , t 2 ) = x(t 1 ) x(t 2 ) , where · refers to an average over an ensemble of 1.8 × 10 8 trajectories that is evaluated for 200 t 1 5 × 10 5 and different time increments h = t 2 − t 1 ∈ {10, 50, 100, 500, 1000, 5000}. The results for the infinite-horizon case, R = 0.1, are shown in Fig. 5, where time t 1 is marked by color and h by symbols, as indicated in the figure caption. As for the LLg, the data accurately collapse to a line when the reduced correlation function φ(t 1 , t 2 )/φ(t 1 , t 1 ) − 1 is plotted as a function of h/t 1 . The mean-square displacement grows with exponent η = 1, and the prediction Eq. (11a) is provided by the thick solid black line. The parameter-free prediction provides values that are smaller by about one order of magnitude than the observed data.
The ratio of the data and the prediction is provided in the upper left inset: The deviation for h/t 1 1 amounts to leading order to a logarithmic law. Encountering logarithmic corrections to scaling is not unexpected for the LG, where they also arise in the mean-square displacement [54,59,60]. Likewise, we see here numerical evidence for a logarithmic contribution, ln(h/t 1 ), to the asymptotic scaling function. On top of that there are noticeable finite-time corrections. The data for different time increments h bend away from the logarithmic line, taking a minimum towards right for t 1 5000, and they increase for larger values of h/t 1 . Similarly to the finite-time corrections of the FnD, the asymptotic scaling is approached from above.
For h/t 1 1 the FnD prediction, that is shown by the thick solid line in the main panel, looks qualitatively different than the numerical data. Rather than following the η = 1 case of Eq. (11a) the numerical data are better described by a power law. A guide to the eye is provided by the thin black line that shows the dependence 1.75 (h/t 1 ) 0.69 . The ratio of the data and this power law is shown in the lower right inset. It reveals that the fit refers only to the final decade of the h = 5000 dependence. Much longer simulations are required, therefore, to make qualified statements about the large h/t 1 asymptotics of the correlation function. In particular, when the correlation decays for vastly different times t 1 and t 2 we would again expect that the scaling function must eventually approach −1 in the large-time limit, 1 t 1 t 2 . We conclude that the important qualitative prediction, Eq. (4), of the FnD model holds also for the LG: In the large-time limit the ratio of the displacement autocorrelation function and the mean square displacement is a function of only h/t 1 . However, in this case the FnD prediction of the scaling function, Eq. (11), provides only a rough idea of the form of the scaling function. It misses logarithmic terms in the small h/t 1 limit, and for large h/t 1 it is off even to leading order.

E. Polygonal billiard channels (PBC)
For PBCs we have numerically computed the autocorrelation function φ(t 1 , t 2 ) for the polygonal channel with infinite horizon H = 1.27 and with finite horizon H = 1.07. The results are shown in the left (H = 1.27) and right (H = 1.07) panel of Fig. 6. We show data for the same combinations of h and t 1 as adopted for the LG and also use the same color coding and symbols. The scaling form provides a data collapse also for these billiards.
For both channels the mean-square displacement scales with η = 1, i.e., with the same exponent as the LG. Hence, the FnD model suggests the same scaling function for these three types of vastly different billiards. However, the data differ for the different billiards.
The data for the infinite-horizon case, H = 1.27, deviate by a factor of 3.5 ± 1 from the FnD prediction, as indicated by the horizontal dashed gray line in the upper left inset in the left panel. Moreover, for small h/t 1 the infinite horizon PBC also shows logarithmic corrections to the FnD prediction, analogous to the correction observed for the LG with infinite horizon. Also for large h/t 1 the data do not follow the logarithmic law predicted by Eq. (11), but they are better described by a power law with an exponent 0.57 that is shown by a thin solid black line in the main panel. The ratio of the data and this power law is shown in the lower right inset. The data for different times t 1 collapse considerably better in this case than in the LG, and they all nicely approach the asymptotic power law. However, also in this case much longer simulations are required to make qualified statements about the large h/t 1 asymptotics of the correlation function.
For H = 1.07 the data lie closer to the FnD prediction. For h/t 1 10 −2 there is a fixed factor of 1.5 between the data and the FnD prediction. This is indicated by a horizontal dashed gray line in the upper left inset, which shows the ratio of the data and the FnD prediction. For the PBC with a finite horizon the FnD prediction therefore provides the right scaling of the data, and it is off in the prefactor of the scaling by 50%. This indicates that the subleading logarithmic contribution of the scaling function for h/t 1 1 is indeed a property of the infinite horizon.
In the range 0.002 < t 1 < 0.2 the offset of the finitehorizon PBC increases to a factor of 3.5 that is indicated by another dashed line in the inset. For still larger h/t 1 the scaling functions of both PBCs coincide. This suggests that the decay of correlations for h/t 1 1 follows the same dependence for PBC with finite and infinite horizon.
We conclude that the important qualitative prediction, Eq. (4), of the FnD model holds also for the PBC: In the large-time limit the ratio of the displacement autocorrelation function and the mean-square displacement is a function of only h/t 1 . In this case the FnD prediction of the scaling function, Eq. (11), provides an accurate description of the small h/t 1 dependence of the billiards with finite horizon, up to a constant factor of 1.5. This agreement is remarkable because the FnD dynamics does not even describe the mean-square displacement of the PBC with finite horizon. This suggests that for h/t 1 1 the scaling form, Eq. (11), is a generic feature of strong anomalous diffusion.

V. DEVIATIONS FROM SCALING
The scaling form, Eq. (4), of the displacement correlation function, Eq. (3), provides a data collapse for all dynamics. For the LLg with small parameter values β LLg we find quantitative agreement between our numerical data and the parameter-free prediction of the FnD model over the full range of the dimensionless time difference h/t 1 , that is varied over 14 orders of magnitude for β LLg = 0.1 (uppermost curve in Fig. 4), for 13 orders of magnitude for β LLg = 0.3 (middle curve) and 11 orders for β LLg = 0.6 (bottom curve). This data closely follows the FnD prediction, showing systematic deviations from scaling only at intermediate values of h/t 1 . We attributed these deviations to corrections to the scaling prediction for finite values of t 1 . For the FnD model such deviations have been established in Eq. (10). These corrections for nonasymptotic regimes are nonuniversal, as may be expected. For the FnD model finite-time data are larger than the asymptotic values, for the LLg we find smaller values. For small values of h/t 1 the numerical data agree quantitatively with the FnD prediction, and the quantitative differences are small over the full range of h/t 1 : The differences increase from 10% for β LLg = 0.1 to 60% for β LLg = 0.6.
For the LW the scaling prediction, Eq. (4), provides a data collapse without appreciable finite-time corrections. For h/t 1 0.01 the FnD model even provides a quantitative description without adjustable parameters. The upper two curves in Fig. 3 demonstrate that this even applies in the ballistic regime. For large h/t 1 the LW data lie below the FnD prediction, with increasing deviations for larger β LW . Comparison with the scaling function for the LW, Eq. (25), reveals that the scaling form, Eq. (4), holds over the full range of h/t 1 , and that the scaling functions of the LW and of FnD match for small h/t 1 .
For the LG with infinite horizon the data show a collapse with a scaling function that follows the trend predicted by FnD (upper curve in Fig. 5). However, there is no quantitative agreement of the scaling functions in this case. For small h/t 1 we identified logarithmic corrections that we expect to persist in this system (dashed line in the upper left inset). Such corrections are in accordance with the logarithmic contribution to the scaling of the statistics of displacements [54]. However, they take the qualitatively new form of a log(h/t 1 ) dependence. For large h/t 1 the FnD model predicts a logarithmic dependence, while the data are better described by a power law. The correlation φ(t 1 , t 2 ) does not decay to the ensemble average if such a power law describes the h/t 1 asymptotics. After all, due to the left-right symmetry of the LG the steadystate equilibrium value of x(t 1 )x(t 2 ) must vanish in the limit 1 t 1 t 2 . The very long transients of the LG [59,60] make it impossible for us to explore this crossover.
The PBC with infinite horizon also enjoys a faithful data collapse with very small scatter for data evaluated at finite times. All data lie in a band that differs only by a factor of 3.5 ± 1 from the FnD prediction (left panel in Fig. 6). Still, the deviations are systematic, of exactly the same type as for the LG. For small h/t 1 the PBC with infinite horizon shows logarithmic corrections. For large h/t 1 the data are better described by a power law than by the logarithmic dependence predicted by the FnD model. The data for the PBC with finite horizon also show a faithful data collapse, with noticeable deviations from the master curve only for the smaller time increment h = 10. For h/t 1 10 −2 the data follow the FnD prediction up to a factor 1.5. Subsequently, they cross over to the dependence also observed for the PBC with infinite horizon. This agreement for large values of h/t 1 indicates a common mechanism for the decay of correlations in the PBC; this is characterized by the autocorrelation function for t 2 t 1 1. It is remarkable that the FnD model does not only suggest a scaling form, Eq. (4), for the displacement correlation function that is followed by all models investigated here, but that there even is quantitative agreement between the predicted scaling function and the numerical data for the LW, LLg, and the PBC with finite horizon. We suggest that this can be explained based on the Buckingham-Pi theorem [62,63]. To this end we consider the following choices of scales and dimensionless parameters: FnD: The time scale is set by the minimum flight time t M = b 1/ξ , Eq. (5b) and below, while space and time scales are related by the unit velocity of the particles. There are no dimensionless groups in this model. SM: The length scale is set by the box size, while space and time scales are related by the unit velocity of the particles. There are no dimensionless groups in this model. LW: The length scale is set by the minimum length of trajectory segments, while space and time scales are related by the unit velocity of the particles. There are no dimensionless groups in this model.
LLg: The length scale is set by the minimum distance of scatterers, while space and time scales are related by the unit velocity of the particles. There are no dimensionless groups in this model.
LG: The length scale is set by the scatterer separation while space and time scales are related by the unit velocity of the particles. The ratio of scatterer distance and scatterer radius provides a dimensionless parameter /R (cf. the inset of Fig. 1).
PBC: The length scale is set by the periodicity 2 x of the channel, while space and time scales are related by the unit velocity of the particles. Three dimensionless groups are required to fully characterize the channel, y t / x, y b / x, and H/ x (cf. the inset of Fig. 2).
For small h/t 1 the FnD model provides a prediction, Eq. (11), on the exponents, while the prefactors of the power laws will in general depend on the dimensionless groups of the models [62,63]. From this perspective, we expect full agreement of the FnD, SM, LW, and the LLg, while the prefactors of the power laws are expected to depend on the system geometry for the LG and PBCs. The scaling behavior for large h/t 1 is not expected to be universal. It will depend on the relaxation to the nonequilibrium distribution and may suffer from breaking of ergodicity and persistent transients.

VI. CONCLUSIONS
The moments of the displacement of many systems that show strongly anomalous superdiffusive transport are dominated by ballistic trajectories or light fronts [18,37,38]. Often they lead to a two-piece linear scaling of the exponent γ (p) of the temporal growth of the pth moment of the probability distribution of the displacement, Eq. (2). Here, we presented the very simple fly-and-die (FnD) model that arguably is the simplest dynamics in this class of systems. Due to the simplicity of its dynamics one can derive analytical expressions for the exponents γ (p) and the displacement autocorrelation function φ(t 1 , t 2 ). Based on the FnD model we have proven that the correlation function follows a universal scaling, Eq. (4), for a range of systems showing strong anomalous diffusion. For all systems where the second moment is governed by a light front the ratio φ(t 1 , t 2 )/φ(t 1 , t 1 ) is a function 1 + C(h/t 1 ) of a single argument, h/t 1 = (t 2 − t 1 )/t 1 with t 2 > t 1 . We thus relate the scalings of the position-position displacement correlation function and of the mean-square displacement, connecting the most widely studied and best characterized property of these systems to a feature that is not known for most strongly anomalous systems. The prediction applies for all systems mentioned in the Introduction.
In the second part of the paper we underpinned our claim by discussing representative models from five classes of widely studied systems: (i) the slicer map, a deterministic, time-discrete dynamical system, (ii) Lévy walks, where strong anomalous diffusion arises from following a route of straight line segments whose length are sampled from a distribution with a power-law tail [16,36], (iii) the Lévy-Lorentz gas, a random-walk model featuring quenched disorder [29], (iv) the Lorentz gas, a chaotic billiard with infinite horizon [22], and (v) polygonal billiard channels, where trajectories only separate when they hit different walls of the polygon [21]. For the billiards we considered systems with finite and infinite horizon.
The FnD model makes two predictions that are robust in the sense that they apply universally to all systems: (1) When the offset value (1 − ν) p c in Eq. (2) is fixed for one of the moments in the large p regime, then the scaling for all higher moments follows without further adjustable parameters. Microscopic details of the dynamics are reflected solely by different values of ν. Specifically, ν = 0 for the FnD dynamics and the SM, ν = 1/2 for the billiard systems, and ν depends on the parameter β LW for the LW and β LLg for the LLg [cf. Eqs. (16) and (18)].
(2) For the FnD dynamics we derived a scaling form, Eq. (11), of the displacement correlation function, Eq. (3). It applies also to the SM, whose correlation function can also be calculated analytically. Due to the peculiarities of the FnD and the SM one cannot expect that the functional form of Eq. (11) is universal. However, it suggests that in general the autocorrelation functions admits a data collapse of the form Eq. (4). This is indeed true for the LW that features the predicted scaling (cf. Fig. 3), but with a different functional form, Eqs. (25). In addition to the LW we tested the data collapse with numerical data for the LLg model with data spanning to the least eleven orders of magnitude in the reduced dimensionless time (Fig. 4), and more than six orders of magnitude for the LG with open horizon and PBCs with open and closed horizon (Fig. 6). For all investigated models the suggested representation of the data provides a data collapse. For the LW and LLg the FnD prediction even provides a parameter-free prediction of the prefactors of the asymptotic scaling laws for small h/t 1 . For the billiards the FnD prediction faithfully describes the trends and the order of magnitude of the correlation function. In Sec. V we related the offset to the presence of dimensionless parameters in the billiards that affect the prefactors of the predicted power-law dependences. Special care is needed in situations where the crossover, p c , of the branches of the spectrum γ (p) arises at the exponent p = p c = 2 of the mean-square displacement. This is the case for billiards with infinite horizon. For the LG and PBC with infinite horizon we observed logarithmic corrections to the scaling prediction of the FnD dynamics that emerge for small h/t 1 . Moreover, the case η = η(2) = 1 is a critical dimension where the FnD takes a logarithmic dependence for large h/t 1 . All considered billiard systems fall into this class, but their scaling in this range can better be fitted by a power law than by a logarithmic dependence. More analytical and numerical work will be needed to fully understand these dependencies. They will be addressed in forthcoming work that focuses on the intriguing features of billiard systems, rather than addressing universal features of correlation functions.
In Eq. (10) we also provided the small t 1 corrections to scaling for the FnD model. These corrections take a nonuniversal, system-specific form. For the FnD model the asymptotic theory is approached from above. The LW data show no noticeable finite time corrections. The LLg data show deviations towards smaller values. For the LG there are noticeable deviations towards larger values. The PBC show deviations towards smaller values only for the smallest considered time difference, h = 10.
We conclude that the FnD model captures the essential features of strong anomalous transport with two underlying scales. Only the presence of a light front and a much slower dynamics are kept, and all surplus features are removed. In the spirit of normal forms in mathematical models our model only accounts for a single long jump. Its simplicity entails a straightforward analytical solution. The present paper scrutinized the resulting predictions for the moments of the displacement and the displacement autocorrelation function. Vastly different dynamics follow the FnD predictions for the high moments of the displacement, Eq. (1) for large p where γ (p) has slope one, and the long-time asymptotics of the displacement autocorrelation function, Eq. (3) for times 1 t 2 2 t 1 2t 2 . Moreover, the FnD establishes a scaling form for the correlations, Eq. (4). In the longtime regime, 1 t 1 t 2 , it provides a data collapse for all investigated models, with model-specific scaling functions for t 2 t 1 . In the opposite limit t 2 − t 1 t 1 t 2 the moments and the correlations are governed by the light front. As a consequence, the parameter-free FnD prediction is followed quantitatively, except that prefactors of scaling laws may differ: Length and time scales shift if and only if a model has dimensionless groups not addressed by the normal form. This finding applies to all systems with strong anomalous diffusion where the second moment is governed by a light front, and it thus calls for experimental verification in a vast range of different physical settings.