Large-deviations approach to thermalization: the case of harmonic chains with conservative noise

We investigate the possibility of characterizing the different thermalization pathways through a large-deviation approach. Specifically, we consider clean, disordered and quasi-periodic harmonic chains under energy and momentum-conserving noise. For their associated master equations, describing the dynamics of normal modes energies, we compute the fluctuations of activity and dynamical entropy in the corresponding biased ensembles. First-order dynamical phase transition are found that originates from different activity regions in action space. At the transitions, the steady-state in the biased ensembles changes from extended to localized, yielding a kind of condensation in normal-modes space. For the disordered and quasi-periodic models, we argue that the phase-diagram has a critical point at a finite value of the disorder or potential strength.


Introduction
Thermalization, the process by which a system reaches a state of thermal equilibrium with its surroundings, stands as one of the fundamental phenomena governing the behavior of diverse physical systems.From quantum particles in condensed matter systems to the dynamics of stars in astrophysical environments, the concept of thermalization underpins our understanding of equilibrium states and relaxation processes across various scales.Particularly fascinating is the key role of relaxation timescales on the very possibility of observing signature of the past events and to store and retrieve memory of them [1].
On general grounds, the problem is usually approached from two main sides.In the context of nonlinear dynamics, emphasis is on how energy spreads from a particular point or region in phase space towards the largest portion of it where energy equipartition among all degrees of freedom occurs.This is the familiar approach pioneered by Fermi, Pasta, Ulam and Tsingou, which has a long-standing history [2] but is still actively investigated.Closeness to integrability usually affects strongly the equilibration in a non-trivial way [3,4,5,6,7,8,9].The situation is even more intriguing for systems with quenched disorder, where and interplay between Anderson localization and chaotic diffusion is at work [10,11,12,13,14].
The second approach is more statistical, in the sense that it rests on some form of stochastic dynamics and on how relatively unlikely, low-entropy states evolve towards maximal entropy one.This á la Boltzmann approach involves equations for distributions derived under suitable assumptions.For instance, in the case of nonlinear lattices this leads to Boltzmann-Peierls [15] and wave turbulence equations [16] .For one-dimensional chains this turned out to be useful to understand anomalous transport [17,18,19,20] and the effect of nonlinear resonances on (pre)thermalization [21,22].
To calculate the probability of these trajectories and how atypical ones arise dynamically one can rely on large deviation theory [23,24,25].In the last decades, such approach has rapidly grown as an effective tool to study fluctuations in statistical mechanics and stochastic processes.For dynamics, it provides a mean to compute the probability of untypical trajectories and value of observables that depart significantly from their average, typical ones.Dynamical phase transitions, namely non-analiticities of the large-deviation functions, may occur as a signal of qualitative changes in the statistical properties.
Our aim here is to explore the possibility to use this approach to explore the thermalization dynamics of many-body classical systems.The idea would be to quantify and characterize regions in phase space that yield possibly slow, or at least sizeably different from average, relaxation pathways.In principle, this should provide insight on the different phase-space trajectories followed during equilibration.
A similar idea has been discussed in [26] (see also the account in [27]).There it has been shown that the large deviation of an observable may have singularities even in the context of the Gaussian model, a paradigmatic non-interacting system, before and after an instantaneous temperature quench.The associated phenomenon of condensation of fluctuations is thus a characteristic feature of the thermalization process.Remarkably, condensation originates from the duality between large deviation events in the given system and typical events in a "biased" one.We also mention a conceptually related work on large deviations for wave-turbulence equations [28].We mostly have in mind the relaxation of isolated systems: as a matter of facts, largedeviation approaches have been reported for nonequilibrium steady states close to what we will deal here, like harmonic networks in contact with thermal baths (see e.g [29,30,31] and references therein).
Since dealing with genuine non-linear forces is a formidable task, it is convenient to adopt an intermediate "mesoscopic" approach whereby nonintegrable interactions are replaced by a random noise.This idea has proved to be useful, for instance, to study the Lyapunov time-scale [32,33] and the onset of diffusive dynamics in almostintegrable systems [34].Also, such approach reproduces many features of nonlinear lattices in steady nonequilibrium states [35,36].This semplification may overlook some of the subtle correlations of the deterministic dynamics but can provide also quantitative results, as for instance the scaling law of the Lyapunov exponent with the energy density for the Fermi-Pasta-Ulam-Tsingou model [33].In this spirit, we focus on a specific model: an harmonic network with momentum and energy conserving noise that we studied in [37].Such model has a the remarkable advantage that, under suitable approximations, it allows to write explicitly a linear master equation describing energy spreading in action space.We will exploit this feature to compute the large deviation properties exactly, through standard methods.
The paper is organized as follows.In Section 2 we present our model system: an harmonic network subject to momentum and energy conserving noise, which is the main object of study.We first briefly recall the kinetic approach presented in [37] that allows to write a master equation for diffusion in action space which is the starting point of the subsequent analysis.In Section 3 we briefly recall the methods used to compute the cumulant generating function of the main observables, the dynamical activity and the dynamical entropy.The first results concerns the simplest case of a chain with translational invariance discussed in Section 4. In the second part we examine the thermalization and large deviation properties of two paradigmatic models, the disordered (Section 5) and quasi-periodic chain (Section 6) that correspond, respectively, to the famous Anderson and Aubry-André -Harper models in the context of electronic transport.

