航空榫连结构接触有限元分析及影响因素探究

任浩杰，何法江，陈志雄

（上海工程技术大学 航空运输学院，上海 201620）

0 引 言

接触问题是指两物体之间由较小的接触面传递作用力，普遍存在于工程应用中，如：过盈装配、传动齿轮和连杆齿形配合等［1］。航空发动机榫连接结构就是典型的接触工作条件下的结构，被广泛应用在发动机风扇、压气机等部件的叶盘之间。研究表明，榫头和榫槽在接触面上的相互作用复杂，局部存在较大的接触应力，容易导致损伤发生破坏，需要计算接触应力，确定其危险区域［2－3］。但是接触问题是一种复杂的边界条件非线性问题，存在各种因素影响，工程应用中很难通过解析获得接触应力。目前主要采用有限元法求解复杂接触问题，有限元法将一个复杂的系统分解为大量的小单元，通过对这些相互作用的单元进行求解，去逼近无限未知量的真实系统。当前国内对榫连接结构的数值计算集中在沿接触面的应力平均值，对高应力梯度区域分析较少且精度偏低，不能精准确定危险失效区域［4］。然而，直接提高精度会导致有限元计算消耗大量的时间与资源，因此需要确定合适建模方式，使之能够高效、高质量求解，并且要严格保证计算收敛，才能得到准确的接触应力分布，进而确定危险失效区域。

本文首先从典型的赫兹接触问题着手，对比有限元解和解析解，以此确定接触问题的有限元法分析流程；其次，对压气机燕尾榫连接结构的接触问题进行分析，探究网格密度对接触应力求解精度的影响，确定危险失效区域，及进一步分析接触应力的影响因素。

1 赫兹接触问题

由于接触问题的非线性，通常采用有限元法进行求解，但在处理和计算中存在计算量过大、收敛困难和精度控制难等问题［5］。赫兹接触理论是接触应力分析的基础，可以通过有限元法对典型的赫兹接触问题进行求解，对比理论解析解，判断有限元分析结果是否合理，并以此确定采用有限元法求解接触问题的流程［6］。以两半径为50 mm、宽100 mm、泊松比0.3、弹性模量200 GPa 的平行圆柱体的接触为例，假设载荷20 KN，分别基于赫兹公式和ANSYS软件计算接触半宽及最大接触应力。算例的接触形式及应力分布，如图1 所示。根据赫兹接触中两个圆柱体的接触应力计算公式可以推导接触半宽a公式（1）和最大接触应力σHmax公式（2）。

图1 接触形式及应力分布Fig.1 Contact pattern and stress distribution

由式（1）和式（2）计算得到a＝0.240 7 mm，σHmax＝528.992 3 MPa。考虑到模型的对称性，对模型进行简化，取模型的对称部分进行有限元仿真分析如图2 所示，得到有限元解a＝0.247 3 mm，σHmax＝528.37 MPa，与根据赫兹接触理论所得的解析解吻合。

图2 两圆柱体接触的有限元模型Fig.2 Finite element model of contact between two cylinders

算例中接触区域网格尺寸设置为0.01 mm。然而，在后续的零部件有限元分析中过小的网格尺寸会导致计算量巨大，因此需要确定合适的网格尺寸范围。本文以接触半宽为参照，接触部分网格尺寸／接触半宽（d／a）设置为0.02、0.05、0.1、0.2、0.5、1、2、3、4、5，得到10 种密度的网格模型，如图3 所示，求解得到最大接触应力随网格密度的变化，如图4所示。结果表明要得到可靠的最大接触应力，d／a应在0.5 以下。

图3 不同网格密度的有限元模型Fig.3 Finite element models with different mesh densities

图4 最大接触应力随网格密度变化Fig.4 Maximum contact stress varies with mesh density

2 燕尾榫连接结构的有限元分析

在明确了有限元法求解接触问题的流程和难点后，对某型号航空发动机燕尾榫连接结构进行了有限元分析。

2.1 结构模型及材料

燕尾榫连接结构模型如图5 所示，榫槽和榫头材料均为TC11，化学成分为Ti 基，常温下机械性能为：σ－1＝540 MPa，σb＝1 130 MPa，σs＝866 MPa，E＝120 GPa，μ＝0.33；TC11 钛合金材料拉伸的应力—应变曲线，如图6 所示［7］。

图5 燕尾榫连接结构模型Fig.5 Dovetail joint structure model

图6 TC11 拉伸的应力—应变曲线Fig.6 Stress－strain curve of TC11 tension

2.2 有限元计算

接触区域的应力应变非常复杂，为了更准确分析接触应力，本文在有限元模型的接触区域细化网格，远离接触区域网格设置较粗糙，以减少计算量如图7 所示。根据燕尾榫连接结构的实际工况设置载荷和约束，轮盘上施加15 KN 均布拉力模拟周向力，在榫头上端面施加10 KN 均布拉力模拟叶片产生的离心力，轮盘施加位移约束，摩擦系数0.5［8］。初始状态下接触面间距为0，接触力为0。计算结果表明，接触区域存在高应力梯度，该工况下接触应力最大值位于接触区域下边缘点，如图8、图9 所示。

