# A Fractal Contact Model for Rough Surfaces considering the Variation of Critical Asperity Levels

Authors: Yong Liu; Haodong Guo; Zhiying Chen; Kunying Ding; Min Dan; Baichun Li; Fangchao Yan
Journal: Advances in Materials Science and Engineering
 (2022)

---

## Abstract

A contact model for rough surfaces based on the fractal theory is proposed in the present work. Firstly, the deformation of the material is divided into four stages: elastic deformation, the first elastoplastic deformation, the second elastoplastic deformation, and full plastic deformation. And the variation of material hardness is considered when analyzing the contact characteristics of a single asperity within the first and second elastoplastic deformation stages. Secondly, the size distribution function of contact spots at different frequency levels is derived. And the expressions of asperity critical frequency levels are rederived. Lastly, the feasibility and credibility of the proposed model are verified by comparison with other contact models and experimental data. The results show that when the variation of the material hardness is considered, the contact area of a single asperity in the first elastoplastic deformation stage becomes larger, while the contact area of a single asperity in the second elastoplastic deformation stage becomes smaller. Moreover, the critical asperity frequency levels of the rough surface are not constant, but the variables are related to the total real contact area of the rough surface and decrease as the real contact area increases. The proposed model is a modification and improvement of the existing fractal contact models, which can lead to a more accurate relationship between the contact load and the total real contact area of the rough surface.


---

## Body

## 1. Introduction

The aero-engine external piping system is mainly used for the transmission of fuel, lubricating oil, hydraulic oil, and air and is an important part of the external accessory device [1–3]. Hundreds of pipelines are installed on an aero-engine. Due to their heat resistance and corrosion resistance, tube connectors in the form of metal-to-metal seals (i.e., no dedicated sealing component used) are often used for the connection between pipelines [4]. Tube connectors are the weakest link of the pipeline system’s sealing performance, and tube connector sealing failure has become one of the pipeline system’s main failure modes [4–6]. Once the metal-to-metal seals are out of work, leakage will be formed. Both stability and reliability of the aero-engine will be affected.The metal-to-metal seal is performed by a direct-metal/metal-tight contact of rough surfaces as shown in Figure1. Although it is simple in structure, the sealing behavior of a metal-to-metal seal is affected by a variety of factors [7–11], among which the surface topography, which usually has a microstructure given by machining processes, is thought to be one of the most important factors [12–14]. Clearly, engineering materials are known to have rough surfaces, and the full control of surface topography at all scales during manufacturing processes is still out of reach [9]. When two rough surfaces come into contact, the topography of the surface leads to imperfect contact and makes the real contact area only a fraction of the nominal area [4, 15]. Noncontact areas communicate with each other to form leakage channels. Thus, the research on leakage mechanism and sealing performance of metal-to-metal seals involves the resulting geometry of contact between the rough surfaces. Accurately characterizing the contact state and extracting the relationship between the contact load and the real contact area of the two sealing surfaces are necessary for the study of the leakage mechanism and sealing performance of metal-to-metal seals [6, 16–18].Figure 1
A schematic showing the flow of sealing medium through the sealing interface of metal-to-metal seal.The research on the contact problem of rough surfaces started as early as the Coulomb’s friction proposed by the French engineer Coulomb [19] in 1781. Hertz [20] gave an analytical solution to the contact problem of frictional elastomers in 1882, which opened up the study of modern contact mechanics. Subsequently, researchers have proposed a variety of contact models to describe the contact behavior of rough surfaces. These models mainly include numerical contact models, statistical contact models, and fractal contact models [15, 21–25]. The numerical contact models generally use digital technology (i.e., SEM, AFM, and STM) to obtain the specific parameters of the surface topography and use the finite element method to simulate the contact behavior of two rough surfaces. The advantage of numerical contact models is that they can obtain simulation results that are closer to the actual situation. But they might require, in general, a large and dense grid, and the computational efficiency might be unacceptable [18]. Statistical contact models use statistical parameters to characterize the rough surfaces. The GW contact model proposed by Greenwood [26] in 1966 is a typical representative of the statistical contact models, and subsequent scholars have improved the GW model from different aspects [27, 28]. The advantage of statistical contact models is that their expressions are simple and clear, which greatly simplifies the derivation of the contact equation between rough surfaces and is conducive to rapid contact analysis between rough surfaces. However, the statistical contact models all simplified the characterizations of the surface topography to various degrees, which consequently caused a large error between the calculated and the actual results. Engineering surfaces are found to have fractal characteristics [29, 30]. When fractal parameters instead of statistic parameters are adopted to characterize rough surfaces, there is no need to make too many assumptions about the surface topography, and the fractal characteristics of rough surfaces can be preserved [31–35]. As a result, although there are some skeptical opinions on the fractal approaches [36], contact models based on fractal theory have received more and more attention from researchers. Bhushan and Bhushan [37] proposed one of the first fractal contact models (MB model) based on the fractal theory and Weierstrass-Mandelbrot (WM) function in 1991. Komvopoulos and Komvopoulos [38, 39] believed that the truncated area of an asperity should not be equal to its real contact area. They presented an improved contact model (WK model) in which the deformation mode of an asperity is divided into three stages: elastic deformation, elastoplastic deformation, and full plastic deformation. Considering the elastoplastic deformation of asperities, Komvopoulos and Komvopoulos [40] established a 3D fractal contact model (YK model) for rough surfaces, which is suitable for both isotropic and anisotropic surfaces. The MB model predicts that the asperity first deforms plastically and then elastically during the loading process. This conclusion is contrary to classical contact mechanics and contradicts people’s intuitive feelings. However, it has never been challenged until Etsion and Etsion [41] published their work (ME model) in 2007. In the ME model, the concept of critical asperity frequency level was proposed for the first time. Subsequently, the ME model was developed into a complete contact model by Huang and Huang (MH model) [42]. Yuan et al. [43] proposed a revised contact model (YC model) based on the MB model and the ME model. Based on the YC model, Yuan et al. have successively proposed a loading-unloading contact model for rough surfaces [35, 44], as well as a normal contact stiffness model for joint surfaces [45].The YC model has solved several deficiencies of the MB model and played a significant role in promoting the research of contact modeling based on the fractal theory. However, the elastic critical frequency level of the YC model is derived on the condition that the asperity height is not greater than itself critical interference. The consequence of this condition in the YC model is that the variation of the critical frequency levels with the total real contact area of the rough surface is not taken into account. In addition, most fractal contact models do not consider the variation of material hardness with deformation. The variation of material hardness has a vital influence on the mechanical properties of the asperities on the rough surface [46], and it hence affects the contact characteristics of the entire rough surface.Therefore, the main purpose of this paper is to propose a fractal contact model considering the variation of the critical asperity levels as well as the variation of the material hardness to extract the relationship between the contact load and the total real contact area of the rough surface. The present work will lay the foundation for the subsequent analysis of the sealing performance of static metal seals.


## 2. Theory and Method

### 2.1. The WM Function

An engineering surface profile often appears random, multiscale, and disordered. The mathematical properties of such a profile are as follows: it is continuous, nondifferentiable, and statistically self-affine [37]. It is found that the WM function satisfies all these properties and is often used to create 2D rough surface profiles. The WM function is given as [37](1)zx=GD−1∑n=nmin∞cos2πγnxγ2−Dn,1<D<2,γ>1,where zx represents the height of the profile along with the x direction. D and G are fractal parameters, representing fractal dimension and characteristic length scale (roughness parameter), respectively. γn determines the frequency spectrum of the surface roughness, and n indicates the frequency level of the asperities. Generally, in order to meet the requirements for high spectral density and phase randomization, γ is taken as 1.5. nmin indicates the low cut-off frequency of the profile.


