基于颗粒堆积算法的堆石料压实密度预测研究

沈超敏，邓 刚，刘斯宏，严 俊，毛航宇，王柳江

(1.河海大学 水利水电学院，江苏 南京 210098；2.中国水利水电科学研究院 水利部水工程建设与安全重点实验室，北京 100048)

1 研究背景

堆石料是由不同粒径的岩石颗粒相互填充形成的颗粒材料，因为其优良的工程特性，在土石坝、堤防、基础等工程中广泛应用。与砂土、粉土等材料定义了明确的粒径区间不同，堆石料的粒径跨度较大，属于典型的宽粒径颗粒材料。以土石坝工程为例，筑坝堆石料的最大粒径可达1000 mm，而粒径较小的颗粒可以低于0.075 mm[1-4]。对于粒径跨度大的颗粒集合体，存在以下两个问题值得探讨：(1)堆石料爆破过程的不确定性和施工填筑过程对粒径的离析作用难免造成堆石料在实际工程中级配的变异性。有研究表明，级配的变化会影响堆石料的变形模量、峰值强度和剪胀特性等[5-7]，对工程的变形与稳定造成不确定性隐患。随着我国在不良地质条件下筑坝、软岩筑坝等工程的增多，堆石料级配的变异性问题进一步凸显。(2)堆石料的缩尺效应。由于工程现场的堆石料最大粒径往远超出室内土力学试验允许的范畴，因此需要将堆石料的级配进行缩尺后开展相关试验。如何建立缩尺后与原型级配堆石料之间工程特性的联系是广泛讨论的难题[8-10]。

堆石料的级配变异性与缩尺效应从本质上属于同一个问题，即如何建立堆石料级配(包括粒径)与其物理力学性质的定量联系。目前已有一些研究针对这个问题进行了一些重要的尝试，且普遍认为颗粒材料的压实密度(最大干密度)是级配与物理力学性质(临界状态位置[11]、压缩模量[12-13]等)的桥梁，因此，有必要在此之前建立级配与堆石料最大干密度之间的定量联系。

目前，已经有部分研究通过室内振动压实和现场大型密度坑试验研究了初始级配对堆石料最大干密度的影响规律，并基于试验结果提出了最大干密度的经验预测方法。例如，朱俊高等[14]通过室内干密度试验研究了堆石料的级配对密度的影响，提出最大干密度与土料级配参数以及缩尺比例间的关系；吴二鲁等[15]采用上述连续级配函数，建立了级配函数参数与最大干密度的统一关系；朱晟等[10]通过大型现场密度试验，研究了不同的初始级配和缩尺方法对堆石料干密度的影响，并评估了缩尺效应对模量的高估程度。对于级配不连续的堆石料也有部分学者开展了研究，如于际都等[16]通过试验研究了间断级配堆石料的压实特性，并基于试验结果建立了间断级配堆石料压实干密度的预测模型；文喜南等[17]在已有颗粒堆积模型的基础上假定了单一粒径的粗粒土颗粒的最优孔隙比与粒径呈幂函数关系，并采用参数优化算法反演模型参数，预测了不同工程粗粒土的最小孔隙比。尽管这些试验结果均具有重要的启示意义，但这些经验公式一般仅适用于某种特定的级配函数，研究结果的可拓展性有待进一步明确。

本文将提出一种可以快速预测给定级配堆石料密度的计算方法，并通过连续、间断级配的堆石料试验结果验证该算法的合理性。最后，基于该预测算法讨论堆石料的粒径截断效应以及缩尺对堆石料最大干密度的影响。

2 颗粒最大干密度快速预测方法

