Hamiltonian Truncation Study of Supersymmetric Quantum Mechanics: S-Matrix and Metastable States

We implement the Rayleigh-Ritz method in supersymmetric quantum mechanics with flat directions, and extract the S-matrix and metastable resonances. The effectiveness of the method is demonstrated in two strongly coupled systems: an N=1 toy supermembrane model, and an N=4 model with a U(1) gauge multiplet and a charged chiral multiplet.

1. We separate the Hamiltonian into a free part H 0 and an interaction potential V , and introduce an IR cutoff at a sufficiently large distance L.
2. We work with a basis of wave functions that diagonalize the free Hamiltonian H 0 in a "box" of size L, and truncate the basis to a finite set by keeping eigenfunctions up to truncation energy Λ. We then evaluate the matrix elements of the full Hamiltonian H 0 + V on this truncated basis, and diagonalize it numerically to find the energy levels. The separation of free and interaction Hamiltonian is such that the interaction potential V is sufficiently smooth, which leads to small mixing of low frequency modes with high frequency modes, and fast convergence of the energy levels with increasing Λ.
3. The energy levels E n of the truncated Hamiltonian for a scattering spectrum are dense in the limit of large L. We subtract from the number of states n up to energy E n a universal IR contribution governed by the asymptotic scattering wave function, and obtain a "renormalized number of states" n(E n ). We refer the collection of points (E n , n(E n )) as the spectral set. 4. We collect the spectral sets for different and sufficiently large values of IR cutoff L. Provided that the truncation energy Λ is sufficiently high, over a finite range of energy E of interest, we will find that the spectral set lies on the union of k smooth curves, k being the number of asymptotic regions, or the effective dimension of the S-matrix at a given energy. In the k = 1 case, where the S-matrix is a single scattering phase, the curve traced out by the spectral set determines the scattering phase φ(E) as a function of energy. In the k > 1 case, by varying the details of the IR cutoff on different asymptotic regions, one can determine the full U (k) S-matrix, up to a relative U (1) k−1 phase ambiguity of the asymptotic scattering wave functions. 5. A metastable state corresponds to a jump of the scattering phase φ(E) by 2π. The spectrum of metastable states and their decay width are determined from the peaks of dφ(E)/dE. Our procedure also allows for the explicit determination of the wave function of metastable resonances.
Our method is explained in more detail in the next section. The convergence with truncation energy is discussed through one dimensional examples in section 3. The application to the toy supermembrane model and the N = 4 SQM will be presented in section 4 and 5 respectively. We conclude with some prospectives on the Rayleigh-Ritz approach to holographic models in section 6.

