Spiral bevel gears are widely used in various fields such as automobiles, aviation, and mining machinery due to their high coincidence degree, strong load-carrying capacity, high transmission ratio, and high transmission efficiency. The manufacturing accuracy of spiral bevel gears is directly related to the errors of numerical control machines. Machine tool errors include geometric errors, thermal errors, and servo control errors. Among them, machine tool geometric errors are repetitive and stable, and can be accurately compensated through the numerical control system. Therefore, studying the characteristics of machine tool geometric errors and determining the key geometric error items is of great importance for tooth surface deviation compensation and the allocation and maintenance of machine tool motion accuracy.

Sensitivity analysis is a method for studying the influence of changes in input factors in a system on the output results. The analysis results are represented by sensitivity coefficients, and those with a large sensitivity coefficient indicate a higher correlation between the input factor and the output result. Sensitivity analysis is an important method for studying the relationship between machine tool geometric errors and machining accuracy. In the aspect of sensitivity analysis of machine tool geometric errors, scholars at home and abroad have conducted in-depth research.

There are mainly two types of sensitivity analysis methods: local sensitivity analysis methods represented by differential method, difference method, and perturbation method, and global sensitivity analysis methods represented by regression analysis method, variance method, and screening method. Through sensitivity analysis, key input factors and non-key input factors can be determined. In actual analysis and calculation, non-key input factors can be ignored first, which can significantly reduce the complexity of the system and the difficulty of analysis, and greatly reduce the calculation amount and processing difficulty of relevant data.

**Local Sensitivity Analysis Method**

Local sensitivity analysis is a kind of single-factor analysis method. By making a small change to a single input factor each time while keeping all other factors unchanged, the sensitivity coefficient is determined based on the differential of the output result with respect to the input parameter or the change of the output result caused by the change of a single input factor. The concept of local sensitivity is clear and the calculation is simple, which is mainly suitable for linear models and models with weak nonlinearity. The analysis process is shown in Figure 1.

[Insert Figure 1: Local Sensitivity Analysis Process]

**Sobol Global Sensitivity Analysis**

The global sensitivity analysis method not only considers the range and distribution law of each parameter, but also takes all parameters into account during the analysis and calculation, considering the influence of the mutual coupling of each input parameter on the result. However, the global sensitivity analysis method generally requires a certain scale of sampling to obtain the input parameters of the model, and then performs sensitivity analysis. Especially for the case with more parameters, the calculation amount will be relatively large.

Taking the Sobol method as an example, it is a global sensitivity method based on variance. This method mainly calculates the sensitivity coefficient of the input parameter by decomposing the variance of the model into sub-items.

Let n represent the number of input parameters of the model, and let xi represent the i-th input parameter, where i = (1, 2,…, n). Use K = f(X) to represent the multivariate function corresponding to the model, where X = (x1, x2,…, xn). f(X) is decomposed as follows:

f(X) = f0 + Σi^n fi(xi) + Σi<j^n fi,j(xi, xj) +… + f1,…,n(x1, x2,…, xn)

where f0 is the expected value of the output result f(X) calculated from the input parameters, which is a constant; fi(xi) represents the output result under the sole action of the input parameter xi; fi,j(xi, xj) represents the output result under the joint action of the input parameters xi and xj; f1,…,n(x1, x2,…, xn) represents the output result under the joint action of all input parameters; and the same applies to other orders.

The total variance D can be expressed in the following form:

D = ∫0^1…∫0^1 f^2(X)dX – f0^2

The variance component can be obtained as follows:

Di1,…,is = ∫0^1…∫0^1 fi1,…,is^2(xi1,…,xis)dxi1…dxis (3)

where 1 ≤ s ≤ n.

According to the Sobol method, by squaring and then integrating both sides of equation (1), we get:

D = Σi=1^n Di + Σi<j^n Di,j +… + D1,2,…,n

The sensitivity coefficient can be obtained through the following formula:

Si1,…,ik = Di1,…,is / D (5)

Since the variance is non-negative, it can be known that Si2,…,is is non-negative. According to equations (4) and (5), we get:

Σi=1^n Si + Σi<j^n Si,j +… + S1,2,…,n = 1

where Si is the first-order sensitivity coefficient corresponding to the variable xi, and Si,j (i < j) is the second-order sensitivity coefficient corresponding to the mutual coupling effect of the variables xi and xj. Similarly, the sensitivity coefficient of each order can be obtained. The global sensitivity coefficient Stot^i corresponding to the variable xi can be obtained by summing the sensitivity coefficients related to the variable xi, and the calculation method is as follows:

Stot^i = Si + Σi=k^s Si,k +… + S1,2,…,n

where the positive integer k takes values in the range 1 ≤ k ≤ n, and k ≠ i.

