von Mises率无关弹塑性本构理论及UMAT源代码解析

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码

1 本构理论

1.1 率形式

对于各向同性线弹性材料,其本构方程为:

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图1

式中假设了应变张量可以分解为弹性应变和塑性应变两部分:

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图2

因此塑性本构的关键在于计算塑性应变的演化。对于率无关弹塑性的本构理论,需要确定以下三个部分:

(1):屈服条件

(2):流动法则

(3):硬化法则

在此采用的是 von Mises 屈服条件

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图3

式中后继屈服应力是等效塑性应变的函数:

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图4

流动法则为:

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图5

式中流动方向的表达式为:

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图6

硬化法则为:

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图7

1.2 Return-mapping算法

上述的本构方程均为率形式。在增量步中,给定增量应变:

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图8

首先假设该增量应变全为弹性应变,计算试验状态下的一些物理量:

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图9

试验状态下的应力

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图10

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图11

试验状态下的屈服函数值:

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图12

利用该试验屈服函数值来判断在该增量步下是否发生了塑性屈服。如果:

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图13

则说明试验状态即为真实状态,即可进行更新:

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图14

反之则需要进行塑性更正,即需要计算塑性乘子的增量,利用以下非线性方程组进行计算:

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图15

可以将该非线性方程组简化至一个非线性方程,过程如下,将该方程组中的第一式分解为球量和偏量两部分:

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图16

因此可以计算应力为:

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图17

将上式中的第二式整理得到:

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图18

可以得到两个张量的方向相同:

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图19

因此偏应力可以用试验状态的信息表示出来:

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图20

代入到最后一个一致性方程中可得:

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图21

即可利用牛顿迭代法对上述非线性方程进行求解,得到塑性乘子增量。求解得到塑性乘子增量之后,即可更新:

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图22

也可以更新塑性应变:

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图23



1.3 一致性切线刚度矩阵

umat除了要求更新应力以及状态变量之外,还需要更新算法的一致性切线刚度模量。当没有发生塑性屈服时,一致性切线刚度矩阵即为弹性矩阵。当发生塑性屈服时,根据应力更新的公式:

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图24

可以计算:

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图25

式中:

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图26

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图27

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图28

式中:

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图29

为硬化曲线的梯度,各向同性硬化曲线使用数据点拟合为分段线性函数:

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图30

因此最终可得到一致性切线刚度矩阵的表达式为:

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图31



2 算法描述

整个 return-mapping 算法框图如下:

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图32



牛顿迭代法算法框图如下:

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图33

3 UMAT代码

umat使用F90编写,其主程序如下:

include "plastic_iso_pack.f90"


subroutine UMAT(stress,statev,ddsdde,sse,spd,scd,                       &
                rpl,ddsddt,drplde,drpldt,                               &
                stran,dstran,time,dtime,temp,dtemp,predef,dpred,cmname, &
                ndi,nshr,ntens,nstatv,props,nprops,coords,drot,pnewdt,  &
                celent,dfgrd0,dfgrd1,noel,npt,layer,kspt,kstep,kinc)


    use plastic_iso_pack


    include 'aba_param.inc'


    character*80 cmname
    dimension    stress(ntens),statev(nstatv),ddsdde(ntens,ntens),   &
                ddsddt(ntens),drplde(ntens),                        &
                stran(ntens),dstran(ntens),time(2),                 &
                predef(1),dpred(1),props(nprops),coords(3),         &
                drot(3,3),dfgrd0(3,3),dfgrd1(3,3)


    !*******************************************************************************
    ! 材料参数
    integer :: n_pt                ! 硬化曲线的数据点个数
    real(8) :: hard_data(int((nprops-2)/2),2)   ! 硬化曲线的数据点表格
    real(8) :: E,nu
    real(8) :: mu,kappa
    


    ! 从props数组中读取材料参数
    n_pt = int((nprops-2)/2)
    E = props(1)
    nu = props(2)
    j=3
    do i = 1, n_pt
        hard_data(i,1) = props(j)
        hard_data(i,2) = props(j+1)
        j = j + 2
    enddo


    ! 更新应力,状态变量以及一致性切线刚度模量
    mu = E / (1.0 + nu) / 2.0
    kappa = E / (1.0 - 2.0*nu) / 3.0
    call plastic_iso_umat(mu,kappa, n_pt,hard_data, stran,dstran, stress,statev,ddsdde)


    return
end subroutine UMAT

umat所需的支持函数编写于单独的module中(plastic_iso_pack.f90)用于调用,该方式可以使得主程序十分简洁。

4 测试

4.1 一个单元单轴加卸载试验

对一个单元进行单轴的加卸载试验,绘制应力应变滞回曲线。Abaqus设置的材料参数为,弹性部分的性质:

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图34

塑性部分的性质:

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图35

umat设置相同的材料参数为:

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图36

以及需要设置状态变量的个数为8。二者计算的结果对比如下:

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图37

二者等效塑性应变的演化对比图如下:

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图38

二者塑性耗散的演化对比图为:

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图39



4.2 带孔板拉伸实验

利用相同的材料性质计算一带孔板的拉伸实验。Abaqus计算的结果如下:

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图40

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图41

使用umat计算的结果如下:

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图42

各向同性硬化von Mises率无关弹塑性本构理论以及umat源代码的图43

可看到二者的结果是完全相同的。


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

QR Code
微信扫一扫,欢迎咨询~

联系我们
武汉格发信息技术有限公司
湖北省武汉市经开区科技园西路6号103孵化器
电话:155-2731-8020 座机:027-59821821
邮件:tanzw@gofarlic.com
Copyright © 2023 Gofarsoft Co.,Ltd. 保留所有权利
遇到许可问题?该如何解决!?
评估许可证实际采购量? 
不清楚软件许可证使用数据? 
收到软件厂商律师函!?  
想要少购买点许可证,节省费用? 
收到软件厂商侵权通告!?  
有正版license,但许可证不够用,需要新购? 
联系方式 155-2731-8020
预留信息,一起解决您的问题
* 姓名:
* 手机:

* 公司名称:

姓名不为空

手机不正确

公司不为空