Chord diagrams, exact correlators in spin glasses and black hole bulk reconstruction

The exact 2-point function of certain physically motivated operators in SYK-like spin glass models is computed, bypassing the Schwinger-Dyson equations. The models possess an IR low energy conformal window, but our results are exact at all time scales. The main tool developed is a new approach to the combinatorics of chord diagrams, allowing to rewrite the spin glass system using an auxiliary Hilbert space, and Hamiltonian, built on the space of open chord diagrams. We argue the latter can be interpreted as the bulk description and that it reduces to the Schwarzian action in the low energy limit.

The SYK model at large N turns out to be solvable, in the IR and other limits, and exhibits several interesting features like conformal regime [1] and maximal chaos exponent [2][3][4]6].
In this work, we will be interested in a particular class of spin glass models introduced in [35], which are close relatives of the SYK model, and derive a formula for the exact 2-pt function of certain operators. The model is the following. Consider n sites with a spin 1 2 degree of freedom on each. Denote the Pauli matrices acting on site i = 1, 2, . . . , n by σ (a) i , with a = 1, 2, 3. Given an integer p, we define a random Hamiltonian H (p) as follows. Let e = (i 1 , ....i p ) be a vector of length p of distinct integers defining a subset of the n sites, and let a = (a 1 , ...a p ) be a second vector of length p, with entries being either 1,2 or 3. Denoting the pair (a, e) by J, we define σ J = σ (a,e) = σ (a 1 ) (1.1) and the spin glass Hamiltonian is where the sum runs over all possible J's, and α J are independent Gaussian variables with zero mean and unit standard deviation (we will drop the superscript p from now on). The relevant parameter controlling the asymptotic density of states is [35] q = e −λ with λ = 4 3 p 2 n , (1.3) and the exact asymptotic density of states of the model (1.2) was computed in [35] in the limit λ fixed, n → ∞ (1.4) We will refer to this as the λ-scaling limit. We will be interested in the limit of λ → 0, where the distribution of eigenvalues approaches a Gaussian distribution (point-wise) and hence we will refer to these models as "Almost Gaussian" spin glass models. The spin glass model (1.2) is quite similar to the SYK model. Apart from replacing Majorana fermions with Pauli matrices, the more critical difference is that in the λ-scaling limit (1.4) p is scaled with √ n, leaving λ as a parameter, whereas in SYK, p is held fixed as n → ∞, while scaling the energies properly to obtain a solution of the model. However, the λ-scaling model with Majorana fermions was discussed in [18] where it was shown to have a low energy limit λ → 0, E → 0, (dubbed "double scaled SYK" model) where the density of states is that of the Schwarzian theory. Hence, models (1.2) can just as well be used to study the physics of AdS 2 . In this work, we discuss the full Almost Gaussian model and use the "double scaled limit" to check our results. The main results in the paper are • We motivate why random operator observables are relevant for black hole physics, i.e., not just random Hamiltonians. This is done in section 2, where we also survey existing results and state the new result on the 2-pt function.
• A new method of computing the distribution of eigenvalues of the Hamiltonian (1.2) in the λ-scaling limit. The new method relies on the reduction in [35] of the spin glass Hamiltonian to chord diagrams but then takes a different route in evaluating the latter. This is done in section 3. Section 4 is an analysis of the λ → 0 limit which parallels Appendix B of [18] in our notation.
• The derivation in section 3 relies on an auxiliary Hilbert space and a Hamiltonian acting on it, which we denote by T . This new Hamiltonian is equivalent 1 to the full Hamiltonian of the spin glass in that whenever the unitary operator e iHt appears, acting on the original Hilbert space, it can be replaced by e iT t acting on the auxiliary Hilbert space. In section 5 we suggest that this is the analogue of the bulk Hamiltonian and show in what limit it reduces to the Schwarzian effective action in its Liouville form.
• In section 6 we compute the exact time dependent 2-pt function of an additional random operator of length p ′ ∼ √ n. This can be reduced to another chord partition function in which one chord is marked. We use the technique developed in section 3 to evaluate it.

Motivation, Setup and Summary of results
We will analyse the spin glass Hamiltonian model (1.1)-(1.2). However, we will probe it using a random operator. The latter will be of a similar statistical type as the Hamiltonian, i.e. it will be defined by the same equation (1.1)-(1.2) but with • a different length parameter p ′ = p, and • a new set of independently drawn coefficients.
In subsection 2.1 we motivate this specific choice of operator. The rest of the section is an "executive summary" of the setup of the model and known results in 2.2, and a summary of the new results in 2.3.

