Symmetry deduction from spectral fluctuations in complex quantum systems

The spectral fluctuations of complex quantum systems are known to be consistent with that from random matrices, but only for desymmetrized spectra. This impedes the analysis of experimentally measured or computed spectra. We show that for a spectrum composed of k>0 independent sequences, its k-th order spacing ratio distribution is identical to its nearest neighbor counterpart with modified Dyson index $\beta$ = k. The fluctuation character and symmetry structure of any arbitrary sequence of levels can be inferred from its higher order fluctuations without desymmetrization. This is shown for random matrices and also verified using measured levels of Ta nucleus, quantum billiards and spin chains.

The spectral fluctuations in complex quantum systems are analyzed using the theoretical framework of random matrix theory (RMT) in many areas of physics [1][2][3][4]. These include few-body systems studied in quantum chaos [5] to interacting many-body systems in condensed matter [6], nuclear [7] and atomic physics [8]. These fluctuations carry signatures of the distinct phases observed in physical systems, viz., integrable or chaotic limit of the underlying classical system [9], metallic or insulating phase [10], localized or thermal phase of many-body systems [11], low-lying shell model or mixing regime of nuclear spectra [12,13].
Beginning with Wigner's surmise [14] in the context of nuclear spectra, the general consensus now is that the spectral fluctuations of complex quantum systems, in suitable limit, display level repulsion consistent with that of an appropriately chosen ensemble of random matrices. For the special case of quantum chaotic systems, the Bohigas-Giannoni-Schmidt conjecture encapsulates this connection between the spectral fluctuations in RMT and chaotic phase of physical systems [15]. This has been amply verified in experiments [16], simulations [17] and derives some theoretical support based on semiclassical techniques [18].
Discrete symmetries of the system, i.e, invariance of the potential under parity, reflection, rotations, are crucial in realizing this connection between spectral fluctuations and dynamical phases. In the presence of symmetries, the Hilbert space of the system splits into invariant subspaces or the Hamiltonian matrix H becomes block diagonal, i.e., H = H 1 ⊕ H 2 ⊕ . . . H m , with each block H i , i = 1, 2...m characterized by good quantum numbers corresponding to the respective symmetries [5]. To compute any measure of spectral fluctuation, all the discrete levels must be drawn from the same subspace. If symmetries are ignored and the levels from different blocks are superposed, the genuine correlation between levels (that might have produced level repulsion) is masked by the near-degeneracies resulting in level clustering. This is misleading since this is also a spectral signature of integrable systems [19].
This implies that the level correlations are sensitive to the presence or absence of symmetries. It is then reasonable to expect the fluctuations of composite spectra, superposed from many independent blocks, to contain information about the symmetry structure of the entire system. However, any measure based on the nearest neighbor (NN) fluctuations, such as the popular NN level spacing distribution, will always tend to the Poissonion limit (level clustering) due to the superposition of noninteracting blocks [20]. In this work, rigorous numerical evidence is presented to show that the higher order level spacing ratio not only identifies the true fluctuation character, viz, level clustering or repulsion, but also allows us to deduce quantitative information about the symmetry structure of the composite Hamiltonian matrix H.
This result effectively obviates the need for symmetry decomposition of quantum systems and also allows any arbitrary sequence of experimentally observed levels, whose symmetry structure is unknown, to be analyzed. This is of considerable interest in RMT as well [21]. Let G be a random matrix such that G = G 1 ⊕ G 2 ⊕ . . . G m , a superposition of m blocks each of which is a Gaussian random matrix. If an arbitrary sequence of eigenvalues of G is given, the fluctuation properties and the block structure of G can be inferred from its higher order fluctuation statistics. The proposed method is straightforward, involving only the calculation of spacing ratios, in contrast to the complicated and approximate methods based on two-level cluster function that deduce m from any composite spectrum [1,26].
Consider a sequence of eigenvalues E i , i = 1, 2, . . . N of a quantum operator or a random matrix. Spectral fluctuations are relatively easier to analyze if spacing ratios defined as r i = Ei+2−Ei+1 Ei+1−Ei , i = 1, 2 · · · N − 2 are used instead of the spacings [22]. This is because the spacing ratios are independent of the local density of states and hence do not require spectral unfolding. For the case of the random matrix ensembles with Dyson index β = 1, 2 and 4, corresponding respectively to the Gaussian orthogonal, unitary and symplectic ensemble, the distribution of spacing ratios is given by [23] For β = 1 these RMT models are applicable to Hamiltonians with time-reversal invariance (TRI), which will be the main focus of this paper. For integrable systems, the ratio distribution becomes P P (r) = 1/(1 + r) 2 . As motivation, in Fig. 1, the numerically computed distribution of NN spacing ratios P (r) is shown for circular (integrable) [24] and stadium (chaotic) [25] billiards. The integrable billiards ( Fig. 1(a)) expectedly agrees with P P (r). Note that stadium billiard has C 2v point group symmetry with four irreducible representations (irreps). If the spectra from each irrep is analyzed separately, by BGS conjecture, an agreement with P (r, 1) of GOE is observed ( Fig. 1(c)). However, in Fig. 1(b), the spectra from all the irreps is superposed, and hence the ratio distribution is closer to P P (r) with pronounced deviation from P (r, 1). In such cases, as we demonstrate below, the true character of spectral fluctuations and the number m of independent spectra superposed can all be inferred using only the higher order spacing ratio distributions without apriori knowledge of its symmetry structure.
To this end, we consider the non-overlapping k-th order spacing ratio, defined as In what follows, the spectra from m independent blocks are superposed, and its distribution of k-th order spacing ratios is denoted by P k (r, β, m). We consider only the case β = 1. For the special case involving NN ratios, we denote P 1 (r, β, 1) = P (r, β). The motivation for considering higher order fluctuation statistics arises from a seminal result conjectured in Ref. [27] and proved by Gunson [28] for the case of circular ensembles of RMT. If two independent spectra from the circular orthogonal ensemble (COE) are superposed, and upon integrating out every alternate eigenvalue, the joint probability distribution of the remaining eigenvalues follow circular unitary ensemble (CUE) statistics. In terms of higher order measures, this result states that the second order statistics of two superposed COE spectra converges to NN statistics of CUE. This is reflected in the distribution of spacings and spacing ratios as well. In the limit of large matrix dimensions, this result holds for Gaussian ensembles too yielding P 2 (r, 1, 2) = P (r, 2) for two superposed spectra. This may be generalized for the superposition of m GOE spectra as implying that its k-th order spacing ratio distribution converges to NN statistics P (r, β ′ ) with β ′ = k. Equation 3 is the main result of the paper. In contrast to this, irrespective of how many uncorrelated spectra are superposed corresponding to integrable systems, the k-th order spacing ratio distribution can be obtained (details in supplementary information [29]) as For k = 1, this reduces to 1 (1+r) 2 , the correct limit for the NN spacing ratio for uncorrelated spectra. We note that Eq. 3 is reminiscent of a scaling relation for higher order spacing ratio distributions reported recently in Ref. [30].
For the superposition of m = 2 to 5 independent GOE spectra, validity of Eq. 3 is verified in Fig. 2 . In this figure, an excellent agreement is seen between the histograms obtained from the computed eigenvalues of GOE matrices and the solid line representing P (r, β ′ = k). For uncorrelated eigenvalues, a similar agreement with Eq. 4 is observed. In order to independently obtain a best quantitative estimate for β ′ in Eq. 3 for a given superposition of m spectra, we compute In this, I m obs (r, 1, m) and I(r, β ′ ) represents the cumulative distribution functions corresponding respectively to the observed histogram P k (r, 1, m) and the postulated function P (r, β ′ ). If the minima of D(β ′ ) occurs at, say, β ′ = β 0 , then β 0 is the best estimate consistent with the observed data. As seen in the insets of Fig. 2, the minima in D(β) coincides with the value of m, the number of superposed spectra.
A complete picture is revealed in Fig. 3 for a superposition of m = 4 independent GOE spectra, where the computed histogram for the k-th order ratio is shown for k = 2 to 7. Based on Eq. 3, we expect it to be consistent with P (r, β ′ = 4). For each k, P k (r, 1, 4) is matched against the corresponding P (r, β ′ ), and D(β ′ ) is calculated. Both visually and quantitatively (the minima of D(β ′ ) in Fig. 3(e)), best agreement is observed for k = 4, verifying the main result in Eq. 3. Significantly, for the superposed spectra, Eqs. 3-4 can be used to infer the correct nature of spectral fluctuations (level repulsion or clustering) and also to determine the number of superposed independent blocks for a random matrix or the number of diagonal blocks (irreps) in the Hamiltonian matrix of complex quantum system. Further, this result will be applied to chaotic systems possessing dif- for the billiards family computed by ignoring their symmetries. This corresponds to superposition of spectra from (a) k = 2, (b) k = 3 and (c) k = 4 irreps. The higher order distributions are best described by P (r, β ′ ) with β ′ = k as dictated by Eq. 3. The insets display D(β ′ ) and its minima corresponds to the correct number of irreps in the system. Also shown as inset is the shape of billiards with an arbitrarily chosen chaotic eigenstate to highlight its symmetry.
ferent symmetries, notably billiards and spin chains, and most importantly to the experimentally measured data of nuclear resonances. First we consider quantum billiards, in which a free particle is confined in a cavity defined by a variety of boundaries [31]. They are popular models in Hamiltonian chaos and mesoscopic physics [32] and variants have been experimentally realized as well [33]. In these systems, modifying the boundary or shape of the billiard changes its symmetry properties and also drives it from integrability to chaos. Their eigenspectrum is obtained by solving the Helmholtz equation with Dirichlet boundary conditions. For a billiard whose boundary is parameterized by r(φ) = r 0 (1 + ǫ cos φ), as ǫ varies from 0 to 1, the system transitions from integrable to chaotic dynamics. For ǫ = 0, a circular billiard shown in Fig.  1(a) is obtained. This is an integrable system and its higher order spacings are in agreement with Eq. 4 (See Ref. [29]). For ǫ = 1, the so-called cardioid billiard is obtained [34], possessing two irreps due to reflection symmetry about the horizontal axis. Thus, eigenlevels obtained disregarding symmetry would correspond to a superposition of two GOE spectra. As anticipated by Eq. 3, its second order spacing ratio distribution P 2 (r, 1, 2) is consistent with P (r, 2) ( Fig. 4(a)). A billiard with three irreps, similar in shape to one that has been experimentally realized [35], is obtained by parametrizing its boundary as r(φ) = r 0 (1 + 0.3 cos(3φ)). This model, with symmetries ignored and after removing degeneracies arising from the two-dimensional irreps, corresponds to a superposition of three chaotic spectra and the best fit for P 3 (r, 1, 3) is provided by P (r, 3) (Fig. 4(b)). A chaotic billiard with four irreps is the well-studied Bunimovich stadium billiard [36] (shown in Fig. 4(c)). This has reflection symmetry about both x and y axes and, in accordance with Eq. 3, displays the best correspondence for P k (r, 1, 4) with P (r, β ′ ) for k = β ′ = 4 (

4(c)
). For all of these cases, insets in Fig. 4 show that the minima of D(β ′ ) corresponds to β ′ = k, where k is the number of irreps. Thus, information about the fluctuation property and irreps can be obtained from higher order fluctuation statistics. Next, a spin-1/2 chain with the Hamiltonian [37] is considered, where L is the number of sites, J xy and J z are the NN coupling strengths in three directions (coupling along x and y being the same), and J ′ xy and J ′ z are the next NN coupling strengths. This system is integrable for η = 0 (as shown in Fig. 1 in Ref. [29]), and chaotic for η 0.2. The total spin in the z-direction, S z , is conserved and the Hamiltonian is block diagonal in the S z basis, with each block corresponding to a given value of S z . However, there still exist other symmetries, which would lead to the form of P (r) appearing to be integrable in this subspace (not shown here). For odd number of sites (L odd ), on computing the higher order spacing ratios and comparing with corresponding P (r, β ′ ), k = β ′ = 2 has the best fit ( Fig. 5(a)). For even number of sites (L even ), however, the best correspondence is for k = β ′ = 4 ( Fig. 5(b)). This is because for L odd or L even , the parity operator (with eigenvalues ±1) commutes with H, leading to two invariant subspaces in a given S z block. For L even , an additional rotational symmetry exists (with eigenvalues ±1) for the corresponding operator giving rise to four irreps. The other parameters used in Figs. 5(a,b) are J xy = J ′ xy = 1.0, J z = J ′ z = 0.5, with L even = 14 and L odd = 15.
Even for systems whose Hamiltonian is not welldefined or unknown as in the case of complex nuclei, experimentally observed nuclear resonance data can be analyzed to characterize its fluctuation statistics and find its number of irreps. We consider a sequence of experi- mentally observed neutron resonances for T a 181 nucleus [38] whose nearest neighbor spacing distribution is discussed in Ref. [13], and it does not match the Wigner surmise. On calculating higher order ratio distributions, remarkably, Eq. 3 holds good for k = 2, and this is further confirmed by the minima of D(β ′ ) for β ′ = 2 in Fig.  6. This indicates that two independent symmetry sectors might be present, and the resonances drawn from each symmetry sector displays level repulsion. This is indeed the case, as confirmed in Refs. [13,38], that this measured sequence consists of a superposition of levels having angular momentum J = 3 and 4. When symmetry decomposed, they are in broad agreement with Wigner surmise. Clearly, for an arbitrary sequence of measured levels, higher order spacing ratios based on Eq. 3 can unambiguously identify the true fluctuation character and the number of symmetry sectors. Typically in experiments and sometimes in simulations as well, the problem of missing levels is encountered [39] leading to incorrect identification of universality class of fluctuations and number of irreps. We tested the robustness of Eq. 3 to randomly missing levels. Randomly chosen entries in a given sequence of levels were deleted from a superposition of GOE spectra. Upon computing D(β ′ ) in each case (details in [29]), it was observed that Eq. 3 holds good even if up to 20 to 30% of the levels are removed as higher order spacings are largely insensitive to randomly missing levels. This is a practical advantage of analyzing higher order statistics.
To summarize, quantum systems must be symmetry decomposed in order to reveal its spectral fluctuation characteristics, level clustering or repulsion, without ambiguity. This implies that the fluctuations carry symmetry information though it was not possible to extract it from NN fluctuation statistics. In this work, it is demonstrated that the higher order spacing ratio distributions can reveal, apart from the fluctuation characteristics, quantitative information about symmetry structure. For a superposition of k independent spectra, the central result relates the k-th order spacing ratio distribution for matrices with Dyson index β = 1 to the corresponding NN statistics with β ′ = k. This powerful relation can be used to determine the number of irreps (or diagonal blocks) present in a Hamiltonian matrix and can be exploited to analyze any arbitrary sequence of experimentally measured or computed levels, even if the system's Hamiltonian and symmetry structure are unknown. Further, the higher order ratio distribution has been derived for uncorrelated eigenvalues, which may be used as a test of integrability. These results have been demonstrated using disparate physical systems like quantum billiards, spin chains and experimentally measured nuclear resonances.

