Imaging conductive nano-domains induced by Gd intercalation in epitaxial bilayer graphene

Abstract

We report nano-infrared (IR) imaging and spectroscopy of epitaxial bilayer graphene (BLG) on silicon carbide (SiC) partially intercalated with gadolinium (Gd). Gd intercalation produces a high density of nanoscale conducting domains that exhibit pronounced IR enhancement at frequencies above the SiC phonon resonance and pronounced amplitude suppression at the resonance. Both effects originate from the increased local optical conductivity induced by Gd. Quantitative modeling of the nano-IR spectra shows that the conductivity of intercalated regions is enhanced by more than a factor of two relative to pristine BLG. This enhancement is attributed to the electronic decoupling of the graphene layers combined with substantial charge transfer from the intercalated atoms. These results demonstrate that controlled metal intercalation enables spatially resolved tuning of the electronic and optical responses of wafer-scale graphene, providing a versatile platform for graphene-based optoelectronic and nanophotonic applications.

---

1 Introduction

Epitaxial growth of graphene and few-layer graphene (FLG) on silicon carbide (SiC) offers a uniquely controlled and reproducible pathway for synthesizing high-quality graphene over wafer-scale areas.^1–4 In this approach, high-temperature thermal decomposition of SiC induces the preferential sublimation of silicon atoms, leaving behind a carbon-rich surface that reconstructs into atomically ordered graphene layers. This self-organized growth process yields graphene films with excellent structural uniformity, low defect densities, and long-range crystalline order.³ Beyond structural quality, epitaxial graphene on SiC offers several intrinsic advantages for device integration. The wide-bandgap, semiconducting nature of SiC provides a mechanically and chemically robust substrate compatible with conventional semiconductor processing, enabling direct integration of graphene into scalable electronic and optoelectronic architectures.⁴ Moreover, the electronic properties of epitaxial graphene—including carrier density, Fermi-level position, and interlayer coupling—can be systematically tuned through growth conditions, substrate termination, and post-growth treatments.^1–4 Collectively, these attributes establish epitaxial graphene on SiC as a versatile, wafer-scale materials platform for developing next-generation electronic, photonic, and plasmonic devices.^5–7

Metal intercalation has emerged as a powerful method to tune the properties of epitaxial graphene on SiC.⁸ By inserting metal atoms at the graphene–substrate interface or between adjacent graphene layers, intercalation can electronically decouple graphene from SiC and modify the interlayer coupling of FLG. As a result, carrier density, Fermi-level position, and band structure can be precisely tailored. Depending on the choice of the intercalant, additional functionalities—such as magnetism, superconductivity, and topologically nontrivial phases—can be introduced,^9–16 significantly expanding the capabilities of graphene-based heterostructures for spintronics, nanoelectronics, and emerging quantum technologies.

Among the various intercalants studied, gadolinium (Gd) stands out due to its large magnetic moment and strong interaction between localized 4f states and itinerant electrons. Intercalating Gd beneath epitaxial graphene may induce magnetic ordering and spin polarization while preserving the high mobility of graphene. The presence of Gd also modifies interlayer coupling and graphene doping via charge transfer, creating spatial variations in carrier density and work function.^17,18 These effects provide opportunities to investigate spin-dependent transport, magneto-optical phenomena, and nanoscale electronic inhomogeneity in graphene.

To date, research on Gd-intercalated epitaxial graphene has relied largely on X-ray photoelectron spectroscopy, low-energy electron diffraction, scanning tunneling microscopy (STM), and angle-resolved photoemission spectroscopy,^18–20 which offer valuable insights into atomic configuration and band-structure evolution. However, how Gd intercalation modifies the local optical response and infrared dynamics at the nanoscale remains largely unexplored. A deeper understanding of how Gd intercalation modifies optical responses and infrared dynamics is essential for advancing its use in optoelectronic and photonic devices and for revealing how nanoscale variations in intercalation drive local electronic inhomogeneities.

