# Surface states of gapped electron systems and semi-metals

## Abstract

With a generic lattice model for electrons occupying a semi-infinite crystal with a hard surface, we study the eigenstates of the system with a bulk band gap (or with nodal points in the gap). The exact solution to the wave functions of scattering states is obtained. With the knowledge on the scattering states, we derive the criterion for the existence of the surface states. It is shown that a surface state is composed by the evanescent waves associated with a given energy band. There may exist the in-band as well as out-band surface states. All the surface states are classified by the continuous energy bands. For electrons in a system with time-reversal symmetry, we rigorously prove the correspondence between the change of Kramers degeneracy of the surface states and the bulk time-reversal invariant. The theory is applicable to systems of (topological) insulators, superconductors, and semi-metals. As examples, we solve the edge states of electrons with/without the spin-orbit interactions in graphene with a hard zigzag edge and that in -wave superconductor with a (1,1,0) surface.

## I introduction

Surface states (SSs) of electron systems (1); (2) can exist in crystals of superconductors (3); (4); (5), semi-metals (6); (7); (8), and topological insulators (TI) (9); (10); (11); (12); (13); (14); (15). The existence of SSs in interface between metals or superconductors leads to sizable electronic tunneling (16); (17); (18). The TI is characterized by the existence of the conducting SSs and an insulating bulk gap generated by the spin-orbit interactions (SOI). The materials existing SSs have prospective applications in electronic/spintronic devises. In particular, the TI can be used to conduct spins on the surface due to spin-Hall effect so that there is no electric resistance and no energy cost (9); (10); (11); (19); (20); (21); (22); (23); (24).

The existence of SSs in TI is considered as topologically protected. There is a correspondence between the existence of SSs and the topological property of the bulk states. For classification of the TI on the basis of the topological invariant (25); (26), Kane et. al. (19); (27); (28); (29) introduce the time reversal (TR) polarization and define the topological invariant for the bulk states. In analogous to Laughlin’s construction for the quantum Hall effect on a cylinder that the change of magnetic flux threading the cylinder can transfer electrons from one end to another through the cylinder (30), the momentum along a crystal axis parallel to the surface of a TI can be considered as the same role as the magnetic flux and its change leads to electron transfer from one surface to another. The transfer of electrons stems from the change of the number of SSs with the change of the magnetic flux (31); (32). In a TI, there are two type states for electrons because of the Kramers degeneracy. The invariant describes the change of the TR polarization between two invariant momenta. Therefore, the invariant is considered as the change of the Kramers degeneracy of SSs between two invariant momenta. A strong TI corresponds to an odd number change of the Kramers degeneracy of SSs. The invariant index introduced by Fu and Kane (29) has been used as a standard criterion for classification of the TIs with inversion symmetry.

In fact, the invariant is a characteristic of the bulk system rather than depending on the conditions of the crystal surface. After all, the edge states in the quantum Hall effect are the edge Landau states, they are certainly different from the SSs of a TI in the absence of the magnetic field. How the SSs form at the surface of a TI and a direct rigorous proof of the correspondence between the SSs and the bulk invariant are not given in the original theory of Kane et. al. (19); (27); (28); (29).

The SSs and the bulk-boundary correspondence have been studied by many theoretical works using various concrete models including lattice (33); (34); (35); (36); (37); (38); (39); (40) and continuous models (41); (42); (43); (44); (45). Most of the studies on the bulk-boundary correspondence are carried out in the context of quantum Hall effect (33); (46); (47). Qi et. al. have shown that the bulk topological quantum order characterized by a nonvanishing Chern number corresponds to the existence of the gapless edge states under certain twisted boundary conditions that allow tunneling between edges in two-dimensional insulators (33). In some of the existing works, the surface states are analytically solved by imposing special boundary conditions (36); (39). The mid-gap (or zero-energy) edge states are studied for the Dirac fermions and topological superconductors (35); (38). For the continuous model, a set of the parameters for the boundary conditions needs to be determined by the bulk system and the crystal potential near the surface (44). It is shown that the boundary conditions strongly affect the spectrum of the surface states and even the existence of the states near the zero momentum parallel to the surface. By the Green’s function theory for the TIs, it is shown that in the presence of interactions between electrons the existence of edge states at the boundary of two topological insulators with different topological invariants is not definitely determined (48). With a concrete model, besides the TRS, the system may have the additional parity or charge inversion symmetry. Certainly, the solution to a concrete model has model-dependence. A question still remains as how to rigorously prove the bulk-surface correspondence based on a generic model.

