Trigonal quasicrystalline states in \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$30^\circ$$\end{document}30∘ rotated double moiré superlattices

We study the lattice configuration and electronic structure of a double moiré superlattice, which is composed of a graphene layer encapsulated by two other layers in a way such that the two hexagonal moiré patterns are arranged in a dodecagonal quasicrystalline configuration. We show that there are between 0 and 4 such configurations depending on the lattice mismatch between graphene and the encapsulating layer. We then reveal the resonant interaction, which is distinct from the conventional 2-, 3-, 4-wave mixing of moiré superlattices, that brings together and hybridizes twelve degenerate Bloch states of monolayer graphene. These states do not fully satisfy the dodecagonal quasicrystalline rotational symmetry due to the symmetry of the wave vectors involved. Instead, their wave functions exhibit trigonal quasicrystalline order, which lacks inversion symmetry, at the energies much closer to the charge neutrality point of graphene.

A special case occurs when two hexagonal lattices are overlapped at θ = 30 • (Fig. 1a). In this instance the atomic arrangement is mapped on to a quasicrystalline lattice, which is ordered but not periodic, with a 12-fold rotational symmetry [20][21][22][23][24][25][26][27][28] . Owing to the momentum mismatch 21 , quasicrystalline twisted bilayer graphene exhibits the electronic structure of almost decoupled bilayer graphene at most energy ranges. Nevertheless, it also hosts unique electronic states which satisfy the 12-fold rotational symmetry 22,27,28 . Such quasicrystalline states arise from the resonant interaction between the states at specific wave vectors via the rotational symmetry of the quasicrystal as well as the translational symmetry of the constituent atomic layers 22, 28 . The red and blue hexagons in Fig. 1b show the first Brillouin zones of the two lattices. The numbered points and dashed lines show the wave vectors of the constituent monolayer states and the interlayer interaction which form the quasicrystalline resonant states. Such quasicrystalline states exhibit a wave amplitude distributed selectively on a limited number of sites in a characteristic 12-fold rotationally symmetric pattern (Fig. 1c). These states, however, appear at the energies (about ±1.7 eV) -far from the charge neutrality point of graphene. Similar quasicrystalline resonant states also arise in any bilayer stacked in a quasicrystalline configuration if all the dominant interlayer interactions occur between the atomic orbitals that have the same magnetic quantum number 28 . Thus, even transition metal dichalcogenides (TMDC) or square lattices can show the quasicrystalline states.
Recently, rapid progress has been made in stacking more than two incommensurate atomic layers and a number of studies have investigated the effects of multiple moiré superlattice potentials on the electronic structure. The most notable example among them is a double moiré system, which is composed of a graphene layer encapsulated by hBN layers (BN/G/BN) 29, 30 . The lattice mismatch between graphene and hBN results in a hexagonal moiré superlattice potential with a superlattice period L M i ( i = 1, 2 ) that can be as long as 14 nm (Fig. 2a). Such a long period [which results in short superlattice reciprocal lattice vectors, Fig. 2b] carves the graphene electronic structures into superlattice bands with an energy scale much smaller than that of pristine graphene (Fig. 2c) 13 . Recently, Leconte and Jung show that BN/G/BN at specific configurations can host two hexagonal moiré patterns overlaid at a twist angle of 30 • , and claimed that the system hosts quasicrystalline electronic structure 31    The wave vectors of the twelve monolayer states C n ( n = 0, 1, . . . , 11 ) which hybridize to quasicrystalline resonant states. The red and blue hexagons, and the red and blue arrows [ a * i and ã * i ( i = 1, 2 )] represent the first Brillouin zones and the reciprocal lattice vectors of each layer. Due to the symmetry, these twelve states are all degenerate in energy, and the dashed lines show that these twelve states interact by the reciprocal lattice vectors of the two layers. Note that these twelve states are centered around the Ŵ point. (c) Local density of states of the quasicrystalline resonant states. The area of the circle is proportional to the squared wave amplitude, and red and blue circles represent the states in the upper and the lower layers, respectively.  13 . The black and red hexagons represent the unit cells of graphene and hBN, respectively, and θ shows the relative orientation. The green hexagons and the vectors L M i ( i = 1, 2 ) represent the unit cells and superlattice vectors of the moiré superlattice, respectively. Here, the lattice constant of hBN is drawn 15% larger than that of graphene to enhance the visibility of the moiré pattern (actual difference is about 1.79%). Inset shows the lattice configuration of graphene; the black and white circles represent the A and B sublattices, a i and τ X show the primitive lattice vectors and the coordinates of the sublattices, respectively. (b) Superlattice Brillouin zone (blue hexagon) near the Dirac point (the region surrounded by blue lines in the inset) and the reciprocal lattice vectors G M i of graphene on hBN. X and Y show the Brillouin zone corners where mini Dirac point appear, and φ shows the relative orientation of G M 1 to the reciprocal lattice vector a * 1 of pristine graphene. Inset shows the first Brillouin zone of graphene, where the black and white circles represent the three equivalent Dirac points, K and K ′ , respectively. (c) The band dispersion of the first two bands in the conduction and valence bands of graphene on hBN with θ = 0 • , which show the band opening at the primary and the mini Dirac points. www.nature.com/scientificreports/ interactions reconstruct the band dispersion of pristine graphene at these wave vectors. Compared to the resonant states of quasicrystalline twisted bilayer graphene where the actual atomic lattices are arranged in a dodecagonal configuration 22,27,28 , the resonant states of BN/G/BN appear at the energies much closer to the charge neutrality point of graphene. However, their wave functions show the quasicrystalline order with a three-fold rotational symmetry rather than fully satisfying the 12-fold rotational symmetry of the double moiré pattern.

