Interference of the Bloch phase in layered materials with stacking shift

# Interference of the Bloch phase in layered materials with stacking shift

## Abstract

In periodic systems, electronic wave functions of the eigenstates exhibit the periodically modulated Bloch phases and are characterized by their wave numbers . We theoretically address the effects of the Bloch phase in general layered materials with stacking shift. When the interlayer shift and the Bloch wave vector satisfy certain conditions, interlayer transitions of electrons are prohibited by the interference of the Bloch phase. We specify the manifolds in the space where the hybridization of the Bloch states between the layers is suppressed in accord with the stacking shift. These manifolds, named stacking-adapted interference manifolds (SAIM), are obviously applicable to general layered materials regardless of detailed atomic configuration within the unit cell. We demonstrate the robustness and usefulness of the SAIM with first-principles calculations for layered boron nitride, transition-metal dichalcogenide, graphite, and black phosphorus. We also apply the SAIM to general three-dimensional crystals to derive special k-point paths for the respective Bravais lattices, along which the Bloch-phase interference strongly suppresses the band dispersion. Our theory provides a general and novel view on the anisotropic electronic kinetics intrinsic to the periodic-lattice structure.

1

## I Introduction

When an electron moves in the condensed matter, interplay of its wave property and atomic configuration often yields nontrivial interference effects, which impose constraints on the electronic dynamics. A celebrated example of such interference effects is the flat band Tasaki-review-PTP1998 (). Despite nonzero hopping integrals, spatially confined electronic eigenstates emerge in exotic structures such as delta-chain, decorated square, and kagome lattices [Fig. 1(a)]. There, possible intriguing phenomena due to electronic correlation have been discussed extensively: The flat-band ferromagnetism Lieb-PRL1989 (); Mielke-JPhysA-1 (); Mielke-JPhysA-2 (); Tasaki-PRL1992 (); Mielke-Tasaki-MathPhys1993 (), fractional quantum hall effect Bergholtz-Liu-review2013 (), etc. Theoretical formulation of the interference effects thus provide us with guiding principles to control the kinetics of electrons, which should facilitate exploration of novel electronic phenomena emerging from the competition of kinetic and local-correlation effects. In this work, we discuss a general effect of interference of the Bloch wave function governing the one-body characteristics of electrons.

The electronic wave function of the eigenstate in the periodic potential has the form Seitz-grouptheory1936 ()

 ψnk(r)=exp(ik⋅r)unk(r), (1)

with being the crystal wavenumber defined in the Brillouin zone (BZ). denotes the lattice-periodic part which satisfies for any lattice vectors . Since the plane-wave part, or the Bloch phase [Fig. 1(b)], does only concern the periodic structure of the system, this should be utilizable to devise interference effects applicable to a wide range of materials.

In this study, we theoretically address such possible interference effects in general layered structures. Theoretical outcomes of our study are directly applicable to recently synthesized atomically thin materials such as graphene Geim-Novoselov-graphene2004 (), boron nitride (BN) Geim-Novoselov-thinlayer (), transition-metal dichalcogenides (TMD) Geim-Novoselov-thinlayer (); Park-Lee-1Tprime-MoTe2 (); Zhou-Dresselhaus-1Tprime-MoTe2 (); Naylor-Johnson-1Tprime-MoTe2-monolayer () and phosphorene Li-Zhang-phosphorene-transistor-nnano (); Liu-Tomanek-phosphorene-mobility-acsnano (). Toward development of novel nanoelectronics devices, numerous experimental and theoretical investigations have been carried out on their electronic properties depending on the number and stacking geometry of the layers CastroNeto-graphene-review2009 (); Xu-Heinz-review-NPhys2014 (); Qian-review-2DMater (). Generally the neighboring layers interact with each other through various mechanisms; electrostatic force, the van der Waals force, electronic interlayer hybridization and charge redistribution Tomanek-DMC-phosphorene (), etc. Here we focus on the hybridization from the viewpoint of the interference of the wave functions. When two layers are placed along with each other, their Bloch states hybridize between the layers. Because of the phase of the wave function, this hybridization is prohibited for some combination of and the relative geometry of the layers [Fig. 1(c)]. This effect should impose k-point and stacking-geometry dependent constraints for the interlayer interaction of the electronic states. As a matter of fact, in graphite, it has been known that the interlayer hybridization at special points of the BZ is much affected by the stacking pattern McClure-model-graphite1957 (); Slonczewski-Weiss-model-graphite1958 (); McClure-rhombo-graphite1969 (); Charlier-Gonze-multilayer-graphite1994 (). By carefully examining electronic band structures in the literature, one can find the quasi-two-dimensional states due to this effect in a variety of works not only on graphite Kresse-graphite-BN-band1996 (); Coehoorn-Wold-TMDC-bandstruct-PES (); Dawson-Bullett-TMTe2 (); McClure-model-graphite1957 (); Slonczewski-Weiss-model-graphite1958 (); McClure-rhombo-graphite1969 (); Charlier-Gonze-multilayer-graphite1994 (); Partoens-multilayer2006 (); MacDonald-graphene-chiraldecomp (); Chiu-Shyu-r-graphite (); Liu-Shen-BN-stack (); Aoki-Amawashi-SSComm (). In addition, it has recently been experimentally demonstrated that this interference can be utilized for controlling the interlayer motion of electrons and holes by the stacking geometry Suzuki2014 (); Akashi2015 (). However, these features have been discussed in distinct material-specific contexts. Our aim in this study is to establish a general theory which not only relates the above-mentioned phenomena but also applies to any layered materials. Such theory should give us a unified insight into how the electronic states in the layered structures are controllable, which will also facilitate efficient nanoelectronics materials design.

In this paper, we present a theory how the Bloch phase governs the possible/impossible interlayer hybridization, which respects only the space group and stacking shift of the layers. For theoretical simplicity, we concentrate on the stacking shift and leave the stacking twist beyond our scope. Still, we will see below that the applicability of our theory is broad.

We specifically formulate special manifolds in the two-dimensional reciprocal space as functions of the stacking shift [Fig. 1 (c)]—the Bloch states in those manifolds cannot hybridize between the layers. These manifolds are determined regardless of the cell periodic part of the wavefunction and summarized as the stacking-adapted interference manifolds (SAIM). The characters of the Bloch states in any of the SAIM can be simply and sensitively controlled by the stacking shift between the layers via the change of the interlayer hybridization. To demonstrate this, we calculate the electronic states in layered BN, TMD and graphite from the first principles and relate their properties to the interlayer stacking shift. We show that, once the electronic properties within the SAIM is known for the monolayer form of these materials, we can predict their stacking dependences with only the information of the interlayer shift vectors. Our results show the usefulness of the SAIM for the control of the electronic states in the layered materials.