Motivation -random observables and factorization
2.1.1 Why random operator probes ?
Since black holes share some properties with chaotic systems [36][37][38][39], they can be thought of as described by a suitable random Hamiltonian. In particular, for AdS black holes, we might want to think about some core of states in the spectrum governed by a random Hamiltonian, describing the near horizon black hole physics, dressed by a "structured" non-random Hamiltonian describing excitations well separated from the horizon. In this picture, one needs to specify the statistical class of the random Hamiltonian. This is precisely what the SYK model achieves, as the relevant class for nearly-AdS 2 spacetimes. The next step is to probe the black hole (BH) using the available bulk probes, such as single trace operators or their analogues. The Hamiltonian in quantum mechanics, or the local energy-momentum tensor in higher dimensions, is one such operator. Probing with the full Hamiltonian does not provide any more information beyond the partition function, but the local energy momentum operator does. In practice, it is another massless field for which we can put sources on the boundary. Just as the full Hamiltonian is indistinguishable from a random operator when acting on the BH states, we can expect that the local energy momentum operator will also be effectively described by some random (local) operator acting on the Hilbert space of the states of the black hole.
But the local energy-momentum tensor is just one of a tower of single trace operators with which we can probe the system. In N = 4 SYM we can use its primary tr(X 2 ) to probe the black hole, or we can just as well use any other of the tr(X n ) operators. If the former is a random operator on the states of the black hole, why should we not expect that all single trace operators be of a similar nature ? We would like to suggest that the relevant probes appearing in General Relativity are random operators on the BH states 2 . The main issue would then be from what ensemble these operators are drawn. If we have some idea about the statistical ensemble of the Hamiltonian, we can try and guess what is the ensemble for the other single trace operators.
Another way to phrase the argument is that the SYK model is dual to AdS 2 in an appropriate large N and energy regimes. But there are other models which realize the same universality class (for example, the one discussed in this paper is based on different spin matrices). So there may be many ways of defining the statistics of the random Hamiltonian which give rise to the same physics -some may be similar to SYK and others may be different. Focusing on the computation of specific operators used to define a specific realization, such as χ i in the SYK model, certainly yields the maximal amount of information about the model but it may not be universal enough throughout the different models. Rather, motivated by the fact that the local energy-momentum tensor is a one "single trace operator" out of many, we would like to suggest that useful probes are random operators appropriately made out of the basic constituents of the theory, just like the Hamiltonian is. The statistical class of these random operators may be more universal throughout the different ways of building models (as we will see in our case).
Yet another argument is the following. In the SYK model, the Hamiltonian is a sum of finite rank polynomials of the χ fields with random couplings. Viewing the χ's as the analogues of the single trace operators in higher dimensions (which is anyhow problematic since they live in SO(N) representations) implies that the Hamiltonian in the black hole regime can be written as a sum of polynomials of single trace operators. This seems to be a very strong assumption for the higher dimensional AdS/CFT dualities. A weaker assumption is that both the local energy momentum tensor and all other single trace operators can be written using some other operators which act on the BH states, which are just used to define the statistical class of the operator and probes. These operators need not be asymptotic observables outside the black hole but rather they just need to be a rich enough set to allow for the correct definition of the statistical class of the observables. This is somewhat against the usual application of the AdS/CFT correspondence where, in this context, the SYK model is taken to be the microscopic theory which defines all of spacetime. In this approach, one is committed to all the operators defined in the model. However, in practice if one is interested in the AdS 2 part, one glues it to an external region in order to break conformal invariance (and the gluing might eventually vary if, for example, one thinks about an AdS 2 near horizon of an object in higher dimensions). It is not clear to what extent the full SYK provides an extension which has an adequate gravity dual outside the AdS 2 region, and even if it does, it is not clear whether it is universal. The right probes on AdS 2 are determined just as much by this outside-of-AdS 2 region since the probe must be defined on the boundary. This means the choice of right probes in the AdS 2 region, within a given model, might be ambiguous in general.

What random operator probes ?
Having argued that random operators are suitable probes, with ensembles related to the one from which the Hamiltonian is taken, in this subsection we would like to discuss another constraint on the ensemble from which probes are drawn, originating from requiring factorization of correlation functions. We will see that it again points us in the direction of almost Gaussian random operators, similar to H (p) .
Within the AdS/CFT correspondence correlation functions of single trace operators factorize at leading order when evaluated in the ground state or in any other state well described by a semiclassical background. This is usually taken to include black holes, although this assertion is on less solid footing there, as the detailed quantum state of the black hole may matter (and surely does over long enough time scales). So the extent to which correlation functions do not factorize will teach us about the role of the quantum state of the black hole, and may also teach us about deviations from the standard Einstein-Hilbert low energy effective action.
In the field theory side of the AdS/CFT correspondence, factorization is a consequence of the large N limit when evaluated around the ground state. Around the black hole background it implies a non-trivial constraint on the statistics of probe operators [40]. Consider a microcanonical ensemble with a small enough energy spread, and consider the 4-pt function 3 Tr where the trace is over the states in this energy band. Since black holes are strongly mixing systems, one might have expected that M -when acting in this energy band -will be described by one of the ordinary random matrix ensembles. An example of this is the often assumed strong form of ETH where the most straightforward interpretation of this formula is as a statement about statistics of the matrix elements 4 . This relation actually comes about by a minimal set of assumptions -that A) only pairwise contraction of the operators matter -after all, we would like to obtain factorization, and B) that all the states in the microcanonical energy band are equivalent and hence the statistics should have full unitary invariance in this energy band -this is also an assumption often made in statistical physics. However, under these circumstances correlation functions do not factorize properly. The ansatz in (2.2) is the same as drawing the operator M from a distribution with measure where N is the number of states in the energy band. So we only need to compute a Gaussian integral. With this measure, in the large N limit, the 4-pt function (2.1) receives only one (planar) contribution. However, factorization implies that there are two contractions. It seems difficult to remedy this within the ordinary ensembles (for example, by changing the measure to e −N V (M,M † ) for a more general V ).
Since there are restrictions to implementing factorization in the simplest ensembles, it is interesting to find additional examples in which correlators factorize. More precisely we would like them to almost factorize -the deviation from exact factorization is then interpreted as bulk interactions. At the level of a single operator, the most naive indicator of factorization -neglecting for the moment the issue of time dependence -is that for an hermitian operator M, where E() is the statistical average over the ensemble from which the operator is drawn (and A is a factor set by the normalization of the operator). The ensembles in [35] are precisely of this type. For any operator of the form (1.1)-(1.2), the distribution of eigenvalues approaches the Gaussian one in the limit λ → 0, so all operators with p √ n ≪ 1 will be approximately Gaussian. If the Hamiltonian has a specific (small) λ, then operators for all other values of (small) λ are in qualitatively a similar statistical class and approximately factorize. We will use them as our probes.

