压实作业是沥青路面施工的关键环节，研究表明国内沥青路面的早期病害和压实直接相关.数值仿真是揭示沥青混合料压实机理的重要方法，文献[1]通过室内试验测试了沥青混合料流变性能，利用有限元法建立了沥青路面压实模型.文献[2-3]利用MTS(material testing system)研究热态沥青混合料室内压实过程的变形特性，并采用有限元软件建立二维沥青路面的压实模型.文献[4-6]采用旋转压实确定了压实参数，建立热力学沥青混合料压实模型，利用现场压实试验证明模型的合理性，研究结果表明，该有限元模型可预测压实过程材料总体宏观的变化趋势.基于有限元建立的路面压实模型，是将三相非均匀的沥青混合料假定为连续均质材料体系，压实过程集料发生较大的变形假定为小变形.因此，采用该方法分析沥青混合料的压实过程存在一定的不足，不能有效地分析材料内部体系的细观结构行为和深层次揭示沥青路面的压实机理.

离散元法能反映非均匀材料的组成特性，擅长模拟非均质材料和大变形过程材料的细观行为. Liu等^[7-10]一直从事该领域的研究，如不规则集料的生成^[7]、集料特性对力学性能的影响^[8]、车-路作用机理^[9-10]等.有关沥青混合料细观研究集中于压实成型后的宏观力学性能与细观结构的关系，对沥青路面的压实尤其是现场压实过程研究较少.文献[11-12]建立了二维压实模型且将集料颗粒视为理想的圆盘单元，忽略了集料的形态特征，与真实集料有较大的差异.本文利用DEM建立考虑集料形态特征的沥青路面三维压实模型，从集料运动位移与角度、接触力、能量等方面研究压实过程中材料的细观行为.

1 试验 1.1 试件设计

以SMA-13为研究对象，沥青胶浆是由细集料(小于2.36 mm的集料)、矿粉和沥青胶结料组成，利用比表面积法确定沥青含量，其级配见表 1.本文采用SBS改性沥青，密度为1.036 g/cm³，SMA-13的沥青含量(质量分数)为6%，根据比表面积法计算得出其沥青胶浆的沥青用量为26.4%，由于沥青胶浆中沥青含量较高且集料级配偏细，其空隙率可认为是零.

采用静压成型圆柱体(R=50 mm，h=150 mm)沥青胶浆试件.该试件不含粗集料且沥青含量高，在高温时变形量大，成型过程中需控制好温度.在制备时试模内壁均匀抹好润滑油，避免因内壁与试件粘结导致试件端部成凹形.在实际压实过程中温度变化使得本构关系及其参数发生改变，鉴于问题的复杂性，压实过程的温度保持恒定.文献[11, 13]利用细观Burger’s模型表征沥青混合料的细观接触模型，并取得了满意的效果.从理论上说，宏观Burger’s模型并不反映材料的塑性变形行为，然而，要更加准确描绘热态沥青混合料的现场压实细观行为必须首先开发宏观的黏弹塑本构模型，其次，通过C++编程开发细观的黏弹塑接触模型，这对编程能力的要求极高；同时细观Burger’s接触模型的参数高达8个，若建立细观黏弹塑接触模型的参数会更多，如何确定材料的参数是一个难题^[14].为此，本文将利用Burger’s近似表示沥青路面压实过程中的黏弹塑行为.

1.2 动态模量测定

为了获得材料的动态模量和确定Burger’s模型的细观参数，利用UTM(universal testing machine)测试不同温度(5、20、30 ℃)和不同频率下(0.01、0.1、0.2、0.5、1、2、5、10、20、25 Hz)的沥青胶浆试件，动态模量测试结果如图 1所示.压实过程中的温度较高，参数难以凭借简单的试验手段获取.通过既有的参考温度根据温度-时间等效原理建立动态模量主曲线，从而获得其他温度和频率下的材料参数^[15].考虑到Burger’s模型在数值仿真过程运行的效率低下，基于温度-频率等效原理的移位因子能显著提升运行效率并具有良好的稳定性和可靠性^[15-16].通过反复计算，最终选取移位因子为1.0×10⁷.利用回归分析获得了Burger’s模型的参数分别为E₁=14.256 MPa，E₂=10.246 MPa，η₁=629.387 MPa·s，η₂=1.817 MPa·s.

2 沥青路面压实模型

PFC(particle flow code)具有比较高的计算效率和方便用户编程，从而成为离散元法重要的分析软件之一.轮胎压路机与压实路面的接触界面复杂，且对混合料施加荷载存在竖直和水平方向作用力.

