Universality in Statistical Measures of Trajectories in Classical Billiard Systems

For classical billiards, we suggest that a matrix of action or length of trajectories in conjunction with statistical measures, level spacing distribution and spectral rigidity, can be used to distinguish chaotic from integrable systems. As examples of 2D chaotic billiards, we considered the Bunimovich stadium billiard and the Sinai billiard. In the level spacing distribution and spectral rigidity, we found GOE behaviour consistent with predictions from random matrix theory. We studied transport properties and computed a diffusion coefficient. For the Sinai billiard, we found normal diffusion, while the stadium billiard showed anomalous diffusion behaviour. As example of a 2D integrable billiard, we considered the rectangular billiard. We found very rigid behaviour with strongly correlated spectra similar to a Dirac comb. These findings present numerical evidence for universality in level spacing fluctuations to hold in classically integrable systems and in classically fully chaotic systems.


Introduction
The idea to model apparently disordered spectra, like those of heavy nuclei, using random matrices was suggested in the mid-50's by Wigner, and then formalized in the early 60's in the work of Dyson and Mehta [1]- [6].They showed that random matrices of Gaussian orthogonal ensembles (GOE) generate a Wigner-type nearestneighbour level spacing (NNS) distribution [6].In a seminal paper, Bohigas, Giannoni and Schmit (BGS) formulated a conjecture [7] stating that time-reversal invariant quantum systems with classically fully chaotic (ergodic) counterpart have universality properties given by random matrix theory (RMT).Experiments in nuclear physics, for example, have shown that spectra originating from different heavy nuclei give the same Wignerian energy level spacing distribution [8].Universality properties in quantum chaos of bound systems, i.e. quantum systems with a fully chaotic classical counterpart, have now been demonstrated in many experiments, computational models and in theoretical studies [6] [9]- [14].Theoretical support of the BGS conjecture came from the semiclassical theory of spectral rigidity by Berry [15] [16], who showed that universal behaviour in the energy level statistics is due to long classical orbits.Sieber and Richter [17] investigated the role of correlations between classical orbits.The semiclassical theory has been further developed by Müller et al. [18] [19].In ref. [19] they presented the "core of a proof" of the BGS conjecture, which, by proving arguments previously used by Berry [15], show that in the semi-classical limit the periodic classical orbits determine the universal fluctuations of quantum energy levels.Further refinements have been made by Keating and Müller [20].
The study of classical strongly chaotic systems (Anosov systems) has revealed that central limit theorems (CLT) hold [21]- [26].This has been proven for the 2D periodic Lorentz gas with finite horizon.The first step of a proof was given by Bunimovich and Sinai [21] and was completed by Bunimovich, Sinai and Chernov [23].At macroscopic times, such deterministic chaotic system converges to Brownian motion, i.e. behaves like a random system [21] [23] [24] [26].This is also supported by the existence of an average diffusion coefficient [21] [23] [27], showing the diffusive character of such chaotic system.Do classical fully chaotic systems also exhibit universality properties?This question was addressed by Argaman et al. [28], who showed that there is universal behaviour in 2 points correlation functions of actions of periodic orbits of classically chaotic systems.As examples they considered the deformed cat-map and the baker-map.This has been elaborated further in a number of studies by Dittes et al. [29], Aurich and Sieber [30], Cohen et al. [31], Tanner [32], Sano [33], Primack and Smilansky [34], Sieber and Richter [17] and Smilansky and Verdene [35].Argaman et al. started from the assumption that spectral fluctuations of chaotic quantum systems follow the predictions of RMT and they derived a universal expression for classical correlation functions of periodic orbits via Gutzwiller's semi-classical trace formula.They concluded "The real challenge, though, is to find out whether these action correlations can be explained on a completely classical level".
An answer was proposed by Laprise et al. [36] [37] who found universal behaviour in classical 2D billiards by looking at fluctuations in spectra of classical action/length matrices from billiard trajectories.They showed that one could distinguish chaos from integrability in classical systems using RMT and an analogue of the BGSconjecture.In particular, they considered the Limaçon/Robnik family of billiards, which interpolates between the chaotic cardiod billiard and the integrable circular billiard.For the cardioid billiard, the level spacing distribution ( ) P s and spectral rigidity ( ) 3 L ∆ were found to be consistent with the GOE behaviour predicted by RMT.For the interpolating case close to the circle, the behaviour approached a Poissonian distribution.The circular billiard itself was found to be very rigid and strongly correlated and yielded ( ) ( ) − .The jump in behaviour at the transition to the circle is associated with the corresponding change in the symmetry group.
This article extends the results of reference [36] in the following directions: 1) We consider the 2D rectangular billiard as another example of an integrable billiard.Compared to the circular billiard, this billiard has lower symmetry (no group property).Nevertheless, it displays strong spectral correlation and rigidity like the circular billiard.2) We present numerical studies for other chaotic 2D billiards: the Sinai-billiard and the Bunimovich stadium billiard.3) In order to understand the observed universal behaviour in chaotic billiards, we present arguments linking such behaviour to CLTs, diffusive and random walk behaviour.In particular, we present a mathematically rigorous result on the distribution of length of trajectories.
The answers we found can be summarized as follows: For the 2D Sinai billiard and the 2D Bunimovich stadium billiard the level spacing distribution ( ) P s and spectral rigidity ( )

