FLAC-命令流及其解释大全

FLAC-命令流及其解释⼤全
1. FLAC3D是有限元程序吗？
答：不是！是有限差分法。

2. 最先需要掌握的命令有哪些？
答：需要掌握gen, ini, app, plo, solve等建模、初始条件、边界条件、后处理和求解的命令。

3. 怎样看模型的样⼦？
答：plo blo gro可以看到不同的group的颜⾊分布
4. 怎样看模型的边界情况？
答：plo gpfix red sk
5. 怎样看模型的体⼒分布？
答：plo fap red sk
6. 怎样看模型的云图？
答：位移：plo con dis (xdis, ydis, zdis)
应⼒：plo con sz (sy, sx, sxy, syz, sxz)
7. 怎样看模型的⽮量图？
答：plo dis (xdis, ydis, zdis)
8. 怎样看模型有多少单元、节点？
答：plo info
9. 怎样输出模型的后处理图？
答：File/Print type/Jpg file，然后选择File/Print，将保存格式选择为jpe⽂件
10. 怎样调⽤⼀个⽂件？
答：File/call或者call命令
10. 如何施加⾯⼒？
答：app nstress
11. 如何调整视图的⼤⼩、⾓度？
答：综合使⽤x, y, z, m, Shift键，配合使⽤Ctrl+R，Ctrl+Z等快捷键
12. 如何进⾏边界约束？
答：fix x ran （约束的是速度，在初始情况下约束等效于位移约束）
13. 如何知道每个单元的ID？
答：⽤⿏标双击单元的表⾯，可以知道单元的ID和坐标
14. 如何进⾏切⽚？
答：plo set plane ori (点坐标) norm (法向⽮量)
plo con sz plane (显⽰z⽅向应⼒的切⽚)
15. 如何保存计算结果？
答：save +⽂件名
16. 如何调⽤已保存的结果？
答：rest +⽂件名；或者File / Restor
17. 如何暂停计算？
答：Esc
18. 如何在程序中进⾏暂停，并可恢复计算？
答：在命令中加⼊pause命令，⽤continue进⾏继续
19. 如何跳过某个计算步？
答：在计算中按空格键跳过本次计算，⾃动进⼊下⼀步
20. Fish是什么东西？
答：是FLAC3D的内置语⾔，可以⽤来进⾏参数化模型、完成命令本⾝不能进⾏的功能
21. Fish是否⼀定要学？
答：可以不⽤，需要的时候查Mannual获得需要的变量就可以了
22. FLAC3D允许的命令⽂件格式有哪些？
答：⽆所谓，只要是⽂本⽂件，什么后缀都可以
23. 如何调⽤⼀些可选模块？
答：config dyn (fluid, creep, cppudm)
后注：这个⼯作很繁琐，需要的时间很多，希望⼴⼤⽹友能够将⾃⼰曾经遇到的常见问题在后续跟贴，也为了将这个FAQ进⾏很好的充实。

24 .如何在圆柱体四周如何施加约束条件？
可以⽤fix ... ran cylinder end1 end2 radius r1 cylinder end1 end2 radius r2 not，其中r2
25. 如何能把⼀个PLOT的图像数据导出来以便⽤其他软件绘图？
回答：⽤set log on 命令，把数据导出来，转到excel⾥处理⼀下，然后⽤surfer或者什么作图软件绘制就⾏了。

26.关于dd、dip的理解：
dip-direction angle, measured in the global xy-plane clockwise from the positive y-axis 这个是从y轴到x顺时针⽅向所转过的⾓度；dip dip angle, measured in the negative z-direction from the global xy-plane 这个是xy平⾯与z轴负⽅向的夹⾓。

dip表⽰对称参照⾯与xy平⾯的夹⾓，对称参照⾯与xy平⾯的夹⾓在xy平⾯的投影是⼀条射线，dd表⽰y轴正向顺时针到那条射线的夹⾓。

