New analytical model of elastic-plastic contact for three-dimensional rough surfaces considering interaction of asperities

The contact calculation of three-dimensional real rough surfaces is the frontier field of tribology and surface science. In this study, we consider the interaction and elastic-plastic deformation characteristics of asperities and further, propose an analytical contact calculation method for rough surfaces considering the interaction of asperities. Based on the watershed algorithm, the rough surface is segmented and the asperities are reconstructed into ellipsoids. According to the height relationship between the asperities, the definition of the deformation reference height of the matrix between each couple of asperities is provided. Subsequently, the calculation formula of the substrate deformation is provided according to the local contact pressure considering the elastic-plastic deformation of the asperity, and the contact state under a specific load is determined using the iterative correction method. The results correspond with those of finite element numerical calculation and the study reveals the following: (1) compared with the results obtained without considering the asperity interaction, contact area, distance, and stiffness will be reduced by 6.6%, 19.6%, and 49.5%, respectively, when the influence of asperity interaction is considered; (2) the interaction of the asperities has the greatest influence on the surface contact distance and stiffness. Under the same load, the existence of asperity interaction will reduce the contact distance, area, and stiffness; (3) considering the interaction of the asperities, the higher asperity will bear more load, but it will simultaneously reduce the contact of the surrounding area and increase that of the distant area. The calculation method proposed in this study has the advantages of high calculation efficiency and accuracy, thus, providing the calculation basis and method for subsequent studies on service performance of rough surfaces, such as the calculation of contact stiffness and fatigue performance analysis of rough surfaces.


