Analysis of ContactMechanical Characteristics of Flexible Parts in Harmonic Gear Reducer

Harmonic gear reducer is widely used in industrial robots, aerospace, optics, and other high-end fields. ,e failure of harmonic gear reducer is mainly caused by the damage of flexible bearing and flexspline of thin-walled vulnerable components. To study the contact mechanical characteristics of flexible components such as flexible bearing and flexspline in harmonic gear reducer, the contact mechanical model of flexible bearing, vibration differential equation of flexspline, and finite element model of each component in harmonic gear reducer were established. Based on the established model of harmonic gear reducer, the influence of the length of flexspline cylinder and the thickness of cylinder bottom on the stress of flexspline is discussed, respectively, and the motion characteristics of flexible bearing are studied. At the same time, the spatial distribution of the displacement of the flexspline and the axial vibration response of the flexspline are studied. ,e correctness of the model established in this paper is verified by experiments. ,e results show that the increase of cylinder length can improve the stress of flexspline in harmonic gear reducer; the wall thickness of cylinder bottommainly affects the stress at the bottom of flexspline but has little effect on the stress of gear ring and smooth cylinder. Along the axis direction of the flexspline, the radial displacement, circumferential displacement, and angular displacement increase linearly with the increase of the axial distance between the cylinder and the bottom. When the excitation frequency is high, the vibration mode of flexspline shell is mainly axial vibration. ,e research results will provide a theoretical reference for the optimal design of harmonic gear reducer and improving the service life of flexible parts.


Introduction
Gear transmission has the advantages of high transmission accuracy and compact structure, which has been widely used in the field of mechanical transmission [1,2]. Harmonic gear reducer has the advantages of large transmission ratio, light weight, small volume, etc. It has the ability to work normally in vacuum state, strong radiation, and large temperature difference environment. Due to its unique advantages, harmonic gear reducer has been widely used in robots [3] and aerospace [4]. In the transmission process of harmonic gear reducer, the flexible bearing and flexspline belong to thin-walled parts and bear variable cycle load for a long time.
Budzik et al. [5] combined the finite element method with the contact element to analyze the influence of different torque values on the stress of flexspline. Rheaume et al. [6] carried out finite element analysis on the torsional stiffness of harmonic gear transmission system and studied the meshing stiffness between circular spline and flexible gear tooth. Folega [7] used the boundary element method to calculate the strength of the teeth of the flexible gear ring in the process of harmonic gear transmission and analyzed the design characteristics of the flexible gear and the influence of different tooth shapes on the stress value of the bottom surface of the ring gear. Li et al. [8] established an optimization model of the wave generator profile and defined the objective function to minimize the maximum circumferential stress of the flexspline assembled with the wave generator of harmonic gear reducer. Li [9] respectively carried out finite element simulation analysis on cylindrical flexspline, cup-shaped flexspline, and hat-shaped flexspline, studied the mechanical characteristics of different flexible wheels, and verified the simulation results with reference to the experiment. Pacana et al. [10] established four kinds of wave generator models: double roller, four roller, cam, and disc. Under the two working conditions of no-load and load, the stress value at the tooth of flexspline of harmonic drive was calculated numerically, and the influence law of the type of wave generator on the stress distribution of flexible gear tooth was studied. Huangfu et al. [11] established the dynamic model for spalled gear pairs, and the meshing characteristics are analyzed by using the contact analysis method of the loaded tooth surface. Routh et al. [12] established the hydrodynamic lubrication equation between the conical gap of cam and flexspline assembly surface and analyzed the cone phenomenon and lubrication mechanism in harmonic drive. Mo et al. [13,14] studied the influence of flexible support stiffness and assembly error phase on load sharing characteristics in gear transmission system. Shi et al. [15] proposed a model considering the temperature-related clearance to analyze the dynamic characteristics of ceramic bearing system. Wang et al. [16] developed a computeraided manufacturing system for cosine gears based on the conjugate theory of digital surfaces and gears. Adams et al. [17] established the model of flexible bearing and studied the vibration frequency of flexible bearing. Walha et al. [18] studied the stress and strain of flexible bearing in harmonic gear reducer by using finite element software. Shah et al. [19] predicted the fatigue life of bearing under dynamic loading and obtained the influence law of the structural size change of bearing inner and outer rings on contact stress and fatigue life. Chen et al. [20] established an analytical finite element model of gear meshing stiffness, considering crack propagation path, and analyzed the influence of crack propagation path on meshing stiffness. Ostapski [21] analyzed the stress distribution of flexible bearing under symmetrical radial load and the influence of assembly error on its fatigue damage. León et al. [22] analyzed the influence of the geometric structure of flexible gear teeth on the performance of harmonic drive. Gravagno et al. [23] studied the influence of the shape of wave generator on the motion error in the process of harmonic drive and quantitatively evaluated it.
From the abovementioned references, it indicates that many models of harmonic gear reducer are established and analyzed, but the comprehensive analysis of the contact mechanical properties of the flexible bearing and the flexspline with thin-wall vulnerable components is rare. Based on the established contact mechanics model and finite element model, the motion characteristics of the flexible bearing and the dynamic response of the flexspline are analyzed. e paper is outlined as follows. After this introduction, Section 2 introduces the mathematical model of the flexible bearing and the flexspline in harmonic gear reducer. e finite element model of flexspline, circular spline, and wave generator in harmonic gear reducer is introduced in Section 3. To verify the correctness of the model, the relevant experiments are carried out in Section 4. In Section 5, the contact mechanical characteristics of flexible parts in harmonic gear reducer are analyzed. e conclusions of present study are summarized in the last section.

