Statistical Analysis for Transmission Error of Gear System with Mechanical and Thermal Deformation Uncertainties

: We establish a robust algorithm to analyze the inﬂuence of system uncertainties on the transmission error of a spur gear pair under 2D simpliﬁcation. The algorithm provides a way of generating smooth cutter proﬁles with machining uncertainties and measuring the thermal deformation through the uncertainties in material properties. Then, it produces realizations of gear tooth proﬁles based on the analytical method for accuracy and computational efﬁciency. Numerical investigations show the statistical analysis on the tooth contact analysis by comparing steel and plastic gears. It is worthwhile remarking that the plastic gear is susceptible to the geometric error caused by thermal deformation. Moreover, although the impact of thermal deformation on steel gear may seem slim, it can have a noticeable inﬂuence when it exists with mechanical uncertainties together.


Introduction
The gear system is a complex elastic mechanical system that is widely used for power transmission in industries. Various gear failures and environmental factors can cause the increase in vibration and noise as described, for example, in [1]. As they affect the performance and durability of the entire system, it is important to precisely measure the noise and vibration produced by gears; see [2], for example, which measures vibrations and noise of a chipper with an internal combustion drive to identify the scale and location of origin of vibration. The transmission error (TE) is a well-known representative criterion for such measurement. Thus, accurate prediction of the TE is indispensable for precise and reliable transmission system design. However, it is not easy because a typical gear system has numerous uncertainties, such as lubrication, material nonlinearity, contact, geometrical tolerance, and nonlinear input loading condition [3,4].
While many studies, e.g., [5][6][7][8], on geometric measurement have succeeded in enhancing the accuracy of geometric specification, gear manufacturing, or wear evaluation, it is still inevitable to admit the existence of geometric errors in gears induced by diverse uncertain factors. Uncertainty quantification (UQ) measures the possible uncertainties and estimates their relations with one another or the effects on the entire system. UQ involves specifying appropriate uncertain factors, quantifying the uncertainties, estimating their effects on the quantity of interest (QoI), and evaluating the sensitivity of the QoI to the uncertain factors. Ultimately, it facilitates adjustment of sensitive uncertain parameters and sophistication of the model design.
Many studies have applied UQ to gear systems in order to examine the influence of uncertain parameters on the TE. The finite element method (FEM) serves as a powerful tool for calculating the tooth deflection and mesh stiffness of the gear system, which has various uncertainties. For example, Li [9] quantitatively investigated the effects of machining errors, pitch errors, and tooth modifications on the load sharing ratio and TE by using specially developed FEM software. Lin and He [10] comprehensively generated a gear model with machining and pitch errors and analyzed the TE using FEM software. Specifically, the helix error and the error in the cutter profile were determined using the sine function. With all the errors considered, the three-dimensional geometry of the tooth surface was generated via GEMTE software. Then, by varying the degree of the uncertainties, the static TE was investigated to analyze the interactions between the gear uncertainties and their effects on the TE.
For reliable analysis of the tooth surface contact, a highly refined FEM mesh is required due to the nonlinearity of contact problems. Numerous studies have developed model reduction algorithms to improve computational efficiency while retaining accuracy (e.g., [11]). Still, for uncertainty quantification, non-intrusive methods such as Monte Carlo simulations require a considerable number of simulations. Thus, despite its advantages, FEM could be inappropriate for the gear design under various uncertainties.
By contrast, an analytical method can rapidly calculate the TE via numerical kinematic simulation [12,13]. This method is computationally efficient and provides an analytical function of the gear mesh stiffness consisting of the design parameters of the gear [14,15]. For example, Wang [16] generated a gear model with tooth modification and analyzed the TE using the analytical method. The modified tooth profile was obtained from the modified cutter profile and used for modeling the tooth surface.
However, the uncertainties in thermal deformation have not been studied sufficiently in the case of both FEM and analytical approaches. Owing to various scientific and industrial advancements, gears can currently operate under heavy loads and high speeds. Under such conditions, the frictional heat on the tooth surface increases the temperature near the contact point, which leads to thermal deformation that can deform the gear tooth and affect the TE of the system significantly. In particular, recent years have witnessed the growing use of polymer-based gears [17][18][19][20][21]. As typical polymer composites have greater thermal sensitivity than steel materials, it is important to examine their thermal deformation in the context of gear design.
In this study, we establish a comprehensive algorithm for UQ of a spur gear pair by considering both mechanical and thermal uncertainties under 2D simplification. We assume that the spur gear has a constant profile in the longitudinal direction and focus on the contact analysis of tooth surface with uncertainties using 2D simplification. The main contributions of this study can be summarized as follows. First, we present a method for realizing smooth cutter profiles with machining uncertainties. It generates a gear tooth profile machined with an arbitrary cutter shape, which is analyzed using the analytical method. Second, we estimate the thermal deformation on the basis of the uncertainties in the gear material properties. To calculate the TE, the displacement due to thermal deformation estimated using the FEM is considered for tooth contact analysis. The influence of the uncertainties in the gear material properties on the TE are statistically analyzed using Monte Carlo simulation. Moreover, we compare the results for steel and plastic spur gear pairs. Accordingly, we clarify that in plastic systems, the thermal deformation has a dominant effect on the machining and pitch errors, whereas the opposite effects are observed in the case of steel. Finally, we analyze the combined influence of three types of uncertainties, i.e., machining errors, pitch errors, and thermal deformation, on the TE. It is worth noting that although thermal deformation alone does not change the TE of steel significantly, it can produce a noticeable difference when it exists with the other type of uncertainties.
The remainder of this paper is organized as follows. Section 2 reviews the analytical method for gear design and TE analysis. Section 3 introduces the proposed method for de-signing a spur gear pair with mechanical uncertainties and thermal deformation. Section 4 describes the TE estimation via loaded tooth contact analysis and presents an algorithm for Monte Carlo simulation. Section 5 summarizes and analyzes the UQ results and discusses the effects of each uncertainty on the TE. Finally, Section 6 concludes the paper.

