含孔洞节理岩体损伤破坏过程的颗粒流数值模拟

含孔洞节理岩体损伤破坏过程的颗粒流数值模拟
张敏思王述红 侯佳男 郭牡丹 杨勇
（东北大学资源与土木工程学院，沈阳 110004）
摘要：基于颗粒流理论，引入颗粒接触连接本构模型，建立了岩体颗粒流数值模型。

通过含孔洞花岗岩试样损伤破坏试验和Particle Flow Code（PFC)离散元程序，对带有圆形空洞的花岗岩的单轴压缩试验进行了数值模拟，模拟试验观察到的空洞周边扰动区的损伤破坏过程。

与室内试验结果对比分析，PFC2D基本再现了脆性岩石的局部受压破坏现象，并且在破坏过程中可以得到应力-应变关系曲线，及其裂缝扩展过程图，进而通过声发射技术来监测微裂缝的演化过程。

通过改变计算模型中颗粒单元的性质，给出了在不同颗粒单元参数时试样的宏观性质，并且损伤破坏试验验证了本模拟的正确性，其研究结果对岩体本构关系的完善有一定的应用价值。

关键词：颗粒流；圆形空洞；数值模拟；声发射
Particle Flow Code Numerical Simulation Of Failure Around a
Circular Opening In Jointed Rock Mass
ZHANG Minsi, WANG Shuhong, HOU Jianan,
GUO Mudan, YANG Yong
（School of Resource and Civil Engineering, Northeastern University, Shenyang 110004, Liaoning,
China）
Abstract: Based on particle flow code（PFC）theory , rock mass numerical models are constituted , in which different contact bond constitutive relation are introduced. The Uniaxial compression tests of failure around a circular opening in Granite sample are simulated by PFC2D respectively, and also by the lab test of failure around a circular opening in Granite sample. The failure process of peripheral disturbance around a circular opening is also simulated. Compared with the indoor test results, the local pressure brittle rock damage phenomenon is represented by PFC2D, and also the stress-strain curve and the crack propagation graph can be obtained through the failure process, the tiny cracks of the evolution process is monitored through the acoustic emission. Through the changes in the calculation model of particle nature of the unit, the macro properties of different parameter in the sample unit are given, more over, the correctness of the simulation is validated by the failure and simulation experiment and there has been certain practical value for the perfect constitutive relation of of rockmass.
Keywords: particle flow code, a circular opening, numerical simulation, the acoustic emission
基金项目：国家高技术研究发展计划（2007AA06Z108）；国家重点基础研究发展计划（2007CB209405）；辽宁省自然科学基金（2009）
作者简介：张敏思，女，山东淄博人，东北大学硕士研究生，主要从事地下结构方面的学习和研究工作。

Email:mscccathy@；王述红，男，江苏泰州人，东北大学副教授，主要从事岩土力学与工程方面的教学和研究工作。

1 引言
节理裂隙岩体在物理、力学特征上较一般工程材料具有更显著的各向异性、非线性及非连续性，其变形与破坏性质十分复杂。

不论从完善岩体力学理论的角度，还是从满足大量岩体工程实践需要的角度，都十分有必要深入开展节理裂隙岩体破坏机理的研究。

国内外对节理裂隙岩体破坏机理研究的常用方法主要有断裂力学方法、损伤力学方法、试验方法、数值方法（如离散元方法）等。

离散单元法由Cundall[1]于1971年首先提出，最早由王咏嘉和剑万禧于1986年介绍给国内的学者。

至今，模拟地下空洞周边的损伤区已经有了一些成功的尝试，Potyondy and Autio [2] and Cundall [3]使用这个程序已经预测了Underground Research Laboratory (URL) in Manitoba, Canada的裂纹区。

2002年，Potyondy[2]模拟了在压缩荷载作用下英闪长岩圆形试验空洞周边的损伤区。

目前，随着离散单元法思想的发展，把其运用到节理岩体的破坏过程仍然有待进一步完善，根据工程实际及理论知识，用数值模拟方式进一步的研究岩体的损伤及其演化，进而预测岩体的宏观力学响应(包括强度、变形与破坏)、分析岩体断裂失效先兆、评价工程岩体的稳定性具有非常重要的意义。

本文提出基于颗粒流理论，引入颗粒接触连接本构关系，建立了含孔洞花岗岩试样模型，进行了数值模拟，模拟试验观察到的空洞周边扰动区的损伤破坏过程，并与室内试验结果[4]对比分析。