Introduction
The surface of mechanical parts is rough and uneven on the micro level, making the real contact area of the surface far smaller than the nominal contact area. Owing to the development of mechanical design with high speed, heavy load, and precision, the influence of the micromorphology of the workpiece has become increasingly significant [1][2][3].
To study the influence of rough surface topography on contact calculation, the characterization and modeling of rough surface topography should be considered first. Based on the statistical theory, Greenwood and Williamson first assumed that the micromorphology comprised a series of discrete hemispherical asperities with equal curvature radii, and their heights satisfied the Gaussian distribution.
Through the characterization of the relevant statistical parameters, the rough surface morphology modeling is completed. Combined with Hertz contact theory, the GW model, which is the earliest rough surface contact model [4], is proposed. The GW model provides an effective method for rough surface contact calculation. Because of its high efficiency and simplicity, the GW model has been extensively used and further improved. For example, different elastic-plastic models [5][6][7], curved surface contact models [8][9][10], asperity shapes [11,12], and surface height distribution functions [13] have been developed. Some researchers are not limited to statistical modeling methods, and have proposed some new theories and research methods. Considering the multi-scale characteristics of rough surface micro morphology, Majumdar proposed the characterization of rough surface morphology according to fractal parameters [14] and modeling according to fractal theory. Combined with the mechanical deformation model of the asperity, the fractal contact model of rough surfaces was proposed [15,16]. To avoid the uncertainty of parameter modeling, Patir et al. [17] relied on the height distribution function and spatial autocorrelation function to perform numerical modeling of rough surface topography, and the relevant deterministic modeling methods have been constantly improved to effectively achieve the numerical reconstruction modeling of deterministic surfaces [18][19][20][21]. For the deterministic surface, to realize the contact calculation and analysis promptly and stably, the researchers simplified and reconstructed the asperities on the deterministic surface, and proposed the relevant rough surface contact calculation model by combining with the mechanical deformation model of the asperity [22,23].
Generally, researchers have used different methods to model the rough surface, and obtained the contact calculation method of the rough surfaces by combining with the mechanical model of the asperity. However, in the related contact calculation, researchers regard the asperity as an independent individual for contact calculation, ignoring the interaction between the asperities. During the actual contact, although the contact area on the rough surface is always discrete, that is, the asperities are spatially independent of each other, there is still interaction between the asperities. For the same surface, different asperities share the same substrate such that the load applied on the asperity will not only lead to the deformation of the asperity itself but also affect the deformation of the substrate in the surrounding area, and subsequently affect the contact of other asperities in the surrounding area. Several relevant studies have been conducted. Based on the statistical model, Zhao et al. [24] calculated the displacement of the substrate under the uniform load according to Love formula. The downward displacement of the average height of the substrate was calculated from the total load and superimposed on the deformation of each asperity. Thereafter, the contact state was obtained according to the calculation method of the statistical model. Ciavarella et al. [25] calculated the average height change of asperities according to the deformation of Hertz pressure in the elastic half space, and further studied the influence of the interaction of asperities on the calculation results of the GW model. Vakis et al. [26] also conducted similar studies, and the calculation results indicated that the density and height distribution of asperities would affect the interaction between the asperities. Furthermore, Chandrasekar et al. [27] studied the relationship between the value of the substrate deformation and the distance, considered the influence of the density of the asperity, and proposed the definition of the substrate. The research indicated that the greater the height of the asperity, the further the contact position is from the substrate, and the less the influence on the substrate deformation. Therefore, the substrate cannot be regarded as the same height, which can further improve the accuracy of calculation. However, limited by the statistical model, the height of the substrate is defined in a www.Springer.com/journal/40544 | Friction unified way, which has some errors with the complex distribution of the actual asperities.
All of the above-mentioned models are based on the simplified calculation model for studying the interaction between the asperities, which is helpful in solving the average contact state of the rough surface efficiently. The related study mainly obtains the relationship between the contact load and distance, but it is unable to obtain the effect of the interaction between asperities on the real contact state of the rough surface, such as the change of the specific contact area. Owing to development of computer technology, the finite element theory [28,29] provides a new method for rough surface contact calculation. According to deterministic modeling, the contact state of rough surfaces can be calculated more accurately, and the interaction between asperities can be considered simultaneously. However, owing to the high-frequency and multi-scale characteristics of micro topography, the computational grid model is very large, the solution speed is low, and the iteration is difficult to converge, making it difficult to meet the efficient and stable computing requirements [30,31].
Based on the analysis of the above research status and aimed at the determined surface topography, we propose a rough surface contact calculation model considering the interaction of asperities. Based on the watershed algorithm, the rough surface is segmented and the asperities are reconstructed. Thereafter, the reconstruction model of the determined surface is obtained. Different substrate heights are determined from the relationship between the asperities. According to the local contact pressure, the substrate deformation around the asperity is calculated to consider the interaction between the asperities. The pressure distribution under the load is calculated and determined using the iterative method, and the specific contact state is obtained. The calculation results of the model are highly consistent with the finite element simulation results, with high efficiency and stability. Based on the calculation model, the influence of the interaction of asperities on the surface contact state is further studied.

