循环荷载下圆砾的宏观动力特性及其细观模拟

Macroscopic Dynamic Characteristics and Microscopic Simulation of Round Gravel Under Cyclic Loading

摘要

为探讨地铁循环荷载作用下地基圆砾的宏观动力特性演化规律及其内部细观结构的动力响应特征，利用DYNTTS大型三轴循环试验系统，开展一系列饱和不排水动三轴试验；同时，基于室内试验建立离散元（DEM）不排水动三轴模型，分析不同相对密度、围压及动应力幅值下圆砾的宏观动力特性演化特征及内部细观参量响应规律。结果表明：当围压较小（σ₃=100 kPa）且相对密度较低（D_r=0.3）时，圆砾的累积塑性应变过大（高于0.61%），列车运行存在严重的安全隐患；当围压较高（σ₃=200、300 kPa）时，累积塑性应变较小（低于0.14%），列车运行相对安全。相对密度和围压越大，试样的累积塑性应变越小，回弹模量越大，试样内部耗损的能量越低，故可考虑采用高压注浆等挤压措施提高土体相对密实度以减小地铁地基的沉降。循环加载过程中，配位数及各向异性均呈周期性变化，配位数的变化方向与循环荷载的变化方向相反，而各向异性的变化方向与循环荷载的变化方向相同。动应力幅值越大，配位数越小，各向异性越强；相对密度及围压越大，配位数越大，各向异性则越弱。强接触主要由轴向应力控制，其对承载力起着主要支撑作用。研究成果可为圆砾地基地铁线路的设计及安全运营提供参考。

Abstract

A series of saturated undrained dynamic triaxial tests are conducted using the DYNTTS large-scale triaxial cyclic test system, and a discrete element (DEM) undrained dynamic triaxial model is established based on indoor experiments to explore the evolution of macroscopic dynamic properties and the internal microstructural dynamic response of round gravel under subway cyclic loading. The macroscopic dynamic characteristics of the round gravel and the internal microscopic parameter response law under different relative densities, confining pressures, and dynamic stress amplitudes are analyzed. The results show that when the confining pressure is small (σ₃=100 kPa) and the relative density is low (D_r=0.3), the accumulated plastic strain of the round gravel is too large (greater than 0.61%), posing a severe safety hazard for train operation. When the confining pressure is higher (σ₃=200, 300 kPa), the cumulative plastic strain is smaller (less than 0.14%), indicating relatively safe train operation. The greater the relative density and confining pressure, the lower the cumulative plastic strain of the specimen, the higher the resilient modulus, and the lower the energy consumption inside the specimen. Therefore, extrusion measures such as high-pressure grouting are considered to improve the relative compactness of the soil to reduce the settlement of the subway foundation. During the cyclic loading process, both the coordination number and anisotropy change periodically. The direction of change in coordination number is opposite to that of cyclic loading, while the direction of change in anisotropy aligns with cyclic loading. With the increase in relative density and confining pressure, the coordination number becomes larger and the anisotropy becomes weaker. With the increase in dynamic stress amplitude, the coordination number decreases and the anisotropy increases. The strong contact is mainly controlled by the axial stresses, which play a significant role in bearing capacity. The research results provide a reference for the design and safe operation of subway lines on round gravel foundations.

Graphical abstract

关键词

循环荷载 / 圆砾 / 动力特性 / 离散元

Key words

cyclic loading / round gravel / dynamic characteristics / discrete element

本刊网刊

圆砾是自然界中较为常见的粗粒土^[1‒2]，其颗粒形状多为球形和亚球形，表面光滑，颗粒级配差，具有强渗透性。随着南宁市多条地铁线路的修建，出现了大量穿越圆砾地层的地下工程。在地铁循环荷载的长期作用之下，圆砾地基极易发生过量的累积沉降，从而导致路基不均匀下沉或因沉降超标而严重威胁行车安全。因此，探究圆砾在地铁循环荷载作用下的宏观动力特性及其内部细观结构的演化特征十分有意义。

