Skip to main content

ORIGINAL RESEARCH article

Front. Carbon, 05 November 2024
Sec. Graphite-ene

Quasicrystalline 30° twisted bilayer graphene: fractal patterns and electronic localization properties

  • 1Instituto de Física, Universidade Federal do Rio de Janeiro, Rio de Janeiro, RJ, Brazil
  • 2Instituto de Física, Universidade Federal Fluminense, Niterói, RJ, Brazil

The recently synthesized 30° twisted bilayer graphene (30°-TBG) systems are unique quasicrystal systems possessing dodecagonal symmetry with graphene’s relativistic properties. We employ a real-space numerical atomistic framework that respects both the dodecagonal rotational symmetry and the massless Dirac nature of the electrons to describe the local density of states of the system. The approach we employ is very efficiency for systems with very large unit cells and does not rely on periodic boundary conditions. These features allow us to address a broad class of multilayer two-dimensional crystal with incommensurate configurations, particularly TBGs. Our results reveal that the 30°-TBG electronic spectrum consist of extended states together with a set of localized wave functions. The localized states exhibit fractal patterns consistent with the quasicrystal tiling.

1 Introduction

The electronic properties of bilayer systems when twisted at specific angles, have been a subject of intense research interest. In particular, twisted bilayer graphene (TBG) near to the so-called magic angles show a plethora of phenomena, such as superconductivity (Cao et al., 2018b; Yankowitz et al., 2019; Lu et al., 2019), magnetism (Sharpe et al., 2019; Lin et al., 2022; Jaoui et al., 2022; Vidarte and Lewenkopf, 2024), and strong electronic correlation effects (Cao et al., 2018a; Lu et al., 2019). Some of these interesting phenomena are also present in twisted multilayer graphene (Brzhezinskaya et al., 2021; Kononenko et al., 2022). Recently, the discovery of dodecagonal quasicrystals (Yao et al., 2018; Ahn et al., 2018) in 30°-twisted bilayer graphene (30°-TBG) has further attracted considerable attention to TBGs, as it offers a unique platform to study the effects arising from long-range order without translational symmetry. Unlike conventional quasicrystals, where all atoms are intrinsically located within a quasiperiodic order (Shechtman et al., 1984; Levine and Steinhardt, 1984), a 30°-TBG is viewed as an extrinsic quasicrystal because its quasiperiodicity arises from the interlayer coupling between two graphene monolayers (Moon et al., 2019).

Dodecagonal 30°-TBG quasicrystals have been successfully synthesized on various substrates, including SiC (Ahn et al., 2018; Bocquet et al., 2020; Fukaya et al., 2021), Pt (Yao et al., 2018), Cu (Yan et al., 2019; Pezzini et al., 2020; Deng et al., 2020; Liu et al., 2021), and a Cu-Ni alloy (Takesaki et al., 2016), using methods such as chemical vapor deposition and carbon segregation. The 30°-TBGs are experimentally identified by low-energy electron diffraction (LEED) (Ahn et al., 2018; Yao et al., 2018; Bocquet et al., 2020; Takesaki et al., 2016), high-energy positron diffraction (HEPD) (Fukaya et al., 2021), transmission electron microscopy (TEM) (Ahn et al., 2018; Deng et al., 2020), scanning tunneling microscopy (STM) (Yan et al., 2019), and Raman spectroscopy (Yao et al., 2018; Pezzini et al., 2020; Liu et al., 2021). Furthermore, magnetotransport measurements have confirmed the low-energy Dirac fermionic nature of these systems (Yan et al., 2019; Pezzini et al., 2020).

Various techniques have been employed to characterize dodecagonal quasicrystalline properties of 30°-TBGs. Angle-resolved photoemission spectroscopy (ARPES) measurements have revealed multiple Dirac cones with 12-fold rotational symmetry in 30°-TBGs (Yao et al., 2018; Ahn et al., 2018), a feature that is highly sensitive to deviations from the 30° twist angle (Ahn et al., 2018). The emergence of mirrored Dirac cones within the Brillouin zone of each graphene layer and a bandgap at the zone boundary, attributed to generalized Umklapp scattering between the two layers, characterizes their coupling (Yao et al., 2018; Ahn et al., 2018). STS measurements (Yan et al., 2019) reveal a suppression of density of state (DOS) consistent with the bandgap observed in previous ARPES studies, suggesting the same origin (Yan et al., 2019). Additionally, time- and angle-resolved photoemission spectroscopy measurements show an unbalanced electron distribution in the replica Dirac cone bands (Suzuki et al., 2019).

The dodecagonal 30°-TBG quasicrystal is a direct realization of the Stampfli construction process (Stampfli, 1986), which involves equilateral triangles, squares, and rhombuses as fundamental units arranged in a dodecagonal pattern. False-colored TEM images of 30°-TBGs (Ahn et al., 2018; Deng et al., 2020) clearly depict the Stampfli tiles in various orientations.

Atomic-resolved STM images of 30°-TBGs (Yan et al., 2019) have offered insights into the unique electron density distribution in real space. A striking feature is the flower-like contrast observed on the “coronene” patches of graphene, which exhibits higher symmetry than adjacent regions. These “flowers” align with the typical patterns of the dodecagonal quasicrystal tiling elements.

The lack of translation invariance in 30°-TBGs poses significant challenges for numerical calculations of its electronic properties. To address this, several strategies have been employed. One approach involves a series of periodic approximations with minimal lattice mismatch between the two layers to effectively reproduce the electronic and optical properties of 30°-TBG within a finite unit cell involving more than 107 atoms (Yu et al., 2019). Another proposed effective method is a k-space tight-binding model (Moon et al., 2019; Yu et al., 2020b). The model reveals wave functions with characteristic fractal inflations of the dodecagonal quasicrystal tiling. Moreover, an inverse participation ratio (IPR) analysis shows that system 12-fold symmetric resonant states are localized (Moon et al., 2019), a hallmark of quasicrystals (Niu and Nori, 1986; Kohmoto et al., 1987). These states are energetically far from the charge neutrality point, but it has been theoretically shown that by applying pressure or an electric field they can be pushed towards the Fermi energy (Yu et al., 2020a).

