Generalized Two-Temperature Model for Coupled ... - ACS Publications

Aug 4, 2017 - number of graphene sheets or flakes and the orientation.4,8,9 To enable rational design of graphene-based composites ..... to integrate ...
0 downloads 9 Views 1MB Size
Subscriber access provided by BOSTON UNIV

Communication

Generalized two-temperature model for coupled phonons in nanosized graphene Meng An, Qichen Song, Xiaoxiang Yu, Han Meng, Dengke Ma, Ruiyang Li, Zelin Jin, Baoling Huang, and Nuo Yang Nano Lett., Just Accepted Manuscript • DOI: 10.1021/acs.nanolett.7b02926 • Publication Date (Web): 04 Aug 2017 Downloaded from http://pubs.acs.org on August 5, 2017

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a free service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are accessible to all readers and citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

Nano Letters is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 27

Nano Letters

Page 1 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

1

Generalized two-temperature model for coupled phonons in nanosized

2

graphene

3

Meng An1,2,#, Qichen Song1,2,#,†, Xiaoxiang Yu1,2, Han Meng1,2, Dengke Ma1,2, Ruiyang Li1,2,

4

Zelin Jin1,2, Baoling Huang3,4, and Nuo Yang1,2*

5

1

6

Wuhan 430074, People's Republic of China

7

2

8

Huazhong University of Science and Technology (HUST), Wuhan 430074, People's Republic

9

of China

State Key Laboratory of Coal Combustion, Huazhong University of Science and Technology,

Nano Interface Center for Energy (NICE), School of Energy and Power Engineering,

10

3

11

Science and Technology, Clear Water Bay, Kowloon, Hong Kong

12

4

13

Shenzhen, 518057, China

Department of Mechanical and Aerospace Engineering, The Hong Kong University of

The Hong Kong University of Science and Technology Shenzhen Research Institute,

14 15 16 17 18 19

# M.A. and Q.S. contributed equally to this work.

20

† Current address: Department of Mechanical Engineering, Massachusetts Institute of Technology, 77

21

Massachusetts Avenue, Cambridge, MA 02139, USA.

22

* To whom correspondence should be addressed. E-mail: (N.Y.) [email protected]

ACS Paragon Plus Environment

Nano Letters

Page 2 of 27

Page 2 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

1

Abstract

2

The design of graphene-based composite with high thermal conductivity requires a

3

comprehensive understanding of phonon coupling in nanosized graphene. We extended the

4

two-temperature model to coupled groups of phonon. The study give new physical quantities,

5

the phonon-phonon coupling factor and length, to characterize the couplings quantitatively.

6

Besides, our proposed coupling length has an obvious dependence on system size. Our

7

studies can not only observe the nonequilibrium between different groups of phonon, but

8

explain theoretically the thermal resistance inside nanosized graphene.

9

KEYWORDS: Phonons, two temperature model, coupling strength, nanosized graphene,

10

molecular dynamics simulations

11

ACS Paragon Plus Environment

Page 3 of 27

Nano Letters

Page 3 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

1

Introduction

2

Owing to the superior thermal conductivity 3000-5000 W/m-K of graphene1, the graphene-

3

based composites2-5 have been widely applied in thermal managements of advanced

4

electronics6, optoelectronics, the photovoltaic solar cell7 and Li-ion battery4. For example, it

5

is found that the thermal conductivity of graphene-based composites at relatively low filler

6

contents can reach up 5.1 W/m-K2. Besides, it has been demonstrated that thermal

7

conductivity of graphene-based composite can be further enhanced by changing the layer

8

number of graphene sheets or flakes and the orientation4, 8, 9. To enable rational design of

9

graphene-based composites with higher thermal conductivity, a fundamental and

10

comprehensive understanding of thermal transport in such materials is essential.

11 12

Extensive studies have concluded that two main thermal resistances play central roles in

13

determining the thermal conductivity of graphene-based composites. The first one is the high

14

interfacial thermal resistance between graphene and matrix materials2, 10. Recent investigation

15

showed that the interfacial thermal resistance could be overcome by alignment arrangement

16

of graphene in composites6. The second one is the phonon coupling thermal resistance

17

between low-frequency out-of-plane (OP) phonons and high-frequency in-plane (IP)

18

phonons11, which is not well studied and responsible for the thermal resistance inside

19

graphene. The heat energy is mostly transferred from matrix materials to OP phonons of

20

graphene due to the strong coupling between the matrix materials and the low-frequency OP

21

phonons (shown in Fig. S1). In graphene-based composite samples, the size of graphene

22

flakes/particles

23

nanoflakes/particles, the OP phonons contribute less to the thermal conductivity than IP

24

phonons because the OP phonons have longer mean free path and wave length than IP

is

nanoscale2,

6

.

The

studies12,

13

indicated

ACS Paragon Plus Environment

that,

in

graphene

