斜齿轮有限元接触分析参数化建模方法

尚振国，高天一，王华

（1.大连海洋大学 机械与动力工程学院，辽宁 大连 116023；2.东北财经大学 实验教学中心，辽宁 大连 116025）

0 引 言

宽斜齿轮传动在大功率齿轮传动装置中的应用日益广泛，虽然相关标准［1-2］提供了斜齿轮强度计算方法，但是这些方法不能提供应力分布、载荷分配等详细啮合状态信息。有限元接触分析法是研究齿轮啮合状态的一种有效方法，建立斜齿轮有限元接触分析模型需要解决数据存储、斜齿轮网格自动生成、啮合齿对接触区网格细化等问题。国内、外已有很多学者研究了斜齿轮有限元接触分析建模方法［3-4］，这些方法采用自编程序实现数据存储、网格自动划分，采用全齿面网格细化［5-6］来满足接触分析对网格密度的要求，导致有限元模型数据结构繁杂，节点数量庞大，严重影响计算效率。本文将斜齿轮三维有限元模型数据信息导入SQL Server 数据库，利用SQL Server数据库强大的数据处理功能实现有限元数据的快速、准确管理；采用接触区局部网格细化代替全齿面网格细化，所以可以采用较小的单元长度而不显著增加单元和节点数量，能够有效地缩减有限元模型节点规模，提高计算效率。这种方法利用VB 编程语言全部实现了程序化和参数化，在输入齿轮基本参数后即可生成三维有限元接触分析模型数据文件，用于斜齿轮啮合状态的数值仿真分析计算。

1 参数化模型的信息结构

一个完整的齿轮箱参数化模型信息结构应该包括几何结构和有限元分析两方面信息，如图1 所示。

图1 方案信息图

根据上述信息结构，在SQL Server 数据库中建立数据表，定义相应的数据项，并建立表间的关联关系。

建模的基本过程为：1）参数生成一个原始的直齿轮，并进行有限元网格划分。2）通过ANSYS 提供的APDL 语言以命令流的方式提取节点编号、坐标及单元与节点对应关系等信息，以文本文件形式保存。再通过VB 编程读入这些文本文件信息，导入SQL Server 数据库。3）在SQL Server 数据库中绕齿轮回转轴线对节点坐标进行旋转变换，生成斜齿轮网格。再根据齿轮修形参数或制造、安装误差批量修改节点坐标，以便进行修形分析或误差影响分析。4）将修改后节点坐标导入ANSYS，根据原单元、节点间关系重建模型，生成最终斜齿轮有限元模型。

2 原始直齿轮有限元网格的生成

建立原始直齿轮有限元模型的基本步骤如下：

1）利用APDL 样条曲线命令BSPLIN 根据渐开线曲线方程生成端面渐开线齿廓，并根据齿根过渡曲线生成对应的齿根圆角，然后再用画圆弧命令LARC 生成齿顶圆和齿根圆，共同构成一个端面齿廓线。

由于斜齿轮齿向修形需要综合考虑轴变形的影响，所以建立端面模型时还应生成所有轴径的圆周线，以便拖拉（Extrude）生成各个阶梯轴段的有限元网格。同时，为了后续细化接触区网格的需要，将轮齿分割为接触侧和非接触侧两部分。

2）为了减小有限元网格的规模，根据斜齿轮啮合重合度确定最多同时参与啮合的轮齿个数，根据需要保留的轮齿个数，利用旋转复制命令生成完整的端面齿廓，如图2 所示。

图2 完整的端面齿廓

3）指定各轮廓线网格划分密度，应用自由四边形方式进行网格划分，并拉伸生成齿轮体网格和各阶梯轴段网格，如图3 所示。

图3 齿轮轴网格

采用上述方法划分的网格绝大多数为六面体网格，只含有少量的退化单元。可以从总体上减少单元、节点数量，提高计算精度和效率，而且单元和节点的排列具有一定的规律性，便于后续的轮齿接触区网格细化及计算结果后处理。

3 斜齿轮有限元模型的生成

3.1 标准斜齿轮有限元网格的生成

由于完全采用拉伸方式生成直齿轮有限元网格，因此所有节点沿轴向方向形成规则的“片层”，在齿轮部分，每一层都相当于一个齿轮端面，因此将所有节点进行分层，同时绕齿轮回转轴线旋转相应角度，将修改后节点坐标导入ANSYS，根据原单元、节点间关系重建模型，得到标准斜齿轮有限元网格。APDL 命令流如下，程序中文本文件ncord 和enlink 分别保存原始直齿轮网格节点坐标和单元与节点对应关系。

*CREATE，FF

*DIM，JNSMM，，2，1，1，

*VREAD，JNSMM(1，1)，JNSM，TXT，，，IJK，2，1(f12.0)

*END

/INPUT，FF

*dim，Nndt，table，JNSMM(1，1)，3，1，Knum，CorData

*tread，Nndt，ncord，txt，，!读入节点坐标

*dim，Erdt，table，JNSMM(2，1)，8，1，Knum，CorData

*tread，Erdt，enlink，txt，，!读入单元节点号

*do，i，1，JNSMM(1，1)

N，i，Nndt(i，1，1)，Nndt(i，2，1)，Nndt(i，3，1)

*enddo !创建节点

*do，i，1，JNSMM(2，1)