Methods
Hexagonal moiré superlattices stacked at 30 • . We consider a trilayer system composed of graphene sandwiched by hBN. Both graphene and hBN are two-dimensional honeycomb lattice whose unit cell comprises of two (A and B) sublattices. Graphene has carbon atoms in both sublattices, while hBN has nitrogen atom on A site and boron atom on B site. The lattice constant of hBN, ã ≈ 0.2504 nm 32 , is slightly larger than that of graphene, a ≈ 0.246 nm , and we use a constant interlayer distance of d = 0.322 nm between the adjacent two layers 33 . Here, we do not consider the lattice relaxation between graphene and hBN 34,35 , since the effects of such relaxation on the electronic structures is an order of a few meV. Nevertheless, our effective theory that respects the lattice symmetry is able to properly describe both the gap at the primary Dirac point and the asymmetric gap opening at the two inequivalent mini Dirac points 13 , as well as the orbital magnetism of the structure 19 .
We define the atomic structure of the double moiré superlattices by starting from a nonrotated arrangement, where the hexagon center of the three layers share the same in-plane position (x, y) = (0, 0) , and the A-B bonds are parallel to each other. We choose a 1 = a(1, 0) and a 2 = a(1/2, 2 ) ], where the upper and lower signs are for the top and bottom layers, respectively, represent the coordinates of the nitrogen and boron atoms in the unit cell (Note that we defined the sublattice coordinates τ X ( X = A, B ) and τX(l) ( X = N, B , l = t, b ) differently from those in our previous work 13 , to make the points with the highest rotational symmetry, i.e., the hexagonal center, as the center of the system. However, both definitions keep the interaction matrices (Eq. 9) the same.). We define the reciprocal lattice vectors a * i and ã * i for graphene and hBN, respectively, so as to satisfy a i · a * j =ã i ·ã * j = 2πδ ij . We then rotate the top and bottom hBN layers with respect to graphene by arbitrary angles θ (t) and θ (b) around the origin, respectively. From now on, we use "BN/G/BN" for this configuration only. Due to the symmetry of the lattice, 0 ≤ θ (l) ≤ 30 • ( l ∈ t, b ) spans all the independent configurations. Figure 3a shows the moiré interference patterns which arise from the lattice mismatch between the top hBN and graphene (left side), and also that from the bottom hBN and graphene (right side), respectively, and  36 . We plot |L M i | and |G M i | against θ in Fig. 3c in red and blue lines, respectively. Now, we will find the configuration where the unit cells of the two hexagonal moiré superlattices have the same size and are overlaid with a relative twist angle of 30 • . In such a configuration, the overlaid two hexagonal superlattices are mapped onto a 12-fold rotationally symmetric quasicrystalline lattice without any translational symmetry, as first shown by Stampfli 37 . From Eq. (1), the former and the latter conditions give Figure 3d shows φ as a function of θ for various ε . The red, green, blue lines correspond to ε > 0 , ε = 0 , ε < 0 , respectively, and the thick black line corresponds to hBN. The two hexagonal moiré superlattices form a dodecagonal quasicrystalline configuration at θ where the line and the dashed horizontal lines cross. If the lattice constant of the top and bottom layers is the same as that of the middle, graphene layer, i.e., ε = 0 , then the two hexagonal moiré patterns cannot have a relative twist angle of 30 • . On the other hands, the systems with ε < 0 , 0 < ε < 0.0353 , and ε ≥ 0.0353 have three, four, and two θ which satisfy the conditions. When the top and bottom layers are hBN, i.e., ε ≈ 0.0179 , the four angles are , are arranged in 12-fold rotational symmetry (Fig. 4a), just like the reciprocal lattice vectors in quasicrystalline twisted bilayer graphene that give rise to the resonant states ( Fig. 1b) 22 .
It should be noted that, however, although the overlap of the two moiré interference patterns are mapped onto a quasicrystalline tiling with 12-fold rotational symmetry, the actual lattice structure belongs to the symmetry group D 3 ; it is invariant under C 3 rotation about the axis perpendicular to the xy-plane and under three C 2 rotation about the axes in the plane, but lacks inversion symmetry. If we replace the top and bottom hBN layers by a material having the same types of atoms in both sublattices, then the lattice has the symmetry group D 6 which is still lower than the 12-fold rotational symmetry. (1) hBN ) represent the Hamiltonian for the intrinsic monolayer graphene and the top (bottom) hBN, respectively, U (t) ( U (b) ) is for the interlayer coupling between the graphene and the top (bottom) hBN. However, since the hBN electronic bands are far from the charge neutrality point of graphene, we can project the total Hamiltonian on to the Bloch bases of graphene p z orbitals at each sublattice, where |R X � is the atomic orbital at the site R X = n 1 a 1 + n 2 a 2 + τ X ( n i ∈ Z , X = A, B ), k is the two-dimensional Bloch wave vectors and N = S tot /S is the number of the graphene unit cells with an area S = ( √ 3/2)a 2 in the total system area S tot . Then, the Hamiltonian near the Dirac point K ξ = −ξ(2a * 1 + a * 2 )/3 of graphene, where ξ = ±1 for K and K ′ , respectively, is reduced to a 2 × 2 form 13 (Note that the theoretical model can be easily expanded to the systems with the encapsulating layers other than hBN; e.g., it is straightforward to expand the model to the atomic layers of which energy bands are close to the charge neutrality point of graphene by explicitly using the Bloch bases for those layers.),  where L = n 1 a 1 + n 2 a 2 , τ X ′ X = τ X ′ − τ X , and is the transfer integral between two p z orbitals at a relative vector R 5,38 , V 0 ppπ ≈ −3.38 eV (Note that the value of V 0 ppπ ≈ −3.38 eV used in this work is different from that used in the previous works ( V 0 ppπ ≈ −2.7 eV ) on the twisted bilayer graphene 6,13,39 and graphene on hexagonal boron nitride 13 . In this work, we scaled V 0 ppπ by a factor of 1.25 to compensate the deviation of the Fermi velocity of a pristine graphene, which affects the entire energy scale, due to the summation over sites in the hopping range.), V 0 ppσ ≈ 0.48 eV , and r 0 ≈ 0.0453 nm 6,39 . The effective potentials by hBN to graphene, V (l) hBN , are explicitly written as 13 (Note that we defined the sublattice coordinates τ X ( X = A, B ) and τX(l) ( X = N, B , l = t, b ) differently from those in our previous work 13 , to make the points with the highest rotational symmetry, i.e., the hexagonal center, as the center of the system. However, both definitions keep the interaction matrices (Eq. 9) the same.) where we truncated much weaker terms O(u 4 0 ) which are associated with longer momentum difference. Here, and (5) is φ rotated from a * 1 (Eq. 1). where ω = e 2πi/3 , and u 0 ≈ −t(K ξ ) ≈ 0.152 eV is the in-plane Fourier transformation of the transfer integral between two p z orbitals (Eq. 7) at q near the Dirac point 5, 38 . We will discuss more about u 0 later. By using V C = 0 , V N = −1.40 eV , and V B = 3.34 eV , as the on-site potential of carbon, nitrogen, and boron atoms, respectively 40 , we get V 0 ≈ 0.0289 eV , V 1 ≈ 0.0210 eV , and ψ ≈ −0.29 (rad). If we replace the top and bottom hBN layers by a material having the same types of atoms in both sublattices, i.e., if V N = V B , then ψ ≡ π/3 (mod π) . The symmetry of such a structure increases to D 6 , and the reduced Hamiltonian H gains the inversion symmetry It is straightforward to show that the reduced Hamiltonian H spans the subspace for any k 0 in the momentum space. To investigate the electronic structures near k 0 , for any practical calculation, we only need a limited number of bases around k 0 inside a certain cut-off circle k c , because the interaction with the states far from k 0 is very weak due to multiple scattering. We can, then, obtain the energy eigenvalues at all the wave vectors in Eq. (13) by diagonalizing the Hamiltonian matrix within the finite bases, and the quasiband dispersion of the system by plotting the energy levels against k 0 . Here the wave number k 0 works like the crystal momentum for the periodic system, and so it can be called the quasicrystal momentum.