In this work, we performed systematic nano-infrared (IR) imaging and spectroscopy studies of Gd-intercalated epitaxial bilayer graphene (BLG) on SiC by using a scattering-type scanning near-field optical microscope (s-SNOM). The s-SNOM was built based on a tapping-mode atomic force microscope (AFM), so we can collect both topography and IR images of the sample simultaneously. By coupling the s-SNOM with various IR sources, we performed both nano-IR imaging and spectroscopy on the sample with a spatial resolution of ∼25 nm, which enables visualization and quantitative mapping of nanoscale domains with distinct IR responses. Prior to the nano-IR measurements, we characterized the surface morphology at different stages of growth using STM. By combining nano-IR, AFM, and STM techniques, we obtained a comprehensive understanding of the nanoscale surface features and electronic properties of Gd-intercalated BLG. All nano-IR and AFM measurements were performed at room temperature (293 K). For STM imaging, the large-area (Fig. 1a and b) and high-resolution (Fig. 1c and d) scans were performed at 293 K and 77 K, respectively.

Fig. 1 (a) STM image of epitaxial BLG after deposition of 0.5 ML of Gd. (b) STM image of the same sample after annealing at 850 °C for 20 minutes. The blue arrows mark some big remnant clusters after annealing. The green arrows mark Gd-intercalated domains. Inset plots the line profile extracted along the red dashed line. (c) A zoomed-in STM scan at the boundary between pristine BLG and intercalated domains. (d) A high-resolution STM image on top of a flat Gd-intercalated domain. Scan parameters of STM: (a) 2 V, 0.5 nA; (b) 3 V, 0.5 nA; (c) 0.2 V, 1 nA; (d) −0.3 V, 0.5 nA.

2 Results and discussion

2.1 Sample preparation and STM characterization

The samples that we studied in this work are epitaxial BLG with partial intercalation, so it is convenient to compare in situ the intercalated and non-intercalated sample regions. The initial BLG samples were grown uniformly on the Si face of 6H-SiC(0001) wafers after annealing the wafers at ∼1570 K under ultra-high vacuum conditions (p < 3 × 10⁻¹⁰ mbar). This is a well-established growth method that is introduced in detail in previous works.^18,21,22 The thickness and uniformity of the BLG samples were then verified by STM before Gd intercalation.¹⁸ Metal intercalation is a kinetic process in which surface-deposited atoms overcome energy barriers to transfer from the top surface to the buried interface. Key parameters that define the transfer are the initial deposited amount of metal atoms, the annealing temperature, and the time spent at this temperature. Gd intercalation is very similar to Dy intercalation, with roughly the same annealing temperatures and intercalated phases.¹⁹ In Fig. 1a, we present the STM topography after the deposition of a nominal 0.5 monolayer (ML) of Gd at room temperature, which shows fractal 3-d islands that are non-crystalline.^23,24 The stepper terrain of the sample surface is also clearly seen, which is a signature feature of the topography of the epitaxial graphene grown on the Si-face of SiC.

After deposition, the sample was annealed at 850 °C for 20 minutes. Fig. 1b presents an STM image acquired after annealing, showing that most of the fractal clusters observed in Fig. 1a have significantly decreased in size. Only a small number of large clusters remain (e.g., those marked by blue arrows), and these features are stable for many days even under ambient conditions. In addition, several relatively flat domains are visible in Fig. 1b (e.g., those marked by green arrows). These domains correspond to Gd-intercalated regions, where Gd atoms occupy the interlayer spacing between graphene layers during the annealing process. Near the centers of these domains, additional step-like features are observed, suggesting the formation of a second Gd layer beneath the surface. As shown by the line profile, the thickness of each Gd layer is approximately 0.25 nm.

These distinct topographic features are likely the result of a defect-assisted intercalation mechanism.^25,26 In this process, Gd atoms initially penetrate the graphene layers through vacancy defects and subsequently spread laterally to form the first intercalated layer. Continued incorporation of Gd atoms through the same defect can lead to the formation of an additional intercalated layer. Fig. 1(c) shows a zoomed-in STM image at the boundary of an intercalated domain, where pristine bilayer graphene (BLG) exhibiting the characteristic 6 × 6 moiré pattern (center right) is surrounded by intercalated regions with increased surface roughness. Similar intercalation-induced surface morphologies have been reported previously.^17,18 Finally, Fig. 1d presents a higher-resolution STM image acquired on a flat intercalated domain, in which the graphene lattice is clearly resolved, confirming that the Gd atoms reside beneath the top graphene layer.

