闽东地区含砾花岗岩残积土的细观剪切特性研究

张宏虎，白 伟，孙明祥，邓 涛

(1.中铁南方投资集团有限公司, 广东 深圳 518052；2.福州大学 土木工程学院, 福建 福州 350108)

花岗岩残积土由花岗岩母岩受长期物理、化学风化作用形成，在我国华东、华南地区有着广泛分布。花岗岩残积土的细观结构复杂、矿物成分多样，尤其是砾石成分对其工程性质有重要影响[1-3]。许多学者通过室内(或现场)试验、数值模拟等手段，对花岗岩残积土的剪切特征展开了细致研究。

胡屏等[4]通过对福州地区花岗岩残积土进行不同粗颗粒含量下的室内大型直剪试验，发现粗颗粒含量的增加，花岗岩残积土的黏聚力呈先减小后增大的趋势。曾庆有[5]、龙志东等[6]分别通过现场大型直剪试验和室内直剪试验，研究花岗岩残积土中粒径成分对抗剪强度参数的影响，表明填料内摩擦角随着粗颗粒含量以及压实度呈近似线性增长。徐晓宇等[7]考虑剪切速率、加载级数和破坏值等因素，对广州地区花岗岩残积土开展室内直剪试验，表明当剪切速率为0.8 mm/min时最能反映原状土的强度特征。Li[8]通过室内大型直剪试验研究混合土和砂砾石的剪切特性，探讨含石量和粗粒形状对抗剪强度的影响。

Graziani等[9]对土石混合体进行了二维直剪试验的颗粒流数值模拟，分析含石量、块石形状和围压等因素对剪切特性的影响规律。徐文杰等[10]通过三维扫描技术建立块石三维形态数据库，进而开展颗粒流数值模拟分析土石混合体中块石效应。杨升等[11]在PFC3D中构建哑铃形和椭圆形两种簇颗粒，探讨颗粒形状对抗剪强度的影响，认为含簇颗粒试样的抗剪强度和剪胀性显著增大。Xu等[12]通过数字图像处理技术生成土石混合体真实模型，探讨内部块石空间分布特征对剪切带形成的影响。赵仕威等[13]和金磊等[14]通过改进随机生成簇颗粒的建模方法，研究颗粒棱角度对剪切强度和剪胀性的影响，以及接触力在剪切过程中的演化规律。张强等[15]通过不规则块石三维重构，建立土石混合体数值模型，探讨块石不规则形状和空间分布对其力学特性的影响。

综合前人研究基础，本文首先通过数字图像处理技术，研究闽东花岗岩残积土砾石的颗粒形状特征，并在PFC中对砾石颗粒进行二维重构；进而对现场大型直剪试验的全过程展开细致的数值模拟，探讨含砾花岗岩残积土的细观剪切特性。

1 砾石颗粒形状特征及颗粒重构

离散元数值模拟中，簇颗粒形状对含砾花岗岩残积土的剪切特性有重要影响。以往采用概化几何形状随机生成的簇颗粒，不能真实的反映闽东地区花岗岩残积土砾石颗粒的形状特征。本文通过数字图像处理技术对花岗岩残积土砾石颗粒的真实形状特征进行统计，构建花岗岩残积土砾石颗粒样本库，进而遴选典型颗粒形状，在PFC中对真实砾石颗粒进行簇颗粒重构[16]。具体步骤如下：

Step1：通过数码相机对砾石颗粒的图像进行采集，部分真实颗粒如图1所示。

图1 现场砾石颗粒

Step2：通过图像处理技术得到颗粒二维轮廓(如图2所示)，并将其转化为一维线轮廓的周期函数，进而采用傅里叶级数对轮廓周期函数进行拟合，得到该颗粒的傅里叶描述子f。

图2 砾石颗粒轮廓曲线

Step3：基于傅里叶描述子f结合有限维向量理论，计算颗粒形状相似度Da，并遴选Da值最大的25种颗粒作为砾石颗粒典型代表。

Step4：在PFC平台下，根据25种典型颗粒的傅里叶描述子f，反算得到颗粒二维轮廓，导入PFC构建25种典型簇颗粒，如图3所示。

2 大型直剪试验颗粒离散元数值模拟

2.1 现场大型直剪试验

