Critical statistics at the mobility edge of QCD Dirac spectra
Abstract
We examine statistical fluctuation of eigenvalues from the nearedge bulk of QCD Dirac spectra above the critical temperature. For completeness we start by reviewing on the spectral property of Anderson tightbinding Hamiltonians as described by nonlinear models and random matrices, and on the scaleinvariant intermediate spectral statistics at the mobility edge. By fitting the level spacing distributions, deformed random matrix ensembles which model multifractality of the wave functions typical of the Anderson localization transition, are shown to provide an excellent effective description for such a critical statistics.
Next we carry over the above strategy for the Anderson Hamiltonians to the Dirac spectra. For the staggered Dirac operators of QCD with 2+1 flavors of dynamical quarks at the physical point and of SU(2) quenched gauge theory, we identify the precise location of the mobility edge as the scaleinvariant fixed point of the level spacing distribution. The eigenvalues around the mobility edge are shown to obey critical statistics described by the aforementioned deformed random matrix ensembles of unitary and symplectic classes. The bestfitting deformation parameter for QCD at the physical point turns out to be consistent with the Anderson Hamiltonian in the unitary class.
Finally, we propose a method of locating the mobility edge at the origin of QCD Dirac spectrum around the critical temperature, by the use of individual eigenvalue distributions of deformed chiral random matrices.
Critical statistics at the mobility edge of QCD Dirac spectra
\FullConference31st International Symposium on Lattice Field Theory  LATTICE 2013
July 29  August 3, 2013
Mainz, Germany
1 Introduction
In this article we shall discuss an intricate parallelism between the two lattice models introduced independently by two Nobel laureates, P.W. Anderson and the late K.G. Wilson, in order to account for physical systems that are apparently unrelated, at least superficially: doped insulator and strong interaction. The bridge between the two proves to be a universal framework for quantum energy levels proposed by another laureate E.P. Wigner, called random matrix theory (RMT) [1].
Historically, it is the statistical fluctuation of lattice Dirac eigenvalues at the ‘hard edge’ of the spectrum, i.e., near the origin, that has been most extensively investigated, because it allows one to directly access the effective chiral Lagrangian that governs the lowenergy, nonperturbative regime of the theory, with a help of chiral RMT [2]. For one instance, measurement of the pion decay constant via the statistical response of small Dirac eigenvalues under the imaginary chemical potential [3, 4] provides a powerful alternative to the conventional method using the temporal correlation of axialvector current and pion operators. On the other hand, the local fluctuation of QCD Dirac eigenvalues in the spectral bulk has attracted lesser attention, as it is not directly connected to the chiral Lagrangian. Halasz and Verbaarschot [5] nevertheless went on to investigate the bulk spectral correlations of the staggered Dirac spectrum from dynamical simulation. They concluded that both short and longrange correlations of Dirac eigenvalues in the scaling regime (and less surprisingly in the strong coupling regime) are perfectly described by the classical RMT of WignerDyson classes. It indicates that entire Dirac eigenstates are extended throughout the whole lattice sites just as in the Anderson Hamiltonian with subcritical randomness. This finding seems to have lead them to presume a possibility of Anderson localization in QCD at or above , as spelled out in the very final sentences of their paper which read (quote):
“A final point of interest is the fate of level correlations during the chiral phase transition. From solid state physics we know a delocalization transition is associated with a transition in the level statistics which raises the hope that such phenomena can be seen in QCD as well.”
This statement sounded rather daring,
because of an obvious and essential difference between
the Anderson tightbinding Hamitonian and the lattice Dirac operator
(beyond whether the disorder sits on or off the diagonal of the matrix):
randomness in Anderson Hamiltonian are mutually independent,
whereas stochastic gauge field variables in QCD are strongly correlated with their neighbors.
To phrase the issue more specifically:
although in their ‘ordinary phases’,
the spectral statistics of both systems allow for a description
in terms of nonlinear models (NLMs)
[6, 7] that reflect the global symmetries of operators in concern,
and the above mentioned difference between the microscopic theories could be irrelevant,
the effect of QCD temperature on the spectral model is unpredictable
as it is not restricted by any symmetry argument.
The purpose of this article, nevertheless, is to draw a definite and affirmative conclusion to the above statement,
from largescale lattice simulations.
We hope our conclusions of identifying the QCD phase transition as
Anderson localization transition at the zero virtuality,
and of establishing the presence of localized states in the hightemperature phase,
provide a novel viewpoint of the issue,
especially in the advent of quarkgluon plasma formation by the heavy ion collision.
In Sect. 2 we start from basic fact on the level statistics of Anderson Hamiltonians and the RMT. We overview the critical statistics at the mobility edge in Sect.3 and its effective description in terms of deformed RM ensembles in Sect. 4. In Sect. 5 we examine the large ensembles of staggered Dirac spectra of hightemperature QCD with 2+1 quarks at the physical point, on lattice sizes up to . Through fitting the local level statistics to the deformed RMs, we confirm the existence of the scaleinvariant mobility edge in the Dirac spectra. In Sect. 6 we propose a strategy of locating the mobility edge at the origin of Dirac spectrum by the use of individual eigenvalues.
2 Anderson Hamiltonians and random matrices
Anderson tightbinding Hamiltonian on a dimensional lattice with or without a magnetic field is defined as [8, 9, 10],
(1) 
and with spinorbit coupling as [11],
(2) 
Here are i.i.d. random variables on the lattice sites , modelling the impurities in the crystal, and the constants in the hopping terms parameterize the strength of the external magnetic field and the spinorbit coupling. The Hamiltonian matrices (1) at and (2) at satisfy realsymmetric and quaternionselfdual conditions, respectively, whereas (1) at satisfies no such (pseudo)reality condition and thus is merely complexHermitian. These three cases are said to belong to the orthogonal, symplectic, and unitary classes, and are assigned the Dyson indices . These oneparticle Hamiltonians are tailored to model the properties of energy eigenvalues and eigenstates of electrons in a disordered metal. For the site energies taken from the Gaussian distribution, one can analytically perform the ensemble averaging of the fold product of characteristic polynomials . After introducing the HubbardStratonovich field and taking the thermodynamic limit to integrate out its heavy components, one can derive a NLM of the universal form [6],
(3) 
Here is the soft component of the HubbardStratonovich field representing the embedding of the NambuGoldstone manifold associated with the symmetry class (quaternionic, complex, and realGrassmannian manifold for ) in , , the diffusion constant that depends on the randomness, and the mean level spacing of an eigenvalue window in the spectrum from which the cluster of s are taken. Correlation function of densities of states, , follows from through the replica trick . If we choose this window away from the band edges such that the eigenvalues are populated densely enough and the typical difference of eigenvalues are much smaller than the Thouless energy defined as , the path integral is dominated by its zero mode,
(4) 
This zerodimensional NLM would also have followed had we started out from the ensemble of evensimpler dense matrices distributed according to the Gaussian measure [12]. For these RM ensembles, known as GOE, GUE, GSE for , the distribution of spacings of adjacent levels normalized (unfolded) by the mean level spacing, , are analytically calculable in limit of large matrix size [13, 14]. For unitary ensembles it is generally expressed in terms of Fredholm determinant, of the integration kernel (the local asymptotic form of in (4)), over the Hilbert space of functions on an interval . Those for orthogonal and symplectic ensembles are similarly expressed in terms of matrix valued kernels, and are known to be related to the corresponding [15, 16]. In the case of bulk correlation of GUE, the pertinent kernel is the sine kernel,
(5) 
which comprises of trigonometric functions originated from (4). The asymptotics of analytically obtained as the above takes the ‘Wigner surmised’ form
(6) 
The level repulsion (for small ) and the rigidity (for large ) in (6) are consequences of extended eigenfunctions that are highly likely for a randomlygenerated dense matrix. Thus the numericallyobserved perfect agreement between the level spacing distributions (LSDs) of Anderson Hamiltonians and random matrices at relatively small randomness and in the spectral bulk (exemplified in Fig.1) is well understood, signifying the extendedness of the singleelectron wave functions in a weakly disordered metal.
3 Critical statistics at the mobility edge
The reduction to the zerodimensional NLM (4), i.e., to the universality class of RMT, might not take place for various reasons: if one simply set the eigenvalue window too close to the band edge, the typical difference of s could be of the same order or larger than the Thouless energy. Alternatively, increasing the disorder strength would diminish the diffusion constant, and accordingly the Thouless energy below the mean level spacing of the energy window in concern. In either case the condition is violated and one has to deal with the whole path integral nonperturbatively, which appeared to be an unfeasible task at the time of perception of the NLMs (3). Accordingly, a perturbative analysis using the expansion from two dimensions was applied and suggested that the functions for the dimensionless conductance (or for the symplectic case) [17],
(7) 
is likely to possess an IRunstable fixed point that separates the metallic and insulating regimes for the orthogonal and unitary classes in three or larger dimensions, and the fixed point for the symplectic class persists even in two dimensions. Once the presence of a fixed point is assumed, the energy spectrum close to the thermodynamic limit is expected to split clearly into the ‘metallic’ region in the band center the ‘insulating’ region at the band edges, depending on the mean level spacing. The energy levels should obey the statistics of random matrices and all eigenfunctions are extended in the former, whereas in the latter all eigenfunctions are localized and accordingly the energy levels should have no correlation (Poisson statistics). As the disorder strength is increased, the boundary of two regions, called mobility edge, move toward the band center and disappear alongside the extended eigenstates between two edges, leading to the metalinsulator transition [18, 19] Being an unstable boundary region separating WignerDyson and Poisson statistics, the mobility edge is expected to exhibit an intermediate level statistics associated with fractal wave functions, corresponding to the NLM precisely at the IRunstable fixed point. Although the width of the mobility edge (in the physical unit) shrinks under an increment of the lattice size, such ‘critical’ statistics [20] should be stable and depend only on the fixed point value of the conductance (which in turn depends on the dimensionality ), and possibly on the boundary condition and the aspect ratio of the lattice [21]. It should otherwise be universal in a sense that it originates from finetuning of a single relevant coupling constant (conductance) and all other irrelevant couplings should play no role [22]. These expectations have been verified numerically on the lattice [23, 24, 25]. For that purpose it is customary to use the uniform onsite randomness of width for a practical reason that the mean level density has a plateau in the band center, which attains maximal efficiency in the spectral averaging (as compared to Gaussian randomness that is theoretically easier to deal with). In Fig.2 we plot the LSDs from the central plateaux of the spectra of Anderson Hamiltonians at and for the orthogonal case [23] and and for the unitary case [24], on cubic lattices of various sizes, with periodic boundary conditions on all sides.
We observe that the LSDs from these eigenvalue windows are indeed scale invariant, and are in between and a hybrid of WignerDyson and Poisson,
(8) 
This observation signifies that the whole band center at these finetuned values of parameters indeed belong to the mobility edge corresponding to the presumed IRunstable fixed point of the NLM. We note that scaleinvariant critical statistics is also observed for the twodimensional Anderson Hamiltonian in the symplectic class [26], but not in other two classes, in accordance with the perturbatively functions (7).
The relationship between the anomalous behavior of the energy eigenfunctions and the fluctuation of the energy levels is widely accepted as follows: if the eigenfunctions from an energy window are multifractal [27], i.e., the inverse participation ratio scales as
(9) 
with in between localized () and extended () states, then two such eigenfunctions overlap only sparsely, for . Since only distant levels become less repulsive and less rigid, it modifies the tail of the LSD to quasiPoissonian, whereas the small behavior is not much affected, leading to (8).
4 Deformed random matrices
Once the existence of critical statistics is established, an immediate challenge is to derive analytically its statistical distributions, such as the LSD and the twolevel correlation function. Since the NLM (3) that originates from the microscopic theory can be solvable only for very exceptional quasi1D cases for which DuistermaatHeckman localization theorem is applicable [28], one natural path is to find a solvable effective model based on the symmetry and universality arguments. One might wonder how the Anderson localization transition, for which the dimension of the system is crucial, could possible be modeled by some effective RM ensemble, which obviously carries no information on the dimensionality. An answer to this frequently asked question is that it is the series of fractal dimensions in (9) that dictates the level statistics (such as level repulsion and rigidity) through the overlap of eigenfunctions, and not the actual dimension of the system. Accordingly, an effective model that possesses the crucial fractal property (9) and that reduces to the classical random matrices when the deformation is turned off, might as well reproduce the critical statistics. In this spirit, three seemingly different ensembles of random matrices have been proposed: (I) StieltjesWigert random matrices with a probability measure [29] (termed as ‘nonclassical ensemble’ in [30]), (II) powerlaw banded random matrices [31], and (III) 1D free fermions at finite temperature [32]. Instead of going deeply into each model, we merely mention the following key properties and refer the details to the original articles:

Ensemble I is a manifestly invariant ensemble. Through an unusual unfolding that results from the mean eigenvalue density and the cumulative density , any level correlator for the unitary class is expressed as in terms of a deformed spectral kernel
(10) for a range of the deformation parameter for which the level correlation is approximately translationinvariant.

Ensemble II has a fixed preferred basis in a sense that the probability measure of the matrix contains an extra factor that favors aligned to (in addition to the conventional Gaussian weight ), and possesses the property (9). Being Gaussian, it is directly mapped to a 1D NLM. For the unitary class, the connected twolevel correlator derived from the supersymmetric version of the NLM (3) takes the form for , to the order [33]. There the function originates from the spectral determinant of the diffusion operator in (3),
(11) This connection relates the phenomenologically introduced parameter in Ensemble I to a physical quantity (conductance) by .

Ensemble III is a modified version of Ensemble B, obtained by replacing the constant for an integration over a group manifold ( in the original version, recovering the invariance of GUE). Being equivalent to a system of free fermions at a temperature related to , the probability distribution of particles’ loci (i.e. of eigenvalues of random matrices) is readily determined, leading again to the same connected twolevel correlator .
Entrusting that these universality among three different ensembles [34] to be an indication of uniqueness of multifractal deformation of the classical random matrices, one of the authors (SMN) computed the LSDs of Ensemble I in three symmetry classes [35, 36] (Fig.3). Ensemble I is most suited for analytical computation of eigenvalue correlation, because established techniques for invariant RM ensembles (such as TracyWidom method of evaluating the Fredholm determinant [14], and the relationships between and [15]) are applicable, only with a replacement of the kernel from (5) to (10). For e.g. the unitary class, the LSD is expressed as in terms of the diagonal resolvent , which satisfies a transcendental equation of Painlevé VI type [35],
(12) 
under the boundary condition
Eq. (12) is a natural extension of the Painlevé V type equation [13] derived for the resolvent of the sine kernel (5) at . One of the authors (SMN) then confirmed that the LSDs at the mobility edge are well fitted, with a single tunable parameter , to these analytic formulas from the deformed kernel (blue curves in Fig.2) for three symmetry classes [36]. The dof of the fitting in the range (with binsize ) is as small as 0.23 (orthogonal) and 0.16 (unitary) for the case of and . This extremely precise matching not merely justifies the validity of the deformed random matrices as an effective model of critical statistics a posteriori, but can even be used as a criterion for an eigenvalue window in the spectrum of a disordered spectra to belong to its mobility edge that separates extended and localized states. Thus we move on to apply this criterion to the Dirac spectrum of QCD at the physical point, in the high temperature phase.
5 QCD Dirac spectra above T
The nature of Dirac eigenstates associated with small eigenvalues in the chirally symmetric phase (in which ) has been debated since the appearance of Ref.[5], partly motivated by the understanding of the fate of when SU is restored. The observation that lowlying eigenvalues at neither obey the Airy statistics of the ‘soft’ band edge of classical random matrices [37, 38], nor the statistics of random matrices at the (multi)critical edge [39, 40], left open the issue of characterizing these eigenvalues and associated eigenstates.
An important step was undertaken by GarcíaGarcía and Osborn who claimed that, right at the temperature of chiral symmetry restoration, the LSD from the spectral window near the origin becomes stable under the increment of spatial size and takes an intermediate form in between Poisson and WignerDyson statistics [41]. From this finding they speculated that Anderson localization is the microscopic mechanism for the chiral phase transition. However, they could not explicitly confirm a transition in the Dirac spectrum from Poisson to WignerDyson statistics. In retrospect, the reason for that was partly the lack of enough statistics that made it necessary to average spectral properties over spectral windows too wide. Another reason was that around the density of localized modes is very low and one needs large volumes to observe clear Poisson statistics.
In order to circumvent those issues, subsequently some of the authors performed simulations on much larger lattices and at temperatures sufficiently above [42, 43, 44]. By choosing the spectral window small enough for the spectral averaging to be justifiable, clear signs of the mobility edge were observed from the scale invariance of the LSD, even at small quark masses. The location of the mobility edge at a physical energy scale (well above the light quark mass) is seen to be stable under the increment of the spatial lattice size, indicating that a finite fraction of Dirac eigenstates per unit volume are localized in the thermodynamic limit.
[fm]  