S-matrix and the density of states
Let us consider a quantum mechanical system with k one-dimensional asymptotic regions where the dispersion relation takes the form E = p 2 , p being the asymptotic momentum. The in and out states will be denoted |E, i in and |E, i out respectively, i = 1, . . . , k. They are related by the S-matrix We now introduce an IR cutoff by placing hard walls at distance x = L i in the i-th asymptotic region, so that the spectrum is discretized. The quantization condition amounts to demanding that the k × k matrix admits a zero eigenvalue, or equivalently, det B(E) = 0.
In the simplest k = 1 case, the S-matrix is a single scattering phase e iφ(E) . The asymptotic quantization condition is Given the spectrum {E n }, we can extract the scattering phase from where the renormalized number of states is given by In practice, we only need to take L to be greater than the effective range of interaction, and collect the spectral set E n , 2π(n + 1 2 ) − 2L E n , n = 0, 1, 2, . . . If there is a metastable state of energy E * and decay width , we expect the scattering phase to behave as for E close to E * . This leads to a peak in the derivative of the scattering phase, (2.10) Now let us generalize this prescription to the case of several asymptotic regions, namely k > 1. To begin with, take L i = L, and denote by e iφ j (E) (j = 1, . . . , k) the eigenvalues of S ij (E). The asymptotic quantization condition can be written as φ in (E n ) + 2L E n = 2π(n + 1 2 ), n = 0, 1, 2, . . . , i n ∈ {1, . . . , k}. (2.11) In other words, now the spectral set (2.8) lies on the union of k curves that are the graphs of the functions φ i (E), i = 1, . . . , k. Numerically, it is again useful to combine the sets (2.8) for different values of L (provided that they are larger than the effective interaction range).
A slight modification of this prescription allows for extracting the full U (k) S-matrix, up to conjugation by a diagonal unitary matrix (since the choice of phase for each asymptotic wave function is a priori ambiguous). Let us take L i = L + y i , with y i finite while taking the large L limit. Now the set of points (2.8) lie on the k curves defined by the k eigenvalues of the matrix as functions of E.
As a nontrivial example, consider the k = 2 case. Up to conjugation by a diagonal U (2) matrix, we can write the S-matrix as is the overall phase, and α(E) and θ(E) are to be determined from the spectral set (2.8). Now take y 1 = y/2, y 2 = −y/2, so that (2.14) The eigenvalues of S are exp i φ 1 (E)+φ 2 (E) 2 ± iγ(E; y) , with cos γ(E; y) = cos θ(E) cos(α(E) + y √ E ). (2.15) γ(E; y) is determined from the difference between the two curves traced out by the spectral set. By maximizing | cos γ(E; y)| with respect to y, we then determine both θ(E) and α(E).
3 Quantum scattering in one dimension from Rayleigh-Ritz method: dependence on truncation energy To gain some intuition, let us consider the scattering problem in one dimension, with the Hamiltonian H = −∂ 2 x + V (x) for x > 0, where the potential V (x) vanishes sufficiently fast in the x → ∞ limit. We also impose the boundary condition that the wave function vanishes at x = 0. In the Hamiltonian truncation approach, we place a hard wall at x = L, and work with the truncated basis of standing waves in the empty box of length L, ψ n (x) = 2 L sin nπx L , n = 1, 2, . . . , N. (3.1) We shall choose N such that the energy at the truncation level is much bigger than the scattering energy E. The matrix elements of the Hamiltonian on this basis are given by One then proceeds to diagonalize the N ×N truncated Hamiltonian matrix H (N ) ≡ (H nm ) 1≤n,m≤N , and read off the scattering φ(E) from the renormalized number of states up to energy E as in (2.5).
The efficiency of the Hamiltonian truncation method relies on the convergence rate of the eigenvalues of the truncated Hamiltonian H (N ) with increasing N or Λ. This may be Assuming the latter has bounded k-th order derivative in x, the off-diagonal H nm decays with n at least as fast as (n/L) −k−1 , and its contribution to the eigenvalues of the Hamiltonian is bounded by (n/L) −2k−4 . The error due to truncation at N may be estimated by summing over n ≥ N , giving a result that scales like N −2k−3 ∼ Λ −k− 3 2 . If the potential V (x) is smooth, we expect the eigenvalues of the truncated Hamiltonian converge with increasing truncation energy Λ faster than any inverse power of Λ.
As an example, consider the delta function potential, The matrix elements of the Hamiltonian on the basis (3.1) are given by The eigenvalues of the truncated Hamiltonian H (N ) are given by solutions of In this case, the error due to the level truncation can be estimated from The error in the energy levels is suppressed by E Λ , agreeing with the expectation of Λ −k− 3 2 for k = −1. This is also confirmed by numerics in Figure 1(a).
If we consider a rectangular potential V (x), which is bounded but its first order derivative in x unbounded, the energy levels are expected to converge with truncation energy as Λ − 3 2 . The numerical results are shown in Figure 1(b). On the other hand, if the potential V (x) is smooth, we expect faster-than-inverse-power convergence with the truncation energy. An example of the form Figure 1(c), and indeed the convergence with truncation energy Λ appears to be exponential.
The lesson here is that, in order to apply the Hamiltonian truncation method efficiently, we should split the Hamiltonian into a "free" part and an interaction "potential", H = H 0 + V , such that while the basis functions diagonalize H 0 , the interaction potential V is sufficiently smooth so that its matrix elements V nm decays sufficiently fast in the n → ∞ limit for fixed m.

The toy supermembrane model
The toy supermembrane model, introduced in [42] (see also [20]), is one of the simplest nontrivial SQM that admits flat directions and a gapless spectrum of scattering states. It is based on a Hilbert space of two-component wave functions in two variables x, y, with a single Hermitian supercharge Q, given by The Hamiltonian is The model admits a discrete symmetry of the dihedral group, which is represented projectively by the generators where P x and P y are parity transform in x and y respectively, and R is the rotation by π 2 on the plane: x → y, y → −x. They obey