L ∆
are consistent with predictions of RMT (GOE behaviour), i.e. show universal behaviour.This behaviour is statistically the same as the one observed in quantum chaos (obtained from energy level spacing distributions).The implication of these findings is that RMT not only represents well the statistical fluctuation properties of the energy spectrum of chaotic quantum systems, but also those of the length spectrum of chaotic classical systems.Moreover, statistical fluctuations obtained from spectra of action/length matrices clearly distinguish chaotic from integrable systems.

Length and Action of Trajectories
In classical systems, chaos information is encoded in trajectories.According to the Alekseev-Brudno theorem [38] the temporal length t is related via the Kolmogorov-Sinai K entropy to the information ( ) ( ) denote the action over ( ) denote the length of the trajectory ( ) . For all pairs of boundary points , i j x x X ∈ , we compute a classical trajectory traj , i j x , connecting those points.We suggest the construction of an action matrix and a length matrix, traj , where i and j are respectively the indices of the final and initial boundary points of the trajectory.Both matrices are viable for statistical analysis of classical chaos.In the case of billiard systems, we consider trajectories where the billiard particle moves with constant velocity u and constant kinetic energy E. Then the action and the length matrix are essentially equivalent, .

Numerical Calculation of Length Matrices
When solving for two boundary points i x and j x in a billiard the solution is generally not unique; the number of possible trajectories varies with the number reb N of rebounds.We have thus classified trajectories-and corresponding length matrices-according to reb N .For reb 15 N = the number of trajectories is quite large ( ) 5 traj 10 N ∼ . In order to limit computational cost, we considered smaller subsets of length matrices, in the following way.For a given number of rebounds reb N and a given pair of boundary points { } , i j x x , we determined the starting angles ν α (measured from the normal to the billiard wall at the boundary point) which corresponds to a trajectory ν .That search was carried out in the  .We repeated this for all combinations of boundary points.We then checked that this procedure did not spoil the statistical properties we aimed to measure by considering different cutoff traj n and different ordering scheme for the set of trajectories.

Integrable Billiard
If one considers integrable quantum systems and analyzes them in terms of the NNS distribution of energy levels and spectral rigidity, then in most cases one finds a Poissonian distribution ( ) P s which is reflected also in the behaviour of ( ) . However, there are examples of integrable quantum systems in 2D where the levels are correlated and the level spacing distribution is not given by a Poissonian.Such cases are called non-generic.Berry and Tabor [39] noted as example uncoupled quantum oscillators in 2D.Casati et al. [40] and later Seligman and Verbaarschot [41] showed that also the integrable 2D quantum rectangular incommensurate billiard does not give an uncorrelated Poissonian level spacing distribution.In a recent study of 2D quantum harmonic oscillators Chakrabarti and Hu [42] found in the case of uncoupled oscillators a level spacing distribution displaying a ( ) .This behaviour is similar to the rigidity of a picket fence (Dirac comb) of equally spaced levels with ( ) They conclude that the 2D quantum oscillator system is highly correlated at short and long- range, is regular and very rigid.
Laprise et al. [36] considered the classical integrable circular billiard and constructed a length matrix from classical trajectories between boundary points located evenly on the billiard wall.They found saturation in the spectral rigidity for large L at

