Boundary Conditions of Weyl Semimetals

We find that generic boundary conditions of Weyl semimetal is dictated by only a single real parameter, in the continuum limit. We determine how the energy dispersions (the Fermi arcs) and the wave functions of edge states depend on this parameter. Lattice models are found to be consistent with our generic observation. Furthermore, the enhanced parameter space of the boundary condition is shown to support a novel topological number.


I. INTRODUCTION
One of the important aspects of topological phases [1,2] is that they bridge condensed matter physics and particle physics. Familiar and important concepts in particle physics, such as topological charges, quantum anomalies and relativity, are completely translated to condensed matter physics, and their experimental realization flourishes and verifies the concept quite nontrivially. The interplay between the condensed matter physics and the particle physics is further expected to provide developments of these fields.
From the theoretical side, topological phases are classified by dimensions and discrete symmetries [3,4]. The symmetries are important for classification, as they are robust against deformations. In particular, at microscopic levels materials have lots of ways for their Hamiltonian to be deformed. However, if one is based on the topological properties which are determined by the discrete symmetries, the classification and the resultant phenomena are robust. This works even at a continuum limit where detailed structure of lattice Hamiltonian disappears. The outcome of the topological phases is the existence of gapless edge mode, as a consequence of the renowned bulk-edge correspondence [5][6][7].
Particle physics Lagrangians are written based on symmetries, too. Normally, particle physics Lagrangians have all possible terms allowed by the required symmetries, in the continuum limit. Then a natural question arises: what are all possible boundary conditions? This question is particularly important in view of the bulkedge correspondence, because one typically uses open boundary conditions in the study of topological phases.
Although to answer this question at the level of listing all possible lattice Hamiltonians with boundaries is far beyond our knowledge, we could approach the answer in the continuum limit. The continuum Hamiltonian is much easier to deal with, and we can ask how general the Hamiltonian at our hand is. Furthermore, it is expected to capture the generic feature of all lattice systems which share the same discrete symmetries as the continuum Hamiltonian.
In this paper, we focus on the 3D Weyl semimetal, which has been recently observed in experiments [8][9][10] after theoretical predictions [11][12][13][14][15][16][17], and study the most generic boundary condition in the continuum limit. We will see how a generic surface term in the Lagrangian can affect the edge states of a Weyl fermion, in a Hamiltonian language -in other words, we will see how a generic modification of the boundary conditions give change in edge dispersions. Surely the existence of edge states for a given Weyl semimetal is explained by the topological number of bulk theory, but how the boundary conditions are related with the edge states, or so-called Fermi arcs, have not been studied well in the literature. [18] To look at genericity under a given symmetry is a particle theory standpoint, which could serve as another new interplay between the two fields.
Our study is divided into two parts, the study of the continuum theory and a verification by lattice models. On the continuum theory side, we start with the standard Hamiltonian with a single Weyl cone, and consider the most general boundary condition, with resultant edge states. In particular, we are careful with allowed number of parameters. We also study the dimensional reduction to 2D, at which topological insulators of class A can be studied. On the lattice model side, a graphene in 2D and a square lattice in 3D are considered.
Our main findings are summarized as follows: • The generic boundary condition of the 3D Weyl semimetal is dictated only by a single real parameter. At low energy, all edge states are labeled by the parameter.
• The parameter describes how the edge dispersion is attached to the Weyl node.
• Results of the lattice models are perfectly consistent with the continuum theory in the vicinity of the Weyl point.
In the course of our study, we also find novel properties of the Weyl semimetals: • Our results of the edge states are all consistent with the bulk-edge correspondence. In particular we propose how to count the edge states for the Weyl semimetals.
• A new topological structure emerges in the enhanced parameter space spanned by the boundary parameter and the two conserved momentum. The edge state carries the topological charge.
Our paper is organized as follows. In Sec. II, we study the 3D Weyl semimetals and their generic boundary conditions in the continuum limit. We derive the relations between energy dispersions, wave functions of edge states and the boundary conditions. We study also a dimensional reduction to 2D topological insulators. In Sec. III, we investigate lattice models and show the consistency with the results obtained in the continuum limit. In Sec. IV we check the bulk-edge correspondence. We further propose how to count the edge states of the Weyl semimetals. In Sec. V we find a new topological structure of the edge states in the enhanced parameter space of the boundary conditions. In the last section, we summarize our paper. In Appendix, we provide two examples of generalization: treatment of generic momentumdependent boundary conditions, and the case with two parallel boundaries, in the continuum limit of the model of the Weyl semimetals.

