Research Article | Open Access

Miao-Ling Lin, Min Feng, Jiang-Bin Wu, Fei-Rong Ran, Tao Chen, Wei-Xia Luo, Heng Wu, Wen-Peng Han, Xin Zhang, Xue-Lu Liu, Yang Xu, Hai Li, Yu-Fang Wang, Ping-Heng Tan, "Intralayer Phonons in Multilayer Graphene Moiré Superlattices", *Research*, vol. 2022, Article ID 9819373, 11 pages, 2022. https://doi.org/10.34133/2022/9819373

# Intralayer Phonons in Multilayer Graphene Moiré Superlattices

#### Abstract

Moiré pattern in twisted multilayers (tMLs) induces many emergent phenomena by subtle variation of atomic registry to modulate quasiparticles and their interactions, such as superconductivity, moiré excitons, and moiré phonons. The periodic superlattice potential introduced by moiré pattern also underlies patterned interlayer coupling at the interface of tMLs. Although this arising patterned interfacial coupling is much weaker than in-plane atomic interactions, it is crucial in moiré systems, as captured by the renormalized interlayer phonons in twisted bilayer transitional metal dichalcogenides. Here, we determine the quantitative relationship between the lattice dynamics of intralayer out-of-plane optical (ZO) phonons and patterned interfacial coupling in multilayer graphene moiré superlattices (MLG-MS) by the proposed perturbation model, which is previously challenging for MLGs due to their out-of-phase displacements of adjacent atoms in one atomic plane. We unveil that patterned interfacial coupling introduces profound modulations on Davydov components of nonfolded ZO phonon that are localized within the AB-stacked constituents, while the coupling results in layer-extended vibrations with symmetry of moiré pattern for moiré ZO phonons. Our work brings further degrees of freedom to engineer moiré physics according to the modulations imprinted on the phonon frequency and wavefunction.

#### 1. Introduction

The misorientation or lattice mismatch between the constituents of van der Waals (vdW) heterostructures (vdWHs) enables the generation of moiré superlattices, which impose nanoscale periodic modulations on the electronic band [1, 2] and optical selection rules [3], and lead to correlated electronic phases [4, 5], superconductivity [6, 7], and moiré excitons [3, 8–10]. Particularly, the interfacial coupling is also modulated by the nanopatterned periodic potential of moiré superlattices, leading to patterned interfacial coupling (PIC) and thus a platform of versatile phonon tunability [11, 12]. Recently, the roles of phonons in many emerging quantum phenomena of vdWHs have been emphasized [13–16]. The control of phonon and its interaction with electronic phases in moiré structures is considered as potential key point for many exotic properties absent in the constituents, such as insulator-metal transition and ferromagnetism [11, 12, 17].

In the light of the above long-pursued goals, the modulation of moiré superlattices and the corresponding PIC on phonon properties have attracted much attention [11, 12, 18, 19]. Previous studies have identified moiré phonons in twisted bilayer transitional metal dichalcogenides (TMDs) [18] that are activated by phonon folding effect in moiré superlattices, whose vibrations are engineered by atomic registry. Meanwhile, the interlayer phonons in reconstructed moiré superlattices of twisted bilayer are renormalized owing to the strong coupling with moiré acoustic phonons [12]. However, in moiré superlattices of graphene or hexagonal boron nitride (hBN) whose monolayer is comprised of only one atomic plane, the imprints of PIC on their intralayer phonons are inherently complex because the displacements of adjacent atoms in one atomic plane are out-of-phase, and those in the adjacent atomic planes vary from local to local due to twist-angle dependent interlayer registries. Recent nano-Raman spectroscopy [11] and nanoinfrared imaging [19] revealed localized lattice dynamics of nonfolded in-plane intralayer phonons in low-angle twisted bilayer graphene (tBLG) and twisted hBN, respectively. In large-angle moiré superlattices beyond the spatial limit of nanotechnologies, the evolution of intralayer phonons upon modulation of PIC remains appealing, especially in multilayer graphene (MLG) moiré superlattices (MLG-MS) with additional and tunable degrees of freedom in quantum manipulation and device characteristics [5, 20–22] owing to the coexistence of PIC and natural interlayer coupling. And two general open questions stand out: (a) How do moiré superlattices modulate the lattice dynamics of intralayer phonons? (b) How does PIC affect the vibrations related to intralayer phonons in the graphene layers that are not adjacent to the moiré interface?

Herein, we reveal the modulation of PIC on the frequency and atomic displacements of intralayer phonons by the resonance Raman spectroscopy of MLG-MS containing both twisted and AB-stacked interfaces. This is manifested in the renormalized Davydov splitting for both intralayer out-of-plane optical (ZO) phonons at the unfolded Brillouin zone (BZ) center, *i.e.*, nonfolded ZO phonons and moiré ZO (mZO) phonons folded from the off-center phonons of untwisted constituents. The renormalized Davydov components of nonfolded ZO phonons are localized within the AB-stacked constituents due to the zero perturbation from PIC, which is rationalized by a perturbation model (PM) and confirmed by the force constant method (FCM). However, the PIC results in vibrations extended to all the stacking layers for the Davydov components of mZO phonons, where the vibration pattern shows similar symmetry to moiré superlattices. This work enhances the understanding about the effect of PIC on phonons in moiré superlattices, pointing to potential opportunities for phonon manipulation.