2.1 颗粒材料降维映射概念颗粒堆积的问题一直是数学和理论物理学界关心的问题。近年来，颗粒物理领域提出了降维映射理论[18]，用于解决球体的密实堆积问题。由于降维映射理论本身并非局限于球形颗粒，因此在水利土木工程中也可以将这个概念拓展到预测一般形态颗粒的最密实堆积。颗粒材料降维的概念的核心是指将符合某一分布P3D(d)的三维颗粒堆积映射成一维的线段。此处，三维的粒径分布P3D(d)的定义是指粒径位于d到d+dd的颗粒数目概率。为了方便理解，图1以球形颗粒为例给出了三维颗粒堆积示意图，假设有一条直线穿越该颗粒堆积体，则该直线被这些颗粒的外表面分割成若干条线段。其中位于各个颗粒内部的线段长度记为L，这些线段的长度分布函数P1D(L)可以与P3D(d)建立如下定量关系[18]：

图1 三维颗粒与一维线段的映射关系

(1)

式中dM为最大的颗粒粒径。

颗粒材料的降维概念建立了三维颗粒与一维线段之间的定量联系，因此，三维颗粒堆积形态的差异均能通过一维线段的排列来表征。如两个颗粒的距离改变对应着这两个颗粒“切割”的线段之间距离的改变。

然而，若要仅从一维的线段排列反映三维颗粒的堆积形态则需要给定这些一维线段间的间隙，换言之，需要给定这些线段之间的势函数。需要说明的是，堆积模拟只要得出堆积的结果，并不需要完全重现颗粒堆积的动力学过程。为此，势函数的具体形式只需要保证颗粒之间不发生重叠，且超过一定间距之后颗粒之间不再具有相互作用即可。为此，本文采用以下最简单的形式：

(2)

式中：V为两条线段间的势能，由于三维的颗粒是硬球，因此V取值仅为零或者无穷大；h为任意两个线段Li和Lj的邻近端的距离；f为颗粒的自由体积参数。式(2)的势函数中可以通过参数f的取值确定一维线段之间的最小间距，当三维颗粒的形状从理想的球形变为不规则形状时，对应的一维线段之间的最小间距增大。Shen等[19]进一步通过模拟发现通过改变参数f的取值，可以反映不同形状系数的颗粒材料堆积。表1给出了一些文献中出现的典型颗粒材料对应的f取值[20]，针对堆石料的参数f取

表1 文献中典型颗粒材料的参数f取值

值还有待进一步明确。式(2)定义的颗粒之间的势函数规定了线段之间的最小间距仅与线段的长度和颗粒的形状有关，并没有考虑不同方向上的颗粒最小间距的差异。对于特殊的颗粒(如片状颗粒、非凸颗粒)的堆积，由于堆积体会存在显著的各向异性，因此式(2)的势函数需要改进后才能考虑各向异性显著的堆积体。

2.2 颗粒材料堆积算法图2给出了计算颗粒材料最小孔隙比emin的算法流程，可以分为三个步骤，每个步骤可以用一层映射表示。

图2 堆石料最大干密度预测算法流程

第一层映射F1将颗粒的级配曲线的函数M(d)映射成颗粒数目概率分布函数P3D(d)。不难根据级配的定义建立两者的联系[21]

(3)

式中：V(d)为粒径为d的颗粒体积；dm和dM分别为试样的最小和最大粒径。

第二层映射F2通过式(1)将粒径数目的概率分布函数P3D(d)映射为一维线段的长度分布函数P1D(L)。并根据概率分布函数P1D(L)生成符合该分布的线段集合L1，L2，…，Ln。本文生成符合概率密度P1D(L)的一维线段集合时采用了舍选算法，其优势是可以数值上生成符合任意概率密度分布的随机变量。为了使生成的线段数量足够具有代表性，已生成的线段集合L1，L2，…，Ln的频数统计与目标分布P1D(L)对比，直到每个长度区间内的差异小于1%，则停止生成新的线段。

第三层映射F3将符合P1D(L)分布的线段进行最优排列，线段的间距由式(2)的势函数确定。当线段集生成后，本文运用较高效的“贪婪算法”将P1D(L)映射到最小孔隙比emin空间。算法的主要优点是将线段按照长度逆序插入，从而避免计算新插入的线段与所有已存在线段的势能。该算法可以总结为：

