斜齿面齿轮齿面接触问题的有限元-线性规划解法

宋相男，靳广虎

(南京航空航天大学 直升机传动技术实验室， 江苏 南京 210016)

0 引言

本质上而言，面齿轮传动是渐开线圆柱齿轮与锥齿轮相啮合的齿轮传动[1]，在直升机主减速器的传动系统、铣床主轴的传动系统以及自行车的无链条传动系统中均有很好的应用。其特有的结构和啮合特性使其具备诸多优势，成为21世纪直升机传动系统研究的重点。

LITVIN.F.L研究了面齿轮齿面成形的理论方法及齿面接触应力的分析[2-3]；沈云波等分析了斜齿面齿轮齿宽限制条件[4-5]；方宗德、WU S H及YE.S.Y等人基于齿轮齿面的柔度矩阵对齿面接触应力进行了计算[6-9]；T.F.Conry等基于线性规划的思想提出了弹性接触问题的解法[10-12]；唐进元基于有限元的思想计算了螺旋锥齿轮的啮合刚度[13]。

本文提出基于有限元-线性规划法对斜齿面齿轮齿面接触问题进行计算的方法，基于齿面接触区域的载荷和变形量计算了齿轮副的啮合刚度及最大压应力。

1 斜齿面齿轮的齿面成形及接触轨迹

1.1 齿面数学模型

斜齿面齿轮齿面十分复杂，无法直接将其齿面方程表示出来，需借助其加工刀具齿面方程，并结合齿轮啮合方程间接对其齿面加以描述。

1) 插齿加工坐标系

根据斜齿面齿轮插齿加工原理，建立如图1所示加工坐标系。其中，ssc(o-xsc,ysc,zsc)为刀具的固定坐标系，ss(o-xs,ys,zs)是与刀具固联的转动坐标系，sfc(o-xfc,yfc,zfc)为面齿轮的固定坐标系，sf(o-xf,yf,zf)是与面齿轮固联的转动坐标系，φs和φf分别为刀具与面齿轮的转角。

2) 刀具的齿面方程

图2为刀具端面齿廓及坐标系，该坐标系与图1中的坐标系对应。其中，θs为刀具渐开线齿廓上的角度参数；θsc为刀具基圆上齿槽宽所对应圆心角之半,其计算公式为

θsc=2π/Ns-(tanαst-αst)

(1)

式中：Ns为刀具齿轮的齿数；αst表示刀具端面压力角。

图1 斜齿面齿轮插齿加工坐标系

图2 刀具齿面坐标系

渐开线齿廓方程为：

(2)

式中：rbs为刀具基圆半径；“±”分别表示渐开线齿廓Ⅰ和Ⅱ。

渐开线斜齿轮刀具的齿面由端面齿廓沿z轴作螺旋运动扫略形成，故刀具齿面1(即齿廓Ⅰ螺旋而成的齿面)方程表示为：

(3)

同理，齿面2方程为：

(4)

式中：γs为端面齿廓沿z轴螺旋的角度；ps为螺旋参数，其计算公式为

ps=2πrbs/tanβb

(5)

式中βb为刀具基圆螺旋角。

3) 斜齿面齿轮的齿面方程

根据包络理论和齿轮啮合原理，则得斜齿面齿轮的齿面方程为：

(6)

式中：Mf,s是坐标系ss到sf的坐标转换矩阵；f(θs,γs,φs)是齿面啮合方程。Mf,s表示为：

其中：v(s, f)为刀具与面齿轮啮合点处的相对运动速度；n表示啮合点处法向量。

1.2 齿面啮合轨迹

采用齿数N1比刀具齿数Ns少1～3个，且其他参数与刀具一致的渐开线斜齿圆柱齿轮与面齿轮啮合，实现点接触传动。因齿面参数相同，故斜齿轮齿面方程与刀具相类似，其齿面方程用r1(θ1,γ1)表示。

啮合传动时，两轮齿面连续相切，则在固定坐标系中，两轮齿面接触点具有相同的坐标和法向量，即：