II. GENERIC BOUNDARIES OF 3D WEYL SEMIMETALS AT CONTINUUM
We start with a generic Weyl semimetal in 3 spatial dimensions. Near the Weyl point, the Hamiltonian in the continuum limit is generically given by Here the Weyl point sits at the origin of the 3-dimensional momentum space. We are interested in the most generic boundary conditions and the spectra, of this system. The states are energy eigenstates, subject to the following boundary condition: where we put a boundary at x 3 = 0, so the material is in the spatial region x 3 ≥ 0. A real constant is the energy eigenvalue, and M is a generic 2 × 2 complex constant matrix. (See App. A for generic momentum-dependent boundary conditions.) The Hamiltonian can effectively describe a two-band system that has a degeneracy at the Weyl point with a definite chirality. It captures the topological nature around the Weyl point, and is equivalent to the Hamiltonian of a Weyl fermion in 3 spatial dimensions. The boundary condition needs to be of the form (3) since at the boundary two components of ψ are related to each other through the boundary condition. The matrix M reflects the arbitrariness in the choice of the boundary condition; One can think of infinitely many kinds of boundaries, starting from just a slicing of the material, to putting various chemical layers on top of the boundary surfaces, such as hydrogen/nitrogen termination or oxidization.
We solve the equations above, which give energy dispersion relations and wave functions. We find that the most generic boundary condition (3) is parameterized only by a single real constant θ + ∈ [0, π).

A. Parameterizing generic boundary condition
In this subsection we explore all the requirements for the boundary condition matrix M in (3) and parametrize it. It turns out that M is parameterized only by two real parameters, and furthermore, the boundary condition (3) needs only one of them. These result from a selfconjugacy of the Hamiltonian (1) and the requirement that M needs to have the eigenvalue −1: the combination M + 1 needs to have a zero eigenvalue, for a nontrivial solution satisfying the boundary condition (3) to exist.
with complex coefficients a 0 = A 0 + iB 0 and a i = A i + iB i ∈ C, where i = 1, 2, 3. We need our Hamiltonian (1) to be self-conjugate, which gives a constraint on the boundary condition. See also [19]. The self-conjugacy condition of Hamiltonian is for arbitrary normalizable ψ 1 and ψ 2 . When we explicitly write the two inner product above as an integration, we find a surface difference between the right hand side and the left hand side, which must vanish: Applying the boundary condition (3) onto this equation gives which is satisfied for any choice of ψ 1 and ψ 2 only when This is required by the self-conjugacy of the Hamiltonian. Equation (7) removes four real d.o.f of M : and results in So we are left with the boundary condition matrix with four real parameters,

Eigenvalues of M
The boundary condition (3) can be regarded as an eigen equation for matrix M . Substituting equation (9) into the determinant of (3) we get The boundary condition requires M to have a real eigenvalue −1, which gives two constraints: As a result, the generic boundary condition matrix should be written as Consequently, M is parametrized by two real parameters. We can choose for parameterizing the matrix, with which the boundary condition can be written as sech χ + i tanh χ cos θ − i sin θ cos θ + i sin θ sech χ − i tanh χ ψ Defining cos θ = sech χ and sin θ = tanh χ and changing variables: we find that the equation (17) becomes Noting a relation the boundary condition is recast to the following simple form It is surprising that actually for the boundary condition we need only a single real parameter θ + . So we conclude that the most generic boundary condition is just dictated by a single real parameter. This equation (20) tells us further that, at the boundary, two components of the fermion need to have the identical magnitude, and the relative phase between them is determined by θ + . This is true for the edge modes as well as the bulk modes.
For our later purpose, we determine here the range of the parameter θ + . First, from the definition of θ and θ , we find that they live on a region 0 < θ ≤ 2π, 0 < θ ≤ π .
Note that sin θ = sech χ ≥ 0. Then one notices that the region can be equally covered by Resultantly, the smallest necessary region for θ + is given by θ + = 0 and θ + = π are the same configuration up to some adjustment of θ − which does not appear in the boundary condition itself.

B. Lagrangian formulation
Since Lagrangian formulation sometimes works easier, here we present an equivalent Lagrangian formulation of what we have seen in terms of the Hamiltonian. In fact we find that the condition (7) shows up naturally in the Lagrangian formulation. Let us describe a generic and consistent boundary condition of a Weyl semimetal in 1+3 spacetime dimensions. Metric convention is chosen as η µν = diag(+, −, −, −) µν .
The bulk Lagrangian (for a right-handed Weyl fermion) is written as where σ µ = (1 2 , σ 1 , σ 2 , σ 3 ). The Dirac equation is which can be rewritten as where i = 1, 2, 3. So the Hamiltonian is i∂ 0 = H, which is the standard Hamiltonian of the Weyl semimetal near the Weyl cone. Let us introduce a surface term in the Lagrangian, for deriving the boundary condition. The total action is The first term is the Weyl Lagrangian. The second integral is with a Hermitian matrix N . A variation ψ → ψ + δψ and ψ † → ψ † + δψ † provides equations at the surface x 3 = 0 as For this to be valid for arbitrary δψ and δψ, we find at the boundary x 3 = 0. These two equations are complex-conjugate to each other. If we write N = iσ 3 M , then (M + 1)ψ The Hermiticity condition N = N † is now written by M as which is found to be equivalent to the Hamiltonian conjugacy constraint on M , (7). Note that this condition follows from the Hermiticity of N , that is, the Hermiticity of the surface Lagrangian. So, in the end, we found that the boundary condition is dictated by a boundary "mass" term with a Hermitian matrix N , The most generic boundary condition is given by a constant Hermitian matrix N , and that is physically natural.