In principle, a surface state is a result of the superposition of the evanescent waves at the surface. By numerical calculation, one expands the SS in terms of all the evanescent waves with the coefficients determined by the Schrödinger equation with boundary condition (49); (50). At a given energy, there may not exist a SS or may exist a number of SSs. Thus, a criterion for the existence of the SSs is desirable. A SS may be composed by part of the evanescent waves rather than from all of them. If one knows and thereby chooses the only possible evanescent waves that result in the SS, the numerical calculation will be very effective. To resolve all these problems requires a fundamental theory on the formation of the SSs.

In this paper, with a generic lattice model, we exactly solve the scattering states of an electron system occupying a semi-infinite crystal with a hard surface. With the result, we derive the criterion for the existence of the SSs and determine the eigenvalues and wave functions of the SSs. We generalize the solution to cases of all possible complicated band structures of the bulk energy. We also analyze the classification of the SSs. For the electron systems with the TRS in TIs, in different from the analogous to the quantum hall effect, we rigorously prove the correspondence between the Kramers degeneracy of the SSs below the Fermi energy and the TR polarization introduced by Kane et. al.. As examples, we solve the edge states of electrons with/without the SOI in graphene for a number of cases and in a -wave superconductor.

## Ii Bulk States

We consider an electron system occupying a semi-infinite lattice with a hard surface. The axis of the coordinates is set as along the inner normal direction of the surface. The unit cell of the lattice contains atoms (or orbitals for an electron). The Hamiltonian of the system is given by

where with creating an electron of spin on th atom of th unit cell, is a matrix, and the sum runs over the unit cells in the semi-infinite space and . Since the momentum parallel to the surface, , is a good quantum number, we will work in the space of but real space along the axis. In this space, the Hamiltonian then reads with . We suppose that the electron hopping is confined within a rang: . For brevity of description, we hereafter may suppress the argument unless it leads to confusion it will be explicitly written out again. We will use the units in which .

First, we study the bulk eigenstates of electrons in the infinite lattice. The bulk states cab be used as the basis for investigating the eigenstates of the system with a surface. For the bulk states, we work in the momentum space. The wave function and the energy of the bulk states are determined by

(1) |

The transpose of the wave function is expressed as

(2) | |||||

where ’s and ’s are the components and is the normalization constant given by

Here, we analyze the property of the wave function. Define the variable with as the momentum along axis. Since is the Fourier transform of with , it can be written as where is a matrix; each element of is a polynomial of . The highest order of the polynomial among these elements is . The energy levels are determined by

(3) |

which is equivalent to

(4) |

In the plane, an energy level as function of defined on the circle . The energy bands correspond to -planes. With analytically continuation, the energy can be defined as function of in the plane. Clearly, is a analytical function of except at because of being analytical. The components ’s and ’s in the wave function can be defined analytical as well in the plane. On the other hand, the normalization constant is not analytical because it is the square root of sum of squares of norms of these components. At the point , the components ’s and ’s may be singular. We note that the exponents of the poles at can be changed by multiplying a power of to these components. By separating out the non-analytical part, the wave function defined in Eq. (2) is an analytical function of except at .

## Iii Preliminary Theory of Eigenstates of Electrons in Semi-Infinite Space

Studying the eigenstates in the semi-infinite space is essentially a one-dimensional problem. We proceed the analysis of eigenstates of electrons from the simpler case to more general cases. In this section, we present wave function of scattering states of the system in a single band without any degeneracy.

Incoming and outgoing waves. The wave function of a scattering state in which the electrons freely move in the system is superposed by the incoming and outgoing and evanescent waves. An incoming wave to the surface with wavenumber is defined as the wave with negative velocity . Accordingly, an outgoing wave is defined as that of positive velocity. The evanescent waves are those waves decaying with increasing the distance from the surface.

Since the energy is a periodic continuous function of the momentum, we suppose there are maximums (and also minimums) in the energy curve within the period and all the maximums (minimums) have the same value . Therefore, we have . Such an energy band is shown in the top panel of Fig. 1. The number of the incoming waves is the same as the outgoing waves. The total number of the incoming and outgoing waves at a given energy is , which is the number of the independent plane waves of bulk states degenerated at the same energy . Under the mapping, the incoming and outgoing wavenumbers at energy are mapped to points on the unit circle in the -plane.

