Metal-to-insulator transition and electron-hole puddle formationin disordered graphene nanoribbons

Metal-to-insulator transition and electron-hole puddle formation
in disordered graphene nanoribbons

Gerald Schubert    Holger Fehske Institut für Physik, Ernst-Moritz-Arndt-Universität Greifswald, 17487 Greifswald, Germany
July 20, 2019

The experimentally observed metal-to-insulator transition in hydrogenated graphene is numerically confirmed for actual sized graphene samples and realistic impurity concentrations. The eigenstates of our tight-binding model with substitutional disorder corroborate the formation of electron-hole-puddles with characteristic length scales comparable to the ones found in experiments. The puddles cause charge inhomogeneities and tend to suppress Anderson localization. Even though, monitoring the charge carrier quantum dynamics and performing a finite-size scaling of the local density of states distribution, we find strong evidence for the existence of localized states in graphene nanoribbons with short-range but also correlated long-range disorder.


The experimental observation of a disorder-induced metal-to-insulator transition in graphene upon hydrogenation BMESHKR09 () has triggered a vivid debate on the nature of this transition. For high concentrations of hydrogen several mechanisms of gap opening have been discussed, such as full to transition, localization of areas, or erasing of midgap states  ENMMBHFBKGN09 () 111Note that we focus exclusively on transport properties and do not address other interesting issues such as defect-induced magnetism or enhanced chemical activity (D. W. Boukhvalov and M. I. Katsnelson, J. Phys. Chem. C 113, 14176 (2009); O. V. Yazyev, Rep. Prog. Phys 73, 056501 (2010)). Graphene, on the other hand, is a truly two-dimensional system, and the one-parameter scaling theory predicts that at zero temperature any finite amount of disorder should lead to Anderson localization (AL) AALR79 (). Otherwise, the existence of a scaling function might be questionable since the Fermi wavelength diverges near the charge neutrality point and there is no spatial scale on which the conductivity is much larger than  NM07 (). So far it seems that AL has not been seen in disordered graphene down to temperatures of liquid-helium CGPNG09 (). This surprising result has been attributed to strong charge carrier density fluctuations that break up the sample into electron-hole puddles ACFD08 (). Within these puddles the local chemical potential deviates enough from the charge neutrality point to allow for electron or hole conductivity. Mesoscopic transport is then determined by activated (variable-range) hopping or leakage between the puddles DAHR11 (). If the formation of electron-hole puddles is suppressed, however, AL might be observed. This has been reported by quite recent experiments in double-layer graphene heterostructures PZJMNCFWTGG11p ().

Previous theoretical work on disordered graphene strongly emphasizes the difference between short- and long-range scatterers. While the former applies to the case of hydrogenation, the latter rather describes the effect of charged impurities in the substrate KCTSMI10 (). Within the Dirac approximation, only short-range impurities cause intervalley scattering, and thus may lead to AL SA02 (). The presence of long-range impurities alone gives rise to intravalley scattering which is not sufficient to localize the charge carriers BTBB07 (); *NKR07. Another factor is the edge geometry of the graphene nanoribbons (GNRs) that determines the universality class of disordered samples as long as the phase coherence length exceeds the system size WTS07 (). Going beyond the Dirac approximation and describing graphene by a tight-binding model, it is natural to ask whether the scattering range is still decisive. Deviations from the idealized linear dispersion, a finite lattice spacing, and the trigonal lattice symmetry, which breaks the rotational symmetry of the Dirac cones, call for a careful numerical analysis of the localization properties within the tight-binding description XX07 (); *BC10; *CCRFP10.

In this work, we prove by unbiased numerics that experimentally relevant concentrations of hydrogen may induce a metal-to-insulator transition in actual-size graphene samples. We show that the single-particle wavefunctions of our disorder model are also localized for correlated long-range disorder. Even for potential fluctuations on an atomistic scale there is strong evidence for electron-hole puddle formation on an intrinsic scale of some ten nm, in agreement with recent experimental observations. DBMLL09 (); *ZBGZC09; *XSBJDWTJL11. In contrast to previous studies on temperature dependent transport in disordered graphene within the semiclassical Boltzmann approach DHL11p (), we restrict ourselves in the following to strictly zero temperature and adopt a purely quantum point of view.

Figure 1: (Color online) Cartoon of the substitutional disorder model describing hydrogenated graphene.

We consider a tight-binding Hamiltonian on the honeycomb lattice with sites, where the operators () create (annihilate) an electron in a Wannier state centered at site , and denotes the nearest-neighbor transfer integral. The landscape of onsite potentials results from the superposition of contributions of randomly distributed Gaussian impurities at positions  RTB07 (): By choice of the range of the individual impurity potentials can be continuously tuned from short-ranged to long-ranged. For and we recover the Anderson model on a GNR ssf09 (). Assuming a fixed for all impurities, the limit results in the binary alloy model in which only distinct sites acquire a finite onsite potential. Vacancies correspond to sites with , leading to a quantum site-percolation scenario SF08 (). The presence of adsorbed hydrogen atoms alters the hybridization of carbon atoms from to , partially removing the corresponding orbital from the -band. We model the yet finite probability of finding electrons at the adsorbant site by a finite value of the disorder strength (see Fig. 1).

