二维中心裂纹板的应力强度因子分析
1. ABAQUS裂纹体的建立
在ABAQUS中，根据裂纹尖端的形状可以建立两种类型的裂纹体：尖锐（Sharp）型裂纹和钝形（Blunt）裂纹，如图1所示。对于尖锐型裂纹，裂纹尖端存在奇异性；而钝形裂纹尖端可以看做是一个具有给定缺口半径的缺口，因此裂纹尖端不存在奇异性，可以按照常规有限元的建模方式来建立。在ABAQUS中，两种类型的裂纹均可以进行应力强度因子分析。
图1 尖锐型裂纹和钝形裂纹
对于两种类型的裂纹，由于裂纹尖端均存在应力集中，在裂纹尖端将产生非常高的应力梯度，因此划分裂纹体时通常需要对裂纹尖端的网格进行细化以获得精确的应力值。需要注意的是，对于尖锐型裂纹，由于裂纹尖端存在奇异性，因此细化裂纹尖端的网格并不能使得裂纹尖端的应力和应变值收敛，减小网格尺寸只会使得裂纹尖端的应力值增大；而对于钝形裂纹，由于裂纹尖端存在钝形缺口，因此细化裂纹尖端的网格将最终得到收敛的缺口应力值。
在采用两种类型的裂纹体计算J积分时，由于J积分的实质是能量释放率，因此在进行线弹性有限元分析时，采用非常粗糙的网格也能够获得精确的J积分值，即使在这种情况下裂纹尖端的局部应力应变场并不是十分精确。而在进行弹塑性断裂力学分析时，通常需要细化裂纹尖端区域来获得精确的J积分值。
本文仅考虑尖锐型裂纹的建立，在ABAQUS中尖锐型裂纹通常也被称为seam裂纹，seam裂纹可以看做是部分或完全插入到模型中的代表裂纹的一条边线，如图2所示。seam裂纹通常由单元边线构成，在这些单元边线上的节点会自动复 制节点，位于这些单元边线上的单元不会共享节点，从而实现了裂纹面的分离。
图2 seam裂纹
如图2中所示，通常裂纹尖端会采用一圈三角形单元进行划分，在三角形单元的外围通常还会围绕多层四边形单元进行过渡。需要注意的是，在ABAQUS中当采用围线积分（contour integral）来计算应力强度因子时，在围线积分的区域只能为四边形或六面体单元，虽然裂纹尖端的网格为三角形单元，但其实际上是退化的四边形单元（degenerated quad），后面将会介绍如何划分裂尖网格。
为了在ABAQUS中定义裂纹，首先需要指定裂纹面（crack front）以及裂纹尖端（crack tip），对于通过ABAQUS创建的二维部件实例，裂纹前沿可以指定为几何点，几何边线以及几何面，裂纹尖端可以指定为几何点，而对于外部导入的二维网格（orphan mesh），裂纹前沿可以指定为节点，单元边和单元面，裂纹尖端尖端可以指定为节点。除此之外还需要指定裂纹面的法向矢量方向或者裂纹扩展的方向，在ABAQUS中裂纹扩展的方向也被称为q向量，该向量将用于围线积分的计算。裂纹的定义如图3所示。
图3 ABAQUS裂纹定义界面
为了在裂尖单元中引入奇异性，需要对单元节点进行特殊的处理。如图4所示，对于8节点的四边形单元（二阶单元，具有中间节点），首先ABAQUS会将四边形单元的其中一条边压缩，假设该单元边由节点a, b和c构成，压缩之后节点a, b和c将合并共同构成裂纹尖端，随后与裂纹尖端相连的两条单元边上的中间节点将会被移动到距离裂纹尖端1/4处的位置。
图4 二阶四边形等参单元奇异性的构建
如果节点a, b和c能够自由移动，并且中间节点移动到了1/4处，如图4中的情况(1)所示，则在退化四边形单元任意位置处的应变具有如下奇异性：
如果节点a, b和c被约束在了一起，并且中间节点移动到了1/4处，如图4中的情况(2)所示，则单元具有如下奇异性：
而由线弹性断裂力学可知，裂纹尖端的应力应变场同样具有r-1/2的奇异性，例如对于I型裂纹，其裂纹尖端的应力场可表示为：
如果节点a, b和c能够自由移动，但中间节点没有移动，则单元具有如下奇异性：
此时应变的奇异性适合于理想塑性的情况。
在ABAQUS，单元奇异性的设置如图5所示。
图5 单元奇异性设置界面
对于中间节点的调整，可以通过改变Midside node parameter的取值来实现，Collapsed element side, single node表示四边形的单元边被压缩，并且裂纹尖端具有相同的节点号；Collapsed element side, duplicate nodes表示虽然四边形的单元边被压缩，但裂纹尖端具有不同的节点号，两种情况如图6所示。
图6 裂纹尖端退化单元控制
从图6中可以看出，虽然两种单元都是退化的四边形单元，均具有4个节点，但裂尖的节点有些不同。Single node对应了图6中右边的情况，此时裂尖实际上只有一个节点，Duplicate nodes对应了图6中左边的情况，裂尖具有两个节点，只是节点坐标是重合的。
如果在退化单元控制中选择No degeneracy，则四边形的单元边将不会压缩，此时裂纹尖端的网格不同于图2中网格，将完全由非退化的四边形单元构成，如图7所示。
图7 裂纹尖端采用非退化四边形网格
如图7所示，如果单元靠近裂纹尖端的中间节点被移动到了1/4处，则在这些单元边上应变具有r-1/2的奇异性，但在单元内部不具有奇异性，这种情况要优于没有奇异性的情况，但没有单元退化的情况好。
需要注意的是，如果模型采用的是外部导入的网格（orphan mesh），此时无法通过图5所示的界面来控制单元奇异性，需要通过Mesh模块的Mesh Edit工具中的Adjust midside来手动调整裂纹尖端附近的中间节点。
如果在裂纹尖端使用了4节点的退化四边形单元（一阶单元），并且裂纹尖端的节点被约束到一起，如图8所示，则单元具有如下奇异性：
图8 一阶退化四边形单元
此时单元的奇异性为r-1，没有达到裂纹尖端应力场的奇异性。
为了在裂纹尖端划分出退化的四边形单元，需要在模型中裂纹尖端区域分割出一块环形区域，并将该区域指定为以四面体为主（quad-dominated）的映射（sweep）网格。
2. 中心裂纹板应力强度因子分析
下面采用ABAQUS对一有限宽中心裂纹板进行应力强度因子分析，中心裂纹板的尺寸如图9所示，板的半宽W=100mm，半高H=200mm，板厚B为1mm，半裂纹a的长度为10mm。该板在远端承受有均匀拉伸应力σ=30MPa。由于模型存在对称性，因此取1/2模型进行分析。
图9 中心裂纹板尺寸
采用ABAQUS建立中心裂纹板的几何模型，由于裂纹尖端需要划分退化的四边形单元，因此需要在几何模型中通过面分割工具在裂纹尖端切割处一块圆形区域，用于划分退化的四边形单元，而圆形区域外围还需要划分一块环形区域，用于划分用于过渡的四边形单元。需要特别注意的是，由于裂纹尖端的网格本质上是退化的四边形单元，因此不能划分为三角形单元，因此在Mesh Control中需要将该区域设置为四边形单元为主（Quad-dominated）的扫略（Sweep）划分方式，不能指定为三角形单元（Tri）。
图10 裂纹尖端的网格控制
对于裂纹尖端区域外围的环形区域，主要起到单元过渡的作用，因此划分方式选择结构化（Structured）的四边形（Quad）单元。对于其他区域，可采用四边形单元或以四边形单元为主的自由划分方式，最终得到裂纹区域的网格控制如图11所示。
图11 裂纹区域网格控制
在ABAQUS的Mesh模块中，不同颜色的区域代表不同的网格划分方式，黄色代表扫略划分方式，绿色代表结构化的映射划分方式，粉色代表自由划分方式。为了保证裂纹区域的网格正常划分，还需要通过Seed Edges工具调整裂纹区域的网格数量，如图12所示。
图12 裂纹区域网格布种
最终可以得到裂纹区域的网格如图13所示。
图13 裂纹区域网格
再次注意，虽然通过前面的划分方法最终在裂纹尖端区域得到了三角形的网格，但其实际上是退化的四边形单元，三角形单元是不能用于围线积分计算的。这一点从单元类型的选择也可以看出，本分析中假设中心裂纹板处于平面应力状态，采用带有缩减积分的8节点四边形单元CPS8R进行计算，虽然此时ABAQUS自动检测处模型中存在三角形单元，并推荐使用6节点的平面应力单元CPS6M，但从后面的分析可知，该单元类型实际上并没有被用于裂纹尖端区域。
在完成裂纹体的网格划分后下面开始创建裂纹，在Interaction模块中选择Special→Crack→ Assign Seam创建Seam裂纹，如前面所提到的，由于裂纹面上的网格是连续的，因此将裂纹面指定为Seam裂纹主要起到了网格分离的作用，如图14所示。
图14 创建Seam裂纹
图14中的红色边线即为中心裂纹的裂纹面。需要注意的是，Seam裂纹只能指定为几何边线或几何面，在外部导入的网格（Orphan mesh）中不能创建Seam裂纹，因此必须保证裂纹面上的节点是分离的。
通过Special→ Crack→ Create选择Contourintegral（围线积分），将裂纹前沿（Crack front）指定为裂纹尖端的几何点，如图15所示，并使用q向量指定裂纹扩展方向，本例中裂纹沿X轴方向扩展，因此可以将q向量设置为沿X轴的单位向量。
图15 指定Crack front
本文中使用的8节点的二阶四边形单元，为了在裂纹尖端产生r-1/2的奇异性，根据前面的讨论，需要将四边形单元的一条边压缩，将中间节点移动到1/4处，并且需要保证裂纹尖端的三个节点被完全约束到一起，如图4中的情况（2）所示，因此在Singularity中将Midside node parameter设置为0.25，并选择Collapse element side, single node。
为了在结果文件中输出应力强度因子，需要在Step模块中创建时间历程输出（History output），并指定围线积分的数量和输出数据的类型，如图16所示。
图16 创建应力强度因子历程数据输出
其他关于边界条件和载荷的设置与常规的静力分析流程完全相同，因此按静力分析的流程进行前处理即可。最终计算得到裂纹尖端的应力分布如图17所示。
提取裂纹尖端的应力强度因子如图18所示。
图18 中心裂纹应力强度因子
从图18中可以看出，由于围线积分1是在裂纹尖端处评估的，因此存在较大的误差，在剔除掉围线积分1之后，剩余围线积分计算得到的平均值为244.713MPa‧mm1/2，与中心裂纹板的解析值243.6MPa‧mm1/2相比误差仅为0.46%，因此可以认为通过ABAQUS的围线积分计算得到的应力强度因子是非常精确的。
为了验证裂纹尖端的单元为压缩的四边形单元，下面将前面分析的模型文件输出为inp文件，并导入到ABAQUS中，此时模型中的部件将仅由孤立网格（orphan mesh）构成，通过Mesh模块的Query工具查询裂纹尖端的网格，可以得到裂纹尖端的单元类型为8节点带缩减积分的平面应力单元CPS8R，因此可以确定裂纹尖端为压缩的四边形单元，而非三角形单元，其节点构成如图19所示。
从图19可以看出，裂纹尖端为具有8节点的退化的四边形单元，而非三角形单元，并且中间节点被移动到了1/4位置处，位于裂纹尖端的三个节点被约束在了一起，这样的节点布置方式能够保证整个单元都具有r-1/2的奇异性。
