Contact Model of Revolute Joint with Clearance Based on Fractal Theory

The contact stiffness of a mechanical bonding surface is an important parameter in determining the normal and radial contact force. To improve the calculation accuracy of the contact force model, the surface roughness of the bonding surface and the energy loss that necessarily occurs during the impact process should be considered comprehensively. To study the normal contact force of a revolute joint with clearance more accurately in the case of dry friction, a nonlinear stiffness coefficient model considering the surface roughness was established based on fractal theory, which considers the elastic, elastic-plastic, and plastic deformations of the asperities of the contact surface during the contact process. On this basis, a modified nonlinear spring damping model was established based on the Lankarani–Nikravesh contact force model. The laws influencing the surface roughness, recovery coefficient, initial velocity, and clearance size on the impact force were revealed, and were compared with the Lankarani–Nikravesh model and a hybrid model using MATLAB. The maximum impact force was obtained using a modified contact force model under different initial velocities, different clearances, and different degrees of surface roughness, and the calculated results were then compared with the experiment results. This study indicates that the modified model can be used more widely than other models, and is suitable for both large and small clearances. In particular, the modified model is more accurate when calculating the contact force of a revolute joint with a small clearance.


Introduction
In recent years, many scholars domestically and abroad have conducted a number of studies on the contact impact problem of dynamics with clearance. A common contact theory is the Hertz contact theory [1]; however, this theory can only describe the elastic deformation of the material during the contact process, and does not consider the energy dissipation. A contact force model with nonlinear damping (Hunt-Crossley model) was thus proposed by Hunt et al. [2], which considers the dissipation of energy during the contact process. Yan et al. [3] introduced the application of a clearance model in the analysis of nonlinear characteristics, a performance evaluation, and a reliability evaluation of the mechanical system dynamics. Han et al. [4] used the joint contact impact force as a parameter to achieve a high-precision lifetime prediction of a mechanism. Lankarani et al. [5] determined the relationship between the damping factor and restitution coefficient by analyzing the energy loss in a system before and after the contact process, and then put forward the Lankarani-Nikravesh continuous contact force model (L-N model). The L-N model not only considers the dissipation of energy during the contact process, but also reflects the influence of the material properties, geometric characteristics, and motion states of the contact body during the contact process. Based on the L-N model and improved elastic foundation contact model [6], a hybrid model was established by Bai et al. [7], through which the nonlinear stiffness coefficient was presented. Wang et al. [8] established an improved nonlinear contact impact force model by considering the axial dimensions, material coefficients, and energy loss comprehensively. However, the previous modeling methods were all established at a macro level, and do not reflect the microscopic characteristics Page 2 of 11 Li et al. Chin. J. Mech. Eng. (2018) 31:109 of the contact surface, such as the influence of the micro topography and surface roughness on the contact stiffness.
There are numerous micro-convex bodies on an actual motion pair surface, all of which have a particular degree of roughness. A rough surface morphology has an important influence on the contact stiffness coefficient [9,10], and affects the contact impact force during an actual contact collision. Therefore, drawing on research methods regarding the contact problems of the surface morphology, a contact model of a revolute joint with clearance when considering the surface morphology can be established, which is significant to a precision analysis and dynamic design of the mechanical system dynamics.
In fact, many asperities appear on a rough contact surface, and a rough surface morphology has an important influence on [11]. Scholars have conducted a number of studies on theoretical contact force models over the years, including the classic GW model [12], WA model [13], and MB model [14]. Furthermore, taking account of the influences of the surface asperities grade, the deformation characteristic, and the contact area, Yuan et al. [15,16] established and analyzed the fractal contact mode between rough surfaces based on the fractal theory. Pan et al. [17,18] derived a theoretical three-dimensional mechanical joint surface based on the W-M function of the three-dimensional joint surface topography. Liu et al. [19] established a fractal model of the tangential contact stiffness of two spherical surfaces when considering the friction factors. Wang et al. [20] established a deterministic model for the elastoplastic contact of rough surfaces under combined loads, and obtained the actual contact area, pressure distribution, and plastic deformation of the micro-convex peaks. A series of studies were carried out on the model of the contact stiffness coefficient of the contact surface. However, there have been relatively few studies establishing the organic connections between the macroscopic levels and microscopic states when considering a rough surface.
In this paper, the contact stiffness coefficient was derived when considering the surface roughness of the contact surface, and a modified contact force model was obtained. Based on the fractal theory, a modified normal contact stiffness coefficient model was also established, which considers the contact surface roughness, as well as the elastic, elastic-plastic, and plastic deformations of the asperities during the contact process. Meanwhile, a modified contact force model was established by combining the L-N model expression. Finally, an experimental system was developed, and the correctness of the model was verified through an experimental study.

