Symmetries and spectral statistics in chaotic conformal field theories II: Maass cusp forms and arithmetic chaos

We continue the study of random matrix universality in two-dimensional conformal field theories. This is facilitated by expanding the spectral form factor in a basis of modular invariant eigenfunctions of the Laplacian on the fundamental domain. The focus of this paper is on the discrete part of the spectrum, which consists of the Maass cusp forms. Both their eigenvalues and Fourier coefficients are sporadic discrete numbers with interesting statistical properties and relations to analytic number theory; this is referred to as `arithmetic chaos'. We show that the near-extremal spectral form factor at late times is only sensitive to a statistical average over these erratic features. Nevertheless, complete information about their statistical distributions is encoded in the spectral form factor if all its spin sectors exhibit universal random matrix eigenvalue repulsion (a `linear ramp'). We `bootstrap' the spectral correlations between the cusp form basis functions that correspond to a universal linear ramp and show that they are unique up to theory-dependent subleading corrections. The statistical treatment of cusp forms provides a natural avenue to fix the subleading corrections in a minimal way, which we observe leads to the same correlations as those described by the [torus]$\times$[interval] wormhole amplitude in AdS${}_3$ gravity.


Introduction
The importance of chaos for conformal field theories and the AdS/CFT correspondence has become increasingly apparent over the years.Quantum chaos is often formulated as a statement about the statistics of the spectrum of energy eigenvalues.Energy levels that are sufficiently close together are expected to show the same statistics as those of an appropriate random matrix ensemble, namely eigenvalue repulsion; the probability of energy levels being nearby decreases as they get closer.This leads to a linear ramp in the spectral form factor, which is the averaged product of partition functions, at late times.Holographic conformal field theories possess a dense spectrum for large enough energies for any spin, and thus could possibly display random matrix universality.
In theories with symmetries, only the parts of the spectrum that are independent of the symmetries can display random matrix universality.In particular, the spectrum of conformal field theories in two dimensions is subject to translation invariance, Virasoro symmetry, and modular invariance.We can remove the consequences of translation invariance and Virasoro symmetry by focusing on conformal primary operators in fixed spin sectors.This leaves the question of modular invariance, which relates primaries of different energy and spin.
In [1], we began investigating the relationship between quantum chaos in two-dimensional CFTs and modular invariance.Motivated by the pure gravity wormhole amplitude found by Cotler and Jensen [2,3] (see also [4]), we argued that random matrix statistics for the "nearextremal" part of the dense spectrum and the corresponding late time linear ramp is an independent feature of each spin sector separately.This is a non-trivial statement because the exact spectrum is fully determined by only the spectrum of spin zero primaries and those of a single non-zero spin.The focus of this analysis was on CFTs where the ramp is encoded solely in the continuous part of the basis of modular invariant functions.There exists a discrete part as well, the Maass cusp forms, that can also encode the ramp.
The cusp forms are interesting objects in their own right: • Cusp forms arise as bound states for a particle moving in the fundamental domain of SL(2, Z).This is a classically chaotic system, however due to its highly symmetric structure it does not obey random matrix statistics.
• Instead, their spectrum of eigenvalues R ± n is bounded from below and is Poisson distributed, i.e corresponds to independent draws from a known distribution. 1 Their Fourier coefficients a (n,±) m for prime spin m are bounded by ±2 and are also Poisson distributed independently drawn from known distributions.As a consequence of Hecke relations the Fourier coefficients for non-prime (composite) spins are polynomials of those with prime spins.This implies that the distributions of non-prime spin Fourier coefficients are also determined by the distributions for prime spins.
We sometimes refer to the collection of eigenvalues and Fourier coefficients as cusp form data. Together their statistical properties are sometimes referred to as arithmetic chaos [5][6][7][8].
A connection between arithmetic chaos and quantum chaos in CFTs is an intriguing possibility (first hinted at in [9]), particularly since quantum chaos in the wormhole amplitude is encoded solely in the cusp forms for all non-zero spins [10].On the one hand, the fact that objects linked to chaos appear in the spectral decomposition of CFTs suggests that the two (very different) types of chaos may be linked in some way.On the other hand, arithmetic chaos is a property of the general modular invariant basis functions, not of the actual CFT A depiction of the statistical approximation to cusp form data.The "spectral staircase" and the erratically distributed Fourier coefficients are replaced with their average values.We will justify this 'statistical' coarse-graining in the time regime that is relevant for random matrix universality.It should be distinguished from 'microcanonical' coarse-graining, which is always required to discuss correlations in the CFT spectrum.
spectrum, so it is also present in integrable CFTs.We intend to clarify the relation in this work.

Summary of results
In this paper, we extend the techniques developed for the continuous SL(2, Z) spectrum in [1] to the cusp forms and identify the relationship between arithmetic and quantum chaos.We show that taking the near-extremal, late time limit automatically implements a statistical averaging over the cusp form data, in particular over their sporadic eigenvalues and erratic Fourier coefficients.On the one hand, this means that quantum chaos in 2d CFT depends only on statistical features of arithmetic chaos, not the detailed structure of the cusp form data. On the other hand, it is remarkable that full information about (i.e., all statistical moments of) the distributions of the cusp form data is encoded in the spectral form factor, assuming it exhibits a linear ramp in all spin sectors.We find the universal form these statistically averaged correlations in the cusp form sector must take to produce a ramp.As with the continuous sector, modular invariance does not spoil the independence of random matrix universality in each separate spin sector, since the statistical averaging proceeds differently in each spin sector; a linear ramp must be imposed as a separate assumption for each spin sector in order to fully and consistently determine the correlations in cusp form expansion coefficients.
By demanding that there be a linear ramp in every spin sector, we are able to "bootstrap" the exact cusp form correlations whose statistical averaging produces random matrix statistics independently in every spin sector.These correlations depend on all moments of the distributions of the Fourier coefficients for prime spin, and are related to well studied numbertheoretic objects.In fact, these correlations are essentially universal and unique under some mild assumptions.The gravitational wormhole amplitude exhibits the same universality, while at the same time having the minimal subleading corrections (in the late time limit) to make it consistent with modular invariance [10]; under some related minimality assumptions our construction reproduces it exactly.
Our presentation is somewhat pedagogical.For the result on cusp form correlations encoding a linear ramp in all spin sectors, see (3.20) and (3.23).We derive this result by investigating statistical features of the sum over 'arithmetically chaotic' cusp forms.We discuss the connection with Euclidean wormholes in section 4.

Outline
The paper is organized as follows.In section 2 we review the setup of [1], introducing the fluctuating part of the partition function and decomposing it in the complete basis of modular invariant functions.In section 3 we analyze how random matrix statistics appears in the cusp forms, and demonstrate its reliance on only arithmetic chaos.We derive an expression whose statistical average produces a ramp in each spin sector and show that it is unique under mild assumptions.In section 4 we then show that this expression exactly matches a calculation in AdS 3 gravity.In the discussion, we put forth some preliminary results on the spectral decomposition of the self-correlations in the spectrum , and how it differs from eigenvalue repulsion.
Conventions are collected in appendix A. We review statistical features of cusp forms in appendix B and discuss some important mathematical properties in appendix C. Appendix D concerns the imprint of linear ramps in a given spin sector onto other spin sectors.

Spectral decomposition of the ramp
We start by reviewing the SL(2, Z) spectral decomposition of the linear ramp and recollecting the work of [1].In the next section we extend this analysis to the cusp forms, and in particularly use statistical properties thereof (dubbed "arithmetic chaos") to simplify the analysis.