2.2 Nano-IR imaging

Following the STM characterizations, we performed s-SNOM imaging studies of the same sample at ambient conditions. The measurement setup is illustrated in Fig. 2a, where the s-SNOM tip is excited by a continuous-wave CO₂ laser that covers discrete wavelengths from 9.2 to 10.7 μm. Fig. 2 plots the AFM topography (Fig. 2b), AFM phase (Fig. 2c), and nano-IR imaging data (Fig. 2d) simultaneously taken on a typical sample area. The AFM phase plotted in Fig. 2c reflects the mechanical stiffness of the sample.²⁷ The nano-IR image was taken at an excitation laser wavelength (λ₀) of 9.4 μm. As discussed in the following paragraphs, this wavelength is away from the strong phonon resonance of SiC centered at 11.5 μm, so the IR responses are mainly from the graphene layers and the Gd-intercalated domains. The near-field IR amplitude (s) measured by s-SNOM is roughly proportional to the average reflection coefficient of all evanescent fields.²⁸ Regions with higher optical conductivity typically exhibit greater reflection coefficients, thus yielding a stronger IR amplitude signal.

Fig. 2 (a) Schematic of s-SNOM study of Gd-intercalated epitaxial BLG sample. (b) and (c) The AFM topography and phase images of a typical sample surface of our Gd-intercalated BLG sample. (d) The nano-IR amplitude image of the sample taken simultaneously with the AFM images. The blue and green arrows in panels (b–d) indicate Gd clusters and Gd-intercalated domains, respectively.

From the AFM topography image (Fig. 2b), one can see again the stepper terrain as well as the big Gd clusters, like those on the STM image (Fig. 1b). The Gd clusters are also seen in the AFM phase image (Fig. 2c), where the majority part appears to be brownish, corresponding to “near-pristine” BLG not affected by the Gd intercalation process. For convenience, we set the phase value (ψ) of the majority region to be zero in Fig. 2c. The Gd clusters (e.g., a pair of such particles marked with blue arrows) appear to have a very high AFM phase (ψ ∼ 0.3 rad.), indicating they are stiffer than the relatively softer pristine BLG. From the nano-IR image (Fig. 2d), we observe that these clusters appear to be darker (i.e., lower IR amplitude), signaling that they are less conductive than BLG. Therefore, Gd clusters are likely not fully crystalline and hence not very metallic.

Another type of feature in the AFM phase image (Fig. 2c) corresponds to patches of domains with sizes varying from tens to hundreds of nanometers. For convenience, we mark a typical one with green arrows in all data images in Fig. 2. These domains have an intermediate phase signal of ψ ∼ 0.1 rad. Indicating they are slightly stiffer than graphene layers. They also appear to be slightly higher (0.3–0.5 nm) than adjacent BLG regions from the AFM topography image (seen more clearly in the zoomed-in images in Fig. 3). Based on these characteristics, we conclude that these domains are Gd-intercalated regions like those seen previously in the STM images (see Fig. 1b and c). From the nano-IR image (Fig. 2d), one can see that these domains appear to have higher IR signals than the rest of the sample areas, so they are more conductive. In addition, we found that the IR contrast of these bright domains varies dramatically from one to another (see Fig. 2d), indicating that they exhibit different conductivities. The variation of IR signals and conductivities is likely due to the differences in the amount of intercalated Gd atoms within each domain.

Fig. 3 (a and b) Zoomed-in AFM topography and AFM phase images of the sample region marked by the white dashed square in Fig. 2b. (c–f) Nano-IR imaging data of the sample when varying the incident light from 9.2 to 10.7 μm. (g) Line profiles extracted along the dashed lines in panels c–f. The profiles are displaced vertically for clarity.

To better visualize all these features, we performed higher-resolution imaging of the sample area marked with a dashed square in Fig. 2b at various excitation wavelengths (λ₀). The imaging results are shown in Fig. 3, where we plot the AFM topography (Fig. 3a), AFM phase (Fig. 3b), and nano-IR images (Fig. 3c–f). Here, one can see more clearly Gd-intercalated domains in the AFM and nano-IR images. Again, we pay attention to the large domain at the lower part of the images (marked with a green arrow). Upon close inspection of the nano-IR images, we find that the distribution of the IR signals within the domain is non-uniform. The center of the domain appears to be brighter than the edges. This might be due to a denser distribution of Gd atoms at the center of the domain, where the graphene properties are more dramatically altered.

