在之前的文章里曾经多次提到过Hashin准则，这是目前区分失效模式的判据中应用最广泛的判据之一，已被Abaqus、Ansys、MSC等大型商业软件所集成。无论中文还是外文有关采用Hashin准则进行复合材料渐进失效分析的文章也是铺天盖地、数不胜数，Hashin于1980年发表的一篇单向纤维增强复合材料失效准则的文章被引用了3790次。在提出该理论时，本来是用于预测单向复合材料失效行为的，然鹅，目前大家基本都在将其应用于层压板的失效预测。

关于Hashin准则的描述以及在WWFE中的表现，之前已经撰文描述过，此处不再赘述，感兴趣的可以点击下方链接了解详情。

聊一聊世界复合材料失效运动会(WWFE)——搞复材失效而不知WWFE你就out了

复合材料失效理论知多少？(一)

本文主要讲解一下Abaqus中使用Hashin失效判据以及基于能量的演化判据进行渐进失效分析时各种参数和变量的定义和来由。有一些读者对这两者的组合使用的非常熟练，但并不了解损伤演化过程中失效判据和临界应变能释放率是如何控制损伤扩展的，希望通过本文能帮助读者对复合材料渐进失效分析有进一步的认识。

01  二维应力应变关系

在Abaqus中Hashin失效准则适用于常规壳单元或者连续壳单元，对应的材料模型为Lamina，应力应变关系如下：

其中，各参数满足稳定性约束条件：

写成刚度矩阵的表达式如下：

简写为：

其中，C0为无损状态下的弹性矩阵。

欲表征材料的应力应变关系，在CAE中仅需要输入弹性常数即可。

图1

02 二维Hashin准则表达式

在Abaqus中引入的二维的Hashin失效准则，共分为纵向拉伸(纤维方向拉伸)、纵向压缩(纤维方向压缩)、横向拉伸(垂直纤维方向拉伸)、横向压缩(垂直纤维方向压缩)等四种失效模式。

上述四个公式用于判断的是材料点的损伤起始，当Fft达到1的时候，表示该材料点纤维方向拉伸损伤刚刚开始，小于1的时候，材料纤维方向无拉伸损伤。其他同理。绝大部分论文仅引用了这四个公式就去做渐进失效分析了，是不完整的，这仅仅能预测损伤的起始，还缺乏损伤起始以后的刚度退化。

与损伤起始相关的是6个强度值，定义如下，6个参数分别是：纵向拉伸强度XT，纵向压缩强度XC，横向拉伸强度YT，横向压缩强度YC，纵向剪切强度SL，横向剪切强度ST。

图2

03 基于能量的损伤演化

