黄 伟，全 威，龚智伟，何继锟，王 敏

（广西大学机械工程学院，广西南宁530004）

大型工程机械在作业时，工作装置起到承受周期性变幅载荷的作用，在剧烈的冲击振动作用下，加快了结构件的磨损，最终使得应力较大的区域产生裂纹后会迅速扩展，导致各部件发生疲劳破坏。为了准确地预测机械结构件裂纹的位置和扩展寿命，以数值模拟与基础试验相结合的疲劳裂纹分析方 法［1-2］早 已 被 广 泛 提 及。Rusiński等［3］和Djurdjevic等［4］采用有限元法分别对矿山机械已开裂动臂和斗轮挖掘机含裂纹臂架进行了数值模拟。最后对比两种结果，后者对水平和垂直方向进行额外的加载模拟，数值分析后表明，水平力是引发初始裂纹萌生的主要因素。樊俊铃等［5］通过在ABAQUS中建立内聚单元模型来研究Ι型裂纹的扩展过程，再采用数值模拟的方法去验证裂尖温度、裂纹扩展速率和应力强度因子幅在裂纹扩展过程中相互之间的函数关系。上官文斌等［6］进行了在变幅加载工况下隔振器橡胶材料的裂纹扩展实验，采用幂函数插值法来计算裂纹扩展速率，从而确立橡胶材料的裂纹扩展模型。

本文以某国产装载机动臂为研究对象，综合考虑动臂的承载情况，通过有限元法在动臂上建立不同的初始裂纹模型，以研究在裂纹扩展过程中应力强度因子值的变化规律。考虑到裂纹扩展快慢的变化，采用数值拟合的方法来建立裂纹扩展速率与深度的函数关系，进而提出了基于Paris公式的裂纹扩展寿命预测方法。

1 动臂模型的建立与仿真

1.1 建立装载机工作装置的刚柔耦合模型

首先建立装载机工作装置各部件的三维模型，进行整体装配后，再以x_t格式将其导入进ADAMS中。此时铲斗与动臂处于初始启动位置如图1所示。

图1 工作装置的刚柔耦合模型电路图Fig.1 Rigid flexible coupling model of working device

给定各部件的密度和质量，同时确保质心位置、各部件间的约束关系与实际情况一致。在不考虑运输过程中保持静止的时间后，历时20 s。通过STEP函数［7］在铲斗齿尖中心处施加3个外力用来模拟实际工况下铲斗的受力情况，分别为物料重力Fig、铲斗插入阻力Fin、转斗阻力Ftn；同时施加两个驱动力，分别控制着转斗油缸和动臂油缸，用来实现铲斗从水平插入物料到最终卸料复位的各阶段性工况。其函数方程为：STEP（time，起始时间，0，终止时间，数量值）。

建立动臂柔性体时，划分网格选取soild187单元。设置动臂铰接孔中心节点为外部节点时，首先用质量单元mass21建立4个质量单元，然后设置在同一个中心平面上，分别对应8个铰接孔的刚性区域如图2所示。在生成MNF模态中性文件后将其替换ADAMS中原有的刚性体动臂，完成动臂的柔性化处理。

图2 动臂柔性化模型Fig.2 Flexible model of boom

1.2 获取各铰接点的载荷历程

将动臂经过柔性化处理后，再对从铲斗水平铲入到最终卸料复位的整个过程进行模拟仿真。然后分别提取动臂与铲斗、摇臂铰接处的载荷历程，如图3和图4所示。

图3 动臂与铲斗铰接点处载荷历程Fig.3 Load history at the joint of boom and bucket

图4 动臂与摇臂铰接点处载荷历程Fig.4 Load history at the joint of boom and rocker

2 动臂危险点应力时程的建立

在整个仿真过程中，由于动臂的位置是随时间变化的，不同位置和时刻的受力大小及方向都是不一样的，所以所得到的载荷谱数据必须经过换算才能施加到Workbench中。已知ADAMS与Workbench中所施加的载荷两者对应的关系为

式中：F、F'分别为Workbench和ADAMS对应方向的力；α为相应力的夹角。