The scattering problem
For a given energy E > 0, there are four scattering in-states |E, i in , and four out-states |E, i out , related by The index i = 1, . . . , 4 label the four asymptotic regions x → ∓∞ (y → 0), and y → ∓∞ (x → 0). The asymptotic states form a representation R of the dihedral symmetry group. Explicitly, the asymptotic scattering wave function of |E, 1 in takes the form and has only outgoing waves in the other three asymptotic regions, namely x → +∞ and y → ±∞. Likewise, the asymptotic wave functions of |E, i in , i = 2, 3, 4, obey and each of them has only outgoing waves in the other three asymptotic regions. The asymptotic wave functions for the out states can be constructed similarly. The dihedral symmetry generators are represented on these asymptotic states by the 4 × 4 matrices The supercharge Q is represented by a matrix ρ in (Q) on the basis of in-states and by ρ out (Q) on the basis of out-states, where The S-matrix lies in a singlet representation contained in R * ⊗ R. Further demanding time reversal symmetry and supersymmetry fixes S(E) up to a phase, namely

Hamiltonian truncation
To proceed, we place hard walls at x = ±L and at y = ±L, and work with the basis of wave functions for each component of ψ(x, y). We can simplify our task slightly by restricting to the U = +1 sector, which amounts to working with the basis , n odd. In the U = +1 sector, the S-matrix has two eigenvalues, ±e iφ(E) . This corresponds to a special case of (2.11) with k = 2, and we can extract φ(E) from the sequence of energy In Figure 2 we plot the density of the wave functions of a scattering state at generic energy and of a metastable state with very small decay width. The latter expectedly has localized support near the origin. The convergence of our numerical results with the truncation energy is exponentially fast, as confirmed in Figure 3.
Our numerical results of scattering phase are shown in Figure 5 and 6. Note importantly that while the scattering phase φ(E) is a priori defined modulo 2π, an abrupt increase of the  ?(E)  phase by 2π as computed from the "renormalized number of states" in energy is physical and signals a metastable or stable bound state. In principle, the metastable states can be detected from poles in the analytically continued S-matrix at complex values of E with negative imaginary parts. When a pole is close to the real axis, it can be detected as a Breit-Wigner peak in dφ(E)/dE.
Going to higher energies is more demanding numerically, as it requires increasing the truncation energy Λ as well as increasing the IR cutoff L. In the energy range where our numerical results are reliable, we find a curious scaling behavior of the renormalized number of states,n(E) ∼ E 1.88 at high energies. Note that a naive semi-classical quantization would suggest a lower bound ∼ E 3 2 , which is also the result one would obtain using a 1-loop truncated Schwinger-Dyson equation [30].  Up to energies where the Hamiltonian truncation approximation is reliable, we observe a power law E p with p ≈ 1.88.

Metastable states
The metastable resonance energies and widths can be read off by fitting the peaks of the dφ(E)/dE plot to a Breit-Wigner distribution (2.10), with the results shown in Figure 7.  excitation level of the y-oscillator is governed by the effective Hamiltonian H eff = −∂ 2 x + (2n + 1)|x| + σ 3 x (4.14) in the large x region, where n = 0, 1, 2, . . .. One component of the n = 0 wave function corresponds to the asymptotic scattering state. The other modes are subject to a confining effective potential V (x) = 2m|x|, m = 1, 2, . . ., giving rise to potential metastable states.
We may try to estimate the metastable resonance energies and their degeneracies using semi-classical Bohr-Sommerfeld quantization for this effective Hamiltonian. The result is E m, = (3πm ) 2 3 , where = 1, 2, . . . is the oscillation period in the effective potential in units of 2π , and there is a further 4-fold degeneracy due to the two components of the wave function and the two asymptotic regions of large positive and negative x. This rough estimate is compared to the direct enumeration of resonance peaks in dφ(E)/dE, as shown in Figure 8 (where only the U = +1 supermultiplets are counted), with qualitative agreement at moderate energies. 1

An N = 4 supersymmetric model
Now we turn to the scattering problem in a simple but nontrivial N = 4 SQM, namely that of a U (1) gauge multiplet coupled to a charged chiral multiplet [7,43], with vanishing Fayet-Iliopoulos parameter. Here we follow the convention of [43]. The Hilbert space consists of 16-component wave functions on R 3 × C, parameterized by coordinates x ∈ R 3 associated with the gauge multiplet, and (z,z) ∈ C associated with the chiral multiplet. The internal degrees of freedom of the wave function comes from quantization of gauginos λ α ,λ β from the gauge multiplet and fermions ψ α ,ψ β from the chiral multiplet. The fermions obey anticommutation relations There are four supercharges, given by Here our convention for the antisymmetric tensor is αβ = − αβ , 12 = 1. The supersymmetry algebra takes the form where G is the U (1) gauge rotation generator G =z∂z − z∂ z − ψ α ∂ ψα + 1, (5.4) and the Hamiltonian H is given by (5.5) The wave functions are restricted to be invariant under G.

