Explorations of Non-Perturbative JT Gravity and Supergravity

Some recently proposed definitions of Jackiw-Teitelboim gravity and supergravities in terms of combinations of minimal string models are explored, with a focus on physics beyond the perturbative expansion in spacetime topology. While this formally involves solving infinite order non-linear differential equations, it is shown that the physics can be extracted to arbitrarily high accuracy in a simple controlled truncation scheme, using a combination of analytical and numerical methods. The non-perturbative spectral densities are explicitly computed and exhibited. The full spectral form factors, involving crucial non-perturbative contributions from wormhole geometries, are also computed and displayed, showing the non-perturbative details of the characteristic `slope', `dip', `ramp' and `plateau' features. It is emphasized that results of this kind can most likely be readily extracted for other types of JT gravity using the same methods. Preliminary results also suggest that a new well-defined non-perturbative completion of ordinary JT gravity using the Hermitian matrix models may exist.


I. INTRODUCTION
There are many reasons to study Jackiw-Teitelboim (JT) gravity [1,2]. One of them is the fact that it is a theory of a two dimensional quantum gravity, where the spacetime is allowed to split and join, changing its topology (characterized by Euler characteristic χ=2−2g−b−c where g counts handles, b boundaries, and c crosscaps). In the full theory the partition function Z(β) is a sum over the contributions from all topologies as well as a non-perturbative part that is not captured by the perturbative expansion in topology: Here, Z χ (β) stands for the contribution to the partition function from surfaces of Euler characteristic χ. It comes with a factor e χS0 , as S 0 is a coupling which multiplies the Einstein-Hilbert action in the model. (Although χ=1 for the (leading) disc order quantities, the subscript 0 will be widely used at leading order henceforth. So the disc level partition function is Z 0 , spectral density is ρ 0 , etc.) The focus of this paper will be on characterizing the full partition function of the theory, including the full non-perturbative physics, by making explicit aspects of the double-scaled matrix model definitions suggested in refs. [3,4], which should be considered companion papers to this one. The beautiful work of refs. [5,6] in defining double-scaled matrix models of (various kinds of) JT gravity is intrinsically perturbative in spirit, since * johnson1@usc.edu they use recursion relations connecting different topologies, and the work in refs. [3,4] is intended as a complementary construction (using minimal strings) that allows more direct access to non-perturbative quantities. The output of this paper will be the first explicit computation of the full spectral densities (and hence the partition functions, by Laplace transform), and explorations of several important phenomena that depend crucially on being able to compute non-perturbative physics. An example of the latter is the 2-point "spectral form factor" shown in figure 1, a quantity that helps in diagnosing universal aspects of quantum chaotic be-arXiv:2006.10959v2 [hep-th] 8 Jul 2020 haviour [7,8]. It was computed using the methods of this paper. This is the first time this quantity (and others like it to be presented later) has been computed fully in JT gravity or supergravity for generic values of β and S 0 , and so some time will be spent unpacking the techniques and results 1 . The late time "plateau" feature of the curve, and the transition to it from the "ramp" behaviour, are intrinsically non-perturbative features of wide interest. There are important non-perturbative effects that show up in the slope part too, in some cases, as will be demonstrated. They can sometimes be dramatic, as will be seen in the supergravity examples presented.  Figure 2: The full spectral density, computed using the methods of this paper, for the (2,2) model of JT supergravity. The dashed blue line is the disc level result of equation (5). Here =1.
Another example of this paper's results is given in figure 2. It is the full spectral density (the thicker line, actually made out of data dots) of the model ρ SJT (E) with the classical result (see equation (5)) plotted as a dashed line for comparison. By Laplace transform, this function defines the full non-perturbative partition function for the supergravity theory, and is computed here explicitly for the first time. This JT supergravity model is in fact the (α, β)=(2, 2) matrix model in the Altland-Zirnbauer [10] classification scheme, or the case Γ= 1 2 in the notation of ref. [4]. The result for the companion (0, 2) case (Γ=− 1 2 ) will be displayed later (see figure 8 on page 8). As can be seen in figure 2, for the (2,2) case non-perturbative effects entirely erase the characteristic classical peak in the spectrum at low energy, which dramatically alters the "slope" part of the spectral form factor as compared to the analogous result for the (0,2) case where a peak persists in the full spectrum.
While ordinary JT gravity is important and interesting (and results will be presented for it), a good deal of attention will be given to these two particular models of JT supergravity. They are of particular interest because the non-perturbative physics is more dramatic, in a sense. It was observed in ref. [6] (and confirmed to be manifest in the minimal model construction of ref. [4]) that beyond the first one or two leading orders of perturbation theory (depending upon the quantity being computed) the entire topological perturbative series vanishes. Therefore the non-perturbative effects uncovered in these models (as will be done here) are placed more in stark relief than other JT gravity systems.
Having shown examples of the key results, the job of the rest of the paper is to explain how to get them, and then to interpret them. The results follow from the nonperturbative construction, proposed in refs. [3,4], of JT gravity and supergravity in terms of minimal string models (of a special type). The basic idea, building on suggestions in refs. [5,11], is to reinterpret the JT system as an infinite set of minimal models (non-linearly) coupled together in a particular way, or equivalently (as explained in ref. [4]) by turning on an infinite set of operators in the minimal string model obtained by taking the k→∞ limit 2 . Since the full information about the kth minimal string model in question (see section III for a quick review of the essentials) involves solving an order 2k+1 highly non-linear ordinary differential equation (ODE), this way of defining JT gravity or supergravity involves solving an infinite order differential equation. This might seem rather daunting, or even formal, but from a pragmatic point of view it is rather straightforward to implement an approximation scheme that allows computation of an answer to a specific concrete question, to whatever accuracy is desired. The point is that the contribution to the model of successively higher orders of derivatives in the ODE grows smaller with increasing k, and so there is a point at which truncating the ODE and solving a finite order equation will give access to the full spectrum all the way up to a given desired energy, to some required accuracy. In other words, this is hardly any different from computing Feynman diagrams up to some sufficiently high order for some field theory problem (except that here the formalism is computing nonperturbative physics, and moreover the series is convergent, not asymptotic.) An outline of the paper is as follows: Section II is a brief summary of some of the (now standard) key ideas in the study of JT gravity that will be used in this paper. It is entirely optional for those who know the subject well, but serves to set context, notation, and (perhaps) some motivation. The deconstruction in terms of minimal models will be lightly explained in Section III.
Refs. [3,4] should be consulted for further details, and the non-perturbative explorations of key toy models presented there. The main task of this paper is to show how to extract non-perturbative results for the full JT (super)gravities. In particular, this section will explain how (using the supergravity examples) the truncation scheme of the previous paragraph works. Section IV will solve the full quantum mechanical system to yield the nonperturbative spectral density (and hence the partition function), for the supergravity cases. Then Section V turns to the non-perturbative spectral form factor for the supergravities, explaining how it is computed and then displaying several results.
Section VI then discusses the analogous construction and results for a non-perturbative definition of ordinary JT gravity obtained (as first presented in ref. [3]) by embedding it into a larger framework that it matches perturbatively (at high energy) but which supplies it with nonperturbatively well-behaved low energy physics. Section VII contains a proposal and brief preliminary discussion of an alternative non-perturbative definition of ordinary JT gravity that may be more natural than that in the previous section. For it to work, a particular kind of solution to the string equations should exist that has seems to have not been considered in this context before, and sample truncated solutions are displayed.
Since most of the results of this paper come from numerically unpacking the highly non-linear system of equations (and also using computer algebra to help unpack them), some Appendices are included with some (it is hoped) helpful technical notes and suggestions about the methods employed, for the reader interested in computing these or other results using this formalism. Appendix A presents a numerical study of the spectral form factor of the Airy model (the double-scaled Gaussian Hermitian matrix model) and compares the results to the known exact expressions, showing how the effects of the truncation to a numerical system are extremely well controlled. This serves as a demonstration of the trustworthiness of the numerical results obtained for the JT gravity and supergravity models in the main body of the paper. Appendix B 1 describes aspects of solving high order differential equations numerically, and Appendix B 2 describes how to solve for the energies and eigenfunctions needed to build the spectrum and spectral form factor. Appendix C lists some important quantities needed in the body of the paper (the Gel'fand-Dikii differential polynomials) and a recursion relation for getting the higher order expressions.
There are some brief closing remarks in the final section, VIII, with thoughts about the potential application of these methods to other systems.

II. JT GRAVITY LIGHTNING TOUR
Although it is a 2D theory of quantum gravity, by virtue of a coupling to a scalar, the dynamics of JT gravity is all on the 1D spacetime boundary. (A good review of much of this is ref. [14].) The boundary can change its shape while keeping its total length fixed to be the inverse temperature β=1/T , the period of Euclidean time. Meanwhile, the bulk spacetime has constant negative curvature (the Ricci scalar R=−2). So the theory is locally AdS 2 , and the leading spacetime (disc topology, i.e., no handles or crosscaps, one boundary) is often called "nearly-AdS 2 " [15][16][17][18], in the sense that, e.g. in Poincaré coordinates, the boundary is not a fixed circle an infinite distance away, but instead a finite loop of length β that is allowed to change its shape. See figure 3. At this order the dynamics of the loop is controlled by a Schwarzian action [16], and the result is: related to the disc order spectral density ρ JT 0 (E) by a Laplace transform. There is a JT supergravity generalization of this result [6,19]: defining a disc order spectral density ρ SJT 0 (E). In each case, the densities are given by: (Henceforth the redefinition √ 2ρ SJT 0 →ρ SJT 0 will be done, to adapt JT conventions of ref. [6] to the matrix model normalization to be used here.) The coupling e −S0 will be denoted in what follows, and indeed it will be the of a key quantum-mechanical system to appear shortly. One interpretation of S 0 is that it is simply the leading (T =0, disc topology) contribution to the entropy. For the ordinary JT case: This leads to a second reason (beyond the one mentioned in the introduction) to study JT gravity. It is a model of the low-temperature (near-extremal) dynamics of certain higher dimensional black holes and branes (see e.g. refs. [20][21][22][23][24]). For example, the metric of a charged black hole in d=4 is well known to become AdS 2 ×S 2 at T =0, and the area A of the two-sphere, S 2 sets the T =0 entropy: A=4S 0 . Turning on a small temperature replaces AdS 2 by "nearly-AdS 2 ", and the horizon area and hence the entropy gets corrections. The JT gravity model captures the dynamics of these corrections. (The dynamical scalar represents the deviation of the area away from extremality.) The 2D dynamics can be thought of as containing black holes in its own right as well, worth studying in their own terms. These are, at leading order, the disc geometries already described. A third reason for studying JT gravity is that it is a low energy holographic dual, in a certain sense [15,17,18,25] of a class of 1D quantum systems that exhibit quantum chaos, such as the Sachdev-Ye-Kitaev (SYK) model [16,26,27]. A key diagnostic of the quantum chaotic behaviour of the system is the 2-point "spectral form factor" Z(β−it)Z(β+it) , which exhibits certain key universal features [7,8,28,29]. Starting out at Z(β) 2 , it decays down a "slope" to a "dip" at during the first epoch of time t, rises along a "ramp" at intermediate times, before levelling off to a "plateau" at late times at a value given by Z(2β) . (See all these features in figure 1, but recall that it is not an SYK spectral form factor, but a gravity one; See below).
The timescales over which these features manifest are important, especially the time to when the plateau sets in, as it gives a measure of how long correlations take to wash away. No single SYK dual cleanly exhibits the universal behaviour individually. There are wild oscillations in the spectral form factor at intermediate and late times 3 . Instead, these features emerge as the timeaveraged behaviour, as can be seen by averaging over an ensemble of models [30]. An important idea in quantum chaos is the notion that random matrix ensembles should capture the universal features seen in the averaged behaviour of a chaotic system (for a review see ref. [7]). This led to the suggestion of refs. [30,31] that a random matrix description of averaged SYK could be available. On the other hand, random matrix models are known to describe, in a "double-scaling" limit [32][33][34][35], the sum over surfaces of a 2D quantum gravity, so this is another way of seeing that there ought to be a dual gravitational description of SYK-like models. This was shown to be more than a coincidence of ideas in ref. [5], where JT gravity was demonstrated to be explicitly consistent with -order by order in the topological expansion-the properties of a double-scaled matrix model. Ref. [6] furnished several more examples and a classification of the 3 In the phraseology of the moment, these later eras are "difficult times" for an SYK model. possibilities in terms of the ten standard random matrix ensembles. So the JT gravity dual (or supergravity dual, for the appropriate generalization of SYK [6,19,[36][37][38][39][40][41]) performs the ensemble average directly. The early time behaviour is controlled by the disconnected diagram constructed of two discs (a pair of AdS 2 black holes), plus corrections, while the later ramp and plateau features come from the cylinder diagram (an AdS 2 wormhole) [42] plus corrections. See figure 4. These amplitudes do not fluctuate chaotically in time, but have smooth behaviour to be expected from geometric objects in a theory of gravity. This can be seen already in the leading computation for the cylinder diagram [5,43], which yields a simple t dependence: Z(β−it)Z(β+it) ∼β −1 β 2 + t 2 gives a rise for the initial part of the ramp behaviour in a regime that would already be beset by fluctuations in any given SYK model. For t β, assuming the transition to the plateau has not yet occurred, this yields a linear rise. In this paper it will be observed that nonperturbative effects can, depending upon the value of , take over rapidly to generate the ramp, and so the linear part is hardly visible at moderate β.
As already mentioned, the plateau in the spectral form factor (and the transition to it from the ramp) is a result of perturbative and, especially, non-perturbative corrections to the leading cylinder contribution. The purpose of this paper is to focus on unpacking the non-perturbative definitions of refs. [3,4] in order to explicitly uncover such effects. Figure 1, already shown above, is a sample of the work reported on in this paper. It is the full spectral form factor for a particular model of JT supergravity. It will be discussed more fully in section V. Now, on to the computations.

III. CONSTRUCTING JT (SUPER)GRAVITY FROM MINIMAL STRINGS
The key ingredients are certain double scaled matrix models that have been used in the past to study certain kinds of "minimal" string theories. (See e.g. refs. [43,44] for reviews.) The details of the string theory constructions do not matter here. The most important fact to know is that some of the models (a subset of the "onecut" matrix models) can be described in terms of an associated 1D quantum mechanics problem [45,46], with Hamiltonian: where the potential u(x) satisfies a non-linear ordinary differential equation (ODE) called a "string equation".
The key task is to build the correct u(x) for the problem in hand. Once it is known, the full spectral density can be extracted by simply solving the spectrum of H and evaluating the fully non-perturbative ρ(E), using an expression given in the next section. It is useful to note that in the limit where just the disc-level physics is kept, the spectral density at this order can be written as a simple integral involving the leading part of the potential, u 0 (x)= lim →0 u(x): where f (u 0 )=−∂x/∂u 0 (x). Turning back to the ingredients, the minimal models will be labelled by an integer index, k. As mentioned before, the models will be combined together to yield the JT (super)gravity. There is a parameter, t k , that will be used to turn on the kth model in the mix. The model is turned on if t k is non-zero. The minimal models in question can be obtained [47][48][49][50][51] 4 by taking the doublescaling limit of models of a complex N ×N matrix M , with a potential V (M † M ) (see also footnote 5). The string equation that needs to be solved is: where the constant Γ will be discussed shortly and Here,R k [u] is the kth order polynomial in u(x) and its x-derivatives defined by Gel'fand and Dikii [53]. They have a purely polynomial in u(x) piece, which is u(x) k , a purely derivative linear piece, u(x) x-differentiated 2k−2 times, and then non-linear mixed terms involving u(x) and its x-derivatives. Here, they are normalized so that the coefficient of u k is unity. The first three are: where a prime denotes an x-derivative times a factor of . It will transpire thatR 4 ,R 5 ,R 6 andR 7 will be used in this paper too, but since they are rather lengthy, some are listed in Appendix C, along with methods for generating others if needed. The boundary condition that ensures good nonperturbative behaviour is, for each model, Note the presence of =e −S0 =1 in the string equation (and the various quantities that make it up). It is very useful for separating the classical parts from the rest, by sending →0, or equivalently, dropping derivatives. For the study of non-perturbative physics, solutions u(x) of the equation will be extracted for =1. Several results presented in figures to come will be for this value, because it allows the non-perturbative effects to be writ large in the results (for spectra, etc.), and therefore seen easily. When instructive to do so, comparison to results with dialled down will be discussed. It is interesting that it is in fact more difficult to solve the string equations for smaller . This is because when derivatives have smaller coefficients, they are allowed to fluctuate more, contributing to the sensitivity when solving these highly non-linear equations, as will be discussed later. (This increased difficulty to get smaller results has the character of a sort of strong/weak coupling duality, in fact.) Turn now to the constant Γ in the string equation (9). With it present, the matrix model is in the (2Γ+1, 2) class in the (α, β) Altland-Zirnbauer classification of matrix ensembles 5 . The two choices Γ=± 1 2 will mostly be considered in this paper, and the two JT supergravity models discussed here are will be labelled (2, 2) and (0, 2).
A particular (say, the mth) minimal model can be studied by setting all the t k to zero except for k=m. As discussed in the previous two papers [3,4], the m=1 case in particular was important as it models the low energy tail of the eigenvalue distribution very well. To get the full behaviour, all of the t k must be turned on in a particular combination. For example, in the case of JT supergravity, the combination (derived in ref. [4]) is: 5 Double scaling means that in the matrix model [47][48][49][50][51] of the complex matrix M , the size N is taken to infinity while couplings in the potential V (M † M ) are tuned to certain critical [54] values. Diagonalizing M turns this into a problem involving its eigenvalues λ i (i=1 · · · N ) at a cost of introducing a van der Monde determinant J= i<j (λ 2 i − λ 2 j ) 2 for the Jacobian. The effective Dyson gas problem for the λ 2 i can be thought of as existing on the positive real line, with a wall at zero. The constant Γ in equation (9) can be thought of as arising from the coefficient of a logarithmic term in the potential of the model (see e.g. ref. [55]), and as such, results in an extra factor λ 2Γ i to the effective integration measure over the ith eigenvalue, giving i λ 2Γ+1 i dλ i . With the factor J included, the model is seen to be in the (2Γ+1, 2) class in the (α, β) Altland-Zirnbauer classification of matrix ensembles, defined for α=0, 1, 2. Actually Γ can be more general integers or half-integers that just these values. So all the infinite models are turned on and the string equation becomes a highly complicated object. But the purpose of this paper is to show that physics can be readily extracted nonetheless.
Here is the reason why. First, look at the disc level. The string equation is expression (9) with the three parts with 2 in front of them removed, and so the solutions are either u 0 =0 or corresponding to the asymptotics given in equation (12). For the second choice, the x<0 regime, the combina- (8), yielding the part of the spectral density expanded in positive powers of E. The simple E − 1 2 part comes from the u 0 =0 behaviour. Integrating f (u 0 ) with respect to u 0 , or simply by looking at equation (14), the explicit potential that gives JT supergravity on the disc is given by the equation where I 0 (s) is the zeroth modified Bessel function of s. This is a remarkably simple form. The issue of tractability becomes the simple issue of how well this potential can be approximated by truncating to a finite number of t k s. The answer boils down to what maximum energy scale E one wants to know the spectrum up to, and to what accuracy. As an example, the full classical potential (15) is plotted in figure 5 alongside two truncations. The first truncation contains just t 1 and t 2 : and it is clear that it is a good approximation for energies up to approximately E∼0.1 after which it begins to deviate considerably. The next example truncation adds t 3 and t 4 : and for energies up to order E∼0.5 it serves rather well. Further improvements come by adding higher orders. The next issue to appreciate is how much the solution changes when all the non-perturbative corrections are included. For all k a shallow well can develop in the central region (slightly to the right of x=0). Crucially, moving away from that region, the deviation of the solution from the disc level behaviour rapidly dwindles, as it matches on to the asymptotic behaviour. The same is true for the coupled solution. Moreover as can be seen from equation (13), the form of the t k as k grows is such that good approximations at the disc level can be found by adding only a small number of minimal models, for a given needed accuracy. For larger Es, the solution becomes hard to distinguish from the classical result, and in that case the exact classical potential can be used, to a good approximation.
The highest truncation levels chosen for the purposes of this paper was to keep all the minimal models up to k=6 (this section) and k=7 (for section VI), although very good results were obtained for lower order truncations too. Since R 7 [k] has the twelfth order derivative of u(x) in it, the string equation (9) is a 14th order differential equation. (It is 12th order for the k=6 truncation.) In general, it is easier to take a derivative of the string equation, whereupon an overall factor of R can be divided out, reducing some of the non-linearity somewhat, at the expense of an increase in the order. For the boundary conditions in question, a 15th (or 13th) order differential equation is not too hard to solve numerically, with care. Some suggestions and notes are given in appendix B 1, for those who wish to carry out their own computations using this framework. Part of the full nonperturbative potential u(x) for the truncation to k=6 is displayed in figure 6. This is for the case Γ=+ 1 2 , i.e., the (2,2) JT supergravity. (It was solved between x=−200 and x=+200.) Notice that it approaches the classical solution and agrees with it rather well up (and beyond, it turns out) to E∼1.3, and so the full solution, out to beyond the x −100 shown, can be used to capture the spectrum of (2,2) JT supergravity with good accuracy (see more in appendix B 2 on how to do this). It is (relatively) easy to do better, if desired, but little visible change was noticed in going to higher orders, in exchange for accessing only a slightly larger maximum energy. Figure 7 shows the solution for Γ=− 1 2 , which will be used to study the properties of the (0,2) JT supergravity.
This solution has a well in the interior (as is quite typical of these solutions), and so slightly more rapid changes  take place there. Rather than use numerical methods to solve for this directly (which are inevitably more sensitive to error in this case), a handy solution-generating technique derived in ref. [56] was used 6 , that allowed it to be generated from the Γ=+ 1 2 solution already found. See Appendix B 1 for more on this.

A. General Remarks
The next step is to solve the full spectral problem for the Hamiltonian H, given the potential u(x) found in the previous section by solving the (truncated) string equation. The relation between the spectrum of H and the JT partition function is as follows. From the minimal string perspective, the JT (super)gravity partition function is simply [5,11] the expectation value of a "macroscopic loop" of length β. The technology for working this out was derived long ago in ref. [45]. (Ref. [43] unpacks the formalism in a useful review.) It is the trace of the exponentiated H, with a projection P inserted: where the operator P≡ µ −∞ dx |x x|, and the upper limit µ will be discussed shortly. Inserting a complete set of states: yields where is the spectral density.
To understand the upper limit µ, return to the leading order expression for the density given in equation (8), and change variables to x. The solid black line in figure 5 is a reminder of the behaviour of u 0 (x). The part of the integral where u 0 (x)=0 begins at x=0 and extends to positive x to some value denoted µ, while (at a given E) the lowest value for x is set by the turning point where E=u(x), at a position denoted −|x 0 | and so: In the full quantum expression (21), there are contributions to the spectral density from the whole integration range down to x=−∞. Intuitively, this is because wavefunctions penetrate to the left beyond the classical turning point E=u(x). (Often, in the classical expression (22), the lower limit is also written as x=−∞ with the understanding that only the real part of the integral contributes to ρ 0 ) 7 . The role of µ becomes clear from focussing on just the contribution from the u 0 (x)=0 part in equation 22. It adds µ/(π √ E) to the classical spectral density. Physically, µ controls how much the classical spectral density is "piled up" against the natural wall at E=0 that stops the energies in the complex matrix models from going negative (see footnote 5). Vanishing µ corresponds to just touching the wall. In the normalization of equation (5) (without the √ 2), the value of µ is 1, and this will be used for much of the rest of this paper. However, it is useful for later to write a more general classical spectral density with a different µ:

B. Computation
Now to the matter of computing the full spectral density (21). Just as in refs. [3,4], a matrix Numerov method [57] was used to solve for the spectrum of H. Conceptually it is a simple problem in finding eigenvalues and eigenfunctions for a Schrödinger operator in 1D. The wavefunctions are free and oscillatory to the far right (x>0), begin to feel the presence of the potential as they move further into the interior (it starts as a 1/x 2 dependence), and then once they hit the potential there is an exponential decay to the left (x<0). As a guide to extracting them with good accuracy, some extra suggestions and notes for interested readers are given in appendix B 2. The same normalization method as the one used in ref. [3] was used for the resulting eigenfunctions.
The key point explained there is that in the far x>0 region, wavefunctions are known to asymptote to a simple form involving the Bessel function of order Γ, where the normalization can be analytically chosen to yield the correct contribution to the disc level spectral density.
The outcome of the numerical spectrum solving was approximately 1000 accurate normalized wavefunctions and their energies, for use in constructing the density in this section, and the spectral form factor in the next. Using them, the spectral density can be constructed using a simple trapezoidal integration to implement equation (21) and the result (for µ=1) is shown in figure 2 for 7 Of course, the integrand of (21) is not written in terms of E and the full u(x), but in terms of the wavefunctions. The precise connection between the two expressions is via the diagonal of the resolvent, (H−E) −1 , of H, denoted R(x). It can be built out of ψ(x, E)ψ * (x, E), using standard Green function methods. Gel'fand and Dikii [53] derived an equation . The leading order solution comes again from dropping derivatives, and is R(x)=±(2 u(x)−E) −1 . In the normalization of this paper ρ 0 =(2/π )Im R(x)dx, yielding the form (22).  √ E divergence present at disc level is erased entirely by non-perturbative effects in the Γ= 1 2 (2,2) case, but not in the Γ=− 1 2 (0,2) case. As discussed in refs. [4,6], for Γ=± 1 2 , there are no perturbative corrections to the spectral density beyond the disc, so all the differences here are due to non-perturbative physics, which makes these JT supergravity cases particularly interesting to study when investigating the results of non-perturbative physics. This stark difference will be reflected in comparisons of the spectral form factor, to be studied in the next section.
Satisfyingly, what has appeared here for the full (2,2) and (0,2) cases are actually grown-up versions of what was discussed for two baby (Bessel) models in ref. [6], and so this constitutes a nice consistency check of the methods of this paper. The Bessel special cases are models of the very low energy tip of the spectral density-the part that (classically) pushes up against the E=0 wall, as mentioned earlier. Generalizing to include arbitrary positive µ (see also refs. [3,4]), the full non-perturbative Bessel spectral density is (after a change of conventions from those papers): where (as a reminder) α=2Γ+1 takes the values 0 and 2. The leading part is the disc contribution and all other perturbative contributions vanish. The oscillating term is the full non-perturbative contribution. At E=0 for α=2, (Γ= 1 2 ) there is the aforementioned cancellation between the divergent disc contribution and the E→0 piece of the oscillating non-perturbative part. This is what is now seen to be present (as it ought to be) in the full JT supergravity spectral density computed explicitly above.

C. A Special Formula, and a Tale of Instantons
Actually, it is possible to go considerably further than the Bessel tail comparison. The structure of the Bessel models above (along with various other technical features of the matrix model recursion relations) led Stanford and Witten to suggest (see appendix E of ref. [6]) that more generally, the JT supergravity spectral density for the (0,2) and (2,2) cases should closely approximate the form 8 : . (25) (Note that their conventions were adapted to match those of the current paper.) Again, the perturbative terms beyond the disc vanish. The argument of the sinusoidal non-perturbative piece is determined in terms of the disc contribution 9 . Crucially, they note that this does not include possible instanton corrections, non-perturbative physics of a form not accounted for by the sinusoidal modulation (hence the use of the ≈ sign in the above). While the expression above was written for the specific value µ=1 (µ is not a variable parameter in ref. [6]), other values of µ (see below equation (21) for its explanation) are quite readily incorporated by this form. This is already obvious for the simple Bessel prototype (24), but further evidence comes from inserting into it the expression (23) for ρ 0 (E, µ) written earlier for the disc level JT supergravity density for other µ. Plotting this analytical result against the numerical results for the spectral density obtained using this paper's methods yields a remarkable agreement, as shown in figures 9 and 10 for Γ= 1 2 and Γ=− 1 2 , for the values µ=1, 5, and 10. As a reminder, the red lines are each made of ∼1800 points individually computed using the truncation and numerical methods outlined, while the blue dots are samples of the analytical formula (25). It is very striking how well they agree.
The areas of disagreement suggest interesting physics, in fact. For smaller µ (e.g., see the case of µ=1 in each figure) it is easy to see disagreement. It is especially quite visible at smaller energies (although at E=0 they agree as they must, since this is covered by the Bessel cases above). The clear pattern in the discrepancy is suggestive of physics and not numerical inaccuracies. 8 The author thanks Douglas Stanford for asking a helpful question about this issue after the 1 st version of this manuscript appeared. 9 In fact, it is easy to guess generalizations of this formula for other half-integer Γ, which seem to give sensible physics using the present methods. It would be interesting to test them.
For example, the fact that the disagreement rapidly disappears as µ increases is striking. Nothing about the numerical methods should show this sort of systematic dependence. As a further example, the form of the disagreement suggests that there is a change in the amplitude of the sinusoidal modulation, again a pattern that is hard to reproduce with mere numerical systematics. A natural guess is that there are further multiplicative factors in the non-perturbative parts that are not present in the analytic formula. They must be small, and of an instanton form (as already anticipated in remarks in ref. [6]). In the kth minimal model the height of the effective potential barrier for one eigenvalue in the background of all the others yields the action of the instanton effects [58] (they are "ZZ-branes" in the modern parlance [59]). Alternatively, it can be computed using a WKB analysis of the string equation itself. This was all worked out long ago for the complex matrix models in ref. [50] (see ref. [60] for a review and more contemporary presentation), and the dependence of the action on µ and is ( 4k 2k+1 )µ 1+ 1 2k / . Therefore a dependence exp(−2µ/ ) should be expected for instanton corrections to JT supergravity constructed out of minimal models in the manner of this paper.
To test this, a rough estimate of the dependence was done by measuring the size of the deviations of the type seen in the figures (9 and 10) for some sample cases, as µ varied from 0 to 5. The deviations fall rapidly in that range (as is evident from the figures), and a logarithmic plot of them against µ yielded a straight line of negative unit slope to good (better than 1%) accuracy. While this can be made more precise, it is already a strong indication that the instanton effects of the expected form are present and accounted for. It will be interesting to explore this further, but this will be left for later work.

D. Going to Weak Coupling
Before moving on it is worth showing the effect of going to weak coupling by reducing . The spectra shown so far are for =1, and working in this regime has the advantage of making much of the non-perturbative physics readily visible. Nevertheless, it is useful to be able to gain access to weaker coupling, not the least because certain dual systems of interest (such as SYK-type models, or large black holes (where S 0 is large)) may be in that regime.
It is possible to reduce the value of , but it is at the expense of not being able to access as high energies in the spectrum for a given truncation, unless there is a compensation in terms of increased computational effort (working on a finer grid and producing more wavefunctions). Additionally, solving the string equation takes more care, as already mentioned, because the derivatives are more in play. Nevertheless, progress is possible. For example, dialling down to e.g., = 1 5 , the string equations can be coaxed to produce the solutions shown (just as close-ups, for brevity) in figure 11. These should be compared to what was shown in the insets of figures 6 and 7. For the (2,2) case, the curve turns the corner more sharply (showing the increased role of the derivatives), while for the (0,2) case, the well is far deeper and narrower (these two features come together to continue to ensure there are no bound states).
It is easy to imagine how these developments continue in order to get to the classical limit →0 (the red dotted line in the figures). The (2,2) curve turns the corner ever more sharply, and for the (0,2) curve the well deepens and narrows and eventually is entirely squeezed away. The techniques already described in the previous two sections can be carried out to study the spectral density, and the result in figure 12 displays the results (red dots) for the Γ= 1 2 or (2,2) case, for the cases µ=1 and 5. (The µ=10 case, which rises to 508 vertically, was omitted for clarity.) Again, the blue dots are samples of the analytic formula (25) with (23) input for ρ 0 (E, µ) (now = 1 5 is inserted), showing again how well this procedure works. Notice that, in accord with the fact that the coupling is weaker, the exp(−µ/ ) instanton effects described in the previous section, which would appear as deviations of the red curves from the blue dots, are too small to seen here.  Figure 12: Two examples (µ=1, lower; µ=5, upper) of the full spectral density for the (2,2) model of JT supergravity, for reduced to 1 5 . The (blue) dots are points from the equation (25). In each case the dashed line is the disc result of equation (23). c.f. the =1 cases of figure 9.
Comparing to the =1 curves in figure 9, the quantum undulations are smaller, for generic E, and the curves stray less from the classical curve, beginning marked deviations only at lower energies compared to the =1 situation, until the curve falls quickly to zero. This is as it should be.

A. General Remarks
As mentioned in the brief review of section II, the spectral form factor is derived from the two point function of Z(β). This has two parts, a disconnected piece Z(β) Z(β ) and a connected piece Z(β)Z(β ) . In the old matrix model language, the connected piece is the connected correlator of two "macroscopic loops", and this is readily written down as [43,45]: = Tr(e −(β+β )H ) − Tr(e −βH Pe −β H P) where P is discussed below equation (18) and The quantity µ, discussed in the previous section, will be set to unity henceforth, corresponding to having the disc level spectral density given in equation (5) (without the overall √ 2).

B. A Phase transtition
Set β=β for a while. Generically, the disconnected part (corresponding to two black holes) and the connected part (a wormhole) are worth studying in their own right as distinct sectors of the quantum gravity that compete for dominance [61] as a function of β. The disconnected part, being the square of the partition function, rapidly decreases with increasing β while the connected part increases. At some point β cr there is a transition, and the connected diagram becomes more dominant. This is also true in the complete (not just perturbative) theory discussed here. This is all nicely under control in the current definitions of JT supergravity. The spectrum has been computed in the previous section and so all the elements in equation (26) are readily computable to the desired accuracy. Figure 13 shows a plot of the (log 10 of the) connected and disconnected pieceswith all perturbative and non-perturbative contributions included-as a function of β, for the (2,2) JT supergravity case, showing the transition at β cr =13.56. A similar computation (yielding a similar graph, omitted) shows that β cr =28.30 for the (0,2) JT supergravity case.
Notably, for the (0, 2) case the amplitudes are over an order of magnitude larger. This is a striking effect attributable entirely to non-perturbative effects. The (2,2) JT supergravity has, as mentioned in the previous section IV, non-perturbative effects that cancel the 1/ √ E behaviour at low energy, coming from the leading disc amplitude. The (0,2) version does not cancel this away, and so while to larger E the two models' spectra are roughly similar (see figures 2 and 8), there is an enhancement at low E that means much larger contributions to the partition function for any fixed β. This marked difference between the non-perturbative physics will make a major appearance in the temporal behaviour of the spectral form factor too, studied next.

C. Time Dependence
The above computation served as a useful guide for what to expect for the spectral form factor, which tracks the correlation function over time. This is done by putting β→β+it and β →β−it, with β fixed, studying the dependence on t. The fixed β can be above or below β cr . The disconnected part will start out as the squared partition function and then decrease with t. This is the "slope" behaviour of the spectral form factor. On the other hand, looking at equation (26), it can be seen that the connected part has a t-independent part, Z(2β), from which is subtracted a positive piece which gets small at large t, with only significant contributions from energies that are close to each other. The value of Z(2β) therefore sets the height of the universal "plateau" feature and the approach to it, the "ramp" has its size set by how rapidly the energy correlations die away at large t. The "dip" region is formed by the process of handing over from the decreasing disconnected part to the increasing connected part. The formalism here allows, using the complete package of almost 1000 good wavefunctions and energies, for this all to be computed to good accuracy for JT supergravity (and for a non-perturbatively well-behaved definition of JT gravity in section VI) for the first time. Using again the JT supergravity example, figure 14 shows the disconnected contribution to the spectral density function, for β=50>β cr , with log 10 axes. The slope behaviour is quite evident. Since the axes are logarithmic, it is easy to see from the figure that the slope of it is roughly −3, suggesting a fall-off of ∼t −3 . This might seem surprising since it is similar to the perturbative fall off rate for ordinary JT gravity (see ref. [30]). The reason for this faster rate is clear, and again attributable to non-perturbative physics. The slope's falloff is controlled by the behaviour of the endpoint of the spectral density. While at the disc level for JT supergravity ρ(E)∼1/ √ E there, producing a t −1 fall-off [62], the fact that in the (2,2) case non-perturbative corrections remove this 1/ √ E behaviour results in the faster fall off more usually associated with the ordinary JT case (and Hermitian matrix models). There are far fewer states in the vicinity of the endpoint. This reasoning predicts that for the (0,2) supergravity, the slope should be closer to a t −1 fall-off. Indeed, this is clear from the behaviour of the disconnected piece for (0,2) supergravity, shown in figure 15, for β=50. There, the slope of the linear part shows t −1 fall off. Continuing in line with these expectations is a computation of the same quantity for the =1/5 case mentioned earlier. The spectrum was shown in figure 12, and from there it is natural to guess that the fall-off would be even faster since the non-perturbative effects have scooped away even more states near the endpoint. A check showed that this is indeed correct, although another figure will not be presented to display the result, to save repetition.
Moving to the connected contribution's time dependence, the (2,2) case is shown in figure 16 for the same value of β=50. The ramp and plateau structures, and the transition between them, are visible. Strikingly, the rise to the plateau is very short-lived, the ramp regime rising only a small (on the logarithmic scale) amount before  transitioning to the plateau.
The beginning shape of the ramp is already anticipated in the perturbative answer, ∼β −1 β 2 + t 2 , long known for two-macroscopic-loop correlators [43] (appropriately continued to yield the t-dependence [5]), but there are strong non-perturbative corrections such that before the long-time linear part can manifest, the other effects turn the ramp over into the plateau.
A similar story is told, initially, by the ramp shape seen for the (0,2) case, shown in figure 17 (again for β=50). However, there's a new feature. The rise is indeed slow, but it is remarkably slow. After almost two orders of magnitude more time has elapsed, as compared to the (2,2) case, the saturation to the plateau has still not quite completed. This is a rather novel feature of this case, and worth further investigation. The origin of this physics is most likely again to be attributed to the peculiar pileup of states in the vicinity of E=0 that this model has. There's an endless supply of closely spaced low-lying states contributing to the part of the form factor that subtracts from the saturation value Z(2β) (see equation (26)). At longer and longer times there are even more low lying states to contribute, and still closely spaced, maintaining their effect of slowing the saturation.
Turning back to the ramp itself, note that this was (for both (2,2) and (0,2)) for the case =1. At smaller values of there is more time for the ramp to develop, with an increased rate of rise before the turnover. Note however that it is only for extremely small that the linear behaviour of the ramp has a chance to appear.
There are two important lessons here. The first is that associating the ramp with linear behaviour (as is sometimes done in the literature) is maybe not the most accurate descriptor. The second is that non-perturbative effects can enhance the appearance of the ramp in the JT supergravity case (in the full spectral form factor made by taking the sum of disconnected and connected parts, even though perturbative expectations might have suggested a reduction [62]. The potential reduction of the ramp feature is based on the idea that a slow t −1 rate of fall might not give the ramp time to develop before the plateau sets in. In fact, non-perturbative effects are seen here to produce a rapid fall (sometimes faster even than the perturbative bosonic t −3 ) giving plenty of time to develop a sharp "dip", a clear ramp, and a smart turnover into the plateau for the (2,2) case.
The sum of the connected and disconnected pieces gives, for (2,2) supergravity, the classic saxophone shape 10 known from studying spectral density functions in a wide range of contexts. It is displayed in figure 18. This is for =1, and the case of =1/5 was already pre- sented in figure 1, at β=35. The latter, being at smaller β and , is larger overall and develops a wider variation. The linear part of the ramp would be even more visible for lower values of these parameters.
For (0, 2) JT supergravity, figure 19 shows the resulting =1 full spectral form factor, again at β=50. It is, as to be anticipated, almost two orders of magnitude larger than for the (2,2) case, because of non-perturbative effects (already discussed). There is (on the logarithmic  scale) a ramp-to-plateau transition (although it is, from the discussion above, much slower than for (2,2)). Also visible is the slower slope-to-dip time (due to its slower decay rate), as already discussed.

D. Temperature Dependence
It is of interest to see how the spectral form factor evolves as a function of temperature. The results for a series of increasing temperatures, β=50, 46, 42, 38, 34 and 30 are presented for the disconnected part of the spectral density in figure 20. The highest temperature curve is at the top. Strikingly, the curves soon merge into each other and follow the fall-off already discussed, regardless of the starting temperature. In figure 21, there is a series of the full spectral density, for the same set of temperatures. Crucially, for comparison purposes, the curves are all uniformly scaled (on the vertical axis) to have the same initial height as the highest temperature (β=30) case, which is the lowermost curve. Therefore, in this figure, relative slopes should not be taken literally. This scaling allows for ready access to some of the more meaningful comparisons to be made, such as the relative size of the curves: Higher temperature (smaller β) gives a vertically larger curve: Higher temperature "shakes" the system up more, resulting in a wider amplitude of deviation from the initial value before it settles down. Interestingly though, the dip time increases slightly with higher temperature, although not dramatically. The rapid ramp time also changes very slowly with β.
Heading toward smaller (but still moderate) values of β (in the region of β β crit ) there are small modulations in  both the disconnected and connected parts of the form factors, accumulating (in the latter) near the cross-over from ramp to plateau. The combined result of these higher temperature structures is that there is a damped wobble as the ramp merges into the plateau. An example of the full spectral form factor showing this feature is given in figure 22 for the case β=14. Whether this is interesting physics or not is not clear 11 . The value of the 11 It is reminiscent of features of an interesting exact expression derived in ref. [63] in the context of an SYK model with source terms.
temperature at which this can be seen seems comfortably below the highest energy allowed by the truncation. Some other other fascinating structures become apparent in the very high temperature regime. How useful they are for the physics in question is debatable since this whole context (the Schwarzian, the connection to black holes, SYK, etc.,) is in a low energy limit. Moreover, high temperature also begins to go beyond the energies for which the truncation of the string equation remains reliable. However it is interesting to observe the features anyway, and could well be instructive for understanding JT models with a cutoff [64][65][66][67]. Looking at the disconnected part of the spectral form factor for the (2,2) supergravity case, a series of dips evolve, becoming more pronounced toward higher temperature (smaller β). See figure 23. They have a clear pattern and structure and are consistent with observations made for large N random matrix systems (even without double scaling). (They are also analytically obtainable in the exact Airy example reviewed in Appendix A.) In fact, the disconnected function is beginning to resemble the form J 1 (t) 2 /t 2 that has been derived analytically for the infinite temperature case. After taking the logarithm, the zeroes of the Bessel function J 1 become the dips in the logarithmic plot. It would be interesting to show that this (or a variation thereof) analytical form emerges in this JT supergravity context as well. Following the numerics to smaller β seems to confirm this (a zero in the connected function also appears, in some examples), although eventually numerical inaccuracies begin to overwhelm the results, presumably because the correct physics needs to include contributions from energies that lie beyond the cutoff on the spectrum up to which the truncated equations are valid.

VI. JT GRAVITY: NON-PERTURBATIVE DEFINITION A
This section presents results analogous to those shown in earlier sections, for the non-perturbative completion of ordinary JT gravity presented in ref. [3]. It might seem odd to have studied the JT supergravity examples first, leaving this case for last, but there is good reason. The non-perturbative physics of this case is more subtle. JT gravity was shown, in ref. [5], to be perturbatively (in the topological expansion) equivalent to a double-scaled Hermitian matrix model, i.e. classified in the Gaussian Unitary ensemble (β=2 in the Dyson-Wigner series). On general grounds, such double-scaled Hermitian matrix models are known to sometimes have non-perturbative (in topology) instabilities, and so it is possible that the JT gravity definition inherits them. More specifically, thinking about the model in terms of constituent minimal models, as in section III, it is made up of an interpolating family of minimal models that have the x→−∞ boundary condition in equation (12) in both directions, and instead solve the string equation R=0. (Recall that R is given in equation (10)). These are the (2k−1, 2) bosonic minimal string models. For k even, these models are non-perturbatively unstable, as has been known for some time [45,[68][69][70]. From the point of view of the spectrum, all the models, when non-perturbative effects are taken into account, have contributions from arbitrarily negative energy sectors. Even though exponentially suppressed, for even k the effective potential turns downward for states at these energies, signalling the system's wish to tunnel to an entirely new solution that is quite different from the one around which perturbation theory was developed. From the perspective of this paper (solving string equations non-perturbatively), this means that for each of those (k even) models, there simply are no real smooth solutions of the kth equation with those conditions [45,69,71,72]. Since JT gravity is made up of, in equal measure, even and odd k models, this strongly suggests that it inherits these problems, as already noted in ref. [5]. (However, see section VII for a possible evasion.) The route that ref. [3] took to supply a nonperturbative definition of JT gravity was to embed it into a larger problem. The minimal models used in previous sections for JT supergravity, also indexed by k, have the same x→−∞ boundary condition as the bosonic minimal models, and in fact when Γ=0 they have identical perturbation theory as solutions to the differential equation (9). Put differently they solve R=0 perturbatively at large −x. This means that if used to construct a JT gravity model, they will yield the same physics at high energies E, but yield different physics as lower energies that is untroubled by the stability issues. The combination of models needed (that will give the Schwarzian spectral density (4)) at high E is as follows: (This relation was first derived in ref. [11], but with different normalization.) In fact it is possible to integrate the f (u 0 ) that results from this combination to find the explicit classical potential u 0 (x) that yields the Schwarzian density, through: the analogue of the case (15) for the JT supergravity found earlier. This non-perturbative taming of JT gravity can be explored extensively along the same lines as done for the JT supergravity models. The truncation scheme works along the same lines, and so there is no need to retread the ideas again. A solution to the string equation (9) with the bosonic JT combination of minimal models (28) was found with the first 7 models turned on, constituting a very good truncation where energies up to E∼1.5 can be trusted. A combination of numerical and analytic methods were used to find the solution to the 15th order differential equation. See appendix B 1 for some tips on how this was done. Figure 24 shows the solution, with the classical (disc level) potential that gives the JT spectral density (4) displayed as a dashed line. Again there is a small well (not deep enough to support bound states) in the interior, and then u(x) settles to zero to the right. These are features shared by the JT supergravity models, as should be expected since their components are being used here as a non-perturbative low energy "regulator" in a sense, removing the leakage to negative energy. Remarkably, the non-perturbative effects will (as happened for the (2,2) JT supergravity case) remove all signs of the classical 1/ √ E contribution to the low energy that would make for a striking and unsatisfactory departure from the physics of JT gravity. The methods of section IV then allowed for the spectrum to be computed (here, µ=0 is used for the upper limit in equation (21)), and figure 25 displays, for the first time, the spectral density of a non-perturbative completion of JT gravity with all perturbative and non-perturbative corrections included (up to this energy). It is clear from the figure that the non-perturbative ripples have already begun to die away and merge into the classical smooth region, showing that this truncation has captured the key physics that is affected by non-perturbative contributions. Notice that the spectrum naturally is bounded below by E=0, as already shown explicitly in preliminary studies in ref. [3]. Note that ref. [3] generalized the construction by turning on the parameter µ that was discussed in the supergravity context in section IV. This is straightforward to do, and there were no additional insights to be gained, and so results are not presented here. The tail of the resulting distribution, not surprisingly, resembles the tails already displayed in ref. [3]. Of course, with the spectrum in hand (approximately 1800 normalized wavefunctions and their energies) the next natural step is to compute the spectral form factor, using the methods of section V. The correlator of two boundaries is readily computed, and the phase transition where the disconnected part (two black holes) hands over to the connected part (wormhole) happens at β cr . (A figure similar to figure 13 is omitted here to avoid repetition.) For β=50, the disconnected, connected, and combined spectral form factor are shown in figures 26, 27, and 28, respectively.
The most striking feature overall is in the disconnected portion of the form factor, controlling the initial slope. As might be expected from the absence of the 1/ √ E low energy behaviour of the regulating models (happily removed by non-perturbative effects), the time dependence of the slope is not t −1 (as it is for the (0,2) JT supergravity, but nor is it the t −3 expected from the classical low energy physics to be read off from the β −3/2 dependence of the partition function (2). Rather, it interpolates be-   tween them, and is ∼t −2 , to the nearest integer. From what was learned from the supergravity cases of last section, the origin of this is clear. The endpoint of the distribution has a different structure (see figure 25), with some non-zero ρ(0) at the end. There are far fewer states than for the (0,2) supergravity case but more than the (2,2) case, and hence the fall-off rate (at least for this value of ) is between that of those two cases. Again, as seen in all the models studied in this paper (and also the special Airy model recalled in Appendix A) non-perturbative effects hasten the transition from dip to ramp to plateau such that the linear part of the classical contribution to the ramp that emerges at long times simply does not have time to develop, for moderate values of .

VII. JT GRAVITY: NON-PERTURBATIVE DEFINITION B
It is worth revisiting the idea mentioned at the beginning of the previous section that the presence of the nonperturbatively unstable even k members of the (2k−1, 2) minimal model family renders unstable the JT gravity from which they are built. This is what motivated the non-perturbative definition of ref. [3] in terms of other, more manifestly stable minimal models. There seems to be another possibility, however. It arises naturally in the context of this paper, and so it will be briefly discussed here, but details of further explorations will be presented elsewhere, if successful.
Simply put, perhaps the presence of the odd k minimal models helps the system avoid the pathologies of the lone even k models. From the perspective presented in ref. [4] that the construction of JT gravity from minimal models is equivalent to being in a particular model and then turning on operators that mixes in the presence of all the others, the view could be taken that the particular model is an odd k model (k can be taken to infinity at the end), and then the (k−1) lower models are turned on. The key here is that they are turned on a finite amount. The Hermitian matrix model string equation R=0 should have a well-defined solution in such a case. This is because the pathologies representing instabilities of the mth (even) model appear in this language only if the coupling t m is sent to infinity, performing the RG flow to land entirely in that model. In that case the solution develops increasingly steep oscillatory components (as t m grows) that herald the appearance of problematic double poles, as noticed early on in ref. [69]. However, the JT prescription only requires the couplings to be certain finite values (28). So, there could well be the right kind of solution present after all, probably a unique one, and it would supply the u(x) needed to define the full non-perturbative JT gravity model.
From the perspective of this paper, all that is needed is to find such a solution at a given truncation, and the truncation should be at an odd k for a solution to be found. Remarkably, a search for such a truncated solution eventually (after some careful coaxing: R=0 is a simpler equation than equation (9), but for the boundary conditions needed it is a skittish creature) yielded solutions of just the right sort. A truncation containing just t 1 , t 2 , and t 3 is presented in figure 29. In both directions it asymptotes to the solution of the cubic equation Meanwhile a solution (but for =3.26; see below) for the truncation up to t 5 is presented in figure 30. It asymptotes to the quintic equation: Notice that in this latter case 12 , the solution quickly approaches (for x<0) the exact disc level potential (29) (plotted as a (red) dotted line) that yields the Schwarzian 12 The solution for u(x) that was sought in this case was for =1, to match the other examples presented, and to be in the physically sensible regime (recall that =e −S 0 ≤1, since the leading entropy S 0 is supposed to be positive). However, this solution was for =3. 26. Getting a solution with all the way down to unity or below has so far proven rather difficult (the equations are easier to solve when the derivatives have larger coefficients), but the evolution of the shape as was reduced toward one seemed to be converging nicely, so it may be possible to complete the job with more care. density (4). So this could be an excellent starting point for solving the non-perturbative spectrum of JT gravity to good accuracy, and computing the spectral form factor. It will have a tail penetrating to negative energy (in this way it is a grown-up version of the k=1 Airy model reviewed in Appendix A), but this will be just fine if the effective potential for one eigenvalue in this region prevents tunnelling, as happens in pure odd k models [70], and so may well be likely here.
It would be very attractive to find a non-perturbative solution of the JT gravity Hermitian matrix model definition that uses the string equations of double scaled Hermitian matrix models after all, and this seems a very promising (and exciting) possibility. Further exploration of this proposal is left for future work.

VIII. CLOSING REMARKS
The purpose of this paper was to explicitly uncover and examine the non-perturbative physics for JT gravity and supergravity that is accessible if they are formulated using minimal model building blocks. This construction is not a simple large k limit of a minimal model, but a more refined affair involving coupling them together in a particular combination, as suggested perturbatively in ref. [11], and extended to non-perturbative physics in refs. [3,4]. The principal applications that demonstrated the facility of the technique was the explicit computation, for the first time, of the full non-perturbative spectral densities of various JT gravity and supergravity models, and the use of these spectral densities to compute the spectral form factor in each case, showing how the non-perturbative effects affect the shape (sometimes dramatically) of this important diagnostic quantity. Having explicit access to the non-perturbative features in this manner turned out to be rather instructive, as extensively discussed in the body of the paper.
Techniques to allow such non-perturbative properties to be extracted, in a consistent and well-defined scheme (for generic values of β and ), have not been presented in the literature before, and it is hoped that these methods and results will go some way to helping uncover more of the fascinating web of interconnections between geometry, quantum mechanics, gravity, and chaos that seems to be emerging from these studies.
It is likely that other models of JT gravity can be "deconstructed" in terms of minimal models in a way analogous to what was done here, and thereby be given a non-perturbative definition, not just in principle, but (as shown here) in useful accessible terms. This includes the possible alternative non-perturbative definition of ordinary JT gravity (using a smooth non-perturbative interpolating solution of the original Hermitian matrix model string equation) suggested in section VII. It could possibly also encompass some of the new kinds of matrix model descriptions of JT gravity black holes mentioned recently in ref. [73]. Perhaps the results of explorations along these lines will be reported soon. The correlator of two boundaries can be computed exactly using properties of the Airy functions out of which the wavefunction is built 13 . The disconnected part is simply the square of the partition function, which can be evaluated by Laplace transform, remembering to include negative energies to incorporate the full nonperturbative spectrum: while implementing equation (26) yields the connected piece to be: These are exact expressions, i.e., both perturbative and non-perturbative parts are incorporated (see refs. [43,74] and the review in the Appendix of ref. [11]). Nevertheless, it is instructive to pretend that only a finite number of wavefunctions are known, (only numerically), for a discrete set of energies up to some maximum energy. The question is then how well the exact expressions can be reproduced. This is the situation of the body of the paper, resulting from the controlled truncation of the infinite order string equation.
The answer to the question, reassuringly, is that a great deal of the important physics is accessible. To show this, a set of 1000 of the wavefunctions (A1) were discretized on the same size grid used in the body of the paper (x is broken up into 20,000 points), for a range of energies −20≥E≥20, and the same code that performed the numerical implementation of the expressions given in equations (20) and (26) were carried out for this exact model. The result presented in figure 32 shows the crossover between the two portions of the correlator as a function of β, and figures 33, 34 and 35 show the disconnected, As might be expected, significant deviations from the dashed lines occur when β becomes too small, including patterns of zeros representing finite size effects. This signals temperatures that excite higher energies that are not included in the numerical scheme (but are in the exact expressions) and hence the results deviate. As long as such extremes are avoided (depending upon the truncation energy chosen), the numerical results are very reliable. This is a good controlled model of the truncation scheme used in the body of the paper, showing that the results obtained are robust.  ical solutions, to help the interested reader learn how to extract useful information for themselves. First, Maple was used in this case (although MATLAB works well too, as probably would other programs). The dsolve routine was used, with an error tolerance of 10 −5 −10 −10 , depending upon the equation being solved. (In fact when these equations were first solved [48,75] it was for the cases k=1, 2, and 3, including cases where the models were non-trivially coupled [76]. Back then, their solution was found by writing a program in FORTRAN that called the routine D02RAF, part of the NAG libraries.) As mentioned in the text, for various reasons, it makes sense to take an additional derivative of the equation. This reduced the non-linearity somewhat, at the expense of increasing the order, which is a small price to pay. This is because the first derivative results in an overall factor of R, which can be divided out since it will not vanish for the solutions of interest. A derivative explicitly removes Γ from the equation, however. Now, the only knowledge the system has of the desired choice of Γ is through subleading (in the small expansion) terms in the boundary conditions, which need to be solved for with all t k present. This is nicely organized on x>0 boundary because the presence of the t k come in one by one at successively higher terms in the 1/x 2 expansion (see ref. [4]). It is less nice to do analytically on the x<0 boundary. There, all the t k contribute at the next order and solving for the order Γ correction requires solving a kth order polynomial. By beyond the k=4 truncation this becomes unpleasant at best. However, if the system is solved on a large enough region, with a small enough discretization, terms beyond the leading left boundary condition can be safely ignored, and a good approximation to the the solution found anyway. (If needed, however, a recursive code for solving for the subleading corrections to the boundary conditions numerically can be employed.) This works well for Γ= 1 2 , while for Γ=− 1 2 a different technique was used, because the numerical approach was less stable due to a more complicated well shape appearing in the interior, which is hard to control in a 13th order differential equation. Ref. [56] noticed that Γ can be changed by an integer using a special "Backlünd" transformation, and actually derived an analytic expression expression showing how to build the new u(x) at Γ±1 from the old u(x) (and its derivatives) at Γ. Here it is: (where R≡R(u Γ ) and a sign has been switched to match the current conventions). So once u(x) for the case of Γ= 1 2 was found using Maple, the output of dsolve contains all the derivatives of u(x) needed to construct u(x) at Γ=− 1 2 , giving the well structure seen in figure 7. In fact, a similar story held for the case of Γ=0 used in section VI. There is also a well structure in the interior, which is hard to solve for numerically when at high order and with boundary conditions far from the structure itself. Experience from studying positive integer Γ suggested that Γ=1 would be smoother to solve for and that was indeed the case. From there, transformation (B1) was used to construct the desired Γ=0 solution.

Suggestions for solving the spectrum
As mentioned in section IV, in order to solve for the solutions to the eigenvalue problem, {E, ψ(E, x)}. a matrix Numerov method [57] was used, as it was in refs. [3,4]. This simply puts the system into a box, and turns the problem into a large matrix diagonalization problem, for a given input potential u(x). This was done using MATLAB. A key point is that it is desirable to have a large number of eigenvalues in the energy range from zero to the chosen highest energy (determined by the level of the truncation of the string equation). So two choices were made to ensure a good set of solutions. The first was to use a large grid, so a grid of 20000×20000 was used. The second was to use a large box. As stated, the spectrum solving method is essentially putting the system into a box, and a portion of the output eigenvalues and eigenfunctions will be affected by the edges of the box. Those should be discarded, and the larger the box, the more useable eigenstates will be available in the reliable energy window. Since, as already observed in section III, the solution for u(x) becomes similar to the disc level behaviour far away enough from the central region, the box can be easily made larger by connecting the solution (solved numerically out to −200≤x≤+200) to a wider region (e.g., −2645≤x≤+2645 for the (2,2) and (0,2) models) where just the exact disc solution u 0 (x) is used. (A smooth (enough) transition between the two solutions was done at x<−100, corresponding to energies well above the cutoff determined by good matching for the truncation, so this does not affect the physics.) Appendix C: Gel'fand- Dikii Polynomials In case they are needed, here are some of the higher order Gel'Fand-Dikii polynomials. In the following equation a prime denotes an x-derivative times a factor of . For high numbers of derivatives, instead, a notation u (n) is used for n primes. The first five are listed here: It will transpire thatR 6 [u] andR 7 [u] will be needed as well, in order to get the required level of accuracy for the quantities computed in this paper. They are rather lengthy quantities, so it is not clear if there is much value in listing them here. Instead, they (and higher oder ones) can be easily computed using the recursion relation: and the requirement that they vanish when u does. More useful therefore are the following three lines of Maple code, which can be iterated, with obvious adjustments, to get to arbitrarily high order: