Relaxation timescales and prethermalisation in d-dimensional long-range quantum spin models
RRelaxation timescales and prethermalisation in d -dimensional long-range quantumspin models Michael Kastner
1, 2, ∗ and Mauritz van den Worm
1, 2 Institute of Theoretical Physics, University of Stellenbosch, Stellenbosch 7600, South Africa National Institute for Theoretical Physics (NITheP), Stellenbosch 7600, South Africa
We report analytic results for the correlation functions of long-range quantum Ising models inarbitrary dimension. In particular, we focus on the long-time evolution and the relevant timescaleson which correlations relax to their equilibrium values. By deriving upper bounds on the correlationfunctions in the large-system limit, we prove that a wide separation of timescales, accompanied bya pronounced prethermalisation plateau, occurs for sufficiently long-ranged interactions.
Empirically it is well established that, after a suf-ficiently long time, most physical many-body systems,whether isolated or coupled to an environment, will equi-librate. In many cases the equilibrium is well describedby a Gibbs state, and this observation is at the basisof equilibrium statistical mechanics. An understandingof the microscopic process leading to thermalisation is,however, still incomplete. Recent experiments with coldatoms, ions, and molecules [1–8] have sparked a revivalof interest in questions related to the foundations of equi-librium statistical mechanics. Substantial progress, oftenbased on typicality techniques, has been made on thetheoretical side in the past few years. Results in vari-ous physical settings have been reported, proving thatequilibration and/or thermalisation takes place for typi-cal quantum systems of sufficient size [9–14].While these results establish that equilibration and/orthermalisation will happen eventually, the time scale ofsuch a relaxation process remains unspecified. Only veryrecently has it become apparent that typicality tech-niques can also be applied for analysing the time scalesof equilibration. The results of these efforts, while pio-neering a promising approach, are not yet fully satisfac-tory, as they either over- [15] or—in a different setting—underestimate [16] the timescales by many orders of mag-nitude. These results indicate that physically realisticmodels or observables are not “typical” in the mathe-matical sense.In this paper we approach the problem of equilibrationtime scales from a different angle, reporting the results ofa model study. This work is a continuation and extensionof a previous paper [17] where exact, analytical expres-sions for equal-time two-point correlation functions havebeen computed for long-range interacting Ising models ina longitudinal magnetic field. While these results are verygeneral and exact, the long-time asymptotics relevant forthe relaxation to equilibrium is not at all obvious fromthe analytical expressions. Motivated by recent ion-trapexperiments where the Ising model with long-range in-teractions can be realised, we have derived in [17] upperbounds on the time-evolution of the various spin–spin ∗ [email protected] correlation functions of the two-dimensional long-rangeIsing model on a triangular lattice. For the long-rangeIsing model on this specific lattice, the asymptotic long-time behaviour can be read off easily from the upperbounds. In the present paper we report generalisationsof these upper bounds to arbitrary regular lattices anddimensionality.We consider long-ranged coupling constants J i,j ∝| i − j | − α decaying like a power law with the distance | i − j | between lattice sites i and j . The exponent α ≥ α = 0 to nearest-neighbourinteractions in the limit α → ∞ . The upper bounds de-rived in this paper are stretched or compressed exponen-tials in time for all values of α . On a d -dimensional latticeand for α < d/
2, we find that some of the spin–spin cor-relation functions relax to their equilibrium values in atwo-step process, governed by two widely separate timescales, while single-spin expectation values relax alreadyon the faster of the two time scales. This kind of be-haviour, characterised by a long-lived quasi-stationarystate in which only some of the expectation values havealready relaxed to their equilibrium value goes under thename of “prethermalisation” and has been discussed ex-tensively in the past few years [18–25].The stretched or compressed exponentials that upper-bound the spin–spin correlation functions do not onlydepend on the exponent α , but they do so in a nonan-alytic, transition-like manner: the long-time asymptoticbehaviour of the spin–spin correlation functions switchesfrom one kind of behaviour to a different functional format the values α = d/ α = d/ −
1. The first of thesethreshold values was already discussed in [17] for the tri-angular lattice in two dimensions, and it is related to theoccurrence of widely separated time scales and prether-malisation. The second threshold value, at α = d/ − α -value) only for lattice dimensions d ≥ I. LONG-RANGE ISING MODEL
Consider a lattice Λ consisting of N sites, to each ofwhich is assigned a spin-1 / a r X i v : . [ c ond - m a t . s t a t - m ec h ] O c t is modelled by a Hilbert space H i = C , and the com-posite system is described by the Ising-type Hamiltonian H (cid:96) = − (cid:88) i ∈ Λ (cid:88) j ∈ Λ \{ i } J i,j σ zi σ zj − h (cid:88) i ∈ Λ σ zi , (1)on the tensor product Hilbert space H = (cid:78) i ∈ Λ H i . Theparameter h ∈ R denotes the magnitude of a homoge-neous external magnetic field in the z direction, and σ zi is the z -Pauli operator acting on the i th component ofthe tensor product space H . At this point the couplingstrengths J i,j between lattice sites i, j ∈ Λ are arbitrary,but we will specialise to power law-decaying long-rangeinteractions at a later stage.As initial states we choose density operators ρ that arediagonal matrices in the σ x tensor-product eigenbasis, ρ = 12 N (cid:18) + (cid:88) i σ xi (cid:18) s xi + (cid:88) j>i σ xj (cid:18) s xxij + (cid:88) k>j σ xk (cid:18) s xxxijk + (cid:88) l>k · · · (cid:19)(cid:19)(cid:19)(cid:19) , (2)where denotes the identity operator on H . The indices i , j , k , l in (2) are summed over the lattice sites. Thischoice of ρ has been exploited previously [17, 26–28], asit leads to particularly simple calculations and results,although generalisations are possible. Starting from ρ ,exact analytic expressions have been reported in [17] forthe time-evolution of the various spin–spin correlationfunctions, e.g. (cid:104) σ xi σ xj (cid:105) ( t ) := Tr (cid:0) e i H (cid:96) t σ xi σ xj e − i H (cid:96) t ρ (cid:1) = (cid:104) σ xi σ xj (cid:105) (0) (cid:0) P − i,j + P + i,j (cid:1) (3) with P ± i,j = 12 (cid:89) k ∈ Λ \{ i,j } cos [2 ( J k,i ± J k,j ) t ] , (4)where we have set h = 0 for simplicity. Other correlationfunctions either behave similarly (like (cid:104) σ yi σ yj (cid:105) ), or simpler(like (cid:104) σ yi σ zj (cid:105) ), or are constant in time (like (cid:104) σ zi σ zj (cid:105) ); see[17] for details. In the following we restrict the presen-tation to the xx -correlation function given in (3), butsimilar techniques can be applied to other correlationfunctions. Upper bounds on one-point functions like (cid:104) σ xi (cid:105) have been reported in the Supplemental Material accom-panying Ref. [29]. II. UPPER BOUND ON THE CORRELATIONS
The expressions in (4) are quasi-periodic in t , and itis therefore precluded that (cid:104) σ xi σ xj (cid:105) ( t ) converges in thelong-time limit for any finite lattice. Only in the thermo-dynamic limit of an infinite number of lattice sites do wehave a chance of observing convergence towards an equi-librium value. To derive such a result, we consider regular d -dimensional lattices. Without loss of generality we con-sider the lattice constants normalised such that there ison average one lattice site per unit (hyper)volume in thelimit of large lattice size. We choose coupling constantsdecaying like a power-law with the Euclidean distance | i − j | between lattice sites i and j , J i,j = J | i − j | α . (5)Without loss of generality we set J = 1. Under theseconditions and in the limit of large lattice size, we obtainthe bounds (cid:12)(cid:12) P − i,j (cid:12)(cid:12) ≤ P − i,j := 12 exp (cid:20) − C − α,d (cid:16) α | i − j | tπ (cid:17) N − α +1) /d (cid:21) for 0 ≤ α < d/ − , exp (cid:20) − C − α,d (cid:16) α | i − j | tπ (cid:17) d/ ( α +1) (cid:21) for α > d/ − , (6a) (cid:12)(cid:12) P + i,j (cid:12)(cid:12) ≤ P + i,j := 12 exp (cid:104) − C + α,d (cid:0) tπ (cid:1) N − α/d (cid:105) for 0 ≤ α < d/ , exp (cid:104) − C + α,d (cid:0) tπ (cid:1) d/α (cid:105) for α > d/ , (6b)with positive constants C + α,d and C − α,d as defined in (23).The proof of the inequalities (6a) and (6b) is postponedto Sec. III.Depending on the value of α , the bounds (6a) and (6b)decay like stretched or compressed exponentials. By nu-merically evaluating the exact expressions (4), we find that the functional form of the bound agrees well, al-though the numerical constants in the bound overesti-mate, as expected, the exact values. From Eqs. (6a) and(6b) one can read off that there are three different regimesof α -values, each with a different relaxation or scaling be-haviour: FIG. 1. Plots of the bound P − i,j + P + i,j on the normalised correlator (cid:104) σ xi σ xj (cid:105) ( t ) / (cid:104) σ xi σ xj (cid:105) (0) as a function of time. The examplesare for lattice dimension d = 3, and the α -values are chosen such that they furnish examples for the three different regimes ofrelaxation behaviour, as discussed in the text. ≤ α < d/ − (cid:12)(cid:12) P − i,j (cid:12)(cid:12) and (cid:12)(cid:12) P + i,j (cid:12)(cid:12) both decay like aGaussian in time, and both do so on time scales thatare N -dependent. The two time scales of relaxation arewidely separated, with (cid:12)(cid:12) P − i,j (cid:12)(cid:12) decaying much slower than (cid:12)(cid:12) P + i,j (cid:12)(cid:12) . The form of the resulting upper bound on (cid:104) σ xi σ xj (cid:105) is shown in Fig. 1a. This regime occurs for positive α only in lattice dimensions d ≥ d/ − < α < d/
2: Again, relaxation takes place in atwo-step process with widely separated time scales. Thefast process described by (cid:12)(cid:12) P + i,j (cid:12)(cid:12) still decays like a Gaussianin time, on a time scale that is N -dependent. The slowtime scale corresponding to (cid:12)(cid:12) P − i,j (cid:12)(cid:12) is independent of thesystem size, with a decay in the form of a compressedexponential. The form of the resulting upper bound on (cid:104) σ xi σ xj (cid:105) is shown in Fig. 1b. α > d/
2: Both terms (cid:12)(cid:12) P − i,j (cid:12)(cid:12) and (cid:12)(cid:12) P + i,j (cid:12)(cid:12) decay to zerolike stretched or compressed exponentials. Relaxationtakes place in a single step, as both relevant time scalesare very similar and independent of N . The form of theresulting upper bound on (cid:104) σ xi σ xj (cid:105) is shown in Fig. 1c.Figs. 2 and 3 show further graphical representationsof the bound, highlighting in particular the qualitativechanges that occur upon variation of the exponent α . FIG. 2. For fixed instances of time, the bound P − i,j + P + i,j isshown as a function of the exponent α . The example is fordimension d = 3 and a lattice of N = 10 × ×
10 sites.The different shaded regions correspond to the three rangesof α -values discussed in the text. On the basis of these results, physical properties ofthe model—including dephasing dynamics, prethermal-isation, and others—can be discussed similarly to thetwo-dimensional case reported in [17], and the reader isreferred to this reference for details.
III. PROOF OF EQS. (6a) AND (6b)
The starting point for the derivation is the exact ex-pression (4), where P ± i,j is given as a product (over alllattice sites) of cosine terms. Since | cos x | ≤
1, we canupper bound the absolute value of this quantity by aproduct over a subset of lattice sites, P ± i,j ≤ (cid:89) k ∈ Λ \ g ± i,j ( t ) cos [2 ( J k,i ± J k,j ) t ] , (7)where we have defined g ± i,j ( t ) := (cid:110) k ∈ Λ : | J k,i ± J k,j ) t | ≥ π (cid:111) . (8) FIG. 3. Contour plot of the bound P − i,j + P + i,j as a functionof t and α for a three-dimensional lattice of N = 10 × × α -valuesdiscussed in the text. FIG. 4. Sketch of the regions g ± i,j ( t ) and h i,j and of the chosencoordinate system as used for the proof in Sec. III. This subset is chosen such that, for all k ∈ Λ \ g ± i,j ( t ), wecan make use of the inequality | cos( πx ) | ≤ − x ≤ exp (cid:0) − x (cid:1) , (9)valid for all | πx | <
2, to write P ± i,j ≤
12 exp (cid:34) − (cid:18) tπ (cid:19) (cid:88) k ∈ Λ \ g ± i,j ( t ) ( J k,i ± J k,j ) (cid:35) . (10)We restrict the k -summation even further by excludingthe hyperslab h i,j sketched in Fig. 4 [30]. The occurrenceof Euclidean distances in the couplings J i,k and J j,k thensuggests to parametrise the lattice sites k ∈ Λ by hyper-spherical coordinates, k ( r, φ , . . . , φ d − ) = r cos φ r sin φ cos φ ... r sin φ · · · sin φ d − cos φ d − r sin φ · · · sin φ d − sin φ d − , (11)with the origin of the coordinate system placed at latticesite i and the z -axis chosen along the line connecting i and j . The couplings can then be written as J i,k = r − α , J j,k = (cid:0) r + 2 rδ cos φ + δ (cid:1) − α/ , (12)where δ = | i − j | denotes the distance between sites i and j . It is then convenient to further restrict the k -summation in (10) by excluding the hyperslab sketchedin Fig. 4. Exploiting also the reflection symmetry of the problem, we arrive at the bound P ± i,j ≤
12 exp (cid:32) − t π (cid:88) k ∈ h i,j \ g ± i,j ( t ) ( J k,i ± J k,j ) (cid:33) (13)with the k -summation restricted to the lattice sites in thehalf plane h i,j = (cid:110) k ( r, φ , . . . , φ d − ) ∈ Λ : 0 ≤ φ ≤ π (cid:111) . (14)For large lattices, we can bound the sum in (13) by anintegral, (cid:88) k ∈ h i,j \ g ± i,j ( t ) ( J k,i ± J k,j ) ≥ π K ( d ) (cid:90) N /d R ± ( t ) d r r d − × (cid:90) π/ d φ cos φ sin d − φ ( J k,i ± J k,j ) , (15)where the prefactor K ( d ) = π for d = 2 , (cid:81) d − m =2 √ π Γ ( ( − m + d ) ) Γ ( (1 − m + d ) ) for d ≥ , (16)originates from the integrations over φ , . . . , φ d − .The lower limit R ± of the r -integration still needs tobe determined such that the region g ± i,j is excluded. I.e.,we need to determine R ± such that | t ( J k,i ± J k,j ) | < π r ≥ R ± ( t ). We are interested in the long-timeasymptotic behaviour, and in this limit large values of r are required to satisfy the above inequality. Hence wecan assume that r is much larger than δ and expand J k,i ± J k,j = 1 r α ± (cid:112) δ + 2 rδ cos φ + r α/ ∼ r − α ± r − α (cid:18) αδ cos φ r (cid:19) (18)to leading order in the small parameter δ/r , yielding J k,i + J k,j ∼ r − α , (19a) J k,i − J k,j ∼ − αδ cos φ r α +1 . (19b)Inserting these asymptotic expressions into (18), we ob-tain R + ( t ) ∼ (cid:18) tπ (cid:19) /α , R − ( t ) ∼ (cid:18) αδtπ (cid:19) / (1+ α ) , (20)valid for sufficiently large t .For similar reasons we can insert the expansions (19a)and (19b) into the integrand of (15). The integrationsbecome elementary in this case, yielding (cid:88) k ∈ h i,j \ g ± i,j ( t ) ( J i,k + J k,j ) ≥ π K ( d )( d − d − α ) (cid:34) N − α/d − (cid:18) tπ (cid:19) d/α − (cid:35) , (21a) (cid:88) k ∈ h i,j \ g ± i,j ( t ) ( J i,k − J k,j ) ≥ π K ( d ) α δ ( d − d − α − (cid:34) N − α +1) /d − (cid:18) αδtπ (cid:19) d/ ( α +1) − (cid:35) , (21b)in the limit of large N and t . Interestingly, the bound in (21a) is independent of the distance δ between the latticesites.Depending on the sign of the exponents 1 − α/d and 1 − α + 1) /d in the N -terms, either the first or the secondterm in the square brackets of (21a) and (21b) will give the dominant contribution in the limit of large system size N . As a result, the asymptotic behaviour of the bounds is different for different ranges of α , (cid:88) k ∈ h i,j \ g ± i,j ( t ) ( J i,k + J k,j ) ≥ π K ( d )( d − d − α ) N − α/d for 0 ≤ α < d/ , − (cid:0) tπ (cid:1) d/α − for α > d/ , (22a) (cid:88) k ∈ h i,j \ g ± i,j ( t ) ( J i,k − J k,j ) ≥ π K ( d ) α δ ( d − d − α − N − α +1) /d for 0 ≤ α < d/ − , − (cid:0) αδtπ (cid:1) d/ ( α +1) − for α > d/ − . (22b)Inserting these expressions into the inequality (10) and defining the positive constants C + α,d = π K ( d )( d − | d − α | , C − α,d = 2 π K ( d )( d − | d − α − | , (23)the main results (6a) and (6b) of the paper follow. [1] T. Kinoshita, T. Wenger, and D. S. Weiss, Nature ,900 (2006).[2] S. Hofferberth, I. Lesanovsky, B. Fischer, T. Schumm,and J. Schmiedmayer, Nature , 324 (2007).[3] J. Struck, C. ¨Olschl¨ager, R. L. Targat, P. Soltan-Panahi,A. Eckardt, M. Lewenstein, P. Windpassinger, andK. Sengstock, Science , 996 (2011).[4] R. Islam, E. E. Edwards, K. Kim, S. Korenblit, C. Noh,H. Carmichael, G.-D. Lin, L.-M. Duan, C.-C. J. Wang,J. Freericks, and C. Monroe, Nat. Commun. , 377(2011).[5] S. Trotzky, Y.-A. Chen, A. Flesch, I. P. McCulloch,U. Schollw¨ock, J. Eisert, and I. Bloch, Nature Phys. , 325 (2012).[6] J. W. Britton, B. C. Sawyer, A. C. Keith, C.-C. J.Wang, J. K. Freericks, H. Uys, M. J. Biercuk, and J. J.Bollinger, Nature , 489 (2012).[7] R. Islam, C. Senko, W. C. Campbell, S. Korenblit,J. Smith, A. Lee, E. E. Edwards, C.-C. J. Wang, J. K.Freericks, and C. Monroe, Science , 583 (2013).[8] B. Yan, S. A. Moses, B. Gadway, J. P. Covey, K. R. A.Hazzard, A. M. Rey, D. S. Jin, and J. Ye, Nature ,521 (2013).[9] P. Reimann, Phys. Rev. Lett. , 190403 (2008).[10] N. Linden, S. Popescu, A. J. Short, and A. Winter, Phys.Rev. E , 061103 (2009). [11] S. Goldstein, J. L. Lebowitz, C. Mastrodonato, R. Tu-mulka, and N. Zangh`ı, Phys. Rev. E , 011109 (2010).[12] P. Reimann, New J. Phys. , 055027 (2010).[13] P. Reimann and M. Kastner, New J. Phys. , 043020(2012).[14] C. Ududec, N. Wiebe, and J. Emerson, Phys. Rev. Lett , 080403 (2013).[15] A. J. Short and T. C. Farrelly, New J. Phys. , 013063(2012).[16] S. Goldstein, T. Hara, and H. Tasaki, New J. Phys. ,045002 (2015).[17] M. van den Worm, B. C. Sawyer, J. J. Bollinger, andM. Kastner, New J. Phys. , 083007 (2013).[18] J. Berges, S. Bors´anyi, and C. Wetterich, Phys. Rev.Lett. , 142002 (2004).[19] M. Moeckel and S. Kehrein, Phys. Rev. Lett. , 175702(2008).[20] M. Eckstein, M. Kollar, and P. Werner, Phys. Rev. Lett. , 056403 (2009).[21] M. Gring, M. Kuhnert, T. Langen, T. Kitagawa,B. Rauer, M. Schreitl, I. Mazets, D. A. Smith, E. Demler,and J. Schmiedmayer, Science , 1318 (2012).[22] J. Marino and A. Silva, Phys. Rev. B , 060408 (2012).[23] D. A. Smith, M. Gring, T. Langen, M. Kuhnert,B. Rauer, R. Geiger, T. Kitagawa, I. Mazets, E. Demler,and J. Schmiedmayer, New J. Phys. , 075011 (2013). [24] Z.-X. Gong and L. M. Duan, New J. Phys. , 113051(2013).[25] M. Marcuzzi, J. Marino, A. Gambassi, and A. Silva,Phys. Rev. Lett. , 197203 (2013).[26] G. G. Emch, J. Math. Phys. , 1198 (1966). [27] M. Kastner, Phys. Rev. Lett , 130601 (2011).[28] M. Kastner, Cent. Eur. J. Phys. , 637 (2012).[29] R. Bachelard and M. Kastner, Phys. Rev. Lett.110