基于LS-D YN A的正交面齿轮动态接触分析

胥国祥,张以都

（虚拟现实技术与系统国家重点实验室，北京航空航天大学，北京100191）

面齿轮传动(Face Gear Drive)是一种圆柱齿轮与圆锥齿轮(面齿轮)相啮合的新型齿轮传动，是采用具有渐开线齿面的刀具经范成加工而成，具有许多独特的优点和几何现象，比如轴向误差不敏感，重合度高，传动振动小等[1~2]。

在齿轮传动的动态接触分析中，传统的齿轮接触分析(TCA)只考虑了轮齿齿面啮合，使接触呈现出一种不连续的接触过程[3]。在全部啮合过程中，轮齿的接触非常复杂，比如轮齿啮合过程中的边缘接触，多齿对轮齿啮合有重合度的问题。因此，原有的TCA方法是不完善的。由于显式动力分析软件LS-DYNA的发展，使模拟齿轮传动动态啮合的仿真成为了可能。

本文采用MATLAB、CATIA软件，建立了正交面齿轮的三维几何模型，并利用LSDYNA对面齿轮啮合过程中的接触状态，进行了动力学仿真分析。

1 面齿轮齿面的基本方程[1,2,4]

图1 正交面齿轮加工图

面齿轮是用直齿渐开线齿轮刀具经范成加工而成。现采用4个坐标系（见图1）。与刀具s和齿轮2一同转动的两个坐标系Ss和S2；刀具s和齿轮2初始位置的2个固定坐标Sso和S20。其中坐标原点Os和O2都在O点；Zso和Zs同轴，为刀具的转动轴线；Z20和Z2同轴，为被加工面齿轮的转动轴线。用和分别表示刀具的转角和齿轮2的转角，在以下讨论中取ZS和Z2的夹角为90°。刀具渐开线齿面方程为：

式中，rbs为刀具渐开线的基圆半径；μs为刀具齿面上一点的轴向参数；θs为刀具渐开线上一点的角度参数。μs和θs为渐开线齿面的参数，式（1）中的上下符号分别对应于刀具齿槽两侧渐开线β-β和γ-γ。θ0s由θ0s=π/2NS-invα确定，其中Ns为刀具齿数，α为刀具压力角。

图2 渐开线齿轮参数图

式中：

M2S是坐标变换矩阵，表示从刀具坐标系向面齿轮坐标系的转换 φ2=q2sφs，q2s=Ns/N2= ω2/ωs；方程 f(μs，θs，φs)由啮合条件ns.vs,2=0（vs,2为啮合点的相对速度）确定。由此得到正交面齿轮的齿面方程：

经推导得到过渡曲面方程的分量形式：

2 面齿轮几何建模

面齿轮的齿面形状比较复杂，直接由面齿轮的方程式（1）和式（2）构造完整的面齿轮有较大困难。为此，通常采用离散点曲面拟合的方法形成面齿轮的齿面。

对式（1）来说，令公式中的 y2为某一常数Yi(i=1，2，3…，n)，即

以面齿轮参数——直齿轮模数5；直齿轮总齿数18；压力角20 deg；面齿轮总齿数100；齿宽70 mm的面齿轮传动为例，将该对齿轮的数字点导入CATIA中，最终可得三维图如图3（CAITA V5 R16）[7]。

图3 面齿轮三维图

面齿轮五齿装配模型如图4。

图4 面齿轮五齿装配图

3 面齿轮有限元分析

3.1 选择单元

将装配图模型导入ANSYS中，并选择合适的单元，在此选用solid164和shell163两种单元。

在ANSYS/LS-DYNA中，solid164单元不具有旋转自由度，不能施加转速和转矩使其转动进行接触分析。定义齿轮内圈表面为shell163单元，并将其定义为刚性体，就可以进行施加转动负载，以进行动力学接触仿真分析。

3.2 设置实常数、材料属性

Shell163单元，壳厚为1，选择S/R Huges-Liu算法，以消除某种沙漏模态。

分析采用单位：kg，mm，s，rad

材料密度：7.83×10-6kg/mm3；弹性模量：2.07×103kPa；泊松比：0.3；

对shell163刚性体单元，设置材料属性时需要设置其平移和旋转约束。

3.3 定义局部坐标系

分别在两齿轮轴线处定义两个局部坐标系，因对刚体加载默认为刚体质心方向，而并非其轴线方向，故需通过局部坐标系来进行其质心位置的修改及后面的加载。使用EDLCS，输入三个点的坐标，进行局部坐标系的定义。

3.4 划分网格

对模型进行网格划分。利用ANSYS中的Mesh Tool命令，设置相应的单元类型、材料常数，对齿轮及其内圈面进行网格划分。划分后的有限元模型如图5。

图5 面齿轮五齿有限元模型

3.5 定义PAR T

PART的定义是具有唯一的TYPE，REAL和MAT号组合的一组单元。通过EDPART命令创建4个PART（如表1所示）。

表1 创建的4个PAR T

3.6 定义刚体质心

由于所建模型为五齿模型，内圈定义为刚体，其默认质心并非在其轴线上，故需重新定义其质心位置，在此之前需定义质心坐标数组，及刚体在自定义质心下的转动惯量数组，再采用EDIPART命令进行刚体质心设置。

EDIPART,2,ADD,CVECT1,TM,1,IVECT1,,11%直齿轮内圈刚体质心定义%

EDIPART,4,ADD,CVECT2,TM,1,IVECT2,,12

%面齿轮内圈刚体质心定义%

3.7 定义接触

在ANSYS/LS-DYNA程序中，没有接触单元，只要定义可能接触的接触表面，它们之间的接触类型，以及与接触有关的一些参数，在程序计算过程中就能保证接触界面之间不发生穿透，并在接触界面相对运动时考虑摩擦力的作用[8,9]。

