## Dirac fermions in quasicrystalline graphene

Quasicrystal lattices, which can have rotational order but lack translational symmetry, can be used to explore electronic properties of materials between crystals and disordered solids. Ahn *et al.* grew graphene bilayers rotated exactly 30° that have 12-fold rotational order. Electron diffraction and microscopy confirmed the formation of quasicrystals, and angle-resolved photoemission spectroscopy revealed anomalous interlayer electronic coupling that was quasi-periodic. The millimeter-scale layers can potentially be transferred to other substrates.

*Science*, this issue p. 782

## Abstract

Quantum states of quasiparticles in solids are dictated by symmetry. We have experimentally demonstrated quantum states of Dirac electrons in a two-dimensional quasicrystal without translational symmetry. A dodecagonal quasicrystalline order was realized by epitaxial growth of twisted bilayer graphene rotated exactly 30°. We grew the graphene quasicrystal up to a millimeter scale on a silicon carbide surface while maintaining the single rotation angle over an entire sample and successfully isolated the quasicrystal from a substrate, demonstrating its structural and chemical stability under ambient conditions. Multiple Dirac cones replicated with the 12-fold rotational symmetry were observed in angle-resolved photoemission spectra, which revealed anomalous strong interlayer coupling with quasi-periodicity. Our study provides a way to explore physical properties of relativistic fermions with controllable quasicrystalline orders.

Quasicrystals, which can have quasi-periodic orders (such as rotational symmetry) without spatial periodicity, have been used to study quantum states between the limits of periodic order and disorder (*1*–*15*). The study of the influence of quasi-periodic order has focused on extended wave functions of ordered states, pseudogaps, and fine structures of density of states (*16*–*28*). However, these studies have been limited to nonrelativistic fermions. We have found experimentally that a relativistic Dirac fermion quasicrystal can be realized when the Dirac electrons in a single-layer graphene are incommensurately modulated by another single-layer graphene that is rotated by exactly 30°. Such a twisted bilayer graphene (TBG) can form a 12-fold quasicrystalline order through interactions between the two layers.

This TBG quasicrystal is a direct material realization of the algorithmic quasicrystal construction process first proposed by Stampfli, with two ideal hexagonal grids rotated 30° with respect to each other (*29*, *30*). We observed that the graphene quasicrystal has distinct Dirac cone replicas in the Brillouin zone. Furthermore, it shows characteristic electronic structures that differ from the periodic TBG crystals with rotation angles near 30°. Angle-resolved photoemission spectroscopy (ARPES) revealed many replicas of Dirac cones with 12-fold rotational symmetry, indicating large enhancement of interlayer coupling with the quasicrystalline order. We show that the graphene quasicrystal can be used to explore quantum states of Dirac electrons in quasi-periodic order, just as a single-layer graphene has served as the quintessential material for studying Dirac fermions in periodic crystals. Furthermore, it has been predicted theoretically that four-dimensional (4D) quantum Hall effects of Dirac fermions can be studied in the 2D quasicrystal (*31*). Thus, the experimental achievement of the quasicrystal of Dirac fermions in two dimensions can be a starting point to approach higher-dimensional physics of quasiparticles in solids.

We grew graphene quasicrystals on the Si face of the 4H-SiC (0001) surface [see figs. S1 and S2 in supplementary materials (SM)]. The upper- and lower-layer graphene of the TBG had rotational angles of 0° and 30° (hereafter referred to as R0° and R30°), respectively, relative to the orientation of the SiC (0001) surface. Usually, when graphene is thermally grown on a Si face of 4H- or 6H-SiC (0001), the graphene assumes the R30° configuration (*32*–*35*). In contrast, by using borazine gas in a vacuum system, a hexagonal boron nitride (h-BN) layer with R0° can be grown epitaxially at 1050°C, where the typical π and σ energy bands of h-BN were observed in the previous ARPES experiments (fig. S1) (*36*). Thus, we first grew a h-BN layer with R0° by using borazine gas at 1050°C and then heated the system at a higher temperature of 1600°C to replace h-BN with a graphene layer while maintaining the same angle of 0° (fig. S1) (*36*). A buffer layer exists between graphene with R0° and bulk SiC, where the buffer layer has the same atomic structure as single-layer graphene but is insulating because some carbon atoms bond to underlying Si atoms (*33*). This procedure enabled the formation of the upper layer of the graphene quasicrystal (*36*). Finally, when the sample was heated further at 1600°C, a graphene layer with R30° grew between the first graphene layer with R0° and the SiC (0001) substrate, resulting in TBG with an exact rotational angle of 30° (fig. S2, B, D, and F). In this case, a buffer layer exists between the graphene quasicrystal and bulk SiC (fig. S2, B, D, and F). Another method can be used to grow the graphene quasicrystal. When hydrogen atoms are intercalated between the buffer layer underlying graphene with R0° and bulk SiC, the buffer layer changes into graphene with R30° (fig. S2, C, E, and G) (*33*). The hydrogen intercalation results in the same graphene quasicrystal on a hydrogen-terminated SiC surface (*33*). In this case, there is no buffer layer between the graphene quasicrystal and bulk SiC (fig. S2, C, E, and G). The nonexistence of boron and nitrogen atoms in the graphene quasicrystal was proved by x-ray photoemission spectroscopy (fig. S3). The Si 2p and C 1s spectra in fig. S3 suggest that the buffer layer in this system is the same as the buffer layer of epitaxial graphene with R30° on SiC (*32*, *33*). The present epitaxial growth has merit over a manual transfer (*37*–*39*) in making the quasicrystal because of the inevitable errors in twisting angles as well as the local distortions during the transfer process.

