Anomalies, absence of local equilibrium and universality in 1-d particles systems

One dimensional systems are under intense investigation, both from theoretical and experimental points of view, since they have rather peculiar characteristics which are of both conceptual and technological interest. We analyze the dependence of the behaviour of one dimensional, time reversal invariant, nonequilibrium systems on the parameters defining their microscopic dynamics. In particular, we consider chains of identical oscillators interacting via hard core elastic collisions and harmonic potentials, driven by boundary Nos\'e-Hoover thermostats. Their behaviour mirrors qualitatively that of stochastically driven systems, showing that anomalous properties are typical of physics in one dimension. Chaos, by itslef, does not lead to standard behaviour, since it does not guarantee local thermodynamic equilibrium. A linear relation is found between density fluctuations and temperature profiles. This link and the temporal asymmetry of fluctuations of the main observables are robust against modifications of thermostat parameters and against perturbations of the dynamics.


Introduction
Understanding the fluctuation properties of nonequilibrium phenomena is one of the major tasks of modern statistical physics [1,2,3,4,5,6,7,8,9]. To this purpose, various generalizations of Onsager-Machlup theory [10] have been proposed, like those of Refs. [8,9]. These papers have also motivated various works on deterministic nonequilibrium particle systems [11,12,13,14,15] meant to investigate the relationship between stochastic and deterministic models of nonequilibrium physics. One basic tenet of classical statistical mechanics, indeed, maintains that the stochastic description is but a reduced (mesoscopic) representation of the deterministic (microscopic) description of a given system of interest. Nevertheless, the relation between the two representations is far from fully understood. For instance, nonequilibrium steady states require thermostats, and one expects the state of a macroscopic system to be unaffected by the details of the thermostatting mechanism. Thus, even in mathematical models of nonequilibrium steady states, one would like stochastic and deterministic thermostats to lead to practically equivalent representations of the same phenomenon. 1 However, different theoretical models may capture different aspects of a physical process, and a complete equivalence should not be expected, except in some limiting situation, as postulated by various equivalence principles [17,18,19,20,21,22,23]. The presence of local thermodynamic equilibrium fosters the equivalence. Indeed, the thermodynamic relations are so largely independent of the details of the microscopic dynamics, hence so generally valid, precisely because of local equilibrium, cf. Section 4 below, Refs. [2,24,25,26,27,28] and references therein. In case local equilibrium cannot be established, one may still wish to organize the different cases in homogeneous groups, or universality classes [29,30], to help our understanding of nonequilibrium phenomena.
It is, then, interesting to point out the conditions under which the asymmetries of fluctuations characterizing the (intrinsically irreversible) stochastic evolutions [8,31] are present in the dynamics of time reversal invariant particle systems. These asymmetries may indeed be experimentally and numerically observed, and are thought to be responsible for the irreversibility of macroscopic phenomena [31]. In Ref. [11], the fluctuations of the current of the nonequilibrium Lorentz gas, with large numbers N of noninteracting particles, were found to be time-symmetric, despite the clearly irreversible behaviour of the system. As explained in Ref. [14], this was due to the lack of interactions among the particles, hence to the lack of correlations among them [15], revealed by the large system limit. 2 On the other hand, the absence of interactions prevents the onset of local thermodynamic equilibrium, despite the convergence to a given steady state. In Refs. [12,13,14,15], it is therefore argued that systems of properly interacting particles should typically have asymmetric fluctuation paths, as predicted for stochastic systems. More precisely, Ref. [14] states that temporally symmetric fluctuation paths of nonequilibrium time reversal invariant deterministic dynamics require the untypical condition of vanishing correlations; generically, fluctuation paths should then be asymmetric. Hence, deterministic reversible and stochastic interacting particle systems are equivalent, i.e. belong to the same universality class, from this point of view.
In the present paper, we continue the investigation of these issues, by looking at variations of the nonequilibrium Fermi-Pasta-Ulam (FPU) model introduced by Lepri, Livi and Politi [33]. This model consists of N anharmonic oscillators of same mass m and positions x j , j = 1, ..., N, interacting via the nearest neighbours FPU-β potential [34]: while the oscillators with j = 1 and j = N further interact with deterministic thermostats known as Nosè-Hoover "thermal baths", and with still walls. The internal energy of the system is then given by: The equations of motion take the form where q 0 = q N+1 = 0 represent the walls and m is the mass of one particle, which we take to be unitary, so that velocitieṡ q i and momenta p i represent the same quantity. Moreover, the interparticle forces are obtained differentiating the potential, F(q) = −V ′ (q), and the Nosè-Hoover thermostats, at "temperatures" T r and T ℓ with response times θ r and θ ℓ , are implemented by the variables ξ ℓ , ξ r , which obeẏ The nonequilibrium FPU-β model is time reversal invariant, but dissipative, in the sense that the time average of the divergence of the equations of motion, −( ξ r + ξ ℓ ), is negative. Therefore, the sum of the Lyapunov exponents is negative, the phase space volumes contract and the system approaches in time a nonequilibrium steady state, characterized by a singular phase space probability distribution.
In the theoretical calculations, the motion of particles is commonly referred to the quantity a, which is the mechanical equilibrium distance between two nearest neighbours [25]. Then, the mean displacement of the position of particle j from its equilibrium position ja, is often assumed to be small, in order to define microscopically the physical observables. It will be pointed out that this is not appropriate, except for particularly small driving forces.
A quantity of interest is the "local virial" [33]. For particles in the bulk (i = 2, ..., N − 1), this is the time average of the product of the displacement of each particle, times the net force acting on it, q i (F i+1 − F i ) , where we have set F i = F(q i − q i−1 ) for simplicity. If the anharmonic part of V is replaced by hard core elastic collisions, the interaction between particles is purely harmonic, and one has F i − F i+1 = k 2 (q i−1 − 2q i + q i+1 ). Therefore, the local virial of particles in the bulk is expressed by ) .
An important feature of this model is that its long wavelength Fourier modes, representing the slow relaxational dynamics, may be considered as approximately independent of the short wavelength modes, representing the fast dynamics [25], i.e. representing some sort of noise added to the slow dynamics. 3 In the slow dynamics, one may further separate a conservative non vanishing harmonic part, which cannot contribute to transport phenomena, from a mode interaction part [25,36,37]. From this standpoint, the nonequilibrium FPU system is similar to the nonequilibrium Lorentz gas, as argued in [14], but the presence of nonvanishing correlations, even if decaying in time, provides a mechanism for asymmetric fluctuation paths, within the corresponding decorrelation time scales, Refs. [13,14]. 4 In this paper, we consider different values of k 2 and a different kind of anharmonicity: that provided by hard core elastic collisions. We also consider the limiting cases in which either the harmonic or the anharmonic interactions vanish. These limiting cases share indeed some peculiarities of noninteracting particle systems, and constitute good candidates for the study of the onset of temporally symmetric fluctuations [14,15]. The fact that all particles are identical, further means that there is no disorder in our chains, hence no reasons for chaos, mixing or ergodicity, except those that may be contributed by the boundaries. We will see that nonequlibrium boundary conditions make correlations develop even in the bulk of such particles chains, and lead to asymmetric fluctuations. The results presented here lead to the following main conclusions: a) The kinetic temperature profiles are linearly related to the deviations of the mean positions of the particles from their equilibrium values, hence to the deviations from uniformity of the density profiles (Section 2): This relation and the temporal asymmetry of fluctuations are robust against modifications of the parameters of the thermostats and of the interaction potentials, besides being robust against stochastic perturbations of the deterministic dynamics (Appendix B).
b) Comparison between our results for deterministic, time reversal invariant, dynamics and those found in the literature for chains with stochastic thermostats, shows that the behaviours concerning the two kinds of thermostats are only qualitatively similar. The two kinds of thermostats do not belong to the same universality class, except in a coarse sense; equivalence is then only expected in the equilibrium limit.
Other relevant conclusions include: c) The deviations of the mean positions of the particles from their equilibrium values are large, except in the purely harmonic case which, however, does not sustain a temperature gradient. Hence the microscopic definitions of nonequilibrium thermodynamic observables, based on small deviations, are generally impaired (Section 2). d) As in stochastically thermostatted systems, investigated by other authors, the properties of steady states depend substantially on the microscopic details of the dynamics, showing that the absence of genuine local thermodynamic equilibrium is typical of physics in 1-d (Section 2). e) As in the stochastic case, temporal asymmetries are ubiquitous in deterministic, time reversal invariant, nonequilibrium particle systems, as long as correlations among particles are present. Such asymmetries are not restricted to large deviations in the large system limit and do not require local thermodynamic equilibrium (Section 3).
f) At variance with common expectations, Nosé-Hoover thermostats alter the bulk behaviour, even when the bulk resembles a noninteracting particle system, for both small and large temperature gradients. This is a manifestation of the nonlocality generically observed in nonequilibrium steady states [8,38,39,40] (Section 4). g) Chaos by itself is not sufficient to establish standard behaviour in 1-d systems, besides not even being necessary [33,41,42,43] (Section 5).
In Section 2, we consider hard spheres on a line, connected in pairs by harmonic forces. In Section 3, the temporal symmetries of the fluctuations of this model are analyzed. In Section 4, the limiting cases of purely harmonic forces and purely elastic collisions are considered. In section 5, the role of chaos is investigated. Section 6 is devoted to concluding remarks. Appendix A quantifies the deviations of the temperature profiles from a theoretical curve. Appendix B gives results about stochatstically perturbed dynamics.

