{'fraction_non_alphanumeric': 0.042281463112213265, 'fraction_numerical': 0.020148791072535647, 'mean_word_length': 4.31499176276771, 'pattern_counts': {'":': 0, '<': 7, '<?xml version=': 0, '>': 4, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Logarithmic scaling of Lyapunov exponents in disordered chiral two-dimensional lattices

We analyze the scaling behavior of the two smallest Lyapunov exponents for electrons propagating on two-dimensional lattices with energies within a very narrow interval around the chiral critical point at E=0 in the presence of a perpendicular random magnetic flux. By a numerical analysis of the energy and size dependence we confirm that the two smallest Lyapunov exponents are functions of a single parameter. The latter is given by ln L/ln xi(E), which is the ratio of the logarithm of the system width L to the logarithm of the correlation length xi(E). Close to the chiral critical point and energy |E|<<E_0, we find a logarithmically divergent energy dependence lnxi(E)proporitonal to |\ln(E_0/|E|)|^{1/2}, where E_0 is a characteristic energy scale. Our data are in agreement with the theoretical prediction of M. Fabrizio and C. Castelliani [Nucl.\Phys.B 583, 542 (2000)] and resolve an inconsistency of previous numerical work.


I. INTRODUCTION
The numerical determination of transport parameters for electrons propagating in disordered two-dimensional systems with chiral symmetry still remains an important unsolved problem. The situation can be represented by a single-band tight-binding model defined on bipartite lattices subjected to purely off-diagonal disorder like a random-magnetic flux with zero mean or real random hopping terms. The latter belongs to the chiral orthogonal universality class while the former is chiral unitary. Due to the chiral symmetry, the model exhibits metallic behavior only at energy E = 0. 1 It is therefore of considerable interest to study the critical properties of the model in the vicinity of the critical point and to investigate its universality.
The traditional finite-size-scaling analysis of disorder driven metal-insulator transitions, i.e., continuous quantum phase transitions at zero temperature, is based on two assumptions. [2][3][4][5] where ξ(E) is the energy dependent correlation length.
Here, vicinity means that both the system size L and ξ(E) are already larger than any other typical length of the model and ξ > L. (ii) At the critical energy E c , the correlation length ξ(E) diverges as with a universal critical exponent ν. Relation (2) was confirmed in a multitude of numerical work on disordered systems in spatial dimension 2 ≤ d ≤ 5 and various physical symmetries. 6 As scaling variables, for example, the localization length, 7 the smallest Lyapunov exponent, 8 the two-terminal conductance, 9,10 the energy level spacings, [11][12][13] and the inverse participation ratio 14 were successfully used.
Recently, the analysis of numerical data obtained for the energy dependence of the two-terminal conductance g on a bricklayer lattice, 15 which represents a generic lattice model for graphene, has led to a power-law energy dependence of the correlation length ξ(E) ∝ |E| −ν , where the critical exponent ν is close to 1/3. Although this outcome is in agreement with previous numerical results [16][17][18][19][20] for square and hexagonal lattices, it is at variance with the Harris criterion 21 which states that ν > 1/d, where d = 2 is the euclidian dimension of the system. More importantly, all numerical data obtained to the present date do not agree with theoretical predictions, 1,22 according to which the correlation length depends logarithmically on the energy, where A is related to the longitudinal conductivity and E 0 is assumed to be of the order of the energy band width. 1 A possible explanation of this disagreement between theory and numerical experiments may be that the energies investigated in the numerical studies, down to 10 −10 so far 20 (in units of the hopping energy), are not sufficiently small in comparison to the unspecified parameter E 0 introduced in the theory. Thus, it could be that the energy interval |E| ≪ E 0 , where the scaling holds, was not reached in previous numerical studies. A second obstacle is the vanishing of the density of states ρ(E) which occurs at E = 0 for hexagonal and bricklayer lattices in the presence of random-magnetic-flux disorder. 23 This behavior persists even in strongly disordered chiral systems so that the two-terminal conductance, which nevertheless turns out to be finite ∼ e 2 /h at the Diracpoint in graphene, 15,[24][25][26] is not a suitable scaling variable for numerical studies. Therefore, it is expedient to investigate instead the smallest Lyapunov exponents, which are associated with the localization length and are not directly affected by the vanishing density of states.
In this paper, we analyze the scaling behavior of both the two-dimensional bricklayer and square-lattice model with random-magnetic-flux disorder. Using the transfermatrix method for quasi-1d systems, [3][4][5]8 we calculate the two smallest Lyapunov exponents z 1 and z 2 for energies very close to E = 0, not achieved in previous work. Lyapunov exponents are more suitable than the conductance since they are less sensitive to the energy dependence of the density of states. Another reason for using Lyapunov exponents is that the analysis of the conductance is far more time consuming, which is crucial since quadruple precision is necessary in our case when energies smaller than 10 −16 are considered. We show that in the vicinity of the critical point our numerical data for the Lyapunov exponents lead to the relation where L is the width of the system. We prove that the correlation length ξ(E) depends logarithmically on the energy (see Eq. (3)) in agreement with the predictions by Fabrizio and Castelliani. 22 Thus, our results resolve the previous discrepancy between analytical theory and numerical calculations.

II. THE MODEL AND METHOD
We study a single-band tight-binding Hamiltonian defined on a two-dimensional square lattice with nearestneighbor hopping and random-flux disorder, which is introduced by complex phase factors in the transfer terms, where c † x,y and c x,y denote creation and annihilation operators of a fermionic particle at site (x, y), respectively. For bricklayer lattices, the prime at the first sum in (5) indicates that only the transfers along every other vertical bond are included. In this way, the square lattice is transformed into a bricklayer where the coordination number is reduced to three nearest-neighbor sites. The bricklayer lattice has the same topology as the honeycomb lattice of graphene and the Hamiltonian (5) possesses the same eigenvalues ±ε i . The phases, which are chosen to be associated only with the vertical bonds in the y-direction, θ x,y;x,y+a = θ x+2a,y;x+2a,y+a − 2πe h φ x,y , are defined by the random magnetic flux φ x,y , which is uniformly distributed −f /2 ≤ φ x,y ≤ f /2 with zero mean and disorder strength 0 ≤ f /(h/e) ≤ 1. The random magnetic flux is pointing perpendicular to the 2d lattice and periodic boundary conditions are applied in the y direction. In contrast to diagonal disorder, this random flux preserves the chiral symmetry for both the square and bricklayer lattices. We fix the units of energy and length scales by the nearest-neighbor hopping energy V = 1 and the lattice constant a = 1, respectively. The disorder strength is taken to be f = 0.5 h/e for the bricklayer and f = 1.0 h/e for the square lattice.
We use the transfer matrix method 6 and collect numerical data for the two smallest Lyapunov exponents z 1 (E, L) and z 2 (E, L). For the system width L and length L x ≫ L, we calculate the transfer matrix M = Lx i M i , and extract the two smallest Lyapunov exponents. The relative uncertainty ǫ(E, L) of our data is 2 · 10 −3 for larger widths L = 192 and L = 160, and decreases down to 10 −4 for the smallest L = 8. This requires the length of the quasi-1d systems L x to be in the range ∼ 10 8 to 10 9 . Since we expect that scaling occurs only in the vicinity of the E = 0 critical point, we consider energies as small as possible, down to the point of |E| = 10 −34 , at least for L ≤ 64. This requires to perform the calculations with quadruple numerical precision.
The specific symmetry of the model provides us with an independent test of the accuracy of our data. Due to the chiral symmetry, the spectrum of Lyapunov exponents must be degenerate at the band center for all L, Deviations from E = 0 remove this degeneracy but the average value, [z 1 (E) + z 2 (E)]/2, equals to z 1 (E = 0) for small values of E.

III. BRICKLAYER LATTICE
The energy dependence of the two Lyapunov exponents is plotted in Fig. 1 for various system widths L of the bricklayer. Our data confirm that z 1 and z 2 are analytical functions of the variable 1/| ln |E||. Therefore, we approximate their energy dependence by the Taylor expansion To estimate the energy E 0 , we first minimize the expression where F (E, L) = β 0 + β 1 (ln L) α1 χ + β 2 (ln L) α2 χ 2 , with respect to parameters α, β, and E 0 . We found that X possesses a minimum when 2.2 < α 1 < 2.5, 3.8 < α 2 < 4.1, and 0.1 < E 0 . It was not possible to obtain a better estimation of the critical parameters from this procedure since small variation of E 0 can be compensated by small change of β 2 and α 2 . Instead, in a more accurate analysis, we fit our numerical data for z 1 and z 2 to the quadratic expansion (7). Fig. 2 shows the L-dependence of the coefficients c 1,2 and d 1,2 . Our data confirm the assumed logarithmic behavior of all coefficients occurring in the expansion and where α 1 and α 2 are close to the anticipated values 2 and 4, respectively. Fig. 3 shows another test of the scaling of the Lyapunov exponents. Following the conventional scaling method, 5 we re-scaled the horizontal axis for the data shown in Fig. 1 by the parameter C(L): χ → χC(L). The such obtained C(L) gives us directly the required scaling behavior as shown in the Fig. 3. The data for both z 1 and z 2 scale to one universal curve. The inset to Fig. 3 confirms the expected power-law relation C(L) ∝ (ln L) 2 .
To obtain a quantitative estimation of the energy E 0 , we repeated the scaling analysis shown in Fig. 3 for various E 0 . Although we recovered the scaling behavior similar to that shown in Fig. 1 (data not shown), the Ldependence of the parameter C(L) depends on the choice of E 0 . As shown in Fig. 4, C(L) ∝ (ln L) κ with the exponent κ decreasing when E 0 increases, converging to κ ≈ 2 for E 0 > ∼ 1. We conclude that the energy E 0 is of the order of unity in our bricklayer model.
Finally, we plot in Fig. 5 the two smallest Lyapunov exponents as a function of a single parameter (ln L) 2 / ln(E 0 /|E|) with E 0 = 0.8. To reduce the finite size corrections, we subtract from the data the values z 1 (E = 0.L) and z 2 (E = 0, L), respectively. All data collapse onto a single curve.

IV. SQUARE LATTICE
Another possibility to check for logarithmic scaling at chiral quantum critical points is the numerical analysis of a simple square lattice. [16][17][18][19][20] We calculated z 1 (E, L) and z 2 (E, L) for L even and Dirichlet boundary conditions in the transverse direction. We found that the scaling analysis is more difficult than for the bricklayer. First, the finite size effects are more pronounced (see lower inset in Fig. 6). We can eliminate them, at least partially, by subtracting the value z(E = 0, L) from z(E, L). 30 Second, in the unperturbed model the van Hove singularity, which appears at E = 0 compared with E = ±1 for the bricklayer, may spoil the scaling analysis. More importantly, following the same procedure as for the bricklayer, we found that the function X given by Eq. (9) possesses a minimum only for small E 0 ∼ 10 −4 although the energy band widths are about the same for both lattices. Since the energy E must be much smaller than E 0 , we had to restrict our analysis to energies |E| ≤ 10 −20 . Fortunately, in such a narrow energy interval, we can neglect the quadratic term in the Taylor expansion (7). As shown in Fig. 6, both z 1 and z 2 are linear functions of χ when L ≤ 32. This enables us to estimate the exponent α 1 from the analysis of the size dependence of the slope, c 1 (L) ∝ (ln L) α1 . This analysis is independent on both the choice of E 0 and finite size effects, provided that the latter do not depend on the energy.
However, the fit turns out to be rather unstable to small changes of the data ensemble. First, the interval of χ is very narrow and almost all data points are accumulated in the right part of this interval. Therefore, the resulting fit is very sensitive to the exact value of z 1 (E = 0). Second, although we calculated our data with high accuracy, Fig. 6 shows that this is still not sufficient for a perfect determination of the slope. To check the accuracy of α 1 , we tested various data ensembles and found that α 1 varies between 1.7 and 2.1. Nevertheless, our data for the square lattice are compatible with a logarithmic scaling relation.

V. CONCLUSIONS
We analyzed the scaling behavior of the two smallest Lyapunov exponents z 1 and z 2 in disordered twodimensional chiral systems defined on a bricklayer and on a square lattice. We found that both z 1 and z 2 follow a logarithmic scaling relation as considered by Sittler and Hinrichsen 27 with a correlation length proposed by Fabrizio and Castelliani. 22 According to Ref. 27, the physical origin of logarithmic scaling is associated with multifractality and local scaling invariance. To the best of our knowledge, the results presented above are the first numerical confirmation of the scaling relation (4), in which the scaling variable is given by the ratio of the logarithm of the system width L to the logarithm of the correlation length ξ, instead of the ratio L/ξ as applied usually. This scaling is accompanied by the logarithmic energy dependence of the correlation length ln ξ(E) ∝ ln(E 0 /|E|) valid for |E| ≪ E 0 . Our results also solve the contradiction between previous numerical work, which apparently did not reach the true scaling regime, and the Harris criterion.
Two methods of the scaling analysis were used. Both confirm that the logarithmic scaling is observable only for very small values of the energy close to the chiral quantum critical point at E = 0. In order to resolve the Lyapunov exponents for energies down to |E| = 10 −34 , the implementation of quadruple precision in the numerical algorithms was necessary. This probably explains why logarithmic scaling was not observed in previous numerical work. [15][16][17][18][19][20] The question arises whether the same logarithmic scaling analysis can be performed also for the two-terminal conductance. At present this seems not possible with our available computing power. In our previous work, 15 we found g(E, L) = g 0 ln(Ẽ ⋆ (L)/|E|) for |E| > E ⋆ = E ⋆ /const., but did not observe any energy and system size dependence of the ensemble averaged conductance g c ≃ 4/π e 2 /h as long as the energy remains smaller than a certain value E ⋆ ∝ L −2 . This size dependent energy interval coincides with the recently observed depression in the density of states. 23 The observation of a tiny logarithmic energy dependence of the conductance, if present at all, would require a far more accurate numerical determination of the ensemble averaged mean conductance.
As shown analytically in Ref. 22, the logarithmic energy dependence of the correlation length near E = 0 is accompanied by a divergence of the density of states ρ(E) ∝ E −1 exp[−(4A ln E 0 /E) 1/2 ]. Such a relation is, however, not found in recent numerical work on a unitary chiral lattice model, where the density of states decreases to zero when E → 0. 23 Also, a different divergency exponent of the density of states ρ(E) ∝ E −1 exp[−1/2(c| ln E/E 0 |) 2/3 ] was derived analytically for the chiral orthogonal model. 28,29 This difference shows also up in the energy dependence of the correlation length. It would be very interesting to see, whether this subtle difference can also be observed in numerical studies on a bricklayer lattice with real random hopping disorder, which belongs to the chiral orthogonal symmetry class.