Set up of the model and summary of known results
The model discussed in [35] is defined in equations (1.1) to (1.4). One of the main results in that paper is that the asymptotic distribution of eigenvalues, in the limit n → ∞ λ fixed, (2.5) 1−e −λ , and vanishes outside this region. The proof proceeds by computing the moments in the following steps: 1) For the first step one needs to define what are chord diagrams. Consider L = 2n dots on a circle -a chord diagram is a pairing of these dots into n pairs. We draw a line connecting each paired dots. i.e., a total of n lines. Denote a specific chord diagram by π. We then denote by k(π) the number of crossings of lines (when we draw the diagram such that each pair of lines intersects at most once). An example of a chord diagram is shown in Figure 1 with n = 8 and with number of crossings k = 2.
In the first step one shows that where the sum is over all the chord diagrams. The expression on the RHS is called the chord partition function, and q = e −λ was defined before in (1.3) in terms of the parameters of the spin glass. For example, the contribution to the sum by the chord diagram shown in Figure 1 would be e −2λ . Chord diagrams were also used in [23] for computing 1/N corrections in the SYK model. In section 3.1 we review this step of the proof in more details since that part of the proof will not change. Furthermore, we will also need to slightly temper with it when computing the 2-pt function.
2) In step 2, one uses the Riordan and Touchard formulae [42,43] and the results of [44] to show that (2.8) are the moments of the distribution (2.6) and further give an explicit formulae for the moments as

Summary of new results
In this paper we discuss a new proof for the value of m L and the energy eigenvalue distribution of the spin glass. We use this to compute the exact two point function for random operators, in the limit λ fixed, n → ∞. Denoting a new random operator by M , it has the form (1.1)-(1.2) (with new randomly chosen coefficients, uncorrelated with those of the Hamiltonian as mentioned in the beginning of section 2) but with a new parameter length parameter p ′ ∝ √ n. More precisely we show that pp ′ n and (a, q) ∞ is the q-Pochammer symbol (see (A.2)). To prove this one evaluates (2.11) We show that the relevant chord diagram which computes this two point function is a chord diagram in which one of the lines is marked, and intersections with this chord are assigned a different weight. An example of a marked chord diagram is given in Figure 2. More precisely : • Given 2n+2 points on a circle, two specific points are connected. This is the "marked" chord. The thick line in Figure 2 connecting the red dots represents the marked chord.
• Between the special points at the ends of the marked chord there are k 1 regular points on one side, and k 2 regular points on the other side (k 1 + k 2 = 2n).
• These remaining 2n points are paired. These will be called "regular" chords.
• Intersection between regular chords is assigned weight q, and the intersection between the regular and marked chords is assigned weightq.
• The marked chord partition function is defined as a sum over pairings of the 2n regular points, with k 1 and k 2 fixed and with weights as above, i.e., where k regular (k marked ) is the number of regular-regular (regular-marked) intersections. For example the 1-marked chord diagram in Figure 2 contributes qq to the m 1,5 .
• Similar to [35] we show that and evaluate the right hand side to obtain (2.10) The evaluation of the various chord partition functions in this work relies on an auxiliary Hilbert space space where there is a natural Hamiltonian whose action is equivalent, in a sense that will be made precise below, to the one of the full Hamiltonian acting on the spin glass Hilbert space. We interpret this auxiliary structure as the bulk dual to the spin glass. Furthermore, we suggest how it reduces to the Schwarzian action in its Liouville form at low energies 5 .

A new derivation of eigenvalue distribution
Given a random Hamiltonian as in (1.1), the authors in [35] compute the asymptotic distribution of eigenvalues in the λ-limit (1.4), by evaluating the moments and by finding the unique distribution compatible with them. E stands for an ensemble average. In section 3.1 we review how [35] reduces the moments (3.1) to evaluating the chord partition function (2.8). [35] then uses the results of [44], and also the formulae of Touchard and Riordan [42,43], to show that the moments (2.8) arise from the distribution given in (2.6), and to give the explicit formulae (2.9) for the moments. Our proof, in section 3.2, replaces this second step, as well as generalizes it to other Chord diagrams, as the one that will appear in the exact 2-pt function in section 6.

From spin glasses to chord diagrams
Given the Hamiltonian (1.2), the computation of the moments proceeds by evaluating the ensemble average of the α J 's first. Any non-vanishing contribution requires at least two insertions of each α J . Moreover, Lemma (4) in [35] shows that the dominant contribution, in the λ-scaling limit (1.4), is when J 1 , ...J L appear exactly in pairs, with higher multiplicities being subleading in n in the large n limit. This gives us the basic structure of chord diagrams where pairing in the chord is defined by having the same J on two different nodes as in Figure 1. Summing over all the relevant values of J amounts to summing over all chord diagrams (i.e all possible pairings of J's), and then sum over all value of J (i.e., both e and a) for each chord. Given a chord diagram we therefore need to evaluate what is the weight that is associated with it, i.e., where there are only L/2 independent J's and the pairing is determined by the chord diagram. The obstruction to immediate evaluation is that σ a i for the same site index i can appear in different J's. However, [35] shows that with probability 1, in the λ-scaling limit, each node can appear in at most two of the chords, enabling the evaluation of the weight.
More precisely, define the intersection of J's by the intersection of the site index, i.e.
(3.4) [35] shows that, for a given J i and J j the size of the overlap is Poisson distributed, and that there is, with probability 1, no triple intersections. I.e., we can assume This statement is summarized in lemma (9) there, and subsequent discussion. Given two sets a and b of integers drawn out of the set {1, 2, ...n} (without repetition in each set), we can think about it as |a| independent processes in which the overlap between the sets increasing by 1 with probability |b| n (in the limit n → ∞). This is a Poisson distribution with mean size of overlap 3λ 4 = |a||b| n . Recall that we scale the size of the set with √ n so this remains finite in the limit n → ∞. The average size of an overlap with an additional index set -say c -is the latter times |c| n → 0, so with probability 1, triple overlaps are empty.
The interplay of chord intersections and overlap of the index sets is the key for evaluating the weight of each chord diagram. As we sum over the J's of the chords, then if two chords do not intersect they will contribute whereas if the chords intersect they give a factor proportional to If there is a non-trivial overlap J 1 J 2 = 0, then these factors will be different. So there is some "penalty" that we pay for each intersection. More precisely, given a chord diagram (i.e. a pairing π), recall that k(π) is the number of pairwise chord intersection. Each chord intersection has a Poisson distributed overlap of sites. Each overlap is independent of the overlap of the other chord intersection. Each overlap of sites (for a given intersection) comes with a factor relative to 1 when the ordering is (aabb) which originates from an overlap in a pair of non-intersecting chords. Therefore, the size, m, of each overlap is Poisson distributed with expectation value 3λ 4 and comes with a weight (− 1 3 ) m . The expectation value of the weight for each chord intersection is therefore e −λ , and the total weight associated with each chord diagram is e −λk(π) . Hence, one finally obtains (2.8).