长期以来，为了解粗粒土在循环荷载作用下的动力特性，有学者利用大型动三轴仪对粗粒土开展了大量的研究，并取得了丰硕成果^[3‒4]。龙尧等^[5]为探讨碎石路基在循环列车荷载作用下的动力特性，开展了一系列的室内动三轴试验，分析了不同循环应力比及围压对碎石路基累积变形的影响，发现循环应力比越大，其软化指数越高；冷伍明等^[6]利用大型动三轴仪对重载铁路路基粗粒土进行了一系列试验，探讨了粗粒土在不同围压、动应力幅值及含水率条件下的累积变形特性，发现增大围压或减小动应力有助于增强土体的稳定性；聂如松等^[7]利用大型动三轴仪研究了粗粒土填料在重载列车荷载作用下的回弹特性，结果发现回弹模量随着围压及动应力幅值的增大而增大，随着含水率的增大而减小；Duong等^[8]通过对3种含水量和4种细粒含量下制备的道砟试样进行循环三轴试验发现，试样在饱和状态下，细粒含量越高，永久轴向变形越大，而在细粒低含量情况下，该趋势相反；Leng等^[9]利用定制的大型三轴设备，对重载铁路路基粗粒土填料进行一系列动三轴试验，分析了含水率、围压和动偏应力对试样弹性模量和累积轴向应变的影响，并提出弹性模量和累积轴向应变的计算公式；Yang等^[10]对重载铁路路基上不同污染程度的碎石进行了一系列动态三轴试验，结果表明，随着污染细颗粒含量的增加，试样的累积塑性应变近似线性增加。此外，在循环荷载作用下，分析材料内部细观结构的演化特征对于深入认识试样的宏观动力特性具有重要意义。近年来，随着计算机的发展，诸多研究人员利用离散元对材料内部细观结构进行了研究。孙逸飞等^[11]为研究循环荷载作用下散粒体级配对其变形特性的影响，建立相应的2维循环离散元模型进行模拟，结果表明散粒体的变形随不均匀系数的增加而增加；张伏光等^[12]采用3维离散元程序，对胶结砂土进行3维循环不排水离散元模拟，证实了离散元模拟能够较好地反映胶结砂土的宏观动力学性质；Nguyen等^[13]利用离散元方法对粗粒土进行一系列不同动应力振幅和加载频率的循环三轴试验，结果表明加载频率对应变累积影响很小，随着应变的累积，配位数略有增加；Nie等^[14]采用离散元方法对循环荷载下砾石‒砂混合料进行模拟，从细观角度描述了循环荷载下接触力的各向异性特性。综上所述，众多学者以含水率、细料含量、围压、固结比和动应力幅值等为主要影响因素，对粗粒土填料的宏、细观动力特性进行研究，取得了丰硕的成果。但是，上述研究的粗粒土大都为重载铁路和高速铁路下的道砟及碎石填料，与之相比，圆砾颗粒形状多为球形和亚球形，表面光滑且级配较差。同时，相较于重载和高速铁路荷载，地铁荷载具有振幅低、运次高及持时长等特性。因此，有必要针对地铁循环荷载作用下铁路地基圆砾的宏观动力特性演化规律及内部细观结构的动力响应开展进一步的研究。

基于此，以南宁地铁地基圆砾作为研究对象，利用DYNTTS大型三轴循环试验系统，进行一系列饱和不排水动三轴试验，探索了在地铁循环荷载作用下，相对密度、围压及动应力幅值对圆砾动力特性的演化规律。同时，基于室内试验，建立了离散元（DEM）不排水动三轴模型，分析了循环荷载下试样内部细观参量的变化，从细观角度探讨了不同相对密度、围压及动应力幅值下圆砾的动力响应机制。

1 试 验

1.1 试验仪器及用土

试验采用粗粒土DYNTTS大型三轴循环试验系统，如图1所示。该系统主要包括轴向加载系统、压力罩、围压/反压控制器、孔隙水压力传感器和数据采集系统等。该系统允许的施加最大围压为1.0 MPa，最大轴向应力为60 kN，最大加载频率为2 Hz，最大轴向行程为75 mm。

试验所用圆砾取自南宁地铁某施工现场，土样及其颗粒级配曲线如图2所示。

1.2 试验方案

根据南宁地铁施工现场的勘察报告及相关文献资料，南宁地区圆砾层大都呈稍密至中密状态，局部密实，故对相对密度D_r=0.3、0.5和0.7的圆砾土进行试验。设置3种围压分别为100、200和300 kPa。南宁地铁隧道大都位于常水位以下，故试样设为饱和状态。地铁荷载下地基土所受动应力大致为10～40 kPa^[15]，故试验设定4级动应力幅值σ{}_{\mathrm{d}}^{\mathrm{a}\mathrm{m}\mathrm{p}}=10、20、30和40 kPa，其中σ{}_{\mathrm{d}}^{\mathrm{a}\mathrm{m}\mathrm{p}}=σ{}_{\mathrm{d}}^{\mathrm{m}\mathrm{a}\mathrm{x}}‒σ{}_{\mathrm{d}}^{\mathrm{m}\mathrm{i}\mathrm{n}}，σ{}_{\mathrm{d}}^{\mathrm{m}\mathrm{a}\mathrm{x}}和σ{}_{\mathrm{d}}^{\mathrm{m}\mathrm{i}\mathrm{n}}分别为循环荷载的最大值和最小值。具体试验方案见表1。

1.3 试验加载条件设置

考虑地铁荷载作用具有瞬时性且地铁隧道埋藏较深，故试验采用固结不排水条件。采用正弦波模拟实际运营条件下列车动荷载对地基的作用^[16]，如图3所示。为降低试验工作量并减少由不同试样间差异引起的结果偏差，采用多级加载法对试样进行加载^[17‒18]，每级动应力幅值下加载次数设定为10 000次。加载频率按照公式f=V/L计算^[19]，其中，南宁地铁速度取V=72 km/h，地铁每节车厢长L=20 m，计算得到加载频率为1 Hz。

1.4 试验过程

砾类土的原状土样较难取得及保存，故均采用重塑土样。按照《铁路工程土工试验规程》（TB 10102—2010）^[20]的要求，采用分层振捣法制备试样，共分10 层。即：首先，根据试验方案称取每层圆砾质量并搅拌均匀，而后分层击实至相应高度，最终制成直径为300 mm、高度为600 mm的试样。试样制备完毕后，从仪器底部连续3 h充入二氧化碳，然后连续12 h充入无气水以溶解试样内部的二氧化碳，最后采用分级加压的方式对试样进行反压饱和。反压共分4级，每级增加50 kPa并维持4 h，最后一级反压增至200 kPa。饱和结束后，若孔压系数B大于0.95，说明试样基本饱和，可对试样进行加载。试验过程如图4所示。