因而本文仅建立钢轮压路机的压实模型.压路机整个系统尤为复杂，难以利用DEM建立仿真模型.同时，当前利用PFC研究材料的细观行为由于计算效率的限制，一般通过缩小尺寸或将有限元与离散元进行耦合，对关键部位采用离散元建模，而其他部位采用有限元模拟.从计算效率和研究重点考虑，本文并未建立整个路面结构，仅通过等效原则建立0.15 m×0.05 m×0.05 m的微型压实模型.

2.1 虚拟沥青路面生成

为了简化虚拟试件，根据现有研究成果^[7-9]，沥青混合料视为由粗集料(大于2.36 mm的集料)、沥青胶浆(小于2.36 mm的集料、矿粉及沥青胶结料)和空隙组成的三相体系.粗集料通过不规则的多面体模拟，详细算法参考前期的研究成果^[17-18].沥青胶浆当作连续颗粒体系；沥青混合料内部的空隙分布规律复杂，与压实成型过程紧密相关，在不研究空隙率的特性时，已有的研究成果将空隙的分布通过随机删除一定数量的沥青胶浆实现^{[14, 16]}.

2.2 接触模型与参数

在PFC数值仿真中，不同的接触模型表征材料内部不同组分的细观行为.根据沥青混合料的组成特点，沥青胶浆内部、沥青胶浆与粗集料间均为黏弹塑接触模型，粗集料间为线性刚度模型和滑动模型，限于篇幅，接触模型的基本原理可参考文献[19].集料参数确定比较复杂，不易验证，以往的研究通过经验选取.本文集料的弹性模量为55.5 GPa，泊松比为0.25，集料间的摩擦系数为0.5.沥青胶浆的参数已通过试验确定，根据宏观、细观参数之间的关系，并通过不断地试算最终可得到模型参数^{[15, 19]}.

2.3 压实作用力

钢轮振动压实过程对沥青路面的作用力视为机械自重G和激振力F₀的分力F的合力，即

$ P = G + F, $

(1)

式中F=F₀sin ωt，ω振动频率，ωt=π/2+2nπ时，F最大, n=0, 1, 2.

钢轮压路机与压实沥青路面的接触特性复杂，与压实温度、压实机械参数及初始压实度等诸多因素相关.随着压实度不断提高，两者接触面积逐渐缩小.为简化压实模型，仅研究单个钢轮的压实作用，将压路机和压实路面接触形状简化为0.05 m×0.03 m的矩形，压路机对路面的压实作用视为均布荷载, 有

$ p = P/(WB). $

(2)

式中：p、P、W分别为压路机的压实应力、作用力和钢轮宽度；B为接地尺寸，是压实轮直径D和轮阻角ζ的函数，B=0.5D/sin ζ.

根据文献[20]对钢轮振动压路机压实应力的研究，压路机对压实沥青路面的作用力为

$ p = 13\left( {G + {F_0}} \right)/\left( {DW} \right). $

(3)

根据PFC3D的基本原理，墙单元无法直接施加荷载，仅通过合适的速度导致离散元模型内部颗粒距离减小，并使得颗粒间发生接触，从而产生反作用力.离散元颗粒单元可以精确控制荷载作用，因此本文采用该方式加载.由于加载的颗粒单元数较多，为了使这些颗粒不发生分离从而影响加载过程，需要对其赋予较大的粘结力.

压实模型中离散元颗粒单元施加荷载的过程如下：1)通过式(3)计算不同压实工况时压路机施加的应力. 2)压实区域为0.05 m×0.03 m，施加荷载的颗粒单元直径为2 mm，单个面积约为4 mm²，共375个颗粒，如图 2所示；基于应力等效原理，根据上述过程所计算的应力，可得到单个离散元颗粒单元所要施加的力.

2.4 碾压时间与边界条件

在沥青混合料的压实过程中，沥青的性质与黏土具有部分类似的特性.部分道路工作者尝试利用固结模型CFF(consolidation fluid flow)揭示热态沥青混合料的压实机理.根据该模型可知：一定的压实功与足够的持续时间是保证有效压实的必备条件；在压实应力一定时，压实遍数的增加，压实功作用于被压实材料的时间越长，即通过控制压路机的碾压时间与压实遍数是等效的.在该压实模型中，通过控制压实时间实现压路机的碾压遍数.当碾压速率为v_d，钢轮压路机和沥青路面的接触长度为l_d，压路机完成一次压实作业的时间为

$ {t_{\rm{d}}} = {l_{\rm{d}}}/{v_{\rm{d}}}, $

(4)

由现行《公路沥青路面施工技术规范》^[21]查得压路机的压实作业的行驶速率为2~6 km/h.当不考虑压实施工过程中的停留时间时，按照式(4)计算得出完成一次碾压的时间为0.018~0.054 s.

