Modeling and analysis of nonlinear axial force generated by automotive drive-shaft systems based on fractal theory

: The nonlinear generated axial force (NGAF) of an automotive drive-shaft system is mainly induced by the tripod joint in the drive-shaft system, which will lead to vibration and noise of a vehicle. On the basis of fractal theory, this paper conducts modeling and analysis of the NGAF from the micro level to study the NGAF more effectively. In order to describe the micro contact and friction states between rollers and tracks inside the tripod joint more accurately, a correction method for correcting the distribution function of asperities between rollers and tracks is proposed. Based on the proposed correction method and a calculation model of the NGAF, a fractal model of the NGAF considering rough surface characteristics is established. The fractal model is mainly related to the fractal parameters, the material parameters, the operating conditions and the correction coefficient of the distribution function of asperities. Using Sobol' global sensitivity analysis method, the first order and overall global sensitivities of the fractal model of the NGAF are subsequently analyzed to determine the influence of the factors on the NGAF. The effectiveness of the fractal model and the analytical method for the NGAF is verified by the experiment and calculation examples. The results show that modeling and analysis of the NGAF considering rough surface characteristics are essential to reveal the relationships between the NGAF and the fractal parameters, the material parameters and the correction coefficient from the micro level, which is helpful for the design and analysis of the NGAF of the drive-shaft system.


Introduction
As an important part of automotive transmission systems, the drive-shaft system can transfer the torque from the transmission to the wheel, and ensure the constant speed of the input and output ends at an articulation angle. In general, due to the friction and motion characteristics inside the tripod joint of driveshaft systems, a periodic axial force, namely, the nonlinear generated axial force (NGAF), will be generated in the operation of drive-shaft systems [1][2][3]. During the driving process of a vehicle, the NGAF will cause vibration and noise [2][3][4][5][6], thus affecting the comfort of the vehicle. In order to reduce the NGAF of driveshaft systems, many researchers have studied the NGAF of drive-shaft systems, as well the researches on the NGAF mainly focus on the experimental researches [1][2][3][4] and the modeling analyses [5][6][7][8].
Since the NGAF is mainly induced by the friction inside the tripod joint, Lee et al. [1][2] developed a special test bench for testing the friction characteristics between rollers and tracks inside the tripod joint under different operation conditions, as well the friction characteristics were modeled and analyzed. In the study, the relationship between the friction characteristics inside the tripod joint and the NGAF was revealed.
However, these friction characteristics cannot reflect the influence of the micro contact and friction characteristics between rollers and tracks inside the tripod joint on the NGAF. In order to study and analyze the NGAF more accurately, Sa et al. [3][4] established a test bench for measuring the NGAF of drive-shaft systems, and studied the influence of the input torque, the type of the tripod joint, the articulation angle and the grease on the NGAF from a macro perspective.
Serveto et al. [5] established an analytic model and a multi-body dynamics model of the NGAF of drive-shaft systems. The two models were mainly related to the input torque, the radius of the tripod joint, the input torque and the Coulomb friction coefficient. In the study, the Coulomb friction coefficient was obtained through minimization of errors between the computed responses and the measured NGAF.
Through the kinematic analysis of the tripod joint, Jo et al. [6] established a calculation model for the NGAF of drive-shaft systems with a tripod joint. In the model, two friction models (a pure sliding friction model and a friction model that considers both sliding and rolling friction) were used to characterize the friction.
The friction coefficients in the friction models were identified by the measurement. Lim et al. [7][8] established multi-body dynamic models of drive-shaft systems with a tripod joint to analyze the NGAF and the dynamic force inside the tripod joint. In the multi-body dynamic models, the Coulomb friction model was used to characterize the friction, while the impact function (which is mainly related to the contact stiffness, the force exponent and the penetration displacement of contact pairs) was used to calculate the normal contact force. However, neither the Coulomb friction model nor the impact function could reflect the micro friction and contact characteristics between contact pairs.
In the aforementioned studies, neither the experimental research nor the modeling analysis of the NGAF considered the influence of micro contact and friction characteristics between the contact pairs inside the tripod joint on the NGAF. In addition, the friction coefficient, which has great influence on the NGAF, was generally identified through minimization of errors between the computed responses and the measured data from a macro perspective. The friction coefficient obtained by the measurement identification, however, cannot directly reflect the influence of the micro characteristics of rough surfaces on the NGAF. Analyzing the influence of the micro characteristics of the rough surface on the NGAF allow us to find the root cause of the NGAF, which is beneficial to the design and optimization of the NGAF of the drive-shaft system. Thus, in order to study the NGAF more effectively, it is essential to consider the influence of the rough surface characteristics between rollers and tracks inside the tripod joint on the NGAF.
Since fractal parameters (which are used to characterize rough surface profiles) of fractal theory do not depend on the sampling length and the instrument resolution [9,10], the fractal theory has a good applicability in the study of mechanical rough surfaces, and thus the fractal theory is usually used to analyze the micro characteristics of contact pairs in practical engineering [11][12][13][14][15]. In the fractal theory analysis, the contact between two rough surfaces is generally regarded as the contact between a flat plane and a rough plane, as well the distribution function of the contact asperities is used to describe the distribution of the contact asperities on the rough plane [16,17]. The contact between rollers and tracks inside the tripod joint, however, is a contact between a spherical surface and a cylindrical surface. Thus, it is necessary to modify the distribution function of the contact asperities so as to obtain the distribution of the contact asperities between rollers and tracks more accurately. In addition, the sensitivity analysis of the NGAF can be performed to obtain the influence of various factors on the NGAF, thus providing guidance for the design and optimization of the NGAF of drive-shaft systems. The primary goals and the contributions of the study are summarized below.
Based on fractal theory and a calculation model of the NGAF, a fractal model of the NGAF considering rough surface characteristics is established. The fractal model of the NGAF is mainly related to the fractal parameters (namely, the fractal dimension and the characteristic length scale), the yield strength, the elastic modulus, the Poisson's ratio, the articulation angle, the shaft angular position of the drive-shaft system and the correction coefficient of the distribution function of contact asperities between rollers and tracks. To describe the contact and friction states between rollers and tracks inside the tripod joint more accurately, a correction method for the distribution function of contact asperities between rollers and tracks is proposed.
Using Sobol' global sensitivity analysis method, the first order global sensitivity and the overall global sensitivity of the fractal model of the NGAF are subsequently analyzed to determine the influence of various factors on the NGAF.
The remaining paper is organized in four different sections. In Section 2, a calculation model of the NGAF is introduced. The fractal model of the NGAF of the drive-shaft system is proposed in Section 3 together with the correction method for the distribution function of contact asperities between rollers and tracks. In Section 4, the measurement method for the NGAF of a drive-shaft system is briefly described together with the model verification, as well the analysis of effects of the factors on the NGAF and the sensitivity analysis of the fractal model of the NGAF are performed through numerical examples. The key conclusions of the study are summarized in Section 5.