### 2.2. The Contact of a Single Asperity

The idea of contact modeling based on the fractal theory is to obtain the “contact area-contact load” relationship of the entire rough surface through integration based on the “contact area-contact load” relationship of a single asperity. Therefore, the first step is to obtain the “contact area-contact load” relationship of a single asperity.


#### 2.2.1. The Existing Elastic Microcontacts

According to the WM function, the profile of the fractal asperity at a certain leveln before deformation is [37](2)znx=GD−1ln2−Dcosπxln,−ln2<x<ln2,where ln represents the length scale of the asperity of level n; that is, the base diameter of the asperity of level n, ln=1/γn.The asperity height is [37](3)δn=zn0=GD−1ln2−D.The contact between two rough surfaces can be simplified to an equivalent rough surface in contact with a rigid flat surface [42]. Figure 2 shows the relationship between the geometric parameters when a fractal asperity is in contact with a rigid smooth surface. According to Figure 2, the interference ωn of the asperity with a length scale ln is(4)ωn=δn−d=GD−1ln2−D−d,where d is the separation distance between the rigid flat surface and the equivalent rough surface. If the surface height zx,y follows the Gaussian distribution, separation distance d can be obtained by the relationship [37](5)ArAa=12π∫d/σ∞exp−x2/2dx=12erfcd2σ,where Ar, Aaare the total real contact area and the nominal contact area, respectively. erfc(x) represents the complementary error function.Figure 2
The contact between a fractal asperity at leveln and a rigid smooth plane, where ln, lnr, and lnt are the base diameter, the real contact diameter, and the truncation diameter, respectively [43].According to the Hertz theory [37], the contact area and contact load of the asperity within the elastic deformation stage can be obtained as(6)an=πRnωn,(7)fne=43ERn1/2ωn3/2,where Rn is the radius of curvature at the asperity summit, Rn=1/d2zn/dx2x=0=lnD/π2GD−1. E is the equivalent elastic modulus, 1/E=1−ν12/E1+1−ν22/E2, and E1, E2, ν1, and ν2 are elastic modulus and Poisson’s ratios of the two rough surfaces. According to equations (6) and (7), the relationship between contact load and real contact area of a single asperity within the elastic deformation stage is(8)fne=4Eπ1/2GD−13lnDan3/2.


#### 2.2.2. The Existing Condition for Initial Yield

The internal stress of the asperity increases with the increase of the contact load or contact interference. An initial yield point will eventually be generated inside the asperity due to excessive stress. The interference corresponding to the initial yield point is termed as the critical interferenceωnec, and is given as [47](9)ωnec=πKH2E2Rn.where K is the hardness factor of the softer material and is given by K=0.454+0.41ν. H denotes the hardness of the softer material and is given by H=2.8Y, and Y is the yield strength of the softer material. The elastic critical contact area of the asperity anec can be obtained as [43](10)anec=πRnωnec=1πKH2E2lnDGD−12.According to equations (7), (9), and (10), we can get the expression of the elastic critical contact load of the asperity fnec, which is(11)fnec=23KHanec.


#### 2.2.3. The Revised Elastoplastic Microcontacts

As the load or interference of the asperity increases, the plastic part inside the asperity will gradually expand to the contact surface, and an annular plastic part will be formed on the contact surface, while the rest of the contact surface surrounded by the plastic region remains elastic deformation. This transition stage, with the transition interference ratio of1<ωn/ωnec≤6(anec<an≤anepc), is known as the first elastoplastic deformation stage. The first elastoplastic critical contact area anepc and the first elastoplastic contact load fnep1 of the asperity are derived by Yuan et al. [43] based on the ME model, which are(12)anepc=7.1197anec,(13)fnep1=23KH×1.1282anec−0.2544an1.2544.In the first elastoplastic deformation stage, the hardness of the material will change with the deformation rather than remaining a constant value [25]. According to equations (12)∼(13), we assume that the material hardness within the first elastoplastic deformation stage HG1an satisfies the following relationship:(14)HG1an=c11Yananecc12,anec<an≤anepc,where c11, c12are the parameters need to be solved.Equation (14) should satisfy two boundary conditions:(15)HG1anec=pemanec,(16)HG1anepc=pepm1anepc,where peman is the average contact pressure of the asperity in the elastic deformation stage, which is given by peman=fne/an. Hence, pemanec can be obtained as pemanec=2/3KH asperity in the first elastoplastic deformation stage and is given by pepm1an=fnep1/an. Substituting equation (14) and equation pemanec=KH into equation (15) yields c11Y=2/3KH. Thus, the parameter c11 can be obtained as(17)c11=2815K.Substituting equations (12)–(14) and (17) into equation (16) yields(18)23KH×7.1197c12=23KH×1.1282anec−0.25447.1197anec0.2544.And parameterc12 can be obtained as(19)c12=ln1.1282×7.11970.2544ln7.1197.Therefore, the contact load of a single asperity in the first elastoplastic deformation stage is revised as(20)fnep1′=HG1anep1·anep1=2815KYanec−c12anc12+1.If the contact load or interference further increases, the plastic part of the asperity gradually expands to envelop the shrinking elastic core. According to Etsion and Etsion [48], the transition interference ratio of this stage is 6<ωn/ωnec≤110(anepc<an≤anpc). Some researchers named this stage the second elastoplastic deformation stage. The second elastoplastic critical contact area anpc and the second elastoplastic contact load fnep2 of the asperity are also given by Yuan et al. [43]:(21)anpc=205.3827anec,(22)fnep2=23KH×1.4988anec−0.1021an1.1021.The hardness of the material will also change with the deformation in the second elastoplastic deformation stage [25]. Similarly, we assume that the material hardness within the second elastoplastic deformation stage HG2an satisfies the following relationship:(23)HG2an=c21Yananecc22,anepc<an≤anpc,where c21 and c22 are the two parameters need to be solved.And equation (23) should also satisfy two boundary conditions:(24)HG2anepc=pem1anepc,(25)HG2anpc=pepm2anpc,where pepm2an is the average contact pressure of the asperity in the second elastoplastic deformation stage and is given by pepm2an=fnep2/an.Substituting equations (12), (13), and (23) into equation (24) yields(26)c21Y7.1197c22=23KH×1.1282×7.11970.2544.Substituting equations (21), (22), and (23) into equation (25) yields(27)c21Y205.3827c22=23KH×1.4988×205.38270.1021.According to equations (26) and (27), parameter c22 can be expressed as(28)c22=ln1.4988+0.1021ln205.3827−ln1.1282−0.2544ln7.1197ln205.3827−ln7.1197.And parameterc21 can be obtained as(29)c21=K×2815×1.1282×7.11970.2544−c22.Therefore, the contact load of a single asperity in the second elastoplastic deformation stage is revised as(30)fnep2′=HG2anep2·anep2=K×2815×1.1282×7.11970.2544−c22Yanec−c22an1+c22.


#### 2.2.4. The Existing Plastic Microcontacts

Whenωn/ωnec>110(an>anpc), the asperity undergoes full plastic deformation. The contact area and contact load of the asperity in this deformation stage are [43](31)an=2πRnωn,fnp=Han.For the sake of clarity, Figure3 shows the deformation law of a single asperity with the fractal asperity frequency level n, as well as the relationship between each critical contact area.Figure 3
The diagram of the deformation law of a single asperity and the relationship between each critical contact area.


### 2.3. The Revised Size Distribution Function of Contact Spots

