俞茂宏的双剪强度理论及后来的统一强度理论是强度理论的大概括,在岩土工程中最常用的Mohr-Coulomb模型等都是统一强度理论的某一特例,因此,该理论得到了大发展。然而,ABAQUS并没有自带统一强度理论的材料本构,因此,要想使用需要自己进行二次开发。目前,网络上关于该理论的二次开发文献数量屈指可数,有的大多也就是介绍了理论部分,在代码编写方面几乎为0。花了十几天时间详细阅读了文献后,本人成功开发了基于ABAQUS的统一强度理论的UMAT本构,经大量模型验证程序可行,特记录一下,也为后来人作部分参考。
1 理论部分
理论部分是最关键的,这部分网上资料众多,这里强烈推荐该理论之父俞茂宏先生所著《Computational Plasticity》,该书专门介绍统一强度理论的,为英文版。下面简单介绍下该理论。
F为统一强度理论的屈服函数表达式,若F>0,则材料发生塑性变形;若F<0,则仍处于弹性阶段。当a与b取不同值时,统一强度理论将退化为下面各种理论。
这也是“统一”的含义。
塑性势函数的一阶导形式与屈服函数一阶导相同,因为后面要使用完全隐式向后欧拉方法,所以需要势函数的二阶导。
同时,对于奇异点,需要进行消除。
当b=1时
当b≠1时
关于理论部分仅做展示,文献众多,看官可自行查阅。
2 程序编写
这里最关键的便是UMAT的编写了,毕竟理论公式都会推但想利用ABAQUS计算还得按照其格式进行编写子程序才行。
UMAT的格式,FORTRAN语言的语法等这里不做介绍,可参见我之前发过的一些小文。下面说下将该理论实现代码的部分过程。
首先是计算流程
当未屈服时,材料处于弹性阶段,此时代码极其简单,只需填写弹性的刚度矩阵作为雅克比矩阵,然后更新应力和状态变量即可。
当F>0时,材料发生屈服,进入塑性阶段。塑性阶段需要判断屈服条件,计算屈服函数和势函数的导函数,进而采用向后欧拉方法将应力拉回屈服面。
计算屈服函数的一阶导作为流动矢量。
代码量很大不一一展示。
3 模型展示
当主应力系数b=0时,统一强度理论退化为mohr-coulomb模型,abaqus自带这种模型。
将长方体模型一端固定,一端加相同压力,观察结果对比
UMAT计算应力结果
自带mohr-coulomb本构计算应力结果
由结果可见,计算精度十分高,验证了模型的正确性。