Dirac Spectrum in Piecewise Constant One-Dimensional Potentials

We study the electronic states of graphene in piecewise constant potentials using the continuum Dirac equation appropriate at low energies, and a transfer matrix method. For superlattice potentials, we identify patterns of induced Dirac points which are present throughout the band structure, and verify for the special case of a particle-hole symmetric potential their presence at zero energy. We also consider the cases of a single trench and a p-n junction embedded in neutral graphene, which are shown to support confined states. An analysis of conductance across these structures demonstrates that these confined states create quantum interference effects which evidence their presence.

We study the electronic states of graphene in piecewise constant potentials using the continuum Dirac equation appropriate at low energies, and a transfer matrix method. For superlattice potentials, we identify patterns of induced Dirac points which are present throughout the band structure, and verify for the special case of a particle-hole symmetric potential their presence at zero energy. We also consider the cases of a single trench and a p-n junction embedded in neutral graphene, which are shown to support confined states. An analysis of conductance across these structures demonstrates that these confined states create quantum interference effects which evidence their presence.

I. INTRODUCTION
Graphene, a single layer of carbon atoms laid out in a honeycomb lattice, is one of the most interesting electronic systems to become available the last few years 1,2 . It differs from conventional two dimensional electron gas (2DEG) systems in that the low energy physics is governed by a Dirac Hamiltonian rather than the more common form for semiconductors, characterized by an effective mass and a band gap. Because of its unusual properties, the possibility of exploiting it for nanoscale applications has become an area of intense investigation.
Nanophysics in graphene may be investigated in systems with very small physical size scales, such as quantum dots 3 or quantum rings [4][5][6][7][8][9] . More commonly, manipulation at the nanoscale is accomplished by the application of electric fields via gate geometries, subjecting the system to potentials varying on a very short length scale. Much recent work has focused on p-n junctions [10][11][12][13][14] and p-n-p junctions 15,16 , whose behavior in graphene is distinct from that of conventional 2DEG's because of the absence of a gap in the spectrum. Very recently, studies of graphene in superlattice potentials have demonstrated the possibility of "band structure engineering" of graphene 17,18 . In particular, for one-dimensional superlattices 19,20 , effective potentials may be designed by which the number of Dirac points at the Fermi energy can be artificially manipulated 21,22,[26][27][28][29] .
In this work, we study graphene systems in the presence of one-dimensional potentials that are piecewise constant. Our approach models graphene in terms of the Dirac equation, which is valid for potentials that vary slowly on the lattice scale, and for which intervalley scattering may be ignored. The allowed energy states for such systems may be found by a transfer matrix method described below. For a one-dimensional superlattice, such problems are reminiscent of dynamical systems with periodically varying parameters. We find that the resulting band structure generically supports a collection of Dirac points at a variety of energies. Dirac points generally appear at energies satisfying ε ≡ Ed/ v F = nπ (E = electron energy, d = superlattice period, v F = Fermi velocity) for n an integer, at q = 0, where q is the momentum perpendicular to the superlattice axis. For particle-hole symmetric potentials, several Dirac points may appear at ε = 0 for non-vanishing values of q [21][22][23][26][27][28][29] . Generically, other Dirac points are present for ε, q = 0.
We also study geometries with a single p or n region, or with a p-n junction region, embedded in a neutral graphene background. We refer to these as a single trench geometry and an embedded p-n junction, respectively. The trench geometry supports localized solutions, at energies for which the electron states are evanescent in the neutral graphene for a given k y = 0. These states disperse as a function of k y such that a finite number of them crosses zero energy. The embedded p-n junction supports localized solutions as well, which may be interpreted in terms of a pair of overlaid trench spectra. If the two trenches are particle-hole conjugates, the individual trench spectra cross zero energy at the same value of k y ; however, because of tunneling between the trenches these become avoided crossings rather than Dirac points. At large k y the gaps at the avoided crossings may be quite small.
Finally we consider transport across these systems within the Landauer formalism, for which the transfer matrix method is easily adapted to yield the conductance of the system. We demonstrate that the conductance has distinctive signatures indicative of the confined states which may be contained at non-vanishing k y for a given embedded structure.
This article is organized as follows. In Section II we express the problem of a Dirac particle in a periodic potential in terms of a path-ordered product, and find the band structure of the simplest case of alternating positive and negative potential regions. Section III develops the transfer matrix method for handling more general piecewise constant potentials. In Section IV we apply the method to compute the spectrum of localized solutions in trench and embedded p-n junction geometries. Section V applies the results of the transfer matrix calculations to compute the conductance of these system. Finally, we conclude with a summary in Section VI.

II. PERIODIC POTENTIAL
In the vicinity of a Dirac point at wavevector K, we write the wavevector k = K + ∆k, and the two component wavefunction as ψ ≡ e iK·x u v . Assuming a scalar potential V (x) which is independent of the coordinate y, the wavefunction φ = u v may be chosen an eigenstate of the Hamiltonian H = e −iK·x H e iK·x : where q = ∆k y . The eigenvalue equation for φ may then be recast in matrix form as ∂ and The solution is where P denotes the path ordering operator, which places smaller values of x to the right. We presume V (x) = V (x+d) is periodic, in which case the Bloch condition becomes where θ x is the Bloch phase. This leads to the condition where Note that M † = −σ x M σ x , where σ x is the Pauli matrix, and that Tr M = 0. Since det(Λ) = 1, the eigenvalues of Λ are given by ζ = 1 2 T ± √ T 2 − 4, where T = Tr (Λ). The energy bands, where the Bloch condition can be satisfied, correspond to the regime |T | ≤ 2. For |T | > 2 the wavefunctions grow exponentially and are not normalizable. The mathematical structure here is common to a variety of problems involving linear equations with periodic coefficients 24 , such as the classical mechanics setting of a pendulum with periodically varying length.
For a continuous periodic function V (x) the path ordered exponential is not analytically tractable. However, many of the essential features are reproduced by considering a piecewise constant V (x) of the form We We then have with α j = d j κ 2 j − q 2 . The resulting trace is (10) The Bloch condition is satisfied when |T | ≤ 2. It is convenient to define the dimensionless quantities where is the average potential. We also define the dimensionless wavevector θ y = q (d 1 + d 2 ). When V is constant, one finds T = 2 cos ε 2 − θ 2 y , and the Bloch condition is satisfied provided |E − V | > v F |q|. This, of course, is the condition that the energy E lies within the Dirac cone.
Next, consider V 1 = V 2 , but with q = 0. The trace is then T = 2 cos ε, and there are allowed states at all energies. The Bloch condition is satisfied marginally when ε = nπ, corresponding to the touching of two energy bands which disperse with q, i.e. to the vanishing of a band gap as q → 0. Indeed this result holds for arbi- , M (x ′ ) = 0 and the path ordering is superfluous. Now consider the case d 1 = d 2 . The spectrum then exhibits a particle-hole symmetry relative to E = V . When E = V , i.e. ε = 0, we have Clearly T ≥ 2, with the marginal condition T = 2 pertaining when where n is an integer, as well as a solution at θ y = 0. Since q is real, there is a maximum allowed value for n: where the brackets indicate the greatest integer function. Thus, two new Dirac points open up each time ω increases by 2π, and there are 2n max +1 Dirac points, including the one at θ y = 0. Similar results have previously been obtained by other methods 21,22 . In Figs. 1 and 2 we plot the contours of the function T (ε, θ y ) for ω = 13 2 π for d 1 = d 2 and d 1 = 2 3 d 2 , respectively. Note the particle hole symmetry present in the former case, which is broken when d 1 = d 2 . Note also the band touchings at ε = nπ for θ y = 0 in both cases. In Fig. 3 we again plot the band structure as in Fig. 1, but showing the high energy structure as well. At energy scales |ε| ≫ 1, the spectrum may be understood as the intersection of the two Dirac cones |ε + ω| > θ y and |ε − ω| > θ y .
As a limiting case of the model, consider i.e. a Dirac comb with a uniform compensating background. The trace in this case is given by where α = (ε − Ω 0 ) 2 − θ 2 y . This limit is distinguished by the fact that all the Dirac points of T (ε, θ y ) lie along the line θ y = 0, as shown in Fig. 4.

III. TRANSFER MATRICES AND LANDAUER CONDUCTANCE THROUGH POTENTIAL SLABS
Consider a piecewise constant potential V (x) as depicted in Fig. 5 in which V = V j over a slab of width d j , with j = 1, . . . , N . We presume the slabs are bounded by pristine regions of graphene with V = 0, and we label these outer regions 0 and N + 1. Wavevectors k j for each of the slabs may be defined by the relations where each k j ≥ 0 is nonnegative. Here, σ j = ±1 is a binary variable in each slab which indicates whether the energy E is above (n-type) or below (p-type) the Dirac point, which occurs at E = V j . We take σ 0 = σ N +1 and In each region j, the solution to the Dirac equation, is either propagating or evanescent. Due to translation invariance in the y-direction, we are free to fix the wavevector q along this axis; all solutions include a factor e iqy . The local Fermi wavevector k j in each region is related to the local density by n j = k 2 j /4π, and the sign of the carrier charge is −σ j . If k j > q, the wavefunction in slab j propagates as a plane wave along the x-direction as well. The solution is of the form Matching the y-dependence gives us the condition To solve for the wavefunction everywhere we must match the value of each component of ψ(x) at the boundaries between slabs. Accordingly, from Eq. 20, we define the matrix which we shall invoke presently. Free propagation across a slab is described by the transfer matrix If k j < q, the solution in slab j has an evanescent (pure exponential) behavior along the x-axis: where ϕ j ∈ − π 2 , π 2 . Momentum conservation along y requires For the evanescent slabs, we define and free propagation across an empty slab is described by The equations guaranteeing continuity of the wavefunction at the slab boundaries are Here, the superscripts R and L refer to the right and left boundaries of the slab. They are related, within a given slab, by the free propagation transfer matrix, We can now compute the transfer matrix M for a general configuration of slabs, which is defined by the rela- In terms of the previously defined quantities, the total transfer matrix is The recurring matrices Q j = M j K j M −1 j may be parameterized by θ j and χ j = kd cos θ j for a propagating region, and by ϕ j and ω j = kd cot ϕ j for an evanescent region, as detailed in the Appendix. As is well-known, the transfer matrix relates data on the left to data on the right of a given region. Incoming and outgoing flux amplitudes are related by the S-matrix, Appealing to the unitarity of S, we then have As discussed in more detail below, the dimensionless 2terminal Landauer conductance, in units of e 2 /h per spin channel, is given by

IV. LOCALIZED SOLUTIONS OF EMBEDDED GEOMETRIES
We first consider situations where the outer regions contain evanescent waves, so the wavefunction is localized within the slabs. In this case, we must have B 0 = 0 and A N +1 = 0, if cot ϕ > 0, which requires M 11 = 0; for the case cot ϕ < 0, we must have A 0 = 0 and B N +1 = 0, leading to M 22 = 0. This condition defines the one-dimensional dispersion E(q), where q = k 0 csc ϕ 0 = k j sin θ j .
(37) The condition M 22 = 0 is obtained by reversing the sign of χ 1 or ω 1 . We now define the dimensionless potential ν ≡ V d/ v F , energy ε ≡ σ 0 k 0 d = ν + σ 1 k 1 d, and wavevector θ y = qd. The eigenvalue condition is then F (ε, ν, θ y ) = 0, where and This form for F holds for either sign of q and is extended to the regime (ε − ν) 2 < θ 2 y via analytic continuation. If we search for eigenvalues at ε = 0, we obtain the condition a graphical analysis of which shows a new root emerging each time ud passes through nπ. Thus, there are ν/π solutions with ε = 0 which lie at finite values of θ y . In Fig. 6, for example, we have ν = 14 and thus there are [14/π] = 4 bands which cross ε = 0 at finite θ y .

B. Embedded junction
Consider next the case of an embedded junction with a double step, When V 1 and V 2 are of opposite sign, this corresponds to a p-n junction embedded in neutral graphene. We Applying the transfer matrix formalism, we once again arrive at an eigenvalue condition F (ε, ν 1 , ν 2 , x 1 , θ y ) = 0, this time with where and where we continue to express energies and wavevectors as dimensionless quantities. As with the single trench case, the expression for F is analytically continued to the regimes where (ε − ν j ) 2 < θ 2 y . Generally speaking, when two semi-infinite undoped regions are separated by a sequence of doped slabs, in addition to the Dirac cones in the bulk one obtains localized subbands with energies lying within the gap between the upper and lower Dirac cones. The wavefunctions for these subband states are localized within the slabs. In Fig. 6, we plot the Dirac cones and subband dispersions for the case of a single strip of width d with potential  Figs. 7 and 8, we investigate the case of a p-n junction embedded in neutral graphene, where both p and n regions have the same width d, and again we choose V = ±14 v F /d for purposes of illustration. Apparently E = 0 always lies within a gap between subbands. Note that the parameters of the p-n junction are such that the spectum is particle-hole symmetric.
The subband structure of a general p-n junction can be understood in terms of the simple diagram in Fig.  9. In the figure, we sketch the three Dirac cones for the undoped bulk region (black lines), the p-type slab (blue lines, V 1 d 1 / v F = 12), and the n-type slab (red lines, V 2 d 2 / v F = −7). For the purpose of generality, we consider a system which has no particle-hole symmetry. Outside the black hatching, there are no bulk states. In regions I and III, there are localized states living along the p-trench. Similarly, in regions I and II, there are localized states living along the n-trench. In region I, there can be resonant tunneling between the p and n trenches, and the localized states live in both trenches. In region IV, which lies outside all three Dirac cones, there can be no states. In Fig. 10, we plot the subband structure for this p-n junction, with equal thicknesses of the p and n regions, superimposing on those curves the results for the individual single barrier problems. A comparison with the sketch in Fig. 9 reveals the basic physics of the energy spectrum. One may see the confined states associated with each of the two trenches, one set dispersing up and the other down. These confined states admix where they are degenerate for the individual trenches, leading to anticrossings. In the case where the two trenches are equal and opposite (e.g., Figs. 7 and 8), this leads to a set of anticrossings around zero energy, rather than the extra Dirac points one finds in the superlattice case.