Quasicrystalline resonant interaction by two moiré superlattice potentials.
In addition to the typical 2-and 3-wave interaction by each moiré superlattice ( Fig. 2c) 13 , the 12-fold rotational symmetry of the wave vectors which couple the monolayer states (Eq. 2) as well as the translational symmetries of the two moiré superlattices enables a unique interaction between twelve degenerate monolayer states. Such a resonant coupling occurs at the twelve symmetric points shown in Fig. 4b, where θ n = 5π 12 + 7nπ 6 − φ i ( n = 0, 1, 2, · · · , 11 ). While the twelve waves which constitute the resonant states in quasicrystalline twisted bilayer graphene are centered around the Ŵ point 28 , the waves involved in the resonant interaction in BN/G/BN are centered around the Dirac point K ξ . These twelve states are degenerate if we ignore the small trigonal warping in this low energy regime. We see that the states at C i strongly interact with the states at C i−1 and C i+1 by the reciprocal lattice vectors of the top and bottom moiré superlattices. The interaction with states at any other k can be safely neglected since the interaction strength is much less or the two states are not degenerate in most cases. Hence, these states form one-dimensional monatomic chain with twelve sites and two pseudospins whose interaction is described by the moiré potential (Eq. 5).
It should be noted that this is not the only resonant coupling in this system. As shown in previous work on quasicrystalline twisted bilayer graphene (e.g., Appendix A in Ref. 28 ), there are more sets of states, with different wave numbers, that show the resonant interaction between the constituent monolayer states. However, the set in Fig. 4b is associated with the strongest interaction |t(q)| and, hence, gives the largest energy separation between the hybridized states.