由于路面路缘石或被压实材料周边材料的挤压，沥青混合料的现场压实边界存在一定的约束.为了合理地模拟沥青路面三维压实模型的边界条件，采用刚度较小的墙单元布置在离散元模型的四周，通过赋予刚度较大的墙单元仿真底部的约束.同时，鉴于压实模型不考虑钢轮压路机产生的水平作用力，需限制加载颗粒单元在水平方向的运动趋势.

3 沥青路面压实特性分析

沥青混合料中的集料运动与初始压实度(摊铺机作业后的压实度)、路面厚度、压实设备参数和材料级配等因素有关.研究表明沥青混合料摊铺后的初始压实度在90%左右，本文SMA-13的初始空隙率假定为14.5%.限于篇幅，本文仅探讨钢轮压路机静压过程，其行驶速度为6 km/h，压实轮宽和直径分别为2.0、1.25 m，荷载为77 kN.根据上述加载颗粒应力计算步骤，可得单个颗粒施加力为1.6 N，从路面厚度、接触力变化以及能量演化机制等方面研究路面细观压实行为.

3.1 路面厚度及集料位移

为了定量分析路面压实厚度演化，选取了有代表性的若干施加荷载的颗粒为研究对象. 图 3为路面厚度的变化情况，路面厚度压缩值随着碾压遍数的增大而增加，而厚度变化差值不断缩小.在初压阶段混合料内部空隙率大，大部分集料间存在一定的空隙，易于压实.然而，随着材料的不断密实，混合料空隙率减小，集料间逐渐形成骨架结构，具有较强的抵抗外界荷载的能力而难以压实，故路面厚度变化差降低.

沥青混合料属于非均质多相复合材料，压实过程中内部集料运动情况存在差异性.通过PFC3D内置语言fish遍历整个模型的集料颗粒，提取x、y、z坐标的累计位移，经编程后处理形成的位移云图，如图 4所示，括号内数值表示在此范围内集料的个数.颗粒单元位移在不同方向呈现非连续、不对称分布. x方向：集料颗粒位移主要为负值，且集中分布在-1.2~0.96 mm，颗粒位移在压实区域无明显的规律. y方向：颗粒位移总体上近似以模型中心线为对称轴呈对称分布，部分位于压实区域附近的集料向上运动，该方向颗粒的位移主要集中在-0.32~0.56 mm. z方向：集料位移集中在-3.5~0.62 mm，荷载区域的集料(占模型集料总比例的21%)总体呈现“V”字形，且颗粒在非荷载区域的位移较小，并表现出非均匀分布的特性.位于荷载区域的右侧，大量集料(占模型集料总比例的25%)存在与荷载方向相反的运动趋势，产生该现象的主要原因是该模型的底部以“墙”作为边界，荷载通过颗粒体系传递至边界时，颗粒无法穿越“墙”，约束边界对集料施加反作用力，且虚拟路面表面无约束，因此颗粒朝上运动，伴随着局部隆起现象.这与实际路面施工时在压实位置周边出现少量的混合料隆起相吻合，也与文献[22]通过数值模拟反映集料的隆起现象一致.根据颗粒位移的分布情况，发现竖直方向颗粒位移大于水平方向，表明颗粒位移受荷载方向的影响.

为了深入分析路面细观压实特性，通过各坐标方向的位移分量可得到集料总位移.基于PFC3D的可视化技术，通过虚拟分割获得集料在竖直与水平方向的位移矢量图.集料方向以箭头示意，位移大小以线长示意，如图 5所示.从位移方面分析，水平位移小于竖直位移，进一步说明集料位移与压实荷载方向正相关.在xz竖直平面：压实区域内集料运动向下，左右侧边界的约束导致集料向上运动.因此，在过渡区形成了类似的“涡流”结构.在xy水平平面：一些集料从加载位置朝两侧发生复杂的运动，也形成“涡流”结构.文献[23]研究了振动压实时沥青混合料内部存在较多“涡”结构，且“涡”结构发生从大到小的演化，本文模拟沥青混合料的压实过程，集料的运动情况不如振动压实激烈，但是产生的现象基本一致.

3.2 接触力演化

在压实应力作用下，混合料内部的空隙因集料-集料、集料-沥青胶浆及沥青胶浆-沥青胶浆的距离的减小而逐渐降低，且材料体系逐渐形成较稳定的骨架结构.从细观结构分析，混合料体系形成稳定的空间结构的本质是材料各组分产生一定的接触力.通过提取不同压实遍数时材料内部的接触力，如图 6所示，由图可知，线段宽度和力大小成正比例，红、黑分别代表接触拉力和压力.从定性角度分析，在不同压实遍数情况下，压实区域接触力的分布密度及其数值大小比非压实区域大，且混合料体系形成的接触力以压为主，以拉为辅.随着压实遍数的增加，混合料内部的接触力不断增大，并逐渐形成较稳定的骨架结构，导致集料的运动位移与角度的增量减小.因此，从不同压实遍数时接触力演化可以合理地解释集料的运动变得微弱.

