Diffusion on edges of insulating graphene with intravalley and intervalley scattering

Diffusion on edges of insulating graphene with intravalley and intervalley scattering


Band gap engineering in graphene may open the routes towards transistor devices in which electric current can be switched off and on at will. One may, however, ask if a semiconducting band gap alone is sufficient to quench the current in graphene. In this paper we demonstrate that despite a bulk band gap graphene can still have metallic conductance along the sample edges (provided that they are shorter than the localization length). We find this for single-layer graphene with a zigzag-type boundary which hosts gapless propagating edge states even in the presence of a bulk band gap. By generating inter-valley scattering, sample disorder reduces the edge conductance. However, for weak scattering a metallic regime emerges with the diffusive conductance per spin, where is the transport mean-free path due to the inter-valley scattering and is the edge length. We also take intra-valley scattering by smooth disorder (e.g. by remote ionized impurities in the substrate) into account. Albeit contributing to the elastic quasiparticle life-time, the intra-valley scattering has no effect on .

I Introduction

Unique electronic properties of single-atomic-layer graphene Neto09 () stem from its two-dimensional (2D) semi-metallic energy spectrum with gapless conical conduction and valence bands. Engineering of a semiconductor-type band gap is expected to provide another desirable means of controlling electric current in graphene, laying the basis for electronic applications Novoselov07 (). Several mechanisms of the gap generation in the single-layer graphene have been discussed in literature Novoselov07 (); Zhou07 (); BN_Gio07 (); Strain_Moh09 (); Strain_Pereira09 (); Strain_Guinea10 (), including breaking the sublattice symmetry on a hexagonal boron-nitride substrate (see e.g. Ref. BN_Gio07, ) and using mechanical strain (see e.g. Refs. Strain_Moh09, ; Strain_Pereira09, ; Strain_Guinea10, ).

Although the opening of the band gap could have a desirable effect on transport in the 2D bulk of the material, it is, generally, not sufficient to control electric conduction near sample boundaries. Boundaries of graphene are natural extended defects that can host unusual electronic states such as edge states appearing on a zigzag boundary Fujita96 (). Various manifestations of such edge states in transport properties and spectroscopy of graphene have been discussed in recent years (e.g. Refs. Waka00, ; Peres06, ; Brey06, ; Sasaki06, ; Koba05, ; Niimi06, ; GT07, ; Akhmerov08, ; Burset08, ; Gusynin08, ; Yazyev08, ; Evaldson08, ; Wimmer08, ; Mucciolo09, ; Basko09, ; GT09a, ; Burset09, ; GT09b, ; Girit09, ; GT09c, ; Volkov09, ; Viana09, ; Ratnikov10, ; Herrera10, ; Qiao11, ; Rozhkov11, ; Roeder11, ; Gunlycke12, ). Most essential for our present discussion is the finding that the edge states remain conducting even when the 2D bulk turns into a band insulator, e.g. due to a staggered potential breaking the sublattice symmetry (see e.g. Refs. Volkov09, ; Qiao11, ). This is indeed expected since the edge states reside on one of the graphene sublattices, and, therefore, the influence of the staggered potential is reduced to an energy shift without dramatic changes in the edge-state dispersion. Thus, the edge states present a potential obstacle for the realization of the band insulator regime in graphene, providing pathways for leakage current. It is of both theoretical and practical interest to identify the factors that may help to reduce the edge conduction. As one of such factors, in this paper we theoretically consider structural disorder involving both smooth potential fluctuations, which couple states within the same graphene valley (intra-valley scattering), and atomically sharp defects generating inter-valley scattering.

Earlier, the influence of disorder on the edge transport was studied numerically for the conventional semi-metallic state of graphene (see e.g. Refs. Waka00, ; Wimmer08, ; Evaldson08, ; Mucciolo09, ; Waka09, ). Interestingly, the edge transport cannot be quenched by usual potential disorder, e.g. by smooth potential fluctuations due to remote ionized impurities. Only atomically sharp defects suppress the edge transport by mixing counter-propagating edge channels through inter-valley scattering. In this paper we extend these findings to the edge transport in the insulating graphene with the staggered potential, . Instead of using numerical approaches, we perform analytic diagrammatic calculations of the edge conductance, explicitly proving that the edge-transport mean-free path is limited only by scattering between graphene’s two valleys, and , and unaffected by smooth potential disorder. Our calculations indicate that the undoped insulating graphene can host diffusive metallic edge states with the conductance given per spin by