Nano Letters

Page 4 of 27

Page 4 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

1

phonons and are suppressed more by the size confinement. So, there is an indirect heat

2

transport inside graphene nanoflakes/particles: the heat will be transferred through the

3

coupling between OP and in-plane (IP) phonons. However, the coupling between OP and IP

4

is relatively weak14, 15. As a result, that leads to a large thermal resistance inside graphene. In

5

this work, we mainly focus on the phonon mode coupling resistance between IP and OP,

6

which determines the heat flux inside nanosized graphene.

7 8

The two-temperature model (TTM) has been successfully applied to investigate the coupling

9

between different energy carriers in thermal nonequilibrium16, 17. In TTM, different energy

10

carriers are considered as two interacting subsystems, such as the electron-phonon18, 19 and

11

phonon-magnon16, 20. Liao et al. have generalized TTM for the coupled phonon-magnon

12

diffusion and predicted a magnon cooling effect16. Besides, some reports show that TTM is

13

implemented with molecular dynamics which is used to study electron-phonon coupled in

14

metal/semiconductor systems19, 21, 22.

15 16

Here, we extended TTM to investigate the coupling between different phonon groups. This

17

method is demonstrated on nanosized graphene, a representative two-dimensional material

18

with the highest known thermal conductivity23,

19

graphene into two groups, namely in-plane (IP) and out-of-plane (OP) phonon group. Based

20

on the temperature profiles of different phonon groups, we calculated the ph-ph coupling

21

factor Gio, comparable with strength of electron-phonon coupling19 and found that Gio is not

22

sensitive to system size. We successfully demonstrated the poor coupling between IP and OP

23

phonon groups in nanosized graphene25-28.

24

Theoretical analysis

24

. We separate phonons in nanosized

ACS Paragon Plus Environment

Page 5 of 27

Nano Letters

Page 5 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

1 2

In this letter, TTM is used to investigate the ph-ph coupling between OP and IP phonon

3

groups. As shown in the configuration (Fig. 1), the entire system is divided into two regions,

4

namely, the left region (l) (where atoms have both OP and IP vibrations) and the right region

5

(r) (only IP vibration). TH (TC) is the prescribed temperature of hot (cold) bath at the left

6

(right) end. As observed in the left region, the temperature difference between OP and IP

7

phonon groups increases and reaches a maximum value, θmax, at the interface. This

8

phenomenon results from the poor coupling between OP and IP phonon groups, which is a

9

necessary condition for TTM.

10 11

The Boltzmann transport equation (BTE) for IP and OP phonons in graphene that determines

12

the phonon distribution function can be written separately as:

13

∂f i  ∂f  + F ⋅ ∇p f i + v ⋅ ∇ r f i =  i  ∂t  ∂t  c

(1a)

14

∂f o  ∂f  + F ⋅ ∇p fo + v ⋅ ∇r f o =  o  ∂t  ∂t c

(1b)

15

where the subscripts o and i denote OP and IP phonons, respectively. When the system

16

reaches steady state without applying external force, the Eq. (1a) and Eq. (1b) can be

17

simplified as:

18

 ∂f  v ⋅ ∇r fi =  i   ∂t  c

(2a)

19

 ∂f  v ⋅ ∇r f o =  o   ∂t c

(2b)

ACS Paragon Plus Environment

Nano Letters

Page 6 of 27

Page 6 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

1

where the right-hand sides are the collision term, which includes three kinds of phonon

2

scattering processes: IP-IP scatterings, OP-OP scatterings and scatterings between IP and OP

3

phonons.

4 5

Previous studies29, 30 have demonstrated that the reflection symmetry in graphene forbids odd

6

number of OP phonons (ZA and ZO) in a three-phonon scattering process, which leads to the

7

relatively weak coupling between OP and IP phonons. Moreover, we also observed the weak

8

couplings by MD simulation results (shown in SI and Fig. S4). The weak coupling has a

9

negligible effect on the distribution and phase space of phonons. IP and OP phonons will be

10

treated as two subsystems. Then, the collision term can be written as:

11

 ∂f   ∂f  v ⋅ ∇r fi =  i  +  i   ∂t io  ∂t ii

(3a)

12

 ∂f   ∂f  v ⋅ ∇r fo =  o  +  o   ∂t io  ∂t oo

(3b)

13

where the subscript ii, oo, io denote IP-IP scatterings, OP-OP scatterings and IP-OP

14

scatterings, respectively. When the relaxation time approximation (RTA) is adopted for IP-IP

15

scatterings, OP-OP scatterings, the Eq. (3a) and Eq. (3b) can write:

16

 ∂f   f - f  v ⋅ ∇r fi =  i  +  i i,0   ∂t io  τ ii ii

(4a)

17

 ∂f   f - f  v ⋅ ∇r f o =  o  +  o o,0   ∂t io  τ oo oo

(4b)

