Evaluation and Optimization of the Oil Jet Lubrication Performance for Orthogonal Face Gear Drive: Modelling, Simulation and Experimental Validation

: The oil jet lubrication performance of a high-speed and heavy-load gear drive is signiﬁcantly inﬂuenced and determined by the oil jet nozzle layout, as there is extremely limited meshing clearance for the impinging oil stream and an inevitable blocking e ﬀ ect by the rotating gears. A novel mathematical model for calculating the impingement depth of lubrication oil jetting on an orthogonal face gear surface has been developed based on meshing face gear theory and the oil jet lubrication process, and this model contains comprehensive design parameters for the jet nozzle layout and face gear pair. Computational ﬂuid dynamic (CFD) numerical simulations for the oil jet lubrication of an orthogonal face gear pair under di ﬀ erent nozzle layout parameters show that a greater mathematically calculated jet impingement depth results in a greater oil volume fraction and oil pressure distribution. The inﬂuences of the jet nozzle layout parameters on the lubrication performance have been analyzed and optimized. The relationship between the measured tooth surface temperature from the experiments and the corresponding calculated impingement depth shows that a lower temperature appears in a situation with a greater impingement depth. Good agreement between the mathematical model with the numerical simulation and the experiment validates the e ﬀ ectiveness and accuracy of the method for evaluating the face gear oil jet lubrication performance when using the impingement depth mathematical model.


Introduction
Large amounts of energy will be lost in overcoming the internal friction between the meshing gears, and this will seriously reduce the mechanical efficiency and performance of the transmission system, especially for aeronautical gears under high-speed and heavy-load operation conditions. Face gear drives are a new type of transmission in which the involute cylindrical gear meshes with the face gear. Face gear drives have been widely used in aeronautical transmission systems, as they offer the advantages of a large contact ratio, good power splitting effect, compact structure, insensitivity to installation errors, etc. [1][2][3][4]. Due to the high-speed and heavy-load conditions in aeronautical transmission applications, as well as the space limitations for installation and working, face gear drives will inevitably produce a large amount of heat and result in energy loss. An effective way to solve these problems is to lubricate the face gear drivers with an oil jet. If the lubrication is insufficient, gear meshing will occur under starved or even dry operation and ultimately lead to tooth surface scuffing, pitting and failure damage. Therefore, it is of great importance to investigate the lubrication performance of face gear drives and optimize the jet nozzle layout.
Currently, three primary methods have been developed to investigate the oil jet lubrication performance of transmission gears. One method is to use numerical simulations or experiments. For example, hydrodynamic models of spur gears were established to investigate the variation of oil volume and pressure through CFD simulation to optimize the jet parameters [5,6]. High-speed infrared testers and infrared sensors have been adopted to measure the instantaneous tooth temperature to study the influence of the speed, load and jet location on the lubrication performance [7,8]. Ouyang et al. [9] proposed a novel model of a spur gear pair based on friction dynamics theory and studied the lubrication performance in the high-speed condition by CFD simulation. Massini et al. [10] exploited a novel rotating test rig to analyze oil jet lubrication through high-speed visualizations. Zhao et al. [11] proposed a simplified model based on CFD method for investigating the ripple source of gear pumps. Moreover, to analyze the transient temperature behavior of spiral bevel gears, Gan et al. [12] proposed a method combining the mixed elastohydrodynamic lubrication with the finite element method and conducted the thermal analysis using the CFD method. Another method is to calculate the oil film thickness and film pressure based on a theoretical equation. For example, Poulios et al. [13] used a quadratic B-spline basis function to obtain a more accurate oil film thickness and pressure field. Furthermore, Wang et al. [14] established formulas for the contact path of the face gear loaded tooth contact analysis (LTCA) and the dimensionless face gear isothermal elastohydrodynamic lubrication (EHL) to obtain the variation in the lubricating oil film thickness and oil film pressure. Following Wang's work, Liu et al. [15] established a comprehensive mechanical efficiency model of a helical gear pair and evaluated the tribological performance in terms of the average film thickness, the friction coefficient, the mechanical power loss, etc. A deterministic model combining contact mechanics with pure extrusion lubrication was established by Xu et al. [16], so as to investigate the extrusion effect of the oil film under the transient boundary lubrication condition. Ahmed et al. [17] presented an automatic locally adaptive finite element solver for fully coupled EHL point contact problems to significantly improve the accuracy of the elastic deformation solution. A calculating model for the full tooth surface flash temperature distribution for face gear drives was established to optimize the anti-scuffing tooth modification schedule for face gear drives based on the Blok flash temperature formula [18], and this model improved the anti-scuffing capacity of face gear drives. Thiagarajan et al. [19] proposed a mixed-thermoelastohydrodynamic (TEHD) model for investigating the effects of surface roughness, fluid structure and thermal interaction on the mixed lubrication in the regions of low film thickness. The final method is to analyze the impingement characteristic of lubricating oil on the tooth surface. Previous research on this aspect is as follows, but only spur and helical gears were studied: Akin et al. [20,21] deduced the formulas for the oil jet impingement depth on spur and helical gears and studied the influence of the offset distance, jet angle and jet velocity on the impingement depth. Ambrose et al. [22] used the lattice Boltzmann method (LBM) to study the oil impingement on a spur gear and analyzed the effect of the oil feed delivery rate on the spreading of the oil jet on the tooth and the splashing profiles.
However, as far as the relevant studies are concerned, it could take a long time to analyze and evaluate the oil jet lubrication performance of gears using CFD simulations or experiments. There are still a few assumptions in solving the equations of the oil film thickness and pressure, so the results are not accurate enough. In addition, until now, only the oil jet impingement depth on the tooth surface of the spur or helical gear pair has been studied. Furthermore, there are relatively few studies on the influence of the jet nozzle layout parameters on the lubrication performance and the optimization of the nozzle layout. Therefore, first, a novel mathematical model for the impingement depth calculation for the orthogonal face gear is established based on face gear meshing theory [23], the involute function of the spur gear [24], and the spatial positional relationship of the face gear tooth surface. Subsequently, the oil volume fraction and oil pressure are obtained by CFD simulation to validate the effectiveness   The derivation steps for the mathematical model of the impingement depth are as follows: 1.
Based on the meshing theory of a face gear, the space coordinate system of the orthogonal face gear is defined.

