Introduction

Hexagonal boron nitride (h-BN) is a structural analog of graphene, consisting of honeycomb lattices of alternating sp2-hybridized B and N atoms. Its atomic flatness, free of dangling bonds and wide bandgap make it a unique material for various applications, such as the interface dielectric in two dimensional electronics1,2. In this respect, high-quality h-BN crystals are essential to prevent charge trapping and scattering caused by defects, thus ensuring optimal performance of electronic devices. The synthesis of single-crystal h-BN is a promising way to eliminate the generation of grain boundaries and defects. However, comparing to graphene, the scalable growth of single-crystal h-BN is much more challenging due to the binary composition of h-BN, which allows for spontaneous formation of nuclei with opposite orientations on most transition metal substrates3,4,5,6. These nuclei have slightly different formation energies, making unidirectional nucleation overwhelmingly difficult. Recently, large-scale oriented nucleation and growth have been achieved by utilizing inequivalent adsorption sites on low-symmetry substrate or lateral docking to step edges, which lift the energy degeneracy of the two oppositely aligned nuclei and result in mono-orientation domains that enable seamless coalescence to form wafer-scale single crystalline h-BN monolayer7,8,9.

For transition metal surfaces with high lattice symmetry, such as (111) and (0001) facets, it has been found that the chemical vapor deposition method shows promise to achieve high-quality h-BN films with minimal density of domain boundaries under certain conditions. For example, monolayer h-BN with macroscopic domain sizes could be grown on Ru(0001) at low precursor pressures due to sparse nucleation, which reduces the density of defects between the merged islands10. In the case of temperature-programmed growth (TPG), the two-step process consisting of low-temperature adsorption of molecular precursors on the substrate and subsequent annealing at high temperature gives rise to an extended h-BN monolayer on Ir(111)11. In situ characterizations by photoelectron diffraction spectroscopy and low energy electron diffraction (LEED) suggest that the scalable growth is facilitated by the formation of domains with highly controlled orientation at the initial stages of nucleation as a result of the interplay between thermal energy and binding energy discrepancies of different types of nuclei. The favored adsorption configuration was identified to be NtopBfcc, referring to N sitting on top of Ir atom and B on fcc site respectively11. However, this assignment was debated by a later study that proposed NtopBhcp, i.e., N sitting on top site and B on hcp site, to be the energetically preferred registry5. The result is consistent with early density functional theory (DFT) calculations suggesting that the NtopBhcp configuration is the lower energy one for h-BN clusters of all sizes12. In addition, there is also a debate regarding the edge termination of h-BN nuclei. Although most results believe that the N-type zigzag edge is more stable on transition metal substrates4,7,8, several works suggest that the B-zigzag edge termination is energetically favored5,13. So far, there are few experiments that could directly identify the edge composition of h-BN clusters formed in the initial nucleation stage.

The discrepancies in observations regarding the nucleation pathway, orientation, and structural evolution of h-BN nuclei on metallic surfaces indicate that a comprehensive understanding of these processes has yet to be achieved. Given the complexity of these mechanisms, an atomic-scale structural elucidation of h-BN during the early stage of nucleation is of paramount importance. The non-contact atomic force microscopy (nc-AFM) enables elemental contrast imaging on monolayer h-BN by combining Kelvin probe force microscopy experiments and nc-AFM image simulations based on DFT calculations14. Therefore, nc-AFM would be a powerful tool for the atomic-scale investigation of h-BN nuclei.

In this study, we investigate the nucleation process of h-BN on Ir(111) using scanning tunneling microscopy (STM) and nc-AFM corroborated by DFT calculations. The choice of Ir(111) as a substrate is based on its moderate reactivity, in contrast to strongly interacting surfaces like Ru(0001) and Rh(111) or weakly interacting surfaces like Cu(111), making it a representative system for studying h-BN growth. We experimentally determine the atomic configuration and precise adsorption sites of h-BN clusters on Ir(111) via element-sensitive structural characterization. Our observations reveal that unidirectional aligned h-BN nanoplates dominate the nucleation process without showing any dependence on size. Based on these findings, we propose an orientation pinning mechanism that explains the nucleation of h-BN on Ir(111). These results provide atomic-scale insights into the early stages of h-BN nucleation and can aid in developing new strategies for growing high-quality h-BN films on metallic substrates.