Rack Cutters and Gear Surfaces
Among the various methods for manufacturing gears, we assume that a rack cutter is used to manufacture the gears and analyze them using an analytical approach [22][23][24]. Figure 1 shows the profile of a basic rack cutter and the relationship between the coordinate systems C c [X c , Y c , Z c ] and C g [X g , Y g , Z g ] of the rack cutter and the spur gear, respectively. The rack design parameters h ap , h f p , p f p are addendum, dedendum and fillet radius, respectively. The coordinate systems C 1 [X 1 , Y 1 , Z 1 ] and C 2 [X 2 , Y 2 , Z 2 ] are the reference coordinate systems of the manufacturing process. Here, u c is the x-coordinate of the position vector r c , d c is the distance between the rack cutter and the gear, and α is the pressure angle of the rack cutter profile. In the coordinate system C c of the rack cutter, the position vector r c of the rack cutter surface can be represented as where θ c denotes the z-coordinate of r c . When the rack cutter translates along O c x c by −φ c , the gear blank rotates about the axis O g z g by φ g at the same time. The relationship between these movements is given by (2) Figure 1. Relationship between the coordinate systems C c and C g of the rack cutter and the spur gear, respectively. Figure 2 shows the trace of the rack cutter profile in the coordinate system C g of the spur gear. The position vector r g of the gear in C g [X g , Y g , Z g ] can be written in terms of the family of rack cutter surfaces r c as The position r g of the gear is described in a three-dimensional volume about u c , φ c , and θ c . Meanwhile, the gear surface is an envelope of the family of rack cutter surfaces, and its profile is an involute curve. The equation of meshing satisfied by the envelope of the family of rack cutter surfaces is given by where r * g is a vector obtained by eliminating the last element of r g . It should be noted that Equation (4) can be rewritten as an explicit function of u c in terms of φ c : which allows us to represent the position r g of the gear in a two-dimensional space in terms of φ c and θ c , which is the manufactured gear surface that we can analyze.

