Non-ergodic delocalized phase in Anderson model on Bethe lattice and regular graph.

# Non-ergodic delocalized phase in Anderson model on Bethe lattice and regular graph.

V. E. Kravtsov Abdus Salam International Center for Theoretical Physics, Strada Costiera 11, 34151 Trieste, Italy L. D. Landau Institute for Theoretical Physics, Chernogolovka, 142432, Moscow region, Russia    B. L. Altshuler Physics Department, Columbia University, 538 West 120th Street, New York, New York 10027, USA    L. B. Ioffe LPTHE - CNRS - UPMC, 4 place Jussieu Paris, 75252, France National Research University Higher School of Economics, Moscow, Russia L. D. Landau Institute for Theoretical Physics, Chernogolovka, 142432, Moscow region, Russia
###### Abstract

We develop a novel analytical approach to the problem of single particle localization in infinite dimensional spaces such as Bethe lattice and random regular graph models. The key ingredient of the approach is the notion of the inverted order thermodynamic limit (IOTL) in which the coupling to the environment goes to zero before the system size goes to infinity. Using IOTL and Replica Symmetry Breaking (RSB) formalism we derive analytical expressions for the fractal dimension that distinguishes between the extended ergodic, , and extended non-ergodic (multifractal), states on the Bethe lattice and random regular graphs with the branching number . We also employ RSB formalism to derive the analytical expression for the typical imaginary part of self-energy in the non-ergodic phase close to the Anderson transition in the conventional thermodynamic limit. We prove the existence of an extended non-ergodic phase in a broad range of disorder strength and energy and establish the phase diagrams of the models as a function of disorder and energy. The results of the analytical theory are compared with large-scale population dynamics and with the exact diagonalization of Anderson model on random regular graphs. We discuss the consequences of these results for the many body localization.

## I Introduction

Recent progress in understanding the dynamical processes of mesoscopic and macroscopic isolated disordered quantum many-body systems is based on the concept of Many-Body Localization (MBL) BAA (); Ogan-Huse (); Ogan-Huse1 (): the many-body eigenstates can be localized in the Hilbert space in a way similar to the conventional real space Anderson localization Anderson () of a single quantum particle in a quenched disorder. Depending on the temperature (total energy) or other tunable parameters the system can find itself either in many-body localized or in the many-body extended phase. In the former case the system cannot be described in terms of conventional Statistical Mechanics: the notion of the thermal equilibrium loses its meaning, as not all positions in the Hilbert space are reachable.

There are reasons to believe that the violation of the conventional thermodynamics does not disappear with the Anderson transition from the localized to the extended state LevPino (); PinoKrav (): in a finite range of the parameters one expects the appearance of a non-ergodic extended phase for which the conventional theory is inapplicable.

In a many body problem the number of the states connected with a given one in the -th order of the perturbation theory in the interaction increases exponentially or faster with .AGKL () Similar situation takes place in the problem of single particle localization on hierarchical lattices such as the Bethe lattice (BL) or random regular graphs of connectivity where the number of sites at a given distance increases exponentially with distance, . One thus believes that these problems might be viewed as toy models of the many body localization. The rapid growth of the number of sites with distance from a given site is very important feature of these problems that distinguishes them from the one-particle Anderson problem in a finite-dimensional space in which the number of sites growth as a power-law. In particular, the slow growth of the number of sites with distance implies that it cannot compensate the exponential decay of the tunneling amplitude with distance. Thus, in this case the resonances either appear at short distances or not at all. This is the reason why in finite dimensional localization all extended quantum states are ergodic and the ergodicity is violated only at the critical point of Anderson transition, which is manifested by the multifractality of the critical quantum states Mirlin-rev ().

Recent numerical studies Biroli (); Our-BL (); AltshulerCuevasIoffeKravtsov2016 () of the Anderson problem on a random regular graph (RRG), which is known to be almost indistinguishable from the Bethe lattice at short length scales, brought up strong evidence in favor of the existence of the non-ergodic phase: the eigenfunctions were found to be multifractal with the fractal dimensions depending on disorder. It was also suggested that the transition (referred to below as ergodic transition) from the extended ergodic (EE) to the non-ergodic extended (NEE) phases is a true transition as evidenced by the jump in the fractal dimensions rather than a crossover AltshulerCuevasIoffeKravtsov2016 ().

Existence of NEE phase and the transition from NEE to EE states has been recently proven KravtsovKhaymovichCuevas2015 (); RP-Bir () for an apparently different model, the random matrix theory with the special diagonal, suggested in 1960 by Rosenzweig and Porter (RP) RPort () and generalized in Ref.KravtsovKhaymovichCuevas2015 (). As we explain below, the property that unifies both models is the self-consistent equations for the Green’s function suggested for the Bethe lattice by Abou-Chakra, Thouless and Anderson AbouChacAnd (). These equation are valid for the Bethe lattice with any connectivity due to its loopless, tree structure. However, being a kind of self-consistent theory, these equations are also valid exactly for the RP model due to its infinite connectivity in the thermodynamic limit.