E，Erdt(i，1，1)，Erdt(i，2，1)，Erdt(i，3，1)，Erdt(i，4，1)，Erdt(i，5，1)，Erdt(i，6，1)，Erdt(i，7，1)，Erdt(i，8，1)

*enddo !创建单元

3.2 修形斜齿轮有限元网格的生成

随着大功率风电齿轮箱向轻型化和环保化方向发展，对斜齿轮传动的动态性能如振动和噪声等，提出了越来越高的要求。齿轮修形技术作为补偿制造、安装误差和改善啮合性能的有效手段，近年来在风电齿轮箱中的应用越来越受到人们的重视。应用本文方法可以方便地生成修形斜齿轮有限元网格，即根据修形曲面方程［7-8］，在SQL Server 数据库中批量修改节点坐标，再导入ANSYS生成有限元模型，以便进行接触区分布、传递误差影响等修形效果仿真分析。

4 接触区六面体网格细化

由于轮齿接触区的宽度很窄，通常只有1～2 mm，为了准确地反映接触区内应力分布情况，接触区内有限元网格密度应该较大，其单元宽度应约为赫兹接触宽度的1/10，因此需要对轮齿接触区进行网格细化，如图4 所示。如何保持细化后的单元仍为六面体单元，而且新产生的单元节点与周围单元节点必须一一对应相连，这是轮齿接触区网格细化的一个难点问题。本文采用分级剖分法，图5 所示为一个六面体单元朝向接触区的面按3×3 剖分示意图，首先沿原六面体单元的一条边建立一组“核”单元，将“核”单元与原单元沿对角连线，形成一个虚拟对角面，后续生成的细化单元均在虚拟对角面上对接。其它份数剖分方法与之类似。

图4 接触区网格

图5 网格细化原理示意图

5 应用实例

表1 齿轮传动参数

以某型风电齿轮箱输出级外啮合斜齿轮为例，建立轮齿有限元接触分析模型，分析修形前后轮齿接触区及接触应力分布情况。齿轮传动参数如表1 所示，该齿轮副的总重合度为4.44，因此齿轮交替处于4 齿啮合区和5 齿啮合区。取刚进入5 齿啮合区位置进行分析，理论上1 号啮合齿对为点接触，接触线长度为零，如图6 所示。

图6 啮合区示意图

图7 和图8 分别为小齿轮和大齿轮轮齿接触区和应力分布图。其中（a）所示为未修形齿轮轮齿接触区和应力分布情况，两齿轮齿顶和齿端均存在边缘接触现象［9-10］，应力较大，并且最大应力出现于小齿轮5 号轮齿上，说明接触区偏向转矩输入端。（b）所示为修形后接触区和应力分布情况，载荷主要由中间3 对齿承担，接触区呈较为理想的椭圆形分布，不但消除了边缘接触现象，而且最大应力明显下降。按赫兹公式计算的接触应力理论值为703.5 MPa，有限元法计算的接触区中间部位接触应力约为675 MPa，两者相近，说明前述有限元建模方法是正确的。

图7 主动齿轮接触面和应力

图8 被动齿轮接触面和应力

6 结 论

采用ANSYS 与SQL Server 数据库相结合技术建立斜齿轮参数化三维整体有限元模型，应用分级剖分法实现了轮齿接触区的局部网格细化，解决了斜齿轮有限元分析快速建模、缩减网格节点规模、提高计算效率等问题，为应用有限元非线性接触分析技术模拟斜齿轮的啮合过程提供了一种新方法。

［1］ ISO 6336-1993 Calculation of load capacity of spur and helical gears［S］.

［2］ GB/T3480-1997 渐开线圆柱齿轮承载能力计算方法［S］.

［3］ 刘辉，吴昌林，杨叔子.参数化啮合斜齿轮三维有限元网格的自动生成［J］.华中理工大学学报，1997，25（4）：13-15.

［4］ 顾守丰，连小珉，颜磊，等.斜齿轮轮齿三维有限元网格自动生成及细化［J］.清华大学学报：自然科学版，1996，36（8）：77-82.

［5］ 丁能根.斜齿轮三维有限元网格和接触单元的自动生成［J］.合肥工业大学学报：自然科学版，2003，26（5）：1094-1097.

［6］ 袁丹青，康仕彬，丛小青，等.直线共轭内啮合齿轮泵齿轮接触强度的有限元分析［J］.排灌机械，2009，27（2）：111-114.

［7］ Litvin F L，Haofenget P，Lagutin S A，et al.Helical and spur gear drive with double crowned pinion tooth surfaces and conjugated gear tooth surfaces：USA，6205879［P］.2001.

［8］ Litvin F L，Gonzalez-perez I，Fuentes A，et al.Topology of Modified Surfaces of Involute Helical Gears with Line Contact Developed for Improvement of Bearing Contact，Reduction of Transmission Errors，and Stress Analysis［J］.Mathematical and Computer Modelling，2005（42）：1063-1078.

［9］ Litvin F L，Fuentes A，Gonzalez-Perez I，et al. Modified involute helical gears：computerized design，simulation of meshing and stress analysis［J］.Comput.Methods Appl.Mech.Engrg.，2003，192：3619-3655

［10］ Chen Yicheng，Tsay C B.Stress analysis of a helical gear set with localized bearing contact［J］.Finite Elements in Analysis and Design，2002，38：707-723.