**Monte Carlo Estimation**

In the process of solving the sensitivity coefficient using the Sobol method, it involves the solution of multiple integrals, which is usually very difficult for complex models. Therefore, the Monte Carlo method is often used to approximately simulate the solution of multiple integrals. The general calculation method is as follows: perform two independent samplings on the input parameters to obtain two independent sampling matrices E and F, as follows:

E = [x11 x12… x1n

x21 x22… x2n

…………

xh1 xk2… xhn]

F = [x11′ x12’… x1n’

x21′ x22’… x2n’

…………

xk1′ xk2’… xkn’]

where k is the number of samples for each sampling.

In order to solve the first-order sensitivity coefficient and the global sensitivity coefficient, based on the matrices E and F, the matrix Er’ is constructed, where i = 1, 2,…, n. The matrix Er is obtained by replacing the i-th column of the matrix E with the i-th column of the matrix F, and the rest is exactly the same as the matrix E.

Er = [x11 x12… x1i’… x1n

x21 x22… x2i’… x2n

………………

xk1 xk2… xki’… xkn]

Therefore, each row of input parameters can separately solve for a model output result. The first-order sensitivity coefficient S and the global sensitivity coefficient Si^tot can be approximately calculated through the following formulas:

Si = (1/k) Σi=1^k f(F)i(f(Er)k – f(E)k) / D (11)

Si^tot ≈ (1/2k) Σk=1^k (f(E)h – f(Er)h)^2 / D

where f(E)s is the output result obtained by substituting the h-th row of the matrix E into the model; f(F)s is the output result obtained by substituting the h-th row of the matrix F into the model; f(Er)i is the output result obtained by substituting the h-th row of the matrix Er into the model; D is the total variance of the calculation model.

**Spiral Bevel Gear Machine Tool Processing Model**

**Machine Tool Geometric Error Classification**

The structure of the spiral bevel gear numerical control machine tool is shown in Figure 2. According to the set numerical control program, the A, B, X, Y, and Z axes move together to process the spiral bevel gear in the Gleason system. The c-axis drives the cutter head to rotate and does not affect the tooth surface generation process. The coordinate transformation relationship of the machine tool is shown in Figure 3.

For the five motion axes related to the tooth surface generation, each axis has geometric errors due to manufacturing, assembly, and other factors, which directly lead to the generation of tooth surface deviations. Each axis has 6 geometric errors, including 3 linear errors and 3 angular errors. For the entire machine tool, 30 geometric errors need to be considered [13].

[Insert Figure 2: Structure Diagram of Spiral Bevel Gear CNC Gear Milling Machine]

[Insert Figure 3: The Machining Coordinate of CNC Machining Tool]

To facilitate analysis, the geometric error sequence is numbered as shown in Table 1, where ε represents the angular error and δ represents the linear error. The error variable is distinguished by two subscripts: the first subscript represents the direction of the error, and the second subscript represents the axis where the error is located. For example, εay represents the angular error in the α direction of the Y-axis, and εyx represents the linear error in the y direction of the X-axis.

[Insert Table 1: Geometric Error and Corresponding Number of Five-Axis CNC Gear Machine]

**Tooth Surface Processing Process**

The generating motion of the workpiece gear is determined by the five axes of X, Y, Z, A, and B. Therefore, the homogeneous transformation matrix from the A-axis to the Y-axis combined with the tool equation can obtain the tooth surface equation rb under ideal conditions [14].

rz = MA · MB · MZ · MX · MY · rt(u, θ)

where M (q = X, Y, Z, A, B) is the motion transformation matrix corresponding to each axis; rt is the expression obtained from the tool equation; u and θ are the parameters of the cutter head. Due to space limitations, the specific process can refer to literature [15].

If the machine tool geometric error is considered, the actual tooth surface equation is as follows:

rg^e = MA · MA^e · MB · MB^e · MZ · MZ^e · MX · MX^e · MY · MY^e · rt(u, θ)

where Mq (q = X, Y, Z, A, B) is the geometric error matrix corresponding to each axis.

**Sensitivity Analysis Example**

**Sampling Calculation**

Taking the gear shown in Table 2 as an example, the theoretical tooth surface of the spiral bevel gear can be obtained according to equation (13). Considering the influence of the machine tool geometric error in actual processing, the actual tooth surface of the spiral bevel gear can be calculated according to equation (14).

[Insert Table 2: Geometrical and Machining Parameters]

In order to facilitate analysis and calculation, the tooth surface is represented by a discrete 15 × 9 lattice. Combining equations (13) and (14), the theoretical tooth surface lattice and the actual tooth surface lattice can be obtained, and the deviation value Kf of the corresponding tooth surface point can be solved through the coordinate values of the corresponding points in the lattice, where f represents the number of the point in the lattice. K represents the tooth surface deviation, which is the output result of the calculation model and is used to measure the magnitude of the tooth surface deviation, and can be obtained through the following formula:

K = Σj=1^13 |Kj| (15)

Given that the variation of each linear parameter is +0.01 mm and the variation of each angular parameter is +27″, taking the tooth surface deviation K as the model output result, the local sensitivity analysis result is shown in Figure 4.

[Insert Figure 4: Local Sensitivity Analysis Result]

Referring to the variation of the input parameters in the local sensitivity analysis, the range of the linear error is given as 0 – 10 μm, the range of the angular error is 0″ – 27″, and the geometric error parameters follow a uniform distribution. Write a sampling program to randomly sample within the given range of the geometric error to generate the machine tool geometric error parameters. Combining equations (13) and (14), the corresponding theoretical tooth surface and actual tooth surface can be calculated, and according to the two kinds of tooth surface lattices, the corresponding tooth surface deviation K can be solved.

Through the Sobol global sensitivity analysis method, the global sensitivity coefficients of 30 machine tool geometric errors can be solved, and arranged in descending order of the global sensitivity coefficient.

**Discussion on Input Parameter Situation**

When using the local sensitivity analysis method for sensitivity analysis, it is necessary to ensure that the variation of the input parameters of the same type is the same, so that the output results are comparable. The selection criteria are independent of the actual range and distribution law of the input parameters, so changes in the actual input parameter value range will not affect the local sensitivity analysis results.

In the previous text, the selection of the input parameter range of the global sensitivity analysis method referred to the characteristics of the local sensitivity analysis method to ensure that the range and distribution law of the same type of parameters are the same. However, in the real situation, the value range of the input parameters is generally not the same, and the change law is often different. Therefore, considering the possible actual situation, taking the geometric error with the number as an example, let the value interval of this geometric error change from 0 – 10 μm to 0 – 20 μm, while keeping other conditions unchanged, the sensitivity analysis result can be obtained.

It can be found that for the global sensitivity analysis method, after changing the value range of the input parameters, the result of the sensitivity analysis changes significantly. Among them, the global sensitivity coefficient of the geometric error item with the number, that is, δs1, changes significantly, from the original 0.023 to 0.154. The global sensitivity coefficients of other items only have minor changes, and the relevant proportions have no significant changes.

**Result Analysis**

By comparing the calculation results of the two sensitivity analysis methods before and after the change in the value range, the following characteristics of the two sensitivity analysis methods for the tooth surface machining deviation calculation model of the spiral bevel gear can be found:

(1) In the ideal case where the value range and distribution law of the input parameters are completely the same, the calculation amount of the local sensitivity analysis method is smaller, and the preliminary analysis results can be obtained quickly. At the same time, the distribution law of the analysis results is the same as that of the global sensitivity analysis results, which has certain reference value.

(2) The Sobol global sensitivity analysis method considers the influence of the value range and distribution law of the input parameters on the analysis results, which is more scientific, the conclusions obtained are more credible, and the applicable range is wider. In addition, most global sensitivity analyses can be used to study the quantitative relationship between the input parameters and the output results and the influence of the interaction between the parameters on the output results. However, compared with the local sensitivity analysis method, the calculation amount increases significantly, and more time is required.

(3) For models with unclear distribution range of input parameters, less available time, and the same or similar distribution of linear or weakly nonlinear models, the local sensitivity analysis can be considered; for nonlinear models with complex distribution range of input parameters and high accuracy requirements, the global sensitivity analysis method can be considered.

**Conclusion**

By applying the two sensitivity analysis methods to the sensitivity analysis of the tooth surface machining deviation model of the spiral bevel gear, the following can be found:

(1) In the ideal case where the value range and distribution law of the input parameters are completely the same, the local sensitivity analysis method has a smaller calculation amount and can obtain preliminary analysis results more quickly. At the same time, the distribution law of the analysis results is the same as that of the global sensitivity analysis results, which has certain reference value.

(2) The Sobol global sensitivity analysis method considers the influence of the value range and distribution law of the input parameters on the analysis results, which is more scientific, and the conclusions obtained are more credible, with a wider range of application. In addition, most global sensitivity analyses can be used to study the quantitative relationship between the input parameters and the output results and the influence of the interaction between the parameters on the output results. However, compared with the local sensitivity analysis method, the calculation amount is significantly increased, and more time is required.

(3) For models with unclear distribution range of input parameters, less available time, and the same or similar distribution of linear or weakly nonlinear models, the local sensitivity analysis can be considered; for nonlinear models with complex distribution range of input parameters and high accuracy requirements, the global sensitivity analysis method can be considered.