LETTER pubs.acs.org/NanoLett
On Mapping Subangstrom Electron Clouds with Force Microscopy C. Alan Wright and Santiago D. Solares* Department of Mechanical Engineering, University of Maryland, College Park, Maryland 20742, United States
bS Supporting Information ABSTRACT: In 2004 Hembacher et al. (Science 2004, 305, 380383) reported simultaneous higher-harmonics atomic force mocroscopy (AFM)/scanning tunneling microscopy (STM) images acquired while scanning a graphite surface with a tungsten tip. They interpreted the observed subatomic features in the AFM images as the signature of lobes of increased electron density at the tungsten tip apex. Although these intriguing images have stirred controversy, an in-depth theoretical feasibility study has not yet been produced. Here we report on the development of a method for simulating higher harmonics AFM images and its application to the same system. Our calculations suggest that four lobes of increased electron density are expected to be present at a W(001) tip apex atom and that the corresponding higher harmonics AFM images of graphite can exhibit 4-fold symmetry features. Despite these promising results, open questions remain since the calculated amplitudes of the higher harmonics generated by the short-range forces are on the order of hundredths of picometers, leading to very small corrugations in the theoretical images. Additionally, the complex, intermittent nature of the tipsample interaction, which causes constant readjustment of the tip and sample orbitals as the tip approaches and retracts from the surface, prevents a direct quantitative connection between the electron density and the AFM image features. KEYWORDS: Atomic force microscopy, subatomic imaging, ultrahigh vacuum, frequency modulation, noncontact AFM n 2004, Hembacher et al.1 reported high-resolution results of a graphite surface characterized with a tungsten tip using simultaneous atomic force microscopy (AFM) and scanning tunneling microscopy (STM). The AFM images were acquired by mapping the higher harmonics response of the oscillating cantilever and exhibited features that fit within the footprint of a single tungsten atom, which the authors interpreted as the signature of lobes of increased electron density at the tip apex. Although these results were not the first to claim subatomic resolution,2 they were the first to be acquired with higher harmonics AFM imaging and the first in which the probe and sample were of different chemical identity. The previous claim of subatomic resolution using AFM2 was based on an experiment which, according to the authors’ interpretation, involved the detachment of a small cluster of atoms from the surface and its attachment to the AFM tip.2,3 Furthermore, the sample was a Si(111) surface, whose adatoms are much farther apart than the atoms in a typical crystal, and the imaging method was standard frequency modulation (FMAFM).4 We focus here on the more challenging problem of simulating the characterization of the graphite/tungsten system using higher-harmonics AFM. Despite the controversy that has surrounded the images of Hembacher et al. and although it has been 7 years since their publication, an in-depth theoretical feasibility study has yet to be produced. Although this paper is neither aimed at settling the controversy nor evaluating the correctness of the experimental approach and results of Hembacher et al., it highlights promising results and lays down a theoretical foundation for the further study of high-resolution higher harmonics AFM imaging.
I
r 2011 American Chemical Society
The origin of the forces in AFM and the length scales over which they act generally necessitate the use of simulation to interpret atomic-scale experimental images. The most common procedure used to model ultrahigh vacuum (UHV) systems involves the use of density functional theory (DFT) to simulate a model AFM tip (usually considering only a few atoms of the apex)57 above a model sample surface in a three-dimensional grid of points.3,5,811 For a specific (x, y) position, the calculated energies (or forces) in the z direction are either numerically differentiated or fit to analytical functions, which are subsequently used to derive the short-range tipsample force curves necessary to simulate an atomic resolution image (if the DFT forces are used instead of the energies, then the analytical function itself becomes the tipsample force curve). This approach has been successful in qualitatively explaining experimental atomic resolution AFM images, in particular with regards to the mechanisms responsible for image contrast.5,11,12 It has also been used to explain what the structure of the tip apex atoms must have been in order to produce a given experimental image, since the small scale of the system makes it impossible to know the exact tip structure during the experiment.3,7,13 DFT simulations for AFM have become more prevalent as computers have become faster, codes more efficient, and exchange-correlation functionals more accurate. However, they have yet to be incorporated into the simulation of AFM images produced with Received: September 5, 2011 Revised: October 11, 2011 Published: October 18, 2011 5026
dx.doi.org/10.1021/nl2030773 | Nano Lett. 2011, 11, 5026–5033
Nano Letters
LETTER
higher harmonics techniques, which have recently emerged at the forefront of scanning probe microscopy.1,14,15 The dynamic nature of FM-AFM only permits the indirect measurement of the tipsample forces, Fts, through the cantilever’s frequency shift vs z position (or z position at constant frequency shift). To recover quantitative values for Fts, an inversion procedure must be applied post hoc.1621 Because this additional step is not experimentally ideal, various methods have been proposed to recover Fts in real time.14,15,20,22 In 2000, using a first-order perturbative approach, D€urig demonstrated that Fts can be reconstructed over the range of the cantilever oscillation by tracking the amplitudes and phases of the higher harmonics of the tip trajectory.22 For a cantilever oscillation given by ψðtÞ ¼
∞
∑ an cosðnωtÞ n¼1
ð1Þ
the relationship between the harmonic amplitudes and the tipsample force is given as an ¼
2 1 Z πk 1 n2
1 1
Tn ðuÞ Fts ðz0 þ a1 uÞpffiffiffiffiffiffiffiffiffiffiffiffi du 1 u2
ð2Þ
for n 6¼ 1, where an is the amplitude of the nth harmonic, k is the cantilever stiffness, and Tn(u) is the nth Chebychev polynomial of the first kind.22,23 It is important to note that this relationship is most useful for cantilever oscillations that are on the order of the short-range interaction. In this case the tipsample force reconstruction requires only a few harmonics because the contribution of each harmonic to the short-range force reconstruction is more significant. When the fundamental oscillation is much larger than the short-range interaction, the low-order harmonic amplitudes become proportional to the frequency shift within D€urig’s derivation, such that mapping those harmonics does not necessarily provide information on the short-range tipsample forces. The theoretical effective frequency (and thus the frequency shift) is readily obtained from the fundamental oscillation amplitude, a1, via eq 2 as sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 1 2 Z 1 u ω ¼ ω0 1 Fts ðz0 þ a1 ð1 þ uÞÞpffiffiffiffiffiffiffiffiffiffiffiffi du k a1 π 1 1 u2 ð3Þ By integrating eq 2 by parts n times Hembacher et al. showed that, an ¼
Z 2 1 An 2 πk 1 n 1 3 3 3 5 3 ::: 3 ð2n 1Þ
dn Fts ðz0 þ a1 uÞ ð1 u2 Þn 1=2 du 1 dzn 1
ð4Þ demonstrating that an correlates to a convolution of dn/dzn[Fts(z0 + a1u)] with a bell-shaped weight function (1 u2)n1/2.1,23 The advantage of this development lies in the fact that higher force gradients have a stronger dependence on tipsample distance.24 Hence, contributions to the higher force derivatives should come primarily from the front-most tip atom (note that the above mathematical derivations require that Fts be continuous and differentiable). The coupling of the higher harmonics amplitudes to higher gradients of Fts makes the latter experimentally accessible and thus implies that mapping the harmonics should provide increased spatial resolution with respect to conventional FM-AFM.1,23,25
Figure 1. The simultaneously recorded tunneling current (left column) and higher harmonics (center column) images from Hembacher et al. (see Figure 2 in ref 1). Each row represents a separate experiment, with gentle tipsurface contacts between experiments. While the tunneling current images show a single maximum within the diameter of a tungsten atom, the higher harmonics images display multiple maxima within the same diameter. The authors interpret these subatomic features as the footprint of bonding lobes of increased electron density at the foremost tungsten tip atom. The right column shows the WignerSeitz cells of body-centered cubic tungsten parallel to [110], [111], and [001], exhibiting 2-, 3-, and 4-fold symmetry, respectively (from top to bottom). Reprinted with permission from ref 1. Copyright 2004 AAAS.
Hembacher et al. performed their 2004 experiments using constant height simultaneous STM/FM-AFM to image the surface with a tungsten tip attached to a qPlus sensor with separate electrodes for measuring tunneling current and cantilever deflection.1,2 The high stiffness of the sensor (∼1800 N/m) enabled a stable 3 Å oscillation amplitude, and its piezoelectric nature provided enhanced sensitivity in the cantilever deflection signal.23 Imaging was conducted in cryogenic UHV, and graphite was chosen as the surface material in order to exploit the small size of carbon as a probe atom—in these experiments, the roles of tip and sample are viewed as reversed, i.e., the sample atoms are imaging the tip according to the reciprocity principle posited by Chen.26 The higher harmonics were collected as a root-meansquare sum after passing the deflection signal through an eighthorder elliptical high-pass filter to eliminate the fundamental frequency (see Supporting Information). The justification for collecting the harmonics this way is given in ref 23 as well as in note 22 in ref 1. The essence of the argument is that (simulated) images produced by individual higher harmonics are similar, so collecting an rms sum of all higher harmonics will not change the nature of the image. Furthermore, this method offers ease of implementation and increases the signal-to-noise ratio. The simultaneously recorded tunneling current and higher harmonics images of Hembacher et al. are shown in Figure 1 with 5027
dx.doi.org/10.1021/nl2030773 |Nano Lett. 2011, 11, 5026–5033
Nano Letters
Figure 2. The three-layer W(001) tip structure (A) and graphite top layer (B) studied in this work. The model surface of graphite consisted of a nine-carbon-ring top layer edge-terminated with either oxygen (red) or hydrogen (white) depending on which was needed to maintain the sp2 hybridization of the carbon atoms. This model also contains a fourcarbon-ring second layer which has been excluded here for clarity. The second layer is included in order to reproduce the difference in electronic states between α- and β-surface carbon atoms. The red circles indicate the grid points for which tipsample force curves were calculated.
each row of images corresponding to a separate experiment (gentle tipsample contacts were made between experiments). The tunneling current images (panels A, C, and E of Figure 1) for each experiment exhibit a single maximum within the diameter of a tungsten atom (∼274 pm), while the higher harmonics images (panels B, D, and F of Figure 1) show multiple maxima within the same diameter. As already stated, the authors interpret these subatomic features as the signature of the electron density at the foremost tungsten tip atom. While both tunneling current and higher harmonics amplitude should exhibit steep distance dependence, tunneling current is proportional to the density of states (DOS) of the tipsample system at or above the Fermi level (as well as the voltage applied between tip and sample),27 so only the most loosely bound electrons will tunnel. For this reason, the STM images cannot resolve the spatial variations in the total charge density.1,28 In fact, it is known that STM can only image the β atoms of the graphite lattice because the α atoms do not have a sufficiently high density of states (DOS) at the Fermi level.28 Hembacher et al. attribute the different symmetries of the subatomic features found in each experiment (see panels B, D, and F of Figure 1) to the crystallographic orientation of the tungsten atoms at the tip apex. WignerSeitz cells of bodycentered cubic tungsten indicate that 2-fold symmetry corresponds to the [110] direction, 3-fold to [111], and 4-fold to [001], as also illustrated in Figure 1 (right column). For the case of [001], the existence of four lobes of increased charge density has previously been demonstrated in periodic surface systems by plane-wave quantum mechanics calculations of the W(001) surface.29,30 These lobes are thought to be the result of covalent-like bonding in the bulk, which is a common feature of transition metals with partially occupied d-shells.1,31 To summarize, the important features of the 2004 experiments that must be considered in a theoretical treatment are as
LETTER
Figure 3. (A) A single isosurface (isovalue = +0.2 e/Å3) of the change in electron density for the model tip structure used. The view is at an angle from below. The apex atom is indicated by a red square. Four lobes of increased charge density are visible. (B) Comparison of the total density (F) to the change in density (dF) after a self-consistent field calculation for the tip. The F = 2.0 e/Å3 isosurface of the total density is shown in gray, along with the dF = +0.2 e/Å3 isosurface shown in orange.
follows: (1) the crystallographic plane of tungsten responsible for image contrast and the presence of bulklike surface states in a nonbulk system such as an AFM tip, (2) the effect of simultaneous STM/AFM on the electronic states of the system, (3) the cantilever dynamics, and (4) the effect of filtering the cantilever trajectory. Each of these features can be the subject of a comprehensive study in itself. We focus here on the first of them, with the specific objective of reproducing the 4-fold symmetry features of Figure 1F. Ultimately, to study the effect of the cantilever dynamics and control systems on the imaging process, continuum simulation (and thus multiscale modeling) is necessary. However, the feasibility of the interpretation must first be tested in the most ideal imaging scenario. To this end, we can avoid the additional uncertainty of the cantilever dynamics and the filtering of a simulated cantilever trajectory by calculating the frequency shift and amplitudes of the higher harmonics directly from the fundamental theory of D€urig (see eq 2), which by the argument of Hembacher et al. should reproduce the subatomic contrast seen experimentally. The advantage of calculating the purely theoretical frequency shift and higher harmonics amplitudes this way is 2-fold. First, as mentioned above, it allows one to study the most ideal imaging scenario in which the higher harmonics amplitudes can be calculated aside from any artifacts that could be introduced through the experimental processing of the cantilever deflection signal (one of the major points of critique that has been brought up against the work of Hembacher et al. is the hypothesis that the 4-fold symmetry features could be the result of filtering artifacts). Consequently, if the subatomic contrast cannot be reproduced in this manner, then one can readily conclude that there is no true connection between the electron density at the apex tip atom and the subatomic features of the higher harmonics images. Second, if the higher harmonics images do reflect the apex atom electron density, then the calculated ideal amplitudes establish a baseline for evaluating the effects of the cantilever dynamics and filtering which would be introduced in a full multiscale simulation. 5028
dx.doi.org/10.1021/nl2030773 |Nano Lett. 2011, 11, 5026–5033
Nano Letters The W(001) tip and graphite surface models used are shown in Figure 2. The DFT code SeqQuest developed by Shultz at Sandia National Laboratory32 was used to run all calculations (the computational methods are described in detail in the Supporting Information). Recalling that quantum mechanics calculations of the W(001) surface reveal four lobes of increased charge density above the top layer atoms, the first question we wish to answer is, Are similar lobes present in nonbulk systems such as AFM tips? Plotted in Figure 3A is a single isosurface (isovalue = +0.2 e/Å3) of the change in electron density, dF, for the tip structure investigated.33 Four lobes of increased charge density are indeed visible at the apex atom of the structure (indicated by the red square). While the presence of the four lobes of increased charge density in a nonbulk tip structure is a promising result, it is important to note that the value of the density increase is small in comparison to the total electron density in the same region of space. Figure 3B illustrates the F = 2.0 e/Å3 isosurface of the total electron density at the apex tip atom (gray), along with the dF = +0.2 e/Å3 isosurface (orange), which is an order of magnitude smaller. Note that only single isosurfaces are plotted in Figure 3. For this specific case, density changes in the interior of the tip structure are not visible because none occur at the selected isovalue (dF = +0.2 e/Å3). Plotted in Figure 4 are two cross sections through the entire change in electron density of the model tip structure, showing all isovalues. Figure 4B indicates that four lobes are also present below the center atom in the top layer of the tip, which is to be expected since this is one of the most bulklike atoms in the system. In order to calculate the force curves necessary to simulate the frequency shift and higher harmonics images, the tip model was placed above the surface at each red grid point shown in Figure 2B. Data points were calculated from 6.00 Å down to 1.00 Å (in steps of 0.25 Å) above the surface, with the z distance defined as the distance between the centers of the apex tungsten tip atom and the surface atoms, prior to relaxation. Over each surface grid point, the tungsten tip began 6.00 Å above the surface. The system was allowed to fully relax, and then the tip atoms were shifted downward by 0.25 Å before the system was allowed to relax again.34 This process resulted in 100 relaxed energy curves each consisting of 21 data points, an example of which is shown in Figure 5A, which is the curve calculated above the solid red grid circle of Figure 2B (directly above one of the β carbon atoms). In order to generate force curves, each energy curve was fit to a function that was then analytically differentiated with respect to z: F(z) = 3E(z). While this seems like a
LETTER
straightforward procedure, finding the optimum function for curve fitting is a critical step. An ill-fitting function undermines the precision of the DFT simulations, which are at the core of this method. In this particular case, the standard three-parameter Morse potential is one such ill-fitting function. One can alter the Morse potential with correction terms to achieve a more accurate fit to the data, but these modified functions either do not offer a satisfactory improvement (e.g., the four-parameter Levine potential35) or offer an improved fit but do not have wellbehaved slopes in regions away from the 21 data points.36 This is uniquely relevant to higher harmonics AFM simulation because we are interested in not only the first energy gradient (to obtain the force) but also the higher energy gradients that
Figure 5. (A) Example final relaxed energy (FRE) (attojoules) versus distance curve. Such curves were collected on all locations (red circles) shown in Figure 2B. Each point on the graph is the result of a separate simulation. (B) Example force curve developed through DFT and curve fitting (red trace) and the reconstruction of the force curve using six harmonics in D€urig’s equations (blue dotted trace). These curves correspond to the location marked with the solid red circle in Figure 2B (directly above one of the β carbon atoms).
Figure 4. (A) Side view of the change in electron density (dF) for the model tip structure, indicating the location of a cross section taken directly below the topmost tip layer (B) and directly below the apex atom (C). Cross section B reveals four bonding lobes below the top layer’s center atom and cross section C reveals four smaller lobes below the apex tip atom, with all isovalues of dF included. 5029
dx.doi.org/10.1021/nl2030773 |Nano Lett. 2011, 11, 5026–5033
Nano Letters
LETTER
Figure 6. Simulated frequency shift (A) and higher harmonics, Vhh (B), images calculated using D€urig’s equations. The grid of Figure 2B has been repeated multiple times in both the x and y directions and then cropped to produce rectangular images. Here the distance of closest approach is 1.88 Å. A red hexagon has been overlaid to indicate the positions of the surface carbon atoms. While the graphite lattice is clearly visible in the frequency shift image, 4-fold symmetry features similar to those observed experimentally appear in the simulated higher harmonics image. (C, D) A close-up side-byside comparison of the experimental (C) and simulated (D) Vhh images revealing 4-fold symmetry features.
couple to the higher harmonics of the cantilever motion. The issue of curve fitting the DFT data is of less importance in the simulation of standard frequency shift images, because slight deviations from the data do not result in drastic qualitative differences between simulated images but instead only alter the quantitative frequency shift values. To ensure continuous, realistic force curves, we used an in-depth, two-stage fitting method which is discussed in detail in the Supporting Information but summarized here. First, the DFT energy data were accurately fit to ninth order polynomials. Because such functions still have illbehaved slopes in regions outside of the data points, the force curves produced through differentiation of this initial fit were not realistic outside of the region of the interaction well. However, one can take new data points from this initial force curve in the accurately represented interaction well and fit them to a standard force curve model. In essence, this second stage of the process saves the region of importance (the interaction well) and truncates the surrounding regions, replacing them with realistic force model behavior. For this second stage fit we used a Levinetype force function.36 Finally, before simulating the frequency shift and higher harmonics images, we determine the number of
harmonics necessary to reconstruct the tipsample force with the 3 Å oscillation amplitude used in the experiment. Figure 5B shows the force curve calculated from the energy data of Figure 5A (red) as well as the reconstruction (dotted blue) calculated from D€urig’s equations using the first six harmonics of the 3 Å cantilever oscillation. The three vertical grid lines represent the cantilever base position (chosen to be 5 Å) and the two end points of the oscillation. Clearly the force curve is well reproduced over the range of oscillation. We now proceed from D€urig’s equations by simulating experimental frequency shift and higher harmonics images. We begin by focusing on the short-range force curves developed through the two-stage fitting method described above, disregarding the long-range van der Waals forces. Figure 6 shows simulated frequency shift (A) and higher harmonics, Vhh (B), images calculated using a 3 Å fundamental oscillation amplitude, with the distance of closest tipsample approach set at 1.88 Å (simulated images at varying closest approach distances are provided in the Supporting Information for not only the tip discussed above but also for a “blunt” tip, whose apex is terminated with four atoms). The gray color map was chosen 5030
dx.doi.org/10.1021/nl2030773 |Nano Lett. 2011, 11, 5026–5033
Nano Letters
LETTER
Figure 7. Illustration of the force cross sections taken relative to (A) the frequency shift and (B) the higher harmonics images of Figure 6. The force corrugations through the carbon ring centers (C) reveal both hole sites and bonds, with the hole sites exhibiting the maximum attractive force. The force cross section through the surface β atoms (D) shows corrugations consistent with the two lobes above each atom.
simply because it provided the best visual contrast. A red hexagon has been overlaid to indicate the positions of the surface carbon atoms. To produce the higher harmonics image, the amplitudes of the first seven harmonics (based on the force reconstruction exercise of Figure 5) as calculated from D€urig’s theory (eq 2) were entered into Hembacher’s equations S-5 and S-6 (found in the Supporting Information).36 While the graphite lattice is clearly visible in the simulated frequency shift image, four-lobed features similar to those observed experimentally appear in the simulated higher harmonics image (see also Figure S-17 in the Supporting Information, which was constructed using a blunt tip and in which the number of features in the higher harmonics image appears to multiply with the number of apex atoms). A direct comparison with the 4-fold experimental higher harmonics image is instructive. Figure 6 also shows close-up views of both the experimental (C) and simulated (D) Vhh images revealing four-lobed features. In both cases, the features appear within a 200 pm diameter. The white cross in the experimental image was placed at a maximum in the tunneling current image, presumed to occur directly above a surface carbon atom.1 In a similar fashion, we have placed a white cross at the exact position of the nearest surface carbon atom in our simulation. Clearly the simulated features are not centered above a single carbon atom. From Figure 6B we see that the brightest of the features actually occurs at a hole position on the surface. The fact that the simulated features are not centered above a carbon atom is not necessarily problematic, but it suggests that the
assumption by Hembacher et al. that the carbon atom position is the same in both tunneling current and higher harmonics images is not guaranteed to be accurate. The experimental observation of offsets between the maxima of simultaneously acquired data channels is not uncommon.37 Note that a comparison between the simultaneously acquired tunneling current and higher harmonics images of panels A and B of Figure 1 as well as panels C and D of Figure 1 reveal that even in these instances the subatomic features are not centered about a maximum in the tunneling current. Another clear difference between the experimental and simulated higher harmonics images is the symmetry of the 4-fold features. While the lobes near the white cross on the experimental image (Figure 6C) appear to be approximately equal in weight, those in the simulated image do not. A possible explanation for this is the influence of the pinned hydrogen and oxygen atoms on the edge of the model surface. It is not unreasonable to conclude that the lack of exact bonding symmetry in the surface model would lead to a lack of symmetry in the calculated force curves and hence in the simulated higher harmonics images. Additionally, close inspection of the experimental image shows that the lobe symmetry is not uniform throughout the surface. Figure 7 shows vertical cross sections of the calculated tipsample forces, which offer further insight into our simulated images. Panels A and B of Figure 7 indicate the location of the cross sections relative to the simulated frequency shift and higher harmonics images of panels A and B of Figure 6, respectively. 5031
dx.doi.org/10.1021/nl2030773 |Nano Lett. 2011, 11, 5026–5033
Nano Letters The force map through the carbon ring centers (Figure 7C) reveals both hole sites (see white cross in Figure 7A and white dashed line in Figure 7C) and bond sites (see blue cross in Figure 7A and blue dashed line in Figure 7C). The hole sites exhibit the greatest attractive force, similar to previous experimental studies using other types of tips (see for example ref 37). Note that the force corrugations are also consistent with the lobe features in Figure 7B. For example, the force cross section through the surface β atoms (Figure 7D) shows two lobes above each atom. The horizontal variation in the forces through the regions resembling the bonding lobes is within 1 nN. The similarity between the images in panels C and D of Figure 6 is a promising result that supports the concept of mapping subatomic features via higher harmonics AFM under ideal conditions. However, it is also important to qualify the meaning of the observed features. While it is reasonable to view them as the signature of the unpaired bonding lobes of the tip apex, they do not correspond to a static imprint of their shape. This can be realized if one considers that for a fundamental amplitude of 3 Å, the tip oscillation spans a total vertical distance of 6 Å, which is significantly greater than the size of the bonding lobes. In this case, it may be more proper to qualify the observed features as a dynamic signature representing the convolution of the cantilever dynamics and parameters with the intermittent short-range forces that emerge between the probe and sample orbitals. In fact, our calculations suggest that the observed subatomic features can change shape and even appear or disappear in a nonmonotonic fashion as a function of the cantilever height, which is not surprising if one considers the evolution of the lobes at the tip apex as the vertical tip position changes (see Supporting Information). Thus, we cannot yet establish a quantitative connection between the electron density of the probe and sample and the higher harmonics images, even under ideal conditions. Additionally, there remain other aspects of the method that require further study, such as the effect of the STM bias voltage, which is expected to induce changes in the geometry of the apex lobes by generating attractive forces between the orbitals of the tip and the sample; the filtering procedure that is used to separate the small higher harmonics signals from the fundamental oscillation (the calculated amplitudes of the first three higher harmonics under the conditions evaluated, in particular zero bias voltage and disregard of longrange attractive forces, were below 0.01% of the fundamental amplitude, in the range of hundredths of picometers38); and the effect of the tip layers directly above the apex atom, whose geometry has a direct effect on the bulklike nature of the apex atom and which our calculations show affects the size of the lobes (see Supporting Information). In summary, we have presented the results of DFT simulations which confirmed the existence of four lobes of increased charge density at the apex atom of a nonbulk W(001) tip, which had previously been observed in quantum mechanics calculations of the W(001) surface.29,30 These bonding lobes have been proposed to be responsible for experimentally observed subatomic, 4-fold symmetry features produced via higher harmonics AFM imaging.1 To test the feasibility of this claim, we first developed a method for simulating higher harmonics AFM imaging. Using the method, we showed that it is possible to obtain four-lobed features in a simulated higher harmonics image that are qualitatively similar to those observed experimentally. Important questions not addressed in this work still remain open, such as the challenges involved in filtering and processing subpicometer
LETTER
higher harmonics amplitudes and the effect of the bias voltage in simultaneous AFM/STM imaging, which could have a significant effect on the magnitude of the short-range tipsample forces and the corresponding harmonic amplitudes. Furthermore a direct quantitative connection between the images and the probe and sample electron density cannot yet be established due to the complex, intermittent nature of the tipsample interaction.
’ ASSOCIATED CONTENT
bS
Supporting Information. Detailed computational methods, additional simulated images for various tipsample distances, images simulated using a blunt tip structure, and details on the experimental harmonics detection. This material is available free of charge via the Internet at http://pubs.acs.org.
’ AUTHOR INFORMATION Corresponding Author
*E-mail:
[email protected].
’ ACKNOWLEDGMENT The authors gratefully acknowledge US NSF award No. CMMI-0841840. ’ REFERENCES (1) Hembacher, S.; Giessibl, F. J.; Mannhart, J. Science 2004, 305, 380–383. (2) Giessibl, F. J. Appl. Phys. Lett. 2000, 76, 1470–1472. (3) Huang, M.; Cuma, M.; Liu, F. Phys. Rev. Lett. 2003, 90, 256101. (4) Albrecht, T. R.; Grutter, P.; Horne, D.; Rugar, D. J. Appl. Phys. 1991, 69, 668–673. (5) Perez, R.; Payne, M. C.; Stich, I.; Terakura, K. Phys. Rev. Lett. 1997, 78, 678–681. (6) Caciuc, V.; Holscher, H.; Blugel, S.; Fuchs, H. Phys. Rev. B 2006, 74, 165318. (7) Pou, P.; Ghasemi, S. A.; Jelinek, P.; Lenosky, T.; Goedecker, S.; Perez, R. Nanotechnology 2009, 20, 264015. (8) Sasaki, N.; Aizawa, H.; Tsukada, M. Appl. Surf. Sci. 2000, 157, 367–372. (9) Dieska, P.; Stich, I.; Perez, R. Nanotechnology 2004, 15, S55. (10) Caciuc, V.; Holscher, H.; Blugel, S. Phys. Rev. B 2005, 72, 035423. (11) Caciuc, V.; Holscher, H. Nanotechnology 2009, 20, 264006. (12) Sugimoto, Y.; Pou, P.; Custance, O.; Jelinek, P.; Morita, S.; Perez, R.; Abe, M. Phys. Rev. B 2006, 73, 205329–9. (13) Giessibl, F. J.; Hembacher, S.; Bielefeldt, H.; Mannhart, J. Science 2000, 289, 422–425. (14) Sahin, O.; Magonov, S.; Su, C.; Quate, C. F.; Solgaard, O. Nat. Nanotechnol. 2007, 2, 507–514. (15) Solares, S. D.; Holscher, H. Nanotechnology 2010, 21, 075702. (16) Giessibl, F. J. Phys. Rev. B 1997, 56, 16010–16015. (17) Gotsmann, B.; Anczykowski, B.; Seidel, C.; Fuchs, H. Appl. Surf. Sci. 1999, 140, 314–319. (18) Giessibl, F. J. Appl. Phys. Lett. 2001, 78, 123–125. (19) D€urig, U. Appl. Phys. Lett. 1999, 75, 433–435. (20) D€urig, U. Appl. Phys. Lett. 2000, 76, 1203–1205. (21) Chawla, G.; Solares, S. D. Meas. Sci. Technol. 2009, 20, 015501. (22) D€urig, U. New J. Phys. 2000, 2, 5.1–5.12. (23) Giessibl, F. J. Surf. Interface Anal. 2006, 38, 1696–1701. (24) If one considers a force that varies as 1/zp, doubling the distance reduces the force by 1/2p, but the nth derivative of the force with respect to z will be 1/2p+n times as small. 5032
dx.doi.org/10.1021/nl2030773 |Nano Lett. 2011, 11, 5026–5033
Nano Letters
LETTER
(25) Giessibl, F. J. Mater. Today 2005, 8, 32–41. (26) Chen, C. J. Introduction to Scanning Tunneling Microscopy; Oxford University Press: New York, 1993. (27) Tersoff, J.; Hamann, D. R. Phys. Rev. B 1985, 31, 805. (28) Hembacher, S.; Giessibl, F. J.; Mannhart, J.; Quate, C. F. Proc. Natl. Acad. Sci. U.S.A. 2003, 100, 12539–12542. (29) Posternak, M.; Krakauer, H.; Freeman, A. J.; Koelling, D. D. Phys. Rev. B 1980, 21, 5601–5612. (30) Mattheiss, L. F.; Hamann, D. R. Phys. Rev. B 1984, 29, 5372–5381. (31) Pauling, L. The Nature of the Chemical Bond, 2nd ed.; Cornell University Press: Ithaca, NY, 1948. (32) Schultz, P. A., unpublished. http://dft.sandia.gov/Quest/. For a description of the method see: Feibelman, P. J. Phys. Rev. B 1987, 35, 26262646. (33) In DFT, the initial electron density F0 in the tipsample system is a superposition of the individual atom electron densities. During each self-consistent field (SCF) calculation for the system, the change in electron density dF is calculated. The total electron density for the system is then the sum of the initial density and the change. (34) Each SCF calculation was converged to within 0.0020 Ry (=0.027 eV). The geometry relaxation was governed by a force convergence criterion of 0.0010 Ry/bohr. During each relaxation, the top layer of the tungsten tip remained fixed, as did the oxygen and hydrogen atoms on the edges of the surface structure. (35) Levine, I. N. J. Chem. Phys. 1966, 45, 827–828. (36) See the Supporting Information for details. (37) Baykara, M. Z.; Schwendemann, T. C.; Altman, E. I.; Schwarz, U. D. Adv. Mater. 2010, 22, 2838–2853. (38) The first three harmonic amplitudes obtained for the force curve shown in Figure 4 are only 0.0049, 0.0098, and 0.0052% of the fundamental amplitude, respectively, which when taken as the only three harmonics in the Vhh sum (see eq S-5 in the Supporting Information) would correspond to 0.029% of the fundamental amplitude signal, after considering the piezoelectric sensitivity factor (eq S-6, Supporting Information).
5033
dx.doi.org/10.1021/nl2030773 |Nano Lett. 2011, 11, 5026–5033