2 试验结果与分析

2.1 累积塑性应变

累积塑性应变是决定地基劣化程度的关键性指标，其大小反映了地基的沉降水平。因此，本文绘制了在地铁循环荷载作用下，不同相对密度、围压及动应力幅值条件下圆砾的累积塑性应变ε_p与加载次数N之间的关系曲线，如图5所示。由图5可知，不同围压条件下，不同相对密度试样的ε_p‒N曲线随动应力幅值的不同，表现出不同的发展形态。根据Werkmeister等^[21]的安定理论（式（1）），将地铁循环荷载作用下圆砾的ε_p‒N曲线分成3种类型：1）塑性安定型（Rang A）。循环荷载下累积塑性应变在加载初始阶段快速增长，随后逐渐变缓直至相对稳定，且累积塑性应变较小，此时试样处于弹性状态。2）塑性蠕变型（Rang B）。循环荷载下累积塑性应变在加载初始阶段快速增长，随着加载继续，应变增长速率变缓但一直处于持续增长的状态，此时试样处于塑性蠕变状态。3）增量破坏型（Rang C）。循环荷载下累积塑性应变急剧增长，在很小的振次范围内，试样已达到破坏状态。

（1）

式中，ε{}_{\mathrm{p}}^{\mathrm{5}\mathrm{000}}和ε{}_{\mathrm{p}}^{\mathrm{3}\mathrm{000}}分别为试样在加载次数为5 000和 3 000时所对应的累积塑性应变。

参考上述分类标准，对试样的累积塑性应变形态进行划分，如图5所示。由图5（a）可知，在低围压σ₃=100 kPa下，D_r=0.3的试样在σ{}_{\mathrm{d}}^{\mathrm{a}\mathrm{m}\mathrm{p}}=40 kPa的作用下，累积塑性应变较大（约0.64%），并有急剧增加的趋势，此时试样处于增量破坏状态，这表明对于埋深较浅（σ₃=100 kPa）且相对密度较小（D_r=0.3）的地铁隧道，地基极有可能产生不均匀沉降，严重威胁列车的行车安全，故应对圆砾地基的强度进行提升。

由图5（b）和（c）可知，在高围压（σ₃=200、300 kPa）下，试样的累积塑性应变曲线发展形态均表现为塑性安定型，这主要是由于高围压限制了试样内部颗粒间的滑移，提高了试样承受外荷载的能力。这表明对于埋深较深的地铁隧道，列车的运行相对安全，不需要对圆砾土地基进行处理。

为更加直观地了解相对密度、围压及动应力幅值对圆砾累积塑性应变的影响，图6给出了试样在4个加载阶段结束后的最终累积塑性应变ε_p与相对密度D_r之间的关系曲线。

由图6可知：相对密度和围压越大，累积塑性应变越小；动应力幅值越大，累积塑性应变越大。以σ₃=100 kPa、σ{}_{\mathrm{d}}^{\mathrm{a}\mathrm{m}\mathrm{p}}=40 kPa为例，当相对密度D_r为0.3、0.5、0.7时，累积塑性应变分别为0.61%、0.16%和0.14%，相较于D_r=0.3，D_r=0.5和0.7的塑性应变分别降低了73.77%和77.05%，可见相对密度的增大能够显著降低试样的累积塑性应变。这表明通过采用高压注浆等挤压措施提高圆砾的相对密实度可有效地降低列车运行期间隧道的沉降。

2.2 滞回曲线

为探索不同动态变形行为试样滞回曲线的演变和特征，图7为当D_r=0.3、σ₃=100 kPa时，试样在σ{}_{\mathrm{d}}^{\mathrm{a}\mathrm{m}\mathrm{p}}=10 kPa（塑性安定型）、30 kPa（塑性蠕变型）和40 kPa（增量破坏型）条件下的典型应力‒应变滞回曲线（为便于观察，仅取特定加载次数下的滞回曲线）。

从图7（a）可以看出，当σ{}_{\mathrm{d}}^{\mathrm{a}\mathrm{m}\mathrm{p}}=10 kPa（塑性安定型），在加载次数较小时，滞回曲线间隔较大，随着加载次数的增加，滞回曲线趋于密集，试样的动应变增加缓慢甚至停止增长。从图7（b）可以看出，当σ{}_{\mathrm{d}}^{\mathrm{a}\mathrm{m}\mathrm{p}}=30 kPa（塑性蠕变型）时，初期的滞回曲线间隔较大，随着加载次数的增加，滞回曲线密度减小并逐渐维持恒定，试样的动应变缓慢增加。由图7（c）可知，当σ{}_{\mathrm{d}}^{\mathrm{a}\mathrm{m}\mathrm{p}}=40 kPa（增量破坏型）时，随着加载次数的增加，滞回曲线间隔快速变大，试样的动应变急剧增加。

为更清晰地了解相对密度及围压对试样滞回曲线的影响，图8给出了当N=10 000时，在不同相对密度（σ₃=100 kPa、σ{}_{\mathrm{d}}^{\mathrm{a}\mathrm{m}\mathrm{p}}=30 kPa）及围压（D_r=0.3、σ{}_{\mathrm{d}}^{\mathrm{a}\mathrm{m}\mathrm{p}}=30 kPa）下试样滞回曲线的演化规律。

