土石混合体单轴压缩试验的三维颗粒流数值模拟

金　磊，曾亚武，叶　阳



土石混合体单轴压缩试验的三维颗粒流数值模拟

金磊，曾亚武，叶阳

（武汉大学土木建筑工程学院，武汉430072）

为探究不同含石量土石混合体在单轴压缩荷载下的变形破坏机理，在不规则颗粒三维离散元精细建模技术的基础上，提出了一种不规则块石三维离散元模型随机生成技术，并建立了符合宏观统计规律的土石混合体三维颗粒流数值模型。进行了低、中、高不同含石量（0，10%，30%，50%，70%，90%）的土石混合体单轴压缩试验颗粒流数值模拟，获得了各种含石量的土石混合体在单轴压缩荷载下的应力应变特征和变形破坏特征，并对其变形破坏的细观机理进行了深入的分析探讨。结果表明：随着含石量的增加，土石混合体试样单轴抗压强度和弹性模量随之减小；中低含石量（0，10%，30%，50%）土石混合体的应力-应变曲线表现出典型的应变软化特征，而高含石量（70%，90%）土石混合体在达到峰值应力后表现出一定程度的塑性流动特征；中低含石量土石混合体试样中形成的微裂纹主要是张拉裂纹，宏观上表现出来的是拉裂破坏，而高含石量土石混合体试样中微裂纹很少，主要是剪切裂纹，宏观上表现出来的是剪切破坏。

土石混合体；不规则颗粒；三维颗粒流；单轴压缩试验；离散元模型

doi：10．11988/ckyyb．20150423

1　研究背景

土石混合体是指第四纪以来形成的，具有一定工程尺度、强度较高的块石、细粒土体及孔隙构成且具有一定含石量的极端不均匀松散岩土介质系统［1］。土石混合体在自然界中分布非常广泛，地球上绝大多数滑坡体都是由土石混合体构成，研究土石混合体变形破坏机理对于进一步完善土石混合体力学理论体系和防治滑坡等地质灾害具有重要的意义。

针对土石混合体这种特殊地质体，国内外许多学者对其物理力学特性、渗透特性和变形破坏机制等进行了较为系统的研究，取得了一系列研究成果。有关土石混合体的研究主要包括土石混合体原位地质调查［2-3］、大型原位实验［4-5］、室内试验［6-7］和数值模拟［8-10］。在土石混合体单轴压缩试验研究方面，WANG Yu等［11］进行了土石混合体室内重塑样单轴压缩试验，并采用CT扫描技术分析了加载过程中裂纹发展特征；Mohammad Afifipour等［12］用砾石和水泥制备了土石混合体的重塑样，进行了高含石量土石混合体的单轴压缩试验，分析了其力学特性和破坏特征；廖秋林等［13］提出了一种基于数码图像的土石混合体结构模型自动生成方法，利用有限元分析了土石混合体在单轴压缩下应力场的结构效应；李世海等［14］基于一种可变形块体离散元建立了土石混合体的三维离散元随机模型，研究了单轴压缩下土石混合体内部应力场分布与含石量和块石大小的关系。

由于天然状态下的土石混合体大多是弱胶结甚至是无黏结的松散岩土体，现场获得未扰动的原状样将变得极其困难。一方面，室内重塑样可以在一定程度上还原土石混合体的某些结构组成特征，但对于单轴压缩试验已有成果都是采用胶结比较强、碎块石粒度较小的试样，与实际的土石混合体有一定的差别。另一方面，土石混合体具有典型的非均质和非连续性，在力学性质上表现出强烈的各向异性，采用传统的建立在宏观连续介质上的力学模型和相应的分析方法，难以合理描述和反映其内部结构的复杂性，也难以从细观水平上揭示其变形破坏机理。离散元法则能较好解决上述难题，基于基本球体的颗粒离散元所采用的接触关系较之块体离散元有了很大的简化，接触搜索耗时大大降低，计算效率显著提高，有着明显的优势。