闽东地区花岗岩残积土分布广、堆积厚度大，以福州滨海快线A1标段为例，盾构穿越花岗岩残积土地层中含有大量砾石颗粒，通过现场土样的筛分试验得到其天然级配曲线如图4所示。其砾石含量(过5 mm筛后的残留颗粒)达到29.7%，其不均匀系数Cu=5.57，曲率系数Cc=1.12，级配良好。同时通过室内密度试验，可得砾石颗粒(粒径大于5 mm)的密度为ρ=2 400 kg/m3，土颗粒平均密度ρ=2 000 kg/m3。

图3 重构后的砾石簇颗粒

图4 花岗岩残积土的级配曲线

根据前人完成的现场大型直剪试验[5]，得到竖向压力0 kPa、160 kPa、320 kPa、560 kPa下，花岗岩残积土的剪应变与剪应力关系曲线，如图5所示。

图5 剪应力-剪应变曲线(现场直剪试验)

根据各曲线的峰值点(A、B、C、D)绘制抗剪强度包线，得到其抗剪强度τ与竖向应力σ的关系式为：τ=0.49σ+6.4；从而得到其抗剪强度指标为：内摩擦角φ=26.1°，黏聚力c=6.4 kPa。

2.2 数值模型建立

数值模拟中剪切盒由上下2面墙体与左右4面墙体构成，其整体尺寸为0.5 m×0.5 m(与现场直剪试验一致)。考虑到数值模型的计算效率，忽略了粒径0.5 mm以下的细颗粒，得到修正级配如图4所示。

在剪切盒范围内首先生成粒径为10 mm～5 mm的簇颗粒，其质量百分比为29.7%；再通过FISH语言生成粒径为5.0 mm～0.5 mm的球颗粒进行填充，最终得到含砾花岗岩残积土大型直剪试验的数值模型，其中簇颗粒180个，球颗粒10 882个，共计11 062个颗粒，如图6所示。

图6 直剪试验的数值模型

上下墙体施加竖向伺服应力分别为160 kPa、320 kPa、560 kPa，待试样达到初始平衡后，保持下盒不动，对上盒墙体施加0.01 mm/step的水平速度，以模拟土样的剪切过程，最终剪切位移达到30 mm时终止试验。

2.3 细观参数标定

数值模拟中所有颗粒之间(包括球颗粒Ball和簇颗粒Clump)的接触，所有颗粒与墙体之间的接触，均采用线性接触本构。以含砾花岗岩残积土大型现场直剪试验结果为参考，采用逼近法对数值模拟中所需的细观参数进行标定，其结果汇总于表1。

图7给出了不同竖向应力下，含砾花岗岩残积土剪切过程中的应力-应变曲线，可以认为数值模拟结果与现场直剪试验结果吻合较好。

(数值模拟与现场直剪试验)

3 含砾花岗岩残积土的细观剪切特性

从物理试验中仅能得到土体的宏观剪切特性(尤其是抗剪强度指标)，而基于颗粒离散元的数值模拟，则可对含砾花岗岩残积土的细观剪切特性展开深入探讨。

3.1 颗粒位移变化分析

以第二组数值模拟(竖向应力320 kPa)为例，剪切终止时，其球颗粒和簇颗粒的水平位移分布如图8所示。试样上下两端的球/簇颗粒，其水平位移分布较为集中(以刚体位移为主)，而试样中部(剪切带范围内)的球/簇颗粒，其水平位移分布较为分散。

为更准确计算剪切带宽度，本文重点关注球/簇颗粒与墙体之间的相对水平位移。以25 mm为一个单位高程，计算该高程范围内球/簇颗粒的平均相对水平位移。仍以第二组数值模拟(竖向应力320 kPa)为例，剪切终止时，其球颗粒和簇颗粒的平均相对水平位移分布如图9所示。可以看出，球颗粒平均相对水平位移沿高程大致呈正态分布，因此可用正态曲线对其平均相对水平位移进行拟合，如图9中虚线所示。

由于下剪切盒保持不动，上剪切盒施加水平剪切位移，因此正态拟合曲线显示的剪切带中心并不在剪切面位置，而是向上偏移u约为9.80 mm。进一步地，本文将剪切带中心上下一个标准差范围定义为剪切带宽度；由此可准确计算含砾花岗岩残积土的剪切带宽度w约为61.38 mm。

图8 颗粒水平位移分布

图9 颗粒平均相对水平位移分布