由图8可知，随着相对密度及围压增加，滞回曲线倾斜程度逐渐增大，试样回弹模量增大，经历了明显的硬化。此外，滞回曲线面积大小可表示每次循环加载过程中试样的能量损耗，由图8还可知，滞回曲线的面积随相对密度及围压的增加逐渐减小，这表明试样内部耗损的能量降低，试样承受动荷载的能力增强。

2.3 回弹模量

回弹模量是能够反映基层和路基刚度的有效指标^[16,22]，如图9所示。

回弹模量M_r定义为^[23]：

M_{\mathrm{r}}=\frac{{\sigma }_{\mathrm{d}}^{\mathrm{a}\mathrm{m}\mathrm{p}}}{{\epsilon }_{\mathrm{1}\mathrm{e}}}（2）

式中，ε_1e为弹性应变。

提取循环荷载下各个滞回曲线的回弹模量，得到不同围压及动应力幅值条件下3种相对密度圆砾的回弹模量M_r与加载次数N之间的关系曲线，如图10所示。从图10可以看出，回弹模量随加载次数的变化趋势与试样的动态变形行为息息相关。当试样处于塑性安定型时，回弹模量表现为在加载初期急剧降低，而后随加载次数的增加逐渐稳定。这主要是因为在加载初期，试样内部颗粒间发生相对滑移，导致试样弹性应变快速增加，进而使试样回弹模量显著减少，随着加载继续，颗粒滑移受限，回弹模量趋于稳定。当试样处于塑性蠕变型时，回弹模量在加载初期也表现出急剧降低，而后随着加载的继续，回弹模量缓慢降低直至稳定。当试样处于增量破坏时，回弹模量随加载次数的增加急剧下降，这是由试样发生破坏，弹性模量急剧增加所致。由图10还可以看出，相对密度及围压越高，回弹模量越大，如：当σ₃=100 kPa、σ{}_{\mathrm{d}}^{\mathrm{a}\mathrm{m}\mathrm{p}}=40 kPa， D_r为0.3、0.5、0.7时，回弹模量分别为48.84、101.14和112.38 MPa，相较于D_r=0.3，D_r=0.5和0.7的回弹模量分别增加了107.1%和131.1%。这说明相对密度和围压的增加能够提高圆砾土的回弹模量和刚度，使其抵抗动荷载的能力得到加强，相应累积应变减小。

3 圆砾动三轴DEM模拟

3.1 数值模型的建立

3.1.1 试样生成

若按圆砾实际级配生成数值试样，生成的颗粒数目高达百万之多，极其耗时且计算机内存无法实现计算。为提高计算效率，将粒径1 mm以下的颗粒用1～2 mm替换；同时，由级配曲线图2可知，粒径为40～60 mm的颗粒极少（占比约4.7%），因此采用等量替换法将40～60 mm的颗粒用5～40 mm的颗粒进行替换。为保证生成试样的均匀性，采用GM法^[24]生成试样，具体过程如图11所示。生成初始试样的高度与宽度分别为400 mm和200 mm。

3.1.2 细观参数标定

参考前人所做粗粒土DEM动三轴模拟，本文DEM本构模型采用线性接触模型^[13,25]。以室内动三轴试验所得累积塑性应变‒振次关系曲线为基础，对模型细观参数进行反复调整标定，结果如图12所示。从图12可以看出，室内试验与数值模拟所得累积塑性应变在数值上基本一致，但由于数值模拟是2D的，而室内试验是3D的，二者结果会存在一定偏差。数值模型具体细观参数见表2。需注意的是，考虑到PFC在进行循环加载时极其耗时，且室内试验每个加载阶段循环加载至1 000次之后累积塑性应变增长幅度极小（前1 000振次所得应变占总应变量90%以上），故在数值试验中，当每个加载阶段分别加载1 000次及应变增幅小于0.005%/100次时，停止数值试验。

3.1.3 试样加载

加载前，对四侧墙体施加围压σ₃对试样进行固结。完成后，采用应力控制法^[26]加载试样，即动应力以正弦变化作用于顶部和底部墙体。加载过程中，通过体积不变法^[27]模拟不排水的过程，如图13所示。

图13中，dx和dy分别为动应力作用下试样的轴向应变和纵向应变，A和y分别为试样的面积和初始高度。此外，MiDi^[28]建议试样在加载过程中应满足式（3），使试样保持准静态状态：

I_{\mathrm{m}}=\dot{\epsilon }d\sqrt[]{\frac{\rho }{p}}\le \mathrm{2.5}\times {\mathrm{10}}^{-\mathrm{3}}（3）

式中，I_m为惯性参数，\dot{\epsilon }为剪切应变率，d为最大粒径，ρ为粒子密度，p为平均应力。

3.2 DEM试验结果分析

3.2.1 配位数

配位数是试样中颗粒的平均接触数，反映了试样内部颗粒的接触程度，配位数越大，说明颗粒体越密实，结构越稳定。根据Thornton^[29]的研究，在离散元中，接触数目少于2的颗粒对试样体系受力基本无贡献。因此本文采用有效配位数（本文简称配位数）描述颗粒间的平均接触数目，其计算公式如下：