Now we examine the λ₀-dependence of IR images (Fig. 3c–f), where one can see that the domain contrast clearly changes with varying IR wavelengths. For quantitative analysis, we plot in Fig. 3g the line profiles extracted across the domain (along the dashed lines in Fig. 3c–f). These profiles are normalized to the background signal of non-intercalated BLG. From Fig. 3g, one can see that the maximum IR contrast increases from ∼70% to ∼130% from 9.2 to 9.6 μm and then drops to about 40% at 10.7 μm. The change in IR contrast among different IR wavelengths, as discussed in detail below, is due to the optical phonon on SiC that peaks around 11.5 μm.

2.3 Nano-IR spectroscopy

A more comprehensive view of the λ₀-dependent IR contrast can be seen in broadband nano-IR spectroscopy data of Gd-intercalated BLG. To perform nano-IR spectroscopy, we excited the s-SNOM with a broadband IR pulsed laser covering 7–13 μm. The IR pulses are generated through a difference frequency generation setup. In Fig. 4a and b, we show the AFM topography and AFM phase image of a sample area where multiple Gd-intercalated domains are clearly seen in the phase image. We then perform a line scan across a random Gd-intercalated domain (along the blue dashed line in Fig. 4a and b). At every pixel along the line, we collect a broadband spectrum. The full linescan data are presented in Fig. 4c and d as two hyperspectral images covering the wavelength ranges of 10.5–12.5 µm (Fig. 4c) and 7.5–10.5 µm (Fig. 4d), respectively. In both images, the horizontal axis represents the tip position, while the vertical axis corresponds to the IR wavelength. The two wavelength regions are shown separately because of the large contrast in IR signal intensity (note the different color scale settings), which arises from the strong SiC phonon resonance near 11.5 µm.

Fig. 4 (a and b) AFM topography and AFM phase images of a typical sample area of Gd-intercalated BLG. (c and d) Broadband hyperspectral maps s(x, ω) covering 10.5–12.5 µm and 7.5–10.5 µm, respectively, which were taken by line-scanning the s-SNOM probe across a bright domain. (e and f) Line profiles across the bright domain at various wavelengths extracted from the hyperspectral maps in panels c and d. The red dashed lines are a guide to the eye.

From Fig. 4c and d, one can spot the Gd-intercalated domain (centered at around x = 800 nm) due to its different IR contrast. At the phonon resonance (Fig. 4c), the domain region appears to be darker compared to neighboring pristine BLG regions. Away from the phonon resonance (Fig. 4d), on the other hand, the domain is brighter. The λ₀-dependent IR contrast of the domain can be visualized more clearly in the horizontal line profiles (Fig. 4e and f) extracted directly from the hyperspectral images. These profiles serve the same purpose as those of Fig. 3g but were taken with a different approach and cover a broader wavelength region. From both the hyperspectral images and profiles, one can see that the Gd-intercalated domain has a positive IR contrast relative to pristine BLG off the SiC phonon (e.g., λ₀ < 10.5 μm), which is consistent with the results shown in Fig. 3g. Nevertheless, when it is on top of the phonon (e.g., λ₀ > 11 μm), the domain shows a negative IR contrast.

The unique IR contrasts can be seen more clearly in Fig. 5a, where we plot together the nano-IR spectra of the Gd-intercalated BLG domain, pristine BLG region, and a bare SiC substrate. As shown in the spectra (Fig. 5a), the SiC phonon resonance, which peaked at 11.5 μm, is the dominant feature within the spectroscopic window. At shorter wavelengths below the phonon resonance, Gd-intercalated BLG has higher signals than pristine BLG and the bare SiC substrate. On top of the phonon resonance, in contrast, the signal of SiC is higher than that of the pristine and intercalated BLG. The distinct IR contrasts on and off the phonon resonance both originate from the high conductivity of the two types of BLG. Away from the phonon, the graphene layers are the main contributors to the s-SNOM signal, so higher conductivity leads to a higher IR amplitude signal. Atop the phonon, SiC is the dominant source of the s-SNOM signal, and graphene layers quench the SiC signal due to the electronic screening effect. Enhanced conductivity of the graphene layers due to Gd intercalation results in more effective screening and consequently a weaker phonon resonance.