An important boost for the interest to single-particle localization problem on BL comes from the recent work Biroli-time-dep () that proposes a mapping of dynamical correlations function in the full many body problem onto single particle correlation functions on RRG and studied them numerically. The results can be applied to the spin correlation function Alet () and the time-dependent even-odd site imbalance Bloch (); Bloch-Science (). The power-law time dependence of the single particle correlation functions on RRG with the exponents that continuously depend on disorder implies similar power-law time-evolution of the corresponding correlation functions in the many body problems Alet (); Bloch-Science (); Bloch (). This behavior is indeed observed Bloch () in a system of interacting cold atoms in a disordered optical lattice.

In this paper we develop an analytical approach to the non-ergodic phase of the Anderson model on the large-connectivity Bethe lattice and random regular graphs which is based on the replica symmetry breaking, the preliminary short version of this paper can be found inArXiv16 ().

The main result of the work is the behavior of the fractal dimension of the wave functions summarized in Fig. 1. Fig. 1a shows the analytical (RSB) and population dynamics results for the fractal dimension for the Bethe lattice. The former predicts that the dimension is a smooth function which varies from to as varies in the range of proving existence of non-ergodic extended phase (NEE). This phase terminates at the Anderson localization transition at large disorder and at the ergodic transition at low disorder. The results of population dynamics coincide with RSB theory for , and corroborate the existence of the non-ergodic extended phase. However, at small the population dynamics predicts gradual crossover to as whilst RSB predicts the transition at . We discuss the origin of the discrepancy and the region of the validity of these results in section XVI. Briefly, we expect that in the bulk of a large Bethe lattice obtained in population dynamics is valid for all in the limit of infinite size , so that the non-ergodic phase survives to the lowest . For RRG the ergodic transition from NEE to EE phase might happen at sufficiently small disorder at which RRG and BL are no longer equivalent.

Fig. 1b shows the dependence of the typical imaginary part of a single-site Green’s function at the band center in the Anderson thermodynamic limit AbouChacAnd () as a function of disorder. Everywhere in the domain of NEE phase it is smaller than averaged over disorder , it approaches it only at small , in ergodic or almost ergodic phase. Furthermore, it becomes exponentially small near the Anderson localization transition.

To verify the results of the analytical theory developed in this paper, we further develop the population dynamics (PD) method by exploiting, in addition to the standard equilibrium PD, a new inflationary PD formalism introduced previously in AltshulerCuevasIoffeKravtsov2016 (). This formalism corresponds to the unusual ("inverted-order") thermodynamic limit (IOTL) in which the bare energy level width prior to the system size , which allows to compute the fractal dimension of a single wave function. The agreement with analytical result appears to be very good for which is not very small: , see Fig. 1a. In a separate computation we verified the results of the analytical theory for the critical behavior of in the conventional thermodynamic limit ( first) in the vicinity of Anderson transition. The results of population dynamics and analytical theory are shown in Fig. 1b. They unambiguously show that as approach from the metallic side of the Anderson transition.

The plan of the remainder of the paper is the following. In section II we define the support set of random wave functions and give a definition of the NEE phase in terms of the scaling of the support set volume with the total volume. In section III we review the behavior of the typical local density of states in the conventional Anderson Thermodynamic Limit (ATL) and in the Inverted Order Thermodynamic Limit (IOTL) which allows to distinguish between the EE and the NEE phases. Section IV formulates the models while section V gives the basic equations for the Green’s functions in these models. In section VI we describe the new method of Inflationary Population Dynamics and derive a relationship between the increment of exponential inflation of the typical imaginary part of Green’s function and the fractal dimension . We derive the basic equations for the one-step Replica Symmetry Breaking in sections VII-IX. In section X we use the basic symmetry of the problem to derive a new algebraic equation for the critical disorder at the localization transition on the Bethe lattice and Random Regular Graphs. This simple equation considerably improves the accuracy of compared to the classical result of Ref.AbouChacAnd () for the small branching numbers . In section XI we derive, within the one-step replica symmetry breaking method, the analytical results for the fractal dimension as a function of disorder at the branching number and compare them with the results of inflationary population dynamics and exact diagonalization on random regular graphs. In section XII we apply the large- approximate solution for to the Generalized Rosenzweig-Porter random matrix ensemble and re-derive the dependence of on the control parameter in the NEE phase that shows a continuous transition to the EE phase discovered earlier in KravtsovKhaymovichCuevas2015 (). In section XIII we derive, in the framework of one-step replica symmetry breaking, the dependence of the typical imaginary part of Green’s function on disorder strength in the Anderson thermodynamic limit and compare it with equilibrium population dynamics numerics in section XIV. In addition, in section XIV we present the results of population dynamics numerics for the correlation function of at different energies and relate it with the noise in interacting spin systems. In section XV we reformulate the condition for the localization and ergodic transitions in terms of the Lyapunov exponents and find the corresponding expressions within the one-step RSB. These analytical results are compared with the population dynamics for the Lyapunov exponents in section XVI. The obtained behavior of Lyapunov exponent allows us to estimate the contribution of large loops present in RRG and obtain the upper bound for the applicability of the analytical and population dynamic results to RRG in section XVI. The phase diagram in the energy-disorder plane for the Bethe lattice and random regular graphs is presented and discussed in section XVII. Section XVIII compares the main results of this paper with the results of other workers on the critical behavior of and on the existence of non-ergodic phases in finite lattices and many body systems. The main results of the paper are summarized in Conclusion, section XIX.

