# Monte Carlo Methods in Scientific Computing Heat Heat Conduction Conduction in in OneOneDimensional Dimensional Systems Systems:: molecular molecular dynamics dynamics and and mode-coupling mode-coupling theory theory Jian-Sheng Wang National University of Singapore 1 Outline

Brief review of 1D heat conduction Introducing a chain model Nonequilibrium molecular dynamics results Projection formulism and mode-coupling theory Conclusion 2 Fourier Law of Heat Conduction Fourier proposed the law of heat conduction in materials as J = T where J is heat current density, is thermal conductivity, and T is temperature. Fourier, Jean Baptiste Joseph, Baron (1768 1830)

3 Normal & Anomalous Heat Transport TL J TH 3D bulk systems obey Fourier law (insulating crystal: Peierls theory of Umklapp scattering process of phonons; gas: kinetic theory, = cvl ) In 1D systems, variety of results are obtained and still controversial. See S Lepri et al, Phys Rep 377 (2003) 1, for a review. 4

Heat Conduction in One-Dimensional Systems 1D harmonic chain, (Rieder, Lebowitz & Lieb, 1967) diverges if momentum is conserved (Prosen & Campbell, 2000) Fermi-Pasta-Ulam model, 2/5 (Lepri et al, 1998) Fluctuating hydrodynamics + Renormalization group, 1/3 (Narayan & Ramaswamy 2002) 5 Approaches to Heat Transport Equilibrium molecular dynamics using linear response theory (Green-Kubo formula) Nonequilibrium steady state (computer) experiment Laudauer formula in quantum regime

6 Ballistic Heat Transport at Low Temperature Laudauer formula for heat current 1 2 I n( ) | t ( ) | d 2 teik ' x eikx re ikx scatter 7 Carbon Nanotube Heat conductivity of Carbon nanotubes at T = 300K by nonequilibrium

molecular dynamics. From S Maruyama, Microscale Thermophysics Engineering, 7 (2003) 41. See also G Zhang and B Li, cond-mat/0403393. 8 Carbon Nanotubes Thermal conductance A of carbon nanotube of length L, determined from equilibrium molecular dynamics with GreenKubo formula, periodic boundary conditions, Tersoff potential. Z Yao, J-S Wang, B Li, and G-R Liu, cond-mat/0402616.

9 Fermi-Pasta-Ulam model A Hamiltonian system with pi 2 H ( p, x) V ( xi 1 xi ) i 1 2m 1 2 2 3 V ( z ) m ( z a ) ( z a) ( z a) 4 2 3 4 N

A strictly one-dimensional model. 10 A Chain Model for Heat Conduction ri = (xi,yi) TL m TH i pi2 1 2 H (p, r ) K r ri 1 ri a 2 i 2m K cos( i )

i Transverse degrees of freedom introduced 11 Nonequilibrium Molecular Dynamics Nos-Hoover thermostats at the ends at temperature TL and TH Compute steady-state heat current: j =(1/N)i d (i ri)/dt, where i is local energy associated with particle i Define thermal conductance by = (TH-TL)/(Na) N is number of particles, a is lattice spacing. 12 Nos-Hoover Dynamics fi L pi , if i N w dp i fi , if N N w i N w

dt fi H pi , if i N N w d L,H dt 1 1 2 1 k BTL , H N w 2 i p iw m 13

Defining Microscopic Heat Current Let the energy density be (r, t ) i (r ri ) i J 0 then J satisfies t A possible choice for total current is N d ( i ri ) Nj J (r , t )dV dt i 1 V 14 Expression of j for the chain model

mji ri (pi pi 1 ) G (i ) ri 1 (pi pi 1 ) G (i 1) ri 1 pi H (i 2, i 1, i 1) ri p i H (i 1, i 1, i ) pi i 1 G (i ) K r (| ri | a)n i 4 H (i, j , k ) K (ni n k cos j ) / | rk | 2 p 1 i K r (| ri 1 | a) 2 (| ri | a) 2 K cos(i ) i 4 2m ri ri 1 ri , ni ri / | ri | 15 Temperature Profile Temperature of i-th particle computed from kBTi= for parameter

set E with N =64 (plus), 256 (dash), 1024 (line). 16 Conductance vs Size N Model parameters (K, TL, TH): slope=1/ 3 Set F (1, 5, 7), B (1, 0.2, 0.4), E (0.3, 0.3, 0.5), H (0, 0.3, 0.5), J (0.05, 0.1, 0.2) , slope=2/5 ln N

m=1, a=2, Kr=1. From J-S Wang & B Li, Phys Rev Lett 92 (2004) 074302. 17 Additional MD data Parameters (K, TL, TH, ), set L(25,1,1.5,0.2) G(10,0.2,0.4,0) K(0.5,1.2,2,0.4) I(0.1,0.3,0.5,0.2) C(0.1,0.2,0.4,0) From J-S Wang and B Li, PRE, 70, 021204 (2004). 18

Mode-Coupling Theory for Heat Conduction Use Fourier components as basic variables Derive equations relating the correlation functions of the variables with the damping of the modes, and the damping of the modes to the square of the correlation functions Evoke Green-Kubo formula to relate correlation function with thermal conductivity 19 Basic Variables (work in Fourier space) m i 2 kj / N Q ( x ja

) e , j N j k m i 2 kj / N Q y e , j N j