Tooth Contact Analysis (TCA)
We generate the manufactured gear surface, which can be represented in 3D. However, we assume that spur gear has a constant profile in the longitudinal direction to simplify the tooth contact. We focus on the contact analysis of spur gear pair under 2D simplification. A contact between two rotating gears must be defined to calculate the transmission error. To analyze the contact between the gears, we consider the relationship between the coordinate systems of the spur gear pair and their local coordinates C S1 [X S1 , Y S1 , Z S1 ] and C S2 [X S2 , Y S2 , Z S2 ]. Figure 3 shows the local coordinates of the spur gear pair, where ψ is the rotation angle of gear and d g is the distance between the two gears. The vectors r gi (i = 1, 2) of the (non-rotating) gear teeth surfaces in each local coordinate system can be represented in terms of the parameters satisfying Equation (5). When the gears rotate, the position vectors S i (i = 1, 2) and the normal vectors N i (i = 1, 2) of the gear surfaces in C si (i = 1, 2) are given by and respectively. Here, S * i is a vector obtained by eliminating the last element of S i . When the two teeth surfaces S 1 and S 2 are in contact, two kinematic conditions must be satisfied at the contact point: The first condition implies that the position vectors are equal, and the second condition implies that the normal vectors are parallel to each other. To simplify Equation (8), we assume that there is no variation on gear surface in the z-axis, the longitudinal direction of gear. Note that the parameter θ c , the z-coordinate of the position vector r c , can be eliminated by following the 2D simplification. Accordingly, three equations P i (i = 1, 2, 3) can be derived from Equation (8): The Newton-Raphson method [25,26] is used to solve Equation (9) and obtain the four unknown parameters (u 1 , u 2 , ψ 1 , ψ 2 ). In the i-th iteration, the method yields a solution X i satisfying which converges to X i −→ X as i → ∞. If one of the unknowns, say ψ 1 , is determined, we can find the solution X for the parameters (u 1 , u 2 , ψ 2 ) from the Newton-Raphson iterations of Equation (10). The TE can be calculated as where n 1 and n 2 are the number of teeth of the driving and driven gears, respectively, and r b2 is the base circle radius of the driven gear.

Loaded Tooth Contact Analysis (LTCA)
It is difficult to predict the performance of a spur gear pair using the TE (see Equation (11)), because the TE under the unloaded condition is affected by only the tooth geometry; the influence of the elastic deflection is ignored. Loaded tooth contact analysis is used to consider the elastic deflection of the tooth, and the transmission error ∆ψ L under the loaded condition is expressed as where ∆ψ U is the transmission error under the unloaded condition and δ E is the elastic deflection of the tooth. The elastic deflection δ E can be calculated using the mesh stiffness of the gear system, which is determined by the bending deflection, fillet-foundation deflection, and contact deflection [14]. Details of the mesh stiffness are described in Appendix A.

Mathematical Modeling of Gears with Uncertainties
Owing to the uncertainties in the manufacturing process, there are differences between the theoretical and manufactured models. In particular, complex errors exist in the tooth surface of the manufactured gear compared to the theoretical tooth surface [10]. As shown in Figure 4, such complex errors include multiple errors in the gear design parameters, each of which is related to mechanical uncertainties and thermal deformation. However, since we consider the tooth contact of spur gear pair under 2D simplification, helix deviation, longitudinal crowning, and assembly misalignments are not considered in this paper. In this section, we introduce the different types of uncertainties in a spur gear pair and describe how we modify the geometry of the gear to consider these uncertainties under 2D simplification.

Mechanical Uncertainties
Mechanical uncertainties imply all the random occurrences in the manufacturing process. The pitch error (PE) and machining tool error (ME) are considered as typical mechanical uncertainties. These errors are considered as factors to estimate the gear precision grade in ISO 1328-1 [27].
Pitch errors occur because of the difference in the interval between the rack cutter teeth. Because these are the only factors that affect the interval of the gear teeth, the uncertainties e p for the pitch errors can be measured by the movement of the rack cutter: for some E p > 0. We consider the machining error in the involute curve of gear to focus on the tooth contact in the gear flank surface because tooth contact usually occurs in the gear flank surface. Meanwhile, machining tool errors arise from differences between the theoretical and real profiles. Although the profile of the rack cutter surface is considered straight theoretically, it is an uneven curve in practice. To consider the uncertainties, n random points {(x 1 , y 1 ) . . . (x n , y n )} are generated within a given finite interval based on the following theoretical profile: where [L a , L b ] is the range of O p X p and denotes the size of the machining tool error on O p Y p . We consider 2 as the randomness range (see Figure 5). To obtain a smooth and continuous curve for the rack cutter profile, we interpolate the random points using Lagrange polynomials, which gives f (x) as follows  Figure 5 shows the random points on the rack cutter profile and their interpolation curve f . Here, are the coordinate systems of the rack cutter profile and the rack cutter, respectively. Further, O p X p coincides with the theoretical rack cutter profile. By considering an additional coordinate system C p besides the original coordinate C c , we can generate discrete random points and adjust their range (i.e., 2 ) in the normal direction of the cutter surface. The rack cutter profile f , where the machining tool uncertainties are considered, is represented by r p in C p . The rack cutter profile in the coordinate system C c is denoted by r c , and it can be written in terms of r p as follows: where u p and θ p are the x-and z-coordinates of the position vector r p , respectively.