The paper contains six Appendices which provide the details of the computations and proofs.

## Ii Support set of random wave functions

A fundamental concept that distinguishes non-ergodic extended states from the ergodic ones is the support set of wave functions. Suppose that wave function amplitudes are ordered:

 |ψa(1)|2≥|ψa(2)|2≥...|ψa(N)|2

and obey the normalization condition:

 N∑i=1|ψa(i)|2=1. (1)

In order to define the support set we introduce the integer valued function that gives the number of sites needed to get the normalization condition with a prescribed accuracy :

 Mϵ∑i=1|ψa(i)|2 ≤1−ϵ (2) but (3) Mϵ+1∑i=1|ψa(i)|2 >1−ϵ.

The manifold of sites contributing to the sum in the left hand side of (2) constitutes a support set of the wave function and the number is the support set volume Our-sset (). The wave function is localized if is finite for any fixed in the limit . It is extended and ergodic (EE) if is finite in this limit, and it is extended, non-ergodic if while .

A special class of extended non-ergodic states are multifractal states for which , where . Because the probability to find the particle in a given state on a site of the support set is almost unity, the typical value of the wave function on the support set sites is . In contrast, the typical value of the wave function on a generic site is much smaller; it is controlled by an exponent : . Qualitatively, the sites of the support set are in resonance with each other while the sites outside the support set have very different energies and are connected to the support set only in high orders of perturbation theory which makes their wave function very small.

One can show Our-sset () that the exponent coincides with the fractal dimension . The latter is determined by the "Shannon entropy" which leading term in the limit is :

 ⟨∑i|ψa(i)|2ln(|ψa(i)|−2)⟩=D1lnN. (4)

## Iii Distribution of LDoS

As argued in AbouChacAnd () the information on the character of wave functions can be extracted from the probability distribution function (PDF) of the generalized local density of states (LDoS) :

 ρi(E)=1π∑a|⟨i|a⟩|2η(E−Ea)2+η2, (5)

where is the Green’s function and is the broadening of energy levels.

For localized wave functions is singular: for all it vanishes in the limit . In contrast, in case of the extended ergodic wave functions the result for the distribution depends on the order of limits. If limit is taken first, before , one gets a stable non-singular AbouChacAnd () with the typical of the order of the averaged value :

 ρtyp≡exp[⟨lnρ⟩]∼⟨ρ⟩. (6)

In the following we shall refer to this order of limits as Anderson thermodynamic limit (ATL).

In this paper we show that the non-ergodic extended states on BL are characterized by a non-singular but extremely broad distribution of , such that is non-zero but parametrically smaller than in ATL:

 0<ρtyp≪⟨ρ⟩. (7)

It becomes exponentially small

 ρtyp∼exp[−a/(1−W/Wc)] (8)

at the disorder strength approaching the Anderson localization transition. Notice that (7) does not by itself prove the existence of a distinct non-ergodic phase because smoothly changes as is decreased and becomes at .

The way to unambiguously characterize the type of extended wave functions is to consider the opposite order of limits when after (inverted-order thermodynamic limit, IOTL). In this limit becomes smaller than the mean level spacing , so that the typical is dominated by the state closest in energy to the observation energy . In this regime can be estimated from (5):

 ρtyp∼ηδ2|ψ|2typ. (9)

For ergodic states , so that one gets

 ρ(erg)typ∼ηN. (10)

For multifractal wave functions the ’main body’ of the wave function is located on its fractal support set, so its typical amplitude is very small; it is characterized by a non-trivial exponent . For this class of states we obtain:

 ρ(mf)typ∼ηN2−α0, (11)

where the exponent .

The exponent is in fact equal to the anomalous dimension, . There are two ways to prove it. One is to use the physical arguments of Ref.AltshulerCuevasIoffeKravtsov2016 () that discussed the crossover from the linear behavior of at small to independent at large and argued that it should occur when which is the distance between the levels in the support set. Because the crosses over from linear function (11) to the constant at these exponents should be equal. Another argument relies on Mirlin-Fyodorov symmetry of multifractal spectrum MF-sym () which gives (see Appendix A):

 2−α0=D1=D, (12)