Continuous-time master equation
Let us first consider a general master equation for the probabilities P ν of a state, ν = 1, . . ., N with transition rates R ν,ν Eq. ( 1) can be rewritten as where we introduced the escape rate from state ν as and 1/r ν the mean residence time in state ν.The calculation of the large-deviation properties for (1) is well established and will be briefly reviewed in the next Section.Before this, we now specify the physical context that we will deal with.Let us consider the following quadratic Hamiltonian, where (q i , p i ), i = 1, . . .N are the canonically-conjugated coordinates ad momenta (for simplicity, we will deal henceforth with the equal-mass case, and set m i = 1).The coupling matrix Φ is semi-positive definite and symmetric so that its eigenvalues are real and non-negative.The equations of motion of the isolated network are In order to have momentum conservation it must be i Φ ij = 0.In this general formulation, the model encompasses several different setups.It includes the standard case of periodic lattices, where Φ is the familiar nearest-neighbor Laplacian matrix, but also the case of disordered structures with random and/or sparse connections as in the case of elastic networks [38,39].It can be regarded as an idealized model of, say, a large atomic cluster or a protein in its native state.
The normal modes' coordinates (Q ν , P ν ) of the network have eigenfrequencies ω ν and are defined by the transformations The χ ν l , ν = 1, . . ., N , are orthonormal and can be taken to be real for the time being.The main observables are the mode energies The classic problem amounts to ask how a (possibly weak) nonlinear interaction potential added to H brings the harmonic modes to equipartition.For conservative dynamics, upon expressing the Hamiltonian in the familiar action-angle variables, a nonintegrable perturbation defines a network of couplings among the unperturbed actions [40,41].The connectivity of such a network will determine relaxation and ergodic properties.In general, one can distinguish such networks depending on whether the number of groups of actions linked by the perturbation scales intensively or extensively on the number of degrees of freedom N .These are termed short or long-range networks, respectively, depending on whether the coupling range if roughly constant or increases proportionally to N [40,41].
In the nonlinear case the coupling may involve three or more actions, resulting in complicated higher-order interaction on a hyper-graph structure.In view of the difficulty of the problem, it is of interest to formulate a hybrid dynamics for Hamiltonian (4) defined by adding on top of the linear deterministic dynamics a conservative stochastic process [42] defined as follows.Consider the network state at time t.At a later time t + τ there occur a "collision" event whereby a couple of particles (m, n), m = n is randomly selected according to the joint probability w n,m and their momenta are exchanged, (p n , p m ) → (p m , p n ).This move conserves energy and momentum.Physically, it can be interpreted as a perfectly elastic rebound as it would occur for a infinite square-well pair potential.The w n,m are usually assigned to describe local interactions (e.g.nearest-neighbors) but may include long-range terms interaction across the lattice as well.The intervals between subsequent collision times τ are also taken as random variables with some preassigned distribution with given, finite, average τ .
In the kinetic approximation [37] the (angle-averaged) energy E ν = E ν of each normal mode obeys the equation This can be regarded as an averaged collision operator coarse-grained over a time scale which is longer than the typical collision time.We work with finite γ = 1 N τ , that corresponds to choosing a finite collision probability per site (as usually employed in the simulations [36,43,44,45]).The form of equation (7) guarantees that ν E ν remains constant as it should.One major advantage is that, the simple form of the collision rule allows us to compute explicitly the transition rates in terms of the eigenvectors [37].First, define V as the a unit vector of components Note that V = V (n,m) is a random vector as the indexes n, m are chosen at random with the prescribed rule.We refer the reader to [37] for the details and just mention the formulas needed henceforth.It is found that where the overline is a shorthand notation for the average over w n,m of the random collisions.Note that R ν,ν ′ = R ν ′ ,ν with the escape rate from state ν as For a generic coupling matrix Φ, such that there is no decoupling among different subsets they are all nonvanishing.They are of the form of the well-known Fermi golden rule, involving the squared amplitudes of the eigenmodes.Formally, the quantities P ν = E ν / E ν (comprised between 0 and 1) satisfy (7).A physical interpretation would be to imagine an large ensemble of quasi-particles with a given distribution among modes, each P ν being the instantaneous fraction of them in mode ν.Altogether, the dynamics is a random walk of quasiparticles in action space.Transitions are determined by the eigenstate structure of the coupling matrix Φ, and not assigned a-priori.
Following the usual arguments of Markov processes, the system is ergodic and approaches microcanonical equilibrium where P ν = 1/N (equipartition) according to the properties of the "collision operator" defined by ( 1) and ( 9).So the thermalization problem reduces to computing its N eigenvalues µ ν , ν = 1, . . ., N .whose absolute values give the spectrum of relaxation rates.The first eigenvalue is µ 1 = 0 and corresponds, as said, to the steady state of energy equipartition.On physical grounds, µ 1 must be non-degenerate since we expect the dynamics to be ergodic under the above hypotheses.All the others µ ν must be strictly negative [46].In particular, the long-time relaxation rate is controlled by spectral gap |µ 2 − µ 1 | = |µ 2 |, µ 2 being referred to as Fiedler eigenvalue in the context of diffusion on graphs.For a finite network, the scaling with N of µ 2 can be expressed in terms of the scaling of the spectral density.More precisely, denoting by ρ(µ) the integrated (cumulative) density of eigenvalues (i.e. the fraction of eigenvalues less than µ), if then d µ is the spectral dimension of the action network.We emphasize that this is a distinct quantity with the standard spectral dimension of the Laplacian matrix Φ entering in the spectral problem for the Hamiltonian (4).The solution of master equation with initial condition P ν (0), is fairly straightforward.Upon expanding in the basis of the eigenvectors ψ (ν) of the matrix W (the stochastic matrix is symmetric so there is no need to distinguish left and right eigenvectors), one can express the solution as This can be easily solved numerically to compute for the observable of interest.