2.
The initial position parameters (t 0 = 0) of the pinion and the face gear have been established separately using an involute function of the spur gear and the spatial position of the face gear tooth surface. After a period of time (t 1 = t), the geometrical positions of the face gear pair and jet stream are calculated.

3.
Since the meshing time of face gear is equal to the time of jet stream reaching the tooth surface, the impingement depth on the pinion and the face gear can be obtained by solving the simultaneous equations.

Definitions of Nozzle Layout Parameters
For the orthogonal face gear drive, the tooth contact type is point contact, and the meshing process can be regarded as the involute spur gear meshing with the rack at different shaft cross sections [25]. The space coordinate system is set up as shown in Figure 2. At this time, there is a circle on the pinion pitch cylinder rolling purely with a circle along the direction of the face gear radius. Point O 0 is the gear coordinate origin, representing the center of the gear-locating surface. The x-axis represents the intersection line between the symmetry surface and the locating surface of the face gear. The z-axis represents the axis of the face gear. The y-axis can be determined by the right-hand rule. Point O represents the intersection point of the face gear axis and the pinion axis. Furthermore, points O p and O g represent the centers of the surfaces through impingement points and perpendicular to the pinion axis and the face gear axis, respectively.
As shown, the impingement depth is directly relevant to the jet orientation parameters x H , y L , z V , θ and φ. x H , y L , z V denote the nozzle exit position; the parameter θ denotes the angle between the jet stream and the z-axis, which is always restricted to π/2 ≤ θ ≤ π/2, and the parameter φ denotes the angle between the jet stream projection line on the xO g y plane and the x-axis. Since the pinion is an involute spur gear with a symmetrical structure, this paper focuses on the case that the jet stream is parallel to the pinion shaft cross section, that is φ = π/2. Additionally, Σ denotes the shaft angle of the face gear pair, and A g is the distance from the pinion axis to the xO 0 y plane.  Figure 2. Description of nozzle layout position and orientation parameters for the orthogonal face gear.