在本算例中，根据齿轮模型，选用表面-表面接触(STS)，定义静态摩擦系数(Static Friction Coefficient)为0.15、动态摩擦系数(Dynamic Friction Coefficient)为0.1，其余采用默认设置。在接触面和目标面分别选取主动轮part1和从动轮part3。

3.8 施加负载

荷载施加方式：在主动轮上施加匀速的角速度，在从动轮上施加恒定的扭矩。

定义荷载时间数组：

*DIM,TIME,ARRAY,4,1,1,,,

*DIM,OMIGA,ARRAY,4,1,1,,,

*DIM,L IJU,ARRAY,4,1,1,,,

其中TIME、OMIGA和L IJU分别为时间、角速度和阻力矩的数组。本算例中，角速度为100 rad/s，总转矩为-2.0×109(mN·mm)：

EDLOAD,ADD,RBOX,11,2,TIME,OMIGA,,,,,

（对小齿轮内圈刚体施加角速度载荷）

EDLOAD,ADD,MX,0,CM1,TIME,LIJU,,,,, （对面齿轮体内圈的节点组施加转拒载荷）

3.9 分析求解

（1）设置求解过程的控制参数。根据小齿轮转速100 rad/s，确定计算时间为0.01 s。结果输出文件的输出步数1 000步。时间历程文件的输出步数100步。

（2）修改K文件、求解。修改关键字：

*MAT_RIGID

1 0.783E-05 0.207E+09 0.300000 0.0 0.0 0.0

-1.00 11 111011

＄小齿轮刚体在局部坐标系下的约束＄

*DEFINE_VECTOR

11 0.0000-300.0000 -5.0000 0.0000-240.0000-5.0000

＄小齿轮内圈转动向量设置＄

*DEFINE_VECTOR

12 0.0000 0.0000-46.9917 0.0000 0.0000 -70.0000

＄面齿轮内圈转动向量设置＄

4 有限元分析结果

使用LS-PREPOST后处理器来处理最终计算结果。

图6 面齿轮在不同时刻的应力云图

面齿轮最大接触应力随小齿轮模数的变化情况（面齿轮齿数为100，小齿轮转速为100 rad/s，压力角20 deg），如图7。

图7 面齿轮最大接触应力随小齿轮模数变化曲线

面齿轮最大接触应力随面齿轮齿数的变化情况（小齿轮模数为5，小齿轮转速为100 rad/s,压力角20 deg），如图8所示。

图8 面齿轮最大接触应力随面齿轮齿数变化曲线

面齿轮最大接触应力随小齿轮转速的变化情况（面齿轮齿数为100，小齿轮模数为5,压力角20deg），如图9所示。

图9 面齿轮最大接触应力随小齿轮齿速变化曲线

由上面各图可以看出：

（1）由图6可以看出，正交面齿轮从啮合开始到退出啮合过程中，轮齿的最大应力位于轮齿中间偏齿顶位置，主要是由于齿顶处曲率半径小。

（2）由图6亦可看出，正交面齿轮在啮合过程中，大多数时间为两对齿参加啮合，其重合度接近2，故与普通圆柱齿轮传动相比，其承载能力相对提高，承载也更加平稳。

（3）本算例中，采用直齿轮齿数比刀具齿数少2，直齿轮和面齿轮传动为点接触传动。在啮合过程中，由图6可看出面齿轮传动已从点接触逐渐转化为面接触。

（4）从图7可以看出，面齿轮的最大接触应力随模数的增加而增加，近似为正比于的关系。

（5）从图8可以看出，随着面齿轮齿数的增加，面齿轮的最大接触应力随之减小，近似反比于面齿轮齿数，这主要是由于随着齿数的增加，重合度有所增加，从而使得最大接触应力有所减小。

（6）从图9可以看出，随着小齿轮转速的升高，面齿轮的最大接触应力有所下降，近似为反比于的关系。

5 结束语

本文详细介绍了用ANSYS/LS-DYNA做齿轮接触仿真分析的步骤，提出了分析过程中应注意的事项，为齿轮动力学分析提供了新的方法。

利用MATLAB和CATIA建立了面齿轮的三维实体模型；利用ANSYS/LS-DYNA计算了面齿轮五齿对在整个啮合过程中齿面的接触压力分布情况；得到了面齿轮最大接触应力随小齿轮模数、面齿轮齿数、小齿轮转速的变化趋势，对面齿轮的分析具有重要意义。

[1]朱如鹏，潘升才.面齿轮传动的研究现状与发展[J].南京航空航天大学学报，1997，29（3）：357-362.

[2]朱如鹏，潘升才.正交面齿轮传动中齿宽设计的研究[J].南京航空航天大学学报，1997，29（3）：357-362.

[3]李润方，王建军.齿轮系统动力学[M].北京：科学出版社，1996.

[4]杨连顺.正交面齿轮传动的啮合特性与强度的计算机辅助分析[D].南京：南京航空航天大学，2001.

[5]薛定宇，陈阳泉.高等应用数学问题的Matlab求解[M].北京：清华大学出版社，2005.

[6]精英科技，郝红伟.Matlab6实例教程[M].北京：中国电力出版社，2001.

[7]李苏红，潘志刚，孟祥宝，等.CATIA V5实体造型与工程图设计[M].北京：科学出版社，2008.

[8]尚晓红，苏建宇，等.ANSYS LS-DYNA动力分析方法与工程实例[M].北京：中国水利水电出版社，2006.

[9]时党勇，李裕春，张胜民.基于ANSYS/LS-DYNA8.1进行显式动力分析[M].北京：清华大学出版社，2005.