and consequently:

 limη→0ρtypη∼ND. (13)

We conclude that behavior of in the IOTL gives directly the support set dimension .

The scaling behavior of the typical value of the wave function amplitude is much easier to determine numerically in the exact diagonalization of finite graphs than the values of the wave function dimensions . In particular, its extrapolation to the graphs of infinite size is much more reliable than that of for . The reason for that is that for broad distributions the average of are controlled by the distribution tail which is more sensitive to the finite size effects and insufficient statitstics than the distribution main body.

The behavior of the distribution function of at in the IOTL does not contain useful information. The reason for it is the presence of the factor in the definition of that implies that at large . This power law dependence simply reflects the probability to find the state very close to the given energy. It serves, however, a useful check on the consistency of the analytic approximations developed below.

## Iv The model

We are considering the Anderson model on a graph with sites:

 ^H=N∑iεi|i⟩⟨i|+N∑i,j=1tij|i⟩⟨j| (14)

Here labels sites of the graph and is connectivity matrix of this graph: equals to if the sites and are connected,otherwise . This class of models is characterized by the on-site disorder: are random on-site energies uniformly distributed in the interval . For the random regular graph (RRG) all sites are statistically equivalent and each of them has neighbors, while the Cayley tree is a directed, hierarchical graph: each site is connected to neighbors of the previous generation and to one site on the next generation. The common feature of both graphs is the local tree structure. The difference is that the Cayley tree is loop-less, while RRG contains loops. Whilst the number of small loops is only for the whole graph, so they cannot have an effect on its properties, the long loops might be more dangerous. A typical random path starting from a site comes back to this site in steps, so a typical loop has the length equal to the graph diameter . Another important feature is that the finite Cayley tree is statistically inhomogeneous: it has a root and a boundary where a finite fraction of states is located.

## V Abou-Chakra-Thouless-Anderson equations

For a general lattice one can write self-consistent equations for the two point Green’s functions . In the absence of the loops, it is possible to derive self-consistent equations for the single site Green’s functions, and where the latter denotes single site Green’s function with the bond removed:

 Gi→k =1E−ϵi−∑j≠kGj→i (15) Gi =1E−ϵi−∑jGj→i (16)

For the Bethe lattice one can introduce the notion of generations: each site of a given generation, , is connected to ancestors (generation ) and descendant (generation ) and focus only on the Green’s functions in which is descendant of : where

 G(ℓ+1)i(E)=1E−ϵi−∑j(i)G(ℓ)j(E) (17)

where are ancestors of .

The equations Eq.(17) are under-determined: the pole-like singularities in the right hand side of this equation might be regularized by adding an infinitesimal imaginary part to similar to Eq.(5). In what follows we will mostly assume the imaginary part always added to in Eq.(17).

At sufficiently small the recursion Eq.(17) might become unstable with respect to addition of . This instability signals the Anderson transition (AT) point and persists everywhere in the extended (both EE and NEE) phases.

One can use the recursive relation Eq.(17) to find the stationary distribution of . This approach was first employed by the authors of the seminal paper AbouChacAnd () who used it to prove the existence of the localized phase on Bethe lattice and to determine the critical disorder of the AT. Recently we have generalized it to identify the non-ergodic phase on Bethe lattice AltshulerCuevasIoffeKravtsov2016 ().

## Vi Inflationary population dynamics and fractal dimension D1

The recursive procedure Eq.(17) is the basis for the recursive algorithm known as population dynamics (PD) pop-dyn () that can be used in two versions. In the linear version we assume infinitesimally small and . In this regime the typical imaginary part increases exponentially with the number of recursion steps in Eq.(17) but remains proportional to :

 ρtyp(ℓ)∝ηeΛℓ, (18)

where is the corresponding increment.

This non-stationary ("inflationary") regime of PD holds at the initial stage of iteration for any sufficiently small . It should be contrasted with the conventional stationary regime when the number of iterations is sufficiently large for the non-linear in terms to become relevant, and the distributions of and reach their equilibrium which is independent of as . This regime exactly corresponds to the Anderson thermodynamic limit.

At an infinitesimal the stationary regime is reached only in the limit. At the same time, the number of generations in a finite graph is limited by the size of the graph. For instance it is equal to for a finite Cayley tree with sites. A part of RRG corresponding to generation with a common ancestor is equivalent to the tree for . Thus, the exponential growth should persist up to . At larger distances the loops have to be taken into account. Consider an iteration of the equations along the loop that corresponds to the positive exponent : . The recursion in this case would predict the infinite growth of the which is impossible for a single finite loop. Clearly, in order to get the correct result one has to stop the recursion after one turn. For short loops the number of loops of length is Bollobas () , so that the probability that a site belongs to a loop of length is . Thus the loops with have vanishingly small probability. A typical loop has the length , where is the graph diameter Bollobas (). Therefore in order to avoid spurious feedback onto itself, the recursion typically has to stop after a number of steps equal to the graph diameter :

 ℓt(N)=d≈lnN/lnK. (19)