For later use, here, we consider the property of the incoming and outgoing waves under the analytical continuation with changing their real wavenumbers to complex wavenumbers. Consider the energy to be analytically continued to complex energy with as a infinitesimal small quantity. Then, the momentum is changed to with determined by

(5) |

For a given , because of the dependence of velocity , the signs of for incoming and outgoing waves are opposite. Therefore, under the analytical continuation in -plane, with the points of incoming waves moving to inside (outside) the unit circle, the points of outgoing waves go outside (inside) it.

To briefly express the wave function of the scattering state, we hereafter denote the lattice coordinate simply as , and as . In addition, we take the wave function of plane wave with wavenumber in real space as

(6) |

with and . By attaching to the denominate, the wave functions satisfy the normalization condition

(7) | |||||

where is the total number of the lattice sites along direction. On the other hand, since the velocity of the evanescent wave of complex wavenumber with is zero, we take its wave function as

(8) | |||||

where . The evanescent wave functions satisfy the following equation

(9) |

Since the contribution to the summation comes from the sites close to the surface, the quantity is negligible small. On the other hand, the magnitude for the plane wave is divergent in the limit .

Scattering states. Now, consider a wave with wavenumber is incoming to the surface. It can be reflected by the surface to all the degenerated outgoing waves of wavenumbers ’s. Besides these waves, there are evanescent waves by the reflection. To see this, we consider the solutions to Eq. (4) for a given energy . Note that the order of the polynomial in Eq. (4) is . Besides the waves with real wavenumbers in the band under consideration, there are additional solutions of complex wavenumbers with equal number of positive and negative imaginary parts. These waves can be considered as the analytical continuation of incoming and outgoing waves of other bands. All the waves of wavenumbers with positive imaginary parts are the evanescent waves. The rest waves of wavenumbers with negative imaginary parts are growing waves, not satisfying the boundary condition at . These growing waves should be ignored in the scattering waves.

Wave function of scattering states. With the above discussion and the definitions, we are ready to express the wave function of the scattering state of the incoming wave with wavenumber reflected to all the outgoing waves ’s and the evanescent waves ’s all degenerated at energy . For brevity, we will denote and simply as and , respectively. Up to a normalization constant, the wave function of this state is given by

(10) | |||||

where and are constants depending on the energy .

The above analysis can also be applied to the process that a single outgoing wave comes from all the incoming waves and the evanescent waves under the surface reflections. The wave function is given by the similar formula as Eq. (10) where the wavenumbers and ’s are understood as the outgoing (with ) and incoming () waves, respectively. Hereafter, we allow the wavenumber variable in the whole range .

For the sake of description, we collect the wave functions in matrix forms as

where means . Now, the formula (10) can be written in the compact form

(11) |

where is a matrix with dimension and is a matrix with dimension .

All the incoming and outgoing and evanescent waves in the wave function given by Eq. (11) satisfy the Schrödinger equation in the bulk region (),

(12) |

The matrices and should be determined by the boundary condition at the surface.

Boundary condition. Close to the surface, for , the Schrödinger equation is different from Eq. (12). Now, it is given by

(13) | |||||

By extending the definition of wave function given by Eq. (10) to , we write the left hand side of Eq. (13) as

(14) | |||||

Note that the first term in the right hand side of Eq. (14) equals . From Eqs. (13) and (14), we have

(15) | |||||

This is the boundary condition for the wave function . We emphasize that from Eq. (15) one cannot assert the vanishing of the wave function at sites .

Remark. The surface of the lattice system is defined as the truncation of the hopping rather than an infinitive potential barrier. The electron cannot go out the lattice because of the vanishing of the hopping to a site outside the lattice. In real space, the wave function of an electron is defined within the lattice . Equation (15) says how the wave function would be if an electron went out the surface.

To proceed, we define the matrix

(16) | |||||

and matrices

(17) | |||||

(18) | |||||

(19) |

The dimensions of and are , and , respectively. The matrices and their product are analytical functions of . The matrix has been written in terms of the product of this analytical part and the non-analytical factor matrix . Similarly, we can define the corresponding matrices for the evanescent waves. For example, we define

