Dissipative Processes in Tribology / D. Dowson el al. (Editors) 0 1994 Elsevier Science B.V. AU rights reserved.
173
A molecularly-basedmodel of sliding friction J. L. Streator
G. W. Woodruff School of Mechanical Engineering Georgia Institute of Technology Atlanta, GA 30332-0405 Kinetic friction is well known to be a dissipative process. The key to understanding friction, then, is to discover the mechanisms by which mechanical energy is dissipated when one body slides upon another. Some recent reports in the literature have discussed the source of dissipation in the context of several models of atomic interaction, including the Frenkel-Kontorova(FK) and the Independent Oscillator (10) model. In the current work, aspects of the FK and I 0 models are incorporated to model the interaction between a solid asperity and an adsorbed interfacial film. As predicted by the I 0 and FK models, calculations indicate that the existence of an appreciable dissipative component of friction depends on the occurrence of an instability in the configurations of the molecules. It is also found that this dissipative component is influenced by whether or not there exist mechanisms for energy to continually propagate away from the interface. When such avenues exist, the friction force is found to be significantly greater than under conditions where the energy can return to the interface. To provide additional insight into the mechanism of dissipation, an analogy is made between the discrete, molecular model and that of a string on an elastic foundation. The conditions governing the propagation of waves with the continuous model suggest that significant dissipation requires a frictional instability. 1. INTRODUCTION Kinetic friction has long been recognized as a dissipative process, but the mechanism by which mechanical energy is lost remains a point of investigation [l]. It is known that the kinetic friction between two bodies results from interactions that occur on the atomic level. It follows, therefore, that the elucidation of frictional dissipation should result from consideration of atomic and/or molecular interactions. Probably the first atomic account of kinetic friction was due to Tomlinson [2]. Tomlinson proposed that when an atom of one body moves tangentially past an atom of another body, the lateral component of force may be different between approach and separation. The non-symmetry in this lateral force arises because, at a certain point, the atoms attain a position of unstable equilibrium. This condition leads to the sudden motion of the atoms toward the nearest stable equilibrium configuration. Tomlinson assumes that this "jump" phenomenon irreversibly transforms all of the stored energy into molecular kinetic energy or heat, thereby accounting for
mechanical losses. Thus the mechanism involves the nonadiabatic motion of atoms. Tomlinson's mechanism appears in the context of the Independent Oscillator (10) model [3, 41, shown schematically in Fig. l(a). Here the atoms of one surface are modeled by masses which are independent of one another but are harmonically coupled to a rigid base. The influence of the opposing surface is modeled by a sinusoidally varying force which translates quasistatically in a direction parallel to the interface. Although drawn vertically, both the interaction force and springs act tangential to the interface. Above a certain magnitude of interaction force the model atoms are displaced by the force in the direction of sliding until they reach a point of unstable equilibrium. Further translation of the upper surface causes the atoms to suddenly spring back. Assuming that this release of energy is irreversible, the instability in the I 0 model accounts for the Occurrence of friction. If, however, the interaction force is low, it can be shown that the molecules remain in positions of stable equilibrium throughout the sliding process. In this case the average tangential force vanishes so
174
L
k
/// Figure 1 : Previous models of dissipation. (a) Independent Oscillator (10) model. (b) Frenkel-Kontorova (FK) model. that the dissipcltive component of friction force is zero [3.4]. Similar effects can be seen in the FrenkelKontorova (FK) [S] model, shown schematically in figure I@). In this model, the atoms are modeled by a set of harmonically coupled masses. In the classical treatment, the system is analyzed for static configurations which lead to instability in the positions of the model atoms. It has been shown that this instability occurs at a critical value of the interaction force between the two surfaces. In the quasi-static analysis, the occurrence of this instability has heen identified as the source of friction. Since the masses of the surface are mutually coupled, the system behavior is influenced by their average spacing as compared to the wavelength of the translating force [ 6 ] , which is taken to be the lattice spacing of the upper surface. It has been found that when the ratio a/b is irrational, yielding incommensurate lattices, the friction force is considerably less than when a/b is a rational number [3.4].
In both the (quasistatic) FK and I 0 models. the Occurrence of friction depends on the existence of unstable equilibrium positions during the translation of one body over another which leads to sudden displacement or "plucking" of the atoms. In contrast. Hirano and Shinjo (1990) [7] performed calculations based on the Morse potentials for selected metals and found that, for these materials, such an instability should not occur. In a separate investigation [8] the authors analyzed a one-dimensional FK model with kinetic energy terms and found that friction could arise under dynamic conditions in cases where the quasistatic analysis would predict no instabilities of the atoms. The authors associated the Occurrence of friction with the transformation of macroscopic kinetic energy into internal motions of the body. Sokoloff (1992) [9] studied the interaction of a translating sinusoidal force with a layered crystalline lattice, modeled by a collection of springs and masses. In this investigation, the author emphasized the importance of internal damping on the generation of friction force for any finite-sized crystal. The Same author also compared frictional forces between commensurate and incommensurate sliding and found the former values to be 12 orders of magnitude greater than the latter [ 101. In the present work. we investigate a simple model of surface interaction which incorporates aspects of the both the 10 and FK models. I n contrast to previous studies, however, we consider directly the role that energy propagation plays i n influencing the kinetic friction. The motivation for this investigation comes from the observation that. when two macroscopic bodies slide. there is an appreciable amount of energy which propagates away from the interface in the form of acoustic waves. The question arises how this mechanism of energy loss affects the magnitude of friction. 2. MOLECULAR MODEL The physical system to be modeled is schematically shown in Fig. 2a. The upper surface is a hard asperity which slides against an adsorbed film of an opposing surface. The mathematical model is shown in Fig. 2b. Since the molecules of the film are anchored to the surface, the masses of
175
vo 0
(a)
*
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
into the contacting bodies, supporting structures. etc. 2.1 Governing Equation
With reference to Fig. 2b. the equation of motion for each of the masses of the model is written
with
Q(x, ,I ) = [H (X, - x, - V,I)- H ( X, - (X, + L ) - V,I )] x
Figure 2: Model of interface. (a) Asperity sliding on adsorbed film. (b) Mathematical model. the model (Fig. 2b) are connected to springs which resist absolute displacement relative to the base. The set of molecules are also assumed to prefer a certain equilibrium spacing which is modeled by the coupling between the masses. Since the upper surface is assumed to be that of a hard solid, it is modeled by a rigid sinusoidal force. All of the couplings are harmonic in nature. Following Shinjo and Hirano (1993) in their application of the FK model [8], we integrate the dynamical equations of motion for each of the masses. To directly investigate the role of energy propagation. the extent of the adsorbed film is taken either to be ( I ) just large enough to cover the contact length over the sliding distance, or (2) effectively unlimited in each direction. In practice, this latter condition means that the surface extent is sufficiently large that, over the duration of the sliding. elastic waves propagate outward and do not return. Physically, this feature is used to model what happens when real 3D bodies experience sliding Contact. In this case, acoustic waves generated at the Contacting interfaces propagate away from the interface region
where m is the mass, k and ke are the spring constants (Fig. 2b), ui is the displacement of the ith mass from its unforced static equilibrium position, xi is the absolute position of the ith mass, Q(xi.1) is the interaction force on the ith mass. Qo is the force amplitude, V, is the sliding velocity, xo is the initial position of the left end of the contact region. L is the length of the contact, a is the wavelength of the periodic force. The function Hc.>is the unit step function. From consideration of the arguments of the step function. it is seen that the interaction force translates from left to right with speed Vo. We can nondimensionalize the governing equation by selecting appropriate length. time and force scales. The parameter b, which gives the equilibrium spacing between the molecules when there is no interaction force (Fig. 2a) is chosen as the length scale, b/Vo is chosen as the time scale, and kb is chosen as the force scale. We also make the following definitions ,41
=by,
t=--'S b
v,
co = b
g
The parameter co has the dimensions of velocity and is nominally the velocity of sound for the medium [ 111. Using the above definitions we have
176
=i'
2.2 Simulation Procedure Equation (2) governs the motion of the ith particle. To study the frictional force we integrate equation (2) for all of the particles in the system using the classical 4th-order Runge-Kutta method [12]. The integration is begun from a position of static equilibrium, which itself is determined by iterative relaxation from an assumed configuration. Starting from this initial static equilibrium position, we sum the interaction force. Q. over all of the particles at each time step. This summation provides the friction force exerted by the translating force field upon the lower surface. The time step was selected to be 0.02 b/co As a check on computational accuracy. we compared the power input to the stored energy of the system. The rate of work performed on each particle is calculated by forming the product of the interaction force on the particle and the particle velocity. After summing over all of the particles in the system at each time step. the power was Integrated in time using Simpson's Rule to provide the total work done by the interaction force at a given instant. The total energy of the system was calculated by summing the kinetic energies of all of the system particles along with the potential energies stored in each of the springs. It was found that there was good agreement between the calculated work done and the energy stored in the system. At the end of the simulation. the relative difference in the two quantities was consistently found to be less than 0.1% .
3.0 RESULTS AND DISCUSSION The governing equation (2) was integrated for a number of parameter combinations. The following values are applicable to all of the calculations to be discussed: a/b = 1, Vdco= 0.02, k& = 0.8. The ratio a/b defines whether or not the interface is commensurate, and the ratio Vdco is the ratio of sliding speed to the nominal acoustic wave speed. With the value of this ratio set to
d
W
v3
y
d b rn I J
l1 -
OS3
-I
Qo/kb=0.05 L h - 5 Qo/kb = 0.1 Qo/kb=0.2
0.2 0.1
d
z 0 F
u
U
0.0 -0.1
-0.2
0
1
2
TIME (Vot/b)
3
Figure 3: Frictional stress vs. time for three force amplitudes. 0.02, the sliding process is slow compared with the speed of relaxation in the interface. 3.1 Effect of Interaction Force Figure 3 shows the (dimensionless) frictional stress vs. (dimensionless) time for three values of the interaction force amplitude and for a (dimensionless) contact length of L/b = 5 . The extent of the model surface film was taken to be unlimited in each direction. The frictional stress is per defined to be the friction force, F, dimensionless contact length L/b, normalized by kb. Thus frictional stress has the form F/kL. Since there is one particle per b spacing. the frictional stress gives a measure of the force per particle in the conlact region. Similarly, one unit of dimensionless time corresponds to the force field translating the equilibrium distance. b, between the molecules. As observed in Fig. 3. For Q&b = 0.05, the frictional stress is nearly sinusoidal being slightly distorted due to some displacement of the masses in the direction of sliding. Moreover, the
177
n
d
w
*
-
Qo/kb=0.05 ~ b - 5 - Qo/kb=O.l 0.03 - - Q$b ~ 0 . 2
W
u g 0.02 zw
u
p
2
2
I
-
0.01 -
0.00
I
I
Figure 4: Kinetic energy within contact for three force amplitudes. amplitude of the frictional stress (0.05) is equal to the amplitude of the normalized interaction force. Since the force profile is nearly symmetrical with respect to the origin, the average frictional stress over the sliding distance is quite small (.OOOOl). When the normalized force amplitude, Q&b, is doubled to 0.1, the friction trace shows a much greater deviation from sinusoidal behavior. Once the frictional stress reaches its peak, there is a faster drop to the minimum value after which some high-frequency oscillations are observed. The average frictional stress remains small but increases to 0.001. When the interaction force is increased to 0.2 kb. there is a sudden drop in the friction force followed by substantial oscillations in the force. In addition, the average frictional stress increases to 0.03. The degree of excitation in the interface can be measured by the amount of kinetic energy generated. Fig. 4 shows the amount of kinetic energy generated in the contact as a function of time. The energy is normalized by kb2 and the dimensionless sliding length L/b. For the case of
Q@b = 0.2, substantial kinetic energy is produced after each drop in the force (see Fig. 3). When the interaction force is reduced to 0.1 kb, the amount of kinetic energy generated is more than order of magnitude less. When Q$b = 0.05. the amount of kinetic energy produced is negligible, as indicated by the horizontal line at the origin. The differences in kinetic energy observed in Fig. 4 compared with the frictional profiles of Fig. 3, indicate the role of frictional instabilities in exciting the interface.
3.2 Effect of Surface Extent Figure 5a shows the frictional stress vs. time for the case of Q&b = 0.5 and L/b = 0.5. Here the surface is unlimited in extent (in each direction). In other words, a sufficient number of surface particles is selected to insure that. during the duration of the simulation, no waves originating from the interface can return to the interface upon reflection from the surface boundaries. As observed in the figure. the friction force demonstrates the plucking instability identified in Fig. 3. Again there is a sudden drop in the friction force followed by substantial oscillations. Fig. 5b shows the (normalized) energy delivered to the surroundings-that is the amount of energy, potential plus kinetic. possessed by the particles ourside of the contact region. The energy is normalized in such a manner that the graph provides the amount of energy lost per particle in the contact region. The figure shows that the energy in the surroundings accumulates with time. As the contact region translates to the right (see Fig. 2b), the particle at the left edge of the contact eventually leaves the contact region while, at the same time, a particle enters the contact from the right. If the particle leaving the contact has been excited by plucking. then its departure will generally result in a greater loss of energy from the contact than is gained from the particle which enters. The key observation here is that once the energy is lost from the contact, it does not return: the energy accumulated outside of the contact increases almost monotonically in time. This steady increase in the energy of the surroundings indicates that energy is continually propagated away from the interface. We compare the foregoing results to those of Fig. 6. In Fig. 6a, frictional stress is computed for
178 the same conditions as those in Fig. Sa. except thal thc surface extent is limited to 10 particles. This number of particles is just enough to cover the contact region for the duration of the simulation (4 units of time). As observed in the figure. the friction reaches an initial peak as in Fig. 5a, but the "stick-slip'' cycle is not repeated to the same degree. Instead. following the initial slip. the force essentially oscillates between positive and negative values. with an average near zero. In fact if we compare the average friction forces in Fig. Sa and Fig. 6a, after time = 1 (i.e.. after the initial stickslip). we find a substantial difference: Fig. 5a corresponds to an average frictional stress of 0.194. while Fig. 6a gives a negative value of -0.03 Clearly, then. the surface extent has a profound effect on the friction force. Consideration of the energy explains the foregoing friction behavior. Figure 6b shows the energy accumulated in the surroundings for the sliding conditions of Fig. 63. In this figure it is observed that the energy of the surroundings initially rises but afterwards does not appreciably increase. This result contrasts that of Fig. Sb i n which the energy lost to the surroundings steadily increases with time. When the surface is unlimited in extent. the energy is continually lost by propagation to the particles surrounding the contact, thus providing a source of dissipation. On the other hand. when the extent of the system is limited. energy does not continually flow out of the contact and dissipative effects are small. 3.3 Effect of Contact Length The role of energy propagation can also be seen when the contact length is varied. Figure 7a shows the amount of energy, per particle in thc contact region. that goes to the surroundings for several values of the contact length: L/b = 5. 10.20, and 40. The surface is of infinire extent. As the contact length increases. the amount of energy lost (per particle in the contact region) decreases, indicating that there is proportionately more energy delivered to the surroundings when the contact length is small than when it is large. This result is expected since, for a shorter contact, a greater percentage of the particles are at the boundary of the contact region. The longer contact length also affects the nature of the friction force. as observed in Fig. 7b.
-0.5
0
3
1
4
TIME~V 0t/b)
Figure 5a: Frictional stress vs. time for Q J k b = 0.5 and L/b = 5. n
stf!
0.6
0.5
Pi
d
3
rn 0
b
cl
0.4
0.3
0.2
$(
u
g 0.1 z W
I
0
1
I
I
I
4
T I M E ~ 0Vt/b; Figure 5b: Energy to surroundings for QJkb = 0.5 and L/b = 5 .
179
This figure provides the friction trace for the same interaction force amplitude as i n Fig. Sa. but with a dimensionless contact length. L/b = 40 instead of L/h = 5 . Initially the friction peaks. as in Fig. Sa. but after the sudden drop, the friction force remains low. Since there is a greater proportion of energy which remains in the contact, the amount of dissipation is less and the average friction is small. Figure 8 shows the time average frictional stress as a function of contact length for several amplitudes of the interaction force. The average is taken over the duration of the sliding which is 3 units of time. In each case the frictional stress is observed to decrease with increasing contact length. The consistency of the trend provides further support for the previous conclusion--namely that the magnitude of the average frictional stress is influenced by the tendency for energy to propagate away from the contact. With a smaller contact region, the energy is more readily lost to the surroundings and the average frictional stress is greater.
3.4 Continuous Approximation Inspection of the governing equation (2) reveals that an analogy can be made between the model of Fig. 2b and that of a string on an elastic foundation. In this regard, we introduce u(x.1) as the transverse string displacement and make the following associations:
as-
0.3
b m A
< z 0 i?
surface extent: 10 particles R
-O.l
2 W
-0.3 -0.5
I
I
I
I
0
1
2
3
TIME(V0t/b)
Figure 6a: Frictional stress vs. time for finite surface extent.
0.5
m
3
0.4 L
0
b 0.3 b
-
-
0.2
-
-
8 A
(3)
2 Q#b=0.5 L/b=5
0.1
at-
Using the above definitions, equation (2) becomes
at?
Y
d
We also let
a h c, ?+ a% pu = (I
0.5
m
a4
ii, e 7
--
3b
0.7
0.0
0
1
3
4
TIME t V 0 t/b) Figure 6b: Energy to surroundings for finite surface extent.
I80
Equation (3) is recognized as the governing equation for a string on an elastic foundation. There are two characteristic homogenous solutions to this equation, given by [ 131:
sw
d
I
I
1
0.4
d
3
cn
0.3
0 where w is the frequency of vibration, y is the
wavenumber. and 7 is magnitude of the wavenumber when the wavenumber is imaginary. Equation (4a) represents a propagating solution and applies when w2 > p. Conversely. equation (4b) represents a non-propagating solution and applies when w2 5 p. Recalling that 0 = the conditions on w indicate that energy will propagate away from the contact region when frequencies are excited that exceed some critical value characteristic of the system. This result indicates that continuous energy loss to the surroundings requires some excitation at the interface. We expect therefore that a plucking instability should be associated with the Occurrence of dissipation and should increase the friction force. We observed this effect in Fig. 3. When Q&b = 0.05 and there is no plucking. the average frictional stress is very low (.ooOOl). On the other hand. with Q&b = 0.2 there was a clear instability and the average frictional stress was markedly higher (.03). Thus a fourfold increase in the interaction force increased the average frictional stress by more than three orders of magnitude.
urn.
3.5 Relation to Other Mechanisms It was shown in previous sections that the magnitude of the dissipative component of friction is influenced by the efficiency with which energy is propagated away from the contacting interface. We now discuss this result in light of other mechanisms of dissipation identified in similar studies. Sokoloff (1992) [HI finds that internal damping within a (finite) crystal lattice is critical for the Occurrence of frictional losses at the interface. In fact, for a purely elastic model. the dissipative component of friction is found to be virtually zero. When damping is present, the elastic waves which emanate from the contact return to the interface
EE-
rn 0 cl
*u
0.2
d 0.1
W
zW 0.0 0
T~ME(V 0t i )
3
Figure 7a: Energy to surroundings for variable contact length.
m
cn
0.3
cn
<
0.1
z F1 -O.l
0
u
U
-0.3
v t
0
1
TIME (V0t i )
3
Figure 7b: Frictional stress vs. time for QJkb = 0.5 and L/b = 40.
181
diminished in energy. In effect some of energy has propagated away from the interface and has not returned. Damping, therefore, can be viewed as a means of facilitating energy propagation away from the contact. Shinjo and Hirano (1993) [9]. studied the FK model with kinetic energy terms. In their work, a set of masses was given an initial velocity and then allowed to move under the influence of a harmonic potential fixed in space. For certain combinations of interaction force and initial velocity. a friction force arose which caused the velocity of the center of mass to decrease to zero. Since the system in question conserved energy with all of the particles comprising the interface, their results demonstrate that dissipation can exist even when no energy is propagated away from the interface. This mechanism, therefore, appears to contradict our observation that energy propagation is primarily responsible for the dissipative component of friction. We may resolve this apparent discrepancy by comparing the macroscopic work done to the net work done. In the study of Shinjo and Hirano. the total energy of the system is conserved. Since any net work done on the system must change the energy of the system, the net rate of work done on the system by the interaction force is identically zero for all time. On the other hand, there is a finite friction force which opposes the translation of the center of mass. This friction force does macroscopic work at a rate equal to -FV,where F is the friction force and V is the velocity of the center of mass. Again, since the net rate of work remains zero, the loss in macroscopic translation energy is exactly balanced by the increase in the internal kinetic and potential energies [9]. Analogous effects are found in our calculations which correspond to the case of a force translating over a stationary surface (Fig. 2). Fig. 9a shows both the macroscopic work and net work done by the interaction force for the operating conditions of Fig. 5a. Each of the work terms is normalized in the usual manner. The macroscopic work done is just the time integration of the product of the friction force and velocity. and defines the degree of dissipation. The net work done is the time integration of the rate of work done on each particle of the system by the interaction force. This
3
0.3
W
m m W
p?
b
r A 0.2 Cll
6 z 0
F: 0.1 u
Ed >
6 0.0
,
0
I
10
1
,
20
I
30
I
I
40
CONTACT LENGTH (L/b) Figure 8: Average frictional stress vs. contact length for several force amplitudes. difference between the two forms of work, therefore, reveals the component of friction which is not associated with energy changes in the system. Hence, there are two components of dissipation: one component of dissipation depends on the flow of energy away from the contact to the rest of the system, while the other occurs even when no net energy is lost. Our concern here is upon which of these mechanisms is predominant. In Fig. 9a. it is seen that the bulk of the frictional loss (i.e.. the macroscopic work) is associated with the increase in the energy of the system (i.e., the net work). The same result is observed in Fig. 9b, where the extent of the surface is limited to 10 particles. Hence, while frictional losses do not necessitate the flow of energy from the contact region, they are greatly enhanced by it. 4.0 CONCLUSIONS A simple model of frictional contact was presented which focuses on the energetics of sliding at the molecular scale. The model. which
182
macroscopic work
0.9
2 0.8 9
w
W
-
0.6
-
n 0.4 d 0
*
-
0.7
z0 0.5
contains aspects of the I 0 and FK models. was used simulate the frictional interactions belween a rigid asperity and an adsorbed film. The rigid asperity was modeled as a translating sinusoidal potential. while the adsorbed film was modeled by a set of harmonically coupled masses which were attached by springs to a rigid base. The frictional behavior of the system was determined by numerical integration of the equations of motion for each of the model molecules. Three key observations werc made:
-
0.3
-
0.2
I I
0
I
I
I
1
4
3
1
T I M E ~ 0Vt/b)
Figure 9a: Macroscopic work vs. total work done on system. 0.5
~
-
0.3
z n 0.2
0
1
These three observations can hc attributed to the same mechanism: energy propagation away from the region of conhct. First, an instability in the motions of the molecules is needed to excire the molecules in the interface. Second. for significant dissipation to occur, an avenue is required through which energy can be lost from the contact. Third. when the contact length is large, there is less availability for energy to propagate away from the interface, yielding lower frictional stress. Although the results presented here correspond to a limited number of cases investigated by simple model, they provide evidence that frictional losses depend, to a great extent. on the way in which energy propagates away from the contact.
surface extent: 10 particles
d 0
*
0.1
0 .o
0
I
I
1
2
1
TIME(V 0t/b)
The average frictional stress was found to hc much higher when the surface was unlimited in extent than when it was finite. The average frictional stress decreased as thc lengrh of contact increased.
Q0/kb = 0.5 L/b=5
9 W
1
- net work - - macroscopic work
2 0.4
e W
~- 1
- 1
An appreciable dissipative component of friction was found to occur only when the interaction force was of sufficient magnitude. This effect corresponded with the onset of an instability in the positions of the model particles.
1
3
Figure 9b: Work comparison for surface of finite extent.
4
ACKNOWLEDGMENT The author would like to thank the Narional Science Foundation for support of his work.
183
REFERENCES 1. Tabor, D., in Fundamentals of Friction: Macroscopic and Microscopic Processes, I. L. Singer and H. M. Pollock, eds., NATO AS1 Series, vol. 220, Kluwer Academic Publishers, The Netherlands, p. 3 (1992).
2. Tomlinson, G. A., Phil. Mag., vol. 7, n. 46, p 905 (1929). 3. McClelland, G. M., in Adhesion and Friction, M. Grunze and H. J. Kreuzer, eds., Springer Series in Surface Science, vol. 17, Springer Verlag, p. 1 (1989). 4. McClelland, G. M., and Glosli, I. N., in Fundamentals of Friction: Macroscopic and Microscopic Processes, I. L. Singer and H. M. Pollock, eds., NATO AS1 Series, vol. 220, Kluwer Academic Publishers, p. 405 (1992). 5. Frenkel Y. and Kontorova T. ,Zh. Eksp. Teor. Fiz., vol. 8, p. 1340 (1938).
6. Peyrard M., and Aubrey, S., J. Phys. C: Solid State Phys., vol. 16, p. 1593 (1983).
7. Hirano, M. and Shinjo, K., Phys. Rev. B. , vol. 41, n. 17, p. 11837 (1990). 8. Shinjo, K., and Hirano, M., Surface Science, vol. 283, p. 473 (1993). 9. Sokoloff, J. B., J. Appl. Phys., vol. 72, n. 4, p. 1262 (1992). 10. Sokoloff, J. B., Phys. Rev. B., vol. 42, p. 760 (1990). 11. Tabor, D., Gases, Liquids, and Solids, 3rd edition, Cambridge University Press, p. 181 (1990). 12. Gear, W. C., Numerical Initial Value Problems in Ordinary Differential Equations, PrenticeHall, New Jersey, p. 47 (1971). 13. Graff, K. F, Wave Motion in Elastic Solids, Ohio State University Press, p. 51 (1975).