Study on gear contact stiffness and backlash of harmonic drive based on fractal theory

Contact stiffness and backlash in the harmonic drive significantly impact a robot's positioning accuracy and vibration characteristics. The height of the harmonic drive tooth pair is typically less than 1 mm, making the measurement and modeling of backlash and contact stiffness inherently complex. This paper proposes a contact stiffness and backlash model by establishing a correlation between fractal parameters and tooth contact load. To obtain the fractal roughness parameters of the real machined tooth surface, a combination of a noncontact optical profiler and the RMS method is employed. Subsequently, the study explores the influence of rough tooth surface and contact force fractal parameters on contact stiffness and gear backlash. The results demonstrate the substantial impact of surface topography parameters and contact force on contact stiffness and backlash. Specifically, an increase in the fractal dimension correlates with a reduction in gear backlash and contact stiffness. Conversely, the fractal roughness parameter exhibits the opposite effect. Notably, an increase in contact force enhances contact stiffness.


Introduction
Harmonic drives (HDs) find extensive application in robot joints owing to their elevated transmission ratios, superior transmission accuracy, and compact dimensions. 1 The pivotal constituents of a harmonic reducer encompass the circular spline (CS), flexspline (FS), and the wave generator. 2During the rotation of the wave generator, numerous tooth pairs on the CS and FS mesh and drive concurrently.The contact stiffness and backlash in the tooth pair model constitute an integral facet of dynamic research, and their precision directly impacts the positioning accuracy and vibration characteristics of the HD. 3 Traditionally, the gear meshing behavior is conceptualized as an interaction of elastomeric bodies.Furthermore, gear contact is simplified to incorporate spring and damping elements.Hence, the investigation into the contact stiffness and backlash of HDs serves as an exemplification of a dynamic problem.
Fractal contact model can be used to describe contact between machined surfaces.][6] Mandelbrot and Benoit 7 first proposed fractal concept, whose characteristics can be observed at the engineering surface.Furthermore, the author presented definitive overview of the origins of his ideas and their new applications.On this basis, Majumdar and Bhushan 8 further analyzed the contact theory of two-dimensional rough surface, more commonly known as the MB model.Yan and Komevpoulos 9 demonstrated fractal characteristics of an engineering surface.Furthermore, the authors described the variation of contact load and stress area with surface roughness parameters in a bolted connection, thus obtaining size independent contact model.Qi et al. 10 investigated normal contact stiffness fractal model between two spheroidal joint surfaces considering friction.Chen et al. 11 modified the backlash equation using the Weierstrass-Mandelbrot (W-M) function of fractal theory, combined it with the traditional gear torsion model, and proposed an improved nonlinear dynamic model of gear pairs with micro characteristics of tooth surfaces.According to the presented investigations, it can be observed that fractal theory is constantly being refined and upgraded when modeling rough surface contacts.
In addition, many studies have been done on the tribological contact of gears.Chen et al. 12 employed fractal theory to express spur gear backlash.Dynamic performance and backlash were described using a fractal dimension parameter.Zhao et al 13 established an improved model for spur gear meshing characteristics analysis with consideration of fractal surface contact and friction.The author analyzed the effect of fractal parameters on mesh characteristics.Liu et al. 14,15 proposed mathematical model of gear surface fractal parameters, finite element simulation load, and stiffness model.The contact load calculated by the finite element method is brought into the MB contact model to obtain a time-varying meshing stiffness model of involute gears and torsional stiffness model of cycloid gears.Based on the conducted literature survey, it is concluded that fractal theory can be effectively used to model gear contact.Tang et al. 16 studied the effects of tooth surface roughness and contact geometry on lubrication characteristics from a microscopic point of view, and the comprehensive performance can be effectively improved by parameter control.Zhao et al. 17 considered the effect of roughness on the normal contact stiffness of gears and adopted the fractal method to a contact model suitable for the contact of gear pair is established, and the influence of friction coefficient on the contact surface is studied on the basis of this model.Yu et al. 18 proposed a general expression for the coefficient of friction between fractal surfaces in hybrid lubrication, which provides a new idea for the study of gear tribological engineering.According to the energy dissipation mechanism of tangential contact damping, Xiong et al. 19 established a tangential contact damping model of tooth surface considering variable friction under dry friction condition by introducing the correction coefficient of tooth surface. in consideration of the effect of tooth surface micromorphology, Zhao et al. 20 established a coupling analysis model for the wear characteristics of involute gear mixed lubrication.
Current harmonic reducer stiffness studies have mainly focused on torsional stiffness.Dhaouadi et al. 21determined that harmonic gear hysteretic stiffness is mainly determined by friction behavior of multi tooth contact.Furthermore, the authors observed that flexible deformation of FS can only affect nonlinearity of the stiffness curve.Therefore, load torque is expressed as superposition of nonlinear stiffness and friction terms.Chen et al. 22 investigated distribution of meshing forces and loading backlashes under varied transmission torques.Ma et al. 23 studied causes of HD gear time-varying stiffness, number of meshing teeth, and meshing length under different torques.However, torsional stiffness model can only represent meshing performance of multiple teeth, not the meshing characteristics of a single tooth pair.However, there is little literature on contact stiffness analysis that considers multi tooth meshing and surface morphology.The purpose of this study is to reveal the influence of surface morphology on the contact stiffness of harmonic reducers and better guide product production and manufacturing.
Contact between CS and FS tooth profiles of harmonic reducer can be equivalent to the contact model of two rough tooth surfaces.In this paper, a mathematic HD contact stiffness and backlash model based on fractal theory is proposed.In "Introduction" section, relevant theoretical research of fractal theory employed in machined tooth surfaces is introduced.Furthermore, the reason why fractal theory is used in this paper is explained.In "Mathematical modeling" section, gear contact stiffness and backlash model is established based on fractal parameters of machined tooth surface and contact theory.Then, graphics and tables are employed for case analysis in "Experimental setup" section.Effect of the instantaneous load and fractal parameters on contact stiffness and backlash of surfaces is analyzed in "Results and discussion" section.Lastly, conclusions are drawn in "Conclusions" section.