Mathematical Models
Harmonic gear reducer has the characteristics of simple structure and high reduction ratio. e structure of B 3 -80 harmonic gear reducer is shown in Figure 1.
In Figure 1, harmonic gear reducer is composed of three key components, namely, wave generator, flexspline, and circular spline. e harmonic gear reducer transmits power or motion through the deformation wave generated by the controllable deformation flexible parts.

Contact Mechanics Model of Flexible Bearing.
With the change of periodic displacement excitation, the contact state between the ball and the inner and outer raceways changes, as well as the contact stiffness and contact force. To analyze the contact mechanical characteristics of the flexible bearing, the contact stiffness of the flexible bearing shown in Figure 2 is calculated. e structural parameters of flexible bearing are shown in Table 1.
According to Hertz contact theory, the ball is defined as the contact I, and the inner and outer rings of the bearing are the contact body II. e plane parallel to the radial plane of the bearing and passing through the normal of the contact surface between the ball and the raceway is defined as the main plane 1, and the plane parallel to the axial plane of the bearing and passing through the geometric center of the ball is the main plane 2. e principal curvature of the ball in contact with the inner ring of the bearing can be expressed as follows: where d w is the ball diameter, d i is the inner diameter of the inner ring, r i is the curvature radius of the inner ring groove, and ρ ij is the curvature of the contact surface of each object. e principal curvature of the ball in contact with the outer ring of the bearing can be expressed as follows: where d o is the outer diameter of the outer ring and r o is the curvature radius of the outer ring groove. e curvature sum of contact pairs can be expressed as e curvature difference of contact pairs can be expressed as follows: 2 Shock and Vibration e equivalent elastic modulus is defined as follows: where E 1 and E 2 are the elastic modulus of ball and bearing ring material, respectively, and μ1 and μ 2 are Poisson's ratio of ball and bearing ring material, respectively.
According to Hertz contact theory, the contact elastic deformation between ball and bearing raceway obeys Hooke's theorem. Compared with the curvature radius of contact surface, the size of contact area is very small. e expressions of contact ellipse area size and contact deformation can be expressed as follows: δ � 4.5Γ 3 ρ π 2 k 2 E * 2 ε where k is the contact ellipse parameter, Γ is the first kind of total elliptic integral, ε is the second type of total elliptic integral, and Q is the load value. By using the least square method, the calculation formulas of k, Γ, and ε can be expressed as follows [24]: where e contact stiffness between the ball and the inner and outer raceways of the bearing can be expressed as follows: When the load is Q � 1500 N, it can be calculated from equations (6) to (9) By substituting the geometric parameters of inner and outer race raceways of flexible bearing into equation (12), the curvature radii of raceways of inner and outer race at different polar angles can be obtained as follows: , where R i is the curvature radius at the inner ring raceway and r o is the curvature radius at the outer ring raceway.