Large-deviations approach
In this section we briefly review the methods employed in the following.The material is standard in the related literature, and we mostly follow reference [47].

Activity
One quantity studied in the literature is the activity, which counts the number of jumps in a trajectory i.e. the number K of configuration changes during a time t.
As the name implies, it is a measure of how active a given trajectory is: an active trajectory has many changes of configuration, an inactive one just few or none.
The relevant operator to study is the tilted generator which for the activity is defined by the N × N matrix with elements [48,49] ( ) Note that his operator has −r ν on the diagonal and the first terms appears only in off-diagonal entries (see for example eqs.( 44) and (71) in Ref. [47]).
Following the general strategy, the large-deviation statistics is encoded in the leading eigenvalue and (right and left) eigenvectors of the tilted operator (also called scaled cumulant generating function).More precisely, we denote by λ K (s, N ) the largest eigenvalue of W K (s) along with its the right eigenvector defined by W K (s)φ K = λ K φ K .Upon normalization, φ K can be considered as a steady probability distribution for the "biased" dynamics (sometimes termed s-ensemble).As it is known [47] for s = 0 the eigenvector yields the density in the steady state, which in our case is the equipartition which is homogeneous.For s < 0, it probed the regime in which the mean activity K/t of histories is typically larger than in the steady state.They correspond to explored configurations where the density is larger than the steady state distribution.On the other hand, at s > 0 histories with smaller K/t are favored.
As an order parameter we may thus consider (see eq. ( 177) in [47]).
that, in our context, measure the occupancy of the steady state and thus of the distribution of modes that are involved in the biased dynamics.