要从细观水平上探究土石混合体的变形破坏机理，一个关键问题就是要建立能够体现其组成结构特征的细观结构数值模型。目前土石混合体细观结构模型的生成方法主要包括基于数字图像处理的模型生成技术和基于统计分析层次上的随机模型生成技术。前者是研究原位实际情况的一种有效的手段和方法，但是其也存在一些不足，例如很难对所有研究区域进行拍照，很难找到为某种研究目的而需要的理想剖面［2］；而且图像处理得到的仅是露头处的局部的平面情形。后者是根据土石混合体结构组成在宏观层次上的统计规律性，建立其相应的细观结构模型，更适合于开展机理性研究，如研究不同含石量对土石混合体物理力学特性和变形破坏机理的影响等。

鉴于上述理由，建立土石混合体三维离散元随机模型，采用颗粒离散元进行数值模拟和分析是从细观水平上研究土石混合体变形破坏机理的一个比较合适的方法。而从已有成果可以看出，目前对土石混合体单轴压缩试验的颗粒离散元模拟还未见报道，其在单轴压缩荷载下变形破坏的细观机理还不为人所知。为此，本文提出了一种不规则块石三维离散元模型随机生成技术；根据土石混合体中块石形状、大小、位置和方向的统计特征，建立一定程度上符合宏观统计规律的土石混合体三维离散元随机模型；对土石混合体单轴压缩试验进行颗粒流数值模拟，重点分析各含石量（0，10%，30%，50%，70%，90%）土石混合体的应力应变关系特征、破坏特征和微裂纹特征，从多个方面、多个尺度上对各含石量的土石混合体的变形破坏机理进行深入的分析，获得了一些有意义的结论。

2　不规则块石数值模型

土石混合体中块石颗粒形状和大小在一定程度上反映了其形成历史，土石混合体的宏观力学行为取决于细观水平上颗粒之间的相互作用，而颗粒形状是影响颗粒间相互作用的一个非常重要的因素。已有颗粒材料的研究结果表明：颗粒形状越不规则，自然堆积状态下孔隙比越大、整体刚度越小、临界摩擦角越大、对应力状态的变化越敏感［15］。因此，要准确模拟土石混合体的力学性质和运动特性，必须首先实现不规则块石的模拟。

2．1不规则块石随机几何模型

目前对随机形状块石的模拟多见于混凝土中粗骨料的研究，已有的建模方法基本上可分为2种：一种是通过将随机生成的三角形或多边形或具有简单数学表达的曲面作为基，然后按一定算法进行延拓生成骨料，并通过凹凸判断来控制骨料为凸形，最后通过冲突判断来实现骨料的无重叠投放［16］；另外一种方法则是将基本几何体用多个平面进行随机切割来获得随机多面体模型［17］。已有的这些随机块体建模方法的算法都比较复杂，效率也都不高，很难生成高含石量的模型，而且随机块体都是凸体，而实际的土石混合体滑坡堆积体中的块石体更多表现的是凹凸并存，因此有必要提出一种新的适用于模拟土石混合体中块石形态特征的建模方法。

目前国际上对不规则颗粒形态的描述主要是分3个尺度来进行。第1个是颗粒尺寸尺度上的，可以定义为球度，它反映了颗粒的总体形状，表示颗粒三向尺寸（长、宽、高）之间的比例；第2个是比颗粒尺寸小一个数量级尺度上的，可以定义为凹凸度，它反映了颗粒表面的宏观凹凸程度；第3个是细观尺度上的颗粒表面的粗糙度。为了更好地模拟土石混合体中块石体的随机几何形态，本文提出了一种基于APDL（ANSYS参数化设计语言）参数化建模的随机块石体建模新技术。现将该技术的主要原理和步骤概述如下：

