新疆某工程粗粒料颗粒破碎大型三轴试验数值模拟

仝 帆，唐新军

(新疆农业大学 水利与土木工程学院， 新疆 乌鲁木齐 830052)

粗粒料在自然界分布广泛，储量丰富，并具有压实性好，透水性强，抗剪强度高，地震荷载作用下不易产生液化等优良工程特性，在工程中得到广泛应用[1]。通过人工爆破开采出来的粗粒料，颗粒比较脆弱，在较高的应力条件下容易破碎。随着我国筑坝技术的不断进步，涌现出一批200 m级的高坝，如已建成的233 m高的水布垭面板堆石坝，在建的261.5 m高的糯扎渡心墙堆石坝等，根据这些已建的高堆石坝变形观测结果，大坝的变形较大，其主要原因就在于粗粒土在高压下重新分布以及颗粒破碎所引起的孔隙压缩[2]。

当土体沿某一剪切面发生破坏时，必须要破坏在此之前互相咬合的颗粒状态，这种破坏一般有两种形式，一是使颗粒从原来的位置上抬起跨越相邻颗粒，二是将颗粒在尖角处剪断后移动，两种破坏形式都会表现为体积的胀大，在高应力下土体颗粒不断发生破碎，影响堆石料坝体的峰值强度、内摩擦角、孔隙比、渗透系数，削弱了其剪胀性、增加其流变性[3]，因此通过研究颗粒破碎来探索分析大坝在运行过程中的变形、沉降是非常有意义的。

国内外诸多学者对颗粒破碎开展过一定研究[3-12]，物理实验获取骨料破碎率的通常做法是，在室内三轴试验前后，应用筛分颗粒的方法来确定各粒组的质量分数差值，以此来表征颗粒破碎，由于室内实验的局限性，以往的实验多注重于颗粒破碎程度的影响分析，本文将应用PFC3D颗粒流数值模拟的方法，记录轴向应变各阶段下的颗粒破碎率，捕捉破碎颗粒绘制破碎颗粒位置云图，探究簇颗粒内异性受力黏结键变化，以期从细观角度揭示颗粒破碎的演变过程。

1 颗粒流模拟破碎颗粒的基本过程

1.1 数值模拟数据来源

本文数据来源于袁铁柱等[13]对新疆某工程筑坝材料的三轴试验研究，试样以片麻岩为主，几何形状极其不规则，为无黏性粗粒土。图1为实验后主应力差与轴向应变的关系曲线，在不同围压下应力-应变关系曲线都表现出应变软化特性，颗粒破碎率也随着围压的增加呈现一定变化趋势，颗粒破碎率见表1。

表1 室内试验破碎率

图1不同围压应力-应变关系曲线

1.2 破碎簇颗粒的建立

本文以不规则Clump颗粒为模板，按照等体积替换原则，以球形Ball颗粒替换Pebble颗粒，建立可破碎簇颗粒Cluster，在簇颗粒内部设置平行黏结键，并依据黏结键的断裂来产生破碎颗粒，可破碎簇颗粒形态见图2。

图2破碎簇颗粒形态

1.3 数值模拟过程

数值模拟过程可以分为以下几个阶段：

(1) 在高600 mm，直径300 mm的压力室内，控制最大粒径为60 mm，最小粒径为5 mm，生成不可破碎骨料颗粒Clump。

(2) 应用Multiply命令，使簇颗粒按照中心位置进行缩进，以此避免簇颗粒跑出加载室墙体，让簇颗粒体在加载室内进行循环调整，使得之前相互堆叠的簇颗粒在循环调整的过程中逐渐弹开，降低簇颗粒内部的最大不平衡力。

(3) 用刚性球体颗粒等体积替换簇颗粒中的每一个颗粒，形成可破碎的簇颗粒Cluster，删除原来的不可破碎颗粒(Clump)，将平行黏结键设置在簇颗粒内部。

(4) 通过对程序语言的调用，在加载室墙体上施加一个稳定的围压，保证数值试验的正常运行。

(5) 在加载板上施加竖向荷载，数值试验开始。

数值试验模型由压力室墙体、加载板、簇颗粒三大部分构成，其中，压力室墙体的刚度应小于加载板的刚度，以便模拟室内试验中橡皮桶的效果，按照上述步骤建模完成后，试样模型如图3所示。

图3数值模型试样图

1.4 接触模型的设置以及细观参数的选取

本文将平行黏结模型设置在簇颗粒内部，通过控制平行黏结模型的法向、切向强度来表征簇颗粒的坚硬程度。同时，将线性滚动模型设置在簇颗粒与簇颗粒之间的接触中，以此增加簇颗粒与簇颗粒之间的抗转动能力[14-16]，弥补了以往数值试验中通过设置较大颗粒表面的摩擦系数来获得更大滑动摩擦力的不足。

通过对试验样本分别赋予不同的细观参数进行数值试验，来获得试验样本的宏观力学特性，对比这些宏观参数与实际工程中的相应参数值，来确定合适的细观参数，这一过程也叫做细观参数的标定，通过大量的数值试验后发现，应用不同细观参数来模拟不同围压下的应力-应变试验较为合理，经过调整，最终选取的细观参数见表2。

表2 细观参数汇总表

2 数值模拟结果分析

2.1 数值试验与室内应力-应变试验结果对比

选用表2内细观参数进行模拟，得到不同围压下室内试验与数值试验应力-应变曲线对比如图4所示。

由图4可知，数值模拟试验与室内试验的应力-应变曲线具有基本一致的发展趋势，表现出应变软化特性，可以较为真实的反应室内物理试验。

