三维裂纹

*create,fract,mac
!******创建文件名为fract的宏文件
/nopr
nsel,all
!******选择所有的节点
*get,n,node,,num,max
!******将节点的最大序号赋予n
cmsel,s,cracktip
!******选择一个名字为cracktip的组分
esln
!******选择连接在已选节点上的单元
*get,elmax,elem,,num,max
!******将单元的最大序号赋予elmax
*do,iel,1,elmax
!******循环开始，从1到elmax
elmi=iel
*if,elmi,le,0,exit
!******??????序号不可能小于0，这两语句是否多余?
*get,eltype,elem,elmi,attr,type
!******将第elmi个单元选择情况
*if,eltype,ne,arg1,cycle
!******当eltype不等于arg1时，跳出循环
n3=nelem(elmi,3)
*if,nsel(n3),le,0,cycle
n7=nelem(elmi,7)
*if,nsel(n7),le,0,cycle
n1=nelem(elmi,1)
n2=nelem(elmi,2)
n5=nelem(elmi,5)
n6=nelem(elmi,6)
x3=0.75*nx(n3)
y3=0.75*ny(n3)
z3=0.75*nz(n3)
x=0.25*nx(n2)+x3
y=0.25*ny(n2)+y3
z=0.25*nz(n2)+z3
n=n+1
n10=n
n,n10,x,y,z
x=0.25*nx(n1)+x3
y=0.25*ny(n1)+y3
z=0.25*nz(n1)+z3
n=n+1
n12=n
n,n12,x,y,z
x7=0.75*nx(n7)
y7=0.75*ny(n7)
z7=0.75*nz(n7)
x=0.25*nx(n6)+x7
y=0.25*ny(n6)+y7
z=0.25*nz(n6)+z7
n=n+1
n14=n
n,n14,x,y,z
x=0.25*nx(n5)+x7
y=0.25*ny(n5)+y7
z=0.25*nz(n5)+z7
n=n+1
n16=n
n,n16,x,y,z
n4=n3
n8=n7
nsel,all
type,3
en,elmi,n1,n2,n3,n4,n5,n6,n7,n8
emore,0,n10,0,n12,0,n14,0,n16
emore,
*enddo
cmsel,u,cracktip
nummrg,node
nsel,all
esel,all
/gopr
*end
!到此宏已经定义完毕
/prep7
*afun,deg
!使用的角度数
innerradius=0.5
outerradius=1
!外半径和内半径分别为0.5和1
scaler=0.015
baseheight=0.51
layerheight=0.18
layeramount=35
rotationangle=6
rotationtimes=90/rotationangle
height=layeramount*layerheight+baseheight
smrt,off
/title,fracture mechanics stress intensity - crack in a finite width plate
/com, ********crack in 3-Ddimensions using solid45 and solid95
antype,static
et,1,solid45
et,2,solid45
et,3,solid95
!定义单元属性
mp,ex,1,2e4
mp,nuxy,1,0.3
!定义材料属性
local,33,1,innerradius
!定义一个局部坐标系
csys,33
n,1
ngen,9,20,1
n,11,scaler
n,171,scaler,180
fill,11,171,7,31,20
!用7个节点填充11节点和171之间的节点，节点序号间隔20
local,33,0,innerradius
!定义新的坐标系
csys,33
fill,1,11,9,2,1,9,20,3
!填充节点
n,15,outerradius-innerradius
!定义15号节点
n,75,outerradius-innerradius,baseheight
!定义75号节点
fill,15,75,2,35,20
!填充
n,155,-innerradius,baseheight
!定义155号节点
ngen,2,2000,155
!此命令需进一步弄清楚====================================
csys,0
n,135,0.6,baseheight
!定义节点
csys,33
fill,75,135,2,95,20,,,0.5
!填充
n,172,-innerradius
ngen,2,2000,172
fill,155,172,5,177,-1,,,0.15
ngen,2,2000,173,177,1
csys,33
fill,11,15,3,,,7,20,3
fill,171,172,5,180,7,,,2
fill,135,155,5,186,7
fill,180,186,5,181,1,5,7,6
!填充完毕
local,44,1,,,,,-90,,0.75
!定义新的44号坐标系
csys,44
ngen,2,2000,all,,,,rotationangle
!下面开始生成单元

