航空齿轮的疲劳寿命分析*

刘 怡，陈玲萍

(湖南工程学院，湖南 湘潭 411100)

0 引言

一般而言，S-N曲线在对数形式下通常表现为三级线性关系，分别是低周疲劳线性区域、高周疲劳线性区域和无限疲劳寿命线性区域。前两者在S-N曲线上表现为斜线段，后者在S-N曲线上表现为水平线段。

对于某型直升机主减速器分扭构型圆柱齿轮而言，其接触疲劳寿命特性也可以用S-N曲线来描述。通常而言，其表达形式有如下3种：

1)Stromeyer方程

(1)

式中：S为交变载荷；S∞为材料的疲劳极限；Sa为应力幅；N为疲劳寿命；A和α为疲劳曲线的形状参数。

2)幂函数形式

(2)

式中：m和C为疲劳曲线的形状参数。

3)三参数对数形式

lgN=C-mlg(S-S∞)

(3)

在公式(1)-公式(3)中，{A，α，S∞}，{m，C}，{m，C，S∞}均为疲劳曲线的形状参数，由材料原件或结构件的疲劳试验获得。

在以上各类型表达方式之中，本文选用幂函数的对数形式的表达式，即对公式(2)进行对数变换，可以得到：

lgN=lgC-mlgSa

(4)

式(4)可以看成两元一次的线性关系表达式，因此选用该表达式来描述某型直升机主减速器分扭构型圆柱齿轮的疲劳特性，具有两方面的优点：一是由于作为两参数的表达形式，在结构上比三参数的表达形式更为精简；二是由于作为两元一次线性方程的表达形式，便于后期针对数据进行拟合处理。

可以看到，一旦确定了齿轮的接触疲劳S-N曲线，并且得到了某工况下的齿轮接触应力幅值，我们即可依据该确定的S-N曲线得到该齿轮在该工况下的疲劳寿命。齿轮的接触疲劳S-N曲线可以通过试验的方法得到，本文中不加以展开，于是重点工作就转变成齿轮接触应力幅值的求解。

传统的齿轮接触应力幅值通常取交变循环载荷作用下接触应力的平均值，与之相比不同的是，本文将齿轮接触应力幅值用最大接触应力值来代替。这种替换方法与传统方法相比较为保守，但对于严苛的直升机主减速器的可靠性要求，具有较高的保证。因而更适用于直升机主减速器中齿轮的疲劳寿命分析。

分扭构型圆柱齿轮传动系统广泛存在于带有双发动机的军用直升机主减速器中，其目的是为了适应具有高转速、重载荷和多工况综合作用下的直升机主减速器的工作环境，提高各部件的可靠性。基于这种特点，有必要对其最重要的二、三级传动链的齿轮(分扭级直齿轮和并车级人字齿轮)进行疲劳强度的校核和疲劳寿命的分析。对于齿轮的疲劳寿命关联性因素，齿轮的应力分析最能反映齿轮的寿命状态。有限元法是一种计算复杂力学模型的有效方法，由于ABAQUS适合解决非线性接触问题，故本文选用ABAQUS对直齿轮传动系统和人字齿传动系统进行最大接触应力的有限元分析，从而给出了一种利用软件求解最大接触应力和最大弯曲应力的方法，为后续对齿轮进行疲劳寿命预测提供依据。

1 有限元分析

1.1 原理阐述

根据动量守恒原理(静态问题)，连续介质在接触力学的表现可以得到以下控制方程：

σij+bi=ρai

(5)

其中，σij是应力张量，bi是体积力张量，ρ是密度，ai是加速度。由公式(1)可以看出，应力-应变关系和应变-位移关系对于连续介质力学的研究是非常有必要的。针对接触问题，根据其运动学条件和动力学条件，考虑位移、应力、应变等参数，可以得到以下约束方程：

g(X)=0，X∈Γc

(6)

(7)

(8)

其中，Γc是发生接触的边界条件，g(X)是点集X在这两个几何体之间的距离，t是接触力，n是外法线向量。公式(6)表示的是两接触实体接触区域之间不存在缝隙，公式(7)表示的是两接触体之间的法向接触力只能是压力，公式(8)表示的是切向力与接触力之间的表达关系。

在有限元方法中，处理接触问题主要有两种计算方法：变分等式方法和变分不等式方法。以位移张量u为自变量，则变分不等式可表示为：

B(u，v)+J(u，w)-J(u，u)≥f(w-u)

(9)

其中，B(u，v)、J(u，v)以及f(v)的表达式分别如下所示，

(10)

(11)

(12)

在以上计算公式中，σij为应力张量，εij为应变张量，bi为体积力张量，ti为接触力张量，Ω为两接触体所占体积，Γt为位移施加的边界，Γc为载荷施加的边界。并且，u和w要满足以下两个条件：

a)u=w=0，位于Γt边界上；

b)nw-g≤0，位于Γc边界上。

