巨石混合体力学性质的数值试验分析

刘 蕾，闫金凯，陈 亮，倪天翔

(1.长安大学地质工程与测绘学院，陕西西安 710054;2.中国地质环境监测院，北京 100081)

0 引言

巨石混合体是指由尺寸较大的块石构成的堆积体，它在我国西南山区广泛分布，通常是滑坡、崩塌等地质灾害的产物，在山区中易形成较为平坦的地形地貌，常被用作为城镇建设的场地［1］。一些国家重大工程的建设也经常与巨石混合体有关，如中国科学院国家天文台“500米口径球面射电望远镜”(FAST)工程建设中就遇到巨石混合体边坡开挖稳定性的问题［2］。由于巨石堆积的结构特征，对巨石混合体边(滑)坡的稳定性评价方法目前尚存在争议。

要对巨石混合体边(滑)坡的稳定性进行准确的评价，首先需要明确巨石混合体的力学性质及参数取值。构成巨石混合体的块石尺寸一般为一米至数米，块石之间互相镶嵌、咬合，导致其宏观力学性质与一般的土体或岩体有很大的差异。

目前对于巨石混合体的力学性质及参数取值方面的研究成果极少，但与之结构特征相近的一些地质体或工程材料的研究方法可作为参考。对于松散堆积体力学性质的研究方法主要有力学试验法［3-8］和数值模拟法［9-11］。力学试验法主要包括野外原位直剪试验、室内大型三轴压缩试验，对于巨石混合体来说，由于其块石尺寸较大，在目前的试验仪器条件下几乎无法进行原位直剪试验;而大型压缩三轴仪的直径一般为300mm，要求试样的最大颗粒直径为60mm，所以利用室内大型三轴压缩试验对巨石混合体的力学性质进行研究时，需进行较大比例的缩尺，试验所得的结果能否真实反应巨石混合体的力学性质缺少相应的验证。数值模拟法不受尺寸的限制，可用于巨石混合体力学性质的研究，但需充分考虑块石形状、块石尺寸、粒径级配、块石的堆积形态等多种因素。

本文通过构建巨石混合体随机块体模型，采用数值模拟的方法对巨石混合体进行双轴压缩数值试验，对巨石混合体的力学性质进行研究。

1 巨石混合体块体模型随机生成技术

本节利用蒙特卡罗随机模拟方法，建立巨石混合体块体的随机生成技术，为后续巨石混合体的力学数值试验提供基础。

本方法采用蒙特卡罗法模拟块石的随机分布，通过抽样的方法获取块石级配、尺寸、形状，运用随机数值模拟技术来生成块体的空间分布形式，进一步构成各种形态的块体单元。

随机块体的生成需要考虑以下几个方面:

①块体形状及边数;

②长轴倾角;

③粒径大小及各尺寸的百分比(即粒径级配);

④块体位置坐标。

基于以上几点考虑，本方法的基本思路是:首先在空间中随机生成直径大小不同的圆，再在圆上截取不同数量、不同位置的点作为块体的角点，连接各角点后即为生成的块体模型。

(1)块体的随机分布

考虑到块体在空间上是随机均匀分布的，因此采用均匀分布模型来生成每个圆的X、Y值，其中，X为块体多边形外接圆的圆心横坐标;Y为多边形外接圆的圆心纵坐标。这样，在圆内生成多边形后就实现了块体在空间上的随机分布。

(2)块体的生成

由于巨石混合体受岩性、原始岩体的节理裂隙分布情况、成因及搬运过程的影响，块体的形状及边数也多种多样。因此在生成块体时，需根据所研究巨石混合体的块体几何特征确定要生成的多边形块体的边数N，然后将各个圆均匀分成N份，在每份圆周上随机生成一个角点，连接各角点生成边数为N的随机块体模型。例如，要生成的块体为四边形，则将圆周分成四份，在每份圆周上随机生成一个角点(A、B、C、D)，则连接四个角点可以构成一个四边形(图1)。由于角点选取的随机性，保证了四边形的形状是随机的，同时四边形的长轴倾角(体现为长轴与X方向的夹角)也是随机的。

图1 四边形随机块体生成示意图Fig.1 Schematic diagram of random quadrilateral block

(3)块体相互侵入关系的判别