129
As shown, the impingement depth is directly relevant to the jet orientation parameters between the jet stream and the z-axis, which is always restricted to   / 2 / 2 π θ π , and the 132 parameter  denotes the angle between the jet stream projection line on the g xO y plane and the 133 x-axis. Since the pinion is an involute spur gear with a symmetrical structure, this paper focuses on 134 the case that the jet stream is parallel to the pinion shaft cross section, that is   / 2 π . Additionally, 135  denotes the shaft angle of the face gear pair, and g A is the distance from the pinion axis to the 136 0 xO y plane.

138
The oil jet orientation parameters , , , x y z θ are known initially, and the distance g A and

139
the oil jet velocity j V are also given. Thus, the process for calculating the impingement depth of the 140 pinion is as follows.

141
At the initial moment (t0 = 0): the position parameters of the face gear pair and the jet stream are 142 as illustrated in Figure 3.

Mathematical Model for the Pinion
The oil jet orientation parameters x H , y L , z V , θ are known initially, and the distance A g and the oil jet velocity V j are also given. Thus, the process for calculating the impingement depth of the pinion is as follows.
At the initial moment (t 0 = 0): the position parameters of the face gear pair and the jet stream are as illustrated in Figure 3.

129
As shown, the impingement depth is directly relevant to the jet orientation parameters 130 , , , x y z θ and  .
, , x y z denote the nozzle exit position; the parameter  denotes the angle 131 between the jet stream and the z-axis, which is always restricted to   / 2 / 2 π θ π , and the 132 parameter  denotes the angle between the jet stream projection line on the g xO y plane and the 133 x-axis. Since the pinion is an involute spur gear with a symmetrical structure, this paper focuses on 134 the case that the jet stream is parallel to the pinion shaft cross section, that is   / 2 π . Additionally, 135  denotes the shaft angle of the face gear pair, and g A is the distance from the pinion axis to the      Figure 3, the projection length of the line from the point O g to the face gear addendum on the surface O g can be expressed as

According to
Substituting the parameter φ = π/2, Equation (1) can be written as where z V denotes the z coordinate value of the point where the jet stream passes through the face gear addendum. Let H be the distance from the face gear addendum to the xO 0 y plane; thus, Energies 2019, 12, 1935 5 of 23 It is known that the jet streamline is parallel to the pinion shaft cross-section; that is, the jet stream is perpendicular to the x-axis. According to Figure It is known that the jet streamline is parallel to the pinion shaft cross-section; that is, the jet 151 stream is perpendicular to the x-axis. According to Figure

156
According to the geometric relationship in Figure 4a, the following equation is obtained: Equations (2), (3) and (4) can be rewritten in a combined form:  According to the geometric relationship in Figure 4a, the following equation is obtained: Equations (2), (3) and (4) can be rewritten in a combined form: where x H and y L represent the x and y coordinate values of the point where the jet stream passes through the face gear addendum, respectively. The following equations can be obtained from Figure 4b: By solving Equations (6), (7), (8) and (9), r α and r z can be expressed as where r z denotes the distance from the intersection point of the jet stream with line O p O 0 to the point O p ; λ denotes the angle between the line r z and the line r α . θ p1 and θ g1 represent the position parameters of the pinion and the face gear at the initial time (t 0 = 0), and their relationship can be deduced by the rotation angle relationship between the face gear and the gear shaper cutter during the machining process. The rotation angles ϕ g and ϕ s of the face gear and the gear shaper cutter satisfy the transmission ratio [26]: where N s and N g denote the numbers of teeth on the shaper and the face gear, respectively. To avoid the interference between the shaper cutter and the edge of the face gear, the face gear drive is changed from an instantaneous line contact to a point contact drive. In this case, the number of teeth on the pinion will be 1-3 teeth less than on the gear shaper cutter [27][28][29]. Figure 5 illustrates an imaginary internal tangency of the shaper cutter and the pinion [30,31].
gear and the gear shaper cutter during the machining process. The rotation angles g φ and s φ of 167 the face gear and the gear shaper cutter satisfy the transmission ratio [26]: where s N and g N denote the numbers of teeth on the shaper and the face gear, respectively.

