裂纹扩展的扩展有限元(xfem)模拟实例详解

基于ABAQUS 扩展有限元的裂纹模拟

化工过程机械622080706010 李建

1 引言

1.1 ABAQUS 断裂力学问题模拟方法

在abaqus中求解断裂问题有两种方法（途径）：一种是基于经典断裂力学的模型；一种是基于损伤力学的模型。

断裂力学模型就是基于线弹性断裂力学及其基础上发展的弹塑性断裂力学等。如果不考虑裂纹的扩展，abaqus可采用seam型裂纹来分析(也可以不建seam，如notch型裂纹)，这就是基于断裂力学的方法。这种方法可以计算裂纹的应力强度因子，J积分及T-应力等。

损伤力学模型是指基于损伤力学发展而来的方法，单元在达到失效的条件后，刚度不断折减，并可能达到完全失效，最后形成断裂带。这两个模型是为解决不同的问题而提出来的，当然他们所处理的问题也有交叉的地方。

1.2 ABAQUS 裂纹扩展数值模拟方法

考虑模拟裂纹扩展，目前abaqus有两种技术：一种是基于debond的技术（包括VCCT）；一种是基于cohesive技术。

debond即节点松绑，或者称为节点释放，当满足一定得释放条件后（COD 等，目前abaqus提供了5种断裂准则），节点释放即裂纹扩展，采用这种方法时也可以计算出围线积分。

cohesive有人把它译为粘聚区模型，或带屈曲模型，多用于模拟film、裂纹扩展及复合材料层间开裂等。cohesive模型属于损伤力学模型，最先由Barenblatt 引入，使用拉伸-张开法则（traction-separation law）来模拟原子晶格的减聚力。这样就避免了裂纹尖端的奇异性。Cohesive 模型与有限元方法结合首先被用于混凝土计算和模拟，后来也被引入金属及复合材料。Cohesive界面单元要服从cohesive 分离法则，法则范围可包括粘塑性、粘弹性、破裂、纤维断裂、动力学失效及循环载荷失效等行为。

此外，从abaqus6.9版本开始还引入了扩展有限元法（XFEM），它既可以模拟静态裂纹，计算应力强度因子和J积分等参量，也可以模拟裂纹的开裂过程。被誉为最具有前途的裂纹数值模拟方法。本文将利用abaqus6.9版本中的扩展有限元法功能模拟常见的Ⅰ型裂纹的扩展。

2 Ⅰ型裂纹的扩展有限元分析

本文针对断裂力学中的平面Ⅰ型裂纹扩展问题用abaqus中的扩展有限元方法进行数值模拟，获得了裂纹扩展的整个过程，裂尖单元的应力变化曲线，以及裂纹尖端塑性区的形状。在此基础上绘制裂纹扩展的能量历史曲线变化趋势图。

2.1 平面裂纹的几何模型

几何模型的尺寸参数如图1所示，其中a=1.5m，b=3m，L=10m，厚度为1m。上下两端分别承受25.32MPa的拉力。

图1 裂纹的几何示意图

2.2 有限元模型

有限元程序采用大型通用ABAQUS6.9软件，选用8节点六面体减缩单元（C3D8R）。网格划分的模型如图2所示。

图2 网格图

2.3 材料性能

在有限元分析中假定材料为理想线弹性的，弹性模量E为2.1×105MPa，泊松比为0.3。本文采用的是基于损伤力学演化的失效准则。具体的参数设置如下。损伤判据为最大主应力失效准则作为损伤起始的判据，最大主应力为84.4MPa。损伤演化选取基于能量的、线性软化的、混合模式的指数损伤演化规律，有关参数为G1C= G2C= G3C=42200N/m，α=1。

2.4 边界条件和初始条件

对于含有裂纹的平板，我们仅仅需要约束住它的刚体位移，保证在在平板两个断面施加应力载荷时，平板不会出现意外的刚体运动。设置裂纹类型为扩展有限元裂纹，扩展区域是整个平板，扩展路径为任意路径。由于计算裂纹扩展实际上是一个大变形问题，所以分析步骤的几何非线性一定要打开。由于裂纹扩展本身是一个强烈的非连续问题，它将导致求解过程的迭代有可能出现不收敛的情况，另外，求解的增量步也会要求很小，这会导致求解时间很长。因此非常有必要对求解过程做一些参数控制，以避免迭代不收敛导致的求解失败的情况的出现。

图3 裂纹体及其扩展区域图4 载荷及边界条件

3 结果分析

3.1 静态裂纹的应力强度因子及J积分的验证

计算应力强度因子及J积分时，需要设置裂纹不能扩展，从而计算静态裂纹的应力强度因子，同时还要在历史变量输出中做相关的设置。另外分析步也需要将几何非线性去除，因为裂纹没有扩展。由此计算得到了裂尖在25.32MPa载荷下的Ⅰ型应力强度因子。同时，我们根据断裂力学理论中关于此模型的理论解如公式（1），计算理论的应力强度因子。最后得到的结果列于表1。

⎪

⎭⎫

⎝⎛=I b a F a K πσ (1)

上式中a ，b 分别是裂纹体的几何尺寸，F 为关于a 和b 比值的函数，可以查表得到，本文中a 与b 的比值为0.5，查表得到F 的函数值为1.50。

表1 应力强度因子的对比表

本文计算值 理论计算值 )(2/1m MPa K ⋅I

58.69

58.28

由此可以计算相对误差为：1.06%，此误差显然属于5%的允许误差范围之内。所以本文计算得到的应力强度因子是可信的。

此外，本文还利用此模型计算了静态裂纹的J 积分值，由于材料是理想线弹性的，所以J 积分与应力强度因子之间存在这样关系，如公式（2）。

E

K J 2

I = (2)

本文在这里列出J 积分的变化趋势图,图中对比了公式（2）的理论解以及有限元结果。从图中可以看出，两者是吻合的，说明了有限元模拟是正确的。J 积分随加载的变化趋势图如图5所示。红色实线表示的是理论结果，黑色点表示有限元结果。

0.0

2.0x10

3

4.0x103

6.0x103

8.0x103

1.0x104

1.2x104

1.4x10

4

1.6x104

1.8x104

J 积分/N m -1

Step time

图5 J 积分历史曲线图