UDEC模拟实例与解析

^`
UDEC 实例翻译与命令解析
翻译：珠穆朗玛
1 地震诱发地层坍塌 Seismic-Induced Groundfall
1.1 问题描述
本例展示使用 UDEC 模拟分析地震诱发地层坍塌的一类的问题,模型见图 1.1,该模型基 于加拿大安大略省萨德伯里市鹰桥公司弗雷则矿 34-1-554 切割断面的一个剖面图的结构和 尺寸. 用二维平面应变模型代表垂直于超采轴向方向的平面效应,超采面高 5m,宽 10m.
假定两个连续节理交叉平面分析:一个角度为 45 度,另一个为-9 度,两者节理间距均为 5m,为了演示的目的,一个近似垂直的“虚拟节理”也被添加到块体内开挖面顶部以增强不稳 定性。

围岩参数来自试验室平均测试数值，假定岩石块体参数如下：
假定块体仅具有弹性行为，
节理假定符合库伦滑动准则，选择典型的教课书数值作为节 理参数，如下：
初始应力状态按各向同性估计为24Mpa（假定垂直荷载由覆盖深度大约800m 的岩层产生）。

1.2 UDEC 分析
UDEC 模拟顺序分三个阶段,首先,模型在初始应力状态下进行无超采固结.其次,进行开挖并且模型循环至平衡状态.本阶段超采面周围的应力分布见图1.2.超采正上方和下方的块体滑动后稳定.在第三阶段.估计了两个不同的峰值速度的地震事件.对所有地震模拟,在问题域的外周边界引入粘滞边界用以消除波的反射.从而模拟有限的岩体,地震事件用施加到模型顶部y 方向的正弦应力波表现.应力波被叠加到已存在的初始地应力上.在第一个模拟中,施加1.25Mpa 的峰值应力,应当注意的是,由于粘滞边界条件实际是在模型顶部, 施加的有效影响应力应该是1.25 MPa/2, or 0.625 MPa.0.02 秒后的开挖面拱顶的应力分布见图1.3,两点的位移被监测,1 点位于开挖面的左角,点2 位于拱顶块体的右角, 图1.4 的位移时间曲线显示两点本质上是弹性反应
.
本例关心的问题是在模型顶部施加的速度和计算速度的对比,下面的公式可以用以估计施加的波速.
使用这个方程,施加的最大波速大概是0.04m/sec,图1.5 显示的峰值波速小于0.06m/sec. 估计的波速和监测波速的不同在于使用的围岩模量.而是没有考虑节理变形的相等变形模量.
在第二个案例中,施加应力波峰值12.5 Mpa(有效应力6.25Mpa).0.02 秒后的开挖拱顶应力分布见图1.6.该图显示出拱顶岩体不受力,表面该块体已经松散并正在下落.对于关心的问题,后来三个时间的几何体和应力分布见图1.8 至图1.10.
在问题的顶部预测的波速(从上面的方程)是0.4m/sec.从模型中计算的波速见图1.11,再次,由于使用的是原岩弹性模量而不是岩体的变形模量导致预测和监测的波速之间的差异.
1.3 节包含了该模型的数据列表,该列表包含了一个FISH 函数(show)被用来创建坍塌的动画文件,每隔0.02 秒俘获一个显示的图片.通过改变FISH 参数time_int 可以改变动画帧的
间隔.视图的总数也可以通过改变snap_shot 的数值进行改变.为了显示80 帧的显示图片而创建的该电影文件需要大概13MB 的硬盘空间.
1.3 数据文件列表
Example 1.1 SEISMIC.DAT
title
SEISMIC INDUCED ROOF COLLAPSE 地震诱发拱顶坍塌
;
round 0.01
; define original boundary of modeled region 定义模型区域的原始边界
block -25,-20 -25,20 25,20 25,-20
; generate joint pattern over entire original region 在整个原始区域生成节理形态
jregion id 1 -25,-25 -25,25 25,25 25,-25
jset 45,0 200,0 0,0 5.0,0 (0,0) range jreg 1
jset -9,0 200,0 0,0 5.0,0 (0,0) range jreg 1
; put in joints needed for the later excavation 为了后面开挖而设置的节理
crack -5.01,-2.51 5.01,-2.51
crack -5.01, 2.51 5.01, 2.51
crack -5,-2.5 -5,2.5
crack 5,-2.5 5,2.5
crack 2.25,2.5 1.93,5.0
; generate fdef zones and assign joint properties (mat=1 & jmat=1;default) 生成单元和设置节理参数generate edge 9.0 range -30,30 -30,30
prop mat=1 d=0.00300 k=39060 g=31780
prop jmat=1 jkn=20000 jks=20000
prop jmat=1 jf=30.0
; apply boundary conditions and initial conditions to 在地应力下施加边界条件和初始条件
; consolidate model under field stresses
bound stress=-24.0, 0.0, -24.0 ygrad=-0.3 0 -0.3
insitu stress=-24.0, 0.0, -24.0 ygrad=-0.3 0 -0.3
bound yvel 0.0 range -26,26 -21,-19
grav 0.0 -10.0
; track the x-displacement, and y-displacement over time 追踪位移
hist solve
hist xdis=0,7 ydis=0,7 type 1
solve rat 1e-5
; save consolidated state
save seismic1.sav
; make excavation
delete range -5,5 -2.5,2.5
solve rat 1e-5
; save excavated state
save seismic2.sav
;
rest seismic2.sav
; apply seismic load from top (peak velocity=0.04 m/sec)
;
; set up nonreflecting boundary
bound mat=1
bound xvisc range -26 -23 -21 21
bound xvisc range 23 26 -21 21
bound xvisc yvisc range -26 26 -21 -19
bound xvisc yvisc range -26 26 19 21
; apply sinusoidal stress wave
bound stress 0 0 -1.25 yhist=cos(100.0,0.0195) range -26 26 19 21 ;
reset time hist disp rot
hist ydis (-4.48,2.57)
hist ydis (0,2.57) yvel (0,2.57) yvel (4,2.57) yvel(-4.48,2.57)
hist yvel (0,20) yvel (25,10) yvel (25,-10) yvel (0,-20)
hist yvel (-25,-10) yvel (-25,10)
hist sxx (25,10) sxx (25,-10) sxx (-25,-10) sxx (-25,10)
hist syy (0,20)
;
damp 0.1 1.0 mass
; 0.02 sec.
cyc time 0.02
save seismic3.sav
;
rest seismic2.sav
; apply seismic load from top (peak velocity=0.4 m/sec)
; set up nonreflecting boundary
bound mat=1
bound xvisc range -26 -23 -21 21
bound xvisc range 23 26 -21 21
bound xvisc yvisc range -26 26 -21 -19
bound xvisc yvisc range -26 26 19 21
; apply sinusoidal stress wave
bound stress 0 0 -12.5 yhist=cos(100.0,0.0195) range -26 26 19 21 reset time hist disp
hist ydis (-4.48,2.57)
hist ydis (0,2.57) yvel (0,2.57) yvel (4,2.57) yvel(-4.48,2.57)
hist yvel (0,20) yvel (25,10) yvel (25,-10) yvel (0,-20)
hist yvel (-25,-10) yvel (-25,10)
hist sxx (25,10) sxx (25,-10) sxx (-25,-10) sxx (-25,10)
hist syy (0,20)
;
damp 0.1 1.0 mass
save seismov.sav
;
; 0.02 sec.
cyc time 0.02
save seismic4.sav
; 0.25 sec.
cyc time 0.23
save seismic5.sav
; 0.50 sec.
cyc time 0.25
save seismic6.sav
; 0.75 sec
cyc time 0.25
save seismic7.sav
;
rest seismov.sav
; make a movie of the groundfall
;
wind -12 12 -12 12
set ovtol 0.05
plot block vel max 2.0 blue stress max 50
movie on
movie file = seismic.dcx
movie step 1000
step 40000
3 隧道支护荷载Tunnel Support Loading
3.1 问题陈述
本例模拟展示了UDEC 在检查衬砌隧道方面的应用,着重强调了荷载在混凝土衬砌中的发展,本例也解释了模拟连续建造操作中独立阶段的模拟程序.
隧道系统的理想几何体见图3.1.系统包含在海床下大约70m(中线)深度,中线间距12m 的两个隧道, 初始水位在隧道中线上方110m 处.服务隧道直径5.24m,衬砌厚度37cm.主隧道直径8.22m,衬砌厚度46cm.服务隧道先于主隧道开挖和衬砌.随后设置主隧道衬砌,水位上升增加到100m.
施工顺序是:
(1)开挖服务隧道excavation of the service tunnel;
(2)衬砌服务隧道lining of the service tunnel; (3)开
挖主隧道excavation of the main tunnel; (4)衬砌主
隧道lining of the main tunnel; and (5)升高水位
raising of the water level.
分析的目的是评价每个施工阶段服务隧道和主隧道支护状况.
本例的材料参数见下:
岩体——开挖隧道的围岩参数为：
弹性模量elastic modulus 0.89 GPa
泊松比Poisson’s ratio 0.35
单轴抗压强度uniaxial compressive strength 3.5 MPa
粘聚力cohesion 1 MPa
密度density 1340 kg/m3
混凝土衬砌——弹性模量为24 GPa ，泊松比为0.19. 假定衬砌为线弹性材料。

3.2 UDEC 分析
本问题创建的UDEC 模型见图3.2，隧道中心线在Y=-70 的位置，注意，模型边界距隧
道开挖面很近，本模型为了演示的目的只是提供一个快速计算，实际计算时应采用大的模型。

底边和侧边界采用滚轴固定，海底之上海水的重量按30m 的水头施加在模型表面上，假定隧道为线性防水衬砌，因此不需要进行瞬态地下水流分析。

忽略孔隙水压力，而将岩体单位重量设置为浮容重。

垂直对水平应力比率假定为静水压力0.3。

对本例而言，上面列的5 个施工顺序按三个阶段进行模拟，主隧道开挖和服务隧道衬砌按一个瞬时的活动进行模拟。

主隧道衬砌和水位上升也假定瞬时产生，这些活动可以通过线性（梯度）降低隧道周围反力或者设定产生一些松弛后安装支护而被分开模拟。

在模拟的第一阶段，在重力应力被初始化到块体中后，服务隧道被开挖，UDEC 循环直到达到平衡状态，弹性位移结果见图3.3.
在第二阶段, 衬砌服务隧道和开挖主隧道.采用16 个梁单元来模拟服务隧道混凝土衬砌, 图3.4 和图3.5 显示了主隧道开挖后的位移和主应力分布状况.注意,整个服务隧道向主隧道转移.
在第三个阶段,一个额外的荷载被施加在模型的顶部用来模拟增加的100m 的水的重量, 使用STRUCT apply pressure 命令来施加静水压力荷载(施加到隧道衬砌上),用以代表隧道中心线上210m 的水位.图3.6 和图3.7 显示了施加额外荷载后相应的轴力和弯矩分布,服务隧道最大轴力为5.59×106N,主隧道最大轴力为9.06×106N.
3.3 数据文件列表
; tunnel support loading
;
block 0 -30 60 -30 60 -90 0 -90 round 0.1
crack 0 -70 60 -70
crack 30 0 30 -90
crack 42 0 42 -90
tun 30 -70 4.11 12
tun 42 -70 2.62 8
tun 30 -70 5.5 12
tun 42 -70 5.5 12
gen edge 2.0Example 3.1 TUNNEL.DAT
海床面
0,-30
30,0
海水面
42,0
60,-30
; initial stress state
bound stress 1.02e5 0 1.02e5 ygrad 1.34e4 0 1.34e4 ;bound stress sxxo sxyo syyo ygrad sxxy sxyy syyy ; sxx = sxxo + (sxxx·x) + (sxxy·y)0,-70
0,-9030,-90 42,-90
60,-70
90
insit stres 1.02e5 0 1.02e5 szz 1.02e5 ygra 1.34e4 0 1.34e4 zgra 0 1.34e4 gravity 0 -10
; rock properties
prop mat=1 d=1340
zone model mohr
zone shear=.33e9 bulk=.99e9 coh=1e6 fric=30.0
;prop mat=1 d=1340 g=.33e9 k=.99e9 coh=1e6 fric=30.0
;change con=3
; elastic joint properties
prop jmat=1 jkn=1e9 jks=1e9 jcoh=1e10 jtens=1e10
; cycle to initial equilibrium
hist solve type 1
solve rat 1e-5
save tun1.sav
; excavate service tunnel
del 40 44 -72 -68
bound -1 1 -91 0 xvel=0.0
bound -1 90 -91 -89 yvel=0.0
bound 59 61 -91 0 xvel=0.0
; histories around tunnel 1
hist ydis 42 -67 sxx 42 -67
hist ydis 42 -73 sxx 42 -73
hist xdis 39 -70 syy 39 -70
hist xdis 45 -70 syy 39 -70
; histories around tunnel 2
hist ydis 30 -65.0 sxx 30 -65.0
hist ydis 30 -75.0 sxx 30 -75.0
hist xdis 25.0 -70 syy 25.0 -70
hist xdis 35.0 -70 syy 35.0 -70
reset disp jdisp
solve rat 1e-5
save tun2.sav
; line service tunnel
struct gen xc=42 yc=-70 npoint=16 mat=5 thick=0.37 fang=-11.25 theta 360 prop mat=5 st_d=2400 st_ymod=24.0e9 st_prat=0.20 st_yield=1e10
prop mat=5 if_kn=1e8 if_ks=1e7 if_coh=1.0e10
; excavate main tunnel
del 28 32 -72 -68
reset disp
solve rat 1e-5
save tun3.sav
; line main tunnel
struc gen xc=30 yc=-70 npoint=8 mat=5 thick=0.46 fang=22.5 theta -360
; add additional load representing raised water level
bound stress 0.0 0.0 -1.0e6 range -1 91 -31 -29
; add hydrostatic loads to tunnel liners
struct apply press 0.0 2.06e6
reset disp
solve rat 1e-5
save tun4.sav
ret
4 流体和动态荷载下的重力坝分析
4.1 问题描述
本例展示了一个坐落在节理岩体基础上的100m 高混凝土重力坝建造问题,节理平均间距是50m,节理方向为20 度和-70 度, 主要研究了两个荷载条件,首先是研究包含流体在岩体裂隙内渗流的水库蓄水效应分析.其次,为了研究潜在地震类型的荷载,一个动态波被施加在模型的基础上进行分析.
4.2 UDEC 分析
本问题的UDEC 模型见图4.1。

模型理想化的展示了一个坐落在节理岩体基础上的重力坝，目的是展示这种类型问题的推荐求解过程。

数据文件见4.3 节，按下面的顺序进行分析。

第一步：重力荷载——空的水库
假定初始应力状态下的有效应力比率为σH/σV = 0.69，水位假定在y = 0，由大坝重量产生的初始应力状态和未蓄水的水库见图4.1。

注意，使用INSITU 命令设定的应力是总应力，被加载在块体上。

对节理裂隙，UDEC 计算有效应力并将域压力设置到静水压力。

第二步：水库蓄水本阶段，假定水位上升到大坝的顶部，在大坝上游侧和岩石基础处施
加静水压力作用，
在侧边界由于施加在大坝上的荷载产生的水平反作用力被认为是滚轴支撑。

对于流体渗流假定下面的条件：
1. 假定沿着模型底部和侧边的节理裂隙连接渗透性为零。

2. 大坝上游的岩体面，通过使用BOUND pp 命令将水头固定在100m（0.98Mpa）。

下游
面，水头设置为零。

3. 大坝和岩体基础交界面假定有较低的渗透性。

4. 使用稳定流体算法（SET flow steady）
选定的第二步结果见图4.2 至图4.6，蓄水后的位移结果见图4.2。

大坝顶部的X 方向和Y 方向位移历史见图4.3。

后来的图形显示蓄水后模型处于平衡状态。

图4.4 的流体速率显示了大多数流体被集中在大坝基础下面的裂隙内。

图4.5 显示了沿着节理（at x = -33.42, y = -30.37）的法向和剪切位移。

正法向位移表明在本阶段节理张开。

在位置x = -22.1, y = -26.3 处沿着裂隙的流体压力历史见图4.6。

^`
^`
第三步：动态荷载
^`
本阶段中，一个垂直传播的正弦波（frequency = 5 Hz）被施加在模型基础上10 秒种。