Entropy production
Another observable of interest is associated to dynamical complexity, Q + = ln P rob(history).In the dynamical systems literature the corresponding largedeviation function is known as the topological pressure.It can be computed in a similar way through the operator W + (s) (see e.g.Eq. ( 36) in Ref. [47]) ) and its the maximal eigenvalue is the topological pressure, large deviation of the entropy production.Following the same prescription as above, we seek for the leading eigenvalue λ + (s, N ) with right eigenvector φ + , W + φ + = λ + φ + .An order parameter ρ + (s) can be defined as in (14) with φ +,ν replacing φ K,ν .The average Q + /t is the Kolmogorov-Sinai entropy h KS = dλ + /ds(0); another important quantity is the topological entropy h top = λ + (1) which is the growth rate of the number of possible histories with time.

Translationally-invariant chains
Let us start from one-dimensional chains with translational invariance i.e. the case where Φ is a circulant matrix.The simplest example would be the standard nearestneighbor coupling, but the results apply to any circulant Φ, including the case of harmonic long-range interactions [50].In fact, in all those instances, the eigenvectors are the familiar lattice Fourier modes and the vector V can be given explicitly where k ν = 2πν N (we use the labeling ν = −N/2 + 1 . . .N/2 throughout this Section).Thus the transition rates in formula ( 9) are known exactly once w n,m is given.Furthermore, for any choice of Φ, the rates only depend on w n,m .For definiteness, we will focus on the case of nearest-neighbor random exchange, namely m = n + 1.This is the mostly studied choice in the literature both for short (nearest-neighbor) [42,51] than for long-range interactions [52].In this case, the transition rates can be computed explicitly since [37] |V (n,m) 16) Note that the escape rate r ν is different for each state ν.Another interesting property is that the off-diagonal terms are of order 1/N , and thus the eigenvalues of the matrix W are well approximated by the diagonal elements and the νth eigenvector is localized on ±ν with all other components being small of order 1/N .The relaxation rate of a generic non-equilibrium initial condition occurs, at long enough times, with a rate |µ 2 | ≈ γ(π/N ) 2 .Thus, the matrices W K and W + are also known analytically: for instance the off-diagonal elements of W + are Let us first discuss the case of the activity, In Figure 1a we report the numerically computed λ K (s) for different sizes.It is shown that that large deviation function has a discontinuity in the derivative at s = 0 corresponding to a dynamical first-order phase transition.Indeed, increasing s leads to a sudden jump in the typical activity dλ K /ds, which corresponds to a dramatic change in the number of configurations explored by histories with activity rate K/t (see the inset of Figure 1a) .As far as the system size scaling is concerned, for s < 0 λ K is N − independent, while for positive s it vanishes as 1/N 2 (see inset in Fig. 1a).This finite-size analysis provide a convincing evidence that the discontinuity persists in the thermodynamics limit.
Taking into account the fact that plane waves come in pairs with opposite wavenumber, it is appropriate to adapt the definition of the order parameter (14) as (and the same for ρ + (s)).This gives indication about the typical wavenumbers which are mostly active along a biased trajectory.Figure 1b shows that ρ K has a discontinuous jump at the transition point s = 0, indicating that different subset of degrees of freedom are involved in the biased dynamics.This is confirmed by inspecting the structure of the leading eigenvector φ K , which is markedly different it the two regimes, see Figures 1c: • at s = 0 the eigenvector φ K is uniform, as expected since it must correspond to the steady state of energy equipartition; • for s > 0, φ K has two large components at ±k 1 plus a uniform background of size O(1/N ); moreover, this shape it is independent of s; • for s < 0, φ K is very small at the band center k ν = 0 were it attains a minimum and increasingly concentrates toward zone-boundary Fourier modes k ν ≈ ±π upon decreasing s.
Thus, trajectories in the inactive phase s > 0 consists mostly the longest-wavelength mode which remains excited, undergoing very little exchange with others.On the other hand, the active phase s < 0 is mostly peaked on short-wavelengths and involve rapid energy exchange among a finite fraction of them.This gives an overall picture of action space, which is consistent with the numerical simulation of the microscopic dynamics [37].The remarkable fact that the transition occurs at s = 0 means that the model dynamics displays "phase coexistence" between fast and slow relaxing modes.To complement the analysis, we also performed a calculation of the large deviation of the dynamical entropy, Figure 2 reports the same type of analysis as above.There is still a similar transition at s = 0, separating two regimes: one "regular" with zero KS entropy and one "chaotic" with finite (perhaps N dependent) entropy.So the natural dynamics present a coexistence of these two.Since the standard KS entropy h KS is given by the the slope at the origin, from the data in the inset of Figure 2a it is seen that for s → 0 + , h KS ≈ 2γ log N , meaning that it subextensive in the system size.It is worth noticing that such scaling occurs in other instances, like for instance the infinite-range Ising model (see Section 6 in [47]).This may be related to the fact that the coupling among modes is of mean-field type, see again eqs.( 16).This feature is different from genuinely chaotic dynamical systems, where typically the KS entropy is extensive, h KS ∼ N .
At this point, we remark that there is indeed a similarity with kinetically constrained models, as Fredrickson-Andersen dynamics and the like, which have been thoroughly studied as toy models of glasses [48,25].Also there, the large-deviation function is singular at s = 0 and the order parameter has a first-order jump.The interpretation is that there are two phases, an active one for s < 0 and an inactive one for s > 0. Physical dynamics take place at s = 0, where the two dynamic phases coexist.It our case the same interpretation applies in action rather than in physical space.

