# Finite thermal conductivity in 1d lattices

## A The rotator model

The simplest example of a classical-spin 1d model with nearest neighbour interactions lies in the class (1) with $`V(x)=1\mathrm{cos}(x)`$ and $`U=0`$. This model can be read also as a chain of $`N`$ coupled pendula, where the $`p_i`$’s and the $`q_i`$’s represent action-angle variables, respectively. It has been extensively studied  as an example of chaotic dynamical system  that becomes integrable both in the small and high energy limits, when it reduces to a harmonic chain and free-rotators, respectively. In the two integrable limits, the relaxation to equilibrium slows down very rapidly for most of the observables of thermodynamic interest (e.g., the specific heat) . As as consequence, the equivalence between ensemble and time averages is established over accessible time scales only inside a limited interval of the energy density $`\epsilon `$. Here we shall discuss heat conduction for values of the energy density corresponding to strongly chaotic behaviour.

## B Numerical analysis of the thermal conductivity

The most natural and direct way to determine $`\kappa `$ consists in simulating a real experiment, by coupling the left and right extrema of the chain with two thermal baths at temperatures $`T_L>T_R`$, respectively. In our simulations we have used Nosé-Hoover models of thermostats , both because they can be easily implemented (integrating the resulting equations with a standard algorithm) and because of the smaller finite-size effects (due to the unavoidable contact resistance).

With this setting, a non-equilibrium stationary state sets in characterized by a non-vanishing heat flux $`J`$ :

$$J=\frac{1}{N}\underset{i}{}j_i=\frac{1}{N}\underset{i}{}\frac{p_i}{2}\left(f_{i+1}+f_i\right)$$
(2)   
where $`f_i=V(q_{i+1}q_i)/q_i=\mathrm{sin}(q_{i+1}q_i)`$ is the interaction force and $`j_i`$ is the local flux at site $`i`$. The total heat flux $`J`$ has to be averaged over a sufficiently long time span to get rid of fluctuations and to ensure the convergence to the stationary regime. This can be tested by monitoring the average heat flux and looking at the scale of its fluctuations. As a result, we have verified that $`210^6`$ time units are sufficient to guarantee a few percents of fluctuations in the worst cases.

The thermal conductivity is determined by assuming the Fourier law, i.e. from the relation $`J=\kappa T`$, where $`T`$ denotes the imposed thermal gradient. The simulations have been performed for $`T_L=0.55`$, $`T_R=0.35`$, and chain lengths ranging from $`N=32`$ to 1024 with fixed boundary conditions. The equations of motion have been integrated with a 4th-order Runge-Kutta algorithm and a time step $`\mathrm{\Delta }t=0.01`$. The results, reported in Fig. 1 clearly reveal a convergence to a $`\kappa `$-value approximately equal to 7 (see the circles). The dotted line represents the best data fit with the function $`a+b/N`$: the agreement is very good, showing that finite-size corrections to $`\kappa `$ are of the order $`O(1/N)`$, as it should be expected because of the thermal contacts. However, more important than the numerical value of the conductivity is its finiteness in spite of the momentum conservation.

In order to test independently the correctness of our results, we have performed direct microcanonical simulations, which allow determining the thermal conductivity through the Green-Kubo formula :

$$\kappa =\frac{1}{T^2}_0^{\mathrm{}}C_J(t)𝑑t$$
(3)   
where $`C_J(t)=NJ(t)J(0)`$ is the flux autocorrelation function at equilibrium and $`T`$ is the temperature. A correct application of the above formula requires fixing the energy density $`\epsilon `$ in such a way that the kinetic temperature (defined as $`T=p^2`$, in agreement with the virial theorem) is close to the average value of the temperature in the previous simulations. The choice $`\epsilon =0.5`$ turns out to be reasonable, as it corresponds to $`T0.46`$. In the absence of thermal baths, the equations of motion are symplectic, so that we have now preferred to use a 6th order McLachlan-Atela integration scheme  (with periodic boundary conditions).

The correlation function has been computed by exploiting the Wiener-Khinchin theorem, i.e. by anti-transforming the Fourier power spectrum. The result of the time integration is almost independent of $`N`$ for $`N>128`$. The gray region in Fig. 1 corresponds to the expected value of $`\kappa `$ taking into account statistical fluctuations. There is not only a clear confirmation of a finite conductivity, but the numerical value obtained with this technique is in close agreement with the direct estimates.

## C Dynamics in the mode space

In order to clarify the difference between the dynamics of the present model and that of FPU-type systems, we have investigated the evolution of the low-frequency Fourier modes. In Fig. 2 we have reported the power spectra of some long-wavelength modes. For the sake of comparison, the same quantities are reported for a diatomic FPU-chain, that is characterized by an anomalous transport. At variance with the FPU model, in the rotators there is no sharp peak (which is a signal of an effective propagation of correlations ). Quite differently, the low-frequency part of the spectrum is described very well by a Lorentzian with halfwidth $`\gamma =Dk^2`$ ($`D4.3`$). This represents an independent proof that energy diffuses, as one expects whenever the Fourier law is established.