169
To avoid the interference between the shaper cutter and the edge of the face gear, the face gear 170 drive is changed from an instantaneous line contact to a point contact drive. In this case, the number 171 of teeth on the pinion will be 1-3 teeth less than on the gear shaper cutter [27][28][29].  O p and O s denote the centers of the pinion and the shaper cutter shaft sections, respectively. Their rotation angles satisfy the following equation: Combining Equation (12) with Equation (14), the relationship between the rotation angle ϕ p of the pinion and ϕ g of the face gear can be expressed as: According to Figure 3 and Equation (15), the initial position parameter θ p1 of the pinion is expressed by where invα p1 is the involute function of the spur gear, representing the spread angle at the intersection point between the pinion pitch circle and the involute; and the pressure angle α p1 on pitch circle of the pinion is expressed as Obviously, from Figure 4, the initial position parameter θ g1 of the face gear is At the moment (t 1 = t), the position parameters of the face gear pair and the jet stream are illustrated in Figure 6.
At the moment (  t t 1 ), the position parameters of the face gear pair and the jet stream are 186 illustrated in Figure 6. Figure 6. Illustration of the impingement depth on the pinion (t1 = t).
189 Figure 6. Illustration of the impingement depth on the pinion (t 1 = t).
As the flowing time of the jet steam is equal to the rotation time of the pinion, which is rotating from the angle θ p1 at the initial time t 0 to the angle θ p2 at time t 1 , the impingement depth d p can be calculated as where where L p represents the impingement distance, and r a represents the addendum radius of the pinion.
where ω p is the angular velocity of the pinion. As can be seen in Figure 6, the position parameter of the pinion at t 1 = t is where invα p2 is the involute function of the spur gear, denoting the spread angle at the impingement point M of the involute shown in Figure 7; α p2 denotes the pressure angle at the impingement point on the volute. Their relationship can be obtained by

195
where p ω is the angular velocity of the pinion.

196
As can be seen in Figure 6, the position parameter of the pinion at 1  t t is   The following equations can be obtained from Figures 6 and 7: where r i denotes the radius of the pinion at the impingement point, and r b denotes the base circle radius of the pinion. Equation (21) can be rewritten as Figure 8 illustrates the projection of the jet stream on the surface O p ; according to the geometric relationship, the following equation can be obtained: Energies 2019, 12, x FOR PEER REVIEW 9 of 26 where i r denotes the radius of the pinion at the impingement point, and b r denotes the base circle 204 radius of the pinion.

205
Equation (21) can be rewritten as This can be simplified by substituting Equation (7) into Equation (27): where ' Δh denotes to the projection of Δh on the surface p O . Hence, the relationship between

'
Δh and Δh is Moreover, Equation (19) can be reformulated as By substituting Equations (3), (10), (11), (20) and (26) into Equations (28) and (29), the 215 mathematical model of the impingement depth on the pinion can be established as follows: This can be simplified by substituting Equation (7) into Equation (27): where ∆h denotes to the projection of ∆h on the surface O p . Hence, the relationship between ∆h and ∆h is ∆h = ∆h Moreover, Equation (19) can be reformulated as By substituting Equations (3), (10), (11), (20) and (26) into Equations (28) and (29), the mathematical model of the impingement depth on the pinion can be established as follows: where

Mathematical Model for the Face Gear
The deduction process of the mathematical model of impingement depth for the face gear is approximately the same as that for the pinion; therefore, only the main derivation steps are presented in this paper. The position parameters of the face gear pair and the jet stream at the initial moment (t 0 = 0 ) and the moment (t 1 = t) are as illustrated in Figure 9a,b, respectively.

218
The deduction process of the mathematical model of impingement depth for the face gear is 219 approximately the same as that for the pinion; therefore, only the main derivation steps are 220 presented in this paper. The position parameters of the face gear pair and the jet stream at the initial 221 moment ( 0 0  t ) and the moment ( 1  t t ) are as illustrated in Figure 9a,b, respectively.

224
Additionally, projections of the jet stream on different surfaces at different moments are shown 225 in Figure 10.

