尖点应力计算的有效性方法

张凤姣，宋少云

(武汉轻工大学机械工程学院，湖北武汉 430023)

随着各国经济飞速的发展，有限元仿真技术开始在各行各业中得到广泛应用［1］。尤其在汽车行业，越来越多的研究者开始对汽车的结构及其零部件进行有限元仿真［2-4］，以考察其力学性能，为设计师们对其结构的设计和优化提供理论依据。

在进行结构分析时，应力集中是一种很棘手的问题［5-6］。所谓应力集中是指受力构件由于外界因素或自身因素几何形状、外形尺寸发生突变而引起局部范围内应力显著增大的现象。笔者在用有限元软件对汽车手刹进行强度分析时，发现对于某些零件，当网格不断细分时，部分尖点会出现应力值无限增大，不会趋于某个值，无法得到确切的应力值。

实际上，有限元理论表明，在尖点(形状突变处)处应力集中系数会无限增大，无法计算真实的应力解［7］。同时，通过有限元仿真的实践也表明，当网格无限细分时，尖点应力的确会无限增大，得不到收敛解。因此，使用有限元软件无法得到尖点处的真实应力值。

为了能够正确地计算尖点处的真实应力值，本文提出了三种计算尖点处应力的方法:圆周节点应力平均法，线段应力外推法和不同倒角应力外推法。并用上述三种方法对一个算例进行计算，把计算结果与精确解相比较，从而考察方法的有效性。

1 算例及精确解

为了能清晰地阐述本文提出的三种方法，先给出一个共同的算例(如图1所示)。该图是一个7字形截面构件，其左边沿和下边沿均固定，而在右边上施加100 MPa的均布载荷，现在要计算倒角处(倒角半径为2 mm)的米塞斯应力值。

图1 7字型构件几何尺寸图

对于这种问题，有限元软件是可以计算其精确应力的。为了得到该点的精确应力，首先对该几何体进行粗糙的网格划分，固定左、下端面并在右端面施加均布载荷，然后进行静力学仿真。在仿真结束后，记录下倒角处的米塞斯应力。

然后对该倒角处进行网格局部细分，连续加密8次，并记录该点的应力值。这样计算7次后，得到该点应力的8个值。取仿真次数为横坐标，而取每次仿真得到的米塞斯应力为纵坐标，得到的曲线如图2所示。从该图可以看出，在第五次仿真后，该点应力趋于收敛。第8次计算结果为120.31MPa，可以认为它是精确解。

图2 倒角处点的米塞斯应力值

2 计算方法

对于上述问题，有限元软件可以得到拐角处的精确解。但是当拐角处倒角半径为0时(尖点)，用有限元软件计算，会发现结果会一直增大，得不到收敛值。但是就实际情况而言，该点的值是客观存在的，而且一定是一个收敛值。那么如何通过计算的方式得到该值呢?

有限元软件对于尖点的应力值求解是发散的，但是对于尖点一定距离处的点的应力计算却是收敛的，而尖点的应力与周围点的应力应该是连续的。基于这种考虑，可以根据尖点附近点的应力来推算尖点处点的应力值。

而从周围点来估算尖点的应力，有很多方法。本文提出三种方法来处理这种问题。

第一种方法:圆周应力平均法。

此方法是以倒角处的点为圆心，分别以1、1.5、2、2.5为半径作圆，分别细分网格6次求得应力曲线图，并考察这些圆周上的应力值，最后平均，求得倒角处应力值。

第二种方法:线段外推法。

此方法是以道交处的点为中心，画6条偏角不同的线段，分别为 30°、45°、60°、90°、120°、150°。对这些线段分别细分网格7次求得应力曲线图，根据这些曲线的走势，选取最收敛的那条曲线进行插值外推中心点的应力值，即所求的倒角处的点应力值。

第三种方法:变倒角半径外推法。

此方法是将尖点处倒不同的角，倒角半径分别为0.5 mm、0.6 mm、0.7 mm、0.8 mm、1 mm、1.2 mm;6种情况下计算尖点处的应力值。将得到的6组数据进行分析，选取各组中最收敛的那条曲线进行三次样条曲线插值外推中心点的应力值，再与精确值作比较，最后得所求的尖点应力值。

通过三种分析方法所得的结果进行比较，最后可以选择一种最好的办法来求解尖点的应力值。

2.1 圆周应力平均法

以倒角处点为圆心，如图3所示。首先以半径为1 mm的圆周开始，绘制出四个圆，半径大小依次增大0.5 mm，并在圆周上6次细分网格，其中半径为0的点即为本文要估算的点。

图3 尖点圆心处四个不同半径的圆

应用有限元软件对上述所绘的4个圆周进行有限元分析，并且依次细分网格6次，得到的应力变化曲线如图4—图7所示(图中横坐标S代表的是圆周上所有节点距离第一个节点的圆弧长度，B.C.D等字母表示的是划分不同网格的次数)。

图4 半径为1 mm圆周上的点应力变化曲线

图5 半径为1.5 mm圆周上的点的应力变化曲线

图6 半径为2 mm圆周上的点的应力变化曲线

