Ergodicity Breaking Transition in Finite Disordered Spin Chains
EErgodicity Breaking Transition in Finite Disordered Spin Chains
Jan ˇSuntajs, Janez Bonˇca,
2, 1
Tomaˇz Prosen, and Lev Vidmar
1, 2 Department of Theoretical Physics, J. Stefan Institute, SI-1000 Ljubljana, Slovenia Department of Physics, Faculty of Mathematics and Physics,University of Ljubljana, SI-1000 Ljubljana, Slovenia
We study disorder-induced ergodicity breaking transition in high-energy eigenstates of interactingspin-1/2 chains. Using exact diagonalization we introduce a cost function approach to quantitativelycompare different scenarios for the eigenstate transition. We study ergodicity indicators such as theeigenstate entanglement entropy and the spectral level spacing ratio, and we consistently find thatan (infinite-order) Kosterlitz-Thouless transition yields a lower cost function when compared to afinite-order transition. Interestingly, we observe that the transition point in finite systems exhibitsnearly thermal properties, i.e., ergodicity indicators at the transition are close to the random matrixtheory predictions.
Introduction.
Generic quantum many-body systemsare expected to be quantum ergodic, implying that long-time averages of local observables after perturbations arethermal [1] and the system satisfies eigenstate thermal-ization hypothesis (ETH) [1–7]. Exceptions from thatgeneric behavior are currently under active investiga-tion. A well-established nonergodic behavior in transla-tionally invariant systems occurs at integrable points [8–16], for which eigenstates do not obey ETH [17–25]and the long-time averages of observables after quantumquenches are described by the generalized Gibbs ensem-ble [23, 26]. Recently, nonergodic properties of certaineigenstates (termed many-body scars) were observed insome translationally invariant models away from inte-grable points [27–34]. For disordered quantum many-body systems in one dimension that are the focus of thisLetter, it is proposed that (almost) all eigenstates be-come nonergodic at large enough disorder due to localiza-tion [35, 36]. This scenario predicts an eigenstate quan-tum phase transition from an ergodic to a nonergodicphase [36–39], the latter named many-body localization(MBL) [40–44].MBL and the corresponding transition were widelystudied by means of numerical approaches [35, 45–73]and phenomenological theories [74–85], as well as exper-imentally [86–92]. Moreover, under certain assumptionsin rigorous approaches, there are arguments about ex-istence of the MBL phase [93, 94]. The nature of thetransition, however, remains less clear. Numerically, thetransition was mainly studied within the framework ofpower-law divergence of correlation length [45, 47, 48, 53], ξ = 1 | W − W ∗ | ν , (1)where W ∗ is the critical disorder and ν is the criticalexponent. The main concern with most of the numericalresults in systems with uncorrelated disorder is that theysuggest ν ∼
1, which violates the Harris bound ν ≥ ν > ξ KT = exp (cid:40) b ± (cid:112) | W − W ∗ | (cid:41) , (2)where b − ( b + ) are nonuniversal parameters below (above)the transition. Hence, there is currently a gap betweenpredictions of exact numerical and phenomenological ap-proaches and as a consequence, the nature of the transi-tion remains an open problem.A new perspective in understanding of the ergodicitybreaking transition was recently obtained by calculatingthe spectral form factor [100], whose finite-size depen-dence was interpreted as a linear drift of the transitionpoint with system size, W ∗ ∝ L , for system sizes wherenumerical diagonalization of the full Hamiltonian matrixis accessible. This result raised the question whethersuch linear drift is an asymptotic feature suggesting thatthe transition to MBL is a crossover, or a preasymp-totic feature consistent with a phase transition at somevery large value of disorder in the thermodynamic limit.Subsequent work [101–103] mostly argued in favor of thesecond option. It is therefore an urgent need to intro-duce new unbiased numerical measures to characterizethe transition, which should also provide a benchmarkfor subsequent phenomenological studies.The goal of this Letter is to quantitatively comparedifferent scenarios of the ergodicity breaking transition.We introduce a cost function approach to describe thequality of the finite-size data collapse of ergodicity indi-cators as functions of L/ξ (i.e., the system size L dividedby the correlation length ξ ). This approach enables us toextract the most optimal form of the correlation lengthand to locate the disorder transition point in finite sys-tems. For the numerically accessible system sizes, ourresults consistently exhibit two main features: the cor-relation length ξ follows the Kosterlitz-Thouless behav-ior (2), and ergodicity indicators at the transition arevery close to random-matrix theory predictions. Model and Methodology.
We study interacting spin-1/2Hamiltonians with on-site disorder on a one-dimensional a r X i v : . [ c ond - m a t . d i s - nn ] A p r periodic lattice with L sites,ˆ H = (cid:88) j =1 J j L (cid:88) (cid:96) =1 (cid:16) ˆ s x(cid:96) ˆ s x(cid:96) + j + ˆ s y(cid:96) ˆ s y(cid:96) + j + ∆ j ˆ s z(cid:96) ˆ s z(cid:96) + j (cid:17) + L (cid:88) (cid:96) =1 w (cid:96) ˆ s z(cid:96) , (3)where ˆ s α(cid:96) ( α = x, y, z ) are spin-1/2 operators at site (cid:96) .We consider the total spin projection s z = 0 sector andset J ≡ W is generated by independent and identically dis-tributed local magnetic fields, with values w (cid:96) ∈ [ − W, W ]drawn from a uniform distribution.We study two disordered models, the J - J model(∆ = ∆ = 0 . J = 1), and the Heisenberg model(∆ = 1 and J = 0). For a given disorder distribution { w (cid:96) } , we calculate exact eigenstates around the centerof the spectrum using shift and invert diagonalizationmethod [104].We focus on two widely studied ergodicity indicatorsthat characterize properties of Hamiltonian eigenstatesand eigenvalues: the eigenstate entanglement entropy S and the spectral level spacing ratio r , respectively. In thecontext of Anderson localization, statistics of Hamilto-nian eigenvalues and the corresponding scaling solutionsas functions of L/ξ represented one of the main numer-ical approaches to detect the transition point [105, 106].Recently, ergodicity indicators S and r have been ex-tensively studied in the context of disordered interactingspin chains [35, 37, 45–57, 107–109].First, we calculate the von Neumann entanglemententropy S α = − Tr { ˆ ρ A ln(ˆ ρ A ) } in an eigenstate | α (cid:105) ,where the subsystem A consists of the first L/ ρ A = Tr L − A { ˆ ρ } is the trace over the re-maining L/ ρ = | α (cid:105)(cid:104) α | . Since we studythe total s z = 0 sector, we divide S α by the corre-sponding random-matrix theory (RMT) result S RMT =( L/
2) ln(2) + (1 / / / − /
2, which includes O (1) contributions and hence minimizes finite-size ef-fects [110]. To calculate the level spacing ratio r , wefirst calculate r α = min { δ α , δ α − } / max { δ α , δ α − } for aneigenstate | α (cid:105) , where δ α = E α +1 − E α is the energylevel spacing [35]. We then obtain S and r by averaging S α /S RMT and r α , respectively, over eigenstates aroundthe center of the spectrum and over different realiza-tions of the disorder distribution { w (cid:96) } [111]. Results for S ( W ) and r ( W ) for different system sizes L are shownfor both models in the insets of Fig. 1 and 2, respec-tively. In the ergodic (small W ) regime, S ≈ r ≈ r GOE ≈ . W in finite systems, S → r → r Poisson = 2 ln(2) − ≈ . S ( W, L ) and r ( W, L ) as functions of
L/ξ . Specifically,we want to establish an unbiased, quantitative measureof the quality of the data collapse for different functionalforms of the correlation length ξ [such as ξ and ξ KT from Eqs. (1)-(2)], and the critical disorder W ∗ , whichis included in the expression for ξ . This is achieved byintroducing the cost function for a quantity X ∈ { S, r } − L/ξ KT . . . S (a) J - J model L = 12 L = 14 L = 16 L = 18 L = 20 W . . . S − L/ξ KT . . . S (b) Heisenberg model L = 12 L = 14 L = 16 L = 18 L = 20 W . . . S FIG. 1. Eigenstate entanglement entropy S for differentsystems sizes L , calculated (a) in the J - J model and (b)in the Heisenberg model. Insets show S as a function of dis-order W . In the main panels, we plot S as a function of L/ξ KT [ − L/ξ KT if W < W ∗ ], where ξ KT is a KT correla-tion length (2), assuming b − = b + ≡ b and the transitionpoint ansatz W ∗ = w + w L . The optimal parameters b , w and w in ξ KT are obtained by minimizing the cost func-tion C S ( ξ KT ) in Eq. (4). The number of data points includedin the minimization procedure is N p = 287 in panel (a) and N p = 225 in panel (b). See Fig. 3 and [111] for details. that consists of N p values at different W and L , C X = (cid:80) N p − j =1 | X j +1 − X j | max { X j } − min { X j } − . (4)In Eq. (4), we sort all N p values of X j according tonondecreasing values of sign[ W - W ∗ ] L/ξ . In the case ofan ideal data collapse, this implies (cid:80) j | X j +1 − X j | =max { X j } − min { X j } and therefore C X = 0. For thelarge data sets studied here, the cost function is alwayspositive, C X >
0. Our goal is to find, for given functionalforms of the correlation length ξ ( W, L ) and the criticaldisorder W ∗ ( L ), the optimal values of fitting parame-ters that minimize C X . We apply the cost function mini-mization algorithm to results for S ( W, L ) and r ( W, L ) at L = 12 , , , ,
20 in the J - J model and the Heisen-berg model (see Supplemental Material [111] for details). Results.
We now describe our main results, valid forboth ergodicity indicators S and r , and for both investi-gated models. For simplicity, we consider a single fitting − L/ξ KT . . . r (a) J - J model L = 12 L = 14 L = 16 L = 18 L = 20 W . . . r − L/ξ KT . . . r (b) Heisenberg model L = 12 L = 14 L = 16 L = 18 L = 20 W . . . r FIG. 2. Level spacing ratio r for different systems sizes L ,calculated (a) in the J - J model and (b) in the Heisenbergmodel. Insets show S as a function of disorder W . In themain panels, we plot r as a function of L/ξ KT [ − L/ξ KT if W < W ∗ ], where ξ KT is a KT correlation length (2), assuming b − = b + ≡ b and the transition point ansatz W ∗ = w + w L .The optimal parameters b , w and w in ξ KT are obtainedby minimizing the cost function C S ( ξ KT ) in Eq. (4). Thenumber of data points included in the minimization procedureis N p = 285 in panel (a) and N p = 175 in panel (b). See Fig. 3and [111] for details. parameter b − = b + ≡ b in ξ KT in Eq. (2). The scenariowith b − (cid:54) = b + and the optimal values of b are discussedin [111].(i) For the simplest functional forms of the transitionpoint we consider two fitting functions W ∗ = w and W ∗ = w + w L , with free parameters w and w . If weonly consider an L -independent function W ∗ = w , thedata collapse (quantified in terms of the cost function)is better as a function of L/ξ than L/ξ KT [left columnsin Tables I and II]. However, the data collapse becomesmuch better (i.e., the cost function becomes much lower)for the KT transition if W ∗ is allowed to increase linearlywith L [central columns in Tables I and II]. The latterstatement holds true also in the special case of zero offset, W ∗ = w L .As an example, we show in Fig. 1 the entanglemententropy S for both models using W ∗ = w + w L , andplot results in the main panels as functions of L/ξ KT . InFig. 2, analogous results are shown for the level spacing
10 12 14 16 18 20 L W ∗ J - J modelHeisenberg model r for W ∗ = w + w LS for W ∗ = w + w Lr for W ∗ = w ∗ ( L ) S for W ∗ = w ∗ ( L ) FIG. 3. Ergodicity breaking transition point W ∗ as a func-tion of system size L , for the J - J model (upper part) andthe Heisenberg model (lower part). Results for W ∗ are ob-tained from the best data collapse using a KT correlationlength ξ KT from Eq. (2), with b as a free parameter (assum-ing b − = b + ). Lines are results for a transition point ansatz W ∗ = w + w L (with free parameters w and w ), sym-bols are results for a function-independent transition point W ∗ = w ∗ ( L ) [with five free parameters for five different sys-tems sizes L = 12 , , , , S (circles and solid lines) and r (squares and dashed lines).In the J - J model we get ( w , w ) = (0 . , .
23) for S and (0 . , .
25) for r , while in the Heisenberg model we get( w , w ) = (0 . , .
10) for S and (0 . , .
09) for r . Therefore,both ergodicity indicators S and r have very similar leadingterm w within the same model, while the subleading term w may be different. ratio r . In both figures, the scaling collapses appear tobe excellent.Results at L/ξ KT ≈ S at the critical pointwere previously discussed in Refs. [52, 82]. The secondobservation is that r ( L/ξ KT ≈
0) and S ( L/ξ KT ≈
0) arevery close to the RMT predictions. Hence, their valuesare close to the thermal values.(ii) The results from (i) are robust towards more gen-eral functional forms of the transition point W ∗ . First,we tested the fitting functions W ∗ = w + w /L and W ∗ = w + w / ln( L ), which are two parameter fits to W ∗ , but always yielded substantially higher cost func-tion when compared to the function W ∗ = w + w L (see [111] for details). Then, we tested a general ansatz W ∗ = w ∗ ( L ) using independent fitting values of w ∗ forevery L [i.e., five different w ∗ for five different systemssizes L = 12 , , , , w ∗ ( L ) afterminimization are rather accurately described by the func-tion W ∗ = w + w L for both ergodicity indicators S and r , as shown in Fig. 3.We compare the values of the cost functions for thedisordered J - J model and the disordered Heisenbergmodel in Tables I and II, respectively. In Figs. S1-S4of [111] we also visually compare the scaling collapses of S and r with and without the drift in the functional de-pendence of the transition point W ∗ . Note that the costfunctions C X [ ξ KT ] for the general model W ∗ = w ∗ ( L )take only slighly lower values than for the linear one W ∗ = w + w L .We interpret our results as evidence that a transi-tion with a KT correlation length ξ KT is more favorablethan a transition with a power-law correlation length ξ .A particularly suggestive evidence supporting the latterstatement is that cost functions using ξ KT with a two-parameter transition point function W ∗ = w + w L [cen-tral columns in Tables I and II] are substantially lowerthan cost functions using ξ with a five-parameter func-tion W ∗ = w ∗ ( L ) [right columns in Tables I and II].Moreover, for the numerically available system sizesthe best scaling collapse of the transition using a KTcorrelation length ξ KT consistently exhibits a linear driftwith system size, W ∗ ∝ L , as suggested by Fig. 3. How-ever, an analogous statement cannot be made for thetransition using ξ . In the latter case, using W ∗ = w ∗ ( L )[right column in Tables I and II], we find degenerate solu-tions (i.e., almost identical values of cost functions) withvery different functional forms of W ∗ .The observed linear scaling W ∗ ∝ w L opens a ques-tion about its fate in the thermodynamic limit. While thesymbols in Fig. 3 show no tendency towards approachinga horizontal line, we can also not exclude scenarios where W ∗ ( L ) ∝ w L represents a small size behavior that even-tually saturates to a finite critical point W ∗∞ in the ther-modynamic limit. As a quantitative estimate of the lat-ter scenario we consider the functional form W ∗ ( L ) = W ∗∞ tanh ( L/L ), where L represents a characteristiclength scale. This form reproduces the linear L depen-dence at L (cid:28) L since W ∗ ( L ) ≈ W ∗∞ L/L + O ( L /L ),where w = W ∗∞ /L . Then, by requiring the leadingterm to be much larger than the subleading term at L = 20, one could estimate a lower bound for W ∗∞ . Us-ing the cost function minimization approach, we estimate L (cid:38)
50 (see [111] for details), and hence W ∗∞ (cid:38) w = 0 .
10) and W ∗∞ (cid:38)
12 in the J - J model (using w = 0 . Discussion.
Intriguingly, our analysis predicts nearlythermal (RMT-like) properties of the transition point infinite systems. This is a consequence of the transitiontaking place at relatively weak disorder: for L ≈
20, itoccurs at W ∗ ≈ W ∗ ≈ J - J model. These values are lower than thoseusually considered in the MBL literature, which mostlyfollowed the initial proposal [35] for the transition pointhaving insulator-like (Poisson) statistics. Indeed, our sce-nario for the transition differs from previous numericalstudies of the same ergodicity indicators in the Heisen-berg model [47, 48], which explored the optimal data col-lapse as a function of L/ξ and obtained W ∗ ≈ . − . TABLE I. Cost function C X , see Eq. (4), in the J - J model.Values of C X are shown for two ergodicity indicators X ∈{ S, r } , using correlation lengths ξ and ξ KT from Eqs. (1)-(2). Columns denote different functional forms of W ∗ usedin ξ and ξ KT . Results are obtained using the data points inFigs. 1(a) and 2(a) for W > . W ∗ = w W ∗ = w + w L W ∗ = w ∗ ( L ) C S [ ξ KT ] 3.99 0.34 0.29 C S [ ξ ] 2.80 1.81 1.71 C r [ ξ KT ] 5.46 1.01 0.92 C r [ ξ ] 4.31 2.71 2.57TABLE II. Cost function C X , see Eq. (4), in the Heisenbergmodel. Values of C X are shown for two ergodicity indica-tors X ∈ { S, r } , using correlation lengths ξ and ξ KT fromEqs. (1)-(2). Columns denote different functional forms of W ∗ used in ξ and ξ KT . Results are obtained using the datapoints in Figs. 1(b) and 2(b) for W > . W ∗ = w W ∗ = w + w L W ∗ = w ∗ ( L ) C S [ ξ KT ] 2.51 0.46 0.29 C S [ ξ ] 1.80 0.94 0.77 C r [ ξ KT ] 2.84 0.51 0.46 C r [ ξ ] 2.16 1.08 1.01 for comparable system sizes.The outcome of our scaling analysis can be consideredas a first step in unifying exact numerical calculationswith RG approaches that also predict a KT-like transi-tion [98, 99]. One of the next goals is to better under-stand the character of the transition point: while our re-sults suggest nearly thermal properties, the RG schemespredict a vanishing density of thermal blocks at the tran-sition [98, 99].Our analysis applies to disorder averages of ergodicityindicators. When the disorder is increased, fluctuationsof ergodicity indicators (at least in finite systems) maybecome anomalous, which has been observed in fluctua-tions of the entanglement entropy [45, 50–53] and in dis-tributions of other observables [66, 113–116]. It remainsopen how the KT character of the transition of averagedquantities is related to other statistical properties of themodel. Conclusions.
In this Letter we proposed an unbiased,quantitative approach based on cost function minimiza-tion to test the nature of a disorder driven transitionin finite quantum spin chains. We argued that certainkey ergodicity indicators exhibit clear signatures of theergodicity breakdown and remarkable finite-size scalingproperties, which do not violate the Harris bound [95–97]and are consistent with the KT quantum phase transi-tion.Moreover, the cost function minimization approach fornumerically accessible system sizes results in two fea-tures: a nearly thermal character of the transition pointand its linear drift with the system size. In fact, bothproperties emerge simultaneously with the KT characterof the transition. The main open question is how closeare the numerical results to the true asymptotic regime,and what is the relation between the observed propertiesof the transition point and the KT nature of the transi-tion when the thermodynamic limit is approached. More work is needed to clarify this.
ACKNOWLEDGMENTS
We acknowledge discussions with F. Heidrich-Meisner.This work is supported by the Slovenian Research Agency(ARRS), Research core fundings No. P1-0044 (J.ˇS., J.B.and L.V.), No. P1-0402 (T.P.) and No. J1-1696 (L.V.),and by the European Research Council (ERC) under Ad-vanced Grant 694544 – OMNES (T.P.). [1] L. D’Alessio, Y. Kafri, A. Polkovnikov, and M. Rigol,From quantum chaos and eigenstate thermalization tostatistical mechanics and thermodynamics, Adv. Phys. , 239 (2016).[2] J. M. Deutsch, Quantum statistical mechanics in aclosed system, Phys. Rev. A , 2046 (1991).[3] M. Srednicki, Chaos and quantum thermalization, Phys.Rev. E , 888 (1994).[4] M. Rigol, V. Dunjko, and M. Olshanii, Thermalizationand its mechanism for generic isolated quantum sys-tems, Nature (London) , 854 (2008).[5] J. Eisert, M. Friesdorf, and C. Gogolin, Quantum many-body systems out of equilibrium, Nat. Phys. , 124(2015).[6] J. M. Deutsch, Eigenstate thermalization hypothesis,Rep. Prog. Phys. , 082001 (2018).[7] T. Mori, T. N. Ikeda, E. Kaminishi, and M. Ueda, Ther-malization and prethermalization in isolated quantumsystems: a theoretical overview, J. Phys. B , 112001(2018).[8] T. Kinoshita, T. Wenger, and S. D. Weiss, A quantumNewton’s cradle, Nature (London) , 900 (2006).[9] Y. Tang, W. Kao, K.-Y. Li, S. Seo, K. Mallayya,M. Rigol, S. Gopalakrishnan, and B. L. Lev, Thermal-ization near integrability in a dipolar quantum Newton’scradle, Phys. Rev. X , 021030 (2018).[10] J.-S. Caux and F. H. L. Essler, Time evolution of lo-cal observables after quenching to an integrable model,Phys. Rev. Lett. , 257203 (2013).[11] B. Wouters, J. De Nardis, M. Brockmann, D. Fioretto,M. Rigol, and J.-S. Caux, Quenching the anisotropicHeisenberg chain: Exact solution and generalized Gibbsensemble predictions, Phys. Rev. Lett. , 117202(2014).[12] B. Pozsgay, M. Mesty´an, M. A. Werner, M. Kormos,G. Zar´and, and G. Tak´acs, Correlations after quantumquenches in the XXZ spin chain: Failure of the gener-alized Gibbs ensemble, Phys. Rev. Lett. , 117203(2014).[13] E. Ilievski, J. De Nardis, B. Wouters, J.-S. Caux,F. H. L. Essler, and T. Prosen, Complete GeneralizedGibbs Ensembles in an Interacting Theory, Phys. Rev.Lett. , 157201 (2015).[14] P. Calabrese, F. H. L. Essler, and G. Mussardo, Intro-duction to ‘quantum integrability in out of equilibriumsystems’, J. Stat. Mech. (2016), 064001. [15] F. H. L. Essler and M. Fagotti, Quench dynamics andrelaxation in isolated integrable quantum spin chains,J. Stat. Mech. (2016), 064002.[16] M. A. Cazalilla and M.-C. Chung, Quantum quenchesin the Luttinger model and its close relatives, J. Stat.Mech. (2016), 064004.[17] M. Rigol, Breakdown of thermalization in finite one-dimensional systems, Phys. Rev. Lett. , 100403(2009).[18] A. C. Cassidy, C. W. Clark, and M. Rigol, Generalizedthermalization in an integrable lattice system, Phys.Rev. Lett. , 140405 (2011).[19] R. Steinigeweg, J. Herbrych, and P. Prelovˇsek, Eigen-state thermalization within isolated spin-chain systems,Phys. Rev. E , 012118 (2013).[20] T. N. Ikeda, Y. Watanabe, and M. Ueda, Finite-sizescaling analysis of the eigenstate thermalization hypoth-esis in a one-dimensional interacting Bose gas, Phys.Rev. E , 012125 (2013).[21] V. Alba, Eigenstate thermalization hypothesis and in-tegrability in quantum spin chains, Phys. Rev. B ,155123 (2015).[22] M. Rigol, Fundamental asymmetry in quenches betweenintegrable and nonintegrable systems, Phys. Rev. Lett. , 100601 (2016).[23] L. Vidmar and M. Rigol, Generalized Gibbs ensemble inintegrable lattice models, J. Stat. Mech. (2016), 064007.[24] T. LeBlond, K. Mallayya, L. Vidmar, and M. Rigol,Entanglement and matrix elements of observables in in-teracting integrable systems, Phys. Rev. E , 062134(2019).[25] M. Mierzejewski and L. Vidmar, Quantitative impactof integrals of motion on the eigenstate thermalizationhypothesis, Phys. Rev. Lett. , 040603 (2020).[26] M. Rigol, V. Dunjko, V. Yurovsky, and M. Olshanii,Relaxation in a completely integrable many-body quan-tum system: An Ab Initio study of the dynamics ofthe highly excited states of 1D lattice hard-core bosons,Phys. Rev. Lett. , 050405 (2007).[27] N. Shiraishi and T. Mori, Systematic construction ofcounterexamples to the eigenstate thermalization hy-pothesis, Phys. Rev. Lett. , 030601 (2017).[28] R. Mondaini, K. Mallayya, L. F. Santos, and M. Rigol,Comment on “Systematic construction of counterexam-ples to the eigenstate thermalization hypothesis”, Phys.Rev. Lett. , 038901 (2018). [29] S. Moudgalya, N. Regnault, and B. A. Bernevig, En-tanglement of exact excited states of Affleck-Kennedy-Lieb-Tasaki models: Exact results, many-body scars,and violation of the strong eigenstate thermalization hy-pothesis, Phys. Rev. B , 235156 (2018).[30] C. J. Turner, A. A. Michailidis, D. A. Abanin, M. Ser-byn, and Z. Papi´c, Quantum scarred eigenstates in aRydberg atom chain: Entanglement, breakdown of ther-malization, and stability to perturbations, Phys. Rev. B , 155134 (2018).[31] C. J. Turner, A. A. Michailidis, D. A. Abanin, M. Ser-byn, and Z. Papi´c, Weak ergodicity breaking from quan-tum many-body scars, Nat. Phys. , 745 (2018).[32] C.-J. Lin and O. I. Motrunich, Exact quantum many-body scar states in the Rydberg-blockaded atom chain,Phys. Rev. Lett. , 173401 (2019).[33] T. Iadecola and M. ˇZnidariˇc, Exact localized and ballis-tic eigenstates in disordered chaotic spin ladders and theFermi-Hubbard model, Phys. Rev. Lett. , 036403(2019).[34] M. Schecter and T. Iadecola, Weak ergodicity breakingand quantum many-body scars in spin-1 XY magnets,Phys. Rev. Lett. , 147201 (2019).[35] V. Oganesyan and D. A. Huse, Localization of inter-acting fermions at high temperature, Phys. Rev. B ,155111 (2007).[36] A. Pal and D. A. Huse, Many-body localization phasetransition, Phys. Rev. B , 174411 (2010).[37] B. Bauer and C. Nayak, Area laws in a many-body lo-calized state and its implications for topological order,J. Stat. Mech. , P09005 (2013).[38] D. A. Huse, R. Nandkishore, V. Oganesyan, A. Pal,and S. L. Sondhi, Localization-protected quantum or-der, Phys. Rev. B , 014206 (2013).[39] D. Pekker, G. Refael, E. Altman, E. Demler, andV. Oganesyan, Hilbert-glass transition: New universal-ity of temperature-tuned many-body dynamical quan-tum criticality, Phys. Rev. X , 011052 (2014).[40] D. Basko, I. Aleiner, and B. Altshuler, Metal-insulatortransition in a weakly interacting many-electron systemwith localized single-particle states, Ann. Phys. ,1126 (2006).[41] I. V. Gornyi, A. D. Mirlin, and D. G. Polyakov, Inter-acting electrons in disordered wires: Anderson localiza-tion and low-T transport, Phys. Rev. Lett. , 206603(2005).[42] R. Nandkishore and D. A. Huse, Many-body-localization and thermalization in quantum statisticalmechanics, Ann. Rev. Cond. Mat. Phys. , 15 (2015).[43] E. Altman and R. Vosk, Universal dynamics and renor-malization in many-body-localized systems, Ann. Rev.Cond. Mat. Phys. , 383 (2015).[44] D. A. Abanin, E. Altman, I. Bloch, and M. Serbyn, Col-loquium: Many-body localization, thermalization, andentanglement, Rev. Mod. Phys. , 021001 (2019).[45] J. A. Kj¨all, J. H. Bardarson, and F. Pollmann, Many-body localization in a disordered quantum Ising chain,Phys. Rev. Lett. , 107204 (2014).[46] T. Devakul and R. R. P. Singh, Early breakdown ofarea-law entanglement at the many-body delocalizationtransition, Phys. Rev. Lett. , 187201 (2015).[47] D. J. Luitz, N. Laflorencie, and F. Alet, Many-bodylocalization edge in the random-field Heisenberg chain,Phys. Rev. B , 081103 (2015). [48] C. L. Bertrand and A. M. Garc´ıa-Garc´ıa, AnomalousThouless energy and critical statistics on the metal-lic side of the many-body localization transition, Phys.Rev. B , 144201 (2016).[49] M. Serbyn and J. E. Moore, Spectral statistics acrossthe many-body localization transition, Phys. Rev. B ,041424 (2016).[50] D. J. Luitz, Long tail distributions near the many-bodylocalization transition, Phys. Rev. B , 134201 (2016).[51] X. Yu, D. J. Luitz, and B. K. Clark, Bimodal entangle-ment entropy distribution in the many-body localizationtransition, Phys. Rev. B , 184202 (2016).[52] V. Khemani, S. P. Lim, D. N. Sheng, and D. A. Huse,Critical properties of the many-body localization tran-sition, Phys. Rev. X , 021013 (2017).[53] V. Khemani, D. N. Sheng, and D. A. Huse, Two univer-sality classes for the many-body localization transition,Phys. Rev. Lett. , 075702 (2017).[54] S.-H. Lin, B. Sbierski, F. Dorfner, C. Karrasch, andF. Heidrich-Meisner, Many-body localization of spinlessfermions with attractive interactions in one dimension,SciPost Phys. , 002 (2018).[55] W. Buijsman, V. Cheianov, and V. Gritsev, Randommatrix ensemble for the level statistics of many-bodylocalization, Phys. Rev. Lett. , 180601 (2019).[56] P. Sierant and J. Zakrzewski, Level statistics across themany-body localization transition, Phys. Rev. B ,104205 (2019).[57] P. Sierant and J. Zakrzewski, Model of level statisticsfor disordered interacting quantum many-body systems,Phys. Rev. B , 104201 (2020).[58] T. C. Berkelbach and D. R. Reichman, Conductivity ofdisordered quantum lattice models at infinite tempera-ture: Many-body localization, Phys. Rev. B , 224429(2010).[59] O. S. Bariˇsi´c and P. Prelovˇsek, Conductivity in a dis-ordered one-dimensional system of interacting fermions,Phys. Rev. B , 161106 (2010).[60] E. Khatami, M. Rigol, A. Rela˜no, and A. M. Garc´ıa-Garc´ıa, Quantum quenches in disordered systems: Ap-proach to thermal equilibrium without a typical relax-ation time, Phys. Rev. E , 050102 (2012).[61] R. Modak and S. Mukerjee, Many-body localization inthe presence of a single-particle mobility edge, Phys.Rev. Lett. , 230401 (2015).[62] S. Bera, H. Schomerus, F. Heidrich-Meisner, and J. H.Bardarson, Many-body localization characterized froma one-particle perspective, Phys. Rev. Lett. , 046603(2015).[63] R. Mondaini and M. Rigol, Many-body localizationand thermalization in disordered Hubbard chains, Phys.Rev. A , 041601 (2015).[64] E. J. Torres-Herrera and L. F. Santos, Dynamics atthe many-body localization transition, Phys. Rev. B ,014208 (2015).[65] M. Serbyn, Z. Papi´c, and D. A. Abanin, Criterion formany-body localization-delocalization phase transition,Phys. Rev. X , 041047 (2015).[66] M. Serbyn, Z. Papi´c, and D. A. Abanin, Thouless energyand multifractality across the many-body localizationtransition, Phys. Rev. B , 104201 (2017).[67] O. S. Bariˇsi´c, J. Kokalj, I. Balog, and P. Prelovˇsek,Dynamical conductivity and its fluctuations along thecrossover to many-body localization, Phys. Rev. B , , 016804(2017).[69] E. V. H. Doggen, F. Schindler, K. S. Tikhonov, A. D.Mirlin, T. Neupert, D. G. Polyakov, and I. V. Gornyi,Many-body localization and delocalization in largequantum chains, Phys. Rev. B , 174202 (2018).[70] N. Mac´e, F. Alet, and N. Laflorencie, Multifractal scal-ings across the many-body localization transition, Phys.Rev. Lett. , 180601 (2019).[71] Z. Lenarˇciˇc, O. Alberton, A. Rosch, and E. Altman,Critical behavior near the many-body localization tran-sition in driven open systems, arXiv:1910.01548.[72] G. Guidici, F. M. Surace, J. E. Ebot, A. Scardicchio,and M. Dalmonte, Breakdown of ergodicity in disor-dered U(1) lattice gauge theories, arXiv:1912.09403.[73] A. L. Corps, R. A. Molina, and A. Rela˜no, Thermal-ization and Thouless energy in the ergodic side of themany-body localization transition, arXiv:2001.03990.[74] M. Serbyn, Z. Papi´c, and D. A. Abanin, Local conserva-tion laws and the structure of the many-body localizedstates, Phys. Rev. Lett. , 127201 (2013).[75] D. A. Huse, R. Nandkishore, and V. Oganesyan,Phenomenology of fully many-body-localized systems,Phys. Rev. B , 174202 (2014).[76] V. Ros, M. M¨uller, and A. Scardicchio, Integrals of mo-tion in the many-body localized phase, Nuc. Phys. B , 420 (2015).[77] A. Chandran, A. Pal, C. R. Laumann, and A. Scardic-chio, Many-body localization beyond eigenstates in alldimensions, Phys. Rev. B , 144203 (2016).[78] D. Pekker, B. K. Clark, V. Oganesyan, and G. Refael,Fixed points of Wegner-Wilson flows and many-bodylocalization, Phys. Rev. Lett. , 075701 (2017).[79] R. Vosk, D. A. Huse, and E. Altman, Theory of themany-body localization transition in one-dimensionalsystems, Phys. Rev. X , 031032 (2015).[80] A. C. Potter, R. Vasseur, and S. A. Parameswaran,Universal properties of many-body delocalization tran-sitions, Phys. Rev. X , 031033 (2015).[81] L. Zhang, B. Zhao, T. Devakul, and D. A. Huse, Many-body localization phase transition: A simplified strong-randomness approximate renormalization group, Phys.Rev. B , 224201 (2016).[82] P. T. Dumitrescu, R. Vasseur, and A. C. Potter, Scalingtheory of entanglement at the many-body localizationtransition, Phys. Rev. Lett. , 110604 (2017).[83] W. De Roeck and F. Huveneers, Stability and insta-bility towards delocalization in many-body localizationsystems, Phys. Rev. B , 155129 (2017).[84] T. Thiery, F. Huveneers, M. M¨uller, and W. De Roeck,Many-body delocalization as a quantum avalanche,Phys. Rev. Lett. , 140601 (2018).[85] S. Balasubramanian, Y. Liao, and V. Galitski, Many-body localization landscape, Phys. Rev. B , 014201(2020).[86] M. Schreiber, S. S. Hodgman, P. Bordia, H. P. L¨uschen,M. H. Fischer, R. Vosk, E. Altman, U. Schneider, andI. Bloch, Observation of many-body localization of in-teracting fermions in a quasi-random optical lattice, Sci-ence , 842 (2015). [87] S. S. Kondov, W. R. McGehee, W. Xu, and B. DeMarco,Disorder-induced localization in a strongly correlatedatomic Hubbard gas, Phys. Rev. Lett. , 083002(2015).[88] J. Smith, A. Lee, P. Richerme, B. Neyenhuis, P. W.Hess, P. Hauke, M. Heyl, D. A. Huse, and C. Mon-roe, Many-body localization in a quantum simulatorwith programmable random disorder, Nat. Phys. ,907 (2016).[89] J.-y. Choi, S. Hild, J. Zeiher, P. Schauß, A. Rubio-Abadal, T. Yefsah, V. Khemani, D. A. Huse, I. Bloch,and C. Gross, Exploring the many-body localizationtransition in two dimensions, Science , 1547 (2016).[90] A. Lukin, M. Rispoli, R. Schittko, M. E. Tai,A. M. Kaufman, S. Choi, V. Khemani, J. L´eonard,and M. Greiner, Probing entanglement in a many-body–localized system, Science , 256 (2019).[91] M. Rispoli, A. Lukin, R. Schittko, S. Kim, M. E. Tai,J. Lonard, and M. Greiner, Quantum critical behaviourat the many-body localization transition, Nature ,385 (2019).[92] B. Chiaro et al , Growth and preservation ofentanglement in a many-body localized system,arxiv:1910.06024.[93] J. Z. Imbrie, On many-body localization for quantumspin chains, J. Stat. Phys. , 998 (2016).[94] J. Z. Imbrie, Diagonalization and many-body localiza-tion for a disordered quantum spin chain, Phys. Rev.Lett. , 027201 (2016).[95] A. B. Harris, Effect of random defects on the criticalbehaviour of ising models, J. Phys. C , 1671 (1974).[96] J. T. Chayes, L. Chayes, D. S. Fisher, and T. Spencer,Finite-size scaling and correlation lengths for disorderedsystems, Phys. Rev. Lett. , 2999 (1986).[97] A. Chandran, C. R. Laumann, and V. Oganesyan, Fi-nite size scaling bounds on many-body localized phasetransitions, arxiv:1509.04285 (2015).[98] P. T. Dumitrescu, A. Goremykina, S. A. Parameswaran,M. Serbyn, and R. Vasseur, Kosterlitz-Thouless scalingat many-body localization phase transitions, Phys. Rev.B , 094205 (2019).[99] A. Goremykina, R. Vasseur, and M. Serbyn, Analyti-cally solvable renormalization group for the many-bodylocalization transition, Phys. Rev. Lett. , 040601(2019).[100] J. ˇSuntajs, J. Bonˇca, T. Prosen, and L. Vidmar,Quantum chaos challenges many-body localization,arXiv:1905.06345.[101] D. A. Abanin, J. H. Bardarson, G. De Tomasi,S. Gopalakrishnan, V. Khemani, S. A. Parameswaran,F. Pollmann, A. C. Potter, M. Serbyn, and R. Vasseur,Distinguishing localization from chaos: challenges infinite-size systems, arXiv:1911.04501.[102] P. Sierant, D. Delande, and J. Zakrzewski, Thoulesstime analysis of anderson and many-body localizationtransitions, arXiv:1911.06221.[103] R. K. Panda, A. Scardicchio, M. Schulz, S. R. Taylor,and M. ˇZnidariˇc, Can we study the many-body locali-sation transition?, EPL , 67003 (2020).[104] F. Pietracaprina, N. Mac´e, D. J. Luitz, and F. Alet,Shift-invert diagonalization of large many-body localiz-ing spin chains, SciPost Phys. , 045 (2018).[105] B. I. Shklovskii, B. Shapiro, B. R. Sears, P. Lambrian-ides, and H. B. Shore, Statistics of spectra of disordered systems near the metal-insulator transition, Phys. Rev.B , 11487 (1993).[106] I. K. Zharekeshev and B. Kramer, Asymptotics of uni-versal probability of neighboring level spacings at theAnderson transition, Phys. Rev. Lett. , 717 (1997).[107] T. Grover, Certain general constraints on the many-body localization transition, arxiv:1405.1471.[108] J. Richter, D. Schubert, and R. Steinigeweg, Decay ofspin-spin correlations in disordered quantum and classi-cal spin chains, Phys. Rev. Research , 013130 (2020).[109] Y. Bar Lev, G. Cohen, and D. R. Reichman, Absence ofdiffusion in an interacting system of spinless fermions ona one-dimensional disordered lattice, Phys. Rev. Lett. , 100601 (2015).[110] L. Vidmar and M. Rigol, Entanglement entropy ofeigenstates of quantum chaotic Hamiltonians, Phys.Rev. Lett. , 220603 (2017).[111] See Supplemental Material for details about the costfunction minimization procedure, qualitative compari-son of different scenarios for the transition point, lowerbound estimate for the critical point, and the com-parison between different forms of the KT correlationlength.[112] Y. Y. Atas, E. Bogomolny, O. Giraud, and G. Roux,Distribution of the ratio of consecutive level spacings inrandom matrix ensembles, Phys. Rev. Lett. , 084101(2013).[113] S. Sahu, S. Xu, and B. Swingle, Scrambling dynam-ics across a thermalization-localization quantum phasetransition, Phys. Rev. Lett. , 165902 (2019).[114] L. Colmenarez, P. A. McClarty, M. Haque, and D. J.Luitz, Statistics of correlations functions in the randomHeisenberg chain, SciPost Phys. , 64 (2019).[115] M. Schulz, S. R. Taylor, A. Scardicchio, andM. ˇZnidariˇc, Phenomenology of anomalous transportin disordered one-dimensional systems, J. Stat. Mech. , 023107 (2020).[116] M. Mierzejewski, M. Sroda, J. Herbrych, andP. Prelovˇsek, Resistivity and its fluctuations in dis-ordered many-body systems: from chains to planes,arxiv:2003.00495 (2020). Supplemental Material:Ergodicity Breaking Transition in Finite Disordered Spin Chains
Jan ˇSuntajs , Janez Bonˇca , , Tomaˇz Prosen and Lev Vidmar , Department of Theoretical Physics, J. Stefan Institute, SI-1000 Ljubljana, Slovenia Department of Physics, Faculty of Mathematics and Physics, University of Ljubljana, SI-1000 Ljubljana, Slovenia
S1. DETAILS ABOUT THE COST FUNCTIONMINIMIZATION PROCEDURE
Here we provide more details about the minimizationprocedure of the cost function C X introduced in Eq. (4)of the main text.We first implement numerical exact diagonalization tocalculate eigenstates and eigenvalues of the Hamiltonianunder investigation. For each Hamiltonian realizationˆ H µ (i.e., for fixed model parameters and for randomlygenerated disorder potentials), we target N eig eigenstatesclosest to the mean energy ¯ E µ = Tr { ˆ H µ } / D , where D =is the Hilbert space dimension in the s z = 0 sector. Weset N eig = 100 in Fig. 1 and N eig = 500 in Fig. 2. Resultsare further averaged over N sample different realizations ofdisorder. For the entanglement entropy S in Fig. 1(a),we use N sample ≥ L ≤
18 and N sample ≥ L = 20, while in Fig. 1(b), we use N sample ≥ L ≤
18 and N sample ≥
100 for L = 20. For the levelspacing ratio r in Fig. 2(a), we use N sample ≥
450 for L ≤
18 and N sample ≥
350 for L = 20, while in Fig. 2(b),we use N sample ≥ L ≤
18 and N sample ≥
100 for L = 20.For a given functional form of the correlation length ξ and the critical disorder W ∗ [which is included in thefunctional form of ξ , see Eqs. (1) and (2)], we thensort numerical values of S and r at different W and L ∈ { , , , , } according to nondecreasing L/ξ .In Fig. 1(a) the data included in the minimization pro-cedure were for 0 . ≤ W ≤
15 (except for for L = 20,where 0 . ≤ W ≤
13) and ∆ W = 0 .
25. In Fig. 1(b) thedata included were for 0 . ≤ W ≤ . . ≤ W ≤
15 (except for L = 20, where0 . ≤ W ≤ .
5) and ∆ W = 0 .
25. In Fig. 2(b) the dataincluded were for 0 . ≤ W ≤ .
75 and ∆ = 0 . scipy to find the optimal set of free parameters of ξ and W ∗ that minimize the cost function. In each realiza-tion, we employ a population size 10 and allow for up to10 iterations. We use the relative tolerance of conver-gence 10 − and employ 10 realizations of the algorithmto verify the precision of fitted parameters of the optimalsolution.In the main text (see Tables I and II) we compared costfunctions with the functional form of the critical disor-der W ∗ = w [one free parameter], W ∗ = w + w L [twofree parameters] and W ∗ = w ∗ ( L ) [five free parametersfor five different L ]. Here we complement those resultsby focusing on two-parameter functional forms W ∗ = TABLE S1. Cost function C X , see Eq. (4), in the J - J model. Values of C X are shown for two ergodicity indica-tors X ∈ { S, r } , using correlation lengths ξ and ξ KT fromEqs. (1)-(2). Columns denote different functional forms of W ∗ used in ξ and ξ KT . Results are obtained using the datapoints in Figs. 1(a) and 2(a) for W > . W ∗ w + w L w + w ln( L ) w + w L C S [ ξ KT ] 1.60 1.40 0.34 C S [ ξ ] 2.60 2.50 1.81 C r [ ξ KT ] 2.73 2.38 1.01 C r [ ξ ] 3.72 3.27 2.71TABLE S2. Cost function C X , see Eq. (4), in the Heisenbergmodel. Values of C X are shown for two ergodicity indica-tors X ∈ { S, r } , using correlation lengths ξ and ξ KT fromEqs. (1)-(2). Columns denote different functional forms of W ∗ used in ξ and ξ KT . Results are obtained using the datapoints in Figs. 1(b) and 2(b) for W > . W ∗ w + w L w + w ln( L ) w + w L C S [ ξ KT ] 1.31 1.08 0.46 C S [ ξ ] 1.33 1.36 0.94 C r [ ξ KT ] 1.58 1.37 0.51 C r [ ξ ] 1.73 1.31 1.08 w + w /L and W ∗ = w + w / ln( L ). These two func-tions imply finite critical disorder W ∗ = w in the ther-modynamic limit L → ∞ . Results for the optimal costfunction are compared to the results for the functionalform with a linear drift with L , i.e., W ∗ = w + w L , seeTables S1 and S2. The main result is that the functionalform of the critical disorder W ∗ always yields a lowercost function if a linear drift with L is allowed. More-over, in the case of a linear drift the solution using theKT correlation length ξ KT is always substantially better.In other cases, there is no considerable difference in costfunctions between solutions using ξ KT or ξ . S2. QUALITATIVE COMPARISON
Tables I-II in the main text show quantitative compar-ison of the cost functions for different functional forms ofthe correlation length ξ and critical disorder W ∗ . Herewe complement these results by showing qualitative (vi-sual) comparison of the best data collapses. We focus on0functional forms of the critical disorder with the lineardrift W ∗ ∝ L and in the absence thereof. In the first case,we show the best data collapse as functions of L/ξ KT us-ing W ∗ = w + w L , see Figs. S1(a)-S4(a). In the secondcase, we show the best data collapse as functions of L/ξ using W ∗ = w , see Figs. S1(b)-S4(b). The results agreewith expectations that the better data collapse in termsof the cost function also yields a visually more convincingdata collapse. S3. LOWER BOUND ESTIMATE FOR THECRITICAL POINT
In the main text we discussed scenarios for the large- L dependence of the transition point W ∗ ( L ). We arguedthat our analysis may provide an estimate for the lowerbound of W ∗ if one assumes that the deviation from thelinear drift W ∗ ∝ w L emerges at system sizes that areonly slightly larger that the maximal system size studiedhere, L = 20.For a quantitative analysis we use the fitting func-tion for the transition point W ∗ ( L ) = W ∗∞ tanh( L/L ).We apply the cost function minimization algorithm for S using ξ KT with two free parameters W ∗∞ and b (i.e., b − = b + ), while we fix L . The corresponding cost func-tion C S as a function of L is shown for both models inFig. S5. We define L as the lower bound for L by re-quiring that C S ( L (cid:38) L ) is independent of L . Whilethe extraction of such a lower bound is less ambiguousfor the J - J model than for the Heisenberg model, weassume for both models L ≈
50. We use this value toestimate the lower bounds for the critical point listed inthe main text.
S4. DIFFERENT FORMS OF THE KTCORRELATION LENGTH
In the main text, we studied the scaling collapses of S and r using the KT correlation length ξ KT from Eq. (2)with the identical parameter b = b − = b + on both sidesof the transition. In the scaling analysis of S in the maintext, we get b = 4 .
87 [ b = 3 .
21] in Fig. 1(a) [1(b)], andin the scaling analysis of r , we get b = 3 .
07 [ b = 1 . b − (cid:54) = b + , i.e., b − and b + are independentfree fitting parameters.We focus on functional forms of the critical disorder W ∗ that yield the lowest cost functions, i.e., W ∗ = w + w L and W ∗ = w ∗ ( L ) [see Tables I and II in themain text]. Results for the cost functions are listed inTable S3 for the J - J model and in Table S4 for theHeisenberg model. The first prominent feature is thatthe parameter b + remains very close to the value of b , while the parameter b − may strongly depart from thisvalue. This is related to the property of the ergodicityindicators S and r being nearly a constant below thetransition, and hence being less sensitive to the choice of b − (which then appears as quite irrelevant parameter).The second prominent feature is that the values of pa-rameters w , w and w ∗ ( L ) in the functional forms of thecritical disorder W ∗ remain essentially unchanged. Thisis shown in Fig. S6, where the results from both scenar-ios b − = b + and b − (cid:54) = b + exhibit fairly good agreement.Hence, we expect that our main results remain robustagainst the choice of relation between b − and b + .Finally, we comment on the values of the parameters b (or b − , b + ) in the KT correlation length ξ KT for differentergodicity indicators within the same model, listed in Ta-bles S3 and S4. These values may suggest that b (or b − , b + ) are not identical for S and r within the same model.Nevertheless, we refrain from making any speculationsabout their asymptotic values. TABLE S3. Cost function C X , see Eq. (4), and the param-eters b − , b + of the correlation length ξ KT , see Eq. (2), in the J - J model. Values of C X are shown for two ergodicity in-dicators X ∈ { S, r } and for different scenarios b − = b + and b − (cid:54) = b + . Columns denote different functional forms of W ∗ used in ξ KT . Results are obtained using the data points inFigs. 1(a) and 2(a) for W > .
5. The corresponding values of W ∗ are shown in Fig. S6. W ∗ = w + w L W ∗ = w ∗ ( L )Case b = b − = b + : C S [ ξ KT ], b b − (cid:54) = b + : C S [ ξ KT ] , b − , b + b = b − = b + : C r [ ξ KT ], b b − (cid:54) = b + : C r [ ξ KT ] , b − , b + C X for the Heisenberg modelusing the data points in Figs. 1(b) and 2(b). Parameters arethe same as in Table S3. W ∗ = w + w L W ∗ = w ∗ ( L )Case b = b − = b + : C S [ ξ KT ], b b − (cid:54) = b + : C S [ ξ KT ] , b − , b + b = b − = b + : C r [ ξ KT ], b b − (cid:54) = b + : C r [ ξ KT ] , b − , b + L/ξ KT . . . S (a) W ∗ = w + w L L = 12 L = 14 L = 16 L = 18 L = 20 W . . . S −
100 0 100
L/ξ . . . S (b) J − J model W ∗ = w L = 12 L = 14 L = 16 L = 18 L = 20 FIG. S1. Eigenstate entanglement entropy S for different systems sizes L in the J - J model. We plot results as a functionof L/ξ , using ξ = ξ KT in (a) [assuming b − = b + ≡ b in Eq. (2), as in Fig. 1(a)], and ξ = ξ in (b). We use the transition pointansatz W ∗ = w + w L in (a) and W ∗ = w in (b). The inset shows results as a function of disorder W . Values of the costfunction are C S ( ξ KT ) = 0 .
34 in (a) and C S ( ξ ) = 2 .
80 in (b).
L/ξ KT . . . S (a) W ∗ = w + w L L = 12 L = 14 L = 16 L = 18 L = 20 W . . . S − −
20 0 20 40
L/ξ . . . S (b) Heisenberg model W ∗ = w L = 12 L = 14 L = 16 L = 18 L = 20 FIG. S2. Eigenstate entanglement entropy S for different systems sizes L in the Heisenberg model. We plot results as afunction of L/ξ , using ξ = ξ KT in (a) [assuming b − = b + ≡ b in Eq. (2), as in Fig. 1(b)], and ξ = ξ in (b). We use thetransition point ansatz W ∗ = w + w L in (a) and W ∗ = w in (b). The inset shows results as a function of disorder W . Valuesof the cost function are C S ( ξ KT ) = 0 .
46 in (a) and C S ( ξ ) = 1 .
80 in (b). − L/ξ KT . . . r (a) W ∗ = w + w L L = 12 L = 14 L = 16 L = 18 L = 20 W . . . r −
100 0 100
L/ξ . . . r (b) J − J model W ∗ = w L = 12 L = 14 L = 16 L = 18 L = 20 FIG. S3. Level spacing ratio r for different systems sizes L in the J - J model. We plot results as a function of L/ξ , using ξ = ξ KT in (a) [assuming b − = b + ≡ b in Eq. (2), as in Fig. 2(a)], and ξ = ξ in (b). We use the transition point ansatz W ∗ = w + w L in (a) and W ∗ = w in (b). The inset shows results as a function of disorder W . Values of the cost functionare C r ( ξ KT ) = 1 .
01 in (a) and C r ( ξ ) = 4 .
31 in (b). − L/ξ KT . . . r (a) W ∗ = w + w L L = 12 L = 14 L = 16 L = 18 L = 20 W . . . r − −
20 0 20 40
L/ξ . . . r (b) Heisenberg model W ∗ = w L = 12 L = 14 L = 16 L = 18 L = 20 FIG. S4. Level spacing ratio r for different systems sizes L in the Heisenberg model. We plot results as a function of L/ξ ,using ξ = ξ KT in (a) [assuming b − = b + ≡ b in Eq. (2), as in Fig. 2(b)], and ξ = ξ in (b). We use the transition point ansatz W ∗ = w + w L in (a) and W ∗ = w in (b). The inset shows results as a function of disorder W . Values of the cost functionare C r ( ξ KT ) = 0 .
51 in (a) and C r ( ξ ) = 2 .
16 in (b).
20 30 40 50 60 70 80 L . . . . C S J - J model (a)
20 30 40 50 60 70 80 L . . . C S Heisenberg model (b)
20 40 60 80 L W ∗ ∞
20 40 60 80 L W ∗ ∞ FIG. S5. Cost functions C S for S in the J - J model (a) and in the Heisenberg model (b), as a function of L . We use the KTcorrelation length ξ KT and the fitting function for the transition point W ∗ ( L ) = W ∗∞ tanh( L/L ). Horizontal lines representthe value w if one uses the fitting function W ∗ ( L ) = w L , and the arrows sketch the onset of roughly L -independent costfunctions C S . Insets: the optimal values of W ∗∞ as a function of L . Dashed lines represent w L .
10 12 14 16 18 20 L W ∗ f o r S Solid and dashed lines: b − = b + Dashed-dotted lines: b − = b + J - J modelHeisenberg model (a)
10 12 14 16 18 20 L W ∗ f o r r Filled symbols: b − = b + Open symbols: b − = b + J - J modelHeisenberg model (b) FIG. S6. Ergodicity breaking transition point W ∗ as a function of system size L , for the J - J model (upper part) and theHeisenberg model (lower part). Results for W ∗ are obtained from the best data collapse using a KT correlation length ξ KT fromEq. (2). Solid and dashed lines, and filled symbols, are identical to the ones in Fig. 3 of the main text (i.e., assuming b − = b + in ξ KT ). Dashed-dotted lines and open symbols are obtained by taking independent parameters b − (cid:54) = b + . All lines are results fora transition point ansatz W ∗ = w + w L (with free parameters w and w ), all symbols are results for a function-independenttransition point W ∗ = w ∗ ( L ) [with five free parameters for five different systems sizes L = 12 , , , , S in panel (a) and the level spacing ratio rr