为了定量分析压实过程中热态沥青混合料接触力的分布规律，通过fish语言编程提取集料-集料的接触力，如图 7所示.

压实过程集料之间的接触力以接触压力为主，压实初期接触压力在(0, 5 N)时的比例超过80%，在第4遍压实时，尽管接触压力的主要分布范围(73%)与第1遍相同，可是在其他区间的接触力比例大于第1遍.根据数据统计分析结果，第1、4遍压实接触力的平均值分别为2.6、4.3 N；中位数分别为0.8、2.1 N.

PFC3D通过测量球的应力平均化过程计算应力张量.为了不影响计算结果，球的直径宜大于材料最大公称粒径的2倍，且球心与模型的边界尺寸宜大于球半径.如图 8所示，球直径30 mm，球心距模型上下边界25 mm，m1、m2、m3分别表示球1、球2、球3，x、y、z方向的应力分别为S₁₁、S₂₂、S₃₃.

总体而言，压实区域m1在x、y、z方向的应力大于m2、m3的应力，且在z方向的应力平均值高达1.1×10⁵ Pa，明显大于x、y方向的应力平均值及其峰值，如图 9所示.由此表明压实区域的应力值与压实荷载的方向正相关.非压实区域m2、m3的应力主要分布在(2×10⁴ Pa, 4×10⁴ Pa)内，由于混合料内部的非均质性，在x、y和z方向的应力值存在波动性和变异性，但总体差异较小，表明非压实区域的应力值与压实荷载的方向无明显关联.上述分析结果与接触力的分布(图 6)情况基本一致，因而利用连续力学的应力平均化过程研究压实混合料的应力具有一定的科学性和合理性.

3.3 能量演化

压路机克服沥青混合料内部阻力做功并转化为其他能量形式，外力做功E_b、动能E_k和应变能E_c计算公式^[19]分别为

$ {E_{\rm{b}}} = \sum\limits_{{N_{\rm{p}}}} {\left( {\left( {m{g_i} + {F_i}} \right)\Delta {U_i} + {M_i}\Delta {\theta _i}} \right)} , $

(5)

$ {E_{\rm{k}}} = \frac{1}{2}\sum\limits_{{N_{\rm{p}}}} {\sum\limits_{i = 1}^6 {{{m'}_i}{v_i}{{^\prime }^2}, } } $

(6)

$ {E_{\rm{c}}} = \frac{1}{2}\sum\limits_{{N_{\rm{c}}}} {\left( {{{\left| {F_i^{\rm{n}}} \right|}^2}/{k^{\rm{n}}} + {{\left| {F_i^{\rm{s}}} \right|}^2}/{k^{\rm{s}}}} \right)} $

(7)

式中：m为质量，g_i为重力加速度，F_i为外力，ΔU_i为位移增量，M_i为外力矩，Δθ_i为旋转角增量，N_p为集料数，m_i^′为广义质量，v_i^′为广义速度，N_c接触数，F_iⁿ、kⁿ分别为法向接触力与刚度，F_i^s、k^s分别为切向接触力与刚度.

在压实初期，E_b、E_c的增量斜率大，随着压实的进行增量斜率不断减小；由于外界荷载刚施加作用力，部分颗粒的运动剧烈，E_k存在异常，然而后期逐渐稳定，如图 10所示.初始阶段，材料体系内部空隙大，在压实荷载的作用下集料运动位移增加，集料的距离减小，集料间的接触数增加.当压实进行到一定阶段，混合料内部集料间的接触数增加，并形成一定的骨架结构，因而集料的运动变弱，逐渐达到稳定值并趋向于平衡状态，压实模型内部能量转换趋于稳定.

4 结论

1) 利用基于离散元的PFC3D分析了SMA-13压实过程的细观特性，发现沥青混合料内部各组分的位移是不连续且不对称的，不同区域的集料运动规律不同，其中非压实区域部分混合料发生隆起的现象；集料的运动、应力受到压实荷载的位置和方向的影响.集料在压实荷载方向上的位移大于水平方向上的位移.

2) 在垂直方向上，压实区域中的颗粒向下移动，非压实区域颗粒向上移动；在水平方向上，大多数颗粒从加载位置向两侧移动和旋转；在压实与非压实过渡区形成“涡状”结构.

3) 刚进行压实时，外力做功和应变能增加较快，后期逐渐变缓；在压实初期，不稳定载荷导致颗粒速度和异常动能，进入后期阶段趋于稳定；外力做功、应变能及动能可较好地反映路面压实模型的状态.