Fig. 5 (a and b) Experimental and modeled nano-IR spectra of the BLG (red), Gd-intercalated BLG (blue), and bare SiC substrate (black). Here, E_F of BLG is set to be 0.34 eV, and the conductivity enhancement factor η of Gd-intercalated BLG is set to be 2.4. (c) Modeled nano-IR spectra with various settings of η. (d) Modeled nano-IR spectra of a 2D sample with different settings of charge scattering rate.

2.4 Quantitative near-field modeling

The spectroscopic features shown in the experimental spectra (Fig. 5a) are well captured by our modeling results (Fig. 5b). Here, the computation was done by using a numerical model built on Comsol Multiphysics. In the model, the s-SNOM tip is approximated as an elongated conducting spheroid with a radius of curvature of 25 nm at the tip apex. The near-field signal is computed by evaluating the total radiating dipole moments of the tip. The graphene sample (both BLG and Gd-intercalated BLG) is modeled as a thin conducting film with an effective thickness of 1 nm. While this setting is close to the real thickness of BLG (2 graphene layers plus a buffer layer) and slightly smaller than that of Gd-intercalated BLG, the computation results do not change significantly (<2% difference) if we choose an effective thickness of 2 or 3 nm. A detailed introduction of the model can be found in a previous paper.²⁹

The key modeling parameter of the BLG sample (with or without Gd-intercalation) is the optical conductivity σ(E), where E is the IR energy. Note that the Fermi energy (E_F) of epitaxial BLG is high (>0.3 eV) (see ref. 10), the interband transitions (only enabled at energies above 2E_F) are far above our energy range (0.1–0.15 eV). Therefore, we only consider the intraband transitions when formulating the optical conductivity. For pristine BLG, the Drude conductivity can be written as:³⁰

Here, γ ≈ 0.39 eV is the interlayer hopping energy of BLG, and Γ is the charge scattering rate. By modeling the experimental spectra, we estimate the E_F of BLG to be 0.34 ± 0.06 eV, which is consistent with previous ARPES results.¹⁸

Based on the BLG modeling result, we can now evaluate the conductivity enhancement by Gd intercalation. Here, we set the optical conductivity of Gd-intercalated BLG in the model to be

σ_Gd-BLG = ησ_BLG, (2)

where η is the enhancement factor and σ_BLG is the conductivity of pristine BLG at E_F = 0.34 eV (see previous paragraph). In Fig. 5b, the modeling spectrum of Gd-intercalated BLG (blue curve) is calculated with η = 2.4, which shows a reasonable agreement with the experimental spectrum in Fig. 5a. Results with other settings of η are shown in Fig. 5c, which indicates that the resonance spectrum is sensitive to Drude conductivity. By comparison, the nano-IR spectra close to the phonon resonance of SiC show only minor variations with respect to the parameter Γ (see Fig. 5d). Note that the enhancement factor η = 2.4 is measured from one particular domain. Considering that there is a large distribution of domains with various sizes and IR signals (see Fig. 2 and 3), the values of η could vary from domain to domain.

The sizable enhancement of conductivity is a direct signature of the significant change of electronic properties of BLG due to Gd intercalation. According to the previous study,²¹ the Gd atoms are likely intercalated between the top two graphene layers when annealed at a relatively low temperature (≤1000 °C). In other words, BLG is decoupled and becomes two isolated SLG due to Gd intercalation (see Fig. 2a). In addition, the intercalated Gd atoms can also significantly dope the adjacent graphene layers due to charge transfer.^21,31 Therefore, the conductivity enhancement observed in our work is likely due to the contribution from both the decoupling and doping of the top two graphene layers by Gd.

For example, in the case of η = 2.4, the two decoupled graphene layers contribute to a total of 2.4 times of σ_BLG(E_F = 0.34 eV). Assuming that the two SLG layers are equally doped by Gd atoms, they each contribute 1.2 times of σ_BLG. Similar to BLG, the interband transitions only occur at very high energies (≥2E_F), so we only need to consider the Drude component of the optical conductivity, as written below:

Based on eqn (3), we estimate that the E_F of each SLG layer is about 0.56 eV, indicating that the graphene layers are highly doped. Note that the estimated value of E_F is only for the measured domain. The spatial variation of IR signals across different domains seen in Fig. 2d indicates that the E_F of decoupled graphene layers also varies. According to previous DFT calculations,²¹ E_F of the decoupled graphene layers due to Gd intercalation can reach as high as 1.6 eV when the Gd atoms are densely distributed between the two top graphene layers.

3 Conclusions

In summary, we have performed a comprehensive nano-IR imaging and spectroscopy study of epitaxial BLG partially intercalated by Gd atoms. We found through nano-IR imaging a large density of Gd-intercalated graphene nano-domains with different sizes and IR contrasts. These domains enhance the IR signals away from the SiC phonon resonance and quench the signal atop the resonance. Both types of responses originate from the increased conductivity of BLG by Gd intercalation. With quantitative modeling of the nano-IR spectra, we estimate the optical conductivity of Gd-intercalated nano-domains to be as high as 2.4 times that of the pristine BLG regions. The enhancement of conductivity is attributed to the decoupling and extra doping of the graphene layers by the intercalated Gd atoms. Our work demonstrates an effective and convenient method for fabricating graphene-based conducting nanodomains, which can be controlled by varying the amount of metal deposition, annealing temperature, and annealing time.²¹ The spatial distribution of nanodomains across the wafer can be further manipulated by introducing shadow masks during metal deposition.^32,33 Our findings highlight the potential of controlled rare-earth intercalation to locally tune the electronic and optical properties of graphene, paving the way for a range of advanced applications such as spatially programmable metasurfaces, reconfigurable photonic circuits, and tunable optoelectronic devices.³⁴

Conflicts of interest

The authors have no conflicts to disclose.

Data availability

The data that support the findings of the study are available from the corresponding author upon reasonable request.

Acknowledgements

Work done at Ames Lab was supported by the U.S. Department of Energy, Office of Basic Energy Science, Division of Materials Sciences and Engineering. Ames Laboratory is operated for the U.S. Department of Energy by Iowa State University under Contract No. DE-AC02-07CH11358. The experimental apparatus and numerical modeling method were partially supported by the National Science Foundation under Grant No. DMR1945560.

References