下面的边界条件被采用。

1. 假定底部边界在水平（剪切）方向为无反射边界且固定垂直方向，动态输入以剪切
应力历史的形式被施加。

2. 类似的条件施加在自由区域的的底部。

3. 基础一边的节点与自由区域的计算相关联。

自由区域被离散化成20 个单元，并设定
为弹塑性行为。

在自由区域使用FFIELD init 命令初始化当前静力水平以匹配模型
中的网格。

4. 在动态阶段，由渗流阶段产生的静力作用仍旧被施加块体上。

因此，在无动态荷载
作用时，块体仍处于平衡状态。

在本阶段，渗流计算被关闭(SET flow off)。

这是近似假定当10 秒动态荷载时无流体产生。

图4.7 至4.9 显示了本阶段在1.5 秒后动态荷载作用下的分析结果。

每个荷载下大坝顶部X 和Y 方向的位移历史见图4.7。

在大坝下方节理的剪切和方向位移（图4.9）显示了剪切位移的累积情况。

注意，在图4.10 中，沿着大坝基础下的第一个节理产生了大量滑移；上游70 度节理张开（无有效应力），这个边产生了较大位移，连同累计的位移显示在图4.8 和4.9.表明大坝很可能破坏。

^`
4.3 数据文件列表
Example 4.1 DAM.DAT
;
; --- dam --- discontinuous joints : 20 and -70 deg. 大坝---不连续节理20和-70度
; --- insitu stresses (k=0.5) 初始应力(K=0.5) ? 与文中叙述不符,应为0.69否
; --- free-field (20 nodes) --- applied only in phase 3 自由区域(20节点),仅在第三阶段应用;第一阶段--- 重力荷载;固定X方向边界
; --- phase 1 --- gravity loading; x-fixed boundaries
;
title
(DAM1) IN-SITU STRESSES ; GRAVITY
;
; rock blocks ; joints (no cohesion) 岩块节理无黏聚力
prop mat=1 d=0.00265 k=33333 g=20000
prop mat=1 jkn=1000 jks=1000 jf=30.0 ;节理参数为啥不用prop jmat命令呢? 下同
; dam and foundation joint 大坝和基础节理
prop mat=2 d=0.00240 k=16667 g=12500
prop mat=2 jkn=1000 jks=1000 jf=30.0 jcoh=2 jtens=2
; above y=-150: free field (low cohesion) ; joints (no cohesion)
; y=-150上部自由区域（低黏聚力）节理无黏聚力
prop mat=3 d=0.00265 k=33333 g=20000 f=30.0 coh=0.2
prop mat=3 jkn=1000 jks=1000 jf=30.0
; below y=-150: free field (cohesion) ; joints (cohesion)
;y=-150下的自由区域（黏聚力），节理（黏聚力）
prop mat=4 d=0.00265 k=33333 g=20000 f=30.0 coh=2
prop mat=4 jkn=1000 jks=1000 jf=30.0 jcoh=2
; viscous boundaries (equivalent elastic properties) ; joints (cohesion)
;粘结边界（等价弹性参数），节理（黏聚力）
prop mat=5 d=0.00265 k=11680 g=11111
prop mat=5 jkn=1000 jks=1000 jf=30.0 jcoh=2
;
round 0.5
; set minimum edge length 设置最小边长
set edge 8.0
; set minimum contact length 设置最小连接长
set clemin=5.0
;
block -200,-200 -200,100 200,100 200,-200 ;模型块体
;
; structure : gravity dam 结构重力坝
crack -201,0 201,0 ;根据几何形体对块体进行分割
crack -40,-1 -40,101
crack 40,0 -40,100
delete range -200,-40 0,100 ;删除多余部分,形成坝体
delete range 40,200 0,100
;
jregion id=1 -200,-200 -200,0 200,0 200,-200 ;节理生成区域-基础区jset 20,0 800,0 0,0 50,0 ( 50,0) range jreg 1 ;20度
jset -70,0 50,0 50,0 50,0 (-50.99, 16.45) range jreg 1 ;-70度
jset -70,0 50,0 50,0 50,0 (-10.40,-21.98) range jreg 1 ; -70度
;
change mat=1 range -200,200 -200,0
change mat=2 range -40,40 0,100
generate edge 60 range -200,200 -200,0 ;坝下基础网格长度60
generate edge 30 range -40,40 0,100 ;坝体网格长度30
;
; all joints 所有节理
change jmat=1 range -210,210 -210,10
; cohesion below y=-150
change jmat=4 range -210,210 -210,-150
; foundation joint 基础下节理
change jmat=2 range -41,41 -1,1 ang=-5,5
;
; boundary conditions : lateral : x-fixed ; bottom : y-fixed 边界条件
bound xvel 0 range -201,-199 -201,1 ;左侧
bound xvel 0 range 199,201 -201,1 ;右侧 bound yvel=0 range -201,201 -201,-199 ;底部
;为初始域压力设置流体密度
; set fluid density (for insitu domain pressures) fluid dens 0.001 bulk 0.0
set flow off
cra ck
-200,100
-40,100
200,100
删除区域1
删除区域2
-200,0
m at=2
50,0 40,0
200,0
X -FIX
y =-150
ma t=1
X-FIX
-200,-200
;
-FIX
200,-200
; set gravity ; set in-situ stresses (total) grav 0 -9.8
insitu -210,210 -210,0 ygrad= 0.017885, 0, 0.025970 & zgrad 0 0.017885 ywtab=0 ;
hist type=2
hist unbal solve_ratio
hist xdis=-40,100 ydis=-40,100 ;
; damp auto solve rat 1e-5 ;
save dam1.sav ;
; --------------------------------------------------------------------- ;
; --- phase 2 --- water loads and flow ; title
(DAM2) WATER LOADS AND FLOW ;
reset dis jdis hist ;
; flow properties ; above y=-150
prop mat=1 jperm=3.0e8 azero=0.001 ares=0.0005
prop mat=3 jperm=0.0 azero=0.001 ares=0.0005
;
; below y=-150
prop mat=4 jperm=3.0e8 azero=0.001 ares=0.0005
prop mat=5 jperm=0.0 azero=0.001 ares=0.0005
;
; foundation joint
prop mat=2 jperm=3.0e8 azero=0.0002 ares=0.0001
;
; set max. aperture
set caprat=2.0
;
; set lateral and bottom boundary contacts to zero permeability change jmat=3 range -201,-199 -150,1
change jmat=3 range 199,201 -150,1
change jmat=5 range -201,-199 -201,-150
change jmat=5 range 199,201 -201,-150
change jmat=5 range -201,201 -201,-199
;在大坝上游固定水头(蓄水后);下游为零
; fix head upstream of dam (full reservoir) ; downstream is zero
bound pp=0.980 range -201,-39 -1,1
;在大坝上游施加垂直水荷载
; apply vertical water load upstream of dam
bound stress= 0, 0, -0.980 range -201,-39 -1,1
;对大坝施加水平荷载
; apply horizontal load to dam
bound stress=-0.980,0,0 ygrad 0.0098,0,0 range -40.1,-39.9 -0.1,100.1 ;固定侧边界(水平)
; fix lateral boundaries (horizontally)
bound xvel=0 range -201,-199 -201,1 ;左侧
bound xvel=0 range 199,201 -201,1 ;右侧
bound yfree range -201,-199 -199,-1 ;左侧
bound yfree range 199,201 -199,1 ;右侧
;
; monitor contact variables at x=-33.42 y=-30.37
hist solve_ratio type=1
hist xdis(-40,100) ydis(-40,100)
hist sdis(-33.42,-30.37) ndis(-33.42,-30.37)
; (flow rate thru contact)
hist flowrate(-33.42,-30.37)
; (domain head) x=-22.14,y=-26.27
hist pp -22.1 -26.3
;
; fluid properties (nonzero bulk modulus)
fluid dens 0.001 bulk 2000.0
;
; switch on fast flow logic
set flow steady
;
solve rat 1e-5
;
save dam2.sav
;
; --- phase 3 --- dynamic loading: shear wave at base of model
; free-field lateral boundaries
; with constant domain pressures ; no fluid flow by setting:
; fluid dens 0.001 bulk 0.0
; dynamic domain pressure is also possible by setting:
; fluid dens 0.001 bulk 2000.0
; the dynamic solution is less conservative and predicts less movement. ;
title
(DAM3) DYNAMIC LOADING
;
; switch off flow
fluid dens 0.001 bulk 0.0
set flow off
;
; generate free-field (20 nodes) ; cohesion below y=-150 ; fixed bottom ffield gen yrange (-200,0) np 20
ffield change mat=3 cons=3 range y -150,0
ffield change mat=4 cons=3 range y -200,-150
;
; initialize FF stresses (same as insitu stresses)
ffield ini sxx 0 0.017885
ffield ini szz 0 0.017885
ffield ini syy 0 0.025970
ffield ini pp 0 -0.0098
;
; fix bottom
ffield base xvel=0
ffield base yvel=0
;
; cycle with FF not applied to grid
; to verify that FF stresses are in equilibrium
reset time hist
hist n 1
hist ffyd -100 1
hist ffsxx -100 1
cycle 100
save dam2a.sav
;
reset time hist
;
; apply dynamic boundary conditions : viscous boundaries and free-field bound mat=5
bound ff range -201,-199 -201,1
bound ff range 199,201 -201,1
bound xvisc range -201,201 -201,-199
;
; amplitude of shear wave : 0.2 MPa ; freq. = 5 hz.
bound stress 0, 0.4, 0 hist=sin(5.0,10.0) range -201,201 -201,-199
;
; fix y-vel at bottom
bound yvel=0 range -201,201 -201,-199
;
; bound. cond. at base of free-field
ffield base sxy=0.4 hist=sin(5.0,10.0)
ffield base xvisc
ffield base yvel=0.0
;
hist ty=3
hist xvel(-40,100) yvel(-40,100)
hist xdis(-40,100) ydis(-40,100)
hist sdis(-33.42,-30.37) ndis(-33.42,-30.37)
hist sstr(-33.42,-30.37) srat(-33.42,-30.37)
hist ffxdis 0,1 ffxdis 0,2
hist ffxvel 0,1 ffyvel 0,1 ffxvel -200,2 ffyvel -200,2
hist xdis -200 0 ydis -200 0
;
; dynamic damping
damp 0.10 2.0 mass
;
; no mass scaling
mscale off
;
cycle time 1.5
;
label hist 1
X Velocity
label hist 2
Y Velocity
label hist 3
X Displacement
label hist 4
Y Displacement
label hist 5
Joint Shear Displacement
label hist 6
Joint Normal Displacement
save dam3.sav
pl hist 1 2
pl hist 3 4
pl hist 5 6
pl bl disp shear stress
;
ret
5 水泥注浆模拟Cement Grouting Simulation
5.1 问题描述
宾汉姆流体模型作为一个合适的水泥基注浆模型被广泛接受(见,例如Littlejohn (1982), Hassler et al. (1987) and Lombardi (1985)),本模拟展示了在UDEC 中宾汉姆流体的应用.
问题几何体展示了一个包含一个圆柱形钻孔(1.2m 直径)的规则节理岩体横断面.岩体假定服从平面应力状态(σxx = 0.2 MPa and σyy = 0.1 MPa).通过保持钻孔内的指定压力来模拟浆液注入过程,压力增量为2000 Pa,每一步均进行流体条件的检查.岩体水力参数和不连续体及浆液参数如下:
^`
5.2 UDEC 分析
UDEC 模型和初始应力状态见图5.1, 通过FLUID cohw 命令设定宾汉姆流体模型的屈服极限强度临界压力梯度.施加在钻孔上的2000Pa 流体压力增量通过FISH 函数press hole 进行设定,见数据文件5.4 节.
直到钻孔内的压力增加到8000pa 时,才产生不稳定流.相对钻孔的不同接触部位的流体速率历史见图5.2.在孔压10000Pa,12000 Pa,16000 Pa 时的流体显示见图5.3,5.4 图5.5. 相应孔内注浆压力的变化见图5.6 到5.8.
5.3 参考文献
Hassler, L., H. Stille and U. Hakansson. “Si m ulation of Grouting in Jointed Rock,” in Proceedings of the 6th International Congress on Rock Mechanics (Montreal, 1987), Vol. 2, pp. 943-946. Rotterdam: A. A. Balkema, 1987.
Littlejohn, G. S. “Design of C ement Based Grouts,” in Proceedings of the Conference on Grouting
in Geotechnical Engineering (New Orleans, February 1982), pp. 35-49. Wallace Hayward Baker,
Ed. New York: ASCE, 1982.
Lom bardi, G. “The Role of Cohesion in Ce ment Grouting of Rock,” in Transactions of the 15th International Conference on Large Dams (Paris, 1985), pp. 235-261 (1985).
^`
^`
^`
^`
5.4 数据文件
;宾汉姆流体定律下的注浆模拟
; grouting simulation with ; bingham fluid flow logic ;标题
title
GROUTING SIMULATION
;问题的几何模型
; problem geometry
round 0.01 ;倒圆角Example 5.1 GROUT.DAT
0,10 10,10
5,9
5,5
block 0 0 0 10 10 10 10 0 ; 块体模型crack 0 5 10 5 ;水平分割线
crack 5 0 5 10 ;垂直分割线
tunnel 5 5 0.6 12 ;钻孔
jset 0 0 10 0 0 0 2 0 ;水平节理
jset 90 0 10 0 0 0 2 0 ;垂直节理0,5
9,5
10,5
gen edge 1.25 ;网格划分
;——显示设置——
;set back iw
;set color red
;plot bl zones
0,010,0
^`
;材料参数
; material properties
; intact rock 围岩
prop mat 1 k 1e10 g 3e9 d 3000
; discontinuities 不连续体
prop jmat 1 jkn 1e10 jks 1e10 jfric 45
prop jmat 1 jperm 1e8 azero 0.0001 ares 0.00005
; insitu stress and boundary conditons 初始应力和边界条件
insitu stress -2e5 0 -1e5 nodis
;nodis 作用在于禁止计算节理法向方向的位移(例如,节理初始时是封闭的)
bound stress -2e5 0 -1e5
; conditions during solution 求解条件
; create hole 钻孔
delete range 4.5 5.5 4.5 5.5 删除钻孔区域
;plot bl zones stress ;图5.1
; cycle to equilibrium 循环到平衡
hist unbal ;记录不平衡力
solve rat 1e-5 ;求解速率
save gr0.sav ;保存
;plot hist 1
; use steady state fluid logic 采用稳定流体规则
set flow steady
reset hist ;重置记录
hist n=1000 ;设置记录步=1000
hist pp 5,5 ;记录点(5,5)的孔隙压力
hist flowrate 9,5 flowrate 5,9 ;记录点(9,5)和(5,9)的流速
; grout properties (density - cohesion) 浆液参数(密度和粘聚力)
fluid den 1000 coh .1
; pressurize hole in .2e4 increments 压力增量0.2e4
;FISH函数定义钻孔增加的压力
def press_hole ;函数名
loop i (1,8) ;定义循环步
nam_t = ’gr’+string(i)+’.sav’;保存的文件名
p_hole = float(i) * 0.2e4 ;以0.2e4增量的孔压
command ;udec命令开始
pfix pres p_hole range 4.9,5.1 4.9,5.1 ;在钻孔周边施加压力
cy 1000 ;循环1000步
save nam_t ;以定义的文件名保存文件endcommand ; 命令结束
endloop ;循环结束
end ;函数定义完毕
press_hole ;调用定义的fish函数
;plot hist 2 3 vs 1 ;图5.2
^`
;plot bound pp ;
ret ;恢复
re gro.sav
plot bl stress ;图5.1
re gr5.sav
plot bound pp ;图5.6
plot bound flow ;图5.3
re gr6.sav
plot bound flow ;图5.4
plot bound pp ;图5.7
5.5 数据文件中命令的解释
7 隧道涌水模拟Inflow into a Tunnel
7.1 问题描述
在岩土力学工程实践中固体变形和地下水流动在大多数状况下是相互依赖的. (1)地下水对岩体中的固体部分施加压力. (2)而岩体变形控制地下水流的条件(例如地
下水通过处的节理内水力裂缝的改变).忽视两种进程的耦合作用(固体变形和地下水流)有时候会导致对岩体力学扰动(如隧道开挖,抽水,建造大坝)产生影响的不正确预测.
本例分析了节理岩体中地下水位下的一个简单的隧道开挖问题(见图7.1).材料参数(例如:初始节理水力裂缝,节理刚度,整个块体集合的强度)和对系统的扰动(例如, 由于隧道开挖引起总应力和的孔隙水压力的改变)的选择需要完全耦合的加以考虑.
地下水排入隧道开挖面内引起地下水位的重大降落和并形成相对于隧道半径上部一定距离的地下水井面.因此,不饱和流动规则被用到本例中以提供一个精确的自由流动模拟. 7.2 UDEC 分析
模型是二维的,两个连续的正交(水平和垂直)设置的节理建立岩体块,见图7.1.两个方向的节理间距均是 10m.本例假定的参数如下:
注意,水的体积模量是2.0 Gpa.采用低的数值(0.2 GPa)有两个原因:(1) 气体溶解在地下水中增加了地下水的可压缩性,(2)当使用可压缩选项时,较低数值的体积模量可以加速收敛达到数值计算的平衡状态.
初始应力是各项同性的,并且作为岩石密度和重力的函数从顶部到底部而增长.初始的水力孔隙压力被定义为地下水位面,在这里是模型的顶面.模型中固体和地下水的远场边界条件均在初始条件下达到平衡.节理中的有效应力、节理平面内总的块体应力和节理内的初始孔隙水压力的法向分量被初始化为不同的变量.因此，变形和节理的水力裂缝作为初始有效应力和节理法向刚度的函数被计算。