(1)将线段按照逆序标记为L1≥L2≥L3≥…≥Ln，并将线段将从长到短依此插入线段集；(2)设定线段间的间隙集{gi}，间隙数量与线段的数目相同。随着第一条线段L1的插入，第一个间隙的长度为g1=fL1。为了插入线段Lj，首先找到并删除已有的线段集中线段间的最大间隙gmax，则新插入的线段会带来两段孔隙，即fLi和max[gmax-(1+f)Lj，fLj]。

按照上述的最小孔隙比算法，最终颗粒体系的最小孔隙比可以按照以下公式计算：

(4)

相应的堆石料最大干密度则可以根据最小孔隙比建立联系

(5)

式中Gs为堆石料的颗粒比重。

3 堆积算法的验证与应用

为了验证该堆积算法在预测不同级配堆石料最大干密度的合理性，本文选用了两组文献中的试验数据进行验证。

第一组用于验证的数据来自文献[16]的间断级配堆石料最大干密度试验，试验的级配通过在连续的基准级配基础上剔除某一粒径组的颗粒制备，其中连续的基准级配曲线由分形分布函数描述[22]

P=(d/dM)3-D

(6)

式中D为级配的分形维数，试验中分别考虑了D=1.9，2.1，2.3和2.5四种不同的基准级配，且针对每个基准级配分别考虑了5～<20 mm，2～<5 mm，0.5～<2 mm，0.25～<0.5 mm以及0.075～<0.25 mm粒径组缺失情况下的最大干密度。

模型预测前，采用某一级配的堆石料最大干密度试验数据标定参数f，本次标定选取的试样为D=2.3，缺失粒径组为0.075～<0.25 mm的试样，标定确定的参数f取值为1.13。标定完成后，该材料的其余试样均采用该参数进行计算预测。计算过程采用了常规的处理器i9-12900H，预测一个试样最大干密度的计算机耗时平均为0.76 s。相较于一些文献中的数值模拟方法计算堆积密度(如离散元等)，计算效率提升非常显著。

图3比较了上述间断级配堆石料最大干密度的试验数据与模型预测结果。从图中可以发现，即使在堆石料的级配曲线存在间断的情况下，本文提出的算法依然能够较好地预测给定级配堆石料的最大干密度。值得注意的是，在密度较低的区间内，图3中的预测和实测密度的偏差不可忽视。可能的原因是文献中的试验并未给出小于0.075 mm的粒径分布，因此在模型预测过程中并未考虑小于0.075 mm的颗粒影响。

图3 间断级配堆石料最大干密度试验和模型预测值对比

图4进一步分析了基础级配的分形维数以及间断位置对堆石料的最大干密度的影响。从图4(a)的试验结果看，无论是连续的级配还是缺失某粒径组的间断级配堆石料，堆石料的最大干密度随着分形维数D的增加而增大。图4(b)的模型预测结果中除了5～<20 mm间断级配中略有差异外，基本较好地预测了这一规律。此外，试验结果表明，对于不同基准级配的堆石料，粒径间断的位置对堆石料影响并没有一致性的规律，而图4(b)的模型预测却能够较好地重现上述试验结果，表明本文提出的堆石料最大干密度计算方法在预测间断级配堆石料堆积特性方面效果显著。

图4 级配分形和间断位置对堆石料最大干密度的影响

第二组验证试样采用了吴二鲁等[15]测试的连续级配砂砾石料进行密度预测研究。这些级配曲线统一由以下级配函数描述

(7)

式中：dM为最大粒径，分别为60、40和20 mm；b和m为函数中的参数，通过对b和m分别取不同的值可以得到不同的级配曲线：当dM=60 mm时，通过改变b和m的值，得到了16组级配曲线，可以研究级配曲线的形态对最大干密度的影响；当b和m取值一定时，改变dM的值，得到了8组级配曲线，可以研究缩尺效应对最大干密度的影响。

图5中对比了采用本文提出的预测方法得到的上述24组连续级配堆石料的最大干密度与文献中的实测结果。计算中参数f取值经标定为0.933。从图中可以看出，本文提出的最大干密度预测方法可以快速预测不同连续级配的堆石料最大干密度，且在预测缩尺试样的最大干密度方面同样效果较好。

