Phonon-limited transport coefficients in extrinsic graphene Enrique Mun˜oz Instituto de F´ısica, Pontificia Universidad Cat´olica de Valpara´ıso PO Box 4059, Valpara´ıso, Chile. (Dated: January 6, 2012) The effect of electron-phonon scattering processes over the thermoelectric properties of extrinsic graphenewasstudied. Electricalandthermalresistivity,aswellasthethermopower,werecalculated within the Bloch theory approximations. Analytical expressions for the different transport coeffi- 2 cients were obtained from a variational solution of the Boltzmann equation. The phonon-limited 1 electrical resistivityρe−ph showsalineardependenceathightemperatures, andfollows ρe−ph∼T4 0 at low temperatures, in agreement with experiments and theory previously reported in the litera- 2 ture. The phonon-limited thermal resistivity at low temperatures exhibits a ∼T dependence, and achieves a nearly constant value at high temperatures. The predicted Seebeck coefficient at very n low temperatures is Q(T)∼−π2k2T/(3eE ), which shows a n−1/2 dependence with the density a B F J ofcarriers, inagreementwith experimentalevidence. Ourresultssuggest thatthermoelectricprop- erties can be controlled by adjusting the Bloch-Gru¨neisen temperature through its dependence on 5 theextrinsic carrier density in graphene. ] l PACSnumbers: 65.80.Ck,72.80.Vp,63.22.Rc l a h - Recently,EfetovandKim1reportedexperimentalmea- 200 s surements of the phonon-limited electrical resistivity of e m graphene samples, by achieving extremely high carrier 150 densities by means of an electrolytic gate, in order to . t minimize the effect of other scattering mechanisms and a 100 m to be away from the neutrality (Dirac) point. They ob- ) servedthat,atlowtemperatures,theelectricalresistivity Ω - ( d displays a T4 dependence, whereas at higher temper- ρ 50 ∼ 1 [con abagmtnlaeouscdereteadesslrelptiidhethretofhiontrehlolaeoninonwtradseylrryaaabsgciylstieniHooffennweacrabttnesh∼tgewaTaneBndeltdnoruceeDmhnleadkcts,ltharSieponoanprrasympg3aarr,n4eo2,de.cmewlTsoehsnhneigectsih,tlwauanitdnteeihdr-- -500 n (11642100....186837.85566 m -2) ρΩ ()(T=0)1150500000 0.1 0.2n-1 (100.-317 m-20).4 0.5 v nal acoustic phonons by the deformation potential ap- 0 50 100 150 200 250 300 7 proximation. Moreover, they calculated the electrical T (K) 5 resistivity2 within the relaxation time approximation, 0 1 which is valid for elastic scattering processes, a condi- FIG. 1: (Color online) Electrical resistivity calculated from . tion not strictly satisfied in electron-phonon scattering. Eq. (31)(lines),comparedwithexperimentaldata(symbols) 1 reported in Ref.(1) at different carrier densities n in units of 0 In this paper, phonon-limited transport coefficients in 1017m−2. Inset: Electrical resistivity at zero temperature, 2 fromdatainRef. (1)(circles),atdifferentcarrierdensitiesn. graphenearestudied,byconsideringthepresenceoftem- 1 Also shown is the linear fit (dashed line) ρ = 14.276+ v: ptieornastuarseinanRdevf.ol(t2a)g,ethgaratdiise,nBtsl.ocUhntdheerorsyimainladrdaesfsourmmpa-- 266.17×1017n−1 (Ω) explained in themain(Tt=e0x)t. i X tionpotential approximationfor the interactionbetween ar eelxepcrtersosniosnasndforlotnhgeitutrdainnaspl-oarctoucsoteicffipciheonntsonasr,eaonbatlyatiniceadl graphene samples at the neutrality point6. However, at high carrier densities the Fermi level is well above the fromavariationalsolutionoftheBoltzmannequation3,4. Dirac point, and the Boltzmann equation then provides This method, at the lowest-order approximation, pro- a correct description. vides equivalent results to the relaxation time approxi- mation and, at higher orders, converges to solutions of the Boltzmann equation which are exact to an arbitrary precision, within the limitations of the set of variational I. ELECTRON-PHONON INTERACTION functions chosen. As has been pointed out before5, the semiclassicalapproachbasedontheBoltzmannequation As is wellknownfrom tight-binding calculations,elec- cannot capture the subtleties of electronic transport in tronic transport in graphene displays relativistic fea- graphene right at the neutrality point. This limitation tures. This is because pristine graphene, free of impu- can, at least in part, be attributed to the inhomogene- rities and defects, behaves as a metal only in the so- ity in the charge distribution experimentally observedin called Dirac points where the conduction and valence 2 bandstouch. Thereexiststwonon-equivalentsuchpoints Notice the presence of the factor cos(θ/2) = cos([ϕ k − within the first Brillouin zone, defined by the vectors ϕk′]/2) which arises from the inner product of the two K = 2π 1/√3, 1 , with a the direct lattice param- pseudo-spinor functions u∗ u defined in Eq. (2). This 1,2 ± a k′· k eter. In the vicinity of each of these points, or val- factor, which would be absent from the inner product (cid:2) (cid:3) leys, it is shown that the dispersion relation is linear, of non-relativistic scalar Schr¨odinger wave-functions, is E (δk) = ~v δk, for δk = k K . In each val- a fingerprint of the relativistic features of electrons in η F η ± | | − ley K , the envelope electronic functions related to this graphene. The more obvious consequence of its presence η dispersion relation can thus be obtained as the pseudo- thesuppressionofbackscattering(θ =π),aphenomenon spinor solutions of an effective two-dimensional Dirac’s observed in graphene and directly related to so called equation, Kleintunneling6–8 in the contextof relativistic quantum mechanics. HˆKψk(±) =vFαˆ·pˆψk(±) =Ek(±)ψk(±). (1) We shallneglect umklappprocesses,andthus we have G = 0 for normal processes. In this case, scattering Here αˆ = (σˆ ,σˆ ) are the Pauli matrices, we have re- x y eventssatisfyenergyand”crystalmomentum”conserva- cfdooerrfirnpeeosdpsiotknivd≡einδ(gnkeefgiogaretenivavece)hcthvoearllsilc,eiynt,yoaranmrdaetliEhzeek(d±en)be=yrgt±yhee~ivtgoFetnkav.laaTlurheeaes twtiiooitnnh,,oǫθkn′teh=oebǫstkca±aint~tsωetrqhineagncdoannkdg′ilt=ei.oknS±qin2qc=.eFtkrh′o2em+rektlh2ee−val2anktttkes′rccoaestqt(ueθar)--, A of the sample, are given by the pseudo-spinors ing events occur within a thin layer near the Fermi sur- face, then k k′ k , which then implies sin(θ/2) = ψk(±)(r)= √12A e−eiiϕϕkk//22 eik·r ≡uk(±)e√ikA·r (2) qti/o2nkFre.strIinctsa∼tDheeb∼ryaenFgmeoodfellaatptipcreoxmimomateinotna,|itnhvisolvceodn|diin- (cid:18)± (cid:19) electron-phonon scattering, q min 2k ,q , with q with ϕk =arctan(ky/kx). theradiusoftheDebyecircle(≤foratw{o-dFimeDn}sionalsysD- In the language of second quantization, the conduc- tem). In geometrical terms, the most restrictive condi- tion electrons field for k in the vicinity of the K = 1,2 tion depends on whether the radius of the Fermi sur- 2π 1/√3, 1 valley is ± a face (a circle centered on each Dirac point for graphene) is larger or smaller than the Debye circle. When writ- (cid:2) (cid:3)ψˆ (r)= ψ(+)(r)cˆ (3) η k kση tenintermsofthefrequencyofthelongitudinal-acoustic Xk,σ phonon modes involved ωq = vsq, the former condition reads ~ω /k min Θ ,Θ , with Θ the Debye with cˆ Fermionic operators [cˆ ,cˆ† ] = q B ≤ { BG D} D kση kση k′σ′η′ + temperature, and ΘBG = 2(vs/vF)EF/kB the Bloch- δk,k′δσ,σ′δη,η′, that destroy (create) a conduction elec- Gru¨neisen temperature. In the experimental setup pre- tron with momentum ~k in the vicinity of the Kη val- sented in Ref.(1), the carrier densities involved are such ley, with spin component σ = , . We model the that the Fermi surface is contained within the Debye {↑ ↓} electron-phonon interaction by the deformation poten- circle, and therefore the Bloch-Gru¨neisen rather than tial approximation2, which is a reasonable one when the Debye temperature imposes the characteristic en- the Fermi surface possesses spherical symmetry3,4. The ergy scale for phonon scatterers. Notice that the analy- two-dimensional version of this is a Fermi ”circle” in sis of the scattering angle above also yields the relation graphene, which corresponds to the intersection of the cos(θ/2) = 1 (q/(2k ))2. Therefore, after Eq. (6) F Dirac cone and the constantenergy plane defined by the − the matrix element for the electron-phonon interaction Fermi energy at ǫ = E . In the deformation potential p k F Hˆint =i M aˆ aˆ† cˆ† cˆ is approximation,the operatorrepresentinglattice dilation e−ph k,q,η,σ q q− −q k+q,ση kση in the second quantization language is P (cid:16) (cid:17) 1/2 M = i(~/2ρ ω )1/2Dq 1 (q/2k )2 . (6) q m q F ∆ˆ(r)=i (2Aρ ω )−1/2(q eˆ ) aˆ†eiq·r aˆ e−iq·r (4) − − m q · L q − q Here, D is the deformation po(cid:16)tential1,3,4 co(cid:17)upling con- q X (cid:0) (cid:1) stant for graphene, and ρ =7.6 10−7 kg/m2 the sur- m Here, ωq =vsq is the frequency for longitudinal acoustic face mass density. × phononsingraphene,aˆ†isthecreationoperatorforlongi- q tudinalacousticphononmodes,andeˆ isthecorrespond- L ing polarization vector for these modes. The electron- II. BOLTZMANN EQUATION phonon interaction Hamiltonian, in the deformation po- tential approximation,is obtained from We introduce, in the usual form, the deviation of the distribution function from equilibrium by the expression Hˆein−tph =i d2r ψˆη†(r)D∆ˆ(r)ψˆη(r) (5) fk = fk0 − χk∂fk0/∂ǫk. Let us consider the electron- Z η=1,2 phonon interaction as the sole scattering mechanism. X Therefore, after Eqs.(5,6) we should consider processes ~q2D2 =i cos(θ/2) aˆ aˆ† cˆ† cˆ of the form k + q ⇄ k′. Under Bloch’s approxima- s2ρmωq q− −q k+q,ση kση tion, which is to assume that the phonon system is in k,Xq,η,σ (cid:16) (cid:17) 3 quasi-equilibrium, the linearized form of the Boltzmann 6 equation when both an external electric field and a tem- perature gradient are imposed can be expressed as 5 350 −vk·∇T∂∂fTk0 −vk·E(−e)∂∂fǫkk0 = -1W K) 4 -1K)235000 Here, gηgkηB=gTσgσkX′,=qh2{χakre−thχek′s}pPinkkq′a−nd{χvakl′le−yχdke}gPenkke′qriaci(e7s), -1κµ[] ( e-ph23 n (1011.73 6m -2)-1~κ (W m 11205050000 respectively,whilethetransitionrateforelectron-phonon 2.86 0 1 4.65 0 50 100 150 200 250 300 scattering processes is 10.8 T (K) Pkkq′ =(2π/~)δG,k′−k−q|Mq|2δ(ǫk−ǫk′ +~ωq) 00 50 100 150 200 250 300 T (K) ×n0qfk0(1−fk0′) (8) FIG. 2: (Color online) Phonon-limited thermal resistivity Hwietrhefmk0oimstehnetuFmermki,-Dainradcdni0qstrthibeutBioonsefodriesltercitbruotniosntatfoesr s[κitei−esphn]−,1in, cuanlcituslaotfed10f1r7omm−2E.qS.(h3o4w)natindtiffheereinnstetcairsritehredceonr-- phonon states with momentum q. All vectors are two- responding thermal conductivity, including the contribution dimensional, with G a vector of the reciprocal lattice. fromCoulombelasticscatteringwithimpurities,andnormal- Forextrinsicgraphene,the FermilevelE =~k v is ized with a nominal packing thicknessof 3.4 ˚A. F F F located above the Dirac point, and depends on the den- sityofcarriersnbytherelationk =√πn. InRef.(1),an F experimentalmethodispresentedwhichallowstocontrol L ) T, which implies κ = L L−1L L . If, the carrier density, and hence the Bloch-Gru¨neisen tem- TT ∇ TE EE ET − TT on the other hand, one considers the case T = 0, then perature. AsdiscussedinRefs.(1,9),thisinturnbecomes ∇ the relation J = L E is obtained. Thus, it is con- EE a wayto controlthe phonon-limitedelectricalresistivity. cluded that the electrical resistivity is ρ = L−1. The In this work, we will show that also the electronic com- variational method3,4 is based on the principleEoEf maxi- ponent of the thermal conductivity, as well as the ther- mal entropy production10, and provides an iterative but mopower(Seebeckcoefficient)whichdeterminesthermo- virtually exact procedure to solve the Boltzmann equa- electric effects, can be controlled in a similar way. We tion,andthereforetocalculatethetransportcoefficients, shallusethevariationalmethodtocalculatethephonon- within the limitations of the set of variational functions limited transportcoefficients which, in contrastwith the chosen. The principle is based on equating the macro- relaxation time approximation, it has the advantage of scopic entropy production rate, as expressed in terms avoiding the assumption of quasi-elasticity in the scat- of the macroscopic currents Eq.(10), with the entropy tering processes. production rate due to microscopic scattering events10, ∂ S = S˙ , as defined from the linearized Boltzmann t scatt equation scattering terms, III. COUPLED CHARGE AND HEAT TRANSPORT S˙scatt = kgηgTσ2 {χk−χk′}2Pkkq′ (11) For a system with thermal and electrostatic poten- B k,q,k′ X tial gradients, the macroscopic entropy production rate is givenin terms of the heat U and chargeJ fluxes by3,4 The electric current J and the heat current U are given by the expressions ∂ S =J E/T +U (1/T). (9) t · ·∇ The Onsager’s theorem10 indicates that the fluxes and J = g g ( e)v χ ∂fk0 η σ k k generalized potentials are then linearly coupled, − ∂ǫk k X ∂f0 J=LEEE+LET∇T U = gηgσ vk(ǫk−EF)χk∂ǫk (12) U=LTEE+LTT∇T (10) Xk k The thermopower Q can be expressed in terms of the The technique is based on expanding the distribution coefficients L by considering the case J = 0 in Eq. αβ function deviation from equilibrium χ = α φ (k) (10),toobtainE=−LETL−E1E∇T. Hence,oneconcludes in terms of a set of variational functikons of tihei foirm tthaiantedQfr=om−tLhEeTseLc−Eo1End. eTqhueattiohner,mUal=con(dLuctiLv−it1yLis ob- φfuin(kct)io=ns(ǫφk −(k)E=F)ik−1(uˆuˆa·nkd).φF(okr)a=mkiniuˆm(aǫlPsetEof)t,rwiael − TE EE ET− 1 · 2 · k − F 4 obtain the charge and heat currents Eq.(10), ∂f0 LEE = Ji[P−1]ijJj J =( e)g g v φ (k) i − η σ k i ∂ǫ Xij k Xk L = T−1 J [P−1] U ∂f0 ET − i ij j U =g g v (ǫ E )φ (k) (13) ij i η σ k k F i X − ∂ǫ k k L = U [P−1] J X TE i ij j ij X Here, g = g = 2 are the valley and spin degeneracies η σ L = T−1 U [P−1] U (19) in graphene, respectively. The coefficients P = P are TT i ij j ij ji − defined as follows Xij Therefore, we proceed to calculate the currents and P = (gηgσ)2 φ (k) φ (k′) φ (k) φ (k′) k′ Pij coefficients. Let us first consider the current J1. ij kBT k,k′,q{ i − i }{ j − j }Pkq The group velocity in Eq.(13) is given by vk = kˆvF, X with kˆ = k/k = (cos(ϕ),sin(ϕ)). The discrete sums (14) are treated in the usual quasi-continuum limit as two- dimensional integrals, A(2π)−2 d2k, and the In terms of the variational solution, after Eq.(12) the k → integral is calculated by using the change of variables electric and heat currents are written as d2k =dϕǫdǫ/(~v )2, P R F J = αiJi J = A(−e)vF +πdϕkˆkˆ uˆ dǫǫ2∂fǫ0 Xi 1 π2(~vF)3 Z−π · Z ∂ǫ U = αiUi (15) Ae EF 2 =uˆ (20) Xi π~ (cid:18)~vF(cid:19) whereas the macroscopic ∂ S and microscopic S˙ en- For the current U1, we consider the quasi-continuum t scatt limitforthesumoverwave-vectorsinEq.(13),withdirec- tropy generation rates in Eqs.(9,11) are expressed by tions kˆ = k/k = (cos(ϕ),sin(ϕ)). Using the integration variables d2k=dϕdǫǫ/(~v )2, we have ∂ S = α T−1J E T−2U T , F t i i i S˙scatt = TXi−1 (cid:2) αiαjP·ij. − ·∇ (cid:3) (16) U1 = π2A(~vvFF)3 Z−ππdϕkˆkˆ·uˆZ dǫǫ2(ǫ−EF)(cid:18)∂∂fǫǫ0(cid:19) Xi,j = uˆ2πA(kBT)2 EF (21) − 3 ~ (~v )2 F The variational problem is to find the set of coefficients α that maximize S˙ , subject to the constraint Similarly, for the current U2 we obtain j scatt { } ∂tS =S˙scatt, which is enforced by a Lagrange multiplier U = AvF π dϕkˆkˆ uˆ dǫǫ2(ǫ E )2 ∂fǫ0 λ, 2 π2(~v )3 · − F ∂ǫ F Z−π Z (cid:18) (cid:19) πA(k T)2 E 2 B F = uˆ (22) δαδ αiαjPTij −λ αi JiT·E −Ui· ∇TT2 =0 − 3 ~ (cid:18)~vF(cid:19) n i,j i (cid:18) (cid:19) Finally, it is straightforwardto verify by inspection that X X (17) J2 =( e)U1. − Let us now consider the coefficient P . From Eq.(8), 11 we first notice that the transition probability rate k′ The solution to this constrained variational problem is Pkq given by λ=1/2, and imposes the condition q = k′ k. We define the directions kˆ′ = (cos(ϕ′),sin(ϕ′))−and kˆ = (cos(ϕ′ θ),sin(ϕ′ θ)), with θ the scattering angle. As befor−e, T αi = [P−1]ij Jj E Uj ∇ (18) in the qu−asi-continuum limit for the sums over wave- · − · T j (cid:18) (cid:19) vectors, we choose the integration variables d2kd2k′ = X dϕ′dθdǫ′ǫ′dǫǫ/(~v )4, to obtain F Here,[P−1]istheinverseofthematrix[P] =P whose ij ij (k T)−1A2 elements are the coefficients Pij defined by Eq.(14). Af- P11 = (Bπ~v )4 dǫǫ dǫ′ǫ′ ter substituting this solution into the expressionsfor the F Z Z π π electric and heat currents, we obtain explicit analytical dθ dϕ′(q uˆ)2 k′ (23) formulas for the transport coefficient tensors defined in × · Pkq Z−π Z−π 5 In Appendix A it is shown that the angular integral At last, we consider the coefficient P . After our defini- 12 dϕ′(q uˆ)2 = πq2, with q = 2k sin(θ/2). Hence, sub- tion Eq.(14), we have F · stituting the explicit expressionfor the transitionproba- bRility rate, Eq.(23) reduces to P = (kBT)−1A2 π dθ π dϕ′ dǫǫ dǫ′ǫ′ 12 (π~v )4 2(k T)−1A2 π F Z−π Z−π Z Z π2B~(~vF)4 Z−πdθq2n0q|Mq|2 ×[(q·uˆ)2(ǫ−EF)+~ωq(k′·uˆ)(q·uˆ)]Pkkq′ (29) AsshowninAppendix A,the angularintegralsaregiven dǫǫ(ǫ+~ω )f0(1 f0 ) (24) × q ǫ − ǫ+~ωq by dϕ′(q uˆ)2 = πq2, and dϕ′(q uˆ)(k′ uˆ) = Z (π/2)[q2+k·′2 k2]. Therefore, using the·proced·ure de- Theenergyintegralisevaluatedusingthetechniquepre- R − R scribed previously, the integral in Eq.(29) is calculated sented in Appendix B. Using the change of variables and expressed in terms of the (z) functions defined in dθ =kF−1[1−(q/(2kF))2]−1/2dq, we obtain Appendix C, Jp (k T)−1A2k ~D2 2kF v 2 ~D2A2E2(k T)4 v T P = B F dq q4+4 s P = F B (T/Θ ) s 11 π2ρ (~v )2 v 12 π2ρ (~v )3(~v )5 J4 BG − v Θ m F Z0 (cid:18) F(cid:19) m s F (cid:20) F (cid:18) BG(cid:19) T 2 π2q4 + q6 [1 q 2]1/2 ×J5(T/ΘBG)+ 8π32 ΘT 2 J4(T/ΘBG)+ 23 ×(cid:18)ΘBG(cid:19) (cid:20) 3 2 (cid:21) −(cid:18)2kF(cid:19) ! (cid:18) BG(cid:19) (cid:26) 2 (1 e−~ωq/(kBT))−1(e~ωq/(kBT) 1)−1 (25) vs (T/Θ )+ 1 (T/Θ ) + 2 (30) × − − ×(cid:18)vF(cid:19) J4 BG 4π2J6 BG ) 3 Eq.(25)iswrittenintermsofthefunctions (z)defined Jp v 3 T 3 in Appendix C, s (T/Θ ) 2π2 (T/Θ ) 7 BG 5 BG ×(cid:18)vF(cid:19) (cid:18)ΘBG(cid:19) J − J # ~D2A2E (k T)4 v 2 (cid:8) (cid:9) F B s P = (T/Θ )+4 11 π2ρm(~vs)5(~vF)3 "J4 BG (cid:18)vF(cid:19) IV. ELECTRICAL RESISTIVITY T 2 π2 1 (T/Θ )+ (T/Θ ) (26) 4 BG 6 BG ×(cid:18)ΘBG(cid:19) (cid:26) 3 J 2J (cid:27)# As discussed in Section III, the electrical resistivity is obtained by setting T = 0 in Eq. (10), to obtain Let us now consider the coefficient P22, which after sim- ρ = L−EE1. By direct su∇bstitution of the coefficients and ilar manipulations as before, can be expressed by the currentsobtainedinSectionIII,weobtainafterEqs.(19) integral form that the leading term which defines the electrical resis- tivity due to electron-phonon scattering processes is (k T)−1A2 π π P = B dθ dϕ′ dǫǫ dǫ′ǫ′ 22 (π~vF)4 Z−π Z−π Z Z ρ (T)=ρ T 4 (Θ /T)+ 4π2 vs 2 ×[(q·uˆ)(ǫ−EF)+~ωq(k′·uˆ)]2Pkkq′ (27) e−ph 0"(cid:18)ΘBG(cid:19) J4 BG 3 (cid:18)vF(cid:19) 6 After expanding the square in Eq.(27), it is shown in T 3 (Θ /T)+ (Θ /T) (31) uAˆ)p2pe=ndπixq2A, thdϕat′(kth′e uˆa)ng=ulaπrki′2n,teagnradls gdiϕve′(,q duˆϕ)′((kq′ · ×(cid:18)ΘBG(cid:19) (cid:26)J4 BG 2π2J6 BG (cid:27)# uˆ) = (π/2)[q2 +k′2 ·k2]. Then, after similar aR·lgebraic· Here, ρ = 8D2k /(e2ρ v v2) is a coefficient with R − R 0 F m s F procedures as described previously, Eq.(27) reduces to dimensions of a two-dimensional resistivity (Ω), and Θ = 2(v /v )E /k is the Bloch-Gru¨neisen temper- ~D2A2E3(k T)4 4π2 BG s F F B P = F B 1+ (1+12 ature. The functions p(ΘBG/T) and their asymptotic 22 π2ρm(~vs)3(~vF)5 (cid:20)(cid:26) 3 properties are definedJin Appendix C. Eq. (31) depends v 2 T 2 2 on the dimensionless parameter T/ΘBG, a feature pre- s (T/Θ ) 1 6 viously obtained2 within the relaxationtime approxima- 4 BG ×(cid:18)vF(cid:19) !)(cid:18)ΘBG(cid:19) J − 3{ − tion. The functions p(ΘBG/T), defined in Appendix J v 2 T 2 16π2 v 2 C, have the property p( ) = p!ζ(p), with ζ(p) the s s J ∞ 6(T/ΘBG)+ Riemann zeta function. Therefore, the low-temperature ×(cid:18)vF(cid:19) )(cid:18)ΘBG(cid:19) J 3 (cid:18)vF(cid:19) (T ΘBG) behavior of the electron-phonon contribu- ≪ 4 2 tion to the electrical resistivity is T v 1 s (T/Θ )+ 1 ×(cid:18)ΘBG(cid:19) ((cid:18)vF(cid:19) J6 BG 10π2 { ρe−ph(T) =ρ04!ζ(4)(T/ΘBG)4+o(T/ΘBG)6(32) v 2 −2(cid:18)vFs(cid:19) )J8(T/ΘBG))# (28) Awet fihnigdh ttheempaesryamtuprteoti(cTlim≫it ρΘe−BpGh)(,T)fro=m E{1q/.(331−) 6 1/10... ρ T/Θ = (π/8)ρ T/Θ . This behavior is cientsandcurrentscalculatedinSectionIII,andneglect- 0 BG 0 BG in agree}ment with experiments1 and earlier theoretical ing terms proportional to (v /v )2 10−4, we obtain s F ∼ resultsbasedontherelaxationtimeapproximation2,thus that the leading contribution to the thermal resistivity supporting the application of the variational method in due to electron-phonon scattering is given by this case. In Fig.(1) we compare the theoretical predic- tion with the experimental data reported in Ref. (1). [κ (T)]−1 = ρ0 (T/Θ )4+ 3 (T/Θ )2 The experimental data are fitted to the expression e−ph T BG 4π2 BG L0 (cid:26)(cid:20) (cid:21) ρ(T)=ρ(T=0)+ρe−ph(T), (33) (Θ /T) 1 (T/Θ )4 (Θ /T) (34) ×J4 BG − 2π2 BG J6 BG withρe−ph(T)definedbyourtheoryEq.(31),andρ(T=0) (cid:27) being the value of the resistivity at zero temperature. Here, =π2k2/(3e2)istheLorenznumberforthefree This parameter can be interpreted as the (temperature L0 B electron gas. At very low temperatures (T Θ ), Eq. independent) contribution of impurity scattering to the ≪ BG (34) predicts the asymptotic limit total electrical resistivity. The impurity contribution to the electrical resistivity has been discussed5,11, and 3 ρ 3 T T it is a sample-dependent property, since it is propor- [κ (T)]−1 = 0 4!ζ(4) +o (35) e−ph Θ 4π2 Θ Θ tional to the surface concentration of impurities nimp. L0 BG BG (cid:18) BG(cid:19) It has been shown5,11 that, for a short range scatter- Itisremarkablethat,accordingtoEq. (35),thephonon- ing potential V (r) = V δ(r), the in-plane resistivity scatt 0 limited thermal resistivity is linear at very low tempera- is ρ n V2, whereas for a long range Coulomb potesncattital1∝1 ρsimcaptt0∝ nimpn−1, this last case showing an ticuarlelys,pinrecdoincttreadstawndithextpheertimypeinctaall∼lyTob2sbeerhveadvioinrtnhoeromreatl-, inversedependencyonthecarrierconcentrationnaswell. three-dimensionalmetals3,4. Thethermalresistivitycon- Therefore,itisexpectedthatthe zerotemperatureresis- tribution due to elastic scatteringwith Coulombimpuri- tivity exhibits adependence onthe carrierconcentration tieshasbeencalculated5,andcorrespondsto[κ ]−1 = of the form ρ = a + bn−1, with the value of the scatt (T=0) 3hT−1u2/(2π2k2E2),withu2 =n Ze2/(16ǫ2ǫ2). The coefficients a and b depending on details of the sample, 0 B F 0 imp 0 total thermal resistivity can be estimated from the ex- particularlyon the concentrationanddistribution ofim- pression [κ(T)]−1 [κ (T)]−1 + [κ (T)]−1. In purities and defects. We verify that this equation is in ∼ e−ph scatt Fig.(2) inset, we compare the value of the total three- excellent agreement with the data reported in Ref. (1), dimensional thermal conductivity at room temperature, as shown in the inset of Fig. (1). We extracted values normalized by a nominal packing ”thickness” of 3.4 ˚A a=14.276(Ω)andb=266.76 1017(Ωm2),withlinear × for the graphene layer. For numerical evaluation, we as- regressioncoefficient r=0.999. sumed Z = 1 for the impurities, and from the exper- As seen in Fig. (1), Eq.(33) provides an excellent fit imental values of the zero temperature electrical resis- to the experimental data. The temperature-dependent tivity, we extracted an average impurity concentration electron-phonon contribution to the electrical resistivity n = 1.3 1015 m−1. For the experimental system is fitted with just two free parameters: The energy pa- imp × reported in Ref.(1), we estimated a dielectric constant rameter D = 23.5 0.5 eV in the deformation poten- tial, and the speed±of sound for longitudinal phonons ǫ = 3.1 representing the average between the SiO2 sub- strate and the PEO polymer electrolyte. As seen in the in graphene v = 24.0 0.6 Km/s. It is remarkable s ± inset of Fig.(2), the total normalized thermal conductiv- that both values are in excellent agreement with inde- ityatroomtemperatureisontheorderofκ 400Wm−1 pendent estimations reported in the literature, particu- K−1. We can compare this with the therm∼al conductiv- larly the speed of sound12,13. It is also relevantto notice ity due to the phonon system, where the experimental that the ratio (v /v )2 10−4 is negligibly small, and s F ∼ value14,15, in excellent agreement with a theory previ- hencetermsproportionaltothisfactorcaninpracticebe ously reported by the present author12, is about12,14,15 neglected in the analytical expressions obtained for the κ˜ = 4300 Wm−1K−1 at room temperature. These electrical resistivity, as well as in other transport coeffi- ph results suggest that, in most of the temperature range, cients discussed along this work. This also explains the the phonon contribution to the thermal conductivity in quantitative success when applying the relaxation time graphene dominates over the electronic contribution, as approximationtoelectron-phononscatteringasreported the present author has pointed out elsewhere12. inprevioustheoreticalwork2,eventhoughthisscattering mechanism is not strictly quasi-elastic. VI. THERMOPOWER V. THERMAL RESISTIVITY AsdiscussedinSectionIII,bysettingJ=0inEq.(10), we obtain that the thermopower (Seebeck coefficient) Q The thermal conductivity is obtained from the gen- is given by Q = L /L . By direct substitution of eral expression Eq. (10) by setting J = 0, which as − ET EE Eq.(19), we obtain that the leading contribution to the discussed in Section III leads to the expression κ = L L−1L L . Bydirectsubstitutionofthecoeffi- thermopower is given by the expression TE EE ET− TT 7 1+ 4π2 T 2 (Θ /T)+ vs T (T/Θ ) 2 T 2 (Θ /T) π2k2T 3 ΘBG J4 BG vF ΘBG J5 BG − ΘBG J6 BG Q(T)= B (cid:20) (cid:16) (cid:17) (cid:21) (cid:16) (cid:17) (cid:16) (cid:17) (36) −3e EF 1+ 4π2 T 2 (Θ /T) 2 T 2 (T/Θ ) 3 ΘBG J4 BG − 3 ΘBG J6 BG (cid:20) (cid:16) (cid:17) (cid:21) (cid:16) (cid:17) 0 VII. CONCLUSIONS We presented a semi-classical theory to calculate -5 the phonon-limited transport coefficients in extrinsic ) 60 graphene at high carrier densities. Our theory is based 1 -K 50 n (1014 m-2) on a variational solution of the Boltzmann equation, µV -10 40 1.5 -1K) and provides explicit analytical expressions for the var- Q ( 2300 12..81 µQ (V n (1017 m-2) iiotyu,sthtrearnmsaplorrtescisoteivffiitcyienatnsd, tshuechrmaospoewleecrt.ricWalersehsioswtievd- -15 1.36 that our analytical results are in excellent agreement 10 2.86 0 4.65 with experimental data arising from two independent 0 50 100 150 200 250 300 10.8 groups, particularly concerning the electrical resistivity1 T(K) -20 and thermopower16. 0 50 100 150 200 250 T (K) Our theoretical results suggest that, in principle, it is possible to control not only the electrical resistivity, but FIG. 3: (Color online) Thermopower at different electronic alsothe thermalresistivity andthermopowerby control- densitiesn,inunitsof1017m−2,calculatedafterEq.(36)with ling the extrinsic carrier concentration in graphene. v =24 Km/s as inferred from experimental data in Ref.(1). s Inset: Experimental data (symbols) for hole thermopower in Ref.(16), compared with Eq.(36) at the corresponding hole densities (no fitting parameters besides the value of v = 24 s VIII. ACKNOWLEDGEMENTS Km/s are used). The author wish to thank FONDECYT grant No 11100064for financial support. Therefore, we find that the thermopower at very low temperatures (T Θ ) becomes BG ≪ Appendix A: Evaluation of angular integrals π2k2T v 5!ζ(5) T Q(T) B 1+ s +o(T3)(37) We consider the elementary integralsoverangularori- ∼−3e E v 4!ζ(4) Θ F (cid:20) F (cid:18) BG(cid:19)(cid:21) entations of the wave-vectorsused in the main text. Let us define the unit vector uˆ =(cos(α),sin(α)), with α an arbitrary (but fixed) direction along the thermal gradi- This result is in agreement with experimental data re- ent or the applied electric field. As in the main text, ported in the literature16. In particular, it has been we define the scattering angle θ by kˆ kˆ′ = cos(θ), observed16 that the thermopower in graphene depends with kˆ, kˆ′ the directions of the initial a·nd final wave- on the carrier density as n−1/2, and shows a linear in T vector after the electron-phonon scattering event takes dependence at very low temperatures. Since the Fermi level is E = ~v √πn, it is clear that Eq.(36) obtained place. Then, we write kˆ′ = (cos(ϕ′),sin(ϕ′)) and kˆ = F F (cos(ϕ′ +θ),sin(ϕ′ +θ)), respectively. Using this nota- from our theory correctly reproduces this feature. The tion, we consider the following integrals: contribution to the total thermopower due to phonon- drag effects has been discussed in Ref.(17), where it is shown that it displays a T3 dependence at low temper- +πdϕ′[kˆ uˆ]2 atures. Clearly then, the phonon drag is a negligible · Z−π contribution to the total thermopower as compared to +π thediffusioncomponentcalculatedinthiswork,whichin = dϕ′[cos(α)cos(ϕ′+θ)+sin(α)sin(ϕ′+θ)]2 agreement with the experimental data16 reproduces the Z−π correctlinearinT dependence atverylowtemperatures. =π[(cos(α))2+(sin(α))2]=π (A1) 8 +πdϕ′[kˆ′ uˆ]2 ~ω/(kBT), we use the identity · Z−π +π = dϕ′[cos(α)cos(ϕ′)+sin(α)sin(ϕ′)]2 +∞ F(x) dx= 1 e−z −1 Z−π (1+ex)(1+e−x−z) − =π[(cos(α))2+(sin(α))2]=π (A2) Z−∞ +∞ (cid:0) ∂f0 (cid:1) [G(x) G(x z)] dx (B1) × − − − ∂x Z−∞ (cid:18) (cid:19) +π dϕ′[kˆ uˆ][kˆ′ uˆ] · · Z−π with G(x) = xF(x′)dx′. This identity is straightfor- +π 0 = dϕ′[cos(α)cos(ϕ′+θ)+sin(α)sin(ϕ′+θ)] wardtoproveusingintegrationbyparts,andnotingthat Z−π f0(x)−f0(x+Rz)=(1−e−z)(1+ex)−1(1+e−x−z)−1. [cos(α)cos(ϕ′)+sin(α)sin(ϕ′)] The right-hand side of the equation can be further × evaluatedusingtheSommerfeldexpansionfortheFermi- +π = dϕ′ (cos(α))2cos(ϕ′+θ)cos(ϕ′) Dirac distribution function at low temperatures, leading Z−π to the formula +(sin(α))2s(cid:2)in(ϕ′+θ)sin(ϕ′)+sin(α)cos(α) ×sin(2ϕ′+θ)]=πcos(θ) (A3) +∞ F(x) dx Directapplicationofthese basicidentities leadsto the Z−∞ (1+ex)(1+e−x−z) following results =(1 e−z)−1 −zF(x)dx+ π2[F′(0) F′( z)] − − 6 − − +π +π (cid:26) Z0 (cid:27) dϕ′[k uˆ][k′ uˆ]=kk′ dϕ′[kˆ uˆ][kˆ′ uˆ] (B2) · · · · Z−π Z−π =πkk′cos(θ) (A4) which is valid up to o((k T/E )4). B F +π +π dϕ′[q uˆ]2 = dϕ′[(k′ k) uˆ]2 · − · Z−π Z−π =π(k2+k′2 kk′cos(θ))=π(k′ k)2 − − =πq2 (A5) Appendix C: Definition and properties of the functions J (z) p +π +π dϕ′[q uˆ][k′ uˆ]= dϕ′[(k′ k) uˆ][k′ uˆ] We have defined the functions · · − · · Z−π Z−π q2 k′2 k2 =π(k′2 2kk′cos(θ))=π + − (A6) − 2 2 z xp 1 (x/z)2 (cid:18) (cid:19) (z)= − dx (C1) Jp (1 e−x)(ex 1) Z0 −p − Appendix B: Evaluation of energy integrals at low temperature The asymptotic limit of these functions, for p > 0 an integer, is given by ( ) = p!ζ(p), with ζ(p) the Rie- p J ∞ In several calculations throughout this work, after mannZetafunction. Inparticular,forp=4and6ithas the change of variables x = (ǫ E )/(k T) and z = the values ζ(4)=π4/90 and ζ(6)=π6/945. F B − 1 D. K. Efetov and P. Kim, Phys. Rev. Lett. 105, 256805 Phys. Rev.B 76, 073412 (2007). (2010). 6 N. M. R.Peres, Rev.Mod. Phys.82, 2673 (2010). 2 E. H. Hwangand S. DasSarma, Phys.Rev.B 77, 115449 7 A. H.Castro, F. Guinea, N.M. R.Peres, K. S. Novoselov (2008). and A. K. Geim, Rev.Mod. Phys.81, 109 (2009). 3 J. M. Ziman, Electrons and Phonons (University Press, 8 C. W. Beenakker, Rev.Mod. Phys.80, 1337 (2008). 1960). 9 M. S. Fuhrer,Physics 3, 106 (2010). 4 O. Madelung, Introduction to Solid-State Theory 10 L. Onsager, Phys.Rev.37, 405 (1931). (Springer, 1978). 11 T. Stauber, N. M. R. Peres and F. Guinea, Phys. Rev. B 5 N.M.R.Peres,J.M.B.LopesdosSantosandT.Stauber, 76, 205423 (2007). 9 12 E. Mun˜oz, J.Lu and B. I.Yakobson,NanoLett.10, 1652 15 A. A. Balandin, S. Ghosh, I. Calizo, D. Teweldebrhan, F. (2010). Miao and C. N.Lau, Nano Lett. 8, 902 (2008). 13 R. Saito, G. Dresselhaus and M. S. Dresselhaus Physical 16 P. Wei, W. Bao, Y. Pu, C. N. Lau and J. Shi Phys. Rev. Properties of Carbon Nanotubes (Icp,1998). Lett. 102, 166808 (2009). 14 S. Ghosh, I. Calizo, D. Teweldebrhan, E. P. Pokatilov, D. 17 S. S.Kubakaddi,Phys.Rev. B 79, 075417 (2009). L.Nika,A.A.Balandin, W.Bao, F.Miao andC.N.Lau, Appl.Phys. Lett. 92, 151911 (2008).