Mathematical modeling
Owing to their superior transmission accuracy and robust stability, harmonic reducers find extensive application in industrial robot joint systems.A standard harmonic reducer is composed of three components: a FS with a cylindrical cup, a CS featuring an inner tooth rim, and a wave generator equipped with a flexible bearing.During the rotation of the harmonic reducer, 30% of all teeth engage in contact simultaneously.Notably, the height of each contact tooth pair is limited to less than 1 mm, with the majority measuring approximately 0.5 mm in height, as depicted in Figure 1.This characteristic significantly contributes to the intricacy of measuring tooth profile and roughness parameters.

Rough contact tooth surface modeling
The contact model of a single tooth pair in a harmonic reducer is depicted in Figure 2. Contrary to the appearance of smooth contact regions, these areas actually comprise numerous asperities.These spherical asperities represent micro features on the surface resulting from the gear machining process.Applying the M-B contact theory, the interaction of two rough contact surfaces can be likened to one rough tooth surface engaging with another that is smooth.Similarly, the contact model involving two asperities can be simplified as the contact of one asperity with a plane.The actual contact area and deformation of each asperity are contingent on its size and the normal force acting upon it.
Mathematical properties of machined tooth surfaces satisfy two-dimensional W-M function 15 : Parameters employed in equation ( 1) are shown in Table 1.For a single asperity, there are four important parameters: curvature radius R, normal deformation δ, real contact area a, and truncated area a ′ .Elastic or plastic deformation of an asperity can be distinguished by analyzing the relevant parameters.The corresponding relationship between geometric parameters of deformed asperity can be expressed as 24 : (2) where the radius of truncated asperity can be described as r ′2 = 2Rδ.The corresponding maximum elastic deformation δ c and contact area a ′ c are key indexes to evaluate whether the plastic deformation occurs or not, which can be given by: where parameter E is the equivalent elastic modulus of two contact tooth pairs.H is the hardness value of the softer tooth surface.Parameter k can be obtained by an empirical formula related to the Poisson ratio k = 0.454 + 0.41ν.
Three-dimensional distribution function model of asperity can then be defined as: where a ′ l indicates the maximum cross-sectional area of the tooth surface asperity, and ψ is related to the fractal dimension which can be obtained as:

Gear contact stiffness modeling
The classical torsional vibration model of a gear pair is illustrated in Figure 3.The time-varying stiffness (K m ) is tantamount to the summation of the contact stiffness of all tooth pairs along with the meshing line.In the context of a HD, where multiple tooth pairs engage in simultaneous meshing, the utilization of a single meshing line is not feasible.Therefore, this article only establishes a contact stiffness model for a single tooth pair, without considering damping effects.
During gear tooth meshing process, contact force is generated at the top of numerous asperities.According to different deformation stages, it can be divided into elastic contact force f e = 4/3ER 1/2 δ 3/2 and plastic contact force f p = Ha ′ .Total elastic normal load F e and plastic normal load F p can be obtained by integrating a distribution function of the asperity in equation (6): Total normal load can be written as: Micro asperities on the contact tooth surface produce contact force regardless of elastic or plastic deformation, but the contact stiffness only exists in the plastic deformation stage.Elastic contact stiffness k e is a reciprocal value of elastic contact load f e to normal deformation δ: In a similar manner, normal contact stiffness K ne is the integral of the contact stiffness of a single asperity k e and the distribution function n(a ′ ): According to a difference in maximum cross-sectional area a ′ of asperities under different load conditions, the normal contact stiffness K m of the tooth surface can be written as:

Gear backlash modeling
Contact surface of machined gears produces clearance without considering profile modification.Therefore, gear backlash is defined by normal clearance distribution function, which has already been previously verified. 12unction f (x)is the backlash function, which is usually written according to equation (14).Figure 4 shows the structure and backlash of a pair of meshing gears: According to Chen et al., 12 backlash b represented by W-M function significantly affects dynamic characteristics of gears.In this paper, W-M function is employed to describe the backlash function.Backlash b described by two-dimensional W-M function can be written as: Parameters presented in equation ( 15) are already defined according to equation (1).

Fractal roughness and fractal dimension parameter identification method
As per equations(1)-( 15), the determination of the fractal roughness parameter (G) and fractal dimension (D) is imperative to acquire the contact stiffness (K N ) and backlash (b).The fractal parameters characterizing the rough tooth surface can be computed using surface data obtained through a noncontact optical profiler.Numerous literatures have presented the calculation method for fractal dimension (D) and G. Notably, in accordance with Zhang et al. 25 and Feng et al., 26 employing the roughness length (RMS) methods yields a more accurate characteristic parameter for the tooth surface.Consequently, the experimental research in this paper is conducted based on this model.
(1) According to Feng et al., 26 the structure-function S(ξ) in RMS method can be represented as: where <> represents the mean value, and z(x) is the same as WM function in equation (1).If equation ( 1) is combined with equation ( 16), it be expressed as: Fractal roughness G = <G s >, parameter C is equal to: where Γ is the gamma equation of two-dimensional profile fractal dimension D s .The relationship between a two-dimensional fractal dimension D s and threedimensional fractal dimension D is: D = <D s >+1.
According to equation ( 17), logarithm of S(ξ) and ξ can be used to obtain D s and G s via least-square fitting.Then, threedimensional fractal dimensions D and G can be obtained from two-dimensional fractal dimensions D s and D s . 4

Experimental setup
Gear hobbing and shaping experiments are conducted to validate the proposed gear contact model, and the topography characteristics of the machined tooth surface were analyzed.This section provides a detailed discussion of the experimental setup and the tooth surface detection process.The YK3610 hobbing machine tool (manufactured in Ningjiang, China) is employed for gear profile hobbing tests, as illustrated in Figure 5.A 40CrNiMaA gear with 200 teeth, featuring a 0.25-module double circular arc profile and a 10.5-mm face width, is machined.The parameters of the FS tooth surface are collected using a noncontact optical profiler (ST-400).
Similarly, a gear profile shaping machine tool (YGK512, manufactured in Yichang, China) is utilized for hobbing tests.Before tooth profile scanning, the CS inner ring gear was separated by wire cutting, as depicted in Figure 6.
The details of the machined gear workpiece and the noncontact optical profiler are listed in Tables 2 and 3, respectively.Based on the results obtained from optical instrument scanning of the tooth profile, tooth root arc, and tooth top arc, the data are more comprehensive than those of the tooth profile alone.Since the same tool under the same working conditions processes both the tooth root and tooth profile, the roughness parameter of the root arc surface is utilized to represent the roughness parameters of the entire tooth profile, as depicted in Figure 7.The region extracted from the rough tooth surface data is henceforth referred to as the "Interested tooth surface."Parameters related to tooth surface roughness, specifically S q and S a , can be directly extracted using the ST400.The RMS method curve should be calculated and plotted according to the equation in "Fractal roughness and fractal dimension parameter identification method" section.
Final experimental results and test are listed in Table 4.