Reconstruction modeling of ellipsoid asperities on rough surfaces
To calculate the influence of the interaction of asperities on the contact state simply and efficiently, researchers usually use statistical parameters to evaluate the calculation, such as the density of the asperity and height distribution function of the asperity [4][5][6][7][8]. The statistical parameters can reflect the characteristic information of the surface micro morphology and significantly simplify the calculation process and improve the calculation efficiency. However, owing to the multi-scale and complexity of the micro morphology, the surfaces with the same characteristic parameters will also present different micro-morphologies. Only the overall contact parameters, such as total contact area, total contact load, and average contact distance, can be obtained through calculation using statistical parameters; thus, it is difficult to analyze the specific contact of different surfaces, such as specific contact pressure distribution and contact area distribution. For the specific calculation of the determined surface, it still needs to rely on the finite element or numerical method, and the calculation efficiency is low.
In addition, to consider the interaction between the asperities, the distance information and height information between the asperities should be obtained. For a certain surface, it is necessary to obtain the relevant information of all the asperities involved in the contact. Therefore, to simplify the calculation process, each contact asperity is simplified and reconstructed into a regular shape, and then its curvature radius, height, position, and other parameters are calculated for further contact calculation.
To precisely describe the complex morphology characteristics of three-dimensional (3-D) rough surfaces and accurately calculate the relevant parameters of the surface asperity, the surface asperities are defined, divided, and subsequently reconstructed into regular ellipsoid asperities based on the reconstruction method proposed by Wen [32].
The contact between rough surfaces can be simplified to the contact between a rough surface and a rigid plane. For the rough surface, the asperities higher than the nominal contact plane | https://mc03.manuscriptcentral.com/friction will participate in the contact. Moreover, the rough surface has a multi-scale nature. Under different scales, different asperity parameters can be obtained. Therefore, the influence of the contact plane must be considered in the definition and division of asperity. Based on the watershed algorithm and contact plane, the definition and division of asperity can be achieved. For the rough surface morphology shown in Fig. 1, it is assumed that water is injected from each local minimum to form a catchment basin and the boundary of the catchment basin forms a watershed. If the minimum and maximum values are interchanged, each asperity will exactly correspond to each catchment basin. With the watershed as the boundary, the region division of the asperity can be achieved. Meanwhile, with the nominal contact plane as the datum plane, only the region higher than the nominal plane is calculated; thus, the definition and division of the asperities can be effectively achieved [32]. The watershed algorithm is a mature algorithm in image science, and the specific implementation process can be observed [33].
After the definition and division of all the asperities are completed, the discrete points in each asperity region are reconstructed, irregular morphology is reconstructed into ellipsoid asperities, and curvature radius, height, center position, and other parameters of all the asperities are calculated according to the reconstructed morphology.
If the coordinates of discrete points on the surface of a single micro convex region are (x, y, z), the ellipsoid equation of the asperity is as follows: The curvature radius R and ellipticity e corresponding to the asperity are calculated using the following formulas: Through the reconstruction modeling of the surface asperities, the division of the asperities can be completed, and the geometric parameters and position parameters of all the asperities can be obtained for the subsequent contact calculation and analysis. In the traditional parameter calculation, the three-point and five-point methods or parameters (such as m 2 /m 4 ) are used. Only the points near the vertex or peak point are used to calculate the parameters of the asperity; thus, it is very easy to be affected by the sampling and resolution [34]. The reconstruction method can reduce the influence of different sampling intervals on the calculation results.

Contact calculation model considering interaction of asperities
Several rough surface contact models simplify the contact between rough surfaces to the contact problem of a high number of asperities [4][5][6][7][8]. However, the assumption that the asperities are independent of each other ignores the influence of their interaction. Although the asperities are independent in space, different asperities share the same substrate; therefore, the contact load will not only lead to the deformation of the asperity itself but also affect the deformation of the substrate in the surrounding area, and subsequently affect the contact of other asperities in the surrounding area. Therefore, in this study, the rough surface is considered to comprise a substrate www.Springer.com/journal/40544 | Friction and a series of discrete asperities. The local pressure distribution at the contacted asperity is calculated using the mechanical deformation model of a single asperity, and the substrate deformation caused by the local pressure is further calculated. According to the substrate deformation, the interaction between the asperities is considered.

Elastoplastic contact model of a single asperity
According to the reconstruction of the surface asperity, the irregular asperity shape is reconstructed into a regular ellipsoid; thereafter, the contact deformation of each asperity can be calculated analytically according to the ellipsoid elastic-plastic deformation model. For a single asperity, the elastic-plastic contact model of a single ellipsoid asperity is obtained based on the elastic-plastic contact deformation theory [11]. The contact deformation between an ellipsoid and a rigid plane is illustrated in Fig. 2, where a and b denote the long and short axis radii of the contact area, respectively (a>b). ω denotes the normal deformation of the asperity, which is determined using the following formula: where z denotes the height of the asperity and d the nominal surface distance between the two surfaces. The contact load F and contact area A of the asperity can be calculated using the normal deformation as follows: where ω 1 and ω 2 denote the critical interferences of the initial yield and fully plastic deformation points, respectively; E denotes the elastic modulus; H denotes the hardness of the material; K denotes the average contact pressure coefficient, which is related to Poisson's ratio v; R denotes the asperity radius; K(e) and E(e) denote the complete elliptic integrals of the first and second types; α and λ denote the intermediate coefficients as follows: Notably, although the model in this study adopts the setting that hardness is 2.8 times the yield strength, the real situation is much more complicated than the setting. The "hardness" is not only related to material properties or geometry but also to the 222 Friction 10(2): 217-231 (2022) | https://mc03.manuscriptcentral.com/friction contact scale [35][36][37]. To reduce the influence of related variables, we adopt the simplest and unified setting in this study.

