ABAQUS子程序.docx

1、ABAQUS子程序ABAQUS用户子程序当用到某个用户子程序时，用户所关心的主要有两方面：一是ABAQUS提供的用户子程序的接口参数。有些参数是ABAQUS传到用户子程序中的，例如SUBROUTINE DLOAD中的KSTEP，KINC，COORDS；有些是需要用户自己定义的，例如F。二是ABAQUS何时调用该用户子程序，对于不同的用户子程序ABAQUS调用的时间是不同的。有些是在每个STEP的开始，有的是STEP结尾，有的是在每个INCREMENT的开始等等。当ABAQUS调用用户子程序是，都会把当前的STEP和INCREMENT利用用户子程序的两个实参KSTEP和KINC传给用户子程序，用

2、户可编个小程序把它们输出到外部文件中，这样对ABAQUS何时调用该用户子程序就会有更深的了解。（子程序中很重要的就是要知道由abaqus提供的那些参量的意义，如下）首先介绍几个子程序：一SUBROUTINE DLOAD(F,KSTEP,KINC,TIME,NOEL,NPT,LAYER,KSPT,COORDS, JLTYP,SNAME) 参数： 1 F为用户定义的是每个积分点所作用的荷载的大小； 2 KSTEP,KINC为ABAQUS传到用户子程序当前的STEP和INCREMENT值； 3 TIME(1),TIME(2)为当前STEP TIME和INCREMENT TIME的值； 4 NOEL,

3、NPT为积分点所在单元的编号和积分点的编号； 5 COORDS为当前积分点的坐标； 6 除F外，所有参数的值都是ABAQUS传到用户子程序中的。 功能： 1 荷载可以被定义为积分点坐标、时间、单元编号和单元节点编号的函数。 2 用户可以从其他程序的结果文件中进行相关操作来定义积分点F的大小。 例1：这个例子在每个积分点施加的荷载不仅是坐标的函数，而且是随STEP变化而变化的。 SUBROUTINE DLOAD(P,KSTEP,KINC,TIME,NOEL,NPT,LAYER,KSPT,COORDS, 1 JLTYP,SNAME) INCLUDE ABA_PARAM.INC C DIMENSIO

4、N TIME(2),COORDS(3) CHARACTER*80 SNAME PARAMETER (PLOAD=100.E4) IF (KSTEP.EQ.1) THEN ！当STEP=1时的荷载大小 P=PLOAD ELSE IF (KSTEP.EQ.2) THEN ！当STEP2时的荷载大小 P=COORDS(1)*PLOAD ！施加在积分点的荷载P是坐标的函数 ELSE IF (KSTEP.EQ.3) THEN ！当STEP=3时的荷载大小 P=COORDS(1)*2*PLOAD ELSE IF (KSTEP.EQ.4) THEN ！当STEP=4时的荷载大小 P=COORDS(1)*3*

5、PLOAD ELSE IF (KSTEP.EQ.5) THEN ！当STEP=5时的荷载大小 P=COORDS(1)*4*PLOAD END IF RETURN END UMAT 子程序具有强大的功能，使用UMAT 子程序：(1) 可以定义材料的本构关系，使用ABAQUS 材料库中没有包含的材料进行计算，扩充程序功能。(2) 几乎可以用于力学行为分析的任何分析过程，几乎可以把用户材料属性赋予ABAQUS 中的任何单元；(3) 必须在UMAT 中提供材料本构模型的雅可比（Jacobian）矩阵，即应力增量对应变增量的变化率。(4) 可以和用户子程序“USDFLD”联合使用，通过“USDFLD”重

6、新定义单元每一物质点上传递到UMAT 中场变量的数值。由于主程序与UMAT 之间存在数据传递，甚至共用一些变量，因此必须遵守有关UMAT 的书写格式，UMAT 中常用的变量在文件开头予以定义，通常格式为：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,

7、DFGRD0,DFGRD1,NOEL,NPT,LAYER,KSPT,KSTEP,KINC)INCLUDE ABA_PARAM.INCCHARACTER*80 CMNAMEDIMENSION 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),DFGRD0(3,3),DFGRD1(3,3)user coding