228
Similarly, the position parameters of the pinion and the face gear at the initial time ( 0 0  t ) and Additionally, projections of the jet stream on different surfaces at different moments are shown in Figure 10.

218
The deduction process of the mathematical model of impingement depth for the face gear is 219 approximately the same as that for the pinion; therefore, only the main derivation steps are 220 presented in this paper. The position parameters of the face gear pair and the jet stream at the initial 221 moment ( 0 0  t ) and the moment ( 1  t t ) are as illustrated in Figure 9a,b, respectively.

224
Additionally, projections of the jet stream on different surfaces at different moments are shown 225 in Figure 10.  Similarly, the position parameters of the pinion and the face gear at the initial time (t 0 = 0) and the face gear at the moment (t 1 = t) are denoted by θ p4 , θ g2 and θ g3 , respectively, which can be calculated using the following expressions where A g denotes the distance from the pinion axis to the xO 0 y plane; Z V denotes the z coordinate value of the point where the jet stream passes through the pinion addendum; α ak4 and α k4 denote the pressure angles at the intersection point of the addendum circle and the pitch circle with the involute, respectively; L g denotes the impingement distance; and X, Y, Z, θ , φ are the jet orientation parameters. According to the definition of the impingement depth, there is an angle between the line from the impingement point to the face gear addendum and the z-axis. The impingement depth on the face gear is assumed to be equal to its projection on the surface O g . Hence, similar to Equation (19), the impingement depth d g on the face gear can be calculated as where R i is the inner radius of the face gear, L g is the impingement distance, which satisfies the expression L g = Y L − Z V − Z V tan θ − ∆h , while ∆h is defined as ∆h = ∆h sin(π − θ ), which is the projection of ∆h on the surface O g , as shown in Figure 10a,c. Therefore, the mathematical model of the impingement depth d g on the face gear can be established as where where V j denotes the jet velocity; ω g denotes the angular velocity of the face gear.

CFD Numerical Simulations
In the case of the given parameters, such as the gear structures, parameters and working conditions, the oil jet lubrication performance can be judged by the oil volume fraction and oil pressure distribution in the meshing area [32][33][34], and a greater oil volume fraction and pressure is commonly recognized as providing better jet lubrication performance. The Fluent program based on the CFD method was adopted to simulate the distribution of oil-air, two-phase flow in the meshing area.

CFD Model and Main Settings
For the face gear drive, the generated heat distribution to the pinion is more than the heat distributed to the face gear at the same time; furthermore, the heat transfer coefficient of the pinion is less than that of the face gear, so the pinion is more prone to damage during the meshing process. Hence, this paper focuses on the jet lubrication performance of the pinion. The separation method [35,36] is used to moderately increase the distance between the pinion and the face gear. The parameters of the face gear pair are listed in Table 1. To ensure the reliability and accuracy of the transient simulation results, grid independence tests were preferentially performed. The distributions of the oil volume fraction and oil pressure at 0.005 s on a specific plane parallel to the coordinate plane xOz were taken for the grid independence tests. Table 2 lists the number of mesh elements in six cases. In the simulations, the number of the mesh elements was ensured as the only independent variable. As seen from simulation results illustrated in Figure 11, the trends of the oil pressure and the oil volume fraction become insignificant with the increasing number of mesh elements. When the mesh elements reach Case 4, the oil pressure and oil volume fraction can be considered to be stable. Therefore, the total mesh elements in all subsequent simulations were controlled at approximately 3.2 million.   According to the grid independence test results, the computation domain was divided into tetrahedral meshes with approximately 3.2 million mesh elements, where the defeature size and curvature size of the engaged teeth were 0.1 mm and 0.4 mm, respectively. Figure 12 presents the final mesh model, where the maximum skewness and the mesh quality are 0.842 and 0.283, respectively. During the progress of oil jet lubrication, the fluid in the gearbox changes from the initially only air into oil-air, two-phase mixture flow. The viscosity and density of the lubrication oil are set as 1.98 mm 2 /s and 959.4 kg/m 3 , respectively.
It is known that the two-phase flow distribution is time-varying, so the pressure-solver and transient-state were adopted in the paper. The VOF multiphase flow model was applied to simulate the oil-air flow in the meshing area [37,38]. Considering the swirling effects generated by high-speed rotating gears, the RNG k − ε turbulence model with higher precision was used. The lower circular surface of the nozzle was set as the velocity inlet, and the velocity and hydraulic diameter were set to 50 m/s and 0.0014 m, respectively. Moreover, non-slip boundary conditions were used for all walls. The dynamic mesh was adopted to simulate the real rotation of the face gear drive. Figure 13 demonstrates the distributions of the oil volume fraction at different moments. To obtain more accurate results, the second-order windward methods of the momentum, the turbulent kinetic energy, and the turbulent dissipation rate were adopted for spatial discretization. Furthermore, the standard SIMPLE algorithm was used for pressure-velocity coupling.