C. Generic edge modes and dispersions
Since the Weyl fermion possesses a topological number, we expect the existence of the edge-localized modes when we introduce a boundary x 3 = 0. In this subsection we look for edge mode solution of the energy eigenvalue problem. With the most generic boundary condition (20) we obtained, we find a dispersion relation and the wave function of the edge mode, which are completely specified by p 1 , p 2 and θ + .

Solving eigenstate equation
Now we look for edge mode solution to eigenvalue equation (2). With an explicit two-component notation the eigenstate equation (2) can be written as This equation can be reorganized into two independent second-order differential equations: We look for the modes localized at the boundary. For the edge modes, we need then the corresponding solutions required by the normalizability are where ξ 0 and η 0 have no dependence on x 3 . These are the edge modes, and in the following we determine the dispersion (p 1 , p 2 ) and the relation between the components ξ 0 and η 0 .

Dispersion relation
We combine the results from eigenvalue equation (2) and boundary condition (3) for edge eigen modes. Substituting equations (37) and (38) into equation (35), we get one independent equation: Writing the boundary condition (20) and equation (39) together, we have This matrix equation contains the essence of the eigenvalue equation (2) and the boundary condition (3) for edge eigen modes. The vanishing determinant condition of (40) gives We move the term to the right, square both sides and cancel p 1 − ip 2 , finding the relation: This is the dispersion relation of the edge states. It is linear with respect to p 1 and p 2 , and speed of light is now anisotropic. Substituting (42) back into equation (41) we also find We write above two equations in a compact way: Interestingly, (44) shows that what the boundary does is only rotating the momenta (p 1 , p 2 ) into ( , α), the energy and the inverse of edge mode decay width (penetration depth). For fixed p 1 and p 2 , we can regard the pair ( , α) as a vector rotating around the origin by 2θ + . When the absolute value of becomes large, α becomes small, then the penetration depth is large. On the other hand, when the absolute value of becomes small, α becomes large and then the penetration depth is small. This coincides with the intuition that the wave function penetration measured from the location of the boundary increases for larger energy of the edge mode.
Plotting the dispersion relation, we actually see in Fig. 1 that the edge dispersion is rotated against the (p 1 , p 2 ) axes by the change of the boundary parameter θ + .

Wave function of edge modes
Let us finally write the wave function of the edge states. We have already used up most of the information and are left with normalization condition only, with which we can determine the wave function completely. Substituting (38) to the normalization condition we obtain a constraint With the second equation of (40), the boundary condition, we can see that the two components should have the same magnitude with a difference of their phases. Combined with (46), they are determined up to an irrelevant overall phase: So the general edge mode wave function is Note that the edge modes exist only in a limited region of the momentum space, since we need to require α > 0. The linear inequality α > 0 specifies a half of the momentum space, only in which the dispersion exists, see Fig. 1.
In the limit α = 0, that is, on the line p 1 sin 2θ + − p 2 cos 2θ + = 0 in the momentum space, the edge mode approaches a non-normalizable mode, which is a constant wave function in the x 3 space. It corresponds to p 3 = 0 bulk mode, whose dispersion is = ± p 2 1 + p 2 2 . In fact, the edge dispersion (42) is identical to that under the  3 4 π, π/2, π/4, 0, −π/4, −π/2, − 3 4 π, respectively. condition α = 0. Therefore we have a consistent picture for any value of θ + : when the edge mode approaches a non-normalizable state in the momentum space, it is consistently and continuously absorbed into the bulk modes. In Fig. 1, we find explicitly that the edge dispersion surface has its boundary on the bulk dispersion surface.
We would like to make one comment about how we could modify the range of α in different setups. If we introduce two boundaries which are parallel to each other, then the condition of the positivity of α does not apply, as the wave functions are normalizable even for a negative value of α. We demonstrate the calculation in App. B.
In summary, we find that the dispersion of the edge state is attached to the bulk Weyl cone in such a way that (i) the edge dispersion is tangential to the Weyl cone, and (ii) the edge dispersion ends at the touching line on the Weyl cone.