1. W. Norimatsu and M. Kusunoki, Phys. Chem. Chem. Phys., 2014, 16, 3501–3511.
2. J. Hass, W. A. de Heer and E. H. Conrad, J. Phys.: Condens. Matter, 2008, 20, 323202.
3. G. R. Yazdi, T. Lakimov and R. Yakimova, Crystal, 2016, 6, 53.
4. N. Mishra, J. Boeckl, N. Motta and F. Iacopi, Phys. Status Solidi A, 2016, 9, 2277–2289.
5. J. Zhao, P. Ji, Y. Li, R. Li, K. Zhang, H. Tian, K. Yu, B. Bian, L. Hao, X. Xiao, W. Griffin, N. Dudeck, R. Moro, L. Ma and W. A. de Heer, Nature, 2024, 625, 60–65.
6. N. N. Nguyen, H. Lee, H. C. Lee and K. Cho, Acc. Chem. Res., 2022, 55, 673–684.
7. X. Wang, C. Tang, R. Li, X. Zhu, W. Jin, L. Hao, R. Sun, J. Liu, Y. Ma and L. Ma, ACS Appl. Nano Mater., 2025, 8, 15799–15807.
8. D. Yang, F. Ma, X. Bian, Q. Xia, K. Xu and T. Hu, J. Phys.: Condens. Matter, 2024, 36, 173003.
9. N. A. Anderson, M. Hupalo, D. Keavney, M. Tringides and D. Vaknin, J. Magn. Magn. Mater., 2019, 474, 666–670.
10. T. P. Kaloni, M. Upadhyay Kahaly and U. Schwingenschlögl, J. Mater. Chem., 2011, 21, 18681–18685.
11. R. Decker, J. Brede, N. Atodiresei, V. Caciuc, S. Blügel and R. Wiesendanger, Phys. Rev. B: Condens. Matter Mater. Phys., 2013, 87, 041403(R).
12. A. P. Durajski, K. M. Skoczylas and R. Szczȩśniak, Phys. Chem. Chem. Phys., 2019, 21, 5925–5931.
13. T. Huempfner, F. Otto, R. Forker, P. Müller and T. Fritz, Adv. Mater. Interfaces, 2023, 10, 2300014.
14. S. Ichinokura, K. Sugawara, A. Takayama, T. Takahashi and S. Hasegawa, ACS Nano, 2016, 10, 2761–2765.
15. Y. Li, P. Tang, P. Chen, J. Wu, B.-L. Gu, Y. Fang, S. B. Zhang and W. Duan, Phys. Rev. B: Condens. Matter Mater. Phys., 2013, 87, 245127.
16. J. Zhou and P. Jena, Phys. Rev. B, 2017, 95, 081102(R).
17. M. Kim, M. Hupalo, M. C. Tringides, B. Schrunk, A. Kaminski, K.-M. Ho and C.-Z. Wang, J. Phys. Chem. C, 2020, 124, 28132–28138.
18. M. Kolmer, B. Schrunk, M. Hupalo, J. Hall, S. Chen, J. Zhang, C.-Z. Wang, A. Kaminski and M. C. Tringides, J. Phys. Chem. C, 2022, 126, 6863–6873.
19. S. Chen, Y. Han, M. Kolmer, J. Hall, M. Hupalo, J. W. Evans and M. C. Tringides, Phys. Rev. B, 2023, 107, 045408.
20. S. A. Gorovikov, A. M. Shikin, G. V. Prudnikova, V. K. Adamchuk, S. L. Molodtsov, C. Laubschat and A. M. Ionov, Surf. Sci., 2001, 474, 98–106.
21. M. Hupalo, E. H. Conrad and M. C. Tringides, Phys. Rev. B, 2009, 80, 041401(R).
22. I. Forbeaux, J. M. Themlin and J. M. Debever, Phys. Rev. B: Condens. Matter Mater. Phys., 1998, 58, 16396–16406.
23. X. Liu, C.-Z. Wang, M. Hupalo, Y.-X. Yao, M. C. Tringides, W.-C. Lu and K.-M. Ho, Phys. Rev. B: Condens. Matter Mater. Phys., 2010, 82, 245408.
24. M. Hupalo, S. Binz and M. C. Tringides, J. Phys.: Condens. Matter, 2011, 23, 045005.
25. Y. Liu, X. Liu, C.-Z. Wang, Y. Han, J. W. Evans, A. Lii-Rosales, M. C. Tringides and P. A. Thiel, J. Phys. Chem. C, 2021, 125, 6954–6962.
26. N. Nayir, M. Y. Sengul, A. L. Costine, P. Reinke, S. Rajabpour, A. Bansal, A. Kozhakhmetov, J. Robinson, J. M. Redwing and A. van Duin, Carbon, 2022, 190, 276–290.
27. S. N. Magonov, V. Elings and M.-H. Whangbo, Surf. Sci., 1997, 375, L385–L391.
28. Z. Fei, G. O. Andreev, W. Bao, L. M. Zhang, A. S. McLeod, C. Wang, M. K. Stewart, Z. Zhao, G. Dominguez, M. Thiemens, M. M. Fogler, M. J. Tauber, A. H. Castro-Neto, C. N. Lau, F. Keilmann and D. N. Basov, Nano Lett., 2011, 11, 4701–4705.
29. Y. Luan, L. McDermott, F. Hu and Z. Fei, Phys. Rev. Appl., 2020, 13, 034020.
30. E. J. Nicol and J. P. Carbotte, Phys. Rev. B: Condens. Matter Mater. Phys., 2008, 77, 155409.
31. M. Alemani, A. Barfuss, B. Geng, C. Girit, P. Reisenauer, M. F. Crommie, F. Wang, A. Zettl and F. Hellman, Phys. Rev. B: Condens. Matter Mater. Phys., 2012, 86, 075433.
32. M. Apanius, P. B. Kaul and A. R. Abramson, Sens. Actuators, A, 2007, 140, 168–175.
33. C. G. Courcimault and M. G. Allen, Nanotechnology, 2004, 15, S528–S533.
34. P. A. Huidobro, M. Kraft, S. A. Maier and J. B. Pendry, ACS Nano, 2016, 10, 5499–5506.

---