图7 半径为2.5 mm圆周上的点的应力变化曲线

由图4—图7可以清楚看出4个不同半径的圆周上所有节点的应力变化情况。根据有限元原理，应用平均值法可以从各圆周上的应力值求得倒角处点的4种应力结果，然后将得出的这些应力与之前的精确值相比较得出相对误差，如表1所示。

表1 圆周推算的应力值

2.2 线段应力外推法

如图8所示，以计算的拐点为圆心，在直径为5mm和直径为20 mm的圆环区域内依次取6条线段，即从 30°开始，分别以 30°、45°、60°、90°、120°、150°为偏角，取6条线段。并且依次7次细分网格，最后读取计算结果。结果如图9—图14所示(图中的S代表的是线段上所有节点到起始第一个点的距离，B.C.D等字母表示的是划分不同网格的次数)。

图8 6条偏角不同的线段

图9 30°线段上的点的应力变化曲线

图10 45°线段上的点的应力变化曲线

图11 60°线段上的点的应力变化曲线

图12 90°线段上的点的应力变化曲线

图13 120°线段上的点的应力变化曲线

图14 150°线段上的点的应力变化曲线

根据这些曲线进行分析，发现30°偏角时，第四次网格细分的时候，应力收敛，通过第四次的曲线进行插值，最后用外推法可以求得尖点处的应力值。同理可得，45°、60、90°、120°、150°的时候均是第四次或者第三次就已经应力收敛，相应的可以插值求得应力值，并与精确值相比较，得到相对误差，其结果如表2所示。

表2 线段推算的应力值

2.3 变倒角半径应力外推法

对图1尖点处(图中倒角R2处)进行6次不同半径的倒角处理，倒角半径大小依次为0.5 mm、0.6 mm、0.7 mm、0.8 mm、1 mm、1.2 mm。

在取6次不同倒角的情况下，均可以通过网格细分得到尖点的应力值，经仿真，其结果如表3所示。

表3 不同倒角推算的应力值

根据表3给出的数据，应用曲线插值尖点应力，并绘制米塞斯应力值与尖点处进行倒角半径值的曲线关系图，如图15所示(图中横坐标代表的是尖点处倒不同倒角的半径大小，纵坐标是米塞斯应力值)。

图15 尖点处不同倒角时的拐点应力值

根据六组将尖点进行不同倒角时的拐点米塞斯应力值，利用上图给出的数据，根据三次样条曲线插值尖点的应力值，得出其值为154.90 MPa。

3 结果分析及讨论

3.1 从表1的数据结果可以得知，当半径不断增大时，试验所得的应力值越来越偏离精确值，但是当半径为1.5 mm的时候，应力值却最接近精确值。半径由小及大远离和接近1.5 mm时，试验值越来越偏离精确值，且误差越来越大。

3.2 从表2的数据结果可以得知，当角度越大时，试验所得的应力值数据越来越接近精确值，当偏角为90°的时候，应力值却最接近精确值，误差最小。但是当偏角大于90°且不断增大的时候，试验值越来越偏离精确值，且误差越来越大。可以看出，沿尖点纵向方向，尖点应力变化梯度确实是应力变化最大的地方。

3.3 从表3的数据结果可以得知，倒角越小越接近真实值，但是当倒角一定小的时候，应力值却无限大，反而计算误差很大。当用样条曲线插值的时候，发现0半径即尖点处的应力值为154.90 Mpa。

3.4 综合比较上述三种方法的计算结果，发现用第一种方法中半径为1.5 mm的圆周节点应力值推算的尖点应力值是最接近真实值的;用第二种方法中90°偏角线段的时候同样也是可以推算的，只是误差稍微大点。当然，第三种方法误差最大，此法不宜。

4 结论

4.1 就此例子而言，第一种方法中，当半径为1.5 mm的圆时，用圆周上节点的应力值去加权平均，得到的应力值是最接近尖点真实值的，误差最小。

4.2 可以通过对倒角附近一定范围的节点的收敛应力，用加权平均或者外推的方法来估算拐角处点的应力值。

4.3 一般情况下，越接近尖点就越能得到收敛的应力值。

［1］王岩.干式切削加工有限元仿真技术研究［D］.大连:大连交通大学，2010.

［2］何大为.基于CVT的混联式混合动力汽车控制与仿真研究［D］.长沙:湖南大学，2008.

［3］亓文果，金先龙，张晓云，等.汽车车身碰撞性能的有限元仿真与改进［J］.上海大学交通学报，2005，39(09):1452-1457.

［4］张瑞妍.汽车起动机减速轴精密成形过程仿真及应用研究［D］.镇江:江苏大学，2009.

［5］于瀚翔，于兰峰.基于有限元法的桁架节点搭接接头应力集中分析［J］.机械设计与制造，2013(3):14-16.

［6］李杨，邓金根，蔚宝华.套管内壁球形腐蚀坑洞的应力集中效应研究［J］.石油机械，2012，40(9):73-77.

［7］唐洪祥，管毓辉.孔口应力集中问题的Cosserat连续体有限元分析［J］.东南大学学报，2013，43(04):849-855.