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