Establishment of Modified Contact Force Model
When two surfaces consisting of asperities are close to each other, the contact only occurs at the top of the asperity owing to the asperities on the contact surface of the actual motion joint. During a study on the surface contact, the contact of two surfaces is usually simplified using the contact between a rough surface and a rigid plane [21].

Contact Deformation of Asperity
Based on fractal theory, the distribution function of asperity on contact surface n a ′ [22] is given as where parameter D is the fractal dimension, ϕ is the domain extension factor for the asperity size distribution, a ′ is the area of a contact spot, a ′ l is the area of the largest contact spot, and A r is the real area of contact.
According to the results in Refs. [23,24], the critical deformation δ c and the elastic and plastic deformations demarcated by the critical area a ′ c are given by where parameter G is a characteristic length of scale of the surface; K is the correlation coefficient [25], which is usually 2.8; H is the hardness of softer materials; and γ is a constant, which is usually 1.5. The composite elastic modulus E * is given by where parameter v is the Poisson ratio. Three deformation stages may occur on the asperity of a rough surface during the contact process, namely, elastic deformation, elastic-plastic deformation, and complete plastic deformation stages.
If δ ≤ δ c , the asperity is an inelastic deformation. Therefore, the elastic load P e on such an asperity, and the contact area of asperity a , are given by the Hertz theory [24] as follows: The range of elastic plastic deformation zone is derived as δ c ≤ δ ≤ 110δ c [26,27]. This stage can be divided into a e = π Rδ, P e = 4 3 ER 1 / 2 δ 3 / 2 . two stages: the first stage of elastic-plastic deformation ( δ c ≤ δ ≤ 6δ c ) and the second stage of elastic-plastic deformation ( 6δ c ≤ δ ≤ 110δ c ). In the elastic-plastic deformation stage, the normal load and contact area of asperity can be written as If δ > 110δ c , the asperity is a plastic deformation. Therefore, the normal load P p and the contact area of asperity a p are given by

Modified Contact Stiffness Model
The ratio of contact force increment to contact distance variation is defined as the contact stiffness, and therefore the stiffness of the asperity under various stages of deformation can be obtained as follows: The total stiffness between the contact surfaces only occurs during the elastic deformation stage, that is, the sum of the stiffness of the three stages of the elastic deformation, the first stage of the elastic-plastic deformation, and the second stage of the elastic-plastic deformation. Therefore, the normal total stiffness between the contact surfaces K n is given by It can be seen from the formula that the contact stiffness K n varies with the actual contact area A r . In fact, the actual contact area A r should be less than the nominal contact area A. Thus, an area proportion coefficient of less than 1 is found, and the relationship between the actual contact area and the nominal contact area can be described as where is assumed as an invariable constant during the contact process.
According to the improved elastic foundation model [7], a pin is equivalent to a rigid wedge. The analytical model for the revolute joints with clearance is illustrated in Figure 1.
The relationship between semi-angle contact ε and the penetration displacement u can be expressed as follows: where parameter ε is the semi-angle contact, u is the penetration displacement, and R 1 and R 2 refer to the radii of the pin and the hole in the plate, respectively. Here, R = R 2 − R 1 represents the radial clearance of the joints.
Because R is much less than R 1 and R 2 , the relationship between l and ε can be expressed as where the parameter l is the length of the chord. If the axial length is b, the nominal contact area A can be expressed as The relationship between the modified stiffness coefficient and the compression depth u can be obtained through formulae (8) and (11) as In addition, the fractal dimensions D and G [22] are given as

Modified Contact Force Model
The L-N model [28,29] considers the energy dissipation during the contact process, the expression of which is as follows: where parameter K is the stiffness coefficient, and D n is the damping coefficient.
In this paper, based on the L-N model expression, the fractal contact force model is established by replacing (15) F n = Ku n + D nu , the stiffness coefficient K with the modified stiffness coefficient K n , and is written in the following form:

Numerical Simulation of Modified Contact Force Model
Using the modified model proposed in this paper, a one-shot contact was analyzed.

Simulation Parameters of the Mechanism
The material and size parameters of the shaft and sleeve were taken as follows: E = 207 GPa, ν = 0.3, σ ѵ = 335 MPa, m = 1 kg, b = 0.01 m, and R 2 = 0.01 m. Because the L-N and hybrid models [6] do not take the roughness of the contact surface into account, roughness Ra in the modified model has a smaller value of 0.05 for unifying the parameters.
With the increase in contact load, the real contact area tends to be a fixed value; however, this value is less than the nominal contact area (approximately 5-10% of the nominal contact area) [30][31][32][33][34]. Therefore, the area proportion coefficient λ in this paper can be valued at 0.05.