Structures of a drive-shaft system and a tripod joint
A drive-shaft system is mainly composed of a ball joint, a tripod joint and a drive-shaft, as shown in Fig. 1(a). The ball joint is connected to the wheel, while the tripod joint is connected to the transmission.
In a drive-shaft system, the NGAF is mainly induced by the tripod joint. The structure of a tripod joint is shown in Fig. 1(b). A tripod joint generally comprises rollers (1), needles (2), a tripod (4) and a housing (6).
There are three tracks (5) parallel to the input shaft for the rollers to roll and slide, as well there is a certain amount of grease in the housing. The tripod has three trunnions (4), 120 degrees apart, which are used for assembling needles and rollers.

Main kinematics results of the drive-shaft system with a tripod joint
The NGAF is related to the kinematics of the drive-shaft system with a tripod joint, and the kinematics have been described in [5,18]. Briefly, the schematic diagram of the kinematic model of a drive-shaft system is shown in Fig. 2 where r is the pitch circle radius (PCR) of the tripod joint; e is the offset, which is defined as the distance between point C and axis z1; φi, z and e are expressed as [5]: where Dd is the distance between point A and axis z1; L is the length of the drive-shaft. Taking roller C1 in the Fig. 2 as the research object, its movement on the track is shown in Fig. 3. The axis system of the track (i0, j0) and that of the roller (i1, j1) correspond to the instantaneous position, and Co is the origin of both the axis systems. In addition, Vc is the velocity of roller C1 in direction j0, as well 1 is the deflection angle of roller C1. During the operation of the drive-shaft system, the roller slides in direction i1 and rolls in direction j1, and thus there are sliding friction and rolling friction between the roller and the track in direction i1 and direction j0, respectively. Since the sliding friction coefficient between the roller and the track is much greater than the rolling friction coefficient [2,5], only the sliding friction between the roller and the track was considered in this paper.
In addition, the relationship among the deflection angle λ1, the articulation angle δ and the shaft angular position φ1 is obtained as [18]:

A calculation model of the NGAF of the drive-shaft system
The schematic diagram of the normal force of rollers acting on tracks is shown in Fig. 4.
are the normal forces of rollers Ci acting on the tracks, as well T is the input torque. In an ideal condition, the offset e in Fig. 2 can be ignored, and according to Eq. (1), it is obtained: Thus, the normal forces (i=1,2,3) of rollers Ci acting on the tracks are: where Fn is the normal force between the roller and the track in the ideal condition. According to Fig. 3, when roller C1 slides along direction i1, the sliding friction force is obtained as： where μ is the sliding friction coefficient; and sign(·) denotes the sign of the velocity 1 .
And then the component of the friction force 1 in the direction j0 is: Thus, combining Eqs. (3), (4), (7) and (8), the component of the friction force 1 along the direction j0, namely, the NGAF between roller C1 and the track, can be expressed as: In general, the input angle φ is equal to the shaft angular position θ, and according to Eq. (2), φ1 = θ.
Eq. (9) thus can be obtained as: Similarly, the NGAFs Fi (i=2,3) between rollers Ci and the tracks can be obtained as follows: where Ai is: The total NGAF of the drive-shaft system, thus, is obtained: It is seen from Eqs. (12) and (13) that the magnitude of the NGAF mainly depends on the friction coefficient μ, the normal force Fn, the articulation angle δ and the shaft angular position θ of the drive-shaft system. In addition, the NGAF changes periodically with the shaft angular position. When the shaft angular position is 90 o , the magnitude of the AFF reaches the maximum [18]. In practical engineering, most attention is paid to the maximum magnitude of the AFF, and thus the shaft angular position is set as 90 o in the following analysis.

A fractal model of the NGAF of the drive-shaft system
In Eq. (13), the normal force and the friction coefficient are macro influencing factors, which cannot directly reflect the influence of the micro contact and friction characteristics between rollers and tracks on the NGAF. In order to study the NGAF more deeply and effectively, fractal models for calculating the normal force and the friction coefficient considering rough surface characteristics are first proposed based on fractal theory, and then the fractal model for calculating the NGAF is proposed by combining Eq. (13) and the fractal models of the normal force and the friction coefficient. In addition, since the contact between rollers and tracks is a contact between a spherical surface and a cylindrical surface, a correction method for the distribution function of contact asperities between rollers and tracks is proposed to describe the contact state between rollers and tracks more accurately.

Fractal function, and deformation characteristics of asperities
The two-dimensional surface profile of the mechanical rough surface can be characterized by the Weierstrass-Mandelbrot (W-M) fractal function [16,19]: The meanings of symbols in Eq. (14) are shown in Table 1. According to the research in [10,16], = 1.5 is a suitable value for high spectral density and phase randomization. In addition, D and G are important parameters affecting the roughness and rough surface characteristics of a rough surface. For a single contact spot on the rough surface, the schematic diagram of its deformation is shown in Fig. 5. In Fig. 5, σ is the deformation displacement of the asperity; R is the radius of the asperity base; ra is the radius of the actual contact area; rn is the radius of the nominal contact area. The deformation displacement σ depends on the values of the peak and trough of the surface profile, and σ can be expressed as follows [16,19]: According to Fig. 5, the radius R of the asperity can be expressed as: Since R is much greater than σ during the deformation of the asperity, according to Eq. (16), it is obtained: The nominal contact area an can be expressed as follows: Combining Eqs. (15), (17) and (18), it is obtained: In general, the deformation of the asperity can be divided into two stages: the elastic deformation stage and the plastic deformation stage. When the contact area of the asperity is greater than the critical contact area, the asperity is in the elastic deformation; when the contact area of the asperity is smaller than the critical contact area, the asperity is in the plastic deformation stage [16]. The critical contact area can be expressed as [16]: where, K is the ratio of the hardness H of the softer material to its yield strength Y; q is related to the Poisson's ratio of the softer material; and E * is the equivalent elastic modulus.
where E1 and E2 are the elastic moduli of the two contact materials, respectively, as well v1 and v2 are the Poisson's ratios of the softer and harder contact materials, respectively.
In the elastic deformation stage, according to Hertz theory, the normal contact force Fe and the contact area ae of the asperity can be expressed as [10]: Substituting Eqs. (15) and (19) into Eq. (22), it is obtained: In the plastic deformation stage, the normal contact force Fp and the contact area ap of the asperity can be expressed as [10]:

Correction method for the distribution function of asperities between the roller and the track
In general, the contact of two rough surfaces can be regarded as the contact between a flat plane and a rough plane. In the contact process, the relationship between the number n of asperities and the contact area a of asperities, that is, the distribution function n(a) of asperities, can be expressed as [20,21]: where, al is the nominal maximum contact area of asperities, Ar is the actual contact area of the two contact surfaces; and ψ is the correction factor for the ratio of the actual contact area Ar to the maximum contact area al of asperities [21].
However, since the outer surface of the roller and the inner surface of the track are a spherical surface and a cylindrical surface, respectively, the contact between the roller and the track can be regarded as a point contact between a spherical surface and a cylindrical surface, as shown in Fig. 6. In the figure, O and Fn are the contact point and the normal force between the roller and the track, respectively; R1 and R2 are the radius of the track cross section and the radius of the roller, respectively; α and Lt are the arc angle of the track cross section and the length of the track, respectively. It is seen from Fig. 6 that the contact between the roller and the track inside the tripod joint is the contact between two curved surfaces, and thus Eq. (25) is not applicable to the contact between the roller and the track. To reflect the distribution of asperities between the roller and the track more accurately, it is necessary to improve Eq. (25). For the contact between two rough curved surfaces, the improved distribution function * ( ) of asperities can be expressed as [22,23]: where χ is the correction coefficient of the improved distribution function * ( ) of asperities, and it can be expressed as [22,23]: where, St is the theoretical contact area between two contact surfaces; Kh is the comprehensive curvature coefficient; and is the sum of the surface areas of the two contact surfaces: where S1 and S2 are the surface areas of the two contact surfaces, respectively. In addition, "+" is for the external contact, while "-" is for the internal contact [23].

Comprehensive curvature coefficient between the roller and the track
The schematic diagram of the static contact between two curved surfaces is shown in Fig. 7 And then the comprehensive curvature coefficient between two curved surfaces is obtained: It is seen from the Fig. 6 that the contact between the roller and the track can be regarded as an internal contact between a spherical surface and a cylindrical surface. According to Eq. (29) and Fig. 6, the equations of the spherical surface and the cylindrical surface near the contact point O are obtained, respectively: Substituting Eqs. (30), (31) and (33)

Correction coefficient for the distribution function of asperities between the roller and the track
Combining Eqs. (35) and (36), the theoretical contact area between the track and the roller is obtained: According to Eq. (28) and Fig. 6, the sum of the surface areas of the spherical surface and the cylindrical surface is: where the outer surface of the roller is regarded as a complete spherical surface.
In addition, since the materials of the roller and the track are basically the same, the elastic modulus and Poisson's ratio of the roller are generally regarded as equal to the elastic modulus and Poisson's ratio of the track in the analysis [25]: According to Eqs. (21) and (39)

Fractal model of the friction coefficient between the roller and the track
The friction coefficient μ in Eq. (13) is another important factor affecting the NGAF, which is related to the normal and tangential contact forces between the roller and the track. When a single asperity yields, the tangential contact force of the asperity can be expressed as [26]: According to Eqs.

A fractal model of the NGAF of the drive-shaft system
Combining Eqs.

Experiments and numerical analysis examples 4.1 Identification of the fractal dimension and the characteristic length scale
The fractal dimension D and the characteristic length scale G are important parameters of the fractal model of the NGAF. Before the validation of the fractal model of the NGAF, the fractal dimension D and the characteristic length scale G need to be identified. In general, the fractal dimension and the characteristic length scale can be identified by combining the measured rough surface profile and the structure-function method [27][28][29]. However, this identification process is complicated. In order to identify the fractal parameters between the roller and the track more effectively, this paper first identifies the fractal dimension through the approximate relationship between the fractal dimension and the roughness, and then identifies the characteristic length scale through minimization of errors between the measured friction coefficient and the friction coefficient calculated by the fractal model of the friction coefficient (Eq. (48)).
According to the literature [30], the fractal dimension D of the rough surface can be approximated by the following formula:  Table 2, which are provided by the manufacturer [25]. Thus, when μ=0.31, D=1.54 and other input parameters of the fractal model are shown in Table 2, the characteristic length scale G can be identified as 5.1×10 -10 m.