Evaluation of the Chord partition function
In this subsection, we will provide a alternative derivation of the chord partition function reproducing the expression for v(E|q) in (2.6). The proof is rather compact, generalizes to more complicated chord partition functions, such as the ones discussed in section 6 6 and suggests a bulk interpretation that we develop in section 5. The evaluation of (2.8) is based on a "hopscotch" recursion relation satisfied by the concept of a partial, or open, chord partition function as follows: • 2 −n Tr(H L ) involves L points in a chord diagram, as indicated in figure 1. Choose one point, i.e. choose one of the H factors, to be the first and begin moving clockwise in the chord diagram. Each time one reaches an extra point and hop over it, we shall refer to it as "a step". As we go along, denote the number of such steps by i, i.e. the number of H factors that were hopped over. In step 1 we hopped over the factor of H that we chose to be the first.
Here k p (π) refers to the number of chord intersections to the left, in the past of our "hopscotch" process. It is convenient to think about v (i) as a column vector and l as its index.
Given this set-up, one can write down a recursion relation for v i . At each step, one can either close a chord (as in Figure 3) or start a new one (as in Figure 4) at the point one is hopping over. If one starts a new line, l changes to l + 1 and the new line enters at the bottom. If one closes a line, it can be either of the l open chords with height between 1 and l. If one closes the line at height p, it crosses (p − 1) lines on its way down. This crossing generates a weight q p−1 when evaluating its contribution to the partial chord partition function. Altogether, the vector of such partition functions satisfies the following recursion relation v with initial condition v (0) l = δ l,0 . The latter can be rewritten in terms of an (L+1)×(L+1) transfer matrix T (L) propagating the partial chord partition function forward with matrix elements (indices running from 0 to L, l 1 (l 2 ) is the row (column) index) describing a matrix with 1's and η l 's on the diagonal below and above the main diagonal, respectively, i.e.
To compute the chord partition function (2.8), define the vector of length L + 1, and then The initial condition v (0) l = δ l,0 dictates the use of the initial state |0 L . Ensuring our procedure counts only chord diagrams that close by the time we reach the L-th point, such that we are computing the usual chord partition function in which all lines are paired, determines the final state.
Notice that we are computing the trace of H L in the original 2 n dimensional Hilbert space, using some auxiliary space based on partial chord diagrams. We shall develop a "bulk" interpretation for the latter in section 5.
Next, given some fixed L, one can always consider a larger L'-sized Hilbert space This allows us to take L ′ → ∞, keeping L fixed. In this infinite dimensional Hilbert space, one can define Hence, T is the infinite dimensional extension of (3.13). This provides an auxiliary Hilbert space and a single matrix T in which one can evaluate all traces as The problem of computing the moments (2.8) reduces to the problem of computing the eigenvalues α of the operator T and expanding the vector |0 in terms of these eigenvectors |α , in an expression where Spec(T ) is the set of eigenvalues, ρ(α) is its density and ψ 0 (α) ≡ 0|α is the overlap of |0 with the |α eigenvector of T . Fortunately, Spec(T ) and the density are very easy to compute and the overlap is given by specific q-Hermite polynomials, as we will see below.
In the notation of the spin glass model, comparing the L dependence in the original moment (3.1) with α L in equation (3.19), suggests the identification where E is the energy of the system, properly interpreted. The asymptotic distribution of the energies should then be identified as

21)
A short example. It is worth while carrying out the procedure above in an explicit, low L case, and compare the result with (2.9). For example m 4 (q) = 2 + q, which can obtained also from the three chord diagrams in Figure 5. In our approach we start with v (0) , act on it 4 times with T (4) (or T ), and project on v (0) . Keeping track of chord histories give the following partial chord partition functions: The symmetric form of the transfer matrix T. The matrix T in (3.13) is not Hermitian, but one can conjugate it to a symmetric version by defining a new matrixT where P is a diagonal matrix with entries (P 0 , P 1 , P 2 . . . ) satisfying where (a; q) l is the q-Pochammer symbol (see (A.1)).T has matrix elements Thus, it is manifestly symmetric, and has the same original moments (3.18) since We will switch between the two transfer matrix descriptions depending on which is more convenient at each stage.

The spectrum of T
Obtaining the spectrum of T is straightforward. The matrix T asymptotes, down the diagonal, to a matrix with 1 one diagonal below the main diagonal and η ∞ = 1 1−q one diagonal above the main diagonal. i.e.
We can think about the eigenvalue problem of T as a scattering problem with the distance along the diagonal playing the role of position. In this interpretation, infinity is captured by the asymptotic form of the operator T far down the diagonal. Hence, this is a scattering problem on the half line with δT acting as a scatterer close to the origin. Indeed, up to an overall rescaling, by conjugating the matrix T asymp , and adding the identity matrix with an appropriate weight, we can bring it to the form with -2 on the diagonal and 1 on the diagonals below and above the main. It is then an approximation to the 2nd derivative operator, making the asymptotic behaviour more familiar in the continuum limit. This interpretation is elaborated in section 5, where the connection between this eigenvalue problem and the Liouville equation is described. However, as far as the spectrum and its density is concerned, the details of the scatterer are not important as both can be read from the behaviour at infinity 7 . So the spectrum of T is the same as that of T asymp which is a Toeplitz tridiagonal matrix (i.e with constant elements one diagonal above and below the main diagonal [45]), for which there is a simple formula for the eigenvalues, which in this case is as was found in [35]. This formula gives us both the spectrum of the T and the density of states on it. Denote θ = sπ n+1 , then in the limit n → ∞ it covers the interval [0, π] with uniform distribution, i.e., inserting a complete set of energy eigenstates is simply done by the replacement (3.31)

