离散元颗粒集合的应力路径相关性研究

2019-09-17灿汪儒鸿彭国园宋强辉6

长江科学院院报 2019年9期

李灿汪儒鸿彭国园宋强辉6

(1.陆军勤务学院军事设施系，重庆 401331； 2.陆军勤务学院岩土力学与地质环境保护重庆市重点实验室,重庆 401331； 3.重庆科技学院，重庆 401331； 4.中国人民解放军95133部队，武汉 430000； 5.西北核技术研究所，西安 710024； 6.重庆地质矿产研究院，重庆 400042)

1 研究背景

天然岩土材料的室内试验存在材料的不均匀性、现场取样对试样的扰动以及试验重复性低等问题。最关键的是，室内试验无法获得加卸载条件下岩土材料内部的微观信息[1]。因此，近年来借助数值模拟方法来模拟岩土材料的研究越来越多。比如连续数值分析方法中的有限单元法、边单元法、有限差分法等；非连续介质分析方法中的块体离散元法、颗粒离散元法等。其中，颗粒离散元法由于不受变形量限制，可方便地处理非连续介质力学问题，体现多相介质的不同物理关系，在模拟岩土材料问题上得到了广泛的应用[2]。比如Wang等[3]对水泥胶结砂土进行了室内三轴压缩排水试验和双轴压缩试验的离散元模拟，分析了人工胶结砂土特性的内在机理以及胶结对剪胀和强度的影响；周健等[4]运用PFC(颗粒流)程序对砂性土坡及黏性土坡分别进行了数值模拟，分析细观参数对土坡破坏型式的影响；蒋明镜等[5-7]利用 PFC2D及NS2D等离散元软件搭建了多种人工胶结砂土的二维数值模型，对胶结砂土的剪切特性和应变局部化进行了分析；崔铁军等[8]利用PFC3D模拟研究了露天矿边坡内不同高度、埋深、装药量的单孔爆破过程。以上研究均较好地利用颗粒离散元方法模拟了部分岩土材料的主要力学性能。

然而，目前应用颗粒离散元方法来展开岩土材料的研究均只考虑了单一应力路径，对数值试样力学特性的分析还停留在模拟常规三轴压缩试验的水平，不能很好地反映不同应力路径对岩土材料在变形、破坏等方面的影响。王伟等[9]、刘恩龙等[10]、汪斌等[11]的室内试验结果表明：不同应力路径对岩土材料力学特性的影响非常大。除此以外，涉及到边坡开挖、基坑施工等的工程都存在复杂的应力环境，基于单一应力路径开发的岩土材料数值模型不一定能合理反映其真实的力学响应，限制了数值模型的实际应用价值。针对以上问题，通过PFC2D自定义Fish语言开发出了能模拟不同应力路径的双轴压缩程序，对离散元颗粒集合进行了大量模拟加卸荷应力路径的双轴压缩试验，得到了相应结论。

2 PFC方法模拟岩土材料简述

PFC为颗粒离散元程序，其基本原理来自分子动力学，是从微观角度研究介质力学特性和行为的工具，它的基本构成为圆盘和圆球颗粒，然后利用边界墙(Wall)约束。计算时不需要给材料参数定义宏观本构关系及对应的参数，只是采用局部接触来反映宏观问题，因此只需要定义颗粒和黏结的几何和力学参数。在模拟岩土材料时，对生成的离散元颗粒集合进行大量数值试验(单轴、双轴、三轴试验等)，根据数值试验的结果不断调整颗粒的微观参数，直至其力学特性与真实材料的力学特性相匹配。模拟岩土材料与真实材料相匹配的过程相对复杂，且已有大量研究，故不作为本文阐述重点，这里主要讨论应用离散元方法生成的颗粒集合是否也同样具有应力路径相关性，从而验证以其模拟岩土材料的可靠性。

图1 颗粒集合双轴压缩模型Fig.1 Model of dual-axial compression on particle sets

3 不同应力路径的双轴压缩试验

下面将详细介绍考虑不同应力路径的离散元颗粒集合平面应变双轴压缩试验。

3.1 二维试验模型搭建

离散元颗粒集合采

用蒋明镜等[12]提出的分层欠压法制得。首先，将颗粒间摩擦系数设置为1.0，颗粒法向及切向接触刚度均设为1.5×108N/m，即将颗粒本身视为刚体。采用分层欠压法制备出初始平面孔隙比e0=0.28、宽600 mm、高1 200 mm的试样，如图1所示。

在图1模型中，离散元颗粒若采用真实岩土颗粒粒径将导致其最终数量远超出程序计算上限。这里参考文献[13]—文献[15]模拟岩土材料时普遍采取的颗粒粒径范围，将试样颗粒最大粒径设为10.0 mm，最小粒径为7.5 mm，不均匀系数为1.3，颗粒密度为2 000 kg/m3。最后，为使该颗粒集合的力学特性贴近真实岩土材料，还需要在颗粒间所有的接触点位置形成胶结，其力学特性由接触黏结模型参数控制。由于本文主要探讨离散元颗粒集合的应力路径相关性，对颗粒微观参数与宏观力学参数的对应关系不加以分析，仅展示颗粒主要微观参数的选取，如表1所示。