定义损伤演化参数的界面如下，输入四种失效模式对应的断裂能:纵向拉伸断裂能Gftc,纵向压缩断裂能Gfcc,横向拉伸断裂能Gmtc,横向压缩断裂能Gmcc。(断裂能，即临界应变能释放率，是裂纹扩展单位面积所需要的能量，一般单位取为N/m，J/m²，N/mm或mJ/mm²。

有关纤维方向断裂韧性测试的文章可点击以下链接：

复合材料纤维方向断裂韧性的测量方法

图3

单纯使用ABAQUS内嵌的渐进失效模型进行失效分析时，了解这些输入参数的物理意义就足够了。

如果进一步深究或者自己去写渐进损伤分析的程序，还会存在一些问题，比如上述断裂能是如何起作用的呢？尤其是四种失效模式均为多应力状态下的复合判据，损伤状态变量又是如何求解的呢？损伤状态变量又是如何控制刚度的衰减的呢？

04 损伤状态变量的求解

损伤起始以后，材料的刚度将会发生逐渐退化，此时开始进入损伤演化阶段，材料的刚度退化程度用损伤状态变量d来表示。

含损伤的材料应力应变关系为：

其中：

上述公式中，

df,dm,ds的表达式又是这样子的：

即纤维方向的损伤取决于纤维方向拉伸损伤状态dft和纤维方向压缩损伤状态dfc，二者只能发生其中的一种，具体发生哪种由沿纤维方向的应力σ11决定，当σ11>=0时，潜在发生的是纤维方向的拉伸损伤，取df=dft，反之，发生纤维方向的压缩损伤，取df=dfc。df的值介于0~1之间，0代表材料完好，1代表彻底失效。

同理，可以得到dm的值。

另外，无论纤维和基体发生拉伸还是压缩损伤，材料均失去承受剪切载荷的能力，因此剪切损伤变量ds并非独立的损伤状态变量，ds取决于独立损伤状态变量dft, dfc, dmt, dmc，只要四个损伤状态中有一个不为0，ds则不为0，四个状态变量有一个达到1，ds必然为1。

那独立的损伤状态dft,dfc,dmt,dmc又是如何求解的呢？

05 独立损伤状态变量的求解

帮助文档以及其他能查到的资料中都给出了下列损伤状态变量d的统一表达式和原理图，却并没有给出详细的解释。

图4

d即为独立的损伤状态变量，要想理解损伤状态计算公式，首先需要引入两个变量定义，

那么，

注意上标0以及红框中标记的限定词“刚好达到1时”，这是理论上的数值，即某一种判据恰好判断因子等于1的一瞬间所记录下的等效位移，记为

该变量代表的是某个材料点某种失效模式下的损伤起始位移，对应的是上图4三角形顶点的横坐标。该损伤起始位移并非一个恒定不变的材料常数，而是随着材料点的应力状态而变化的，不同的单元有可能得到不同的损伤起始位移值。这一点也是很多人不解的地方。在自己写子程序的时候，正确理解这一点很重要。

类似地，当某种失效判据判断因子恰好达到1时，记录下该时刻的等效应力，该等效应力即为当前材料点在当前应力状态下，某种失效模式所对应的等效强度，也就是上图中三角形顶点对应的纵坐标。同样，等效强度也不是恒定不变的材料常数。

具体到不同失效模式下，等效位移和等效应力的表达式也是有所不同的，在Abaqus中，纤维拉伸、纤维压缩、基体拉伸、基体压缩状态下四种等效位移和等效应力的表达式如下：

其中，Lc为单元的特征长度，一阶面单元中特征长度是面积的平方根，一阶体单元中，特征长度是体积的立方根。

根据上述公式，再加上四种失效判据的公式，随着载荷的不断增加，就可以得到每个材料点每一种失效模式下的损伤起始位移和等效强度，损伤起始位移表达式如下：

等效强度表达式如下：

至此，每种失效模式下，损伤演化部分的损伤起始位移和等效强度就求解出来了，也就是本构三角形顶点对应的横坐标和纵坐标已知了。未知的是三角形底边的长，对应的也就是每个材料点上彻底失效时的位移。

图5

每种失效模式彻底失效时对应的失效位移计算公式很简单：

其中Gc为临界应变能释放率，也就是在CAE中输入的损伤演化部分的断裂能(前面提到的Gftc,Gfcc,Gmtc,Gmcc)，也是三角形的面积。底边的长度即为两倍三角形的面积除以三角形的高，高就是等效强度，已在上一步求解得出。此时三角形中的各参数都已经可以求解出，就可以带入独立损伤状态变量d的计算公式了。

图6

以纤维拉伸失效为例，

同样的方式，四种独立的损伤状态变量就都可以求解出来了。

四种独立的损伤状态变量对应的断裂能在图3中输入。

文中多次提到失效判据判断因子恰好等于的一刻，这仅仅是理论解，在实际的程序计算中，几乎不可能使判断因子恰好等于1，有可能大于1，这种情况在自己写子程序时也是必须要考虑的，以后再专门写文章讲解这个问题吧。

本文完结，文中后半部分某些公式为自己编写，无参考文献，表述不当之处，敬请批评指正。觉得有用，欢迎转发！

Abaqus基础课程10讲

Abaqus/ATOM非参优化(Tosca)

Abaqus非线性和工程实例进阶17讲Isight-Abaqus参数优化和材料拟合课程如何提高虚拟仿真的影响力？内心坚韧，眼中有光