Measurements and model verification
To validate the fractal model of the NGAF of the drive-shaft system, a test bench is designed for measuring the NGAF of the drive-shaft system. The pictorial illustration and the schematic diagram of the test bench are shown in Figs. 8(a) and 8(b), respectively. The detailed measurement method for the NGAF of a drive-shaft system have been introduced in [31]. In brief, the ball joint is connected to a control motor, while the tripod joint is connected to a hydraulic actuator, as shown in Fig. 8(b). During the measurement, the hydraulic actuator applies an input torque, while the control motor applies a rotation speed. By adjusting the relative position of the hydraulic actuator end and the control motor end, a specified articulation angle can be determined. At the hydraulic actuator end, four force sensors are evenly installed to measure the NGAF of the drive-shaft system, as shown in Fig. 8(b).  Taking a drive-shaft system as the test object, the dimension parameters, the material parameters and the input torque of the drive-shaft system are shown in Table 2. Through the identification in Section 4.1, the fractal dimension and the characteristic length scale are obtained as 1.54 and 5.1×10 -10 m, respectively.
Substituting the parameters shown in Table 2  can be obtained on the basis of the test bench of the drive-shaft system, and the measured results are also shown in Fig. 9. During the measurement, the input torque is set as 100Nm, as well the rotation speed is set as 200rpm [25]. In addition, since the influence of grease is not considered in the proposed fractal model of the NGAF of the drive-shaft system, the grease inside the tripod joint is cleaned. It is seen from Fig. 9 that the calculated NGAF is basically consistent with the measured NGAF, indicating the effectiveness of the NGAF fractal model of the drive-shaft system.