SL(2, Z) spectral theory
Beginning from the full CFT partition function on a torus with modular parameter τ = x+iy, Z(x, y), in [1,9] the authors introduce a fluctuating partition function Z P (x, y) by a series of steps intended to account for the symmetries of the problem: • First, the partition function is divided by that of a single non-compact boson, Z 0 = 1/(y 1/2 |η(x + iy)| 2 ), to remove all Virasoro descendants in a modular invariant fashion.
• Then one realizes that the 'censored' part of the spectrum (i.e., states with h or h ≤ c−1 24 , equivalently E ≤ 2π m − 1  12 ≡ E m ) is not typically chaotic.Therefore, one subtracts off this part of the spectrum.In addition, one also removes the part of the dense spectrum (h, h > c− 1  24 , equivalently E > E m ) that is determined from the censored spectrum by symmetries (such as modular S-transformations); together, these two parts are called the modular completion of the censored spectrum, Z C (x, y).
The resulting "fluctuating" partition function Z P (x, y) is the object that can display quantum chaos.Finally, we write this object in terms of a decomposition into sectors with definite spin: This is not quite an ordinary partition function: the density of states it describes corresponds to fluctuations around the average density of states.
To understand why this is, we have to explain the process of modular completion.There are different ways of performing the modular completion of the censored spectrum ρ C (E), which are all modular invariant and do not introduce new censored states.We focus on the kind introduced in [11], where the modular completion of each censored state results in a continuous density of states in the dense part of the spectrum 2 .For example, the modular completion of the vacuum gives a continuous density of states for the dense spectrum that includes the well-known Cardy formula for the leading average density of states at high energies.Other censored states give additional contributions that are subleading at high energy and large central charge.
In effect, the modular completion defines our coarse-graining procedure: where 'D' refers to the dense part of the full spectrum. 3With this prescription, subtracting the modular completion of the censored spectrum amounts to eliminating the latter while also removing the average density of states from the heavy spectrum.Explicitly, is the fluctuating density of states corresponding to the fluctuating partition function, where we include the normalization factors from the non-compact boson.Thus, (2.1) is the modular-invariant object that can display quantum chaos.It is 'fluctuating' in the sense that ρ P (E) has a vanishing microcanonical average (in particular it has both positive and negative contributions).
The linear ramp: We are interested in the universal correlations that this fluctuating spectrum exhibits due to quantum chaos.In particular, a quantum chaotic CFT will have a universal asymptotic contribution to the variance of Z P , describing eigenvalue repulsion.This is often called the 'linear ramp' and corresponds to analytically continuing y 1 → β + iT and y 2 → β − iT and taking the large T limit of the spectral form factor.A linear ramp is captured in y i variables by the following limiting behavior: where '. ..' denotes subleading terms in the large y i limit.This should be read as a statement about each spin sector separately.We choose here (and henceforth) a normalization for the ramp corresponding to the GOE universality class, which matches the discussion in [10,13].The normalization would be different for other universality classes, in particular it would differ by a factor 1 2 for GUE as in [3] and our previous work [1]. 4 .

