Universal hard-edge statistics of non-Hermitian random matrices

Random matrix theory is a powerful tool for understanding spectral correlations inherent in quantum chaotic systems. Despite diverse applications of non-Hermitian random matrix theory, the role of symmetry remains to be fully established. Here, we comprehensively investigate the impact of symmetry on the level statistics around the spectral origin -- hard-edge statistics -- and expand the classification of spectral statistics to encompass all the 38 symmetry classes of non-Hermitian random matrices. Within this classification, we discern 28 symmetry classes characterized by distinct hard-edge statistics from the level statistics in the bulk of spectra, which are further categorized into two groups, namely the Altland-Zirnbauer$_0$ classification and beyond. We introduce and elucidate quantitative measures capturing the universal hard-edge statistics for all the symmetry classes. Furthermore, through extensive numerical calculations, we study various open quantum systems in different symmetry classes, including quadratic and many-body Lindbladians, as well as non-Hermitian Hamiltonians. We show that these systems manifest the same hard-edge statistics as random matrices and that their ensemble-average spectral distributions around the origin exhibit emergent symmetry conforming to the random-matrix behavior. Our results establish a comprehensive understanding of non-Hermitian random matrix theory and are useful in detecting quantum chaos or its absence in open quantum systems.


I. INTRODUCTION
Symmetry plays a pivotal role in random matrix theory and is crucial for comprehending spectra of complex systems [1].Dyson's threefold way initially classifies Hermitian random matrices according to time-reversal symmetry (TRS) [2].The bulk of spectra is characterized by the threefold universal level statistics, termed the Wigner-Dyson statistics [3,4].These universal level statistics appear in the following renowned random matrix ensembles: Gaussian unitary ensemble (lacking TRS), Gaussian orthogonal ensemble (possessing TRS with sign +1), and Gaussian symplectic ensemble (possessing TRS with sign −1).Extensive theoretical and experimental investigations have demonstrated that the spectral correlations of chaotic closed quantum systems align with these universal random-matrix statistics with corresponding symmetry [5][6][7][8][9].Subsequently, random matrix ensembles in the chiral and Bogoliubov-de Gennes (BdG) symmetry classes, which incorporate chiral symmetry [10,11] and particle-hole symmetry [12], were introduced.These symmetries dictate universal properties in the proximity of the symmetry-preserving point E = 0 in the energy spectrum.Together with TRS, chiral and particlehole symmetries form the comprehensive classification of Hermitian random matrices, known as the Altland-Zirnbauer (AZ) tenfold classification [12].Importantly, this well-established Hermitian random matrix theory has found broad applications across diverse fields of physics.In condensed matter physics, it elucidates universal transport phenomena of disordered metals and su-perconductors [13,14] and also serves as a powerful tool to detect the Anderson transitions [15] and many-bodylocalization transitions [16,17].In nuclear physics, random matrices provide an effective model of chiral symmetry breaking [18].
While researchers historically focused on Hermitian operators, non-Hermitian operators are equally important in various physical systems.Prime examples include chiral Dirac operators with chemical potential in quantum chromodynamics (QCD) [19][20][21] and scattering matrices in quantum transport phenomena [13,14].Non-Hermitian matrices appear naturally in stochastic processes and network systems, having important applications in biological systems [22][23][24][25][26][27][28][29].Realistic physical systems are inevitably coupled with their surrounding environments, and non-Hermitian operators also play pivotal roles in the dynamics of such open physical systems [30].For the Markovian environment, the dynamics of open quantum systems is described by the Lindblad master equation [31,32], where the generator of the dynamics, known as the Lindbladian, is a non-Hermitian operator.The time evolution of classical synthetic materials such as optical and photonic systems with energy gain or loss [33][34][35][36], along with quantum systems under continuous measurement [37][38][39][40], can also be characterized by non-Hermitian Hamiltonians.Thanks to advances in state-of-the-art experiments and theories, the past years have witnessed remarkable development in non-Hermitian physics [41][42][43][44].
Non-Hermitian random matrix theory is useful for understanding the spectral correlations and dynamics of open quantum systems.Universal dynamical features of open quantum chaotic systems were investigated using non-Hermitian random matrices [73,75,76,80].It was conjectured that the statistics of complex energy levels in non-integrable open quantum systems follow those of non-Hermitian random matrices, which has been numerically verified in several models [71,[77][78][79]82].Level correlations of non-Hermitian random matrices in the bulk of complex spectra (i.e., away from special points, lines, and edges of complex spectra) belong to the threefold universality classes [72,77,107].These universality classes depend only on an extension of TRS for non-Hermitian operators, time-reversal symmetry † [46], while the influence of other defining symmetries appears elusive.Given the prevalence of almost all 38 symmetry classes in realistic physical systems (e.g., see Refs.[94,97,98]), characterizing the unique universality classes of level correlations in each symmetry class and exploring the impact of symmetry stands out as one of the crucial challenges in non-Hermitian random matrix theory.Despite its importance and some discrete works on level statistics beyond the threefold way, especially those on non-Hermitian but real random matrices [108][109][110][111][112], a comprehensive understanding of the universality classes within the 38-fold symmetry classification is still lacking.
In this work, we demonstrate that the level statistics around the spectral origin faithfully reveal the impact of all defining symmetries.In Hermitian random matrix theory, such level statistics around the spectral origin are known as the hard-edge statistics.While the spectral origin of complex spectra may look like a hole or node, we find a rich analogy of the level statistics with the Hermitian case.Thus, we call the level statistics around the spectral origin the hard-edge statistics even in non-Hermitian random matrix theory.We comprehensively explore possible universality classes of the hard-edge statistics of non-Hermitian random matrices in all 38 symmetry classes.We identify 28 symmetry classes where the level statistics around the spectral origin differ from those in the bulk of the spectra or around the real or imaginary axis by analyzing their defining symmetries.This approach is akin to that in Hermitian random matrix theory, where the hard-edge statistics distinct from the bulk statistics appear in the chiral and BdG symmetry classes but not in the Wigner-Dyson symmetry classes [1,12,14].By contrast, non-Hermiticity enriches the symmetry classification and leads to more diverse universal spectral correlations.Specifically, we categorize these 28 symmetry classes into two groups and investigate all of them, summarized in Tables II and III.Notably, similar to the Hermitian counterparts, only a few eigenvalues closest to the spectral origin contribute to the hard-edge statistics.Consequently, analyzing the hard-edge statistics and diagnosing symmetry classes requires studying ensembles of operators.
The first group, summarized in Table II, comprises seven symmetry classes, where the spectral origin is the only high-symmetry point in the complex plane, and the hard-edge statistics exhibit spectral U(1) rotation symmetry.These seven symmetry classes, together with the threefold classes based on time-reversal symmetry † , form the tenfold classification, which we dub Altland-Zirnbauer 0 (AZ 0 ) classification.We employ the density of complex eigenvalues and the distributions of the complex level ratios around the spectral origin to characterize the universality classes of the hard-edge statistics.We numerically obtain seven universal ratio distributions and summarize their properties in Table II.In a similar spirit to the Wigner surmise, we calculate these distributions analytically for small non-Hermitian random matrices.We show that the asymptotic behavior observed in large matrices can be accurately estimated from that in small ones.Notably, distributions of the complex eigenvalue with the minimum modulus were investigated previously [94,113,114].By contrast, we introduce a new statistical quantity that captures the hard-edge statistics-ratio of the complex eigenvalue with the minimum modulus to that with the second minimum modulus, which we dub the complex level ratio.We clarify the difference between the two approaches and discuss the potential advantages of the complex level ratio.
The second group, summarized in Table III, comprises the remaining 21 symmetry classes which have no spectral U(1) rotation symmetry around the origin.We investigate the distributions of the eigenvalues z min with the smallest modulus.We numerically demonstrate that in certain symmetry classes, the distributions exhibit deltafunction peaks on the real and/or imaginary axes, in contrast to the preceding seven symmetry classes in the AZ 0 classification.We find that the probability of z min being real or purely imaginary is characteristic and universal in each symmetry class, summarizing these universal values for all these 21 symmetry classes in Table III.In particular, we focus on three representative symmetry classes (i.e., classes BDI, CII, and AII + S + ) and carefully examine their universal distributions of z min .The obtained distributions are qualitatively different when z min is complex, real, or purely imaginary, and in certain classes, the distributions exhibit the point group symmetry D 4 .We provide theoretical explanations for our numerical findings.
To further demonstrate the universality of our newly found hard-edge statistics, we investigate diverse open quantum systems, including quadratic and many-body Lindbladians, as well as non-Hermitian Hamiltonians.For the seven symmetry classes in the AZ 0 classification, we construct corresponding physical models.An ensemble of a physical model often has statistical symmetry (ensemble symmetry), while symmetry of a single realization of the model determines the relevant symmetry class.Notably, the ensemble symmetry of physical models is generally lower than the ensemble symmetry of random matrices in the same symmetry class.Nevertheless, the hard-edge statistics of non-Hermitian random matrices in the AZ 0 classification manifest themselves in these physical systems in the seven symmetry classes, signaling emergent spectral U(1) symmetry.For the 21 symmetry classes beyond the AZ 0 classification, we construct models in the three representative symmetry classes.We show that they also exhibit the hard-edge statistics and emergent spectral D 4 symmetry consistent with non-Hermitian random matrices.
The rest of this work is organized as follows.In Sec.II, we begin with the symmetry classification of non-Hermitian random matrix theory.We identify the symmetry classes exhibiting the unique hard-edge statistics.In Secs.III and IV, we investigate the hard-edge statistics of non-Hermitian random matrices in the Gaussian ensembles within and beyond the AZ 0 classification, respectively.Sections V and VI are respectively dedicated to open quantum physical models within and beyond the AZ 0 classification.We uncover that the hard-edge statistics of these physical models coincide with non-Hermitian random matrices, demonstrating that our characterization of the hard-edge statistics efficiently captures the chaotic behavior or integrability of open quantum systems.Section VII is devoted to conclusions and discussions.