Results and discussions

Preparation and structural characterizations of h-BN nuclei

The synthesis scheme of h-BN on Ir(111) is illustrated in Fig. 1a. Following the TPG method described in ref. 11, vaporized precursor15,16 was deposited on the substrate at room temperature with a subsequent annealing at 700 K. Large-scale STM images show that nanoscale motifs of various sizes sparsely decorate the Ir(111) surface (Supplementary Fig. 1). They are uniformly distributed on terraces, with no clear tendency for step attachment. Close-up STM and corresponding nc-AFM images of h-BN nanoclusters are shown in Fig. 1b-f. Although no structural details are resolved in STM topographies, the nc-AFM images acquired using metallic tip reveal triangular shape with honeycomb lattices and zigzag-type edges.

Fig. 1: Preparation and structural characterization of h-BN nuclei on Ir(111).
figure 1

a Schematic illustration of the two-step temperature-programmed growth of h-BN nanoplates on Ir(111). The precursors, borazine (B3N3H6) and aminoborane (BH2·NH2), are generated from thermal decomposition of ammonia borane (BH3·NH3). b–f STM and the corresponding nc-AFM images of h-BN nuclei for (b) n = 3 (n denotes the number of hexagonal rings on the edge of the h-BN nuclei), c n = 4, d n = 5, e n = 6, and f n = 7. STM parameters: V = 0.5 V, I = 10 pA; nc-AFM parameters: V = 0 V, amplitude = 100 pm. Scale bars: b–f, 5 Å. The STM and nc-AFM images are acquired using metallic tip.

Full size image

The smallest h-BN nucleus observed in our experiment maintains the triangular shape and possesses three hexagonal rings on the edge (referred to n = 3), as displayed in Fig. 1b. Because the B:N ratio of the triangular nuclei is deviated from the stoichiometric ratio of the precursor, we suggest that the building blocks in the h-BN nuclei formation are elemental atoms of B and N, which is also supported by the previously reported X-ray photoelectron spectroscopy (XPS) results that the atomic fragments are observed at temperatures higher than 330 K due to the break-up of B-N bonds in molecular precursors17. In contrast to polycyclic aromatic hydrocarbons with hydrogen-terminated edge18, the less visible atoms at the peripheries of h-BN nanoclusters suggests that edge atoms are bonded to the underlying Ir(111) rather than being hydrogenated, similar to the cases of graphene nuclei grown on metal substrates in UHV system19. The force spectroscopy ∆f(z) conducted on a central atom and an edge atom in the h-BN nucleus with n = 4 reveals their adsorption heights difference of ~ 22 pm, which reasonably agrees well with the DFT calculations (Supplementary Fig. 2).

Despite of size variation, the h-BN nanoclusters align their edges with \(\left[1\bar{1}0\right]\) or equivalent crystallographic directions of Ir(111), indicating a strict registry between h-BN crystallization and the substrate lattices. Supplementary Fig. 1 presents a large-scale STM image, while Supplementary Fig. 3 provides additional nc-AFM images, both of which further support this conclusion. A close examination of 2097 nanoclusters in these images reveals that the majority of the nanoclusters exhibit a consistent orientation, whereas a smaller fraction orient in the opposite direction, i.e., rotated by 60°. We have not found any significant correlation between the preferred orientation of h-BN nuclei and their size or the nucleation density. The preferential orientation of h-BN nuclei observed in our experiments is consistent with the early study, where LEED patterns obtained on TPG-prepared h-BN domains reveal a favored crystalline orientation11.

