Table Of ContentBoundary conditions forDiracfermions onaterminated honeycomb lattice
A. R. Akhmerov and C. W. J. Beenakker
Instituut-Lorentz, Universiteit Leiden, P.O.Box 9506, 2300 RALeiden, TheNetherlands
(Dated:October,2007)
WederivetheboundaryconditionfortheDiracequationcorrespondingtoatight-bindingmodelonatwo-
dimensionalhoneycomblatticeterminatedalonganarbitarydirection.Zigzagboundaryconditionsresultgener-
icallyonce theboundary isnot parallel tothebonds. Sinceahoneycomb stripwithzigzagedges isgapless,
thisimpliesthatconfinementbylatticeterminationdoesnotingeneralproduceaninsulatingnanoribbon. We
8
considertheopeningofagapinagraphenenanoribbonbyastaggeredpotentialattheedgeandderivethecor-
0
respondingboundaryconditionfortheDiracequation.Weanalyzetheedgestatesinananoribbonforarbitrary
0
boundaryconditionsandidentifyaclassofpropagatingedgestatesthatcomplementtheknownlocalizededge
2
statesatazigzagboundary.
n
a PACSnumbers:73.21.Hb,73.22.Dj,73.22.-f,73.63.Bd
J
1
1 I. INTRODUCTION dencecharacteristicofinsulatingbehaviorrequiresthespecial
armchairorientation(ϕamultipleof60◦),atwhichthedecay
l] Theelectronicpropertiesof graphenecanbe describedby ratef(ϕ)vanishes.
al a differenceequation (representinga tight-bindingmodelon Confinement by a mass term in the Dirac equation does
h a honeycomb lattice) or by a differential equation (the two- produceanexcitationgapregardlessoftheorientationofthe
s- dimensionalDiracequation)[1,2]. Thetwodescriptionsare boundary. We show how the infinite-mass boundary condi-
e equivalent at large length scales and low energies, provided tion of Ref. [8] can be approached starting from the zigzag
m the Dirac equation is supplemented by boundary conditions boundary condition, by introducing a local potential differ-
. consistentwiththetight-bindingmodel.Theseboundarycon- enceonthetwo sublatticesin thetight-bindingmodel. Such
t
a ditionsdependon a variety of microscopicproperties, deter- a staggered potential follows from atomistic calculations[3]
m minedbyatomisticcalculations[3]. andmaywellbetheoriginoftheinsulatingbehaviorobserved
experimentallyingraphenenanoribbons[9,10].
- For a general theoretical description, it is useful to know
d whatboundaryconditionsonthe Dirac equationare allowed Theoutlineofthispaperisasfollows.InSec.IIweformu-
n by the basic physical principles of current conservation and late,followingRefs.[4,5],thegeneralboundaryconditionof
o
(presence or absence of) time reversal symmetry — inde- theDiracequationonwhichouranalysisisbased. InSec.III
c
[ pendently of any specific microscopic input. This problem we derivefrom the tight-bindingmodelthe boundarycondi-
was solved in Refs. [4, 5]. The general boundary condi- tioncorrespondingtoanarbitrarydirectionoflatticetermina-
3
tion depends on one mixing angle Λ (which vanishes if the tion.InSec.IVweanalyzetheeffectofastaggeredboundary
v
3 boundarydoesnotbreaktimereversalsymmetry),onethree- potentialon the boundarycondition. In Sec. V we calculate
2 dimensionalunitvectornperpendiculartothenormaltothe the dispersion relation for a graphene nanoribbon with arbi-
7 boundary, and one three-dimensional unit vector ν on the trary boundary conditions. We identify dispersive (= propa-
2 Blochsphereofvalleyisospins. Altogether,fourrealparame- gating)edgestateswhichgeneralizetheknowndispersionless
0. tersfixtheboundarycondition. (= localized)edge states at a zigzag boundary[11]. The ex-
1 Inthepresentpaperweinvestigatehowtheboundarycondi- ponentialdependenceof the gap∆ onthe nanoribbonwidth
7 tiondependsonthecrystallographicorientationofthebound- iscalculatedinSec.VIbothanalyticallyandnumerically.We
0 ary. As the orientation is incremented by 30◦ the boundary concludeinSec.VII.
:
v configurationswitchesfromarmchair(paralleltoone-thirdof
i thecarbon-carbonbonds)tozigzag(perpendiculartoanother
X
one-thirdofthebonds).Theboundaryconditionsforthearm- II. GENERALBOUNDARYCONDITION
r chair and zigzag orientations are known [6]. Here we show
a
that the boundarycondition for intermediate orientationsre- Thelong-wavelengthandlow-energyelectronicexcitations
mainsofthezigzagform,sothatthearmchairboundarycon- ingraphenearedescribedbytheDiracequation
ditionisonlyreachedforadiscretesetoforientations.
Sincethezigzagboundaryconditiondoesnotopenupagap HΨ=εΨ (2.1)
in the excitation spectrum [6], the implication of our result
(notnoticedinearlierstudies[7])isthataterminatedhoney- withHamiltonian
comblatticeofarbitraryorientationismetallicratherthanin-
sulating.Wepresenttight-bindingmodelcalculationstoshow H =vτ (σ p) (2.2)
0
⊗ ·
that, indeed, the gap∆ exp[ f(ϕ)W/a] in a nanoribbon
∝ −
atcrystallographicorientationϕvanishesexponentiallywhen actingona four-componentspinorwave functionΨ. Here v
its width W becomes large compared to the lattice constant is the Fermi velocity and p = i~ is the momentum op-
− ∇
a,characteristicofmetallicbehavior. The∆ 1/W depen- erator. Matricesτ ,σ are Paulimatricesin valleyspace and
i i
∝
2
sublatticespace,respectively(withunitmatricesτ ,σ ). The III. LATTICETERMINATIONBOUNDARY
0 0
currentoperatorinthedirectionnisn J =vτ (σ n).
0
· ⊗ ·
The Hamiltonian H is written in the valley isotropic rep- Thehoneycomblatticeofacarbonmonolayerisatriangu-
resentation of Ref. [5]. The alternative representationH′ = larlattice(latticeconstanta)withtwoatomsperunitcell,re-
vτz (σ p)ofRef. [4]is obtainedbytheunitarytransfor- ferredtoasAandBatoms(seeFig.1a). TheAandBatoms
⊗ ·
mation separately form two triangular sublattices. The A atoms are
connectedonlytoBatoms,andviceversa. Thetight-binding
H′ =UHU†, U = 1(τ +τ ) σ +1(τ τ ) σ . (2.3)
2 0 z ⊗ 0 2 0− z ⊗ z equationsonthehoneycomblatticearegivenby
As described in Ref. [4], the general energy-independent
εψ (r)=t[ψ (r)+ψ (r R )+ψ (r R )],
boundaryconditionhastheformofalocallinearrestrictionon A B B 1 B 2
− − (3.1)
thecomponentsofthespinorwavefunctionattheboundary: εψB(r)=t[ψA(r)+ψA(r+R1)+ψA(r+R2)].
Ψ=MΨ. (2.4) Here t is the hopping energy, ψ (r) and ψ (r) are the
A B
electron wave functions on A and B atoms belonging to
The 4 4 matrix M has eigenvalue 1 in a two-dimensional
× the same unit cell at a discrete coordinate r, while R1 =
subspace containing Ψ, and without loss of generality we
(a√3/2, a/2), R = (a√3/2,a/2) are lattice vectors as
2
mayassumethatM haseigenvalue 1intheorthogonaltwo- −
showninFig.1a.
−
dimensionalsubspace. ThismeansthatM maybechosenas
Regardlessofhowthelatticeisterminated,Eq.(3.1)hasthe
aHermitianandunitarymatrix,
electron-holesymmetryψ ψ ,ε ε. Forthelong-
B B
→ − → −
M =M†, M2 =1. (2.5) wavelength Dirac Hamiltonian (2.2) this symmetry is trans-
latedintotheanticommutationrelation
Therequirementofabsenceofcurrentnormaltothebound-
ary, Hσz τz +σz τzH =0. (3.2)
⊗ ⊗
Ψn J Ψ =0, (2.6) Electron-holesymmetryfurtherrestrictstheboundarymatrix
B
h | · | i
M inEq.(2.10)totwoclasses: zigzag-like(ν = zˆ,n=zˆ)
with nB a unit vector normal to the boundary and pointing and armchair-like (ν = n = 0). In this sect±ion we will
z z
outwards,isequivalenttotherequirementofanticommutation showthatthe zigzag-likeboundaryconditionappliesgeneri-
ofthematrixM withthecurrentoperator, callytoanarbitraryorientationofthelatticetermination.The
armchair-likeboundaryconditionis only reachedfor special
M,n J =0. (2.7)
{ B · } orientations.
ThatEq.(2.7)impliesEq.(2.6)followsfrom Ψn J Ψ =
B
h | · | i
ΨM(n J)M Ψ = Ψn J Ψ . The converse is
B B
hpro|veninA·pp.A.| i −h | · | i A. Characterizationoftheboundary
Wearenowfacedwiththeproblemofdeterminingthemost
general4 4matrixM thatsatisfiesEqs.(2.5)and(2.7). Ref. Aterminatedhoneycomblatticeconsistsofsiteswiththree
×
[4]obtainedtwofamiliesoftwo-parametersolutionsandtwo neighborsintheinteriorandsiteswithonlyoneortwoneigh-
more families of three-parameter solutions. These solutions bors at the boundary. The absent neighboring sites are in-
are subsets of the single four-parameter family of solutions dicated by open circles in Fig. 1 and the dangling bonds by
obtainedinRef.[5], thin line segments. The tight-binding model demands that
the wave function vanishes on the set of absent sites, so the
M =sinΛτ (n σ)+cosΛ(ν τ) (n σ), (2.8)
0 1 2 firststepinouranalysisisthecharacterizationofthisset. We
⊗ · · ⊗ ·
assume that the absent sites form a one-dimensional super-
whereν,n ,n arethree-dimensionalunitvectors,suchthat
1 2 lattice, consisting ofa supercellof N emptysites, translated
n and n are mutually orthogonal and also orthogonal to
1 2 overmultiplesofasuperlatticevectorT. Sincetheboundary
n . Aproofthat(2.8)isindeedthemostgeneralsolutionis
B
superlattice is part of the honeycomb lattice, we may write
giveninApp.A. OnecanalsocheckthatthesolutionsofRef.
[4]aresubsetsofM′ =UMU†. T = nR1 +mR2 withnandm non-negativeintegers. For
example, in Fig. 1 we have n = 1, m = 4. Without loss
In this work we will restrict ourselvesto boundarycondi-
ofgenerality,andforlaterconvenience,wemayassumethat
tionsthatdonotbreaktimereversalsymmetry. Thetimere-
m n=0(modulo3).
versaloperatorinthevalleyisotropicrepresentationis −
The angle ϕ between T and the armchair orientation (the
T = (τ σ ) , (2.9) x-axisinFig.1)isgivenby
y y
− ⊗ C
with the operator of complex conjugation. The boundary 1 n m π π
C ϕ=arctan − , ϕ . (3.3)
condition preserves time reversal symmetry if M commutes (cid:18)√3n+m(cid:19) −6 ≤ ≤ 6
withT. ThisimpliesthatthemixingangleΛ = 0,sothatM
isrestrictedtoathree-parameterfamily, The armchair orientation correspondsto ϕ = 0, while ϕ =
π/6correspondstothezigzagorientation. (Becauseofthe
M =(ν τ) (n σ), n nB. (2.10) ±π/3periodicityweonlyneedtoconsider ϕ π/6.)
· ⊗ · ⊥ | |≤
3
with ~k = p T. While the continuous quantum number
·
k (0,2π) describes the propagation along the boundary,
∈
a second (discrete) quantum number λ describes how these
boundarymodesdecayawayfromtheboundary. Weselectλ
bydemandingthattheBlochwave(3.4)isalsoasolutionof
ψ(r+R )=λψ(r). (3.5)
3
ThelatticevectorR = R R hasanonzerocomponent
3 1 2
−
acosϕ > a√3/2 perpendicular to T. We need λ 1 to
| | ≤
preventψ(r)fromdivergingintheinteriorofthelattice. The
decaylengthl inthedirectionperpendiculartoT isgiven
decay
by
acosϕ
l = − . (3.6)
decay
ln λ
| |
FIG. 1: (a) Honeycomb latice constructed from a unit cell (grey TheboundarymodessatisfyingEqs.(3.4)and(3.5)arecal-
rhombus) containing twoatoms(labeledAandB), translatedover culatedin App.B fromthe tight-bindingmodel. In the low-
latticevectors R1 and R2. Panels b,c,d show three different peri- energy regime of interest (energies ε small compared to t)
odicboundarieswiththesameperiodT =nR1+mR2.Atomson there is an independentset of modeson each sublattice. On
the boundary (connected by thick solidlines) have dangling bonds
sublatticeAthequantumnumbersλandkarerelatedby
(thindottedlinesegments)toemptyneighboringsites(opencircles).
The number N of missing sites and N′ of dangling bonds per pe- ( 1 λ)m+n =exp(ik)λn (3.7a)
riod is ≥ n+m. Panel d shows a minimal boundary, for which − −
N =N′ =n+m.
andonsublatticeBtheyarerelatedby
( 1 λ)m+n =exp(ik)λm. (3.7b)
The number N of empty sites per period T can be arbi- − −
trarily large, but it cannot be smaller than n + m. Like- For a given k there are roots λ of Eq. (3.7a) having
A p
wise, the number N′ of dangling bonds per period cannot absolutevalue 1,withcNorrespondingboundarymodesψp.
≤
be smaller than n + m. We call the boundary minimal if We sort these modes according to their decay lengths from
N = N′ = n+m. For example, the boundary in Fig. 1d short to long, ldecay(λp) ldecay(λp+1), or λp λp+1 .
≤ | | ≤ | |
′ ThewavefunctiononsublatticeAisasuperpositionofthese
isminimal(N = N = 5), whilethe boundariesin Figs.1b
′ ′ modes
and1carenotminimal(N =7,N =9andN =5,N =7,
respectively). In what follows we will restrict our consider- NA
ations to minimal boundaries, both for reasons of analytical ψ(A) = α ψ , (3.8)
p p
simplicity[12]andforphysicalreasons(itisnaturaltoexpect Xp=1
thattheminimalboundaryisenergeticallymostfavorablefor
agivenorientation). withcoefficientsα suchthatψ(A) vanishesontheN miss-
p A
′
We conclude this subsection with a property of minimal ingAsites. Similarlythereare rootsλ ofEq.(3.7b)with
boundaries that we will need later on. The N empty sites λ′ 1, λ′ λ′ . ThecNorBrespondinpgboundarymodes
per period can be divided into NA empty sites on sublattice |forpm|≤thew|avpe|≤fun|ctpi+on1|onsublatticeB,
AandN emptysitesonsublatticeB. Aminimalboundary
B
isconstructedfromntranslationsoverR ,eachcontributing NB
oneemptyAsite,andmtranslationsover1R2,eachcontribut- ψ(B) = α′pψp′, (3.9)
ing one empty B site. Hence, N = n and N = m for a Xp=1
A B
minimalboundary.
withα′ suchthatψ(B)vanishesontheN missingBsites.
p B
B. Boundarymodes C. Derivationoftheboundarycondition
Theboundarybreaksthetwo-dimensionaltranslationalin- To derivetheboundaryconditionfortheDirac equationit
varianceoverR1andR2,butaone-dimensionaltranslational is sufficient to consider the boundary modes in the k 0
invariance over T = nR1 + mR2 remains. The quasimo- limit. The characteristic equations (3.7) for k = 0→each
mentum p along the boundary is therefore a good quantum have a pair of solutions λ = exp( 2iπ/3) that do not
±
±
number.ThecorrespondingBlochstatesatisfies depend on n and m. Since λ = 1, these modes do
±
| |
not decay as one moves away from the boundary. The cor-
ψ(r+T)=exp(ik)ψ(r), (3.4) responding eigenstate exp( iK r) is a plane wave with
± ·
4
wave vector K = (4/3)πR /a2. One readily checks that D. Precisionoftheboundarycondition
3
this Bloch state also satisfies Eq. (3.4) with k = 0 [since
K T =2π(n m)/3=0(modulo2π)]. Ataperfectzigzagorarmchairedgethefourcomponentsof
· −
Thewavefunctions(3.8)and(3.9)onsublatticesAandB theDiracspinorΨaresufficienttomeettheboundarycondi-
inthelimitk 0taketheform tion.Neartheboundarieswithlargerperiodandmorecompli-
→
catedstructurethewavefunction(3.10)alsonecessarilycon-
NA−2 ′
ψ(A) =Ψ eiK·r+Ψ e−iK·r+ α ψ , (3.10a) tainsseveralboundarymodesψp,ψpthatdecayawayfromthe
1 4 p p boundary.Thedecaylengthδoftheslowestdecayingmodeis
Xp=1
the distance at which the boundaryis indistinguishablefrom
ψ(B) =Ψ eiK·r+Ψ e−iK·r+NB−2α′ψ′. (3.10b) aperfectarmchairorzigzagedge.Atdistancessmallerthanδ
2 3 p p theboundaryconditionbreaksdown.
Xp=1 Inthecaseofanarmchair-likeboundary(withn= m),all
′
thecoefficientsα andα inEqs. (3.10)mustbenonzeroto
The four amplitudes(Ψ , iΨ , iΨ , Ψ ) Ψ form the p p
1 − 2 3 − 4 ≡ satisfy theboundarycondition. Themaximaldecaylengthδ
four-componentspinorΨintheDiracequation(2.1). There-
isthenequaltothedecaylengthoftheboundarymodeψ
maining 2and 2termsdescribedecayingboundary n−1
NA− NB− whichhas the largest λ. Itcan be estimated fromthe char-
modesofthetight-bindingmodelthatarenotincludedinthe | |
acteristic equations(3.7) that δ T . Hence the larger the
Diracequation. ≈ | |
period of an armchair-like boundary, the larger the distance
We are now ready to determine what restriction on Ψ is
from the boundary at which the boundary condition breaks
imposedby the boundaryconditionon ψ(A) and ψ(B). This
down.
restriction is the required boundary condition for the Dirac
Forthezigzag-likeboundarythesituationis different. On
equation.InApp.Bwecalculatethat,fork =0,
onesublatticetherearemoreboundarymodesthanconditions
imposed by the presence of the boundary and on the other
=n (n m)/3+1, (3.11)
NA − − sublatticetherearelessboundarymodesthanconditions. Let
B =m (m n)/3+1, (3.12) usassume thatsublattice A hasmore modesthan conditions
N − −
(which happens if n < m). The quickest decaying set of
sothat A+ B =n+m+2isthetotalnumberofunknown boundarymodessufficienttosatisfythetight-bindingbound-
N N
amplitudesin Eqs. (3.8) and (3.9). These have to be chosen
ary condition contains n modes ψ with p n. The dis-
p
such that ψ(A) and ψ(B) vanish on NA and NB lattice sites tance δ from the boundary within which the≤boundary con-
respectively. For the minimal boundaryunder consideration
dition breaks down is then equal to the decay length of the
we haveNA = n equationstodetermine A unknownsand slowestdecayingmodeψ inthissetandisgivenby
N n
N =mequationstodetermine unknowns.
B B
N
Three cases can be distinguished [in each case n m = δ =l (λ )= acosϕ/ln λ . (3.13)
decay n n
− − | |
0(modulo3)]:
[SeeEq.(3.6).]
1. If n > m then n and m+2, so Ψ = As derivedin App. B for the case of large periods T
Ψ4 =0,whileΨN2Aan≤dΨ3areNunBde≥termined. 1 a, the quantum number λn satisfies the following sys|tem| ≫of
equations:
2. If n < m then n and m+2, so Ψ =
NB ≤ NA ≥ 2 1+λ m+n = λ n, (3.14a)
Ψ =0,whileΨ andΨ areundetermined. n n
3 1 4 | | | |
n n
arg(1+λ ) arg( λ )= π. (3.14b)
n n
3. If n = m then = n + 1 and = m + 1, so − n+m − n+m
A B
N N
Ψ = Ψ and Ψ = Ψ .
| 1| | 4| | 2| | 3| Thesolutionλnofthisequationandhencethedecaylengthδ
donotdependonthelength T oftheperiod,butonlyonthe
Ineachcasetheboundaryconditionisofthecanonicalform | |
ration/(n+m)=(1 √3tanϕ)/2,whichisafunctionof
Ψ=(ν τ) (n σ)Ψwith
−
· ⊗ · theangleϕbetweenT andthearmchairorientation[seeEq.
1. ν = zˆ, n = zˆ if n > m (zigzag-type boundary (3.3)]. Inthecasen >mwhensublatticeB hasmoremodes
condit−ion). thanconditions,thelargestdecaylengthδfollowsuponinter-
changingnandm.
2. ν = zˆ,n = zˆifn < m(zigzag-typeboundarycondi- AsseenfromFig.2, theresultingdistanceδ withinwhich
tion). the zigzag-type boundary condition breaks down is zero for
thezigzagorientation(ϕ = π/6)andtendsto infinityasthe
3. ν zˆ=0,n zˆ=0ifn=m(armchair-typeboundary orientationoftheboundaryapproachesthearmchairorienta-
co·ndition). · tion (ϕ = 0). (For finite periodsthe divergenceis cut off at
δ T a.)Theincreaseofδnearthearmchairorientation
∼| |≫
Weconcludethattheboundaryconditionisofzigzag-typefor isratherslow: Forϕ & 0.1thezigzag-typeboundarycondi-
anyorientationT oftheboundary,unlessT isparalleltothe tionremainspreciseonthescaleofafewunitcellsawayfrom
bonds[sothatn=mandϕ=0(moduloπ/3)]. theboundary.
5
Thedensityofedgestatesismaximalρ = 1/3aforaperfect
zigzagedgeanditdecreasescontinuouslywhentheboundary
orientation ϕ approaches the armchair one. Eq. (3.16) ex-
plainsthenumericaldataofRef.[11],providingananalytical
formulaforthedensityofedgestates.
IV. STAGGEREDBOUNDARYPOTENTIAL
The electron-hole symmetry (3.2), which restricts the
boundary condition to being either of zigzag-type or of
armchair-type,isbrokenbyanelectrostaticpotential.Herewe
consider,motivatedbyRef. [3], the effectofastaggeredpo-
tentialatthezigzagboundary.Weshowthattheeffectofthis
potentialistochangetheboundaryconditioninacontinuous
wayfromΨ= τ σ ΨtoΨ= τ (σ [zˆ n ])Ψ.The
z z z B
± ⊗ ± ⊗ · ×
first boundary condition is of zigzag-type, while the second
boundary condition is produced by an infinitely large mass
termattheboundary[8].
FIG.2: Dependenceontheorientationϕofthedistanceδfromthe
The staggered potential consists of a potential V = +µ,
boundary within which the zigzag-type boundary condition breaks A
down. Thecurveiscalculatedfromformula(3.14)validinthelimit VB = µ on the A-sites and B-sites in a total of 2N rows
−
|T|≫aoflargeperiods. Theboundaryconditionbecomesprecise closest to the zigzag edge parallel to the y-axis (see Fig. 3).
uponapproachingthezigzagorientationϕ=π/6. Since this potential does not mix the valleys, the boundary
conditionnearazigzagedgewithstaggeredpotentialhasthe
form
Although the presented derivation is only valid for peri-
odic boundaries and low energies, such that the wavelength Ψ= τ (σ cosθ+σ sinθ)Ψ, (4.1)
z z y
− ⊗
ismuchlargerthanthelength T oftheboundaryperiod,we
arguethattheseconditionsma|yb|erelaxed. Indeed,sincethe inaccordwiththegeneralboundarycondition(2.10). Forθ =
boundaryconditionislocal,itcannotdependonthestructure 0,πwehavethezigzagboundaryconditionandforθ = π/2
±
oftheboundaryfaraway,hencetheperiodicityofthebound- wehavetheinfinite-massboundarycondition.
ary cannotinfluence the boundarycondition. It can also not To calculate the angle θ we substitute Eq. (3.10) into the
dependonthewavelengthoncethewavelengthislargerthan tight-binding equation (3.1) (including the staggered poten-
thetypicalsizeofaboundaryfeature(ratherthanthelengthof tialattheleft-handside)andsearchforasolutioninthelimit
theperiod).Sinceformostboundariesbothδandthescaleof ε = 0. Theboundaryconditionispreciseforthezigzagori-
theboundaryroughnessareoftheorderofseveralunitcells, entation, so we may set αp = α′p = 0. It is sufficient to
weconcludethatthezigzagboundaryconditionisingeneral consider a single valley, so we also set Ψ3 = Ψ4 = 0. The
agoodapproximation. remaining nonzero components are Ψ1eiK·r ψA(i)eiKy
andΨ eiK·r ψ (i)eiKy,whereiintheargu≡mentofψ
2 B A,B
≡
numbersthe unit cell away from the edge and we have used
E. Densityofedgestatesnearazigzag-likeboundary that K points in the y-direction. The resulting difference
equationsare
A zigzag boundaryis known to support a band of disper-
µψ (i)=t[ψ (i) ψ (i 1)], i=1,2,...N, (4.2a)
A B B
sionless states [11], which are localized within several unit − − −
µψ (i)=t[ψ (i) ψ (i+1)], i=0,1,2,...N 1,
cellsneartheboundary. Wecalculatethe1Ddensityofthese B A A
− −
(4.2b)
edgestatesnearanarbitraryzigzag-likeboundary. Againas-
suming that the sublattice A hasmore boundarymodesthan ψ (0)=0. (4.2c)
A
conditions(n < m), for each k there are (k) N lin-
A A
earlyindependentstates(3.8),satisfyingthNebound−arycondi- For the Ψ1,Ψ2 components of the Dirac spinor Ψ the
boundarycondition(4.1)isequivalentto
tion. For k = 0 the number of boundary modes is equal to
6
=n (m n)/3,sothatforeachkthereare
NA − − ψA(N)/ψB(N)= tan(θ/2). (4.3)
−
N = (k) n=(m n)/3 (3.15)
states NA − − SubstitutingthesolutionofEq.(4.2)intoEq.(4.3)gives
edgestates. Thenumberoftheedgestatesforthecasewhen
1+sinh(κ)sinh(κ+2Nµ/t)
n > magainfollowsuponinterchangingnandm. Theden- cosθ = , (4.4)
sityρofedgestatesperunitlengthisgivenby cosh(κ)cosh(κ+2Nµ/t)
Nstates m n 2 withsinhκ=µ/2t.Eq.(4.4)isexactforN 1,butitisac-
ρ= = | − | = sinϕ. (3.16) ≫
T 3a√n2+nm+m2 3a| | curatewithin2%foranyN.Thedependenceoftheparameter
| |
6
θoftheboundaryconditiononthestaggeredpotentialstrength Inthissectionweconsiderthemostgeneralboundarycon-
µ is shown in Fig. 4 forvariousvaluesof N. The boundary dition(2.10),constrainedonlybytime-reversalsymmetry.We
conditionisclosesttotheinfinitemassforµ/t 1/N,while donotrequirethattheboundaryispurelyaterminationofthe
∼
theregimesµ/t 1/N orµ/t 1correspondtoa zigzag lattice,butallowforarbitrarylocalelectricfieldsandstrained
≪ ≫
boundarycondition. bonds.TheconclusionofSec.III,thattheboundarycondition
iseitherzigzag-likeorarmchair-like,doesnotapplytherefore
totheanalysisgiveninthissection.
The general solution of the Dirac equation (2.1) in the
nanoribbonhastheformΨ(x,y)=Ψ (x)eiky. Weimpose
n,k
thegeneralboundarycondition(2.10),
Ψ(0,y)=(ν τ) (n σ)Ψ(0,y), (5.1a)
1 1
· ⊗ ·
Ψ(W,y)=(ν τ) (n σ)Ψ(W,y), (5.1b)
2 2
· ⊗ ·
withthree-dimensionalunitvectorsν , n , restrictedbyn
i i i
·
xˆ = 0 (i = 1,2). (Thereisnorestrictiononthe ν .) Valley
i
isotropyoftheDiracHamiltonian(2.2)impliesthatthespec-
FIG.3: ZigzagboundarywithV = +µontheA-sites(filleddots) trumdoesnotdependonν andν separatelybutonlyonthe
1 2
and V = −µ on the B-sites (empty dots). The staggered poten- angleγ between them. Thespectrumdepends, therefore,on
tialextendsover2N rowsofatomsnearesttothezigzagedge. The threeparameters: Theangleγ andtheanglesθ ,θ between
1 2
integericountsthenumberofunitcellsawayfromtheedge. thez-axisandthevectorsn ,n .
1 2
The Dirac equation HΨ = εΨ has two plane wave solu-
tionsΨ exp(iky+iqx)foragivenεandk,corresponding
∝
tothetwo(realorimaginary)transversewavenumbersqthat
solve (~v)2(k2 +q2) = ε2. Each of these two plane waves
has a twofold valley degeneracy, so there are four indepen-
dentsolutionsintotal. Sincethewavefunctioninaribbonisa
linearcombinationofthesefourwaves,andsinceeachofthe
Eqs.(5.1a,5.1b)hasatwo-dimensionalkernel,theseequations
providefourlinearlyindependentequationstodeterminefour
unknowns.TheconditionthatEq.(5.1)hasnonzerosolutions
gives an implicit equation for the dispersion relation of the
nanoribbon:
cosθ cosθ (cosω cos2Ω)+cosωsinθ sinθ sin2Ω
1 2 1 2
−
sinΩ[sinΩcosγ+sinωsin(θ θ )]=0, (5.2)
1 2
− −
whereω2 =4W2[(ε/~v)2 k2]andcosΩ=~vk/ε.
−
For θ = θ = 0 and γ = π Eq. (5.2) reproduces the
1 2
transcendental equation of Ref. 6 for the dispersion relation
ofazigzagribbon.Inthecaseθ =θ =π/2ofanarmchair-
1 2
FIG.4: Plotoftheparameterθintheboundarycondition(4.1)ata
likenanoribbon,Eq.(5.2)simplifiesto
zigzagedgewiththestaggeredpotentialofFig.3.Thecurvesarecal-
culatedfromEq.(4.4). Thevaluesθ = 0andθ = π/2correspond,
cosω =cosγ. (5.3)
respectively,tothezigzagandinfinite-massboundaryconditions.
This is the only case when the transverse wave function
Ψ (x) is independent of the longitudinal wave number k.
n,k
InFig.5weplotthedispersionrelationsforseveraldifferent
V. DISPERSIONRELATIONOFANANORIBBON boundaryconditions.
The low energy modes of a nanoribbon with ε < ~v k
| | | |
Agraphenenanoribbonisacarbonmonolayerconfinedto [seepanelsa-dofFig.5]haveimaginarytransversemomen-
alongandnarrowstrip. Theenergyspectrumε (k)ofthen- tum since q2 = (ε/~v)2 k2 < 0. If q becomes larger
n
− | |
thtransversemodeisafunctionofthewavenumberk along than the ribbon width W, the corresponding wave function
thestrip. Thisdispersionrelationisnonlinearbecauseofthe becomeslocalizedattheedgesofthenanoribbonanddecays
confinement,which also may openup a gapin the spectrum inthebulk.Thedispersionrelation(5.2)forsuchanedgestate
aroundzeroenergy. We calculatethe dependenceofthe dis- simplifiestoε=~v k sinθ forthestatelocalizednearx=0
1
persionrelationonthe boundaryconditionsatthetwo edges and ε = ~v k sin|θ| for the state localized near x = W.
2
− | |
x=0andx=W ofthenanoribbon(takenalongthey-axis). These dispersive edge states with velocity vsinθ generalize
7
opposite zigzag edges have the same staggered potential, so
thattheboundaryconditionis
Ψ(0,y)=+τ (σ cosθ+σ sinθ)Ψ(0,y), (5.4a)
z z y
⊗
Ψ(W,y)= τ (σ cosθ+σ sinθ)Ψ(W,y). (5.4b)
z z y
− ⊗
Thedependenceofθontheparametersµ,N ofthestaggered
potentialisgivenbyEq.(4.4). Thisboundaryconditioncor-
respondstoγ = π,θ = θ = θ, sothatithasagapforany
1 2
nonzeroθ. AsshowninFig.6,∆(θ)increasesmonotonically
with θ from the zigzag limit ∆(0) = 0 to the infinite-mass
limit∆(π/2)=π~v/W.
3
]
W2
/
v
h¯
[
∆ 1
0
0 π/4 π/2
θ
FIG. 6: Dependence of the band gap ∆ on the parameter θ in the
staggeredpotentialboundarycondition(5.4).
VI. BANDGAPOFATERMINATEDHONEYCOMB
LATTICE
Inthissectionwe returntothecase ofaboundaryformed
purelybyterminationofthelattice.Ananoribbonwithzigzag
boundaryconditionhaszerobandgapaccordingtotheDirac
equation(Fig. 5a). Accordingto the tight-bindingequations
FIG.5: Dispersionrelationofnanoribbonswithdifferentboundary
conditions. Thelarge-wavenumberasymptotes|ε| = ~v|k|ofbulk thereis a nonzerogap∆, which howevervanishesexponen-
statesareshownbydashedlines. Modesthatdonotapproachthese tially with increasing width W of the nanoribbon. We esti-
asymptotesareedgestateswithdispersion|ε| = ~v|ksinθi|. The matethedecayrateof∆(W)asfollows.
zigzag ribbon with γ = π and θ1 = θ2 = 0 (a) exhibits disper- Thelow energystates in a zigzag-typenanoribbonare the
sionlessedgestatesatzeroenergy[11]. Ifθ1 orθ2 arenonzero(b, hybridizedzeroenergyedgestatesattheoppositeboundaries.
c)theedgestatesacquirelineardispersionandifsinθ1sinθ2 > 0 Theenergyεofsuchstatesmaybeestimatedfromtheoverlap
(c) aband gap opens. If γ isunequal to0or π (d) thevalleys are betweenthe edgestates localized atthe oppositeedges, ε =
mixedwhichmakesallthelevelcrossingsavoidedandopensaband (~v/W)exp( W/l ). In a perfectzigzag ribbon there
gap. Armchair-likeribbonswithθ1 = θ2 = π/2(e,f)aretheonly a±reedgestatesw−ithl decay=0(andε= 0),sothatthereisno
ribbonshavingnoedgestates. decay
bandgap. Fora ribbonwith a morecomplicatededgeshape
the decay length of an edge state is limited by δ, the length
within which the boundarycondition breaks down (see Sec.
theknown[11]dispersionlessedgestatesata zigzagbound-
III.D). This length scale provides the analytical estimate of
ary(withsinθ =0).
thebandgapinazigzag-likeribbon:
Inspectionofthedispersionrelation(5.2)givesthefollow-
ing condition for the presence of a gap in the spectrum of
theDirac equationwitharbitraryboundarycondition: Either ∆ ~ve−W/δ, (6.1)
the valleys should be mixed (γ = 0,π) or the edge states ∼ W
6
atoppositeboundariesshouldhaveenergiesofoppositesign
withδgivenbyEqs.(3.13)and(3.14).
(sinθ sinθ >0forγ =πorsinθ sinθ <0forγ =0).
1 2 1 2 Thebandgapofanarmchair-likeribbonis
Asanexample,wecalculatethebandgapforthestaggered
potentialboundaryconditionofSec.IV. We assumethatthe ∆=(~v/W)arccos(cosγ) (6.2)
8
2
[see Eq. (5.3) and panels e,f of Fig. 5]. Addinganotherrow
ofatomsincreasesthenanoribbonwidthbyonehalfofaunit
cellandincreasesγ byK R = 4π/3,sotheproduct∆W
3
·
insucharibbonisanoscillatoryfunctionofW withaperiod
of1.5unitcells.
To test these analytical estimates, we have calculated
∆(W) numerically for various orientations and configura-
)
ϕ
tions of boundaries. As seen from Fig. 7, in ribbons with 1
(
f
a non-armchair boundary the gap decays exponentially
∝
exp[ f(ϕ)W/a] as a functionof W. Nanoribbonswith the
−
sameorientationϕbutdifferentperiod T havethesamede-
| |
cay rate f. As seen in Fig. 8, the decay rate obtained nu-
merically agrees well with the analytical estimate f = a/δ
following from Eq. (6.1) (with δ given as a function of ϕ in
Fig. 2). The numerical results of Fig. 7 are consistent with
earlierstudiesofthe orientationdependenceofthebandgap 0 π/6
innanoribbons[7],buttheexponentialdecreaseofthegapfor ϕ
non-armchairribbonswasnotnoticedinthosestudies.
For completenesswe show in Fig. 9 ournumericalresults
FIG. 8: Dependence of the gap decay rate on the orientation ϕ of
forthebandgapinanarmchair-likenanoribbon(ϕ=0). We
the boundary (defined in the inset of Fig. 2). The dots are the fits
see that the gap oscillates with a period of 1.5 unit cells, in tonumericalresultsofthetight-bindingequations,thesolidcurveis
agreementwithEq.(6.2). theanalyticalestimate(6.1).
FIG.7:Dependenceofthebandgap∆ofzigzag-likenanoribbonson
thewidthW.Thecurvesintheleftpanelarecalculatednumerically FIG.9: Dependenceofthebandgap∆onthewidthW foranarm-
fromthetight-bindingequations.Therightpanelshowsthestructure chair ribbon (dashed line) and for a ribbon witha boundary of the
oftheboundary,repeatedperiodicallyalongbothedges. sameorientationbutwithalargerperiod(solidline).Thecurvesare
calculatednumericallyfromthetight-bindingequations.
VII. CONCLUSION ∆ (~v/W)exp( W/δ)inananoribbonofwidthW. We
≈ −
havetestedouranalyticalresultsfor∆withthenumericalso-
In summary, we have demonstrated that the zigzag-type lutionofthetight-bindingequationsandfindgoodagreement.
boundary condition Ψ = τ σ Ψ applies generically While the lattice termination by itself can only produce
z z
± ⊗
to a terminated honeycomb lattice. The boundary condition zigzag or armchair-type boundary conditions, other types of
switchesfromtheplus-signtotheminus-signatthearmchair boundaryconditionscanbereachedbybreakingtheelectron-
orientationϕ=0(moduloπ/3),whentheboundaryisparal- holesymmetryofthe tight-bindingequations. We havecon-
lelto1/3ofallthecarbon-carbonbonds(seeFig.10). sideredtheeffectofastaggeredpotentialatazigzagboundary
The distance δ from the edge within which the boundary (producedforexamplebyedgemagnetization[3]), andhave
conditionbreaksdownis minimal(= 0) atthe zigzagorien- calculatedthe correspondingboundarycondition. Itinterpo-
tation ϕ = π/6 (moduloπ/3) and maximal at the armchair latessmoothlybetweenthezigzagandinfinite-massboundary
orientation.Thisisthelengthscalethatgovernsthebandgap conditions,openingupagapinthespectrumthatdependson
9
thestrengthandrangeofthestaggeredpotential. Using only the Hermiticity of M, we have the 16-parameter
We have calculated the dispersion relation for arbitrary representation
boundaryconditionsandfoundthattheedgestateswhichare
3
dispersionlessatazigzagedgeacquireadispersionformore
M = (τ σ )c , (A2)
general boundary conditions. Such propagating edge states i⊗ j ij
exist, for example, near a zigzagedge with staggeredpoten- iX,j=0
tial. with real coefficients c . Anticommutationwith the current
ij
Ourdiscoverythatthezigzagboundaryconditionisgeneric operatorbringsthisdowntothe8-parameterform
explainsthefindingsofseveralcomputersimulations[11,13,
14] in which behavior characteristic of a zigzag edge was 3
observed at non-zigzag orientations. It also implies that the M = τi (ni σ), (A3)
⊗ ·
mechanismofgapopeningatazigzagedgeofRef.[3](pro- Xi=0
duction of a staggered potential by magnetization) applies
where the n ’s are three-dimensional vectors orthogonal to
i
generically to any ϕ = 0. This may explain why the band n . The absence of off-diagonalterms in M2 requires that
6 B
gap measurementsof Ref. [10] producedresults that did not
thevectorsn , n , n aremultiplesofaunitvectorn˜ which
1 2 3
dependonthecrystallographicorientationofthenanoribbon.
isorthogonalton . ThematrixM maynowberewrittenas
0
M =τ (n σ)+(ν˜ τ) (n˜ σ). (A4)
0 0
⊗ · · ⊗ ·
TheequalityM2 =1furtherdemandsn2+ν˜2 =1,leading
0
tothe4-parameterrepresentation(2.8)afterredefinitionofthe
vectors.
APPENDIXB:DERIVATIONOFTHEBOUNDARYMODES
We derive the characteristic equation (3.7) from the tight-
FIG.10:Thesetwographeneflakes(orquantumdots)bothhavethe
samezigzag-typeboundarycondition: Ψ = ±τz ⊗σzΨ. Thesign binding equation (3.1) and the definitions of the boundary
switchesbetween+and−whenthetangenttotheboundaryhasan modes(3.4) and(3.5). In thelow energylimitε/t a/T
≪ | |
anglewiththex-axiswhichisamultipleof60◦. wemaysetε 0inEq.(3.1),soitsplitsintotwodecoupled
→
sets of equationsfor the wave function on sublattices A and
B:
ψ (r)+ψ (r R )+ψ (r R )=0, (B1a)
B B 1 B 2
Acknowledgments − −
ψ (r)+ψ (r+R )+ψ (r+R )=0. (B1b)
A A 1 A 2
ThisresearchwassupportedbytheDutchScienceFounda- SubstitutingR byR +R intheseequationsandusingthe
1 2 3
tionNWO/FOM.WeacknowledgehelpfuldiscussionswithI. definition(3.5)ofλweexpressψ(r+R )throughψ(r),
2
Adagideli,J.H.Bardarson,Ya. B.Bazaliy,andI.Snyman.
ψ (r+R )= (1+λ)−1ψ (r), (B2a)
B 2 B
−
ψ (r+R )= (1+λ)ψ (r). (B2b)
A 2 A
−
APPENDIXA:DERIVATIONOFTHEGENERAL
BOUNDARYCONDITION(2.8) Eqs.(3.5)and(B2)togetherallowtofindtheboundarymode
withagivenvalueofλonthewholelattice:
We first show that the anticommutation relation (2.7) fol-
ψ (r+pR +qR )=λq( 1 λ)−pψ (r), (B3a)
lows from the current conservation requirement (2.6). The B 2 3 − − B
ψ (r+pR +qR )=λq( 1 λ)pψ (r), (B3b)
current operator in the basis of eigenvectors of M has the A 2 3 A
− −
blockform
withpandqarbitraryintegers.Substitutingψ(r+T)intoEq.
(3.4) from Eq. (B3) and using T = (n+m)R +nR we
X Y 1 0 2 3
nB ·J =(cid:18)Y† Z(cid:19), M =(cid:18)0 1(cid:19). (A1) arriveatthecharacteristicequation(3.7).
− We now find the rootsof the Eq. (3.7) for a givenk. It is
The HermitiansubblockX acts in the two-dimensionalsub- sufficientto analyze the equationfor sublattice A only since
spaceofeigenvectorsofM witheigenvalue1. Toensurethat the calculationfor sublattice B is the same after interchang-
Ψn J Ψ = 0foranyΨinthissubspaceitisnecessary ing n and m. The analysis of Eq. (3.7a) simplifies in polar
B
hand| suffi·cie|ntithat X = 0. The identity (n J)2 = 1 is coordinates,
B
equivalenttoYY† =1andZ =0,hence M,n· J =0.
{ B· } 1+λm+n = λn (B4)
We now show that the most general 4 4 matrix M that | | | |
× (m+n)arg( 1 λ) k narg(λ)=2πl, (B5)
satisfies Eqs. (2.5) and (2.7) has the 4-parameterform (2.8).
− − − −
10
with l = 0, 1, 2.... The curve defined by Eq. (B4) is tween λ and λ∗, the incrementof the left-hand side of Eq.
a contour on±the±complex plane around the point λ = 1 (B5)betnweenλ∗nandλ mustbeequalto2π(n 1) 2πn
− n n − ≈
which crosses points λ = 1/2 i√3/2 (see Fig. 11). (for T a),whichimmediatelyleadstoEq.(3.14)forλ .
± n
− ± | |≫
Theleft-handsideofEq.(B5)isamonotonicfunctionofthe
positiononthiscontour.Ifitincreasesby2π∆lontheinterval
Im(λ)
betweentworootsoftheequation,thenthereare∆l 1roots
−
inside this interval. For k = 0 both λ and λ are rootsof
− +
thecharacteristicequation. Sointhiscasethenumber of
A
N
roots lying inside the unit circle can be calculated from the
incrementof the left-hand side of Eq. (B5) between λ and
−
λ+: Re(λ)
1 0
−
1 2π 2π n m
= (n+m) +n 1=n − 1.
A
N 2π (cid:20) 3 3 (cid:21)− − 3 −
(B6)
Similarly,onsublatticeB,wehave(uponinterchangingnand
m),
FIG.11:Plotofthesolutionsofthecharacteristicequations(B4,B5)
m n
=m − 1. (B7) forn = 5, m = 11, andk = 0. Thedotsaretheroots, thesolid
B
N − 3 − curveisthecontourdescribedbyEq.(B4),andthedashedcirclesare
unitcircleswithcentersat0and−1.
The same method can be applied to calculate λ . Since
n
therearen 1 rootsonthecontourdefinedbyEq.(B4)be-
−
[1] P.R.Wallace,Phys.Rev.71,622(1947). arXiv:cond-mat/0701599.
[2] D.P.DiVincenzoandE.J.Mele,Phys.Rev.B29,1685(1984). [10] M.Y.Han,B.Oezyilmaz, Y.Zhang,andPh.Kim,Phys.Rev.
[3] Y.-W.Son,M.L.Cohen,andS.G.Louie,Phys.Rev.Lett.97, Lett.98,206805(2007).
216803(2006). [11] K.Nakada,M.Fujita,G.Dresselhaus,andM.S.Dresselhaus,
[4] E.McCannandV.I.Fal’ko,J.Phys.Condens.Matter16,2371 Phys.Rev.B54,17954(1996).
(2004). [12] ThemethoddescribedinSec.IIIcanbegeneralizedtobound-
[5] A.R.AkhmerovandC.W.J.Beenakker, Phys.Rev.Lett.98, arieswithN′ >n+msuchasthe“stronglydisorderedzigzag
157003(2007). boundary” of I. Martin and Ya. M. Blanter, arXiv:0705.0532.
[6] L.BreyandH.A.Fertig,Phys.Rev.B73,235411(2006). Forthesenon-minimalboundariesthezigzagboundarycondi-
[7] M.Ezawa,Phys.Rev.B73,045432(2006);PhysicaStatusSo- tionisstillgeneric.
lidi(c)4,489(2007). [13] A.RycerzandC.W.J.Beenakker,arXiv:0709.3397.
[8] M.V.BerryandR.J.Mondragon,Proc.R.Soc.LondonA412, [14] A.Rycerz,arXiv:0710.2859.
53(1987).
[9] Z. Chen, Y.-M. Lin, M. J. Rooks, and Ph. Avouris,