D. Reduction to 2D
It is important that the analysis given above can be consistently translated to 2D topological insulator of class A. It is just a dimensional reduction of the pre-vious Hamiltonian from 3D to 2D, by a replacement of one of the momenta -p 2 with a constant mass parameter m. This means that we can study most generic boundary condition of the class A topological insulator in the continuum limit, and its consequence in the edge dispersions.
By the dimensional reduction, the Hamiltonian of the 2D gapped fermion is given as The analysis of the boundary condition we had before for the Weyl semimetals does not change, since it is just a renaming of p 2 . So it is identical to our previous (20): The dispersion relation and the inverse decay width α are given simply by a replacement of p 2 with m: The same is applied for the edge mode wave function: Fig. 2 shows the bulk and the edge dispersions for various choices of the boundary parameter θ + . Since our procedure is just replacing the momentum p 2 by a constant m, it amounts to choosing a plane of constant p 2 in the 3D Weyl semimetal dispersion given in Fig. 1. Taking a cross-section, we find that the 3D Weyl dispersion and the edge dispersion reduce to dispersions of the gapped bulk and the linear edge modes in 2D. It is interesting that the rotation 2θ + in the (p 1 , p 2 ) plane for the 3D Weyl semimetals can inherits its nature in the 2D topological gapped system in an nontrivial manner. The form of the edge dispersion, as a function of p 1 , looks quite nontrivial in Fig. 2. For some special choice of the value θ + = π/2, the edge dispersion eventually disappear. For some other values of θ + , the edge dispersion becomes a flat band.
By taking a massless limit m = 0 for the bulk system, the edge dispersion is simply given by The existence condition of the edge state is p 1 sin 2θ + (= α) > 0. So the edge dispersion, which is a half line, emanates linearly from the Dirac point of the graphene by the slope cos 2θ + , where the parameter θ + can range 0 ≤ θ + < π.

III. LATTICE MODELS
The effective model study shown above exhibits an interesting behavior of the edge state depending on the boundary condition. Let us then show how such an argument on the boundary condition is realized in lattice models with tight-binding Hamiltonians.

A. Boundary condition for discretized model
In the effective continuum theory the boundary condition requires some conditions due to self-conjugacy of the Hamiltonian. Following this argument, we consider the boundary condition with the discretized lattice model.
First of all, we cannot directly apply the continuum theory argument to the lattice model because this argument relies on the integral by parts: We need to replace the differential operator with a difference operator which does not satisfy the Leibniz rule. We have to be careful about dealing with the boundary of the discrete lattice system.
To demonstrate how the self-conjugacy characterizes the boundary condition, we consider a discrete model defined on a finite one-dimensional lattice labeled by n = 1, . . . , N . The self-conjugate operator we consider here is H = −iσ∇ where σ is a Hermitian matrix to be taken as a Pauli matrix, and the difference operator is defined This difference operator reduces to the differential operator in the continuum limit, so that the operator becomes the standard Dirac Hamiltonian H → −iσ∂ x in the limit. Since they are related to each other, i∇ † ψ n+1 = −i∇ψ n , this is locally self-conjugate. However, as discussed before, we need to take care of the boundary: The discrete Dirac Hamiltonian is self-conjugate up to the boundary term (56) where we introduced auxiliary fields ψ 0 and ψ N +1 . The second line shows the surface term in this case, and the self-conjugacy of the Hamiltonian requires that this part should vanish We have two possibilities to solve this condition. The first is the periodic boundary condition ψ n = ψ n+N for ∀n ∈ {1, . . . , N }. Then these two terms cancel each other. The second is the situation that the both two terms vanish independently, which corresponds to the open boundary condition. Let us focus on the first term ψ † 0 (iσ)ψ 1 , and apply the boundary condition which is analogous to that considered in continuum theory We remark that this boundary condition is assigned to both ψ n=0 and ψ n=1 , although the former one is just an auxiliary field. If this matrix M satisfies the first term vanishes due to the same argument given in the continuum theory shown in Sec. II A. We also apply a similar condition to the opposite boundary n = N, N + 1 with a matrix M which is not necessarily the same as M as long as satisfying the condition (59).

