Journal of Nanomaterials

Journal of Nanomaterials / 2010 / Article
Special Issue


View this Special Issue

Research Article | Open Access

Volume 2010 |Article ID 537657 |

Hengji Zhang, Geunsik Lee, Alexandre F. Fonseca, Tammie L. Borders, Kyeongjae Cho, "Isotope Effect on the Thermal Conductivity of Graphene", Journal of Nanomaterials, vol. 2010, Article ID 537657, 5 pages, 2010.

Isotope Effect on the Thermal Conductivity of Graphene

Academic Editor: Rakesh Joshi
Received01 Jun 2010
Accepted07 Jul 2010
Published15 Aug 2010


The thermal conductivity (TC) of isolated graphene with different concentrations of isotope (C13) is studied with equilibrium molecular dynamics method at 300 K. In the limit of pure C12 or C13 graphene, TC of graphene in zigzag and armchair directions are ~630 W/mK and ~1000W/mK, respectively. We find that the TC of graphene can be maximally reduced by ~80%, in both armchair and zigzag directions, when a random distribution of C12 and C13 is assumed at different doping concentrations. Therefore, our simulation results suggest an effective way to tune the TC of graphene without changing its atomic and electronic structure, thus yielding a promising application for nanoelectronics and thermoelectricity of graphene-based nano device.

1. Introduction

Since it was fabricated in 2004 [1], graphene, a monolayer of sp2-bonded network of carbon atoms, has attracted much attention for its unique electronic properties [2]. Meanwhile, both recent theory and experiment studies [3, 4] have revealed that isolated graphene has shown an unusual high thermal transport capability, which is of great importance in thermal management of nanoelectronics. Different from conventional metallic materials, thermal energy carriers for graphene are mainly in the form of phonon vibrations [35], and phonon contribution to TC is approximately 50 times larger than electron contribution at room temperature [3], which suggests that electron thermal transport is negligible in our case. Furthermore, the electronic contribution to TC would be independent of isotope effect as is electronically identical to . For these reasons, we will study only the phonon contribution to TC in this paper.

In addition to utilizing its high TC, another possible application of graphene has been investigated for thermoelectric energy conversion [6], where low TC but high electric conductivity for graphene is required for obtaining high thermoelectric efficiency. Such efficiency is expressed as , where σ is the electric conductivity, S is the Seebeck coefficient, and is the thermal conductivity. To achieve high ZT for graphene, a general scheme is to minimize while keeping and S less changed. One practical method is to dope the graphene with the stable isotope since electronic structure of graphene is unchanged in this doping. In a pure crystal, one without defects or dislocations, phonon scattering in the presence of different isotopes has been strongly correlated with changes in thermal conductivity. Similar to the observation for isotope effects on TC of Ge, diamond, and boron nitride nanotubes [7, 8], it is interesting to check how effectively isotope-doping method will reduce TC of Graphene.

Modeling of thermal transport can be achieved by using Boltzmann transport equation (BTE) [9, 10], or molecular dynamics (MD) simulations [11]. In BTE method, the single mode relaxation time (SMRT) approximation is a commonly used technique involving the assignment of the relaxation time to different phonon scattering mechanisms. The relaxation time can be either fitted to the experimental TC value [12] or determined from MD simulations [10]. In MD simulation approach, TC can be predicted from either nonequilibrium MD [13], where a temperature is applied across the simulation cell, or equilibrium MD [14], where the so-called Green-Kubo method is used to compute TC from heat current fluctuations. One advantage for MD method is that there is no assumption needed for phonon interactions. As long as the phonon dispersion and anharmonicity of the potential are accurate, MD method provides a robust way to accurately compute thermal transport.

In Section 2, we introduce Green-Kubo method and describe the simulation procedures to compute TC. In Section 3, we show the isotope effects on thermal transport of graphene, and discuss the reason for TC reduction and its possible application for thermoelectric application. Section 4 presents a summary and conclusions.

2. Green-Kubo MD Simulation Method