Disordered chain
Let us now turn to the simplest non-homogeneous lattice: the harmonic chain with disorder in the pinning potential with σ i being i.i.d.variables with uniform distributions in [0, w], gauged by the disorder strength parameter w.Periodic boundary conditions are assumed.As it is well known, the eigenstates are the exponentially-localized Anderson modes, whose localization length decreases with w [53,54].Several variants in this type of model with conservative noise have been considered earlier [55,56,57].We consider the standard case of nearest-neighbor collisions and individual realizations of the quenched disorder {σ i }.
Before discussing large deviations, we illustrate some features that can be inferred from the kinetic equations (see also [37] for a comparison with numerics).The localization level affect drastically the connectivity in action space as measured by the transition rates.Indeed, the R ν,ν ′ are all non-vanishing, but can be exponentially small due to the small overlap of the corresponding eigenvectors, see again Eq. ( 16).In Figure 3a we plot the spectral gap as a function of disorder.There is a size-dependent disorder value, separating fast relaxation regime from a slower one.Qualitatively, this correspond to a change in the structure of the matrix W from fully connected (as in the ordered case) to almost-banded.Above such value, the spectral gap vanishes as 1/N 2 , see the inset in Figure 3a.
A key indicator that has been often employed to characterize thermalization is the non-equilibrium (Shannon) spectral entropy [58] For large t, S approaches its equipartition value log N .In general S will depend on the initial conditions.In Figure 3b,c we report the evolution of S(t) for two values of the disorder strength and different sizes obtained by numerical solution of the master equation, starting from the ground-state excitation P ν (0) = δ ν,ν0 .
(ν 0 = 1 in the calculations).There is a crossover from an initial exponential to an intermediate logarithmic growth.The initial regime can be accounted for by a "mean-field" approximation [59]: assuming that energy of the initially excited mode is evenly redistributed among the other for ν = ν 0 yielding the approximation where E ν0 (t) ≈ (1 − 1/N ) exp(−µ ν0 t) + 1/N .This approximation works well for longrange action networks i.e. for weak disorder.The second regime is readily understood as a one-dimensional diffusive process, where the number of excited modes grows as a √ t.This is in agreement with the known result that for diffusion on graphs the Shannon entropy grows as ds 2 log t, if the spectrum of the corresponding Laplacian operator shows scaling with a finite spectral dimension d s [59].It also fits with the observed wavepacket spreading law in nonlinear disordered chains in the strongly chaotic regime [12,60] where indeed microscopic chaos can be effectively approximated by a stochastic process.
Let us now turn to the calculation of the large deviation properties.In Figure 4a we report data for the activity.There is again a signature of a first-order dynamical transition as in the translationally-invariant case and thus of an active and an inactive phase.However, the critical point is now shifted to a nonzero value s = s c .As for the finite-size analysis, note that the curve collapse upon dividing by N , signaling extensivity, at variance with the ordered case.Moreover, the data in Figure 4a show that for N > 10 3 the curves collapse very accurately also close to s c .Both the activity dλ K (s)/ds (inset of Figure 4a) than the order parameter ρ K (s) (Figure 4b) display the same qualitative behavior, with a step discontinuity at s = s c .Moreover, we observed that s c increase upon increasing the disorder strength w.The shift of the critical point implies that the unbiased dynamics at s = 0 has a finite activity rate.The structure of the eigenvector φ K displays a qualitative change at s = s c from a rather homogeneous profile to a localized one (Figure 4c).
In Figure 5, we report the calculation for the entropy.The data confirms that the dynamical transition is shifted to a finite s where dλ + (s)/ds and the order parameter ρ + (s) undergo a finite jump for large enough N .A noteworthy difference with the respect to the ordered case is that the eigenvalue is proportional to N in the whole s range.Thus also the slope at the origin is proportional to N , meaning that in this case there is a finite Kolmogorov-Sinai entropy density h KS .This feature suggest a form of more homogeneous "chaos" in presence of heterogeneity.Finally, the shape of the eigenvector φ + displays a qualitative change at s = s c from a rather homogeneous profile to a localized one (Figure 4c), similar to what observed for the activity.
What is the dynamical phase diagram of the model?We have seen that the phase transition persist upon increasing w.On the other hand, for large disorder the  analysis in [37] and the data in Figure 3c suggest that, for fixed N , the thermalization basically amounts to a diffusive process on a network with finite spectral dimension (equal to one).In this limit, the problem should thus amount to a continuous-time random walk and we do not expect any dynamical phase transition.Indeed, for a such a process λ + (s) should be roughly Poissonian with no singularity [47].
Based on the above limiting cases, we may argue that the transition line s c (w) in the (s, w) plane should terminate at a critical point.In Figure 6a-d we show λ K (s) and λ + (s) divided by N for two different values of the disorder strength.The data confirm the expectation: for w = 10 the curves have smooth derivatives and there is no sign of criticality.In Figure 6e we draw a schematic phase diagram for the activity based on the above consideration.A similar diagram should be valid for the entropy.

