Slow Quantum Thermalization and Many-Body Revivals from Mixed Phase Space

The relaxation of few-body quantum systems can strongly depend on the initial state when the system’s semiclassical phase space is mixed; i.e., regions of chaotic motion coexist with regular islands. In recent years, there has been much effort to understand the process of thermalization in strongly interacting quantum systems that often lack an obvious semiclassical limit. The time-dependent variational principle (TDVP) allows one to systematically derive an effective classical (nonlinear) dynamical system by projecting unitary many-body dynamics onto a manifold of weakly entangled variational states. We demonstrate that such dynamical systems generally possess mixed phase space. When TDVP errors are small, the mixed phase space leaves a footprint on the exact dynamics of the quantum model. For example, when the system is initialized in a state belonging to a stable periodic orbit or the surrounding regular region, it exhibits persistent many-body quantum revivals. As a proof of principle, we identify new types of “quantum many-body scars,” i.e., initial states that lead to long-time oscillations in a model of interacting Rydberg atoms in one and two dimensions. Intriguingly, the initial states that give rise to most robust revivals are typically entangled states. On the other hand, even when TDVP errors are large, as in the thermalizing tilted-field Ising model, initializing the system in a regular region of phase space leads to a surprising slowdown of thermalization. Our work establishes TDVP as a method for identifying interacting quantum systems with anomalous dynamics in arbitrary dimensions. Moreover, the mixed phase space classical variational equations allow one to find slowly thermalizing initial conditions in interacting models. Our results shed light on a link between classical and quantum chaos, pointing toward possible extensions of the classical Kolmogorov-Arnold-Moser theorem to quantum systems.

Popular Summary

Intuition suggests that generic, strongly interacting quantum systems establish thermal equilibrium even in isolation. In the process of thermalization, quantum information encoded in the system becomes irrecoverable. Hence, it is highly desirable to understand how thermal equilibrium is established and if it can be delayed. In this work, we provide a framework that allows one to find the slowly thermalizing states in strongly interacting quantum systems.

Using mathematical techniques to map the quantum dynamics of a strongly interacting system into a classical nonlinear dynamical system, we show that the resulting dynamical system generally has a mixed phase space where regions with regular motion coexist with chaotic islands. In addition, we characterize the error encountered in approximating quantum dynamics by a classical dynamical system.

When the error is small, initializing the quantum system in the regular region of phase space can give rise to spectacular many-body revivals, in which the system periodically returns to its original state. This generalizes the mechanism of quantum many-body scars—revivals that happen when there is an unstable periodic trajectory. However, even when the error is large, we show that our framework can find states with parametrically slower thermalization.

This work provides a practical tool to search for slowly thermalizing states of strongly interacting quantum systems. Moreover, our work puts the intriguing connection between classical and quantum chaos in many-body systems on a firm basis.

Figure 1

(a) The state of a quantum system that is short-range entangled and translationally invariant with a unit cell of size $K$ can be described by an MPS with the same size of unit cell. (b) The wave function of the square lattice is approximated by a TTS in the form of a Cayley tree with connectivity $C=4$. The state with $K=2$ unit cell is described by a tree where two different tensors alternate.

Figure 2

(a) Poincaré sections at ${\theta }_2=0$ for a three-site unit cell reveal mixed phase space with large regular regions. A stable periodic orbit that gives rise to robust revivals is denoted by star, while diamond denotes a longer orbit that does not give rise to revivals in exact dynamics. The dashed line represents a particular cut through the regular region, discussed in Sec. 2c. (b) Variational TDVP prediction for the local density of the Rydberg excitation (top) and the entanglement entropy (bottom) for the shortest-period trajectory labeled by star on the Poincaré section in (a). The inset shows the unit-cell choice and entanglement cuts. (c) Quantum dynamics for $|{\mathbb{Z}}_3⟩$ initial state and two entangled initial states marked in (a) obtained from iTEBD with $\chi =900$. Top: The decay of fidelity per spin has regular minima that correspond to revivals of fidelity in a finite-size system. Bottom: Entanglement entropy. The stable “star” trajectory features a better fidelity revival, as well as strongly suppressed entanglement growth. The growth of entanglement limits the iTEBD simulation to shorter times for $|{\mathbb{Z}}_3⟩$, diamond states.