Analysis of Contact Process under Different Restitution Coefficients
Assume that the size of the clearance is 0.5 mm, the initial contact velocity is 1 m/s, and the restitution coefficient is different. The simulation results of the modified model are shown in Figure 2. Figure 2(a) shows the characteristics of damping hysteresis during the contact process. The damping of the material will cause a loss of energy, and the restitution coefficient will be smaller, whereas the loss of energy during the contact process will be larger. The modified contact force model presented in this paper can better describe the loss of energy during the contact process when the restitution coefficient is small. Furthermore, Figure 2(b) and (c) also conform to the change in contact force and deformation during actual contact. Therefore, it is possible to describe the contact process with different restitution coefficients.

Analysis of Contact Process at Different Initial Velocities
Assume that the clearance size is 0.5 mm, the restitution coefficient is 0.9, and the initial contact velocity differs. The simulation results of the modified model are as shown in Figure 3. It can be seen from Figure 3 that, with an increase in the initial contact velocity, the contact force and deformation increase rapidly, and the contact time is shortened, which is consistent with the actual contact process. In addition, when the contact force is at maximum, the deformation also reaches the maximum value, which is consistent with the theoretical formula.

Analysis of Contact Process for Different Amounts of Roughness
Assume that the clearance size is 0.5 mm, the restitution coefficient is 0.9, and the initial contact velocity is 1 m/s. The simulation results of the modified model when choosing different amounts of contact surface roughness are as shown in Figure 4. As Figure 4 shows, the roughness of the contact surface is smaller, the contact force is larger, and the time to reach the maximum contact force is shorter. Moreover, the roughness is smaller, as is the maximum deformation; in addition, the time to reach the maximum deformation is also shorter, the reason for which is that a smaller amount of contact surface roughness leads to a larger fractal dimension, and a larger stiffness coefficient. When the initial conditions are the same, the stiffness is larger, the contact force is greater, and the deformation is smaller.

Comparative Analysis of Models
In this paper, a single contact simulation of the L-N model and the hybrid model established in Ref. [6] was carried out, and the results were compared with those of the modified model.

Comparative Analysis of Models under Different Clearance Sizes
Assume that the restitution coefficient is 0.9, the initial contact velocity is 1 m/s, and the clearance size is different. The contact force of the three models varies with the amount of deformation, as shown in Figure 5. As shown in Figure 5, the modified model shows a slight difference with the L-N and hybrid models with a large clearance, whereas the difference between the modified model and the hybrid model is greater than that of the L-N model when the clearance decreases, the reason for which is that the L-N model is based on the Hertz contact theory, and is suitable for the case of noncoordinated contact with a large clearance. For a small clearance, the error in the calculation result is large. The modified model proposed in this paper and the hybrid model are suitable for both large and small clearances, and their calculation results are more accurate. In addition, when the clearance is small, the contact becomes closely conformal contact. Because the modified model described in this paper is based on the fractal theory and micro-surface topography, it is more accurate when calculating a small-clearance contact.

Comparative Analysis of Models under Different Restitution Coefficients
Assume that the clearance size is 0.5 mm, the initial contact velocity is 1 m/s, and the restitution coefficient is different. The contact force of the three models varies with the amount of deformation, as shown in Figure 6. As shown in Figure 6, when the restitution coefficient is large, the calculated results of the three models are relatively similar, whereas when the restitution coefficient is small, the calculated results of the modified model in this paper are similar to those of the hybrid model, which is quite different from the L-N model. The damping force terms in the modified model and the hybrid model can reflect the loss of energy better during the contact process with a smaller restitution coefficient, thereby improving the calculation results. The modified contact force model has a broader range of application.

Experiment Verification
A crash test rig was designed, the results of which prove the correctness of the theoretical analysis.

Test Rig
The test rig is shown in Figure 7, which consists of a slider crank mechanism, a simulation contact part, and the data acquisition system. By controlling the speed of the motor, the slider crank mechanism can make contact block 1 obtain different initial velocities. The simulation contact part can simulate the contact between the shaft and sleeve when designing contact blocks 1 and 2 as convex and concave. The data acquisition system can transfer the signal of the test rig to the data acquisition instrument, and display it in a computer in real time.
During the experiment, contact block 1 was driven by the motor, and impacted contact block 2, which was fixed on the left end. As shown in Figure 8, the maximum contact force during the contact process was measured using a c = 1 mm  acquisition instrument and the charge amplifier, and displayed in real time in the computer. The initial contact velocity measured through the experiment is substituted into the modified model for a numerical calculation, and the other parameters of the modified model are as given in Table 1.

Experimental Results and Analysis
The experiment and numerical results are shown in Tables 2, 3, 4 for different speeds, amounts of surface roughness, and clearances. The experiment data are consistent with the theoretical data as shown in the following tables, which illustrates that the modified normal contact force model proposed in this article is correct. The causes of the error may be from (1) vibrations and