In contrast to prior findings with CO-functionalized tips14, h-BN lattices consistently displayed as honeycomb patterns in nc-AFM images acquired using a metal tip (Supplementary Fig. 4)20. An intriguing feature of the nc-AFM images is the apparently larger honeycomb rings in the interior of the h-BN clusters, which is contrary to that of PAHs. In the latter case, imaging distortion is typically observed on the peripheral moieties due to the deflection of CO on the tip apex21. In this study, the enlarged rings have been consistently observed in nc-AFM images acquired using metallic tip. The Δf line profile along the symmetric axis of a h-BN nanoplate with n = 4 (Fig. 2a) reveals that apparent size of the central hexagonal ring of BN is measured to be 3.42 Å (Fig. 2b), which is ~14% larger than the value of bulk h-BN. Although such magnification of bond length might occur for nc-AFM imaging by utilizing CO-functionalized tip, it is unexpected for metallic tips, which typically have rigid geometric structures22. It is worth mentioning that nc-AFM image records the frequency shift, i.e., force gradient between tip and sample. Different from I(z) spectroscopy in STM mode, the Δf(z) spectroscopy in nc-AFM mode is non-monotonic, that is, the slope is either positive or negative on the opposite side of the minimum of Δf(z) curve, as illustrated in the inset of Fig. 2c. Therefore, the constant-height frequency shift imaging does not present the accurate height and location of the adsorbent. The minimum of Δf(z) spectroscopy (Δfmin) provides information about the position where the Pauli repulsion begins to significantly contribute to the total interactions. Thus, the Δfmin can provide valuable information regarding the targeted sample, and may serve as a criterion for estimating the atomic positions and elemental composition of h-BN nuclei.

Fig. 2: Determination of accurate locations of B and N atoms in h-BN nuclei.
figure 2

a High-resolution nc-AFM image of h-BN nanoplate with n = 4. Scale bar: 3 Å. b Frequency shift (Δf) line profile recorded along the white arrow in a at z = 0 pm. c Line map of Δf(z) along the same trace in (a). Data are recorded in constant-height AFM mode shown for various tip-sample separation, displayed in steps of 5 pm. Inset: Schematic illustration of the structural model and the extraction of Δfmin for each force curve. d The Δfmin trajectory along the same arrow. The value of Δfmin is extracted from the line map of (c), and the long-range force has been subtracted from Δf(z) curve before the determination of Δfmin for each location.

Full size image

We acquired two-dimensional Δf map of the h-BN nanoplate in the x-z plane, with the x-axis aligned along the direction indicated by the white arrow in Fig. 2a, and the z-axis perpendicular to the h-BN nanoplate (Fig. 2c). The Δfmin for each position is extracted and plotted as a function of x position. Before that, the long-range forces have been eliminated by excluding the force curve recorded on vicinal Ir(111) surface. The atomic locations in h-BN nanoplate are clearly distinguished as peaks in Δfmin(x) curve (Fig. 2d). Furthermore, the atom located on the vertex of h-BN triangle is also recognized (marked by blue arrow in Fig. 2d), which is invisible in constant-height profile. In this case, the precise location of each atom is determined and the size of the central BN ring is corrected to be 3.1 Å.

Elemental identification of h-BN via force spectroscopy

To distinguish the B and N atoms in h-BN, we conduct atomic force spectroscopy measurements using metallic tip. To minimize the influence of the surrounding environment on the tip-sample interaction, we select the central region of h-BN nanoplates with n = 3, 4, 5 for the force spectroscopy measurements and extract the long-range force contributions from the substrate. As shown in Fig. 3a, when the tip is positioned far from the sample, the force spectra measured on different atoms are indistinguishable, because the van der Waals attractive force is dominant between tip and sample. As the tip approaches the surface, short-range Pauli repulsive force and chemical bonding force become predominant, leading to the observable distinctions in the force spectra obtained from adjacent atoms within h-BN nuclei.