In addition to the direct application, we believe that the SAIM is a useful concept to derive various non-trivial electronic properties in general periodic structures. We exemplify this point by decomposing general three-dimensional lattices into stacked layers. It is consequently revealed that special k-point paths in the BZ are present where the band dispersion is anomalously suppressed regardless of the atomic configuration within the unit cell. We show the robustness of these paths by examining the first-principles calculations for semiconducting lithium hydride (LiH) and sodium chloride (NaCl) with rocksalt-type crystal structure. These paths, which we name the Bloch-phase induced flat-band paths (BIFP), indicate hidden anisotropic characteristics of the Bloch states intrinsic to the three-dimensional lattice structure.

This paper is organized as follows. In Sec. II, we discuss the Bloch-phase interference in stacked one-dimensional chain as a simple example for later generalization. The general theory for arbitrary layers and derivation of the SAIM are presented in Sec. III. The first-principles electronic structure calculations for BN, TMD, graphite and black phosphorus and their interpretation based on the SAIM are presented in Sec. IV. The consideration on the general three-dimensional lattice is exhibited in Sec. V, where the BIFP is introduced through an analysis on the tight-binding model on the face-centered cubic (FCC) lattice and the first-principles band dispersions along the BIFP are examined for LiH and NaCl. Miscellaneous topics are discussed in Sec. VI. Section VII is devoted to summary and conclusions of the present work. The full lists of the SAIM and BIFP are provided in Appendices B and E.

## Ii Introductory example: one-dimensional chain

We exemplify the present concept with an -wave single orbital tight-binding model on a one-dimensional chain [Fig. 2 (a)] with periodic boundary condition. The primitive lattice vector is denoted by . The electronic eigenstates of this system is characterized by the one-dimensional Bloch wavenumber .

Next, we put an identical chain in parallel with the original one with finite shift [Fig. 2(b)]. Since the primitive lattice vector of the total system remains , inter-chain hybridization occurs only between the Bloch states having the same . Although there seems to be no general reason for prohibiting such hybridization, it is not the case for : When , the hybridization between the states of the two systems is suppressed. This is illustrated in Fig. 2(b): the hopping amplitudes from the two neighboring sites of the other chain cancel with each other due to their opposite phases of the wave function.

Let us discuss this fact more formally. For the two-chain system with shift vector , we define the Bloch eigenstates for the th chain as (=1, 2). The hybridization amplitude between the Bloch states of the respective chains is then defined by , where denotes the total Hamiltonian including the hopping between the sites belonging to different chains [Fig. 2(b)]. Here and hereafter, we use the terms “hopping” and “hybridization” to describe the transition between the atomic (or molecular) orbitals and Bloch states, respectively. The individual chain is invariant under mirror reflection with respect to planes perpendicular to the chain . Note that there are two inequivalent types of planes which cannot be connected with the lattice vector: ones crossing the sites (plane 1) and ones crossing between the sites (plane 2) [See Fig. 2(b)]. When , plane 1 for the one chain and plane 2 for the other chain coincide with each other and therefore the total Hamiltonian retains the global mirror symmetry [Fig. 2(b)]. We then get

 tπ(a/2) = ⟨1,π|H(a/2)|2,π⟩ (2) = ⟨1,π|mmH(a/2)mm|2,π⟩ = [⟨1,π|m][mH(a/2)m][m|2,π⟩] = −tπ(a/2),

where we have used and and . The last line forces to be zero.

Thus, the inter-chain hybridization between the Bloch states are prohibited when . The key factors in the proof are the following: (i) There are two types of mirror planes which cannot be connected with the lattice vector, (ii) the state is the eigenstate of the mirror reflection, whose eigenvalue depends on the planes, and (iii) when , the common mirror planes correspond to plane 1 and plane 2 for the first and second chains, respectively. We next extend the present discussion for general two-dimensional crystals.

## Iii General theory on interlayer hybridization

We discuss the two-dimensional periodic system based on its classification by the net, which is a line graph connecting the points generated by the two-dimensional primitive translation vectors (see Fig. 3(a) for example) and therefore the two-dimensional variant of the Bravais lattice footnote-net-origin (). Every periodic two-dimensional crystal structure corresponds to either of the five nets: square, hexagonal, diamond, rectangular, or oblique nets. The definitions of the nets are provided in Appendix B.1 for convenience. Note that the concept of net is also applicable to when the atoms are arranged within a nonzero range in the direction perpendicular to the net. We also employ the “in-plane” symmetry operation ; namely, the operation on the two-dimensional space spanned by the net []. In the discussions below, the set of the in-plane symmetry is limited to the orthogonal operation, which keep the vector inner product invariant. There is a one-to-one correspondence between the possible set of compatible with the crystal structure and the corresponding net; for example, only the rotation by () is compatible with the oblique net.

Let us next consider a two-dimensional periodic system which has a certain in-plane symmetry . We do not address the one-dimensional case explicitly since it is trivial footnote-1d (). For some combinations of and net, there are two or more possible reference axis/planes in the unit cell, which cannot be connected to each other by the primitive translation vectors. Define the fractional translation connecting these references by , where represents the shift vector. Labeling the symmetry operations with respect to different references by , these elements are related by

 S2 = Tτ1S1T−1τ1, S3 = Tτ2S1T−1τ2, ⋯ , Sn = Tτn−1S1T−1τn−1. (3)

An example of the hexagonal net with = is illustrated in Fig 3(a). This net is formed by the lattice vectors and (red or dark arrows) and there are three possible rotation axes (filled triangles) in the unit cell. These axes are connected by the vectors and (cyan or light arrows).

Here we introduce an electronic one-body Hamiltonian for the two-dimensional periodic layer, which is characterized by real-space ionic potential. Regardless of the form of the potential, the corresponding eigenstates, the Bloch states, take the form , where denotes the cell-periodic part composed of the atomic orbitals. Note that the wave functions introduced in Eq. (1) are given by and , respectively. Later in this subsection, we focus only one non-degenerate band and omit the band index ( and ). When is an eigenfunction of ,

 Sjψk(r) = exp[i(φpwj(r)+φcellj(r))]ψk(r) ≡ exp(iφj)ψk(r), (j=1,2,…,n)

where denotes the contribution of the plane-wave part () and that of the cell-periodic part (). The following formula is then derived (see Appendix A for detail)

 φpwj(r)−φpw1(r) = (k−Sk)⋅τj−1, (5) φcellj(r)−φcell1(r) = 0, (6)

and consequently we get

 Sjψk(r) = Tτj−1S1T−1τj−1ψk(r) (7) = exp[i(k−Sk)⋅τj−1+iφ1]ψk(r).

In the BZ, we then define the stacking-adapted interference manifolds (SAIM) by the following:

 k∈Kin(S;τ)⇔ (k−Sk)⋅τ≠0 mod 2π (8) and  k−Sk=n1b1+n2b2(∃n1,n2∈Z), (9)

where and are the reciprocal vectors. The condition Eq. (9) is required by the assumption that is an eigenstate of . With being a vector connecting different references of , is, if nonempty, located at the boundary of the BZ inui (). Figure 3(b) illustrates [=] for the hexagonal net.

