螺旋锥齿轮齿面加载性能多目标优化
王琪 , 周驰 , 桂良进 , 范子杰     
清华大学 汽车工程系, 汽车安全与节能国家重点实验室, 北京 100084
摘要:针对螺旋锥齿轮齿面设计需进行多次试算的问题,提出了齿面加载性能多目标优化问题数学模型的建立及求解方法。建立了包含最大接触应力、加载传动误差、加载接触区范围和齿根弯曲应力的齿面加载性能多目标优化问题的数学模型;建立了考虑齿根弯曲应力计算的加载接触分析模型,用于计算优化问题中的目标函数和约束函数;采用Kriging代理模型结合多目标遗传算法对该优化问题进行求解;通过对驱动桥螺旋锥齿轮副进行多目标优化,实现了接触区完全位于理想接触区且无边缘接触,最大接触应力降低11.7%,加载传动误差降低27.9%,大小轮最大齿根弯曲应力分别降低2.0%和12.6%。通过对优化后的齿面进行驱动桥加载接触印迹台架试验,对接触区进行了验证,结果显示优化得到的接触区与试验结果吻合良好,验证了优化方法的可行性和正确性。
关键词螺旋锥齿轮    齿面优化    齿面加载接触分析    加载接触性能    
Multi-objective optimization of loaded spiral bevel and hypoid gears
WANG Qi, ZHOU Chi, GUI Liangjin, FAN Zijie     
State Key Laboratory of Automotive Safety and Energy, Department of Automotive Engineering, Tsinghua University, Beijing 100084, China
Abstract: The loaded stresses in the tooth flank of spiral bevel and hypoid gears were optimized without a trial-and-error procedure in a multi-objective optimization model. The multi-objective optimization model was developed considering the maximum contact pressure, load transmission error, load contact pattern and the maximum bending stress. A semi-analytical loaded tooth contact analysis model for the root bending stress was developed to calculate the objectives and constraints. Kriging surrogate models of the objectives and constraints together with a multi-objective genetic algorithm were used to solve the optimization model. The model was then used to optimize a spiral bevel gear pair for a drive axle. The loaded contact pattern was completely in the ideal zone with the maximum contact pressure decreased by 11.7%, the loaded transmission error decreased by 27.9% and the maximum bending stresses in the wheel and pinion decreased by 2.0% and 12.6%. The optimal design was installed in a drive axle for loaded contact tests that showed that the actual loaded contact pattern coincided well with the predicted loaded contact pattern, which verified the accuracy of the optimization approach.
Key words: spiral bevel and hypoid gear     tooth flank optimization     loaded tooth contact analysis     loaded contact performance    

螺旋锥齿轮具有传动平稳、承载能力高等优点,广泛用于汽车和航空领域。这类齿轮要求受载后接触区不出现边缘接触,加载传动误差峰峰值(loaded transmission error, LTE)、最大接触应力和最大齿根弯曲应力较小。为获得性能良好的齿面,对齿面加载性能进行优化设计具有重要意义。

对于齿面设计,Litvin[1]提出了局部综合法,通过控制参考点的局部综合参数实现对齿面的控制,也有学者直接采用齿面的修形面进行齿面设计[2-3]。在这些设计方法的基础上,为了得到理想接触区,Achtmann等[4]确立了空载接触椭圆的理想位置,对空载接触椭圆进行了优化;文[5-7]以局部综合法的部分控制参数为设计变量,研究了各参数对接触迹线的影响,对空载接触迹线进行了优化设计;刘光磊等[8-9]以变性法加工中的变性系数为变量,对空载传动误差曲线对称性进行了优化设计。

在齿面空载性能优化基础上,有学者对加载性能进行了相关研究。Simon[10-11]研究了刀具参数对最大接触应力的影响,并根据其分析结果找到了合适的刀具参数值,但未考虑加载接触区位置、齿根弯曲应力等性能;Artoni等[12]和Gabiccini等[13]以修形面的参数为设计变量,对齿面加载接触区域面积进行了优化,但其优化模型仅考虑了接触区面积,未考虑其他性能指标;Artoni等[14]建立了以加载传动误差为目标并控制接触区位置的齿面优化模型,但其模型中没有考虑齿根弯曲应力。Artoni等[15]在建立的单目标优化模型[14]基础上,建立了以最大接触应力和加载传动误差等为目标、以接触区域为约束的多目标优化模型,通过加权方法将多目标问题转换为单目标问题,并采用不同初始点进行多次单目标优化问题计算得到Pareto最优解集,但其Pareto解集的获得需进行多次单目标优化,且优化模型中也未考虑齿根弯曲应力。

综上可以看出,现有研究对于螺旋锥齿轮齿面加载性能的优化存在如下问题:1)对于齿面加载性能的多目标优化问题研究较少;2)在研究中有部分加载性能未考虑,如齿根弯曲应力。存在这些问题的主要原因在于,进行多目标优化问题求解时一般需要进行多次齿面加载接触分析(loaded tooth contact analysis, LTCA),而LTCA常采用半解析方法计算[16-17],需利用轮齿粗糙的有限元模型计算弯曲变形,但采用粗糙的有限元模型无法准确计算齿根弯曲应力[18],若采用精细的有限元模型则由于整体计算量大,不适于优化过程中多次计算。