C_{\mathrm{n}}=\frac{\mathrm{2}{N}_{\mathrm{c}}-N_{\mathrm{b}\mathrm{1}}}{N_{\mathrm{b}}-N_{\mathrm{b}\mathrm{1}}-N_{\mathrm{b}\mathrm{0}}}（4）

式中，C_n为有效配位数，N_c和N_b分别为总的接触数和总的颗粒数，N_b1和N_b0分别为仅有一个接触的颗粒数和无接触的颗粒数。

以D_r=0.5、σ₃=200 kPa为例，不同动应力幅值下试样配位数C_n随加载次数N的变化曲线如图14所示。由图14可知，随着动应力幅值的增大，配位数越来越小，其上下变化幅度越来越大，这表明动应力越大，试样越疏松。这主要是由于动应力幅值增加，试样内部超孔隙水压力增加，有效应力就越低，试样内部颗粒接触就越弱。此外，在一个正弦加载周期（a—b—c）中，配位数呈反正弦变化（A—B—C），即配位数的变化方向与循环荷载的变化方向相反。

为了便于观察循环荷载下不同相对密度及围压时试样配位数的变化规律，取每个循环加载中配位数的平均值进行分析。不同相对密度及围压下试样的配位数C_n随加载次数N的变化曲线如图15所示。由图15可知：相对密度及围压越大，配位数越大，颗粒间接触愈紧密。此外，相同动应力幅值下，相对密度和围压越小，配位数的下降幅度越大。以σ{}_{\mathrm{d}}^{\mathrm{a}\mathrm{m}\mathrm{p}}=40 kPa为例，当相对密度分别为0.7、0.5、0.3时，配位数的下降幅度分别为0.014、0.022和0.063。这表明相对密度及围压的增加可提高试样的密实性；同时，也从细观层面解释了增加相对密度及围压可提高圆砾强度和抵抗变形能力的原因。

3.2.2 组构的演化

为便于分析试样细观特性，选取试样的4种加载状态，如图16所示。

诸多研究表明，细观组构在控制颗粒材料宏观力学特性上具有重要作用^[30‒31]。颗粒间法向接触力是细观组构的重要参量，其分布（大小和方向）可较好地反映材料的各向异性。为更好地描述颗粒间法向接触力的分布情况，Rothenburg等^[32]提出颗粒间法向接触力分布可采用式（5）所示的傅里叶级数进行拟合：

f_{\mathrm{n}}(\theta )=f_{\mathrm{0}}[\mathrm{1}+a_{\mathrm{n}}\mathrm{c}\mathrm{o}\mathrm{s}\mathrm{2}(\theta -{\theta }_{\mathrm{n}})]（5）

式中：f_n(θ)为颗粒间法向接触力的分布函数，θ为接触法向与水平方向的夹角；f₀为颗粒间所有接触的平均法向接触力；a_n为傅里叶系数（各向异性系数），该值的大小反映了各向异性的复杂程度；θ_n为颗粒间法向接触力各向异性的主方向（与水平方向的夹角）。

以D_r=0.5、σ₃=200 kPa为例，在不同动应力幅值下，颗粒间平均法向接触力分布如图17所示。图17中，颗粒间的法向接触力以10º为一个统计区间，并通过傅里叶函数（式（5））进行拟合，结果表明整体拟合效果良好（虚线为拟合曲线）。

由图17可知，在同一动应力幅值下，随着循环加载状态1变化至状态4，试样的法向接触力分布形状的变化过程为圆形—纺锤形—圆形—纺锤形，相应试样的各向异性变化过程为各向同性—各向异性—各向同性—各向异性。试样加载时，法向接触力分布的主方向集中在约90º方向，与最大主应力方向一致。同时，可以发现，随着动应力幅值的增加，法向接触力越大，试样越易被破坏。此外，动应力幅值越大，各向异性系数a_n越大，即试样的各向异性程度越强，如：在加载状态4时，相较于σ{}_{\mathrm{d}}^{\mathrm{a}\mathrm{m}\mathrm{p}}=10 kPa，σ{}_{\mathrm{d}}^{\mathrm{a}\mathrm{m}\mathrm{p}}=20、30和40 kPa的各向异性程度分别增加了104.72%、215.12%和339.56%。

为进一步地了解循环荷载作用下不同相对密度试样平均法向接触力的分布情况，图18给出了当σ₃=200 kPa、σ{}_{\mathrm{d}}^{\mathrm{a}\mathrm{m}\mathrm{p}}=40 kPa时不同相对密度试样平均法向接触力的分布。

由图18可知：循环荷载在最小值（状态1和状态3）时，平均法向接触力分布呈圆形，表现出各向同性的特征；而循环荷载在最大值（状态2和状态4）时，法向接触力分布呈纺锤形，表现出各向异性特征。同时可发现，在同一加载状态下，相对密度越大，垂直向的平均法向接触力越小，这意味着较小的法向接触应力即可承受同等动应力的作用，表明试样越不易发生破坏，在宏观上表现为相对密度越大，累积应变越小。此外，在同一加载状态下，相对密度越大，各向异性系数a_n越小，即试样的各向异性越弱，如：在加载状态4时，相较于D_r=0.3，D_r=0.5和0.7时试样的各向异性分别下降了14.2%和22.64%。这表明相对密度的增加能降低试样的各向异性，但降低幅度相对较小。

