Multi-objective optimization of gearbox based on panel acoustic participation and response surface methodology

Abstract

Aiming at noise problem of two-stage gearbox, this paper uses the panel acoustic contribution analysis (PACA) and response surface methodology to conduct structural-borne acoustics analysis and multi-objective optimization of gearbox. Using this optimization, optimal layout and parameters of the added ribs can be determined. The dynamic model of a spur gear-shaft-bearing system is built in consideration of the time-varying mesh stiffness, the time-varying bearing stiffness, and the flexibility of the shaft, while dynamic load of the bearing is obtained through solving the model. The dynamic load of bearings is considered as excitation to analyze the vibration and noise radiation characteristics of housing. Furthermore, acoustic transfer vector, modal acoustic contribution and PACA are used to select the region with the most acoustic contribution. The variable ribs structure model (VRSM) is established in this region. Finally, the VRSM are adopted to optimize the housing. The optimization model based on the peak sound pressure level (PSPL) and the ribs mass (R_m) is established by using response surface method. Through multi-objective optimization, the variable in VRSM can be modified or deleted to reduce the radiated noise of the gearbox. The optimization results demonstrate that the PSPL after optimization is decreased by 4.88 dB at frequency of 3130 Hz. At the same mass, the noise reduction effect of VRSM is 12% higher than that of standard rib. Therefore, the optimization effect of the housing structure is obvious. PACA-VRSM-RSM- NSGA-II technology provides a new research scheme for reducing the vibration and noise of the gearbox.

Introduction

Gear system was widely used in automotive, aerospace, and power generation industries due to its compact structure and high transmission efficiency. The vibration of gearbox is mostly caused by the meshing action of gear teeth, in which it passes through gears, axles and bearings, and finally transfers to housing, which causing forced vibration of housing and radiating noise into the surrounding environment. The finite element/boundary element method (FEM/BEM) is the prevailing method applied to the numerical analysis for the gearbox vibration and noise.^1–3 Kubur et al. ⁴ proposed a dynamic model of a multi-shaft helical gear reduction unit formed by N flexible shafts, and studied the influence of bearing stiffness, shaft size, and other parameters on the vibration of the gearbox. Later, Wei et al. ⁵ established the coupled nonlinear dynamics model of a wind turbine planetary gear transmission system, and proposed a lightweight and reliability constrained optimization design method according to the sequential quadratic programming algorithm. Subsequently, Chang et al. ⁶ proposed a comprehensive coupled gearbox system dynamics model of a gear-shaft-bearing-shell system, which improved the accuracy of bearing response and thus more accurately estimated the radiation noise of the gearbox. In order to solve the noise problem generated by operating gearboxes, Han et al. ⁷ proposed a computational method for predicting gearbox noise from steady dynamic response to acoustic radiation calculation.

Based on the gear dynamics model, many scholars have provided different schemes to reduce the radiated noise of the gearbox. To reduce the periodic change of meshing stiffness, Dogruer et al. ⁸ proposed a nonlinear controller which adjusts the torque acting on the driving gear to minimize the adverse effect of time-varying mesh stiffness. And a spur gearbox was designed which does not emit much noise and vibration. Jolivet et al. ⁹ addressed the issue of the impact tooth flank finishing on noise and vibration. Ghosh et al. ¹⁰ used the genetic algorithm to optimize the geometrical shape of the gear, reducing the noise of the tractor gearbox. In addition, the surface ribbing of the shell is one of the main methods of noise reduction. Many of structures contain an important account of ribs. Moyne et al. ¹¹ discussed the influence of adding ribs to an automobile gearbox shell on the radiation noise. The influence of ribs was theoretically discussed and three different physical effects were distinguished: the “vibrating,” the “obstacle,” and the “source” effects. Furthermore, panel acoustic contribution analysis (PACA) is an important method for optimizing the shell structure, which can determine the approximate position of the ribs plate distribution. Based on the dynamics model of a gear system, Zhou et al. ¹² added ribs using the acoustic contributions of the surface of the gearbox, and designed two gear stiffness improvement schemes to reduce the radiated noise of the gearbox. Although PACA can identify the panel area of the greatest acoustic contribution, it cannot determine the exact location of the ribs in this area, nor the geometry of the added ribs, which requires further calculations. Wang et al. ¹³ proposed an optimal ribs layout design method to reduce radiation noise based on topology optimization and acoustic contribution analysis, verified the accuracy of the result by modal test and acoustic measurement. However, the results of topology optimization are complex, diverse and difficult to achieve. Compared to topology optimization, traditional structure optimization is often costly and inefficient. Artificial intelligence algorithms can significantly reduce calculation time and material cost. Wang et al. ¹⁴ proposed a multi-objective optimization design based on the response surface method and PACA, and pasted a damping layer on the panel surface to reduce radiation noise. Zhou et al. ¹⁵ adopted the free damping structure and new stiffener structure to optimize the rear end cover of a double planetary gear housing. The optimization model based on the vibration acceleration of the rear end cover surface was established by applying K-S function and response surface method.

However, there has been no serious attempt to determine the optimum layout and parameters of ribs. To address this, our investigation adds acoustic transfer vector (ATV) and modal acoustic contribution (MAC) to help identify the regions with the greatest acoustic contribution, enriching the screening process of the PACA method. Based on ATV-MAC-PACA, the variable ribs structure model (VRSM) as design variables in which utilized ribs width, thickness, number, and direction is proposed. Through multi-objective optimization, the variable in VRSM can be modified or deleted to reduce the radiated noise of the housing. PACA-VRSM-RSM- NSGA-II technology provides a new research scheme for the layout of gearbox ribs.