针对现有研究中存在的问题,本文首先建立了考虑齿根弯曲应力的螺旋锥齿轮齿面多目标优化问题的数学模型;其次,为计算加载性能指标,在LTCA方法中加入了基于载荷分布和重新加密网格的齿根弯曲应力计算;然后,利用Kriging代理模型结合多目标遗传算法(non-dominated sorting genetic algorithm Ⅱ, NSGA-Ⅱ)进行了多目标优化问题的求解;最后,以驱动桥中的螺旋锥齿轮副为例进行了优化,并取优化后的齿面进行了加工与加载接触印迹台架试验。

1 齿面加载性能多目标优化模型的建立

本文齿面设计采用基于修形面的设计方法[2],对于所研究的格里森(Gleason)收缩齿螺旋锥齿轮,修形面是指小轮点接触齿面与小轮完全共轭齿面的差曲面。为了实现齿面设计,通常会给定修形面控制参数[2],包括:1)空载传动误差曲线二次项系数a,如式(1)所示,用于控制空载时的传动误差;2)预设接触迹线斜率ks,如图 1所示,图 1将小轮修形面上的点放入二维投影面对应的网格中,投影面内以-1~1范围表示,斜率ks可以根据图 1中坐标系计算得到,它用于控制空载接触迹线;3)预设接触椭圆半宽b,如图 1所示,用于控制接触椭圆半宽。

图 1 小轮投影面内接触迹线示意图

另外,由于在齿轮实际啮合时会受到系统变形产生的错位量的影响,因此在设计修形面时考虑了图 2所示的系统错位量[19]。考虑错位量后,接触中心取图 1A点以保证空载传动误差曲线的对称性,最终得到图 3所示的小轮修形面示意图。

$ {\phi _2} - \phi _2^0 = \frac{1}{{{i_{12}}}}\left( {{\phi _1} - \phi _1^0} \right) + a{\left( {{\phi _1} - \phi _1^0} \right)^2}. $ (1)
图 2 螺旋锥齿轮错位量示意图

图 3 小轮修形面示意图

式中:ϕ1ϕ2分别为小轮和大轮啮合时的转角,ϕ10ϕ20分别为小轮和大轮参考点位置转角,i12为理论传动比。a为预设抛物线二次项系数,它通过式(2)可预控空载传动误差(unloaded transmission error, UTE)。式(2)中z1为小轮齿数。

$ {\rm{UTE}} = a{\left( {\frac{{2{\rm{ \mathsf{ π} }}/{z_1}}}{2}} \right)^2}. $ (2)

计算得到修形面后(以Γ表示),利用齿面误差反求方法[20-22]得到实现修形面Γ对应的机床加工参数。为计算修形面Γ对应的加工参数,将Γ看作误差曲面。通常齿面方程位置矢量r和法矢量n