B. 3D Weyl semimetal model
We then incorporate the boundary condition to the Weyl semimetal model on a lattice. We consider the Hamiltonian defined on a 3D lattice, where n is a three-dimensional vector n = (n 1 , n 2 , n 3 ) ∈ Z 3 , and the operator is given by In the Fourier basis, the Bloch Hamiltonian is obtained as which exhibits Weyl points at .
The parameter c tunes the Weyl point positions. The energy band spectrum is drawn in Fig. 3 with p 3 = 0 and c = 1. We see two Weyl points at (p 1 , p 2 ) = (±π/2, 0) at this section. We introduce the boundary to this model. Suppose the lattice is defined on the region n 3 ≥ 1, and impose the boundary condition where the matrix M satisfies Then the situation is completely parallel with the continuum theory studied in the previous section. The matrix M is parametrized by two parameters, θ + and θ − , and the boundary condition is rephrased in terms of these parameters, which is equivalent to Thus it depends only on the parameter θ + in the end. We consider the spectrum and wave function behavior of the edge state under the boundary condition (64). For the eigenvalue equation the Hamiltonian has a matrix form in the partial Fourier basis, where the off-diagonal element is given by which behaves as ∆(p) ∼ p 1 ± ip 2 in the vicinity of the Weyl points. The sign ± depends on the chirality of the Weyl points.
Assuming that the wavefunction is given by with the real parameter |β| ≤ 1, the eigenvalue equation (68) is equivalent to To obtain a non-trivial solution to this zero mode equation, we asign the condition det D = 0 which yields There are two solutions for β ≥ 0 and β ≤ 0. We remark that these two possibilities correspond to the doublers at p 3 = 0 and π in the momentum space.
Then, together with the boundary condition (67), the zero mode equation (72) gives Since β ∈ R, we obtain (p) = − cos 2θ + Re ∆(p) − sin 2θ + Im ∆(p) , (76) α(p) = sin 2θ + Re ∆(p) − cos 2θ + Im ∆(p) , which is rewritten as a matrix form where we defineα(p) := 1 2 β −1 − β , and from (70), the Comparing with the continuum theory (44), now the situation is parallel under the replacemtnt (p 1 , p 2 , α(p)) −→ (Re ∆(p), Im ∆(p),α(p)) . Fig. 4 shows the boundary parameter dependence of the bulk and edge state dispersions. The edge state spectrum has a support only where the normalizability conditioin is satisfied |β| ≤ 1. As mentioned before, there are two solutions corresponding to positive and negative β. We focus on the positive solution in the following. When we change the parameter θ + , the edge state spectrum rotates around the Weyl points. The orientation, that is, how the edge state spectrum winds, depends on the chirality of the Weyl nodes. This result is consistent with the continuum theory analysis in particular in the vicinity of the Weyl points. To see the parameter dependence more explicitly, let us take the constant energy section of the spectrum, which yields the Fermi arc, shown in Figs. 5 and 6. This shows that the parameter characterizing the boundary condition θ + plays a role of the rotation angle of the Fermi arc, as studied in continuum theory. In the present case of the lattice models, the Fermi arc ends on the Weyl points and have a finite support in the momentum space. Such a behavior of the Fermi arc has been experimentally observed, for example, in the transition metal pnictide family [20].
Let us comment on the Fermi arc behavior corresponding to the negative β solution. In this case the Fermi arc appears in the regionα(p) < 0, which is complement tõ α(p) > 0 for the positive β solution. In addition, the rotation orientation for β < 0, depending on the boundary condition θ + , is the same as that for β > 0, as shown in Fig. 7, and thus the winding number is also the same for both cases. This implies that these two contributions from the positive and negative β solutions are not canceled with each other.

C. Reduction to 2D class A system
As mentioned in Sec. II D for continuum theory, the 3D Weyl semimetal system is translated to 2D class A system using the dimensional reduction. Let us try this reduction to apply the systematic study on a lattice, and see how the topological edge state behaves under the generic boundary condition.
Replacing the momentum p 1 or p 2 with the constant mass parameter m in the 3D Hamiltonian (83), we obtain the 2D class A system This is the dimensional reduction along the p 1,2direction. Putting the boundary at n 3 = 1, we can apply the same self-conjugacy argument to obtain the boundary condition as the 3D system (66). Thus the edge state spectrum is given by Figs. 8 and 9 show the bulk and edge state dispersions depending on the boundary condition parameter θ + , which is again consistent with the continuum theory. Such a dependence of the boundary condition was recently predicted to be observed in monolayer silicene/germanene/stanene nanoribbons [21].

IV. THE BULK-EDGE CORRESPONDENCE
In this section, we study the relation between the bulk and the edge states, and check the bulk-edge correspondence. We will see that for the most generic boundary conditions, the bulk-edge correspondence for the 2D topological phase works perfectly. The bulk-edge correspondence [5][6][7] for topological insulators is well-known, while that for 3D Weyl semimetals has been understood in a way through a dimensional reduction to 2D. We claim here that the bulk-edge correspondence for the 3D Weyl semimetals can be defined as follows: the topological number of the bulk counts the chirality of the Weyl fermions, while the topological number of the edge is defined through the orientation of the Fermi arcs attached to the Weyl cone. This is based on our analysis on most general boundary conditions of the Weyl semimetals. In the following, first we investigate the case of the two dimensions, then later we discuss the case of the three dimensions, the Weyl semimetals.

A. 2D topological phases with the most general boundary conditions
For the bulk-edge correspondence for the 2D topological insulator of class A, the formula is given by where k is the bulk topological number (the TKNN number), and n ± counts the number of right/left moving edge states.
In two dimensions, both sides of (85) should be understood as the difference under sign flip of the parameter m for each cone. The topological number k is defined as the difference of the TKNN number when m changes its sign: As for the topological number of the gapless edge states, in our single fermion problem, we choose it as the sign of When it has a plus sign, we denote is as n + = 1, and when it has a minus sign, we denote it as n − = 1.
In this way, for all possible values of θ + , we show the bulk-edge correspondence. This means that the correspondence is true for any consistent boundary conditions in 2 dimensions.

