Energy levels of hybrid monolayer-bilayer graphene quantum dots

PHYSICAL REVIEW B 93, 165410 (2016) Energy levels of hybrid monolayer-bilayer graphene quantum dots M. Mirzakhani,1,2,* M. Zarenia,1,† S. A. Ketabi,2,‡ D. R. da Costa,3 and F. M. Peeters1,3,§ 1 Department of Physics, University of Antwerp, Groenenborgerlaan 171, B-2020 Antwerp, Belgium 2 School of Physics, Damghan University, P. O. Box: 36716-41167, Damghan, Iran 3 Departamento de Fı́sica, Universidade Federal do Ceará, 60455-900 Fortaleza, Ceará, Brazil (Received 28 January 2016; published 8 April 2016) Often real samples of graphene consist of islands of both monolayer and bilayer graphene. Bound states in such hybrid quantum dots are investigated for (i) a circular single-layer graphene quantum dot surrounded by an infinite bilayer graphene sheet and (ii) a circular bilayer graphene quantum dot surrounded by an infinite single-layer graphene. Using the continuum model and applying zigzag boundary conditions at the single-layer–bilayer graphene interface, we obtain analytical results for the energy levels and the corresponding wave spinors. Their dependence on perpendicular magnetic and electric fields are studied for both types of quantum dots. The energy levels exhibit characteristics of interface states, and we find anticrossings and closing of the energy gap in the presence of a bias potential. DOI: 10.1103/PhysRevB.93.165410 I. INTRODUCTION Quantum dots (QDs) in both single-layer and bilayer graphene have been the subject of intensive research during the last few years, owning to their unique electronic and optical properties [1–14]. Single-layer graphene (SLG) QDs are small flakes cut out from graphene in which carrier confinement is due to the quantum size effect. The electronic and optical properties of such QDs depend on the shape and edges of the dot. For example, in the presence of zigzag edges, the energy spectrum of SLG QDs exhibits zero-energy levels, while with armchair edges the spectrum displays an energy gap [2,8,16]. Electrostatic confinement of electrons in integrable graphene QDs was also proposed in which the effect of edges is no longer important [15]. Bilayer graphene (BLG) consists of two van der Waals (vdW) coupled layers of SLG and has a gapless and paraboliclike spectrum at low energies [17]. Unlike SLG, an external electric field, realized by external gate potentials, can induce a tunable band gap in the energy spectrum of BLG. Engendering this gap using nanostructured gates led to the realization of electrostatic defined BLG QDs [5,18]. In such BLG QDs, the confinement is due to the electrostatic potentials and therefore the effect of edges are no longer important. Such QDs were recently realized by two different experimental groups [19,20]. Most graphene samples exfoliated from graphite consist of islands of one or few layers of graphene [21–23]. In these samples, the influence of junctions between different graphene regions play a significant role in their transport and electronic properties. It was shown that SLG-BLG hybrid systems exhibit unusual transport properties due to the different quantum Hall (QH) states in the SLG and BLG regions [24]. An unconventional Landau quantization was recently observed at the interface of such hybrid systems [25]. The Landau levels of an infinite SLG-BLG system were theoretically studied for both zigzag and armchair boundary conditions at the 1D interface [26]. Here, we investigate a very different geometry and demonstrate that carriers can be confined in SLG and BLG islands in a hybrid QD-like structure made of SLG-BLG junctions. This novel type of QDs can be realized by the (accidental) nanostructuring of one of the graphene layers in bilayer graphene. For convenience, we will restrict ourselves to circular QDs. This will allow us to present analytical results, which we will compare with a pure numerical approach. We propose the following two types of hybrid QDs: SLG-infinite BLG—a circular SLG QD surrounded by an infinite BLG sheet [see Fig. 1(a)] and BLG-infinite SLG—a circular BLG QD surrounded by an infinite SLG sheet [see Fig. 1(b)]. Taking the circular geometry with radius R for the QDs, we employ the continuum model, i.e., solving the Dirac-Weyl equation, and obtain analytical results for the energy levels and corresponding wave functions. We study the effect of both perpendicular electric and magnetic fields on the energy levels. For zero-magnetic field, we demonstrate that SLG-infinite BLG, in contrast to the BLG-infinite SLG QD, exhibit confined states in the presence of an external electric field. If such a circular QD is cut out of BLG, one will have both armchair and zigzag edges. However, in order to obtain analytical results for the energy levels and to observe features brought by the zigzag edges in the spectrum, we will implement the zigzag boundary condition at the SLG-BLG junction in our continuum approach. Breaking the inversion symmetry due to the interface and breaking the time reversal symmetry with a magnetic field, the two Dirac valleys K and K ′ , should be studied separately. The paper is organized as follows. In Sec. II, we consider SLG-infinite BLG QDs in the (A) absence and (B) presence of a perpendicular magnetic field. In both cases, the effect of an external electric field is studied. Section III concerns numerical results for BLG-infinite SLG QDs. We conclude the manuscript in Sec. IV. II. SLG-INFINITE BLG QUANTUM DOTS * [email protected] † [email protected] ‡ [email protected] § [email protected] 2469-9950/2016/93(16)/165410(11) A. Zero magnetic field First, we investigate the energy levels of a circular SLG QD embedded in infinite BLG [see Fig. 1(a)]. This system can be 165410-1 ©2016 American Physical Society PHYSICAL REVIEW B 93, 165410 (2016) MIRZAKHANI, ZARENIA, KETABI, DA COSTA, AND PEETERS 2R form [28] 2R (a)  τ (ρ,ϕ) = eimϕ (b) A1 B1 A2 B2 R R B B A A FIG. 1. Schematic pictures of the proposed circular SLG-BLG hybrid QDs with radius R. (a) SLG-infinite BLG QD: circular SLG dot surrounded by an infinite BLG. (b) BLG-infinite SLG QD: circular BLG dot surrounded by an infinite SLG. The upper pictures show a side view of the systems. considered as an infinite BLG sheet where a circle of radius R is cut out from its upper graphene layer. So we assume that one layer of BLG, containing A1 and B1 sublattices, seamlessly continues to the SLG with A and B sublattices, while the other graphene layer composed of A2 and B2 sublattices is sharply cut at the boundary r = R. We obtain the corresponding Hamiltonian in both SLG and BLG regions and by implementing zigzag boundary conditions to one of the graphene layers at the SLG-BLG interface, we calculate the energy levels. The dynamics of carriers in the honeycomb lattice of covalent-bond carbon atoms of single layer graphene can be described by the following Hamiltonian, which in zero magnetic field is given by [27] H = vF p · σ + U1 I, (1) where vF ≈ 106 m/s is the Fermi velocity, p = (px ,py ) is the two-dimensional momentum operator, σ denotes the Pauli matrices, and U1 is the potential applied to SLG. We assume that the carriers are confined in a circular area of radius R, with zigzag boundary. In polar coordinates and dimensionless units, the Hamiltonian (1) reduces to the form   u1 π+ H= , (2) π− u1 with the momentum operator π± = −ie±iτ ϕ   ∂ iτ ∂ , ± ∂ρ ρ ∂ϕ (3) where the dimensionless variables are ρ = r/R and u1 = U1 R/vF . r and ϕ are the radial and azimuthal coordinates of the cylindrical coordinate system, respectively. The two valleys are labeled by the quantum number τ , which is τ = +1 for the K valley and τ = −1 for the K ′ valley. The Schrödinger equation becomes H(ρ,ϕ) = ε(ρ,ϕ), (4) where the carrier energy E, is written in dimensionless units as ε = ER/vF . The two-component wave function has the   φAτ (ρ) , ie−iτ ϕ φBτ (ρ) (5) where m = 0,±1,±2, . . . denotes the angular momentum label. The components φA and φB correspond to different sublattices A and B, respectively. Solving Eq. (4), the radial dependence of the spinor components is described by   d mτ − 1 τ − φB (ρ) = (ε − u1 )φAτ (ρ), dρ ρ   d τm τ φA (ρ) = −(ε − u1 )φBτ (ρ). + (6) dρ ρ Decoupling the above equations, we arrive at the Bessel differential equation for φAτ : d 2 φAτ (ρ) dφAτ (ρ) + ρ + [(ε − u1 )2 ρ 2 − m2 ]φAτ (ρ) = 0, dρ 2 ρ (7) with the solution ρ2 φAτ (ρ) = C τ Jm (aρ), (8) τ where a = ε − u1 and C is the normalization constant. The second component of the wave function can be obtained from Eq. (6) as φBτ (ρ) = −τ C τ Jm−τ (aρ). Thus the wave function becomes   C τ Jm (aρ) τ imϕ .  (ρ,ϕ) = e ieiτ ϕ τ C τ Jm+τ (aρ) (9) (10) The BLG region can be described in terms of four sublattices, labeled A1, B1, for the lower layer and A2, B2, for the upper layer [see Fig. 1(a)]. The A1 and B2 sites are coupled via a nearest-neighbor interlayer hopping term t ≈ 0.4 eV. The BLG Hamiltonian in the vicinity of the K point, is given by (in dimensionless units) [18,29] HK = H0K + ( u/2)σz , (11) with u0 ⎜π− K H0 = ⎝ ′ t 0 ⎛ π+ u0 0 0 t′ 0 u0 π+ ⎞ 0 0⎟ , π− ⎠ u0 (12) where t ′ = tR/vF , u0 = (u1 + u2 )/2, u = u1 − u2 , and u1,2 = U1,2 R/vF , with U1 and U2 the potentials at the two layers. The operator σz is defined as   I 0 , (13) σz = 0 −I where I is the 2 × 2 identity matrix. The Hamiltonian at the K ′ point is obtained by interchanging π+ and π− in Eq. (12). The eigenstates of Hamiltonian (11) are four-component spinors [30] ⎛ K ⎞ φA1 (ρ)eimϕ ⎜iφ K (ρ)ei(m−1)ϕ ⎟ ⎜ ⎟ K (14) (ρ,ϕ) = ⎜ B1K ⎟, ⎝ φB2 (ρ)eimϕ ⎠ 165410-2 K iφA2 (ρ)ei(m+1)ϕ PHYSICAL REVIEW B 93, 165410 (2016) ENERGY LEVELS OF HYBRID MONOLAYER-BILAYER . . . where m is the angular momentum label. Solving the Schrödinger equation, the radial dependence of the spinor components are described by   d m−1 K K K φB1 (ρ) = (α − δ)φA1 − (ρ) − t ′ φB2 , dρ ρ   m K d K φ (ρ) = −(α − δ)φB1 + (ρ), dρ ρ A1   m+1 K d K K φA2 (ρ) = (α + δ)φB2 + (ρ) − t ′ φA1 (ρ), dρ ρ   m K d K φ (ρ) = −(α + δ)φA2 − (ρ), (15) dρ ρ B2 where τ = ±1, distinguishes the boundary conditions for the two valleys. The above conditions lead to a system of equations from which we obtain the eigenvalues. For the K point, with the help of the wave functions (10), (14), (18), and (18), we arrive at ⎞ ⎛ K⎞ ⎛ C −Jm (a) Km (γ+ ) Km (γ− ) K M K ⎝C1 ⎠ = ⎝Jm−1 (a) b+ Km−1 (γ+ ) b− Km−1 (γ− )⎠ 0 c+ Km (γ+ ) c− Km (γ− ) C2K ⎛ K⎞ C K ⎝ C (22) × 1 ⎠ = 0, C2K where α = ε − u0 and δ = (u1 − u2 )/2. These equations can K : be decoupled and we obtain for φA1  2  d 1 d m2 K K + (ρ) = γ±2 φA1 (ρ), (16) − 2 φA1 2 dρ ρ dρ ρ where b± = γ± /(α − δ) and c± = ((α − δ)2 + γ±2 )/(α − δ)t ′ . The corresponding calculations for the K ′ point leads to ⎛ ′⎞ ⎛ ⎞ CK −Jm (a) c+ Km (γ+ ) c− Km (γ− ) ′ ′⎜ ⎟ M K ⎝C1K ⎠ = ⎝−Jm+1 (a) d+ Km+1 (γ+ ) d− Km+1 (γ− )⎠ ′ 0 Km (γ+ ) Km (γ− ) C2K ⎛ ′⎞ CK ⎜C K ′ ⎟ × ⎝ 1 ⎠ = 0, (23) K′ C2 where the potential-dependent eigenvalues are γ± = {−(α 2 + δ 2 ) ± [(α 2 − δ 2 )t ′2 + 4α 2 δ 2 ]1/2 }1/2 . (17) The differential equation (16) is the known modified Bessel equation. Here we choose the modified Bessel function of the second kind Km (γ± ), as the appropriate solutions vanishing at r → ∞. Thus we have K φA1 (ρ) = C1K Km (γ+ ρ) + C2K Km (γ− ρ). (18) Using Eqs. (15), we obtain the other spinor components: 1 C K γ+ Km−1 (γ+ ρ) + C2K γ− Km−1 (γ− ρ) , α−δ 1 1 K φB2 C K ((α − δ)2 + γ+2 )Km (γ+ ρ) (ρ) = (α − δ)t ′ 1  + C2K ((α − δ)2 + γ−2 Km (γ− ρ) , K φB1 (ρ) = K (ρ) = φA2 (α 2 1 C K γ+ ((α − δ)2 + γ+2 )Km+1 (γ+ ρ) − δ 2 )t ′ 1 + C2K γ− ((α − δ)2 + γ−2 )Km+1 (γ− ρ) , (19) where CjK (j = 1,2) are the normalization constants. The wave function for the K ′ valley can be written as ⎞ ⎛ K′ φA1 (ρ)eimϕ ⎜iφ K ′ (ρ)ei(m+1)ϕ ⎟ ⎟ ⎜ K′ (20) (ρ,ϕ) = ⎜ B1K ′ ⎟. ⎝ φB2 (ρ)eimϕ ⎠ ′ K iφA2 (ρ)ei(m−1)ϕ Solving the Schrödinger equation (4) for the K ′ valley and comparing with the differential equations (15), we find K′ K′ K′ K′ K K K K (φA1 ,φB1 ,φB2 ,φA2 ) = (φB2 ,φA2 ,φA1 ,φB1 ). Now, we apply zigzag boundary conditions [26] at the SLGBLG interface. These conditions yield Aτ (ρ,ϕ) = τ A1 (ρ,ϕ)|ρ=1 , τ B1 (ρ,ϕ)|ρ=1 , 0= τ B2 (ρ,ϕ)|ρ=1 , Bτ (ρ,ϕ) = (21) where d± = c± γ± /(α + δ). The nonzero eigenenergies are the ′ solutions of det |M K | = 0 and det |M K | = 0. The zero-energy states can be investigated separately by solving Eqs. (6) and (15) in the case of ǫ = 0 and U1 = U2 = 0. Applying the boundary conditions (21), one finds zero-energy states at the K (K ′ ) valley only for m  0 (m  0). Returning to the nonzero eigenenergies for √ the unbiased case U1 = U2 = 0, we find γ± = (−ε2 ± ε2 t ′2 )1/2 which is pure imaginary when |ε| > t ′ . In the interval |ε| < t ′ , γ+ is real while γ− is pure imaginary. Requiring det |M K | = ′ det |M K | = 0, we find a continuum energy band with no discrete levels and thus no confined states in the QD. In the presence of bias, γ+ = γ−∗ when |ε| < u1 . For this interval we can select the real (imaginary) part of the modified Bessel functions Km (γ+ ) (Km (γ− )) as our solutions. Figure 2 shows the energy levels with the angular momenta m = 0,±1,±2,±3 as a function of the dot radius R, for the biased potential U1 = −U2 = 0.1 eV. An energy gap appears between the conduction and valence bands. Notice that the band gap is  2 given by [17] g = |△U |t/ △U + t 2 (see the solid black horizontal lines in Fig. 2), coming from the Mexican-hat shaped low-energy dispersion in pristine BLG. For U ≪ t, the band gap is g ≈ U . For both valleys, the number of energy levels increases, as the dot radius increases and the band gap decreases to zero. In both cases, each set of energy levels are approximately equally spaced for fixed m. The energy spectrum corresponding to the K valley [Fig. 2(a)] exhibits the symmetry EK (m) = EK (−m), which does not hold for the K ′ valley levels. This is different from SLG [28] and BLG QD [5] flakes in which the K and K ′ energy levels are degenerate at zero magnetic field. This difference is due to the zigzag boundary condition applied to the SLG-BLG interface which removes the layer symmetry and thus the valley symmetry in BLG. The energy spectrum in Fig. 2(b) shows groups of energy 165410-3 PHYSICAL REVIEW B 93, 165410 (2016) MIRZAKHANI, ZARENIA, KETABI, DA COSTA, AND PEETERS carriers are confined in both SLG QD and at the SLG-BLG interface. The set of energy levels with E ∼ + g /2 are interface states that are predominantly confined along the zigzag edge of the second layer with low density in the dot area of the first layer [see Fig. 3(c)]. (a) E (eV) m=0 m=±1 m=±2 m=±3 1 2 B. Nonzero magnetic field In the presence of a perpendicular magnetic field B = B êz , one needs to replace the canonical momentum p by the gaugeinvariant kinetic momentum p + eA(r) in the Hamiltonians (2) and (12). A(r) = (0,Br/2,0) is the vector potential in the symmetric gauge. In the case of √ B = 0, it is more convenient to use the magnetic length lB = /eB as the unit of length, facilitating the interpretation of our numerical results.√This results in the following dimensionless quantities, ξ = r/ 2lB , ε = ElB /vF , and u1 = U1 lB /vF . In these units, π± in Eqs. (2) and (12) takes the form   iτ ∂ 1 ∂ ± ∓ τξ , (24) π± = −ie±iτ ϕ √ ξ ∂ϕ 2 ∂ξ 0.10 3 (b) 4 E (eV) 0.05 0.00 0.05 0.10 0 10 20 30 40 50 60 70 R (nm) FIG. 2. Energy levels of SLG-infinite BLG QD for m = 0,±1,±2,±3 as a function of the dot radius R, in the presence of bias U1 = −U2 = 0.1 eV at the (a) K and (b) K ′ valleys. The energy levels of the K valley (a) satisfy EK (m) = EK (−m). The solid (dashed) curves are for m > 0 (m < 0). levels with E ∼ + g /2 and m > 0. These levels correspond to states that are mainly confined at the interface [see Fig. 3(c)]. Figures 3(a)–3(d) show the electron densities in each layer corresponding to the points labeled by (1)–(4) in Figs. 2(a) and 2(b). The energy levels at the K valley [see Figs. 3(a) and 3(b) for the levels (1) and (2)] show that the confinement is mostly in the SLG QD, while for the energy levels corresponding to the K ′ valley [see Figs. 3(c) and 3(d)] the (a) (b) and the radial Dirac-Weyl equation for the spinor components of the wave function (5) becomes   1 ∂ (τ m − 1) − − τ ξ φBτ (ξ ) = (ε − u1 )φAτ (ξ ), √ ξ 2 ∂ξ   τm 1 ∂ + + τ ξ φAτ (ξ ) = −(ε − u1 )φBτ (ξ ). (25) √ ξ 2 ∂ξ After decoupling the above equations, we obtain     ξ2 τ 1 d2 1 d m2 φ (ξ ) − + − 2 + (m − τ ) + 2 dξ 2 ξ dξ ξ 2 A = (ε − u1 )2 φAτ (ξ ). (26) 2 Using the ansatz φAτ (ξ ) = ξ |m| e−ξ /2 f (ξ 2 ), Eq. (26) yields the confluent hypergeometric ordinary differential equation ξ̃ df (ξ̃ ) d 2 f (ξ̃ ) + (|m| + 1 − ξ̃ ) 2 d ξ̃ d ξ̃  1 − |m| + m − τ + 1 − (ε − u1 )2 f (ξ̃ ) = 0, (27) 2 where ξ̃ → ξ 2 . The solutions are the confluent hypergeometric functions of the first kind 1 F˜1 (a,b,ξ 2 ) with   a = 21 |m| + m − τ + 1 − (ε − u1 )2 , b = |m| + 1. (28) Then Aτ becomes 0.054 Aτ (ξ,ϕ) = C τ eimϕ ξ |m| e−ξ 0.026 (c) (d) 0.082 0.054 FIG. 3. Probability density corresponding to the states indicated by (1), (2), (3), and (4) in the energy spectrum of Fig. 2 for R = 30 nm. Blue solid curves refer to layer 1, and red dashed curves denote layer 2. 2 /2 ˜ 1 F1 (a,b,ξ 2 ), (29) τ where C is the normalization constant. The second component of the wave function is extracted from the second equation of Eq. (25): ei(m−τ )ϕ 2 Bτ (ξ,ϕ) = −iC τ √ ξ |m| e−ξ /2 2(ε − u1 )  × 2ξ a 1 F˜1 (a + 1,b + 1,ξ 2 ) 165410-4 +  |m| + τ m + (τ − 1)ξ ξ   2 ˜ F (a,b,ξ ) . (30) 1 1 PHYSICAL REVIEW B 93, 165410 (2016) ENERGY LEVELS OF HYBRID MONOLAYER-BILAYER . . . Using Eqs. (4), (11), (12), (14), and (24) in the new dimensionless units the radial dependence of the spinor components in BLG are described by   1 d (m − 1) K − − ξ φB1 (ξ ) √ ξ 2 dξ K K = (α − δ)φA1 (ξ ) − t ′ φB2 (ξ ),   1 d m K K + + ξ φA1 (ξ ) = −(α − δ)φB1 (ξ ), √ dξ ξ 2   1 d (m + 1) K (ξ ) + + ξ φA2 √ ξ 2 dξ K K = (α + δ)φB2 (ξ ) − t ′ φA1 (ξ ),   1 d m K K − − ξ φB2 (ξ ) = −(α + δ)φA2 (ξ ), √ ξ 2 dξ (31) where α = ε − u0 , δ = (u1 − u2 )/2, u0 = (u1 + u2 )/2, u1,2 = U1,2 lB /vF , ε = ElB /vF , and t ′ = tlB /vF . U1 and U2 are the potentials on the two different graphene layers. K These equations can be decoupled to obtain for φA1     1 d2 ξ2 K m2 1 d − + m + φ (ξ ) − + 2 dξ 2 ξ dξ ξ2 2 A1 = K γ± (ε)φA1 (ξ ), with the matrix elements   m11 = − 1 F̃1 a; b; ξR2 ,   m12 = U a+ ,b,ξR2 ,   m13 = U a− ,b,ξR2 ,    1 m21 = 2ξR a1 F̃1 a + 1; b + 1; ξR2 ε − u1    (|m| + m) F̃1 a; b; ξR2 , + ξR 1    1 2ξR a+ U a+ + 1; b + 1; ξR2 m22 = α−δ   (|m| + m)  − U a+ ; b; ξR2 , ξR    1 2ξR a− U a− + 1; b + 1; ξR2 m23 = α−δ   (|m| + m)  − U a− ; b; ξR2 , ξR m31 = 0,   m32 = 4ξR2 a+ (a+ + 1)U a+ + 2,b + 2,ξR2   + 4 ξR2 − |m| − 1 a+ U a+ + 1,b + 1,ξR2   + 2[(α − δ)2 − |m| − m]U a+ ,b,ξR2 ,   m33 = 4ξR2 a− (a− + 1)U a− + 2,b + 2,ξR2   + 4 ξR2 − |m| − 1 a− U a− + 1,b + 1,ξR2   (37) + 2[(α − δ)2 − |m| − m]U a− ,b,ξR2 . (32) where γ± = α 2 + δ 2 ± [(α 2 − δ 2 )t ′2 + (1 − 2αδ)2 ]1/2 . (33) 2 K (ξ ) = ξ |m| e−ξ /2 f (ξ 2 ), similar for the Using the ansatz φA1 SLG region, we arrive at the confluent hypergeometric ordinary differential equation with solutions 1 F˜1 (a± ,b,ξ 2 ) and U (a± ,b,ξ 2 ), where a± = 21 (|m| + m + 1 − γ± ), b = |m| + 1. (34) Now for the BLG region, we need to take the confluent hypergeometric functions of the second kind U (a,b,x) which K decays exponentially for r → ∞. Then φA1 becomes K φA1 (ξ ) = ξ |m| e−ξ 2 /2 C1K U (a+ ,b,ξ 2 ) + C2K U (a− ,b,ξ 2 ) , (35) One can similarly obtain the corresponding matrix for the K ′ K′ . The nonzero-energy levels are obtained from the valley MZZ K |=0 condition det |MZZ In order to find the solutions of the zero-energy states, one can solve Eqs. (25) and(31) in the case of ǫ = 0 (U1 = U2 = 0). Applying the boundary conditions (21) at the interface shows that zero-energy states exist for all angular momenta m, at the K valley, and only for m < 0 at the K ′ valley. 1. Unbiased system where CjK (j = 1,2) are the normalization constants. The other spinor components of the wave function can be obtained using K and employing the properties of Eq. (31) by inserting φA1 the confluent hypergeometric function. The wave function K′ K′ K′ K′ at the K ′ point can be obtained from (φA1 ,φB1 ,φB2 ,φA2 )= K K K K (φB2 ,φA2 ,φA1 ,φB1 ). Applying the boundary conditions (21) for the K point at √ the interface ξ = ξR = R/( 2lB ), we arrive at ⎛ K⎞ C K ⎝C K ⎠ MZZ = 0, 1 C2K where K MZZ ⎛ m11 = ⎝ m21 m31 m12 m22 m32 ⎞ m13 m23 ⎠, m33 (36) Figure 4 shows the energy levels as a function of magnetic field at the K (solid curves) and K ′ (dashed curves) valleys when U1 = U2 = 0. The energy levels are shown for the angular momenta, m = −1 (blue), m = 0 (green), and m = 1 (red) with QD radius R = 30 nm. As mentioned before, the SLG-infinite BLG QD in the absence of unbiased magnetic field displays a continuum energy band which is also apparent in Fig. 4. In this case, there are many degenerate zero-energy states at zero magnetic field corresponding to all angular momenta for both K and K ′ valleys. Increasing the magnetic field, zero-energy degenerate levels for each m is lifted due to the breaking of the time reversal symmetry. Breaking of the inversion symmetry due to the interface removes the degeneracy of the K and K ′ valleys. The spectrum shows anticrossings, which is due to the influence of the SLG-BLG interface. At high magnetic fields, the energy levels merge into 165410-5 PHYSICAL REVIEW B 93, 165410 (2016) MIRZAKHANI, ZARENIA, KETABI, DA COSTA, AND PEETERS 3.0 (a) R 2.0 B 1.5 m=-1 c m=0 a 0 5 10 K K’ b 1.3 10 0.0 3.0 1.4 m=1 0.5 15 15 20 B (T) 25 2.5 30 35 2.0 1.5 1.0 0.5 K 0.0 0 √ the Landau levels (LLs) of SLG, (ε = ± 2n, n is the LLs’ number [26]). In a strong magnetic field, the carriers become localized at the center of the SLG dot and will not be influenced by effects due to the dot interface. It should be mentioned that electron-hole symmetry exists for both valleys, because all K K′ [Eq. (36)] and MZZ are even function matrix elements in MZZ of the dimensionless energy ε, when U1 = U2 = 0. An enlargement around a particular anticrossing point is given in the inset of Fig. 4. The wave functions as well as the probability densities for the points labeled by (a), (b), and (c) in Fig. 4 are shown in Fig. 5. The upper panels show the wave functions of layer 1 which are continuous in both SLG and BLG. The wave functions of layer 2 are plotted in the middle panels. Point (a) corresponds to confinement in both BLG regions and the SLG-BLG interface, while for the energy state (c) the electrons are confined inside the SLG QD. Right K '(a) B1 K' A (b) (c) K' A1 K' B K' A2 K' B2 | K' 2 Layer 2 | K' 2 Layer1 | r (nm) r (nm) (c) (b) U=0 FIG. 4. Energy spectrum of unbiased SLG-infinite BLG QD (i.e., U1 = U2 = 0) as a function of magnetic field for m = −1 (blue), m = 0, (green) and m = 1 (red) with the dot radius R = 30 nm at the K (solid curves) and K ′ (dashed curves) valleys. The inset is an enlargement of the black square box for a particular anticrossing of energy levels in the K ′ valley. Wave Function (a. u.) R B2 A2 A 1.0 | B1 A1 1.5 (ћvF/lB) (ћvF/lB) 2.5 r (nm) FIG. 5. The wave functions corresponding to the points (a), (b), and (c) of the inset in Fig. 4. The upper panel shows the wave functions for layer 1, the middle panel displays the wave functions for layer 2, and the lowest panel shows the density of the bound states in the two layers for the K ′ valley. 15 B (T) 30 K’ 0 15 B (T) 30 FIG. 6. (a) Schematic pictures of the terminated systems, SLG QD (left) and BLG antidot (right). Lower panel displays the energy spectrum of SLG-infinite BLG QD (black solid curves) and the terminated systems, SLG QD (blue dashed curves) and BLG antidot (red dashed curves) as a function of magnetic field for the two valleys (b) K and (c) K ′ . The dot radius is R = 30 nm and m = 1. at the anticrossing [i.e., point (b)] one finds confinement in both SLG and BLG. Plateaulike features appear in the energy spectrum, that can be understood when comparing the energy levels of terminated SLG QD and BLG antidot. Consider two terminated systems, SLG QD and BLG antidot with zigzag edges as shown in Fig. 6(a). The boundary condition is Aτ = 0 for SLG QD and τ τ A2 = 0 for BLG antidot. Lower panels of Figs. 6(b) B1 = and 6(c) show the energy spectrum of SLG-infinite BLG QD (black solid curves) and those of the terminated systems (dashed curves) at the two valleys K and K ′ for m = 1. The spectrum of SLG-infinite BLG QD resembles that of the terminated systems with an energy gap opened at every crossing point. Such energy gaps at the crossing points can be interpreted as due to the hybridization between SLG QD and BLG states. The energy levels of the lowest bound states for the unbiased system are shown in Fig. 7 as a function of the dot radius for B = 10 T , and m = −1, 0, 1 at the K (solid curves) and K ′ (dashed curves) valleys. For R → 0, the energy levels correspond to the LLs of unbiased bilayer graphene given by [31]   t4 1 t2 2 ′ + (2n + 1)E0 ± + (2n′ + 1)t 2 E02 + E04 , ε=± E0 2 4 (38) where E0 = vF / lB , n′ = n + (|m| + m)/2, and n = 0,1,2, . . ., with m > 0 and m < 0 states being degenerate. 165410-6 PHYSICAL REVIEW B 93, 165410 (2016) ENERGY LEVELS OF HYBRID MONOLAYER-BILAYER . . . 3 (a) (ћvF/lB) (ћvF/lB) 2 SLG LLs + u 0 1 m=1 m=0 K K’ m=-1 m=0 m=-1 m=1 1 2 U=0 3 0 3 R (nm) FIG. 7. Energy states of unbiased SLG-infinite BLG QD (U1 = U2 = 0) as a function of the dot radius with B = 10 T, for m = −1,0,1 at the valleys K (solid curves) and K ′ (dashed curves). 2. Biased system The dependence of the spectrum on the magnetic field, at the two valleys K and K ′ for a biased system with U1 = −U2 = 0.1 eV, are shown respectively in Figs. 8(a) and 8(b). The results are presented for m = −1 (blue), m = 0 (green), and m = 1 (red). The energy levels show a band gap between the conduction and valence bands at the K and K ′ points. For small magnetic field (B → 0), this band gap exhibits a divergence when expressed in the units of vF / lB . Applying a gate potential breaks the degeneracy of the lowest-energy states, and the spectrum becomes strongly dependent on m. However, as the magnetic field increases, the magnetic confinement becomes important, as seen by the lifting of the degeneracy of the states, and the energy levels approach the LLs of SLG √ (see black dashed curves, i.e., ε = ± 2n + u1 ). Furthermore, both K and K ′ spectra show electron-hole asymmetry because of the breaking of inversion symmetry due to the presence of the external gate potentials. The behavior of the electron states in both valleys is qualitatively similar, but the hole states display different behavior which for the K valley are nearly degenerate. Results for the spectrum of localized states as a function of R are shown in Fig. 9 for the biased case, U1 = −U2 = 0.1 eV. The energy levels are plotted for angular momentum labels, m = −1 (blue), m = 0 (green), and m = 1 (red) with B = 10 T at the valleys (a) K and (b) K ′ . When R → 0, for both valleys, energy levels coincide with the LLs of biased BLG as it should be. The LLs of biased BLG is determined by the 20 30 40 20 B (T) 30 40 (b) 2 (ћvF/lB) Increasing R, the LLs of BLG split for different valleys (breaking of inversion symmetry) as well as for different angular momenta (breaking of time reversal symmetry) and approach the LLs of SLG. The energy levels corresponding to the K valley in Fig. 7, demonstrate the merging of the nth BLG LLs at R → 0 to the nth SLG LLs at larger R. For the K ′ energy levels, the nth LLs of BLG approach the (n + 1)th LLs of SLG. Similar analysis is also available for the terminated systems, illustrated in Fig. 6, exhibiting the appearance of plateau feature in the spectrum. 10 m=-1 m=0 1 SLG LLs + u 0 m=1 1 2 3 0 10 FIG. 8. Energy spectrum of biased SLG-infinite BLG QD as a function of magnetic field for m = −1,0,1 at the (a) K and (b) K ′ valleys with the dot radius R = 30 nm. The black dashed curves depict the SLG LLs +u1 . The applied bias is U1 = −U2 = 0.1 eV. equation [31] [(α + δ)2 − 2(n′ + 1)][(α − δ)2 − 2n′ ] − (α 2 − δ 2 )t ′2 = 0. (39) For small R, the lowest energy levels become nearly degenerate, forming two energy bands around ε = ±u1 = ±1.23. As R increases, degeneracy of the energy levels is lifted and finally they connect to the LLs of SLG subject to the external potential u1 . Figure 10 shows the spectrum for terminated systems (dashed curves) and biased SLG-infinite BLG QD (black solid curves) as a function of R at the valleys K and K ′ for m = 1. The resemblance of the two different spectra are also evident for the biased case. III. BLG-INFINITE SLG QUANTUM DOTS A. Zero magnetic field In this section we choose the inverse of the previous QD system in which a BLG QD is surrounded by an infinite SLG [see Fig. 1(b)]. This system can be considered as an infinite BLG sheet in which a circle of radius R from its upper layer is left and the other part is removed. Similar to the previous section the Hamiltonian is solved for both parts of the system and then we choose the appropriate wave functions to satisfy the extreme conditions when r → 0 for bilayer dot and r → ∞ for SLG. Here, the modified Bessel function of the first kind Im , 165410-7 PHYSICAL REVIEW B 93, 165410 (2016) MIRZAKHANI, ZARENIA, KETABI, DA COSTA, AND PEETERS (ћvF/lB) (ћvF/lB) (a) m=-1 m=0 m=1 2.5 2.0 SLG LLs + u m=-1 m=1 m=0 4 BLG LLs 1.5 K K’ U=0 1.0 0.5 0.0 0 0.1 0.2 B (T) (b) 3 FIG. 11. Energy spectrum of unbiased BLG-infinite SLG QD (i.e., U1 = U2 = 0) as a function of magnetic field for R = 30 nm and m = −1,0,1 at the two valleys, K (solid curves) and K ′ (dashed curves). Black dot-dashed curves are the LLs of unbiased infinite BLG. The inset shows an enlargement of the energy levels for low magnetic fields. (ћvF/lB) 2 1 m=-1 0 m=0 m=1 1 SLG LLs + u in the presence of bias. Thus there are no bound states in BLG-infinite SLG QDs for B = 0. 2 3 0 10 20 30 40 R (nm) 50 60 70 B. Nonzero magnetic field FIG. 9. Energy spectrum of biased SLG-infinite BLG QD (U1 = −U2 = 0.1 eV) as a function of dot radius with B = 10 T and for m = −1,0,1 at the valleys K (a) and K ′ (b). Black dashed lines are the LLs of SLG +u1 . and the Bessel function of the second kind Ym , are the solutions for regions of the BLG QD and infinite SLG, respectively. As discussed before (see Sec. II A), there is no unique linear combination of real or imaginary parts of Im (γ± ) and Ym (γ± ) from which unique discrete energies can be obtained even 4 (a) (b) 3 (ћvF/lB) 2 1 K 0 K’ 1 2 3 0 35 R (nm) 70 0 35 R (nm) 70 FIG. 10. Energy spectrum of SLG-infinite BLG QD (black solid curves) and the terminated systems, SLG QD (blue dashed curves) and BLG antidot (red dashed curves) as a function of dot radius R, for the two valleys (a) K and (b) K ′ . The magnetic field is B = 10 T and m = 1. In the presence of a magnetic field, the calculations are similar to those presented in Sec. II B. To avoid repetition, we limit ourselves here to the numerical results. Similar analysis for the case of zero-energy levels as in previous section, one finds that zero-energy states exist only for m  0 at the K valley, and for all momenta at the K ′ valley. 1. Unbiased system The dependence of the spectrum on magnetic field, for R = 30 nm and m = −1 (blue), m = 0 (green), and m = 1 (red) at the two valleys K (solid curves) and K ′ (dashed curves) is shown in Fig. 11. The inset shows an enlargement of the states for small magnetic fields. As we see, for B = 0, the energy levels coincide with the LLs of SLG. The energy levels are nearly degenerate for m > 0 and m < 0 at very low magnetic fields (0–0.1 T). As the magnetic field increases, the states merge to form the LLs of an unbiased infinite BLG sheet [black dashed curves, see Eq. (38)], which indicates that the carriers become strongly localized at the center of the dot. The energy levels in Fig. 11 show that the LLs of SLG connect to LLs of BLG with the same Landau level indices for the K valley. However, for the K ′ valley, the SLG and BLG LLs connect by nSL → nBL − 1. The spectrum of terminated systems [in this case, BLG QD (red dashed) and SLG antidot (blue dashed)], and the BLG-infinite SGL QD (black solid curves) are shown in Fig. 12 at the two valleys K and K ′ for m = −1. The oscillatory feature in the energy levels can be understood in relation to the terminated systems as explained in Fig. 6. In Fig. 13, we show the dependence on R of the low energy spectrum, at the two valleys K (solid curves) and K ′ (dashed curves) for m = −1 (blue), m = 0 (green), and m = 1 (red) with B = 10 T. The R → 0 limit corresponds to the case of 165410-8 PHYSICAL REVIEW B 93, 165410 (2016) ENERGY LEVELS OF HYBRID MONOLAYER-BILAYER . . . 3 (b) (a) (a) (ћvF/lB) (ћvF/lB) 2 1 m=-1 m=0 m=1 0 Biased BLG LLs 1 2 K B (T) 3 0 K’ SLG sheet, and the spectrum agrees with the LLs of SLG with the states being degenerate for all angular momenta and for both valleys. When R increases, angular momentum as well as valley-degenerate levels split, showing flat plateau features for certain ranges of dot radius, and eventually merge into the LLs of unbiased BLG. The way of the LLs of SLG and BLG are connected is similar to those discussed for Fig. 11. 2. Biased system (ћvF/lB) Figure 14 displays the energy levels for the biased case with U1 = −U2 = 0.1 eV and R = 30 nm as a function of magnetic field for (a) K and (b) K ′ valleys. As in the unbiased case, the spectrum is plotted for m = −1 (blue), m = 0 (green), and m = 1 (red). The energy levels show a band gap between the conduction and valence bands. For very small values of B, the spectrum becomes degenerate and forms a continuum band. However, as the magnetic field increases, the degeneracy of the levels is lifted for each angular momentum m=1 K’ K’ m=-1 20 30 40 (b) B (T) R (nm) FIG. 13. Energy levels of unbiased BLG-infinite SLG QD (U1 = U2 = 0) as a function of dot radius R, for m = −1,0,1 at the two valleys, K (solid curves) and K ′ (dashed curves). The magnetic field is B = 10 T. (ћvF/lB) FIG. 12. Energy levels of BLG-infinite SLG QD (black solid curves) and the terminated systems, BLG QD (blue dashed curves) and SLG antidot (red dashed curves) as a function of magnetic field at the two valleys K (a) and K ′ (b) for m = −1. The dot radius is R = 30 nm. m=0 10 m=-1 m=0 Biased BLG LLs m=1 B (T) FIG. 14. Energy spectrum of biased BLG-infinite SLG QD (i.e., U1 = −U2 = 0.1 eV) as a function of the magnetic field for R = 30 nm and m = −1,0,1 at the two valleys, K (a) and K ′ (b). Black dashed curves are the LLs of a biased infinite BLG. (with the large shift for low lying states) and eventually approaching the LLs of biased BLG (blacked dashed curves). For each value of m, the hole energy levels show anticrossings when the energy levels approach the LLs of biased BLG. There are some other features, e.g., electron-hole asymmetry, degeneracy lifting between the states of the different valleys, and different behavior of the hole states in the K valley, similar to SLG-infinite QDs. The spectrum as a function of R for the biased case with U1 = −U2 = 0.1 eV, is plotted in Fig. 15 for the valleys (a) K and (b) K ′ . The other parameters are the same as for the unbiased case shown in Fig. 13. For R = 0, the spectrum corresponds to the LLs of the SLG sheet, being degenerate for all m. With increasing dot radius, the degeneracy of the levels for different m is lifted, and the levels connect to different LLs of BLG. The results show that the low-energy LLs of SLG converge to the low-energy LLs of BLG, as the dot radius increases, and they form two bands around these energies. It is also seen that the hole states, in the case of the K valley are approximately degenerate. In graphene QDs, it is possible to define a scaling factor for the maximum of the magnetic field by setting the dot size R equal to the cyclotron radius at the Fermi energy, i.e., R = lB2 kF [32]. Having E = vF kF and E = (vF kF )2 /t for the low energy dispersion of monolayer and bilayer graphene, respectively, one can obtain kF and consequently the scaling factor with respect to the magnetic field and size 165410-9 PHYSICAL REVIEW B 93, 165410 (2016) MIRZAKHANI, ZARENIA, KETABI, DA COSTA, AND PEETERS (ћvF/lB) (a) m=1 m=0 m=-1 (ћvF/lB) (b) m=0 m=-1 m=1 R (nm) FIG. 15. Energy spectrum of biased BLG-infinite SLG QD as a function of dot radius, for m = −1,0,1 at the two valleys, K (a) and K ′ (b). The magnetic field is B = 10 T and U1 = −U2 = 0.1 eV. The black dashed lines show the SLG LLs +u1 . of the confinement area, i.e., √RSLG = E/(evF B) for SLGinfinite BLG QD and RBLG = Et/(evF B) for BLG-infinite SLG QD. IV. CONCLUSION Using the continuum model, i.e., solving the Dirac-Weyl equation, we obtained analytical results for the energy levels and corresponding wave functions for two new types of QDs in hybrid SLG-BLG systems: SLG-infinite BLG QDs and BLG-infinite SLG QDs. We implemented the zigzag boundary condition at the SLG-BLG junction in order to observe features brought by the zigzag edges in the spectrum. Both Dirac valleys were investigated separately because of the breaking of inversion symmetry due to the presence of the SLG-BLG interface. [1] A. V. Rozhkov, G. Giavaras, Y. P. Bliokh, V. Freilikher, and F. Nori, Phys. Rep. 503, 77 (2011). [2] Z. Z. Zhang, K. Chang, and F. M. Peeters, Phys. Rev. B 77, 235411 (2008). [3] D. Bischoff, A. Varlet, P. Simonet, M. Eich, H. C. Overweg, T. Ihn, and K. Ensslin, Appl. Phys. Rev. 2, 031301 (2015). [4] B. Trauzettel, D. V. Bulaev, D. Loss, and G. Burkard, Nat. Phys. 3, 192 (2007). No bound states are observed in both types of QDs when B = 0 and when no bias potential is applied to the graphene layers. In biased SLG-infinite BLG QDs with B = 0, bound states and the opening of an energy gap between the electron and hole energy levels are found, and this gap closes when increasing the dot radius. The energy spectrum corresponding to the K valley exhibits the symmetry EK (m) = EK (−m), which no longer holds for the K ′ valley levels. In the presence of a perpendicular magnetic field without bias, degeneracy of energy levels is lifted. Increasing the magnetic field, the spectrum exhibits anticrossings and the levels merge into the LLs of SLG which indicates that the carriers become strongly localized at the center of the SLG QD. The anticrossings are due to the interplay of confined states in both SLG and BLG region and at the SLG-BLG interface. In the presence of an external electric field, our results showed that the degeneracy of the lowest-energy states are lifted at low magnetic fields and the opening of an energy gap is found. Furthermore, both the K and K ′ spectrum show electron-hole asymmetry because of the breaking of the inversion symmetry by the external potential. For BLG-infinite SLG QDs, we found that there are no bound states with or without external electric field when B = 0. Energy spectrum of unbiased BLG-infinite SLG QDs for B = 0 showed that the degeneracy of the states is lifted, and for high magnetic fields, the energy levels merge into the LLs of unbiased BLG. The biased case showed an opening of energy gap. Dependence of spectrum on dot radius in the SLG-infinite BLG QDs (BLG-infinite SLG QDs) shows that in the R → 0 limit, the energy levels correspond to the unbiased LLs of BLG (SLG). Increasing R, the LLs of BLG (SLG) split for different angular momenta and approach the LLs of SLG (BLG). For the case of biased SLG-infinite BLG QDs with B = 0, the energy gap for R → 0 closes with increasing R. The energy levels of the proposed new QDs can be investigated experimentally by STM, which measures the local density of states (LDOS). From the LDOS, one can obtain information about the position of the bound states and the localized electron/hole distribution. ACKNOWLEDGMENT This work was supported by the Fonds Wetenschappelijk Onderzoek (FWO)-CNPq project between Flanders and Brazil and the Brazilian Science Without Borders program. [5] J. M. Pereira Jr, F. M. Peeters, P. Vasilopoulos, R. N. Costa Filho, and G. A. Farias, Phys. Rev. B 79, 195403 (2009). [6] P. Recher and B. Trauzettel, Nanotechnol. 21, 302001 (2010). [7] A. V. Rozhkov and F. Nori, Phys. Rev. B 81, 155401 (2010). [8] M. Zarenia, A. Chaves, G. A. Farias, and F. M. Peeters, Phys. Rev. B 84, 245403 (2011). [9] A. D. Güçlü, P. Potasz, and P. Hawrylak, Phys. Rev. B 84, 035425 (2011). 165410-10 PHYSICAL REVIEW B 93, 165410 (2016) ENERGY LEVELS OF HYBRID MONOLAYER-BILAYER . . . [10] D. P. Zebrowski, E. Wach, and B. Szafran, Phys. Rev. B 88, 165405 (2013). [11] M. Zarenia, B. Partoens, T. Chakraborty, and F. M. Peeters, Phys. Rev. B 88, 245432 (2013). [12] A. Müller, B. Kaestner, F. Hohls, T. Weimann, K. Pierz, and H. W. Schumacher, J. Appl. Phys. 115, 233710 (2014). [13] D. R. da Costa, M. Zarenia, A. Chaves, G. A. Farias, and F. M. Peeters, Carbon 78, 392 (2014). [14] D. Bischoff, M. Eich, O. Zilberberg, C. Rössler, T. Ihn, and K. Ensslin, Nano Lett. 15, 6003 (2015). [15] J. H. Bardarson, M. Titov, and P. W. Brouwer, Phys. Rev. Lett. 102, 226803 (2009). [16] M. Ezawa, Phys. Rev. B 76, 245415 (2007). [17] E. McCann, D. S. L. Abergel, and V. I. Falko, Eur. Phys. J. Spec. Top. 148, 91 (2007). [18] J. M. Pereira, P. Vasilopoulos, and F. M. Peeters, Nano Lett. 7, 946 (2007). [19] M. T. Allen, J. Martin, and A. Yacoby, Nat. Commun. 3, 934 (2012). [20] A. M. Goossens, S. C. M. Driessen, T. A. Baart, K. Watanabe, T. Taniguchi, and L. M. K. Vandersypen, Nano Lett. 12, 4656 (2012). [21] Y. Kobayashi, K. I. Fukui, T. Enoki, K. Kusakabe, and Y. Kaburagi, Phys. Rev. B 71, 193406 (2005). [22] Y. Niimi, T. Matsui, H. Kambara, K. Tagami, M. Tsukada, and H. Fukuyama, Phys. Rev. B 73, 085421 (2006). [23] G. M. Rutter, J. N. Crain, N. P. Guisinger, P. N. First, and J. A. Stroscio, J. Vac. Sci. Technology. A 26, 938 (2008). [24] J. Tian, Y. Jiang, I. Childres, H. Cao, J. Hu, and Y. P. Chen, Phys. Rev. B 88, 125410 (2013). [25] W. Yan, S.-Y. Li, L.-J. Yin, J.-B. Qiao, J.-C. Nie, and L. He, arXiv:1502.00785. [26] M. Koshino, T. Nakanishi, and T. Ando, Phys. Rev. B 82, 205436 (2010). [27] P. Recher, B. Trauzettel, A. Rycerz, Y. M. Blanter, C. W. J. Beenakker, and A. F. Morpurgo, Phys. Rev. B 76, 235404 (2007). [28] M. Grüjić, M. Zarenia, A. Chaves, M. Tadić, G. A. Farias, and F. M. Peeters, Phys. Rev. B 84, 205441 (2011). [29] E. McCann and V. I. Falko, Phys. Rev. Lett. 96, 086805 (2006). [30] M. Zarenia, F. M. Peeters, G. A. Farias, and J. M. Pereira, Jr, Nano Lett. 9, 4088 (2009). [31] J. M. Pereira, Jr, F. M. Peeters, and P. Vasilopoulos, Phys. Rev. B 76, 115419 (2007). [32] D. K. Ferry, L. Huang, R. Yang, Y. C. Lai, and R. Akis, J. Phys.: Conf. Ser. 220, 012015 (2010). 165410-11