Quantum vs classical dynamics in a spin-boson system: manifestations of spectral correlations and scarring
D Villasenor, S Pilatowsky-Cameo, M A Bastarrachea-Magnani, S Lerma-Hernandez, L F Santos, J G Hirsch
QQuantum vs classical dynamics in a spin-boson system:manifestations of spectral correlations and scarring
David Villase˜nor, Sa´ul Pilatowsky-Cameo, Miguel A. Bastarrachea-Magnani, Sergio Lerma-Hern´andez, Lea F. Santos, and Jorge G. Hirsch Instituto de Ciencias Nucleares, Universidad Nacional Aut´onoma de M´exico,Apdo. Postal 70-543, C.P. 04510 Cd. Mx., Mexico Department of Physics and Astronomy, Aarhus University,Ny Munkegade, DK-8000 Aarhus C, Denmark Facultad de F´ısica, Universidad Veracruzana,Circuito Aguirre Beltr´an s/n, Xalapa, Veracruz 91000, Mexico Department of Physics, Yeshiva University, New York, New York 10016, USA
We investigate the entire evolution of the Dicke model, which is a two-degree-of-freedominteracting spin-boson model of great experimental interest. Our objects of study are thequantum and classical survival probabilities of initial coherent states and the correspondingclassical evolution of the Wigner distribution in phase space. We show that major aspectsof the system are uncovered by analyzing its long-time dynamics, such as whether the ini-tial state is in a regular or chaotic region, in the vicinity of a separatrix, or yet close to anunstable periodic orbit. We demonstrate that a ratio of two between the quantum and clas-sical asymptotic values of the survival probability is a clear indicator of maximal quantumergodicity. In this case, the quantum survival probability develops a correlation hole, that isnonexistent in its classical version and results in a longer equilibration time for the quantumdynamics. These findings are corroborated by our analytical expressions for the survivalprobability and the equilibration time. a r X i v : . [ c ond - m a t . o t h e r] F e b I. INTRODUCTION
Experimental advances in the studies of isolated quantum systems have resulted in ever longercoherence times [1–4], strongly motivating theoretical and experimental interests in long-timequantum dynamics. Recent works on quantum chaos, which showed that the short-time exponen-tial growth of out-of-time correlators (OTOCs) [5–11] is not a universal signature of chaos, butcan emerge also near critical points [12–14], have offered another motivation to switch the atten-tion to long-time dynamics [13, 15]. At large times, the quantum-classical correspondence breaksdown, and effects that are purely quantum come to light. It is at this point that universal spectralcorrelations can get manifested [16–24].At long times, equilibration eventually takes place. In classical mechanics, the mixing prop-erties of chaotic dynamics have provided a fundamental mechanism to explain the equilibrationprocess and the ergodic properties of physical systems. The passage to the quantum domain entailsnew phenomena such as superpositions, quantum interferences, and the effects of the quantizedspectra [16–24] and quantum scars [25–35]. Even though isolated quantum systems are describedby linear equations, one can still talk about equilibration in the sense of saturation of the dynam-ics. That is, the evolution of the observable reaches a point, where it simply fluctuates around itsasymptotic value, and these fluctuations decrease with the system size [36–40].The present work investigates the entire evolution of the Dicke model, with special emphasis onthe long-time dynamics This two-degree-of-freedom interacting spin-boson model was introducedto explain the collective phenomenon of superradiance [41–44], a phenomenon that has been ex-perimentally studied with cold atoms in optical cavities [45–50]. Depending on the parameters andexcitation energy, the model presents regular and chaotic regions [43, 51–53]. It has been used instudies of nonequilibrium dynamics [22, 54–57] and as a paradigm of ultra-strong coupling regimein several systems [58–60]. Experimentally, the model can be studied by means of cavity assistedRaman transitions [48, 61] and with trapped ions [62, 63].The quantity that we use for our studies is the probability to find the initial state at a later time,the so-called survival probability or return probability. We consider initial coherent states, whichenables a direct comparison between the exact quantum evolution and its classical descriptionobtained with the truncated Wigner approximation (TWA) [64–66]. This comparison allows forthe identification of features in the equilibration process that are purely quantum. This includesquantum fluctuations and the effects of correlations between the eigenvalues. Both are associatedwith the discreteness of the spectrum and therefore manifest themselves at long times.We also compare the results for the survival probability with the classical evolution of theWigner distribution in phase space. With this parallel, we gain a deeper understanding of thedifferent features of the quantum dynamics that emerge at different time scales.We place initial coherent states in the regular and chaotic region and show how they can bedistinguished by analyzing the behavior of the survival probability at long times. The quantum andclassical relaxation times do not coincide in the chaotic regime, but they do in the regular case.We also demonstrate how the ratio between the asymptotic values of the quantum and classicalsurvival probabilities tells us about the proximity of the initial state to the separatrix of the modelor to an unstable periodic orbit. In addition, this ratio can be used to gauge the degree of scarring ofthe initial state. One of our strongest results is to show that a ratio equal to two indicates maximalquantum ergodicity. For this limit, we have an analytical expression for the entire evolution of thesurvival probability and for the equilibration time.This paper is organized as follows. The core of the work is in Sec. IV, where a detailed com-parative study of the classical and quantum evolution of the survival probability is presented forthe regular and chaotic regimes. In preparation to this analysis, Sec. II describes the Dicke Hamil-tonian and its classical limit, as well as the initial coherent states studied. In Sec. III, we introducethe survival probability, its relation with the local density of states (LDoS), and its classical ap-proximation using the TWA. Our conclusions are presented in Sec. V.
II. DICKE MODEL
The Dicke model [41] represents a set of N two-level atoms with atomic transition frequency ω interacting with a single mode of a radiation field with frequency ω . It is described by thefollowing Hamiltonian, ˆ H D = ω ˆ a † ˆ a + ω ˆ J z + 2 γ √N ˆ J x (ˆ a † + ˆ a ) , (1)where (cid:126) = 1 , ˆ a ( ˆ a † ) is the bosonic annihilation (creation) operator of the field mode, ˆ J x,y,z = (cid:80) N k =1 ˆ σ kx,y,z are collective pseudo-spin operators given by the sum of Pauli matrices ˆ σ x,y,z , and γ is the spin-boson interaction strength. When γ reaches a critical value γ c = √ ωω / , a second-order quantum phase transition takes place [42, 43]. There, the system goes from the normal phase( γ < γ c ), where the ground state is characterized by all atoms in their ground state and no photons,to the the superradiant phase ( γ > γ c ), where the ground state has a macroscopic population ofphotons and excited atoms.The eigenvalues j ( j + 1) of the total spin operator ˆ J = ˆ J x + ˆ J y + ˆ J z determine the differentinvariant subspaces. We work with the maximum value j = N / , which defines a symmetricatomic subspace that includes the ground state. The Hamiltonian ˆ H D commutes also with theparity operator ˆΠ = e iπ ˆΛ , where ˆΛ = ˆ a † ˆ a + ˆ J z + j ˆ1 represents the total number of excitations witheigenvalues Λ = n + m + j . Here, n indicates the number of photons and m + j is the number ofexcited atoms, m being the eigenvalue of the operator ˆ J z . A. Classical Limit of Dicke Hamiltonian
The corresponding classical Hamiltonian is obtained using Glauber coherent states for thebosonic sectors [52, 53, 67, 68], | q, p (cid:105) = e − ( j/ ( q + p ) e (cid:104) √ j/ q + ip ) (cid:105) ˆ a † | (cid:105) , (2)and Bloch coherent states for the pseudo-spin sectors, | Q, P (cid:105) = (cid:32) − Z (cid:33) j e [ ( Q + iP ) / √ − Z ] ˆ J + | j, − j (cid:105) , (3)where Z = Q + P . The canonical j -independent variables ( q, p ) and ( Q, P ) are associated withthe photonic and atomic degrees of freedom, respectively. The ket | (cid:105) denotes the photon vacuumand | j, − j (cid:105) , the state with all atoms in their ground state. The rescaled classical Hamiltonian h cl is given by (see Appendix A for details), h cl ≡ (cid:104) q, p ; Q, P | ˆ H D | q, p ; Q, P (cid:105) j (4) = ω q + p ) + ω Z + 2 γQq (cid:114) − Z − ω . The rescaled classical Hamiltonian and its four-dimensional phase space M are independent of j .This is equivalent to working with an effective Planck constant (cid:126) eff = 1 /j [69].Depending on the parameters and energies, the model displays chaotic or regular dynamics.As a case study, we choose ω = ω , j = 100 , and a coupling strength in the superradiant phase γ = 2 γ c . For these parameters, the normalized energy (cid:15) = Eω j (5)of the ground state is (cid:15) GS = − . . The dynamics is regular from (cid:15) GS up to (cid:15) ≈ − . and chaoticabove this point [53]. B. Initial States
Our initial states are the Glauber-Bloch coherent states | q, p ; Q, P (cid:105) . They allow for a directconnection between the quantum states and the points in the classical phase space of h cl , andalso for a relatively simple calculation of the Wigner distribution needed for the evaluation of theclassical dynamics.To select the initial states, we restrict ourselves to the hyperplane p = 0 and solve the second-degree equation h cl ( q, p, Q, P ) = ω (cid:15) in q . This equation has two solutions, q − and q + , with q − ≤ q + . Our initial states are centered at ( q, p, Q, P ) = ( q + , , Q , P ) . Two of them have energyin the regular region, (cid:15) R = − . , and two of them have the energy shell fully covered by chaotictrajectories, (cid:15) C = − . . The criteria for our choices and the specific values of Q and P aredescribed below.
1. Regular Regime: (cid:15) R = − . Poincar´e sections at p = 0 of the classical Hamiltonian at energy (cid:15) R = − . are shown inFig. 1 (a). The phase space splits in three regions of regular trajectories, whose borders are definedby the separatrix marked with large black dots. We choose state I at ( Q , P ) = (1 , and indicateit with a blue dot in Fig. 1 (a). State II is close to the separatrix, ( Q , P ) = (1 . , , and is shownwith a red dot in the same figure.The structure of the phase space reflects the quasi-conserved quantities of the Dicke Hamil-tonian at low energies [70]. They can be identified by means of an adiabatic approximation,where the dynamics separates into two parts, one fast-evolving mode and a slow one, effectivelydecoupling the boson and pseudo-spin dynamics [71]. For the parameters considered here, a quasi-constant of motion is given by the nutation angle of the pseudo-spin, which precesses fast aroundan axis whose direction oscillates slowly in a way dictated by the slow bosonic variables. StateI is at the center of the slow-boson regular region located between . ≤ Q ≤ . , where theboson modes are expressed by very small nutation angles and large amplitudes of the precessionaxis’ oscillations. The rightmost region in Fig. 1 (a) corresponds to modes of the fast pseudo-spindegree of freedom, where the dynamics has maximal nutation angles and very small oscillationsof the precession axis instead.Quasi-integral behaviors are destroyed in non-linear systems due to the growth and consequentoverlap of nonlinear resonances between the integrable modes. These non-linear resonances ariselocally in the phase space [72] and are separated by a separatrix, around which a stochastic layeris formed [73]. Chaos emerges from it, destroying the remnant regular surfaces and opening the IIIIV (a) (b)(c) (d)
FIG. 1. Poincar´e sections ( p = 0 ) for the rescaled classical Hamiltonian h cl at energies (cid:15) R = − . (a) and (cid:15) C = − . (c). In panel (a), large black points mark the separatrix of the regular modes, the blue pointindicates the center of the coherent state I with coordinates ( Q , P ) = (1 , , and the red point indicatesthe state II centered close to the separatrix at ( Q , P ) = (1 . , . The closed curves encircling these twopoints represent the spreading of the coherent state wave function up to e − . In (b) and (d): Participationratio [Eq. (6)] of the coherent states centered in each point of the Poincar´e surface and projected in theHamiltonian eigenbasis. In panel (d), the cyan point indicates the initial coherent state III centered at ( Q , P ) = (1 . , , which has a low PR ( P R = 1066 ) and the red point marks the initial state IV centeredat ( Q , P ) = ( − . , . , which has a high PR ( P R = 5743 ). phase space to the diffusion of trajectories, which start wandering through the whole region ofthe resonance overlap. In the Poincar´e sections of Fig. 1 (a), one can see non-linear resonancesbetween the adiabatic modes described above, specifically in the moon-shaped region of regulartrajectories rotating around the point ( Q, P ) (cid:39) (1 . , . State II is in this region, but very close tothe separatrix.In Fig. 1 (b), we plot the participation ratio P R of the coherent states centered in each point ofthe Poincar´e surface and projected into the energy eigenstates, P R = 1 (cid:80) k | c k | , (6)where c k = (cid:104) E k | q + , Q , P (cid:105) and ˆ H | E k (cid:105) = E k | E k (cid:105) . The participation ratio is strongly correlatedwith the underlying classical dynamics [53, 74]. For instance, the stochastic layer appearing inthe classical dynamics around the separatrix is associated with larger values of P R , which indicatestates that are more delocalized in the energy eigenbasis.
2. Chaotic Regime: (cid:15) C = − . Poincar´e sections for the selected energy (cid:15) C = − . are shown in Fig. 1 (c). They reveal aregion of hard chaos, where chaotic trajectories, all with the same positive Lyapunov exponent,densely fill the whole phase space. From these Poincar´e sections, no particular region can beidentified, but the participation-ratio map in Fig. 1 (d) provides a richer picture, with coherentstates showing different levels of spreading in the energy eigenbasis. To analyze how the structureof the initial states affects the dynamics and equilibration, we select two initial states. State III,indicated with the cyan point in Fig. 1 (d), is located in the region with small values of P R , at ( Q , P ) = (1 . , . State IV, marked with a red point in Fig. 1 (d), is in the region of large valuesof P R , at ( Q , P ) = ( − . , . . These two states are representative of the typical cases usedin the studies of the dynamics of chaotic regions in [28]. III. SURVIVAL PROBABILITY
The survival probability, S P , is the probability of finding an evolved quantum state back in theinitial state | Ψ(0) (cid:105) = (cid:80) k c k | E k (cid:105) , S P ( t ) = |(cid:104) Ψ(0) | Ψ( t ) (cid:105)| = (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:88) k | c k | e − iE k t (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)(cid:12) . (7)By introducing the local density of states (LDoS) or strength function, that is the energy distribu-tion weighted by the components | c k | of the initial state, G ( E ) = (cid:88) k | c k | δ ( E − E k ) , (8)we can also write the survival probability as the squared norm of the Fourier transform of G ( E ) as S P ( t ) = (cid:12)(cid:12)(cid:12)(cid:12)(cid:90) d E G ( E ) e − iEt (cid:12)(cid:12)(cid:12)(cid:12) . (9)The evolution of the survival probability shows different behaviors at different time scales [20,21, 75, 76]. By smoothing the LDoS, one gets insight on how its structure affects the dynamics atdifferent times. The smoothing is done through a finite resolution function, given by ρ T ( E ) = (cid:88) k | c k | Υ T ( E − E k ) , (10)where Υ T ( E − E k ) = ( T /π ) sinc [( E − E k ) T ] , and sinc ( x ) = sin( x ) /x . The time resolution T = π/ ∆ reflects aspects of the LDoS that are of order ∆ in energy. A. Time Scales of the Survival Probability
Figures 2 (a) and (b) and Figs. 3 (a), (b), (d), (e), (g), (h) show the smoothed LDoS, ρ T ( E ) ,of the initial coherent states described in Sec. II B for different time resolutions [Eq. (10)]. InFigs. 2 (c) and (d) and Figs. 3 (j) and (k), we show the infinite-time resolution LDoS, G ( E ) . ω T = . ω T → ∞ FIG. 2. Smoothed local density of states ρ T ( E ) (a,b) and infinite-time resolution LDoS G ( E ) (c,d) forthe coherent states I (a,c) and II (b,d) in the regular region [see Fig. 1 (a)]. In (a) and (b): The timeresolution is ω T = 0 . [ ∆ / ( ω j ) = 0 . ]. Solid blue lines are the smoothed energy profiles ρ T ( E ) ,given by Eq. (10). Red solid lines are the Gaussian energy profiles ρ ( E ) , given by Eq. (11). In (c) and (d): ω T → ∞ [ ∆ / ( ω j ) → ]. Blue dots represent the numerical components of the LDoS. Both profiles arecentered at (cid:15) R = − . with standard deviation σ/ ( ω j ) = 0 . for state I and σ/ ( ω j ) = 0 . for stateII.
1. Survival Probability: Initial Decay
For all cases, up to a time resolution T ∗ that depends on the state, we find a very good Gaussiandistribution for the coherent states profiles, ρ ( E ) = 1 √ πσ e − ( E − E c ) / (2 σ ) . (11)The distributions are centered at the energy E c of the initial state and have width given by theenergy standard deviation σ , which can be calculated numerically or even analytically [77–79].According to Eq. (9), the Gaussian envelope of the LDoS leads to the initial Gaussian decay ofthe survival probability, S P ( t ) = e − σ t , (12)which is consistent with the universal quadratic behavior, S P ( t (cid:28) σ − ) ≈ − σ t , for very shorttimes.In addition to the initial coherent states, we show in Figs. 3 (c), (f), (i), and (l) the smoothedand infinite-time resolution LDoS of a random initial state. Motivated by the high level of delocal-ization of the initial coherent state IV in the chaotic regime, the random state is built with energycomponents generated randomly around the Gaussian distribution of state IV as | c ( r ) k | = r k ρ ( E k ) A ν ( E k ) . (13)Above, ρ ( E k ) is obtained from Eq. (11), ν ( E k ) is given by the equation for the density of states(DoS) provided in the Appendix A, r k are random numbers from an exponential distribution P ( r ) = λe − λr , and A = (cid:80) q r q ρ ( E q ) /ν ( E q ) is a normalization constant. The division by the (g) (h) (i)(d) (e) (f(a) (b) (c))(l)(j) (k) ω T = . ω T = ω T = ω T → ∞ (m) FIG. 3. Smoothed local density of states ρ T ( E ) (a)-(i) and infinite-time resolution LDoS G ( E ) (j)-(l) forthe coherent states III (a,d,g,j) and IV (b,e,h,k) in the chaotic region [see Fig. 1 (d)], and a random state(c,f,i,l) (see text). In (a)-(c): The time resolution is ω T = 0 . [ ∆ / ( ω j ) = 0 . ]. In (d)-(f): ω T = 2 [ ∆ / ( ω j ) = 0 . ]. In (g)-(i): ω T = 3 [ ∆ / ( ω j ) = 0 . ]. Solid blue lines are the smoothed energyprofiles ρ T ( E ) [Eq. (10)], and red solid lines are the Gaussian energy profiles ρ ( E ) [Eq. (11)]. In (j)-(l): ω T → ∞ [ ∆ / ( ω j ) → ]. Blue dots are the numerical components of the LDoS. All profiles are centeredat (cid:15) C = − . with standard deviation σ/ ( ω j ) = 0 . for state III and σ/ ( ω j ) = 0 . for both thestate IV and the random state. Panel (m): Histogram of the numbers r CS k (light gray bars) obtained withEq. (14) for the coherent state IV. The red dashed curve represents the best fit to an exponential distribution P ( r ) = λe − λr , where λ = 0 . . DoS is done to compensate for the different energy densities and to guarantee a smooth envelop-ing distribution ρ ( E ) . The exponential distribution for generating the random numbers is usedbecause if we evaluate the numbers r CS k = A ν ( E k ) | c k | ρ ( E k ) , (14)for the components c k of the coherent state IV in the energy eigenbasis, we obtain the histogramin Fig. 3 (m), which is very well fitted with an exponential distribution.By comparing the smoothed LDoS of the state III in Fig. 3 (d), the state IV in Fig. 3 (e), andthe random state in Fig. 3 (f), all of them with time resolution ω T = 2 , one notices that thesmoothed LDoS of the random state already deviates from a Gaussian, which must affect its short-time dynamics. While the initial coherent states III and IV are expected to remain on the Gaussiandecay at this time scale, the random state should already diverge from it, as we indeed confirm inSec. IV. For an even higher time resolution, as ω T = 3 in Figs. 3 (g)-(i), the smoothed LDoS ofstate III also deviates from the Gaussian, and so does its Gaussian decay (see Sec. IV).
2. Survival Probability: Intermediate Times
The behavior of the survival probability at intermediate times depends on the initial coher-ent state, as one can anticipate from the infinite-time resolution LDoS in Figs. 2 (c) and (d) andFigs. 3 (j) and (k). For the regular coherent state I [Fig. 2 (c)], very few components have nonegligible values and the LDoS is well described by a Gaussian distribution. For the regular co-herent state II located close to the separatrix [Fig. 2 (d)], the number of participating energy levelsis larger, but they are still organized according to a set of different Gaussians with different am-plitudes and centers [79]. In contrast, the components of the chaotic coherent states in Fig. 3have a very different structure. For the coherent state III, with small P R [Fig. 3 (j)], the compo-nents are bunched around some specific energy levels, while the coherent state IV, with a large P R [Fig. 3 (k)], counts with the participation of most components. The consequences of thesedistributions to the evolution of the survival probability are discussed in Sec. IV.
3. Survival Probability: Asymptotic Values
The asymptotic value of the survival probability, S ∞ P = (cid:104) S P ( t ) (cid:105) t →∞ = lim t →∞ t (cid:90) t d t (cid:48) S P ( t (cid:48) ) , (15)can be derived from S P ( t ) = (cid:88) k (cid:54) = l | c l | | c k | e − i ( E k − E l ) t + (cid:88) k | c k | . (16)In the absence of energy degeneracies, the first term on the right-hand-side of the equation cancelsout on average, so S ∞ P = (cid:88) k | c k | , (17)which is the inverse of the participation ratio shown in Eq. (6). However, when the the energylevels have degeneracies of degree d k (that is, | E k , m (cid:105) with m = 1 , .., d k ), the first term in Eq. (16)contributes with additional terms to the asymptotic value, which is now given by S ∞ P = (cid:88) E k d k (cid:88) m =1 | c k,m | , (18)where c k,m are the components in the degenerate space of E k . In the superradiant phase of theDicke model, for the energy range going from the ground state to (cid:15) = − , the spontaneous break-ing of the parity symmetry [80] produces an energy spectrum with two-fold degeneracies . For Because of tunneling effects in the energy region with (cid:15) < − , parity partners are not exactly degenerated, but weverified that their energy differences are smaller than the numerical precision of our numerical calculations. (cid:15) C is well above − , the influence of thedegeneracies is rather marginal. B. Classical Limit of the Survival Probability
To find the classical limit of the survival probability, we use the Wigner formalism [81]. Detailsare given in Appendix B and Appendix C. The basic idea is to use the overlap property betweentwo arbitrary quantum states | A (cid:105) and | B (cid:105) to write the survival probability in terms of the Wignerfunction W [82], (cid:12)(cid:12) (cid:104) A | B (cid:105) (cid:12)(cid:12) = (2 π (cid:126) ) d (cid:90) d u W A ( u ) W B ( u ) , (19)where d represents the degrees of freedom of the system. We scale the Wigner function of ourGlauber-Bloch initial coherent states and work with w ( u , t ) = (1 /j ) W ( √ jq, √ jp, Q, P, t ) , where u = ( q, p, Q, P ) is a point in the j -scaled phase space M . This gives S P ( t ) = (cid:18) πj (cid:19) (cid:90) M d u w ( u , w ( u , t ) . (20)To finally obtain the classical limit, we use the TWA. As shown in Appendix C, the short-timedependence of the Wigner function can be written in terms of the Hamiltonian flow ϕ t : M → M .One finds that w ( u , t ) = w ( ϕ − t ( u ) , for short times, so the classical survival probability can bedefined as S P ( t ) = (cid:18) πj (cid:19) (cid:90) M d u w ( u ) w ( ϕ − t ( u )) , (21)where w ( u ) = w ( u , t = 0) . This quantity is numerically constructed through a Monte Carlomethod. Its asymptotic value, S ∞ P = (cid:10) S P ( t ) (cid:11) t →∞ , is obtained with Eq. (15). IV. CLASSICAL AND QUANTUM DYNAMICS
By comparing the behaviors of the quantum and classical survival probabilities for the fourinitial coherent states selected in Sec. II B, we are able:– To understand why S P ( t ) for coherent initial states reaches values close to zero at short times.– To track down the origins of the different long-time behaviors of initial states that have thesame energy.– To identify properties of the survival probability that are purely quantum and which are seenonly at long times.– To explain why the saturation times of the quantum and classical dynamics in the chaoticregime do not coincide. A. Regular Region
We analyze first the quantum S P ( t ) and the classical S P ( t ) for the initial coherent states I andII [Fig. 1 (a)], which are in the regular regime ( (cid:15) R = − . ).1
1. Initial State I: Center of the Regular Mode
An analytical equation is available for the survival probability of the initial coherent state I,located at the center of the slow-boson regular region [79]. It is given by S P ( t ) ≈ ω σ √ π (cid:88) n =1 e − n (cid:32) ω σ + t t D (cid:33) cos( nω t ) + S ∞ P , (22)where the index n denotes the distance between the eigenenergies of the levels with no negligiblecomponents, ω is related with the spacing between neighboring levels, σ is the width of the LDoS, t D = ω / ( σ | e | ) is the decay time, with e being the anharmonicity of the spectrum, and S ∞ P = ω σ √ π . (23)The time t D comes from the Gaussian decay of the slowest component n = 1 .In the bottom panel of Fig. 4, we compare the full quantum evolution of S P ( t ) with the ana-lytical expression in Eq. (22) and the classical S P ( t ) . The latter two show perfect agreement andthey are indistinguishable from the quantum result up to the decay time t D . Soon after this point, S P ( t ) shows quantum fluctuations around S ∞ P that emerge due to the discreteness of the spectrumand which the analytical expression in Eq. (22) and S P ( t ) are unable to reproduce.The remarkable agreement between the quantum and classical results allows for a more intu-itive classical interpretation of the behavior of the survival probability. In the top panels of Fig. 4,we show the classical evolution of the Wigner distribution in the phase space. The top row givesthe distribution projected onto the Q - P plane and the bottom row, the distribution on the q - p plane.Each column corresponds to a specific instant of time, as indicated. The classical survival prob-ability measures the percentage of the Wigner distribution, in color, that is inside of the startingregion occupied at t = 0 (See Appendix D for a detailed explanation). This starting region ismarked with a black outline which is inside a bigger one marking the available phase space. Theparallel between the evolution of the survival probability and of the Wigner distribution in thephase space goes as follows (animations are available in the Supplemental Material (SM) [83]):– At ω t = 0 (first column), the distribution is entirely within the starting region. Concurrently,the value of the survival probability in the bottom panel is one.– As the distribution moves out of the starting region, the survival probability follows a Gaus-sian decay. At ω t = 2 (second column), the distribution is effectively outside the initial regionand the survival probability becomes close to zero .– After a full period, the distribution comes back to the starting region at ω t = 6 . (thirdcolumn) and we have the first revival of S P ( t ) .– After several periods, the distribution further spreads over the classical orbit and the amplitudeof the revivals of the survival probability decreases, as for example at ω t = 100 (fourth column).– Classically, at ω t = 5000 (fifth column), the Wigner distribution ends up filling homoge-neously the region covered by the trajectories of the phase-space points of the initial distribution.At this point the classical S P ( t ) reaches the constant value S ∞ P , while S P ( t ) fluctuates around theasymptotic value S ∞ P ≈ S ∞ P . Note that even though in the Q - P projection the distribution appears to be still partially inside the starting region,this is merely an artifact caused by the projection. In the q - p projection, the distribution is already clearly outsidethe starting region. This should be kept in mind when interpreting these figures, both projections are important andcannot be treated independently. FIG. 4. Top panels: Classical evolution of the Wigner distribution projected onto the Q - P plane (top row)and q - p plane (bottom row) for the initial coherent state I in Fig. 1 (a). The initial distribution (smallblack circle) is inside the available phase space (big black circle). Each column represents an instant oftime, as indicated. (Animations can be found in the SM [83].) Bottom large panel: Quantum survivalprobability (light gray solid line), analytical expression from Eq. (22) (orange dashed line), and classicalsurvival probability (dark gray solid line) for the initial coherent state I. The horizontal red dotted line showsthe asymptotic value S ∞ P [Eq. (23)]. The decay time t D is indicated with a vertical line. The analyticalexpression uses the parameters from the numerical data ( σ, ω , e ) = (4 . , . , − . × − ) , whichgives ω t D = 403 . .
2. Initial State II: Close to the Separatrix
The LDoS of the initial coherent state II presents different Gaussian distributions with differentamplitudes, centers and widths [Fig. 2 (d)]. This indicates that this state activates the two adia-batic modes (both the slow boson mode and the fast pseudo-spin mode), as well as the non-linearresonances between them. This is indeed expected, since state II is very close to the separatrix[Fig. 1 (a)]. An analytical expression can also be obtained for this case by generalizing Eq. (22),where in addition to contributions from Gaussians, we also need to take into account interferencesbetween them [79].In the bottom panel of Fig. 5, we show the quantum and classical survival probability for thestate II. The agreement between the two is good, but contrary to Fig. 4, the values of S P ( t ) are3slightly smaller than those for S P ( t ) . This is better seen with the temporal averages of the curves,shown in Fig. 5 in blue for S P ( t ) and in red for S P ( t ) . FIG. 5. Top panels: Classical evolution of the Wigner distribution projected onto the Q - P plane (toprow) and q - p plane (bottom row) for the initial coherent state II from Fig. 1 (a). The initial distribution(small black circle) is inside the available phase space (big black circle). Each column represents a time,as indicated. (Animations can be found in the SM [83].) Bottom panel: Quantum survival probability(light gray solid line), its temporal average (blue solid line), classical survival probability (dark gray solidline), and its temporal average (red solid line) for the initial coherent state II. The temporal averages arecomputed for temporal windows of constant size in the logarithmic scale. The horizontal blue dotted lineis the quantum asymptotic value S ∞ P and the red dashed line corresponds to the classical asymptotic value S ∞ P . To better understand the small differences between S P ( t ) and S P ( t ) , we calculate their asymp-totic values for initial coherent states ( q, p, Q, P ) = ( q + , , Q i , , where Q i covers all possiblepoints of the selected Poincar´e surface in Fig. 1 (a). We plot the values of S ∞ P and S ∞ P in Fig. 6.One sees that S ∞ P is smaller than S ∞ P near the separatrix regions ( Q = 0 . , . , and . ),while the two values get closer near the centers of the adiabatic-modes regions ( Q = 1 . , . ) andof the center of the non-linear resonances ( Q = 1 . ). Since a coherent state located in the sepa-ratrix has a Wigner distribution defined in regions classically not connected by the trajectories, weattribute the small difference between S P ( t ) and S P ( t ) to a dynamic tunneling effect that takesplace in the quantum regime [84], but is absent in the classical limit.4 FIG. 6. Quantum asymptotic value S ∞ P (blue solid line) and classical asymptotic value S ∞ P (red solid line)for a set of initial coherent states in the same regular Poincar´e section of Fig. 1 (a) with P = 0 and Q covering the whole allowed interval. The vertical lines marked with S indicate the separatrix at differentcrossing points, the one marked with C I indicates the center of an adiabatic-mode region and the one markedwith C II indicates the center of the non-linear resonances. Even though the quantum-classical agreement for state II is not exact, the comparison betweenthe results for the survival probability in the bottom panel and the evolution of the Wigner dis-tribution in the top panels is similar to that presented in Fig. 4 and explains specific features of S P ( t ) , such as the revival at ω t = 6 . . At very long times, such as ω t = 5000 in the figure, theWigner distribution fills homogeneously the classical region that is covered by the trajectories ofthe initial distribution. Since this region is larger than in the case of state I, the asymptotic valuein Fig. 5 is one order of magnitude smaller than in Fig. 4. B. Chaotic Region
We now study the quantum S P ( t ) and the classical S P ( t ) for the initial coherent states III andIV [Fig. 1 (d)], which are in the chaotic regime ( (cid:15) R = − . ). We start the analysis with state III,which has a smaller number of contributing energy eigenbasis (smaller P R ) than state IV.
1. Initial State III: Low P R Despite being in the chaotic region, initial coherent states such as state III count with a relativelysmall number of contributing energy eigenstates, as suggested by the low value of P R in Fig. 1 (d).As we show below, these states lead to large recurrences at intermediate times and, contrary to theregular case, the equilibration time for the classical and quantum dynamics no longer coincide.The equilibration time for the classical dynamics is now longer than the equilibration time for thequantum S P ( t ) .The quantum and classical survival probabilities of the initial coherent state III are shown inthe bottom panel of Fig. 7. There is excellent agreement from ω t = 0 up to times beyond thedecaying recurrences. Initially both curves decay on a Gaussian to values close to zero and thenrevivals appear. These two features can be well understood by studying the classical evolutionof the Wigner distribution shown in the top panels of Fig. 7. Analogously to the discussionsabout Fig. 4 and Fig. 5, as the originally localized Wigner distribution at ω t = 0 (first column)5 FIG. 7. Top panels: Classical evolution of the Wigner distribution projected onto the Q - P plane (top row)and q - p plane (bottom row) for the initial coherent state III from Fig. 1 (d). The initial distribution (smallblack circle) is inside the available phase space. Each column represents a time, as indicated. (Animationscan be found in the SM [83].) The black arrows indicate an enhancement of the distribution around anunstable periodic orbit (see text). Bottom panel: Quantum survival probability (light gray solid line), itstemporal average (blue solid line), classical survival probability (dark gray solid line), and its temporalaverage (red solid line) for the initial coherent state III. The temporal averages are computed for temporalwindows of constant size in the logarithmic scale. The horizontal blue dotted line is the quantum asymptoticvalue S ∞ P and the horizontal red dashed line corresponds to the classical asymptotic value S ∞ P . The insetshows the classical approximation (TWA) of the survival probability in a semi-log scale, which makes clearthat the oscillations are periodic. moves outside the starting region at ω t = 2 (second column), the survival probability becomeseffectively zero. The recurrences are connected with the return of the distribution to the startingregion, such as at ω t = 5 . . (third column). These recurrences are periodic, as confirmed withthe inset in the bottom panel of Fig. 7.Soon after the last revival, the quantum-classical correspondence breaks down. At this point thequantum S P ( t ) reaches an asymptotic value, around which one finds large quantum fluctuations,while the classical S P ( t ) continues decreasing and equilibrates at a longer time. The asymptoticvalue of the classical survival probability is obtained using the TWA and the ergodic hypothesis6(see Appendix E), S ∞ P = 12 √ πσν c . (24)The asymptotic result S ∞ P of the quantum survival probability is more than twice this value.Large S ∞ P , or equivalently small P R , in the chaotic region is usually associated with the phe-nomenon of scarring, which is indeed the case here. The classical Wigner distribution at ω t = 400 (fifth column in the top panel of Fig. 7) is enhanced in a small closed region indicated with blackarrows in the figure. This reveals the presence of unstable periodic orbits of relatively short period.These orbits are responsible for the short-time periodic revivals of the quantum S P ( t ) and classical S P ( t ) .For longer times, the periodic orbits produce opposing effects in the quantum and classicalregimes. Classically, we observe a slow decay of S P ( t ) towards its asymptotic value and thus along classical equilibration time if compared to the quantum case. Recurrences imply that the dy-namics revisits part of the phase space that was initially covered instead of exploring new regions,which slows down the full spread over the phase space [28, 29]. In the quantum evolution, the scar-ring decreases the equilibration time. As discussed in Ref. [27] and as can be seen in Fig. 3 (g),the proximity of an initial state to unstable periodic orbits of short period produces a particularstructure in the distribution of the energy eigenbasis components. Many of these components arevery small, while the large ones are organized in bunches separate in energy by ∆ E s ≈ π/τ ,where τ is the period of the classical unstable periodic orbit. These highly populated eigenstatesare scarred , meaning that they are concentrated in the phase space around the unstable periodicorbits [28, 29]. Because of this concentration, the phase space available for the quantum evolu-tion of the initial coherent state is effectively shrunk, resulting in an smaller equilibration time for S P ( t ) than for S P ( t ) .
2. Initial State IV: High P R The behavior of the survival probability presented in the previous subsection is not general.Most initial coherent states in the chaotic region are similar to state IV, being highly delocalizedin the energy eigenbasis. In fact, as discussed in Sec. III A, the energy distribution of state IV iscomparable to that of a random state. This latter state is very useful, because one can derive ananalytical expression for its survival probability. The expression for the average over an ensembleof initial random states is given by (see [22] and Appendix F for details) S ( r ) P ( t ) = 1 − S ( r ) , ∞ P η − (cid:34) ηe − σ t − b (cid:18) Dt π (cid:19)(cid:35) + S ( r ) , ∞ P , (25)where η = 2 √ πσν c and ν c = ν ( E c ) is the DoS [see Eq. (A4) in Appendix A] evaluated inthe center of the energy profile E c , σ is the width of the LDoS, the b function is the Gaussianorthogonal ensemble (GOE) two-level form factor studied in random matrix theory [85], and D =2 /ν c is the mean level spacing of the correlated eigenvalues. By comparing the quantum survivalprobability with S ( r ) P ( t ) and the classical S P ( t ) , we can explain the different behaviors of S P ( t ) found at different time scales.In the bottom panel of Fig. 8, we show the quantum S P ( t ) of state IV (light grey solid line),its time average (blue solid line), the classical S P ( t ) (dark grey solid line), its time average (redsolid line), and the time average of the analytical expression S ( r ) P ( t ) (orange solid line). At short7times, S P ( t ) and S P ( t ) overlap. There is also perfect agreement between the Gaussian decay ofthe temporal average of S P ( t ) , S P ( t ) and S ( r ) P ( t ) , as expected from the great similarity betweenthe smoothed LDoS at low time resolution ( ω T = 0 . ) of the coherent and the random state[Fig. 3 (b) and (c)]. However, the curves for S P ( t ) and S ( r ) P ( t ) diverge after ω t ∼ . . One seesthat S P ( t ) remains on the Gaussian decay, while S ( r ) P ( t ) reaches a plateau. This divergence canalso be understood from the smoothed LDoS, but now at higher time resolution ( ω T = 2 , ) [seeFigs. 3 (e) and (f), Figs. 3 (h) and (i)]: while the smoothed LDoS for the state IV is still a goodGaussian, that of the random state shows deviations. This discrepancy implies that the componentsof the coherent state IV are not exactly random.The quantum and classical survival probabilities continue their Gaussian decay for ω t > . and reach values close to zero, which can once again be understood from the analysis of theclassical evolution of the Wigner distribution shown in the top panels of Fig. 8. As the initiallylocalized Wigner distribution at ω t = 0 (first column) moves out of the starting region at ω t = 2 and ω t = 8 . (second and third columns), the survival probability becomes effectively zero.Beyond these times and contrary to what one sees for the regular regime and for state III, thequantum-classical correspondence breaks down before the quantum equilibration, as we explainnext.At ω t ≈ (marked with a vertical line in the large bottom panel of Fig. 8), the classical sur-vival probability attains its asymptotic equilibration value S ∞ P and the Wigner distribution coversthe entire energy shell, as seen in the top panels of Fig. 8 for ω t = 31 . and ω t = 400 (fourthand fifth columns). In contrast, the quantum S P ( t ) continues to raise, leading to a longer equi-libration time. For ω t > , the temporal averaged S P ( t ) coincides again extremely well withthe analytical expression for S ( r ) P ( t ) . During this ramp toward equilibration, the dynamics is con-trolled by the b function, which reflects the correlations between the energy levels in the chaoticregime. For the chaotic Dicke model, these correlations are equivalent to those in the spectrumof GOE random matrices [43, 52, 68, 86], which justifies the use in Eq. (25) of the same form ofthe b function from random matrix theory. At this time scale the dynamics becomes universal.The region of the ramp is referred to in Fig. 8 as “GOE correlations” and is commonly known ascorrelation hole [16–22, 87, 88]. This is a quantum effect associated with the discreteness of thespectrum, which does not appear for the classical S P ( t ) .In contrast to the scarred state III, the quantum survival probability of the state IV reaches itsasymptotic value at a time t r that is longer than the classical relaxation time, and is given by thesame expression obtained with the analytical equation for S ( r ) P ( t ) (see Ref. [22] and Appendix F), t r = πν c (cid:112) δ S P , (26)where δ S P is a small parameter indicating that the values of S P ( t > t r ) are already within thequantum fluctuations around S ∞ P . This time is marked with a vertical line in Fig. 8.The asymptotic value of the quantum survival probability agrees with the saturation point of S ( r ) P ( t ) (see Appendix F), S ∞ P = S ( r ) , ∞ P = (cid:104) r (cid:105)(cid:104) r (cid:105) √ πσν c = 1 √ πσν c , (27)where (cid:104) r (cid:105) and (cid:104) r (cid:105) are the first and second moments of the exponential distribution used to generatethe random numbers in Eq. (13), and the ratio (cid:104) r (cid:105)(cid:104) r (cid:105) = 2 (28)8 FIG. 8. Top panels: Classical evolution of the Wigner distribution projected onto the Q - P plane (top row)and q - p plane (bottom row) for the initial coherent state IV from Fig. 1 (d). The initial distribution (smallblack circle) is inside the available phase space. Each column represents a time, as indicated. (Animationscan be found in the SM [83].) Bottom panel: Quantum survival probability (light gray solid line), itstemporal average (blue solid line), classical survival probability (dark gray solid line), and its temporalaverage (red solid line) for the initial coherent state IV, as well as the analytical expression for randomstates given in Eq. (27) (orange solid line). The temporal averages are computed for temporal windows ofconstant size in the logarithmic scale. The horizontal blue dotted line is the quantum asymptotic value S ∞ P [Eq. (27)] and the horizontal red dashed line corresponds to the classical asymptotic value S ∞ P [Eq. (24)].The vertical line at ω t ≈ marks the correlation hole and beginning of the ramp of S P ( t ) towardssaturation. The relaxation time t r [Eq. (26)], marked with a vertical line, is calculated using δ S P = 0 . , so ω t r = 345 . . is determined by the fluctuations of the energy components of the LDoS of the random state withrespect to the exact Gaussian envelope for state IV. One sees that due to this ratio, the value of S ∞ P is twice as large as the asymptotic value of the classical survival probability, which is givenby Eq. (24). If the fluctuations of the energy components of the LDoS of the random state wereabsent and (cid:104) r (cid:105) / (cid:104) r (cid:105) was one, S ∞ P would coincide with the classical result. The origin of thesefluctuations is rooted in a remaining structure present even in strongly chaotic eigenstates (theso-called nodal structure [27, 29, 89]), which effectively limits the ergodicity of the quantumevolution [90] when compared to the classical dynamics.9We close this discussion with an interesting observation. The minimum value of S ( r ) P ( t ) ,reached right after the Gaussian decay, coincides exactly with S ∞ P . That is, before the ramptowards equilibration, S ( r ) P ( t ) stabilizes at the value associated with the ergodicity of the classicalevolution. Why the random state is able to reach this greatest level of spreading to later contractto the value S ∞ P of maximal quantum ergodicity is an open question to us.
3. Quantum Asymptotic Values and Unstable Periodic Orbits
The purpose of this subsection is to show numerically that most initial coherent states in thechaotic region are indeed marginally affected by unstable periodic orbits, being well described bythe behavior of the survival probability reported in Fig. 8.In the central panel of Fig. 9, we show the ratio R = S ( r ) , ∞ P S ∞ P = 2 S ∞ P S ∞ P (29)for a mesh of initial coherent states distributed over the same Poincar´e surface as in Fig. 1 (c). Inthe equation above, S ( r ) , ∞ P is given by the analytical expression in Eq. (27), which is the asymptoticvalue of an initial coherent state without the influence of unstable periodic orbits. We see that formost of the initial coherent states the ratio R is very close to one (light color). There are few ratiosthat are smaller than one (dark color) and therefore affected by unstable periodic orbits. The ratio R is twice the factor F introduced in [28] to gauge the fraction of phase space explored by thequantum state.Around the central panel of Fig. 9, several plots of the quantum (blue solid line) and classical(red solid line) survival probabilities averaged over temporal windows are shown together with theanalytical expression in Eq. (25) for the random initial state (orange solid line). Panels (a), (b),(d)-(g) make it clear that a ratio R ≈ leads to a generic behavior of the quantum S P ( t ) , welldescribed by S ( r ) P ( t ) . In all these panels, there appears a ramp towards equilibration associatedwith the presence of correlated eigenvalues. These random-like coherent states contrast with thosewhere ratio R < , which are affected by quantum scarring, as in Fig. 9 (c) and (h). These lattercases lead to recurrences in the evolution of the survival probability and a strong influence of theunstable periodic orbits in the quantum dynamics. V. CONCLUSIONS
Our comprehensive study of the dynamics and equilibration process of the Dicke model showsthat the quantum survival probability of coherent states agrees extremely well with its classicallimit obtained with the truncated Wigner approximation up to the point where either the classicalor the quantum survival probability equilibrates. It is at this point, at very long-times, that one canidentify major aspects of the model, such as those itemized below. – Features of the equilibration process that are purely quantum and caused by the dis-creteness of the energy spectrum.
One example is the quantum fluctuations around the asymp-totic value of the quantum survival probability, which are present in both the regular and chaoticregime. The other is the onset of the correlation hole caused by the correlations between theeigenvalues, which is a universal behavior emerging in the chaotic regime for non-scarred initialstates.0 (a) (b) (c)(e)(d) (g) (h)(f ) . . . . . . . . . . R . FIG. 9. Central panel: Map of the ratio R [Eq. (29)] for a mesh of initial coherent states in the chaoticregime. Panels (a)-(h) show the quantum (blue solid line) and classical (red solid line) temporal averagedsurvival probability and the analytical expression (orange solid line) for random states [Eq. (25)]. Thehorizontal dashed line represents the quantum asymptotic value S ∞ P . – Whether the initial state is located in the regular or chaotic regime. The equilibration timeof the classical and quantum dynamics coincides in the regular regime, but differ in the chaoticregion. In the latter case, correlations between the eigenvalues or the proximity to a periodic orbitseparates the curves of the classical and quantum survival probability, one saturating earlier thanthe other. – Whether the initial state in the regular regime is close to the separatrix.
For initialstates located in the border of disconnected classical phase-space regions, the dynamical tunnelingbetween these regions leads to an asymptotic value of the quantum survival probability slightlylower than the classical value, while away from the separatrix the asymptotic values coincide. – The degree of scarring of the initial coherent states in the chaotic regime.
For initial co-herent states affected by quantum scarring, the ratio between the asymptotic value of the quantumsurvival probability and the asymptotic value of the classical survival probability is larger than 2and the equilibration of the quantum dynamics happens before the classical one. The larger theratio is, the larger the degree of scarring. – The onset of maximal quantum ergodicity.
In the chaotic regime, for initial coherentstates unaffected by quantum scarring, the quantum evolution of the survival probability at longtimes coincides with that for random states. This analogy allows us to derive analytical results1for S ∞ P and for the relaxation time. We find that the ratio between the asymptotic values of thequantum and classical survival probabilities equals 2 and classical equilibrium takes place beforethe saturation of the quantum dynamics. The fact that this ratio is not 1 indicates that the quantumevolution is more restricted than the classical one due to remaining structures of the quantumstates.This work makes evident the great advantage of having access to the TWA in the phase space,which can be used to explain various features of the quantum evolution of the survival probability.This includes time intervals where the values of S P ( t ) are close to zero, revivals, and the presenceof periodic orbits. ACKNOWLEDGEMENTS
We acknowledge the support of the Computation Center - ICN, in particular to Enrique Pala-cios, Luciano D´ıaz, and Eduardo Murrieta, and valuable conversations with Jonathan Torresand Jorge Ch´avez-Carlos. We acknowledge financial support from Mexican CONACyT projectCB2015-01/255702, DGAPA- UNAM project IN109417. LFS is supported by the NSF grantNo. DMR-1936006. LFS and JGH thank the hospitality of the Aspen Center for Physics and theSimons Center for Geometry and Physics at Stony Brook University, where some of the researchfor this paper was performed.
Appendix A: Classical Limit of the Dicke Hamiltonian
To obtain the classical Hamiltonian in Eq. (4), we use the Glauber [Eq. (2)] and Bloch [Eq. (3)]coherent states expressed in terms of the general parameters α, z ∈ C as a tensor product of theform, | α, z (cid:105) = | α (cid:105) ⊗ | z (cid:105) = e −| α | / (1 + | z | ) j ∞ (cid:88) n =0 j (cid:88) m = − j (cid:115)(cid:18) jj + m (cid:19) α n z j + m √ n ! | n (cid:105) ⊗ | j, m (cid:105) , (A1)and take the expectation value of the Dicke Hamiltonian ˆ H D [67], H cl = (cid:104) α, z | ˆ H D | α, z (cid:105) = ω | α | − jω − | z | | z | + γ (cid:112) j z + z ∗ | z | ( α + α ∗ ) . (A2)Considering the harmonic oscillator α = (cid:114) j q + ip ) and the Bloch sphere z = (cid:114) j z − j z e − iφ parameters in terms of canonical variables ( q, p ) and ( φ, j z ) , we obtain H cl = j (cid:20) ω q + p ) + ω j z + 2 γ (cid:112) − j z q cos( φ ) (cid:21) . (A3)To finally get Eq. (4), we perform a canonical transformation to the atomic variables ( j z , φ ) → ( Q, P ) , where Q = (cid:112) j z ) cos( φ ) and P = − (cid:112) j z ) sin( φ ) satisfy the Poisson bracket { Q, P } = 1 , and rescale the overall classical Hamiltonian H cl by j .2To obtain the representation of coherent states given by Eq. (2) and Eq. (3), we express thecoherent state parameters α and z in terms of the canonical variables ( q, p, Q, P ) , such that α = (cid:114) j q + ip ) and z = 1 √ − Z ( Q + iP ) with Z = Q + P .With the classical Hamiltonian H cl , a semi-classical approximation to the DoS can be found byintegrating the phase space volume available for a given energy. It reads [68, 91] ν ( (cid:15) ) = 2 jω π (cid:82) y + y − d y cos − (cid:18)(cid:113) y − (cid:15) )¯ γ (1 − y ) (cid:19) , (cid:15) gs ≤ (cid:15) < − (cid:15) + π (cid:82) y + (cid:15) d y cos − (cid:18)(cid:113) y − (cid:15) )¯ γ (1 − y ) (cid:19) , | (cid:15) | ≤ , (cid:15) > , (A4)where y ± = − ¯ γ − (cid:16) ¯ γ − ∓ (cid:112) (cid:15) − (cid:15) ) (cid:17) and ¯ γ = γ/γ c . The ground state energy is given by (cid:15) gs = − for the normal phase, and by (cid:15) gs = (cid:15) = − (cid:0) ¯ γ + ¯ γ − (cid:1) for the superradiant phase. InFig. 10, we compare the DoS obtained numerically with Eq. (A4), showing excellent agreement.The figure also shows the energies selected for our studies throughout this paper. ϵ R ϵ C - - - - ϵ ω ν / j FIG. 10. Density of states evaluated numerically (blue dots) with bin size ∆ (cid:15) = 0 . and analytical expres-sion (A4) (red solid line). Hamiltonian parameters: ω = ω = 1 , ¯ γ = 2 , and j = 100 . Vertical blackdotted and dashed lines indicate respectively the energies (cid:15) R = − . and (cid:15) C = − . selected for our stud-ies. A truncated Hilbert space was employed using the basis of Refs. [92, 93], ensuring 30825 convergedeigenenergies, which range from the ground state energy (cid:15) gs = − . to the truncation at (cid:15) T = 0 . . Appendix B: Wigner Function of Glauber and Bloch Coherent States
As we see from Appendix A, the coherent states for the Dicke model are the product of Glauberand Bloch coherent states [Eq. (A1)]. Thus, the associated Wigner function is the product ofthe Wigner functions associated to the Glauber and Bloch states. For a Glauber state | α (cid:105) with α ( q, p ) = (cid:113) j ( q + ip ) , the Wigner function is given by a normal distribution w q ,p ( q, p ) = jπ e − jd , (B1)3with d = (cid:113) ( q − q ) + ( p − p ) . For a Bloch coherent state, the Wigner function in variables ( θ, φ ) may be written as a sum of Legendre polynomials P k ( x ) [94] w θ ,φ ( θ, φ ) = (2 j )!4 π j (cid:88) k =0 (cid:115) (2 k + 1)(2 j − k )!(2 j + k + 1)! P k (cos Θ) , (B2)where Θ is the angle between ( θ, φ ) and ( θ , φ ) obtained from cos Θ = cos θ cos θ + sin θ sin θ cos( φ − φ ) . (B3)For large j values this is very well approximated by a normal distribution on the Bloch sphere w θ ,φ ( θ, φ ) ≈ jπ e − j Θ . (B4)This approximation is very useful, because sampling from normal distributions is easy and com-putationally cheap. The complete coherent Wigner function is just the product of the Wignerfunctions given by Eq. (B1) and Eq. (B4), w u ( u ) = (cid:18) jπ (cid:19) e − j ( d +Θ ) . (B5) Appendix C: Truncated Wigner Aproximation and Monte Carlo Method
The temporal evolution of the Wigner function is governed by the so-called Moyal equation [95] ∂w∂t ( u , t ) = (cid:8) w ( u , t ) , h cl ( u ) (cid:9) M . (C1)Here, { A, B } M = 2 (cid:126) A sin (cid:20) (cid:126) (cid:16) ←− ∂ q −→ ∂ p − ←− ∂ p −→ ∂ q (cid:17)(cid:21) B represents the Moyal bracket. Taylor expandingthe sine, one may write { A, B } M = { A, B } + O (cid:0) (cid:126) (cid:1) where {· , ·} is the Poisson bracket. So ∂w∂t ( u , t ) = (cid:8) w ( u , t ) , h cl ( u ) (cid:9) + O (cid:0) j − (cid:1) . (C2)If we ignore the j − order terms in this equation, we are left with the classical Liouville equa-tion. This is known as the truncated Wigner approximation (TWA) and yields the correct quantumevolution for small times.Within this approximation, the Wigner function remains constant along classical trajectories inphase space, so the time dependence of w ( u , t ) may be written in terms of the Hamiltonian flow ϕ t : M → M . This function describes the time evolution of an initial condition u ∈ M by u ( t ) = ϕ t ( u ) , (C3) Formally, the Hamiltonian that should be used in the Moyal equation is the Weyl transform of the quantum Hamil-tonian. Our Hamiltonian was not obtained by Weyl transformation but by the calculation of the expectation valuewith coherent states. The difference between both Hamiltonians turns out to be equal to the constant energy ω + ω .Because only the derivatives of the Hamiltonian appear in the equation, this constant number makes no actualdifference. ϕ = Id, ϕ − t = ( ϕ t ) − , and ϕ t + t = ϕ t ◦ ϕ t .Staying constant along classical trajectories means that for any pair of times t and t w ( ϕ t ( u ) , t ) = w ( ϕ t ( u ) , t ) . (C4)In particular, taking t = 0 and performing an adequate change of variables, w ( u , t ) = w ( ϕ − t ( u ) , . (C5)Inserting Eq. (C5) into Eq. (20), we find Eq. (21).Note that Eq. (21) may be written as S P ( t ) = (cid:68) w (cid:0) ϕ − t ( u ) (cid:1)(cid:69) w . (C6)Provided that w is everywhere positive, as it is for coherent states [see Eq. (B5)], this expectedvalue may be efficiently approximated using a Monte Carlo method, S P ( t ) ≈ M M (cid:88) i =1 w ( ϕ − t ( u i )) , (C7)where the points u i ∈ M are randomly sampled from the initial distribution w , and M is asufficiently large, albeit computationally accessible, integer. Appendix D: Intuitive Interpretation of the Classical Survival Probability
An intuitive understanding of the classical survival probability may be obtained by consideringa simple measurable set
S ⊆ M (a sphere, for example). Let S be its indicator function .Consider the normalized distribution ρ S = S /V S in lieu of the Wigner function. Here V S = (cid:82) d u S ( u ) is the volume of S . In this case, to retain normalization, our prefactor becomes /V S instead of (2 π/j ) , and then, from Eq. (21), S P ( t ) = 1 V S (cid:90) M d u S ( u ) S ( ϕ − t ( u )) . (D1)Since S ( ϕ − t ( u )) = ϕ t ( S ) ( u ) and using that the product of the indicator functions of two sets isthe indicator function of their intersection, we get S P ( t ) = 1 V S (cid:90) M d u S∩ ϕ t ( S ) ( u ) = V S∩ ϕ t ( S ) V S . (D2)This result is easy to interpret: S P ( t ) is the percentage of the volume of S that is found backinside S after time t . In other words, it is the probability that a point sampled from S is back in S after time t [25, 96]. For any set
X ⊆ M , its indicator function X : M → { , } equals 1 for points in X and 0 for points outside. Appendix E: Asymptotic Value of the Classical Survival Probability
To obtain an expression for S ∞ P , we assume ergodicity. This is a reasonable assumption forchaotic behaviors. If the Hamiltonian flow ϕ t is ergodic over the energy shells in phase space,then, for any real function f ( u ) of phase space, temporal averages in composition with the floware equal to space averages over the corresponding energy shell, that is, for any fixed point u ∈ M with energy E = H cl ( u ) , (cid:10) f ( ϕ t ( u )) (cid:11) t →∞ = (cid:104) f (cid:105) E = j (2 π ) ν ( E ) (cid:90) M d v δ ( H cl ( v ) − E ) f ( v ) , (E1)where j − (2 π ) ν ( E ) is the volume of the energy shell for E/j in M , and ν ( E ) is given byEq. (A4).It is then straightforward to calculate (cid:10) S P ( t ) (cid:11) t →∞ = (cid:18) πj (cid:19) (cid:90) M d u w ( u ) (cid:10) w ( ϕ − t ( u )) (cid:11) t →∞ = (cid:18) πj (cid:19) (cid:90) M d u w ( u ) (cid:104) w (cid:105) H cl ( u ) = (cid:18) πj (cid:19) (cid:90) ∞ E gs d E (cid:104) w (cid:105) E ν ( E ) . (E2)The last equality is obtained by using (cid:82) ∞ E gs d E δ ( E − H cl ( u )) inside the integral, a change inthe integration order, and a substitution of the value of (cid:104) w (cid:105) E . In the above, E gs represents theground state energy. The classical energy distribution for the state associated to w is ρ cl ( E ) = (cid:90) M d u δ ( E − H cl ( u )) w ( u )=(2 π ) j − ν ( E ) (cid:104) w (cid:105) E , (E3)therefore (cid:10) S P ( t ) (cid:11) t →∞ = (cid:90) ∞ E gs d E ρ cl ( E ) ν ( E ) . (E4)This result is exact, but we may approximate ν ( E ) [see Eq. (A4)] by the value at the center of theclassical energy distribution ( ν ( E ) ≈ ν ( E c ) = ν c ) to get (cid:10) S P ( t ) (cid:11) t →∞ = 1 ν c (cid:90) ∞ E gs d E ρ cl ( E ) . (E5)By further approximating ρ cl with the Gaussian distribution given by Eq. (11), we obtain Eq. (24). Appendix F: Analytical Expression of the Survival Probability for a Random Ensemble with aGaussian Energy Profile
An analytical expression for the survival probability averaged over an ensemble of randominitial states was derived in [22]. The properties of the ensemble are as follows. Its members are6constrained to have a smooth LDoS ρ ( E ) and their energy components are given by | c ( r ) k | = r k ρ ( E k ) A ν ( E k ) , (F1)where ν ( E ) is the DoS, the numbers r k are randomly generated from a given probability distribu-tion P ( r ) , and A = (cid:80) q r q ρ ( E q ) /ν ( E q ) is a normalization constant. This gives, S ( r ) P ( t ) = 1 − S ( r ) , ∞ P η − (cid:34) ηS bc P ( t ) − b (cid:18) Dt π (cid:19)(cid:35) + S ( r ) , ∞ P , (F2)where the initial decay of the survival probability is dictated by S bc P ( t ) = (cid:12)(cid:12)(cid:12)(cid:12)(cid:90) d Eρ ( E ) e − iEt (cid:12)(cid:12)(cid:12)(cid:12) , (F3)the effective dimension of the ensemble η = ν c (cid:82) d E ρ ( E ) , (F4)where ν c = ν ( E c ) is the density of states evaluated in the center of the energy profile, and (cid:104) r n (cid:105) are the n -th moments of the distribution P ( r ) , and the asymptotic value of the survival probabilitycorresponds to S ( r ) , ∞ P = (cid:104) r (cid:105)(cid:104) r (cid:105) η . (F5)The function b is the two-level form factor of the GOE [85] b (¯ t ) =[1 − t + ¯ t ln(2¯ t + 1)]Θ(1 − ¯ t )++ (cid:34) ¯ t ln (cid:18) t + 12¯ t − (cid:19) − (cid:35) Θ(¯ t − , (F6)where Θ is the Heaviside step function. The factor D in the argument of b is the mean levelspacing of the correlated eigenvalues.For the random ensemble that we consider in this paper, ρ ( E ) has a Gaussian energy profile[see Eq. (11)] and the random numbers r k are generated from the exponential distribution P ( r ) = λe − λr with (cid:104) r n (cid:105) = n ! /λ n in Fig. 3, which implies that S bc P ( t ) = e − σ t , (F7) η = 2 √ πσν c , (F8) S ( r ) , ∞ P = 2 η = 1 √ πσν c . (F9)Because the correlations in the spectrum appears only for energy levels in the same paritysector, the mean level spacing of correlated eigenvalues is D = ( D + + D − ) / , where the meanlevel spacing for each parity sector, D ± = 1 /ν ± , is given by the respective density of states. Thesedensities are, in turn, given by ν ± = ν c / , with ν c the density of states of the whole spectrum,yielding D = 2 /ν c . From the previous results for S ( r ) , ∞ P , η , S bc P , and D , we get Eq. (25).7To derive the relaxation time t r of the ensemble-averaged survival probability, we consider theasymptotic form of b , which grows toward saturation following a power-law behavior b (cid:18) tπν c (cid:19) → π ν c t for tπν c (cid:29) . (F10)At this temporal scale, the contribution of the initial decay S bc P is negligible and the asymptoticform of (F2) is given by S ( r ) P ( t ) → − η b (cid:18) tπν c (cid:19) + S ( r ) , ∞ P = S ( r ) , ∞ P (cid:32) − π ν c t (cid:33) , (F11)where in the last step we have used η = 2 /S ( r ) , ∞ P . We define the relaxation time according to S ( r ) P ( t r ) = (1 − δ S P ) S ( r ) , ∞ P , (F12)where δ S P is a small parameter determining the point where S ( r ) P ( t ) is already within the fluctua-tions around the asymptotic value. By substituting Eq. (F11) in Eq. (F12), we obtain Eq. (26). [1] T. Kinoshita, T. Wenger, and D. S. Weiss, “A quantum Newton’s cradle,” Nature , 900 (2006).[2] J. Simon, W. S. Bakr, R. Ma, M. E. Tai, P. M. Preiss, and M. Greiner, “Quantum simulation ofantiferromagnetic spin chains in an optical lattice,” Nature (London) , 307–312 (2011).[3] Michael Schreiber, Sean S. Hodgman, Pranjal Bordia, Henrik P. L¨uschen, Mark H. Fischer, RonenVosk, Ehud Altman, Ulrich Schneider, and Immanuel Bloch, “Observation of many-body localizationof interacting fermions in a quasirandom optical lattice,” Science , 842–845 (2015).[4] Adam M. Kaufman, Alexander Lukin M. Eric Tai, Matthew Rispoli, Robert Schittko, Philipp M.Preiss, and Markus Greiner, “Quantum thermalization through entanglement in an isolated many-body system,” Science , 794 (2016).[5] Juan Maldacena and Douglas Stanford, “Remarks on the Sachdev-Ye-Kitaev model,” Phys. Rev. D , 106002 (2016).[6] Efim B. Rozenbaum, Sriram Ganeshan, and Victor Galitski, “Lyapunov exponent and out-of-time-ordered correlator’s growth rate in a chaotic system,” Phys. Rev. Lett. , 086801 (2017).[7] Efim B. Rozenbaum, Sriram Ganeshan, and Victor Galitski, “Universal level statistics of the out-of-time-ordered operator,” Phys. Rev. B , 035112 (2019).[8] Koji Hashimoto, Keiju Murata, and Ryosuke Yoshii, “Out-of-time-order correlators in quantum me-chanics,” J. High Energy Phys. , 138 (2017).[9] Ignacio Garc´ıa-Mata, Marcos Saraceno, Rodolfo A. Jalabert, Augusto J. Roncaglia, and Diego A.Wisniacki, “Chaos signatures in the short and long time behavior of the out-of-time ordered correla-tor,” Phys. Rev. Lett. , 210601 (2018).[10] Rodolfo A. Jalabert, Ignacio Garc´ıa-Mata, and Diego A. Wisniacki, “Semiclassical theory of out-of-time-order correlators for low-dimensional classically chaotic systems,” Phys. Rev. E , 062218(2018).[11] Jorge Ch´avez-Carlos, B. L´opez-del Carpio, Miguel A. Bastarrachea-Magnani, Pavel Str´ansk´y, SergioLerma-Hern´andez, Lea F. Santos, and Jorge G. Hirsch, “Quantum and classical Lyapunov exponentsin atom-field interaction systems,” Phys. Rev. Lett. , 024101 (2019). [12] Silvia Pappalardi, Angelo Russomanno, Bojan ˇZunkoviˇc, Fernando Iemini, Alessandro Silva, andRosario Fazio, “Scrambling and entanglement spreading in long-range spin chains,” Phys. Rev. B ,134303 (2018).[13] Quirin Hummel, Benjamin Geiger, Juan Diego Urbina, and Klaus Richter, “Reversible quantum in-formation spreading in many-body systems near criticality,” Phys. Rev. Lett. , 160401 (2019).[14] Sa´ul Pilatowsky-Cameo, Jorge Ch´avez-Carlos, Miguel A. Bastarrachea-Magnani, Pavel Str´ansk´y, Ser-gio Lerma-Hern´andez, Lea F. Santos, and Jorge G. Hirsch, “Positive quantum Lyapunov exponentsin experimental systems with a regular classical limit,” Phys. Rev. E , 010202(R) (2020).[15] Emiliano M. Fortes, Ignacio Garc´ıa-Mata, Rodolfo A. Jalabert, and Diego A. Wisniacki, “Gaug-ing classical and quantum integrability through out-of-time-ordered correlators,” Phys. Rev. E ,042201 (2019).[16] Luc Leviandier, Maurice Lombardi, R´emi Jost, and Jean Paul Pique, “Fourier transform: A tool tomeasure statistical level properties in very complex spectra,” Phys. Rev. Lett. , 2449–2452 (1986).[17] Joshua Wilkie and Paul Brumer, “Time-dependent manifestations of quantum chaos,” Phys. Rev. Lett. , 1185–1188 (1991).[18] Y. Alhassid and R. D. Levine, “Spectral autocorrelation function in the statistical theory of energylevels,” Phys. Rev. A , 4650–4653 (1992).[19] E. J. Torres-Herrera and L. F. Santos, “Dynamical manifestations of quantum chaos: Correlation holeand bulge,” Phil. Trans. R. Soc. A , 20160434 (2017).[20] E. J. Torres-Herrera, Antonio M. Garc´ıa-Garc´ıa, and Lea F. Santos, “Generic dynamical featuresof quenched interacting quantum systems: Survival probability, density imbalance, and out-of-time-ordered correlator,” Phys. Rev. B , 060303 (2018).[21] Mauro Schiulaz, E. Jonathan Torres-Herrera, and Lea F. Santos, “Thouless and relaxation time scalesin many-body quantum systems,” Phys. Rev. B , 174313 (2019).[22] S. Lerma-Hern´andez, D. Villase˜nor, M. A. Bastarrachea-Magnani, E. J. Torres-Herrera, L. F. Santos,and J. G. Hirsch, “Dynamical signatures of quantum chaos and relaxation time scales in a spin-bosonsystem,” Phys. Rev. E , 012218 (2019).[23] Jordan S. Cotler, Guy Gur-Ari, Masanori Hanada, Joseph Polchinski, Phil Saad, Stephen H. Shenker,Douglas Stanford, Alexandre Streicher, and Masaki Tezuka, “Black holes and random matrices,” J.High Energy Phys. , 118 (2017).[24] Tokiro Numasawa, “Late time quantum chaos of pure states in random matrices and in the sachdev-ye-kitaev model,” Phys. Rev. D , 126017 (2019).[25] Philip Pechukas, ““quantum chaos” in the irregular spectrum,” Chem. Phys. Lett. , 553–557 (1982).[26] E B Stechel and E J Heller, “Quantum ergodicity and spectral chaos,” Ann. Rev. Phys. Chem. ,563–589 (1984).[27] Eric J. Heller, “Bound-state eigenfunctions of classically chaotic Hamiltonian systems: Scars of peri-odic orbits,” Phys. Rev. Lett. , 1515–1518 (1984).[28] Eric J. Heller, “Quantum localization and the rate of exploration of phase space,” Phys. Rev. A ,1360–1370 (1987).[29] E. J. Heller, “Wavepacket dynamics and quantum chaology,” in Les Houches Summer School 1991 onChaos and Quantum Physics , edited by M.-J. Giannoni, A. Voros, and J. Zinn Justin (Springer, 1991).[30] W. E. Bies, L. Kaplan, and E. J. Heller, “Scarring effects on tunneling in chaotic double-well poten-tials,” Phys. Rev. E , 016204 (2001).[31] M. V. Berry, “Regular and irregular semiclassical wavefunctions,” J. Phys. A , 2083–2091 (1977).[32] Kenneth G. Kay, “Toward a comprehensive semiclassical ergodic theory,” J. Chem. Phys. , 3026–3050 (1983). [33] I C Percival, “Regular and irregular spectra,” J. Phys. B , L229–L232 (1973).[34] C. J. Turner, A. A. Michailidis, D. A. Abanin, M. Serbyn, and Z. Papif´a, “Weak ergodicity breakingfrom quantum many-body scars,” Nat. Phys. , 745–749 (2018).[35] Wen Wei Ho, Soonwon Choi, Hannes Pichler, and Mikhail D. Lukin, “Periodic orbits, entanglement,and quantum many-body scars in constrained models: Matrix product state approach,” Phys. Rev. Lett. , 040603 (2019).[36] Peter Reimann, “Foundation of statistical mechanics under experimentally realistic conditions,” Phys.Rev. Lett. , 190403 (2008).[37] A. J. Short, “Equilibration of quantum systems and subsystems,” New J. Phys. , 053009 (2011).[38] A. J. Short and T. C. Farrelly, “Quantum equilibration in finite time,” New J. Phys. , 013063 (2012).[39] Pablo R. Zangara, Axel D. Dente, E. J. Torres-Herrera, Horacio M. Pastawski, An´ıbal Iucci, andLea F. Santos, “Time fluctuations in isolated quantum systems of interacting particles,” Phys. Rev. E , 032913 (2013).[40] Kai He, Lea F. Santos, Tod M. Wright, and Marcos Rigol, “Single-particle and many-body analysesof a quasiperiodic integrable system after a quench,” Phys. Rev. A , 063637 (2013).[41] R. H. Dicke, “Coherence in spontaneous radiation processes,” Phys. Rev. , 99 (1954).[42] Klaus Hepp and Elliott H Lieb, “On the superradiant phase transition for molecules in a quantizedradiation field: the Dicke maser model,” Ann. Phys. (N.Y.) , 360 – 404 (1973); Klaus Hepp andElliott H. Lieb, “Equilibrium statistical mechanics of matter interacting with the quantized radiationfield,” Phys. Rev. A , 2517–2525 (1973); Y. K. Wang and F. T. Hioe, “Phase transition in the Dickemodel of superradiance,” Phys. Rev. A , 831–836 (1973).[43] Clive Emary and Tobias Brandes, “Chaos and the quantum phase transition in the Dicke model,” Phys.Rev. E , 066203 (2003).[44] Barry M. Garraway, “The Dicke model in quantum optics: Dicke model revisited,” Philos. Trans.Royal Soc. A , 1137 (2011).[45] Kristian Baumann, Christine Guerlin, Ferdinand Brennecke, and Tilman Esslinger, “Dicke quantumphase transition with a superfluid gas in an optical cavity,” Nature (London) , 1301 (2010).[46] K. Baumann, R. Mottl, F. Brennecke, and T. Esslinger, “Exploring symmetry breaking at the Dickequantum phase transition,” Phys. Rev. Lett. , 140402 (2011).[47] Helmut Ritsch, Peter Domokos, Ferdinand Brennecke, and Tilman Esslinger, “Cold atoms in cavity-generated dynamical optical potentials,” Rev. Mod. Phys. , 553–601 (2013).[48] Markus P. Baden, Kyle J. Arnold, Arne L. Grimsmo, Scott Parkins, and Murray D. Barrett, “Re-alization of the Dicke model using cavity-assisted Raman transitions,” Phys. Rev. Lett. , 020408(2014).[49] J. Klinder, H. Keßler, M. Reza Bakhtiari, M. Thorwart, and A. Hemmerich, “Observation of a super-radiant mott insulator in the Dicke-Hubbard model,” Phys. Rev. Lett. , 230403 (2015).[50] Alicia J. Koll´ar, Alexander T. Papageorge, Varun D. Vaidya, Yudan Guo, Jonathan Keeling, andBenjamin L. Lev, “Supermode-density-wave-polariton condensation with a Bose-Einstein condensatein a multimode cavity,” Nat. Comm. , 14386 (2017).[51] C.H Lewenkopf, M.C Nemes, V Marvulle, M.P Pato, and W.F Wreszinski, “Level statistics transitionsin the spin-boson model,” Phys. Lett. A , 113 – 116 (1991).[52] Clive Emary and Tobias Brandes, “Quantum chaos triggered by precursors of a quantum phase transi-tion: The Dicke model,” Phys. Rev. Lett. , 044101 (2003); M. A. Bastarrachea-Magnani, S. Lerma-Hern´andez, and J. G. Hirsch, “Comparative quantum and semiclassical analysis of atom-field systems.ii. Chaos and regularity,” Phys. Rev. A , 032102 (2014); Miguel Angel Bastarrachea-Magnani,Baldemar L´opez del Carpio, Sergio Lerma-Hern´andez, and Jorge G Hirsch, “Chaos in the Dicke model: quantum and semiclassical analysis,” Phys. Scr. , 068015 (2015).[53] M. A. Bastarrachea-Magnani, B. L´opez-del-Carpio, J. Ch´avez-Carlos, S. Lerma-Hern´andez, and J. G.Hirsch, “Delocalization and quantum chaos in atom-field systems,” Phys. Rev. E , 022215 (2016);J. Ch´avez-Carlos, M. A. Bastarrachea-Magnani, S. Lerma-Hern´andez, and J. G. Hirsch, “Classicalchaos in atom-field systems,” Phys. Rev. E , 022209 (2016).[54] P. P´erez-Fern´andez, P. Cejnar, J. M. Arias, J. Dukelsky, J. E. Garc´ıa-Ramos, and A. Rela˜no, “Quantumquench influenced by an excited-state phase transition,” Phys. Rev. A , 033802 (2011).[55] Alexander Altland and Fritz Haake, “Quantum chaos and effective thermalization,” Phys. Rev. Lett. , 073601 (2012).[56] Michal Kloc, Pavel Str´ansk´y, and Pavel Cejnar, “Quantum quench dynamics in Dicke superradiancemodels,” Phys. Rev. A , 013836 (2018).[57] Peter Kirton, Mor M. Roses, Jonathan Keeling, and Emanuele G. Dalla Torre, “Introduction to theDicke model: From equilibrium to nonequilibrium, and vice versa,” Advanced Quantum Technologies , 1800043 (2019).[58] Daniele De Bernardis, Tuomas Jaako, and Peter Rabl, “Cavity quantum electrodynamics in the non-perturbative regime,” Phys. Rev. A , 043820 (2018).[59] Anton Frisk Kockum, Adam Miranowicz, Simone De Liberato, Salvatore Savasta, and Franco Nori,“Ultrastrong coupling between light and matter,” Nature Reviews Physics , 19–40 (2019).[60] P. Forn-D´ıaz, L. Lamata, E. Rico, J. Kono, and E. Solano, “Ultrastrong coupling regimes of light-matter interaction,” Rev. Mod. Phys. , 025005 (2019).[61] Zhiqiang Zhang, Chern Hui Lee, Ravi Kumar, K. J. Arnold, Stuart J. Masson, A. L. Grimsmo, A. S.Parkins, and M. D. Barrett, “Dicke-model simulation via cavity-assisted Raman transitions,” Phys.Rev. A , 043858 (2018).[62] J Cohn, A Safavi-Naini, R J Lewis-Swan, J G Bohnet, M G¨arttner, K A Gilmore, J E Jordan, A M Rey,J J Bollinger, and J K Freericks, “Bang-bang shortcut to adiabaticity in the Dicke model as realizedin a penning trap experiment,” New J. Phys. , 055013 (2018).[63] A. Safavi-Naini, R. J. Lewis-Swan, J. G. Bohnet, M. G¨arttner, K. A. Gilmore, J. E. Jordan, J. Cohn,J. K. Freericks, A. M. Rey, and J. J. Bollinger, “Verification of a many-ion simulator of the Dickemodel through slow quenches across a phase transition,” Phys. Rev. Lett. , 040503 (2018).[64] M. J. Steel and M. J. Collett, “Quantum state of two trapped Bose-Einstein condensates with a Joseph-son coupling,” Phys. Rev. A , 2920–2930 (1998).[65] Anatoli Polkovnikov, “Phase space representation of quantum dynamics,” Ann. Phys. , 1790 –1852 (2010).[66] J. Schachenmayer, A. Pikovski, and A. M. Rey, “Many-body quantum spin dynamics with MonteCarlo trajectories on a discrete phase space,” Phys. Rev. X , 011022 (2015).[67] M. A. M. de Aguiar, K. Furuya, C. H. Lewenkopf, and M. C. Nemes, “Particle-spin coupling in achaotic system: Localization-delocalization in the husimi distributions,” EPL (Europhys. Lett.) ,125 (1991); M.A.M de Aguiar, K Furuya, C.H Lewenkopf, and M.C Nemes, “Chaos in a spin-bosonsystem: Classical analysis,” Annals of Physics , 291 – 312 (1992).[68] M. A. Bastarrachea-Magnani, S. Lerma-Hern´andez, and J. G. Hirsch, “Comparative quantum andsemiclassical analysis of atom-field systems. I. Density of states and excited-state quantum phasetransitions,” Phys. Rev. A , 032101 (2014).[69] A. D. Ribeiro, M. A. M. de Aguiar, and A. F. R. de Toledo Piza, “The semiclassical coherent statepropagator for systems with spin,” J. Phys. A , 3085 (2006).[70] Armando Rela˜no, Carlos Esebbag, and Jorge Dukelsky, “Excited-state quantum phase transitions inthe two-spin elliptic gaudin model,” Phys. Rev. E , 052110 (2016). [71] M A Bastarrachea-Magnani, A Rela˜no, S Lerma-Hern´andez, B L´opez del Carpio, J Ch´avez-Carlos,and J G Hirsch, “Adiabatic invariants for the regular region of the Dicke model,” J. Phys. A , 144002(2017).[72] L. E. Reichl and W. M. Zheng, “Nonlinear resonance and chaos in conservative systems,” in Directionsin Chaos Volume 1 (1987) pp. 17–90.[73] D. A. Usikov G. M. Zaslavsky, R. Z. Sagdeev and A. A. Chernikov,
Weak Chaos and Quasi-RegularPatterns (Cambridge University Press, Bristol, England, 1991).[74] M. A. Bastarrachea-Magnani, B. L´opez-del-Carpio, J. Ch´avez-Carlos, S. Lerma-Hern´andez, and J. G.Hirsch, “Regularity and chaos in cavity qed,” Phys. Scr. , 054003 (2017).[75] Marco T´avora, E. J. Torres-Herrera, and Lea F. Santos, “Inevitable power-law behavior of isolatedmany-body quantum systems and how it anticipates thermalization,” Phys. Rev. A , 041603 (2016).[76] Marco T´avora, E. J. Torres-Herrera, and Lea F. Santos, “Power-law decay exponents: A dynamicalcriterion for predicting thermalization,” Phys. Rev. A , 013604 (2017).[77] E. J. Torres-Herrera, M Vyas, and Lea F. Santos, “General features of the relaxation dynamics ofinteracting quantum systems,” New J. Phys. , 063010 (2014).[78] John Schliemann, “Coherent quantum dynamics: What fluctuations can tell,” Phys. Rev. A , 022108(2015).[79] Sergio Lerma-Hern´andez, Jorge Ch´avez-Carlos, Miguel A. Bastarrachea-Magnani, Lea F. Santos,and Jorge G. Hirsch, “Analytical description of the survival probability of coherent states in regularregimes,” J. Phys. A , 475302 (2018).[80] Ricardo Puebla, Armando Rela˜no, and Joaqu´ın Retamosa, “Excited-state phase transition leading tosymmetry-breaking steady states in the Dicke model,” Phys. Rev. A , 023819 (2013).[81] E. Wigner, “On the quantum correction for thermodynamic equilibrium,” Phys. Rev. , 749–759(1932).[82] William Case, “Wigner functions and Weyl transforms for pedestrians,” Am. J. Phys. (2008),10.1119/1.2957889.[83] See Supplemental Material..[84] S. Keshavamurthy and P. Schlagheck, Dynamical Tunneling. Theory and Expriment (CRC Press, BocaRaton, 2011).[85] M. L. Mehta,
Random Matrices (Academic Press, Boston, 1991).[86] Wouter Buijsman, Vladimir Gritsev, and Rudolf Sprik, “Nonergodicity in the anisotropic Dickemodel,” Phys. Rev. Lett. , 080601 (2017).[87] H. Alt, H.-D. Gr¨af, T. Guhr, H. L. Harney, R. Hofferbert, H. Rehfeld, A. Richter, and P. Schardt,“Correlation-hole method for the spectra of superconducting microwave billiards,” Phys. Rev. E ,6674–6683 (1997).[88] T. Gorin and T. H. Seligman, “Signatures of the correlation hole in total and partial cross sections,”Phys. Rev. E , 026214 (2002).[89] S. Nonnenmacher, “Anatomy of quantum chaotic eigenstates,” S´eminaire Poincar´e XIV , 177–220(2010).[90] E. B. Stechel, “Quantum ergodicity and a quantum measure algebra,” J. Chem. Phys. , 364–371(1985).[91] Tobias Brandes, “Excited-state quantum phase transitions in Dicke superradiance models,” Phys. Rev.E , 032133 (2013).[92] Miguel A Bastarrachea-Magnani and Jorge G Hirsch, “Efficient basis for the Dicke model: I. Theoryand convergence in energy,” Phys. Scripta , 014005 (2014). [93] Jorge G Hirsch and Miguel A Bastarrachea-Magnani, “Efficient basis for the Dicke model: II. Wavefunction convergence and excited states,” Phys. Scripta , 014018 (2014).[94] Popo Yang, Iv´an F Valtierra, Andrei B Klimov, Shin-Tza Wu, Ray-Kuang Lee, Luis L S´anchez-Soto,and Gerd Leuchs, “The wigner flow on the sphere,” Phys. Scr. , 044001 (2019).[95] Andrei B Klimov, Jos´e Luis Romero, and Hubert de Guise, “Generalized SU(2) covariant wignerfunctions and some of their applications,” J. Phys. A , 323001 (2017).[96] T. Gorin, Tomaz Prosen, Thomas H. Seligman, and Marko ˇZnidariˇc, “Dynamics of Loschmidt echoesand fidelity decay,” Phys. Rep.435