A random matrix model with localization and ergodic transitions

Motivated by the problem of Many-Body Localization and the recent numerical results for the level and eigenfunction statistics on the random regular graphs, a generalization of the Rosenzweig-Porter random matrix model is suggested that possesses two localization transitions as the parameter $\gamma$ of the model varies from 0 to $\infty$. One of them is the Anderson transition from the localized to the extended states that happens at $\gamma=2$. The other one at $\gamma=1$ is the transition from the extended non-ergodic (multifractal) states to the extended ergodic states similar to the eigenstates of the Gaussian Orthogonal Ensemble. We computed the two-level spectral correlation function, the spectrum of multifractality $f(\alpha)$ and the wave function overlap which all show the transitions at $\gamma=1$ and $\gamma=2$.


Introduction
Motivated by the problem of many-body (MB) Localization [1] and the applicability of the Boltzmann's statistics in interacting disordered media [2], there was recently a revival of interest to the Anderson localization (AL) problem on hierarchical lattices such as the Bethe lattice (BL) or the random regular graph (RRG). Due to hierarchical structure of the Fock space connected by the two-body interaction, statistics of random wave functions in such models is an important playground for MB localization. In particular, the non-ergodic extended phase on disordered hierarchical lattices could model a breakdown of conventional Boltzmann statistics in interacting MB systems and an emergence of a phase of a 'bad metal' [3] or unconventional fluid phases [4] in systems of interacting particles.
However, even for the one-particle AL existence of such a phase in a finite interval of disorder strengths is a highly non-trivial issue.
According to earlier studies [5,6] there is only one transition in such models at a disorder strength W W AT = which is the AL transition (AT) that separates the localized and ergodic extended states. However, recent numerical studies [7] of level statistics on RRG seem to indicate on the second transition at W W W ET AT = < which is identified as the transition between the ergodic and non-ergodic extended states (ET). Subsequent studies [8,9] raise doubts about the existence of the second transition on RRG. Numerical results of [8,9] indicate on the non-ergodic states on RRG in a wide range of disorder strengths down to very low disorder W W 5 17.5 AT = »  , while in [9] it is demonstrated how an apparent non-ergodic behavior for the intermediate matrix sizes N in Levy random matrix (RM) ensemble evolves into the ergodic one at larger Nʼs. Complexity of RRG and the controversy associated with existence of the ergodic transition at W W ET = necessitate a search for a simpler model in which such a transition may occur.
Inspired by the success of Wigner-Dyson RM theory [10] which predictions are relevant in such seemingly different fields of physics as nuclear physics and nano-and mesoscipic physics, our goal is to search for a RM model that would be able to give a simple and universal description of all the three phases: good metal, MB insulator and 'bad metal', which are relevant in the problem of MB localization. An important heuristic argument to construct such a model is that RRG with disordered on-site energies i e is essentially a two-step disorder ensemble. The disorder of the first level is the structural disorder due to the random structure of RRG where each of N sites of the graph is connected with the fixed number K 1 + of other sites in a random manner.
An ensemble of tight-binding models on such graphs with deterministic on-site energies i e and hopping integrals is believed to be equivalent to the Gaussian RM ensemble [11]. The disorder of the second level is produced by randomization of i e fluctuating independently around zero with the distribution function p ( ) e .
For numerical calculations this distribution is often taken in the form p O N 0 ( ) l~and γ is the main control parameter of the problem.
One can estimate the strength of disorder required for the AL transition as corresponding to the typical fluctuation of diagonal matrix element equal to the typical off-diagonal matrix element times the coordination number K. For the coordination number K Ñ (each site is connected with any other one) this results in N 1 s~, or N 1 2 s~. However, this estimation does not take into account a random, sign-alternating character of the off-diagonal matrix elements. It is likely that for sign-alternating hopping there is another relevant coordination number Ñ with the critical scaling N 1 s~. As we show below it corresponds to the ergodic transition. For technical reasons the most significant progress in the analytical studies of the model was achieved [13,14,16] for the 'unitary' RP (URP) ensemble. The conclusion was that at 2 g = the spectral formfactor (two-level correlation function) is neither of the Wigner-Dyson nor of the Poisson form [13,15,16] which is typical for the AL transition point. In contrast, at 1 g = the level statistics was found to be GUE [14].
The papers [13,14,16] have a status of classic keynote papers in the field.
The Dyson ideas of the Brownian motion of energy levels first applied to the RP ensemble in [13] were developed in the series of works [17,18]. It was shown that the possible transitions in the level statistics are associated with the fixed points of parameter We conclude that the only fixed point of Λ is possible at 2 g = , and no transition at 1 g = can be obtained from the results of [17,18].
In this paper by a more sophisticated analysis of the two-level correlations and the eigenfunctions statistics we show that the above extension of the RP model indeed contains not one but two transitions. One of them at 2 g = corresponds to the transition from the extended to the localized states. However, the extended states emerging at 2 g < are not ergodic: their support set contains infinitely many N D 1 sites in the N  ¥ limit, which, however, is a zero fraction of all sites, since D 1 1 < . Such non-ergodic extended states on RRG are recently discussed in [8]. With further decrease of γ the second transition at 1 g = happens which is a transition from the non-ergodic extended states to the ergodic extended states with D 1 1 = similar to the eigenstates of the GOE.
We prove this statement in three steps. As the first step we use the perturbative arguments to compute the statistics of wave function amplitude r o 2 | ( )| y in a certain observation point r o . We obtain a drastic change of the character of this distribution at 1 g = and 2 g = which is summarized in figure 1. This result is fully confirmed by a numerical diagonalization of the Hamiltonian (see figures 2-4). It is also confirmed by the numerical   analysis of the moments of random wave functions which determine their Shannon entropy and the support set dimension D 1 (see figure 5). Then we compute numerically the overlap of amplitudes for two different wave functions with the energy difference ω and find the scaling with N of the Thouless energy E N z Th~-which exponent z changes abruptly at 1 g = and 2 g = (see figure 6). Finally, we perform a rigorous calculation of the spectral form-factor which also shows the transition at 1 g = and 2 g = (see figure 7). In the last section we compare the corresponding results for our model and for the RRG and demonstrate their similarity. It allows us to unify both models in a special universality class of random hierarchical models which differs from the one realized in localization transition points of two-and three-dimensional Anderson models. Further details concerning this model can be found in [20].