下面先介绍一下实验过程和结果，然后用所建立的模型进行数值分析，最后用颗粒流程序PFC2D来模拟这个试验，将模型校准以后再进行带有孔洞试件的数值模拟。

结果显示，PFC2D再现了脆性岩石的局部受压破坏现象，数值试验的破坏区与试验室观察到破坏区非常吻合。

2 试验过程
2.1 试验模型
使用取自河北曲阳的带有空洞的长方体花岗岩试件上进行试验。

试件的加载方向的长度为100mm，截面为100mm*100mm（图1）。

侧向尺寸为100mm*200mm，孔洞尺寸20mm。

使用的花岗岩力学参数为：弹性模量取73.8MPa（国际岩石力学实验建议），体积模量取43.9GPa,泊松比取0.22，密度取2600kg/m3。

在整个试验过程中移动垂直驱动器来保持侧向位移变化率为0.02mm/min。

a）岩样正面
a）Positive figure of rock
sample
b）岩样右侧面
b）Right side figure of rock
sample
图1 试验模型
岩石类材料的破坏过程实际上是材料内部微裂纹的损伤演化过程，随着荷载的增加，微裂纹慢慢的延伸和扩展，最终形成宏观裂缝导致材料破坏[5]。

岩石渐进破裂失稳过程中，岩石应变、声发射及其损伤程度之间的关系为：N D
ε∝∝。

岩石属于脆性材料，其损伤属于脆塑性损伤[6-8]。

用声发射次数表示岩石
损伤程度，基于岩石声发射次数的本构方程可以表示为：
00
exp[()]m N E N σε=⋅−[9] 此本构方程反映了某种岩石在达到一定应变时岩石的损伤情况，但m 值需要实验来确定，要直观确定岩石的损伤是比较困难的问题。

2.2 试验方法及结果
应用加拿大ESG 公司生产的声发射系统，建立在到时残差模型，在最小二乘法的基础上，开发单纯形定位方法[10]，进而反演声发射事件模式。

由于花岗岩颗粒坚硬，尺寸较大，变形也较大，在加载的初期没有明显的声发射现象，然后随着荷载的增加，定位点沿着裂纹扩张的方向（岩样的两端）不断的增加。

应变（荷载）较小的阶段，内部发生的破裂事件是较少的，测得的声发射要比实际产生的少，随着应变（荷载）的增大，声发射数的增加，岩石的损伤就有了明显的变化。

把花岗岩的参数带入到损伤本构模型中得到的应力-应变曲线与试验结果相比较，得到图2。

应力 /M P a 应变
图2 理论值与试验值的对比
由于声发射信号本身复杂多变，声发射数本身产生的误差，导致理论曲线与试验曲线存在少许偏差[11]。

45MPa 55MPa 65MPa 75.5MPa
图3 声发射正面定位图
45MPa 55MPa 65MPa 75.5MPa
图4 声发射侧面定位图
图5 花岗岩岩样的破坏形式
图3及图4显示了施加不同荷载后，花岗岩试样的整个破坏过程的声发射定位结果。

从破坏的形态上观察，岩样是由于纵向的（加载方向）裂缝引起的破坏，但是与完整岩样有明显区别的是，裂缝在孔洞两侧交叉、汇集。

由定位图中可以看出，孔洞的四周是声发射密集的地方。

图5是试验结束时事件的照片，它显示了最后的破坏组成。

可以观察到花岗岩的破坏形态非常散碎。

我们可以明显的看到缺口的损伤，
3 含孔洞岩石单轴试验数值模拟
3.1 颗粒流程序 PFC 计算思想
PFC2D (Particle Flow Code in 2 Dimension) 即二维颗粒流程序，是通过离散单元方法来模拟圆形颗粒介质的运动及其颗粒间的相互作用[1, 12]。

本次先模仿具有不同颗粒的全边界问题，然后把本构行为植入到模型中去，这样，本构反应就由大量的遵从简单本构关系的颗粒间的相互反应体现出来了[13-14]。

需要注意的是，PFC2D 中涉及的n K ，s K ，fric
等细观力学参量表征的是颗粒间的力学行为。

通过程序自带的FISH 语言，PFC2D 程序在细观方面可以监测试验过程中试样内部各点的接触力、应力–应变及微裂隙数量等物理量的变化，宏观方面可以观察到试样表面裂隙的发育情况、试样整体的位移及试样宏观的应力–应变关系等。