Recent theoretical predictions for 30°-TBG systems are rich and diverse. Numerical simulations suggest that the perfect superlubricity characterized by a scale-invariant sliding force related to the contact area, can be achieved for geometric sequences of dodecagonal quasicrystal tiling elements (Koren and Duerig, 2016). The system’s quasiperiodicity and weak interlayer coupling can originate quantum oscillations with spiral Fermi surfaces (Spurrier and Cooper, 2019). We note that the nature of interlayer hybridization selection rules governing interlayer coupling has been thoroughly investigated, via symmetry and group representation theory (Yu et al., 2022). Unlike the only allowed equivalent hybridization in D6h untwisted graphene bilayer, D6 twisted graphene bilayer permits equivalent and mixed hybridizations. Furthermore, D6d dodecagonal quasicrystal in 30°-TBG allows for both equivalent and nonequivalent hybridizations (Yu et al., 2022). Additionally, the presence of an exchange field and a Rashba spin-orbit coupling have been proposed to induce a nontrivial topological phase in 30°-TBG systems (Li and Wang, 2020).

In this study, we reexamine the electronic properties of 30°-TBGs, with a focus on the spatial distribution of the localized states, using the Haydock-Heine-Kelly (HHK) recursion technique (Haydock et al., 1972; Haydock et al., 1975; Haydock, 1980; Vidarte and Lewenkopf, 2022). This method relies on the nearsightedness of electronic matter (see, for instance, Ref (Prodan and Kohn, 2005)) and, as a real-space approach, unlike the standard techniques (Martin, 2004), eliminates the need for periodic boundary conditions, making it ideal for analyzing strongly disordered, amorphous, and quasicrystalline electronic systems. We demonstrate that the HHK method outperforms the previously proposed approaches for 30°-TBG systems. By developing a geometric model to interpret the LDOS results, we are able to uncover the spatial distribution of the LDOS suppression observed in Ref (Yao et al., 2018). and gain insights into the fractal nature of the localized electronic states.

2 Theory and methods

2.1 Lattice structure

We define the dodecagonal quasicrystal structure for a 30°-TBG system by initiating from AA-stacked bilayer graphene (Castro Neto et al., 2009) and subsequently rotating the top layer by 30° around the shared hexagonal center of both layers. The primitive lattice vectors of the bottom layer are given by a1b=3a0êx and a2b=3a0/2êx+3êy, while the top layer ones are defined by rotating the bottom layer vectors, namely, ait=R(π/6)aib. Here, the interlayer spacing between graphene layers is d0=3.35 Å and the carbon-carbon bond length is a0=1.42 Å (Castro Neto et al., 2009).

The obtained 30°-TBG system exhibits a dodecagonal quasicrystal structure with D6d point group symmetry, consistent with experimental observations (Yu et al., 2022). The system’s symmetry is primarily defined by six reflection planes σd and six two-fold rotation axes C2. Carbon atoms equidistant from the twist axis are equivalent and related through the application of these symmetry operators.

The 30°-TBG lattice can be spatially mapped onto a quasicrystal lattice model constructed by Stampfli tiles (Stampfli, 1986) discussed below. Here, we discuss the spatial distribution of lattice sites and propose a geometrical model based on two distinct Stampfli quasicrystal lattice construction methods.

Figure 1A illustrates the first Stampfli construction method obtained from the grid method (Stampfli, 1986; Ahn et al., 2018). The quasicrystal lattice is composed of equilateral triangles, squares and rhombuses (with 30° acute angle), arranged with different orientations to maintain a 12-fold rotational symmetry pattern, while lacking translational symmetry. The inflation rule of Stampfli tiles involved rotating the lattice by 15° around its perfect rotational symmetry center and scaling side lengths by a factor 2+3 (Stampfli, 1986). This process generates a new lattice whose vertices coincide with the ones of the original lattice.

Figure 1
www.frontiersin.org

Figure 1. Stampfli quasicrystal construction methods. (A) First construction method used a grid with 12-fold rotational symmetry consisting of equilateral triangles, squares and rhombuses. (B) Second construction method as a compound tessellation consisting of equilateral triangles and squares. The red and black circles indicate the dodecagonal centers. (C) Two dodecagons of identical size for both Stampfli construction methods.

Figure 1B presents the second construction method as a compound tessellation (Koren and Duerig, 2016; Stampfli, 1986). The dodecagonal compound tessellation is decomposed into equilateral triangles and squares, as depicted in Figure 1C. The center of each dodecagon corresponds to a vertex of a larger, similar dodecagonal element, with side lengths scaled by a factor 2+3 (Stampfli, 1986).

In the next section, we use the properties of the Stampfli construction methods to discuss electronic emergent localization features of 30°-TBG systems. More specifically, we employ the third inflation of Stampfli tiles, seed multiplied by a factor 2+33/2, as used in the first construction method. This ensures that the internal structures of each tile kind are similar (Ahn et al., 2018).

Let us introduce a geometric model that identifies the main locations of enhanced concentration of lattice sites. Figure 2 highlights five distinct radial distances denoted as r1, r2, r3, c1 and c2, where the Stampfli tile vertices serve as central points of emergent geometric site concentration. Central points formed by a single kind of Stampli tiles, coinciding with the main twist axis perpendicular to the plane, exhibit a perfect 12-fold symmetry. In contrast, central points formed by the vertices of different Stampfli tiles, corresponding to carbon atom sites at the same radial distances, are characterized by imperfect rings. Most imperfect rings exhibit an imbalance in emergent state measures, as detailed in Section 3. The imbalance primarily originates from the distinct symmetry groups of the Stampfli tiles: rhombuses, triangles and squares possess twofold, threefold and fourfold symmetry, respectively, see Figure 2A. Near-perfect rings are readily identified by the dodecagonal centers used in the second construction method, red circles in Figure 1B, where the center point regions exhibit pronounced 12-fold symmetry and balanced emergent state measures.

Figure 2
www.frontiersin.org

Figure 2. Geometric model: The five different highlighted-color, marked as r1, r2, r3, c1, c2, indicate an enhanced concentration of carbon atom sites at five different radial distances, using as references (A) the Stampfli tile vertices and (B) the dodecagonal centers obtained by the first and second construction methods, respectively. Red and blue lines correspond to the top and bottom graphene layers, respectively.