Statistics of eigenfunction amplitudes
As the off-diagonal matrix elements in equation (1) are small, one can employ the perturbation theory for computing the distribution function of the amplitudes x N r o 2 | ( )| y = . The first order perturbation theory gives: where the maximum of r n ( ) y is supposed to be at r r n = .
. Although it is hard to prove such a convergence rigorously, a plausible argument in its favor is that the effective coordination number of oscillatory contributions is N rather than N. The corresponding criterion of convergence is N N W Consider the regular part of the characteristic function Q e N r That is why it is only the expansion of equation (3) at small ζ what matters for P(x) at 1 g > . Thus we obtain for the regular part of the function P(x): There are two normalization conditions for P(x): the normalization of probability equation (5) and the normalization of the wave function equation (6): A caution, however, should be taken: by normalization r 1 on any lattice site cannot exceed 1, and therefore x N  . One can see that the above estimation for x max is valid only for In order to compensate for the deficiency of normalization in equation (6) one has to assume a singular part of One can see that for 2 g > equation (6) is dominated by the singular term, and A N 1 = -. This corresponds to the strongly localized wave functions. The mechanism of emergence of the singular term at the AL transition at 2 g = is somewhat similar to the Bose-condensation, where the singular term also appears because of the deficiency of normalization of the Bose-Einstein distribution.
One can express the distribution function equation (4) through the spectrum of fractal dimensions [8,21]: Using equations (4)-(7) one obtains: The upper cutoff max a g = corresponds to the lower cutoff x min . The lower cutoff min a depends on γ. In the localized region 2 g > , figure 1(a), 0 min a = . At the AL transition point 2 g = the function f ( ) a has the same triangular shape as at W W AT = on RRG, figure 1(b). In the region of the extended non-ergodic states 1 2 g < < , figure 1(c), 2 0 min a g = -> . It is remarkable that in the entire region 1 2   g the symmetry [21,22] Finally, at 1 g = the two limits min a and max a collapse in one point 1 a = which marks the transition point to the ergodic state, figure 1(d).
Note that f ( ) a for 2 g > (see figure 1(a)) has a singular peak at 0 a = which corresponds to the singular . This singular f ( ) a is not a limit of any convex function. However, one may easily see that all the moments N N q 2 q | | y á ñ~thave the same exponents q t as for the 'convex' f ( ) a shown by the dashed line in figure 1(a): f ( ) a a g = for 0 a g < < . Such a triangular f ( ) a with the slope smaller than 1/2 also holds in the localized phase on RRG [8].
The numerical calculation of f ( ) a which involve the rectification and extrapolation procedures described in [8], fully confirms the above results. In figure 2 we present the results of this calculation for N=2 8 -2 14 and the extrapolated f ( ) a (shown by a solid red line) for 3 g = which perfectly coincides with the prediction of our perturbative analysis above. The similar coincidence was obtained for 1.5 g = , while for 0.75 g = the distribution function P(x) is practically indistinguishable from the Porter-Thomas distribution of the GOE.
3. The support set dimension D 1 By calculating the Legendre transform q t [21] of f ( ) a (shown in figure 1) one finds that in the intermediate for q 1 2 > are the same and equal to: Thus the support set of a typical wave function is a fractal containing N N D 2 1 = g sites. As N 2  ¥ g in the limit N  ¥ it is an extended state. However the support set contains a fraction of all sites F N 1 = g tending to zero in this limit. Thus it is a non-ergodic state. In the localized phase 2 g > (including the critical point 2 g = ) we obtained: One can see that the fractal dimensions D 0 q = only at q 1 g > , while they are non-zero and negative for q 0 1 g < < . This is not exactly the behavior of the typical Anderson insulator where all fractal dimensions with q 0 > are equal to zero. The behavior similar to equation (10) were found in certain two-dimensional random Dirac models [25][26][27][28]. Such a quasi-localized phase is referred to as the frozen phase and the corresponding transition is known as the freezing transition. In such a phase a typical wave function amplitude has several sharp peaks separated by valleys where  figure 5. However, for 1.5 g < at our systems sizes the support set dimension D 1 approaches the values expected in N  ¥ limit (dashed line in figure 5). This fact implies that for 1.5 g < we reached N N c  and thus it may serve as numerical evidence of convergence and existence of non-ergodic extended phase in the thermodynamic limit.
We also introduced the N 1 corrections to the fit with its magnitude C N 1 being a measure of the global curvature of the x x ln á ñ versus N ln dependence. Remarkably, C N 1 changes sign at both the transition points 1 g = and 2 g = (though the positive C N 1 is very small for 2 g > ). We also checked that it changes sign at the localization transition point of the 3D Anderson model (not shown). We believe that the changing of sign of C N 1 is a convenient way to identify the points of both localization and ergodic transitions.