where the transport mean-free path is expressed in terms of the Fourier transform of the intervalley-disorder correlation function, , and the edge-state velocity , and are the edge-state Fermi points relative to the and valleys, respectively [ is the distance between the source and drain, assumed much larger than , but smaller than the typical size of the system that exhibits full localization Waka09 ()]. The specifics of the insulating graphene lies in the fact that the Fermi points are shifted with respect to and points by the staggered potential . Given furthermore the low edge-state velocity , Peres06 (); Sasaki06 (); Akhmerov08 () a metallic diffusive regime under weak scattering condition emerges without any doping of the material. This is in stark contrast with the conventional metallic transport which occurs when the Fermi level is pushed into a conduction or valence band.

The subsequent sections give a complete account of our theoretical approach: In Sec. II we introduce the model for the edge states in disorder-free insulating graphene with the staggered potential and calculate the edge-state Green’s functions. In Sec. III we introduce the model of disorder, calculate the disorder-averaged Green’s functions, the renormalized edge velocity and, finally, the edge conductance from Kubo formula. Section IV summarizes our results.

Ii Edge states in insulating single-layer graphene

ii.1 Boundary problem

We begin by analyzing the edge states in disorder-free graphene described by the effective four-band Hamiltonian:


where Pauli matrices and represent the valley and sublattice degrees of freedom, respectively [throughout the paper products of - and - matrices should be understood as direct products], is the bulk Fermi velocity determined by the nearest-neighbor-hopping energy and lattice constant, and is the staggered (e.g. substrate-induced BN_Gio07 ()) sublattice potential. Equation (2) adopts the following convention for the basis states:


where and label the sublattices and valleys, respectively.

Following our previous work on the edge states in semi-metallic graphene GT07 (); GT09a (); GT09b (); GT09c () we will work with the Green’s function, , defined by the equation


where energy includes an infinitesimal imaginary part , with () for the retarded (advanced) Green’s function, and is a unit matrix composed of the unit matrices in valley () and sublattice () spaces. Equation (4) will be solved in a semispace with a single edge at described by the boundary condition Akhmerov08 (); McCann04 () (see also Appendix A):


It involves two unit vectors , orthogonal to each other and to the vector normal to the boundary , ensuring the vanishing of the particle current normal to the edge Akhmerov08 (); McCann04 (). The vector components or serve to parametrize the boundary types considered below in Sec. II.3 and II.4.

The only restrictions on the boundary condition (5) are the time-reversal symmetry and the absence of the intervalley coupling. Thus, Eq.  (5) can be seen as a generalized continuum model for zigzag-type edges. It cannot be applied to the armchair edges because the latter couples the and valleys. However, extended armchair edges are unlikely to occur because they are less stable than zigzag edges (see, e.g., Ref. Girit09, ). As to the short-length armchair edges, they can be treated as a special type of the boundary defects causing inter-valley scattering, which is considered later in Sec. III. Therefore, Eq. (5) is a good starting point for analyzing defect-free graphene. It is also easy to verify that Eq. (5) ensures vanishing of the normal component of the particle current .

ii.2 Green’s function solution

The Green’s function is block-diagonal in valley space,


where are matrices which, in the basis defined by Eq. (3), have the following structures:


Let us first calculate the matrix elements of . Expanding

where is the edge length, and writing Eq. (4) in components, it is straightforward to express the off-diagonal elements and in terms of the diagonal ones as follows


where and satisfy the equations:


The boundary conditions for Eqs. (10) and (11) follow from Eqs. (5), (7) and (8) as


We seek the solutions to Eqs. (10) and (11) in the form of the linear combinations:

where the first terms are the Green’s functions of the unterminated sublattices, while the second ones are the decaying solutions of the corresponding homogeneous equations. The coefficients are obtained from boundary conditions (13) and (15) with the following results:


The first terms in Eqs. (17) and (18) vanish at the boundary and do not have poles within the gap, , implying that the edge states are entirely described by the second terms. The latter have a pole within the gap at . Assuming that the energy is close to this pole, we can neglect the first terms in Eqs. (17) and (18) and find a compact expression for the matrix [see Eq. (9)]:


The expression for can be obtained from Eq. (19) by replacing , , and :


The poles in Eqs. (19) and (20) yield the edge spectrum in valleys and under the condition that the arguments of the Heaviside (theta) functions in the numerators are positive. These restrictions are, in turn, enforced by the positiveness of the inverse decay length, . Consequently, the edge spectrum is given by


These edge states can also be viewed as the solution of a strictly 1D problem described by Green’s functions (19) and (20) integrated over the transverse coordinate :


Below we analyze these results for two commonly considered confinement types, viz. the zigzag edge and the mass confinement.

ii.3 Zigzag edge

Figure 1: Edge states (k) and within bulk band gap in valleys and , respectively: (a) zigzag boundary [see Eqs. (24) and (25) with ] and (b) mass confinement [see e.g. Eq. (29)]. Edge momentum is measured in units of .

In our parametrization of the boundary condition, the zigzag-type edge corresponds to the limit . In this case the edge-state spectrum, Eqs. (21) and (22), reduces to


where is the edge-state velocity,


It vanishes for , which corresponds to the flat edge-state band. In what follows we keep the velocity finite (but small compared with the bulk velocity ), so that the edge states remain dispersive. We also note that the staggered sublattice potential shifts the edge states by energy such that they cross the mid-gap energy at finite wave-vectors [defined in Eq. (1)] with respect to the and points (Fig. 1a). The energy shift reflects the fact that the zigzag edge states reside on one of the sublattices. This is seen from the matrix structure of the edge Green’s function [Eqs. (7) and (23)] which for has only two diagonal nonzero matrix elements. To simplify the model, from now on we will work with the effective 1D Green’s function: 1


where and are the projector matrices in the valley and sublattice spaces, respectively:


matrix elements and are given by


Equation (28) describes the left-moving () and right-moving () states originating from valleys and , respectively. Comparing Eqs. (7) and (8) with (26), we see that both left- and right-movers reside on the same sublattice (the case of sublattice corresponds to the boundary condition with ), as we mentioned in the introduction.

ii.4 Mass confinement

This confinement type Berry87 () is realized in the limit , resulting in the edge-state spectrum,


These edge states have the velocity equal to the bulk one, , and cross the mid-gap energy at points and . Unlike the zigzag edge, there is no energy shift due to the staggered potential because in this case the edge states propagate on two sublattices. This is again seen from the matrix structure of the edge Green’s function [Eqs. (7) and (23)] which for has both diagonal and off-diagonal matrix elements in each valley:


where the valley projectors are defined in Eq. (27) and matrix elements are given by


Equations (30) and (31) recover the corresponding results of Ref. GT10, for the edge states in 2D topological insulators. Their transport properties are reviewed in detail in Ref. GT12, .

Iii Edge conductance

Figure 2: Diagrammatic representations of (a) Kubo formula (32) for the conductance, (b) Dyson equation (36) for disorder-averaged Green’s functions (thick lines) and (c) equation (III.3) for renormalized current vertex in ladder approximation. Thin and dashed lines correspond to the unperturbed Green’s function and disorder correlator, respectively.

In the rest of the paper we focus specifically on transport properties of the edge states on a zigzag-type boundary. In the metallic regime the conductance of the edge states is given by Kubo formula:


where are the disorder-averaged retarded and advanced edge Green’s functions, and are the bare and renormalized edge-state velocity operators, and is the length of the edge (e.g. the distance between the source and drain contacts). To evaluate Eq. (32) we need to calculate first and in the presence of disorder.

iii.1 Potential and inter-valley coupling disorder

We assume that edge disorder can be described by equation