Hard core and harmonic potentials
In the present case, the harmonic part of the potential of Eq.(1) is retained, while the term proportional to β is replaced by hard core elastic collisions between particles of radius r = 1. This model shares some features of non-interacting particles systems. Indeed, isolated one-dimensional systems of elastic particles of equal mass and size exchange their momenta in such a way that the overall motion is equivalent to that of a system of non-interacting particles, in which each particle preserves its momentum [44]. Moreover, in isolated systems, the harmonic potential should not constitute an important source of correlations, since the oscillations it entails amount to independent normal modes. Nevertheless, two features distinguish our model from systems of non-interacting particles: • the elongation of the harmonic springs discriminates the case in which particles bounce back at collisions, from the case in which they pass through one another; • the presence of the thermostats at the boundaries of the chain spoils the normal modes decomposition of the purely harmonic motion.
The question, which we have investigated numerically, is how crucially these facts affect the behaviour of the system. To verify that the steady state is indeed achieved, we have computed the local "temperatures", p 2 i , and the virial profiles, at different times along a simulation of length t max . We have found that the temperature profile rapidly approaches its asymptotic form (a quarter of our t max suffices for good accuracy), while the virial expression takes much longer to settle down, cf. Figure 1 for the case with N = 400, k 2 = 1, θ ℓ = θ r = 1.
Similar profiles have been found in numerous one-dimensional nonequilibrium chains of oscillators, with anharmonic interaction potentials, like the FPU-β chains, [25,45,46]. In particular, assuming that dynamical chaos practically amounts to some degree of stochasticity in the evolution, our model should behave similarly to that of Ref. [27], which has stochastic boundary thermostats, harmonic interactions and stochastic exchanges of momenta between nearest neighbours, meant to mimic (momentum and energy preserving) hard core interactions. In a suitable continuum limit, the corresponding energy profile takes the form [27]: where ζ is the Riemann ζ-function, and the space variable is rescaled so that the left end of the chain corresponds to x = −1 while the right end corresponds to x = 1. However, as shown by the left panel of Figure 2, and deduced from a comparison of our Figure 1 with Figure 2 of [27], the temperature profiles obtained by Lepri, Mejía-Monasterio and Politi are not completely equivalent to ours. In the first place, the slope of profile (6) has square root divergences at the boundaries, which are less steep than ours. Secondly, the profiles of Ref. [27] are odd with respect to the centre of the chain, even at finite N, while our profiles are not. One further difference is the simple k 2 dependence of the finite−N temperature profiles of [27], as opposed to the strongly irregular dependence of our profiles, cf. central panel of Figure 2. To investigate this question, we have analyzed the case with k 2 = 1, T ℓ = 320, T r = 20, θ ℓ = θ r = 1, for N = 100, 200, 250, 300, 400. Differently from [27], we found no evidence of the N −1/3 rate of convergence towards the asymptotic profile T (x) and, in most of the chain, the difference between our profiles and T (x) was observed to grow with N. Therefore, as shown in Appendix A, either a scaling regime towards T (x) sets in at N larger than 400, or our asymptotic profiles will not correspond to those of [27]. On the other hand, Appendix A shows that our profiles converge rapidly to a given asymptotic profile, which is not odd. This may be related to the fact that, differently from the stochastic case, the energy exchange between Nosé-Hoover boundary thermostats and 1-d (or almost 1-d) systems of particles is more efficient at the hot side than at the cold side, cf. Ref. [47] and Section 4 below. 5 5 In Ref. [47] it was observed that small temperature gradients allow the cold thermostat to behave as "thermodynamically" as the hot thermostat. The boundary thermostats of 1-d chains act only on momentum degrees of freedom, making it difficult, especially at large temperature gradients, for local equipartition to be established, except in the presence of mass diffusion [48]. However, the difference between T and the profiles obtained by Nosé-Hoover boundary thermostats persists at small gradients, as further evidenced in Section 4. i profile reaches rather rapidly its asymptotic form (the curves for the four simulation lengths are indistinguishable), while the virial expression takes longer to converge. Because of hard core interactions, kinetic temperature and virial do not coincide [49].
Despite these facts, the temperature profiles concerning nonequilibrium chains of oscillators typically enjoy the same qualitative behaviour, consisting of steep curves at the boundaries, due to contact resistance, interpolated by almost straight lines in the bulk. The nature of this behaviour is far from obvious, cf. the discussion in Ref. [27,46]. In any event, the strong and irregular dependence of the steady states of both deterministic and stochastic 1-d systems, on microscopic details of the dynamics, as well as the violations of Fourier Law, reveal the absence of genuine local thermodynamic equilibrium and of diffusion, see e.g. Refs. [16,25,26,28,45].
For certain observables, this dependence on the microscopic mechanisms, hence the lack of local thermodynamic equilibrium, persists in the large N limit, [16,27,28]. Similarly, completely different situations are determined by the boundary conditions, such as the divergence of the conductivity of disordered harmonic chains with free boundaries, in the thermodynamic limit, as opposed to its vanishing trend in chains with fixed boundaries [28,30]. Moreover, the scaling behavior of the conductivity with system size crucially depends on the spectral properties of the heat baths [50]. Rieder, Lebowitz and Lieb interpreted this kind of results as "unphysical" [16]. Indeed, were 1-d chains of oscillators to represent macroscopic objects enjoying thermodynamic properties, the thermostat dependence of the heat fluxes of [16], or of our temperature profiles, would be unrealistic. However, recent theoretical and technological developments allow different interpretations. The study of nearly 1-d systems, especially with non-macroscopic numbers of elementary constituents, indicates that standard thermodynamic relations typically fail to describe their behaviour, although a complete understanding of their properties is currently missing, see e.g. [30,46,51,52] and references therein. Thus, a possible absence of local thermodynamic equilibrium does not need to be unphysical. Indeed, local thermodynamic equilibrium requires a sufficiently fast decay of correlations of all relevant observables, which would lead to diffusive transport. But this is seldom afforded by (quasi-) 1-d systems [26,53], 6 which include modern technological artifacts as well as structures found in nature. It is also known that 1-d systems are affected by peculiarities such as the cumulative O(N) fluctuations, about the particles equilibrium positions [54]. These frustrate the microscopic definitions of observables based on the assumption of small fluctuations and, more importantly, are at odds with the properties of solids. Furthermore, we have found that temperature gradients induce large displacements of the average equilibrium positions, making the distribution of matter inhomogeneus. For three different values of k 2 , the left panel of Figure 3 portrays these average deviations, q i = x i − ai , for i = 1, ..., N = 100. The right panel of the figure portrays the quantity x i+1 − x i , which is related to the inverse of the density of particles. One observes that nonequilibrium boundary conditions lead to a shift of all particles towards the cold side of the chain, hence to a density gradient.
Denoting by · k 2 the averages computed with elastic constant equal to k 2 , closer examination shows that: • all q i k 2 are positive and, in the bulk, are approximated by a parabola with maximum at i > N/2 for k 2 0, while the maximum of q i 0 occurs at i = N/2; • q N k 2 = q N k2 for all pairs k 2 ,k 2 , while q 1 k 2 = q 1 k2 only if both k 2 ,k 2 0; • the dependence of q i k 2 on k 2 is not monotonic, like the temperature profile is not.
Therefore, the systems with pure hard sphere interactions, whose bulk corresponds to noninteracting particles, behave differently from the systems with both hard core and harmonic interactions. Indeed, while the temperature profiles at k 2 = 0 do not substantially differ from those at k 2 0, the net energy transfer with k 2 = 0 is singular, since it is equivalent to that of a single particle, bouncing back and forth between the hot and cold walls, independently of the value of N. Recalling that a = 5, in our calculations, these plots show that the steady state (average) particle distribution is not uniform and that the deviations of particles positions from their equilibrium values are large and correlated to the kinetic temperature profiles. Indeed, x i+1 − x i which, apart from an unessential additive constant, is the discretized derivative of q i , qualitatively approximates the typical temperature field (compare e.g. the right panel of Figure 3 with Figures 1 and  2), hence a linear relation between the two quantities may be surmised: The good agreement between T i and the kinetic temperature profile p 2 i is demonstrated by Fig.4.  This relation between temperature profiles and rescaled average displacements of neighbouring particles is robust against modifications of all the parameters of the dynamics, including thermostat parameters and interaction potentials, besides being robust against (small and large) stochastic perturbations of the dynamics, cf. Appendix B. Apart from the temporal asymmetries of the fluctuations of the main observables, considered in the next section, this is the only result which does not show a delicate dependence on the details of the microscopic dynamics.

