Correlators for the Wigner-Smith time-delay matrix of chaotic cavities

We study the Wigner-Smith time-delay matrix $Q$ of a ballistic quantum dot supporting $N$ scattering channels. We compute the $v$-point correlators of the power traces $\mathrm{Tr} Q^{\kappa}$ for arbitrary $v\geq1$ at leading order for large $N$ using techniques from the random matrix theory approach to quantum chromodynamics. We conjecture that the cumulants of the $\mathrm{Tr} Q^{\kappa}$'s are integer-valued at leading order in $N$ and include a MATHEMATICA code that computes their generating functions recursively.


Introduction
The proper delay times τ 1 , . . . , τ N are the eigenvalues of the Wigner-Smith time-delay matrix Q, which is defined in term of the scattering matrix S(E) by [19,40,46] (1) The Wigner-Smith matrix provides an overall description of the dynamical aspects of a scattering process in terms of the phase shift between a scattered wavefunction and a freely propagating one. Ballistic quantum dots provide a physical realization of systems whose scattering matrix S(E) is effectively modelled by a random unitary matrix [5], where N is the total number of scattering channels. The joint probability law of the rescaled inverse delay times λ i = (N τ i ) −1 was derived in [11] when the scattering matrix S(E) belongs to the Dyson's circular ensembles [18] and was generalised in [31] for systems whose symmetries belong to the classification introduced by Altland and Zirnbauer [2,48]. It is known in Random Matrix Theory (RMT) as the eigenvalue distribution of a Wishart-Laguerre matrix ensemble (Eq. (31)-(32) below). This ensemble belongs to the class of one-cut β-ensembles, where β is the Dyson index. Without loss of generality, we will consider systems such that β = 1, 2 or 4 according to the internal symmetries of the problem (see the reviews [7,25,42] for an extensive discussion).
In this paper we study the power traces of the time-delay matrix Q, (2) T κ = N κ−1 TrQ κ , κ ≥ 0, and we derive new formulae for their joint cumulants in the limit N → ∞.
Traces of the Wigner-Smith matrix for chaotic cavities have been investigated intensively in the past. Many results are available on the Wigner time-delay T 1 = TrQ, both for finite and large N [22,28,32,32,41,43]. A few results for higher-powers T κ = TrQ κ are also available. The covariances of T 1 and T 2 have been computed at leading order in N in [24]. Soon after one of the authors [14] considered the joint distribution of all the power traces T κ (κ ≥ 0) and computed the averages and the full covariance structure at leading order in N . Recently, exact finite-N formulae for averages of arbitrary products of power traces T κ (κ ≥ 1) have been derived for β = 2 [34]. The regime of non-ideal couplings with the external leads, of obvious importance in physical applications, has also been investigated [37], but is not considered in this paper.
In deriving our results we follow ideas similar to those discussed in [15]. Here, however, we prefer to use the convenient language and formalism of the so-called 'loop calculus' [33] developed in the RMT approach to quantum chromodynamics. The functional derivative method discussed in [15] can be concretely realized in this framework by means of the 'loop insertion operator'. In addition, the derivation of our results is simplified by using a suitable 'set of coordinates' first introduced by Ambjørn et al. in [4]. We are not aware of any general result on the leading order cumulants of random Wigner-Smith matrices (3) so far; this is somewhat astonishing since the method of loop equations -well established in a different area of physics -can be applied directly.
The aim of this paper is thus threefold. First, we present a number of new results on the Wigner-Smith matrix of ballistic quantum dots. We believe that the ideas discussed in this work may be useful in related problems; using techniques similar to those used here it is possible, for instance, to obtain results for the moments of transmission matrices (see also recent applications of loop equation formalism in Gaussian and circular β-ensembles [47]). Hence, our second aim is to popularize the loop calculus in the random matrix approach of quantum transport in chaotic cavities. Third, we provide strong evidences of an underlying combinatorial problem whose interpretation might be found through careful considerations of the semiclassical approaches to quantum chaos.
The plan of this paper is as follows. We start in Section 2 with some preliminaries on the one-cut β-ensembles. In Section 3 the iterative scheme for the v-point Green functions is presented. Section 4 contains the main result (Theorem 1) and its proof. In Section 5 we discuss further results (Theorem 2) and a conjecture. The paper is complemented by three appendices where we collect a few properties of cumulants (Appendix A), we present the derivation of some identities involved in the iterative scheme for the Green's functions (Appendix B), and we include a MATHEMATICA code to implement in a systematic way the generating functions (Appendix C).