where the first diagonal term accounts for smooth random potential fluctuations, e.g. due to remote ionized impurities in the substrate. The potential is characterized by the correlation function


with averaging over disorder configurations (e.g. over impurity coordinates). The second term in Eq. (33) describes atomically sharp defects (lattice vacancies, short-length armchair edges etc.) which couple the valleys (and, generally, the sublattices). For this disorder type we use the correlation function


where the Kronecker symbols imply completely uncorrelated valley and sublattice disorder components. The calculations below are done for spatially isotropic disorder with and .

iii.2 Disorder-averaged Green’s function

Within the standard self-consistent Born approximation 2 the disorder-averaged Green’s functions can be obtained from the Dyson equation (see diagram in Fig. 2b) which (with suppressed superstripts for brevity) is given by


where is the self-energy:


and and are the Fourier transforms of the correlation functions and [see Eqs. (34) and (35)]. We seek the solution to Eq. (36) in the form of the projector expansion:


with two unknown scalar functions and . The ansatz (38) is valid only for the zigzag-type edge [Cf. Eq. (26)]. Inserting this into Eq. (37) we have


Here the first two terms include the potential scattering and inter-valley coupling between the left and right movers, and . This coupling originates from those disorder terms which swap both the valleys and sublattices. The other two terms in Eq. (39) result from the disorder which swaps either the valleys or the sublattices. The latter has no effect on because upon inserting Eq. (39) into Eq. (36) the corresponding products of projectors vanish: . Notice that the summation over the valley and sublattice indices in Eq. (39) yields the factor of 4.

Inserting Eq. (39) into Eq. (36) and collecting coefficients at we obtain algebraic equations for :


where are scalar functions given by


Like in conventional metals (see e.g. Ref. Rammer, ), Eqs. (40) can now be solved using the sharpness of the Green’s functions near Fermi points , which yields the finite quasiparticle life-time :


where the sign is positive or negative for the retarded or advanced functions, respectively. The potential and inter-valley scattering mechanisms give additive contributions to the spectral broadening:


where is the density of states per valley and spin. Since potential disorder cannot cause backscattering, the corresponding scattering rate involves the correlation function at zero momentum transfer , which corresponds to forwardscattering. In contrast, inter-valley scattering occurs between the Fermi points (see Fig. 1), so that the corresponding scattering rate involves the correlation function with finite momentum transfer .

iii.3 Vertex renormalization

We demonstrate below the interplay of the scattering times and [Eqs. (43) and (44)] in the disorder-renormalized velocity , which is one of the central result of this paper. In order to calculate the renormalized velocity in the conductance formula (32) we consider the vertex equation in the usual ladder approximation (see e.g. Ref. Rammer, and diagram in Fig. 2c):


Like the edge Green’s function (26) the bare edge velocity matrix (46) has only two diagonal elements and corresponding to two counter-prapagating channels from different valleys. We seek the solution to Eq. (III.3) in the form of the projector expansion:


with four unknown scalar functions , , and . Inserting Eqs. (38), (46) and (47) into Eq. (III.3) and collecting the coefficients at the projectors, we find that can be expressed through by means of


and satisfy the following equations:


Here again the integral over can be calculated using the sharpness of the Green’s functions at [see Eqs. (42) - (43)], which yields the following result:


In parallel, we perform the -integration in the conductance formula (32), obtaining


The required sum of the renormalized velocities at Fermi points is obtained from Eq. (50) as


yielding, finally, the edge conductance:


Iv Results and conclusions

We have demonstrated that in single-layer graphene with insulating bulk zigzag-like edges provide pathways for metallic conduction. Both intra- and inter-valley scattering have been taken into account. Although both scattering mechanisms contribute to the elastic broadening of the spectrum (43), only the intervalley scattering time (44) enters the edge conductance (53). The transport mean-free path can then be identified as . We emphasize that Eq. (51) holds under weak scattering condition , enforced by the staggered potential in zigzag-type terminated graphene. Let us discuss qualitatively the dependence of the edge conductance on the staggered potential and disorder strength, assuming a Gaussian correlation function for the inter-valley disorder,


