Quantitative Subsurface Atomic Structure Fingerprint for 2D Materials

Jun 5, 2016 - Interfaces and subsurface layers are critical for the performance of devices made of 2D materials and heterostructures. Facile, nondestr...
2 downloads 0 Views 4MB Size
Download PDF
Quantitative Subsurface Atomic Structure Fingerprint for 2D Materials and Heterostructures by First-Principles-Calibrated Contact-Resonance Atomic Force Microscopy Qing Tu,†,§ Björn Lange,†,§ Zehra Parlak,† Joao Marcelo J. Lopes,‡ Volker Blum,*,† and Stefan Zauscher*,† †

Department of Mechanical Engineering and Materials Science, Duke University, Durham, North Carolina 27708, United States Paul-Drude-Institut für Festkörperelektronik, Hausvogteiplatz 5-7, D-10117 Berlin, Germany

‡

S Supporting Information *

ABSTRACT: Interfaces and subsurface layers are critical for the performance of devices made of 2D materials and heterostructures. Facile, nondestructive, and quantitative ways to characterize the structure of atomically thin, layered materials are thus essential to ensure control of the resultant properties. Here, we show that contact-resonance atomic force microscopywhich is exquisitely sensitive to stiffness changes that arise from even a single atomic layer of a van der Waals-adhered materialis a powerful experimental tool to address this challenge. A combined density functional theory and continuum modeling approach is introduced that yields sub-surface-sensitive, nanomechanical fingerprints associated with specific, well-defined structure models of individual surface domains. Where such models are known, this information can be correlated with experimentally obtained contact-resonance frequency maps to reveal the (sub)surface structure of different domains on the sample. KEYWORDS: ab initio calculations, contact-resonance atomic force microscopy, elastic properties, surfaces and interfaces, 2D materials and heterostructures ronments (e.g., low-energy electron diffraction (LEED),5,6,19−21 low-energy electron microscopy (LEEM),6,19−22 grazing incidence X-ray diffraction,23 X-ray standing waves13), or provide only area integrated information (e.g., X-ray photoelectron spectroscopy5,20 and Raman spectroscopy15,16,20). Several atomic force microscopy (AFM)-based nondestructive techniques have been developed to indirectly probe interfaces and/or subsurface structures. For instance, scattering-type scanning near-field optical microscopy24 provides insights into lateral defect structures by measuring the plasmonic response of a surface layer. However, this technique demands the existence of high momenta plasmons and is therefore not suitable for a large range of material. Tip-enhanced Raman spectroscopy (TERS) can provide nanoscale chemical analysis of a sample surface,25−27 but the stringent prerequisites for the tips and sample substrates limit the application of TERS.27 Peakforce tapping AFM has been used to detect the commensurate−

T

wo-dimensional van der Waals (vdW)-bonded materials, including graphene, hexagonal boron nitride (hBN), and transition metal dichalcogenides (MX2), possess extraordinary electrical, optical, and mechanical properties, holding great promise for high-performance device applications.1 Harnessing the unique properties of 2D materials and heterostructures for nanoscale device applications requires precise control over the interface to the growth substrate,2−6 the layer number,7−9 the lateral arrangement of different domains,10,11 and the nature of subsurface structure elements.12 Interfaces between these 2D materials are critical for device operation,1,12,13 for example, in vertical heterostructures, and interactions through interfaces and/or intentional or unintentional intercalations6 can chemically change or mechanically distort the 2D material layers, altering their electrical,13,14 optical,14−16 and mechanical17 properties. It is thus critical to obtain precise information about the identity and the in- and out-of-plane arrangement of material domains in layered 2D material heterostructures.18 Many existing surface and interface characterization techniques are destructive (e.g., cross-sectional transmission electron microscopy19), require restrictive experimental envi© 2016 American Chemical Society

Received: April 8, 2016 Accepted: June 4, 2016 Published: June 5, 2016 6491

DOI: 10.1021/acsnano.6b02402 ACS Nano 2016, 10, 6491−6500

Article

www.acsnano.org

Article

ACS Nano

Figure 1. Overview of our first-principles-informed CR-AFM approach. In CR-AFM (a), the sample is placed on an ultrasonic actuator, and the AFM cantilever, in contact with the sample, is driven and maintained at resonance by dual AC resonance tracking while imaging the sample.31 Regions on the sample with different mechanical properties (i.e., I and II) have different CR frequencies ( f I and f II) (b). In addition to the resonance frequency, the applied forces (F⊥I and F⊥II) while imaging the different regions are recorded. The cantilever dynamics is modeled by an Euler−Bernoulli beam (c) with a linear spring k* located at the tip position (L1). The tip−sample contact is approximated by a Hertzian contact model (d), i.e., an elastic sphere with radius Rt in contact with an elastic half-plane. With the experimentally measured frequencies and forces, the calibrated relative tip positions (L1/L), and invoking a Hertzian contact model, the Euler−Bernoulli beam equation is used to derive the reduced contact modulus ratio (EII*/EI*)exp. Alternatively, the AFM tip can also be considered as an ultrasonic transducer with a size equal to the contact area (radius a), and the generated acoustic waves (dark arrows) can propagate into the layered materials and across the interfaces (transmission and reflection (light arrows)) (d). This situation can be described by a continuum model, for which the stiffness matrix for each material layer (Ckl,i(I,II)) in each of the different regions on the sample is derived from an atomistic structure model (e) via first-principles calculations (f). Using the experimentally measured frequencies, forces, and the AFM tip radius, the continuum model will generate a reduced modulus ratio (E*II /E*I )th, characteristic of the atomistic structure models used. By comparing the theoretically derived with the experimentally derived reduced modulus ratio, one can then test the validity of a particular structure model, revealing detailed subsurface/interface information.

istic calculations and a continuum mechanics model to provide quantitative, atomic-layer-resolved nanomechanical information. This approach yields a nanomechanical, sub-surfacesensitive “fingerprint” of the atomic structure of an interface, demonstrated here for heterostructures involving graphene and other 2D materials. Fast continuum simulations, informed by databases30 of elastic properties of 2D and 3D materials, will enable the rapid identification of experimentally observed structure fingerprints.

incommensurate transition in graphene−h-BN heterostructures by analyzing qualitative out-of-plane Young’s modulus contrast patterns.17 Recently, contact-resonance atomic force microscopy (CR-AFM) in liquid has been used to characterize the out-of-plane mechanical properties of 2D Ti3C2Tx (MXene)based electrodes in situ during alkaline cation intercalation and extraction.28 In this case, the film studied was micrometers thick and the measured mechanical property changes reflect the cumulative effects of the interfaces in the film. Here, we report on an innovative approach that synergistically combines CR-AFM29a nondestructive, dynamic scanning probe microscopy method that can be employed in a wide range of environmental conditionswith first-principles atom-

RESULTS AND DISCUSSION First-Principles-Informed CR-AFM. Figure 1 summarizes our combined approach. CR-AFM is based on conventional 6492

