In search of a many-body mobility edge with matrix product states in a Generalized Aubry-André model with interactions
IIn search of many-body mobility edge with matrix product states in a GeneralizedAubry-Andr´e model with interactions
Nicholas Pomata, Sriram Ganeshan,
2, 3 and Tzu-Chieh Wei C. N. Yang Institute for Theoretical Physics and Department of Physics and Astronomy,State University of New York at Stony Brook, NY 11794-3840, United States Physics Department, City College of New York, NY 10031 CUNY Graduate Center, New York NY, 10031
We investigate the possibility of a many-body mobility edge in the Generalized Aubry-Andr´e(GAA) model with interactions using the Shift-Invert Matrix Product States (SIMPS) algorithm(Phys. Rev. Lett. 118, 017201 (2017)). The non-interacting GAA model is a one-dimensionalquasiperiodic model with a self-duality induced mobility edge. The advantage of SIMPS is that ittargets many-body states in an energy-resolved fashion and does not require all many-body statesto be localized for convergence, which allows us to test if the interacting GAA model manifests amany-body mobility edge. Our analysis indicates that the targeted states in the presence of thesingle particle mobility edge match neither ‘MBL-like’ fully-converged localized states nor the fullydelocalized case where SIMPS fails to converge. An entanglement-scaling analysis as a function ofthe finite bond dimension indicates that the many-body states in the vicinity of a single-particlemobility edge behave closer to how delocalized states manifest within the SIMPS method.
Isolated quantum systems are conjectured to equili-brate at the level of a single eigenstate via subsystemthermalization in the absence of a bath. This conjectureis known as the Eigenstate Thermalization Hypothesis(ETH) [1, 2]. Over the past decade, many-body local-ization (MBL) has emerged as a candidate phase whichmaximally violates ETH, where all the eigenstates fail toequilibrate at the subsystem level. Currently, there is aconsensus that MBL exists in one dimension with shortrange interactions [3–5]. Experiments indicate the exis-tence of MBL in a number of platforms [6–8] , althoughhow to unambiguously quantify MBL in an experimentis still work in progress.A natural question then arises as to whether MBLalways represents the most generic violation to ETH,where all eigenstates are non-thermal, or there can becases where only part of the many body spectrum willbe localized. This many-body mobility edge scenario wasoriginally presented in the works of Basko, Aleiner andAltshuler [3] where they found a possible many-body de-localization phase transition at finite temperature. Nu-merical works have observed evidence for a many bodymobility edge [9–17], although finite size effects plaguethe reliability of these results. However, the works of DeRoeck et al. [18] claim to exclude the possibility of any mobility edge using avalanche arguments. More recently,experiments have shown evidence for many-body mobil-ity edge in a shallow lattice limit of the Aubre-Andremodel [19, 20]. It is an open question if the experimentalobservation of a non-ergodic phase is an indication of amore robust violation to ETH or simply a finite-size andfinite-time effect. The question of the presence or absenceof many-body mobility edges remains unresolved.In this letter, we consider the fate of many-body local-ization in the presence of a single particle mobility edge atlarge system sizes. We consider the interacting version of the generalized Aubry-Andr´e (GAA) model, which pos-sesses a mobility edge protected by self-duality in its sin-gle particle spectrum. Machine learning methods have in-dicated the existence of a non-ergodic metal in the centerof the many body spectrum of this model [21]. Recent ex-periments have realized the bosonic version of the GAAmodel in the synthetic lattices of laser-coupled atomicmomentum modes and studied the influence of weak in-teractions on the mobility edge [22]. In order to addressthis question at large system sizes, we use the energy-targeting Shift-Invert Matrix Product States (SIMPS)algorithm of Yu, Pekker, and Clark [23]. We show thatthe SIMPS method should be capable of identifying amany-body mobility edge due to its energy-targeting na-ture. We benchmark the properties of the targeted ma-trix product state in the mobility-edge regime to thatof the convergent fully-localized regime and the fully-delocalized regime where the algorithm is expected tofail to converge. The bond-dimension dependence of en-tanglement scaling with subsystem size shows evidenceof delocalization when compared with how delocalizedstates manifest in the fully delocalized regime of the GAAmodel within the SIMPS method.
Model —
We consider the GAA model [24], with a self-duality protected mobility edge. This model (when usingopen boundary conditions, as we do) has the Hamiltonian H = t L − (cid:88) n =1 ( ψ † n ψ n +1 + ψ † n +1 ψ n )+ 2 λ L (cid:88) n =1 cos(2 πbn + φ )1 − α cos(2 πbn + φ ) ψ † n ψ n , (1)which becomes the standard Aubry-Andr´e model when α = 0, with the phase φ determining a family of “disorder a r X i v : . [ c ond - m a t . d i s - nn ] D ec FIG. 1. The spectrum and inverse participation ratio (IPR)for the noninteracting generalized Aubry-Andr´e model at α = 0 . λ and any t >
0, obtained by diagonal-izing single-particle Hamiltonians. We use a length N = 987chain, selecting a Fibonacci number F and using periodicboundary conditions with b approximated by F /F to min-imize boundary effects. The distinction between delocalized(IPR ∼ /N ) and localized (IPR > /N , increasing to 1)behavior on either side of the self-dual line αE = 2( t − λ )(green) is clear. Disorder strengths λ studied herein in theinteracting case are marked in black. realizations”. We also use the standard choice of b as theinverse golden ratio √ . When t > | α | <
1, thisis determined to be self-dual for energies αE = 2( t − λ ) . (2)As shown in Fig. 1, this condition fixes a mobility edgefor nonzero α at E = α ( t − λ ). Later works consideringan interacting version of this model [25–27], constructedwith the simple addition of a four-fermion term H = H + V (cid:88) n ψ † n ψ n ψ † n +1 ψ n +1 , (3)have analyzed it with exact diagonalization.The main goal of this work is to expand the systemsize substantially using the SIMPS method [23]. Thetradeoff for large system sizes is that the finite bond di-mension cuts off the entanglement of the state. If themany-body state obtained by the SIMPS is localized thenthe entanglement does not scale with the cut size and isunaffected by the increasing bond dimension. However,for a delocalized case the entanglement scales as a func-tion of the cut size and can only become a volume lawif the bond dimension diverges as 2 L/ (for system size L ). However, for a bond dimension χ , the single-cutentanglement entropy of an MPS is limited to precisely log χ , since the bond dimension is exactly the Schmidtrank across a given cut. Thus we trade system-size lim-itations suffered in exact diagonalization methods withfinite-entanglement limitations due to limited bond di-mension. The advantage of this tradeoff is that we canbenchmark states against fully localized and fully delo-calized systems in terms of how their properties scalewith the bond dimension while making finite-size effectsnegligible. The SIMPS method —
Friesdorf et al. [28] have shownthat Matrix Product States (MPS) can efficiently rep-resent excited eigenstates of localized systems and aretherefore an effective means of non-perturbatively ana-lyzing localized systems at large system sizes. To extractMPS approximations of eigenstates, we use the SIMPSalgorithm, which we outline in detail in the Section Aof the supplemental material. SIMPS and other MPSalgorithms can only attempt to diagonalize a system un-der the assumption that it is localized, meaning thatthey will otherwise give “false positives” of relatively low-entanglement states which are not approximations of anyeigenstates and are instead linear combinations of stateswith similar entropies. Indeed, if we assume a typical en-ergy spacing, at a typical energy E , of s ∼ − L E (where L is the system size), an equal combination of n adja-cent states would have energy variance on the order of∆ E /E ∼ n − L (see (A1) of the supplemental mate-rial). At a system size L = 64, such a combination couldstill have energy error below machine precision if it con-sisted of 2 distinct eigenstates. However, our intuitionand the benchmarks we use suggest this is not the case,e.g. a low-energy-error superposition like that would stillhave high entropy and thus could not be replicated as anMPS.We note that, in the similar case of MPS approxima-tions to critical ground states, there exist well-establishedscaling relations [29, 30]. These relations include theasymptotic behavior of the correlation length with re-spect to (a) the single-cut entanglement entropy, (b) thebond dimension, and (c) the energy error relative to thetrue ground state. Such a relation, applied to excitedstates of disordered ergodic systems, would be necessaryin order to distinguish with any certainty the phases wehope to observe. In the absence of such an asymptoticdescription, we attempt to extract empirical relationshipsand benchmarks with respect to fully localized and fullydelocalized cases which can help separate localized andergodic phases. Model Parameters —
The systems we will be primarilyconsidering will have an interaction strength
V /t = 1and a mobility-edge parameter α = 0 .
3. We enforceparticle-number conservation as a global U (1) symme-try in order to restrict to half-filling. As shown in Fig-ures S3, S4, and S7 of the supplemental material, we FIG. 2. A violin plot representing the distribution of energyerror for various disorder strengths; the area of a shape ina given region approximates the frequency of samples withinthe 1188 sample states, and lines are placed at the medians.In all cases, the error substantially decreases with bond di-mension; however, the errors at λ = 0 . λ = 3 . λ = 1 . λ = 1 .
5, meanwhile, only ap-proach either scale at high bond dimension, when the λ = 1 . have demonstrated to our satisfaction that the qualityof states obtained via SIMPS does not vary significantlywith system size for sufficiently large system size (roughly L > L = 64. Wechoose sample “disorder strengths” λ of 0 . , . , . , . φ = πn/
6. In each system, for each of the bonddimensions 10, 14, 20, 25, and 30, we sample 99 targetenergies equally spaced between 0 . ε max < ε < . ε max [31] and again between 0 . ε max < ε < . ε max . In par-ticular this means that for each value used of the “dis-order” strength λ and the bond dimension χ , we have asample size of 1188 states. Benchmark I: SIMPS applied to a fully-localized many–body system —
We begin by applying SIMPS to thedisorder strength λ = 3 .
5, that is, we tune the system to
FIG. 3. Single-cut entropy by cut location (i.e. distance fromthe nearest endpoint) for strong disorder, λ = 3 .
5, at L =64. The entropy has evidently saturated by χ = 30 for thedifferent cuts and energies, and it remains low moving intothe bulk, as is generally expected for the area-law behavior oflocalized states.FIG. 4. Single-cut entropy by cut location for weak disorder, λ = 0 .
5, at L = 64. While there are signs of convergence(to volume-law behavior) very close to the boundary, as seenin greater detail in Fig. 6, in general the entropy is high andfailing to converge. Also notably, it peaks near the bound-ary, at a location l ∼ be far into the region corresponding to single-particle lo-calization. In Fig. 2, we see that the states we find giventhese parameters have very low energy variance ∆ E . Infact, for high bond dimensions χ >
20, ∆ E appears tosaturate at about 10 − , of order comparable to the tol-erance of the subroutines of our SIMPS implementation.We then consider entanglement entropy, as in Fig. 3. Wesee, first, that as a function of bond dimension the en-tropy has also largely saturated by χ = 30 (in fact, inFigures S10 and S11 of the supplemental materials we seethat the entire entropy distribution has largely convergedwith respect to bond dimensions); indeed, the movementwe see before that point is likely attributable to a reduc-tion of bias against higher-entropy states. Moreover, weobserve that neither the single-cut entropy nor the bond-dimension corrections to it grow significantly as we moveinto the bulk of the system, ruling out the possibility thatour evidence of localization can be viewed primarily as afinite-size effect. Benchmark II: SIMPS applied to a fully delocalized sys-tem —
The SIMPS algorithm is supposed to fail withany finite bond dimension for the fully delocalized casedue to volume law entanglement scaling. However, wecan quantify this failure in the form of energy varianceand entanglement scaling with bond dimension. Withinour model, the “disorder” strength λ = 0 .
5, using theabove remaining parameters, corresponds to full delocal-ization in the single-particle case. We find in Fig. 2 thatthe energy errors are very large, eclipsing the values thatwould be predicted by na¨ıvely combining the cutoff errorand density of states. The implication of this is promis-ing: even as the system size becomes large, the algorithmcannot produce pseudo-eigenstates of the delocalized sys-tem which exploit tight energy spacings to exhibit smallenergy error.In Fig. 4, we find rapid growth of entanglement entropyas we move up to 5 sites into the bulk of the system; no-tably, while a failure to converge is apparent away fromthe boundary, near the boundary we see convergence tosomething resembling a volume law. The full entangle-ment distribution is given in Fig. S10 of the supplemen-tal material. Further away from the boundary, we seethe entanglement entropy fall again before settling intoasymptotic behavior; this is evidently an artifact of finiteentanglement given that the peak moves into the bulk aswe increase the bond dimension.
Mobility edge case from SIMPS as compared to bench-marks I (localized) and II (delocalized) —
Now wepresent the main result of this letter. In the presence ofthe single-particle mobility edge, the localization prop-erties of the many-body interacting states can in prin-ciple have four outcomes: (1) all many-body states arelocalized; (2) the many body spectrum has a mobilityedge; (3) all many-body states are delocalized; and (4)the exotic case of a spectrum containing non-ergodic ex-tended states. Even though the SIMPS algorithm can-not unambiguously discriminate among all of these fourscenarios, it can locate the existence of localized statesin the many-body spectrum in an energy-resolved way.Thus we can address the question of whether the many-body spectrum contains any localized states when thesingle-particle spectrum possess a mobility edge withinnumerically accessible bond dimensions.We consider the disorder strengths λ = 1 .
2, for which afull single-particle band is delocalized, and λ = 1 .
5, whichis fully localized (but with longer localization length inbands closer to the critical line of the mobility edge thanin benchmark I, λ = 3 . FIG. 5. Single-cut entropy by cut location for disorderstrengths λ = 1 . λ = 1 . χ = 30 is apparent only in the lower energyrange with λ = 1 .
5. In the higher energy range, meanwhile,we observe for both λ convergence towards the boundary toapparently volume-law behavior as well as a clear peak in theentropy at l ∼
5, as is also seen with λ = 0 . n = 1 may be easier to represent at low Schmidtrank than those observed at n = 2. We also compare Fig. 3,which has a similar dip that does not precede growth into thebulk. entropy scaling up to (cid:96) (cid:39)
4, are more delocalized thanthose near the edge of the band, where we only see a hintof entropy scaling with cut size.In Fig. 6, we plot entanglement entropy versus bonddimension at small cut sizes ( l = 2 , , ,
5) averaged overthe energy windows selected near the band edge and cen-ter. The full entanglement distribution for this case isexplored in greater detail in Fig. S10 of the supplemen-tal material. Within the bond dimensions we have used,we do not observe saturation of entanglement entropy incontrast to the benchmark I. In the absence of such asaturation, we can still use the dependence of entangle-ment on both cut size (i.e., for single cuts, the distancebetween the cut and the boundary) and bond dimension.For λ = 1 . FIG. 6. Average single-cut entanglement entropy, togetherwith the standard error of the mean, as a function of thebond dimension χ for several small cuts in the various casesconsidered. We expect to see logarithmic growth with even-tual convergence either to a constant value (in the localized,area-law case) or a value proportional to the distance from theboundary (in the delocalized, volume-law case). The asymp-totic behavior generally anticipated in the bulk is logarithmicgrowth if the entropy is volume-law and convergence if it isarea-law; near the boundary, in the volume-law case, we ex-pect convergence to a value determined by that volume law. possibility of entanglement scaling convergence for largerbond dimensions. Meanwhile, when λ = 1 .
5, we observeentanglement scaling similar to the delocalized case, atthe band center. However, the scaling at the band edgefor λ = 1 . Conclusions —
We have observed what appears to bea compelling distinction between thermalized and local-ized behavior in an interacting quasiperiodic system at areasonably large system size L = 64. In particular, wefind that we can extract “good” eigenstates with low en-tropy when the strength of the quasiperiodic “disorder”is high; conversely, when it is low, we only find eigen-states of poor quality (as measured by the energy error)whose entanglement entropy moreover increases substan-tially with bond dimension in accordance with a volumelaw. When we compare these two cases with intermedi-ate cases selected for the possibility of seeing a mobilityedge, we find evidence that the many body spectrum ap-pears delocalized for λ = 1 . Acknowledgements —
The authors would like to thankJed Pixley and Bryan Clark for helpful discussions. Weadditionally thank Rutgers University for their hospital-ity during the workshop “Quasiperiodicity and Fractalityin Quantum Statistical Physics”, where part of this workwas completed. S.G. was supported by the National Sci-ence Foundation under Grant OMA-1936351. T.-C.W.and N.P. were supported by the National Science Foun-dation under Grant No. PHY 1915165. [1] J. M. Deutsch.
Quantum statistical mechanics in a closedsystem . Physical Review A, , 2046 (1991).[2] M. Srednicki. Chaos and quantum thermalization . Phys-ical Review E, , 888 (1994).[3] D. M. Basko, I. L. Aleiner, and B. L. Altshuler. Metal–insulator transition in a weakly interacting many-electronsystem with localized single-particle states . Annals ofphysics, , 1126–1205 (2006).[4] J. Z. Imbrie.
On Many-Body Localization for QuantumSpin Chains . Journal of Statistical Physics, , 998–1048 (2016).[5] J. Z. Imbrie.
Diagonalization and many-body localizationfor a disordered quantum spin chain . Physical reviewletters, , 027201 (2016).[6] 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 inter-acting fermions in a quasirandom optical lattice . Science, , 842–845 (2015).[7] P. Bordia, H. P. L¨uschen, S. S. Hodgman, M. Schreiber,I. Bloch, and U. Schneider.
Coupling identical one-dimensional many-body localized systems . Physical re-view letters, , 140401 (2016).[8] B. Chiaro, C. Neill, A. Bohrdt, M. Filippone, F. Arute,K. Arya, R. Babbush, D. Bacon, J. Bardin, R. Barends,et al.
Growth and preservation of entanglementin a many-body localized system . arXiv preprintarXiv:1910.06024 (2019).[9] J. A. Kj¨all, J. H. Bardarson, and F. Pollmann.
Many-bodylocalization in a disordered quantum Ising chain . Physicalreview letters, , 107204 (2014).[10] X. Li, S. Ganeshan, J. Pixley, and S. D. Sarma.
Many-body localization and quantum nonergodicity in a modelwith a single-particle mobility edge . Physical review let-ters, , 186601 (2015).[11] R. Modak and S. Mukerjee.
Many-body localization inthe presence of a single-particle mobility edge . Physical review letters, , 230401 (2015).[12] M. Serbyn, Z. Papi´c, and D. A. Abanin.
Criterionfor many-body localization-delocalization phase transi-tion . Physical Review X, , 041047 (2015).[13] D. J. Luitz, N. Laflorencie, and F. Alet. Many-body local-ization edge in the random-field Heisenberg chain . Phys-ical Review B, , 081103 (2015).[14] S. D. Geraedts, R. N. Bhatt, and R. Nandkishore. Emer-gent local integrals of motion without a complete set oflocalized eigenstates . Physical Review B, , 064204(2017).[15] M. Goihl, J. Eisert, and C. Krumnow. Exploration of thestability of many-body localized systems in the presenceof a small bath . Physical Review B, , 195145 (2019).[16] R. Yao and J. Zakrzewski. Many-body localization inBose-Hubbard model: evidence for the mobility edge .arXiv preprint arXiv:2002.00381 (2020).[17] P. Brighi, D. A. Abanin, and M. Serbyn.
Stability ofmobility edges in disordered interacting systems . PhysicalReview B, , 060202 (2020).[18] W. De Roeck, F. Huveneers, M. M¨uller, and M. Schiulaz.
Absence of many-body mobility edges . Phys. Rev. B, ,014203 (2016).[19] H. P. L¨uschen, S. Scherg, T. Kohlert, M. Schreiber,P. Bordia, X. Li, S. D. Sarma, and I. Bloch. Single-particle mobility edge in a one-dimensional quasiperi-odic optical lattice . Physical review letters, , 160404(2018).[20] T. Kohlert, S. Scherg, X. Li, H. P. L¨uschen, S. D. Sarma,I. Bloch, and M. Aidelsburger.
Observation of many-bodylocalization in a one-dimensional system with a single-particle mobility edge . Physical Review Letters, ,170403 (2019).[21] Y.-T. Hsu, X. Li, D.-L. Deng, and S. D. Sarma.
Machinelearning many-body localization: Search for the elusivenonergodic metal . Physical Review Letters, , 245701(2018).[22] F. A. An, K. Padavi´c, E. J. Meier, S. Hegde, S. Ganeshan,J. Pixley, S. Vishveshwara, and B. Gadway.
Observationof tunable mobility edges in generalized Aubry-Andr \ ’ { e } lattices . arXiv preprint arXiv:2007.01393 (2020).[23] X. Yu, D. Pekker, and B. K. Clark. Finding Matrix Prod-uct State Representations of Highly Excited Eigenstatesof Many-Body Localized Hamiltonians . Phys. Rev. Lett., , 017201 (2017).[24] S. Ganeshan, J. H. Pixley, and S. Das Sarma.
NearestNeighbor Tight Binding Models with an Exact MobilityEdge in One Dimension . Phys. Rev. Lett., , 146601(2015).[25] X. Li, S. Ganeshan, J. H. Pixley, and S. Das Sarma.
Many-Body Localization and Quantum Nonergodicity ina Model with a Single-Particle Mobility Edge . Phys. Rev.Lett., , 186601 (2015).[26] X. Li, J. H. Pixley, D.-L. Deng, S. Ganeshan, and S. DasSarma.
Quantum nonergodicity and fermion localizationin a system with a single-particle mobility edge . Phys.Rev. B, , 184204 (2016).[27] D.-L. Deng, S. Ganeshan, X. Li, R. Modak, S. Mukerjee,and J. H. Pixley. Many-body localization in incommen-surate models with a mobility edge . Annalen der Physik, , 1600399 (2017).[28] M. Friesdorf, A. H. Werner, W. Brown, V. B. Scholz, andJ. Eisert.
Many-Body Localization Implies that Eigenvec-tors are Matrix-Product States . Phys. Rev. Lett., , 170505 (2015).[29] F. Pollmann, S. Mukerjee, A. M. Turner, and J. E.Moore.
Theory of Finite-Entanglement Scaling at One-Dimensional Quantum Critical Points . Phys. Rev. Lett., , 255701 (2009).[30] V. Stojevic, J. Haegeman, I. P. McCulloch, L. Taglia-cozzo, and F. Verstraete.
Conformal data from finite en-tanglement scaling . Phys. Rev. B, , 035120 (2015).[31] Here we use ε to signify the energy density above theground state, such that the spectrum lies between 0 and ε max ; in general we can determine the ground state en-ergy and the maximal energy with relative ease to deter-mine where this range lies.[32] B. Clark. personal communication.[33] M. Hauru and G. Vidal. Uhlmann fidelities from tensornetworks . Phys. Rev. A, , 042316 (2018).[34] S. Bera and A. Lakshminarayan. Local entanglementstructure across a many-body localization transition .Phys. Rev. B, , 134204 (2016).[35] C. G. West and T.-C. Wei. Global and short-range en-tanglement properties in excited, many-body localized spinchains . arXiv e-prints, arXiv:1809.04689 (2018). = (cid:104) ψ i + | O † O | ψ i + (cid:105) (cid:104) ψ i + | O | ψ i (cid:105) FIG. S1. The tensor equation underlying the SIMPS algorithm developed in Yu, Pekker, & Clark[23] : for O = H − E ,update | ψ i +1 (cid:105) by solving for the light-blue tensor T (and then decompose T using SVD). The purpose of this is to optimize (cid:104) ψ i +1 | O † | ψ i (cid:105) subject to the constraint (cid:104) ψ i +1 | O † O | ψ i +1 (cid:105) = 1, but with a Lagrange multiplier, absorbed into T , that can besolved for with normalization. Supplemental Material
A. The SIMPS algorithm
In order to perform the SIMPS algorithm, as withother DMRG-based MPS algorithm, we begin by express-ing the Hamiltonian as an automaton-style matrix prod-uct operator, formed in this case from the Jordan-Wignertransform of the fermionic Hamiltonian, expressed interms of operator-valued matrices as O n = σ − σ + √ V σ z W n σ z − E σ z σ + σ − σ z √ V σ z ,W n = 2 λ cos(2 πbn + ψ )1 − α cos(2 πbn + φ ) , with boundary vectors v L = (0 , , , ,
1) and v R =(1 , , , , | ψ (cid:105) and a target energy E , incorporated into the MPOas above.2. Given an iteration | ψ i (cid:105) , optimize the next iteration | ψ i +1 (cid:105) as follows:3. | ψ i +1 (cid:105) may be initialized randomly, but we have,following a suggestion by Clark [32], initialized itwith | ψ i (cid:105)
4. Site-by-site, optimize | ψ i +1 (cid:105) to satisfy ( H − E ) | ψ i +1 (cid:105) = ψ i : that is, apply the shifted and in-verted Hamiltonian.5. To do so, we represent this equation as the maxi-mization of (cid:104) ψ i +1 | ( H − E ) | ψ i (cid:105) , subject to the con- FIG. S2. A replication of the test run by Yu, Pekker, andClark to eliminate entropy bias of the SIMPS algorithm fromconsideration. The disordered Heisenberg model, with dis-order strength W = 8, is analyzed for 101 disorder realiza-tions on a 10-site spin chain. Each of the 2 L eigenstatesare extracted via exact diagonalization; then SIMPS, withbond dimension χ = 12, is used for 2 L equally-spaced tar-get energies, with post-processing to remove duplicates (if |(cid:104) ψ i | ψ j (cid:105)| > .
3, we exclude the state with greater ∆ E ). Thisproduces 103,424 states via exact diagonalization and 45,917states via SIMPS. The latter is about three times as manystates as in the original test, which explains the difference innoise in that case. As with the original visualization, SIMPSyields significantly fewer states at nearly all entropy ranges,the exceptions being the particularly small ones ( S ≤ . S = log 2, with the difference infrequency being visibly greater at, e.g., S ∼ S ∼ . straint (cid:104) ψ i +1 | ( H − E ) | ψ i +1 (cid:105) = (cid:107) ψ i (cid:107) which willbe uniquely satisfied by ( H − E ) − ψ i . FIG. S3. A violin plot displaying the distribution of energy er-rors, in the systems under consideration within the main textbut for a greater selection of system sizes: for χ ∈ { , , } ,the systems are analyzed at L = 128, and for χ ∈ { , } , thesystems are analyzed at L = 256. Note that the the L = 64data is contained within Fig. 2 of the main text. This plotdemonstrates that, with these parameters, system size doesnot affect the quality of states within the sizes considered.
6. This is done, site-by-site, by solving the diagram-matic equation in Fig. S1 for individual tensors. Wenote that we find it preferable to update two sitesat once (i.e. the tensor being optimized consistsof the contraction of tensors at two sites), espe-cially when enforcing charge/fermion-number con-servation, in order to speed up convergence. Theresulting two-site tensor is then split via SVD toupdate the MPS.7. This may be repeated until | ψ i +1 (cid:105) has converged;alternatively, when initializing ψ i +1 with ψ i , veryfew sweeps (optimizing the tensors at each site)may be conducted per iteration, as the goal of con-vergence is the eigenvalue equation ψ ∝ ( H − E ) ψ ,which should be more accurate after each sweep.8. Repeat until the energy has converged, or until amaximum number of iterations has been reached. FIG. S4. A violin plot displaying the distribution of energyerrors at several system sizes & λ , given bond dimension χ =10, for tests sampling the full spectra of systems In the original work of Yu, Pekker, and Clark [23], theauthors claim that SIMPS is “sampling states at a givenentanglement with the same frequency as ED and hencethere is no systemic bias.” This would be quite remark-able, given the general expectation that entropy may di-verge approaching a transition - in particular, for anybond dimension χ there should be truly localized stateswith entanglement entropy at some cut in excess of themaximum log χ . (In fact, in a good approximation of aphysical state, it may be expected that the entropy ceil-ing should be even less than that absolute maximum, as,for example, was shown for infinite MPS approximationsfor ground states of critical spin chains by Pollmann etal. [29].) Although they acknowledge a “failure of SIMPSto find high-quality eigenstates in [the] near-ergodic andergodic regime”, they do not explain why there shouldbe a hard boundary between regimes near to and farfrom ergodicity. Moreover, in the data they provide asevidence for this claim (Fig. S2), the divergence of theproportion of SIMPS states from that of ED states athigher entanglement entropy seems apparent (if small),and likely statistically significant. To confirm statisticalsignificance, we replicate the test they use to producethese data as faithfully as possible, yielding data thatclearly replicates the major features of this figure, par- FIG. S5. The distribution of energy errors at several(smaller) bond dimensions, given system size L = 128,for tests sampling the full spectra of systems with λ ∈{ . , . , . , . , . , . , . } ticularly a divergence between sampling rates at higherentropies, in Fig. S2.In addition to favoring true low-entropy eigenstates,we have noted that the SIMPS algorithm will produce“false” eigenstates when no low-entropy eigenstates areavailable, as is evidenced by the fact that the algorithmproduces any states at all within the presumed ergodicregime. To attempt to constrain the false eigenstates weobserve, we may try to approximate a worst-case scenarioby supposing that there exist n consecutive eigenstates,of some separation s : that is, taking the crudest possibleapproximation, the energies take the form E k = E + ks .Then the energy variance would be∆ E ≡ (cid:104) H (cid:105) − (cid:104) H (cid:105) = (cid:104) ( H − E ) (cid:105) − (cid:104) H − E (cid:105) = n − (cid:88) k =0 s k n − (cid:32) n − (cid:88) k =0 skn (cid:33) = n ( n − s . (A1)We do not presuppose the order of magnitude of n , as thepossible entropy reduction from such a combination ishighly dependent on the nature of the eigenstates them-selves. We may, however, presume that the worst-caseenergy spacing s is of the order 2 − L , such that FIG. S6. Average single-cut entanglement entropy, togetherwith the standard error of the mean, as a function of the bonddimension χ for several larger cuts (as in Fig. 6 of the maintext but with the subsequent four cuts). ∆ E ∼ n L L . (A2) B. Additional datasets
In addition to the dataset described and referenced inthe main text, we have two additional datasets that wewill reference on occasion in these Supplemental Mate-rials. As in the main text, each uses the Hamiltoniandefined by (1) and (3), with α = 0 . t = V = 1, and b = √ , and with protection of U (1) symmetry to re-strict to half-filling. • The exploratory trials: For disorder strengths to0.5, 1.0, 1.5, 2.0, 2.5, 3.0, and 3.5 (that is, λ = n/ n from 1 to 7), we take 24 disorder realiza-tions (that is, values of φ in (1)). With bonddimension χ = 10, we analyze systems of size L = 16 , , , χ = 6 and χ = 14 to systems of size L = 128. In each case, we select 400 target energiesthat encompass the full energy spectrum (notingthat this means we will see, and reject, a numberof copies of the states with lowest and highest en-ergy). We have used this dataset to produce Fig-ures S4, S5, and S7. In the latter, we also havea subset of those conditions, namely χ = 10 and χ = 14 for L = 128, with at λ = 1 . FIG. S7. Entropy histograms for various values of λ , given various bond dimensions and system sizes; this is useful todemonstrate that the effect of system size is minimal for L ≥
64. Vertical lines mark S = ln χ , the maximum possible entropyat a given bond dimension. • In Fig. S3, we display results from trials using thesame parameters as those in the main text, exceptwith (for χ = 10, χ = 14, and χ = 20) L = 128in addition to L = 64, as well as (for χ = 10 and χ = 14) L = 256. C. Energy error
We have presented our primary results on the energyerror in Fig. 2; we show additional results for larger sys-tem sizes in Fig. S3. We here note briefly how it is calcu-lated. In particular, it emerges fairly easily from SIMPScalculation: contracting the transfer matrices used to ob-tain the LHS of Fig. S1 gives (cid:104) ( H − E ) (cid:105) , and (cid:104) H − E (cid:105) (which does not come pre-calculated) is easily computedthrough the simpler contraction of the RHS of the same(replacing | ψ i (cid:105) with | ψ i +1 (cid:105) .)In Figs. S4 and S5, we introduce an additional setof simulations, in which we sampled the entire spec-trum for a substantial number of “disorder” strengths λ , with bond dimensions among { , , } and lengths(for χ = 10; otherwise L = 128) in { , , , } , toexamine how our results scale with bond dimension andwith system size. D. Single-cut entanglement entropy
By keeping the MPS in (bi)canonical form, we areable to extract entanglement entropies directly from theSchmidt coefficients which are stored as part of theansatz. We have explored how average entropies, at givendistances from the boundary, scale with the bond di-mension χ in Figs. 6, 3, 4, and 5 of the main text. Inparticular, in Fig. 6, we examined the entropy scalingat several specific cuts; we repeat this further into thebulk in Fig. S6, and then go deeper by examining thecorresponding entropy distributions in each case in Fig-ures S10 and S11, respectively. In Fig. S7, we take adifferent approach and examine the entropy distributionfor all cuts at various λ , L , and χ .One feature that is clearly visible in many of thesehistograms is a peak at S = ln 2, corresponding to two-site resonances. We confirm that this is not a numericalartifact using exact diagonalization in Fig. S8.1 FIG. S8. A histogram of all single-cut entropies, for all states,of 12-site systems with λ = 2, given 1920 “disorder realiza-tions” (i.e. phases φ ). This demonstrates a very clear reso-nance peak at S = ln 2, marked by the dotted line. E. Energy wandering
Another quantity we may use to diagnose the good-ness of states is the so-called “energy wandering”, thedifference between the energy of a state and the targetenergy used to obtain it. The idea behind using this isto determine whether or not approximate eigenstates ofadequate quality are sufficiently common. In Fig. S9 wecompare the distribution of values of (cid:107) E − E (cid:107) with thatof ∆ E . F. Uhlmann fidelities
Inspired by, and using methods based on, [33], we com-pute Uhlmann fidelities: if the reduced density matrix ofan eigenstate ψ i on a segment A is ρ i,A = tr A | ψ i (cid:105)(cid:104) ψ i | ,then the Uhlmann fidelity between ψ i and ψ j on A is F = tr (cid:113) √ ρ i,A ρ j,A √ ρ i,A . (F1)Ideally, • In the localized case, the distribution of these quan-tities will be determined by so-called “l-bits”: if ψ i and ψ j differ on an l-bit whose support is within A , then F = 0; if they agree on all l-bits mostlysupported within A , then F ∼
1; and intermedi-ate values will only occur when there are l-bits onwhich ψ i and ψ j differ that have significant supportboth inside and outside of A . • In the fully ergodic regime, where ρ i,A should be fully determined by the energy (astr A exp( − β ( E i ) H ), with β ( E i ) the inverse-temperature corresponding to E i ), we expect aless discrete distribution of F , with values con-tinuously dependent on the energy difference andstochastically dependent on the choice of region A . FIG. S9. Scatter plot of energy error versus energy wander-ing, to attempt to determine whether a relation exists be-tween quality of states & how far from the target energy thealgorithm must “wander” to find it
In Fig. S12, we examine the distribution of Uhlmann fi-delities in various systems considered, for various sizes ofregion A . We see that the typical behavior, in the local-ized case, is a bimodal distribution with one narrow peakat 0 corresponding to cases differing on l-bits supportedin A and another narrow peak at 1 corresponding to casesagreeing on l-bits overlapping A , with the former shrink-ing and the latter growing both as the size of A increases(so that F = 1 would require agreement on more l-bits)and as we move to the higher-energy band, which we ex-pect to have larger localization length. In the delocalizedcase, meanwhile, we see a broad unimodal distributionwhose peak, in addition to lowering as the size of A andthe energy of the band increase, raises as the bond dimen-sion increases, suggesting an increase in similarity as theaccuracy improves. (It is not truly unimodal however; asmall peak at F = 0 which narrows with increasing bonddimension suggests that the pseudo-eigenstates obtainedin this case do sometimes have features which resemblel-bits.In analyzing the λ = 1 . λ = 1 . • For the most part, the distribution in the higher-energy band is closer to a unimodal distribution likethe one seen in the delocalized case; in (a) a secondpeak close to F = 1 in the width-3 distributionexists but grows less distinct with increasing bonddimension • The distribution in the lower-energy band is moreconsistently bimodal, although with lower maximaat nonzero fidelity. • In particular, as bond dimension increases the dis-tributions in the higher-energy band appear to con-verge towards a unimodal distribution; this seemsto be the case, though it is less clear, for the lower-energy band when λ = 1 .
2. However, for the lower-energy band with λ = 1 .
5, we see apparent con-vergence toward a bimodal distribution in (a) and(b).
G. Localization Lengths
We follow [34, 35] in using two measures of entangle-ment between two qubits, negativity and concurrence, toattempt to estimate the localization length, which shoulddiverge approaching a localization transition or mobilityedge. In particular, given concurrence values C i,j or neg-ativity values N i,j between two sites of a given state, wefit the nontrivial values to (see eqs. 21 and 22 of [35]) C i,i ± n = k C ± i exp( − n/ζ C ± i ) N i,i ± n = k N ± i exp( − n/ζ N ± i ) . In Fig. S13, we take, for each state, an average of allthese 1 /ζ C ± i and, separately, 1 /ζ N ± i , weighted by σ /ζ .While this confirms some basic expectations – the local-ization lengths of eigenstates with λ = 3 . λ = 1 . λ = 1 . λ = 0 .
5. We must therefore conclude that we will not beable to perform much meaningful analysis on this data.3 (a) (cid:96) = 1(b) (cid:96) = 2(c) (cid:96) = 3(d) (cid:96) = 4(e) (cid:96) = 5
FIG. S10. Histograms of single-cut entanglement entropy by cut position, corresponding to the data in Fig. 6 of the main text (a) (cid:96) = 6(b) (cid:96) = 7(c) (cid:96) = 8(d) (cid:96) = 9(e)half-cut FIG. S11. Histograms of single-cut entanglement entropy by cut position, corresponding to the data in Fig. S6 of the maintext, as well as (e) histograms of half-cut entanglement entropy (i.e. (cid:96) = 32 for L = 64) (a)width 3(b)width 5(c)width 8 FIG. S12. Violin plots of Uhlmann fidelities among systems with various disorder strengths and at different energy ranges,for segments of width (a) 3, (b) 5, and (c) 8. All fidelities calculated are between states with the same disorder strength anddisorder sample , bond dimension, and energy range. All segments of the given length within the 64-site chain are considered.6