In the MB model, the size distribution function of contact spots is defined as [37](32)na=D2alD2a−D+22,0<a<al,1<D<2,where al denotes the largest contact area of asperity. And the total real contact area is [37](33)Ar=∫0alnaada=D2−Dal.According to equation (2), the period of asperity at frequency level n is Tn=2π/2πγn, and the period of asperity at frequency level n+1 is Tn+1=2π/2πγn+1. Therefore, Tn+1/Tn=1/γ. Hence, we can obtain the relationship between the size distribution function of the asperity at frequency level n and the size distribution function of the asperity at frequency level n+1, which is nn+1a=γnna. Letting the size distribution function of the asperity at the initial frequency level nmin be nnmina=Mna, then the contact spot size distribution function at each frequency level can be obtained as(34)nna=Mγn−nminna,nmin≤n≤nmax,where M is a parameter that needs to be solved. According to equations (33) and (34), the total real contact area of asperities at each frequency level is(35)Anr=∫0anlnnaada=MD2−Dγn−nminanl,where anl is the largest contact area of contact spots at frequency level n. According to the geometry relationship shown in Figure 2, we assume that the largest contact spot of frequency level nmin is the largest contact spot of the entire rough surface, i.e., anminl=al. And we assume that the relationship between the largest contact spots of two adjacent frequency levels is an+1l/anl=1/γ2. Then, the largest contact area of the asperity at any frequency level would be anl=1/γ2n−nminal. Thus equation (35) can be reexpressed as(36)Anr=Mγn−nminD2−Dal=Mγn−nminAr.And the relationship betweenanl and Ar can be obtained as(37)anl=2−DD1γn−nmin2Ar.The parameterM can be determined according to the relationship Anminr+Anmin+1r+⋯+Anmaxr=Ar. Substituting equation (36) into it yields M=1−1/γ/1−1/γnmax−nmin+1.


### 2.4. The Revised Critical Asperity Frequency Levels

For a certain frequency leveln, the relationship between the interference ωn and the asperity height δn is ωn≤δn [46]. If the height of an asperity is less than its elastic critical interference, inelastic deformation occurs. Thus, the condition of elastic deformation for an asperity given by the YC model is δn≤ωnec [43]. However, according to Etsion and Etsion [41], the condition for elastic deformation of an asperity should be ωn≤ωnec, which indicates that the asperity occurs elastic deformation when the interference is not greater than its critical interference. Etsion and Etison [41] further derived the inequality used to solve the elastic critical length scale lec(or the elastic critical asperity frequency level nec), which is 1−cosπlnt/2ln≤KH/2E⋅ln/GD−12. Unfortunately, we cannot get the analytical solution of lec according to this inequality, nor can we get an analytical solution of nec according to the relationship ln=1/γn.Now, we reanalyze the conditions for judging the elastic deformation of the asperity. At a mean surface separation distanced, if the largest interference of the asperity at frequency level n is not greater than its critical interference, all of the contact asperities at frequency level n deform elastically. Thus, the condition we give to judge the elastic deformation of the asperity is(38)ωnl≤ωnec,where ωnl is the largest interference of the contact spot at frequency level n. According to equation (6), equation (39) can be re-written as(39)anl≤anec.Substituting equations (10) and (37) into the above inequality, the elastic critical frequency level can be obtained from the following equation:(40)2−DD1γnec−nmin2Ar=1πKH2E21γ2necD1G2D−1.According to the above inequality, the elastic critical frequency levelnec is not only related to the material parameters H,E and the topography parameters D,G but also related to the total real contact area of the two rough surfaces Ar. Similarly, the first elastoplastic critical frequency level nepc and the second elastoplastic critical frequency level npc can be obtained from the following two equations:(41)2−DD1γnepc−nmin2Ar=7.11971πKH2E21γ2nepcD1G2D−1,2−DD1γnpc−nmin2Ar=205.38271πKH2E21γ2npcD1G2D−1.According to equations (41)∼(43), the equations used to calculate the critical frequency levels nec,nepc,npc are not only related to the material parameters H,E and topography parameters D,G, but also related to the total real contact area of the rough surface Ar. As the given parameter Ar changes, so do the calculation results of the critical frequency levels nec,nepc,npc. That is, the critical frequency levels nec,nepc,npc variate with the change of the Ar value. However, in the equations for calculating the critical frequency levels given by the YC model, the critical frequency levels are only related to the material parameters and topography parameters, and the influence of the total real contact area Ar on the critical frequency levels is not considered.


### 2.5. The Revised Real Contact Area and Contact Load of the Rough Surface

When the asperity frequency level ranges fromnmin to nec, i.e., nmin≤n≤nec, only elastic deformation takes place. And the real contact area Ar1 and contact load Fr1 are obtained as(42)Ar1=Are=∑n=nminnec∫0anlMγn−nminnaada,Fr1=Fre=∑n=nminnec∫0anlfneMγn−nminnada.When the asperity frequency level belongs tonec<n≤nepc, elastic or the first elastoplastic deformation can take place. The real contact area Ar2 is obtained as(43)Ar2=Are+Arep1=∑n=nec+1nepc∫0anecMγn−nminnaada+∑n=nec+1nepc∫anecanlMγn−nminnaada.And the contact loadFr2 is obtained as(44)Fr2=Fre+Frep1=∑n=nec+1nepc∫0anecfneMγn−nminnada+∑n=nec+1nepc∫anecanlfnep1′Mγn−nminnada.When the asperity frequency level isnepc<n≤npc, elastic deformation, the first elastoplastic deformation, or the second elastoplastic deformation will occur. The real contact area Ar3 is evaluated as(45)Ar3=Are+Arep1+Arep2=∑n=nepc+1npc∫0anecMγn−nminnaada+∑n=nepc+1npc∫anecanepcMγn−nminnaada+∑n=nepc+1npc∫anepcanlMγn−nminnaada,.And the contact loadFr3 is evaluated as(46)Fr3=Fre+Frep1+Frep2=∑n=nepc+1npc∫0anecfneMγn−nminnada+∑n=nepc+1npc∫anecanepcfnep1′Mγn−nminnada+∑n=nepc+1npc∫anepcanlfnep2′Mγn−nminnada,.When the asperity frequency level belongs tonpc<n≤nmax, all of the four types of deformation can take place. The real contact area Ar4 is calculated as(47)Ar4=Are+Arep1+Arep2+Arp=∑n=npc+1nmax∫0anecMγn−nminnaada+∑n=npc+1nmax∫anecanepcMγn−nminnaada,+∑n=npc+1nmax∫anepcanpcMγn−nminnaada+∑n=npc+1nmax∫anpcanlMγn−nminnaada.And the contact loadFr4 is estimated as(48)Fr4=Fre+Frep1+Frep2+Frp=∑n=npc+1nmax∫0anecfneMγn−nminnada+∑n=npc+1nmax∫anecanepcfnep1′Mγn−nminnada+∑n=npc+1nmax∫anepcanpcfnep2′Mγn−nminnada+∑n=npc+1nmax∫anpcanlfnpMγn−nminnada.The detailed calculation results of the real contact area and contact load for the rough surface are provided in Appendix A. The total real contact area and contact load of all asperity levels are estimated by(49)Ar=Ar1+Ar2+Ar3+Ar4.Fr=Fr1+Fr2+Fr3+Fr4.The nondimensional forms ofAr and Fr are defined as [43](50)Ar∗=ArAa,Fr∗=FrAaE,where the nominal contact area Aa is given by Aa=L2, L=1/γnmin. For the sake of clarity, the flowchart of the solution procedure for extracting the relationship between the total real contact area and contact load is shown in Figure 4.Figure 4
Flowchart for extracting the relationship between the total real contact area and contact load.


