Marc中用户子程序UEALASTOMER的使用 1 概述 Marc在2005r1版增加了通用的橡胶类材料本构模型定义用户子程序UEALASTOMER,通过该用户子程序,用户很容易将自己研究得到的本构模型耦合到Marc软件之中,或者将自己的一些材料数据编入到Marc程序中。下面就该用户子程序的使用方法做一个简单介绍。 2用户子程序介绍 该用户子程序主要是让用户定义自己的超弹性本构模型。可以通过模型定义选项:FOAM, MOONEY, ARRUDBOYCE, GENT和OGDEN来激活。该用户子程序必须要和大应变选项LARGER STRAIN, 2 一起使用。对于一些不支持该大应变选项的单元类型如rebar单元,是不能采用该用户子程序的。 对于可压缩泡沫材料,有四种类型的应变能函数可以用UEALASTOMER用户子程序定义。具体的类型由模型定义选项中FOAM中的第三行第4个数据来确定: 1)iflag=1,基于不变量模型; 2)iflag=2,基于主伸长模型; 3)iflag=3,基于不变量模型并分为体积量和偏量部分; 4)iflag=4,基于主伸长模型并分为体积量和偏量部分; 对于泡沫模型( iflag=1,2,3,4),Marc采用传统的位移单元,不考虑体积约束。对于不可压缩和近似不可压缩橡胶类材料,橡胶模型可以采用iflag=5或6。如果采用了FEATURE,3402 选项,体积应变能也可以通过iflag=7来定义。 因为橡胶类材料都是近似不可压缩的,因而将应变能函数分为体积量部分和偏量部分在数值上效率更高一些。一种混合表达式,将静水压作为独立变量,在Marc中用于克服由于体积约束带来的数值困难。在混合模式中,采用压力和体积应变成线性的假设。因而,只有应变能函数中的偏量部分才需要在用户子程序中定义。 然而,有时压力和体积应变之间不是线性关系。此时,就有了由压力、体积率和位移构成的三场表达式。目前FEATURE,3402参数选项用于激活三场表达式。此时,应变能函数中的体积项和偏量项都可以在用户子程序UEALASTOMER中定义。 基于不变量的橡胶模型通过模型数据MOONEY, ARRUDBOYCE或GENT模型定义选项在用户子程序UEALASTOMER中定义。 5)iflag=5,基于不变量模型,仅有偏量部分; 基于主伸长的橡胶模型通过在模型数据OGDEN的第三行第3个数据设为3后在用户子程序UEALASTOMER中定义。 6)iflag=6,基于主伸长模型,仅有偏量部分; 值得注意的是,当iflag=5或6时只有偏量部分在UEALASTOMER中定义,体积项在Marc内计算。为此,必须在模型定义选项中定义体积模量。 对于不使用用户子程序时及材料数据直接输入时,如果输入数据中没有定义体积模量,体积模量由程序自动设为初始剪切模量的5000倍。但如果采用了用户子程序UEALASTOMER,不能确定初始的剪切模量,因而体积模量必须直接定义。 7)iflag=7,应变能的体积项,基于主伸长和不变量的模型均可以用。 橡胶类模型可以用于体积可压缩和体积不可压缩材料。 对于用户自己研发的材料本构模型,在编用户子程序前,一般要先推导出相关的应变能函数的导数和偏导数,然后进行用户子程序的编程。参见下表: file:///C:\Users\user\AppData\Local\Temp\msohtmlclip1\01\clip_image002.giffile:///C:\Users\user\AppData\Local\Temp\msohtmlclip1\01\clip_image004.gif
例如,对于很多用户比较熟悉的Mooney模型:W=C1*(I1-3)+C2*(I2-3),可推导出W1 =C1,W2 =C2。 3 应用实例 下面举一个O型圈受压的例子来说明该用户子程序的用法,该例子选自Marc用户手册E卷7.23节。材料实际为OGDEN模型,参数如下: file:///C:\Users\user\AppData\Local\Temp\msohtmlclip1\01\clip_image006.gif 看起来是两项式,实际上第一项的系数µ为0,因此实际是一项式,但编程时时按照两项式处理。 iflag=2,因此是基于主伸长的泡沫本构模型,因而需要在用户子程序中定义W1、W2、W3、W11、W22、W33、W12、W23、W31。得到所需9项的表达式即可编用户子程序,在Marc的安装目录下有各子程序的模板文件,拷贝到工作目录中进行修改,如下所示: subroutine uelastomer(iflag,m,nn,mats,be,x1,x2,x3,detft, $ enerd,w1,w2,w3,w11,w22,w33,w12,w23,w31, $ dudj,du2dj,dt,dtdl,iarray,array) c implicit logical (a-z) c ** Start of generated type statements ** real*8 aaa, alfa, array, bbb, be, beta, ccc, cmu, det, detbeta real*8 detft, dt, dtdl, du2dj, dudj, enerd, enerdn integer i, iarray, iflag, m, mats, nn, nser real*8 pa_alfa, pa_beta, pa_cmu, w1, w11, w12, w2, w22, w23, w3 real*8 w31, w33, x1, x1alfa, x2, x2alfa, x3, x3alfa c ** End of generated type statements ** dimension m(2),be(6),dt(*),dtdl(*) c c iflag=2, 类型为2 c c defined material parameters,定义模型参数 c dimension pa_cmu(10),pa_alfa(10),pa_beta(10) nser=2 pa_cmu(1)=0.d0 pa_alfa(1)=2.d0 pa_beta(1)=-1.d0 pa_cmu(2)=-3.2d1 pa_alfa(2)=-2.d0 pa_beta(2)=1.d0 c det=x1*x2*x3 c enerd=0.0d0 w1=0.0d0 w2=0.0d0 w3=0.0d0 w11=0.0d0 w22=0.0d0 w33=0.0d0 w12=0.0d0 w23=0.0d0 w31=0.0d0 do i=1,nser c c extract foam parameters,提取材料参数 c cmu=pa_cmu(i) alfa=pa_alfa(i) beta=pa_beta(i) c c calculate stress/tangent coefficients,计算应力/切线系数 c x1alfa=x1**alfa x2alfa=x2**alfa x3alfa=x3**alfa detbeta=det**beta c if(abs(alfa).gt.1.d-20.and.abs(beta).gt.1.d-20) $ enerdn=(x1alfa+x2alfa+x3alfa-3.d0)/alfa $ +(1.d0-detbeta)/beta enerd=enerd+enerdn*cmu c w1=w1+cmu*(x1alfa-detbeta)/x1 w2=w2+cmu*(x2alfa-detbeta)/x2 w3=w3+cmu*(x3alfa-detbeta)/x3 c aaa=alfa-1.d0 bbb=(beta-1.d0)*detbeta w11=w11+cmu*(aaa*x1alfa-bbb)/x1/x1 w22=w22+cmu*(aaa*x2alfa-bbb)/x2/x2 w33=w33+cmu*(aaa*x3alfa-bbb)/x3/x3 c ccc=cmu*beta w12=w12-ccc*detbeta/x1/x2 w23=w23-ccc*detbeta/x2/x3 w31=w31-ccc*detbeta/x3/x1 c enddo return end 定义材料特性时选用下图所示的选项: file:///C:\Users\user\AppData\Local\Temp\msohtmlclip1\01\clip_image008.jpg file:///C:\Users\user\AppData\Local\Temp\msohtmlclip1\01\clip_image010.jpg 图 1 经典界面和新界面材料特性定义菜单 在运行作业的窗口中指定用户子程序文件,如下图所示,其它操作与直接输入材料参数是相同。 file:///C:\Users\user\AppData\Local\Temp\msohtmlclip1\01\clip_image011.giffile:///C:\Users\user\AppData\Local\Temp\msohtmlclip1\01\clip_image012.giffile:///C:\Users\user\AppData\Local\Temp\msohtmlclip1\01\clip_image014.jpg file:///C:\Users\user\AppData\Local\Temp\msohtmlclip1\01\clip_image016.jpg
图 2 经典界面和新界面运行选用用户子程序文件菜单
以下是初始模型图以及分析得到的O型圈在受压过程中位移云图和最后的接触力分布图。
图3 变形前模型图
图4 对压量50%时的位移云图 图5 对压量100%时的位移云图 图6 对压量100%时的接触力分布云图
4 参考信息 已经安装Marc2013.1的用户,可以在Marc的安装路径以下位置找到相关文件:D:\MSC.Software\Marc\2013.1.0\marc2013.1\demo 适用版本:Marc 2005及Marc2005以后版本
|