#### 2. Results

##### 2.1. Enhanced Davydov Components of ZO Modes in MLG-MS

Moiré phonons linked with the reciprocal unit vectors of the moiré superlattices within the BZ interior of untwisted constituents can be activated by the phonon folding effect in the Raman spectroscopy of twisted bilayers [18, 23–27], where rich resonance mechanisms of moiré phonon modes were reported [18, 26, 28]. Similar phenomena are also expected in MLG-MS, which are assembled by -layer graphene (LG) and LG and are denoted as (See Methods, Supplementary Figure S1 and Figure S2). If , LG (LG) is AB-stacked. The Raman spectrum of in Figure 1(a) is resonantly excited by 1.96 eV, and its mode intensity () is ~20 times stronger than that of 1LG, similar to the giant enhancement of G mode in BLG [28–31]. Closer inspection shows that the G band exhibits obvious doubling, which is common in ( or 1) and is connected to a resonant process determined by twist angle [30]. The relative Raman intensity between the two subpeaks is sensitive to the excitation energy. The characteristic R peak [23] below the mode originates from the moiré phonon associated with the transverse optical (TO) phonon branch, which is denoted as mTO phonon here. This peak is observed at ~1510 , which corresponds to a twist angle () of 11.3° based on the -dependent mTO phonon frequency [26, 32]. Figures 1(b) and 1(c) show the stacking schematics of with (also in Supplementary Figure S3) and the corresponding reciprocal lattices, respectively. The lengths of basic wave vectors of moiré superlattices with are defined by and [18, 23]. The phonons at () points determined by within the BZ of the constituents, moiré phonons, are folded back to the BZ center of moiré superlattices. Because shows symmetry [33], its in-plane and out-of-plane vibrations become Raman-active and are assigned to and irreducible representations, respectively. Thus, the moiré phonons are expected to be observed in Raman spectra. Besides the mTO or R mode, additional Raman peaks are also observed in in the spectral range of 400–900 , when compared with those in 1LG and 3LG. Based on the phonon dispersion [34] and the deduced , the extra peaks in at ~408 cm^{-1} and ~602 cm^{-1} are assigned to the moiré phonons related to transverse and longitudinal acoustic phonon branches, and they are denoted as mTA and mLA, respectively. Furthermore, the modes at ~844 cm^{-1} and ~870 cm^{-1} are, respectively, moiré (*i.e.*, mZO) phonon and nonfolded phonons linked with the ZO phonon branch. The mTA and mLA modes are observed in parallel (VV) and crossed (HV) configurations, whereas mZO and nonfolded ZO phonons appear only in the VV configuration. This further confirms the above assignments. In contrast to one Raman peak related to each moiré phonon and no feature linked with nonfolded ZO mode in tBLG [26, 32], three peaks with a width of ~1 cm^{-1} each are clearly resolved for both mZO and nonfolded ZO modes in .

Fundamentally, the optical modes in LG undergo Davydov splitting due to the weak interlayer coupling and split up into Davydov components [35]. For example, the silent ZO mode () in 1LG is expected to split up into two silent and one Raman-active modes for 3LG, although none of them are observed in Raman spectroscopy due to Raman-inactivity or weak electron-phonon coupling. Considering the comparable interfacial and interlayer breathing coupling in MLG-MS [33], four Davydov components are expected for ZO phonons in , inconsistent with the three Davydov-split peaks of nonfolded ZO and mZO modes observed in experiments. Notably, the mZO and ZO phonon modes can only be resonantly detected in a narrow excitation energy range.

To understand the observed Davydov splitting of the nonfolded ZO and mZO phonons in MLG-MS and their modulation by interlayer coupling at AB-stacked and twisted interfaces, we measure the Raman spectra of other under the resonant excitation conditions, as plotted in Figure 2(a). As the commonly used laser in the visible region cannot match the van Hove singularities in electronic joint density of states of all the optically allowed transitions in with close to 0° and 30° [26], in this work, we only present the results in with . Similar to (Figure 1(a)), the G band in ( or ) consists of two subbands. The Raman peaks at ~870 cm^{-1} originate from nonfolded ZO phonon, while the modes with peak positions varying with samples are mZO phonon modes, e.g., peaks at ~844 cm^{-1} in and t(1 + 3)LG with and peaks at ~855 cm^{-1} in and with . Except for , the nonfolded ZO and mZO modes in are both activated and exhibit Davydov splitting, although (mZO) is much weaker than (ZO). Interestingly, the number of observed Davydov components of nonfolded ZO modes in is directly related to the number of layers of its constituents, but not equal to its total number (*i.e.*, ). For example, we observe three and two Davydov components of nonfolded ZO modes in and , respectively, rather than four Davydov components. This suggests that the interfacial coupling in is too weak to induce Davydov components of nonfolded ZO phonons. Furthermore, the frequency difference between the two adjacent Davydov components stays constant (2.0 cm^{-1}) in with different (Figure 2(b)).