DOI: 10.1021/acsnano.6b02402 ACS Nano 2016, 10, 6491−6500

Article

ACS Nano

Figure 2. CR-AFM exemplified for epitaxial graphene grown on 6H-SiC(0001). (a) Atomistic models of (left to right) ZLG, MLG, and BLG. The blue carbon layer, G0, indicates the buffer layer. The first (G1) and second (G2) vdW-bonded graphene layers are shown in light green and red, respectively. Silicon atoms at the interface are shown in brown and carbon atoms in gray. (b,c) CR-AFM height image and corresponding first flexural mode contact-resonant frequency mapping over a terrace on the ZLG sample. (d,e) Raman spectra of the ZLG and the MLG samples, respectively. (f) 2D peak of the MLG sample, fitted by a single Lorentzian peak (fwhm ∼45 ± 2 cm−1). (g) AFM height image and (h) corresponding CR-AFM first flexural mode contact-resonance frequency mapping on the MLG sample. The blue arrows indicate MLG growth fingers extending into the BLG region of the sample. (i) High-resolution CR-AFM frequency map of the region shown in red in (h). Other types of fingers (referred to as region “O” finger) with much higher CR frequency than that of the MLG domains are identified, and one of them is circled by a white dashed line. (j) Average height cross-sectional profiles (direction: from left to right) for the regions indicated in (i). For averaged f and F⊥ of the domains identified in (c) and (i), see SI Tables S1 and S2.

k* can be related to the cantilever flexural resonance frequencies, f, the primary observable of CR-AFM, by solving the beam dynamics analytically (see Supporting Information (SI), Section I). For a given flexural vibration mode, a higher contact stiffness will result in a higher resonance frequency.32 The contact stiffness thus reflects the elastic properties of (i) each atomic layer in the sample, of (ii) each interface, and of (iii) the probe tip. The Hertzian contact model (Figure 1c), widely accepted to describe the AFM tip−sample contact for elastic materials,32−34 connects k* to the reduced elastic modulus of the tip−sample contact, E*, via

contact-mode AFM but uses an additional actuatortypically mounted beneath the samplethat drives the coupled cantilever−sample system (Figure 1a) at resonance. CR-AFM is thus able to map structure-dependent differences via the contact-resonance frequency f of different regions (Figure 1b), reflecting nanoscale changes in the elastic and interfacial properties of thin, layered materials systems.29,32 We show that these measured changes can be quantitatively compared to the theoretically expected response from domains of different structure. Remarkably, our approach is sensitive to changes of the structure of atomically thin layers and their sequence in the investigated region. In CR-AFM, the AFM cantilever can be considered as a Euler−Bernoulli cantilever beam (Figure 1c). The tip−sample contact is represented by an elastic spring,29 whose spring constant reflects the effective tip−sample contact stiffness k*

k* =

dF dδ

F⊥

k* =

3

6F ⊥R tE*2

(2)

(see Methods for details). Rt is the effective radius of the assumed spherical tip in contact with a flat surface. In practice, both Rt and F⊥ are determined experimentally along with k*. Importantly, eq 2 shows that the dependence on Rt can be eliminated by focusing on ratios of the reduced elastic moduli, EII*/EI*, of two different surface regions I and II:

(1)

where δ is the indentation depth of the probe and F⊥ is the total force applied normal to the surface. F⊥ includes both the set point force applied through the AFM cantilever and the local interaction between the AFM cantilever tip and the sample surface determined by force−distance “pull-off” curves.

E II* = E I* 6493

⎛ k * ⎞3 F ⊥ ⎜ II ⎟ I⊥ ⎝ kI* ⎠ FII

(3) DOI: 10.1021/acsnano.6b02402 ACS Nano 2016, 10, 6491−6500

Article

ACS Nano

remnants left behind during the anisotropic growth of a new graphene layer. Thus, G0-terminated “fingers” should extend into newly formed MLG domains, and G1-terminated fingers should extend into newly formed BLG domains. Indeed, Figure 2h (blue arrows) reveals some features that are consistent with MLG fingers extending into BLG domains because these fingers have contact-resonance frequencies similar to those of MLG domains. However, CR-AFM reveals another group of “fingers” (labeled “O” in Figure 2i) that also extends into the BLG area. These fingers exhibit CR frequencies considerably higher than those in the BLG and MLG regions (SI Table S2). AFM height measurements show that, on this sample, their height difference to BLG is smaller than the step height between BLG and MLG (Figure 2j). Additionally, after storing the sample in air for several months, regions “O” show height changes that are not observed for MLG and BLG (SI Figure S4). Thus, the fingers labeled “O” cannot be MLG. The experimentally derived reduced elastic moduli ratio (eq 3) for the ZLG to MLG regions on the ZLG sample in Figure 2c is (EZLG * /EMLG * )exp = 1.082 ± 0.024. For the MLG sample (Figure 2i), (E*BLG/E*MLG)exp = 0.928 ± 0.012 and (E*O/E*MLG)exp = 1.140 ± 0.028 for the BLG to MLG and the “O” to MLG regions, respectively. AFM tapping-mode phase images (SI Figure S2b,d) also provide layer-number-dependent contrast and qualitatively corroborate the observed stiffness differences.19,40 However, in contrast to the CR-AFM procedure summarized in Figure 1, it is not possible to deconvolute the phase contrast into quantitative, mechanical property information on the single atomic layer level.19,33,41 Theoretical Derivation of Reduced Moduli Ratios. To deconvolute the experimentally measured, reduced elastic moduli and, where possible, to unequivocally assign a layer structure with single-atomic-layer resolution, we calculate elastic modulus ratios (EII*/EI*)th by a continuum mechanics modeling approach.34 This approach combines Hertzian contact theory with the mechanical impedance of thin, layered materials derived from acoustic wave propagation. The AFM tip−sample contact (Figure 1c) is considered to be a finitesized transducer (contact radius a) that modulates the layered material surface with ultrasonic pressure fields. When a ≪ λ0, where λ0 is the shortest wavelength in the probed layered solid at the actuator frequency, the radiation loss is very small. Thus, the imaginary part of the effective mechanical radiation impedance Im(Zeff), which reflects the elastic system response, dominates. Im(Zeff) is related to a lossless spring element, equivalent to the sample surface stiffness ks*, which reflects the elastic properties of the layered sample (see Methods and SI Section VII):34