T Pk Qk , A ( Pk , Pk , Qk , Qk ) k k 1, 2, , N 20 Equation of Motion for A A H H LA, L t q p p q Formal solution:

tL A(t , p, q) A( p (t ), q(t )) e A(0, p, q) 21 Projection Operator & Equation Define 1 PX X , A A, A A We have P 2 P Apply P and 1P to the equation of

motion, we get two coupled equations. Solving them, we get t dA(t ) iA(t ) (t s ) A( s )ds Rt dt 0 22 Projection Method (Zwanzig and Mori) Equation for dynamical correlation function: t G (t ) (t )G ( )d iG (t ) 0 where G(t) is correlation matrix of normalmode Canonical coordinates (Pk,Qk). is related to the correlation of random force.

23 Definitions G (t ) A(t ), A(0) 1 A, A 1 0 (t ) Rt , R A, A 1 i A , A

A, A PX X , A 1 A, A A, P 2 P Rt et (1 P ) L (1 P )LA dA(t ) LA(t ) dt L is Liouville operator 24 Correlation function equation and its solution (in FourierLaplace space)

t G (t ) (t )G ( )d iG (t ) 0 izt G[ z ] e G (t )dt Define 0 the equation can be solved as G[ z ] 1 i ( z ) [ z ] in particular g k [ z ]

0 Qk (t )Qk (0)* | Qk |2 iz [ z] e izt dt 2 2 k , , (k ) z izk [ z ] 25 Small Oscillation Effective Hamiltonian 1 2 2

H eff ( P, Q) Pk (k ) Qk 2 k , k p q 0 2 (3) v Q Q

Q v Q k , p,q k p q k , p,q k QpQq Equations of motion P H , Q H k k Qk Pk 26 Equation of Motion of Modes 2

Qk (k ) Qk .. Q Q .. Q Q k' k '' k' k' k ' k ''k ( ) 2 Q Q k

k k .. Qk' Qk'' k ' k ''k 27 Determine Effective Hamiltonian Model Parameters from MD 2 k |Q | vk , p ,q vk(3),p ,q

1 1 2 O(v ), 2 (k ) k BT QkQ pQq 2 | Qk |2 | Q p |2 | Qq |2 QkQp Qq 6 | Qk |2 | Qp |2 | Qq |2 28 Mode-Coupling Approximation (t) RQQ (t) g(t)g(t) [mean-field type]

29 Full Mode-Coupling Equations k (t ) K p ,q g p (t ) g q (t ) K p,q g p (t ) g q (t ), p q k k (t ) K p(3),q g p (t ) g q (t ),

p q k p q k iz k [ z] gk [ z] 2 , , 2 z (k ) iz k [ z ] g k [ z ] is Fourier-Laplace transform of g k (t ) 30 Damping Function [z]z] Molecular Dynamics

Mode-Coupling Theory From J-S Wang & B Li, PRE 70, 021204 (2004). 31 Correlation Functions Correlation function g(t) for the slowest longitudinal and transverse modes. Black line: mode-coupling, red dash: MD. N = 256. g(t) e-tcos(t) 32 Decay or Damping Rate

Decay rate of the mode vs mode index k. p = 2k/(Na) is lattice momentum. N = 1024. slope=3/2 longitudinal slope=2 transverse Symbols are from MD, lines from mode-coupling theory. Straight lines have slopes 3/2 and 2, respectively. 33

Mode-Coupling Theory in the Continuum Limit 1 (t ) 2 1 (t ) 2 / a 2 (3) 2 dq K g

( t ) K g ( t ) q q /a / a dq K

g ( t ) g q q (t ) /a p, (t ) p 2 v, (t ) 34 Asymptotic Solution The mode-coupling equations predict, for large system size N, and small z : 2 k 2 p [ z ] c p , p Na 1

2 p [ z ] bz p , = 2 If there is no transverse coupling, = z(-1/3)p2 (Result of Lepri). 35 Mode-Coupling [z]z]/p2 || slope = 1/2 slope = 0 At parameter set B. Blue dash : asymptotic analytical result, red line : Full theory on N =1024, solid line : N limit theory

36 Green-Kubo Formula 1 k BT 2 aN J (t ) J (0) dt , 0 k J bk Qk Pk * , bk ik p k , J (t ) J (0) |bk |2 g QQ ,k (t ) g PP ,k (t ) g QP ,k (t ) 2

k , k , 2 b 2 1/(2 ) g ( t ) t k k k

37 Green-Kubo Integrand Parameter set B. Red circle: molecular dynamics, solid line: modecoupling theory (N = 1024), blue line: asymptotic slope of 2/3. 38 N with Periodic Boundary Condition slope=1/ 2 Molecular dynamics Mode-coupling

from GreenKubo formula on finite systems with periodic boundary conditions, for parameter set B (Kr=1, K=1, T=0.3) 39 Relation between Exponent in and If mode decay with z-p2, then With periodic B.C. thermal conductance N 1- With open B.C. N 1-1/(2-) Mode coupling theory gives =1/2 with transverse motion, and =1/3 for strictly 1D system. 40

Conclusion Quantitative agreement between modecoupling theory and molecular dynamics is achieved Molecular dynamics and mode-coupling theory support 1/3 power-law divergence for thermal conduction in 1D models with transverse motion, 2/5 law if there are no transverse degrees of freedom. 41