!=====================================
csys,0
e,2,22,1,1,2002,2022,2001,2001
egen,8,20,-1
e,2,3,23,22,2002,2003,2023,2022
egen,8,20,-1
egen,9,1,-8
!=====================================
type,2
emodif,1
*repeat,8,1
nummrg,node
NSEL,S,LOC,X,0
NSEL,R,LOC,Y,0
cm,cracktip,node
/nerr,0
fract,2
!处理裂尖的单元，
csys,0
type,1
e,208,2208,172,172,209,2209,173,173
egen,6,1,-1
!与ngen命令要好好的研究=========================
e,213,2213,177,177,214,2214,155,155
e,201,2201,2208,208,202,2202,2209,209
egen,6,1,-1
e,194,2194,2201,201,195,2195,2202,202
egen,6,1,-1
e,187,2187,2194,194,188,2188,2195,195
egen,6,1,-1
e,180,2180,2187,187,181,2181,2188,188
egen,6,1,-1
e,171,2171,2180,180,151,2151,2181,181
e,151,2151,2181,181,131,2131,2182,182
e,131,2131,2182,182,132,2132,2183,183
egen,4,1,-1
e,11,12,32,31,2011,2012,2032,2031
egen,4,1,-1
e,31,32,52,51,2031,2032,2052,2051
egen,4,1,-1
e,51,52,72,71,2051,2052,2072,2071
egen,4,1,-1
e,71,72,92,91,2071,2072,2092,2091
egen,4,1,-1
e,91,92,112,111,2091,2092,2112,2111
egen,4,1,-1
e,111,2111,2131,131,112,2112,2132,132
egen,4,1,-1
!生成了契模型单元
csys,44
egen,rotationtimes,2000,all,,,,,,,,,rotationangle
!现在不明白产生的裂纹是圆形的还是椭圆形的
!旋转的度数为6度，旋转角度为90度，旋转次数为15
!旋转结果不是一个圆的原因应该在坐标系44这里
!44号坐标系0.75为椭圆比例，所以生成的裂纹为椭圆形裂纹
!=========================================================
csys,0
n,216,1.2
!********定义216号节点
fill,15,216,1,215
n,10216,1.2,0,-0.5
fill,10015,10216,1,10215
fill,216,10216,4,2216,2000
csys,44
fill,215,10215,4,2215,2000
csys,0
n,218,1.2,0,-1.2
fill,10216,218,1,217
n,30218,0,0,-1.2
fill,218,30218,9,12218,2000
n,30215,0,0,-0.85
fill,30215,30218,2,30216,1
csys,44
fill,217,30217,9,12217,2000
fill,10216,30216,9,12216,2000
fill,10215,30215,9,12215,2000
csys,0
n,220,1.2,0.17
fill,35,220,1,219
n,10220,1.2,0.17,-0.5
fill,10035,10220,1,10219
fill,220,10220,4,2220,2000
csys,44
fill,219,10219,4,2219,2000
csys,0
n,222,1.2,0.17,-1.2
fill,10220,222,1,221
n,30222,0,0.17,-1.2
fill,222,30222,9,12222,2000
n,30219,0,0.17,-0.85
fill,30219,30222,2,30220,1
csys,44
fill,221,30221,9,12221,2000
fill,10220,30220,9,12220,2000
fill,10219,30219,9,12219,2000
e,15,215,2215,2015,35,219,2219,2035
egen,15,2000,-1
e,215,216,2216,2215,219,220,2220,2219
egen,15,2000,-1
e,10216,217,12217,12216,10220,221,12221,12220
e,12216,12217,14217,14216,12220,12221,14221,14220
egen,9,2000,-1
e,217,218,12218,12217,221,222,12222,12221
e,12217,12218,14218,14217,12221,12222,14222,14221
egen,9,2000,-1
csys,0
egen,3,20,2101,2115,1,,,,,,0,0.17
egen,3,100,2116,2150,1,,,,,,0,0.17
nsel,s,loc,x,1.2
ngen,2,10,all,,,0.2
e,216,226,2226,2216,220,230,2230,2220
egen,5,2000,-1
e,220,230,2230,2220,320,330,2330,2320
egen,5,2000,-1
e,320,330,2330,2320,420,430,2430,2420
egen,5,2000,-1
e,10