Effects of factors on the NGAF
According to the analysis in Section 3.5, the main factors influencing the NGAF include the fractal parameters (namely, the fractal dimension and the characteristic length scale), the material parameters (namely, the elastic modulus, the Poisson's ratio and the yield strength), the operating conditions (namely, the articulation angle and the shaft angular position) and the correction coefficient. Since the effects of the operating conditions on the NGAF have been analyzed in [18,31], this section only analyzes the effects of the fractal parameters, the material parameters and the correction coefficient on the NGAF. The initial analysis values of the input parameters of the NGAF fractal model (Eq. (49)) are shown in Table 3.
In Table 3, the size parameters and material parameters related to the tripod joint are the design parameters of an existing tripod joint, which are supplied by the manufacturer [25]; the initial analysis value of the input torque remains consistent with the test condition in Section 4.2; and 6 o is the design value of the articulation angle of the drive-shaft system [25]. According to the analysis in Section 2.3, the magnitude of the AFF reaches the maximum when the shaft angular position is 90 o , and thus the initial analysis value of the shaft angular position is set as 90 o . In addition, the correction coefficient is calculated by substituting the size parameters and material parameters related to the tripod joint as well as the input torque into Eq.
(41), as well the fractal dimension and the characteristic length scale are obtained by identifying the rough surface characteristics between the rollers and tracks inside the tripod joint, as described in Section 4.1.
It should be pointed out that the correction coefficient, the fractal dimension and the characteristic length scale are important parameters that characterize the real rough surface characteristics between the roller and the track. The correction coefficient corrects the contact shape and contact state between the roller and the track by modifying the distribution of asperities between the roller and the track, as well the fractal dimension and the characteristic length scale directly reflect the roughness and the rough surface profiles between the roller and the track. (1) Effects of the fractal dimension and the characteristic length scale on the NGAF When magnitudes of the characteristic length scale are 1×10 -10 m, 5×10 -10 m and 10×10 -10 m, respectively, and other parameters are shown in Table 3, the relationship between the NGAF and the fractal dimension is shown in Fig. 10(a). It is seen that with the increase of the fractal dimension, the NGAF first increases and then decreases. When the fractal dimension is about 1.63, the NGAF reaches the maximum value. According to the analysis in [10], with the increase of the fractal dimension D, the critical contact area decreases, and then a large number of asperities in the elastic stage increase; when the fractal dimension increases to a specific value, the asperities in the elastic stage increases to the limit; as the fractal dimension continues to increase, the percentage of asperities in the elastic stage decreases since a large number of asperities smaller than the critical contact area increase.
Thus When fractal dimensions are 1.55, 1.65, and 1.75, respectively, and other parameters are shown in Table 3, the relationship between the NGAF and the characteristic length scale is shown in Fig. 10(b). It is seen that with the increase of the characteristic length scale, the NGAF decreases. This is because with the increase of the characteristic length scale, the radius of asperities decreases and the percentage of asperities in the elastic stage decreases [10], and thus the tangential contact force between the roller and the track will decreases. The NGAF thus can be reduced by reducing the characteristic length scale that affects the profile of the rough surface.
(2) Effects of the yield strength, the Poisson's ratio and the elastic modulus on the NGAF (3) Effects of the correction coefficient on the NGAF Figure 12 illustrates the variation in the NGAF with the variation in the correction coefficient for two different fractal dimensions (1.55 and 1.65). The results suggest that the NGAF decreases with the increase of the correction coefficient. Under a certain load, when the correction coefficient increases, it means that the contact area between the roller and the track increases, which lead to an increase in the number of asperities and then a decrease in the load borne by each asperity. If the load borne by each asperity decreases, its contact area is likely to be smaller than the critical contact area, and thus more asperities are in the plastic deformation stage. Thus, with the increase of the correction coefficient, the tangential contact force between the roller and the track decreases, resulting in the decrease of the NGAF. In summary, in addition to the yield strength, other influencing factors, namely, the fractal parameters, the Poisson's ratio, the elastic modulus and the correction coefficient, all have a nonlinear effect on the NGAF. Moreover, the influencing factors (including the fractal parameters, the material parameters and the correction coefficient) affect the magnitude of the NGAF by affecting the percentage of contact asperities between the roller and the track in the elastic stage. According to the above analyses, it is seen that the relationships between the NGAF and the influencing factors can be studied and revealed from the micro level on the basis of the proposed fractal model of the NGAF, which will provide an important reference for reducing the NGAF and the design of the drive-shaft system.

Sensitivity analysis of the fractal model of the NGAF
To further determine effects of the factors in Section 4.3 on the NGAF, the sensitivity analysis of the fractal model (Eq. (49)) of the NGAF is carried out by using Sobol' global sensitivity method [32]. The main ideal of Sobol' global sensitivity method is to decompose the model into the sum of 2 n incremental terms, and then calculate the deviation and total variance of the model response through sampling, so as to obtain the first order global sensitivity and the overall global sensitivity of the model. In the sensitivity analysis, the value ranges of the fractal dimension, the characteristic length scale, the elastic modulus, the Poisson's ratio, the yield strength, the articulation angle, the shaft angular position and the correction coefficient are shown in Table 4, as well other input parameters of the fractal model (Eq. (49)) are shown in Table 3.
When 1,000,000 random samples are sampled for the influencing factors in Table 4, the sensitivity analysis results are shown in Table 5. It is seen that comparing with the fractal dimension, the shaft angular position, the Poisson's ratio, the yield strength and the correction coefficient, the first order global sensitivity coefficients of the characteristic length scale, the elastic modulus and the articulation angle are larger, and thus these factors have a greater influence on the NGAF. In addition, for the characteristic length scale, the elastic modulus and the articulation angle, the absolute differences between their first order global sensitivity coefficients and their overall global sensitivity coefficients are large, which indicates that there is a large interaction among these influencing factors on the NGAF. According to the overall global sensitivity, influences of the considered factors on the NGAF are ranked from top to bottom as follows: the articulation angle, the characteristic length scale, the elastic modulus, the yield strength, the fractal dimension, the shaft angular position, the Poisson's ratio and the correction coefficient.

Conclusion
The method for correcting the distribution of asperities between the roller and the track as well as the fractal model of the NGAF proposed in this paper are effective and useful for revealing the relationships between the NGAF and the influencing factors (including the fractal parameters, the material parameters and the correction coefficient) from the micro level, which will provide an important reference for reducing the NGAF and the design of the drive-shaft system.
The influencing factors affect the magnitude of the NGAF strongly by affecting the percentage of contact asperities between the roller and the track in the elastic stage. With the increase of the fractal dimension, the NGAF first increases and then decreases. By improving the machining precision of contact surfaces between the roller and the track, the NGAF of a drive-shaft system can be reduced when the fractal dimension is large. The results suggest decrease in the NGAF with increasing the characteristic length scale, and thus the NGAF can be reduced by reducing the characteristic length scale. The NGAF increases with the increase of the yield strength, while the NGAF decreases with the increase of the Poisson's ratio and the elastic modulus. Thus, the NGAF also can be reduced effectively by reducing the yield strength and increasing the Poisson's ratio and the elastic modulus. By changing the material and size parameters of the roller and the track, the correction coefficient can be increased, and then the NGAF can be reduced.

Replication of results
The results of the above numerical analysis can be obtained through the Matlab code. The relevant Matlab code can be seen in the supplementary material.