Fig. 3: Identification of B and N atoms in h-BN nanoplates using force spectroscopy.
figure 3

a ∆f(z) spectra measured on B and N atoms around the central region of the h-BN nanoplates with n = 3, 4, 5. b The corresponding short-range force curves calculated from (a). The minimums of force spectra are marked by triangular symbols. The contribution of long-range force has been subtracted from ∆f(z) spectra. Scale bars: 5 Å.

Full size image

The force spectra are separated into two branches around their ∆fmin, as plotted by pink and blue symbols in Fig. 3a. For h-BN nanoplate with n = 3, the Δfmin on the central atom (pink curve) is less negative than that of the neighboring atoms (blue curves). This character is reversed in the case of h-BN nanoplate with n = 5. The force curves recorded on the central ring of h-BN with n = 4 demonstrate alternative Δfmin. Considering the higher electron density of N atoms, it is expected that the N atom would generate a stronger chemical bonding attraction with the metallic tip apex, compared to that for B atoms. Therefore, we tentatively assign the blue curves with larger attractive force to N atoms, and the pink curves to B atoms. The force measured on the B atoms exhibits a maximum attractive value of ~ 0.70 ± 0.01 nN, while reaches ~ 0.85 ± 0.04 nN on N atoms (Fig. 3b). Comparatively, the force obtained using a metal tip is approximately one order of magnitude larger than the force obtained using a CO-tip23 highlighting the enhanced chemical reactivity of the metallic tip apex.

DFT calculations are carried out to elucidate the underlying mechanism for the elemental identification capability of metallic nc-AFM probes. The interaction energies (Eint) between metallic tip and h-BN are calculated for different values of z, which is indicated in Fig. 4c. The values of Eint are further fitted by Morse potential to get the ∆f versus z curves by quadratic differentiating (see Methods for details). The calculated ∆f(z) curve (Fig. 4a) on N atom exhibits more negative value of ∆fmin than that on B, agreeing well with the experimental assignments. The projected density of states (PDOS) analysis demonstrates that as the metallic tip approaches the h-BN clusters, there is a prominent increasing of PDOS near the Fermi energy (EF) for the N atom, while that of B atom remains unaltered (Fig. 4b). The calculated charge density difference reveals that the metallic tip induces a larger electron accumulation on N atom than that on B atom as well (Fig. 4c, d). The increased PDOS of N atom is mainly due to the hybridization of its pz orbital with d orbitals of metallic tip, resulting a larger attraction than B atom. Different metallic tip apexes were employed in our calculations and they consistently produced similar results regarding the elemental identification of h-BN (Supplementary Fig. 5). Thus, the short-range chemical bonding force play a significant role in the interactions between metallic tip and h-BN nanoplates, enabling elemental identification for metallic tip via force spectroscopy. By accurately assigning the B and N atoms, the edge termination of h-BN nanoplates shown in Fig. 3 can be concluded to N-zigzag type.

Fig. 4: DFT calculated interactions between metallic tip and h-BN.
figure 4

a DFT calculation derived ∆f(z) curve on B and N atom of h-BN nanoplate with n = 4 by ultilizing a W tip. Inset: the calculated interaction energies (Eint) on B and N (hollow stars), and the fitting curves (solid lines) by Morse potential. b PDOS of W-tip, and differential PDOS of B atom and N atom during the tip approaching (The differential PDOS of B and N shown here is calculated by PDOSz=7.65 Å – PDOSz=8.95 Å). c, d Charge density difference of h-BN nanoplate with z = 8.95 Å and z = 7.65 Å for tip on (c) N atom and d B atom. The isosurface energy is 0.012 e bohr−3, where blue and red mean accumulation and depletion of charge densities, respectively.

Full size image