284
It is known that the two-phase flow distribution is time-varying, so the pressure-solver and 285 transient-state were adopted in the paper. The VOF multiphase flow model was applied to simulate 286 the oil-air flow in the meshing area [37,38]. Considering the swirling effects generated by high-speed rotating gears, the RNG k ε  turbulence model with higher precision was used. The lower circular 288 surface of the nozzle was set as the velocity inlet, and the velocity and hydraulic diameter were set to 289 50 m/s and 0.0014 m, respectively. Moreover, non-slip boundary conditions were used for all walls.

290
The dynamic mesh was adopted to simulate the real rotation of the face gear drive. Figure 13 291 demonstrates the distributions of the oil volume fraction at different moments. To obtain more 292 accurate results, the second-order windward methods of the momentum, the turbulent kinetic 293 energy, and the turbulent dissipation rate were adopted for spatial discretization. Furthermore, the 294 standard SIMPLE algorithm was used for pressure-velocity coupling.  For the determination of the time step size, an initial value of 5 × 10 −6 s was obtained by the minimum edge length of the nozzle divided by the fluid velocity. With respect to the solving efficiency and accuracy, it was necessary to perform time step independence tests in the six cases of 5 × 10 −6 s, 2 × 10 −5 s, 3.5 × 10 −5 s, 5 × 10 −5 s, 6.5 × 10 −5 s and 8 × 10 −5 s. Similar to the grid independence tests, the distributions of the oil pressure and oil volume fraction on the specific plane were used for evaluation. As shown in Figure 14, the changes of the oil pressure and oil volume fraction were within an allowable error range before 5 × 10 −5 s. Subsequently, as the time step increased, the oil pressure and oil volume fraction were significantly affected. Thus, the time step was determined to be 5 × 10 −5 s. , , , , and . Similar to the grid 301 independence tests, the distributions of the oil pressure and oil volume fraction on the specific plane 302 were used for evaluation. As shown in Figure 14, the changes of the oil pressure and oil volume 303 fraction were within an allowable error range before

308
The simulation results should be obtained in a convergence state. In this paper, whether the  The simulation results should be obtained in a convergence state. In this paper, whether the calculation reached convergence was mainly judged by the following three aspects: Firstly, all residuals were set to 0.001, including the Continuity Equation, the Momentum Equations, the Turbulent Kinetic Energy Equation and the Turbulent Kinetic Energy Dissipation Rate Equation. Secondly, when the residual values all dropped below 0.001, the residual was adjusted to continue the simulations. Meanwhile, the oil pressure and oil volume fraction in the meshing area were monitored. Finally, when the monitored variables almost exhibited no further change with the increase of the iteration step, the mass and momentum data would be conserved. If the error was within the allowable range, the simulation could be considered to reach the convergence.
According to Equations (31) and (32), the impingement depth on the pinion surface is directly related to the nozzle layout parameters x H , y L , z V , θ. Since the nozzle layout parameters are restricted by the gearbox space, the initial nozzle position is determined as (70, 52, 32, 105 • ). Each parameter of the initial position is investigated to verify the impingement depth mathematical model and to optimize the nozzle layout.