where is the disorder correlation length and is the root-mean-square amplitude of the disorder. From Eqs. (1) and (54) we have


We see that the edge conductance exponentially increases with . The reason is that the inter-valley scattering involves the finite momentum transfer between the Fermi points , with the scattering probability reducing with . On the other hand, gets suppressed algebraically as with increasing root-mean-square disorder amplitude , which may help in practice to reduce the edge conductance.

G.T. acknowledges the hospitality of the Max Planck Institute for the Physics of Complex Systems in Dresden where this work was initiated. M.H. thanks the DFG for support in the Emmy Noether Programme and through Forscherguppe 760.

Appendix A Derivation of boundary condition (5)

We begin by deriving the boundary condition (5) for the Green’s function at the sample edge . To be concrete we consider the retarded Green’s function in real space and time . It is a matrix in valley and sublattice space, with matrix elements standardly expressed through the annihilation and creation field operators as


where (and independently ) runs over the index set and of the basis states introduced in Eq. (3), the double brackets denote averaging with the equilibrium statistical operator, and is the Heaviside function.

It is obvious from Eq. (56) that the boundary condition for the Green’s function is just the same as for . Indeed, the boundary condition for (derived earlier in Refs. Akhmerov08, ; McCann04, ) can be written as


where are the elements of the matrix


which is defined in Eq. (5) and main text. Inserting Eq. (57) into Eq. (56) at , we obtain the boundary condition for the matrix elements of the Green’s function:


which along with Eq. (58) yields Eq. (5).

Appendix B Derivation of boundary conditions (13) and (15)

In order to derive these equations we start with the boundary condition for the upper block of the Green’s function, [see Eqs. (5), (7) and (8)], which has explicit form


Therefore, for the upper diagonal element we have


Expanding in plane waves and using the relation [see Eq. (9)],