II. SYMMETRY CLASSIFICATION OF NON-HERMITIAN RANDOM MATRIX THEORY
We start with the AZ classification of Hermitian operators H.It consists of two types of anti-unitary symmetry, time-reversal symmetry (TRS): (1) particle-hole symmetry (PHS): and their combined unitary symmetry chiral symmetry (CS): where T + , C − , and Γ are all unitary operators.Non-Hermiticity H ̸ = H † ramifies the definitions of TRS, PHS, and CS, and enriches the symmetry classification [51].Replacing H in their definitions by H † results in two types of additional anti-unitary symmetry, time-reversal symmetry † (TRS † ) : particle-hole symmetry † (PHS † ) : and additional unitary symmetry, sublattice symmetry (SLS): The combination of CS and SLS gives another unitary symmetry pseudo-Hermiticity (pH): These four symmetries, together with the three original symmetries in the AZ classification, constitute the 38-fold symmetry classification [45,46] . Thus, the particle-hole operation generally transforms an eigenvalue z to another eigenvalue −z, and hence the complex spectrum is symmetric about the origin.Notably, while PHS, in general, creates opposite-sign pairs (z, −z) of complex eigenvalues, it imposes a special constraint on the eigenvalue invariant under the particle-hole operation, i.e., the zero eigenvalue z = 0.If a zero mode exists, its left eigenvector coincides with v T C † − .While exact zero modes typically do not appear in generic matrices, the level interaction around the spectral origin should be influenced by PHS.Similarly, SLS relates an eigenvalue z to −z.On the other hand, TRS, PHS † , CS, and pH accompany complex or Hermitian conjugation: while TRS and pH relate an eigenvalue z to z * , PHS † and CS relate an eigenvalue z to −z * .In contrast to these symmetries, TRS † does not lead to such pairs of eigenvalues between TABLE I. Altland-Zirnbauer † (AZ † ) classification of non-Hermitian random matrices based on time-reversal symmetry (TRS), time-reversal symmetry † (TRS † ), and pseudo-Hermiticity (pH).If H respects TRS or pH, iH respects particle-hole symmetry † (PHS † ) or chiral symmetry (CS), respectively.If H belongs to the symmetry class in the column "Class (H)", iH belongs to the corresponding class in "Class (iH)".In the nomenclature of symmetry classes, η represents pH, and the subscript of η± denotes whether pH commutes (+) or anti-commutes (−) with TRS † .In the columns "TRS (PHS † )" and "TRS † ", the entry "±1" specifies the sign of the anti-unitary symmetry, and "0" means its absence.The influence of symmetries on the level statistics around the real axis is summarized [112].If there is a soft gap in the density of complex eigenvalues, ρc(x+iy), around the real axis y = 0, its small-|y| behavior is listed.The column "Real eigenvalues" indicates whether random matrices have a finite number of real eigenvalues ( √ ) or not (×); for the five symmetry classes with " √ ", the number N real of real eigenvalues after the ensemble average is subextensive, N real ∝ √ N , with N being the dimension of the matrix.

