Quasicrystalline electronic states in twisted bilayers and the effects of interlayer and sublattice symmetries
EElectronic states in two-dimensional bilayer quasicrystals and the effects of interlayerand sublattice symmetries
J. A. Crosse
1, 2 and Pilkyung Moon
1, 2, 3, 4, ∗ Arts and Sciences, NYU Shanghai, Shanghai, China NYU-ECNU Institute of Physics at NYU Shanghai, Shanghai, China Department of Physics, New York University, New York, USA State Key Laboratory of Precision Spectroscopy,East China Normal University, Shanghai, China (Dated: June 12, 2020)We study the electronic structure of quasicrystals composed of incommensurate stacks of atomiclayers. We consider two systems: a pair of square lattices with a relative twist angle of θ = 45 ◦ anda pair of hexagonal lattices with a relative twist angle of θ = 30 ◦ , with various interlayer interactionstrengths. This constitutes every two-dimensional bilayer quasi-crystal system. We investigate theresonant coupling governing the quasicrystalline order in each quasicrystal symmetry, and calculatethe quasi-band dispersion. We find that some quasicrystalline states, which are usually obscuredby additional weakly coupled states, are more prominent, i.e., ”exposed”, in the systems withstrong interlayer interaction. We also show that we can switch the states between quasicrystallineconfiguration and its layer components, by turning on and off the interlayer symmetry. On the otherhand, hexagonal lattices with sublattice potential asymmetry, e.g., transition metal dichalcogenideof hexagonal Boron Nitride, give the quasicrystalline states quite similar to those in the system inthe absence of sublattice asymmetry. I. INTRODUCTION
When two honeycomb lattices are overlapped on topof the other at a twist angle θ = 30 ◦ , the atomic ar-rangement is mapped on to a quasicrystalline lattice,which is ordered but not periodic, with a 12-fold rota-tional symmetry . Recently, it has been demonstratedthat bilayer graphene with a precise rotation angle of 30 ◦ exhibits the atomic structures satisfying the quasicrys-talline tiling as well as the spectrum respecting the 12-fold rotational symmetry . Similar structures have alsobeen realized by growing bilayer graphene on top of theNi or Cu surface , and also by a transfer method .The conventional moir´e effective theory, which is basedon the period of the moir´e pattern arising from the inter-ference between the lattice periods, cannot describe theelectronic structures of such quasicrystals composed ofincommensurate stack of atomic layers (hereafter ”vdW-QCs”) since the rotational symmetry of quasicrystalsdoes not commute with translation. In our previouswork, we developed a momentum-space tight-bindingmodel which can describe the electronic structures ofatomic layers stacked at any configuration without rely-ing on the moir´e periodicity . This model enabled us toreveal the quasi-band dispersion and the emergence of theelectronic states having the quasicrystalline order in thevdW-QC composed of two graphene layers stacked at 30 ◦ by fully respecting the rotational symmetry of quasicrys-tals as well as the translational symmetry of constituentlayers. While conventional quasicrystals can be viewedas intrinsic quasicrystals where all the atomic sites areintrinsically arranged in the quasiperiodic order, vdW-QCs are regarded as extrinsic quasicrystals , in that theyare composed of a pair of perfect crystals having indepen-dent periodicities, and the quasiperiodic nature appears only in the perturbational coupling between the two sub-systems. Thus, vdW-QCs provide a unique opportunityto design quasicrystalline states by using atomic layerswith various symmetries and interaction strength.In this paper, we numerically investigate the electronicstructures of vdW-QCs for every possible rotational sym-metry in two-dimensional space. Since a periodic two-dimensional atomic layer can have 2-, 4-, 6-fold rotationalsymmetry, we can make only 8-fold [octagonal, Fig. 1(a)]or 12-fold [dodecagonal, Fig. 1(d)] vdW-QCs by stackingtwo square lattices at 45 ◦ or by stacking two hexagonallattices at 30 ◦ , respectively. We first find the resonantcondition, which gives quasicrystalline order to the elec-tronic states, in each system, and calculate the quasi-band dispersion and density of states (DOS) for variousinterlayer interaction strength. We identify the featureswhich arise from the quasicrystalline order as opposed tothose arising from the interaction common to any other θ in the spectrum of vdW-QCs. In addition, we show thatsome quasicrystalline states, which are usually obscuredby additional weakly coupled states, are more prominentin vdW-QCs with strong interlayer interaction.We also investigate the effects of lifting both inter-layer and sublattice symmetry on the electronic struc-ture. Since the quasicrystalline order arises from the res-onant interaction between the states in both layers, inter-layer potential asymmetry results in a dramatic changein the electronic structure. In addition, we can switch be-tween states respecting quasicrystal symmetry and thosesatisfying only half the symmetry by turning on and offthe interlayer symmetry. On the other hand, we showthat sublattice potential asymmetry, such as that in tran-sition metal dichalcogenide, does not make a dramaticdifference - it results only in a constant energy shift.We also analytically interpret the mixing between the a r X i v : . [ c ond - m a t . m e s - h a ll ] J un quasicrystalline states, which may influence the physicalproperties (such as optical selection rules), and changesto the location of the band edge via interlayer or sublat-tice potential asymmetry.The paper is organized as follows. In Sec. II, we presentthe atomic structures and tight-binding model for vdW-QCs, and introduce the dual tight-binding approach inthe momentum space. In Sec. III A, we find the res-onant interaction which gives the quasicrystalline elec-tronic states in vdW-QCs. In Sec. III B and III C, we de-rive the minimal Hamiltonian and calculate the band dis-persion and wave functions of octagonal and dodecagonalvdW-QCs, respectively. We also investigate the effects ofvarious interlayer interaction strengths, the features aris-ing from 2-wave mixing, and the effects of the lifting ofinterlayer and sublattice potential asymmetry. A briefconclusion is given in Sec. IV. II. THEORETICAL METHODSA. Atomic structure and Brillouin zones ofquasicrystalline twisted bilayers
We define the atomic structure of the octagonal vdW-QCs by starting from two perfectly overlapping squarelattices and rotating the layer 2 around the center of thesquare by θ = 45 ◦ [Fig. 1(a)]. We set xy coordinates par-allel to the layers and z axis perpendicular to the plane.The system belongs to the symmetry group D d , and it isinvariant under an improper rotation R ( π/ M z , where R ( θ ) is the rotation by an angle θ around z axis, and M z is the mirror reflection with respect to xy plane.The primitive lattice vectors of layer 1 are taken as a = a (1 ,
0) and a = a (0 , a is the lattice con-stant, and those of the layer 2 as ˜ a i = R ( π/ a i . In thispaper, we model the square lattices by a minimal, oneorbital model with one sublattice site. Then, the atomicpositions are given by R X = n a + n a + τ X (layer 1) , R ˜ X = ˜ n ˜ a + ˜ n ˜ a + τ ˜ X (layer 2) , (1)where n i and ˜ n i are integers, X ( ˜ X ) denotes the sublat-tice site (only one in this case) of layer 1 (layer 2) of whichposition in the unit cell is defined by τ X = ( a/ , a/ τ ˜ X = R ( π/ τ X + d e z ]. Here, d is the interlayer spacingbetween the two layers and e z is the unit vector normalto the layer. The reciprocal lattice vectors of layer 1 aregiven by a ∗ = (2 π/a )(1 ,
0) and a ∗ = (2 π/a )(0 , a ∗ i = R ( π/ a ∗ i .Likewise, we define the atomic structure of the do-decagonal vdW-QCs by starting from two perfectly over-lapping honeycomb lattices (i.e., AA-stacked bilayers)and rotating the layer 2 around the center of the hexagonby θ = 30 ◦ [Fig. 1(d)]. The system belongs to the sym-metry group D d , and it is invariant under an improperrotation R ( π/ M z . The primitive lattice vectors of layer 1 are taken as a = a (1 ,
0) and a = a (1 / , √ / a is the lattice constant, and those of the layer2 as ˜ a i = R ( π/ a i . The atomic positions are given byEq. (1), where X = A, B ( ˜ X = ˜ A, ˜ B ) denotes the sub-lattice site of layer 1 (2), and τ X and τ ˜ X are the sub-lattice positions in the unit cell, defined by τ A = − τ , τ B = τ , τ ˜ A = − R ( π/ τ + d e z , τ ˜ B = R ( π/ τ + d e z with τ = (0 , a/ √ d is the interlayer spac-ing between the two layers. The reciprocal lattice vec-tors of layer 1 are given by a ∗ = (2 π/a )(1 , − / √
3) and a ∗ = (2 π/a )(0 , / √ a ∗ i = R ( π/ a ∗ i .The red and blue squares (hexagons) in Fig. 1(c)[Fig. 1(f)] show the Brillouin zones of layer 1 and 2 inoctagonal (dodecagonal) vdW-QCs, respectively. B. Tight-binding model for van der Waals bilayers
We model both systems by the tight-binding modelwith a single orbital of arbitrary atomic species. Al-though we use p z orbitals in this paper, just like themodel of graphene and hexagonal Boron Nitride, but itcan be any other orbital optimal for each system; e.g., intransition metal dichalcogenides, we can use d -orbitals ,and so on. The Hamiltonian is spanned by the Blochbases of each sublattice, | k , X (cid:105) = 1 √ N (cid:88) R X e i k · R X | R X (cid:105) (layer 1) , | ˜ k , ˜ X (cid:105) = 1 √ N (cid:88) R ˜ X e i ˜ k · R ˜ X | R ˜ X (cid:105) (layer 2) , (2)where | R X (cid:105) is the atomic p z orbital at the site R X , k and˜ k are the two-dimensional Bloch wave vectors and N = S/S tot is the number of the unit cells with an area S [ S = a for square lattices and S = ( √ / a for honeycomblattices] in the total system area S tot .We use a Slater-Koster parametrization for thetransfer integral between any two p z orbitals, − T ( R ) = V ppπ (cid:34) − (cid:18) R · e z | R | (cid:19) (cid:35) + V ppσ (cid:18) R · e z | R | (cid:19) , (3)where R is the relative vector between two atoms, and V ppπ = V ppπ e − ( | R |− a ) /δ ,V ppσ = V ppσ e − ( | R |− a ) /δ , (4)for square lattices and V ppπ = V ppπ e − ( | R |− a/ √ /δ ,V ppσ = V ppσ e − ( | R |− . a ) /δ , (5)for honeycomb lattices so that the first-nearest neighborintralayer coupling becomes V ppπ ( < δ = δ = 0 . a . x [ a ] -8 0 8-4 4 y [ a ] -8-4048 (a) x [ a ] -8 0 8-4 4 y [ a ] -8-4048 (d) k x [1/ a ] -4 -2 0 2 4 k y [ / a ] -4-2024 (c) k x [1/ a ] -15 -5 0 10 15-10 5 k y [ / a ] -15-1001015-55 (b) k x [1/ a ] -20 -10 0 10 20 k y [ / a ] -20-1001020 (e) k x [1/ a ] -4 -2 0 2 4 k y [ / a ] -4-2024 (f) ~**a ~a **a a ~**a ~ a **a MC
02 34 567 89 1011 1
KK~
FIG. 1. (a) Real-space lattice structures of octagonal vdW-QCs. Red and blue squares (circles) represent the unit cells(atomic sites) of layer 1 and 2, respectively. (b) Dual tight-binding lattice in the momentum space for octagonal vdW-QCs (seetext). Red and blue squares show the extended Brillouin zones of layer 1 and 2, respectively. The number n represents theposition of Q n ( n = 0 , , , · · · , k for layer 1, and blue ones represent the inverted wave numbers ˆ k − ˜ k for layer 2, where ˆ k is taken as Q here. (c) The wave numbers C n ( n = 0 , , , · · · ,
7) at the cross points between the first Brillouin zones of thetwo lattices, which are the original positions of k (layer 1) and ˜ k (layer 2) associated with Q n . The dashed line indicates theconnection in the 8-ring Hamiltonian as in (b). Due to the symmetry, these 8 wave numbers are all degenerate in energy. (d),(e), (f) Plots similar to (a), (b), (c) for dodecagonal vdW-QCs. Here, n for Q n and C n runs from 0 to 11, and ˆ k = . The total tight-binding Hamiltonian is expressed as H = H + H + U + H V + H ∆ (6)where H and H are the Hamiltonian for the intrinsicsquare or honeycomb lattices of layer 1 and 2, respec-tively, U is for the interlayer coupling, H V is for the inter-layer potential asymmetry, and H ∆ is for the sublatticepotential asymmetry. Note that the square lattices withthe minimal, one orbital model considered in this workdoes not have the H ∆ term as there is only one sublatticein this case. The intralayer matrix elements of layer 1 aregiven by (cid:104) k (cid:48) , X (cid:48) |H + H V + H ∆ | k , X (cid:105) = [ h X,X (cid:48) ( k ) + { V / s ∆ ∆ / } δ X,X (cid:48) ] δ k (cid:48) , k ,h X,X (cid:48) ( k ) = (cid:88) L − T ( L + τ X (cid:48) X ) e − i k · ( L + τ X (cid:48) X ) , (7) where L = n a + n a , τ X (cid:48) X = τ X (cid:48) − τ X , V is themagnitude of the interlayer potential asymmetry, s ∆ is+1 and − X = A and B , respectively, and ∆ is themagnitude of the sublattice potential asymmetry. Sim-ilarly, the matrix for H is given by replacing k with R ( − π/ k in a square lattice and R ( − π/ k in a hon-eycomb lattice, V / − V /
2, and s ∆ by +1 and − A and ˜ B , respectively. The interlayermatrix element between layer 1 and 2 is written as (cid:104) ˜ k , ˜ X | U | k , X (cid:105) = − (cid:88) G , ˜ G t ( k + G ) e − i G · τ X + i ˜ G · τ ˜ X δ k + G , ˜ k + ˜ G , (8)where G = m a ∗ + m a ∗ and ˜ G = ˜ m ˜ a ∗ + ˜ m ˜ a ∗ ( m , m , ˜ m , ˜ m ∈ Z ) run over all the reciprocal pointsof layer 1 and 2, respectively. Here t ( q ) = 1 S (cid:90) T ( r + z ˜ XX e z ) e − i q · r d r (9)is the in-plane Fourier transform of the transfer integral,where z ˜ XX = ( τ ˜ X − τ X ) · e z . For H and H , we consideronly the nearest neighbor interactions in octagonal vdW-QCs to keep the symmetric and simple cosine bands, andall interactions of | R | within 3 a in dodecagonal vdW-QCsto make it consistent with previous works . C. Effects of using different atomic layers andinterlayer interaction
Different vdW-QCs exhibit different electronic struc-tures, since different atomic layers have different inter-layer interaction. In addition, we can also tune the inter-layer interaction in a given vdW-QC with the interlayerdistance d , e.g., by applying an external pressure or in-tercalation. However, if two different vdW-QCs have thesame sublattice configuration and orbitals, the interac-tions [Eq. (8)] in the two systems have the same phaseand the same direction q / | q | dependence, and differ onlyin the magnitude of the interaction t ( q ). In most vander Waals systems with d > a , t ( q ) scales with V ppσ and exponentially decays with d . Thus, we can inves-tigate the electronic structures of various vdW-QCs bysimply scaling t ( q ). In Sec. III B and III C, we will firstinvestigate the electronic structures of a vdW-QC with t ( q ) obtained at a specific combination of ( V ppσ , d ), theninvestigate the effects of using different materials and in-teraction strength by tuning t ( q ). D. Dual tight-binding model in momentum space
In real-space, since quasicrystals do not have period-icity, we need infinitely many bases to solve Eq. (6). Al-though some conventional approximations with a finitenumber of bases, such as a periodic approximant or afinite-size model, can give an energy spectrum quite sim-ilar to the actual spectrum, the resulting wave functionslose their long-range quasicrystalline nature.Instead, we can solve Eq. (6) rigorously by using atight-binding model in k space, which is the dual coun-terpart of the original tight-binding Hamiltonian in thereal space. Equation (8) shows that the interlayer inter-action occurs between the states satisfying the general-ized Umklapp scattering condition k + G = ˜ k + ˜ G . It isstraightforward to show that the entire subspace spannedby H + H + U from a layer 1’s Bloch state at ˆ k is given by {| k , X (cid:105) | k = ˆ k + ˜ G , ∀ ˜ G } and {| ˜ k , ˜ X (cid:105) | ˜ k = ˆ k + G , ∀ G } .According to Eq. (8), the interaction strength between k = ˆ k + ˜ G and ˜ k = ˆ k + G is given by t ( q ) where q = k + G = ˜ k + ˜ G = k + ˜ k − ˆ k . Then, the interactionstrength can be visualized by the diagram Figs. 1(b) and (e), where all the layer 2’s wave points ˜ k are inverted toˆ k − ˜ k , and overlapped with the layer 1’s wave points k . Inthe map, the quantity | q | = | k + ˜ k − ˆ k | is the geometricaldistance between two points, and the interaction takesplace only between the points located in close distance,since t ( q ) decays in large q . Therefore, in practical cal-culation, we only need a limited number of states aroundˆ k inside a certain cut-off circle k c . Hence, we can obtainthe energy eigenvalues at ˆ k by diagonalizing the Hamil-tonian matrix within a finite set of bases. If the k pointsare viewed as “sites”, the whole system can be recog-nized as a tight-binding lattice in k space, which is dualto the original Hamiltonian in the real space. This en-ables us to calculate the electronic structures of almostevery possible stack of atomic layers without relying onmoir´e periodicity. The method as well as its validity areexplained in detail in Ref. .Below, we will first study the eigenstates of H + H + U by using the dual-tight binding method, then investigatethe effects of H V and H ∆ in terms of those eigenstates. III. RESULTS AND DISCUSSIONA. Resonant states respecting the rotationalsymmetry of quasicrystals
The rotational symmetry of the quasicrystal as well asthe translational symmetries of the constituent atomiclayers [accordingly, the generalized Umklapp scatteringcondition Eq. (8)] reveal the most dominant interaction,which comes from the resonance between degeneratestates, in each vdW-QC. In octagonal vdW-QCs, we seethat the eight symmetric points Q n = [ − π/a, − π/a ] +( a ∗ / √ nπ/ , sin( nπ/ n = 0 , , , · · · ,
7) forma circular chain in the dual tight-binding lattice withˆ k = Q . The chain has a radius of a ∗ / √ ≡ | a ∗ i | / √ √ π/a , and is indicated by the dashed ring in Fig. 1(b).Noting that the layer 2’s wave points are inverted, thesepoints are associated with layer 1’s Bloch wave numbers k = Q n for even n ’s and layer 2’s ˜ k = Q − Q n for odd n ’s. Figure 1(c) shows the original positions of k (layer 1)and ˜ k (layer 2) associated with Q n in the first Brillouinzone, C n = √ a ∗ sin 22 . ◦ [cos(3 π/ nπ/ , sin(3 π/ nπ/ n = 0 , , , · · · , t = t ( | C n | ). The interaction to otherneighboring states in the dual tight-binding lattice canbe safely neglected since the interaction strength is muchless than t and the two states are not degenerate in mostcases.Likewise, in dodecagonal vdW-QCs, wesee that the twelve symmetric points Q n = a ∗ [cos( nπ/ , sin( nπ/ n = 0 , , , · · · ,
11) form acircular chain in the dual tight-binding lattice withˆ k = and the radius is a ∗ ≡ | a ∗ i | = 4 π/ ( √ a ) [dashedring in Fig. 1(e)]. These points are associated withlayer 1’s Bloch wave numbers k = Q n for even n ’sand layer 2’s ˜ k = − Q n for odd n ’s, and Fig. 1(f)shows the original positions of k (layer 1) and ˜ k (layer 2) associated with Q n in the first Brillouin zone, C n = 2 a ∗ sin 15 ◦ [cos(5 π/ nπ/ , sin(5 π/ nπ/ n = 0 , , , · · · , t = t ( | C n | ).It should be noted that these states are not the only setof states which show the resonant coupling in each sys-tem. As we shown in Appendix A, there are more setsof states, with different wave numbers, showing the res-onant interaction respecting the rotational symmetry ofthe quasicrystals. However, the sets in Figs. 1(b) and (e)give the strongest interaction, i.e., largest energy separa-tion, since these states form the rings with the shortestdistance between neighboring states in the dual tight-binding lattices. B. Octagonal quasicrystal
1. Hamiltonian
In octagonal vdW-QCs, the strongest resonant inter-action occurs at ˆ k = Q . Thus, by replacing ˆ k with Q + k , we can express the Hamiltonian H = H ring + H V , (10)in the vicinity of k = , in the bases of( | k (0) (cid:105) , | k (1) (cid:105) , · · · , | k (7) (cid:105) ), where | k ( n ) (cid:105) is | k + C n , X (cid:105) for even n (layer 1) and | k + C n , ˜ X (cid:105) for odd n (layer2). Here, H ring ( k ) = H (0) − t − t − t H (1) − t − t H (2) − t . . . . . . . . . − t H (6) − t − t − t H (7) , (11)is the Hamiltonian matrix of the resonant ring inthe absence of interlayer potential asymmetry, where H ( n ) ( k ) = h X,X [ R ( − nπ/ k + C ], and we neglect the k dependence of the interlayer matrix element t ( q ). Thediagonal elements H ( n ) represent monolayer’s Hamilto-nian at k = k + C n for even n and ˜ k = k + C n for odd n . Note that h X,X in H ( n ) is same for any n , and the de-pendence of the diagonal elements on n solely comes from R ( − nπ/ k in the argument of h X,X . Consequently,the ring Hamiltonian H ring is obviously symmetric un-der rotation by a single span of the ring (i.e., moving C n to C n +1 ), which actually corresponds to the operation[ R ( π/ M z ] (225 ◦ rotation and swapping layer 1 and 2)in the original system. In addition, H ring has a particle-hole symmetry with respect to the energy E = h , where h = h X,X ( C ) ≈ − V ppπ (cos √ π + 1), up to the firstorder to k (Appendix B). H V is the Hamiltonian representing the interlayer po-tential asymmetry, H V = V σ z σ z σ z σ z , (12)where V ( ≥
0) represents the difference in the electro-static energies between the two layers, and σ i is the Paulimatrix. With H V , the Hamiltonian H , which was origi-nally in the form of one-dimensional monatomic chain inthe dual-tight binding lattice, becomes that of diatomicchain with alternating on-site potential.
2. Band structures and wave functions
Figure 2(a) shows the band structures of the octagonalvdW-QCs near C n , in the absence of interlayer potentialasymmetry (i.e., V = 0), plotted as a function of k , andFig. 2(b) shows its closer view near k = . We choose( V ppσ , d ) = ( − . V ppπ , . a ), which gives the interac-tion strength between the neighboring sites of the circu-lar chain in the dual-tight binding lattice t = 0 . V ppπ .The eight parabolic bands are arranged on a circle witha radius ∆ k = (2 − √ π/a by the Umklapp scattering[Eq. (8)], and they are strongly hybridized near k = .As a result, the originally degenerate eight states of thesquare lattices split into different energies, and exhibitthe characteristic dispersion including parabolic band-bottoms at this high energy regime and the frilled bandedge, which is flat up to the first order to k .At k = , H ring can be analytically diagonalized toobtain a set of energies E m = h − t cos q m , (13)which have the energy span of 4 t , where q m = (5 π/ m with m = 0 , ± , ± , ± , v m =( µ − m , µ − m , µ − m , · · · , µ m ) / √ µ m = e iq m ) is the coeffi-cient to the Bloch bases | k (0) (cid:105) , | k (1) (cid:105) , · · · , | k (7) (cid:105) . Here thestates with m = ± s ( s = 1 , ,
3) form twofold doublets,and belong to two-dimensional E s irreducible represen-tation of D d point group, while the m = 0 and 4 arenon-degenerate, and belong to A and B , respectively.If we disregard the z -position difference, the index m canbe regarded as quantized angular momentum. The factthat there are 8 unique values for m as well as the fact m = 0, 4 y [ a ] m = ±2 y [ a ] m = ±1, ±3 y [ a ] x [ a ] 8-8 40-4 x [ a ]8-8 40-4 x [ a ] ene r g y [ | V pp π | ] k [1/ a ] ene r g y [ | V pp π | ] k [1/ a ] (a) (b)(c) ±2±14 m = 0±3 FIG. 2. (a) Electronic structures of octagonal vdW-QCscalculated by the 8-ring effective model. Blue and red arrowsshow the band opening by the interlayer and intralayer 2-wave mixing, respectively (see Sec. III B 3). (b) Detailed bandstructures near k = C n [the region encircled by the blackdashed rectangle in (a)] with index m indicating quantizedangular momentum in 8-fold rotational symmetry. (c) LDOSat k = C n characterized by m , where the area of the circleis proportional to the squared wave amplitude, and red andblue circles represent the states in the upper and the lowerlayers, respectively. that the eigenvalue of R ( π/ M z is given by e iπm/ arethe evidence that the quasicrystalline electronic statesrespect an 8-fold rotational symmetry.The 8-wave resonant coupling also gives rise to a char-acteristic pattern in the wave function. Figure 2(c) showsthe wave functions at k = where the hybridization isthe most prominent. We can see that the wave amplitudeis distributed on a limited number of sites in a 8-fold ro-tationally symmetric pattern.
3. Effects of interlayer interaction strength and other kindsof interaction
As discussed in Sec. II C, we can calculate the qua-sicrystalline states of various octagonal vdW-QCs, whichare either composed of other materials or different in-terlayer distance d , by simply scaling the magnitude ofthe interlayer interaction t ( q ). Figure 3(a) shows theband structures of an octagonal vdW-QC with t ( q ) 2times larger than that in Fig. 2. Although the inter-action strength varies with q , hereafter we label each ene r g y [ | V pp π | ] ene r g y [ | V pp π | ] k [1/ a ] (a) (b) ±2±14 m = 0±3 k y [ / a ] k y [ / a ] k y [ a ] k x [ a ] 4-4 20-2 k x [ a ]4-4 20-2 k x [ a ] (c) |t | [| V ppπ |] α β min max FIG. 3. (a) Electronic structures of octagonal vdW-QCs withthe interlayer interaction | t | t . The white dashed line corresponds to the DOS forthe system considered in Fig. 2. The white arrows show thepseudogaps by the interlayer and intralayer 2-wave interac-tion. (c) Three representative interaction residing in octago-nal vdW-QCs; quasicrystalline interaction by 8-wave mixing(left), interlayer 2-wave mixing (middle), and intralayer 2-wave mixing (right). Red and blue contours show the Fermisurfaces of layer 1 and 2, respectively, black dashed circlesshow the wave numbers where the interaction occurs, and theblue dashed line shows the reciprocal lattice vector of layer 2which mediates the interaction between the states in layer 1. system with t ( q ) at q showing the strongest quasicrys-talline interaction, t [= t ( C n )]. The stronger t makesthe 8-waves interact over a much wider area in the Bril-louin zone, and the energy spacing between the quasicrys-talline states larger. Accordingly, the m = ± t , where the white dashedline corresponds to the DOS for the system consideredin Fig. 2. The large DOS observed at lower energiesin the systems with a large | t | reflects the flat bandsarising from the 2-wave mixing at Γ. While the eightBloch states centered at C n are sufficient to fully de-scribe the resonant interaction governing the quasicrys-talline states, some other minor interactions are not cap-tured by these bases. Thus, in order to calculate the DOSin wider energy range, we used more (32 waves) bases tocalculate the DOS in Fig. 3. We also plot the band edgesof the quasicrystalline states [Eq. (13)] by black dashedlines. As | t | increases, the energy spacing between theedges increases and the height of the DOS peaks alsogrows rapidly. It should be noted that some band edges(e.g., m = ±
2) lead to a series of characteristic spikypeaks in DOS and dips (pseudogaps) in between, whileother edges are buried in the DOS of weakly coupledstates. Thus, quasicrystalline features, such as local den-sity of states (LDOS) with 8-fold rotational symmetryand relevant physical properties, are most prominent atthe energies where the band edges coincide with the spikypeaks in DOS. As changing t does not break the symme-try of the Hamiltonian, it neither changes the symmetrynor the degeneracy of quasicrystalline states.In addition to the features from the quasicrystalline8-wave mixing, Fig. 3(b) shows the peaks and pseudo-gaps associated with other kinds of interaction. We plotthe wave numbers associated with these interactions inFig. 3(c), together with the Fermi surfaces. The middlepanel shows the 2-wave mixing between the states in dif-ferent layers , which occurs when the Fermi surfaces ofthe two layers meet, while the right panel shows the 2-wave mixing between the states in the same layer assistedby the potential of the opposite layer . Blue and redarrows in Figs. 2(a) and 3(a) show the band opening bythe interlayer and intralayer 2-wave mixing, respectively,whose size also increases with | t | . The interlayer interac-tion strength t ( q ) involved in the interlayer and intralayer2-wave mixing is 0.472 and 1.49 times the interactionstrength t for the 8-wave interaction. However, the in-tralayer mixing exhibits a band opening smaller thanthe interlayer mixing partly due to the two successiveinterlayer interaction and partly due to the energy differ-ence between the states in opposite layers. At t → E = cos √ π + 1( ≈ . − √ π ( ≈ . − π/ √ ≈ .
42) in unit of | V ppπ | . It should be noted that the states and band open-ing arising from these three mixing are continuously con-nected to each other in the Brillouin zone [Figs. 3(a) and(c)]. Unlike the quasicrystalline 8-wave interaction, boththe 2-wave mixing processes, labelled α and β in Fig. 3(c),can occur in bilayer square lattices stacked at any rota-tion angle θ , i.e., at usual moir´e superlattices. However,it is straightforward to show that α and β occur at differ-ent energies in the systems with θ other than 45 ◦ . Andat small θ , the interlayer mixing appears at the energiescloser to the band maximum at M than the intralayermixing, while the order is reversed as θ approaches 45 ◦ .
4. Effects of interlayer potential asymmetry
Figure 4(a) shows the band dispersion near C n of oc-tagonal vdW-QCs under three different interlayer asym-metric potential, V = 0, 0 . | V ppπ | , 0 . | V ppπ | . Since (d) s = 0 x [ a ]8-840-4 y [ a ] y [ a ] s = 2 x [ a ]8-840-4 y [ a ] y [ a ] s = ±1 x [ a ]8-840-4 y [ a ] y [ a ] | t | [ | V pp π | ] V [| V ppπ |] deg r ee o f m i x i ng [ % ] (c) V = 0 V = 0.1 V = 0.2 (a) E s − h [ | V pp π | ] k [1/ a ] V [| V ppπ |] (b) V /2 m s FIG. 4. (a) Band dispersion near C n of octagonal vdW-QCsunder three different interlayer potential asymmetry, V = 0,0 . | V ppπ | , 0 . | V ppπ | . (b) Band edges at C n with various V .Indices m and s show the angular momentum of the pristinequasicrystalline states with 8-fold rotational symmetry andthat of 4-fold rotational symmetry under interlayer potentialasymmetry. Dashed arrows show the interaction between theconstituent quasicrystalline states by H V . (c) Degree of mix-ing between the constituent states in Ψ s = ± (see text) withrespect to the interlayer interaction strength | t | and inter-layer potential asymmetry V . (d) Plots similar to Fig. 2(c)for V (cid:54) = 0. The top and bottom panels show the LDOS of theupper and lower bands, respectively. Eq. (10) satisfies Σ (cid:48)− ( H− h I )Σ (cid:48) = − ( H− h I ) at k = for Σ (cid:48) = diag( iσ y , iσ y , iσ y , iσ y ) and I is an 8 × H has a particle-hole symmetry with respectto the energy E = h . As V increases, however, thestates with m = ± k = lose the degeneracy, andall the band edges move away from E m = ± (= h ).We can obtain further insight on the effects of the inter-layer potential asymmetry from the analytic expressionof the energies at k = . The interlayer potential asym-metry couples the eigenstates of H ring ( k = ) that haveangular momenta that differ by ± (cid:104) v m (cid:48) |H V | v m (cid:105) = (cid:26) − V / m − m (cid:48) ≡ , , (14)since the diagonal elements of H V work as a staggeredpotential with 1/4 period of the ring in the dual tight-binding lattice. Thus, the Hamiltonian matrix in thebases of the quasicrystalline states is reduced to four 2 × H m,m (cid:48) = (cid:18) E m − V / − V / E m (cid:48) (cid:19) , (15)for ( m, m (cid:48) ) = (0 , , (1 , − , (2 , − , (3 , − v m , originate from the resonant inter-action between the degenerate states | k ( n ) (cid:105) in the two lay-ers, the interlayer potential asymmetry breaks the 8-foldrotational symmetry of the states by lifting the degen-eracy of | k ( n ) (cid:105) . This reduces the allowed angular quan-tum numbers to s (= 0 , ± , ≡ m ≡ m (cid:48) (mod 4), whichindicates a 4-fold rotational symmetry. We obtain thefollowing energies and wave functions E s = h ± (cid:113) t cos (5 πs/
4) + V / , Ψ s = c m v m + c m (cid:48) v m (cid:48) , (16)where ( c m , c m (cid:48) ) is (sin( φ/ , − cos( φ/ φ/ , sin( φ/ φ = tan − ( V / (4 t cos q m )), and plot E s against V inFig. 4(b). The states with s = 1 and s = − v ∗ m = v − m symmetry of thewave functions. At small V , the interlayer interaction | t | suppresses the energy shift of s = 0 , ± . On the other hand, the two stateswith s = 2 ( m = ±
2) are composed of two degeneratequasicrystalline states, m = 2 and m = −
2. Thus, theirband edges shift as much as the applied bias in oppositedirections and are not affected by the interlayer interac-tion t . As V increases, the overall energy span of theseresonant states increases, while the energy spacing be-tween the adjacent states decreases.The dashed arrows in Fig. 4(b) show the interactionbetween quasicrystalline states between v m and v m (cid:48) by H V , and Fig. 4(c) shows the degree of mixing in Ψ s = ± ,which we defined as (1 − || c m | − | c m (cid:48) | | ) ×
100 [%]. Sys-tems with | t | < V exhibit stronger mixing, which willinfluence the transition behavior, such as the optical se-lection rule, in vdW-QCs. The states with s = 0 exhibita similar, but slightly weaker, mixing owing to the largerenergy difference between E m and E m (cid:48) in Ψ s =0 . How-ever, the states with s = 2 are special in that the con-stituent states v and v − are always fully mixed, i.e., c m (cid:48) = − c m for the upper band and c m (cid:48) = c m for the lowerband, regardless of the values of t and V . Again, thisis due to the degeneracy between the constituent states m = 2 and m (cid:48) = − s = 0 , ± , V = 0 . | V ppπ | in the top and bottompanels in Fig. 4(d), respectively. Due to the interlayerpotential asymmetry, the wave functions Ψ s are more orless spatially polarized to either layer. And the strongerthe mixing, the more the wave functions are layer polar-ized; for example, Ψ s = ± exhibit more polarization than Ψ s =0 . This is becauseΨ s = 12 ( c m + c m (cid:48) )( v m + v m (cid:48) ) + 12 ( c m − c m (cid:48) )( v m − v m (cid:48) ) , (17)where ( c m , c m (cid:48) ) ∈ R , and v m + v m (cid:48) and v m − v m (cid:48) areperfectly polarized to layer 2 and 1, respectively, since µ m (cid:48) = − µ m . Thus, as the mixing becomes stronger, theupper bands ( c m (cid:48) ≈ − c m ) consist mostly of | k ( n ) (cid:105) witheven n (i.e., layer 1) while the lower bands ( c m (cid:48) ≈ c m )consist mostly of | k ( n ) (cid:105) with odd n (i.e., layer 2). Again,the states with s = 2 are special in that their wave func-tions Ψ s =2 are perfectly polarized to either layer regard-less of the values of t and V because the constituentstates v and v − are always fully mixed. This is similarto the case of an one-dimensional diatomic chain whosesublattices stop completely at the acoustic and opticalmodes.The C n in Fig. 1(c), which are associated with | k ( n ) (cid:105) ,remain the same since the interlayer potential asymmetrydoes not change the Umklapp scattering paths. Thus, theLDOS profile of each layer-polarized state, which is asso-ciated with C n for even n (layer 1) or C n for odd n (layer2), is exactly consistent with the profile of each layerin the absence of the potential asymmetry [Fig. 2(c)].Therefore, we can switch between the quasicrystallinestates and their layer components by applying an elec-tric field. C. Dodecagonal quasicrystal
1. Hamiltonian
In dodecagonal vdW-QCs, the strongest resonant in-teraction occurs at ˆ k = . Thus, by replacing ˆ k with k ,we can express the Hamiltonian of the resonant ring H ring in the absence of the interlayer and sublattice potentialasymmetry by a 24 ×
24 matrix H ring ( k ) = H (0) W † WW H (1) W † W H (2) W † . . . . . . . . . W H (10) W † W † W H (11) , (18) H ( n ) ( k ) = (cid:32) h ( n ) AA h ( n ) AB h ( n ) BA h ( n ) BB (cid:33) , W = − t (cid:18) ω ω ∗ (cid:19) , (19)in the bases of ( | k (0) (cid:105) , | ˜ k (1) (cid:105) , | k (2) (cid:105) , | ˜ k (3) (cid:105) , · · · , | k (11) (cid:105) ).Here, k ( n ) = k + Q n for even n (layer 1) and ˜ k ( n ) = k − Q n for odd n (layer 2), where | k ( n ) (cid:105) and | ˜ k ( n ) (cid:105) are( | k ( n ) , X (cid:105) , | k ( n ) , X (cid:48) (cid:105) ) and ( | ˜ k ( n ) , ˜ X (cid:105) , | ˜ k ( n ) , ˜ X (cid:48) (cid:105) ) with thesublattices X and X (cid:48) are arranged in the order of ( A, B )or ( ˜ A, ˜ B ) for n ≡ , B, A ) or ( ˜ B, ˜ A )for n ≡ ,
2. And h ( n ) X (cid:48) X ( k ) = h X (cid:48) X [ R ( − nπ/ k + Q ], ω = e πi/ , and we neglect the k dependence of the in-terlayer matrix element t ( q ).In the given bases order, the Hamiltonian representingthe interlayer and sublattice potential asymmetry are ex-pressed by H V = V I , − I , I , − I , I , − I , I , − I , I , − I , I , − I ) , (20)and H ∆ = ∆2 diag( σ z , − σ z , − σ z , σ z , σ z , − σ z , − σ z , σ z ,σ z , − σ z , − σ z , σ z ) , (21)respectively, where V represents the difference in the elec-trostatic energies between the two layers, ∆ is the differ-ence between the on-site potentials between two sublat-tices, and I is a 2 × H = H ring + H V + H ∆ . (22)
2. Band structures and wave functions
We plot the band structures near C n of the dodecago-nal vdW-QCs with ( V ppσ , d ) = ( − . V ppπ , . a ) and V = 0, ∆ = 0 in Fig. 5(a), and their closer view in (b).The twelve Dirac cones are arranged on a circle with aradius ∆ k = 4(2 − √ π/ (3 a ), and they are stronglyhybridized near k = with t = 0 . V ppπ to ex-hibit the characteristic dispersion including flat bandbottoms, the Mexican-hat edges, and the frilled bandedges. We can get the electronic structures of variousdodecagonal vdW-QCs by using the proper V ppπ ; e.g., V ppπ = − .
38 eV gives the spectrum of vdW-QC com-posed of two graphene layers, which is known as qua-sicrystalline twisted bilayer graphene.At k = , H ring can be analytically diagonalized toobtain a set of energies (neglecting the constant energy) E ± m = t cos q m ± (cid:113) t sin q m + ( h − t cos q m ) , (23)where h = h AB ( Q ) = h BA ( Q ) = − . V ppπ , and q m = (7 π/ m with m = − , − , · · · , , ± in Eq. (23), respectively.The energy scaling in conduction band is, however, muchsmaller than that in the valence band since the wavefunction of the conduction band of intrinsic graphene,having the same phases between the sublattices, sup-presses the interlayer interaction by a factor of 3. The in-dex m is a quantized angular momentum respecting the 12-fold rotational symmetry. The states with m = ± s ( s = 1 , , , ,
5) form twofold doublets, while the m = 0and 6 are non-degenerate. Note that the interaction re-sponsible for the formation of the quasicrystalline statesonly weakly affects the spectrum at energies away fromthe momentum matching conditions; e.g., in a quasicrys-talline twisted bilayer graphene there is no meaningfulchange on the Fermi velocity , since E ± m are far from theDirac point. Further detailed analysis on H ring can befound in Ref. .Figure 5(c) shows the LDOS of the quasicrys-talline states, where we can see that the wave am-plitude distribute selectively on a limited number ofsites in a characteristic 12-fold rotationally symmet-ric pattern. The wave functions for E ± m are v ± m =(1 / √ µ − m , µ − m , µ − m , · · · , µ m ) (cid:78) ( c ± m, , c ± m, ) ( µ m = e iq m ), where ( c + m, , c + m, ) = (sin( φ m / , cos( φ m / c − m, , c − m, ) = (cos( φ m / , − sin( φ m / φ m = tan − [( h − t cos q m ) / ( √ t sin q m )]. Since the Hamiltonian has asymmetry of Σ (cid:48)− H ring Σ = H ∗ ring Σ = diag( σ x , σ x , σ x , σ x , σ x , σ x ) (24)at k = , the states with angular momentum m and − m are degenerate and v ±− m = σ x ( v ± m ) ∗ , and it is straight-forward to show that their LDOS profiles are exactlythe same to each other. Figure 5(c) also shows that thestates with ± m exhibit LDOS profiles which look similarto those of 6 ∓ m ; v ± m clearly shows that the states with m = 0 and 6, and also the states with m = 3 and − m = ± ± m = ± ±
4) become different as | t | grows.Likewise, the ± m states in the conduction band exhibitLDOS profiles which look similar to the 3 ∓ m ones inthe valence band in the systems with a small | t /V ppπ | .
3. Effects of interlayer interaction strength and other kindsof interaction
Figure 6(a) shows the valence band structures of adodecagonal vdW-QC with a interlayer interaction t ( q )that is 2 times larger than the one in Fig. 5. The en-ergy spacing between the quasicrystalline states becomeslarger, and the flat band area of m = 6 ( m = ±
4) statein the valence band in Fig. 6(a) is approximately 2-times(5-times) as large as that in Fig. 5(a), and it is 28-times(70-times) bigger than the flat band area of magic-angletwisted bilayer graphene. As a greater number of theelectronic states are involved in the flat bands, we expectto see stronger electron-electron interacting effect. Thedensity map in Fig. 6(b) shows the DOS of dodecagonalvdW-QCs with various t calculated by using 182-wavebases. The white dashed line corresponds to the DOS0 ene r g y [ | V pp π | ] CB: m = ±3 / VB: m = 0, 6 y [ a ] x [ a ] CB: m = ±2, ±4 / VB: m = ±1, ±5 y [ a ] x [ a ] -0.2 k [1/ a ] -1.5 -1.0 -0.5 0 0.5 1.0 1.5-0.4-0.60.8-0.8 (a) (b)(c) ±5±36 m = 0 ene r g y [ | V pp π | ] k [1/ a ] ±2±4±1 m = 0±1±36±4±2±5 CB: m = 0, 6 / VB: m = ±3 y [ a ] x [ a ] CB: m = ±1, ±5 / VB: m = ±2, ±4 y [ a ] x [ a ] FIG. 5. Plots similar to Fig. 2 for dodecagonal vdW-QCscalculated by the 12-ring effective model. Blue and red arrowsin (a) show the band opening by the interlayer and intralayer2-wave mixing, respectively. The top and bottom panels in(b) show the quasicrystalline states in the conduction bandand valence band, respectively. for the system considered in Fig. 5, and the black dashedlines show the band edges of the quasicrystalline states[Eq. (23)].The systems with larger | t | exhibit higher DOS peaksowing to the increase of the flat band area in the k space.Not every quasicrystalline state leads to a DOS display-ing spiky peaks interspersed with pseudogaps, so qua-sicrystalline features would be most prominent at the min max k y [ / a ] k y [ / a ] k y [ a ] k x [ a ] 4-4 20-2 k x [ a ]4-4 20-2 k x [ a ] (c) k [1/ a ] -1.5 -1.0 -0.5 0 0.5 1.0 1.5 (a) (b) ene r g y [ | V pp π | ] ±5 ±36 m = 0 ±2 ±4±1 m = 0±1±36±4±2±5 |t | [| V ppπ |] α β ≈ FIG. 6. (a) Electronic structures of dodecagonal vdW-QCs with the interlayer interaction | t | t . The white dashed line correspondsto the DOS for the system considered in Fig. 5. (c) Plotssimilar to Fig. 3(c) in dodecagonal vdW-QCs. energies where the band edges coincide with the spikypeaks in DOS, especially at the m = 6 , ± , ± m = 6 , ± , ± t neither changes thesymmetry nor degeneracy of quasicrystalline states. Inmost practical parameter ranges, a system with larger | t | exhibits a larger energy spacing between the qua-sicrystalline states in both the conduction and valencebands. Note that, in the systems with extremely stronginterlayer interaction ( | t | > h / | t | increases(Appendix C). Such a condition, however, is hard to beachieved in the most practical systems. Thus, hereafter,we will consider the systems with | t | < h / E = − . | V ppπ | when t = 0, show to theband edges and pseudogaps arising from intralayer 2-wave mixing . We plot the band opening associatedthe interlayer and intralayer mixing as the blue and redarrows in Figs. 5(a) and 6(a), respectively, and visualizethese interactions in Fig. 6(c). The band opening by the12-wave mixing in the conduction band is much smallerthan that in the valence band for the same reason as theband opening via the 12-wave mixing (Sec. III C 2). Al-though the interlayer interaction strength t ( q ) involvedin the intralayer 2-wave mixing is about 1.50 times theinteraction strength in the interlayer 2-wave mixing, theintralayer mixing exhibits smaller band opening since itinvolves two successive interlayer interaction between thestates with different energies. At sufficiently large | t | ( > . | V ppπ | ), however, the intralayer interaction givesband opening throughout the entire Brillouin zone inthe valence band except in the vicinity of the quasicrys-talline states. This makes the quasicrystalline stateswith m = 6 , ± t → E = 0. Again, the states and band opening arisingfrom these three different mixings are continuously con-nected to each other in the Brillouin zone [Figs. 6(a) and(c)]. The 2-wave mixing can occur in bilayer honeycomblattices stacked at any rotation angle θ , and α and β oc-cur at different energies when θ (cid:54) = 30 ◦ . At small θ , theinterlayer mixing appears at energies closer to the bandedge at K (and K (cid:48) ) than the intralayer mixing, while theorder is reversed as θ approaches 30 ◦ .
4. Effects of interlayer potential asymmetry
Figures 7(a) and (d) show the dispersion in the con-duction band and valence band near C n of dodecagonalvdW-QCs under three different interlayer asymmetric po-tential, V = 0, 0 . | V ppπ | , 0 . | V ppπ | . Again, the stateswith m = ± k = lose their degeneracy, and all theband edges, in both the conduction band and valenceband, move away from E m = ± as V increases.At k = , the interlayer potential asymmetry couplesthe eigenstates of H ring , v bm ( b is + for conduction bandand − for valence band), having their angular momentadiffer by ± (cid:104) v b (cid:48) m (cid:48) |H V | v bm (cid:105) = (cid:26) − V ( c b (cid:48) m (cid:48) , c bm, + c b (cid:48) m (cid:48) , c bm, ) ( m − m (cid:48) ≡ , , (25)since the diagonal elements of H V work as a staggeredpotential with 1/6 period of the ring in the dual tight-binding lattice. Thus, the Hamiltonian matrix is reducedto 4 × H m,m (cid:48) = (cid:18) C m − V R ( φ m − φ m (cid:48) ) − V R − ( φ m − φ m (cid:48) ) C m (cid:48) (cid:19) , C m = (cid:18) E − m E + m (cid:19) , (26) in the bases of ( v − m , v + m , v − m (cid:48) , v + m (cid:48) ) for ( m, m (cid:48) ) =(0 , , (1 , − , (2 , − , (3 , − , (4 , − , (5 , − R ( φ ) is a rotation matrix. The electronic stateslose the 12-fold rotational symmetry, since H V liftsthe degeneracy of | k ( n ) (cid:105) in different layers, andare characterized by the angular quantum number s (= 0 , ± , ± , ≡ m ≡ m (cid:48) (mod 6) for a 6-fold rota-tional symmetry. We obtain two band edges E + s in theconduction band and another two E − s in the valenceband by diagonalizing Eq. (26). These are plottedagainst V in Figs. 7(b) and (e), for the conduction andvalence bands, respectively. It is straightforward toshow that the states with a quantum number s and − s ( s = 1 ,
2) are degenerate, since the reduced Hamiltoniansatisfies ˜Σ − ˜ H m,m (cid:48) ˜Σ = ˜ H − m, − m (cid:48) , ˜Σ = diag(1 , − , , − . (27)In most practical cases, Eq. (26) can be further reducedto two 2 × φ m − φ m (cid:48) ) / ≈ φ m ≈ φ m (cid:48) ≈ ◦ ).The dashed arrows in Figs. 7(b) and (e) show the in-teraction between quasicrystalline states between v ± m and v ± m (cid:48) by H V . The states in the conduction band exhibitlarger mixing between the constituent quasicrystallinestates than those in the valence band, due to the smallerenergy difference (not shown). And, similar to the octag-onal vdW-QCs [Fig. 4(c)], materials with weaker | t | un-der larger V experience larger energy shift, mixing, andaccordingly larger spatial layer-polarization because ofEq. (17) (Appendix D). Again, the states with s = 3 arespecial in that the constituent states v ± m =3 and v ± m = − are always fully mixed regardless of the values of t and V , due to the degeneracy between E ± m =3 and E ± m = − .We plot the LDOS of the states in the conduction bandand valence band with s = 0 , ± , ± , V = 0 . | V ppπ | in Figs. 7(c) and (f), respectively, where the top and bot-tom panels in each figure show the LDOS of the upperand lower bands, respectively. Again, the stronger themixing, the more the wave functions are layer polarized,and the wave functions Ψ ± s =3 are mostly polarized to ei-ther layer even at very weak V . The LDOS profile ofeach layer-polarized state is exactly consistent with theprofile of each layer in the absence of the potential asym-metry [Fig. 5(c)], since H V does not change the Umklappscattering paths.In dodecagonal vdW-QCs with sublattice symmetry(i.e., ∆ = 0), an interlayer potential asymmetry does notopen a gap at the Dirac point. This is because the co-existence of the time reversal symmetry and the in-place C (180 ◦ ) rotation symmetry requires vanishing of theBerry curvature at any nondegenerate point the the en-ergy band , and this guarantees the robustness of bandtouching points in two-dimensional systems . Just like2 (c) s = 0 s = ±1 s = ±2 s = 3 -8 -4 0 4 8 -8 -4 0 4 8 -8 -4 0 4 8 -8 -4 0 4 8 x [ a ] x [ a ] x [ a ] x [ a ] (f) s = 0 s = ±1 s = ±2 s = 3 -8 -4 0 4 8 -8 -4 0 4 8 -8 -4 0 4 8 -8 -4 0 4 8 x [ a ] x [ a ] x [ a ] x [ a ] k [1/ a ] -0.4 -0.2 0 0.2 0.40.650.750.700.600.55 E s [ | V pp π | ] -0.6-0.7-0.8-0.9 E s [ | V pp π | ] k [1/ a ] -0.4 -0.2 0 0.2 0.4 (a)(d) V [| V ppπ |] V [| V ppπ |] (b)(e) m ±16±4±30±5±2 s ±13±20±13±20 m ±50±2±36±1±4 s ±13±20±13±20 +− FIG. 7. (a) Conduction band dispersion near C n of dodecagonal vdW-QCs under three different interlayer potential asymmetry, V = 0, 0 . | V ppπ | , 0 . | V ppπ | . (b) Conduction band edges at C n with various V . Indices m and s show the angular momentumof the pristine quasicrystalline states with 12-fold rotational symmetry and that of 6-fold rotational symmetry under interlayerpotential asymmetry. Dashed arrows show the interaction between the constituent quasicrystalline states by H V . (c) Plotssimilar to Fig. 5(c) for V (cid:54) = 0 in the conduction band. The top and bottom panels show the LDOS of the upper and lowerbands, respectively. (d), (e), (f): Plots similar to (a), (b), (c) for valence band. twisted bilayer graphenes with any rotation angle , do-decagonal vdW-QCs composed of two honeycomb latticeswith ∆ = 0 has the C symmetry, even in the presenceof interlayer potential asymmetry because C does notflip the layers. Thus, the Dirac points of dodecagonalvdW-QCs with ∆ = 0 are protected even in the presenceof interlayer potential asymmetry.
5. Effects of sublattice potential asymmetry
We plot the conduction band and valence band near C n of a dodecagonal vdW-QC with sublattice potentialasymmetry of ∆ = − . V ppπ in Figs. 8(b) and (e), re-spectively, and plot the bands in the absence the asym-metry in (a) and (d) as a reference. As usual in mosttransition-metal dichalcogenides, ∆ (cid:54) = 0 in the currentmodel makes a band opening as large as ∆ at the energyrange centered at E = 0. Unlike the interlayer potentialasymmetry, however, we can see that breaking the sub-lattice symmetry does not make dramatic change to the band structures near the quasicrystalline states.The potential which breaks the sublattice symmetry, H ∆ , couples the eigenstates of H ring , v bm , whose angularmomenta differ by ± (cid:104) v b (cid:48) m (cid:48) |H ∆ | v bm (cid:105) = (1 − i )∆4 ( c b (cid:48) m (cid:48) , c bm, − c b (cid:48) m, c bm, ) ( m − m (cid:48) ≡ − , (1+ i )∆4 ( c b (cid:48) m (cid:48) , c bm, − c b (cid:48) m, c bm, ) ( m − m (cid:48) ≡ , , (28)since the diagonal elements of H ∆ work as a poten-tial with 1/3 period of the ring in the dual tight-binding lattice. Thus, H ∆ couples quasicrystallinestates with four different m , ( m , m (cid:48) , m , m (cid:48) ) =(0 , , , − , (1 , , − , − , ( − , , , − k [1/ a ] -0.4 -0.2 0 0.2 0.4 k [1/ a ] -0.4 -0.2 0 0.2 0.4 (b)(e) k [1/ a ] -0.4 -0.2 0 0.2 0.4 k [1/ a ] -0.4 -0.2 0 0.2 0.4 E s [ | V pp π | ] E s [ | V pp π | ] -0.8-0.7-0.6-0.5-1.0-0.9 (a)(d) Δ [| V ppπ |] Δ [| V ppπ |] (f)(c) m ±50±2±36±1±4 m ±16±4±30±5±2 FIG. 8. (a) and (b): Conduction band dispersion near C n of a dodecagonal vdW-QC with sublattice potential asymmetryof (a) ∆ = 0 and (b) − . V ppπ , respectively. (c) Band edges at C n with various ∆. Indices m show the angular momentumof the pristine quasicrystalline states with 12-fold rotational symmetry. (d), (e), (f): Plots similar to (a), (b), (c) for valenceband. tonian matrix is reduced to 8 × H m ,m (cid:48) ,m ,m (cid:48) = C m D m ,m (cid:48) D † m (cid:48) ,m D † m ,m (cid:48) C m (cid:48) D m (cid:48) ,m D † m (cid:48) ,m C m D m ,m (cid:48) D m (cid:48) ,m D † m ,m (cid:48) C m (cid:48) , D m a ,m b = (1 + i )∆4 (cid:18) cos ( φ m a + φ m b ) sin ( φ m a + φ m b )sin ( φ m a + φ m b ) − cos ( φ m a + φ m b ) (cid:19) (29)in the bases of ( v − m , v + m , v − m , v + m , v − m (cid:48) , v + m (cid:48) , v − m (cid:48) , v + m (cid:48) ).The electronic states lose the 12-fold rotational symme-try, and are characterized by the angular quantum num-ber s (= 0 , ± ≡ m i , (mod 3) ( m i = m , m , m (cid:48) , m (cid:48) ) fora 3-fold rotational symmetry. Again, the states with aquantum number s = 1 and s = − × E ± s in the unitof | V ppπ | plotted against ∆. Unlike the interlayer po-tential asymmetry, ∆ merely shifts the energies slightlyaway from the Dirac point and does not make dramaticchange to the quasicrystalline states, which is consistentwith the band structures in Figs. 8(b) and (f). This isbecause, | c bm, | ≈ | c bm, | in most practical systems with | t | < h / φ m ≈ ◦ ) at this high energy regime. Thus,the potential, which has the opposite sign between thesublattices, is almost cancelled in the intraband interac-tion ( b = b (cid:48) ), due to the phase cancellation. On the otherhand, although materials with much higher | t | have fi-nite contribution from the sublattice phases in the in-traband interaction, the overall interaction is still very4weak since the energy difference between the quasicrys-talline states increases as | t | grows. In any system, theinterband interaction ( b (cid:54) = b (cid:48) ) is always negligible dueto the large energy difference. Thus, dodecagonal vdW-QCs composed of transition metal dichalcogenides willalso exhibit the quasicrystalline states and DOS analo-gous to the quasicrystals composed of bilayer graphene.If the two layers have different ∆, other than a simplesign difference, the degeneracy of the m = ± IV. CONCLUSIONS
We investigated the electronic structures of quasicrys-tals composed of the incommensurate stack of atomiclayers (vdW-QCs) with various rotational symmetry andvarious interlayer interaction strength. We found the res-onant coupling of the monolayer states which gives qua-sicrystalline order in each quasicrystal symmetry, andcalculate the quasi-band dispersion and wave functionsrespecting the rotational symmetry of quasicrystal. Al-though the quasicrystalline states coexist in energy withweakly coupled states (e.g., the states from 2-wave mix-ing) in general, we showed that some quasicrystallinestates, which are usually obscured by additional weaklycoupled states, are more prominent in quasicrystals withstrong interlayer interaction. Besides, we showed that wecan switch the states between quasicrystalline configura-tion and its layer components, and also mix the statesby lifting the layer degeneracy, e.g., by applying elec-tric field. We analytically showed that the states in themiddle are always fully mixed and 100% layer polarizedregardless of the magnitude of the interlayer asymme-try. On the other hand, hexagonal lattices with sublat-tice potential asymmetry, e.g., transition metal dichalco-genide or hexagonal Boron Nitride, give the quasicrys-talline states quite similar to those in the system in theabsence of sublattice asymmetry. This is the first theoret-ical work which investigated the formation of quasicrys-talline states in vdW-QCs for every possible rotationalsymmetry and sublattice symmetry in two-dimensionalsystem, which will lead to extended exploration of richquasicrystal physics in designer quasicrystals.
Note added.
During the completion of this work, webecame aware of recent theoretical works on the pres-sure and electric field dependence of quasicrystalline elec-tronic states in 30 ◦ twisted bilayer graphene . ACKNOWLEDGMENTS
This work was supported by Science and Tech-nology Commission of Shanghai Municipality grantno. 19ZR1436400, and NYU-ECNU Institute of Physicsat NYU Shanghai. This research was carried out on theHigh Performance Computing resources at NYU Shang-hai.
Appendix A: Quasicrystalline states with weakerinteraction
As mentioned in Sec. III A, the sets of the waves inFigs. 1(b) and (e) are not the only set of states whichshow the resonant coupling in each system. We can findmore sets of states, with different wave numbers, showingthe resonant interaction respecting the rotational sym-metry of the quasicrystals. For example, the eight states k (red) and ˜ k (cid:48) (blue) in Fig. 9(a) also form a circularchain in the dual tight-binding lattice. Note that ˆ k forthese states (= ) is different from that for the statesin Fig. 1(b). These states are mapped to k (red) and˜ k (blue) in the first Brillouin zone, experience a reso-nant interaction, and form quasicrystalline states. Figure9(c) shows the band structures near the quasicrystallinestates arising from these eight states. It should be notedthat, however, the strength of the resonant interaction, | t ( q ) | , for the states in Fig. 9(a) is much weaker thanthat for the states in Fig. 1(b). This is because | t ( q ) | decays fast as | q | grows, and the former states have thechain with a longer segment length (= | q | ). Thus, theband opening in Fig. 9(c) is much smaller than that inFigs. 2(a) and (b). Dodecagonal vdW-QCs also havemore sets of states showing the resonant interaction. Inmost systems, however, such states can be mostly ne-glected since their interaction strengths are very weak,and they are also mixed with other types of interaction(e.g., 2-wave mixing). Thus, the sets in Figs. 1(b) and (e)give the strongest interaction, i.e., largest energy separa-tion and clear quasicrystalline order, since these statesform the rings with the shortest distance between neigh-boring states in the dual tight-binding lattices. Appendix B: Particle-hole symmetry of theHamiltonian of octagonal quasicrystal
By considering only the nearest neighbor pairs in theintralayer interaction, the ring Hamiltonian H ring of theoctagonal vdW-QCs [Eq. (11)], up to the first order to k = ( k ,x , k ,y ), can be transformed to5 ene r g y [ | V pp π | ] k [1/ a ] -0.03 -0.02 -0.01 0 0.01 0.02 0.031.101.021.041.061.08 (c) k x [1/ a ] -15 -5 0 10 15-10 5 k y [ / a ] -15-1001015-55 (a) k x [1/ a ] -4 -2 0 2 4 k y [ / a ] -4-2024 (b) FIG. 9. (a) and (b) Plots similar to Figs. 1(b) and (c) for the next strongest quasicrystalline interaction in octagonal vdW-QCs,where ˆ k = . (c) Electronic structures near the second dominant quasicrystalline states of octagonal vdW-QCs calculated bythe 8-ring effective model. U − H ring ( k ) U = H (cid:48) ring ( k ) = H (cid:48) ( − C ∗ CC H (cid:48) ( − C ∗ C H (cid:48) ( − C ∗ . . . . . . . . . C H (cid:48) (3) C ∗ C ∗ C H (cid:48) (4) (B1)with a transformation matrix U = ( v − , v − , v − , · · · , v , v ) , (B2)where v m = √ ( µ − m , µ − m , · · · , µ m , µ m ) T ( µ m = e i πm/ )is the eigenstate of the quasicrystalline state with aquantized angular momentum of m , H (cid:48) ( m ) = h − t cos(5 πm/ C ( k ) = sin( √ π ) aV ppπ ( k ,x − ik ,y ). Then, it is straightforward to show that H (cid:48) ring has a particle-hole symmetry with respect to the energy E = h , Σ − ( H (cid:48) ring − h I )Σ = − ( H (cid:48) ring − h I ) , Σ = σ z σ z σ z σ z , (B3)where I is an 8 × k ) is an eigenstate of H (cid:48) ring with an en-ergy of h + E , then Σ − Ψ( k ) is an eigenstate of energy h − E . Appendix C: Band edges of quasicrystalline states ofdodecagonal vdW-QCs with various t Figure 10 shows the band edges of quasicrystallinestates of dodecagonal vdW-QCs [Eq. (23)] with various interlayer interaction strength t up to vary strong in-teraction regime. We can see that the energy spacingbetween the quasicrystalline states increases as | t | in-creases in most practical interaction strength, i.e., | t | Figures 7(b) and (e) show that the states with s = 0determine the energy span of the resonant states in thepresence of the interlayer potential asymmetry, in boththe conduction band and valence band. Equation (26)shows the coupling between the quasicrystalline stateswith angular momentum m by interlayer potential asym-metry. In most practical cases, the matrix can be furtherreduced to two 2 × H ± m,m (cid:48) = (cid:18) E ± m − V / − V / E ± m (cid:48) (cid:19) , for the conduction band ( ˜ H + m,m (cid:48) ) and valence band( ˜ H − m,m (cid:48) ), since the interaction between the state in theconduction band and the state in valence band is almostnegligible due to the large energy difference. For the cou-6 FIG. 10. Band edges of quasicrystalline states of dodecago-nal vdW-QCs [Eq. (23)] with various interlayer interactionstrength t . The black dashed line corresponds to the bandedges for the system considered in Fig. 5, and the numbersshow the quantized angular momentum m . pled states with s = 0 (i.e., m = 0 and m (cid:48) = 6), the in-teraction between the conduction band and valence band is completely forbidden due to the sublattice symmetry.Then, we get E − s =0 = − h ± (cid:112) t + V / E + s =0 = h ± (cid:112) t + V / | h | > | t | .Thus, the states in the conduction band exhibit smallerenergy span than those in the valence band.The wave functions of the higher energy statesin both the conduction band and valence bands areΨ ± s =0 = sin( ˜ φ/ v ± m =0 + cos( ˜ φ/ v ± m =6 , and the lowerenergy states are Ψ ± s =0 = cos( ˜ φ/ v ± m =0 − sin( ˜ φ/ v ± m =6 ,where ˜ φ is tan − ( − V / (6 t )) for valence band andtan − ( V / (2 t )) for conduction band. As ˜ φ becomesclose to 90 ◦ , i.e., in materials with smaller | t | and | V | >> | t | , Ψ ± s =0 becomes (1 / √ v ± + v ± ) forthe upper state and (1 / √ v ± − v ± ) for the lowerstate. Since v ± = (1 / √ , , , , · · · , (cid:78) (1 , ± v ± = (1 / √ − , , − , , · · · , (cid:78) (1 , ± ± s =0 becomes polarized to either layer, i.e., the state is mostlycomposed of the Bloch bases with n of even (layer 1) orodd (layer 2) numbers, in the systems with small inter-layer interaction strength | t . On the other hand, thestates with s = 3 are always 100% polarized to eitherlayer. ∗ Corresponding author: [email protected] P. Stampfli, Helv. Phys. Acta , 1260 (1986). S. J. Ahn, P. Moon, T.-H. Kim, H.-W. Kim, H.-C. Shin,E. H. Kim, H. W. Cha, S.-J. Kahng, P. Kim, M. Koshino,Y.-W. Son, C.-W. Yang, and J. R. Ahn, Science , 782(2018). T. Suzuki, T. Iimori, S. J. Ahn, Y. Zhao, M. Watanabe,J. Xu, M. Fujisawa, T. Kanai, N. Ishii, J. Itatani, et al. ,ACS Nano , 11981 (2019). Y. Takesaki, K. Kawahara, H. Hibino, S. Okada, M. Tsuji,and H. Ago, Chemistry of Materials , 4583 (2016). W. Yao, E. Wang, C. Bao, Y. Zhang, K. Zhang, K. Bao,C. K. Chan, C. Chen, J. Avila, M. C. Asensio, et al. , Proc.Natl. Acad. Sci. , 6928 (2018). X.-D. Chen, W. Xin, W.-S. Jiang, Z.-B. Liu, Y. Chen, andJ.-G. Tian, Adv. Mater. , 2563 (2016). S. Pezzini, V. Miseikis, G. Piccinini, S. Forti, S. Pace,R. Engelke, F. Rossella, K. Watanabe, T. Taniguchi,P. Kim, and C. Coletti, Nano Lett. , 3313 (2020). P. Moon, M. Koshino, and Y.-W. Son, Phys. Rev. B ,165430 (2019). D. Xiao, G.-B. Liu, W. Feng, X. Xu, and W. Yao, Phys.Rev. Lett. , 196802 (2012). J. C. Slater and G. F. Koster, Phys. Rev. , 1498 (1954). P. Moon and M. Koshino, Phys. Rev. B , 195458 (2012). G. Trambly de Laissardi`ere, D. Mayou, and L. Magaud, Nano Lett. , 804 (2010). E. J. Mele, Phys. Rev. B , 161405 (2010). R. Bistritzer and A. H. MacDonald, Proc. Natl. Acad. Sci. , 12233 (2011). M. Koshino, New J. Phys. , 015014 (2015). P. Moon and M. Koshino, Phys. Rev. B , 205404 (2013). M. Koshino, P. Moon, and Y.-W. Son, Phys. Rev. B ,035405 (2015). P. Moon, Y.-W. Son, and M. Koshino, Phys. Rev. B ,155427 (2014). Note that the ratio between V ppσ and V ppπ of the hexagonallattices used in this work is different from that used in theprevious works ( V ppσ ≈ . 48 eV and V ppπ ≈ − . and graphene on hexag-onal boron nitride . In this work, we scaled V ppπ by a fac-tor of 1 . 25, while keeping V ppσ ≈ . 48 eV, to compensatethe deviation of the Fermi velocity of a pristine graphenedue to the summation over sites in the hopping range andmake the band dispersion and the energies of the van Hovesingularities consistent with the experimental results . M. Koshino, Phys. Rev. B , 115409 (2013). G. Yu, M. I. Katsnelson, and S. Yuan, Pressure and electricfield dependence of quasicrystalline electronic states in 30 ◦ twisted bilayer graphene, arXiv preprint arXiv:2003.11879(2020). P. Moon and M. Koshino, Phys. Rev. B90