To efficiently identify the edge termination of h-BN nuclei in batches, large-scale STM imaging combined with STM simulations are performed. We find that the corners of all h-BN nuclei appear as local depressions in low-bias STM images, as seen in the inset of Fig. 5a. In view of typical edge types (B- or N- zigzag) and atomic registrations (NtopBhcp or NtopBfcc), the STM morphologies of a h-BN nucleus (n = 4) in the four representative scenarios are simulated and presented in Fig. 5b, c. Apparently, the N-edge nuclei in NtopBhcp and NtopBfcc configurations both exhibit truncated triangular characteristics. In contrast, the two nuclei with B-edge both have the regular triangle shape. The good agreement between simulations and experimental observations further consolidates that the h-BN nuclei on Ir(111) are N-terminated at the edges. Besides, DFT calculations reveal that the N atom-terminated edges are energetically favorable than the B-edge in all cases, more pronounced for larger h-BN nuclei (Supplementary Fig. 6).

Fig. 5: Experimental and simulated STM images of h-BN nuclei on Ir(111).
figure 5

a Large-scale STM image of h-BN nanoplates with diverse sizes (V = − 20 mV, I = 0.1 nA). Inset: magnified STM image of h-BN nucleus with n = 4. The white arrows indicate depression contrast in the vertices of triangular h-BN nanoplate. b, c The DFT-optimized models of triangular h-BN nanoplates on Ir(111) with (b) N edge and (c) B edge, and their corresponding STM simulations (V = − 20 mV). The adsorption sites of NtopBhcp and NtopBfcc, represent that N atoms are on top sites and B atoms are on hcp sites and fcc sites of Ir(111), respectively.

Full size image

Orientation-pinning mechanism of h-BN nucleation

Comparing the two orientations of h-BN nanoplates with NtopBhcp and NtopBfcc registrations, we found that the formation energies of the nuclei exhibit a size-dependency. Specifically, the NtopBhcp configuration is energetically more favored than NtopBfcc for nuclei with sizes n = 2 and 3, while for larger nuclei (n = 4,5), NtopBfcc is the more stable registration (Fig. 6). To experimentally determine the adsorption sites of h-BN nuclei on the Ir(111) substrate, we used the “adjusted constant height AFM” method24. This involves scanning a CO-functionalized tip at a small distance from the substrate to extract the stacking sequence of the atoms at step edges, allowing us to distinguish hcp and fcc hollow sites on the terrace (Supplementary Fig. 7). Using this approach, we confirmed that the majority of h-BN nuclei takes the NtopBhcp registry regardless of their size. This seems to conflict with our calculations in Fig. 6 showing that the NtopBhcp configuration prevails for the initial nuclei (n = 2,3), while the NtopBfcc registry is favored as nuclei grow.

Fig. 6: The orientation pinning effect of h-BN nucleation on Ir(111).
figure 6

DFT calculated formation energy of h-BN nuclei with variable size and atomic registration.

Full size image

The observed NtopBhcp registry of h-BN nuclei in our experiment can be explained by considering the rotation energy barrier between the two possible orientations of the h-BN clusters. For h-BN nuclei of n = 4 and larger, the NtopBfcc registry is thermodynamically favored. However, for h-BN nuclei of n = 3, the NtopBhcp configuration is adopted, indicating that the growing nucleus must overcome the rotation barrier from Ntop Bhcp to NtopBfcc. Our calculations suggest that the rotation barrier for h-BN nucleus (n = 4) is approximately 4.4 eV, which is too high to be overcome at the growth temperature. This orientation pinning phenomenon, where the initial clusters formed at the early stages of nucleation determine the orientation distributions of domains in the subsequent growth process, plays a crucial role in the synthesis of single-crystalline h-BN films7,8. The estimated energy difference between the NtopBhcp and NtopBfcc orientations, as determined from the statistical analysis of experimental data on their relative ratio, is approximately 0.144 eV (for details see the discussion in supplementary Fig. 1). This value agrees well with the calculated formation energy difference for the n = 2 h-BN nucleus (0.136 eV as illustrated in Fig. 6), suggesting that the orientation pinning phenomenon occurs in the very early stage of h-BN nucleation. Therefore, the optimal procedure for producing large-area single-crystal h-BN film should involve two steps. The first step is the low-temperature production of mono-oriented h-BN nuclei, which are energetically favored under thermodynamic control. The second step is the growth of nuclei into domains. Because the orientation of the nuclei has been pinned, a relatively higher temperature can assist the thermal decomposition of the precursor and enhance the growth rate of h-BN domains in the second step, provided that the concentrations of the decomposition products are not sufficient to form new nuclei of the critical size.