Quasi-periodic chain
To further test the approach and the above results, let us consider the chain with a quasi-periodic on-site potential Here, ∆ > 0 is strength of the local force and the parameter β measures the degree of incommensurability of the potential with the lattice period and ϕ in an arbitrary phase.The eigenstate problem can be recast in the form (where ǫ = ω 2 + 1) which is readily recognized as the celebrated Aubry-André-Harper model [61,62,63].A well-known feature of the model is that it displays a localization transition at ∆ = 2, where the eigenmodes change from extended to exponentially localized.For a finite, periodic chain of length N a typical choice is to set β as a rational approximant, for instance the ratio of consecutive Fibonacci numbers, such that for large N β approaches the golden mean.
To compare with the disordered chain, in Figure 7a we report the spectral gap as a function of ∆.For ∆ < 2 the eigenvectors are extended and the spectral gap remains finite and independent of the size.Instead, as soon as localization sets in for ∆ > 2, |µ 1 | vanishes as 1/N 2 (see the inset of Figure 7a) and thermalization slows down considerably.The typical growth of the spectral entropy follow the same two basic routes as above: S(t) grows exponentially or logarithmically in time for the two cases respectively (see Figures 7b-c) and saturates at the equipartition value.
Thus, as far as average properties are concerned, the quasiperiodic lattice behaves very much as the disordered chain, the main difference being that the 'closing' of the  spectral gap (and the change from exponential to logarithmic growth of S) occurs sharply at ∆ = 2, independently of the lattice size.
As for the large deviation properties, we checked that upon increasing the strength of the potential ∆ the behavior is definitely similar to the disordered chain.For It is also of interest to study the behavior of the the large deviation functions across the localization transition ∆ = 2.In particular the numerical solution shows that the dynamical phase transition does not occur around this value.We just discuss the case of the entropy, as the activity displays similar features.For instance, in Figure 8a,b we report the scaled λ + for two values of ∆ below and above the transition.There is no sign of discontinuity neither in λ + nor in the order parameter ρ + , as shown in Figures 8c,d which only shows a smooth monotonic growth in the extended phase ∆ < 2 (Figure 8c) and a maximum around s = 0 in the localized phase ∆ > 2 (Figure 8d).Another difference is in the eigenvectors profiles that do not display a sharp change from extended to localized as in the disordered chain (compare Figure 8e-h with the lower panels of Figure 5).Actually, for ∆ < 2, φ + changes gradually from a peaked distribution for negative s to a uniform one for s > 0, albeit with some structure (Figure 8e-f).Instead, for ∆ > 2, φ + remains extended and fluctuates wildly (Figure 8g-h), displaying a fluctuating pattern which may reflect the underlying structure of the eigenfunctions.

