Table Of ContentNumerical analysis of electronic conductivity in graphene with
resonant adsorbates: comparison of monolayer and Bernal bilayer
Ahmed Missaoui1, Jouda Jemaa Khabthani2, Nejm-Eddine Jaidane1, Didier Mayou3,4 and Guy Trambly de
Laissardi`ere5
7
1 1 LaboratoiredeSpectroscopieAtomiqueMol´eculaireetApplications,D´epartementdePhysique,Facult´edesSciencesdeTunis,
0 Universit´e Tunis El Manar, Campus Universitaire 1060 Tunis, Tunisia
2 2 LaboratoiredePhysiquedelamati`erecondens´ee,D´epartementdePhysique,Facult´edesSciencesdeTunis,Universit´eTunis
El Manar, Campus Universitaire 1060 Tunis, Tunisia
n
3 Universit´e Grenoble Alpes, Inst NEEL, 38042 Grenoble, France
a
J 4 CNRS, Inst NEEL, 38042 Grenoble, France
5 Laboratoire de Physique th´eorique et Mod´elisation, CNRS and Universit´e de Cergy-Pontoise, 95302 Cergy-Pontoise, France
7
2
January 27, 2017
]
l
l Abstract. Wedescribetheelectronicconductivity,asafunctionoftheFermienergy,intheBernalbilayer
a
h graphene (BLG) in presence of a random distribution of vacancies that simulate resonant adsorbates. We
- compareittomonolayer(MLG)withthesamedefectconcentrations.Thesetransportpropertiesarerelated
s to the values of fundamental length scales such as the elastic mean free path L , the localization length
e e
ξ and the inelastic mean free path L . Usually the later, which reflect the effect of inelastic scattering by
m i
phonons,stronglydependsontemperatureT.InBLGanadditionalcharacteristicdistancel existswhich
1
. isthetypicaltravelingdistancebetweentwointerlayerhoppingevents.Wefindthatwhentheconcentration
t
a of defects is smaller than 1%–2%, one has l ≤L (cid:28)ξ and the BLG has transport properties that differ
1 e
m
from those of the MLG independently of L (T). Whereas for larger concentration of defects L <l (cid:28)ξ,
i e 1
- and depending on Li(T), the transport in the BLG can be equivalent (or not) to that of two decoupled
d MLG. We compare two tight-binding model Hamiltonians with and without hopping beyond the nearest
n
neighbors.
o
c
PACS. 7 2.15.Lh – 7 2.15.Rn – 7 3.20.Hb – 7 2.80.Vp – 7 3.23.-b
[
1
v 1 Introduction tonian that takes only into account the hopping between
6
thenearestneighborsorbitals. Yetsomestudiesshowthe
1
Graphene consists of a monolayer (MLG) carbon atoms, importance of hopping beyond nearest neighbors on elec-
2
8 with sp2 hybridization, forming a 2D honeycomb lattice, tronic structure and transport properties [4,14,24,25,32].
0 withtwoequivalentatoms–atomAandatomB–inaunit In this paper, we present electronic properties of MLG
. cell. Linear dispersion relation of the p electron states andBLGobtainedbytwotight-binding(TB)models:the
1 z
close to the Fermi energy induces many fascinating trans- standard model with nearest neighbor only (TB1) and a
0
7 port properties which give rise to potential device ap- TBmodelincludingtheeffectofthehoppingbeyondnear-
1 plications [1,2,3,4]. Few-layer graphene also present un- est neighbors (TB2). This second model, which is more
: usualproperties.InparticulartheBernalbilayergraphene realistic, predicts some differences in the transport at en-
v
(BLG) with AB stacking, as in graphite, breaks the atom ergies close to the resonant energy of scatters.
i
X A / atom B symmetry and leads to quadratic dispersion We consider local defects, such as adsorbates or va-
r relation [5,6,7,8,9,10,11,12,13]. cancies, that are resonant scatters. Local defects tend to
a Electronictransportissensitivetostaticdefectswhich scatter electrons in an isotropic way for each valley and
are for example screened charged impurities, or local de- lead also to strong intervalley scattering. The T matrix
fects like vacancies or adsorbates, (hydrogen, adatoms of a local defect usually depends strongly on the energy.
or admolecules, chemically bound to one carbon atom of In the case of simple vacancies or adsorbates (atoms or
the surface of graphene layer). Theoretical studies of the molecules) that create one covalent bound with a carbon
effects introduced by the adsorbates on the conductivity atomsofMLG(BLG),theT matrixdivergesattheenergy
has been done for MLG (Refs. [14,15,16,17,18,19,20,21, E (with TB1 model E =E =0). For this reason,
MG MG D
22,23,24,25] and Refs. therein), and for BLG [26,27,28, theses scatters are called resonant scatters. The adsor-
29,30,31,32,12].MostofthemconsiderastandardHamil- bate is simulated by a simple vacancy in the plane of p
z
2 A. Missaoui et al.: Numerical analysis of conductivity in bilayer graphene
orbital as usually done [4,12,20,23,32]. Indeed the cova-
lent bonding between the adsorbate and the carbon atom A1 B1 A2 (cid:1)0 B2
ofgraphenetowhichitislinked,eliminatesthepz orbital A2
from the relevant energy window. The scatterers are dis- er B2 (cid:1)
y 1
tributed randomly in both planes and with the same con- p la (cid:1)0
centration in both planes. We consider here that the up o
T
B1 A1
and down spins are degenerate i.e. we deal with a para-
magnetic state. Indeed the existence and the effect of a
magnetic state for various adsorbates or vacancies is still
debated [33,34]. Let us emphasize that in the case of a Bottom layer
magnetic state the up and down spins give two different
Fig. 1. Planar view of the crystal structure of AB stacked
contributions to the conductivity but the individual con-
bilayer graphene. Atoms A and B on the lower layer are
1 1
tribution of each spin can be analyzed from the results
shown as blue circles. A , B on the upper layer are orange
2 2
discussed here.
circles.
We first determined the density of states (DOS) in
disordered MLG and BLG in presence of static scatter-
ers (vacancies) with various concentrations c from 0.5% 2 Electronic structure
to 10%. Elastic mean free path L , which depends on
e
the distribution of scatters and on energy E, is also com- 2.1 Tight binding Hamiltonian models
puted. From diffusive properties of wave packet in the
structure,theelectricalconductivityσ iscomputedversus BLG can be considered as two coupled MLGs with the
E andtheinelasticmeanfreepathL .L duetoelectron- toplayershiftedbyacarbonbondfromthebottomlayer.
i i
phonons interaction or magnetic fied. Roughly speaking, Consequently, BLG consists of four atoms in its unit cell,
largeL correspondtolowtemperaturelimitandsmallL two carbons A , B from the unit cell in the bottom layer
i i 1 1
to room temperature. The numerical method used takes and A , B in the the top layer where A sits at the top
2 2 2
into account all quantum effects. We show that difference of A (Fig. 1). We used a tight binding scheme (TB) [36].
1
between MLG and BLG is explained by considering the Only p orbitals are taken into account since we studied
z
average distance l over which a charge carrier travels in the electronic properties around the Fermi energy level.
1
a layer between two interlayer hoppings [12,35]. As ex- Interlayer interactions are not restricted to the ppσ terms
plained in this paper, l (cid:39)2nm and it is almost indepen- but ppπ terms have also to be introduced. The Hamilto-
1
dentontheconcentrationofdefects.IndeedBLGhassome nian has the form:
similar properties to MLG for small L and large c, i.e.
when Le (cid:39) Li < l1. But for small c valuies, Le (cid:39) Li ≥ l1, Hˆ =(cid:88)(cid:15)i|i(cid:105)(cid:104)i|+(cid:88)tij|i(cid:105)(cid:104)j|, (1)
theeffectsofinterlayerhoppingaffecttheelectronicprop- i (i,j)
erties of BLG with respect to MLG case. In this case, the
conductivityoftheBLGvarieswithclikeinusualmetals where i is the orbital located at ri and the sum runs over
where DOS is finite whereas MLG behaves like a semi- all neighboring i, j sites. The energy on the site is taken
metals with Dirac electrons. Moreover, different regimes equal to zero for first nearest neighbor model. tij is the
oftransportinBLGarefounddependingonthevaluesof hopping element matrix between site i and site j, com-
the energy, like in MLG [20,23,25]. puted from the Slater-Koster parameters,
Finally the localization length ξ is computed and we t = n2V (r ) + (1−n2)V (r ), (2)
ij c ppσ ij c ppπ ij
study the localization regime (which can be observed ex-
perimentally in a very low temperature regime, i.e. large where V and V depend from the distance r and
ppσ ppπ ij
L ). For the concentrations studied here, l < ξ, which n is the cosines direction along the Oz axis. It is either
i 1 c
means that in this localization regime the coupling be- equal to zero or to a constant because the two graphene
tween the two layers plays always a significant role. Thus layershavebeenkeptflatinourmodel.TheSlaterKoster
behavior of the BLG and MLG are different and localiza- parametersareexponentiallydecayingfunctionofthedis-
tion length is larger in the BLG. tance:
(cid:16) (cid:16) r (cid:17)(cid:17)
In Sect. 2, the two tight binding (TB) models used V (r )=−γ exp q 1− ij , (3)
are described and the corresponding density of states are ppπ ij 0 π a
(cid:18) (cid:18) (cid:19)(cid:19)
discussed. Transport properties are presented in Sect. 3. r
V (r )=γ exp q 1− ij . (4)
We first describe rapidly the computational method and ppσ ij 1 σ a
1
therelevantlengthstoanalyzeconductivityinBLG(Sect.
3.1). Then elastic scattering length L (Sect. 3.2) and mi- Itallows,accordingtothevalueofq ,totakeintoaccount
e π
croscopic conductivity σ (Sect. 3.3) are presented. Fi- both first neighbors or second neighbors. a is the nearest
M
nallyquantuminterferencecorrectionstotheconductivity neighbor distance within a layer, a = 1.418˚A, and a is
1
and localization length ξ are analyzed (Sect. 3.4). Sect. 4 the interlayer distance, a = 3.349˚A. First neighbors in-
1
is devoted to concluding remarks. teractioninaplaneischaracterizedbythecommonlyused
A. Missaoui et al.: Numerical analysis of conductivity in bilayer graphene 3
value γ =2.7eVandthe secondneighborsinteraction γ(cid:48) 0.2
issetto00.1γ0.Theratioqπ/ainEq.3isfixedbythevalu0e ((aa)) ❍∆ cc == 01.5
of the γ(cid:48). The interlayer coupling between two p orbitals ∇ c = 5
0 z ✩
in π configuration is γ1, γ1 = 0.48eV, and it is fixed to ★ ▲ cc == 1 00.5
obtain a good fit with ab initio calculation around Dirac m)] ● c = 1
energy in AA stacking and AB bernal stacking [36]. It is ato ★▼ cc == 510
wfoorrVthpptπoqaaσnnod=teVptqphaπσa:t=welochgao((cid:48)oγ−s0e/aγt0(cid:48)h)e s=am2e.2d1e8c˚Aa−y1c,oefficie(n5)t n[states/(eV 0.1 ★ ▼✩ ∆ ✩ ∇★
1 ∇
● ▼
▲ ❍
with a(cid:48) = 2.456˚A the distance between second neighbors ❍ ∆
●
in a plane. All p orbitals have the same on-site energy
z
(cid:15) in two planes. In order to obtain a Dirac energy E 0
i D -1.2 -0.8 -0.4 0 0.4 0.8 1.2
equal to zero in MLG, one fixes (cid:15)i equal to −0.78eV for Energy [eV]
TB model with hopping beyond first neighbors. This is 0.2
necessary because hopping beyond first neighbors breaks (b) ∆ c = 0.5
❍
the electron/hole symmetry and then shifts ED value [4, ❒ cc == 13
21].WehaveusedthismodelHamiltonianinourprevious ∇ c = 5
✩
c = 10
wgroarpkh[e3n6e,3b7i,la3y8e]rt.oInstuodrdyetrhteoelaencatrlyozneictshteruecffteucrte orofthaotepd- m)] ★ ▲● cc == 01.5
pthinegsibmepyolensdtTthBemfirosdt-enle(iTghBb1o)rwdiitshtafinrcsets-,newieghcboonrsihdoeprpailnsgo eV ato 0.1 ▼ ★■▼ ccc === 3510 ▲
omnolydAesalnedxdepsctlarhiiebneesddamainbeotvhpeea.rianmtreotdeursctitohnanwethceoncsoidmeprltehteatTrBes2- n [states/( ● ✩ ■ ❒❍ ★ ✩▼
∆
onant adsorbates –simple atoms or molecules such as H, ∇ ● ❒ ∇
OH,CH3–createacovalentbondwithsomecarbonatoms ❍■
of the BLG. To simulate this covalent bond, we assume ∆ ▲
thatthepz orbitalofthecarbon,thatisjustbelowthead- 0
-1.2 -0.8 -0.4 0 0.4 0.8 1.2
sorbate,isremoved.Inourcalculationsthemono-vacancies
Energy [eV]
are distributed at random between the both planes with
a finite concentration c. Fig. 2. (color online) Total density of states (DOS) n ver-
sus energy E, for concentrations c [%] of resonant adsorbates
(vacancies): (dashed line) TB1 (first neighbor hopping only),
2.2 Density of states (solid line) TB2 (with hopping beyond first neighbor); (a) for
MLG, (b) for BLG. The black dashed curve is pristine MLG
Figure 2(b) shows the total density of states (total DOS) (BLG), i.e. c=0, calculated by TB1.
n(E) in BLG for different concentrations c = 0.5% to
10% of defects randomly distributed. For comparison, the
DOSs of MLG with same TB Hamiltonian [23,24] are
alsoshown(Fig.2(a)).Withoutdefects,c=0,BLGDOS
E (cid:39) −0.2eV to E (cid:39) −0.3eV when defect concen-
andMLGDOSdifferforenergiesE suchas−γ <E <γ MG MG
1 1 tration increases from c=0.5% to 10%, like in monolayer
[12].Forsmallcconcentrations,c<∼1%,thisdistinction
graphene (Fig. 2(a)) [4,14,24,25]. In the following we
isstillobserved.Butremarkably,MLGandBLGhavevery
distinguish three cases according to energy values:
similar DOS, for concentration of defects c larger than
∼1%.
With TB1 model (first-neighbor hopping only), states
occur at energy E = 0. This is reminiscent of the (i) Sufficiently large energies, BLG DOS and MLG DOS
MG
midgap state produced by asymmetry between the num- are similar and they are not strongly modified by the
ber of atoms A and B in monolayer graphene [41,14]. In presence of resonant defects.
agreement with previous findings for MLG [4,14,16,20] (ii) Energies in the “plateau” due to vacancies but not in
and BLG [12,29,32], for large values of energies E the the midgap states E (cid:54)=EMG.
DOS is weakly affected by the presence of disorder. Fi- (iii) Energies in the midgap states, i.e. E = EMG = 0 for
nallyneartheDiracenergythereisanintermediateregime TB1(withfirst-neighborhoppingsonly)andE (cid:39)EMG
where the pseudo-gap is filled (“plateau”). for TB2 (with hopping beyond first neighbors).
WiththeTB2modelincludingallneighbors,themidgap
state is no longer at E = 0, but it becomes a broad
peak at negative energy E . E value varies from Thesethreecasescorrespondtodifferenttransportregimes.
MG MG
4 A. Missaoui et al.: Numerical analysis of conductivity in bilayer graphene
√
3 Transport properties V (MLG) = 2V (MLG) [38]. According to the renor-
0 B
malization theory [40] in 2D systems with static defects,
3.1 Computational method and relevant lengths diffusivity D always goes to zero at very large τ . At each
i
energy, the microscopic diffusivity D (microscopic con-
M
IntheframeworkoftheKubo-Greenwoodformulaforelec- ductivity σ ) is defined as the maximum value of D(τ )
M i
tronic transport properties, the quantum diffusion coef- (σ(τ )). We compute also the elastic mean-free path L
i e
ficient D (diffusivity) and conductivity σ are computed along the x-axis, from the relation [23],
by using the polynomial expansion method, developed by
Mayou,Khanna,RocheandTriozon[42,43,44,45,46].This 1 (cid:26)L2(E,τ )(cid:27) 2D (E)
L (E)= Max i i = M . (9)
numericalapproachallowsveryefficientcalculationsbyre- e V (E) τi τ V (E)
0 i 0
cursion (Lanczos algorithm) in real-space which take into
account all quantum effects. It has been used to study L is the average distance between two elastic scattering
e
quantumtransportindisorderedgraphene,chemicallydopedevents. At each energy, the elastic scattering times τ is
e
grapheneandbilayer(seeRefs.[18,19,20,22,23,24,32]and deduced from L by L (E)=V (E)τ (E).
e e 0 e
Refs. therein). Our calculations are performed on sample In our calculations τ and L are considered as ad-
i i
containing up to a few 107 carbon atoms, which corre- justable parameters. Roughly speaking, when L (cid:28) L
i e
sponds to typical sizes of about one micrometer square (τ (cid:28) τ ) the inelastic disorder dominates; it should cor-
i e
and allows to study systems with elastic mean-free length respond to very high temperatures. When L (cid:39) L (τ (cid:39)
i e i
of the order of few hundred nanometers. τ ), the conductivity is equal to microscopic conductivity,
e
Elastic scattering events are taken into account in the which should correspond to high temperature cases, typ-
Hamiltonian,buteffectsofinelasticscatteringbyphonons ically room temperature. And when L (cid:29) L (τ (cid:29) τ ),
i e i e
attemperatureT arenotincludedintheHamiltonian.To quantumlocalizationwilldominatestransportproperties;
consider the inelastic processes, we introduce an inelas- thisisthelowtemperaturelimit.Inthefollowingwethere-
tic scattering time τ (T) beyond which the propagation fore discuss the two important cases: L (cid:39) L (τ (cid:39) τ )
i i e i e
becomes diffusive due to the destruction of coherence by and L (cid:28)L (τ (cid:28)τ ).
i e i e
these inelastic processes. The effect of a magnetic field on Figure3showsthevariationoftheconductivityσ and
the electron propagation is not included directly in the theinelasticmeanfreepathL versusτ forenergiescorre-
i i
TB model, but a magnetic field B can have also a similar sponding to the three previous cases (see Sect. 2.2). The
incoherent dephasing effect. This dephasing effect occurs first case (i) (E = 1.5eV in Fig. 3) corresponding to a
on a length L (B) such that the flux of the magnetic field Boltzmanntransport:forlargevaluesofτ ,theconductiv-
i i
enclosed in the disk of radius L (B) is equal to the flux ity σ is almost constant as expected in a diffusive regime.
i
quantum h/e, i.e. L (B) (cid:39) (cid:112)h/eB. We treat these two This regime corresponds to energies for which the DOS is
i
dephasing effects in a phenomenological way through a weaklyaffectedbyscatters.Thethirdcase(iii)(E =EMG
Relaxation Time Approximation (RTA) as described here in Fig. 3) is determined by the transport of the midgap
after. In the RTA, the conductivity along the x-axis is states which are localized states. The latter case (ii) is
given by, [23] an intermediary regime between the previous two: for τi
closed to the elastic scattering time τ , there is a diffu-
e
σ(E ,τ )=e2n(E )D(E ,τ ), (6) sive behavior where the σ(τ ) reaches a maximum, σ ;
F i F F i i M
L2(E ,τ ) for larger values of τi, τi (cid:29) τe, σ(τi) decreases progres-
D(E ,τ )= i F i , (7) sively as expected in localization regime due to Anderson
F i 2τ
i localization in 2D [40].
In BLG another relevant time is the average traveling
where E is the Fermi energy, n(E ) is the density of
F F
timet betweentwointerlayerhoppingsofthechargecar-
states (DOS) and L is the inelastic mean-free path con- 1
i
riers, which is associated to an average traveling distance
ductivity along the x-axis. L (E ,τ ) is the typical dis-
i F i
l in a layer between two interlayer hoppings [12,35]. In
tance of propagation during the time interval τ for elec- 1
i
perfect BLG typical values of t and l can be easily es-
trons at energy E. 1 1
timated : t = ¯h/Γ (cid:39) 2×10−15s, where Γ (cid:39) 0.4eV is
WecomputethedistanceL ,thediffusivityD andthe 1 1 1
i
theinterlayerhoppingparametersoftheHamiltonian,and
conductivity σ at all inelastic scattering times τ and all
i
l (cid:39)V t (cid:39)2nm [35] where V is the velocity in MLG,
energies E for model Hamiltonian that includes inelastic 1 m 1 m
V (cid:39)V (cid:39)106ms−1. When there is elastic disorder such
scatters distributed randomly in the super-cell. At short m 0
that τ <t the value of t can be modified. A simple ar-
times τ –i.e. τ lower than elastic scattering time τ – the e 1 1
i i e
gumentmaybegivenasfollows:ABlochstateoftheMLG
propagation is ballistic and the conductivity σ increases
is still coupled to Bloch states of the other layer by the
when τ increases (Fig. 3),
i
sameintensity,typicallyΓ ,butthesestatesarenolonger
1
Li(E,τi)(cid:39)V0(E)τi when τi (cid:28)τe, (8) eigenstates and have a typical lifetime τe. Because of that
they have a spectral width W (cid:39) ¯h/τ . From the Fermi
e
where V (E) is a velocity at the energy E and short time Golden rule the typical time needed to jump from one
0
t. In crystals, V ≥ V where V is the Boltzmann ve- layer to the other will be such that h¯/t(cid:48) (cid:39)Γ2/W. There-
0 B B 1 1
locity (intra-band velocity) [39]. In BLG and MLG, V fore the new value of the interlayer hopping time t(cid:48) will
0 1
and V have the same order of magnitude: V (BLG) = be larger and will be such that t(cid:48) (cid:39)t (¯h/Γ τ ). Since the
B 0 1 1 1 e
A. Missaoui et al.: Numerical analysis of conductivity in bilayer graphene 5
1000
(a) ▲
10 (a) TB1 ❁ ❁ 100 10 ▲ ∆ ● ∆
❉ ❖ ❍
●
]0 ❁ 10 m] ❍ ❚
σ [G 1 ❖ ❉ ❉ ❖ 1 L [ni L [nm]e 1 ❚
0.1 ❁ ❖ E = EMG = 0 eV ▼ ▼ ∇
❖ ❉ E = 0.2 eV 0.1 ∇
❉ ❁ E = 1.5 eV ★✩ ✩ ★
0.0110-17 10-16 10-15 10-14 10-13 10-12 10-11 0.01 ★ ✩
τ [s]
i
100 1000 0.1
❁ -1.2 -0.8 -0.4 0 0.4 0.8 1.2
Energy [eV]
(b) TB2
10 ❁ ❉ 100 (b)
❖
10 ∆
]0 ❁ ❖ ❉ 10 m] ▲ ∆ ▲▲▲▲▲▲ ●
G 1 n
σ [ ❖ 1 L [i ● ❍ ❍ ❚
❚
❖
0.1 ❁❉ ❖❁❉ EEE == = 01 E..25M eeGVV= -0.22 eV 0.1 L [nm]e 1 ❒■ ❒ ■ ▼
∇
0.01 0.01 ▼ ∇
10-17 10-16 10-15 10-14 10-13 10-12 10-11 ✩ ★ ✩
τ [s]
i ★
✩
Fig.3. (coloronline)Conductivityσ(solidline)andinelastic
scatteringlengthL (dashedline)inBLGversusinelasticscat-
i 0.1
teringtimeτ forconcentrationc=1%ofresonantadsorbates -1.2 -0.8 -0.4 0 0.4 0.8 1.2
i
(vacancies), and for 3 energy values: (a) TB1 (first neighbor Energy [eV]
hoppingonly),(b)TB2(withhoppingbeyondfirstneighbors).
Fig. 4. (color online) Elastic mean free path L (E) versus
G =2e2/h. e
0 energy E (a) for MLG, (b) for BLG: (dashed line) TB1 (first
neighborhoppingonly),(solidline)TB2(withhoppingbeyond
first neighbor). Concentration of vacancies: triangle up c =
0.5%,circlec=1%,rectanglec=2%,squarec=3%,triangle
propagationisdiffusiveonthetimescalet>τ withdiffu-
e down c=5%, star c=10%.
sioncoefficientD (cid:39)V2τ ,thenewlengthl(cid:48) inpresenceof
0 e 1
defects is obtained from the relation, l(cid:48)2/t(cid:48) (cid:39) V2τ , and
1 1 0 e
thus l(cid:48) (cid:39) l (cid:39) 2nm depends weakly on disorder (i.e. l(cid:48)
1 1 1
almost independent on τ ). MLG. This influences strongly the localization regime
e
Asshowninthefollowing,electronicpropertiesofdis- as we discuss in Sect. 3.4.
orderedBLGdependonthevaluesofthelengthL ,ξ and
e
l which are characteristic of the BLG and of the amount
1
of elastic scatters:
3.2 Elastic mean-free path
– For low concentration of defects, c<c =1%-2%, and
l
for E (cid:54)= E , we have, l ≤ L (cid:28) ξ, and thus elec-
MG 1 e
tronic properties of BLG are influenced by interlayer Theelasticmean-freepathL (Eq.(9))alongthex-axisas
e
hoppings for every L ≥L values. afunctionoftheE isshownfigure4fordifferentvaluesof
i e
– Forlargerconcentrationofdefects,c>c ,oneobtains, vacancyconcentrationscinMLGandBLGwithbothTB
l
L(cid:32) <l (cid:28)ξ.Whentheeffectofquantuminterferences models.Itdependsontheenergyevenintheintermediate
e 1
on conductivity is small, i.e. when L (cid:39)L <l , BLG regimeandittakesafiniteandnon-zerovalueforE =E
e i 1 D
behavesastwodecoupledMLG.ForL(cid:32) <l ≤L (cid:28)ξ, butstayscomparabletothedistancedbetweenadsorbates
e 1 i √
the coupling between the two planes can affect the (vacancies) defined by, d(cid:39)1/ n , where n is the adsor-
a a
propagation of charge carriers before inelastic scatter- batesdensity.Numericalresults(Fig.4)showthatL val-
e
ing makes the propagation diffusive (i.e. on the length ues in BLG and MLG are close to each other. Moreover,
scaleL ).Inthisregime,quantumcorrectionstotrans- L (E) < l (cid:39) 2nm (i.e. τ < t ) for c > c (cid:39) 1%–2%;
i e 1 e 1 l
portarenotthesameintheBLGandintwodecoupled whereas for smaller c, L (E)≥l (i.e. τ ≥t ) for c<c .
e 1 e 1 l
6 A. Missaoui et al.: Numerical analysis of conductivity in bilayer graphene
-0.4 -0.2 0 0.2 -0.6 -0.3 0 0.3 0.6
10
20
(a) ∆ 1.5 ▲ ◗ 2 ✚ ✚
● ▲ (a) ❍ ∆ ∆
◗
8 ❍ ▲ 1 ❚ ★▼ ● 16 ✚ ❚ ❍
[G]M0 6 ● 0.5 ★▼ ★ ❚▼ σ [G]M0128 ❚ ❍ ∆ 1.51 ❚
σ
4 ❚ ▲ ● ❚ ✚∆ ❍
❍
∆ ▼ 4 ❚
2 ∇
▼ ★ 0
★ ✩ ∇ ✩ -1.2 -0.8 -0.4 0 0.4 0.8 1.2
Energy [eV]
0
-1.2 -0.8 -0.4 0 0.4 0.8 1.2 -0.6 -0.3 0 0.3 0.6
Energy [eV] 20
-0.4 -0.2 0 0.2 2
(b) ❚
20
▲ 16 ❚
(b) ●
2.4 ▲ 1.5 ❒ ❒
★
16 ❍ ▲▲▲▲▲▲ 1.6 ▼❚ ■ ▼■ ❚ ▼● ■ [G]M012 1 ∇✩ ∇ ✩
12 σ 8 ❚
G]0 ❚ ● 0.8 ❒
σ [M 8 ∆ ● ❍❚ 4 ∇ ❒ ❚ ∇
✩ ✩
❒ ∆ ❒
■ ▲ ■ ▼ 0 -1.2 -0.8 -0.4 0 0.4 0.8 1.2
4 ▼
∇ ∇ ★ Energy [eV]
★ ✩ ✩
Fig. 6. (color online) Microscopic conductivity σ in BLG
M
0 for TB1 (first neighbor hopping only): (a) c < c and σ
-1.2 -0.8 -0.4 0 0.4 0.8 1.2 l M
with a constant minimum value (plateau) which decreases as
Energy [eV]
c increases, (b) c>c and σ reaches a minimum value inde-
l M
Fig. 5. (color online) Microscopic conductivity σM versus pendentofcvalue(seetext).Concentrationofvacancies:semi
energy E for (a) MLG and (b) BLG: (dashed line) TB1 (first circle c = 0.05%, cross c = 0.25%, triangle up c = 0.5%, cir-
neighborhoppingonly),(solidline)TB2(withhoppingbeyond cle c = 1%, rectangle c = 2%, square c = 3%, triangle down
first neighbor). Concentration of vacancies: triangle up c = c = 5%, star c = 10%. Insert shows enlarged curves around
0.5%,circlec=1%,rectanglec=2%,squarec=3%,triangle E =0. G =2e2/h.
D 0
down c = 5%, star c = 10%. Insert shows enlarged curves
around E for TB2. G =2e2/h.
D 0
σ is calculated with the Bloch-Boltzmann approach [4,
B
12]. In this regime the conductivity decreases with the
3.3 Microscopic conductivity concentration of defects.
Intheintermediateenergyvaluesregime(ii),thesemi-
AsshownFig.3,thestaticscatteringeventsperturbstronglyclassical approach fails and the behavior depends on c.
thewavepacketpropagationandamaximumvalueofthe From Fig.6(a), forsmallcvalues,typicallyc<c (cid:39)1%–
l
conductivity σ(τ ), called microscopic conductivity σ , 2% (i.e. L (cid:39) L ≥ l ), σ seems to reach a constant
i M i e 1 M
σ (E) = e2n(E)D (E), is reached. σ calculated from minimum value (“plateau”), but this minimum σ value
M M M M
both TB models are shown in Figs. 5 and 6 for different decreasesascincreases. Thisconcentrationdependenceis
concentrationscofvacanciesinMLGandBLG.According specific to BLG and is not observed in MLG. For larger
to the renormalization theory [40], this value is obtained c (Fig. 6(b)), c > c (i.e. L (cid:39) L < l ), σ reaches
l i e 1 M
whentheinelasticmeanfreepathL andtheelasticmean a minimum value independent on c value: σ (cid:39) 1.2G
i M 0
free path L are comparable, L (cid:39)L , which corresponds where G = 2e2/h. This values for BLG is [29,30] two
e e i 0
to τ (cid:39) τ . As L and τ decrease when the temperature times the universal value of the conductivity, ∼4e2/(πh),
i e i i
T increases, the microscopic conductivity is a good es- expected in presence of resonant scatters in MLG [16,18,
timation of the high temperature conductivity (or room 20,22,23,25].
temperature conductivity). Results with TB2 model (including hopping beyond
For energies corresponding the to Boltzmann regime, nearest neighbors) (fig. 5) show that a plateau of the
i.e. regime (i) described in Sec. 2.2, σ (cid:39) σ , where microscopic conductivity near the Dirac energy exists in
M B
A. Missaoui et al.: Numerical analysis of conductivity in bilayer graphene 7
1 0.8
(a) EE == 00..0150 (b) E = 0.28 105 c = 0.5
E = 0.14 E = 0.32 c = 1
0.8 EE == 00..1284 EE == 00..346 0.6 104 cc == 35
G]00.6 EEEE ==== 0---000.1...1120050 EEE === 00-0..13.6 0.4 [nm] ξ110023 cc == 15
σ[
0.4 10
0.2 10 0.1 0.2 0.3 0.4
0.2
Energy [eV]
Fig. 8. (color online) Localization length versus energy for
0 0 differentconcentrationsc[%]ofvacancies:(dashedline)MLG,
0.1 1 10 100 100010000 0.1 1 10 100 1000
L[nm] L [nm] (solid line) BLG. This results are computed with TB1 model
i i
2 1.6 (first neighbor hopping only).
E = -0.06
(c) E = -0.06 (d) E = 0.20
E = -0.10 E = 0.20
1.6 E = -0.10 E = 0.28
EE == --11..4400 EE == 00..2382 1.2 localization effects that are a consequence of quantum in-
E = 0.18 E = 0.32 terferences.Inthatcase,forallvacancyconcentrationval-
G]01.2 EEE === 000...122800 EE == 00..3366 0.8 duieffsecr,enLti l(cid:29)ocal1liz(ai.tei.onτib(cid:29)ehta1v)i.oTrhtehraenforMeLBGL.GAsshosuhlodwhnavine
σ [
Fig. 7, σ(L ) follows the linear variation with the loga-
0.8 i
rithm of the inelastic mean free path L , like in the case
i
0.4 of monolayer graphene [22,23],
0.4
(cid:18) (cid:19)
L
σ(E,L )=σ (E)−αG log i , (10)
i 0 0 L (E)
0 0 e
0.1 1 10 100 100010000 0.1 1 10 100 1000
Li [nm] Li [nm] where G0 = 2e2/h, Le(E) is defined by Eq. (9), and σ0
valuesareintherangeofσ values.Forlowconcentration
M
ofdefectsc=0.5%and1%onecanestimateα(cid:39)0.26,and
Fig. 7. (coloronline)Conductivityσ(L )versustheinelastic
i for larger c, c = 3% and 5%, α (cid:39) 0.32. These values are
scatteringlengthL in(a)(b)MLGand(c)(d)BLG,forcon-
i closetotheresultfoundinmonolayergraphenewithsame
centration (a) (c) 1% and (b) (d) 5%, at different energies E
computational method [23] and close to the prediction of
[eV] in the plateau of σ (E): (dashed line) TB1 (first neigh-
M the perturbation theory of 2D Anderson localization for
borhoppingonly),(solidline)TB2(withhoppingbeyondfirst
which α=1/π [40]. As in the monolayer case, this linear
neighbor). G =2e2/h.
0
variation of σ with logL is found for both models, with
i
nearest neighbor hopping only and with hopping beyond
MLG and BLG, but is not symmetric due to the symme- first nearest neighbors.
try breaking electron-hole. Nevertheless in this case, for We finally define the localization ξ length from the
c > cl, σM values are still close to the universal conduc- expression (10) by extrapolation of σ(Li) curves (Fig. 7)
tivity plateau; whereas for c<cl, the minimum value of when : σ(Li =ξ)=0, and then
σM decreases as c increases like with TB1 model. (cid:18)σ (E)(cid:19)
Forenergiesinthemidgapstates(iii)withTB1model ξ(E)=L (E)exp 0 . (11)
e αG
(first neighbor hopping only), an anomalous behavior of 0
the conductivity is obtained and there is a peak of σ
M The ξ values for energies in the plateau of σ (i.e. case
at E = E = 0. With TB2 (including hopping beyond M
MG (ii) described in Sect. 2.2) are shown figure 8. For large
nearestneighbors),thisanomalousbehaviorisstillslightly
concentrations of defects (c > c ), ξ in MLG and BLG
presentatE (cid:39)E ,butthechangeintheconductivityis l
MG is almost independent on the energy. Moreover ξ(BLG) is
rathersmall.Thus,asinmonolayergraphene[22],conduc-
always larger than ξ(MLG). This difference results from
tion by the midgap states is very specific to TB1 model.
the fact that σ (cid:39)2σ (Sect. 3.3) and then, for
0,BLG 0,MLG
c>c ,
l
(cid:18) ξ(E) (cid:19) (cid:18) ξ(E) (cid:19)2
3.4 Quantum localization regime (cid:39) , (12)
L (E) L (E)
e BLG e MLG
Inthissectionweconsiderthecaseoflargeinelasticmean- whereL aresimilarinMLGandBLG(Sec.3.2).Itisthus
e
free path, L (cid:29) L (i.e. τ (cid:29) τ ), which should corre- amultilayereffectonquantuminterferencesthatmodifies
i e i e
sponds to the low temperature limit. The conductivity in the 2D behavior with respect to MLG cases. For low de-
the intermediate energy case (ii) (see Sect. 2.2) presents fects concentration (c < c ), interlayer coupling modifies
l
8 A. Missaoui et al.: Numerical analysis of conductivity in bilayer graphene
alsoquantuminterferences. Therefore,foreveryresonant Acknowledgment
scatterer concentration, quantum interference corrections
to the conductivity in BLG and in MLG are not similar. The authors wish to thank C. Berger, W. A. de Heer,
L. Magaud, P. Mallet and J.-Y. Veuillen for fruitful dis-
cussions.Thenumericalcalculationshavebeenperformed
at Institut N´eel, Grenoble, and at the Centre de Cal-
4 Conclusion culs (CDC), Universit´e de Cergy-Pontoise. This work was
supported by the Tunisian French Cooperation Project
(Grant No. CMCU 15G1306) and the project ANR-15-
Toconcludewehavestudiednumericallythequantumdif-
CE24-0017.
fusionofchargecarriersinmonolayergrapheneandBernal
bilayer graphene in the presence of local defects. These
defects are simulated by simple vacancies randomly dis-
Author contribution statement
tributed in the structure. Among the fundamental length
scales in the MLG and BLG there are the elastic mean
Numerical calculations have been performed by A. Mis-
free path L , the localization length ξ and the inelastic
e
saoui and G. Trambly de Laissardi`ere. All authors took
mean-free path L which in real systems depends on tem-
i
part in drafting of the manuscript as well as in analysis
perature or on magnetic field. For the bilayer, there is
and interpretation of the results.
an additional length scale which is the typical distance
l (cid:39)2nm over which an electron travels in a plane before
1
hoppingtotheotherplane.WehavecomparedtheBilayer
and monolayer transport properties for identical concen- References
trations of vacancies. We show that these properties can
be either similar or different depending on the compari- 1. P. R. Wallace, Phys. Rev. 71, 622 (1947)
son between l and the other three length scales L , L , 2. C. Berger, Z. Song, X. Li, X. Wu, N. Brown, C. Naud, D.
1 e i
ξ. This relation explains essentially the numerical results Mayou, T. Li, J. Hass, A. N. Marchenkov, E. H. Conrad, P.
detailed in this paper. N. First, W. A. de Heer, Science 312, 1191 (2006)
3. M.I.Katsnelson,K.S.Novoselov,A.K.Geim,NaturePhys
Our results show that for strong concentration of de-
2, 9 (2006)
fects, c > c (cid:39) 1%–2%, the bilayer graphene could be
l
4. A. H. Castro Neto, F. Guinea, N. M. R. Peres, K. S.
equivalent to two independent disordered monolayers of
Novoselov, A. K. Geim, Rev. Mod. Phys 81, 109 (2009)
graphene,becausetheelasticmeanfreepathL issmaller
e 5. S. Latil, L. Henrard, Phys. Rev. Lett 97, 036803 (2006)
thantheaveragedistancel .Thereforeforc>c ,theuni-
1 l 6. E. V. Castro, K. S. Novoselov, S. V. Morozov, N. M. R.
versal aspects of the conductivity are present in bilayer,
Peres, J. M. B. Lopes dos Santos, J. Nilsson, F. Guinea, A.
as in monolayer graphene, with (TB1) or without (TB2)
K. Geim, A. H. Castro Neto, Phys. Rev. Lett 99, 216802
the hopping beyond nearest neighbors.
(2007)
In the high temperature limit, i.e when inelastic scat- 7. F. Zhang, B. Sahu, H. Min, A. H. MacDonald, Phys. Rev.
tering length Li is small, Li (cid:39) Le, the conductivity in B 82, 035409 (2010)
bilayer is almost equal to two times the universal min- 8. F. Varchon, P. Mallet, J.-Y. Veuillen, L. Magaud, Phys.
imum plateau of microscopic conductivity in monolayer Rev. B 77, 235412 (2008)
graphene (except for the Dirac energy with TB1 model 9. T. Ohta, A. Bostwick, T. Seyller, K. Horn, E. Rotenberg,
that takes only into account nearest neighbor hopping). Science 313, 951 (2006)
For smaller c, c < c , the BLG should be considered like 10. I. Brihuega, P. Mallet, C. Bena, S. Bose, C. Michaelis, L.
l
a usual metal: with static defects the minimum of micro- Vitali,F.Varchon,L.Magaud,K.Kern,J.Y.Veuillen,Phys.
scopic conductivity of BLG increases when the c values Rev. Lett. 101, 206802 (2008)
decreases. For the parameters studied here, the local- 11. S.DasSarma,S.Adam,E.H.Hwang,E.RossiRev.Mod.
ization length ξ is larger than the traveling distance l Phys 83, 407 (2011)
1
12. E. McCann and M. Koshino, Rep. Prog. Phys. 76 056503
between two interlayer hopping; and therefore the BLG
(2013)
and MLG have different localization lengths at the same
13. S. Ulstrup, J. C. Johannsen, F. Cilento, J. A. Miwa,
concentration, even if they have similar elastic-mean-free
A. Crepaldi, M. Zacchigna, C. Cacho, R. Chapman, E.
paths. The localization length is the largest in the BLG.
Springate, S. Mammadov, F. Fromm, C. Raidel, T. Seyller,
In the limit, i.e L (cid:29) L , (which is relevant at low tem-
i e F.Parmigiani,M.Grioni,P.D.C.King,P.Hofmann1,Phys.
perature) and for all c values, the conductivity follows a
Rev. Lett. 112, 257401 (2014)
linearvariationwiththelogarithmofL inMLGandBLG
i 14. V.M.Pereira,J.M.B.L.dosSantos,A.H.C.Neto,Phys.
andforbothTB1models(nearestneighborhoppingonly)
Rev. B 77, (2008) 115109
and TB2 (with hopping beyond nearest neighbors), ex-
15. J. P. Robinson, H. Schomerus, L. Oroszla´ny, and V. I.
cepted for energyin the midgap states forTB1. This is in Falko, Phys. Rev. Lett. 101, 196803 (2008)
good agreement with two-dimensional Anderson localiza- 16. S.Yuan,H.DeRaedt,M.I.Katsnelson,Phys.Rev.B82,
tion and consistent with the expected universal behavior 115448 (2010)
of conductivity of a two dimensional disordered system 17. S.Yuan,T.O.Wehling,A.I.Lichtenstein,M.I.Katsnel-
[40]. son, Phys. Rev. Lett. 109, 156601 (2012)
A. Missaoui et al.: Numerical analysis of conductivity in bilayer graphene 9
18. A.Lherbier,S.M.-M.Dubois,X.Declerck,Y.-M.Niquet,
S. Roche, J.-C. Charlier, Phys. Rev. B 86, 075402 (2012)
19. S. Roche, N. Leconte, F. Ortmann, A. Lherbier, D. Sori-
ano,andJ.-C.Charlier,SolidStatesComm.152,1404(2012)
20. A.Cresti,F.Ortmann,T.Louvet,D.VanTuan,S.Roche,
Phys. Rev. Lett. 110, 196601 (2013)
21. A. Kretinin, G. L. Yu, R. Jalil, Y. Cao, F. Withers, A.
Mishchenko,M.I.Katsnelson,K.S.Novoselov,A.K.Geim,
F. Guinea, Phys. Rev. B 88, 165427 (2013)
22. G. Trambly de Laissardi`ere, D. Mayou, Modern Physics
Lett. B 25, 1019 (2011)
23. G. Trambly de Laissardi`ere, D. Mayou, Phys. Rev. Lett.
111, 146601 (2013)
24. G. Trambly de Laissardi`ere, D. Mayou, Adv. Nat. Sci.:
Nanosci. Nanotechnol. 5, 015007 (2014)
25. P.-L.Zhao,S.Yuan,M.I.Katsnelson,H.DeRaedt,Phys.
Rev. B 92, 045437 (2015)
26. M. Koshino, T. Ando, Phys. Rev. B 73, 245403 (2006)
27. S. Adam, S. Das Sarma, Phys. Rev. B 77, 115436 (2008)
28. M. Koshino, New J. Phys. 11, 095010 (2009)
29. S.Yuan,H.DeRaedt,M.I.Katsnelson,Phys.Rev.B82,
235409 (2010)
30. J.W.Gonza´lez,H.Santos,M.Pacheco,L.Chico,L.Brey,
Phys. Rev. B 81, 195406 (2010)
31. A. Ferreira, J. Viana-Gomes, J. Nilsson, E. R. Mucciolo,
N.M.R.Peres,A.H.CastroNeto,Phys.Rev.B83,165402
(2011)
32. D. Van Tuan, S. Roche, Phys. Rev. B 93, 041403 (2016)
33. R. R. Nair, M. Sepioni, I.-L. Tsai, O. Lehtinen, J.
Keinonen, A. V. Krasheninnikov, T. Thomson, A. K. Geim,
and I. V. Grigorieva, Nature Phys. 8, 199 (2012).
34. W.L.Scopel,W.S.Paz,J.C.C.Freitas,Sol.Stat.Comm.
240, 5 (2016)
35. I. Snyman, C. W. J. Beenakker, Phys. Rev. B 75, 045322
(2007)
36. G. Trambly de Laissardi`ere, D. Mayou, L. Magaud, Phys.
Rev. B 86, 125413 (2012)
37. G. Trambly de Laissardi`ere, D. Mayou, L. Magaud, Nano
Lett. 10, 804 (2010)
38. G. Trambly de Laissardi`ere, O. F. Namarvar, D. Mayou,
L. Magaud, Phys. Rev. B 93, 235135 (2016)
39. G.TramblydeLaissardi`ere,D.Mayou,C.R.Physique15,
70 (2014)
40. P. A. Lee, T. V. Ramakrishnan, Rev. Mod. Phys. 57, 287
(1985)
41. P. W. Brouwer, E. Racine, A. Furusaki, Y. Hatsugai, Y.
Morita, and C. Mudry, Phys. Rev. B 66, 014204 (2002)
42. D. Mayou, Europhysics Lett. 6, 549 (1988)
43. D. Mayou and S.N. Khanna, J. Phys. I France 5, 1199
(1995)
44. S.RocheandD.Mayou,Phys.Rev.Lett.79,2518(1997)
45. S. Roche and D. Mayou Phys. Rev. B 60, 322 (1999)
46. F. Triozon, J. Vidal, R. Mosseri, D. Mayou, Phys. Rev. B
65, 220202 (2002)