（1）随机球度基面生成。将实际块石体三向尺寸中最大的作为建模时x向长度，根据另外2个方向的长度与x向长度之比得到y向长度和z向长度（规定y向长度比z向长度大）。根据得到的x，y，z向尺寸即可建立一个中心在原点的长方体，删除长方体，保留组成长方体的6个面，即可得到规定球度的基面。应用时可根据实际统计得到的块石三向尺寸的特征，设定三向尺寸的随机数变化范围，由此可得到符合实际统计规律的随机球度的基面。

（2）随机凹凸面生成。将第1步得到的球度基面划分网格，单元形状宜选择三角形壳单元，单元尺寸可根据颗粒表面凹凸点个数确定，如较少（4个左右）则可将各边的分段个数设为1，以此类推。划分网格后记录各三角形面单元的组成节点信息，得到各个节点的球坐标，然后对各节点球坐标中的第1个分量（该节点到原点的距离）进行随机比例的缩小变换，随机数的取值范围可根据实际颗粒表面凹凸程度试算对比来确定。接下来根据前面记录的三角形面单元与其组成节点的关系以及随机变换后的新节点坐标重新生成面，即可得到需要的随机凹凸度的曲面。

（3）随机几何模型生成。根据变换后得到的随机曲面围成的空间即可建立随机几何体模型，将单元类型改为四面体单元，对生成的几何体进行体网格划分，即可得到满足后续离散元建模要求的几何模型。

根据已有土石混合体中块石形态特征的统计调查并经试算对比，本文在块石随机模型的生成过程中进行了一定的简化，将块石颗粒y和z向尺寸设为相同，y和z向长度与x向长度比值设为在0．3～0．8的范围内随机取值，初始随机基面网格划分尺寸设为每个边的分段数为2，各节点球坐标第1分量缩小比例在0．4～1．0的范围内随机取值。在上述参数取值情形下，由第1步得到的某个随机球度的基面（x向长度为40 mm）如图1（a）所示，由第2，3步得到的随机几何模型如图1（b）所示。

图1　典型不规则块石几何模型Fig．1　Typical geometric model of irregularly shaped rock block

从以上建模方法和结果可以看出，本文提出的建模方法简单实用，只需要经过3步操作即可获得不规则块石的随机几何模型。至于不规则颗粒表面的粗糙度，几何模型中暂不用考虑，因为一方面可以通过数值模拟时对颗粒设定不同的细观摩擦系数来表示，另一方面由基本球体填充后建立的离散元数值模型的表面是不平整的，可近似模拟粗糙度。

2．2不规则块石离散元模型及其密度优化

根据2．1节方法所得到的不规则块石几何模型，基于不规则颗粒三维离散元精细建模技术［18］，可以建立不规则块石的三维离散元模型。其主要步骤归纳如下：

（1）针对所获得的已划分网格的不规则块石几何模型，导出模型的表面节点坐标和内部节点坐标信息。

（2）对于任意一个内部节点，计算其到每个表面节点的距离，只记录并存储其中的最小值。对所有内部节点都进行一次计算，即可得到所有内部节点到表面各自对应的最小距离值。然后以内部节点为球心，以对应的最小距离值为半径生成球体。所有的球体组合就能刻画出不规则块石的外形。

（3）根据研究的需要，引入最小球体半径Rlim和有效空间系数K来对基本填充球体数量进行优化，以节省计算资源。

图1（b）中建立的不规则块石几何模型划分网格后有2 314个表面节点和6 671个内部节点，则初始填充球个数为6 671个，如图2（a）所示，这样生成颗粒集合体时球体数量将会非常多，造成计算效率低下，因此需要引入控制参数进行优化。取最小球体半径Rlim为2 mm，有效空间系数K为1．10，则优化后的离散元模型如图2（b）所示，该离散元模型由544个球体构成。由图可见，优化后的模型填充球体数大大减少，而模型形态特征基本相同，可以大大节省计算资源。

图2　典型不规则块石三维离散元模型Fig．2　Typical 3-D DEM model of irregularly shaped rock block

