Spectrum of π electrons in bilayer graphene nanoribbons and nanotubes: an analytical approach J. Ruseckas and G. Juzeliu¯nas ∗ Institute of Theoretical Physics and Astronomy, Vilnius University A. Goˇstauto 12, LT-01108 Vilnius, Lithuania I. V. Zozoulenko Solid State Electronics, ITN, Link¨oping University, 601 74 Nork¨oping,Sweden 1 We present an analytical description of π electrons of a finite size bilayer graphene within a 1 framework of the tight-binding model. The bilayered structures considered here are characterized 0 by a rectangular geometry and have a finite size in one or both directions with armchair- and 2 zigzag-shaped edges. We provide an exact analytical description of the spectrum of π electrons in n the zigzag and armchair bilayer graphene nanoribbons and nanotubes. We analyze the dispersion a relations, thedensity of states, and theconductance quantization. J 8 PACSnumbers: 73.22.Pr 1 ] l l a h - s e m . t a m - d n o c [ 2 v 3 7 6 1 . 0 1 0 1 : v i X r a 2 I. INTRODUCTION Sinceitsisolationin2004,graphene—asinglesheetofcarbonatomsarrangedinahoneycomblattice—hasattracted an enormous attention because of its highly unusual electronic and transport properties that are strikingly different from those of conventional semiconductor-based two-dimensional electronic systems (for a review see Refs. 1–4). It has been immediately realizedthe significance andthe potential impactof this new materialfor electronics. This far, it has been demonstrated that the graphene has the highest carrier mobility at room temperature in comparison to any known material5. However, graphene is a semimetal with no gap and zero density of states at the Fermi energy. This makes it difficult to utilize it in electronic devices such as field effect transistor (FET) requiring a large on/off current ratio. The energy gap can be opened in a bilayer graphene by applying a gate voltage between the layers6. Thisgate-inducedbandgapwasdemonstratedbyOostingaet al.7, andthe on/offcurrentratioofaround100atroom temperature for a dual-gate bilayer graphene FET was reported by the IBM8. Another way to introduce the gap is to pattern graphene into nanoribbons9,10. The conductance of graphene nanoribbons (GNRs) with lithographically etched edges indeed revealed the gap in the transport measurements11,12. This gap has been subsequently understood as the edge-disorder-inducedtransportgap13–15 rather than the intrinsic energy gap expected in ideal GNRs due to the confinement9 or electron interactions and edge effects10 During last years the great progress has been achieved in fabrication and patterning of the GNRs with ultrasmooth and/or atomically controlled edges. This includes e.g. a controlled formation of edges by Joule heating16, unzipping carbon nanotubes to form nanoribbons17, chemical route to produce nanoribbons with ultrasmooth edges18 and atomically precise bottom-up fabrication of GNRs19. All these advances in nanoribbons fabrication will hopefully enable not before long the electronic measurement in near-perfect nanoribbons free from the edge or bulk disorder defects. An important insight into electronic properties of graphene and GNRs can be obtained from exact analytical approaches. The analytic calculations for the electronic structure of the GNRs have been reported in Refs. 20– 24. The electronic structure of the bilayer graphene was addressed in Refs. 25–29 where the analytical results were presented (both exact and perturbative). We are not however aware of analytical treatment of bilayer GNRs (Note that a numericalstudy of the magnetobandstructure of the GNRs was reported in Ref. 30 and 31, and the analytical and numericaltreatment of the edge states in the bi- and N-layer graphene and GNRs was presented in Refs. 29 and 32). The purpose of the present study is to provide an exact analytical description of the spectrum of π electrons in the zigzag and armchair bilayer nanoribbons and nanotubes including the dispersion relations, the density of states, and the conductance quantization. The paper is organized as follows: In order to illustrate our method, in Sec. II we present known analytical results for a simpler system, monolayer graphene of the finite size. Subsequently in Sec. III we derive the main analytical expressions for the energy spectrum of finite-size structures of bilayer graphene. These expressions are used in Sec. IV to analyze the energy spectrum of various bilayer graphene structures near the Fermi energy. Finally, Sec. V summarizes our findings. II. SINGLE LAYER GRAPHENE Analytical expressions for the π electron spectrum in GNRs and graphene nanotubes (GNTs), based on tight- binding model, were provided in Ref. 23. In this section we will rederive the same expressions in an analytically simpler way. Our method more clearly shows the connection between solutions for the infinite sheet of graphene and for the finite-size sheet. In addition, simpler method will allow us to derive later on analytical expressions of the π electron spectrum for more complex systems, bilayer GNRs and GNTs. A. Electron spectrum in infinite sheet of graphene First we will consider π electron spectrum in an infinite sheet of graphene. Hexagonal structure of graphene is shown in Fig 1a. The structure of the graphene can be viewed as a hexagonal lattice with a basis of two atoms per unitcell. TheCartesiancomponentsofthelatticevectorsa anda area(3/2,√3/2)anda(3/2, √3/2),respectively. 1 2 − Here a 1.42˚A is the carbon-carbon distance1. The three nearest-neighbor vectors are given by δ =a(1/2,√3/2), 1 ≈ δ =a(1/2, √3/2), and δ =a( 1,0). The tight-binding Hamiltonian for electrons in graphene has the form 2 1 − − Hgr =−t (a†ibj +b†jai), (1) i,j hXi 3 FIG. 1. (Color online) (a) Honeycomb lattice structure of graphene, made out of two interpenetrating triangular lattices. a1 and a2 are the lattice unit vectors, and δi, i = 1,2,3 are the nearest-neighbor vectors. (b) Indication of labels of carbon atomsin therectangularunitcell. (c)Brillouin zonesforhexagonalunitcell(solid hexagon) andrectangularunitcell(dashed rectangle). TheDiracpointsareindicatedbysolidcirclesforthehexagonalunitcellandhollowcirclesfortherectangularunit cell. where the operators a and b annihilate an electron on sublattice A at site RA and on sublattice B at site RB, i i i i respectively. The parameter t is the nearest-neighbor hopping energy (t 2.8eV). From now on we will write all ≈ energies in the units of the hoping integral t, therefore we will set t = 1. Let us label the elementary cells of the lattice with twonumbers p andq. Thenthe atomsinthe sublattices AandB arepositionedatRA =pa +qa and p,q 1 2 RB =δ +pa +qa , respectively. p,q 1 1 2 The π electron wave function satisfies the Schr¨odinger equation, HΨ=EΨ. (2) We search for the eigenvectors of the Hamiltonian (1) in the form of the plane waves (Bloch states) by taking the probability amplitudes to find an atom in the sites RA and RB of the sublattices A and B as p,q p,q ψpA,q =cAeik·RAp,q, ψpB,q =cBeik·RBp,q. (3) Thus Eq. (2) yields the eiganvalue equations for the coefficients cA and cB EcA =cBφ˜(k), (4) − EcB =cAφ˜( k), (5) − − where φ˜(k) eik·δ1 +eik·δ2 +eik·δ3. (6) ≡ From Eqs. (4) and (5) we get the eigenenergies and the corresponding coefficients determining the eigenvectors φ˜(k) E(k)=s φ˜(k) , cA = , cB =1, (7) 1| | −E(k) where s = 1. In the anticipation of the rectangular geometry we introduce dimensionless Cartesian components of 1 ± the wave vector κ=3ak , ξ =√3ak (8) x y instead of the wave vector components k and k . Then using the coordinates of the vectors δ we have x y j ξ φ˜(k)=e−iκ3 +2eiκ6 cos (9) 2 (cid:18) (cid:19) and the expression for the eigenenergies becomes1 ξ ξ κ E(k)=s 1+4cos2 +4cos cos . (10) 1 s 2 2 2 (cid:18) (cid:19) (cid:18) (cid:19) (cid:16) (cid:17) 4 For satisfying boundary conditions it is useful to adopt a larger unit cell characterized the same geometry as the whole sheetofgraphene. Since we areinterestedinconfigurationsofthe graphenewith rectangulargeometry,we will use a rectangular unit cell, as it has been done in Ref. 23. Such unit cell has four atoms labeled with symbols l , λ , ρ , r , as it is shown in Fig. 1b. The atoms with labels l and ρ belong to the sublattice A , the atoms with labels λ and r belong to the sublattice B. The position of the unit cell is indicated with two numbers n and m. The first Brillouin zone corresponding to the rectangular unit cell contains the values of the wave vectors κ, ξ in the intervals π κ<π , π ξ <π. We search for the eigenvectors having the form of plane waves, − ≤ − ≤ ψ =c eiξm+iκn, (11) m,n,α α where α=l,ρ,λ,r. This solution can be obtained from Eq. (3) using the equalities cr =cB, cρ =cAe−ik·δ1, cλ =cBe−ik·a1, cl =cAe−i2akx. (12) The Brillouin zones corresponding to hexagonal and rectangular unit cells are shown in Fig. 1c. Compared to the area of the Brillouin zone of the hexagonal unit cell, the area of the Brillouin zone of the rectangular unit cell is two times smaller. Smaller Brillouin zone leads to the appearance of additional dispersion branches. Those dispersion branches can be taken into account by using two values of the wave vector κ in Eqs. (10) and (7), the one with π κ < π and another obtained replacing κ by 2π +κ. Using Eqs. (7), (12) we obtain the coefficients of the − ≤ eigenvectors cr =1, cρ = e−iξ2 φ(κ,ξ) , (13) − E(κ,ξ) φ(κ,ξ) cl = s3e−iκ2 , cλ =s3e−i21(κ+ξ), (14) − E(κ,ξ) where ξ φ(κ,ξ)=s3e−iκ2 +2cos (15) 2 (cid:18) (cid:19) and s = 1 indicates the dispersion branches that appear due to the smaller Brillouin zone. The equation for the 3 ± energy now becomes ξ ξ κ E(κ,ξ)=s 1+4cos2 +s 4cos cos . (16) 1 3 s 2 2 2 (cid:18) (cid:19) (cid:18) (cid:19) (cid:16) (cid:17) This equation has been obtained in23. Zero energy points of the graphene honeycomb lattice with dispersionrelation (10) are at the points K = (2π,2π/3) and K = (2π, 2π/3), where coordinates are given in (κ,ξ) space. K points ′ − correspond to the corners of the first Brillouin zone. Using the Brillouin zone corresponding to the rectangular unit cell, the zero energy points have coordinates 0, 2π and the number of these points is only two, as it is shown in ± 3 Fig. 1c. (cid:0) (cid:1) Since we will consider finite-size graphene sheets, evanescent solutions become important. Solution exponentially decreasing or increasing in the x-direction can be obtained by taking κ = iκ in Eqs. (13), (14) and (16), whereas | | solutionexponentially decreasingorincreasinginthe y-directioncanbe obtainedby takingξ =iξ . The dependency | | oftheenergyonκwhenξ =0isshowninFig.2. Weseethatthebrancheswithrealandimaginaryκdonotintersect at κ >0. | | B. Electron spectrum in various single layer graphene structures From the boundary conditions we get restrictions on the possible values of the wave vectors κ, ξ. We will consider thestructuresofgraphenethathaveasetofN rectangularunitcellsinthex(armchair)directionandasetof +1/2 N rectangularunitcellsinthey (zigzag)direction,sothatthereare hexagonsalongthey axis. Notethatrectangular N unit cell shown in Fig. 1b extends over the whole hexagon in the y direction, whereas it extends over more that one hexagon in x direction. Using periodic boundary condition, corresponding to the graphene torus, we get that the possible values of the wave vectors κ, ξ are 2π 1 ξ = j, j = N , N +1,..., N − (17) j − 2 − 2 2 N (cid:22) (cid:23) (cid:22) (cid:23) (cid:22) (cid:23) 2π N N N 1 κ = ν, ν = , +1,..., − (18) ν N − 2 − 2 2 (cid:22) (cid:23) (cid:22) (cid:23) (cid:22) (cid:23) 5 4 2 E 0 -2 -4 0 0.5 1 1.5 2 2.5 3 |κ| FIG.2. (Color online) Left: dispersion branchesofgraphenefor rectangular unitcell, calculated according toEq.(16). Right: dispersion branchesfor ξ=0, showing propagating solutions (red solid) and evanescent solutions (green dashed). Here denotes the integer part of a number. Thus the spectrum of graphene torus is given by Eq. (16) replacing κ ⌊·⌋ and ξ by κ and ξ . ν j For graphene armchair nanotubes one has the periodic boundary condition in the x direction and the requirement ψ = ψ = ψ = ψ = 0 for the y direction. Since the energy (16), does not depend on the sign 0,n,r 0,n,l +1,n,l +1,n,r of wave vector ξ, wNe will searchNfor the eigenvectors of the Hamiltonian (1) as a superposition of periodic solutions Eq. (11) with ξ and ξ, − ψm,n,α =acα(ξ,κν)eiξm+iκνn+bcα( ξ,κν)e−iξm+iκνn, (19) − where κ is given by Eq. (18) and ξ needs to be determined. From the boundary conditions we get a system of two ν equations for the coefficients a and b ac (ξ,κ )+bc ( ξ,κ )=0, (20) r,l ν r,l ν − aeiξ( +1)c (ξ,κ )+be iξ( +1)c ( ξ,κ )=0. (21) N r,l ν − N r,l ν − This system of equations has non-zero solutions only when the determinant is zero. From Eqs. (13), (14) it follows that the coefficients c (ξ,κ) do not depend on the sign of ξ and we get the condition sin(ξ( +1))=0 or r,l N πj ξ = , j =1,..., (22) +1 N N Additionally, there are two N-fold degenerate levels corresponding to ξ = π with energies E = 1. The states of ± those levels have zero wave function amplitudes at the l and r sites. For graphene zigzag nanotubes one has the periodic boundary condition in the y direction and the condition ψ =ψ =0 for the x direction. Similarly as for the armchair nanotubes, the energy (16), does not depend m,0,r m,N+1,l on the signof wave vector κ, and we searchfor the eigenvectorsof the Hamiltonian (1) as a superposition of periodic solutions Eq. (11) with κ and κ, − ψm,n,α =acα(ξj,κ)eiξjm+iκn+bcα(ξj, κ)eiξjm−iκn, (23) − where ξ is given by Eq. (17) and κ needs to be determined. From the boundary conditions we get a system of two j equations for the coefficients a and b ac (ξ ,κ)+bc (ξ , κ)=0, (24) r j r j − ac (ξ ,κ)eiκ(N+1)+bc (ξ , κ)e iκ(N+1) =0. (25) l j l j − − Using Eqs. (13), (14) we obtain that non-zero solutions are possible when sin(κN) ξ j = s 2cos . (26) sin κ N + 1 − 3 2 2 (cid:18) (cid:19) (cid:0) (cid:0) (cid:1)(cid:1) 6 The possible values of wave vector κ should obey this equation. The same condition has been obtained in Ref. 23. Equation(26) allowsfor the imaginaryvalues of wavevectorκ. The imaginaryvalues appear when ξc < ξ <π and j s = 1, where the critical value ξc = 2arccos N/(2N +1) of the wave vector ξ is obtained from Eq.|(26|) setting 3 − κ = 0. In the limit N from the condition (26) with imaginary κ and Eq. (16) follows that E = 0: edge states → ∞ (cid:0) (cid:1) near zigzag edges in the semi-infinite system have zero energy. For N sheet of graphene open boundary conditions in the y direction are the same as for armchair nanotubes ×N and in the x direction are the same as for zigzag nanotubes. Since the resulting conditions for the wave vectors κ, ξ arenotcoupled,theeigenvectoroftheHamiltonian(1)isasuperpositionoffourperiodicsolutionshavingallpossible combinations of the signs of κ and ξ and the possible values of the wave vectors are given by Eqs. (22) and (26). In addition there are two N-fold degenerate levels corresponding to ξ =π with energies E = 1. ± III. BILAYER GRAPHENE Now we will consider the spectrum of π electrons in bilayer graphene. The tight-binding Hamiltonian for electrons in bilayer graphene has the form Hbi =V (a†j,2aj,2+b†j,2bj,2−a†j,1aj,1−b†j,1bj,1) j X −t (a†i,pbj,p+b†j,pai,p)−t⊥ (a†j,1aj,2+a†j,2aj,1), (27) i,j ,p j hXi X where the operators a and b annihilate an electron on sublattice A at site RAp and on sublattice B at site i,p i,p p i p RBp, respectively. The index p=1,2numbers the layersin the bilayersystem. Inthe Hamiltonian(27) weneglected i the terms corresponding to the hopping between atom B and atom B , with the hopping energy γ , and the terms 1 2 3 corresponding to the hopping between atom A (A ) and and atom B (B ) with the hopping energy γ . Neglect 1 2 2 1 4 of those hopping terms leads to the minimal model of bilayer graphene28. The parameter t (t 0.4eV) is the hopping energy between atom A and atom A while V is half the shift in the electrochemical⊥pot⊥en≈tial between the 1 2 two layers. Similarly as for the monolayer graphene, we will express all the energies in the units of t. A. Electron spectrum in infinite sheet of bilayer graphene We will proceed similarly as in the previous Section and will analyze an infinite system at first. The atoms in the sublattices A and A are positioned at RA1,2 = pa +qa , in the sublattice B the atoms are positioned at 1 2 p,q 1 2 1 RpB,1q = δ1+pa1+qa2 and in the sublattice B2 the atoms are positioned at RBp,2q = −δ1+pa1+qa2. We search for the eigenvectors of the Hamiltonian (27) in the form of the plane waves. The probability amplitudes to find an atom in the sites RA1,2 and RB1,2 of the sublattices A and B are p,q p,q j j ψpA,1q,2 =cA1,2eik·RAp,1q,2 , ψpB,1q,2 =cB1,2eik·RBp,1q,2 . (28) The coefficients cAp and cBp obey the eigenvalue equations EcA1 =VcA1 +cB1φ˜(k)+γcA2, (29) − EcB1 =VcB1 +cA1φ˜( k), (30) − − EcA2 = VcA2 +cB2φ˜( k)+γcA1, (31) − − − EcB2 = VcB2 +cA2φ˜(k). (32) − − Here energy E, potential V and interaction between layers γ t /t are in the units of the hoping integral t. Using the nearest-neighbor hopping energy t 2.8eV and the hop≡pin⊥g energy between two layers t 0.4eV one gets γ 0.14. When V =0, the π electron sp≈ectrum is determined by the equation ⊥ ≈ ≈ γ γ2 E(k)=s s + + φ˜(k)2 , (33) 1 2 2 r 4 | | ! 7 where s ,s = 1. The coefficients of the eigenvector are 1 2 ± E(k) cA1 = , cB1 =1, (34) −φ˜( k) − E(k) φ˜(k) cA2 =s s , cB2 = s s . 1 2φ˜( k) − 1 2φ˜( k) − − When V =0 the spectrum is 6 γ2 γ4 E(k)=s +V2+ φ˜(k)2+s + φ˜(k)2(4V2+γ2) (35) 1 2 s 2 | | 4 | | r and the coefficients of the eigenvector are E(k)+V cA1 = , cB1 =1, (36) − φ˜( k) − E(k) V φ˜(k) cA2 = − f(k), cB2 = f(k), φ˜( k) −φ˜( k) − − where the function (E(k)+V)2 φ˜(k)2 f(k)= −| | (37) γ(E(k) V) − describes the contribution of the second sheet of graphene to the eigenvector. Finite-size bilayer graphene sheets can be in AB-α or AB-β stacking, as is shown in Fig. 3a,b. Similarly as for graphene monolayer, we will use rectangular unit cells, one shifted with respect to the other, in each layer of bilayer graphene. However, the position of rectangular cells are different for different stacking types. Rectangular unit cells have eight atoms with labels l , λ , ρ , r and l , λ , ρ , r , as is shown in Fig. 3c,d. For the AB-α stacking the 1 1 1 1 2 2 2 2 atoms with labels l , ρ belong to the sublattice A , atoms λ , r to the sublattice B , atoms l , ρ to the sublattice 1 1 1 1 1 1 2 2 A andatomsλ , r to the sublattice B . Forthe AB-β stackingthe atoms withlabels l , ρ belongto the sublattice 2 2 2 2 1 1 B , atoms λ , r to the sublattice A , atoms l , ρ to the sublattice B and atoms λ , r to the sublattice A . 1 1 1 1 2 2 2 2 2 2 We search for the solutions of the form ψ =c eiξm+iκn (38) m,n,αp αp where α=l,ρ,λ,r is the label of atoms and p=1,2 is the number of the layer. For the AB-α stacking this solution can be obtained from Eq. (28) using the equalities cr1 =cB1, cρ1 =cA1e−ik·δ1, cλ1 =cB1e−ik·a1, cl1 =cA1e−i2akx, (39) cr2 =cB2e−iakx, cρ2 =cA2e−ik·δ1, cλ2 =cB2eik·δ2, cl2 =cA2eiakx (40) whereas for the AB-β stacking the coefficients are cr1 =(cA1)∗, cρ1 =(cB1)∗e−ik·δ1, cλ1 =(cA1)∗e−ik·a1, cl1 =(cB1)∗e−i2akx, (41) cr2 =(cA2)∗e−ik·a2, cρ2 =(cB2)∗e−iakx, cλ2 =(cA2)∗, cl2 =(cB2)∗eik·δ1. (42) Similarlyasformonolayergraphene,totakeintoaccountthesmallerBrillouinzoneweneedtwodispersionbranches: one with κ and one with 2π +κ . Using Eq. (34) or Eq. (36) we obtain the coefficients of the eigenvectors. The expressions for the coefficients are presented in Appendix A. The expression for the energy becomes γ2 γ4 E(κ,ξ)=s +V2+ φ(κ,ξ)2+s + φ(κ,ξ)2(4V2+γ2) (43) 1 2 s 2 | | 4 | | r which reduces to γ γ2 E(κ,ξ)=s s + + φ(κ,ξ)2 (44) 1 2 2 r 4 | | ! 8 FIG. 3. (Color online) Upper part: sublatices A1 , A2 , B1 , B2 on bilayer graphene in AB-α stacking (a) and AB-β stacking (b). Lowerpart: indicationoflabelsofcarbonatomsusedinthedescriptionoftheπelectronspectrumforthebilayergraphene with AB-αstacking (c) and AB-β stacking (d). for V =0. Here ξ ξ κ φ(κ,ξ)2 =1+4cos2 +s 4cos cos (45) 3 | | 2 2 2 (cid:18) (cid:19) (cid:18) (cid:19) (cid:16) (cid:17) and s = 1 indicates the dispersion branches that appear due to the smaller Brillouin zone. 3 ± Inadditiontothepropagatingwaves,forfinite-sizebilayergraphenesheetsevanescentsolutionsbecomeimportant. Solution exponentially decreasing or increasing in the x-direction can be obtained by taking κ = iκ. Solution | | exponentially decreasingorincreasingin the y-directioncanbe obtainedby takingξ =iξ . In additionto the purely | | imaginaryξ therearesolutions,correspondingtos = 1,havingcomplexvaluesofξ. Thedependency oftheenergy 3 − onthe wavevectorκ whenthe wavevectorξ isconstantandonthe wavevectorξ whenthe wavevectorκis constant is shown in Fig. 4. We see that now, in contrast to the graphene monolayer, the branches with real and imaginary κ can have the same energy. B. Electron spectrum in various bilayer graphene structures We will consider the structures of bilayer graphene that have a set of N rectangular unit cells in the x (armchair) direction and a set of +1/2 rectangular unit cells in the y (zigzag) direction, so that there are hexagons along N N the y axis. Note that the rectangular unit cells shown in Figs. 3c and 3d extend over the whole hexagon in the y direction,whereastheyextendovermorethatonehexagoninxdirection. Inprinciple,inthecaseofbilayergraphene nanotubes the numbers N or in for the inner and outer cylinders are different. However, for simplicity we will N consider them as the same, which is a good approximation for sufficiently large tubes when N or . →∞ N →∞ Similarly as forgraphene monolayer,fromthe boundaryconditions we get restrictionsonthe possible values of the wave vectors κ, ξ. Using periodic boundary condition, corresponding to the bilayer graphene torus, we get that the 9 4 4 2 2 E 0 E 0 -2 -2 -4 -4 0 0.5 1 1.5 2 2.5 3 0 0.5 1 1.5 2 2.5 3 |κ| Re ξ+Im ξ FIG.4. (Coloronline)Dispersion branchesofbilayergraphene: dependencyoftheenergyonthewavevectorκwhenthewave vector ξ is constant (ξ=0) (left) and on thewave vector ξ when thewave vector κ is constant (κ=1.0) (right). Propagating solutions are shown with red solid line, evanescent solutions with green dashed line, and evanescent oscillating solutions with complex valueof thewave vectorξ are shown with blue dotted line. Inorder to show thestructureof thedispersion branches more clearly, thevalueof theparameter γ is set sufficiently large, γ =0.5. possible values of the wave vectors κ, ξ are given by Eqs. (17), (18). For bilayer graphene armchair nanotubes one has the periodic boundary condition in the x direction and the condition ψ =ψ =ψ =ψ =0 (46) 0,n,rp 0,n,lp N+1,n,lp N+1,n,rp for the y direction. Here p=1,2 is the number of the layer. This condition is the same for both the AB-α and AB-β stackings. For bilayergraphene with AB-α stackingthe coefficients c (ξ,κ) do not depend onthe signof ξ and we rp,lp get the same conditions (18), (22) for the wave vectors κ, ξ, as for the monolayer graphene armchair tubes. For bilayer graphene with AB-β stacking the coefficients c (ξ,κ) depend on the sign of ξ, and condition for the rp,lp possiblevaluesofthewavevectorξ ismuchmorecomplicated. Thereareeightboundaryconditionsinthey direction. In bilayer graphene there are four eigenstates with different wave vectors along y direction, ξ(1), ξ(2), ξ(3) and ξ(4), having the same energy: E(κ,ξ(1)) = E(κ,ξ(2)) = E(κ,ξ(3)) = E(κ,ξ(4)), as is evident from Fig. 4. Two or four of the wave vectors ξ(1), ξ(2), ξ(3), ξ(4) can be imaginary or complex numbers. Since the energy does not depend on the sign of ξ, we can form a wave function from superposition of eight waves. From the boundary conditions (46) resulting resulting setof linearequations canhavenonzerosolutiononly if 8 8 determinantis zero. Analytical form × of this condition in is too large and too complicated to be useful. Forbilayergraphenezigzagnanotubesonehastheperiodicboundaryconditioninthey directionandthecondition ψ =ψ =ψ =ψ =0 (47) m,0,r1 m,N+1,l1 m,0,l2 m,N+1,r2 for the x direction. Here p=1,2 is the number of the layer. This condition is the same for both the AB-α and AB-β stackings. Inthe bilayergraphenethere aretwoeigenstateswithwavevectorsalongx direction,κ(1) andκ(2),having differentabsolute values but correspondingthe same energy: E(κ(1),ξ)=E(κ(2),ξ). One orboth ofthe wavevectors κ(1) , κ(2) can be imaginary. The energy can be equal only if the signs s , s obey the condition 1 2 s(2)s(2) = s(1)s(1) (48) 1 2 − 1 2 When the bias potential is zero, V =0, from the equality of the energy we can express κ(2): κ(2) κ(1) γ s(2)cos =s(1)cos +s(1)s(1) E(κ(1),ξ) (49) 3 (cid:18) 2 (cid:19) 3 (cid:18) 2 (cid:19) 1 2 2cos ξ 2 (cid:16) (cid:17) When V =0 then 6 κ(2) κ(1) γ s(2)cos =s(1)cos 4E2V2+γ2(E2 V2) (50) 3 (cid:18) 2 (cid:19) 3 (cid:18) 2 (cid:19)± 2cos ξ − 2 p (cid:16) (cid:17) There are four boundary conditions in the x direction. Since the energy does not depend on the sign of κ, we can form a wave function from superposition of four waves. From the boundary conditions (47) resulting set of linear 10 equations can have nonzero solution only if 4 4 determinant is zero. The possible values of the wave vector ξ is × given by Eq. (17), and the conditions for the possible values of the wave vector κ are given in the Appendix B. For N sheet of bilayer graphene open boundary conditions in the y direction are the same as for armchair ×N nanotubes, Eq. (46) and in the x direction are the same as for zigzag nanotubes, Eq. (47). For AB-α stacking, the conditions for the possible values of the wave vectors κ,ξ are combination of the conditions for zigzag and armchair bilayergraphenetubes. Specifically,whenV =0,theconditionsaregivenbyEqs.(22)and(B1)or(B2). WhenV =0 6 thenthe conditionsaregivenby Eqs.(22) and(B3). For AB-β stackingitis impossible toseparateconditionsfor the wave vector ξ from the conditions for the wave vector κ. The resulting expressions are very large and complicated. C. Summary of the possible values of wave vectors For structures of bilayer graphene, the energy spectrum is completely determined by Eq. (44) or (43) with appro- priate expressions for wave vectors κ and ξ. Equations presentedin Appendix B make one quantum number dependent onthe other. This dependence appears because of zigzag-shaped edges. For structures where zigzag edges do not exist or their effect can be disregarded the wave vector κ can be replaced by a continuous variable. ν Thus, the possible values of wave vectors for various structures are as follows: ForthearmchairbilayergrapheneribbonofinfinitelengthwithAB-αstacking,thewavevectorsaredetermined • by πj 0 κ π, ξ = , j =1,..., (51) j ≤ ≤ +1 N N For the armchair bilayer graphene ribbon of infinite length with AB-β staking we have 0 κ π and the • ≤ ≤ equation for the possible values of ξ is complicated. For the zigzag bilayer graphene ribbon of infinite length we have 0 ξ π, the conditions for the possible • ≤ ≤ values of κ, given in Appendix B, are different for AB-α and AB-β stackings. ForthezigzagbilayercarbontubeofinfinitelengthwithAB-αorAB-βstacking,thewavevectorsaredetermined • by 2π 1 0 κ π, ξ = j, j = N , N +1,..., N − (52) j ≤ ≤ − 2 − 2 2 N (cid:22) (cid:23) (cid:22) (cid:23) (cid:22) (cid:23) For the armchair bilayer carbon tube of infinite length with AB-α or AB-β stacking: • 2π N N N 1 κ = ν, ν = , +1,..., − , 0 ξ π. (53) ν N − 2 − 2 2 ≤ ≤ (cid:22) (cid:23) (cid:22) (cid:23) (cid:22) (cid:23) Takingintoaccounttherangesofthepossiblevaluesofthewavevectors,zeroenergypointsforvariousstructures with bias potential V =0 are as follows: For zigzag bilayer carbon tube zero energy points are 0,2π , 0, 2π . • 3 − 3 For armchair bilayer carbon tube, zero energy point is(cid:0) 0,2π(cid:1) .(cid:0) (cid:1) • 3 The dispersion of armchair bilayer graphene ribbon has(cid:0)only(cid:1)one zero-energy point 0,2π • 3 For zigzag bilayer graphene ribbon dispersion this point cannot be shown in the rea(cid:0)l plan(cid:1)e. • IV. BAND STRUCTURE NEAR THE FERMI ENERGY In this Section only a part of the spectrum with smallest absolute value of the energy is in focus. This part corresponds to s = 1, s = 1. In order to obtain an approximate expression for the energy spectrum near the 2 3 − −