Calculation model of substrate deformation
In the deformation model of asperity, there is the assumption that the substrate of asperity does not deform, which is not consistent with the actual situation. In actual contact, under the action of a load, the asperity deforms directly, and the substrate will also deform as it is not rigid, thus, affecting the change of the surrounding area morphology. As illustrated in Fig. 3, there are three asperities, and the middle asperity is the highest. Therefore, when the rough surface contacts a rigid plane, the middle asperity is the first to contact and deform. Simultaneously, owing to the load effect, the substrate will also deform; thus, the height of the two lower asperities in the surrounding area will be reduced, affecting the subsequent contact situation.
According to Hertz contact theory, under Hertz contact pressure, the deformation formula of an elastic half space substrate is as follows: where p 0 denotes the maximum Hertz contact pressure, a the half contact width, and r the horizontal distance. According to the above theory, Ciavarella et al. [25] proposed the relevant calculation model of asperity interaction. However, there are two deficiencies in the relevant theories improved in this study. First, the above calculation method ignores the change of pressure distribution along the height direction. It can be observed from the calculation formula that p 0 denotes the maximum pressure at the contact position and the calculation target is the deformation of the substrate. The contact position and substrate are not at the same height. The maximum pressure will decrease with the increase in the height difference, and the influence on the deformation of the substrate will decrease. Chandrasekar et al. [27] also indicated this problem and proposed the related problems of substrate definition. Considering the influence of substrate height, the contact pressure at the substrate position is determined by the distance between the contact surface height and the substrate height, whereas the maximum contact pressure p z at the substrate height is determined by the maximum pressure at the contact surface position p0 and the distance h as follows: Limited by the requirement of the statistical model to use a unified benchmark for calculation, Chandrasekar et al. and the subsequent researchers defined the lowest point between all the asperities as the substrate height. However, for the deterministic surface, the substrate height between the asperities is not exactly the same. For two asperities, considering the interaction, the height of the substrate should be the lowest point between the two asperities. As illustrated in Fig. 4, when the influence of asperity 2 on the asperities in the surrounding area is considered, the substrate height is clearly not the same for the different asperities on the left and right.
In addition, for the substrate, it is reasonable that the deformation is in the elastic range, but for the contact of the asperity, it is very easy to produce elastic-plastic deformation. The deformation of the substrate is calculated using the pressure calculated from the contact model of asperity. Accordingly, it is unreasonable to calculate the deformation of the substrate by Hertz contact pressure and Hertz www.Springer.com/journal/40544 | Friction contact half width. Therefore, in this study, we used the maximum contact pressure and contact half width of the elastic-plastic contact model of asperity for calculation.
To compare the influence of different substrate heights and the elastic-plastic deformation of the asperity on the calculation results, the profile illustrated in Fig. 4 is selected for calculation and comparison. In Fig. 4, the heights of asperities 1 and 3 are equal. The highest value of asperity 2 is 12 μm, whereas the radius of curvature is 4 μm. The height of the lowest point between asperities 1 and 2 is -4 μm, and the height of the lowest point between asperities 2 and 3 is -12 μm. According to the morphology data, an isotropic elastic-plastic contact model is established by importing the data into the commercial finite element software ABAQUS. As illustrated in Fig. 5, the upper surface is simplified as a rigid plane and the lower surface is the target surface. A fixed constraint is applied at the bottom of the target surface substrate, and a symmetrical constraint is applied to both sides of the substrate. Four-node linear element (CPE4) is used to mesh the substrate, and the mesh near the contact surface is refined. Based on the finite element elastic-plastic model, the contact between the profile and a rigid plane is calculated. The actual contact only occurs between the plane and the highest asperity 2. The displacement of the top points of asperities 1 and 3 is recorded. In addition, let the substrate height between asperities 1 and 2 be -4 μm, and the substrate height between asperities 2 and 3 be -12 μm. According to formulas (5)- (14), the displacement of asperities 1 and 3 under the influence of asperity 2 is calculated. The results of pure elastic Hertz model are also calculated, and compared as illustrated in Fig. 6.
It can be observed from the Hertz results in the figure that if the Hertz contact pressure and the contact half width are used, the influence of the interaction of asperities will be significantly underestimated. The results of the finite element model indicate that different substrate heights can affect the value of asperity interaction. For the model presented in Fig. 3, the other parameters of asperities 1 and 3 are similar, only the height of the substrate is different. When asperity 2 is forced to deform, the deformation of asperity 1 under the interaction of asperities is larger than that of  | https://mc03.manuscriptcentral.com/friction asperity 3. As the load increases, the gap between the two increases, that is, the influence of different substrate heights increases, and the model in this study exhibits a similar trend. In the finite element results presented in Fig. 6, the maximum difference of deformation between the two asperities is 15.8%. For the model in this study, the calculation result is 16.7%. As the results of the finite element model are constrained by boundary conditions, mesh generation, and other factors, it is difficult for the analytical formula to describe the deformation accurately. Consequently, there is a specific error between the two; however, generally, the results of this study and the finite element model exhibit good consistency.
Therefore, according to the elastic-plastic contact model of asperity, the contact pressure and half contact width are calculated, and the influence of different substrate heights is considered. For the calculation of the interaction between two asperities, the height of the substrate is determined by the lowest point between the two asperities to improve the accuracy of the calculation of the interaction. Considering the horizontal distance between the asperities r and the vertical distance between the height of the contact surface and height of the substrate h, the formula for calculating the substrate deformation under the influence of asperity interaction is modified as follows: where i denotes the serial number of the asperity deformed under a load and n the serial number of the asperity affected by the interaction.

