粗粒土三轴试验颗粒流细观参数敏感性分析

沈筠，莘子健

（1.安徽省交控建设管理有限公司，安徽 合肥 230000;2.中国地质大学（武汉）工程学院，湖北 武汉 430074）

0 引言

数值模拟计算是岩土工程领域重要的研究方法，岩土工程数值计算方法总体上可以分为两大类。一类是基于连续介质力学理论的方法，如有限单元法（FEM）和快速拉格朗日法（FLAC）等。从宏观意义上说，岩体可以视为连续介质。但从本质上讲，岩土材料都是由离散的、尺寸不一、形状各异的颗粒或块体组成的，非连续性是岩土材料的一种重要基本特性。例如，土就是松散颗粒的堆积物，同样，天然岩体也是由被结构面切割而成的大小不一、形态各异的岩石块体所组成。散粒岩土材料的力学特性有着重要的工程应用，如泥砂的沉淀，土堤、土(岩）坡、铁路道砟等的稳定性研究，散粒岩土材料的力学特性研究是岩土力学中最基本的也是最重要的问题之一。

连续介质力学理论把岩土材料作为一个整体来考虑，研究的重点放在建立粒子集合的本构关系，从粒子集合整体的角度研究散粒体介质的力学行为。可以预见，基于连续介质力学理论的方法不能体现颗粒间的复杂相互作用及高度非线性行为，也不能真实刻画散体材料的流动变形特征。

第二类岩土工程数值计算方法，即基于非连续介质力学理论的方法，如离散单元法PFC、UDEC、3DEC等，充分考虑了岩土材料的非连续性，可以有效解决上述问题，因此也越来越受到重视，成国文等通过三维空间随机点生成与颗粒绑定技术，利用颗粒流软件PFC3D，结合工程实际建立了不同含石率的颗粒流三维模型，对模型进行了双轴模拟试验。刘先珊等建立基于平行粘结接触的三维颗粒流数值模型（PFC3D）用以模拟剪切试验中砂岩结构的破坏机理；李识博等采用三维颗粒流软件（PFC3D）对陇西地区相同含水率、不同围压和相同围压、不同含水率条件下的黄土三轴固结不排水剪切试验进行了数值模拟；赵金凤等针对土石混合体的细观材料特性，分别采用球形颗粒单元和非规则组合颗粒单元模拟土体和块石材料，对其在不同含石量和颗粒粘结强度下的直剪试验过程进行离散元分析；毛海涛等采用三维颗粒流软件（PFC3D）并结合控制变量法与莫尔-库伦破坏准则，对不同含水率及围压条件下的非饱和紫色土三轴固结不排水试验进行了数值模拟。

但是颗粒离散元法也有其不足，采用颗粒流方法分析问题时，复杂模型的建立相当繁琐，除此之外细观参数的标定也较为困难且无明显关系可循，一般情况下均采用试错法进行参数标定，有些甚至没有严格的参数标定过程。但是试错法较为盲目，如果能基于一些定性关系，再进行参数标定，将会使得参数标定的效率大大提高。因此，本文基于典型的颗粒离散元软件PFC就室内粗粒土三轴试验建立了离散元模型，并对其细观参数进行了敏感性分析，得到了一些定性关系。

2 三轴试验的颗粒流数值模拟

2.1 PFC基本原理

PFC中的颗粒为刚体，但在力学关系上允许重叠，以模拟颗粒之间的接触力，通过显示差分算法，在颗粒间交替运用牛顿运动定律和力-位移定律进行计算求解。在PFC中，材料的本构特性是通过接触本构模型来模拟的。PFC5.0中提供了10种内嵌的接触模型，常用的接触本构模型有线性接触模型、接触黏结模型和平行粘结模型。当模拟土体时，一般选用线性接触模型，这类模型的特点是抗压、抗剪，但不抗拉；当模拟岩石、混凝土等胶结程度较高、可以承受弯矩等荷载时，一般采用平行粘结模型。本文选择的接触模型为线性接触模型来对粗粒土的三轴试验进行模拟和细观参数的敏感性分析。

2.2 数值模型的生成

三轴模型试样为直径250mm、高度500mm的圆柱样。土颗粒使用三维球体表示，粒径分布范围为6mm～9mm，土颗粒密度为2.7g/cm，孔隙率为0.36，采用的初始细观参数如表1所示，最终生成颗粒8924个，建立的试样模型如图1、图2所示。

初始细观参数 表1

图1三轴试验模型

图2 三轴试样模型

2.3 细观参数敏感性分析

2.3.1 围压的影响

在不同围压下进行模拟得到的应力-应变曲线如图3，由图可知，随着围压的增大，峰值应力也随之明显增大，围压0.4MPa时为0.64MPa，围压0.8MPa时为1.12MPa，围 压1.2MPa时 为1.53MPa。围压为0.4MPa时，应力-应变关系为应变软化型，随着围压的增加，当围压为0.8和1.2MPa时转变为应力硬化型，但是在线性变形段内，应力-应变曲线的斜率较为接近，即弹性模量相差不大，且达到峰值时对应的应变值近似相同。

图3 不同围压下的应力-应变曲线

2.3.2 有效模量的影响

通过设置不同的Emod即有效模量进行模拟，得到了相应的应力应变曲线如图4，由图可知，随着有效模量的增加，应力-应变曲线的弹性模量也随着明显增加，应力-应变关系均为应变硬化型，但达到峰值应力时对应的应变值相差较大，当有效模量为10kPa时，在应变为24%时曲线取得峰值应力，当有效模量为100kPa，在应变18%时取得峰值应力，有效模量为1000kPa时，在应变11%时取得峰值应力，但是不同曲线的峰值应力较为接近，相差不大。

图4 不同有效模量下的应力-应变曲线

2.3.3 摩擦系数的影响

通过设置不同的摩擦系数Fric进行模拟，得到相应的应力-应变曲线如图5，由图可知，随着摩擦系数的增加，应力-应变曲线的峰值应力也随之增加，但是应力-应变关系都为应变硬化型，且曲线的弹性模量近乎相同。

图5 不同摩擦系数下的应力-应变曲线

2.3.4 刚度比的影响

通过设置不同的刚度比（Kratio）进行模拟，得到的应力-应变曲线、泊松比-轴应变曲线如图6、图7，由图6可知，随着刚度比的增大，弹性模量有了一定的增加，当刚度比由0.5增加到1.0时，应力-应变曲线变化不大，刚度比增加到2.0时，峰值强度由1.06MPa增加到1.16MPa。由图7可知，随着刚度比的增加，泊松比也随之增加，在轴应变为4%时，不同刚度比下的泊松比差值最大，在刚度比为0.5、1.0、2.0时分别为0.37、0.44、0.5，随着轴应变的增加，刚度比为0.5和1.0下的泊松比逐渐接近，而刚度比为2.0时的泊松比则保持稳定。

图6 不同刚度比下的应力-应变曲线

图7 不同刚度比下的泊松比-轴应变曲线

2.3.5 颗粒形状的影响

PFC中默认生成的颗粒为圆盘或者球体，并不能完全反映真实颗粒的形状，会对模拟的可靠性造成影响。通过编写相关的Fish代码，将生成的球体颗粒（ball）全部替换为不规则团簇（clump），单个块体和替换后的试样模型如图8所示。将替换后的试样进行三轴试验模拟，并与替换前的应力-应变曲线进行对比，如图9所示，由图可知，将球体颗粒替换为不规则团簇后，试样的峰值应力有了很大的提高，由1.0MPa提高到了1.9MPa，但是在不规则团簇试样在曲线的后半段出现了上下起伏的波动现象，是因为不规则的团簇在剪切过程中出现错动、翻滚，其位置随着剪切过程一直发生着调整。

图8 改进的试样模型

图9 不同颗粒形状下的应力-应变曲线

3 结论

本文基于PFC3D对室内粗粒土的三轴剪切试验进行了数值模拟，并对细观参数进行了敏感性分析，定性探讨了宏细观参数之间的联系，得出了以下结论：

①围压对三轴试样剪切的应力-应变影响较大，随着围压的增大，峰值应力也随之增大，而弹性模量变化不大，且在相同的轴应变处取得峰值应力；

②有效模量对应力-应变曲线的弹性模量影响较大，随着有效模量的增加，弹性模量也随之增加，且峰值应力有了一定的增加，但是取得峰值应力时对应的轴向应变相差较大，有效模量越大，越早达到峰值应力；

③摩擦系数对应力-应变曲线的峰值应力影响较大，随着摩擦系数的增加，试样的峰值应力随之增加，而弹性模量变化不大；

④刚度比对应力-应变曲线影响不大，随着刚度比增加，弹性模量和峰值强度略微增加，而对泊松比-轴应变曲线影响较大，在轴应变较小时，泊松比差距较大，随着轴向应变的增加，泊松比差值越来越小；

⑤土颗粒的形状对应力-应变曲线影响较大，当土颗粒形状由球体替换为不规则块体时，应力-应变曲线的峰值应力剧烈增加，弹性模量也随之增大，但是由于替换了全部的球体颗粒，使得应力-应变曲线的后半段出现了波动，是因为不规则块体随着剪切进行发生的翻滚、跳动。