##### 2.2. Localized Davydov Components of Nonfolded ZO Modes

Davydov splitting of intralayer phonons induced by interlayer coupling in multilayer and twisted TMDs is well represented by the vdW model [36–38], where the interlayer coupling can be treated as an overall force due to the in-phase atomic displacements in one atomic plane. The frequency differences () between each Davydov component and the lowest-frequency one are closely related to the frequency of layer-breathing (LB) phonons. We tried to extend this vdW model to the nonfolded ZO phonons in both with and without interfacial LB coupling considered (Supplementary Figure S4), while obvious discrepancy is present between the calculated and experimental results. This suggests that the vdW model is not valid in understanding the Davydov splitting of nonfolded ZO phonons in MLG-MS, where the adjacent atoms in one atomic plane are out-of-phase and interlayer coupling cannot be assumed as a residual restoring force to all atoms. Looking for insights into the Davydov splitting of nonfolded ZO phonons, the Davydov components are derived from the linear superpositions of ZO phonons in each graphene layer. The influence from weak interlayer coupling can be assumed as a perturbation to such equivalent ZO phonons, which is not only associated with periodic potential of PIC but also the atomic displacements of nonfolded ZO phonon. In this case, we develop a model from perspective of perturbation to understand the Davydov splitting of nonfolded ZO phonons in , where the unperturbed mode in 1LG is treated as the basis, and the perturbation from interlayer (interfacial) coupling is (), as schematized in Figure 2(c).

Taking as an example, the reduced Hamiltonian is
where is the energy of an unperturbed ZO phonon, *i.e.*, in 1LG. is set as 870.9 to match the nonfolded ZO mode frequencies to the experimental ones in and other . is introduced to follow the change of the on-site energy in the middle layer due to the new geometry, similar to the on-site energy variation due to environmental changes in the electronic tight-binding model [39]. Considering that the number of observed Davydov components of nonfolded ZO modes in is related to and rather than (Figure 2(a)), the impact of interfacial coupling on the nonfolded ZO modes is negligible. Thus, we assume , and then, we solve the secular equation det. According to the observed splitting frequency of 3.16 in , and are deduced as 1.1 and 2.3 , respectively. We find the positive and negative have reasonable physical significance and can reproduce all the experimental results, as discussed below. The corresponding eigenvectors are , , and . The vibration corresponding to is localized in 1LG (), which cannot be observed due to its Raman inactivity in 1LG. The latter two eigenvectors correspond to the two observed Davydov components in , *i.e.*, (), and frequency increases with . and modes are associated with the out-of-phase and in-phase superpositions of the two unperturbed ZO modes, respectively.

Similar analysis with can be extended to understand the Davydov splitting of nonfolded ZO phonon in other . The corresponding calculated eigenvectors (*i.e.*, atomic displacements) are shown in Figures 2(d)–2(f)). In , two equivalent modes with their vibrations localized in the top or bottom 1LG are predicted by the PM. Thus, no ZO mode is activated in , in line with the experimental result (Figure 2(a)). In addition, there are two Davydov components of nonfolded ZO phonons in , each doubly degenerate (Supplementary Figure S5). For , four eigenvectors are obtained, among which () localizes in the 1LG constituent and the other three localize within the 3LG constituents. With the abovedetermined values of (1.1 ) and (-2.3 ), between and components and between and components are, respectively, 2.1 and 4.1 , which agree well with the experimental data, as elucidated in Figure 2(b). In particular, the PM also reproduces the five Davydov components of nonfolded ZO modes in , among which the two highest-frequency components exhibit almost the same frequency. In this case, four Davydov components are resolved in Raman spectroscopy. These results indicate the validity of PM to understand the localized Davydov components of nonfolded ZO phonons within the constituents of , including the assumption of negligible perturbation from PIC (*i.e.*, ). Only one is used in the PM for the middle layers in MLG-MS, no matter they are adjacent to the twisted or AB-stacked interface. Notably, and are closely linked with both the PIC and atomic displacements of nonfolded ZO phonon in one atomic plane. Although the moiré pattern in MLG-MS largely reduces the perturbation from interfacial coupling on nonfolded ZO phonons due to the planar locally mismatched periodicity of the charge density variation [33], it keeps the on-site energy variation in the two interfacial layers of MLG-MS the same as that in AB stacking, which may be ascribed to the comparable LB coupling at AB-stacked and twisted interfaces. Thus, the PIC in MLG-MS significantly modulates the lattice dynamics of nonfolded ZO phonons.

##### 2.3. Renormalized Davydov Components of Nonfolded ZO Modes by PIC