Figure 3

(a) A rapid drop in the largest component of the Fourier spectrum of a local observable, $n_1(t)$, signals the onset of chaotic behavior. (b) The entropy of a region of size $L$ at fixed time $t_0=15$ remains low when the quantum dynamics is launched near the center of the regular region (we use the iTEBD algorithm to simulate quantum dynamics). When the border of the regular region is approached, the value of entanglement rapidly increases with changing the initial state, defined by the value of $c$. When $c$ approaches one, the dynamics becomes insensitive to the changes in the initial state due to rapid thermalization. (c) The derivative of entropy along the trajectory is most sensitive around the boundary of the regular region.

Figure 4

(a) TDVP dynamics on a tree (solid lines) compares favorably with exact dynamics on a $4\times 4$ lattice (dashed lines). Different colors denote different sites in the unit cell. (b) Scaling of the period of the isolated periodic orbit as a function of the connectivity of the TTS. Stars indicate exact results for chain, honeycomb, and square lattices.

Figure 5

(a) Poincaré section for ${\mu }_z=0.225$ at $({\theta }_1^{*}=1.25\pi ,{\dot{\theta }}_1<0)$ and $⟨H⟩=0$ reveals mixed phase space with large regular regions. The Poincaré section is obtained by initializing the system on a grid of points with the step $\delta {\theta }_2=\delta {\varphi }_2=8\times {10}^{-3}$ and numerically integrating the evolution up to a time $t=3000$ using the ndsolve routine in Mathematica. (b) Comparison between TDVP predictions for the dynamics of local observables (${\mu }_z=0.225$) on the periodic trajectory (solid lines) and the exact dynamics of a $4\times 4$ lattice with periodic boundary conditions. (c) Bloch sphere visualization of the individual spin dynamics on the TDVP trajectory. The trajectory on the big meridian corresponds to the case ${\mu }_z=0$. Nonzero chemical potential ${\mu }_z=0.225$, 0.65 causes a progressive tilting of the trajectory and the development of a knot around the south pole.

Figure 6

(a) TDVP quantum leakage presents an upper bound on $f_T$. Despite not being tight, this bound captures the qualitative effect of the applied perturbation. The data illustrate that $f_T$ is system size independent and therefore is well defined for ${\mu }_3\in [-0.6,0.4]$. The inset shows that the Floquet exponent $\lambda$ becomes nonzero only for ${\mu }_3>0.255$, suggesting that there is no direct relation between $\lambda$ and leakage. Finally, for ${\mu }_3\ge 0.45$, the periodic trajectory disappears.

Figure 7

Bipartite entanglement dynamics for different MPS initial states shows a strong dependence of entanglement growth on the initial state. The inset reveals that the amount of entanglement at time $t_f=5$ is correlated with the MPS leakage averaged over TDVP dynamics, ${\mathrm{\Gamma }}_{t_f}/t_f^2$, where the periodic trajectory is shown by star. Using initial conditions on the periodic trajectory results in slow entanglement dynamics (thick dashed line), but there are initial conditions that give even slower entanglement growth. The data are obtained with iTEBD with a maximum truncation error of approximately ${10}^{-10}$.

Figure 8

