岩土工程数值计算方法

岩

土

工

程

数

值

计

算

方

法

报

告

学院：土木与环境工程学院

姓名：xxxxxx

学号：xxxxxxxx

三维有限差分稳定性分析

一、FLAC3D基本原理

FLAC(Fast Lagrangian Analysis of Continua)是由美国Itasca 咨询公司研究开发的显式有限差分程序，可用于工程力学计算，模拟岩石、土等材料的力学行为。由于其采用了显式拉格朗日算法及混合离散划分单元技术，使得该程序能较好地模拟地质材料在达到强度极限或屈服极限时发生的破坏和塑性流动，分析渐进破坏和失稳，特别适用于模拟大变形。材料通过单元和区域表示，根据计算对象的形状构成相应的网格。每个单元在外载和边界约束条件下，按照给出的本构关系产生力学响应。FLAC 软件主要是为岩土工程稳定性分析开发的岩石力学计算程序，它包括了反映地质材料力学效应的特殊计算功能，能够计算地质类材料的高度非线性(包括应变硬化/软化)、不可逆剪切破坏和压密、粘弹(蠕变)、空隙介质的应力—渗流耦合及动力学行为等。

FLAC 提供了多种材料本构模型：各向同性弹性模型、横观各向同性弹性模型、摩尔-库仑塑性模型、应变硬化/软化塑性模型、德鲁克-普拉格塑性模型、遍布节理模型、双屈服塑性模型、霍克-布朗模型、空单元模型等。另外，程序设有界面单元，可以模拟断层、节理和摩擦边界的滑动、张开和闭和行为。支护结构，如砌衬、锚杆、支架等与围岩的相互作用也可以在FLAC 中进行模拟。同时，用户可根据自己的需要在FLAC 中创建自己的本构模型，进行各种特殊修正和补充。

FLAC 采用显式算法来获得模型全部运动方程的时间步长解，从而可以追踪材料的渐进破坏和跨落，这对研究开采的时间效应和空间效应是非常重要的。此外，程序允许输入多种材料类型，亦可在计算过程中改变某个局部的材料参数，增强了程序使用的灵活性，用来提供采动区域的跨落过程和开采中的充填过程。FLAC 具有强大的后处理功能，用户可以直接在屏幕上绘制图形，或以文件形式创建和输出打印多种形式的图形。使用者还可以根据需要，将若干个变量合并在同一幅图形中进行研究分析。

FLAC3D 的计算过程如图1.1 所示，在每个时步内，首先，根据高斯定律，由节点速度求出新的应变速度；根据本构方程计算出各时步内单元新的应力；然后，由运动方程计算新的结点速度和位移。

在各个时步内进行循环计算，当计算达到平衡状态时，不平衡力趋近于0，如果有塑性流动产生，则不平衡力趋近某一数值。

图1.1 FLAC 3D

循环求解过程 二、节理化模型

节理化模型是各向异性塑性模型，它包括了包含在摩尔—库仑体内特殊方向上的弱面。根据应力状态、弱面走向以及模型体和弱面的材料特性的不同，屈服可能发生在模型体内，或者发生在弱面上，或者在两个部位同时发生。

这种模型在FLAC 中的实现方法是首先判别总体破坏，同时应用到和FLAC 中摩尔—库仑模型中相同的相关塑性修正，然后对更新的应力在弱面上产生的破坏进行分析，同时对这些应力分别进行进一步的校正。弱面内的破坏准则存在于包含了拉应力路径的摩尔—库仑屈服条件的局部形式中，与局部剪切流动法则不相关联而与局部拉应力流动法则相关联。

图1.2显示了整体坐标系（x ，y ）和局部坐标系(x ′，y ′)下存在于摩尔—库仑体内的软弱面。

图1.2 整体坐标系中沿θ角方向的软弱面

为了简化这部分的符号，定义σij ，它对应于由各阶段总体破坏的塑性修正的应用引起的应力分量，这种总应力分解成局部应力后可表示为：

θσθθσθσσ2221221111

sin cos sin 2cos ++='

θσθθσθσσ2221221122

cos cos sin 2sin +-=' 3333σσ='

)sin (cos cos sin )(2212221112

θθσθθσσσ-+--=' 式中：θ——节理夹角（从x 坐标轴逆时针方向算起）。

依照约定，由τ表示弱面上切向引力分量的大小，相应的应变变量为γ，可以看到： 12

12e τσγ'='=

有了这两个符号，弹性增量的法则的局部表达式可以表示为：

11

1112223322

1332112233

13321122()()()2e e e e e e e e e e

e e e e e e e e e G σαασαασαατγ??????????????''''=++''''=++''''=++=

式中：α1=K ＋4G ／3 ，α2=K －2G ／3 ，上标e 代表“弹性部分”。

软弱面的破坏准则可以在图1.3所示的22σ'，τ平面内表示出来。

图1.3 FLAC 中的软弱面破坏准则

将摩尔—库仑破坏准则定义为s f ＝0，从点A 到点B 的局部破坏包络线定义为：

j j s C f +'--=φστtan 22

将拉应力破坏准则定义为t f ＝0，从点B 到点C 的局部破坏包络线定义为：

22

σσ'-=t

j t f 式中：j ?，j c ，t j σ——分别为弱面的摩擦角、粘聚力和抗拉强度。 对于摩擦角不为零的弱面，抗拉强度的最大值定义如下：

j

j t

j C φσtan max ,=

剪切和张拉势函数s g 和t g 对应于不相关联的流动法则，剪胀角j ψ对应于相

关联的流动法则，它们分别表示为：

j s g ψστtan 22

'--= 22

σ'-=t g 破坏准则边界附近的流动法则用FLAC 中摩尔—库仑模型中所述的方法定

义。这里，由函数h(22σ'，τ)=0表示(22σ'，τ)平面中s f ＝0和t f ＝0曲线的对

角线（图1.4），此函数形式为：

)(22

t

j p j p j h σσαττ-'--=

图1.4 节理化模型中用于定义弱面流动法则的区域

其中p j τ和p j α是两个常量，定义如下：

t j j j p j C σφτtan -=

j j p j φφαtan tan 12-+=