Rectangular Billiard
As example of an integrable classical billiard we consider the 2D rectangular billiard, shown in Figure 1.
The shape is determined by the parameters a, b, which were chosen to be 5 a = and 1 b = .The boundary , located in the upper right quarter of the billiard wall, are distributed regularly, with perimeter spacing given by ( ) ( ) For a given pair of boundary points, we found that the behaviour of the number of trajectories versus the number of rebounds is linear (not shown).The error behaviour of trajectories as a function of the number of rebounds has been obtained by taking into account traj 20 N = trajectories corresponding to different initial conditions.The error was stable in the regime reb 0 100 N < < and did not increase beyond 16 10 − .This is in contrast to an exponential increase found in chaotic billiards.
The symmetry of the rectangular billiard is mirror symmetry under reflection about the x-and y-axes (with origin at centre of rectangle).The symmetry shows up in the shape of trajectories.For example, a trajectory (1) going from starting point N x to endpoint 1 N x − bouncing once from the lower boundary wall has the same shape as the trajectory (2) going from starting point x − bouncing once from the lower boundary wall.Trajectory ( 2) is obtained by a translation of trajectory (1) in x-direction by the amount of perimeter spacing, Equation (6).However, such discrete translations do not form a group (they would form a group if one would consider the rectangular billiard with periodic boundary conditions at all billiard boundary walls, which would map the billiard on a torus, and the symmetry group would then be a rotation group).The length matrix is not exactly a circular matrix either, however, it shares with circular matrices the property of having a number of different matrix elements equal to or less than the rank of the matrix.This property implies strong correlations between length matrix elements, which translates to strong correlations of eigenvalues of the length matrix.To sum up, although there is no exact symmetry group, there are residues of a "broken symmetry  .This is consistent also with the NNS level distribution ( ) P s of the length matrix shown in Figure 2(c) where one observes ( ) ( ) indicating the behaviour of a Dirac comb.The rectangular billiard turns out to be highly correlated at short and long range, to be very rigid and to be regular.This is possibly evidence for universal behaviour in the integrable case.Comparing the behaviour of the rectangular billiard with the circular billiard [36], we observe that they differ in their symmetry properties.In the circular billiard hopping from one boundary point to its neighbour stands for a group operation.The corresponding operation in the rectangular billiard has no group property.However, the resulting strong correlation and spectral rigidity are found to be very similar for both billiards.

Sinai Billiard
For general closed 2D billiards, the mean free path length free Λ in between two collisions is given by the billiard geometry via [43] free π , where Q stands for the billiard area and Q ∂ for the length of the wall.Let us now consider a particle of mass m moving in the 2D Sinai billiard (see Figure 3(a)).The classical Sinai billiard system is known to be fully chaotic [44].The billiard is symmetric under mirror operation about x-and y-axes, and about the diagonals of angles π 4 and 3π 4 , respectively, passing through the centre.In order not to mix different symmetry classes, we consider the billiard with boundary points 1 , , N x x  located in one eighth of the exterior billiard wall (see Figure 3(a)).It is important to make sure that the chaotic behaviour comes from the dynamics of the system and not from random location of boundary points.Thus we have located the boundary points in a regular way, with perimeter spacing The rule of dynamics is free motion in the interior region and elastic specular collision at the central disc and the exterior square wall.We have classified trajectories using the scheme presented in Section 2. A global For reb 48 N > saturation is reached (the relative error is in the order of one).Such exponential behaviour (see Figure 3(b)) is related to the largest positive Lyapunov exponent of the system and thus represents a global characteristic of a chaotic system.Like the exponential behaviour of number of trajectories, Equation ( 9), also the exponential error behaviour distinguishes between chaotic billiards and integrable billiards.
We found a relative error of about 13 10 − for reb 15 N ≤ .After unfolding the spectra we retain 8 significant digits.We used the technique of Gaussian broadening [13] [45] to unfold the raw spectrum.In this method there is a free parameter, which we tuned to reproduce the auto-correlation coefficient ( 0.27 C = − for Wigner distribution [8] [46]).For a given trajectory index ν we computed a level spacing distribution.Afterwards we superimposed the spectra corresponding to different trajectory indices ν .For reb 15 N = the number of trajectories is quite large ( ) 5 traj 10 N ∼ . In order to limit computational cost, we have restricted ourselves by considering smaller subsets of length matrices.It turns out that traj 20 n = gives a reasonably good statistics.Note that, in this way, we don't consider all trajectories and have to order them arbitrarily as explained in Section 2. Therefore, we verified that the resulting distributions are independent of the particular subset of trajectories and ordering used within statistical errors.The resulting level spacing distribution ( ) P s , and are consistent with GOE behaviour also for ( ) . As in the quantum Sinai billiard, the level spacing distribution has been shown to give a Wigner distribution [7].
The histogram of the length matrix elements itself is shown in Figure 4(a).The distribution looks close to a Gaussian.Determining if it is a pure Gaussian, is a question physically relevant for the following reason: If the distribution ( ) P Λ is a Gaussian, then the matrix elements obey GOE statistics.Then RMT [6] implies that the level fluctuation statistics  P Λ .The random walk model is mathematically a Markov chain, i.e., the length of each piece of straight path in between subsequent collisions is given by a random number.Two subsequent random numbers are statistically independent.On the other hand, the chaotic Sinai billiard is a deterministic system, and the length of two subsequent pieces of straight trajectory are not independent.Hence it is plausible that the distribution ( ) P Λ of length matrix elements for the Sinai billiard is not given by an exact Gaussian.
Mathematical note.The BGS-conjecture does not state that the matrix elements of a quantum Hamiltonian must be distributed like a GOE ensemble.The conjecture rather only says that the statistical fluctuations of the eigenvalue spacings obtained from the quantum Hamiltonian are the same as those from a GOE ensemble, giving a Wignerian distribution.In other words, it is possible that the matrix elements of the quantum Hamiltonian be distributed quite differently from a Gaussian and that its level spacing distribution be nevertheless Wignerian.Such a situation, where the distribution of matrix elements is not GOE, but the level fluctuation statistics is GOE, occurs in nuclear physics.An example is the distribution of the Hamiltonian matrix elements obtained from nuclear shell model calculations [12].In this model, there are vanishing Hamiltonian matrix elements.This implies that the number of independent matrix elements is much smaller than in a random matrix of the same size.However, Brody et al. [47] showed that the 2-body residual interaction in the shell model yields matrix elements of random character following a Gaussian distribution.In particular, they showed that spectral fluctuation properties from such ensembles with orthogonal symmetry are identical to those from GOE.That implies that GOE is meaningful to predict spectral fluctuation properties of nuclei governed by 2-body interactions, though the Hamiltanian does not follow a Gaussian distribution.

