;-----------------------------------------------------
; ---- Excavation and Support for a Shallow Tunnel ---
;-----------------------------------------------------
new ;新建项目
set fish autocreate off
title 'Excavation and Support for a Shallow Tunnel' ;定义题目
; generate primitive components of grid
; concrete liner - upper tunnel
gen zon cshell p0 0 0 0 p1 7 0 0 p2 0 51 0 p3 0 0 5.5 &
dim 5 5 5 5 size 2 51 10
group zone 'concrete liner'
;
; upper tunnel
gen zon cylinder p0 0 0 0 p1 5 0 0 p2 0 51 0 p3 0 0 5 &
size 5 51 10
group zone tunnel range group 'concrete liner' not
;
; lower tunnel & liner
gen zone brick p0 0 0 -4.5 p1 add 7 0 0 p2 add 0 51 0 p3 add 0 0 4.5 &
size 7 51 3
;
; surrounding rock (8 primitives)
gen zon radcyl p0 0 0 0 p1 27 0 0 p2 0 51 0 p3 0 0 25 &
dim 7 5.5 7 5.5 size 5 51 10 8 rat 1 1 1 1.3
;
gen zone brick p0 7 0 -4.5 p1 27 0 -15 p2 add 0 51 0 p3 7 0 0 &
p4 27 51 -15 p5 7 51 0 p6 27 0 0 p7 27 51 0 &
size 8 51 3 ratio 1.3 1 1
;
gen zone brick p0 0 0 -15 p1 add 27 0 0 p2 add 0 51 0 p3 0 0 -4.5 &
p4 27 51 -15 p5 0 51 -4.5 p6 7 0 -4.5 p7 7 51 -4.5 &
size 7 51 8 rat 1 1 0.7692307692307692
;
gen zon brick p0 0 0 25 p1 add 27 0 0 p2 add 0 51 0 p3 add 0 0 10 &
size 5 51 2
;
gen zon bric p0 27 0 25 p1 add 17 0 0 p2 add 0 51 0 p3 add 0 0 10 &
size 2 51 2 rat 2 1 1
;
gen zon bric p0 27 0 -15 p1 add 17 0 0 p2 add 0 51 0 p3 add 0 0 40 &
size 2 51 8 rat 2 1 1
;
gen zon bric p0 27 0 -40 p1 add 17 0 0 p2 add 0 51 0 p3 add 0 0 25 &
size 2 51 2 rat 2 1 0.5
;
gen zon bric p0 0 0 -40 p1 add 27 0 0 p2 add 0 51 0 p3 add 0 0 25 &
size 7 51 2 rat 1 1 0.5
;
; assign names to groups of zones
group zone rock range group 'concrete liner' not group tunnel not
;
; assign Mohr-Coulomb material model
model mech mohr
pro bulk 50e6 she 18e6 fric 20 coh 25e3 ten 0 dil 0 range z 25 35
pro bulk 4e8 she 1.5e8 fric 20 coh 50e3 ten 5e3 dil 3 range z -50 25
; assign boundary conditions ;施加边界条件，后面可以直接修改为具体的数字
fix x range x -.1 .1
fix x range x 43.9 44.1
fix z range z -40.1 -39.9
fix y range y -.1 .1
fix y range y 50.9 51.1
; assign initial stress state ;初始应力状态
set grav 0 0 -10
ini density 2200
ini szz -770e3 grad 0 0 22000
ini sxx -770e3 grad 0 0 22000
ini syy -385e3 grad 0 0 11000 ;; monitor variables in model ;hist add unbal ;hist add gp zdisp 0 0 5.5
hist add gp xdisp 7 0 0
hist add gp zdisp 0 0 0
施加初始应力
模型中变量的监控监控不平衡力
hist add gp zdisp 0 0 35
hist add gp zdisp 0 30 5.5
hist add gp xdisp 7 30 0
hist add gp zdisp 0 30 0
hist add gp zdisp 0 30 35
hist add gp zdisp 0 12 35
hist add gp zdisp 0 18 35
hist add gp zdisp 0 24 35
hist add gp zdisp 0 36 35
hist add gp zdisp 5 30 35
hist add gp zdisp 10 30 35 ;;
监控以上这些点的x及z方向位移变化
sav geom1
;
def conc_parm ;定义支护参数，parm即parameter，参数的意思
global bmc = 20.7e9 ;定义体积模量为全局变量，b代表bulk，m代表modulus，c代表concrete
global smc = 12.6e9 ;定义剪切模量为全局变量，s代表shear，m代表modulus，c代表concrete
end
@conc_parm
;
; define the locations of cable patterns 1, 2 and 3
;
def cab_parm ;定义锚杆参数
global x_b = get_array(4,3) ;定义数组（4,3），即锚杆的位置
global z_b = get_array(4,3)
global y0 = -3 ;将锚杆的位置以数组的形式表示出来（x,y,z）
x_b(1,1) = 0.8
x_b(2,1) = 2.1
x_b(3,1) = 3.5
x_b(4,1) = 5.5
z_b(1,1) = 5.5
z_b(2,1) = 2.4
z_b(3,1) = 4.7
z_b(4,1) = 1.5
x_b(1,2) = 0.8
x_b(2,2) = 0.8
x_b(3,2) = 3.5
x_b(4,2) = 5.5
z_b(1,2) = 0.6
z_b(2,2) = 4.0
z_b(3,2) = 2.4
z_b(4,2) = 0.6
x_b(1,3) = 0.8
x_b(2,3) = 2.6
x_b(3,3) = 5.0
x_b(4,3) = 3.5
z_b(1,3) = 2.4
z_b(2,3) = 4.0
z_b(3,3) = 3.0
z_b(4,3) = 0.6 ;对每根锚杆的位置进行赋值
end
def inip(iidx) ;global x1 = x_b(1,iidx)
global x2 = x_b(2,iidx)
global x3 = x_b(3,iidx)
global x4 = x_b(4,iidx)
global z1 = z_b(1,iidx)
global z2 = z_b(2,iidx)
global z3 = z_b(3,iidx)
定义初始锚杆位置
global z4 = z_b(4,iidx)
end
@cab_parm
;
; install initial cables ;;
def ins_cab initial cables
global iidx ;global cab_seg ;global cab_seg_m ;loop iidx (1,3)
inip(iidx)
安装初始锚杆
定义要初始安装的锚杆，其中ins_cab表示install 定义全局变量iidx
定义锚杆划分单元数
锚杆长度（有待进一步确认）
;
cab_seg = cab_seg_m-3*(3-iidx)
global y1 = 0.
global y2 = float(cab_seg)
command
sel cable id @iidx begin @x1 @y1 @z1 end @x1 @y2 @z1 nseg @cab_seg
sel cable id @iidx begin @x2 @y1 @z2 end @x2 @y2 @z2 nseg @cab_seg
sel cable id @iidx begin @x3 @y1 @z3 end @x3 @y2 @z3 nseg @cab_seg
sel cable id @iidx begin @x4 @y1 @z4 end @x4 @y2 @z4 nseg @cab_seg
sel cable pro emod 45e9 xcarea 1.57e-3 gr_per 1.0 &
yten 25e4 gr_k 17.5e6 gr_c 20e4 range id @iidx ;施作初始锚杆
end_command
end_loop
end
set @cab_seg_m 15
@ins_cab
; install pre-support concrete ;预支护
;
sel shell id 10 group rock range cyl end1 0 0 -1.5 end2 0 1 -1.5 rad 7.4 &
cyl end1 0 0 -1.5 end2 0 1 -1.5 rad 6.7 not &
z -0.1 6
sel shell prop isotropic 10.5e9,0.25 thickness 0.3 density 2500
def monit
global ipt_surf = gp_near(0,30,35) ;地表
global ipt_crown = gp_near(0,30,5.5) ;拱顶
global ipt_spring = gp_near(7,30,0)
end
@monit
save m_ini
table 1 name 'ground surface at tunnel center line'
table 2 name 'tunnel crown' ;隧道拱顶
table 3 name 'tunnel sidewall' ;定义表格的名字，隧道边墙
;
; FISH function to control excavation and support sequence
def excav
y0 = y0 3
local cut_i = y0/3 1
global cut
loop cut (cut_i,16)
local cut_cur = cut
local ii = out(' EXCAVATION STEP ' string(cut))
y0 = 3*(cut-1)
y1 = y0 3
global yp0 = y0 1
global yp1 = y1 1
global ys0 = yp0-3
global ys1 = yp1-3
global yc0 = y0-3
global yc1 = y1-3
global id_ = 10
; id_ = 10*(cut 1) ; use if shells unconnected
command
; install pre support concrete
sel shell id @id_ group rock &
range cyl end1 0 @yp0 -1.5 end2 0 @yp1 -1.5 rad 7.4 &
cyl end1 0 @yp0 -1.5 end2 0 @yp1 -1.5 rad 6.7 not &
z -0.1 6
sel shell prop isotropic 10.5e9,0.25 thickness 0.3 density 2500 &
ran y @yp0 @yp1
; excavate next cut
model mech null range group tunnel y @y0 @y1
model mech null range group 'concrete liner' y @y0 @y1
; delete-cables in the excavated area
sel delete cable range id 1 y @y0 @y1
sel delete cable range id 2 y @y0 @y1
sel delete cable range id 3 y @y0 @y1
end_command
local cut_1 = cut-1
iidx=int(cut_1-3*(cut_1/3)) 1
y2=min(y1 15,51)
inip(iidx)
ii = out(' CABLE BOLT PATTERN ' string(iidx))
command
sel delete cable range id @iidx
; install new cables
sel cable id @iidx begin @x1 @y1 @z1 end @x1 @cab_seg_m
sel cable id @iidx begin @x2 @y1 @z2 end @x2 @cab_seg_m
sel cable id @iidx begin @x3 @y1 @z3 end @x3 @cab_seg_m
sel cable id @iidx begin @x4 @y1 @z4 end @x4 @cab_seg_m
sel cable pro emod 45e9 xcarea 1.57e-3 gr_per 1.0 &
yten 25e4 gr_k 17.5e6 gr_c 20e4 ran id @iidx
; shotcrete
@y2 @z1 @y2 @z2 @y2 @z3 @y2 @z4 nseg nseg nseg nseg
sel shell prop isotropic 10.5e9,0.25 thickness 0.5 density 2500 &
ran y @ys0 @ys1
end_command
if cut > 1 then
command
; concrete liner
model mech el range group 'concrete liner' y @yc0 @yc1
prop bulk @bmc sh @smc range group 'concrete liner' y @yc0 @yc1
end_command
end_if
command
step 3000
end_command
; store displacements in tables ;将位移储存在表格中
xtable(1,cut) = 3.0 * cut
ytable(1,cut) = gp_zdisp(ipt_surf)
xtable(2,cut) = 3.0 * cut
ytable(2,cut) = gp_zdisp(ipt_crown)
xtable(3,cut) = 3.0 * cut
ytable(3,cut) = gp_zdisp(ipt_spring)
command
save m1
end_command
if cut=5 then
command
save m1_15
end_command
end_if
if cut=9 then
command
save m1_27
end_command
end_if
if cut=10 then
command
save m1_30
end_command
end_if
end_loop
end
@excav
return
