非饱和重塑与结构性黄土平面应变试验三维离散元模拟

蒋明镜，王优群，卢国文，张 鹏

(1.天津大学 建筑工程学院土木系， 天津 300072；2.天津大学 水利工程仿真与安全国家重点试验室， 天津 300072；3.同济大学 土木工程防灾国家重点试验室， 上海 200092；4.同济大学 土木工程学院地下建筑与工程系， 上海 200092)

天然黄土是一种主要以粉粒尺度构成的特殊的结构性土[1]，因其沉积特点与所处气候环境黄土多以非饱和状态存在。同时，黄土特有的大孔隙架空结构和颗粒间的胶结作用使其宏、微观力学性质与重塑土有较大差异，当该结构遭到扰动时将呈现出复杂的应力状态，进而诱发各种工程地质灾害。

在黄土地区的实际工程中，如深基坑坑壁土体[2]、高填方路堤堤身[3]、隧道围岩体[4]等近似于平面应变问题的工程屡见不鲜，开展平面应变条件下非饱和结构性黄土的宏微观力学特性研究具有重要的理论意义与应用价值。然而，目前许多平面应变问题的工程参数设计仍采用常规三轴试验条件下测试的力学指标，并未考虑中主应力变化引起应力路径改变对土体强度的影响，应用于平面应变问题分析具有一定的局限性。基于此，张玉等[5-7]通过平面应变状态的真三轴室内试验，揭示平面应变条件下土的力学性质指标和变形、强度准则。但由于天然结构性黄土原位取样难度大，可重复性较差、真三轴试验装设计复杂且成本高，从而制约了结构性黄土宏微观力学性质的研究和发展。随着微观测试技术的发展，诸多学者通过电镜扫描[8]、CT扫描[9]、压汞测试[10]等观测手段分析了黄土微观结构特征，对宏观力学特性受控于其微观特性影响有了更深入地认识，但较难观测试样微观结构及其演化过程。

Cundall等[11]于1979年提出了适合于颗粒体介质的颗粒离散元法(Distinct Element Method，DEM)，建议从土体的微观形成及胶结机理出发研究土体宏观力学性质及变形。该方法具有成本低，可重复性好，能够实时监测试样颗粒、接触、粒间接触、胶结破坏等微观信息演化等特点。近年来，离散单元法逐步应用于非饱和土的宏微观特性的分析。熊超[12]考虑了粒间毛细力影响，通过修改PFC软件内置接触刚度模型将毛细力作用植入到离散元模型中，对非饱和土进行三轴压缩模拟，探究了不同微观参数对宏观力学行为的影响。蒋明镜等[13]更进一步地在含抗转动和抗扭转完整接触模型的基础上引入范德华力和毛细力以及不可恢复的化学胶结作用，建立了可以全面考虑含水率-孔隙比-吸力耦合作用的非饱和结构性黄土三维微观接触模型。该模型考虑全面，可以较好地反映非饱和黄土的宏微观力学状态，已成功应用于黄土的侧限压缩[14]、等向压缩[15]试验。

基于此，本文引入该非饱和结构性黄土三维接触模型，开展了非饱和重塑与结构性黄土的三维DEM分析，模拟常含水率平面应变试验，分析了压缩曲线与微观特征。可为建立考虑微观因素的平面应变条件下非饱和结构性黄土本构理论提供参考。

1 试验部分

1.1 模型介绍及参数确定

合理的离散元微观接触模型是准确模拟非饱和结构性黄土的关键。采用团队近期提出的非饱和结构性黄土三维接触模型，该模型在颗粒间胶结接触力学特性室内接触试验和数值模拟基础上，提出了胶结计算刚度和强度公式，考虑了颗粒吸引力和不可恢复的化学胶结作用。该模型以三维含抗转动、抗扭转完整接触模型为基础，在不考虑化学胶结作用时可用于非饱和重塑黄土离散元分析，考虑化学胶结作用时用于非饱和结构性黄土离散元分析。表1为DEM模拟拟采用的材料参数，关于模型的具体介绍以及参数标定方法详见文献[13,16]。

表1 胶结接触模型微观参数

1.2 离散元试样制备