Supermultiplets
The model admits the symmetry SU (2) J × U (1) R , where the SU (2) J rotates the R 3 as well as the fermions, whereas the U (1) R rotates the fermions only. Their generators are The supercharges Q α transforms as a doublet under SU (2) J and carry R-charge −1. Likewise, Q α carry R-charge +1. There is also a Z 2 symmetry taking R to −R.
The Hilbert space splits into sectors of R-charge R = 0, ±1, and ±2. The R = 2, SU (2) spin j states are necessarily annihilated byQ α , and form a multiplet with R = 1 spin j ± 1 2 and R = 0 spin j states. The R = 2 sector does not admit scattering states, and thus they belong to supermultiplets that consist of only normalizable energy eigenstates, which we refer to as R = 2 multiplets. Likewise, there is another multiplet consisting of stable R = −2 and R = 0 states of spin j, and R = −1 states of spin j ± 1 2 . All scattering states must lie in supermultiplets that consist of R = ±1 states of spin j and R = 0 states of spin j ± 1 2 . We refer to these supermultiplets as R = 1 multiplets. The asymptotic region is given by the simultaneous limit of r = | x| → ∞ andr = √ 2|z| → 0.

Stable states from the R = 2 multiplets
It suffices to examine the R = 2 sector wave function, whose r,r dependence is governed by the reduced Hamiltonian where [43] H 0 = − 1 2r 2 ∂ r r 2 ∂ r − 1 2r ∂rr∂r + j(j + 1) The Hamiltonian (5.7) has a discrete spectrum and thus we can adopt the standard Rayleigh-Ritz method, imposing an IR cutoff r < L andr <L and demanding that the wave function vanishes at r = L and atr =L. Note that the wave function need not vanish at r = 0 or at r = 0.
Denote by f j,n (r) = a j,n √ r J j ( 2ω j,n r) Note that here j is an integer, and so the r-dependence of the basis wave functions can be expressed in terms of elementary functions, while ther-dependence is expressed through the Bessel function J 1 .
We then diagonalize the matrix n, m|H|k, = δ nm δ k ω j+ 1 2 ,n with a suitable level truncation on n, m, k, . Numerical results are shown in Figure 9.

The R = 1 sector
To compute the S-matrix of scattering states and metastable resonances, we will study the R = 1, spin j sector using the Hamiltonian truncation method. In the spectrum we will find scattering states that are in a supermultiplet with R = 0 and R = −1 states, as well as a discrete set of stable states that are in a supermultiplet with R = 0 and R = 2 or R = −2 states. By comparison with the spectrum of R = 2 sector, we will be able to remove the stable states from the spectrum of the R = 1 sector, and identify the density of scattering states (subtracting IR divergent contribution) and thereby their scattering matrix.
The reduced Hamiltonian of the R = ±1 sector can be written in the form  We will then diagonalize the matrix n, m, s|H j,|R|=1 |k, , t = δ nm δ k δ st ω j+δ s2 +δ s4 ,n L 2 + ω δ s3 +δ s4 ,m L 2 +  Note that despite that the wave function has 4 components, there is only one scattering state at a given energy E. Therefore, the S-matrix is determined by a single scattering phase. The asymptotic wave function takes the form The number of states only has a linear divergence in L. After subtracting off the IR divergence, we obtain the renormalized number of states In Figure 11 we plotted the density of the wave functions in the (r, r) plane of a scattering state at a generic energy and of a metastable state with very small decay width.