Spectral decomposition:
It is useful to expand the fluctuating partition functions in a complete basis of normalizable modular invariant functions on the fundamental domain F [9] (see also [10,[14][15][16]).Such eigenfunctions consist of a continuous spectrum of Eisenstein series E s (y) with s ∈ 1 2 + iR, and a discrete spectrum of Maass cusp forms ν n,± (y): in addition to a constant function, ∆ F ν 0 (τ ) = 0. Note that there are both even (+) and odd (−) cusp forms, so there are two sets of eigenvalues {R ± n }.These are randomly distributed, which we will quantify later.After expanding the fluctuating part of the partition function in this modular invariant basis, the expansion coefficients are then unconstrained by modular invariance and their statistical properties are a good diagnostic of chaos.To write this, we refine the decomposition (2.1): where the spectrum consists of the following pieces:5 spin 0, continuous: spin > 0, continuous: with vol(F) = π 3 and the norm of cusp forms refers to the Petersson norm.The modular invariant expansion coefficients are {z 0 , z n,± , z 1 2 +iα }, and we are interested in their variance and how it encodes the linear ramp.Using the explicit basis functions for m > 0, where Λ(s ).We also trivially obtain a decomposition into bases of Bessel functions (both with continuous order, K iα , and sporadic discrete order, K iR ± n ) by defining the spin-dependent spectral overlap coefficients: (2.10) The fact that {z 0 , z n,± , z 1 2 +iα } are independent of spin leads to spectral determinacy [9]: full knowledge of Z m P, cont./disc.,,±(y) for only m = 0 and a single non-zero spin determines the partition function for every other spin. 7That the coefficients must be independent of spin will prove to be important.

Linear ramp from correlations in spectral overlap coefficients
We wish to discuss how the ramp (2.4) translates into specific universal correlations between the coefficients of the spectral decomposition.This discussion should a priori be had for each spin sector individually.
we give the translation: (2.8) 7 For partition functions that arise as Poincaré series, further constraints follow [10].However, we do not assume this here.
Ramp for spin 0: For spins m 1 = m 2 = 0, the ramp is encoded in Z 0 P (y 1 ) Z 0 P (y 2 ) , and in particular it is fully determined by correlations in the overlap coefficients with Eisenstein series [1]:8 (2.11) where terms subleading in the large y i limit (denoted as 'const.')are required to obtain a modular invariant expression.While the correlations in z 1 2 +iα will in general contain more information, the above should be understood as the universal contribution that is due to a ramp in the spin 0 sector. 9amp for non-zero spins: For non-zero spins, decomposing the ramp using (2.6) and noting that the leading term shown in (2.4) is even in spin, it is clear that the ramp could apriori be encoded in cross-correlations between any of the terms in (2.6) (subject to producing the correct parity).We will now discuss the form of the correlations that can encode a ramp in a single spin sector.However we will later stress the limitations of this approach when asking for linear ramps in more than one spin sector simultaneously, as these are not independent of each other and additional consistency conditions must be imposed.
Using the spectral decomposition of the partition function, we can write its even and odd parts for spin m ≥ 1 in a basis of Bessel functions, where each individual term is still modular invariant by construction: (2.12) Let us now briefly review how the ramp could be encoded the Eisenstein series correlations (see [1]).For the continuous part of the spectral decomposition, we can use the orthogonality of Bessel functions to invert the α-integral in (2.12): This allows us to translate the universal expression for RMT eigenvalue repulsion, (2.4), into an expression for the correlations of z m (α) coefficients by performing two correlated y-integrals of the form (2.13): This shows how a ramp for specific spin m can be encoded in the coefficients of the Eisenstein series, and the required correlations are again diagonal in the spectral parameter. 10It is straightforward to verify this result explicitly by transforming (2.14) back to y-variables, which reproduces (2.4).We can equivalently write (2.14) as: (2.15) The fact that this relation depends explicitly on spin might be understood as follows: the existence of a ramp in each spin sector gives partial information about the correlation of the modular invariant coefficients in different regimes, roughly organized by scale of oscillation as function of α.The different regimes are spin-dependent, so (2.15) is to be understood as being valid only in the regime informed by the asymptotic form of the spin-m partition function. 11Nevertheless, after discussing a similar situation for the cusp forms, we will return to this issue in section 5 and point out some subtleties.
3 Ramp from cusp forms -the statistical approximation In this section we extend the above analysis to the cusp forms.The main tool we use is a continuous approximation to the sum over the cusp forms, which utilizes statistical information, known as "arithmetic chaos".We introduce the approximation in the context of a simple ansatz that yields ramps in a single spin sector from cusp form sums, and show how this is reproduced to good accuracy by the statistical approximation.This sets the stage for an improved ansatz, discussed in the next section, where we "assemble" ramps for all spin sectors simultaneously.We would like to explore the type of correlations in the overlap coefficients, ⟨z n 1 ,± z n 2 ,± ⟩ ramp , which yield a linear ramp through a sum over cusp forms: 10 Since z m (α) is even in α by definition, we refer to the presence of the symmetrized sum of delta-functions in (2.14) as diagonal. 11To organize the information conveniently and discuss the relationship between all the statements implied by the ramp in different spin sectors, ref. [1]  where we normalized such that a (n,±) m=1 = 1.The Fourier coefficients are distributed according to a Wigner semi-circle for prime spins m → ∞.Studying the nearest neighbor spacings reveals that both the eigenvalues and the Fourier coefficients (for fixed spin) are Poisson distributed -a fact we shall refer to as arithmetic chaos; see Appendix B for details and plots.In a sense, arithmetic chaos is more akin to an integrable rather than a chaotic structure.One of our goals is to elucidate the relationship between this randomness in the structure of the Maass cusp form expansion and the genuine quantum chaos described by the linear ramp in the spectral form factor.To reproduce the ramp from a sum over cusp forms, we will have to address this interplay.
A central ingredient in our analysis is a certain continuum approximation to the discrete sum over cusp forms; relatedly, we will argue that all cusp form data can be replaced with its statistical average, which we will explain in turn.Before giving details, let us summarize the steps we will follow: 1.To find ⟨z n 1 ,± z n 2 ,± ⟩ such that the sum (3.1) yields a ramp, we first note that the linearly increasing density of eigenvalues R ± n allows us to approximate the sum by an integral over a continuous eigenvalue density.We will argue that this approximation becomes arbitrarily good for large y i .Equivalently, we can think of the large y i limit as implementing a statistical averaging over eigenvalues.
2. While less obvious, we will show that the large y i limit also acts as a statistical averaging over the Fourier coefficients a (n,±) m . Since they appear squared in the spectral form factor, the cusp form sum is effectively only sensitive to their statistical variance.Thanks to certain Hecke relations, the information contained in the variances of Fourier coefficients for all spins m is equivalent to the information contained in the full distribution of those with prime m.

Using these statistical properties, we illustrate what kind of correlations ⟨z
can yield a ramp in a given spin sector.We then show how to get a ramp in every spin sector in a very constrained way.The correlations thus obtained come with a certain amount of freedom.We show that fixing this freedom in the simplest possible way leads to a result that matches the pure gravity wormhole amplitude [3].
Throughout this section we make extensive use of a database of 5832 even and 6282 odd Maass cusp forms (corresponding to eigenvalues R ± n < 400), computed in [17] (see also [18] for a subset).We also assume the non-degeneracy of cusp forms, which is a widely believed but unproven conjecture.

Statistical treatment of the sum over eigenvalues R ± n
A ramp can be encoded in the coefficients of Maass cusp forms; to extract this, we need to invert the discrete part of (2.12).This requires an appropriate regularization of the integrals over Bessel functions to ensure their orthogonality in the discrete solution space.We avoid this technical point for the moment by working with an approximate continuous representation.This will allow us to derive the solution.We will see that this representation utilizes many of the statistical properties of the cusp form, thus connecting arithmetic chaos to the expansion of the ramp in the cusp forms.
To start we define the density of cusp forms by µ ± (R), defined by where z m R,± is a smooth function of R such that z m R ± n ,± ≡ z m n,± .We will justify by construction that this is consistent with a ramp.
The asymptotic density of cusp forms can be approximated by a continuous function, using the 'Weyl law' (see for example [19,20]): (3.5) The lower cutoff r ± > 0 in (3.4) is chosen appropriately such as to avoid over-counting of the constant cusp form.We review this approximation and various other statistical properties of the Maass cusp forms in appendix B, see in particular figure 5.For the purpose of our analysis, the smooth approximation to the density of eigenvalues R ± n sometimes allows us to replace sums by integrals: n is performed up to some cutoff such that convergence is achieved within the displayable accuracy.The plots show that the sum converges to the ramp linear y/(2π) up to an mdependent constant that is subleading as y → ∞.
which one might expect to hold for sufficiently smooth functions f .Clearly the approximation is better for functions f with support at larger values of R, since the eigenvalue density increases linearly with R; thus, more precisely, for any ε > 0 and sufficiently smooth functions f , there is a sufficiently large n 0 such that Working with this continuous approximation, a calculation identical to (2.14) gives: One can immediately see that this approximate continuous expression translates into the following correlations for the discrete coefficients: or, equivalently, the cusp form sum (3.1) encoding a linear ramp should be of the form We quantify how the sum over cusp forms indexed by n depends on the terms with small n (and hence small R ± n ).We plot the ramp in the spectral form factor computed using only values of n for which R ± n > R min and normalize it by the complete result.Asymptotically as y → ∞ this ratio converges to 1, no matter how many low-lying values of R ± n we exclude.We show the cases of even (left) and odd (right) cusp forms separately (they are almost indistinguishable).Solid lines correspond to spin m = 1, dashed lines to m = 2.This is approximate in the following sense.To evaluate the sum we can proceed in two ways: (i) analytically, we can approximate the sum by an integral as in (3.6), which in turn recovers the exact ramp in every spin sector (by construction).The approximation is then due to replacing the sum by an integral.This approximation becomes increasingly good for y i → ∞ because the support of the Bessel functions shifts to larger values of R ± n where these are more dense.(ii) Numerically, we can confirm directly that the discrete sum (3.10) does also reproduce the ramp up to an error (a subleading constant shift) that goes to zero as y i → ∞. 12 Figure 2 illustrates the result (both for even and odd parity cusp forms).We see that the numerical evaluation of the Maass cusp form sum asymptotes to the expected linear ramp for large values of y i (we only show the case y 1 = y 2 ≡ y, but other cross sections of the (y 1 , y 2 ) plane were checked similarly).

Dominance of large R ±
n as y i → ∞: We wish to comment more on the convergence of the sum over cusp forms indexed by n, (3.10), and the continuous approximation.As functions of y, the Bessel functions K iR ± n (2πmy) have strong oscillations for 0 < 2πmy ≲ R ± n , with amplitude of order 2π/R ± n e −πR ± n /2 , and subsequently decay exponentially like 1/(4my) e −2πmy , independent of R ± n .First, the exponential decay implies that the sum over n converges and can thus be truncated in numerical evaluation.More non-trivially, the sum is dominated by terms with increasingly large values of R ± n .We verify this numerically in figure 3: we compare the ramp computed using all relevant terms in the sum with the partial result obtained by dropping all terms with 0 < R ± n < R min .We observe that the ratio of these two quantities approaches 1 as y → ∞, for any choice of R min .Another way to present the same result would be as follows: any partial sum over only low-lying R ± n converges to a finite constant (times the usual e −2πm(y 1 +y 2 ) ) as y i → ∞, as the Bessel functions become independent of R ± n .This is therefore subleading to the ramp: lim for any n max .Effectively the sum over n is dominated by a window R min ≲ R ± n ≲ R max where both R min and R max increase with y.Thus using the continuous approximation, i.e., treating the eigenvalues R ± n statistically for large y, is valid. 13o summarize, the sum (3.10) can be (roughly) split into three pieces, which qualitatively contribute as follows to the spectral form factor: exponentially small (3.12)

Statistical treatment of the Fourier coefficients a (n,±) m
Let us return to the sum over cusp forms, (3.1).We wish to address the following question: what form of correlations ⟨z n 1 ,± z n 2 ,± ⟩ yields the ramp (3.10)?Naively, it seems that we have already answered this question in (3.9).However, that expression, taken literally, would via (2.10) give a different, spin-dependent form of ⟨z n 1 ,± z n 2 ,± ⟩ for every spin, which clearly cannot be correct.So how is (3.10) consistent with spin-independent correlations ⟨z n 1 ,± z n 2 ,± ⟩?
To resolve this conundrum, we take a detour to discuss properties of the Fourier coefficients a (n,±) m of the cusp forms.
The Fourier coefficients are erratic, see (3.2).What does this mean for the validity of our continuous approximation to the eigenvalues?We have already seen that the sum over n can be split into three pieces, see (3.12).Both R min and R max increase indefinitely as y i → ∞, so the density of R ± n , which are of interest to regime (2), increases as well.Summing over an increasingly dense set of R ± n acts as a statistical coarse-graining over the remaining ndependent terms in the sum.In particular, the product of the Fourier coefficients appearing in the sum and the correlations ⟨z n 1 ,± z n 2 ,± ⟩ get averaged over.We therefore expect to be able to replace the discrete erratic Fourier coefficients by their statistical distribution.

Distribution of Fourier coefficients:
The statistical distribution of the Fourier coefficients is a well-known topic of mathematical research, and we review it in some detail in appendix B. Let us only point out the most crucial aspects.First, the asymptotic distribution of a (n,±) m for fixed prime spins m ≡ p is well known [21]: For large prime spins, this approaches a Wigner semicircle (2π) −1 √ 4 − x 2 .Another notable feature is that the distribution suggests that |a (n,±) p | < 2 for all n, a property known as the Ramanujan-Petersson conjecture [21].We are interested in moments of these distributions.Since the sum (3.1) features the squares of Fourier coefficients, a statistical feature of particular interest is their variance which has the following exact value for prime spins: 14 See (B.5) for higher moments.We use the notation (• • • ) to denote statistical averaging (over n).This is independent of the microcanonical averaging, denoted by ⟨• • • ⟩, which we always use to discuss correlations in the coarse-grained CFT spectrum.
For non-prime spins m, the variances N ± m are determined by the distributions for prime spins.Importantly, not only the variances of the distributions for prime spins, but also their higher moments are needed.The reason is that the Fourier coefficients themselves are determined as non-linear polynomials of those for prime spins by a certain Hecke algebra, see (B.9) for some examples.Statistical averaging over such polynomials requires knowledge 14 It is interesting to note that since the Fourier coefficients for prime spins are Poisson distributed, as shown in appendix B, the variance in (3.14) already implies delta-functions in spin and eigenvalue index, p1, p2 prime: a This suggests that arithmetic chaos is linked to the diagonal approximation in the periodic orbit picture of [10].The delta-function in the eigenvalue indices persists even for non-prime spins and is therefore tied to the effective statistical averaging implemented by the correlated cusp form sums (3.1).Note, however, that we will later average over summands involving higher moments of Fourier coefficients, which complicates the picture.
of higher moments of the prime distributions.In summary, the following three pieces of information are equivalent: where values for prime m are printed in boldface (see tables 1 and 2 for more details).
Statistical averaging in the spectral form factor: Whenever the statistical averaging over n applies to our cusp form sum over n, it means that we can replace discrete erratic expressions by their statistical average.This amounts to a significant simplification for evaluating sums such as (3.1): for large y i the exact squared Fourier coefficients (which oscillate erratically) can be replaced by their mean value, i.e., the variance of their distribution (3.13), thus 'forgetting' about the detailed sporadic values and only keeping track of statistical information.This explains how it was possible that the correlations ⟨z n 1 ,± z n 2 ,± ⟩ that would follow from (3.9) would depend on spin in such an extremely fine tuned way as to cancel all erratic Fourier coefficients a : the correlations ⟨z n 1 ,± z n 2 ,± ⟩ do not actually need to cancel the Fourier coefficients exactly, but only on average.As we will see, this is indeed possible in a spin-independent way.
Focusing on a single spin sector, the fact that the Fourier coefficients only need to cancel on average means we would expect to reproduce the linear ramp in the spin m sector from This corresponds to correlations of the form 15 15 The second approximation, ≈ 24 π 2 , is valid asymptotically for very large n, i.e., for very large yi.For all numerical results in this paper, this approximation is not good enough and is not used.We compute the Maass cusp form sum using the variance of the Fourier coefficients instead of their exact values in (3.18).For large y increasingly many Fourier coefficients contribute to the sum over n, which means that their square can be increasingly well approximated by their variance.We therefore reproduce the linear ramp asymptotically (up to a subleading constant shift), c.f. figure 2. The left (right) shows the case of even (odd) parity cusp forms.In the odd case the ramps for different spins lie almost on top of each other.Insets show larger values of y.
We can check the validity of this claim numerically by computing the sum (3.18) and comparing it with the true form of the ramp.As can be seen in figure 4, for large y the correct linear ramp is approached, again up to a constant which is subleading for y → ∞.
Evidently, (3.19) still depends on the spin m via the normalization N ± m , albeit much more weakly than had we tried to cancel the erratic Fourier coefficients in (3.18) exactly (term by term).It is therefore still not a good candidate for correlations ⟨z n,± z n,± ⟩ that yield linear ramps independent of the choice of spin.And indeed, correlations of the form (3.19) only yield a ramp with the correct slope in the spin m sector.In other spin sectors m ′ , we would need a similar form of correlations, but with a different normalization 1/N ± m ′ .We will remedy this situation in the following subsection.

Ramps in all spin sectors: number theory and uniqueness
As we have seen, (3.19)only encodes the ramp in the spin m superselection sector, but it 'contaminates' the slope of any putative ramp in other spin sectors.The basic assumption of quantum chaos, however, would be a linear ramp with the correct slope in all spin sectors.To achieve this, let us now take the statistical averaging one step further and improve the naive ansatz (3.19) such that it works on average for every spin sector, i.e., in a spin-independent way.We wish to write: with a spin-independent function f (n,±) such that The effective averaging over n will then guarantee that in the limit y i → ∞, we recover the ramp for all spins m: e −2πm(y 1 +y 2 ) + . . .
where we replaced a (n,±) m 2 f (n,±) by its average according to (3.21) in the second line and then simply applied (3.10).We denote subleading terms by '. ..'.
We will refer to f (n,±) as the arithmetic kernel associated with the cusp form ν n,± .This name is inspired by the fact that any function satisfying (3.21) must obviously depend on all Fourier coefficients for all spins in a very fine-tuned way such that it produces just the right normalization for the ramp in every spin sector.It must, in a sense, encode all the information loosely referred to as arithmetic chaos, such as Hecke relations (B.8) and the statistical distribution of Fourier coefficients (3.13).Note that the ansatz (3.20) assumes diagonality in n i .We will justify by construction that this is a consistent assumption.Note further that the condition (3.21) really only needs to hold asymptotically as a statement about the average over terms in the spectral form factor with large n. 16Deviations for small n will only affect subleading terms in the late-time spectral form factor.We fix this ambiguity in the minmal way, i.e., by imposing (3.21) as an average over all n as written.
Given all the information encoded in 'arithmetic chaos', it is remarkable that such a function exists.We will now first write down this function, then explain why it works, and then derive it, showing that it is essentially unique (under the above assumptions).The arithmetic kernel satisfying (3.21) is given by Let us first confirm that this function satisfies (3.21).We do this in three steps: 1.If m ≡ p is prime: Since the Fourier coefficients for prime spins are independently distributed, we only need to know the second and fourth moments of the distributions (3.13), which are easy to calculate.We immediately find: where every factor is individually 1 due to the following statistical facts: 2. If m = p k is a prime power: For prime power spins, we can analogously show that every factor in an expression similar to (3.24) is 1.For the first factor (p ′ = p) we need some more non-trivial facts about the Fourier coefficients, which follow from the Hecke multiplicativity rules (B.8).The required properties are (see appendix C and in particular Lemma 4): (3.27) The property (3.21) follows factor by factor.
The arithmetic kernel f (n,±) has a deep number theoretical meaning in terms of Hecke L-functions.We elaborate on these fascinating mathematical properties in appendix C. We can also derive f (n,±) from physical requirements, i.e., by merely imposing (3.21) in all spin sectors.We sketch the derivation below, delegating details to appendix C.3.

Uniqueness of the arithmetic kernel:
We will now derive the arithmetic kernel (3.23) by arguing that the requirement (3.21) fixes it uniquely (within an ansatz class).First recall that Fourier coefficients have multiplicative properties due to them being eigenvalues of Hecke operators.In particular, if the spin has a prime factor decomposition as in (3.27), since a for prime p.Since we already assumed diagonality in eigenvalues R ± n i in (3.20), odd powers of Fourier coefficients will average to zero, and we discard them in our ansatz.Such terms would not affect the construction of the universal ramp, but they would change the subleading behavior of the late time spectral form factor. Discarding odd powers in the ansatz can thus be viewed as a minimality assumption about the ansatz.It would be interesting to constrain such ambiguities further, using input from the off-diagonal sector.
The condition (3.21) yields an infinite linear system constraining the parameters c p,r in terms of moments of the distribution of Fourier coefficients.After some investigation (see appendix C.3), this system can be written as follows: Making extensive use of (i) Hecke relations and (ii) all moments of the distributions of prime Fourier coefficiens, the solution of this system for m = p k is unique: Using (3.27), the condition (3.21) for all m is then solved by While we have made simplifying assumptions in the derivation of this kernel (see the discussion after (3.21)), its uniqueness within a large class of possibilities is remarkable.We show in the next subsection that the structure of the result (3.20), (3.23) is more than just a mathematical curiosity; it has a number theoretical interpretation and its simplicity is in fact intimately tied to a calculation in AdS 3 pure gravity.

Matching universal correlations to the AdS 3 wormhole
Our 'bottom-up' construction of the spectral overlap coefficients encoding the linear ramp was based on minimal assumptions about quantum chaos in all spin sectors and consistency with the symmetries of CFTs.We also assumed a certain minimality in the ansatz for the arithemtic kernel f (n,±) , which then allowed us to fully determine it.In this section we compare this 'minimally consistent' arithmetic kernel with the wormhole amplitude found in AdS 3 pure gravity, which also exhibits such linear ramps.We find detailed agreement.Demanding universal eigenvalue repulsion (i.e., a linear ramp) in every spin sector of the CFT, and assuming that for m > 0 this property is encoded in the cusp form sector alone, we constructed the following form of spectral correlations as the simplest consistent possibility: By virtue of being spin-independent, these correlations provide a manifestly modular invariant encoding of a linear ramp in all spin sectors.(Of course, the 'bare' asymptotic ramp is not modular invariant by itself, so the subleading corrections produced by (4.1) are important.)Let us now turn to gravity.The spectral decomposition of the T 2 ×I wormhole amplitude in AdS 3 pure gravity [2,3] was given in [10], and provides an explicit example of a modular invariant spectral form factor that contains a ramp in the large y i limit. 17In our notation it corresponds to the following non-zero variances: 18 where the cusp form norms are computed with respect to the Petersson inner product (see appendix C for more details, and figure 9 for concrete values).The second line is meant to indicate that both the even and odd correlations as indicated give a ramp with correct normalization.In a CFT with parity symmetry, the even and odd ramps describe chaos in different parity superselection sectors.Now compare our result (4.1) with (4.2).The continuous part of the correlations, which encodes the spin 0 ramp, matches immediately (which is by construction).More interestingly, the discrete correlations, which we constructed by imposing quantum chaotic universality consistently across spin sectors, also match the gravity result.To see this, we need an important fact from arithmetic number theory, which is derived and explained in appendix C. The central observation is that our arithmetic kernel f (n,±) is a particular meromorphic symmetric square L-function L (n,±) ν×ν (s) evaluated at s = 1.For every single cusp form, this function provides a generalization of the Riemann zeta-function that encodes all the statistical properties and Hecke relations between different spin Fourier coefficients.The precise statement is: 3) The intermediate steps in this equation are reviewed in appendix C.This establishes equality of, on the one hand, the correlations found from demanding a 'bare' linear ramp in all 17 We thank Scott Collier for private conversation on this result. 18To compare with [10], note that To compare with [3], note that we introduced an additional factor of 2 in the wormhole amplitude to match the GOE universality class, c.f., [10,13].
spin sectors (taking into account the mechanism of statistical averaging over cusp forms and constructing a minimal spin-independent arithmetic kernel) in (4.1), and, on the other hand, the pure gravity result, (4.2).
It is interesting to note that the spin-0 ramp, encoded in the Eisenstein sector, can similarly be expressed in terms of a suitable L-function: where is the meromorphic continuation of a sum over Fourier coefficients.See appendix C.4 for more details.
To summarize, we have found that a linear ramp in all spin sectors m ≥ 1 is encoded in the following sum over cusp forms in the near-extremal limit: e −2π|m|(y 1 +y 2 ) + . . .(for all m) (4.5) The agreement of the wormhole amplitude with the 'minimal' realization of quantum chaos across spin sectors was called the MaxRMT principle in [10].It amounts to the statement that the gravity amplitude is the minimal modular completion of a spectral form factor exhibiting linear ramps.More precisely, ref. [10] shows that the wormhole amplitude is the minimal extension of the 'bare' ramp, after imposing 'diagonal' and 'Hecke' projections onto correlated eigenvalues and eigenfunctions in the spectral decomposition of the spectral form factor. 19 Our investigation similarly imposed some minimality requirements: the main assumptions were the realization of quantum chaos in all spin sectors and modular invariance; we argued that these assumptions required a spin-independent form of the arithmetic kernel and then constructed the simplest consistent kernel from an ansatz (3.20) by solving the statistical constraints.By discarding from the ansatz any terms that would be invisible to our statistical condition (3.21), we fixed it fully and recovered the wormhole amplitude.This simplicity provides a statistical perspective based on arithmetic chaos on the MaxRMT principle of [10].