2  0  2.60    16, 24, 32, 48  4  2.6  3k  256 
3  3.75  .125  24, 28, , 48  4  394MeV  7k40k  2561k 
The present study is based on two sets of lattice simulations, one quenched with the gauge group SU(2) and a dynamical SU(3) simulation with 2+1 flavors of stout smeared staggered quarks at the physical point. The parameters of the simulations are summarized in Table 1. More details of the quark action, scale setting and quark masses for the dynamical simulation can be found in Refs. [45, 46]. In the following we show that around the mobility edge in the spectrum, the unfolded level spacing distribution is well described by a deformed random matrix model. For determining the local spectral statistics throughout the spectrum, the eigenvalue windows are chosen as small as possible while attaining sufficient statistics. We tried spectral averaging over one (i.e. no spectral averaging) to two level spacings for the smallest lattice ( for SU(2) quenched), and of six to twelve level spacings for the largest lattice ( for SU(3) dynamical), around each designated point in the spectrum. The spacings are unfolded by their mean value within each window, . We plot the probability density of with a binsize of 0.05 for the interval , and find a bestfitting LSD of deformed random matrix ensembles (in symplectic class for SU(2) quenched and in unitary class for SU(3) dynamical) for these 80 values by varying the deformation parameter . Arbitrariness in the choice of the binsize could be avoided by fitting the cumulative LSD to the deformed random matrices, but in order to achieve acute sensitivity to the deformation parameter we employed the LSD itself for fitting. In Fig. 4 we exhibit sample plots of LSDs from Dirac eigenvalue windows centered at some s, and the bestfitting LSDs of deformed random matrices.
In the vicinity of the mobility edge (to be identified below), the results from deformed random matrices fit quite nicely to the Dirac data, with /dof of around 1. As the LSD changes through the spectrum so does the deformation parameter corresponding to the deformed random matrix ensemble providing the best fit to the data. In Fig. 5 we plot the bestfitting deformation parameter versus the center of each eigenvalue window . It turns out that the enveloping curves of for each lattice size are rather insensitive to the number of the level spacings within the window, so for the case of SU(3) dynamical we show in Fig.5 left a plot with a fixed number of LSs as: two for , four for , and six for .
From these figures one clearly sees that there are two regions in the spectrum. For small , the deformation parameter, , increases with the system size and is expected to go to its Poisson limit as . In this part of the spectrum, modes are localized. For large the deformation parameter decreases as the system becomes bigger and is expected to go to its WignerDyson statistics limit as . This corresponds to delocalized modes. In between, there is a fixed point where is independent of the lattice size. We identify this point with the mobility edge separating localized and delocalized eigenmodes in the thermodynamic limit.
By minimizing the variance of among various , the location of this mobility edge is determined as for SU(2) quenched and for SU(3) dynamical. The latter translates to MeV in physical unit. Furthermore, LSDs precisely from tiny eigenvalue windows including are scale invariant and are all juxtaposed right on the top of the predictions from the deformed RMs of symplectic class at , and of unitary class at (Fig. 6). /dof for the latter are . We noticed, rather unexpectedly, that the latter value of the deformation parameter at the mobility edge of the SU(3) Dirac spectrum is consistent with the value that corresponds to the mobility edge of the Anderson Hamiltonian in the unitary class on an isotropic 3D lattice (Fig. 2, right).
Finally we examine the shape parameters of the LSDs. A customary choice in the study of Anderson Hamiltonians is the variance of and the area up to the crossing point () of Wigner and Poisson distributions, . In Fig.7 we plot the shape parameters of LSDs from various eigenvalue windows, on lattices. One observes that as the location of the window moves form the origin to the bulk, shape parameters universally align on a specific curve connecting the Poisson and WignerDyson limits, regardless of the spatial size. Although this curve slightly deviates from the one from the deformed unitary random matrices (red curve), these two cross at one point (see the inset of Fig.7), precisely corresponding to the mobility edge. This findings may provide a further support on our claim that the mobility edge in the Dirac spectrum in the high temperature phase survives the thermodynamic limit, and its spectral fluctuation is described by the deformed random matrices characterizing the fixed point of the NLM from Anderson Hamiltonians.
6 Discussion: Critical statistics near the origin at T
We showed that the level spacing statistics across the Andersontype transition in the QCD staggeredDirac spectrum can be described by the same deformed RM ensemble that describes the transition in the Anderson model. Together with the compatibility of the critical exponents of these two models [47], our result lends further support to the idea that these two transitions indeed belong to the same universality class. Further confirmation of our findings in the case of GinspargWilson Dirac operators [42], close to the physical point, would be an important next step.
So far we discussed the transition in the spectrum at a fixed temperature above . To answer the original question of “whether chiral symmetry restoration is driven by the Anderson localization of the quasizero modes of the Dirac operator”, one should approach the pseudocritical temperature of QCD from above, carefully monitoring the level statistics within each fine spectral window. By extrapolating the dependence of the mobility edge on the temperature, two of us (TGK, FP) found that the mobility edge goes to zero around MeV [44], which is compatible with the pseudocritical temperature of the QCD transition.
It would be highly desirable to study how the spectral statistics changes as the mobility edge goes to zero at . In particular, it would be important to verify that the scale invariant spectral statistics is still described by the corresponding deformed random matrix model. There is, however, an additional complication here. It is only in the bulk of the spectrum that there is no difference between the chiral and the nonchiral version of the given matrix model. However, as the temperature is lowered and shifts to the spectrum edge, one has to use the chiral version of the matrix model to describe the spectral statistics.
For this purpose, GarcíaGarcía and Verbaarschot had already tailored a chiral version of the deformed RM ensemble of type III [48]. They have obtained an exact form of the twolevel correlation function , and in the approximation of keeping and discarding in the asymptotic limit and , it takes the form
(13) 
This ‘deformed Bessel’ kernel is the chiral counterpart of the ‘deformed sine’ kernel in (10) and reduces to this nonchiral deformed kernel in the limit with fixed. A practical problem of using the microscopic level density for fitting the Dirac spectrum is that becomes rather structureless at a finite deformation parameter (see Fig.8 below) [49]. On the other hand, the level number variance (an integral transform of the twolevel correlation function ) for large is sensitive to the parameter , but one would need extremely large lattices for the window at the origin containing dozens of eigenvalues to be uniformly fitted with a single parameter .
To circumvent this practical problem, we propose an alternative strategy of using the individual eigenvalue distributions [50] that have characteristic peaky shapes of their own and respond sensitively to the deformation, instead of the spectral density that comprises of these peaks, . One could in principle apply TracyWidom method to (13) to obtain a closed analytic equation of Painlevé type for analogous to (12), but for the actual fitting purpose it is sufficient to evaluate the Fredholm determinant and the resolvents by the Nyströmtype approximation [51, 4]. Here we employ the Gaussian quadrature of order and exhibit the distributions of the five smallest eigenvalues as the deformation parameter is varied in the range (Fig.8).
Computational details, including applications to orthogonal and symplectic ensembles, will appear elsewhere. These distributions have WignerDyson and Poisson asymptotics on different sides of the peaks,
(14) 
Because of this very characteristic shapes,
they are unmistakable candidates for fitting the eigenvalues from the mobility edge occurring at the origin.
After verifying that eigenvalues in the mobility edge of the chiral Anderson Hamiltonian occurring
around the origin is described by these distributions,
and including the effect of quark masses in the regime [52, 50],
we consider it a challenging but feasible task to compare the distributions of
each of smallest QCD Dirac eigenvalues at and around and to answer to the point
quoted in the Introduction.
Acknowledgements: It is our great pleasure to thank the organizers of LATTICE2013 for the chance of presenting a work that was in progress and unpublished at the time of the conference, in the plenary session. SMN thanks E. Itou and J.J.M. Verbaarschot for useful communications.
References
 M. L. Mehta, Random Matrices, 3rd ed. (Academic Press, New York, 2004).
 E. V. Shuryak and J. J. M. Verbaarschot, Nucl. Phys. A 560, 306 (1993).
 P. H. Damgaard, U. M. Heller, K. Splittorff, and B. Svetitsky, Phys. Rev. D 72, 091501 (2005).
 S. M. Nishigaki, Prog. Theor. Phys. 128, 1283 (2012); Phys. Rev. D 86, 114505 (2012).
 M. A. Halasz and J. J. M. Verbaarschot Phys. Rev. Lett. 74, 3920 (1995).
 F. Wegner, Z. Phys. B 36, 209 (1980).
 K. B. Efetov, Adv. Phys. 32, 53 (1983).
 P. W. Anderson, Phys. Rev. 109, 1492 (1958).
 E. Hofstetter and M. Schreiber, Phys. Rev. B 49, 14726 (1994).
 E. Hofstetter and M. Schreiber, Phys. Rev. Lett. 73, 3137 (1994).
 S. N. Evangelou and D. E. Katsanos, J. Stat. Phys. 85, 525 (1996).
 J .J. M. Verbaarschot, H. A. Weidenmüller, and M. R. Zirnbauer, Phys. Rep. 129, 7367 (1985).
 M. Jimbo, T. Miwa, Y. Môri, and M. Sato, Physica D 1, 80 (1980).
 C. A. Tracy and H. Widom, Comm. Math. Phys. 163, 33 (1994).
 C. A. Tracy and H. Widom, Comm. Math. Phys. 177, 727 (1996).
 P. Desrosiers and P. J. Forrester, Nonlinearity 19, 1643 (2006).
 F. Wegner, Nucl. Phys. B 316, 663 (1989).
 B. Kramer and A. MacKinnon, Rep. Prog. Phys. 56, 1469 (1993).
 A. D. Mirlin, F. Evers, I. V. Gornyi, and P. M. Ostrovsky, Int. J. Mod. Phys. B 24, 1557 (2010).
 B. I. Shkhlovskii, B. Shapiro, B. R. Sears, P. Lambrianides, and H. B. Shore, Phys. Rev. B 47, 11487 (1993).
 L. Schweitzer and H. Potempa, Physica A 266, 486 (1999); J. Phys. Cond. Matter 10, L431 (1998).
 B. Kramer and A. MacKinnon, Phys. Rev. Lett. 47, 1546 (1981).
 S. N. Evangelou, Phys. Rev. B 49, 16805 (1994); I. K. Zharekeshev and B. Kramer, Phys. Rev. B 51, 17239 (1995); I. Varga, E. Hofstetter, M. Schreiber, and J. Pipek, Phys. Rev. B 52, 7783 (1995).

