{'fraction_non_alphanumeric': 0.038070603664978685, 'fraction_numerical': 0.011803708696127364, 'mean_word_length': 4.32389449185415, 'pattern_counts': {'":': 0, '<': 13, '<?xml version=': 0, '>': 9, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 2, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
The renormalization flow of the hierarchical Anderson model at weak disorder

We study the flow of the renormalized model parameters obtained from a sequence of simple transformations of the 1D Anderson model with long-range hierarchical hopping. Combining numerical results with a perturbative approach for the flow equations, we identify three qualitatively different regimes at weak disorder. For a sufficiently fast decay of the hopping energy, the Cauchy distribution is the only stable fixed-point of the flow equations, whereas for sufficiently slowly decaying hopping energy the renormalized parameters flow to a delta peak fixed-point distribution. In an intermediate range of the hopping decay, both fixed-point distributions are stable and the stationary solution is determined by the initial configuration of the random parameters. We present results for the critical decay of the hopping energy separating the different regimes.


I. INTRODUCTION
The localization of a quantum particle in the presence of a random potential remains a very active topic in condensed matter physics. 1 The prototypical model of localization is the Anderson tight-binding model with short-ranged, nearest-neighbour hoppings on an hypercubic lattice. 2 Such model undergoes a transition between extended and localized wavefunctions for sufficiently high spatial dimensions. 3 In the context of random matrices, a popular class of models is represented by the Wigner ensemble 4 , where the matrix elements are Gaussian distributed random variables and the fullyconnected infinite-range character of the hopping energies prevents the wavefunctions to become localized. To study an intermediate situation between these two fundamental models, one requires a hopping energy decaying slowly as a function of the intersite distance. The hierarchical Anderson model (HAM), introduced originally by Bovier, 5 is a tight-binding model with onsite disorder and hopping energies organized in an hierarchical block structure. The hopping energy falls off as a power-law for large intersite distances, allowing to interpolate smoothly between models with short-range and infinite-range hopping energy. Hierarchical models have a long tradition in statistical physics, which goes back to Dyson,6 and they constitute an approximate route to study the behavior of models defined in terms of the standard short-range Laplacian on the hypercubic lattice, such as the classical random walk 7,8 and interacting spin systems. 9,10 Besides that, hierarchical models are conveniently designed such that they preserve their structure under renormalization transformations, 5,9,11 being amenable to an exact and thorough analysis.
Contrary to the rigorous results established for the density of states (DOS), 8,[12][13][14][15] less work has been devoted to the study of the nature of the eigenstates of the HAM. In a recent paper, 16 the authors have shown the existence of an extended phase for a sufficiently slow decay of the hopping energy, in contrast to a previous conjecture stating that all states should be localized. 17 The results are based on a renormalization procedure for the resolvent matrix, which allows to compute numerically the inverse participation ratio (IPR) for extremely large system sizes. In the present work we address the problem of identifying and analyzing the qualitative change in the fixed-point distribution of the flow equations corresponding to the onset of a localization transition in the HAM. The stationary solution of the flow equations has been studied so far only in the strong disorder regime, 17 where the Cauchy distribution is the only fixed-point and all eigenstates are localized.
Here we complement the work presented in Ref. [16] by studying the flow equations for the renormalized random potentials (RRP) resulting from the consecutive elimination of the degrees of freedom of the resolvent matrix via a simple change of integration variables. We focus on the stability of the fixed-point distribution of the RRP at the band edge of the pure spectrum, when a small amount of on-site disorder is introduced. The motivation for studying this specific situation is twofold. First, the problem concerns the survival of the band edge extended wavefunction of the pure model in the presence of weak disorder. Second, this is the energy range where the integrated DOS of the HAM has a similar behaviour as that exhibited by short-range systems in finite dimensions, and we expect that our work provides further insights on the behavior of the latter.
By means of the numerical solution of the flow equations, combined with a perturbative approach, we show that three distinct regimes emerge at weak disorder. The distribution of the RRP flows to a Cauchy distribution fixed-point provided the hopping energy decays sufficiently fast. For a sufficiently slow decay of the hopping energy, the fluctuations of the RRP vanish exponentially and the flow converges to a delta peak distribution. In an intermediate region of the hopping energy decay, our numerical results suggest that both fixed-point distributions are stable and the asymptotic behaviour depends on the specific microscopic configuration of the on-site disorder. We quantify the basin of attraction of both solutions by computing numerically the fraction of the flow that has evolved to a delta peak distribution. As it will be explained later, our conclusions are valid for finite, but very large system sizes.
The paper is organized as follows. We define the model in section II. The renormalization procedure used in deriving the flow equations is explained in section III, while a perturbative expansion of these equations is presented in section IV. The discussion of the numerical results, guided by the outcome of the perturbative approach, is left for section V. Finally, we present some final remarks in the last section.

II. THE HIERARCHICAL ANDERSON MODEL
Tight-binding models constitute the simplest lattice models to study the diffusion of a quantum particle in the presence of a spatially random potential. 2 We consider a one-dimensional chain of unity lattice spacing composed of L = 2 N sites i = 1, . . . , L, with random potentials {ε i } i=1,...,L drawn from a distribution p(ε). At this stage there is no need to specify p(ε) and we keep the model definitions as general as possible. In the hierarchical Anderson model the kinetic energy is given in terms of a hierarchical Laplacian. 5 Inspired by the original work of Dyson, 6 we define the Hamiltonian as follows where |i is the canonical site basis. The hierarchy of hopping energies has a total number of N levels, where p = 1 and p = N denote, respectively, the lowest and the highest level of the hierarchy. At each level the system is divided into 2 N −p distinct blocks, each of which contains 2 p sites. The hopping between any two sites within a single block of level p has energy JV p , while the hopping between sites in two different blocks is determined by levels higher in the hierarchy and has energy t p = J N n=p V n , where J sets the scale of energy. A schematic representation of this hierarchical block structure of the kinetic energy is presented in Ref. [16].
Distinctly from the case of ultrametric random matrices, 18 where the hierarchical structure is encoded in the choice of variances for the Gaussian distributed hoppings between the sites, here {V p } p=1,...,N are nonrandom parameters. We choose them to decay as a function of the level index according to V p = 2 −α(p−1) , where α > 1 controls the speed of the decay. This restriction on α ensures that, in the absence of disorder, the support of the DOS is bounded for L → ∞ (see Eq. (4)). For N ≫ 1, the magnitude of the hopping energy between two sites separated by a distance of O(L) scales as O(1/L α ), exhibiting the same long-distance behaviour as a tight-binding model with size L and hopping energy decaying as a power α of the intersite distance. [19][20][21][22][23] For p(ε) = δ(ε), the eigenvalues and eigenvectors of the Hamiltonian (1) can be computed analytically 5,17 and the average DOS reads where The average DOS is a series of Dirac delta peaks, which may be interpreted as arising from flat bands. Each peak in ρ pure (E) corresponds to a level of the hierarchy and the factor 2 −p comes from the degeneracy induced by the symmetry between the blocks at each level. The delta peaks accumulate at the upper spectral edge where α > 1 ensures that E pure ∞ = lim p→∞ E pure p < ∞. The IPR of a normalized eigenstate |ψ is defined as In the pure model the IPR of the eigenstate at E = E pure ∞ scales as I = 1/L, corresponding to an extended wavefunction.
The integrated density of states of the pure HAM is 8 where Therefore, close to the upper spectral edge E pure ∞ , the integrated DOS exhibits the asymptotic behaviour 8,12,14,15 The number d s is the spectral dimension 12,14,15 and its definition is motivated by noting that the same band edge asymptotics of the integrated DOS is observed in the case of the short-range Laplacian on an hypercubic lattice, for which the spectral and the spatial dimension coincide. The integrated DOS of the HAM also presents the same band edge asymptotics as that found in the pure one-dimensional tight-binding model with power-law decaying hopping energy, with an exponent in the range 1 < α < 2. 24 Therefore, the integrated DOS and the IPR of models with long-range hopping energies exhibit, in the neighbourhood of the upper spectral edge, the same behaviour as that found in the short-range Laplacian on spatial dimension D, as long as α is chosen such that d s = D. 8 Consistent with that, the HAM undergoes a localization transition close to E pure ∞ , 16 and this is the interesting region to study the flow of the renormalized parameters.

III. THE RENORMALIZATION FLOW EQUATIONS
In this section we discuss the main ideas involved in the derivation of the equations describing the flow of the RRP. The central object of our approach is the resolvent matrix The resolvent elements in the site basis can be represented in terms of Gaussian integrals according to where We have introduced the shorthand notation x 1,...,A ≡ x 1 , . . . , x A to represent sets of variables. The local parameters include the random potentials, while W (N ) encodes the hierarchical hopping contribution The essential idea consists in obtaining a recursion relation between the resolvent of a system with 2 N sites and the resolvent of a system with 2 N −1 sites, with renormalized model parameters. The change of integration variables (7) allows to calculate explicitly the integrals over ..,2 N −1 , halving the number of degrees of freedom. The function S (N −1) , following from this integration, has the same formal structure as Eq. (8), reflecting the invariance of the Hamiltonian under a renormalization transformation. After applying this change of variables ℓ times in a consecutive way, we obtain an expression for the resolvent elements G which is formally the same as Eq. (7), but depends on the function The renormalized parameters fulfill the recurrence equations 17 where i = 1, . . . , 2 N −ℓ and ℓ = 1, . . . , N . The initial values {µ ..,2 N and J (0) are the parameters of the resolvent in the original model. Equations (10) and (11) hold for a single realization of the random Hamiltonian H N with a finite size L = 2 N and random potentials drawn from an arbitrary distribution p(ε).
This procedure further provides a set of recursion relations for the resolvent matrix elements. After performing ℓ = N changes of integration variables, we end up with a single site resolvent characterized by the renormalized parameter µ (N ) 1 . This is the initial condition for the iteration of the resolvent recurrence equations from ℓ = N to ℓ = 1, which finally restores {G ii } obtained from this procedure, we refer the reader to Ref. [16].
Here we study the flow of the distribution P (ℓ) (µ) of the random variables {µ i } are interpreted as renormalized random potentials, the distinction between localized and extended states should be accompanied by a qualitative change of the fixed-point distribution P (∞) (µ) in the limit η → 0. Throughout the rest of the paper we work directly at η = 0, such that {µ (ℓ) i } are real variables and the statistical properties derived from Eq. (10) are valid for a finite system size L. In spite of that, we will be interested in the behaviour of P (∞) (µ) when L becomes very large, which eventually leads to strong fluctuations of {µ (ℓ) i } due to the presence of arbitrarily small denominators in Eq. (10). These unbounded fluctuations are suppressed by any nonzero value of η, affecting the stability of the different fixed-point distributions in a decisive way. An analogous approach has been used in the context of Levy random matrices, 25,26 where the resolvent matrix elements are calculated directly at η = 0.
The random potentials {ε i } enter solely in the initial distribution P (0) (µ) and they constitute the unique source of randomness in the flow of P (ℓ) (µ). We expect that a fixed-point distribution P (∞) (µ) is attained for finite values of ℓ, provided L is sufficiently large.
The distribution P (ℓ) (µ) can be computed analytically in two limiting situations. In the pure model, where p(ε) = δ(ε), it is easy to show that By setting E = E pure ∞ and taking the limit ℓ → ∞ we obtain the fixed-point distribution P The second solvable case is represented by a Cauchy distribution characterized by a scale parameter γ > 0 and a divergent variance. In this case, one has 17 Setting once again E = E pure ∞ , in the ℓ → ∞ limit we obtain the fixed-point distribution .
The convergence towards the stationary solution P (∞) p (µ) is naturally interpreted as a signature of the extended phase, since the RRP do not fluctuate from site to site. Besides that, the band edge wavefunction corresponding to P (∞) p (µ) uniformly spreads throughout the whole system. On the other hand, the strong fluctuations of the RRP, due to the Cauchy distribution P (∞) c (µ) and its divergent variance, are characteristic of the localized phase. We point out that spectral localization has been proven in the whole range of parameters when {ε i } i=1,...,L are Cauchy distributed random variables. 8 The study of the pure model or of the initial Cauchy distribution is less interesting, since the extended or localized fixed-points are stable in the whole parameter space, depending whether we choose a distribution p(ε) with zero or infinite variance, respectively. The choice of a distribution p(ε) with a finite variance will eventually lead to a competition for stability among P (∞) c (µ) and P (∞) p (µ). For a distribution p(ε) with a finite variance, spectral localization has been proven for α > 3/2, 12 while numerical results for the average IPR support the presence of extended states in the same range of α. 16

IV. WEAK DISORDER EXPANSION
In order to perform an expansion of Eq. (10) in powers of the disorder strength W , we rescale the random potentials as ε i → W ε i and assume that they are drawn from a distribution with ε i ε = 0 and ε i ε j ε = δ ij . We assume that W/J ≪ 1 and, for the initial iteration steps, we expand Eq. (10) up to order O(W 2 ), from which we derive the following expression for arbitrary ℓ From Eq. (14) one can compute the average and the standard deviation in which we have retained terms up to O(W 2 ). By calculating (µ The behavior of m ∞ (E) determines whether the perturbative expansion is convergent or not. One immediately notes that m ∞ (E) diverges whenever we choose E at one of the energies of the pure spectrum. This situation is trivial in the sense that the eigenstates at E = E pure p (p < ∞) are localized for arbitrary weak disorder. 16 The extended eigenstate which may remain stable for W > 0 is located at E = E pure ∞ . In this case, the behaviour of m ℓ (E pure ∞ ) for ℓ → ∞ depends on α according to .
From the perturbation expansion it follows that the delta peak fixed-point distribution becomes unstable for arbitrary weak disorder as long as α ≥ 2.
Up to now we have been disregarding the conditions of validity of the perturbative approach. Let us have a closer look on this issue by making the following change of variables which allows us to rewrite Eq. (10) as follows .
(19) From Eq. (19) it is more straightforward to understand why perturbation might fail. The expansion of Eq. (19) up to order O(W 2 ) is a good approximation throughout the whole renormalization flux provided that |ν | for a certain i and ℓ, a small denominator arises in Eq. (19), and the approximation given by Eq. (14) breaks down for ν (ℓ+1) i . This resonance-like effect yields RRP with anomalous large magnitudes and we expect that the variance of their distribution will exhibit an abrupt increase.
Although the failure of the perturbative results depends crucially on the fluctuations of the RRP, we can estimate the value of ℓ at which the perturbation breaks down for E = E pure ∞ . Let us assume that W/J ≪ 1 and the flow evolves according to perturbation in the first iteration steps, since E pure For 1 < α < 3 2 , there is no positive value of ℓ * which solves Eq. (20), since the right hand side diverges as a function of ℓ * and W/J ≪ 1 by construction. For α > 3 2 , there is always a value of ℓ * for which Eq. (20) is fulfilled, since the right hand side vanishes exponentially for increasing ℓ * . This value is given by For fixed W/J ≪ 1, ℓ * → ∞ as α approaches 3/2 from above. Equation (21) predicts that the perturbation expansion does not break down for α < 3/2, such that the delta peak distribution is the only stationary solution.

V. NUMERICAL RESULTS
In this section we discuss the numerical results for the evolution of P (ℓ) (µ) in connection with the perturbative approach of the previous section. The fact that Eq. (10) is defined for a finite system size L represents a serious numerical restriction, since the total number of iteration steps is limited by L. In order to overcome this issue, a different route is followed in the numerical calculation of P (ℓ) (µ). The RRP µ (ℓ) 1 , . . . , µ (ℓ) 2 N −ℓ , at a given layer ℓ, are statistically independent random variables and their distribution P (ℓ) (µ) depends only upon P (ℓ−1) (µ). This allows us to implement a population dynamics approach, which consists in parametrizing the distribution P (ℓ) (µ) by a large number N of stochastic variables representing instances of µ. To update P (ℓ) (µ), we choose at random two variables from the pool representing P (ℓ−1) (µ), which are used to update, according to Eq. (10), a single variable extracted at random from the pool of layer ℓ. This updating rule is repeated until P (ℓ) (µ) reaches a stationary form. One expects that the statistical properties of the RRP converge to a well-defined limit for large enough N . We remark that sample to sample fluctuations may arise in the population dynamics algorithm due to finite values of N . In this sense, the population size N plays an analogous role as L in finite size calculations of Eq. (10). For detailed discussions of the population dynamics algorithm in the context of spin-glasses and random matrices, we refer the reader to Refs. [27] and [28], respectively.
In all numerical results presented in this section, p(ε) is a Gaussian distribution with mean zero and standard deviation W . The initial values of µ  from the population dynamics method, while the solid lines are Gaussian distributions with mean zero and standard deviations for different values of ℓ, given by Eq. (17). As can be seen, the agreement between the numerical and the perturbation results is excellent for this value of α, where the delta peak is the only stable fixed-point distribution. Figure 1 illustrates the typical flow in the extended phase: the initial Gaussian distribution P (0) (µ) shrinks exponentially to a delta peak, characterizing the absence of fluctuations and the spatial homogeneity of the RRP.
In figure 2 we show the flow of P (ℓ) (µ) for W = 10 −2 , α = 2.25 and relatively large values of ℓ. The perturbative approach breaks down and P (ℓ) (µ) evolves to a Cauchy fixed-point distribution, as can be noticed from the comparison between the population dynamics data (symbols) and a Cauchy distribution obtained from a fitting of the data for ℓ = 15 (solid line). This is the only stable fixed-point distribution for this choice of α, E and W . The presence of large, scale-free fluctuations in the RRP typically yields localized eigenstates.
In order to clarify the breaking mechanism of the perturbative approach, figure 3 exhibits the standard deviation ∆ ℓ of P (ℓ) (µ) for α = 1.25 and α = 2.25, corresponding to the data in figures 1 and 2, respectively. For α = 1.25, ∆ ℓ vanishes exponentially as a function of ℓ according to Eq. (17). For α = 2.25, the flow of ∆ ℓ is described by Eq. (17) up to a certain ℓ, at which the presence of small denominators in Eq. (10) leads to an abrupt increase of ∆ ℓ by many orders of magnitude. This mechanism is responsible for the emergence of strong fluctuations in the RRP, driving the system to the Cauchy fixed-point distribution. In fact, the erratic behaviour of ∆ ℓ for α = 2.25 and large values of ℓ is a signature that P (ℓ) (µ) has evolved to a Cauchy distribution.
For further larger values of ℓ, we eventually found that the Cauchy distribution usually becomes unstable and the parameters {µ (ℓ) i } flow back to a Gaussian distribution, until they finally reach the delta peak distribution for ℓ → ∞. This effect is clearly illustrated in figure 4, where we present the standard deviation of P (ℓ) (µ) up to ℓ = 50, for α = 2.25, W = 10 −2 and E = E pure ∞ . For intermediate values of ℓ the standard deviation exhibits the erratic behaviour typical of the regime where P (ℓ) (µ) evolves to a Cauchy fixed-point distribution. However, for ℓ ≥ ℓ c the standard deviation presents once more the decay ∆ ℓ ∝ 2 − ℓ 2 , reflecting the Gaussian behaviour of P (ℓ) (µ). In order to probe the effect of the population size N on the stability of the Cauchy fixed-point, we have computed the average of ℓ c over a certain number of independent runs of the population dynamics algorithm. The outcome for α = 2.25, as a function of N , is displayed in figure 5. The data show that the mean value ℓ c diverges as a logarithmic function of N , strongly indicating that the second Gaussian regime for larger ℓ is just an artifact of the finite values of N , and the Cauchy distribution is the only stable solution for N → ∞ and large values of α.
In the population dynamics method, the sample to sample fluctuations of the initial configuration µ   figure 4) as a function of the population size N for α = 2.25 and E = E pure ∞ . The distribution p(ε) has a Gaussian form, with mean zero and standard deviation W = 10 −2 . The average ℓc is computed using 50 independent runs of the population dynamics algorithm. The solid line is the best fit ℓc = a + b ln N of the data, with parameters a = 11.1(1.3) and b = 2.14(9).
proportional to the fraction of initial configurations that flow to P (∞) (µ). From a total of S independent runs of the population dynamics algorithm, let us define F S as the fraction of runs in which the standard deviation of P (ℓ) (µ) is given by Eq. (17). We also define F S , i.e., the average value of F S over different sets of samples of fixed size S. The quantity F S provides a measure of the size of the basin of attraction of the fixed-point delta peak distribution.
We have computed the average fraction F S over five independent sets, each one containing S = 50 samples. The behavior of F S as a function of N is displayed in figure 6, for W = 10 −2 , E = E pure ∞ and different values of α. As it can be seen, in the regime of large N we have that F S = 0 for α = 2, whereas F S ≃ 1 for α = 1.5. We have checked that F S → 1 for fixed α = 1.5 and decreasing W . For α = 1.6 and α = 1.7, the fraction F S approaches a value 0 < F S < 1 when N ≫ 1. The numerical results on figure 6 strongly suggest that, for a certain interval of values of α, both the delta peak and the Cauchy distribution are stable fixed-point distributions and the asymptotic behaviour depends fundamentally on the initial configuration of the RRP.

VI. FINAL REMARKS
We have studied analytically and numerically the flow of the distribution of the renormalized random potentials (RRP) in the hierarchical Anderson model (HAM), characterized by a hopping energy decaying as a power-law with exponent α. More specifically, we have focused on the stability of the fixed-point distribution of the flow at the upper spectral edge of the pure model, when a small amount of on-site disorder is added to the system. For large values of α (short-range hopping), the RRP flow to a Cauchy fixed-point distribution, independently of their initial configuration. This is consistent with the localization of all eigenstates in low-dimensional tight-binding models with short-range hoppings. 3 For small values of α (long-range hopping), the fluctuations of the RRP vanish exponentially and the flow converges to a delta peak distribution. This is somehow consistent with the Wigner ensemble of random matrices, 4 where the fully-connected infinite-range hoppings delocalize all eigenvectors. In an intermediate range of α we have found that the delta peak and the Cauchy distribution are both stable fixedpoints, and the asymptotic flow depends on the specific realization of the on-site disorder.
Although Eq. (21) implies that the perturbative approach for the flow equations breaks down for 3/2 < α < 2, we have found numerically that the RRP flow either to a delta peak or to a Cauchy distribution in this range of α, depending on the initial configuration of the random parameters. We point out that Eq. (21) has been derived under a very crude assumption, namely that all RRP are of the same order of magnitude in the initial steps of the flow, which amounts to neglect spatial fluctuations. In spite of that, numerical and analytical results seem to agree that for α < 3/2 (d s > 4) the delta peak is the only fixed-point distribution, while for α > 2 (d s < 2) the RRP always flow to the Cauchy fixed-point distribution.
Rigorous results have shown that for α > 3/2 the HAM spectrum contains solely a pure-point contribution, 14 whereas numerical results for the inverse participation ratio support the existence of extended eigenstates in the range 3/2 ≤ α 2, 16 which coincides with the regime where the delta peak and the Cauchy distribution coexist as stationary solutions of the flow equations. Overall, these results may indicate the presence of a mixed phase in the HAM, exhibiting features of localized and extended states. The study of the spatial decay of the wavefunctions and of the level-spacing distribution could provide valuable information about the physical properties in this intermediate regime of α. Analogous examples of mixed behaviour of localized and extended features have been reported in the study of Levy random matrices 25,29 and, more recently, in the Anderson model on the Bethe lattice. 30 It would be interesting to investigate whether such unusual behaviour observed in the HAM is present close to the band edge of high-dimensional tight-binding models with short-range hoppings.