VB环境下ANSYS二次开发：超弹性材料模型算法

ANSYS现有版本虽包括了金属、橡胶、Drucker-Prager、混凝土等众多材料模型，但仍无法满足工程计算需要，为了弥补这一不足，ANSYS为用户提供了开发材料模型的接口，即UPFs。通过修改、编译连接相关用户子程序，可以得到各种符合用户需要的材料模型。ANSYS的TB，HYPER命令给用户提供了各种不可压缩和可压缩的超弹性材料模型，比如：Polynomid Form模型、Mooney-Rivlin模型、Neo-Hookean模型、Yeoh模型、Arruda-Boyce模型、Gent模型、Ogden模型、Hyperfoam模型以及Blatz-Ko模型等。但是对于需要使用另外模型的用户，则需要UserHyper用户子程序来编写自己的超弹性材料模型。

UserHyper用户子程序介绍

用户可以使用如下命令调用用户定义的超弹性材料模型：

TB,HYPER,,,,USER

可以使用所有支持超弹性材料的单元。

UserHyper用户子程序输入输出变量说明如下：

*deck.UserHyper USERDISTRIB parallel gal

subroutine UscrHyper(

＆ prophy,irtcomp,nprophy,invar,

＆ potential,plnvDer)

c*********************************************************

c 输入变量

c =================

c prophy (dp,ar(*),i) 材料特性数组

c npoophy (int,sc,i) 材料常数个数

c ivar dp,ar(3) 不变量

c 输出变量

c =================

c incomp (log,sc,i) 可压缩或不可压缩

c potential dp,sc 应变能

c pInvDer dp,ar(10) 应变能wtr对il,i2,j的导数

c l-应变能wtr对il的导数

c 2-应变能wtr对i2的导数

c 3-应变能wtr对il的二阶导数

c 4-应变能wtr对il和i2的导数

c 5-应变能wtr对i2的导数

c 6-应变能wtr对il和j的导数

c 7-应变能wtr对i2和j的导数

c 8-应变能wtr对j的导数

c 9-应变能wtr对j的二阶导数

模型算法

超弹性理论中一个很重要的量是应变能函数，UserHyper用户子程序的应变能函数适用于基于三个应变不变量的率无关且各向同性超弹性模型。如果用户想要定义率相关、非各向同性或非弹性材料模型，则需要利用前面介绍的usermat用户子程序。ANSYS安装文件夹中的UserHyper.F子程序基于Arruda-Boyce模型。这里将要介绍的是一种ANSYS没有的超弹性材料模型。