(20) | |||||

(21) |

With these definitions, equation (15) reads

(22) |

By denoting

equation (22) reads

(23) |

Since the outgoing (or incoming) wave functions as the diagonal elements in and the evanescent wave functions in are all independent, the columns of matrix are linear independent vectors. Therefore, matrix is invertible. We thus have

(24) |

which is the solution to the constants and .

Here, we go further to analyze the solution. Recall that is constructed from the real waves of bulk states. Since the diagonal elements (column vectors) of are independent wave functions, the columns of are thereby independent vectors. Therefore, the rank of is . According to the algebra theory, there exist independent row vectors in . The other row vectors depend linearly on these independent vectors. These vectors can be eliminated by a row transformation . Therefore, can be transformed to a matrix with a block of square matrix of rank and the rest block as zero matrix,

(25) |

Similarly, we have

(26) |

Under the transform , the matrix is transformed to with . Denote as

(27) |

Its inverse is given by

(28) |

Write in the block form,

(29) |

where the dimension of is , and the dimensions of other three matrices are automatically determined. With these notations, the solution to reads

(30) | |||||

which means

(31) | |||||

(32) |

Actually, the matrix is obtained from the matrix by the row transformation . Since is analytical function of , the matrix is also analytical. The matrix can be considered as a function of one momentum, for example, since all the momenta are associated with the same energy . In , only the analytical components ’s and ’s in the wave functions appearing in Eq. (2) are involved. In the -plane under the mapping , we can suppose that is not the (only possible) pole of because it can be gauged away by multiplying a power of to the components ’s and ’s as noted earlier. Therefore, the matrix is an analytical function of inside the unit circle .

Surface States. Now, we consider the analytical continuation of the scattering wave function in the -plane. By changing the variable from the unit circle to the inside, the energy changes to the outside of the energy band. As discussed after Eq. (5), when goes into the unit circle, the exponential with varies to outside the unit circle, giving rise to the waves with growing amplitudes as increasing in the wave function given by Eq. (11). However, at a point, , there may exist such a nonzero vector (of components) that

(33) |

We therefore have . Applying Eq. (33) to Eq. (11), we get a wave function

(34) |

where all the components are evanescent waves with complex wavenumbers all with positive imaginary parts. The components of the vector are the amplitudes of the corresponding evanescent waves in . Thus, the wave function describes a surface state.

From Eq. (34), we get a conclusion that the SS consists of evanescent waves only belonging to the band under consideration. There is no any other evanescent wave from other bands in the components of the SS. This is different from the scattering states in which the evanescent waves from other bands are mixed in. This is an interesting result.

In other words, the above conclusion can be expressed as that all the SSs of the system are classified by the continuous bands.

Criterion for the existence of SS. The existence of SS is determined by Eq. (33). We denote as . Since all the elements in the diagonal of are nonzero numbers, the nonzero vector corresponds to the nonzero vector . The existing of the nonzero vector means that the determinant of vanishes

(35) |

This is the criterion for the existence of the SSs. This is another central result of this work.

The matrix here takes the role of the generalized Jost function (matrix) in the theory of eigenstates of electrons in a field of central force (51); (52). In that case, the existence of bound states is determined by the zeros of the determinant of the generalized Jost function. All the bound states are classified with the angular momentum. However, the relationship between and the scattering matrix in the present case is not the same as that in the central force problem; there are the additional evanescent waves from other bands in the present case.

To count the number of the SSs of a band, we need to find out all the zeros of . When arrives at , the other variables simultaneously arrive at inside the unit circle. When arrives at any one of these points, vanishes. Therefore, all these points are the zeros of . All these zeros correspond to a single SS. To count the number of zeros of , one usually performs contour integral in the complex plane.

We note that close to some kind of zeros , may vary as with an integer. We need to distinguish these irregular zeros from the regular ones in . For regular zeros , with as positive integer when close to . ( can be an integer larger than unity for the case when the regular zeros can group into multiple ones.) If there is only one kind of such irregular zeros, can be written as with and as the regular and irregular parts, respectively. For counting the SSs with irregular zeros of by performing contour integral, the contour should surround each of these irregular zero points times. Under these considerations, the number of the surface states of the energy band is then given by