Temporal symmetry of fluctuation paths
In the literature, various notions of fluctuation path have been investigated. Given the substantial equivalence of the results based on such different notions, we adopt the first definition of fluctuation-relaxation path of Ref. [14], denoted by FR1. 7 Denote by M the phase space of our system, by S t : M → M the time evolution operator, and by X : M → I R an observable of interest. Consider an initial phase Γ ∈ M , in the support of the steady state phase space probability distribution, and denote by X t the quantity X(S t Γ), under the assumption that all, but a set of vanishing probability, such initial conditions enjoy the same statistics. Choose a fluctuation value T (X). The FR1 fluctuation path is defined as follows: Definition. Assume that Xt = T (X). Then, for any t 0 , τ > 0, the FR1 fluctuation path of duration 2t 0 based att is the curve As in Ref. [14], the symmetry properties of the fluctuation paths of the observable X are here assessed setting a threshold T (X) = I E(X) + 3σ(X), three standard deviations above the mean I E(X), while the width of the observation time interval is taken to be 2t 0 = 2. Then, all such time intervals are translated by −t, so that the centers of all fluctuations coincide at t = 0. The first fluctuating observable we consider is the local dimensionless density of particles in the center of the chain, defined in a box of size L, ideally with 1 ≪ L ≪ N: where ρ 0 = L/(L − 1)a is the equilibrium density.
where τ labels a bin, and C τ (t) = 1 if t belongs to the τ-th bin, while C τ (t) = 0 if it does not. The quantity m(τ) is the number of integration time steps which make a bin. 8 Figure 5 represents the normalized average paths, R τ = (ρ τ − min ρ τ )/(max ρ τ − min ρ τ ), which start at zero in the first bin and reach 1, their maximum value, in the central bin. The result is that the average path is temporally asymmetric, with growth steeper than relaxation, and that, differently from the observables considered in Section 2, it is not appreciably affected by k 2 . 7 See [13,14] for a discussion of the ambiguities in the definitions of fluctuation paths in deterministic dynamics. Despite such ambiguities, the different notions of fluctuation path have led to analogous conclusions and can be used interchangeably to asses the symmetry properties of nonequilibrium fluctuations. Furthermore, the different notions coincide in the large N limit, cf. [12,13,14,15]. 8 In our simulations, the integration time step equals 10 −3 for small N, and 10 −4 for large N, so that the numerical errors are of the same order.  This is consistent with Refs. [14,15], in which asymmetric paths are expected to be typical of nonequilibrium steady states of interacting particle systems. In Figure 5, the asymmetry seems to grow with N, as in [14], but a comparison of cases with different N is frustrated by the unclear dependence on N of the many parameters entering the definition of the asymmetry, like T (X) and t 0 . For instance, Figure 5 seems to imply that average growth (τ < 0) and relaxation (τ > 0) tend to become linear as N increases, but this may be due to the fact that the observation time t 0 , which ought to be connected with the correlations decay rate, increases with N. The right panel of Figure 5 would then only illustrate a smaller fraction of the average path than that of the left panel. As a matter of fact, Figure 6 shows that the instantaneous asymmetry, δ τ = ρ τ − ρ −τ , almost saturates with growing τ, but it saturates at larger values of τ for larger N. On the other hand, it is not clear how t 0 should be modified with N; it is only obvious that ρ τ − ρ −τ must decrease back to 0, for τ sufficiently long that correlations have decayed [15]. Unfortunately, that time is too long to allow us to collect a statistically relevant sample of fluctuation paths.  Let us introduce the asymmetry ν, as a normalized cumulative difference between fluctuation and relaxation trajectories, The probability that this quantity be positive, P(ν) > 0, is reported in Table 1 for N = 100, where one observes that, from this point of view, the k 2 = 0 case does not differ substantially from the k 2 > 0 cases. The above figures and table suggest that nonequilibrium steady states of identical hard spheres on a line do not behave like independent particles, although in equilibrium, when the thermostats are removed, they do. Away from equilibrium, particles appear to develop space correlations, which reach the bulk, eventually connecting them to the boundaries, even in the absence of springs. Therefore, the observed asymmetries reveal that a form of nonlocality (long range correlations) [38,39,40,55] is common in nonequilibrium states of both deterministic and stochastic models.
Another observable of interest is the heat flux. To define the heat flux in the presence of hard core collisions, we adopt the method of planes, developed by Todd, Daivis and Evans [56], according to which the fluctuating heat flow J q , through a plane transversal to the medium, can be decomposed as the sum of two terms: The term J K q represents the kinetic part of the heat flux, which is expressed by: where i denotes a particle crossing the plane at (discrete) times t i,m i , c xi is the component along direction x of the particle's velocity with respect to the streaming velocity, A is the area of the plane (if the system is 3 dimensional), and is the internal energy concerning particle i. Moreover, v i is the velocity of the particle in the laboratory frame, u(x i ) is the streaming velocity at position x i , and φ i j is the interaction potential between particles i and j. The term J U q represents the potential part of the heat flux, which is defined by: where F i j is the force that particle j exerts on particle i. Figure 7 shows the average fluctuations of J U q for the cases with N = 100, T ℓ = 320, T r = 20, and for different values of k 2 . Analogously to the density average path, the average J U q fluctuation rises faster than it relaxes, and is not appreciably affected by the value of k 2 . Results similar to those reported have been obtained for all different parameters choices, even with k 2 = 0, when the bulk behaviour should have closely resembled that of non-interacting particles. This is consistent with the discussion of the single file behaviour of [53], according to which the order of particles introduces persistent correlations, even in the presence of positive Lyapunov exponents. The results appear to be robust and independent of the time reversal parity of the observable at hand, since local density and heat flux have opposite parity.
The results of this section definitely prove that asymmetric fluctuation paths are typical of deterministic, time reversal invariant nonequilibrium models.

