Control over band structure and tunneling in Bilayer Graphene induced by velocity engineering
CControl over band structure and tunneling in Bilayer Grapheneinduced by velocity engineering
Hosein Cheraghchi ∗ and Fatemeh Adinehvand School of Physics, Damghan University, 36716-41167, Damghan, Iran (Dated: November 9, 2018)The band structure and transport properties of massive Dirac Fermions in bilayer graphene withvelocity modulation in space are investigated in presence of the previously created band gap. Itis pointed out that the velocity engineering is considered as a factor to control the band gap ofsymmetry-broken bilayer graphene. The band gap is direct and independent of velocity value ifvelocity modulated in two layers is set up equally. Otherwise, in the case of interlayer asymmetricvelocity, not only the band gap is indirect, but also the electron-hole symmetry fails. This band gapis controllable by the ratio of the velocity modulated in the upper layer to the velocity modulatedin the lower layer. In more detail, the shift of momentum from the conduction band edge to thevalence band edge can be engineered by the gate bias and velocity ratio. A transfer matrix methodis also elaborated to calculate four-band coherent conductance through a velocity barrier possiblysubjected to a gate bias. Electronic transport depends on the ratio of velocity modulated insidethe barrier to the one for surrounding regions. As a result, a quantum version of total internalreflection is observed for enough thick velocity barriers. Moreover, a transport gap originating fromthe applied gate bias is engineered by modulating velocity of the carriers in the upper and lowerlayers.
PACS numbers: 72.80.Vp,73.22.Pr,73.23.Ad,73.63.-b
I. INTRODUCTION
Charge carriers in monolayer graphene at low energies,near the neutrality point, are described by Dirac fermionswith a velocity that is independent of wavelength . Thisunique property proposes an analogous between Diracfermions and electromagnetic or mechanical waves in op-tics and acoustics. Furthermore, this brings several un-usual electronic properties such as anomalous integer and fractional quantum Hall effects, electronic focus-ing by means of a rectangular potential barrier (Veselagolensing) , Klein tunneling and minimal conductivity .Spatial modulation of wave velocity has been originallystudied in optics, acoustics and recently in photonic crys-tals . The idea can be also applied for Dirac fermionwaves by defining a velocity barrier as the region in whichthe Fermi velocity differs from the one in the surround-ing background. In analogous with optics, some opticalrules are expected to be valid for massless Dirac fermionwaves propagating in monolayer graphene sheets .There are several ways to engineer the Fermi veloc-ity ( v F ) by means of a control over the electron-electroninteraction in graphene. Enhancement in the electron-electron interaction induces an increase in the Fermivelocity . Furthermore, an enhancement in v F whichis logarithmic in the carrier concentration n has beenestablished in experiments and also described by therenormalization group theory . Modifications of cur-vature of graphene sheet , periodic potentials and di-electric screening are some of propositions for en-gineering v F via the electron-electron interaction. The v F of graphene is inversely proportional to the dielec-tric constant of the environment embedding graphenesheet . Structures with velocity modulation in space FIG. 1: a) Schematic diagram of bilayer graphene junctionin AB stacking with velocity modulation in space. b) Energyband structure of bilayer graphene for different regions withdifferent velocities. At the same time which vertically gatebias δ is present, velocity may be experimentally modulatedin each layer of BLG. The ratio of velocity modulated in theupper layer ( v u ) to the lower layer ( v d ) controls the feature ofthe spectrum as well as tunneling through a velocity barrier. can be also made by application of appropriate doping or placing a grounded metal plane as a screening planeclose to graphene . In the presence of the screeningplanes, speed of carriers is smaller than the speed at iso-lated graphene sheet. Recently, in a 2d electron gas, anartificial graphene has been proposed by modulating aperiodic potential of honeycomb symmetry . Electronsin artificial graphene sheets behave like massless Dirac a r X i v : . [ c ond - m a t . m e s - h a ll ] N ov fermions with a tunable Fermi velocity.The electronic properties of monolayer graphene sheetswith spatial modulation of the Fermi velocity have beeninvestigated in literature . However, the elec-tronic properties in bilayer graphene (BLG) with an in-terlayer asymmetric velocity have not been elucidated indetail so far. There are numbers of different experimentsin which a controllable direct band gap is observed ingated bilayer graphene . However, the amount ofcurrent in the off-state still remains high . Thisoff-current has been attributed to several sources suchas edge states , the presence of disorder , coexistenceof massive and massless Dirac fermions in twisted AA-stacking bilayer graphene grown on SiC . Strain isthe other known factor which controls the band gap inBLG .In this work, we point out that the velocity modifi-cation in symmetry-broken BLG, as an inevitable ex-perimental factor, is able to control the band gap. Inthe absence of the gate bias δ = 0, symmetric or asym-metric velocity modulation in two layers is not able tocreate a gap in the band structure. However, the previ-ously created gap δ (cid:54) = 0 can be controlled by the ratioof modulated velocity in the upper layer to the lowerlayer η . The band gap is direct if velocity of itinerantquasi-particles in each layer is set up equally. This gap isindependent of velocity, while the momentum attributedto the band gap is inversely proportional to the velocity.On the other hand, the band gap is indirect for non-equalvelocities modulated in layers. In this case, the bandstructure and subsequently the band gap are controlledby η . The Shift of momentum from the conduction bandedge to the valence band edge depends on the gate biasand velocity ratio. Moreover, the electron-hole symme-try fails when η (cid:54) = 1. This kind of control over the bandstructure which is induced by different velocity modula-tion in each layer, opens up the possibility of new deviceapplications in nanoelectronics. More importantly, in aBLG under application of gate bias, experiments have tobe care about the transition of direct to indirect bandgap. This transition can be induced by modification ofvelocity in layers originating from several experimentalrequirements such as coating a metallic gate electrode,changing carrier concentration by using application of agate voltage, strain and etc.To manifest such a control over the gap, we develop atransfer matrix approach to investigate transport proper-ties through the velocity barrier subjected to a gate biasin BLG. A schematic diagram of the proposed system ispresented in Fig.1 which indicates simultaneous velocityand electrostatic junction. The proposed method is basedon a four-band Hamiltonian for AB stacking . As aresult, similar to monolayer graphene , a total in-ternal reflection occurs for Dirac fermion waves hitting ona thick barrier at the angles of incidence greater than acritical angle. Moreover, it is observed that the trans-port gap depends on the velocity ratio η at large gatebias. This gap is induced by application of a symmetry breaking factor in the barrier region.We organize this paper as the following: In section II,we present four-band Hamiltonian and a general formulafor deriving the spectrum in the presence of velocity mod-ulation in addition to vertically applied gate bias. Thenin section III, we switch to calculate transport propertiesthough a velocity junction possibly subjected to an ex-ternal gate bias in generic form. Finally, the last sectionincludes the results. II. HAMILTONIAN AND BAND STRUCTUREIN PRESENCE OF INTERLAYER ASYMMETRY
The four-band Hamiltonian of bilayer graphene closeto the Dirac point (i.e say the valley of K point) for ABstacking is described as the follow: H = (cid:18) − i (cid:126) v u ( σ. ∇ ) † + V u I FF − i (cid:126) v d ( σ. ∇ ) + V d I (cid:19) (1)where F = (cid:18) t
00 0 (cid:19) , − i (cid:126) vσ. ∇ = (cid:18) π † π (cid:19) and I is the unit matrix. Here, π = − i (cid:126) v ( ∂ x − k y ), t = 390 meV is the coupling energy between the layers. V u = V + δ and V d = V − δ describe an asymmet-ric factor which can be applied by a vertically gate biasor doping. This interlayer asymmetry emerges as a dif-ference between on-site energies belonging to each layer.Another interlayer asymmetry can be induced by differ-ent modulation in the velocity of itinerant quasi-particlesin the upper and lower layers, v u = ξ u v F and v d = ξ d v F respectively. v F is the commonly Fermi velocity used forgraphene. V is the gate voltage applied on both layerssetting up to zero. 2 δ is the potential difference betweenthe upper and lower layers induced by a gate bias ordoping. The eigen function of the above Hamiltonian is written as Ψ = (cid:0) ψ uA ψ uB ψ dB ψ dA (cid:1) (cid:62) . By solvingthe eigenvalue equation of H Ψ = E Ψ, band structurecan be calculated in the gapless case or in the presenceof previously applied gate bias.At the same time which vertically gate bias is present,velocity may be experimentally modulated in each layerof BLG. In a gaped BLG, we will show that there ispossibility for engineering the previously created gap byusing a velocity modulation in each layer. In the presenceof a gate bias accompanied with an interlayer asymmetryin velocity, the BLG’s spectrum can be extracted fromthe following equation which presents k ( E ) . k ( E ) = [ a ± (cid:112) a − b ] /v u a ( E, η, δ ) = [ η ( E − δ ) + ( E + δ ) ] / b ( E, η, δ ) = η ( E − δ )( E − δ − t ) η = ξ u /ξ d (2) FIG. 2: Bulk band structure of bilayer graphene for sev-eral values of velocity modulated in the lower layer ξ d whilevelocity in the upper layer ( ξ u = 1) is fixed. Here, ξ = v/v F . If the gate voltage is nonzero V (cid:54) = 0, functions of a and b in the above equation depend on ε = E − V instead of E .Based on the velocity ratio ( η ) and in the presence of agate bias (2 δ ), we will indicate that BLG has two differentbehaviors. For η = 1, BLG behaves as a semiconductorwith direct band gap, while for η (cid:54) = 1, it behaves as asemiconductor with indirect band gap. In the case of δ = 0, independent of η , there is no gap in BLG. Bulkband structure calculated by the above equation for δ = 0is shown in Fig.(2) and for δ (cid:54) = 0 in Fig.(3) and Fig.(4a).To investigate the behavior of the energy gap E gap , onecan simply derive the following conditions to emerge theextermum points of E ( k ). Based on Eq.(2), there are twoconditions to satisfy the extermum condition ∂E/∂k = 0: (cid:40) b = 0 ∀ k = 0 b = a ∀ k = ± v − u √ a (3)An immediate result from Eq.3 is that the energy gap E gap ( η ) = E c ( η ) − E v ( η ) depends on the velocity ratio η , where E c and E v are the conduction and valence bandedges; however, the momentum attributed to the con-duction k c and valence k v band edges depend on bothvariables of ξ u and ξ d . The condition b = 0 results infour eigenvalues of Hamiltonian at the Dirac point cal-culated as E = ± δ and ± ( t + δ ). These eigenvalues andconsequently the energy gap appeared at k = 0 are in-dependent of the velocity ratio η . The condition b = a leads to the energy gap at the k-points derived by thefollowing equation: k c/v ( ξ u , ξ d ) = ± v − u (cid:112) a ( E c/v , η, δ ). A. Gapless band structure in presence of interlayersymmetric potential
Let us first concentrate on the gapless case with noexternal gate biasing ( δ = 0) which conserves chiral sym-metry. Based on Eq.(2), the four band spectrum for a FIG. 3: Bulk band structure of bilayer graphene for the caseof the same velocity modulated in the upper and lower layerswhenever a band gap is previously created by application ofa vertically gate bias. Note that in this case, velocity in twolayers are equal to each other but possibly can be differentfrom v F . In other word, ξ = ξ u = ξ d . BLG with a tunable velocity in each layer ( v u (cid:54) = v d ) canbe derived as the following , E = ± (cid:114) ϕ ( k ) + ( − (cid:15) (cid:113) ϕ ( k ) − v u v d k ϕ ( k ) = (( v u + v d ) k + t ) / . (4)where (cid:15) = 1 and 2 are attributed to the low and highenergy bands, respectively. In the case of (cid:15) = 1, thereis no band gap at the Fermi Dirac point ( k = 0). Thewhole spectrum is robust against the exchange of v u by v d . This robustness can also be derived by exchanging of η → /η in Eq.(2). In this case, the only real solutionsfor the extermum points derived by Eq.(3), are E = 0and ± t which emerge at k = 0.The chiral symmetry is conserved even though quasi-particles have different velocities in each layer. In thiscase, modulation of velocities in each layer just changesthe effective mass of quasi-particles. Fig.(2) shows theenergy bands of BLG with different velocities in eachlayer. The band structure is symmetric and behaves asa parabolic form. As a conclusion, without any applica-tion of potential difference, only interlayer asymmetry invelocity is not able to break the electron-hole symmetry. B. Band structure in presence of interlayerasymmetry in potential but symmetry in velocity:Direct band gap
In the presence of an interlayer asymmetric factor suchas an external gate bias ( δ (cid:54) = 0) and in the special case ofthe same velocities setting up on each layer v u = v d = v ( η = 1), the four band spectrum is described as : E = ( vk ) + δ + t / − (cid:15) (cid:112) ( vk ) (4 δ + t ) + t / k (5)As shown in Fig.(3), the low energy band (cid:15) = 1 dis-plays a Mexican hat shape. Despite turning external gatebias on, the band structure still remains symmetric giv-ing rise the electron-hole symmetry. The functions of a ( E, δ ), b ( E, δ ) defined in Eq.(2) are independent of η .Therefore, the band gap is independent of the velocitywhich is modulated in layers. The requisite condition forderiving the band gap ( b = a in Eq.(3)) results in asymmetric solution for the conduction and valence bandedges, E c = − E v = tδ √ δ + t ∀ k (cid:54) = 0 (6)So the band gap is written as E gap = 2 E c . At k = 0,the gap is fixed to the value 2 δ . Because a ( E v ) = a ( E c ),one can conclude that the momentum of the conductionand valence band edges emerge at the same point k c = k v = k gap from the center of valley. k gap = ± δv F (cid:114) t + 2 δ t + 4 δ ξ (7)Consequently, the band gap is direct and the momen-tum attributed to the gap is inversely proportional to thevelocity ξ . For the limit of low external gate bias δ (cid:28) t ,the band gap tends to the gate bias E gap → δ . However,for large potential differences δ (cid:29) t , the band gap tendsto saturate at the interlayer hopping energy E sat.gap → t .For both limits, the momentum attributed to the bandgap behaves as k gap ∝ δ/v . For the case of slower ve-locity ξ <
1, the effective mass at the conduction andvalence band edges is heavier than the effective mass forthe one with faster velocity ξ > C. Band structure in presence of interlayerasymmetry in potential and velocity: Indirect bandgap
In this case, interlayer asymmetry is applied on bothof electrostatic potential and also velocity of itinerantquasi-particles ( v u (cid:54) = v d ). In the case of η (cid:54) = 1 and δ (cid:54) = 0,there is an asymmetry between the conduction and va-lence bands of the spectrum giving rise the electron-holeasymmetry . Consequently, the conduction and valenceband edges are appeared at asymmetric energy pointsmeasuring from the band center E = 0. As a result,the momentum attributed to the conduction and valenceband edges emerges at different points, k c (cid:54) = k v . There-fore, the band gap is indirect . The shift of momen-tum from the conduction band edge to the valence bandedge (∆ k = k c − k v ) depends on velocity in each layer.Although the band gap just depends on the velocity ra- FIG. 4: a) Bulk band structure of bilayer graphene for aninterlayer asymmetry in velocity whenever a band gap is pre-viously created by application of a vertically gate bias. Notethat velocity modulation is ξ u for the upper layer and fixedfor the lower layer ξ d = 1. For more detail, the electro-static potential on the upper and lower layers are set (+ δ )and ( − δ ), respectively. Here the potential difference is set2 δ = 400 meV . b) The momentum of the conduction k c andvalence k v band edges, and also the momentum shift from theconduction band edge to the valence band edge ∆ k = k c − k v as a function of ξ u . c) The conduction E c and valence E v band edges and also the energy gap E gap in terms of ξ u (orhere η ). tio η , however, in a fixed velocity ratio, the whole featureof the spectrum is sensitive to both values of velocity at-tributed to the upper and lower layers. Let us set velocityof the lower layer to be fixed as ξ d = 1 while ξ u is tun-able. The lower and upper layers are characterized bythe electrostatic potential of − δ and δ , respectively.The asymmetric band structure is represented in FIG. 5: a) The energy gap in terms of the potential differencebetween the upper and lower layer for several velocity ratio.The inset figure shows the saturated band gap in terms of thevelocity ratio for δ = 4000 meV . The numerical calculationshown in the inset figure confirms the analytical derivation ofEq.8. b) The momentum shift from the conduction band edgeto the valence band edge in terms of the potential differencebetween the upper and lower layer for several velocity ratio. Fig.(4a) for three values of velocity of the upper layer ξ u . Although the band structure is asymmetric, its formpreserves the ’ Mexican hat ’ shape. In the appendix A,we have provided a comparison between the electron-holeasymmetry arising from the full Hamiltonian of BLG andthe dominant Hamiltonian which is considered in thiswork.Fig.(4b) shows the momentum attributed to the con-duction k c and valence k v band edges and also their mo-mentum shift ∆ k in terms of ξ u . As it is obviously ob-served, both of k c and k v decreases with ξ u . Moreover,their curves intersect each other at η = 1 which resultsin the direct band gap. However, for all values of η (cid:54) = 1,the band gap is indirect. For ξ u <
1, the momentumshift of k c away from the Dirac point is larger than themomentum shift of k v .By finding roots of Eq.(3), the conduction and valenceband edges are computed in terms of system parameters.Fig.(4c) indicates dependence of E c , E v and also E gap on the velocity ratio η . The curves related to E c and E v never intersect each other. In all ranges of η , E c > E v . SoBLG always behaves as a semiconductor, not metal norsemi-metal. The energy gap has a maximum at η = 1 inwhich the gap is direct. A sharp variation of E gap with η is seen for the range of η <
1. Parameters of η cr.c and η cr.v are those critical velocity ratios in which E c or E v crossthe band center E = 0. The curvature width of function E gap ( η ) is measured by ∆ η cr. = η cr.v − η cr.c . The critical velocity ratio for the valence and conduction band edgesis derived as the following form: η cr.v/c = 1 + 2( t/δ ) [1 ± (cid:112) δ/t ) ]. In both limits of δ (cid:28) t and δ (cid:29) t , thewidth of the peak which emerges in E gap ( η ), tends to∆ η cr. → t/δ . As a conclusion, for large gate bias δ ,there is a sharp variation in the energy gap as a functionof the velocity ratio. In large velocity ratio η → ∞ , theasymptotic solution of Eq.(3) for the conduction bandedge is E c → δ . In this limit, the momentum attributedto the conduction band edge behaves as a power law with v u ; k c → δ/v u . In the opposite limit of η →
0, theasymptotic solution for the valence band edge is E v →− δ . So, the momentum attributed to the valence bandedge tends to the constant; k v → δ/v F .Although the energy gap increases with the externalgate bias, as it is shown in Fig.(5a), the energy gap iscontrollable by means of the velocity ratio in large δ . Infact, for δ (cid:29) t , the band gap saturates with the gate volt-age at the value which is proportional to the interlayercoupling ( t ). In this limit, by applying the approximationof ( | E − δ |(cid:29) t ) in Eq.3, one can analytically derivethat the saturated band gap at δ (cid:29) t behaves with thevelocity ratio as the following form; E sat.gap ( η ) = 2 √ ηη + 1 t. (8)In the special case of η = 1, the band gap saturatesat E sat.gap ( η = 1) → t . As shown in the inset Fig.(5a),numerical calculations completely confirm this analyticalderivation. The momentum shift ∆ k , which measureshow much the gap is indirect, can be manipulated byusing the gate bias. Fig.(5b) represents the momentumshift from k c to k v in respect of the gate bias for sev-eral values of ξ u . This momentum shift from k c to k v increases with the gate bias. If we transform the velocityratio as η → /η , in the spectrum feature, the conductionband will be exchanged with the valence band. Further-more, based on Eq.(3), the band gap is robust againsttransformation of η → /η .In addition to the direct measurements of the spec-trum, the dependence of the energy gap on the velocityratio can be manifested in transport properties througha velocity junction. III. TRANSPORT PROPERTIES ACROSSNON-UNIFORM POTENTIAL AND VELOCITYJUNCTIONS
Let us consider a BLG sheet in which the velocity ofitinerant quasi-particles in the upper and lower layersvaries in space; representing as v u ( −→ r ) and v d ( −→ r ). Weassume that variation of velocity is smooth on the scaleof the lattice constant. In this section, we outline the ap-proach used to investigate transport properties througha barrier of velocity and potential. A. Current Density Operator
First, by using the continuity equation, we derive thecurrent density operator. The continuity equation is asthe following, ∇ .j = − ∂ t ρ (9)where ρ = Ψ † Ψ is the charge and j is the current den-sity operator. By using the Schroedinger equation, di-vergence of the current density operator is written as ∇ .j = [( H Ψ) † Ψ − Ψ † ( H Ψ)] /i (cid:126) (10)By substitution of H from Eq.1 and two componentspinor as Ψ = (cid:18) ψ u ψ d (cid:19) in the above equation, we have i (cid:126) ∇ .j = − (cid:18) ( − i (cid:126) v u ( σ. ∇ ) † + δ ) ψ u + F ψ d F ψ u + ( − i (cid:126) v d ( σ. ∇ ) − δ ) ψ d (cid:19) † (cid:18) ψ u ψ d (cid:19) + (cid:0) ψ † u ψ † d (cid:1) (cid:18) ( − i (cid:126) v u ( σ. ∇ ) † + δ ) ψ u + F ψ d F ψ u + ( − i (cid:126) v d ( σ. ∇ ) − δ ) ψ d (cid:19) (11)After simplification, it is derived that interestingly, thecurrent density operator is independent of the gate bias δ and also the hopping matrix F . ∇ .j = (cid:104) v u ∇ . ( ψ † u σ † ψ u ) + v d ∇ . ( ψ † d σψ d ) (cid:105) (12)Therefore, current density operator for a BLG sheet ispresented as, j = (cid:18) ψ u ψ d (cid:19) † (cid:18) v u σ † v d σ (cid:19) (cid:18) ψ u ψ d (cid:19) . (13)Finally, the current density in the i’th region can bewritten in the following compact form. j i = Φ † i ΣΦ i (14)where the auxiliary spinor is defined as Φ i = (cid:101) v i Ψ i andΣ = (cid:18) σ † σ (cid:19) , (cid:101) v i = (cid:18)(cid:112) v iu (cid:112) v id (cid:19) . B. Transfer Matrix Method
We assume a plane wave solution for the four-bandHamiltonian. So the wave function in each region witha constant potential is written as the following matrixproduct, Ψ( x ) = P ( x ) ∗ A , where P ( x ) and A are theplane wave and coefficient matrices, respectively. Detailof matrices P and A are accessible in appendix B and also Refs.(35,38). The local current density in terms ofmatrices P ( x ) and A in each region reads as the followingform, j i = A † i P † i (cid:101) v † i Σ (cid:101) v i P i A i (15)where the auxiliary spinor in Eq.14 has been replaced byΦ i = (cid:101) v i P i A i . The continuity equation of −→∇ . −→ j ( −→ r ) = 0leads to the boundary matching condition at interfacesof a junction. On the other word, conservation of thecurrent density results in the continuity of the auxiliaryspinor Φ i on the boundaries of the barrier junction.Φ = Φ = ⇒ (cid:101) v Ψ = (cid:101) v Ψ Referring to the schematic cartoon shown in Fig.1, weconsider a simultaneous barrier of velocity, v ( x ) = v u = v d = v F I : x < , III : x > w v u , v d II : 0 < x < w (16)and electrostatic potential. At the same time, the bar-rier can be subjected to a gate bias. V ( x ) = V u = V d = V I : x < , III : x > w (cid:26) V u = V (cid:48) + δ/ V d = V (cid:48) − δ/ II : 0 < x < w(17)By applying continuity of the auxiliary spinor on theboundaries of the barrier, one can connect the coefficientmatrix related to the last region A to the coefficientmatrix for the first region A . A = M A M = P − (0)˜ v − ˜ v P (0) P − ( w )˜ v − ˜ v P ( w ) (18)where M is the transfer matrix. We assume that theenergy range of incidence particles in the first region islimited to the range of 0 < ε < t . Consequently,the wave numbers α (1)+ and α (3)+ which are defined in theappendix.B, are real while α (1) − and α (3) − are imaginary.In this range of energy, coefficient matrices in the firstand third regions are proposed as the following form. A = (cid:0) r e g (cid:1) (cid:62) , A = (cid:0) t e d (cid:1) (cid:62) For the first region, e g is the coefficient of growing evanes-cent state and r is the coefficient of reflection. In the lastregion, t is the transmission coefficient and e d is the co-efficient of decaying evanescent state. By rearrangementof Eq. 18, the coefficient of transmission is derived as afunction of the transfer matrix elements as the following; t = [ M − M M /M ] − . (19)The transmission probability of particles through a bar-rier is defined as the ratio of out-flowing current to in-flowing current. T = J out J in (20)where J out is the out-flowing current in the last regionand J in is the in-flowing current incidence from the firstregion. By using Eq.15, the transmission probability canbe represented as the following form. T = (cid:0) t (cid:1) P † ˜ v † Σ˜ v P t (cid:0) (cid:1) P † ˜ v † Σ˜ v P (21)The conductance is calculated by using Landauer for-malism in the linear regime. Transport is coherent andis calculated at zero temperature. Conductance is pro-portional to angularly averaged transmission projectedalong the current direction. G = 2 G (cid:90) π/ T ( E, ϕ ) cos( ϕ ) dϕ (22)where G = e mvw/ (cid:126) . C. Transport across a single velocity barrier
The behavior of a beam produced by Dirac fermionswhenever hit on the barrier region, is similar to the be-havior of an optical beam passing through dielectric ma-terials. In the subsequent sections, we will show that aquantum mechanical version of well-known laws in geo-metrical optics can be also applied on the propagation ofDirac fermions in BLG.
Case i : Let us consider tunneling through a singlevelocity-induced sharp barrier. For a pure velocity bar-rier, type of quasi-particles inside and outside of the bar-rier is the same for all ranges of energy. For normalincidence θ = 0 and in absence of any gate bias, trans-mission coefficient for a velocity barrier with unity veloc-ity ratio η = 1 inside and outside of the barrier can beanalytically calculated as t = e iα w [cos( α w ) − iS sin( α w )] − (23)where S = 12 ( ω α ω α + ω α ω α ) FIG. 6: a) Transmission probability as a function of inci-dence angle θ for the case of the same velocity modulatedin both layers. The velocity ratio ξ u = ξ d = v /v F is setto 0 . , . , . , .
0. Estimated critical angle for the velocityratio 1 . . θ cr. ≈
56 and 42. b) Transmissionprobability in terms of the velocity ratio for several incidenceangles. we consider a thick velocity barrier with the width w = 100 nm for parts (a) and (b). c) Conductance as a func-tion of the width for the velocity ratio equal to 0 . , . , . E = 10 meV and α = [ ω + ω t/v ] / , α = [ ω + ω t/v ] / are thewave vectors along the x-axis direction outside and insidethe velocity barrier, respectively. Here, scaled energy ineach region is defined as ω = E/v and ω = E/v . Re-placing defined parameters in S , results in S = 1. There-fore, transmission probability is derived as the followingform, T = | t | = 1[cos ( α w ) + sin ( α w )] = 1 . (24)As a result, independence of all barrier parameters, trans-mission at the normal incidence is always perfect. Thisbehavior is similar to what we expect from the standardKlein tunneling. This transparency at the normal in-cidence will be demonstrated numerically in Fig.6. Atarbitrary incidence angle, the wave vectors along the x-axis direction in the regions I and II can be representedas the following. α = (cid:115) v ( E + tE ) − k y , α = (cid:115) v ( E + tE ) − k y (25)Suppose that the velocity outside the barrier v is setto be v F . Conservation of the energy E and the com-ponent k y of the wave vector across the barrier leads tothe following compact form for the wave vector inside thevelocity barrier. α = k (cid:114) ξ − sin θ (26)where ξ = ξ u = ξ d = v /v F is the velocity of quasi-particles inside the barrier scaled by v F . θ is the in-cidence angle of quasi-particles which hit on the barrierfrom the region I. For the range of ξ >
1, a look at Eq.26obviously demonstrates that there are some evanecsentmodes in the barrier region (in which α is imaginary) ifonly the incidence angle θ is greater than a critical anglewhich is defined as, θ cr. = arcsin(1 /ξ ) . (27)In analogous with optics, the total internal reflection(TIR) emerges when a Dirac fermion wave hits from adenser medium (region I) on a rarer medium (the bar-rier region II). This behavior is interpreted as ξ > . To demonstrate such a crit-ical angle in BLG, we plot transmission probability T ( θ )as a function of the incidence angle in Fig.6 for severalvalues of velocity. For ξ > θ > θ cr. , transmis-sion is negligible for enough thick barriers. We have alsochecked that variation of transmission around the criti-cal TIR angle is more sharp for the multiple structure ofvelocity barriers in compared with the single velocity bar-rier. Furthermore, as indicated in Fig.(6b), transmissionprobability shows a sharp change in behavior at ξ = 1.In the case of η (cid:54) = 1, the larger velocity modulated in theupper or lower layer, the smaller critical angle emerges.The critical angle just depends on ξ . So this property ismore appropriate for designing a waveguide based on theBLG substrates .As a conclusion for Eq.26, for the range of ξ <
1, thewave vector inside the barrier α is real which gives risethe propagating modes. Consequently, some resonancestates are expected to emerge. The resonance states obeythe following resonance condition, α ( θ , ξ, E ) w = nπ ,where n is the resonant order. As seen in Fig.(6a,b), the FIG. 7: a) A 3D contour-plot of transmission probability as afunction of incidence angle θ and energy ε = E − V (cid:48) for thevelocity ratio a) ξ u = ξ d = 1 . > ξ u = ξ d = 0 . < ξ > ξ <
1. The gate bias is set to δ = 40 meV . Thegate potential applying on both layers is set to V (cid:48) = 80 meV for inside the barrier and V = 0 for outside the barrier. Thebarrier width is w = 50 nm . velocity barrier is transparent against the propagation ofDirac fermionic waves at the resonance states. The reso-nance states emerges at the special values of the incidenceangle, the barrier width and those velocities belonging tothe range of ξ <
1. To distinguish the propagating fromthe evanescent modes, we study conductance as a func-tion of the barrier width in Fig.(6c) for several values ofvelocity. For ξ <
1, conductance has an oscillatory be-havior with the barrier width originating from the prop-agating modes. on the other hand, conductance dropssharply to zero for ξ >
D. Transport across velocity barrier in presence ofa gate bias
Case ii : In this case, velocity of carriers changes (still η = 1) in the barrier region where a perpendicular gatebias is simultaneously applied ( δ (cid:54) = 0). To clarify trans-port properties of the mentioned system, it is worth tolook at the wave number inside the barrier. The wavenumber along the x-axis direction is presented as FIG. 8: Conductance in terms of Fermi energy for severalvalues of the velocity ratio η . Transport gap depends on thevelocity ratio. The gate bias is equal to δ = 125 meV . Thegate potential applied on both layers is V (cid:48) = 140 meV forinside the barrier and V = 0 for outside the barrier. Thebarrier width is w = 30 nm . α = k (cid:114) µλξ − sin θ (28)where µ = ε + δ + (cid:113) ε δ − t ( δ − ε ) λ = ε + tε Following procedure of the previous section gives a criti-cal incidence angle as θ cr. ( E ) = arcsin( (cid:114) µλ ξ ). Note thatin this formula, the critical angle depends on the Fermienergy. Therefore, in the presence of electrostatic gatepotentials, this set-up is not proposed as an appropriatematerial for designing waveguide. However, a definitionof the critical angle for such system is useful to interpretthe behavior of transmission.A 3 D contour plot of transmission in terms of the in-cidence angle and energy is indicated in Fig.7 for twovelocity values: a) ξ = 1 . > ξ = 0 . <
1. For ξ = 1 . ε > ξ = 0 .
5. As a conclusion, in addition to the parameters µ and λ which are energy dependence, the velocity ξ stillcan play an important role to tune transporting modes.To manifest such a property, we study conductance asa function of energy in Fig.(7c) for ξ = 0 . .
5. Itis interesting that conductance for the velocity ξ = 0 . ξ = 1 .
5. Moreover,conductance has an oscillatory behavior with the Fermienergy if ξ <
1. However, it behaves smoothly with theenergy if ξ > Case iii : In the last case, in addition to the gate bi-asing ( δ (cid:54) = 0), we modulate velocity in layers not to beequal to each other η (cid:54) = 1. Since the band gap of the bar-rier portion depends on the velocity ratio η , we expect tomanifest this property by concentrating on the transportgap. Fig.8 represents conductance in terms of energy forseveral values of the velocity ratio η . What is novel isthat the transport gap appeared in conductance dependson the velocity ratio η . The behavior of the conductionand valence band edges with the velocity ratio is in goodagreement with those shown in Fig.(4c).Referring to Fig.(5a), dependence of the band gap onthe velocity ratio is strong when the gate bias is large. Soone can observe that the transport gap depends on thevelocity ratio. The transport gap is remarkable when athick velocity barrier is manipulated in the presence of alarge gate bias δ . Maximum band gap emerges at η = 1. IV. CONCLUSION
In the presence of a previously applied gate bias, theelectronic band structure of bilayer graphene is investi-gated when quasi-particles have different Fermi velocityin each layer. We address that the velocity engineering isone of the inevitable experimental factors which affectsthe transport gap in the broken-symmetry BLG.In absence of any electrostatic potential, only the mod-ulation of velocity in layers does not cause to open a bandgap. In other words, the chiral symmetry conserves forpurely velocity modulation δ = 0 while this symmetrywill break when a gate bias is subsequently applied onBLG. It should be noted that in the presence of a gatebias δ (cid:54) = 0, the electron-hole symmetry preserves when-ever the same velocity is modulated in both layers; η = 1.In addition, the band structure keeps its ’ Mexican hat ’shape with a direct band gap. Moreover, the band gapis independent of velocity value. The maximum value ofthe band gap occurs at η = 1. The momentum attributedto the band gap is inversely proportional to the velocity.In a generic case, non-equal velocities in two layers( η (cid:54) = 1) result in the transition of the direct-to-indirectband gap. The band gap depends on the velocity ratio η and has a peak at η = 1. Interestingly, the electron-holesymmetry fails, however the band structure still keepsits ’ Mexican hat ’ shape. The shift of momentum fromthe conduction band edge to the valence band edge isincreased with the gate bias.In the second part, we elaborate a transfer matrixmethod to calculate coherent tunneling through a veloc-ity barrier possibly subjected to a gate potential. In anal-ogous with optics, we propose a total internal reflectionangle θ cr. so that transmission becomes sharply negligiblefor the incidence angles larger than θ cr. . The transportgap which is induced by application of the gate bias inthe barrier region, depends on the velocity ratio.0 FIG. 9: The electron-hole asymmetric factor as a functionof momentum. The gate bias is equal to δ = 400 meV . V. ACKNOWLEDGEMENT
We highly acknowledge R. Asgari for his useful com-ments during improvement of this work. One of the au-thors, (H.C), thanks the institute for research in funda-mental sciences (IPM) and also the international centerfor theoretical physics (ICTP) for their hospitality andsupport during a visit in which part of this work wasdone. We should also thank M. Barbier for his commentsin the four-band tunneling.
Appendix A: The electron-hole asymmetry
To measure the electron-hole (e-h) asymmetry, we de-fine the e-h asymmetric factor as the following; ( | E c | − | E v | ) / | E c | . For the case of equal velocities modulatedin both layers η = 1, the e-h asymmetric factor is zero forthe studied Hamiltonian shown in Eq.1. However, as seenin Fig.9, this asymmetric factor increases with the mo-mentum very faster than a linear behavior . This factorreaches to the value of 2 in the special momentum. It isinteresting that by application of the transformation of η → /η , the e-h asymmetric factor behaves as r → /r .In addition to the velocity modulation, the e-h asym-metry is also originated from the inter-layer coupling ( γ )between A − A B − B . At the first orderapproximation, we have not considered such term in thedominant Hamiltonian shown in Eq.1. In fact, the mostimportant terms which affect the main feature of theband structure are γ and γ = t . Here, γ is the intra-layer hopping between A − B A − B v F in the tight-binding approximationand γ = t is the inter-layer coupling between A − B γ , behaves as 4 γ /γ . The well-established values for the hopping parameters are equal to γ ≈ . eV and γ ≈ eV . So the e-h asymmetric factor originating from γ is in order of magnitude 0 . . As a conclusion,in the presence of the previously created band gap, thee-h asymmetry arising from the velocity engineering isa dominant factor in compared with the e-h asymmetrycaused by parameter γ . Appendix B: Wave Function
The eigenfunction of four band Hamiltonian of Eq.1 isdefined with the following spinor.Ψ( x ) = P ( x ) A (B1)where coefficient matrix is written as A = (cid:0) U A U B D B D A (cid:1) (cid:62) and plane wave matrix is presented as P ( x ) = e iα + x e − iα + x e iα − x e − iα − x f ++ e iα + x f − + e − iα + x f + − e iα − x f −− e − iα − x s + e iα + x s + e − iα + x s − e iα − x s − e − iα − x g ++ s + e iα + x g − + s + e − iα + x g + − s − e iα − x g −− s − e − iα − x (B2) f ± + = v u ± α + − ik y ε − δ , f ±− = v u ± α − − ik y ε − δg ± + = v d ± α + + ik y ε + δ g ±− = v d ± α − + ik y ε + δs ± = ( ε − δ ) − v u [( α ± ) + k y ] t ( ε − δ ) , ε i = E − V i (B3)where α + and α − are the wave vectors along the currentdirection ( x ) which is defined as α ± = (cid:113) a ( ε, η, δ ) − v u k y ± (cid:112) a ( ε, η, δ ) − b ( ε, η, δ ) /v u (B4). If the gate voltage turns on, a and b defined in Eq.2are function of ε in stead of E . ∗ Electronic address: [email protected] A. H. Castro Neto, et al. , Rev. Mod. Phys. , 109 (2009). V. P. Gusynin and S. G. Sharapov, Phys. Rev. Lett. ,146801 (2005). X. Du, I. Skachko, F. Duerr, A. Luican, and E. Y. Andrei,Nature , 192 (2009); K. I. Bolotin, F. Ghahari, M. D.Shulman, H. L. Stormer, and P. Kim, Nature , 196(2009). V. V. Cheianov, V. Falko, and B. L. Altshuler, Science , 1252 (2007). J. R. Williams, L. DiCarlo, and C. M. Marcus, Science , 638 (2007). N. Stander, B. Huard, and D. Goldhaber-Gordon, Phys.Rev. Lett. , 026807 (2009) . K. S. Novoselov, A. K. Geim, S. V. Morozov, D. Jiang, M.I. Katsnelson, I. V. Grigorieva, S. V. Dubonos, and A. A. Firsov, Nature , 197 (2005); ibid:
Nat. Phys. , 177(2006). J. D. Joannopoulos, et al. , Photonic Crystal: Molding theFlow of Light, Princeton University Press, Princeton NJ,2008. A. Concha and Z. Tesanovic, Phys. Rev. B. , 033413(2010). P. E. Trevisanutto, C. Giorgetti, L. Reining, M. Ladisa, V.Olevano, Phys. Rev. Lett. , 226405 (2008). D. C. Elias, et al. , Nat. Phys. , 701 (2011). X. Du, I. Skachko, et al. ,Nat. Nano. , 491 (2008). C. H. Park, et al. , Nat. Phys. , 213 (2008). D. A. Siegel, et al. , Natl. Acad. Sci. USA , 11365(2011). C. Jang, et al. , Phys. Rev. Lett. , 146805 (2008). H. Choongyu, et al. , Nat. Scientific Reports , 590 (2012). C. Attaccalite and A. Rubio, Phys. Status. Solidi. B ,2523 (2009); M. Polini, et al. , Solid. State. Commu. ,58, (2007); G. Borghi, et al. , Solid. State. Commu. ,1117, (2009). A. Raoux, M. Polini, et al. , Phys. Rev. B , 073407(2010); G. Borghi, et al. , solid state commun. , 1117(2009). C. H. Park, S. G. Louie, Nano Lett. , 1793 (2009); M.Gibertini, et al. , Phys. Rev. B , 241406 (2009). J. H. Yuan, et al. , J. Appl. Phys. , 103706 (2011); Jian-Hui Yuann, et al. , Physica B , 4214 (2011). P. M. Krstajic, P. Vasilopoulos, J. Phys. Cond. Matt. ,135302 (2011). N. M. R. Peres, J. Phys. Cond. Matt. ,095501 (2009). T. Ohta, et al. , Science , 951 (2006). E. Castro, et al. , Phys. Rev. Lett. , 216802 (2007). J. Oostinga, et al. , Nature Mater. , 151 (2008). Y. Zhang, et al. , Nature, , 820 (2009). E. Mccann, M. Koshino, Rep. Prog. Phys. 76,056503(2013). I. Martin, et al. , Nature Phys. 7, 38-42 (2010); R. Weitz, et al. , Science 330, 812-816 (2010). P. Hofmann, Nature Mater. , 874875 (2013). J. Li, et al. , Nature Phys. 7 38 (2011). E. V. Castro, et al. , Phys. Rev. Lett. , 036802 (2011);E. Rossi, S. Das Sarma, Phys. Rev. Lett. , 155502(2011). K. S. Kim, et al. , Nature Mater. , 887892 (2013). S-M. Choi, S-H. Jhi,and Y-W. Son , Nano. Lett. , 3486(2010). B. Verberck, B. Partoens, F. M. Peeters, and B. Trauzettel,Phys. Rev. B. , 125403 (2012). M. Barbier, P. Vasilopoulos, F. M. Peeters and J. M.Pereira, Phys. Rev. B. , 155402 (2009) B. Van Duppen, F. M. Peeters, Phys. Rev. B. , 205427(2013). M. Mucha-Kruczynski, et al. , Semicond. Sci. Technol. ,033001(2010). H. cheraghchi, F. Adinehvand, J. Phys. Cond. Matt.24