We determined the crystal orientation of the graphene layers from low-energy electron diffraction (LEED) and transmission electron microscopy (TEM) measurements (Fig. 1 and figs. S4 and S5). The upper and lower layers of the TBG can be clearly distinguished from the LEED pattern (Fig. 1A) because the intensity of the 1×1 LEED pattern of the upper-layer graphene (the red hexagon in Fig. 1A) is higher than that of the lower-layer graphene (the blue hexagon in Fig. 1A). The LEED spots other than the 1×1 LEED patterns come from the quasi-periodic order induced by the interlayer interaction. As shown in Fig. 1B, our observation is fully consistent with the simulated pattern from an atomic structure model of graphene quasicrystal (see also fig. S6). By using the 1-mm by 1-mm electron beam, we also confirmed that the LEED pattern is uniform throughout the entire sample size of up to 4 mm by 7 mm, where the sample size is limited only by our experimental setup and can be scaled up to a wafer scale (fig. S5). Furthermore, the coherent length of the graphene quasicrystal was similar to that of typical epitaxial graphene with R30° grown on SiC (fig. S4).

The graphene quasicrystal can be spatially mapped onto a quasicrystal lattice model constructed by dodecagonal compound tessellations (Fig. 1, D and E, and figs. S7, D and E; S8; and S9) (*29*, *30*). Squares, rhombuses, and equilateral triangles with different orientations can fill the entire space with a 12-fold rotationally symmetric pattern without translational symmetry. Because the Stampfli tiles have a fractal structure with self-similarity, the same pattern emerges at a larger scale with an irrational scaling factor. For the graphene quasicrystal, the Stampfli tiles have the scaling factor of (Fig. 1, D and E, and figs. S7 to S9) (*29*). As shown in the atomic model (Fig. 1, C and D), the graphene quasicrystal results in Stampfli tiles such as equilateral triangles and rhombuses. In the false-colored TEM image of the graphene quasicrystal transferred from a SiC wafer to a TEM grid (Fig. 1F and fig. S7D), the Stampfli tiles, including squares (blue), rhombuses (red), and equilateral triangles (green) with different orientations, are clearly observed. The LEED pattern and TEM image clearly indicate that the TBG with a rotational angle of 30° has a quasicrystalline order with 12-fold rotational symmetry. The detailed TEM images match the Stampfli tiling very well, with a few minute differences between them (figs. S8 to S12). The slight differences originate mainly from different internal structures of Stampfli tiles due to the different radial distances from the center (figs. S10 to S12). In addition, because our TEM image was taken for our transferred sample on a TEM grid, the height of the freestanding graphene could vary in space slightly, thus defocusing the image locally. The successful TEM experiments show that the graphene quasicrystal can be isolated from a substrate and is chemically and structurally stable at room temperature in air.

The characteristic electronic structures were determined from ARPES measurements. Figure 2A shows the constant energy map of ARPES spectra of the graphene quasicrystal, and Fig. 2C shows the enlarged view in the vicinity of the zero crystal momentum denoted by the Γ point (see also fig. S13). Several replicas of main Dirac cones respecting the 12-fold rotational symmetry were observed throughout the momentum space (Fig. 2), in sharp contrast to the usual TBGs with an arbitrary rotational angle other than 30° (*40*, *41*). Considering the large photon beam size of 3 mm by 3 mm in our ARPES system, the energy-momentum dispersions throughout the entire momentum space ensure the uniform sample quality, as in the LEED experiments. Like the LEED patterns, the intensities of the photoemission spectra of the lower graphene layer are weaker than those of the upper graphene layer (Fig. 2A) because of attenuated intensities for photoemitted electrons from the lower layer.