(7)

依据上述理论分析，利用matlab编写程序，实现斜齿面齿轮齿面及接触轨迹的可视化，如图3。

图3 斜齿面齿轮齿面及接触轨迹

2 斜齿面齿轮齿面接触问题计算理论

2.1 面齿轮弹性接触问题的基本理论

图4 斜齿面齿轮接触情况

取区域Ⅱ包含实际接触区域Ⅰ，区域Ⅱ上各点的弹性力为Pi，区域Ⅱ上点是否在实际接触区域Ⅰ上的判别条件为：

1) 在接触区域Ⅰ:

(8)

2) 不在接触区域Ⅰ:

(9)

式中：δ为两轮齿的整体接近量，δ=δ1+δf。

各点的接触弹性变形可表示为：

ωi=ai·pi

(10)

式中ai表示各点的齿面柔度系数。

设区域Ⅱ上共有N个节点，则各点弹性力之和为F，即

(11)

由式(8)-式(11)可得到区域Ⅱ上各节点的初始距离hi、节点力pi、弹性变形ωi和轮齿的整体接近量δ满足：

(12)

记区域Ⅱ上斜齿轮和面齿轮齿面上各节点的柔度系数之和为一N×N阶矩阵A

(13)

则式(12)可改写为：

(14)

式中：P=[P1,P2,…,Pn]T，H=[h1,h2,…,hn]T为N×1阶列向量。

2.2 两齿面对应节点间初始距离hi的求解

(15)

式中：t为坐标参数，其值为任意实数。

(16)

3 齿面载荷分布的线性规划算法

3.1 面齿轮弹性接触问题的数学模型

将变形协调条件式(8)和式(9)改写为

(17)

引入松弛变量Yi(Yi≥0)，则上式可写为

(18)

(19)

将式(18)代入式(14)可获得斜齿面齿轮齿面弹性接触问题的数学模型：

(20)

式中：{e}为N×1阶单位向量；[I]为N×N阶单位矩阵；P=[P1,P2,…,Pn]T、H=[h1,h2,…,hn]T为N×1阶列向量。

该模型中需求解的未知量为(P,δ,Y)，Pi≥0，δ≥0，Yi≥0，且各未知量均需满足接触条件式(19)。

3.2 改进的单纯形法计算齿面载荷

参照线性规划的一般形式，将斜齿面齿轮的弹性接触问题改写为：

在约束条件：

(21)

下,求目标函数

(22)

的最小值。

式中：Zi为人工变量，其值非负。

该问题不同于一般线性规划问题之处在于其附加有接触条件式(19)，将该问题表示如表1所示。

表1 面齿轮弹性接触问题的线性规划表示

表1第一行元素为待求的各未知量，往后每行代表一个约束条件，共有N+2个约束条件，最右边一列的值为等式约束式(21)的右边常数项。

将表1的第N+2约束行元素减去前N+1约束行元素之和，并将结果放置于原第N+2约束行，获得改进后的形式如表2所示。

表2 面齿轮弹性接触问题的改进线性规划表示

表2中，di、Zd可分别表示为：

(23)

由表2可知，Zi的系数矩阵为一单位阵，故令Zi为基变量，令(P,δ,Y)为非基变量，取值为0，得初始基本可行解。因未知量数远大于方程数，故有多组基本可行解，但仅有一组为所需解。由线性规划的求解过程可知，得到一组基本可行解后，需不停地进行转轴运算，使各变量不断进入或退出基变量，找出一组使目标函数取最小值的可行解。

任何一个非基变量xs增加一正增量Δxs，将使目标函数变化为Zd+dsΔxs，为使目标函数减小最快，应使ds为di中的最小负值。找出变量xs所在列的所有正系数ais，其对应最右边的常数项为bi，作比值bi/ais，找出所有比值的最小者。假设在r行，则以r行s列元素为转轴元素进行转轴运算。此时，原非基变量xs进入基本可行解，而相应的Zr则退出基本可行解。