## D Temperature dependence of the thermal conductivity

The most natural question arising from these results concerns the reason for the striking difference with other symmetric models such as the FPU-$`\beta `$ system. As long as each $`(q_{i+1}q_i)`$ remains confined to the same valley of the potential, there cannot be any qualitative difference with the models previously studied in the literature. Jumps thorugh the barrier, however, appear to act as localized random kicks that contribute to scattering the low-frequency modes and thus to a finite conductivity. If this intuition is correct, one should find analogies between the temperature dependence of the conductivity and the average escape time from the potential well. To this aim, we have computed $`\kappa `$ for different temperature values by performing microcanonical simulations with various energy densities. From the data reported in Fig. 3, one can notice a divergence for $`T0`$ of the type $`\kappa \mathrm{exp}(\alpha /T)`$ with $`\alpha 1.2`$. An even more convincing evidence of this behavior is provided by the temperature dependence of the average escape time (see the triangles in Fig. 3) with an exponent $`\alpha 2`$. The latter behaviour can be explained by assuming that the jumps are the results of activation processes. Accordingly, the probability of their occurrence is proportional to $`\mathrm{exp}(\mathrm{\Delta }V/T)`$, where $`\mathrm{\Delta }V`$ is the barrier height and the Boltzmann constant is fixed equal to 1 (as implicitely done throughout the Letter). Since $`\mathrm{\Delta }V=2`$, the whole interpretation is consistent. Moreover, in the absence of jumps, the dependence of the conductivity on the length should be the same as in FPU-systems, i.e. $`\kappa N^{2/5}`$. Therefore, a low-frequency mode travelling along the chain should experience a conductivity of the order of $`\overline{N}^{2/5}`$, where $`\overline{N}`$ is the average separation between jumps. Under the assumption of a uniform distribution of phase jumps, their spatial separation is of the same order of their time separation, so that we can expect that $`\kappa \mathrm{exp}[2\mathrm{\Delta }V/(5T)]`$. On the one hand, this heuristic argument explains why and how such jumps contribute to a normal transport. On the other hand, the numerical disagreement between the observed and the expected value of the exponent $`\alpha `$ (1.2 vs. 0.8) indicates that our analysis needs refinements. In fact, we should, e.g., notice that, in the low-energy limit, nonlinearities become negligible, implying that deviations from the asymptotic law $`\kappa N^{2/5}`$ should become relevant.

## E Further checks

In order to test the conjecture that jumps are responsible for a normal heat transport, we have decided to investigate some other models. First, we have considered a double-well potential $`V(x)=x^2/2+x^4/4`$ (the same as in FPU with a different sign for the harmonic term). The results of the direct simulations are reported in Fig.1 (see triangles) for a temperature corresponding again to a quarter of the barrier height. The finiteness of the conductivity ans its numerical value is confirmed by the computation of $`\kappa `$ through the Green-Kubo formula (see the light-grey shaded region).

Finally, we have considered an asymmetric version of the rotator model, namely $`V(x)=A\mathrm{cos}(x)+0.4\mathrm{sin}(2x)`$, where $`A`$ is fixed in such a way that the minimum of the potential energy is zero, and the temperature again corresponds to one quarter of the barrier height. In this case too, the conductivity is finite, confirming our empirical idea that the jumps are responsible for breaking the coherence of the energy flux and, in turn, for the finite conductivity. However, from the point of view of Ref. , this result is quite unexpected, since, in view of the model asymmetry, one expects that the pressure $`\varphi =f_i/N`$ is non-zero. Nonetheless, microcanonical simulations show that although the distribution of forces is definitely asymmetric, their average value is numerically 0. This can be understood by noticing that in view of the boundedness of the potential, the system cannot sustain any compression. Therefore, no contradiction exists with Ref. .

In conclusion, in this Letter we have reported about the first evidence of normal heat transport in 1d systems with momentum conservation. Such a behaviour appears to be connected with jumps between neighbouring potential valleys. From the dynamical point of view, it is natural to ask what are the peculiar properties of such jumps that make them so different from other types of nonlinear fluctuations that may occur in single-well type of potentials. The only clearly distinctive feature that we have found is the “hyperbolic” type of behaviour in the vicinity of a maximum of the potential which has to be confronted with the typical “ellyptic” character of the oscillations around the minima. We hope to be able to understand in the future whether this is truly the reason for the difference in FPU and rotator systems.

We thank S. Lepri for several profitable discussions. One of us (AP) thanks Gat Omri for some remarks about our last comments. ISI in Torino is acknowledged for the kind hospitality during the Workshop on ”Complexity and Chaos” 1998, where part of this work was performed.