FACTA UNIVERSITATIS Series: Mechanical Engineering Vol. 17, N o 2, 2019, pp. 207 - 215 https://doi.org/10.22190/FUME190404026D © 2019 by University of Niš, Serbia | Creative Commons License: CC BY-NC-ND Original scientific paper  VERIFICATION OF RABINOWICZ’ CRITERION BY DIRECT MOLECULAR DYNAMICS MODELING Andrey I. Dmitriev 1,2 , Anton Yu. Nikonov 1,2 , Werner Österle 3 , Bai Cheng Jim 4 1 Institute of Strength Physics and Materials Science SB RAS, Tomsk, Russia 2 Tomsk State University, Tomsk, Russia 3 Federal Institute for Materials Research and Testing, Berlin, Germany 4 Institute for Composite Materials, Kaiserslautern, Germany Abstract. In the paper we use direct molecular dynamics modeling to validate the criterion for formation of wear debris proposed by E. Rabinowicz in 1958. A conventional molecular dynamics using a classical Tersoff’s potential was applied to simulate the sliding behavior within a thin film corresponding to a tribofilm formed from silica nano-particles in amorphous-like state. The simulation was carried out by varying the initial temperature and the spatial size of the simulated crystallite. The results show the change in sliding behavior of silica-based tribofilm depending on the temperature and the size parameter of the system under consideration. Thus increasing the temperature provides smooth sliding while at moderate conditions wear process can occur via debris formation. Our estimations show good correlation between predicted critical size of the simulated system and calculated energetic characteristics. Key Words: Molecular Dynamics, Adhesive Wear, Wear Debris, Critical Size, Surface Microasperities, Thermal Conditions, Shear Resistance Force 1. INTRODUCTION The problem of investigating relationships between surface roughness of solids and friction mechanisms is of great practical and scientific interest. As mentioned in [1–3], E. Rabinowicz was one of the first to explore this issue. In 1958 he proposed an original criterion that determines the contact mode of two friction surfaces by means of estimating contact stresses [4]. As long as the stresses in the vicinity of contacting microasperities are Received April 04, 2019 / Accepted June 18, 2019 Corresponding author: Andrey I. Dmitriev Institute of Strength Physics and Materials Science SB RAS, 634055, pr. Akademicheski 2/4, Tomsk, Russia. E-mail: dmitr@ispms.ru 208 A.I. DMITIREV, A.YU. NIKONOV, W. ÖSTERLE, B.C. JIM less than the yield strength of a material, it is plastically deformed. In the case the stresses exceed the yield strength, cracking of the surface layer of the material and formation of wear debris can occur. Thus, according to this criterion, the formation of wear debris and transition to adhesive wear are possible only when the contact area is larger than a critical value. The criterion proposed by E. Rabinowicz can be written as follows 2 0 2G w D    , (1) where σ 2 0 – shear strength, G – shear modulus and ∆w – increment of surface energy. Although this assumption can be verified numerically, there were some difficulties caused by the dominance of the methods of continual description of simulated medium. The limitation of continual models lies in their inability to explicitly simulate (without loss of accuracy) wear process because of complicated description of the formation of debonded fragments and new contact spots. Therefore, the particle method and its various modifications seem to be substantially more effective [5–7]. In terms of this discrete approach the processes of multiple cracking and rebinding of the medium can be explicitly simulated. Molecular dynamic simulation is one of the most widespread methods of discrete description. Most commonly, atomistic simulation of adhesive wear predicts continuous smoothening of surfaces rather than sliding of contacting bodies with formation of wear debris. First of all, this is due to the limitations on spatial dimensions of the simulated system. In this regard, it is worth to emphasize the studies [8, 9] done by the group of Molinari, which have firstly predicted the formation of wear particles induced by the damage of the surface layer of a material during adhesive interaction of surface microasperities within the framework of empirical molecular dynamics. This result was obtained by elaborating simple pair potentials with tuned inelastic properties, which then could be associated with macroscopic behavior of actual materials. The authors have proposed a simple analytic expression that is completely similar to Eq. (1) and that predicts the transition to the adhesive wear mechanism in terms of single asperity size: * 2 ( ) j w d G     (2) where σ 2 j like σ 2 0 in Eq. (1) characterizes shear strength, while 𝜆 is the coefficient related to the surface roughness geometry. In order to study the relation between surface roughness parameters and wear mechanisms a single contact model elaborated within the framework of empirical molecular dynamics was used in [8, 9], which did not imply a tight link to the scale of a simulated sample. At the same time, worth mentioning is the study of the authors of the present paper [10], where MD simulation of a tribofilm consisting of amorphous-like silicon dioxide particles resulted in the prediction of the conditions of wear mode transition of two contacting bodies from a laminar sliding regime to an abrasive wear regime with the formation of wear debris. In contrast to empirical molecular dynamics with phenomenological pair interaction, the Tersoff three-particle interatomic potential [11] adapted for calculations of different Si-O compounds like α-quartz and SiO2 glass [12] was used there. Moreover, in [10] it was considered a continuous tribofilm with Verification of Rabinowicz Criterion by Direct Molecular Dynamics Modeling 209 periodical boundary conditions rather than single contacting microasperities with various dimensions, and the tribofilm temperature was used as a criterion for transition between the friction regimes. When the kinetic temperature of the sample was above the glass- transition temperature of SiO2, the sliding regime was unstable and accompanied by the formation of an agglomerate of silicon and oxygen atoms. When simulating shear deformation of the tribolayer at the temperature corresponding to extreme conditions of a thermal burst in a contact patch, the sliding turns into the laminar regime. Considering equations (1) and (2), it can be supposed that the system temperature specified in the simulation determines energetic characteristics of the simulated medium, in particular the ∆w and σ 2 j(σ 2 0) values. However, this assumption requires detailed verification. Thus, the objective of this work is to study the sliding behavior of silica-based tribofilm in order to understand the influence of temperature of the atomic system and size parameter of the simulated sample on adhesive wear mechanisms. 2. NUMERICAL MODEL The initial structure of amorphous silica sample was created by the following algorithm. At the beginning the initial structure of alpha-quartz (SiO2) was formed by bonding Si–O tetrahedra together. The bulk sample in a shape of the parallelepiped with the following geometry: 15.0 x 11.6 x 8.3 nm along x (Lx), y (Ly) and z (Lz) directions respectively was used as a basic specimen (see Fig. 1a). (a) (b) Fig. 1 (a) Initial structure of an alpha-quartz sample before heating. (b) Resulting structure of the amorphous silica sample and a loading scheme for sliding simulation. Hereafter the atoms colored in yellow are used to indicate the features of deformation along the sample. The following zones are indicated: A – atoms under loading (loaded layers), B – crystallite layers, C – melted/amorphous area. Purple color denotes O atoms, green – Si atoms 210 A.I. DMITIREV, A.YU. NIKONOV, W. ÖSTERLE, B.C. JIM To study the influence of size parameter on sliding regime of the simulated sample three other specimens where the width (Lx) was equal to 10.0, 22.5 and 30.0 nm were generated as well. Thus, the width of the narrowest specimen was corresponded to 0.67 width of the base sample, and the number of atoms in it did not exceed 110 thousand. The width of two other samples corresponded to width of the base sample multiplied by 1.5 and 2.0, which contained 270 and 360 thousand atoms, respectively. Each specimen was divided into three zones: A, B and C, as shown in Fig. 1. To prepare the amorphous layer the central zone C in each specimen was subjected to heating up to a temperature of 6000 K within the framework of a microcanonical ensemble NVE until this fragment was completely melted. Height (Ly) of zone C in all samples was identical and equal to 10 nm. After 10 ps at 6000 K the melted fragment was quenched to the required temperature during 170 ps. At the end of this procedure, we obtained a block with amorphous silica interlayer as shown in Fig. 1b. For generation of positions and velocities the modeled sample was considered as an NVE ensemble that maintained the number of particles N, the occupied volume V and the energy of the system E. The integration of Nose-Hoover style non-Hamiltonian equations of motion was utilized with a time step Δt of 0.5 fs. To imitate the extension into a tribofilm, periodic boundary conditions were assigned along x and z direction. All simulations were performed using LAMMPS [13], while for the visualization of the simulation results and structure analysis the visualization tool OVITO [14] was used. The calculations were performed on the SKIF Cyberia supercomputer resource of Tomsk State University. The modeled structure was located between the two loaded layers and subjected to a sliding loading with constant velocities (V) +15 and –15 m/s applied to the upper and lower layers, respectively as shown schematically in Fig. 1b. Thus the total sliding velocity was 30 m/s. Despite the rather large values, these are the characteristic velocities used in MD simulation, since the calculation times are usually only a few nanoseconds. Simultaneously with shear loading a normal force (F) of 1.3 nN was applied along +y- and –y-directions yielding a total contact pressure of about 350 MPa. The temperature of the whole specimen was kept constant within NVT ensemble and was varied in the range from 300K (ambient temperature conditions) to 1100 K (for sliding simulation under flash temperature conditions at the local contact) [15]. The given temperature is achieved by using the velocity rescaling method during the MD simulation process from the energy balance described by Eq. (3): 2 1 3 2 2 N i i b i m v k T N  (3) where N is the atom number, kb is the Boltzmann constant, mi and vi are the i th atom mass and velocity, respectively. Since the mechanical properties of amorphous material is closely related to its density, the resulting density was kept constant at 70 at./nm³ (2.25g/cm3) for all tests. The integration of Nose-Hoover style non-Hamiltonian equations of motion was utilized with a time step Δt of 0.5 fs. Note that due to the action of periodic boundary conditions along the direction of applied load, we have not a spatially independent tribofilm, but a spatially limited fragment that regularly repeats in this direction. However, within this fragment, the implementation of various physical processes and mechanism of deformation are not limited. This means that the variable width of the fragment being modeled can be that critical spatial parameter determining the character of the wear process. Verification of Rabinowicz Criterion by Direct Molecular Dynamics Modeling 211 3. RESULTS AND DISCUSSION As mentioned above, the performed earlier sliding simulation using the Tersoff potential showed a pronounced effect of the system temperature on the ability of a silica tribofilm to exhibit smooth sliding [9]. The latter was obtained at the high temperature only, whereas at 300 K a stick-slip type of sliding was observed. Following the goal of the paper, we simulate the sliding of amorphous silica sample with the width of 15.0 nm at few different temperatures: 300 K, 400 K, 500 K, 600 K, 700 K, and 1100 K. (a) (b) (c) (d) (e) Fig. 2 Structures of the central fragment of different silica specimens after 1 ns of sliding: (a) Lx = 10.0 nm, T = 300 K; (b) Lx = 15.0 nm, T = 500 K; (c) Lx = 22.5 nm, T = 500 K; (d) Lx = 30.0 nm, T = 500 K. (e) Trajectories of atoms of the specimen with Lx = 15.0 nm after 1.0 ns of sliding at 500 K. Arrow marks the direction of atoms rotation According to the simulation results the critical temperature (Tc) above which the transition from unstable behavior with debris formation to smooth sliding conditions is taking place is about 500 K. For two other specimens with Lx=22.5 nm and Lx=30.0 nm, temperature Tc is also close to 500 K, while for the narrowest specimen with the width of only 10.0 nm the critical temperature is much lower and is about 300 K. This means that the size of the 212 A.I. DMITIREV, A.YU. NIKONOV, W. ÖSTERLE, B.C. JIM fragment being modeled is so small that debris formation in it is possible only under conditions of increasing viscosity (adhesive properties) of the modeled system in comparison with large-sized samples. Fig. 2 denotes the resulting structures of all considered specimens at the stage of sliding where wear debris forms in the amorphous layer. According to the trajectories of atoms (see Fig. 2e) these wear particles seem to adopt the function of nano-sized rollers. The same samples were also studied at elevated temperatures. Fig. 3 denotes the resulting structures of the specimens at the conditions of smooth sliding. In comparison with the data presented in Fig. 2, only the temperature was increased by 100 K, while all other parameters remained the same. It is well seen that the velocity accommodation between the upper and lower part of the specimen in that case occurs almost entirely within the amorphous layer. In contrast to the conditions of the critical temperature, no damaged regions are observed and shearing takes place homogeneously within the layer. (a) (b) (c) (d) (e) Fig. 3 Structures of the central fragment of different silica specimens after 1 ns of sliding: (a) Lx = 10.0 nm, T = 400 K; (b) Lx = 15.0 nm, T = 600 K; (c) Lx = 22.5 nm, T = 600 K; (d) Lx = 30.0 nm, T = 600 K. (e) Trajectories of atoms of the specimen with Lx = 15.0 nm after 2.0 ns of sliding at 700 K Arrows show the velocity distribution of atoms Instead of the movement of aggregates of linked atoms, now single atoms move as they do in liquid films – in opposite directions along the interface between the upper and lower part and a neutral layer in the middle. This is well indicated by arrows in Fig. 3e, where atoms trajectories at steady state motion are shown. Comparing Figs. 2 and 3 it can be Verification of Rabinowicz Criterion by Direct Molecular Dynamics Modeling 213 seen that increase in temperature above the critical value (~500K in our case for most samples) leads to change in the sliding regime from unstable with possible formation of wear debris to smooth sliding. Moreover, it was found that the critical temperature decreases rapidly with decreasing sample width of less than 15 nm, which corresponds to the width of the basic sample in our calculations. Thus, as we assumed, the temperature of the system acts as a controlled parameter that ensures a change in the adhesive properties of the model material. Increasing the temperature leads to a decrease in adhesion forces, and hence to a change in the sliding regime. At the same time, as mentioned early in the Introduction, there is a critical length scale of the system, which can change the conditions for transitions of the sliding regime. Fig. 4 demonstrates the dependence of the sample width on the Tc value. Fig. 4 Dependence of critical conditions for transitions of sliding regime on sample width and temperature conditions In order to verify Eqs. (1) and (2) we calculated the values of the resistance force to tangential motion within the amorphous interlayer acting on the atoms belonging to the loaded layers for the base sample at various temperatures. In addition, the estimation of the corresponding change in specific surface energy was made. The latter is calculated as the difference between the energy of a crystallite with a free surface and the same crystallite with periodic boundary conditions divided by the area of the free surface. Calculated results are presented in Fig. 5 and summarized in Table 1. According to the data, as the temperature rises, both energy characteristics fall. So the peak value for the shear resistance force at 300 K is ~1590 eV/Å and at 400 K is about 1500 eV/Å. At the same time, the specific surface energy decreases from 3.93 J/m 2 to 3.74 J/m 2 . In this case, the shear modulus, as clearly seen from Fig. 5, remains practically unchanged. Substituting the obtained values into Eqs. (1) or (2), we can find that the ratio of the critical size parameter d * (D) for temperatures 300K and 400K is about 0.93, while the estimation according to data in Fig. 4 gives a close value of ~ 0.83. Thus, the results of direct MD modeling indicate that the criterion proposed by E. Rabinowicz can be used to estimate the critical size parameter leading to the formation of wear debris in the model of tribofilm as well. 214 A.I. DMITIREV, A.YU. NIKONOV, W. ÖSTERLE, B.C. JIM Fig. 5 Time dependence of the resistance force to tangential motion within the silica amorphous interlayer in the base sample at various thermal conditions Table 1 Energetic parameters calculated for basic sample at various thermal conditions Temperature, K Specific surface energy, J/m 2 Resistance force, eV/Å 300 3,9349 1590 400 3,7791 1500 500 3,5737 1420 600 3,3939 1300 700 3,3493 1150 1100 3,1029 950 5. CONCLUSION Direct molecular dynamics modeling of the shear deformation of a silica sample containing an amorphous interlayer has demonstrated the effect of temperature and size of the system on the character of its sliding while other loading conditions remained constant. It was found that the conditions for changing the regime of sliding for the modeled tribofilm can be qualitatively and quantitatively described by the criterion proposed by E. Rabinowicz [4], namely the existence of critical size of a single roughness leading to friction mode switching from smooth sliding to abrasive wear. Calculations showed that the ratio of two critical sizes of the system, obtained by direct MD simulation and the ratio of two critical length parameters calculated on the basis of Eqs. (1) or (2), give similar values. This confirms the validity of the Rabinowicz’s criterion not only for the macroscopic systems but also for atomic objects. As limitations of the proposed model of the tribofilm, note that the use of periodic boundary conditions along the loading direction has certain effects on the conditions of wear debris formation and therefore has to be studied additionally in details. Verification of Rabinowicz Criterion by Direct Molecular Dynamics Modeling 215 Acknowledgements: Investigations have been carried out with the financial support from Russian Science Foundation for Basic Research Grant No 18-508-12054 and the Fundamental Research Program of the State Academies of Sciences for 2013-2020, line of research III.23.2.4. REFERENCES 1. Popov, V.L., 2018, Adhesive wear: generalized Rabinowicz’ criteria, Facta Universitatis-Series Mechanical Engineering, 16(1), pp. 29-39. 2. Popov, V.L., Pohrt, R., 2018, Adhesive wear and particle emission: Numerical approach based on asperity- free formulation of Rabinowicz criterion, Friction, 6(3), pp. 260-273. 3. Popov, V.L., Popova, E., 2018, 60 year of Rabinowicz’ criterion for adhesive wear, Friction, 6(3), pp. 341-348. 4. Rabinowicz, E., 1958, The effect of size on the looseness of wear fragments, Wear, 2, pp. 4–8. 5. Cheng, H., Shuku, T., Thoeni, K., Temppone, P., Luding, S., Magnanimo, V., 2019, An iterative Bayesian filtering framework for fast and automated calibration of DEM models, Computer Methods in Applied Mechanics and Engineering, 350, pp. 268-294. 6. Österle, W., Dmitriev, A.I., Kloss H., 2012, Does ultra-mild wear play any role for dry friction applications, such as automotive braking?, Faraday Discussions, 156, pp. 159-171. 7. Dmitriev, A.I., Nikonov, A.Y., Österle, W., 2017, Molecular dynamics sliding simulations of amorphous Ni, Ni-P and nanocrystalline Ni films, Computational Material Science, 129, pp. 231-238. 8. Aghababaei, R., Warner, D.H., Molinari, J.-F., 2016, Critical length scale controls adhesive wear mechanisms, Nature Communications, 7, pp. 11816/1-11816/8. 9. Molinari, J.-F., Aghababaei, R., Brink, T., Frerot, L., Milanese, E., 2018, Adhesive wear mechanisms uncovered by atomistic simulations, Friction, 6(3), pp. 245-259. 10. Dmitriev, A.I., Nikonov, A.Yu., Österle, W., 2016, MD Sliding simulations of amorphous tribofilms consisting of either SiO2 or carbon, Lubricants, 4(3), pp. 1-24. 11. Tersoff, J., 1988, New empirical approach for the structure and energy of covalent systems, Physical Review B, 37, pp. 6991-7000. 12. Munetoh, S., Motooka, T., Moriguchi, K., Shintani, A., 2007, Interatomic potential for Si–O systems using Tersoff parameterization, Computational Materials Science, 39(2), 334-339. 13. Plimpton, S., 1995, Fast parallel algorithms for short-range molecular dynamics, Journal of Computational Physics, 117(1), pp. 1-19. 14. Stukowski, A., 2010, Visualization and analysis of atomistic simulation data with OVITO–the Open Visualization Tool, Modelling and Simulation in Materials Science and Engineering, 18, pp. 015012/1- 015012/7. 15. Dmitriev, A.I., Österle, W., 2014, Modelling the sliding behaviour of tribofilms forming during automotive braking: impact of loading parameters and property range of constituents, Tribology Letter, 53, pp. 337–351.