Stochastic mechanism of energy dissipation in noncontact atomic force microscopy studied using molecular dynamics with Langevin boundary conditions T. Trevethan and L. Kantorovich Department of Physics, Kings College London, Strand, London WC2R 2LS, United Kingdom (Received 4 February 2004 revised manuscript received 11 May 2004 published 14 September 2004) Based on the stochastic friction force theory of energy dissipation in non-contact atomic force microscopy (NC-AFM), we performed realistic molecular dynamics (MD) simulations on the MgO s001d surface to complement previous studies which only considered low frequency phonons in detail. We employed and calibrated Langevin boundary conditions to reduce effects due to the finite system size and thus to mimic an infinite lattice. The calculated dissipation energies are many orders of magnitude smaller than those observed experimentally and are similar to those calculated previously using simple analytical models. These findings suggest that this mechanism is not responsible for the observed energy dissipation. DOI: 10.1103/PhysRevB.70.115411 PACS number(s): 68.37.Ps I. INTRODUCTION Non-contact atomic force microscopy (NC-AFM) is de- veloping into an important tool in surface and nanoscience, which has now achieved atomic resolution images of many insulating and conducting surfaces.1 In NC-AFM operating in frequency modulation (FM) mode, a tip with an atomi- cally sharp nanoasperity attached to a cantilever is oscillated at its resonant frequency above a surface and scanned later- ally, with a feedback mechanism regulating the cantilever height so as to maintain a constant frequency shift.2 In addi- tion to the frequency signal, the excitation signal which is used to maintain the constant amplitude of the cantilever oscillations, can also be recorded during the scan. This exci- tation signal which is directly related to the damping of the oscillating cantilever also exhibits strong atomic periodicity, suggesting that energy is lost to the surface via the chemical interaction of the nanotip with atoms in the surface.3 These so-called dissipation images are however less well under- stood than the constant frequency shift topography images, and the mechanism of their formation is still the subject of debate, since in the non-contact regime the tip and surface are usually considered to remain distinctly separated throughout the oscillation cycle. Several mechanisms have been proposed to explain atomic resolution dissipation images and these are reviewed in Ref. 4. Only two mechanisms however are considered likely to be responsible for the observed dissipation: adhe- sion hysteresis and the stochastic friction force mechanism. The adhesion hysteresis mechanism relies on a strong and reversible structural change in the surface and/or tip leading to the tip following different potential energy surfaces on approach and retraction. We have recently shown5 how re- versible tip induced surface instabilities can lead via this mechanism to dissipated energies comparable to those ob- served experimentally, however, the theory of this mecha- nism is not yet complete. In the stochastic friction force mechanism (or the ���Brown- ian motion��� mechanism), the tip experiences a frictional force due to the thermal fluctuations of the atoms in the surface and tip, similar to that of a massive Brownian par- ticle immersed in a fluid of much lighter particles. This mechanism was first proposed in Ref. 6 and its theory was later developed using classical nonequilibrium statistical theory7���9 and then mixed classical and quantum ensembles 10���12 it is now considered to be complete. In this theory, the tip, at a given position Q above the surface, ex- periences a frictional force proportional to its velocity and acting in the opposite direction, -jsQdQ. �� The friction coef- ficient, jsQd, arises due to the random fluctuations of the force between the surface and tip caused by the thermal mo- tion of the atoms in the surface and tip. According to the theory, the friction coefficient at a given tip position Q above the surface is proportional to the time integral of the en- semble average autocorrelation function of the vertical force fluctuations (about the average tip force), DFstd, acting on the tip at that point: jsQd = 1 kBT E��� 0 kDFs0d D Fstdleqdt. s1d Several studies have attempted to calculate the friction coefficient based on this theory using different methods, and in all cases the dissipated energies predicted are at least four orders of magnitude smaller than those observed experimen- tally above plane crystal surfaces in the non-contact regime. In Refs. 8 and 12 the friction is calculated analytically for the basic model of a single atom tip interacting with a single atom in the surface within the Debye model. Although for- mally the entire phonon frequency spectrum (up to the De- bye frequency vD) was taken into account, this model treat- ment is only exact near v =0 when the phonon density of states is correct. In Ref. 13 the friction is calculated for a realistic system, with a multiatom tip interacting with several atoms in the surface, however again only phonons near v =0 were considered. The conclusion drawn from these cal- culations, based on the dissipated energies predicted, is that this model is unable to account for the experimental obser- vations and therefore atomic scale dissipation in NC-AFM is not due to this mechanism. To validate this conclusion, the contribution of all phonons to the calculated friction coeffi- cient and thus the dissipation energy within a realistic atom- PHYSICAL REVIEW B 70, 115411 (2004) 1098-0121/2004/70(11)/115411(7)/$22.50 ��2004 The American Physical Society 70 115411-1