Figure 2: (Color online) Spatial distribution of the normalized LDOS at for the binary alloy model with potential difference and impurity concentration (left column) and (right column). The disorder configurations shown in the top panel were used for both BC. Data obtained by exact diagonalization (ED), GNR sample size , corresponding to atoms.

Experimental results by Bostwick et al. BMESHKR09 () suggest a metal-to-insulator transition in graphene for a hydrogen coverage as low as . In Fig. 2 we contrast the spatial distribution of the local density of states (LDOS) at the Dirac point energy for a hydrogen coverage slightly below and above this threshold. For zigzag boundaries the well known edge states persist even in the presence of weak disorder. Impurities on the sublattice with high LDOS values near the GNR edges drastically reduce these values. On the other sublattice they do not have any effect. In the bulk of the ribbon, the LDOS is slightly enhanced as compared to the ordered case. Positive interference traps the wavefunction on sites in between the impurities. For periodic boundary conditions (PBC) the spatial LDOS distribution is clearly distinct for both impurity concentrations: We observe only slight local perturbations of the perfectly extended state for low impurity concentrations but a clearly localized state at hydrogen coverage. Measurements on a sample of this size therefore yield metallic (insulating) behavior for coverages (). The observed metallic character seems to be merely a consequence of a finite localization length that exceeds the system size. Note that the calculated state characteristics and experimental results agree qualitatively for PBC only. This underlines that the observed localization properties are intrinsic to short-range disordered bulk graphene. Edge effects arise on top, but are to a certain extent irrelevant in experiments, especially if mobilities are measured using a multiterminal Hall geometry KSZXH09 ().

Figure 3: (Color online) Upper panel: Distribution of the LDOS at experimentally relevant energies for the binary alloy model with different impurity concentrations and PBC. The sample width . Normalization of the distribution to directly relates its position to the height of the maximum and its width. Inset: Magnification of the averaged DOS with indications of the energies for which the LDOS distributions are shown. Lower panel: Finite-size scaling of the LDOS distribution. Data obtained by the kernel polynomial method with resolution adapted to the level-spacing, (for details see WWAF06 (); ssbfv10 ()).

In order to assert that AL takes place in infinite GNRs, we analyze the distribution of the normalized LDOS in Fig. 3. We restrict ourselves to three characteristic energies which are shown in the inset together with the averaged density of states . The behavior of the LDOS distribution upon finite-size scaling (lower panel) is a powerful criterion to detect AL for different kinds of disordered systems ssbfv10 (), even in presence of interactions BAF04 (); *SBH11. Extended states are characterized by an being independent of the system size. Otherwise sensitivity of the distribution to the system size indicates localization, which we indeed observe for binary alloy disordered GNRs for all energies and both impurity concentrations.

For a given state, the shape of the LDOS distribution and the extent of its shifting depend on ; the more pronounced the shift and the more asymmetric , the shorter is . Larger impurity concentrations enhance localization, as can be seen from the asymmetric shape of for . The persisting size dependence for proves localization also for such a weak randomness. Here the high sensibility of the LDOS distribution to the ratio of and system size is of vital importance. It allows to detect localization also in the case of weak disorder for which distinctly exceeds the system size and consequently is concentrated around unity.

Figure 4: (Color online) Time evolution of the local particle density on disordered zigzag GNRs after, at , finite (ordered) leads were attached to the left and right of the GNR. is normalized to the actual mean particle density on the GNR. The color scale is identical to Fig. 2. Sample dimensions are , corresponding to lattice sites. Left and middle column: Binary alloy model with . Right column: Gaussian correlated disorder model with , where is the inter-carbon distance. Here the potential is normalized to . The insets show magnifications of in the quasistationary regime together with the corresponding potential landscape.
Figure 5: (Color online) Quasistationary local particle density summed over the transverse direction. Each subpanel shows how changing one control parameter influences the baseline case (red dashed line: binary disorder, , , zigzag BC) while keeping the others fixed.