Thermal Deformation Uncertainties
In this study, we assume that either mechanical uncertainties or thermal deformation can influence geometric errors in gears. While mechanical uncertainties have been considered as an important factor affecting the performance of gear pairs (or systems), thermal deformation has been relatively less studied and ignored in most previous works [9,10]. Relative significance of the two types of uncertainties depends on the material properties. For example, in steel gears, thermal deformation barely changes geometric errors, and the impact of mechanical uncertainties is dominant. On the other hand, thermal deformation can be more influential than mechanical uncertainties. Figure 6 shows the geometric error induced from the machining error (black) and from thermal deformation (red), respectively, when a plastic gear operates at 200 RPM under 20 Nm. We use the machining error described in Section 3.1 when = 4. The implementation of thermal deformation is described in the remainder of this section. Here, we consider deterministic thermal deformation, i.e., η = 0; see Equation (17) below. Figure 6 exhibits that thermal deformation can lead to higher geometric error than machining error in the case of a plastic gear, which is unlikely to happen in steel gears. . Geometric error induced from the machining error (black) and from thermal deformation (red), respectively, when a plastic gear operates at 200 RPM under 20 Nm. The parameters for machining error and thermal uncertainties are set as = 4 and η = 0, respectively. It shows that thermal deformation can lead to higher geometric error than machining error in the case of a plastic gear, which is unlikely to happen in steel gears.
Friction between the gear meshing surfaces makes a gear pair operate and continuously generates heat on the tooth surface. The frictional heat increases the temperature of the gear pair, resulting in thermal deformation of the gear tooth. To estimate the thermal deformation, we measure the change in displacement of a gear tooth at the temporal steady-state when heat transfer is complete. By performing FE analysis at the steady-state, we estimate the transition of displacement from the initial state. The displacement change depends on the thermal expansion coefficient and steady-state temperature. To clarify, we consider the frictional heat as the only heat source and assume that the operating conditions such as torque, speed, lubrication, and ambient temperature are given a priori following the design process of the entire system. The heat convection parameters used in thermal analysis are constant values that are determined according to the lubrication condition and operation environment. Thus, when the frictional heat is constant, material properties are the only factors that determine the steady-state temperature and the displacement change.
In this study, we assume randomnesses in five essential material properties, i.e., Young's modulus, thermal expansion coefficient, density, conductivity, and specific heat capacity, to investigate the uncertain thermal deformation. All material properties are assumed to follow uniform distributions within the ranges thermal expansion coe f f icient, µ * 3 : density, µ * 4 : conductivity, µ * 5 : speci f ic heat capacity).
Here, µ i represents each material property with uncertainty, and η denotes the size of the uncertainties of the material properties.
Frictional heat is related to the contact stress and relative sliding speed, which depends on where the contact occurs. Each position of the gear tooth has a different temperature when the gear pair reaches a steady state. Thus, the gear has a different shape compared to the involute curve after thermal deformation. Rolling friction, friction caused by deformation, and sliding friction are the three main aspects of the frictional heat generated between the gear meshing surfaces [28]. Rolling friction and friction caused by deformation have a negligible effect on heat generation. Therefore, only sliding friction is considered to predict the frictional heat. The sliding friction is determined by the contact stress, relative sliding speed, and friction coefficient between the gear meshing surfaces. The frictional heat at each contact point is used for thermal analysis based on the FEM [29].
The thermal deformation in the steady state is measured by the nodal displacements of the tooth surface in the coordinate system C c . Here, U x and U y denote the xand ydirectional nodal displacements, respectively. The nodal displacements U x and U y can be expressed in terms of u p , the x-coordinate in C p , using the relationship between u c and u p in Equation (16).
To construct continuous displacement curves of the thermal deformation, we interpolate the discrete nodal displacements U x and U y . Note that the profile of the gear surface consists of (n + 1) nodes, whose nodal positions are denoted by u pi in the coordinate system C p , i = 0, 1, · · · , n. We denote the nodal displacements U x and U y at u p = u pi by U xi and U yi , respectively. Figure 7 (black diamonds) shows the discrete nodal displacements U xi and U yi , i = 0, 1, · · · , n. By using the monomial polynomials u j p , j = 0, 1, · · · , n, the displacement curves of thermal deformation in the xand y-directions, denoted by H x and H y , respectively, can be expressed as where  Figure 7 (red circles) shows the interpolation curves for the thermal deformations, i.e., H x (u p ) and H y (u p ).
To represent the gear geometry with thermal deformation, the continuous curves H x and H y should be considered in the mathematical method for designing cylindrical gears. The position vector r T , which expresses the gear geometry with thermal deformation, is represented as where M T is a transformation matrix used to consider the effect of thermal deformation.