对其他两组数值模拟结果也进行类似的细观统计分析，得到不同竖向应力下，含砾花岗岩残积土的剪切带偏移量及剪切带宽度，如表2所示。由表2可知，竖向应力越大，剪切终止时，其偏移量u与剪切带宽度w逐渐减小，且减小速率逐渐趋缓；在同一竖向应力下，根据球颗粒与簇颗粒统计得到的剪切带偏移量u与剪切带宽度w几乎相同。

表2 不同竖向应力下剪切带偏移量和剪切带宽度(基于颗粒相对位移统计)

3.2 颗粒旋转变化分析

颗粒的旋转也是反映剪切带特性的重要指标。以25 mm为一个单位高程，计算该高程范围内球/簇颗粒的平均旋转量。仍以第二组数值模拟(竖向应力320 kPa)为例，剪切终止时，其球颗粒和簇颗粒的平均旋转量分布如图10所示。

图10 颗粒平均旋转量分布图

可以看出，球颗粒平均旋转量沿高程大致呈正态分布，因此可用正态曲线对其平均旋转量进行拟合，如图10中虚线所示。同样地，根据3.1节中确定剪切带中心与剪切带宽度的方法，则根据球颗粒和簇颗粒统计得到的剪切带偏移量u分别约为20.01 mm、23.38 mm，剪切带宽度w分别约为89.82 mm、92.13 mm。

需要说明的是，由于簇颗粒具有不规则的几何形状和粗糙的表面，在剪切过程中与其他颗粒的接触作用更强，限制了簇颗粒自身运动，因此簇颗粒的平均旋转量明显小于球颗粒。这也可能是造成二者统计所得结果(剪切带中心偏移和剪切带宽度)不一致的重要原因。

对其余两组数值模拟结果也进行类似的细观统计分析，得到不同竖向应力下，含砾花岗岩残积土的剪切带偏移量及剪切带宽度，如表3所示。对比表2可知，相比基于颗粒相对位移统计结果，基于颗粒旋转量统计得到的剪切带偏移量和剪切带宽度显著偏大。

表3 不同竖向应力下剪切带偏移量和剪切带宽度(基于颗粒旋转量统计)

鉴于使用球颗粒/簇颗粒的平均相对位移统计所得结果的一致性较好，建议使用基于平均相对水平位移的统计结果，来表征含砾花岗岩残积土的剪切带宽度及剪切带中心偏移。

3.3 接触力链分析

数值模拟中的颗粒接触包括Ball-ball接触(细颗粒间)和Ball-clump接触(细颗粒和砾石颗粒间)。若将接触力大于1.5 kN定义为强接触，图11绘制了不同剪切位移下，颗粒间强接触力链的分布示意图(为便于观察，图中仅展示砾石颗粒，细颗粒未画出)。竖向压力伺服完成后，颗粒间强接触力链分布稀疏；当剪切位移行进至15 mm时，强力链沿45°方向基本贯通；当剪切位移行进至30 mm时，强力链几乎遍布全部试样，其优势力链方向仍为45°。

同时，强力链主要集中在砾石颗粒间；而细颗粒则填充试样孔隙，使得各颗粒间咬合更加充分。二者共同作用使得试样在剪切过程中能够发生连续的挤压和变形，利于形成较为稳定的强力链，增强了试样的抗剪强度。上述接触力链分析也从细观机理上，佐证了前人关于砾石成分对花岗岩残积土抗剪强度指标影响的试验结果[2-4]。

(竖向压力320 kPa)

4 结 论

以闽东地区含砾花岗岩残积土为对象，首先通过数字图像技术细致描述了砾石的颗粒形状特征，并在PFC中对砾石颗粒进行二维重构；进而对现场大型直剪试验的全过程展开细致的数值模拟，探讨含砾花岗岩残积土的细观剪切特性。

球颗粒/簇颗粒的平均相对水平位移，沿高程大致呈正态分布。其剪切带中心不在剪切面位置，而是略微向上偏移，上述剪切带中心偏移量和剪切带宽度均随竖向压力的增大而逐渐减小。

球颗粒/簇颗粒的平均转动量，沿高程大致呈正态分布；但基于颗粒旋转量统计得到的剪切带偏移量和剪切带宽度显著偏大，且球颗粒和簇颗粒统计所得结论的一致性较差。因此仍建议使用基于平均相对水平位移的统计结果，来表征含砾花岗岩残积土的剪切带宽度及剪切带中心偏移。

颗粒间接触力链的分析表明，强力链主要集中在砾石颗粒间；而细颗粒则使各颗粒间咬合更加充分，利于形成较为稳定的强力链，增强了试样的抗剪强度。