$ \left\{ \begin{array}{l} \mathit{\boldsymbol{r}} = \mathit{\boldsymbol{r}}\left( {\theta ,\varphi ,{\xi _j}} \right)\\ \mathit{\boldsymbol{n}} = \mathit{\boldsymbol{n}}\left( {\theta ,\varphi ,{\xi _j}} \right) \end{array} \right.,j = 1,2, \cdots ,11. $ (3)

其中:θφ分别为加工时的刀盘转角和摇台转角,ξj为小轮的11个加工参数(基本摇台角、垂直轮位、水平轮位、安装根锥角、径向刀位、滚比、床位、刀倾角、刀转角、刀盘半径、刀具齿形角)。

对于齿面上一点i,定义由加工参数微小变动δξj引起该点齿面法向误差为δRi(式(4)),则加工参数变动δξj对齿面的灵敏度向量为s(式(5)),将所有加工参数的灵敏度向量组集成矩阵,即可得到齿面法向误差灵敏度矩阵S(式(6))。

$ \begin{array}{*{20}{c}} {{\rm{ \mathsf{ δ} }}{\mathit{\boldsymbol{R}}_i} = \left( {{\mathit{\boldsymbol{r}}_i}\left( {{\xi _j} + {\rm{ \mathsf{ δ} }}{\xi _j}} \right) - {\mathit{\boldsymbol{r}}_i}\left( {{\xi _j}} \right)} \right) \cdot {\mathit{\boldsymbol{n}}_i}\left( {{\xi _j}} \right),}\\ {j = 1,2, \cdots ,11;} \end{array} $ (4)
$ \mathit{\boldsymbol{s}}\left( {{\xi _j}} \right) = \frac{{\partial {\mathit{\boldsymbol{R}}_i}}}{{\partial {\xi _j}}}{\rm{ \mathsf{ δ} }}{\xi _j},\;\;\;j = 1,2, \cdots ,11. $ (5)

其中i代表齿面上的第i个点。

$ \mathit{\boldsymbol{S}} = \left[ {\mathit{\boldsymbol{s}}\left( {{\xi _1}} \right),\mathit{\boldsymbol{s}}\left( {{\xi _2}} \right), \cdots ,\mathit{\boldsymbol{s}}\left( {{\xi _{11}}} \right)} \right]. $ (6)

设加工参数组成的调整向量为{δξj}。为了求解叠加修形面Γ后加工参数的变化值,建立加工参数反求问题的数学模型如下:

$ \begin{array}{l} {\rm{find}}\left\{ {{\rm{ \mathsf{ δ} }}{\xi _j}} \right\},j = 1,2, \cdots ,11,\\ \min {\left\| {\mathit{\boldsymbol{S}}\left\{ {{\rm{ \mathsf{ δ} }}{\xi _j}} \right\} + \mathit{\boldsymbol{ \boldsymbol{\varGamma} }}} \right\|_2}. \end{array} $ (7)

采用最小二乘法求解式(7),可得到修形面控制参数对应的齿面加工参数。在修形面设计基础上,建立齿面加载性能多目标优化问题数学模型。

1.1 优化设计变量

以齿面设计中用到的3个控制参数作为优化设计变量:a用于控制空载传动误差曲线,ks用于控制接触迹线,b用于控制接触椭圆长轴。变量的设计空间需人为给定。对于a,其值不能使空载传动误差曲线过于平缓,否则易趋于完全共轭齿面,也不能过陡,否则会导致空载噪声变大。通常推荐控制UTE值在-40~-60 μrad,利用式(2)可推算a的范围。接触迹线斜率按照最大和最小能取到的斜率范围取值。接触椭圆半宽b一般取为0.3~0.4。

1.2 目标函数与约束函数

以加载传动误差峰峰值LTE和齿面最大接触应力σc为目标函数。加载传动误差为螺旋锥齿轮振动噪声的激励源,一般来说,该值越小,振动激励源越小。对于σc,为保证齿面不出现失效,在不出现边缘接触的情况下σc越小越好。

为了准确反映加载接触区位置,将加载接触分析结果放置在大轮投影面内,如图 4所示。黑线为LTCA得到的接触区;蓝线为加载后理想的接触区,理想接触区位置为大小端各向内缩进齿宽的5%,齿顶和齿根各向内缩进全齿高的5%,整个区域占齿宽和齿高的90%。为保证不出现边缘接触,定义

$ \Delta S = {S_{\rm{L}}} - {S_{\rm{c}}}. $ (8)
图 4 (网络版彩图)加载接触区域示意图

式中:SL为LTCA得到的接触区面积;Sc为理想接触区与LTCA得到的接触区重合部分的面积;ΔS为两者差值,即LTCA接触区落在理想接触区以外的面积,即图 4中阴影区域,通过控制该区域小于设定容差εS,可以保证加载接触区域位于理想区域内,且不出现边缘接触。

控制加载接触区的同时,还应控制小轮和大轮最大齿根弯曲应力σpσg不超过许用值σb

最终,根据设计变量、目标及约束函数,齿面加载性能多目标优化问题的数学模型如式(9)所示,其中ksab的下标1、2分别为取值范围上、下限。

$ \begin{array}{*{20}{c}} {\min {\sigma _{\rm{c}}},{\rm{LTE,}}}\\ {{\rm{Vars}}:{k_{\rm{s}}},a,b,}\\ {{\rm{s}}.\;{\rm{t}}.:\Delta S - {\varepsilon _S} \le 0,}\\ {{\sigma _{\rm{p}}} - {\sigma _{\rm{b}}} \le 0,}\\ {{\sigma _{\rm{g}}} - {\sigma _{\rm{b}}} \le 0,}\\ {{k_{\rm{s}}} \in \left[ {{k_{{\rm{s1}}}},{k_{{\rm{s2}}}}} \right],a \in \left[ {{a_{\rm{1}}},{a_{\rm{2}}}} \right],b \in \left[ {{b_{\rm{1}}},{b_{\rm{2}}}} \right].} \end{array} $ (9)
2 齿面加载接触分析(LTCA)方法

在对齿面多目标优化问题式(9)进行求解时,需要进行多次LTCA,本文采用了半解析LTCA计算方法[16]以提高计算效率。为考虑齿根弯曲应力,LTCA计算过程分为两部分,如图 5所示:第1部分是齿面载荷分布计算;第2部分是生成精细有限元网格,利用齿面载荷分布进行齿根弯曲应力计算。

图 5 (网络版彩图)半解析LTCA计算流程图

LTCA齿面载荷分布计算如图 5中黑色虚线框所示。LTCA计算过程中,一个完整的齿啮入到啮出包含多个时刻,每个时刻通常有2对或3对齿发生接触。每个时刻需要同时满足位移协调方程式(10)和转矩协调方程式(11)。位移协调是指每一啮合时刻、每个齿面接触线上的点自身受力后产生的变形协调。LTCA方法中将齿面受载后产生的变形分为接触变形、剪切变形和弯曲变形[17]。轮齿弯曲变形通过建立轮齿粗糙的有限元模型获得,若接触点载荷向量为F,则式(10)中(Rp+ Rg)F即为弯曲变形量;接触变形和剪切变形Sc(F)采用Weber经验公式[17]计算。通过位移协调方程可以保证每个齿各接触点变形协调,而各齿间载荷的分配则通过转矩协调方程实现。转矩协调方程是指在某一啮合时刻,所有发生接触的齿面上,各接触点上的力产生的转矩之和与输入转矩平衡,如式(11)所示。

$ \left( {{\mathit{\boldsymbol{R}}_{\rm{p}}} + {\mathit{\boldsymbol{R}}_{\rm{g}}}} \right)\mathit{\boldsymbol{F}} + {\mathit{\boldsymbol{S}}_{\rm{c}}}\left( \mathit{\boldsymbol{F}} \right) - \mathit{\boldsymbol{Z}} + {\mathit{\boldsymbol{d}}_0} < {\mathit{\boldsymbol{\varepsilon }}_1}, $ (10)
$ \sum {\left( {\mathit{\boldsymbol{F}} \times \mathit{\boldsymbol{r}}} \right) \cdot \mathit{\boldsymbol{p}}} - {T_{\rm{s}}} < {\varepsilon _2}. $ (11)

式中:RpRg分别为小轮和大轮齿面接触点法向柔度矩阵,对于柔度矩阵等的计算可参考文[16],在此不再赘述;F为接触载荷向量;Sc(F)为接触变形和剪切变形;Z为齿轮刚体转动产生的位移,可以等效看成附加的一个传动误差;d0为两齿面接触点的初始间距;r为接触点位置矢量;p为齿轮转动轴线矢量;Ts为总转矩;ε1ε2为收敛容差。

在获得齿面载荷分布后,第2部分计算齿根弯曲应力,如图 5中红色虚线框所示。齿根弯曲应力的计算需要生成精细的轮齿有限元网格,并将LTCA计算第1部分得到的载荷分布施加于精细的有限元网格中,为提高计算效率, 直接将加载后模型导入商用有限元软件Nastran并返回结果。

综上所述,图 5中LTCA计算的主要流程归纳如下:

步骤1  利用修形面获得大小轮齿面点,计算可能的接触点,生成粗糙的有限元模型用于弯曲变形计算。

步骤2  给定啮合初始时刻t0=1,计算各齿面之间的法向间距d0

步骤3  附加转动调整量TE(TE0为迭代初值)。

步骤4  找到发生接触的接触点,对其施加载荷F(F0为迭代初值),若接触线上点的变形不满足位移协调方程(10),则对载荷值进行迭代(ΔF)直至满足;若满足位移协调方程,则转步骤5。

步骤5  计算各齿面接触点上的载荷产生的转矩是否满足转矩协调方程(11),若不满足则对施加的调整量TE进行迭代(ΔTE),得到新的调整量,并转到步骤3;若满足转矩协调方程,则获得该时刻各齿面载荷分布。

步骤6  利用程序生成精细有限元模型,并将步骤5中得到的载荷分布施加到精细模型上。调用Nastran求解该模型,返回该时刻最大齿根弯曲应力值。

步骤7  更新计算时刻至下一啮合位置,重复步骤2—6,直至啮合结束时刻。

3 多目标优化问题求解

从式(9)可以看出,此优化问题在求解过程中需要进行多次LTCA计算,但考虑齿根弯曲应力后的LTCA计算量大大增加,需要耗费大量计算时间。为提高计算效率,本文采用代理模型与经典多目标优化算法相结合的方式求解式(9)。求解过程主要步骤包括试验设计与代理模型的建立、基于代理模型的多目标优化问题求解。

3.1 试验设计与代理模型建立

本文选用空间填充效果较好的拉丁超立方(Latin hypercube sampling, LHS)设计,通过LHS获得试验设计样本点。

建立试验设计样本点后,选择Kriging方法建立代理模型,Kriging方法适用于高度非线性、多峰值问题[23]。Kriging模型利用部分已知信息模拟未知点的信息,它由局部模型和全局偏差组成,响应值与自变量之间关系为

$ y\left( x \right) = f\left( x \right)\beta + z\left( x \right). $ (12)

其中:y(x)为未知Kriging模型;f(x)为回归模型,通常为多项式函数;β为待定的回归系数,z(x)为均值为0、方差为σ2的随机过程。在模型中,f(x)和β为确定性部分,提供了设计空间内的全局近似模拟;z(x)为模型的局部偏差。

Kriging代理模型的建立过程详见文[24],其中主要包括多项式回归模型和相关函数的选择,本文选择二次多项式模型作为回归模型,并选择Gauss相关函数[24]

$ R\left( {{x^i},{x^j}} \right) = \exp \left( { - \sum\limits_{k = 1}^n {{\theta _k}{{\left| {x_k^i - x_k^j} \right|}^2}} } \right). $ (13)

其中:θk为未知参数,表示待估的第k个元素;ij代表试验设计第ij个样本点。

3.2 基于代理模型的多目标优化问题求解

对于多目标优化问题,为了获得其Pareto解集,常用的优化方法为多目标遗传算法NSGA-Ⅱ[25]。采用代理模型进行多目标优化时,由于对复杂问题直接进行试验设计获得的代理模型精度难以满足要求,因此需要对样本点进行补充以提高代理模型在最优解集附近的精度。基于设计点提高期望(expected improvement, EI)指标[24]动态补充样本点方法是Kriging代理模型中最常使用的方法。对目标函数y=f(x),EI指标表达式为式(14)。EI值越大,则表示在设计点x目标函数有可能有更好解,或者在设计点x代理模型精度较差。

$ \begin{array}{*{20}{c}} {{\rm{EI}}\left( x \right) = {I_x}\mathit{\Phi }\left( {\frac{{{I_x}}}{s}} \right) + s\phi \left( {\frac{{{I_x}}}{s}} \right),}\\ {{I_x} = \left\{ \begin{array}{l} {f_{\min }} - \bar y,\;\;\;\;{f_{\min }} - \bar y > 0;\\ 0,\;\;\;\;\;\;\;\;\;\;\;\;\;{f_{\min }} - \bar y \le 0. \end{array} \right.} \end{array} $ (14)

其中:fmin为样本点中该目标函数的最小值,ys分别为x处代理模型预测值和预测均方根误差,Φ是均值为y、标准差为s的累计分布函数,ϕ是均值为y、标准差为s的概率密度函数。Jeong等[26]将EI推广到多目标优化问题,用每个目标函数的EI值作为多目标优化时目标函数的适应度值,寻找各目标函数中具有最大EI值的设计点,并将其补充至样本点。基于EI值,本文采用代理模型和NSGA-Ⅱ求解多目标优化问题的计算流程如图 6所示。

图 6 求解多目标优化问题流程图

步骤1  采用试验设计样本点拟合优化问题目标函数和约束函数的Kriging代理模型。

步骤2  建立式(15)用于最大EI值的求解。利用步骤1得到的代理模型和NSGA-Ⅱ方法,求解式(15)获得各目标函数的最大EI值及其对应的设计点。式(15)中目标函数和约束函数均采用步骤1得到的代理模型计算。

$ \begin{array}{*{20}{c}} {\min \;{{\text{EI}}_{{\sigma _{\rm{c}}}}},{{\text{EI}}_{{\rm{LTE}}}}{\rm{,}}}\\ {{\rm{Vars}}:{k_{\rm{s}}},a,b,}\\ {{\rm{s}}.\;{\rm{t}}.:\Delta S - {\varepsilon _S} \le 0,}\\ {{\sigma _{\rm{p}}} - {\sigma _{\rm{b}}} \le 0,}\\ {{\sigma _{\rm{g}}} - {\sigma _{\rm{b}}} \le 0,}\\ {{k_{\rm{s}}} \in \left[ {{k_{{\rm{s1}}}},{k_{{\rm{s2}}}}} \right],a \in \left[ {{a_{\rm{1}}},{a_{\rm{2}}}} \right],b \in \left[ {{b_{\rm{1}}},{b_{\rm{2}}}} \right].} \end{array} $ (15)

步骤3  判断各目标函数最大EI值(EIσc和EILTE)是否小于各自的设定阈值。若大于设定阈值,则利用LTCA方法计算各最大EI值设计点的目标函数和约束函数,并将该设计点加入样本点集后转步骤1;若小于设定阈值,则样本点补充完成,转步骤4。

步骤4  利用补充后的样本点拟合Kriging代理模型,采用该代理模型和NSGA-Ⅱ方法对优化问题式(9)进行求解,并获得最终Pareto解集。

4 优化算例

将本文方法用于螺旋锥齿轮齿面优化设计,齿轮基本设计参数见表 1。以正车面为例,运行工况为满载输入扭矩1 500.0 N·m,齿轮错位量见表 2图 7a为原始设计LTCA加载接触区,在齿顶和小端脱出齿面边界;图 7b为原始设计接触应力云图,虚线为理想接触区,接触区存在边缘接触,最大接触应力为1 583.6 MPa,出现在齿顶;加载传动误差LTE为28.67 μrad;图 7c7d分别为原始设计大小轮齿根弯曲应力云图,大小轮最大齿根弯曲应力分别为291.9和321.5 MPa。建立齿面加载性能多目标优化问题数学模型如式(16)所示。最大弯曲应力许用值300.0 MPa,容差εS取0.1 mm2

$ \begin{array}{*{20}{c}} {\min {\sigma _{\rm{c}}},{\rm{LTE,}}}\\ {{\rm{Vars}}:{k_{\rm{s}}},a,b,}\\ {{\rm{s}}.\;{\rm{t}}.:\Delta S - 0.1 \le 0,}\\ {{\sigma _{\rm{p}}} - 300 \le 0,}\\ {{\sigma _{\rm{g}}} - 300 \le 0,}\\ {{k_{\rm{s}}} \in \left[ { - 7, - 2.5} \right],a \in \left[ { - 0.0035, - 0.002} \right],}\\ {b \in \left[ {0.33,0.4} \right].} \end{array} $ (16)
表 1 螺旋锥齿轮设计参数表
齿数 23(小轮),35(大轮)
旋向 左(小轮),右(大轮)
中点螺旋角/(°) 35
齿面宽/mm 41(小轮),41(大轮)
模数 6.772
轴交角/(°) 90
平均压力角/(°) 22.5
大轮刀盘名义半径/mm 127
轮齿收缩方式 标准收缩
全齿高系数 1.888
工作齿高系数 1.7

表 2 运行工况齿轮错位量
ΔP/μm ΔW/μm ΔE/μm ΔΣ/mrad
小轮 41.48 -69.38 -58.62 0.029 1
大轮 61.98 24.89 -53.88 -0.023 4
综合 103.46 -44.49 -112.50 0.005 7

图 7 (网络版彩图)原始设计加载性能

对优化问题的设计变量进行试验设计,样本点数取为33[24],采用LTCA方法计算样本点真实响应。建立优化问题式(16)的代理模型,采用图 6流程通过EI指标对样本点进行补充以提高模型精度。设置σc和LTE的EI阈值分别为5和0.1(阈值通常可设定为样本最小值的1%或更小[24]),经过35次补充后,各目标函数最大EI值小于设定阈值,完成样本点补充。EI值变化过程曲线如图 8所示。采用交叉检验方法对代理模型的精度进行检验,得到目标函数σc和LTE的拟合优度R2值分别为0.934和0.956,精度满足要求。

图 8 目标函数EI值变化曲线

获得代理模型后,采用NSGA-Ⅱ多目标优化方法进行优化计算(种群数目100,代数200,交叉概率为0.9,变异概率为0.05),最后计算得到优化结果Pareto解集如图 9所示。取Pareto中间位置的解G1进行LTCA分析和实际生产加工,G1解对应的优化设计变量为a=-0.002 48、b=0.373 8、ks=-4.466 4,对应的小轮凹面加工参数见表 3。LTCA分析计算结果与代理模型计算结果对比见表 4,代理模型与LTCA计算最大误差为2.38%(ΔS绝对误差很小,不再计算相对误差),进一步说明代理模型精度满足要求。最优解G1的加载性能计算结果与原始设计对比见表 5。齿面加载接触区如图 10a所示,全部位于理想接触区内;加载接触应力如图 10b所示,最大值为1 398.7 MPa;大轮与小轮齿根弯曲应力最大时刻分别如图 10c10d所示,对应的最大齿根弯曲应力分别为286.0和281.0 MPa;加载传动误差曲线如图 11所示,LTE值为20.67 μrad。

图 9 优化问题的Pareto解集

表 3 G1解优化设计变量对应的小轮凹面加工参数
基本摇台角/(°) 63.848 67
垂直轮位/mm 0.418 11
水平轮位/mm -1.512 51
安装根锥角/(°) 32.197 86
径向刀位/mm 113.243 37
滚比 1.807 72
床位/mm -0.726 30
刀转角/(°) 202.025 09
刀倾角/(°) 0.100 72
刀盘半径/mm 123.472 99
刀具齿形角/(°) 19.616 70

表 4G1目标函数与约束函数代理模型计算值与LTCA计算值对比
计算模型 σc/MPa LTE/μrad σg/MPa σp/MPa ΔS/mm2
代理模型 1 394.2 20.34 285.9 274.3 0.077
LTCA 1 398.7 20.67 286.0 281.0 0.099
误差/% 0.32 1.60 0.03 2.38

表 5 优化后Pareto解集中G1解与原始设计的计算结果对比
σc/MPa LTE/μrad σg/MPa σp/MPa ΔS/mm2
原始设计 1 583.6 28.67 291.9 321.5 10.15
优化设计 1 398.7 20.67 286.0 281.0 0.099
降低百分比/% 11.7 27.9 2.0 12.6 99.0

图 10 (网络版彩图)解G1优化设计的加载性能

图 11 优化设计得到的加载传动误差曲线与未优化时的对比

为了对接触区进行验证,按解G1进行磨齿加工,并将加工后的齿轮副安装到驱动桥主减速器中,对其进行加载接触台架试验。加工及台架试验整体流程见图 12。首先,按照表 3得到优化设计变量对应的加工参数在Gleason 600G磨齿机进行磨齿加工,并利用齿面误差检测的手段保证齿面精度满足小于10 μm,加工得到齿轮副后将大小轮在滚动检测机上进行空载运动滚检,图 12所示的空载接触区与设计区域吻合良好;其次,将其安装到驱动桥主减速器中;然后,在驱动桥试验台架上进行加载接触印迹试验。驱动桥加载接触印迹试验台架主要组成部分如图 13所示,包括了驱动桥支撑座、轮毂与连接法兰、主轴和左右两轮边的加载电机等。在进行加载接触印迹试验时,首先给定一较小载荷使齿轮的小轮凹面和大轮凸面开始工作;其次,在桥壳开孔处,在大轮齿面上涂抹红丹粉用于后续接触区位置的观测(啮合后不含红丹粉的区域是齿面加载接触区);然后,在台架的操作台逐渐增加扭矩至满载后,调低输入轴转速,通过桥壳开孔拍照记录接触区。进行两次该过程试验,两次得到的接触区一致。将试验后的加载接触区与LTCA分析结果进行对比,试验得到的接触区和LTCA得到的接触区如图 14所示,为便于对比,将LTCA得到的接触区边界点绘制为三维图形。从图 14中可以看出,按优化解加工的齿轮副在满载时试验得到的接触区与优化设计得到的计算结果一致。

图 12 (网络版彩图)螺旋锥齿轮驱动桥加载接触印迹台架试验流程

图 13 (网络版彩图)驱动桥加载接触印迹试验台架

图 14 (网络版彩图)优化后LTCA分析与试验得到的接触区对比

由计算和试验结果可以得出如下结论:

1) 试验得到的齿面空载接触区与理论计算结果一致(图 12),加载接触区与LTCA分析结果一致(图 14),验证了半解析法LTCA计算得到的齿面接触区的正确性。

2) 优化前LTCA计算得到的齿面最大接触应力为1 583.6 MPa,出现在大轮齿顶位置,由于在齿顶存在边缘接触导致产生了齿面应力集中。优化后接触区位于理想接触区内,最大接触应力出现在齿面中部位置,通过图 10b可以看出, 载荷分布比原设计更加合理,最大接触应力比原设计减小11.7%。在齿面载荷分布改善的基础上,加载传动误差峰峰值比原设计减小27.9%。由于加载传动误差是主要噪声的激励源,也是在后续进行动力学分析中的激励的输入值,峰峰值的下降有利于降低齿轮的噪声指标。