Conclusion

In summary, nucleation of h-BN on the Ir(111) surface at the early stage is investigated using a combination of STM, nc-AFM and DFT calculation. We demonstrate the the capability of force spectroscopy to identify individual B and N atoms within h-BN nuclei using a reactive metal tip. Based on that, the precise structures of triangular h-BN nuclei are carefully determined, exhibiting a unidirectional orientation with N-terminal zigzag edges and specific atomic registration of NtopBhcp to Ir(111). These characteristics are primarily governed by the initial thermodynamically favored structure and the subsequent orientation pinning effect during domain expansion. We believe that this approach provides a versatile means to investigate the nucleation and growth of other binary and ternary two dimensional materials, facilitating the synthesis of single crystals and large-scale manufacturing processes.

Methods

Preparation of h-BN on Ir(111)

The experiments were performed in an ultrahigh-vacuum (UHV) chamber with a base pressure of 1 × 10−10 mbar. The Ir(111) surface was carefully cleaned using repeated cycles of argon ion sputtering followed by annealing at 1100 K. The h-BN nanoplates were grown on the Ir(111) by depositing precursor molecules at a partial pressure of 8 × 10−8 mbar for 10 mins at room temperature and annealing at 700 K for 10 mins.

STM/nc-AFM characterization

The STM and nc-AFM measurements were performed in a UHV chamber at 5 K equipped with a ScientaOmicron low-temperature STM/nc-AFM system operated with a Nanonis controller. The bias voltage (V) was applied to the tip, and the sample was grounded. The nc-AFM images were recorded in constant-height mode based on a qPlus sensor design25. Several different nc-AFM tips are used in this work and the parameters for force spectra measurements in Fig. 3 are: n = 3, f0 = 39675 Hz, Q = 66390, A = 100 pm; n = 4, f0 = 40027 Hz, Q = 51415, A = 100 pm; n = 5, f0 = 39675 Hz, Q = 66390, A = 50 pm. The STM and nc-AFM images were processed by the Gwyddion software. To approximately remove the van der Waals contribution of the substrate from the force spectroscopy, Δf(z) spectrum recorded on the bare Ir(111) surface was subtracted from those taken on h-BN nanoplates. The short-range force spectra F(z) were calculated by a deconvolution process of the corresponding Δf(z) using Sader-Jarvis method26.

DFT calculations

Spin-polarized density functional theory (DFT) calculations were performed by using the plane wave-based Vienna ab initio simulation package (VASP)27,28. The generalized gradient approximation (GGA) method with Perdew-Burke-Ernzerhof (PBE) functional was employed to describe the exchange-correlation interaction among electrons29. The van der Waals (vdW) correction with the Grimme approach (DFT-D2) was included in the interaction between h-BN nanosheet and substrates. The energy cutoff for the plane wave-basis expansion was set to 520 eV and the atomic relaxation was continued until the forces acting on atoms were smaller than 0.02 eV Å−1. Ir(111) surface was modeled using an 8 × 8 slab with three Ir atomic layers, of which the bottom layer was fixed and separated by 20 Å of vacuum. The Brillouin zone was sampled with 3 × 3 × 1 Gamma-center k-point mesh, and the electronic states were smeared using the Fermi scheme with a broadening width of 0.2 eV. To estimate the interactions between tip and h-BN/Ir(111) nanosheet, we modeled the tip by a 14-atom W(100) pyramid. Simulated STM topography is based on the Tersoff-Hamann tunneling model with a height of 3 Å above the h-BN/Ir(111) nanosheet for states between 0 and 0.2 eV below the Fermi level.