Discussion
To summarize, we note again that the Euclidean wormhole amplitude (4.2) describes a universal part of the spectral correlations in any individual chaotic CFT, which dominates the late time near-extremal limit.We constructed the same object 'bottom-up' by imposing quantum chaos (in the form of a linear ramp) in every spin sector separately and consistently balancing the imprints ramps in any given spin sector have on the slope of ramps in other spin sectors.
We delineated the way in which a solution can be constructed based on statistical considerations of Maass cusp forms.A crucial role was played by the effective statistical averaging over erratic data defining the modular invariant Maass cusp forms.It is due to this averaging that, on the one hand, all statistical information about 'arithmetic chaos' is encoded in the collection of linear ramps, while, one the other hand, detailed erratic features of cusp forms are washed out and a single spin-independent form of chaotic correlations could be bootstrapped.There is some freedom in the constructed solution, which would affect subleading corrections to the spectral form factor; the match with the gravitational result was established by not making use of any of this freedom, i.e., fixing it in the minmal way.We conclude with some further comments.

Spectral determinacy
It was found in [10] that the spectral decomposition of the AdS 3 wormhole amplitude is such that the correlations in the Eisenstein series coefficients and those in the Maass cusp form coefficients are identical.We found the same result by imposing statistical universalities (quantum chaos) in all spin sectors and implementing them in a minimal way through a sum over cusp forms.This strengthens the spectral determinacy property of general two-dimensional CFTs [9], as in these examples all spin sectors exhibit identical correlations ('strong spectral determinacy' [10]).How is this consistent with one of the basic assumptions of quantum chaos, i.e., the independence of spectral universalities in each symmetry superselection sector?We take the following perspective: even though the statistical approximation required that our result (4.1) for spin m > 0 linear ramps had to be the same for all spins, it nevertheless encodes separate input from all spin sectors.This is manifest when we consider the arithmetic kernel f (n,±) : it contains all squared Fourier coefficients for all spin sectors in a highly fine-tuned way such as to ensure the correct statistical property (3.21) for all spin sectors.For example, had we only imposed the ramp in some particular spin sector, then the naive ansatz (3.19) would have been sufficient.But this would have impacted the slope of the ramp in all other spin sectors.Finding the universal kernel that yields the correct slope for all spins required us to separately assume the existence of a ramp for all spins and input the corresponding information into the construction of f (n,±) in a correlated way. 20e can summarize this as follows: imposing random matrix universality in just one given spin sector leaves a lot of freedom for the choice of the cusp form correlations ⟨z n,± z n,± ⟩, thanks to statistical coarse-graining in the late time limit.It does by no means imply a linear ramp with the correct slope for any other independent spin sector.But imposing random matrix universality in all spin sectors, leads to enough constraints to determine a universal, spin-independent form for the correlations describing the leading order linear ramp.Further, the statistical conditions we investigated naturally led to a 'minimal' solution of this problem, which agrees with the gravity result.