B. 3D Weyl semimetals with the most general boundary conditions
Let us turn to the case of the 3 dimensions, the Weyl semimetals. In three dimensions, the topological number K is defined by the wrapping number of a map b i (p j ) which shows up in the Hamiltonian Our Hamiltonian (1) is given by b i = p i and the Weyl node is at p i = 0. Considering a two-sphere surrounding the Weyl node, we obtain We claim the bulk-edge correspondence for the 3D Weyl semimetal is given by where K is the topological number defined above. We define N andÑ to count the numbers of edge states with independent orientations with respect to the orientation of the bulk dispersion cone, as we will see below.
To discuss the orientation, we have to view the bulk dispersion in the subspace (p 1 , p 2 ) since the edge dispersion lives in that space. First, in the (p 1 , p 2 , ) space, we notice that all constant p 3 slices of the bulk states have the same orientation. Let us make further a slice at a constant positive energy . The cross-section of the bulk dispersion is a circle (see Fig. 11). The orientation of the circle is definite due to the topological number (assuming b 3 = p 3 ).
The constant energy slice of the edge dispersions defines the Fermi arcs. Since generically the edge state dispersions are planes tangent to the bulk dispersions, the Fermi arcs share the same property. The number N counts the number of Fermi arcs which are tangential to the bulk dispersion circle and emanates in a counterclockwise orientation. On the other hand, the numberÑ counts that in a clockwise orientation [22]. Our claim for the bulk-edge correspondence is that this orientation of the bulk circle remains the same for the edge (the Fermi arcs).
Let us check this explicitly for two typical examples. In Fig. 11 we show the Hamiltonian (1) with θ + = 0, and the case of the Hamiltonian H = −p 1 σ 1 + p 2 σ 2 + p 3 σ 3 with θ + = 0. The former case has K = 1 as explained before, while the latter case has K = −1. As we can see in Fig. 11, it is obvious that we have (N,Ñ ) = (1, 0) for the former case, and (N,Ñ ) = (0, 1) for the latter case. So, they are consistent with our claim of the bulk-edge correspondence (92).
All the edge states in Fig. 1 have the same N andÑ according to our definition: (N,Ñ ) = (1, 0). So they are consistent again with (92). The examples in the lattice models we considered are shown to be consistent with the bulk-edge correspondence. Note that Fermi arcs join Weyl nodes, and our counting works for each Weyl node.
To be more precise, each Fermi arc has two end points, and one end has (N,Ñ ) = (1, 0) while the other end has (N,Ñ ) = (0, 1). So the numbers are assigned to each end point of the Fermi arc.

V. NEW TOPOLOGICAL STRUCTURE AND BERRY PHASE
In Sec. II, we saw that the boundary condition of the 3D Weyl semimetals, as well as that of the 2D system, has only a single real parameter θ + . Since this θ + is a new parameter describing the system with a boundary, we can think of it as a coordinate in the theory space. Normally, for a given Hamiltonian, the theory space is spanned by parameters of the Hamiltonian. From a topological viewpoint, the parameters are conserved momenta. The 3D Weyl semimetals are of that category, and the parameter dependence of the Hamiltonian defines the bulk topological charge. Now, once we introduce the boundary, one of the momenta becomes ill-defined and drops off from the list of the parameters. However, there shows up a set of new parameters describing the boundary condition. From the analyses of this paper, we know that the new parameter is only θ + . So, the most general parameter space of the 3D Weyl semimetals with a boundary in the continuum limit is described by (p 1 , p 2 , θ + ).
To look for a novel topological structure of the system with a boundary, we study the wave function of the edge states, which depends only on the three parameters (p 1 , p 2 , θ + ). The nontrivial topological structure can often be detected by a Berry phase in the parameter space. We will find that, for the present case, the only nonvanishing Berry connection is that for θ + , and it provides us with a nontrivial Berry phase along a path in the parameter space.
Before getting to the details, we first note that the important part is just the phase of the wave function, to obtain a nonzero Berry connection. Suppose that the phase of the wave function does not depend on a parameter β. Then we find easily that Under this equality, the Berry connection associated with the parameter β is which means the vanishing of the Berry connection, A β = 0. Now, if we look at our general edge wave function (48), the phase does not depend on p 1 and p 2 . Therefore, we conclude in our generic parameter space.
However the phase of the wave function (48) depends on θ + , and we can calculate the Berry connection as Note here that the Berry connection for the edge state, which has x 3 dependence, needs in its definition the integral over x 3 space so that the connection becomes Hermitian. So, we have a nontrivial Berry connection along the parameter θ + , which dictates the most general boundary condition of the 3D Weyl semimetals.
With the non-vanishing Berry connection at hand, let us study a possible topological charge. The range of the parameter θ + is, as was analyzed earlier, the period 0 < θ + ≤ π. The point θ + = 0 is identical with the point θ + = π, so it describes a circle [23]. We consider a path going around this circle once, with some dependence on p 1 and p 2 . Let us calculate the Berry phase φ B along this path in the parameter space p 1 , p 2 and θ + . Using (95), we obtain This means that the edge state has a new topological charge, and its value is 1/2. In this calculation we considered a closed path in the (p 1 , p 2 , θ + ) space. Let us check whether the path exist or not. For the edge state to exist, we need to satisfy the normalizability condition α > 0. This amounts to a nontrivial relation α = p 1 sin 2θ + − p 2 cos 2θ + > 0.
In changing θ + from 0 to π, it is necessary to choose a path in the (p 1 , p 2 ) space so that the above inequality is satisfied. An example of such a path is given by (p 1 , p 2 ) = c (sin 2θ + , − cos 2θ + ) with a positive constant c. One may expect that the 2D case should have a similar topological number. Unfortunately, this is not the case. Since α has to be positive, a constant m gives a constraint on θ + : This has no solution for p 1 for a given m and all possible θ + . For example, for any positive m, at θ + = 0, there is no p 1 satisfying the inequality. In the same manner, for any negative m, at θ + = π/2, there is no p 1 . So, the constancy of m does not allow any path going from θ + = 0 to θ + = π.
If we consider a special case of m = 0, then we can find a path satisfying (98). An explicit example is p 1 = c sin 2θ + with a positive constant c.
We conclude that the edge states of the 3D Weyl semimetals acquire a new topological charge in the space of parameters of the boundary conditions. The 2D gapped case eliminates the topological charge, but 2D gapless systems can have the topological charge. The topological charge winds the space of θ + , the only parameter dictating the most generic boundary conditions.