Harmonic chain in a nonequilibrium steady state
In this section we consider the model studied by Rieder, Lebowitz and Lieb [16], with the stochastic boundary thermostats replaced by Nosé-Hoover thermostats, which make the dynamics time reversal invariant. The oscillators are thus coupled by the potential energy V of Eq.(1), with β = 0. In equilibrium, this system is equivalent to a set of independent particles, the normal modes.
If applicable, the local version of the virial theorem now implies [33]: Clearly, in the absence of thermostats, the local virial theorem does apply to each independent mode of oscillation, but in general, this has to be checked. As in the case of stochastic baths considered by [16], we find that Eq. (15) holds and that the kinetic temperature profile is characterized by jumps at the boundaries, and by approximately flat profiles in the bulk. As in [16], harmonic chains do not sustain a kinetic temperature gradient but, differently from the case of [16], the approximate equipartition of energy in the bulk does not correspond to the mean of the left and right temperatures: it stays closer to the higher temperature, cf. Fig.8. The shorter θ, the stronger the coupling between particle and thermostat. Results at short θ = θ ℓ = θ r are reported in the left half of the left panel, while results at large θ are reported in the right half. After some oscillation, decreasing θ makes the profile settle around the value 240, which is much closer to T ℓ = 320 than to T r = 20. Note the apparent discontinuity of the profile between θ = 0.5 and θ = 0.3, which leads to the flat bulk profile at (T ℓ + T r )/2. The right panel portrays the profiles with θ ℓ = 1 and decreasing θ r : the system equilibrates with the hot thermostat independently of θ r . In both panels N = 100. Some dependence of the bulk behaviour on the characteristics of the thermostats is not surprising: the choice of the parameters determines the efficiency of the energy exchange between thermostats and thermostatted system. This, in turn, may affect the behaviour of a non-macroscopic system since, by definition, this does not enjoy any local thermodynamic equilibrium. In our case, the relevant parameters consist of the target temperatures T ℓ and T r , of the thermostat characteristic times θ ℓ and θ r , and of the elastic constant k 2 . We have studied the behaviour of the profile in the θ = θ ℓ = θ r → 0 limit, the limit of strong coupling between thermostats and thermostatted particles, as well as in weak coupling cases. Figure 8 (left panel), shows that the kinetic temperature profile settles closer to T ℓ than to T r , at all values of θ, except at θ = 0.3. The situation does not change if θ r is shorter than θ ℓ , so that the cold thermostat is more strongly coupled to the system than the hot thermostat, as illustrated by the right panel of Fig.8. For θ = 0.3, the deterministic time reversal invariant Nosé-Hoover thermostats behave most similarly to the irreversible stochastic ones of [16], as far as the kinetic temperature profile is concerned. Note, in particular, the flatness of the θ = 0.3 profile, and the temperature overshoot on the cold side, predicted by the theory of [16]. The equivalence with the results for stochastic thermostats of [16] is not complete, however, because the temperature does not overshoot in the hot side.  i.e. weak coupling, it settles within [T r , T ℓ ], while for low θ, i.e. strong coupling, it reaches even lower than T r . The profile depends quite irregularly on θ. For θ = 0.7 and θ = 1 the bulk temperature settles close to (T ℓ + T r )/2. This is consistent with the observations of Ref. [47], on a system of rotating disks and pointlike particles, thermostatted by Nosé-Hoover mechanisms acting on the boundary disks. There, the hotter side reaches local equilibrium more easily than the colder one, even if not necessarily at the target temperature of the thermostat. The thermostat closer to an equilibrium state (as it should be) seems, then, more efficient in driving the bulk of the chain than the other thermostat. In Ref. [47], it was further observed that small gradients lead to standard behaviour of the cold thermostat. Differently, in our case, the temperature profiles remain quite irregular even under relatively small gradients, and only for a limited set of θ's do they lie close to (T ℓ + T r )/2. The profiles worsen in the small θ limit, i.e. for strong couplings, when they even brim over the interval [T r , T ℓ ], cf. Fig.9. The strong dependence on the microscopic parameters is confirmed by the departure at large θ of the bulk profiles from the mean of the boundary values, cf. the case θ = 2 in the left panel of Fig.9.
To further assess the relation between our chains and those of [16] we have produced the histograms of p i , and have verified that they are not Gaussian, although to different degrees for the different i's. Indeed, the distributions of p 1 and p N/2 differ much less from a normal distribution than that of p N (N is the particle interacting with the cold bath at temperature T ℓ ), as indicated by normal distribution plots not reported here. Also, local Gaussian distributions are better and better approximated at both ends of the chain, as in [47], if the temperature gradient is reduced. Increasing k 2 , without changing the thermostats properties, is another route towards approximately Gaussians local distributions. This confirms that closer to equilibrium it is easier to obtain the equivalence between deterministic and stochastic thermostats, as generally expected [26], although the temperature profiles demonstrate that complete equivalence is not to be expected.
The peculiar behaviour of the purely harmonic nonequilibrium chains is evidenced also by quantities such as x i+1 − x i and q i , whose profiles are much more irregular than those of Figure 3, cf. Figure 10. Here, differently from Figure 3