Verification and Optimization of Parameter x H
According to the initial position, four groups of the parameter x H are set preferentially as 65 mm, 70 mm, 75 mm and 80 mm, respectively, while the other three parameters remain the same. Based on the mathematical model established above, the values of the impingement depths corresponding to the four groups are obtained by the implicit function in the program MATLAB, as given in Table 3. It can be seen that a maximum impingement depth can be obtained with x H = 75 mm. Table 3. Values of parameter x H and the corresponding calculated impingement depths.

Number
x H (mm) Accordingly, the CFD models using the parameters provided in Table 3 are established. Figure 15 illustrates the distributions of the oil volume fractions in the meshing area in the convergence state of the calculation. Obviously, Figure 15c exhibits superior characteristics on the pinion surface, as the oil distribution is more concentrated and uniform with relatively little oil flowing out. Accordingly, the CFD models using the parameters provided in Table 3 are established. Figure   331 15 illustrates the distributions of the oil volume fractions in the meshing area in the convergence 332 state of the calculation. Obviously, Figure 15c exhibits superior characteristics on the pinion surface, 333 as the oil distribution is more concentrated and uniform with relatively little oil flowing out.   In the meshing cycle, the trends of the oil volume fraction and oil pressure generally increase, and both can reach the maximum values in the case of x H = 75 mm; consequently, a better lubrication performance can be achieved. Mutual verification was achieved between the mathematical model and the numerical simulation. The oil volume fractions and pressures were observed to fluctuate at different moments; this result was due to the fact that the jet stream was blocked by the rotating gear and could not enter the meshing area smoothly. Meshing time/(1e-03s) Furthermore, the parameter 75 mm H x  may be a point close to the optimal or peak value, and the optimal value could be more accurate, as the solution interval is smaller. Figure 17 illustrates  Furthermore, the parameter x H = 75 mm may be a point close to the optimal or peak value, and the optimal value could be more accurate, as the solution interval is smaller. Figure 17 illustrates the calculated impingement depth within the range of x H ∈ [70 mm, 80 mm] with an interval length of only 1 mm. The maximum impingement depth is observed to appear when x H = 71 mm, which is regarded as the optimal value for generating a better lubrication performance.

Verification and Optimization of Parameter y L
Similarly, four groups of the parameter y L are set preferentially as 37 mm, 42 mm, 47 mm and 52 mm, respectively, while the parameter x H is set at the optimal value of 71 mm, as determined previously, while the other two parameters remain the same. The corresponding impingement depths are calculated as shown in Table 4, where negative values indicate that the jet stream could not reach the tooth surface. The maximum impingement depth could be obtained when y L = 42 mm. Furthermore, the parameter 75 mm H x  may be a point close to the optimal or peak value, 345 and the optimal value could be more accurate, as the solution interval is smaller. Figure 17 illustrates   Accordingly, CFD models were established with the different parameters provided in Table 4. Figure 18 illustrates the distribution of the oil volume fractions in the meshing area; obviously, Figure 18b,c exhibits superior distribution characteristics. Accordingly, CFD models were established with the different parameters provided in Table 4.    Figure 19 further presents that both the oil volume fraction and the oil pressure can achieve maximum values in the case of y L = 42 mm. A good agreement between the mathematical model and the CFD numerical simulation can be achieved. However, moving the nozzle closer to the meshing point does not improve the lubrication performance, as can be seen from the figure; when the y L decreases to 37 mm, the oil volume fraction and oil pressure both decrease significantly. This condition is due to the noticeable accumulation of oil in the meshing area, leading to more heat accumulation and a poor cooling effect.    x  , as given in Table 5. Furthermore, Figure 20 illustrates the calculated impingement depth within a range of y L ∈ [37 mm, 47 mm] with a smaller interval length of 1 mm. The maximum impingement depth occurs when y L = 40 mm, which is regarded as an optimal value.   x  , as given in Table 5.

Verification and Optimization of Parameter z V
When the value of parameter z V changes, while the jet angle parameter θ remains unchanged, the jet stream may be prevented by the gear teeth, leading to unexpected lubrication performance. According to the above research, the parameter z V is set as 27 mm, 32 mm, 37 mm and 42 mm in the case of x H = 71 mm and x H = 40 mm, as given in Table 5. Multiple calculated impingement depths appear to be negative values, indicating that the jet stream could not successfully jet onto the tooth surface. Therefore, there is no need to compare the oil volume fraction and oil pressure distributions, since only the model number 2 has a positive impingement depth. Using a method similar to that employed in Sections 3.2 and 3.3, an appreciable maximum impingement depth when z V = 32 mm could be obtained, as seen in Figure 21.      Table 6. The maximum impingement depth 395 could be obtained when the jet angle was  105 .