Eigensystem of T matrix
The previous asymptotic discussion suggests to parametrise the eigenvalues of the matrix T as E(µ) ≡ 2µ √ 1−q . Let v (µ) be the corresponding eigenvector 8 . This allows us to write (3.11) together with the recursion relation (3.10) as Just for this subsection, we will allow l = −1 and define v For the wave functions, or form factors, we will be more specific below. Also, reading the spectrum and density from infinity also assumes that there are no bound states near the origin. 8 v (µ) ∝ ψ(α) of section 3. For now the normalization is different though.
Using the full range θ is dictated by the the discussion of the spectrum in section 3.2.1. At this point, it is more useful to switch toT , since we need to conjugate the form factor. The eigenvectors of symmetric transfer matrixT are just P v µ . In components, the eigenvectors arê where N (µ, q) is a normalization which is fixed by the requirement that the states are delta function normalized in θ, which gives (see Appendix B) N (µ, q) = √ (q;q)∞|(e 2iθ ;q)∞| √ 2π . With this one can easily write down the matrix element 9 The moments of the distribution eq(3.19) can then be computed to be where we have defined the distribution Below we will show that this formula is the same distribution given in (2.6).
Matching to the result in [35]. Recall that [35] obtained the moments m L (q) as moments of distribution v(E|q) given in (2.6), i.e

Switching to angular variables via
followed by a change of variables as to obtain that (3.40) matches our result (3.38). 9 Recall that the density of states ρ(θ) is uniform.

The q → 1 limit of the distribution
Our interest in approximately factorized correlators, suggests to work in the regime λ → 0 or q → 1. The analysis is most easily done by arranging the Pochhammer symbols into Jacobi Theta functions and performing modular transformations. The results in this section are similar to Appendix B in [18] and [24] after a suitable substitution that takes us between our model and the SYK model in the same scaling as above.
To study the λ → 0 limit of the distribution Ψ(θ, q) (3.39), it is convenient to rewrite it in terms of Jacobi Theta functions. Using (A.4), Using the modular transformation (A.5) we rewrite it as (4.2) and the λ → 0 limit becomes The last equality follows since the exponential overcomes the hyperbolic cosine factor for m ≥ 1 given that θ ≤ π. Plugging the above λ → 0 expansion in (4.2), one gets where in the last step we used θ ≥ 0. This determines the dominant contribution to the distribution (3.39) to be Notice this function is symmetric under E → −E, and vanishes at the edges E max = −E min = 2 √ λ , which correspond to θ = 0 and θ = π, respectively, since E = 2 cos(θ) √ 1−q . The distribution (4.5) has several interesting regimes: • The λ → 0 with E fixed regime. As highlighted in [35], pointwise, which is the Gaussian limit of the distribution (4.5).
• The other interesting behaviour is close to edges. Setting ϕ = π − θ, we begin with the limit ϕ = π − θ ∝ λ (4.7) (where by this we also include ϕ λ ≫ 1 fixed, λ → 0). In this case the distribution becomes The quadratic term in the exponential can also be neglected in this regime, giving rise to the density of states of the Schwarzian theory (∝ sinh(2π (E − E min )/λ 3 2 )) after recalling that near the edge • Note that we can actually neglect the quadratic term in the exponential already at λ ≪ ϕ ≪ √ λ, and extend the Schwarzian regime. This points to another simplification of the spectrum which actually covers the bulk of the spectrum at λ ≪ ϕ ≪ π−λ. In this range we can also expand the second sinh, and obtain that the distribution is just a gaussian in (ϕ − π/2). The center of this range includes the Gaussian-inenergy distribution, and its edges overlaps with the Schwarzian distribution. It would be interesting to find a symmetry argument for this entire range.

The canonical ensemble in the q → 1 limit
Similarly one can analyze the canonical partition function in the limit q → 1. Using the variable ϕ ≡ π − θ as before We can treat most of the spectrum using the discussion in bullet 3 above, leaving out only a very low temperature regime where φ ∼ λ. We will prefer however to split the discussion according to first two bullets, i.e., to a high temperature phase and a Schwarzian phase which then splits into a low temperature and a very low temperature phase. Both of the latter are obtained from the Schwarzian density of states and go smoothly into each other, and we make this division mainly for the sake of the discussion of the 2-pt function in section 6, for which the difference between these regimes is more meaningful.
• High Temperature phase (when λ − 1 2 ≫ β) : localizing in the region |ϕ − π 2 | ≪ 1, reduces the partition function to a Gaussian around E = 0, and the partition function can then be written as It is clear that the gaussian cuts off the integral if ϕ deviates from π 2 . To evaluate Z(β), we set π 2 − ϕ ≡ x. The limit above translates to x ≪ 1, in which case we approximate the cosine by x 2 to obtain Z(β) ∼ e β 2 2 and the integral is supported at x = β √ λ 2 with a width of order √ λ.
• Low Temperature phase (when λ − 3 2 ≫ β ≫ λ − 1 2 ) : generically, one expects ϕ ≪ 1, but the thermodynamic behaviour of the system is sensitive to how small ϕ is compared to λ, due to the argument in the sinh factor in (4.8). Consider the regime λ ≪ ϕ ≪ 1, where the distribution is approximated by Expanding the Boltzmann factor, the partition function reduces to 10 Notice the integral is mainly supported near ϕ = π β √ λ with a width of λ 1 4 √ β . The consistency with the assumption λ ≪ ϕ ≪ 1 requires β ≪ λ − 3 2 . This regime, along with the next one, are part of the conformal low energy limit of the theory.
• Very Low Temperature phase (when β ≫ λ − 3 2 ) : consider the regime ϕ ≪ λ. After linearising both the sin and sinh factors, the distribution (4.8) simplifies to 2λ − 2ϕ 2 λ ϕ 2 (4.14) The Boltzmann factor in the partition function cuts off the integral around ϕ ∼ λ 1 4 √ β . This is consistent with our regime ϕ ≪ λ, since β ≫ λ − 3 2 . The partition function can then be evaluated as where in evaluating the integral, we have replaced the upper limit by ∞.
Relation to previous work. The low energy behaviour identified in (4.8) is the one discussed in appendix B in [18] and in [24]. To make the comparison with [18] easier, notice the density of states (3.41) can be written as (4.16) 10 The term − 2ϕ 2 λ in the exponent is negligible compared to −βϕ 2 √ λ due to β √ λ ≫ 1.
where µ = cos θ and N = 1−q 2k+1 . This matches equation (81) in [18] by identifying their parameters a, λ s , J as a = µ , J = √ λ e λ/8 , λ s = λ 2 (4.17) whereas both energies are the same 11 . The second matching follows from the observation that our variances equal unity, as in [35], whereas our normalisation was TrH 2 /TrI = 1 (see equation (80) in [18]). The third matching is due to the Majorana nature of the fermions in the SYK model. Thus the density of states in [35] is exactly the same as the doubled scaled SYK, up to these identifications. The further triple scaled limit, isolating the Schwarzian action in SYK, corresponds to the low energy behaviour captured by the density of states (4.8) in our set-up.