3) 优化前齿面接触区在齿顶和小端存在脱出齿面的情况,且存在边缘接触产生的应力集中;优化后齿面接触区齿顶、齿根、小端和大端均未脱出且无边缘接触。理想接触区外的面积减少99%,加载接触区保持在设定的边界内。

4) 优化前齿轮最大齿根弯曲应力云图如图 7c7d所示,小轮最大齿根弯曲应力大于大轮;优化后大小轮的最大齿根弯曲应力均比原设计减小,小轮最大齿根弯曲应力降低后略低于大轮。大轮最大齿根弯曲应力降低2%,小轮降低12.6%。由于小轮齿数比大轮少,因此小轮齿根弯曲应力下降较多有利于增加小轮的寿命。

5) 采用本文齿面加载性能多目标优化模型得到了齿面优化问题的Pareto解集,其中的典型解G1各性能指标达到了优化效果。

5 结论

本文针对螺旋锥齿轮齿面加载性能优化问题,建立了齿面加载性能多目标优化问题的数学模型,该模型增加了加载分析中的齿根弯曲应力性能指标。为计算优化问题中的目标函数和约束函数,建立了考虑齿根弯曲应力计算的半解析齿面加载接触分析模型。此外,针对优化模型求解效率的问题,在齿面多目标优化问题中引入了Kriging代理模型,采用基于代理模型结合多目标遗传算法的方式对齿面多目标优化问题进行了求解。算例及试验结果表明:本文提出的齿面多目标优化方法能够较准确地获得齿面优化设计的Pareto解集,同时兼具较高的计算效率。