Differential Equation of Flexspline
Vibration. e dynamic characteristics of flexspline have great influence on the whole transmission system. e cylindrical cup-shaped flexspline can be regarded as a rotating structure composed of a smooth cylindrical shell and a ring gear, and the tooth size of the ring gear part is very small, which can be transformed into a cylindrical shell by the principle of inertia equivalence [25]. A microelement is intercepted on the curved surface of flexspline shell of harmonic gear reducer, as shown in Figure 3.
In Figure 3, the x direction of the coordinate system is the diameter direction of the flexspline shell, the y direction of the coordinate system is the circular tangent direction of the flexspline shell, and the z direction of the coordinate system is the axis direction of the flexspline shell, N z , N θ , N θz , and N zθ are membrane forces per unit length, Q θ and Q z are transverse shear forces per unit length, M z and M θ are bending moments per unit length, M zθ and M θz are torque per unit length, q x and q θ are external load components per unit area along x, y, and z directions, and g 1 , g 2 , and g 3 are inertial force components per unit area along x, y, and z directions.
Under the \joint action of film force, transverse shear force, bending moment, torque, external load component, and inertial force component, the force of microelement of harmonic gear reducer flexor shell should be kept in balance.
rough the formula derivation of the mathematical model, the force equilibrium condition of microelement segment can be obtained as follows: where R is the radius of curvature of the neutral surface. e inertia force components of unit area in x, y, and z directions at any point on the curved surface of flexspline shell of harmonic gear reducer meet the following formula: where u, v, and w are the axial, tangential, and normal dislocations of any point on the neutral surface of the microelement, ρ is the density of 30 CrMnSi of alloy structural steel, and δ is the thickness of the smooth cylinder wall of the flexible wheel. Substituting equation (15) into equation (14), we can get the following equation: e stress and displacement of the cylindrical cupshaped flexspline shell satisfy the following relations: where E and µ are the elastic modulus and Poisson's ratio of alloy structural steel 30 CrMnSi of flexspline material, respectively.
e differential equations of shell-free vibration shown in equation (18) can be established by simultaneous equations (16) and (17).

Shock and Vibration
Under the free vertical simply supported boundary conditions at both ends of the flexible shell with length L, the corresponding modes of vibration of the three displacement components u, v, and w can be expressed as follows [26]: where m and n represent the axial half wave number and circumferential full wave number of the corresponding vibration modes of the three displacement components of the shell, respectively. By combining equations (18) and (19), the homogeneous linear equations of the displacement component corresponding to the free vibration of the flexspline shell can be written as e elements a ij (i, j � 1, 2, 3) in coefficient matrix A are as follows:

Shock and Vibration
where Ω 2 is the dimensionless frequency coefficient, and its relationship with circular frequency ω can be expressed as follows: e vector x can be expressed as follows: e column vector b can be expressed as follows: To obtain the nonzero solution of equation (20), the determinant of coefficient matrix A must be zero. After sorting out, the unary cubic equation of dimensionless frequency coefficient can be written as e coefficients a, b, and c of the equation can be expressed as follows: where the a i (i � 1, 2, 3), b i (i � 1, 2, 3, 4), and c i (i � 1, 2, 3, 4, 5) defined in equation (27) are expressed in "Appendix A." According to equation (19), the corresponding displacement modes of the three roots of equation (26) can be obtained as follows: where the B, C, and D defined in equation (28) are expressed in "Appendix B."

Finite Element Models
Firstly, the complete structure of flexspline, circular spline, and wave generator is designed, and then the finite element model is established according to their structural dimensions.