The reasoning above neglects the constructive interference between different paths leading to the same point, we shall discuss the consistency of this assumption below in Section XVI. Here we notice that the prescription Eq.(19) gives a correct result for the fractal dimensions in the Rosenzweig-Porter random matrix ensemble earlier obtained in KravtsovKhaymovichCuevas2015 (). Note that this random matrix theory can be mapped on a graph where every site is connected to any other site directly and thus the corresponding graph diameter is , despite , so that there are plenty of loops on this graph!

Combining (18) and (19) one obtains for a finite RRG:

 ρtyp∼ηNΛ/lnK. (20)

Now, comparing Eq.(20) with Eq.(11) we see that the inflationary PD corresponds to the inverted-order thermodynamic limit (IOTL), and using Eq.(12) we obtain AltshulerCuevasIoffeKravtsov2016 () for the fractal dimension :

 D1=ΛlnK. (21)

## Vii Large connectivity approximation

The increment in Eq.(18) was computed numerically for in AltshulerCuevasIoffeKravtsov2016 () using a specially designed non-stationary population dynamics algorithm. The results unambiguously show that is a continuous function of disorder that vanishes and changes sign at the AT point and grows almost linearly as decreases below . At low the function flattens out, so that never exceeds . This is the simplest test for consistency of Eq.(21), as cannot exceed 1. We discuss the behavior at small in more detail in section XI.

In this section we derive an analytic expression for and in the case of large connectivity . Linearizing the right hand side of (17) in we obtain:

 IG(ℓ)i(E)=∑j(i)IG(ℓ−1)j(E−ϵi−RΣ(ℓ)i)2 (22)

The general method for the solution of the equations of this type was developed in Derrida () that employed mapping to traveling wave problem. More compact solution uses replica approach and one step replica symmetry breaking (see e.g. IoffeMezard (); FeigelmanIoffeMezard ()).

We begin with the expression for - (and -) dependent , where

 Z(E)=∑Pℓ∏k=11(E−ϵ(k)P)2. (23)

In Eq.(23) determines a path that goes from an initial point in generation 1 to a point in generation (), and is the random on-site energy renormalized by the real part of self energy on this path in the -th generation.

The function has a meaning of the free energy of a polymer on the Bethe lattice Derrida ()with unusual disordered site energies . In order to compute it we use the replica method:

 Λ(E) =limn→0 1nℓ(¯¯¯¯¯¯¯Zn−1).

Where can be written as:

 ¯¯¯¯¯¯¯Zn=¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯∑P1,...Pnn∏a=1ℓ∏k=11(E−ϵ(k)Pa)2. (24)

We refer to as entries and the product of entries in (24) as a path. Without replica symmetry breaking (RSB) there would be pathes contributing to , each path containing entries. The one step replica symmetry breaking implies that the main contribution is given by paths where these entries are grouped into groups of identical entries each, considering the contribution of different groups as statistically independent (the "RSB ansatz"). The RSB solution for the increment is obtained by minimization of with respect to :

 Λ(E)=minmΛ(E,m)≡Λ(E,m0). (25)

One step RSB gives exact result for the free energy, , (24) for any distribution function of entries provided that these entries are not correlated between different sites. In order to understand the reason for this we define the ’free energy’ of individual path by

 fP=1ℓℓ∑k=1ln(E−ϵ(k)Pa)2 (26)

The partition sum, is controlled by the path that corresponds to the minimal energy, . The number of paths with energies larger than grows exponentially, , so that main contribution to comes from the thermodynamically large number of paths with . The large number of paths implies that is self averaging and thus can be computed straightforwardly. On the other hand, all these paths have the same as the optimal one in thermodynamic limit which allows one to extract the value of from this computation. These arguments implicitly assume that the two paths are either fully correlated (if they pass through the same point) or completely uncorrelated (if they pass through different points). In this situation, the one step RSB approximation gives exact results to this problem (see also IoffeMezard (); FeigelmanIoffeMezard ()).

In the presence of non-local correlations between entries caused by the correlations in , the configurational entropy acquires a non-trivial dependence on the distance between paths. In this case one probably needs to introduce re-weighting factors in order to get the self averaging partition function, which is equivalent to the full RSB scheme instead of the one-step RSB ansatz. The development of such scheme is beyond the scope of this paper and our abilities.

Formally, the next step is averaging with respect to random on-site energies:

 Λ(E,m) = limn→01n⎡⎢⎣(K∫F(ϵ)dϵ|E−ϵ|2m)n/m−1⎤⎥⎦ (27) = 1mln(K~Im),

where

 ~Im=∫F(ϵ)dϵ|E−ϵ|2m (28)