Contact model of a rough surface
The contact between two rough surfaces can be replaced by the contact between an equivalent rough surface and a rigid smooth surface such that only the contact between a rough surface and a rigid smooth surface is considered. According to the reconstruction method of surface asperity, the curvature radius, height, position, and other parameters of all asperities can be obtained by reconstructing the determined surface topography. Without considering the interaction of asperities, the normal deformation of a single asperity can be determined according to the average surface distance d and the geometric parameters of the asperity. Thereafter, the contact state of the asperity can be calculated according to the mechanics model of a single asperity, and the contact state of the entire surface is obtained by calculating the contact state of all asperities and superimposing them.
However, considering the interaction of asperities, all the asperities involved in the contact will affect the asperities in the surrounding area by reducing the height of the substrate. Therefore, for each asperity in contact, its maximum height z requires correction. The amount of correction is the superposition of the deformation of the substrate under the action of all other contacted asperities as follows: where z n ' denotes the modified height of the asperity and s the total number of asperities in contact. The normal deformation is calculated according to the modified height of the asperity; thus, the influence of the interaction of asperities on the calculation results can be included. Considering that the calculation of the substrate deformation requires specific pressure distribution and that the change of the substrate deformation will also affect the specific pressure distribution, an iterative correction method is adopted to calculate the accurate contact state. The specific calculation process is as follows: (1) Input the surface mechanical parameters (E, H, v), provide the total load F, set the error threshold e, and assume the average height of the surface to be the reference surface height.
(2) According to the reference surface, the asperities on the rough surface are reconstructed, and the height z, radius of curvature R, along with position coordinates of the center point (x, y) of all contacted asperities are obtained.
www.Springer.com/journal/40544 | Friction (3) According to the nominal surface distance d, the curvature radius R, and the height of the asperity, the load on a single asperity F i is calculated according to formulas (5)- (7), and the total load F t is obtained through superposition. Adjust the value of d to make the total load equal to F. Assume the height of the reference surface to be the height of the contact surface d.
(4) According to the load F i , position of the center point of all the contacted asperities, and height of the contact surface, the correction of the maximum height of all contacted asperities is obtained using formulas (13)- (15).
(5) Correct the height of the asperities, repeat steps 2 and 3, and assume the height of the new contact surface to be d 1 . If d 1 -d < e, the final calculation result will be obtained; otherwise, repeat steps 4 and 5.