Illustration of the calculation of the one-point correlation function for $C=3$. The gray colored bonds denote a vectorization of the bonds pointing toward the same direction when the contractions are performed. Graphical definitions of (a) the local tensor $A$; (b) the gauge symmetry [Eq. (b5)], where a hatched circle is introduced as a shorthand notation for a Kronecker delta in virtual index space; (c) the local transfer matrix [Eq. (b6)]; (d) the contraction of the local operator $O_k$ by the local tensor; and (e) the right dominant eigenvalue of the transfer matrix. (f) illustrates the calculation of Eq. (b8) in three steps. In the first step, we show a part of the tree which includes the site $k$, where the local operator acts. The gray colored vertices imply contracted tensors in those sites. In the second step, the gauge symmetry (b) is used, and the whole tree is reduced to a semi-infinite one-dimensional tensor network (the green line is a guide to the eye). In the final step, the infinite product of local transfer matrices is replaced by the right dominant eigenvector of the matrix; see Eq. (b8).

Figure 9

We derive the general bound on the norm of the vector $|\phi (t)⟩$ that represents the mismatch between the exact wave function and its TDVP “image” after the TDVP evolution returns to its initial state. We note that the projection of the exact evolution onto the TDVP manifold does not coincide with the TDVP trajectory at finite times.

Figure 10

The Poincaré section $({\theta }_2=0,{\dot{\theta }}_2<0)$ for the deformed PXP model with ${\mu }_3=-0.3$ shows an overall increase of the regular islands.

Figure 11

TDVP analysis of trajectories in the deformed PXP model with ${\mu }_3=0.3$. (a) The Poincaré section $({\theta }_2=0,{\dot{\theta }}_2<0)$ reveals a smaller size of regular islands in phase space. Black symbols star, diamond, up-pointing triangle, circle, square, and down-pointing triangle show the locations of some periodic orbits. For each of the orbits, in (b) we show the value of the Floquet exponent, quantum leakage ${\mathrm{\Gamma }}_P$, and the dynamics of local observables $n_{1,2,3}$. The physical period $P$ is defined as the time after which the local observables return to their initial values and can be smaller than the period of the orbit for the MPS parameters.

Figure 12

The projection of the Poincaré section onto the ($\theta ,\varphi$) plane in the vicinity of the periodic orbit reveals a regular region of phase space.

Figure 13

Dynamics of entanglement of subsystems of size $l$ at one boundary, $S_l=S(l)/2$, in TFIM for “fast” (solid lines) and “slow” (dashed lines) initial states defined in the text. Different colors correspond to different subsystem sizes that range from $l=1$ (magenta) to $l=9$ (red). Black lines correspond to the entanglement entropy of a half-infinite subsystem. The data are obtained with iTEBD with truncation error $\epsilon <{10}^{-5}$. The inset shows the same data for $l=5,\dots ,9$ with rescaled axes, highlighting that the difference in entanglement spreading cannot be explained by the constant time shift and persists for largest subsystems.

Figure 14

(a) Comparing the expectation values of magnetization per spin between different system sizes suggests the applicability of ETH. The inset confirms that the fluctuations of magnetization between adjacent eigenstates in the energy window $0.19L+[-1,1]$ are suppressed, in accordance with the ETH prediction, as $1/\sqrt{D}$, where $D$ is the Hilbert space dimension. (b) Bipartite entanglement entropy of the eigenstates for $L=20$ shows no outliers near the middle of the spectrum. All data are obtained for the zero-momentum inversion-symmetric sector of TFIM with periodic boundary conditions using exact diagonalization.

Figure 15

The deformation ${\mu }_2=0.24$ and ${\mu }_z=0.4$ effectively destroys the eigenstates with low bipartite entanglement entropy in the pure PXP model. The data are obtained for the $L=28$ chain with periodic boundary conditions in the zero-momentum inversion-symmetric sector.

Figure 16

Entanglement growth in the deformed PXP model [Eq. (d1)] strongly depends on initial conditions. All initial conditions are chosen to have the same energy density $⟨H⟩/L=0$ and the same entanglement entropy. The light blue line corresponds to the periodic trajectory. The inset shows the correlation between the value of entanglement at late time and quantum leakage ${\mathrm{\Gamma }}_{t_f}/t_f^2$, $t_f=10$.