The Green-Kubo formula [14] derived from linear response theory can express thermal conductivity tensor in terms of equilibrium heat current-current autocorrelation in the form where is the system volume defined as the area of graphene multiplied with van der Waals thickness (3.35 Angstrom), kB is the Boltzmann constant, T is the system temperature, and τm is the time required to be longer than the time for current-current correlations to decay to zero [11]. is denoted as the heat flux in α or β direction, its expression is commonly defined as where , , and are the position, velocity and site energy of atom , respectively. is the potential energy at site . Then, we use Hardy’s definition [15], where is the force exerted on atom i by atom j. One merit of Hardy’s definition is that it is independent of pair-wise or many-body potential formulas. Based on the above equations, heat current at each step is recorded on disk as a quantity defined by atom positions, velocities, and interatomic forces, which can be extracted from atomistic simulations. The last step is to compute TC tensor with the known heat current in (1).

In the equilibrium MD simulations, we used the second generation REBO carbon potential for its accuracy in describing bond strength and anharmonicity of carbon materials [16]. A unit cell of graphene, with a periodic boundary condition as shown in Figure 1, has been thermalized to 300 K for 200 ps with Berendson thermostats. Afterwards, the heat current of graphene is recorded every 2 fs in nine microcanonical ensemble simulations with uncorrelated initial conditions. Each microcanonical simulation needs to run up to 10 ns to obtain converged TC value. Compared with nonequilibrium molecular dynamics (NEMD), the Green-Kubo method is indeed less computationally efficient; however, it is free of troubles such as finite size and boundary effects, which are commonly unavoidable in NEMD.

3. Results and Discussion

We computed the TC of mass defect-free graphene as 630 W/mk and 1000 W/mk in armchair and zigzag direction, respectively. This is lower than the reported experimental data [4, 5]. The reason has to do with the discrepancies on phonon dispersion between experiment measurement and theoretical data predicted with original REBO potential [18]. Although the absolute TC of graphene is overall underestimated in our simulation, the relative difference of TC caused by isotope mass defect is still physically meaningful, and the normalized TC reduction seen in Figure 4 properly reflects the isotope effect on TC of graphene.

One primary test done before studying isotope effects on thermal transport is the convergence test for graphene with different unit cell periodic boundary lengths. As shown in Figure 2, we have tested three cases with the boundary length of 1.6, 3.4, and 5.1 nm respectively, and the converged result suggests that the majority contribution from different phonon vibration modes are well included in our simulations. Then, to be efficient, we choose the 1.6 nm unit cell as the structure for the following simulations.

The first thing we learned from the simulations is that TC of pure or graphene in zigzag direction is 58% greater than that in armchair direction. This can be conceptually explained from the acoustic phonon dispersion and Grüneisen parameter of graphene calculated with REBO potential as shown in Figure 3. In the SMRT approximation, the BTE method can express TC as a summation of each phonon mode contribution, and where and represent phonon modes and wave momentum, and is the specific heat of the phonons, which is constant in the classical case. is the group velocity for mode , is the relaxation time, and T is the temperature. From Figure 3(a), we notice that both out-of-plane acoustic (ZA) and transverse acoustic (TA) modes in zigzag direction have apparently larger group velocity than that in armchair direction, longitudinal acoustic (LA) modes group velocity in both directions is quite similar. Then, from Figure 3(b), the Grüneisen parameter, which is plotted as a function of phonon mode and momentum , shows a smaller difference in zigzag and armchair directions for each vibration mode. As the Grüneisen parameters for LA, TA, and ZA are similar, it is reasonable to assume that the relaxation times for phonon vibration in zigzag and armchair direction are the same. Based on the above analysis, we conclude that the difference in the TC between zigzag and armchair directions mainly comes from their different group velocities.

In the isotope effect study, we generated a wide range of graphene samples with different C13 concentrations randomly distributed, since it is a more realistic possible configuration after the synthesis of graphene. Recently, Mingo et al. have proposed a possible method to generate isotope clusters in graphene, and theoretically demonstrated TC reduction using nonequilibrium Green’s function method [19]. Figure 4 shows the normalized graphene TC values as a function of concentration in armchair and zigzag directions, and the maximum TC reduction for graphene can be made between 25% and 75% of the atomic concentrations.