Results and discussions
To verify the correctness of the method proposed in this study, the contact calculation and verification analysis are performed based on measured data of rough surface topography, and the results are compared with those of Zhao's [24], KE [6], Chandrasekar's [27], and the finite element models. Notably, the KE model is discontinuous in the critical stage, but the calculation results still have high accuracy in the small deformation range [8,38].
Three groups of ground workpiece surfaces were measured using the white light interferometer Wyko NT9100; the sampling interval was 0.99 μm, total area of sampling area was A n = L x L y = 0.1 0.1 = 0.01 mm 2 , and surface morphology measured was as illustrated in Fig. 7.
According to the surface morphology, the surface asperity parameters are calculated as summarized in Table 1, and the surface mechanical parameters of the workpiece are listed in Table 2. The statistical parameters provided in Table 1 are for the calculation requirements of a statistical model; however, the calculation of the model in this study does not require statistical parameters.
First, the influence of the interaction of asperities on the contact distance is studied. Based on the      In the above models, the KE model does not consider the influence of asperity interactions. Zhao's model is based on a statistical model, whereas Chandrasekar's model is based on KE model, but both of them consider the influence of asperity interactions. It can be observed from Figs. 8-10 that compared with the model without the consideration of the interaction of asperities, under the same load, the dimensionless contact distance predicted by Chandrasekar's model and the model in this paper decreases, whereas the prediction result of Zhao's model is the opposite when the interaction of asperities is considered. As the load increases, the influence of asperity interactions on the prediction result of contact distance also increases. It can be observed from Fig. 8 that the contact distance between the statistical (KE and Chandrasekar) and deterministic models (finite element and current models) at the initial contact position is significantly different. This is because the statistical model usually assumes that the surface height distribution obeys a certain probability density function, such as Gaussian distribution; thus, the statistical model assumes that the two surfaces start to contact at the position where the contact distance is 3σ. However, for the data of surface 1, the maximum surface height is approximately 4σ, and the contact occurs at this position, indicating that the statistical model may have a large error in the calculation and prediction of some surfaces. Generally, the calculation results of the current and finite element models correspond, verifying the accuracy of the calculation method in this our study.
To make a more detailed comparison between the calculation results of the current and finite element models, and simultaneously study the influence of the interaction of asperities on the contact calculation results based on the data of surface 2, the contact calculation is performed according to the finite element and simplified models without considering the interaction of asperities and the model in this paper. The correlation rules of dimensionless contact load (F/A n /E), dimensionless contact distance (d/σ), and dimensionless contact area (A/A n ) are illustrated in Figs. 11 and 12. Meanwhile, the contact distance is mainly related to the contact stiffness in the www.Springer.com/journal/40544 | Friction engineering performance parameters. Therefore, the correlation rules of the dimensionless contact load and dimensionless contact stiffness (σK/A n /E) are illustrated in Fig. 13.
It can be observed from Figs. 11-13 that the calculation results considering the interaction of   asperities are closer to the prediction results of the finite element model. Based on the results of the finite element model, under the same load, the maximum deviation of the calculation results for contact area, distance, and stiffness are 10.7%, 22.0%, and 93%, respectively, without considering the interaction of asperities. When the interaction of asperities is considered, the maximum deviation of the calculation results for contact area, distance, and stiffness are only 4.2%, 4.8%, and 14.2%, respectively. Compared with the calculation results without considering the interaction of asperities, under the same load, the interaction of asperities reduces the contact area, clearance, and stiffness. In the case of a low load, the results of the three models are almost the same. As the load increases, the influence of the interaction of asperities increases significantly. When the dimensionless load is 5× 10 -3 , the predicted contact area, clearance, and stiffness are reduced by 6.6%, 19.6%, and 49.5%, respectively, owing to the interaction of asperities.
The results indicate that the interaction of asperities has little effect on the contact area, but has a more significant effect on the contact distance and stiffness.
The total contact load, area, and distance represent the change trend of the contact state of the entire surface. Therefore, the above results mainly compare and study the influence of the interaction of asperities on the contact state of the entire surface. In this paper, the specific contact area of the rough surface is further studied and analyzed. As KE, Zhao's, and Chandrasekar's models are calculation models based on statistical characteristic parameters, it is impossible to obtain the specific contact area. Therefore, the calculation results of the finite element model and the model in this paper are mainly compared. Simultaneously, to compare the influence of the interaction of asperities on the calculation results, according to the model in this paper, the calculation results without considering the interaction of asperities are also included in the comparison. Based on the data of surface 2, when the dimensionless load is 0.0005, 0.001, 0.005, and 0.01 (low to high loads), the contact stress nephograms on the contact area obtained using the three calculation methods are presented in Figs. 14-16.    As can be observed, Fig. 14 presents the calculation result without considering the interaction of asperities, Fig. 15 presents the result after considering the interaction of asperities, and Fig. 16 shows the calculation result of the finite element model. Generally, the contact area calculated by the current model is very close to that of the finite element model, further verifying the reliability of the calculation method proposed in this study. For the grinding surface 2 illustrated in Fig. 7(b), there are dense peaks in the three areas on the left, middle, and right, and the height distribution is low on the left and high on the right. During the contact process from low to high loads, the surface starts to contact from the highest peak on the right side and gradually expands to the peaks in the middle and left areas. The influence of the interaction of asperities is studied by contrast. It can be observed from (a) and (b) in Figs. 14 and 15 that under a low load (0.0005 or 0.001), the surface mainly contacts at the highest peak on the right side and the middle part, and the effect of the interaction of asperities is not significant. The contact area calculated by the two models is almost the same. Owing to the increase of the load, the interaction of asperities has an effect on the distribution of the contact area. As illustrated in (c) and (d) in Figs. 14 and 15, considering the interaction of asperities, the contact area in the middle decreases, www.Springer.com/journal/40544 | Friction whereas the contact area on the left increases; that is, under the interaction of asperities, the contact area near the high asperity decreases, whereas the contact area further from the high asperity increases. This is because the high asperity bears a large load, resulting in a significant downward movement of the substrate in the surrounding area and a drop in the height of the asperities. The area further from the high asperity is less affected; thus, there is little change in the height of the asperities, increasing the contact probability. It can be observed from the calculation results of the finite element model that the calculation results under the consideration of the interaction of asperities and the finite element results are more consistent, proving the accuracy of the calculation method in this paper.

Conclusions
1) In this study, we propose a new analytical model of elastic-plastic contact for 3D rough surfaces considering the interaction of asperities and the elastic-plastic deformation of the asperities. A new definition of the height of substrate is proposed, thus, improving the accuracy of predicting the interaction effect of asperities. The calculation method proposed in this study has the advantages of high calculation efficiency and accuracy.
2) The existence of asperity interaction has the most significant influence on the surface contact distance and stiffness. Under the same load, the existence of asperity interaction will reduce the contact distance, stiffness, and area. For the specific contact area, higher asperities will bear more load when the interaction of asperities is considered; however, the contact of the surrounding area will be reduced while the contact of the far area is increased. Owing to the increase of the load, the influence of the asperity interaction will increase significantly.
3) In the calculation of the rough surface contact, neglecting the interaction of asperities will result in errors. The results indicate that the predicted contact area, distance, and stiffness decrease by 6.6%, 19.6%, and 49.5%, respectively, owing to the interaction of asperities. 4) Using the method proposed in this study, contact calculation and analysis of rough surfaces can be performed accurately and efficiently, providing a scientific basis for the research work in the fields of contact stiffness calculation and fatigue performance analysis of rough surfaces.