Deriving chaos
A first-principles, bottom-up derivation of chaos in holographic CFTs still eludes us.While we now understand the relationship between quantum chaos and modular invariance better, quantum chaos is still a basic assumption that we show is consistent with other features of the 2d CFTs.This is contrasted with the wormhole amplitude in gravity, which can be derived from first principles.Some standard properties of holographic CFTs might be sufficient for such a derivation, in particular the assumptions that yield a dense spectrum above the extremal limit (large central charge, no conserved currents, and a twist gap).A promising path towards this would be the construction of an Efetov sigma model as in [23,24], similar to how chaos is derived in the SYK model [25].It would be fascinating to see if such an approach can be adopted using recent discussions of random matrix ensembles for 2d CFT operator data and OPE coefficients, which furnish approximate solutions to the bootstrap equations [26][27][28].

The plateau
In chaotic quantum mechanics the universal form of eigenvalue correlations is expected to take the random matrix form for sufficiently close energy levels, depending on the universality class (see, e.g., [29]).For the GUE universality class, this is The first term describes the disconnected part, the third the famous sine-kernel which gives rise to the ramp in the time domain.We now wish to discuss the second term, i.e., the tautological "self-correlations", to provide comparison with the chaotic case -eigenvalue repulsion and the ramp -discussed before.In quantum chaotic systems this term gives rise to the eventual plateau for sufficiently long times (or equivalently for sufficiently close eigenvalues), but this term is even more universal as it also exists in integrable systems.Systems with Poissonian statistics are completely described by the first two terms in (5.1), up to non-universal terms at early times. 21hile the ramp appears to a natural object in the spectral decomposition, and can be described as analogous to the "diagonal approximation" [10] in a periodic orbit expansion (i.e., the correlations in spectral eigenvalues, α and R ± n , are diagonal), we will see that the plateau is perhaps less natural.This is consistent with the analogy with the semi-classical periodic orbits, for which the plateau is non-perturbative.Note also that in gravity calculations, the plateau is much more difficult to obtain than the ramp; in JT gravity it arises from an infinite sum of wormhole geometries [30,31].We will now offer a few comments on the spectral decomposition of the plateau, leaving a full analysis for future work.
Self-correlations: First, we comment on the expected height of the plateau in a quantum chaotic system, and see how this is reproduced in our language with the fluctuating partition function Z m P (y).Recall from (2.2) that the density of states for the fluctuating partition function is just the density of states for the dense spectrum minus its average, ρ P (E) = ρ D (E) − ⟨ρ D (E)⟩.This means that the second and third terms in (5.1) have to do with the correlation of the fluctuating partition function.Thus the height of the plateau we expect from considering the fluctuating partition function is just that of a standard partition function (multiplied by √ y 1 y 2 e π 6 (y 1 +y 2 ) from the definition of Z P ).Focusing on the second term in (5.1): where ⟨ρ m D (E 1 )⟩ is the average density of spin m Virasoro primaries. 22Note that the plateau coefficient is given by Z m P,D , which is not the modular invariant, fluctuating, dense partition function.It is just the standard partition function for the dense primaries of spin m; in particular it is not modular invariant.
By taking y i → ∞, we can estimate the plateau height:23 Comparing (5.4) to the ramp, we see that the ramp and plateau become equal to each other when √ y 1 y 2 ∼ T = ⟨ρ m D (E m )⟩ ≡ ∆(E m ) −1 , i.e., when at times of order the inverse mean-level spacing at threshold energy, as expected from general considerations.