Lyapunov exponents: the harmonic and the hard core cases
In this section we consider the Lyapunov exponents of the two limiting cases of purely harmonic and purely hard core interactions. These exponents have been computed in tangent space, with the usual algorithm devised by Benettin, Galgani, Giorgilli and Strelcyn [57]. The phase space contraction rate χ, i.e. the dynamical dissipation, can be separately computed as the negative of the average of the divergence of the vector field, neglecting the instantaneous elastic collisions, which do not contribute to the variations of the phase space volumes [58]. In our case, χ = ( ξ ℓ + ξ r ) and χ ≥ 0, where equality characterizes equilibrium states, while χ > 0 characterizes nonequilibrium states. The motion of purely harmonic chains, in the absence of thermostats, is hamiltonian and fully integrable, hence all Lyapunov exponents vanish. Table 2 shows the values of the largest Lyapunov exponent, λ 1 , for several different kinetic temperatures at the left end of a chain of 100 particles, while the right end is subjected to a thermostat with T r = 20. It is interesting to note that even at equilibrium, i.e. for T ℓ = T r , the presence of the thermostats makes chaotic the dynamics. Indeed, λ 1 is higher at lower temperature differences, since higher dissipations imply smaller attractors and more orderly states. It also appears that λ 1 saturates as a function of the temperature difference, while the dissipation keeps increasing significantly. As a function of the relaxation times of the thermostats, the largest Lyapunov exponent decreases monotonically. For instance, the largest exponent for N = 100 and T ℓ = 320 decreases from the value λ 1 = 0.00779 reported in Table 2 for θ = 1 to λ 1 = 0.000444 for θ = 10. This is a consequence of the fact that the motion of single oscillators is integrable and that large θ practically decouples the thermostats from the system of interest. Table 3 reports the computed values of λ 1 and χ, at fixed temperature difference and varying chain length. The Lyapunov exponent decreases with N, hence with the temperature gradient, while the dissipation χ is practically constant.  Table 3. Largest Lyapunov exponent and dissipation for the purely harmonic case, with k 2 = 1, and fixed temperature difference: T r = 20, T ℓ = 320. The global gradient decreases with N, like λ 1 , while χ is practically constant.
The picture is completed by the behaviour with N at fixed temperature gradient, cf. Table 4. 10 The data indicate that two competing effects contribute to the values of the largest Lyapunov exponent: a) decreasing the temperature gradient, the system approaches the equilibrium state, which enjoys the largest value of λ 1 afforded by a chain of given length; b) increasing N reduces the impact of the thermostats on the bulk of the chain, which then better approximates an isolated, fully integrable system, with vanishing Lyapunov exponents.
As λ 1 does not grow indefinitely with (T ℓ − T r ), increasing N leads to a decrease of λ 1 . The dissipation χ , on the other hand, is more directly related to the value of (T ℓ − T r ), and grows with it.  Table 4. Largest Lyapunov exponent and dissipation for the purely harmonic case, with k 2 = 1, T r = 20 and fixed global temperature gradient (T ℓ − T r )/N = 3.
An interesting feature of the purely harmonic case, with Nosé-Hoover thermostats at different temperatures, is that non-chaotic dynamics can at least be realized in short chains. For instance, in the case with N = 3, T ℓ = 80, k 2 = 1 and T r = 20, all Lyapunov exponents are negative, except for one vanishing exponent, and the motion is periodic. As a matter of fact, even some cases with larger N might share these features, but their largest exponents appear to be only marginally positive, hence hard to distinguish from vanishing values. Differently, the case with N = 4 has clearly positive exponents.
Consider now the case in which particles interact via hard core elastic collisions only, which, in the absence of thermostats, would be the other limiting, fully integrable case. The largest exponents are commonly found to be positive and behave as in the purely harmonic case, cf. Table 5.  Table 5. Largest Lyapunov exponent and dynamical dissipation for hard particles with k 2 = 0, N = 100 and T r = 20. The behaviours of both λ 1 and χ do not differ substantially from those of purely harmonic cases.
The comparison between hard core and purely harmonic cases illustrates one interesting fact, consistent with the observations of Ref. [41]. Both kinds of systems have positive Lyapunov exponents, but the purely harmonic case shows a quite peculiar bulk equipartition of energy, while the chains of purely hard core particles verify more standard nonequilibrium conditions. The randomness entailed by hard collisions, which occur at practically random times and positions, as assumed in [27], breaks correlations among particles more efficiently than generic chaotic mechanisms do and favours standard behaviour, as argued also in [41]. However, the other peculiarities noted in the previous sections still do not allow us to speak of genuine local thermodynamic equilibrium for chains of hard particles. That chaos does not suffice for standard behaviour, in FPU systems, has been observed in the past [33,43].