Bulk reconstruction
In section 3 we presented a new derivation of the density of energies in v(E|q) in the λ-scaling limit (1.4) keeping β finite and where E() on the left hand side is the average over the ensemble of Hamiltonians. The range of integration on the right hand side is the spectrum of the random Hamiltonian, and its "randomness" now hides in the 1/n corrections which are neglected in this limit. Loosely, one can hope that for a specific realization of the Hamiltonian H, one can write with probability 1 (or 1 − O(1/n)) on the space of random Hamiltonians, in the large n limit. In this case one is dealing with a specific Hamiltonian on the left hand side. This single Hamiltonian realization corresponds to the boundary field theory Hamiltonian in the AdS/CFT written in terms of the fundamental field theory objects -in our case the spin operators. The operator is random and only in the n → ∞ limit its spectrum converges to anything universal. In this section, we suggest that the operator T (orT ) is the bulk Hamiltonian, i.e. the analogue of the bulk Hamiltonian for the near-AdS background -whose low energy limit is given by the Schwarzian action -but extended to the full model. Recall that the parameter E appearing in the right hand side of (5.1) and (5.2) can be reinterpreted as the energy of 11 Our normalizations are different from [18] since their distribution ρs(E) integrates to 2 N 2 whereas our v(E|q) integrates to 1 the field theory Hamiltonian, but it is also the eigenvalue of the operator T (orT ) which acts on the (altogether different) Hilbert space of weights of open chord lines. Whereas the spectrum of H changes from realization to realization, the matrixT is fixed. There is no contradiction since we work in the limit n → ∞, fixed λ limit, where the spectrum of H is universal.
Furthermore, the operatorT can be used as the Hamiltonian not only for the partition function, but for a much broader set of computations. It should be clear that the insertion of any finite polynomial of H in expectation values involving density matrices of the form for any analytic weight function f (E), can be turned into the insertion of the same polynomial withT as its argument, following the procedure described in section 3. In other words, the insertion of e −itH in expectation values involving (5.3) can be exactly replaced by e −itT , while the density matrix itself is mapped into the density matrix (as an operator in the Hilbert space defined on the chord diagram side) Having two different Hamiltonians, acting on different Hilbert spaces but propagating the system in exactly the same way, supports the dual interpretation we suggest forT . This means that we can access a large set of weights on the energy eigenstates as long as the function is smoother than the energy spacing (actually smoother than 1 n for the entire energy band). This is not in contradiction with what we know about the bulk Hamiltonian (anything which extends the low energy effective action), since it is not clear that it should be able to capture states whose support on close by energy states is rapidly varying 12 .
Phrased differently we regard E, when used as the eigenvalue ofT , as a parameter which scans over the allowed energy range only after taking the limit n → ∞. It is not the discrete spectrum of energies of the finite n system. It should be viewed as a coarse grained version of the latter, very much like the energy measured in gravity is a coarse grained version of the discrete set of energies of the field theory (when defined on a compact space). Going from the eigenvalues ofT to the eigenvalues of H at finite n is an interesting problem, and it is similar to seeing -in General Relativity -the discreteness in energies of a black hole.
The above discussion, together with the behaviour of the partition function in the low temperature regime, suggests the low energy physics for q → 1 should be governed by the Schwarzian action (in the gravity dual), as in the SYK model. In the following, we derive this connection by matching the continuum limit of the equation determining the spectrum ofT with Liouville quantum mechanics 13 , which can be written as the Schwarzian action, as discussed in [46,47] 14 . 12 Unless, for example, one believes in the microstate program in its strongest form where one can choose a specific energy eigenstate in the most extreme case. 13 We would like to thank D. Bagrets for a discussion of this point. 14 See [48] for a 2d CFT perspective on this matter.
To take the continuum limit, it is convenient to define the matrixT ≡ ST S −1 where S is a diagonal matrix with entries S ii = (−1) i . Notice that solving for the eigenvalues of theT matrix still resembles a scattering problem on the half line, with the index i of the vector measuring the distance from the origin, just like it did for theT , T matrices. The asymptotic form of theT matrix is In the continuum limit, the above matrix includes the second derivative operator. To make this more precise, define φ = log(q)i (5.6) Using the form of theT operator in (3.25), its continuum limit equals Notice the potential term comes from the expansion 1−q i+1 (3.12), which is accurate since i is large and q → 1, from below.
The eigenvalue problem then reduces to the quantum mechanical eigenvalue problem This is equivalent to the Liouville form of the Schwarzian action in equation 32 in [47], given by after a constant shift of φ. In [47] M was the scale M = N log N 64J √ π (for the SYK model with quartic interactions). For us it is set by | log(q)| −2 ∼ λ −2 .
The prescription in [46,47] (and in [48] for 2D case) requires that, in the path integral, we sum over trajectories that begin and end in the strong coupling region φ → ∞. This is in qualitative agreement with our prescription since we place the state v 0 as initial and final states. Recall that v 0 = (1, 0, 0, 0....), i.e., only the i = 0 term is turned on, which where the term q i is the largest. In terms of φ, e φ is largest which is indeed the analogue of the Liouville strong coupling region. The models are of course not exactly the same since the model in [46] captures the low energy and theT matrix captures the full dynamics.
This also gives an interpretation of the index i via its relation to φ. φ(t) measures where the AdS 2 space is glued to whatever non-universal UV we have (the leading effect being the Schwarzian action), i.e., φ(t) parametrizes the length of the AdS 2 throat. We see that in the full model the size of AdS 2 is actually quantized, giving rise to a minimal size AdS which corresponds to the state v 0 .
It is worth reiterating that the density of states for the H Hamiltonian is different than the density of states for the matrixT , even in the large n limit. Rather the density of states in the former is related to the density of states in the latter H by equation (3.21) or, equivalently, by which means that we have to put a specific initial and final states forT in order to compute the partition function. It is tempting to interpret this in Minkowski space as a computation with initial and final states at the past and future singularities of the black hole.
6 The two point function