Bunimovich Stadium Billiard
Let us consider the 2D Bunimovich stadium billiard with semi-axes a and b (Figure 5).The billiard is known to be fully chaotic [48] [49].The billiard is symmetric under mirror operation about x-and y-axes.In order not to mix different symmetry classes, we consider the billiard with boundary points 1 , , N x x  located in one quarter of the billiard wall.The quarter perimeter has the length ( ) . The nodes n x are equidistantly distributed on the quarter perimeter of the stadium, with perimeter spacing given by The rule of dynamics is free motion in the interior region and elastic specular collision at the exterior square wall.For a given pair of boundary points, the number of trajectories traj N connecting these points increases on average with the number of rebounds exponentially, according to ( ) Such behaviour is qualitatively similar to that found in the Sinai billiard.We also measured the numerical error following the method used in the Sinai billiard.Likewise, we found a regime of exponential behaviour followed by a regime of saturation (not shown).On average the exponential increase is given by ( ) reb reb exp 14.3 , 0.247 0.005, 58.
We have classified trajectories by the the number of rebounds reb N , which was kept fixed.In this way, we generate an ensemble of length matrices { } In order to see if the observed Wigner distribution in the level spacing distribution depends on the statistical method of averaging over several trajectories (i.e.length matrices ν Λ ) we have tested an alternative, where we averaged over different stadium geometries.All stadium billiards are fully chaotic, and if the Wigner behaviour is universal, this should show up for any stadium geometry.Superimposing spectra corresponding to different stadium shapes should improve statistics.In their paper on quantum chaos in nuclear physics, Bohigas et al. [50] have analyzed spectra of a series of different heavy nuclei (but having the same quantum numbers) and showed that their level fluctuation statistics agrees with the prediction of RMT.These different nuclei correspond to different potentials.In analogy to that we considered here length spectra corresponding to different stadium shapes.The stadium wall represents a curve where the potential jumps from zero to infinity.Hence different stadium shapes correspond to different potentials.By superimposing spectra from different stadium shapes we obtained results displayed in Figure 7. Within statistical errors the results are equivalent to those obtained by superimposing spectra from different trajectories (Figure 6).

