来自 UDEC 6.0 Manual →EXAMPLE APPLICATIONS →Seismic-Induced Groundfall

1.问题描述

此示例用来模拟地震引起的隧道顶部块体掉落问题，以说明使用UDEC分析此类动力问题的方法。

涉及的问题：(1)交错节理组的设置； (2)粘性边界(自由场边界)的施加；(3)应力波的输入

模型如下图所示，简化为平面应变问题。 中部开挖的部分高5 m、宽10 m。

假定有两组节理与分析平面相交：一个方向为 45º，另一个方向为 -9º。 两组的关节间距均为 5 m。 出于演示目的，还向开挖部分顶部的砌块添加了一条垂直的“人造” 节理，以增强动力作用下的不稳定性。

根据完整岩石的实验室测试的平均值，假定了岩石的几种材料特性并假定块体为弹性体：

选择的节理模型为库仑滑移节理模型，并选择节理的力学参数如下：

地应力状态认为是各向同性，估计值为24 MPa（相当于约800m深度处岩层产生的垂直荷载）。

2.UDEC分析

UDEC 建模分为三个阶段。

首先，在地应力作用下，对没有开挖的模型运行至平衡状态。 然后开挖并将模型运行到平衡状态。 下图一幅图显示了此阶段开挖部分周围的应力分布。

最后，分析在两个具有不同峰值速度的地震作用下模型的响应。

对于所有地震这种动力问题的模拟，一般在研究域的外围引入粘性边界（或设置自由场边界），以消除应力波反射，从而模拟处于大地的无限大区域。 地震事件由施加在模型顶部Y向的正弦应力波表示，施加的应力波叠加在现有的地应力上。

在第一个模拟中，施加了1.25 MPa（applide stress=1.25 × cos( 2π 100t) ）的峰值应力。 注意，由于只在模型顶部粘性边界上施加应力波，“有效”的施加应力为1.25 MPa的一半（ 0.625 MPa）。 下图显示了0.02秒后的开挖的顶板中的应力分布。

监测了两个点的位移，点1位于开挖的左上角； 点2位于顶板块体的右上角。 监测点的位移与时间的关系图（如下图）基本上显示了弹性响应。

将估计的施加速度与模型顶部的计算速度进行比较。 以下方程式可用于估计施加的速度：

使用该方程式，发现施加的最大速度约为 0.04 m/sec。 下图显示了Y方向约 0.06 m/sec 的峰值速度。 估计速度与模拟测得速度之间产生差异的原因是由于使用了完整岩石的模量而不是等效变形模量（考虑了节理变形）导致的。

在第二个模拟中，施加了峰值应力为 12.5 MPa（applied stress = 12.5 × cos(2π 100t) ），有效应力为 6.25MPa的应力波。 下图显示了0.02秒后的开挖顶板中的应力分布。 该图显示了顶板块体没有受力，表明已经松动。

监测点位移与时间的关系图也表明顶板块体已经松动并且正在下降。 后面的几张图显示了正在下降的块体。

同样，顶部的预测速度（根据上述方程式）为 0.4 m/sec。 模型计算的速度如下图。 预测的速度和模拟测量的速度之间产生差异的原因同上。

命令末尾使用MOVIE命令生成图片，这将生成一系列名为 “SEISMIC.DCX” 的图片文件

3.命令说明

查看完整的命令流请点击下方的文件：
seismic.dat

只说一些重点命令的含义，其他的应该一看就明白。

group zone 'mat1' ;单元分组，命名为“mat1”

;赋予块体弹性体属性
zone model elastic density 0.003 bulk 3.906E4 shear 3.178E4 range group 'mat1'

group joint 'jmat1' ;节理分组并命名

;赋予库伦滑移节理模型
joint model area jks 2E4 jkn 2E4 jfriction 30 range group 'jmat1'

; new contact default，设置新接触的默认属性
set jcondf joint model area jks=2E4 jkn=2E4 jfriction=30

;施加边界条件和初始条件，方式为Y方向的应力梯度
boundary stress -24.0,0.0,-24.0 xgrad 0.0,0.0,0.0 ygrad -0.03,0.0,-0.03
insitu stress -24.0,0.0,-24.0 xgrad 0.0,0.0,0.0 ygrad -0.03,0.0,-0.03

;固定模型底部
boundary yvelocity 0 range -26,26 -21,-19
set gravity=0 -10

;设置监测点记录位移并求解平衡
history xdisplace 0.0,7.0
history ydisplace 0.0,7.0
solve ratio 1.0E-5

后面为开挖的命令

delete range -5,5 -2.5,2.5 ;开挖，即删除中间块体
solve ratio 1.0E-5

模型四周设置黏性边界条件(自由场边界)

boundary xvisc ff_bulk=39060.0 ff_shear=31780.0 ff_density=0.0030 range &
 -26,-23 -21,21
boundary xvisc ff_bulk=39060.0 ff_shear=31780.0 ff_density=0.0030 range &
 23,26 -21,21
boundary yvisc ff_bulk=39060.0 ff_shear=31780.0 ff_density=0.0030 range &
 -26,26 -21,-19
boundary yvisc ff_bulk=39060.0 ff_shear=31780.0 ff_density=0.0030 range &
 -26,26 19,21

再在模型顶部施加应力波，yhistory为乘子。可作为乘号理解：applide stress=1.25 × cos( 2π 100t)

boundary stress 0.0,0.0,-1.25 yhistory= cosine(100,1.95E-2) range -26,26 19,21
;cosine(100,1.95e-2)中 100 指频率，1.95e-2 指施加的时间

如果想生成顶板掉落的动画。可以使用下面的命令

restore 'seismic4.sav'
;make a movie of the groundfall
;
wind -12 12 -12 12
set ovtol 0.05
plot block vel max 5.0 blue stress max 50 ;设置输出位移和应力图
movie on
movie file = seismic.dcx ;设置文件名称
movie step 1000 ;设置每隔1000步记录一张图
step 40000

生成动画的movie命令不太好用，图片格式只能指定为.dcx。。输出的图也不清晰

如果要生成动画，可以尝试在循环体中使用set plot 等命令，不过一般我也没输出过图片，以后等我搞清楚了再写吧。