本次数值模拟的创新之处是在颗粒装配产生初始应力时采用接触刚度模型；加载过程中采用滑动模型和平行连接模型，以更好的模拟花岗岩材料。

3.2 数值模拟
以室内试验为基础，进行细观颗粒流模拟。

PFC 试样是根据给定的粒径大小，粒径比及孔隙率，按照随机分布规律排列的一些颗粒构成，颗粒通过摩擦忽略不计的四道墙体来约束（见图6）。

首先按照实际试
样的尺寸，建立0.1m×0.2m 基本模型，在生成PFC 试样时设定颗粒试样是由不同半径的颗粒单元组成，颗粒半径的分布采用从max R 到min R 的正态分布，经过大量PFC 试样的模拟试算，选定花岗岩PFC 试样粒径min R =1.75mm 到max R / min R =1.66。

PFC 模型中主要的控制参数见表1。

表1 花岗岩PFC 模型基本参数
颗粒比重/ kg/m 3
最小颗粒半径/mm 粒径放大系数颗粒接触模量/GPa 花岗岩 2630 1.75 1.66 21
a ）墙体
b ）墙体+颗粒（未平衡）
c ）无孔洞颗粒试样（平衡）
d ）含孔洞颗粒试样（平衡）
图6 颗粒模型图
模型建好以后，为了与室内单轴试验曲线特征相匹配，需要进行一系列的PFC 数值模拟试验，通过反复调整PFC 模型的输入参数，直到数值试验结果与实际物理模型试验结果基本一致。

PFC2D 程序里的参数有：颗粒接触模量c E 、颗粒接触法向和切向刚度n K 和s K 、颗粒的法向连接强度和切向强度n T 和s T 、颗粒间的摩擦系数μ等。

需要说明的是，为了确定在程序校准中需要使用的微观力学参数，本次数值分析中使用宏观力学参数中的杨氏模量、泊松比和单轴压缩强度。

经过反复调整，单轴试验数值模拟各参数见表2。

图7为数值模拟应力应变曲线。

表2 单轴压缩试验数值模拟参数
接触模量
/Ec a MPa 法向与切向刚度比率s n K K /
摩擦系数μ 法向连接强度a MPa 切向连接强度a MPa 花岗岩颗粒体 21.0 2.0 0.50 64.0 64.0
应力 /M P a 应变
图7 单轴试验数值模拟应力-应变曲线
在数值模拟的最后，通过测绘PFC2D 试件的裂缝可以辨别出试件的损伤区。

在破坏过程中随着裂纹的出现，我们观察到裂纹是由压应力比较大的空洞边缘开始发展的 [15-19]。

持续加载过程中在缺口处会形成更多的损伤。

缺口的进一步屈服刺激着试件两边破裂区附近裂纹的出现。

最后试件表现出一个宏观连续的裂缝区。

含孔洞岩石含孔洞岩石满足岩石脆性破坏四阶段，即裂隙闭合阶段；线弹性变形阶段；裂隙稳定扩展阶段；裂隙失稳扩展阶段。

图8 中分别显示了单轴压缩模拟试验中微裂缝的变化，以及受压粒子的运动方向。

其中受压后微裂缝的变化与试验得到的声发射定位图相对应，孔洞岩样在达到70%峰值强度之后，孔洞的水平方向出现裂缝，说明有明显的声发射现象出现。

图9是数值模拟试验最终试样图，它显示了最后的破坏组成。

与室内试验相比，同样我们可以明显的看到缺口的损伤，在应力集中洞的侧边界也形成了明显的两条缺口。

a ）PFC 试件受压产生
微裂缝 b ）PFC 试件受压粒子的 运动方向
图8 单轴压缩条件下产生的微裂缝及粒子运动方向
图9 PFC 试样的破坏形式
4 结果及讨论
通过试验及数值模拟试件损伤过程，刚开始试件受力后， 由于在初始变形阶段受力较小， 一般看不出明显变化。

但随着应力的增加，内部结构发生明显变化， 或有新的细观裂纹产生，于是可以很清晰的分析试样的破坏趋势。

从而得到其应力-应变曲线。

图10为数值模拟曲线与室内试验结果对比，可见PFC2D 试样的应力-应变曲线和实际试验测试结果基本吻合。