用该方法求解斜齿面齿轮的齿面接触问题时，有附加接触条件式(19)，故对进入基本可行解的变量有限制条件。假设某次进入基本可行解的变量选为Ps，则需检查与之对应的Ys是否为基变量。若否，Ps可自由进入基本可行解；若是，Ps进入后Ys退出，则此次Ps仍可进入，否则找出di中除ds以外的最小负值所在列的变量进入基本可行解；同理，若某次进入基本可行解的变量选为Ys，则应对相应的Ps进行检查。

重复上述运算过程，直至得到满足接触条件的基本可行解。实际计算表明，在2×(N+1)个循环内便可以得到唯一可行解。

3.3 齿轮副综合啮合刚度的求解

斜齿面齿轮的齿厚在齿宽方向上不断变化，故目前还未见斜齿面齿轮副啮合刚度的解析求解法。基于以上求解，可获得接触区域上各节点的载荷和变形，从而求得单对齿啮合刚度：

(24)

式中δu为单对轮齿的综合弹性变形。

本文中，δu定义为整片接触区域上所有对应节点间弹性变形的平均值，即

(25)

式中m为实际接触的节点个数。

多对齿接触时，各接触齿对间可看作是并联弹簧，故将同时啮合轮齿的单对齿啮合刚度进行分段叠加，即可得啮合副的综合啮合刚度Km：

(26)

式中s为同时啮合的轮齿对数。

3.4 实例计算

斜齿轮的法面模数mn=4,齿数N1=27，压力角α=22.5°，螺旋角β=10°，刀具齿数Ns=28，面齿轮齿数Nf=131。

1) 初始距离hi

图5所示为接触点附近两齿面对应点间的初始距离图。由图可知，斜齿面齿轮内径靠近齿根和外径靠近齿顶处，两齿面上对应点间的初始距离较大，而在接触点附近，两齿面上对应点间初始距离较小，该趋势与两齿面形状吻合。

图5 接触区域对应点间初始距离

2) 接触区域载荷分布

图6所示为斜齿面齿轮齿面接触区域上的载荷分布，该接触区域上最大载荷出现在初始接触点(椭圆中心)处，接触区域的载荷在空间内呈半椭球体分布，这与Hertz接触理论是吻合的。

图6 接触变形区域载荷分布

3) 啮合刚度

根据式(22)和式(23)求得单对齿从啮入到啮出过程中，啮合刚度计算值及傅里叶拟合曲线如图7。

图7 单齿啮合刚度曲线

将单对齿啮合刚度曲线向左、右移动角度Δφ可得临近齿对的单齿啮合刚度曲线如图8。Δφ为斜齿圆柱齿轮转过一个周期的转角与重合度的比值。

图8 各啮合齿对的啮合刚度

将各单齿对的啮合刚度按式(24)进行分段叠加，得综合啮合刚度，如图9。

图9 综合啮合刚度

4) 最大接触应力计算结果与验证

采用本文方法和Hertz理论法分别计算单对齿啮合时接触区域的最大接触应力，结果如图10所示。其中，施加的载荷F=2500N。由图可知，啮合过程中，两种方法所得接触应力的最大误差为7.87%，最小误差接近0。

图10 斜齿面齿轮齿面最大接触应力

4 结语

本文主要做了如下工作：

1)通过包络理论和齿轮啮合原理推导了斜齿面齿轮的齿面方程，并基于matlab实现其齿面及接触轨迹的可视化。

2)建立适合斜齿面齿轮接触问题的线性规划模型，对单纯形算法加以改进，求得齿面载荷分布和各接触位置的弹性变形，给出了基于有限元和线性规划求解斜齿面齿轮副啮合刚度的方法。

3) 通过与传统的Hertz接触理论计算出的接触应力进行对比，得到的误差均在8%以内，证明了本文方法的有效性。

本文为斜齿面齿轮齿面接触问题提出了一种行之有效的计算方法，给出了一种计算斜齿面齿轮副啮合刚度的途径，为斜齿面齿轮传动中接触强度的分析提供了一定的理论基础。