直剪试验下砂土剪切带形成的细观变形机制

潘远阳，魏玉峰，李远征，杨 豪，雷 壮

(成都理工大学 地质灾害防治与地质环境保护国家重点实验室，四川 成都 610059)

砂土是自然界中较为典型的岩土体材料，广泛应用于工程建设领域中，具有离散性、各向异性、非均匀性等特殊性质。砂土介质在外荷载作用下的力学特性及变形演化机制，一直以来都是岩土力学领域研究的热点，其中关于剪切带的研究取得了一系列丰富的成果。目前关于剪切带的理论研究多偏向于三轴试验，且大部分研究中未对剪切带的形成及厚度等问题进行重点分析。Vardoulakis等[1]基于X射线对砂土双轴试验中剪切带的厚度进行了深入研究，发现剪切带的厚度受试样平均粒径的影响。Tatsuoka等[2]基于摄影测量技术对平面应变压缩试验的应变场进行了研究，得到了类似结论；同时发现在峰值应力后开始发生明显的应变局部化现象，最终贯通形成剪切带，且峰值后应变主要集中于剪切带内，剪切带外则发生轻微的变化。Potts等[3-4]采用弹塑性模型，基于有限元法对砂土直剪试验的力学特征及剪切带的应变局部化问题进行了分析。但有限元法无法揭示颗粒的受力及变形演化规律，也不能反映颗粒孔隙率等变化情况，而离散元法在一定程度上可以解决这些问题。Masson等[5]运用离散元法对松散及密实试样的剪切力学行为进行了宏、细观分析研究。Gutierrez等[6]基于离散元法从微观力学角度研究了粒状材料剪切过程中的尺寸效应及其对试样宏观性质的影响。郭聚坤等[7]基于粗糙度对结构物细砂界面的剪切特性进行了研究，发现剪切面的形态与结构物的粗糙度有关。Oda等[8]在研究中发现，颗粒材料微观变形机制的主要影响因素是颗粒间的转动，而不是颗粒间的滑动。孙其诚等[9]从接触力与接触网络结构方面，对颗粒介质剪切过程中的力链进行了深入研究。可见，利用离散元法对砂土介质的细观力学行为进行研究，可行性良好，也是对物理试验的有效补充。

本文基于离散元法对密实砂土直剪试验进行数值模拟，主要对剪切过程中的核心区域（剪切带）进行分析，从试样特定状态的受力与变形宏细观特征入手，重点分析直剪试验中颗粒试样的局部变形及剪切带的形成演化特征；从颗粒的接触、相互作用及其运动过程角度分析试样整体的细观力学特性；并结合颗粒试样的配位数、平均转动量及孔隙变化率等细观参数分析剪切带内外区域的差异。

1 直剪试验数值模拟

1.1 数值模型

数值试验中剪切盒的尺寸为80 mm×60 mm，颗粒形状为圆形，粒径介于0.6～0.8 mm内均匀分布，平均粒径d50为0.7 mm，墙体摩擦系数为0，试样初始孔隙率为0.16，颗粒共计9 683个。模型中采用6道刚性墙体模拟物理试验的剪切盒，为防止颗粒在剪切过程中溢出，在上下剪切盒中间生成2道水平刚性墙体。颗粒间的接触为接触刚度模型，主要细观参数有颗粒间的接触刚度、颗粒与墙体间的接触刚度、孔隙率及密度等，参考Wang[10]及蒋明镜等[11]的数值试验参数并进行了适当调整，具体取值见表1。

为保证模型有较好的收敛性，本试验中颗粒间及颗粒与墙体间的法向接触刚度与切向刚度比值取1.4。颗粒试样生成后，先通过伺服控制系统使试样处于固结稳定状态，随后固定下剪切盒，且保持法向应力恒定，对上盒施加1个水平向右的速度进行剪切（模拟中施加的速度等效于力的作用），以模拟密砂物理直剪试验，数值模型如图1所示。

为分析试样在外荷载作用下的剪切变形发展情况，本文在室内物理试验的基础上，对离散元模型试样分别施加50，100，200和300 kPa共4个等级的法向荷载，以进行数值试验研究；通过反复调试数值模型中的细观参数，直至数值模拟试验结果与物理试验结果基本一致(见图2)，进而保证后文数值模拟试验成果的可靠性。从物理试验结果可发现，各级正应力下的剪应力差距接近倍数关系主要由试样间的黏聚力、内摩擦角及仪器误差所引起。

图 1 数值模型Fig. 1 Numerical model