In Fig. 4 we contrast the quantum dynamics of a particle injected into zigzag GNRs with binary or Gaussian correlated disorder and impurity concentrations of and . As initial state we prepare an exact eigenstate of the ordered infinite graphene lattice in the lead left to the sample (see Fig. 1). After bringing the lead in contact with the sample, we let the system evolve in time by solving the time-dependent Schrödinger equation by a Chebyshev expansion technique SF08 (); FSSWFB09 (); *LSYD11. Due to the coupling with the disordered GNR, is not an eigenstate of the overall system but comprises admixtures of other states, mainly from the vicinity of . The snapshot at , where , confirms the intuition that spreading is faster the lower the impurity concentration is (see Fig. 4). At all states have reached quasistationarity, and we can extract their characteristic features. For the state spans the whole sample. The inset reveals its puddle-like structure with density fluctuations of two orders of magnitude on length scales of . At an impurity concentration of the particle density is reduced about two orders of magnitude between left and right edge of the GNR, providing a direct visualization of AL. Note that the local structure of remains puddle-like, but the spatial extent of the puddles is substantially reduced below . From the similarity of both local structures one may argue that for larger systems also states for will be localized. Note that correlated disorder results in a markedly smoother potential landscape. In addition, the electron-hole puddles are superimposed by a coarse-grained filamentary structure.

Integrating the local particle density over the transverse ribbon direction allows for a more quantitative analysis of the localization properties (see Fig. 5). For a fixed configuration of impurity positions we vary the relevant control parameters and extract from fitting the quasistationary density to an exponential decay. Removing part of the impurities results in a larger [Fig. 5(a)]. For the considered ribbon width the difference between zigzag- and PBC is marginal, while armchair edges drastically reduce the transmission [Fig. 5(b)]. This reflects the mismatch of preferred transport direction and ribbon axis. Varying the incident particle energy [Fig. 5(c)] we observe an enhanced transmission for , which might be attributed to resonant (localized) states. Note that at the position of the chemical potential of graphene on a SiC substrate, , is even shorter than at the Dirac point . For correlated long-range disorder [Fig. 5(d)], increases with , yielding extended states if . For the states are still clearly localized. This disagrees with results in the literature obtained within the Dirac approximation BTBB07 (); *NKR07, where the valley quantum number is conserved and localization is suppressed in absence of intervalley scattering SA02 (). We attribute this difference to the lattice discreteness and the breaking of the rotational symmetry of the Dirac cones by the trigonal symmetry of the honeycomb lattice. Moreover, even for a narrow banded initial state the dynamics is influenced by states from the whole energy spectrum where the graphene dispersion significantly deviates from the linear approximation. If all these aspects were taken into account, long-ranged disorder may cause localization within a tight-binding description in accordance with one-parameter scaling ZHBWXL09 ().

Figure 6: (Color online) Left: Finite-size scaling of the LDOS distributions for the Gaussian correlated disorder model with parameters matching Fig. 4. Right: CDF of the LDOS for various potential ranges .

Whether AL really occurs for correlated disorder can be proven by performing a finite-size scaling of the LDOS distribution. In doing so we find evidence for localization for from the shifting of the LDOS distribution (see left panel of Fig. 6). A remarkable feature of is the salient tail that develops for small values of on increasing . In the right panel of Fig. 6 this results in a kinked cumulated distribution function (CDF) instead of the approximate straight CDF for uncorrelated disorder 222 We would like to stress that there is numerical evidence against a strict log-normal shape of the LDOS distribution,for which the CDF would be an exact straight line, see Rodriguez et al. (2009), arXiv:1107.5736v1..

Figure 7: (Color online) Spatial distribution of the normalized LDOS at for the disorder configuration shown in the left panel with PBC. System sizes and color mapping are the same as in Fig. 2. Data obtained by ED.

So far we considered strong disorder () and weak randomness (), for which it is tempting to relate average puddle size and distance between the impurities. Interestingly, electron-hole puddles also arise for weak disorder strength and strong randomness, corresponding to a disorder landscape varying on an atomistic scale without any correlations. Such a modeling might be regarded as an attempt to capture the effect of the buffer layer forming between epitaxially grown graphene and its SiC substrate EBHJKLMORRRSWWS09 (). In this setup we refrain from considering a correlated potential landscape in order not to impose any a priori correlations in the LDOS. For weak disorder strength, , the LDOS at nevertheless becomes puddle-like with a characteristic scale of - (see middle panel of Fig. 7). Note that the choice of BC has no qualitative impact on the LDOS in this limit of the binary alloy model. Obviously, any subtle differences in the state characteristics induced by BC are masked by the randomness of the potential landscape. Increasing the potential difference to , drops below the system size and we observe clearly localized states.

To conclude, hydrogenated graphene behaves at zero temperature as a ’normal’ two-dimensional disordered system concerning AL, provided the extensions of ultra-high-quality samples become very large. If the localization length noticeably exceeds the system size, the sample nevertheless shows metallic behavior. We find that also certain long-range correlations in the disorder landscape yield localized single-particle wavefunctions. Most notably, we show that disorder-induced electron-hole puddles may arise for both disorder types. The intrinsic scale of the puddle-like structures in the eigenstates is not simply set by the distance between impurities, but results from subtle quantum interference effects. Even for atomic scale fluctuations of the disorder potential they might exceed 5 nm, which is in the range of experimentally measured values. The presence of electron-hole puddles, leading to intra- and inter-puddle transport, drives the system away from the metal-to-insulator transition, thereby masking AL PZJMNCFWTGG11p (). This resolves the puzzle why AL is so hard to detect in disordered graphene and GNRs.