Generalities about the one-cut ensembles
Consider an ensemble of N × N random matrices where β = 1, 2 or 4 if φ belongs to the set of real symmetric, complex hermitian or quaternion self-dual matrices, respectively. This ensemble is invariant under the adjoint action of the group O(N ), U (N ) or Sp(N ) for β = 1, 2 or 4, respectively, and hence it is usually referred to as an invariant ensemble. We assume that the real potential V is independent of N and β and such that its derivative V ′ (the force field) is a rational function. For our purposes it is convenient to consider the parametric family of potentials where it is understood that only finitely many t i 's are nonzero. The nonzero couplings are assumed to be such that the partition function Z =´dφe −βN TrV (φ) is finite. Of course, there are many admissible choices for the couplings (hence for V ). One of the simplest is the quadratic potential V (x) = x 2 /4 (that is t 2 = 1/2 and t i = 0 for i = 2), corresponding to the celebrated Gaussian ensembles GOE, GUE and GSE. Here we anticipate that in our problem we will set the couplings t 0 = −t 1 = −1/2 and t i = 0 for i ≥ 2. Note that necessarily t 0 ≤ 0 in (5) and, whenever t 0 = 0 it is implicitly assumed that the probability measure (4) is supported on the set of positive definite matrices, φ ≥ 0. Averages of observables are defined in the usual way Observables that depend uniquely on the spectrum {λ i } (i = 1, . . . , N ) of φ are particularly meaningful. In this case, (6) can be reduced to an integral over the joint probability distribution of the eigenvalues of φ, which is given by For all v ≥ 1, let us define the v-point connected Green function of φ Note that G v (z 1 , . . . , z v ) is a nonrandom complex function of v variables and is well-defined for ℑz i > 0. Of course, G v (z 1 , . . . , z v ) is symmetric in its variables. We denote the leading order term in N of (8) by (9) G v,0 (z 1 , . . . , z v ) = lim N →∞ The indices 'v, 0' denote that G v,0 is the generating function of the leading order in N (the g = 0 term in the 'genus expansion' 1/N g ) of the v-th order cumulants. The existence of a 1/N expansion has been proved in [20] for polynomial potentials and β = 2, and later for a general class of one-cut matrix models and every β > 0 in [9]. Note that, in general, G v,0 (z 1 , . . . , z v ) is not continuous. For large N the eigenvalues of φ concentrate on a subset of R that corresponds to the singularities of G v,0 . For instance, the discontinuity of the 1-point Green function G 1,0 (z) = lim N →∞ on the real line is proportional to the limit density of the eigenvalues of φ.
We now assume that the potential V is such that the singularities of G 1,0 (z) consist of a finite single cut [A, B], i.e. the limit density of the eigenvalues is concentrated on a single interval. Under this one-cut assumption we have [3] (10) Hereafter the integration contour C turns anticlockwise around the cut [A, B] without enclosing the point ω = z. The real edges of the cut A and B are determined by the normalization condition G 1,0 (z) ∼ 1/z as z → ∞:

Iterative scheme for the correlators
The following ideas have been developed in the RMT approach to quantum chromodynamics and 2D quantum gravity (see [1,3,4,12,13,17,26,33,45]). We include them to make the paper self-contained. In the following, we compute the v-point correlators using the functional derivative method outlined below. The precise asymptotic meaning of these formal computations can be found in the mathematical papers [9,20].
3.1. Functional derivative method. The whole family of Green's functions (at leading order in N ) can be obtained from G 1,0 (z) using an iterative procedure involving functional derivatives. Indeed, we can write the variational formula, where the functional derivative operator is defined as Hence, for all v > 1:

3.2.
Total derivative formula. The functional derivative method can be implemented by using the 'loop-insertion operator' and a suitable set of variables. A natural choice would be to work with the couplings t i 's that parametrize the potential V (x) in (5). However, it was realized in [4] that explicit calculations are easier by working in terms of the edges A and B and the coordinates M ℓ and J ℓ defined as Note that, from (11), The loop-insertion operator [3,4,20] is defined as The action of this operator has an interesting interpretation in the theory of symmetric forms on Riemann surfaces (see, for instance, [13]). Using The loop-insertion operator and the new set of variables provide a handy realization of the functional derivative in terms of a total derivative formula: Therefore, in order to apply the functional derivative δ/δV (z), we need to compute the variation of A, B, M ℓ and J ℓ with respect to the external potential V (z): and For completeness we present the derivation of these formulae in Appendix B.

3.3.
Multi-point Green functions. The 1-point Green function in the one-cut regime is given by (10) supplemented by (11). A first functional differentiation provides the 2-point correlator (24) Note that this result depends on the potential V (x) only through the edges A and B of the cut (determined by (11)). This is called macroscopic universality of the 2-point correlators [6,10]. G 2,0 (z 1 , z 2 ) has no explicit dependence on V (x). As a consequence of the functional derivative identity (12) the v-point Green functions (v ≥ 2) have no explicit dependence on the potential. However, the universality of the correlators decreases as v increases. For instance, the 3-point correlator is given by Note that δA δV (z) and δB δV (z) contain M 1 and J 1 (see (20)- (21)). Hence G 3,0 depends on V through A, B, M 1 and J 1 (it is 'less universal' than G 2,0 ). A further application of the functional derivative gives the 4-point correlator It is easy to see that in general, the v-point Green function (v ≥ 2) depends on the edges A and B and only on the first v − 2 elements of the sequences {M ℓ } and {J ℓ }.
For general v > 2 we can write The recurrence relation (27) supplemented by the initial datum G 2,0 (z 1 , z 2 ) in (45) provides the full family of v-point connected correlators at leading order in N for the class of one-cut matrix ensembles.
Remark 1. We emphasise that this route is independent of the specific model (the potential V (x)) and is valid within the class of one-cut ensembles. The leading order of the v-point Green function has the same structure determined by the one-cut property. Equations (10)-(11)-(24)- (27) give the generic solution of the multi-point correlators G v,0 as a function depending parametrically on A, B, M ℓ and J ℓ . This general structure is applied to a specific model by inserting the model-dependent parameters A, B, M ℓ and J ℓ determined by (11)-(15)- (16). (Note that the values of the edges A and B as well as {M ℓ } and {J ℓ } are functions of the couplings t i of the potential.) For v ≥ 3 the multi-point Green functions can be written as is a symmetric polynomial of degree at most 2v − 5 in each variable. Note also that at leading order in N the Dyson index β plays no special role other than being a multiplicative parameter of the multi-point Green functions G v,0 . (See [21] for additional remarks.) (3):

Main Result
The indices 'v, 0' denote that F v,0 is the generating function of the of the v-th order cumulants at leading order in N (the 'genus' g = 0 term) . Here we provide explicit expressions for F v,0 (z 1 , . . . , z v ). The joint distribution of the inverse delay times This distribution is of the form (7) with external potential which belongs to the general class (5) with the choice t 0 = −1/2, t 1 = 1/2 and t i = 0 for i ≥ 2.
Since V (x) in (32) is strictly convex, for large N the density of inverse delay times concentrates on a single interval (the one-cut assumption is then satisfied). Hence, the joint law of the inverse delay times λ i 's fits in the iterative scheme described in Section 3. Before stating the result we need the following Lemma that specifies the actual values of the edges A, B and the variables M ℓ , J ℓ for the Wigner-Smith inverse matrix. In the following we will use the lower-case symbols a, b, m 1 , j 1 , etc. to denote these model-dependent values. Similarly, we denote the v-point Green function of the inverse delay times and its leading order in N by In (33) the connected average refers to the joint distribution of the inverse delay times (31). Lemma 1. The leading order 1-point Green function of the inverse delay times is Let us denote with V given in (32) and a, b in (36). Then Remark 2. Equations (35)-(36) mean that, as N → ∞, the density of inverse delay times is the Marčenko-Pastur distribution [30] supported on the interval [a, b]: We can now state our central result.