2.2 颗粒破碎演化过程

每种围压作用下进行数值模拟得到的试样达到土体剪切破坏时的破碎率与室内试验得到的破碎率对比见表3，表中的断裂键数代表发生破坏的连接数目，破碎率指的是断裂的连接数目与总的连接数目的比值。由表中数据可以看出，数值模拟得到的破碎率与室内试验得到的破碎率误差较小，表2所提供的平行黏结强度可以真实反映实际材料宏观破碎特性的要求，选择参数是合理的。

表3 数值模拟与室内试验破碎率对比表

为了进一步阐述可破碎性土在荷载作用下的破碎演化规律，研究了在围压400 kPa、600 kPa、800 kPa下破碎率随轴向应变改变而发生的变化规律，见图5。

图4数值试验与室内试验结果对比

图5中，曲线为颗粒破碎率随轴向应变的变化曲线，点画线为标示线用以辅助说明颗粒破碎率随轴向应变的变化趋势，从图5中看出，破碎率随轴向应变的变化分为三个阶段：第一阶段a—b，簇颗粒所受轴向荷载较小，破碎率速度增长曲线平缓；第二阶段b—c，轴向荷载不断加大，应力集中现象明显，破碎率速度增长曲线较陡；第三阶段c—d，破碎率速度增长曲线逐渐平缓原因有两点，首先颗粒大量破碎，填充空隙，提高土体密实度抑制颗粒破碎。其次，颗粒破碎有优先级，大颗粒相比小颗粒更容易发生破碎，即在加载过程中，在同样的应力条件下，粒径较大的簇颗粒首先发生破碎，其次才是粒径较小的簇颗粒破碎，使得在应变软化阶段，达到临界状态后，样本破碎率增长缓慢。

图5破碎率和轴向应变关系

2.3 剪切位移场及追踪破碎颗粒变化规律

通过记录颗粒在剪切过程中颗粒的位移场，表达在加载阶段剪切带逐步形成的过程；追踪破碎颗粒位置，直观反应剪切过程中随轴向荷载的改变，颗粒破碎情况。图6揭示的是在不同围压下，同一轴向应变时刻，破碎颗粒位置以及簇颗粒位移场随轴向应变的变化。

加载桶内黑色短线表示破碎颗粒位置，对比不同围压在同一轴向应变时刻的簇颗粒位移场及破碎颗粒位置图发现，随着围压的增加，簇颗粒的定向位移趋势逐渐明显，颗粒破碎数量逐渐增多，随着围压的增加，较为明显的剪切带出现的越早并且剪切带宽度越窄。对比同一围压在不同轴向应变时刻的簇颗粒位移场和颗粒破碎发现，图6(a′)为加载初期颗粒位移场，颗粒位移方向，大小毫无规律，伴有少量的颗粒破碎，破碎颗粒多集中在加载桶壁处。图6(b′)为峰值强度之前的位移场，加载桶内颗粒开始定向移动，位移场呈轴对称分布，破碎颗粒数量逐渐增多，开始有在剪切带附近破碎颗粒集中的趋势。图6(c′)～图6(e′)为峰值强度出现后的位移场，该阶段颗粒定向移动趋势明显，破坏颗粒咬合平衡体系，破碎颗粒的数量也逐渐增多，对比颗粒位移场图6(c′)和图6(e′)可看出剪切带已经形成，大致呈45°倾角，这与Marsal R J等[17]的结论：堆石料在平面应变试验中大致沿45°角的平面发生破坏是一致的。

图6剪切位移场及破碎颗粒的形成过程

2.4 颗粒平行黏结键拉压平衡分析

数值模拟结果发现，未破碎颗粒簇内部受剪应力黏结键个数与受压应力黏结键个数随着轴应变的变化呈轴对称分布，不同围压下，黏结键个数变化如图7所示。

分析图7可知，随着轴向应变的变化，受压应力的黏结键个数逐渐减少，而受拉应力的黏结键个数逐渐增多，从簇颗粒运动趋势的细观角度分析其原因，来揭露簇颗粒破碎的内在原因。在加载初期，孔隙所构成“塌陷区”居多，受荷载作用，簇颗粒产生定向移动的趋势，运动到“塌陷区”上部的簇颗粒无法获得支撑，在轴向荷载作用下簇颗粒由受压应力转为拉应力，达到极限拉应力时，簇颗粒产生破坏。加载后期，大量破碎颗粒将孔隙逐步填充完毕，“塌陷区”数量明显减少，使簇颗粒破碎率增长缓慢。从簇颗粒在荷载作用下受力黏结键的拉压特性变化也可分析出，在荷载作用下颗粒破碎以受拉应力的剪切破坏为主。

图7不同围压下黏结键个数和轴向应变关系

3 结 论

(1) 峰值强度出现之前，颗粒位移场呈轴对称分布，少量颗粒发生破碎，峰值强度出现后，颗粒的定向位移趋势明显，剪切带附近的颗粒破碎量大，且颗粒破碎是由骨料内部向外部拓展的过程。

(2) 随着围压得增加，颗粒破碎率逐渐增加；对比同一围压下，不同轴向应变的颗粒破碎率变化发现，颗粒破碎率曲线随轴向应变的变化分为加载初期的缓慢增长期、应变硬化阶段快速增长期、以及应变软化阶段的平稳增长期三个阶段。

(3) 分析加载前后破碎颗粒内部黏结键特性变化，表明颗粒破碎以剪切破坏为主。