Statistical Analysis for Gear System with Uncertainties
For statistical analysis, we generated spur gear pairs with mechanical uncertainties and thermal deformation and calculated the TE of a spur gear pair using LTCA. The deterministic LTCA is described in Section 2.2. The transmission error ∆ψ L under the loaded condition is the sum of the transmission errors ∆ψ U under the unloaded condition and elastic deflection δ E , where both ∆ψ U and δ E are affected by uncertainties. Further, ∆ψ U is calculated using the TCA of the spur gear pair, which has a profile with uncertainties, and δ E is calculated using the gear profile and material properties that have uncertainties.
To measure and analyze the responses of a spur gear pair with comprehensive uncertainties, we performed Monte Carlo (MC) simulation [30][31][32]. For this purpose, we varied the main uncertain parameters e a (Equation (13)), (Equation (14)), and η (Equation (17)), and generated each MC sample by following the Latin hypercube sampling (LHS) method [33,34]. The detailed algorithm for the MC simulation for generating spur gear pairs is described in Section 4.2.

LTCA for Gear System with Uncertainties
To calculate the transmission error ∆ψ U under the unloaded condition, we considered the gear profile with uncertainties. The equation of meshing satisfied by the envelope of the family of rack cutter surfaces cannot be rearranged as an explicit function because of the uncertainties of the rack cutter surface. Therefore, we need to consider the equation of meshing as an additional condition along with the kinematic conditions (8) satisfied at the contact points, which results in three kinematic conditions: where Here, S ei , N ei , r pi , and r Ti (i = 1, 2) are the position vector, normal vector of the gear surface with uncertainties in C si , position vector of the gear surface in C p , and position vector of the gear surface with thermal deformation in C g , respectively. The subscripts i = 1, 2 denote driving and driven gears, respectively. Further, S * e1 , S * e2 , and r * p are vectors obtained by eliminating the last elements of S e1 , S e2 , and r p , respectively. To simplify Equation (21), we assume that there is no variation on gear surface in the z-axis, the longitudianl direction of gear. Note that parameter θ p , the z-coordinate of the position vector r p , can be eliminated as following the 2D simplification. From Equation (21), five equations P ei (i = 1, 2, 3, 4, 5) are derived: Using the Newton-Raphson method, we can obtain the six unknown parameters (u p1 , φ 1 , ψ 1 , u p2 , φ 2 , and ψ 2 ) that satisfy Equation (23).
The elastic deflection δ E is affected by the gear profile, material properties, and applied load. However, to focus on gear profile and material properties related to uncertainties, we exclude the load sharing effect in this study. The mesh stiffness changes according to the contact zone, the shape of the profile, and elastic modulus. Thus, we consider the uncertainties in the gear profiles and elastic modulus when calculating the bending deflection δ e b , fillet-foundation deflection δ e f , and contact deflection δ e h , all of which together determine the elastic deflection δ E .
The gear profile with uncertainties consists of (n + 1) points, and the ith point is expressed as r i g , as shown in Figure 8. A gear tooth is regarded as a nonuniform cantilever beam to calculate the bending deflection. The position of r i g is used to measure the area of the tooth cross section S i and the area moment of inertia from the nonuniform cantilever beam. The gear profile changed by the uncertainties affects the contact zone, which can be obtained using Equation (23) of the gear system. The contact zone changes the area where the force is applied and thus determines the bending deflection. The three deflections can be written as where E e , A e i , and I e i are the elastic modulus, area of the tooth cross section S i , and area moment of inertia, respectively, each with uncertainties. The detailed notations are described in Appendix A.
The mesh stiffness K e meshing of the spur gear pair with uncertainties is expressed as where k e b1 and k e b2 are the bending stiffness of the driving and driven gears, respectively, with uncertainties. k e f 1 and k e f 2 are the fillet-foundation stiffnesses of the driving and driven gears, respectively, with uncertainties.