参考文献
[1] LITVIN F L. Gear geometry and applied theory[M]. Englewood Cliffs, USA: Prentice-Hall, 1994.
[2] SHIH Y P. A novel ease-off flank modification methodology for spiral bevel and hypoid gears[J]. Mechanism and Machine Theory, 2010, 45(8): 1108–1124. DOI:10.1016/j.mechmachtheory.2010.03.010
[3] 曹雪梅, 邓效忠, 聂少武. 基于共轭齿面修正的航空弧齿锥齿轮高阶传动误差齿面拓扑结构设计[J]. 航空动力学报, 2015, 30(1): 195–200.
CAO X M, DENG X Z, NIE S W. Ease-off flank topography design for aviation spiral bevel gears with higher-order transmission errors by modification of conjugate flank[J]. Journal of Aerospace Power, 2015, 30(1): 195–200. (in Chinese)
[4] ACHTMANN J, BAÄR G. Optimized bearing ellipses of hypoid gears[J]. Journal of Mechanical Design, 2004, 125(4): 739–745.
[5] 田行斌, 方宗德. 弧齿锥齿轮加工参数的全局优化设计[J]. 航空动力学报, 2000, 15(1): 75–77.
TIAN X B, FANG Z D. Overall optimization on machine-tool settings of spiral bevel gears[J]. Journal of Aerospace Power, 2000, 15(1): 75–77. (in Chinese)
[6] 周彦伟, 方宗德, 邓效忠. 航空弧齿锥齿轮齿面接触应力计算与优化[J]. 航空动力学报, 2002, 17(3): 373–376.
ZHOU Y W, FANG Z D, DENG X Z. Calculation and optimization for tooth contact stress of spiral bevel gears[J]. Journal of Aerospace Power, 2002, 17(3): 373–376. (in Chinese)
[7] 张金良, 李少华, 方宗德, 等. 准双曲面齿轮的齿面优化设计[J]. 机械科学与技术, 2007, 26(3): 366–368.
ZHANG J L, LI S H, FANG Z D, et al. Hypoid gear machining parameter optimization[J]. Mechanical Science and Technology, 2007, 26(3): 366–368. (in Chinese)
[8] 刘光磊, 樊红卫, 谷霁红, 等. 弧齿锥齿轮传动误差曲线优化的半变性法[J]. 航空动力学报, 2010, 25(8): 1888–1893.
LIU G L, FAN H W, GU J H, et al. Half-modified roll for the optimization of transmission error function of spiral bevel gears[J]. Journal of Aerospace Power, 2010, 25(8): 1888–1893. (in Chinese)
[9] 刘光磊, 刘则良. 弧齿锥齿轮传动误差曲线的双层优化法[J]. 农业机械学报, 2012, 43(7): 223–227, 222.
LIU G L, LIU Z L. Dual-layer optimization for transmission errors of spiral bevel gears[J]. Transactions of the Chinese Society for Agricultural Machinery, 2012, 43(7): 223–227, 222. (in Chinese)
[10] SIMON V. Optimal machine tool setting for hypoid gears improving load distribution[J]. Journal of Mechanical Design, 2001, 123(4): 577–582. DOI:10.1115/1.1414129
[11] SIMON V. Design of face-hobbed spiral bevel gears with reduced maximum tooth contact pressure and transmission errors[J]. Chinese Journal of Aeronautics, 2013, 26(3): 777–790. DOI:10.1016/j.cja.2013.05.005
[12] ARTONI A, BRACCI A, GABICCINI M, et al. Optimization of the loaded contact pattern in hypoid gears by automatic topography modification[J]. Journal of Mechanical Design, 2008, 131(1): 011008.
[13] GABICCINI M, BRACCI A, GUIGGIANI M. Robust optimization of the loaded contact pattern in hypoid gears with uncertain misalignments[J]. Journal of Mechanical Design, 2010, 132(4): 041010. DOI:10.1115/1.4001485
[14] ARTONI A, KOLIVAND M, KAHRAMAN A. An ease-off based optimization of the loaded transmission error of hypoid gears[J]. Journal of Mechanical Design, 2009, 132(1): 011010.
[15] ARTONI A, GABICCINI M, GUIGGIANI M, et al. Multi-objective ease-off optimization of hypoid gears for their efficiency, noise and durability performances[C]//Proceedings of the International Design Engineering Technical Conferences and Computers and Information in Engineering Conference. Washington DC, USA: ASME, 2011: DETC2011-47715.
[16] 王琪, 周驰, 桂良进, 等. 准双曲面齿轮LTCA弯曲变形约束方法[J]. 航空动力学报, 2017, 32(11): 2800–2807.
WANG Q, ZHOU C, GUI L J, et al. Constraint method of calculating tooth bending deformation of hypoid gears in LTCA[J]. Journal of Aerospace Power, 2017, 32(11): 2800–2807. (in Chinese)
[17] WILCOX L E, CHIMNER T D, NOWELL G C. Improved finite element model for calculating stresses in bevel and hypoid gear teeth[R]. Alexandria, USA: AGMA, 1997.
[18] HOTAIT M A, KAHRAMAN A, NISHINO T. An investigation of root stresses of hypoid gears with misalignments[J]. Journal of Mechanical Design, 2011, 133(7): 071006. DOI:10.1115/1.4004224
[19] 周驰, 彭钱磊, 丁炜琦, 等. 汽车驱动桥主减速器支承刚度的有限元分析[J]. 汽车工程, 2016, 38(8): 981–988.
ZHOU C, PENG Q L, DING W Q, et al. Finite element analysis on the support stiffness of final drive in vehicle drive axle[J]. Automotive Engineering, 2016, 38(8): 981–988. (in Chinese)
[20] GOSSELIN C, NONAKA T, SHIONO Y, et al. Identification of the machine settings of real hypoid gear tooth surfaces[J]. Journal of Mechanical Design, 1998, 120(3): 429–440. DOI:10.1115/1.2829170
[21] LIN C Y, TSAY C B, FONG Z H. Computer-aided manufacturing of spiral bevel and hypoid gears with minimum surface-deviation[J]. Mechanism and Machine Theory, 1998, 33(6): 785–803. DOI:10.1016/S0094-114X(97)00101-8
[22] ARTONI A, GABICCINI M, GUIGGIANI M. Nonlinear identification of machine settings for flank form modifications in hypoid gears[J]. Journal of Mechanical Design, 2008, 130(11): 112602. DOI:10.1115/1.2976454
[23] GIUNTA A A, WATSON L T. A comparison of approximation modeling techniques: Polynomial versus interpolating model[R]. Blacksburg, USA: AIAA, 1998.
[24] JONES D R, SCHONLAU M, WELCH W J. Efficient global optimization of expensive black-box functions[J]. Journal of Global Optimization, 1998, 13(4): 455–492. DOI:10.1023/A:1008306431147
[25] DEB K, PRATAP A, AGARWAL S, et al. A fast and elitist multiobjective genetic algorithm:NSGA-Ⅱ[J]. IEEE Transactions on Evolutionary Computation, 2002, 6(2): 182–197. DOI:10.1109/4235.996017
[26] JEONG S, OBAYASHI S. Efficient global optimization (EGO) for multi-objective problem and data mining[C]//Proceedings of 2005 IEEE Congress on Evolutionary Computation. Edinburgh, UK: IEEE, 2005, 3: 2138-2145.