为探究地铁循环荷载作用下，不同围压下试样平均法向接触力分布情况，图19给出了当D_r=0.5、σ{}_{\mathrm{d}}^{\mathrm{a}\mathrm{m}\mathrm{p}}=40 kPa时不同围压下试样平均法向接触力的分布。由图19可知，在同一加载状态下，围压对平均法向接触力影响较大，即：围压越大，法向接触力越大，宏观表现为围压越大，试样承受荷载的能力越强。此外，同一加载状态下，围压越大，各向异性系数a_n越小，即试样的各向异性程度越弱，如：加载状态为状态4时，相较于σ₃=100 kPa，σ₃=200和300 kPa时的各向异性程度分别降低了32.92%和59.5%，这表明围压的增加能够降低试样的各向异性。同时，与相对密度对试样各向异性的影响程度对比可知，围压对试样各向异性的降低效果明显高于相对密度。

为定量研究各向异性对材料力学特性的影响，Satake^[33]提出组构及组构张量的概念对试样的细观结构进行了分析。组构张量定义为：

{\mathit{\Phi }}_{\mathit{i}\mathit{j}}=\frac{\mathrm{1}}{N_c}{\displaystyle \sum _{c=\mathrm{1}}^{N_c}}{n}_i^c{n}_j^c（6）

式中，N_c 为颗粒系统内总的接触数目，n{}_i^c和n{}_j^c分别为第c个接触的单位法向向量在i和j方向上的投影长度。在此基础上，参考广义偏应力的概念，常用偏组构Φ_d=Φ₁‒Φ₃定量地描述接触的各向异性程度，其中Φ₁、Φ₃分别为组构张量的最大、最小主值。

不同动应力幅值下，试样（D_r=0.5、σ₃=200 kPa）的偏组构随加载次数的变化曲线如图20所示。由图20可知：随着动应力幅值的增大，偏组构越来越大，其上下变化幅度也越来越大，这表明试样的各向异性程度越强。此外，在一个正弦加载周期中（a—b—c），偏组构的变化方向呈正弦变化（A—B—C），即偏组构的变化方向与循环荷载的变化方向相同。

为便于观察循环荷载下不同相对密度及围压下试样偏组构变化规律，取每个循环中偏组构的平均值进行分析。不同相对密度及围压下试样的偏组构随加载次数的变化如图21所示。由图21可知：相对密度及围压越大，偏组构越小，试样的各向异性越弱；在同一动应力幅值下，相对密度和围压越大，偏组构变化幅度越小，表明增大相对密度及围压能降低试样各向异性程度。根据Radjai等^[34]提出的方法，可用平均接触力将试样整体接触分为强接触和弱接触。

不同加载状态下，试样（D_r=0.5、σ₃=200 kPa）整体接触、强接触和弱接触的法向分布（接触法向）如图22所示。

由图22可知：状态1时，强接触沿轴向分布较弱，弱接触分布较均匀，整体接触呈各向同性；状态2时，强接触沿轴向分布十分明显，弱接触沿水平方向发展，整体接触沿轴向发展；状态3和4时，接触分布分别与状态1和2类似，但接触沿轴向分布程度略强于状态1和2。这表明强接触主要由轴向应力控制，其对承载力起主要支撑作用，而弱接触主要由围压控制。

4 结 论

以南宁地铁地基圆砾为研究对象，利用DYNTTS大型三轴循环试验系统，开展了一系列饱和不排水动三轴试验，分析了相对密度、围压及动应力幅值对圆砾的动力特性及累积塑性应变的影响；同时，基于圆砾室内试验，利用PFC2D内嵌fish语言，编写了相应的DEM循环双轴不排水模型程序，从细观角度揭示了不同相对密度、围压及动应力幅值下圆砾的动力响应机制。主要结论如下：

1）对于埋深较浅（σ₃=100 kPa）且地基相对密度较小（D_r=0.3）的地铁隧道，圆砾的累积塑性应变较大（高于0.61%），并有急剧上升的趋势，严重威胁列车行车安全，应对圆砾地基的强度进行提升。当地铁隧道埋深较深（σ₃=200、300 kPa）时，累积塑性应变较小（低于0.14%），列车的运行相对安全，不需要对圆砾地基进行处理。

2）增加围压和相对密度，能够显著降低试样的累积塑性应变（降低幅度约77.05%），提高试样的回弹模量（最大提高幅度约131.1%），表明高压注浆等挤压措施对提高土体相对密实度以减小列车运行期间隧道的沉降十分有效。

3）加载阶段，配位数的变化方向与循环荷载的变化方向相反，呈反正弦变化。随着动应力幅值的增大，配位数越来越小，颗粒间接触越疏松。相对密度及围压越大，配位数越大，试样颗粒间接触越紧密，在宏观上体现为试样承受动荷载能力越强。

4）在循环加载过程中，各向异性的变化方向与循环荷载的变化方向相同，呈正弦变化。随着动应力幅值增加，试样的各向异性程度增大；随着相对密度及围压增加，各向异性程度减小。强接触对承载力起着主要支撑作用，其主要由轴向应力控制，而弱接触主要由围压控制。

参考文献