A few values of lim
Proof of Theorem 1. To prove the Theorem we first establish the identity (42) relating the generating functions F v,0 (z 1 , . . . , z v ) and the leading order multipoint Green functions g v,0 (z 1 , . . . , z v ). For a family of linear statistics X (j) on the inverse delay times λ i (here f j (x) denotes a generic function) the identity follows by Cauchy's integral formula with g v defined in (33). In (52), the integration contours enclose anticlockwise the inverse delay times λ i and no singularities of Let us specialize (53) to the Wigner-Smith time-delay moments The integral can be evaluated as the sum of residues. The possible singular sets are z i = 0 and z i = ∞. For v ≥ 2, the integrand has no residue at z i = ∞ and hence the only contribution comes from the residue at z i = 0: where ǫ is sufficiently small and the factor (−1) v takes into account the fact that z i = 0 is 'outside' the contour C. The right-hand side of (55) is a generic coefficients of the power series expansion of z 1 · · · z v g v,0 (z 1 , . . . , z v ) at z 1 = · · · = z v = 0, and this proves the identity (42) for v ≥ 2. A similar argument applies to v = 1, where the additional δ v,1 term in F 1,0 is due to the fact that g 1,0 is not analytic at z = 0 and z = ∞. As already discussed, the Green functions of the inverse delay times λ i 's are amenable to the iterative scheme presented in Section 3. Hence, assuming Lemma 1, the proof is complete.
The rest of this section is devoted to proving Lemma 1.
Proof of Lemma 1. The lemma can be easily proved by using residue calculus. We first compute the edges a and b. For the potential (32) of the inverse delay times V ′ (x) = (x − 1)/(2x), (11) reads We find the condition for the edges a + b = 6, ab = 1, and hence (36). The 1-point Green function is and hence (35). We now compute {m ℓ } and {j ℓ }. For ℓ ≥ 1 we find Inserting the values of a = 3− √ 8 and b = 3+ √ 8 we eventually find (40).