## 2.1. The WM Function

An engineering surface profile often appears random, multiscale, and disordered. The mathematical properties of such a profile are as follows: it is continuous, nondifferentiable, and statistically self-affine [37]. It is found that the WM function satisfies all these properties and is often used to create 2D rough surface profiles. The WM function is given as [37](1)zx=GD−1∑n=nmin∞cos2πγnxγ2−Dn,1<D<2,γ>1,where zx represents the height of the profile along with the x direction. D and G are fractal parameters, representing fractal dimension and characteristic length scale (roughness parameter), respectively. γn determines the frequency spectrum of the surface roughness, and n indicates the frequency level of the asperities. Generally, in order to meet the requirements for high spectral density and phase randomization, γ is taken as 1.5. nmin indicates the low cut-off frequency of the profile.


## 2.2. The Contact of a Single Asperity

The idea of contact modeling based on the fractal theory is to obtain the “contact area-contact load” relationship of the entire rough surface through integration based on the “contact area-contact load” relationship of a single asperity. Therefore, the first step is to obtain the “contact area-contact load” relationship of a single asperity.


### 2.2.1. The Existing Elastic Microcontacts

According to the WM function, the profile of the fractal asperity at a certain leveln before deformation is [37](2)znx=GD−1ln2−Dcosπxln,−ln2<x<ln2,where ln represents the length scale of the asperity of level n; that is, the base diameter of the asperity of level n, ln=1/γn.The asperity height is [37](3)δn=zn0=GD−1ln2−D.The contact between two rough surfaces can be simplified to an equivalent rough surface in contact with a rigid flat surface [42]. Figure 2 shows the relationship between the geometric parameters when a fractal asperity is in contact with a rigid smooth surface. According to Figure 2, the interference ωn of the asperity with a length scale ln is(4)ωn=δn−d=GD−1ln2−D−d,where d is the separation distance between the rigid flat surface and the equivalent rough surface. If the surface height zx,y follows the Gaussian distribution, separation distance d can be obtained by the relationship [37](5)ArAa=12π∫d/σ∞exp−x2/2dx=12erfcd2σ,where Ar, Aaare the total real contact area and the nominal contact area, respectively. erfc(x) represents the complementary error function.Figure 2
The contact between a fractal asperity at leveln and a rigid smooth plane, where ln, lnr, and lnt are the base diameter, the real contact diameter, and the truncation diameter, respectively [43].According to the Hertz theory [37], the contact area and contact load of the asperity within the elastic deformation stage can be obtained as(6)an=πRnωn,(7)fne=43ERn1/2ωn3/2,where Rn is the radius of curvature at the asperity summit, Rn=1/d2zn/dx2x=0=lnD/π2GD−1. E is the equivalent elastic modulus, 1/E=1−ν12/E1+1−ν22/E2, and E1, E2, ν1, and ν2 are elastic modulus and Poisson’s ratios of the two rough surfaces. According to equations (6) and (7), the relationship between contact load and real contact area of a single asperity within the elastic deformation stage is(8)fne=4Eπ1/2GD−13lnDan3/2.


### 2.2.2. The Existing Condition for Initial Yield

The internal stress of the asperity increases with the increase of the contact load or contact interference. An initial yield point will eventually be generated inside the asperity due to excessive stress. The interference corresponding to the initial yield point is termed as the critical interferenceωnec, and is given as [47](9)ωnec=πKH2E2Rn.where K is the hardness factor of the softer material and is given by K=0.454+0.41ν. H denotes the hardness of the softer material and is given by H=2.8Y, and Y is the yield strength of the softer material. The elastic critical contact area of the asperity anec can be obtained as [43](10)anec=πRnωnec=1πKH2E2lnDGD−12.According to equations (7), (9), and (10), we can get the expression of the elastic critical contact load of the asperity fnec, which is(11)fnec=23KHanec.


### 2.2.3. The Revised Elastoplastic Microcontacts

As the load or interference of the asperity increases, the plastic part inside the asperity will gradually expand to the contact surface, and an annular plastic part will be formed on the contact surface, while the rest of the contact surface surrounded by the plastic region remains elastic deformation. This transition stage, with the transition interference ratio of1<ωn/ωnec≤6(anec<an≤anepc), is known as the first elastoplastic deformation stage. The first elastoplastic critical contact area anepc and the first elastoplastic contact load fnep1 of the asperity are derived by Yuan et al. [43] based on the ME model, which are(12)anepc=7.1197anec,(13)fnep1=23KH×1.1282anec−0.2544an1.2544.In the first elastoplastic deformation stage, the hardness of the material will change with the deformation rather than remaining a constant value [25]. According to equations (12)∼(13), we assume that the material hardness within the first elastoplastic deformation stage HG1an satisfies the following relationship:(14)HG1an=c11Yananecc12,anec<an≤anepc,where c11, c12are the parameters need to be solved.Equation (14) should satisfy two boundary conditions:(15)HG1anec=pemanec,(16)HG1anepc=pepm1anepc,where peman is the average contact pressure of the asperity in the elastic deformation stage, which is given by peman=fne/an. Hence, pemanec can be obtained as pemanec=2/3KH asperity in the first elastoplastic deformation stage and is given by pepm1an=fnep1/an. Substituting equation (14) and equation pemanec=KH into equation (15) yields c11Y=2/3KH. Thus, the parameter c11 can be obtained as(17)c11=2815K.Substituting equations (12)–(14) and (17) into equation (16) yields(18)23KH×7.1197c12=23KH×1.1282anec−0.25447.1197anec0.2544.And parameterc12 can be obtained as(19)c12=ln1.1282×7.11970.2544ln7.1197.Therefore, the contact load of a single asperity in the first elastoplastic deformation stage is revised as(20)fnep1′=HG1anep1·anep1=2815KYanec−c12anc12+1.If the contact load or interference further increases, the plastic part of the asperity gradually expands to envelop the shrinking elastic core. According to Etsion and Etsion [48], the transition interference ratio of this stage is 6<ωn/ωnec≤110(anepc<an≤anpc). Some researchers named this stage the second elastoplastic deformation stage. The second elastoplastic critical contact area anpc and the second elastoplastic contact load fnep2 of the asperity are also given by Yuan et al. [43]:(21)anpc=205.3827anec,(22)fnep2=23KH×1.4988anec−0.1021an1.1021.The hardness of the material will also change with the deformation in the second elastoplastic deformation stage [25]. Similarly, we assume that the material hardness within the second elastoplastic deformation stage HG2an satisfies the following relationship:(23)HG2an=c21Yananecc22,anepc<an≤anpc,where c21 and c22 are the two parameters need to be solved.And equation (23) should also satisfy two boundary conditions:(24)HG2anepc=pem1anepc,(25)HG2anpc=pepm2anpc,where pepm2an is the average contact pressure of the asperity in the second elastoplastic deformation stage and is given by pepm2an=fnep2/an.Substituting equations (12), (13), and (23) into equation (24) yields(26)c21Y7.1197c22=23KH×1.1282×7.11970.2544.Substituting equations (21), (22), and (23) into equation (25) yields(27)c21Y205.3827c22=23KH×1.4988×205.38270.1021.According to equations (26) and (27), parameter c22 can be expressed as(28)c22=ln1.4988+0.1021ln205.3827−ln1.1282−0.2544ln7.1197ln205.3827−ln7.1197.And parameterc21 can be obtained as(29)c21=K×2815×1.1282×7.11970.2544−c22.Therefore, the contact load of a single asperity in the second elastoplastic deformation stage is revised as(30)fnep2′=HG2anep2·anep2=K×2815×1.1282×7.11970.2544−c22Yanec−c22an1+c22.