We acknowledge financial support by DFG through the graphene priority program SPP 1459, KONWIHR-II as well as granting of computing time on HLRB Munich.


  • (1) A. Bostwick et al.,  Phys. Rev. Lett. 103, 056404 (2009)
  • (2) D. C. Elias et al., Science 323, 610 (2009).
  • (3) Note that we focus exclusively on transport properties and do not address other interesting issues such as defect-induced magnetism or enhanced chemical activity (see, e.g., D. W. Boukhvalov and M. I. Katsnelson, J. Phys. Chem. C 113, 14176 (2009); O. V. Yazyev, Rep. Prog. Phys 73, 056501 (2010)).
  • (4) E. Abrahams et al., Phys. Rev. Lett. 42, 673 (1979)
  • (5) K. Nomura and A. H. MacDonald, Phys. Rev. Lett. 98, 076602 (2007)
  • (6) A. H. Castro Neto et al., Rev. Mod. Phys. 81, 109 (2009)
  • (7) S. Adam et al., Phys. Rev. Lett. 101, 046404 (2008)
  • (8) S. Das Sarma et al., Rev. Mod. Phys. 83, 407 (2011)
  • (9) L. A. Ponomarenko et al., DOI: 10.1038/NPHYS2114
  • (10) J. Katoch et al., Phys. Rev. B 82, 081417 (2010)
  • (11) H. Suzuura and T. Ando, Phys. Rev. Lett. 89, 266603 (2002); J. Wurm et al. (2011) arXiv:1111.5969v1
  • (12) J. H. Bardarson et al., Phys. Rev. Lett. 99, 106801 (2007)
  • (13) K. Nomura, M. Koshino, and S. Ryu, Phys. Rev. Lett. 99, 146806 (2007)
  • (14) K. Wakabayashi, Y. Takane, and M. Sigrist, Phys. Rev. Lett. 99, 036601 (2007)
  • (15) S. Xiong and Y. Xiong, Phys. Rev. B. 76, 214204 (2007);
  • (16) J. Bang and K. J. Chang, ibid. 81, 193412 (2010);
  • (17) A. Chaves et al., ibid. 82, 205430 (2010);
  • (18) A. Deshpande et al., Phys. Rev. B 79, 205411 (2009)
  • (19) Y. Zhang et al., Nature Phys. 5, 722 (2009)
  • (20) J. Xue et al., Nature Mat. 10, 282 (2011)
  • (21) S. Das Sarma, E. H. Hwang, and Q. Li (2011), arXiv:1109.0988v1
  • (22) A. Rycerz, J. Tworzydło and C. W. J. Beenakker, EPL 79, 57003 (2007).
  • (23) G. Schubert et al., Phys. Rev. B 79, 235116 (2009)
  • (24) G. Schubert and H. Fehske, Phys. Rev. B 77, 245130 (2008)
  • (25) J. W. Klos et al., Phys. Rev. B 80, 245432 (2009)
  • (26) A. Weiße et al., Rev. Mod. Phys. 78, 275 (2006).
  • (27) G. Schubert et al., Phys. Rev. B 81, 155106 (2010)
  • (28) F. X. Bronold et al., Philos. Mag. 84, 673 (2004)
  • (29) D. Semmler et al., Phys. Rev. B 84, 115113 (2011)
  • (30) H. Fehske et al., Phys. Lett. A 373, 2182 (2009)
  • (31) Q. Liang et al., J. Phys. C 23, 345502 (2011)
  • (32) Y.-Y. Zhang et al., Phys. Rev. Lett. 102, 106401 (2009)
  • (33) Note that there is numerical evidence against a strict log-normal LDOS distribution, i.e. the CDF should slightly deviate from straight line, see Rodriguez et al. (2009), arXiv:1107.5736v1.
  • (34) K. V. Emtsev et al., Nature Mat. 8, 203 (2009)
Comments 0
Request Comment
You are adding the first comment!
How to quickly get a good reply:
  • Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
  • Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
  • Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minimum 40 characters and the title a minimum of 5 characters
Add comment
Loading ...
This is a comment super asjknd jkasnjk adsnkj
The feedback must be of minumum 40 characters
The feedback must be of minumum 40 characters

You are asking your first question!
How to quickly get a good answer:
  • Keep your question short and to the point
  • Check for grammar or spelling errors.
  • Phrase it like a question
Test description