In this equation is the box-shaped on-site energy distribution function, and the contribution of is completely neglected which can be justified at very large . Indeed, as it has been shown in AbouChacAnd (), in the limit the relevant range of disorder potential corresponds to so, the real part of the self-energy in the denominator of Eq.(22) can be neglected.

Then the increment found from Eq.(27) takes the form:

 Λ=2ln(Wc(E)W). (29)

In Eq.(29) the critical disorder close to the middle of the band is defined as:

 lnWc2=12m0lnK1−2m0. (30)

where is found from the minimization condition Eq.(25):

 2m01−2m0=lnK1−2m0 (31)

Combining (30,31) to exclude we get an equation for :

 Kln(Wc2)=Wc2e, (32)

which is exactly the "upper bound" equation for of Ref.AbouChacAnd (). At large one obtains with logarithmic accuracy

 Wc≈2eKln(eK) (33)

in agreement with AbouChacAnd ().

## Viii RSB parameter m0 and Abou-Chakra-Thouless-Anderson exponent β

The parameter found from Eq.(31) (which is independent of in the leading approximation) has a special physical meaning, it is related to the power-law dependence of the distribution function. To establish this correspondence consider the behavior of the moments in IOTL.

The value of computed in section VII describes the exponential growth of the typical value of . One can use the same method to determine that governs the growth of the moments . Repeating the arguments of section VII we get provided that is sufficiently small, , so that analytic continuation to this value of has the same structure as continuation to employed in section VII . For the solution disappears which describes the fact that higher moments of grow faster with than its typical value or even diverge indicating the absence of linear response (). This behavior of the moments implies that the distribution function of acquires a stationary form at with the power law tail with the lower cutoff and the upper cutoff that grows with . Indeed, for this distribution all moments are finite for and diverge for . The same conclusion can be obtained by solving the equation for the evolution of the distribution function directly (see FeigelmanIoffeMezard () and Appendix B).

In particular, at the Anderson transition determines the power of (or ) in the power-law distribution function in both IOTL and ATL limits (see Appendix A):

 P(ρ)∝1ρ1+m0, (34)

Thus the RSB parameter at is identical to the exponent introduced in Ref.AbouChacAnd (). The work AbouChacAnd () also shows that in agreement with the general arguments of section III the exponent at the transition is

 β=1/2 (35)

The same result, Eq.(35), follows from the duality (96),(97) for a linear .

On the other hand, one obtains from Eq.(31):

 m0≈1/2−1/(2lnK). (36)

which coincides with the exact result within the accuracy of the approximation that neglected the effects of the real part of the Green’s function. In the next section we discuss how one can take into account these effects and develop better approximation.

## Ix Minimal account for the real part of self-energy

Fairly large errors of the “upper-limit” value of (32) can be traced back to the inaccurate value of (36) at that differs from the exact result, , by This difference is due to the complete neglect of the real part of the self energy in the denominator of (22) and the resulting logarithmic divergence of the average at . In order to improve the accuracy of the analytic theory we introduce the effective distribution function of the real part of instead of the distribution of the on-site energies . Because in the resulting model the entries in (24) remain uncorrelated, for this effective distribution the one-step RSB ansatz is treated exactly and leads to (27). This will allow us to restore the exact value of at the transition point and dramatically reduce the error in the value of

We emphasize that function is an effective distribution which takes into account correlations of different entries caused by correlations of in a given path in Eq.(23) when is not neglected. Indeed, if is anomalously large, the Green’s function at the affected sites should be anomalously small. This effect leads (see Appendix D for detailed derivation) to the symmetry of the PDF of the product along a path of length , which is equivalent to the symmetry of :

 Feff(ϵ+E)=Feff(ϵ−1+E). (37)

Notice that the introduction of does not solve the problem of non-local correlations between different paths which may invalidate the one-step RSB ansatz.

The simplest approximation for the effective distribution function that is close to the original distribution but obeys the symmetry Eq.(37) at is

 Feff(ϵ)=Feff(1/ϵ)=θ(|ϵ|−2/W)θ(W/2−|ϵ|)W−4/W, (38)

Thus the minimal account of is equivalent to imposing the symmetry (37) which eliminates small . Physically, it describes the level repulsion from the state at energy . We will see below that this is a crucial step with many implications. For instance it allows for the replica-symmetric solution which corresponds to . Notice that in the absence of the gap at low implied by the distribution (38) this solution did not exist as always diverges at .

Now the critical disorder and are found from the solution of the system of equations:

 ~Im =K−1, (39a) ∂~Im∂m =0, (39b)

where

 ~Im=∫Feff(ϵ+E)dϵ|ϵ|2m.

On can easily see that the symmetry Eq.(37) results in:

 ~Im =~I1−m, (40a) ∂~Im =−∂~I1−m, (40b)