图 2 数值试验与物理试验结果对比Fig. 2 Comparison of numerical and physical test results

1.2 应力应变关系

由于在常规直剪试验中，剪切盒前后两侧边界不发生相对运动，即剪切过程中试样的体积只沿竖直方向发生膨胀或压缩，试验中可根据试样整体的纵向位移变化反映其剪胀性。利用数值试验可记录剪切过程中的剪应力、剪切位移及体积变化情况，从而导出试样在100 kPa法向力作用下的应力-应变曲线（图3）。

从图3可看出：BC段为应变硬化阶段，C点为峰值状态，C点后试样处于应变软化阶段，即峰值应力后试样的变形模式失去原来的均匀性而呈现软化特性，这将引起应变局部化现象。其中，试样的峰值强度为50 kPa，剪应力与法向应力的比值为5.0，残余强度约28 kPa，残余应力比为2.8。由剪应力-剪切位移及体积变化曲线可知，试样在剪切过程中经历了A～F共6个时刻，基于此可对试样的宏、细观力学行为及结构演化进行研究，表2为6个时刻对应的详细特征。其中A点为初始剪切状态，B点为对应剪缩体积应变率为零时刻，C为应力峰值点，该点体积应变率最大；F点为残余状态，体积应变率也为零。由体积变化曲线知，试样在剪切过程中先发生体积剪缩后膨胀，整体表现出剪胀性，结果与赵联桢等[12]对冻砂土接触面循环剪切性能研究中出现的剪胀现象相似。从细观角度看，试样的宏观变形与剪切过程中颗粒的滑移和转动有密切联系；在峰值强度时颗粒试样完成第一次排列，随后在达到残余强度时完成颗粒的重排列，此时颗粒的运动对试样的剪胀性影响较小，试样体积趋于稳定。

图 3 应力-位移及体积变化曲线Fig. 3 Curves of stress-displacement and volume change

表 2 剪切状态特征Tab. 2 Shear state feature

2 剪切带形成的宏观特征

2.1 试样局部变形

在岩土体的力学性质中，非连续、不均匀及局部变形问题是重点研究的对象。剪切试验中由于试样内部某点率先达到峰值状态，随后软化而产生较大应变，这将引起应变局部化现象[13]。为研究颗粒试样的变形及发展情况，在模拟试验前，对剪切盒中的试样进行了标记，即设置了9个相同宽度的纵向条块，监测到试样的变形演化情况如图4所示。

由图4可看出，未剪切前试样均匀分布，纵列无明显变形，主要为纯剪切模式。剪切开始后，靠近剪切盒左侧边界的条块开始产生轻微裂纹，其余条块无明显变形；当剪切应力达到峰值状态时，试样中的裂纹进一步扩展，且在盒体中部出现了局部变形现象，各条块中部均表现为凹向左侧，主要是盒体左侧边界受到伺服系统施加的速度所致。峰值状态后，试样变形表现出非线性特征，随着剪应力的继续作用，各条块的变形由左至右逐渐发育，且变形局部化现象越来越显著。剪切结束后，变形主要集中于上下剪切盒之间的剪切面内。剪切面左侧边界处变形最大，向右逐渐递减。综上可知，试样在纯剪切阶段整体变形较小，峰值状态时整体变形明显，峰值状态后变形以非线性应力应变为主导，由左侧边界向右逐渐发育，最后在剪切面形成狭长的条带状变形区域。

图 4 试样局部变形Fig. 4 Sample local deformation

2.2 剪切带形成的形态特征

2.2.1 剪切带形状 由于直剪试验的剪切面限定在上下盒之间的平面，故剪切带的分布也主要集中于此区域。剪切过程中，试样发生局部化变形，将会在水平剪切面内形成狭长的条带状区域（即剪切带），为观测剪切带的形成与演化过程，数值试验中对剪切盒内的颗粒试样进行了动态监测，结果如图5所示。

分析图5可见，未剪切前试样均匀分布，且无明显变形现象。剪切开始后，靠近剪切盒中部的区域开始产生轻微裂纹，且逐步扩展；当剪应力达到峰值状态C点时，裂纹进一步扩展并在盒体中部贯通形成狭长的带状破坏区域（即剪切带）；且剪切带表现为凹向盒体下侧，主要是由于试样在剪切过程中的相互挤压和错动，剪切带内的颗粒以滚动和翻转的形式跨越盒体下部的试样所致。峰值状态后，颗粒试样处于软化状态，随着剪切应力的增大，剪切带的宽度有所增加；随后在法向荷载的约束下，剪切带的宽度略微变窄，最后逐渐趋于稳定发展，直至试验结束。

