format: text

text:
# Statistical-mechanical formulation of Lyapunov exponents

## I Introduction

Lyapunov exponents are an important tool for the characterization of dynamical systems. Their very definition has a strong statistical mechanical flavor, as it involves ‘extensivity in time’ of certain quantities regardless of the initial (‘border’) conditions. It is then natural to seek to express them in an explicitly statistical mechanical way, in terms of a partition function.

In this paper, building upon previous work, mainly of Graham Gozzi and Reuter (1994) and Gozzi Graham (1988), we show how the study of Lyapunov exponents, as well as Ruelle’s angle, can be cast as a quantum many-body problem — in fact a rather standard one, except that the Hamiltonian is in general non-hermitian. As we shall see, this does not bring a miraculous solution to all calculational problems, but is serves two purposes:

* Because the problem is formulated as a standard quantum many-body one, all the tools developed in that wider context are available. Some, but not all of them have already been used as approximation schemes for Lyapunov exponents, but others such as Hartree-Fock, dynamic mean-field theory, and the renormalization group seem promising.
* More important, general theoretical results are re-expressed in a language that is often familiar. A typical example is when one asks whether finite-dimensional systems have a limit Lyapunov density function $`\rho (\lambda )`$: In this setting the question becomes whether a quantum finite-dimensional system has an extensive free energy for all chemical potentials. Although this does not in itself prove the existence of $`\rho (\lambda )`$ in the thermodynamic limit, it renders it intuitive and acceptable, at least up to the level of rigor of theoretical physics.

Consider a general dynamical system:

$$\dot{x}_i=f_i(𝐱,\eta );i=1,\mathrm{},N$$
(1)
and, in particular, the version with additive noise:

$$\dot{x}_i=f_i(𝐱)+\eta _i$$
(2)
where $`\eta _i`$ is a Gaussian white variable with variance $`2T_i`$. In the limit of zero noise we have a standard dynamical system and for a particular form of $`f`$ a Hamiltonian system. Clearly, an alternative way to study (2) is to go to the Fokker-Planck (or ‘Kramers’, or ‘Liouville’, depending on the context) description of evolution of probabilities:

$$\dot{P}(𝐱,t)=H_{FP}P(𝐱,t)$$
(3)
where

$$H_{FP}\frac{}{x_i}\left(T\frac{}{x_i}+f_i\right)$$
(4)
Here and in what follows summation of repeated indices is assumed, unless otherwise stated. The Fokker-Planck operator acts on the space of functions of the coordinates $`𝐱`$: it resembles a Schrödinger operator, although it is in general non-Hermitian. The noiseless limit is subtle, and is the subject of ergodic theory.

Introduce now two sets of fermion $`a_i^{}`$, $`b_i^{}`$ and boson $`\alpha _i^{}`$, $`\beta _i^{}`$ creation operators ($`i=1,\mathrm{},N`$); and the corresponding vacuum, defined by:

$$a_i|=b_i|=\alpha _i|=\beta _i|=0;i$$
(5)
It will turn out that all the information we search is obtained directly from the following generalization of the Fokker-Planck operator:

$$H_S=H_{FP}+V_{kl}(𝐱)(a_k^{}a_l+b_k^{}b_l+\alpha _k^{}\alpha _l+\beta _k^{}\beta _l)$$
(6)
where we denote:

$$V_{kl}(𝐱)\frac{f_k}{x_l}$$
(7)
$`H_S`$ acts on the product space of functions of the coordinates and of the number of fermions and bosons of each type. (Clearly, $`H_S`$ coincides with $`H_{FP}`$ when restricted to the zero fermion and boson subspace.) This is the many-body system mentioned above.

If the $`x_i`$ are lattice variables, and the $`f_i`$ are short-range interactions, then the system (6) defines a quantum (non hermitian) theory, also having short-range interactions. If the system is instead off-lattice, the $`x_i`$ describe the position of the particles and the $`(a_i,b_i,\alpha _i,\beta _i)`$ play the role of ‘spin’ degrees of freedom carried by the quantum particles: if the $`f_i`$ are short-ranged, then both the direct and the spin-spin interaction are also short range.

An alternative strategy, that leads to the generalized Lyapunov exponents, is based on replicas of fermions $`a_i^\gamma `$, $`b_i^\gamma `$ with $`\gamma =1,\mathrm{},q`$ and the operator:

$$H_q=H_{FP}+V_{kl}(𝐱)(a_k^\gamma a_l^\gamma +b_k^\gamma b_l^\gamma )$$
(8)
The Lyapunov exponents are defined as follows. The separation of two infinitesimally close trajectories (evolving under the same noise, see Arnold Arnold and Wihstutz (1986)) $`𝐱(𝐭)`$ and $`𝐱(𝐭)+𝐲(𝐭)`$ is described by the evolution of vectors in the tangent space:

$$y_k(t)=U_{ki}(𝐱_𝐨,t)y_i(0)$$
(9)
where $`𝐱_𝐨`$ is the initial condition and $`U`$ is defined as the solution of the linear equation:

$$\dot{U}_{ki}(t)=V_{kj}(x,t)U_{ji}(t)$$
(10)
The $`N`$ Lyapunov exponents $`\lambda _i`$ measure the rate of growth of the volume in the tangent space. We construct

$$A(t)U(t)U^{}(t)$$
(11)
and consider its eigenvalues $`𝒜_1(t)𝒜_2(t)\mathrm{}𝒜_N(t)`$ (a set per initial condition and/or noise realization). The Lyapunov exponents are

$$\lambda _i\underset{t\mathrm{}}{lim}\frac{1}{t}\mathrm{log}𝒜_i$$
(12)
Throughout this paper averages $``$ are over the noise and/or the initial condition. Note that the $`\lambda _i`$ are automatically labeled in a decreasing order $`\lambda _1\lambda _2\mathrm{}\lambda _N`$. Their existence has been proved under very general conditions by Oseledec Oseledec (1968).

One can in fact consider the generalized Lyapunov exponents $`\lambda _i^q`$ (GLE), of interest as a measure of the intermittency (see Benzi et al. Benzi et al. (1985)) defined for example as:

$$\lambda _i^{2q}\underset{t\mathrm{}}{lim}\frac{1}{t}\mathrm{log}𝒜_i^q$$
(13)
(a more standard definition will be given below in (17)).

Adopting the nomenclature of disordered systems, we shall call the true Lyapunov exponents ‘quenched’, and the $`q=1,2,\mathrm{}`$ set the ‘annealed’ exponents. The quenched exponents can be formally obtained as the analytic continuation $`q0`$ (the replica trick), but we shall also consider a direct evaluation here.

Years ago, Graham Graham (1988) constructed a representation of the annealed Lyapunov exponents using fermions valid for systems with one degree of freedom. Later, Gozzi Gozzi and Reuter (1994) derived a generalization to several degrees of freedom. His approach differs from ours in two respects: first, his expression, if averaged over initial conditions or noise, yields the annealed quantities. Second, and also important, he considers the eigenvalues of $`U`$ rather than $`UU^{}`$, a different (though potentially interesting) quantity. In the present treatment we derive an expression for the quenched exponents directly by introducing auxiliary fermions and bosons: the technique of supersymmetry Efetov (1893, 1997)). We also give expressions for the generalized exponents for even $`q`$. In all cases we consider the operator $`UU^{}`$, rather than $`U`$: this is the reason why we need to introduce the fermions and bosons in pairs $`(a_i,b_i)`$ and $`(\alpha _i,\beta _i)`$ (cfr. (6,8)).

Lyapunov exponents appear naturally in the context of localization problems for quantum disordered systems: one considers exponents associated with the product of transfer matrices, which play the role of the evolution operators. (see e.g. Derrida et al. (1987); Weigt and Monasson (1996)) Thus, the Green-function formalism has been treated with supersymmetric (SUSY) techniques to extract correlations (see Balents and Fisher Balents and Fisher (1997), where there is also a detailed analysis of the supersymmetry group and representation, and Markoš Markoš (1989) for a path-integral formalism). In these works supersymmetry was used in order to obtain the Green function, and then the quenched quantities were obtained via the replica trick.

This paper is organized as follows: In section II we present the formalism for the expression of usual and generalized Lyapunov exponents. For Hamiltonian systems we show that the Lyapunov pairing rule follows trivially from a particle-hole symmetry, and we also discuss the computation of Ruelle’s angle. In section III we discuss the statistical mechanics of Lyapunov exponents of macroscopic systems, in particular the existence of a Lyapunov density function. We then present two examples of application: in IV to problems of Random Matrices, and in section V to a Hamiltonian Mean Field model (HMF).

## II The formalism

### II.1 Basic quantities

The average expansion of $`p`$-dimensional volumes evolving with (2) in the space of phases can be expressed, in terms of $`U`$ as:

$$R_p(x_o,t)=\underset{i_1,\mathrm{},i_p}{}det\left[(\xi _{i_1}^{},\mathrm{}..,\xi _{i_p}^{})^{}\left[U^{}U\right](x_o,t)(\xi _{i_1},\mathrm{}..,\xi _{i_p})\right]$$
(14)
where $`\xi _{i_1}\mathrm{}\xi _{i_p}`$ are all sets of $`p`$ orthonormal vectors. Introducing the cumulative Lyapunov exponents $`\mathrm{\Lambda }_i`$ as

$$\mathrm{\Lambda }_i\lambda _1+\mathrm{}+\lambda _i$$
(15)
one has

$$\mathrm{\Lambda }_p=\underset{t\mathrm{}}{lim}\frac{1}{2t}\mathrm{ln}R_p(x_o,t)$$
(16)
The expression (16) is not suited to be written as an integral of an exponential (a partition function) because of the logarithm. This is the usual quenched versus annealed problem in disordered systems: we are interested in the average of the logarithm and not the logarithm of the average. Below we shall overcome it by means of supersymmetry, but we shall first also give expressions for the annealed quantities.

The generalized Lyapunov exponent (GLE) (see Benzi et al. (1985)) can be used to measure intermittency (i.e. rare trajectories having unusual Lyapunov exponents) and are relevant quantities in the so-called “thermodynamic formalism” of chaotic systems. They are defined via the moments of the $`R_p(x_o,t)`$:

$$\mathrm{\Lambda }_p^{2q}\lambda _1^{2q}+\mathrm{}+\lambda _p^{2q}=\underset{t\mathrm{}}{lim}\frac{1}{t}\mathrm{ln}(R_p(x_o,t)^q)$$
(17)
To the extent that the moments $`R_p(x_o,t)^q`$ are sufficient to reconstruct the distribution law of $`R_p(x_o,t)`$ and the averages for all real $`q`$, one can use them to find the quenched average (16). In such cases the quenched quantities can be extracted from (17) using the replica trick:

$$\mathrm{\Lambda }_p=\frac{d\mathrm{\Lambda }_p^q}{dq}(q=0)$$
(18)
### II.2 Generalized exponents

Let us first obtain an expression for the annealed exponents $`\lambda _i^{2q}`$ for integer $`q`$. This calculation is the closest to the construction of Gozzi Gozzi and Reuter (1994), with the important difference that we consider the limit eigenvalues of the matrix $`A=UU^{}`$ and not of $`U`$ itself.

Introducing $`p`$ pairs of fermions $`a_i`$ and $`b_i`$ one can write:

$$R_p(x_o,t)=|\psi _p𝒯\left[e^{_0^tH_1(x(t^{}))𝑑t^{}}\right]\psi _p^{}|$$
(19)
where $`𝒯`$ denotes time-order and

$`H_1(t)`$ $``$ $`V(𝐱)_{ij}(a_i^{}a_j+b_i^{}b_j)`$ (20)
$`\psi _p`$ $``$ $`{\displaystyle \frac{1}{p!}}{\displaystyle \underset{i_1,\mathrm{},i_p}{}}a_{i_1}\mathrm{}a_{i_p}b_{i_1}\mathrm{}b_{i_p}`$ (21)
The explicit time-dependence of $`H_1`$, is given by the evolution of $`𝐱`$ Eq. (2) via $`V_{ij}(𝐱)`$. More generally, one can write, in terms of $`q`$ replicas of the set of fermions $`a_i^\gamma `$ and $`b_i^\gamma `$:

$$R_p^q(x_o,t)=|\psi _p^q𝒯\left[e^{_0^tH_1^q(x(t^{}))𝑑t^{}}\right]\psi _p^q|$$
(22)
where

$`H_1^q(t)`$ $``$ $`V(𝐱)_{ij}(a_i^\gamma a_j^\gamma +b_i^\gamma b_j^\gamma )`$ (23)
$`\psi _p^q`$ $``$ $`\mathrm{\Pi }_{\gamma =1}^q{\displaystyle \underset{i_1,\mathrm{},i_p}{}}a_{i_1}^\gamma \mathrm{}a_{i_p}^\gamma b_{i_1}^\gamma \mathrm{}b_{i_p}^\gamma `$ (24)
(see appendix A). This is not the final expression, since we have not yet imposed that $`𝐱`$ evolves according to the equation of motion.

The use of two sets ($`a`$’s and $`b`$’s) of fermions are needed in order to follow the rate of growth of the norm of vectors in the tangent space and not the evolution of the vectors themselves: we are interested in the eigenvalues of $`U^{}U`$ and not in those of $`U`$. (This is obviously not a problem in dimension one Graham (1988)). In dimensions larger than one, Gozzi Gozzi and Reuter (1994) introduced one family of fermions and thus studied the eigenvalues of $`U`$: these do not in general coincide with the Lyapunov exponents (see references Gordhirsch et al. (1987) and Crisanti et al. (1993) for a discussion), although they surely give relevant information. In section IV we will explicitly show in an example how the $`a`$ and $`b`$-fermions interfere in a non-trivial way.

We now write an expression for the average over noise and/or initial conditions of (22) using the information that $`𝐱(𝐱_𝐨,𝐭)`$ evolves according to (2), and the probability density follows (4). One way to do this is to express the weight of each trajectory as:

$$R_p^q(𝐱_𝐨,t)=Dx\delta (x𝐱(𝐱_𝐨,\eta ,t^{}))|\psi _p^q\left[e^{_0^tH_1^q(x(t^{}))𝑑t^{}}\right]\psi _p^q|_\eta $$
(25)
where $`Dx`$ means the flat functional integral over trajectories and the delta-function imposes that $`x(t)`$ satisfies the equation of motion. This passage is just the standard textbook exercise of going from a Langevin to a Fokker-Planck description (see for example Zinn-Justin (1996)). The result is that the probability evolves through the Fokker-Planck equation, and we have:

$$R_p^q(𝐱_𝐨,t)=1||\psi _p^q\left[e^{tH_R}\right]|\psi _p^q||𝐱_𝐨$$
(26)
where $`H_q=H_{FP}H_1^q`$ is given in (8). The ket $`1|`$ is the flat measure $`1|𝐱=1`$ and $`𝐱|𝐱_𝐨=\delta (𝐱𝐱_𝐨)`$ is the distribution associated to $`𝐱_𝐨`$. Note that the time-ordering is automatic, as the evolution of $`𝐱`$ is taken care of by $`H_{FP}`$. In the limit $`t\mathrm{}`$ the logarithm of $`R_p^q(𝐱_𝐨,t)`$ is given by the smallest eigenvalue of $`H_q`$ in the subspace not orthogonal to the vectors $`1||\psi _p^q`$ and $`\psi _p^q||𝐱_𝐨`$.

We can now describe a practical algorithm for computing the GLE’s. First we have to identify the smallest invariant subspace containing the states appearing to the right and to the left of (26). Clearly, the fermion numbers $`N_a^\gamma =_ia_i^\gamma a_i^\gamma `$ and $`N_b^\gamma =_ib_i^\gamma b_i^\gamma `$ are conserved, and the calculation of $`R_p^q(𝐱_𝐨,t)`$ involves working in the subspace $`N_b^\gamma =N_a^\gamma =p`$. Furthermore, because $`H_q`$ and $`\psi _p^q`$ commute with the operators $`𝒫^\gamma `$ defined by:

$$𝒫^\gamma a_i^\gamma 𝒫^\gamma =b_i^\gamma ;𝒫^\gamma b_i^\gamma 𝒫^\gamma =a_i^\gamma $$
(27)
we should look into the subspace of eigenfunctions symmetric under exchange of $`a`$’s and $`b`$’s, having eigenvalue one under the $`𝒫^\gamma `$: here is where two families of fermions get mixed. Because the quantity (26) is by construction positive, one can see that the eigenvalue having the smallest real part within this subspace has zero imaginary part.

#### II.2.1 Variational calculations

A variational approach can immediately be implemented for the calculation of the generalized (and, with an act of faith, the usual) Lyapunov exponents. Since one is looking for the lowest eigenvalue of $`H_q`$ within a given subspace, one can use a variational trial function $`|\phi `$

$$\mathrm{\Lambda }_i^{2q}\underset{\phi }{\mathrm{min}}\text{Re}\frac{\phi |H_q|\phi }{\phi |\phi }$$
(28)
If the family of variational functions is parametrized for every $`i,q`$, one obtains an approximation which is an explicit function of $`q`$: one can then envisage computing also an approximation for the usual Lyapunov exponents via (17).

### II.3 Ordinary (quenched) Lyapunov exponents

In this section we construct an expression for the (quenched) Lyapunov exponents using a supersymmetry formalism. For the purposes of clarity, we shall do so in two steps: first a naïve calculation that does not take into account the convergence of the sums, but is closer to the standard supersymmetry treatments in other contexts. Next, we use the Borel transform technique to work more properly.

#### II.3.1 Careless calculation

We introduce a set of fermion and boson operators as in (5) and (6), and define the number operators:

$`N_a=a_k^{}a_k`$ ; $`N_b=b_k^{}b_k`$
$`N_\alpha =\alpha _k^{}\alpha _k`$ ; $`N_\beta =\beta _k^{}\beta _k`$
$`N_{bos}={\displaystyle \frac{N_\alpha +N_\beta }{2}}`$ ; $`N_{fer}={\displaystyle \frac{N_a+N_b}{2}}`$
$`\overline{N}=`$ $`N_{fer}`$ $`+N_{bos}`$ (29)
which will commute with all other operators. We shall also need:

$$fa_k^{}b_k^{};\overline{f}\alpha _k^{}\beta _k^{}$$
(30)
Let us now introduce, for a given trajectory $`𝐱(𝐭)`$, the quantity:

$$H_1^S(t)=V_{kl}(𝐱)(a_k^{}a_l+b_k^{}b_l+\alpha _k^{}\alpha _l+\beta _k^{}\beta _l)$$
(31)
With these notations we compute $`Z(\mu ,\overline{\mu })`$ as:

$$|e^{f^{}+\overline{f}^{}}\left(𝒯e^{t(H_1^S\overline{\mu }N_{bos}\mu N_{fer})}\right)e^{f+\overline{f}}(1)^{N_{Bos}}|$$
(32)
an one can easily show that (Appendix B):

$$Z(\mu ,\overline{\mu })=\frac{det[1+e^{\mu t}A]}{det[1+e^{\overline{\mu }t}A]}$$
(33)
where $`A(t)`$ is associated with the trajectory $`𝐱`$. This function will generate all the Lyapunov exponents for the trajectory as:

$$G(\mu )\underset{t\mathrm{}}{lim}\frac{1}{t}\frac{Z(\mu ,\overline{\mu })}{\mu }|_{\overline{\mu }=\mu }=\underset{t\mathrm{}}{lim}\underset{j=1}{\overset{N}{}}\frac{e^{\mu t}𝒜_j}{1+e^{\mu t}𝒜_j}$$
(34)
since for large $`t`$, $`G(\mu )`$ is a ladder with steps at the values of $`\mu `$ corresponding to the $`\lambda _i`$’s (i.e. it is the integral of the Lyapunov distribution function). Just as in the previous subsection, we wish to calculate $`G(\mu )`$. Again, this is directly done in the Fokker-Planck formalism by the quantity:

$$G(\mu )\underset{t\mathrm{}}{lim}\frac{1}{t}\frac{Z(\mu ,\overline{\mu })}{\mu }|_{\overline{\mu }=\mu }$$
(35)
with $`Z(\mu ,\overline{\mu })`$ written as:

$$\varphi _L|e^{t(H_S+\mu N_{fer}+\overline{\mu }N_{bos})}(1)^{N_{Bos}}|\varphi _R$$
(36)
where $`H_S=H_{FP}H_1^S`$ is given in (6). The left and right eigenvectors are

$`|\varphi _R`$ $`=`$ $`{\displaystyle 𝑑𝐱_𝐨P(𝐱_𝐨)e^{f+\overline{f}}||𝐱_𝐨}`$
$`\varphi _L|`$ $`=`$ $`||e^{f^{}+\overline{f}^{}}`$ (37)
where $`P(𝐱_𝐨)`$ is the initial condition distribution. We can also write $`G(\mu )`$ as an expectation value:

$$G(\mu )=\varphi _L|e^{t\{H_S+\mu \overline{N}\}}N_{fer}(1)^{N_{bos}}|\varphi _R$$
(38)
Equation (33) can be expressed as a formal series in powers of $`\mathrm{exp}(t\mu )`$ and $`\mathrm{exp}(t\overline{\mu })`$. This same series is reproduced by (32), or, in averaged version, (36). Each term of the form $`\mathrm{exp}[t(n_1\mu +n_2\overline{\mu })]`$ corresponds to an expectation value in the subspace of $`n_1`$ fermions and $`n_2`$ bosons. Clearly, (33) has a very small convergence radius: this is because the number of bosons, unlike the number of fermions, is unlimited. In other words, expressions (33), (32) and (36) have only a formal meaning. One can still work with them if at the end of a derivation one can resum exactly the series, in which case one has in fact performed the analytic continuation in $`\mu `$ and $`\overline{\mu }`$.

#### II.3.2 Borel transform

Let us now give a more proper construction to the supersymmetric expression of the Lyapunov exponents. It will turn out that the formalism that emerges is not much more complicated than the one in the previous paragraphs.

Briefly, the Borel transform technique consists of going from the formal series:

$$h(y)a_0+a_1y+a_2y^2+\mathrm{}$$
(39)
to its convergent transform

$$h^B(y)=a_0+\frac{a_1}{1!}y+\frac{a_2}{2!}y^2+\mathrm{}$$
(40)
which can then be inverted. In our case, we shall take

$$h(y)a_0+a_1y+a_2y^2+\mathrm{}=\frac{det[1+ye^{\mu t}A]}{det[1+ye^{\overline{\mu }t}A]}$$
(41)
(cfr. Eq. (33)), and define as the Borel-transformed partition function:

$$Z^B(\mu ,\overline{\mu })=h^B(y)|_{y=1}$$
(42)
that is, we are dividing by $`(n_1+n_2)!`$ the term having $`n_1`$ fermions and $`n_2`$ bosons. Repeating the construction above, it is easy to see that:

$$Z^B(\mu ,\overline{\mu })=\varphi _L^B|e^{t(H_S+\mu N_{fer}+\overline{\mu }N_{bos})}(1)^{N_{Bos}}|\varphi _R^B$$
(43)
which is just like (36) except that the left and right eigenvectors are

$`|\varphi _R^B`$ $`=`$ $`{\displaystyle 𝑑𝐱_𝐨P(𝐱_𝐨)g(f+\overline{f})||𝐱_𝐨}`$
$`\varphi _L^B|`$ $`=`$ $`{\displaystyle }||g(f^{}+\overline{f}^{})`$ (44)
where the function $`g(x)`$ is:

$$g(x)=\underset{n=0}{\overset{\mathrm{}}{}}\frac{x^n}{n!\sqrt{n!}}$$
(45)
We can now define as before:

$$G^B(\mu )\underset{t\mathrm{}}{lim}\frac{1}{t}\frac{Z^B(\mu ,\overline{\mu })}{\mu }|_{\overline{\mu }=\mu }$$
(46)
which we can also write as an expectation value:

$$G^B(\mu )=\varphi _L^B|e^{t\{H_S+\mu \overline{N}\}}N_{fer}(1)^{N_{bos}}|\varphi _R^B$$
(47)
Interestingly enough, we can retrieve the information directly (without the need to antitransform) from $`G^B(\mu )`$, since one can show (Appendix C) that:

$$G^B(\mu )=\underset{j=1}{\overset{N}{}}\theta (\lambda _j\mu )$$
(48)
again gives the same ladder function.