The explanation of the almost “parabolic” shape of TC in Figure 4 can be found in an earlier important relation derived by Klemens [20]. As it is commonly accepted, the longer the phonon mean free path is, the larger TC would be. Klemens found that, for isotope scattering, the mean free path is proportional to , where and is the temperature. Ci and Mi represent the concentration and the mass of the constituent isotope atoms, respectively. Thus, the mean free path directly has to do with the g factor, which is the mass variation of isotope atoms. In our simulation, factor reaches its maximum for around 50% of C13, so there we have the minimum phonon mean free path and the minimum TC. As atomic concentration approaches to zero or 100%, becomes smaller, phonon free path is larger, and TC increases. It will not get to infinity because of the Umklapp and other phonon scattering mechanisms.

Among various methods that modulate the TC of graphene, the isotope-doping method provides an efficient way to improve thermoelectric efficiency, since isotope atoms do not change the electronic structure of graphene, the electric conductivity, and the Seebeck coefficient remains the same after the doping. As our simulations suggest, in armchair and zigzag directions can be reduced by up to 80%, which means an improvement of the ZT by a factor of 5 since ZT is inversely proportional to TC. Recently, two important simulation studies were carried out on ZT of armchair [21] and zigzag [22] graphene nanoribbons (GNRs). In both works, the authors looked for defect and disordered structures that minimize TC while keep electric conductivity less changed, thus maximizing the ZT. Compared with this, the isotope-doping method appears to be another attractive method that would promote ZT even further for GNRs when it is combined with those methods [21, 22]. We believe our results can motivate new experiments to check the efficiency of isotope-doping method for thermoelectric application.

4. Conclusion

In summary, we have studied isotope effect on TC of isolated graphene with equilibrium MD methods. Our simulation results suggest that TC of graphene can be effectively reduced by up to 80% in armchair and zigzag directions for isotope concentrations as low as 25%. The phenomenon that mass defect can reduce TC is explained with the relation between phonon mean free path and mass variation of the isotope mixtures [20]. Finally, our study shows that doping graphene with carbon isotope C13 could be a practical way to reduce TC without changing its electric property, thus promoting thermoelectric coefficient.