Concluding remarks
We have analyzed deterministic chains of oscillators with harmonic forces and elastic collisions, to asses the equivalence of deterministic time reversal invariant and stochastic models of thermostats. Close to equilibrium, it is easier to obtain equivalent behaviours than far from equilibrium, although the peculiarities of 1-d dynamics prevent complete equivalence, even in the large N limit. Only from a qualitative standpoint, various properties are common to the different models, like the form of the temperature profiles, the temporal asymmetries of flucutations and the nonlocality of nonequilibrium steady states [3,8,9,38,40]. Similarly, the strong dependence of the bulk behaviour on microscopic details of the dynamics, such as the boundary conditions, which is enhanced by the growth of the dissipation, reveals qualitative similarities between the different kinds of thermostatted evolutions. This qualitative similarity is robust against changes of the interaction potentials. For instance, the overall behaviour of systems made of purely harmonic identical oscillators differs from that implied by hard core and other non-linear forces, but this happens for both stochastic and deterministic thermostats.
Therefore, from relatively coarse, qualitative, viewpoints may we state that the different models of thermostat are equivalent. This weak equivalence suffices to conclude that anomalous behaviours are in fact typical of physics in (quasi) 1-dimension, but does not suffice to predict the kind of anomalies that one should expect in different situations. As a matter of fact, Appendix B shows that numerous parameters affect the temperature profiles, including the elastic constant k 2 and any perturbation of the interaction forces. In particular, the presence of different degrees of stochasticity produces different behaviours, and only in the equilibrium limit may the equivalence of the different models be established.
Independently of the nature of the thermostats, "realistic" particle interactions, or mechanisms preventing the conservation of the most obvious dynamical quantities do not lead to the onset of local thermodynamic equilibrium [26,41,46]. One reason is that certain correlations (e.g. particles order) persist in time, violating the molecular chaos hypothesis of kinetic theory, in 1-d or quasi-1-d systems. Thus, some kind of anomalous behaviour is to be generically expected, as a consequence of the low dimensionality of the dynamics, rather than of the peculiarities of thermostats. Indeed, in 3 dimensions, the same thermostatting mechanisms do not lead to anomalous behaviour. Even in cases of apparently normal transport in 1-d, this is not as robust, with respect to variations of the microscopic parameters, as it is in macroscopic thermodynamic phenomena [26,41], indicating that genuine local thermodynamic equilibrium does not hold in those cases. Experimental evidence, although presently scanty, confirms this picture, cf. e.g. [46,51,52,53,60,61,62]. Further study is desired, to understand the relation between the peculiarities of the 1-d models so far considered in the literature and those of real systems, not in local thermodynamic equilibrium.
We found the linear relation (5) between the particles mean displacements from their equilibrium positions and the temperature profiles. Apart from the temporal asymmetries of the fluctuations of the main observables, considered in Section 3, relation (5) is the only result which does not show a delicate dependence on the details of the microscopic dynamics. It is robust against all modifications of the dynamics, which we have investigated, and must therefore play an important role in the identification of universality classes and in the equivalence of thermostats in 1-d systems Our analysis of systems with purely harmonic and purely hard core interactions indicates that chaos, or a generic source of randomness, per se, does not suffice to establish standard nonequilibrium steady states. On the other hand, it has been observed by various authors that dynamical chaos is not even necessary for that, see e.g. [41,42].