Spectral decomposition:
We can now analyze how the plateau appears in the Eisenstein series; the cusp forms come with new technical issues, and we relegate their discussion to appendix D.5.For the spin 0 case, we plug (5.4) into the usual integral transform (similar to (2.13)): (5.5) The most interesting feature of this expression is that it is not diagonal in α i .24  2We use the average density from [32], given by For the spin m i > 0 case, we find similarly:25 (5.7) Again, the correlations for the plateau in any spin sector are not diagonal (there is no deltafunction imposing α 1 = ±α 2 ).From the perspective of [10], this means the plateau does not come from the diagonal approximation analogous to the semi-classical periodic orbits, as one would expect.
Similar to our discussion of the ramp, we can ask if (5.7) should be improved by imposing a plateau consistently across all spin sectors.We leave such an analysis to the future, but discuss the question of the imprint of a plateau in a given spin sector onto other spin sectors, using numerical evidence, in appendix D.4.
Such a procedure is necessary as the Fourier transform of |ξ| only makes sense as a distribution i.e. when integrated against test functions, and the resulting distribution cannot be defined without some method of regularizing the singularity at α1 ± α2 = 0.The first method is by subtracting the first two terms in the Taylor series so that the singularity becomes removable; the second is to use integration by parts and discard boundary terms, which makes the singularity integrable.We illustrate this in figure 5.
The Fourier coefficients {a (n,±) p } for fixed prime spin p and ordered by increasing corresponding eigenvalue R ± n are equidistributed according to the distributions [19,21] which approaches the Wigner semi-circle 1 2π √ 4 − x 2 as p → ∞.Equidistribution means that averages over all cusp forms for a given spin can be replaced with averages over the distribution, i.e., lim This is illustrated in figure 6.
We now investigate the nearest-neighbour level spacing, both for the eigenvalues R ± n and for the Fourier coefficients a . This provides a more numerically tractable way of analyzing the correlations than the density of states two-point function.After "unfolding" the spectrum, 26 we calculate the distribution of the difference between nearest-neighbour levels: 26 Unfolding the spectrum corresponds to replacing each member , the spectrum has constant unit density in x variables [34].An integrable spectrum is distributed according to Poissonian statistics, P P (s) = e −s , which means level attraction: P P (s) → 1 as s → 0. Chaotic spectra, on the other hand, are distributed according to the ensemble with appropriate symmetries, e.g., the Gaussian orthogonal ensemble with P GOE (s) = 1 2 πse −πs 2 /4 , which exhibits level repulsion: P GOE (s) → 0 as s → 0.
The eigenvalues of the cusp forms are known to obey Poissonian statistics [5], and we find that the Fourier coefficients for prime spin obey the same, shown in figure 7; hence, both the eigenvalues and Fourier coefficients are distributed randomly but not chaotically.Effectively, for any given spin m, the Fourier coefficients for different n are independent random variables.
We can equivalently characterize the distributions (B.2) through their moments.For the distributions of prime spin Fourier coefficients, (B.2), the odd moments are zero and the even moments are: 27p prime: a (B.5) For example: (n,±) p 6 = (p + 1)(5p 2 + 4p + 1) p 3 (B.6)These distributions for Fourier coefficients are specifically for prime spins.All Fourier coefficients for non-prime (composite) spins are fully determined in terms of these by Hecke relations because Maas cusp forms are eigenfunctions of Hecke operators: where This implies immediately: for example, if p, p ′ are prime we get the important multiplicative relation: a + δ pp ′ (see Lemma 1 for more relations of this type).The Hecke relations allow us to construct the non-prime Fourier coefficients from the prime ones.This in turn implies that the variances ('normalization factors' N ± m ) of the distributions of Fourier coefficients for non-prime spins follow from higher moments of the prime distributions.We give a few examples: We give exact analytical and numerical values for some of these moments in tables 1 and 2. The composite spins m = 4 and m = 8 are special cases of a general result, see Lemma 4. The numerical values for the second moments are within a few percent of the theoretical values.This error increases for higher moments due to the limited number of cusp forms available numerically.The numerical results for odd forms are consistently slightly better because we have more of them available.More generally, computing just the variances N ± m for all non-prime m requires knowledge of all moments of the distributions of prime coefficients.Since the distributions are bounded, their moments determine the distributions fully.In other words, knowledge of all variances N ± m is equivalent to complete knowledge of all the prime distributions (B.2).We show the variance of the Fourier coefficients for a large number of spins in figure 8.   Table 2: Numerical values of the moments of the distributions of Fourier coefficients, computed using the Fourier coefficients for cusp forms with eigenvalue R ± n < 400 (separately for even and odd parity).