### 2.2.4. The Existing Plastic Microcontacts

Whenωn/ωnec>110(an>anpc), the asperity undergoes full plastic deformation. The contact area and contact load of the asperity in this deformation stage are [43](31)an=2πRnωn,fnp=Han.For the sake of clarity, Figure3 shows the deformation law of a single asperity with the fractal asperity frequency level n, as well as the relationship between each critical contact area.Figure 3
The diagram of the deformation law of a single asperity and the relationship between each critical contact area.


## 2.2.1. The Existing Elastic Microcontacts

According to the WM function, the profile of the fractal asperity at a certain leveln before deformation is [37](2)znx=GD−1ln2−Dcosπxln,−ln2<x<ln2,where ln represents the length scale of the asperity of level n; that is, the base diameter of the asperity of level n, ln=1/γn.The asperity height is [37](3)δn=zn0=GD−1ln2−D.The contact between two rough surfaces can be simplified to an equivalent rough surface in contact with a rigid flat surface [42]. Figure 2 shows the relationship between the geometric parameters when a fractal asperity is in contact with a rigid smooth surface. According to Figure 2, the interference ωn of the asperity with a length scale ln is(4)ωn=δn−d=GD−1ln2−D−d,where d is the separation distance between the rigid flat surface and the equivalent rough surface. If the surface height zx,y follows the Gaussian distribution, separation distance d can be obtained by the relationship [37](5)ArAa=12π∫d/σ∞exp−x2/2dx=12erfcd2σ,where Ar, Aaare the total real contact area and the nominal contact area, respectively. erfc(x) represents the complementary error function.Figure 2
The contact between a fractal asperity at leveln and a rigid smooth plane, where ln, lnr, and lnt are the base diameter, the real contact diameter, and the truncation diameter, respectively [43].According to the Hertz theory [37], the contact area and contact load of the asperity within the elastic deformation stage can be obtained as(6)an=πRnωn,(7)fne=43ERn1/2ωn3/2,where Rn is the radius of curvature at the asperity summit, Rn=1/d2zn/dx2x=0=lnD/π2GD−1. E is the equivalent elastic modulus, 1/E=1−ν12/E1+1−ν22/E2, and E1, E2, ν1, and ν2 are elastic modulus and Poisson’s ratios of the two rough surfaces. According to equations (6) and (7), the relationship between contact load and real contact area of a single asperity within the elastic deformation stage is(8)fne=4Eπ1/2GD−13lnDan3/2.


## 2.2.2. The Existing Condition for Initial Yield

The internal stress of the asperity increases with the increase of the contact load or contact interference. An initial yield point will eventually be generated inside the asperity due to excessive stress. The interference corresponding to the initial yield point is termed as the critical interferenceωnec, and is given as [47](9)ωnec=πKH2E2Rn.where K is the hardness factor of the softer material and is given by K=0.454+0.41ν. H denotes the hardness of the softer material and is given by H=2.8Y, and Y is the yield strength of the softer material. The elastic critical contact area of the asperity anec can be obtained as [43](10)anec=πRnωnec=1πKH2E2lnDGD−12.According to equations (7), (9), and (10), we can get the expression of the elastic critical contact load of the asperity fnec, which is(11)fnec=23KHanec.


## 2.2.3. The Revised Elastoplastic Microcontacts

As the load or interference of the asperity increases, the plastic part inside the asperity will gradually expand to the contact surface, and an annular plastic part will be formed on the contact surface, while the rest of the contact surface surrounded by the plastic region remains elastic deformation. This transition stage, with the transition interference ratio of1<ωn/ωnec≤6(anec<an≤anepc), is known as the first elastoplastic deformation stage. The first elastoplastic critical contact area anepc and the first elastoplastic contact load fnep1 of the asperity are derived by Yuan et al. [43] based on the ME model, which are(12)anepc=7.1197anec,(13)fnep1=23KH×1.1282anec−0.2544an1.2544.In the first elastoplastic deformation stage, the hardness of the material will change with the deformation rather than remaining a constant value [25]. According to equations (12)∼(13), we assume that the material hardness within the first elastoplastic deformation stage HG1an satisfies the following relationship:(14)HG1an=c11Yananecc12,anec<an≤anepc,where c11, c12are the parameters need to be solved.Equation (14) should satisfy two boundary conditions:(15)HG1anec=pemanec,(16)HG1anepc=pepm1anepc,where peman is the average contact pressure of the asperity in the elastic deformation stage, which is given by peman=fne/an. Hence, pemanec can be obtained as pemanec=2/3KH asperity in the first elastoplastic deformation stage and is given by pepm1an=fnep1/an. Substituting equation (14) and equation pemanec=KH into equation (15) yields c11Y=2/3KH. Thus, the parameter c11 can be obtained as(17)c11=2815K.Substituting equations (12)–(14) and (17) into equation (16) yields(18)23KH×7.1197c12=23KH×1.1282anec−0.25447.1197anec0.2544.And parameterc12 can be obtained as(19)c12=ln1.1282×7.11970.2544ln7.1197.Therefore, the contact load of a single asperity in the first elastoplastic deformation stage is revised as(20)fnep1′=HG1anep1·anep1=2815KYanec−c12anc12+1.If the contact load or interference further increases, the plastic part of the asperity gradually expands to envelop the shrinking elastic core. According to Etsion and Etsion [48], the transition interference ratio of this stage is 6<ωn/ωnec≤110(anepc<an≤anpc). Some researchers named this stage the second elastoplastic deformation stage. The second elastoplastic critical contact area anpc and the second elastoplastic contact load fnep2 of the asperity are also given by Yuan et al. [43]:(21)anpc=205.3827anec,(22)fnep2=23KH×1.4988anec−0.1021an1.1021.The hardness of the material will also change with the deformation in the second elastoplastic deformation stage [25]. Similarly, we assume that the material hardness within the second elastoplastic deformation stage HG2an satisfies the following relationship:(23)HG2an=c21Yananecc22,anepc<an≤anpc,where c21 and c22 are the two parameters need to be solved.And equation (23) should also satisfy two boundary conditions:(24)HG2anepc=pem1anepc,(25)HG2anpc=pepm2anpc,where pepm2an is the average contact pressure of the asperity in the second elastoplastic deformation stage and is given by pepm2an=fnep2/an.Substituting equations (12), (13), and (23) into equation (24) yields(26)c21Y7.1197c22=23KH×1.1282×7.11970.2544.Substituting equations (21), (22), and (23) into equation (25) yields(27)c21Y205.3827c22=23KH×1.4988×205.38270.1021.According to equations (26) and (27), parameter c22 can be expressed as(28)c22=ln1.4988+0.1021ln205.3827−ln1.1282−0.2544ln7.1197ln205.3827−ln7.1197.And parameterc21 can be obtained as(29)c21=K×2815×1.1282×7.11970.2544−c22.Therefore, the contact load of a single asperity in the second elastoplastic deformation stage is revised as(30)fnep2′=HG2anep2·anep2=K×2815×1.1282×7.11970.2544−c22Yanec−c22an1+c22.


## 2.2.4. The Existing Plastic Microcontacts