VI. SUMMARY AND DISCUSSION
In this paper, we have studied the most general boundary conditions of the 3D Weyl semimetals in the continuum limit around the Weyl point. The boundary conditions are shown to be dictated by a single real parameter θ + which takes a value on a circle in the range 0 < 2θ + ≤ 2π. The edge state wave functions and their dispersion relations are obtained, and we find that the dispersion plane as a function of the remaining conserved momenta (p 1 , p 2 ) terminates at the bulk dispersion cone, as a tangential plane to it. The parameter θ + corresponds to the rotation angle of the edge dispersion plane relative to the bulk Weyl cone.
We build lattice Hamiltonian with a parameter at the boundary of the 3D square lattice, which reproduces the θ + dependence of the edge dispersion relations. Introduction of a boundary mass term at the edge of the lattice leads to various shape of the Fermi arcs joining Weyl nodes with opposite chiralities.
Through a dimensional reduction, the system becomes a 2D topological insulator of class A. The bulk-edge correspondence is found to be consistent for all values of θ + , meaning that for any generic boundary condition the bulk-edge correspondence works. We propose how to count the edge modes of 3D Weyl semimetals so that it becomes consistent with the correspondence between the number of the edge states and the bulk topological charges.
Furthermore, we discover a new topological number for the edge states of the 3D Weyl semimetals. The topological charge is associated with a Berry phase along the path parameterized by θ + , the new parameter dictating the whole boundary conditions. Various values of the new parameter θ + can be realized in experiments. For example, a hydrogen termination of graphene and related materials has been studied [21], and the dispersion relations obtained from the microscopic lattice Hamiltonians exhibit a behavior which we have generically studied in this paper. To what extent our θ + could be realized in experiments is one of the interesting future issues.
The meaning of the newly found topological charge needs to be verified in more details. The topological charge owned by the edge modes was studied in [24] for 4D topological insulators, inspired by a connection to superstring theory [25]. The edge topological charge would indicate existence of edge-of-edge states. For the present case, the topological charge is obtained by the change in θ + , that is, the boundary condition itself. Therefore, to have such an edge-of-edge state, one needs to change θ + as a function of the location on the boundary. It would be quite interesting if there exists such a topologically protected edge-of-edge state for 2D and 3D gapless systems.
In this paper our motivation came from a particletheoretical viewpoint to explore all possible boundary conditions in the continuum limit of the Weyl semimetals. As we summarized above, our approach turned out to give fruitful results in condensed matter physics. Bridges between condensed matter physics and particle physics will play crucial roles more in coming advances in physics.
Appendix A: Momentum-dependent generic boundary condition In this appendix, we study generic momentumdependent boundary conditions. We shall prove that the two properties of the edge state dispersions described at the end of Sec.II C are not modified even under a generic momentum dependence in the boundary conditions.
In the analysis so far, we have considered only a constant matrix M . This is because basically we are dealing with the small momentum limit where the Weyl point is approximated by a relativistic dispersion relation. Since the bulk Hamiltonian is linear in the momenta, We could allow a linear dependence also for the boundary Lagrangian, which means M linear in p. Since the boundary is at x 3 = 0, good quantum numbers are only p 1 and p 2 , so in M let us allow a linear dependence in p 1 and p 2 . Since these are just parameters, the constraint equation for M is obtained in the same manner, to have A 2 1 + A 2 2 − B 2 3 = 1. Although A 1 , A 2 and B 3 are linear in p 1 and p 2 , this constraint equation is quadratic. This means that the linear approximation of the boundary matrix M breaks down. One needs higher order terms in p 1 and p 2 to have a consistent boundary condition.
Let us formally continue the study of having a consistent momentum-dependent boundary condition and allow higher order terms of the momenta in M . Since we have shown above that generic boundary condition is completely dictated by the parameter θ + , this means that the θ + will depend on the momenta, θ + (p 1 , p 2 ). Thus, we have a function-parameter family of the boundary condition, whose edge state has a dispersion = −p 1 cos(2θ + (p 1 , p 2 )) − p 2 sin(2θ + (p 1 , p 2 )). (A1) Let us show that this generalized boundary condition also shares the same properties as that for the constant θ + . First, let us check that the edge dispersion intersects with the bulk dispersion. The latter is a disk for a given energy , Supposing that the boundary of the disk shares a point with the edge dispersion, we denote it as (p 1 , p 2 ) = (cos a, sin a).
The existence of a means that the edge dispersion intersects with the bulk dispersion. Substituting this expression to the edge state dispersion (A1), we find an equation for a, cos(a − 2θ + ( cos a, sin a)) = −1.
Since the left hand side of this equation is a periodic function of a, we can show that this equation always have odd number of solutions for a in the period 0 ≤ a < 2π. Thus the edge dispersion always intersects with the bulk dispersion.
Next, let us show that the edge dispersion is always tangential to the Weyl cone. At the intersection point, we have a solution a = a 0 of the equation (A4). At a = a 0 , the tangential line of the bulk dispersion at a constant has a slope − cot a 0 in the (p 1 , p 2 ) space, since the Weyl cone at a slice of constant is just the disk. Let us show that the edge dispersion (A1) also has the same slope at the intersection. Differentiating the dispersion relation (A1) with respect to p 2 , we find dp 1 dp 2 = − sin 2θ + cos 2θ + + sin 2θ + cos 2θ + p 1 − p 2 d(2θ + ) dp 2 . (A6) By using the intersection condition (A4) and (A3), we can simplify this equation to dp 2 dp 1 = − cot a 0 (A7) which shows nothing but the slope − cot a 0 in the (p 1 , p 2 ) space. This completes the proof that the edge dispersion is tangential to the bulk dispersion. Finally, we show that the intersection point is in fact the point where the edge dispersion ends. To show it, we notice that the termination condition of the edge dispersion is α = 0, since the existence of the edge is certified by α > 0. Now, for the momentum-dependent θ + , we have the same formula α = p 1 sin (2θ + (p 1 , p 2 )) − p 2 cos (2θ + (p 1 , p 2 )) . This means that the edge dispersion is actually terminated at the intersection point, thus the edge is absorbed into the bulk.
In this manner, we can show that, even when a completely generic momentum dependence is allowed in the boundary condition, the edge dispersion keeps the properties that it is tangential to the bulk Weyl cone and it is terminated there.