基于基本球体重叠构建不规则块石的建模方法不可避免地会出现模型密度分布不均匀的现象，为此本文采用平均密度的方法使得模型的质量与实际块石的质量相等，用公式可以表示为：

式中：ρm为优化后各球体赋值密度；Vclump为模型体积；ρ0为颗粒实际密度；Vi，overlap为每个球体中重叠的体积。

3　土石混合体数值模型

3．1颗粒流数值平台（PFC3D）概述

三维颗粒流PFC3D是根据颗粒离散元理论，采用显示差分算法，交替运用力-位移定律和牛顿运动定律来模拟颗粒之间的相互作用和运动过程。PFC3D提供了3种基本的接触模型，它们分别是接触刚度模型、滑移模型和黏结模型。其中黏结模型分为接触黏结和平行黏结［19］。根据土石混合体的结构组成特征，本文选取线性接触刚度模型、滑移模型和接触黏结模型来进行模拟。

PFC3D中提供了clump即“团粒”以供不规则颗粒的模拟，它由基本的球体颗粒聚集而成，组成颗粒间相对位置保持固定，组成颗粒间的接触计算在循环的过程中被跳过，因此大大节约了计算时间。clump边缘颗粒与其外部颗粒产生相互作用，但是不管是否有接触力作用在其上面，clump不会解体［19］。土石混合体中的块石的强度和刚度较细粒土体基质要大得多，因此采用clump来模拟不规则块石是比较合适的。

3．2土石混合体数值试样

天然形成的土石混合体中块石颗粒粒径变化范围较大，为了尽可能使试样的块石粒径与实际一致，并考虑到计算效率，本文所建立的数值试样是直径为200 mm、高度为500 mm的大型圆柱形试样。根据Medley［20］研究，土/石阈值可定义为

dS/RT=0．05LC。（3）式中：dS/RT为土/石阈值；LC为土石混合体的工程特征尺度，对于三轴试验试样，可取为试样直径，因此本文的土/石阈值可定为10 mm。鉴于本文的研究目的是探究土石混合体的变形破坏机理，为了提高计算效率，本文对颗粒级配进行了一定程度的简化，土体颗粒直接简化为球体，粒径取值范围为8～10 mm。对于块石，其粒径简化为在20～40 mm范围内服从正态分布，这里所说的粒径并不是块石的实际粒径，而是指与该粒径的球体体积相同的当量粒径。采用第1节所述的土石混合体中块石的形态特征参数取值范围和随机建模方法建立20～40 mm范围内块石的随机离散元模型库，本次研究所建立的块石颗粒库包括100个粒径的随机颗粒，基本满足统计规律的要求。

建立随机块石颗粒库以后，就可以着手建立土石混合体数值试样了。首先要将块石投放到随机的位置上，这里提出一个简单实用的方法：先利用PFC3D中generate命令生成的球体不会重叠的特点，在比试样较高一些的空间中生成数目较大、粒径也较大的球体（为了不让投放的块石重叠），则这些球体的球心坐标是随机的，将这些坐标输出即可获得足够的空间随机点坐标；选取第一个随机点的坐标，随机选择颗粒库中的颗粒球体坐标文件，利用VC++语言编制程序将选择的颗粒平移到选取的相应的随机点上，同时进行随机旋转，这样就投放好了第1个块石颗粒；同理继续投放，直到块石颗粒的总体积达到需要的含石量；投放完成以后，在PFC3D平台上试运行块石颗粒文件，查看有无接触，若无接触即满足要求，若存在接触则需重新设定随机坐标的空间高度或者球体直径大小，重新投放。