基于图1所示粒径级配曲线以及表1所示接触模型参数进行制样，颗粒中值粒径d50=2×10-5m；颗粒密度ρ=2 710 kg/m3；局部阻尼系数为0.7；含水率为7.1%。制样步骤如下：

(1) 采用分层欠压法[17]分6层生成含有40 842个颗粒的松散无胶结长方体试样用以模拟重塑黄土，试样宽高比为2∶1∶1，与室内试验的尺寸相同，目标孔隙比e0=1.2。考虑到黄土粒间范德华力对黄土沉积过程的影响，制样中引入范德华力的作用，用以制成松散均匀的大孔隙试样，取范德华力σvan=4 kPa。

(2) 在侧限条件下对该无胶结试样进行预压，竖向应力为12.5 kPa，用以模拟地表下1 m左右的土体单元应力状态。

(3) 在无胶结试样内部符合胶结条件区域生成胶结，用以模拟结构性黄土。

图1 离散元模拟颗粒级配曲线

1.3 实验步骤

DEM模拟步骤与室内试验步骤保持一致，步骤如下：

(1) 固结阶段。对制备完成的试样进行等向加载至不同围压，该过程由伺服系统的六面刚性墙完成。

(2) 剪切阶段。固定中主应力方向的刚性墙，保持中主应力方向应变为零，用以模拟平面应变条件。保持伺服围压与含水率不变，使大主应力方向刚性墙相向运动，应变率为5%/min，实时记录各项指标变化。

2 结果分析

2.1 应力应变响应

图2给出了非饱和黄土离散元试样平面应变试验的应力应变关系规律，本文的目的旨在从宏微观尺度探究其力学特性的内在机理，并不追求与特定黄土室内试验结果完全吻合。图中可以看出，重塑与结构性土试样均表现出压硬性，即峰值偏应力随围压(σ3)的增大而增大；非饱和重塑土的应力应变关系在低围压时呈现弱软化或剪切滑移的理想塑性表现，随围压升高软化现象逐渐明显，这可能与高围压条件下随着竖向应力的升高导致颗粒破碎有关。而结构性黄土试样应力应变关系均为应变软化，且围压较小时，软化较明显。低围压时，土的抗剪强度由压应力作用下的颗粒摩擦和土组构与胶结的抗剪应力提供，而高围压时，土的抗剪强度由压应力作用下的颗粒摩擦和土组构与胶结破损的抗剪应力提供，说明高围压导致黄土原生结构损伤，土体颗粒重新排列，次生结构在固结阶段已经开始形成，导致在高围压下结构性黄土的表现更趋近于重塑黄土；而低围压情况下土体结构性未遭到破坏，在峰值应力过后因胶结持续破坏会导致较重塑土更大的应变软化现象。变软化现象在常围压下出现与剪胀有关，对于结构性土，除剪胀外还与胶结力等结构性丧失有关[18]。与此同时，受结构性影响，各级围压下结构性黄土峰值偏应力及弹性模量大于重塑试样。

图2 非饱和黄土平面应变试验应力应变关系离散元模拟结果

2.2 胶结破坏分析

颗粒间胶结作用为非饱和结构性土区别于重塑的主要微观特征，粒间胶结的存在以及在受荷增湿过程中的破坏造成了结构性土特有的宏观力学表现。图3给出了非饱和结构性黄土离散元试样平面应变试验胶结破坏数与轴向应变的关系。图中可以看出，胶结破坏在剪切初期几乎未发生破坏，在剪切过程中在阈值轴向应变2%左右开始出现胶结破坏并迅速增多，随着破坏数量接近完全破坏时，破坏速率才减慢最终接近平稳。对比应力应变关系可以发现，胶结开始破坏均发生在峰值应力出现之前，可以推断出胶结力对于结构性土的峰值应力有很大的影响。且对比不同围压发现：围压越大，胶结破坏数越多，增长速度越快。

图3 非饱和结构性黄土离散元平面应变试验胶结破坏数与轴向应变的关系

2.3 能量转化