Verification and Optimization of Jet Angle θ
Based on the above research, the maximum impingement depth and the best lubrication performance could be obtained in the case of x H = 71 mm, y L = 40 mm and z V = 32 mm. The jet angles were set preferentially as 100 • , 102.5 • , 105 • , and 107.5 • , respectively, and the corresponding impingement depths are presented in Table 6. The maximum impingement depth could be obtained when the jet angle was 105 • . Table 6. Different jet angles θ and the corresponding impingement depth.

Number
x H (mm)    Furthermore, Figure 24 illustrates that the calculated impingement depth within the range of θ ∈ [102.5 • , 107.5 • ] with the smaller interval length of only 0.5 degrees. It can be seen that the maximum impingement depth appears when θ = 105 • , which is regarded as the optimal value.  417 Figure 25 shows the whole experimental system, including a gearbox equipped with a nozzle layout 418 adjustment device, a lubrication oil supply system, a driving motor, a loading motor, torque 419 transducers, a thermal infrared imager, and a high-speed camera. The main experimental 420 parameters are shown in Table 7.  From the analysis and optimization above, an optimal combination of jet nozzle layout parameters that can generate a greater impingement depth, a greater oil volume fraction and a greater oil pressure for this face gear pair are determined and exhibit a better lubrication performance.

Experiments
To validate the feasibility and reliability of the theoretical method for evaluating the lubrication performance using the impingement depth mathematical model, oil jet lubrication experiments were performed under different nozzle layouts, which was an important prerequisite for deriving the novel impingement depth mathematical model for face gears and conducting numerical simulations. Figure 25 shows the whole experimental system, including a gearbox equipped with a nozzle layout adjustment device, a lubrication oil supply system, a driving motor, a loading motor, torque transducers, a thermal infrared imager, and a high-speed camera. The main experimental parameters are shown in Table 7.  417 Figure 25 shows the whole experimental system, including a gearbox equipped with a nozzle layout 418 adjustment device, a lubrication oil supply system, a driving motor, a loading motor, torque 419 transducers, a thermal infrared imager, and a high-speed camera. The main experimental 420 parameters are shown in Table 7.  Different nozzle layout parameters, including the jet angle and the jet offset distance, as shown in Table 8, were designed and arranged in the experiments. Figure 26 shows the jet oil states photographed by a high-speed camera. Furthermore, a thermal infrared imager was used to capture pictures of the tooth surface temperatures, as listed in Table 9. The measured temperatures for one model test are shown in Figure 27. Different nozzle layout parameters, including the jet angle and the jet offset distance, as shown 425 in Table 8, were designed and arranged in the experiments. Figure 26 shows the jet oil states 426 photographed by a high-speed camera. Furthermore, a thermal infrared imager was used to capture 427 pictures of the tooth surface temperatures, as listed in Table 9. The measured temperatures for one 428 model test are shown in Figure 27.     437 Figure 28 and Table 9 Table 9 summarize the measured tooth surface temperatures corresponding to the calculated impingement depths. 437 Figure 28 and Table 9 summarize the measured tooth surface temperatures corresponding to the 438 calculated impingement depths.  As shown, with an increase in the calculated impingement depth, the tooth surface temperature 444 is reduced, and consequently, better lubrication and cooling performance can be achieved, which 445 can validate the effectiveness and accuracy of the method used to evaluate the oil jet lubrication 446 performance using the impingement depth mathematical model. 447 Figure 28. Trends of the tooth surface temperatures with the corresponding calculated impingement depths.

429
As shown, with an increase in the calculated impingement depth, the tooth surface temperature is reduced, and consequently, better lubrication and cooling performance can be achieved, which can validate the effectiveness and accuracy of the method used to evaluate the oil jet lubrication performance using the impingement depth mathematical model.