Appendix A
To assess the relevance of our results for the large N limit, and to compare with the results of Refs. [27,28] we have produced temperature profiles for hard core collisions and k 2 = 1, with T ℓ = 320, T r = 20, θ ℓ = θ r = 1, and growing N. Figure 11 shows these profiles, plotted as functions of the rescaled variable x = i/N. One observes a scaling similar to that of Fig.1 in Ref. [43], where the asymptotic profile is well approximated, on the scale of the figure, at moderately large values of N.  The difference ∆ N (y) = T (y) − T N (y) between the analytic profile T of Ref. [27] and our profiles, T N , with normalized variable y ∈ [−1, 1], for finite N, does not scale with N as the difference between asymptotic and finite N profiles of [27].
If the scaling was confirmed, we would have had where g(y) depends on the details of the model, but not on N. As a matter of fact, the absolute value of ∆ N (y) initially grows with N, rather than decreasing, e.g. in the center of chain. This does not prevent the scaling to set in at still larger values of N, although this looks unlikely, given the observed convergence of our profiles to a different asymptotic shape.

Appendix B
We consider a stochastic perturbation of the model of Section 2. If the dynamics of Section 2 leads two particles to collide at a given instant of time t, they will actually collide with probability 1 − p, p ∈ [0, 1], at time t. Taking p = 0 yields the dynamics of Section 2, while p = 1 yields the purely harmonic dynamics of Section 4, because no collision takes place. For a system of N = 100 particles, we find that the temperature profiles strongly depend on the values of the pair of parameters p and k 2 , confirming that nonequilibrium chains of oscillators may hardly be part of a unique universality class, except in the equilibrium limit, i.e. for temperatures at the boundaries very close to each other. The states of these chains are indeed too far from local equilibrium, and sensitive to many details of the microscopic dynamics, in general.
In particular, Fig. 13 shows the temperature profiles for different values of p and k 2 = 0.1, for two temperature differences (T ℓ −T r ). The dependence on p is quite strong although, in both cases, the profiles settle around the unperturbed ones, when p decreases. The cases with lower temperature difference more rapidly collapse on a unique shape, different from the theoretical T of [27]. and T r = 20. The right panel reports the case with T ℓ = 120 and T r = 20. The continuous line represents the temperature profile T (x) of eq. (6). The lower temperature jump leads to more rapid convergence to the unperturbed profiles, as p → 0. Figure 14, shows that increasing the rigidity of the chains, reduces the effect of the parameter p. Comparing Figs. 13 and 14 leads to the conclusion that, even for a fixed kind of thermostats, the equivalence of the different microscopic dynamics requires various parameters to be adjusted, in general. In our case, k 2 plays an important role, which it did not in Refs. [27,28]. Taking of k 2 = 2 does not yield profiles which approximate the theoretical one of [27] more closely than k 2 = 1.  The only result which is robust against all the modifications of the dynamics, which we have investigated, is the validity of the linear relation (5), which must then play an important role in the identification of universality classes and in the equivalence of thermostats in 1-d systems. Figures 15 and 16   For a given choice of the parameters of the deterministic model, the values of the parameters β 1 and β 2 of Eq.(5) do not depend on p, as long as p < 0.9. They change discontinuosly for p ≥ 0.9, i.e. close to the purely harmonic chains, which are quite peculiar models, compared to the others.  (5) is confirmed in all instances as discussed in Section 2.