In summary, we use two complementary lattice construction methods to study electronic localization properties of 30°-TBGs. The first construction method effectively identifies all emergent localization measures, particularly those with exhibiting imbalances. The second construction method makes it easier to identify emergent localization measures with pronounced 12-fold symmetry.

2.2 Model Hamiltonian

The low-energy electronic structure of an incommensurate 30°-TBG system is effectively described by the single-orbital tight-binding Hamiltonian (Castro Neto et al., 2009; Reich et al., 2002).

H=i>jtijcicj+H. c,(1)

where i and j represent the lattice points and tij is the transfer integral between the Wannier pz electronic orbitals centered at the carbon atoms at the sites i and j. The key parameter in this model is the transfer integral tijt(Rij), with Rij=RiRj, where Ri and Rj are lattice vectors associated with the sites i and j. The hopping term tij is cast as (Trambly de Laissardière et al., 2010; Moon and Koshino, 2013; Uryu, 2004)

tR=VppπR1RezR2+VppσRRezR2,(2)

with

VppπR=Vppπ0expRa0δandVppσR=Vppσ0expRd0δ,(3)

where Vppπ0=2.7 eV is the transfer integral between the nearest-neighbor atoms that belong to the same graphene layer, Vppσ0=0.48 eV is the interlayer transfer integral between atoms located at different graphene sheets, and δ=0.45 Å is a characteristic decay length of the transfer integral (Uryu, 2004; Trambly de Laissardière et al., 2010). For R>5 Å the transfer integral is exponentially small and can be safely neglected. The intra- and interlayer transfer integrals have been obtained from a fit of the dispersions of graphene monolayer and graphene AB-stacked bilayer calculated by density functional theory methods (Trambly de Laissardière et al., 2010).

Variations in the interlayer coupling strength, Vppσ0, result in a monotonically shift of the energy peaks corresponding to the largest spectral weight of the atomic sites with the strongest interlayer interaction. These states are located at energies far from the Fermi level. On the other hand, variations in the intralayer hopping parameter, Vppπ0, primarily affect the Fermi velocity. The combined influence of both the intralayer and interlayer hopping parameters determines the energy position of the density of states suppression arising from the intersection of two Dirac cones, as discussed in Ref. (Yao et al., 2018).

The influence of an external magnetic field B is accounted for by the Peierls substitution (Peierls, 1933; Cresti, 2021). This approach modifies the hopping integrals between lattice sites to account for the magnetic field effect on electron motion, namely,

tijtijexpieRiRjdrAr,(4)

where is the Planck constant and e is the electron charge. We choose the vector potential as A(r)=0,Bx,0, that gives a constant magnetic field B=×A=Bêz perpendicular to TBG.

2.3 Numerical method

We compute the local spectral function of 30°-TBG quasicrystals, described by Equations 14, using the Haydock-Heine-Kelly (HHK) recursion technique (Haydock et al., 1972; Haydock et al., 1975; Haydock, 1980; Vidarte and Lewenkopf, 2022). The HHK method is an O(N) approach which consists of a recursive procedure (Haydock, 1980) that transforms an arbitrary sparse Hamiltonian matrix in a tridiagonal form. Next, one evaluates the diagonal Green’s function, Giir(ϵ)Gr(Ri,Ri,ϵ), by a continued fraction expansion, that is very amenable for numerical calculation. A detailed presentation of the method and its implementation can be found, for instance, in Refs. (Haydock, 1980; Vidarte and Lewenkopf, 2022).

The local density of states (LDOS) at any site j can be written as

ρjϵ=1πImGjjrϵ1πlimη0+ImGjjϵ+iη.(5)

In practice, a finite η serves as a convenient regularization parameter (Vidarte and Lewenkopf, 2022; Haydock et al., 1975). The latter is related to physical processes when addressing realistic systems, since it accounts for the imaginary part of the Green’s function self-energy due, for instance, to disorder effects that are ubiquitous in graphene systems (Mucciolo and Lewenkopf, 2010).

As noted in the introduction, the HHK method is an order N real-space approach that does not require periodic boundary conditions periodic boundary conditions or finite system sizes. Here, we increase the number of considered lattice sites until the LDOS converges, typically with approximately N107 atomic sites. Although, in principle, the calculation must be repeated for each atomic site, the geometric model introduced above, combined with the self-similarity of the electronic spectral functions for different lattice inflations (see below), indicates that only a finite number of sites need to be considered. This feature is essential for the computational efficiency of the method.

3 Results

Using the HHK method we compute the LDOS given by Equation 5 and confirm that the equivalent sites at identical radial distances share the same local spectral function due to the most relevant symmetry operators, as discussed in Section 2.1. Consequently, our calculations require the computation of only one twenty-forth of the total number of carbon atoms, within a specified circular region of the 30°-TBG dodecagonal quasicrystal.

Figure 3 presents a numerical simulation of the electron density distribution at the charge neutrality point, ϵCNPdϵρi(ϵ), in a real-space image of the 30°-TBG system. The real-space image reveals a homogeneous electron density distribution within each graphene layer. When the two layers are considered together, a non-periodic dodecagonal pattern emerges. This finding is in line with experimentally obtained atomic resolution STM images (Deng et al., 2020; Yan et al., 2019).

Figure 3
www.frontiersin.org

Figure 3. Electron density distribution in a real-space image of the 30°-TBG system. Black lines represent the third inflation of Stampfli tiles using the first construction method, while the white lines correspond to the second method. The areas of the red and blue disks are proportional to the local charge density in the top and bottom layers, respectively. The purple regions correspond to overlapping blue and red disks. The green and sky-blue circles indicate the flower-like structure. The figures on the right side correspond to a zoom-in of the restructure on the left side.

The inset in Figure 3 reveals an intriguing feature of the space-resolved electron density distribution of a 30°-TBG: the flower-like structures. These structures indicate local-symmetry regions associated with 12-atom rings. Collectively, the flower-like structures exhibit a characteristic dodecagonal quasicrystal pattern with strict 12-fold rotational symmetry. Within each flower-like structure, the local symmetry regions also possess 12-fold rotational symmetry. In the Stampfli rhombuses, specifically on the smaller diagonal, two flower-like structures are linked by two carbon atoms located in the second ring closest to the center points.