（在本例总，初始节理变形可以被约束并且零应力下的裂缝宽度ao 成为初始应力状态下的裂缝宽度），模型的初始应力状态见图7.2，而初始孔
隙水压力分布见图7.3. 由于几何体不规则需要一些计算步。

在固体和地下水模型中建立和初始化应力状态后，隧道被开挖。

隧道的开挖在固体模型（隧道边界上的总残余应力被降低为零）和地下水模型（隧道边界上的孔隙水压力被降低至零）内引起了一个扰动。

但是，岩体内的力学进程时标和节理内的流体进程时标具有不同的数量级，本例中，模型中固体的反应需要很少的时间尺度，反应是瞬时的。

因此，隧道开挖引起模型的第一步反应是不排水变形：当节理内孔隙水压力是固体模型变形和地下水的体积模量（只是没有流体）的函数，固体模型产生变形。

在不排水时孔隙水压力和模型变形见图7.4 和图7.5。

^`
首先，不排水阶段的模拟达到力学平衡后，第二步是模拟由于地下水进入隧道引起的排水变形（如固结），这个进程是地下水消耗控制下的时间效应，本阶段模拟真实流动时间直到达到稳定流态和固体达到平衡状态。

在UDEC 中有三个选项可以模拟耦合流动的固结（由于地下水进入隧道排水带来的变形的时间效应）效应。

（1）瞬时的可压缩流（SET flow compressible）；（2）稳定流(SET flow steady);和（3）试验算法的瞬态流(SET flow transient)。

如果稳定流计算是单独的路径（例如线性模型），那么在稳定状态条件三个选项必须产生相同的结果。

稳定流选项不必考虑不饱和流——孔隙水压力确定地下水水位，但是，在本问题中，从地下水位上的不饱和单元渗透在极限状态下必须变成零。

这样，三种模式应该在极限状态下收敛。

这可以从图7.7、图7.8、图7.9 的结果加以证实，这些结果从可压缩和瞬态流选项获得，获得的途径是使用稳定流选项（It can be verified from Figures 7.7, 7.8 and 7.9 that the results obtained from the compressible and transient flow options approach those obtained using the steady-flow option.），记录点的位置见图7.6。

需要注意的是, 对于稳定流和瞬态选项的计算时间是一个远快过可压缩流选项一个数量级的计算时间（It must be noted that the calculation times for the steadyflow and the transient options are an order of magnitude
^`
faster than that for the compressible-flow option.）看来这个实验的瞬态逻辑是解决这个问题的最好方法, 因为对于稳定流选项计算时间
几乎是相同的,而它模拟模型的完全瞬态响应。