图 5 剪切带形成及演化Fig. 5 Formation and evolution of shear bands

2.2.2 倾角 剪切带的倾角是指剪切带内最大竖向应变区域与水平剪切面的夹角，即直剪试验中剪切带与剪切面之间的夹角。根据数值试验中的模拟结果进行测量，发现剪切带与剪切面的夹角较小；峰值状态前剪切带未贯通且主要集中于剪切面左右两侧，此时剪切变形区域与剪切面的倾角约为10°；峰值时刻剪切带逐步贯通且在剪切盒中部趋于近水平状态，剪切带与剪切面的夹角约为13°；峰值状态后剪切带与剪切面的夹角逐渐增大，最终稳定在11°左右。

2.2.3 剪切带厚度 剪切带的厚度是指剪切面中最大竖向变形区域所包括范围的距离，剪切带厚度体现颗粒试样对外荷载的响应。为了揭示剪切带厚度的影响，数值模拟试验中通过分别改变颗粒粒径、刚度、孔隙率及颗粒摩擦系数等参数的大小，观测剪切带厚度的变化情况。研究结果表明，影响剪切带厚度的主要因素有平均粒径d50及颗粒摩擦系数f。

由图6可发现，当平均粒径d50恒定，改变颗粒摩擦系数的大小时，可发现随着颗粒摩擦系数的增大，剪切带厚度逐渐变窄，主要是由于剪切过程中颗粒间的咬合力和摩擦阻力增大，颗粒试样的运动将受到约束，从而导致颗粒的滑移和翻转幅度减弱。当颗粒摩擦系数恒定，改变平均粒径d50的大小时，可发现随着平均粒径的增大，剪切带厚度在初始阶段急剧增大，当平均粒径增大到原来的2倍时，剪切带的厚度达到最大值，此时剪切带的厚度约为21d50；随后改变平均粒径的大小，发现剪切带的厚度逐渐变小，且在平均粒径增大到原来的5倍时，剪切带的厚度达到最小值，此时剪切带的厚度约为12d50；随后增大平均粒径，剪切带的厚度处于平稳发展的趋势。

图 6 剪切带厚度与平均粒径及摩擦系数的关系Fig. 6 Relationship between shear band thickness and average particle size and friction coefficient

综上可知，剪切带的厚度受平均粒径及颗粒摩擦系数的影响较为显著，且当颗粒摩擦系数恒定时，剪切带的厚度约为（12～21）d50。

3 剪切带形成及演化细观特征

3.1 力链网络演化特征

颗粒物质力学研究结果表明，土颗粒通常以密集排布为主，相邻颗粒间接触形成许多传递外荷载的准直线型链状结构，称为力链。孙其诚等[14]认为力链是颗粒接触力传递的路径，其中传递较大力、与该准直线路径对应的若干颗粒组成的链状结构，称为强力链；弱力链则传递较小力；并提出了强力链FC的大小判据（FC大于平均接触力F）和角度判据，Z为平均配位数）。而力链网络是指土体颗粒在外荷载作用下相互挤压形成的网状结构，即颗粒接触网络；力链网络是颗粒骨架体系在细观层面上的表现形式，其复杂动力学响应机制决定颗粒体系的宏观力学特性，是颗粒物质力学领域重点研究的问题。

直剪试验中剪切带内外颗粒试样的变形或破坏均反映在力链尺度上，图7为砂土试样在100 kPa法向力作用下不同剪切时刻对应的力链网络结构变化。试样在初始固结稳定后，强力链较均匀且主要以环状形式沿盒体竖直方向分布；剪切进行后，随着对左侧墙体施加一定的速度，力链网络的结构开始发生变化，即开始向剪切带中部聚集；从剪切开始到峰值应力阶段，力链分布主要集中于剪切盒左上角和右下角区域，且呈线状分布。峰值应力状态时，力链网络的变化更加集中，主要表现为盒体左上角和右下角区域的力链向剪切盒中部贯通，且沿剪切盒呈对角线分布，此时强力链聚集并起骨架作用；峰值应力阶段后，随着剪切位移的增加，力链网络的分布继续变化；当剪切到D阶段后，强力链密集程度逐渐变疏散，且沿剪切带由初始的贯通状态向周围渐变扩散。这进一步体现出强力链的演化与剪切带的形成及受力有关。