单元试验通过应力控制或应变控制对试样做功，在实验过程中，墙体做的功转化为颗粒弹性能、胶结弹性能、塑性耗散能、阻尼能等。图4、图5给出了非饱和重塑黄土和结构性黄土离散元试样平面应变试验墙体做功和能量变化。图中可以发现，重塑土中墙体做功绝大部分转化为塑性耗散能，塑性耗能表现为试样内颗粒之间发生的滑动、转动和扭动，宏观上表现为塑性变形。颗粒弹性能在剪切过程中不断增大，而达到峰值应力附近保持稳定，说明在峰值应力后土颗粒本身将不再承担额外剪切荷载，荷载将由压应力作用下的颗粒摩擦的抗剪应力提供。随着围压的增加颗粒排列更为紧密，塑性耗能、阻尼耗能及颗粒弹性能均增大，且颗粒性能峰值随围压的增加而增大。

图4 不同围压下非饱和重塑黄土能量变化与轴向应变的关系

图5 不同围压下非饱和结构性黄土能量变化与轴向应变的关系

结构性土中墙体做功绝大部分转化为塑性耗散能结论不变。结构性土中，颗粒弹性能在剪切初期迅速增大到达峰值后有一定的下降；考虑到胶结作用，胶结弹性能受围压的影响较小，随着剪切的进行与胶结破坏数的增加，胶结弹性承担的耗能会逐渐降低。当胶结开始破坏时，胶结破坏耗能迅速增加，而后逐渐趋于平缓，由此可知胶结破坏在受剪过程中为持续性破坏，围压增大，胶结破坏耗能也随着增加。

2.4 力学配位数

力学配位数是试样中扣除接触数量小于1的平均每个颗粒拥有的接触数量，与土体的力学性质关系密切，具体计算方法参见文献[19]。图6给出了非饱和黄土离散元试样平面应变试验力学配位数与轴向应变的关系。图中可以看出，重塑黄土的力学配位数在剪切初期可能由于体缩有一定的上升，随后随着球应力的增加，基本达到稳定值后保持不变。而结构性土力学配位数在剪切初期有一定上升，随着剪切的进行，力学配位数有一定的下降，而后达到稳定。且围压越大，力学配位数也增多。结构性黄土的力学配位数初始状态下比重塑土力学配位数大。在小围压条件下，结构性黄土在胶结大量破坏后其力学配位数会快速靠近重塑黄土。

图6 非饱和黄土平面应变试验离散元模拟力学配位数与轴向应变的关系

2.5 接触偏组构

随机排土体的接触偏组构是反映土体颗粒间接触的各向异性程度的重要微观变量，接触偏组构越大，粒间接触方向的各向异性越强，当材料为各向同性时，偏组构为零[20]。图7给出了非饱和黄土离散元试样平面应变试验接触偏组构与轴向应变的关系。从图中可以看出，围压对于接触偏组构影响较小。重塑样接触偏组构随轴向应变的增加出现先快速增加，达到峰值后缓慢下降的趋势。而结构性土在剪切初期接触偏组构有轻微下降，当胶结开始大量破坏时，接触偏组构开始迅速增加达到峰值后开始趋于平稳。对比可以发现，结构性土与重塑土接触偏组构在达到平稳段后趋于一致，由此说明胶结会阻碍各向异性的发展。

图7 非饱和黄土平面应变试验离散元模拟接触偏组构与轴向应变的关系

3 结 论

(1) 重塑与结构性土试样均表现出压硬性，体变均呈轻微剪胀；非饱和重塑土的应力应变关系围压较小时呈应变硬化，围压较大时，出现应变软化；而结构性黄土试样应力应变关系均为应变软化，且围压较小时，软化较明显，结构性表现充分。

(2) 胶结在剪切初期基本不发生破坏，当轴向应变到达阈值之后开始大量破坏，对比应力应变关系可知，在不同围压下，胶结开始破坏均发生在峰值应力出现之前，可以推断出胶结力对于结构性土的峰值应力有很大的影响。

(3) 重塑土与结构性黄土试样中墙体做功绝大部分转化为塑性耗散能。结构性土中胶结弹性能受围压的影响较小，随着剪切的进行与胶结破坏数的增加，胶结弹性能会逐渐降低。胶结破坏在受剪过程中为持续性破坏，胶结破坏耗能会随着围压增大而增加。

(4) 重塑与结构性黄土的力学配位数在剪切初期有一定的上升，且围压越大，力学配位数也增多。围压对于接触偏组构影响较小。