Universality in Chaotic Billiards
The numerical experiments with chaotic billiards investigated above led us to the following observations: 1) For a given pair of boundary points reb N , the number of classical trajectories traj N shows exponential behaviour as function of reb N .2) For a given pair of boundary points, the numerical error ε shows exponential behaviour as function of the number of rebounds reb N .
3) The NNS level fluctuation distribution ( ) P s obtained from eigenvalues of the length matrix Λ shows a Wignerian distribution.This is consistent with GOE behaviour predicted by RMT. 4) The spectral rigidity ( ) obtained from the spectrum of length matrix Λ is also found consistent with GOE behaviour predicted by RMT. 5) The distribution ( ) Λ is found to be close to a Gaussian.
The leading Gaussian behaviour in the distribtion ( ) P Λ is possibly universal.How can we understand such behaviour of ( ) P Λ in chaotic billiards?The Sinai billiard is equivalent to a system where the billiard ball moves in an open system of equal spherical discs located on a 2D rectangular regular grid, known as the 2D (diluted) Lorentz gas model.The billiard ball alternates between free motion and collisions with the circular discs.There are two types of 2D periodic Lorentz gases (or Sinai billiards on a torus): One has a horizon, where free paths between collisions are bounded (the scatterers are sufficiently dense to block every direction of motion.The other type has an infinite horizon, where the particle can move indefinitely without collision with any disc.The Sinai billiard investigated above (see Figure 3(a)) belongs to this class.The CLT, the existence of finite diffusion coefficient and the convergence to Brownian motion were proven only for the periodic Lorentz gas with finite horizon [21] [23].In the case of infinite horizon, the CLT has not been proven.Moreover, in this case the diffusion coefficient is infinite and there is no convergence to Brownian motion [51].
Actually, for the periodic Lorentz gas with finite horizon it can be shown rigorously that the distribution of length of trajectories becomes a Gaussian distribution in the limit of many bounces.This holds when the initial points of trajectories are distributed randomly on the billiard boundary.Then Chernov and Markarian [52] prove the following result for all z −∞ < < ∞ .Here ( ) x τ denotes the return function and ( ) ( ) ( ) is the time of the n-th collision. 2f σ is the variance of the random variable depending on the observable f, which here is the return function ( ) is the collision map and M is the collision space (phase space of billiard map).Thus the distribution of the travel time until the n-th collision obeys the CLT, and hence converges in distribution to a normal, i.e.Gaussian distribution in the limit where the number of collisions goes to infinity.In the case where the billiard particle moves with constant speed u, this means also the length of trajectory n n ut Λ = obeys the CLT, and the length distribution ( ) n P Λ tends to a normal Gaussian distribution for n → ∞ .
This result seems to support our numerical findings of a (near) Gaussian distribution of length of trajectories.However, the scenario where the above mathematical result holds and the scenario of our numerical study differ in two respects, namely in the distribution of boundary points and in the horizon of billiard.For the purpose of statistical analysis in terms of RMT we are interested in the distribution ( ) n P Λ , which corresponds to the case where the billiard particle starts from and arrives at boundary positions taken from a discrete set of boundary points, which is different from an initial random set of boundary points.Nevertheless, it seems plausible that a similar mathematical result may hold also for the discrete set of initial and final boundary points, such that ( ) tends to a Gaussian.Second, the mathematical result holds for the case of finite horizon, while the Sinai billiard in our numerical study has an infinite horizon.At a first glance, it looks puzzling in absence of strong mathematical results (for infinite horizon) that our numerical results nevertheless show (near) Gaussian behaviour for the length distribution.A possible explanation is this: The choice of the discrete set of regularly distributed boundary points selects preferentially trajectories which often bounce at the central disc, while trajectories with infinite horizon (no bounces at disc) are avoided.This is supported also by numerical evidence that the distribution ( ) n P Λ of the Sinai billiard (Figure 4(a)) and that of the rectangular billiard viewed as Sinai billiard with central disc of radius zero (Figure 2(a)) are quite different.
In the case of the stadium billiard, Bálint and Gouëzel [53] have shown that a limit theorem also holds.They proved that the limit distribution ( ) n → ∞ is Gaussian, however, the scaling factor is not given by the standard expression n , but rather by log n n .As consequence, there is no finite standard diffusion coefficient and no convergence to Brownian motion.Bálint and Gouëzel prove that the distribution of length of trajectories, tends to a gaussian distribution in the limit where the number of collisions n goes to infinity.Like in the Lorentz gas, this corresponds to a random distribution of initial points on the boundary.This is consistent with our numerical results of the stadium billiard which show in ( ) P Λ small deviations from a Gaussian comparable in magnitude with those in the Sinai model.They are likely due to the small number of bounces.
Let us suppose that ( ) n P Λ does tend to a Gaussian.As a consequence, in this limit (and after suitable normalization) the length matrix Λ itself will become a GOE matrix [11].This is sufficient to guarantee in this limit that the level fluctuations ( ) P s and spectral rigidity ( ) 3 L ∆ obtained from the matrix Λ show GOE behaviour [6].This does establish universal behaviour in the limit of many bounces in the case of the rectangular Lorentz gas/Sinai billiard as well as the Bunimovich stadium billiard.Now we want to address the following questions: 1) Concerning universality observed in classical chaotic billiards, what are the underlying physical principles?We will give a heuristic description-not a rigorous derivation-of the physical principles leading to the phenomenon of universality.Let us consider chaotic billiards in the regime of macroscopic times, i.e., when the billiard particle does a large number of bounces.Consider a trajectory carrying out reb N bounces.Because the system is a classically deterministic system, each segment of trajectory between two subsequent bounces depends on and is completely determined by the preceeding segment of trajectory.However, if one computes the correlator between trajectory segment n and trajectory segment n k + , such correlation tends to zero, when k becomes large (see below equations ( 16) and ( 17), and reference [52]).This effect is due to the dynamics of classical chaos (positive Lyapunov exponent).Thus two segments of the same trajectory, which are sufficiently distant in terms of travel time or number of intermediate bounces, will become statistically independent.This makes the system similar to a random walk system, where any two segments of a trajectory are statistically independent.To see this, consider trajectories of the chaotic billiard composed of N segments ( 1 N − bounces), and choose N k  .Then we can consider all segments j i, with i j k − > to be statistically independent.Consequently, for N k  almost all segments are statistically independent from each other and hence the system will statistically behave similar to a random walk system.This means in the regime of a large number of rebounds that the system looses its memory of the detailed underlying laws of physics.It becomes universal.
Such universality manifests itself in the statistical observable of distribution of length of trajecories ( ) P Λ which for all chaotic billiards studied here shows a near-Gaussian character (like the random walk system).
2) What is the physical significance of such universality?If one considers classically chaotic billiards in the regime of macroscopic times, where Λ denotes the matrix of length of trajectories, and if one considers as statistical observable the distribution of length of trajectories ( ) P Λ , then ( ) P Λ will asymptotically approach a Gaussian and Λ will asymptotically become a random matrix from a Gaussian Orthogonal Ensemble (GOE).Hence the level spacing distribution of the length matrix Λ will be described by GOE random matrix statistics.However, one should be careful.The matrix of length Λ has not lost all information on the physics of the underlying system.For instance, the variance of length 2 ∆Λ turns out to behave differently in the Sinai billiard (on a in the finite horizon regime, scaling like n, in contrast to the stadium billiard, where it scales like ( ) n .Statistical behaviour depends on the particular observable considered.Thus using the length marix Λ , which via ( ) P Λ and ( ) P s displays universal behaviour, we can compute physical observables, e.g., being related to 2 ∆Λ .Such observable yields transport properties, in particular the diffusion coefficient.Results of transport properties are presented in Section 6.
3) Are there connections between universality and physical quantities which are easily observed in real systems?In particular, is such universality related to transport properties of the classical system?The answer is yes, and we will show in Section 6 how a transport coefficient can be obtained from the length matrix Λ in the regime of universality.
Moreover, the behaviour of the chaotic billiard systems when approaching the regime of universality is characterized by laws specific for the particular billiard system.As example consider the decay laws of correlation functions.Consider a billiard trajectory and consider as observable f the segment of trajectory from n to 1 n + , while g denotes the segment from n k + to 1 n k + + .Then the correlation function between observables f and g for many billiards shows either an exponential fall-off behaviour [52] ( ) , e , or polynomial fall off behaviour ( ) where A and B are constants, 0 k ≠ and the exponent 0 α > is specific for the particular billiard.For example, the triangular Lorentz gas as well as the Bunimovich stadium billiard obey exponential falloff behaviour [52].
Here we suggest for chaotic billiard systems that the approach towards universality (i.e., increasing the number of rebounds) contains further physical information characteristic for the system.In particular, we expect that the distribution ( ) for reb N → ∞ , where the asymptotic approach is characterized by an exponent characteristic for the particular billiard.
4) Are those universality properties related to thermodynamical observables?If we consider the chaotic billiards consisting of a single particle moving in a rigid environment of scatterers, then it does not make sense to talk about thermodynamics.Thus for systems considered in this work the answer is no.However, if one considers billiards of many particles, then thermodynamics (as a function of temperature) will influence the dynamics.We shall defer the study of such effects on universality to future work.
5) Is such universality related to spectral statistics of the corresponding quantum system, i.e., what is the relation between universality in classical chaotic systems and universality in chaotic quantum systems, as defined via the Bohigas-Giannoni-Schmit conjecture [7]?This is a very interesting question, for which we do not have an answer.The finding that universality properties exist in both, classical and quantum chaos, may be a hint that such a relation actually exists.On the other hand, for universality in quantum chaos in the semiclassical regime, periodic orbits play a crucial role (via Gutzwillers trace formula).In contrast, in our study of chaotic classical billiards universality is captured in the length matrix of bouncing, zig-zag, non-periodic trajectories.In order to address such question further, we suggest to consider the path integral relation (in Euclidean, i.e. imaginary time) which relates the quantum Hamltonian H to the classical Euclidean action Eucl S .For billiards, recall that the length of a trajectory is essentially equivalent to its action.Thus spectral properties of the quantum Hamiltonian may indeed be related to the spectrum of the classical action.