Scattering states from the R = 1 multiplets
The stable states in the R = 2 sector with spin j are in the same supermultiplet as R = 1 states with spin j ± 1 2 . Therefore, the spectrum of scattering states of R = 1 and spin j can be obtained by removing from the R = 1 spectrum stable states in the same supermultiplet as R = 2, spin j ± 1 2 states. The resulting scattering phase φ(E), and the metastable resonances from dφ(E)/dE, are shown in Figure 12. Figure 13 shows the metastable resonance energies and widths obtained by fitting to a Breit-Wigner distribution.
As j increases, the effective potential pushes the wave function away from the origin, and it becomes important to correct for the finite L effects in order to obtain accurate scattering phases at low energies. This is achieved by comparison with the free Hamiltonian with an effective rotational potential j(j+1) 2r 2 at the same IR cutoff L, and subtract off the contribution ton(E) of this system (which admits a trivial S-matrix). The scattering phases (up to constant shifts) for j up to 21 2 are shown in Figure 14. The results are in agreement with the Coulomb branch effective Hamiltonian [44] which is valid in the limit of fixed E and large j, that gives a trivial scattering phase (when the FI parameter is set to zero). The widths of the metastable resonances generally become narrower as j increases, as the wave functions are pushed to large values of r, and the excited modes inr direction become more stable against decaying to the asymptotic states. n(E) j=1/2 j=3/2 j=5/2 j=7/2 j=9/2 j=11/2 j=13/2 j=15/2 j=17/2 j=19/2 j=21/2 Figure 14: Numerical results ofn(E) for the R = 1 multiplets with j = 1 2 , 3 2 , . . . , 21 2 as indicated in the legend, after removing the stable states in R = 2 multiplets from the R = 1 sector, and correcting a finite L effect due to a rotational effective potential j(j + 1)/2r 2 that affects low energy scattering states. Here we have taken the IR cutoffs L = 25,L = 5 and truncation energy Λ = 200.

Discussion
While the Rayleigh-Ritz method is very general and can be applied to any quantum system, supersymmetry provides interesting classes of models with flat directions and nontrivial spectra of metastable states, and in particular precise holographic models that allow for accessing the semi-classical gravity regime in the bulk. To develop tools for analyzing the real time dynamics of strongly coupled SQM and the unitary evolution of black hole microstates is the main motivation for this work.
In applying the Rayleigh-Ritz method to the numerical study of the S-matrix, two key observations were made in this paper. Firstly, we do not need to take the IR cutoff L to be excessively large; on the other hand, it is extremely important to sample over a set of moderately large L values, so that the spectral set {(E n , n)} trace out the curves of scattering phases. This is crucial for determining the metastable resonances and their widths to high accuracy. Secondly, by varying the IR cutoffs in various asymptotic regions, we can extract not just a single scattering phase but the entire S-matrix. The second observation is not necessary for the models considered in this paper, but will be useful for more general models that have several asymptotic regions unrelated by symmetries, such as SQMs that admit noncompact Higgs and Coulomb branches.
Ultimately, we would like to implement our method in the BFSS matrix quantum mechanics [10] and extract the spectrum of black hole microstates in the holographic dual [12][13][14][15]. This is substantially more complicated than the toy supermembrane model and the N = 4 SQM considered here, but they share a number of common features. While in principle the large N limit is required for accessing the semi-classical gravity regime in the bulk, there are indications that this may not be entirely necessary [34][35][36], and we may extract relevant physics even for N = 2 or 3. The SU (2) BFSS MQM involves a priori 2 24 components of the wave function that depend on 27 bosonic coordinates. However, there are only 3 SO(9) × SU (2) invariants formed out of the bosonic coordinates, which will play the role of r,r variables in the N = 4 SQM analyzed in this paper. Furthermore, the 2 24 component fermionic wave function reduces to a few hundred irreducible representations of SO(9) × SU (2) [45,46], and the direct diagonalization of a truncated Hamiltonian could be manageable when restricted to a sector with fixed SO(9) angular momentum. We hope to report on this in the near future.
As discussed for the 1D models in section 3.1, it is important to split the Hamiltonian such that the interaction term is smooth in r andr. One could use a basis that does not leave a smooth interaction term, but the convergence in truncation energy becomes slower, requiring significantly more memory resources and time (especially in the diagonalization). To avoid this, we used a basis of Bessel functions which despite requiring numerical integrations for the matrix elements of the Hamiltonian, leads to exponential convergence in truncation energy. We illustrate this with a sample calculation at small IR cutoff lengths L andL in Figure 15.
We parallel processed the computation of matrix elements using the 'parfor' function of MATLAB's Parallel Computing Toolbox. A typical run of our code used a pool of 30 workers to distribute the process.   sin( nπr L ) sin( mπr L ) adapted to a free Hamiltonian that does not include the rotational potential. While the spectrum computed using the Bessel basis converges exponentially fast, the results computed using the simpler basis have very slow convergence with truncation energy Λ (and the convergence becomes worse with increasing L).
The second step is to diagonalize the matrix. For that end we used MATLAB's function eig() that automatically implements multithreading of available cores, and hence no explicit parallelization is needed. A few typical runs and the parameters of the numerics are shown in Table 1 Table 1: Sample parameters of a typical run for the N = 4 model in the R = 1, j = 1/2 sector with L =L = 10, and for the toy supermembrane model with L = 20, computed with a pool of 30 workers.