块石投放结束以后，可以通过“半径膨胀法”生成符合级配要求的土体颗粒。这里根据所需的含石量（块石体积与土石总体积的比值）计算土颗粒和块石颗粒体积时，试样孔隙率假设为0，这样得到的初始试样是高度＞500 mm颗粒之间非常疏松的颗粒集合体。然后，让颗粒在较大的重力加速度（本文中取4 g）下自由下落到基本达到平衡，将位置高于500 mm的颗粒删除，循环至平衡，即可得到一个较为密实、块石分布与实际较为相似的土石混合体单轴压缩数值试样。本文研究建立了含石量为0，10%，30%，50%，70%，90%的土石混合体数值试样。这里给出含石量为0%，50%和90%的试样的颗粒流模型，如图3所示，土颗粒为蓝色，块石由不同颜色的clumps表示。

图3　土石混合体数值试样Fig．3　Numerical specimens of soil-rock mixture

3．3模型参数与试验方法

参考土石混合体相关的试验结果和其他数值模拟成果［1，8-10］，经过参数试算，本文选取了一组具有代表性的颗粒细观参数。其中，土颗粒之间的接触黏结法向强度和切向强度设为20 N，土颗粒与块石颗粒之间的接触黏结法向强度和切向强度设为2 N，块石与块石颗粒间无接触黏结。上下边界墙的法向刚度为1×109N/m，切向刚度和摩擦系数均为0，其它参数见表1。

表1　土石混合体细观力学参数Table 1　Meso-mechanical parameters of soil-rock mixture

采用PFC3D进行无侧限单轴压缩试验的数值模拟时，把它看成是围压很小、侧墙的刚度也很小（不影响其力学行为）的特殊的三轴压缩试验。经试算，这里围压设为100 Pa，侧墙的刚度设为1 000 N/m。试验时，首先使试样在0．4 MPa的围压下进行各向等压固结，然后按上述参数值施加颗粒间的黏结，以模拟实际边坡中一定深度处的土石混合体的初始状态；再利用伺服控制程序对试样施加100 Pa的围压，这时试样在黏结的拉力和颗粒间压力作用下达到自平衡状态，随后即可准备加载。模型的上下边界墙作为加载板，加载速度为0．001 m/s，这里需要指出颗粒流静力学数值模拟中的速度与真实的速度大小是不一样的［21］，采用应变控制加载方式。在加载过程中对轴向应力、轴向应变、微裂纹总数、剪切裂纹数、张拉裂纹数等量进行了采样和记录，以便后续进行分析。

4　数值试验结果分析

4．1应力-应变特征

根据数值试验结果，绘制各含石量土石混合体单轴压缩的应力-应变关系曲线，如图4所示。从图4可以看出，随着含石量的增加，土石混合体试样单轴抗压强度和弹性模量逐渐降低，而破坏应变逐渐增大，延性增强，这与Mohammad Afifipour等［12］总结的室内试验结果是一致的。与胶结土体（含石量为0）相比，含石量为10%的土石混合体单轴抗压强度和弹性模量降低不明显，含石量为30%时则开始有明显降低。低含石量（0，10%，30%）的土石混合体在单轴压缩荷载下表现出明显的应变软化特征：应力-应变曲线呈下凹形的初始压密阶段；压密后的弹性变形阶段，应力-应变曲线基本呈直线；稳定破裂发展阶段，应力-应变曲线开始出现波动；不稳定的破裂发展阶段，应力迅速降低；强度丧失和完全破坏阶段，应力维持在较低水平波动。中等含石量（50%）的土石混合体与低含石量的试样相比，单轴抗压强度和弹性模量均显著降低，其应力-应变曲线也表现出一定的应变软化特征。高含石量（70%，90%）的土石混合体单轴抗压强度基本维持在很小值（2 kPa左右），达到峰值应力后应力-应变曲线表现出一定程度的塑性流动特征。

图4　不同含石量试样应力-应变关系曲线Fig．4　Stress-strain relationship curves under different rock block proportions

图5　不同含石量试样破坏形态和颗粒位移Fig．5　Failure mode and particle displacement of specimens with various rock block proportions

4．2破坏特征

图5给出了各含石量土石混合体试样单轴压缩试验的破坏形态，土颗粒为浅灰色球体，块石为不同颜色（灰度）的clumps，图中黑色箭线表示中间纵截面上的颗粒位移矢量。