where . Eq.(40) implies that , so that is an exact solution to the equations (39b). The critical disorder is then found from the first equation

 ~I1/2(Wc)=K−1. (41)

## X Improved large-K approximation for Wc

As we have seen, the Abou-Chakra-Thouless-Anderson "upper bound" (32) for has an accuracy of . The symmetry (37) allows one to take into account all terms and obtain a new estimate for the critical disorder, which accuracy is at least .

Computing using Eq.(38) one reduces Eq.(41) to

 2Kln(Wc2)=Wc2−2Wc. (42)

The results of solution of this algebraic equation for different connectivity is summarized in Table I.

One can notice an excellent agreement with numerics even for the minimal . The results for large are expected to be even more accurate.

We conclude that the correct symmetry improves at lot the large-K approximation and leads to an extremely simple and powerful formula for which accuracy exceeds by far any approximations to the exact Abou-Chacra-Thouless-Anderson theory known so far.

## Xi Analytical results for D(w) and m(W) at the band center E=0.

The results of Section IX allows one to get the analytical results for anomalous dimension . Plugging (38) into (27) we compute the increment . Finally we use Eq.(21) to convert it into fractal dimension . The resulting prediction of the RSB theory for the behavior of at is displayed in Fig. 4. In this figure we also compare the RSB result for with the results of the population dynamics and the results of the direct numerical diagonalization for finite RRG. The latter were obtained by using the data of Ref.Our-BL () for the distribution of for RRG of moderate sizes . In more detail, we have computed the finite-size spectrum of fractal dimensions defined by

 f(α,N) =ln[NP(ln|ψenv|2)]/lnN, (43) α =−ln|ψenv|2/lnN

where is wave function envelope and is the corresponding probability density, extrapolated it to and found from the maximum point of the extrapolated . More details can be found in Appendix F. Here we only note that the spectrum of fractal dimensions translates into

 P(lnZ)=Aexp[lnNf(lnZlnN)] (44)

where . The form (44) of the distibution is a very general one that corresponds to the fractality of the physical quantity described by variable For linear function it is reduced to the power law. More generally it describes the crossover from the power law to a Gaussian-like behavior. This distribution function appear in a variety of problems, including classical onesKravtsovPekola2016 ().

We notice excellent agreement between the results obtained by population dynamics and the data of the direct diagonalization away from Anderson transition, . At larger population dynamics and direct diagonalization results deviate from each other due to the rapidly growing correlation volume as approaches As discussed in Section XIII the correlation volume at , so the deviations of the results of direct diagonalization from the infinite size limit at are only to be expected. The results of the RSB theory and population dynamics are in a very reasonable quantitative agreement with each other at .

At low the population dynamics and RSB give qualitatively different predictions. The former predicts gradual crossover to that follows the scaling behavior with whilst RSB predicts the ergodic transition at . Both approaches correspond to infinite sizes, the difference between them is due to incomplete account of the effects of the real part of the Green’s function in the analytic RSB solution. Generally, one expects that analytic solution is exact at large at which all effects of the real part of the Green’s function are small. This can be verified by computing the dependence of . Both the general arguments of section III and population dynamics (see Fig.16) predicts the distribution function at This translates into for all (section VIII). For RSB solution at but it deviates from it at . These deviations are small in a wide range of :

 m0−12≈3ln(2Kln(W2)W2−2W)2ln2(W2) (45)

confirming the accuracy of RSB approach at large . The Fig. 5 displays dependence that confirms that it stays close to in a wide range of Strong deviations of from for at exactly correspond to the deviations of the RSB and population dynamics results shown in Fig. 4.

We now discuss in more detail the predictions of the RSB solution for low where its accuracy is uncertain. As decreases below , increases monotonically from and it reaches at

 (WE/2)W2E+4W2E−4 =e√K, (46) WE ≈2e√K, K≫1 (47)

At this point the RSB solution terminates (see Fig. 6), because only are allowed in RSB solution. This proves, within one-step RSB, existence of the ergodic transition from the non-ergodic extended (multifractal) phase described by the RSB solution to the extended ergodic phase described by the replica symmetric (RS) solution with . Existence of such a RS solution and the fact that at is a consequence of the symmetry Eq.(38). Indeed, at (and ) we have:

 Λ(m=1)=ln(K∫Feff(ϵ)dϵϵ2). (48)

Because of the symmetry of , changing the variables of integration converts the integral in Eq.(48) into the normalization integral for the effective distribution function . Then we immediately obtain from Eq.(21) that the RS solution corresponds to , i.e. to the ergodic extended phase.

In Appendix C we prove that the existence of termination point of the RSB solution at a non-zero such that , , and is a generic feature of the theory. It occurs at any function obeying the symmetry Eq.(37) and decreasing sufficiently fast at large and small , e.g faster than at .