3.1.1 力链光密度 为研究剪切带内外区域的强力链分布差异，运用专业图像分析软件Image pro-plus对力链演化图进行图像处理与分析，主要根据图中强弱力链的颜色深浅进行识别，提取出各阶段的强力链网状线条图像（图8），进而统计出力链的光密度（也称光强度）分布。光密度是生物细胞图像处理中常用的指标，光密度的定义为吸收光的物质的光学密度，其大小直接与染色物质的量成正比。平均光密度，其中：为 累积光密度；为被测区域面积。

图 7 剪切带力链演化Fig. 7 Force chain evolution

图 8 强力链图像提取Fig. 8 Strong chain image extraction

图 9 峰值时刻力链光密度Fig. 9 Optical density diagram of force chain at peak time

随后基于该软件自带的测量功能对强力链进行数据统计与分析，以剪切面为水平面，力链与水平剪切面的夹角为力链角度，提取出峰值时刻强力链的光密度数据。进而以力链角度为纵坐标，力链光密度为横坐标建立直角坐标系，得到峰值时刻强力链的光密度统计结果（图9）。

在图9（a）中，统计数据点呈纵向条带状分布。颜色较深的数据点代表弱力链，其灰度值集中于90～105，且角度分布范围介于0°～180°，这表明弱力链在水平面内均存在；而颜色较浅的数据点则代表强力链，其灰度值多集中于220～250，由图9可发现强力链的密度远远高于弱力链，约为弱力链的2～3倍，且图中部分数据点集中于90°附近呈水平直线状分布，这与峰值时刻部分强力链沿剪切盒呈竖向分布有关；而颜色介于黑色与灰色之间的数据点则代表平均光密度，图中平均光密度的灰度值集中于125～150。

图9（b）中，黑色数据点代表弱力链，红色数据点代表强力链，蓝色数据点则为强弱力链的平均值；3组数据点在0°～180°范围内均有分布，这表明强弱力链在剪切盒平面内均有分布。从光密度数据点的分布差异可发现，红色数据点的密集程度较高且多集中在1.7左右，表示强力链的分布密度比较集中；而黑色数据点较为稀疏且多集中在0.8～1.2区间内，表示弱力链的分布较松散；这均符合峰值时刻的力链网络分布情况。

3.1.2 强力链面积占比 同力链光密度处理方法一样，运用Image pro-plus软件对各阶段的强力链进行网状线条提取，可得到强力链图像。随后通过软件自带的测量功能进行强力链的面积统计及计算。文中主要对初始剪切A点、峰值时刻C点及剪切F点强力链的数据进行多次随机统计，并基于多次测量数据算出平均值，最终得到散点统计图（图10）。其中强力链面积占比，其中：为强力链面积；A1为被测区域面积。

由图10可看出，初始剪切A点的强力链面积占比为39%，主要是初始剪切时刻颗粒试样受法向荷载的作用处于固结状态，强力链多为竖向均匀分布。而峰值时刻C点的强力链面积占比约为52%，大于初始时刻，主要是由于剪切试验开始后剪应力的作用加强了颗粒的受力，从而使得盒体内的强力链增多；而剪切F点的强力链面积占比最小，约为19%，主要是剪切试验接近结束时，法向荷载的约束作用越来越小，强力链的分布达到最小值。

图 10 各阶段强力链面积占比Fig. 10 Strength chain area ratios at each stage

3.2 颗粒转动情况

颗粒转动指标是表征直剪过程中剪切带形成及演化的重要细观物理量。Oda等[15]在双轴试验研究中发现颗粒试样的体积变化主要反映在剪切带区域，且这种体积变化与颗粒体系的旋转紧密相关。文中基于数值试验对剪切盒内的颗粒转动量进行了监测，得到如图11所示的颗粒转动演化情况。

图 11 剪切带颗粒转动演化Fig. 11 Rotation diagram of shear band particles

由图11中的颜色分区可发现，颗粒试样的转动存在显著的差异性。首先体现在剪切初始阶段，试样的转动主要从剪切盒上下两侧开始变化，且有向剪切盒中部区域靠拢的趋势；其次是当剪应力达到峰值状态C点时，颗粒的转动发生剧烈的变化，转动区域明显缩小并向剪切盒中部聚集，如图11（c）所示；随后这种趋势更加强烈，最后在剪应力及法向荷载的作用下，颗粒转动范围逐步缩小并集中于剪切带区域，此时可测得剪切带厚度约（13～19）d50，满足前文中关于剪切带厚度的研究结果。

基于以上对颗粒试样转动过程的数值模拟，试验过程中对剪切带内外区域的转动角度进行统计，可得到如图12所示的关系曲线。