从图5可见，低含石量（0，10%，30%）土石混合体试样的破坏形态和颗粒位移趋势基本是相同的。试样上部发生明显的侧向膨胀，上部颗粒基本上向两边近水平运动，边缘有局部颗粒脱落，表现出明显的类似于岩石的压致拉裂特征，这与已有土石混合体室内重塑样单轴压缩的破坏形态基本是一致的。中等含石量（50%）土石混合体试样上部和下部颗粒均有明显的侧向膨胀和运动，试样局部边缘出现颗粒脱落，也表现出一定程度的拉裂特征。而高含石量（70%，90%）土石混合体试样没有出现明显的局部突出变形，呈整体式变形，上部颗粒基本沿45°（相对于水平方向）分别向左下和右下方运动，下部颗粒则基本沿45°分别向左上和右上方运动。

4．3细观机理分析

表2是3．3节所述试验过程中设置了初始黏结后的试样在加上较小围压（初始扰动）前后黏结数量变化情况，其中BNSS1为试样扰动前土颗粒间黏结数、BNSR1为试样扰动前土颗粒与块石颗粒间黏结数、BNSS2为试样扰动后土颗粒间黏结数、BNSR2为试样扰动后土颗粒与块石颗粒间黏结数。可以看出各含石量的土石混合体在加上较小围压达到自平衡的过程中，土颗粒间的黏结数变化不大，而土颗粒与块石颗粒间黏结数变化非常大，加上较小围压后即自平衡后，土颗粒与块石颗粒间的接触黏结基本已全部破坏。这是因为土颗粒与块石颗粒之间的黏结强度比土颗粒间小得多，而不规则块石周围又容易出现应力集中，所以即使在很小的扰动下，也基本被破坏，加载时可不考虑。

表2　初始扰动前后试样黏结数Table 2　Number of binding between particles in specimens before and after initial disturbance

图6 含石量为30%和70%试样微裂纹数与轴向应变关系曲线Fig．6　Relationship between micro-crack number and axial strain with rock block proportion of 30%，70%，respectively

图7　不同含石量试样微裂纹分布Fig．7　Distribution of micro-cracks in specimen of different rock block proportions

从表2可以看出，随着含石量的增加，土颗粒间黏结数不断减少，加载前初始试样中总黏结数也不断减少。中低含石量（0，10%，30%，50%）土石混合体试样中的黏结数很多，而高含石量（70%，90%）土石混合体试样中的黏结数较少。在单轴压缩荷载下，由于没有侧向的有效约束，高含石量土石混合体试样中块石之间无法像三轴压缩试验中那样有效地接触、咬合，所以其强度主要取决于颗粒间的黏结。因此，随含石量增加，试样中黏结数越少，整体胶结程度越弱，试样单轴抗压强度和弹性模量不断降低。

从图5可以看出，中低含石量的土石混合体中块石并没有直接接触，而是分散在黏结的土体中，试样整体上呈弱胶结的状态，所以其应力-应变曲线表现出典型的应变软化特征。而高含石量的土石混合体试样中不规则块石相互接触，其间填充的土颗粒很少，孔隙率较大，基本无胶结，所以其应力-应变曲线在峰值后表现出一定程度的塑性流动特征。

将加载过程中监测的微裂纹（颗粒间黏结破裂）总数、剪切破坏形成的微裂纹数和张拉破坏形成的微裂纹数随轴向应变的发展绘成曲线。结果表明，中低含石量（0，10%，30%，50%）的土石混合体微裂纹发展规律基本一致，高含石量（70%，90%）的土石混合体微裂纹发展规律也基本一致。现分别给出具有代表性的含石量为30%和70%的土石混合体试样的微裂纹数与轴向应变的关系曲线，如图6所示。并给出各试样最终的微裂纹分布图，如图7所示，其中拉裂纹用黑色的八边形表示，剪切裂纹用深灰色的八边形表示。