V. TWO-PROBE LANDAUER CONDUCTANCE
In this section, we consider transmission through a slab geometry for the more general case where the graphene leads are doped. Consider the geometry of Fig. 5. Define a vertical (ŷ-directed) surface Σ somewhere to the right of the N th slab. The current across that surface will be a sum of three contributions. First, there will be carriers injected from the left side, with chemical potential µ L , moving at angle θ ∈ − 1 2 π, 1 2 π , which make it past the slabs with probability T (E, θ). These are all right-movers. Second, there will be carriers from the right lead, with chemical potential µ R , moving at an angle θ ∈ 1 2 π, 3 2 π . Because of this range of θ, all these are left-movers. Finally, with probability R ′ (E, θ), each of these left-moving carriers reflects off the barrier region, scattering into a state with θ ′ = π − θ. Putting this all together, where L Σ is the length of the surface and N (E) is the density of states, and g = 4 accounts for spin and valley degeneracy. J(E, θ) is given by We taken Σ =x. Unitarity of the S-matrix yields T (E, θ) + R(E, π + θ) = 1. Expanding in the difference µ L − µ R , we obtain at T = 0 where V = (µ L − µ R )/e is the voltage drop. The average transmission coefficient is defined as The carrier density at zero temperature is n = gk 2 F /4π, so the Fermi energy is E F = v F k F = v F 4πn/g. Thus, Defining the effective channel number, the two-terminal Landauer conductance becomes Expressed in a convenient set of units, the resistance is wheren is the dimensionless carrier density (in the leads) in units of 10 12 cm −2 andL y is the length of the barrier region in units of microns.
For a single trench of width d 1 and carrier density n 1 embedded between graphene leads of carrier density n 0 ,
Here we have assumed k 1 < q so the trench contains evanescent solutions, but by a simple analytic continuation our result also holds for the k 1 > q case.
In Fig. 11 we plot the average transmission coefficient T (E F ) as a function of the density n 1 for a single trench with d 1 = 100 nm, for several different values of the lead density n 0 . There are three features worthy of note. First, there is a nonzero transmission even when n 1 = 0, due to the presence of evanescent solutions within the trench 25 . Second, there is (obviously) perfect transmission when n 1 = n 0 . Third, the transmission falls off rapidly as n 1 → 0, and is smooth for 0 < n 1 < n 0 , which is a regime where the carriers in both regions are of the same sign, yet total internal reflection can occur at the interface, with the critical angle given by sin θ c 0 = n 1 /n 0 1/2 . Finally, when n 1 > n 0 or when the carriers in the leads and the trench are of opposite sign, one observes quantum oscillations in the transmission, due to constructive and destructive interference effects from scattering at the walls of the trench.
It is interesting to compare Fig. 11 with results in Ref. 11, which reports the measured conductance through an electrostatically defined trench potential in a wide graphene ribbon as a function of the gate voltage defining the trench. As in our results, the experiments indicate a rapidly falling conductance when n 1 has the same sign as n 0 , as n 1 approaches zero, and after passing through zero one observes oscillations due to quantum interference. A prominent difference is the absence in the experimental results of a rising background around which these oscillations occur. In our model, this appears because the density of states in the trench increases as |n 1 | increases. It is likely that the difference between our theoretical results and those observed in Ref. 11 are due to the relatively smooth potential induced by the remote gates in the experiment. In particular this can yield relatively large "depletion" regions between the trench and the graphene leads which may dominate the resistance of the device, and overwhelm the increase in conductance related to the trench density of states. In principle, one could model the device as a series of piecewise constant potentials and use the transfer matrix method we have developed to obtain a more quantitatively accurate model.
In Fig. 12 we plot T (E F ) as a function of n 1 for an embedded p-n junction with d 1 = d 2 = 50 nm and n 2 = −n 1 . For such a symmetric junction the average transmission is symmetric under n 1 → −n 1 . The steplike features present in the transmission here represent scattering effects due to the localized solutions and the avoided crossings in their dispersions, as discussed in the previous section. Thus the bound states support an (in principle) measurable signature in conductance through the system.

VI. SUMMARY
In this paper we have studied the spectra of and transport properties through graphene in the presence of several piecewise constant potentials. We demonstrated that these may be investigated in a straightforward manner using a transfer matrix approach. For the case of superlattice potentials, we demonstrated that exact solutions of the Dirac equation generically contain Dirac points centered at zero energy for particle-hole symmetric potentials 21,22 , and that these are special cases of Dirac points that appear throughout in the band structure of superlattice potentials. We also considered the cases of a single trench and an embedded p-n junction in neutral graphene, both of which may support confined states at large enough values of k y , the momentum along the translationally invariant direction. An analysis of conductance through these structures shows quantum interference effects which may be interpreted as signatures of these bound states, in principle allowing these unique features of graphene nanostructures to be seen experimentally.
Note added: While this manuscript was in the final stages of preparation, a manuscript by Barbier, Vasilopoulos, and Peeters appeared 23 which also considers the effect of a superlattice on monolayer graphene and the emergence of additional Dirac points. We are grateful to F. M. Peeters for subsequently alerting us to the substantial previous work which he and his collaborators have contributed on the subject of Dirac subband spectra and transport in the presence of one-dimensional potentials [26][27][28] . We are also grateful to Y. P. Bliokh for drawing our attention to Ref. 29, which also addresses subband formation in the presence of a superlattice structure.