216,10226,227,217,10220,10230,231,221
e,10220,10230,231,221,10320,10330,331,321
e,10320,10330,331,321,10420,10430,431,421
e,217,227,228,218,221,231,232,222
e,221,231,232,222,321,331,332,322
e,321,331,332,322,421,431,432,422
nsel,s,loc,x,1.4
esln,s,0,all
egen,15,300000,all,,,,,,,,0.2
nsel,s,loc,y,baseheight
ngen,2,60000,all,,,,layerheight
csys,44
e,214,2214,155,155,60214,62214,60155,60155
egen,rotationtimes,2000,-1
e,207,2207,2214,214,60207,62207,62214,60214
egen,rotationtimes,2000,-1
e,200,2200,2207,207,60200,62200,62207,60207
egen,rotationtimes,2000,-1
e,193,2193,2200,200,60193,62193,62200,60200
egen,rotationtimes,2000,-1
e,186,2186,2193,193,60186,62186,62193,60193
egen,rotationtimes,2000,-1
e,135,2135,2186,186,60135,62135,62186,60186
egen,rotationtimes,2000,-1
e,115,2115,2135,135,60115,62115,62135,60135
egen,rotationtimes,2000,-1
e,95,2095,2115,115,60095,62095,62115,60115
egen,rotationtimes,2000,-1
e,75,2075,2095,95,60075,62075,62095,60095
egen,rotationtimes,2000,-1
e,75,259,2259,2075,60075,60259,62259,62075
egen,15,2000,-1
e,259,420,2420,2259,60259,60420,62420,62259
egen,15,2000,-1
e,10420,421,12421,12420,70420,60421,72421,72420
e,12420,12421,14421,14420,72420,72421,74421,74420
egen,9,2000,-1
e,421,422,12422,12421,60421,60422,72422,72421
e,12421,12422,14422,14421,72421,72422,74422,74421
egen,9,2000,-1
csys,0
e,420,430,2430,2420,60420,60430,62430,62420
egen,5,2000,-1
e,10420,10430,431,421,70420,70430,60431,60421
e,421,431,432,422,60421,60431,60432,60422
e,430,300430,302430,2430,60430,360430,362430,62430,60430
egen,5,2000,-1
e,10430,310430,300431,431,70430,370430,360431,60431
e,431,300431,300432,432,60431,360431,360432,60432
egen,14,300000,-8
nummrg,all
numcmp,all
nsel,s,loc,y,baseheight+layerheight
esln,s,0,all
egen,layeramount,400000,all,,,,,,,,,layerheight
!生成20层这么多，每层的高度为0.18
nummrg,node
nsel,all
esel,all
/nerr,defa
/output
outpr,all
nsel,all
nsel,s,loc,x,0
dsym,symm,x
d,all,ux
nsel,all
nsel,s,loc,y,0
csys,44
nsel,r,loc,x,innerradius,10
dsym,symm,y
d,all,uy
csys,0
nsel,all
nsel,s,loc,y,height
sf,all,pres,-1000
nsel,all
esel,all
finish
/output,scratch
/solu
solve
finish
/output
/post1
c*****in post1 determine ki (stress intensity factor) using kcalc!*****
csys,33
path,ki1,3,,48
ppath,1,1
ppath,2,111
ppath,3,112
kcalc,,,1
!624.58,626.27,630.26,637.98,647.78,641.90,482.85