Neither population dynamics nor RSB theory precludes the ’first order’ jump in for finite graphs in which loops become important at large scales. In RSB theory gives a limit of stability of the non-ergodic extended phase. The actual ergodic transition may occur before this limit is reached, as the replica symmetric solution exists in the entire region . In this case it should be a first order transition at similar to the one observed in AltshulerCuevasIoffeKravtsov2016 (). We estimate the effect of the loops of large sizes in finite RRG in Section XVI and conclude that they might become relevant at in agreement with the transition observed AltshulerCuevasIoffeKravtsov2016 () in the data for the largest graphs.

## Xii Application to Rosenzweig-Porter model

The generalized Rosenzweig-Porter random matrix model (GRP) is probably the simplest model in which both localization and ergodic transitions happen, with the non-ergodic extended phase existing in between KravtsovKhaymovichCuevas2015 (). It plays the same role for the field of quantum non-ergodicity as the random energy model for classical spin glasses.

The model is formally defined RPort (); KravtsovKhaymovichCuevas2015 () as a Hermitian matrix with random Gaussian entries independently fluctuating about zero with the variance , and , where is an -independent number. By changing the energy scale one may define , where and . Thus the GRP model corresponds to . The AT critical point in the limit corresponds to (see Ref.KravtsovKhaymovichCuevas2015 () and references therein) and thus .

We apply Eq.(29) to GRP, as it should be valid for any graph with connectivity . The GRP can be mapped on a graph where each site is connected with any other site directly and thus in this model and the graph diameter . Thus we immediately obtain from Eq.(29):

 Λ=2ln(NNγ/2)=(2−γ)lnN. (49)

Now, Eq.(18) terminated at and (11,12) give

 ρtyp ∼ηeΛ∼ηND1 (50) D1(γ) =2−γ (51)

This result coincides with the fractal dimensions (valid for all ) for the GRP obtained in Ref.KravtsovKhaymovichCuevas2015 () from completely different arguments. Note that minimizing is in the entire region of non-ergodic extended states in the limit . This implies that the exponent in the power-law dependence (34) is for all values of , in agreement with general expectations and the results of works KravtsovKhaymovichCuevas2015 (); RP-Bir (). We conclude that the exact result of the RSB theory in this case is associated with the value of the exponent in the entire region of non-ergodic states.

## Xiii RSB results for ρtyp

The typical local density of states in Anderson thermodynamic limit, is an important characteristic of a strongly disordered system. Physically, it characterizes the inverse escape time, i.e. the time needed for a particle to leave a vicinity of a given site. This time is finite in delocalized regime but becomes infinite as . We discuss the relation between and physical properties in the end of this section.

In this section we compute as a function of disorder using the one-step RSB and compare it with the results of the population dynamics. Our goal is to obtain a stationary distribution of for RRG where all sites are statistically equivalent. Note that stationarity of the probability distribution, i.e. its independence of and by no means implies homogeneity of for a particular realization of disorder.

It will be more convenient for us to solve the equations for the imaginary part of the self-energy related to the imaginary part of the Green’s function by . As we show below close to the critical point the typical becomes exponentially small in . This strong dependence on is the same for and that differ from each other only by factors and . Below we shall focus on the strong exponential dependence of these quantities and ignore the difference between them.

The recursion equation for follows directly from (17):

 S(ℓ+1)i=∑j(i)S(ℓ)jϵ2j+(S(ℓ)j)2, (52)

where , and is a sum of independent random on the ancestors sites.

The power law distribution of , which is a general property one-step RSB solution (see Appendices B,E), implies that the contribution to the moment comes from a wide region of . Indeed, for this and only this moment the integral is logarithmically divergent. The wide distribution of individual terms in the sum (52) implies that in this sum one term is much larger than others, so that the -th power of the sum is equal to sum of the powers. This allows us to write the closed equation for :

 ⟨ρm0⟩=K⟨ρm0(ϵ2+ρ2)m0⟩ϵ, (53)

where denotes averaging with the distribution function approximated by Eq.(38).

Eq. (53) can be rewritten in terms of the distribution function

 ⟨ρm⟩=K∫dρρmP0(ρ)Ξ(ρ;W,m), (54)

where

 Ξ(ρ;W,m)=∫dϵ(ϵ2+ρ2)mFeff(ϵ). (55)

The averaging over and in the same generation are independent, because on the tree depends only on in the previous generations.

We now use the definition valid for any moment to arrive at the equation:

 ∫dρρm0P0(ρ)[Ξ(ρ;W,m0))−K−1]=0, (56)

where is taken from the solution of (25).

Equality (56) is an implicit equation for . In order to make it explicit we note that is the stationary distribution function at In the spirit of Ginzburg-Landau theory we assume that this function does not change significantly when appears below at large so that

 P0(ρ)∝1/ρ1+m0(W),ρ≳ρtyp (57)

One can also show that in the vicinity of the Anderson transition point the curvature of as a function of is very small: decreases faster than as