Flexspline Model.
e inner diameter of flexspline is 80 mm, the nominal transmission ratio is 85, the modulus is 0.5 mm, the allowable torque of low-speed shaft is 90 N·m, the allowable radial load of low-speed shaft is 3000 N, and the weight is 9.2 kg. e involute tooth profile with stable motion, uniform force, and mature processing technology is selected. e pressure angle α � 20°, addendum height coefficient h * a � 1.0, tip clearance coefficient c * � 0.35, and radial deformation coefficient w * � 1.0. e diameter of index circle of flexspline is 84 mm, the diameter of root circle is 83 mm, and the diameter of addendum circle is 84.6 mm. e length of flexspline cylinder is 70 mm, the width of gear ring is 12 mm, the wall thickness of gear ring is 1 mm, the wall thickness of smooth cylinder is 0.8 mm, and the wall thickness of cylinder bottom is 0.8 mm. e finite element model of flexspline as shown in Figure 4 is established by bottom-up method.
In Figure 4, 30 CrMnSi alloy structural steel is selected as the flexspline material. Its elastic modulus is 204 GPa, Poisson's ratio is 0.29, and density is 6691 kg/m 3 .

Circular Spline Model.
In the process of the harmonic gear reducer transmission, the flexspline and the circular spline tooth mesh to complete the transmission of the motion and the force. erefore, the modulus, pressure angle, addendum coefficient, and clearance coefficient of the circular spline are the same as those of the flexspline. e diameter of index circle of circular spline is 85 mm, the diameter of root circle is 86.7 mm, and the diameter of addendum circle is 84.35 mm. e wall thickness of gear ring is 14 mm, and the thickness of annular circular spline is 7 mm. e finite element model of circular spline established by bottom-up method is shown in Figure 5.
In Figure 5, the material of circular spline is 45 steel with Poisson's ratio of 0.3, elastic modulus of 210 GPa, and density of 7850 kg/m 3 .

Wave Generator Model.
e wave generator is composed of cam and flexible bearing. It is the component that produces deformation wave in the transmission process of harmonic gear reducer. e structure of mechanical cosine cam wave generator is shown in Figure 6. e polar coordinate equation of cam profile is as follows: where θ is the polar angle and ρ(θ) is the polar diameter. According to the cam polar coordinate equation, the plane profile is generated, and then the cosine cam 3D model is formed by stretching. Finally, the model is divided into left and right half cams along the short axis direction, and the two half cams are meshed by mapping mesh generation method. rough calculation, the maximum polar diameter of cosine cam is 30.5 mm, and the minimum polar diameter is 29.5 mm. e finite element model of cosine cam is shown in Figure 7.
In Figure 7, the material of cam is 45 steel with Poisson's ratio of 0.3, elastic modulus of 210 GPa, and density of 7850 kg/m 3 .
Flexible bearing is another part of the wave generator. Like ordinary bearing, it is composed of inner and outer rings, cage, and ball. However, the thickness of inner and outer rings is much thinner than that of ordinary bearing, which belongs to thin-walled bearing. e finite element model of flexible bearing is shown in Figure 8.

Shock and Vibration
In Figure 8, ZGCr15 is selected as the material of the inner and outer rings of the flexible bearing, with Poisson's ratio 0.3, the elastic modulus 207 GPa, and the density 7800 kg/m 3 . Nylon 66 was selected as the cage material, with Poisson's ratio of 0.28, elastic modulus of 8.3 GPa, and density of 1150 kg/m 3 . GCr1 is chosen as the material of ball, whose Poisson's ratio is 0.29, elastic modulus is 217 GPa, and density is 7850 kg/m 3 .