The authors are grateful to support from Lockheed Martin. G. Lee acknowledges support from SWAN, and A. F. Fonseca acknowledges partial support from the Brazilian agency CNPq. We thank Rodney Ruoff for suggesting to study isotope effect in graphene. The authors also appreciate discussions with Davide Donadio and Giulia Galli.


  1. K. S. Novoselov, A. K. Geim, S. V. Morozov et al., “Electric field in atomically thin carbon films,” Science, vol. 306, no. 5696, pp. 666–669, 2004. View at: Publisher Site | Google Scholar
  2. A. K. Geim and K. S. Novoselov, “The rise of graphene,” Nature Materials, vol. 6, no. 3, pp. 183–191, 2007. View at: Publisher Site | Google Scholar
  3. K. Saito, J. Nakamura, and A. Natori, “Ballistic thermal conductance of a graphene sheet,” Physical Review B, vol. 76, no. 11, Article ID 115409, 2007. View at: Publisher Site | Google Scholar
  4. J. H. Seol, I. Jo, A. L. Moore et al., “Two-dimensional phonon transport in supported graphene,” Science, vol. 328, no. 5975, pp. 213–216, 2010. View at: Publisher Site | Google Scholar
  5. A. A. Balandin, S. Ghosh, W. Bao et al., “Superior thermal conductivity of single-layer graphene,” Nano Letters, vol. 8, no. 3, pp. 902–907, 2008. View at: Publisher Site | Google Scholar
  6. D. Dragoman and M. Dragoman, “Giant thermoelectric effect in graphene,” Applied Physics Letters, vol. 91, no. 20, Article ID 203116, 2007. View at: Publisher Site | Google Scholar
  7. D. T. Morelli, J. P. Heremans, and G. A. Slack, “Estimation of the isotope effect on the lattice thermal conductivity of group IV and group III-V semiconductors,” Physical Review B, vol. 66, no. 19, Article ID 195304, 9 pages, 2002. View at: Google Scholar
  8. C. W. Chang, A. M. Fennimore, A. Afanasiev et al., “Isotope effect on the thermal conductivity of boron nitride nanotubes,” Physical Review Letters, vol. 97, no. 8, Article ID 085901, 2006. View at: Publisher Site | Google Scholar
  9. A. J. H. McGaughey and M. Kaviany, “Quantitative validation of the Boltzmann transport equation phonon thermal conductivity model under the single-mode relaxation time approximation,” Physical Review B, vol. 69, no. 9, Article ID 094303, 12 pages, 2004. View at: Google Scholar
  10. D. Donadio and G. Galli, “Thermal conductivity of isolated and interacting carbon nanotubes: comparing results from molecular dynamics and the Boltzmann transport equation,” Physical Review Letters, vol. 99, no. 25, Article ID 255502, 2007. View at: Publisher Site | Google Scholar
  11. P. K. Schelling, S. R. Phillpot, and P. Keblinski, “Comparison of atomic-level simulation methods for computing thermal conductivity,” Physical Review B, vol. 65, no. 14, Article ID 144306, 12 pages, 2002. View at: Google Scholar
  12. S. E. Lemehov, V. Sobolev, and P. Van Uffelen, “Modelling thermal conductivity and self-irradiation effects in mixed oxide fuels,” Journal of Nuclear Materials, vol. 320, no. 1-2, pp. 66–76, 2003. View at: Publisher Site | Google Scholar
  13. J. Hu, X. Ruan, and Y. P. Chen, “Thermal conductivity and thermal rectification in graphene nanoribbons: a molecular dynamics study,” Nano Letters, vol. 9, no. 7, pp. 2730–2735, 2009. View at: Publisher Site | Google Scholar
  14. A. S. Henry and G. Chen, “Spectral phonon transport properties of silicon based on molecular dynamics simulations and lattice dynamics,” Journal of Computational and Theoretical Nanoscience, vol. 5, no. 2, pp. 141–152, 2008. View at: Google Scholar
  15. R. J. Hardy, “Energy-flux operator for a lattice,” Physical Review, vol. 132, no. 1, pp. 168–177, 1963. View at: Publisher Site | Google Scholar
  16. D. W. Brenner, O. A. Shenderova, J. A. Harrison, S. J. Stuart, B. Ni, and S. B. Sinnott, “A second-generation reactive empirical bond order (REBO) potential energy expression for hydrocarbons,” Journal of Physics Condensed Matter, vol. 14, no. 4, pp. 783–802, 2002. View at: Publisher Site | Google Scholar
  17. The PHONON package of D. Alfe, 2008,
  18. L. Lindsay and D. A. Brodio, “Optimized Tersoff and Brenner empirical potential parameters for lattice dynamics and phonon thermal transport in carbon nanotubes and graphene,” Physical Review B, vol. 81, Article ID 205441, 6 pages, 2010. View at: Publisher Site | Google Scholar
  19. N. Mingo, K. Esfarjani, D. A. Broido, and D. A. Stewart, “Cluster scattering effects on phonon conduction in graphene,” Physical Review B, vol. 81, no. 4, Article ID 045408, 6 pages, 2010. View at: Publisher Site | Google Scholar
  20. P. G. Klemens, “The scattering of low-frequency lattice waves by static imperfections,” Proceedings of the Physical Society A, vol. 68, no. 12, pp. 1113–1128, 1955. View at: Publisher Site | Google Scholar
  21. Y. Ouyang and J. Guo, “A theoretical study on thermoelectric properties of graphene nanoribbons,” Applied Physics Letters, vol. 94, no. 26, Article ID 263107, 2009. View at: Publisher Site | Google Scholar
  22. H. Sevinçli and G. Cuniberti, “Enhanced thermoelectric figure of merit in edge-disordered zigzag graphene nanoribbons,” Physical Review B, vol. 81, Article ID 113401, 2010. View at: Google Scholar

Copyright © 2010 Hengji Zhang et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

More related articles

 PDF Download Citation Citation
 Download other formatsMore
 Order printed copiesOrder

Related articles

We are committed to sharing findings related to COVID-19 as quickly as possible. We will be providing unlimited waivers of publication charges for accepted research articles as well as case reports and case series related to COVID-19. Review articles are excluded from this waiver policy. Sign up here as a reviewer to help fast-track new submissions.