8、 to define DDSDDE, STRESS, STATEV, SSE, SPD, SCDand, if necessary, RPL, DDSDDT, DRPLDE, DRPLDT, PNEWDTRETURNENDUMAT 中的应力矩阵、应变矩阵以及矩阵DDSDDE ，DDSDDT ，DRPLDE 等，都是直接分量存储在前，剪切分量存储在后。直接分量有NDI 个，剪切分量有NSHR 个。各分量之间的顺序根据单元自由度的不同有一些差异，所以编写UMAT 时要考虑到所使用单元的类别。下面对UMAT 中用到的一些变量进行说明：DDSDDE（ NTENS， NTENS）是一个NTENS 维的方

9、阵，称作雅可比矩阵， ， 是应力的增量， 是应变的增量， DDSDDE（ I，J）表示增量步结束时第J 个应变分量的改变引起的第I 个应力分量的变化。通常雅可比是一个对称矩阵，除非在“*USER MATERIAL”语句中加入了“UNSYMM”参数。STRESS （NTENS）应力张量矩阵，对应NDI 个直接分量和NSHR 个剪切分量。在增量步的开始，应力张量矩阵中的数值通过UMAT 和主程序之间的接口传递到UMAT 中，在增量步的结束UMAT 将对应力张量矩阵更新。对于包含刚体转动的有限应变问题，一个增量步调用UMAT 之前就已经对应力张量的进行了刚体转动，因此在UMAT 中只需处理应力张量的

10、共旋部分。UMAT 中应力张量的度量为柯西（真实）应力。STATEV （NSTATEV）用于存储状态变量的矩阵，在增量步开始时将数值传递到UMAT 中。也可在子程序USDFLD或UEXPAN 中先更新数据，然后增量步开始时将更新后的数据传递到UMAT 中。在增量步的结束必须更新状态变量矩阵中的数据。和应力张量矩阵不同的是：对于有限应变问题，除了材料本构行为引起的数据更新以外，状态变量矩阵中的任何矢量或者张量都必须通过旋转来考虑材料的刚体运动。状态变量矩阵的维数，等于关键字“*DEPVAR”定义的数值。状态变量矩阵的维数通过ABAQUS 输入文件中的关键字“*DEPVAR”定义，关键字下面数据行

11、的数值即为状态变量矩阵的维数。材料常数的个数，等于关键字“*USER MATERIAL”中“CONSTANTS”常数设定的值。PROPS （NPROPS）材料常数矩阵，矩阵中元素的数值对应于关键字“*USER MATERIAL”下面的数据行。SSE ， SPD ， SCD分别定义每一增量步的弹性应变能，塑性耗散和蠕变耗散。它们对计算结果没有影响，仅仅作为能量输出。其他变量：STRAN( NTENS) ：应变矩阵；DSTRAN( NTENS) ：应变增量矩阵；DTIME ：增量步的时间增量；NDI ：直接应力分量的个数；NSHR ：剪切应力分量的个数；NTENS ：总应力分量的个数， NTENS

12、 NDI NSHR = + 。使用UMAT 时需要注意单元的沙漏控制刚度和横向剪切刚度。通常减缩积分单元的沙漏控制刚度和板、壳、梁单元的横向剪切刚度是通过材料属性中的弹性性质定义的。这些刚度基于材料初始剪切模量的值，通常在材料定义中通过“*ELASTIC”选项定义。但是使用UMAT 的时候，ABAQUS 对程序输入文件进行预处理的时候得不到剪切模量的数值。所以这时候用户必须使用“*HOURGLASS STIFFNESS” 选项来定义具有沙漏模式的单元的沙漏控制刚度， 使用“*TRANSVERSE SHEAR STIFFNESS”选项来定义板、壳、梁单元的横向剪切刚度。几个关于子程序的问题及相应

13、解答Q: 本人在用umat作本构模型时， *static, 1,500,0.000001,0.1 此时要求的增量步很多，即每次增量要很小， *static 1,500 时，在弹性向塑性过度时，出现错误，增量过大，出现尖点.? A: YOU CAN TRY AS FOLLOWS: *STEP,EXTRAPOLATION=NO,INC=2000000 *STATIC 0.001,500.0,0.00001,0.1。 Q: 在abaqus中，如果采用umat,利用自己的本构，如何让abaqus明白这种材料的弹塑性应 变，也就是说，如何让程序返回弹性应变与塑性应变，好在output中输出，我曾想用最笨