Experimental Verification
To verify the correctness of the finite element model of harmonic gear reducer established in this paper, the strain gauge tester is used as the data collector to study the deformation characteristics of flexspline of harmonic gear reducer. In the experiment, the motor, the input shaft, the coupling, the connecting shaft, the harmonic gear reducer, the output shaft, the deep groove ball bearing, the washer, and the nut are connected in turn. BFH1K-1.5EB resistance strain gauge sensor, XD-37GB555 DC motor, and B 3 -80 harmonic gear reducer are selected. e motor speed was set to 6 r/min. e strain is measured by resistance strain gauge, and the signal is extracted by strain gauge tester. e resistance strain gauge is pasted in the middle of the circumferential outer surface of the flexspline of the harmonic gear reducer. e experimental platform of flexspline deformation characteristics in harmonic gear reducer is shown in Figure 9.
In Figure 9, the load is generated by washer, nut, and TLA type preset torque wrench. e nut is pretightened by the torque wrench, and the load is generated by the friction between the circular surface of the washer and the side of the bearing bracket. e load can be adjusted by changing the value of the preset torque wrench; when the washer and nut are removed, the output shaft idles to simulate the no-load condition.     Shock and Vibration e strain value corresponding to the flexspline rotation in harmonic gear reducer is extracted, and the experimental data are compared with the theoretical value obtained by finite element contact mechanics model simulation. e comparison results are shown in Figure 10.
In Figure 10, the simulation results of the finite element contact mechanics model established in this paper are consistent with the experimental results, and the relative error is small. erefore, the correctness of the finite element contact mechanics model of harmonic gear reducer established in this paper is verified, and the model can be used to predict the deformation of flexspline of harmonic gear reducer.

Sensitivity Analysis of Flexspline Stress Characteristics.
To analyze the stress sensitivity of flexspline in harmonic gear reducer, the influence of cylinder length and cylinder bottom wall thickness on flexspline stress is discussed.

Influence of Cylinder Length on Flexspline Stress.
e length L of the cylinder determines the geometric structure of the flexspline in the axial direction. It is one of the key structural parameters of the cylindrical cup flexspline and has a significant effect on the stress characteristics of the flexspline. erefore, the stress characteristics of flexspline are analyzed in the range of L � 50∼85 mm. e stress nephogram of flexspline under different cylinder lengths is shown in Figure 11.
In Figure 11, when the length L of the flexspline increases from 50 mm to 85 mm, the maximum equivalent stress decreases by 34.3%. erefore, the increase of cylinder length L can effectively improve the stress condition of flexspline.
e flexspline is divided into three parts: gear ring, cylinder, and bottom of flexspline, and the stress influence curve of flexspline tube length as shown in Figure 12 can be obtained.
In Figure 12, compared with the bottom part of the ring gear, the stress of the ring gear is greater. In the range of L � 50∼70 mm, the maximum equivalent stress of the ring, cylinder, and bottom of the flexspline decreases obviously, from 683 MPa to 514 MPa, 682 MPa to 365 MPa, and 142 MPa to 98.1 MPa. In the range of L � 70∼85 mm, the decrease trend of the maximum equivalent stress tends to be gentle, from 514 MPa to 488 MPa, 406 MPa to 303 MPa, and 98.1 MPa to 79.2 MPa. It shows that the stress concentration of flexspline decreases with the increase of cylinder length.

Influence of Bottom Wall ickness on Flexspline Stress.
To analyze the influence of the bottom wall thicknessδ″ on the stress of the flexspline, the stress characteristics of the flexspline were analyzed in the range of wall thicknessδ ″ � 0.80∼1.15 mm. e stress nephogram of flexspline under different wall thicknesses of cylinder bottom is shown in Figure 13.
In Figure 13, with the increase of the wall thickness at the bottom of the cylinder, the maximum equivalent stress of the flexspline is about 506 MPa. It shows that the change of the wall thickness at the bottom of the cylinder has little effect on the stress of the flexspline. e flexspline is divided into three parts: gear ring, cylinder, and flexspline bottom, and the stress curve of

Influence of Flexible Ring on Bearing Contact Stiffness.
To explore the influence of flexible ring on contact stiffness of bearing in harmonic gear reducer, the influence of ring stiffness and flexibility characteristics on contact stiffness was analyzed with ordinary rigid bearing as the comparison object. e difference of contact rigidity between ball and inner and outer rings of rigid bearing and flexible bearing at different polar angles is shown in Figure 15.
In Figure 15, in the vicinity of 0°, 90°, 180°, and 270°, the contact rigidity difference between the ball bearing of rigid ring and flexible ring and the channel is small. e difference values of the inner ring contact stiffness △K i were 2.39 × 10 8 N/m, 2.16 × 10 8 N/m, 2.39 × 10 8 N/m, and 2.16 × 10 8 N/m, respectively.
e difference values of the outer ring contact rigidity △K o were 2.04 × 10 8 N/m, 1.73 × 10 8 N/m, 2.04 × 10 8 N/m, and 1.73 × 10 8 N/m. Because the polar angles of 0°, 90°, 180°, and 270°are the positions of the major and minor axes of the cosine cam wave generator, the value of the contact stiffness between the ball and the inner and outer rings is closer.

