15, CHAOS 015119 2005
Heat conduction in the FrenkelKontorova model
Bambi Hu
Department of Physics, Centre for Nonlinear Studies, and The Beijing-Hong Kong-Singapore Joint Centre for Nonlinear and Complex Systems (Hong Kong), Hong Kong Baptist University, Kowloon Tong, Hong Kong, China and Department of Physics, University of Houston, Houston, Texas 77204-5005
Lei Yanga
Department of Physics, Centre for Nonlinear Studies, and The Beijing-Hong Kong-Singapore Joint Centre for Nonlinear and Complex Systems (Hong Kong), Hong Kong Baptist University, Kowloon Tong, Hong Kong, China
Received 28 October 2004; accepted 10 January 2005; published online 28 March 2005 Heat conduction is an old yet important problem. Since Fourier introduced the law bearing his name almost 200 years ago, a first-principle derivation of this simple law from statistical mechanics is still lacking. Worse still, the validity of this law in low dimensions, and the necessary and sufficient conditions for its validity are far from clear. In this paper we will review recent works on heat conduction in a simple nonintegrable model called the FrenkelKontorova model. The thermal conductivity of this model has been found to be finite. We will study the dependence of the thermal conductivity on the temperature and other parameters of the model such as the strength and the periodicity of the external potential. We will also discuss other related problems such as phase transitions and finite-size effects. The study of heat conduction is not only of theoretical interest but also of practical interest. We will show various recent designs of thermal rectifiers and thermal diodes by coupling nonlinear chains together. The study of heat conduction in low dimensions is also important to the understanding of the thermal properties of carbon nanotubes. 2005 American Institute of Physics. DOI: 10.1063/1.1862552 Heat conduction is an old yet important problem. Since Fourier introduced the law bearing his name almost 200 years ago, a first-principle derivation of this simple law from statistical mechanics is still lacking. Worse still, the validity of this law in low dimensions has recently been called into question. In this work we will review recent works on heat conduction in a simple model called the FrenkelKontorova model. We will study the thermal conductivity, its dependence on the temperature and various parameters of the model, phase transitions, and finite-size effects. As potential applications, we will show various designs of thermal rectifiers and diodes. This study will be of relevance to carbon nanotubes.
I. INTRODUCTION
The Fourier law describing heat conduction states that the heat flux J is proportional to the temperature gradient T and the coefficient of proportionality is the thermal conductivity : J=- T. 1
It has been almost 200 years since Fourier proposed this phenomenological law that a first-principle derivation from statistical mechanics is still lacking. Worse still, even the validity of the Fourier law in low dimensions is far from clear. In integrable systems, Lebowitz et al.1 proved that a temperature gradient cannot be established and thus the Fourier law
a
Electronic mail: [e-mail]
is not valid. However, in nonintegrable systems, the situation is far more complex. In some nonintegrable systems, the Fourier law is obeyed; yet in other nonintegrable systems, the Fourier law is not obeyed. Indeed it would be very desirable to find out the complete set of necessary and sufficient conditions for the validity of the Fourier law. In this paper, we will study heat conduction in a nonintegrable model called the FrenkelKontorova FK model.2 The FK model was first proposed by Frenkel and Kontorova3 in 1938 to study surface phenomena. Since then it has found applications in a wide variety of physical systems, such as adsorbed monolayers, Josephson junctions, charge density waves, magnetic spirals, tribology, and DNA denaturation. Despite its deceptively simple form, the FK model exhibits very rich and complex behaviors. Its heat-conduction property413 is a case in point. The study of heat conduction is not purely of theoretical interest; it also has practical implications. The study of electric currents has led to the invention of electric rectifiers, diodes, and transistors. However, virtually no study of thermal rectifiers, diodes, or transistors has been made. Very recently Terraneo, Peyrard, and Casati12 pointed out the possibility of the design of a thermal diode by coupling two nonlinear chains. However, the gain of this thermal diode was only a factor of 2 or 3. Using instead two FK chains, Li, Wang, and Casati13 have been able to raise the gain to a factor of 100. We will present another design of controlled thermal conductors.10 In this model the conductivity can change by a factor of 100 00 and thus one can control the
2005 American Institute of Physics
1054-1500/2005/15 1 /015119/9/$22.50
15, 015119-1
015119-2
B. Hu and L. Yang
Chaos 15, 015119 2005
FIG. 1. Normalnormal plot of N : a is for the FK model, for the FPU model. Here, T = Th + Tl / 2.
= 3.0, b is
FIG. 2. a Loglog plot of the heat flux correlation function CJ t = J t J 0 ; b plot of CJ t . N = 4096, T = 0.5.
system from a good thermal conductor to a good thermal insulator. The study of heat conduction is also of importance to the study of carbon nanotubes. Although most studies on carbon nanotubes have been focused on their electric conduction properties, the thermal conduction properties of carbon nanotubes are equally important. So the study of heat conduction in low dimensions is also of practical interest insofar as nanotechnology is concerned.
II. THERMAL CONDUCTIVITY OF THE FRENKELKONTOROVA MODEL
Heat conduction in classical systems has attracted much attention recently. This problem relates macroscopic physical phenomena to their microscopic statistical properties see Ref. 14 for a recent review . A general Hamiltonian for the study of heat conduction is as follows:
N
H=
i=1
p2 i + U qi+1,qi + V qi , 2m
2
the temperature difference at the two ends of a chain, j the heat flux, and J = jN the total heat flux. The thermal conductivity is infinite if system 2 is harmonic1 or integrable. In this case, the energy carriers phonons or solitons are not scattered and they propagate ballistically. When a temperature difference T is applied to the two ends of a chain, a constant heat flux j independent of N sets in. The conductivity is found to be proportional to N. For nonintegrable models, Prosen and Campbell15 have proved a theorem classifying them into two categories according to the conservation or nonconservation of momentum. In nonintegrable models without an external potential, momentum is conserved and there exist long wavelength modes. Ballistic transport in the infinite wavelength limit may happen and the thermal conductivity diverges with the system size as N . In nonintegrable models with an external potential, translational invariance is broken and momentum is not conserved. The thermal conductivity is finite. However, there are also exceptions to this theorem.16 In the FrenkelKontorova model,
N
where N is the total number of the particles, m = 1 the mass of particles, pi the momentum of the ith particle, qi its displacement, U the interparticle potential, and V the external potential. is the strength of the interparticle potential, here = 1. is the strength of the external potential. The thermal conductivity is defined as = J / T = jN / T, where T is
H=
i=1
p2 i + qi+1 - qi - l0 2 + 1 - cos qi 2m 2 2 a
,
3
where 2a is the periodicity of the external potential, l0 = 2.0 the natural length of the spring, l the atomic mean distance,
015119-3
Heat conduction in FK model
Chaos 15, 015119 2005
FIG. 4. Loglog plot of
N.
ture for -N0 i 0 and Ti = Tl low temperature for N i N + N0. We set N0 = 8, = 0.1, and Th - Tl / Tl 0.1. To obtain a steady state, the total integration time is typically 108 109 units. Here, the total heat flux J = jN is the onedimensional version of the general expression of the total heat flux.18 We have checked that this is sufficient for the system to reach a steady state since the temperature profile in the central region is linear and the time-averaged heat flux ji is independent of the site i.
FIG. 3. Two temperature profiles.
and l / 2a the winding number. When l / 2a is a rational number, the FK model is in the commensurate phase. When l / 2a is an irrational number, the FK model is in the incommensurate phase. The ground state of the incommensurate FK model undergoes a phase transition by the breaking of analyticity.17 Most studies of heat conduction in the FK model consider the commensurate phase. The thermal conductivity is found to be finite. In the commensurate phase, translational invariance is broken and an acoustic branch does not exist in the phonon spectrum. Thus it is reasonable to expect a finite thermal conductivity. In our simulations, we use a fixed boundary condition and the Langevin thermostat. We consider a chain of N0 + N + N0 oscillators. The central N oscillators follow the Hamiltonian equations of motion while the outer 2N0 ones satisfy qi = - qi - qi + i ,
where = U + V is the total potential, i the white Gaussian noise, i t k t = 2 Ti ik t - t . Ti = Th high tempera-
FIG. 5. a Loglog plot of
e
N . b Plot of
e
T.
015119-4
B. Hu and L. Yang
Chaos 15, 015119 2005
FIG. 6. a Loglog plot of e N and b N . b Loglog plot of . varies from 0.1 to 3.0. Th = 1.0, Tl = 0.9. b
e
and
FIG. 7. Loglog plot of N and varies from 128 to 8192. a , b b l / 2a for N = 2048.
b
. Th = 1.0, Tl = 0.9, = 2.0, and N N and e N for different l / 2a; c
In Fig. 1, we plot N . For comparison, we have also plotted N for the FPU model with the following potential: U qi+1,qi = V qi = 0. The thermal conductivity is found to be finite in the FK model and divergent in the FPU model. We have also performed direct microcanonical simulations. The thermal conductivity is given by the GreenKubo formula 1 NT2
t
1 2
qi+1 - qi - l0 2 +
1 4
qi+1 - qi - l0 4 , 4
In Fig. 2, CJ t shows an exponential decay. This result agrees with that in Ref. 6. Thus the thermal conductivity is finite for mid-values of the parameters in the commensurate phase.
III. TEMPERATURE DEPENDENCE OF THE THERMAL CONDUCTIVITY
= lim lim
t N
CJ t dt,
0
where CJ t = J t J 0 is the total heat flux correlation function. We consider a finite chain N = 4096 with a periodic boundary condition. First, the chain is embedded in a Langevin heat bath with a temperature T = 0.5. After the system has reached thermal equilibrium, the thermostat is removed and then the decay of the heat flux correlation function of the isolated chain is measured. To increase the accuracy of the correlation function, its averaged value is taken over 3000 different realizations of the initial thermalization of the chain.
Since heat conduction in the FK model is normal, it is interesting to study its temperature dependence. When the temperature tends to zero, the FK model tends to a harmonic chain with a sinusoidal external potential. In this case the thermal conductivity should diverge.14 Then it is natural to ask whether there is a nonequilibrium phase transition in the low-temperature regime.11 When the system size is finite, there is a boundary jump of the temperature profile. It is nontrivial19 that the conductivity is influenced by the boundary jump. By investigating the temperature profiles, we have made two observations: 1 the boundary jump decreases when the system size N increases; 2 the boundary jump increases when the temperature of the heat bath decreases. In Fig. 3, two temperature profiles are shown. For an estimate of the finite-size effect of the boundary jump, we can calculate the effective conductivity e by using the gradient averaged only over the central linear region. We find that the temperature profile in the re1 4 gion 5 N , 5 N is usually linear. We take this part as the cen-
015119-5
Heat conduction in FK model
Chaos 15, 015119 2005
FIG. 8. Frequency spectra of particle vibrations for different a.
FIG. 9. Frequency spectra of particle vibrations for different . N = 256.
tral linear region. Here = b, the bulk conductivity. We can estimate the finite-size effect of the boundary jump by comparing e and b. For = 3.0, we plot e N and b N at T = 0.75, 1.5, 3.0 for N = 32 2, 048. In Fig. 4, for T = 3.0, e N and b N overlap when N 256. The numerical results show clearly e N = b N = constant. For T = 1.5, the difference between 512. The results e N and b N becomes smaller when N show e N = b N = constant also. For T = 0.75, the difference between e and b remains noticeable until N = 2048. The numerical results show that e and b increase as N increases so it is hard to draw a clear conclusion. But it seems reasonable to expect that e and b will approach a constant as N 5000. These results indicate that the finitesize effect becomes larger as T decreases. One needs to simulate a very long chain to obtain a confirmed conductivity when T is very small. Hence, it is important to take into account the finite-size effect of the boundary jump. If we only calculate b for N 512, then one is tempted to conclude that there is a phase transition7 in the regime 0.75,3.0 , as shown in the dash window in Fig. 4. In the FK model, as the temperature decreases, the system approaches the harmonic limit of the external potential. The thermal conductivity increases. On the other hand, when the temperature becomes very high, each particle can jump freely between the minima of the external potential. The thermal conductivity increases also. So there should exist a temperature at which the thermal conductivity reaches a minimum.6 This is shown in Fig. 5. Figure 5 a shows e N for = T 1.5 to 75.0. e N becomes flat when N 1024. e T is plotted in Fig. 5 b for N = 1024. The minimum of thermal conductivity is found at T 4.5. When T 4.5, e T is a monotonically decreasing function. When T 4.5, e T is a monotonically increasing function. Calculating both e and b gives us an estimate of the finite-size effect of the boundary jump. The finite-size effect becomes larger as T decreases. One therefore needs to simulate a very long chain to confirm the results. So from numerical simulations alone it is very difficult to obtain a clear conclusion on the existence or nonexistence of a phase
transition.7 Furthermore, our simulations show that there is no phase transition in very large regimes of the temperature. We believe that the scattering of the energy carrier is related to the thermal fluctuation of the heat bath, and the probability of the scattering may be proportional to e-1/T, as in the rotator model.9 So, even if a phase transition exists, it probably will occur only at an extremely low temperature.
IV. PARAMETER DEPENDENCE OF THE THERMAL CONDUCTIVITY AND CONTROLLED THERMAL CONDUCTORS
The thermal conductivity depends on various parameters of the FK model. In particular, the strength9 and the periodicity10 of the external potential can affect the thermal conductivity. This information can be used to design a controlled thermal conductor. It is clear the strength of the external potential can change the thermal conductivity at a fixed temperature. There are two limiting cases. As the strength of the external potential 0, the FK model becomes a purely harmonic system. Its thermal conductivity is infinite. As , the FK model becomes a chain of independent oscillators. It is a 0.1, 3 , the simulathermal insulator. In the mid-range, tion results are plotted in Fig. 6. In Fig. 6 a , for = 1, e N and b N coincide when N 128. The numerical results show that e N = b N = constant. For = 0.5, the difference between e N and 512. The results show e N b N decreases when N = b N = constant also. For = 0.1, the difference between e and b remains noticeable until N = 2048. e and b increase as N increases but the slope flattens. It seems a reasonable conclusion that e and b will approach a constant as N . These results indicate that the finite-size effect of the boundary jump becomes larger for a smaller . One needs to simulate a very long chain to obtain a confirmed thermal conductivity when is very small. The relation between and is shown in Fig. 6 b . It is expected that the thermal conductivity is a monotonically decreasing function of . = 0.1 / = 3.0 5000. It means that one Here the ratio can adjust to change the thermal conductivity. In our simulations, a phase transition cannot be found for 0.1, 3 .
015119-6
B. Hu and L. Yang
Chaos 15, 015119 2005
FIG. 10. Sketch of a diode model.
We now study how the thermal conductivity depends on the periodicity of the external potential 2a as a is varied from 0.5 to 8. First consider two limiting cases when and T take mid-values. When a , the relative atomic displacement q / a = qi+1 - qi / a tends to zero. The FK model approaches its continuum limit, the sine-Gordon chain. The sine-Gordon chain is an integrable system. Its energy and momentum are conserved. Thus the thermal conductivity should be infinite. On the other hand, when a 0, a particle needs to overcome an infinite potential well. Thus the thermal conductivity should be zero. We therefore expect the thermal conductivity would change from zero to infinity when a changes from zero to infinity. The simulation results are plotted in Figs. 7 a 7 c . In Fig. 7 a , b N is shown for a = 8, 4, 8 / 3, 2, 8 / 5, 4 / 3, 1, 1 / 2, and 34/ 21. For a 8 / 3, b N is nearly a constant. For a = 4, b N increases as N increases, but the slope flattens. In Fig. 7 b , the finite-size effect of the boundary jump is shown for a = 8 , 4 , 2. For a = 2, b N and e N coincide for all N and they are nearly flat. For a = 4, b N and
FIG. 12. Loglog plot of the total positive heat flux J+ N , the total negative heat flux J- N , and their ratio = J+ N / J- N on the diode model.
FIG. 11. Temperature profiles of the diode model.
1024 and approaches a constant e N coincide when N when N 1024. For a = 8, b N and e N do not overlap even for our longest simulation N = 8192. It means that the finite-size effect from the boundary jump exists for N = 128 to 8192. It seems that the finite-size effect may disappear when N 20 000. For a larger a, we need a larger-size simulation to avoid the finite-size effect. So it is hard to say whether a phase transition exists by numerical simulations alone. The simulation results suggest that, for a = 1 / 2 to 4, the thermal conductivity is finite. The finite-size effect becomes larger when a increases. For a = 8, although the finitesize effect exists, the slope of b N and e N has a small decrease as N increases. Thus the thermal conductivity may very likely be finite as N . In Fig. 7 c , the relation between and 1 / a is shown. The conductivity is a monotonically decreasing function of 1 / a. It is reasonable to assume that the soft modes including the zero-mode play an important role in heat conduction.20 We have studied the frequency spectra of the particle vibrations as we vary a and . In Fig. 8, the frequency spectra move to the low frequencies as a increases. So there are more soft modes and they have longer mean free path. Thus the system has a larger thermal conductivity. This agrees with Fig. 7 c . In Fig. 9, the frequency spectra move to the high frequencies as increases. This agrees with Fig. 6 b . By adjusting the strength or the periodicity of the external potential, we can design a controlled thermal conductor.
015119-7
Heat conduction in FK model
Chaos 15, 015119 2005
FIG. 13. Loglog plot of the positive heat flux j + N , the negative heat flux j- N .
The thermal conductivity of the controlled heat conductors can change nearly 10 000 times. It means that one can adjust or a to change the FK model from a thermal conductor to a thermal insulator. In practical applications, changing the periodicity of the external potential may be a more economical way than changing the strength of the external potential.
V. THERMAL DIODE AND RECTIFIER
FIG. 14. Frequency spectra of particle vibrations in the diode model. N = 256.
Thermal rectification can occur, for example, in a layer of fluid. Heating a layer of fluid from below, convection occurs and the current can transfer heat from below to above. However, if the fluid is heated from above, convection cannot occur.12 The thermal current from above to below is much small than before. Terraneo et al.12 suggested the possibility of controlling heat flow in a nonlinear chain connected to two thermostats at different temperatures. A strongly nonlinear chain is sandwiched between two weakly anharmonic chains. The phonon bands of the left and right chains do not change significantly when the orientation of the gradient is reversed. In the central part, its phonon band undergoes a big change when the direction of the gradient is reversed. Recently, Ref. 13 presented another design of a thermal diode. The overlap of the phonon bands of the two parts leads to the rectifying effect. They show that the ratio of the heat fluxes in opposite directions is 100 when the system size is about 100. Here we show a different design of a thermal diode. Two FK chains are coupled by a harmonic interparticle potential with strength k. The two FK chains have the same interpar1 ticle potential U qi+1 , qi = 2 qi+1 - qi - l0 2, = 1.0 but different external potentials Left: V qi = 1 1 - cos qi , 2 aL =
L,
profiles with two heat baths in the original and reversed directions. Here N = 128, k = 0.1, Th = 0.135, Tl = 0.035. We have also studied the size dependence of the thermal diode. The total positive and negative heat flux are shown in 250. When Fig. 12. In Fig. 12 a , k = 0.05. When N = 64, N = 2048, 5. In Fig. 12 b , k = 0.1. When N = 64, 200. 2. In both cases, decreases by increasWhen N = 2048, ing N. In Fig. 13, the positive and negative heat flux are shown. We study the frequency spectra of the particle vibration of the diode model in Fig. 14. The solid line shows the frequency spectra of the left part, the dot line shows the frequency spectra of the right part. In Fig. 14 a , the two frequency bands have a large overlap. Thus the system is a thermal conductor. In Fig. 14 b , the two frequency bands have nearly no overlap. Thus the system is essentially a thermal insulator. It is also possible to design a thermal rectifier by connecting two FK chain directly:
N/2
H=
i=1 N
p2 1 1 i + qi+1 - qi - l0 2 + 1 - cos qi 2 2 2 aL p2 1 1 i + qi+1 - qi - l0 2 + 1 - cos qi 2 2 2 aR , 6
5 1 Right: V qi = 1 - cos qi , 2 aR =
R,
+
i=N/2
where aL = 0.5, L = 0.5 and aR = 2.0, R = 2.0. A sketch of the model is shown in Fig. 10. Figure 11 shows the temperature
where aL = 16, aR = 8 / 3. Here, Th = 1.0 and Tl = 0.1. The heat flux ratio in the two opposite directions is just about 3.5.
015119-8
B. Hu and L. Yang
Chaos 15, 015119 2005
temperature profiles of the heat baths in the original and reversed directions. The temperatures of the heat baths are Th = 2.0 and Tl = 0.5. The heat flux ratio between the two opposite directions is about 2.5. These designs of the thermal diode and the thermal rectifier are based on the mismatch of the left part and right part of the phonon bands at different temperatures.12 In model 5 , the ratio of the heat fluxes in opposite directions is about 200.
VI. SUMMARY
FIG. 15. Temperature profiles of the rectifier model 6 . N = 128.
Figure 15 shows the temperature profiles of the two heat baths in the original and reversed directions. We will also show the design of a thermal rectifier similar to that of Ref. 12. The thermal rectifier consists of three parts. Each part has the same Hamiltonian but different periodicity. The left part is a weakly nonlinear chain, a = 16; the central part is a strongly nonlinear chain, a = 8 / 3; the right part is a weakly nonlinear chain, a = 8. Figure 16 shows two
In the paper, we have reviewed some recent works on heat conduction in the FrenkelKontorova model. The thermal conductivity of the FK model in the commensurate phase is finite. The total momentum of the system is not conserved because translational invariance is broken by the external potentials. We have also studied the dependence of the thermal conductivity on the temperature and other parameters of the model. In particular, we have investigated how the strength and the periodicity of the external potential affect the thermal conductivity. We have also discussed the phase transition problem and the finite-size effect. Our simulations show that there is no phase transition in very large regimes of the parameters. We suggest that calculating both the effective conductivity and the bulk conductivity will provide an estimate of the finitesize effect of the boundary jump. For potential applications, we have presented various designs of controlled thermal conductors. By adjusting the periodicity of the external potential, the thermal conductivity can change nearly 10 000 times. Therefore one can control the system from a good thermal conductor to a thermal insulator. We have also presented one design of a thermal diode and two designs of a thermal rectifier. By coupling two or three FK chains with different periodicity, we show how the thermal diode and the thermal rectifier work.
ACKNOWLEDGMENTS
We would like to thank Dr. D. Campbell, Dr. G. Casati, Dr. D. He, Dr. B. Li, Dr. Y. Zhang, and members of the Centre for Nonlinear Studies for useful discussions. This work was supported in part by grants from the Hong Kong Research Grants Council RGC and the Hong Kong Baptist University Faculty Research Grant FRG .
Z. Rieder, J. L. Lebowitz, and E. Lieb, J. Math. Phys. 8, 1073 1967 ; R. J. Rubin and W. L. Greer, ibid. 12, 1686 1971 . O. M. Braun and Y. S. Kivshar, The Frenkel-Kontorova Model: Concepts, Methods, and Applications Springer, Berlin, 2003 . 3 Ya. Frenkel and T. Kontorova, Zh. Eksp. Teor. Fiz. 8, 89 1938 . 4 M. J. Gillan and R. W. Holloway, J. Phys. C 18, 5705 1985 . 5 B. Hu, B. Li, and H. Zhao, Phys. Rev. E 57, 2992 1998 . 6 G. P. Tsironis, A. R. Bishop, A. V. Savin, and A. V. Zolotaryuk, Phys. Rev. E 60, 6610 1999 . 7 A. V. Savin and O. V. Gendelman, Phys. Rev. E 67, 041205 2003 . 8 L. Yang and P. Grassberger, arXiv: cond-mat/0306173 2003 . 9 L. Yang, contributed talk given in Dynamics Days Asia-Pacific: The Third International Conference on Nonlinear Science, Singapore, 2004. 10 B. Hu and L. Yang, preprint CNS-04-9 of the Centre for Nonlinear Studies, Hong Kong Baptist University, submitted for publication August, 2004 . 11 B. Hu and L. Yang, "Phase Diagram of the Frenkel-Kontorova Model,"
2 1
FIG. 16. Temperature profiles of the rectifier model 7 . N = 128.
015119-9
Heat conduction in FK model
Chaos 15, 015119 2005 Giardina, R. Livi, A. Politi, and M. Vassalli, ibid. 84, 2144 2000 . S. Aubry, Springer Ser. Solid-State Sci. 8, 264 1978 . 18 R. Kubo, M. Toda, and N. Hashitsume, Statistical Physics II, in Springer Series in Solid State Sciences Vol. 31 Springer, Berlin, 1991 . 19 K. Aoki and D. Kusnezov, Phys. Rev. Lett. 86, 4029 2001 . 20 P. Grassberger, W. Nalder, and L. Yang, Phys. Rev. Lett. 89, 180601 2002 .
17
preprint CNS-04-7 of the Centre for Nonlinear Studies, Hong Kong Baptist University, 2004. 12 M. Terraneo, M. Peyrard, and G. Casati, Phys. Rev. Lett. 88, 094302 2002 . 13 B. Li, L. Wang and G. Casati, Phys. Rev. Lett. 93, 184301 2004 . 14 S. Lepri, R. Livi and A. Politi, Phys. Rep. 377, 1 2003 . 15 T. Prosen and D. K. Campbell, Phys. Rev. Lett. 84, 2857 2000 . 16 O. V. Gendelman and A. V. Savin, Phys. Rev. Lett. 84, 2381 2000 ; C.