Late Time Quantum Chaos of pure states in the SYK model

In this letter, we study the return amplitude, which is the overlap between the initial state and the time evolved state, in the Sachdev-Ye-Kitaev (SYK) model. Initial states are taken to be product states in a spin basis. We numerically study the return amplitude by exactly diagonalizing the Hamiltonian. We also derive the analytic expression for the return amplitude in random matrix theory. The SYK results agree with the random matrix expectation. We also study the time evolution under the different Hamiltonian that describes the traversable wormholes in projected black holes in the context of holography. The time evolution now depends on the choice of initial product states. The results are again explained by random matrix theory. In the symplectic ensemble cases, we observed an interesting pattern of the return amplitude where they show the second dip, ramp and plateau like behavior.

In this letter, we study the time dependence of pure states in the SYK model. Here we consider a special class of pure states that were first considered in [12,13] and studied further in [14][15][16]. They are simultaneous eigenstates of spin operators that are constructed from Majorana fermions. One interesting physical interpretation of these states is that they are states after projection measurements of maximally entangled states by these spin operators [12,17]. In this context, we can interpret our setup as time evolution after projection measurements. We expect that the time evolution starts flipping these spins under the SYK Hamiltonian and we get a more general superpositions of these product states. In this paper, we consider the return amplitude, which is the square of the fidelity and used in the similar setup in conformal field theory [18,19]. We can also consider the time evolution under deformed Hamiltonians after the measurements. Here we consider a deformation proposed by [12]. They are interpreted as a traversable wormhole protocol to see the inside of black holes [12,20]. This deformed Hamiltonian can also be seen as deformation of an integrable Hamiltonian with degenerate spectrum by the chaotic SYK Hamiltonian.