变分等式的方法一般被广大的商业软件包所广泛采用，其表达式为：

(13)

其中：上标t和0分别表示的是时间变量为t时刻和初始时刻。文中对应力指标的分析，选用ABAQUS有限元分析软件，计算方法选择变分等式进行数值解析。

由于有限元法的应用范围很广，可以计算更复杂的问题，如冲击和断裂等，所以选择有效的算法来计算接触分析就显得尤为重要。对于复杂的问题，预测何时何地发生接触是非常困难和繁琐的任务。有时候，两实体已经发生了接触。而有时，在当前步骤或者更远的下一步会发生接触，所以接触的状态不能提前预知。在有限元中，对于接触的分析可以通过将其描述成为求解接触区域内的位移场，从而在接触边界条件的约束下，将系统势能达到最小值。边界条件的约束算法常用的是：拉格朗日乘子法和罚函数法。

拉格朗日乘子法的系统矩阵的表达式为：

(14)

其中，*K为刚度矩阵，tFc为拉格朗日的乘子向量-未知接触力向量，U为节点位移的向量，Q为在接触面的约束矩阵，*R为在节点接触部位的载荷向量，τP为节点外载荷向量。

罚函数法的原理是将接触非线性问题转化为材料非线性问题。其表达公式为：

(*K+Kα)U=*R-τPa

(15)

其中，*K是由惩罚参数得到的，而τPa由间隙函数得到的。此时如果两实体相互接触，抗力会阻止它们的边界进行相互渗透。

在拉格朗日乘子法中，拉格朗日乘数为了满足运动条件而产生一种抗力。而罚函数法中，惩罚参数就像一个弹簧的加入，以此来抵抗两实体进行接触。拉格朗日乘子法更准确，但总体刚度矩阵的大小有明显的增大。收敛的刑罚方法的精度则取决于所选的惩罚参数。如果参数太大，由于系统矩阵是病态的，计算结果可能会偏离正确值。但是如果太小，接触对可能会出现不合理的渗透。在实际分析中，导致非线性的原因是不同的，因此求解接触问题时需根据不同的工况进行选择。本文为了尽可能准确的计算出接触应力和主应力的分布情况和大小，因此在计算过程中主要采用拉格朗日乘子法。

1.2 分扭级和并车级齿轮几何建模

对齿轮接触分析的第一步就是对其进行几何建模，一般的有限元软件能建立简单的几何模型，由于齿轮结构的几何形状较为复杂，建模时参数流程设计较为繁杂，且容易出错，因此考虑采用常用的三维软件SolidWorks对分扭级和并车级齿轮进行几何建模。其中关于SolidWorks开发Geartrax插件能用于精确齿轮的自动设计和齿轮副的设计，能对直齿轮、斜齿轮、锥齿轮、链轮蜗轮蜗杆、带轮进行设计。Geartrax直齿轮生成插件如图1所示。

对于渐开线圆柱齿轮几何建模，只需在Geartrax插件中输入齿轮对的模数、齿数、齿数比、齿宽等就能在SolidWorks中分别生成大小齿轮。并在SolidWorks中将生成的齿轮装配，图2为两级传动系统单独拆分出来的示意图，图3、图4分别为分扭直齿轮几何模型和并车人字齿几何模型。

在ABAQUS中直接导入已建立好的齿轮模型，并设置相应的材料参数。本文所用分扭级和并车级的齿轮材料均为9310钢，其基本力学性能如表1所示。

表1 齿轮材料的基本力学性能

对于分扭级传动系统(直齿轮)可以直接在ABAQUS中通过六线法[1]对单齿进行切分后直接划分出结构化的网格模型，然后通过阵列得到全齿的模型。由于大轮轮齿较多，为了简化计算过程，本文对大轮只选用部分齿(两个大轮各选用9齿)进行计算。推算出的有限元模型见图5、图6，有限元单元划分成8节六面体的网格单元。其中共包含单元总数1035264个，节点总数1182236个。

Litvin等做过很多关于齿轮接触分析的研究[2-5]，其所用的齿轮接触分析大多为三齿或者五齿模型，因此上述采用部分齿来代替全齿进行接触分析的方法是可行的。

在要分析的存在接触的轮齿之间通过拉格朗日乘子法建立接触关系。以齿轮的中心点作为模型建立的参考点，并在参考点上对应点位设计新的坐标系。结合齿轮内圈点位的布局情况，构建齿轮内圈相关的耦合约束模型。非完整齿轮的两侧面也应与参考点一起进行耦合。选取各参考点上的坐标系作为相应耦合约束中的坐标系。

1.3 边界条件的设定

为确保齿轮在初始设定值呈收敛特征，文中对于设计的两种模型，分成4个分析步骤对齿轮系统的接触程度进行分析。

(1)分扭级传动系统的边界条件

分扭级由一个小轮和两个大轮组成，在有限元分析时通过给小轮施加转角，大轮施加扭矩的方式进行分析。