To understand the origin of the 12-fold Dirac cone replicas, we performed a numerical simulation of the ARPES spectra by using a single-particle tight-binding model that has been successfully used for describing various electronic properties of TBGs of general angles (*42*, *43*) (see SM for the detailed methods). Figure 3A shows the constant energy map at the Fermi energy (i.e., 0.3 eV above the Dirac point) of simulated ARPES spectra. The red (blue) contours indicate that those Dirac cones originate mainly from inter- and intralayer scatterings of the wave functions of the upper (lower) graphene layer. In stacked crystals, a state in the upper layer is scattered to every state in the lower layer, which satisfies the generalized Umklapp scattering condition of **k + G = k′ + G′**. Here, **k** and **G** (**k′** and **G′**) are the in-plane wave vector and reciprocal vector of the upper (lower) layer, respectively (*44*, *45*).

The coupling strength between and is proportional to *t*(**k** + **G**) = *t*(**k′** + **G′**), where *t*(**q**) is the in-plane Fourier transform of the interlayer transfer integral and where **q** is the equivalent wave vector common to both layers (*44*). Figure 3, B to D, shows some examples of the Umklapp scattering from the Dirac point (**k** = **K**) of the upper layer, which shows that the scattering with different vectors **G** is mapped to distinct wave vectors **k′** in the lower layer. We can also check this finding given that, though the states have the same wave function at the lattice point in the upper layer regardless of **G**, the states which satisfy **k′** + **G′** = **k** + **G** give distinct wave functions in the lower layer for different **G** (see SM for a detailed explanation). Each scattering event makes a replica of the Dirac cone (originally at **k** = **K**) at **k′**. Because **G** and **G′** are neither identical nor commensurate in graphene quasicrystal, there are many distinct Dirac cone replicas in the momentum space. Likewise, the Dirac cone at **k′** of the lower layer is replicated at **k** of the upper layer, satisfying **k + G = k′ + G′**.

From these calculations, the observed positions of Dirac cones in our ARPES experiment are reproduced in our simulations as shown in Fig. 3, E and F. The observed positions, however, cannot be reproduced if the angle deviates from 30° very slightly (fig. S13, D and E). For example, if the twisted angle is 29.9576°, the TBG has a translational symmetry with a periodic supercell that is 1351 times the graphene unit cell, and the computed positions of the Dirac cone replicas in this quasicrystal approximant were not consistent with the experimental ones (fig. S13E).

Close inspection of ARPES intensities of Dirac cone replicas highlights the important role of higher-order Umklapp scatterings in realizing 12-fold Dirac cones in graphene quasicrystals. Within a single Umklapp scattering process, the scattering strength is proportional to the wave vector **q** component of the in-plane Fourier transform of the interlayer transfer integral *t*(**q**). Figure 3, B to D, shows the examples of the scattering with the shortest three vectors , and Fig. 3E shows the entire map of the Dirac cone replicas up to the shortest eight vectors . The numbers in the circles indicate the order of the length of (see SM and fig. S19) (*44*, *45*). For a single Umklapp scattering, the scattering intensity should decrease rapidly with as shown in Fig. 3G because *t*(**q**) decays exponentially with large , as demonstrated in fig. S19C. However, because the number of Umklapp scattering events which give the same also increases very rapidly with increasing , the intensities of higher-order Dirac cone replicas increase markedly compared with single-scattering intensities, as shown in Fig. 3G (e.g., intensities are enhanced by factors of 1.15 × 10^{5} and 7.40 × 10^{11} with = 6.14 and 9.01 Å^{−1}, respectively). Thus, the resulting intensities for larger become similar (Fig. 3G).

In our experiment, we observed a similar trend in ARPES intensities (Fig. 3G). The first two intensities of Dirac cone replicas decrease exponentially, and the remaining Dirac cones show almost constant intensities. Although a general trend of intensities in our experiment agrees well with theoretical expectations, the absolute relative intensities do not match the computed ones. The discrepancies between experimental and simulated intensities of the Dirac cone replicas may originate from hidden scattering paths for large Umklapp processes, such as with the surface reconstruction of the SiC substrate (*32*, *34*, *46*–*48*), or impurities or localized states that may not be properly accounted for in the simulation on the basis of single-particle pictures.