Finally, we stack an identical layer with shift and label the Bloch states of the th layer as [See Fig. 3(c), e.g.]. A short summary of the consequence is the following: The hybridization between the Bloch states of different layers are canceled for due to the interference of the Bloch phase.

The total one-body Hamiltonian has the form . Here, and denote the Hamiltonian for the first and second layers in isolated forms and and are the layer-diagonal potential terms due to the presence of the neighboring layers. represents all the hopping terms between the layers. Since the lattice translation vectors of the total and individual systems are common, regardless of , the total Hamiltonian and the interlayer hopping term are block-diagonalized for each . Using as a basis set, the matrix representation of is

 H(τ)=(⟨1,k|H(τ)|1,k⟩⟨1,k|H(τ)|2,k⟩⟨2,k|H(τ)|1,k⟩⟨2,k|H(τ)|2,k⟩). (10)

When , the total Hamiltonian is invariant under global symmetry operation , which corresponds to and for the first and second layer, respectively. When , for the non-diagonal element, we get

 ⟨1,k|H(τ)|2,k⟩ = ⟨1,k|h12(τ)|2,k⟩ (11) = ⟨1,k|S−1Sh12(τ)S−1S|2,k⟩ = [⟨1,k|S−11]h12(τ)[T−τS1T−1−τ|2,k⟩] = e−i(k−Sk)⋅τ⟨1,k|h12(τ)|2,k⟩.

Here we have used Eqs. (III.1) and (7). For , and the only allowed value for is zero. The prefactor comes solely from the plane-wave part of the wave function (see Eqs.(5) and (6)) and therefore represents the interference of the Bloch phase comment-v12 ().

The manifold is uniquely determined according to the kind of the net, symmetry element and the shift vector . For example, we show for the rectangular net in Fig. 4. The lattice and reciprocal vectors and the corresponding BZ are shown in panel (a). The possible symmetry elements compatible with this net are the two-fold rotation () and mirror reflections with respect to the and planes ( and , respectively). When the whole crystal is symmetric, there are four inequivalent axes in the unit cell and there are three vectors connecting them [Fig. 4(b)]. Elementary calculations give us the following:

 Kin(C2;τ1) = {k|  k=b1/2,  b1/2+b2/2}, (12) Kin(C2;τ2) = {k|  k=b1/2,  b2/2}, (13) Kin(C2;τ3) = {k|  k=b2/2,  b1/2+b2/2}. (14)

When the crystal is symmetric, on the other hand, there are two inequivalent mirror planes in the unit cell [Fig. 4(c)]. The vector connecting them can be taken arbitrarily with respect to the component parallel to the planes: with . Regardless of , we get

 Kin(mx;τ) = {k|k=b1/2+tb2} (−1/2≤t≤1/2). (15)

A similar discussion applies to . With ,

 Kin(my;τ) = {k|k=tb1+b2/2} (−1/2≤t≤1/2). (16)

The complete list of non-empty for all the nets is available in Appendix B.

We thus introduced a concept of SAIM for the respective nets and established its relation to interlayer hybridization of the Bloch states, which can be summarized as follows. Suppose a single-band monolayer system which is invariant under in-plane symmetry operation . When an identical layer is stacked side by side, a theorem on interlayer hybridization is derived. Namely,

 ⟨1,k|H(τ)|2,k⟩=0 (k∈Kin(S;τ)) (17)

where and denote the in-plane shift between the adjacent layers and the total Hamiltonian, respectively. If there are two or more symmetry references for within the unit cell, is non-empty with set to connect the different symmetry references. Note that this theorem does not concern the character of the cell-periodic parts of the Bloch states.

#### Layers related by symmetry operation

Frequently, we find multilayer structures where the neighboring layers are not exactly identical but related with each other by some symmetry operation ; namely, inversion, reflection, etc.. Extending the present theory to such cases is in principle possible if the corresponding symmetry operation does not change the shape and orientation of the two-dimensional BZ. For example, three-dimensional hexagonal BN consists of two kinds of BN layers which are the inverse of the others. Although the atomic configurations are different for the two adjacent layers, their BZs are common. We here address an extension to such systems to enhance the applicability of the present theory.

First, we need to redefine the “shift” between the layers. Namely, we set the positions of the nets for the two layers so that their nodes correspond to the centers of the cell-periodic part of the Bloch states . The in-plane shift is defined with these nets. Next, we see if is transformed by by either of the following:

 kΣ→ k (18) kΣ→ −k. (19)

In the former case, the character of for the first and second layers are the same, whereas in the latter for the second layer is complex conjugate of the first since they are related by , with being the time reversal. Equation (17) is obviously applicable if either of the following is satisfied: (I) is transformed by Eq. (18), or (II) is transformed by Eq. (19) and is a real representation of . Otherwise, whether the interlayer hybridization is zero or nonzero depends on the specific character of . In summary, the extension is achieved with only two additional steps for the cell-periodic part: (a) Redefine the in-plane shift between the layers based on the positions of the orbital and (b) examine if it is real representation of .

When we apply the present results to realistic models based on the crystal structure of the material, step (b) above requires a burdensome calculation to specify the character of . To facilitate efficient applications, we here exemplify situations where Eq. (17) is obviously applicable: When is the mirror plane parallel to the layers; when neither nor is compatible with the net; is formed by the atomic orbital with zero orbital angular momentum with respect to the axis perpendicular to the layers (, , , ); etc.

### iii.2 Multilayer and bulk systems

When a number of identical or symmetry-related layers are stacked so that the shape of their corresponding two-dimensional BZs are the same, the total Hamiltonian remains block-diagonalized with respect to the crystal wavenumber defined for the original layers. By arranging the interlayer shift vectors, desired combinations of hybridizations among the layers can be suppressed for .