In principle, the PM with the abovedetermined values of and is also applicable to estimate the Davydov components of ZO phonons in AB-stacked MLG. For example, two Davydov components at 869.8 and 872 are predicted in 2LG (Figure 3(a)), while three Davydov components are anticipated in 3LG (Figure 3(b)). One would intuitively expect that the Davydov splitting in is the same as those in AB-stacked LG and LG due to the negligible effect of interfacial coupling on nonfolded ZO phonons. However, the case is different due to the on-site energy variation, as exemplified by the comparison between 2LG, 3LG, and . In (Figure 3(c)), the reduced Hamiltonian can be divided into two nonzero submatrices (as divided by the dotted lines) related to the 2LG and 3LG constituents due to , whereas these two submatrices are different from those in pristine 2LG and 3LG due to the presence of in the interfacial layers of the 2LG and 3LG constituents. Thus, the resulting frequencies and atomic displacements of Davydov components related to nonfolded ZO phonon in (Figure 3(c)) are distinguished from those in pristine 2LG (Figure 3(a)) and 3LG (Figure 3(b)). In addition, when the top 2LG and bottom 3LG are assembled in AB stacking at the twisted interface and a 5LG is formed, the perturbation from the interfacial coupling is . Thus, the ZO phonon splits up into five Davydov components, which stems from the vibrations of the five stacking layers, as elucidated in Figure 3(d). Their frequencies also show obvious distinction from those in . These contrasts imply profound modulations of the PIC on the lattice dynamics of nonfolded ZO phonons in MLG-MS.

More imprints for the PIC on nonfolded ZO phonons can be found in the atomic displacements calculated by the FCM by considering the Born-von Karman model with the interactions along the radial and tangential directions (see Materials and Methods) [40, 41]. For calculation simplicity, we take with as an example. A supercell of the atomic structure for 2LG and the vibration patterns of are plotted (Figures 4(a) and 4(b)). in 2LG exhibits similar atomic displacements to except for in-phase superposition of the two nonfolded ZO phonons in each layer. For , the calculation presents two Davydov components in the 2LG constituent due to the interlayer coupling and one mode with vibration localized in the 1LG constituent, as shown in Figures 4(c)–4(e)). This confirms that the interfacial coupling is too weak to couple the three graphene layers for nonfolded ZO phonons. Additionally, for the and modes in , the vibration amplitude ratios of the superposition atoms between the two bottom layers are 1.5:(1) and 1 : 1.5, respectively, distinct from those in pristine 2LG (Figure 4(b)). This underpins the modulation of moiré superlattices on atomic displacements of nonfolded ZO phonons in MLG-MS. The above analysis based on the PM and FCM is universally valid to estimate the Davydov components of intralayer phonons and the effects of interlayer/interfacial coupling on them in other vdWHs.

We notice that the localized and renormalized Davydov components linked with the nonfolded ZO phonon are ubiquitously observed in MLG-MS under resonance excitation, as demonstrated by the Raman spectra of , , and with various in Figure 5. For , and components with a splitting frequency of are present. As for , three peaks (*i.e.*, , , and components) are observed when , while two distinct Davydov components are observed when , whose splitting frequency is quite close to the estimated between and components and between and components by PM, respectively. Four Davydov components are observed in , where each frequency difference between two given Davydov components varies slightly with twist angle. Despite the Davydov splitting in MLG-MS is independent of , the resonant excitation energy and the relative Raman intensity of Davydov components are closely associated with , pointing to the twist-angle dependent electronic transition and electron-phonon coupling. It is also worth noting that the reduced symmetry by twisting and resonance excitation are prerequisites for the observations of Davydov components of ZO phonons in MLG-MS.

##### 2.4. Layer-Extended Davydov Components of mZO Modes Modulated by Moiré Pattern

We further look for insights into the effects of the PIC in moiré pattern on the mZO phonons. The maximum number of observed Davydov components linked with mZO phonons in can be up to . For specific , the number of observed Davydov components is dependent on and laser excitation, as shown in Figures 6(a) and 6(b). We first assume that the observed splitting peaks in only correspond to Davydov components of mZO phonon in its constituents due to the weak impact from interfacial coupling. Following this hypothesis and the phonon dispersion calculated by FCM (Figure 6(c)), we should observe three peaks of mZO phonons in a , in which the one corresponding to the mZO phonon from 1LG constituent shows frequency between the two Davydov components of mZO phonons in 2LG. However, by comparing the mZO peaks in with at the same of in Figure 6(a), we find that the frequency of the mZO peak in approaches the highest-frequency components of mZO mode in . This implies the frequencies of the Davydov components of mZO peaks in MLG-MS are altered by the interfacial coupling. When providing further insights into the frequency differences between each Davydov entity and the lowest-frequency one, as shown in Figures 6(d) and 6(e), we find that each stays constant for with ranging from 7.9 to 16.4. And it is also the case in . This indicates that the perturbations from the interlayer and interfacial couplings vary slightly with large in MLG-MS. We also extended the proposed PM to understand the Davydov components of mZO phonons in and to estimate the perturbations from interlayer and interfacial couplings. We found the calculated results (dashed lines) based on PM roughly agree with the experimental results when and . Hence, the perturbations from interlayer coupling and PIC on mZO phonons in are equivalent; thus, the atomic displacements of mZO Davydov components extend to all the stacking layers in MLG-MS. This is distinct from the case of nonfolded ZO phonons.

