通用离散元用户指导
（U D E C 3.1）
北京科技大学
2004.9
目 录
1 引 言 ............................................................ 1 1.1 总 论 ........................................................... 1 1.2 与其他方法的比较 ................................................ 2 1.3 一般特性 ........................................................ 2 1.4 应用领域 ........................................................ 3 2 开始启动........................................................... 4 2.1 安装和启动程序 .................................................. 4 2.1.7 内存赋值 .................................................... 4 2.1.9 运行UDEC .................................................... 5 2.1.10 安装测试程序 ............................................... 5 2.2 简单演示－通用命令的应用 ......................................... 5 2.3 概念与术语 ...................................................... 6 2.4 UDEC模型：初始块体的划分 ........................................ 8 2.5 命令语法 ........................................................ 9 2.6 UDEC应用基础 .................................................... 10 2.6.1 块体划分 .................................................... 10 2.6.2 指定材料模型 ................................................ 16 2.6.2.1 块体模型 ................................................ 16 2.6.2.2 节理模型 ................................................ 17 2.6.3 施加边界条件和初始条件 ....................................... 19 2.6.4 迭代为初始平衡 .............................................. 21 2.6.5 进行改变和分析 .............................................. 24 2.6.6 保存或恢复计算状态 ........................................... 25 2.6.7 简单分析的总结 .............................................. 25 2.8 系统单位 ........................................................ 26 3 用UDEC求解问题 .................................................... 27 3.1 一般性研究 ...................................................... 27 3.1.1 第1步：定义分析模型的对象 ................................... 28 3.1.2 第2步：产生物理系统的概念图形 ............................... 28 3.1.3 第3步：建造和运行简单的理想模型 ............................. 28 3.1.4 第4步：综合特定问题的数据 ................................... 29 3.1.5 第5步：准备一系列详细的运行模型 ............................. 29 3.1.6 第6步：进行模型计算 ......................................... 29
－1－
3.1.7 第7步：提供结果和解释 ....................................... 30 3.2 产生模型 ....................................................... 30 3.2.1 确定UDEC模型合适的计算范围 .................................. 30 3.2.2 产生节理 .................................................... 32
3.2.2.1 统计节理组生成器 ....................................... 32 3.2.2.2 VORONOI多边形生成器 ................................... 34 3.2.2.3 例子 .................................................. 34 3.2.3 产生内部边界形状............................................. 35 3.3 变形块体和刚体的选择 ............................................ 38 3.4 边界条件 ....................................................... 42 3.4.1 应力边界 .................................................... 42
3.4.1.1 施加应力梯度 ........................................... 43 3.4.1.2 改变边界应力 ........................................... 44 3.4.1.3 打印和绘图 ............................................ 44 3.4.1.4 提示和建议 ............................................ 45 3.4.2 位移边界 .................................................... 46 3.4.3 真实边界－选择合理类型 ....................................... 46 3.4.4 人工边界 .................................................... 46
3.4.4.1 对称轴 ................................................ 46 3.4.4.2 截取边界 .............................................. 46 3.4.4.3 边界元边界 ............................................ 49
3.5 初始条件 ....................................................... 50 3.5.1 在均匀介质中的均匀应力：无重力 ............................... 50 3.5.2 无节理介质中具有梯度变化的应力：均匀材料 ...................... 51 3.5.3 无节理介质中具有梯度变化的应力：非均匀材料 .................... 51 3.5.4 具有非均匀单元的密实模型 ..................................... 52 3.5.5 随模型变化的初始应力 ......................................... 53 3.5.6 节理化介质的应力............................................. 54 3.5.7 绘制应力等值线图............................................. 55 3.6 加载与施工模拟 ................................................. 57 3.7 选择本构模型 ................................................... 62 3.7.1 变形块体材料模型............................................. 63 3.7.2 节理材料模型 ................................................ 64 3.7.3 合理模型的选择 .............................................. 65 3.8 材料性质 ....................................................... 71
－2－
3.8.1 岩块性质 .................................................... 71
3.8.1.1 质量密度 .............................................. 71 3.8.1.2 基本变形性质 ........................................... 71 3.8.1.3 基本强度性质 ........................................... 72 3.8.1.4 峰后效应 .............................................. 73 3.8.1.5 现场性质参数的外延 ..................................... 77 3.8.2 节理性质 .................................................... 80 3.9 提示和建议 ..................................................... 81 3.9.1 节理几何形状的选择 ........................................... 81 3.9.2 设计模型 .................................................... 81 3.9.3 检查模型运行时间............................................. 82 3.9.4 对允许时间的影响............................................. 82 3.9.5 单元密度的考虑 .............................................. 83 3.9.6 检查模型响应 ................................................ 83 3.9.7 检查块体接触 ................................................ 83 3.9.8 应用体积模量和剪切模量 ....................................... 83 3.9.9 选择阻尼 .................................................... 84 3.9.10 给块体和节理模型指定模型和赋值 .............................. 84 3.9.11 避免圆角误差 ............................................... 85 3.9.12 接触嵌入 ................................................... 85 3.9.13 非联结块体 ................................................. 86 3.9.14 初始化变量 ................................................. 86 3.9.15 确定坍塌荷载 ............................................... 86 3.9.16 确定安全系数 ............................................... 86 3.10 解 释 ......................................................... 88 3.10.1 不平衡力 ................................................... 88 3.10.2 块体/网格结点的速度 ......................................... 88 3.10.3 块体破坏的塑性指标 .......................................... 89 3.11 模拟方法 ...................................................... 89 3.11.1 有限数据系统模拟 ............................................ 89 3.11.2 混沌系统的模拟 ............................................. 90 3.11.3 局部化、物理的不稳定性和应力路径 ............................ 91
－3－
1 引 言
1.1 总 论
通用离散元程序（UDEC，Universal Distinct Element Code）是一个处理不连续介质的二维离散元程序。UDEC用于模拟非连续介质（如岩体中的节理裂隙等）承受静载或动载作用下的响应。非连续介质是通过离散的块体集合体加以表示。不连续面处理为块体间的边界面，允许块体沿不连续面发生较大位移和转动。块体可以是刚体或变形体。变形块体被划分成有限个单元网格，且每一单元根据给定的“应力－应变”准则，表现为线性或非线性特性。不连续面发生法向和切向的相对运动也由线性或非线性“力－位移”的关系控制。在UDEC中，为完整块体和不连续面开发了几种材料特性模型，用来模拟不连续地质界面可能显现的典型特性。UDEC是基于“拉格朗日”算法很好地模拟块体系统的变形和大位移。
UDEC包含了功能强大的程序语言FISH函数。借助于FISH函数，用户可以编写自己的功能函数，扩展UDEC的应用功能。FISH函数为简化分析，适应特殊要求的UDEC的用户，提供了一个强有力的工具。
UDEC采用的离散单元法理论由Cundall（1971）首次提出，至今已经过了20多年的发展。在1985年，Cundall博士和Itasca公司在IBM系列兼容微机上开发了UDEC工程计算应用程序。该软件为建立数以千块模型的高速计算而设计。基于浮点运算速度的优势和低成本的内置RAM，用UDEC程序可大大地提高了计算大规模问题的能力。例如，在具有4MB RAM的微机上，UDEC能够求解2500个刚体（或1000个具有8个自由度变形体）的模型。该模型的求解速度大约为每分钟200次。在RAM确定的情况下，其计算速度是与模型的块体数量成线性关系。
对于典型的模型，约1500个刚体（或500个变形体）或更少，在UDEC中采用的显式解法，大约需要2000～4000计算步可以获得问题的解。例如，一个500个变形块体的模型，计算4000步大约需要6 min。因此，典型的工程问题用UDEC计算仅需几十分钟或几个小时。
UDEC是一个命令驱动（而不是菜单驱动）的计算程序。尽管菜单驱动程序易于初次学习，但在UDEC中所提供的命令驱动结构具有如下优点：
1、输入的“语言”是基于可识别的文字命令，使你易于识别每一个命令的作用（例如BOUNDARY命令，是指施加模型的边界条件）。
2、工程模拟通常是按照系列施工顺序构成 ―― 即，构造原岩应力，施加作用的荷载、开挖隧道、安装支护等。一系列（从文件或键盘上）输入命令完全对应于实际的施工顺序。
3、根据文本编辑器，很容易对UDEC数据文件进行编辑和修改。几个数据文件能
－1－
相互连接，进行多个问题求解，这对于进行参数的灵敏度分析是十分有用的。
4、命令驱动结构允许用户开发前后处理程序，控制UDEC必要的输入/ 输出。用户可以为一系列UDEC的模拟，编写节理模拟函数，产生特定的节理结构。可采用FISH程序语言，并插入到输入的文件中，使计算很容易实现。
1.2 与其他方法的比较
对于UDEC程序，一个共同的问题是，UDEC是一个有限元程序还是离散元程序？他们的主要区别是什么？UDEC程序与其他程序有何关系？为回答上述问题，下面将给予解释。
许多有限元、边界单元和拉格朗日有限差分程序都具有“界面单元”或“节理单元”，使程序能够模拟问题中的不连续面，扩大程序的应用范围。然而，他们的公式在一个或多个方面通常受到限制：首先，当考虑很多相互切割的节理就可能打乱系统的逻辑关系；其次，不可能自动识别新的接触面进行自动考虑；第三，计算公式可能有小位移和无转动条件限制，所以通常适用连续介质的程序。
术语“离散单元法”（Discrete element method）意味着： （a）允许离散块体发生有限的位移和转动，包括完全脱离； （b）在计算过程中，自动识别新的接触面。
在不连续介质中，如果没有第一个属性，程序不可能产生某些重要的机理。如果没有第二个特性，程序将限制在事先已知的相互作用的有限块体数。离散元法（Distinct element method）是由Cundall和Strack（1979）采用变形接触和显式、时间域的初始运动方程（而不是变换，块体方程）提出的特殊的离散单元法程序。
离散单元法的计算机程序主要有以下四类：
1、Distinct Element Programs － 该类程序采用显式时间步直接进行运动方程的求解。块体可以是刚体或变形体（通过细分成单元）；接触面是可变形的。UDEC就属此类。
2、Modal Methods － 该类方法类似于刚体离散单元法，但对于变形体采用模型叠加技术。
3、Discontinuous Deformation Analysis － 接触是刚体，块体可以是刚体或变形体。通过迭代算法可以获得非嵌入条件；块体变形性基于应变模型的叠加。
4、Momentum-Exchange Methods－接触面和块体都是刚体：块体接触面在瞬时碰撞的过程中惯性矩发生交换，可以表征滑动和摩擦特性。
1.3 一般特性
UDEC主要用于岩石边坡的渐进破坏研究及评价岩体的节理、裂隙、断层、层面对地下工程和岩石基础的影响。UDEC对研究不连续特征的潜在破坏模型是十分理想的工具。
－2－
当地质结构特征明显且易于明确描述的情况适宜使用该程序进行分析。UDEC开发了人工或自动节理生成器，用以模拟产生岩体中一组或多组不连续面。在模型中，可以产生变化范围较大的节理模式。屏幕绘图工具允许用户随时观看节理模型。在最后确定所选择的节理模型前，能容易进行调整与修改。
也可以获得不同的节理材料特性。基本模型是指定节理弹性刚度、摩擦角、粘结力、张拉强度和剪胀特性的库仑滑动准则。对该模型的改进包括随着位移的发展而粘结力和张拉强度的降低弱化。在此还可获得一个比较复杂的模拟连续屈服的节理模型，用以模拟弱化为累积塑性剪切位移函数的连续变化特性。作为一个选择模型，还可获得Barton－Bandis节理模型。节理模型和性质参数也可分别赋给单一节理或节理组。应当注意，即使地质图上所显示的节理为直线段，节理的几何粗糙度也可以通过节理材料模型加以表征。
UDEC的块体可以是刚体或变形体。对于变形块体，开发了包括用于开挖模拟的空模型（null）、应变硬化/软化的剪切屈服破坏模型以及非线性不可逆的剪切破坏和压缩模型。因此，块体能被用来模拟回填、土体介质以及完整岩石。
UDEC的基本公式假设为二维平面应变模型。此条件涉及断面保持为定值，并在平行于该断面的平面上作用荷载的无限长结构。所以，非连续面也被假设为平面特性。另外，UDEC提供了一个平面应力问题的选择。对于平面应变分析，如果在垂直于平面方向的应力?zz，为最大或最小主应力，在垂直于平面方向， 块体可能出现塑性屈服，
UDEC的显式求解算法允许进行动态或静态分析。对于动态计算，用户指定的速度或应力波可作为外部的边界条件或者内部激励直接输入到模型中。一个简单的动态波型库也可以获取。UDEC为动力分析设计了自由边界条件。
在静态分析中，包括了应力（力）和固定位移（速度为零）两种边界条件。边界条件在不同的位置可以是不同的。同时，在UDEC中还可以获得边界元边界，用于模拟无限弹性边界。也可以获得半平面解用来描述自由面效应。
UDEC还能够模拟通过模型中的孔隙和不连续面的流体流动。在此认为块体是不可渗透的。岩体的渗透率取决于节理的力学变形。也能够进行力学－流体全耦合分析。反过来，节理水压也将影响力学特性。流体被处理为平行板的粘性流。
程序中的结构单元可用于模拟岩体加固和工程表面支护。加固包括端部锚固、全长锚索和锚杆。表面支护模拟诸如喷射混凝土、混凝土衬砌和其他形式的隧道支护。 UDEC包含一个强有力的程序语言，FISH，能够使用户定义新的变量和函数。FISH是一个编辑器。通过UDEC数据文件进入程序被翻译并储存在内存中。
1.4 应用领域
UDEC最初是为节理岩石边坡的稳定性分析开发的。对于块体不连续公式和运动方
－3－
程（包括惯性项）采用显式时间步求解方法，便于块状岩体边坡的渐进破坏分析和大变形运动研究。
UDEC常用于采矿工程，已经进行了深部地下采矿洞室的静态与动态分析。洞室围岩破坏诱发的断裂、滑移是用UDEC分析研究的实例之一。通过在模型的边界施加动应力或速度波研究爆破影响。地震诱发的断层滑移也通过采用连续屈服节理模型进行了研究。结构单元已经用于模拟全长岩锚和喷射混凝土的各种岩体加固系统。
UDEC还应用于地下结构和深部高辐射废料的储存研究领域。通过应用热模型，UDEC已经应用于模拟与核废料相关的热荷载效应。
UDEC在作为一个计算设计工具，仍受到一定的限制。然而，程序较适用于研究节理效应的潜在破坏机理。节理岩体特性是一个“有限数据系统”－即，在很大程度上内部结构和应力状态是未知和不可知的。因此，建立一个完备的节理模型是不可能的。而且，UDEC是一个二维程序，除了特殊情况外，不可能表征具有三维结构的节理模型。不过，应用UDEC程序，可以从现象学的角度研究节理岩体地下工程开挖响应。该方法可加深岩石力学设计中各种不同现象的相互影响的理解。采用这种方法，工程师能够通过识别地下工程可能产生不可接受的变形或加载导致的破坏机理，从而揭示工程所潜在的诸多问题。
值得注意的是，UDEC程序对于模拟颗粒流动或动态分析火山喷发是不适宜的。对于该类研究，可以采用PFC2D程序。
2 开始启动
2.1 安装和启动程序
本节为首次使用UDEC的用户提供指导。如果你熟悉该程序仅仅是偶尔使用，你会发现本节尤其是2.6节对于改变你原有印象是有帮助的。UDEC程序共有65个主命令，有接近400个关键词。 2.1.7 内存赋值
UDEC自动调节内存大小达到8MB。可以通过下列命令查询、改变内存值： Uedc m Uedc 14 Print mem
如果更多的内存可以获得，其内存能够通过应用环绕磁盘文件获得额外内存。表2.2给出了最大块体数与所需内存的关系。
－4－
表2.年2 RAM与最大块体单元
RAM（MB） 2 4 8 16 最大刚性块体数 400 2500 7500 15000 最大变形块体数* 300 1000 3000 7000 * 假设每块体8个自由度。块体最大数随自由度的增加而减少。 2.1.9 运行UDEC call file.dat 2.1.10 安装测试程序
有三个简单的数据文件，test1.dat 、test2.dat、test3.dat用于程序测试。
2.2 简单演示－通用命令的应用
Block (0,0) (0,20) (20,20) (20,0) ; 产生一个块体 plot block ; 显示该块体 划分初始块体成小块体。 Crack (0,2) (20,8) Crack (5,3) (5,20) Crack (5,12) (20,18)
固定最下和最左块体，使之不可移动的命令如下： fix rang 0,20 0,5 fix rang 0,5 0,20
该命令固定形心处在0 prop mat=1 dens=2000 prop jmat=1 jkn=1.33e7 jks=1.33e7 jfric=20.0 对于该问题，所有的块体密度被指定为2000kg/m3。所有的节理切向刚度和法向刚度分别被指定为1.33e7，节理面的摩擦角为20o。下面将会发现，不同节理和块体可以赋予不同性质参数。 其次，在x和y方向的重力加速度可以通过如下命令予以赋值： set gravity 0，-10.0 为了吸收振动能量，引入阻尼命令 damp local －5－ 上述命令是UDEC的缺省阻尼条件，因此，DAMP local 实际上并不需要。我们在此仅仅是强调这是静态分析。 对于该点，问题很容易被执行。正如在后面能看到的，通过观测特定点的岩体运动有助于进行工程特性判断。在该问题中，我们监测模型右角点y方向的速度，记录该运动所采用的命令是： hist yvel (20,20) type 1 关键词type 是在屏幕上以指定的间隔显示其值。 Step 100 ; 迭代次数 在计算过程中，当前的循环数，计算时间、最大不平衡力，在点（20，20）的y方向速度以每间隔10次显示在屏幕上。 Plot hist 1 Title HEAD> ASIMPLE SLOPE STABILITY EXAMPLE EQUILIBRIUM STAGE Plot block Save slope.sav 通过最左边的块体来研究边坡的特性： delete range 0,5 0,20 命令delete 将删除形心位于0 Cycle 1000 Plot block velocity 2.3 概念与术语 UDEC所涉及的一些术语大部分与其他应力分析程序类似。在UDEC模型中采用一些特殊的术语来描述不连续面特征。按分类给出如下的基本定义。图2.6给出所给出的术语定义。 UDEC MODEL －UDEC模型：是用户为模拟实际的物理模型建立的。当称之为UDEC模型，就意味着为数值求解定义的求解条件的一系列命令。 BLOCK － 块体：是离散单元计算的基本单元体。通过切割一个块体成多个小的块体产生UDEC模型。每一块可能与其他块体分离或通过界面力与其他块体相互作用的独立块体。 CONTACT － 接触：每一块体通过点接触与相邻块体连接。接触可以认为是施加外力到每一块体的边界条件。 DISCONTINUITY － 不连续面：是分离岩体成离散部分的地质特征。不连续面包括岩体中的节理、裂隙、断层和其他不连续特征。 －6－ 图2.6 UDEC模型的例子 DOMAIN － 区域：是指块体间的空洞或空间。DOMAIN在UDEC模型中被处理为实体。每一个DOMAIN是由两个或多个接触面确定的封闭区域。外DOMAIN是指围绕UDEC模型的区域。 ZONE－单元：是由有限个单元组成的变形块体。在每一单元计算力学变化和温度变化。在UDEC采用三角单元。 GRIDPOINT －结点（或节点）：节点包括有限单元的角点。每一单元涉及三个节点。一对x和y坐标定义每一个节点。因此确定了有限单元的精确位置。另一节点的术语是node。 MODEL BOUNDARY – 模型的边界：是一个UDEC模型的周边。边界与模型的外区域一致。内边界（即模型内的孔洞）也是模型的边界。每一内边界通过内区域定义。 BOUNDARY CONDITION － 边界条件：是约束或控制模型的边界（即对于力学问题固定位移或外力）。 INITIAL CONDITION － 初始条件：模型受扰动（开挖）或加载（支护）前的原岩应力状态。 NULL BLOCK－开挖块：表示模型中的空域（即材料不存在）。空块体可在后来加以改变，例如，模拟回填（但一旦块体从模型中删除，就不可能恢复）。 STRUCTURAL ELEMENT － 结构单元：用来表征结构（如隧道衬砌、锚杆和锚索）与岩体的相互作用的一维单元。结构单元也可以具有材料非线性。在大应变模型中可以出现几何非线性。 STEP －求解（或迭代）：尽管一个大的问题需要上万次计算才能达到稳定解，但一般典型问题的求解需要2000～4000次循环，可以获得系统的平衡或稳态流。 －7－ STATIC SOLUTION － 静态解：当模型中动能的变化速率接近可以忽略的情况时，UDEC就认为达到了静态或拟静态解。 UNBALANCED FORCE － 不平衡力：表示当静力分析中的力所处于的不平衡状态（即节理开始滑动或塑性流动）。 DYNAMIC SOLUTION － 动力解：尽管系统的缺省为静态求解过程，但可以进行动态分析。对于动态分析，全运动方程（包括惯性项）可被求解。动能的消耗产生直接影响。 2.4 UDEC模型：初始块体的划分 UDEC模型首先生成整个计算范围的单一块体。然后，通过用地质结构特征（如断层、节理裂隙等）和工程结构（如地下洞室与隧道等）作为边界，切割该块体成小的块体来考虑模型特征。 模型的所有块体都是通过块体质心和角点的坐标（x和y）确定。块体接触面以及变形块体的节点也通过他们的坐标位置确定。产生模型包括由端点坐标（x，y坐标）所定义的线段（splits）切割模型块体。 UDEC模型所有的条目（块体、角点、接触面、空区、节点和单元）都是通过位于主数组中的地址编号，由UDEC自动的、唯一识别和确定。这些编码号也可以用作特殊的单元。编码系统并不是顺序编码，所以用户必须通过绘图或打印加以识别。 例如，图2.7说明一个UDEC模型块体在x和y方向皆为10个单位（比如10m）。模型通过一水平不连续面（x=0,y=5 to x=10,y=5）划分成两个块体。这两个块体具有编号为2和118。块体通过位于块体角点之间的接触面连接。接触号是223和260。内部区域由两个接触面产生并由内部区域号297所识别。 显示在图2.7中的模型由列在下表中的命令产生。 block (0,0) (0,10) (10,10) (10,0) Crack (0,5) (10,5) plot hold block num cont num dnum 两块体的每一块可通过产生有限单元形成变形体。图2.8给出了上部块体划分为8个单元和下部块体划分为4个单元的单元号。在两块体间产生了一个新的接触面（序号为606）。位于块体棱上的任何节点总会产生。新的接触606对应于上部块体的棱产生的新的角点。 Gen quad 11,6 range 0,10 0,5 Gen quad 10 Plot hold zone num cont num －8－ 图2.7 UDEC模型块体被划分成两个刚体 图2.8 包含两个变形块体的UDEC模型 2.5 命令语法 UDEC中所有命令都是面向单词，并由主要命令单词和随后的一个或多个关键词或值构成。某些命令接受开关，即关键词修改命令的作用。每一命令都具有下列格式： COMMAND keyword value … < keyword value … > 在此，位于 < >内的参数为选择参数。而位于（ ）表示可以该参数为任意给定的值。 －9－ 命令可依次写在命令行中。可能你已注意到，命令关键词仅前面几个字母为黑体。实际输入时仅输入这些黑体字母就可由系统识别。 2.6 UDEC应用基础 UDEC是基于命令驱动格式。命令单词控制程序的运行。本节将提供给新用户一些基本命令。为了建立一个UDEC模型进行模拟，必须考虑计算问题的基本成分： （1）由切割产生几何问题，由此建立离散单元块体模型； （2）本构特性和材料性质； （3）边界条件和初始条件。 块体模型定义问题的几何体。本构特性及所涉及的材料参数反映模型在受到干扰后的力学响应。边界条件和初始应力定义了原岩状态，即在未受到扰动（开挖、支护、爆破等）前的应力和位移状态。 在UDEC中定义了这些条件后，可以进行改变（即开挖材料或改变边界条件），从而计算产生模型响应。像UDEC一类显式求解技术所获得问题的实际解与传统的隐式求解方法的结果有所不同。UDEC采用的是显式时间步求解代数方程，其解是在一系列计算迭代后才获得。在UDEC中计算迭代步数可以通过用户控制。用户必须确定所进行的求解迭代步数是否达到了实际问题的解。 图2.9给出了采用UDCE进行静态分析的求解一般过程。由于这求解程序符合实际物理模型的生产工序和实际条件，因此其计算过程是方便的。采用上述过程进行简单的应力分析的UDEC基本命令将叙述如下。 2.6.1 块体划分 UDEC模型是通过切割初始的UDEC块体成小的块体代表模型的实际边界。采用下述命令，建立模型块体。 Block x1,y1 x2,y2 x3,y3 在此，（x1，y1），（x2，y2），（x3，y3）… 是定义块体角点的坐标对。角点必须按顺时针方向排列。角点应当与物理模型的边界条件一致。块体有很多角点，但通常从4角点块体做起。 在UDEC中所有块体都有“圆角”，其目的在于避免块体悬挂在有棱角的节点上。由于块体悬挂引起应力集中。然而，圆角值存在与模型有关的上限值。对于变形块体，最大圆角长度应当不超过块体平均棱长的1％。圆角长度可以如下命令加以改变： round d 在此，d是圆角距离（缺省值是d=0.5）。模型中的所有圆角长度都是相同的。 建议在block命令前指定圆角长度。在block 命令后，键入plot block命令，就能够显示圆角的效果。 －10－ 建立模型 （1）生成模型块体，切割块体产生计算模型的几何体； （2）定义本构模型和材料参数； （3）指定边界条件和初始条件。 迭代计算使之平衡（模拟未扰动前状态） 检查模型响应 进行变化分析 例如： ? 地下开挖 ? 改变边界条件 迭代计算 重复其他变化进行响应分析 检查模型响应 图2.9 静态分析的一般过程 UDEC有几个命令用于产生计算模型的几何体。生成地质结构（即节理）的两个主要命令如下： Crack Jset Crack 命令用于产生块体中单一直线特征的裂缝。裂缝由端点坐标（x1，y1）和（x2，y2）所确定。 Jset 命令则是自动节理组生成器。根据所给定的特征参数（即倾角、迹长、岩桥长度、间距和空间位置）产生一组裂缝。 Crack 和Jset 两个命令用于产生UDEC块体中的地质不连续面，即节理并不一定完全将岩块切割成分离两个块体。然而，UDEC需要连续断裂（即所有断裂都必须切割块体）。由crack或Jset命令所产生的不连续面位置将被储存。为产生块体内连续断裂，可采用crack命令产生的断裂。刚性块体在计算过程中，或变形块体在单元划分时，没有连接形成完整块体的裂缝将被删除。 －11－ 下面的例子说明用Crack和Jset命令切割块体。这两个命令的详细描述将在命令表的第1节给出。节理生成器将在3.2.2节给予详细解释。 Example 2.4 产生简单的UDEC模型 round 0.1 block (0,0) (0,10) (10,10) (10,0) Crack (0,5) (10,5) 在这个最简单的模型中，切割块体涉及选择位置和指定裂缝。通过键入这些命令，就可以产生10×10个单位的块体，然后劈裂成两个块体。CRACK命令产生一个连续、水平的贯通模型的裂缝。注意圆角的长度指定为0.1。 通过键入如下命令，可以产生一个槽口： crack 2.5,10 5.0, 7.5 crack 5.0, 7.5 7.5,10 通过键入以下命令，就显示出包括块体地址号的块体图形。 Plot block num 通过采用DELETE 命令，能从模型中删除一个块体。例如，为了删除槽口块体，键入如下命令： delete range block 368 或 delete range 4.5 ,5.5 8,10 4.5< x <5.5 和8 < y < 10的范围必须包含被删除块体的形心。注意，当对模型进行某些操作时，采用坐标范围是比较明智的。与问题相关的地址号有时发生变化。 对于JSET命令的参数需要4组数据对参数值。每一数据对中的第一个值是均值，而第二个是对应于均值的最大均方差（相对于均匀概率分布形式）。第一组数据对是节理迹线与x-坐标轴的正方向的夹角。第二对数据节理迹线长度；第三组数据是不连续节理的岩桥长度；第四组数据是节理间距。还有一些选择参数，可以用于产生一组比较复杂的节理模式。JSET命令的一般应用在3.2.2节加以讨论。采用JSET命令产生两组节理组的应用在例2.5中得到说明。 例2.5 两组连续节理组的产生 new round 0.01 block (0,0) (0,20) (20,20) (20,0) jset (45,0) (5,0.5) (0.5,0) (2,0) jset (-10,0) (5,0.5) (0.2,0) (1.5,0) 在上例所生成的节理图如图2.11所示。第一个JSET命令产生一组与x轴方向夹角为45o具有间距为2个单位的连续节理。第二个JSET命令产生与x方向夹角为 -10o、间 －12－ 距为1.5个单位的连续节理。圆角长度的选择可能影响节理组的产状。节理的位置可能由于当棱长小于2倍的圆角长度不能产生块体而可能改变节理的位置。 如果在JSET命令前，增大圆角的长度（比如说0.1），模型中某些节理的位置将发生改变。用JSET命令产生节理可能涉及某些试错法。第3节给出了进行产生节理过程的建议。 图2.11 产生的两组连续节理模型 当产生大小悬殊的块体时，建议从模型中删除较小块体，以提高模型的计算效率。在例2.5中，块体尺寸的变化范围从1.751×10 Print max 键入如下命令，删除极小块体： delete range area 3e-2 所有面积小于3×10-2的块体都从模型中删除。通常，将小于最大块体的1％左右小块被删除后对计算结果的影响并不显著。 最后，注意到NEW命令用在第二个例子，以便允许开始一个新的模型。当切割块体（尤其当采用JSET）时，一个重要问题是综合考虑块体数与计算速度的协调。计算速度与模型的块体数（或变形体单元数）成函数关系。根据经验，模型具有大约1200刚体（或具有8自由度的500变形体）进行2000～4000迭代步就能获得静态问题的解。对于90MHz的微机，对于500个变形块体模型允许4000步大约需要10分钟。根据你的计算机的计算速度，可以估算出一个模型所需的计算时间。 通过切割UDEC块体形成工程结构形状，这须在进行工程开挖前实施。通常采用三个命令来产生形状： crack tunnel arc －3 到3.679，可以由以下命令查找: －13－ 前面已经给予介绍CRACK 命令。TUNNEL命令产生圆形形状。该圆由用户指定的裂缝段数构成。ARC命令由用户指定的角度，产生弧形断裂模型。可以结合这些命令产生各种形状的UDEC块体。例2.6给出的命令产生断层切割一个圆形隧道的模型。 EXAMPLE 2.6 断层切割一个圆形隧道 New Round 0.1 Block -10,-10 -10,10 10,10 10,-10 Tunn 0,0 2 16 Crack -5,10 5,-10 Plo hold block num 所生成的模型如图2.12所示。圆形隧道的圆形坐标（0,0）、半径为2和划分成16个裂缝段。由于隧道全部处于块体内部，所以仅用TUNNEL命令不能产生独立的块体。必须采用CRACK切割模型块体的边从而产生新的块体。如果用户运行仅用TUNNEL命令所产生的模型，则隧道裂缝在运行前被删除。通过引入CRACK命令，连接隧道裂缝延伸到模型外边界从而形成连续的裂缝，因此，形成有隧道和断层构成的块体（如图2.12所示）。 应当注意，裂缝并不贯穿隧道的周边。如果TUNNEL命令先给出，随后的CRACK或JSET命令并不贯穿隧道。首先应用TUNNEL是较为方便的，因为隧道开挖仅涉及删除一个块体，即 delete range block 1920 或 delete range -1,1 -1,1 将模拟圆形隧道的开挖。 图2.12 断层切割圆形隧道 －14－ 例2.7给出了一条断层切割一个马蹄形隧道： New Round 0.1 Block -10,-10 -10,15 10,15 10,-10 arc 0,5 2,5 180 8 Crack -2,0 -2,5 Crack -2,0 2,0 Crack 2,0 2,5 Crack -5,15 5,-10 隧道的形状如图2.13所示。隧道顶弧的圆心在（0，5），起始点在（2，5）和180o的圆心角，逆时针画圆，划分成8段。前三个CRACK命令产生隧道的边墙和底板，后一个CRACK产生一条切割隧道的断层。开挖隧道也可以通过删除块体编号或块体形心位置来实现。另一选择就是指定包含该隧道块体 window -2,2 0,7 delete range window 图2.13 断层与马蹄形隧道相交 2.6.2 指定材料模型 2.6.2.1 块体模型 一旦完成块体切割，必须对所有的块体和不连续面指定材料特性。缺省为所有的块体皆为刚体。在多数分析中，块体应为变形体。仅仅在应力水平较低或岩块材料具有高强度和低变形的情况才能够应用刚性块体的假设。 －15－ 块体的变形特征通过以下命令定义： gen edge v 或 gen quad v GEN 命令激活三角形网格有限单元自动生成器。命令GEN edge v 将作用于任意形状的块体。其v值定义三角形单元的最大边长，即v值越小，块体中的单元越小。应当注意的是：具有高的边长比值的块体并不能产生单元，其极限的比重近似为1:10。 通过 Plot zone 检查模型单元。 采用命令GEN quad v，指定模型为塑性材料模型的单元。该类型的单元提供了对于塑性问题的精确解。然而，GEN quad 命令可能对某些形状的块体不起作用。在此情况下，应当采用GEN edge 。 在UDEC中为变形块体（单元）开发了7种材料模型。对大部分用户，最常用的三种模型如下： change cons=0 ; null model change cons=1; elastic model change cons=3; Mohr-coulomb model CHANGE 命令改变块体为指定的变形块体。Cons=0意味着模型块体材料被移出或开挖。这允许用户改变块体在以后的某些阶段返回为弹性或弹塑性材料。如果块体被删除，则以后计算阶段不可再恢复。 Cons=1 改变块体为各向同性弹性特性；而cons=3则改变块体为摩尔－库仑模型，考虑塑性特性。缺省值为所有变形体则自动改变为cons=1。 块体改变为cons=1和cons=3必须提供PROPERTY mat 命令给块体赋予材料参数值。注意性质参数不要赋给特定的块体，而是赋给材料号。材料参数可以赋值给多达50种材料号。然后，材料号再赋给具有CHANGE mat 命令的块体。 对于弹性模型，需要的性质为： （1）密度 （2）体积模量 （3）剪切模量 注意：体积模量K、剪切模量G与杨氏模量E、泊松比?之间的关系如下： K＝或 E＝EE，G＝ （31－2?）（21??）9KG3K－2G，?＝ （23K＋G）3K＋G －16－ 对于摩尔－库仑塑性模型，需要的性质为： （1）密度 （2）体积模量 （3）剪切模量 （4）内摩擦角 （5）粘聚力 （6）剪胀角； （7）抗拉强度 如果上述参数没有赋值，系统自动赋零值。 对于UDEC程序，对于上述两个模型，密度、体积模量和剪切模量必须赋予正值。 2.6.2.2 节理模型 除了给块体赋予材料模型外，还应对模型中的所有不连续面（即接触面）赋予材料模型。对于不连续面，有四种本构模型。在UDEC中开发了四种节理本构模型。但对于大部分模型分析，最适宜的模型有库仑滑动模型（完全弹塑性），可通过如下模量赋予不连续结构面： change jcons=2 所有不连续结构面的缺省模型是Jcons=2。 节理材料模型也是通过PROPERTY mat 命令赋予材料性质参数。如同块体，参数并不是直接赋给不连续面，而是材料号。材料号是通过CHANGE mat 命令赋给节理。对于库仑滑动模型，所需要的参数是： （1）法向刚度 （2）切向刚度 （3）内摩擦角 （4）粘聚力 （5）剪胀角； （6）抗拉强度 如果所有的参数没有给赋值，他们的缺省值为零。在UDEC中，必须给节理的法向和切向刚度赋值，并为正值。例2.8演示了材料模型的应用。 例2.8 指定材料模型与性质参数 New Round 0.1 Block -10,-10 -10,10 10,10 10,-10 tunnel 0,0 2,16 jset -70,0 40,0 0,0 40,0 -1,-1 －17－ jset -50,0 40,0 0,0 3,0 0, 2 gen edge 2.0 change jmat=2 range angle -51,-49 change jmat=5 range angle -71,-69 pro mat=1 d=2500 b=1.5e9 s=0.6e9 pro jmat=1 jkn=2e9 jks=2e9 jcoh=1e10 jten=1e10 pro jmat=2 jkn=2e9 jks=1e9 jfr=45 pro jmat=5 jkn=2e9 jks=1e9 jfr=10 change cons=0 range -1,1 -1,1 在上面的例子中，一条与x轴方向成－70o单一断层切割圆形隧道。还有倾角为－50o和间距为3m的节理组。块体是变形的和具有最大边长是2m的三角形单元。图2.14显示了隧道、节理和单元。块体是弹性的，其性质参数是通过PROP mat=1赋值。50o的节理组是通过性质参数号2赋值。70o的断层参数是由材料参数号5赋值。jmat=2，angle=-51, -49 仅将性质参数赋予-51o和-49o之间。jmat=5，angle=-71, -69仅将其值赋予范围角位于-71o和-69o之间。所有其他不连续面被赋予材料号1，这意味着将具有高粘聚力和抗拉强度的粘结效应。这些不连续面是“虚拟的”节理和对应于后来将隧道开挖。虚拟节理不发生滑动或张开。可用以下命令检查材料号。 Plot block mat 隧道块体被改变为cons=0是模拟开挖。 Plot zone可得到下图。 图2.14 圆形隧道与70o断层和50o节理组构成的模型 2.6.3 施加边界条件和初始条件 －18－ 在完成所有块体切割（节理切割）和变形单元划分之后，应施加边界条件和初始条件。施加力学边界条件通常采用BOUNDARY命令。该命令用来指定力、应力和速度（位移）边界条件。边界力和应力能够施加到刚体和变形体的边界上，但速度（位移）边界仅适用于变形块体（见施加刚性块体的边界命令FIX，FREE和LOAD）。表2.4提供了边界条件命令的总结和效果。 BOUNDARY xload 和yload命令施加x-和y-方向的分力到边界角点。BOUNDARY stress 施加应力张量到边界上。BOUNDARY xvel 和yvel 在选择的边界结点施加x-和y-方向的速度。 注意的是：应用BOUNDARY命令所产生的条件或约束将不发生改变（除非用户再次改变）。 表2.4 边界条件命令总结 命 令 BOUNDARY Stress Xload Yload Xvel Yvel FIX FREE LOAD Xload Yload 效 果 施加总应力到刚体或变形体块体的边界上 施加刚体或变形体边界的x方向的荷载 施加刚体或变形体边界的y方向的荷载 施加变形体边界的x方向的速度（位移） 施加变形体边界的y方向的速度（位移） 固定刚体边界的速度（位移） 释放刚体的速度（位移） 施加x方向的荷载到刚体的边界 施加y方向的荷载到刚体的边界 初始应力条件能够被指定到所有的变形单元和所有刚体或变形体之间节理的法向应力和剪切应力。INSITU命令用来初始化应力。采用该命令，可以赋值初始应力。边界的初始条件应用于例题2.8，其命令应用见例2.9。 例2.9 施加的边界条件和初始条件 boundary stress -10e6,0,0 range -11,-9 -10,10 boundary stress -10e6,0,0 range 9, 11 -10,10 boundary stress 0,0，-5e6 range - 10,10 9,11 boundary yvel =0.0 range -10,10 -11, -9 insitu stress -10e6,0,-5e6 szz -4.8e6 10MPa压应力被作用到模型的左、右边界的x方向上。5MPa的压应力（负号为压）被施加到上部边界的y方向。底部边界的y方向的运动被固定。固定边界的位移对 －19－ 于考虑重力的情况下是尤其重要的。注意：应力边界影响所有的自由度。因此，在施加速度边界条件前的同一边界，应当施加应力边界条件，否则，所指定的速度约束将不起作用。而且，四个BOUND命令的每一个都应给出x和y的坐标范围。应当提醒的是为了确保通过BOUND命令所产生的响应完全位于指定的范围。键入如下命令： Print bound 和 Plot bound xcond Plot bound ycond 来检查边界条件。 命令INSITU在x方向始化所有的应力为-10MPa和在y方向初始化应力为-5MPa。在平面之外的z方向也给予初值?zz=-4.8MPa。对于弹性块体分析，z方向的应力未初始化并不影响平面应变问题的解。然而，对于塑性分析，z方向的应力可能影响到破坏状态，因此，应当慎重选择?zz的应力初始化。 2.6.4 迭代为初始平衡 UDEC模型在进行开挖模拟前必须进行初始状态的平衡计算。施加合适的边界条件和初始条件，使模型与初始的平衡状态相吻合。然而，对于复杂的几何形状和多介质材料的情况，在给定的边界条件和初始条件下，进行计算获得平衡是十分必要的。其计算可采用STEP（或CYCLE或SOLVE）命令。借助于STEP命令，为达到模型的平衡，用户指定循环步进行计算。当每一刚体形心的结点力 或变形体的结点力接近零，模型就处于平衡状态。当激活STEP命令后，最大的结点力矢量（称之为不平衡力）由UDEC进行监测，并在屏幕上显示。在此用户能够估计模型何时达到平衡状态。 对于任何模型的数值分析，不平衡力不可能完全达到零。但当最大的结点不平衡力与初始所施加的总的力比较相对较小时，就可认为模型达到平衡状态。例如，如果最大不平衡力从最初的1MN降低到100N，此时（最大不平衡力与初始的不平衡力之比为0.01％），则认为模型达到平衡。 采用UDEC进行数值分析判断模型平衡是一个重要的问题。用户必须确定模型在何时达到平衡状态（即问题的解）。在UDEC中设置有一些特征，用于支持这种决策。如记录最大不平衡力历史是其中之一： hist unbal 除此之外，还有速度历史（即某一结点的速度或位移），命令如下： hist xvel 5, 5 hist ydisp 0, 11 第一个是记录位移坐标（x=5, y=5）附近结点x方向的速度，而第二个是记录接近坐标（x=0, y=11）位置处y方向的位移。在进行数百次或数千次迭代后，这些历史记录将绘图和显示其平衡条件。 －20－ 图2.15所示的模型为节理岩体中开挖一矩形区域。在此岩体中存在一组节理和一条断层。模型受10MPa的静水应力场和重力作用。重力是用以下命令： set grav 0.0 , -9.81 第一个是x方向的加速度，第二个值为y方向的加速度为9.81m/sec2（向下作用）。 当考虑重力所引起的应力变化不大可以忽略重力的作用。在例2.10中，尽管考虑重力，有助于识别围绕洞室周围的松散块体受重力的作用，但是考虑重力引起的应力变化小于0.5MPa，相对于10MPa的原岩应力可以忽略不予考虑。 ro 0.1 bl -10,-10 -10,10 10,10 10,-10 cr -2,-2 -2,2 cr -2,2 2,2 cr 2,2 2,-2 cr 2,-2 -2,-2 jset 70,0 40,0 0,0 40,0 -2,0 jset -50,0 40,0 0,0 3,0 1,2 gen edge 2.0 chan jmat=2 range angle -51,-49 chan jmat=5 range angle 69,71 prop mat=1 d=2500 b=1.5e9 s=.6e9 prop jmat=1 jkn=2e9 jks=2e9 jcoh=1e10 jten=1e10 prop jmat=2 jkn=2e9 jks=1e9 jfr=45 prop jmat=5 jkn=2e9 jks=1e9 jfr=5 bound stress 0,0,-10e6 range -10,10 9,11 bound xvel=0.0 range -11,-9 -10,10 bound xvel=0.0 range 9,11 -10,10 bound yvel=0.0 range -10,10 -11,-9 insitu stress -10e6,0,-10e6 szz -10.0e6 set grav 0.0 -9.81 hist unbal hist ydis 0,2 ; solve for 10 step 700 save fall1.sav 如果位于开挖体顶板的块体被分离，则由于重力作用将落到洞室。这将在后面的2.6.5节中详细说明。应当注意，当重力应力与原岩应力具有相同的量级，则用INSITU命令施加应力梯度（应力随高度的变化）以加速初始平衡状态的收敛。 －21－ 图2.15 节理岩体中的矩形开挖体 初始不平衡力近似为2.0MPa。在进行700步后，降至约为10N。通过绘制两个历史可以发现，最大的不平衡力已接近零，位移接近2.7×10-3m常量（图2.16和2.17）。 图2.16 最大不平衡力历史 图2.17 在（0，2）位置的y位移历史 －22－ 模型的不平衡力在接触力和块体角点与结点的力略有差异，这与圆形角点有关。圆形在块体角点产生很小的“洞”。为不平衡力的平衡增加某些迭代步是必要的。 设置不平衡力的值： solve force = f 在此，f 为用户定义的不平衡力值。 检查不连续面的破坏条件是十分重要的。对于这个问题，初始x发现和y方向的应力分量是相同的，因此，模型中的节理不能滑动，这可以用以下命令加以验证： Plot bou slip 所显示的图形是由模型外边界和满足库仑滑动准则的节理组成。对于选择的模型条件，有可能出现在原岩应力状态下节理发生滑动。例如，改变水平应力分量?zz到-5MPa并回到例2.10。则断层将沿着整个长度上发生滑移。如何识别初始应力在平衡计算过程中节理发生的滑移，用户应当重新评估所选择的原岩应力参数和不连续面强度。沿着节理长度滑移的模型表明该模型并非可靠。 在模拟开挖前确保模型处于平衡状态是十分重要的。通过记录几种历史以考察最大的不平衡力的衰减。如果所进行的计算步超过模型达到平衡所需的计算步，并不会影响计算结果。然而，如果不充分的计算步将影响模型的计算结果。 UDEC计算可在任何时间通过按 2.6.5 进行改变和分析 UDEC允许在求解过程中的任意部位改变模型条件。这些变化可能具有以下形式： （1）开挖材料； （2）增加或删除边界荷载或应力； （3）固定或释放边界结点的速度（位移）； （4）改变材料模型或块体和变形体的性质参数。 可以用DELETE 命令或CHANGE cons=0命令模拟材料开挖。用BOUNDARY xload，yload 或stress 命令施加荷载和应力。通过采用BOUNDARY xvel 或yvel命令固定边界角点。通过BOUNDARY xfree 和yfree 命令移去边界约束。用CHANGE命令改变变形块体和不连续面的材料模型。而用PROPERTY命令可改变材料性质参数。 很显然，几种命令可以重复应用，进行各种模型的改变。例如，从初始平衡状态，应用这些命令继续例2.10获得例2.11。 例2.11 开挖隧道和监测其响应 rest fall1.sav delete -2,2 -2,2 －23－ reset disp reset hist hist unbal hist ydis 0,2 step 2000 plot blo stress disp save fall2.sav 由于采用DELETE命令，模拟开挖矩形洞室，导致模型应力的变化。结点位移与历史记录被重新设定，仅由开挖所引起位移变化被监测。建议在块体被删除后重新设定历史位置。 在矩形洞室开挖后产生很高的不平衡力，因此需要进行计算使之获得重新平衡。然而，在此情况下，没有观测到不平衡力接近很小的值，而处于其值为0.017MN的常值。进而y位移历史记录也显示在计算2000步后，在位置（0，2）的位移仍向下运动。开挖顶板上的块体已经从围岩脱离和掉落到洞室内。这由图2.18清楚地看出。不平衡力不可能接近于零，因为顶板块体应自由下落。 如果预计模型变化将导致破坏（即力的平衡条件不能获得）就不要用SOLVE命令求解。 图2.18 洞室顶板块体发生冒落 2.6.6 保存或恢复计算状态 当进行分步计算时，另外两个命令SAVE和RESTORE是有用的。在一个阶段的结 －24－ 尾（即初始平衡），采用如下命令，可以保存模型状态。 Save file.sav 式中，file.sav 是一个用户定义的文件名。扩展名.sav定义这个文件是一个保存文件。这个文件可以采用如下命令进行恢复： rest file.sav 2.6.7 简单分析的总结 在表2.5中给出了本节所介绍的主要命令。更常见的情况求解问题需要从进行一个简单的问题计算开始。 表2.5 简单问题分析的基本命令 功 能 产生块体模型 ROUND BLOCK 切割块体 CRACK JSET TUNNEL ARC 块体和节理的材料模型和参数 GEN CHANGE PROPERTY 边界条件和初始条件 BOUNDARY INSITU 初始平衡（具有重力） DAMP local SET gravity STEP SOLVE 模型变化 DELETE CHANGE PROPERTY BOUNDARY CABLE 监测模型响应 HISTORY PLOT 保存或恢复当前状态 SAVE RESTORE 命 令 －25－ 2.8 系统单位 表2.6 系统单位－力学参数 SI Length Density Force Stress Gravity m kg/m3 N Pa m/sec2 m 103kg/m3 kN kPa m/sec2 M 106kg/m3 MN MPa m/sec2 cm 106g/m3 Mdynes Bar cm/sec2 3 用UDEC求解问题 本章为应用UDEC求解岩石力学工程问题提供指导。在3.1节给出地质力学分析阶段的建议。在3.2节通过例3.10 明确了模型准则和求解过程必须考虑的问题，涉及的问题如下： （1）产生模型（3.2节）； （2）刚体或变形块体的选择（3.3节）； （3）边界和初始条件（3.4和3.8节） （4）加载和模拟顺序（3.6节）； （5）块体和节理模型和材料参数的选择（3.7和3.8节）； （6）改进模型效率的方法（3.9节）； （7）计算结果的解释（3.10节）。 最后，地质力学领域的模拟原理参见3.11节。在该领域进行模型分析的新手可能希望首先咨询该节。地质力学模拟方法与其他工程领域，如结构工程存在很大的不同。进行地质力学分析始终记住这一点是十分重要的。 3.1 一般性研究 模拟地质工程的过程涉及一些特殊的考虑，其设计方法也与其他人工材料结构不同。在岩土体上建造结构或在其中开挖分析与设计，必然面对相对少的现场数据以及材料的变形和强度性质参数存在较大变化的情况。获得岩土工程现场完整的现场资料是不可能的。例如，原岩应力、材料性质和不连续面特性等信息仅是部分的。 由于为设计预测所输入的必要信息是有限的，所以，地质力学数值模型主要用于理解影响系统特征的力学机理。一旦掌握了系统的特性，然后，为工程设计过程探索一些简单的计算。 面对地质工程研究总是缺少满意的数据以及缺乏对材料性质的充分理解，而在其他 －26－ 领域，在具有充足的数据的情况下，应用UDEC可直接用于工程设计。当获得合理的数据应用此程序时总能获得合理的结果。应当认识到，如图3.1所示的应用的过度阶段。 典型情况 复杂的地质条件； 简单的地质条件 不可获取的数据； 投入资金 无试验经费 进行现场研究 数 据 研 究 无 资料充分 机 理 通过参数研究 预 测 研 究 现场特性分类 （直接用以设计） 图3.1 模型研究图谱 UDEC程序可用于模型特性的预测（如图3.1的右边），或仅作为“数值试验”来测试一些设想（图3.1的左边）。正是现场资料（和资金）而不是程序决定了应用情况。如果具有高质量的和足够的数据，UDEC就能够给出好的预测。 由于大部分的UDEC分析是处于较少数据的情况下进行的，所以，本节将探讨类似于试验研究的数值模型研究技术。数值模型决不应当被认为是一个一端接受信息而另一端输出结果的“黑箱”。为了获得可以合理的解释，必须十分注意准备数值“样本”和多样本“试验”。表3.1 列出了进行成功的数值模拟试验的建议步骤。下面分别讨论： 表3.1 地质力学问题的数值分析步骤 第 1步 第2步 第3步 第4步 第5步 第6步 第7步 定义模型分析的对象 产生模型系统的概念图形 建造和允许简单的理想模型 搜集模型所需的计算数据 准备一系列用于分析的详细模型 进行模型计算 提供结果和解释 3.1.1 第1步：定义分析模型的对象 一个分析模型所了解内容与深入程度常常取决于分析的目的。例如，如果是为解释系统的特性所提出的两种相互冲突机理的决策，此时可建造一个较粗糙的模型，用于两种机理的研究。如果试图涉及存在于实际模型中的复杂条件，然而，如果可能对模型的响应产生微不足道的影响或与模型计算的目的毫不相关的计算特征可以被忽略。 3.1.2 第2步：产生物理系统的概念图形 －27－ 重要的是构思出实际问题的图形，便于初步估计在所施加的条件下，预测系统的基本特性。当准备这个图形时，应当回答几个问题。例如，该系统是否稳定？主要力学响应是线性还是非线性？是否存在可能影响特性的不连续面？是否存在地下水的影响？实际的系统物理结构是否还存在其他几何问题？ 这些考虑将表征了诸如模型的几何形状、块体材料模型、边界条件以及初始平衡条件等数值模型的总体特征。这将决定是否采用三维模型或二维模型？ 3.1.3 第3步：建造和运行简单的理想模型 当为数值分析理想化一个物理系统时，较有效的方法是在构筑详细的模型之前，首先建造和运行一个简单的测试模型。为产生数据和对问题的理解，应在尽可能早的阶段产生这样的一个简单模型。其结果可用于进一步系统的概念图形。在简单模型运行后可能需要重复第2步。 简单的模型能揭示一些问题，以便在进行深入分析之前加以修正。例如，所选择的材料模型是否能够代表所期望的特性？边界条件是否影响模型的响应？基于简单模型的计算结果能够有助于指导对分析起重要影响作用的数据研究方案。 3.1.4 第4步：综合特定问题的数据 对于一个模型分析所需的数据类型包括： ? 详细的几何参数（即地下洞室形状、地表形态、坝形状，岩石或土体结构）； ? 地质结构的位置（即断层、层理，节理组等）； ? 材料特性（即弹性或塑性性质，峰后特性）； ? 初始条件（即原岩应力状态，孔隙压力，饱和度）； ? 外部加载（即爆破荷载、洞壁压力）。 由于分析所涉及的条件（尤其应力状态、变形和强度性质）存在很大程度的不确定性，为研究必须选择参数的合理变化范围。基于简单模型的计算结果（第3步）常常能够有助于确定变化范围。 3.1.5 第5步：准备一系列详细的运行模型 通常数值分析用于研究不同的破坏机理、研究一定范围变化的计算参数的系列分析。当为计算准备一系列计算模型时，应考虑如下一些方面的问题： （1） 每一个计算需要花费多少时间？如果模型运行的时间过长，为达到有用的结 论所需要获得足够的信息可能是困难的。为缩短计算时间，可以考虑在多个计算机上，运行参数变化的计算。 （2） 应考虑保存所需要的模型在计算过程中的中间状态，以便每一参数的变化不 必重复计算。例如，如果分析几个加载或卸载阶段，用户应当能够返回运行已经进行的任一阶段，以便改变一个参数后从那一阶段继续计算。 －28－ （3） 在模型中是否设置足够的监测位置（历史记录），为进行清楚地解释模型计 算结果和不同计算参数的比较分析提供足够的信息？在模型中设置几个参数变化的监测点，对计算过程监测是有帮助的。尤其模型中的最大不平衡力应当被监测，以便检查在分析的每一阶段的平衡或破坏。 3.1.6 第6步：进行模型计算 在进行一系列模型分析之前，最好首先选择一个或两个模型进行详细地分析。这些运行应当随时被中断，确保达到预期的效果。一旦能够确信模型的计算是正确的，几个模型数据文件被联系在一起，进行一系列模型的连续分析。在连续运行的任何时间，应有可能中断计算，查看结果，然后继续或修改模型。 3.1.7 第7步：提供结果和解释 求解的最后一步是为进行清楚地解释分析提供计算结果。最好是通过直接在屏幕上显示或输出的图形结果。图形结果应当提供便于进行计算与现场观测结果的比较方式。图形应当能够清楚地分析所感兴趣的区域，例如应力集中位置、模型中稳定与不稳定区域。模型中任何变量的数值也能够容易获得,为详细解释模型的响应。 为有效地进行地质工程问题分析，我们建议了上节介绍应遵循的七个步骤。下面章节将论述在UDEC的应用中为模型研究中涉及到的这些步骤中的每一方面所涉及的特定问题。 3.2 产生模型 UDEC程序在产生几何模型的方式与传统的数值分析程序有所不同。首先产生计算范围的单一块体。然后，这个块体被切割成小的块体。模型中块体的边界是地质结构面或工程结构（如开挖体边界）。这种切割处理被称之为节理生产的几何体。然而，“节理”代表物理模型中的实际地质结构和人造结构边界或将被移去或在以后连续的计算步中改变材料。对于后者，节理是虚拟的，其存在不应影响模型的计算结果。虚拟节理的表征将在3.2.3节中讨论。 3.2.1 确定UDEC模型合适的计算范围 UDEC几何模型必须具有足够大的范围，在感兴趣的区域内，包含主要的地质结构特征，由此代表真实的实际的物理问题。考虑的方面如下： （1） 处于何处的地质结构（即断层、节理和层面）应详细描述？ （2） 模型边界的位置对模型的影响程度如何？ （3） 如果应用变形块体，在关心的区域，何种密度的单元可满足问题的精度？ 上述三个方面决定了实际分析UDEC模型的规模。 －29－ 如果仅有很少的地质结构（即两个或三个切割断层或遍布空间的节理组），这些可能通过CRACK或JSET命令单独输入。记住UDEC是一个二维程序。除了特殊情况，忽略三维效应。如果地质结构不能用垂直于分析平面的二维特征加以表征，则可能需要采用三维分析（例如ITASCA程序3DEC）。 单个独立特征的断裂可用一种或两种方式进行输入。一种是用CRACK命令，他给出断裂的两个端点；另一种是用JSET命令，给出断裂的倾角和断裂通过的位置。例如，或者 Crack (0,0) (10,10) 或 Jset (45,0) (20,0) (0,0) (100,0) (5,5) 都可用来定义一个倾角为45o，且通过坐标点（x=5, y=5）的节理。为了在UDEC计算中得到认可，节理必须是连续的（即完全劈裂块体成两个部分）。然而，节理可以由分段连续、具有不同角点的分段构成。CRACK和JSET产生不连续节理段（SPLIT命令与CRACK具有相同的形式，但他不能产生不连续节理）。例如，燕尾状节理可以用例3.1中的命令生成。 例3.1 燕尾状节理 Round 0.1 Block 0,0 0,10 10,10 10,0 cr 0,5 2.5,6 cr 2.5,6 5,5 cr 5,5 7.5,6 cr 7.5,6 10,5 在上述例子中，几个CRACK命令的顺序没有限制。如果被后来产生的裂缝交切，内部裂缝被保存在一个临时性的文件中和后来被应用。在节理生成过程中，任何内部或部分贯穿裂缝，在模型运行过程中被删除。当块体被赋予可变形的（GEN），内部或部分裂缝也自动被删除。也可采用JDELETE命令实施人工删除。 对于块体作指定的圆角长度可能局部影响节理的产生。最小块体棱长定义为圆角长度的两倍。因此，节理段可能背离这个准则。例如，例3.2 显示了指定一个裂缝的端点位置处于角点圆角长度的两倍，该裂缝通过角点的位置坐标发生错位。UDEC并对任何这样错位不会给出任何警告，所以用户应当通过采用PLOT block 命令随时进行检查。 例3.2 圆角长度对产生裂缝的影响 New Round 0.2 ; Rou 0.1 －30－ Block 0,0 0,10 10,10 10,0 cr 0.3,0 9.7,10 如果圆角长度小到0.1，则裂缝将处在指定的位置。 SET edge 命令运行用户人为定义的最小块体棱长。基于这个命令，用户为精确求解而避免块体具有小的棱长，或反之具有很大的边长比，能够设定一个较小的圆角长度。例如，如果设置命令SET edge 0.4 和ROUND 0.1，则不会产生小于0.4的块体棱长和块体的圆角长度将为0.1。这些命令必须在BLOCK命令之前给出。 模型中用于描述地质特征的节理数（即块体数）存在一个限制。这涉及模型的范围和块体的单元数（如果采用变形单元）。实际的限制依赖于如同表2.2所示的可利用计算机内存。在进行节理生成时必须考虑这个限制。根据经验，最后的方式总是从较少节理开始，然后，如果有必要再逐渐增加节理来达到预期的效果。应当避免试图建立一个复杂节理模型的诱惑，关于此问题将在3.11节中进一步讨论。 节理自动生成器可以在UDEC中得到，可根据实际测量的参数（即节理倾角、间距、长度和岩桥长度等）产生节理组。JSET命令也起到节理生成器的作用。如在第3.2.2节的例子所描述的那样。在UDEC模型还可获得一种特殊的生成Voroni形状多边形的节理生成器。在3.2.2节给出该程序的应用实例。高级用户还可编写自己的节理生成器。这可通过将节理生成器所定义节理段端点坐标（x，y）对组成的表列来实现。FISH函数可以完成这个自动处理过程。这个列表能够被UDEC直接读取。另外，图形数字化仪也可以用来产生节理端点的坐标对。 请记住，在UDEC内的节理是作为直线段显示。许多线段可能需要用不规则的节理形态来拟合。模拟者必须决定用UDEC模拟实际节理模式的水平。几何的不规则对节理的影响也可以通过节理材料模型性质（即变化节理面的性质参数）加以考虑。 如上所述，模型边界必须具有足够的远，以致模型对边界不产生影响。第3.4节论述了边界的影响结果。一般地，对于单一地下开挖工程，边界离开挖边界的距离应当大于开挖跨度的5倍左右。 然而，合适的距离取决于分析的目的。如果分析目的主要用于考虑破坏，然而，模型边界可以靠近一些。如果关注的是位移（变形），则距离边界的距离需要增加。 借助于经验估计边界对模型的影响是重要的。从一个粗糙的模型和矩形边界开始，分别研究固定边界和自由边界条件研究改变边界的距离的影响。考察感兴趣的模型区域内的应力或位移，随边界距离变化所发生的变化，来评价边界对结果的影响。参见第3.4节所给出的边界影响研究实例。 一旦完成块体切割和确定模型的边界位置，下一步就是考虑应采用块体单元的大小与网格密度。较密的网格单元应当处在高应力区或高梯度变形区（即在开挖区附近）。为了高精度，单元形状尺寸之比（即三角边与高之比）也应尽可能接近于1。对于5:1 －31－ 的情况可能是不精确的。同时建议相邻块体单元的大小不应有较大的突变。合理的精度是相邻两单元面积之比不应当超过4:1。 3.2.2 产生节理 UDEC提供了两个节理生成器：一个统计节理生成器，由传统的岩石力学参数所定义的参数产生节理；另一个是Voronoi分块式节理生成器，用于产生随机尺寸的多边形块体。统计节理生成器是采用JSET命令和岩体不连续面的几何参数统计特征值。Voronoi生成器是采用VORONOI命令和划分块体为随机大小和形状的子块。该生成器将在下节加以描述。 3.2.2.1 统计节理组生成器 JSET节理生成器是根据所选定的统计参数生成节理模式。基于这样的节理模式，特殊的几何参数对力学特性产生的影响可以进行定量描述。同时，在获得现场观测的节理模式情况下，采用人工模式与观测模式相匹配方法决定生成的参数。 一个节理组可以通过8个生成参数表述：4个几何参数的均值和4个随机参数的均方差。由JSET命令给出的参数如下： JSET am、ad、tm、td、gm、gd、sm、sd?x0、y0??ad0? 在此，a－节理与x轴的夹角； t－节理段迹线长度； g－两节理段间的长度（即岩桥长度）； s－垂直于节理迹线的间距。 对于上述的每一对值，前一个带有下标m的是均值；第二个带有下标d的是均方差。图3.2给出了参数的说明。 图3.2 节理组参数 －32－ 最后的三个参数是选择参数： x0、y0 － 节理起始坐标（整体坐标轴）。节理将从点（x0、y0）开始产生节 理，将充满由选择参数range 所定义的整个范围。 ad0 － 所有节理与所给定的节理迹线方向的偏差。 通过定义的一个限制区域（range），生成的节理能够限制在所选择的模型区域。见理论与背景第1.1.3节中的range 关键词的描述。在大部分情况下，通过jregion n 关键词定义各种区域，在此n 是指定JREGION 命令的参考序号id。 JREGION 定义了一个凸多边形区域来限制节理组的生成范围。该命令参数如下： JREGION id n x1 y1 x2 y2 x3 y3 x4 y4 每一个节理区域是通过id序号识别。区域的坐标按顺时针方向定义了节理产生的边界。如果给出选择的关键词delete ，在此之前由JSET、VORONOI或CRACK命令生成的节理将全部被删除。这可避免当指定多个节理区域情况下，被相邻区域产生的节理所切割。 其他相对于JSET命令（即 mat n）的关键词range也能用来限制块体被切割。如果对JSET命令没有指定范围，将产生遍及整个区域的节理。 JSET 生成器也可能通过设定很大的tm和sm值用来产生贯穿UDEC块体的一条节理。请记住，JSET能够产生不连续节理。不连续节理完全处在块体里边不能看见。当一个变形块体开始生成单元或当一个刚体模型开始运行时，所有的节理段将被删除。 3.2.2.2 VORONOI多边形生成器 VORONOI生成器产生随机大小的多边形块体。UDEC模型中的一个或多个块体能够在细分成任意大小的Voronoi子块。该节理生成器对模拟裂缝扩展是有用的。当Voronoi块体间的节理强度被超过时将发生断裂。 VORONOI命令具有下列形式： VORONOI edge l < iterations n> 对于Voronoi多边形指定平均棱长。多边形具有随机尺寸，但具有平均棱长l。Voronoi块体的大小尺寸能够通过增加迭代次数变得均匀缺省值n=5。也能指定圆角长度。圆角长度，v 必须至少小于块体棱长l的20倍。 当应用VORONOI命令时，多边形区域应当略大于被细分的块体区域。这将抑制边界影响。Range 关键词在产生Voronoi块体区域的用法，与JSET命令具有相同的方式。 Voronoi算法根据随机分布点从多边形区域开始。然后允许内部点移动。迭代过程运动到这些点。迭代步越高，点间距越均匀。接下来，所有点产生三角形。最后，通过做具有公共边所有三角形的垂直平分线，生成Voronoi多边形。多边形在所指定区域被其边界所截取。 3.2.2.3 例子 －33－ 为了说明节理的生成过程，下面给出了产生节理的几个例子。 例3.3 四组规则节理组 new ro 0.01 blo 0 0 0 20 20 20 20 0 jset 0 0 50 0 0 0 3 0 jset 90 0 50 0 0 0 3.5 0 jset 30 0 50 0 0 0 4 0 jset 50 0 50 0 0 0 6 0 ret 3.2.3 产生内部边界形状 当完成一个计算区域的UDEC建模后，还必须定义与实际物理模型的边界条件相一致的块体边界条件。这些可能是代表地下开挖，或孔洞或诸如土坝的人工结构或天然边坡表面特征的外部边界。如果实际物理问题具有复杂的边界，重要的问题是评估简化模型对对计算结果的影响程度，即简化后的简单模型是否能够重新给出重要的模型机理。 在模拟所有的物理边界（包括后来的模型中增加或开挖的区域）在求解之前必须加以定义。在以后的分析中增加结构形状必须在事先定义，并通过CHANGE cons=0或ZONE model null 移动，只有当在出现的阶段（如分析筑坝阶段）再出现。 用下列命令实施边界形状的产生： CRACK JSET TUNNEL ARC 每一个命令切割块体成一条或多条裂缝，集合组成所要求的形状。CRACK和JSET命令产生节理的直线段。TUNNEL和ARC命令分别形成圆形线段和弧段。 在大部分情况下，在人造特征之前应当产生模型的自然特征（如节理、层面、地表形状、坡面等）。例如，在例3.9中首先产生节理结构，然后形成圆形开挖体。 这个问题的模型显示在图3.11中。有一组倾角为110o和间距0.5的节理组和三组三条独立节理。每一条的倾角为5o。圆形隧道的圆心为（0，0）、半径为2和边界由32段组成。由JSET产生的节理贯穿隧道周长。 注意，在UDEC中，正倾角是从x轴的正轴方向顺时针到节理进行测量的角度。 因为隧道是圆形的，所以隧道的开挖可以直接用DELETE annulus命令删除。 Delete range annulus 0，0 0，2 另外，如果TUNNEL命令在JSET命令前给出，则节理并不贯穿隧道周边（如图 －34－ 3.12）。因此给出坐标范围也可以删除隧道，即 delete range -0.1,0.1 -0.1,0.1 该命令仅适用于TUNNEL命令。 如果TUNNEL命令被首先给出，应当检查模型以确保隧道内部无节理化区域的存在并不影响初始应力状态。比较在模型中，无隧道和有隧道边界两种情况的应力分布，确信应力状态不受边界的影响。 图3.11 在生成节理组后产生隧道 图3.12 在节理组生成前产生隧道 不考虑天然节理与内部边界节理（虚拟的）的产生顺序。但以后必须考虑沿隧道周边的虚拟节理并不影响模型响应。周围块体间的虚拟节理必须被粘合在一起，以阻止滑 －35－ 动或张开。这可通过设定较高的粘聚力和抗拉强度值来实现。剪切刚度和法向刚度也必须赋予足够高的值以阻止两块体产生不同的节理位移。 正是给节理赋予很高的刚度以阻止虚拟节理发生运动。然而，由于UDEC的计算步取决于刚度，因此，当给予很好的刚度值，系统的响应和收敛速度将十分缓慢。所以，建议最低的刚度应与小裂隙变形分析所采用的刚度相一致。一个好的经验是节理刚度，kn和ks，取10倍于相邻单元的最高刚度等效值。单元在法线方向的等效刚度（产生单位位移所需要的应力）是： ?(K?43G?max?? （3.1） ??zmin?在式中，K、G － 分别是体积模量和剪切模量； ?zmin － 相邻单元在垂直方向的最小宽度（图3.13）。 max??符号表示对于与节理相邻的所有单元的最大值（即相邻节理可能存在几种材料）。 图3.13 单元尺寸用于刚度计算 通过采用CHANGE 命令能够改变虚拟节理的材料性质（见2.6.2节）。例如，在上述的例子中，隧道节理性质由下列命令进行改变： change jmat=2 range ann 0,0 1.99, 2.01 在圆心为（0,0）、半径r1=1.99和r2=2.01之间的环状节理将具有节理性质号为2。该变化将由以下命令进行查看： plot mat joint 虚拟节理的性质参数通过PROPERTY jmat=2 命令予以赋值。 现论述建立模型的最后一个问题。对于连续介质分析时，通常利用开挖体的对称条 －36－ 件，以利于减小模型尺寸。但对称条件并不易应用于不连续介质，因为，不连续特征的存在，限制模型的对称性，除非在特殊情况。 3.3 变形块体和刚体的选择 不连续介质分析的一个重要问题，是采用刚性块体还是采用表征完整岩石特征的变形块体的决策。本节将探讨选择刚体或变形体的问题。如果需要进行变形分析，由几个不同的模型可以用于模拟块体的变形特征。这些在3.7节讨论。 如在理论与背景的第1节所述，早期的离散元程序假设块体是刚性的。然而，对于考虑块体的变形特性的重要性已经得到共识，尤其研究地下工程的稳定性分析和深埋结构的动力响应更为重要。在离散元分析中考虑块体的变形特征的一个最明显的理由是表征受侧向约束的“泊松比效应”。 岩石力学问题通常对岩体泊松比的选择是十分敏感的。这是因为节理和完整岩石对压力是敏感的：其破坏准则是侧向压力的函数（即摩尔－库仑准则）。获得节理岩体的真实泊松特性，对于具有实际意义的数值模拟是十分关键的。 岩体有效泊松比是由两个部分组成：（1）由节理位移产生；（2）由完整岩石的弹性性质引起。除非在浅埋条件或低侧向应力水平，完整岩石的压缩性对整个岩体的压缩性起到重要作用。因此，完整岩石泊松比对节理岩体泊松比产生重要影响。 严格地讲，单一的泊松比?，仅仅定义弹性材料。最多仅适用于极少节理产生各向弹性特性。所以，为讨论各向异性的岩体材料，定义“泊松效应”是十分必要的。 泊松效应将定义为：当在垂直方向施加荷载的条件，不允许在水平方向产生应变（位移）时水平应力与垂直应力之比。 假设平面应变条件。作为一个例子，对于各向同性弹性材料的泊松效应为： ?xx? (3.2) ??yy1??推导垂直节理模型的泊松效应如图3.14所示。如果该节理用刚性块体模拟，施加的垂直应力将根本不产生水平应力。这显然不实际的，这是因为忽略了由完整岩石的泊松比产生的水平应力。 节理和岩块以串联方式相互作用。换句话说，作用在节理上和岩块上的应力是相等的。节理岩体的总应变是节理应变和岩块应变之和。岩体弹性性质可以通过叠加节理和岩块的变形得到： ??xx?rock?Cjointing???=C?yy?????xx???? (3.3) ?yy? －37－ 图3.14 岩体内具有水平和垂直节理的泊松效应模型 如果以各向同性弹性材料模拟完整岩块，其变形矩阵为： C节理变形矩阵为： rock1???E-???1?? ??? ? (3.4) 1-???Cjointing?1? 0?sk?n? (3.5) ??1??0 ?skn???式中，s － 节理间距，kn － 节理的法线刚度。 如果方程式（3.3）中的?xx=0，则 total?xxC12?total (3.6) ?yyC11式中，Ctotal?Crock?Cjointing。 因此，节理岩体的总泊松效应为： ?xx?(1??) (3.7) ??yyE/(skn)?(1??)(1??)方程(3.7)给出了以E/(skn)作为变量的函数（如图3.15所示）。而且此图是为了检验此公式所进行UDEC模拟的结果。其比值E/(skn)是与节理刚度相关的完整岩块刚度的度量。对于E/(skn)的低值，岩体泊松效应主要受控于完整岩块的弹性性质。对于 E/(skn)的高值，泊松效应受控于节理。 现在考虑变化的不同倾角岩体的泊松效应。泊松效应是节理产状和弹性参数的函 －38－ 数。考虑如图3.16所示的特例。岩体包含两组等间距，倾角为与x轴的夹角?的节理组。节理弹性性质是由法线刚度和切向刚度组成。完整岩块假设为完全刚性的。 图3.15 垂直节理岩体的泊松效应（完整岩石的泊松比 ??0.3） 图3.16 与x方向夹角为?、间距为s的节理岩体泊松效应模型 对于这种节理模型的泊松效应为 ?xxcos2? [(kn/ks)?1] (3.8) ??yysin2? ?cos2? (kn/ks)上述公式中几个不同的?值的图形如图3.17所示。在此图也给出了UDEC数值模拟结果。UDEC模拟的结果与方程（3.8）十分接近。 －39－ 图3.17 具有不同节理角的节理岩体（块体是刚体）的泊松效应 方程（3.8）表明了在数值分析中选用合理的节理剪切刚度是重要的。剪切刚度与法线刚度之比值戏剧性的影响岩体泊松效应。如果剪切刚度等于法线刚度，则泊松效应为零。对于从2.0到10.0比较合理的kn/ks比值，泊松效应是相当的高，达到0.9。 其次，完整岩块的弹性性质对其影响将在?＝45o的情况下加以研究。在下面的分析中，节理为垂直分布模式，处理完整岩石为各向同性弹性材料。岩体总弹性性质将由叠加节理和岩石的变形获得： 具有倾角为45o，等间距的两组节理的为： Cjointing?12sknks?ks?kn ks?kn??k?k k?k? nsn??s因此，岩体总泊松效应为： ?xx[?(1??)]/E?(kn?ks)/(2sknks) (3.9) ??yy[?(1??)(1??)]/E?(kn?ks)/(2sknks)方程（3.9）中对于泊松比?=0.2和几个不同的E/(skn)比值的情况下，给出的图形如图3.18所示。对于较低的E/(skn)比值，岩体的泊松效应受控于完整岩石的弹性性质。对于较高的E/(skn)比值，泊松效应受控于节理。 －40－ 图3.18 倾角?＝45、等间距两组节理岩体的泊松效应（泊松比?=0.2） o3.4 边界条件 在数值模型中，边界条件是由表征模型边界的场变量（即应力或位移）组成。边界分为两类：实际的和人工的。实际边界是模型中存在的物理对象（即隧道界面或地面）。人工界面并不真正存在，但为了封闭所选择的单元数（即块体）必须引入。能够施加的边界类型是类似的。这些条件将首先讨论。然后（在第3.4.4节），对涉及到边界的位置和人工边界的选择以及对计算结果的影响提出某些建议。 力学边界具有两类：位移和应力。自由边界是应力边界的特例。两类力学边界在3.4.1和3.4.2中讨论。第三种类型，“边界元边界在3.4.4中论述。用于动力分析的粘滞边界和自由场边界在理论和背景中的第4节中讨论。 3.4.1 应力边界 UDEC模型的缺省边界是无约束的自由边界。力或应力可以通过BOUNDARY命令施加到任意整个边界或部分边界上。用stress关键词可以指定平面应力张量（?xx、?xy、?yy）的每一个单独应力分量。例如，命令 Boundary stress 0,-1e6,-2e6 range 0,10 -1,1 将施加?xx?0、?xy??10、和?yy??2?10到位于坐标窗口0< x <10，-1 Print boundary 每一外边界的角点将以表的形式列出指定的边界值。边界在模型计算过程中可能运动，所以，用户必须检查坐标窗口是否足够大，使窗口包含所有的边界角点。 另外，边界条件也可在两边界的角点间施加。作用的边界从第一个边界角点地址顺 66 －41－ 时针到第二个边界角点地址施加。角点的地址可以用下面的命令进行查看： Print boundary state 角点地址以标题COR和列的方式显示。最好给出坐标范围而不是角点地址，因为，应用的地址与模型相关，即如果改变命令排序，地址可能发生改变。 在UDEC中按照一般符合惯例，压应力为负号。而且，UDEC实际上施加应力分量作为力或者是产生作用到给定的边界平面上的应力张量的摩擦力。摩擦力被分成两个分量：永久的和瞬时的。动力分析采用随时间变化的永久摩擦力是常荷载和瞬时荷载。采用xload和yload 关键词也可将单个x方向和y方向的力的分量施加边界上。用LOAD命令也可指定荷载于刚性块体上。刚性块体是施加到块体的形心上。 3.4.1.1 施加应力梯度 BOUNDARY命令可以增加关键词xgrad和ygrad，此关键词允许应力或力在指定的边界上按线性变化。在上面的每一个关键词后面的关键词，给出了在x或y方向变化的应力分量： Xgrad sxxx axyx syyx Ygrad sxxy sxyy syyy 应力从坐标原点（x=0，y=0）沿边界按距离线性变化： 0?xx??xx?(sxxx)x?(sxxy)y 0?xy??xy?(sxyx)x?(sxyy)y 0?yy??yy?(syyx)x?(syyy)y 0式中，?xx，?xy，?yy是原点的应力分量。 00上式运算过程由下面的例子加以解释： boundary stress 0,0,-10e6 ygrad 0,0, 1e5 range -0.1,0.1 -100,0 0 原点的应力是?xx=0，?xy=0，?yy=-10×106。 00在y方向的应力变化?yy的计算式为 ?yy = -10×106＋（105）y ?yy在y=-100的值为-20×106。在边界上，y的变量是从原点开始线性变化。 3.4.1.2 改变边界应力 正如上述讨论，瞬时荷载可用history关键词，为动力分析施加边界条件。对于静力分析，在UDEC模拟过程中，改变边界的应力值也是必要的。例如，作用到基础的荷载可能发生改变。为了考虑应力或荷载突然变化的影响，需要给出新的BOUNDARY命令，对原来边界角点施加改变后的应力或荷载。 －42－ 在这种情况下，新的边界应力或荷载值将被叠加到原有的值上。如果应力被移动，当前的值应当以相反的符号给出。如果一瞬时荷载被改变（即用history 关键词指定荷载），用相同的history 所给出新的荷载被叠加到存在的荷载上。然而，用不同的history所给出的一个新的瞬时荷载，将取代老的荷载。 3.4.1.3 打印和绘图 边界应力和荷载可以通过命令PRINT bound 和 PLOT bound xcond 或 PLOT bound ycond 进行验证。PRINT bound 命令列出每一角点所指定的边界角点地址和对应的当前的值。一旦BOUNDARY 命令被发出，整个模型的外边界生成边界角点列表。选择的关键词可以用PRINT bound 命令来检查边界的不同条件。例如 print bound force 列出永久力（fx,fy）和在当前荷载步增加的荷载增量（fxi，fyi）的表列。如果施加瞬时荷载（用BOUND … hist命令），总荷载为永久荷载与当前循环步的瞬时荷载之和。 需要施加边界条件的边界角点地址能够通过以下命令输出显示。 Print bound state PLOT bound xcond 或 ycond 命令通过符号表示施加在x-或y-方向边界条件的类型。 3.4.1.4 提示和建议 本节将探讨施加应力边界条件所存在的一些琢磨不定的困难。应用UDEC程序，可将应力施加于一个没有位移约束物体的边界上（不同于很多的有限元程序，在此需要某些约束）。物体将完全按照实际物体所作用的方式产生反作用，即如果边界应力没有平衡，则整个物体将产生运动。例3.10说明这个效应。 图3.19给出了所生成的图形。施加的?11引起水平应力作用到物体上。由于该物体是倾斜的，所以力导致运动引起物体的转动。相似但影响有所不同的是从一个物体中开挖应力边界附近的材料。物体最初在重力作用下平衡，但移去材料减小重量，整个物体开始向上运动（例3.11）。 在运行中这个数据文件（例3.11）所遇到的困难能够通过固定底部边界予以消除。第3.4.4节包含涉及这样的人工边界的位置的信息。 最后，应力边界影响所有的自由度。所以速度边界条件必须在应力边界条件之后给出。如果应力边界在速度边界之后施加，给出速度的响应将会消失。例3.12说明此问题。 －43－ 图3.19 力作用于倾斜物体产生的旋转位移 3.4.2 位移边界 UDEC模型不能直接控制位移。事实上，他们对计算过程并不起作用。为了施加已知的位移到边界上，有必要固定边界或对于给定的迭代步，给出边界的位移速率（采用COMMARY命令）。如果要求的位移是D，对应于时间为增量T（计算步），施加的速度为V（即D=VT），在此，T??tN，?t是时间步，N是迭代步数（循环数）。实际上，为减小模型系统的波动，V应当赋予小值而N应大值。 BOUNDARY命令用来沿着边界（而不是与x或y坐标轴一致）固定在x或y方向（BOUND xvel 或yvel）或法线或切线方向（BOUND nvel 或svel）的变形块体的结点速度。刚性块体的速度可以用FIX命令。如果FIX命令放在INITIAL命令之前，速度能够用用户所选定值加以固定。也可用FISH函数加以改变。 随时间变化的速度历史能够通过BOUND… hist 命令对刚体或变形体予以施加。这个history 关键词必须出现在指定速度历史边界BOUND xvel 或BOUND yvel的同一行。历史也可通过FISH函数施加。正如在3.4.1.4节所讨论的，速度边界应当在应力边界之后给出。变形块体的固定速度边界条件能够通过BOUND xfree 或BOUND yfree命令得以释放。对应刚体，则用FREE命令。 通过在指定的速度关键词之后增加关键词gvel，速度也允许在指定的边界范围按线性变化。Gvel关键词后有六个参数，描述了速度分量在x或y方向的变化： gvel vx0 vy0 vxx vxy vyx vyy 从坐标原点（0，0），速度随距离线性变化： vx=vx0+(vxx)x+(vxy)y vy=vy0+(vyx)x+(vyy)y (3.11) －44－ 在式中，vx0和vy0是在原点的速度分量。 3.4.3 真实边界－选择合理类型 准确地识别计算模型中某些特殊表面边界条件的类型有时是很困难的。例如，在一个三轴试验模型中，通过台板施加的荷载是应力边界？或台板作为刚体而视为位移边界？当然，包括台板的整个试验机能够一并模拟，但可能是非常耗时的。记住如果存在很大差异的刚度，UDEC花很长时间才能收敛。一般情况，如果施加荷载的物体与样本相比非常“刚”（比如说，是弱者的20倍）。则可以作为应力控制边界加以模拟。很清楚，作用在物体表面上的流体压力属于后一类。位于节理岩体上的基础，常常可移去常值的速度而作为刚性边界加以处理，用于研究岩体的破坏荷载。这个研究还有另一优点－它非常容易地控制试验和获得较理想的荷载或位移图形。人所共知的刚性试验机远比一般试验机稳定。 3.4.4 人工边界 人工边界分为两类：对称轴和截取边界线。 3.4.4.1 对称轴 当几何形状和荷载是关于一轴或多轴对称时，就可能考虑其特点进行模型的简化。例如，如果一结构的对称轴是垂直的，则对称轴上的水平位移均为零。所以，我们能够取其垂线作为边界，边界上的所有结点，用BOUND xvel=0设置x方向的位移为零。如果在对称轴上的速度为零，也可采用此命令进行设置为零。在此情况下，y方向的位移并不受到影响。类似的情况是水平轴为结构的对称轴。命令BOUND nvel=0能够用来设置与坐标轴成一定角度的对称轴。 正如在3.2.3节所讨论，不连续面的存在使得对称性的利用变得困难。当采用对称性，应当十分小心地考虑节理产状的影响。 3.4.4.2 截取边界 当模拟无限大物体（即地下隧道）或非常大物体时，由于内存或计算时间的限制，其计算模型不可能涉及整个物体。人工边界应置于所关心区域的足够远，以致不影响模型效应。知道这些人工边界应置多远以及对模型的计算应力和位移产生多大的误差是有用的。基于包含两个隧道的弹性介质模型进行一系列数值试验。该模型仅包含用于产生两个隧道的虚拟节理。图3.20和图3.21分别显示最小和最小的模型。 在所有的计算中，块体中的单元尺寸是相同的，以致消除单元离散的影响。考虑两个物理量：在较大跨度隧道顶板中心点的垂直位移，uy，和两个隧道间柱中点的垂直应力，?yy。原岩应力比是2:1（垂直与水平之比）。对应每一矩形边界的几何体，进行两中边界条件的计算：长应力边界和零位移边界。另外，还进行了UDEC边界元边界的计算。例3.13包含上述三种边界的数据文件。 －45－ 图3.20 包含两个隧道的小模型 图3.21 包含两个隧道的大模型 图3.22用无量纲的形式显示了计算结果。观测的位移和应力相对于他们的无影响值进行正则化处理。“边界尺寸”的值是平均宽度和高度，“物体尺寸”是贯穿两隧道的平均距离。从图3.22中给予我们几点提示： －46－