Results and discussion
By using the Bloch bases ( |k (0) � , |k (1) � , · · · , |k (11) � ) near the twelve wave vectors k (n) = C n + k 0 , where |k (n) � is (|k (n) , A�, η|k (n) , B�) with η = ω ξ ×floor(n/4) , we can express the Hamiltonian of the resonant ring H ξ ring by a 24 × 24 matrix Here (11) t(q) = 1 S T(r + zX X e z )e −iq·r dr  (15) is symmetric under rotation by four span of the ring (i.e., moving C n to C n+4 ), which actually corresponds to the [R(7π/6)] 4 ( 120 • rotation) of the entire system. This means that the resonant states of the BN/G/BN have a three-fold rotational symmetry, which is much lower than the 12-fold rotational symmetry of the double moiré patterns of the system. This Hamiltonian cannot obtain a six-fold rotational symmetry, even if we replace the top and bottom hBN layers by a material having the same types of atoms in both sublattices (e.g., V N = V B ), since (i) atomic structure lacks 12-fold rotational symmetry and (ii) the 12 wave vectors involved are centered at K ξ which has three-fold rotational symmetry (Fig. 2b). Thus, the relevant terms cannot be gauged out by a similarity transformation. There is, however, an exception in the systems with V N = V B ; the wave functions at k 0 = 0 , and only at this k 0 , show a six-fold rotational symmetry ( Fig. 6d-f).
The interlayer interaction strength u 0 (Eq. 10) deviates from −t(K ξ ) as the distance between k (n) and K ξ increases. This effect becomes obvious in the system with a longer |G M,(l) i | , e.g., in BN/G/BN with θ 4 . However, what is more important is that, the u 0 associated with the interaction between the neighboring Bloch states |k (n) � ( n = 0, 1, · · · , 11 ) [dashed lines in Fig. 4b] are not the same. Since the Fourier transformation of the transfer integral between p z orbitals, t(q) (Eq. 11), are isotropic along the in-plane direction, t(q) depends only on |q| . If any electronic state is mainly comprised of three monolayer states of which waves vectors are arranged in a threefold rotationally symmetric way around K ξ , e.g., K ξ + k , K ξ + R(2π/3)k , K ξ + R(4π/3)k , then the |q| associated with the interaction between the monolayer states are identical since we can freely choose one among the three equivalent K ξ , i.e., K ξ , R(2π/3)K ξ , R(4π/3)K ξ , for each state. This is the case that happens at the mini Dirac point (the hexagonal Brillouin zone corners) of graphene on hBN 13 . On the contrary, t(q) for each interaction of the resonant coupling in BN/G/BN are not identical, since there are more than three states involved. As a result, u 2 0 in Eq. (10) varies ±2, 3% for θ 1 , θ 2 , θ 3 , and ±11% for θ 4 . Nevertheless, since the twelve wave vectors are arranged in a 12-fold rotationally symmetric way, equally spaced triplets (i.e., n ∈ {0, 4, 8} , n ∈ {1, 5, 9} , n ∈ {2, 6, 10} , n ∈ {3, 7, 11} ) satisfy the three-fold rotational symmetry with K ξ . As a result, the oscillation of u 2 0 is consistent with the rotational symmetry of H ξ ring and does not reduce the three-fold symmetry further. Hereafter, we will consider the structures with θ 1 , θ 2 , θ 3 only, and use u 0 = −t(K ξ ) ≈ 0.152 eV , i.e., ignore the variation of u 0 .
Since a hBN layer lacks the inversion symmetry, it is natural to ask whether the band structure changes if we rotate one of the hBN layers by 180 • (We label this structure BN/G/NB). The BN/G/NB also belongs to the symmetry group D 3 , while the three C 2 axes are 60 • rotated from those of BN/G/BN. It is well known that twisted double bilayer graphene [41][42][43][44] and BN/G/BN at general angles 31 show the change of electronic structures with respect to such a change. The 180 • rotation of hBN corresponds to the swap of the boron and nitrogen atoms. Thus, we can get the effective potential Eq. (8) of the moiré superlattice from such a layer by replacing ψ by −ψ + 2π/3 (Eq. 10), while keeping that of the other moiré superlattice unchanged. However, we can reduce the Hamiltonian of BN/G/NB to that of BN/G/BN (Eq. 15) by a similarity transformation that multiplies e 2iξψ ω −ξ ( e −2iξψ ω ξ ) to the Bloch bases |k (n) � with n ≡ 2, 3 (mod 4) for θ 1 , θ 3 , θ 4 ( θ 2 ). As a result, the resonant states are invariant with respect to the replacement of BN to NB. Figure 5a shows the valence band structures of BN/G/BN at θ 1 near the resonant states at C n plotted as a function of k 0 . The twelve Dirac cones are arranged on a circle with a radius �k = 2|G M,(l) i | sin(π/12) and they are strongly hybridized near k 0 = 0 to exhibit the characteristic dispersion including parabolic bottoms, a frilled band edge, and new Dirac points at −0.164 eV and −0.203 eV . We have similar resonant states also in the conduction band, while the energy spacing between the resonant states is much smaller than in the valence band, just like the cases of graphene on hBN 13 and quasicrystalline twisted bilayer graphene 22 . As predicted by the symmetry of H ξ ring , the band structures exhibit three-fold rotational symmetry around k 0 = 0 , as we can clearly see from the energy contours in Fig. 5b. The structures with the other angles, θ i ( i = 2, 3, 4 ), also show similar band dispersion, except that the resonant states are formed at the energies far from the charge neutrality point of graphene, since they have longer |K ξ − C n | . The energy splitting between the resonant states in BN/G/BN, O(|V 1 |) , is much smaller than that in the quasicrystalline twisted bilayer graphene, O(|u 0 |) , since the interaction here involves a second order scattering through the hBN layer. Nevertheless, the resonant states of BN/G/ BN appear at the energies much closer to the charge neutrality point of graphene than those of quasicrystalline twisted bilayer graphene (about ±1.7 eV) , since the wave vectors responsible for the interaction in BN/G/ BN ( O(|G M,(l) i |) ) are much shorter than those in quasicrystalline twisted bilayer graphene ( O(|a * i |) ). Thus, the resonant states of BN/G/BN appear at much smaller, experimentally feasible, electron densities.
At k 0 = 0 , we can reduce H ξ ring to an 8 × 8 form, where χ is e 2iξψ ω −ξ , and K stands for complex conjugation. Thus, the states with (m, m ′ ) = (0, −ξ) form twofold doublets, and belong to two-dimensional E irreducible representation of D 3 point group, while the states m = ξ is non-degenerate, and belong to either of A 1 or A 2 . Figure 6 shows the quasicrystalline wave functions of the third resonant state ( m = ξ ) in the hole side at k 0 = 0 on the graphene lattice. Figure 6a-c show the wave functions in a system in the absence of the inversion symmetry, i.e., V N = V B and ψ = π/3 (mod π) , such as graphene encapsulated by hBN layers, and Fig. 6d-f show those in a system with the inversion symmetry, i.e., V N = V B and ψ ≡ π/3 (mod π) , such as graphene encapsulated by a material having the same types of atoms in both sublattices. In both systems, the wave amplitude show the quasicrystalline order which is distributed on a limited number of sites in a pattern which is incompatible with the periodicity. Unlike the rotational symmetry of the double moiré pattern (12-fold), however, the wave functions in the absence of the inversion symmetry show a three-fold rotational symmetry, while those in the presence of the inversion symmetry show a six-fold rotational symmetry, which slightly resembles the probability distribution of the system with a true 12-fold rotational symmetry 22 . It should be noted that, however, the wave functions of the system with the inversion symmetry lose the six-fold rotational symmetry at k 0 = 0 . Figure 6a and d, b and e, c and f show the wave functions at θ 1 ,θ 2 ,θ 3 , respectively. The length scale of the patterns is much larger than that of the quasicrystalline wave functions in quasicrystalline twisted bilayer graphene 22  | is almost proportional to θ (Fig. 3c). Such an electron distribution would be prominent at the Finally, we discuss several factors that may affect the results in this work. (i) Effects of the theoretical model and parameters: In this work, we described the quasicrystalline resonant states by the effective model in the framework of a tight-binding Hamiltonian. Then, it is natural to ask how more accurate models, such as firstprinciple methods, affect the results. In our effective model, we describe the electronic states by using a few monolayer Bloch states which contribute significantly to the states. In first-principle methods, we get fully coupled wave functions of the states, and may analyze the contributions of each monolayer Bloch state by the projection of the wave functions to the Bloch states. The difference between the two models mainly comes from (1) the interaction between the resonant chain (Eq. 15) and the ignored monolayer Bloch states, and (2) the inaccurate parameters used in the tight-binding model (Eq. 7). Since (1) originate from the truncation of the bases, we can easily minimize the difference from (1) by using more bases in our model. However, since the energies of such truncated Bloch states are usually very different from the energy of the states which form the resonant chain, they make little change to the band dispersion. In quasicrystalline twisted bilayer graphene, for example, increasing the number of bases from the minimum 12-182 merely gives the replicas of the same bands at other wave vectors and make almost no change to the band dispersion ( Fig. 6 in 22 ). For (2), on the other hand, we agree that the parameters used in this work might not be accurate. The difference in the interlayer interaction strength u 0 will make a difference in the interaction strength within the resonant chain (Eq. 15), which affects the energy scale of the band dispersion, e.g., the energy spacing between the bands in Fig. 5a. The inaccuracy in V B and V N will make a difference in ψ (Eq. 10), which affects the slope and the degree of anisotropy of the band dispersion. Nevertheless, the wave functions still keep the D 3 symmetry, and the quasicrystalline states of any similar systems, which have the same lattice symmetry and the same type of the dominant atomic orbitals, are described by the universal form Eq. (15) with slightly different parameters.
(ii) Effects of using TMDC as encapsulating materials: We can get similar quasicrystalline resonant states in the system composed of a monolayer graphene encapsulated by two TMDC layers at suitable angle. Since the metal atom in a TMDC layer is far from the interface, the dominant interlayer interaction comes from the orbital hybridization between the p z orbital of the chalcogen atom in TMDC and the p z orbital of the carbon atom in www.nature.com/scientificreports/ graphene. The resonant interaction in this work couples many monolayer states at various points in the Brillouin zone. Thus, to make the magnitude of every interaction which forms the resonant chain identical, so that every Bloch states can contribute with an equal footing, t(q) (Eq. 11) should be isotropic, or has a rotational symmetry the same as C n . Since t(q) between the atomic orbitals with the same magnetic quantum number is isotropic 28 , the graphene encapsulated by TMDC layers will also exhibit the resonant states respecting the quasicrystalline rotational symmetry. The difference between the lattice symmetry of the hBN layer (honeycomb lattice) and the chalcogen layer (hexagonal or triangular lattice) makes a phase difference between the A-to-A and A-to-B interaction in Eq. (9).
(iii) Effects of lattice relaxation: Incommensurately stacked atomic layers tend to exhibit local lattice relaxation, especially when the moiré pattern is sufficiently large 30,34,35 . Thus, our systems, especially the structure with an angle of θ 1 , will also exhibit such a relaxation. In this work, however, we did not consider such effects, since the lattice relaxation mainly influences the band opening at the main Dirac point and Brillouin zone boundaries. In more detail, the lattice relaxation changes the parameters used (e.g., the u 0 for A-to-A and A-to-B in Eq. (9) would become different), so there will be some differences in the details of the band dispersion [see (i) above]. However, such a difference does not change the overall symmetry and the universal equation describing the energies and wave functions of the resonant states.
(iv) Tuning conductivity by gating: It is well-known that the conventional quasicrystals are poor conductors at the energies where the quasicrystalline states are dominant (i.e., at pseudogaps). And the van der Waals quasicrystals would be also likely poor conductors at the energies where the quasicrystalline states are dominant (i.e., where there is less decoupled states), as we can expect from the large band curvatures and the wave functions distributed selectively on a limited number of sites. Thus, if we tune the Fermi energies near the quasicrystalline states by gating, which is possible only when the states are within a feasible range from the charge neutrality point as in the current system, then we will be able to tune conductivity significantly by gating.