Appendix B: Two parallel boundaries
In the main part of this paper, we study the general case of a single flat surface as a boundary. For realistic materials, two parallel boundaries are typical, and in this appendix we analyze the Weyl semimetal with two parallel boundaries, x 3 = 0 and x 3 = L, in the continuum limit.
Each boundary can have the parameter θ + , so the number of the parameters grows as one introduces many boundaries. In this subsection, just for simplicity, we consider the case with identical values of θ + for the two boundaries.
Since the bulk Hamiltonian is not altered, the energy eigen equation (35) does not change, thus a generic solution is with (37). Note here that we need to include another mode exp[+αx 3 ] because the allowed region of x 3 is a finite period 0 ≤ x 2 ≤ L. Now, (35) leads to On the other hand, the boundary condition at x 3 = 0 and x 3 = L leads to (1 e −2iθ+ ) ξ 1 + ξ 2 η 1 + η 2 = 0, (1 e −2iθ+ ) e −αL ξ 1 + e αL ξ 2 e −αL η 1 + e αL η 2 = 0 (B5) which are equivalent to (1 e −2iθ+ ) So, we can solve the wave function as if two boundaries are just a set of copies of a single boundary, thanks to our assumption that the same boundary conditions are shared by the two. Together with (B2) and (B3), we have zero-mode equations However this equation has no solution for generic p 1 and p 2 . Therefore, either (B7) or (B8) is solved by a trivial vanishing solution. When (B7) has a nontrivial solution, ξ 2 = η 2 = 0, then the wave function and the energy eigenvalue is completely identical to the previous case of a single boundary: = −p 1 cos 2θ + − p 2 sin 2θ + , (B10) α = p 1 sin 2θ + − p 2 cos 2θ + .
(B11) This is the edge state dispersion, which turn out to be identical to the single boundary case, (44). On the other hand, when (B8) has a nontrivial solution, then ξ 1 = η 1 = 0 and we obtain another edge state = −p 1 cos 2θ + − p 2 sin 2θ + , (B12) α = −p 1 sin 2θ + + p 2 cos 2θ + . (B13) This expression differs by just a sign of α, compared to the previous one. Furthermore, since the front factor is exp[αx 3 ] which differs only by the sign of α compared to the previous one, this solution turns out to be identical to the previous one. So, we conclude that we have a single edge state given by (B10) and (B11). The only difference compared to the single boundary case is that now there is no restriction α > 0. In fact, depending on the momenta (p 1 , p 2 ), α can be positive or negative, and depending on it, whether the wave function localizes at x 3 = 0 or x 3 = L will be determined. Thus the obtained edge state, although written in a unified form, represents both the modes, one localized at x 3 = 0 and the other at x 3 = L.