C Hecke relations, cusp form norms, and L-functions
In this appendix we provide some mathematical details relating to the norms of cusp forms and their relation to objects in analytic number theory.In order to be pedagogical, we provide step-by-step proofs of some important statements.Most of the general definitions and results can be found in the literature, see, for example, [35][36][37][38].m , while the composite primes have much more complicated behaviour (which is determined by the Hecke relations and higher moments of the distribution of prime Fourier coefficients); in this range, all variances are O(1), but as m grows, the maximum possible value grows as well.

C.1 Hecke relations and L-functions
We first note some useful relations between Fourier coefficients of prime power spins.Lemma 1.Let p, p 1 , . . .p r be distinct primes.Then:  Theorem 1.The norms of the cusp forms satisfy: ν×ν (1) . (C.11) Proof: (See, e.g., references [40,41].)We compute the norm using the Rankin-Selberg trick, i.e., note that a constant expression can be computed as the residue at s = 1 with an Eisenstein series, which in turn allows for unfolding of the fundamental domain: To evaluate the residue, we note that the sum in the last line is related to a Rankin-Selberg zeta function and has an Euler product formula ( [37], Lemma 3. Taking the residue at s = 1 of both sides and using Res s=1 ζ(s) = 1 ν×ν (1) . (C.15)

C.2 Statistical averages and moments of L-functions
Having reviewed some basic facts about the cusp form norms and the distribution of their Fourier coefficients, we can now state some of the crucial properties that hold after statistical averaging over n.Let us first state the following useful Lemma 4. Let k ∈ N and p be prime.Then the average over n yields the following results: (n,±)

16)
Proof: To prove (i), we start with Lemma 1(iv) and evaluate its average using the moments of the distributions for prime Fourier coefficients (B.5): where Θ(n) = 1 if n > 0 and vanishes otherwise.The second result, (ii), can be proven in a similar fashion, using Lemma 1(ii) to simplify.To prove (iii), we use Lemma 1(ii) to calculate as follows: Finally, to prove (iv), we use (B.5) and the series representation of the hypergeometric function: Corollary: For any k ∈ N and p prime, (C.20) Proof: Follows immediately from Lemma 4(i) and (iii).
Finally, the central property needed in our analysis of the gravity amplitude concerns the interplay of the moments of distributions of Fourier coefficients and the cusp form norms: Theorem 2. Let m ≥ 1 be any integer spin.Then the statistical averaging over different cusp forms indexed by n yields: Proof: The first equality follows from Theorem 1.To prove the second equality, write the L-function in terms of its Euler product: where we applied the Corollary of Lemma 4 factor by factor after decomposing a (n,±) m into factors of Fourier coefficients of prime powers (Lemma 1(i)): In the main text we verified that the arithmetic kernel (3.23) has the required properties to produce a ramp in all spin sectors.We also outlined how to derive it, but provide more details here.The derivation essentially also shows that it is unique, up to modifications, which are invisible to our averaging condition (3.21).
Construction of f (n,±) : The goal is to find a function f (n,±) which satisfies (3.21).Since the Fourier coefficients have multiplicativity properties determined by Hecke relations, it is convenient to begin by decomposing a (n,±) m into coefficients with prime-power index: where p i are distinct primes.Let us therefore first find a function f (n,±) p which is fine tuned to Fourier coefficients with prime-power index p k such that: It is important to note that a (n,±) p k is fully determined by powers of a (n,±) p , see Lemma 1(iv).
In order for a condition such as (C.25) to hold, the function f .This is captured by an ansatz of the following form: We only need even powers of the Fourier coefficients because any odd powers will have vanishing expectation value.We also do not need any Fourier coefficients with spin other than p because these are distributed independent of a (n,±) p , so they can be absorbed into c p,r as far as the averaged (C.25) is concerned.The condition (C.25) then amounts to an infinite number of constraints on c p,r .For example: for all k.This can be written as an infinite list of equations labelled by k, which we call This can be thought of as an (infinite dimensional) triangular matrix acting on the vector of moments of Fourier coefficients.If cp,0 ̸ = 0, we can form a linear combination which cancels all terms but one: This contradicts the assumption, so we must have cp,0 = 0. Next, if cp,1 ̸ = 0 we could form a similar linear combination which is again contradictory and thus implies cp,1 = 0. Continuing this way, we must have cp,k = 0 for all k.Therefore, there does not exist any solution different from c p,k .
To summarize, we have shown that solves (C.25) and is unique as far as our ansatz is concerned. 29This function will give a ramp in all spin sectors of the form m = p k .From the multiplicative property of the Fourier 29 There are ways to modify f (n,±) p that are invisible to the statistical averaging.For example, one can add odd powers of a (n,±) p with arbitrary coefficients, as these will vanish in the evaluation of (3.21).See main text for more comments.
coefficients, (C.24), it is then clear how to construct the function that will yield a linear ramp for all spins m = p k 1 1 • • • p kr r ; indeed, we simply construct it as , where we used ζ(2) = p prime (1 − p −2 ) −1 = π 2 6 and we used the symmetric square L-function associated with the cusp form ν n,± , see Lemma 3. It is a meromorphic function in s with a potential pole at s = 1 [36] (in our case there is no pole, so we can simply evaluate at s = 1). 30 This completes our derivation of the arithmetic kernel f (n,±) .

C.4 L-function for Eisenstein series
It is natural to define the following continuous family of L-functions for the Eisenstein series E 1 2 +iα (x, y) in terms of their Fourier coefficients (c.f., (2.9)): 31 Proof: We expand the zeta-functions formally in the domain where they converge: Note in particular: which is similar to the symmetric square L-function for cusp forms evaluated at s = 1.See [22] for more details on this connection. 30Note that there would be a pole at s = 1 if there was some prime Fourier coefficient with a (n,±) p = ±2.It is unproven but widely believed to be true that such a Fourier coefficient does not exist (Ramanijan-Petersson conjecture) [21]. 31The factor 1 2 is unconventional, but will make the following discussion more convenient.

