请问对圆柱壳施加轴压载荷，无论加载有多大，提取的载荷位移曲线总是线性的呢，也考虑了材料的退化，刚度不是应该变化吗，谢谢
下面是所采用的命令流
*creat,reaction,mac !生成reaction宏文件，将如下命令加入，以便input
*cfopen,reaction,txt !在默认路径下生产reaction.txt文件
!*cfopen,d:\file1,txt
*vwrite
(' 位移载荷大小/mm 支座反力/N')
*vwrite,DDX(1),FFZ(1)
(2f15.2)
*cfclos
*end
Finish
/clear
!/title, pprogressive failure analysis of carbon fiber/epoxy composite laminates
!using CDM(40MPa hydrogen storage tank)
/uis,msgpop,3
/prep7
et,1,95
et,2,64
mp,ex,1,70e3
mp,nuxy,1,0.3
tb,biso,1,1,2
tbdata,,246,600
*afun,deg
E1=181e3
E2=10.3e3
v12=0.28
v23=0.49
!!!定义约束反力数组
*dim,FFZ,array,100
!!!定义临时位移数组，后续定义DDX的实际长度即为数据长度
*dim,DDX,array,100
*dim,ls,array,1000
*dim,offsy,array,5000
*dim,offsz,array,5000
*dim,offsyz,array,5000
*dim,offstrainz,array,5000
*dim,offstrainsita,array,5000
*dim,offstrainzsita,array,5000
*dim,ons1,array,5000
*dim,ons2,array,5000
*dim,ons6,array,5000
*dim,onstrain1,array,5000
*dim,onstrain2,array,5000
*dim,onstrain6,array,5000
*dim,H1,array,5000
*dim,H2,array,5000
*dim,H6,array,5000
*dim,H12,array,5000
*dim,Tsaiwu,array,5000
*dim,Y1,array,5000
*dim,Y2,array,5000
*dim,Y6,array,5000
F22=1.0/(298*298)
F66=1.0/(778*778)
F11=1.0/(2150*2150)
F12= -0.5*sqrt(F11*F22)
LFN=0
FFN=0
MCN=0
SSN=0
*dim,MLFN,array,100
*dim,MFFN,array,100
*dim,MMCN,array,100
*dim,MSSN,array,100
*dim,D1,array,5000
*dim,D2,array,5000
*dim,D6,array,5000
*dim,angle,array,10
*dim,m,array,10
*dim,n,array,10
*dim,te,array,6,6,5000
*dim,ts,array,6,6,5000
*dim,tst,array,6,6,5000
*dim,unit,array,6,6,5000
*dim,s,array,6,6,5000
*dim,c,array,6,6,5000
*dim,cp,array,6,6,5000
*dim,midlpf,array,6,6,5000
angle(1)=90
angle(2)=-90
angle(3)=18.9
angle(4)=-18.9
angle(5)=90
angle(6)=-90
angle(7)=28.6
angle(8)=-28.6
angle(9)=90
angle(10)=-90
*do,i,1,10
*vfun,m(i),cos,angle(i)
*vfun,n(i),sin,angle(i)
*enddo
*do,i,1,10 ! 转轴矩阵
*do,j,(i-1)*500+1,i*500
te(1,1,j)=m(i)*m(i),n(i)*n(i),0,0,0,-2*m(i)*n(i)
te(1,2,j)=n(i)*n(i),m(i)*m(i),0,0,0,2*m(i)*n(i)
te(1,3,j)=0,0,1,0,0,0
te(1,4,j)=0,0,0,m(i),n(i),0
te(1,5,j)=0,0,0,-n(i),m(i),0
te(1,6,j)=m(i)*n(i),-m(i)*n(i),0,0,0,m(i)*m(i)-n(i)*n(i)
*enddo
*enddo
*do,i,1,5000
unit(1,1,i)=1,0,0,0,0,0
unit(1,2,i)=0,1,0,0,0,0
unit(1,3,i)=0,0,1,0,0,0
unit(1,4,i)=0,0,0,1,0,0
unit(1,5,i)=0,0,0,0,1,0
unit(1,6,i)=0,0,0,0,0,1
*enddo
*do,i,1,5000
*moper,tst(1,1,i),te(1,1,i),solv,unit(1,1,i) ! 转轴矩阵
*enddo
*do,i,1,5000
*mfun,ts(1,1,i),tran,tst(1,1,i) ! 矩阵转置
*enddo
*do,i,1,5000 ! 编写刚度矩阵
s(1,1,i)=1.0/E1,-v12/E1,-v12/E1,0,0,0
s(1,2,i)=-v12/E1,1.0/E2,-v23/E1,0,0,0
s(1,3,i)=-v12/E1,-v23/E1,1.0/E2,0,0,0
s(1,4,i)=0,0,0,2*(1+v23)/E2,0,0
s(1,5,i)=0,0,0,0,2*(1+v12)/E1,0
s(1,6,i)=0,0,0,0,0,2*(1+v12)/E1
*enddo
*do,i,1,5000
*moper,c(1,1,i),s(1,1,i),solv,unit(1,1,i)
*enddo
*create, material
*do,j,1,5000 ! 引入破坏变量
c(1,1,j)=(1-D1(j))*c(1,1,j)
c(2,2,j)=(1-D2(j))*c(2,2,j)
c(6,6,j)=(1-D6(j))*c(6,6,j)
*enddo
*do,j,1,5000
*moper,midlpf(1,1,j),ts(1,1,j),mult,c(1,1,j) ! 矩阵相乘
*enddo
*do,j,1,5000
*moper,cp(1,1,j),midlpf(1,1,j),mult,tst(1,1,j)
*enddo
*do,j,1,5000
TB,ANEL,j+1,1,21,0
TBDATA,,cp(1,1,j),cp(1,2,j),cp(1,3,j),0,0,cp(1,6,j)
TBDATA,,cp(2,2,j),cp(2,3,j),0,0,cp(2,6,j),cp(3,3,j)
TBDATA,,0,0,cp(3,6,j),cp(4,4,j),cp(4,5,j),0
TBDATA,,cp(5,5,j),0,cp(6,6,j)
*enddo
*end
*use,material
csys,1
!cyl4,0,0,44,0,45.8,30,160
*do,i,1,10
cyl4,0,0,45.8+(i-1)*0.42,0,45.8+i*0.42,30,160
*enddo
vglue,all
numcmp,all
!LESIZE,12,,,50
!LESIZE,1,,,2
!LESIZE,3,,,2
!LESIZE,6,,,2
!LESIZE,8,,,2
!LESIZE,2,,,10
!LESIZE,4,,,10
!LESIZE,5,,,10
!LESIZE,7,,,10
LSEL,S,LENGTH,, 0.42 !选择所有长度为T_Ply的线段
LESIZE,all, , , 1, , , , , !分成1份
LSEL,S,LENGTH,, 160 !选择所有长度为T_Ply的线段
LESIZE,all, , , 50, , , , , !分成1份
allsel
LESIZE,2,,,10
LESIZE,4,,,10
LESIZE,5,,,10
LESIZE,7,,,10
*do,i,1,9
LESIZE,13+(i-1)*4,,,10
*enddo
*do,i,1,9
LESIZE,14+(i-1)*4,,,10
*enddo
alls
!type,1
!mat,1
!vmesh,1
type,2
*do,i,1,10
mat,i+1
vmesh,i
*enddo
nummrg,node,,,,low
nummrg,elem,,,,low
*do,i,1,5000
emodif,i,mat,i+1
*enddo
/pnum,mat,1
eplot
Finish
/solu
allsel
da,18,symm
da,6,symm
*do,i,1,8
da,(i-1)*4+22,symm
*enddo
da,5,symm
*do,i,1,9
da,(i-1)*4+19,symm
*enddo
allsel
asel,s,loc,z,0 !施加z约束
da,all,uz,0
allsel
!d,1000,uz,0
rescontrol,define,all,1,1
sfa,37,1,pres,-500
DDX(1)=-500
allsel
rescontrol,define,all,1,1
nsubst,5
cnvtol,f,,0.1,,0.1
solve
finish
!!!提取支座反力：
/post1
nsel,s,loc,z,0
nplot
fsum
*get,FFZ(1),fsum,0,item,FZ
Allsel
eplot
*do,i,2,50
/post1
rsys,1
!etable,sequ,s,eqv
etable,sy,s,y
etable,sz,s,z
etable,syz,s,yz
etable,strainy,epel,y
etable,strainz,epel,z
etable,strainyz,epel,yz
! *do,j,1,1000
! *get,ls(j),elem,j,etab,sequ
!*enddo
*do,j,1,5000
*get,offsy(j),elem,j,etab,sy
*get,offsz(j),elem,j,etab,sz
*get,offsyz(j),elem,j,etab,syz
*enddo
*do,j,1,5000
*get,offstrainsita(j),elem,j,etab,strainy
*get,offstrainz(j),elem,j,etab,strainz
*get,offstrainzsita(j),elem,j,etab,strainyz
*enddo
*do,j,1,10
*do,k,(j-1)*500+1,j*500
ons1(k)=m(j)*m(j)*offsz(k)+n(j)*n(j)*offsy(k)-2*m(j)*n(j)*offsyz(k)
ons2(k)=n(j)*n(j)*offsz(k)+m(j)*m(j)*offsy(k)+2*m(j)*n(j)*offsyz(k)
ons6(k)=m(j)*n(j)*offsz(k)-m(j)*n(j)*offsy(k)+(m(j)*m(j)-n(j)*n(j))*offsyz(k)
*enddo
*enddo
*do,j,1,10
*do,k,(j-1)*500+1,j*500
onstrain1(k)=m(j)*m(j)*offstrainz(k)+n(j)*n(j)*offstrainsita(k)-m(j)*n(j)*offstrainzsita(k)
onstrain2(k)=n(j)*n(j)*offstrainz(k)+m(j)*m(j)*offstrainsita(k)+m(j)*n(j)*offstrainzsita(k)
onstrain6(k)=2*m(j)*n(j)*(offstrainz(k)-offstrainsita(k))+(m(j)*m(j)-n(j)*n(j))*offstrainzsita(k)
*enddo
*enddo
*do,j,1,5000
Y1(j)=abs(0.5*ons1(j)*onstrain1(j))
Y2(j)=abs(0.5*ons2(j)*onstrain2(j))
Y6(j)=abs(0.5*ons6(j)*onstrain6(j))
*enddo
*do,j,1,5000
H1(j)=F11*ons1(j)*ons1(j)
H2(j)=F22*ons2(j)*ons2(j)
H6(j)=F66*ons6(j)*ons6(j)
H12(j)=2*F12*ons1(j)*ons2(j)
Tsaiwu(j)=H1(j)+H2(j)+H6(j)+H12(j)
*enddo
parsav,all,parameter,txt
/solu
antype,,rest
parres,new,parameter,txt
!*do,j,1,1000
! *if,ls(j),gt,324,then
! ekill,j
!LFN=LFN+1
!*endif
! *enddo
*do,j,1,5000
*if,Tsaiwu(j),gt,1,then
*if,H1(j),gt,H2(j),and,H1(j),gt,H6(j),then
D1(j)=1-exp(-0.01*Y1(j))
FFN=FFN+1
*elseif,H2(j),gt,H1(j),and,H2(j),gt,H6(j),then
D2(j)=Y2(j)/40
MCN=MCN+1
*else
D6(j)=Y6(j)/15
SSN=SSN+1
*endif
*endif
*enddo
!MLFN(i-1)=LFN
MFFN(i-1)=FFN
MMCN(i-1)=MCN
MSSN(i-1)=SSN
/prep7
*do,j,1,5000
TBDE,ANEL,j+1
*enddo
*use,material
/solu
allsel,all
sfa,37,1,pres,-( 500+i*200)
DDX(i)= -(500+i*200)
allsel
nsubst,5
cnvtol,f,,0.1,,0.1
ncnv,2
solve
!!!提取约束反力
/post1
SET,LAST
nsel,s,loc,z,0
nplot
fsum
*get,FFZ(i),fsum,0,item,FZ
allsel
eplot
*enddo
reaction
!!!输出最后一个载荷步的图像各层图片
i=i+1
/post1
rsys,1
!etable,sequ,s,eqv
etable,sy,s,y
etable,sz,s,z
etable,syz,s,yz
etable,strainy,epel,y
etable,strainz,epel,z
etable,strainyz,epel,yz
! *do,j,1,1000
! *get,ls(j),elem,j,etab,sequ
!*enddo
*do,j,1,5000
*get,offsy(j),elem,j,etab,sy
*get,offsz(j),elem,j,etab,sz
*get,offsyz(j),elem,j,etab,syz
*enddo
*do,j,1,5000
*get,offstrainsita(j),elem,j,etab,strainy
*get,offstrainz(j),elem,j,etab,strainz
*get,offstrainzsita(j),elem,j,etab,strainyz
*enddo
*do,j,1,10
*do,k,(j-1)*500+1,j*500
ons1(k)=m(j)*m(j)*offsz(k)+n(j)*n(j)*offsy(k)-2*m(j)*n(j)*offsyz(k)
ons2(k)=n(j)*n(j)*offsz(k)+m(j)*m(j)*offsy(k)+2*m(j)*n(j)*offsyz(k)
ons6(k)=m(j)*n(j)*offsz(k)-m(j)*n(j)*offsy(k)+(m(j)*m(j)-n(j)*n(j))*offsyz(k)
*enddo
*enddo
*do,j,1,10
*do,k,(j-1)*500+1,j*500
onstrain1(k)=m(j)*m(j)*offstrainz(k)+n(j)*n(j)*offstrainsita(k)-m(j)*n(j)*offstrainzsita(k)
onstrain2(k)=n(j)*n(j)*offstrainz(k)+m(j)*m(j)*offstrainsita(k)+m(j)*n(j)*offstrainzsita(k)
onstrain6(k)=2*m(j)*n(j)*(offstrainz(k)-offstrainsita(k))+(m(j)*m(j)-n(j)*n(j))*offstrainzsita(k)
*enddo
*enddo
*do,j,1,5000
Y1(j)=abs(0.5*ons1(j)*onstrain1(j))
Y2(j)=abs(0.5*ons2(j)*onstrain2(j))
Y6(j)=abs(0.5*ons6(j)*onstrain6(j))
*enddo
*do,j,1,5000
H1(j)=F11*ons1(j)*ons1(j)
H2(j)=F22*ons2(j)*ons2(j)
H6(j)=F66*ons6(j)*ons6(j)
H12(j)=2*F12*ons1(j)*ons2(j)
Tsaiwu(j)=H1(j)+H2(j)+H6(j)+H12(j)
*enddo
|
IP卡
狗仔卡
发表于 2017-2-20 20:45:06
提升卡
沉默卡
喧嚣卡
变色卡
抢沙发
千斤顶
显身卡