14、 地方法，在uvarm中定义输出，利用getvrm获取材料点的值，但无法获取增量应力，材料常 数等，研究了帮助中的例子，umatmst3.inp,umatmst3.for,他采用mises J2 流动理论，我 在output history 显示他已进入塑性状态，但他的PE仍然为0！? A: 用uvar( )勉强成功 。 Q: 偶在umat中调用求主应力函数 CALL SPRINC(STRESS,PS,LSTR,NDI,NSHR) 后，存储主应力得数组PS中 各个主应力排列顺序是什么？ PS1PS2PS3 ? PS1PS2PS3PS2 ? A: 第二个 。 个人觉得： umat实现自己的本构没

15、有固定的方法，对于不同的本构有可能必须采用不同的方法。这要 靠自己不断地摸索。有可能一种方法对于简单加载问题还行，但有可能对于复杂问题并不 收敛。最重要一点，就是umat中采用的算法必须consistent.再就是ddsdde必须正确，（如 果采用back_Euler 方法等一些算法，ddsdde错误有时不影响结果（对于简单加载问题没有 影响，能收敛，），但对于复杂问题不收敛。 uptonow，你这个算法对于Mises,hill,J2，J2d等一类的屈服函数是正确的，但具体的本构 还要灵活运用，这我也正学习，正在摸索。 有时，umat需要很强的有限元基础，并且对采用的本构要很熟悉，不要在一颗树

16、上吊死才 好。首先要确认自己的umat没有错误，如果没有，但就是不收敛（在不断减小加载步长的 情况下，当然最好对步长不敏感，特别是对于粘弹性，粘塑性，内变量一类的材料，有的 本构取决于背应力的计算）。 那就应该考虑换一种算法。 一点体会，请大家探讨。Q: abaqus-uamt的老问题，缺少 ABA_PARAM.INC文件 ? A: 在cvf6.5调试时，显示 缺少 ABA_PARAM.INC文件！ 这个没有任何关系的，这个错误将在ABAQUS调用UMAT的时候自动会找到，仅仅有这个错误 将没有任何影响的。也就是说，ABAQUS中调用的时候，实际并不存在这个错误。FT，忘了 说一句了，你把AB

17、A_PARAM.INC.dp或ABA_PARAM.INC.sp拷到你的程序工作空间后 ，应该将 把ABA_PARAM.INC.dp或ABA_PARAM.INC.sp的后缀.sp或.dp去掉，即将ABA_PARAM.IN C.dp或 ABA_PARAM.INC.sp改名为ABA_PARAM.INC。 呵呵，他的意思是在Visual Fortran中调试其子程序，我觉得这是一个好办法，我当时也 是这麽办的，毕竟在ABAQUS中调试是非常麻烦的，只有当你的UMAT没有语法或者明显的逻 辑错误，你在ABAQUS中调试才能事半功倍。Q：uvarm可以输出到哪里？.odb可以么？ 另外那个strav？就是

18、自己定义用于umat的那个数组里的数可以输出到.odb里么？ A：在umat中，statev是不能用在output中的，statev只是作为一个解的 状态变量，说来惭愧，我是在umat中定义peeq,(peeq)的求解一般在 弹塑性力学书上有（等效塑性应变），用write()写入一个临时文件， （切记：这个文件unit号不要与abaqus中的重合，因为他有一些系统默认 的文件号，） 然后在uvarm中读取，以uvarm输出，因为uvarm可以以odb的形式输出， 支持output,field,output,history Q：在本版看了一个一维固结的例子，其中含有用户子程序，如下 SUBROU

19、TINE UFIELD(FIELD,KFIELD,NSECPT,KSTEP,KINC,TIME,NODE, COORDS,TEMP,DTEMP) INCLUDE ABA_PARAM.INC DIMENSION FIELD(NSECPT),TIME(2),COORDS(3),TEMP(NSECPT), DTEMP(NSECPT) KFIELD=1 FIELD(1)=COORDS(2) RETURN END 我在VF6.5中进行调试，提示找不到ABA_PARAM.INC 请问大侠这如何解决。 还有我因为是初次接触用户子程序，我查阅了本版所有的相关贴子，都讲的不太详细，我将问题总结一下，大侠们能不能