Next, suppose an infinite number of the layers are stacked so that a certain satisfies for all the shift vectors between the adjacent layers. The eigenstates formed by are then labeled by the crystal wavenumber perpendicular to the layers . However, their energies are almost independent of because of the absence of the interlayer hybridization; namely, the flat band dispersion is seen along from . Note that there is still a possibility of interlayer hybridization between distant layers. For example, when the hexagonal nets are stacked with uniform shift ([Fig. 3(a)], the -nearest neighbor layers () are regarded as unshifted from the original layer. The hybridization amplitudes among , and can then be nonzero, respectively. In summary, the band dispersion in direction tends to be much flatter through than through other points because of the canceled interlayer hybridizations.

### iii.3 Multiband and degeneracy

In the above discussions, we have concentrated on the intra-band hybridization by limiting ourselves to the single band case. For general mutiband systems, interlayer hybridizations between different orbitals are possible even for , if their orbital characters are different. Nevertheless, we can ignore such inter-orbital hybridizations when the respective bands in focus are non-degenerate and well separated from other bands in energy at in the monolayer system. This point is elaborated in Appendix C. Finally, the following rule of thumb is established: In general multilayer systems, the Bloch states at for each layer are well protected from interlayer hybridizations if only their energy eigenvalues are nondegenerate and well separated from other states. In Sec. IV, we demonstrate that this seemingly naive rule of thumb indeed applies to a wide range of multilayer systems.

### iii.4 Implication of stacking-adapted interference manifolds

We have seen that the interlayer hybridization is generally canceled by the Bloch phase for within . This is valid even in the multiband case for the intraband hybridization. These results imply a notable fact: When we examine the electronic structure of layered systems with changing the interlayer shift, the Bloch states around with any and always exhibit noticeable variation. This special feature provides us with a guiding principle to explore the shift-dependent electronic property of layered systems, that we should examine the electronic states around . For the other points, the shift-vector dependence crucially differs by systems and its general features are probably difficult to derive. For example, in the single-band system formed by -wave atomic orbitals, the interlayer hybridization obviously varies little by changing the interlayer shift for far from . This is reasonably inferred since the overlap integral of the cell-periodic parts of the Bloch functions in the neighboring layers have positive definite values regardless of the shift. In this case, one finds the appreciable shift-vector dependence only in the vicinity of with any .

## Iv Applications to layered semiconductors

In the previous section, we have presented a theory relating the interlayer hybridization of the Bloch states and the in-plane shift between adjacent layers based on a general modeling of layered materials. The significant consequence is that whether the hybridization is canceled or not can be discussed only by examining the relative shift of adjacent layers. We next show that our theory indeed applies to real materials. First-principles electronic structure calculation is employed for this purpose, where the inter-site hopping is treated without modeling. We take for example the layered BN, graphite, group-VI TMD and black phosphorus. The structures of the former three in the monolayer form commonly correspond to the hexagonal net. They form various multilayer structures with shifts so that the axes for the neighboring layers are shared, and then corresponds to the corners of the hexagonal BZ, and [Fig. 3 (b)]. On the other hand, black phosphorus monolayer (phosphorene) forms in the rectangular net structure. According to the in-plane symmetry and interlayer shift, possibly corresponds to the either edges of the BZ(Fig. 4).

We first examine the electronic states in bilayer BN with various stacking geometries, showing that whether the interlayer hybridization is suppressed or not is correctly predicted by simply seeing how the two layers are shifted from each other. Next, anomalously flat band dispersion predicted from the present theory is examined in the bulk BN. We also address molybdenum disulfide MoS. Interlayer hybridization in multilayer MoS has been recently discussed in the literature Liu-review2015 (); Akashi2015 (), which is actually an example of the present general theorem. Fi, we work on layered graphite. The two-fold degeneracy at the and points in the monolayer form defies straightforward application of the present theorem as discussed in Sec. III.3. Nevertheless, we see that a minimal extension enables us simple geometric analysis on the interlayer hybridization. The possible interlayer hybridizations and the band structure along the line (along the direction passing through the point) for these systems have been discussed before in various preceding works, though in separate contexts. Our theory rederives results consistent with them from a unified viewpoint. Finally, through the comparison of the electronic structures in the black phosphorus with the experimentally observed shifted stacking and pathological unshifted stacking, we see that our theory is well applicable to the structure other than the hexagonal one. The detailed conditions for the first-principles calculations are given in Appendix  D.

### iv.1 Boron nitrides

#### Bilayer

Here we investigate the bilayers of BN with three types of stacking labeled , and , which are depicted in Fig. 5(c). corresponds to the -BN bilayer; type corresponds to the neighboring two layers in the bulk rhombohedral BN; corresponds to the -BN bilayer with shift, where nitrogen atoms in the upper layer are placed over those in the lower. The former two have been experimentally realized Warner-Buchner-BNbilayer-topography2010 (), whereas the latter stacking has been studied only theoretically and found to be unstable Ooi-Adams-variousBN2006 (); Constantinescu-Kuc-BNstacking-PRL2013 ().

The crystal and electronic band structures of the monolayer BN are displayed in Fig. 5(a). The crystal structure belongs to the hexagonal net, where bold dashed arrows indicate the primitive lattice vectors. The monolayer is -symmetric and therefore there are two additional axes in the unit cell indicated by and . The valence-top and conduction-bottom states at the and points are well separated in energy from other states, which justifies the neglect of inter-orbital hybridization between the layers.

The above pieces of information are enough for analyzing the interlayer hybridization in the stacking because the two layers are related only by translation. The in-plane shift between the two layers yields and therefore interlayer hybridization amplitudes must be respectively zero for the conduction-bottom and valence-top states at the and points. This is confirmed by observing the wave functions of the -point valence-top and conduction-bottom states: Due to the inequivalence of the electrostatic potentials for the states in the upper and lower layers, small band splitting is seen in the valence-top and conduction-bottom states. As clearly seen in Fig. 5(d), the calculated wave functions for these states (1–4) are localized around each layer, indicating the zero hybridization between them.

The analysis for the and stackings requires additional information for the monolayer band structure because the second layer is an inverse of the first. Namely, the conduction-bottom and valence-top states at the and points are formed by the B- and N- orbitals, respectively. The present cases then correspond to case (II) in Sec. III.1, which allows us to rely again on the simple geometry-based argument regardless of the specific orbital characters.

Introducing the hexagonal nets for the B and N sites, respectively, we see that in the stacking the B-net (N-net) of the upper layer is shifted from that of the lower layer by (). This shift yields zero hybridizations at the and points for the both states. This is immediately confirmed by the calculated band structure: The and -point valence-top and conduction-bottom states are doubly degenerate [Fig. 5(c)]. Since the degeneracy is lifted when either the potential difference between the layers or the interlayer hybridization is nonzero, the degeneracy indicates the zero interlayer hybridization amplitudes. In the case, on the other hand, the B-nets are shifted from each other whereas the N-nets are not. Then, only the conduction-bottom states should be protected from the interlayer hybridization. This is also confirmed with the calculated band structure displayed in panel (c), where only the valence-top states split. The wave functions of the valence-top states are equally distributed to the two layers, which indicates the formation of the bonding/antibonding states due to the interlayer hybridization. The conduction-bottom states are doubly degenerate, indicating the zero interlayer hybridization.

#### Bulk

We have thus demonstrated the relevance of our theory relating the stacking geometry and interlayer hybridization with the bilayer systems. Next, we proceed to the infinite-layer case. Although the band structure in this system has been studied Kresse-graphite-BN-band1996 (), we revisit it in view of the present theory, deriving which hybridization processes dominate the bandwidth for various stacking geometries.

We address the three types of structure: (i) the rhombohedral BN, where the monolayers are infinitely stacked with constant shift [Fig. 6(b)]; (ii) the hexagonal BN, where the layers are -stacked (See Figs 5(c) and 6(e)); (iii) the AA BN, where the layers are stacked without shift or inversion [Fig. 6(h)]. For the - and -point states, the rules of interlayer hybridization established in Sec. IV.1.1 straightforwardly yield the complete list of possible/impossible hybridizations for all the three stackings. In the rhombohedral BN, hybridizations to nearest-neighbor and next-nearest-neighbor layers are prohibited because of the relative shifts by and , respectively. The most relevant hybridization allowed is therefore those between the third-nearest neighbors [Fig. 6(b)]. In the hexagonal BN, on the other hand, the hybridizations between the next-nearest neighbor layers are allowed because they are unshifted [Fig. 6(e)]. Finally, in the AA BN, hybridization between the neighboring layers is obviously allowed [Fig. 6(h)]. We depict the geometrically allowed hybridizations for the three stackings by arrows in Fig. 6(c), (f), and (h). Therefore, it is anticipated that the bandwidth along the line has the smallest value in the rhombohedral structure, followed in order by the hexagonal and AA. The first-principles band structures indeed show this trend [Fig. 6(d), (g), (j)].

We also display the calculated band structures in the path in a broader energy range in Fig. 7. Apart from the Brillouin-Zone folding with the hexagonal stacking, the band structures for the three stackings are very similar in the path. Significant differences are seen only in the path, which indicates that the phase interference effect due to the stacking geometry is relevant only around the -point states. Interestingly, in the path, the other valence bands are also flat in the rhombohedral and hexagonal cases, demonstrating the generality of the present hybridization rule.

### iv.2 Transition-metal dichalcogenides

The TMD (=Mo, W; =S, Se, Te) mono and multilayer have been recently much in focus for its valley-contrasted phenomena Xu-Heinz-review-NPhys2014 (). Many features of TMD are common to the BN in their monolayer form: Seen from above, the former structure is identical to the latter. The only difference is that the atoms are placed in a non-zero range with respect to the direction in TMD (Figs. 8 (a), (b), and (c)). The valence-band top and conduction-band bottom are located at the and points of the hexagonal BZ. Most studies on the TMD multilayer concern either the 2H or 3R stacking depicted in Figs. 8 (d) and (e). In prior to this work, for MoS, some of the present authors have found that the interlayer hybridization of the conduction-band bottom states is zero with both the 2H and 3R stackings, whereas that for the valence-band top states is zero with only the 3R stacking Suzuki2014 () (Figs. 8 (d) and (e)). This has been proved to be the consequence of the phase interference Liu-review2015 (); Akashi2015 (). We revisit this topic to see that the previous results are re-derived in a simplified fashion with the present theory.

First, the monolayer structure belongs to the hexagonal net. The crystal is symmetric and therefore there are two possible shift vectors and [Fig. 8 (c)] which yield . The valence-top and conduction-bottom states are non-degenerate and well separated from other bands. These pieces of information are enough for the 3R stacking. We straightforwardly get to the consequence that the interlayer hybridization is prohibited for both the valence-top and conduction-bottom states. For the 2H stacking, where the upper layer is inverted from the lower, we need additional information: The cell-periodic parts of the valence-band top and conduction-band bottom Bloch states are formed by Mo orbital with angular momentum and orbital with , respectively Xiao-MoS2-PRL (). The case (II) in Sec. III.1 applies only for the conduction-bottom. We next introduce a hexagonal net for the Mo sites. In the 2H stacking, the Mo-nets are shifted from each other by [Fig. 8(e)]. The argument based on the geometry immediately yields zero hybridization for the conduction-bottom states comment-valence (). The present theory thus gives a simple way of understanding the possible/impossible interlayer hybridizations of the valley states of the TMD.

### iv.3 Graphite

In the band structure of graphene monolayer, the celebrated Dirac cone emerges at the and points. Since the states at the Dirac cone are two-fold degenerate, interlayer hybridizations among all of the degenerate states inevitably become relevant. Still, we can apply the analysis based on the layer geometries. The key property is that the Dirac-cone states are formed by the C- orbitals at the inequivalent C sites [C and C in Fig. 9(a)]. We then introduce the auxiliary hexagonal nets for the C and C sites, respectively [Fig. 9(a)]. Once we know that the two orbitals both belong to the same representation of , we can later forget about the specific character of the orbitals and focus only on the geometric relations between the C- and C-nets. The analysis below corresponds to a simplified counterpart of Charlier et al. (Ref. Charlier-Gonze-multilayer-graphite1994, )

Let us consider the stacked bilayer of graphene with shift as depicted in Fig. 9(b). Extracting the nets of the two layers, we see that the C-net of the upper layer is stacked onto the C-net of the lower layer without shift, whereas the shift vectors for the other combinations are either or . Therefore, the only allowed interlayer hybridization for the - and -point states is that between the C state of the upper and C state of the lower [Fig. 9(c)]. Note that for the bilayer stacked without shift, we straightforwardly get to the rule that the interlayer hybridizations across the C and C orbitals are prohibited [Fig. 9(d)].

On the basis of the rules derived above, we predict the band structures of bulk graphite in the path with three types of stacking: Rhombohedral (with shift , , ), Bernal (, , ), and AA (without shift) stacking [Fig. 10(a), (d), (g)].

Let us begin with the rhombohedral stacking. We first depict the allowed hybridizations between the neighboring layers by thick arrows in Fig. 10(b) based on Fig. 9(c). These hybridizations yield localized bonding and antibonding states. A dispersive band should not be formed, however, since there is no channel passing through the whole crystal. Next, we add the allowed hybridizations between the next-nearest neighbor layers as thin arrows based on Fig. 9(c). With these hybridizations, channels through the crystal are formed. Therefore, there should emerge bands whose bandwidths are governed by the second-nearest neighbor hybridization. The prediction is confirmed by the first-principle calculation as depicted in Fig. 10(c). The inter-band splitting reflects the bonding/antibonding nature dominated by the nearest-neighbor hybridization.

We also give a similar discussion for the Bernal stacking. Here, with only the nearest-neighbor hybridizations, a connected channel through the crystal is formed [Fig. 10(e)]. For the second-nearest-neighbor layers, the relative shift is zero and therefore the rule depicted in Fig. 9(d) applies. The remaining disconnected states become connected by the second-nearest-neighbor hybridizations and form channels through the crystal. Hence, there should emerge two types of bands, the former of which has larger bandwidth dominated by the nearest-neighbor hybridizations and the latter of which is flatter and dominated by the next-nearest hybridizations. The result of first-principles calculation is shown in Fig. 10(f).

Finally, we consider the AA case. According to the rule of Fig. 9(d), all the interlayer hybridizations are allowed among the C and C states, respectively. This yields two distinct channels connected by the nearest-neighbor hybridizations [Fig. 10(h)]. These channels represent emergence of two-fold degenerate bands dominated by the nearest-neighbor hybridizations, which is consistent with the first-principles calculation presented in Fig. 10(i).

Comparing the band structures, one would notice similarities in the energy scales. The more dispersive bands in the Bernal stacking has almost the same bandwidth as that of the band in the AA. The splitting of the rhombohedral bands is approximately half the bandwidth of the AA bands. On the other hand, the widths of the bands in the rhombohedral stacking and the flatter band in the Bernal stacking are almost the same. It is now clear that the former and latter features are dominated by the nearest-neighbor and next-nearest neighbor interlayer hybridizations, respectively.

### iv.4 Black phosphorus

The recent successful exfoliation of black phosphorus monolayer (black phosphorene) has revealed its remarkable property as a two-dimensional semiconductor Li-Zhang-phosphorene-transistor-nnano (); Liu-Tomanek-phosphorene-mobility-acsnano (), which has added a new member to the possible building blocks of two-dimensional semiconductor devices. Numerous studies on its electronic structures in the monolayer and multilayer forms have been reported so far (for a review, see Refs. Morita-review1986, ; Liu-Ye-review2015, ; Ling-Dresselhaus-review2015, ; Carvalho-CastroNeto-review2016, ). We reanalyze them from the viewpoint of the present theory. The black phosphorene forms in rectangular-net structure (Fig. 11(a)), which is a marked feature compared with the systems addressed in the previous subsections, and therefore this should give an example if our theory on the interlayer hybridization is valid in a net other than the hexagonal one.

As revealed in previous studies, the black phosphorene exhibits a direct gap at the point, whereas at the Brillouin-zone edge, where is located (Fig. 11(b)), all the bands are doubly degenerate (Fig. 11(c)), which are forced by the presence of the two-fold screw and glide mirror symmetries (Fig. 11(a)) Herring-nonsymmolph (); Heine-textbook (). As in the case of graphene, for such degenerate bands, we cannot derive general selection rules for the interlayer hybridization only from the stacking shift. Nevertheless, we find below that the two statements derived from our theory are valid: (i) The interlayer hybridization of the states within the manifold tends to be much suppressed (Sec. III.3), and (ii) the stacking-shift dependence of the electronic band structure is especially appreciable around (Sec. III.4).

We calculated the band structure of the bulk black phosphorus; namely, black phosphorene in layered bulk forms. The familiar type of bulk black phosphorus [AB stacking (Fig. 11(d))] has a structure where the inequivalent mirror planes are shared with the neighboring layers. In this case, the interlayer shift vector and (wavy line in Fig. 11(b)) can therefore be employed. The calculations were then carried out for the AB and pathological AA stacking cases [; Fig. 11(e)]. The band structures for the Z––Y–T and U–X–S–R paths (Fig. 11(f)) are shown in Fig. 11(g). Remarkably, the bands along the paths –Z and X–U are quite similar. On the other hand, in the paths Y–T and S–R, which are across , the similarity is apparently lower and the dispersions tend to be smaller in the AB case. This feature represents the marked contrast in the hybridization properties in the two types of stacking and provides another example demonstrating the general applicability of the present theory.

## V Application to general three-dimensional crystals

In the previous sections, we have considered what happens on the electronic state when multiple layers are stacked, in a unified view of the interference effect on the hybridization. We have demonstrated that our theory enables us analysis based only on the stacking geometry for points within . Below, we address an inverse problem: What can be derived by considering decomposition of three-dimensional crystals into stacked layers? Generally, three-dimensional periodic crystals are classified by the 14 Bravais lattices. Some of the Bravais lattices are decomposed into layers stacked so that the in-plane shift vector for the neighboring layers makes non-empty. Consequently, we obtain maps of -point paths for respective Bravais lattices, where the electronic band dispersions tend to be anomalously flat. This is due to the cancellation of the hybridization between the neighboring layers originating from the Bloch phase. We thereby propose a concept of the Bloch-phase induced flat-band paths (BIFP) for each Bravais lattice. Our theory reveals hidden -dependent anisotropy of the Bloch states in the three-dimensional lattice structure.

First, we consider a single-band tight-binding model on the FCC lattice [Fig. 12(a)]. We demonstrate decomposition of the lattice into the layers with shift and derive the -point paths where the band has little dispersion. This prediction is confirmed with first-principle calculations of the band structures in the NaCl-type semiconductors, whose structures belong to FCC.

### v.1 FCC tight-binding model

We here consider the single-band tight-binding model on the FCC lattice

 H = H1+H2+H3+⋯ (20) = t1∑⟨ij⟩1c†icj+t2∑⟨ij⟩2c†icj+t3∑⟨ij⟩3c†icj+… = ∑k[ϵ1(k)+ϵ2(k)+ϵ3(k)+…]~c†k~ck,

where denotes the th nearest-neighbor intersite hopping terms and denotes the summation over the th nearest pairs of sites. For example, the distances between the nearest-neighbor sites, next-nearest-neighbor sites, third nearest-neighbor sites are , , , , respectively, with being the lattice parameter [Fig. 12(a)]. We define () as creation (annihilation) operator for the local orbital belonging to a one-dimensional representation of the point group , which corresponds to the highest symmetry compatible with the FCC lattice. denotes the th nearest neighbor hopping amplitude. Regardless of the specific values of , the Hamiltonian is diagonalized with respect to the crystal wavenumber defined within the FCC BZ [Fig. 12(b)]. denotes the energy dispersion originating from .

Seen from the (1 1 1) direction, the FCC lattice structure can be regarded as hexagonal nets stacked with shift, where the common global axis corresponds to the inequivalent axes for the neighboring layers [Fig. 12(c)]. The “intralayer” primitive lattice vectors are introduced as depicted in solid arrows in panel (c), and the corresponding two-dimensional BZ for the individual nets is determined to be a green plane in panel (d). Then, let us construct partial Bloch sums for each layer by summing up the local orbitals within each layer. Here, labels the index of the layer and is the wavenumber defined for the two-dimensional BZ. The partial Bloch sums can be hybridized between the layers via inter-site hopping terms connecting the sites in different layers. The resulting eigenstates are labeled by the component of the wavenumber and show energy dispersion along this direction. However, between any pairs of layers shifted by [Fig. 12(c)], hybridizations are perfectly canceled for because of Eq. (17). This occurs when is at the “ or points” (namely, and points of the two-dimensional hexagonal BZ). Namely, we get

 ⟨ΨL,K2D|H|ΨL′,K2D⟩=0 (21) (L′−L=3n+1,3n+2;n=0,±1,±2,…). (22)

Furthermore, it is easily found in the FCC crystal structure that up to the fifth-nearest hopping terms cannot connect the sites in the third- or further-nearest neighbor layers. In conjunction with Eq. (22), this fact yield

 ⟨ΨL,K2D|Hl|ΨL′,K2D⟩=0(∀L≠L′;l=1,2,3,4,5). (23)

This means that the energy dispersions along through the and points are not affected by these short-range hopping terms and dominated by . Namely,

 ϵl(k)=const.(l=1,2,3,4,5) (24)

for . Thus, the band dispersion along these paths is anomalously flat even if the near-site hopping amplitudes are large. Discussions for the symmetrically equivalent directions such as (1 1 -1) also yield different paths where the band is anomalously flat.

Next, the FCC lattice can be decomposed into the stacked square nets if seen from the (1 0 0) direction. With this view, inequivalent axes, axes, and planes of the odd-numbered and even-numbered layers coincide [Fig. 12(f)]. The intralayer primitive lattice vectors and the corresponding BZ are defined as depicted in panels (f) and (g). Interlayer hybridizations are then prohibited for with and , which represent all the points at the BZ edge (red or dark solid lines in panels (g) and (h); see Appendix B). This prohibition yields

 ϵ1(k)=const. (25)

for .

Seen from the (1 1 0) direction, on the other hand, the lattice can be decomposed into the stacked rectangular nets (panel (i)). With this stacking, the global axis and mirror planes and correspond to inequivalent symmetry reference axis/planes for the neighboring layers, respectively. The intralayer primitive lattice vectors and the corresponding BZ are defined as depicted in panels (i) and (j). According to the results in Fig. 4, interlayer hybridizations are prohibited for all the Bloch states at the BZ edge (panels (j) and (k)). In this case, however, even can generate the hopping between the second-nearest-neighbor layers, which are not shifted from each other. Then, cannot be constant along the paths through the BZ edges toward the (1 1 0) direction.

Applying the present theory to the tight-binding model on the FCC lattice, we found -point paths where some near-site hopping terms does not contribute to the band dispersion. Significantly, in the paths along the (1 1 1) and (1 0 0) direction, the nearest-neighbor hopping terms do not give any dispersion. Here we note two key points in the above discussions: They respect only the kind of Bravais lattice and the symmetry compatible to it; second, although we have addressed a single-band model, the cancellation of the interlayer hybridization at is obviously valid in multiband systems, at least for intra-orbital ones (Sec. III.3). Anomalously flat dispersions are therefore expected to be frequently realized along the same paths regardless of the atomic configuration in the unit cell. We hence define the BIFP as these paths for the FCC lattice. The general emergence of flat bands along the BIFP is later exemplified by the first-principles calculations for solid materials with the NaCl-type structure.

### v.2 NaCl-type crystalline materials

We here investigate the band structures of NaCl-type semiconductors: LiH and NaCl. Their structure does not show clear layered character as depicted in Fig. 13(a). Nevertheless, below we observe many anomalously flat dispersions along the BIFP. This represent that the corresponding Bloch states have very anisotropic effective masses along the BIFP despite the three-dimensional crystal structure, which is due to the interference of the Bloch phase induced by the stacking with shift.

We calculated the band structures of LiH and NaCl from the first principles for the paths drawn in Figs. 13(b) and (c). The detailed conditions for the first-principles calculations are available in Appendix  D. Since the NaCl structure is -symmetric seen from (1 1 1) and -, -, - and -symmetric seen from (1 0 0), the , , and paths belong to the BIFP, where the band structure is expected to show little dispersions. Also, any paths along the (1 0 0) direction through the edge of the 2D square BZ belong to the BIFP. In order to examine this, we introduce the and points at the halfway in the and paths, respectively, and take the path.

For LiH [Fig. 13(d)], we see that the dispersion of the valence band formed by the H- orbitals is much less in the , , and paths than in the and paths, which verifies our expectations. On the other hand, for the conduction bands, although we see relatively flat dispersions for some of the bands, the dispersions are generally appreciable even in the former paths.

A possible origin of the dispersions is obviously the hybridization between distant layers because the wave functions of the conduction states are more delocalized than those of the valence states. Another origin is the inter-orbital interlayer hybridizations. The conduction bands are formed by Li- and orbitals, which are approximately degenerate in energy. Since the in-plane Bloch sums formed by them corresponds to various representations of the in-plane symmetry, many types of hybridizations are possible even between the neighboring layers, though we do not go further about this point. Note that for the (1 0 0) direction, bands are flatter in the path compared with in the and paths. This is probably due to the energy splitting between the Li conduction states: In the isolated monolayer form, they well split at the point because of the low symmetry, which reduce the impact of the inter-band hybridizations in the multilayer form (See Appendix C).

For NaCl [Fig. 13(e)], we also see a number of bands with small dispersion in the BIFP. However, the dispersions of the valence bands are more significant compared with the LiH case. This is due to the degeneracy of the Cl- states forming the valence bands. Note that we again find that the dispersions tend to be flat especially in the path for the both valence and conduction bands. Similarly to the LiH case, this is probably due to the lifted degeneracy of the Cl- states induced by the low symmetry.

We have thus demonstrated that band dispersions indeed tend to be small in the BIFP. Although the magnitude of the flattening is essentially affected by the band degeneracy, we assert that similar flat dispersions should appear in the same paths for general FCC crystals.

### v.3 General three-dimensional lattices

It is worth mentioning again that the above discussion for the FCC system is always possible regardless of the structure within the unit cell. The key property of the FCC lattice is that it can be decomposed into layers so that their in-plane shift vector makes nonzero in the 2D BZ for these layers. As a matter of fact, 11 of the 14 Bravais lattices have this property (Table 1). Since the Brillouin-Zone shape is determined by the Bravais lattice, for these lattices, we can also define the BIFP and formulate one-to-one correspondences for the BZ and the BIFP. We provide in Appendix E the representative BIFP derived from the decomposition. Note that we cannot give the complete list of BIFP for some Bravais lattices (see Appendix F).

In Sec. V.1, we have seen that the nearest-neighbor hopping does affect the band structure along the BIFP with respect to the (1 1 0) direction. This represents that the probability of observing the flat bands along each BIFP depends on the lattice geometry. We can characterize this probability by referring to the distance between the neighboring layers in the decomposition scheme. For the BIFP yielded by the in-plane symmetry , we define the threshold distance as

 dth≡dth(S)={2clayer  (S≠C3)3clayer  (S=C3). (26)

Namely, corresponds to the shortest distance between the “unshifted” layers—between which interlayer hybridization is not canceled. The following statement is consequently valid. Take an intra-orbital inter-site hopping between a certain pair of sites: If the distance between the sites is shorter than , such hopping does not contribute to the band dispersion along the BIFP. Equantion (26) clarifies the advantage of the BIFPs yielded by , where we should find flat dispersions especially frequently. The paths for bulk -BN, 3R-MoS and -graphite in Sec. IV and those along the (1 1 1)-BIFP in the NaCl-type crystals in Sec. V.2 correspond to this case.

We occasionally encounter the cases that, in addition to the BIFP related to the Bravais lattice, specific atomic configurations within the unit cell and lattice parameters result in extra -point paths where near-site hoppings do not affect the band dispersions. We have in fact seen a typical example in Sec. IV.1.2 as -BN. Although the primitive hexagonal lattice cannot generally yield the BIFP along the (0 0 1) direction (see Appendix E), the hopping terms between the nearest-neighbor layers do not affect the dispersion along the path. Construction of the BIFP for such systems can be done by introducing the site-specific nets as demonstrated in Sec. IV. The one-to-one correspondences between each space group (230 in total) and the BIFP could be formulated, though we do not address this issue in this paper.

## Vi Discussions

### vi.1 Relevance to experiments

In the previous sections, we have demonstrated that there are -points where the interlayer hybridizations are prohibited depending on the stacking shift and that much affects the electronic structures in multilayered systems. This characteristic is especially utilizable for van der Waals-coupled multilayer materials because they can adopt various stacking patterns thanks to the weak interlayer binding. In particular, when the valence-band top and/or conduction-band bottom are located at any of in the monolayer form, one can demonstrate stacking-dependent characteristics of electronic states with low-energy probes. It is expected to be often the case because includes some of the high-symmetry special points, where the band extrema are located. In fact, some of the authors have shown that the dimensionality of the excitons in MoS is controlled by the stacking geometry through the reflectivity measurement Akashi2015 (). Since there has recently been a notable progress in the experimental technique to control the stacking shift Kim-graphene-pret-a-porte-PRB (); Jiang-Wu-MoS2-fold-NNano (), our theory should have a broad potential applicability.

Even in the general three-dimensional crystals, the Bloch-phase interference can serve significant effects if the stacked-bilayer structure is formed locally. For example, BiS-based superconductors Mizuguchi-Bi4O4S3-PRBR2012 (); Mizuguchi-LaOFBiS2-JPSJ2012 () having BiS bilayers show abrupt increase of the superconducting transition temperature  Wolowiec-Yazici-pressure-PRB2013 () concomitantly with the pressure-induced structural phase transition from the tetragonal to monoclinic structure Tomita-Mizuguchi-pressure-structure-JPSJ2014 (); Guo-Yuan-EuBiS2F-pressure-structure-PRB2015 (). It has been pointed out that the strong interlayer hybridization at the conduction-band bottom is switched on by this transition Ochi-BiS2-JPSJ (), which could be crucial for enhancing the . This switching is due to the change of the interlayer shift, which breaks the interference condition. For the low-dimensional Bloch states emerging along the BIFP, however, it will be difficult to observe. Such states usually emerge in the middle of the valence/conduction bands because the BIFP usually do not pass through the special points of the three-dimensional BZ. In order to observe such high energy states, heavy doping or a wavenumber- and energy-resolved electronic probe is necessary. Still, we believe that our maps of the BIFP are useful for electronic materials search and design because of their general applicability. With a broad range of candidate materials, it could be possible to find suitable ones where we can detect the low-dimensional Bloch states. In fact, in LiH, the energy of the flat band in the path is near the Fermi energy, where holes could be injected with chemical or field-effect doping.

### vi.2 Remark for tight-binding modeling

A general remark on the development of tight-binding models for layered materials is derived from our theory. When one theoretically studies a layered system, only a small number of inter-site hopping matrix elements are usually retained in the model, particularly for those between the sites of different layers. However, when the layers are shifted from each other, such modeling will suffer from underestimation of interlayer transfer of electrons in the states within the SAIM, unless the hopping terms between “unshifted” layers are included. Take for example the case of the graphite with the Bernal or rhombohedral stackings (Sec. IV.3). If only the hopping matrix elements up to the nearest-neighbor layers are retained in the modeling, the electrons in the -point states never propagate across the whole layers (See Fig. 10 (e)). To avoid such too idealized situations, the inter-site hopping between the second-nearest neighbor layers ( in the Partoens-Peeters modeling Partoens-multilayer2006 ()) is needed for the both stackings.

### vi.3 Possible extensions

In this work, we have addressed only the interlayer transitions of electrons induced by the intrinsic Hamiltonian. The important fact is that the interlayer transitions concerning electronic states in the SAIM is subject to the stacking geometry. This is also extended to more general interlayer transitions via external perturbations such as incident light. The only difference is that the operator does not always transform according to the identity representation of the space group as does the Hamiltonian. Reconsideration of optical absorption spectra for few-layer graphene Mak-Heinz-fewlayer-graphene-PRL2010 () from the present viewpoint could give a new insight. Also, an extension to phononic systems is possible considering that the displacement vector of the phonon mode obeys the Bloch theorem. By expanding the total Born-Oppenheimer energy with respect to the mode displacement vectors for the isolated layers, we can derive the rule for which the second derivatives become zero. Our maps for the SAIM and BIFP will be useful when one considers such extensions.

## Vii Summary and Conclusion

We have addressed the Bloch-phase interference on the electronic motion in layered crystals and presented a general theory on the interlayer hybridizations in view of the stacking geometry. The key result is the development of the SAIM in the two-dimensional BZ, displaying which Bloch states are prevented from the interlayer hybridization depending on the stacking shift. We have demonstrated that our theory is useful for studying the electronic states in multilayer systems with applications to BN, graphite, TMD and black phosphorus. Also, we have shown that the prevention of the interlayer hybridization is robust even in general three-dimensional crystals and it induces strongly anisotropic electronic bands where effects of some inter-site hopping amplitudes are canceled. The theory is based on the classification of the periodic lattice structure and therefore in principle applicable to every material with periodic structure, which gives us a simple view on understanding and controlling the electronic states.

###### Acknowledgements.
This work was supported by MEXT Element Strategy Initiative to Form Core Research Center in Japan and JSPS KAKENHI Grant Number 15K20940 (R. A.) from Japan Society for the Promotion of Science (JSPS). We thank Ryotaro Arita, Shinji Tsuneyuki, Masayuki Ochi, Peter Maksym, Susumu Saito, and Ryuji Suzuki for careful reading of the manuscript and valuable comments.

## Appendix A Derivation of Eqs. (5) and (6)

We first derive Eq. (5). Suppose the origin of is located on the reference point of . We then trivially get

 S1eik⋅r=eik⋅S−1r=ei(Sk)⋅r (27)

and

 φpw1(r)=ei(Sk−k)⋅r. (28)

Equation (5) is also obtained straightforwardly:

 Sjeik⋅r ≡ Tτj−1S1T−1τj−1eik⋅r (29) = e−i(Sk−k)⋅τj−1ei(Sk)⋅r,

where we used for general function . This form is translated to Eq. (5).

Before going to Eq. (6), we introduce the plane-wave expansion of the Bloch function

 ψk(r)=eik⋅r[∑GeiG⋅rck,G], (30)

where the factor in the bracket corresponds to the cell-periodic part and is the expansion coefficient. is the reciprocal lattice vector. To obtain Eq. (6), we here derive constraints between , and imposed by the symmetry. Since is an eigenfunction of , the following expression holds: . Transforming this equation into the plane-wave expansion form,

 ∑Gei(SG+Sk)⋅rck,G=∑Gei(G+k)⋅reiφ1ck,G. (31)

Dividing the both sides by and redefining the vector , we get

 ∑GeiG⋅rck,S−1(G−G0)=∑GeiG⋅reiφ1ck,G. (33)

Here, a reciprocal lattice vector is given by , which is required for to be an eigenfunction of . The uniqueness of the Fourier expansion yields

 eiφ1=ck,S−1(G−G0)ck,G  (∀G). (34)

Applying the similar discussion to ,

 eiφjck,Gck,S−1(G−G0)ei(S−1G0)⋅τj−1 =e−i(G−S−1G)⋅τj−1  (∀G). (35)

Since the left-hand side is independent of because of Eq. (34), the right-hand side must also be independent of ; namely,

 exp[−i(G−S−1G)⋅τj−1]=const.=