Temperature dependence of thermal conductivities of coupled rotator lattice and the momentum diffusion in standard map

0.22 MB·English

by  Yunyun Li

Temperature dependence of thermal conductivities of coupled rotator lattice and the
momentum diffusion in standard map
Yunyun Li,1 Nianbei Li,1,∗ and Baowen Li1,2,3,†
1Center for Phononics and Thermal Energy Science and School of Physics Science and Engineering,
Tongji University, 200092 Shanghai, People’s Republic of China
2Department of Physics and Centre for Computational Science and Engineering,
National University of Singapore, Singapore 117546, Republic of Singapore
5 3NUS Graduate School for Integrative Sciences and Engineering, Singapore 117456, Republic of Singapore
1
0 Incontrarytoother1Dmomentum-conservinglatticessuchastheFermi-Pasta-Ulamβ (FPU-β)
2 lattice,the1Dcoupledrotatorlatticeisanotableexceptionwhichconservestotalmomentumwhile
n exhibitsnormalheatconductionbehavior. Thetemperaturebehaviorofthethermalconductivities
a of 1D coupled rotator lattice had been studied in previous works trying to reveal the underlying
J physical mechanism for normal heat conduction. However, two different temperature behaviors
of thermal conductivities have been claimed for the same coupled rotator lattice. These different
9
2 temperature behaviors also intrigue the debate whether there is a phase transition of thermal con-
ductivities as thefunction of temperature. In this work, we will revisit the temperature dependent

] thermalconductivitiesforthe1Dcoupledrotatorlattice. Wefindthatthetemperaturedependence
h follows a power law behavior which is different with the previously found temperature behaviors.
c Our results also support the claim that there is no phase transition for 1D coupled rotator lattice.
e
We also give some discussion about the similarity of diffusion behaviors between the 1D coupled
m
rotator lattice and the single kicked rotator also called theChirikov standard map.
-
t
a PACSnumbers: 05.60.-k,44.10.+i,05.45.-a
t
s
.
t I. INTRODUCTION pendence of thermal conductivity of κ(T) ∝ e−T/A was
a
m claimed where A is a positive constant. There it was
also argued that a possible phase transition at temper-
- The exploration of underlying mechanism for anoma-
d ature around T ∼ 0.2 − 0.3 where heat conduction is
lousandnormalheatconductioninlowdimensionalsys-
n normal above this temperature while anomalous below
tems represents a huge challenge in the area of statisti-
o this temperature exists [62]. Although this phase transi-
c cal physics [1–4]. After enormous efforts for more than
tion statement was challenged as a finite size effect [63],
[ onedecadefromnumericalsimulations[5–44],theoretical
a later work supported this phase transition conjecture
predictions[45–57]andexperimentalobservations[58–60],
1 after deriving a similar temperature dependent thermal
v there is still no consensus for the exact physical mecha- conductivity as κ(T)∝e−T/A [50].
9 nism causing anomalous heat conduction. It is believed
1 that momentum conservation plays an important role Most recently, the 1D coupled rotator lattice re-
4 in determining the actual heat conduction behavior. In attracts much attention due to the new finding of simul-
7 general,1D non-integrablelattices with momentum con- taneouslyexistingnormaldiffusionofmomentumaswell
0 servingpropertyshouldhaveanomalousheatconduction as heat energy [64]. It is argued that the normal be-
1. where the thermal conductivity κ diverges with the lat- haviormight be due to the reducednumber ofconserved
0 tice size N as κ ∝ Nα where 0< α < 1 [1–4]. However, quantitieswhichisuniqueforperiodicinteractionpoten-
5 the 1Dcoupledrotatorlattice is a wellknownexception. tials [65, 66]. But the new debate is that whether the
1 It exhibits normal heat conduction behavior despite its stretch is conservedor not in 1D coupled rotator lattice.
:
v momentum conserving nature [61, 62]. The investigation of 1D coupled rotator lattice will be
i the key to unravel the true mechanism behind the con-
X The normal heat conduction in 1D coupled rotator
nection between momentum conservation and normal or
r lattice was discovered via numerical simulations by two anomalous heat conduction. Therefore, it is the right
a groups independently [61, 62]. In order to understand
time to revisit the temperature dependence of thermal
the underlying mechanism, the temperature dependence
conductivity of 1D coupled rotatorlattice as a first step.
of thermal conductivity has been studied in detail in
bothworks. InRef. [61], the temperature dependence of In this work, we will revisit the temperature depen-
thermal conductivity of coupled rotator lattice has been dence of thermal conductivities for the 1D coupled rota-
found to be like κ(T) ∝ e∆V/T where ∆V is a posi- tor lattice. We find that the temperature dependence is
tive constant. However, in Ref. [62], a temperature de- neither κ(T) ∝ e∆V/T nor κ(T) ∝ e−T/A as previously
claimed [61, 62]. The actual temperature dependence is
a power-law dependence as κ(T) ∝ T−3.2. The possible
connectionwiththemomentumdiffusionofsinglekicked
∗Electronicaddress: [email protected] rotator or the Chirikov standard map has also been dis-
†Electronicaddress: [email protected] cussed. In order to determine whether there is a phase
2
transition, we also present the temperature dependent
thermal conductance at different system sizes. All the 0.11
thermal conductances for different sizes collapse to the T= 0.1
same value at low temperatures while approach to the Ti
power-law behavior as κ(T) ∝ T−3.3 at high tempera- 0.10
tures. However, the crossover temperature decreases as
thesystemsizeincreases. Thisfactindicatesthatthereis
(a)
no phase transitionbetween normaland anomalousheat 0.09
0.0 0.5 1.0
conduction. In thermodynamical limit, the heat conduc- i/ N
0.9
tionisnormalforalltemperaturesexceptthetrivialzero N=50
T= 0.8 N=100
temperature point. In Sec. II the lattice model will be
N=200
introduced. Numerical results and discussions will be Ti N=400
0.8 N=800
presented in Sec. III and the results will be summarized
N=1600
in Sec. IV.
(b)
0.7
0.0 0.5 1.0
II. MODEL i/N
FIG. 1: (color online). Temperature profiles for the 1D cou-
pled rotator lattice at (a) T = 0.1 and (b) T = 0.8. The
The Hamiltonian for the 1D coupled rotator lattice is
lattice sizes are N = 50,100,200,400,800 and 1600 and the
defined as the following [61, 62]:
color ruleof thelines arethesamefor both (a) and(b). The
firstandlast atom areputintocontactwith aLangevinheat
N p2 bathwithtemperaturesetasT =T(1±∆)where∆=0.1
H = i +K[1−cos(q −q )] (1) L/R
(cid:20) 2 i+1 i (cid:21) here.
Xi=1
whereq andp denotethedisplacementfromequilibrium
i i This is because the kinetic energy is proportional to the
and momentum for i-th atom, respectively. The mass of
temperatureas p2 =T andthepotentialenergyinEq.
theatommandtheBoltzmannconstantk hasbeenset i
B (1) is confined b(cid:10)y t(cid:11)he cosine function.
intounity. TheparameterK withenergydimensionrep-
Withoutgettingintonumerics,wecangetaqualitative
resentsthecouplingstrengthoftheinter-atompotential.
picture of the thermal conductivity of 1D coupled rota-
Therefore, the system temperature T can be rescaled by
tor lattice as the function of temperature. The thermal
T/K and we can also set K =1 for simplicity [67]. The
conductivity κ(T) will diverge in the low temperature
equations of motion for i-th atom are
limit approaching to the harmonic limit and will decay
to zero in the high temperature limit as approaching to
q˙ = p
i i the unconnected N free particles. As a general picture,
p˙i = Ksin(qi+1−qi)−Ksin(qi−qi−1) (2) the thermal conductivity κ(T) will decrease as the tem-
perature increases.
At low temperature limit, the displacements are small
In non-equilibrium numerical simulations, we put the
values so that |q −q | ≪ 1. The Hamiltonian can be
i+1 i firstandlastatomintocontactwithLangevinheatbath.
expanded into the Taylor series as:
The temperatures for left and right heat bath are set as
T =T(1±∆) where T denotes the averagetempera-
H = N p2i + (qi+1−qi)2 tuLr/eRand ∆ = 0.1 gives rise to the temperature gradient
(cid:20) 2 2 along the lattice. The fixed boundary conditions with
Xi=1
q =q =0 are also been applied.
(q −q )4 (q −q )6 0 N+1
i+1 i i+1 i
− + −... (3)
4! 6! (cid:21)
III. RESULTS AND DISCUSSIONS
This Hamiltonian will approach to the integrable Har-
monic lattice only at zero temperature T = 0. The first
Before we discuss the results of thermal conductivi-
stable nonlinear potential term will be the sextic poten-
ties, we first show the temperature profiles. In Fig. 1,
tial.
the temperature profiles at two different temperatures
On the other hand, at infinite high temperature limit
T = 0.1 and T = 0.8 are plotted. The lattice sizes are
T =∞,the rotatorlattice willapproachtoanotherinte-
chosen as N = 50,100,200,400,800and 1600. For rela-
grable system consisting of N independent and isolated
tive low temperature at T =0.1, the temperature jumps
free particles as
at two boundaries are obvious. However, the tempera-
ture jumps are reduced with the increase of lattice size
N p2
H = i (4) N as can be seen in Fig. 1(a). In Fig. 1(b) where the
Xi=1 2 temperature is relatively high at T = 0.8, all the tem-
3

100 N=1600
- 3.2 10
~ T
10
- 3.2
~ T
N=50 1
N=100
1 N=200
N=400
N=800
N=1600 0.1
0.1
0.01 0.1 1 1 2
T
T
FIG. 2: (color online). Thermal conductivities κ(T) as FIG. 3: (color online). Thermal conductivities κ(T) as the
the function of temperature for different lattice sizes N = function of temperature for lattice size N = 1600. The data
50,100,200,400,800 and 1600. The straight line is propor- aretakenfromFig. 2andcanbebestfittedtobeapower-law
tional to T−3.2 which describes the temperature behavior of dependenceas κ(T)∝T−3.2.
κ(T) at high temperatureregion very well.
more than two orders of magnitudes for the κ(T) value
peratureprofilescollapseto the samestraightline asthe as shown in Fig. 3. This also explains the poor fitting
temperature jumps are very small for all lattice sizes. of the κ(T) ∝ e∆V/T dependence in Ref. [61] and the
Inordertoobtainthetemperaturedependenceofther- κ(T)∝e−T/A dependence in Ref. [62].
malconductivities,weneedfirstdefinethewayhowκ(T) Thispower-lawdependenceofκ(T)∝T−3.2cannotbe
can be calculated numerically. We notice that the tem- explained by the effective phonon theory which is able
perature profiles all collapse to the same straight line if to predict the temperature dependent thermal conduc-
the temperature jumps can be ignored at high tempera- tivities for other typical 1D nonlinear lattices such as
tures or long lattice sizes. It is appropriate to define the FPU-β lattice and generalized nonlinear Klein-Gordon
thermal conductivity κ(T) as: lattices [68–75]. According to the effective phonon the-
ory,thethermalconductivityforlowtemperaturerotator
JN
lattice with Hamiltonian of Eq. (3) can be derived as
κ(T)= (5)
T −T
L R
1
κ(T)∝ ∝T−2 (6)
whereJ =hJiiistheaverageheatfluxalongthelatticein ǫ
the stationary state and the local heat flux J is defined
i where ǫ is the nonlinearity strength with the following
as J = −p Ksin(q −q ) via the energy continuity
i i i+1 i temperature dependence
equation.
InFig. 2, the thermalconductivities κ(T)asthe func- (q −q )6
tion of temperature are plotted for different lattice sizes ǫ∝ i+1 i ∝T2 (7)
(cid:10)h(q −q )2i(cid:11)
N =50,100,200,400,800and1600. Athightemperature i+1 i
region, all the thermal conductivities κ(T) for different at low temperature region. Here we consider the sextic
lattice sizes collapse together indicating the saturation potential term as the lowest nonlinear term because the
of thermal conductivities as the increase of lattice size. dynamicsgovernedbythenegativequarticpotentialterm
This is characteristic for lattices with normal heat con- is unstable. Therefore, the actual temperature behavior
duction. As the temperature decreases,the thermalcon- ofκ(T)∝T−3.2 is steeper thanthe prediction ofκ(T)∝
ductivities first increases and then becomes flat. This is T−2 from effective phonon theory.
because the phonon mean free paths are getting longer Although it is difficult to unravel the exact physical
as the temperature reduces and the ballistic regime will mechanism behind the power-law dependence of ther-
be approached if the phonon mean free paths are longer mal conductivities for coupled rotator lattice, it is very
thanthelattice size. Theboundaryjumps willdominate helpful tolook intothe transportpropertiesofthe single
the temperature profiles as in Fig. 1(a) and the defini- kicked rotator (the Chirikov standard map) [76]. As the
tion of κ(T) of Eq. (5) will no longer reflect the actual name indicates, the coupled rotator lattice is a kind of
thermal conductivities. connected kicked rotators. The equations of motion for
At high temperatures, it is clearly seen that the ther- the single kicked rotator are
mal conductivities follows a power-law dependence as
κ(T)∝T−3.2. For the longestsize we consideredhere as pn+1−pn = Ksin(qn)
N =1600,this power-lawbehaviorcanbe best fitted for q −q = p (8)
n+1 n n+1
4

-1
10
10
-2
DP 1 - 3.2 10
~ T
N=50
-3 N=100
10 N=200 - 3.2
0.1 N=400 ~ T
N=800
-4 N=1600
10
0.01
1 2 3 0.01 0.1 1
T T
FIG. 4: (color online). Momentum diffusion constant DP as FIG. 5: (color online). Heat conductance σ as the
the function of temperature for 1D coupled rotator lattice. function of temperature T for different lattice size N =
ThenumericaldataareobtainedviaequilibriumMDsimula- 50,100,200,400,800 and 1600. All the other parameters are
tions asin Ref. [64]. Thesolid line of T−3.2 is guided for the thesame as in Fig. 2. Theleft-pointing arrow representsthe
eye. trend of decreasing crossover temperature as the lattice size
increases.
where q and p denote the coordinate and momentum
n n
after n-th kick. As a consistence check, we plot the momentum diffu-
The variation of momentum p is unbounded and can sion constant D for the coupled rotator lattice as the
P
be characterized by a normal diffusion [76, 77] functionoftemperatureinFig. 4. Thesametemperature
∆p2(t) ∼Dt (9) dependence of DP ∝T−3.2 has been obtained.
Therefore, the prediction of κ(T) ∝ T−3 with anal-
(cid:10) (cid:11)
where D is the diffusion constant. This is similar to the ogy to the single kicked rotator at high temperature or
normal momentum diffusion for 1D coupled rotator lat- low K region is close to the numerical observation of
tice. Accordingto Ref. [76,77],the diffusionconstantD κ(T)∝T−3.2. Thisindicatestheremightbesomedeeper
depends on the coupling strength K as connectionbetweenthediffusionbehaviorofsinglekicked
rotator and 1D coupled rotator lattice. In addition, the
D ∝ (K−1.2)3, 1.2<K <4 (10)
analogy analysis also predicts a κ(T)∝T−2 behavior at
D ∝ K2, K >4 (11) low temperature region which is the same as the predic-
tionfromtheeffectivephonontheoryofEq. (6). Therea-
where 1.2 is the chaos threshold for the kicked rotator.
sonwhythistemperaturebehaviorcannotbeobservedin
Ifoneassumethetransportpropertiesofcoupledrota-
numerical simulations might be due to the severe finite
tor lattice are the same as that of single kicked rotator,
size effect at low temperature region as can be seen in
one would expect that the energy diffusion constant D
E
Fig. 1(a).
as well as the momentum diffusion constant D of cou-
P
pled rotator lattice should also follow the same depen- In the final part we will discuss the issue about the
dence as that of single kicked rotator as in Eq. (10). To possiblephasetransitionbetweennormalandanomalous
translatetheK dependenceintoT dependence,wenotice heat conduction for 1D coupled rotator lattice. As we
that the parameter K in single kicked rotator of Eq. (8) haveshowninFig. 3,theactualtemperaturedependence
plays the same role as that in the coupled rotatorlattice ofκ(T)isapower-lawdependenceasκ(T)∝T−3.2. This
asin Eq. (2). As wehavediscussedabove,the dynamics indicatesthe analyticderivationofκ(T)∝e−T/A inRef.
ofcoupledrotatorlatticecanberescaledwithT/K. The [50] originally for a lattice with pinned on-site potential
K dependence should be inversely proportionalto the T can not be extended to 1D coupled rotator lattice.
dependence. LowK valueregioncorrespondstothehigh On the other hand, the claim that the transition tem-
T regionandviceverse. Andfornormalheatconduction, perature is around(0.2−0.3)in Ref. [62] has been chal-
the thermal conductivity κ is proportionalto the energy lenged in Ref. [63] with numerical simulations of longer
diffusion constant D . This will finally give rise to the sizes. The possible phase transitioncouldbe a finite size
E
prediction oftemperature dependent thermal conductiv- effect. Here we reconfirm the finite size effect by giving
ity κ(T)ofcoupledrotatorlattice fromthe analogywith a more illustrative picture in Fig. 5. The heat conduc-
single kicked rotator: tance σ ≡ κ(T)/N as the function of temperature has
been plotted for different lattice sizes. At high temper-
κ(T) ∝ T−2, low T (12)
ature regions, all the heat conductances follow the same
κ(T) ∝ T−3, high T (13) temperature behavior of σ ∝ T−3.2. As the tempera-
5
ture decreases,the phononmeanfree pathwillovercome connectionwiththesinglekickedrotatororthe Chirikov
the lattice size and the heat conductance will saturate standard map has been discussed where a κ(T) ∝ T−3
to a value determined by the harmonic limit of coupled dependence can be implied. Our results also recon-
rotator lattice. As can be seen from Fig. 5, the heat firmthatthepreviouslyclaimedpossiblephasetransition
conductance for short lattice will first bend and become should be a finite size effect.
flat as the temperature decreases. This crossover tem-
perature decreases as the lattice size increases smoothly
and no sign of phase transition can be observed. It is V. ACKNOWLEDGMENTS
also noticed that the crossover temperature happens to
be around(0.2−0.3)for lattice sizesof afew thousands. The numericalcalculations were carriedout at Shang-
hai Supercomputer Center. This work has been sup-
portedbytheNSFChinawithgrantNo. 11334007(Y.L.,
IV. SUMMARY
N.L., B.L.), the NSF China with grant No. 11347216
(Y.L), Tongji University under Grant No. 2013KJ025
We have systematically investigated the temperature (Y.L), the NSF China with Grant No. 11205114(N.L.),
dependence of thermal conductivities of 1D coupled ro- the Program for New Century Excellent Talents of the
tator lattice. The actual temperature dependence is a Ministry of Education of China with Grant No. NCET-
power-lawdependence ofκ(T)∝T−3.2 whichisdifferent 12-0409 (N.L.) and the Shanghai Rising-Star Program
with the observations of previous studies. The possible with grant No. 13QA1403600(N.L.).
[1] S.Lepri,R.Livi,andA.Politi,Phys.Rev.Lett.78,1896 [27] H. Zhao, Phys. Rev.Lett. 96, 140602 (2006)
(1997) [28] D. He, S. Buyukdagli, and B. Hu, Phys. Rev. E 78,
[2] S.Lepri,R.Livi,andA.Politi,Phys.Rep.377,1(2003) 061103 (2008)
[3] A.Dhar, Adv.Phys.57, 457 (2008) [29] Z.-G.Shao,L.Yang,W.-R.Zhong,D.-H.He,andB.Hu,
[4] S. Liu, X. F. Xu, R. G. Xie, G. Zhang and B. Li, Eur. Phys. Rev.E 78, 061130 (2008)
Phys.J. B 85, 337 (2013) [30] Y. Dubi and M. Di Ventra, Phys. Rev. E 79, 042101
[5] B.Hu,B.Li,andH.Zhao,Phys.Rev.E57,2992(1998) (2009)
[6] P.Tong,B.Li,andB.Hu,Phys.Rev.B59,8639(1999) [31] A. Henry and G. Chen, Phys.Rev. B 79, 144305 (2009)
[7] T. Hatano, Phys. Rev.E 59, R1 (1999) [32] K. Saito and A. Dhar, Phys. Rev. Lett. 104, 040601
[8] G. P. Tsironis, A. R. Bishop, A. V. Savin, and A. V. (2010)
Zolotaryuk, Phys.Rev.E 60, 6610 (1999) [33] L.Wang,D.He,andB.Hu,Phys.Rev.Lett.105,160601
[9] A. Sarmiento, R. Reigada, A. H. Romero, and K. Lin- (2010)
denberg,Phys. Rev.E 60, 5317 (1999) [34] N.Yang,G. Zhang,and B. Li,NanoToday 5, 85 (2010)
[10] A.Dhar and D.Dhar, Phys. Rev.Lett. 82, 480 (1999) [35] L. Wang and T. Wang, EPL 93, 54002 (2011)
[11] D. Alonso, R. Artuso, G. Casati, and I. Guarneri, Phys. [36] L. Wang, B. Hu, and B. Li, Phys. Rev. E 86, 040101
Rev.Lett. 82, 1859 (1999) (2012)
[12] B.Hu,B.Li,andH.Zhao,Phys.Rev.E61,3828(2000) [37] D. Xiong, J. Wang, Y. Zhang, and H. Zhao, Phys. Rev.
[13] K.AokiandD.Kusnezov,Phys.Lett.A265,250(2000) E 85, 020102(R) (2012)
[14] B.Li,H.Zhao,andB.Hu,Phys.Rev.Lett.86,63(2001) [38] G. T. Landi and M. J. de Oliveira, Phys. Rev. E 87,
[15] A.Dhar, Phys.Rev. Lett.86, 5882 (2001) 052126 (2013)
[16] K. Aoki and D. Kusnezov, Phys. Rev. Lett. 86, 4029 [39] D. Xiong, Y. Zhang, and H. Zhao, Phys. Rev. E 90,
(2001) 022117 (2014)
[17] Y.Zhang and H. Zhao, Phys. Rev.E 66, 026106 (2002) [40] S.G.Das,A.Dhar,K.Saito,C.B.Mendl,andH.Spohn,
[18] B. Li, L.Wang, and B. Hu,Phys.Rev.Lett. 88, 223901 Phys. Rev.E 90, 012124 (2014)
(2002) [41] C. B. Mendl and H. Spohn, Phys. Rev. E 90, 012147
[19] K.Saito, EPL 61, 34 (2003) (2014)
[20] A. V. Savin and O. V. Gendelman, Phys. Rev. E 67, [42] A.V.SavinandY.A.Kosevich,Phys.Rev.E89,032102
041205 (2003) (2014)
[21] S.Lepri,R.Livi,andA.Politi, Phys.Rev.E68,067102 [43] S. Liu, J. Liu, P. Hanggi, C. Wu, and B. Li, Phys. Rev.
(2003) B 90, 174304 (2014)
[22] D.Segal,A.Nitzan,andP.H¨anggi,J.Chem.Phys.119, [44] L. Wang, L. Xu, and H. Zhao, Phys. Rev. E 91, 012110
6840 (2003) (2015)
[23] O.V.Gendelman and A.V.Savin,Phys.Rev.Lett. 92, [45] S. Lepri, R. Livi, and A. Politi, Europhy. Lett. 43, 271
074301 (2004) (1998)
[24] B. Li, J. Wang, L. Wang, and G. Zhang, Chaos 15, [46] S. Lepri, Phys. Rev.E 58, 7165 (1998)
015121 (2005) [47] O. Narayan and S. Ramaswamy, Phys. Rev. Lett. 89,
[25] G.ZhangandB. Li,J. Chem.Phys.123, 114714 (2005) 200601 (2002)
[26] H. Zhao, Z. Wen, Y. Zhang, and D. Zheng, Phys. Rev. [48] J.-S.WangandB.Li,Phys.Rev.Lett.92,074302(2004)
Lett.94, 025507 (2005) [49] P. Cipriani, S. Denisov, and A. Politi, Phys. Rev. Lett.
6
94, 244301 (2005) Rev. Lett.84, 2144 (2000).
[50] E. Pereira and R. Falcao, Phys. Rev. Lett. 96, 100601 [62] O. V.Gendelman and A. V.Savin,Phys.Rev.Lett.84,
(2006) 2381 (2000).
[51] L. Delfini, S. Lepri, R. Livi, and A. Politi, Phys. Rev. E [63] L. Yangand B. Hu,Phys.Rev.Lett. 94, 219404 (2005).
73, 060201 (2006) [64] Y.Li,S.Liu,N.Li,P.Hanggi,andB.Li,arXiv:1407.1161
[52] G. Basile, C. Bernardin, and S. Olla, Phys. Rev. Lett. [65] S. G. Das and A. Dhar, arXiv:1411.5247v2
96, 204303 (2006) [66] H. Spohn,arXiv:1411.3907
[53] H.van Beijeren, Phys. Rev.Lett. 108, 180601 (2012) [67] N.Li,J.Ren,L.Wang,G.Zhang,P.Hanggi, andB.Li,
[54] C.B.MendlandH.Spohn,Phys.Rev.Lett.111,230601 Rev. Mod. Phys. 84, 1045 (2012).
(2013) [68] N.Li,P.Tong, andB.Li, Europhys.Lett.7549(2006).
[55] E. Pereira, R. Falcao, and H.C.F. Lemos, Phys. Rev. E [69] N. Li and B. Li, Europhys.Lett. 78, 34001 (2007).
87, 032158, (2013) [70] N. Li and B. Li, Phys.Rev.E 76, 011108 (2007).
[56] H.Spohn,J. Stat.Phys. 154, 1191 (2014) [71] N. Li and B. Li, J. Phys. Soc. Japan 78, 044001 (2009).
[57] S. Liu, P. H¨anggi, N. Li, J. Ren, and B. Li, Phys. Rev. [72] N.Li, B. Li, andS.Flach, Phys.Rev.Lett.105, 054102
Lett.112, 040601 (2014) (2010).
[58] C. W. Chang, D. Okawa, H. Garcia, A. Majumdar, and [73] N. Li and B. Li, AIPAdvances2, 041408 (2012).
A.Zettl, Phys. Rev.Lett. 101, 075903 (2008) [74] N. Li and B. Li, Phys.Rev.E 87, 042125 (2013).
[59] X. Xu, L. F. C. Pereira, Y. Wang, J. Wu, K. Zhang, X. [75] L. Yang, N. Li, and B. Li, Phys. Rev. E 90, 062122
Zhao, S. Bae, C. T. Bui, R. Xie, J. T. L. Thong, B. H. (2014).
Hong, K. P. Loh, D. Donadio, B. Li, and B. Ozyilmaz, [76] B. V.Chirikov, Phys. Rep.52, 263 (1979).
Nat.Commun. 5, 3689 (2014) [77] R. S. MacKay, J. D. Meiss, I. C. Percival, Phys. D 13,
[60] T.Meier, F.Menges, P.Nirmalraj, H.H¨olscher,H.Riel, 55 (1984).
and B. Gotsmann, Phys.Rev.Lett. 113, 060801 (2014)
[61] C. Giardina, R. Livi, A. Politi, and M. Vassalli, Phys.