分析步1：在齿面接触部位存在齿轮间隙处，作为齿轮运转的初始值，通过对大齿轮施加较小的转动作用力，而此时小齿轮加固定边界条件，这样两个齿轮才能够完全接触，降低齿轮的间隙，实现在初始状态下齿轮迭代收敛。

分析步2：保持小轮固定约束，释放大轮参考节点上的旋转位移约束，给两个大轮施加初始扭矩，但该扭矩值不宜过大。

分析步3：将小齿轮的边界条件进行固定，增加大齿轮的负载扭矩力，使其达到设定值的50%。

分析步4：保持大轮的扭矩值和边界条件不变，释放小轮的旋转自由度，将小齿轮沿着旋转轴绕自身旋转至设计的角度，从而能够模拟小齿轮与接触面的齿合作用，降低大齿轮负载扭矩过高带动齿轮运转的工况值[6]。

(2)并车级传动系统的边界条件

由于并车级是一个大轮和两个小轮组成的系统，因此在对其进行有限元分析时载荷和边界条件的施加方式与分扭级存在一定的差别，各分析步的具体施加方式如下：

分析步1：在齿面接触部位存在齿轮间隙处，作为齿轮运转的初始值，通过对小齿轮自由转度上施加较小的作用力，而此时大齿轮加固定边界条件。

分析步2：将大齿轮按固定条件进行约束，释放小齿轮在参考点的旋转位移约束，并对小齿轮增加初始扭矩值，但扭矩值的选取不能过高。

分析步3：将大齿轮持续固定约束，调整小齿轮的负载扭矩值，使得该扭矩值达到设定的最高额定值。

分析步4：设定小齿轮的扭矩值、边界条件恒定状态，调整大齿轮的旋转自由度，将齿轮沿着模型的轴线部位旋转，从而能够模拟大齿轮与接触面的齿合作用，降低大齿轮负载扭矩过高带动齿轮运转的工况值。

2 分析与总结

待ABAQUS完成上述分析后，在ABAQUS后处理中可以直接读取得到分扭级和并车级的最大接触应力和最大弯曲应力的变化情况。有限元方法对模型的接触应力分布，解析后的结果见图7、图8。

根据有限元分析结果可知，分扭级的重合度小于2，并车级的重合度小于3。因此，分扭级从第二齿进入啮合时刻开始分析(在3.11时刻第二次进入啮合)，且由于分扭级的网格划分的较密，计算量大，最终只计算到3.7时刻时即终止计算。从上述图形可知，该计算结果中包含两对完整的轮齿的接触过程，已能说明分扭级的齿根弯曲应力和接触应力的变化情况。并车级从第三齿进入啮合的时刻开始分析(在3.23时刻第三齿进入啮合)，由于并车级重合度较高，让其计算到4.0时刻后再退出有限元计算。

结合上述边界条件可知，分扭级时小轮施加转角，相当于小轮作主动轮，轮齿从齿根啮合到齿顶，所以在单齿啮合区其弯曲应力逐渐增大。同理，分扭级大轮在单齿啮合区弯曲应力逐渐减小。图6的分析结果与上述结论完全一致，由此也可说明整个分析过程是正确的。

由图7和图8可知，两种情况下小轮的齿根弯曲应力要大于大轮，在分扭级大小轮的接触应力差别不大且两个大轮的接触应力和弯曲应力基本一致。并车级大轮的接触应力大于小轮且两个小轮的接触应力和弯曲应力基本一致。齿轮在基础部位的疲劳值，在文中选取应力达到最大值的疲劳状态进行计算。

根据图7和图8可知分扭级和并车级的最大接触应力和齿根位置的最大主应力与通过AGMA和HB标准计算给出的结果比较如表2所示。

由表2可知，有限元分析结果和理论结果基本一致。但是并车级的接触应力与理论值偏差较大。主要是因为该对齿轮存在一定的边缘接触(如图9(b)所示)。

在ABAQUS中可以观察到两种情况下不同时刻小轮上的接触应力和齿根应力的分布情况，分别如图9、图10所示。

由图9(a)可以看出，分扭级齿轮系统在工作时齿宽方向的接触较为均匀。由图9(b)可知，并车级齿轮在啮合过程中靠近轮齿端面和齿顶的地方均存在一定的边缘接触。由于边缘接触易引起应力集中(图9(b))，从而大大影响轮齿的接触疲劳寿命，因此这种现象应尽量避免。在工程中常通过修形的方法来避免边缘接触。

通过上述分析可知，本文选用ABAQUS对直齿轮传动系统和人字齿传动系统进行最大接触应力的有限元分析，计算得到了两种传动系统在啮合过程中的最大接触应力和最大弯曲应力的变化情况。后续可以根据得到的接触应力和弯曲应力与S-N曲线进行计算从而对齿轮进行疲劳寿命预测。