并不简单的弹塑性本构子程序

初学材料力学就知道最常见的金属一般都是弹塑性的。

所谓弹塑性，就是把材料性能划分成了两个阶段，前面的阶段是弹性，比较好理解，载荷与变形线性变化。后面塑性，就是指材料继续变形，但是载荷不往上走了，或者即便走也变慢了。而且即便完全卸载，第二个阶段的变形仍然会保留。

材料如此，人亦如此，过度消耗是补不回来的。弹塑性材料有屈服强度这个概念，就是指进入塑性后，本来向上的曲线开始低头了，所谓之“屈服”。

只要做结构强度方向，弹塑性几乎是个天天都能听到的词，以至于我对它毫无“敬畏之心”，总觉得这个玩意很简单。尤其是我研究生开始做复合材料力学以后，就觉得复合材料比金属高端多了。我们材料是各向异性的，刚度矩阵更复杂，我们还有蔡吴、蔡希尔、哈辛一堆“高级”失效准则，材料还可以分层失效，写到论文里面更好看，更别提失效因子、渐进损伤，总之就是牛掰。

后来我第一次写弹塑性本构的时候，懵了。这玩意比我想象的要复杂的多。首先弹塑性这个问题并不简单，要想解释清楚它，需要从材料微观层面，了解晶体位错等等现象。甚至于到2011年，写这些问题的综述还能发一篇Nature。归根结底，我们并未完全研究透材料的弹塑性行为，以及相关的强度、韧性问题。

即便就本构层面而言，弹塑性光一个塑性流动方向要想写出来就不容易，网上能看到一大堆公式，各种导数偏导数。

问题是在UAMT/VUMAT里面是很难做这种偏导的，包括迭代数值计算，不是完全不能，而是写出来大概率各种报错，还不好调试找原因。在子程序里面，最稳妥的就是写加减乘除。

那时候写弹塑性本构，对我理解子程序以及ABAQUS逻辑，起到了非常重要的作用。我的体会是，学写子程序，应该先写弹性，接着就写弹塑性，这样才能打好基础。像我当时属于是回头补课。

在写弹塑性本构之前，我对塑性流动是干嘛使的没有直观概念。写的时候我才明白，由于只能先算出来等效塑性应变，没有流动方向的话，就无法把它转换到各个应变分量，不知道应变分量就无法计算应力。这玩意从数学上讲，是一个转换公式。

我们目前重工业上大部分的结构材料还是金属，尽管ABAQUS中有自带的JC模型，但是如果要模拟更复杂的情况，学会写弹塑性本构就十分必要。

本期就给一个弹塑性VUMAT拉伸失效的案例，结合单元删除技术，模拟结构破坏过程。

本构模型

采用经典老演员JC模型描述本案例的弹塑性本构：

为了模拟结构破坏，采用如下准则判断单元完全失效，满足其一即可：

（1）材料Mises应力达到极限值；

（2）材料极限应变达到极限值。

子程序结构

子程序的基本结构如下：

1.初始化准备工作 

程序首先进行初始化准备工作，读入材料的弹性参数、强度参数、硬化参数以及应变率相关参数，然后构建弹性刚度矩阵，为后续计算奠定基础。

2.进入材料点循环

接下来进入材料点循环，对每个积分点逐一进行计算。对于每个材料点，程序首先读取上一步的状态变量，包括累积的等效塑性应变、应力状态以及背应力等内部变量。

3.失效判断

程序随后进行失效判断，检查材料是否满足失效准则。判断依据包括两个方面：一是等效塑性应变是否超过极限应变阈值，二是等效应力是否达到破坏强度。一旦满足任一失效条件，程序将材料标记为失效状态，并大幅降低其刚度以模拟材料的承载能力丧失。

4.本构响应计算阶段

在本构响应计算阶段，程序考虑了应变率效应和材料硬化特性，更新当前的屈服应力。同时计算应力偏量，得到米塞斯等效应力和塑性流动方向，这些是判断材料是否屈服的关键参数。

5.弹塑性判别

然后进行弹塑性判别。将当前等效应力与更新后的屈服应力进行比较：

• 若未达到屈服，材料表现为弹性响应，应变增量全部转化为弹性应变，应力通过弹性刚度矩阵直接计算得到。
• 若超过屈服，材料进入塑性状态，此时需要计算塑性应变增量。程序根据屈服准则和流动法则，确定塑性应变的大小和方向，将总应变增量分解为弹性部分和塑性部分，然后用弹性应变计算更新后的应力。

6.状态变量更新环节

完成应力应变更新后，程序进入状态变量更新环节，更新等效塑性应变、塑性应变率、背应力等内部变量，这些变量将传递到下一时间步继续使用。程序还计算并存储材料的内能，用于能量平衡检查。

7.循环处理

处理完当前材料点后，程序继续循环处理下一个材料点，直至所有积分点计算完毕，最后返回主程序。 

测试模型

在ABAQUS的PART模块，用回转体方法，创建一个哑铃形试验件，根据对称性，建立1/4模型。

底部固定，顶部拉伸，内部面设置对称条件。

为了加速计算进程，设置质量缩放。

拉伸破坏过程

随着载荷增加，试验件中间位置出现出现“颈缩”，变形达到一定程度，试验件断开。