[1] Ma Shaokun, Duan Zhibo, Huang Zhen,et al.Study on the stability of shield tunnel face in clay and clay-gravel stratum through large-scale physical model tests with transparent soil[J].Tunnelling and Underground Space Technology,2022,119:104199. doi:10.1016/j.tust.2021.104199

[2] Tang Kaishun, Xie Xiongyao, Yang Lei.Research on mech-anical characteristics of gravel soil based on large-scale triaxial tests[J].Chinese Journal of Underground Space an-d Engineering,2014,10(3):580‒585.

[3] 唐开顺,谢雄耀,杨磊.圆砾土大型三轴试验力学特性研究[J].地下空间与工程学报,2014,10(3):580‒585.

[4] Yang Haihua, Liu Liang, Liu Hanlong,et al.Tests for dyna-mic modulus and damping ratio of high polymer cementitious Gobi soil[J].Journal of Vibration and Shock,2023,42(3):12‒20. doi:10.13465/j.cnki.jvs.2023.03.002

[5] 杨海华,刘亮,刘汉龙,等.高聚物胶凝戈壁土的动模量及阻尼比试验研究[J].振动与冲击,2023,42(3):12‒20. doi:10.13465/j.cnki.jvs.2023.03.002

[6] Cui Yujun.Mechanical behaviour of coarse grains/fines m-ixture under monotonic and cyclic loadings[J].Transportation Geotechnics,2018,17:91‒97. doi:10.1016/j.trgeo.2018.09.016

[7] Long Yao, Zhang Jiasheng, Ding Jianyuan,et al.The cyclic load experiments and an accumulated deformation model for coarse-grained soil filling[J].Journal of Vibration and Shock,2017,36(16):128‒133.

[8] 龙尧,张家生,丁建源,等.粗粒土路基循环荷载试验及累积变形模型研究[J].振动与冲击,2017,36(16):128‒133.

[9] Leng Wuming, Zhou Wenquan, Nie Rusong,et al.Analysis of dynamic characteristics and accumulative deformation of coarse-grained soil filling of heavy-haul railway[J].Rock and Soil Mechanics,2016,37(3):728‒736. doi:10.16285/j.rsm.2016.03.015

[10] 冷伍明,周文权,聂如松,等.重载铁路粗粒土填料动力特性及累积变形分析[J].岩土力学,2016,37(3):728‒736. doi:10.16285/j.rsm.2016.03.015

[11] Nie Rusong, Sun Baoli, Cheng Longhu,et al.Resilient characteristics test of coarse-grained soil filler for heavy haul railway subgrade[J].Journal of the China Railway Society,2022,44(1):96‒104.

[12] 聂如松,孙宝莉,程龙虎,等.重载铁路路基粗粒土填料回弹特性试验研究[J].铁道学报,2022,44(1):96‒104.

[13] Duong T V, Tang A M, Cui Yujun,et al.Effects of fines and water contents on the mechanical behavior of interlayer soil in ancient railway sub-structure[J].Soils and Foundations,2013,53(6):868‒878. doi:10.1016/j.sandf.2013.10.006

[14] Leng Wuming, Xiao Yuanjie, Nie Rusong,et al.Investigating strength and deformation characteristics of heavy-haul railway embankment materials using large-scale undrained cyclic triaxial tests[J].International Journal of Geomechanics,2017,17(9):04017074. doi:10.1061/(asce)gm.1943-5622.0000956

[15] Yang Zhihao, Yue Zurun, Tai Bowen.Investigation of the deformation and strength properties of fouled graded macadam materials in heavy-haul railway subgrade beds[J].Construction and Building Materials,2021,273:121778. doi:10.1016/j.conbuildmat.2020.121778

[16] Sun Yifei, Song Shunxiang, Gao Yufeng.Discrete element simulation of the grading dependence of the cyclic behavior of granular aggregates[J].Journal of China University of Mining & Technology,2018,47(4):874‒878.

[17] 孙逸飞,宋顺翔,高玉峰.级配散粒体循环荷载下变形特性的离散元模拟[J].中国矿业大学学报,2018,47(4):874‒878.

[18] Zhang Fuguang, Nie Zhuochen, Chen Mengfei,et al.DEM analysis of macro-and micro-mechanical behaviors of cemented sand subjected to undrained cyclic loading[J].Chinese Journal of Geotechnical Engineering,2021,43(3):456‒464. doi:10.11779/CJGE202103008

[19] 张伏光,聂卓琛,陈孟飞,等.不排水循环荷载条件下胶结砂土宏微观力学性质离散元模拟研究[J].岩土工程学报,2021,43(3):456‒464. doi:10.11779/CJGE202103008

[20] Nguyen N S, François S, Degrande G.Discrete modeling of strain accumulation in granular soils under low amplitude cyclic loading[J].Computers and Geotechnics,2014,62:232‒243. doi:10.1016/j.compgeo.2014.07.015

[21] Nie Zhihong, Qi Qun, Wang Xiang,et al.DEM investigation of strain behaviour and force chain evolution of gravel‒sand mixtures subjected to cyclic loading[J].Particuology,2022,68:13‒28. doi:10.1016/j.partic.2021.10.006

[22] Tang Yiqun, Wang Yanling, Huang Yu,et al.Dynamic str-ength and dynamic stress‒strain relation of silt soil under traffic loading[J].Journal of Tongji University,2004,32(6):701‒704.