图5 连续级配堆石料最大干密度试验和模型预测值对比

在验证了提出的预测方法适用于文献[15]中的砂砾石料后，本文进一步探索该算法在土石料级配优化方面的应用前景。图6为实测和基于本文算法预测得到的参数b和m对砂砾石料最大干密度的影响规律。图中的离散点为实测或预测数据点的位置。从图6(a)的试验结果可以看出，堆石料的最大干密度在b≈0.1，m≈0.6时出现极大值，且最大干密度的极值区域在图中呈“马鞍状”分布；远离该区域的最大干密度出现不同程度的减小。需要说明的是，云图中除了实测点之外的结果均是插值得到的，由于图6(a)中参数b在-1～0区间内的数据点较少，因此图中插值得到的结果可信度并不高。为了更准确地了解级配的变异性对该砂砾石料最大干密度的影响规律，图6(b)采用堆积算法进行了更加精细的预测，数据点数量相较于实测点的16个增加至36个。从预测的结果看，参数m和b对最大干密度的影响规律基本与实测结果一致，但在m∈(0.4，0.6)，b∈(-1，-0.2)区间内得到了实测数据缺失部分的密度分布规律：在b趋向于-0.6，m趋向于0.4值附近，密度存在极大值。

图6 参数b和m对连续级配砂砾石料最大干密度的影响

4 关于细料截断以及缩尺效应的讨论

4.1 截断误差的影响本文提出的堆积算法是基于已知的级配曲线预测堆石料的最小孔隙比，并换算成最大干密度。然而，实际的堆石料级配信息往往并不完整，小于筛分试验中最小筛径的颗粒粒径分布属于未知信息。这种粒径截断必然造成模型在预测过程中存在误差。显然小于最小筛径的细料占比越多，模型的预测的误差的可能性也会越大。此外，堆石料由原型级配到室内级配的缩尺的过程也伴随着改变截断粒径百分比的问题，有必要研究粒径截断在堆石料密度的缩尺过程中扮演的角色。

文献[15]中最小的筛分粒径均为0.5 mm，小于0.5 mm的颗粒粒径分布并未试验获得。因此采用本文提出的预测模型仅能考虑粗料(粒径大于等于0.5 mm)形成的最密实堆积，细料的影响并未考虑。因为文献[15]中的不同的试样来源于同一料源，且这部分细料没有进一步筛分，因此可以假定各试样中小于0.5 mm的颗粒粒径分布相同。

假定某堆石料密实试样，剔除了细料后粗料(试样Ⅰ)先进行最大干密度试验，得到的最小孔隙比为

(8)

式中Va1和Vs1分别为粗料密堆积体中孔隙与固相的体积。作为对比，不剔除细料的堆石料(试样Ⅱ)也进行了最大干密度试验。试样Ⅱ相比于试样Ⅰ，增加的细料含量造成的影响包含两个方面：一是对孔隙体积的影响，另一个是对固相体积的影响。试样Ⅱ的孔隙比写作

(9)

式中Va2和Vs2分别为试样Ⅱ相较于试样Ⅰ，增加细料造成的孔隙和固相的体积增量。

固相部分的体积Vs1和Vs2是固定不变的，可以根据级配曲线中截断粒径对应的百分比得知；对于孔隙部分的体积，Va1可以根据本文提出的堆积算法得到，而Va2的值并不固定，由于堆积的模式不同，当细料填充粗料的孔隙时，Va2的取值甚至可能为负值。

假定增加细料会同时改变孔隙和固相的体积，但不同试样中细料的孔隙与固相的体积比值为一恒定值，即

Va2/Vs2=e2=const

(10)

该假定的主要依据是试验中的砂砾石料来源于同一料源，且宽级配的堆石料试样中，细料部分的堆积与粗骨料的堆积独立性较高，因此小于最小筛径的颗粒级配虽然未知，却基本相同。

结合式(8)—(10)，不难建立emin与粗料最小孔隙比e′min之间的关系

emin-e′min=θe2

(11)

式中θ为截断粒径对应的细料含量百分比。