在生成块体时，要保证各块体之间不会发生相互侵入的状况。本方法中由于多边形在圆形内部生成，因此只要保证圆与圆之间不发生相互侵入即可。如图2所示，圆与圆相互侵入关系的判别准则为:d≥r1+r2，式中，d为两个圆心之间的距离，r1、r2分别为两个圆的半径。即任意两个圆的圆心之间的距离不小于这两个圆的半径之和，这样可以保证任何两个圆之间都不发生侵入，圆内生成的多边形相互之间也不会发生侵入(图2)。

图2 随机圆形侵入关系示意图Fig.2 Discrimination of random round separation

根据上述方法得出的随机块体之间必定存在较大空隙，块体之间为完全不接触的关系，可将块体模型转入相应的离散元计算软件，使块体在重力作用下自然堆积，所得的模型即为考虑块体形状、尺寸、空间分布、堆积形态等多因素在内的巨石混合体随机块体模型。

2 基于随机块体的巨石混合体双轴压缩数值试验

2.1 计算模型的建立

本次数值试验以某巨石混合体为例建立计算模型。根据现场调查的结果，巨石混合体中块石的粒径级配曲线如图3所示。在生成块体模型时，将块石尺寸概化为4m、3m、2m、1m四个代表粒径，各粒径对应的百分含量如表1所示。块体形状采用了边数分别为4、5、6的三种多边形，各多边形的数量及位置采用随机均匀分布的方法生成。

图3 巨石混合体粒径级配曲线Fig.3 Particle size distribution curves of large rock mixture

表1 随机生成块体粒径及百分含量表Table 1 Particle size and percentage of random block

利用随机块体生成软件进行块体的生成，然后导入计算软件中进行块体的堆积，最后生成的巨石混合体双轴压缩试样如图4所示，模型尺寸为20m×40m。

图4 巨石混合体双轴压缩模型Fig.4 Biaxial compression model of large rock mixture

2.2 计算条件

本次数值试验利用CDEM软件进行计算，将试样中块体边界设为接触面，块体材料采用莫尔-库仑模型，试样中接触面及块体内部网格划分情况如图5所示。块体材料参数及接触面参数如表2、表3所示。

图5 接触面及网格划分Fig.5 Contact surface and grid

表2 块体物理力学参数表Table 2 Recommended physical and mechanical parameters of the block

表3 接触面计算参数表Table 3 Recommended physical and mechanical parameters of the surface

计算时，将模型的底部固定，两侧面施加不同程度的围压(围压分别为 0.2MPa、0.5MPa、1MPa、2MPa、3MPa)，对模型顶部施加竖向和水平向的位移约束，使试样在围压作用下进行压缩，待稳定后，在模型顶部施加竖直向下的均匀速率使模型压缩，竖向速率为1×10-6m/s。模型的边界条件如图6所示。压缩过程中，在接近模型顶部的位置选取一条水平向的测线，在测线上等间距选取若干个监测点(本次试验选取的监测点为20个)，分别监测测线上各个部位的位移及应力。

图6 模型边界条件示意图Fig.6 Boundary conditions of the model

3 计算结果与分析

3.1 变形与破坏规律

图7为围压为2MPa时模型的竖向位移云图。从图中可知，模型在压缩过程中体现出明显的上部位移大、下部位移小的特征。模型受压后，首先是块体之间被压密，块体发生滑动、转动变形，由于侧向围压的作用，块体之间的变形受到一定的限制，到一定程度后模型完全被压密。此后随着压缩的继续进行，块体应力不断增大，逐渐发生块体内部的变形破坏。从图8所示的模型塑性剪应变云图来看，发生变形破坏的块体主要位于模型中部。由于巨石混合体块石堆积的不均匀性，导致双轴压缩情况下剪切破碎带不像均质连续材料那样规则，但总体上剪切破坏带大致仍按X型分布。

从图9所示的不同围压下巨石混合体模型的破坏情况来看，巨石混合体在双轴压缩情况下存在较明显的剪胀现象，并且随围压的增大，剪胀程度呈减小的趋势。其原因在于，在低围压条件下巨石混合体的块石之间接触相对松散，块石强度较块石之间的咬合力和摩阻力要大得多，低围压条件下剪切过程中主要是克服块石之间的咬合力和摩阻力。块石沿着相互之间的接触面滑移时，必然引起模型体积的增大，从而表现出剪胀特性。而随着围压的增大，块石接触越加紧密，咬合力和摩阻力也相应增大，与块石自身强度的差值相对减小，剪切过程中不仅仅是克服块石之间的咬合力和摩阻力，还伴随着块石自身的变形，从而表现出剪胀程度减小。