Whenωn/ωnec>110(an>anpc), the asperity undergoes full plastic deformation. The contact area and contact load of the asperity in this deformation stage are [43](31)an=2πRnωn,fnp=Han.For the sake of clarity, Figure3 shows the deformation law of a single asperity with the fractal asperity frequency level n, as well as the relationship between each critical contact area.Figure 3
The diagram of the deformation law of a single asperity and the relationship between each critical contact area.


## 2.3. The Revised Size Distribution Function of Contact Spots

In the MB model, the size distribution function of contact spots is defined as [37](32)na=D2alD2a−D+22,0<a<al,1<D<2,where al denotes the largest contact area of asperity. And the total real contact area is [37](33)Ar=∫0alnaada=D2−Dal.According to equation (2), the period of asperity at frequency level n is Tn=2π/2πγn, and the period of asperity at frequency level n+1 is Tn+1=2π/2πγn+1. Therefore, Tn+1/Tn=1/γ. Hence, we can obtain the relationship between the size distribution function of the asperity at frequency level n and the size distribution function of the asperity at frequency level n+1, which is nn+1a=γnna. Letting the size distribution function of the asperity at the initial frequency level nmin be nnmina=Mna, then the contact spot size distribution function at each frequency level can be obtained as(34)nna=Mγn−nminna,nmin≤n≤nmax,where M is a parameter that needs to be solved. According to equations (33) and (34), the total real contact area of asperities at each frequency level is(35)Anr=∫0anlnnaada=MD2−Dγn−nminanl,where anl is the largest contact area of contact spots at frequency level n. According to the geometry relationship shown in Figure 2, we assume that the largest contact spot of frequency level nmin is the largest contact spot of the entire rough surface, i.e., anminl=al. And we assume that the relationship between the largest contact spots of two adjacent frequency levels is an+1l/anl=1/γ2. Then, the largest contact area of the asperity at any frequency level would be anl=1/γ2n−nminal. Thus equation (35) can be reexpressed as(36)Anr=Mγn−nminD2−Dal=Mγn−nminAr.And the relationship betweenanl and Ar can be obtained as(37)anl=2−DD1γn−nmin2Ar.The parameterM can be determined according to the relationship Anminr+Anmin+1r+⋯+Anmaxr=Ar. Substituting equation (36) into it yields M=1−1/γ/1−1/γnmax−nmin+1.


## 2.4. The Revised Critical Asperity Frequency Levels

For a certain frequency leveln, the relationship between the interference ωn and the asperity height δn is ωn≤δn [46]. If the height of an asperity is less than its elastic critical interference, inelastic deformation occurs. Thus, the condition of elastic deformation for an asperity given by the YC model is δn≤ωnec [43]. However, according to Etsion and Etsion [41], the condition for elastic deformation of an asperity should be ωn≤ωnec, which indicates that the asperity occurs elastic deformation when the interference is not greater than its critical interference. Etsion and Etison [41] further derived the inequality used to solve the elastic critical length scale lec(or the elastic critical asperity frequency level nec), which is 1−cosπlnt/2ln≤KH/2E⋅ln/GD−12. Unfortunately, we cannot get the analytical solution of lec according to this inequality, nor can we get an analytical solution of nec according to the relationship ln=1/γn.Now, we reanalyze the conditions for judging the elastic deformation of the asperity. At a mean surface separation distanced, if the largest interference of the asperity at frequency level n is not greater than its critical interference, all of the contact asperities at frequency level n deform elastically. Thus, the condition we give to judge the elastic deformation of the asperity is(38)ωnl≤ωnec,where ωnl is the largest interference of the contact spot at frequency level n. According to equation (6), equation (39) can be re-written as(39)anl≤anec.Substituting equations (10) and (37) into the above inequality, the elastic critical frequency level can be obtained from the following equation:(40)2−DD1γnec−nmin2Ar=1πKH2E21γ2necD1G2D−1.According to the above inequality, the elastic critical frequency levelnec is not only related to the material parameters H,E and the topography parameters D,G but also related to the total real contact area of the two rough surfaces Ar. Similarly, the first elastoplastic critical frequency level nepc and the second elastoplastic critical frequency level npc can be obtained from the following two equations:(41)2−DD1γnepc−nmin2Ar=7.11971πKH2E21γ2nepcD1G2D−1,2−DD1γnpc−nmin2Ar=205.38271πKH2E21γ2npcD1G2D−1.According to equations (41)∼(43), the equations used to calculate the critical frequency levels nec,nepc,npc are not only related to the material parameters H,E and topography parameters D,G, but also related to the total real contact area of the rough surface Ar. As the given parameter Ar changes, so do the calculation results of the critical frequency levels nec,nepc,npc. That is, the critical frequency levels nec,nepc,npc variate with the change of the Ar value. However, in the equations for calculating the critical frequency levels given by the YC model, the critical frequency levels are only related to the material parameters and topography parameters, and the influence of the total real contact area Ar on the critical frequency levels is not considered.


## 2.5. The Revised Real Contact Area and Contact Load of the Rough Surface

When the asperity frequency level ranges fromnmin to nec, i.e., nmin≤n≤nec, only elastic deformation takes place. And the real contact area Ar1 and contact load Fr1 are obtained as(42)Ar1=Are=∑n=nminnec∫0anlMγn−nminnaada,Fr1=Fre=∑n=nminnec∫0anlfneMγn−nminnada.When the asperity frequency level belongs tonec<n≤nepc, elastic or the first elastoplastic deformation can take place. The real contact area Ar2 is obtained as(43)Ar2=Are+Arep1=∑n=nec+1nepc∫0anecMγn−nminnaada+∑n=nec+1nepc∫anecanlMγn−nminnaada.And the contact loadFr2 is obtained as(44)Fr2=Fre+Frep1=∑n=nec+1nepc∫0anecfneMγn−nminnada+∑n=nec+1nepc∫anecanlfnep1′Mγn−nminnada.When the asperity frequency level isnepc<n≤npc, elastic deformation, the first elastoplastic deformation, or the second elastoplastic deformation will occur. The real contact area Ar3 is evaluated as(45)Ar3=Are+Arep1+Arep2=∑n=nepc+1npc∫0anecMγn−nminnaada+∑n=nepc+1npc∫anecanepcMγn−nminnaada+∑n=nepc+1npc∫anepcanlMγn−nminnaada,.And the contact loadFr3 is evaluated as(46)Fr3=Fre+Frep1+Frep2=∑n=nepc+1npc∫0anecfneMγn−nminnada+∑n=nepc+1npc∫anecanepcfnep1′Mγn−nminnada+∑n=nepc+1npc∫anepcanlfnep2′Mγn−nminnada,.When the asperity frequency level belongs tonpc<n≤nmax, all of the four types of deformation can take place. The real contact area Ar4 is calculated as(47)Ar4=Are+Arep1+Arep2+Arp=∑n=npc+1nmax∫0anecMγn−nminnaada+∑n=npc+1nmax∫anecanepcMγn−nminnaada,+∑n=npc+1nmax∫anepcanpcMγn−nminnaada+∑n=npc+1nmax∫anpcanlMγn−nminnaada.And the contact loadFr4 is estimated as(48)Fr4=Fre+Frep1+Frep2+Frp=∑n=npc+1nmax∫0anecfneMγn−nminnada+∑n=npc+1nmax∫anecanepcfnep1′Mγn−nminnada+∑n=npc+1nmax∫anepcanpcfnep2′Mγn−nminnada+∑n=npc+1nmax∫anpcanlfnpMγn−nminnada.The detailed calculation results of the real contact area and contact load for the rough surface are provided in Appendix A. The total real contact area and contact load of all asperity levels are estimated by(49)Ar=Ar1+Ar2+Ar3+Ar4.Fr=Fr1+Fr2+Fr3+Fr4.The nondimensional forms ofAr and Fr are defined as [43](50)Ar∗=ArAa,Fr∗=FrAaE,where the nominal contact area Aa is given by Aa=L2, L=1/γnmin. For the sake of clarity, the flowchart of the solution procedure for extracting the relationship between the total real contact area and contact load is shown in Figure 4.Figure 4
Flowchart for extracting the relationship between the total real contact area and contact load.