Monte Carlo (MC) Simulation
Monte Carlo (MC) simulation was used to evaluate the effects of the respective uncertainties on the transmission error. The uncertainties in the pitch, machining, and material properties are characterized by the parameters E p , , and η, respectively. The parameters indicate the range or magnitude of the respective uncertainties. To estimate the influence of these uncertainties on the spur gear pair, we generated gear model samples under a randomly given pair (E a , , η).
Specifically, for the pitch uncertainty, we randomly choose N values of e i p ∈ [−E p , E p ], (i = 1, 2, · · · , N), for a given uncertainty parameter E p > 0. From Equation (13), the movement of the rack cutter φ g is measured for each e i p . Second, to obtain the rack cutter profile samples in which machining uncertainty is considered, we adopt the Latin hypercube sampling (LHS) method [33,34]. It is known that LHS exhibits a faster rate of convergence than a simple random sampling method. For a given machining uncertainty parameter > 0, the LHS operates by dividing the uncertainty ranges into non-overlapping subspaces of equal size and obtains a random point from each subspace. Figure 9 shows an example of the sampling for the machining uncertainty. Here, we suppose that = 1 and divide the X, Y axes into m sub (= 8) subspaces of equal size. Out of m 2 sub subspaces in R 2 , we choose m sub subspaces such that their rows and columns are not overlapping with each other. In each subspace, a point is randomly selected. Then, as described in Section 3.1 and Figure 5, we obtain a rack cutter profile sample by interpolating these points. For a fixed > 0, we can generate N different sets of points from the LHS subspaces, which leads to N rack cutter profile samples. In the numerical tests, we use a fixed number of subspaces for the LHS sampling m sub = 8; see Algorithms 1 and 2 and Section 5. By interpolating these points, we obtain a rack cutter profile sample; see Section 3.1 and Figure 5.
Finally, to calculate the geometric error of thermal deformation, uncertainties in the material properties are considered. Five material properties (µ 1 : Young's modulus, µ 2 : thermal expansion coefficient, µ 3 : density, µ 4 : conductivity, µ 5 : specific heat capacity) generated randomly within a given range [−µ i (1 − η), µ i (1 + η)], (i = 1, 2, · · · , N) were used for the thermal analysis. We can obtain data on the geometric error generated by thermal deformation through thermal analysis using FEM. Two types of materials (steel and plastic) were considered to investigate the influence of thermal deformation on the transmission error.
In summary, from a given pair (E a , , η), we have N samples of gear models in which all three types of uncertainties are considered. The pitch and machining uncertainties are related to the N rack cutter samples and their movements. Through the rack cutters, the corresponding N gear profiles are generated, where the displacement due to thermal deformation is affected by the material property.
After obtaining the N gear model samples, we combine them with each other and generate (N(N − 1))/2 gear systems with uncertainties. Then, LTCA is conducted to calculate the transmission error under the loaded condition. Depending on which gears are combined and the related uncertainties, the spur gear pair exhibits different transmission errors. Detailed algorithms for obtaining (N(N − 1))/2 spur gear pair samples for given uncertainty parameters (E a , , η) are described in Algorithms 1 and 2.
Based on the spur gear pair samples with uncertainties, we performed statistical analysis of the corresponding transmission errors. Randomly select e i p ∈ [−E p , E p ].

Numerical Studies
The influences of the mechanical uncertainties and thermal deformation on the transmission errors of a spur gear pair were examined using the proposed method with numerical examples. The tooth contact analysis of numerical examples is under 2D simplification, and we focus on studying the influence of uncertainties on 2D spur gear pairs. We considered three numerical examples for investigating the influences of the uncertainties and transmission errors with a comprehensive error. Table 1 summarizes the design parameters and material properties of the spur gear pair used in the numerical example. There are two cases of gear pairs, i.e., steel and plastic. For lubrication, we consider its effect on the steel gear pairs only since, in general, plastic gear pairs do not use lubrication due to material characteristics. Young's modulus (MPa) 206,000 2150 Thermal expansion coefficient (10 −6 / • C) 11.5 85 Density (T/mm 3 ) 7.83 × 10 −9 1.13 × 10 −9 Conductivity (W/(m × K)) 50 0.29 Heat capacity (mJ/(T × K)) 4.58 × 10 8 1.67 × 10 9

Effect of Machining Error (ME) on Transmission Error (TE)
Thermal deformation and pitch error are not considered as uncertainties in this example. To study the influence of the machining error on the transmission error, the mean and standard deviation of the transmission error are analyzed under unloaded and loaded conditions.
The mean and standard deviation of the transmission error under the unloaded condition are shown in Figure 10. The performance of the spur gear pair under the unloaded condition is affected by the shape of the tooth profile. When = 0, the tooth profile is an ideal involute curve, and the spur gear pair has no transmission error. As increases, the fluctuation of the transmission error increases, which can be confirmed as shown in Figure 10. Figure 11 shows the mean and standard deviation of the transmission error under the loaded condition. The plastic spur gear pair has larger fluctuations than the steel spur gear pair in both the mean and the standard deviation of the transmission error because plastic has lower stiffness than steel. Furthermore, there is a difference in the tendency of the transmission error. The transmission error of the steel spur gear pair is similar to the transmission error under the unloaded condition, whereas that of the plastic spur gear pair is not. This implies that plastic is susceptible to applied force. As increases, the standard deviation of the peak-to-peak transmission error (PPTE) increases. From Figure 12, the range of the PPTE of the designed spur gear pair can be predicted. The machining error under 2D simplification can be compared with total profile deviation, which is one of the factors in ISO precision grade. It is possible to calculate the level of machining error required to satisfy the transmission error of the desired gear quality.   There is no thermal deformation or pitch error (η = 0, e p = 0).

Effect of Thermal Deformation (TD) on Transmission Error (TE)
In this example, LTCA is conducted on a spur gear pair with thermal deformation under unloaded and loaded conditions to study the transmission error. To compare the influence of thermal deformation on the transmission error between the spur gear pairs with different materials, TCA was conducted on steel and plastic, and the mechanical uncertainty was not considered.
The transmission error of the spur gear pair under the unloaded condition was generated by the geometric error caused by thermal deformation. The difference in transmission error between steel and plastic gear pairs results from the difference in material properties. From Figure 13, it is found that the plastic spur gear pair is more susceptible to thermal deformation than the steel spur gear pair. When η increases, the increase in the standard deviation of the plastic spur gear pair is larger than that of the steel spur gear pair. This difference becomes larger under the loaded condition. The standard deviation of the steel spur gear pair under the loaded condition is similar to that under the unloaded condition, whereas that of the plastic spur gear pair increases by approximately two times, as shown in Figure 14. The PPTE generated by thermal deformation is similar regardless of the size of η. As η increases, the standard deviation of the PPTE of the steel and plastic spur gear pairs tends to increase; however, it is not as large as the standard deviation of the machining error, as shown in Figure 15. The more precise the spur gear pairs with smaller machining errors, the greater is the transmission error generated by thermal deformation. When we design a spur gear pair with a material that is susceptible to thermal deformation, considering the range of the transmission error expected during operation will be useful for a safer and more precise design.

Effect of Comprehensive Error on Transmission Error (TE)
A spur gear pair with a comprehensive error consisting of a machining error, thermal deformation, and a pitch error is generated, and LTCA is conducted on plastic and steel. Figures 16-19 show the transmission errors of the spur gear pair as one of the samples. Samples of the spur gear pair have a machining error ( = 0.5), thermal deformation (η = 0.0125), and a pitch error e p . The pitch error of the driving gear is e p = −0.001 deg in Figures 16 and 17 and e p = 0.001 deg in Figures 18 and 19. The driven gear has no pitch error.
Transmission errors with a comprehensive error under unloaded and loaded conditions are calculated to be larger than the transmission errors of the three uncertainties. In the case of a steel gear under the unloaded condition, the transmission error with thermal deformation is the smallest among the three uncertainties. Although the transmission error caused by only thermal deformation is small, thermal deformation coupled with mechanical uncertainties causes a significant transmission error. When thermal deformation is coupled with other uncertainties, it has a more considerable impact on the transmission error than magnitude. From Figures 16 and 17, it is found that the effect of the interaction of thermal deformation is different according to the material properties. The effect of the interaction of each uncertainty is clearly shown in Figure 20. There is a difference between the comprehensive uncertainty (red line) and the linear sum of each uncertainty (blue line). Note that the effect of the interaction of each uncertainty leads to the difference between the geometric error of the comprehensive uncertainty and the linear sum of each uncertainty. In the case of a plastic gear under the unloaded condition, the tendency of the transmission error with the comprehensive error is similar to that of thermal deformation. From this fact, it is found that thermal deformation has the greatest effect on the transmission error among the three uncertainties. In the case of a plastic gear under the loaded condition, all transmission errors with uncertainty have similar tendencies owing to elastic deflection. The gear tooth, which has a pitch error, causes an error in the tooth contact pattern. When the gear tooth contact is late or early, the difference in the contact pattern appears as a transmission error. From comparison between Figures 16 and 17 (or Figures 18 and 19), it is found that the magnitude of pitch error affects the major tendency of transmission error, not a detailed waveform.

Conclusions
We proposed an analytical method that calculates the TE of a spur gear pair by considering various materials and statistically analyzed the range of the TE caused by various uncertainties. In contrast to the conventional approach, which considers only the mechanical uncertainty, we additionally considered the thermal deformation for tooth contact analysis. This approach is appropriate in that thermal deformation can generate a large geometric error according to the material. To demonstrate the importance of the geometric error caused by thermal deformation, the magnitudes of geometric errors caused by mechanical uncertainty and thermal deformation were compared in the case of a plastic gear. Through various numerical examples, we investigated the tendency and influence of each uncertainty on the TE. We believe that considering the uncertainties can lead to more precise gear design in terms of predicting the performance range and managing the gear quality. Although steel gears are not susceptible to thermal deformation, we demonstrated that the TE of a spur gear pair with a comprehensive error is greater than that of a spur gear pair with mechanical uncertainties due to thermal deformation. In this study, we only investigated a numerical method that calculates the TE of a spur gear pair with uncertainties. However, this approach should be validated experimentally, and the influence of thermal deformation on the TE should be examined using the proposed method. Furthermore, the proposed method can be extended to other gear systems, such as helical, worm, and harmonic drive systems, since the approach is based on the kinematic relationship between the rack and the gear. Studying spur gear pairs that are susceptible to uncertainties in the tooth geometry is an important direction for future research. Data Availability Statement: Data available on request due to privacy/ethical restrictions. The data presented in this study are available on request from the corresponding authors. The data are not publicly available due to their containing information that could compromise the privacy of research participants.

Conflicts of Interest:
The authors declare no conflict of interest. The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript, or in the decision to publish the results.

Appendix A. Mesh Stiffness
Here, as in [35], a tooth is regarded as a nonuniform cantilever beam to calculate the bending deflection δ b : 1/I i =(1/I i + 1/I i+1 )/2, Here, F is the applied force, α m is the operating pressure angle, G is the shear modulus, s h is a shear factor, e i and d i are illustrated in Figure A1, I i is the area moment of inertia, A i is the area of the tooth cross section S i , E is Young's modulus, and ν is Poisson's ratio.
The equation for the fillet-foundation deflection δ f is based on Muskhelishvili's theory, which was enhanced by Sainsot [36,37]. By assuming a linear distribution of the normal stress and a constant shear stress at the tooth root, the gear-body-induced tooth deflection can be calculated analytically as Here, W is the tooth width, while u f , S f , h f i , and θ f are given in Figure A1. Further, A i , B i , C i , D i , E i , and F i are given in Table A1, while X * i represents the coefficients L * i , M * i , P * i , and Q * i . Figure A1. Nonuniform cantilever beam gear model. From the results of [38], the contact stiffness is based on Hertzian 2D theory. According to Hertzian law, the two elastic bodies in contact can be approximated by two cylinders on the 2D plane. The contact stiffness is constant in longitudinal direction and along the entire line of action. The contact deflection δ h is given by The mesh stiffness is calculated by superimposing the various stiffnesses, namely bending stiffness, fillet-foundation stiffness, and contact stiffness, as