The central point of this paper is that, by our approach, the reduced elastic moduli, measured for different surface regions on a given sample, can be related quantitatively to the specific surface and subsurface atomic structure(s) of these regions. As shown in Figure 1, E*II /E*I can be inferred either purely from experimental data (via the Euler−Bernoulli beam model) or independently from theory by postulating surface and subsurface atomic structure models for different regions (Figure 1e). The elastic properties of each atomic layer are calculated from first principles and passed into a continuum mechanical model (see Methods) that predicts the reduced elastic modulus for the AFM tip−sample system. The measured ratios (EII*/EI*)exp can be compared to first-principles-based ratios (EII*/EI*)th to establish the consistency between experiment and the postulated structure models. The result is a quantitative, nanomechanical “fingerprint” of the surface and subsurface atomic structure of the different regions on a sample. CR-AFM of Epitaxial Graphene on SiC. Epitaxial graphene films, grown by high-temperature Si sublimation on SiC,3,4,20,23,35 are excellent models to demonstrate the power of our approach. The structure and bonding of each plane are well understood in terms of strain-free (6√3 × 6√3)-R30° interface models (Figure 2a)21,35,36 and also amenable to first-principles calculations of elastic properties. The partially σbonded, honeycomb c plane, G0, is referred to as the “buffer layer” (not yet graphene). We call structure models (including the SiC substrate) terminated by the buffer layer G0 “zero-layer graphene” (ZLG) (Figure 2a, left). Monolayer graphene (MLG) (Figure 2a, middle) and bilayer graphene (BLG) (Figure 2a, right) structure models are defined by adding vdWbonded graphene layers G1 and G2, respectively, to the SiC substrate and the “buffer layer” G0. Figure 2b−j shows AFM images and Raman spectra from two different samples grown at conditions chosen to yield predominantly ZLG and MLG domains (see Methods). Their surfaces consist of micrometersized terraces arising from the underlying SiC (Figure 2g and SI Figure S2). Raman spectroscopy confirms the predominant graphene layer number on the samples. The “ZLG sample” (Figure 2d) shows two broad Raman peaks around 1330 and 1575 cm−1 and one smaller peak around 1480 cm−1. These peaks and the absence of a 2D peak are consistent with ZLG as the predominant domain.15 CRAFM imaging additionally shows some MLG, softer than ZLG, along the terrace edges (Figure 2c and discussion below). The “MLG sample” (Figure 2e) shows a sharp Raman G peak, a broad D peak, and a 2D peak that can be fitted by a single Lorentzian distribution function (fwhm ∼45 ± 2 cm−1, Figure 2f), without apparent shoulders. These features are consistent with G1 coverage.4,16 Some bilayer graphene occurs at terrace edges (Figure 2h), in agreement with literature reports.19,20 Figure 2c,h,i shows that CR-AFM imaging yields an astonishingly clear, mechanical contrast between different regions on the samples. On the ZLG sample (Figure 2c), the CR frequencies for the first flexural mode in the region labeled “ZLG” are higher than those in the region labeled “MLG”. On the MLG sample (Figure 2h,i), the region labeled “MLG” has CR frequencies higher than those in the region labeled “BLG”. We show below that these assignments agree quantitatively with results derived theoretically for the corresponding structure models. In addition to large MLG and BLG regions, some finger-like features and small islands are visible on the MLG sample (Figure 2h). In earlier work,38,39 such features were explained as

ks* = Im(2fπ 2a 2Zeff )

(4)

where f is the frequency of the acoustic wave, that is, the CR frequency. Importantly, the continuum model requires as input the atomic-layer-resolved elastic constants of the assumed atomic structure of a given sample region. We calculate in-plane and out-of-plane stiffness tensor elements, cij, from atomistic layer models using density functional theory (DFT) based on the Tkatchenko-Scheffler42 (TS) vdW-corrected Perdew−Burke− Ernzerhof43 generalized gradient approximation (PBE+TS) and linear response theory44 (SI Sections III and IV and Tables S5− S7). The in-plane cij for hexagonal planes (c11, c12, c22, c13, c44, c55, and c66) are calculated for equivalent planes in the vdW6494

DOI: 10.1021/acsnano.6b02402 ACS Nano 2016, 10, 6491−6500

Article

ACS Nano

