Two-dimensional Dirac particles in a Pöschl-Teller waveguide

We obtain exact solutions to the two-dimensional (2D) Dirac equation for the one-dimensional Pöschl-Teller potential which contains an asymmetry term. The eigenfunctions are expressed in terms of Heun confluent functions, while the eigenvalues are determined via the solutions of a simple transcendental equation. For the symmetric case, the eigenfunctions of the supercritical states are expressed as spheroidal wave functions, and approximate analytical expressions are obtained for the corresponding eigenvalues. A universal condition for any square integrable symmetric potential is obtained for the minimum strength of the potential required to hold a bound state of zero energy. Applications for smooth electron waveguides in 2D Dirac-Weyl systems are discussed.

The Pöschl-Teller potential 1 plays an important role in many fields of physics 2 from modeling diatomic molecules and quantum many-body systems [3][4][5] , to applications in astrophysics 6,7 , optical waveguides 8 and quantum wells 9, 10 , through to Bose-Einstein and Fermionic condensates 11,12 , and supersymmetric quantum mechanics 13 . For the one-dimensional Schrödinger equation, the hyperbolic symmetric form can be solved in terms of associated Legendre polynomials and the eigenvalues are known explicitly 1,14 . We consider an analogous relativistic problem, that of a two-dimensional Dirac particle, confined by a one-dimensional Pöschl-Teller potential. Several solutions have been obtained for the Dirac equation with central Pöschl-Teller potentials [15][16][17][18][19][20][21] , and the hyperbolic-secant potential is also known to admit analytic solutions for both the one and two-body one-dimensional Dirac problems [22][23][24] . Modified Pöschl-Teller potential potentials have also been employed in numerical simulations of potential barriers in bilayer graphene 25 .
With the recent explosion of research in Dirac materials 26 there has been a renewed interest in quasi-relativistic phenomena considered in condensed matter systems of different dimensionalities. This is due to the fact that the same equations which govern Dirac fermions in relativity, map directly to the equations of motion describing the quasi-particles in systems such as graphene 27 , carbon nanotubes 28 , topological insulators [29][30][31] transition metal dichalcogenides 32 and 3D Weyl semimetals 33 . Massless Dirac particles are notoriously difficult to confine; however, it has been demonstrated that certain types of one-dimensional electrostatic waveguides in graphene, possess zero energy-modes which are truly confined within the waveguide 22,23 , and that the number of these zero energy-modes is equal to the number of supercritical states (i.e. bound states whose energy, E = −M, where M is the particle's effective mass). Transmission resonances and supercriticality of Dirac particles through one-dimensional potentials have been studied extensively [34][35][36][37][38][39][40][41][42][43][44][45] . The majority of studies model top-gated carbon-based nanostructures using abrupt potentials [46][47][48][49][50][51][52][53][54][55][56] . However, experimental potential profiles vary smoothly over many lattice constants, with even the smallest of gate generated potentials being far from square 57 . There is also some controversy concerning waveguides which are defined by smoothly decaying, square integrable functions, which decay at large distances as 1/x n , where n > 1. Numerical experiments imply that such waveguides always contain a zero-energy mode 58 , whereas analysis based upon relativistic Levinson theorem, says there is a minimum potential strength required to observe a zero-energy mode 45,[59][60][61][62] . Our result supports the latter, and we demonstrate this through a simple analysis of supercritical states of zero energy.
In 2D Dirac materials, the low-energy spectrum of the charge carriers can be described by a Dirac Hamiltonian 63 of the form , σ x,y,z are the Pauli spin matrices. v plays the role of the speed of light and k z is proportional to the particle in-plane effective mass. In graphene, the charge carriers are massless, k z = 0, and the dispersion is linear, v = v F ≈ 10 6 m/s is the Fermi velocity and s K has the value of ±1 for the K and K′ valley respectively 64 . For narrow-gap nanotubes and certain types of graphene nanoribbons 65, 66 the operator k y can by replaced by the number k y = E g /(2ħv F ), which in the absence of applied field is fixed by geometry, where E g is the value of the bandgap. For nanotubes, E g can be controlled by applying a magnetic field along the nanotube axis 24,[67][68][69][70] . When k z is finite, Eq. (1) can be used as a simple model for silicene 71 or Weyl semimetal 26,72 . It has also been proposed that when graphene is subjected to a periodic potential on the lattice scale, for example, graphene on top of a lattice-matched hexagonal boron nitride 73 the Dirac fermions acquire mass with E g = 2ħv F k z being of the order of 53 meV.
In what follows we shall consider a particle described by the Hamiltonian (1) subjected to a one-dimensional potential U(x), which varies on a scale much larger than the lattice constant of the corresponding Dirac material, therefore allowing us to neglect inter-valley scattering for the case of graphene. We shall also set s K = 1, and note that the other valley's wave function can be obtained by exchanging k y → −k y . The Hamiltonian acts on the two-component Dirac wavefunction ik y A B y to yield the coupled first-order differential equations and L is a constant. V = UL/ħv F and the charge carrier energy, ε, have been scaled such that E = εL/ħv F . The charge carriers propagate along the y-direction with wave vector k y = Δ/L, which is measured relative to the Dirac point, Ψ A (x) and Ψ B (x) are the wavefunctions associated with the A and B sublattices of graphene and finally M = k z L represents an effective mass in dimensionless units.
In what follows we consider the relativistic quasi-one-dimensional Pöschl-Teller potential problem which can be applied to describe e.g. graphene waveguides. We obtain the exact energy eigenfunctions for this potential and formulate a method for calculating the eigenvalues of the bound states. We then analyze the energy-spectrum of the symmetric Pöschl-Teller potential and obtain expressions for the eigenvalues of the supercritical states. By analyzing the zero-energy supercritical states we obtain a universal threshold condition for the minimum potential strength required for a potential to possess a zero-energy mode, for any square-integrable potential. We also show that the eigenfunctions in the non-relativistic limit restore the one-dimensional Schrödinger equation solutions. Finally, we analyze the eigenvalue spectrum for the modified Pöschl-Teller potential which includes an asymmetry term.

Relativistic one-dimensional Pöschl-Teller problem
In this section we consider the potential which is a linear combination of the symmetric Pöschl-Teller potential with an additional term which enables the introduction of asymmetry 74 . This potential belongs to the class of quantum models, which are quasi-exactly solvable 23,[75][76][77][78][79][80][81] , where only some of the eigenfunctions and eigenvalues are found explicitly. The depth of the well is given by -(a − b) 2 /4a, and the potential width is characterized by the parameter L, which was introduced after Eq. (4). For the case of b = 0, the potential transforms into the symmetric Pöschl-Teller potential, while if a = 0, the potential is a smooth potential step, which can be used to model a p-n junction 82,83 . The symmetric and asymmetric forms of the potential are plotted in Fig. 1. (3,4) to be reduced to a single second-order differential equation in where s = −(−1) j and j = 1, 2 correspond to the spinor components Ψ 1 and Ψ 2 respectively. By making the natural change of variable = +  Z x [1 tanh( )]/2 and using the wave function of the form allows Eq. (6) to be reduced to the Heun confluent equation 84 in variable Z: 2 , s α,β,γ = ±1 and the parameters α, β, γ, δ s and η s are: This same method of reducing a system of coupled first-order differential equations to the Heun confluent equation has been exploited to solve various generalisations of the quantum Rabi model 85 , and notably the quasi-exact solutions of the Pöschl-Teller family potentials and Rabi systems are closely related 86 . In some instances, the resulting Heun confluent functions can be terminated as a finite polynomial 23,87 , allowing particular eigenvalues to be obtained exactly, providing the parameters obey special relations. When this is not the case, the energy spectrum can be obtained fully via the Wronskian method 23,85,88,89 , which is the method we shall utilize. Equation (7) has regular singularities at Z = 0 and 1, and an irregular one at Z = ∞ which is outside the domain of  x . The solutions to Eq. (7) are given by , therefore, as expected from the general theory of second order differential equations, the full solution is given by a linear combination of just two linearly independent functions where β is not an integer. It should be noted that when a = 0 and α = δ = 0, the Heun confluent functions appearing in Eq. (9) reduces to a Gauss hypergeometric functions and for the case of a massless particle Eq. (9) reduce down to the solutions obtained in ref. 82. If, however, β is an integer then H(α, −β, γ, δ, η s , Z) is divergent and B i has to be set to zero, and the second linearly independent solution can be constructed from a series expansion about the point 1 − Z. The solutions to the Heun confluent equation thus far have been given as power series expansions about the point Z = 0. However, these power series rapidly diverge as Z approaches the second singularity; therefore, at Z = 1 we must construct solutions as power series expansions in variable 1 − Z: The constants C i and D i are found by matching the two power series expansions and their derivatives at Z 0 where 0 < Z 0 < 1. For bound states, we require that ∆ > E 2 2 and ∆ > − E b ( ) 2 2 . These conditions ensure that the bound states are inside the effective bandgap (which accounts for the motion along the y-axis). As x → −∞, Z → 0 and as x → ∞, Z → 1, therefore we may write the asymptotic expressions of Ψ j as Therefore, for bound states, B i (A i ) is zero for s β = 1 (s β = −1) and D i (C i ) is zero for s γ = 1 (s γ = −1). Clearly the choice of s β and s γ is arbitrary, therefore, from hereon in we set both to 1 unless otherwise stated. In this instance, the energy eigenvalues are found from the condition:

Symmetric Pöschl-Teller potential solutions.
In general, relating Ψ 1 to Ψ 2 is non-trivial, since neither a known expression exists which connects Heun confluent functions about two different singular points for arbitrary parameters, nor is there a general expression relating the derivative of the confluent Heun function to other confluent Heun functions, though particular instances have been obtained 91,92 . However, for the symmetric Pöschl-Teller potential (i.e. b = 0) one can obtain the relation: ( , , , , , ) In pristine graphene, only the symmetric form of Eq. (5) will contain non-leaky modes at zero energy. Non-zero-energy modes will have a finite lifetime since they can always couple to continuum states outside of the waveguide, whereas zero-energy modes are fully confined since the density of states vanishes at zero energy outside of the well. Asymmetric forms of Eq. (5) never contain truly bound modes since the density of states cannot vanish on both sides of the potential simultaneously. Notably, this is somewhat counterintuitive as for the Schrödinger problem a symmetric potential always contains a bound state, which can be removed by asymmetry. The emergence of bound states for a relativistic problem with an infinitely wide barrier is a manifestation of the Klein tunneling phenomenon 22,23 .
We shall now consider the symmetric form of Eq. (5) for massless particles. Accordingly, we set b = 0 and M = 0, and in this instance the symmetrized real functions 22,23 , are given by where Φ = H(α, β, γ, δ 1 , η 1 , Z) and θ β = −E tan / . By employing the identity Eq. (14), the derivatives appearing in Eq. (13) can be expressed in terms of Heun confluent functions. It immediately follows that at the origin Ψ Ψ + Ψ Ψ = * * 0 I II I II , which in terms of the functions Ψ A and Ψ B yields: . Therefore, the condition from which the eigenvalues of the spectrum are determined, Eq. (17), can be written as . This condition can be understood in terms of parity. In principle, one can construct from these functions odd and even solutions. However, since the even modes of Ψ I occur at the same energies as the odd modes of Ψ II and vice versa, one can obtain the eigenvalues when the symmetrized functions Ψ I or Ψ II are zero at the origin. The functions Ψ I and Ψ II are related to the earlier introduced functions Ψ 1 and Ψ 2 by Using Eqs (19,20) together with Eqs (9),(10) allows Ψ I and Ψ II to be expressed explicitly in terms of Heun functions.
Eq. (18) is formally the same as Eq. (13) but computationally faster. In Fig. 2 we plot the numerically obtained solutions of Eq. (18) for the potential defined by a = 24 and b = 0. The dashed lines represent the boundary at which the bound states merge with the continuum which occurs at the energies E = ±Δ and E + a/4 = Δ. For the potential of strength a = 24 we find that there are four zero-energy bound modes, occurring at Δ = 0.597, 2.276, 3.817 and 5.282. Their normalized wavefunctions are shown in Fig. 3.
As mentioned previously, the number of zero-energy modes is equal to the number of supercritical states (neglecting the spin and valley degrees of freedom). For the symmetric Pöschl-Teller potential, the eigenvalues of these supercritical states can be determined approximately, via simple analytic expressions. Moving to the symmetric basis, (Ψ I , Ψ II ) T , allows the pair of coupled first order differential equations, Eq. (3) and Eq. (4) to be reduced to a single second order differential equation in Ψ II : For supercritical states, E = −Δ, Eq. (21) transforms into the differential equation for the angular prolate spheroidal wave functions 93 :   where N = 1, 2, 3, … and q = 2N − 1 93 . Keeping only the terms of expansion shown in Eq. (23) yields the following eigenvalues: , it was found that there was a minimum potential strength of V 0 = 1/2, required to observe a zero-energy mode 22 . According to the Landauer formula, the conductance along the waveguide when the Fermi level is set to the Dirac point is 4ne 2 /h, where n is the number of zero-energy modes. The existence of a threshold in the potential strength needed for the waveguide to contain a zero-energy mode allowed us to suggest that such waveguides could be used as switchable devices. However, later numerical calculations utilizing a variable phase method implied that power-decaying potentials always possess a bound mode 58 . This result cast serious doubt in the validity of employing exponentially decaying potentials as a suitable model for graphene waveguides, since realistic potential profiles decay a power of distance rather than exponentially. Notably, the threshold potential strength at which the first zero-energy mode appears can be x x For even modes of Ψ I , C 2 = 0, whereas, for odd modes of Ψ I , C 1 = 0. In the absence of the potential, when E = Δ = 0 the two first order differential equations in Ψ I and Ψ II decouple, and Eq. (21) reduces to a first order differential equation. As E = −Δ → 0 and x → ±∞ (where the potential is zero), Eq. (25) and Eq. (26) are required to be linearly dependent 94 and the Wronskian of the solutions Ψ I and Ψ II is zero 95 . Consequently, the bound modes satisfy the condition: |Ψ I (±∞)| 2 = |Ψ II (±∞)| 2 . Therefore, the threshold potential strength at which the first zero-energy mode appears is found by the condition 0 Therefore, for any square-integrable potential, the threshold for the appearance of the first bound state of zero energy is only a function of the integrated potential. Notably, this is the same result obtained as relativistic Levinson theorem 45,[59][60][61][62] . For the Pöschl-Teller potential, Eq. (27) yields a = π, which agrees with Eq. (18). Finally, it should be noted that in the non-relativistic limit, Eq. (10) restores the well known results 14 for the bound state functions of the Schrödinger equation for the Pöschl-Teller potential. In the limit that α and δ are much smaller than β and γ (i.e. large Δ): ( ( , , , , , )) ( , 1 ;1 ; 1 ), and 2 F 1 is the Gauss hypergeometric function. Substituting E = E SE + Δ, b = 0 and s β = −s γ = 1 into Eq. (10), results in the non-relativistic bound state functions  Modified Pöschl-Teller potential solutions. We shall now consider the case of finite b, which represents a smooth asymmetric waveguide. Previously considered asymmetric waveguides varied abruptly on the same scale as the graphene lattice constant 54,[96][97][98][99] . We shall now consider more realistic smooth asymmetric waveguides, which fit closer to experimentally achievable potential profiles. In Fig. 4, we plot the energy spectrum for the potential defined by the parameters a = 24 and b = 2. The introduction of the asymmetry term b reduces the number of modes at E = 0, which are now quasi-bound modes for the massless case (Fig. 4a), since they can couple to continuum states outside of the waveguide. Naturally, for massive Dirac fermions full confinement is possible across a range of energies. In Fig. 5, we show a schematic diagram of the dispersion of a gapped Dirac material, subjected to the modified Pöschl-Teller potential defined by a = 24 and b = 2. For a particle of mass M = 2, it can be seen that for the energy range E = 0 to E = 2 there are no continuum states outside of the well, therefore in that range all bound solutions will be non-leaky. The corresponding energy spectrum of confined states is shown in Fig. 4b.

Conclusions
We have analyzed the behavior of quasi-relativistic two-dimensional particles subjected to a modified Pöschl-Teller potential. Our results have direct applications to electronic waveguides in 2D Dirac materials. For the symmetric Pöschl-Teller potential, explicit forms were obtained for the eigenvalues of supercritical states. A universal expression, for any symmetric potential, was obtained for the critical potential strength required to observe the first zero-energy state. The well known results for the Pöschl-Teller potential are recovered in the non-relativistic limit.