式(11)表明，若剔除材料中的细料进行密度预测，则预测密度的误差与材料中实际的细料含量成正比。图7给出了砂砾石料的最大干密度预测的误差与材料细料含量占比的关系。可以从图中看出，预测密度的误差与细料含量θ基本成正比，结果符合式(11)的预测结果。且可以根据图7中拟合直线的斜率确定e2的值，为e2=0.4042。

图7 细料含量与最小孔隙比预测误差的关系

将上述e2的取值代入式(11)，可以修正砂砾石料最大干密度(见图8)。可以看出，相较于不考虑细料含量的影响，修正后的预测结果整体上更合理。值得注意的是，图5中未修正的预测结果和修正之后的虽然有差异，但两者都能较好地预测不同级配砂砾石料的最大干密度。未修正的预测模型基本能够预测最大干密度的原因主要是：尽管试验中试样的细料含量有所差异，但总体上仍以粗料为主；其次，本试验中的细料源自天然的砂砾石土，具有一定的自相似性。然而，若细料含量较高的非连续级配试样，则细料的截断误差不可忽视。

图8 砂砾石料最大干密度试验和修正后模型预测值对比

4.2 密度尺寸效应堆石料制样密度的缩尺效应研究是研究堆石料力学行为无法回避的基本问题之一。堆石料的缩尺方法分为剔除法、等量替代法、相似法和混合法四种，每种缩尺方法均有特定的工程使用范畴[23]。因为不同工程的堆石料级配曲线的差异和选用缩尺方法的不同，很难从堆石料制样密度的缩尺效应试验数据中得出统一结论[24]。堆石料密度出现细料截断的典型场景之一是堆石料缩尺后制样[25]。以图9中的相似级配法为例，在室内制备堆石料的缩尺试样时，小于筛分粒径的颗粒不再进一步筛分。因此，缩尺前后的试样中，小于筛径的细料并未进行缩尺。

图9 缩尺制样过程中的粒径截断

显然，若试样中的细料部分级配未知，堆石料最大干密度的缩尺将无法准确预测。针对文献[15]中的缩尺试样，图10给出了试验与模型预测的砂砾石料最大干密度的缩尺效应。其中，图10(b)的模型预测中采用式(10)的修正方法考虑了细料截断的影响。考虑到原型和缩尺后的堆石料级配曲线存在差异，但是两者的堆石料颗粒是一致的，因此在模型预测时，保持室内和现场试验时参数f的取值一致，仅通过改变输入的级配曲线得到原型和缩尺后的最大干密度。从图中可见，堆石料的最大干密度均随着缩尺比例(原型尺寸/缩尺后尺寸)的增大而减小，且在明确细料截断影响的前提下，模型能够较好地预测给定级配堆石料的缩尺效应。

图10 试验与模型预测的砂砾石料最大干密度缩尺效应对比

然而，式(10)中参数e2在实际情况下是未知量。为了进一步明确细料、粗料的粒径分布在堆石料密度的缩尺效应中扮演的角色，本文采用堆积算法进行了堆石料密度缩尺效应的敏感性分析。首先采用了b=-0.2，m=0.4的粗料级配曲线，细料部分采用了0.49～0.5 mm之间的均匀分布细料(记为L0.49-0.5)、0.1～0.11 mm之间的均匀分布细料(记为L0.1-0.11)以及0.075～0.5 mm之间的宽级配均匀分布细料(记为L0.075-0.5)。从图11(a)可以看出，无论细料部分的级配如何变化，砂砾石料的最大干密度总是随着缩尺比例的增大而减小。从颗粒的堆积角度看，堆石料的缩尺造成了两个影响：一是粗骨架孔隙尺寸的降低，二是细料的含量增加。随着缩尺比例的增加，降低粗料的孔隙尺寸后，总体上不利于填充细料，导致缩尺后的试样密度降低；而增加细料含量对密度的影响却存在一定的不确定性。例如，对于粗料中L0.49-0.5的细料，最大粒径为40 mm的试样和60 mm的试样最大干密度就不存在显著差异。

