Characterizing the dynamical phase diagram of the Dicke model via classical and quantum probes
R. J. Lewis-Swan, S. R. Muleady, D. Barberena, J. J. Bollinger, A. M. Rey
CCharacterizing the dynamical phase diagram of the Dicke model via classical andquantum probes
R. J. Lewis-Swan,
1, 2
S. R. Muleady,
3, 4, ∗ D. Barberena,
3, 4, ∗ J. J. Bollinger, and A. M. Rey
3, 4 Homer L. Dodge Department of Physics and Astronomy,The University of Oklahoma, Norman, Oklahoma 73019, USA Center for Quantum Research and Technology, The University of Oklahoma, Norman, Oklahoma 73019, USA JILA, NIST, Department of Physics, University of Colorado, Boulder, CO 80309, USA Center for Theory of Quantum Matter, University of Colorado, Boulder, CO 80309, USA National Institute of Standards and Technology, Boulder, Colorado 80305, USA (Dated: February 5, 2021)We theoretically study the dynamical phase diagram of the Dicke model in both classical andquantum limits using large, experimentally relevant system sizes. Our analysis elucidates that themodel features dynamical critical points that are distinct from previously investigated excited-stateequilibrium transitions. Moreover, our numerical calculations demonstrate that mean-field featuresof the dynamics remain valid in the exact quantum dynamics, but we also find that in regimeswhere quantum effects dominate signatures of the dynamical phases and chaos can persist in purelyquantum metrics such as entanglement and correlations. Our predictions can be verified in currentquantum simulators of the Dicke model including arrays of trapped ions.
Introduction:
Advances in atomic, molecular and optical(AMO) quantum simulators are driving a surge in theinvestigation of dynamical phase transitions (DPTs) andassociated non-equilibrium phases of matter [1–4]. In aclosed system, a DPT is a critical point separating dis-tinct dynamical behaviours that emerge after a suddenquench of a control parameter [5–11], and can be definedvia non-analytic behaviour in a time-averaged order pa-rameter [12–14].To date, the majority of experimental investigations inthis direction have been tailored towards integrable mod-els, featuring effective infinite-range interactions betweenspins, which admit analytical treatments [15–18]. Richernon-integrable models have been pursued in trappedion systems [19], but the associated complexity of thequantum dynamics limited the theoretical analysis ofthe DPT to small system sizes and prevented a rigor-ous scaling analysis. Hence, it is highly desirable tofind and study DPTs in non-integrable models featuringnovel non-equilibrium phenomena that are both imple-mentable in tunable quantum simulators and theoreti-cally tractable under controllable approximations.We advance this direction by studying a DPT in theiconic Dicke model [20–22], which describes the collectivecoupling of many spins to a single harmonic oscillator.The model is attractive as it features an array of phenom-ena, such as non-integrability [23–25], chaos [26–28] andequilibrium quantum phase transitions (QPTs) in bothground and excited states [29–34], but involves only apair of disparate degrees of freedom such that the modelremains amenable to analytic and numerical treatments.Moreover, the model is already studied in state-of-the-artAMO quantum simulators, including trapped-ion arrays[35, 36] and cavity-QED [37–40]. Here, we investigate theDPT in analytically tractable spin- and boson-dominatedlimits, as well as a non-integrable regime where near- resonant coupling of spin and bosons leads to chaoticdynamical phases also seen in other non-integrable sys-tems [5, 6]. While recent works have linked DPTs to co-existing excited-state QPTs [17, 33, 41], we determine theDPT is a clearly distinct non-equilibrium phenomenonthat explicitly manifests at different parameter values.By studying large, experimentally relevant, system sizesusing efficient numerical methods we are able to demon-strate that mean-field features of the dynamics remainvalid in the exact quantum dynamics. Conversely, inregimes where quantum effects dominate we find that sig-natures of the DPT and chaos persist in purely quantummetrics such as entanglement and correlations.
Model:
The Dicke Hamiltonian for N spin-1 / H D = 2 g √ N (cid:0) ˆ a + ˆ a † (cid:1) ˆ S z + δ ˆ a † ˆ a + Ω ˆ S x . (1)Here, ˆ a (ˆ a † ) is the bosonic annihilation (creation) opera-tor of an oscillator with frequency δ , ˆ S α = 1 / (cid:80) Nj =1 ˆ σ αj are collective spin operators for α = x, y, z and ˆ σ αj Pauli matrices of the j th spin. The spins are subjectto a transverse field of strength Ω and the spin-bosoncoupling is characterized by g . The Hamiltonian ex-hibits a Z parity symmetry associated with the operatorˆΠ = e iπ ( ˆ S x +ˆ a † ˆ a + N/ such that [ ˆ H D , ˆΠ] = 0.The equilibrium phase diagram of the Dicke model fea-tures a ground-state QPT at a critical coupling g QPT = √ δ Ω / g (cid:29) g QPT ) phase, with degenerate ground-states | ψ s (cid:105) = √ [ | − ( N/ z (cid:105) ⊗| α s (cid:105) ±| ( N/ z (cid:105) ⊗| − α s (cid:105) ] associated withdifferent parity sectors, and a normal ( g (cid:28) g QPT ) phasewith | ψ n (cid:105) ≈ | ( N/ x (cid:105) ⊗ | (cid:105) [35]. Here, we have definedcollective spin states via ˆ S x,y,z | m x,y,z (cid:105) = m x,y,z | m x,y,z (cid:105) and |± α s (cid:105) is the bosonic coherent state for α s = g √ N /δ . a r X i v : . [ qu a n t - ph ] F e b Additionally, an excited-state QPT (EQPT) exists in thespectrum of the superradiant phase, g > g
QPT , definedby a critical energy E c = − Ω N/ E < E c . Dynamical phase diagram:
We study the dynamicalphase diagram that arises after a quench of the trans-verse field. Concretely, the system is initialized in theground-state of ˆ H D at fixed g and δ with Ω = 0, suchthat | ψ (0) (cid:105) = | ( − N/ z (cid:105) ⊗ | α s (cid:105) , and the transverse fieldis then quenched to a final value Ω (cid:54) = 0. To garnerinsight into the dynamics we first study the classicalmodel (Figs. 1-3) before probing the role of quantumfluctuations. The classical limit of the Dicke model isequivalent to solving the Heisenberg equations of mo-tion for operators under a mean-field approximation,wherein expectation values are factorized according to (cid:104) ˆ O ( t ) ˆ O ( t ) (cid:105) = (cid:104) ˆ O ( t ) (cid:105)(cid:104) ˆ O ( t ) (cid:105) [42]. For brevity we adoptthe notation O ≡ (cid:104) ˆ O ( t ) (cid:105) herein.It is convenient to study dynamics in terms of two vari-ables: ˜ g = 2 g/ √ δ Ω ≡ g/g QPT and δ/ Ω. The formercharacterizes the effective strength of the spin-boson in-teraction relative to the single-particle terms, while thelatter describes the relative energy scales of spin andbosonic excitations and loosely expresses the relative im-portance of each degree of freedom to the dynamics. Inthe limit δ/ Ω (cid:29) g ∼ H eff = ( χ/N ) ˆ S z + Ω ˆ S x with theboson-mediated interaction characterized by χ ≡ g /δ [43]. For δ/ Ω (cid:28) g ∼ δ/ Ω ∼ δ/ Ω (cid:29) δ/ Ω (cid:28) δ/ Ω ∼ g is shownin the time-traces of Fig. 1 for SDR and BDR. In bothlimits we identify that the dynamics can be character-ized as either trapped or untrapped . The former occurswhen the spin-boson interaction dominates the Hamil-tonian and leads to a locking of the spins and bosonsclose to their initial configuration. In the SDR this hasbeen interpreted as a self-generated detuning ∝ (cid:104) ˆ S z (cid:105) ˆ S z that locks out rotations due to the transverse field [18].Conversely, the untrapped dynamics are characterized bylarge coherent oscillations in S z and X = (cid:104) ˆ a + ˆ a † (cid:105) dom-inated by the single-particle terms of the Hamiltonian.We classify the dynamical phases and identify a DPTusing a pair of interchangeable time-averaged orderparameters S z = lim T →∞ (1 /T ) (cid:82) T S z ( t ) dt and X = -1-0.50 FIG. 1. Dynamical phase diagram for time-averaged orderparameter S z (center) and typical time traces (surrounding)for the Dicke model. Parameters of time-traces are (clock-wise from top left): ( δ/ Ω , ˜ g ) = (0 . , δ/ Ω , ˜ g ) = (4 , δ/ Ω , ˜ g ) = (4 ,
1) and ( δ/ Ω , ˜ g ) = (0 . , S z isobtained up to a maximum time gT = 10 . lim T →∞ (1 /T ) (cid:82) T X ( t ) dt . The trapped phase is definedby non-zero S z (cid:54) = 0 and X (cid:54) = 0, while in the untrappedphase S z = X = 0. The phases are separated by a crit-ical coupling [main panel of Fig. 1]: i) ˜ g s DPT ≈ √ g b DPT ≈ / in theboson-dominated regime. The dynamical phases in theLMG model have recently been observed in a cavity-QEDquantum simulator [18].In the non-integrable RR, δ/ Ω ∼
1, the dynamicalphase diagram is more complex. The main panel of Fig. 1illustrates that the time-averaged order parameter be-comes noisy in the untrapped phase. In this regime typ-ical time-traces [Fig. 2(a)] feature erratic oscillations inboth spin and boson observables and there exist shortperiods where the system abruptly becomes re-trapped.This behaviour signals a chaotic dynamical phase [5, 6]that arises due to the known chaos of the Dicke modelfor ˜ g > δ ∼ Ω [24, 26, 27, 29].We confirm this by computing the Lyapunov exponent λ L [27] in Fig. 2(b), as a function of ˜ g and δ/ Ω for the ini-tial condition corresponding to | ψ (0) (cid:105) . We find regions ofchaos in parameter space, signalled by λ L > S z (Fig. 1).The dynamical phase diagram of SDR and BDR is cap-tured by an effective model involving only the dominant FIG. 2. (a) Typical time-trace near the transition in the res-onant region, ˜ g ≈ .
299 and δ/ Ω = 1. (b) Characterizationof phase-space via Lyapunov exponent λ L . Star indicates pa-rameter regime of time-trace. degree of freedom. Specifically, the mean-field dynamicsare reduced to an equivalent picture of a classical parti-cle with co-ordinate ξ = S z (spin-dominated) or ξ = X (boson-dominated) confined within a 1D potential V ( ξ )that depends only on ˜ g , δ/ Ω and the initial state. Nearthe DPT V ( ξ ) is well approximated by a double well po-tential [Fig. 3(a)].The initial position and velocity of the particle, ξ (0)and ˙ ξ (0), and the relative height of the central maximum, V (0), characterizes the dynamics of the system [6]. For˜ g > ˜ g DPT and the particle initially located in one well, ξ (0) (cid:54) = 0, the particle has insufficient mechanical energyto overcome the barrier and remains confined. On theother hand, for ˜ g < ˜ g DPT the particle has sufficient en-ergy to pass over the barrier and traverses freely betweenboth wells. The former scenario describes trapped dy-namics, ξ (cid:54) = 0, whilst the latter describes the untrappedphase, ξ = 0. The critical point is defined as the con-dition for which the particle first surpasses the centralbarrier and we find ˜ g s DPT = √ g b DPT = 3 / , respec-tively, in agreement with Fig. 1.Notably, the location of the DPT in the SDR coincideswith the EQPT of the Dicke Hamiltonian [31, 32, 34],e.g., the energy of the initial state matches the EQPTcritical energy, E ≡ (cid:104) ψ | ˆ H | ψ (cid:105)| ˜ g s DPT = E c . In fact,prior work has suggested that the DPT could be inter-preted through the lens of the EQPT [17, 33, 41] andthe critical point is determined via this energy relation.However, we explicitly find this is not valid in the BDR,where (cid:104) ψ | ˆ H | ψ (cid:105)| ˜ g b DPT > E c [42], explicitly demonstrat-ing that the DPT is generally a distinct non-equilibriumphenomenon.The potential model is further evidenced by inspect-ing typical trajectories of trapped and untrapped casesin the classical phase-space, as shown in Fig. 3(b) for thethree regimes. In BDR and SDR we observe well definedorbits in the dominant degree of freedom that are cen-tered around fixed point(s) in phase-space correspondingto the minima of the potential well(s) [42]. The tra-jectories of the complementary slaved observables showsimilar behaviour, but tend to densely fill out the avail-able phase-space due to fast micromotion on top of theslower orbits arising from the enslavement to the domi- -1 0 1123 -1 0 1123 -1 0 FIG. 3. (a) In the integrable limits the DPT can be de-scribed as a particle (green markers) in a 1D potential withco-ordinate ξ = S z ( δ/ Ω (cid:29)
1) or ξ = X ( δ/ Ω (cid:28) g (cid:28) ˜ g DPT (untrapped phase, left) the particle has suf-ficient energy to traverse both wells of the potential, whereasfor large ˜ g (cid:29) ˜ g DPT (trapped phase, right) the particle re-mains energetically confined to a single well. (b) Typicaltrajectories in phase-space for trapped (red) and untrapped(blue) phases at δ/ Ω = (0 . , ,
10) (top to bottom). For thespin phase-space we use co-ordinates r = 1 + 2 S x /N and ϕ = arctan( S y /S z ) (c) Dynamical phase diagram as a func-tion of initial state. At small detuning δ/ Ω = 0 . α (top panel), and at large detuning δ/ Ω = 4we vary the tipping angle θ relative to the south pole of theBloch sphere ( S z = − N/
2) (bottom panel). nant degree of freedom (see also Fig. 1). In the resonantregime the potential description breaks down as no singledegree of freedom dominates. Consequently, we observea lack of clear orbits in phase-space, although the erratictrajectories do fill out distinct volumes of phase-space ineach dynamical phase.DPTs in isolated systems are intrinsically dependenton initial conditions [18]. We demonstrate this by prob-ing initial states of the form | ψ (0) (cid:105) = | ( − N/ n (cid:105) ⊗ | α (cid:105) where α ∈ R and we define ˆ S n | m n (cid:105) = m n | m n (cid:105) forˆ S n = ˆS · n with n = (0 , sin( θ ) , cos( θ )) defined by thetipping angle θ from the south pole of the collectiveBloch sphere. In Fig. 3(c) we probe the dependence of˜ g DPT on the initial amplitude α in the BDR, while inthe SDR we vary the initial tipping angle θ . We ob-serve good agreement with analytic predictions for ˜ g DPT based on the potential model [42]. Small deviations areobserved for θ ≈ ± π/ α/α s →
0, which correspondto initializing a particle near the central maximum of thepotential such that even a small off-resonant exchangeof energy with the complementary degree of freedom isenough for the particle to overcome the barrier.
FIG. 4. Signatures of DPT in quantum dynamics (seeRef. [42] for details). (a) Time traces of quantum dynam-ics in the resonant regime, ˜ g = 1 . δ/ Ω = 1, for N = 40(light gray), 200 (gray), and 600 (black). Mean-field solu-tion (green) is shown for reference. (b) Time-averaged en-tanglement entropy S vN = (1 /T ) (cid:82) T + t t dtS vN ( t ) for S vN = − Tr[ˆ ρ s log(ˆ ρ s )] where ˆ ρ s is the reduced density matrix of thespins. (c) Temporal fluctuations of entanglement entropy,(∆ S vN ) = (1 /T ) (cid:82) T + t t dt [ S vN ( t ) − S vN ] , normalized by S vN . (d) Growth rate γ of quantum fluctuations, obtainedby an empirical fit to (cid:104) (∆ ˆ S z ) (cid:105) ∝ (1 − e − γt ). Data in (b)-(d)is obtained for N = 600. Time averages in (b) and (c) arecomputed over a window t ∈ [ t , T + t ] for gt , gT = 150, soas to reduce transient effects in steady state estimates. Grayregions correspond to excluded data. Quantum dynamics:
Using an efficient exact diagonal-ization method [42] we are able to simulate the quantumdynamics and observe mean-field signatures of the DPTthat survive quantum noise. For example, in Fig. 4(a)we probe a typical time-trace of (cid:104) ˆ S z (cid:105) in the resonantregime [same as Fig. 2(a)] for experimentally realistic N = 40 , ,
600 and find that it is possible to graduallyobserve the signature re-trapping dynamics with increas-ing N , before quantum fluctuations dephase (cid:104) ˆ S z (cid:105) → S vN = − Tr[ˆ ρ s log(ˆ ρ s )] where ˆ ρ s is the reduceddensity matrix of the spins. For δ/ Ω (cid:46) S vN = (1 /T ) (cid:82) T + t t dtS vN ( t ) demar-cates the trapped and untrapped phases consistent withthe order parameter S z , with S vNtrapped (cid:28) S vNuntrapped .The entanglement decreases with increasing δ/ Ω as thebosons become eliminated, although we expect that, e.g.,bipartite entanglement between the spins will retain in-dications of the DPT. Panel (c) also demonstrates thatthe temporal fluctuations of the entanglement vary withthe underlying integrability of the Hamiltonian [46, 47].Comparing to the Lyapunov exponent in Fig. 2, we see a correlation between regions of chaotic (non-integrable)dynamics and suppressed temporal fluctuations of S vN ( t )(relative to S vN ).Panel (d) evidences that the critical region of the DPTcan be diagnosed by the rapid growth of quantum fluctu-ations. By empirically fitting (cid:104) (∆ ˆ S z ) (cid:105) = (cid:104) ˆ S z (cid:105) − (cid:104) ˆ S z (cid:105) ∝ (1 − e − γt ) we find that the rate at which quantum fluctua-tions buildup, characterized by γ , is largest for ˜ g ∼ ˜ g DPT .This is consistent with the effective potential description,as near the critical coupling the mean-field trajectoryspends an increasing amount of time probing the cen-tral barrier (e.g., unstable fixed point) [6, 42], motivat-ing the expectation that quantum fluctuations will growexponentially. Nevertheless, we also note that a similarlyrapid growth of fluctuations is observed away from theDPT for 1 (cid:46) δ/ Ω (cid:46) Experimental realization:
The full dynamical phase di-agram of the Dicke model can be studied in a range ofcurrent state-of-the-art AMO quantum simulators, butmost readily in arrays of trapped ions. In particular, theDicke model was recently realized in a 2D Penning trapconfiguration [35, 36] where a spin-1 / g can be in principle controlled viathe applied laser power. Moreover, δ can be varied bycontrolling the detuning of the lasers relative to the fre-quency of the targeted center-of-mass mode [48]. Consid-ering the simulator reported in Refs. [35, 36], we predictit should be possible to probe regimes 0 . (cid:46) δ/ Ω (cid:46) g ∼ N ≈
200 [49] are possible in the 2D geometry,which is sufficient to access the signatures of the DPT asin Figs. 1 and 4 [42]. Looking ahead, implementing theDicke model in 3D ion crystals [50] would open a path toeven larger system sizes that are beyond the capabilityof numerical methods.
Conclusion:
We have studied a DPT in the Dicke modeland unique features arising from non-integrable andchaotic regimes of parameter space. Our numerical studyindicates signatures of the DPT survive in the quantumdynamics and are accessible in current state-of-the-artAMO quantum simulators based on, e.g., trapped-ion ar-rays. This can motivate future investigations connectingto quantum phenomena such as, e.g., information scram-bling [27, 51]. Additionally, studying entanglement dy-namics associated with DPTs in collective systems willopen new directions for the generation of metrologicallyuseful states for quantum-enhanced sensors [52] or fre-quency and time standards [53].
Acknowledgements:
We acknowledge helpful discussionswith Jamir Marino and Matthew Affolter. This work issupported by the AFOSR grant FA9550-18-1-0319, bythe DARPA and ARO grant W911NF-16-1-0576, theARO single investigator award W911NF-19-1-0210, theNSF PHY1820885, NSF JILA-PFC PHY-1734006 andNSF QLCI-2016244 grants, and by NIST. ∗ These two authors contributed equally[1] M. Heyl, A. Polkovnikov, and S. Kehrein, Phys. Rev.Lett. , 135704 (2013).[2] M. Schir´o and M. Fabrizio, Phys. Rev. Lett. , 076401(2010).[3] B. Sciolla and G. Biroli, Journal of Statistical Mechanics:Theory and Experiment , P11003 (2011).[4] B. ˇZunkoviˇc, M. Heyl, M. Knap, and A. Silva, Phys.Rev. Lett. , 130601 (2018).[5] A. Lerose, J. Marino, B. ˇZunkoviˇc, A. Gambassi, andA. Silva, Phys. Rev. Lett. , 130603 (2018).[6] A. Lerose, B. ˇZunkoviˇc, J. Marino, A. Gambassi, andA. Silva, Phys. Rev. B , 045128 (2019).[7] M. Schir´o and M. Fabrizio, Phys. Rev. B , 165105(2011).[8] F. Peronaci, M. Schir´o, and M. Capone, Phys. Rev. Lett. , 257001 (2015).[9] B. Sciolla and G. Biroli, Phys. Rev. B , 201110 (2013).[10] A. Chiocchetta, M. Tavora, A. Gambassi, and A. Mitra,Phys. Rev. B , 220302 (2015).[11] A. Chiocchetta, A. Gambassi, S. Diehl, and J. Marino,Phys. Rev. Lett. , 135701 (2017).[12] M. Eckstein, M. Kollar, and P. Werner, Phys. Rev. Lett. , 056403 (2009).[13] A. Gambassi and P. Calabrese, EPL , 66007 (2011).[14] P. Smacchia, M. Knap, E. Demler, and A. Silva, Phys.Rev. B , 205136 (2015).[15] S. Smale, P. He, B. A. Olsen, K. G. Jackson, H. Sharum,S. Trotzky, J. Marino, A. M. Rey, and J. H. Thywissen,Science Advances (2019), 10.1126/sciadv.aax1568.[16] H.-X. Yang, T. Tian, Y.-B. Yang, L.-Y. Qiu, H.-Y. Liang,A.-J. Chu, C. B. Da˘g, Y. Xu, Y. Liu, and L.-M. Duan,Phys. Rev. A , 013622 (2019).[17] T. Tian, H.-X. Yang, L.-Y. Qiu, H.-Y. Liang, Y.-B. Yang,Y. Xu, and L.-M. Duan, Phys. Rev. Lett. , 043001(2020).[18] J. A. Muniz, D. Barberena, R. J. Lewis-Swan, D. J.Young, J. R. K. Cline, A. M. Rey, and J. K. Thompson,Nature , 602 (2020).[19] J. Zhang, G. Pagano, P. W. Hess, A. Kyprianidis,P. Becker, H. Kaplan, A. V. Gorshkov, Z.-X. Gong, andC. Monroe, Nature , 601 (2017).[20] R. H. Dicke, Phys. Rev. , 99 (1954). [21] B. M. Garraway, Philos. Trans. Royal Soc. A , 1137(2011).[22] P. Kirton, M. M. Roses, J. Keeling, and E. G.Dalla Torre, Advanced Quantum Technologies , 1800043(2019).[23] A. Altland and F. Haake, New Journal of Physics ,073011 (2012).[24] A. Altland and F. Haake, Phys. Rev. Lett. , 073601(2012).[25] W. Buijsman, V. Gritsev, and R. Sprik, Phys. Rev. Lett. , 080601 (2017).[26] J. Ch´avez-Carlos, M. A. Bastarrachea-Magnani,S. Lerma-Hern´andez, and J. G. Hirsch, Phys. Rev. E , 022209 (2016).[27] R. Lewis-Swan, A. Safavi-Naini, J. J. Bollinger, andA. M. Rey, Nature communications , 1 (2019).[28] D. Villase˜nor, S. Pilatowsky-Cameo, M. A. Bastarrachea-Magnani, S. Lerma-Hern´andez, L. F. Santos, and J. G.Hirsch, New Journal of Physics , 063036 (2020).[29] C. Emary and T. Brandes, Phys. Rev. Lett. , 044101(2003).[30] C. Emary and T. Brandes, Phys. Rev. E , 066203(2003).[31] P. P´erez-Fern´andez, A. Rela˜no, J. M. Arias, P. Cejnar,J. Dukelsky, and J. E. Garc´ıa-Ramos, Phys. Rev. E ,046208 (2011).[32] T. Brandes, Phys. Rev. E , 032133 (2013).[33] R. Puebla and A. Rela˜no, EPL (Europhysics Letters) , 50007 (2014).[34] M. A. Bastarrachea-Magnani, S. Lerma-Hern´andez, andJ. G. Hirsch, Phys. Rev. A , 032101 (2014).[35] 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, Phys. Rev.Lett. , 040503 (2018).[36] 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, New Journal ofPhysics , 055013 (2018).[37] K. Baumann, C. Guerlin, F. Brennecke, andT. Esslinger, Nature , 1301 (2010).[38] J. Klinder, H. Keßler, M. Wolke, L. Mathey, and A. Hem-merich, PNAS , 3290 (2015).[39] Z. Zhang, C. H. Lee, R. Kumar, K. J. Arnold, S. J. Mas-son, A. L. Grimsmo, A. S. Parkins, and M. D. Barrett,Phys. Rev. A , 043858 (2018).[40] R. M. Kroeze, Y. Guo, V. D. Vaidya, J. Keeling, andB. L. Lev, Phys. Rev. Lett. , 163601 (2018).[41] R. Puebla, Phys. Rev. B , 220302 (2020).[42] “See Supplemental Material at [URL will be inserted bypublisher].”.[43] A. Safavi-Naini, R. J. Lewis-Swan, J. G. Bohnet,M. Garttner, K. A. Gilmore, E. Jordan, J. Cohn, J. K.Freericks, A. M. Rey, and J. J. Bollinger, “Exploring adi-abatic quantum dynamics of the dicke model in a trappedion quantum simulator,” (2017), arXiv:1711.07392.[44] “Manuscript in preparation.”.[45] S. Strogatz, Nonlinear Dynamics and Chaos: With Appli-cations to Physics, Biology, Chemistry, and Engineering ,Studies in Nonlinearity (Avalon Publishing, 2014).[46] H. B. Kaplan, L. Guo, W. L. Tan, A. De, F. Marquardt,G. Pagano, and C. Monroe, Phys. Rev. Lett. , 120605(2020).[47] T. Kiendl and F. Marquardt, Phys. Rev. Lett. , , 489 (2012).[49] J. G. Bohnet, B. C. Sawyer, J. W. Britton, M. L. Wall,A. M. Rey, M. Foss-Feig, and J. J. Bollinger, Science , 1297 (2016).[50] A. Mortensen, E. Nielsen, T. Matthey, and M. Drewsen,Phys. Rev. Lett. , 103001 (2006).[51] M. K. Joshi, A. Elben, B. Vermersch, T. Brydges,C. Maier, P. Zoller, R. Blatt, and C. F. Roos, Phys.Rev. Lett. , 240505 (2020).[52] K. A. Gilmore, J. G. Bohnet, B. C. Sawyer, J. W. Britton,and J. J. Bollinger, Phys. Rev. Lett. , 263602 (2017).[53] S. M. Brewer, J.-S. Chen, A. M. Hankin, E. R. Clements,C. W. Chou, D. J. Wineland, D. B. Hume, and D. R. Leibrandt, Phys. Rev. Lett. , 033201 (2019).[54] R. B. Sidje, ACM Trans. Math. Softw. , 130 (1998).[55] M. G¨arttner, J. G. Bohnet, M. Safavi-Naini, M. L. Wall,J. J. Bollinger, and A. M. Rey, Nat. Phys. , 781 (2017).[56] Fast oscillations in S z at frequency Ω are overlaid ontop of the slower oscillation envelope. While it is pos-sible, but perhaps workload intensive, to temporally re-solve these in the current experiment, we further pointout that constructing a reliable time-average of S z (andthus distinguishing the phases) does not require one todo so. Importantly, one must only ensure that data is nottaken at, e.g., multiples of the period of these oscillationsto avoid pathological aliasing/sampling problems in thedata.[57] K. Gilmore and M. Affolter, private communication. Supplemental Material: Characterizing the dynamical phase diagram of the Dickemodel via classical and quantum probes
CLASSICAL LIMIT OF DICKE MODEL: MEAN-FIELD EQUATIONS OF MOTION
The mean field equations of motion (EOM) can be derived from the Heisenberg equations of motion by replacingoperators ˆ a and ( ˆ S x , ˆ S y , ˆ S z ) with the classical variables β/ √ N and ( s x , s y , s z )/ N , respectively:˙ β = − iδβ − igs z ˙ s x = − g ( β + β ∗ ) s y ˙ s y = 2 g ( β + β ∗ ) s x − Ω s z ˙ s z = Ω s y . (S1)These should be combined with the initial conditions β (0) = gβ /δ and [ s x (0) , s y (0) , s z (0)] = [sin θ , , − cos θ ], where β is a quantity of order 1 and θ is the initial tipping angle of the spins with respect to the − z direction in the xz plane. Since the initial condition for β is parameter dependent, we define ˜ β = δβ/g and ˜ β (0) = β so that theequations of motion take a more insightful form:˙˜ β = − iδ ( ˜ β + s z )˙ s x = − Ω (cid:20) ˜ g β + ˜ β ∗ ) s y (cid:21) ˙ s y = Ω (cid:20) ˜ g β + ˜ β ∗ ) s x − s z (cid:21) ˙ s z = Ω s y , (S2)where ˜ g = g/g QP T = 2 g/ √ δ Ω is a normalized coupling constant. This form of the EOM allows us to unambiguouslyanalyze the spin- ( δ/ Ω) → ∞ ) and boson-dominated ( δ/ Ω →
0) limits at fixed ˜ g . ANALYTIC TREATMENT OF DPT IN INTEGRABLE REGIMES
In this section we present a detailed outline of the classical particle in a potential description of the DPT in theintegrable spin- and boson-dominated limits. We derive each case separately, drawing from the above mean-field EOM[Eq. (S2)], and present the dependence of the critical point of the DPT on the initial state.
Spin-dominated regime: δ/ Ω (cid:29) First, we address the limit where the spins are the dominant degree of freedom. To derive an effective model ofthe dynamics we take the limit δ/ Ω → ∞ but keeping Ω fixed and formally eliminate the bosonic degree of freedomso we are left with a simpler pure spin description. This is achieved by observing from Eq. (S2) that the bosons willevolve on a much faster timescale (1 /δ ) than the spins (1 / Ω). Thus, the dynamics of the bosons is composed of twodifferent contributions: a slow drift created by the slow time evolution of s z plus very fast oscillations about it. Theslow dynamics is characterized by the locally averaged value ˜ β , which is obtained by setting ˙˜ β = 0 in Eq. (S2):˜ β av = − s z . (S3)As a first approximation, the spin degree of freedom will only be sensitive to the dynamics of ˜ β av , such that thespins evolve under a closed set of equations of motion (EOM):˙ s x = Ω˜ g s z s y ˙ s y = − Ω (cid:0) ˜ g s z s x + s z (cid:1) ˙ s z = Ω s y . (S4)These equations can equivalently be derived from a classical Hamiltonian h = Ω (cid:18) − ˜ g s z + s x (cid:19) , (S5)with the standard angular momentum Poisson brackets: { s i , s j } = (cid:15) ijk s k .The energy h is a conserved quantity defined with respect to the initial condition of the spins, [ s x (0) , s y (0) , s z (0)] =[sin θ , , − cos θ ]: h = Ω (cid:18) − ˜ g θ + sin θ (cid:19) . (S6)Similarly, the total spin length, s x + s y + s z = 1, is also conserved under the dynamics. Together, this allows us toeliminate s x (using conservation of h ) and s y (using the EOM), to obtain a single equation for s z :( ˙ s z )2 + Ω (cid:26)(cid:20) ˜ g s z − cos θ ) + sin θ (cid:21) + s z − (cid:27)(cid:124) (cid:123)(cid:122) (cid:125) V ( s z ) = 0 . (S7)This differential equation describes the motion of a classical particle in a one-dimensional potential V ( s z ). Thedynamics of the particle is dictated by the form of the potential and the initial condition: The particle is initializedat s z (0) = − cos θ with potential energy V ( s z (0)) = 0 and kinetic energy ˙ s z (0) / g →
0, the particle canfreely move from s z = − cos θ to s z = cos θ and back. This describes the untrapped regime. As ˜ g is increased, V ( s z )develops a maximum at s z = 0. The height of this central barrier increases until it matches the initial mechanicalenergy of the particle for a critical coupling ˜ g s DPT ( θ ) = √ (cid:115) θ cos θ . (S8)For ˜ g ≥ ˜ g s DPT ( θ ) the particle is confined to s z ≤ θ = 0, which is the primary focus of the main text, ˜ g s DPT = √ Boson-dominated regime: δ/ Ω (cid:28) Here, we address the limit where the bosons are the dominant degree of freedom. To derive an effective model of thedynamics we take the limit δ/ Ω → δ fixed. Complementary to the prior analysis, the spin degree of freedomnow has a very short evolution timescale (1 / Ω) as compared to the bosons (1 /δ ). Hence, from the perspective of thespins, the ˜ β term in their equations of motion [see Eq. (S2)] is static. Hence, the spin dynamics can be understood - - - - - - - - - FIG. S1. Effective potential for the spin degree of freedom when δ (cid:29) Ω and θ = 0 (left) or θ = π/ s z = − cos( θ ) and then moves to the right. When ˜ g < ˜ g s DPT (0) (blue) the central barrier is below the initial “energy”and so s z oscillates between − cos( θ ) and cos( θ ). When ˜ g > ˜ g s DPT (blue) the particle cannot overcome the central barrier andso s z < θ in this instance), and in the spin-dominated regime this solely introduces the dependenceof the critical coupling ˜ g s DPT ( θ ) on the initial state. as a fast precession with respect to a “static” ˜ β -dependent axis. Furthermore, the slow time evolution of ˜ β inducesadiabatic evolution of the spins on top of their fast precession, with instantaneous axis of rotation given by n t ≡ (cid:113) g ˜ β R (cid:0) , , ˜ g ˜ β R (cid:1) , (S9)where ˜ β R = ( ˜ β + ˜ β ∗ ) / β and n t owes its time dependence to ˜ β R .The drastic difference in time-scales between the spin and bosonic degrees of freedom means that the bosons areonly sensitive to the local average of the spin vector (e.g., we ignore the rapid oscillations on 1 / Ω timescales) s : s av = ( n t =0 · s ) n t , (S10)with s ≡ s (0) the initial condition. For simplicity of analysis in the following, we will take θ = 0 (relevant for themain text). Then, the boson equation of motion becomes˙˜ β = − iδ (cid:32) ˜ β − ˜ g β (cid:112) g β ˜ g ˜ β R (cid:113) g ˜ β R (cid:33) , (S11)which can be recast in terms of an effective potential by eliminating the imaginary part of ˜ β :˙˜ β R δ (cid:32) ˜ β R − β (cid:113) g ˜ β R (cid:112) g β + β − β (cid:19)(cid:124) (cid:123)(cid:122) (cid:125) V ( ˜ β R ) = 0 . (S12)The analysis of this differential equation follows identically to that of the spin-dominated regime. For ˜ g → ± β , which corresponds to the untrapped regime. As ˜ g increases, a maximum developsat ˜ β R = 0 and eventually acts as a barrier that keeps β > g b DPT ( β ) = (cid:20) − β β (2 − β ) (cid:21) / . (S13)For β = 1, we recover ˜ g DPT = 3 / as per the result quoted in the main text. An important difference with the spinpotential is that in this case the form of the effective boson potential depends on the initial state of the spins as wellas the initial state of the bosons, whereas the spin potential is insensitive to the initial state of the bosons. While ourresult here, Eq. (S13), specializes to θ = 0, we emphasize that the calculation is easily generalized for arbitrary θ . COMPARISON TO EQPT
The Dicke model has been intensely studied with regards to the presence of a transition in the excited-state spectrumof the Hamiltonian (EQPT). Not only was the nature of the EQPT subject to initial debate [31–33], although it isnow accepted as a true quantum phase transition, but also there was initially speculation as to whether the EQPTmight be connected to various non-equilibrium features of the model such as the emergence of chaos in the energyspectrum [26, 31].The EQPT occurs for g > g
QPT (e.g., when the coupling is large enough that the ground-state corresponds to thesuperradiant phase) at a critical energy E c = − Ω N/
2. It is signaled by a non-analyticity in the density of states, inparticular a logarithmic divergence in the first derivative [32], as well as singularities in (cid:104) ˆ S x (cid:105) and (cid:104) ˆ a † ˆ a (cid:105) computed withrespect to the energy eigenstates [31]. The EQPT can be explicitly connected to the classical model in the large- N limit by noting that it arises due to the emergence of a saddle point in the classical phase-space [32, 34].The DPT we have studied appears to be quantitatively connected to the EQPT in a certain limit. In particular,it is straightforward to show that the DPT in the spin-dominated regime coincides exactly with the EQPT: We findthat the energy of the initial state, E = (cid:104) ψ (0) | ˆ H D | ψ (0) (cid:105) = − g N/δ , coincides with the critical energy E c = − Ω N/ g DPT = (cid:112) δ Ω /
2. Indeed, going a step further, it is straightforward to prove that the DPT in the spin-dominatedregime always precisely occurs when E = E c for any initial state (e.g, arbitrary θ , as covered in the previoustreatment of the DPT by the particle in a potential model). This close correspondence is consistent with prior workin the literature [33], which suggested that EQPTs might be accessible through a type of DPT driven by the breakingof the parity symmetry in the energy spectrum. Moreover, there has been a recent experimental demonstration [17]connecting a DPT to an EQPT in a spinor Bose-Einstein condensate.Nevertheless, we highlight that this coincidence is not universal in the Dicke model. Specifically, the boson-dominated regime provides conclusive evidence that the DPT and EQPT are, in general, disconnected. To make thisabundantly clear: In the boson-dominated regime the energy of the prototypical initial state evaluated at the DPT isgiven by E | ˜ g b DPT = − g N/δ | g =3 / √ δ Ω / = −√ N/ > E c . NUMERICAL SIMULATION OF QUANTUM DYNAMICS
We efficiently simulate the quantum dynamics of the Dicke model by solving the time-dependent Schr¨odingerequation using a Krylov-subspace projection method [54]. Using the symmetries of the system, we are able tomake our numerical calculation tractable by: i) representing the spin degree of freedom in the fully symmetricsubspace of Dicke states | S, m z (cid:105) with S = N/
2, defined such that ˆ S | S, m z (cid:105) = S ( S + 1) | S, m z (cid:105) for ˆ S = (cid:80) α = x,y,z ˆ S α and ˆ S z | S, m z (cid:105) = m z | S, m z (cid:105) , and ii) representing the boson degree of freedom in a truncated Fock basis | n (cid:105) where n ∈ [0 , n max ] for some maximum occupation n max . While it is straightforward to simulate the evolution using a fixed n max , this may require extensive trial runs and benchmarking to ensure convergence of all relevant observables.To avoid this issue, and to efficiently simulate large system sizes to long times over a wide range of parameterspace, we utilize a bosonic Hilbert space that dynamically grows or shrinks throughout the evolution. Initially, we set n max by bounding the error in the initial wavefunction normalization 1 − |(cid:104) ψ (0) | ψ (0) (cid:105)| ≤ (cid:15) = 10 − ; by checking forconvergence as (cid:15) is varied, we have verified that this bound is sufficient to generate relative errors in the spin-bosonentanglement entropy dynamics on the order of 10 − at long times, and even smaller errors for relevant observables.Throughout the evolution, we periodically check (on a time interval τ ) the normalization of our wavefunction whenprojected onto states with boson occupations of n max and n max −
1. If the error introduced by the projection ontoeither of these sets of states exceeds (cid:15) at any point in this interval, we increase n max by a fixed amount ∆ n andrecompute the evolution over this interval, repeating until the error tolerance is satisfied. In general, we find that theshort-time dynamics requires a much larger value of n max than the initial state; in contrast, however, the late-timedynamics typically require a much smaller value of n max . To speed up the computation of the late-time dynamics,we also allow for a decrease in the value of n max , likewise checking the projection onto all boson occupations between n max − ∆ n and n max . If the projection onto this set of states is less than (cid:15) at any point over some interval τ , wereduce n max by ∆ n , and truncate away the corresponding states from our Hilbert space.For all simulations in Fig. 4(b)-(d), we use N = 600 and evolve to a fixed time gt = 300, initializing the system inthe state | ( − N/ z (cid:105) ⊗ | α s (cid:105) with α s = ( g/δ ) √ N as described in the main text. In Fig. S2, we show typical dynamics ofthe entanglement entropy, S vN = − Tr[ˆ ρ s log(ˆ ρ s )] where ˆ ρ s is the reduced density matrix of the spins, magnetization (cid:104) ˆ S z (cid:105) , and variance (cid:104) (∆ ˆ S z ) (cid:105) = (cid:104) ˆ S z (cid:105) − (cid:104) ˆ S z (cid:105) for various parameter regimes. To produce Fig. 4 of the main manuscript,we simulate the Dicke model for a range of 51 values of ˜ g and 39 values of δ/ Ω, equally-spaced between the bounds
FIG. S2. Typical time-traces of the entanglement entropy S vN (top row), magnetization (cid:104) ˆ S z (cid:105) (middle row), and variance (cid:104) (∆ ˆ S z ) (cid:105) (bottom row) for the quantum Dicke model with N = 600, up to time gt = 300. Columns correspond to differentvalues of δ/ Ω, with δ/ Ω = 0 . . . g = 1 .
1, corresponding to the untrapped phase, whilered curves are for ˜ g = 1 .
6, corresponding to the trapped phase. Black curves are selected for ˜ g occurring at the DPT for thecorresponding values of δ/ Ω: for δ/ Ω = 0 .
5, 1 .
2, and 4 .
0, we have ˜ g = 1 .
32, 1 .
29, and 1 .
43, respectively. The dashed linesplotted alongside the entropy and variance are the corresponding fits to the function A (1 − e − γt ) (used to extract the timescale γ for the variance in Fig. 4d of the main text), while the dotted lines plotted alongside the variance are the corresponding fitsto the translated logistic function (see text). shown on the plot. For a small set of parameters, typically with δ/ Ω (cid:46) . g (cid:38) . A (1 − e − γt )with free fitting parameters A (steady state value) and γ , shown as dashed lines in Fig. S2. We generally find that theentanglement entropy S vN ( t ) is consistent with such a functional form, with the exception of transient oscillations onfast time-scales. For some parameters we find that the early time dynamics of the spin variance (cid:104) (∆ ˆ S z ) (cid:105) is capturedmore accurately by a translated logistic function, a/ (1 + e − c ( t − t ) ) − b , than an exponential function; here, a providesthe overall scale of the function to be consistent with the observed steady state, c is the logistic growth rate, t is thefunction midpoint, and b is an offset to allow for the fact that (cid:104) (∆ ˆ S z ) (cid:105) = 0 at t = 0. Fitting to this function providesan alternative estimate for the relaxation timescale of the dynamics via γ logistic = c + 1 /t , where we account for boththe logistic growth timescale as well as the delay in the onset of this growth.In Fig. S3, we plot various notions of the relaxation timescales for the dynamics, including the value of γ obtainedby an exponential fit to the spin variance, which is plotted in Fig. 4(d) of the main text. Despite the complexnature of the dynamics over the wide range of behaviors we observe, we find that these timescales – obtained byfits to relatively simple functions – are generally consistent with each other, and thus serve as reliably proxies insystematically comparing the dynamics over the considered parameter space.Lastly, in Fig. S4, we also provide results for N = 200 and N = 40, analogous to those for N = 600 shown in Fig. 4of the main manuscript, with the addition of the steady-state magnetization (cid:104) ˆ S z (cid:105) . For N = 40, the timescale γ isextracted via an empirical fit of the entanglement entropy S vN ( t ) ∝ (1 − e − γt ), as opposed to the spin fluctuations (cid:104) (∆ ˆ S z ) (cid:105) , which is used for N = 200 and N = 600; the spin fluctuations for smaller systems exhibit noisy oscillationsthat prevent a consistent estimate of a relaxation timescale. While signatures of the DPT are still evident in the FIG. S3. Comparison of relaxation timescales for the quantum dynamics with N = 600. Panels (a) and (b) correspond totimescales γ obtained through fits of the entropy S vN and the variance (cid:104) (∆ ˆ S z ) (cid:105) , respectively, to the exponential function A (1 − e − γt ); panel (b) is identical to Fig. 4(d) of the main text. Panel (c) corresponds to the timescale γ logistic obtainedthrough a fit to a shifted logistic function (see text). single-body observables for systems as small as N = 40, features of both the DPT and classical chaos are overshadowedby the associated quantum noise in this system. However, for N = 200, which is within the capabilities of currenttrapped-ion platforms, signatures of classical chaos and the DPT feature much more prominently. EXPERIMENTAL REALIZATION IN A TRAPPED ION QUANTUM SIMULATOR
In the main text we discuss how the dynamical phases we report can be realized in an implemenation of the Dickemodel in a 2D trapped ion array formed in a Penning trap experiment [35, 36]. Here, we support these commentswith a brief comparison of the achievable coherent time-scales relative to decoherence and technical limitations.First, we establish the time-scales of the dynamics we are most interested in. To approach this problem concretely,we will focus on the limiting cases of the boson- and spin-dominated for which we can identify the relevant energy-scales from our prior analytic treatment in the mean-field limit. We expect quantum corrections to predominantly leadto, e.g., damping of oscillations due to the quantum noise. The resonant regime interpolates between the integrablespin and boson regimes, and thus should be suitably estimated from our results here.The key requirement of the experiment will be to distinguish the untrapped from trapped phase. In currentexperiments measurement of the phonon modes is difficult, and so we focus on the oscillations in the spin degree offreedom, S z (note that both S z and X are always equally suitable observables to diagnose the DPT). In particular,we must be able to observe clear oscillations in S z when ˜ g < ˜ g DPT . This implies that the frequency of the oscillationsmust be large compared to relevant decoherence. In the 2D Penning trap of Refs. [35, 36, 43, 49, 55] the dominantsource of decoherence is single-particle elastic dephasing of the spins at a rate Γ el . In comparison, the frequency scaleof the relevant oscillations in S z in the boson-dominated regime is given by δ [56], whereas in the spin-dominated(LMG) regime the oscillation frequency is Ω. Thus we require δ/ Γ el (cid:38) δ/ Ω (cid:28) / Γ el (cid:38) δ/ Ω (cid:28) g = 1, as this is the characteristic scale of the DPT in both the spin- and boson-dominated regimes.Given a fixed ratio η = δ/ Ω we then have from the condition ˜ g = 1 that Ω = 2 g/ √ η or δ = 2 g √ η . Finally, this allowsus to rewrite the condition on our oscillation frequency scales as δ Γ el = 2 g √ η Γ el for η (cid:28) , (S14)ΩΓ el = 2 g √ η Γ el for η (cid:29) . (S15)To connect to the experimental system we consider a pair of accessible spin-boson coupling strengths, 2 g/ (2 π ) =(1 . , .
7) kHz [57]. These values are similar to that reported in Ref. [35]. The decoherence rate Γ el is correlated withthe spin-boson coupling strength, as it arises due to spontaneous emission from the applied off-resonant ODF lasers[49], and for the above values of g we correspondingly have Γ el ≈ (350 , − . Thus, from an example calculationwith 2 g/ (2 π ) = 3 . δ/ Γ el = Ω / Γ el ≈
13 at η = 0 . η = 10 respectively. Both values demonstrate that,in principle, the dynamical phases should be accessible across the full range of η in the current trapped ion simulator. FIG. S4. Comparison of quantum dynamics for N = 200 (a-d) and N = 40 (e-h). Panels (a)/(e) corresponds to the averagemagnetization, S z = (1 /T ) (cid:82) T + t t dt (cid:104) ˆ S z ( t ) (cid:105) . (b)/(f) Time-averaged entanglement entropy S vN = (1 /T ) (cid:82) T + t t dtS vN ( t ) for S vN = − Tr[ˆ ρ s log(ˆ ρ s )] where ˆ ρ s is the reduced density matrix of the spins, analogous to Fig. 4(b) of the main manuscript. (c)/(g)Temporal fluctuations of entanglement entropy, (∆ S vN ) = (1 /T ) (cid:82) T + t t dt [ S vN ( t ) − S vN ] , normalized by S vN , analogous toFig. 4(c) of the main manuscript. (d)/(h) Relaxation timescale γ of the dynamics. For panel (d), this is obtained by anempirical fit to (cid:104) (∆ ˆ S z ) (cid:105) ∝ (1 − e − γt ), analogous to Fig. 4(d) of the main manuscript; for panel (h), this is obtained by anempirical fit to S vN ∝ (1 − e − γt ). Time averages in (a)-(c), (e)-(g) are computed over a window t ∈ [ t , T + t ] for gt , gT = 150,so as to reduce transient effects in steady state estimates. A more systematic and careful study is presented in Fig. S5, where we also indicate the role of other technicallimitations. In particular, the available range of detunings δ is limited by: i) coupling to other normal modes of thecrystal and ii) fluctuations of the COM frequency related to impurity ions in the crystal and thermal occupation ofthe in-plane motional modes [57]. The former places an approximate limit of δ/ (2 π ) (cid:46)
10 kHz while the latter boundsthe detuning from below as δ/ (2 π ) (cid:38)
50 Hz. One might expect COM fluctuations could be impactful in regimesof η (cid:28) δ is favoured, while coupling to other modes will be relevant for η (cid:29) δ is typicallylarger. Examining the predicted frequency scales in Fig. S5, we find that COM fluctuations should not be an issuebut coupling to other modes can become problematic for η (cid:29) g is set too large. In particular, one should favour2 g/ (2 π ) = (1 .
88) kHz to ensure δ/ (2 π ) <
10 kHz at the largest values of η . For completeness, the transverse fieldtakes values in the range 0 . (cid:46) Ω / (2 π ) (cid:46)
12 kHz in both plots, which should be readily accessible. -1 -1 FIG. S5. Comparison of coherent and dissipative time-scales for experimental system based on Ref. [35]. Solid lines indicate δ/ Γ el for boson-dominated regime ( δ/ Ω (cid:28) / Γ el for spin-dominated regime ( δ/ Ω (cid:29) δ/ Γ el : Upper red shaded area indicates detunings that may introduce couplingsto other normal modes of the crystal, δ/ (2 π ) >
10 kHz, while lower shaded area indicates detunings below the limit imposedby COM fluctuations δ/ (2 π ) <<