The formation energy (Eƒ) of h-BN nanosheet on Ir(111) is defined as the following equation:

$${E}_{f}={E}_{h-{{{{{\rm{BN}}}}}}/{{{{{\rm{Ir}}}}}}(111)}-{E}_{{{{{{\rm{Ir}}}}}}\left(111\right)}-{{{{{\rm{n}}}}}}{E}_{{{{{{\rm{N}}}}}}}-{{{{{\rm{m}}}}}}{E}_{{{{{{\rm{B}}}}}}}$$

where \({E}_{h-{{{{{\rm{BN}}}}}}/{{{{{\rm{Ir}}}}}}(111)}\) and \({E}_{{{{{{\rm{Ir}}}}}}\left(111\right)}\) are the total energies of h-BN/Ir(111) nanosheet and Ir(111) slab. \({E}_{{{{{{\rm{N}}}}}}}\) and \({E}_{{{{{{\rm{B}}}}}}}\) represent the energies of N and B atoms in N2H4 molecular and BN monolayer, respectively, as shown in the following equations:

$${E}_{{{{{{\rm{N}}}}}}}={E}_{{{{{{{\rm{N}}}}}}}_{2}{{{{{{\rm{H}}}}}}}_{4},({{{{{\rm{g}}}}}})}-2{E}_{{{{{{{\rm{H}}}}}}}_{2},({{{{{\rm{g}}}}}})}$$
$${E}_{{{{{{\rm{B}}}}}}}={E}_{{{{{{\rm{BN}}}}}},{{{{{\rm{ML}}}}}}}-{E}_{{{{{{\rm{N}}}}}}}$$

where \({E}_{{{{{{{\rm{N}}}}}}}_{2}{{{{{{\rm{H}}}}}}}_{4}}\) and \({E}_{{{{{{{\rm{H}}}}}}}_{2}}\) are the total energies of N2H4 and H2 gas molecules. \({E}_{{{{{{\rm{BN}}}}}},{{{{{\rm{ML}}}}}}}\) is the total energy of BN monolayer. n and m represent the numbers of N and B atoms in h-BN nanosheet.

To obtain the force spectra acting on the B/N sites, we calculated the total energies of the W-tip and h-BN nanoplates on Ir(111) for different distances z. The interaction energies Eint were then acquired by subtracting the total energies of the W-tip and the h-BN nanoplates adsorbed on three layers Ir-atoms from the total energy of the combined system. These interaction energies were fitted by a Morse potential by the following equation30:

$$V\left(z\right)={V}_{0}\left\{{\left[1-\exp (-2b\frac{z-c}{c})\right]}^{2}-1\right\}+d$$

V(z) is the interaction energy as a function of the “tip-surface” distance z, which is defined as the distance before relaxation. Parameters b, c, and V0 are related to the strength and extent of the two-body chemical interaction. d is a constant related to the zero-point energy and can be chosen arbitrarily as it is irrelevant to the calculation of force spectra. Thus, the force gradient can be obtained analytically from the derivative of the fitted Morse potential of the above Eq. The force-distance curve \(F\left(z\right)\) was extracted by differentiating the Morse fit function with respect to z:

$$F\left(z\right)=-\frac{\partial {E\left(z\right)}_{{{{{\mathrm{int}}}}}}}{\partial z}$$

The frequency shift-distance curve Δf(z) shown in Fig. 4a is extracted by differentiating the force curve with respect to z:

$$\varDelta f=-\frac{{f}_{0}}{2{k}_{0}}\frac{\partial F\left(z\right)}{\partial z}$$