20、详细的讲解一下， 1子程序格式(程序后缀是.f; .f90; .for;.obj?) 2CAE中如何调用，command下如何调用？ 3若有多个子程序同时存在，如何处理 4我对VF不是很熟，是否可以用VC，C编写子程序？A: 若要在vf中调试，那么应该根据需要把SITE文件夹中的ABA_PARAM_DP.INC（双精度）或ABA_PARAM_SP.INC（单精度）拷到相应的位置，并改名为ABA_PARAM.INC即可。 1。我试过，.for格是应该是不可以的，至少6.2和6.3版本应该是不行，其他的没用过，没有发言权。 在Abaqus中，运行abaqus j=jobname user=user

21、name时，默认的用户子程序后缀名是.for（.f,.f90应该都不行的，手册上也有讲过），只有在username.for文件没有找到的情况下，才会去搜索username.obj，如果两者都没有，就会报错误信息。 如果username包括扩展名for或obj，那么就根据各自的扩展名ABAQUS会自动选择进行操作。 2。cae中在creat job的job manager中的general中可以指定子程序；command下用命令:abaqus j=jobname user=userfilename （无后缀）； 3。将其写在一个文件中即可，然后用一个总的子程序调用（具体参见手册） 4。据说6.4

22、的将可以，6.3的你可以尝试着将VC，C程序编译为obj文件，没试过。在你的工作目录下应该已经存在ufield.obj和uvarm.obj这两个文件（这两个文件应该是你分别单独调试ufield.FOR和uvarm.FOR时自动编译生成的，你可以将他们删掉试试看），但是由于你的FOR文件中已经有了UVARM和UFIELD这两个subroutine，显然会造成重复定义，请查实。 Q: 假定采用mises屈服准则。 1 在调用UMAT之前，ABAQUS传递给UMAT本次增量开始时的应力sigma(0)，总应变E,应变增 量delta(E)。 状态变量保存：弹性应变，塑性应变，等效塑性应变。 （请问，

23、状态变量保存的弹性应变塑性应变是否等于ABAQUS传递给UMAT的总应变？ ？） 2，然后在UMAT中利用上述的ABAQUS传递的量和状态变量得到DDSDDE矩阵，然后返回给 ABAQUS，ABAQUS根据delta(sigma)=ddsdde*delta(E)，并且得到本次增量结束时的应力 sigma(1)=sigma(0)+delta(sigma)=sigma(0)+ddsdde*delta(E) 3 然后更新本次增量结束时的状态变量：弹性应变，塑性应变，等效塑性应变以供下次调 用UMAT 请问手册上UMAT必须更新应力，可是根据上述我的理解好像是ABAQUS根据UMAT提供的本 次增量的

24、DDSDDE在ABAQUS中更新，请问到底是怎么回事？谢谢！ A: 1 在调用UMAT之前，ABAQUS传递给UMAT本次增量开始时的应力sigma(0)，总应变E,应. : 量delta(E)。 : 状态变量保存：弹性应变，塑性应变，等效塑性应变。 2，然后在UMAT中利用上述的ABAQUS传递的量和状态变量得到DDSDDE矩阵，然后返回给 : ABAQUS，ABAQUS根据delta(sigma)=ddsdde*delta(E)，并且得到本次增量结束时的应力 3 然后更新本次增量结束时的状态变量：弹性应变，塑性应变，等效塑性应变以供下以上有些answers本人并没有亲自证实，如有问题请及时

25、更正！abaqus学习经验第一章 ABAQUS 简介1 (pp7)在开始 程序 ABAQUS 6.5-1ABAQUS COMMAND，DOS 提示符下输入命令Abaqus fetch job = 可以提取想要的算例input 文件。第二章 ABAQUS 基本使用方法2(pp15)快捷键：Ctrl+Alt+左键来缩放模型；Ctrl+Alt+中键来平移模型；Ctrl+Alt+右键来旋转模型。(pp16)ABAQUS/CAE 不会自动保存模型数据，用户应当每隔一段时间自己保存模型以避免意外丢失。3(pp17)平面应力问题的截面属性类型是Solid（实心体）而不是Shell（壳）。ABAQUS/CAE