Clear signatures of the modulation from the moiré pattern on the mZO Davydov components are found in the atomic displacements calculated by the FCM, as demonstrated in Figure 6(f) (also in Supplementary Figure S6 and Figure S7). The calculations demonstrate that the mZO Davydov components stem from the vibrations of the three stacking graphene layers in , in consistent with the predictions from PM. More interestingly, the vibration patterns at each layer show symmetry. However, the corresponding unfolded ZO() phonons in 3LG exhibit distinctly different vibration patterns, without any signature of symmetry, as shown in Figure 6(g). This further confirms the crucial role of the moiré pattern in the modulation for both frequency and lattice dynamics of mZO phonons in MLG-MS. The excitation-energy dependent mZO Davydov components observed in with different (Figures 6(a) and 6(b)) imply the -dependent modulation of the moiré pattern on electronic properties, moiré phonons, and the corresponding electron-phonon coupling.

#### 3. Discussion

We have observed the Davydov splitting of ZO and mZO modes renormalized by PIC in MLG-MS using resonance Raman spectroscopy. The Davydov components of nonfolded ZO phonons are localized within the multilayer constituents due to the negligible perturbation from interfacial coupling, which is well represented by PM and FCM. In contrast, all the stacking layers are coupled by interfacial coupling for the vibrations of mZO phonons, and the corresponding atomic displacements show similar symmetry to the moiré pattern. The distinct difference of the vibration patterns between the Davydov components of nonfolded ZO and mZO phonon modes underpins the crucial role of PIC in the modulation of lattice dynamics for intralayer phonons in MLG-MS. This also provides potential routes to separately control the properties of nonfolded and moiré intralayer phonons. Furthermore, the rich and overlapped modulations of moiré superlattices on varied quasiparticles and the related interactions, e.g., electrons, excitons, and intralayer/interlayer phonons, will broaden the prospects for device applications of moiré materials. Remarkably, this work demonstrates that common Raman spectroscopy can provide information on the renormalized lattice dynamics in nanoscale moiré landscape of large-angle MLG-MS beyond the limit of nanotechnologies, which can be extended to other multilayer moiré superlattices. Our findings provide new insights into the rich and complex phonon physics in moiré structures, encouraging more theoretical and experimental works on this topic.

#### 4. Materials and Methods

##### 4.1. Sample Preparations

MLG flakes are mechanically exfoliated from the highly oriented pyrolytic graphite onto a wafer. During the exfoliation, an -layer graphene (LG) may be folded onto an LG randomly and thus form the [42], such as the () in Supplementary Figure S1(a). Alternatively, an LG on one substrate can also be transferred onto an LG on another substrate to form the . [43] The samples and with the same twisted angle () are prepared in this way (Supplementary Figure S1(b)). After transferring, the samples were annealed in flowing H/N 300 C gas for two hours to remove the residues. In addition, the studied () and () in Figure 4 were directly grown by the low-pressure chemical vapor deposition growth [44]. The layer number of MLG flakes was identified by the Raman spectra and the optical contrast (Supplementary Figure S2) [45–47]. The twist angle between the two constituents of MLG-MS is determined by comparing the characteristic mTO(R) mode frequency (below the G band) with the phonon dispersion along direction [34]. Additionally, the twist angle is confirmed by the -dependent resonance energies and mTO(R) phonon frequencies [26, 32]. Taking with as an example, Supplementary Figure S3 presents a schematic diagram for its atomic structure, which shows evident moiré pattern.

##### 4.2. Raman Measurements

The Raman spectroscopy was measured under backscattering configuration at room temperature, using the Jobin-Yvon HR800 system equipped with a liquid nitrogen cooled charge-coupled (CCD) detector, a 100× objective len and several gratings. The excitation energies are 1.58 eV and 1.71 eV from a Ti: Sapphire laser, 1.96 eV and 2.09 eV from a He-Ne laser, 1.85 eV from a diode pumped solid-state laser, 1.83 eV, 2.18 eV, 2.34 eV, and 2.41 eV from a Kr laser, and 2.54 eV, 2.60 eV, and 2.71 eV from an Ar laser. The resolution at 2.60 eV is 0.07 cm^{-1} per CCD pixel. Plasma lines are removed from laser signals by using BragGrate Bandpass filters (OptiGrate Corp.). The typical laser power is 1 mW to avoid heating. The acquisition time for each spectrum is 1200 seconds to enhance the signal-to-noise ratio.

##### 4.3. Force Constant Methods

The theoretical calculation is based on the empirical force constant method (FCM) by considering the Born-von Karman model of the lattice dynamics for atomic coupling. By solving the crystal dynamics described by the phonon dispersions of various samples, such as graphene, MLG, and can be obtained. In this equation, is the mass of the atom in the unit cell; and enumerate the Cartesian coordinates , , and ; specifies the displacement of the atom in the unit cell; and (,) is the atomic force constants that is the second derivative of the crystal potential energy with respect to the atomic displacements taken at the equilibrium position, including the intralayer couplings and the interlayer interactions. By applying Born-von Karman model to describe the carbon-carbon intralayer interaction and taking four nearest neighbor atoms into account, the fitting parameters of the force constant can be given by the phonon dispersion of graphene. Then, the radial and tangential interlayer force constants are considered to reproduce the phonon dispersion along direction of bulk graphite [40]. With these fitting parameters for the empirical force-constant model, we can easily access the varied phonon dispersions of and the atomic displacements of all the phonons in the BZ.