10
20
30
405060
70
80
应力 /M P a 应变
图10 PFC 试验结果与实测曲线对比
同样，对照PFC2D 试样模型得出的数值模拟结果和试验观察到的结果，不难发现，其结果很相符，破坏荷载以及裂纹模式都和试验结果也很吻合。

除此之外，在数值模型中微裂纹模式清晰地显示了出现于空洞两侧的相互竞争的裂缝，这和试验室试验中由声发射确定的位置是一致的。

综上，与室内试验结果对比分析，PFC2D 基本再现了脆性岩石的局部受压破坏现象，并且在破坏过程中可以得到其应力应变关系曲线，及其裂缝扩展过程图，通过声发射技术来监测微裂缝的演化过程。

通过改变计算模型中颗粒单元的性质，给出了在不同颗粒单元参数时试样的宏观性质，其结果对岩体的本构关系的完善有一定的应用价值。

参考文献
[1] CUNDALL P A, STRACK O D L. A discrete numerical model for graunlar assemblies [J ]. Geotechnique,
1979, 29(1): 47-65．
[2] POTYONDY D, AUTIO J. Bonded-particle simulations of the insitu failure test at Olkiluoto. in:
ELSWORTH D, TINUCCI J, HEASLEY K, editors. Rock mechanics in the national interest, Proceedings of
the 38th US Rock Mechanics Symposium, 2001,pp. 1553-60．
[3] POTYONDY D, CUNDALL P. The PFC model for rock: predicting rock-mass damage at the underground
research laboratory. Itasca Consulting Group, Inc. Report no. 06819-REP-01200-10061-R00,2001．
[4] 侯佳男．含孔洞岩石损坏破坏过程声发射试验研究[D]．沈阳：东北大学，2009．
[5] 陶纪南，张克利，郑晋峰．岩石破裂过程声发射特征参数的研究[J]．岩石力学与工程学报，1996（15）
增刊：452-455．
[6] 秦跃平，孙文标，王磊．岩石损伤力学模型分析[J]．岩石力学与工程学报，2003，22（5）：702-705．
[7] TANG C. A. Numerical Simulation of rock failure and associated seismicit [J]. Int. J. Rock Mech. Min. Sci.,
1997, 34: 249-262．
[8] 赖勇，张永兴．细观损伤理论确定岩体参数的方法研究[J]．岩石力学与工程学报，2008（27）增刊：
2617-2623．
[9] TANG C. A. Numerical Simulation of Rock Failure and Associated Seismicit [J], Int. J. Rock Mech. Min. Sci.,
1997, 34: 249-262．
[10] 张薇，薛嘉庆．最优化方法[M]．沈阳：东北大学出版社，2004．
[11] 胡新亮，马胜利，高景春，等．相对定位方法在非完整岩体声发射定位中的应用[J]．岩石力学与工
程学报，2004，23（2）：277–283．
[12] CUNDALL P A, STRACK O D L. Particle flow code in 2 Dimensions[A]. Itasca Consulting Group Inc.,
1999．
[13] 刘顺桂，刘海宁，王思敬，等．断续节理直剪试验与PFC2D数值模拟分析[J]．岩石力学与工程学报，
2007，27（9）：1828–1836．
[14] Itasca Consulting Group Inc.. PFC2D (particle flow code in 2 dimensions) theory and background[R].
Minnesota, USA: Itasca Consulting Group Inc., 2002．
[15] HOEK E. Brittle failure of rocks in rock mechanics in engineering practice. In: Stagg KG, Zienkiewicz
OC, editors. Rock mechanics in engineering practice. New York: Wiley; 1968. pp. 99–124．
[16] BRACE W F, PAULDING B, SCHOLZ C. Dilatancy in the fracture of crystalline crocks. J Geophys Res
1966;71(16): 3939–53．
[17] BIENIAWSKI Z T. Mechanism of brittle fracture of rock, parts I, II and III. Int J Rock Mech Min Sci
Geomech Abstr 1967, 4(4): 395–430．
[18] LOCKNER D A, BYERLEE J D, KUKSENKO V, et al. Observation of quasi-static fault growth from
acoustic emissions. In: Evans B, Wong T-f, editors. Fault mechanics and transport properties of rocks.
New York: Academic Press; 1992. pp. 3–31．
[19] MARTIN C D, CHANDLER N A. The progressive fracture of Lac du Bonnet granite. Int J Rock Mech Min
Sci 1994, 31: 643–59．。