表1 接触黏结模型颗粒主要微观参数Table 1 Major microscopic parameters of bonding particles

3.2 双轴压缩试验

3.2.1 围压控制

通过PFC2D自定义Fish语言可实现2种不同应力路径下的双轴压缩试验模拟。首先，使材料在初始围压中等向固结，即将围压与轴压以同等速率加载至相同的数值；然后，通过伺服控制使轴向加载速率恒定为应变5%/min，该加载速率能基本保证试样处于准静态变形条件，且试样内部应力场在应变局部变化之前处于均匀状态。在这一过程中，通过控制围压大小的变化，使材料沿不同的应力路径被轴向加载直至破坏。根据应力路径的不同，可将加载过程中围压的变化情况分为2类：第1类，轴向加载过程中保持围压不变；第2类，轴向加载过程中围压逐渐变小。

除此以外，对预加载模型设置不同的初始围压进行多组试验，初始围压分别为1.5，1.0，0.5 MPa。为便于读者理解，以初始围压1.5 MPa为例，展示双轴压缩过程中2类不同应力路径下的围压及轴压变化情况，如图2所示。

图2 2类不同应力路径下的围压及轴压变化情况Fig.2 Changes in confining pressure and axial pressure under two different stress paths

3.2.2 结果监测

在双轴压缩试验过程中，竖直向应力即轴压定义为σ1，水平向应力即围压定义为σ3。平均应力p=σ1+σ3/2，偏应力q=σ1-σ3。通过监测以上数据，生成试样的偏应力-轴向应变曲线。

4 试验结果的比较

将初始围压大小不同时的双轴压缩测试结果进行数值模拟之间的比较，并将数值模拟结果与刘恩龙等[10]的室内实验结果进行横向比较。

4.1 数值模拟之间的比较

在初始围压大小分别为1.5，1.0，0.5 MPa的条件下，离散元颗粒集合在不同应力路径下的偏应力-轴向应变曲线如图3所示。

图3 2类不同应力路径下偏应力-轴向应变曲线Fig.3 Curves of deviatoric stress versus axial strain under two different stress paths

从图3可以看出：不同初始围压条件下，围压变小的应力路径使离散元颗粒集合能达到的峰值强度均远小于围压不变的应力路径；同时，2类应力路径都表现出围压越大，材料的峰值强度更高，出现应变软化现象对应的轴向应变更小。上述结果说明应力路径对离散元颗粒集合的影响主要体现在材料破坏强度上面，且具有明显的规律性。

4.2 与真实岩土材料试验的横向对比

为更好地将离散元颗粒集合所具有的应力路径相关性与真实岩土材料进行比对，这里将数值模拟的结果与刘恩龙等[10]所做的针对结构性土应力路径相关性的研究进行比对。这是因为结构性土作为一种岩土材料本身就具有大孔隙、粒间胶结的特点，与离散元颗粒集合较为相似。需要说明的是：在离散元颗粒集合中不存在孔隙水，对其进行双轴压缩试验时不需要考虑孔隙水压力的问题。就试验条件而言，相当于常规三轴压缩模拟的固结排水试验。因此，将数值模拟的结果与刘恩龙等[10]室内试验中的固结排水试验结果进行比对。图4为结构性土在固结排水试验中的总应力路径示意图。

图4 结构性土总应力路径示意图[10]Fig.4 Total stress path of structural soil[10]

图4中OA段表示等向固结阶段；AC段表示围压不变而增大竖向应力直至试样破坏的应力路径；AB段表示减小围压的同时保持竖向应力不变直至试样破坏的应力路径。从破坏包络线可以看出：围压不变的应力路径比围压变小的应力路径使得结构性土达到的破坏强度更高。相比之下，离散元颗粒集合在不同应力路径影响下也具有相似的破坏轨迹，进一步验证了离散元颗粒集合也具有与真实岩土材料相似的应力路径相关性。

4.3 颗粒微观现象分析

为进一步体现应力路径对离散元颗粒集合破坏特征的影响，结合颗粒微观接触力链进行分析。在试样内部所有相互接触的颗粒间画一条直线，并用该直线的粗细程度表示两颗粒之间接触力的大小(线段越粗表示颗粒间接触力越大)，由这些粗细不同的线段组成的网格即是接触力链。因此可通过颗粒接触力链来分析材料内部的应力分布情况。图5为离散元颗粒集合在2类应力路径影响下双轴加载过程中的接触力链示意图。其中围压变小的应力路径下颗粒集合内部应力分布不均匀，上下两端受力多，左右两侧受力少；相比之下，围压不变的应力路径下颗粒集合的应力分布更加均匀，材料能更充分地发挥强度特性，因此峰值强度更高。