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

基于ABAQUS扩展有限元的裂纹模拟化工过程机械 622080706010 李建1引言1.1 ABAQUS断裂力学问题模拟方法在abaqu s中求解断裂问题有两种方法（途径）：一种是基于经典断裂力学的 模型；一种是基于损伤力学的模型。断裂力学模型就是基于线弹性断裂力学及其基础上发展的弹塑性断裂力学 等。如果不考虑裂纹的扩展，abaqus可采用seam型裂纹来分析（也可以不建seam， 如notch型裂纹），这就是基于断裂力学的方法。这种方法可以计算裂纹的应力强 度因子，J积分及T-应力等。损伤力学模型是指基于损伤力学发展而来的方法，单元在达到失效的条件 后，刚度不断折减，并可能达到完全失效，最后形成断裂带。这两个模型是为解 决不同的问题而提出来的，当然他们所处理的问题也有交叉的地方。1.2 ABAQUS裂纹扩展数值模拟方法考虑模拟裂纹扩展，目前abaqus有两种技术：一种是基于debond的技术（包 括V CCT）； 一种是基于cohesive技术。debond即节点松绑，或者称为节点释放，当满足一定得释放条件后（COD 等，目前abaqus提供了5种断裂准则），节点释放即裂纹扩展，采用这种方法时也 可以计算出围线积分。cohesive有人把它译为粘聚区模型，或带屈曲模型，多用于模拟film、裂纹 扩展及复合材料层间开裂等。cohesive模型属于损伤力学模型，最先由Barenblatt 引入，使用拉伸-张开法则（traction-separation law）来模拟原子晶格的减聚力。 这样就避免了裂纹尖端的奇异性。Cohesive模型与有限元方法结合首先被用于 混凝土计算和模拟，后来也被引入金属及复合材料。Cohesive界面单元要服从 cohesive分离法则，法则范围可包括粘塑性、粘弹性、破裂、纤维断裂、动力学 失效及循环载荷失效等行为。此外，从abaqus6.9版本开始还引入了扩展有限元法（XFEM），它既可以模 拟静态裂纹，计算应力强度因子和J积分等参量，也可以模拟裂纹的开裂过程。 被誉为最具有前途的裂纹数值模拟方法。本文将利用abaqus6.9版本中的扩展有限 元法功能模拟常见的I型裂纹的扩展。2I型裂纹的扩展有限元分析本文针对断裂力学中的平面1型裂纹扩展问题用abaqus中的扩展有限元方 法进行数值模拟，获得了裂纹扩展的整个过程，裂尖单元的应力变化曲线，以及 裂纹尖端塑性区的形状。在此基础上绘制裂纹扩展的能量历史曲线变化趋势图。2.1平面裂纹的几何模型几何模型的尺寸参数如图1所示，其中i=1.5m, b=3m, L=10m，厚度为1m。上下两端分别承受25.32MPa的拉力。图1裂纹的几何示意图2.2有限元模型有限元程序采用大型通用ABAQUS6.9软件，选用8节点六面体减缩单元（C3D8R）。网格划分的模型如图2所示。2.3材料性能在有限元分析中假定材料为理想线弹性的，弹性模量E为2.1x105MPa,泊松 比为0.3。本文采用的是基于损伤力学演化的失效准则。具体的参数设置如下。 损伤判据为最大主应力失效准则作为损伤起始的判据，最大主应力为84.4MPa。 损伤演化选取基于能量的、线性软化的、混合模式的指数损伤演化规律，有关参 数为 G1C= G2C= G3C=42200N/m，o=1。2.4边界条件和初始条件对于含有裂纹的平板，我们仅仅需要约束住它的刚体位移，保证在在平板两 个断面施加应力载荷时，平板不会出现意外的刚体运动。设置裂纹类型为扩展有 限元裂纹，扩展区域是整个平板，扩展路径为任意路径。由于计算裂纹扩展实际 上是一个大变形问题，所以分析步骤的几何非线性一定要打开。由于裂纹扩展本 身是一个强烈的非连续问题，它将导致求解过程的迭代有可能出现不收敛的情 况，另外，求解的增量步也会要求很小，这会导致求解时间很长。因此非常有必 要对求解过程做一些参数控制，以避免迭代不收敛导致的求解失败的情况的出 现。图3裂纹体及其扩展区域图4载荷及边界条件3结果分析3.1静态裂纹的应力强度因子及J积分的验证计算应力强度因子及J积分时，需要设置裂纹不能扩展，从而计算静态裂纹 的应力强度因子，同时还要在历史变量输出中做相关的设置。另外分析步也需要 将几何非线性去除，因为裂纹没有扩展。由此计算得到了裂尖在25.32MPa载荷 下的I型应力强度因子。同时，我们根据断裂力学理论中关于此模型的理论解如 公式（1），计算理论的应力强度因子。最后得到的结果列于表1。上式中a, b分别是裂纹体的几何尺寸，F为关于a和b比值的函数，可以查表 得到，本文中a与b的比值为0.5，查表得到F的函数值为1.50。表1应力强度因子的对比表本文计算值理论计算值Ki (MPa - mi/2)58.6958.28由此可以计算相对误差为：1.06%，此误差显然属于5%的允许误差范围之内。 所以本文计算得到的应力强度因子是可信的。此外，本文还利用此模型计算了静态裂纹的J积分值，由于材料是理想线弹 性的，所以J积分与应力强度因子之间存在这样关系，如公式（2）。本文在这里列出J积分的变化趋势图，图中对比了公式（2）的理论解以及有 限元结果。从图中可以看出，两者是吻合的，说明了有限元模拟是正确的。J积 分随加载的变化趋势图如图5所示。红色实线表示的是理论结果，黑色点表示有 限元结果。kN分积iRiimiiiiifliStep time=0.7976Step time=0.8708Step time=0.6743Step time=0.9551Step time=0.1143Step time=0.3943TBW如1也仰minIMHI Hl IIStep time=0.9994Step time=1Step time=0.9998图6裂纹扩展过程从上述的裂纹扩展过程的应力分布图，我们可以得到如下几点结论，证明我 们的数值模拟具有一定的正确性。首先，在裂纹尖端出现了应力集中，这是和断 裂力学理论符合的。其次，观察裂纹附近的应力分布，我们可以看到应力分布的 趋势是与理论计算的塑性区的形状大致相同的，理论计算的塑性区形状如图7所 示。图7理论上的塑形区形状图3.3裂尖单元应力变化其次，考察裂尖单元的应力随载荷增加的变化。实际上裂尖单元应力值的具 体大小并没有意义，因为表征断裂韧强度的是应力强度因子和J积分。而单元应 力随载荷增加导致的变化可以帮助我们理解裂尖单元在起裂到完全断裂的整个 过程。观察图8,我们可以大致得到这个裂尖单元参与断裂过程的整个历史。首先， 在应力加载的早些时候，裂尖单元的应力随着载荷的增加而增加，此时裂尖单元 的应力并没有达到损伤判据的临界应力，所以单元没有起裂。随着载荷的继续增 加，应力值继续增加，当到大概0.7976左右时，裂尖的最大主应力达到了最大主 应力损伤判据的临界值，于是裂纹起裂，直至完全裂开，单元的应力奇异性消失， 裂尖单元转变为一个普通的非裂尖单元。这个过程对应于图8中的右边应力增大 后有急剧减小的曲线。之后，由于载荷还没有完全加载完毕，所以裂尖单元在转 变为普通单元之后随着载荷的继续增加，其单元应力又会随之在增加。以上就是一个裂尖单元在整个加载过程中的单元应力历史变化的三个阶段。0.20.40.60.81.0Step time1401201008060402000.0辅力应主大最图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个单元长度。Step time=0.7976I 1 口 I I I iBBrni图9随加载历史扩展的动态裂纹图3.5裂纹扩展的能量历史曲线绘出整个模型的总能量，动能，内能和外力功随计算过程的历史曲线如图10 所示。由图可以看出，总能量和动能在整个过程中并没有发生显著的增加，可以 表明整个计算过程基本是稳定的。而伪应变能在计算过程中有稍微的增长。图10裂纹扩展的能量历史曲线4结论本文的工作是基于abaqus6.9版本的扩展有限元功能计算了 I型裂纹的扩展。 得到如下几点结论。首先，本文计算了1型裂纹的应力强度因子，结构表明计算值与理论值的误 差在5%以内，结果可靠有效。同时计算了裂纹的J积分值，获得了积分的历史 曲线，有限元结果与理论结果吻合。裂纹的扩展过程与理论是吻合的，整个裂纹扩展的计算没有出现不稳定的情 况。裂尖出现了应力的奇异性，裂纹扩展平稳。并通过研究裂尖单元的应力历史 曲线，直观的获得了裂纹在一个单元上扩展的三个典型阶段。其次，对裂纹扩展的过程进行了分析，指出了裂纹扩展的大概几个阶段以及 具体的裂纹扩展过程。最后，对裂纹扩展过程中的能量变化作了简单描述，说明了计算过程的稳定 性，并验证了系统的能量守恒关系。