Design methods and processes

This section analyzes the vibration and noise calculation method of a two-stage gearbox system. And the PACA-VRSM-RSM-NSGA-II method for the gearbox structure optimization is as shown in Figure 1. The entire analysis process is composed of five components, namely the (a) gearbox dynamic excitation calculation, (b) gearbox dynamic response analysis, (c) gearbox radiated noise calculation, (d) determine ribs addition area and VRSM, and (e) Multi-objective optimization.

Figure 1.

Multi-objective optimization calculation process of ribs layout of the two-stage gearbox.

(a) Gearbox dynamic excitation calculation. Analytical model is established by using lumped mass method. In the model, the influence of the TVMS, TVBS, and shaft flexibility is considered. The dynamic load of the bearings is calculated by the Newmark integral method. ¹⁶

(b) Gearbox dynamic response analysis. Modal analysis of the shell finite element model with boundary conditions is carried out. The dynamic load of the system is then added to the corresponding position of the shell finite element model. The modal superposition method is used to solve the dynamic response of the gearbox structure, and the vibration response and surface acceleration of the housing are also obtained.

(c) Gearbox radiated noise calculation. The surface element of the finite element model is extracted, and the shell element is used to fill the holes on the surface of the structure. Finally, the gearbox boundary element model is obtained by remeshing the surface elements of gearbox. And the obtained boundary element model and ATV are used to solve the radiated noise of the gearbox.

(d) Determine ribs addition region and VRSM. Three screening methods, ATV, MAC, and PACA, were used to determine the region with the greatest acoustic contribution. And the variable ribs model is added in this region.

(e) Multi-objective optimization. Eleven structural parameters of the VRSM are used as input variables. PSPL and R_m are used as output variables. Response surface models of PSPL and R_m are then constructed by central composite design (CCD), and R² is used to test the accuracy of the model. When the model accuracy meets the requirements, the obtained response surface objective function obtained is exported to NSGA-II, and the optimal ideal point is determined by analyzing the Pareto solution set. Finally, an experimental comparison is conducted to obtain the ideal VRSM.

Gearbox dynamic excitation calculation

The object of analysis in this paper is a two-stage gear reducer, whose model is as shown in Figure 2, which is a two-stage straight toothed cylinder gear reducer, and the minor chamfering and feature in structure is simplified in the model.

Figure 2.

(a) 3D model of two-stage gearbox, (b) Internal structure of the gearbox.

The specific parameters of the gear are shown in Table 1, and the bearing parameters are shown in Table 2.

Table 1.

Gear system parameters.

Parameters	1st stage driving gear	1st stage driven gear	2nd stage driving gear	2nd stage driven gear
Number of teeth	36	90	29	100
Modulus (mm)	1.5	1.5	1.5	1.5
Pressure angle (o)	20	20	20	20
Tooth width (mm)	12	12	12	12

Table 2.

Bearing parameters.

Diameter of inner race (mm)	Diameter of outer race (mm)	Ball number	Ball diameter (mm)
28.7	46.6	8	8.7
Pitch circle diameter (mm)	Radial internal clearance	Curvature radius (mm)	Goodness fit
37.65	0.5	4.5	0.5172

Gear transmission system dynamic model

The spur gear transmission system model has six degrees of freedom (DOF), namely rotation DOF of driving and driven wheels, the translation DOF along x and y directions (see Figure 3). Then the relative total displacement of meshing elements of spur gear pair δ is defined as

$$\delta =x_p\sin \alpha +y_p\cos \alpha -r_p{\theta }_p-x_g\sin \alpha -y_g\cos \alpha -r_g{\theta }_g-e\left(t\right)$$

(1)

where α is the pressure angle, r _p is the base circle radius of the driving gear, and r _g is the base circle radius of the driven gear.

Figure 3.

Dynamics model of gearbox transmission system.

Considering the influence of time-varying meshing stiffness and the gear meshing error, the motion differential equation of the meshing element of the spur gear can be expressed as