图7 围压2MPa时模型竖向位移云图Fig.7 Vertical displacement of the model when the confining pressure is 2 MPa

图8 围压2MPa时模型塑性剪应变云图Fig.8 Plastic shear strain of the model when the confining pressure is 2 MPa

图9 不同围压条件下模型竖向位移云图Fig.9 Vertical displacement of the model in different confining pressure

3.2 应力－应变关系

图10为主应力差-轴向应变曲线。从图中可知，各曲线存在不同程度的波动现象，其原因可能与剪切过程中巨石混合体块石之间不断发生变形，堆积形态不断发生变化有关。总体上讲，在低荷载作用下巨石混合体的应力应变曲线呈直线状，表现出弹性变形的特点，随着荷载的进一步增大，模型迅速屈服，应力应变曲线下凹。其原因在于，巨石混合体由块石堆积而成，块石之间互相镶嵌、咬合，荷载较小时不足以克服块石之间的咬合力和摩阻力，巨石混合体的变形主要表现为块石自身的弹性变形;而当荷载进一步增大，块石之间就可能产生相对的滑移，从而导致应力-应变曲线迅速下凹。

在相同的轴向应变条件下，巨石混合体的应力与围压有关，围压越高应力越大。在低围压条件下(0.2MPa)巨石混合体的应力-应变曲线有较明显的峰值强度，模型破坏时有一定程度的应力降，表现出应变软化特性;在高围压条件下(围压＞0.5MPa)应力-应变曲线没有明显的峰值强度，模型屈服后随轴向应变的增大持续增大，表现出应变硬化特性。其原因可能在于，块石产生相对滑移后，在荷载作用下与滑移方向上的块石接触，产生新的镶嵌、咬合，因而使得模型在剪切破坏过程中没有明显的峰值强度，表现出应变硬化的特点，而此种现象在围压越高时越明显。

图10 巨石混合体应力—轴向应变曲线Fig.10 Stress－axial strain curves

3.3 抗剪强度

根据各级围压下模型的轴向峰值强度(没有明显峰值强度时，取轴向应变为10%对应的轴向应力值)，可绘制模型在不同围压条件下的莫尔圆(图11)。

图11 莫尔圆和强度包线Fig.11 Mohr's circles and Mohr-Coulomb failure envelopes

由图11可知，较低应力条件下巨石混合体的强度包络线呈直线状，符合莫尔-库仑强度理论;而应力水平较高时，强度包络线向下弯曲，呈下凹形状，说明在较高的应力水平下巨石混合体剪切破坏时并不遵从莫尔-库仑强度理论。对于这种具有非直线形强度包络线的岩土体的抗剪强度及其抗剪强度参数的评价，目前存在多种方法和观点。对于巨石混合体滑坡来说，由于坡体的厚度一般不会太大(通常几米或十几米)，所以巨石混合体通常处于较低的应力水平下。因此，在进行巨石混合体滑坡的研究中，可认为巨石混合体的强度符合莫尔-库仑理论，其抗剪强度参数可按莫尔-库仑理论表达。

通过对包络线直线段的拟合，可得出巨石混合体在较低应力条件下的抗剪强度参数。为分析巨石混合体的抗剪强度参数与块石表面参数的关系，分别取块石表面摩擦角为 35°、30°、25°进行了计算，计算结果如表4所示。

表4 巨石混合体抗剪强度参数表Table 4 Shear strength parameters of huge rock deposits

从表4中可知，计算时认为块体表面不存在粘聚力，而计算所得的巨石混合体的综合抗剪强度参数中C值不为零，这里的C值可以理解为块体之间产生的咬合力。计算所得的巨石混合体内摩擦角比块体表面的摩擦角普遍高10°左右。随着块体表面摩擦角的降低，巨石混合体的C、Φ值均有所减小。

4 结论

本文建立了基于蒙特卡罗法的巨石混合体随机块体生成技术，进行了基于随机块体的巨石混合体双轴压缩数值试验，对巨石混合体的力学性质及参数取值进行了研究，主要结论如下:

(1)在双轴压缩条件下，巨石混合体表现出高围压条件下的应变硬化和低围压条件下的应变软化特征;

(2)巨石混合体在剪切破坏时表现出剪胀特性，在低围压条件下剪胀更明显;

(3)块石堆积的结构特征导致巨石混合体的剪切破坏带具有不规则性，主要位于试样中部，大致呈X型;

(4)巨石混合体的强度包络线在低围压下(＜0.5MPa)为直线，高围压下(＞0.5MPa)呈下凹型。巨石混合体滑坡中的巨石混合体通常处于低应力条件，可按强度包络线的直线段求取抗剪强度参数，认为其符合莫尔-库仑强度理论。

(5)巨石混合体的抗剪强度参数与块体表面的摩擦系数有关，块石表面的摩擦角越大则巨石混合体的内摩擦角越大，通常巨石混合体的内摩擦角比块石表面摩擦角大10°左右。

［1］殷跃平，李廷强，唐军.四川省丹巴县城滑坡失稳及应急加固研究［J］.岩石力学与工程学报，2008，27(5):971-978.YIN Yueping， LI Tingqiang， TANG Jun.Landslide reactivation and emergency stabilization on Danba county town in Sichuan province［J］.Chinese Journal of Rock Mechanics and Engineering，2008，27(5):971 – 978.

［2］闫金凯，马娟，冯春.FAST台址巨石混合体边坡开挖稳定性分析［J］.中国地质灾害与防治学报，2012，23(2):25-29.YAN Jinkai，MA Juan，FENG Chun.Study on stability of large rock mixture slope under excavation in FAST location［J］.The Chinese Journal of Geological Hazard and Control，2012，23(2):25-29.

［3］Marsal R.J.，Large scale testing of rockfill materials，J.Soil Mech.Fdn Engng Div.Am.Soc.Civ.Engrs.93，1967，SM 2:27-34.

［4］Marachi N.D.，Strength and deformation characteristics of rockfill materials，PhD thesis，University of California，1969.

［5］R.J.Marsal， Mechanical properties of rockfill，Embankment dam engineering，casagrande volume，New York，1973:109-200.

［6］Indraratna B.et al.Design for grouted rock bolts based on the convergence control method，Int.Rock Mech.Sci.Geomech.Abstr.27，1986，No.4:268-281.

［7］秦红玉，刘汉龙，高玉峰，等.粗粒料强度和变形的大型三轴试验研究［J］.岩土力学，2004，25(10):1575-1580.QIN Hongyu， LIU Hanlong， GAO Yufeng， etal.Research on strength and deformation behavior of coarse aggregates based on large-scale triaxial tests［J］.Rock and Soil Mechanics，2004，25(10):1575-1580.

［8］涂国祥，黄润秋，邓辉，等.某巨型冰水堆积体强度特性大型常规三轴试验［J］.山地学报，2010，28(2):147-153.TU Guoxiang，HUANG Runqiu，DENG Hui，et al.Study on the strength and deformation behabior of a huge outwash deposits based on large-scale triaxial tests［J］.Journal of Mountain Science，2010，28(2):147 – 153.

［9］李世海，汪远年.三维离散元土石混合体随机计算模型及单向加载试验数值模拟［J］.岩土工程学报，2004，26(2):172-177.LI Shihai，WANG Yuannian.Stochastic model and numerical simulation of uniaxial loading test for rock and soil blending by 3D-DEM［J］.Chinese Journal of Geotechnical Engineering，2004，26(2):172-177.

［10］徐文杰，胡瑞林，岳中琦，等.土石混合体细观结构及力学特性数值模拟研究［J］.岩石力学与工程学报，2007，26(2):300-311.XU Wenjie， HU Ruilin， YUE Zhongqi， etal.Mesostructural characterand numerical simulation of mechanical properties of soil-rock mixtures［J］.Chinese Journal of Rock Mechanics and Engineering，2007，26(2):300-311.

［11］周伟，常晓林，周创兵，等.堆石体应力变形细观模拟的随机散粒体不连续变形模型及其应用［J］.岩石力学与工程学报，2009，28(3):491-499.ZHOU Wei，CHANG Xiaolin，ZHOU Chuangbing，et al.Stochastic granule discontinuous deformation model og rockfill and its application［J］.Chinese Journal of Rock Mechanics and Engineering，2009，28(3):491 – 499.