## 3. Results and Discussion

In order to obtain the contact parameters, the equivalent rough surface parameters [40, 46, 49] used in the present work are shown in Table 1. Let the topography parameters be G=2.5×10−9 m, D=1.5. Let the sampling length and the resolution be L=0.5×10−4 m (L=1/γnmin), Ls=1.5×10−9m (Ls=1/γnmax), respectively. The nominal contact area, the range of asperity frequency levels, and the critical asperity frequency levels, hence, can be calculated as Aa=L2=0.25×10−8m2, nmin=24, nmax=50, nec=32, nepc=37, npc=45. It should be noted that parameters D and G refer to the fractal parameters of the equivalent rough surface, which can be calculated from the topography parameters of the two contact surfaces. The calculation procedure is given in Appendix B.Table 1
Equivalent rough surface parameters.
ParametersValuesEquivalent elastic modulusE7.2×1010 N/㎡Poisson’s ratioν0.17Initial hardnessH5.5×109 N/㎡Profile scale parameterG1×10−9∼1×10−6 mFractal dimensionD1<D<2Figure5 shows the comparison of the proposed model and the MB model on the critical contact areas and the largest contact area at different asperity levels. The value of nondimensional total real contact area Ar∗(Ar∗=Ar/Aa) is 0.4. It is can be seen from Figure 5 that the value of critical contact areas (anec, anepc, anpc) and the largest contact area anl of the present work decrease with the increase of frequency level n. And when n is taken as a specific value, the relationship between anec, anepc and anpc of the present work is anec<anepc<annc. Moreover, in the present work, the largest contact area is smaller than the elastic critical contact area when the value of n belongs to 0,32. The largest contact area is larger than the elastic critical contact area and smaller than the first elastoplastic critical contact area if the asperity level n belongs to 33,37. The largest contact area is larger than the first elastoplastic critical contact area and smaller than the second elastoplastic critical contact area if the asperity level n belongs to 38,45. The largest contact area is larger than the second elastoplastic critical contact area if the asperity level n belongs to 46,50. However, the largest contact area anl and the elastic critical contact area anec of the MB model do not change with the increase of n. And the largest contact area of the MB model is always larger than its elastic critical contact area. In addition, Figure 5 also shows the relationship between the asperity frequency level and the critical contact areas in the YC model. And the relationships between n and critical contact areas in the YC model are the same as those in the proposed model. The analytic expression of anl is not given by the YC model, so we cannot give a comparison of anl between the proposed model and the YC model in Figure 5.Figure 5
Comparison of the proposed model and the MB model on the critical contact areas and the largest contact area at different asperity frequency levels (Ar/Aa=0.4).The relationships between the fractal dimensionD and the critical asperity frequency levels of the proposed model are shown in Figure 6. Figure 6 also shows the changing trend of critical asperity frequency levels with fractal dimension D in the YC model. The characteristic length scale is G=2.5×10−9m, and the nondimensional total real contact area Ar∗ is 0.4. It should be noted that, to facilitate comparison, the critical frequency level is set to be 0 if it is not within the range of nmin and nmax. Figure 6 indicates that the elastic critical level, elastoplastic critical level, and plastic critical level all increase with the increase of fractal dimension in both the proposed model and the YC model. Moreover, the relationship between nec, nepc and npc in the proposed model is the same as that in the YC model, which is nec<nepc<npc. In addition, when D is taken as a specific value, the critical levels calculated by the proposed model are larger than the critical levels of the same type calculated by the YC model. The reason for this difference is that the method of solving the elastic critical level in the proposed model is different from that in the YC model. The conditions for the proposed model and the YC model to solve the elastic critical level are ωnl≤ωnec and δn≤ωnec, respectively. And the asperity interference will not be greater than its height, i.e., ωnl≤δn, so the value of nec calculated by the YC model is smaller than that calculated by the proposed model. It is also due to the same reason that the value of nepc and npc calculated by the YC model are smaller than those calculated by the proposed model.Figure 6
Comparison of the YC model and the proposed model on the critical frequency levels at different fractal dimensionsD(G=2.5×10−9 m, Ar/Aa=0.4).The relationships between nondimensional total real contact areaAr∗ and the critical asperity frequency levels of the proposed model and the YC model are presented in Figure 7. Figure 7 indicates that the critical frequency levels of the proposed model decrease with the increase of Ar∗ value. However, as the Ar∗ value increases, the critical frequency levels of the YC model remain unchanged. According to equations (41)∼(43), the asperity critical frequency levels of the proposed model are related to the total real contact area of the rough surface, and the values of n and the total real contact area Ar are negatively correlated. Therefore, the critical frequency levels in the proposed model decrease with the increase of Ar∗. The influence of the total real contact area on the asperity critical frequency levels is not considered in the YC model, and the change of Ar∗ will consequently not affect the calculation results of the critical frequency levels.Figure 7
Comparison of the YC model and the proposed model on the critical asperity frequency levels with different nondimensional total real contact areas (G=2.5×10−9 m, D=1.5).Figure8 shows the “contact load–contact area” relationships of a single asperity both with and without considering the variation of material hardness within the first elastoplastic deformation stage. The frequency level n is 33, and the value of Ar∗ is 0.4. According to Figure 8, whether the variation of the material hardness is considered or not, the “contact load–contact area” relationship of a single asperity shows the same trend; that is, the contact area increases with the increase of contact load. In addition, the contact area with the consideration of material hardness variation is larger than that without the consideration of material hardness variation with the same contact load. And the difference between them gradually shrinks as the contact load increases.Figure 8
The contact load-contact area relationship of a single asperity in the first elastoplastic deformation stage (n=33, Ar/Aa=0.4).According to the previous analysis, the material hardness within the first elastoplastic deformation stage can be expressed as a function of the contact area of a single asperity. Figure9(a) presents the relationships between the material hardness and the contact area of a single asperity for different values of n. It is can be seen from Figure 9(a) that the material hardness in the first elastoplastic deformation stage increases with the increase of contact area. And when the contact area is taken as a fixed value, the larger the value of n, the larger the material hardness. Figure 9(b) shows that when the contact area is taken as a fixed value, the smaller the value of D, the larger the material hardness. Figure 9(c) demonstrates that when the contact area is taken as a fixed value, the larger the value of G, the larger the material hardness.Figure 9
Relationships between the material hardness and the contact area of a single asperity within the first elastoplastic deformation stage. .td.jin Bool Technology Co. (a)D=1.5,Ar/Aa=0.4. (b)n=33,G=2.5×10−9,Ar/Aa=0.4. (c)n=33,D=1.5,Ar/Aa=0.4.
(a)(b)(c)The “contact load–contact area” relationships of a single asperity with and without the consideration of material hardness variation within the second elastoplastic deformation stage are shown in Figure10. The asperity frequency level is 38, and the value of Ar∗ is 0.4. In the second elastoplastic deformation stage, regardless of whether the variation of material hardness is considered or not, the contact area increases with the increase of contact load. And the real contact area without considering the variation of material hardness is larger than the contact area considering the variation of material hardness for the same contact load. Moreover, as the contact load increases, the difference between the contact areas tends to increase.Figure 10
The contact load-contact area relationship of a single asperity in the second elastoplastic deformation stage (n=38, Ar/Aa=0.4).Figure11 shows the relationship between the material hardness and the contact area of a single asperity for different values of asperity level n, fractal dimension D, and roughness parameter G in the second elastoplastic deformation stage. Similar to the first elastoplastic deformation stage, the material hardness increases as the contact area of a single asperity increases. Moreover, the influences of different values of asperity level, fractal dimension, and roughness parameter on the relationship between the contact area and material hardness are the same as those in the first elastoplastic deformation stage.Figure 11
Relationships between the material hardness and contact area of a single asperity within the second elastoplastic deformation stage. (a)D=1.5,G=2.5×Ar/Aa=0.4. (b)n=38,G=2.5×10−9,Ar/Aa=0.4. (c)n=38,D=1.5,Ar/Aa=0.4.
(a)(b)(c)Figure12(a) indicates the comparison of the proposed model and the MB model on the relationship between nondimensional total contact load Fr∗ and nondimensional total real contact area Ar∗ for different roughness parameters. The comparison of the proposed model and the MB model on the relationship between nondimensional total real contact area Ar∗ and the fraction of elastic contact area for different roughness parameters is shown in Figure 12(b). In order to obtain reasonable comparison results, some surface topography parameters and material parameters are given, such as D=1.5, Y/E=0.01. It is can be seen from Figure 12(a) that, in the MB model, as G/Aa increases, the Fr∗ value required to produce a specific Ar∗ value increases. The same changing trend is also found in the presented model. In addition, the Fr∗ value of the proposed model required to produce a particular Ar∗ value is less than that of the MB model. Figure 12(b) indicates that, in the MB model, the fraction of elastic contact area increases as Ar∗ increases. This is the contradiction between the MB model and classical contact mechanics. Classical contact mechanics believes that as the real contact area increases, the contact load increases, which means that more asperities will undergo inelastic deformation. As a result, the fraction of elastic contact area should decrease. In addition, it is can be obtained from Figure 12(b) that the fraction of elastic contact area in the MB model decreases with the increase of G/Aa. In the proposed model, the fraction of elastic contact area is approximately equal to 1, and it hardly changes with the change of Ar∗ value. It does not change with the change of G/Aa value either.Figure 12
Comparison of the proposed model and the MB model (a) on the relationship between nondimensional total contact load and nondimensional total real contact area and (b) on the relationship between nondimensional total real contact area and the fraction of elastic contact area for different roughness parameters (nmax=24, nmax=50, D=1.5, Y/E=0.01).
(a)(b)Figure13(a) shows the comparison of the present model and the MB model on the relationship between nondimensional total contact load Fr∗ and nondimensional total real contact area Ar∗ for different values of D. Figure 13(a) indicates that, in the MB model, the Fr∗ value required for a particular Ar∗ value decreases with the increase of D when D belongs to 1.1,1.5. When D is between 1.5 and 1.9, as the value of D increases, the Fr∗ value required to produce a specific Ar∗ value increases. However, in the proposed model, the Fr∗ value needed for a specific Ar∗ value decreases as the value of D increases within the global value range of D. Furthermore, it is also can be seen from Figure 13(a) that the Fr∗ value corresponding to the same Ar∗ value in the proposed model is less than that in the MB model when the fractal dimension takes the same value. Figure 13(b) indicates that, in the MB model, the fraction of elastic contact area is 0 when D=1.1, which means that no asperity deforms elastically. And when D=1.5 and D=1.9, the fraction of elastic contact area increases with the increase of Ar∗ value. As mentioned above, this phenomenon is contrary to classical contact mechanics, while in the presented model, when D=1.1, the fraction of elastic contact area decreases as Ar∗ value increases. This trend of change is consistent with classical contact mechanics. When D takes the value of 1.5 and 1.9, the fraction of elastic contact area in the presented model hardly changes with the increase of Ar∗ value and is approximately equal to 1. In addition, the fraction of elastic contact area in the proposed model is greater than that in the MB model when the fractal dimension takes the same value. It should be noted that the YC model does not give the calculation formula of the largest contact area of the asperity anl at frequency level n, which makes it impossible for us to program the YC model. Thus, we do not give a comparison between the proposed model and the YC model regarding the Ar/Aa∼Fr/AaE relationship and the Are/Ar∼Ar/Aa relationship. Nevertheless, the idea of contact modeling on rough surfaces in reference [43] still has a vital reference significance.Figure 13
Comparison of the proposed model and the MB model (a) on the relationship between nondimensional total contact load and nondimensional total real contact area and (b) on the relationship between nondimensional total real contact area and the fraction of elastic contact area for different fractal dimensions (nmax=24, nmax=50, G/Aa=1.5, Y/E=0.01).
(a)(b)Figure14 presents the comparison results of the proposed model and other contact models with experimental data on the relationship between the nondimensional total contact load and the nondimensional total real contact area. In order to ensure the reliability of the comparison results, some parameters used in the present work are the same as those adopted by the MB model and the YC model, i.e., D=1.49, G/Aa=1×10−10, Y/E=0.05. It should be noted that the critical contact areas of the proposed model are scale-dependent. Therefore, the values of the critical contact areas are related to the sampling length and the resolution of the measuring instrument. Different sampling lengths and instrument resolutions will result in different critical contact areas, which will lead to the change in the relationship between nondimensional total contact load Fr∗ and nondimensional total real contact area Ar∗. Figure 14(a) shows the comparison results for nmin=67, nmax=90. At lower loads, the MB model is in good agreement with the experimental data. Under the same load within the lower loads, the results of the proposed model and the YC model are larger than the experimental data, and the proposed model is closer to experimental data than the YC model. At higher loads, the MB model deviates seriously from the experimental results. The explanation given by Majumdar and Bushan is that the interaction of asperities is not considered. However, both the proposed model and the YC model are close to the experimental values at higher loads and are better than the results of the GW model. Figure 14(b) indicates that the contact areas of the proposed model and the YC model are larger than the experimental value for the same load at lower loads, and the YC model is closer to the experimental results. At higher loads, the results of the proposed model and the YC model are less than the experimental values, and the proposed model is in better agreement with the experimental results. In addition, it can be seen from Figures 14(a) and 14(b) that, in the entire loading range, the agreement between the proposed model and the experimental data is better than that of other models.Figure 14
Comparison of the proposed model and other contact models with experiment data. (a)nmin=67,nmax=90. (b)nmin=70,nmax=90.
(a)(b)


## 4. Conclusions

An improved fractal contact model considering the variation of the critical asperity levels as well as the variation of the material hardness is proposed in the present work. The main conclusions are as follows:(1)
The real contact area of a single asperity obtained by considering the variation of material hardness is greater than that without considering the variation of material hardness within the first elastoplastic deformation stage, while in the second elastoplastic deformation stage, the real contact area of a single asperity considering the variation of material hardness is less than that without considering the variation of material hardness.(2)
The size distribution functions of the contact spots at different frequency levels are derived. The expressions of asperity critical frequency levels are rederived. The results show that the critical asperity levels are not constant values, but variable values related to the total real contact area of the rough surface and decrease with the increase of the total real contact area.(3)
The proposed model is a modification and improvement of the existing fractal contact models, which can lead to a more accurate relationship between the contact load and the total real contact area of the rough surface.The proposed model is helpful for the analysis of the sealing performance of static metal seals under different contact pressures in our subsequent studies.



---

*Source: 2985674-2022-09-21.xml*