Contrary to what has been described in the experimental literature (Yan et al., 2019), all flower-like structures are easily identified by the first Stampfli construction model, and not by the second one.

3.1 Local electronic structures

Figure 4 shows the evaluation of “integrated spectral functions” for various Stampfli inflation orders, considering unit structures such as equilateral triangles, squares, and rhombuses. The integrated spectral function of each Stampfli unit structure is computed as the sum of the set of local spectral functions, Fi(ϵ)=jiρj(ϵ), where j labels the carbon atoms sites within the Stampfli tile and i indicates radial distance order. The averaged deviations between the spectral functions for different inflation orders and the reference spectral function, F1(ϵ), are computed as Fi(ϵ)F1(ϵ)2dϵ. This analysis provides insights into the impact of inflation order on the electronic properties of the Stampfli tiles.

Figure 4
www.frontiersin.org

Figure 4. Electronic structure evaluation with different orders of inflation for unit structures such as (A) triangles, (B) squares, and (C) rhombuses. Means correspond to the sum of discrepancies between two spectral functions over the energy band for different radial distances. The error bars quantify the standard deviations from the mean value that are mainly due to the singularities of the valence band.

In Figure 4 error bars quantify the standard deviation from the mean value, primarily due to the peaks associated to the valence band singularities that vary along the radial distance. For a given kind of Stampfli tile, the first inflation shows fluctuations of order 103 and 104 along the radial distance. By the third inflation, standard deviations are negligible for all spectral functions, regardless of radial distances. A small standard deviation indicates greater similarity with the first spectral function. Thus, the third inflation is the smallest order at which the spectral functions of the Stampfli unit structures exhibit typical patterns characteristic of the dodecagonal quasicrystal. By accounting for the normalization factor that relates the spectral functions of the different Stampfli unit structures, we observe an excellent agreement within the numerical precision. Our findings suggest that the electronic structure of the system can be calculated using only a single Stampfli unit structure.

In summary, the crystalline structure and local spectral functions of 30°-TBGs match the Stampfli tiles, with minor discrepancies observed for the first and second inflations. These slight differences primarily stem from variations in the internal structures of Stampfli tiles at different radial distances (Ahn et al., 2018). As the sizes of the Stampfli tiles increase, these differences become increasingly insignificant, revealing the fractal structure of the LDOS.

3.2 Emergent localization

Let us now examine the localized states corresponding to the singularities in the valence band (Moon et al., 2019). Figure 5 presents the 30°-TBG system density of states, DOS(ϵ)=iρi(ϵ), and the inverse participation ratio, IPR(ϵ)=i[ρi(ϵ)]2, based on the quasicrystal lattice model constructed by Stampfli tiles.

Figure 5
www.frontiersin.org

Figure 5. (A) DOS of the dodecagonal quasicrystalline in 30°-TBG. (B) Peaks marked as α, β, γ, τ, μ correspond to the DOS singularities in the valence band. (C) IPR for the same energy window.

Figures 6A–E depict real-space images of emergent localized state of a 30°-TBG system, corresponding to the peaks α, β, γ, τ, and μ of DOS singularities in the valence band, respectively. The sizes of the red (blue) circles are proportional to the LDOS at the top (bottom) layer site at which the circle is centered. The peak energies α, β, γ, τ and μ are approximately 1.163 eV, 0.949 eV, 0.702 eV, 0.842 eV and 0.44 eV, respectively. All emergent states exhibit a non-periodic dodecagonal pattern and can be spatially represented by a quasicrystal lattice model at the third inflation of the Spampfli tiles.

Figure 6
www.frontiersin.org

Figure 6. Real-space images of emergent states localization of a 30°-TBG system. (A–E) correspond to the peaks α, β, γ, τ, μ marked in Figure 5B, respectively. Black lines represent the third inflation of Stampfli tiles used in the first construction method. The areas of the red and blue disks are proportional to the LDOS at a particular energy in the top and bottom layers, respectively. The black circles of (D), (E), and the insets of (A–C) are the zoom-in structures in relation to the central point that coincides with the twist axis, with perfect 12-fold symmetry. The red and green circles of the insets of (A–C) are the zoom-in structure in relation to the central point formed by the vertices of different Stampfli tiles. Bar graphs of (F–E) correspond to the DOS, ρ(ϵi, d) at the energy peaks ϵi with i={α,β,γ,τ,μ}, respectively. See main text for the definitions of ρ(ϵi, d) and ρmax(ϵi, d). The histograms relate the LDOS at different radial distances d={r1,r2,r3,c1,c2}, indicated in Figure 2.

Figures 6F–J present histograms illustrating the spatial distribution of states corresponding to the energy peaks ϵi with i={α,β,γ,τ,μ}, respectively. We analyze the variation of emergent localized states along the spacial radial distances, d={r1,r2,r3,c1,c2}, as indicated in Figure 2, where the vertices of the Stampfli tiles serve as central points. Each histogram shows the average LDOS for the energy peak ϵi, ρ(ϵi,d)=jdρj(ϵi), where j corresponds to the atomic sites located at a distance d. The average LDOS is normalized by its maximum value, ρmax(ϵi,d).

The black circles in the insets of Figures 6A–C and in Figures 6D, E correspond to a zoom-in view of the structures of the central point regions formed by exclusively by Stampfli rhombuses associated with the energy peaks ϵi with i={α,β,γ,τ,μ}, respectively. The distribution of the LDOS in each peak, ρj(ϵi), ehxibits a perfect 12-fold symmetry. In the insets of Figures 6A–C, the red and green circles identify regions corresponding to dodecagonal centers used in the second construction method and the vertices of different Stampfli tiles used in the first construction method, respectively. While the red circles do not locally preserve 12-fold symmetry, their distribution is close to the latter. However, the green circles deviate significantly form 12-fold symmetry, displaying an imbalance due to the different symmetry groups of the Stampfli tiles.

We observe that the average LDOS at specific radial distances, ρ(ϵi,d), in regions with strong imbalances resemble those in the regions with perfect 12-fold symmetry.

Figures 6A, F demonstrate that the α-peak states are highly localized on the second ring of carbon atoms closest to the center points, at a radial distance r2, with ρmax(α,r2)=0.203 a.u. Interestingly, these second-ring atomic sites can be considered as opposite of dimer sites, due to their significant spacing.