从图3和图4中可以看出，在5.5 s时两处铰接点受力最大，此时恰好为铲斗插入物料并要翻转的瞬间。按照式（1）换算后得到动臂与铲斗处F'x为416 000 N，F'z为101 000 N，动臂与摇臂处F'x为655 000 N，F'z为21 976 N。将此时的动臂设好约束关系后进行静力学分析如图5所示。

图5 5.5 s时动臂应力图Fig.5 Boom stress diagram at 5.5 s

从图5可以看出，在动臂与摇臂铰孔处，σmax为178.73 MPa。该位置在应力集中作用下，容易产生疲劳裂纹。因此，本文考虑铰孔内壁的裂纹容易被忽视的实际情况，将此处确定为动臂进行裂纹分析的危险位置。按照此方法将所得的载荷历程分为200步，历时20 s，逐一进行各个时刻的静力学计算，最终得到危险点的应力时程曲线如图6所示。

图6 动臂应力时程图Fig.6 Stress time history of boom

3 动臂裂纹扩展仿真计算

3.1 有限元法确定应力强度因子

多数含裂纹结构件的尺寸大小、裂纹位置以及裂纹本身的几何尺寸或承载情况都十分复杂，相比较传统的解析法，采用有限元法来确定应力强度因子值则更为简单。基于变分原理的有限元法，能将整体结构进行离散化，用插值法去分析离散化后各个单元内的相互力学作用，通过有限矩阵形式去计算结构件的整体行为。常见的用来求解应力强度因子的有限元法包括裂纹张开位移法、应力法、虚裂纹闭合法和虚裂纹扩展法。

根据线弹性断裂力学理论，按照裂纹的受力特点，可以将裂纹的基本类型抽象化分为Ⅰ型（张开型裂纹）、Ⅱ型（滑移型裂纹）和Ⅲ型（撕开型裂纹），而3种裂纹尖端区域位移与应力都含有的奇异性［8］。在求解裂纹尖端应力场或位移场时，在求解精度方面，裂纹尖端的位移法要比应力法更好，其应力强度因子位移法计算公式［9］为

式中：Δυ、Δu、Δω分别为由裂纹分开的两个平面的X、Y、Z方向的相对位移；E为材料的弹性模量；G为剪切模量；r为极坐标半径；v为应变；k为描述疲劳裂纹的基本参数。

在工程结构中，最容易导致疲劳断裂的是Ⅰ型，若结构中同时存在两种及两种以上的复合型裂纹时，出于安全考虑，往往将其看作Ⅰ型裂纹进行研究［10］。因此，在本文中，主要采用Workbench的Fracture工具求解应力强度因子KΙ值，将动臂的裂纹扩展简化为纯Ⅰ型裂纹扩展进行分析。

3.2 裂纹模型的建立

动臂材料采用的是Q345钢，其材料性能如表1所示［11］。

表1 Q345钢的材料属性Tab.1 Material properties of Q345 steel

在动臂与摇臂铰接孔内壁上建立裂纹模型时，为了使结果更可靠，需要对裂纹尖端区域合理地划分网格，裂尖的网格要划分得足够细，以模拟裂纹尖端的特性。以裂纹c=20 mm、a=10 mm为例，通过Fracture工具的三维四面体划分网格后自动生成裂纹，如图7所示。

图7 铰孔处裂纹的有限元模型Fig.7 The finite element model of cracks in r eamed holes

计算后可以发现在初始裂纹裂尖出现了明显的应力集中，周围应力场呈现“蝴蝶状”分布。而应力强度因子最大值出现在裂纹尖端附近，，如图8所示。

图8 初始裂纹的应力强度因子值及应力场Fig.8 Stress intensity factor value and str ess field of initial crack

3.3 确定应力强度因子幅表达式

在计算裂纹扩展寿命时，工程结构件所承受的载荷谱和疲劳裂纹扩展速率是必不可少的。Paris提出应力强度因子不仅是表征裂纹尖端附近应力、应力场的主要参量，而且是控制裂纹扩展速率的主要参量。由此，提出Paris公式：

式中：N为循环周次；C为描述疲劳裂纹扩展性能的基本参数；ΔK=Kmax-Kmin、Kmax和Kmin分别为一个应力循环内的最大应力σmax和最小应力σmin所对应的应力强度因子值。