Figure 3. Calculation of the elastic properties of graphene on SiC(0001). (a) Schematic illustration of the linear spring model. Layers with layer-dependent stiffness tensor elements c33,i are connected by individual springs with spring constant ki. (b) Layer-dependent stress−strain curves for a (6√3 × 6√3)-R30° BLG interface model are shown in the inset. Sn denotes SiC bilayers; G0, G1, G2, and G3 refer to the buffer layer, first, second, and third graphene layer, respectively. The c33,i elements are obtained from fits to the linear part of the stress−strain curves. (c) Calculated layer-dependent stiffness values c33,i for the (6√3 × 6√3)-R30° ZLG, MLG, and BLG interface models (Figure 2a). (d) Reduced modulus ratios (EII*/EI*)th predicted for different atomic structure models compared with experimentally measured modulus ratios (EII*/EI*)exp on the ZLG and MLG samples in Figure 2c,i, respectively. The error bar for (EII*/EI*)th is calculated by considering the uncertainty */EMLG * )exp is compared to the predicted modulus ratio in tip radius Rt calibration (see SI Section I). The experimentally measured ratio (EO * /EMLG * )th, but the actual atomic structure of region O is unknown. (EZLG

Figure 3c compares the layer-dependent c33,i derived from (6√3 × 6√3)-R30° ZLG, MLG, and BLG slab models. Three distinct layer stiffness regimes are apparent. The highest c33,i values in the range of 500−520 GPa are found for the SiC bilayers (Si). The pristine buffer layer (G0 in ZLG) yields c33,i ≈ 270 GPa. When graphene is present on the buffer layer, its c33,i decreases to 231 GPa. As expected, the vdW-bonded graphene planes G1 and G2 exhibit the lowest stiffness c33,i ≈ 60 GPa. Figure 3d compares reduced elastic modulus ratios (E*II / EI*)th, predicted for different atomic structure models, to the experimentally measured ratios (EII*/EI*)exp on the ZLG and MLG samples (Figure 2c,i). For the ZLG sample, (EZLG * / E*MLG)th = 1.088 is in excellent agreement with (E*ZLG/E*MLG)exp = 1.082 ± 0.024, validating our structure assignments of these areas. Likewise, (EBLG * /EMLG * )th = 0.936 agrees nearly perfectly * /EMLG * )exp = 0.928 ± 0.012 on the MLG sample. This with (EBLG leaves the question of the identity of the finger-like structures, labeled “region O” in Figure 2i. The measured ratio (E*O/ EMLG * )exp = 1.140 ± 0.028 is high and could indicate a substrate area not covered by graphene. In Figure 3d, we thus compare * )exp to the predicted ratio (EZLG * /EMLG * )th = 1.094 on (EO*/EMLG the MLG sample. The discrepancy between the experimental ratio and the ratio calculated assuming ZLG suggests that “region O” cannot be assigned unequivocally to a simple ZLG surface area. The growth models in refs 38 and 39 and the

bonded bulk materials (using graphite in the case of graphene). The out-of-plane tensor elements c33,i, schematically indicated in Figure 3a, depend on the atomic plane number i and are obtained from detailed atomic structure models of each region. The experimentally observed (6√3 × 6√3)-R30° periodicity of graphene leads to challenging model sizes (1742 atoms per supercell for the ZLG slab model with six bilayers of SiC plus 338 atoms per supercell for each additional graphene layer), handled by the FHI-aims all-electron code45,46 with “tight” numerical settings (SI Section V and Table S8) and the ELPA eigenvalue solver library.47,48 For continuous media, Hooke’s law relates the stress σ3 and the layer-dependent strain ϵ3,i in the z-direction via the rank-4 stiffness tensor element c33,i σ3 = c33, i ϵ3, i

(5)

As shown in Figure 3b, we calculate the c33,i by applying an external surface stress σ3 and averaging the resulting atomic displacements in each plane to obtain the strain, ϵ3, i =

⟨Δzi⟩ ⟨Δzi0⟩

− 1, where ⟨Δzi⟩ and ⟨Δz0i ⟩ are the laterally

averaged, strained and unstrained layer distances, respectively. Interplanar spring constants ki (Figure 3a) are determined from ki =

c33, iA ⟨Δzi⟩

, where A denotes the loaded area. 6495

DOI: 10.1021/acsnano.6b02402 ACS Nano 2016, 10, 6491−6500

Article

ACS Nano

Figure 4. CR-AFM image on an oxygen-intercalated MLG sample reveals distinct subsurface regions. (a) Height and (b) corresponding CR frequency mapping of the first flexural mode. (c,d) Zoom-in height and frequency images, respectively, of the area shown in the red box in (b). CR-AFM reveals three distinct regions marked “I”, “II”, and “III” in (d). (e) Comparison of reduced moduli ratios using various models to the experimental ratio of region “II” to “III”. The white area reflects the experimental error range. (f) Atomistic structure models of the −OH interface (left) and Si2O5 interface (right). O atoms are shown in red. For domain-averaged f and F⊥ in subfigure (d), see SI Table S3.

occurrence of actual MLG fingers in Figure 2h suggest that “region O” could be related to, for example, oxidized former MLG areas no longer covered by graphene after handling the sample in air. While a detailed search for definite atomic structure models for these areas is outside the scope of this paper, the success of our nanomechanical fingerprint is that we can clearly distinguish “region O” fingers from MLG fingers, and that we can rule out a simple ZLG model for “region O” fingers and similar areas on the MLG sample. CR-AFM of Oxygen-Intercalated Graphene on SiC. Next, we demonstrate the power of our approach for graphene with deliberately modified subsurface interfaces, revealing subsurface features and information on atomic structures. By annealing MLG samples in air, oxygen and/or water vapor can passivate the Si dangling bonds and break Si−C bonds at the buffer layer−SiC interface.4,5,49−52 A pristine MLG region can thus be transformed to oxygen-intercalated, quasi-free-standing bilayer graphene (O-QF-BLG) without a buffer layer underneath, as verified by Raman spectroscopy and angle-resolved photoemission spectroscopy (ARPES).4 Similarly, pristine ZLG can be transformed to O-QF-MLG.52 In Figure 4a−d, we show topographic and corresponding CR-AFM images of a MLG sample that was oxygen-intercalated as described in ref 4. Raman spectroscopy confirms the absence of a buffer layer (SI Figure S5), indicating an essentially defect-free graphene surface layer. Likewise, AFM height images are consistent with homogeneous graphene-covered terraces (Figure 4a,c). In contrast, CR-AFM (Figure 4b,d) reveals a surprising richness of structural features and morphological variability that is not apparent in either AFM height images or Raman spectroscopy. These features are thus very likely located underneath the seemingly intact graphene surface layer. As an example, we focus on the zoomed-in region in Figure 4d. Here, the nanomechanical response helps us to distinguish three different

regions, labeled I, II, and III. Regions II and III are both extended with similar but not identical CR frequencies and are also seen in topographic AFM. In contrast, region I, which is readily apparent in CR-AFM and consists of localized islands ≈100 nm in diameter with a much higher CR frequency, is not identifiable in the topography. Visually, Figure 4c suggests that region I is more abundant within region II than within region III. Although (as shown below) a definitive assignment of all three regions by theory alone remains ambiguous without clear starting points for atomistic structure models, the qualitative revelation of abundant, localized subsurface features, not apparent in AFM height images or Raman, is a spectacular accomplishment of CR-AFM. Considering the growth history of this sample, one would expect the areas near the substrate step edges (region III in Figure 4c) to correspond to BLG that was subsequently intercalated with oxygen. Region II, away from the substrate steps, could be MLG, now also intercalated with oxygen. In a simple scenario (supported by ARPES and Raman4), region II should therefore now correspond to O-QF-BLG, while region III would correspond to oxygen-intercalated, quasi-free-standing trilayer graphene (O-QF-TLG). CR-AFM shows immediately that the reality is more complex and varied, revealing that region II is in detail inhomogeneous and includes large numbers of localized region I defects where oxidation has likely progressed further. However, for a quantitative comparison to theory, we do not have sufficient information about the oxidized structure formed underneath the apparently intact surface graphene planes. As shown in Figure 4e, the experimentally observed reduced moduli ratio (E*II /E*III)exp = 1.044 ± 0.031. Assuming areas II and III to be MLG and BLG supported by a buffer layer G0, we predict (E*II /E*III)th = 1.11 ± 0.01, that is, significantly higher than experiment. As an example of an oxide-supported model, we consider an 6496

DOI: 10.1021/acsnano.6b02402 ACS Nano 2016, 10, 6491−6500

Article

ACS Nano intercalated two-dimensional (√3 × √3)-R30° Si 2 O 3 (“silicate”) layer, connected by a linear Si−O−Si bond to SiC. If not covered by graphene, this model is identical to an experimentally established Si2O5 superstructure on SiC(0001) (Figure 4f).53 However, with 1.09 ± 0.01, the theoretically predicted modulus ratio of BLG supported by silicate (BLGSi2O5) to TLG supported by silicate (TLGSi2O5) is still higher than experiment (see Methods and SI Section VI). As long as the interface underneath the BLG and TLG is kept identical (shown for hypothetical OH-terminated interface models and G0 buffer layer interfaces in Figure 4e), the predicted (EII*/E*III)th values are similarly high. Thus, the interfaces under regions II and III in Figure 4d are likely not identical. Two different reasonable interface combinations for II/III (BLGSi2O5/TLGG0 andBLGOH/TLGSi2O5) that would be compatible with the experimentally observed elastic modulus ratio are shown as examples in Figure 4e. Another simple possibility is suggested by observing that, in our calculations, the silicate adlayer is softer than SiC itself. Thus, the experimental modulus ratio could also be explained by a SiO2 film underneath the O-QF-BLG of region II thicker than that of the SiO2 film underneath the suspected O-QF-TLG of region III. This assignment would be consistent with the observed occurrence of further oxidized region I embedded in region II. The observed enhanced stiffness of region I could then be due to a slightly volume expanded, localized three-dimensional oxidized domain “pushing” against the covering, suspended bilayer graphene from underneath. However, in all cases, further evidence from other methods or from the growth process itself would be required to reveal definitive atomistic structure models for regions II and III. Beyond Graphene. Finally, we show in Figure 5 that CRAFM is able to resolve differences between other well-defined

200 nN, and f = 200 kHz. Figure 5 shows that the different models yield distinctly different ratios (EII*/EI*)th. This sensitivity does not originate from different c33,i values alone. Rather, the elastic moduli also reflect different in-plane cij and distinct interplanar spacings in the different materials that enter the continuum model described above.

CONCLUSIONS In this paper, we show that CR-AFM, combined with firstprinciples theory and modeling, is a powerful experimental tool to investigate the interfaces and subsurface layers of 2D material heterostructures. CR-AFM enables the nondestructive and quantitative characterization of structural, compositional, and/or interfacial bonding changes in lateral and vertical material heterostructures. Remarkably, CR-AFM is exquisitely sensitive to the stiffness changes that arise from even a single atomic layer of a vdW-adhered material. If likely atomistic structure model candidates are available, a combined DFT and continuum modeling approach can be used to analyze the experimentally obtained contact-resonance frequency maps quantitatively to yield sub-surface-sensitive, structural fingerprints of different regions of interest on the sample. CR-AFM is able to detect subsurface defects easily and with high lateral resolution even if these are not apparent in other surfacesensitive techniques such as AFM height imaging or Raman spectroscopy. In combination with first-principles-calculated elastic stiffness parameters of 2D interfaces and atomic layers, the continuum mechanics model can be readily solved without expensive computational resources. For cases where welldefined structure models of individual surface domains are known, our approach thus constitutes a versatile nanomechanical technique to resolve and interpret surface and subsurface structure across the wide field of 2D materials and their heterostructures. METHODS Graphene Growth on SiC(0001). Epitaxial graphene was grown on a 6H-SiC(0001) substrate, cut from a nominally on-axis, 2 in., ndoped SiC wafer, polished on the (0001) face. The substrates were chemically cleaned in n-butylacetate, acetone, and methanol. Hydrogen etching and graphene growth were both performed in an induction heating furnace. Hydrogen etching treatment was carried out at 1673 K for 15 min in a forming gas atmosphere (95 atom % Ar and 5 atom % H) of 900 mbar and a flow rate of 500 sccm. The “ZLG sample” was then produced23 by thermally treating the samples in the same system at a temperature of 1723 K for 15 min in an Ar atmosphere of 900 mbar and a flow rate of 100 sccm. The “MLG sample” was prepared4,37 by exposure to 900 mbar Ar atmosphere also with a flow rate of 500 sccm at 1873 K for 15 min. Oxygen intercalation was achieved4 by thermally treating the monolayer epitaxial graphene on 6H-SiC (0001) for 40 min at 873 K in air with a preceding heating ramp of 50 K/min. Raman spectra of the samples were acquired with a Horiba Jobin Yvon LabRam Aramis Raman microscope equipped with CCD camera and a HeNe laser (λ ∼ 633 nm). During the measurement, the laser spot (∼1 μm) was placed inside the SiC terraces to minimize the signal disturbance from the edges. All Raman data were corrected by subtracting the spectrum from that of a H-etched SiC sample obtained from the same wafer.4,15 AFM Measurements. All AFM measurements were conducted with an Asylum MFP-3D AFM (Oxford Instrument, CA) in dry, ambient environment (20% ≤ RH ≤ 28%, T ∼ 297 ± 1 K). Prior to any measurement, the deflection sensitivity was calibrated by engaging the cantilever on a silicon surface. The spring constant kc of the cantilever was determined from the power spectral density of the thermal noise fluctuations in air54 by fitting the first free resonant peak

Figure 5. Calculated reduced elastic moduli ratios E*th,II/E*th,I for different 2D layer models (II = graphene, h-BN, MoS2, or MoO3) suspended on H-QF-MLG vs the bare H-QF-MLG model surface as reference (area I). Assumed tip parameters: radius Rt = 100 nm; force F⊥ = 200 nN; frequency f = 200 kHz.

2D material domains. Specifically, we consider conceptual (√3 × √3)-R30° interface models for 2D hexagonal boron nitride, molybdenum disulfide, and molybdenum trioxide on hydrogenintercalated, quasi-free-standing monolayer graphene (H-QFMLG) on SiC(0001).13 To mimic typical CR-AFM experimental input parameters (Figure 1), we use Rt = 100 nm, F⊥ = 6497

DOI: 10.1021/acsnano.6b02402 ACS Nano 2016, 10, 6491−6500

Article

ACS Nano to equations for a simple harmonic oscillator55 (see SI Table S4 for calibration details). For CR-AFM mapping, the built-in dual actuation resonance tracking approach of the Asylum MFP-3D AFM was used to track the contact-resonance frequency while scanning the sample surface.31 Ultrasonic actuation was implemented by gluing the sample on a small ultrasonic transducer (see Figure 1a) with a broadband resonance of 2.25 MHz (V133-RM, Olympus NDT). The driving signal sent from the AFM controller ranged between 100 and 200 mV and yielded 10−30 mV of contact-resonance amplitude in air. The local interaction between the AFM cantilever tip and the sample surface was determined from adhesion measurements on the region of interest in the sample. The relative tip position L1/L (see Figure 1c) was determined by measuring the first and second flexural mode resonance frequencies of the cantilever in constant force contact with a silicon surface. The tip radius Rt was determined by measuring the contact stiffness on a silicon wafer using the known bulk material properties of Si56 and applying a Hertzian contact model (see below). The average CR frequency of a region in the image is the mean value of the Gaussian fit to the frequency distribution of the entire region (see SI Figure S3 for instance). Hertzian Contact. Hertzian contact theory describes the elastic contact of a sphere with a plane (Figure 1d) and is often used to approximate an AFM tip indenting a sample surface with a force F⊥. The contact radius, a, and the contact stiffness, k*, are given by32

a=

3

3F ⊥R * , k* = 4E*

3

k* 6F ⊥R*E*2 , E* = 2a

since the lateral cij continue to be taken from unstrained graphene planes. A trilayer graphene-covered variant of the (√3 × √3)-R30° model in the SI (Figure S7) demonstrates the validity range of the linear fitting of c33,i by stress−strain curves. We finally explored a more recent description of dispersion terms than PBE+TS, based on many-body dispersion (PBE+MBD).58 However, this method is untested for the mechanical response properties at issue here and, in its present implementation, could only be applied to small supercell models. The calculations were performed based on the (√3 × √3)-R30° interface model. The agreement of the calculated elastic modulus ratios with the experimental data is somewhat worse than the agreement obtained by the PBE+TS method (SI Figure S6). Continuum Model for Layered Material Systems. Here, we consider the effective mechanical impedance, Zeff, of a layered materials system following the approach outlined in ref 34; for further details, see SI Section VII. Briefly, the imaginary part of the effective impedance can be related to the sample surface stiffness, k*s ks* =

Es* =

1 1 1 = + E* Et* Es*

(8)

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acsnano.6b02402. Detailed description of the methods and materials used in this work (PDF)

AUTHOR INFORMATION Corresponding Authors

*E-mail: [email protected]. *E-mail: [email protected]. Author Contributions §

Q.T. and B.L. contributed equally to this work.

Ei 1 − νi

(11)

S Supporting Information *

Rt and Rs are the curvature radii of the AFM tip and the sample surface, respectively. As Rs → ∞, R* → Rt. Et* and Es* are the reduced moduli of the AFM tip and the sample, respectively, which are related to the mechanical properties of the materials. The reduced modulus of an isotropic material Ei* is Ei* =

ks* 2a

ASSOCIATED CONTENT

where R* and E* are the reduced radius and the reduced modulus, respectively, and are given by (7)

(10)

which yields the reduced moduli of the layered material system

(6)

1 1 1 = + R* Rt Rs

normal force on contanct = Im(2fπ 2a 2Zeff ) displacement

2

Notes

(9)

The authors declare no competing financial interest.

where Ei and νi are the Young’s modulus and the Poisson ratio of the material, respectively. For anisotropic or layered materials, the reduced modulus of the material is a more complicated expression and can be calculated by the method proposed in ref 34 (see below). DFT Calculations of Out-of-Plane Mechanical Properties for Epitaxial Graphene on SiC. Supercells for DFT calculations of (6√3 × 6√3)-R30° models as described in the main text contain six SiC bilayers and up to three layers of graphene. The total number of atoms per supercell is thus 1742, 2080, and 2418 atoms for ZLG, MLG, and BLG, respectively. The carbon side of the slab model was passivated with hydrogen atoms (Figure 1e). To model bulk-like conditions, the two bottom-most SiC bilayer were kept fixed during the ionic relaxation. The atomic displacements in the slab model were then averaged over individual layers. By applying different, fixed external surface stresses σ3 to the atoms of the topmost surface plane and taking finite differences of the resulting displacements, we obtained the layer-dependent, out-of-plane stiffness tensor element c33,i by fitting to the linear range of the stress−strain curves (Figure 3b and eq 5). For computational expediency, we also performed calculations based on a much smaller (√3 × √3)-R30° interface model for conceptual models of the O-intercalated surfaces. In order to add commensurate graphene planes to this model, the graphene planes incur ≈8% strain.57 For MLG and BLG, we verified that the deviations of calculated c33,i to the full (6√3 × 6√3)-R30° models are small

ACKNOWLEDGMENTS V.B., S.Z., Q.T., and B.L. gratefully acknowledge support from Duke University’s Energy Research Seed Fund. This work was financially supported by the NSF through the Research Triangle MRSEC (DMR-11-21107). This research used resources of the Argonne Leadership Computing Facility, which is a DOE Office of Science User Facility supported under Contract DE-AC02-06CH11357. The authors thank Dr. Robert Glass and Cree Inc. for providing SiC wafers free of charge. REFERENCES (1) Ferrari, A. C.; Bonaccorso, F.; Fal’ko, V.; Novoselov, K. S.; Roche, S.; Boggild, P.; Borini, S.; Koppens, F. H. L.; Palermo, V.; Pugno, N.; Garrido, J. A.; Sordan, R.; Bianco, A.; Ballerini, L.; Prato, M.; Lidorikis, E.; Kivioja, J.; Marinelli, C.; Ryhanen, T.; Morpurgo, A.; et al. Science and Technology Roadmap for Graphene, Related TwoDimensional Crystals, and Hybrid Systems. Nanoscale 2015, 7, 4598− 4810. (2) Kim, K. K.; Hsu, A.; Jia, X.; Kim, S. M.; Shi, Y.; Hofmann, M.; Nezich, D.; Rodriguez-Nieva, J. F.; Dresselhaus, M.; Palacios, T.; Kong, J. Synthesis of Monolayer Hexagonal Boron Nitride on Cu Foil Using Chemical Vapor Deposition. Nano Lett. 2012, 12, 161−166. 6498

DOI: 10.1021/acsnano.6b02402 ACS Nano 2016, 10, 6491−6500

Article

ACS Nano (3) Hertel, S.; Waldmann, D.; Jobst, J.; Albert, A.; Albrecht, M.; Reshanov, S.; Schöner, A.; Krieger, M.; Weber, H. B. Tailoring the Graphene/Silicon Carbide Interface for Monolithic Wafer-Scale Electronics. Nat. Commun. 2012, 3, 957−6. (4) Oliveira, M. H., Jr.; Schumann, T.; Fromm, F.; Koch, R.; Ostler, M.; Ramsteiner, M.; Seyller, T.; Lopes, J. M. J.; Riechert, H. Formation of High-Quality Quasi-Free-Standing Bilayer Graphene on SiC(0001) by Oxygen Intercalation upon Annealing in Air. Carbon 2013, 52, 83− 89. (5) Oida, S.; McFeely, F. R.; Hannon, J. B.; Tromp, R. M.; Copel, M.; Chen, Z.; Sun, Y.; Farmer, D. B.; Yurkas, J. Decoupling Graphene from SiC(0001) via Oxidation. Phys. Rev. B: Condens. Matter Mater. Phys. 2010, 82, 041411. (6) Mathieu, C.; Lalmi, B.; Mentes, T. O.; Pallecchi, E.; Locatelli, A.; Latil, S.; Belkhou, R.; Ouerghi, A. Effect of Oxygen Adsorption on the Local Properties of Epitaxial Graphene on SiC(0001). Phys. Rev. B: Condens. Matter Mater. Phys. 2012, 86, 035435. (7) Yu, Y.; Huang, S.-Y.; Li, Y.; Steinmann, S. N.; Yang, W.; Cao, L. Layer-Dependent Electrocatalysis of MoS2 for Hydrogen Evolution. Nano Lett. 2014, 14, 553−538. (8) Zhu, W.; Perebeinos, V.; Freitag, M.; Avouris, P. Carrier Scattering, Mobilities, and Electrostatic Potential in Monolayer, Bilayer, and Trilayer Graphene. Phys. Rev. B: Condens. Matter Mater. Phys. 2009, 80, 235402. (9) Zhang, Y.; Tang, T.-T.; Girit, C.; Hao, Z.; Martin, M. C.; Zettl, A.; Crommie, M. F.; Shen, Y. R.; Wang, F. Direct Observation of a Widely Tunable Bandgap in Bilayer Graphene. Nature 2009, 459, 820−823. (10) Liu, Z.; Ma, L.; Shi, G.; Zhou, W.; Gong, Y.; Lei, S.; Yang, X.; Zhang, J.; Yu, J.; Hackenberg, K. P.; Babakhani, A.; Idrobo, J.-C.; Vajtai, R.; Lou, J.; Ajayan, P. M. In-Plane Heterostructures of Graphene and Hexagonal Boron Nitride with Controlled Domain Sizes. Nat. Nanotechnol. 2013, 8, 119−124. (11) Chen, W.; Chen, S.; Qi, D. C.; Gao, X. Y.; Wee, A. T. S. Surface Transfer P-Type Doping of Epitaxial Graphene. J. Am. Chem. Soc. 2007, 129, 10418−10422. (12) Geim, A. K.; Grigorieva, I. V. Van Der Waals Heterostructures. Nature 2013, 499, 419−425. (13) Sforzini, J.; Nemec, L.; Denig, T.; Stadtmüller, B.; Lee, T. L.; Kumpf, C.; Soubatch, S.; Starke, U.; Rinke, P.; Blum, V.; Bocquet, F. C.; Tautz, F. S. Approaching Truly Freestanding Graphene: The Structure of Hydrogen-Intercalated Graphene on 6H-SiC(0001). Phys. Rev. Lett. 2015, 114, 106804. (14) Kumar, H.; Er, D.; Dong, L.; Li, J.; Shenoy, V. B. Elastic Deformations in 2D Van Der Waals Heterostructures and Their Impact on Optoelectronic Properties: Predictions from a Multiscale Computational Approach. Sci. Rep. 2015, 5, 10872. (15) Fromm, F.; Oliveira, M. H., Jr.; Molina-Sánchez, A.; Hundhausen, M.; Lopes, J. M. J.; Riechert, H.; Wirtz, L.; Seyller, T. Contribution of the Buffer Layer to the Raman Spectrum of Epitaxial Graphene on SiC(0001). New J. Phys. 2013, 15, 043031. (16) Lee, D. S.; Riedl, C.; Krauss, B.; von Klitzing, K.; Starke, U.; Smet, J. H. Raman Spectra of Epitaxial Graphene on SiC and of Epitaxial Graphene Transferred to SiO2. Nano Lett. 2008, 8, 4320− 4325. (17) Woods, C. R.; Britnell, L.; Eckmann, A.; Ma, R. S.; Lu, J. C.; Guo, H. M.; Lin, X.; Yu, G. L.; Cao, Y.; Gorbachev, R. V.; Kretinin, A. V.; Park, J.; Ponomarenko, L. A.; Katsnelson, M. I.; Gornostyrev, Y. N.; Watanabe, K.; Taniguchi, T.; Casiraghi, C.; Gao, H.-J.; Geim, A. K.; et al. Commensurate-Incommensurate Transition in Graphene on Hexagonal Boron Nitride. Nat. Phys. 2014, 10, 451−456. (18) Gong, Y.; Lin, J.; Wang, X.; Shi, G.; Lei, S.; Lin, Z.; Zou, X.; Ye, G.; Vajtai, R.; Yakobson, B. I.; Terrones, H.; Terrones, M.; Tay, B. K.; Lou, J.; Pantelides, S. T.; Liu, Z.; Zhou, W.; Ajayan, P. M. Vertical and In-Plane Heterostructures from WS2/MoS2 Monolayers. Nat. Mater. 2014, 13, 1135−1142. (19) Hibino, H.; Kageshima, H.; Nagase, M. Epitaxial Few-Layer Graphene: Towards Single Crystal Growth. J. Phys. D: Appl. Phys. 2010, 43, 374005.

(20) Emtsev, K. V.; Bostwick, A.; Horn, K.; Jobst, J.; Kellogg, G. L.; Ley, L.; McChesney, J. L.; Ohta, T.; Reshanov, S. A.; Röhrl, J.; Rotenberg, E.; Schmid, A. K.; Waldmann, D.; Weber, H. B.; Seyller, T. Towards Wafer-Size Graphene Layers by Atmospheric Pressure Graphitization of Silicon Carbide. Nat. Mater. 2009, 8, 203−207. (21) Riedl, C.; Coletti, C.; Starke, U. Structural and Electronic Properties of Epitaxial Graphene on SiC(0001): a Review of Growth, Characterization, Transfer Doping and Hydrogen Intercalation. J. Phys. D: Appl. Phys. 2010, 43, 374009−18. (22) Heinz, K.; Starke, U.; Bernhardt, J.; Schardt, J. Surface Structure of Hexagonal SiC Surfaces: Key to Crystal Growth and Interface Formation? Appl. Surf. Sci. 2000, 162, 9−18. (23) Schumann, T.; Dubslaff, M.; Oliveira, M., Jr.; Hanke, M.; Fromm, F.; Seyller, T.; Nemec, L.; Blum, V.; Scheffler, M.; Lopes, J. M. J.; Riechert, H. Structural Investigation of Nanocrystalline Graphene Grown on (6√3 × 6√3) R 30o-Reconstructed SiC Surfaces by Molecular Beam Epitaxy. New J. Phys. 2013, 15, 123034. (24) Fei, Z.; Rodin, A. S.; Andreev, G. O.; Bao, W.; McLeod, A. S.; Wagner, M.; Zhang, L. M.; Zhao, Z.; Thiemens, M.; Dominguez, G.; Fogler, M. M.; Neto, A. H. C.; Lau, C. N.; Keilmann, F.; Basov, D. N. Gate-Tuning of Graphene Plasmons Revealed by Infrared NanoImaging. Nature 2012, 487, 82−85. (25) Hoffmann, G. G.; de With, G.; Loos, J. Micro-Raman and TipEnhanced Raman Spectroscopy of Carbon Allotropes. Macromol. Symp. 2008, 265, 1−11. (26) Saito, Y.; Verma, P.; Masui, K.; Inouye, Y.; Kawata, S. NanoScale Analysis of Graphene Layers by Tip-Enhanced Near-Field Raman Spectroscopy. J. Raman Spectrosc. 2009, 40, 1434−1440. (27) Schmid, T.; Opilik, L.; Blum, C.; Zenobi, R. Nanoscale Chemical Imaging Using Tip-Enhanced Raman Spectroscopy: A Critical Review. Angew. Chem., Int. Ed. 2013, 52, 5940−5954. (28) Come, J.; Xie, Y.; Naguib, M.; Jesse, S.; Kalinin, S. V.; Gogotsi, Y.; Kent, P. R. C.; Balke, N. Nanoscale Elastic Changes in 2D Ti3C2Tx (MXene) Pseudocapacitive Electrodes. Adv. Energy Mater. 2016, 6, 1502990. (29) Rabe, U.; Janser, K.; Arnold, W. Vibrations of Free and SurfaceCoupled Atomic Force Microscope Cantilevers: Theory and Experiment. Rev. Sci. Instrum. 1996, 67, 3281−3293. (30) de Jong, M.; Chen, W.; Angsten, T.; Jain, A.; Notestine, R.; Gamst, A.; Sluiter, M.; Ande, C. K.; van der Zwaag, S.; Plata, J. J.; Toher, C.; Curtarolo, S.; Ceder, G.; Persson, K. A.; Asta, M. Charting the Complete Elastic Properties of Inorganic Crystalline Compounds. Sci. Data 2015, 2, 150009. (31) Rodriguez, B. J.; Callahan, C.; Kalinin, S. V.; Proksch, R. DualFrequency Resonance-Tracking Atomic Force Microscopy. Nanotechnology 2007, 18, 475504. (32) Rabe, U. In Applied Scanning Probe Methods II; Bhushan, B., Fuchs, H., Eds.; NanoScience and Technology; Springer: Berlin, 2006; Chapter on Atomic Force Acoustic Microscopy, pp 37−90. (33) Garcia, R.; Herruzo, E. T. The Emergence of Multifrequency Force Microscopy. Nat. Nanotechnol. 2012, 7, 217−226. (34) Yaralioglu, G. G.; Degertekin, F. L.; Crozier, K. B.; Quate, C. F. Contact Stiffness of Layered Materials for Ultrasonic Atomic Force Microscopy. J. Appl. Phys. 2000, 87, 7491−7496. (35) Nemec, L.; Blum, V.; Rinke, P.; Scheffler, M. Thermodynamic Equilibrium Conditions of Graphene Films on SiC. Phys. Rev. Lett. 2013, 111, 065502. (36) Riedl, C.; Starke, U. Structural Properties of the GrapheneSiC(0001) Interface as a Key for the Preparation of Homogeneous Large-Terrace Graphene Surfaces. Phys. Rev. B: Condens. Matter Mater. Phys. 2007, 76, 245406. (37) Oliveira, M. H., Jr.; Schumann, T.; Ramsteiner, M.; Lopes, J. M. J.; Riechert, H. Influence of the Silicon Carbide Surface Morphology on the Epitaxial Graphene Formation. Appl. Phys. Lett. 2011, 99, 111901. (38) Borovikov, V.; Zangwill, A. Step-Edge Instability during Epitaxial Growth of Graphene from SiC(0001). Phys. Rev. B: Condens. Matter Mater. Phys. 2009, 80, 121406. 6499

DOI: 10.1021/acsnano.6b02402 ACS Nano 2016, 10, 6491−6500

Article

ACS Nano (39) Tetlow, H.; de Boer, J. P.; Ford, I. J.; Vvedensky, D. D.; Coraux, J.; Kantorovich, L. Growth of Epitaxial Graphene: Theory and Experiment. Phys. Rep. 2014, 542, 195−295. (40) Magonov, S.; Elings, V.; Whangbo, M.-H. Phase Imaging and Stiffness in Tapping-Mode Atomic Force Microscopy. Surf. Sci. 1997, 375, L385−L391. (41) Garcia, R.; Perez, R. Dynamic Atomic Force Microscopy Methods. Surf. Sci. Rep. 2002, 47, 197−301. (42) Tkatchenko, A.; Scheffler, M. Accurate Molecular Van Der Waals Interactions from Ground-State Electron Density and FreeAtom Reference Data. Phys. Rev. Lett. 2009, 102, 073005. (43) Perdew, J. P.; Burke, K.; Ernzerhof, M. Generalized Gradient Approximation Made Simple. Phys. Rev. Lett. 1996, 77, 3865−3868. (44) Beckstein, O.; Klepeis, J. E.; Hart, G. L. W.; Pankratov, O. FirstPrinciples Elastic Constants and Electronic Structure of α-Pt2Si and Pt2Si. Phys. Rev. B: Condens. Matter Mater. Phys. 2001, 63, 134112. (45) Blum, V.; Gehrke, R.; Hanke, F.; Havu, P.; Havu, V.; Ren, X.; Reuter, K.; Scheffler, M. Ab initio Molecular Simulations with Numeric Atom-Centered Orbitals. Comput. Phys. Commun. 2009, 180, 2175−2196. (46) Havu, V.; Blum, V.; Havu, P.; Scheffler, M. Efficient O(N) Integration for All-Electron Electronic Structure Calculation Using Numeric Basis Functions. J. Comput. Phys. 2009, 228, 8367−8379. (47) Auckenthaler, T.; Blum, V.; Bungartz, H.-J.; Huckle, T.; Johanni, R.; Krämer, L.; Lang, B.; Lederer, H.; Willems, P. R. Parallel Solution of Partial Symmetric Eigenvalue Problems from Electronic Structure Calculations. Parallel Comput. 2011, 37, 783−794. (48) Marek, A.; Blum, V.; Johanni, R.; Havu, V.; Lang, B.; Auckenthaler, T.; Heinecke, A.; Bungartz, H. J.; Lederer, H. The ELPA Library: Scalable Parallel Eigenvalue Solutions for Electronic Structure Theory and Computational Science. J. Phys.: Condens. Matter 2014, 26, 213201. (49) Bom, N. M.; Oliveira, M. H., Jr.; Soares, G. V.; Radtke, C.; Lopes, J. M. J.; Riechert, H. Synergistic Effect of H2O and O2 on the Decoupling of Epitaxial Monolayer Graphene from SiC(0001) via Thermal Treatments. Carbon 2014, 78, 298−304. (50) Ostler, M. Quasi-Free-Standing Graphene on Silicon Carbide. Ph.D. Thesis, University of Erlangen-Nürnberg, 2014. (51) Ostler, M.; Koch, R. J.; Speck, F.; Fromm, F.; Vita, H.; Hundhausen, M.; Horn, K.; Seyller, T. Decoupling the Graphene Buffer Layer from SiC(0001) via Interface Oxidation. Mater. Sci. Forum 2012, 717−720, 649−652. (52) Ostler, M.; Fromm, F.; Koch, R. J.; Wehrfritz, P.; Speck, F.; Vita, H.; Bottcher, S.; Horn, K.; Seyller, T. Buffer Layer Free Graphene on SiC(0001) via Interface Oxidation in Water Vapor. Carbon 2014, 70, 258−265. (53) Bernhardt, J.; Schardt, J.; Starke, U.; Heinz, K. Epitaxially Ideal Oxide-Semiconductor Interfaces: Silicate Adlayers on Hexagonal (0001) and (0001¯) SiC Surfaces. Appl. Phys. Lett. 1999, 74, 1084− 1086. (54) Hutter, J. L.; Bechhoefer, J. Calibration of Atomic-Force Microscope Tips. Rev. Sci. Instrum. 1993, 64, 1868−1873. (55) Walters, D. A.; Cleveland, J. P.; Thomson, N. H.; Hansma, P. K.; Wendman, M. A.; Gurley, G.; Elings, V. Short Cantilevers for Atomic Force Microscopy. Rev. Sci. Instrum. 1996, 67, 3583−3590. (56) Dolbow, J.; Gosz, M. Effect of Out-of-Plane Properties of a Polyimide Film on the Stress Fields in Microelectronic Structures. Mech. Mater. 1996, 23, 311−321. (57) Mattausch, A.; Pankratov, O. Ab initio Study of Graphene on SiC. Phys. Rev. Lett. 2007, 99, 076802. (58) Ambrosetti, A.; Reilly, A. M.; DiStasio, R. A. J.; Tkatchenko, A. Long-Range Correlation Energy Calculated from Coupled Atomic Response Functions. J. Chem. Phys. 2014, 140, 18A508.

6500

DOI: 10.1021/acsnano.6b02402 ACS Nano 2016, 10, 6491−6500