Figures 6B, G reveal a quasi-homogeneous emergent states distribution for the β-peak, with ρmax(β,c2)=0.119 a.u., indicating less less localization at radial distance c1. Figures 6I, J show that the peaks τ and μ peaks also exhibit quasi-homogeneous emergent states distributions, with ρmax(τ,r1)=0.089 a.u. and ρmax(μ,c1)=0.062 a.u., respectively. While the τ-peak states are less localized at c1, the μ-peak states are less localized at r2 and c2.

Figure 6H illustrates that the γ-peak states are highly localized at a radial distance c1, coinciding with the dimer sites, with ρmax(γ,c1)=0.102 a.u. In Figure 6C, the largest localizations are indicated by the purple disks, which correspond to overlapping red and blue disks representing the top and bottom layers, respectively.

3.3 Suppression of the DOS

We will now focus on the electronic density of states (DOS) suppression reported in Refs. (Yao et al., 2018; Yan et al., 2019). Our goal is to use the geometric model introduced in Section 2.1 to understand the spatial distribution of the electronic states corresponding to the suppressed DOS.

Figures 7A–F show the emergent suppressions in the LDOS of the 30°-TBG systems. These suppressions are quantified by Sp1ρ(ϵi,d)/ρmax(ϵi,d) at the energy ϵi=0.536 eV, that lies between the γ and μ-peak energies. A larger Sp value indicates a greater suppression of the LDOS, suggesting stronger interlayer coupling between the graphene sheets. The minimum and maximum LDOS values correspond to ρmin(ϵi,r2)=0.019 a.u. and ρmax(ϵi,c1)=0.05 a.u., respectively.

Figure 7
www.frontiersin.org

Figure 7. Emergent suppression strength in the LDOS of a 30°-TBG system. The areas of the red and blue disks are proportional to the suppression strengths in the top and bottom layers, respectively. (A–C) Zoom-in of the structures close to the twist axis, in the dodecagonal centers and in the vertices of different Stampfli tiles, respectively. (D–F) represent the internal structures of the Stampfli tiles used in the first construction method. The histogram in (G) relates the LDOS at different radial distances, indicated in Figure 2.

In Figure 7A, the black circle corresponds to a zoomed-in view of the structure close to the twist axis. The real-space distribution of suppression strengths exhibits a perfect 12-fold symmetry. The red circle, Figure 7B, which corresponds to the dodecagonal centers, shows a LDOS distribution that only slightly deviates from 12-fold symmetry. In contrast, the green circle, corresponding to the vertices of different Stampfli tiles, exhibits a much stronger imbalance in the second ring closest to the center, violating the 12-fold symmetry. Figures 7D–F illustrate the distributions of suppression strengths within the internal structures of Stampfli tiles, exhibiting double, triple and quadruple symmetry. Figure 7G reveals that the emergent suppression strengths are highly localized at radial distances r2 and c2.

3.4 Landau levels spectrum

Let us now perform a numerical investigation of the LL spectra in 30°-TBG systems and compare these results with experimental observations (Yan et al., 2019). For that purpose, we compute the DOS, ρ(ϵ), by integrating the LDOS, ρj(ϵ), over the atomic sites within a Stampfli unitary structure.

Figure 8A shows a series of DOS for a pristine 30°-TBG system subjected to different perpendicular magnetic fields B. The LLs are broadened using η=5 meV (Vidarte and Lewenkopf, 2022). The curves are offset in the y-axis for clarity. The spectra consist of a sequence of well-defined LL peaks labeled by the index N with increasing field. Both the LL peak spacing and LL peak heights (related to the LL degeneracy) grow with increasing magnetic fields (Goerbig, 2011).

Figure 8
www.frontiersin.org

Figure 8. (A) DOS (in arbitrary units) as a function energy E (in meV) of a pristine 30°-TBG system for magnetic fields of B=6, 8, 9 and 10 T. (B) LL peak energies (in meV) as a function of |N|B, where N is the LL index. The symbols correspond to the peaks in (A), and the solid line shows the linear fit with the experimental data of Ref. (Yan et al., 2019). (C) Comparison between the LL spectrum versus Sgn(N)(|N|B)1/2 for 30°-TBG and monolayer systems for B=25 T. The inset shows the deviation ΔϵN=ϵN30°TBGϵNG (in meV) as a function the LL index N.

The LL energies, ϵN, are determined by using a smaller regularization parameter, η=0.1 meV, and fitting the computed DOS of each LL peak with a Lorentzian distribution. Figure 8B shows that our results agree with the relation ϵNϵF=Sgn(N)vF2e|N|B, where vF the Fermi velocity, that is unique to massless Dirac fermions. For clarity, we subtract the Fermi energy ϵF, measured at the energy of the charge neutrality point. Fitting both electron and hole branches of the experimental data yields a Fermi velocity of vF=(1.100±0.001)×106 m/s, as reported in Ref. (Yan et al., 2019). In summary, our numerical results show a good agreement with the STS measurements, confirming massless Dirac fermions in 30°-TBG systems.

Let us compare the 30°-TBG LL energies with those of graphene monolayers. Figure 8C shows ϵN for 30°-TBG and graphene monolayer systems as a function of (|N|B)1/2 for 25 T. We verify that for low energies, corresponding to |N|20, the LL energies show, to very good accuracy, a linear dependence with (|N|B)1/2. Thus, we conclude that the electronic properties of the decagonal quasicrystal in 30°-TBG behave like a system with two separated monolayer graphene in the low-energy region.

Finally, let us examine the large-|N| limit, where one expects the LLs energies to deviate from the linear behavior (Goerbig, 2011; Vidarte and Lewenkopf, 2022). In the case of 30°-TBGs, these deviations arise due to long-range hopping terms that fit the nonlinear features of the dispersion relation (Vidarte and Lewenkopf, 2022; Plochocka et al., 2008). The inset of Figure 8C gives the LL energy difference between 30°-TBG and graphene monolayer systems, ΔϵN=ϵN30° TBGϵNG. We observe a notable discrepancy in the valence band of 20 meV for N=15. Nonetheless, the discrepancy is relatively small, reflecting the similarity of the DOS in the 30°-TBG and graphene systems at B=0. At higher energies, near ϵ1 eV where the electronic states are localized (at B=0), we anticipate larger deviations. However, given that the LLs spacing at such energies is 10 meV at B=25 T, experimental detection of these effects is challenging due to the broadening caused by the intrinsic disorder and the significant excitation and/or doping required.

