Equilibration in low-dimensional quantum matrix models

Matrix models play an important role in studies of quantum gravity, being candidates for a formulation of M-theory, but are notoriously difficult to solve. In this work, we present a fresh approach by introducing a novel exact model provably equivalent with low-dimensional bosonic matrix models. In this equivalent model significant local structure becomes apparent and it can serve as a simple toy model for analytical and precise numerical study. We derive a substantial part of the low energy spectrum, find a conserved charge, and are able to derive numerically the Regge trajectories. To exemplify the usefulness of the approach, we address questions of equilibration starting from a non-equilibrium situation, building upon an intuition from quantum information. We finally discuss possible generalizations of the approach.

Matrix models play an important role in studies of quantum gravity, being candidates for a formulation of M-theory, but are notoriously difficult to solve. In this work, we present a fresh approach by introducing a novel exact model provably equivalent with low-dimensional bosonic matrix models. In this equivalent model significant local structure becomes apparent and it can serve as a simple toy model for analytical and precise numerical study. We derive a substantial part of the low energy spectrum, find a conserved charge, and are able to derive numerically the Regge trajectories. To exemplify the usefulness of the approach, we address questions of equilibration starting from a non-equilibrium situation, building upon an intuition from quantum information. We finally discuss possible generalizations of the approach.
Matrix quantum mechanics has received significant attention in recent years, mainly for its suspected connection with quantum gravity. Most prominently, Banks, Fischler, Shenker and Susskind (BFSS) [1] have proposed that supersymmetric matrix quantum mechanics, called matrix theory, gives a formulation of M-theory. This is an exciting proposal about the fundamental degrees of freedom in nature, which may allow us to address questions that are out of reach of semi-classical gravity. Questions of particular interest are on the quantum properties of black holes, such as their microscopic constituents, thermalization, and the process of evaporation [2].
BFSS matrix theory is the maximally supersymmetric Yang-Mills theory in 0 + 1 dimensions with U (N ) gauge symmetry. The bosonic degrees of freedom of this theory are N × N Hermitian matrices X (i) , i = 1, . . . , d. Their diagonal elements denote the position of D0-branes-point-like nonperturbative objects in string theory-and their off-diagonal elements represent the lowest modes of open strings stretched between D0-branes; i is the spatial index, and the dimension of space for the BFSS theory is d = 9. Configurations that are related by U (N ) (under which X (i) transforms in the adjoint representation) are to be identified. The U (1) part (X (i) being proportional to the unit matrix) is decoupled from the rest (and is interpreted as the center of mass position), so we will concentrate on the case of traceless matrices with SU (N ) gauge symmetry. There is no mass term in the action, and the interaction among the X (i) is given by a quartic term.
Unsurprisingly, this is an intricate theory and it has not been solved to date. To start with, it is a highly non-local model, in that every local degree of freedom directly couples to any other degree of freedom. The interaction graph is that of the complete graph. What is more, the model is quartic with strong interactions, resulting in the situation that naïve approximation schemes such as perturbation theory appear doomed to failure.
This does not mean that important features of the model are not known. Relevant developments include: A proof of the existence of a zero-energy supersymmetric bound state [3], and the continuous spectrum above it [4]; a qualitative understanding of scattering amplitudes [5]; perturbative computations of the effective potentials around block diagonal configurations [6]; the study of finite temperature theory by Monte Carlo methods and the comparison with black hole thermodynamics [7]; computations of correlation functions from gaugegravity correspondence [8] and Monte Carlo simulations [9]. Still, key questions are open, essentially originating from the non-local highly interacting nature of the model.
In this work, we present a fresh approach to the problem, which constitutes a major simplification of the problem in low dimensions. We first present a new equivalent formulation in low dimensions, taking full advantage of the gauge symmetry and exhibiting an aspect of locality in the model. This formulation is amenable to analytical study and precise numerical analysis. We expect this to be an interesting toy model for higher-dimensional situations. The mindset is that instead of aiming at, say, simulating the classical version of the model in high dimensions with Monte Carlo methods, we present an exact analysis for low dimensions, but one where most features can be accurately discussed. It turns out that in this reformulation, questions of equilibration [10,11] in non-equilibrium can be precisely posed.
Matrix models and relations to other models. We will study bosonic matrix models defined by the Hamiltonian where a t a , with the N × N traceless Hermitian generators t a of SU (N ), and position operators x a are the conjugate momenta. This Hamiltonian describes 0 + 1 dimensional Yang-Mills theory. The gauge field is not dynamic and can be set to zero by a gauge transformation; its equation of motion imposes the constraint that states under consideration are required to be singlets of SU (N ). In this work we will primarily treat the d = 2, N = 2 case [12]. (1) for N = 2 and d = 2, with the constraint that the states are SU(2) singlets, is equivalent with