从图6、图7可见，各含石量土石混合体张拉裂纹数、剪切裂纹数和微裂纹总数的变化趋势基本一致：加载初期微裂纹数比较少，增长缓慢，峰值前后发展迅速，裂纹数量急剧增加，后期趋于平稳。中低含石量土石混合体试样中形成微裂纹主要是张拉裂纹，主要分布在试样上部，所以宏观上表现出来的是上部拉裂破坏。而高含石量的土石混合体试样微裂纹很少，主要是剪切裂纹，宏观上表现出来的是剪切破坏。

5结论

为探究不同含石量土石混合体在单轴压缩荷载下的变形破坏机理，本文提出了一种不规则块石三维离散元模型随机生成技术，建立了在一定程度上符合宏观统计规律的土石混合体三维离散元随机模型，对土石混合体单轴压缩试验进行了颗粒流数值模拟，深入分析了其变形破坏的细观机理，在理论上取得了一些有意义的结论：

（1）随含石量增加，土石混合体单轴抗压强度和弹性模量随之减小。

（2）中低含石量（0，10%，30%，50%）的土石混合体中不规则块石没有直接接触，而是分散在黏结的土体中，试样整体上呈弱胶结的状态，其应力-应变曲线表现出典型的应变软化特征。而高含石量（70%，90%）的土石混合体不规则块石相互接触，其间填充颗粒很少，孔隙率较大，基本无胶结，其应力-应变曲线在峰值应力后表现出一定程度的塑性流动特征。

（3）中低含石量（0，10%，30%，50%）土石混合体试样中形成微裂纹主要是张拉裂纹，宏观上表现出来的是拉裂破坏；而高含石量（70%，90%）的土石混合体试样微裂纹很少，主要是剪切裂纹，宏观上表现出来的是剪切破坏。

［1］ 徐文杰，胡瑞林，岳中崎．土-石混合体随机细观结构生成系统的研发及其细观结构力学数值试验研究［J］．岩石力学与工程学报，2009，28（8）：1652-1665．

［2］ 油新华．土石混合体的随机结构模型及其应用研究［D］．北京：北方交通大学，2001．

［3］GRAZIANI A，ROSSINI C，ROTONDA T．Characterization and DEM Modeling of Shear Zones at a Large Dam Foundation［J］．International Journal of Geomechanics，2012，12（6）：648-664．

［4］李晓，廖秋林，赫建明．土石混合体力学特性的原位试验研究［J］．岩石力学与工程学报，2007，26（12）：2377-2384．

［5］COLI N，BERRY P，BOLDINI D．In situ Non-conventional Shear Tests for the Mechanical Characterization of a Bimrock［J］．International Journal of Rock Mechanics and Mining Sciences，2011，48（1）：95-102．

［6］ 董云．土石混合料强度特性的试验研究［J］．岩土力学，2007，28（6）：1269-1274．

［7］刘建锋，徐进，高春玉，等．土石混合料干密度和粒度的强度效应研究［J］．岩石力学与工程学报，2007，26（增1）：3304-3310．

［8］金磊，曾亚武，李欢，等．基于不规则颗粒离散元的土石混合体大三轴数值模拟［J］．岩土工程学报，2015，37（5）：829-838．

［9］贾学明，柴贺军，郑颖人．土石混合料大型直剪试验的颗粒离散元细观力学模拟研究［J］．岩土力学，2010，31（9）：2695-2703．

［10］丁秀丽，李耀旭，王新．基于数字图像的土石混合体力学性质的颗粒流模拟［J］．岩石力学与工程学报，2010，29（3）：477-484．

［11］WANG Yu，LI Xiao，ZHANG Bo，et al．Meso-damage Cracking Characteristics Analysis for Rock and Soil Aggregate with CT Test［J］．Science China Technological Sciences，2014，57（7）：1361-1371．

