## Abstract

We discuss the decay rates of chaotic quantum systems coupled to noise. We model both the Hamiltonian and the system-noise coupling by random N×N Hermitian matrices, and study the spectral properties of the resulting Liouvillian superoperator. We consider various random-matrix ensembles, and find that for all of them the asymptotic decay rate remains nonzero in the thermodynamic limit; i.e., the spectrum of the superoperator is gapped as N→∞. For finite N, the probability of finding a very small gap vanishes as P(Δ)∼ΔcN, where c is insensitive to the dissipation strength. A sharp spectral transition takes place as the dissipation strength is increased: for dissipation beyond a critical strength, the slowest-decaying eigenvalues of the Liouvillian correspond to isolated "midgap" states. We give evidence that midgap states exist also for nonrandom system-noise coupling and discuss some experimental implications of the above results.

Original language | American English |
---|---|

Article number | 234103 |

Journal | Physical Review Letters |

Volume | 123 |

Issue number | 23 |

DOIs | |

State | Published - 5 Dec 2019 |

### Bibliographical note

Funding Information:https://orcid.org/0000-0002-0999-2355 Can Tankut 1 Oganesyan Vadim 1,2 Orgad Dror 3 Gopalakrishnan Sarang 1,2 Initiative for the Theoretical Sciences, The Graduate Center, 1 CUNY , New York, New York 10012, USA 2 Department of Physics and Astronomy , College of Staten Island, Staten Island, New York 10314, USA Racah Institute of Physics, 3 The Hebrew University , Jerusalem 91904, Israel 5 December 2019 6 December 2019 123 23 234103 1 April 2019 © 2019 American Physical Society 2019 American Physical Society We discuss the decay rates of chaotic quantum systems coupled to noise. We model both the Hamiltonian and the system-noise coupling by random N × N Hermitian matrices, and study the spectral properties of the resulting Liouvillian superoperator. We consider various random-matrix ensembles, and find that for all of them the asymptotic decay rate remains nonzero in the thermodynamic limit; i.e., the spectrum of the superoperator is gapped as N → ∞ . For finite N , the probability of finding a very small gap vanishes as P ( Δ ) ∼ Δ c N , where c is insensitive to the dissipation strength. A sharp spectral transition takes place as the dissipation strength is increased: for dissipation beyond a critical strength, the slowest-decaying eigenvalues of the Liouvillian correspond to isolated “midgap” states. We give evidence that midgap states exist also for nonrandom system-noise coupling and discuss some experimental implications of the above results. National Science Foundation 10.13039/100000001 DMR-1653271 DMR-1508538 United States-Israel Binational Science Foundation 10.13039/501100001742 2018159 How nonequilibrium systems relax to steady states is a central topic in many-body dynamics. For a system coupled to an external environment, the rate of approach to the steady state depends nontrivially on the system-environment coupling, i.e., the dissipation strength. Weak dissipation enhances decay, but strong dissipation can suppress it through the quantum Zeno effect [1–4] . These phenomena have been extensively studied, both theoretically and experimentally, for specific models [5–15] . Here, in contrast, we discuss them within the generic setting of random matrix theory (RMT) [16] . Historically, RMT was introduced to describe complex dynamical systems for which a microscopic description would be intractable [17,18] . RMT is believed to describe the generic long-time behavior of chaotic quantum systems [19–21] . It predicts universal features in the density of states and level statistics that have been verified numerically and experimentally in many settings [22] . RMT has also been extended to open quantum systems, mainly via scattering theory [23,24] , and has been used to model evolution under effective non-Hermitian Hamiltonians [25] . Given this background, it is notable that, until very recently [26,27] , RMT was not applied to the Lindblad master equation [28–30] , which is the standard framework for describing open quantum systems in fields ranging from quantum optics to mesoscopics. Much is known about the steady states of specific Lindblad equations, and some notions of universality have been developed for these states and the phase transitions between them [13,31] . However, the issue of universality in the dynamics of master equations, e.g., their approach to a steady state, remains largely unexplored (see however Refs. [32–37] ). Here, we use RMT to identify universal dynamical properties of a class of open quantum systems: those coupled to classical white noise, or equivalently to purely dephasing Markovian baths [38–40] . We study such dynamical properties by exploring the spectrum of the Liouvillian “superoperator” L (defined below), which generates the dissipative evolution of the system’s density matrix ρ . The spectrum of L has a zero mode (the steady state), which in the systems we consider is unique and proportional to the identity (i.e., an infinite-temperature state). The other eigenmodes of L have strictly negative real parts, which correspond to decay rates. Thus, the eigenvalue with smallest nonzero real part sets the timescale on which the system asymptotically approaches its steady state. We find that in the N → ∞ limit this approach is always “fast” in the sense that very small decay rates become infinitely improbable. However, the distribution of the slow decay rates evolves nonanalytically with the dissipation strength. While for weak dissipation the asymptotic decay rate resides at the edge of a continuum of low-lying L eigenstates, it, and possibly the next few slowest rates, split from the continuum and become isolated “midgap” states [41] beyond a critical dissipation strength, see Fig. 1 . This transition seems not to have been noticed in previous work, although related phenomena have been observed in the mathematical literature on classical dissipative systems [42,43] . In the following we establish these results, comment on the finite- N scaling of the gap probability distribution, demonstrate that midgap states appear regardless of details of the system-noise coupling and discuss some of their experimental implications. 1 10.1103/PhysRevLett.123.234103.f1 FIG. 1. (a) Distribution of Liouvillian eigenvalues in the complex plane near the origin, obtained from 60 realizations of an N = 100 system with large dissipation, γ = 50 . (b) Probability density of the real part of the eigenvalues. A well-defined crossing point in the finite-size scaling analysis identifies the edge of the continuous spectrum, which is largely due to complex eigenvalues. The purely real eigenvalues dominate the small | λ | regime. Clearly visible is an isolated state at λ ≃ - 0.032 , and possibly another at λ ≃ - 0.057 . (c) Spectral gap and the edge of the continuum as a function of dissipation γ . The gap values obtained from the largest accessible systems N = 160 and by extrapolation from smaller sizes agree and scale as ( γ / 2 ) ± 1 for weak/strong dissipation. For γ ≥ 6 the edge of the continuum deviates from the gap, and the eigenvalues closest to zero become isolated. This transition is also evident from the inset, which shows that the gap between the first two “excited” eigenstates of the Liouvillian shrinks with N for γ < 6 and increases for γ > 6 . Master equation.— We consider the Lindblad equation ∂ t ρ = L ρ ≡ - i [ H , ρ ] + γ ∑ k = 1 n d ( L k ρ L k † - 1 2 { L k † L k , ρ } ) , (1) where H is its Hamiltonian, L k are “jump operators” representing the coupling of the system to its environment, and γ is the dissipation strength. We focus on the case where the Hamiltonian is an N × N random matrix from the Gaussian Orthogonal Ensemble (GOE), whose elements have zero mean and variance ⟨ H i j H k l ⟩ = ( 1 / 2 N ) ( δ i k δ j l + δ i l δ j k ) , and where there is a single jump operator ( n d = 1 ), statistically independent of H and also drawn from the GOE. Other cases will be addressed briefly at the end (and in [44,45] ). Note that we have scaled the distribution such that the N → ∞ spectra of H and L reside in [ - 2 , 2 ] . We concentrate on the large- N limit, but also discuss some nonperturbative features that arise at finite N . General properties.— The right-hand side of Eq. (1) is a linear “superoperator” acting on ρ , allowing us to set up the eigenvalue problem L ρ = λ ρ . The Liouvillian L is a completely positive trace preserving map, and for γ ≥ 0 all eigenvalues λ satisfy Re λ ≤ 0 , with at least one eigenvalue at λ = 0 . For the model considered here, the steady state is unique and given by the infinite temperature thermal state proportional to the identity. Furthermore, the evolution under the generator L does not take states out of the space of density matrices, which requires it to preserve both trace and Hermiticity. It then follows that ( L ρ ) † = L ( ρ † ) and as an N 2 × N 2 matrix its components obey the relation L j i l k = L i j k l * . This “conjugation” symmetry can be used to show that the eigenvalues of L come in complex conjugate pairs [45] , which can collide at an “exceptional point” and become purely real [50] , or vice versa. A nontrivial property of the model is that at least N eigenvalues are necessarily purely real. This property follows from combining the conjugation symmetry noted above with the additional symmetry of the Liouvillian L i j k l = L k l i j which follows when H and L are real symmetric matrices [45] . Weak dissipation.— In the absence of dissipation the dynamics is governed by the operator - i [ H , · ] . The Liouvillian eigenvalues correspond to differences, - i ( E α - E β ) , between Hamiltonian eigenenergies. This spectrum has N eigenstates of the form | α ⟩ ⟨ α | with eigenvalue zero, and N ( N - 1 ) states of the form | α ⟩ ⟨ β | with imaginary nonzero eigenvalues that come in complex conjugate pairs. In the language of nuclear magnetic resonance, the former correspond to “populations” and the latter to “coherences” [51] . Level repulsion in the spectrum of H manifests itself as a suppression of the density of eigenvalues with small nonzero imaginary parts | Im ( λ ) | ≲ 1 / N . When γ is small, we can treat it perturbatively. To first order the populations and coherences decouple. In the N -dimensional subspace of populations one obtains a symmetric classical master equation, in which the off-diagonal terms are γ | L α β | 2 with L α β = ⟨ α | L | β ⟩ , while the diagonal terms are determined by the constraint that each column should sum to zero. Every term in this matrix is sign definite, so one can disorder average the master equation to get a matrix element (i.e., a transition rate) ⟨ | L α β | 2 ⟩ ∼ 1 / ( 2 N ) between any pair of populations. The resulting spectrum in the population subspace contains a unique steady state and N - 1 degenerate states with eigenvalue λ = - γ / 2 . Meanwhile, at first order and up to 1 / N corrections, each coherence | α ⟩ ⟨ β | picks up a (real) perturbative shift of the form - ( γ / 2 ) ∑ η ( | L α η | 2 + | L β η | 2 ) = - γ / 2 . Therefore, N 2 - 1 states rigidly move away from the steady state by an amount - γ / 2 . Since there is a hard gap at first order in γ , we expect it to be robust, provided that higher-order perturbative corrections do not diverge in the large- N limit. We have checked this to order γ 2 [45] . Strong dissipation.— As γ → ∞ the spectrum of L becomes λ a b = - ( γ / 2 ) ( κ a - κ b ) 2 , where κ a are the L eigenvalues. Here too, there are N “diagonal” zero modes of the form | a ⟩ ⟨ a | , which we call “ L -populations.” We dub the remaining eigenstates, of the form | a ⟩ ⟨ b | , “ L -coherences.” In this limit the spectrum is real and gapless, with bandwidth 4 γ and level spacing ∼ γ / N 2 . Next, we perturb in the Hamiltonian, taking N → ∞ at large but finite γ . At first order, an L -coherence | a ⟩ ⟨ b | picks up a purely imaginary O ( 1 / N ) shift τ a b = i ( ⟨ b | H | b ⟩ - ⟨ a | H | a ⟩ ) , with the result λ a b = - γ 2 ( κ a - κ b ) 2 + i τ a b . (2) The L -populations, however, remain exact zero modes to first order. To resolve these degeneracies, one must go to second order in the Hamiltonian where populations and coherences mix. One can disentangle the subspaces via a Schrieffer-Wolff transformation [47] . Each L -population connects to 2 ( N - 1 ) coherences, obtained by changing one (but not both) of its indices. Consequently, any two populations, | a ⟩ ⟨ a | and | b ⟩ ⟨ b | are connected by the second order matrix element φ a b = - | H a b | 2 λ a b - | H a b | 2 λ a b * = 4 γ | H a b | 2 ( κ a - κ b ) 2 ( κ a - κ b ) 4 + ( 2 / γ ) 2 τ a b 2 . (3) One may worry that since many [ O ( N 3 / 4 ) ] of the L -coherences which are coupled to a population have eigenvalues that are smaller than a typical H a b , the coherences are strongly coupled to the populations and there is no way to separate them. However, note that the matrix elements (3) vanish for pairs of L -populations with very similar eigenvalues κ a ≈ κ b . Consequently, coherences with small | λ a b | do not couple strongly to the populations and the required separation between the two sets can be carried out to second order in H [45] . In the subspace of L -populations the spectrum consists of two types of states. The lowest few eigenstates are isolated and remain separate from the continuum even after disorder averaging. Their form and their eigenvalues λ n = - ( 2 / γ ) n , with integer 0 ≤ n ≤ n max , can be obtained analytically [45] , see Fig. 2 for comparison with numerical results. For more negative λ , the peaks corresponding to individual eigenvalues become more closely spaced, broader, and merge into a continuum upon disorder averaging [45] . Concomitantly, in the L -coherences subspace, second-order perturbation theory shifts states near λ = 0 by an amount - c / γ [45] , with c a constant, resulting in an overall gapped spectrum. Numerically, we find c > 2 such that the λ 1 state (with possibly additional n > 1 levels) always appears to the right of the continuum, leading to the large- γ behavior of the gap Δ = 2 / γ , in accord with the quantum Zeno effect. 2 10.1103/PhysRevLett.123.234103.f2 FIG. 2. The first (a) and second (b) “excited” eigenstates of the Liouvillian, obtained from averaging over 60 realizations of an N = 100 system with γ = 50 . Shown are components along the L -populations | κ ⟩ ⟨ κ | , parametrized by their L eigenvalue, κ . The components along the L -coherences are down by 2–3 orders of magnitude and are negligible in (a). Their spectral weight in (b) is accounted for by rescaling the eigenstate. The expected N , γ → ∞ behavior in terms of the Chebyshev polynomials U n ( κ / 2 ) is shown in red. Numerical results.— Figure 1(c) presents numerical results on the distribution of gaps at various values of γ and N . Since the results obtained by finite-size extrapolation agree well with those obtained for the largest manageable systems ( N = 160 ) using the Arnoldi method, we expect that the calculated gap values are close to the thermodynamic limit. Both methods give a transition from Δ = γ / 2 at weak dissipation to Δ = 2 / γ at large dissipation with a maximum around γ ≈ 3.5 . We used finite-size scaling to locate the edge of the continuous spectrum of L . Within our numerical constraints we were able to detect, in the range 0.1 ≤ γ ≤ 100 , sharpening of the density of states (DOS) edge with system size, giving a well-defined crossing point [Fig. 1(b) ]. This is consistent with the conjecture that the DOS jumps discontinuously at the edge of the spectrum. When γ ≲ 6 the jump in the DOS occurs at the gap, but for larger γ the jump is at a larger value of | λ | than the gap. In this regime, it seems that there is at least one isolated state, and possibly multiple states, at decay rates appreciably lower than the edge of the continuous spectrum. These can clearly be seen as isolated peaks in the DOS for the larger system sizes in Fig. 1(b) . To better locate this transition, we computed the γ dependence of the distance between the two smallest (in absolute value) nonzero eigenvalues of L . A clear flow reversal is seen in the inset of Fig. 1(c) , indicating a phase transition at γ c ≈ 6 . For γ ≲ 6 the two eigenvalues approach each other with increasing system size, as expected for continuum states. For larger γ , however, they move away from each other with increasing N , the DOS jump steepens and the isolated eigenvalue becomes more clearly separated from the continuum. Finite-size corrections.— Our results indicate that as N → ∞ the disorder-averaged density of states is zero at sufficiently small | λ | . We now discuss how this hard gap sets in, by considering finite-size systems. We argue on general grounds that the density of states as λ → 0 is sensitive to the symmetry classes of H and L k , as well as to the size of the Hilbert space N and the number of distinct jump operators L k ; however, its λ → 0 shape is γ independent. For the case of n d jump operators, the probability of finding a gap Δ scales as P ( Δ ) ∼ Δ β n d ( N - 1 ) / 2 - 1 , (4) where β = 1 , 2, 4 for orthogonal, unitary, and symplectic ensembles. This prediction is compared with numerics for the orthogonal ensemble (Fig. 3 ), where results for other ensembles are shown in [45] . The figures plot the asymptotic density of eigenvalues, and show that it compares favorably with the scaling behavior, Eq. (4) . When the Hamiltonian and jump operators belong to different ensembles, the appropriate value of β is that corresponding to the lower-symmetry (higher β ) ensemble. 3 10.1103/PhysRevLett.123.234103.f3 FIG. 3. Eigenvalue density close to λ = 0 at γ = 2 for small system sizes, N = 3 , 5, 7, averaged over 5 × 10 5 realizations; the behavior is consistent with Eq. (4) (straight lines). Equation (4) can be understood via a golden-rule argument. For simplicity we first treat the case of a real ( β = 1 ) H and L with n d = 1 . Consider a system initialized in the eigenstate | n ⟩ of the Hamiltonian. In the presence of a Markovian bath (which has an energy-independent density of states) its decay rate is given by Γ n = ∑ m ≠ n | ⟨ n | L | m ⟩ | 2 , where m runs over all other eigenstates of H . There are N - 1 matrix elements in this expression, and (given that H and L are mutually uncorrelated and taken from the GOE) L n m ≡ ⟨ n | L | m ⟩ can be regarded as independent real Gaussian random numbers. For Γ n to be small we need all matrix elements to be small. One can approximate the cumulative probability distribution P ( Γ n ≤ Δ ) ∼ ∏ m P ( L n m 2 ≤ Δ ) ∝ ( Δ ) N - 1 , (5) from which Eq. (4) follows by differentiation. The other cases are simple generalizations. For example, in the unitary ensemble, each matrix element is a complex number, and is only small when its real and imaginary parts are separately small, yielding a factor of two in the exponent. In general, therefore, the disorder-averaged DOS at a given small λ vanishes exponentially in N . Discussion.— In this Letter we have explored the spectral structure of low-lying (i.e., long-lived) states of a dephasing Lindblad master equation with random Hamiltonian and jump operator. While the conditions under which a random-matrix Hamiltonian captures the essential features of a physical system are well understood and generic, the conditions for the system-bath coupling to permit an RMT description are less clear. Hence, we have checked that our main results, i.e., the gap in the L spectrum and the appearance of midgap states for large γ , are robust to details of the coupling. Indeed, we find that L -level repulsion is inessential since a model with uniformly distributed L eigenvalues yields similar results [45] . Moreover, the coupling need not be random at all, as shown by Fig. 4(a) where a gap and midgap states appear for a jump operator with evenly spaced eigenvalues. Rather, the key ingredients are (1) that both H and L have bounded spectra, and (2) that eigenvectors of L look random in the eigenbasis of H and vice versa [45] . 4 10.1103/PhysRevLett.123.234103.f4 FIG. 4. (a) Probability density of the real part of the Liouvillian eigenvalues near the origin for a model whose H is drawn from the GOE and where the eigenvalues of the jump operator, L , are regularly spaced in [ - 2 , 2 ] . The results are for γ = 50 and averaged over 100 realizations. (b) The first two “excited” eigenstates of the Liouvillian with N = 60 . Shown are components along the L -populations | κ ⟩ ⟨ κ | , parametrized by their L eigenvalue, κ . (c) Autocorrelations in a model with GOE H and L , γ = 200 and N = 60 , averaged over 100 realizations. The solid line depicts 1 / t decay. The most direct experimental consequence of our results is the nonanalytic behavior of the approach to steady state as dissipation increases. Furthermore, our analysis implies two outcomes for systems in the strong-dissipation regime. First, the decay of an initial L -population, | κ ⟩ ⟨ κ | , has a specific functional form on long timescales. This can be probed by measuring the observable which determines the system-bath coupling (e.g., the dipole moment of a chaotic quantum dot coupled to electromagnetic noise). The system then collapses to an L -population, which for large γ subsequently evolves via hopping in the subspace of L -populations. Due to the relatively local structure of the hopping (3) , probability spreads ballistically, eventually becoming uniform. This is reflected by the autocorrelations between consecutive measurements outcomes C ( t ) = Tr [ | κ ⟩ ⟨ κ | e L t | κ ⟩ ⟨ κ | ] , which decays as 1 / t for γ / 2 n max < t < γ / 2 , as follows from expanding | κ ⟩ ⟨ κ | in the L eigenmodes [45] , see Fig. 4(c) . Here, n max is the number of states separated from the continuum, which we find grows with N [45] . Second, the late-time distribution of the measurement outcome is determined by the slowest-decaying modes, and its deviation from the distribution in the steady state is given by the distribution of L -populations in the midgap states. We have obtained explicit expressions for the latter within the model studied, illustrating that they correspond to the longest-wavelength modulations of the L -populations, see Fig. 2 . We have also checked that their general form is insensitive to the choice of model, as shown by Fig. 4(b) and in [45] . This work has focused on the GOE because of its simplicity and its direct relation to noisy dynamics. A companion paper [44] finds similar features for ensembles of non-Hermitian L with no symmetries. A natural question is how far these results extend to other random matrix ensembles, as well as to systems with locality constraints [52] , such as band matrices. S. G. acknowledges support from NSF Grant No. DMR-1653271. S. G., V. O., and D. O. acknowledge support by the United States-Israel Binational Science Foundation (Grant No. 2018159). V. O. acknowledges support from NSF Grant No. DMR-1508538. T. C. acknowledges helpful discussions with A. Abanov, J. Feinberg, T. Seligman, W. Tarnowski, and V. K. Varma. D. O. acknowledges useful discussions with O. Agam.

Funding Information:

S.G. acknowledges support from NSF Grant No.DMR-1653271. S.G., V.O., and D.O. acknowledge support by the United States-Israel Binational Science Foundation (Grant No.2018159). V.O. acknowledges support from NSF Grant No.DMR-1508538. T.C. acknowledges helpful discussions with A. Abanov, J. Feinberg, T. Seligman, W. Tarnowski, and V.K. Varma. D.O. acknowledges useful discussions with O. Agam.

Publisher Copyright:

© 2019 American Physical Society.