（稳定流计算仅在稳定状态被更正。

）稳定状态流动的孔隙压力见图7.10，在这个显示中稳定状态的水面被极限状态下的孔
隙水压分布所指示。

节理内的流动速率见图7.11。

^`
7.3 数据文件列表
Example 7.1 INFLOW.DAT ;--------------------------------------------
;
; INFLOW INTO A TUNNEL FROM A FREE SURFACE
; 自由水面水流入隧道
; -------------------------------------------
;
config tflow ;transient fluid flow (for SET flow transient)
title
Inflow into a tunnel
block 0,100 0,200 200,200 200,100
round 0.1
; add the tunnel
tunnel 200. 150. 10 40 ;圆心（200,150）半径10,分段数40
;
; two orthogonal joint sets
;
jset 90 0 200 0 0 0 10 0 5 0
jset 0 0 200 0 0 0 10 0 0 5
;set back iw ;背景色为白色
;set color b ;模型线为黑色
;set box red ;图例色为红色
;set border yellow ;屏幕边界色为黄色
;plot block
gen edge=20 range 0 200 0 200
;plot zones
; material properties
;
prop mat=1 d=2700 bulk=2.e10 g=1.5e10
prop jmat=1 jkn=1e10 jks=1.e10 jcoh=0 jfric=30.0
prop jmat=1 azero=1.e-3 ares=0.5e-3 jperm=3e2
fluid dens=1000 bulk=2e9
set capratio=3. ;最大接触裂隙宽度被限制为最大残余宽度的3倍。

;; Maximum contact hydraulic aperture is limited to r times the maximum residual aperture. (Default is 5.)
; boundary conditions
;
boun stress -5.4e6 0 0 ygrad 2.7e4 0 0 range -0.1 0.1 99.9 200.1 ;左侧应力x=0,y=100~200
boun pp 2e6 pygrad -1.e4 range -0.1 0.1 99.9 200.1 ;左侧水压
boun yvel 0 pp 1e6 range -0.1 200.1 99.9 100.1 ;模型底部限制y方向位移
boun xvel 0 imperm range 199.9 200.1 99.9 200.1 ;模型右侧不可渗透边界，限制X方向位移为0 boun imperm range -0.1 200.1 199.9 200.1 ;模型顶边界为不可渗透边界
不可渗
透边界
,200200,200
不可渗
透边界
str ess
2e6 5.4e6,100200,100
; 模型概念图
; initial conditions
;
gravity 0 -10
insitu stress -5.4e6 0 -5.4e6 ygrad 2.7e4 0 2.7e4 ywtable 200
;
;Insitu ywtable y w ——Domains with y-coordinate below y w are assigned pressure,
; steady state flow only
;
set flow steady
;
; settle under gravity
;
hist unbal
hist ydisp 100,200
;
solve rat 1e-5
;
save inflow1a.sav
pr max
; PRINT max prints out the current maximum unbalanced volume，as well as the actual number of ;iterations performed in the latest fluid step
;---------------------
;
; compressible flow
;
;---------------------
;
rest inflow1a.sav
reset jdisp disp
set flow comp
fluid bulk 200e6 ；0.2Gpa
set flow off
set caprat 3.
;
; remove tunnel block and calculate undrained deformation
;
delete ann 200 150 0 10
reset hist
hist unbal
hist pp 150 150
hist ydisp 200 200
solve rat 1e-5
save influndr.sav
;
; drain groundwater and cycle to steady state
;
reset hist time jdisp disp
set ftime 0.0
hist ncyc 100
hist flowtime
hist unbal
; hist flowrate
hist flow 191 155 flow 191 145 flow 196 158 flow 196 142
hist pp 150 150
hist ydisp 200 200
;
set flow comp
;
fluid bulk 200e6
;
set nfmech 5 ; set nfmech n : maximum number of mechanical steps executed between fluid-flow steps & when nflow is nonzero for SET flow compressible (default: n = 1)
cycle time 1
; cycle time t —t is problem time duration, in seconds, for this increment of cycling
set nfmech 4
cycle time 2
set nfmech 3
cycle time 3
set nfmech 2
cycle time 4
set nfmech 1
cycle time 5
save inflo1c1.sav
set nfmech 1
cycle time 10
save inflo1c2.sav
set nflow 5
; set nflow n——n fluid-flow timesteps are performed during each cycle (default: n =1) for
; SET flow compressible.
cycle time 10
save inflo1c3.sav
set nflow 10
cycle time 20
save inflo1c4.sav
set nflow 50
cycle time 50
save inflo1c5.sav
;-------------------- ----------------------------------------
;
; steady-state flow
;
;-------------------- ----------------------------------------
rest inflow1a.sav
set caprat 3.
;
; remove tunnel block and drain
;
delete ann 200 150 0 10
;
; cycle to steady state
;
reset hist time jdisp disp
set ftime 0.0
hist unbal
; hist flowrate
hist flow 191 155 flow 191 145 flow 196 158 flow 196 142 hist pp 150 150
hist ydisp 200 200
;
set flow steady
cy 5000
save inflow1s.sav
;-------------------- ----------------------------------------
;
; transient flow
;
;-------------------- ----------------------------------------
res influndr.sav
reset hist time jdisp disp
set ftime 0.0
hist ncyc 1
hist flowtime
hist unbal
; hist flowrate
hist flow 191 155 flow 191 145 flow 196 158 flow 196 142 hist pp 150 150
hist ydisp 200 200
;
set flow trans
set maxmech 1000
set voltol 0.0001
;
;
set dtflow 0.01
cycle 10
save inflo1t1.sav
set dtflow 0.03
cycle 8
save inflo1t2.sav
set dtflow 0.1
cycle 10
save inflo1t3.sav
set dtflow 0.3
cycle 10
save inflo1t4.sav
set dtflow 1.
save inflo1t5.sav
return
8 节理边坡渗流分析Flow through a Jointed Rock Slope
8.1 问题描述
岩体节理边坡的稳定性受边坡后面水位所影响。

本例中，水位阶段升高直到边坡失稳，当节理内的流体压力增加（节理内的有效正应力降低）以至于边坡面节理的抗剪强度达到极限时，边坡产生破坏。

包含规则节理边坡的几何模型见图8.1。

水位在坡址按四个阶段逐渐升高，四个阶段的水位高程分别是6m、8m、9m 和 10m。

在每一阶段均进行一个稳定流体的分析。

假定下面的材料参数为节理边坡的岩体参数。