| 查看: 2844 | 回复: 2 | ||
[求助]
pfc3d空心圆柱建模，请教一下问题出在哪？
|
|
;fname: triax_1.DAT SET random def make_walls ; create walls: a cylinder and two plates extend = 0.1 rad_cy = 0.5*width rad_cy1=0.2*width w_stiff= 1e8 _z0 =-extend _z1 = height*(1.0 + extend) command wall type cylinder id=1 kn=w_stiff end1 0 0 _z0 end2 0 0 _z1 rad rad_cy rad_cy wall type cylinder id=1000 kn=w_stiff end1 0 0 _z0 end2 0 0 _z1 rad rad_cy1 rad_cy1 end_command _x0 = -rad_cy*(1.0 + extend) _y0 = -rad_cy*(1.0 + extend) _z0 = 0.0 _x1 = rad_cy*(1.0 + extend) _y1 = -rad_cy*(1.0 + extend) _z1 = 0.0 _x2 = rad_cy*(1.0 + extend) _y2 = rad_cy*(1.0 + extend) _z2 = 0.0 _x3 = -rad_cy*(1.0 + extend) _y3 = rad_cy*(1.0 + extend) _z3 = 0.0 command wall id=5 kn=w_stiff face (_x0,_y0,_z0) (_x1,_y1,_z1) (_x2,_y2,_z2) & (_x3,_y3,_z3) end_command _x0 = -rad_cy*(1.0 + extend) _y0 = -rad_cy*(1.0 + extend) _z0 = height _x1 = -rad_cy*(1.0 + extend) _y1 = rad_cy*(1.0 + extend) _z1 = height _x2 = rad_cy*(1.0 + extend) _y2 = rad_cy*(1.0 + extend) _z2 = height _x3 = rad_cy*(1.0 + extend) _y3 = -rad_cy*(1.0 + extend) _z3 = height command wall id=6 kn=w_stiff face (_x0,_y0,_z0) (_x1,_y1,_z1) (_x2,_y2,_z2) & (_x3,_y3,_z3) end_command end ; ---------------------------------------------------- def assemble ; assemble sample s_stiff=0.0 ; initial stiffnesses n_stiff=1e8 tot_vol = height * pi * rad_cy^2.0-height * pi * rad_cy1^2.0 rbar = 0.5 * (rlo + rhi) num = int((1.0 - poros) * tot_vol / (4.0 / 3.0 * pi * rbar^3)) mult = 1.6 ; initial radius multiplication factor rlo_0 = rlo / mult rhi_0 = rhi / mult command gen id=1,num rad=rlo_0,rhi_0 x=-1.0,1.0 y=-1.0,1.0 z=0.0,height & filter ff_cylinder prop dens=1000 ks=s_stiff kn=n_stiff end_command ii = out(string(num)+'particles were created') sum = 0.0 ; get actual porosity bp = ball_head loop while bp # null sum = sum + 4.0 / 3.0 * pi * b_rad(bp)^3 bp = b_next(bp) end_loop pmeas = 1.0 - sum / tot_vol mult = ((1.0 - poros) / (1.0 - pmeas))^(1.0/3.0) command ini rad mul mult cycle 1000 prop ks=1e8 fric 0.25 cycle 250 end_command end ; ---------------------------------------------------- def cws ; change lateral wall stiffnesses command wall type cylinder id 1 kn=w_stiff wall type cylinder id 1000 kn=w_stiff end_command end ; ---------------------------------------------------- def ff_cylinder ff_cylinder = 0 _brad = fc_arg(0) _bx = fc_arg(1) _by = fc_arg(2) _bz = fc_arg(3) _rad = sqrt(_bx^2 + _by^2) if _rad + _brad > rad_cy then if _rad - _brad < rad_cy1 then ff_cylinder = 1 end_if end_if end ; ---------------------------------------------------- macro zero 'ini xvel 0 yvel 0 zvel 0 xspin 0 yspin 0 zspin 0' SET height=4.0 width=2.0 rlo=0.15 rhi=0.2 poros=0.6 make_walls assemble set w_stiff= 1e8; make lateral wall stiffness=1/10 of ball stiffness cws cyc 500 zero plot create assembly plot set cap size 25 plot set mag 2.5 plot add ball lorange plot add wall red save tt_ass.SAV return |
回复此楼» 猜你喜欢
拟解决的关键科学问题还要不要写
已经有8人回复
-
26申博
已经有3人回复
-
存款400万可以在学校里躺平吗
已经有22人回复
-
最失望的一年
已经有4人回复
-
国自然申请面上模板最新2026版出了吗？
已经有19人回复
-
请教限项目规定
已经有3人回复
-
基金委咋了？2026年的指南还没有出来？
已经有10人回复
-
基金申报
已经有6人回复
-
推荐一本书
已经有13人回复
-
疑惑？
已经有5人回复
2楼2015-12-11 20:58:35
3楼2016-03-13 11:36:07