蠕变中文解释
ANSYS 提供了两个用户徐变方程:USERCR.F 和USERCREEP.F 。其中:
显式徐变用USERCR.F ;前提是C6 = 100
隐式徐变用USERCREEP.F ,前提是TBOPT =100
(1)用户徐变子程序usercr ,用于显式徐变
subroutine usercr (elem,intpt,mat,ncomp,kfirst,kfsteq,e,posn,d,
x proptb,timval,timinc,tem,dtem,toffst,fluen,dfluen,epel,epcrp,
x statev,usvr,delcr)
c
c *** 基本功能: 允许用户写自己的徐变规律。该逻辑仅在C6=100时可用。
c *** 次要功能: 演示用户徐变方程的编写
c *** 注意-本文件包含ANSYS 机要信息***
c *** ansys(r) copyright(c) 2000
c *** ansys, inc.
c
c 输入变量:
c | (译者注) c |
c | 类型: int -整型
c | dp -双精度型 c | 长度: sc -标量
c | ar( , )-数组 c | 目的: in -输入
c | out -输出
c | inout -输入输出
c 变量 (类型,长度,目的) - 描述
c elem (int,sc,in) - 单元号 (标识)
c intpt (int,sc,in) - 单元积分点数
c mat (int,sc,in) - 材料引用号
c ncomp (int,sc,in) - 应力/应变分量数 (1,4 or 6)
c 1 - x
c 4 - x,y,z,xy
c 6 - x,y,z,xy,yz,xz
c kfirst (int,sc,in) - 若是首次则值为1,否则为0
c (对把状态变量初始化为非零值有用) c kfsteq (int,sc,in) - 若是子步中首次平衡迭代则为1,否则为0 c e (dp,sc,in) - 杨氏弹性模量
c posn (dp,sc,in) - 泊松比
c d (dp,ar(ncomp,ncomp),in) - 弹性应力-应变矩阵
c proptb (dp,ar(72),in) - 使用TB 命令输入的材料属性值
c (勿使用PROPTB(13), 因为它用在别处) c timval (dp,sc,in) - 当前时间
c timinc (dp,sc,in) - 这一子步中的时间增量
c tem (dp,sc,in) - 这一子步结束时的温度
c dtem (dp,sc,in) - 这一子步结束时的温度增量
c toffst (dp,sc,in) - 相对于绝对零度的温度偏置
c fluen (dp,sc,in) - 这一子步结束时的流量
c dfluen (dp,sc,in) - 这一子步的流量增量
c
c epel (dp,ar(ncomp),inout) - 弹性应变
c epcrp (dp,ar(ncomp),inout) - 前面子步的徐变应变
c statev (dp,ar(ncomp*5+2),inout) - 前面(已收敛)子步中的状态变量
c 该变量仅用于显式徐变,并且是指一 c 个内部变量,但不是用TB,stat 命令定 c 义的用于隐式徐变(用户徐变)以及用 c 户材料的内部变量
c usvr (dp,ar(nuval,nintp),inout) - 附加状态变量,来自前面的迭代(使用 c NSVR 命令保存)
c
c 输出变量:
c 变量 (类型,长度,目的) - 描述
c epel (dp,ar(ncomp),inout) - 考虑徐变增量后的弹性应变调整
(elastic strain adjusted for creep increment)
c epcrp (dp,ar(ncomp),inout) - 更新后的徐变应变
c statev (dp,ar(ncomp*5+2),inout) - 更新后的各状态变量
c usvr (dp,ar(nuval,nintp),inout) - 更新后的附加状态变量
c delcr (dp,sc,out) - 等效徐变应变增量(用于徐变率计算) c
c FORTRAN 参数 (由用户定义):
c 变量 (类型) - 描述
c nuval (int) - 每个积分点处附加状态变量数目
c nintp (int) - 本子程序使用的单元最大积分点数目 c (最大为14)
c 注意: nuval ×nintp = nstv (用nsvr 命令),不能超过840!
c
c 内部变量:
c 变量 (类型,长度) - 描述
c con (dp,sc) - 临时变量
c del (dp,ar(6)) - 徐变应变增量
c epet (dp,sc) - 等效弹性应变 (徐变前)
c ept (dp,ar(6)) - 总应变
c eptot (dp,sc) - 等效总应变, 弹性 + 徐变
c sigen (dp,sc) - 等效应力 (徐变前)
c temabs (dp,sc) - 绝对温度
#Include "impcom.inc"
external erhandler
c
c 用户定义的FORTRAN 参数
c --- usvr数据的容量
integer nuval,nintp
parameter (nuval=1,nintp=1)
c
c 外部子程序和函数
external egen
double precision egen
external vapb,vapb1,vamb1,vmult
c
c 整型变量
integer elem,intpt,mat,ncomp,kfirst,kfsteq
c
c 双精度变量 (&为续行标识符,原文为x ,为清晰起见这里该为&。译者注)
double precision
& e,posn,d(ncomp,ncomp),proptb(72),timval,timinc,tem,dtem,toffst,
& fluen,dfluen,epel(ncomp),epcrp(ncomp),statev(ncomp*5+2),
& usvr(nuval,nintp),delcr,
& temabs,con,del(6),epet,sigen,eptot,ept(6)
c
c **** 当前版本的ANSYS 程序不累计徐变效应(creep work)的影响,如需要,建议按如下步骤进行:
c 1. 在本子程序中计算需要的徐变效应。
C 这可以在单元中累加实现。可以这样做:
c 在子程序的开头 if (intpt .eq. 1) usvr(1,1) = 0.0,
c 然后在所有积分点累加usvr(1,1) = usvr(1,1) + xxxxx
c 注意这仅仅是徐变效应密度,如果需要总徐变效应,需要乘以体积。这c 将在下面的步骤中阐述。
C 事实上,比这要稍微复杂一些,我们希望这个累加仅对收敛的平衡迭 c 代进行,具体实现方法为:
c if (kfstps .eq. 1) work = 0.0
c if (kfsteq .eq. 1) workold = work
c work = workold + delta(work)
c 这样,在收敛的平衡迭代中,效应将是正确的
c 2. 在userou 中,将数据从 usvr(1,1) 移到udbdat(1).
c userou 仅在所有积分点都出来完毕后才被调用
c 3. 在post1, 在nmisc record(各种各样的不可和记录). 的末尾, c 获取udbdat(1)
c 4. 乘以单元体积并累加.
c 这个过程假设徐变效应密度在不规则单元区域内不迅速变化
c
c 这个子程序演示用基本徐变方程计算徐变应变。更详细的资料参见ANSYS 理论手册
3.3。
c 在大应变分析中,所有应变都是hencky 应变. (即对数应变)
c 所有量都是旋转中立(rotation-neutralized )并且是基于旋转单元坐标系统(不是全局)
c
c ***** 初始化 *****
c
c ---初始检查
c --- 当作者进行修改时,应删去下面的警告
if (intpt.eq.1 .and. kfirst.eq.1)
& call erhandler('usercr',5000,2,
& 'ANSYS,INC.- USERCR 的替代版已使用 -'
& ,0.0d0,' ')
c
c ---- 作者应删除下面一行,它仅仅是为分析者提供的
con = d(1,1)+dfluen+dtem+fluen+statev(1)+usvr(1,1)+
& dble(elem+kfsteq+mat)
c ---检查作者定义的 usvr 容量是否超界
c
if (nuval*nintp.gt.840)
&call erhandler('usercr',5010,4,
& ' USVR 超出允许的最大存储量.'
& ,0.0d0,' ')
c
c --- 初始化徐变应变增量以防没有徐变
delcr = 0.0d0
c --- 没有徐变倘若时间没有变化
if (timinc .le. 0.0d0 .and. timval .le. 0.0d0) go to 999
c --- 没有徐变倘若温度未定义
temabs = tem + toffst
if (temabs .le. 0.0d0) then
c
call erhandler('usercr',5020,3,
& ' 徐变Temperature= %G 应大于绝对零度.'
& ,temabs,' ')
c
go to 999
endif
c
c ***** 定义等效应变和应力 *****
c
c --- 利用函数egen 定义等效应变
epet = egen (ncomp,epel(1),posn)
c --- 若应变为0,则无徐变
if (epet .eq. 0.0d0) go to 999
c --- 定义应力
sigen = e*epet
c
c
c ***** 定义徐变应变率 *****
c
c ***** 用户修改正式开始
c ***** 开始徐变方程的例子 - same as c6 = 0
c
c --- 若c1为0则跳过
if (proptb(1) .eq. 0.0d0) go to 999
c
c ---定义等效总应变
call vapb (epel(1),epcrp(1),ept(1),ncomp)
eptot = egen (ncomp,ept(1),posn)
c
c --- 定义徐变应变率
if (eptot .gt. 0.0d0) then
delcr = delcr + exp( log(proptb(1)) + proptb(2)*log(sigen) +
& proptb(3)*log(eptot) - proptb(4)/temabs)
endif
c
c ***** 徐变方程的例子到此结束
c ***** 用户修改正常结束
c
c ***** 为每一个分量计算徐变应变增量 *****
c
c --- 根据徐变应变率计算徐变应变增量
delcr = delcr*timinc
c
c --- 应用prandtl-reuss 关系为每一分量计算增量
if (ncomp .eq. 1) then
del(1) = delcr*epel(1)/epet
else
con = delcr/epet
con = con/(2.0d0*(1.0d0+posn))
call vmult (epel(4),del(4),ncomp-3,3.0d0*con)
del(1) = con*(2.0d0*epel(1) - epel(2) - epel(3))
del(2) = con*(2.0d0*epel(2) - epel(3) - epel(1))
del(3) = con*(2.0d0*epel(3) - epel(1) - epel(2))
endif
c
c --- 更新应变
call vapb1 (epcrp(1),del(1),ncomp)
call vamb1 (epel(1),del(1),ncomp)
c
999 return
end
================================
(2)用户徐变子程序usercreep ,用于隐式徐变
SUBROUTINE usercreep (impflg, ldstep, isubst, matId , elemId,
& kDInPt, kLayer, kSecPt, nstatv, nprop,
& prop , time , dtime , temp , dtemp ,
& toffst, statev, creqv , pres , seqv ,
& delcr , dcrda)
c************************************************************************* c *** 基本功能 ***
c 定义徐变规律,当TB,CREEP 选项 TBOPT=100时使用.
c 演示如何补充用户徐变方程
c
c 徐变方程为
c dotcreq := k0 * seqv ^ n * creqv ^ m * exp (-b/T)
c 其中
c seqv 等效有效应力(Von-Mises 应力)
c creqv 等效有效徐变
c T 温度
c k0, m, n, b 材料常数
c
c 该模型对应于基本徐变函数 TBOPT = 1
c
c gal 10.01.1998 c
c************************************************************************* c
c | (译者注) c |
c | 类型: in -整型
c | dp -双精度型 c | 长度: sc -标量 c | ar( , )-数组 c | 目的: i -输入
c | o -输出 c
c
c impflg (in ,sc ,i)
c ldstep (in ,sc ,i)
c isubst (in ,sc ,i)
c matId (in ,sc ,i)
c elemId (in ,sc ,i)
c kDInPt (in ,sc ,i)
c kLayer (in ,sc ,i)
c kSecPt (in ,sc ,i)
c nstatv (in ,sc ,i)
c nprop (in ,sc ,i)
c prop (dp ,ar(*),i)
c
c time
c dtime
c temp
c dtemp
c toffst (dp, sc, i)
c seqv (dp ,sc , i)
c creqv (dp ,sc , i)
c pres (dp ,sc , i)
c
c
c 输入输出参数
c statev (dp,ar(*), i/o)
c
c | i/o-输入输出 | l -局部变量 显式/隐式积分标识(目前未使用) 当前荷载步 当前子步 材料号 单元号 材料积分点 层号 断面点(Section point) 状态变量数 材料属性数组长度 材料属性数组 该数组传递所有与TB,CREEP 有关的TBDATA 命令定义的在温度temp 时的徐变常数(不要使用prop(13),因为它被用在别处) 当前时间 当前时间增量 当前温度 当前温度增量 温度相对于绝对零度的偏置值 等效有效应力 等效有效徐变应变 静水压力应力, (Sxx+Syy+Szz)/3 用户定义的在时间't' / 't+dt'时的内部状态变量. 该数组将在时间增量的开始,传递这些变量的值。在该子程序中,如果其中的某些状态变量与徐变行为有关,那么它们必须在时间增量的终点被更新。
c delcr (dp ,sc , o) 徐变应变增量
c dcrda (dp,ar(*), o) 输出数组
c dcrda(1) – 徐变应变增量对有效应力的导出量 c dcrda(2) – 徐变应变增量对徐变应变的导出量 derivitive of incremental creep
c strain to creep strain c
c 局部变量
c c1,c2,c3,c4 (dp, sc, l) 临时变量,存储徐变常数
c con1 (dp ,sc, l) 临时变量
c t (dp ,sc, l) 临时变量
c
c************************************************************************* c
c
#include "impcom.inc"
c
DOUBLE PRECISION ZERO
PARAMETER (ZERO = 0.0d0)
c
c ---参数列表
c
INTEGER ldstep, isubst, matId , elemId,
& kDInPt, kLayer, kSecPt, nstatv,
& impflg, nprop
DOUBLE PRECISION dtime , time , temp , dtemp , toffst,
& creqv , seqv , pres
DOUBLE PRECISION prop(*), dcrda(*), statev(nstatv)
c
c ---局部变量
c
DOUBLE PRECISION c1 , c2 , c3 , c4 ,
& con1 , delcr , t
c
c************************************************************************* c
c ***当应力和徐变应变都为0时,跳过
if (seqv.LE.ZERO.AND.creqv.LE.ZERO) GO TO 990
c *** 施加温度偏置(即转化为绝对温度)
t = temp + toffst
c *** 基本徐变方程
c delcr := c1 * seqv ^ n * creqv ^ m * exp (-b/T) * dtime
c1 = prop(1)
c2 = prop(2)
c3 = prop(3)
c4 = prop(4)
c *** 用户需要确定c4为非零值,温度也非零
con1 = ZERO
if(c4.ne.ZERO .and. t.gt.ZERO) con1 = c4/t
c *** 计算徐变应变增量
if (creqv .le. TINY) creqv = sqrt(TINY)
delcr = ZERO
IF(c1.gt.ZERO) delcr = (exp( log(c1) + c2 * log(seqv) + & c3 * log(creqv) - con1 )) * dtime c *** 徐变应变增量对有效应力的导出量
dcrda(1)= c2 * delcr / seqv
c ***徐变应变增量对有效徐变的导出量
dcrda(2)= c3 * delcr / creqv
c *** 将有效徐变应变写入最后的状态变量以便验证
statev(nstatv) = creqv
990 continue
return
end