istic model of the tip���surface junction is necessary. An important feature of the expression (1) for the friction coefficient is that the tip force autocorrelation function is an equilibrium ensemble averaged property for a given tip po- sition above the surface which can be evaluated via an equi- librium classical molecular dynamics (MD) simulation. By performing a series of equilibrium MD simulations at many different tip positions Q above a surface, the friction coeffi- cient as a function of tip position, jsQd, can be determined within a realistic atomistic model of the system. Using this function and typical experimental parameters of tip oscilla- tion, it is possible to calculate the dissipated energy per tip oscillation cycle as a function of tip closest approach and lateral position. These calculations are the main subject of this article. As will be clear from the following, our MD simulations are in fact complimentary to the calculations per- formed previously.8,13 Note that we assume here that atomic instabilities do not occur at the close approach in other words, we are con- cerned only with the friction in the noncontact regime when local ���soft��� modes are not present. Of course, if the tip does come sufficiently close to the surface, the nonequilibrium theory14 predicts that the friction will only be caused by ���fast��� phonons, whereas ���soft��� phonons are responsible for the adhesion effects due to atomic instabilities. However, treatment of dissipation due to the adhesion hysteresis mechanism is not the subject of the present paper. The plan of the remainder of this paper is as follows: the next section outlines the details of the simulation and the boundary procedure within the MD simulations to eliminate finite size effects the third section details the calculation of the friction coefficient and dissipated energies, which are then compared with experimental values and the analytical models mentioned above the final section consists of a dis- cussion of our results and conclusions. II. SIMULATION MODEL Since the friction coefficient (1) is expressed via the time integral of the fluctuating tip-surface force autocorrelation function, kDFs0d D Fstdleq, we start this section by discussing the calculation of this integral. We shall discuss how this issue is intrinsically related to the boundary conditions im- posed on the finite system of the simulation. A. Treatment of boundary conditions It is well known that particular care must be taken when evaluating the integrals of correlation functions of thermally fluctuating variables in finite systems15 (for a more general discussion of this, so-called ���plateau problem,��� see Refs. 16 and 17), and that these integrals only yield converged results when the thermodynamic limit is taken before the infinite time limit in the integral this is a major issue in systems typically composing of N . 103 ���104 atoms. This means that if the time integral (1) is calculated over a finite time t (i.e., when the upper infinite limit is replaced by t), j does not converge to a finite nonzero value if the limit t������ is taken before the N ������ thermodynamic limit. In other words, the value of j depends on the choice of the integration limit t and on the system size. This type of problem is encountered in studies which have attempted to evaluate the friction coefficient of a Brownian particle in a fluid using MD simulations.15,18,19 In Ref. 19 the equations of motion of particles are modified so that the momentum fluctuations of the system, and hence dynamical correlation functions become nearly independent of system size. It was possible to achieve this by modifying the veloc- ity distribution of the mobile fluid particles crossing the simulation cell boundaries, with a randomly generated distri- bution function, which leads to the Brownian particle force autocorrelation function integral converging to the physical value, which is nearly independent of the system size. The main idea of the boundary conditions discussed above is in randomizing the movement of boundary atoms. Although the idea itself seems to be very attractive, it cannot be implemented directly in the present investigation of a solid in which atoms vibrate about their equilibrium lattice positions in addition to that, our systems do not possess any periodic symmetry. In the present case another approach, known as the Generalized Langevin Equation (GLE) or the ���ghost atom��� approach,20,21 can be applied. In this method, which was first used to model atom���surface scattering, the system is divided into a finite primary region 1, and the rest of the crystal which is named the secondary region 2. Only the motion of the atoms in the primary region is considered explicitly the effect of the atomic vibrations in the infinite secondary region is introduced via additional forces acting on all atoms in the primary region. There are two types of forces (to be referred to as Langevin terms) to be added to the equations of motion of the atoms in the primary region to mimic the effect of an infinite solid: a generalized friction to represent the dissipation of energy from the primary region to the rest of the lattice, and a fluctuating (stochastic) force to represent the energy imparted to the primary zone due to the thermal vibrations of the lattice. The dynamics of the pri- mary region can then be recovered through the correct choice of the forms of the generalized friction and fluctuating force.21 Thus, the equation of motion for every atom in the pri- mary region resembles the Langevin equation with memory.22 For the degree of freedom i it can be written as follows: ] pi ] t = Fi - Et 0 o jP1 Gijst - t8dpjst8ddt8 + f i std, i P 1, s2d where pi is the classical momentum, Fistd is the force due to interaction with other degrees of freedom, Gijstd is the kernel of the generalized friction term, and f i std is the fluctuating (random) force. The exact forms of the generalized friction and fluctuating force are as intractable as the infinite set of equations that describe the atomic vibrations in the infinite secondary re- gion. However, it is possible to make several approximations to allow simplified calculation and acceptably preserve the effects of the secondary region atoms.20,21,23 First, by disre- garding the memory of the friction (assuming the limit of T. TREVETHAN AND L. KANTOROVICH PHYSICAL REVIEW B 70, 115411 (2004) 115411-2