#### Data Availability

All data used to support the findings of this study are available from the corresponding author upon request.

#### Conflicts of Interest

The authors declare that there is no conflict of interest regarding the publication of this article.

#### Authors’ Contributions

P.T. conceived the idea and directed and supervised the project; M.L. and P.T. designed the experiments; F.R., W.L., H.W., W.H., Y.X., and H.L. prepared the samples; M.L. performed the experiments; M.F. and Y.W. calculated the atomic displacements in twisted multilayer graphene by force constant model; P.T. and M.L. analyzed the data with inputs from J.W., T.C., X.Z., and X.L.; and P.T. and M.L. developed the perturbation model and wrote the manuscript with input from all authors.

#### Acknowledgments

This work was supported by the National Natural Science Foundation of China (Grant nos. 12004377, 11874350, and 21571101) and CAS Key Research Program of Frontier Sciences (Grant nos. ZDBS-LY-SLH004 and XDPB22).

#### Supplementary Materials

Figure S1: optical images of . Figure S2: optical contrast and Raman spectra of . Figure S3: moiré pattern in . Figure S4: frequency differences between Davydov components of nonfolded ZO modes based on vdW model. Figure S5: Davydov splitting of nonfolded ZO modes in and . Figure S6: atomic displacements of Davydov components of mZO modes in . Figure S7: atomic displacements of Davydov components of the corresponding unfolded ZO modes in 3LG.* (Supplementary Materials)*

#### References