26、 推荐的建模方法是把整个数值模型（如材料、边界条件、载荷等）都直接定义在几何模型上。载荷类型Pressure 的含义是单位面积上的力，正值表示压力，负值表示拉力。4(pp22)对于应力集中问题，使用二次单元可以提高应力结果的精度。5(pp23)Dismiss 和Cancel 按钮的作用都是关闭当前对话框，其区别在于：前者出现在包含只读数据的对话框中；后者出现在允许作出修改的对话框中，点击Cancel 按钮可关闭对话框，而不保存所修改的内容。6(pp26)每个模型中只能有一个装配件，它是由一个或多个实体组成的，所谓的“实体”（instance）是部件（part）在装配件中的一种映射，一个部件可以

27、对应多个实体。材料和截面属性定义在部件上，相互作用（interaction）、边界条件、载荷等定义在实体上，网格可以定义在部件上或实体上，对求解过程和输出结果的控制参数定义在整个模型上。7(pp26) ABAQUS/CAE 中的部件有两种：几何部件（native part）和网格部件（orphan mesh part）。创建几何部件有两种方法：（1）使用Part 功能模块中的拉伸、旋转、扫掠、倒角和放样等特征来直接创建几何部件。（2）导入已有的CAD 模型文件，方法是：点击主菜单FileImportPart。网格部件不包含特征，只包含节点、单元、 面、集合的信息。创建网格部件有三种方法：（1）

28、导入ODB 文件中的网格。（2）导入INP 文件中的网格。（3）把几何部件转化为网格部件，方法是：进入Mesh 功能模块，点击主菜单MeshCreate Mesh Part。8(pp31)初始分析步只有一个，名称是initial，它不能被编辑、重命名、替换、复制或删除。在初始分析步之后，需要创建一个或多个后续分析步，主要有两大类：（1）通用分析步（general analysisstep）可以用于线性或非线性分析。常用的通用分析步包含以下类型：Static, General: ABAQUS/Standard 静力分析Dynamics, Implicit: ABAQUS/Standard 隐式动

29、力分析Dynamics, Explicit: ABAQUS/ Explicit 显式动态分析（2）线性摄动分析步（linear perturbation step）只能用来分析线性问题。在ABAQUS/Explicit 中不能使用线性摄动分析步。在ABAQUS/Standard 中以下分析类型总是采用线性摄动分析步。Buckle: 线性特征值屈曲。Frequency: 频率提取分析。Modal dynamics: 瞬时模态动态分析。Random response: 随机响应分析。Response spectrum: 反应谱分析。Steady-state dynamics: 稳态动态分析。9（p

30、p33）在静态分析中，如果模型中不含阻尼或与速率相关的材料性质，“时间”就没有实际的物理意义。为方便起见，一般都把分析步时间设为默认的1。每创建一个分析步，ABAQUS/CAE 就会自动生成一个该分析步的输出要求。10 （pp34）自适应网格主要用于ABAQUS/Explicit 以及ABAQUS/Standard 中的表面磨损过程模拟。在一般的ABAQUS/Standard 分析中，尽管也可设定自适应网格，但不会起到明显的作用。Step 功能模块中，主菜单OtherAdaptive Mesh Domain 和OtherAdaptive Mesh Controls 分别设置划分区域和参数。11（pp37）使用主菜单Field 可以定义场变量（包括初始速度场和温度场变量）。有些场变量与分析步有关，也有些仅仅作用于分析的开始阶段。使用主菜单Load Case 可以定义载荷状况。载荷状况由一系列的载荷和边界条件组成，用于静力摄动分析和稳态动力分析。12（pp42）独立实体是对部件的复制，可以直接对独立实体划分网格，而不能对相应的部件划分网格。非独立实体是部件的指针，不能直接对非独立实体划分网格，而只能对相应的部件划分网格。由网格部件创建的实体都是非独立实体。13（pp45）Quad 单元（二维区域内完全使用四边形网格）和Hex 单元（三维区域内完全使用六面体网格）可以用较小的