2.The SYK Model
In the SYK model, we consider even N Majorana fermions ψ i that are obeyed to the anti commutation relation {ψ i , ψ j } = δ ij . The Hamiltonian of the SYK model with q body interactions is given by Here J a1···aq are random couplings with mean J a1···aq J = 0 and variance J 2 a1···aq J = J 2 (q−1)! N q−1 . The q = 4 model is the original SYK model that we mainly focus on in this letter. This system has two important symmetries [6][21]. The first one is the anti unitary symmetry T . This symmetry satisfies T ψ a T −1 = ψ a . The SYK Hamiltonian (1) is invariant under T when q = 0 (mod 4). The other important symmetry is the mod 2 fermion number operator (−1) F with ((−1) F ) 2 = 1. These symmetries can have global anomalies depending on N (mod 8) [21]. These anomalies are the origin to realize all of Gaussian unitary (GUE), orthogonal (GOE) and symplectic (GSE) ensembles in the SYK model [6]. Here we summarize the results [6,22] in the table.I.
Now, we consider the return amplitude. The return amplitude is Here |ψ 0 is a initial state and H is the Hamiltonian of the system. This definition is also applicable to any quantum systems. We take the square of | ψ 0 | e −iHt |ψ 0 |, which is different from the definition in [18,19]. Our choice makes the relation to so called spectral form factor clear, which is studied extensively in the field of quantum chaos and also recently in holography [7,[23][24][25][26]. Because we expect that the late time behavior of the SYK model is governed by random matrices [7], it is good to study first the return amplitude in random matrix theory. Now we want to compute where dH is the Haar measure on the space of L × L Hermitian matrices. We can compute this in random matrix theory using the Haar integrals. The results are . This g(t) is the spectral form factor, which diagnoses the energy level correlations in chaotic systems. We also find that this relation also holds in GSE ensemble by replacing the aevrage to GSE , using the Haar integral [27,28]. In GOE ensemble, the return amplitude depends on the initial state |ψ 0 1 . The equation (3) means that the return amplitude is given essentially by the constant L shift of the spectral form factor. When t = 0, the spectral form factor is simply given by the square of the dimension of the Hilbert space g(0) = L 2 . Under the time evolution, g(t) GUE decreases and hit the minimal value. These regimes are called the slope and the dip [7]. Then, g(t) GUE increases linearly. This linear growth is called the ramp [7] and this reflect the long range eigenvalue correlations in chaotic systems. Finally g(t) GUE saturate the late time value g(∞) GUE = L in sufficiently late time. This is called the plateau [7]. The plateau value generically coincides with the infinite time average: where N E is the degeneracy of each energy level E. This plateau value is much smaller than the initial value g(0) = L 2 but still bigger than O(1). The relation (3) between the return amplitude g p (t) and the spectral form factor g(t) says that the return amplitude also shows the slope, the dip, the ramp and the plateau. The plateau value for the return amplitude is given by g p (t) GUE = 2/(L + 1), which is also much smaller than the initial value g p (0) = 1.
(3) also means that the initial state can |ψ 0 evolves to other states that are orthogonal to |ψ 0 . Let us pick a state |ψ 1 that satisfies ψ 1 |ψ 0 = 0. In the similar way with the return amplitude, we can calculate the overlap | ψ 1 | e −iHt |ψ 0 | 2 . That becomes This amplitude increases under the time evolution from 0, then slightly decreases and finally saturate the late time value 1/(L + 1). Interestingly, this late time value is not equal to the late time value of the return amplitude g p (∞) = 2/(L + 1) but the half of that. On the other hand, we get dU | ψ| U |φ | 2 = 1/L for any choice of |ψ and |φ for Haar random unitary U . This is because the average with lim T →∞ dHe − L 2 TrH 2 F (e −iHt ) is not equivalent to the Haar random unitary average dU F (U ) for functions F on the space of L × L unitary groups [30].
The spectral form factor have a finite temperature generalization g(t; β) = Tr(e −βH−iHt )Tr(e −βH+iHt ). A finite temperature analog of the return amplitude, we consider We can think of this as the return amplitude with the initial state e − β 2 H |ψ 0 . As is the case with the finite temperature spectral form factor [7], when we take the ensemble average , the annealed disorder is not equal to the quenched dis- where E is an ensemble average. In this paper we consider the annealed disorder in which the analytic treatment becomes easy [7,31]. Another motivation to take the annealed version is its similarity with the unnormalized cylinder amplitude in Quantum Field Theories [18,19]. The ensemble average of the finite temperature return amplitude (6) in random matrix theory becomes We call the right hand side of (3) and (7) the shifted spectral form factor. Now we consider the return amplitude in the SYK model. The states we consider is defined as follows. First there are N/2 spin operators in the SYK model defined by S k = −2iψ 2k−1 ψ 2k . These satisfy S 2 k = 1 and eigenvalues of S k are given by s k = ±1. The pure states |B s we consider is defined as common eigenstates of these spin operators: This defines 2 N 2 states, one for each choice of the spins s k . These states form a basis of the SYK Hilbert spaces. By including the Euclidean evolution in the SYK Hamiltonian, we can produce lower energy states: The return amplitude for the SYK model are Currently we do not have any technique to compute analytically (10) for finite N , we numerically computed the finite temperature return amplitude (10)  FIG. 1. These are numerical plots for the N = 14 q = 4 SYK model. We take the disorder average for 1500 samples except for the single sample case. We put β = 1.5. The return amplitude is defined (10) and we choose the state that satisfies S k |B ↑···↑ = |B ↑···↑ . The shifted spectral form factor is the right hand side of (7) for the SYK Hamiltonian. We normalize them so that the initial values become 1. A single sample of gp(t) shows erratic oscillation around the averaged return amplitude at late time and it is not self averaging [32].
observe the slope, the dip, the ramp and the plateau in the return amplitude in the SYK model. The early time decay is almost the same with the spectral form factor. In the large N limit, for any |B s [12] in the leading of 1/N expansion. The early time dependence is captured by the analytic continuation of β → β + it from the leading term in 1/N , we expect the match between them. On the other hand, the ramp and the plateau region we do not expect that because they are non perturbative effects in 1/N expansion [7,30]. The plot shows that they take different value at late time. The late time behavior is expected to be governed by random matrix theory [7]. We also expect this to the return amplitude. To confirm this, we compare the return amplitude with the shifted spectral form factor (7) where the ensemble average g(t; β) GUE is replaced by the SYK coupling average g(t; β) J . We also restrict the Hamiltonian to the fixed (−1) F charge sector in the shifted spectral form factor because only that acts on the state |B s . The plots agree very well and these results also support the random matrix behavior in the late time in the SYK model.

A Deformed Hamiltonian
Next we consider the following "mass term" Hamiltonian [12]: This Hamiltonian is diagonalized in the |B s state basis. Especially, the unique ground state of this Hamiltonian is given by |B s with spin {s k } and energy E (0) 0 = −N/4. By flipping some spins from the ground state |B s , we obtain the whole energy eigenstates. The excited state energy levels are given by There are energy gaps , which is given by , between the bands. Now we consider the Hamiltonian that contains the both of the SYK term and (11): This Hamiltonian was originally proposed to describe the traversable wormhole after projection measurements [12]. We call this H def the deformed Hamiltonian. Here we consider the regime that µ is large and we can treat the SYK term as perturbation. This can be seen as a perturbation of the integrable system with degenerate spectrum by the chaotic Hamiltonian. 2 . We also concentrate on the infinite temperature cases. Because µ is large, exact energy levels from m and form band like structure. By exactly diagonalizing the Hamiltonian, we can study the return amplitude under the deformed Hamiltonian. We show the numerical results in Fig.2. Here we explain the results. First, the ground state of H def is almost given by |B s with m = 0 because the gap between the ground state and  (3) is given by L = dm. We see the kink around the transition time from the ramp to the plateau in the m = 3 case, which reflects GSE statistics [7].
the first excited states is µ, which is sufficiently large and suppresses the mixing with other |B s states. In this case, the return amplitude does not decay and also shows oscillation at early time. The behavior of |B s with m = 0 are described by random matrix theory. According to the perturbation theory of quantum mechanics the first order shift of energy levels are determined by the projection of H SY K onto the degenerate energy levels [36]. Therefore at early time the projection of H SY K on each degenerate energy levels determines the time evolution of |B s . We know the dimension of the projected Hamiltonian from (12). By considering the spectral form factor restricted on the m-th band and its shift by L = d m , we see the very good agreement with the return amplitude. Because the spectral form factor and the return amplitude depend on the level statistics, next we determines the symmetry class of each band. To see this, we need to know the symmetry properties of |B s . These are eigenstates of (−1) F . The anti unitary T flips the spin T S k T −1 = −S k . Therefore, we find that T |B s = e iθs |B −s where e iθs is a phase factor and |B −s is the state that satisfies S k |B −s = −s k |B −s . From this, we find for the projection operator onto the flip number m sector P m = #flip=m |B s B s | satisfies This means that T flips the bands.
Pm − N 4 +m and P i = |ψ i,1 ψ i,1 | + |ψ i,2 ψ i,2 | with two eigenstates of i-th degenerate eigenvalues. We choose the basis that satisfies T |ψ i,1 = |ψ i,2 . From these definitions, we find that T P i T −1 = P i and T Q N . By solving this symmetry constraint in the basis with T = 0 −1 1 0 K , we obtain with a real number x and a complex number a. This means that in a generic Hamiltonian the off diagonal element a enters and the degeneracies are removed. Therefore the degeneracy is removed at the second order of the perturbation. This difference of the order means that the return amplitude (and the spectral form factor) does not see the true degeneracy at early time. We see numerically in fig.2 that they show the first order degeneracy at the first plateau, but after that they show the second slop, dip, ramp, and plateau. The second plateau value is smaller than the first plateau value because (4) means that degeneracies give larger plateau values. To see the level statistics further, we study the distribution of the adjacent gap ratio [6,[37][38][39] for each m-th band. The adjacent gap ratio is defined as r n = En+1−En En−En−1 for an ordered spectrum E n−1 < E n < E n+1 . The distribution of the ratio r n in Poisson statistics is On the other hand, in random matrices, the distribution of the ratio r n becomes [40] For GOE β = 1 and Z β = 27 8 . For GUE β = 2 and 4π 81 √ 3 . For GSE β = 4 and Z β = 4π 729 . r → 0 behavior p(r) ∼ r β represents the level repulsion. We study this gap ratio numerically in N = 16 case and compare with the random matrix case. The m = N 4 = 4th band shows GOE statistics and the m = 3rd band shows GUE statistics. These result agree with our symmetry analysis above. m = N 4 -th bands in N = 4 (mod 8), the degeneracies are The right is the gap ratio r ′′ n = E 2n+2 −E 2n+1 E 2n −E 2n−1 . We take the average over 2000 samples for both cases.
removed at different order of perturbation. Therefore, E 2n − E 2n−1 is of order µ −1 while E 2n+1 − E 2n is of order 1. In this case we expect E 2n+2 − E 2n is determined by the first order perturbation and looks like GSE ensemble when µ is large. We study the gap ratio r ′ n = E2n+2−E2n E2n−E2n−2 numerically and the results show GSE statistics. On the other hand, the gap ratio r ′′ n = E2n+2−E2n+1 E2n−E2n−1 is also an order one quantity because both of the numerator and the denominator are of order µ −1 . These gap E 2n − E 2n−1 come from (15) that is Hermitian. Therefore we expect this ratio shows GUE statistics. We study this gap ratio r ′′ n numerically and they show GUE statistics as expected.
4.Discussion The return amplitude in random matrices is exactly calculated and related to the spectral form factor. Our numerical study also shows that this relation holds in the SYK model. Initially the return amplitude decays but at late time first they grows from the dip and then saturates the plateau value. Because the random matrix behavior is expected to be universal in chaotic systems, we expect that these behaviors are true even in more generic chaotic systems like conformal field theory after a suitable average like a time average. When we deform the SYK Hamiltonian by a mass term, the return amplitude depends on the choice of initial states. If we choose the initial product states to be a ground state of the mass term, the return amplitude does not decay. This deformation prevents the initial state to be scrambled and protects from thermalization. In gravity side, we can interpret this as the disappearance of the black hole horizon and we have access to the black hole interiors. If we flip the m spins from the almost ground state one, return amplitude decays and their behaviors are again explained by random matrix theory. The most interesting case is the m = N 4 cases in N = 4 (mod 8) cases where we see the second dip, ramp and plateau. In early time the return amplitude behave like GSE statistics with 2 degeneracies at each level, but at late time it realizes that the degeneracies are actually removed and finally decays to smaller values. This serves an example of the return amplitude or the spectral form factor with complicated patterns.
Acknowledgement We thank A.Maloney, S. Harrison and T.Takayanagi for useful discussions. We also thank R.Namba and D.Yoshida about discussions for numerical calculations. TN is supported by JSPS fellowships and the Simons Foundation through the It From Qubit collaboration.

A.Explicit Realization of Majorana fermions
In this appendix we give an explicit representation of Majorana fermions. We follow the notation of [30]. We can realize the fermions as the tensor products of the Pauli matrices: Here are the Pauli operators on i-th site and we omit the symbols ⊗ for tensor product in (18). Then, the k-th spin operator S k becomes This confirms that S k is the spin operator that measure the eigenvalues of Z k . In this basis, we can write the |B s state as for the state with S k |B s = s k |B s . We also give the explicit form of anti unitary symmetry T and mod 2 fermion number (−1) F in this basis. The fermion number (−1) F is given by T depends on the N . When N/2 is odd case that corresponds to N = 2, 6 (mod 8) , where K is the antiunitary operator that takes the complex conjugate. When N/2 is even case that corresponds to N = 0, 4 (mod 8) , Using this explicit representation, we can show the symmetry property in Table.I. Though these realization gives a way to see the symmetry property, they are characterized by topological invariants[21] and independent from explicit realization.

B. Haar integrals and the derivation of the return amplitude in Random Matrices
We derive the return amplitude in random matrix theory. The key observation is that the measure is invari This integral have four same unitary matrices. To evaluate this, we need the following integral: Using this integral, we obtain In a similar way, we obtain for orthogonal |ψ 0 and |ψ 1 .
For GSE or GOE cases, the Haar integral for U ∈ U (L) is replaced by the Haar integral for symplectic groups S ∈ Sp(L/2) or orthogonal groups O ∈ O(L). The Haar integral for symplectic groups are given by [27,28] dS where C ij is the antisymmetric invariant. This coupling gives the inner product between the Kramers' pairs. Because this coupling is antisymmetric, the diagonal part v i C ij v j for |ψ 0 = i v i |e i vanishes and they can be ignored in the calculation of the return amplitude. Using this integral, we obtain which is the same expression with the GUE case though the Haar integrals themselves are different.
The Haar integral for orthogonal groups are given by [27,28] dO where C ij is the symmetric coupling and we can choose a basis with C ij = δ ij . Unlike the case of symplectic groups, the diagonal part v i C ij v j for |ψ 0 = i v i |e i does not vanish and the return amplitude depends on states.
C. Symmetry analysis of the Deformed SYK Hamiltonian at 1st order perturbation In this section we study the symmetry class H m in detail that is imposed by T symmetry. Because T relates the m-th band and the ( N 2 − m)-th band, it is sufficient to consider the the symmetry constraint on the following submatrix: where H m = P m H SY K P m and S m = P m H SY K P N 2 −m . Using the symmetry of the SYK Hamiltonian This means S m = 0 for N = 2, 6 (mod 8). For N = 0, 2 (mod 8) cases, T 2 = 1. In these cases, we can choose the basis with

D. Evolution Amplitude in the SYK model
In this section we consider the overlap between time evolved state e −iHt |ψ 0 and states |ψ 1 that is orthogonal to |ψ 0 , which we consider in random matrices in (5), in the SYK model. Especially, we consider the amplitude between |B s with different spins: Here we call this the evolution amplitude. In figure.4, we compare the return amplitude and the evolution amplitude numerically in the SYK model. After the ensemble average, they show that the plateau value in the return amplitude is clearly larger than the evolution amplitude. Even in single sample, the return amplitude looks to oscillate around the larger average value than the the evolution amplitude. As an example of the return amplitude, we consider | B ↑↑↑↑↑↑↑ | e −iH SY K t |B ↑↑↑↑↑↑↑ | 2 , and as an example of the evolution amplitude we consider | B ↓↓↑↑↑↑↑ | e −iH SY K t |B ↑↑↑↑↑↑↑ | 2 . Top: We take the average over 1500 samples for both of the return amplitude and the evolution amplitude. Bottom: The return amplitude and the evolution amplitude with single sample.

E. Time Average of single sample in the SYK model
The spectral form factor is not self averaging [32], and we need some averages to get smooth behaviors. In the SYK model, we take the ensemble average over the coupling J a1···aq . In this appendix, we consider the time average of the return amplitude in a single sample in the SYK model, which is another kind of average. First we consider the infinite time average , which gives the averaged plateau value, in general quantum systems without degeneracy. By decomposing the state |ψ 0 = i c i |E i , we obtain Therefore the time average is not exactly the same with the plateau value in ensemble average (3). If we take the average of (35) over states , it becomes the late time value 2/(L + 1) in the ensemble average. In the SYK model, we consider the following time average: (36) Here we take the time average between 1 2 t < t ′ < 3 2 t around time t, which is taken in [23]. We show the plots of the return amplitude and the spectral form factor with time average in figure.5. Both of the return amplitude and the spectral form factor show the slope, dip, ramp and plateau behavior. This motivate us to expect that the return amplitude in more generic quantum systems like chaotic CFTs shows these structure after the time average. We also compare the shifted spectral form factor (35) where the ensemble average of the spectral form factor is replaced by the time average. As we pointed out, their late time value is not exactly the same, but the behavior shows good agreement on each time. FIG. 5. This is the plot of the return amplitude in the N = 18, q = 4 SYK model with the time average. We choose |B ↑↑↑↑↑↑↑↑↑ as the initial state. We can see the slope, dip, ramp and plateau even in the time average cases. As we mentioned, the plateau value of the time averaged return amplitude have a small deviation from the shifted spectral form factor.