Further results on the cumulants and a conjecture
In the previous section we have obtained the generating functions of all cumulants of power traces of Q. One may be satisfied with the belief that this result tells us everything about the time-delay matrix at leading order in the number of scattering channels N . However, we notice that all the values of α[κ 1 , . . . , κ v ] reported in Table 1 are integer numbers.
In this section we address the questions: are the α[κ 1 , . . . , κ v ]'s integers? and if yes, what do they count? We feel that finding an answer to these questions is both challenging and intriguing. In this section we investigate these issues further.
Lemma 2. Let n be a nonnegative integer. Then (z 2 −6z+1) −n/2 is analytic around z = 0 and the coefficients of the series expansion are nonnegative integers.
Proof. The thesis is a straightforward consequence of the classical identity where P ℓ (t) is the Legendre polynomial of degree ℓ: For t = 3 we write explicitly Hence, for all nonnegative integers n: Theorem 2. Let α[κ 1 , . . . , κ v ] be defined as in (3).
Proof. (i) The first assertion is well-known. (The numbers α[κ] = 1, 2, 6, 22, . . . (κ ≥ 1) enumerate some combinatorial objects, see [23], and they are referred to as the sequence of the 'large Schröder numbers' [38].) (ii) The generating function of α[κ 1 , κ 2 ] is given by (see (48)) and is analytic at (z 1 , z 2 ) = (0, 0). By Lemma 2 the term in the square bracket is analytic about (z 1 , z 2 ) = (0, 0) with integer series coefficients. Note that The products of power series with integer coefficients is again a power series with integer coefficients; therefore, the coefficients α[κ 1 , κ 2 ] of the series expansion of (59) are integers. (Since (59) is analytic, all the coefficients of negative powers of the Cauchy product vanish.) It is possible to perform the Cauchy product in closed form and obtain the explicit expression where This expression shows that α[κ 1 , κ 2 ] is positive.
Remark 4. It is worth remarking that given enough algebraic computation (e.g. by using the algorithm of Appendix C), it is possible to verify that α[κ 1 , κ 2 , . . . , κ v ] are integers for a given large v by the same method.
Theorem 2 and the inspection of the first few values of α[κ 1 , . . . , κ v ] for v > 4 (see Table 1) suggest the following conjecture.
Conjecture. For all v ≥ 1 and κ i ≥ 0 The fact that these cumulants are integer-valued is quite unexpected. Available finite-N formulae [32,34] are cast in a form which is not suitable for a straightforward asymptotic expansion in 1/N , and seem to be of scarce help to solve this conjecture. In addition, it is not clear what (if any) enumeration problem is solved by the family of numbers α[κ 1 , . . . , κ v ] (so far a combinatorial interpretation is available only for α[κ 1 ], see [35]). Some hints may come from semiclassical approaches to quantum transport. In fact, RMT predictions in quantum transport are usually confirmed by independent calculations for β = 1 and β = 2 based on semiclassical considerations. Averages in RMT correspond to sums over pairs of classical trajectories with encounters connecting the leads with the interior of the cavity [8,27,28,44]. In absence of time reversal symmetry (β = 2), the set of classical trajectories can be partitioned according to topological properties. Each class of trajectories has a representative diagram. When computing cumulants (connected averages) one only needs to consider connected diagrams. The contribution c(D) of a diagram D (i.e. the sum over all trajectories in the class represented by D) can be evaluated by counting incoming channels I(D), links L(D) and encounters E(D) inside the cavity containing at most one end point [28]  Assuming the equivalence between RMT and semiclassical methods, the cumulants of the Wigner-Smith matrix C v (T κ 1 , . . . , T κv ) for β = 2 are given by sums over v sets of paired trajectories whose representative diagrams are connected. These sums can be converted into sums over diagrams whose contribution is given by (63). The leading order in N to C v (T κ 1 , . . . , T κv ) comes from diagrams D without loops that can be represented as rooted planar fat trees [8,28,34] (genus-0 maps). From (63) it is clear that this leading order term is integer-valued. In view of the relation (3) between different Dyson indices, this argument shows the plausibility of our conjecture. A more precise mapping between RMT averages and semiclassical sum rules is likely to shed some light on the underlying combinatorics of α[κ 1 , . . . , κ v ]. This semiclassical argument, various finite-N moments [34] for β = 2, and other evidences not reported here [16,29], suggest that the conjecture could be extended to higher order corrections in 1/N . Let us consider a collection of random variables ξ i (i ≥ 1) with finite moments of all order. The v-th cumulant of ξ 1 , . . . , ξ v is defined according to

Acknowledgements
where π denote a generic partition of {1, . . . , v} of size |π|. For instance, the lowest order cumulants are The variables ξ 1 , . . . , ξ v in the definition (64) do not need to be distinct.
The derivative with respect to A is The last derivative requires a bit more work. Starting from we observe that the following decomposition holds with C j = Res The computation of [δJ ℓ /δV (z)] goes along identical lines. Eventually one finds (22)-(23).

Appendix C. Mathematica Code
Below follows a MATHEMATICA code implementing the generating function F v,0 (z 1 , . . . , z v ) of cumulants of the Wigner-Smith matrix. The code provided is tailored for the Wigner-Smith problem but can be easily modified for other problems.
The parameter v determines the highest order in the recursion: the code produces the generating functions F 1,0 , F 2,0 , F 3,0 , . . . , F v,0 up to v = v. The code first implements the Green functions G 1,0 , . . . , G v,0 using the iterative scheme outlined in Section 3. Then the set of the model parameters is assigned and the generating functions are produced.
The parameters β, A, B, M 1 , J 1 , etc. defined by (11) and (15) The next part of the code specifies the model-dependent parameters a, b, {m ℓ } and {j ℓ } denoted by a, b, m, j, respectively, and the one-point function G 1,0 (G [1] in the code). They are hard-coded in the following four lines. This is the only section of the code that depends on the details of the potential V (x).  In the last part of the code the numerical values of a,b,m,j are assigned to the symbols A,B,M1,J1, etc., and the list of generating functions F 1,0 , . . . , F v,0 is constructed according to identity (42).