Overlap of different wave functions
Next we compute numerically the overlap of different wave functions K N r r . Note that for the ergodic wave functions of GOE K 1 ( ) w = is independent of the energy difference ω, as in this case the overlap is always 100%. Our results presented in figure 6 show that for 1 g > the overlap K ( ) w has a plateau at E Th w < which is followed by a fast decrease K . The Thouless energy E Th [29] that separates the GOE-like behavior (plateau) from the system specific behavior (K 2 ( ) w w µ -), depends on N as a power-law N z -. However, the scaling exponent z is different in all the three phases (see figure 6). In the localized phase ( ) , and the plateau extends to entire spectral band-width. The emergence of such a plateau that survives the limit N  ¥ is a signature of the ergodic state [30]. Surprisingly, the overlap function K N 0 as N  ¥ at any fixed ω and 1 g > . This phenomenon of 'repulsion of wave functions' [30] is a peculiar feature of our model. The non-ergodic fractal states at the localization transition points of the two and three-dimensional Anderson models of the Dyson symmetry classes, as well as those of the power-law banded random matrices [21,30,31] show a different behavior. In these models, the Thouless energy for fractal states is proportional to the mean level spacing E Th d and the behavior for E N Th ( ) w > is described by the conventional Chalker's scaling [32,34,35]: Only at a very large energy separations ω of the order of the total spectral band-width, the 'repulsion of wave functions' was observed [30].
A remarkable feature of the present model is that the Thouless energy in the region of extended non-ergodic states is much larger than the mean level spacing: One can interpret this relationship as a non-trivial dynamical scaling exponent For non-interacting systems in two or three-dimensions in the point of Anderson transition z=1 for all Dyson symmetry classes. A non-trivial z is known only in two-dimensional systems described by the Dirac equation with random vector-potential which belong to chiral symmetry classes [25][26][27][28] where the freezing transition is observed.
In terms of the dynamical exponent z the leading power-law term in the Chalker's scaling for E Th w  can be rewritten as [28]: In our model we have: One can consider equation (15) as a particular case of expansion in E 1 Th w  with the leading term exponent 1 m = : in which the coefficient c 0 is zero. We will see below that equations (12) and (16) hold for the RRG too. However in this case the coefficient c 1 0~i s non-zero. Thus one can speak of the special universality class of models with a non-trivial z 1 ¹ and 1 m = which the present model belongs to together with the RRG model.