⾸先应该按照dd的⽅向⼤体确定这个⾯的朝向，dd 指的是从y轴正⽅向按顺时针（clockwise）⽅向转向所要确定⾯的法线⽅向在xy平⾯上的投影的夹⾓，然后再确定dip，dip指的是从xy平⾯转向所要确定的平⾯的⾓．
27. 如何调⽤FISH⽂件？
新⼿调⽤FISH⽂件出错的原因：是：⽂件名为fishcall.fis.txt(因为在记事本⾥这个txt的后缀是默认的）。

修改⽅法是：将fishcall.fis.txt⾥⾯的内容全部复制到新的记事本⾥，再按保存，保存时候将保存类型改为所有⽂件(*,*),这时你便可在⽂件名栏内输⼊fishcall.fis，类型也就正确了。

28. 如何理解Tension-p shear-p tension-n shear-n表⽰的含义？
在计算的循环⾥⾯,每个循环中,每个zone都依据failure criterion处于不同的状态,shear 和tension分别表⽰因受剪和受拉⽽处于塑性状态.n表⽰now,p表⽰previous,即分别表⽰在现在和以前的循环当中处于塑性状态.
29. 绘制相关曲线命令的⼀些解释
1) 在plot hist m vs n的形式⾥,m代表y轴,n代表x轴(不管m,n的正负);
2) "-"表⽰对其值作"mirror",如果⼤家使⽤AUTOCAD的话就明⽩"mirror"的含义了,⽐如上⾯的hist 1 gp ydisp的值是0~100,那么vs -1就变成-100~0.以此类推.
3) 为什么要这么做呢?主要是为了符合岩⼟⼯程的习惯需要,例如桩载荷试验曲线,沉降曲线等都不是画在第⼀象限内.
30. 如何查看剖⾯上任⼀点的位移值？
print gp disp range x() y( ) z( )或
print gp disp range id
plot gp disp range x() y( ) z( )或
plot gp disp range id
31. fix x y z range x -0.1 0.1 ..........................a
fix x range x -0.1 0.1................................b
请问a式与b式有何区别？
a 表⽰在x=0 这个平⾯是固定的
b 表⽰在x=0 这个平⾯只x⽅向不能发⽣运动
32. initial 与apply 有何区别？
initial初始化命令，如初始化计算体的应⼒状态等；
apply边界条件限制命令，如施加边界的⼒、位移等约束等。

initial的应⼒状态会随计算过程的发⽣⽽发⽣改变，⼀般体⼒需要初始化，⽽apply施加的边界条件不会发⽣变化。

33. 如何查看各个时段不平衡⼒的具体数值？
采⽤his来记录计算，包括位移应⼒等命令his unbal
his gp(zone) zdis range (0 0 0) 或者id=？
导出数据命令his write n vs m begin 时步 end 时步file filename.his
n表⽰纪录的id m表⽰时步
要导出不平衡⼒的具体数值
his unbal
step 100000 or solve
his write 1 vs step begin 1 end 1000 file 123.his
使⽤上述命令就可以查看各个时步下的不平衡⼒的具体数值
34. 如何显⽰变形轮廓线的命令？
plo ske magf 10 其中10为放⼤系数
2．FLAC2d命令的FAQ：
1、怎样开始⼊门？
答：初步学习⾸先看说明书User's gurid的Section2 Getting started和Section 3 Problem sloving with FLAC
2、哪些命令需要掌握？
答：grid gen prop interface ini fix set hist step solve pause cont sav 等，在学习的过程中遇到这些命令和不会的命令应该查看说明书的Command Reference，⾥⾯讲的很清楚，特别需要注意的是⼀定要仔细看⼀下Command Reference开头的Optional Range Phrases，这是讲范围的，因为很多命令的施加都需要加到模型的⼀部分，这⾥就告诉你怎样来表⽰这⼀部分。

这⼀原则同样适合与其他5个ITASCA软件。

还有⼀点就是六个软件（FLAC FALC3D UDEC 3DEC PFC PFC3D）的range⼤不相同，要注意区分。

3、FLAC和FLAC3D建模有什么区别？
答：FLAC和FLAC3D基本的建模思想⼤致是⼀样的，不同的是FLAC是⾸先⽤GRID命令把整体的⽹格画出来，然后⽤
MODELL NULL命令把⽹格分成许多部分，再⽤GEN命令把这些部分⼀块⼀块搭接起来。

然后再⽤GEN LINE，GEN ARC等命令画出要开挖的边界，在FLAC中产⽣⽹格的命令只有⼀个GRID，⽽且只能⽤⼀次。

