裂纹扩展的扩展有限元xfem模拟实例详解.doc

基于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。
(1)
上式中a，b分别是裂纹体的几何尺寸，F为关于a和b比值的函数，可以查表得到，本文中a与b的比值为0.5，查表得到F的函数值为1.50。

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

本文计算值
理论计算值

58.69
58.28

由此可以计算相对误差为：1.06%，此误差显然属于5%的允许误差范围之内。所以本文计算得到的应力强度因子是可信的。
此外，本文还利用此模型计算了静态裂纹的J积分值，由于材料是理想线弹性的，所以J积分与应力强度因子之间存在这样关系，如公式（2）。
(2)
本文在这里列出J积分的变化趋势图,图中对比了公式（2）的理论解以及有限元结果。从图中可以看出，两者是吻合的，说明了有限元模拟是正确的。J积分随加载的变化趋势图如图5所示。红色实线表示的是理论结果，黑色点表示有限元结果。

图5 J积分历史曲线图

3.2 裂纹扩展过程展示

Step time=0.1143
Step time=0.3943
Step time=0.6743

Step time=0.7976
Step time=0.8708
Step time=0.9551

Step time=0.9994
Step time=0.9998
Step time=1
图6 裂纹扩展过程

从上述的裂纹扩展过程的应力分布图，我们可以得到如下几点结论，证明我们的数值模拟具有一定的正确性。首先，在裂纹尖端出现了应力集中，这是和断裂力学理论符合的。其次，观察裂纹附近的应力分布，我们可以看到应力分布的趋势是与理论计算的塑性区的形状大致相同的，理论计算的塑性区形状如图7所示。

图7 理论上的塑形区形状图

3.3 裂尖单元应力变化
其次，考察裂尖单元的应力随载荷增加的变化。实际上裂尖单元应力值的具体大小并没有意义，因为表征断裂韧强度的是应力强度因子和J积分。而单元应力随载荷增加导致的变化可以帮助我们理解裂尖单元在起裂到完全断裂的整个过程。
观察图8，我们可以大致得到这个裂尖单元参与断裂过程的整个历史。首先，在应力加载的早些时候，裂尖单元的应力随着载荷的增加而增加，此时裂尖单元的应力并没有达到损伤判据的临界应力，所以单元没有起裂。随着载荷的继续增加，应力值继续增加，当到大概0.7976左右时，裂尖的最大主应力达到了最大主应力损伤判据的临界值，于是裂纹起裂，直至完全裂开，单元的应力奇异性消失，裂尖单元转变为一个普通的非裂尖单元。这个过程对应于图8中的右边应力增大后有急剧减小的曲线。之后，由于载荷还没有完全加载完毕，所以裂尖单元在转变为普通单元之后随着载荷的继续增加，其单元应力又会随之在增加。
以上就是一个裂尖单元在整个加载过程中的单元应力历史变化的三个阶段。

图8 裂纹尖端单元应力历史曲线

3.4 裂纹扩展分析
如图9所示，载荷从零开始不断加载。随着载荷的增加，裂尖处的单元应力也不断增加。当裂尖处的最大主应力值达到临界值时，裂尖处的单元开始失效，裂纹开始穿过单元扩展，时间步大概是0.7976左右，此时裂尖处开始形成粘结裂纹(cohesive crack)。从0-0.7976这个过程可以认为是裂纹孕育期。之后随着裂尖处的能量释放率达到裂纹扩展阻力率GC时，裂尖处的粘结裂纹开始扩展成真实裂纹，裂尖单元的XFEM值达到1。裂尖处的单元损伤值达到临界值，时间步大概是0.9996左右。此时认为这一结构开始失效，裂纹失去平衡，开始失稳扩展，真实裂纹开始形成并不断扩展。从0.7976-0.9996这个过程可以认为是裂纹的萌生过程。0.9996以后裂纹失稳扩展，结构失效。从图中可以看出初始裂纹长度为12个单元距离，到最终加载结束时，真实裂纹长度为15个单元长度，粘结裂纹长度为7个单元长度。

图9 随加载历史扩展的动态裂纹图

3.5 裂纹扩展的能量历史曲线
绘出整个模型的总能量，动能，内能和外力功随计算过程的历史曲线如图10所示。由图可以看出，总能量和动能在整个过程中并没有发生显著的增加，可以表明整个计算过程基本是稳定的。而伪应变能在计算过程中有稍微的增长。

图10 裂纹扩展的能量历史曲线

4 结论
本文的工作是基于abaqus6.9版本的扩展有限元功能计算了Ⅰ型裂纹的扩展。得到如下几点结论。
首先，本文计算了Ⅰ型裂纹的应力强度因子，结构表明计算值与理论值的误差在5%以内，结果可靠有效。同时计算了裂纹的J积分值，获得了J积分的历史曲线，有限元结果与理论结果吻合。
裂纹的扩展过程与理论是吻合的，整个裂纹扩展的计算没有出现不稳定的情况。裂尖出现了应力的奇异性，裂纹扩展平稳。并通过研究裂尖单元的应力历史曲线，直观的获得了裂纹在一个单元上扩展的三个典型阶段。
其次，对裂纹扩展的过程进行了分析，指出了裂纹扩展的大概几个阶段以及具体的裂纹扩展过程。
最后，对裂纹扩展过程中的能量变化作了简单描述，说明了计算过程的稳定性，并验证了系统的能量守恒关系。