Spectral form-factor
Finally, we present the results of a rigorous calculation of the spectral form-factor C t t , e . It follows from equation (17) that for 2 g > , R 1 ( )  t in the N  ¥ limit, which corresponds to completely uncorrelated energy levels and the exact Poisson statistics. Another important feature of equation (17) is that R 0 1 ( ) = for all 1 g > .
In figure 7 we plot the unfolded spectral form-factor R ( ) t for the two phases: (i) the critical phase of the AL transition at 2 g = and (ii) the intermediate phase at 1 2 g < < . One can see that while at 2 g = the function R ( ) t has a non-trivial N  ¥ limit, for 1 2 g < < the limit coincides with that of the GUE, except for the where there is a jump in R ( ) t . This jump is a hallmark of the intermediate phase.
To demonstrate this more clearly we blow up the region of small t by re-scaling the variable corresponds to t t | | -¢ smaller than the inverse total spectral band-width. For 1 g = we have N 1 Thus the value S 0 e 2 2 ( )~p -L at 1 g = is smaller than 1. So, in addition to the specific critical behavior of S(u) at 2 g = (shown by the red curve in figure 7) one obtains yet another critical behavior of S(u) at the ergodic transition 1 g = (shown by the dashed yellow line in figure 7) which is stable in the N  ¥ limit and is characterized by S 0 1 This is how the GUE limit S u 0 ( ) º is reached. d~gappears for n at 1 2 g < < : The level compressibility χ [33,34] is ill-defined in our model because of the jump in the spectral form factor Formally it can take any value from 0 c = to 1 c = depending on the parameter n N 2 g -.
This is in contrast to other RM models with multifractal eigenstates, e.g. PLBRM and Moshe-Neuberger-Shapiro models (see [34] and references therein) where the level compressibility is well defined and takes a definite value 0 1 c < < in the N  ¥ limit.

Comparison with RRG model
In the introduction we mentioned a heuristic relation between the Anderson model on a hierarchical RRG and the RP model which has no apparent hierarchical structure. It is instructive now to compare the main results of this paper with the corresponding results for RRG.
First of all we recall (see figure 1 and the corresponding explanations in the text) that all the moments q 2 | | y ñ in the localized and the AT critical phases of our model have exactly the same q-dependence as in the corresponding phases of the RRG. The N -and g -dependence of the moments is also very similar (see figure 5 with figure 8) to the corresponding N and W-dependences for the random wave functions obtained by the exact diagonalization of the Anderson model on RRG with the branching number K 2 and = N 2000 16000 = -. However, there is an important difference: we found only one point of changing the sign of C N 1 on RRG which corresponds to the known point of the AL transition at W » 17.5.
In figure 9 we demonstrate that in the case of RRG the scaling of the Thouless energy with the system size follows the same equation (12)

= -
(see inset of figure 9). This is very different from the usual scaling K N F N N N ,  . Our conclusion is that the localized and AT critical states are very similar for our model and RRG. The extended states show non-ergodicity in a broad interval of W and g and are characterized by the Thouless energy which in both models is much larger than the mean level spacing. However, the existence or nonexistence of the ergodic transition is more subtle and depends on tiny features of the model. It exists in our model and most probably does not exist on RRG with the branching number K=2. Nonetheless our study largely confirms expectation on the similarity between the RRG and RP models. This allows us to speak on the special class of models with the explicit (RRG) or hidden (RP) hierarchical structure.