$$\{\begin{array}{c}{m}_p{\ddot{x}}_p+c_m\dot{\delta } \sin \alpha +k_m\delta \sin \alpha +f_s\sin \alpha =0\\ m_p{\ddot{y}}_p+c_m\dot{\delta } \cos \alpha +k_m\delta \cos \alpha +f_s\cos \alpha =0\\ I_p{\ddot{\theta }}_p-c_m\dot{\delta }{r}_p-k_m\delta r_p-f_s{r}_p=0\\ m_g{\ddot{x}}_g-c_m\dot{\delta } \sin \alpha -k_m\delta \sin \alpha -f_s\sin \alpha =0\\ m_g{\ddot{y}}_g-c_m\dot{\delta } \cos \alpha -k_m\delta \cos \alpha -f_s\cos \alpha =0\\ I_g{\ddot{\theta }}_g-c_m\dot{\delta }{r}_g-k_m\delta r_g-f_s{r}_g=0\end{array}$$

(2)

where m is the mass of gear, k _m is the meshing stiffness, c _m is the meshing damping, I is the moment of inertia of gear and f _s is the normal impact force.

During the operation of the gear transmission system, the bearing supports the shaft system and transmits vibration from the gears to the housing. ¹⁸ The bearing stiffness matrix can be expressed as

$$K_m^e=\mathrm{diag}\left(k_{bx},k_{by},0\right)$$

(3)

The finite element model of the gearbox (see Figure 4(a)) is established according to the gear meshing, bearing, and shaft elements, with a total of 29 elements and 32 nodes. ¹⁷ The stiffness matrix of each element is superimposed according to the position of the element nodes to obtain the overall stiffness matrix, as shown in Figure 4(b). The generalized displacement vector of the system is

$$\left\{x\right\}=\left\{x_{1x},x_{1y},x_{1\theta },x_{2x},x_{2y},x_{2\theta }\dots \dots x_{nx},x_{ny},x_{n\theta }\right\}$$

(4)

where n denotes the number of nodes.

Figure 4.

Solution of overall stiffness of gearbox. (a) Finite element model of gearbox, (b) Overall stiffness matrix of gearbox.

The dynamics equation of the gear system can be expressed as

$$\left[m\right]\left\{\ddot{x}\left(t\right)\right\}\mathrm{+}\left[c\right]\left\{\dot{x}\left(t\right)\right\}\mathrm{+}\left[k\left(t\right)\right]\left\{x\left(t\right)\right\}=F\left(t\right)$$

(5)

where {x(t)} is the generalized nodal displacement vector, [m], [c], and [k] are the matrix of global mass, damping, and stiffness, respectively, and {[F(t)]} is the external force vector received by the system.

Bearing dynamic load

The input speed is 5400 r/min, and the output torque is 100 Nm. The time domain and frequency domain diagrams of dynamic loads of each bearing are obtained by Newmark integral method, as shown in Figure 5. where different colors represent different bearings. f _m1 (3240 Hz) and f _m2 (1044 Hz) represent the meshing frequency of the first gear pair and the second gear pair, respectively.

Figure 5.

The bearing dynamic load of the gear system. (a) Time response of loads from bearing 1 to bearing 6, (b) Spectrum of loads from bearing 1 to bearing 6.

Experimental verification

The platform used in this experiment is the SQI Wind turbine transmission system test bench (shown in Figure 6(a)). The main devices include the test bed, drive motor, controller, gearbox, and data collector, which can be used to study the dynamic analysis of the gear transmission in this paper. The driving motor is used to adjust the speed automatically, increasing 240 rpm each time, and the acceleration value after each speed change is recorded by the acceleration sensor (shown in Figure 6(b)). The experimental results are compared with the simulation results (see Figure 7). It can be seen that the simulation results are very close to the trend of the experiment, indicating that the simulation results can replace the experiment.

Figure 6.

Transmission system test bench. (a) System equipment, (b) Sensor measurement position.

Figure 7.

Comparison between simulation and experiment.

Gearbox dynamic response analysis

Figure 8 is the finite element model of the gearbox. The material of the plate and ribs are both steels having Young’s modulus E = 2.07 × 10¹¹ Pa, Poisson’s ratio μ = 0.3, and mass density ρ = 7800 kg/m³. In order to apply dynamic load to the gearbox, a mass node is created in the center of the bearing hole and the node is coupled with the node on the inside surface of the bearing hole. Then the dynamic load is applied on the central node. A displacement constraint is applied at the bottom of the gearbox. The finite element model of housing and the first eight natural vibration modes of the gearbox are shown in Figure 8. When the bottom of the gearbox is restrained, the vibration of the lower part of the gearbox is weak, while the vibration amplitude of the top part is obvious.

Figure 8.

Finite element model of housing, (b) Mode of vibration of gearbox.

Gearbox radiated noise calculation

The BEM is obtained by extracting surface elements from the gearbox finite element model and repairing the bearing holes with shell elements. The acoustic fluid is air having mass density ρ0 = 1.2 kg/m3 and sound speed c = 343.4 m/s. Based on the boundary element model, the hemispherical shape acoustic field is established. Four points with the same radius are selected on the acoustic field as the measured field points, namely, field points A, B, C, and D, as shown in Figure 9(a). In the case of small perturbations, the acoustic equation can be considered linear. Thus, a linear relationship can be established between the input (the normal vibration of the structure surface) and the output (the sound pressure at a point in the sound field). According to the relationship shown in Figure 9(b), the sound pressure at a certain field point can be expressed as

$$p\left(\omega \right)=ATV{\left(\omega \right)}^T{v}_n\left(\omega \right)$$

(6)

where ATV is the acoustic transfer vector calculation, and ω is the angular frequency.

Figure 9.

(a) Gearbox field point distribution, (b) ATV.

The boundary element method is used to solve the noise radiation of the gearbox by taking the normal vibration acceleration of the shell’s outer surface as the boundary condition. Figure 10 shows the noise spectrum when the input speed is 5400 rpm, and the torque is 100 N·m. Due to the consideration of shaft flexibility, the overall noise increases first and then decreases. From the sound pressure level curves, it is well noticed that peak is generated on the position of 3 × fm2. The corresponding frequency of the maximum radiated noise is 3 × f _m2 (3132 Hz), which is close to the eighth order natural frequency of housing (3098.35 Hz). At 3 × f _m2 (3132 Hz), the SPL at field point A is the largest, so field point A is selected as the target field point. The radiated noise at point A can be reduced by suppressing the PSPL at 3 × f _m2 .

Figure 10.

Noise spectrum of gearbox.

Determine ribs addition region and VRSM

In this study, the regions with larger ATV can be found by the analysis of ATV. Then, the regions with larger MAC can be obtained by MAC. Finally, panels can be divided on the regions with larger MAC and ATV concurrently and the acoustic contribution of each panel can be obtained by PACA. This process can be found the region with maximum acoustic contribution. Based on ATV-MAC-PACA, the VRSM was established.

Acoustic transfer vector (ATV)

Figure 11 demonstrates the amplitude of the acoustic transfer vector at 3 × f _m2 (3132 Hz) to field point A. It can be seen that the amplitude on the upper part of the gearbox is significantly greater than that on the bottom, indicating that the normal vibration on the upper part of the gearbox has a more obvious impact on the PSPL. In summary, the first screening region is the upper part of the gearbox.

Figure 11.

Amplitude of acoustic transfer vector at field point A at 3 × f _m2 .

Modal acoustic contribution analysis (MACA)

Using the modal superposition method, the displacement response of the multi-free system can be obtained as

$$x\left(\omega \right)={\displaystyle \sum _{k=0}^n\frac{f\left(\omega \right){\phi }_k^T}{a_k\left(i\omega -{\lambda }_k\right)}}{\phi }_k$$

(7)

where k is the order of the mode, $q_k\left(\omega \right)=\frac{f\left(\omega \right){\phi }_k^T}{a_k\left(i\omega -{\lambda }_k\right)}$ is the structural modal participation factor, λ _k is the k-th eigenvalue of the system, {φ} _k denotes the mode shape of the k-th order; {f (ω)} is the load vector, and α _k denotes a constant related to the characteristics of the system.

By projecting {x(ω)} onto the normal direction of the box surface and taking the derivative, the vibration velocity in the normal direction of the gear box surface can be obtained as

$$v_n\left(\omega \right)=i\omega {\displaystyle \sum _{k=1}^n{q}_k\left(\omega \right){\phi }_{jn}}$$

(8)

where {φ} _jn is the component of the mode shape on the normal line of the gearbox surface.

By substituting equations (9) into (10), the sound pressure at any point in the sound field can be obtained as

$$p\left(\omega \right)=ATV{\left(\omega \right)}^{T}i\omega {\displaystyle \sum _{k=1}^n{q}_k\left(\omega \right)}{\phi }_{nk}={\displaystyle \sum _{k=1}^{n}i\omega q_k}\left(\omega \right)ATV{\left(\omega \right)}^T{\phi }_{nk}={\displaystyle \sum _{k=1}^n{P}_{sk}}$$

(9)

where $P_{sk}={\displaystyle \sum _{k=1}^{n}i\omega q_k}\left(\omega \right)ATV{\left(\omega \right)}^T{\phi }_{nk}$ denotes the sound pressure generated by mode K. The projection of P _sk in the direction of the total sound pressure P(ω) is called the modal acoustic contribution.

The acoustic contribution of each mode at field point A is obtained by MACA, as shown in Figure 12. It can be seen that the eighth mode has a large acoustic contribution to the PSPL at field point A. The corresponding main mode shape of this mode is shown in Figure 12. The eighth main mode vibrates obviously at the top of the gearbox. Therefore, the result of the second screening is the top region of the gearbox.

Figure 12.

The modal acoustic contributions at the field point A at the meshing frequency.

Panel acoustic contribution analysis (PACA)

The main panels corresponding to the first eight modes are divided, and PACA analysis was carried out. The acoustical contribution of each panel to the field point A at 3×f _m2 (3132 Hz) was obtained, as shown in Figure 13. It can be seen from the figure that panel 1 and panel 2 have the largest contribution, more than 25 percent, so the screening results through PACA are panel 1 and panel 2.

Figure 13.

Analysis results of acoustical contribution of plate surface at meshing frequency.

In summary, Using the ATV-MAC-PACA method, the central region of the top of the gearbox has the largest acoustic contribution to the field point A at 3×f _m2 (3132 Hz). Therefore, the VRSM can be added at this region to reduce the radiated noise of the gearbox.

Variable ribs structure model (VRSM)

Although the general region of the added ribs can be obtained from the previous screening, the research on the distribution and number of ribs and the specific parameters of each rib requires further design. The ribs layout method adopted in this paper refers to literature. ¹⁹

Figure 14(a) is the rectangular region determined according to the screening results. To study the effect of the layout and parameters of the ribs on reducing the radiated noise of the shell, the ribs model is added to the top center of the shell (note: all model ribs are added to the inner surface of the top of the box body to prevent disturbance of the ATV and improve calculation efficiency). The characteristic of VRSM is characterized by the fact that there is no definite ribs layout and size, and variables in the model can be changed as required. Thus, different VRSMs can be realized according to different optimization results. As shown in Figure 14(b), the model contains 11 design variables; namely, the number of horizontal and vertical ribs (n ₁ and n ₂ ), the width of the ribs (d ₁ , d ₂ , d ₃ , and d ₄ ), and the thickness of the ribs (t ₁ , t ₂ , t ₃ , and t ₄ ). φ is the rotation angle of the ribs relative to the x-axis (shown in Figure 14(c)). In addition, the spacing y ₁ and y ₂ between the ribs are determined by the number of horizontal and vertical ribs, respectively; y ₁ and y ₂ reflect variation the density of rib distribution in the gearbox. It can be seen that through the selection of value range, the rib height and thickness variables can be set to zero, and then the redundant ribs can be deleted in the corresponding optimization.

Figure 14.

3D model of ribs layout design. (a) the position where the ribs are added in the gearbox, (b)VRSM, and (c) variable φ: the rotation Angle of the ribs relative to the x-axis.

Multi-objective optimization

There is a great relationship between the mass of ribs and peak value of radiated noise of the gearbox, and the two are in conflict in essence, so it is very necessary to establish a multi-objective optimization model.

Response surface model (RSM)

Based on statistical and mathematical methods, RSM can build and evaluate a variety of complex problems. In this method, the output variable (system response) is approximately fitted to the specified design space to replace the real response surface. ¹⁸ In order to give consideration to the fitting accuracy and optimization calculation time, the central composite design (CCD) was adopted. In this design, alpha = 0.5, the center point is 1 and the factorial core uses the min-res V method. The relation between the independent variable in the response surface and response Z can be expressed as

$$\stackrel{\mathrm{\Lambda }}{Z}=f\left(n_1,n_2,d_1,d_2,d_3,d_4,t_1,t_2,t_3,t_4,\phi \right)+\epsilon$$

(10)

where ε is the calculation error.

RSM usually adopts quadratic polynomial model to fit the parameter model of the response. Its approximate function can be expressed as

$$\stackrel{\mathrm{\Lambda }}{Z}=c_0+{\displaystyle \sum _{i=1}^n{c}_i}{x}_i+{\displaystyle \sum _{i=1}^n{c}_{ii}}{\left(x_i\right)}^2+i=1nj>icijxixj$$

(11)

where Z denote the estimated response; c ₀ denote the constant term coefficient; c _i denotes the coefficients of the input factor x _i ; c _ii denotes the quadratic coefficients on the independent variable x _i and c _ij denotes the coefficients of the interaction term between the input variable x _i and x _j .

Multi-objective genetic algorithm

Genetic Algorithms is a heuristic algorithm developed on the basis of biological evolutionary theory and Genetic thought. Compared with single objective optimization, multi-objective optimization can better solve problems with multiple objective requirements. NSGA-II proposed by Deb et al. ²⁰ is a multi-objective optimization algorithm based on Pareto optimal solution, that is, a fast non-dominated multi-objective optimization algorithm with elite retention strategy. The three major characteristics of NSGA-II are the proposed congestion degree and congestion comparison operator, and the introduction of non-dominated ranking and elite strategy.

In this paper, the optimal VRSM is solved while considering multiple targets. The objective functions of Rm and PSPL are obtained by a large number of experimental data and RSM. The optimization problem model is as follows

$$\{\begin{array}{l}\mathrm{Minimize}{R}_m=f_1\left(n_1,n_2,d_1,d_2,d_3,d_4,t_1,t_2,t_3,t_4,\phi \right)\\ \mathrm{Minimize}PSPL=f_2\left(n_1,n_2,d_1,d_2,d_3,d_4,t_1,t_2,t_3,t_4,\phi \right)\\ \mathrm{subject to constraints}\\ 2\le n_1,n_2\le 6\\ 0\le d_1,d_2,d_3,d_4\le 5\\ 0\le t_1,t_2,t_3,t_4\le 6\\ 0\le \phi \le 90\end{array}$$

(12)

where f ₁ (·) and f ₂ (·) denote objective functions for R_m and PSPL, respectively.

Results and discussion

As per CCD design, 91 experiments are performed for VRSM. See Appendix 1 for specific simulation experiment data. The regression analysis is an assessment of the accuracy and predictive power of response surface models. When the response surface model is more accurate, the numerical model will be closer to the physical model, and the predicted value will be more accurate. Generally, R² and adjusted R² are used to characterize the accuracy of response surface model, and the mathematical formulations are given by

$$R^2=1-\frac{S_{SE}}{S_{ST}}$$

(13)

$$R_{adj}^2=\frac{S_{SE}\left(N-K-1\right)}{S_{ST}\left(N-1\right)}$$

(14)

$$S_{SE}={\displaystyle \sum _{i=1}^N{\left(Y_i-y_i\right)}^2}$$

(15)

$$S_{ST}={\displaystyle \sum _{i=1}^N{Y}_i^2}-{\left({\displaystyle \sum _{i=1}^N{y}_i}\right)}^{2}N$$

(16)

where S _SE denotes the sum of squares of regression residuals, S _ST denotes the sum of squares of total variation, N denotes the number of experimental sites, K denotes the number of unknown coefficients, Y _i denotes the response vector of the test point, y _i denotes the RSM approximate estimate vector.

Quadratic equation for PSPL

A mathematical model for prediction of PSPL using VRSM as the desiccant material is described in equation (17)

$$\begin{array}{l}PSPL=93.11-1.77n_1+0.86n_2-0.18d_1+0.14d_2-0.20d_3+0.45d_4+0.05t_1+0.14t_2\\ +0.15t_3+0.07t_4+0.01\phi -0.02n_1{d}_3+0.01n_1{t}_2-0.04n_1{t}_3-0.03n_2{d}_1-0.04\times n_2{d}_4\\ -0.03n_2{t}_2-0.01d_1{d}_2+0.02d_1{d}_3-0.04d_1{t}_1-0.01d_1{t}_3-0.01d_2{t}_1-0.04\times d_2{t}_2+0.02d_2{t}_3\\ -0.01d_2{t}_4+0.001d_2\phi +0.01d_3{d}_4+0.01d_3{t}_1-0.08d_3{t}_3+0.02d_3{t}_4-0.02d_4{t}_1-0.02d_4{t}_2\\ +0.04d_4{t}_3-0.07d_4{t}_4+0.01t_1{t}_2+0.01t_2{t}_3+0.02t_3{t}_4-0.0006t_4\phi +0.22n_1^2-0.09n_2^2\\ +0.06d_1^2-0.03d_2^2+0.03d_3^2-0.06d_4^2-0.003t_1^2-0.02t_2^2-0.03t_3^2-0.01t_4^2-0.0001{\phi }^2\end{array}$$

(17)

The response surface model of PSPL is established through the above calculation formula. Before using the mathematical model, the goodness of developed models is assessed by evaluating the coefficient of determination (R²). The values of R² = 0.9832 and adjusted R² = 0.9632 are very close to 1, and the difference between the predicted value and the adjusted value is less than 0.2. Therefore, the mathematical model established for the optimization objective meets the accuracy requirements. As shown in Table 3, when the resulting p-value is less than 0.0001, the model can be considered significantly important. The items exceeding the recognized standard value of 0.05 can then be deleted to improve the model accuracy. Nevertheless, the purpose of this study is to investigate the effects of all parameters on response variables, and as many parameters as possible can be retained for research when the accuracy is satisfied. According to the sum of squares, the linear term has a contribution of 44.22%. The contribution of the cross-interaction term is 36.45%, and the contribution of self-interaction term was 0.11%. Figure 15 shows the two-dimensional contour plot of part of the response surface, in which the contour plot of the PSPL is in good agreement with the experimental data. The contour lines in the figure are divided into five levels, from red to blue, indicating gradual changes from high to low. It can be observed that when the other nine variables are fixed, and the horizontal coordinate and vertical coordinate variables change, PSPL will also change accordingly (Figure 15).

Table 3.

Factors and levels of the central composition design.

Source	PSPL		Rm		Source	PSPL		Rm
	Sum of squares	p-value prob > f	Sum of squares	p-value prob > f		Sum of squares	p-value prob > f	Sum of squares	p-value prob > f
Model	176.62	<0.0001	1,564,182	<0.0001	d2t3	1.13	0.0003	2576	0.0002
n1	4.10	<0.0001	43,291	<0.0001	d3d4	0.34	0.0363	1605	0.0019
n2	2.86	<0.0001	28,604	<0.0001	d3t1	0.66	0.0047	924	0.0143
d1	5.63	<0.0001	46,528	<0.0001	d3t3	17.62	<0.0001	114,466	<0.0001
d2	3.07	<0.0001	28,738	<0.0001	d3t4	0.79	0.0021	4168	<0.0001
d3	13.36	<0.0001	93,237	<0.0001	d4t2	1.48	0.0001	4263	<0.0001
d4	15.23	<0.0001	94,319	<0.0001	d4t3	3.43	<0.0001	9530	<0.0001
t1	2.47	<0.0001	12,380	<0.0001	d4t4	14.21	<0.0001	90,296	<0.0001
t2	5.44	<0.0001	49,436	<0.0001	t1t2	0.38	0.0286	3470	<0.0001
t3	16.10	<0.0001	117,150	<0.0001	t2t3	0.46	0.0162	1018	0.0105
t4	9.83	<0.0001	69,887	<0.0001	t3t4	1.24	0.0002	3643	<0.0001
φ	0.03	0.5288	721	0.0283	t4φ	0.36	0.0334	2410	0.0003
n1d3	0.77	0.0024	5360	<0.0001	n12	0.11	0.2308	23.76	0.6763
n1t3	2.90	0.0358	20,953	<0.0001	n22	0.02	0.6147	4.06	0.8627
n2d1	0.86	0.0000	2904	0.0001	d12	0.02	0.6163	0.43	0.9550
$n_2{d}_4$	2.58	0.0014	682	0.0326	d22	0.004	0.8164	0.72	0.9419
n2t2	0.99	<0.0001	5197	<0.0001	d32	0.01	0.7762	0.60	0.9469
d1d2	0.22	0.0007	2434	<0.0001	d42	0.02	0.6339	0.64	0.9455
d1d3	0.52	0.0885	1591	<0.0001	t12	0.001	0.9666	3.30	0.8762
d1t1	5.66	0.0113	38,861	<0.0001	t22	0.002	0.8589	4.19	0.8607
d1t3	0.26	<0.0001	2501	<0.0001	t32	0.01	0.6928	13.85	0.7498
d2t1	0.58	0.0682	7048	<0.0001	t42	0.002	0.8803	14.97	0.7403
d2t2	4.50	0.0074	30,060	<0.0001	φ2	0.01	0.7333	1.85	0.9071
Model statistics
R2			0.9832	0.9979
Adjusted R2			0.9632	0.9924
Predicted R2			0.9120	0.9445
Adeq precision			29.5684	61.1448

Figure 15.

Contour plots for PSPL in terms of (a) n₁vs. d₃, (b) n₂vs. d₄, (c) n₂vs. t₂, (d) d₁vs. d₃, (e) d₂vs. t₁, (f) d₂vs. t₂, (g) d₂vs. t₃, (h) d₃vs. d₄, (i) d₃vs. t₁, (j) d₄vs. t₂, (k) d₄vs. t₃, (l) d₄vs. t₄, (m) t₁vs. t₂, (n) t₂vs. t₃, (o) t₃vs. t₄, and (p) t₄vs. φ.

Figure 16.

Contour plots for R_m in terms of (a) n₁vs. d₃, (b) n₂vs. d₄, (c) n₂vs. t₂, (d) d₁vs. d₃, (e) d₂vs. t₁, (f) d₂vs. t₂, (g) d₂vs. t₃, (h) d₃vs. d₄, (i) d₃vs. t₁, (j) d₄vs. t₂, (k) d₄vs. t₃, (l) d₄vs. t₄, (m) t₁vs. t₂, (n) t₂vs. t₃, (o) t₃vs. t₄, and (p) t₄vs. φ.

Quadratic equations for R_m

Regression function concerning R_m can be expressed as follows

$$\begin{array}{l}{R}_m=65.13+8.56n_1-4.92n_2-20.41d_1-2.79d_2-1.69d_3-9.42d_4-8.56t_1-10.66t_2\\ -19.93t_3-8.89t_4+0.46\phi +2.88n_1{d}_1+2.33n_1{d}_3+1.25n_1{d}_4+0.61n_1{t}_1+3.63n_1{t}_3+1.20n_1{t}_4\\ -1.8n_2{d}_1+1.36n_2{d}_2+0.79n_2{d}_3+2.43n_2{d}_4+2.14n_2{t}_2+2.08n_2{t}_3-0.65n_2{t}_4+0.07n_2\phi \\ +1.07d_1{d}_2-1.07d_1{d}_3-0.42d_1{d}_4+4.01d_1{t}_1+1.10d_1{t}_3+0.34d_1{t}_4+0.42d_2{d}_3-0.82d_2{d}_4\\ +1.67d_2{t}_1+3.69d_2{t}_2-1.10d_2{t}_3-0.09d_2\phi -0.94d_3{d}_4-0.67d_3{t}_1-0.62d_3{t}_2+6.81d_3{t}_3\\ -1.27d_3{t}_4+1.44d_4{t}_2-2.04d_4{t}_3+7.08d_4{t}_4-0.95t_1{t}_2-0.37t_1{t}_3-0.37t_1{t}_3-0.62t_1{t}_4\\ +0.04t_1\phi -0.67t_2{t}_3+0.07t_2\phi -1.10t_3{t}_4+0.05t_3\phi +0.06t_4\phi -3.29n_1^2-1.36n_2^2-0.28d_1^2\\ -0.37d_2^2-0.34d_3^2-0.34d_4^2+0.54t_1^2+0.61t_2^2+1.12t_3^2+1.16t_4^2-0.002{\phi }^2\end{array}$$

(18)

The above regression function is mathematically verified by the ANOVA test, and the obtained results are shown in Table 3. Based on the sum of squares, it can be seen that the contribution ratio of the linear term is 37.35%, and the contribution ratio of the cross-action term is 22.76%. According to the corresponding calculated statistical error value R ² = 0.9979 and adjusted R ² = 0.9924. The results show that the experimental data are well-adapted to the established model. Figure 16 depicts the contour map of part of R_m, showing the influence of each parameter and the interaction between the parameters on the corresponding response.

Optimization results

Statistical analysis shows that the results obtained by RSM correspond well with the experimental data. Therefore, the regression equation obtained by the response surface method can be used as the objective function of the subsequent multi-objective algorithm. NSGA-II is used to calculate the optimal solution set of PSPL and Rm. NSGA-II has been implemented in this study for a population size of 50, a crossover fraction of 0.8, a Pareto front population fraction of 0.35. Due to the conflict and competition among optimization objectives, multi-objective optimization has no unique solution, and its solution is a set of optimal non-inferior solutions. It is noteworthy that in the present study, NSGA-II has been run 20 times for the case of optimization, and corresponding first ranked Pareto fronts have been selected. Figure 17 shows the Pareto solution set for PSPL and Rm solved by NSGA-II. The coordinates of the points in the Pareto solution set in the figure are derived from Table 4. Table 4 shows the optimal values of the two objective functions and the corresponding values of the 11 parameters. The values of variables in the table are rounded results.

Figure 17.

Pareto solution set calculated with NSGA-II.

Table 4.

The corresponding parameter variables and corresponding values of the output variables in the Pareto solution set.

Variables											Objective function
n1	n2	d1	d2	d3	d4	t1	t2	t3	t4	φ	FSPL	Rm
4	6	3	3	2	5	5	6	1	3	11	87.68	280.24
4	5	5	3	5	5	6	4	6	6	78	86.49	521.49
5	5	4	3	3	5	6	4	6	6	80	86.71	499.77
5	5	4	4	3	5	6	4	6	6	86	86.79	472.89
4	6	3	5	5	4	5	6	6	1	38	86.95	426.18
4	6	3	5	4	4	5	6	6	1	43	86.99	422.70
4	6	3	5	4	5	5	6	6	1	41	87.29	376.95
4	6	2	1	4	5	4	6	5	5	33	87.44	345.47
4	6	1	4	2	5	3	6	3	6	17	87.61	330.93
4	6	1	5	2	5	3	6	3	4	13	87.63	313.50
3	6	5	4	1	5	4	5	1	3	21	87.68	303.13
4	6	4	0	1	5	5	6	0	2	10	88.40	175.97
4	6	3	0	1	5	5	5	1	3	6	88.50	173.18
5	2	1	3	3	0	3	0	6	2	90	88.79	156.86
4	6	4	0	3	5	5	4	0	2	8	88.82	138.81
4	6	4	0	3	5	5	3	0	1	5	89.15	128.64
3	6	3	1	1	5	6	3	0	1	1	89.20	125.56
4	6	4	0	3	5	5	3	0	1	5	89.31	114.61
4	6	4	0	1	5	5	1	0	1	1	89.39	105.53
4	6	4	0	1	5	6	0	0	1	3	89.42	97.15
4	6	4	0	1	5	6	1	0	1	3	89.54	93.40
3	2	2	4	2	0	1	0	5	3	84	90.22	45.18
4	2	1	3	1	0	1	1	6	3	74	90.32	41.68
4	2	1	5	1	1	1	1	5	5	89	90.42	34.21
4	2	1	4	1	0	1	0	5	3	88	90.44	23.00
4	2	1	4	1	0	1	0	4	3	80	90.69	13.43
5	2	3	3	2	0	1	0	1	3	81	90.75	10.85
4	2	1	4	1	0	1	0	4	3	85	90.78	8.55
4	2	1	4	0	0	1	1	2	4	88	90.96	6.28

As seen in Table 4, 4 and 6 occur most frequently among all optimization results for n1 and n2. The highest number of occurrences in d1, d2, d3, and d4 are 1, 4, 1, and 5, respectively. The highest number of occurrences in t1, t2, t3, and t4 are 5, 6, 6, and 1, respectively. The optimal value of φ is close to 90°. Among all variables, such as the number, width, and thickness of ribs and φ, the thickness of the ribs has the highest value and the most obvious effect on noise reduction. Therefore, the primary consideration should be to increase the thickness of the ribs. Figure 17 clearly illustrates that at point A (n ₁ = 4, n ₂ = 5, d ₁ = 5 mm, d ₂ = 3 mm, d ₃ = 5 mm, d ₄ = 5 mm, t ₁ = 6 mm, t ₂ = 4 mm, t ₃ = 6 mm, t ₄ = 6 mm, φ = 78°), PSPL is reduced to the minimum, so point A can be regarded as the optimal layout condition for the most obvious noise reduction of the gearbox after adding variable ribs. R_m drops to its lowest at point B (n ₁ = 4, n ₂ = 2, d ₁ = 1 mm, d ₂ = 4 mm, d ₃ = 0 mm, d ₄ = 0 mm, t ₁ = 1 mm, t ₂ = 1 mm, t ₃ = 2 mm, t ₄ = 4 mm, φ = 88°). Thus, point B is considered to be the optimal layout condition for variable ribs with the lowest mass. In order to verify the rationality of the theoretical points, the parameters of the variable ribs model corresponding to points A and B in the multi-objective optimization results are substituted into the calculation flow of the gearbox radiation noise, and the obtained results are compared with points A and B. The comparison results are provided in Table 5. In the corresponding output variables of the point A, the PSPL and R_m optimization results and simulation values errors are 0.52% and +1.8 g, respectively. In the output variables corresponding to point B, the theoretical and calculated error values of PSPL and R_m are 0.67% and +4.82 g, respectively. Figure 18 shows the spectral comparison of sound pressure levels of the improved front and rear gearboxes, where the noise reduction amplitude of the ideal point A is very obvious compared with that before adding ribs, while the noise reduction amplitude of ideal point B is very small. Researchers can choose the appropriate variable ribs layout according to the specific requirements of the R_m. For this purpose, the relationship between PSPL and R_m is deduced according to the Pareto front as follows

$$R_m=177873.68-3890.95\left(PSPL\right)+22.28{\left(PSPL\right)}^2$$

(19)

The region obtained by PACA method is used as the base area of the rib to establish the standard rib. Adjust the thickness of the standard rib so that the mass of the standard rib is the same as the mass of the VRSM at the optimal point A. The PSPL of the standard rib and VRSM at the field point A is shown in Table 6. At the same mass, the noise reduction effect of VRSM is 12% higher than that of standard rib.

Table 5.

The calculated ideal points are compared with the corresponding experimental data.

Design variable point	Output variable	Estimated	Experimental	Deviation
A, B
A	PSPL	86.49	86.95	0.53%
	Rm	521.49	523.29	+1.8 g
B	PSPL	90.96	91.58	0.67%
	Rm	6.28	11.1	+4.82 g

Figure 18.

radiation noise spectrum of the original structure and each scheme model at field point A.

Table 6.

The PSPL at field point A.

	Standard rib, dB	VRSM, dB
PSPL	87.48	86.95

Conclusions

In this paper, a two-stage gearbox radiated noise model is established, and the response surface method and NSGA-II can modify the variables in the VRSM to reduce the noise of the gearbox. The conclusions of this paper are as follows:

1. ATV and MAC improve the accuracy of the PACA method. The ATV-MAC-PACA method can also be used to add damping materials, shock absorbers, and other mechanisms.

2. The thickness of the ribs is the key variable for noise reduction. The PSPL after optimization is decreased by 4.88 dB at frequency of 3130 Hz. At the same mass, the noise reduction effect of VRSM is 12% higher than that of standard rib.

3. The VRSM should be determined according to different field points and different excitation frequencies. A disadvantage of this method is that acoustic contribution analysis must be performed again after changing the working conditions and considering the field points to achieve a good noise reduction effect.