⽽FLAC3D是直接⽤GEN 命令⼀块⼀块搭接，GEN可以⽤很多次。

4、怎样查看模型？
答：plot grid 可以查看⽹格，plot grid num 可以查看节点号。

5、怎样查看边界条件？
答：plot fix
6、怎样查看预先留设的开挖边界？
答：plot grid mark。

7、怎么查看模型的云图？
答：pl syy (sxx xdis ydis sig1..) fill ，注意，带fill关键字就是查看云图，不带fill 关键字就是查看等直线图。

8、怎样查看模型的⽮量图？
答：pl grid vel (disp) red
9、怎样查看INTERFACE？
答：pl grid iface id red
print iface id ; id 是interface的编号
10、FLAC怎么调整视图？
答：ctrl+z 或shift+z ,然后点击要放⼤部位的左上⾓和右下⾓就可以放⼤某个部位，如果要缩⼩，双击屏幕。

3 FLAC 使⽤步骤
FLAC 使⽤步骤
1.FLAC程式使⽤前准备步骤
步骤1：依⽐例画出所欲分析的资料
在纸上画出地点的位置、地层资料、并简单标⽰距离及深度资料。

步骤2：换算输⼊资料成同⼀单位
将现有地层资料，如 Density, Bulk modulus, Young`s modulus, tension, cohesion, friction Angle 等资料，换算成同⼀单位。

附注：需谨慎检查输⼊资料的单位，如因单位不同⽽造成过⼤或过⼩的值，将会造成 FLAC ⽆法计算，⽽产⽣ERROR 信息。

步骤3：应⽤公式简略计算
应⽤公式或依据经验，简略算出FLAC 输出资料的范围，以作为Debug 及输出资料分析时验证。

步骤4：建⽴x, y 坐标与node i, j 之间的关系
在图上距离及深度的关系，建⽴x, y 坐标系统，再由x, y 坐标系统，转换与⽹格间系，为了便于以后输出资料的分析，故应确实掌握⽹
格的位置及其相对应的x, y 坐标。

建议在敏感区域使⽤较密的⽹格，其它地⽅则使⽤较疏的⽹格，刚开始执⾏程序时，不宜使⽤⽹格太⼤的⽹格数⽬，因尽量使⽹格总数少
于1000，以节省时间。

2.FLAC 3D 程序的编写步骤：
1 Config ________
2 Grid ________
3 Model ________
4求起始的应⼒平衡
(1)建⽴x, y 坐标与⽹格的关系，建议使⽤Gen 指⽰：
Gen x1,y1 x2,y2 ,x3,y3 x4,y4 i=i0,i1 j=j0,j1 详细指令参见使⽤⼿册，FLAC 程式可⾃动产⽣x, y 坐标与⽹格的关系，但由于产⽣的⽹格坐标不易控制，将对其它的⼯作产⽣负⾯影响，故依使⽤前步骤4所建⽴的关系，将⽹格依其疏密程度需要的不同，实际控制⽹格
的坐标。

(2)设定材料性质：prop
(3)设定外⼒：Set Grav, Apply Pressure, ini sxx, Syy
(4)设定边界条件：fix, free
(5)求起始的应⼒平衡：solve
(6)储存：Save
5求⼯程的影响
求出区域内的应⼒分布情况后，再依⼯程的流程及步骤阶段执⾏各⼯程进⾏过程的影响，建议使⽤以下的步骤：
(1)调出起初的应⼒平衡：re_____ .sav
(2)设定新的材料性质：model,prop
(3)设定新的⽀撑性质：struct
(4)设定新的外⼒ apply
(5)设定边界条件 apply,fix
(6)求⼯程时的应⼒平衡 solve, step, cycle
(7)储存 save file.sav
4．Flac3d 的命令
1.在FLAC3D 中⽤到的是体积模量K 和剪切模量G ，它们与杨⽒模量E 与泊松⽐υ之间的关系为：)
1(2,)21(3υυ+=-=E G E K 。

2.在FLAC3D 中边坡的稳定性计算原理仍为折减系数法:k
k c c 折减前折减后折减前折减后??tan tan ,==
,其中属于临界值的K 即为边坡的安全系数。

3.边坡稳定性的判断⽅法之⼀为最⼤不平衡⼒的收敛性，最⼤不平衡⼒的定义为：在每⼀个计算循环中，外⼒通过⽹格节点传递分配到体系各节点时，所有节点的外⼒与内⼒之差中的最⼤值。

最⼤不平衡⼩于某⼀临界值作为边坡的收敛标准。

命令流为：set mech force 。

4.由于FLAC3D 在出图背景是灰⾊，图形为彩⾊，在论⽂打印时若我们需要⿊⽩或者灰度的图形，可以使⽤快捷键CTRL+G 或者在图形状态下执⾏settings/grayscale 命令，将图形转换为灰度显⽰。

5.动画的⽣成，flac3d 可以⽣成动画，⼤致步骤为，⾸先命令要显⽰的图像，设置动画，动画开启，计算。

命令流为：
Ini xd 0 yd 0 zd 0；将各个⽅向的位移清零
Ini xv 0 yv 0 zv 0；将各个⽅向的应⼒清零
App nstress -100e3 ran z 2.9 3.1 x 1 2 y 1 2；在指定体内添加应⼒
Plot set rot 20 0 30；将模型转动⼀定的可视⾓度
Plot con szz ou on magf 10;将竖向应⼒云图放⼤⼗倍以后输出
Plot add hist 1；在动画中加⼊最⼤不平衡⼒，编号为1
Set movie avi step 1 file 1.avi ；定义动画，动画记录的频率为1时步⼀频率Movie start ；动画开始
Solve ；计算
第⼆章命令流按照顺序进⾏
2-1定义⼀个FISH 函数
new
def abc
abc = 25 * 3 + 5
End
print abc
2-2使⽤⼀个变量
new
def abc
hh = 25
abc = hh * 3 + 5
End
Print hh
Print abc
2-3对变量和函数的理解
new
def abc
hh = 25
abc = hh * 3 + 5
End
set abc=0 hh=0
print hh
print abc
print hh
new
def abc
abc = hh * 3 + 5
end
set hh=25
print abc
set abc=0 hh=0
print hh
print abc
print hh
2-4获取变量的历史记录
new
gen zone brick size 1 2 1
model mohr
prop shear=1e8 bulk=2e8 cohes=1e5 tens=1e10
fix x y z range y -0.1 0.1
apply yvel -1e-5 range y 1.9 2.1
plot set rotation 0 0 45
plot block group
def get_ad
ad1 = gp_near(0,2,0)
ad2 = gp_near(1,2,0)
ad3 = gp_near(0,2,1)
ad4 = gp_near(1,2,1)
end
get_ad
def load
load=gp_yfunbal(ad1)+gp_yfunbal(ad2)+gp_yfunbal(ad3)+gp_yfunbal(ad4) end hist load
hist gp ydis 0,2,0
step 1000
plot his 1 vs -2
2-5⽤FISH函数计算体积模量和剪砌模量
new
def derive
s_mod = y_mod / (2.0 * (1.0 + p_ratio))
b_mod = y_mod / (3.0 * (1.0 - 2.0 * p_ratio))
end
set y_mod = 5e8 p_ratio = 0.25
derive
print b_mod
print s_mod
2-6 在FLAC输⼊中使⽤符号变量
New
def derive
s_mod = y_mod / (2.0 * (1.0 + p_ratio))
b_mod = y_mod / (3.0 * (1.0 - 2.0 * p_ratio))
end
set y_mod = 5e8 p_ratio = 0.25
derive
gen zone brick size 2,2,2
model elastic
prop bulk=b_mod shear=s_mod
print zone prop bulk
print zone prop shear
2-7 控制循环
New
def xxx
sum = 0
prod = 1
loop n (1,10)
sum = sum + n
prod = prod * n
end_loop
end
xxx
print sum, prod
new
gen zone brick p0 (0,0,0) p1 (-10,0,0) p2 (0,10,0) p3 (0,0,-10) model elas plot set rotation 0 0 45
plot block group
def install
pnt = zone_head
loop while pnt #null
z_depth = -z_zcen(pnt)
y_mod = y_zero + cc * sqrt(z_depth)
z_prop(pnt, ’shear’) = y_mod / (2.0*(1.0+p_ratio))
z_prop(pnt, ’bulk’) = y_mod / (3.0*(1.0-2.0*p_ratio))
pnt = zone_next(pnt)
end_loop
end
set p_ratio=0.25 y_zero=1e7 cc=1e8
install
2-8 拆分命令⾏
new ;example of a sum of many things
def long_sum
temp = v1 + v2 + v3 + v4 + v5 + v6 + v7 + v8 + v9 + v10 long_sum = temp + v11 + v12 + v13 + v14 + v15 end
2-9 变量类型
new
def haveone
aa = 2
bb = 3.4
cc = ’Have a nice day’
dd = aa * bb
ee = cc + ’, old chap’
end
haveone
print fish
2-10 IF条件语句
new
def abc
if xx > 0 then
abc = 1000
else
abc = -1000
end_if
end
set xx = 10
print abc
set xx = 0
print abc
2-11 索单元⾃动⽣成
new
gen zone brick size 10 3 5
plot set rotation 0 0 45
plot block group
def place_cables
loop n (1,5)
z_d = float(n) - 0.5
command
sel cable beg 0.0,1.5,z_d end 7.0,1.5,z_d nseg 7 end_command end_loop
end
place_cables
plot grid sel geom red
new
gen zone brick size 10 3 5
plot set rotation 15 0 60
plot block group
mod mohr
prop bulk 1e8 shear .3e8 fric 35
prop coh 1e3 tens 1e3
ini dens 1000
set grav 0,0,-10
fix x y z range z -.1 .1
fix y range y -.1 .1
fix y range y 2.9 3.1
fix x range x -.1 .1
fix x range x 9.9 10.1
set large
hist unbal
solve
save cab_str.sav
ini xdis 0 ydis 0 zdis 0
hist gp xdisp 0,1,5
def place_cables
loop n (1,5)
z_d = 5.5 - float(n)
z_t = z_d + 0.5
z_b = z_d - 0.5
command
free x range x -.1,.1 z z_b z_t
solve
sel cable beg 0.0,0.5,z_d end 7.0,0.5,z_d nseg 7
sel cable beg 0.0,1.5,z_d end 7.0,1.5,z_d nseg 7
sel cable beg 0.0,2.5,z_d end 7.0,2.5,z_d nseg 7
sel cable prop emod 2e10 ytension 1e8 xcarea 1.0 &
gr_k 2e10 gr_coh 1e10 gr_per 1.0
end_command
end_loop
end
place_cables
save cab_end.sav
plot sketch sel cable force red
2-12圆形隧道开挖模拟计算
;建⽴模型
gen zon radcyl p0 0 0 0 p1 6 0 0 p2 0 1 0 p3 0 0 6 &
size 4 2 8 4 dim 3 3 3 3 rat 1 1 1 1.2 group outsiderock
gen zone cshell p0 0 0 0 p1 3 0 0 p2 0 1 0 p3 0 0 3 &
size 1 2 8 4 dim 2.7 2.7 2.7 2.7 rat 1 1 1 1 group concretliner fill group insiderock gen zon reflect dip 90 dd 90 orig 0 0 0 gen zon reflect dip 0 dd 0 ori 0 0 0
gen zon brick p0 0 0 6 p1 6 0 6 p2 0 1 6 p3 0 0 13 size 4 2 6 group outsiderock1
gen zon brick p0 0 0 -12 p1 6 0 -12 p2 0 1 -12 p3 0 0 -6 size 4 2 5 group outsiderock2
gen zon brick p0 6 0 0 p1 21 0 0 p2 6 1 0 p3 6 0 6 size 10 2 4 group outsiderock3
gen zon reflect dip 0 dd 0 orig 0 0 0 range group outsiderock3
gen zon brick p0 6 0 6 p1 21 0 6 p2 6 1 6 p3 6 0 13 size 10 2 6 group outsiderock4
gen zon reflect dip 90 dd 90 orig 0 0 0 range x -0.1 6.1 z -6.1 -12.1
gen zon reflect dip 90 dd 90 orig 0 0 0 range x 6.1 21.1 z -12.1 13.1
;绘制模型图
plot block group
plot add axes red
;plot set rotation 0 0 45 ⽤于显⽰三维模型
;设置重⼒
set gravity 0 0 -10
;给定边界条件
fix z range z -12.01,-11.99
fix x range x -21.01,-20.99
fix x range x 20.99,21.01
fix y range y -0.01 0.01
fix y range y 0.99,1.01
;求解⾃重应⼒场
model mohr
ini density 1800 ;围岩的密度
prop bulk=1.47e8 shear=5.6e7 fric=20 coh=5.0e4 tension=1.0e4 ;体积、剪切、摩擦⾓、凝聚⼒、抗拉强度set mech ratio=1e-4
solve
save Gravsol.sav
plot cont zdisp outl on
plot cont szz
;⽑洞开挖计算
initial xdisp=0 ydisp=0 zdisp=0
model null range group insiderock any group concretliner any
plot block group
plot add axes red
set mech ratio=5e-4
solve
save Kaiwsol.sav
plot cont zdisp
plot cont sdisp
plot cont szz
model elas range group concretliner any
plot block group
plot add axes red
ini density 2500 range group concretliner any ;衬砌混凝⼟的密度
prop bulk=16.67e9,shear=12.5e9 range group concretliner any ;衬砌混凝⼟的体积弹模、剪切弹模set mech ratio=1e-4
solve
save zhihusol.sav
plot cont zdisp plot cont sdisp plot cont szz plot cont xzz ;完成计算分析
第四章命令流按照顺序进⾏
4-1数组
new
def afill ;fill matrix with random numbers
array var(4,3)
loop m (1,4)
loop n (1,3)
var(m,n) = urand
endloop
endloop
end
def ashow ;display contents of matrix
loop m (1,4)
hed = ’ ’
msg = ’ ’+string(m)
loop n (1,3)
hed = hed + ’ ’+string(n)
msg = msg + ’ ’+string(var(m,n))
endloop
if m = 1
dum = out(hed)
endif
dum = out(msg)
endloop
4-2函数操作
new
def xxx
aa = 2 * 3
xxx = aa + bb
end
new
def stress_sum
stress_sum = 0.0
pnt = zone head
loop while pnt # null
Stress_sum = stress_sum + z_sxx(pnt) pnt=z_next(pnt)
end_loop
end
new
def stress_sum
sum = 0.0
pnt = zone head
loop while pnt # null
sum = sum + z_sxx(pnt)
pnt=z next(pnt)
end loop
stress_sum = sum
end
new
def h_var_1
ipz = z_near(1,2,3)
H_var_1 = sxx(ipz) + fac
H_var_2 = syy(ipz) + fac
H_var_3 = szz(ipz) + fac
H_var_4 = sxy(ipz) + fac
4-3函数删除与重定义
new
def joe
ii=out(‘This is A Function’)
end
def fred
joe
end
fred
def joe
ii=(‘This is A New Function’)
end
fred
4-4字符串
new
def in_def
xx = in(msg+’(’+’default:’+string(default)+’):’) if type(xx) = 3
in_def = default
else
in_def = xx
endif
end
;
def moduli_data
default = 1.0e9
msg=’Input Young‘s modulus ’
Y_mod = in_def
default = 0.25
msg=’Input Poisson‘s ratio ’
p_ratio = in_def
if p_ratio = 0.5 then
p_ratio = in_def
endif
s_mod = y_mod / (2.0 * (1.0 + p_ratio))
b_mod = y_mod / (3.0 * (1.0 - 2.0 * p_ratio))
end
;
moduli_data
gen zone brick size 2,2,2
model elastic
prop bulk = b_mod shear = s_mod
print p_ratio y_mod b_mod s_mod
pause
print zone prop bulk
pause
print zone prop shear
4-5马蹄形隧道⽹格
new
gen zone radcyl size 5 10 6 12 rat 1 1 1 1.2 &
p0 0,0,0 p1 100,0,0 p2 0,200,0 p3 0,0,100
gen zone radtun size 5 10 5 12 rat 1 1 1 1.2 &
p0 0,0,0 p1 0,0,-100 p2 0,200,0 p3 100,0,0
gen zone reflect dip 90 dd 270 origin 0,0,0
plot set rotation 30 0 30
plot block group
4-6复杂形状⽹格⽣成
new
; 主隧道⽹格
gen zon radcyl p0 15 0 0 p1 23 0 0 p2 15 50 0 p3 15 0 8 & size 4 10 6 4 dim 4 4 4 4 rat 1 1 1 1 fill
gen zon reflect dip 90 dd 90 orig 15 0 0
gen zon reflect dip 0 ori 0 0 0
; 辅助隧道⽹格
gen point id 1 (2.969848,0.0,-0.575736)
p5 0 50 8 p6 7 0 8 p7 7 50 8 p8 point 1 p10 point 2 & size 3 10 6 4 dim 3 3 3 3 rat 1 1 1 1
gen zon radcyl p0 0 0 -1 p1 0 0 -8 p2 0 50 -1 p3 7 0 0 p4 0 50 -8 &
p5 7 50 0 p6 7 0 -8 p7 7 50 -8 p9 point 1 p11 point 2 & size 3 10 6 4 dim 3 3 3 3 rat 1 1 1 1
;衬砌⽹格
sel shell range cyl end1 0 0 -1 end2 0 50 -1 rad 3
; 隧道外围边界⽹格
gen zone radtun p0 7 0 0 p1 50 0 0 p2 7 50 0 p3 15 0 50 p4 50 50 0 & p5 15 50 50 p6 50 0 50 p7 50 50 50 & p8 23 0 0 p9 7 0 8 p10 23 50 0 p11 7 50 8 &
size 6 10 3 10 rat 1 1 1 1.1
gen zone brick p0 0 0 8 p1 7 0 8 p2 0 50 8 p3 0 0 50 &
p4 7 50 8 p5 0 50 50 p6 15 0 50 p7 15 50 50 &
size 3 10 10 rat 1 1 1.1
;对称得出1/2模型
gen zon reflect dip 0 ori 0 0 0 range x 0 23 y 0 50 z 8 50
gen zon reflect dip 0 ori 0 0 0 range x 23 50 y 0 50 z 0 50
;建⽴主隧道和辅助隧道块名称
group service range cyl end1 0 0 -1 end2 0 50 -1 rad 3
group main range cyl end1 15 0 0 end2 15 50 0 rad 4
;对称得出完整模型
gen zon reflect dip 90 dd 270 ori 0 0 0
;⽹格显⽰
plot set rotation 30 0 30
plot block group
4-7⽹格连接
new
gen zone brick size 4 4 4 p0 0,0,0 p1 4,0,0 p2 0,4,0 p3 0,0,2
gen zone brick size 8 8 4 p0 0,0,2 p1 4,0,2 p2 0,4,2 p3 0,0,4
attach face range z 1.9 2.1
model elas
prop bulk 8e9 shear 5e9
fix z range z -.1 .1
fix x range x -.1 .1
fix x range x 3.9 4.1
plot set rotation 30 0 30
plot block group
solve
save att.sav
plot cont zdisp outl on
单⼀⽹格（⼩）程序
new
gen zone brick size 8 8 8 p0 0,0,0 p1 4,0,0 p2 0,4,0 p3 0,0,4
model elas
prop bulk 8e9 shear 5e9
fix z range z -.1 .1
fix x range x -.1 .1
fix x range x 3.9 4.1
fix y range y -.1 .1
fix y range y 3.9 4.1
apply szz -1e6 range z 3.9 4.1 x 0,2 y 0,2
hist unbal
;⽹格显⽰
plot set rotation 30 0 30
plot block group
solve
save noattx.sav
plot cont zdisp outl on
单⼀⽹格（⼤）程序
new
gen zone brick size 4 4 4 p0 0,0,0 p1 4,0,0 p2 0,4,0 p3 0,0,4 model elas prop bulk 8e9 shear 5e9
fix z range z -.1 .1
fix x range x -.1 .1
fix x range x 3.9 4.1
plot set rotation 30 0 30
plot block group
solve
save noattx.sav
plot cont zdisp outl on
4-8⽴⽅体洞⽳⽹格⽣成
new
def parm
rad=4.0
len=10.0
in_size=6
rad_size=10
end
parm
gen zone radbrick edge len size in_size in_size in_size rad_size & rat 1.0 1.0 1.0 1.2 dim rad rad rad
;对称得出完整模型
gen zon reflect dip 0 dd 0 ori 0 0 0
gen zon reflect dip 90 dd 270 ori 0 0 0
;⽹格显⽰
plot set rotation 30 0 30
plot block group
4-9球体洞⽳⽹格⽣成
new
def parm
rad=4.0
len=10.0
in_size=6
rad_size=10
end
p_gp=gp_head
loop while p_gp#null
px=gp_xpos(p_gp)
py=gp_ypos(p_gp)
pz=gp_zpos(p_gp)
dist=sqrt(px*px+py*py+pz*pz)
if dist>0 then
k=rad/dist
ax=px*k
ay=py*k
az=pz*k
maxp=max(px,max(py,pz))
k=len/maxp
bx=px*k
by=py*k
bz=pz*k
u=(maxp-rad)/(len-rad)
gp_xpos(p_gp)=ax+u*(bx-ax)
gp_ypos(p_gp)=ay+u*(by-ay)
gp_zpos(p_gp)=az+u*(bz-az)
end_if
p_gp=gp_next(p_gp)
end_loop
end
gen zone radbrick edge len size in_size in_size in_size rad_size & rat 1.0 1.0 1.0 1.2 dim rad rad rad make_sphere
;对称得出完整模型
gen zon reflect dip 0 dd 0 ori 0 0 0
gen zon reflect dip 90 dd 270 ori 0 0 0
;⽹格显⽰
plot set rotation 15 0 30
plot block group
4-10应⼒边界
p2 (0,4,0) p3 (2,0,3.464)
model elastic
prop bulk 1e8 shear .3e8
apply nstress -1e6 range plane dip 60 dd 270 origin 0.1,0,0 above step 1 ;⽹格显⽰
plot set rotation 30 0 30
plot block group
plot grid red fap green
4-11改变应⼒边界条件
new
gen zone brick size 6 6 6
model elas
prop bulk 1e8 shear 7e7
fix x range x -0.1 0.1
def superstep
loop ns (1,n_steps)
x_stress = stress_inc
command
apply sxx add x_stress range x 5.9,6.1 y 0,6 z 0,2 step 100
end_command
end_loop
end
set n_steps=100 stress_inc=-1e3
hist zone sxx 6,0,0
plot create sxx_hist
plot add hist 1
plot show
Superstep
;⽹格显⽰
plot set rotation 30 0 30
plot block group
plot grid red fap green
new
fix x range x -0.1 0.1
def x_stress
x_stress = stress_inc * step
end
set stress_inc = -1e3
apply sxx 1.0 hist x_stress range x 5.9,6.1 y 0,6 z 0,2 hist zone sxx 6,0,0 hist x_stress
step 100
;⽹格显⽰
plot set rotation 30 0 30
plot grid red fap green
new
gen zone brick size 6 6 6
model elas
prop bulk 1e8 shear 7e7
fix x range x -0.1 0.1
table 1 0,0 100,-1e5
apply sxx 1.0 hist table 1 range x 5.9,6.1 y 0,6 z 0,2
hist zone sxx 6,0,0
step 100;⽹格显⽰
plot set rotation 30 0 30
plot grid red fap green
new
gen zone brick size 6,6,6 p1 6,0,-1
model elas
prop bulk 8e9 shear 5e9
apply sxx -2e6 range x -0.1 0.1
apply sxx -2e6 range x 5.9 6.1
step 500
plot bound fap
plot grid disp
new
set grav 0 0 -10
ini dens 1000
fix x range x -.1 .1
fix x range x 4.9 5.1
fix y range y -.1 .1
fix y range y 4.9 5.1
ini szz -5e4 grad 0 0 -1e4
app szz -5e4 range z -.1 .1
solve
model null range x 1,4 y 1 4 z 3 5
step 100
plot set plane dip 90 dd 180 origin 0,2.5,0
plot add surf plane behind yell
plot add vel plane behind black
plot set rotation 30 0 30
plot block group
plot show
4-12位移边界
new
gen zone brick size 4,4,4 p3 2,0,3.464
model elas
prop bulk 1e8 shear .3e8
macro left_boun ’plane dip 60 dd 270 origin 0.1,0,0 above’macro right_boun ’plane dip 60 dd 270 origin 3.9,0,0 below’apply nvel 0.1 plane dip 60 dd 270 range left_boun。