Detailed shapes of Dirac cone replicas in the graphene quasicrystal also show their characteristic features. They were all electron-doped, and the Dirac point energies of replicas on the upper and lower layers were nearly identical (Fig. 4, A and B). The average values of the Dirac points of the upper- and lower-layer graphene of the graphene quasicrystal measured from different samples with different photon energies are approximately 0.29 and 0.30 eV, respectively, where the overall energy resolution was approximately 60 meV (fig. S15). In typical epitaxial graphene on a SiC (0001) substrate, electrons were transferred from the SiC substrate to graphene (*32*–*34*, *46*). Likewise, we expected a relative shift of the Dirac points between the two layers caused by the effective perpendicular electric field, as shown in recent studies on gated TBGs (*49*, *50*). Unlike gated TBGs (*49*, *50*), however, all Dirac cones in the sample had almost identical Dirac point energies, showing an anomalous doping behavior. Within a single-particle tight-binding approach, the doping dependence shown in Fig. 4 could not be described very well, implying either stronger interlayer coupling or larger interlayer screening. We also could not exclude other dopant sources such as carbon or silicon clusters. Further theoretical and experimental work is needed to resolve the anomalous doping in future.

This anomalous interlayer screening is in accordance with the exponential enhancement of ARPES intensities of high-order Dirac cone replicas (Fig. 3G) because of the higher-order Umklapp scatterings. Regardless of the large variation in ARPES intensities of Dirac cone replicas, their Fermi velocities are almost constant (Fig. 4C; see also fig. S13F) (*51*, *52*). The overall interaction strength between Dirac cones on upper and lower layers can be estimated as , where Δ*E* is the energy difference between and states. In TBGs, Δ*E* for the dominant interaction near the Dirac point is proportional to the rotation angle. Thus, TBGs with small angles exhibit very strong coupling between and because Δ*E* is very small, causing strong distortion of the Dirac cone structures, such as the additional renormalization of Fermi velocity (*42*, *53*). In the graphene quasicrystal, however, Δ*E* is very large, so we expect negligible distortion of the Dirac cone structures. As expected, all of the Dirac cones showed a well-defined linear momentum-energy relation (Fig. 4, B and C). Furthermore, the Fermi velocities of the Dirac cones did not differ as much as would be expected from Umklapp scattering with large values. The energy distribution curves (fig. S14) show that a typical pseudogap feature in epitaxial graphene (*46*, *54*) occurs at the Dirac points for both layers. Having established the quasicrystalline order in our TBG system and shown its consequence in realizing a distinct Dirac electron system, we expect that the TBG rotated by exactly 30° will be an important material for studying various unsolved issues (*26*–*28*, *30*, *31*) in anomalous physical properties of quasicrystals.

## Supplementary Materials

This is an article distributed under the terms of the Science Journals Default License.

## References and Notes

**Acknowledgments:**We thank S.-H. Kang, H. J. Choi, E.-G. Moon, L. A. Wray, A. Kent, M. Yang, and K. Ihm for fruitful discussions.

**Funding:**This work was supported by the National Research Foundation (NRF) of Korea (NRF-2015R1A2A2A01004853 and NRF-2017M2A2A6A01019384). M.K. was supported by JSPS KAKENHI grants JP25107005 and JP15K21722. Y.-W.S was supported by the NRF of Korea (grant 2017R1A5A1014862, SRC program vdWMRC center). P.M. was supported by the NYU Shanghai (start-up funds), the NYU-ECNU Institute of Physics, the NYU Global Seed Grant for Collaborative Research program, and the NSF of China (grant 11550110177, Research Fund for International Young Scientists program). P.K. was supported by the Global Research Laboratory program (2015K1A1A2033332) through the NRF of Korea. C.-W.Y. was supported by the National Research Council of Science and Technology (NST) (CRC-15-06-KIGAM) and the NRF of Korea (NRF-2011-0030058 and NRF-2015R1D1A1A01059653).

**Author contributions:**S.J.A., H.-W.K., H.-C.S., E.H.K., and J.R.A. contributed to the growth of the graphene quasicrystal and ARPES experiments. T.-H.K., H.W.C., and C.-W.Y. contributed to the TEM experiments. P.M., M.K., and Y.-W.S. contributed to the theoretical calculations. P.K. and S.-J.K. contributed to the interpretation of the experiments.

**Competing interests:**None declared.

**Data and materials availability:**All data needed to evaluate the conclusions in the paper are present in the paper or the supplementary materials.