基于FLAC_3D_二次开发的不连续变形模拟方法探讨

图 3 二次开发流程图 Fig. 3 Flow chart of secondary development
屈服强度为 60.62 MPa，屈服强度所对应的应变值约 为 5.5 mm。该模拟结果与文献[23]中的 60.96 MPa 和 5.0 mm（见图 6）是非常接近的。
3 单轴压缩试验验证
为验证程序的正确性，取典型的单轴压缩试件进 行受力分析。 3.1 计算模型 为与文献[23] 中的室内试验结果进行对比分析， 选取的单轴压缩试件，其尺寸为 20 mm×40 mm，弹 性模量为 14 GPa，泊松比为 0.26。利用 FLAC3D 以速 度方式进行动态加载模拟（上下两端自由，仅施加速 度荷载） 。 具体模型尺寸如图 4 所示， 计算参数见表 1。
Discussion of discontinuous deformation simulation method based on secondary development of FLAC3D
ZHENG Wen-hua1, LIU De-jun2
(1. Institute of Foundation Engineering China Academy of Building Research, State Key Laboratory of Building Safety and Built Environment, Beijing 100013, China; 2. Department of Geotechnical Engineering, Tongji University, Shanghai 200092, China)
程模拟分析中具有一定的局限性。 FLAC3D 软件采用的是连续介质快速拉格朗日分 析方法，其能进行非线性大变形问题的模拟分析，且 计算结果容易收敛，在岩土工程的分析计算中有着较 其他数值模拟方法更大的优点， 被广泛地应用于土建、 [11-13] 。已有部分学者利用 FLAC3D 交通、采矿等领域 软件中的置空功能（model null）来进行裂纹开裂过程 的模拟分析[14,15]。但是该方法存在的问题是将单元置 空后，其两侧单元不再具有任何相互作用，即裂纹完 全断开。这与实际工程中裂纹张开过程中相互作用力 逐渐消失以及产生滑移变形的实质是不相符的。 因此， 3D 本文在 FLAC 软件的基础上，探讨通过 FLAC3D 内 部提供的 interface 接触面单元来模拟劈裂裂纹的张开 扩展过程。
接触面由一系列三角形单元组合形成，其中每一 个三角形单元由 3 个节点组成。一般来说，接触面单 元与实体单元的表面黏结在一起；接触面的节点自动 在每一个接触面单元的顶点创建而成。当其他网格的 表面和接触面单元发生接触时，接触面的节点就会发 现这种接触，并且由法向、切向刚度及滑动的相关性 质所定义。每一个接触面的节点都有一个相关联的特 征面积。图 1 为接触单元模型。 1.2 接触面单元的本构模型 接触面单元本构模型的组成部分如图 2 所示。图 kn 为法向刚度， Ts 为抗拉强度， Ss 2 中 ks 为剪切刚度， 为抗剪强度， D 为膨胀角， S 为滑动组件。
图 1 接触单元模型 Fig. 1 Interface model
1.1
接触面单元的组成成分
2 不连续变形模拟程序的二次开发
2.1 FLAC3D 不连续变形模拟思路 莫尔库仑模型可用于模拟岩土体的受力变形过 程，利用其塑性屈服准则可给出塑性区岩土体的分布 范围，但是该模型不能考虑进入塑性区岩土体的强度 折减和损伤效应。本文提出的 FLAC3D 不连续变形模 拟方法中，充分考虑了塑性区岩土体的强度折减及损 伤效应。 在模拟计算过程中，裂纹开裂准则选取最大拉应 力准则[18-20]，即当单元体所受最大拉应力超过其抗拉 强度时，该单元发生破裂，即裂纹扩展。接触面单元 的力学参数在条件允许时应根据实测的损伤区岩土体 的力学参数来设定，可利用下式计算初始损伤[21]： Vpm 2 D0 1 ( ) ， (4) Vpr 式中，Vpm，Vpr分别为初始损伤岩体和无损岩体（完 整岩石）的弹性纵波波速，若无实测值可参照以往学 者的工程经验和研究结果[22]来设定。 在 FLAC3D 中，接触面的法向刚度和切向刚度可 由下列公式计算得出： EE r kn ， (5) s(Er E ) GG r ks ， (6) s (G r G ) 式中，E、G 为岩体的杨氏模量和剪切模量，Er、Gr 为完整岩石的杨氏模量和剪切模量。
0
引
言
随着计算机技术的飞速发展，数值分析方法因其 计算成本低、效率高等自身优势，而被广泛地应用于 岩土工程问题的模拟分析。岩土工程的数值分析方法 可以分为连续变形分析方法和不连续变形分析方法， 主要有：有限元法、有限差分法、离散元法、DDA 法 [1-2] 、流形元法以及无网格法等。其中，常用于进行不 连续变形分析的软件有 UDEC、3DEC[3]、PFC[4]等离 散元方法以及 DDA 不连续变形分析方法。 90 年代末， [5-8] 并在最近几年被广泛用于 出现了扩展有限元方法 ， Abaqus 软件裂纹开裂的模拟分析中[9-10]。 Abaqus 软件通过扩展有限元的方法可以较好地 模拟裂纹扩展开裂 的不连续变形 过程，但是由于 Abaqus 软件采用的是有限元法， 在进行一些复杂问题 的模拟分析时，特别是几何大变形问题时，计算工作 量相当大，且计算结果较难收敛，因此，在实际的工
─────── 收稿日期：2013–03–03
198
岩 土 工 程 学 报
2013 年
1
FLAC3D 不连续变形分析的基础
FLAC3D 软件中提供了 interface 接触面单元，该 单元可用于模拟节理、断层、桩土接触面等物体间的 接触、滑移、张开等问题[16-17]。该接触面单元具有摩 擦角、黏聚力、膨胀角、法向刚度、切向刚度以及拉、 剪黏结强度等参数，能很好地模拟裂纹面间的滑移、 张开直至完全断开的过程。
图 2 接触单元本构模型[17] Fig. 2 Interface constitutive model
1.3
接触面单元的求解过程 FLAC3D 中提供接触面单元类型：当其受外力值
增刊 1
郑文华，等. 基于 FLAC3D 二次开发的不连续变形模拟方法探讨
199
接触面单元的建立可以有效地模拟岩土体的损伤 破裂过程。在接触面建立的开始阶段，接触单元等效 于损伤区岩土体，在进一步的受力过程中会产生相应 的滑移变形；当剪应力大于剪切强度或法向拉应力大 于抗拉强度时， 接触面间的连接会断开， 即裂纹张开。 2.2 二次开发程序的实现过程 在动态循环计算 利用 FLAC3D 中 FISH 语言功能， 过程中加入接触面单元的建模及求解过程，实现图 3 所示的 FLAC3D 不连续变形模拟分析的二次开发。
图 4 试件模型及加载示意图 Fig. 4 Schematic of sample and loading
3.2 计算结果 在单轴压缩过程中，随着竖向荷载的施加，试件 受力过程中会经历弹塑性阶段、损伤阶段、断裂阶段 一系列过程。在断裂阶段，岩土单元强度降低，但仍 保留部分残余强度。在这一过程中，接触面单元的引 入可以很好地实现岩土体损伤以及破裂过程的模拟。 根据模拟过程中记录的试件应力–应变关系曲线（见 图 5）可以看出：在 FLAC3D 模拟计算中，该试件的
Abstract: On the basis of Mohr model, a simulation of discontinuous deformation is studied using the FLAC3D inner interface model and the secondary development platform of the FISH language-compiles program included in FLAC3D software. The problem of dynamic modeling of interface elements in FLAC3D is solved, and uniaxial compression tests and comparative studies are carried out. The analytical results reveal that the discontinuous deformation simulation method based on FLAC3D is effective and feasible. The simulated results are similar to the experimental ones, and the proposed method can effectively catch the developing process of cracking. This study provides an new view for simulation of discontinuous deformation. Key words: discontinuous deformation; crack; secondary development; interface; uniaxial compression
小于黏结强度时，接触面单元为弹性；当剪应力大于 剪切强度或法向拉应力 （有效应力） 大于抗拉强度时， 接触面间的连接断开。同时，当接触面在连接未断开 时，不产生分离变形，仅产生剪切滑移变形。 （1）弹性阶段的控制方程： (1) Fn( t t ) kn un A n A ， ( t t ) (t ) ( t (1/ 2) t ) ， (2) Fsi Fsi ks usi A si A 式中，Fn( t t ) 为 t t 时刻的法向力，Fs(i t t ) 为 t t 时 刻的剪切力。 （2）库仑滑移阶段的控制方程 (3) Fsmax cA tan ( Fn pA) 。 从以上描述中，可以看出 FLAC3D 中的接触面单 元可以用于模拟实际工程中裂纹的张开扩展过程。但 是由于 FLAC3D 中接触面的建立较为复杂，且其动态 接触面建模过程需要利用 FISH 语言进行 FLAC3D 软件 的二次开发。因此，实现 FLAC3D 不连续变形分析模 拟的关键问题是 interface 接触面单元动态建模过程的 二次开发。