E. Hofstetter,
Phys. Rev. B 54, 4552 (1996);
M. Batsch, L. Schweitzer, I. K. Zharekeshev, and B. Kramer, Phys. Rev. Lett. 77, 1552 (1996).  T. Kawarabayashi, T. Ohtsuki, K. Slevin, and Y. Ono, Phys. Rev. Lett. 77, 3593 (1996).
 Y. Asada, K. Slevin, and T. Ohtsuki, Phys. Rev. Lett. 89, 256601 (2002).
 M. Schreiber and H. Grussbach, Phys. Rev. Lett. 67, 607 (1991).
 A. Lamacraft, B. D. Simons, and M. R. Zirnbauer, Phys. Rev. B 70, 075412 (2004).
 K. A. Muttalib, Y. Chen, M. E. H. Ismail, and V. N. Nicopoulos, Phys. Rev. Lett. 71, 471 (1993).
 P. J. Forrester, LogGases and Random Matrices (Princeton Univ. Press, Princeton, 2010).
 A. D. Mirlin, Y. V. Fyodorov, F.M. Dittes, J. Quezada, and T. H. Seligman, Phys. Rev. E 54, 3221 (1996).
 M. Moshe, H. Neuberger, and B. Shapiro, Phys. Rev. Lett. 73, 1497 (1994).
 A. V. Andreev, B. D. Simons, and B. L. Altshuler, J. Math. Phys. 37, 4968 (1996).
 V. E. Kravtsov and K. A. Muttalib, Phys. Rev. Lett. 79, 1913 (1997).
 S. M. Nishigaki, Phys. Rev. E 58, R6915 (1998).
 S. M. Nishigaki, Phys. Rev. E 59, 2853 (1999).
 P. H. Damgaard, U. M. Heller, R. Niclasen, and K. Rummukainen, Nucl. Phys. B 583, 347 (2000).
 F. Farchioni, P. de Forcrand, I. Hip, C. B. Lang, and K. Splittorff, Phys. Rev. D 62, 014503 (2000).
 A. D. Jackson, M. K. Şener, and J. J. M. Verbaarschot, Nucl. Phys. B 479, 707 (1996).
 G. Akemann, P. H. Damgaard, U. Magnea, and S. M. Nishigaki, Nucl. Phys. B 519, 682 (1998).
 A. M. GarcíaGarcía and J. C. Osborn, Phys. Rev. D 75, 034503 (2007).
 T. G. Kovács, Phys. Rev. Lett. 104, 031601 (2010).
 T. G. Kovács and F. Pittler, Phys. Rev. Lett. 105, 192001 (2010).
 T. G. Kovács and F. Pittler, Phys. Rev. D 86, 114515 (2012).
 Y. Aoki, Z. Fodor, S. D. Katz, and K. K. Szabó, JHEP 0601, 089 (2006).
 S. Borsányi, G. Endrődi, Z. Fodor, A. Jakovác, S. D. Katz, S. Krieg, C. Ratti, and K. K. Szabó, JHEP 1011, 077 (2010).
 M. Giordano, T. G. Kovács, and F. Pittler, PoS LATTICE 2013 (2013) 213.
 A. M. GarcíaGarcía and J. J. M. Verbaarschot, Nucl. Phys. B 586, 668 (2000).
 A. M. GarcíaGarcía and K. Takahashi, Nucl. Phys. B 700, 361 (2004).
 P. H. Damgaard and S. M. Nishigaki, Phys. Rev. D 63, 045012 (2001).
 F. Bornemann, Math. Comp. 79, 871 (2010).
 P. H. Damgaard and S. M. Nishigaki, Nucl. Phys. B 518, 495 (1998).