［12］AFIFIPOUR M，MOAREFVAND P．Mechanical Behavior of Bimrocks Having High Rock Block Proportion［J］．International Journal of Rock Mechanics and Mining Sciences，2014，65（1）：40-48．

［13］廖秋林，李晓，朱万成，等．基于数码图像土石混合体结构建模及其力学结构效应的数值分析［J］．岩石力学与工程学报，2010，29（1）：155-162．

［14］李世海，汪远年．三维离散元土石混合体随机计算模型及单向加载试验数值模拟［J］．岩土工程学报，2004，26（2）：172-177．

［15］CHO G C，DODDS J，SANTAMARINA J C．Particle ShapeEffectsonPackingDensity，Stiffness， and Strength：Natural and Crushed Sands［J］．Journal of Geotechnical and Geoenvironmental Engineering，2006，132（5）：591-602．

［16］宋来忠，彭刚，姜袁．混凝土三维随机参数化骨料模型［J］．水利学报，2012，43（1）：91-98．

［17］NASSAUER B，LIEDKE T，KUNA M．Polyhedral Particles for the Discrete Element Method［J］．Granular Matter，2013，15（1）：85-93．

［18］ZENG Ya-wu，JIN Lei，DU Xin，et al．Refined Modeling and Movement Characteristics Analyses of Irregularly Shaped Particles［J］．International Journal for Numerical and Analytical Methods in Geomechanics，2015，39（4）：388-408．

［19］Itasca Consulting Group，Inc．Universal Distinct Element Code Version 3．1［M］．Minneapolis：Itasca Consulting Group，2006．

［20］MEDLEY E W．The Engineering Characterization of Melanges and Similar Block-in-matrix Rocks（Bimrocks）［D］．Berkeley，California：Department of Civil Engineering，University of California，1994．

［21］CHO N，MARTIN C D，SEGO D C．A Clumped Particle Model for Rock［J］．International Journal of Rock Mechanics&Mining Sciences，2007，44（7）：997-1010．

（编辑：刘运飞）

30%，50%）show a typical feature of strain softening；whereas stress-strain curves of soil-rock mixture with high rock block proportions（70%，90%）show some degree of plastic flow characteristics after peak stress．Micro-cracks in the specimens of SRM with medium-low rock block proportions are mainly tensile cracks，resulting in macroscopic tensile rupture；while few micro-cracks develop in the specimens of soil-rock mixture with high rock block proportions，which are mainly shear cracks，resulting in macroscopic shear failure．

Three-Dimensional Particle Flow Simulation of Uniaxial Compression Tests on Soil-rock Mixture

JIN Lei，ZENG Ya-wu，YE Yang
（School of Civil Engineering，Wuhan University，Wuhan430072，China）

In order to study the deformation failure mechanism of soil-rock mixture under uniaxial compression in the presence of different rock block proportions，on the basis of the refined method for constructing discret element method（DEM）model of irregularly shaped particles，a random technique to constructing DEM model of irregularly shaped rock block is proposed and 3-D particle flow models of soil-rock mixture which conform to the macroscopic statistical rule are established．Then，particle flow simulations of uniaxial compression tests are conducted．It is noteworthy that with the increment of rock block proportions，the uniaxial compressive strength and the elastic modulus both decrease．Stress-strain curves of soil-rock mixture with medium-low rock block proportions（0，10%，

soil-rock mixture；irregularly shaped particles；3-D particle flow；uniaxial compression test；DEM model

TU411

A

1001-5485（2016）07-0093-07

2015-05-20；

2015-06-12

国家自然科学基金项目（41272342）

金磊（1989-），男，安徽安庆人，博士研究生，主要从事岩土体变形破坏机理等方面的研究，（电话）18062654787（电子信箱）whujinlei@whu．edu．cn。

曾亚武（1964-），男，湖北安陆人，教授，博士生导师，主要从事岩石力学与工程方面的教学与研究工作，（电话）13808653023（电子信箱）zengyw@whu．edu．cn。