(36) |

where the dependence of has been written out explicitly, and the factor eliminates the multiple counting. We will give an example for the existence of the irregular zeros later.

Equation (36) is different from that for the bound states of a central force. In the latter case, the relevant zeros for the bound states are all regular ones () (51); (52). The reason is that the energy-momentum relation for free electrons in that case is the simple parabolic form rather than the complicated band structures discussed here.

As seen from Eq. (36), in general, the integral

may not be an integer. For the case there is only one kind of irregular zeros, we can infer the exponent of the irregular zeros from the value of integral .

## Iv General Band Structures

In real problems, the energy band structures are generally not the simple one as we discussed above. Here, we extend the above analysis to more general cases with complicated band structures.

(i) Band with non-equal-valued maximums (minimums). Without loss of generality, we consider the case that there is one maximum with a different value from the highest maximum. Such an energy band with is depicted in Fig. 2. Above a critical energy , the number of the incoming waves (also of the outgoing waves) is . The absent incoming and outgoing waves are lost at the wavenumber . Actually, the wavenumbers of the absent incoming and outgoing waves for have become complex and , respectively. They are determined by . By expanding the energy at , we have

To satisfy , we must have and the vanishing of all odd order derivatives of the energy with respective to . For , the wave of is an evanescent wave damping with distance from the surface, while the wave of is growing with distance. As depicted in top panel in Fig. 2, we use the red-dashed vertical line to present the absent waves. Under the mapping , this vertical line is mapped to the line along direction with angle equal to .

For the scattering states of energy with , we need to add the evanescent wave to and to drop the growing wave from since it does not satisfy the boundary condition as the distance becoming infinitive. We then get a similar expression for the wave function as given by Eq. (11).

In-band surface states. On the other hand, there may be the possibility that the evanescent wave evolves to a SS. This is an in-band SS. To get the criterion for the existence of this SS, we still consider the matrix with dimension for . Now, one column of , say the th column, is associated with the evanescent wave. Suppose by row transformation is transformed to

(37) |

The block is associated with the real (incoming or outgoing) waves and thereby . Under analytical continuation as going above from below, if at a certain energy (), we will get an nonzero vector (with the superscript denoting the transpose) satisfying and thus obtain a SS coming from the evanescent wave. Clearly, the criterion for the existence of this SS is still given by Eq. (35).

As shown in Fi. 2, for energy in the regime , the real incoming and outgoing waves are defined respectively on the arcs and on the unit circle in the -plane, while the evanescent wave is defined on the line along the radial direction. When the real wavenumber varies from , the energy varies as . Meanwhile, the point reflecting the evanescent wave moves as in the -plane. The green points and are the possible zeros of giving rise to a single in-band SS.

The contour of the integral in Eq. (36) is along the entire unit circle with a small arc close to point shown in lower panel of Fig. 2. With the small arc close to the point in the contour, all the zeros for the in-band state are included in the integral. Since the evanescent wave defined on the line in composing the matrix and the real wave defined in the arc are associated with the same energy, the small arc close to point is a correspondence to the one close to .

When increasing the energy , the point moves toward . As the energy going above the upper bound of the band, all the waves involved in composing the matrix become evanescent waves and the original zero point moves inside the unit circle on the radial line . Here, the angle of radial line is the momentum of another maximum of the energy band.

For a more general band structure, all the minimums may not equal valued. We can extend the above analysis to this case as well. We use the dash-dot lines to define the corresponding evanescent waves. For the above example, the dash-dot lines in unit circle in Fig. 2 are depicted for the evanescent waves with energy below the lower bound of the band. The angle of each dash-dot line is the momentum of the corresponding minimum of the band.

There may be the cases that the highest minimum energy is larger than the lowest maximum energy. An in-band SS may consist of the evanescent waves defined in both dash and dash-dot lines.

Zeros distribution. From the example, we get the following conclusions. (1) For the in-band SSs, some of the zeros of distribute on the (dashed and/or dash-dot) radial lines in the unit circle and the rest zeros distribute on the unit circle. (2) For the SSs with energy above the top bound of the band, all the zeros distribute on the dashed radial lines. (3) For the SSs with energy below the lower bound of the band, all the zeros distribute on the dash-dot radial lines. On these lines, since the energy is real, the zeros of give rise to physical SSs.