Transport Properties from Length of Trajectories
Above we have shown for the Bunimovich stadium billiard and for the Sinai billiard that they display universality properties via the statistical behaviour of the matrix of length of trajectories.Here we will show that such universal behaviour is related to relevant physical quantities.In particular we will extract transport properties from the length matrix (note the analogy to computation of transport properties in quantum in semi-classical regime via Gutzwiller formula).As examples we consider the stadium billiard.
In systems for which the CLT is verified, the diffusive character is manifested by a linear relation between the time of travel and the variance of position.The diffusion coefficient d (in 2 dimensions) is given by the Einstein relation ( ) where τ denotes the time of travel and ( ) var X  denotes the variance of the position.A similar linear relation between the time of travel and the variance of the length of trajecories can be demonstrated for the periodic triangular Lorentz gas with finite horizon, verifying the CLT (results not shown here): where travel time τ is related to mean trajectory length of n collisions n Λ and velocity u via n u τ = Λ .Based on this approach, we define a diffusion coefficient with respect to the variable Λ : for the Bunmovich stadium billiard, which is a chaotic system with concave repeller/scatterer.We have chosen to analyze such system, and compute its transport properties because it reveals a very interesting non-standard diffusion behaviour.In a numerical modelling study of the stadium billiard, Borgognoni et al. [54] have studied its diffusion behaviour by testing if a linear relation holds between the time of travel and the variance of angular momentum.They found that a linear relation fits quite well their data and computed the diffusion coefficients as ratio of variance of angular momentum and time of travel (or number of bounces).However, in 2006, Bálint and Gouëzel [53] proved rigorously for the stadium billiard that a "non-standard" limit theorem holds, where the width does not follow the usual n as in Equation ( 22), but rather follows an anomal log n n law.As consequence, the system does not converge to Brownian motion in the limit of many bounces, and the standard diffusion coefficient does not exist.
We carried out numerical simulations using trajectories of the length matrix Λ and calculated the total time of travel τ and the variance of length ( ) var Λ .According to reference [53], one expects that the variance should scale like the travel time with a logarithmic correction ( ) ( ) We define the diffusion coefficient by We have done statistical tests of = with an error of 2.8 (max) and  One observes that the statistical error of the fit including the logarithmic correction is smaller than that of the linear fit.An estimate of the anomalous diffusion constant from the logarithmic fit is shown in Figure 8(c).This means 1) numerical consistance with the mathematical result obtained by Bálint and Gouëzel [53], 2) the Bunimovich stadium billiard, although considered a non-generic billiard in quantum chaos due to the "bouncing ball states" of wave solutions, is a very interesting system in classical chaos in the asymptotic regime of many bounces, which displays anomalous diffusion behaviour, 3) the length matrix Λ which generates universal behaviour in ( ) P Λ and ( ) P s , contains infor- mation allowing to distinguish normal diffusion from anomalous diffusion.