The exact 2-pt function
As explained in section 2, we want to compute correlators of random operators M taken from the same universality class as the Hamiltonian (1.2). Hence, these are defined by where J is now a string of p m distinct sites and Pauli matrices. The sum runs over all such possible J's, and m J are independent Gaussian variables with zero mean and unit standard deviation (in particular they are also independent of the coefficients α J in H).
There are two relevant parameters that we will keep fixed in the limit n → ∞. The first is the analogue of λ (see (1.3)) for the random operator (6.1) for any value of β and t. The formalism developed below, based on the set-up in section 3, proceeds by evaluating, and then resumming, expressions of the form This formalism can be extended to compute any n-pt function [49]. The strategy is to reduce the computation to some relevant chord partition function, and then to evaluate it. The identification of the relevant partition function follows the discussion in section 3. The Gaussian integration over the random coefficients of the operators still pairs them. Hence, one can still think in terms of chord diagrams. The only difference is that the Gaussian integral over m J 's pairs the two M insertions, whereas . . . . . . To evaluate the marked chord diagram, we need to modify the "hopscotch" procedure described in section 3.
where W (q) = Diag(1,q,q 2 , ....) (6.7) encodes the intersection of an H-line with an M -line, when the former hops over the latter.
In the next subsection we will evaluate this expression for a special case of q andq. But before we do that, we will perform a quick check on our results above.
A check. Before evaluating (6.11) for a special case of q andq, one can perform a check by taking theq → 1 limit. There should be no cost for the H lines crossing the M lines in this limit. Hence, it must be that To check (6.11) is compatible with this behaviour, notice that nearq → 1, (q 2 ; q) ∞ → 0, due to the first term in the product. Hence (6.11) vanishes, unless θ 1 → θ 2 , since an additional zero in the denominator occurs then 15 . Hence, the integrand in (6.11) behaves like a delta function whose strength is given by 15 Another zero may appear in the denominator when θ1 + θ2 = π but this appears in a co-dimension 2 in the range of integration.
6.2 The q → 1 limit withq = q m The exact 2-pt function of M (6.11) holds for all ranges of time (which are held fixed in the n → ∞ limit). In the remainder of this section we will compute the formula in a specific case, which is the low energy regime where conformal symmetry is expected to appear, as discussed in section 4 and in Appendix B of [18]. We will work in the limit q,q → 1 since we want to work in the limit in which the correlators of each operator separately approximately factorize. However, more and more terms contribute in the Pochammer symbols in this limit, similar to what we had for the partition function, and hence it is important how we take this limit. Since in gravity non-factorization of correlation functions for different operators is governed by the same parameters (e.g. the same 1/N), then the rates of q → 1 andq → 1 should be related.
A particularly simple case to analyze isq = q m with m an integer. This has technical advantages, but it is also physically interesting because it corresponds to p m = m p . (6.15) That is, if the Hamiltonian is made out of a sum of strings of p spin operators (with random coefficients), then the random operator M is made out of a string of m · p spin operators. This is reminiscent of the statement that, say for 4D, N = 4 SYM, the Hamiltonian is a descendant of Tr(X 2 ), yet we can probe the system with low energy fields, which correspond to single trace operators of the form Tr(X n ), n > 2, and their conformal descendants. As discussed in [18] and section 4, our model has a conformal low energy limit. Hence, conformal symmetry should assign a dimension one to the Hamiltonian. If the fundamental fields (in this case the spin operators) can be assigned a specific conformal dimension, and if this conformal dimension is additive in composite operators -as in the SYK model on both counts -then one expects the conformal dimension of M to be m. We will see how our exact formula matches this, up to the existence of mixing with operators of lower dimension when we work at finite temperature. Despite this, our exact formula always has an overlap with an operator of the right dimension.
Before doing the computation we would like to recall an additional formula to which we will compare our result. We will actually be computing the "two sided correlator" This computation is slightly easier than the ordinary thermal correlator. We refer to this as the two sided correlator since, in an eternal black hole in AdS, it is the relevant correlator when there is one operator on each of the boundaries. For an particle of mass M in the BTZ black hole this correlator is (see for example [37] in the Eikonal approximation, with a shock wave there) where l is the AdS 3 radius, M is the mass of the particle and M l is the conformal dimension of the associated operator. This is what is expected from conformal invariance. In the single sided correlator the cosh is replaced by a sinh, to obtain the expected short distance behaviour of 1/t 2M l , and the correct Euclidean time periodicity.