在超弹性理论中，应变能函数被划分为偏部分(deviatoric part)和体积部分(volumetric part)，因此一般情况下应变能函数的输入变量是3个修正的应变不变量(分别对应UserHyper用户子程序的输入变量invar(l)、invar(2)、invar(3））：

数值实施

在超弹性理论中，可以将应变能函数分成两部分处理，即偏应变能部分和体积应变能部分。本节使用的模型偏应变能函数采用Horgan-Saccomandi建议的形式：

UserHyper用户子程序的详细代码如下：

*deck.UserHyper USERDISTRIB parallel gal

subroutine UscrHyper(

＆ prophy,irtcomp,nprophy,invar,

＆ potential,plnvDer)

c*********************************************************

c 输入变量

c =================

c prophy (dp,ar(*),i) 材料特性数组

c npoophy (int,sc,i) 材料常数个数

c ivar dp,ar(3) 不变量

c 输出变量

c =================

c incomp (log,sc,i) 可压缩或不可压缩

c potential dp,sc 应变能

c pInvDer dp,ar(10) 应变能wtr对il,i2,j的导数

c l-应变能wtr对il的导数

c 2-应变能wtr对i2的导数

c 3-应变能wtr对il的二阶导数

c 4-应变能wtr对il和i2的导数

c 5-应变能wtr对i2的导数

c 6-应变能wtr对il和j的导数

c 7-应变能wtr对i2和j的导数

c 8-应变能wtr对j的导数

c 9-应变能wtr对j的二阶导数

c**************************************************************

c---参数

#include “impcom.inc”

DOUBLE PRECISION mu,jm,alpha,bulk

DOUBLE PRECISION il,i2,i3

DOUBLE PRECISION jdcnom, vpotcn

c ----变量列

INTEGER nprophy

DOUBLE PRECISION prophy(*),invar(*),

& potential,pInvDer(*)

LOGICAL incomp

c----局部变量

mu =prophy(l)/2.0d0

jm =prophy(2)

bulk =prophy(3)

alpha=prophy(4)

i1 = invar<l)

i2 = invar(2)

i3 = invar(3)

c-----计算偏应变能

jdenom=(jm*jm**2.0d0-jm**2.0d0*il+jm*i2-1.OdO)

potential =-mu*jm*log(jdenom /((jm-1.OdO)* !偏应变能Wd

&(jm - 1.0d0)**2.0d0))

plnvDer(1)= mu*jm*jm**2.0d0/jdenom !计算导数值

pInvDer(2)=-mu*jm**2.0d0/jdenom

pInvDer(3)= mu*jm*jm**2.0d0*jm**2.0d0/jdenom**2.0d0

pInvDer(4)=-mu*jm**2.0d0*jm**2.0d0/jdenom**2.0d0

pInvDcr(5)= mu*jm*jm**2.0d0/jdenom**2.0d0

c----计算体积应变能

pInvDer(6)= O.OdO

pInvDer(7)= O.OdO

pInvDer(8)= O.OdO

pInvDer(9)= O.OdO

IF(alpha.gt.0.0001d0)THEN !通过参数a判断材料是否可压缩

incomp=.FALSE. !可压缩

vpoten =bulk/alpha**2.0d0*(cosh(alpha*(i3-1.OdO))

!体积皮变能 Wv

& -l.OdO)

potential=potential+vpoten !总应变能

plnvDcr(8)= bulk/alpha**1.OdO*sinh(alpha*(i3-l.OdO))

!计算导数值

pInvDer(9)= bulk*cosh(alpha*(i3-l.OdO))

END IF

RETURN

END

生成并调用宏文件

在ANSYS中，宏是包含一系列ansys命令并且后缀为.MAC或.mac的命令文件。宏文件往往记录一系列频繁使用的ansys命令流，实现某种有限元分析或其他算法功能。宏文件在ansys中可以当作定义的ansys命令进行使用，可以带有宏输入参数，也可以有内部变量，同时在宏内部可以直接引用总体变量。除了执行一系列的ansys命令之外，宏还可以调用GUI函数或把值传递给参数。

利用*USE命令调用宏文件，并向宏文件传递参数：

*USE,Name,ARG1,ARG2,ARG3,ARG4,ARG5,ARG6,ARG7,ARG8,ARG9,AR10,AR11,AR12,AR13,AR14,AR15,AR16,AR17,AR18

其中，Name是宏文件名，ARGI到AR18是宏文件用到的参数值。

APDL实现过程

下面为两个简单的橡胶类材料受力分析的实例，目的是与ANSYS自带的Gent模型比较，以便验证前面建立的用户超弹性模型的正确性。通过模拟单轴拉伸试验考察Horgan-Saccomandi偏应变能函数，通过模拟静水压缩考察Bischoff体积应变能函数。

1．单轴拉伸

建立两个SOLID185单元，边界条件完全相同，只是使用的材料不同，如图8-7所示。命令流( \chp8\userhyper\userhyper_uniaxial.inp)如下：

finisb

/clear

！设置参数

STRE1TCH=4.9 !拉伸量参数

！进入前处理嚣

/prep7

et,l，185 ！设置单元类型

r,l

tb,hyper,1,1,4,user ！l号材料为用户材料

tbdata,,100,25,0,0 !a=0，所以是不可压缩材料

tb,hyper,2,1,3,gent !2号材料为ANSYS标准gent

tbdata,,lOO,25.4-3,0 !gent超弹材料参数

！几何建模

block，，1，，1，，1

block,,2,3,,1,,1

esize,,1 !仅划分一个单元

mat,1$vmesh,1 !划分网格

mat,2$vmesh,2

！设置边界条件

nsel,s,loc,x,0

nsel,s,loc,x,2

d,all,ux

nsel,s,loc,y,0

d,all,uy

nsel,s,loc,z,0

d,all,uz

nsel,s,loc,y,1

nsel,r,loc,x,0,1

cp,next,uy,all

MYPARM1-ndnext(0)

d*ndnext(0),uy,STRETCH-l !施加拉伸约束

nsel,s,loc,y,1

nsel,r,loc,x,2,3

cp,next,uy,all

d*ndnext(0),uy,STRETCH-l ！施加拉伸约束

allsel,all

finish

!进入求解器

/solu

antype,static ！静态分析类型

nlgeom,on ！打开大变形开关

outres,all,all ！输出所有计算结果

nsubst,100,1000,50 ！求解子步数设置

solve

save

finish

！1.定义工作文件名及文件标题

!(1)定义工作文件名：

/FILNAME,shiyan_thermal,0

!(2)定义文件标题

/TITLE,shiyan of experiment

!2.定义单元类型及材料属性

！（1）定义单元类型：

/PREP7

ET,1,185

！（2）设置材料属性

MPTEMP,,,,,,,,

MPTEMP,1,0

MPDATA,EX,1,,2e11

MPDATA,PRXY,1,,0.3

MPTEMP,,,,,,,,

MPTEMP,1,0

MPDATA,DENS,1,,7800

!3.建立几何模型

！（1）显示工作平面

WPSYYLE,,,,,,,,1

!(2)建立模型

BLOCK,0,10,0,20,0,2,

/REPLO

wpoff,3,5,0

/REPLO

CYLIND,1, ,0,10,0,360,

/REPLO

wpoff,0,10,0

CYLIND,1, ,0,10,0,360,

wpoff,0,-5,0

CYLIND,1, ,0,10,0,360,

wpoff,4,0,0

CYLIND,1, ,0,10,0,360,

wpoff,0,5,0

CYLIND,1, ,0,10,0,360,

wpoff,0,-10,0

CYLIND,1, ,0,10,0,360,

！（3）体布尔操作

FLST,3,6,6,ORDE,2

FITEM,3,2

FITEM,3,-7

VSBV, 1,P51X

！4.生成有限元模型

/REPLO

FINISH

/AUX12

FINISH

/PREP7

MSHAPE,1,3D

MSHKEY,0

CM,_Y,VOLU

VSEL, , , , 8

CM,_Y1,VOLU

CHKMSH,'VOLU'

CMSEL,S,_Y

VMESH,_Y1

CMDELE,_Y

CMDELE,_Y1

CMDELE,_Y2

！5.加载以及求解

！（1）加载

//REPLO

VPLOT

FINISH

/SOL

FLST,2,1,5,ORDE,1

FITEM,2,3

/GO

DA,P51X,ALL,0

/REPLO

FLST,2,1,4,ORDE,1

FITEM,2,7

/GO

FLST,2,1,5,ORDE,1

FITEM,2,4

/GO

SFA,P51X,1,PRES,1000

/STATUS,SOLU

SOLVE

FINISH

/POST1

PLDISP,2

/EFACET,1

PLNSOL, S,EQV, 0,1.0

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

返回上级列表

请

联系我们

下一篇: ANSYS DesignModeler 3D曲线特征创建：点文件方式 >