D Effects of chaos across different spin sectors
In this appendix we show (and review) that the existence of a ramp (or plateau) in a given spin sector is generically not sufficient to conclude the existence of a ramp (or plateau) in another spin sector.We previously showed this for the Eisenstein series in [1], and only briefly review those results here.We mainly focus here on extending this result to the Maass cusp form spectrum.
D.1 Signatures of a spin m = 0 ramp at spin (m 1 , m 2 ) Recall the simple form of correlations ⟨z 1 2 +iα 1 z 1 2 +iα 2 ⟩ that correspond to a linear ramp in the m = 0 sector, (2.11).Since these correlations also enter into the spectral form factor for all m > 0, we can ask about their imprint onto higher spin sectors.We previously found in [1] (numerically) for the contribution to the spin (m 1 , m 2 ) spectral form factor due to the existence of a ramp at spin 0:  where "⊃ 0 " means that we only consider the contribution to the left hand side that is implied by the existence of a ramp at spin m = 0.The first few spin-dependent prefactors in this expression are where the notation "⊃ m " means that we consider the contribution to the spectral form factor that is implied by the existence of a ramp in the spin m sector.For m = 1, this expression is particularly simple.Its numerical evaluation gives So, unlike for spin 0, the higher spin ramps do imprint onto the slope of ramps in other spin sectors.This is analogous to the situation with the 'naive' ansatz for the spin m ramp in the cusp form case discussed in the main text, see (3.19).It would be interesting to analyze if the naive ansatz for ramps in the Eisenstein sector can be improved, or if the above analysis hints at a deeper inconsistency with ramps for m > 0 being encoded in Eisenstein series at all.
D.3 Signatures of a spin m ramp at spin (m 1 , m 2 ): Maass cusp forms The numerical analysis of the previous subsection can be generalized straightforwardly to the case of cusp forms.To this end, we adapt the calculation (D.3): let us assume that the spin m spectrum of Maass cusp forms contains a linear ramp.As discussed in the main text, the statistical averaging over cusp form data, which is automatic in the large y i limit, means that there are different choices of correlations which would all yield a linear ramp in some given spin sector.Ultimately, we found (4.1) by demanding a ramp with the correct slope in every spin sector.That is, we demanded that the imprint of any spin sector is the same on any other spin sector.In this appendix we analyze the consequences of working with less fine-tuned spectral correlations that are engineered to only describe RMT statistics in a fixed spin sector.In particular, consider the most naive ansatz, obtained by taking (3.9) and simply dividing out the Fourier coefficients: Such a correlation implies that the spin (m 1 , m 2 ) sector must contain the following term: Z m 1 P,disc.,±(y 1 ) Z m 2 P,disc.,±(y 2 ) ⊃ m n 1 ,n 2 ⟨z n 1 ,± z n 2 ,± ⟩ spin m ramp naive' a (n 1 ,±) m 1 a (n 2 ,±) This can be evaluated numerically, which yields similar results as in the case of Eisenstein series discussed in the previous subsection: Z m 1 P,disc.,±(y 1 ) Z m 2 P,disc.,±(y 2 ) ⊃ m η ± m,m 1 × δ m 1 m 2 1 π y 1 y 2 y 1 + y 2 e −2πm 1 (y 1 +y 2 ) (y i ≫ 1) , (D.8) where the spin-dependent coefficient η ± m,m 1 is generally different from 1.This is illustrated in figure 10.From the curves in that figure, we find the following numerical fit: The fact that the prefactor η ± m,m 1 in (D.8) is not 1 means that the ramp at spin (m 1 , m 2 ) is not fully encoded in the ramp at spin m.Random matrix universality in one spin sector therefore does not imply random matrix universality in a different spin sector -as is consistent with general expectations in the theory of quantum chaos.Instead, one must fine-tune the approximation (D.6) in a way that is informed by cusp form data in all other spin sectors.This is achieved by the arithmetic kernel f (n,±) discussed in the main text.

D.4 Independence of the plateaus
Similar to the case of the ramps analyzed above and in [1], here we discuss the numerical imprint of a plateau in one spin sector onto other sectors.We begin with the imprint of a 33 Computing N ± m 1 /N ± m using the finite number of cusp forms available to us, i.e., using (3.17), yields agreement within ∼ 2%.This shows that a still much larger number of cusp forms is required in order to get very close to the theoretical values for η ± m,m 1 .

Figure 1 :
Figure1: A depiction of the statistical approximation to cusp form data.The "spectral staircase" and the erratically distributed Fourier coefficients are replaced with their average values.We will justify this 'statistical' coarse-graining in the time regime that is relevant for random matrix universality.It should be distinguished from 'microcanonical' coarse-graining, which is always required to discuss correlations in the CFT spectrum.

Figure 2 :
Figure2: Numerical verification of the encoding of a linear ramp in correlations of even (left) and odd (right) Maass cusp forms, according to (3.10), for y 1 = y 2 ≡ y.The summation over n is performed up to some cutoff such that convergence is achieved within the displayable accuracy.The plots show that the sum converges to the ramp linear y/(2π) up to an mdependent constant that is subleading as y → ∞.

2 y 2 )Figure 3 :
Figure3: We quantify how the sum over cusp forms indexed by n depends on the terms with small n (and hence small R ± n ).We plot the ramp in the spectral form factor computed using only values of n for which R ± n > R min and normalize it by the complete result.Asymptotically as y → ∞ this ratio converges to 1, no matter how many low-lying values of R ± n we exclude.We show the cases of even (left) and odd (right) cusp forms separately (they are almost indistinguishable).Solid lines correspond to spin m = 1, dashed lines to m = 2.

Figure 4 :
Figure4: We compute the Maass cusp form sum using the variance of the Fourier coefficients instead of their exact values in(3.18).For large y increasingly many Fourier coefficients contribute to the sum over n, which means that their square can be increasingly well approximated by their variance.We therefore reproduce the linear ramp asymptotically (up to a subleading constant shift), c.f. figure2.The left (right) shows the case of even (odd) parity cusp forms.In the odd case the ramps for different spins lie almost on top of each other.Insets show larger values of y.

26 ) 3 .
Note that these properties encode all information about the distributions(3.13).Arbitrary m: For any general integer m, there is a prime factorization m = p k 1 1 • • • p kr r .The Hecke multiplicativity rules (B.8) imply for different primes p it is useful to first solve the problem (3.21) for prime power spins, m = p k .Consider an ansatz of the form f (n,±) p = r≥0 c p,r a (n,±) p 2r (3.28)

Figure 5 :
Figure 5: Comparison of the exact counting function of discrete eigenvalues of the Laplacian with the Weyl law approximation as given in (B.1).

Figure 6 :
Figure 6: The distribution of the first 5832 even and 6282 odd Fourier coefficients for prime spins, compared to the distribution they are equidistributed with respect to.

Figure 7 :
Figure 7: The distribution of the nearest neighbor spacing of the first 5832 even and 6282 odd eigenvalues and spin 2 Fourier coefficients, compared to the Poissonian expectation.Note that the statistics for all other prime spins is similar.

Figure 8 :
Figure 8: The theoretical variance of the Fourier coefficients as a function of m, up to the 100th prime.The prime coefficients have variance 1 + 1m , while the composite primes have much more complicated behaviour (which is determined by the Hecke relations and higher moments of the distribution of prime Fourier coefficients); in this range, all variances are O(1), but as m grows, the maximum possible value grows as well.

1 )
Proof: (i) and (ii) follow immediately from the Hecke algebra (B.8).(iii) follows by viewing (ii) as a recursion relation and solving it in terms of a (n,±) p .(iv) follows from squaring and simplifying (iii).

Figure 10 :
Figure10: Numerical evaluation of (D.7): we show the imprint of a ramp in the spin m = 1 sector onto the spin (m 1 , m 2 ) sector of the spectral form factor, assuming the naive form of correlations (D.6) which is not engineered to have information about any other spin sector.While the asymptotic contribution has the correct linear y-dependence (for m 1 = m 2 ), it does not have the correct slope to account for all the information encoded in a ramp.

Figure 12 :
Figure 12: Plot of the imprint of the spin 1 plateau on other spin sectors for y 1 = y 2 = y.We again normalize by ⟨ρ D (E m )⟩e −2πm(y 1 +y 2 ) ; at large c, ⟨ρ D (E 1 )⟩/⟨ρ D (E m )⟩ ≈ e −2π(m−1) .The result is a function that grows like √ y; creating a similar plot for a fixed y 2 yields the functional form in (D.14).

Table 1 :
Exact values of the moments of the distributions of Fourier coefficients.For prime m, these follow from (B.5).For composite m, they are constructed from prime moments using Hecke relations.
r a .h.s. is the k-th Catalan number.Using a general recursion relation of the moments of Fourier coefficients (Lemma 4(iv)), we can write the r.h.s. as follows: ± z n 2 ,± spin m ramp naive' ≈ 1 δ m 1 m 2 δ n 1 n 2 .(D.6)