可发现，颗粒试样的转动量具有以下规律：（1）剪切前颗粒试样不发生旋转（为了便于观察，图12中以沿纵坐标的虚线表示）；（2）剪切带内颗粒的转动变化量最大，且沿剪切带上下两侧逐渐递减；（3）剪切试验中颗粒的转动是个持续不断的过程，当剪切达到峰值状态后，颗粒的转动仍继续发展，直至试验结束。

图 12 各阶段颗粒转动量统计Fig. 12 Statistics of particle rotation at each stage

3.3 配位数

颗粒试样的转动必将引起配位数的变化。配位数是指颗粒接触的平均数目，其大小反映了颗粒接触的密实程度，配位数越大，表明颗粒间的接触越充分[16]，更有利于骨架体系承担外力。密砂剪胀现象的实质是颗粒体系在外荷载作用下发生重排列引起的，试验中由于剪切力的作用，颗粒会发生滑移、滚动及重组，这直接影响配位数的变化。为了监测剪切带内外区域的配位数变化情况，将试样按横向进行分区，从下至上依次标记为1，2，3区，可得如图13所示的配位数变化曲线。

由图13可见，剪切过程中不同区域内的配位数变化存在明显差异；整体而言，剪切后各区域的配位数均有不同程度的降低：其中剪切带区域的降低幅度最大，当剪切到1 mm左右后变化减缓，且呈波状曲线递进，当剪切位移约3.3 mm时其配位数达到最小值；而位于剪切带上下两侧的配位数变化幅度相对较小；剪切带下侧区域内呈曲折线降低；而剪切带上侧部分颗粒位于盒体顶部，受挤压作用最小，故在剪切初始阶段，其配位数相对另外两个区域较小，在剪切到1 mm后其变化呈曲折线递进。综上可知，剪切带内试样的配位数减小幅度相对带外两侧更大。

图 13 配位数变化曲线Fig. 13 Coordination curve

图 14 孔隙率变化曲线Fig. 14 Porosity curve

3.4 孔隙率演化

同配位数一样，孔隙率反映了颗粒接触的密实程度，其大小与试样受力密切相关。为监测试验中剪切带内外区域的孔隙率变化情况，剪切开始前对试样进行横向分区，从下至上依次标记为1，2，3区，可得如图14的关系曲线。

由图14可发现，剪切试验开始后，各区域的孔隙率变化表现出明显的差异，其中剪切带区域内的孔隙率变化幅度最大，且曲线呈逐渐增长的趋势，直至试验结束；而剪切带上下两侧区域内颗粒的孔隙率变化幅度相对较小，其中剪切带下侧区域内的孔隙率变化呈折线增长，其增长幅度仅次于剪切带；而剪切带上侧区域则呈波状曲线平稳演化，且其孔隙变化率为负值。故密砂在直剪试验中，颗粒的孔隙率增长从剪切带向上下两侧逐渐发展，且在剪切带区域发生明显的剪胀现象，而剪切带外两侧则发生一定的剪胀和剪缩；总体而言，剪切带区域对试验整体的影响最大。

4 结 语

本文基于离散元数值方法，对砂土介质（以密砂为例）的直剪试验进行了模拟，从细观角度分析了试样剪切带的形成及演化机制，得到以下结论：

（1）离散元数值模拟方法模拟砂土直剪试验的结果与室内物理试验结果较为吻合，可用于研究粒状土体的力学特性。

（2）剪切过程中试样的应力、应变发展具有明显的峰值及峰后软化现象，这导致试样的变形模式从均匀变形向局部变形转变，最终在剪切面区域内发生塑性变形而形成狭长的剪切带，剪切带的厚度受平均粒径及颗粒摩擦系数的直接影响，且当颗粒摩擦系数恒定时剪切带的厚度约为（12～21）d50，反映出剪切带对外荷载的响应较为剧烈。

（3）颗粒体系的力链网络结构随剪切过程不断变化。通过对剪切各阶段的强力链进行网状线条提取分析，发现初始阶段弱力链集中且占主导地位，随着剪切的进行，强弱力链分布发生调整，强力链逐渐增多并集中于剪切带区域，且峰值时刻的强力链最为密集，说明试验中力链形态的分布与剪切带形成及演化具有对应性。

（4）颗粒试样的平均转动量、配位数及孔隙率的变化具有相似性，三者变化幅度中较大的范围均集中于剪切带区域，这表明剪切带是直剪试验的核心区域，对颗粒试样的细观参数起着决定性作用。