4 Conclusion and discussion

We have studied the local spectral properties and the emergence of real-space electronic localization in dodecagonal quasicrystal structures of 30°-TBG systems. The HHK method’s efficiency in calculating local spectral functions without relying on periodic boundary conditions outperforms atomistic calculations (Moon et al., 2019) by orders of magnitude. Furthermore, our findings exhibit excellent agreement with the localized states reported in Ref (Moon et al., 2019). and the experimental density of states observed in Refs. (Yao et al., 2018; Yan et al., 2019).

Our results indicate that the electronic structure of 30°-TBG system can be accurately calculated using only one of Stampfli tiles at the third inflation level, despite the observation that the spectral functions of the three Stampfli tiles are similar but not identical. The spectral functions across the occupied electronic band reveal patterns characteristic of the dodecagonal quasicrystal, influenced by the varying number of carbon atom sites in the internal crystalline structure of the Stampfli tiles.

We developed a simple geometric model, applicable to the two lattice construction methods, to analyze the LDOS features in 30°-TBG systems. The emergent states associated with the singularities in the valence band are all identified with the Stampfli tiles, primarily those with imbalances in the non-perfect rings centered on the polygon vertices. Perfect and near-perfect rings are identified with decagonal centers that exhibit emergent localization measures with pronounced 12-fold symmetry and optimal balance.

We find that the most pronounced effects of interlayer coupling between two graphene layers of the 30°-TBG system are localized in the second ring centered on the polygon vertices (r2) and in the centers of the Stampfly squares (c2). This observation suggests a significant suppression of the LDOS at the intersection between the original and mirrored Dirac cones.

As an outlook, the efficiency of our method enables the analysis of twisted multilayer graphene systems (Brzhezinskaya et al., 2021; Kononenko et al., 2022), both commensurate and incommensurate, which involve a significantly larger number of atomic sites. Of particular interest are the recently observed moiré quasicrystals in twisted trilayer graphene (Hao et al., 2024), which present a promising direction for further development of our work.

Data availability statement

The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.

Author contributions

KV: Conceptualization, Writing–original draft, Writing–review and editing, Formal Analysis, Investigation, Methodology. CL: Conceptualization, Writing–original draft, Writing–review and editing, Supervision.

Funding

The author(s) declare that financial support was received for the research, authorship, and/or publication of this article. This work was partially supported by the Brazilian Institute of Science and Technology (INCT) in Carbon Nanomaterials and the Brazilian agencies CAPES, CNPq, FAPEMIG, and FAPERJ.

Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

References

Ahn, S. J., Moon, P., Kim, T.-H., Kim, H.-W., Shin, H.-C., Kim, E. H., et al. (2018). Dirac electrons in a dodecagonal graphene quasicrystal. Science 361, 782–786. doi:10.1126/science.aar8412

PubMed Abstract | CrossRef Full Text | Google Scholar

Bocquet, F. C., Lin, Y.-R., Franke, M., Samiseresht, N., Parhizkar, S., Soubatch, S., et al. (2020). Surfactant-mediated epitaxial growth of single-layer graphene in an unconventional orientation on SiC. Phys. Rev. Lett. 125, 106102. doi:10.1103/PhysRevLett.125.106102

PubMed Abstract | CrossRef Full Text | Google Scholar

Brzhezinskaya, M., Kononenko, O., Matveev, V., Zotov, A., Khodos, I. I., Levashov, V., et al. (2021). Engineering of numerous moiré superlattices in twisted multilayer graphene for twistronics and straintronics applications. ACS Nano 15, 12358–12366. doi:10.1021/acsnano.1c04286

PubMed Abstract | CrossRef Full Text | Google Scholar

Cao, Y., Fatemi, V., Demir, A., Fang, S., Tomarken, S. L., Luo, J. Y., et al. (2018a). Correlated insulator behaviour at half-filling in magic-angle graphene superlattices. Nature 556, 80–84. doi:10.1038/nature26154

PubMed Abstract | CrossRef Full Text | Google Scholar

Cao, Y., Fatemi, V., Fang, S., Watanabe, K., Taniguchi, T., Kaxiras, E., et al. (2018b). Unconventional superconductivity in magic-angle graphene superlattices. Nature 556, 43–50. doi:10.1038/nature26160

PubMed Abstract | CrossRef Full Text | Google Scholar

Castro Neto, A. H., Guinea, F., Peres, N. M. R., Novoselov, K. S., and Geim, A. K. (2009). The electronic properties of graphene. Rev. Mod. Phys. 81, 109–162. doi:10.1103/RevModPhys.81.109

CrossRef Full Text | Google Scholar

Cresti, A. (2021). Convenient Peierls phase choice for periodic atomistic systems under magnetic field. Phys. Rev. B 103, 045402. doi:10.1103/PhysRevB.103.045402

CrossRef Full Text | Google Scholar

Deng, B., Wang, B., Li, N., Li, R., Wang, Y., Tang, J., et al. (2020). Interlayer decoupling in 30° twisted bilayer graphene quasicrystal. ACS Nano 14, 1656–1664. doi:10.1021/acsnano.9b07091

PubMed Abstract | CrossRef Full Text | Google Scholar

Fukaya, Y., Zhao, Y., Kim, H.-W., Ahn, J. R., Fukidome, H., and Matsuda, I. (2021). Atomic arrangements of quasicrystal bilayer graphene: interlayer distance expansion. Phys. Rev. B 104, L180202. doi:10.1103/PhysRevB.104.L180202

CrossRef Full Text | Google Scholar

Goerbig, M. O. (2011). Electronic properties of graphene in a strong magnetic field. Rev. Mod. Phys. 83, 1193–1243. doi:10.1103/RevModPhys.83.1193

CrossRef Full Text | Google Scholar