The reduced formula
Whenq = q m , the identity allows to write the 2-pt function (6.11) as Notice the finite product can be rewritten, within the integral, as (6.20) Taking the derivatives outside of the integral, allows to write the integrand in terms of ϑ functions (see (A.4)) depending only on q, The q → 1 limit of the Jacobi Theta functions is evaluated as in (4.4), bringing the 2-pt function to the form This is the exact 2-pt function for q =q m in the limit λ → 0. In the next subsections, we study the function I(β, t, q), from which all m > 1 correlators can be extracted, in the low temperature and very low temperature regimes (or long time, and very long time regimes).

Low and very low temperature regimes
Since the integral I(β, t, q) localizes near the edges at low energies, we define φ i = π − θ i . Expanding the integral near φ i ∼ 0, As explained in section 4.1, the low energy (and very low energy) regime satisfies β √ λ ≫ 1. To study the behaviour of the Gaussian factors in the above integral, it is convenient to rescale the integration variable ϕ i ≡ φ i λ together with the time and temperature parameters β ≡ λ 3/2 β,t ≡ λ 3/2 t. The low energy regime is equivalently described by allowing to approximate (6.23) by This integral has two regimes, following a similar discussion for the partition function in section 4.1 : • The low energy-long time regime characterised byβ,t ≪ 1, where the integral receives contributions from the range ϕ i ≫ 1.
• The very low temperature regime, or very long time scale regime, characterised bỹ β ≫ 1 ort ≫ 1, where the integral receives contributions primarily from ϕ i ≪ 1.

Low temperature regime
The low energy-long time regimeβ,t ≪ 1 allows to extend the range of integration to ∞ since the gaussian in the integrand cuts off the integral well before the limits in (6.25). Notice also the integral is supported at large values of ϕ 1 , ϕ 2 , allowing us to approximate three of the sinh functions by their larger exponentials where we changed variables to ϕ = ϕ 1 +ϕ 2 2 , σ = ϕ 2 −ϕ 1 2 . Due to the 1 sinh(2πσ) term, the σ integral receives contributions from finite σ, whereas its limit of integration is ±ϕ, much larger quantities. This means we can trade the σ limits by ±∞. Furthermore, we can also neglect the e −βσ 2 term and the σ 2 term relative to ϕ 2 in the (ϕ 2 − σ 2 ) term. After these approximations, our integral reduces to Sinceβ ≪ 1 we can in any case neglect the ϕ s dependence in the numerator. The integral shows different behaviours depending on the scaling oft : • Whent √β ≪ 1, the ϕ s dependence in the denominator can be neglected and, to leading order inβ, the result is I(β, t, q) = λ where we changed the integration variable to ϕ ≡ ϕs √β + π−t β in the second step. Due to the large t/β, ort/β, (6.31) differs from (6.30) by an additional et 2 β .

Very low temperature
Whenβ = βλ 3/2 ≫ 1, the angles ϕ 1 and ϕ 2 are localized to a range much smaller than 1. This allows to expand the sinh functions in (6.25) to obtain I(β, t, q) = λ 6 2 e where we traded the upper limit with ∞. This may have the following interpretation. This quantity equals where E() is the statistical average and β 1 , β 2 are related to β, t. This means that where E i measures the energy of the state above the ground state. We can interpret this as if the operator M acts as an underlying gaussian random matrix which couples to low energy states with form factors ϕ 2 , i.e. consider a set of random vectors where the sum is up to some energy higher than the scale set by the very low temperature, and c i,α are independent complex Gaussian random variables with mean zero and unit standard deviation. Take M to be a random Gaussian Hermitian matrix in terms of these variables M = α,β |v α M α,β v β | (6.36) where andM are independent complex Gaussian variables with mean zero and standard deviation 1. In this case (6.34) is satisfied.

6.4
The final correlator The evaluation of the exact 2-pt function (6.22) requires to compute the action of the operator D m (t, β) in (6.20) on I(β, t, q) and interpret the result. It is easy to read the results without actually having to worry about the details of D m (t, β).
• The low temperature regime corresponds to the conformal regime, when the fluctuations of the pseudo-Goldstone modes are still small. Assigning the Hamiltonian H the conformal dimension 1, one would expect an operator made of m spin operators to have dimension m.
This is exactly what happens in our formulas. For m = 1, the operator D 1 (t, β) reduces to the identity. Hence, our result (6.30) is the correlator for an operator of dimension 1 i.e., ∼ 1 cosh 2 ( πt β ) . For m > 1, there exists operator mixing, but we can extract the operator content from the correlator as follows. To isolate the conformal dimensions of the participating operators, first insert the operators on the same side, or equivalently take t = iβ/2+t ′ . This turns the cosh into a sinh. Second, take the limit t ′ ≪ β. In this case the leading contribution 16 in D m (t, β) acts on it with ∂ 2m−2 t turning the correlator into 1/t 2m which is the 2-pt function for an operator of dimension m.
• For the very low temperature/long time regime we can compare (6.32) with equation (67) in [47]. Although their discussion is for SYK model with quartic interactions, it is within the Liouville description of the Schwarzian action. Since our spin glass model reproduces the latter in this very low temperature regime, both results should be similar. The finite temperature 2-pt function of a pair of Majorana fermions in the SYK model at long times/low temperatures (in the conventions used in [47]) equals G(τ ) ∼ − M 2 β 1/2 √ J sgn(τ ) τ 3/2 (β − τ ) 3/2 , τ ≫ M ≡ N log N 64 √ π J (6.37) For a 2-pt function of higher dimension operators, the time dependence (at long time and low temperature) remains with the same power, except that the coefficient of power of M in front of the expression increases.
To match with (6.32), one needs to work with Lorentzian time, τ = it and to shift the time imaginary axis by t → t − i β 2 . Altogether, where we analytically continued back to lorentzian time in the last step.