Class (H)
Density ρc(x + iy) around the real axis Real eigenvalues different eigenvalues but imposes constraints on generic individual eigenvalues.
Importantly, the influence of symmetry on the level statistics is contingent upon the invariance of eigenvalues under the symmetry operation.This is also the case for Hermitian random matrices, where TRS is relevant to generic real eigenvalues, but PHS and CS (or equivalently, SLS) are relevant only to the zero eigenvalue.As also discussed above, TRS † constrains generic complex eigenvalues and influences the spectral correlations in the bulk of complex spectra.On the other hand, the other six types of symmetry affect the level statistics only around the symmetry-preserving lines or points.For example, TRS transforms a complex eigenvalue z to another complex eigenvalue z * , as discussed previously.As a result, TRS is respected only by eigenvalues satisfying z = z * , i.e., real eigenvalues z ∈ R, and hence changes the spectral statistics only around the real axis.Similar to TRS, the other three types of symmetry, PHS † , CS, and pH, are also relevant to the level statistics only around the real or imaginary axis.By contrast, the remaining two types of symmetry, PHS and SLS, are not necessarily relevant even around the real or imaginary axis and influence the level statistics only around the spectral origin.
Consistent with this general understanding, numerical investigations confirmed that symmetries other than TRS † do not influence the level statistics within the bulk of spectra (i.e., generic eigenvalues away from the spectral origin, as well as real and imaginary axes) [77].TRS † alone gives rise to the threefold classification, comprising classes A, AI † , and AII † (see the first three rows in Table I or II), and yields the three distinct universal-ity classes of level statistics for generic complex eigenvalues [77].Specially, let z α be a complex eigenvalue, z α;NN be its nearest-neighbor eigenvalue (i.e., |z α;NN − z α | = min β |z α − z β |), and z α;NNN be its next-nearest-neighbor eigenvalue.The distributions of the complex level spacing [71,72,77,107] z α;NN − z α and the complex levelspacing ratio (z α;NN −z α )/(z α;NNN −z α ) [79] were used to characterize the spectral correlations and found to converge to universal level statistics for large random matrices.
Furthermore, TRS, TRS † , and pH give rise to the tenfold classification (see Table I).This tenfold classification is equivalent to the Altland-Zirnbauer † (AZ † ) classification [46] since TRS is equivalent to PHS † [51].In this classification, the three symmetry classes (i.e., classes A, AI † , and AII † ) where real eigenvalues do not have higher symmetry than generic ones are characterized only by TRS † and already included in the previous threefold classification.Consequently, 10 − 3 = 7 distinct universality classes of level statistics on and around the real axis emerge.It was found that these symmetries determine whether random matrices have real eigenvalues, as well as the behavior of the density of complex eigenvalues around the real axis (see a summary in Table I) [55,[107][108][109][110]112].The level spacings or level-spacing ratios of real eigenvalues (if exist) also exhibit universal and characteristic distributions in each symmetry class [112].The study of level statistics on or around the imaginary axis does not yield new universality classes since the real and imaginary axes are transformed into each other via the multiplication by i.
All the seven symmetries in Eqs. ( 1)-( 7) within the 38- fold classification influence the level statistics around the spectral origin.In the tenfold AZ † classification discussed above, the spectral origin does not possess higher symmetry compared to real eigenvalues.In the remaining 38 − 10 = 28 symmetry classes, by contrast, the spectral origin exhibits higher symmetry and hence is special.In this work, we show that the unique level statistics indeed appear around the spectral origin in these 28 symmetry classes.We categorize these 28 symmetry classes into two groups according to whether they involve symmetries associated with complex or Hermitian conjugation (i.e., TRS, PHS † , CS, and pH).
In the first group, symmetry classes involve TRS † , PHS, and SLS, but do not involve symmetries associated with complex or Hermitian conjugation.TRS † , PHS, and SLS lead to the tenfold classification in Table II different from the AZ † classification in Table I.Here, we call this additional tenfold classification Altland-Zirnbauer 0 (AZ 0 ) classification since it is concerned with the spectral origin.Besides the three symmetry classes (i.e., classes A, AI † , and AII † ) that only involve TRS † and are counted before, the seven symmetry classes (i.e., classes AIII † , BDI 0 , CII 0 , D, C, CI 0 , and DIII 0 ) are in the first group.Owing to PHS, SLS, or both of them, which maps an eigenvalue z to another eigenvalue −z, the complex spectra of generic matrices in any of these symmetry classes are invariant under the reflection about the spectral origin and hence respect the point group symmetry C 2 [see a schematic in Fig. 1 (a)].Away from the origin, the spectral statistics depend solely on TRS † since PHS and SLS merely create opposite-sign pairs of eigenvalues, (z, −z).This is also the case for non-zero real or purely imaginary eigenvalues.Around the origin, by contrast, even PHS and SLS change the spectral correlations and lead  1)-( 7)].Timereversal symmetry † (TRS † ), particle-hole symmetry † (PHS † ), and chiral symmetry (CS) form the 10-fold AZ † classification (blue circle; see Table I for details).TRS † , particle-hole symmetry (PHS), and sublattice symmetry (SLS) form the 10fold AZ0 classification (orange circle; see Table II   to unique spectral statistics, as shown below.If a matrix H belongs to one of these seven symmetry classes, e iϕ H also belongs to the same symmetry class for arbitrary ϕ ∈ R (see Table II).
The remaining 28 − 7 = 21 symmetry classes are in the second group and involve at least one of TRS, PHS † , CS, and pH, as summarized in Table III.These symmetries accompany complex or Hermitian conjugation and hence map an eigenvalue z to z * or −z * .In these classes, SLS, PHS, or both of them also exist.Consequently, the real or imaginary axis exhibits higher symmetry than generic points, and the spectral origin exhibits even higher symmetry.Owing to the combination of these symmetries, H does not necessarily belong to the same symmetry class as e iϕ H for generic ϕ ∈ R, in contrast to the AZ 0 symmetry classification.Rather, complex spectra of generic non-Hermitian matrices in these 21 symmetry classes are invariant under the reflections about the real and imaginary axes, and hence respect the point group symmetry D 2 [see a schematic in Fig. 1 (b)].Figure 2 summarizes the aforementioned subgroups of symmetry classes and illustrates their relationship.
Below, we comprehensively study the spectral statistics around the origin-hard-edge statistics-of non-Hermitian random matrices in all the AZ 0 symmetry classes (Sec.III), as well as the other 21 symmetry classes with a particular focus on three typical classes (classes BDI, CI, and AII + S + ) (Sec.IV).Through these analyses, we expand the classification of the universal spectral statistics of non-Hermitian random matrices to encompass 38-fold symmetry classes.
TABLE II.Altland-Zirnbauer0 (AZ0) classification of non-Hermitian random matrices based on time-reversal symmetry † (TRS † ), particle-hole symmetry (PHS), and sublattice symmetry (SLS).The column "Equivalent class" shows the equivalent symmetry classes in Ref. [46].In the columns "TRS † " and "PHS", the entry "±1" specifies the sign of the anti-unitary symmetry, and "0" means its absence.The density ρ(|z|) of the modulus of eigenvalues and the level-ratio distribution pr(r) [Eq.(11)] share the same small-argument behaviors, which are summarized in the columns "ρ(|z|) and pr(|z|)".The mean values ⟨r⟩ and ⟨cos θ⟩ of random matrices are summarized, where re iθ is the level ratio around the spectral origin.The statistics are obtained by diagonalizing 10 7 realizations of 128 × 128 random matrices.The numbers in the parentheses denote standard deviations.

Class
Equivalent class TRS † PHS SLS ρ(|z|) and pr(|z|) FIG. 3. Density of complex eigenvalues of non-Hermitian random matrices in the Gaussian ensemble within the AZ0 classification.While the density is nearly constant away from the spectral origin, it deviates from this constant around the origin, manifesting distinct spectral correlations around the origin due to higher symmetry.The behavior of the density, as well as the concomitant spectral correlations, around the origin varies among different symmetry classes.In each symmetry class, the result is obtained by 10 7 samples of 128 × 128 random matrices.In classes (c) CII0 and (g) DIII0, all eigenvalues are twofold degenerate, and we only count the degenerate pairs once when calculating the density.

III. NON-HERMITIAN RANDOM MATRICES IN THE AZ0 CLASSIFICATION
We consider non-Hermitian random matrices in the seven symmetry classes within the AZ 0 classification (Table II).Specifically, we take the Gaussian ensemble with the probability density function with a positive constant β > 0. We follow the terminology in the AZ classification: classes AIII † , BDI 0 , and CII 0 are categorized as chiral symmetry classes, and classes D, C, CI 0 , and DIII 0 are categorized as BdG classes.The explicit forms of non-Hermitian random matrices in each symmetry class are enumerated in Appendix A. Similar to Hermitian random matrix theory [14], the total number νd E of topologically protected exact zero modes can manifest themselves and influence the level statistics around the origin in some symmetry classes.Here, d E refers to the degrees of degeneracy of generic eigenvalues: d E = 2 in the presence of TRS † with the negative sign −1 (i.e., classes AIII † , CII 0 , and DIII 0 ) and d E = 1 otherwise.Specifically, a non-Hermitian random matrix H in class D (H T = −H) has νd E = 1 zero mode if the matrix dimension of H is odd.A non-Hermitian random matrix H in class DIII 0 (H T = −H and σ y H T σ y = H) has νd E = 2 zero modes if the matrix dimension of H is 4n + 2 (n ∈ Z).For non-Hermitian random matrices in classes AIII † , BDI 0 , and CII 0 , on the other hand, νd E equals the difference between the numbers of positive and negative eigenvalues of the symmetry operator S. While random matrices having such zero modes can be physically relevant, for example, to non-Hermitian superconductors with vortices [117,118], dissipative SYK-type models [119], and QCD with topological charges [19], we focus on random matrices with νd E = 0 in this work.
We numerically diagonalize 10 7 samples of N × N (N = 128) non-Hermitian random matrices in each symmetry class.Let z 1 be the eigenvalue with the smallest modulus.We normalize eigenvalues such that ⟨|z 1 |⟩ = 1, where ⟨•⟩ denotes the average over the disorder ensemble.After the ensemble average, the density of complex eigenvalues does not depend on the phases of eigenvalues but only on their modulus (Fig. 3).It is rotationally invariant with respect to the origin and exhibits spectral U(1) rotation symmetry.This spectral U(1) symmetry arises because the Gaussian ensemble includes e iϕ H with the same probability for arbitrary ϕ ∈ R. The eigenvalues are distributed almost uniformly in a circle except around the origin and boundary.In classes BDI 0 and D, the density of eigenvalues around the origin is non-zero and higher than that in the bulk, revealing the absence of level repulsion at the spectral origin.On the other hand, in classes AIII † , CII 0 , C, CI 0 , and DIII 0 , the density vanishes at the origin.
Notably, the small-|z| behavior of the density of the modulus of eigenvalues z i , characterizes the different strength of the level repulsion in different symmetry classes (Fig. 4).Among the three chiral symmetry classes, the descending order of ρ(|z|) for |z| ≪ 1 is classes BDI 0 , AIII † , and CII 0 [see Fig. (classes CII 0 and C) .
(10) The coefficients of this proportional relationship are given in the caption of Fig. 4. As noted above, the linear decay of ρ(|z|) for classes D and BDI 0 does not originate from the level repulsion but solely from the integral measure along the angular direction.The higher density of eigenvalues, given as ρ(|z|)/π|z|, around the spectral origin compared to that in the bulk is a characteristic feature of these two symmetry classes.By contrast, the cubic decay (with a potential logarithmic correction) for the other symmetry classes indicates the level repulsion around the spectral origin.In a similar spirit to the Wigner surmise, the qualitative behavior of level statistics for large random matrices can be estimated from those of small matrices.In Appendix B, we analytically calculate ρ(|z|) of small random matrices, exhibiting the identical small-|z| behavior observed in large ones.
To further characterize the level correlations around the spectral origin, we introduce a new quantitative measure of the hard-edge statistics for non-Hermitian operators.Owing to PHS, SLS, or both of them, complex eigenvalues must appear in the opposite-sign pairs (+z, −z).Let ±z 1 , ±z 2 , ±z 3 , . . .be eigenvalues (or Kramers pairs of eigenvalues) of H in the ascending order of the modulus (i.e., |z  11)] around the spectral origin for the seven symmetry classes in the AZ0 classification.The error bars represent one standard deviation evaluated by the bootstrap method [116].
introduce the complex level ratio by We characterize the level-ratio statistics by the distribution p r (r) of the module r of the ratio and the distribution p θ (θ) of its angle (phase) θ, where we have r ∈ [0, 1] and θ ∈ [−π/2, π/2) by definition.When the levels are uncorrelated, In Sec.V D, we show that such level-ratio statistics indeed appear in an integrable open quantum system.The complex level ratio is the ratio between two complex eigenvalues, whereas a complex spacing ratio previously used for studies of spectral correlations in the bulk spectra [79,120,121] is a ratio between two spacings of eigenvalues.For Hermitian operators, one can arrange real eigenvalues λ i 's in ascending order and introduce the spacing ratio as min(s i , s i−1 )/ max(s i , s i−1 ) with s i ≡ λ i+1 − λ i [120,121].For non-Hermitian operators, one can select a complex eigenvalue z α , find its nearest-neighboring eigenvalue z α;NN and next-nearestneighboring eigenvalue z α;NNN in the complex plane, and introduce the complex spacing ratio as (z α;NN − z α )/(z α;NNN − z α ) [79].The complex level ratio and the complex spacing ratio share the same aim-to utilize quantities independent of the normalization of the random-matrix ensembles.Since the former involves two eigenvalues while the latter involves three, the two measures quantify the different types of spectral correlations.
For the random-matrix ensemble, both the radial distribution p r (r) and the angular (phase) distribution p θ (θ) take non-trivial forms owing to the level correlations between the level pairs (+z 1 , −z 1 ) and (+z 2 , −z 2 ).From the numerical diagonalization of 128 × 128 random matrices, we obtain the radial distribution p r (r) and angular distribution p θ (θ) of the complex level ratio re iθ (Fig. 5).In each symmetry class, the distributions p r (r) and p θ (θ) show characteristic behavior, significantly distinct from that of uncorrelated levels.The small-r behavior of p r (r) is similar to the small-|z| behavior of the density ρ(|z|) of the modulus of eigenvalues in the same symmetry class.For example, among the three chiral symmetry classes, the descending order of p r (r) for r ≪ 1 is the same as that of ρ(|z|) for |z| ≪ 1 [compare Fig. 4 (a) with Fig. 5 (a)].Fitting p r (r) with r ≪ 1, we find that p r (r) shows the same power-law decay as ρ(|z|) in Eq. ( 10) and Table II.The density ρ(|z|) should be mostly contributed by z 1 for |z 1 | ≪ 1, leading to r ∝ |z 1 |.Consequently, the small-argument behaviors of both ρ(|z|) and p r (r) are mainly determined by the level correlations between the eigenvalue pair (+z 1 , −z 1 ) closest to the spectral origin, which underlies the similarity between ρ(|z|) and p r (r).Moreover, the angular distributions p θ (θ) are unique in different symmetry classes and characterize the different universality classes [Fig.5 (b)].First, we have p θ (θ) = p θ (−θ) since the matrix ensembles are invariant under complex conjugation, transforming θ to −θ.Second, p θ (θ) reaches its minimum and maximum at θ = 0 and π/2, respectively.The pairs (+z 1 , −z 1 ) and (+z 2 , −z 2 ) tend to be oriented perpendicularly to each other rather than in parallel, revealing the presence of level correlations.We also obtain the analytic formulas of p N =4 The mean values ⟨r⟩ and ⟨cos θ⟩ of the complex level ratio re iθ serve as useful quantitative indicators of the hardedge statistics among the different universality classes.These values in each symmetry class are summarized in Table II.Notably, the distribution of r in each symmetry class should converge to a universal form for sufficiently large matrix sizes.To clarify this convergence, we utilize ⟨r⟩ and ⟨cos θ⟩ as functions of the matrix size N and find that they already reach the convergence for N ≥ 128 (see Fig. 6).A similar rapid convergence of the distribution of |z 1 | was also reported in Ref. [94].In Sec V, we also compare the random-matrix behavior of the hardedge statistics with the distributions of physical models, which further provides evidence of the universality.
In the literature, the distributions of the eigenvalue z 1 with the smallest modulus were mainly used to characterize the hard-edge statistics (see Ref. [122] for the Hermitian case and Ref. [94] for the non-Hermitian case, as well as the references therein).Our newly proposed distributions of the complex level ratios in Eq. ( 11) differ from those of z 1 in several aspects.The distribution of re iθ does not only depend on the modulus r but also on the phase θ and hence is more informative than that of z 1 .Studying p r (r) and p θ (θ) does not require normalization while studying the distribution of z 1 requires normalization, for example, by the mean value of |z 1 |.Furthermore, the distribution of re iθ can be more effective in distinguishing between integrable and non-integrable models, as demonstrated in Sec.V D.

IV. NON-HERMITIAN RANDOM MATRICES BEYOND THE AZ0 AND AZ † CLASSIFICATION
We consider non-Hermitian random matrices in the 21 symmetry classes beyond the AZ 0 and AZ † classification.Similar to the previous analysis, we take the Gaussian ensemble in Eq. ( 8).The relevant symmetries of these classes are summarized in Table III.We first study classes BDI, CI, and AII + S + as representative classes (see Sec. IV A), and later discuss all the other classes in a similar manner (see Sec. IV B and Appendix C).In Appendix C, we provide additional characterization of the hard-edge statistics for these general symmetry classes.

A. Classes BDI, CI, and AII + S+
Non-Hermitian random matrices in class BDI (CI) respect TRS with sign +1 and CS commuting (anticommuting) with TRS.On the other hand, non-Hermitian random matrices in class AII + S + respect TRS with sign −1 and SLS commuting with TRS.See Table III, as well as Appendix A for explicit forms of non-Hermitian random matrices in each symmetry class.We numerically diagonalize 10 7 samples of N × N (N = 256) non-Hermitian random matrices in each of the three symmetry classes.In these symmetry classes, the combination of symmetries leads to quartets (z, −z, z * , −z * ) of complex eigenvalues, and hence the spectrum exhibits D 2 symmetry [see Fig. 1 (b)].In all the three symmetry classes, the eigenvalues are distributed almost uniformly in a circle in the complex plane, except around its origin, its circumference, and the real and imaginary axes.This observation is consistent with the symmetries of these random matrices since the spectral origin, as well as the real and imaginary axes, possesses higher symmetries than generic points in the bulk of the spectrum.In classes BDI and CI, a finite number N real/imag of eigenvalues are real or purely imaginary.Here, N real/imag scales as √ N after the ensemble averaging (not shown here), ex-hibiting subextensive behavior similar to that of random matrices in the AZ † classification (Table I).In class AII + S + , by contrast, all eigenvalues are complex [Figs.7 (d)-(f)].We find that the presence of real or purely imaginary eigenvalues is determined by the relevant symmetries on the real and imaginary axes.
To characterize the level statistics around the spectral origin, we first consider the distribution p min (z min ) of the eigenvalue z min closest to the origin.To make z min unique, we consider the distribution for Re z min ≥ 0 and Im z min ≥ 0. We normalize z min such that ⟨|z min |⟩ = 1.Because of the possible presence of real and purely imaginary eigenvalues, the probability density p min (z min = x + iy) of z min is decomposed as p min (z min = x + iy) = p min;r (x)δ(y) + p min;i (y)δ(x) + p min;c (x, y) , (15) where δ(•) is the Dirac delta function, and p min;r (x), p min;i (y), and p min;c (x, y) represent the distributions of z min when it is real, purely imaginary, and complex, respectively.The distribution p min (z min ) is normalized by and  (20) for |z min | ≪ 1.The joint probability distribution of complex eigenvalues for non-Hermitian random matrices in class AII + S + with arbitrary N was obtained previously [123].In Appendix B, we verify that 4 × 4 random matrices capture the same small-|z min | behavior as that of large random matrices.
The distribution of complex z min is qualitatively different from the distribution p min;r/i (|z min |) of real or purely imaginary z min (if exists) even in the same symmetry class.In fact, in class BDI, we have [Fig.7 (a)] In class CI, we have [Figs.7 (b)] They are significantly distinct from the small-|z min | behavior of p (r) min;c (|z min |).We provide an intuitive explanation for the stronger level repulsion for complex z min .As discussed before, all of ±z min and ±z * min are eigenvalues.For complex z min , these four eigenvalues differ from each other and share the smallest modulus.By contrast, for real z min (z min = z * min ) or purely imaginary z min (z min = −z * min ), only the two different eigenvalues share the smallest modulus.For complex z min with |z min | ≪ 1, there are at least three eigenvalues (i.e., −z min and ±z * min ) whose distance to z min is smaller than 2|z min | and hence much smaller than unity.In contrast, for real or purely imaginary z min with |z min | ≪ 1, only at least one eigenvalue satisfies this requirement.This difference leads to the distinct strength of level repulsion.In addition, we analytically calculate the hard-edge statistics for 4 × 4 non-Hermitian random matrices in class BDI in Appendix B, which also show distinct small-|z min | behaviors for complex, real, and purely imaginary z min , consistent with those of large random matrices.
Furthermore, the angular distribution  I [55,110,112].In both classes BDI and CI, non-zero real eigenvalues respect TRS with sign +1 while purely imaginary eigenvalues respect CS, equivalent to pH.The exponents in the power-law decay of p (θ) min;c align with those of ρ c (x + iy).The difference between classes BDI and CI lies in the signs of PHS, not influencing symmetries of non-zero real or purely imaginary eigenvalues.
In class AII + S + , non-zero real eigenvalues respect TRS with sign −1 while non-zero purely imaginary eigenvalues respect PHS † with sign −1, equivalent to TRS with sign −1 [51].The exponent in the power-law de-cay of p (θ) min;c also aligns with that of ρ c (x, y).It is also notable that p (θ) min;c (θ z ) in class AII + S + exhibits higher symmetry than those in classes BDI and CI.Because of the equivalence between TRS and PHS † , if H belongs to class AII + S + , iH also belongs to class AII + S + , both of which occur with the same probability in the Gaussian ensemble.Consequently, the probability of z being the eigenvalue with the smallest modulus is the same as that of iz.In combination with the spectral D 2 symmetry due to TRS and SLS, the distribution of z min should thus exhibit the point group symmetry D 4 after the ensemble average.This spectral D 4 symmetry requires p TABLE III.Twenty-one symmetry classes of non-Hermitian random matrices beyond the AZ † and AZ0 classification.If H belongs to the symmetry class in the column "Class (H)", iH belongs to the corresponding class in "Class (iH)".If H and iH belong to the same symmetry class, the class is in bold.In the nomenclature of symmetry classes, S represents sublattice symmetry (SLS) and η represents pseudo-Hermiticity (pH).The subscript of S± or η± denotes whether SLS or pH commutes (+) or anti-commutes (−) with time-reveal symmetry (TRS), particle-hole symmetry (PHS), or chiral symmetry (CS).The first and second subscripts of S±± denote the commutation (+) or anti-commutation (−) relations between SLS and TRS, and between SLS and PHS, respectively.In the columns "TRS", "TRS † ", "PHS", and "PHS † ", the entry "±1" specifies the sign of the anti-unitary symmetry, and "0" means its absence.The column "[Γ, S]±=0" specifies whether CS (Γ) and SLS (S) commute (+) or anti-commute (−) with each other if both of them exist.The symmetry classes are categorized into different groups by the number (N = 0, 2, 4) of anti-unitary symmetries.zmin denotes the eigenvalue with the smallest modulus.The level statistics around the spectral origin are characterized by ⟨|zmin| 2 ⟩ and the probabilities of zmin being real or purely imaginary.If these probabilities equal each other within the range of statistical errors, they are in bold.The statistics are obtained by diagonalizing 10 7 samples of 256 × 256 random matrices for classes BDI, CI, and AII + S+ and diagonalizing 10 6 samples of 256 × 256 random matrices for the other symmetry classes.
Beyond the AZ † and AZ 0 classification, we now investigate the remaining 18 symmetry classes.As indicators of the universal classes of the level statistics around the spectral origin, we use ⟨|z min | 2 ⟩ and the probability of z min being real or purely imaginary, denoted by Pr(z min ∈ R) in Eq. (17) and Pr(iz min ∈ R) in Eq. ( 18), respectively.Since these indicators already reach the convergence for 256 × 256 non-Hermitian random matrices in the three representative classes (see Fig. 8), we obtain these values in all the remaining 18 classes through the exact diagonalization of random matrices with the same size, as summarized in Table III.We expect these values to be universal and characterize the universality classes of the complex spectral statistics.For the three representative classes, we shortly verify the universality also by comparison with physical models (see Table V in Sec.VI for details); the exhaustive verification for physical models in the other symmetry classes is left for future study.
In each symmetry class, the distributions of z min exhibit characteristic behaviors.Within a subset of symmetry classes, the probability of z min being real or purely imaginary is non-zero, similar to classes BDI and CI.We find that whether Pr(z min ∈ R) is zero depends only on symmetries respected by the real axis (i.e., TRS, TRS † , and pH), regardless of the other symmetries.These sym- metries, comprising the AZ † classification, also determine whether non-Hermitian random matrices in the AZ † classification have real eigenvalues (see Sec. II and a summary in Table I).For example, in the AZ † classification, a random matrix in class AI + η − has TRS with sign +1 and TRS † with sign −1, leading to a subextensive number of real eigenvalues (see Table I).Correspondingly, classes BDI + S −− , CI + S ++ , BDI + S +− , and CI + S −+ have the same signs of TRS and TRS † , and non-Hermitian random matrices in these symmetry classes exhibit Pr(z min ∈ R) ̸ = 0 (see Table III).Similarly, whether Pr(iz min ∈ R) is zero depends only on the symmetries respected by pure imaginary eigenvalues (i.e., PHS † , TRS † , and CS).
Another notable feature is the spectral D 4 symmetry of the level statistics in some symmetry classes where H and iH belong to the same symmetry class.Such symmetry classes are shown in bold font in Table III.When Pr(z min ∈ R) and Pr(iz min ∈ R) overlap within the range of statistical errors, they are also shown in bold font in Table III.These two bold-font items consistently appear in the same rows in Table III due to the spectral D 4 symmetry of the distributions of z min , analogous to our discussion for class AII + S + .
In the remaining 18 symmetry classes, we also find the distinct small-|z min | behavior between the distribution of complex z min and that of real or purely imaginary z min .Moreover, the angular distribution p Thereby, our discussions on the level repulsion in the three representative classes can also be applied to the remaining 18 symmetry classes.A more detailed numerical study of the asymptotic behavior of the distributions, in a similar manner to Eqs. ( 20)- (25), and their theoretical explanations are left for future study.

V. OPEN QUANTUM SYSTEMS IN THE AZ0 CLASSIFICATION
In quantum chaotic systems that respect PHS or SLS and are isolated from the environment, the spectral statistics around zero energy E = 0 can be described by Hermitian random matrices in the same symmetry class (see, for example, Ref. [14] and references therein).Here, we calculate the density of complex eigenvalues and the distribution of the complex level ratios in Eq. ( 11) around the spectral origin of open quantum systems that respect PHS or SLS.We find that they are consistent with the statistics of non-Hermitian random matrices studied in Sec.III.As a comparison, we also study the level statistics of an integrable non-Hermitian many-body Hamiltonian with PHS.

A. Lindblad equation
When a quantum system is coupled to a Markovian environment, the time evolution of its reduced density matrix ρ is described by the Lindblad equation [30][31][32] Here, H is the Hamiltonian of the system that governs the unitary dynamics of ρ, and L µ 's represent dissipators originating from the coupling with the environment.The Lindbladian L is a superoperator that acts on the density matrix.Through the vectorization ρ = ρ ij |i⟩ ⟨j| → |ρ⟩ = ρ ij |i⟩ |j⟩, the density matrix can be mapped to a state in the double Hilbert space.Correspondingly, L becomes a non-Hermitian operator in the double Hilbert space, given as For a generic system with N complex fermions, the dimension of the double Hilbert space is (2 N ) 2 = 4 N .Therefore, L is a 4 N × 4 N operator, and numerical evaluation of its spectrum for larger N costs computation time.
Notably, quadratic Lindbladians, comprising quadratic Hamiltonians and linear dissipators, effectively reduce to single-particle non-Hermitian Hamiltonians [124,125], which facilitates numerical investigations of relatively larger systems.A non-interacting fermionic Hamiltonian and linear dissipators are generally given by in the Nambu basis Ψ . By construction, the 2N × 2N Hermitian matrix H respects PHS H = −τ x H T τ x , where τ x is the Pauli matrix acting on the (c, c † ) space.We define a 2N × 2N matrix by and its particle-hole-symmetric combination Using H and Γ, we can construct a 2N × 2N singleparticle non-Hermitian Hamiltonian From complex eigenvalues λ i 's of this effective singleparticle non-Hermitian Hamiltonian K, the 4 N eigenvalues of L are given as with n i = 0 or 1. Notably, K always satisfies PHS † , In this work, we consider H with the conserved particle number, given by The coupling with the environment is assumed to manifest itself in single-particle loss; the dissipators L µ consist only of annihilation operators, In such a case, the effective single-particle non-Hermitian Hamiltonian K is block-diagonalized as While PHS † in Eq. ( 33) relates the two reduced blocks h−il T l * and −h T −il † l = −(h−il T l * ) * , the reduced block h − il T l * itself does not necessarily respect PHS † .The symmetry of quadratic Lindbladians is determined solely by h − il T l * without intrinsic PHS † .Since h is Hermitian and l T l * is positive semidefinite, complex eigenvalues z i 's of h − il T l * must satisfy the constraint It was argued that quadratic Lindbladians cannot respect TRS, PHS, pH, or SLS since these symmetries violate this non-positivity constraint [95].As discussed previously, TRS or pH (PHS or SLS) requires complex eigenvalues to appear in (z, z * ) [(z, −z)] pairs, which is incompatible with Eq. ( 37).However, if we shift the Lindbladian by its trace, the traceless part is free from such constraints [97,98,126,127].In principle, the diagonal block of the traceless part, given as can respect any of the seven symmetries in Eqs. ( 1)- (7).
Hereafter, we investigate the symmetry classes and the hard-edge statistics of the shifted quadratic Lindbladian, represented by the effective non-Hermitian Hamiltonian h eff .

B. Quadratic Lindbladians
We study three prototypical models of quadratic Lindbladians in the AZ 0 classification.We consider the following non-interacting Hermitian Hamiltonians on the three-dimensional cubic lattice i ∈ Z 3 in the chiral symmetry classes (i.e., classes AIII, BDI, and CII): where the bracket ⟨i, j⟩ represents the nearest neighboring sites.Each site contains only one orbital but can accompany the spin degree of freedom.We assume the periodic boundary conditions and the system size of L 3 .For any hopping amplitude t i,j and L ∈ 2Z, the Hamiltonian respects SLS with S i,j ≡ δ i,j (−1) (ix+iy+iz) .Depending on the choices of t i,j , H realizes models in all three chiral symmetry classes, as follows.
Let us first consider cases without the spin degree of freedom: (i) Class AIII.-Choose t i,j = e iθ i,j , where θ i,j ∈ R is distributed uniformly and independently in [0, 2π).The Hamiltonian realizes the random flux model [128,129] and hence belongs to class CII.
These models were originally proposed for solid-state materials but also found realizations in synthetic materials such as cold atoms and photonic systems [131][132][133][134].The coupling between cold atoms can be designed to respect CS [132], and the hopping t i,j can be controlled by synthetic gauge fields [135,136].
We introduce the single-particle loss by dissipators that consist of a linear combination of annihilation operators [48,51,137,138].For each Hamiltonian H, we choose the following linear dissipators such that the relevant symmetries are respected: (i) Class AIII † .-ForH AIII , we introduce the dissipators as (ii) Class BDI 0 .-ForH BDI , we introduce the dissipators as (iii) Class CII 0 .-ForH CII , we introduce the dissipators as In all the cases, γ ≥ 0 denotes the dissipative coupling strength between the system and environment.The dissipative coupling is assumed to be identical among different sites and maintain translation invariance.The effective non-Hermitian Hamiltonians [Eq.( 38)] associated with these dissipators are denoted by h AIII † eff , h BDI0 eff , and h CII0 eff , respectively.Even in the presence of dissipation, they still respect SLS and TRS † with the same symmetry operators and hence belong to Classes AIII † , BDI 0 , and CII 0 .The explicit forms of the Hamiltonians and their symmetries can be found in Appendix D.
We investigate these effective non-Hermitian Hamiltonians numerically (see the parameters in Appendix D).For each class, we diagonalize 2 × 10 5 samples in the disorder ensemble.We study the statistics of single-particle eigenvalues of the effective non-Hermitian Hamiltonians, rather than many-body eigenvalues of the Lindbladians.Figures 9 (a), (c), and (e) show the distributions of complex eigenvalues for a single realization of each Hamiltonian.The full complex spectra deviate from the circular law of non-Hermitian random matrices [107,139] and are generally non-universal.The range of each spectrum along the imaginary-axis direction is much smaller than that along the real-axis direction because of the small dissipation strength γ.
We investigate the spectral properties around the origin and demonstrate that they conform to the universal hard-edge statistics determined solely by symmetry.We normalize the eigenvalues of the Hamiltonians such that ⟨|z 1 |⟩ = 1, where z 1 denotes the complex eigenvalue with the smallest modulus.After the ensemble average, the density of complex eigenvalues around the spectral origin does not depend on the phases but only on the modulus, and hence exhibits spectral U( 1 IV, we compare ⟨r⟩ and ⟨cos θ⟩ of the effective non-Hermitian Hamiltonians (and also those studied in Sec.V C) and non-Hermitian random matrices.In each symmetry class, the corresponding statistical quantities overlap with each other within two standard deviations.Our results demonstrate the universality of the hard-edge statistics of non-Hermitian random matrices and reveal the chaotic nature of the quadratic Lindbladians.
As also discussed above, the complex spectra around the origin exhibit U(1) rotation symmetry even for the physical models.Notably, this spectral U(1) symmetry of the ensemble-average density does not arise from the original symmetry of the disorder ensemble of the Hamiltonians but from the chaotic behavior.In fact, unlike the Gaussian ensemble of random matrices, e iϕ h X eff (X = AIII † , BDI 0 , CII 0 ) does not belong to the same disorder ensemble of h X eff for generic ϕ ∈ R. Nevertheless, its hard-edge statistics can still be described by random matrices, resulting in the emergence of the spectral U(1) symmetry.Notably, the spectra away from the origin, es-pecially spectral boundaries [see Fig. 9 (a), (c), and (e)], do not respect the U(1) symmetry.Thus, in such systems, an energy threshold E c ≥ 0 should exist, akin to the Thouless energy [140], such that the spectra are well described by the random matrices for |z| ≤ E c .In the study of closed quantum systems, the Thouless energy was shown to satisfy E c ∝ L −2 [141], exhibiting the same order as the mean level spacing of two-dimensional systems.In our open quantum systems studied here, E c is at least larger than 2⟨|z 1 |⟩, while an accurate estimate of E c is left for future study.It should also be noted that the emergent spectral U(1) symmetry, as well as the universality of the hard-edge statistics, is exact only in the infinite-size limit L → ∞.
The eigenstates around the spectral origin of the effective non-Hermitian Hamiltonians h eff correspond to the decay modes in the original Lindbladians.Before the spectral shift in Eq. ( 38), these eigenvalues have non-zero imaginary parts −2γ.In our models, we have γ = 0.1 for classes AIII † and BDI 0 , and γ = 0.4 for class CII 0 .Owing to the weak dissipation, 2γ is much smaller than the bandwidth ≈ 10 (i.e., range of the spectrum along the real axis; see Fig. 9), and the decay time of these eigenstates is estimated as 1/(2γ) (the reduced Planck constant is set to ℏ = 1).Thus, this decay should only have a small effect on the short-time dynamics in the time regime t ≲ 1/(2γ).It is also notable that physically relevant initial states can have a large overlap with the eigenstates around the spectral origin of the shifted Lindbladians [127].In such a case, the short-time dynamics should be predominantly determined by the level statistics and the properties of eigenstates around the origin.

C. Non-Hermitian superconductors
Superconductors with energy gain and loss are one of the most intensively studied systems in the field of non-Hermitian physics [46,101].Superconductors attached to leads or reservoirs are also effectively described by non-Hermitian Hamiltonians [99,100,103].In the following, we study four prototypical models of non-Hermitian superconductors that belong to classes D, C, CI 0 , and DIII 0 in the AZ 0 classification and demonstrate its dissipative quantum chaos described by non-Hermitian random matrices with symmetry.
(i) For symmetry class D, we introduce an energy gain and loss into a two-dimensional p-wave superconductor [142], Here, ϵ r i and ϵ i i are random chemical potential and random energy gain or loss that are distributed uniformly and independently in [−W 1 /2, W 1 /2] and [−W 2 /2, W 2 /2], respectively.Additionally, te iϕ and t represent the hopping amplitudes, and ∆ and i∆ are the p-wave superconducting paring potentials.In the Nambu basis where τ x is the Pauli matrix acting on the particle-hole space.Hence, H D belongs to class D (see Appendix D for further details).
(ii) For symmetry class C, we introduce a random imbalanced s-wave pairing potential into a threedimensional superconductor with conserved spin, with the Nambu basis Ψ i = (c ↑ , c † ↓ ).Here, ϵ i is a random potential distributed uniformly and independently in [−W 1 /2, W 1 /2].t and ϕ characterize the hopping amplitude and phase, respectively, ∆ ∈ R represents a Hermitian pairing potential, and ∆ s i ∈ R represents the random imbalanced s-wave paring potential [143][144][145][146][147] and hence belongs to class C. (iii) For symmetry class CI 0 , we consider a non-Hermitian superconductor similar to H C in Eq. ( 45) but with a different pairing potential.We denote its Hamiltonian as H CI0 ; for further details, refer to Appendix D.
(iv) For symmetry class DIII 0 , we study a non-Hermitian extension of a three-dimensional timereversal-invariant superconductor with the spin-orbit coupling [148,149].We denote its Hamiltonian as H DIII0 ; additional information is provided in Appendix D.
We study single-particle eigenvalues of the BdG Hamiltonians of these four non-Hermitian superconductors numerically with the parameters provided in Appendix D. For each Hamiltonian, we diagonalize 2 × 10 5 disorder realizations.Figures 11 (a), (c), (e), and (g) show the spectral distributions of single realizations.The ensemble average of the density of complex eigenvalues around the origin exhibits emergent spectral U(1) rotation symmetry [Figs.11 (b), (d), (f), and (h)].Importantly, both the density of the modulus (after the normalization such that ⟨|z 1 |⟩ = 1) and the distributions of the complex level ratio z 1 /z 2 in Eq. ( 11) match well with those of non-Hermitian random matrices in the same symmetry class, as shown in Fig. 12 and Table IV.These results further demonstrate that the hard-edge statistics of non-Hermitian random matrices well describe the universal level statistics of open quantum chaotic systems.
The random-matrix behavior of the hard-edge statistics in these models suggests the delocalization of the eigenstates around the spectral origin.Disorder can drive quadratic non-Hermitian Hamiltonians into the Anderson-localized phases [53][54][55][56][57][58][59][60][61][62][63][64][65][66], where eigenvalues should be uncorrelated.For the quadratic Lindbladians h X eff (X = AIII † , BDI 0 , CII 0 ) in Sec.V B, disorder manifests not in on-site potential but in random hopping, resulting in SLS and the concomitant absence of localization [15].By contrast, the non-Hermitian superconductors studied in this subsection include on-site disorder.Hence, we study their weak disorder regimes and anticipate that they have yet to undergo the Anderson localization, supported by their hard-edge statistics.Notably, the finite-size scaling of the hard-edge statistics can capture possible Anderson transitions for stronger disorder.
Additionally, symmetry plays a fundamental role in spectral statistics and Anderson localization.In the models studied in this work, the spectral origin exhibits higher symmetry than generic points in the bulk of spectra.Thus, it is plausible that eigenstates around the spectral origin are extended while those away from the origin are localized, or vice versa.A notable Hermitian example of the former case is found in two-dimensional disordered superconductors in class D [15,150], where the spectral statistics exhibit the random-matrix behavior around the origin but become uncorrelated away from the origin.We still consider such models with extended states only around the origin as chaotic, because the extended states may contribute significantly to chaotic behavior in quantum transport phenomena.

D. Integrable model: non-Hermitian SYK2 model
In the preceding sections, we classify the universal hard-edge statistics in various open quantum chaotic models.In this section, on the other hand, we study an open quantum integrable model and demonstrate that the hard-edge statistics provide a diagnosis for the integrability of non-Hermitian many-body operators.
It was found that H q=2 SYK as a many-body Hamiltonian respects PHS with sign +1 for N ∈ 8Z [94].Additionally, its level statistics in the bulk of the spectrum (i.e., away from the origin or the edge of the spectrum) were found to be accurately described by the Poisson statistics of complex numbers (or equivalently, two-dimensional Poisson statistics).In the presence of integrability, many-body PHS does not necessarily imply that the model belongs to symmetry class D. Generally, the integrability implies a block-diagonalized structure of many-body Hamiltonians, and their relevant symmetry classes should be determined by symmetries of the reduced blocks.Conversely, the hard-edge statistics can also be used as a diagnostic tool for detecting the integrability of many-body Hamiltonians with symmetries.
To demonstrate this in the non-Hermitian SYK 2 model, we numerically diagonalize 2 × 10 5 realizations of H q=2 SYK with N = 24, and compare the statistics of its many-body eigenvalues z i 's [Eq.( 48)] around the spectral origin with the hard-edge statistics from non-Hermitian random matrices in class D (Fig. 13).The radial distribution p r (r) and the angular distribution of the complex level ratio re iθ ≡ z 1 /z 2 are consistent with the complex Poisson statistics in Eq. ( 14) [see Figs.
implying that the density of complex eigenvalues, given by ρ(|z|)/π|z|, is not higher but almost constant around the spectral origin.All these distributions are significantly different from the random-matrix behavior in class D, providing clear evidence of the integrability of the model.We also study the statistics of the eigenvalue z min ≡ z 1 with the smallest modulus [94,122].Its radial distribution satisfies [Fig.13 (d)] The distribution conforms to the complex Poisson statistics and deviates from the random-matrix statistics in class D, indicating that p min (|z min |) also serves as a tool to detect quantum chaos or its absence.However, p min (|z min |) of non-Hermitian random matrices in class D is close to the complex Poisson statistics [Fig.13 (d)].Thus, a practical use of p min (|z min |) for detecting quantum chaos is less efficient than the level ratio statistics, especially for many-body models, where numerical calculations are constrained by limited system size and finitesize effects are substantial.

VI. OPEN QUANTUM SYSTEMS BEYOND THE AZ0 AND AZ † CLASSIFICATION
In this section, we study the hard-edge statistics of open quantum systems beyond the AZ 0 classification.In Sec.VI A, we study the single-particle spectra of noninteracting systems.In Sec.VI B, by contrast, we consider the many-body spectra of interacting Lindbladians.Specifically, we focus on the distributions of the complex eigenvalue min closest to the spectral origin in the three representative symmetry classes, i.e., classes BDI, CI, and AII + S + .We that the numerically obtained distributions of the physical models match well the distributions of non-Hermitian random matrices in symmetry classes.Notably, in classes BDI and CI, the distributions of z show deltafunction peaks on the real and imaginary axes, consistent with the random-matrix behavior.

A. Non-interacting systems
We introduce three non-interacting models in classes BDI, CI, and AII + S + , respectively.For class BDI, we add the following dissipators to the random-hopping model H BDI [i.e., Eq. ( 39) with t i,j ∈ R].For class AII + S + , we add the following dissipators to the SU(2) model H CII [i.e., Eq. ( 39) with t i,j ∈ SU( 2 SYK in Eq. ( 47).(a) Density of eigenvalues in the complex plane, given by ρ(|z|)/π|z|.(b) Radial distributions pr(r) and (c) angular distributions p θ (θ) of the complex level ratios re iθ [Eq.(11)] around the spectral origin.(d) Radial distribution of the eigenvalue zmin with the smallest modulus.The points with the error bars refer to the non-Hermitian SYK model with q = 2, and the solid curves are obtained from non-Hermitian random matrices (RM) in class D. As a non-Hermitian model in class CI, we investigate a three-dimensional non-Hermitian superconductor sim-ilar to Eq. ( 45) but with a different form of imbalanced pairing potentials.We refer to its Hamiltonian as H CI , and the details can be found in Appendix D.
For each symmetry class, we numerically diagonalize 2 × 10 5 samples of Hamiltonians with the parameters provided in Appendix D. For the models in classes BDI and CI, a subextensive number of eigenvalues are real or purely imaginary [Figs.14 (a) and (b)].The probability of z min being real or purely imaginary coincides with that of non-Hermitian random matrices in the same symmetry classes (see Table V).For the model in class AII + S + , by contrast, all eigenvalues are complex, also consistent with the random-matrix behavior [Fig.14 (c)].Furthermore, the distributions of z min , including the radial distributions of real, purely imaginary, and complex z min , and the angular distributions of complex z min , in all these models are well described by non-Hermitian random matrices [Figs.14 (d)-(f)].In Table V, we compare ⟨|z min | 2 ⟩ from non-Hermitian random matrices and that from the non-Hermitian Hamiltonians.The values are consistent among different systems in the same symmetry class.

B. Many-body Lindbladians
As a prototypical physical model with strongly correlated fermions coupled to the environment, we study the Lindbladian L for the SYK model [98,[155][156][157].The Hamiltonian is given by the original SYK model consisting of the q-body random interactions, where γ i 's are Majorana operators, N is the number of Majorana fermions, and J i1i2•••iq 's are independent real Gaussian variables with zero means and ⟨J 2 i1i2•••iq ⟩ = ((q − 1)!/N q ) J 2 .Then, we add the following linear dissipators L m , where K m;i 's are independent complex Gaussian variables with zero means and As discussed in Sec.V A, after the vectorization of the density matrix, the Lindbladian L becomes a non-Hermitian operator in the double Hilbert space.It can be represented with the basis of adjoint fermions [84] (see also Appendix E for details), facilitating analytic analyses.The SYK Lindbladian L respects a unitary symmetry where (−1) F is the total fermion parity operator in the double Hilbert space.For q = 6, the traceless Lindbladian L − TrL/Tr1 belongs to symmetry class BDI for N ∈ 4Z and class CI for N ∈ 4Z + 2 in each sector of total fermion parity (−1) F = ±1 (see Ref. [98] and Appendix E for derivations).We numerically diagonalize 10 5 realizations of the traceless Lindbladians with the following two sets of parameters: (i) Class BDI.-N = M = 12, J = 1, and K = 0.4; (ii) Class CI.-N = M = 10, J = 1, and K = 0.3.The SYK Lindbladians exhibit distinct distributions of z min for different N [Figs.15 (b), (d)], consistent with non-Hermitian random matrices in the corresponding symmetry classes (see Table V).As also discussed in Sec.V B, the hard-edge statistics should be relevant to the short-time dynamics of the SYK Lindbladians.

VII. CONCLUSIONS AND DISCUSSIONS
In 28 out of 38 symmetry classes of non-Hermitian random matrices, we find that the spectral origin in the complex plane respects higher symmetry compared to other points.We comprehensively investigate the hardedge statistics in these 28 symmetry classes.This work, combined with the prior works on the threefold universal statistics in the bulk [77] and the sevenfold universal statistics on the real axis [112], forms a foundational comprehension of non-Hermitian random matrix theory across all the 38 symmetry classes.In the seven symmetry classes of the AZ 0 classification, we introduce complex level ratios around the spectral origin and analyze the universality classes of the hard-edge statistics.We numerically obtain the seven universal ratio distributions, with both radial and angular components exhibiting the characteristic behavior.In 21 symmetry classes beyond the AZ 0 classification, we utilize the distributions of complex eigenvalues z min with the smallest modulus to characterize the hard-edge statistics.In certain symmetry classes, z min can be real or purely imaginary with the universal probability.This probability serves as practically useful indicators for determining universality classes.Additionally, we perform analytic calculations for small random matrices and elucidate the essential features of our numerical results.
The results of this paper also serve as a diagnostic tool of non-integrability or its absence in open quantum systems with symmetry, as well as the transitions between them.We investigate diverse physical models in the seven AZ 0 symmetry classes and the three representative symmetry classes beyond the AZ 0 classification.In contrast to dense random matrices in the Gaussian ensembles, ensembles of physical Hamiltonians or Lindbladians are regarded as ensembles of sparse matrices because of the locality constraint.Nevertheless, the random-matrix hard-edge statistics emerge in both physical Hamiltonians and Lindbladians, showcasing their universality.Notably, in the presence of symmetry, the spectral statistics around the origin can exhibit distinct behavior compared to that in the bulk.In such systems, chaos and integrability cannot be diagnosed solely by the bulk statistics, and thus both bulk and hard-edge statistics are significant.Prime examples in closed quantum systems include twodimensional disordered superconductors, in which only zero modes can be free from the Anderson localization and exhibit quantum chaotic behavior [15,150].Examining the spectral statistics in both the bulk and around the origin also provides insights into fundamental questions regarding the impact of symmetry on the universality classes of chaos-integrability transitions.
Symmetry of ensembles of physical systems is generally lower than symmetry of the Gaussian random-matrix ensembles.However, the hard-edge statistics of physical systems conform to the random-matrix behavior, which implies the emergent spectral U(1) or D 4 symmetry of the ensemble-average distributions of the complex eigenvalues around the spectral origin.It is worth mentioning that the ensemble-average distributions of real eigenvalues of chaotic Hermitian physical systems in the standard symmetry classes (i.e., classes A, AI, and AII) should also be symmetric about the origin within several mean level spacings, provided that the origin is away from the edge of the spectra.Since generic Hamiltonian ensembles of Hermitian systems do not have such spectral inversion symmetry, this implies emergent Z 2 symmetry.Nonetheless, the spectral origin of Hermitian systems in the stan-dard symmetry classes does not possess any higher symmetry than generic eigenvalues, and the Z 2 symmetry arises because the density of eigenvalues should be almost constant within the range of several mean level spacings.By contrast, the emergent symmetry in the non-Hermitian systems studied in this work is distinct and more intricate than the Hermitian case.The spectral origin of these non-Hermitian systems exhibits higher symmetry, and the density around the origin displays a non-constant but D 4 -symmetric (or U(1)-symmetric) distribution within the range of several mean level spacings.
It is interesting to explore alternative characterizations of hard-edge statistics and their applicability to physical systems.For example, in Appendix C, for non-Hermitian random matrices in classes BDI, CI, and AII + S + , we consider the distributions of complex, real, and purely imaginary eigenvalues, separately.Despite slower convergence with respect to the matrix size, these distributions also appear to be universal and align with physical Hamiltonians.Other potential characterizations, including the level number variance within a fixed range around the spectral origin and n-point correlation functions [158], are also left for future study.
The symmetry classification of open quantum systems is not limited to Hamiltonians or Lindbladians but also applies to dynamical generators of the non-Markovian dynamics [159].Investigating and comparing their hardedge statistics with our results can provide insights into these less-explored, but physically relevant, open quantum systems.The applications of non-Hermitian random matrix theory in classical systems, including complex networks with symmetry, also merit future studies.

(B6)
Here, an extra factor of two is due to the assumption |z 1 | ≤ |z 2 |.For r ≪ 1, we have ρ r (r) ∝ −r 3 ln r, which is consistent with the numerical results from large random matrices (see Table II).We can also integrate p joint (|y Since the phase θ of the level ratio z 1 /z 2 equals (ϕ 1 − ϕ 2 )/2, its distribution function is obtained as The angular distribution p N =4 θ reaches its maximum and minimum at θ = ±π/2 and θ = 0, respectively, consistent with the numerical results of large random matrices [see Fig. 10 (c with Notably, w 1 and w 2 are independent complex random variables with the identical probability distribution in the complex plane [77], The radial distribution of w is given as The radial distribution of t ≡ w 1 /w 2 is given as Thus, the probability density of t in the complex plane is The complex level ratio v ≡ re iθ ≡ z 1 /z 2 = (1+ t)/(1−t) is a meromorphic function of t.Then, the probability density of v in the complex plane is given as Integrating p v (v) over θ or r, we obtain the radial and angular distributions of the complex level ratio v as respectively.For r ≪ 1, we have p N =4 r (r) ∝ r, and p N =4 θ reaches its maximum and minimum at θ = ±π/2 and θ = 0, respectively, consistent with numerical results from large random matrices [see Table II and Fig. 10 (c)].for small |z| ≪ 1.This linear decay of ρ (|z|) does not necessarily mean the level repulsion around the origin since it just arises from the two-dimensional integral measure along the angle direction.A similar result was obtained in Ref. [94].
A generic 2 × 2 non-Hermitian matrix H in class BDI 0 (or equivalently, class D + S + ), which respects SLS σ z Hσ z = −H and PHS H T = −H, is again given as Eq.(B19).Consequently, the probability distribution function ρ (|z|) is the same as ρ (|z|) in class D [i.e., Eq. (B21)].However, this is specific to small non-Hermitian matrices.For larger non-Hermitian matrices, ρ (|z|)'s are different between classes D and BDI 0 owing to the many-level effect, as shown in Fig. 4.

b. Class C
A generic 2 × 2 non-Hermitian matrix H in class C, which respects PHS τ y H T τ y = −H, is given as The probability distribution function ρ (|z|) of |z| = a 2 x + a 2 y + a 2 z is equivalent to the level-spacing distribution of 2 × 2 non-Hermitian random matrices in class A, which is analytically obtained as [77] ρ (|z|) = 2C 4  3 |z| with the normalization constant Here, ρ (|z|) is normalized by The cubic decay of ρ (|z|) means the level repulsion around the origin.A similar result was obtained in Ref. [94].

c. Class CII0
No generic 2 × 2 non-Hermitian matrix H is present in class CII 0 (or equivalently, class C + S + ).A generic 4×4 non-Hermitian matrix H in class CII 0 , which respects TRS † σ y H T σ y = H and PHS τ y H T τ y = −H, is given as The eigenvalues are obtained as and exhibit the two-fold degeneracy due to TRS † .Thus, the probability distribution function ρ (|z|) is equivalent to the level-spacing distribution of 2 × 2 non-Hermitian random matrices in class AII † , which is analytically obtained as [77] ρ (|z|) = 2C A generic 2 × 2 non-Hermitian matrix H in class AIII † , which respects SLS σ z Hσ z = −H, is given as The two eigenvalues of H are obtained as Notably, the density of complex eigenvalues, ρ (|z|), is equivalent to the level-spacing statistics of 2 × 2 non-Hermitian matrices in class AI † , which are analytically obtained as [77] ρ where A similar result was obtained in Refs.[94,113].In the previous analyses [77,112], the logarithmic correction arises only in the presence of TRS † with the positive sign +1.By contrast, the logarithmic correction appears in the hard-edge statistics with SLS, even in the absence of TRS † .No generic 2 × 2 non-Hermitian matrix H is present in class DIII 0 (or equivalently, class D + S − ).A generic 4×4 non-Hermitian matrix H in class DIII 0 , which respects TRS † σ y H T σ y = H and SLS (σ z τ z ) H (σ z τ z ) = −H, is given as The modulus of the two eigenvalues is |z| = a 2 x + a 2 y , and hence the probability distribution function ρ (|z|) is equivalent to ρ (|z|) in class AIII † [i.e., Eq. (B32)].
A generic 2 × 2 non-Hermitian matrix H in class CI 0 (or equivalently, class C + S − ), which respects PHS τ y H T τ y = −H and TRS † H T = H, is given as Therefore, w 1 and w 2 are independent random variables with the same probability density function.Here, w 1 and w 2 are either real or purely imaginary, and we require Re w 1 ≥ 0 and Im w 2 ≥ 0. The probability density p w (w 1 = α + iβ) of w 1 in the complex plane is given as [112] p w (w (B38) According to whether w 1 and w 2 are real or purely imaginary, we have three possible cases, as follows.
(i) Both w 1 and w 2 are real.-Inthis case, all the four eigenvalues of H are real.Without loss of generality, we can assume w 1 ≥ w 2 ≥ 0.Then, we have z 2 = w 1 + w 2 > z 1 = w 1 − w 2 ≥ 0, and hence z 1 is the eigenvalue with the smallest modulus.The probability density p min;r (x) of z 1 = x ≥ 0 is given as For non-Hermitian random matrices in some symmetry classes beyond the AZ † and AZ 0 classification (see Table III), a subextensive number of eigenvalues are real or purely imaginary.We also find that in these symmetry classes, the level correlations of complex eigenvalues are qualitatively different from those of the real or purely imaginary eigenvalues (see Fig. 7).To characterize the level statistics in these symmetry classes, besides studying the eigenvalue z min with the smallest modulus (see Sec. IV for details), we here study complex, real, and purely imaginary eigenvalues, separately.Let c min , R min , and I min be the complex, real, and purely imaginary eigenvalues with the smallest modulus, respectively.Owing to the defining symmetries, we can require Re c min , Im c min > 0, R min ≥ 0, and I min ≥ 0 without loss of generality.As in Sec.IV A, we here focus on the three representative symmetry classes, i.e., classes BDI, CI, and AII + S + .
We first investigate the distribution of c min = |c min |e iϕ for non-Hermitian random matrices in classes BDI, CI, and AII + S + .We normalize c min such that ⟨|c min |⟩ = 1.In class AII + S + , all eigenvalues are complex, and hence we have c min = z min .For N ≥ 256, |c min |2 and ⟨cos ϕ⟩ converge to the characteristic values in each symmetry class (see Fig. 17 (C1) This is the same as the small-|z min | behavior of p min;c (|z min |) in the same symmetry classes (Fig. 7).The angular distributions p(ϕ) of c min are also similar to those of complex z min [compare Figs.18 (b), (d) with Fig. 14].
Then, we study the distributions of R min and I min in classes BDI and CI.In each symmetry class, the distributions p(R min ) and p(I min ) share the same small-R min or small-I min behavior as p min:r/i (|z min |).However, p(R min ) and p(I min ) are long-tailed in contrast to p min:r/i (|z min |).The statistics of I min converge much more slowly than the statistics of c min and do not converge even for N = 1000 [see Figs.19 (c) and (d)].Thus, p(|c min |) or p(R min ) should be a more computationally efficient measure of the level statistics than p(I min ).Furthermore, we compare p(|c min |) and p(ϕ) of open quantum systems studied in Sec.VI A with those of non-Hermitian random matrices.We observe that p(|c min |) of open quantum systems in different symmetry classes exhibit distinct behavior consistent with those of non-Hermitian random matrices (Fig. 20).The angular distributions p(ϕ) also match well with those of non-Hermitian random matrices (Fig. 20).These results further demonstrate that the level statistics of random matrices well describe universal properties of quantum chaotic systems with symmetry.In Sec.V B, we consider Hermitian Hamiltonians H AIII , H BDI , and H CII in the form ⟨i,j⟩ c † i t i,j c j + H.c. [Eq.(39)] and their dissipators.Their effective non-Hermitian Hamiltonians in Eq. ( 38) are given as respectively.These effective non-Hermitian Hamiltonians still respect SLS and TRS † with the same symmetry operators.In fact, we have , and h CII0 eff belong to symmetry classes AIII † , BDI 0 , and CII 0 , respectively.

Non-Hermitian superconductors
In Sec.V C, we consider non-Hermitian superconductors in classes D, C, CI 0 , and DIII 0 in the AZ 0 classification.Here, we provide some details of these models.
(i) For symmetry class D, the phase e iϕ of the hopping term in the Hamiltonian [Eq.(43)] can originate from a spatial modulation of the superconducting phase.In fact, with a gauge transformation c i → e iixϕ c i , the hopping phase can be eliminated [142], and the paring terms change as i∆c i+ex c i → i∆e i(2ixϕ+ϕ) c i+ex c i and ∆c i+ey c i → ∆e i(2ixϕ) c i+ey c i .In the Nambu basis Ψ i = (c i , c † i ) T , the Hamiltonian is written as the Nambu basis is given as where σ µ 's are Pauli matrices acting on spin space.In terms of the Nambu basis Ψ i = (c ↑ , c ↓ , c † ↑ , c † ↓ ) T and τ µ 's that act on its particle-hole space, the Hamiltonian reads, where ϵ r i and ϵ i i are a random chemical potential and random energy gain and loss that are distributed uniformly in [−W 1 /2, W 1 /2] and [−W 2 /2, W 2 /2], respectively.This Hamiltonian respects τ x H T DIII0 τ x = −H DIII0 , σ y H T DIII0 σ y = H DIII0 , (D15) and hence belongs to class DIII 0 .
For each Hamiltonian, we diagonalize 2 × 10 5 disorder realizations, with the parameters chosen as follows: Here, the lattice is either square lattice or cubic lattice with the linear dimension L. In Sec.VI A, for symmetry class CI, we consider a three-dimensional non-Hermitian superconductor with imbalanced pairing potentials.In the Nambu basis, the Hamiltonian H CI reads and hence belongs to class CI.We diagonalize 2 × 10 5 realizations with the parameters: t = 2, µ = 0, ∆ x = ∆ y = 1, ∆ z = 1/2, W 1 = 6, W 2 = 3, and L = 12.

FIG. 1 .
FIG. 1. Symmetry of complex spectra (denoted by the colored region) for non-Hermitian random matrices.(a) Seven symmetry classes in the Altland-Zirnbauer0 (AZ0) classification.The complex spectrum is invariant under reflection about the spectral origin but not under reflection about the real or imaginary axis, resulting in the point group symmetry C2.(b) Twenty-one symmetry classes beyond the AZ0 classification.The complex spectrum is invariant under reflections about both real and imaginary axes, indicating the point group symmetry D2.
4 (a)].Among the four BdG symmetry classes, for |z| ≪ 1, ρ(z) in class D is the largest, that in class C is the smallest, and those in classes CI 0 and DIII 0 are close to each other [see Figs. 4 (b) and (c)].Note that the density ρ(|z|) of the modulus differs from the density of eigenvalues in the complex plane by a Jacobian, and hence we have ρ(|z| = 0) = 0 even in the absence of the level repulsion for classes D and BDI 0 .Fitting ρ(|z|), we find the asymptotic behavior of ρ(|z|) for |z| ≪ 1 (also summarized in Table II), and BDI 0 ) ; −|z| 3 ln |z| (classes AIII † , CI 0 , and DIII 0 ) ; |z| 3

5 FIG. 5 .
FIG. 5. (a)Radial distributions pr(r) and (b) angular distributions p θ (θ) of the complex level ratio re iθ [Eq.(11)] around the spectral origin for the seven symmetry classes in the AZ0 classification.The error bars represent one standard deviation evaluated by the bootstrap method[116].

r
(r) and p N =4 θ (θ) for 4 × 4 non-Hermitian random matrices in classes AIII † and D, which show behavior consistent with large ones (see Appendix B).

FIG. 7 .
FIG. 7. (a)-(c) Distributions of the eigenvalue zmin with the smallest modulus of non-Hermitian random matrices in classes BDI, CI, and AII + S+.The blue, orange, and yellow curves represent the radial distributions of zmin when it is complex, real, and purely imaginary, respectively.Inset: angular distribution of complex zmin.The distributions are obtained by diagonalizing 10 7 samples of 256 × 256 random matrices in each symmetry class.(d)-(f) Eigenvalues of a single realization of a 256 × 256 random matrix in each symmetry class.Complex, real, and purely imaginary eigenvalues are denoted by the black, blue, and red points, respectively.(g) Small-|zmin| behavior of the radial distributions of complex zmin in different symmetry classes, with the solid curves being fitting curves.

5 FIG. 8 .
FIG. 8. (a) Probability of the eigenvalue zmin with the smallest modulus being real or purely imaginary as functions of the matrix size N in different symmetry classes.In class AII + S+, zmin is always complex.(b) ⟨|zmin| 2 ⟩ as functions of the matrix size N in different symmetry classes.The statistics are obtained by diagonalization of N sample matrices for each class and size with N sample = 10 7 , 10 6 , 2 × 10 5 for N = 256, N < 256, N > 256, respectively.The error bars represent one standard deviation.The dashed lines denote the mean values obtained by matrices with N = 256.
) rotation symmetry [Figs.9 (b), (d), and (f)].For the Hamiltonians in different symmetry classes, ρ(|z|)'s show distinct behavior consistent with the characteristic forms of non-Hermitian random matrices in the corresponding symmetry classes [Fig.10 (a)].The complex level ratios re iθ ≡ z 1 /z 2 in Eq. (11) around the spectral origin also obey the same distributions of non-Hermitian random matrices in the same symmetry class [Figs.10 (b) and (c)].To further validate this consistency, in Table

4 FIG. 9 .
FIG. 9. Complex eigenvalues of single realizations of the effective non-Hermitian Hamiltonians for the quadratic Lindbladians in classes (a) AIII † , (c) BDI0, and (e) CII0 (see Sec. V B). (b), (d), (f) Ensemble average of the density of complex eigenvalues around the spectral origin of the same Hamiltonians as the left subfigures.

FIG. 10 .
FIG. 10.Hard-edge statistics of the effective non-Hermitian Hamiltonians of the quadratic Lindbladians in classes AIII † , BDI0, and CII0 (see Sec. V B).(a) Density ρ(|z|) of the modulus of complex eigenvalues.(b) Radial distributions pr(r) and (c) angular distributions p θ (θ) of the complex level ratios re iθ [Eq.(11)] around the spectral origin.The colored points with the error bars refer to the Hamiltonians, and the black solid lines are obtained from non-Hermitian random matrices in the same symmetry classes.

3 FIG. 11 .
FIG. 11.Complex eigenvalues of single realizations of the Hamiltonians of non-Hermitian superconductors in classes (a) D, (c) C, (e) CI0, and (g) DIII0(see Sec.V C). (b), (d), (f), (h) Ensemble average of the density of complex eigenvalues around the spectral origin of the same Hamiltonians as the left subfigures.

FIG. 12 .
FIG. 12. Hard-edge statistics of the Hamiltonians of non-Hermitian superconductors in classes D, C, CI0, and DIII0 (see Sec. V C).(a) Density ρ(|z|) of the modulus of eigenvalues.(b) Radial distributions pr(r) and (c) angular distributions p θ (θ) of the complex level ratios re iθ [Eq.(11)] around the spectral origin.The colored points with the error bars refer to the Hamiltonians, and the black solid curves are obtained from non-Hermitian random matrices in the corresponding symmetry classes.

FIG. 13 .
FIG. 13.Hard-edge statistics of the two-body non-Hermitian SYK model H q=2SYK in Eq. (47).(a) Density of eigenvalues in the complex plane, given by ρ(|z|)/π|z|.(b) Radial distributions pr(r) and (c) angular distributions p θ (θ) of the complex level ratios re iθ [Eq.(11)] around the spectral origin.(d) Radial distribution of the eigenvalue zmin with the smallest modulus.The points with the error bars refer to the non-Hermitian SYK model with q = 2, and the solid curves are obtained from non-Hermitian random matrices (RM) in class D.

1 FIG. 14 .
FIG. 14. (a)-(c) Complex eigenvalues of single realizations of the non-Hermitian Hamiltonians (see Sec. VI A).Complex, real, and purely imaginary eigenvalues are denoted by the black, blue, and red points, respectively.(d)-(f) Distributions of the complex eigenvalue zmin with the smallest modulus.The blue, orange, and yellow points with the error bars represent the radial distributions of zmin when it is complex, real, and purely imaginary, respectively.Inset: angular distribution of complex zmin.The solid curves correspond to non-Hermitian random matrices (RM).(a), (d) Effective non-Hermitian Hamiltonian h BDI eff for class BDI.(b), (e) Hamiltonian HCI of the non-Hermitian superconductor for class CI. (c), (f) Effective non-Hermitian Hamiltonian h AII+S + eff for class AII + S+.

FIG. 15 .
FIG. 15.Complex eigenvalues in the symmetry sector (−1) F = −1 of single realizations of the many-body Lindbladians in (a) class BDI (N = 12) and (b) class CI (N = 10).Complex, real, and purely imaginary eigenvalues are denoted by the black, blue, and red points, respectively.The distributions of the complex eigenvalue zmin with the smallest modulus for the shifted many-body Lindbladians L − TrL/Tr1 in (b) class BDI and (d) class CI.The blue, orange, and yellow points with the error bars represent the radial distributions of zmin when it is complex, real, and purely imaginary, respectively.Inset: angular distribution of complex zmin.The solid curves correspond to non-Hermitian random matrices (RM).

1 FIG. 16 .
FIG.16.Distributions of the eigenvalue zmin with the smallest modulus of non-Hermitian random matrices in 18 out of 21 symmetry classes beyond the AZ † and AZ0 classification.The distributions in the three symmetry classes (i.e., classes BDI, CI, and AII + S+) not shown here are found in Fig.7.The blue, orange, and yellow curves represent the radial distributions of zmin when it is complex (C), real (R), and purely imaginary (I), respectively.In some symmetry classes (e.g., class AI + S+), the orange and yellow curves overlap almost perfectly, and the orange curve is hardly seen.Inset: angular distribution of complex zmin.The distributions are obtained by diagonalizing 10 6 samples of 256 × 256 random matrices in each symmetry class.

FIG. 17 .
FIG. 17. Mean values |cmin| 2 and ⟨cos ϕ⟩ as functions of the matrix size N , where cmin = |cmin|e iϕ is the complex eigenvalue with the smallest modulus.For each symmetry class and each matrix size, N sample samples are diagonalized, with N sample = 10 7 , 10 6 , 2 × 10 5 for N = 256, N < 256, N > 256, respectively.The dashed lines denote the mean values obtained by matrices with N = 256.

FIG. 19 .
FIG. 19.Distributions of the complex, real, and purely imaginary eigenvalues with the smallest modulus for 256×256 non-Hermitian random matrices in (a) class BDI and (b) class CI.For a small proportion p of random matrices, real eigenvalues (p ≈ 0.02%) or purely imaginary eigenvalues (p ≈ 3.5%) are absent.We exclude such matrices in the statistics of real or imaginary eigenvalues.(c) Mean values I 2 min and (d) probability that no purely imaginary eigenvalues exist as functions of the matrix size N .

FIG. 20
FIG. 20. (a), (c) Radial and (b), (d) angular distributions of the complex eigenvalue cmin with the smallest modulus for the non-Hermitian Hamiltonians and many-body Lindbladians.The points with the error bars correspond to the quantum systems, and the solid curves are obtained from the random matrices in the same symmetry classes.(a), (b) Effective non-Hermitian Hamiltonian in Eq. (D5) for class BDI and non-Hermitian superconductor in Eq. (D16) for class CI. (c), (d) Many-body Lindbladians of the SYK model in Eqs.(53) and(54) for classes BDI (N = 12) and CI (N = 10), where we use the eigenvalues in the subsector of (−1) F = −1 to obtain the statistics.

TABLE IV .
(11) values ⟨r⟩ and ⟨cos θ⟩ of non-Hermitian random matrices (RM) and Hamiltonians in the AZ0 classification (Secs.V B and V C), where re iθ denotes the complex level ratio [Eq.(11)] around the spectral origin.
and is denoted by H AIII .The Hamiltonian H AIII only respects SLS and hence belongs to class AIII.
Consider next the case with the spin degree of freedom: (iii) Class CII.-Choose 2×2 hopping matrices t i,j that are distributed uniformly and independently on the Harr measure of SU(2) matrices.The Hamiltonian realizes the SU(2) model [130] and is denoted by H CII .In addition to SLS, H CII respects TRS with sign −1 σ y H * CII σ y = H CII

TABLE V .
Mean value ⟨|zmin|⟩ and probability of zmin being real or purely imaginary for non-Hermitian random matrices, the Hamiltonians in Sec.VI A, and the SYK Lindbladians in Sec.VI B. Here, zmin denotes the complex eigenvalue with the smallest modulus.