Theorem 1 (Equivalent model) The Hamiltonian
all other elements being zero.
Proof of equivalence. We now prove that the above models are indeed equivalent. We start from the explicit form of 2×2 traceless Hermitian matrices, with Pauli matrices σ a . Expressing P (i) similarly, we have a set of pairs of canonical coordinates, In what follows, it is convenient to turn to the position representation in radial coordinates, in which the interaction term ( is the angle between 1 and 2 and r 1 , r 2 ∈ (0, ∞). Singlet states are given by linear combinations of ψ Accordingly, the state vector in the position representation takes the form ψ(r 1 , The Hamiltonian in the radial position representation becomes −∂ 2 , where ∆ S is the Laplace operator on the unit sphere. We rewrite the interaction term, using sin 2 (θ 1, where θ 1,2 is the angle on the great circle, and we use the expansion ψ with the Legendre polynomials P ℓ . This relation can be seen, e.g., by setting θ 2 = 0 using the rotational symmetry. The coefficients c ℓ,2,ℓ ′ are related to the Wigner 3-J symbol, The Hamiltonian as an operator acting on ρ takes the form of (2) with 3A ℓ, (3). Discussion of the model. Interestingly, the bandedness of A introduces a notion of locality to the model not apparent in the original form Eq. (1). There is a fast convergence A ℓ,ℓ → 1/2, A ℓ,ℓ±2 → −1/4, also ℓ(ℓ + 1) ≈ ℓ 2 , for growing values of ℓ. Hence, the part acting on l 2 reminds of the harmonic oscillator on a lattice [13]. The even and odd sublattices are decoupled, as there are no non-zero matrix elements coupling ℓ and ℓ ± 1. The dynamics on L 2 (Ê + ) ⊗ L 2 (Ê + ) determines joint effective "mass" and "spring constants" for the two indirectly coupled systems in l 2 , providing indirect interaction between these systems.
The model is invariant under the spatial SO(2) rotation. The generator in the representation of Eq. (4) is (1) , and we have an equivalent operator Q in the representation of Eq. (2) , all other elements being 0. Note that because of the invariance of the model under the action of Q, we have [Q, H] = 0, and the Hilbert space is divided into eigenspaces of Q. The spectral values of Q are q ∈ 2 , the even integers. This angular momentum will become important in the analysis of the model.
Numerical analysis. The aim of our numerical analysis is the low energy regime of the model. We compute the spectrum of H in the different eigenspaces of Q, as well as the spatial extension of the state, X 2 ⊗ ½ + ½ ⊗ X 2 1/2 . For our (banded) Hamiltonian, the operator norm error can easily be bounded, which yields a bound to the spectral values of the true Hamiltonian, and gives us confidence in the significance of our results.
The Hilbert space of this model is infinite dimensional, so we make use of a finite-dimensional approximation. We stress that without the presented reformulation, already the smallest non-trivial instance of the model-N = 2 and d = 2seems inaccessible even on super-computers, and it would be also difficult to implement the gauge constraint numerically. Without the transformation we perform, we are facing a Hilbert space H = L 2 (Ê) ⊗6 without a local structure, the exponent growing rapidly with N and d. Sparse methods allow for approximations of the position and momentum operators on local systems in such dimensions with, say, 10 × 10-matrices, which yields a very crude approximation, unlike our approach. Due to reformulation (2), living on H = L 2 (Ê + ) ⊗ L 2 (Ê + ) ⊗ l 2 we are able to relax the strong dimensional restrictions in a finite-dimensional approximation considerably. We restrict L 2 (Ê + ) to the Hilbert space spanned by the first h 0 odd eigenfunctions of the harmonic os-cillator, and we restrict the dimension of l 2 to some value ℓ 0 . For a good asymptotic overview, we choose h 0 = 107 and ℓ 0 = 156. This choice provides about 800 good eigenstates. The computation of the spectrum in the different eigenspaces of Q is performed by first diagonalizing H and then computing the expectation value of Q 2 for each energy eigenstate. The spectrum of H is non-degenerate and hence eigenstates of H are eigenstates of Q. We double-check the results with a priori knowledge that Q takes integral values and Q 2 has a two-fold degeneracy (from Q ↔ −Q). Eigenstates with too large deviations from this are discarded, and because we pay attention to their order, we will discard all of them after the first discarded energy.
Results. The ground state energy is E GS = 1.05535 . . . . We now use the label q ∈ 2 , the even integers, for the eigenvalue of Q. We enumerate the energies within each Qeigenspace with a parameter n ∈ AE 0 , hence we have states {|q, n } with Q|q, n = q|q, n and H|q, n = E q,n |q, n . The data suggest power laws both in q and n, as does a semiclassical analysis for the model. There is a scaling argument which suggests the dependence E q,0 ∼ |q| 2/3 [14]. One can fit the data for n = 0 with (E q,0 + E 0 ) α = a + b|q|, a relation which qualitatively holds very well for a range of exponents α ≈ 1.5 . . . 2.3 provided the other constants are chosen appropriately. For α = 2 this relation is known as the linear Regge trajectory and corresponds to the behavior of a relativistic string [15], which has been observed in the 2 + 1 dimensional Yang-Mills theory [16]. The best fit we obtain when leaving all constants subject to variation yields α ≈ 1.62(2), E 0 ≈ 1.6(1), a ≈ 2.8(3) and b ≈ 1.0(1). As the figure shows, the mapping (n, q) → (E q,n + E 0 ) α is an affine function of both n and q-the prefactor of n is c ≈ 1.9(1)-up to certain values of n, where we notice almost degenerate pairs. After that, linear growth continues, see Fig. 1. Low values of |q|, especially q = 0, are dominated by the degeneracies and hence the linear growth is mostly hidden, but the onset of degeneracies comes later for larger values of |q|. Note that the growth of (E q,n + E 0 ) α with n stays always sublinear for all values of q. We also consider the spatial extension of the states, X 2 ⊗ ½ + ½ ⊗ X 2 1/2 . The extensions within each eigenspace of Q are affine functions of E, except that around the degenerate pairs states have a much smaller size.
Equilibration in matrix models. Equipped with these insights, we now turn to questions of equilibration. Such considerations complement studies of classical limits of matrix models [17], as well as ones that give rise to bounds to scrambling based on the locality of the models and Lieb-Robinson bounds [18]. In the quantum setting, 'equilibration' refers to the situation that for most times, expectation values take values as if the system was in the time-averaged state ω := 1 T Depicted is the value of (Eq,n + E0) α on the vertical axis with n on the horizontal axis, using α ≈ 1.62. The black line belongs to q = 0 and the others to increasing |q|, from green to red. There is affine growth of (Eq,n+ E0) α over both n and |q| in certain regions. The affine growth over n stops at certain points where we notice almost degenerate pairs of energies. The growth over n is always sublinear. (Small plot: sizes) Depicted is the value of X 2 ⊗½+½⊗X 2 1/2 for fixed q = 3, as an example. This value is approximately an affine function of H , with occasional collapses to much smaller values at the positions of the almost degenerate energies.
average is the effective dimension d −1 eff := n | ψ(0)|q, n | 4 . In fact, for any d S -dimensional subsystem S the expected deviation in trace-norm from the time average is For the above matrix model, we immediately get an estimate for the effective dimension of a random state of a microcanonical energy window, slightly improving on Ref. [10]. Denote with d ∆ the number of spectral values of H q contained in an interval [E, E + ∆] for a given fixed value of q, and let I denote the corresponding index set k ∈ I ⇔ E q,k ∈ [E, E + ∆] in this eigenspace of Q.
Using the Weingarten function calculus for computing moments of entries of Haar-random unitaries U , and making use of the fact that each entry is identically distributed, one finds Obviously, having a large number d ∆ of energy levels of H q within an energy window [E, E + ∆] implies also that there are observables that will take long to equilibrate.

Observation 3 (Equilibration times) Consider the initial state vector
Then the deviation from the infinite time average ω = i∈I |q, i q, i|, normalized by the operator norm of the observable O ≤ 2, fulfills This follows immediately from the fact that the left hand side is bounded from below by i,j∈I which in turn is bounded from below by the right hand side.
We now assume that H q is not degenerate. This is supported by the numerical results, although step-like features in the spectrum exist, which are very close to degeneracies. This does not, however, invalidate the following argument unless a degeneracy is exact. Furthermore, there is strong evidence presented above that within each eigenspace of Q we have sublinear growth of (E q,n + E 0 ) α with n, where α = 1.62 in Fig. 1. Hence for ∆ ≪ E, d ∆ ≥ α∆(E + E 0 ) α−1 /c, where c is the n-proportionality factor in (E q,n +E 0 ) α = a+b|q|+cn, such that, with Observation 2, That is to say, for large initial energies, one expects a strong equilibration and expectation values often take the values close to the ones of the time average. The sublinear growth of (E q,n + E 0 ) α allows to find an upper limit for the gaps in the spectrum. We find that for energies within one eigenspace of Q and above E we have (E q,i+1 − E q,i ) ≤ c/α(E + E 0 ) α−1 , so we can estimate with Observation 3 that So for this specific initial condition, the system will be close to the infinite time average in expectation, but we can make the equilibration time scale as large as we want by shifting the energy up. These statements can only be concluded from the microscopic Hamiltonian once the spectrum as above has been identified, as facilitated by our new normal form.
Paths towards generalizations. To establish a relation between matrix models and black hole physics, we may need d ≥ 3, since gravity in lower dimensions is special. The case of higher N and d is under study. Hamiltonian (1) with general N and d can be mapped to a model which involves d particles in a N 2 − 1 dimensional space spanned by the adjoint representation of SU (N ). Wave functions, which have to be singlets of SU (N ), are parameterized by r i = (2Tr X (i) X (i) ) 1/2 , i = 1, . . . , d, together with the parameters which specify the SU (N ) representation. For N = 2, d = 2, there is only one such parameter ℓ, but for d ≥ 3, singlets are constructed from more than two parties, and we need more parameters. However, there is still a notion of locality in the gauge space, since each term in the interaction Tr X (i) , X (j) 2 can change the SU (N ) quantum number only by a certain amount. We believe this locality provides an important clue for analytical and numerical studies.
Conclusion. We have constructed a model which is equivalent to the N = 2, d = 2 bosonic matrix model, and studied its spectrum as well as its equilibration properties. We found an affine dependence of (E q,n + E 0 ) α , with plausible α ≈ 1.5 . . . 2.3, on |q| and n in certain regions. The dependence on |q| for the lowest n (the leading Regge trajectory) is well-known to be the behavior of a string if α = 2, and so is the spatial extension, which we found to be proportional to the energy (except for certain states). It is remarkable that this kind of non-perturbative behavior, which is usually found with extensive Monte Carlo simulations in lattice gauge theory, is found by straightforward diagonalization in our approach.
Although the 1/N expansion of any gauge theory can be represented as strings [19], it is far from obvious what kind of string theory our model with N = 2 should correspond to. In fact, the power law dependence on n without any degeneracy (up to some n, for fixed q) resembles a system of one or two oscillators, such as [20], rather than strings, which have an infinite number of modes. The states that have small spatial extension, found at larger n with fixed q, may represent composite states. We hope our result can serve as data for the understanding of gauge/gravity correspondence for small N .
Regarding the question of equilibration, we found convincing numerical evidence that for high energies, states will for almost all times look like their long time average. We were able to show, however, that observables and state preparations exist, where the equilibration lasts arbitrarily long.
Matrix quantum mechanics studied in this work is obtained from pure Yang-Mills theory in d + 1 dimensions by a dimensional reduction (i.e. by assuming that the fields depends only on time). The Lagrangian of the d + 1 dimensional theory is given by where the repeated indices µ, ν = 0, 1, . . . , d are summed over. The field strength is Gauge field components A µ are associated with the Lie algebra of the gauge group; for SU (N ), they are represented as traceless Hermitian N × N matrices. The gauge transformation by an SU (N ) element U is given by and the Lagrangian is invariant under this transformation. The Lagrangian of matrix quantum mechanics (Yang-Mills theory in 0 + 1 dimensions) is obtained by setting the spatial derivatives to zero. Writing A 0 = A and A i = X (i) , the field strengths become and we obtain where the repeated indices i, j = 1, . . . , d are summed over. The infinitesimal form of a gauge transformation for U = e −iΛ is In 0 + 1 dimension, the gauge field A is not a dynamical field, since it has no kinetic term. We can set A to zero by a gauge transformation; A only plays the role of Lagrange multiplier, which imposes a constraint (equation of motion w.r.t. A) on the system. We will perform the canonical quantization in the A = 0 gauge. Consider each real number of the matrix elements as a dynamical variable. The trace part represents the center of mass. The corresponding momentum is conserved, and its dynamics is decoupled from the rest. So we will concentrate on the study of the relative motion.
For the explicit analysis, it is convenient to expand the fields in the basis t a , where t a are traceless Hermitian matrices which satisfy and where f abc is the structure constant of SU (N ). For the SU (2) gauge group, studied in the main text, we have t a = σ a /2 and f abc = ǫ abc . Gauge transformation (5) for x (i) a with Λ = a λ a t a is The conjugate momenta for x For SU (2), by using 3 a=1 ǫ abc ǫ ade = δ bd δ ce − δ be δ cd , and redefining x (i) a → 2 1/3 x (i) a , p (i) a → 2 −1/3 p (i) a and H → 2 −4/3 H we get Eq. (4). The equation of motion for A (obtained by varying the Lagrangian by a a and setting a a = 0 afterwards) is Noether's theorem tells us that the V a are the generators of SU (N ) transformations. We impose constraint (7) by requiring that the physical state vectors |ψ are annihilated by the generators V a , V a |ψ = 0.
In other words, |ψ are singlets of SU (N ).

Chew-Frautschi plot
In Fig. 2 we show a Chew-Frautschi plot of the model. An asymptotic affine dependence between (E q,n + E 0 ) 2 and |q| is clearly seen in the large energy / large angular momentum regime.
FIG. 2. Chew-Frautschi plot and size plot (large plot: energies). Depicted is the value of (Eq,n + E0) 2 , similarly to Fig. 1, but fitting the data while fixing α = 2. In this plot, we show graphs belonging to constant values of n, while varying over values of the angular momentum q. An asymptotic affine dependence between (Eq,n + E0) 2 and |q| is clearly seen in the large energy / large angular momentum regime. (Small plot: sizes) Depicted is the value of X 2 ⊗ ½+½⊗X 2 1/2 for fixed n = 11, as an example, over the full range of q. This graph follows, asymptotically for large |q| and n, an affine dependence as well, with a low energy region where the size and energy are smaller than a linear extrapolation of the asymptotics implies. Is resembles the behavior of an initially strongly bound system which is less tightly bound for high energies.