ABAQUS UMAT二次开发初探

初识ABAQUS UMAT二次开发

在这一期本来小编想介绍如何使用UMAT实现弹塑性本构关系，可是后来再百度中发现这些东西已经烂大街了。小编索性就自己臆想了一个本构关系，即应力-应变关系曲线符合二次抛物线函数。这样既有了上升段，也有下降段。通俗的讲就是可以实现随着荷载递增结构内力先增大，其后降低直至失去承载能力。

下面小编就带大家一起看看这种本构关系是如何在ABAQUS中用UMAT实现的。

1模型简介

为了让大家更好的理解UMAT的使用方法。这一期的模型跟上一期一样采用最简单的桁架结构，如图1所示。模型由三根杆件组成，一根竖杆、两根斜杆。竖杆长为30，与斜杆的夹角均为45度，截面面积为1。下面三个节点均施加x、y向约束，上部节点施加一个向上的节点位移荷载0.0017。单元仍旧采用Truss单元，每根杆件划分为1个单元。

图1  有限元模型及边界条件示意图

材料属性设置：在User Material中定义两个参数，一个值是-7244290000，一个值是412000，这两个数分别是应变-应变二次抛物线函数的两个系数a和b。通过变量props传递给UAMT子程序。应力-应变关系曲线如图2所示：

图2  本构关系曲线（应力-应变关系曲线）

分析步和求解器设置：小编在试算的时候发现，如果将荷载定义为节点力，那么使用abaqus/standard Static General程序求解收敛性不稳定，使用Dynamic/implicit的拟静态程序求解收敛性会好很多。本例施加的是节点位移荷载，我们选用Static General即可。在增量设置里选用Fixed类型，最大增量步设置为100步，增量步长设置为0.01。求解器类型选用迭代法，其它保持默认即可。

有关abaqus/standard求解非线性问题的相关介绍或知识，大家可以查看abaqus的帮助文档或者庄茁老师的《基于abaqus的有限元分析和应用》的第八章，介绍的非常详细。

2UMAT子程序设置

本例用到的UMAT相关知识与上一期基本相同。唯一不同的是多使用了一个变量STRAN。STRAN（NTENS）：应变矩阵（或向量），包括NDI个直接分量，NSHR个剪切分量，与应力变量STRESS相对应。整个UAMT子程序如下：

 1    SUBROUTINE UMAT(STRESS,STATEV,DDSDDE,SSE,SPD,SCD, 2    1 RPL,DDSDDT,DRPLDE,DRPLDT, 3    2 STRAN,DSTRAN,TIME,DTIME,TEMP,DTEMP,PREDEF,DPRED,CMNAME, 4    3 NDI,NSHR,NTENS,NSTATV,PROPS,NPROPS,COORDS,DROT,PNEWDT, 5    4 CELENT,DFGRD0,DFGRD1,NOEL,NPT,LAYER,KSPT,JSTEP,KINC) 6C 7    INCLUDE 'ABA_PARAM.INC' 8C 9    CHARACTER*80 CMNAME10    DIMENSION STRESS(NTENS),STATEV(NSTATV),11    1 DDSDDE(NTENS,NTENS),DDSDDT(NTENS),DRPLDE(NTENS),12    2 STRAN(NTENS),DSTRAN(NTENS),TIME(2),PREDEF(1),DPRED(1),13    3 PROPS(NPROPS),COORDS(3),DROT(3,3),DFGRD0(3,3),DFGRD1(3,3),14    4 JSTEP(4)15C16    REAL*8  emod17c ----------------------------------------------------------------18c   umat for nonlinear problems19c   be used for Truss element20c ----------------------------------------------------------------21c   props(1) - a22c   props(2) - b23c ----------------------------------------------------------------24c25    if (ndi.ne.1) then26        write (7,*) 'This umat may only be used for elements with one direct stress components'27        call xit28    endif29c30c   elastic properties31    emod = 2*props(1)*STRAN(1)+props(2)32c33c   elastic stiffness34c35    ddsdde(1, 1)=emod36c37c    calculate stress38c39    stress(1) = stress(1)+ ddsdde(1,1)* dstran(1)40c41    RETURN42    END

代码1~15行是UMAT的固定代码。

第16行定义了一个实数变量emod用于存储剪切刚度。

第25~28行，用于单元类型判断。当单元节点自由度数不等于1时，报错。

第31行计算剪切刚度，即应力-应变关系函数的一阶导函数。不清楚的可以看一下abaqus帮助文档中关于非线性求解器的介绍，或者百度一下牛顿迭代法的相关内容。

第35~42行，与上一期相同，都是用于根据应变增量更新应力。

3计算结果

荷载作用下结构各单元的内力变化示意图：

图3  结构的内力变化示意图

   由图3可知，竖杆的内力变化明显经历了两个过程，从随荷载的增大而增大到随荷载的增大而降低。竖杆的应力随荷载增量的变化曲线，如下图所示。

图4  竖杆的应力随荷载增量的变化曲线

可以从abaqus计算结果中导出竖杆的应力-应变关系曲线，如下图所示。

图5  从abaqus计算结果中导出的竖杆应力-应变曲线

对比图5和图2可知，模型中竖杆的应变-应变曲线与输入的本构关系完全一致，可见UMAT程序正确。

4结语

Abaqus UMAT用户子程序接口功能非常强大，编写一个漂亮的UMAT子程序需要一定的力学基础（尤其是有限元、连续介质力学）。本文旨在让大家对UMAT有一个初步的了解，以及使用UMAT需要用到哪些基本知识，好在今后的学习中有侧重点。后续的推文中，小编还会进一步结合我们有限元理论模块的进度，进一步讲解Abaqus UMAT的使用。

免责声明：本文系网络转载或改编，未找到原创作者，版权归原作者所有。如涉及版权，请联系删