Summary
This paper is about classical chaos occuring widely in nature, for example in astrophysics, meteorology and dynamics of the atmosphere, fluid and ocean dynamics, climate change, chemical reactions, biology, physiology, neuroscience, or medicine.We have suggested to extend random matrix theory, used in chaotic quantum systems, to classically chaotic and integrable systems.We have studied fully chaotic as well as integrable billiards and used a statistical description based on the length of trajectories to discriminate chaotic versus integrable behaviour. Results: 1) In chaotic billiards (stadium and Sinai billiards), the NNS distribution ( ) P s as well as ( ) 3 L ∆ obtained from the length matrix show GOE behaviour, i.e., they are universal.The implication is that RMT not only models spectral statistical fluctuation properties in chaotic quantum systems, but goes beyond and applies as well to chaotic classical systems.
2) The distribution of length matrix elements ( ) P Λ itself is universal to leading order.The difference between the shape of ( ) P Λ and a Gaussian distribution gives a quantitative measure of how much the motion of the billiard ball in a chaotic billiard differs from a random walk.
3) For the integrable rectangular billiard we find a correlation coefficient 1 C ≈ , in spectral rigidity ( )

L ∆
and in the NNS distribution ( ) P s strong evidence for rigid behaviour and strong correlation between neigh- bour eigenvalues.Such behaviour can be understood from the observation that integrability introduces strong correlations in length matrix elements.This proliferates to the spectra.4) In contrast, for integrable quantum systems the NNS distribution generally shows Poissonian behaviour with correlation coefficient 1 C = (there are a few exceptions).Thus from the point of view of level spacing distributions, a marked difference shows up between the classical and the quantum world.In our opinion, such difference is due to the group properties/symmetries of classical trajectories.
Future directions: 1) We plan to do numerical studies to investigate if universality also holds in chaotic potential systems.
2) We hope that our findings may contribute to obtain a unified description of both, quantum and classical chaos, and help understanding why quantum chaos is typically weaker than classical chaos, e.g., via an effective quantum action [55] [56].
3) The global statistical approach to classical chaos proposed here may help to give insight into the problem of ergodicity breaking in Hamiltonian systems (e.g., dense packing of discs in the Lorentz gas model [27]).
) this indicates high correlation at short and long range and very rigid behaviour.