Results and discussion
Tooth surface topography characteristics are discussed in this section.Effect of rough tooth surface fractal parameters and contact force on contact stiffness and gear backlash is studied.

Effect of fractal parameters on backlash b
Simulation of backlash b with different fractal dimensions D is shown in Figure 8.With an increase in fractal dimension, tooth surface becomes smoother.As shown in Table 4, average surface backlash value is used as the evaluation basis.In a similar manner, the smoother the tooth surface, the smaller the contact clearance.Nevertheless, the prevailing manufacturing accuracy for small module gears typically falls within the range of 5-7 levels.In theory, a larger fractal dimension correlates with a smaller K.However, owing to the impact of manufacturing accuracy, Figure 8 demonstrates that a fractal dimension of approximately 2.4 is relatively optimal.In Figure 9, backlash variation with different fractal roughness parameters G is shown.When G changes from 1E-9 to 1E-15, the tooth surface becomes smooth and backlash becomes uniform and stable.According to Table 5, G should be less than 1E-13 to ensure that the backlash reaches micron level.

Effect of fractal parameters on contact stiffness K
Different fractal dimension and roughness parameters also affect the contact stiffness of the tooth surface.Contact load between CS and FS tooth of the harmonic reducer is assumed as 50 N.Considering that the actual fractal dimension D of machined tooth surface is between 2.3 and 2.6, stiffness variation law is discussed in this range.Variation of contact stiffness with tooth surface parameters is shown in Figure 10.According to Figure 10(a), the data are interpolated and fitted, contact stiffness decreases with an increase in fractal dimension D. Coupled with the tooth surface morphology illustrated in Figure 8, it can be deduced that an escalation in micro-convex body density leads to a reduction in backlash and tooth surface contact stiffness.
Figure 10(b) illustrates the variation in contact stiffness (K) as the tooth surface fractal roughness parameter G changes, as detailed in Table 5. Notably, an augmentation in the fractal roughness parameter G corresponds to an enhancement in contact stiffness.However, currently the manufacturing accuracy of small module gears is usually around 5-7 levels.In theory, the larger the fractal dimension, the smaller the K.However, due to the influence of manufacturing accuracy, it can be seen from Figure 8 that a fractal dimension of around 2.4 is relatively suitable.

Effect of contact force on contact stiffness K
As previously mentioned, 30% of all teeth are simultaneously in contact in typical harmonic reducers.The contact position is considered to be on the indexing circle.Contact load of a single tooth pair under rated   working load can be expressed as: where T denotes the load torque, N is the number of contact tooth pairs, and d is diameter of the dividing circle.
Based on the provided data in "Experimental setup" section and Table 1, contact stiffness curves under different contact loads are investigated.Rated operating torque of this type of harmonic reducer is 52 Nm.Contact load under rated condition is 34.67 N. Therefore, contact load range can be set between 10 N and 100 N to evaluate contact stiffness variation, as shown in Figure 11.As illustrated in Figure 11, contact stiffness exhibits a linear increase in tandem with contact load.Specifically, the contact stiffness of the harmonic reducer gear surface reaches 6 × 10 7 under rated conditions.The experimentally acquired fractal parameters, denoted as D and G, are detailed in Table 3.

Conclusions
In this paper, mesh characteristics of harmonic gear contact tooth surfaces are theoretically and experimentally investigated.The following conclusions are drawn: 1. Contact teeth surface of harmonic reducer can be described by M-B contact model.Furthermore, gear contact stiffness and backlash are functions related to fractal dimension D, fractal roughness parameter G, and contact load F.
2. The results show that contact stiffness decreases with an increase in fractal dimension D. In addition, an increase in fractal roughness parameter G and contact force F contributes to improving the contact stiffness.3. Increasing fractal dimension D and decreasing fractal roughness parameter G can reduce the contact backlash.
In addition, this model is not only suitable for modeling the contact stiffness and backlash of HDs but also for identifying meshing parameters of gears with relatively small modules.

Figure 1 .
Figure 1.Core parts and tooth profile of harmonic reducer.

Figure 2 .Figure 3 .
Figure 2. Macroscopic and microscopic structure of tooth surface contact area.

Figure 4 .
Figure 4. Backlash function model of a gear-pair.

Figure 8 .
Figure 8. Simulation of backlash (b) with different values of fractal dimensions (D).

Table 4 .
Experimental and test results.

Table 5 .
Average backlash for different values of D and G.