Analytical expression
For a given sequence of uncorrelated eigenvalues, E 1 ≤ E 2 ≤ · · · E N , the spacings between nearest neighbours is defined as s i = E i+1 − E i , i = 1, 2, · · · N . The distribution of these spacings is of the form P (s) = e −s , and hence distributions of spacings and spacing ratios for integrable quantum systems are termed Poissonian.
The ratios of nearest neighbour spacings for these systems are defined as r i = s i+1 /s i , i = 1, 2, · · · N , and the distribution of these ratios is of the form P (r) = 1/(1 + r) 2 [1].
Ratios of higher order spacings may be defined as To obtain a form for the distribution of r (k) , the higher order spacings may be expressed in terms of nearest neighbour spacings as Then the distribution of s (k) i may be calculated as the distribution of a sum of k random variables s i , each of which is distributed as P (s) = e −s . For simplicity, s (k) i is denoted as z below. The distribution of z is given by Then the distribution of higher order spacing ratios is simply the distribution of the quotient of two random variables, each of which is distributed as Eq. 9. This distribution may be calculated as P (k) P (r) = |z|P (rz)P (z)dz (10) Substituting for P (z) and P (rz) from Eq. 9, P (k) This can be evaluated in terms of the incomplete gamma function Γ(x) as P (k) For k = 1, it reduces to the familiar form For k = 2, and for k = 4, P Comparison of analytical form of P

EFFECT OF MISSING LEVELS
The effect of missing levels in a given sequence of superposed spectra is studied by randomly deleting a fixed percentage of levels, and then calculating higher order spacing ratios, from a superpostion of GOE spectra of dimension N = 40000. In each case, D(β ′ ) is calculated, and the value of β ′ corresponding to the minima of D(β ′ ) corresponds to the best fit. Eq. 4 of the main paper, namely, P k (r, 1, m) = P (r, β ′ ), where β ′ = m = k, the expected value of β ′ is 2 (for k = 2) and 4 (for k = 4) respectively, given by the blue and red dashed lines in Fig 8. It may be observed that assuming even a 10% fluctuation in the numerical evaluation of β ′ , a significant deviation from the predicted β ′ occurs only when about 20% of the levels are missing. A similar behavior was seen for spin chains with 2 and 4 irreps as well (not shown). * E-mail: harshini.t@gmail.com † E-mail: santh@iiserpune.ac.in