- R. Bistritzer and A. H. MacDonald, “Moiré bands in twisted double-layer graphene,”
*Proceedings of the National Academy of Sciences*, vol. 108, no. 30, pp. 12233–12237, 2011. View at: Publisher Site | Google Scholar - H. Yoo, R. Engelke, S. Carr et al., “Atomic and electronic reconstruction at the van der Waals interface in twisted bilayer graphene,”
*Nature Materials*, vol. 18, no. 5, pp. 448–453, 2019. View at: Publisher Site | Google Scholar - H. Yu, G.-B. Liu, J. Tang, X. Xu, and W. Yao, “Moiré excitons: From programmable quantum emitter arrays to spin-orbit–coupled artificial lattices,”
*Science Advances*, vol. 3, no. 11, article e1701696, 2017. View at: Publisher Site | Google Scholar - Y. Cao, V. Fatemi, A. Demir et al., “Correlated insulator behaviour at half-filling in magic-angle graphene superlattices,”
*Nature*, vol. 556, no. 7699, pp. 80–84, 2018. View at: Publisher Site | Google Scholar - X. Liu, Z. Hao, E. Khalaf et al., “Tunable spin-polarized correlated states in twisted double bilayer graphene,”
*Nature*, vol. 583, no. 7815, pp. 221–225, 2020. View at: Publisher Site | Google Scholar - Y. Cao, V. Fatemi, S. Fang et al., “Unconventional superconductivity in magic-angle graphene superlattices,”
*Nature*, vol. 556, no. 7699, pp. 43–50, 2018. View at: Publisher Site | Google Scholar - J. M. Park, Y. Cao, K. Watanabe, T. Taniguchi, and P. Jarillo-Herrero, “Tunable strongly coupled superconductivity in magic-angle twisted trilayer graphene,”
*Nature*, vol. 590, no. 7845, pp. 249–255, 2021. View at: Publisher Site | Google Scholar - K. Tran, G. Moody, F. Wu et al., “Evidence for moiré excitons in van der Waals heterostructures,”
*Nature*, vol. 567, no. 7746, pp. 71–75, 2019. View at: Publisher Site | Google Scholar - C. Jin, E. C. Regan, A. Yan et al., “Observation of moiré excitons in WSe
_{2}/WS_{2}heterostructure superlattices,”*Nature*, vol. 567, no. 7746, pp. 76–80, 2019. View at: Publisher Site | Google Scholar - K. L. Seyler, P. Rivera, H. Yu et al., “Signatures of moiré-trapped valley excitons in MoSe
_{2}/WSe_{2}heterobilayers,”*Nature*, vol. 567, no. 7746, pp. 66–70, 2019. View at: Publisher Site | Google Scholar - A. C. Gadelha, D. A. A. Ohlberg, C. Rabelo et al., “Localization of lattice dynamics in low-angle twisted bilayer graphene,”
*Nature*, vol. 590, no. 7846, pp. 405–409, 2021. View at: Publisher Site | Google Scholar - J. Quan, L. Linhart, M.-L. Lin et al., “Phonon renormalization in reconstructed MoS
_{2}moiré superlattices,”*Nature Materials*, vol. 20, no. 8, pp. 1100–1105, 2021. View at: Publisher Site | Google Scholar - C. Jin, J. Kim, J. Suh et al., “Interlayer electron-phonon coupling in WSe
_{2}/hBN heterostructures,”*Nature Physics*, vol. 13, no. 2, pp. 127–131, 2017. View at: Publisher Site | Google Scholar - K. Chaudhary, M. Tamagnone, M. Rezaee et al., “Engineering phonon polaritons in van der waals heterostructures to enhance in-plane optical anisotropy,”
*Science Advances*, vol. 5, no. 4, article eaau7171, 2019. View at: Publisher Site | Google Scholar - M.-L. Lin, Y. Zhou, J.-B. Wu et al., “Cross-dimensional electron-phonon coupling in van der Waals heterostructures,”
*Nature Communications*, vol. 10, no. 1, p. 2419, 2019. View at: Publisher Site | Google Scholar - P. Merkl, C. K. Yong, M. Liebich et al., “Proximity control of interlayer exciton-phonon hybridization in van der Waals heterostructures,”
*Nature Communications*, vol. 12, no. 1, p. 1719, 2021. View at: Publisher Site | Google Scholar - Y. Tang, L. Li, T. Li et al., “Simulation of Hubbard model physics in WSe
_{2}/WS_{2}moiré superlattices,”*Nature*, vol. 579, no. 7799, pp. 353–358, 2020. View at: Publisher Site | Google Scholar - M.-L. Lin, Q.-H. Tan, J.-B. Wu et al., “Moiré phonons in twisted bilayer MoS
_{2},”*ACS Nano*, vol. 12, no. 8, pp. 8770–8780, 2018. View at: Publisher Site | Google Scholar - S. L. Moore, C. J. Ciccarino, D. Halbertal et al., “Nanoscale lattice dynamics in hexagonal boron nitride moiré superlattices,”
*Nature Communications*, vol. 12, no. 1, p. 5714, 2021. View at: Publisher Site | Google Scholar - G. Chen, A. L. Sharpe, P. Gallagher et al., “Signatures of tunable superconductivity in a trilayer graphene moiré superlattice,”
*Nature*, vol. 572, no. 7768, pp. 215–219, 2019. View at: Publisher Site | Google Scholar - P. Rickhaus, F. K. de Vries, J. Zhu et al., “Correlated electron-hole state in twisted double-bilayer graphene,”
*Science*, vol. 373, no. 6560, pp. 1257–1260, 2021. View at: Publisher Site | Google Scholar - M. He, Y. Li, J. Cai et al., “Symmetry breaking in twisted double bilayer graphene,”
*Nature Physics*, vol. 17, no. 1, pp. 26–30, 2021. View at: Publisher Site | Google Scholar - V. Carozo, C. M. Almeida, E. H. M. Ferreira, L. G. Cançado, C. A. Achete, and A. Jorio, “Raman signature of graphene superlattices,”
*Nano Letters*, vol. 11, no. 11, pp. 4527–4534, 2011. View at: Publisher Site | Google Scholar - R. He, T.-F. Chung, C. Delaney et al., “Observation of low energy Raman modes in twisted bilayer graphene,”
*Nano Letters*, vol. 13, no. 8, pp. 3594–3601, 2013. View at: Publisher Site | Google Scholar - P. Ramnani, M. R. Neupane, S. Ge, A. A. Balandin, R. K. Lake, and A. Mulchandani, “Raman spectra of twisted CVD bilayer graphene,”
*Carbon*, vol. 123, pp. 302–306, 2017. View at: Publisher Site | Google Scholar - G. S. Eliel, M. V. Moutinho, A. C. Gadelha et al., “Intralayer and interlayer electron-phonon interactions in twisted graphene heterostructures,”
*Nature Communications*, vol. 9, no. 1, p. 1221, 2018. View at: Publisher Site | Google Scholar - P. Parzefall, J. Holler, M. Scheuck et al., “Moiré phonons in twisted MoSe
_{2}–WSe_{2}heterobilayers and their correlation with interlayer excitons,”*2D Materials*, vol. 8, no. 3, article 035030, 2021. View at: Publisher Site | Google Scholar - A. Righi, P. Venezuela, H. Chacham et al., “Resonance Raman spectroscopy in twisted bilayer graphene,”
*Solid State Communications*, vol. 175-176, pp. 13–17, 2013. View at: Publisher Site | Google Scholar - V. Carozo, C. M. Almeida, B. Fragneaud et al., “Resonance effects on the raman spectra of graphene superlattices,”
*Physical Review B*, vol. 88, no. 8, article 085401, 2013. View at: Publisher Site | Google Scholar - J.-B. Wu, X. Zhang, M. Ijäs et al., “Resonant Raman spectroscopy of twisted multilayer graphene,”
*Nature communications*, vol. 5, no. 1, p. 5309, 2014. View at: Publisher Site | Google Scholar - M. C. DeCapua, Y.-C. Wu, T. Taniguchi, K. Watanabe, and J. Yan, “Probing the bright exciton state in twisted bilayer graphene via resonant Raman scattering,”
*Applied Physics Letters*, vol. 119, no. 1, article 013103, 2021. View at: Publisher Site | Google Scholar - J. Campos-Delgado, L. G. Cançado, C. A. Achete, A. Jorio, and J. P. Raskin, “Raman scattering study of the phonon dispersion in twisted bilayer graphene,”
*Nano Research*, vol. 6, no. 4, pp. 269–274, 2013. View at: Publisher Site | Google Scholar - J.-B. Wu, Z.-X. Hu, X. Zhang et al., “Interface coupling in twisted multilayer graphene by resonant Raman spectroscopy of layer breathing modes,”
*ACS Nano*, vol. 9, no. 7, pp. 7440–7449, 2015. View at: Publisher Site | Google Scholar - P. Venezuela, M. Lazzeri, and F. Mauri, “Theory of double-resonant raman spectra in graphene: intensity and line shape of defect-induced and two-phonon bands,”
*Physical Review B*, vol. 84, no. 3, article 035433, 2011. View at: Publisher Site | Google Scholar - J.-B. Wu, M.-L. Lin, X. Cong, H.-N. Liu, and P.-H. Tan, “Raman spectroscopy of graphene-based materials and its applications in related devices,”
*Chemical Society Reviews*, vol. 47, no. 5, pp. 1822–1873, 2018. View at: Publisher Site | Google Scholar - Q.-J. Song, Q.-H. Tan, X. Zhang et al., “Physical origin of davydov splitting and resonant Raman spectroscopy of davydov components in multilayer MoTe
_{2},”*Physical Review B*, vol. 93, no. 11, article 115409, 2016. View at: Publisher Site | Google Scholar - Q.-H. Tan, X. Zhang, X.-D. Luo, J. Zhang, and P.-H. Tan, “Layer-number dependent high-frequency vibration modes in few-layer transition metal dichalcogenides induced by interlayer couplings,”
*Journal of Semiconductors*, vol. 38, no. 3, article 031006, 2017. View at: Publisher Site | Google Scholar - Y.-C. Leng, M.-L. Lin, Y. Zhou et al., “Intrinsic effect of interfacial coupling on the high-frequency intralayer modes in twisted multilayer MoTe
_{2},”*Nanoscale*, vol. 13, no. 21, pp. 9732–9739, 2021. View at: Publisher Site | Google Scholar - J. L. Mercer and M. Y. Chou, “Tight-binding model with intra-atomic matrix elements,”
*Physical Review B*, vol. 49, no. 12, pp. 8506–8509, 1994. View at: Publisher Site | Google Scholar - H. Wang, Y. F. Wang, X. W. Cao, M. Feng, and G. X. Lan, “Vibrational properties of graphene and graphene layers,”
*Journal of Raman Specroscopy*, vol. 40, no. 12, pp. 1791–1796, 2009. View at: Publisher Site | Google Scholar - A. I. Cocemasov, D. L. Nika, and A. A. Balandin, “Phonons in twisted bilayer graphene,”
*Physical Review B*, vol. 88, no. 3, article 035428, 2013. View at: Publisher Site | Google Scholar - K. S. Novoselov, D. Jiang, F. Schedin et al., “Two-dimensional atomic crystals,”
*Proceedings of the National Academy of Sciences*, vol. 102, no. 30, pp. 10451–10453, 2005. View at: Publisher Site | Google Scholar - C. R. Dean, A. F. Young, I. Meric et al., “Boron nitride substrates for high-quality graphene electronics,”
*Nature Nanotechnology*, vol. 5, no. 10, pp. 722–726, 2010. View at: Publisher Site | Google Scholar - Y. Song, D. Pan, Y. Cheng, P. Wang, P. Zhao, and H. Wang, “Growth of large graphene single crystal inside a restricted chamber by chemical vapor deposition,”
*Carbon*, vol. 95, pp. 1027–1032, 2015. View at: Publisher Site | Google Scholar - A. C. Ferrari, J. C. Meyer, V. Scardaci et al., “Raman spectrum of graphene and graphene layers,”
*Physical Review Letters*, vol. 97, no. 18, article 187401, 2006. View at: Publisher Site | Google Scholar - W. Zhao, P. Tan, J. Zhang, and J. Liu, “Charge transfer and optical phonon mixing in few-layer graphene chemically doped with sulfuric acid,”
*Physical Review B*, vol. 82, no. 24, article 245423, 2010. View at: Publisher Site | Google Scholar - X.-L. Li, W.-P. Han, J.-B. Wu, X.-F. Qiao, J. Zhang, and P.-H. Tan, “Layer-number dependent optical properties of 2D materials and their application for thickness determination,”
*Advanced Functional Materials*, vol. 27, no. 19, article 1604468, 2017. View at: Publisher Site | Google Scholar

#### Copyright

Copyright © 2022 Miao-Ling Lin et al. Exclusive Licensee Science and Technology Review Publishing House. Distributed under a Creative Commons Attribution License (CC BY 4.0).