#### II.3.3 Symmetries

Although we shall not make here much use of the symmetries in the problem, let us briefly mention them for completeness. Clearly, the Hamiltonian $`H_S`$ is invariant under all the (supersymmetric) transformations rotating simultaneously the $`\alpha _i,\beta _i,a_i,b_i`$ so as to leave the quadratic form $`f+\overline{f}`$ and the fermion-boson part of $`H_S`$ invariant. Then, expectation values can be written in the standard form:

$$O=Tr[(1)^{N_{fer}}𝐂O]$$
(49)
with $`𝐂`$ supersymmetric:

$$𝐂e^{t\{H_S+(\mu +i\pi )(N_{fer}+N_{bos})\}}|\varphi _R^B\varphi _L^B|$$
(50)
Both the original and the Borel transformed versions have the same symmetries.

One can show that supersymmetry is responsible for the fact that $`Z(\mu ,\mu )`$ and $`Z^B(\mu ,\mu )`$ are independent of $`\mu `$: the constancy of the normalization is indeed the underlying reason why we can use the method to obtain quenched averages. (See Ref. Balents and Fisher (1997)) for a discussion in detail of the supersymmetry group and representations).

### II.4 Discussion

The calculation of generalized exponents $`\lambda _i^{2q}`$ is done, as we have seen, by computing the lowest eigenvalue within a subspace of the Hilbert space. This is because the large-$`t`$ limit automatically projects onto the corresponding eigenstate. Approximate and numerical methods for the estimation of the ground state of Schrödinger-like operators abound in the literature, we have already mentioned the variational principle.

If one wishes to extend these results for the quenched exponents, an analytic continuation to $`q0`$ is needed. This is easily done (although with a leap of faith) when an explicit expression for all even $`q`$ is available. Such will be the case in a variational, a perturbative or a mean-field computation.

On the other hand, the supersymmetry method yields the Lyapunov exponents without the need of any continuation. However, there is a price to pay: expressions (38) and (47) involve a sum of terms within subspaces of any number of bosons. In the Borel transformed version (47) this sum is convergent for all finite $`t`$. However, one can see that the largest term corresponds to a boson number of the order of $`\mathrm{exp}(\lambda _1t)`$, and this number grows as we consider larger times. In other words, we can only perform this sum for finite $`t`$, and only then make $`t\mathrm{}`$. Again, this is no problem if an analytic expression is available (perhaps as a result of an approximation), but it does seem problematic to attack a problem numerically this way.

### II.5 Hamiltonian systems: Pairing rule and Ruelle Angle

Damped Hamiltonian systems are a particular case of the dynamics (2), which can be written as

$`\dot{q}`$ $`=`$ $`{\displaystyle \frac{}{p_i}}`$
$`\dot{p}`$ $`=`$ $`{\displaystyle \frac{}{q_i}}\gamma p_i+\gamma ^{\frac{1}{2}}\eta _i`$ (51)
where $`\gamma `$ measures the intensity of the coupling to the bath. The simplectic structure of Hamilton’s equations has consequences.

#### II.5.1 Pairing of exponents

One of the proprieties we can easily infer from our formulation is the pairing-rule for the Lyapunov systems. This pairing rule was proved for a class of dynamical systems (named ‘quasi-hamiltonian’ by Dressler Dressler (1988)); in our formalism it is the result of a particle-hole symmetry. This is most easily seen for the annealed exponents. In the case

$$=\frac{p_i^2}{2}+𝒱(𝐪)$$
(52)
the $`V_{ij}`$ read,

$$V=\left(\begin{array}{cc}0& \mathrm{𝟏}\\ \frac{^2𝒱}{q_iq_j}& \gamma \mathrm{𝟏}\end{array}\right)$$
(53)
Denoting $`x_i=q_i`$ for ($`i=1,\mathrm{},N`$) and $`x_i=p_i`$ for ($`i=N+1,\mathrm{},2N`$), we consider the transformation:

$$\overline{a}_i^{}=ϵ_{ij}a_j\overline{b}_i^{}=ϵ_{ij}b_j$$
(54)
where

$$ϵ=\left(\begin{array}{cc}\mathrm{𝟎}& \mathrm{𝟏}\\ \mathrm{𝟏}& \mathrm{𝟎}\end{array}\right)$$
(55)
(here and in what follows $`\mathrm{𝟎}`$ and $`\mathrm{𝟏}`$ are the null matrix and the identity matrix in the $`N\times N`$ space).

Under this transformation the quantities from (19) transform as:

$`\overline{H_q}H_q2q\gamma (NN_{fer})`$ (56)
$`\overline{\psi _p}|\overline{0}=\psi _{Np}|0`$ (57)
This particle-hole symmetry, a consequence of the symplectic structure of the evolution operator, was noted by Gozzi Gozzi and Reuter (1994). One then gets:

$$\mathrm{\Lambda }_p^{2q}=\mathrm{\Lambda }_{Np}^{2q}+(Np)2q\gamma $$
(58)
which in terms of the exponents becomes

$$\lambda _n^{2q}=2q\gamma \lambda _{Nn+1}^{2q}$$
(59)
Using (18) and (15) we can infer that the Lyapunov exponents are symmetric around $`\frac{\gamma }{2}`$. Note that this result holds for arbitrary $`\eta _i(t)`$, and not only in the case in which it is a white noise.

#### II.5.2 Ruelle’s angle

The difficulty in computing the Lyapunov exponents is due to the fact that, in the tangent space, any instantaneous frame generically turns. In order to define an angle associated with these transformations Ruelle Ruelle (1985) defined a rotation number (or Ruelle frequency -RF- after Dressler Dressler and Lauterborn (1990)) and proved an additive theorem for it. To construct this number we can use a polar decomposition of the evolution operator: $`U=A^{\frac{1}{2}}Q`$ with $`AUU^{}`$ simplectic and symmetric and $`Q`$ simplectic and orthogonal:

$$Q=\left(\begin{array}{cc}X& Y\\ Y& X\end{array}\right)$$
(60)
Clearly, $`TX+iY`$ is unitary; and we can extract a rotation number from $`det(T)=e^{im}`$ (the determinant of an unitary matrix is a pure phase).

The Ruelle angle can be extracted directly from $`U`$ (see Littlejohn (1986)) in the following way <sup>1</sup><sup>1</sup>1The argument presented by Littlejohn Littlejohn (1986) is for closed trajectories but it works the same for open ones.: write the evolution operator as

$$U=\left(\begin{array}{cc}U_A& U_B\\ U_C& U_D\end{array}\right)$$
(61)
Then, $`m(t)`$ will be the argument of $`z(t)=det(U_A+iU_B)`$. We can write:

$$e^{2im(t)}=\frac{z(t)^2}{z(t)z(t)^{}}=\frac{det(U_A+iU_B)^2}{det(U_AU_A^{}+U_BU_B^{})}$$
(62)
where we have used the fact that, since $`U`$ is a real simplectic matrix, $`U_AU_B^{}=U_BU_A^{}`$.

We can use fermions to write the determinant quantities as:

$$det(U_A+iU_B)=|\underset{i=1}{\overset{N}{}}a_{q_i}𝒯e^{tV}\underset{i=1}{\overset{N}{}}(a_{q_i}^{}+ia_{p_i}^{})|$$
(63)
(here we adopt the notation $`a_{q_i}`$ $`a_{p_i}`$ for the fermion corresponding to the $`ith`$ position or momentum degree of freedom ($`i=1,\mathrm{},N`$). We can further use the $`b`$ fermions to obtain the square of this determinant. To compute the denominator of (62), we first note that $`O=U_AU_A^{}+U_BU_B^{}`$ is the top-left block in the $`A=UU^{}`$. To have his determinant in the denominator we can use the relation:

$$\frac{1}{det(1+O)}=|e^{\alpha _{q_i}\beta _{q_i}}Te^{tV}(1)^{N_{bos}}e^{\alpha _{q_i}^{}\beta _{q_i}^{}+\alpha _{p_i}^{}\beta _{p_i}^{}}|$$
(64)
But, since $`O`$ is an $`N\times N`$ block of $`A`$, it will generically have a projection in the direction of the $`N`$ largest eigenvalues, so for large times:

$$det(1+O(t))det(O(t))e^{S_>t}$$
(65)
where $`S_>`$ is the sum of the largest $`N`$ Lyapunov exponents. Putting everything together, we can write

$$e^{2im(t)}=1||e^{\alpha _{q_i}\beta _{q_i}}\underset{k=1}{\overset{N}{}}a_{q_k}b_{q_k}\left[e^{tH_S}\right]\underset{l=1}{\overset{N}{}}(a_{q_l}^{}+ia_{p_l}^{})(b_{q_l}^{}+ib_{p_l}^{})(1)^{N_{bos}}e^{\alpha _{q_i}^{}\beta _{q_i}^{}+\alpha _{p_i}^{}\beta _{p_i}^{}}||𝐱_𝐨$$
(66)
Note that the quantity we have calculated is not quite the quenched averaged angle, but rather the average of the exponential.

## III Macroscopic systems

An active field of research is the information that Lyapunov exponents can provide in extensive systems in thermodynamic limit. The natural object to study is the Lyapunov density function:

$$\rho (\mu )\underset{i=1}{\overset{N}{}}\delta (\lambda _i\mu )$$
(67)
or, better, the cumulative version:

$$C(\mu )_\mu ^{\mathrm{}}\rho (\lambda )𝑑\lambda $$
(68)
In particular $`C(0)`$ is the sum of positive exponents, related to the Kolmogorov-Sinai (KS) entropy. The existence of a thermodynamic limit for the Lyapunov densities has been conjectured (Ruelle (1982); Livi et al. (1986)). A problem immediately arises: the Lyapunov exponents are themselves the result of the limit $`t\mathrm{}`$, and the question as to whether this, and the thermodynamic limits commute is not obvious. Typical examples when they do not is when there are macroscopic motions that take times that diverge with the size. We shall find a clear example in the Hamiltonian mean-field model below (Section V): the particles organize in an anisotropic object which can turn collectively like a rotor: In contact with a bath there is diffusion of the collective angle, but, since the moment of inertia scales with $`N`$, the collective motion is absent if we consider $`N\mathrm{}`$ before $`t\mathrm{}`$. SinaiSinai (1996) has shown that in this last order of limits the densities are well defined in a system of confined particles with pairwise interactions: we shall see that this result is very natural.

In the previous sections we have shown that $`C(\mu )`$ is the large time limit of $`C_t(\mu )`$ with

$$C_t(\mu )=G^B(\mu ,t)=\frac{1}{t}\frac{\mathrm{ln}Z^B(\mu ,\overline{\mu },t)}{\mu }|_{\mu =\overline{\mu }}$$
(69)
(the logarithm in the r.h.s. has no effect, since the normalization is one).

Consider this expression: $`Z^B`$ is a partition function associated with the ‘quantum’ hamiltonian $`H_S`$, where the time plays the role of an inverse temperature, $`\mu `$ of a chemical potential and $`C_t(\mu )`$ of the derivative of a free energy density with respect to the chemical potential (i.e. a particle number per unit volume). The Lyapunov density is hence a form of compressibility.

As mentioned in the introduction, if the original problem is on a lattice and has nearest neighbor interactions, the fermions and bosons are also lattice variables interacting with the nearest neighbor variables (through $`V_{ij}`$). We have then a ‘quantum’ lattice problem with short range interactions. On the other hand, the system could be off-lattice, and the $`x_i`$ be a set of d-dimensional vectors describing the position of the particles interacting via short-range pair forces $`f_i(𝐱)=_jf(x_ix_j)`$. The variables $`(a_i,b_i,\alpha _i,\beta _i)`$ play the role of ‘spin’ degrees of freedom carried by quantum particles, both the direct and the spin-spin interaction are also short range.

All in all, we are asking whether a quantum theory with short range interactions has a good thermodynamic limit with a well defined free-energy density. There is only one non-standard feature if we ask for the $`t\mathrm{}`$ taken before the thermodynamic limit: this is, as we have seen, like asking in a statistical mechanical problem about the zero temperature limit taken before the thermodynamic limit - sometimes a tricky question.

The arguments on extensivity become more subtle if we wish to study the thermodynamic limit of the largest Lyapunov exponent: this is like asking in a particle system not what is the chemical potential needed to create a certain particle density, but rather to create a single particle in the whole system: clearly this is a question of order $`O(1/N)`$. We shall return to this point in section V.

An interesting special case is the behavior of the Lyapunov exponents close to zero in a system with soft modes. In the present context this concerns the properties of a statistical mechanic problem around $`\mu =0`$, i.e. free from external chemical potential. For example, from Eq. (38) we have:

$$\rho (\lambda =0)=\varphi _L^B|e^{tH_S}N_{fer}(N_{fer}+N_{bos})|\varphi _R^B$$
(70)
A very intriguing possibility that immediately comes to mind when working in the present framework is that of studying universality properties in critical points using renormalization-group ideas and techniques.

## IV Random Matrices

In this section we use our formalism to derive some results already obtained for some Random Matrix (RM) models. For brevity we shall only do this in the pure fermion (annealed) case, although the supersymmetric approach can also be applied.

The main lesson we shall obtain is that these systems become, by virtue of the disorder, interacting fermion problems. As such, they can be very well attacked by some of the many methods devised for such cases: Feynman diagrams of course, but also ressumations such as Hartree-Fock.

### IV.1 Weak disorder expansion

The first model where we can show the power of our approach is the one proposed by Derrida & al.Derrida et al. (1987). They study the weak disorder expansion of the quenched Lyapunov exponents for a product of the form:

$$P=\underset{t}{}U_tU_t=B_o+ϵB$$
(71)
where $`B_o`$ is a fixed matrix, $`B`$ is a random matrix, and $`ϵ`$ is a small parameter. We shall study the case with $`B`$ a multi-dimensional, Gaussian white noise, with zero mean (as the finite mean can be safely included in the constant matrix):

$`B_{ij}(t)`$ $`=`$ $`0`$
$`B_{ij}(t_1)B_{kl}(t_2)`$ $`=`$ $`\sigma _{ij,kl}\delta (t_1t_2)`$ (72)
Following Derrida we study the case with the matrix $`B_o`$ having non-degenerate eigenvalues, well separated.

The first step in order to use our formalism is the set-up a continuous-time variant of the problem:

$$U(dt)=1+B_odt+ϵBdt$$
(73)
Next, we exponentiate this expression. Due to the non-continuous character of the random term $`B`$ the correct form of the exponential is:

$$U(dt)=e^{B_odt+Bdt\frac{1}{2}B^2dt^2}$$
(74)
This means that the evolution of $`U`$ will be given by (as in Eq. (10)):

$$V_{ij}(t)=B_{oij}+B_{ij}(t)\frac{1}{2}\underset{k}{}B_{ik}B_{kj}$$
(75)
We can derive now the path-integral form of $`R_1^q`$ as:

$$|\underset{\gamma =1}{\overset{q}{}}a_{i_1}^\gamma \mathrm{}a_{i_n}^\gamma b_{i_1}^\gamma \mathrm{}b_{i_n}^\gamma |e^{tH}|\underset{\gamma =1}{\overset{q}{}}a_{l_1}^\gamma \mathrm{}a_{l_n}^\gamma b_{l_1}^\gamma \mathrm{}b_{l_n}^\gamma |$$
(76)
($`\gamma `$ is the replica index). We used as a Hamiltonian

$$H=V_{ij}(x)\underset{\gamma =1}{\overset{q}{}}(a_i^\gamma a_j^\gamma +b_i^\gamma b_j^\gamma )$$
(77)
The average $`R_1^q`$ will be expressed by integrating the gaussian noise $`B_{ij}`$ in (76) which will transform $`H`$ into

$`H^{}`$ $`=`$ $`H_o+ϵ^2H_I`$ (78)
$`H_o`$ $`=`$ $`B_{oii}{\displaystyle \underset{\gamma =1}{\overset{q}{}}}(a_i^\gamma a_i^\gamma +b_i^\gamma b_i^\gamma )`$ (79)
$`H_I`$ $`=`$ $`{\displaystyle \frac{B_{mn}B_{ij}}{2}}{\displaystyle \underset{\gamma ^{},\gamma =1}{\overset{q}{}}}(a_i^\gamma a_m^\gamma ^{}a_n^\gamma ^{}a_j^\gamma +`$ (80)
$`+b_i^\gamma b_m^\gamma ^{}b_n^\gamma ^{}b_j^\gamma +a_i^\gamma b_m^\gamma ^{}b_n^\gamma ^{}a_j^\gamma )`$
We have used the fact that we can diagonalize $`B_o`$ and work on that base (which will not change the fermionic states). We have also arranged the creators and the destructors of fermions in normal form. The expression (80) can be proved to be the correct one by going back to the Suzuki-Trotter product, performing the averages there, and then reconstructing the continuous version. We are now in the possession of a time-independent Hamiltonian. Note that after integration of the noise we get a non-trivial result only because of the presence of two types of fermions: without them we would have lost the higher moments of the noise. The replica trick is easy to implement for perturbation expansions, since these yield an explicit dependence as a polynomial in $`q`$ of each term of the expansion. This is exactly the program we carry on now. First, the eigenstates of $`H_o`$ are:

$$\psi _o=\underset{\gamma }{\overset{q}{}}\underset{i}{\overset{p}{}}a_i^\gamma b_i^\gamma $$
(81)
and the corresponding eigenvalue is:

$$2q\underset{i}{\overset{p}{}}\epsilon _i$$
(82)
($`\epsilon _i`$ are the eigenvalues of $`B_o`$). The first non-zero order in $`ϵ`$ is in $`ϵ^2`$; in order to compute it we must use first-order perturbation theory to obtain:

$`\psi _o|H_I|\psi _o^{}={\displaystyle \frac{B_{ij}B_{ji}}{2}}(2q)+3q^2{\displaystyle \frac{B_{ii}B_{jj}}{2}}`$ (83)
For the second order perturbation theory (which will give the contribution in $`ϵ^2`$) we must identify the states connected by the perturbation. We can now retain the terms linear in $`q`$, and hence obtain the quenched average (the coefficients of the linear term (see (17)):

$`\epsilon _p`$ $`=`$ $`{\displaystyle \underset{i=1}{\overset{p}{}}}\epsilon _iϵ^2{\displaystyle \frac{B_{ij}B_{ji}}{2}}+`$ (84)
$`+ϵ^4{\displaystyle \underset{i=1}{\overset{p}{}}}{\displaystyle \underset{l=1}{\overset{p}{}}}{\displaystyle \underset{j=1}{\overset{p}{}}}{\displaystyle \underset{m>p}{}}{\displaystyle \frac{B_{ij}B_{mi}B_{jl}B_{lm}}{\epsilon _j\epsilon _m}}`$
$`{\displaystyle \frac{ϵ^4}{2}}{\displaystyle \underset{i=1}{\overset{p}{}}}{\displaystyle \underset{j=1}{\overset{p}{}}}{\displaystyle \underset{n>p}{}}{\displaystyle \underset{m>p}{}}{\displaystyle \frac{B_{im}B_{jn}B_{mi}B_{nj}}{\epsilon _j+\epsilon _i\epsilon _m\epsilon _n}}`$
This expansion is the continuous-time, Gaussian white noise equivalent of the expression obtained by DerridaDerrida et al. (1987).

### IV.2 Parisi-Vulpiani

Another model that can be revisited is the random matrix model introduced by Parisi and Vulpiani Parisi and Vulpiani (198) (see also Crisanti Crisanti et al. (1993)) in order to mimic some systems that show strong chaos. A continuous-time version of this model, in the case of only one spatial degree of freedom ($`N=2`$) is a linearized evolution of the form

$$V(t)=\left(\begin{array}{cc}0& 1\\ \eta (t)& 0\end{array}\right)$$
(85)
where $`\eta `$ is an Gaussian noise with a mean $`r`$ and a deviation $`\sigma `$.

Using this definition we can compute the average of

$$R^q(t)=P(\eta )𝑑\eta \psi ^q|𝒯\left[e^{_0^tH(\eta (t^{}))𝑑t^{}}\right]|\psi ^q$$
(86)
where

$$\psi ^q=\underset{\gamma =1}{\overset{q}{}}(a_q^\gamma b_q^\gamma +a_p^\gamma b_p^q)$$
(87)
and

$$H=\underset{\gamma =1}{\overset{q}{}}(a_q^\gamma a_p^\gamma +b_q^\gamma b_p^\gamma )+\eta (t)\underset{\gamma =1}{\overset{q}{}}(a_p^\gamma a_q^\gamma +b_p^\gamma b_q^\gamma )$$
(88)
This Hamiltonian form of the exponential quantities associated with the Lyapunov exponents can be checked to be true using small-time developments in the Suzuki-Trotter formula, just as in the previous section. We can integrate the noise in (86) and obtain:

$$R(t)^q=\psi ^q|𝒯\left[e^{tH^{}}\right]|\psi ^q$$
(89)
where the averaged Hamiltonian can be again checked in the discretized version:

$`H^{}`$ $`=`$ $`{\displaystyle \underset{\gamma =1}{\overset{q}{}}}(a_q^\gamma a_p^\gamma +b_q^\gamma b_p^\gamma )+r{\displaystyle \underset{\gamma =1}{\overset{q}{}}}(a_p^\gamma a_q^\gamma +b_p^\gamma b_q^\gamma )+`$ (90)
$`+`$ $`{\displaystyle \frac{\sigma ^2}{2}}\left({\displaystyle \underset{\gamma =1}{\overset{q}{}}}(a_p^\gamma a_q^\gamma +b_p^\gamma b_q^\gamma )\right)^2`$
Let us first concentrate on the zero average case $`r=0`$. We make a redefinition of the fermions as:

$`a_qt_1a_q;a_q^{}{\displaystyle \frac{a_q^{}}{t_1}}`$ ; $`a_pt_2a_p;a_p^{}{\displaystyle \frac{a_p^{}}{t_2}}`$
$`{\displaystyle \frac{t_1}{t_2}}`$ $`=`$ $`(\sigma ^2)^{\frac{1}{3}}`$ (91)
In terms of the new fermions, the dependence on $`\sigma `$ becomes explicit:

$$H^{}(\sigma ^2)^{\frac{1}{3}}\left[\underset{\gamma =1}{\overset{q}{}}(a_q^\gamma a_p^\gamma +b_q^\gamma b_p^\gamma )+\frac{1}{2}\left(\underset{\gamma =1}{\overset{q}{}}(a_p^\gamma a_q^\gamma +b_p^\gamma b_q^\gamma )\right)^2\right]$$
(92)
We now have only to diagonalize this time-independent Hamiltonian on subspaces with exactly two fermions of each replica.

We have done this numerically for values of $`q`$ up to 25 within the $`3^q`$-dimensional basis generated by:

$$a_p^\gamma b_p^\gamma ;a_q^\gamma b_q^\gamma ;\frac{1}{\sqrt{2}}(a_p^\gamma b_q^\gamma +a_q^\gamma b_p^\gamma )$$
(93)
(symmetric under Eq. (27)). For small $`q`$ it is easy to diagonalize $`H`$, for example for $`q=1`$ we obtain $`\lambda _1^2=(2\sigma ^2)^{\frac{1}{3}}`$.

Let us note here that the Hamiltonian (90) is, for the case $`q=1`$ the matrix

$$\left(\begin{array}{ccc}0& 0& \sqrt{2}\\ \sigma ^2& 0& \sqrt{2}r\\ \sqrt{2}r& \sqrt{2}& 0\end{array}\right)$$
(94)
essentially the same used by Anteneodo and Vallejos Anteneodo and Vallejos (2002) to derive the same results.

Let us now turn to the case in which the noise has non-zero average. We can transform (90) using transformations (91) with

$$\frac{t_1}{t_2}=\sqrt{r}$$
(95)
The Hamiltonian becomes $`H^{}=\sqrt{r}H^{\prime \prime }`$ with

$`H^{\prime \prime }`$ $`=`$ $`{\displaystyle \underset{\gamma =1}{\overset{q}{}}}(a_q^\gamma a_p^\gamma +b_q^\gamma b_p^\gamma )+{\displaystyle \underset{\gamma =1}{\overset{q}{}}}(a_p^\gamma a_q^\gamma +b_p^\gamma b_q^\gamma )+`$ (96)
$`+`$ $`s^2\left({\displaystyle \underset{\gamma =1}{\overset{q}{}}}(a_p^\gamma a_q^\gamma +b_p^\gamma b_q^\gamma )\right)^2`$
where and $`s^2=\frac{\sigma ^2}{r^{3/2}}`$. There is a crossover (see Lima Lima and Rahibe (1990)) in the Lyapunov exponent dependence between the limits of small and large $`s`$ (see Fig. 1).

Finally, let us give an example on the role played by spatial structures. We consider a one dimensional system with a noise $`\eta `$ (Eq. (85)) random variables with correlations depending on the distance between sites:

$`\eta _{ij}(t)`$ $`=`$ $`0`$
$`\eta _{ij}(t)\eta _{kl}(t^{})`$ $`=`$ $`\sigma ^2(|ij|)(\delta _{ik}\delta _{jl}+\delta _{il}\delta _{jk})\delta (tt^{})`$
We shall compute the annealed Lyapunov exponent ($`\lambda _1^2`$). The Hamiltonian is

$`H={\displaystyle \underset{l=1}{\overset{M}{}}}(a_{q_l}^{}a_{p_l}+b_{q_l}^{}b_{p_l})+{\displaystyle \underset{l,m=1}{\overset{M}{}}}\eta _{lm}(a_{p_l}^{}a_{q_m}+b_{p_l}^{}b_{q_m})`$
After averaging it becomes

$`H={\displaystyle \underset{l=1}{\overset{M}{}}}(a_{q_l}^{}a_{p_l}+b_{q_l}^{}b_{p_l})+{\displaystyle \underset{l=1}{\overset{M}{}}}\sigma ^2(0)(a_{p_l}^{}a_{q_l}+b_{p_l}^{}b_{q_l})^2`$
$`+{\displaystyle \underset{kl=1}{\overset{M}{}}}{\displaystyle \frac{\sigma ^2(|lk|)}{2}}(a_{p_l}^{}a_{q_k}+b_{p_l}^{}b_{q_l}+a_{p_k}^{}a_{q_l}+b_{p_k}^{}b_{q_l}^{})^2`$
The system becomes explicitly translational-invariant. We have to diagonalize this Hamiltonian on the states

$`{\displaystyle \frac{1}{\sqrt{M}}}{\displaystyle \underset{l=1}{\overset{M}{}}}a_{q_l}^{}b_{q_l}^{}`$
$`{\displaystyle \frac{1}{\sqrt{M}}}{\displaystyle \underset{l=1}{\overset{M}{}}}a_{p_l}^{}b_{p_l}^{}`$
$`{\displaystyle \underset{l=1}{\overset{M}{}}}{\displaystyle \frac{1}{\sqrt{2M}}}(a_{q_l}^{}b_{p_l}^{}+a_{p_l}^{}b_{q_l}^{})`$ (100)
where $`M`$ is the total number of sites. Within this subspace, $`H`$ is identical to the matrix (94) with $`r=0`$ and we get:

$$\sigma ^2=M(\underset{i=1}{}\sigma ^2(i)+2\sigma ^2(0))=M\underset{j=1}{}\eta _{ij}\eta _{ij}$$
(101)
again a result obtained before in Refs. Parisi and Vulpiani (198); Anteneodo and Vallejos (2002).

## V Hamiltonian Mean Field system

The Hamiltonian mean field system we consider Antoni and Ruffo (1995) is composed of $`N`$ coupled rotators with a classical hamiltonian:

$$H=\underset{i}{}\frac{p_i^2}{2}+\frac{J}{N}\underset{i,j}{}[1\mathrm{cos}(q_iq_j)]$$
(102)
This system has (in the canonical ensemble) a phase transition at $`T_c=0.5J`$ between a ferromagnetic and a paramagnetic phase Antoni and Ruffo (1995); Latora et al. (1999, 2000)

The equations of motion are, if we allow for noise and dissipation:

$`{\displaystyle \frac{dq_i}{dt}}`$ $`=`$ $`p_i`$
$`{\displaystyle \frac{dp_i}{dt}}`$ $`=`$ $`JM_x\mathrm{sin}q_i+JM_y\mathrm{cos}q_i\gamma p_i+\sqrt{2\gamma T}\eta _i`$
where we have introduced:

$$M_x\frac{1}{N}\underset{i}{}\mathrm{cos}q_i;M_y\frac{1}{N}\underset{i}{}\mathrm{sin}q_i$$
(104)
The matrix $`V_{ij}`$ governing the evolution operator in tangent space is:

$$V=\left(\begin{array}{cc}\mathrm{𝟎}& \mathrm{𝟏}\\ E& \gamma \mathrm{𝟏}\end{array}\right)$$
(105)
with

$`E_{ij}={\displaystyle \frac{J}{N}}\mathrm{cos}(q_iq_j)forij`$
$`E_{ii}={\displaystyle \frac{J}{N}}{\displaystyle \underset{i}{}}\mathrm{cos}(q_iq_j){\displaystyle \frac{J}{N}}=`$
$`=JM_x\mathrm{cos}q_i+JM_y\mathrm{sin}q_i{\displaystyle \frac{J}{N}}`$ (106)
We can use now the function defined in (32) in order to compute the Lyapunov spectrum. The supersymmetry Hamiltonian becomes a sum of single particle terms:

$`H_S`$ $`=`$ $`{\displaystyle \underset{i}{}}H_S^{(i)}`$
$`H_S^{(i)}`$ $``$ $`{\displaystyle \frac{}{p_i}}\left[JM_x\mathrm{sin}q_iJM_y\mathrm{cos}q_i\gamma p_i\right]`$ (107)
$`+T\gamma {\displaystyle \frac{^2}{p_i^2}}p_i{\displaystyle \frac{}{q_i}}J\nu _{q_i}^l\nu _{p_i}^l`$
$`+\gamma \nu _{p_i}^l\nu _{p_i}^l+J(M_x\mathrm{cos}q_i+M_y\mathrm{sin}q_i)\nu _{p_i}^l\nu _{q_i}^l`$
$`J\nu _{p_i}^l\mathrm{cos}q_i𝒞^lJ\nu _{p_i}^l\mathrm{sin}q_i𝒮^l`$
(here no summation over $`i`$ is implied). For compactness we have defined $`\nu ^l`$ for $`l=1,2,3,4`$:

$`\nu _{q_i}^l`$ $``$ $`(a_{q_i},b_{q_i},\alpha _{q_i},\beta _{q_i})`$
$`\nu _{p_i}^l`$ $``$ $`(a_{p_i},b_{p_i},\alpha _{p_i},\beta _{p_i})`$ (108)
and the collective operators:

$$𝒞^l\frac{1}{N}\underset{i}{}\nu _{q_i}^l\mathrm{cos}q_i;𝒮^l\frac{1}{N}\underset{i}{}\nu _{q_i}^l\mathrm{sin}q_i$$
(109)
Not surprisingly, the operator $`H_S`$ is a mean-field (quantum) operator itself, and we can solve the problem with any of the usual methods. For example, introducing explicitly the collective variables using the functional delta-functions:

$`{\displaystyle D[M_x]\delta (NM_x\underset{i}{}\mathrm{cos}q_i)}={\displaystyle D[M_x]D[\widehat{M}_x]e^{\widehat{M}_x(NM_x_i\mathrm{cos}q_i)}}`$
$`{\displaystyle D[M_y]\delta (NM_y\underset{i}{}\mathrm{sin}q_i)}={\displaystyle D[M_y]D[\widehat{M}_y]e^{\widehat{M}_y(NM_y_i\mathrm{sin}q_i)}}`$
$`{\displaystyle D[𝒞^l]\delta (N𝒞^l\underset{i}{}\nu _{q_i}^l\mathrm{cos}q_i)}={\displaystyle D[𝒞^l]D[\widehat{𝒞}^l]e^{\widehat{𝒞}^l(N𝒞^l_i\nu _{q_i}^l\mathrm{cos}q_i)}}`$
$`{\displaystyle D[𝒮^l]\delta (N𝒮^l\underset{i}{}\nu _{q_i}^l\mathrm{sin}q_i)}={\displaystyle D[𝒮^l]D[\widehat{𝒮}^l]e^{\widehat{𝒮}^l(N𝒮^l_i\nu _{q_i}^l\mathrm{sin}q_i)}}`$ (110)
we can write:

$$Z(\mu ,\overline{\mu })=D[M_x]D[\widehat{M}_x]\underset{l}{}D[𝒞^l]D[\widehat{𝒞}^l]D[𝒮^l]D[\widehat{𝒮}^l]e^{N({\scriptscriptstyle 𝑑t(\widehat{M}_xM_x+\widehat{M}_yM_y+\widehat{𝒞}^l𝒞^l+\widehat{𝒮}^l𝒮^l)}W)}$$
(111)
where $`W`$ is the action for a single pair of variables $`(q_i,p_i)`$:

$$e^W=\varphi _L^{(i)}|𝒯e^{t(H_{eff}^{(i)}+\mu N_{fer}^{(i)}+\overline{\mu }N_{bos}^{(i)})}(1)^{N_{bos}^{(i)}}|\varphi _R^{(i)}$$
(112)
(no summation). The left and right eigenvectors are

$`|\varphi _R^{(i)}`$ $`=`$ $`{\displaystyle 𝑑x_iP(x_i)e^{f^{(i)}+\overline{f}^{(i)}}||x_i}`$
$`\varphi _L^{(i)}|`$ $`=`$ $`||e^{f^{(i)}+\overline{f}^{(i)}}`$ (113)
where for simplicity we have assumed that all pairs of variables $`(q_i,p_i)`$ have the same initial distribution, and hence left and right vectors are in product form. The single particle effective Hamiltonian ($`H_{eff}^{(i)}`$) equals:

$$H_S^{(i)}\widehat{M}_x\mathrm{cos}q_i\widehat{M}_y\mathrm{sin}q_i\widehat{𝒞}^l\nu _{q_i}^l\mathrm{cos}q_i\widehat{𝒮}^l\nu _{q_i}^l\mathrm{sin}q_i$$
(114)
which has in principle a time dependence through the collective variables. We wish to evaluate $`<Z(\mu ,\mu ^{})>`$ by saddle-point method in the thermodynamical limit. The saddle point equations for the ordinary collective variables read:

$`\widehat{M_x}`$ $`=`$ $`{\displaystyle \frac{W}{M_x}}M_x={\displaystyle \frac{W}{\widehat{M_x}}}`$
$`\widehat{M_y}`$ $`=`$ $`{\displaystyle \frac{W}{M_y}}M_y={\displaystyle \frac{W}{\widehat{M_y}}}`$ (115)
We also have saddle point equations for the sixteen fermionic and bosonic collective variables $`\widehat{𝒞}^l,𝒞^l,\widehat{𝒮}^l,𝒮^l`$.

We have 20 equations of type (115) which we can solve assuming that the system is at $`t=0`$ already in thermodynamic equilibrium – the vector to the right in (112) is Gibbs-distributed. This implies that the collective variables may be constant in time, an assumption we verify later. First of all, it is easy to see that the saddle-point values of the $`\widehat{𝒞}^l,𝒞^l,\widehat{𝒮}^l`$ and $`𝒮^l`$ vanish (Appendix E). Next, one can see that if the endpoint of the trajectories is left free (i.e. we are not conditioning to a specific arrival point), then causality implies that $`\widehat{M}_x=\widehat{M}_y=0`$ (Appendix E). Under these assumptions, the only variables with non-zero saddle-point values are:

$$M_x=\mathrm{cos}q_i_M;M_y=\mathrm{sin}q_i_M$$
(116)
where the average $`_M`$ is taken with the single particle dynamics:

$`{\displaystyle \frac{dq}{dt}}`$ $`=`$ $`p`$
$`{\displaystyle \frac{dp}{dt}}`$ $`=`$ $`JM\mathrm{sin}q\gamma p+\sqrt{2\gamma T}\eta `$ (117)
Indeed, this yields the equilibrium value of $`M`$, as solution of the equation:

$$M=\frac{I_1(y)}{I_0(y)}y=\frac{MJ}{T}$$
(118)
(where $`I`$ are the Bessel functions, see Antoni and Ruffo (1995)). What we have shown is that the Lyapunov density function of the system is, to leading order in $`N`$ the sum of Lyapunov densities for a single particle system moving according to Eqn. (117).

In the limit of zero coupling to the bath $`\gamma =0`$ the exponents for a single particle are zero, as a consequence of conservation of energy and the pairing rule. Note that the fact that the energy of each particle is conserved separately is an artifact of the large $`N`$-limit: there is a coupling between the fluctuations at the following order in $`1/N`$. In conclusion, in the thermodynamical limit the function $`G(\mu )`$ is a step of height $`2N`$ at zero.

A vanishing largest Lyapunov exponent $`\lambda _i`$ has already been obtained both numerically (Latora et al. (1999, 2000)) and analytically (Anteneodo and Vallejos (2002); Firpo (1998)) in the paramagnetic phase, but not in the ferromagnetic phase. This is not in contradiction with our results: according to our calculation one can still have a vanishing fraction of non-zero Lyapunov exponents, that do not contribute to the density function in the large-$`N`$ limit. If we wish to calculate the largest exponent $`\lambda _1`$ we should take our calculation to the next order.

Figure 2 shows a numerical calculation of the spectrum of Lyapunov exponents for various values of the number of particles $`N`$. The simulations are carried out in the microcanonical ensemble at an energy where the greatest Lyapunov exponent does not show an important variation with $`N`$ ($`E=0.3`$, see Latora et al. (1999)). Even though the computer time requirements do not allow us to get close to the thermodynamical limit, these numerical results show a reasonably good agreement with the hypothesis of a step spectrum.

At finite $`\gamma `$ the Lyapunov distribution function is no longer a step, even in the thermodynamic limit. We can obtain it very easily by solving the one-particle system (117) numerically. The results depend strongly of $`\gamma `$; there is a scaling law in $`\gamma ^{1/3}`$ below the critical temperature and an identical zero $`\lambda _1`$ in the paramagnetic phase (see Fig. 3).

Surprisingly, after rescaling the largest Lyapunov exponent behaves very much like the one computed by Firpo (Firpo (1998)) in the canonical ensemble. There are also qualitative and quantitative resemblances with the one obtained in micro-canonical simulations of a large (but obviously finite) number of particles (see Latora et al. (1999, 2000)).

We can conclude that, at finite $`\gamma `$, $`G(\mu )`$ should have two sigmoids of height $`N`$ located symmetrically around zero.

## VI Conclusions

Writing the Lyapunov exponents as a statistical mechanic object provides a different perspective of the problem, allowing to transpose much of the knowledge and intuition developed in that wider context.

Perhaps the most clear example are the questions related to the thermodynamic limit. We have seen that the existence of a limit Lyapunov density function is the kind of extensivity property that most theoretical physicists would accept without proof — at least for the times that do not diverge with the system size.

As to the different approximation schemes, a case in point is the continuous product of random matrices Section IV. We have seen there that the problem is mapped into a system of interacting fermions. We have treated it as a perturbative expansion of weak interaction reobtaining the results of Ref. Derrida et al. (1987). However, once in the language of interacting particles, one immediately thinks of other, more global approximations; we have already mentioned Hartree-Fock, which can easily be implemented both in the pure fermion or in the supersymmetric formalism. Similarly, the standard mean-field treatment of the model of section V can be extended for the calculation of the exponents. This also suggests that it may be interesting in general to construct a local mean-field approximation for problems with space: this could give a simple analytic handle on the spatial structures involved with each exponent.

A question into which we have not looked in detail is intermittency. The generating function of generalized Lyapunov exponents $`\lambda _i^{2q}`$ involves a distribution function in space that depends on $`q`$. Physically, this arises because we are conditioning the probability of a trajectory to having an unusual value of the exponents. Hence, studying the lowest eigenvector of $`H_q`$ gives us information on the spatial structures responsible for intermittent behavior.

Yet another interesting question is to look at systems at or near criticality, and borrow methods and ideas from the rich theory of critical phenomena to infer results on the behavior of the exponents there.

Acknowledgment

We wish to thank Stefano Ruffo for clarification and encouragement at all stages of the preparation of this work.

## Appendix A

In this appendix we show how the use of fermions transform the quantity $`R`$ defined in (14) can be written as a certain matrix element of the hamiltonian (21) (equation 19)).

Let $`a_i^{}`$ be a set of fermion operators

$$[a_i,a_j^{}]_+=\delta _{ij}$$
(119)
We can encode (14) by writing:

$$O(t)=a_k^{}y_k(t)$$
(120)
We put:

$$\dot{O}(t)=[V_{kl}a_k^{}a_l,O(t)]_{}$$
(121)
It is easy to check that the $`y_i`$ then satisfy (9). The solution to (121) is:

$$O(t)=\left(𝒯e^{tV_{kl}a_k^{}a_l}\right)O(0)\left(𝒯e^{tV_{k^{}l^{}}a_k^{}^{}a_l^{}}\right)^1$$
(122)
which implies:

$$\left(𝒯e^{tV_{kl}a_k^{}a_l}\right)a_i^{}\left(𝒯e^{tV_{k^{}l^{}}a_k^{}^{}a_l^{}}\right)^1=U_{ji}a_j^{}$$
(123)
Let us now add a second set of fermions $`b_k^{}`$, and construct the operator:

$$fa_k^{}b_k^{}$$
(124)
so we can write $`U_{ji}U_{ri}a_j^{}b_r^{}`$ as:

$$\left(𝒯e^{tV_{kl}(a_k^{}a_l+b_k^{}b_l)}\right)f\left(𝒯e^{tV_{k^{}l^{}}(a_k^{}^{}a_l^{}+b_k^{}^{}b_l^{})}\right)^1$$
(125)
(note that the same constructions can be made with bosons instead of fermions.) Denoting $`|`$ the fermion vacuum, Eq. (125) implies:

$`|a_{i_1}\mathrm{}a_{i_p}b_{i_1}\mathrm{}b_{i_p}|𝒯\left[e^{V_{ij}(a_i^{}a_j+b_i^{}b_j)}\right]|a_{l_1}^{}\mathrm{}a_{l_p}^{}b_{l_1}^{}\mathrm{}b_{l_p}^{}|=`$
$`U_{k_1l_1}\mathrm{}U_{k_pl_p}U_{j_1l_1}\mathrm{}U_{j_pl_p}|a_{i_1}\mathrm{}a_{i_p}a_{k_1}^{}\mathrm{}a_{k_p}^{}||b_{i_1}\mathrm{}b_{i_p}b_{j_1}^{}\mathrm{}b_{j_p}^{}|`$ (126)
Summing over $`i_1,\mathrm{},i_p`$ the first part transforms into:

$$U_{k_1l_1}U_{j_1l_1}\mathrm{}U_{k_pl_p}U_{j_pl_p}=(UU^{})_{k_1j_1}\mathrm{}(UU^{})_{k_pj_p}$$
(127)
Denoting with $`\sigma ((j_1,\mathrm{},j_p)(k_1,\mathrm{},k_p))`$ the signature of the permutation between the two sets of indices the second (fermionic) part of (126) equals

$$\sigma ((i_1,\mathrm{},i_p)(k_1,\mathrm{},k_p))\sigma ((j_1,\mathrm{},j_p)(i_1,\mathrm{},i_p))$$
(128)
Using (127) (126) and the properties of the signature it follows that (126) equals

$$(p!)\underset{j_1,k_1,\mathrm{}.j_p,k_p}{}A_{j_1k_1}\mathrm{}A_{j_nk_p}\sigma ((j_1,\mathrm{},j_p)(k_1,\mathrm{},k_p))$$
(129)
(remember that $`A=U^{}U`$). After summing on $`k_1,\mathrm{},k_p`$ the expression (129) is, by definition:

$$(p!)^2det\left[(y_{i_1}^{},\mathrm{}..,y_{i_p}^{})^{}\left[U^{}U\right](y_{i_1},\mathrm{}..,y_{i_p})\right]$$
(130)
## Appendix B

In this appendix we derive formula (32). Denoting $`|`$ the fermion vacuum, Eq. (125) implies:

$$\left(𝒯e^{tV_{kl}(a_k^{}a_l+b_k^{}b_l)}\right)e^f|=e^{A_{jr}a_j^{}b_r^{}}|$$
(131)
and:

$`|e^f^{}\left(𝒯e^{tV_{kl}(a_k^{}a_l+b_k^{}b_l)\mu N_{fer}t}\right)e^f|=`$
$`=|e^{b_ia_i}e^{\mu N_{fer}}e^{\gamma A_{jr}a_j^{}b_r^{}}|=det[1+e^{\mu t}A]`$ (132)
The last equality is most easily checked by performing a rotation of the fermions to diagonalize $`A`$. Then we have:

$`|e^{b_ia_i}e^{\mu N_{fer}t}e^{A_{jr}a_j^{}b_r^{}}|`$ $`=`$ (133)
$`|e^{b_i^{}a_i^{}e^{\mu t}}e^{𝒜_ja_j^{}b_j^{}}|`$ $`=`$
$`{\displaystyle \underset{k}{}}|e^{b_k^{}a_k^{}e^{\mu t}}e^{𝒜_ka_k^{}b_k^{}}|`$
where $`𝒜`$ are the eigenvalues of $`A`$ (no summation on k); in the last expression there is no sum on $`i`$ in the exponents. Each factor can be developed into $`1+e^{\mu t}𝒜_i`$, so the whole expression will be:

$$\underset{i}{}(1+e^{\mu t}𝒜_i)=det[1+e^{\mu t}A]$$
(134)
Consider now the same steps, but now replacing the fermions $`a_k^{}`$ and $`b_k^{}`$ by bosons $`\alpha _k^{}`$ and $`\beta _k^{}`$ as in (30); then we obtain the analogue of (132):

$`|e^{\overline{f}^{}}\left(𝒯e^{tV_{kl}(\alpha _k^{}\alpha _l+\beta _k^{}\beta _l)\mu N_{bos}t}\right)e^{\overline{f}}|=`$
$`=|e^{\beta _i\alpha _ie^{\mu t}}e^{A_{jr}\alpha _j^{}\beta _r^{}}|=det[1e^{\mu t}A]^1`$ (135)
which can again be proved by performing a rotation of the bosons to diagonalize $`A`$. After the same step as in (133) the terms that survive are those with an equal number of the creation and destruction operators and we can use

$$|(\alpha ^n\beta ^n)|=n!$$
(136)
to obtain the analogue of (134):

$$\underset{i}{}(1+e^{\mu t}𝒜_i+e^{2\mu t}𝒜_i^2+\mathrm{})=\underset{i}{}\frac{1}{1e^{\mu t}𝒜_i}$$
(137)
One can change the sign in (137) we can introduce some “negative norm” bosonic states, or as:

$`{\displaystyle \frac{1}{det[1+A]}}=|e^{\overline{f}^{}}\left(𝒯e^{tV_{kl}(\alpha _k^{}\alpha _l+\beta _k^{}\beta _l)}\right)e^{\overline{f}}|`$
$`=|e^{\overline{f}^{}}\left(𝒯e^{tV_{kl}(\alpha _k^{}\alpha _l+\beta _k^{}\beta _l)}\right)(1)^{N_{bos}}e^{\overline{f}}|`$ (138)
## Appendix C

In this appendix we show the ladder structure of the functions $`G`$ defined in (38) and (47). The common ingredient after developing the exponentials and matching the terms allowed by the conservation of bosons and fermions is:

$$|(f^{}+\overline{f}^{})^n\left(𝒯(1)^{N_{bos}}e^{tV_{kl}(a_k^{}a_l+b_k^{}b_l)+tV_{kl}(\alpha _k^{}\alpha _l+\beta _k^{}\beta _l)\mu t(N_{bos}+N_{fer})}\right)N_{fer}(f+\overline{f})^n|$$
(139)
We can commute the “evolution” operator with $`N_{fer}`$ and apply it on the right ket; the action of the operators number of bosons and fermions $`N_{bos}+N_{fer}`$ (on the left ket) can be easily computed so the expression becomes:

$$e^{n\mu t}|\frac{(f^{}\overline{f}^{})^n}{n!}N_F\frac{(𝒜_ia_i^{}b_i^{}+𝒜_i\alpha _i^{}\beta _i^{})^n}{n!}|$$
(140)
The expansion of the powers in terms of individual creators an destructors:

$$e^{n\mu t}(n!)^2\underset{k_1+\mathrm{}+k_N=n}{}|\underset{j}{}\frac{(b_ja_j\beta _j\alpha _j)^{k_j}}{k_j!}N_{fer}\underset{j}{}(𝒜_j)^{k_j}\frac{(a_j^{}b_j^{}+\alpha _j^{}\beta _j^{})^{k_j}}{k_j!}|$$
(141)
Here there is no sum on repeated index inside the parenthesis. Expanding further $`N_{fer}`$ in terms of fermions operators this term can be written as:

$`e^{n\mu t}(n!)^2{\displaystyle \underset{l}{}}{\displaystyle \underset{{\scriptscriptstyle k_i}=n}{}}|{\displaystyle \frac{(b_la_l\beta _l\alpha _l)^{k_l}}{k_l!}}{\displaystyle \frac{a_l^{}a_l+b_l^{}b_l}{2}}𝒜^{k_l}{\displaystyle \frac{(a_l^{}b_l^{}+\alpha _l^{}\beta _l^{})^{k_l}}{k_l!}}{\displaystyle \underset{jl}{}}{\displaystyle \frac{(b_ja_j\beta _j\alpha _j)^{k_j}}{k_j!}}𝒜_j^{k_j}{\displaystyle \frac{(a_j^{}b_j^{}+\alpha _j^{}\beta _j^{})^{k_j}}{k_j!}}|`$ (142)
But for each $`l`$, at $`jl`$ and $`k_j>0`$ a factor of the product above will be:

$`(𝒜_j)^{k_j}({\displaystyle \frac{|\beta _j^{k_j}\alpha _j^{k_j}\alpha _j^{k_j}\beta _j^{k_j}|}{k_j!^2}}{\displaystyle \frac{k_j^2|\beta _j^{k_j1}\alpha _j^{k_j1}\alpha _j^{k_j1}\beta _j^{k_j1}||b_ja_ja_j^{}b_j^{}|}{k_j!^2}}))`$ (143)
and we can easily see that this equals zero. So, inside each term of the sum on l, only that with $`k_l=n`$ and the others $`k_j=0`$ will survive.

The sum (141) simplifies to:

$$e^{n\mu t}(n!)^2\underset{l}{}\frac{𝒜_l^n(1)^{n1}n^2|\beta _l^{n1}\alpha _l^{n1}\alpha _l^{n1}\beta _l^{n1}||b_la_l\frac{a_l^{}a_l+b_l^{}b_l}{2}a_l^{}b_l^{}|}{n!^2}=e^{n\mu t}\underset{l}{}𝒜_l^n(1)^{n1}$$
(144)
Now we can reconstruct the sum on n for the different cases; the Borel construction gives uses the series

$$\underset{j=1}{\overset{N}{}}\underset{n=1}{\overset{\mathrm{}}{}}\frac{(e^{\mu t}𝒜_j)^n}{n!}=\underset{j=1}{\overset{N}{}}[1exp(e^{\mu t}𝒜_j)]$$
(145)
The series from (145) are convergent at each time and for each trajectory (in a well behaving, smooth potential); one can formally sum the general terms (144) without the $`\frac{1}{n!}`$ and obtain:

$$\underset{j=1}{\overset{N}{}}\underset{n=1}{\overset{\mathrm{}}{}}(e^{\mu t}𝒜_j)^n=\underset{j=1}{\overset{N}{}}\frac{e^{\mu t}𝒜_j}{1+e^{\mu t}𝒜_j}$$
(146)
Using the exponential form of $`𝒜_j`$, in the limit $`t\mathrm{}`$, both (145) and (146) will give

$$G(\mu )=G^B(\mu )=\underset{j=1}{\overset{N}{}}\theta (\lambda \mu )$$
(147)
## Appendix D

In this appendix we derive the fourth term in $`ϵ`$ of the weak disorder expansion treated in section (IV.1). All the eigenstates of $`H_o`$ can be obtained by destroying and creating fermions in $`\psi _o`$. Due to the fact that the perturbation contains two creators and two destructors the states connected by the perturbation with $`\psi _o`$ are only those that differ from it by one or two fermions. The “replica composition” of the states connected by the perturbation must be the same (two different replicas are not connected).

One set of states is formed by destroying one fermion from the first $`p`$ and creating one (from the last $`Np`$); both fermions must be in the same replica family:

$`\psi _{b\gamma jm}=b_j^\gamma b_m^\gamma \psi _0`$
$`\psi _{a\gamma jm}=a_j^\gamma a_m^\gamma \psi _0`$
$`jpm>p`$ (148)
These eigenstates will contribute to the Lyapunov exponent by a second order perturbation term:

$`{\displaystyle \underset{\gamma =1}{\overset{q}{}}}{\displaystyle \underset{j=1}{\overset{p}{}}}{\displaystyle \underset{m>p}{}}{\displaystyle \frac{\psi _{b\gamma jm}|V|\psi _0\psi _0|V|\psi _{b\gamma jm}}{\epsilon _j\epsilon _m}}+`$ (149)
$`{\displaystyle \underset{\gamma =1}{\overset{q}{}}}{\displaystyle \underset{j=1}{\overset{p}{}}}{\displaystyle \underset{m>p}{}}{\displaystyle \frac{\psi _{a\gamma jm}|V|\psi _0\psi _0|V|\psi _{a\gamma jm}}{\epsilon _j\epsilon _m}}`$ (150)
The linear part in q of this term is:

$$2q\underset{i=1}{\overset{p}{}}\underset{l=1}{\overset{p}{}}\underset{j=1}{\overset{p}{}}\underset{m>p}{}\frac{B_{ij}B_{mi}B_{jl}B_{lm}}{\epsilon _j\epsilon _m}$$
(151)
The states that have two different fermions are of three types:

$`\psi _{b\gamma in;b\gamma ^{}jm}=b_i^\gamma b_j^\gamma ^{}b_m^\gamma ^{}b_n^\gamma \psi _0`$
$`\psi _{a\gamma in;b\gamma ^{}jm}=a_i^\gamma b_j^\gamma ^{}b_m^\gamma ^{}a_n^\gamma \psi _0`$
$`\psi _{a\gamma in;a\gamma ^{}jm}=a_i^\gamma a_j^\gamma ^{}a_m^\gamma ^{}a_n^\gamma \psi _0`$
$`ipjpm>pn>p`$ (152)
and each of them will contribute with terms of the form

$$\underset{\gamma =1}{\overset{q}{}}\underset{j=1}{\overset{p}{}}\underset{m>p}{}\underset{\gamma ^{}=1}{\overset{q}{}}\underset{i=1}{\overset{p}{}}\underset{n>p}{}\frac{\psi |V|\psi _0\psi _0|V|\psi }{\epsilon _j+\epsilon _i\epsilon _m\epsilon _n}$$
(153)
The linear term in q will be:

$$q\underset{i=1}{\overset{p}{}}\underset{j=1}{\overset{p}{}}\underset{n>p}{}\underset{m>p}{}\frac{B_{im}B_{jn}B_{mi}B_{nj}}{\epsilon _j+\epsilon _i\epsilon _m\epsilon _n}$$
(154)
## Appendix E

In this appendix we study the consistency equations for the collective variables in the HMF model. As we can see in from the equation (34) we need quantities at $`\mu =\overline{\mu }`$; the derivation implied by the definition of the $`G`$ function will be carried out on $`<Z(\mu ,\overline{\mu })>`$ after the saddle point evaluation. So, we must take the derivative of the expression (111). In this point, the derivative reads:

$$\frac{Z^S(\mu ,\overline{\mu })>}{\mu }|_{\overline{\mu }=\mu }=\frac{1}{Z^S}\left[\frac{W(\mu ,\overline{\mu })}{\mu }|_{\overline{\mu }=\mu }+\underset{k}{}\frac{Z^S(\mu ,\overline{\mu })>}{𝒳_k}\frac{𝒳_k}{\mu }|_{\overline{\mu }=\mu }\right]$$
(155)
where we take $`𝒳_k`$, for $`k=1,\mathrm{},20`$, to be the vector of the collective variables, and $`Z^S`$ is the value of $`<Z(\mu ,\overline{\mu })>`$ at the saddle (111). We see that the second part of the rhs of this equation is zero because of the consistency equations (115); in conclusion we are interested in the values of the collective variables at $`\mu =\overline{\mu }`$ and only the direct dependence on the two variables of $`W`$ will be interesting for us. With this assumption we write the consistency equation for the boson and fermion variables as:

$$𝒮^l=\frac{\varphi _L^{(i)}|𝒯e^{t(H_{eff}^{(i)}+\mu \overline{N})}\nu _{q_i}^l\mathrm{sin}q_i(1)^{N_{bos}^{(i)}}|\varphi _R^{(i)}}{e^W}$$
(156)
$$\widehat{𝒮}^l=\frac{\varphi _L^{(i)}|𝒯e^{t(H_{eff}^{(i)}+\mu \overline{N})}\nu _{p_i}^l\mathrm{sin}q_i(1)^{N_{bos}^{(i)}}|\varphi _R^{(i)}}{e^W}$$
(157)
$$𝒞^l=\frac{\varphi _L^{(i)}|𝒯e^{t(H_{eff}^{(i)}+\mu \overline{N})}\nu _{q_i}^l\mathrm{cos}q_i(1)^{N_{bos}^{(i)}}|\varphi _R^{(i)}}{e^W}$$
(158)
$$\widehat{𝒞}^l=\frac{\varphi _L^{(i)}|𝒯e^{t(H_{eff}^{(i)}+\mu \overline{N})}\nu _{p_i}^l\mathrm{cos}q_i(1)^{N_{bos}^{(i)}}|\varphi _R^{(i)}}{e^W}$$
(159)
All those equations contain expectation values of single boson and fermion operators which necessarily vanish. This is normal, as we could expect from the beginning that those variables integrate out, contributing only with a non-exponential prefactor. This prefactor (as expected from the supersymmetry considerations ) will be one in the limit $`\mu =\overline{\mu }`$ but in the process of derivation will give terms of order $`O(\frac{1}{N})`$.

There are four equations left can be divided in two parts. The first parts contains the variables $`\widehat{M}_x`$ and $`\widehat{M}_y:`$

$$\widehat{M}_x=J\frac{\varphi _L^{(i)}|𝒯e^{t(H_{eff}^{(i)}+\mu \overline{N})}(\frac{}{p_i}\mathrm{sin}q_i+\mathrm{cos}q_i\nu _{p_i}^l\nu _{q_i}^l)(1)^{N_{bos}^{(i)}}|\varphi _R^{(i)}}{e^W}$$
(160)
$$\widehat{M}_x=J\frac{\varphi _L^{(i)}|𝒯e^{t(H_{eff}^{(i)}+\mu \overline{N})}(\frac{}{p_i}\mathrm{sin}q_i\mathrm{sin}q_i\nu _{p_i}^l\nu _{q_i}^l)(1)^{N_{bos}^{(i)}}|\varphi _R^{(i)}}{e^W}$$
(161)
The first term of the rhs of these equations is zero (as demanded by the causality) and the second is also zero, due to the mismatch in the fermion-boson operators.

Finally, the last two equations are the only nontrivial ones:

$$M_x=\frac{\varphi _L^{(i)}|𝒯e^{t(H_{eff}^{(i)}+\mu \overline{N})}\mathrm{cos}(q_i)(1)^{N_{bos}^{(i)}}|\varphi _R^{(i)}}{e^W}$$
(162)
$$M_y=\frac{\varphi _L^{(i)}|𝒯e^{t(H_{eff}^{(i)}+\mu \overline{N})}\mathrm{sin}(q_i)(1)^{N_{bos}^{(i)}}|\varphi _R^{(i)}}{e^W}$$
(163)
Now the hamiltonian $`H_{eff}^{(i)}`$ is much simpler:

$$H_{eff}^{(i)}=\frac{}{p_i}\left[JM\mathrm{sin}q_i\gamma p_i\right]+T\gamma \frac{^2}{p_i^2}p_i\frac{}{q_i}J\nu _{q_i}^l\nu _{p_i}^l+\gamma \nu _{p_i}^l\nu _{p_i}^l+JM\mathrm{cos}q_i\nu _{p_i}^l\nu _{q_i}^l$$
(164)
where we used also the rotational symmetry (which define a saddle manifold) on the space of $`M_x`$ and $`M_y`$ which allow us to fix $`M_y=0`$ and $`M_x=M`$; looking back at (111) and (155) we can conclude that

$$G(\mu )=N\underset{t\mathrm{}}{lim}\frac{1}{te^W}\frac{W(\mu ,\overline{\mu })}{\mu }|_{\overline{\mu }=\mu }$$
(165)
Formula (112) combined with the simple expression of $`H_{eff}^{(i)}`$ (164) allow us to infer that

$$G(\mu )=NG^{(i)}(\mu )$$
(166)
where $`G^{(i)}(\mu )`$ characterize the Lyapunov spectrum of a single particle with the dynamics (117); this spectrum contains in fact only two exponents.


text:
# Hanbury Brown Twiss effects in channel mixing normal-superconducting systems

## 1 Introduction

In normal mesoscopic multiterminal conductors, the cross correlation between currents flowing into different terminals are manifestly negative. As shown in Ref. , this is a consequence of the fermionic statistics of the electrons. Recently, such negative correlations were observed experimentally , in an electrical analog  of the Hanbury Brown Twiss experiment  with photons.

The manifestly negative cross correlations were derived under the assumption of noninteracting electrons. In contrast, long range Coulomb interactions between electrons can give rise to positive cross correlations between currents at capacitive contacts . Moreover, in Ref.  it is demonstrated that for inelastic scattering between edge states in a multiterminal conductor in the quantum Hall regime, the current cross correlation can be positive. Another situation where positive cross correlations are possible is when the normal conductor is connected to a superconducting reservoir .

In a normal-superconducting system, electrons can be converted into holes via Andreev reflection at the interface between the normal conductor and the superconductor. The Andreev reflection thus introduces correlations between electrons and holes in the normal conductor, a phenomena known as the proximity effect . The influence of the proximity effect on the current auto-correlations, i.e. the shot noise, in two-terminal normal-superconductor junctions has been studied in Refs. .

In multiterminal conductors, Andreev reflection can lead to positive cross correlations between currents flowing in the contacts to the normal reservoirs . This was originally predicted for single mode conductors . In contrast, a discussion of a very asymmetrical geometry in which one of the contacts is a tunneling tip suggested negative correlations in the presence of channel-mixing scattering . Interestingly, in Ref.  it was shown that in a multiterminal diffusive normal-superconducting system with perfect interfaces between the normal conductor and superconducting reservoir, the cross correlations are manifestly negative in the absence of the proximity effect. Positive correlations can be enforced with the help of ferromagnetic contacts  a case which we will not discuss here further.

Clearly, these very different predictions and their sensitivity to the specific choice of parameters and geometry demand a more systematic investigation of the role of channel-mixing and the proximity effect. To investigate the role of channel mixing Boerlin et al.  have analyzed a geometry with diffusive tunnel junctions. Independently, the authors  present a systematic analysis of the current correlations in a system consisting of a chaotic quantum dot connected via point contacts to one superconducting and two normal reservoirs. Using the scattering approach for normal-superconducting systems  in combination with a random matrix theory  description of the chaotic dot has the advantage of allowing us to analyze both the case with and without proximity effect on the same footing. The generic properties of the model also makes our results qualitatively relevant for multiterminal normal-superconducting structures with random scattering.

## 2 The model

The three terminal dot-superconductor junction is shown in Fig. 1. A chaotic quantum dot (see Ref.  for definition) is connected to two normal reservoirs ($`N_1`$ and $`N_2`$) and one superconducting reservoir ($`S`$) via quantum point contacts. The contact to the normal reservoir is perfect, the contact to the superconducting reservoir has a mode independent transparency $`\mathrm{\Gamma }_S`$.

The contacts support $`N`$ and $`M`$ transverse modes respectively. The conductance of the point contacts are much larger than the conductance quanta $`2e^2/h`$, i.e. $`N,M\mathrm{\Gamma }_S1`$, so Coulomb blockade effects in the dot can be neglected. The two normal reservoirs are held at the same potential $`V`$ and the potential of the superconducting reservoir is zero. The magnetic flux in the dot is $`\mathrm{\Phi }`$.

## 3 Current cross correlations

Due to the random scattering in the dot, the current $`I_i(t)`$ in contact $`i`$ fluctuates around its quantum statistical average $`\overline{I_i}`$. We calculate the zero-frequency spectral density of the current cross-correlations

$$P_{12}=2𝑑t\overline{\mathrm{\Delta }I_1(t)\mathrm{\Delta }I_2(0)},$$
(1)
where $`\mathrm{\Delta }I_j(t)=I_i(t)\overline{I_i}`$. The correlation $`P_{12}`$ can be expressed  in terms of the scattering matrix $`S`$ of the total system, dot and superconductor. To proceed $`S`$ can be expressed  in terms of the scattering matrix $`S_d`$ of the dot and the Andreev reflection amplitude at the contact-superconductor interface.

We consider the limit of zero temperature and a potential $`eV`$ much lower than the inverse dwell time of the dot, $`E_{Th}`$, where the energy dependence of the scattering matrix $`S_d`$ can be neglected. We also assume $`eV\mathrm{\Delta }`$, where $`\mathrm{\Delta }`$ is the superconducting gap. This restricts the transport to an energy interval where no quasiparticles can escape into the superconductor. In this regime, the scattering matrix $`S`$ describes only scattering between the normal reservoirs.

Noting that the current fluctuation is just the sum of the fluctuations of electron and hole currents, the noise power can be conveniently written ,

$$P_{12}=P_{12}^{ee}+P_{12}^{hh}+P_{12}^{eh}+P_{12}^{he}$$
(2)
where $`P_{12}^{\alpha \beta }`$ is the correlation between $`\alpha `$ and $`\beta `$ quasiparticle currents. Unlike the cross correlations in purely normal systems , which are manifestly negative, the cross-correlation $`P_{12}`$ can be positive, because the correlations between different types of quasiparticles, $`P_{12}^{eh}+P_{12}^{he}`$, are positive.

The ensemble averaged correlations $`P_{12}`$ are calculated using the statistical properties of the scattering matrix $`S_d`$ of the dot . The details of these calculations are presented in Refs. , here we just show the result in Fig. 2 for various transparencies $`\mathrm{\Gamma }_S`$ of the dot-superconductor contact.

Two different regimes for the magnetic flux are considered. For a flux much smaller than the flux quanta, $`\mathrm{\Phi }h/e`$, the magnetic flux has no effect on the proximity effect and can be completely neglected. In the opposite regime, a magnetic flux larger than the flux quantum, $`\mathrm{\Phi }h/e`$, effectively breaks the time reversal symmetry in the dot and suppresses the proximity effect .

We see from Fig. 2 that suppressing the proximity effect has a very different effect in the two limits; dominating coupling of the dot to the normal reservoirs, $`M\mathrm{\Gamma }_SN`$, and to the superconducting reservoir $`M\mathrm{\Gamma }_SN`$. In the limit of strong coupling to the superconducting reservoir, $`M\mathrm{\Gamma }_SN`$, suppressing the proximity effect with a magnetic flux in the dot completely suppresses the positive correlations.

In this limit, in the presence of the proximity effect, there is a gap in the spectrum around the Fermi energy in the dot . As a consequence, quasiparticles injected from one normal reservoir are Andreev reflected, effectively direct at the contact-dot interface , with unity probability back to the same reservoir. This is shown in Fig. 3 for a perfect contact to the superconducting reservoir, $`\mathrm{\Gamma }_S=1`$, however, the scattering probabilities (as a function of $`2N/M\mathrm{\Gamma }_S`$) are essentially independent of $`\mathrm{\Gamma }_S`$. For $`2N/M\mathrm{\Gamma }_S0`$, the scattering process is thus deterministic, there is no partition of incoming quasiparticles and hence no noise, $`P_{12}=0`$.

Increasing the coupling to the normal reservoir, the probability of normal reflection as well as cross Andreev reflection, from one reservoir to the other, becomes finite. As is clear from Fig. 3, the normal reflection is the dominant process of the two. On a formal level, we can thus neglect the terms in $`P_{12}`$ containing the cross Andreev reflection amplitude. As seen in Fig. 3, this gives $`P_{12}=2P_0(N/M\mathrm{\Gamma }_S)`$, positive since only terms in $`P_{12}^{eh}+P_{12}^{he}`$ contribute. Here $`P_0=4e^3/h`$.

Suppressing the proximity effect, the gap in the spectrum is suppressed and the scattering processes are strongly modified. This is clearly seen by comparing the two figures in Fig. 3, showing the scattering probabilities for $`\mathrm{\Gamma }_S=1`$. In the limit of $`2N/M0`$, an injected electron, which in the presence of the proximity effect was directly back reflected as hole with unity probability to the same contact, now has a probability $`1/4`$ each to leave the dot in the following four ways: as an electron into contact 1 or 2 or a hole into contact 1 or 2.

As a consequence, the correlations become negative, i.e. there are no positive correlations in the absence of the proximity effect, similar to what was found for a metallic diffusive system . Fig. 2 demonstrates that this is independent on the transparency $`\mathrm{\Gamma }_S`$ of the barrier in the dot-superconductor contact.

In the opposite limit, with a weak coupling of the dot to the superconductor, $`M\mathrm{\Gamma }_SN\mathrm{\Gamma }_N`$, the cross correlations are independent of the proximity effect, as seen in Fig. 2. This is also observed for the different scattering probabilities in Fig. 3. An analytical calculation of the correlations in this limit gives the simple result

$$\frac{P_{12}}{P_0}=\frac{M}{2N}R_{eh}(12R_{eh})$$
(3)
where $`R_{eh}=\mathrm{\Gamma }_S^2/(2\mathrm{\Gamma }_S)^2`$ is the Andreev reflection probability of quasiparticles incident in the dot-superconductor contact. There is a crossover from negative to positive correlations, which occurs for $`R_{eh}=1/2`$, i.e $`\mathrm{\Gamma }_S=2(\sqrt{2}1)0.83`$. In this limit of weak coupling to the superconducting reservoir, where the proximity effect plays no role, it is possible to explain the current correlations with straightforward particle counting arguments.

## 4 Particle counting model

In the limit of weak coupling to the superconductor, a quasiparticle injected from one of the normal reservoir will at the most scatter once at the contact between the dot and the superconductor (with a probability $`R_{eh}`$ to Andreev reflect), before leaving the dot. Noting that an Andreev reflection at the dot-superconducting interface effectively leads to an injection of a pair of particles (electrons or holes), we can equivalently say that in the limit of weak coupling to the superconductor, none of the particles in such an injected pair will ever return to the dot-superconductor contact, but instead directly and independently leave the dot into one of the normal reservoirs, with probability $`T_1=T_2=1/2`$ for each reservoir. This process is shown schematically in Fig. 4.

Based on this, we can, using particle counting statistics arguments, calculate the probability $`P(N_1,N_2,N_p)`$ that $`N_1(N_2)`$ electrons have ended up in the normal reservoir $`1(2)`$, when $`N_p`$ pairs have tried to enter the dot. From the probability $`P(N_1,N_2,N_p)`$, we can calculate (the details are presented in Ref. ) the cross correlations, giving

$`P_{12}=N_1N_2N_1N_2R_{eh}(12R_{eh}).`$ (4)
This is just the same as was found in Eq. (3), using a completely quantum mechanical approach.

As a comparison, considering the injection of single particles (see Fig. 4), corresponding to the situation when the superconducting reservoir is replaced by a normal reservoir, we instead get manifestly negative cross correlations.

From this we can conclude that in the limit of a weak coupling to the superconductor, the current correlations are determined by the Andreev reflection at the dot-superconductor interface together with independent partitioning of particles inside the dot.

## 5 Conclusions

In conclusion, we have investigated the relation between current cross correlations and the proximity effect in a three terminal chaotic-dot superconductor junction. We find that both the sign and magnitude of the correlations can be changed by a weak magnetic field in the dot, suppressing the proximity effect. In the limit of weak coupling of the dot to the superconducting reservoir, the current correlation can be found by particle counting arguments.

We acknowledge discussions with H. Schomerus, Y. Nazarov and W. Belzig. This work was supported by the Swiss National Science Foundation and the program for Materials with Novel Electronic Properties.