Conclusions
This work is an attempt to employ large-deviations methods to get insights on the thermalization of classical many-body systems.Exploiting the mathematical simplicity of conservative noise dynamics, we are able to compute the large deviation functions for two relevant quantities, the activity and dynamical entropy.For both observables, there is evidence of a dynamical phase transitions that correspond to trajectories of different qualitative nature.We thus distinguish different region in action space corresponding to different thermalization pathways.
An appealing way to interpret the transition from extended to localized in the leading eigenvectors, is to describe it as a form of condensation [26,27] occurring in action/mode space.In other words, most trajectories (a macroscopic fraction) in the localized phase concentrate on a subset of action space.The structure of the eigenvector gives information on the actual modes mostly involved in the dynamics.
When formulated in terms of the kinetic equations, the thermalization of such system is basically a random walk on a large network in action space.Every move conserves the total energy and only changes the entropy and the graph is undirected.From this point of view, there is a close connection with the general problem of large deviations of observables of random walks evolving on random graphs.This issue has been investigated recently [64,65].Indeed, sudden changes in degree fluctuations, similar to dynamical phase transitions, are related to localization transition are found that may be connected with what discussed here.

Figure 1 .
Figure 1.Large deviations of activity for the translationally invariant chain with conservative noise, γ = 1.(a) the largest eigenvalue of W K (s) as a function of s, for different chain lengths N ; the inset in an enlargement around s = 0 showing the eigenvalues multiplied by N 2 ; (b) the order parameter ρ K (s); (c) the eigenvectors φ K for N = 256 as a function of the wavenumber kν : left panel s − 0.005, −0.01, −0.0500.1,middle s = 0, right s = 0.1.The data show evidence of a dynamical phase transition at s = 0 at which the activity displays a first-order jump and φ K changes from extended to localized.

