7816
Langmuir 1999, 15, 7816-7821
Molecular Dynamics Simulation of Friction of Hydrocarbon Thin Films Hiroyuki Tamura, Muneo Yoshida, Kenichi Kusakabe, Chung Young-Mo, Ryuji Miura, Momoji Kubo, Kazuo Teraishi, Abhijit Chatterjee, and Akira Miyamoto* Department of Materials Chemistry, Graduate School of Engineering, Tohoku University, Aoba-ku, Sendai 980-8579, Japan Received April 29, 1998. In Final Form: June 3, 1999 Molecular Dynamics (MD) simulations were performed to investigate the dynamic behavior of hydrocarbon molecules under shear conditions. Frictional properties of cyclohexane, n-hexane, and iso-hexane thin films confined between two solid surfaces were calculated. Because the affinity of the solid surfaces in these simulations is strong, slippages occurred at inner parts of the confined films, whereas no slippages were observed at the solid boundaries. The hexagonal closest packing structure was observed for the adsorbed cyclohexane molecular layers. The branched methyl groups in the iso-hexane molecules increase the shear stress between the molecular layers. For the n-hexane monolayer, molecules were observed to roll during the sliding simulations. Rolling of the n-hexane molecules decreased the shear stress.
1. Introduction In machinery and in the automobile, electronics, and many other industries, the control of rheological and tribological phenomena is particularly important to prevent wear and save energy. For example, to increase the capacity of magnetic storage systems, lubrication of the head-disk interface must be improved. To rationalize rheological and tribological mechanisms, the atomic-level origins of friction have been investigated, both experimentally and theoretically. New experimental techniques, such as the surface force apparatus (SFA) and atomic force microscopy (AFM), have been used to study the molecular and atomic-level rheology and tribology. The SFA provides information on the structural and rheological properties of thin molecular films,1-5- such as surface-induced phase transitions, increasing viscosity, and quantized frictional forces. Recently, molecular dynamics (MD) simulations have been widely used to obtain molecular and atomiclevel information about tribological and rheological phenomena; for example, frictional phenomena of two surfaces separated by Lennard-Jones fluids were simulated.6-10 In these simulations, it was observed that for thin molecular films, the molecules ordered into discrete layers between the surfaces, resulting in a decrease in the diffusion rate and an increase in viscosity. The stick-slip motions observed in thin films of the sphere-structured molecules were also investigated using the L-J fluid model.8 Moreover, the behavior of confined chain molecules and their shear properties have been investigated using * Corresponding author. Tel.: +81 22 217 7233. Fax: +81 22 217 7235. E-mail:
[email protected]. (1) Gee, M. L.; McGuiggan, P. M.; Israelachvili, J. N.; Homola, A. M. J. Chem. Phys. 1990, 93, 1895. (2) Granick, S. Science 1991, 253, 1374. (3) Klein, J.; Kumacheva, E. Science 1995, 269, 816. (4) Klein, J.; Kumacheva, E. J. Chem. Phys. 1998, 108, 6996. (5) Klein, J.; Kumacheva, E. J. Chem. Phys. 1998, 108, 7010. (6) Bitsanis, I.; Magda, J. J.; Tirrell, M.; Davis, H. T. J. Chem. Phys. 1987, 87, 1733. (7) Bitsanis, I.; Somers, S. A.; Davis, H. T.; Tirrell, M. J. Chem. Phys. 1990, 93, 3427. (8) Lupkowski, M.; Swol, F. J. Chem. Phys. 1991, 95, 1995. (9) Thompson, P. A.; Robbins, M. O. Phys. Rev. A 1990, 41, 6830. (10) Gao, J.; Luedtke, W. D.; Landman, U. Phys. Rev. Lett. 1997, 79, 705.
the model L-J chain molecules.11-21 Such atomistic information about the tribological and rheological processes is important not only for fundamental scientific understanding but also for the technological development of lubricants. The molecular architecture of a lubricant is believed to have a considerable effect on tribological and rheological phenomena such as frictional properties, reactivities of lubricants, and wear of surfaces. However, the properties of lubricants measured experimentally reflect the combined effect of various parameters such as the orientation, the packing, and the interaction of molecules. To design new lubricant molecules, it is necessary to understand the effects of the molecular architecture on the mechanisms of lubrication. Recently, the Quantitative Structure Property Relationship (QSPR) method has been applied to the design of lubricant molecules.22,23 QSPR provides a relationship between the molecular architecture and the frictional properties of lubricants. MD simulations are also useful in designing lubricants due to its ability to take into account atomic-level dynamic behaviors and energetics. For example, the effect of molecular branching on confined structures and frictional properties has been investigated both experimentally1,24 and theoretically.11,17,18 (11) Gao, J.; Luedtke, W. D.; Landman, U. J. Chem. Phys. 1997, 106, 4309. (12) Dijkstra, M. Europhys. Lett. 1997, 37, 281. (13) Bitsanis, I. A.; Pan, C. J. Chem. Phys. 1993, 99, 5520. (14) Gupta, S. A.; Cochran, H. D.; Westermann-Clark, G. B.; Bitsanis, I. A. J. Chem. Phys. 1994, 100, 8444. (15) Padilla, P.; Toxvaerd, S. J. Chem. Phys. 1994, 101, 1490. (16) Stevens, M. J.; Mondello, M.; Grest, G. S.; Cui, S. T.; Cochran, H. D.; Cummings, P. T. J. Chem. Phys. 1997, 106, 7303. (17) Gupta, S. A.; Cochran, H. D.; Cummings, P. T. J. Chem. Phys. 1997, 107, 10316. (18) Gupta, S. A.; Cochran, H. D.; Cummings, P. T. J. Chem. Phys. 1997, 107, 10327. (19) Manias, E.; Bitsanis, I.; Hadziioannou, G.; Brinke, G. Europhys. Lett. 1996, 33, 371. (20) Manias, E.; Hadziioannou, G.; Bitsanis, I.; Brinke, G. Europhys. Lett. 1993, 24, 99. (21) Manias, E.; Hadziioannou, G.; Brinke, G. Langmuir 1996, 12, 4587. (22) Tsubouchi, T.; Hata, H. Tribology Int. 1994, 27, 183. (23) Tsubouchi, T.; Hata, H. Tribology Int. 1995, 28, 335.
10.1021/la9805084 CCC: $18.00 © 1999 American Chemical Society Published on Web 09/11/1999
Molecular Dynamics Simulation of Friction
Langmuir, Vol. 15, No. 22, 1999 7817
E)
∑b Hb {1 - exp[-a(b - b0)]}2 + ∑θ Hθ (θ - θ0)2 +
∑φ Hφ [1 + cos(nφ - φ0)] + ∑i ∑ j>i
Figure 1. Unit cell of the model used for the simulation.
The friction coefficients of molecules increase as the degree of branching increases due to elastohydrodynamic lubrication (EHL),24 whereas the frictional forces of branched molecules are lower than those of straight molecules in SFA measurements.1 It has been observed that the friction properties of lubricant films depend not only on the molecular structure but also on the film thickness.1-5 In previous nonequilibrium molecular dynamics (NEMD) simulations, the friction coefficient of a confined cyclohexane film was found to be higher than that of benzene.25 In the present study, the dynamics of cyclohexane, n-hexane, and iso-hexane molecules under shear are simulated to investigate the effects of film thickness and molecular architecture on frictional properties. 2. Methods Figure 1 shows the model system of the thin film of hydrocarbon molecules confined between two atomically smooth Fe(001) surfaces, where the periodic boundary conditions are imposed in the x and y directions. Each surface model is composed of 324 atoms and has a unit cell size of 25.8 Å × 25.8 Å. The vibrations of the lowest layer of the lower Fe solid and the highest layer of the upper Fe solid were frozen. The interaction between atoms in the Fe(001) surfaces is expressed by the Morse potential:
E)
∑i ∑ j>i
[
Aij
-
r12 ij
Bij r6ij
+
]
ZiZje2 rij
(2)
Here, each term describes the Morse potential (Hb ) bond energy, a ) force constant, b ) bond length, b0 ) bond length at minimum energy), the angle of bending (Hθ ) force constant, θ ) bending angle, θ0 ) bending angle at minimum energy), the torsion angles (Hφ ) force constant, n ) order of rotation axis, φ ) torsion angle, φ0 ) torsion angle at minimum energy), the van der Waals (Aij ) (Ai‚Aj)1/2, Bij ) (Bi‚Bj)1/2) and the coulomb (Z ) charge, e ) elementary electric charge) potentials, respectively. The Lennard-Jones type potential was adopted for the interaction between Fe atoms and organic molecules. Values of Ai and Bi for Fe atoms were determined previously in reproducing the cohesive energy and the density of an Fe solid.29 The force-field parameters used here are summarized in Figure 2. Lubricant molecules were placed between two solid surfaces, and equilibrium structures were obtained at 300 K under a constant normal force (equivalent to 1 GPa of pressure) loaded on the frozen atoms of both Fe(001) solids. Then, sliding simulations were performed for 50 ps by sliding the frozen layer of the lower Fe(001) solid, with the frozen layer of the upper Fe(001) solid held fixed. During the sliding simulations, a constant normal force, equivalent to 1 GPa, was loaded on the frozen atoms of both Fe(001) solids. Heat generated by sliding was removed by scaling the velocities of Fe atoms at 300 K. The temperature of the confined molecules increased slightly when sliding began. The time step and the sliding velocity during the sliding simulations were 0.5 fs and 100 m/s, respectively. The Verlet algorithm29 was used for calculation of atomic motions; atom positions were stored every 50 fs. The sum of the forces acting on the atoms of the Fe(001) solid in each direction were averaged every 50 fs. The frictional force was considered to be proportional to the area of contact for the molecularly smooth surfaces.1 In this study, shear stress S was defined as Fx ) S × A, where Fx and A are the average lateral force during the sliding simulations and the area of the unit cell, respectively. 3. Results and Discussion
D [exp{-2β(rij - r0)} - 2 exp{-β(rij - r0)}] (1)
where D ) 0.42 eV is the bond energy, β ) 1.98 Å-1 is the force constant, rij is the bond length, and r0 ) 2.84 Å is the bond length at minimum energy, which was determined previously in reproducing the vaporization energy and the elastic constants.26 The interaction between atoms in organic molecules was approximated by the consistent valence force field (CVFF) potential:27 (24) Muraki, M. Tribology Int. 1987, 20, 347. (25) Yamano, H.; Shiota, K.; Miura, R.; Katagiri, M.; Kubo, M.; Stirling, A.; Broclawik, E.; Miyamoto, A.; Tsubouchi, T. Thin Solid Films 1996, 282, 598. (26) Girifalco, L. A.; Weizer, V. G. Phys. Rev. 1959, 114, 687. (27) Dauber-Osguthorpe, P.; Roberts, V. A.; Osguthorpe, D. J.; Wolff, J.; Genest, M.; Hagler, A. T. Proteins: Struct., Funct., Genet. 1988, 4, 31.
3.1. Effects of Film Thickness on Frictional Properties. To investigate the influence of film thickness and molecular density on frictional properties, of the cyclohexane molecules were calculated as a function of the number of molecules between the two solid surfaces. The average film thickness and the average shear stress in 50 ps of the sliding simulations are summarized in Figure 3 a and b, respectively. All cyclohexane molecules maintained the chair conformation during the simulations. When the number of molecules in the unit cell was 20 (the area per molecule is 33.3 Å2) or less, the monolayer structure was obtained. The structure of the cyclohexane monolayers could be regarded as hexagonal closest packing (hcp). The orientation of the cyclohexane rings is almost parallel to the solid surfaces when the number of molecules is 16 (41.6 Å2/molecule) or less; however, at densities (28) Halicioglu, T.; Pound, G. M. Phys. Status Solidi 1975, 30, 619. (29) Swope, W. C.; Andersen, H. C.; Behrens, P. H.; Wilson, K. R. J. Chem. Phys. 1982, 76, 637.
7818
Langmuir, Vol. 15, No. 22, 1999
Tamura et al.
Figure 3. (a) Film thicknesses and (b) average shear stress of the cyclohexane films as a function of the number of molecules between the two solid surfaces.
beyond this value, the orientation of the cyclohexane rings is more disordered. In the cases of 16 and 20 molecules, in which clearer ordered structures are observed, lower friction coefficients are obtained. It was observed both experimentally30 and theoretically31 that for adsorbed noble gas atoms, the static friction of the solid monolayers
was lower than those of the fluid monolayers. It seems that for the cyclohexane monolayers, low friction is obtained when the mobility of the molecules is restricted. The density and velocity profiles of the cyclohexane molecules normal to the wall are shown in Figure 4. The center of mass of the cyclohexane molecules are inhomogeneously distributed normal to the solid surfaces during the sliding simulations, similar to the previous study on L-J oligomers.21 As the affinity of the solid surfaces in the present simulations is strong, slippages occurred at inner parts of the confined films, whereas no significant slippages were observed at the solid boundaries except for the monolayers. Similar behavior was also observed in MD simulations using the L-J oligomers when the wall affinity was strong enough.21 Nine discrete layers were observed when 160 molecules were confined; in this case, the first and second layers from the solid surface were locked, so slippage occurred at the middle part of the film (Figure 4b). The calculated shear stress of bilayers was lower than that of monolayers. A similar tendency was observed with SFA measurements of friction between mica surfaces and spherical molecules.1 The difference between the shear stress of the monolayers and that of bilayers can be explained in terms of the difference between molecular-
(30) Krim, J.; Solina, D. H.; Chiarello, R. Phys. Rev. Lett. 1991, 66, 181.
(31) Cieplak, M.; Smith, E. D.; Robbins, M. O. Science 1994, 265, 1209.
Figure 2. Potential parameters used in the simulation.
Molecular Dynamics Simulation of Friction
Langmuir, Vol. 15, No. 22, 1999 7819
Figure 5. (a) Snapshots of the cyclohexane bilayer (32 molecules) during the sliding simulation. Only the bond lines between carbon atoms are shown. (b) Schematic model of the change in packing between upper and lower cyclohexane layers. Molecules from lower and upper layers are represented by white and dark hexagons, respectively.
Figure 4. Density profile and the velocity profile for the cyclohexane molecules: (a) 4 layer; (b) 9 layer. Solid lines and filled points indicate the distribution along the Z direction for the local density (number of molecules/Å) and the velocity of molecules, respectively.
molecular friction and molecular-solid friction. Low shear stress was obtained when 32 molecules were confined. In this case, a slight boundary slip was observed due to the incommensurate contact between the ordered cyclohexane layers and the Fe(001) lattice. The calculated shear stress of three layers was lower than that of the bilayers, as observed in the SFA measurement,1 however, the difference was very slight in the present simulations. In the cases of more than three layers, shear stress did not depend on film thickness. When film thickness was increased to nine layers, shear stress was similar to that of four layers although the layered structure was more disordered. 3.2. Dynamic Behavior of Cyclohexane Molecules. To obtain insight into the friction between the cyclohexane layers, the change in shear stress and the corresponding change in molecular behavior were investigated. Snapshots and schematic illustrations of the cyclohexane bilayer consisting of 32 molecules from the sliding simulation are shown in Figure 5 a and b, respectively. A clear epitaxial matching between upper and lower cyclohexane layers was observed at this density. A periodic change in the packing between upper and lower molecular layers was observed during the sliding simulations. The change in shear stress and in the van der Waals energy of the cyclohexane bilayer during the sliding simulation is shown in Figure 6 a and b, respectively. The shear stress and van der Waals energy showed an oscillation with the same periodicity, which can be explained as follows. At the valley of the van der Waals energy, the cyclohexane molecules form roughly an hcp structure. As sliding proceeds, the hcp structure of the cyclohexane molecules is distorted, and the van der Waals energy increases. In this process, a peak of shear stress is observed due to the shear resistance at the interface between upper and lower
Figure 6. Change in (a) shear stress and (b) van der Waals energy during the sliding simulation of the cyclohexane bilayer (32 molecules).
molecular layers. When sliding is continued further, the cyclohexane bilayer regains a stable hcp structure, accompanied by a decrease in the van der Waals energy and a valley of shear stress. The above process is repeated during the sliding simulations. The periodicity of the shear stress corresponded to the periodicity of the hcp structure of the cyclohexane layer. Such behavior is caused by the contribution of the in-plane structure to the stress component parallel to the shearing direction.32 For bilayers, particularly those with 32 molecules, accumulation (32) Haugstad, G.; Gladfelter, W. L. Langmuir 1993, 9, 3717.
7820
Langmuir, Vol. 15, No. 22, 1999
Tamura et al.
Table 1. Average Shear Stress (108 N/m2) of Each Model in 50 ps of Sliding Simulations monolayer (16 molecules) bilayer (32 molecules)
n-hexane
iso-hexane
cyclohexane
2.2 1.0
3.2 1.8
3.6 1.8
and release of shear stress were clearly observed. This feature decreased as the film thickness increased. It has been observed that confined cyclohexane thin films can sustain a shear stress up to some yield value.3-5 To investigate this behavior, spring simulations similar to the ones used for the L-J fluid8 were performed on the cyclohexane bilayer (32 molecules). The bottom wall was pulled by a spring connected to a stage that was moving at constant velocity. Initially, the cyclohexane bilayer maintained the equilibrium structure, and shear stress increased as the spring stretched. At shear stresses beyond the yield point, which corresponds to the peak of the shear stress, a slippage occurred at the interface of upper and lower molecular layers. It was observed in the MD simulations using the L-J fluid that stick-slip behavior depended on the stage velocity and the spring constant; in particular, stick-slip behavior was not observed at high stage velocities and/or high spring constants.8 In the present simulations, where the spring velocity and the spring constant were 100 m/s and 80 N/m, respectively, stick-slip behavior was observed due to the shear resistance between the cyclohexane layers. 3.3. Frictional Properties of the Monolayer and Bilayer of n-Hexane and iso-Hexane Molecules. Sliding simulations of n-hexane and iso-hexane molecules were performed to investigate the effect of branching. For similar confined short chain molecules, discontinuous changes in film thickness and shear stress have been observed in SFA measurements.1 In MD simulations using bead-spring models, the layered structure of confined molecules were reproduced for both straight and branched chain molecules.11,17 In the present simulations, layered structures were found for both n-hexane and iso-hexane molecules. In monolayers, orientation of molecules did not change much during the sliding simulations. The average shear stress in 50 ps of sliding simulations is summarized in Table 1. Shear stress of the bilayers was lower than that of the monolayers, an observation that was also made using SFA measurements of friction between mica surfaces and layered molecular films.1 The changes in shear stress and in the van der Waals energy of the n-hexane monolayer during the sliding simulations are shown in Figure 7 a and b, respectively. The periodicity of the shear stress and the van der Waals energy corresponds to the periodicity of the Fe(001) lattice. The same periodicity in shear stress and van der Waals energy was also found during the sliding simulations of the iso-hexane and cyclohexane monolayers. These results suggest that the periodicity of the shear stress found in the monolayers may originate from the periodicity of the lattice. 3.4. Effect of Molecular Rolling. Shear stress of the iso-hexane monolayer was higher than that of the nhexane monolayer. During the sliding simulations of the n-hexane monolayer, the rolling of molecules around the chain axis was observed for the n-hexane molecules, which oriented perpendicular to the sliding direction. However, in the case of the iso-hexane monolayer, molecular rolling was not observed. It seems that the branched methyl groups of the iso-hexane molecule disturbed rolling of the molecules. This rolling seems to decrease the shear stress. To clarify this effect, frictional properties of an anisotropic n-hexane monolayer were investigated. The n-hexane
Figure 7. Change in (a) shear stress and (b) van der Waals energy during the sliding simulation of the n-hexane monolayer.
Figure 8. Equilibrium structure of the anisotropic n-hexane monolayer, constructed to explore effects of rolling along the chain axis (as shown in the illustration under the picture). Average shear stress (108 N/m2) along the x and y directions was 1.5 and 2.8, respectively.
molecules were arranged anisotropically on the solid surfaces as shown in Figure 8. Shear stress along the x direction, which is perpendicular to the carbon chains, and along the y direction, which is parallel to the carbon chains, was calculated. Average shear stress for sliding along the x and y directions was 1.5 and 2.9 (108 N/m2), respectively. As expected, shear stress is lower in the x direction than in the y direction due to rolling of the
Molecular Dynamics Simulation of Friction
Langmuir, Vol. 15, No. 22, 1999 7821
Figure 9. Snapshots of the n-hexane bilayer during the sliding simulation. n-Hexane molecules and Fe atoms are represented by bond lines and white circles, respectively: (a) Before sliding, molecules were bridged; (b) after 10 ps of sliding, the bridged molecules were pulled along the sliding direction; (c) after 15 ps of sliding, the bridged molecule was torn off from one side; (d) after 50 ps of sliding, all molecules were incorporated into either upper or lower layers.
molecules. This effect is of interest in interpreting the frictional properties of the monolayer consisting of chained molecules. In the present simulations, molecular rolling was observed only in the monolayer of the n-hexane molecules. Therefore, the effects of rolling on frictional properties were not apparent in the confined multilayer of chained molecules. It was observed in previous MD simulations that diffusion coefficients of short chained molecules on the wall were higher than in the inner part of the fluid film.15 This behavior seems to be caused by mechanisms other than molecular rolling. 3.5. Effect of Bridging. As in the sliding simulations of the n-hexane and iso-hexane confined films, slip also occurred at the inner part of the confined films, except for in the monolayers. This tendency depends on the affinity of the walls.20,21 The layered structure was enhanced after sliding; molecules at the solid boundaries especially tend to adsorb more parallel to the solid surfaces. A similar tendency has also been observed in NEMD simulations of L-J oligomers.21 In the present simulations, molecular alignment tends to decrease the shear stress. Snapshots of the n-hexane bilayer during the sliding simulations and the corresponding shear stress and van der Waals energy are shown in Figures 9 and 10, respectively. In these simulations, two molecules were bridged between the upper and lower molecular layers before sliding. The van der Waals energy of the n-hexane bilayer initially increased because the bridging molecules were pulled. Shear stress decreased after all molecules were adsorbed parallel to the solid surfaces. When the bridging molecules were torn off from one side, the van der Waals energy decreased, in agreement with a previous study.33 3.6. Effect of Branching. Shear-induced alignments were also observed for iso-hexane molecules. In the present simulations, for both n-hexane and iso-hexane molecules, the rate of bridging affected the shear stress. However, when all molecules at solid boundaries adsorbed parallel to the solid surfaces, the shear stress of the iso-hexane bilayer was higher than that of the n-hexane bilayer. These results suggest that the branched methyl groups in isohexane molecules increased the friction between upper and lower molecular layers. MD simulations using the L-J oligomers model also showed that the branches (33) Manias, E.; Bitsanis, I.; Hadziioannou, G.; Brinke, G. Mol. Phys. 1995, 85, 1017.
Figure 10. Change in (a) shear stress and (b) van der Waals energy of the n-hexane bilayer. Average shear stress (108 N/m2), from 0 to 25 ps and from 25 to 50 ps, is 1.5 and 0.5, respectively.
interlock with the branches of other molecules.17 It has been observed under EHL conditions that the friction coefficients of chain molecules increase as the degree of branching increases.24 On the other hand, SFA measurements showed that the frictional forces of branched alkane molecules are smaller than those of the n-alkane molecules if the chain length is long enough.1 These results indicate that the effect of branching on frictional properties depends on the lubricant chain length, applied shear rate, film thickness, and affinity of solid surfaces. In the present study, using short chain molecules and a strong wall affinity, branching increased the shear stress between molecular layers. 4. Conclusions A comparative study of cyclohexane, n-hexane, and isohexane molecules confined and sheared between two solid surfaces was performed using NEMD simulations. As the affinity of the solid surfaces in these simulations is strong, interlayer slippages occurred at inner parts of the confined films, whereas no slippages were observed at the solid boundaries.21 A hexagonal closest packed structure was observed for the adsorbed cyclohexane molecular layers. The shear stress showed a periodic oscillation for the cyclohexane bilayer that was correlated with the periodicity of the packing structure. When the stable structure of cyclohexane molecules was disturbed by sliding, high shear stress was generated. The branched methyl groups of the iso-hexane molecules resulted in an increase of the shear stress between molecular layers. In confined monolayers, the branched methyl groups of the iso-hexane molecules prevented the molecules from rolling. The investigation of other technologically interesting lubricants, such as traction fluids and perfluoropolyether, are currently in progress. LA9805084