图7 燕尾榫连接结构有限元模型Fig.7 Finite element model of dovetail joint structure

图8 应力分布云图Fig.8 Stress distribution nephogram

图9 接触区域应力分布Fig.9 Stress distribution in the contact area

2.3 计算结果的收敛性分析

根据计算结果已初步确定接触区域下边缘存在高应力梯度，且网格密度对有限元解产生影响。为保证所得到的有限元解有效，对接触区域下边缘依次分别细化为0.048 mm、0.024 mm、0.012 mm、0.008 mm、0.004 mm。求解得到不同网格密度下应力最大值，并计算各种应力最大值之间的误差，并根据公式（3）判断结果得收敛性［9］。即相邻两次模型的解误差在5%内即可认为收敛，不必再继续细化网格。各种应力最大值随网格密度的变化趋势如图10 所示，可以看出最大等效Mises 应力随网格加密而明显增大，但到一定尺寸后变化趋势减小；切向应力σxx、法向应力σyy以及剪切应力τxy同等效Mises应力变化趋势相同。不同模型计算误差分析结果见表1，表明随着网格密度的增大，最大应力值随之增大，计算结果逐渐收敛。网格1 同网格2 以及网格2 同网格3 的计算结果相差较大，且都超过5%，表明该网格密度下的计算结果没有收敛。网格3 和网格4 的计算结果相差4.3%，表明网格3 的计算结果已经达到网格收敛性要求。力分即若要得到准确的接触应布，求解模型的网格密度需大于或等于网格3。

图10 应力最大值随网格密度的变化趋势Fig.10 Trend of maximum stress with grid density

表1 不同尺寸网格模型计算误差分析Tab.1 Calculation error analysis of mesh models of different sizes

2.4 接触区域网格的优化

通过上文的分析计算，已经确定在榫连接接触区域下边缘存在高应力梯度，计算结果的精度受网格密度影响大。高密度网格模型直接导致计算量大幅度增加，需要对接触区域的网格进行优化。由于高应力梯度区域只存在于接触区域下边缘，因此只对接触区域下边缘进行网格加密。改进后结果相差不大，但是计算时间明显减少，对于后续处理更复杂的接触结构具有很大帮助。

3 接触应力影响因素的进一步分析

3.1 载荷对接触应力的影响

发动机在实际工作运转中存在多种工况，与之对应的榫连接结构也存在着不同的工况。不同载荷下的接触应力也会产生变化，为了分析不同工况下接触应力的变化，取10 KN、11 KN、13 KN、15 KN 4种工况载荷，接触面初始间距为0，摩擦系数取0.5，建立相同的模型进行分析。4 种工况下，接触应力都存在明显的峰值，最大接触应力值接近且都位于接触下边缘区，最大值所在区域相差不大。4 种工况下接触应力的分布和变化趋势一致，如图11 所示，接触应力随载荷增大而增大，沿接触面长度急剧增大达到峰值，随后急剧减小最后趋近稳定，4 种载荷下危险失效区域一致。

图11 不同工况下接触应力分布图Fig.11 Contact stress distribution under different working conditions

3.2 摩擦系数对接触应力的影响

摩擦系数也是影响接触应力的重要因素之一。为分析摩擦系数对接触应力的影响，在接触面设置不同摩擦系数的模型进行数值分析。钛合金材料在常温时摩擦系数一般在0.5 左右，因此摩擦系数取0.2、0.4、0.5、0.6、0.8，采用相同网格密度的模型，载荷为10 KN，初始状态接触间距为0。5 种不同的接触面摩擦系数下，接触应力的分布及变化趋势一致，如图12 所示，随着摩擦系数增大，接触面上接触应力反而减小，同时越靠近接触边缘，摩擦系数的影响越明显；随着摩擦系数的改变，最大接触应力出现区域变化不大，危险失效区域依旧在接触面下边缘。

图12 不同摩擦系数下接触应力分布图Fig.12 Contact stress distribution under different friction coefficients

4 结束语

本文对赫兹接触问题进行了有限元分析，明确了有限元分析接触问题的流程和难点，结果表明网格密度影响数值求解精度。要得到可靠结果，网格尺寸和接触半宽之比应在0.5 以下。对航空发动机燕尾榫连接结构接触应力进行分析，着重探究了网格密度对计算结果收敛性的影响及载荷和摩擦系数对接触应力的影响，得出结论：在燕尾榫连接结构接触面存在较高的应力梯度；常规的网格密度在高应力梯度求解结果的精度不够，需要优化高应力梯度位置网格；最大接触应力点位于接触区域下边缘，接触区域下边缘为危险失效区域；随着载荷和接触面摩擦系数的变化，接触应力分布趋势不变，数值随载荷增大而增大，随摩擦系数增大而减小。上述的高应力梯度区域的网格优化、网格收敛性的分析，对有限元法接触分析建模、平衡计算精度和计算量具有参考意义；榫槽在接触区域下边缘位置存在危险，对榫连接结构的设计和应用具有实际意义。