由图6可知,在各个时刻的等效应力都不同,以两个不同时刻的应力差值为应力幅,则会计算出不同的裂纹扩展寿命。因此,取动臂一个循环周期内的应力最大值及最小值,以此为差值的应力幅来计算动臂的裂纹扩展寿命会较小,但更可靠[12]。因此在本文中,对c=20 mm、a=10 mm的初始裂纹取σmax=178.73 MPa、σmin=10.31 MPa时所对应的时刻,求得ΔK为

分别建立半椭圆形长半轴c为15、20、25、30 mm，短半轴（深度）a为2、3、4、5、8、10、15、20、25、30、40、50、60 mm的裂纹扩展模型。按照上述方法计算应力强度因子幅值ΔK，单位，整理后如表2所示。

表2 不同裂纹模型的应力强度因子数值Tab.2 Stress intensity factors of different crack models mm

通过Matlab对上述数据进行3次多项式拟合后，分别得到c为15、20、25、30 mm时的ΔK表达式：

将其绘制如图9所示。

图9 ΔK-a的关系曲线Fig.9 The relation curve ofΔK-a

可以看出，不同的初始裂纹长度，其应力强度因子幅ΔK随着裂纹扩展过程中深度a的增加逐渐增大，而裂纹初始长度c越大，则裂纹应力强度因子幅ΔK增长速度越快。

3.4 疲劳裂纹寿命的预测

基于断裂力学理论，将裂纹扩展过程划分为3个区域阶段。第Ι区为不扩展区，此时ΔK＜ΔKth。第ΙΙ区为裂纹稳定扩展区，是估算疲劳裂纹寿命的主要部分，此时ΔK≥ΔKth。第ΙΙΙ区为裂纹急速扩展区，裂纹扩展速率急速上升而产生瞬间断裂，在工程中可以忽略不计。

当处于第ΙΙ区时，通常采用Paris式（5）来描述，它是疲劳裂纹扩展分析中应用最广泛的裂纹扩展模型，将其变换形式后得

式中：a0为裂纹初始深度，在建立裂纹模型时已确定；ΔK为裂纹尖端应力强度因子幅；C、m为Paris经验公式系数，本文中取C=1.17×10-9，m=4.47［13］；ac为临界裂纹尺寸，由线弹性断裂力学公式所确定：

本文中KIC取，裂纹结构的几何形状因子F取1.12，最后确定ac=68.8 mm。从而计算出初始裂纹长度为20 mm时，裂纹深度a0增长到ac，动臂铰孔内壁的裂纹扩展寿命。

依此方法，分别计算出在不同的初始裂纹长度c时，裂纹深度随着循环次数的变化，将结果通过Matlab拟合后得到动臂危险点处的裂纹扩展曲线如图10所示。

从图10中可以看出，在4种初始裂纹长度不同、裂纹深度扩展到30 mm时，裂纹扩展速率都显著增加，并逐渐趋近于临界裂纹深度，即将进行断裂失效。

图10 裂纹扩展循环次数曲线Fig.10 Curve of crack growth cycles

4 结论

（1）首先通过有限元法确定裂纹模型，然后采用数值拟合的方法绘制裂纹扩展过程中的相关图表，最后基于Paris公式进行裂纹疲劳寿命的计算。与传统的断裂力学解析法相比，该方法简化了分析机械结构件寿命的整个过程，能够针对性地发现结构缺陷并快速分析，其结果更偏于保守。

（2）Paris公式能够准确地解释在单个循环载荷周期下装载机动臂的疲劳裂纹扩展情况。裂纹扩展模拟仿真和基于Paris公式的分析结果都表明，裂纹扩展速率随应力强度因子幅的增大而增大，且随着裂纹的初始长度的增大而增大。

（3）从裂纹扩展循环次数曲线可以看出，不同长度的初始裂纹的循环寿命都在临近30 mm的深度之前趋于缓慢上升的趋势，之后急速增加。因此，当动臂铰孔位置处裂纹扩展到将近30 mm时，应当尽快采取相应的维护措施。