Hao, C.-Y., Zhan, Z., Pantaleón, P. A., He, J.-Q., Zhao, Y.-X., Watanabe, K., et al. (2024). Robust flat bands in twisted trilayer graphene moiré quasicrystals. Nat. Commun. 15, 8437. doi:10.1038/s41467-024-52784-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Haydock, R. (1980). “The recursive solution of the Schrödinger equation,”. Solid State Phys., 35, 215–294. doi:10.1016/S0081-1947(08)60505-6

CrossRef Full Text | Google Scholar

Haydock, R., Heine, V., and Kelly, M. J. (1972). Electronic structure based on the local atomic environment for tight-binding bands. J. Phys. C. Solid State Phys. 5, 2845–2858. doi:10.1088/0022-3719/5/20/004

CrossRef Full Text | Google Scholar

Haydock, R., Heine, V., and Kelly, M. J. (1975). Electronic structure based on the local atomic environment for tight-binding bands: II. J. Phys. C. Solid State Phys. 8, 2591–2605. doi:10.1088/0022-3719/8/16/011

CrossRef Full Text | Google Scholar

Jaoui, A., Das, I., Di Battista, G., Díez-Mérida, J., Lu, X., Watanabe, K., et al. (2022). Quantum critical behaviour in magic-angle twisted bilayer graphene. Nat. Phys. 18, 633–638. doi:10.1038/s41567-022-01556-5

CrossRef Full Text | Google Scholar

Kohmoto, M., Sutherland, B., and Tang, C. (1987). Critical wave functions and a cantor-set spectrum of a one-dimensional quasicrystal model. Phys. Rev. B 35, 1020–1033. doi:10.1103/PhysRevB.35.1020

PubMed Abstract | CrossRef Full Text | Google Scholar

Kononenko, O., Brzhezinskaya, M., Zotov, A., Korepanov, V., Levashov, V., Matveev, V., et al. (2022). Influence of numerous moiré superlattices on transport properties of twisted multilayer graphene. Carbon 194, 52–61. doi:10.1016/j.carbon.2022.03.033

CrossRef Full Text | Google Scholar

Koren, E., and Duerig, U. (2016). Superlubricity in quasicrystalline twisted bilayer graphene. Phys. Rev. B 93, 201404. doi:10.1103/PhysRevB.93.201404

CrossRef Full Text | Google Scholar

Levine, D., and Steinhardt, P. J. (1984). Quasicrystals: a new class of ordered structures. Phys. Rev. Lett. 53, 2477–2480. doi:10.1103/PhysRevLett.53.2477

CrossRef Full Text | Google Scholar

Li, Z., and Wang, Z. F. (2020). Quantum anomalous Hall effect in twisted bilayer graphene quasicrystal. Chin. Phys. B 29, 107101. doi:10.1088/1674-1056/abab77

CrossRef Full Text | Google Scholar

Lin, J.-X., Zhang, Y.-H., Morissette, E., Wang, Z., Liu, S., Rhodes, D., et al. (2022). Spin-orbit–driven ferromagnetism at half moiré filling in magic-angle twisted bilayer graphene. Science 375, 437–441. doi:10.1126/science.abh2889

PubMed Abstract | CrossRef Full Text | Google Scholar

Liu, J., Wang, Z., Ling, D., Wei, D., Lv, W., Kang, X., et al. (2021). Synchronous growth of 30°-twisted bilayer graphene domains with millimeter scale. 2D Mater. 8, 021002. doi:10.1088/2053-1583/abda0e

CrossRef Full Text | Google Scholar

Lu, X., Stepanov, P., Yang, W., Xie, M., Aamir, M. A., Das, I., et al. (2019). Superconductors, orbital magnets and correlated states in magic-angle bilayer graphene. Nature 574, 653–657. doi:10.1038/s41586-019-1695-0

PubMed Abstract | CrossRef Full Text | Google Scholar

Martin, R. (2004). Electronic structure: basic theory and practical methods. Cambridge University Press.

Google Scholar

Moon, P., and Koshino, M. (2013). Optical absorption in twisted bilayer graphene. Phys. Rev. B 87, 205404. doi:10.1103/PhysRevB.87.205404

CrossRef Full Text | Google Scholar

Moon, P., Koshino, M., and Son, Y.-W. (2019). Quasicrystalline electronic states in rotated twisted bilayer graphene. Phys. Rev. B 99, 165430. doi:10.1103/PhysRevB.99.165430

CrossRef Full Text | Google Scholar

Mucciolo, E. R., and Lewenkopf, C. H. (2010). Disorder and electronic transport in graphene. J. Phys. Condens. Matter 22, 273201. doi:10.1088/0953-8984/22/27/273201

PubMed Abstract | CrossRef Full Text | Google Scholar

Niu, Q., and Nori, F. (1986). Renormalization-group study of one-dimensional quasiperiodic systems. Phys. Rev. Lett. 57, 2057–2060. doi:10.1103/PhysRevLett.57.2057

PubMed Abstract | CrossRef Full Text | Google Scholar

Peierls, R. (1933). Zur Theorie des Diamagnetismus von Leitungselektronen. Z. Phys. 80, 763–791. doi:10.1007/BF01342591

CrossRef Full Text | Google Scholar

Pezzini, S., Mišeikis, V., Piccinini, G., Forti, S., Pace, S., Engelke, R., et al. (2020). 30°-Twisted bilayer graphene quasicrystals from chemical vapor deposition. Nano Lett. 20, 3313–3319. doi:10.1021/acs.nanolett.0c00172

PubMed Abstract | CrossRef Full Text | Google Scholar

Plochocka, P., Faugeras, C., Orlita, M., Sadowski, M. L., Martinez, G., Potemski, M., et al. (2008). High-energy limit of massless Dirac fermions in multilayer graphene using magneto-optical transmission spectroscopy. Phys. Rev. Lett. 100, 087401. doi:10.1103/PhysRevLett.100.087401

PubMed Abstract | CrossRef Full Text | Google Scholar

Prodan, E., and Kohn, W. (2005). Nearsightedness of electronic matter. Proc. Natl. Acad. Sci. U.S.A. 102, 11635–11638. doi:10.1073/pnas.0505436102

PubMed Abstract | CrossRef Full Text | Google Scholar

Reich, S., Maultzsch, J., Thomsen, C., and Ordejón, P. (2002). Tight-binding description of graphene. Phys. Rev. B 66, 035412. doi:10.1103/PhysRevB.66.035412