Analysis of Motion Characteristics of Flexible Bearing.
e ball at the spherical coordinates (35, 0) was defined as ball no. 1, and then the node N52521 at the spherical center and the node N53753 at the contact area between the sphere and the raceway were extracted, respectively. e rolling body center and surface node velocity curve of nodes N52521 and N53753 are shown in Figure 16.
In Figure 16, the velocity curve of node N52521 at the center of the ball presents the law of sines and cosines along the direction of x axis and y axis. e phenomenon of fluctuation appeared within a small range of the velocity curve. e fluctuation was caused by the gap of 0.1 mm between the ball and the pocket hole of the cage, and there is a slight collision between the ball and the cage during the ball movement. e joint velocity of the node N53753 on the surface of the ball changes within the range of 0∼13.57 m/s, and the peak and trough of the joint velocity alternate after 0.006 s. e corresponding time of the peak of the velocity represents the contact between the node and the raceway of the bearing inner ring, and the corresponding time of the trough of the velocity represents the contact between the node and the raceway of the bearing outer ring. e motion period of the node is 0.003 s, that is, the rotation period of the ball is 0.003 s.

Analysis of the Spatial Distribution State of the Flexspline
Displacement. At three different cross-sections of the flexspline barrel, the spatial distribution of the radial displacement, circumferential displacement, axial displacement, and angular displacement of the flexspline along with the change of polar angle and the axial distance between the flexspline barrel and the barrel bottom are shown in Figure 17.
In Figure 17, along the direction of the pole angle θ, the radial displacement, circumferential displacement, axial displacement, and angular displacement of the flexspline conform to the law of sines and cosines. Along the direction of the central axis of the flexible cylinder, with the increase of the axial distance Z from the cylinder bottom, the values of radial displacement, circumferential displacement, and angular displacement gradually increase, showing a trend of linear increase. Compared with radial and circumferential displacements, the value of axial displacement is one order of magnitude smaller, and its value does not change with the size of the axial distance Z.

Analysis of Axial Vibration Response of Flexspline.
e circumferential full wave number n is kept unchanged. When the axial half wave number m is taken as different values, the axial displacement vibration response of the flexspline is obtained, and the influence of the axial half wave number m on the axial vibration response of the flexspline is analyzed.
e contour map of the axial displacement of flexspline is shown in Figure 18.
It can be seen from Figures 18(a) to 18(c) that when m � 1, along the axis of the cylinder, there is a peak value and a valley value in the contour map. Along the circumference of the cylinder, there are two peak values, one valley value or one peak value, and two valley values. When i takes 1, 2, and 3, respectively, the proportion relationship of axial displacement amplitude is 1 : −0.416 : 0.017. It can be seen from  Figures 18(d) to 18(f ) that when m � 2, there are two peaks, one valley or one peak, two valleys along the axis, and circumference of the cylinder. When i is 1, 2, and 3, respectively, the proportional relationship of axial displacement amplitude is 1 : −0.101 : −0.006. erefore, the increase of axial half wave number m will enhance the axial displacement vibration response of flexspline.

Conclusions
In this paper, the contact mechanics model of flexible parts in harmonic gear reducer is established. e influence of flexspline cylinder length and flexspline cylinder bottom wall thickness on flexspline stress is analyzed. e motion characteristics of flexible bearing in contact process are studied, and the vibration response of flexspline is analyzed. e conclusions can be summarized as follows: (1) e stress concentration of flexspline decreases with the increase of cylinder length, and the increase of cylinder length can effectively improve the stress condition of flexspline of harmonic gear reducer.

Data Availability
e data used to support the findings of this study are included within the article. e processed data are available from the corresponding author upon request.

Conflicts of Interest
e authors declare that they have no conflicts of interest.