- Email: [email protected]
Molecular dynamics simulation of temperature profile in partially hydrogenated graphene and graphene with grain boundary
Journal of Molecular Graphics and Modelling 62 (2015) 38–42
Contents lists available at ScienceDirect
Journal of Molecular Graphics and Modelling journal homepage: www.elsevier.com/locate/JMGM
Molecular dynamics simulation of temperature profile in partially hydrogenated graphene and graphene with grain boundary Erfan Lotfi a , M. Neek-Amal b,∗ , M. Elahi a a b
Department of Physics, Science and Research Branch, Islamic Azad University, Tehran, Iran Department of Physics, Shahid Rajaee Teacher Training University, Lavizan, Tehran 16788, Iran
a r t i c l e
i n f o
Article history: Received 27 May 2015 Received in revised form 13 August 2015 Accepted 14 August 2015 Available online 24 August 2015 Keywords: Graphene Hydrogenated graphene Grain boundaries Temperature distribution
a b s t r a c t Temperature profile in graphene, graphene with grain boundary and vacancy defects and hydrogenated graphene with different percentage of H-atoms are determined using molecular dynamics simulation. We also obtained the temperature profile in a graphene nanoribbon containing two types of grain boundaries with different misorientation angles, = 21.8◦ and = 32.2◦ . We found that a temperature gap appears in the temperature profile of a graphene nanoribbon with a grain boundary at the middle. Moreover, we found that the temperature profile in the partially hydrogenated graphene varies with the percentage of hydrogens, i.e. the C:H ratio. Our results show that a grain boundary line in the graphene sheet can change the thermal transport through the system which might be useful for controlling thermal flow in nanostructured graphene. © 2015 Elsevier Inc. All rights reserved.
1. Introduction Graphene is a crystalline allotrope of carbon with two dimensional honeycomb lattice structure and sp2 in-plane covalent bonds [1,2], which gained considerable attention in recent decade due to its exceptional physical properties [3–8]. Thermal properties of graphene are particularly interesting because it has the largest thermal conductivity among other materials, i.e. of about 1000 W m−1 K−1 . Therefore, graphene is a candidate for applications of thermal transport and heat management in nanoelectronic devices [9]. Faugeras et al. used local laser excitation and temperature readout from the intensity ratio of Stokes to anti-Stokes Raman scattering signals in order to study the temperature distribution in graphene [10]. They found a radial distribution for the temperature profile and reported thermal conductivity of graphene to be about 630 W m−1 K−1 . On the other hand several theoretical works on the thermal conductivity and thermal transport in graphene and hydrogenated graphene [11,12] confirm high thermal conductivity in the carbon based nano-structures (for a review on thermal conductivity of graphene see [13] and references therein). One of the proposed applications of thermal conductivity of graphene is
∗ Corresponding author at: Department of Physics, Shahid Rajaee Teacher Training University, Lavizan, Tehran 16788, Iran. E-mail address: [email protected] (M. Neek-Amal). http://dx.doi.org/10.1016/j.jmgm.2015.08.007 1093-3263/© 2015 Elsevier Inc. All rights reserved.
to create a driving force with a temperature gradient in order to move objects [14,15]. On the other hand hydrocarbons, i.e. (CH)n , are the simplest organic structures made of merely carbon and hydrogen atoms [16] which have different thermal conductivity than pristine graphene [17]. The latter is due to the sp3 bond formation in the hydrogenated graphene. Experimentally, it has been shown that hydrogenated graphene is obtained reversibly starting from a pristine graphene layer [18]. It is also a material of high interest due to its potential applications in nanoelectronics [19]. We recently found that fully hydrogenated graphene is an un-rippled system in contrast to graphene [20] and has a lower melting temperature than graphene [21]. “In the presence of F/H adatoms the C-bonds in graphene transit from sp2 to sp3 hybridization, which turn the conjugated, graphitic C C bonds into single C C bonds. This re-hybridization turns the lattice structure into an angstrom scale out-of-plane buckled shaped membrane known commonly as chair configuration.” Moreover, graphene with various kinds of defects has been extensively investigated in the past few years [22], e.g. the buckling of graphene nano-ribbons with a grain boundary were studied using atomistic simulations with free and supported boundary conditions [23]. It was found that graphene with vacancy defects is a less stiff material as compared to perfect graphene [24]. In this study, using state of the art molecular dynamics simulations, we obtain the radial profile of temperature in perfect graphene, graphene with vacancy defects and graphene with grain boundary. We present an analytical model which gives
E. Lotfi et al. / Journal of Molecular Graphics and Modelling 62 (2015) 38–42
39
good agreement with our MD results for pristine graphene and hydrogenated graphene. We show that the temperature profile in graphene is strongly affected by the presence of randomly distributed vacancy defects, grain boundaries but slightly varies with the number of hydrogens. We also study the temperature profile along a graphene nanoribbons containing a grain boundary. We found that the well known linear temperature gradient along a system with a grain boundary is divided into two lines which are separated by a temperature gap. The paper is organized as the following. In Section 2 we present our model. In Section 3 we review the used molecular dynamics simulation method. In Section 4 the radial temperature profile is and corresponding results are presented. Finally we summarize the paper in Section 5.
ring. The temperature of the middle rings is calculated in each time step and averaged over 1.22 ns.
2. The model
∇ 2 T (r) + q(r) = 0,
The partially hydrogenated systems are generated by removing randomly H atoms form both sides of graphane. The presence of two high-symmetry directions in graphene, armchair and zigzag, and any (local) misorientation angle close to 0 or 60 can be considered as small-angle grain boundary along these two directions, respectively. The large-angle grain-boundary, LAGBI, LAGBII, which are characterized by = 21.8◦ and = 32◦ has already been suggested in the literature [25]. In fact the two kind of grain boundaries, i.e. LAGBI and LAGBII are typical interfaces between domains of graphene with different crystallographic orientation. Mutual orientations of the two crystalline domains are described by the misorientation angle which for LAGBI is = 21.8◦ and for LAGBII is = 32.2◦ [26]. Top view of the pristine graphene, graphane (chair like configuration of fully hydrogenated graphene), graphene with grain boundary line which are named LAGBI and LAGBII are shown in Fig. 1(a)–(d), respectively. Two typical systems with 20% and 80% hydrogens are shown in Fig. 1(e)–(f). A temperature gradient is introduced in the radial direction by increasing the temperature of a small region (hot spot) in the center of the system, see Fig. 2(a). Such a model simulates the experimental set-up of Ref. [10]. For graphene nanoribbons we apply a temperature gradient along the longitudinal direction by keeping the two ends at different temperatures.
where q(r) ∝ P exp(−|r − r0 |2 / 2 ) is the heat source which causes heat flow into the system (it can be generated by light absorption from a laser beam with power P focused around r0 ) and is the thermal conductivity of the system (in practice corresponds to the estimated radius of the laser spot on the sample). The solution of Eq. (1) is a logarithmic function
3. Computational methods Here we give details about our molecular dynamics simulations. All molecular dynamics simulations were performed using the Large-Scale Atomic/Molecular Massively Parallel Simulator package (LAMMPS) [33]. The simulated sample is a square sheet of graphene with length 35 nm in zig-zag and arm chair direction. In order to apply a temperature gradient we keep the temperature of the central atoms (a circle with 476 atoms) to TH = 310 K and the outer circle to TC = 290 K. To calculate the middle region temperature we divided our sample into several rings (which are shown by different colors in Fig. 2(a)). We used Tersoff potential [28] to model the interaction between carbon atoms in graphene and AIREBO potential [29,30] to describe the interaction between C C and C H atoms in partially hydrogenated graphene and graphane. The temperature of the hot and the cold spots are controlled by Nose Hoover thermostat [31,32] using an NVT ensemble. In order to find steady state condition we perform 100 ps simulation. During equilibration the middle rings are simulated by an NVE ensemble. Moreover, to prevent non-desired crumpling at the edges the other atoms beyond the outer circle shown in Fig. 2(a) are fixed. Because the two hot and cold regions are hold at different temperatures, a steady heat flux flows from the innermost circle to the outermost
4. Radial temperature profile At some initial time, say t = 0, a nonuniform temperature profile T(r, t) is applied to graphene. This temperature profile generates a local heat current J(r, t). Due to energy conservation and basic thermodynamics (continuity equation) one can write q(r, t) = ca (T )
∂T (r, t) = −∇ · J, ∂t
(1)
where ca is the specific heat per unit area. Using Fourier’s law for the steady state the latter equation can be written in a time independent form and yields the radial temperature distribution in the graphene sheet by solving the following radial differential equation
T ≈ T0 −
˛ ln 2
r r0
(2)
−
˛ , 4
(3)
where T0 is the central region temperature, ˛ = (P/(d)) and = 0.5772 is the Euler’s constant with ‘d’ the thickness of membrane [10]. In order to numerically analyze our results, we use a simplified version of Eq. (3) by introducing two fitting parameters A and T1 : T ≈ T1 − A ln(r),
(4)
where T1 is linearly proportional to the temperature in the central region and ‘A’ is determined by thermal conductivity, the total absorbed laser power for given laser spot region and the edge temperature. 4.1. Graphene By applying a temperature gradient of about 20 K between the two hot and cold spots of the system, and after reaching steady state, we calculated the radial distribution of temperature (see Fig. 3). In Fig. 3 we show two typical MD simulations results (symbols) and corresponding fits (solid curve and solid line) for the systems depicted in Fig. 2, i.e. perfect graphene and graphene with 2% vacancy defects. We found that the temperature profile for LAGBI and LAGBII are close to each other and to the one with 2% vacancy defect. There is an obvious difference in the temperature profile between defected graphene (linear fit) and perfect graphene (logarithmic fit). It is seen that the radial temperature profile in the defected graphene does not follow the above theoretical model. Therefore using Eq. (4) is only valid for perfect graphene without any kind of defects such as vacancies and grain boundaries. Since Eq. (4) is based on a continuum model we conclude that a continuum model is not applicable for the heat transfer in defected graphene. The fitting parameters for different studied systems are listed in Table 1. 4.2. Graphane and partially hydrogenated graphene We also studied the temperature profile in graphane and partially hydrogenated graphene with different C:H ratio. Using ‘A’ as a fitting parameter in Eq. (4), we found that the radial temperature profile obeys the above mentioned continuum model
40
E. Lotfi et al. / Journal of Molecular Graphics and Modelling 62 (2015) 38–42
Fig. 1. Perfect graphene (a), fully hydrogenated graphene (graphane) (b), graphene with a grain boundary line with two different misorientation angles = 21.8◦ (LAGBI) (c), and = 32.2◦ (LAGBII) (d). The lines are guides to the eye to see the orientation of the grain boundaries. The panels (e) and (f) show two systems with 20% and 80% hydrogens with respect to fully hydrogenated system.
(Eq. (4)). In Fig. 4, we show two typical temperature profiles for hydrogenated graphene with 20% and 80% randomly distributed hydrogen. Removing hydrogens from fully covered graphene, only alters slightly the logarithmic function of the temperature profile Table 1 Fitting parameters in Eq. (4) for various systems.
Graphene 20%H 40%H 60%H 80%H 100%
A
T1
5.299 11.34 12.60 13.63 14.09 9.53
318.14 317.13 319.44 321.09 321.80 312.37
but the overall change obeys Eq. (4). Since in the defected graphene (presence of vacancies and grain boundaries) some of the sp2 bonds are removed and some dangling bonds are formed (unsaturated bonds), we expect that the continuum model is no longer valid. However the presence of hydrogens in the system does not change the hexagonal network of graphene while introduces some local sp3 bonds with local buckling, i.e. the system can be considered as continuum plate and Eq. (4) is applicable.
4.3. Temperature profile along graphene nanoribbon In order to investigate the temperature profile in graphene nanoribbons with LAGBI and LAGBII grain boundaries, we applied
E. Lotfi et al. / Journal of Molecular Graphics and Modelling 62 (2015) 38–42
41
Fig. 2. The used models to study the temperature distribution in (a) perfect graphene and (b) graphene with 2% vacancies. Different colors refer to the different segments which are used to calculate the local temperatures.
a temperature gradient along the x-direction of a graphene nanoribbon with LAGBI and LAGBII grain boundaries. When heat flow is restricted to one direction along a two-dimensional planar surface, the solution of Fourier heat diffusion equation: d dx
dT (x)
dx
= 0,
(5)
where x is either zig-zag or armchair direction. If the distance between a hot reservoir (kept at T = Thot ) and a cold sink (kept at T = Tcold ) is L within the plane, then Eq. (5) can be solved and results in the following solution for T(x) T (x) = Tcold +
(Thot − Tcold )x . L
(6)
The temperature of the hot and cold ends (see Fig. 5(a)) is kept at TH = 310 K and TL = 290 K, respectively. It is important to notice that we divided the ribbon into several segments (which are shown by different colors in Fig. 5(a)) and the temperature in each rectangle was calculated by averaging over 1.22 ns. The temperature
305
Graphene A+B*Ln(r) 2% vacancies A+B*r
profile along the graphene nanoribbon is shown in Fig. 5(b). For comparison purposes, we simulate perfect nanoribbons (without grain boundary) of graphene with arm-chair and zig-zag direction along the temperature gradient. Surprisingly, there is a clear gap in the temperature profile of LAGBI and LAGBII nanoribbons, i.e. a considerable change in the temperature before and after LAGB lines which act as a heat flow barrier. It is attributed to the fact that when the heat flux flows along the ribbons with grain boundaries, these defects scatter heat (originated from phonon scattering). Furthermore, along the grain boundaries, the possibility of inelastic phonon scattering increases [34]. The temperature gap created in LAGBI is a little larger than that of LAGBII. This effect is promising for controlling and engineering heat-flow in graphene devices which also motivates more experimental and theoretical studies. In conclusion, in the perfect graphene as we see from circular and square dots in Fig. 5, a linear decrease is obtained in agreement with Eq. (6). However in the presence of grain boundary, we obtained two different regions in both sides of the grain boundary line which are separated by a large gap in the midpoint. In each side we obtain a linear decrease. It is interesting that the grain boundary plays the role of hot reservoir kept at T = 200 K (cold sink kept at T = 150 K) for the right (left) hand side region. The slope of all the fitted lines are
305
300
T (K)
20%H A+B*Ln(r) 80%H A+B*Ln(r)
300
T(K)
295
295
290
50
100
150
r (Å) Fig. 3. Radial temperature profile in different studied systems shown in Fig. 2, i.e. perfect graphene and graphene with 2% vacancy. The solid curve is a fit on the MD data for temperature profile of perfect graphene using Eq. (4) (blue curve) and the solid line is a linear fit (gray line) to the MD data for graphene with 2% vacancy. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of the article.)
290
50
100
150
r (Å) Fig. 4. Radial temperature profile in two partially hydrogenated systems with randomly 20% and 80% hydrogens. The solid lines are fits using Eq. (4).
42
E. Lotfi et al. / Journal of Molecular Graphics and Modelling 62 (2015) 38–42
graphene, and graphene with grain boundary were investigated. We discovered a temperature gap along the grain boundary line in a graphene nanoribbon which can alter heat transfer through the system. Such a grain boundary acts as a barrier in front of the heat flux. The radial temperature profile obtained form continuum model (T1 − A ln(r)) is only applicable on the perfect graphene and the hydrogenated graphene. The temperature profile in defected graphene and graphene with grain boundary linearly decrease through the system. Acknowledgement We would like to acknowledge Francois Peeters for his valuable comments. References
Fig. 5. (a) The model which shows the hot (red) and the cold (blue) regions in LAGBII system with TH = 310 K and TC = 290 K. (b) The MD results for the temperature profile in the systems with LAGBI and LAGBII grain boundaries (symbols) and corresponding linear fits (solid lines). The gray rectangle refers to the position of grain boundary. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of the article.)
almost the same while for both sides the length in Eq. (6) should be substituted by L/2. 5. Conclusions In summary, using molecular dynamics simulations temperature profile in graphene, graphane, partially hydrogenated
[1] K.S. Novoselov, A.K. Geim, S.V. Morozov, D. Jiang, Y. Zhang, S.V. Dubonos, I.V. Grigorieva, A.A. Firsov, Science 306 (2004) 666. [2] A.K. Geim, K.S. Novoselov, Nat. Mater. 6 (2007) 186. [3] A.K. Geim, Science 324 (2009) 1530. [4] C. Xu, H. Li, K. Banerjee, IEEE Trans. Electron Devices 55 (2008) 3264. [5] F. Schedin, A.K. Geim, S.V. Morozov, E.W. Hill, P. Blake, M.I. Katsnelson, K.S. Novoselov, Nat. Mater. 6 (2007) 652. [6] C. Gomez-Navarro, R.T. Weitz, A.M. Bittner, M. Scolari, A. Mews, M. Burghard, K. Kern, Nano Lett. 7 (2007) 3499. [7] T.J. Booth, P. Blake, R.R. Nair, D. Jiang, E.W. Hill, U. Bangert, A. Bleloch, M. Gass, K.S. Novoselov, M.I. Katsnelson, A.K. Geim, Nano Lett. 8 (2008) 2442. [8] F. Wang, Y.B. Zhang, C.S. Tian, C. Girit, A. Zettl, M. Crommie, Y.R. Shen, Science 320 (2008) 206. [9] A.A. Balandin, S. Ghosh, W.Z. Bao, I. Calizo, D. Teweldebrhan, F. Miao, C.N. Lau, Nano Lett. 8 (2008) 902. [10] C. Faugeras, B. Faugeras, M. Olita, M. Potemski, R.R. Nair, A.K. Geim, ACS Nano 4 (2010) 1889. [11] H. Shen, Micro Nano Lett. 8 (2013) 740. [12] M. Neek-Amal, F.M. Peeters, Phys. Rev. B 83 (2011) 235437. [13] K.M.F. Shahil, A.A. Balandin, Solid State Commun. 152 (2012) 1331. [14] A. Lohrasebi, M. Neek-Amal, M.R. Ejtehadi, Phys. Rev. E 83 (2011) 042601. [15] A. Barreiro, R. Rurali, E.R. Hernandez, J. Moser, T. Pichler, L. Forro, A. Bachtold, Science 320 (2008) 775. [16] G.A. Olah, A. Molnar, Hydrocarbon Chemistry, Wiley Interscience, Hoboken, NJ, 2003. [17] Q.-X. Pei, Z.-D. Sha, Y.-W. Zhang, Carbon 49 (2011) 4752. [18] D.C. Elias, R.R. Nair, T.M.G. Mohiuddin, S.V. Morozov, P. Blake, M.P. Halsall, A.C. Ferrari, D.W. Boukhvalov, M.I. Katsnelson, A.K. Geim, K.S. Novoselov, Science 323 (2009) 610. [19] M.Z.S. Flores, P.A.S. Autreto, S.B. Legoas, D.S. Galvao, Nanotechnology 20 (2009) 465704; O. Leenaerts, H. Peelaers, A.D. Hernandez-Nieves, B. Partoens, F.M. Peeters, Phys. Rev. B 82 (2010) 195436. [20] S. Costamagna, M. Neek-Amal, J.H. Los, F.M. Peeters, Phys. Rev. B 86 (4) (2012) 041408. [21] S.K. Singh, S. Costamagna, M. Neek-Amal, F.M. Peeters, J. Phys. Chem. C 118 (2014) 4460. [22] A.H. Castro Neto, F. Guinea, N.M.R. Peres, K.S. Novoselov, A.K. Geim, Rev. Mod. Phys. 81 (2009) 109. [23] M. Neek-Amal, F.M. Peeters, Appl. Phys. Lett. 100 (2012) 101905. [24] M. Neek-Amal, F.M. Peeters, Appl. Phys. Lett. 97 (2010) 153118. [25] P. Simonis, C. Goffaux, P.A. Thiry, L.P. Biro, Ph. Lambin, V. Meunier, Surf. Sci. 511 (2002) 319. [26] O.V. Yazyev, S.G. Louie, Phys. Rev. B 81 (2010) 195420. [28] J. Tersoff, Phys. Rev. B 39 (1989) 5566. [29] S.J. Stuart, A.B. Tutein, J.A. Harrison, J. Chem. Phys. 112 (2000) 6472. [30] S. Plimpton, J. Comput. Phys. 117 (1995) 1. [31] S. Nose, J. Chem. Phys. 81 (1984) 511. [32] W.G. Hovver, Phys. Rev. A 31 (1985) 1695. [33] See http://lammps.Sandia.gov for more information about LAMMPS. [34] V. Samvedi, V. Tomar, Nanotechnology 20 (2009) 365701.
Contents lists available at ScienceDirect
Journal of Molecular Graphics and Modelling journal homepage: www.elsevier.com/locate/JMGM
Molecular dynamics simulation of temperature profile in partially hydrogenated graphene and graphene with grain boundary Erfan Lotfi a , M. Neek-Amal b,∗ , M. Elahi a a b
Department of Physics, Science and Research Branch, Islamic Azad University, Tehran, Iran Department of Physics, Shahid Rajaee Teacher Training University, Lavizan, Tehran 16788, Iran
a r t i c l e
i n f o
Article history: Received 27 May 2015 Received in revised form 13 August 2015 Accepted 14 August 2015 Available online 24 August 2015 Keywords: Graphene Hydrogenated graphene Grain boundaries Temperature distribution
a b s t r a c t Temperature profile in graphene, graphene with grain boundary and vacancy defects and hydrogenated graphene with different percentage of H-atoms are determined using molecular dynamics simulation. We also obtained the temperature profile in a graphene nanoribbon containing two types of grain boundaries with different misorientation angles, = 21.8◦ and = 32.2◦ . We found that a temperature gap appears in the temperature profile of a graphene nanoribbon with a grain boundary at the middle. Moreover, we found that the temperature profile in the partially hydrogenated graphene varies with the percentage of hydrogens, i.e. the C:H ratio. Our results show that a grain boundary line in the graphene sheet can change the thermal transport through the system which might be useful for controlling thermal flow in nanostructured graphene. © 2015 Elsevier Inc. All rights reserved.
1. Introduction Graphene is a crystalline allotrope of carbon with two dimensional honeycomb lattice structure and sp2 in-plane covalent bonds [1,2], which gained considerable attention in recent decade due to its exceptional physical properties [3–8]. Thermal properties of graphene are particularly interesting because it has the largest thermal conductivity among other materials, i.e. of about 1000 W m−1 K−1 . Therefore, graphene is a candidate for applications of thermal transport and heat management in nanoelectronic devices [9]. Faugeras et al. used local laser excitation and temperature readout from the intensity ratio of Stokes to anti-Stokes Raman scattering signals in order to study the temperature distribution in graphene [10]. They found a radial distribution for the temperature profile and reported thermal conductivity of graphene to be about 630 W m−1 K−1 . On the other hand several theoretical works on the thermal conductivity and thermal transport in graphene and hydrogenated graphene [11,12] confirm high thermal conductivity in the carbon based nano-structures (for a review on thermal conductivity of graphene see [13] and references therein). One of the proposed applications of thermal conductivity of graphene is
∗ Corresponding author at: Department of Physics, Shahid Rajaee Teacher Training University, Lavizan, Tehran 16788, Iran. E-mail address: [email protected] (M. Neek-Amal). http://dx.doi.org/10.1016/j.jmgm.2015.08.007 1093-3263/© 2015 Elsevier Inc. All rights reserved.
to create a driving force with a temperature gradient in order to move objects [14,15]. On the other hand hydrocarbons, i.e. (CH)n , are the simplest organic structures made of merely carbon and hydrogen atoms [16] which have different thermal conductivity than pristine graphene [17]. The latter is due to the sp3 bond formation in the hydrogenated graphene. Experimentally, it has been shown that hydrogenated graphene is obtained reversibly starting from a pristine graphene layer [18]. It is also a material of high interest due to its potential applications in nanoelectronics [19]. We recently found that fully hydrogenated graphene is an un-rippled system in contrast to graphene [20] and has a lower melting temperature than graphene [21]. “In the presence of F/H adatoms the C-bonds in graphene transit from sp2 to sp3 hybridization, which turn the conjugated, graphitic C C bonds into single C C bonds. This re-hybridization turns the lattice structure into an angstrom scale out-of-plane buckled shaped membrane known commonly as chair configuration.” Moreover, graphene with various kinds of defects has been extensively investigated in the past few years [22], e.g. the buckling of graphene nano-ribbons with a grain boundary were studied using atomistic simulations with free and supported boundary conditions [23]. It was found that graphene with vacancy defects is a less stiff material as compared to perfect graphene [24]. In this study, using state of the art molecular dynamics simulations, we obtain the radial profile of temperature in perfect graphene, graphene with vacancy defects and graphene with grain boundary. We present an analytical model which gives
E. Lotfi et al. / Journal of Molecular Graphics and Modelling 62 (2015) 38–42
39
good agreement with our MD results for pristine graphene and hydrogenated graphene. We show that the temperature profile in graphene is strongly affected by the presence of randomly distributed vacancy defects, grain boundaries but slightly varies with the number of hydrogens. We also study the temperature profile along a graphene nanoribbons containing a grain boundary. We found that the well known linear temperature gradient along a system with a grain boundary is divided into two lines which are separated by a temperature gap. The paper is organized as the following. In Section 2 we present our model. In Section 3 we review the used molecular dynamics simulation method. In Section 4 the radial temperature profile is and corresponding results are presented. Finally we summarize the paper in Section 5.
ring. The temperature of the middle rings is calculated in each time step and averaged over 1.22 ns.
2. The model
∇ 2 T (r) + q(r) = 0,
The partially hydrogenated systems are generated by removing randomly H atoms form both sides of graphane. The presence of two high-symmetry directions in graphene, armchair and zigzag, and any (local) misorientation angle close to 0 or 60 can be considered as small-angle grain boundary along these two directions, respectively. The large-angle grain-boundary, LAGBI, LAGBII, which are characterized by = 21.8◦ and = 32◦ has already been suggested in the literature [25]. In fact the two kind of grain boundaries, i.e. LAGBI and LAGBII are typical interfaces between domains of graphene with different crystallographic orientation. Mutual orientations of the two crystalline domains are described by the misorientation angle which for LAGBI is = 21.8◦ and for LAGBII is = 32.2◦ [26]. Top view of the pristine graphene, graphane (chair like configuration of fully hydrogenated graphene), graphene with grain boundary line which are named LAGBI and LAGBII are shown in Fig. 1(a)–(d), respectively. Two typical systems with 20% and 80% hydrogens are shown in Fig. 1(e)–(f). A temperature gradient is introduced in the radial direction by increasing the temperature of a small region (hot spot) in the center of the system, see Fig. 2(a). Such a model simulates the experimental set-up of Ref. [10]. For graphene nanoribbons we apply a temperature gradient along the longitudinal direction by keeping the two ends at different temperatures.
where q(r) ∝ P exp(−|r − r0 |2 / 2 ) is the heat source which causes heat flow into the system (it can be generated by light absorption from a laser beam with power P focused around r0 ) and is the thermal conductivity of the system (in practice corresponds to the estimated radius of the laser spot on the sample). The solution of Eq. (1) is a logarithmic function
3. Computational methods Here we give details about our molecular dynamics simulations. All molecular dynamics simulations were performed using the Large-Scale Atomic/Molecular Massively Parallel Simulator package (LAMMPS) [33]. The simulated sample is a square sheet of graphene with length 35 nm in zig-zag and arm chair direction. In order to apply a temperature gradient we keep the temperature of the central atoms (a circle with 476 atoms) to TH = 310 K and the outer circle to TC = 290 K. To calculate the middle region temperature we divided our sample into several rings (which are shown by different colors in Fig. 2(a)). We used Tersoff potential [28] to model the interaction between carbon atoms in graphene and AIREBO potential [29,30] to describe the interaction between C C and C H atoms in partially hydrogenated graphene and graphane. The temperature of the hot and the cold spots are controlled by Nose Hoover thermostat [31,32] using an NVT ensemble. In order to find steady state condition we perform 100 ps simulation. During equilibration the middle rings are simulated by an NVE ensemble. Moreover, to prevent non-desired crumpling at the edges the other atoms beyond the outer circle shown in Fig. 2(a) are fixed. Because the two hot and cold regions are hold at different temperatures, a steady heat flux flows from the innermost circle to the outermost
4. Radial temperature profile At some initial time, say t = 0, a nonuniform temperature profile T(r, t) is applied to graphene. This temperature profile generates a local heat current J(r, t). Due to energy conservation and basic thermodynamics (continuity equation) one can write q(r, t) = ca (T )
∂T (r, t) = −∇ · J, ∂t
(1)
where ca is the specific heat per unit area. Using Fourier’s law for the steady state the latter equation can be written in a time independent form and yields the radial temperature distribution in the graphene sheet by solving the following radial differential equation
T ≈ T0 −
˛ ln 2
r r0
(2)
−
˛ , 4
(3)
where T0 is the central region temperature, ˛ = (P/(d)) and = 0.5772 is the Euler’s constant with ‘d’ the thickness of membrane [10]. In order to numerically analyze our results, we use a simplified version of Eq. (3) by introducing two fitting parameters A and T1 : T ≈ T1 − A ln(r),
(4)
where T1 is linearly proportional to the temperature in the central region and ‘A’ is determined by thermal conductivity, the total absorbed laser power for given laser spot region and the edge temperature. 4.1. Graphene By applying a temperature gradient of about 20 K between the two hot and cold spots of the system, and after reaching steady state, we calculated the radial distribution of temperature (see Fig. 3). In Fig. 3 we show two typical MD simulations results (symbols) and corresponding fits (solid curve and solid line) for the systems depicted in Fig. 2, i.e. perfect graphene and graphene with 2% vacancy defects. We found that the temperature profile for LAGBI and LAGBII are close to each other and to the one with 2% vacancy defect. There is an obvious difference in the temperature profile between defected graphene (linear fit) and perfect graphene (logarithmic fit). It is seen that the radial temperature profile in the defected graphene does not follow the above theoretical model. Therefore using Eq. (4) is only valid for perfect graphene without any kind of defects such as vacancies and grain boundaries. Since Eq. (4) is based on a continuum model we conclude that a continuum model is not applicable for the heat transfer in defected graphene. The fitting parameters for different studied systems are listed in Table 1. 4.2. Graphane and partially hydrogenated graphene We also studied the temperature profile in graphane and partially hydrogenated graphene with different C:H ratio. Using ‘A’ as a fitting parameter in Eq. (4), we found that the radial temperature profile obeys the above mentioned continuum model
40
E. Lotfi et al. / Journal of Molecular Graphics and Modelling 62 (2015) 38–42
Fig. 1. Perfect graphene (a), fully hydrogenated graphene (graphane) (b), graphene with a grain boundary line with two different misorientation angles = 21.8◦ (LAGBI) (c), and = 32.2◦ (LAGBII) (d). The lines are guides to the eye to see the orientation of the grain boundaries. The panels (e) and (f) show two systems with 20% and 80% hydrogens with respect to fully hydrogenated system.
(Eq. (4)). In Fig. 4, we show two typical temperature profiles for hydrogenated graphene with 20% and 80% randomly distributed hydrogen. Removing hydrogens from fully covered graphene, only alters slightly the logarithmic function of the temperature profile Table 1 Fitting parameters in Eq. (4) for various systems.
Graphene 20%H 40%H 60%H 80%H 100%
A
T1
5.299 11.34 12.60 13.63 14.09 9.53
318.14 317.13 319.44 321.09 321.80 312.37
but the overall change obeys Eq. (4). Since in the defected graphene (presence of vacancies and grain boundaries) some of the sp2 bonds are removed and some dangling bonds are formed (unsaturated bonds), we expect that the continuum model is no longer valid. However the presence of hydrogens in the system does not change the hexagonal network of graphene while introduces some local sp3 bonds with local buckling, i.e. the system can be considered as continuum plate and Eq. (4) is applicable.
4.3. Temperature profile along graphene nanoribbon In order to investigate the temperature profile in graphene nanoribbons with LAGBI and LAGBII grain boundaries, we applied
E. Lotfi et al. / Journal of Molecular Graphics and Modelling 62 (2015) 38–42
41
Fig. 2. The used models to study the temperature distribution in (a) perfect graphene and (b) graphene with 2% vacancies. Different colors refer to the different segments which are used to calculate the local temperatures.
a temperature gradient along the x-direction of a graphene nanoribbon with LAGBI and LAGBII grain boundaries. When heat flow is restricted to one direction along a two-dimensional planar surface, the solution of Fourier heat diffusion equation: d dx
dT (x)
dx
= 0,
(5)
where x is either zig-zag or armchair direction. If the distance between a hot reservoir (kept at T = Thot ) and a cold sink (kept at T = Tcold ) is L within the plane, then Eq. (5) can be solved and results in the following solution for T(x) T (x) = Tcold +
(Thot − Tcold )x . L
(6)
The temperature of the hot and cold ends (see Fig. 5(a)) is kept at TH = 310 K and TL = 290 K, respectively. It is important to notice that we divided the ribbon into several segments (which are shown by different colors in Fig. 5(a)) and the temperature in each rectangle was calculated by averaging over 1.22 ns. The temperature
305
Graphene A+B*Ln(r) 2% vacancies A+B*r
profile along the graphene nanoribbon is shown in Fig. 5(b). For comparison purposes, we simulate perfect nanoribbons (without grain boundary) of graphene with arm-chair and zig-zag direction along the temperature gradient. Surprisingly, there is a clear gap in the temperature profile of LAGBI and LAGBII nanoribbons, i.e. a considerable change in the temperature before and after LAGB lines which act as a heat flow barrier. It is attributed to the fact that when the heat flux flows along the ribbons with grain boundaries, these defects scatter heat (originated from phonon scattering). Furthermore, along the grain boundaries, the possibility of inelastic phonon scattering increases [34]. The temperature gap created in LAGBI is a little larger than that of LAGBII. This effect is promising for controlling and engineering heat-flow in graphene devices which also motivates more experimental and theoretical studies. In conclusion, in the perfect graphene as we see from circular and square dots in Fig. 5, a linear decrease is obtained in agreement with Eq. (6). However in the presence of grain boundary, we obtained two different regions in both sides of the grain boundary line which are separated by a large gap in the midpoint. In each side we obtain a linear decrease. It is interesting that the grain boundary plays the role of hot reservoir kept at T = 200 K (cold sink kept at T = 150 K) for the right (left) hand side region. The slope of all the fitted lines are
305
300
T (K)
20%H A+B*Ln(r) 80%H A+B*Ln(r)
300
T(K)
295
295
290
50
100
150
r (Å) Fig. 3. Radial temperature profile in different studied systems shown in Fig. 2, i.e. perfect graphene and graphene with 2% vacancy. The solid curve is a fit on the MD data for temperature profile of perfect graphene using Eq. (4) (blue curve) and the solid line is a linear fit (gray line) to the MD data for graphene with 2% vacancy. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of the article.)
290
50
100
150
r (Å) Fig. 4. Radial temperature profile in two partially hydrogenated systems with randomly 20% and 80% hydrogens. The solid lines are fits using Eq. (4).
42
E. Lotfi et al. / Journal of Molecular Graphics and Modelling 62 (2015) 38–42
graphene, and graphene with grain boundary were investigated. We discovered a temperature gap along the grain boundary line in a graphene nanoribbon which can alter heat transfer through the system. Such a grain boundary acts as a barrier in front of the heat flux. The radial temperature profile obtained form continuum model (T1 − A ln(r)) is only applicable on the perfect graphene and the hydrogenated graphene. The temperature profile in defected graphene and graphene with grain boundary linearly decrease through the system. Acknowledgement We would like to acknowledge Francois Peeters for his valuable comments. References
Fig. 5. (a) The model which shows the hot (red) and the cold (blue) regions in LAGBII system with TH = 310 K and TC = 290 K. (b) The MD results for the temperature profile in the systems with LAGBI and LAGBII grain boundaries (symbols) and corresponding linear fits (solid lines). The gray rectangle refers to the position of grain boundary. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of the article.)
almost the same while for both sides the length in Eq. (6) should be substituted by L/2. 5. Conclusions In summary, using molecular dynamics simulations temperature profile in graphene, graphane, partially hydrogenated
[1] K.S. Novoselov, A.K. Geim, S.V. Morozov, D. Jiang, Y. Zhang, S.V. Dubonos, I.V. Grigorieva, A.A. Firsov, Science 306 (2004) 666. [2] A.K. Geim, K.S. Novoselov, Nat. Mater. 6 (2007) 186. [3] A.K. Geim, Science 324 (2009) 1530. [4] C. Xu, H. Li, K. Banerjee, IEEE Trans. Electron Devices 55 (2008) 3264. [5] F. Schedin, A.K. Geim, S.V. Morozov, E.W. Hill, P. Blake, M.I. Katsnelson, K.S. Novoselov, Nat. Mater. 6 (2007) 652. [6] C. Gomez-Navarro, R.T. Weitz, A.M. Bittner, M. Scolari, A. Mews, M. Burghard, K. Kern, Nano Lett. 7 (2007) 3499. [7] T.J. Booth, P. Blake, R.R. Nair, D. Jiang, E.W. Hill, U. Bangert, A. Bleloch, M. Gass, K.S. Novoselov, M.I. Katsnelson, A.K. Geim, Nano Lett. 8 (2008) 2442. [8] F. Wang, Y.B. Zhang, C.S. Tian, C. Girit, A. Zettl, M. Crommie, Y.R. Shen, Science 320 (2008) 206. [9] A.A. Balandin, S. Ghosh, W.Z. Bao, I. Calizo, D. Teweldebrhan, F. Miao, C.N. Lau, Nano Lett. 8 (2008) 902. [10] C. Faugeras, B. Faugeras, M. Olita, M. Potemski, R.R. Nair, A.K. Geim, ACS Nano 4 (2010) 1889. [11] H. Shen, Micro Nano Lett. 8 (2013) 740. [12] M. Neek-Amal, F.M. Peeters, Phys. Rev. B 83 (2011) 235437. [13] K.M.F. Shahil, A.A. Balandin, Solid State Commun. 152 (2012) 1331. [14] A. Lohrasebi, M. Neek-Amal, M.R. Ejtehadi, Phys. Rev. E 83 (2011) 042601. [15] A. Barreiro, R. Rurali, E.R. Hernandez, J. Moser, T. Pichler, L. Forro, A. Bachtold, Science 320 (2008) 775. [16] G.A. Olah, A. Molnar, Hydrocarbon Chemistry, Wiley Interscience, Hoboken, NJ, 2003. [17] Q.-X. Pei, Z.-D. Sha, Y.-W. Zhang, Carbon 49 (2011) 4752. [18] D.C. Elias, R.R. Nair, T.M.G. Mohiuddin, S.V. Morozov, P. Blake, M.P. Halsall, A.C. Ferrari, D.W. Boukhvalov, M.I. Katsnelson, A.K. Geim, K.S. Novoselov, Science 323 (2009) 610. [19] M.Z.S. Flores, P.A.S. Autreto, S.B. Legoas, D.S. Galvao, Nanotechnology 20 (2009) 465704; O. Leenaerts, H. Peelaers, A.D. Hernandez-Nieves, B. Partoens, F.M. Peeters, Phys. Rev. B 82 (2010) 195436. [20] S. Costamagna, M. Neek-Amal, J.H. Los, F.M. Peeters, Phys. Rev. B 86 (4) (2012) 041408. [21] S.K. Singh, S. Costamagna, M. Neek-Amal, F.M. Peeters, J. Phys. Chem. C 118 (2014) 4460. [22] A.H. Castro Neto, F. Guinea, N.M.R. Peres, K.S. Novoselov, A.K. Geim, Rev. Mod. Phys. 81 (2009) 109. [23] M. Neek-Amal, F.M. Peeters, Appl. Phys. Lett. 100 (2012) 101905. [24] M. Neek-Amal, F.M. Peeters, Appl. Phys. Lett. 97 (2010) 153118. [25] P. Simonis, C. Goffaux, P.A. Thiry, L.P. Biro, Ph. Lambin, V. Meunier, Surf. Sci. 511 (2002) 319. [26] O.V. Yazyev, S.G. Louie, Phys. Rev. B 81 (2010) 195420. [28] J. Tersoff, Phys. Rev. B 39 (1989) 5566. [29] S.J. Stuart, A.B. Tutein, J.A. Harrison, J. Chem. Phys. 112 (2000) 6472. [30] S. Plimpton, J. Comput. Phys. 117 (1995) 1. [31] S. Nose, J. Chem. Phys. 81 (1984) 511. [32] W.G. Hovver, Phys. Rev. A 31 (1985) 1695. [33] See http://lammps.Sandia.gov for more information about LAMMPS. [34] V. Samvedi, V. Tomar, Nanotechnology 20 (2009) 365701.