[23] 唐益群,王艳玲,黄雨,等.地铁行车荷载下土体动强度和动应力‒应变关系[J].同济大学学报(自然科学版),2004,32(6):701‒704.

[24] Tang L, Yan M H, Ling X Z,et al.Dynamic behaviours of railway's base course materials subjected to long-term low-level cyclic loading:Experimental study and empirical mo-del[J].Géotechnique,2017,67(6):537‒545. doi:10.1680/jgeot.16.p.152

[25] Wang Jiaquan, Zhu Mengke, Lin Zhinan,et al. Influence of fines content on static and dynamic characteristics of saturated gravelly soil[J].China Civil Engineering Journal,2023,56 (5):112‒121.

[26] 王家全,祝梦柯,林志南,等.细粒含量对饱和砾性土静动力学特性的影响[J].土木工程学报,2023,56(5):112‒121.

[27] Qi Shuai, Cui Yujun, Dupla J C,et al.Investigation of the parallel gradation method based on the response of track-bed materials under cyclic loadings[J].Transportation Geotechnics,2020,24:100360. doi:10.1016/j.trgeo.2020.100360

[28] Indraratna B, Nimbalkar S, Rujikiatkamjorn C.From theory to practice in track geomechanics—Australian perspective for synthetic inclusions[J].Transportation Geotechnics,2014, 1(4):171‒187. doi:10.1016/j.trgeo.2014.07.004

[29] 中华人民共和国铁道部. 铁路工程土工试验规程:TB 10102—2010 [S].北京:中国铁道出版社,2010.

[30] Werkmeister S, Dawson A R, Wellner F.Permanent deformation behaviour of granular materials[J].Road Materials and Pavement Design,2005,6(1):31‒51. doi:10.3166/rmpd.6.31-51

[31] Erlingsson S, Rahman M S.Evaluation of permanent deformation characteristics of unbound granular materials by means of multistage repeated-load triaxial tests[J].Transportation Research Record:Journal of the Transportation Research Board,2013,2369(1):11‒19. doi:10.3141/2369-02

[32] Rollins K M, Evans M D, Diehl N B,et al.Shear modulus and damping relationships for gravels[J].Journal of Geotechnical and Geoenvironmental Engineering,1998,124(5):396‒405. doi:10.1061/(asce)1090-0241(1998)124:5(396)

[33] Duan N, Cheng Y P.A modified method of generating sp-ecimens for a 2D DEM centrifuge model[C]//Proceedings of the Geo-Chicago 2016.Chicago:American Society of Civil Engineers,2016:610‒620. doi:10.1061/9780784480151.060

[34] Wang Yong, Wang Yanli.Liquefaction characteristics of gr-avelly soil under cyclic loading with constant strain amplitude by experimental and numerical investigations[J].Soil Dynamics and Earthquake Engineering,2017,92:388‒396. doi:10.1016/j.soildyn.2016.10.029

[35] Xu X M, Ling D S, Cheng Y P,et al.Correlation between liquefaction resistance and shear wave velocity of granular soils:A micromechanical perspective[J].Géotechnique,2015,65(5):337‒348. doi:10.1680/geot.sip.15.p.022

[36] Wang Rui, Fu Pengcheng, Zhang Jianmin,et al.DEM study of fabric features governing undrained post-liquefaction s-hear deformation of sand[J].Acta Geotechnica,2016,11(6):1321‒1337. doi:10.1007/s11440-016-0499-8

[37] MiDi G.On dense granular flows[J].The European Physical Journal E,2004,14(4):341‒365. doi:10.1140/epje/i2003-10153-0

[38] Thornton C.Numerical simulations of deviatoric shear deformation of granular media[J].Géotechnique,2000,50(1):43‒53. doi:10.1680/geot.2000.50.1.43

[39] Yang Jiming, Zhang Xiaoyong, Zhang Fuyou,et al.Mesoscopic study on bearing characteristics of pile foundation under pile‒soil‒cap combined interaction in sand[J].Rock and Soil Mechanics,2020,41(7):2271‒2282.

[40] 杨济铭,张小勇,张福友,等.砂土中桩‒土‒承台协同作用下桩基承载特性细观研究[J].岩土力学,2020,41(7):2271‒2282.

[41] Liu Yang, Zhang Duo, Wu Shunchuan,et al.DEM investigation on the evolution of fabric under true triaxial conditions in granular materials[J].International Journal of Geomechanics,2020,20(8):04020110. doi:10.1061/(asce)gm.1943-5622.0001740

[42] Rothenburg L, Bathurst R J.Analytical study of induced anisotropy in idealized granular materials[J].Géotechnique,1989,39(4):601‒614. doi:10.1680/geot.1989.39.4.601

[43] Satake M.The role of the characteristic line in static soil behavior[C]//IUTAM Symposium on Deformation and Fa-ilure of Granular Materials.Delft:A Balkema,1982:63‒68.

[44] Radjai F, Jean M, Moreau J J,et al.Force distributions in dense two-dimensional granular systems[J].Physical Rev-iew Letters,1996,77(2):274‒277. doi:10.1103/physrevlett.77.274

基金资助

国家自然科学基金项目(52268062)

广西自然科学基金重点项目(2020GXNSFDA238024)

广西重点研发计划项目(GUIKE AB22080061)

广西交通运输行业重点科技项目（GXJT‒2020‒02‒08）