CrossRef Full Text | Google Scholar

Sharpe, A. L., Fox, E. J., Barnard, A. W., Finney, J., Watanabe, K., Taniguchi, T., et al. (2019). Emergent ferromagnetism near three-quarters filling in twisted bilayer graphene. Science 365, 605–608. doi:10.1126/science.aaw3780

PubMed Abstract | CrossRef Full Text | Google Scholar

Shechtman, D., Blech, I., Gratias, D., and Cahn, J. W. (1984). Metallic phase with long-range orientational order and no translational symmetry. Phys. Rev. Lett. 53, 1951–1953. doi:10.1103/PhysRevLett.53.1951

CrossRef Full Text | Google Scholar

Spurrier, S., and Cooper, N. R. (2019). Theory of quantum oscillations in quasicrystals: quantizing spiral Fermi surfaces. Phys. Rev. B 100, 081405. doi:10.1103/PhysRevB.100.081405

CrossRef Full Text | Google Scholar

Stampfli, P. (1986). A dodecagonal quasiperiodic lattice in two dimensions. Helv. Phys. Acta 59, 1260–1263.

Google Scholar

Suzuki, T., Iimori, T., Ahn, S. J., Zhao, Y., Watanabe, M., Xu, J., et al. (2019). Ultrafast unbalanced electron distributions in quasicrystalline 30° twisted bilayer graphene. ACS Nano 13, 11981–11987. doi:10.1021/acsnano.9b06091

PubMed Abstract | CrossRef Full Text | Google Scholar

Takesaki, Y., Kawahara, K., Hibino, H., Okada, S., Tsuji, M., and Ago, H. (2016). Highly uniform bilayer graphene on epitaxial Cu–Ni(111) alloy. Chem. Mater. 28, 4583–4592. doi:10.1021/acs.chemmater.6b01137

CrossRef Full Text | Google Scholar

Trambly de Laissardière, G., Mayou, D., and Magaud, L. (2010). Localization of Dirac electrons in rotated graphene bilayers. Nano Lett. 10, 804–808. doi:10.1021/nl902948m

PubMed Abstract | CrossRef Full Text | Google Scholar

Uryu, S. (2004). Electronic states and quantum transport in double-wall carbon nanotubes. Phys. Rev. B 69, 075402. doi:10.1103/PhysRevB.69.075402

CrossRef Full Text | Google Scholar

Vidarte, K. J. U., and Lewenkopf, C. (2022). High-energy Landau levels in graphene beyond nearest-neighbor hopping processes: corrections to the effective Dirac Hamiltonian. Phys. Rev. B 106, 155414. doi:10.1103/PhysRevB.106.155414

CrossRef Full Text | Google Scholar

Vidarte, K. J. U., and Lewenkopf, C. (2024). Magnetic properties of low-angle twisted bilayer graphene at three-quarters filling. doi:10.48550/arXiv.2404.08177

CrossRef Full Text | Google Scholar

Yan, C., Ma, D.-L., Qiao, J.-B., Zhong, H.-Y., Yang, L., Li, S.-Y., et al. (2019). Scanning tunneling microscopy study of the quasicrystalline 30° twisted bilayer graphene. 2D Mater. 6, 045041. doi:10.1088/2053-1583/ab3b16

CrossRef Full Text | Google Scholar

Yankowitz, M., Chen, S., Polshyn, H., Zhang, Y., Watanabe, K., Taniguchi, T., et al. (2019). Tuning superconductivity in twisted bilayer graphene. Science 363, 1059–1064. doi:10.1126/science.aav1910

PubMed Abstract | CrossRef Full Text | Google Scholar

Yao, W., Wang, E., Bao, C., Zhang, Y., Zhang, K., Bao, K., et al. (2018). Quasicrystalline 30° twisted bilayer graphene as an incommensurate superlattice with strong interlayer coupling. Proc. Natl. Acad. Sci. U.S.A. 115, 6928–6933. doi:10.1073/pnas.1720865115

PubMed Abstract | CrossRef Full Text | Google Scholar

Yu, G., Katsnelson, M. I., and Yuan, S. (2020a). Pressure and electric field dependence of quasicrystalline electronic states in twisted bilayer graphene. Phys. Rev. B 102, 045113. doi:10.1103/PhysRevB.102.045113

CrossRef Full Text | Google Scholar

Yu, G., Wang, Y., Katsnelson, M. I., Lin, H.-Q., and Yuan, S. (2022). Interlayer hybridization in graphene quasicrystal and other bilayer graphene systems. Phys. Rev. B 105, 125403. doi:10.1103/PhysRevB.105.125403

CrossRef Full Text | Google Scholar

Yu, G., Wu, Z., Zhan, Z., Katsnelson, M. I., and Yuan, S. (2019). Dodecagonal bilayer graphene quasicrystal and its approximants. npj Comput. Mater. 5, 122. doi:10.1038/s41524-019-0258-0

CrossRef Full Text | Google Scholar

Yu, G., Wu, Z., Zhan, Z., Katsnelson, M. I., and Yuan, S. (2020b). Electronic structure of twisted double bilayer graphene. Phys. Rev. B 102, 115123. doi:10.1103/PhysRevB.102.115123

CrossRef Full Text | Google Scholar

Keywords: quasicrystals, twisted bilayer graphene, electronic structure, Landau level spectrum, Stampfli tiles, HHK recursive technique

Citation: Vidarte KJU and Lewenkopf C (2024) Quasicrystalline 30° twisted bilayer graphene: fractal patterns and electronic localization properties. Front. Carbon 3:1496179. doi: 10.3389/frcrb.2024.1496179

Received: 13 September 2024; Accepted: 22 October 2024;
Published: 05 November 2024.

Edited by:

Kaihui Liu, Peking University, China

Reviewed by:

Julian David Correa, Universidad de Medellín, Colombia
Habil. Maria Brzhezinskaya, Helmholtz Association of German Research Centers (HZ), Germany

Copyright © 2024 Vidarte and Lewenkopf. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

*Correspondence: Kevin J. U. Vidarte, am9obm1hcnVydmlAZ21haWwuY29t

Disclaimer: All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article or claim that may be made by its manufacturer is not guaranteed or endorsed by the publisher.