有限元软件-通用二次开发思路.ppt

单击此处编辑母版标题样式,单击此处编辑母版文本样式,第二级,第三级,第四级,第五级,*,主要内容,1.,通用有限元软件二次开发简介,2.,用户自定义材料本构关系,3.,小结,1,通用有限元软件二次开发简介：,ABAQUS,ABAQUS,二次开发工具,用户子程序：,Fortran,VC,自定义载荷，边界，本构关系，后处理,脚本语言：,PYTHON,对,ABAQUS,功能进行全面的用户更新,2,通用有限元软件二次开发简介,：,ABAQUS,ABAQUS,用户子程序,CREEP,：定义时相关，黏塑性行为（蠕变和膨胀）,DFLOW,：定义固结分析中的非均匀空隙流速,DFLUX,：定义热传导中的非均匀热流,DISP,：定义非均布边界条件,DLOAD,：定义非均布载荷,FILM,：定义固结沉降分析中的非均匀渗流系数,FRIC,：定义接触面中的摩擦行为,3,通用有限元软件二次开发简介,：,ABAQUS,ABAQUS,用户子程序,GAPCON,：定义耦合温度,-,位移分析或纯热传导分析,中接触面的导率,GAPELECTR,：定义耦合热电分析中的导率,HARDINI,：定义初始等效塑性应变和初始背应力张量,HETVAL,：提供热传导分析中的内部热生成,MPC,：定义多点约束,ORIENT,：提供局部材料方向或运动耦合约束中的局部方向或惯性作用的局部刚体方向,4,通用有限元软件二次开发简介,：,ABAQUS,ABAQUS,用户子程序,RSURFU,：定义刚体面,SDVINI,：定义初始求解相关的状态变量,SIGINI,：定义初始应力场,UCORR,：定义随机响应载荷的交叉相关性质,UEL,：定义单元,UEXPAN,：定义热应变增量,UEXTERNALDB,：管理用户定义的外部数据库并计算模型无关的历史信息,5,通用有限元软件二次开发简介,：,ABAQUS,ABAQUS,用户子程序,UFIELD,：定义预定义场变量,UFLUID,：定义静水流体单元的流密度,UFLUIDLEAKOFF,：定义空隙压力粘着单元的流体渗漏系数,UGENS,：定义壳界面的力学行为,UHARD,：定义各向同性或混合硬化的屈服面尺寸和硬化参数,UHYPEL,：定义亚弹性应力应变关系,6,通用有限元软件二次开发简介,：,ABAQUS,ABAQUS,用户子程序,UHYPER,：定义超弹性材料,UINTER,：定义接触面的面交互行为,UMASFL,：定义对流,/,扩散分析中的质量流动率条件,UMAT,：定义材料的力学行为,UMATHT,：定义材料热行为,UMESHMOTION,：指定自适应网格中的网格运动约束,UMULLINS,：定义,Mullins,效应材料模型的损伤变量,7,通用有限元软件二次开发简介,：,ABAQUS,ABAQUS,用户子程序,UPOREP,：定义初始流体空隙压力,UPRESS,：指定等效压应力条件,UPSD,：定义随机响应载荷的频率相关性,URDFIL,：读结果文件,USDFLD,：指定材料点的域变量,UTRACLOAD,：指定非均布牵引载荷,UTRS,：定义粘弹性材料的减缩时间平移函数,UVRAM,：单元输出,8,通用有限元软件二次开发简介,：,ABAQUS,ABAQUS,用户子程序,UWAVE,：定义,ABAQUS/AQUA,分析波运动,VOIDRI,：定义初始空穴比,9,通用有限元软件二次开发简介,：,ABAQUS,ABAQUS,用户子程序示例,1,：不同载荷步间改变弹性模量,Input file,*HEADING,用户自定义损伤弹性模型（,USDFLD,）,*,ELEMENT,TYPE=T2D2,ELSET=ONE,1,1,2,*NODE,1,0.,0.,2,10.,0.,*SOLID SECTION,ELSET=ONE,MATERIAL=ELASTIC,1.,10,通用有限元软件二次开发简介,：,ABAQUS,ABAQUS,用户子程序示例,1,*MATERIAL,NAME=ELASTIC,*ELASTIC,DEPENDENCIES=1,2000.,0.3,0.,0.00,1500.,0.3,0.,0.01,1200.,0.3,0.,0.02,1000.,0.3,0.,0.04,*USER DEFINED FIELD,*DEPVAR,1,11,通用有限元软件二次开发简介,：,ABAQUS,ABAQUS,用户子程序示例,1,*BOUNDARY,1,1,2,2,2,*STEP,*STATIC,0.1,1.0,0.0,0.1,*CLOAD,2,1,20.,*END STEP,12,通用有限元软件二次开发简介,：,ABAQUS,ABAQUS,用户子程序示例,1,*STATIC,0.1,1.0,0.0,0.1,*CLOAD,2,1,0.,*END STEP,*STEP,INC=20,*STATIC,0.1,2.0,0.0,0.1,*CLOAD,2,1,40.,*END STEP,13,通用有限元软件二次开发简介,：,ABAQUS,ABAQUS,用户子程序示例,1,SUBROUTINE USDFLD(FIELD,STATEV,PNEWDT,DIRECT,T,CELENT,TIME,DTIME,CMNAME,ORNAME,NFIELD,NSTATV,NOEL,NPT,LAYER,KSPT,KSTEP,KINC,NDI,NSHR,COORD,JMAC,JMATYP,MATLAYO,LACCFLA),INCLUDE ABA_PARAM.INC,CHARACTER*80 CMNAME,ORNAME,CHARACTER*3 FLGRAY(15),DIMENSION FIELD(NFIELD),STATEV(NSTATV),DIRECT(3,3),1 T(33),TIME(2),DIMENSION ARRAY(15),JARRAY(15),JMAC(*),JMATYP(*),1 COORD(*),14,通用有限元软件二次开发简介,：,ABAQUS,ABAQUS,用户子程序示例,1,C Absolute value of current strain:,CALL GETVRM(E,ARRAY,JARRAY,FLGRAY,JRCD,JMAC,JMATYP,MATLAYO,LACCFLA),EPS=ABS(ARRAY(1),C Maximum value of strain up to this point in time:,CALL GETVRM(SDV,ARRAY,JARRAY,FLGRAY,JRCD,JMAC,JMATYP,MATLAYO,LACCFLA),EPSMAX=ARRAY(1),C Use the maximum strain as a field variable,FIELD(1)=MAX(EPS,EPSMAX),C Store the maximum strain as a solution dependent state,C variable,STATEV(1)=FIELD(1),15,通用有限元软件二次开发简介,：,ABAQUS,ABAQUS,用户子程序示例,1,C If error,write comment to.DAT file:,IF(JRCD.NE.0)THEN,WRITE(6,*)REQUEST ERROR IN USDFLD FOR ELEMENT NUMBER,NOEL,INTEGRATION POINT NUMBER,NPT,ENDIF,C,RETURN,END,16,通用有限元软件二次开发简介,：,ABAQUS,ABAQUS,用户子程序示例,2,：定义蠕变模型,17,通用有限元软件二次开发简介,：,ABAQUS,ABAQUS,用户子程序示例,2,：定义蠕变模型,SUBROUTINE CREEP(DECRA,DESWA,STATEV,SERD,EC,ESW,P,QTILD,1 TEMPDTEMP,PREDEF,DPRED,TIME,DTIME,CMNAME,LEXIMP,LEND,2 COORDS,NSTATV,NOEL,NPT,LAYER,KSPT,KSTEP,KINC),C,INCLUDE ABA_PARAM.INC,C,CHARACTER*80 CMNAME,C,DIMENSION DECRA(5),DESWA(5),STATEV(*),PREDEF(*),DPRED(*),1 TIME(2),COORDS(*),EC(2),ESW(2),18,通用有限元软件二次开发简介,：,ABAQUS,ABAQUS,用户子程序示例,2,：定义蠕变模型,A=,SIG0=,AN=,T1=EXP(QTILD/SIG0),T2=EXP(QTILD/SIG0),DECRA(1)=A*(.5*(T1T2)*AN*DTIME,IF(LEXIMP.EQ.1)THEN,DECRA(5)=AN*A*(.5*(T1T2)*(AN1.)*DTIME/,1 SIG0*.5*(T1+T2),END IF,RETURN,END,19,通用有限元软件二次开发简介,：,ANSYS,ANSYS,二次开发工具,APDL:,通过参数化模型来自动完成一些通用性强的任务；,UIDL,：用户界面设计语言，允许用户设计,ANSYS,图形界面；,UPFs,：,Fortran90,函数及例程以扩展或修改程序的功能，包括定义新材料本构，新单元，新的屈服准则，自定义优化算法，将,ANSYS,作为一个子程序来调用等。,20,通用有限元软件二次开发简介,：,ANSYS,ANSYS-APDL *Do,循环,*do !,起始行,!,循环语句块,*,enddo !,结束行,不允许用,label,分支语句*,if,或*,go,命令跳出,do,循环语句；,不允许用,label,将程序跳到另一行，但可以用,if,then,else,来实现；,do,循环结构中，第一次循环后自动禁止命令结果输出；欲得到所有结果输出，在,do,循环结构中使用,/gopr,或,/go,语句；,/clear,命令不会清除,do,循环的堆栈，但它会删除所有的参数。可在,/clear,命令前运行,/parsav,命令来防止。,21,通用有限元软件二次开发简介,：,ANSYS,ANSYS-APDL,*do,i,0,64,1!,采用,do,循环设置移动载荷,*,set,tim,tim+1,time,tim,nsel,all,fdele,all,all,nsel,all,nsel,s,loc,x,(64+i)*0.5!,施加,5,个集中载荷,f,all,fy,-load,nsel,all,nsel,s,loc,x,(64+i-3)*0.5,*enddo,22,通用有限元软件二次开发简介,：,ANSYS,ANSYS-APDL,CYC1=2!,循环,2圈,*,DO,I,1,CYC1,1!DO,循环开始,LSEL,S,LOC,Y,60,TIME,5,NSUBST,5,AUTOTS,0,SFL,ALL,PRES,PEAK1,ALLSEL,LSWRITE,LSEL,S,LOC,Y,60,AUTOTS,0,NSUBST,5,TIME,10,SFL,ALL,PRES,0,ALLSEL,23,通用有限元软件二次开发简介,：,ANSYS,ANSYS-APDL,LSWRITE,LSEL,S,LOC,Y,60,TIME,12.5,NSUBST,5,AUTOTS,0,SFL,ALL,PRES,VALLEY1,ALLSEL,LSWRITE,LSEL,S,LOCY,60,NSUBST,5,AUTOTS,0,TIME,15,SFL,ALL,PRES,0,ALLSEL,LSWRITE,*ENDDO!DO,循环结束,24,通用有限元软件二次开发简介,：,ANSYS,ANSYS-APDL,：,Repeat,n,1,10,360/(nnode-1),0!,创建,nnode,个节点，夹角为,360/nnode,*repeat,nnode,1,0,360/nnode,0!,重复执行上述命令,nnode,次,25,通用有限元软件二次开发简介,：,ANSYS,ANSYS-UIDL,：单行参数输入,*,ASK,，,Par,，,Query,，,DVAL,其中，,Par,为参数名称，用于存储用户输入的参数。,Query,是询问信息，用户可以输入最多包含,54,个字符串的提示信息以方便正确输入参数。,DVAL,是用户用空响应时程序自动赋给该参数的缺省值。用户用空格响应时则表示删除该参数。,*ASK,RADIUS,INPUT THE RADIUS OF CIRCLE,4,26,通用有限元软件二次开发简介,：,ANSYS,ANSYS-UIDL,：多行参数输入,MULTIPRO,START,PROMPT_NUM,*CSET,STRT_LOC,END_LOC,PARAM_NAME,PROMPT_STRING,DEF_VALUE,MULTIPRO,END,其中，,START,为第一个参数，用于标识,MULTIPRO,指令的开始；,P ROMPT_NUM,为一整型数，等于执行“,MULTIPRO”,命令行后的*,CSET,参数输入提示行的数目，至少有一个*,CSET,命令省略了,DEF_VALUE,参数或,DEF_VALUE,为,0,，才必须用到该参数。,27,通用有限元软件二次开发简介,：,ANSYS,ANSYS-UIDL,：多行参数输入,MULTIPRO,START,3,*CSET,1,3,EX_1,YOUNGS MODULUS(MPa),2.06E5,*CSET,4,6,NUXY_1,POISSIONS RATIO,0.3,*CSET,7,9,DENS_1DENSITY(103Kg/mm3),7.85e-9,*CSET,61,62,ENTER THE ATTRIBUTES OF,MATERIAL 1,*CSET,63,64,NOTE:UNIT OF LENGHTH IS MM!,MULTIPRO,END,/PREP7,MP,EX,1,EX_1,MP,NUXY,1,NUXY_1,MP,DENS1,1,DENS_1,28,通用有限元软件二次开发简介,：,ANSYS,ANSYS-UIDL,：多行参数输入,29,通用有限元软件二次开发简介,：,ANSYS,ANSYS-UIDL,：自定义工具栏,C,：,Program FilesAnsys Incv90ANSYSapdlstart90.ans!,视安装路径而定用文本编辑器打开后，在最后一行追加如下代码：,!*,/psearch,e:ansys !,宏文件存放路径,*,ABBR,Config,Kconfig,*ABBR,etype,Ketype,*ABBR,Mp,Kmp,*ABBR,Biso,Kbiso,*ABBR,Chaboche,kchaboche,*ABBR,Model,Kmodel,*ABBR,Load,Kload,*ABBR,Solve,KSolve,*ABBR,SS_center,Kscenter,*ABBR,SS_edge,Kedge,*ABBR,Avi_seqv,Kflash,30,通用有限元软件二次开发简介,：,ANSYS,ANSYS-UIDL,：多行参数输入,ANSYS,高级工程应用实例与二次开发,31,通用有限元软件二次开发简介,：,MARC,MSC.MARC,二次开发工具,用户子程序：,Fortran,自定义载荷，边界，本构关系，后处理,脚本语言：,PYTHON,用户开发新的软件功能，仅调用,MARC,的求解器。,32,通用有限元软件二次开发简介,：,MARC,subroutine ufxord(xord,ncrd,n),implicit real*8(a-h,o-z)dp,dimension xord(ncrd),r=24.d0/3.14159d0,t=xord(2)/r,xord(2)=r*sin(t),xord(3)=r*cos(t),write(6,1)n,(xord(k),k=1,ncrd),1 format(i5,3e13.5),return,end,33,通用有限元软件二次开发简介,：,MARC,34,二用户子程序定义本构模型,ABAQUS,：,Plastic,模型包括了各向同性硬化、随,动硬化混合模型，,Prager,模型，,AF,模,型等,UMAT,ANSYS,：刚塑性、双线性、多段线性模型，,塑性,Chaboche,模型，,ANADE,模型,非,线性塑性等,USERMAT,MARC,：刚塑性、双线性、多段线性模型，粘,塑性,Chaboche,模型等,HYPERLA2,35,本构模型简介,在材料力学行为的数学描述中，通常通过本构方程来表征材料的反应，即在本构方程中给出应力作为物体变形历史的函数。这些本构方程即为本构模型。,不同的材料力学行为需要不同的本构模型来描述，如粘性流体、橡胶和混凝土，它们的本构模型截然不同。,在一维固体力学中，本构关系经常归属于材料的应力,-,应变关系。,36,本构模型简介,描述固体材料力学行为的本构模型大致归类：,弹性,非弹性,率无关弹性,率相关弹性,率无关非弹性,率相关非弹性,线弹性,非线性弹性,粘弹性,亚弹性,超弹性,(,伪弹性或拟塑性,),双线性弹塑性,多线性弹塑性,非线性弹塑性各向同性硬化,随动硬化,蠕变,粘塑性,各向同性硬化,随动硬化,固体材料唯象本构模型,仅蠕变,体积膨胀蠕变,非金属塑性,混凝土,岩土,1,2,3,4,非各向同性硬化,非各向同性硬化,超塑性,37,线弹性本构模型的有限元实现,路径无关,可逆,无能量耗散,单轴状态下：,多轴状态下：,38,线弹性本构模型的有限元实现,其中，,有,81,个独立常数，与应力张量的,9,个分量对应。,应力和应变的对称性要求应力的,6,个独立分量与应变的,6,个独立分量有关。独立常数减少为,36,个。即,66,矩阵。,又由于矩阵的主对称性，独立弹性常数减少为,n(n+1)/2=21,个。,39,线弹性本构模型的有限元实现,由广义胡克定律可知：,40,线弹性本构模型的有限元实现,本构模型的离散：,41,线弹性本构模型的有限元实现,SUBROUTINE UMAT(,STRESS,STATEV,DDSDDE,SSE,SPD,SCD,1 RPL,DDSDDT,DRPLDE,DRPLDT,2,STRAN,DSTRAN,TIME,DTIME,TEMP,DTEMP,PREDEF,DPRED,CMNAME,3 NDI,NSHR,NTENS,NSTATV,PROPS,NPROPS,COORDS,DROT,PNEWDT,4 CELENT,DFGRD0,DFGRD1,NOEL,NPT,LAYER,KSPT,KSTEP,KINC),CINCLUDE ABA_PARAM.INC,C,以,ABAQUS,为例，调用用户材料子程序,UMAT,42,线弹性本构模型的有限元实现,CHARACTER*80 CMNAME,DIMENSION STRESS(NTENS),STATEV(NSTATV),1 DDSDDE(NTENS,NTENS),2 DDSDDT(NTENS),DRPLDE(NTENS),3 STRAN(NTENS),DSTRAN(NTENS),TIME(2),PREDEF(1),DPRED(1),4 PROPS(NPROPS),COORDS(3),DROT(3,3),DFGRD0(3,3),DFGRD1(3,3),DIMENSION DSTRES(6),D(3,3),C,C PROPS(1)E,弹性模量,C PROPS(2)NU,泊松比,C,NTENS,：总的应力或应变分量个数，三维问题为,6,个，平面和轴对称问题为,4,个，一维问题为,1,个。,NTENS=NDI+NSHR,43,线弹性本构模型的有限元实现,Cif(NDI.NE.3)THEN,CWRITE(7,*)THIS UMAT BE USED FOR 3D,CCALL XIT,CENDIF,C,C ELASTIC PROPERTIES,PARAMETER(ONE=1.0D0,TWO=2.0D0,THREE=3.0D0),EMOD=PROPS(1),E,ENU=PROPS(2),v,EBULK3=EMOD/(ONE-TWO*ENU),K,EG2=EMOD/(ONE+ENU),2G,EG=EG2/TWO,G,EG3=THREE*EG,3G,ELAM=(EBULK3-EG2)/THREE,C,计算弹性常数,44,线弹性本构模型的有限元实现,C,C ELASTIC STIFFNESS,C,DO K1=1,NDI,DO K2=1,NDI,DDSDDE(K2,K1)=ELAM,END DO,DDSDDE(K1,K1)=EG2+ELAM,END DO,DO K1=NDI+1,NTENS,DDSDDE(K1,K1)=EG,END DO,C,计算弹性刚度矩阵C=,DDSDDE(ntens,ntens),45,线弹性本构模型的有限元实现,C,C CALCULATE STRESS,C,DO K1=1,NTENS,DO K2=1,NTENS,STRESS(K2)=STRESS(K2)+DDSDDE(K2,K1)*DSTRAN(K1),END DO,END DO,C,RETURN,END,计算当前应力,46,线弹性本构模型的有限元实现,1/8,单元验证（,C3D8,）,*Heading,*Job name:Linearelastic Model name:C3D8,*Preprint,echo=NO,model=NO,history=NO,contact=NO,*PARTS,*Part,name=PART-1-1,*Node,1,0.,0.,0.,2,1.,0.,0.,3,1.,1.,0.,4,0.,1.,0.,5,0.,0.,1.,6,1.,0.,1.,7,1.,1.,1.,8,0.,1.,1.,*Element,type=C3D8,1,1,2,3,4,5,6,7,8,*Elset,elset=ALLE,1,*Solid Section,elset=ALLE,material=ALLE,1.,*End Part,47,线弹性本构模型的有限元实现,1/8,单元验证（,C3D8,）,*ASSEMBLY,*,*Assembly,name=Assembly,*,*Instance,name=PART-1-1,part=PART-1-1,*End Instance,*,*Nset,nset=ALLN,instance=PART-1-1,generate,1,8,1,*Nset,nset=,N1,internal,instance=PART-1-1,4,*Nset,nset=,N2,internal,instance=PART-1-1,5,*Nset,nset=,N3,internal,instance=PART-1-1,1,*Nset,nset=,N4,internal,instance=PART-1-1,4,5,8,*Elset,elset=,E1,internal,instance=PART-1-1,1,*Surface,type=ELEMENT,name=,S1,internal,E1,S4,*End Assembly,48,线弹性本构模型的有限元实现,1/8,单元验证（,C3D8,）,*Amplitude,name=AMP-1,0.,0.,1.,1.,2.,0.,*,*MATERIALS,*,*Material,name=ALLE,*Depvar,0,*User Material,constants=2,2,.,e6,0.3,49,线弹性本构模型的有限元实现,1/8,单元验证（,C3D8,）,*BOUNDARY CONDITIONS,*,*Name:Disp-BC-1 Type:Displacement/Rotation,*Boundary,N1,3,3,*Name:Disp-BC-2 Type:Displacement/Rotation,*Boundary,N2,2,2,*Name:Disp-BC-3 Type:Displacement/Rotation,*Boundary,N3,1,1,*Name:Disp-BC-4 Type:Displacement/Rotation,*Boundary,N3,2,2,*Name:Disp-BC-5 Type:Displacement/Rotation,*Boundary,N3,3,3,*Name:Disp-BC-6 Type:Displacement/Rotation,*Boundary,N4,1,1,50,线弹性本构模型的有限元实现,1/8,单元验证（,C3D8,）,*-,*,*STEP:Step-1,*,*Step,name=Step-1,inc=200,*Static,0.1,2.,0.1,0.1,*,*LOADS,*,*Name:SURFFORCE-1 Type:Pressure,*Dsload,amplitude=AMP-1,S1,P,-200.,*-,51,线弹性本构模型的有限元实现,1/8,单元验证（,C3D8,）,*OUTPUT REQUESTS,*,*Restart,write,frequency=0,*,*FIELD OUTPUT:F-Output-1,*,*Output,field,*Element Output,E,S,*,*HISTORY OUTPUT:H-Output-1,*,*Output,history,variable=PRESELECT,*End Step,52,弹塑性本构模型的有限元实现,首先以双线性弹塑性本构模型为例，阐述弹塑性有限元实现的过程。,应变率分解：,弹性本构方程：,屈服条件：,塑性流动法则：,加卸载条件：,一致性条件：,应力率,-,应变率关系：,连续体弹塑性切线模量：,关联流动法则，一般取53,弹塑性本构模型的有限元实现,径向回退算法：,求当前应力：,54,弹塑性本构模型的有限元实现,径向回退算法：,针对某一具体的屈服面，如：,55,弹塑性本构模型的有限元实现,连续体弹塑性切线模量：,56,弹塑性本构模型的有限元实现,径向回退算法,1,1.,设初始值,2.,第,K,次迭代时,检查屈服条件,收敛,4.,更新塑性应变和内变量,5.,计算连续体弹塑性切线模量,Else,转到第,3,步。,IF,3.,计算塑性参数的增量,k=k+1,，转到第,2,步,57,弹塑性本构模型的有限元实现,径向回退算法,2,1.,设初始值,2.,检查屈服条件,弹性状态,4.,更新塑性应变和内变量,5.,计算连续体弹塑性切线模量,Else,转到第,3,步。,IF,3.,计算塑性参数的增量,收敛,IF,Else,转到第,3,步。,ABAQUS,58,弹塑性本构模型的有限元实现,(,各向同性硬化,),SUBROUTINE UMAT(STRESS,STATEV,DDSDDE,SSE,SPD,SCD,1 RPL,DDSDDT,DRPLDE,DRPLDT,STRAN,DSTRAN,2 TIME,DTIME,TEMP,DTEMP,PREDEF,DPRED,MATERL,NDI,NSHR,NTENS,3 NSTATV,PROPS,NPROPS,COORDS,DROT,PNEWDT,CELENT,4 DFGRD0DFGRD1,NOEL,NPT,KSLAY,KSPT,KSTEP,KINC)C,INCLUDE ABA_PARAM.INCC,CHARACTER*80 MATERL,DIMENSION STRESS(NTENS),STATEV(NSTATV),1 DDSDDE(NTENS,NTENS),DDSDDT(NTENS),DRPLDE(NTENS),2 STRAN(NTENS),DSTRAN(NTENS),TIME(2),PREDEF(1),DPRED(1),3 PROPS(NPROPS),COORDS(3),DROT(3,3),4 DFGRD0(3,3),DFGRD1(3,3),59,弹塑性本构模型的有限元实现,(,各向同性硬化,),C DIMENSION EELAS(6),EPLAS(6),FLOW(6),PARAMETER(ONE=1.0D0,TWO=2.0D0,THREE=3.0D0,SIX=6.0D0),DATA NEWTON,TOLER/10,1.D-6/,C UMAT FOR ISOTROPIC ELASTICITY AND ISOTROPIC PLASTICITYC J2 FLOW THEORY,C CAN NOT BE USED FOR PLANE STRESS,C,C PROPS(1)E,E,C PROPS(2)NU,V,C PROPS(3)SYIELD,y,C PROPS(4)-HARDENING MODULUS,H,C-,60,弹塑性本构模型的有限元实现,C,IF(NDI.NE.3)THEN,WRITE(6,1),1 FORMAT(/,30X,*ERROR-THIS UMAT MAY ONLY BE USED FOR,ELEMENTS WITH THREE DIRECT STRESS COMPONENTS),ENDIF,C,C ELASTIC PROPERTIESC,EMOD=PROPS(1),ENU=PROPS(2),SYIEL0=PROPS(3),HARD=PROPS(4),IF(ENU.GT.0.4999.AND.ENU.LT.0.5001)ENU=0.499,EBULK3=EMOD/(ONE-TWO*ENU),EG2=EMOD/(ONE+ENU),EG=EG2/TWO,EG3=THREE*EG,ELAM=(EBULK3-EG2)/THREE,61,弹塑性本构模型的有限元实现,C ELASTIC STIFFNESS,C,DO 20 K1=1,NTENS,DO 10 K2=1,NTENS,DDSDDE(K2,K1)=0.0,10 CONTINUE,CONTINUE,C,DO 40 K1=1,NDI,DO 30 K2=1,NDI,DDSDDE(K2,K1)=ELAM,CONTINUE,DDSDDE(K1,K1)=EG2+ELAM,40 CONTINUE,DO 50 K1=NDI+1,NTENS,DDSDDE(K1,K1)=EG,50 CONTINUE,62,弹塑性本构模型的有限元实现,C CALCULATE STRESS FROM ELASTIC STRAINSC,DO 70 K1=1,NTENS,DO 60 K2=1,NTENS STRESS(K2)=STRESS(K2)+DDSDDE(K2,K1)*DSTRAN(K1),CONTINUE,70 CONTINUE,C RECOVER ELASTIC AND PLASTIC STRAINS,DO 80 K1=1,NTENS,EELAS(K1)=STATEV(K1)+DSTRANK1)EPLAS(K1)=STATEV(K1+NTENS),CONTINUE,EQPLAS=STATEV(1+2*NTENS),弹性应变,塑性应变,等效塑性应变,63,C MISES STRESSC,SMISES=(STRESS(1)-STRESS(2)*(STRESS(1)-STRESS(2)+,(STRESS(2)-STRESS(3)*(STRESS(2)-STRESS(3)+,(STRESS(3)-STRESS(1)*(STRESS(3)-STRESS(1),DO 90 K1=NDI+1,NTENS,SMISES=SMISES+SIX*STRESS(K1)*STRESS(K1),CONTINUE,SMISES=SQRT(SMISES/TWO),IF(SMISES.GT.(1.0+TOLER1)*SYIEL0)THEN,C FLOW DIRECTIONC,SHYDRO=(STRESS(1)+STRESS(2)+STRESS(3)/THREE,ONESY=ONE/SMISES,DO 110 K1=1,NDI,FLOW(K1)=ONESY*(STRESS(K1)-SHYDRO),CONTINUE,DO 120 K1=NDI+1,NTENS,FLOW(K1)=STRESS(K1)*ONESY,120 CONTINUE,弹塑性本构模型的有限元实现,64,C SOLVE FOR EQUIV STRESS,NEWTON ITERATION,C SYIELD=SYIEL0,DEQPL=0.0,DO 130 KEWTON=1,NEWTON,RHS=SMISES-EG3*DEQPL-SYIELD,DEQPL=DEQPL+RHS/(EG3+HARD),EQPLAS=EQPLAS+DEQPL,SYIELD=SYIEL0+HARD*DEQPL,IF(ABS(RHS).LT.TOLER2*SYIEL0)GOTO 140,CONTINUE,WRITE(6,2)NEWTON 2 FORMAT(/,30X,*WARNING-PLASTICITY ALGORITHM DID NOT,CONVERGE AFTER,I3,ITERATIONS),140 CONTINUE,弹塑性本构模型的有限元实现,65,C CALC STRESS AND UPDATE STRAINSC,DO 150 K1=1,NDI,STRESS(K1)=FLOW(K1)*SYIELD+SHYDRO EPLAS(K1)=EPLAS(K1)+THREE*FLOW(K1)*DEQPL/TWO EELAS(K1)=EELAS(K1)-THREE*FLOW(K1)*DEQPL/TWO,CONTINUE,DO 160 K1=NDI+1,NTENS,STRESS(K1)=FLOW(K1)*SYIELD EPLAS(K1)=EPLASK1)+THREE*FLOW(K1)*DEQPL,EELAS(K1)=EELAS(K1)-THREE*FLOW(K1)*DEQPL,CONTINUE,EQPLAS=EQPLAS+DEQPL,弹塑性本构模型的有限元实现,66,C JACOBIAN,EFFHRD=EG3*HARD/(EG3+HARD),EFFG=EG*SYIELD/SMISES,EFFG2=TWO*EFFG,EFFG3=THREE*EFFG2/TWO,EFFLAM=(EBULK3-EFFG2)/THREE,DO 220 K1=1,NDI,DO 210 K2=1,NDI,DDSDDE(K2,K1)=EFFLAM,CONTINUE,DDSDDE(K1,K1)=EFFG2+EFFLAM,CONTINUE,DO 230 K1=NDI+1,NTENS,DDSDDE(K1,K1)=EFFG,CONTINUE,DO 250 K1=1,NTENS,DO 240 K2=1,NTENS,DDSDDE(K2,K1)=DDSDDE(K2,K1)+FLOW(K2)*FLOW(K1),*(EFFHRD-EFFG3),CONTINUE,CONTINUE,ENDIF,ENDIF,弹塑性本构模型的有限元实现,67,C STORE STRAINS IN STATE VARIABLE ARRAY,C,DO 310 K1=1,NTENS,STATEV(K1)=EELAS(K1)STATEV(K1+NTENS)=EPLAS(K1),CONTINUE,STATEV(,1+2*NTENS,)=EQPLASC,RETURN,END,C-,ABAQUS,的,INP,文件与线性的,INP,类似，只需更改材料常数和状态变量。,弹塑性本构模型的有限元实现,*MATERIALS,*,*Material,name=ALLE,*Depvar,13,*User Material,constants=,4,2,.,e6,0.3,200.,5800.,68,作为有生命力的通用软件，用户可开发性,的友好与否至关重要。通常在完成一项,工作时，需要 几个软件协同工作，扬长,避短