we obtain from Eqs. (61) and (62) a closed boundary condition for , which after elementary algebra yields the boundary condition (13). Repeating step by step the same calculation for in Eq. (60) leads to the boundary condition (15) for the other diagonal matrix element of the Green’s function.


  1. This is justified because of the small edge-state width as , so that it can be made smaller than other length scales in the -direction.
  2. This method cannot, however, be used for edge states in the strong localization regime studied, e.g. in Ref. Waka09, .


  1. A.H. Castro Neto, F. Guinea, N.M. Peres, K.S. Novoselov, and A.K. Geim, Rev. Mod. Phys. 81, 109 (2009).
  2. K. Novoselov, Nat. Materials 6, 720 (2007).
  3. S. Y. Zhou, G.-H. Gweon, A. V. Fedorov, P. N. First, W. A. de Heer, D.-H. Lee, F. Guinea, A. H. Castro Neto, and A. Lanzara, Nat. Materials 6, 770 (2007).
  4. G. Giovannetti, P. A. Khomyakov, G. Brocks, P. J. Kelly, and J. van den Brink, Phys. Rev. B 76, 073103 (2007).
  5. T. M. G. Mohiuddin, A. Lombardo, R. R. Nair, A. Bonetti, G. Savini, R. Jalil, N. Bonini, D.M. Basko, C. Galiotis, N. Marzari, K. S. Novoselov, A. K. Geim, and A. C. Ferrari, Phys. Rev. B, 79, 205433 (2009).
  6. V. M. Pereira, A. H. Castro Neto, and N. M. R. Peres, Phys. Rev. B 80, 045401 (2009).
  7. F. Guinea, M. I. Katsnelson, and A. K. Geim, Nat. Phys. 6, 30 (2010).
  8. M. Fujita, K. Wakabayashi, K. Nakada, and K. Kusakabe, J. Phys. Soc. Jpn. 65, 1920 (1996).
  9. K. Wakabayashi and M. Sigrist, Phys. Rev. Lett. 84, 3390 (2000).
  10. N. M. R. Peres, F. Guinea, and A. H. Castro Neto, Phys. Rev. B 73, 125411 (2006).
  11. L. Brey and H. A. Fertig, Phys. Rev. B 73, 235411 (2006).
  12. K. Sasaki, S. Murakami, and R. Saito, Appl. Phys. Lett. 88, 113110 (2006).
  13. Y. Kobayashi, K. I. Fukui, T. Enoki, K. Kusakabe, and Y. Kaburagi, Phys. Rev. B 71, 193406 (2005).
  14. Y. Niimi, T. Matsui, H. Kambara, K. Tagami, M. Tsukada, and H. Fukuyama, Phys. Rev. B 73, 085421 (2006).
  15. G. Tkachov, Phys. Rev. B 76, 235409 (2007).
  16. A. R. Akhmerov and C. W. J. Beenakker, Phys. Rev. B 77, 085423 (2008).
  17. P. Burset, A. Levy Yeyati, and A. Martin-Rodero, Phys. Rev. B 77, 205425 (2008).
  18. V. P. Gusynin, V. A. Miransky, S. G. Sharapov, and I. A. Shovkovy, Phys. Rev. B 77, 205409 (2008).
  19. O. V. Yazyev and M. I. Katsnelson, Phys. Rev. Lett. 100, 047209 (2008).
  20. M. Evaldsson, I. V. Zozoulenko, H. Xu and T. Heinzel, Phys. Rev. B 78, 161407(R) (2008).
  21. M. Wimmer, I. Adagideli, S. Berber, D. Tomanek, and K. Richter, Phys. Rev. Lett. 100, 177207 (2008).
  22. E. R. Mucciolo, A. H. Castro Neto, and C. H. Lewenkopf, Phys. Rev. B 79, 075407 (2009).
  23. K. Wakabayashi, Y. Takane, M. Yamamoto and M. Sigrist, New J. Phys. 11, 095016 (2009).
  24. D. M. Basko, Phys. Rev. B 79, 205428 (2009).
  25. G. Tkachov, Phys. Rev. B 79, 045429 (2009).
  26. P. Burset, W. Herrera, and A. Levy Yeyati, Phys. Rev. B 80, 041402 (2009).
  27. G. Tkachov and M. Hentschel, Phys. Rev. B 79, 195422 (2009).
  28. C. Ö. Girit, J. C. Meyer, R. Erni, M. D. Rossell, C. Kisielowski, L. Yang, C.-H. Park, M. F. Crommie, M. L. Cohen, S. G. Louie, and A. Zettl, Science 323, 1705 (2009).
  29. G. Tkachov and M. Hentschel, Eur. Phys. J. B 69, 499 (2009).
  30. V. A. Volkov and I. V. Zagorodnev, Low Temp. Phys. 35, 2 (2009).
  31. J. Viana-Gomes, V. M. Pereira, and N. M. R. Peres, Phys. Rev. B 80, 245436 (2009).
  32. P. V. Ratnikov and A. P. Silin, Fiz. Tverd. Tela 52, 1639 (2010) [Phys. Solid State 52, 1763 (2010)].
  33. W. Herrera, P. Burset, and A. Levy Yeyati, J. Phys.: Condens. Matter 22, 275304 (2010).
  34. Z. Qiao, S. A. Yang, B. Wang, Y. Yao, and Q. Niu, Phys. Rev. B 84, 035431 (2011).
  35. A.V. Rozhkov, G. Giavaras, Y. P. Bliokh, V. Freilikher, and F. Nori, Phys. Rep. 503, 77 (2011).
  36. G. Röder, G. Tkachov, and M. Hentschel, Europhys. Lett. 94, 67002 (2011).
  37. D. Gunlycke and C. T. White, J. Vac. Sci. Technol. B 30, 03D112 (2012).
  38. E. McCann and V. I. Fal’ko, J. Phys.: Condens. Matter 16, 2371 (2004).
  39. M. V. Berry and R. J. Mondragon, Proc. R. Soc. Lond. A Math. Phys. Sci. 412, 53 (1987).
  40. G. Tkachov and E. M. Hankiewicz, Phys. Rev. Lett. 104, 166803 (2010).
  41. G. Tkachov and E. M. Hankiewicz, arXiv:1208.1466.
  42. J. Rammer, Quantum Transport Theory (Westview Press, Boulder, CO, 2004).