Figure 2 .
Figure 2. Large deviation function for the dynamical entropy for the translationally invariant chain with conservative noise, γ = 1.(a) The largest eigenvalue of W + (s) as a function of s for different chain lengths N ; the inset in an enlargement around s = 0 showing the eigenvalues divided by ln N and the dashed line has a slope 2γ; (b) the parameter ρ + (s); (c) the eigenvectors φ + for N = 256 as a function of the wavenumber kν : left panel s = −0.1,middle s = 0, right s = 0.005, 0.01, 0.050, 0.1 .

Figure 3 .
Figure 3. Thermalization of the disordered chain with conservative noise γ = 1:(a) the spectral gap as a function of the disorder strength w, the inset shows that µ 2 is proportional to 1/N 2 for w large; (b,c) time evolution of the spectral entropy S(t) for two values of w corresponding to weak and strong disorder respectively.Data from numerical solution of the master equation starting from the groundstate mode ν = 1.The dashed line is the law S(t) = 1 2 log t + const..

Figure 4 .
Figure 4. Large deviations of activity for the disordered chain with conservative noise w = 2, γ = 1: (a) the rescaled largest eigenvalue of W K (s) as a function of s for different chain lengths N ; In the inset we report the negative of the numerical derivative −λ ′ K (s)/N computed by finite differences.(b) the order parameter ρ K (s); (c) the eigenvectors φ K for N = 2048.The data show evidence of a dynamical phase transition at s = sc ≈ 0.0376 at which the activity displays a first-order jump and the leading eigenvector of W K changes from extended to localized.

Figure 5 .
Figure 5. Large deviations of entropy production for the disordered chain with conservative noise w = 2 γ = 1: (a) the largest eigenvalue of W + (s) as a function of s for different chain lengths N ; the dashed line has a slope 2γ.In the inset we report the the numerical derivative λ ′ K (s)/N computed by finite differences.(b) the parameter ρ + (s); (c) the eigenvectors for N = 1024 as a function of the lattice site number left s = −0.1,middle s = 0, right s = 0.1 evidence of transition at s ≈ −0.059 at which the eigenvector changes from extended to localized.

Figure 6 .
Figure 6.Disordered chain with conservative noise: comparison of the scaled large deviation functions for two values of the disorder strength w = 2, 10 (left and middle panels respectively).(a,b) scaled largest eigenvalue λ K and (c,d) λ + .The dashed line in panel (d) is a fit with a Poissonian law yielding 0.38(e 5.0s − 1); (e) schematic phase diagram for the activity in the (s, w) plane where the continuous line is the transition line from active to inactive state.The red dot marks the critical point where the transition disappear.

5 Figure 7 .
Figure 7. Thermalization of the the quasi-periodic chain with with conservative noise starting from the ground-state mode ν = 1, γ = 1; (a) the spectral gap as a function of the potential strength ∆, inset shows that µ 2 is proportional to 1/N 2 for ∆ > 2. (b,c) time evolution of the spectral entropy S(t) for two values ∆ < 2 a and ∆ > 2 respectively.The dashed line is 1 2 log t.In all calculation we set β as the ratio of consecutive Fibonacci numbers and ϕ = π/5.

Figure 8 .
Figure 8. Large deviations of dynamical entropy for the quasiperiodic chain with conservative noise γ = 1 for different sizes and two value of the potential strength, ∆ = 1.5 (left columns) and ∆ = 2.5 (right columns) ; (a,b) the scaled largest eigenvalue of W + (s) as a function of s for different chain lengths N ; (c,d) the order parameter ρ + (s); (e-h) the eigenvectors φ + for N = 1597 as a function of the lattice site number (e,g) s = −0.1 (f,h) s = 0.1.