Figure 1 .Figure 2 .
Figure 1.2D Rectangular billiard.Horizontal length a, vertical length b.Trajectories go from boundary points i x to j x located in the upper right quarter of the billiard wall.

Figure 2 (
b) shows the spectral rigidity ( ) 3 L ∆ of the spectrum of the length matrix.It rapidly saturates to a value Dirac comb.We have also computed the correlation coefficient to obtain 0.996 C =

Figure 3 .
Figure 3. (a) Geometry of 2D Sinai billiard (not to scale).Parameters: radius of sphere 1 R = , half-length of box 5 d = .Trajectories go from boundary points i x to j x located in a segment of one eighth of the wall of the billiard; (b) Relative error as function of number of rebounds for 2D Sinai billiard.
n = trajectories, are shown in Figure 4(b) & Figure 4(c).The results show a Wigner distribution for ( ) statistics (which seems consistent with the numerical results shown in Figure 4(b) & Figure 4(c).We pointed out above that the random walk model gives a Gaussian distribution for the histogram of the length matrix elements ( )

Figure 5 .
Figure 5. 2D stadium billiard.Trajectories go from boundary points i x to jx located in a segment of one quarter of the billiard wall.

Figure 6 (
to all possible trajectory indices ν .For each trajectory index ν we computed a level spacing distribution.Afterwards we superimposed the spectra coming from a smaller subset of trajectory indices traj 1, , n ν =  (like we did for the Sinai billiard).The resulting level spacing distribution ( ) P s of length eigenvalues and the spectral rigidity b) & Figure 6(c).The results are consistent with a Wigner distribution, i.e.GOE behaviour.The histogram of the length matrix elements ( ) P Λ shown in Figure 6(a) looks close to a Gaussian.

Figure 6 .
Figure 6.2D stadium billiard.Average over trajectories.Geometry parameters 2.2 a b = .Boundary points on quarter of wall 40 N = .(a) Distribution of length matrix elements (after subtraction of mean Λ ).reb 30 N = ; (b) Level spacing distribution ( ) P s from length matrix Λ .

2 ∆ΛFigure 8 .
comparing a fit linear in n with a fit linear in Making fits of the form ( )

Figure 8 .
Figure 8. Anomalous diffusion in stadium billiard.(a) Variance of Λ linear fit in n; (b) Variance of Λ fit including logarithmic term ( ) log n n ; (c) Diffusion constant from including logarithmic term.

Figure 8 ( 1 0
Figure 8(a), compared to the fit with logarithmic correction giving 2 characteristic feature of chaos is encoded in the number of classical trajectories.For the Sinai billiard we found that the number traj N of trajectories averaged over all pairs of boundary points increases as function of the number of rebounds reb N like an exponential,