ABAQUS UEL - Embedded crack model 在CST单元中的简单实现

工程实际中经常遇见材料开裂问题。对于受到外力作用的结构，材料的承载力随着结构的变形发生改变，比如延性材料。因此可以采用连续体塑性力学框架对结构进行分析模拟。但是对于脆性材料，连续体塑性力学尽管可以模拟出结构宏观的力-位移关系，但是很难模拟出材料的开裂破坏。

为了采用有限元模拟材料的开裂，cohesive单元常被预设在连续体单元之间，通过cohesive塑性或者损伤本构模拟出材料强度下降的过程。但是，这种方式使得裂缝只能在预设cohesive单元的位置开裂。为了更好的模拟连续体单元的开裂，Embedded crack model假设在单元内部存在开裂面，并通过cohesive zone model定义开裂面的本构关系。开裂面的角度通过最大主应力确定。这种方式极大地方便了网格划分过程。 本文通过3节点CST单元，介绍如何使用UEL子程序在ABAQUS中实现简单的Embedded crack model。

1. Embedded crack model 的概念

假设有一个三角形单元如下，

在节点1的水平位置施加位移0.05，单元的变形为，

对于延展性良好的材料，上述变形发生时材料发生塑性变形；但是对于脆性材料，单元发生开裂，开裂的效果为下图，

开裂面上的坐标系存在法向（n）和切向（s）行为。材料发生开裂时，开裂面的本构是非线性的，而其他连续体部分依然是线弹性的。

相比起将整个单元都定义为非线性，这种将非线性开裂面内嵌在线弹性连续体内的方法更接近脆性材料的性质。

2. 开裂面本构关系

连续体线弹性的本构采用Plane strain relation，在大部分有限元参考书中都有介绍，其刚度矩阵为，

开裂面的本构建立在塑性力学框架中，采用了简单的 tension cut-off 屈服准则。

tn 表示开裂面上法向应力，ft 则是材料的受拉强度，可以表达为，

3. 数值算法公式推导

当开裂发生时，其应变计算为，

和弹性单元不同，这里的应变不再是简单将应变-位移矩阵B乘以节点位移向量U，而是减去了由于开裂导致的多余位移He。另外需要提到的是，H矩阵的表达式是根据节点开裂角度定义的，如下图，

开裂面和连续体之间的应力转换为，

c 和 s 分别是

开裂角度和最大应力的计算为，

通过本构关系计算应力为，

所以，接触面上的应力可以表示为，

把上式写为rate形式，

再联系consistency 条件，

以及associate flow rule

计算得到了塑性因子增量，

把塑性因子增量带回到前面的式子，可以得到切向刚度矩阵，

将以上公式进行整理，得到算法如下，

FE 模拟结果

单元在受拉时，由于材料开裂导致了强度下降。但是在ABAQUS后显示时，单元表现出的是变形而不是开裂的状态，这是CST单元缺陷。

UEL 和 input文件在收费内容中。