ACS Paragon Plus Environment

Page 7 of 27

Nano Letters

Page 7 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

1

where f i , 0 and f o, 0 is the IP and OP phonon distribution function at equilibrium (the Bose-

2

Einstein distribution f = ( exp ( hω k B T ) − 1) , respectively.

3

time for IP and OP phonons, respectively. After multiplying Eq. (4a) and Eq. (4b) by ħω, and

4

integrating over all wavevector q , the last terms on the right-hand sides of Eq. (4a) and Eq.

5

(4b) drop out because they are odd function with respect to all wavevector q

6

described as:

-1

τ ii and τoo are

the relaxation

31

. And they are

7

 ∂f  ∇r ⋅ ∑ hωq vfi = ∑ hωq  i   ∂t io q q

(5a)

8

 ∂f  ∇r ⋅ ∑ hωq vfo = ∑ hωq  o   ∂t io q q

(5b)

9

Denoting J i / o = ∑ hωq vf i / o and ∂Ei / o ∂t = ∑ hωq ( ∂f i / o ∂t )io , they write: q

q

10

∇ ⋅ Ji = −

∂Ei ∂t

(6a)

11

∇ ⋅ Jo = −

∂Eo ∂t

(6b)

12

The scatterings between IP-OP are responsible for the local energy exchange between IP and

13

OP. We use

14

the rate of energy transfer between IP and OP is then described by32:

1 τ io to describe the scattering rate for the coupling between IP and OP. Then

15

∆Ei ci (Ti nter − Ti ) = ∆t τ io

(7a)

16

∆Eo co (Tinter − To ) = ∆t τ io

(7b)

ACS Paragon Plus Environment

Nano Letters

Page 8 of 27

Page 8 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

1

where the intermediate temperature that both channels approach is Tinter =

2

Therefore, the rate of energy transfer is written as:

ciTi + coTo . ci + co

3

∆Ei ci co (To − Ti ) = = Gio (To − Ti ) τ io (ci + co ) ∆t

(8a)

4

∆Eo ci co (Ti − To ) = = Gio (Ti − To ) ∆t τ io (ci + co )

(8b)

Gio =

5

ci co τ io (ci + co )

(8c)

6

where Gio, is defined as the coupling factor. Then, with the heat diffusion equation, the Eq.

7

(7a) and Eq. (7b) can be written as:

8

∇qi = Gio (To − Ti )

(9a)

9

∇qo = Gio (Ti − To )

(9b)

10

We then apply the Fourier Law q = −κ∇ T to rewrite the Eq. (9a) and Eq. (9b):

11

−κi∇2Ti = Gio (To − Ti )

(10a)

12

−κo∇2To = Gio (Ti − To )

(10b)

13

Considering the one-dimensional temperature gradient in graphene, the two-dimensional heat

14

transfer problem through IP, OP and the coupling between OP and IP can be well described

15

by a one-dimensional (1D) model. The governing equations for coupled phonon transport

16

states:

17

d 2To κo 2 = Gio (To − Ti ) dx

ACS Paragon Plus Environment

(11a)

Page 9 of 27

Nano Letters

Page 9 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

1

- κi

d 2Ti = Gio (To − Ti ) dx 2

(11b)

2

Subtracting Eq. (11a) from Eq. (11b), we can obtain d 2θ dx2 − γ 2θ = 0 , where θ = To − Ti and

3

γ = Gio (1 κ i + 1 κ o ) . Based on the boundary condition θ x→−∞ = 0 and θ x=0 = θ max , the

4

temperature difference profile is written as θ = θ max exp(γx ) .

5 6

The validity of TTM requires that the coupling between two carriers is so weak that the

7

coupling has a negligible effect on the distribution and phase space of each carrier. As far as

8

we know, TTM has been successfully used to describe the coupling of electron-phonon18, 19,

9

and phonon-magnon16, 20. For the coupling between IP phonons and OP phonons, both them

10

are bosons, which is similar the phonon-magnon coupling. Both the experimental and

11

simulation work14, 15 have implied the coupling of IP-OP is much weaker than that of either

12

IP-IP or OP-OP. In addition, this weak coupling is confirmed by our MD results (shown in

13

Fig. S2). It is deemed that the IP-OP coupling has a negligible effect on the distribution and

14

phase space of either IP or OP phonon group. Therefore, the IP-OP coupling satisfies the

15

requirement of TTM. We have computed more simulation results (Fig. S6) to show the strong

16

coupling between two in-plane directions (x and y). Because TTM is only applicable to the

17

weak coupling between different energy carriers, the coupling factor between two IP

18

directions should not be resolved using TTM.

19 20

Here, we propose the ph-ph coupling length, lc, shown in Fig. 1, to quantitatively characterize

21

the coupling between them. Its physical meaning is the distance required to equilibrate OP

ACS Paragon Plus Environment

Nano Letters

Page 10 of 27

Page 10 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

1

and IP phonon groups when the temperature difference exists. Specifically, we define such a

2

characteristic length as the distance between the position of θmax and 5%θmax: lc =

3

− ln(5%)

γ



3

γ

=

3 1 1  G io  +   κi κo 

(12)

4 5

Numerical Simulations

6

The ph-ph coupling length of nanosized graphene is numerical calculated by means of TTM-

7

MD, and the detailed non-equilibrium MD simulation setup is shown in Fig. 2(a). The heat

8

source with a higher temperature TH=310K is applied to the atoms in red region and the heat

9

sink with a lower temperature TC=290K is applied to atoms in blue region where only IP

10

phonons exist. All of our simulations are performed by large-scale atomic/molecular

11

massively parallel simulator (LAMMPS) packages33. The fixed (periodic) boundary

12

conditions are used along the x (y) direction. The optimized Tersoff potential34 is applied to

13

describe interatomic interactions, which has successfully reproduced the thermal transport

14

properties of graphene35, 36. The detailed parameters of optimized Tersoff potential are shown

15

in Supporting Information (SI). The velocity Verlet algorithm is adopted to integrate the

16

discrete differential equations of motions. The time step is set as 0.5fs. We relaxed the

17

graphene structure in the isothermal-isobaric (NPT) ensemble. After NPT relaxation,

18

simulations are performed for 2ns to reach a steady state. After that, a time average of the

19

temperature and heat current is performed for 10ns.

20

21

The thermal conductivity is calculated based on the Fourier’s Law37 κ = − J A ⋅∇T , where J

22

denotes the heat current transported from the hot bath to cold bath, A is the cross section area

ACS Paragon Plus Environment

Page 11 of 27

Nano Letters

Page 11 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

1

and ∇T is the temperature gradient along the x direction. The temperature of different

2

phonon groups is computed by Tp ( x ) =

N

∑m v i

i, p

vi , p

NkB , where vi , p is the velocity

i =1

3

vector of phonon vibrating along p direction, i.e. IP or OP, of atom i. N and kB are the number

4

of atoms in the simulation cell and the Boltzmann constant, respectively. Besides, it is

5

calculated that the thermal conductivity of IP (OP) including only the atomic vibration along

6

in-plane (out-of-plane) direction by freezing the other direction. The simulations detail is

7

shown in Supplementary Information (Fig. S3).

8 9

Results and Discussion

10

In our simulations, the width of simulation cell is set as 5.2 nm. The lattice constant (a) and

11

thickness (d) of graphene are 0.143 nm and 0.335 nm, respectively. The size effect could

12

arise if the width is not sufficiently large38, 39. The thermal conductivity reaches a saturated

13

value when the width is larger than 5.2 nm which is consistent with our previous works39, 40.

14

Besides, to examine the accuracy of MD results we calculated thermal conductivity (κ) of

15

nanosized suspended single-layer graphene with a length (L) as 33 nm. The calculated κ is

16

1082 ± 103 W/m-K, which is consistent with previous studies35, 36.

17 18

The temperature profile of TTM-MD is presented in Fig. 2(b), where the system length along

19

x direction is 101 nm. In the left region, there are obvious differences between the

20

temperature curves of IP and OP phonons. It is found that the maximum temperature

21

difference (θmax) is 5.1 K at the interface. In order to extract the ph-ph coupling length (lc),

22

the profile of θ in Fig. 2(c) is fitted exponentially, based on theoretical solution of Eq. (11a

23

and 11b). Then, lc can be obtained based on its definition as the distance between the position

ACS Paragon Plus Environment

Nano Letters

Page 12 of 27

Page 12 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

1

of θmax and 5%θmax. If the system length is not large enough to equilibrate different phonon

2

groups, the coupling length would be larger than system length. That is IP and OP

3

subsystems are under nonequilibrium state. We calculated four simulation cells with different

4

length, and lc are 60 nm, 64 nm, 70 nm and 77 nm corresponding to length is 33 nm, 60 nm,

5

67 nm, and 101 nm, respectively. The profiles of temperature difference with L=33nm, 67nm,

6

70nm is shown in Supplementary Information (Fig. S4).

7 8

To evaluate the ph-ph coupling factor (Gio) based on Eq. (12), we calculated the thermal

9

conductivity of IP/OP phonon groups (κi/κo) including only the atomic vibration along in-

10

plane (out-of-plane) direction (shown in Fig. 3). It is also shown that the summation of κi and

11

κo approximately equals to the thermal conductivity of pristine graphene (blue open circles),

12

instead of an obvious reduction by scatterings. This phenomenon in nanosized graphene

13

arises from the poor ph-ph coupling between OP and IP26, 28, 41 and the decoupled phonon

14

scattering interaction between OP and IP of freezing method28. Moreover, it is noted that our

15

MD results on graphene with nano-confinement shows (in Fig. 3) the relative contribution of

16

OP phonons is around 25%, which agrees with the results from other groups by theoretical

17

prediction42 and MD results13, 43, 44. On the other side, Lindsay et al.30 found that, in an

18

infinite graphene sheet, the OP phonons in graphene dominate the thermal transport because

19

the reflection symmetry in two-dimensional graphene significantly restricts the phase space

20

for phonon-phonon scattering of OP phonons45. This discrepancy may come from three

21

possible reasons. The first is the presence of normal scatterings which is excluded in the

22

iterative three-phonon calculation. The second reason may come from the essential difference

23

between three-phonon calculation and MD simulation. In anharmonic lattice dynamics

24

calculation, the graphene sheet remains a perfect plane with the reflection symmetry perfectly

25

preserved, while in MD simulations the atoms are not at their equilibrium positions and the

ACS Paragon Plus Environment

Page 13 of 27

Nano Letters

Page 13 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

1

graphene is not a plane. As a result, the reflection symmetry may not be well presented in

2

MD simulations. Thirdly, another possible reason is the size effect. That is, Lindsay et.al

3

studied bulk graphene45, while our MD results and other groups calculated the thermal

4

conductivity of finite-sized nanostructured graphene.

5 6

To show that our results are robust, we also have used the AIREBO potential46 to calculate

7

the thermal conductivity of different phonon groups in nanosized graphene (as shown in Fig.

8

S5). It is found that the relative contribution of OP phonons is smaller than IP phonons for

9

both two potentials.

10 11

In low dimensional nanostructures, the thermal conductivity becomes dependent on the

12

characteristic length of system due to the length is comparable to phonon mean free path

13

(MFP)47. Recent studies9, 23, 48-50 showed that the thermal conductivity of graphene increases

14

logarithmically (~logL) with the size of sample, which is taken to fit MD results shown in Fig.

15

3.

16 17

With the κi, κo and lc for several systems with different lengths, we can calculate the ph-ph

18

coupling factor Gio based on Eq. (12) (shown in Fig. 4(a)). To the best of our knowledge, it is

19

the first time to extract the ph-ph coupling factor. It is found that Gio in nanosized graphene is

20

not sensitive to system length, and the average value is 4.26×1017 W/m3-K, which is

21

comparable with the electron-phonon coupling factor (Gep) in metals ranging from

22

5.5×1016W/m3-K to 2.6×1017W/m3-K19. To verify our results, the coupling factor between IP

23

and OP phonon groups is also calculated from anharmonic lattice dynamics and first-

ACS Paragon Plus Environment

Nano Letters

Page 14 of 27

Page 14 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

1

principles (calculation details in SI). The Gio is calculated as 1.07 × 1017 W/m 3 -K , which is the

2

same order of magnitude as the value from MD simulation, 4.26 × 1017 W/m 3 -K .The

3

discrepancy may come from the fact that the reflection symmetry in graphene forbids odd

4

number of OP modes (ZA and ZO) in a three-phonon scattering process29, 30. As a result, the

5

interaction between IP and OP phonon can be weak. In first-principle calculation, this

6

reflection is forced to be preserved, yet in MD simulation and in reality, the existence of

7

boundary can break the reflection symmetry and induces more scattering that is involved with

8

both IP and OP phonons. Besides, the first-principles calculation used here fails to capture

9

higher-order phonon scattering, which is naturally embedded in MD simulation. This may

10

also lead to the underestimated coupling factor from first-principle.

11 12

As shown in Fig. 4(b), we obtained the size dependence of ph-ph coupling length of

13

nanosized graphene based on the size dependence of the thermal conductivity. It is clearly

14

observed that the ph-ph coupling length, lc, diverges when κ ~ logL23. Combining with the

15

definition of lc in Eq. (12), it is easy to understand that the size dependence of thermal

16

conductivity contributes to the size effect of lc.

17 18

The ph-ph coupling factor and length may provide us a new way in understanding heat

19

transfer. For example, the weak coupling means there are less scatterings between IP and OP

20

phonons, which attribute to the ultra-high thermal conductivity of suspended graphene sheet.

21

The coupling factor and length provide a perspective in explanation quantitatively. As we

22

know, the thermal performance of graphene-based composites is determined from the

23

interfacial thermal resistance (named as outer resistance)51, 52. Besides, based on our finding

24

above, the longer ph-ph coupling length means it is difficult to reach the equilibrium state

ACS Paragon Plus Environment

Page 15 of 27

Nano Letters

Page 15 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

1

between IP and OP phonons, which corresponds to a larger inside thermal resistance. Both

2

the outer and inside resistances contribute to the unsatisfying thermal performance of

3

graphene-based composites.

4

Conclusion

5

In summary, we have extended the TTM to investigate the coupled phonons implemented

6

with molecular dynamics. The merit of this work lies in the fact that we applied the

7

techniques widely used in the field of electron-phonon interactions to investigate the ph-ph

8

interactions, utilized the analogy between field-driven electron and phonon, and combined IP

9

and OP phonon groups in a smooth way. The magnitude of ph-ph coupling factor has been

10

estimated in nanosized graphene and is independent of system size. Besides, we proposed the

11

ph-ph coupling length associated with the ph-ph coupling strength. Our studies can not only

12

observe the nonequilibrium between different groups of phonon, but explain theoretically the

13

thermal resistance inside nanosized graphene.

14 15

Acknowledgements

16

This project is supported by the National Natural Science Foundation of China (Grant

17

51576076). We are grateful to Jingtao Lü, Gang Zhang, Lifa Zhang, Jinlong Ma and Shiqian

18

Hu for useful discussions. The authors thank the National Supercomputing Center in Tianjin

19

(NSCC-TJ) and China Scientific Computing Grid (ScGrid) for providing assistance in

20

computations.

21

ACS Paragon Plus Environment

Nano Letters

Page 16 of 27

Page 16 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

1

Figure 1 (Color on-line) Representative temperature profile in TTM for the coupled phonons.

2

The system consists of two regions, namely, the left region (l) (where atoms have both OP

3

and IP vibrations) and the right region (r) (only IP vibration). To,l, Ti,l are temperatures for OP

4

and IP phonon group in the left region. Ti,r is the temperature of IP phonon group in the right

5

region. θmax is the maximum temperature difference between OP and IP phonon group. lc is

6

the ph-ph coupling length, which is defined as the distance between the position of θmax and

7

5%θmax.

8

ACS Paragon Plus Environment

Page 17 of 27

Nano Letters

Page 17 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

1

Figure 2 (Color on-line) (a) Schematic illustration of molecular dynamics simulation setup.

2

The black (green) arrow denotes the vibration along the out-of-plane (in-plane) direction. The

3

blue dot line is the interface of different simulation regions. The atoms in the left region

4

vibrate freely in three directions while those in the right region could only vibrate along in-

5

plane direction. The fixed (periodic) boundary condition is applied to along x (y) direction.

6

Hot bath (red atoms) TH and cold bath (blue atoms) TC are applied. The temperatures of two

7

heat baths are set as TH=310K, TC=290K, respectively. (b) The temperature profile of

8

different phonon groups. (c) The temperature difference distribution between OP and IP

9

phonon group when the system length is set as L=101nm. The fitting line is based on

10

θ = θ max exp ( γ x ) . The red dot lines denote the maximum temperature difference and its five

11

percent, respectively.

12

ACS Paragon Plus Environment

Nano Letters

Page 18 of 27

Page 18 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

1

Figure 3. (Color on-line) The length dependence of thermal conductivity with different

2

phonon groups in nanosized graphene calculated from NEMD. Pristine graphene (Pristine),

3

atoms vibrated along two in-plane directions only (Only IP), atoms vibrated along our-of-

4

plane direction only (only OP), and the summation of the values of only IP and only OP. The

5

system length ranges from 16 nm to 100 nm. The fitting lines are based on κ ~ logL.

6

ACS Paragon Plus Environment

Page 19 of 27

Nano Letters

Page 19 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

1

Figure 4. (Color on-line) (a) The ph-ph coupling factor versus the system lengths calculated

2

by MD and ab initio methods. (b) The length dependence of the ph-ph coupling length from

3

MD simulation, where the thermal conductivity of nanosized graphene follows κ~logL.

4

ACS Paragon Plus Environment

Nano Letters

Page 20 of 27

Page 20 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

References

1 2

1.

Balandin, A. A.; Ghosh, S.; Bao, W.; Calizo, I.; Teweldebrhan, D.; Miao, F.; Lau, C.

3

N. Nano Lett. 2008, 8, (3), 902-907.

4

2.

Shahil, K. M.; Balandin, A. A. Nano Lett. 2012, 12, (2), 861-867.

5

3.

Renteria, J. D.; Nika, D. L.; Balandin, A. A. Appl. Sci. 2014, 4, (4), 525-547.

6

4.

Goli, P.; Legedza, S.; Dhar, A.; Salgado, R.; Renteria, J.; Balandin, A. A. J. Power

7

Sources 2014, 248, 37-43.

8

5.

Nika, D. L.; Balandin, A. A. Rep. Prog. Phys. 2017, 80, (3), 036502.

9

6.

Renteria, J.; Legedza, S.; Salgado, R.; Balandin, M.; Ramirez, S.; Saadah, M.; Kargar,

10

F.; Balandin, A. Materials & Design 2015, 88, 214-221.

11

7.

12

2854-2877.

13

8.

Goyal, V.; Balandin, A. A. Appl. Phys. Lett. 2012, 100, (7), 073113.

14

9.

Balandin, A. A. Nat. Mater. 2011, 10, (8), 569-81.

15

10.

Shahil, K. M. F.; Balandin, A. A. Solid State Commun. 2012, 152, (15), 1331-1340.

16

11.

Wu, X.; Luo, T. J. Appl. Phys. 2014, 115, (1), 014901.

17

12.

Wei, Z.; Yang, J.; Bi, K.; Chen, Y. J. Appl. Phys. 2014, 116, (15), 153503.

18

13.

Qiu, B.; Ruan, X. Appl. Phys. Lett. 2012, 100, (19), 193101.

19

14.

Sullivan, S.; Vallabhaneni, A.; Kholmanov, I.; Ruan, X.; Murthy, J.; Shi, L. Nano

20

Lett. 2017, 17, (3), 2049-2056.

21

15.

22

(12), 125432.

23

16.

Liao, B.; Zhou, J.; Chen, G. Phys. Rev. Lett. 2014, 113, (2), 025902.

24

17.

Majumdar, A.; Reddy, P. Appl. Phys. Lett. 2004, 84, (23), 4768.

25

18.

Chen, G. J. Appl. Phys. 2005, 97, (8), 083707.

Notton, G.; Cristofari, C.; Mattei, M.; Poggi, P. Appl. Therm. Eng. 2005, 25, (17),

Vallabhaneni, A. K.; Singh, D.; Bao, H.; Murthy, J.; Ruan, X. Phys. Rev. B 2016, 93,

ACS Paragon Plus Environment

Page 21 of 27

Nano Letters

Page 21 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

1

19.

Wang, Y.; Ruan, X.; Roy, A. K. Phys. Rev. B 2012, 85, (20), 205311.

2

20.

An, K.; Olsson, K. S.; Weathers, A.; Sullivan, S.; Chen, X.; Li, X.; Marshall, L. G.;

3

Ma, X.; Klimovich, N.; Zhou, J. Phys. Rev. Lett. 2016, 117, (10), 107202.

4

21.

Lu, Z. X.; Wang, Y.; Ruan, X. L. Phys. Rev. B 2016, 93, (6), 064302.

5

22.

Yang, N.; Luo, T.; Esfarjani, K.; Henry, A.; Tian, Z.; Shiomi, J.; Chalopin, Y.; Li, B.;

6

Chen, G. J. Comput. Theor. Nanos. 2015, 12, (2), 168-174.

7

23.

8

Xie, R.; Thong, J. T.; Hong, B. H.; Loh, K. P.; Donadio, D.; Li, B.; Ozyilmaz, B. Nat.

9

Commun. 2014, 5, 3689.

Xu, X.; Pereira, L. F.; Wang, Y.; Wu, J.; Zhang, K.; Zhao, X.; Bae, S.; Tinh Bui, C.;

10

24.

Balandin, A. A.; Ghosh, S.; Bao, W.; Calizo, I.; Teweldebrhan, D.; Miao, F.; Lau, C.

11

N. Nano Lett. 2008, 8, (3), 902-7.

12

25.

Zhang, J.; Wang, Y.; Wang, X. Nanoscale 2013, 5, (23), 11598-603.

13

26.

Zhang, J.; Wang, X.; Xie, H. Phys. Lett. A 2013, 377, (41), 2970-2978.

14

27.

Zhang, J.; Huang, X.; Yue, Y.; Wang, J.; Wang, X. Phys. Rev. B 2011, 84, (23),

15

235416.

16

28.

Zhang, H.; Lee, G.; Cho, K. Phys. Rev. B 2011, 84, (11), 115460.

17

29.

Gu, X.; Wei, Y.; Yin, X.; Li, B.; Yang, R. arXiv preprint arXiv:1705.06156 2017.

18

30.

Lindsay, L.; Broido, D. A.; Mingo, N. Phys. Rev. B 2010, 82, (11), 115427.

19

31.

Chen, G., Nanoscale energy transport and conversion: a parallel treatment of

20

electrons, molecules, phonons, and photons. Oxford University Press: New York, 2005.

21

32.

Sanders, D.; Walton, D. Phys. Rev. B 1977, 15, (3), 1489.

22

33.

Plimpton, S. J. Comput. Phys. 1995, 117, (1), 1-19.

23

34.

Lindsay, L.; Broido, D. A. Phys. Rev. B 2010, 81, (20), 205441.

24

35.

Yang, L.; Chen, J.; Yang, N.; Li, B. Int. J. Heat Mass Transfer 2015, 91, 428-432.

25

36.

Chen, J.; Zhang, G.; Li, B. Nanoscale 2013, 5, (2), 532-6.

ACS Paragon Plus Environment

Nano Letters

Page 22 of 27

Page 22 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

1

37.

Yang, N.; Zhang, G.; Li, B. Nano Lett. 2008, 8, (1), 276-280.

2

38.

Schelling, P. K.; Phillpot, S. R.; Keblinski, P. Phys. Rev. B 2002, 65, (14), 144306.

3

39.

Shiqian, H.; Meng, A.; Nuo, Y.; Baowen, L. Nanotechnology 2016, 27, (26), 265702.

4

40.

Song, Q.; An, M.; Chen, X.; Peng, Z.; Zang, J.; Yang, N. Nanoscale 2016, 8, (32),

5

14943-9.

6

41.

Zhang, J.; Wang, X.; Xie, H. Phys. Lett. A 2013, 377, (9), 721-726.

7

42.

Nika, D. L.; Ghosh, S.; Pokatilov, E. P.; Balandin, A. A. Appl. Phys. Lett. 2009, 94,

8

(20), 203103.

9

43.

Feng, T.; Ruan, X.; Ye, Z.; Cao, B. Phys. Rev. B 2015, 91, (22), 224301.

10

44.

Zou, J.-H.; Cao, B.-Y. Appl. Phys. Lett. 2017, 110, (10), 103106.

11

45.

Lindsay, L.; Li, W.; Carrete, J.; Mingo, N.; Broido, D. A.; Reinecke, T. L. Phys. Rev.

12

B 2014, 89, (15), 155426.

13

46.

Stuart, S. J.; Tutein, A. B.; Harrison, J. A. J. Chem. Phys. 2000, 112, (14), 6472-6486.

14

47.

Yang, N.; Xu, X.; Zhang, G.; Li, B. AIP Advances 2012, 2, (4), 041410.

15

48.

Ghosh, S.; Bao, W.; Nika, D. L.; Subrina, S.; Pokatilov, E. P.; Lau, C. N.; Balandin,

16

A. A. Nat. Mater. 2010, 9, (7), 555-8.

17

49.

Nika, D. L.; Askerov, A. S.; Balandin, A. A. Nano Lett. 2012, 12, (6), 3238-3244.

18

50.

Nika, D. L.; Balandin, A. A. J. Phys.: Condens. Matter 2012, 24, (23), 233203.

19

51.

Liu, Y.; Huang, J.; Yang, B.; Sumpter, B. G.; Qiao, R. Carbon 2014, 75, 169-177.

20

52.

Hu, L.; Desai, T.; Keblinski, P. Phys. Rev. B 2011, 83, (19), 195423.

21

ACS Paragon Plus Environment

Page 23 of 27

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Nano Letters

Figure 1 (Color on-line) Representative temperature profile in TTM for the coupled phonons. The system consists of two regions, namely, the left region (l) (where atoms have both OP and IP vibrations) and the right region (r) (only IP vibration). To,l, Ti,l are temperatures for OP and IP phonon group in the left region. Ti,r is the temperature of IP phonon group in the right region. θmax is the maximum temperature difference between OP and IP phonon group. lc is the ph-ph coupling length, which is defined as the distance between the position of θmax and 5%θmax. 287x190mm (300 x 300 DPI)

ACS Paragon Plus Environment

Nano Letters

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 2 (Color on-line) (a) Schematic illustration of molecular dynamics simulation setup. The black (green) arrow denotes the vibration along the out-of-plane (in-plane) direction. The blue dot line is the interface of different simulation regions. The atoms in the left region vibrate freely in three directions while those in the right region could only vibrate along in-plane direction. The fixed (periodic) boundary condition is applied to along x (y) direction. Hot bath (red atoms) TH and cold bath (blue atoms) TC are applied. The temperatures of two heat baths are set as TH=310 K, TC=290 K, respectively. (b) The temperature profile of different phonon groups. (c) The temperature difference distribution between OP and IP phonon group when the system length is set as L=101 nm. The fitting line is based on θ=θmaxexp(γx). The red dot lines denote the maximum temperature difference and its five percent, respectively. 270x164mm (300 x 300 DPI)

ACS Paragon Plus Environment

Page 24 of 27

Page 25 of 27

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Nano Letters

Figure 3. (Color on-line) The length dependence of thermal conductivity with different phonon groups in nanosized graphene calculated from NEMD. Pristine graphene (Pristine), atoms vibrated along two in-plane directions only (Only IP), atoms vibrated along our-of-plane direction only (only OP), and the summation of the values of only IP and only OP. The system length ranges from 16 nm to 100 nm. The fitting lines are based on κ ~ logL. 174x201mm (300 x 300 DPI)

ACS Paragon Plus Environment

Nano Letters

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 4. (Color on-line) (a) The ph-ph coupling factor versus the system lengths calculated by MD and ab initio methods. (b) The length dependence of the ph-ph coupling length from MD simulation, where the thermal conductivity of graphene follows κ~logL. 287x201mm (300 x 300 DPI)

ACS Paragon Plus Environment

Page 26 of 27

Page 27 of 27

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Nano Letters

TOC Graphic 40x35mm (300 x 300 DPI)

ACS Paragon Plus Environment