图11 砂砾石料最大干密度缩尺效应的敏感性分析

同时，假定细料级配为L0.49-0.5，粗料分别采用b=0.3、m=1，b=0.6、m=0.8，b=-1、m=0.6，b=-0.2、m=0.4四种不同的级配。图11(b)给出了上述四种堆积条件下的粗粒料缩尺规律，可以看出，不同粗料级配下的砂砾石料缩尺规律各不相同。对于b=-1、m=0.6和b=0.3、m=1的粗料级配，缩尺后的试样最大干密度并无显著变化，而b=0.6、m=0.8和b=-0.2、m=0.4的试样则其密度随着堆石料粒径的增大呈增大趋势。

从图11的模拟结果可以看出，堆石料最大干密度的缩尺与粗料和细料部分的级配均有关。缩尺以后的砂砾石料最大干密度以降低为主，但并没有统一的规律。这一模拟结果也印证了不同文献中关于堆石料密度的缩尺规律不一致的试验结果。因此，为了更好地预测原型尺寸堆石料的最大干密度，一个可行的方法是在筛分细料部分的粒径分布，并采用颗粒堆积算法预测原型尺寸的堆石料最大干密度。

5 结论

堆石料的最大干密度是反映其工程力学特性的重要指标。本文基于颗粒材料的降维映射原理，提出了能够实现给定级配堆石料最大干密度的快速预测算法。进一步将该算法运用于堆石料最大干密度的变异性研究，并探讨了细料截断以及缩尺效应对堆石料最大干密度的影响。得到的主要结论有：

(1)本文提出的颗粒堆积算法仅包含一个计算参数f，通过一个堆石料试样标定完参数后，即可以实时计算预测该堆石料在其他给定级配情况下的最大干密度。与文献中连续和间断级配堆石料的试验结果对比发现，颗粒堆积算法不仅具有较好的预测精度，也可以进行比试验更高效的级配优化分析。

(2)由于堆石料细料部分级配属于未知信息，预测堆石料最大干密度会造成不可避免的截断误差，误差的绝对值与该堆石料的细料含量呈正比。对粗料和细料的堆积特性进行分析进一步发现，由于堆石料形成过程造成的自相似性，不同级配的堆石料细料部分的密堆积孔隙比基本一致。

(3)堆积算法可以较好地根据缩尺后试样的最大干密度预测原型级配试样的最大干密度。且可以从堆积角度给出最大干密度缩尺效应的内涵：从颗粒堆积角度看，堆石料的缩尺事实上是粗料部分的缩尺，细料部分仅发生含量的变化。粗料部分缩尺后骨架的孔隙尺寸降低，不利于细料填充，进而导致缩尺后的试样密度降低；而增加细料含量对密度的影响却存在一定的不确定性。采用堆积算法开展的敏感性分析也表明，堆石料最大干密度的缩尺效应并没有统一的规律，与粗料和细料具体的粒径分布都密切相关。

本文提出的堆石料密实堆积算法在预测堆石料的最大干密度中呈现出较好的效果和广泛的应用前景，该模型后续有望用于预测现场堆石料的压实程度、分析堆石料密度的缩尺效应以及开展堆石料级配的优化设计。不过在运用该模型时需要明确该模型的使用范畴：(1)本模型对最大干密度的预测是基于当前级配预测当前级配下的最大干密度，并不考虑堆积过程中的颗粒破碎影响；(2)模型主要针对总体呈现各向同性的颗粒堆积，不能适用于片状、针状等特殊的形态颗粒的堆积。作为一个尚处于初期的计算模型，在以下方面还有完善的空间：(1)本计算模型是基于堆石料的当前级配预测堆石料在当前级配下能达到的最密实堆积，对于软岩等破碎显著的堆石料，后续研究中有必要将颗粒破碎演化方程与堆积算法统一，提出考虑颗粒破碎的堆石料最大干密度预测算法；(2)堆积算法中的参数f与颗粒的总体形状有关，后续研究中有必要进一步建立f与颗粒形态参数之间的定量关系。