Kronig-Penney model on bilayer graphene: spectrum and transmission periodic in the strength of the barriers
KKronig-Penney model on bilayer graphene: spectrum and transmission periodic in thestrength of the barriers
M. Barbier, P. Vasilopoulos, and F. M. Peeters Department of Physics, University of Antwerp, Groenenborgerlaan 171, B-2020 Antwerpen, Belgium Department of Physics, Concordia University, 7141 Sherbrooke Ouest, Montréal, Quebec, Canada H4B 1R6
We show that the transmission through single and double δ -function potential barriers of strength P = V W b / (cid:126) v F in bilayer graphene is periodic in P with period π . For a certain range of P valueswe find states that are bound to the potential barrier and that run along the potential barrier.Similar periodic behaviour is found for the conductance. The spectrum of a periodic succession of δ -function barriers (Kronig-Penney model) in bilayer graphene is periodic in P with period π . For P smaller than a critical value P c , the spectrum exhibits two Dirac points while for P larger than P c an energy gap opens. These results are extended to the case of a superlattice of δ -function barrierswith P alternating in sign between successive barriers; the corresponding spectrum is periodic in Pwith period π . PACS numbers: 71.10.Pm, 73.21.-b, 72.80.Vp
I. INTRODUCTION
Graphene, a one-atom thick layer of carbon atoms, hasbecome a research attraction pole since its experimentaldiscovery in 2004 . Since carriers in graphene behave likerelativistic and chiral massless fermions with a linear-in-wave vector spectrum, many interesting features could betested with this material such as the Klein paradox ,which was recently observed , the anomalous quantumHall effect, etc., see Ref. 5 for two recent reviews. Theeffort to realise this Klein tunnelling through a potentialbarrier also lead to other interesting features, such asresonant tunnelling through double barriers . With thepossibility to fabricate devices with single-layer graphene,bilayer graphene has also been extensively investigatedand been shown to possess extraordinary electronic be-haviour, such as a gapless spectrum, in the absence ofbias, and chiral carriers . Many of these nanostructurescould be given another functionality if based on bilayerinstead of single-layer graphene.The electronic band structure can be modified by theapplication of a periodic potential and/or magnetic bar-riers. Such superlattices (SLs) are commonly used toalter the band structure of nanomaterials. In single-layer graphene already a number of papers relate theirwork to the theoretical understanding of such periodicstructures . Much less experimental and theoreticalwork has been done on bilayer graphene .We will study the spectrum, the transmission, and theconductance of bilayer graphene through an array of po-tential barriers using a simple model: the Kronig-Penney(KP) model , i.e. a one-dimensional periodic successionof δ -function barriers on bilayer graphene. The advantageof such a model system is that, 1) a lot can be done ana-lytically, 2) the system is clearly defined, 3) and it is pos-sible to show a number of exact relations. The present re-search is also motivated by our recent findings for single-layer graphene , where very interesting and unexpectedproperties were found, for instance, that the transmis- sion and energy spectrum are periodic in the strength ofthe δ -function barriers. Surprisingly, we find that for bi-layer graphene similar, but different, properties are foundas function of the strength of the δ -function potentialbarriers. Due to the different electronic spectra close tothe Dirac point, i.e., linear for graphene and quadraticfor bilayer graphene, we find very different transmissionprobabilities through a finite number of barriers and verydifferent energy spectrum, for a superlattice of δ -functionbarriers, between single-layer and bilayer graphene.The paper is organised as follows. In Sec. II we brieflypresent the formalism. In Sec. III we give results for thetransmission and conductance through a single δ -functionbarrier. We dedicate Sec. IV to bound states of a single δ -function barrier and Sec. V to those of two such barri-ers. In Sec. VI we present the spectrum for the KP modeland in Sec. VII that for an extended KP model by con-sidering two δ -function barriers with opposite strength inthe unit cell. Finally, in Sec. VIII we make a summaryand concluding remarks. II. BASIC FORMALISM
We describe the electronic structure of an infinitelylarge flat graphene bilayer by the continuum nearest-neighbour, tight-binding model and consider solutionswith energy and wave vector close to the K ( K (cid:48) ) point.The corresponding four-band Hamiltonian and eigen-states Ψ are H = V v F π t ⊥ v F π † V t ⊥ V v F π † v F π V , ψ = ψ A ψ B ψ B (cid:48) ψ A (cid:48) . (1)with π = p x + ip y ( p x,y = − i (cid:126) ∂ x,y ) and p the mo-mentum operator. We apply one-dimensional potentials V ( x, y ) = V ( x ) and consequently the wave function canbe written as ψ ( x, y ) = ψ ( x ) e ik y y with the momentum in a r X i v : . [ c ond - m a t . m e s - h a ll ] J a n the y-direction a constant of motion. Solving the time-independent Schrödinger equation H ψ = Eψ we obtain,for constant V ( x, y ) = V , the spectrum and the eigen-states. The latter are given by Eq. (A5) in App. A andthe spectrum by Eq. (A2) ε = u + 1 / ± (cid:112) / k ,ε = u − / ± (cid:112) / k , (2)where we used the dimensionless variables, ε = E/t ⊥ , u = V /t ⊥ , x → xt ⊥ / (cid:126) v F , k y → (cid:126) v F k y /t ⊥ and ε (cid:48) = ε − u ; v F = 10 m/s, and t ⊥ = 0 . eV expresses the couplingbetween the two layers.For later purposes we also give the frequently used two-band Hamiltonian H = − v F t ⊥ (cid:18) π † π (cid:19) + V, (3)and the corresponding spectrum E − V = ± ( v F (cid:126) /t ⊥ )( k x + k y ) . (4)As seen, there are qualitative differences between the twospectra (compare Eq. (4) with Eq. (2)) that will be re-flected in those for the transmission and conductance insome of the cases studied. As the approximation of thetwo-band Hamiltonian is only valid for E (cid:28) t ⊥ , we canexpect a qualitative difference with the four-band Hamil-tonian if | E | (cid:38) t ⊥ . III. TRANSMISSION THROUGH A δ -FUNCTION BARRIER We assume | E − V | < t ⊥ outside the barrier such thatwe obtain one pair of localised and one pair of travellingeigenstates in the well regions characterised by wave vec-tors α and β , where α is real and β imaginary, see App.A. Consider an incident wave with wave vector α fromthe left (normalised to unity); part of it will be reflected,with amplitude r , and part of it will be transmitted withamplitude t . Then the transmission is T = | t | . Also,there are growing and decaying evanescent states nearthe barrier, with coefficients e g and e d , respectively. Therelation between the coefficients can be written in theform N t e d = r e g . (5)This leads to a system of linear equations that can bewritten in matrix form = N N N − N N N N N − tre d e g , (6) with N ij the coefficients of the transfer matrix N . De-noting the matrix in Eq. (6) by Q , we can evaluate thecoefficients from (cid:0) t, r, e d , e g , (cid:1) T = Q − (cid:0) , , , (cid:1) T .As a result, to obtain the transmission ampli-tude t it is sufficient to find the matrix element ( Q − ) = [ N − N N /N ] − .We model a δ -function barrier as the limiting case ofa square barrier, with height V and width W b shown inFig. 1, represented by V ( x ) = V Θ( x )Θ( W b − x ) . The b FIG. 1: Schematics of the potential V(x) of a single squarebarrier. transfer matrix N for this δ -function barrier is calculatedin appendix B and the limits V → ∞ and W b → aretaken such that P = V W b / (cid:126) v F is kept constant.The transmission T = | t | for α real and β imaginaryis obtained from the inverse amplitude, t = cos P + iµ sin P + ( α − β ) k y αβε sin P cos P + iν sin P , (7)where µ = ( ε +1 / /α and ν = ( ε − / /β . Contour plotsof the transmission T are shown in Figs. 2(a) and 2(b)for strengths P = 0 . π and P = 0 . π , respectively.The transmission remains invariant under the trans-formations: P → P + nπ, k y → − k y . (8)The first property is in contrast with what is obtained inRef. 3. In the latter work it was found, by using the × Hamiltonian, that the transmission T should be zero for k y ≈ and E < V , while we can see here that for certainstrengths P = nπ there is perfect transmission. The lastproperty is due to the fact that k y only appears squaredin the expression for the transmission. Notice that incontrast to single-layer graphene the transmission for ε ≈ is practically zero. The cone for nonzero transmissionshifts to ε = 1 / − cos P ) with increasing P till P = π . An area with T = 0 appears when α is imaginary,i.e., for ε + ε − k y < (as no propagating states areavailable in this area, we expect bound states to appear).From Figs. 2(a), 2(b) it is apparent that the transmissionin the forward direction, i.e., for k y ≈ , is in generalsmaller than ; accordingly, there is no Klein tunneling.However, for P = nπ , with n an integer, the barrierbecomes perfectly transparent.For P = nπ we have V = (cid:126) v F ( nπ/W b ) . If the elec-tron wave vector is k = nπ/W b its energy equals theheight of the potential barrier and consequently there isa quasi-bound state and thus a resonance . The condi-tion on the wave vector implies W b = nλ/ where λ isthe wavelength. This is the standard condition for Fabry-Perot resonances. Notice though that the invariance ofthe transmission under the change P → P + nπ is notequivalent to the Fabry-Perot resonance condition. FIG. 2: (Color online) Contour plot of the transmission for P = 0 . π in (a) and P = 0 . π in (b). In (b) the boundstate, shown by the red curve, is at positive energy. The whitearea shows the part where α is imaginary. The probabilitydistribution | ψ ( x ) | of the bound state is plotted in (c) forvarious values of k y and in (d) for different values of P . From the transmission we can calculate the conduc-tance G given, at zero temperature, by G/G = (cid:90) π/ − π/ T ( E, φ ) cos φ d φ, (9)where G = (4 e / πh )[ E F + t ⊥ E F ] / / (cid:126) v F ; the angleof incidence φ is determined by tan φ = k y /α . It is notpossible to obtain the conductance analytically, thereforewe evaluate this integral numerically.The conductance is a periodic function of P (since thetransmission is) with period π . Fig. 3 shows a contourplot of the conductance for one period. As seen, the con-ductance has a sharp minimum at ε = 1 / − cos P ) : thisis due to the cone feature in the transmission which shiftsto higher energies with increasing P . Such a sharp min-imum was not present in the conductance of single-layergraphene when applied to the same δ -function potentialbarrier . FIG. 3: (Color online) (a) Contour plot of the conductance G . (b) Slices of G along constant P . IV. BOUND STATES OF A SINGLE δ -FUNCTION BARRIER The bound states here are states that are localised inthe x-direction close to the barrier but are free to movealong the barrier, i.e. in the y-direction. Such boundstates are characterised by the fact that the wave functiondecreases exponentially in the x direction, i.e., the wavevectors α and β are imaginary. This leads to e g e g = N e d e d , (10)which we can write as = Q e d e g e d e g , (11)where the matrix Q is the same as in Eq. (6). In or-der for this homogeneous algebraic set of equations tohave a nontrivial solution, the determinant of Q must bezero. This gives rise to a transcendental equation for thedispersion relation det Q = N N − N N = 0 , (12)which can be written explicitly [cos P + iµ sin P ][cos P + iν sin P ]+ ( α − β ) k y αβε sin P = 0 , (13)This expression is invariant under the transformations P → P + nπ, k y → − k y ,
3) ( ε, P ) → ( − ε, π − P ) . (14)Furthermore, there is one bound state for k y > and π/ < P < π . For P < π/ we can see that there isalso a single bound state for negative energies from thethird property above. From this transcendental formulaone can find the solution for the energy ε as function of k y numerically. We show the bound state by the solidred curve in Fig. 2(b). This state is bound to the po-tential in the x direction but moves as a free particlein the y direction. We have two such states, one thatmoves along the + y direction and one along the − y di-rection. The numerical solution approximates the curve ε = cos P [ − / / k y ) / ] . If one uses the × Hamiltonian one obtains the dispersion relation given inAppendix C by Eq. (C1). By solving this equation onefinds for each value of P two bound states one for pos-itive and one for negative k y . Moreover, for positive P these bands have a hole like behaviour and for negative P an electron like behaviour. Only for small P do theseresults coincide with those from the × Hamiltonian.The wave function ψ ( x ) of such a bound state is char-acterised by the coefficients e g , e g on the left, and e d and e d on the right side of the barrier. We can obtainthe latter coefficients by using Eq. (11), by assuming e g = 1 and afterwards normalising the total probabilityto unity in dimensionless units. The wave function ψ ( x ) to the left and right of the barrier can be determined fromthese coefficients by using Appendix A. In Figs. 2(c), (d)we show the probability distribution | ψ ( x ) | of a boundstate for a single δ -function barrier: in (c) we show itfor several k y values and in (d) for different values of P .One can see that the bound state is localised around thebarrier and is less smeared out with increasing k y . No-tice that the bound state is more strongly confined for P = π/ and that | ψ ( x ) | is invariant under the trans-formation P → π − P . V. TRANSMISSION THROUGH TWO δ -FUNCTION BARRIERS We consider a system of two barriers, separated bya distance L , with strengths P and P , respectively, asshown schematically in Fig. 4. We have L → Lt ⊥ / (cid:126) v F ≡ . L/nm which for L = 10 nm, v F = 10 m/s, and t ⊥ = 0 . eV equals . in dimensionless units. Thewave functions in the different regions are related as fol-lows FIG. 4: A system of two δ -function barriers with strengths P and P placed a distance L apart. ψ (0) = S ψ (0) , ψ (0) = S (cid:48) ψ ( L ) ,ψ (1) = S ψ ( L ) , ψ (0) = S S (cid:48) S ψ ( L ) , (15) where S (cid:48) = GM (1) G − represents a shift from x=0 tox=L and the matrices S and S are equal to the matrix N (cid:48) of Eq. (B3) with P = P and P = P , respectively.Using the transfer matrix N = G − S S (cid:48) S GM ( L ) weobtain the transmission T = | t | .In Fig. 5 the transmission T ( ε, k y ) is shown for par-allel (a), (b) and anti-parallel (b), (c) δ -function barrierswith equal strength, i.e., for | P | = | P | , that are sepa-rated by L = 10 nm, with P = 0 . π in (a) and P = 0 . π in (b). For P = π/ , the transmission amplitude t forparallel barriers equals − t for anti-parallel ones and thetransmission T is the same, as well the formula for thebound states. Hence panel (b) is the same for paralleland anti-parallel barriers. The contour plot of the trans-mission has a very particular structure which is very dif-ferent from the single-barrier case. There are two boundstates for each sign of k y , which are shown in panel (d) forparallel and panel (e) for anti-parallel barriers. For anti-parallel barriers these states have mirror-symmetry withrespect to ε = 0 but for parallel barriers this symmetry isabsent. For parallel barriers the change P → π − P willflip the spectrum of the bound state. The spectrum of thebound states extends into the low-energy transmissionregion and gives rise to a pronounced resonance. Noticethat for certain P values (Figs. 5(a) and 5(d)) the energydispersion for the bound state has a camelback shape forsmall k y , indicating free propagating states along the y direction with velocity opposite to that for larger k y val-ues. Contrasting Fig. 2(b) with Fig. 5(d)-(e)we see thatthe free-particle like spectrum of Fig. 2(b) for the boundstates of a single δ -function barrier is strongly modifiedwhen two δ -function barriers are present.From the transfer matrix we find that the transmissionis invariant under the change P → P + nπ and k y → − k y for parallel barriers, which was also the case of a singlebarrier, cf. Eq. (14). In addition, it is also invariant, foranti-parallel barriers, under the change P → π − P. (16)The conductance G is calculated numerically as in thecase of a single barrier. We show it for (anti-)parallel δ -function barriers of equal strengths in Fig. 6. Thesymmetry G ( P + nπ ) = G ( P ) of the single barrierconductance holds here as well. Further, we see thatfor anti-parallel barriers G has the additional symmetry G ( P ) = G ( π − P ) as the transmission does. VI. KRONIG-PENNEY MODEL
We consider an infinite sequence of equidistant δ -function potential barriers, i.e., a superlattice (SL), withpotential V ( x ) = P (cid:88) n δ ( x − nL ) . (17)As this potential is periodic the wave function should bea Bloch function. Further, we know how to relate the FIG. 5: (Color online) Panels (a), (b), and (c): contour plots of the transmission through two δ -function barriers of equalstrength P = | P | = | P | separated by a distance L = 10 nm. For parallel barriers we took P = 0 . π in (a) and P = 0 . π in(b). For anti-parallel barriers results are given for P = 0 . π in (b) and P = 0 . π in (c). The solid red curves in the whitebackground region is the spectrum for the bound states. Panels (d) and (e) show the dispersion relation of the bound statesfor various strengths | P | , respectively, for parallel and anti-parallel barriers. The thin black curves delimit the region wherebound states are possible.FIG. 6: (Color online) Contour plot of the conductance oftwo δ -function barriers with strength | P | = | P | = P andinter-barrier distance L = 10 nm. Panel (a) is for parallelbarriers and panel (c) for anti-parallel barriers. Panels (b)and (d) show the conductance, along constant P , extractedfrom panels (a) and (c), respectively. coefficients A of the wave function before the barrierwith those ( A ) after it, see Appendix B. The result is ψ ( L ) = e ik x L ψ (0) , A = N A , (18)with k x the Bloch wave vector. From these boundaryconditions we can extract the relation e − ik x L M ( L ) A = N A , (19)with the matrix M ( x ) given by Eq. (A4). The determi-nant of the coefficients in Eq. (19) must be zero, i.e., det[ e − ik x L M ( L ) − N ] = 0 . (20) If k y = 0 , which corresponds to the pure 1D case, onecan easily obtain the dispersion relation because the firsttwo and the last two components of the wave functiondecouple. Two transcendental equations are found cos k x L = cos αL cos P + 12 (cid:16) αε + εα (cid:17) sin αL sin P, (21a) cos k x L = cos βL cos P + 12 (cid:18) βε + εβ (cid:19) sin βL sin P. (21b)Since β is imaginary for < E < t ⊥ , we can write Eq.(21b) as cos k x L = cosh | β | L cos P − | β | + ε | β | ε sinh | β | L sin P, (22)which makes it easier to compare with the spectrum ofthe KP model obtained from the × Hamiltonian, seeEq. (3). The latter is given by the two relations cos k x L = cos κL + ( P/ κ ) sin κL, (23a) cos k x L = cosh κL − ( P/ κ ) sinh κL, (23b)with κ = √ ε . This dispersion relation, which has thesame form as the one for standard electrons, is not peri-odic in P and the difference from that of the four-bandHamiltonian is due to the fact that the former is not validfor high potential barriers. One can also contrast the dis-persion relations (21) and (23) with the correspondingone on single-layer graphene cos k x L = cos λL cos P + sin λL sin P, (24)where λ = E/ ( (cid:126) v F ) . This dispersion relation is also pe-riodic in P .In Fig. 7 we plot slices of the energy spectrum for k y =0 . There is a qualitative difference, between the four-band and the two-band approximation for P = π . Onlywhen P is small does the difference between the two 1Ddispersion relations become small. Therefore, we will nolonger present results from the × Hamiltonian thoughit has been used frequently due to its simplicity. Thepresent results indicate that one should be very carefulwhen using the × Hamiltonian in bilayer graphene.Notice that for P = 0 . π the electron and hole bandsoverlap and cross each other close to | k y | ≈ . π/L ) .That is, this is the spectrum of a semi-metal. Thesecrossing points move to the edge of the Brillouin zone(BZ) for P = π resulting in a zero-gap semiconductor.At the edge of the BZ the spectrum is parabolic for lowenergies.For k y (cid:54) = 0 , the dispersion relation can be written ex-plicitly in the form FIG. 7: (Color online) Slices of the spectrum of a KP SLwith L = 10 nm along k x , for k y = 0 , with P = 0 . π in(a) and (b) and P = π in (c) and (d). The results in (a)and (c) are obtained from the four-band Hamiltonian (1) andthose in (b) and (d) from the two-band one (3). The solid anddashed curves originate, respectively, from Eqs. (21a,23a) and(21b,23b). cos 2 k x L + C cos k x L + C / , (25)where C = − αL + cos βL ) cos P − ( d α + d β ) sin P, (26)and C = (2 + k y /ε ) + (2 − k y /ε ) cos αL cos βL + [( ε − k y ) + ε (2 ε − αL sin βL/ αβε − (cid:110) k y /ε − (2 + k y /ε ) cos αL cos βL + (cid:2) ε − / − k y + k y /ε (cid:3) sin αL sin βL/αβ (cid:111) cos 2 P + [ d α cos βL + d β cos αL ] sin 2 P, (27)with d α = (2 ε + 1) sin αL/α and d β = (2 ε −
1) sin βL/β .The wave vectors α = [ ε + ε − k y ] / and β = [ ε − ε − k y ] / are pure real or imaginary. If β becomes imaginary,the dispersion relation is still real ( β → i | β | and sin βL → i sinh | β | L ). Further, if α becomes imaginary, that is for α → i | α | , the dispersion relation is real. The dispersionrelation has the following invariance properties: E ( k x , k y , P ) = E ( k x , k y , P + 2 nπ ) , (28a) E ( k x , k y , P ) = − E ( π/L − k x , k y , π − P ) , (28b) E ( k x , k y , P ) = E ( k x , − k y , P ) . (28c)In Fig. 8 we show the lowest conduction and highest va- FIG. 8: (Color online) SL Spectrum for L = 10 nm, thelowest conduction and highest valence band for P = 0 . π in(a) and P = 0 . π in (b), are shown.FIG. 9: (Color online) SL spectrum for L = 10 nm . Thedashed blue, solid red, and dash-dotted purple curves are,respectively, for strengths P = 0 . π , P = 0 . π , and P =0 . π . (a) shows the spectrum vs k x for k y = 0 while (b)shows it vs k y for k x at the value where the bands cross. Theposition of the touching points and the size of the energy gapare shown in (c) as a function of P . The dash-dotted, bluecurve and the solid, black curve show k x, and the energyvalue of the touching points, respectively. For P > P c a gapappears; the energies of the conduction band minimum and ofthe valence band maximum are shown by the red and purple,solid curve, respectively. lence bands of the energy spectrum of the KP model for P = 0 . π in (a) and P = 0 . π in (b). The former hastwo touching points which can also be viewed as over-lapping conduction and valence bands as in a semi-metaland the latter exhibits an energy gap. In Fig. 9 slices ofFigs. 8(a), (b) are plotted for k y = 0 . The spectrum ofbilayer graphene has a single touching point at the origin.When the strength P is small, this point shifts away fromzero energy along the k x axis with k y = 0 and splits intotwo points. It is interesting to know when and wherethese touching points emerge. To find out we observethat at the crossing point both relations (21) should befulfilled. If these two relations are subtracted we obtainthe transcendental equation αL − cos βL ) cos P + (1 /
2) ( d α − d β ) sin P, (29)where d α = (2 ε + 1) sin αL/α and d β = (2 ε −
1) sin βL/β .We can solve Eq. (28) numerically for the energy ε . Forsmall P and small L this energy can be approximated by ε = P/L . Afterwards we can put this solution back intoone of the dispersion relations to obtain k x .In Figs. 9(a), (b) we show slices along the k x axisfor k y = 0 and along the k y axis for the k x value of atouching point, k x, . We see that as the touching pointsmove away from the K point, the cross sections show amore linear behaviour in the k y direction. The position ofthe touching points is plotted in Fig. 9(c) as a functionof P . The dash-dotted blue curve corresponds to thevalue of k x, (right y axis), while the energy value of thetouching point is given by the black solid curve. Thistouching point moves to the edge of the BZ which occursfor P = P c ≈ . π . At this point a gap opens (theenergies of the top of the valence band and of the bottomof the conduction band are shown by the lower purple andupper red solid curve, respectively) and increases with P .Because of property 2) in Eq. 28(b) we plot the resultsonly for P < π/ . We draw attention to the fact thatthe dispersion relation differs to large extent for large P from the one that results from the × Hamiltoniangiven in Appendix C. This is already apparent from thefact that the dispersion relation does not exhibit any ofthe periodic in P behaviours given by Eqs. (28a) and(28b).An important question is whether the above period-icities in P still remain approximately valid outside therange of validity of the KP model. To assess that webriefly look at a square-barrier SL with barriers of fi-nite width W b and compare the spectra with those of theKP model. We assume the height of the barrier to be V / (cid:126) v F = P/W b , such that V W b / (cid:126) v F = P . The SL pe-riod we use is nm and the width W b = 0 . L = 2 . nm.For P = π/ the corresponding height is then V ≈ t ⊥ .To fit in the continuum model we require that the poten-tial barriers be smooth over the carbon-carbon distancewhich is a ≈ . nm. In Fig. 10 we show the spectra forthe KP model and this SL. Comparing (a) and (b) wesee that for P small the difference between both modelsis rather small. If we take P = π/ though, this dif-ference becomes large, especially for the first conductionand valence minibands, as shown in panels (c) and (d). The latter energy bands are flat for large k y in the KPmodel, while they diverge from the horizontal line (E=0)for a finite barrier width. From panel (f), which showsthe discrepancy of the SL minibands between the exactones and those obtained from the KP model, we see thatthe spectra with P = 0 . π are closer to the KP modelthan those for P = π/ . Fig. 10(e) demonstrates that theperiodicity of the spectrum in P within the KP model,i.e., its invariance under the change P → P + 2 nπ , ispresent only as a rough approximation away from it. VII. EXTENDED KRONIG-PENNEY MODEL
In this model we replace the single δ -function barrierin the unit cell by two barriers with strengths P and − P . Then the SL potential is given by V ( x ) = P (cid:88) n δ ( x − nL ) − P (cid:88) n δ ( x − ( n +1 / L ) . (30)Here we will restrict ourselves to the important case of P = P . For this potential we can also use Eq. (20)of Sec. VI, with the transfer matrix N replaced by theappropriate one of Sec. V.First, let us consider the spectrum along k y = 0 whichis determined by the transcendental equations cos k x L = cos αL cos P + D α sin P, (31a) cos k x L = cos βL cos P + D β sin P, (31b)with D γ = (cid:2) ( γ + ε ) cos γL − γ + ε (cid:3) / γ ε . It ismore convenient to look at the crossing points becausethe spectrum is symmetric around zero energy. This fol-lows from the form of the potential (its spatial averageis zero) or from the dispersion relation (31a): the change ε → − ε entails α ↔ β and the crossings in the spectrumare easily obtained by taking the limit ε → in one ofthe dispersion relations. This gives the value of k x at thecrossings k x, = ± arccos[1 − ( L /
8) sin P ] /L, (32)and the crossing points are at ( ε, k x , k y ) = (0 , ± k x, , .If the k x, value is not real, then there is no solution atzero energy and a gap arises in the spectrum. From Eq.(31a) we see that for sin P > /L a band gap arises.In Fig. 11 we show the lowest conduction and highestvalence band for (a) P = 0 . π , and (b) P = 0 . π . Ifwe make the correspondence with the KP model of Sec.V we see that this model leads to qualitatively similar(but not identical) spectra shown in Figs. 8(a) and 8(b):one should take P twice as large in the correspondingKP model of Sec. V in order to have a similar spectrum.Here we have the interesting property that the spectrumexhibits mirror symmetry with respect to ε = 0 whichmakes the analysis of the touching points and of the gapeasier.In Fig. 12 we plot the k x value (dash-dotted, bluecurve) of the touching points k x, versus P , if there is FIG. 10: (Color online) Spectrum of a SL with L = 50 nm, (a) and (b) are for P = 0 . π and (c) and (d) are for P = π/ . (a),(c) and (e) are for a rectangular-barrier SL with W b = 0 . L and u = P/W b , while (b) and (d) are for the KP model. (e) showsthe spectrum for u corresponding to P = (1 / π ; the dashed curves show the contours of the spectrum in (c) for P = π/ .(f) Shows the discrepancy of the SL minibands between the exact ones and those obtained from the KP model, averaged over k space (where we used k y L/π = 6 as a cut-off). The conduction (valence) minibands are numbered with positive (negative)integers.FIG. 11: (Color online) The first conduction and valenceminibands for the extended KP model for L = 10 nm with P = 0 . π in (a) and P = 0 . π in (b).FIG. 12: (Color online) Plot of the ± k x, values, for which theminibands touch each other, as a function of P (dash-dotted,blue curve), and the size of the band gap E gap (solid, redcurve). The calculation is done for the extended KP modelwith L = 10 nm. no gap, and the size of the gap E gap (solid, red curve)if there is one. The touching points move toward theBZ boundary with increasing P . Beyond the P value forwhich the boundary is reached, a gap appears betweenthe conduction and valence minibands. VIII. CONCLUSIONS
We investigated the transmission through single anddouble δ -function potential barriers on bilayer grapheneusing the four-band Hamiltonian. The transmission andconductance are found to be periodic functions of thestrength of the barriers P = V W b / (cid:126) v F with period π .The same periodicity was previously obtained for suchbarriers on single-layer graphene . We emphasise thatthe periodicity obtained here implies that the transmis-sion satisfies the relation T ( k x , k y , P ) = T ( k x , k y , P + nπ ) for arbitrary values of k x , k y , P , and integer n . Inprevious theoretical work on graphene and bilayergraphene Fabry-Pérot resonances were studied and T = 1 was found for particular values of α , the elec-tron momentum inside the barrier along the x axis. Fora rectangular barrier of width W and Schrödinger-typeelectrons, Fabry-Pérot resonances occur for αW = nπ and E > V as well as in the case of a quantum wellfor E > , V < . In graphene, because of Klein tun-nelling, the latter condition on energy is not needed. Be-cause α depends on the energy and the potential barrierheight in the combination E − V , any periodicity of T in the energy is equivalent to a periodicity in V if noapproximations are made, e.g., E (cid:28) V , etc. Althoughthis may appear similar to the periodicity in P , thereare fundamental differences. As shown in Ref. 21, theFabry-Pérot resonances are not exactly described by thecondition αW = nπ (see Fig. 3 in Ref. 21) while the pe-riodicity of T in the effective barrier strength is exactly nπ . Furthermore, the Fabry-Pérot resonances are foundfor T = 1 , while the periodicity of T in P is valid for anyvalue of T between and .Further, we studied the spectrum of the KP modeland found it to be periodic in the strength P with period π . In the extended KP model this period reduces to π . This difference is a consequence of the fact that forthe extended SL the unit cell contains two δ -functionbarriers. These periodicities are identical to the onefound earlier in the (extended) KP model on single-layergraphene. We found that the SL conduction and valenceminibands touch each other at two points or that thereis a energy gap between them. In addition, we found asimple relation describing the position of these touchingpoints. None of these periodic behaviours results fromthe two-band Hamiltonian; this clearly indicates that thetwo-band Hamiltonian is an incorrect description of theKP model in bilayer graphene. In general, results de-rived from these two tight-binding Hamiltonians agreewell only for small energies . The precise energy rangesare not explicitly known and may depend on the partic-ular property studied. For the range pertaining to thefour-band Hamiltonian ab-initio results indicate thatit is approximately from − eV to + . eV.The question arises whether the above periodicities in P survive when the potential barriers have a finite width.To assess that we briefly investigated the spectrum of arectangular SL potential with thin barriers and comparedit with that in the KP limit. We showed with some ex-amples that for specific SL parameters the KP model isacceptable in a narrow range of P and only as a rough ap-proximation away from this range. The same conclusionholds for the periodicity of the KP model.The main differences between the results of this workand those of our previous one, Ref. 17, are as follows. Incontrast to monolayer graphene we found here that:1) The conductance for a single δ -function potential bar-rier depends on the Fermi energy and drops almost tozero for certain values of E and P . 2) The KP model(and its extended version) in bilayer graphene can opena band gap; if there is no such gap, two touching pointsappear in the spectrum instead of one. 3) The Dirac linefound in the extended KP model in single-layer grapheneis not found in bilayer graphene. Acknowledgments
This work was supported by IMEC, the FlemishScience Foundation (FWO-Vl), the Belgian SciencePolicy (IAP), and the Canadian NSERC Grant No.OGP0121756.
Appendix A: Eigenvalues and eigenstates for aconstant potential
Starting with the Hamiltonian (1) for a one-dimensional potential V ( x, y ) = V ( x ) , the time-independent Schrödinger equation H ψ = Eψ leads to − i ( ∂ x − k y ) ψ B = ε (cid:48) ψ A − ψ B (cid:48) , − i ( ∂ x + k y ) ψ A = ε (cid:48) ψ B , − i ( ∂ x + k y ) ψ A (cid:48) = ε (cid:48) ψ B (cid:48) − ψ A , − i ( ∂ x − k y ) ψ B (cid:48) = ε (cid:48) ψ A (cid:48) , (A1)The spectrum and the corresponding eigenstates can be obtained, for constant V ( x, y ) = V , by progressiveelimination of the unknowns in Eq. (A1) and solutionof the resulting second-order differential equations. Theresult for the spectrum is ε = u + 1 / ± (cid:112) / k ,ε = u − / ± (cid:112) / k . (A2)The unnormalised eigenstates are given by the columnsof the matrix GM , where G = f α + f α − f β + f β − − − f α − f α + − f β − − f β + , (A3)with f α,β ± = − i ( k y ± i ( α, β )) /ε (cid:48) ; α = [ ε (cid:48) + ε (cid:48) − k y ] / and β = [ ε (cid:48) − ε (cid:48) − k y ] / are the wave vectors. M is givenby M = e iαx e − iαx e iβx
00 0 0 e − iβx . (A4)The wave function in a region of constant potential is alinear combination of the eigenstates and can be written Ψ( x ) = ψ A ψ B ψ B (cid:48) ψ A (cid:48) = GM ABCD . (A5)We can reduce its complexity by the linear transforma-tion Ψ( x ) → R Ψ( x ) where R = 12 − −
11 0 1 00 1 0 1 , (A6)which transforms Ψ( x ) to Ψ( x ) = (1 / ψ A − ψ B (cid:48) , ψ B − ψ A (cid:48) , ψ A + ψ B (cid:48) , ψ B + ψ A (cid:48) ) T . Then the basis functions aregiven by the columns of GM with G = α/ε (cid:48) − α/ε (cid:48) − ik y /ε (cid:48) − ik y /ε (cid:48) − ik y /ε (cid:48) − ik y /ε (cid:48) β/ε (cid:48) − β/ε (cid:48) . (A7)The matrix M is unchanged under the transformation R and the new Ψ( x ) fulfils the same boundary conditionsas the old one. Appendix B: The transfer matrix
We denote the wave function to the left of, inside, andto the right of the barrier by ψ j ( x ) = G j M j A j , with j = , , and , respectively. Further, we have G = G and M = M . The continuity of the wave function at x = 0 and x = W b gives the boundary conditions ψ (0) = ψ (0) and ψ ( W b ) = ψ ( W b ) . In explicit matrix notation thisgives G A = G A and G M ( W b ) A = M ( W b ) G A ,where A = G − G M − ( W b ) G − G M ( W b ) A . Thenthe transfer matrix N can be written as N = G − G M − ( W b ) G − G M ( W b ) . Let us define N (cid:48) = G M − ( W b ) G − , which leads to ψ (0) = N (cid:48) ψ ( W b ) .To treat the case of a δ -function barrier we take thelimits V → ∞ and W b → such that the dimensionlesspotential strength P = V W b / (cid:126) v F is kept constant. Then G and M ( W b ) simplify to G = − − , (B1) M ( W b ) = e iP e − iP e iP
00 0 0 e − iP , (B2)and N (cid:48) becomes N (cid:48) = cos P i sin P i sin P cos P P i sin P i sin P cos P . (B3) Appendix C: Results for the × Hamiltonian
Using the × Hamiltonian (3) instead of the × one can sometimes lead to unexpectedly different re- sults; below we give a few examples. In a slightly mod-ified notation pertinent to the × Hamiltonian we set α = [ − ε + k y ] / , β = [ ε + k y ] / , and use the samedimensionless units as before.Bound states for a single δ -function barrier u ( x ) = P δ ( x ) , without accompanying propagating states, arepossible if k y = 0 or k y > | ε | . In the former case thesingle solution is ε = − sign ( P ) P / . In the latter onethe dispersion relation is ε ( P + 2 α )( P − β ) + 2 P k y ( αβ − k y ) = 0 . (C1)The dispersion relation for the KP model obtainedfrom the × Hamiltonian is cos(2 kL ) + 2 F cos( kL ) + F = 0 , (C2)where F = − cosh( βL ) − cosh( αL ) + P β sinh( βL ) − P α sinh( αL ) ,F = 1 αβε (cid:8) αβ ( ε + k y P / β cosh( βL ) (cid:2) α (2 ε − k y ) cosh( αL ) + ε P sinh( αL ) (cid:3) − P βL ) (cid:2) α ( ε − k y / P sinh( αL ) + 2 ε α cosh( αL ) (cid:3)(cid:27) . (C3) K. S. Novoselov, A. K. Geim, S. V. Morozov, D. Jiang, Y.Zhang, S. V. Dubonos, I. V. Grigorieva, and A. A. Firsov,Science , 666 (2004). O. Klein, Z. Phys. , 157 (1929). M. I. Katsnelson, K. S. Novoselov, and A. K. Geim, NaturePhysics , 620 (2006). N. Stander, B. Huard, and D. Goldhaber-Gordon, Phys.Rev. Lett. , 026807 (2009); A. F. Young and P. Kim,Nature Phys. , 222 (2009). A. H. Castro Neto, F. Guinea, N. M. R. Peres, K. S.Novoselov, and A. K. Geim, Rev. Mod. Phys. , 109(2009); C. W. J. Beenakker, ibid , 1197 (2008). J. M. Pereira Jr., P. Vasilopoulos, and F. M. Peeters, Appl.Phys. Lett. , 132122 (2007). E. McCann, D. S. L. Abergel, and V. I. Fal’ko, Solid StateComm. , 110 (2007). C.-H. Park, L. Yang, Y.-W. Son, M. L. Cohen, and S. G.Louie, Nature Phys. , 213 (2008). S. Ghosh and M. Sharma, J. Phys.: Condens. Matter ,292204 (2009). Y. P. Bliokh, V. Freilikher, S. Savel’ev, and F. Nori, Phys.Rev. B , 075123 (2009). I. Snyman, Phys. Rev. B , 054303 (2009). R. Nasir, K. Sabeeh, and M. Tahir, Phys. Rev. B ,085402 (2010). M. Barbier, F. M. Peeters, P. Vasilopoulos, and J. M.Pereira Jr., Phys. Rev. B , 115446 (2008). C. Bai and X. Zhang, Phys. Rev. B , 075430 (2007). M. Barbier, P. Vasilopoulos, F. M. Peeters, and J. M.Pereira Jr., Phys. Rev. B , 155402 (2009). C. Kittel,
Introduction to Solid State Physics, 5th edn ,(John Wiley & Sons, Inc., 1976). M. Barbier, P. Vasilopoulos, and F. M. Peeters, Phys. Rev.B , 205415 (2009). E. McCann and V. I. Fal’ko, Phys. Rev. Lett. , 086805(2006). A. Matulis and F. M. Peeters, Phys. Rev. B , 115423(2008). A. V. Shytov, M. S. Rudner, and L. S. Levitov, Phys. Rev.Lett. , 156804 (2008); P. G. Silvestrov and K. B. Efe- tov, ibid. , 016802 (2007); F. Young and Philip Kim,Nature Physics , 222 (2009). I. Snyman and C. W. J. Beenakker, Phys. Rev. B ,045322 (2007). M. Ramezani Masir, P. Vasilopoulos, and F. M. Peeters, Phys. Rev. B , 115417 (2010). S. Latil and L. Henrard, Phys. Rev. Lett.97