(ii) Overlapped bands. When there exist overlapping between energy bands, an incoming wave can be reflected to outgoing waves of all overlapped bands. All the overlapped bands should be considered as one band. The upper bound of the band is given by the highest maximum in the overlapped bands, while the lower bound is the lowest minimum. The rank is then the sum of the individual bands. An example of two bands overlapping is shown in Fig. 3. The vertical lines for defining the evanescent waves are also depicted in Fig. 3. Now, the incoming (outgoing) waves include all that of the individual bands. Thus, the formalism given in Sec. III along with the extension discussed in case (i) is applicable to the present case.

To apply the contour integral given by Eq. (36) for the present case, here, we need to clearly define the mapping and the contour in the -plane for energy bands with overlapping.

The mapping. We take of the red band. Then, all other ’s of red band and ’s of blue band in the incoming (outgoing) waves are determined by . We can also take , which is a function of . Under the mapping, the red and blue bands are mapped to -plane and -plane, respectively. The vertical lines are mapped to lines in the corresponding planes. Denote the numbers of maximums of red band as and the blue band as . Their least common multiple is . Now, the closed contour is defined as running times on the unit circle in the red band plane. With this definition, runs times on the unit circle in the -plane.

As shown in Fig. 3, the zeros of are the and points in the -plane and -plane, respectively. The points in -plane are equivalent. They give rise to a single SS. We need to perform the contour integral in one plane, for example, the red plane. Since is a function of , the zeros in -plane are automatically counted. In -plane, since the contour surrounding the zeros times, we need to divide the contour integral by to eliminate the multiple counting of the number of SSs. Of course, instead of the -plane, the integral can be performed in the -plane.

The analysis can be extended to the case of overlapped bands of a number larger than two.

(iii) Degenerated bands. For an energy band with degeneracy , there may exist surface reflections between the degenerated states. This is a special case of bands with overlapping. The wave functions for the scattering states can be still written as Eq. (11) but with the components of understood as

Though the momenta of the components of are the same , they are independent column vectors. By definition, the rank of now is . On the other hand, the least common multiple of the numbers of ’s is its self. The formula for counting the number of SSs is the same as Eq. (36).

So far, we have considered all possible cases of the band structures in real problems.

For latter use, we give an equivalent formula for . We formally write it as

(38) |

where is a functional of matrix . The functional means the process of generating out with row transformation, putting (or for overlapped bands) as the exponent to , and then doing special treatment for the irregular zeros of .

## V System with TRS

Here, we apply the previous results to the electron system with the TRS.

Time-reversal symmetry implies that the electron system is invariant under the operation with the Pauli matrix (operating in the spin space) and the complex conjugation operator. In real space under consideration, the TRS is reflected by .

For the bulk states, by the TRS, , there exists the Kramers degeneracy between the bulk states and with energy . Therefore, the energy bands come in pairs as depicted in Fig. 4. The eigenstates are divided into two types of states and . For the case of only the Kramers degeneracy existing, these eigenstates satisfy the relation (27)

(39) |

with as a phase quantity. There may exist other degeneracy in the states because of other possible symmetries. The TR operator is then reflected by an unitary transformation . Suppose the degeneracy in states of each type is . All the energy bands will be grouped to bands. By writing the wave function in a compact form, the action of the TR operator on the degenerated states is given by

(40) |

where is

and is a unitary matrix. For the case of two type states not degenerated at the same momentum, is an off-diagonal-block matrix

(41) |

where is a unitary matrix, and is the transpose of .

At the TRS points ’s (in the Brillouin zone) that , is an asymmetric matrix because of .

Now, we study the eigenstates of electrons in the system with a surface. We confine ourselves to the case of .

(i) Number difference between the two type SSs. Here, we temporarily suppose the two type states are not degenerated at the same momentum and there is no surface reflection between the waves of different type. We then can calculate the number of the SSs of each type. Type I and type II surface states are determined by matrix and , respectively. The two matrices are related by

(42) | |||||

where is a diagonal (block) matrix given by

Here, is the number of the maximum of energy band for type I (II) states. Without any other degeneracy except the Kramers degeneracy, is given by as defined by Eq. (39) with the index of the band under consideration. For the case there exist degeneracy other than the Kramers degeneracy in the type I states, is a matrix as that appeared in Eq. (