A brief review of thermal transport in mesoscopic systems from nonequilibrium Green's function approach
aa r X i v : . [ c ond - m a t . m e s - h a ll ] F e b A brief review of thermal transport in mesoscopic systems from nonequilibriumGreen’s function approach
Zhizhou Yu, Guohuan Xiong, and Lifa Zhang ∗ NNU-SULI Thermal Energy Research Center (NSTER) & Centerfor Quantum Transport and Thermal Energy Science (CQTES),School of Physics and Technology, Nanjing Normal University, Nanjing 210023, China
With the rapidly increasing integration density and power density in nanoscale electronic devices,the thermal management concerning heat generation and energy harvesting becomes quite crucial.Since phonon is the major heat carrier in semiconductors, thermal transport due to phonons inmesoscopic systems has attracted much attention. In quantum transport studies, the nonequilibriumGreen’s function (NEGF) method is a versatile and powerful tool that has been developed for severaldecades. In this review, we will discuss theoretical investigations of thermal transport using theNEGF approach from two aspects. For the aspect of phonon transport, the phonon NEGF methodis briefly introduced and its applications on thermal transport in mesoscopic systems including one-dimensional atomic chains, multi-terminal systems, and transient phonon transport are discussed.For the aspect of thermoelectric transport, the caloritronic effects in which the charge, spin, andvalley degrees of freedom are manipulated by the temperature gradient are discussed. The time-dependent thermoelectric behavior is also presented in the transient regime within the partitionedscheme based on the NEGF method.
I. INTRODUCTION
As transistor gate lengths are scaled down into the10-nm regime with the rapid development of nanotech-nology, millions of transistors are fabricated within asquare millimeter in the integrated circuit chip . Withthe increasing transistor density in chips, the power den-sity raises rapidly, which becomes the roadblock for thecontinued miniaturization of integrated circuits since theenhanced chip temperature prevents the reliable perfor-mance of integrated circuits. In order to design next-generation devices with low energy consumption, it iscrucial to study the thermal transport in nanostruc-tures to understand heat generation and dissipation.Recently, numerous researchers have proposed varioustheoretical models to study the fundamental physics inthermal transport and carried out experiments on low-dimensional nanomaterials to show their potential appli-cations in thermal engineering .Phonon, the physical quasiparticle representing themechanical vibrations, is responsible for the transmis-sion of heat in solids. Understanding and controllingthe transport properties of phonons provide opportuni-ties to reduce heat consumption and utilize waste heat.Various prototypical phononic devices such as thermaldiodes , thermal transistors , thermal logic gates , andthermal memories have been proposed to manipulatethe heat flow at the nanoscale. Recently, the chirality ofphonons has been observed experimentally in monolayertungsten diselenide . The discovery of chiral phononshas received wide attention in emerging fields such asvalleytronics and topological states . Therefore, ex-ploring the mechanisms of phonon transport and scatter-ing in nanoscale phononic devices is of great importancefor artificially tuning thermal transport properties for fu-ture heat management in electronic devices and specificapplications in phononic devices. Apart from phonon transport, the thermoelectric ef-fect which describes a direct conversion from heat energyto electric energy and vice versa, is another major con-cern in the field of thermal transport due to its poten-tial applications in harvesting and recovering heat. Theperformance of thermoelectric materials at a certain tem-perature is evaluated by the dimensionless figure of merit( ZT ). The big challenge lying behind the thermoelectrictechnology is the improvement of ZT value of thermo-electric materials, namely, simultaneous enhancement inthe electrical conductivity and reduction in the latticethermal conductivity . In the past decades, the ther-moelectric behavior of a series of low-dimensional mate-rials has been theoretically predicted and experimentallystudied, which exhibits huge potential in the applicationof high-performance thermoelectric devices . How-ever, it is still an open question and a long way to searchfor better thermoelectric materials and further improvethe ZT value.The method of nonequilibrium Green’s function(NEGF) is a versatile and powerful tool to study bothelectronic and phononic transport properties in nanoscalematerials. The NEGF method was used to investigatequantum electric transport by Caroli et al. for the firsttime in 1971 . An explicit formula for the transmissioncoefficient and tunneling current was derived in termsof the Green’s function. A Landauer formula for thecurrent through an interacting electron region was de-rived by Meir and Wingreen, which provided a modernframework to study the electronic transport in meso-scopic systems . The general formula of time-dependentelectric current through the interacting and noninteract-ing mesoscopic systems was derived using the KeldyshNEGF technique . Besides the electronic transport, theNEGF method was used to treat the phonon transport insolid junctions by Wang et al. and the formula of ther-mal current due to atomic vibrations was presented interms of Green’s function . Within the NEGF ap-proach, many-body effects in quantum transport suchas electron-phonon and electron-electron interactions canbe included through self-energies without deviating theframework . The NEGF method was also combinedwith the density functional theory (DFT) which is an art-of-the-state technique for modeling and predicting theelectronic transport properties of nanomaterials .In this review, we aim to give a brief summary of theo-retical studies on thermal transport including the phononand thermoelectric transport in mesoscopic systems byusing the NEGF method. In Sec. II, we first introducethe phonon NEGF method and its applications on ther-mal transport. The interfacial thermal transport in one-dimensional atomic chains, phonon transport in multi-terminal systems, and time-dependent phonon transportin the transient regime are discussed. In Sec. III, thebasic concepts of thermoelectricity are introduced. Thedc thermoelectric transport and its application on spinand valley caloritronics are discussed within the lin-ear response theory. The time-dependent thermoelectrictransport in the transient regime within the partitionedscheme was also presented. Finally, a brief conclusionand outlook are given in Sec. IV. II. PHONON TRANSPORTA. NEGF method for phonon transport
Various methods have been used to study the phonontransport, such as molecular dynamics (MD) andBoltzmann transport equation (BTE) method . TheMD method can incorporate nonlinearity. However, itis only valid at high temperatures and becomes not ac-curate at low temperatures due to its classical nature.The BTE method is usually used to study the thermalconductivities for bulk materials and can not be used forsystems without translational invariance. For the meso-scopic system in which quantum effects dominate thephonon transport, NEGF is an effective approach in awhole diffusive to ballistic regime . In this sec-tion, we first give a quick review of the NEGF techniquein phononic systems.
1. phonon current
We consider a nonconducting solid that only the vibra-tional degrees of freedom are treated. The Hamiltonianis given by , H = X α = L,C,R H α + u † L V LC u C + u † C V CR u R , (1)where L, C, R denotes the left lead, central region, andright lead, respectively. H α = 12 ˙ u † α ˙ u α + 12 u † α K α u α , (2) where u α is the column vector consisting of all displace-ment variables in region α and ˙ u α is the correspondingconjugate momentum. K α is the spring constant matrix. V CL = V † LC and V CR = V † RC are the coupling matricesof the central region to the left and right leads, respec-tively. The dynamic matrix for a full linear system canbe written as, K = K L V LC V CL K C V CR V RC K R . (3)The phonon current flow from the left lead to the cen-tral region can be defined as , J L = −h ˙ H L ( t ) i . (4)By using the Heisenberg equation of motion, we can ob-tain J L = h ˙ u † L ( t ) V LC u C ( t ) i . (5)By defining the following lesser Green’s function , G 1] is the Bose-Einsteindistribution function in lead α ( k B = 1 for simplicity)and Ξ( ω ) = Tr( G r Λ L G a Λ R ) , (18)is the phonon transmission coefficient in the form of theCaroli formula . More details of the basic definition andproperties of phonon NEGF can be found in Refs. 39 and40.We define the phonon thermal conductance as , κ ph = lim ∆ T → J ∆ T , (19)where ∆ T is the temperature difference of two leads. Forballistic transport, the phonon conductance can be ex-pressed in the form of Landauer-like formula κ ph = Z + ∞ dω π ω Ξ( ω ) ∂n∂T . (20) 2. Nonlinear systems In the following, we discuss the quantum self-consistentmean-field theory based on the NEGF method to dealwith nonlinear thermal transport. We introduce thequartic interaction term into the Hamiltonian as an ex-ample, which can be given by H n = 14 X ijkl T ijkl u C,i u C,k u C,j u C,l . (21) We can also handle the cubic interaction term for thethermal transport. By applying the equation of motion,the Green’s function with the nonlinearity can be writtenas , ∂ ∂τ G im ( τ, τ ) + X j K C,ij G jm ( τ, τ )+ X jkl T ijkl G jklm ( τ, τ, τ, τ )= − δ ( τ − τ ) δ im − X j Z dτ Π ij ( τ, τ ) G jm ( τ , τ ) , (22)where G ( τ , τ , τ , τ ) = − i h T c u ( τ ) u ( τ ) u ( τ ) u ( τ ) i with T c the time-order operator is the four-point Green’s func-tion. Within the mean-field approximation, the four-point Green’s function can be represented by the two-point Green’s function , − iG ( τ , τ , τ , τ ) ≈ G ( τ , τ ) G ( τ , τ ) + G ( τ , τ ) G ( τ , τ )+ G ( τ , τ ) G ( τ , τ ) . (23)Then we can obtain ∂ ∂τ G im ( τ, τ ) + X j K C,ij G jm ( τ, τ )+3 i X jkl T ijkl G kl (0) G jm ( τ, τ )= − δ ( τ − τ ) δ im − X j Z dτ Π ij ( τ, τ ) G jm ( τ , τ ) . (24)Therefore, we can account the nonlinearity by the follow-ing self-energy,Π n,ij = 3 i X kl T ijkl G kl (0) = 3 X kl T ijkl h u k u l i , (25)where h u k u l i = i Z ∞ dω π G < ( ω ) . (26)We note that this nonlinear self-energy is real and it onlyshifts the frequencies of phonon modes.By introducing the nonlinear self-energy, the retardedGreen’s function with nonlinearity can be then writtenas G r = [( ω + i + ) I − K C − Π r − Π n ] − . (27)With the help of the Keldysh equation, i.e., Eq. (15), theretarded Green’s function can be solved self-consistently.Since we are considering an effectively harmonic problem,the phonon current can be still calculated from Eq. (17). 3. Electron-phonon interaction For thermal transport through the metal-semiconductor interface, energy must transfer betweenelectrons and phonons. Therefore, it is highly desirableto understand the heat dissipation for thermal transportthrough the interface with the electron-phonon interac-tion. The electron-phonon coupling in the central regioncan be described by H ephC = X i ǫ i d † i d i + X ijk M kij d † i d j u k , (28)where d † i ( d i ) is the electron creation (annihilation) op-erator and ǫ i is the electron energy level in the centralregion. M kij is the electron-phonon coupling matrix ele-ment. Since the Landauer formula of phonon current isonly applicable to quasi-ballistic transport, one needs touse the Meir-Wingree formula, i.e., Eq. (13), to calculatethe phonon current of inelastic processes, including theelectron-phonon scattering.The electron-phonon interaction is included as a per-turbation. The full retarded Green’s function withinthe electron-phonon coupling can be obtained from theDyson equation ,¯ G r = G r + G r Π reph ¯ G r , (29)where G r given in Eq. (14) is the bare phonon re-tarded Green’s function without electron-phonon interac-tion. The Keldysh equation for the system with electron-phonon coupling becomes¯ G < = ¯ G r (Π < + Π In thermal transport, the interfacial thermal scatteringbecomes extremely important as the dimension of ther-mal devices shrinks to the nanoscale. In low dimensionalsystem, it was found that the interfaces can dramaticallyaffect the thermal transport . In recent years, in-terfacial thermal transport has been extensively studiedby both classical and quantum approaches. To study theinterfacial thermal transport, the most widely used mod-els are the acoustic mismatch model and the diffusemismatch model . However, both models are lack ofaccuracy in calculating the interfacial thermal resistancesince they neglect the atomic details of actual interfacestructures. The NEGF approach, which is a powerfulmethod to treat nonequilibrium and interacting systems,has been extensively applied to study interfacial ther-mal transport. Moreover, the NEGF method can offer astraightforward way to treat nonlinear systems.One-dimensional atomic chain model has been exten-sively used to study the interfacial thermal transport,which can provide fundamental physical pictures forpractical thermal devices. The one-dimensional atomicchain consists of two semi-infinite leads and a central re-gion, as shown in Fig. 1. The left and right leads are inequilibrium at different temperatures T L and T R , respec-tively. The central region is coupled with the left andright leads by harmonic springs with constant strength k and k , respectively. The left lead, central region,and right lead are all harmonic chains with the springconstant and mass k , m , k , m , k , m , respectively.The total Hamiltonian of the one-dimensional atomicchain can be given by H = X α = L,C,R H α + 12 k ( x L, − x C, ) + 12 k ( x C,N C − x R, ) . (34)Here, H α = N α X i =1 m α ˙ x α,i + N α − X i =1 k α ( x α,i − x α,i +1 ) , (35)where x α,i is the relative displacement of i th atom in part α . N α is the number of atoms in part α . Note that forthe semi-infinite leads, N L and N R are infinite.The simplest model of the one-dimensional atomicchain is the single-junction case, namely, two semi-infinite leads are directly connected by a spring with aconstant strength k . For the Hamiltonian of the single-junction case, we can set k = 0, N C = 0, and replace x C, by x R, in Eq. (34). Within the NEGF approach,the transmission coefficient can be given by ,Ξ( ω ) = B C C | A A − B | . (36) k m , k k m , k m , k ...N C FIG. 1. Schematic of the one-dimensional atomic chainmodel. The central region is coupled with the left and rightleads by harmonic springs with constant strength k and k ,respectively. The left lead, central region, and right lead areall harmonic chains with the spring constant and mass k , m , k , m , k , m , respectively. INTERFACIAL THERMAL TRANSPORT IN ATOMIC JUNCTIONS PHYSICAL REVIEW B , 064303 (2011)FIG. 1. (Color online) Schematic of the 1D atomic chain model.The size of the center part is 8. The left and right regions aretwo semi-infinite harmonic atomic chains at different temperatures,and . The three parts are coupled by harmonic springs with and , all of which are harmonic chains withspring constant as and and , and andrespectively. resulting in the transmission coefficient through the centerpart (1 − | )(1 − | (16)here ij and are determined by Eqs. (12) and ( ), re-spectively; is the number of atoms in the center atomicchain. From this expression, we can find that the transmissioncoefficient oscillates with frequency, and is between theenvelope lines of maximum and minimum transmission,are max (1 − | )(1 − | (1 − | 23 21 for constructive interference and min (1 − | )(1 (1 + | 23 21 for destructive interference. IV. RESULTS AND DISCUSSIONSA. Thermal transport in 1D one-junction chains In Sec. III, we have derived the analytical expressions forthe phonon transmission coefficient for the point-junction andextended-junction (two-point-junction) cases Eqs. (11), (12),and (16) by using the scattering boundary method. Usingthese analytical expressions, we analyze the role of variousparameters on the thermal transport in one- and two- pointjunctions.Figure shows the transmission coefficient as a function offrequency for a different interface spring constant for thepoint-junction model. The maximum frequency at which thetransmission coefficient is above zero is equal to the minimumof 2 /m and 2 /m . In Fig. 2(a), the two semi-infiniteatomic chains have the same mass and spring constant. Whenthe interface coupling equals to that of the chains, thetransmission is equal to one in the whole frequency domain,because of the homogeneity of the chain structure. If increases or decreases, the transmission coefficient decreases.If we set /m /m , the transmission coefficient exhibitssimilar behavior, the only difference is that the transmissioncoefficient changes to the value obtained by Eq. (15). InFig. 2(b), the two semi-infinite atomic chains have differentmasses and spring constants. The transmission decreases withincreased frequency for all the coupling values . Also,it appears that for a given frequency the transmission ismaximized for a value residing between and . FromEq. (11) and Eq. (12), 0, if 0; and ] hasdefinite value 1 − | 1) (1(1 (1 , if = ∞ 2. (Color online) Transmission coefficient vs frequencyfor different interface coupling in one-junction chains.Transmission in one junction connected by the same semi-infiniteatomic chains with , m 0; solid, dashed,and dash-dotted lines correspond to 1, 0.5, 1.0and 2.0, respectively. (b) Transmission in one junction connectedby two different semi-infinite atomic chains with , m, k 5, 1.0, 1.5, 3.0 and 8.0,vely. The maximum transmission concept results in themaximum junction conductance as shown in Fig. . Withthe increasing of , we find that the conductance will firstincrease, then arrive at a maximum value, and then slightlydecrease and at last it will tend to a constant. We find that themaximum transmission or conductance occurs at given by 12 12 (17)That is, when the coupling spring stiffness is equal to theharmonic average of spring connecting atoms in the two semi-infinite chains. In Fig. , we show the thermal conductancevs the ratio of and . For the two semi-infinite chainssame mass = , the maximum conductance FIG. 3. (Color online) Thermal conductance vs interface coupling in the point-junction model. Here , m FIG. 2. Thermal conductance σ as a function of interfacecoupling k in the single-junction model. Here k = 1 . m = 1 . 0. Reproduced with permission from Ref. [51]. Here, A i = ω − k i m i (1 − λ i ) − k m i , (37) B = k √ m m , (38) C i = ωm i q k i m i − ω m i , (39)where λ i = e iq i a i with q i the wave vector and a i theinteratomic spacing.Figure 2 presents the thermal conductance as a func-tion of interfacial coupling k in the single-junctionmodel. It is found that the thermal conductance ini-tially increases with the increasing interfacial coupling k and reaches a maximum value. It then decreasesslightly and finally approaches a constant value. Zhanget al. found that the maximum thermal conductanceoccurs when the interface spring equals the harmonic av-erage of the spring constants in two semi-infinite leads,namely, k satisfies k = k ,m = 2 k k k + k . (40) Besides, the effect of impurity mass and mechanical ad-hesion on phonon transport was investigated by Salton-stall et al. by introducing an impurity mass and variablebonding into the single-junction model . For the case ofinterface mass, it is found the maximum transmission oc-curs when the interface mass equals the arithmetic meanof the mass on either side of the interface. For the case ofthe interface spring, one can maximize the transmissionwhen the interface spring is set to the harmonic mean ofthe spring constants in two semi-infinite leads, namely, k ,m .The single-junction model can be extended to the two-junction model which involves a central part. In the two-junction model, the transmission wave is scattered by twoboundaries, which results in multiple reflections. Thetransmission behavior can be considered as the combina-tion of the transmission in the single-junction model andthe oscillatory behavior due to the multiple scattering.For the two-junction model with homogenous mass andcoupling in the central part, it is found that the phonontransmission oscillates with frequency in the envelopelines of minimum and maximum transmission which canbe determined by the single-junction model . The in-terfacial thermal conductance of two-junction model forvarious mass-graded and coupling-graded materials wasinvestigated by Xiong et al. The optimized homoge-nous coupler , the arithmetic mass-graded and coupling-graded coupler, the geometric mass-graded and coupling-graded coupler, and the coupler with both geometricgraded mass and coupling were studied. Relative tothe optimized homogenous couplers, the mass-graded orcoupling-graded structures were found to be applicable toimprove the interfacial thermal conductance of two leadmaterials with both mismatched impedance and mis-matched cutoff frequencies . For the couplers with bothgeometric graded mass and geometric graded coupling,the interfacial thermal conductance can be maximum en-hanced nearly up to sixfold compared to the optimizedhomogenous case. They also found that the interfacialthermal conductance decreases with the increasing cut-off frequency ratios for all six cases due to the increasingmismatch of the cutoff frequency .In the above, we discuss the interfacial thermal trans-port in one-dimensional atomic chains with only linearcoupling interactions. However, the nonlinear effect atthe interface is another crucial issue for further under-standing the fundamental physical mechanism of phonontransport. Zhang et al. introduced a fourth-order non-linear interaction into the one-dimensional atomic chainmodel and studied the thermal transport through a solid-solid interface . By using the quantum self-consistentmean-field theory based on the NEGF method, theyfound that the nonlinear interaction λ plays a role tomodulate the interfacial linear coupling k and the ef-fective interfacial coupling can be given by k ,eff = k + 3 λ (cid:18) h u i m − h u u i√ m m + h u i m (cid:19) . (41)It was also found that in the weak-interfacial-couplingregime, the interfacial thermal transport is enhanced bythe nonlinearity, while the enhancement vanishes in thestrong-interfacial-coupling regime.The phonon transport with the weak electron-phononinteraction was also studied in one-dimensional atomicchains. Based on the NEGF method. L¨u etal. derived the electrical and energy current of thecoupled electron-phonon system by introducing theelectron-phonon interaction within the adiabatic Born-Oppenheimer approximation . They showed that theself-consistent Born approximation fulfills the electricaland energy current conservation. Zhang et al. studiedthe thermal conductance and thermal rectification acrossthe metal-insulator interface with electron-phonon inter-action by using the NEGF method . They found thethermal conductance has a nonmonotonic behavior asa function of the average temperature of both phononleads. Moreover, by considering the same temperatureof left and right phonon leads and setting k = 0, thephonon contribution in metal was excluded to avoid di-vergence. Figure 3(a) presents that the thermal rectifi-cation changes its sign with the increase of temperatureat a relatively larger electron-phonon interaction. Whilethe thermal rectification remains negative at a very weakelectron-phonon interaction, as shown in Fig. 3(b). Thereverse of thermal rectification can be explained by therelation of thermal currents in the forward and backwarddirections. At a weak electron-phonon interaction, theforward thermal current is smaller than the backwardone, which results in the negative thermal rectification,as presented in Fig. 3(c). When the electron-phonon in-teraction is strong, the forward thermal current becomeslarger than the backward one since more electrons faraway from the Fermi surface contribute to the thermalenergy, leading to the positive thermal rectification.Besides, the interfacial thermal transport was stud-ied across anharmonic systems via the one-dimensionalatomic chain model. He et al. developed a quantumself-consistent approach to renormalize the anharmonicHamiltonian to an effective harmonic one, which was usedto calculate the interfacial phonon transport within theframework of NEGF method . Fang et al. studied theanharmonic phonon transport across interfaces in non-linear one-dimensional lattice chains based on the equi-librium MD simulation. An efficient method to calcu-late the frequency-dependent anharmonic phonon trans-mission coefficients was proposed based on the linear re-sponse theory .Recently, interfacial phonon transports have beenextensively studied across the interfaces basedon various nanostructures such as single-moleculejunctions , self-assembled monolayer interfaces ,one-dimensional nanotube junctions , and two-dimensional heterojunctions . These studies onthermal conductance through actual interfaces confirmsthe general rules obtained from the NEGF method in theone-dimensional atomic chains. Hu et al. investigated J. Phys.: Condens. Matter L Zhang et al e 5. of the metal–insulator interface versus temperature at ep versus the electron–phonon interaction ep at different temperatures. (c) The forward thermal currentversus ep at available density of states (LADOS) at = | = | )ω ) (dashed line). (d1) ((d2)) andto the forward (backward) transport at 25 and 1.0, respectively. (e1)–(e4) The forward (solid lines) and (ε) For all the curves,; from (b) to (e) 4. 5. Discussions Green’s functions in the calculation A, we have calculated the self-energies of thevergence we obtain all the Green’sof in equations ( ) and ( ), thenwe can calculate the heat currents. If we do not consider thein the metal part, we should choose theto be a very tiny positive value, thus the heatto the left part will be very small and negligible; aswn in figure , if 10 flowing intois almost zero. If 0, the iterationto calculate the self-energies for the EPI will be divergent,but if 10 , the iteration is convergent. Thus herea role as a very small onsite potentialis similar as the tiny onsite ineffect [36]. We plot the bares functions s functionsin figure . Without the EPI, the bare Green’s functions ofto the originas shown in figure but the symmetry forbe broken and change greatly if we turn on the EPI, asbe seen in figure However, the Green’s functions fordo not change much. role of Fermi energy of electrons We have set the onsite 0 in the main text, thus the Fermigy of the electron system is 0. However, we can, as shown in figureof the thermal conductance andhave maximum value shifts aHowever, the curves of the thermal conductance andhave the same properties as in the case at 0.ofreversal of rectification still holds for theof nonzero Fermi energy.We have used a simplified one-dimensional model toFor some quasi-one-dimensional metal–insulator interfacebe recast asa one-dimensional model, and then our results of thebe applied. The thermalof an interface system comes from the two materialses and the interface between them. Our study in FIG. 3. (a) Thermal rectification R of the metal-insulator in-terface as a function of temperature T for different tempera-ture gradients with a electron-phonon interaction V ep = 0 . R as a function of electron-phononinteraction V ep at different temperatures. (c) Thermal cur-rent in the forward (solid line) and backward (dashed line)transport a function of electron-phonon interaction V ep at T = 0 . 25. Reproduced with permission from Ref. [46]. the phonon transport across a self-assembled monolayerof alkanethiol molecules sandwiched between gold andsilicon substrates using the MD simulation. They foundthat the transmission coefficients exhibit strong andoscillatory dependence on frequency, which agrees withthe phonon transmission behavior in the two-junctionmodel . The interfacial thermal conductance of par-tially unzipped carbon nanotubes was studied by usingthe NEGF method . The armchair carbon nanotubewas longitudinally unzipped to obtain curved zigzaggraphene nanoribbons in its central part, as shown inFig. 4(c). In Fig. 4(a), Chen et al. presented that thethermal conductance exhibits a linear dependence onthe width of the unzipped graphene nanoribbon region.This can be explained by the enhanced phonon transportchannels of carbon nanotubes with a wider width ofthe unzipped region from the phonon transmission ofpartially unzipped carbon nanotubes (PUCNTs) shownin Fig. 4(b). Such a linear behavior of the thermal FIG. 4. (a) Scaled thermal conductance at 300 K of m -PUCNT( n , n ) as a function of the scaled width m/ n . Thescaled thermal conductance is defined as the ration of ther-mal conductance of m -PUCNT( n , n ) to the thermal conduc-tance of a pristine ( n , n ) carbon nanotube. m is the num-ber of zigzag carbon atom chains in the unzipped part. (b)Phonon transmissions of m -PUCNT(10,10) as a function ofphonon frequency. The highest (lowest) value is representedby red (blue) color. (c) Phonon local density of states of a7-PUCNT(6,6) at ω = 1000 cm − . Reproduced with permis-sion from Ref. [63]. conductance to the width of the unzipped graphenenanoribbon region implies that the key factor determinedthe phonon conduction is the width of the central part. C. Multi-lead systems In Eqs. (13) and (17), the thermal currents of systemswith two leads are derived. These formulas can be usedin the same form for systems with multiple leads whenthere are no interactions between leads. Similar to thetheory of B¨uttiker on the electronic transport in systemswith multiple leads, the thermal current flowing out the α lead can be given by J α = Z + ∞ dω π ω X β = α Ξ βα ( ω )( n α − n β ) , (42)where Ξ βα ( ω ) = Tr( G r Λ α G a Λ β ) , (43)is the transmission coefficient between the α and β leads.The ballistic thermal transport in three-terminal junc-tions was studied by Zhang et al. in which the thermalcurrent of the third lead is set to be zero by adjusting itsbath temperature . The thermal rectification is foundin asymmetric three-terminal junctions due to the inco-herent phonon scattering from the control lead. By in- ference versus magnetic field at temperature45 K. The hexagons and squares correspond to central regions for thea nearest-neighbor coupling. The red dottedto a linear fit from 0 to 40 T. The size of the center region fors is 9 6, the same as the inset (c) in figure 3. Results In the following calculation, we assume a lattice constant Å, and the force constanteV 4. The ratio of the longitudinal and transverse soundto be /v 2. Then the speed of sound for longitudinal acousticis about 4000 m s . As mentioned above, is estimated to be about 3 10 Hz10 eV at 1 T. We set all the couplings between the leads and central region theregion have the same spring constants for simplicity.10), we find that the relative Hall temperature difference is an oddof magnetic field. The Onsager relation, αβ βα , always holds due to theof the conductance. Furthermore, if there is a symmetry operation h that S K S S A S = − αβ αβ . If this relation is true, then there is no PHE in the system. Theseconsistent with a different treatment for bulk systems based on the22].We discuss numerical results in the following. Figure ws the temperature differenceat temperature 45 K for the honeycomb and square latticesFor the honeycomb case, the Hall temperature is odd and linearin the magnetic field between 0 and 40 T, in that range the slope of the curve is 3 10 K Tto the experimental data in [ ]. When the magnetic field is extremely large, it New Journal of Physics FIG. 5. Hall temperature difference R as a function of mag-netic filed B at temperature T = 5 . 45 K. Reproduced withpermission from Ref. [78]. troducing the spin-phonon interaction, the thermal rec-tification can be found in symmetric three-terminal junc-tions with an external magnetic field. The ballistic ther-mal rectification effect was also studied analytically andnumerically in asymmetric three-terminal mesoscopic di-electric systems . The model of three-terminal junctionsis widely extended to study the thermal transport in var-ious two-dimensional nanomaterials . For instance,Ouyang et al. studied the phonon rectification effectof asymmetric three-terminal graphene nanojunctions .They found that the rectification efficiency is strongly de-pendent on the asymmetry of graphene nanojunctions,which can be significantly improved by increasing thewidth difference between left and right leads. Moreover,the mode-dependent phonon transport in three-terminalgraphene nanojunctions was investigated by Gu et al.based on the NEGF method and the acoustic modeswere found to contribute higher transmission coefficientsbetween the zigzag graphene nanoribbon and the thirdlead .In electronic transport, four-terminal devices havebeen extensively used to study the spin Hall effect fortwo-dimensional mesoscopic systems in which a trans-verse charge accumulation is induced by a longitudinalelectric field . Analogous to the electric Hall effect,the phonon Hall effect where a transverse heat flow in di-electrics is induced by a longitudinal temperature differ-ence has been discovered experimentally in 2005 . Usingthe NEGF approach, Zhang et al. studied the phononHall effect for paramagnetic dielectrics in four-terminalnanojunctions. Fig. 5 presents the Hall temperaturedifference for the honeycomb and square lattices withnearest-neighbor couplings under different magnetic filedat the temperature of T = 5 . 45 K. For the honeycomblattice, it is found that the Hall temperature differenceexhibits the linear relation to the magnetic field lesserthan 40 T. The fitted slope is about 3 × − K T − ,which is comparable to the experimental results . Whenthe magnetic field is extremely large, the Hall tempera-ture difference decreases slightly with the increasing mag-netic field. However, the phonon Hall effect can not beobtained in the square lattice with nearest-neighbor cou-plings due to the mirror reflection symmetry of the dy-namic matrix. Once the next-neighbor couplings is con-sidered in the square lattice, the phonon Hall effect canthen be obtained. D. Time-dependent phonon transport in thetransient regime In the past decade, most of the theoretical works onthermal transport focus on the calculation of steady-statephenomena. However, the time-dependent phonon cur-rent in the transient regime is also an important question.Recently, the transient phonon transport was studied inarbitrary harmonic systems connected to phonon bathsby abruptly turning on the coupling between leads withinthe partition scheme based on the NEGF method .Considering a single-junction one-dimensional chainmodel in which the left and right leads are initially un-coupled. Before t = 0, it is assumed that the left andright leads are in thermal equilibrium with temperature T L and T R , respectively. The coupling between the leftand right leads is suddenly switched on at t = 0 by aninterparticle harmonic potential with a spring constant k . The time-dependent phonon current in the transientregime can be expressed as J L ( t ) = k Im (cid:20) ∂G RL,< ( t , t ) ∂t (cid:21) t = t = t . (44)Here, the time-derivative of G RL,< ( t , t ) is given by ∂G RL,< ( t , t ) ∂t = − k Z t dt a G RL,r ( t , t a ) ∂G RL,< ( t a , t ) ∂t − k Z t dt a G RL,< ( t , t a ) ∂G RL,a ( t a , t ) ∂t + k Z t dt a Z t dt b G RL,r ( t , t a ) G RL,< ( t a , t b ) × ∂G RL,a ( t b , t ) ∂t + ∂G RL,< ( t , t ) ∂t , (45)where G RL,< ( t , t ) = − k Z t dt a h g R,r ( t − t a ) g L,< ( t a − t )+ g R,< ( t − t a ) g L,a ( t a − t ) i , (46)and G RL,β ( t , t ) = − k Z t dt a G RL,β ( t , t a ) G RL,β ( t a , t )+ G RL,β ( t , t ) , (47) PHYSICAL REVIEW B , 019902(E) (2011) Transient behavior of heat transport in a thermal switch [Phys. Rev. B 81, 052302 (2010)] Eduardo C. Cuansing and Jian-Sheng Wang (Received 29 November 2010; published 6 January 2011)DOI: 10.1103/PhysRevB.83.019902 PACS number(s): 44 10 i, 63 22 m, 66 70 Lm, 99 10 Cd We discovered an error in our code that affected thedynamical behavior of the energy current. In determiningthe current we numerically calculate several integrals usingthe trapezoidal rule. The integrals in Eq. (7), however, werenumerically calculated incorrectly because of an error (amisplaced division by 2) in the code. This resulted in thecurrent to erroneously decay faster than the correct behavior.The figures for the numerical results in our manuscript shouldtherefore be replaced by Figs. wn in this erratum.Comparing the previously published figures and the cor-rected figures, we find that for long times the current inboth cases do approach the steady-state values calculatedindependently from the Landauer formula. Furthermore, theinitial negative spike in the transient current do also occur inboth cases. The difference lies in how fast the current decaysto the long-time steady-state value. From the corrected figures,the characteristic decay time is about 30 10 s. All of ourother conclusions remain the same. time [10 -14 s] -160-120-80-4004080 I L [ n W ] time [10 -14 s] -160-120-80-4004080 I R [ n W ] first-order all orders first-orderall orders (a) (b) 1. (Color online) Corrected figure replacing the previouslyFig. . Shown are plots of the current flowing out of theand (b) right leads. The (red) lines are the results when only thein the perturbation is used in the calculation. The leftlead has temperature 330 K while the right lead has temperatureK. The interparticle spring constant is 625 eV u)the on-site spring constant is 0625 eV u). -14 [ n W ] [ n W ] (b) 2. (Color online) Corrected figures replacing the previouslyFig. 4. (a) The sum of the currents, , when theaverage temperature between the leads are 10 K (red triangles),K (green squares), and K (blue circles). Theoffsets of the leads are (b) Plots of as functionsof the average temperature at time 12 7[ ] (red triangles),24 6[ ] (green squares), and 38 2[ ] (blue circles), where10 s. c u rr e n t [ n W ] 3. (Color online) Corrected figures replacing the previouslypublished Fig. 5. Shown are plots of the current when the left andhave the same temperature , where 10 K for theand K for the (blue) circles.1098-0121/2011/83(1)/019902(1) ©2011 American Physical Society FIG. 6. Time-dependent phonon current flowing out of the(a) left and (b) right leads in the transient regime. The redlines are the results when only the first-order term in theperturbation is used in the calculation. The temperaturesof the left and right leads are set to be T L = 330 K and T R = 270 K, respectively. Reproduced with permission fromRef. [85]. with β = r, a . The first-order term of Eq. (47) can beexpressed as G RL,β ( t , t ) = − k Z t dt a g R,β ( t − t a ) g L,β ( t a − t ) . (48)In order to calculate the time-dependent phonon cur-rent, the time variable is discretized into a large num-bers of segments. Since the analytic expressions for theequilibrium surface Green’s function g L,γ and g R,γ ( γ = r, a, < ) in Eqs. (46) and (48) have been given in the fre-quency domain , the corresponding time-dependent sur-face Green’s functions can be numerically calculated byFourier transform to obtain G RL,γ ( t , t ). Then one cansolve G RL,r ( t , t a ) and ∂G RL,a ( t a ,t ) ∂t required in Eq. (45)from Eq. (47) by transforming the integral into a sum.Finally, by solving the time-derivative of the Green’sfunction G RL,< ( t , t ) in Eq. (45), the time-dependentphonon current in the transient regime can be calculated.Figure 6 plots the time-dependent phonon current inthe transient regime by setting the temperatures of leftand right leads to be T L = 330 K and T R = 270 K,respectively. Once the coupling between left and rightleads is switched on, the transient currents of both leadsflow in an unexpected direction, namely, flow from thecolder lead to the hotter one. The transient currents thenincrease to positive and gradually approach the steady-state that can be calculated directly from the Landauerformula in the long-time limit. The time-dependent cur-rents exhibit oscillatory behavior and the oscillation fre-quency is comparable to the highest phonon frequenciesavailable in the system.In addition, the transient behavior of time-dependentphonon current can also be studied by the full-countingstatistics of heat transport in harmonic junctions basedon the NEGF technique . Wang et al. derivedthe generating function of energy counting statistics forphononic junctions which can be expressed in terms ofcontour-ordered Green’s function as ,ln Z ( ξ ) = − 12 Tr j,τ ln(1 − G Π A ) . (49)Here, the notation Tr j,τ represents the trace in bothspace index j and contour time τ . G is the Green’s func-tion defined on the Keldysh contour and Π A is obtainedfrom the difference of the original lead self-energy andthe lead energy shifted by the contour time arguments.In the long-time limit, the cumulant generating functionfor large t M can be expressed using Green’s functions inthe frequency domain,ln Z ( ξ ) = − t M Z + ∞−∞ dω π ln det { − G r Π L G a Π R [( e iξω − n L +( e − iξω − n R + ( e iξω − e − iξω − n L n R ] } . (50)This formalism is first given by Saito and Dhar andsatisfies the steady-state fluctuation theorem. Agarwallaet al. then investigated the full counting statistics ofheat transferred in harmonic chains in the presence ofboth temperature gradients and time-dependent driv-ing forces . The cumulant generating function for heattransferred from the leads to the central region was calcu-lated based on the two-time measurement concept usingthe NEGF method. The transient behavior and steady-state fluctuations were studied in atomic chains with dif-ferent initial conditions and the results were generalizedfor systems with multiple heat baths. III. THERMOELECTRIC TRANSPORTA. dc thermoelectric transport Since the observation of the Seebeck effect which re-vealed the interplay between thermal gradient and elec-tric potential, thermoelectricity has attracted much at-tention due to its potential applications in power gener-ation and refrigeration. Recently, the Seebeck effect wasstudied in various nanostructures which provides new op-portunities for designing thermoelectric devices with high ZT values . The Seebeck coefficient has been suc-cessfully measured in molecular junctions by trappingmolecules between two gold electrodes, which offers apromising way to study the fundamental physics in ther-moelectric energy conservation . A significant ZT valueof 0.6 is achieved experimentally at room temperaturein one-dimensional silicon nanowires with rough surfaceswhich exhibits a 100-fold reduction of thermal conduc-tivity due to the efficient phonon scattering .In dc transport, the thermopower is related to theelectric conductance of nanodevices which can be sim-ply modeled by the well-known Landauer-B¨uttiker for-malism within the NEGF approach. Similar to the phonon energy current, the electric current and the elec-tric heat current for spin-degenerate systems can be givenby ( ~ = e = 1 for simplicity) I = Z + ∞−∞ dǫπ T ( ǫ )( f L − f R ) , (51)and I h = Z + ∞−∞ dǫπ ( ǫ − µ ) T ( ǫ )( f L − f R ) . (52)Here, f α ( ǫ ) = 1 / [exp( ǫ − µ α T α )+1] is the Fermi-Dirac distri-bution function with µ α the chemical potential in lead α ( k B = 1 for simplicity). T ( E ) is the electric transmissionspectrum T ( E ) = Tr( G r Γ L G a Γ R ) , (53)where Γ α and G r ( a ) are the electric bandwidth functionof lead α and the electric retarded (advanced) Green’sfunction, respectively.In the linear response, namely, under small bias voltageand small temperature gradient, the electric current canbe linearly expanded I = Z + ∞−∞ dǫπ T ( ǫ ) (cid:20) − ∂f∂ǫ ∆ V − ∂f∂ǫ (cid:18) ǫ − µT (cid:19) ∆ T (cid:21) ≡ L ∆ V + L T ∆ T. (54)Here, ∆ V = V L − V R and ∆ T = T L − T R are the bias dif-ference and temperature difference between the left andright leads, respectively, and L n = − Z + ∞−∞ dǫπ T ( ǫ )( ǫ − µ ) n (cid:18) ∂f∂ǫ (cid:19) . (55)Similarly, the electric heat current can be expanded as, I h = L ∆ V + L T ∆ T. (56)The Seebeck coefficient, also called thermopower,which measures the magnitude of ∆ V to balance theelectric current along the reverse direction due to ∆ T ,is defined as, S = − ∆ V ∆ T (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) I =0 . (57)From Eq. (54), it is easy to obtain, S = − T L L . (58)We can also define the thermal conductance of elec-trons κ el when the electric current is zero. From I = 0and I h = κ el ∆ T , we obtain κ el = 1 T (cid:18) L − L L (cid:19) . (59)0By further defining the electron conductance from Ohm’slaw G = I/ ∆ V = L , we can replace L n in Eqs. (54) and(56) by G , S , and κ el and obtain the following relation (cid:18) II h (cid:19) = (cid:18) G G S G ST κ el + G S T (cid:19) (cid:18) ∆ V ∆ T (cid:19) . (60)Moreover, we can describe the thermoelectric effect bythe figure of merit, ZT , which gives the maximum effi-ciency of energy conservation in thermoelectric devices.It can be calculated by ZT = G S κ el + κ ph T. (61)Based on the NEGF method, the thermoelectric prop-erties were extensively investigated in low-dimensionalnanostructures . Gunst et al. studied the thermo-electric properties of graphene antidot structures by us-ing the π -tight-binding model. They found that the ZT can exceed 0.25 at room temperature and it is highlysensitive to the structure of antidot edges . Chen etal. studied the thermoelectric properties of graphenenanoribbons, junctions, and superlattices . Their find-ings indicate that the thermoelectric behavior is con-trolled by the width of the narrower part of graphenejunctions. Moreover, the thermoelectric transport wasstudied in hybrid graphene and boron nitride nanorib-bons and it was found that the ZT value can be signif-icantly enhanced by periodically embedding hexagonalboron nitride into graphene nanoribbons . Besides, theeffect of electron-phonon coupling and electron-electroninteraction on thermoelectric transport was studied in asingle molecular junction and it was found that ZT canbe enhanced by increasing electron-phonon coupling andCoulomb repulsion .In 2008, the spin Seebeck effect which generates thespin voltage from temperature gradient has been ob-served experimentally in a metallic magnet by Uchidaet al . How to manipulating and control the spin de-grees of freedom in thermal ways has attracted muchattention. Spin caloritronics concerning coupled spin,charge, and energy transport in magnetic structureswas introduced to focus on the relations between spinand heat current . In spin caloritronics, variousnonequilibrium phenomenons driven by thermal gradi-ent have been investigated such as thermal spin trans-fer torque , spin-polarized currents , and purespin currents . Using the first-principles calcula-tion combined with the NEGF method, a strongly spin-polarized current due to temperature difference was ob-tained in magnetized zigzag graphene nanoribbons bybreaking the electron-hole symmetry . The spin cur-rent can be completely polarized by tuning the gate volt-age. Moreover, a pure spin current was generated in atriangulene-based molecular junction on a large scale bychanging the temperature gradient and gate voltage .Apart from the charge and spin degrees of freedom,the valley degree of freedom can be used in valleytron-ics for the application of information processing similar (cid:11)(cid:68)(cid:12)(cid:11)(cid:69)(cid:12) FIG. 7. (a) Schematic diagram of zigzag graphene nanorib-bons with two semi-infinity leads (blue shadow). Two staticgate regions with v g = 0 . v g (orangeshadow) is tunable in the central region. (b) Valley current asa function of v g under different temperature gradients withfixed T R = 0 K. Reproduced with permission from Ref. [120]. to spin used in spintronics . A complete valley po-larized electronic current has been obtained by simplyintroducing the line defect in graphene . The genera-tion of a pure bulk valley current without net charge cur-rent through quantum pumping has also been reported ingraphene by using the well-known Dirac Hamiltonian .Analogous to spin caloritronics, valley caloritronics, acombination of valleytronics and thermoelectrics, hasbeen proposed to generate a valley polarized current ora pure valley current using thermal means .The valley Seebeck effect was first proposed in gatetunable zigzag graphene nanoribbons by Yu et al. usingthe tight-binding model within the NEGF framework ,as shown in Fig. 7(a). From the unique band structure ofzigzag graphene nanoribbons, one can find that the mo-mentum and valley index of electrons in the first subbandare locked together. Therefore, the left- and right-movingelectrons have valley index K and K ′ , respectively. Sinceat given energy the sign of f L − f R determines the di-rection of electron flow and the valley index, the valleycurrent of zigzag graphene nanoribbons can be simplyexpressed as , I v = Z dE π sgn( f L − f R )( f L − f R ) T ( E ) . (62)A pure valley current can be generated by the thermalgradient as well as the external bias. In order to con-1trol the pure valley current, the gate voltage v g appliedin the central region is modulated. Fig. 7(b) presentsthe pure valley current as a function of v g at differenttemperature gradient with T R = 0 K. It is found thereis a threshold gate voltage to open the valley current.Both the threshold gate voltage and on valley currentare proportional to the temperature gradient and thevalley current reaches the maximum value at the neu-tral gate voltage. These behaviors suggest the potentialapplications as a valley field-effect transistor driven bythe temperature gradient.Moreover, the dephasing effect and doping effect onthe valley Seebeck effect in zigzag graphene nanoribbonswere studied . It was found that the dephasing effectonly reduces the magnitude of pure valley current. Whilethe valley polarized current occurs by random doping ofboron and nitrogen atoms and the valley polarizationcan be effectively tuned by the doping concentration.Both the valley polarized current and pure valley currentcan also be obtained in wedge-shaped zigzag graphenenanoribbon junctions . In addition to graphene-basednanostructures, valley and spin thermoelectric transporthas also been investigated in silicene junctions andgroup-IV monolayers . B. time-dependent thermoelectric transport in thetransient regime Besides the static thermoelectric behavior, time-dependent thermoelectric transport is also an importantissue that may provide fundamental insights to under-stand the thermal response of mesoscopic systems. Gen-erally, there are two different schemes to study time-dependent quantum electronic transport. One is thepartition-free scheme (Cini scheme) in which the ini-tial state of the system is assumed to be at equilibriumthat can be described by a thermal density matrix .Then the system can be perturbed by applying a time-dependent voltage bias. Another way is the partitionedscheme (Caroli scheme) which assumes that the two-probe system is disconnected initially and the couplingbetween the scattering region and two leads is treated asthe time-dependent perturbation . In the following,we will discuss the time-dependent thermoelectric trans-port in the transient regime using the NEGF methodwithin the Caroli scheme.Within the Caroli scheme, the leads are assumed to atequilibrium states with the temperature T α and appliedbias V α before t = 0 and the couplings between leadsand the central region are turned on at t = 0. The ex-act solution of the transient electric current that beyondwide-band limit (WBL) can be given by I L ( t ) = Z dǫ π Tr[ A ( ǫ, t )Σ < ( ǫ ) B L ( ǫ, t )+ A ( ǫ, t )Σ In this review, we focus on the thermal transportin mesoscopic systems studied by using the NEGF ap-proach. We first give a brief introduction to the phononNEGF method and the detailed formalism of phononcurrent is presented in terms of phonon Green’s func-tion. Various theoretical investigations on quantumthermal transport in mesoscopic systems are discussed,which covers the interfacial thermal transport in one- dimensional atomic chains, the effect of nonlinearity andelectron-phonon coupling on the interfacial thermal con-ductance, phonon transport in multi-terminal systems,and time-dependent phonon transport in the transientregime. We also introduce the application of the NEGFmethod on the thermoelectric transport within the linearresponse theory. The formalism of the Seebeck coefficientand ZT value in the dc thermoelectric transport are givenand they are extended to the spin and valley caloritron-ics. The time-dependent thermoelectric transport in thetransient regime is further discussed within the Carolischeme.There are still many issues that deserve future inves-tigation in the field of thermal transport. For instance,manipulating phonons in two- and three-dimensional in-terfaces to achieve low interfacial thermal conductance,controlling the chirality of phonon in topological insula-tors, controlling other (quasi) particles such as magnonsand skyrmions by the means of thermal, and discoveringnew materials with low thermal conductance and highelectric conductance for optimized thermoelectric perfor-mance. From the aspect of the development of the NEGFmethod, time-dependent thermal transport, higher-orderfluctuations of thermal current, and the NEGF-DFTframework for phonon transport, are still open to ad-dress. We hope this brief review can inspire more investi-gations on quantum thermal transport and provide help-ful guidance on thermal engineering and applications. ACKNOWLEDGMENTS This work was financially supported by the Na-tional Natural Science Foundation of China (Grants Nos.12074190, 11975125, 11890703, and 11874221). ∗ [email protected] E. Pop, S. Sinha, and K. E. Goodson, Proceedings of theIEEE , 1587 (2006). N. Li, J. Ren, L. Wang, G. Zhang, P. H¨anggi, and B. Li,Rev. Mod. Phys. , 1045 (2012). G. Zhang and Y.-W. Zhang,Chinese Physics B , 034401 (2017). X. Chen, Y. Liu, and W. Duan,Small Methods , 1700343 (2018). D. G. Cahill, P. V. Braun, G. Chen, D. R. Clarke, S. Fan,K. E. Goodson, P. Keblinski, W. P. King, G. D. Mahan,A. Majumdar, H. J. Maris, S. R. Phillpot, E. Pop, andL. Shi, Applied Physics Reviews , 011305 (2014). B. Li, L. Wang, and G. Casati,Phys. Rev. Lett. , 184301 (2004). B. Li, L. Wang, and G. Casati,Applied Physics Letters , 143501 (2006). W. Chung Lo, L. Wang, and B. Li,Journal of the Physical Society of Japan , 054402 (2008). L. Wang and B. Li, Phys. Rev. Lett. , 177208 (2007). L. Wang and B. Li, Phys. Rev. Lett. , 267203 (2008). H. Zhu, J. Yi, M.-Y. Li, J. Xiao, L. Zhang, C.-W.Yang, R. A. Kaindl, L.-J. Li, Y. Wang, and X. Zhang,Science , 579 (2018). J. Lu, C. Qiu, M. Ke, and Z. Liu,Phys. Rev. Lett. , 093901 (2016). J. Lu, C. Qiu, L. Ye, X. Fan, M. Ke, F. Zhang, andZ. Liu, Nature Physics , 369 (2017). Y. Liu, Y. Xu, S.-C. Zhang, and W. Duan,Phys. Rev. B , 064106 (2017). S. Twaha, J. Zhu, Y. Yan, and B. Li,Renewable and Sustainable Energy Reviews , 698 (2016). D. Li, Y. Gong, Y. Chen, J. Lin, Q. Khan,Y. Zhang, Y. Li, H. Zhang, and H. Xie,Nano-Micro Letters , 36 (2020). L.-D. Zhao, S.-H. Lo, Y. Zhang, H. Sun, G. Tan, C. Uher,C. Wolverton, V. P. Dravid, and M. G. Kanatzidis,Nature , 373 (2014). M.-J. Lee, J.-H. Ahn, J. H. Sung, H. Heo, S. G. Jeon,W. Lee, J. Y. Song, K.-H. Hong, B. Choi, S.-H. Lee, andM.-H. Jo, Nature Communications , 12011 (2016). C. Chang, M. Wu, D. He, Y. Pei, C.-F. Wu, X. Wu, H. Yu, F. Zhu, K. Wang, Y. Chen, L. Huang, J.-F. Li, J. He, andL.-D. Zhao, Science , 778 (2018). H. Babaei, J. M. Khodadadi, and S. Sinha,Applied Physics Letters , 193901 (2014). C. Caroli, R. Combescot,P. Nozieres, and D. Saint-James,J. of Physics C: Solid State Physics , 916 (1971). Y. Meir and N. S. Wingreen,Phys. Rev. Lett. , 2512 (1992). A.-P. Jauho, N. S. Wingreen, and Y. Meir,Phys. Rev. B , 5528 (1994). J.-S. Wang, J. Wang, and N. Zeng,Phys. Rev. B , 033408 (2006). J.-S. Wang, N. Zeng, J. Wang, and C. K. Gan,Phys. Rev. E , 061128 (2007). N. Sergueev, D. Roubtsov, and H. Guo,Phys. Rev. Lett. , 146803 (2005). T. Shimazaki and Y. Asai,Phys. Rev. B , 075110 (2008). M. Paulsson, T. Frederiksen, and M. Brandbyge,Phys. Rev. B , 201101 (2005). A. Ferretti, A. Calzolari, R. Di Felice, F. Manghi,M. J. Caldas, M. B. Nardelli, and E. Molinari,Phys. Rev. Lett. , 116802 (2005). K. S. Thygesen and A. Rubio,Phys. Rev. B , 115333 (2008). J. Taylor, H. Guo, and J. Wang,Phys. Rev. B , 245407 (2001). M. Brandbyge, J.-L. Mozos, P. Ordej´on, J. Taylor, andK. Stokbro, Phys. Rev. B , 165401 (2002). Z.-Y. Ong and E. Pop, Phys. Rev. B , 075471 (2011). G. Zhang and B. Li, J. Chem. Phys. , 114714 (2005). G. Zhang and H. Zhang, Nanoscale , 4604 (2011). R. Yang and G. Chen, Phys. Rev. B , 195316 (2004). W. Li, N. Mingo, L. Lindsay, D. A. Broido, D. A. Stewart,and N. A. Katcho, Phys. Rev. B , 195436 (2012). W. Li, J. Carrete, N. A. Katcho, and N. Mingo,Comput. Phys. Commun. , 1747 (2014). J.-S. Wang, J. Wang, and J. T. L¨u,Eur. Phys. J. B. , 381 (2008). H. L. J. T. Jian-Sheng Wang, Bijay Kumar Agarwalla,Front. Phys. , 673 (2014). H. Haug and A.-P. Jauho, Quantum Kinetics in Transportand Optics of Semiconductors (Springer-Verlag, Berlin,1998). N. Mingo and L. Yang, Phys. Rev. B , 245406 (2003). T. Yamamoto and K. Watanabe,Phys. Rev. Lett. , 255503 (2006). L. Zhang, J. Thingna, D. He, J.-S. Wang, and B. Li,EPL (Europhysics Letters) , 64002 (2013). J. T. L¨u and J.-S. Wang, Phys. Rev. B , 165418 (2007). L. Zhang, J.-T. L¨u, J.-S. Wang, and B. Li,J. Phys.: Condens. Matt. , 445801 (2013). K. Gordiz and A. Henry,Scientific Reports , 18361 (2015). J. Chen, J. H. Walther, and P. Koumoutsakos,Advanced Functional Materials , 7539 (2015). W. A. Little, Canadian Journal of Physics , 334 (1959). E. T. Swartz and R. O. Pohl,Rev. Mod. Phys. , 605 (1989). L. Zhang, P. Keblinski, J.-S. Wang, and B. Li,Phys. Rev. B , 064303 (2011). C. B. Saltonstall, C. A. Polanco, J. C. Duda,A. W. Ghosh, P. M. Norris, and P. E. Hopkins, Journal of Applied Physics , 013516 (2013). G. Xiong, J.-S. Wang, D. Ma, and L. Zhang,EPL (Europhysics Letters) , 54007 (2020). B. Chen and L. Zhang,J. Phys.: Condens. Matt. , 125401 (2015). D. He, J. Thingna, J.-S. Wang, and B. Li,Phys. Rev. B , 155411 (2016). J. Fang, X. Qian, C. Y. Zhao, B. Li, and X. Gu,Phys. Rev. E , 022133 (2020). J. C. Kl¨ockner, M. B¨urkle, J. C. Cuevas, and F. Pauly,Phys. Rev. B , 205425 (2016). J. C. Kl¨ockner, R. Siebler, J. C. Cuevas, and F. Pauly,Phys. Rev. B , 245404 (2017). L. Cui, R. Miao, C. Jiang, E. Meyhofer, and P. Reddy,The Journal of Chemical Physics , 092201 (2017). L. Hu, L. Zhang, M. Hu, J.-S. Wang, B. Li, and P. Ke-blinski, Phys. Rev. B , 235427 (2010). J. Lu, K. Yuan, F. Sun, K. Zheng, Z. Zhang,J. Zhu, X. Wang, X. Zhang, Y. Zhuang,Y. Ma, X. Cao, J. Zhang, and D. Tang,ACS Applied Materials & Interfaces , 42708 (2019). H. Fan, M. Wang, D. Han,J. Zhang, J. Zhang, and X. Wang,The Journal of Physical Chemistry C , 16748 (2020). X. Chen, Y. Xu, X. Zou, B.-L. Gu, and W. Duan,Phys. Rev. B , 155438 (2013). W. Zhang, N. Mingo, and T. S. Fisher,Phys. Rev. B , 195429 (2007). Y. Xu, X. Chen, B.-L. Gu, and W. Duan,Applied Physics Letters , 233116 (2009). Y. Xu, X. Chen, J.-S. Wang, B.-L. Gu, and W. Duan,Phys. Rev. B , 195425 (2010). Z. Ding, Q.-X. Pei, J.-W. Jiang, W. Huang, and Y.-W.Zhang, Carbon , 888 (2016). S. Sadasivam, N. Ye, J. P. Feser, J. Charles,K. Miao, T. Kubis, and T. S. Fisher,Phys. Rev. B , 085310 (2017). Z. Zhang, Y. Xie, Q. Peng, and Y. Chen,Nanotechnology , 445703 (2016). Y. Blanter and M. B¨uttiker, Phys. Rep. , 1 (2000). M. B¨uttiker, Phys. Rev. Lett. , 1761 (1986). M. B¨uttiker, Phys. Rev. B , 9375 (1988). L. Zhang, J.-S. Wang, and B. Li,Phys. Rev. B , 100301 (2010). Y. Ming, Z. X. Wang, Z. J. Ding, and H. M. Li,New J. Phys. , 103041 (2010). T. Ouyang, Y. Chen, Y. Xie, X. L. Wei, K. Yang, P. Yang,and J. Zhong, Phys. Rev. B , 245403 (2010). Z.-X. Xie, K.-M. Li, L.-M. Tang, C.-N. Pan, and K.-Q.Chen, Appl. Phys. Lett. , 183110 (2012). Carbon , 818 (2020). L. Zhang, J.-S. Wang, and B. Li,New J. Phys. , 113038 (2009). Y. Xing, Q.-f. Sun, and J. Wang,Phys. Rev. B , 205339 (2006). Y. Xing, Q.-f. Sun, and J. Wang,Phys. Rev. B , 075324 (2007). M. Wei, M. Zhou, B. Wang, and Y. Xing,Phys. Rev. B , 075432 (2020). C. Strohm, G. L. J. A. Rikken, and P. Wyder,Phys. Rev. Lett. , 155901 (2005). E. C. Cuansing and J.-S. Wang,Phys. Rev. B , 052302 (2010). R. Tuovinen, N. S¨akkinen, D. Karlsson, G. Stefanucci, and R. van Leeuwen, Phys. Rev. B , 214301 (2016). E. C. Cuansing and J.-S. Wang,Phys. Rev. B , 019902 (2011). J.-S. Wang, B. K. Agarwalla, and H. Li,Phys. Rev. B , 153412 (2011). B. K. Agarwalla, B. Li, and J.-S. Wang,Phys. Rev. E , 051142 (2012). B. K. Agarwalla, J.-H. Jiang, and D. Segal,Phys. Rev. B , 245418 (2015). K. Saito and A. Dhar,Phys. Rev. Lett. , 180601 (2007). K. Saito and A. Dhar, Phys. Rev. E , 041121 (2011). Y. Dubi and M. Di Ventra,Rev. Mod. Phys. , 131 (2011). A. I. Hochbaum, R. Chen, R. D. Delgado, W. Liang,E. C. Garnett, M. Najarian, A. Majumdar, and P. Yang,Nature , 163 (2008). P. Reddy, S.-Y. Jang, R. A. Segalman, and A. Majumdar,Science , 1568 (2007). T. Gunst, T. Markussen, A.-P. Jauho, and M. Brandbyge,Phys. Rev. B , 155449 (2011). Y. Chen, T. Jayasekera, A. Calzo-lari, K. W. Kim, and M. B. Nardelli,Journal of Physics: Condensed Matter , 372202 (2010). K. Yang, Y. Chen, R. D’Agosta, Y. Xie, J. Zhong, andA. Rubio, Phys. Rev. B , 045425 (2012). Y. Xing, Q.-f. Sun, and J. Wang,Phys. Rev. B , 235411 (2009). M.-M. Wei, Y.-T. Zhang, A.-M. Guo, J.-J. Liu, Y. Xing,and Q.-F. Sun, Phys. Rev. B , 245432 (2016). B. Wang, J. Zhou, R. Yang, and B. Li,New Journal of Physics , 065018 (2014). J. Li, B. Wang, F. Xu, Y. Wei, and J. Wang,Phys. Rev. B , 195426 (2016). B. Zhou, B. Zhou, Y. Yao, G. Zhou, and M. Hu,Journal of Physics: Condensed Matter , 405303 (2017). P. N. Butcher, J. Phys.: Condens. Matt. , 4869 (1990). U. Sivan and Y. Imry, Phys. Rev. B , 551 (1986). G. D. Mahan, Many-Particle Physics (Springer, NewYork, 2000). J. Ren, J.-X. Zhu, J. E. Gubernatis, C. Wang, and B. Li,Phys. Rev. B , 155443 (2012). K. Uchida, S. Takahashi, K. Harii, J. Ieda,W. Koshibae, K. Ando, S. Maekawa, and E. Saitoh,Nature , 778 (2008). G. E. Bauer, A. H. MacDonald, and S. Maekawa,Solid State Communications , 459 (2010), spinCaloritronics. G. E. W. Bauer, E. Saitoh, and B. J. van Wees,Nature Materials , 391 (2012). M. Hatami, G. E. W. Bauer, Q. Zhang, and P. J. Kelly,Phys. Rev. Lett. , 066603 (2007). Z. Zhang, L. Bai, X. Chen, H. Guo, X. L. Fan,D. S. Xue, D. Houssameddine, and C.-M. Hu,Phys. Rev. B , 064414 (2016). M. Zeng, Y. Feng, and G. Liang,Nano Letters , 1369 (2011). X.-Q. Yu, Z.-G. Zhu, G. Su, and A.-P. Jauho,Phys. Rev. Lett. , 246601 (2015). S.-g. Cheng, Y. Xing, Q.-f. Sun, and X. C. Xie,Phys. Rev. B , 045302 (2008). Q. Wang, J. Li, Y. Nie, F. Xu, Y. Yu, and B. Wang,Phys. Chem. Chem. Phys. , 15736 (2018). D. Xiao, W. Yao, and Q. Niu,Phys. Rev. Lett. , 236809 (2007). C. E. Nebel, Nat. Mater. , 690 (2013). A. Rycerz, J. Tworzydlo, and C. W. J. Beenakker,Nat. Phys. , 172 (2009). D. Gunlycke and C. T. White,Phys. Rev. Lett. , 136806 (2011). Y. Jiang, T. Low, K. Chang, M. I. Katsnelson, andF. Guinea, Phys. Rev. Lett. , 046601 (2013). Z. Yu, F. Xu, and J. Wang, Carbon , 451 (2016). L. Zhang, Z. Yu, F. Xu, and J. Wang,Carbon , 183 (2018). X. Chen, L. Zhang, and H. Guo,Phys. Rev. B , 155427 (2015). X. Zhai, W. Gao, X. Cai, D. Fan, Z. Yang, and L. Meng,Phys. Rev. B , 245405 (2016). Z. Ping Niu and S. Dong,Applied Physics Letters , 202401 (2014). X. Zhai, S. Wang, and Y. Zhang,New Journal of Physics , 063007 (2017). G. Stefanucci and C.-O. Almbladh,Phys. Rev. B , 195318 (2004). M. Cini, Phys. Rev. B , 5887 (1980). C. Caroli, R. Combescot, D. Led-erer, P. Nozieres, and D. Saint-James,J. Phys. C: Solid State Phys. , 2598 (1971). Z. Yu, J. Yuan, and J. Wang,Phys. Rev. B , 235433 (2020). B. Wang, J. Wang, and H. Guo,Phys. Rev. Lett. , 398 (1999). J. Chen, M. ShangGuan, and J. Wang,New Journal of Physics , 053034 (2015). X. Chen, J. Yuan, G. Tang, J. Wang,Z. Zhang, C.-M. Hu, and H. Guo,Journal of Physics D: Applied Physics , 274004 (2018). F. G. Eich, A. Principi, M. Di Ventra, and G. Vignale,Phys. Rev. B , 115116 (2014). F. G. Eich, M. Di Ventra, and G. Vignale,Phys. Rev. B , 134309 (2016). ˇC. Lozej and T. Rejec, Phys. Rev. B , 075427 (2018). A. Cr´epieux, F. ˇSimkovic, B. Cambon, and F. Michelini,Phys. Rev. B , 153417 (2011). A. Kara Slimane, P. Reck, and G. Fleury,Phys. Rev. B , 235413 (2020). M. M. Odashima and C. H. Lewenkopf,Phys. Rev. B , 104301 (2017). M. Ridley and R. Tuovinen,J. Low Temp. Phys. , 380 (2018). A.-M. Dar´e and P. Lombardo,Phys. Rev. B , 035303 (2016). Z. Yu, L. Zhang, Y. Xing, and J. Wang,Phys. Rev. B , 115428 (2014). Z. Yu, G.-M. Tang, and J. Wang,Phys. Rev. B , 195419 (2016). H. Li, B. K. Agarwalla, and J.-S. Wang,Phys. Rev. B , 165425 (2012). M. Ridley, M. Galperin, E. Gull, and G. Cohen,Phys. Rev. B , 165127 (2019). G. Tang, J. Thingna, and J. Wang,Phys. Rev. B , 155430 (2018). G. Tang, X. Chen, J. Ren, and J. Wang,Phys. Rev. B97