水工沥青混凝土三轴试验的三维细观模拟

刘 璇，唐新军，王建祥

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

水工沥青混凝土三轴试验的三维细观模拟

刘 璇，唐新军，王建祥

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

基于颗粒离散元理论，利用沥青混凝土三轴试验进行三维数值模拟，对水工沥青混凝土的细观性能研究具有重要意义。通过数值模拟发现：根据沥青混凝土的材料组成及特性，选择合适的材料单元间的接触本构模型以及合理的细观力学参数，使用颗粒流PFC3D软件可近似模拟沥青混凝土的室内三轴试验，并可再现试件剪切破坏的发展过程；数值模拟试验所得到的应力-应变曲线及体积应变-轴向应变曲线与室内三轴试验结果基本一致。

颗粒流；沥青混凝土；细观模拟；应力-应变关系

沥青混凝土因其稳定和优良的防渗性能，近年来已被普遍应用在水利工程中。沥青混凝土由粗、细骨料、沥青和孔隙组成，其骨料颗粒数较多、形状各异、分布随机，为刚性或拟刚性材料；沥青是典型的黏弹性材料。长期以来，我国主要依据室内试验[1-2]及数值模拟方法[3-7]结合传统土力学理论进行沥青混凝土的性质研究，且主要集中在宏观性能方面。近年来，对水工筑坝材料细观性质的研究越来越受到国内外学者的重视[8-12]。

离散元法在20世纪70年代由Cundall首先提出[13]。按照单元类型可分为以块体为基本单元的块体离散元法，以及以圆盘或圆球为基本单元的颗粒离散元法。其中，颗粒离散元法专用于研究分散颗粒介质材料的细观力学特性。周里群等[14]人依据颗粒离散元法，生成二维沥青混凝土试件，完成了道路沥青混凝土切削仿真问题的模拟。武利强等[15]人使用二维颗粒流PFC2D软件探究了两种不同成型制样方式对沥青混凝土三轴试验结果的影响。以往对沥青混凝土的模拟，多从二维角度进行平面建模分析，本文使用三维颗粒流PFC3D(Particle Flow Code 3D)软件，对沥青混凝土三轴试验进行三维模拟，从细观角度解释沥青混凝土的应力-应变及变形特性，分析剪切破坏的发展过程，为沥青混凝土细观性质的研究提供参考。

1 颗粒流软件的基本原理

1.1 颗粒流基本计算原理

颗粒流PFC3D/2D软件以颗粒离散元理论为基础，将材料介质视为刚性颗粒群，各颗粒单元之间通过接触独立作用。

1.1.1 力-位移之间的关系

在软件运算过程中，各接触单元间通过力-位移定律实现接触力的传递。依据力-位移定律，各接触单元间的力矢量Fi可分解为：

(1)

法向接触力分量为：

(2)

式中：Kn为接触的法向刚度；Un为接触的法向相对位移；ni为法向单位向量。

切向接触力以增量形式计算：

(3)

式中：Ks为接触的切向刚度；ΔUs为每个时间步长的位移增量。

1.1.2 运动法则

每个接触单元的运动方式由其承受的合力与合力矩确定，可使用接触单元上一点的平移和颗粒的转动来表示：

平动方程：Fi=m(xi+gi)

(4)

式中：Fi为作用于颗粒单元上的合力；m为颗粒质量；gi为体积力加速度。

转动方程：Mi=Hi

(5)

式中：Mi为作用于颗粒单元上的合力矩；Hi为角动量。

1.2 接触本构模型选择

沥青混合料的材料特性通过选择合适的接触本构模型进行设置。两个接触集料单元之间的相对运动使用滑动模型进行模拟；集料单元间的弹性关系使用接触刚度模型加以描述。黏结模型限定了集料单元间黏结材料承受的法向力和剪力合力的最大值，沥青混凝土中的沥青属于黏弹性材料，根据胶浆理论[16]，选择使用平行黏结模型来模拟沥青胶浆[17]的力学特性。

平行黏结模型可模拟两个集料单元间介质材料的胶结状态，且可与滑动模型共同作用，一个平行黏结键被视为一系列分散在接触点横断面上的弹簧。这些有着法向和切向刚度的弹簧在接触单元之间传递力及力矩，两个接触单元受到力和力矩的作用会产生相对运动；这种作用在两接触单元上的力和力矩，与黏结材料的最大正应力和剪应力有关。当这些应力超过其相应的黏结强度，即导致平行黏结键断裂，力与力矩随即从断裂的平行黏结模型中移除。

2 数值模拟试验过程

2.1 室内三轴试验及数据

对某水库心墙用沥青混凝土进行室内三轴试验。沥青混凝土级配组成见表1，沥青用量用油石比表示，为6.6%；孔隙率为0.9%；按照设计配合比制成Φ100 mm×200 mm的沥青混凝土试件，密度为2.33 g/cm3；将制备好的试件在8℃下恒温3 h以上，在0.2 MPa、0.4 MPa、0.6 MPa三个围压下进行三轴试验；试验过程均在8℃±0.5℃的环境中进行。

表1 沥青混凝土集料级配

2.2 数值模拟试件的生成

在颗粒流PFC3D软件运算过程中，颗粒数目过多会导致运算缓慢，为提高软件计算效率，将沥青砂浆组颗粒最小直径控制为1 mm；按粒径将试件颗粒分为粗骨料组19 mm～9.5 mm、9.5 mm～4.75 mm、4.75 mm～2.36 mm，沥青砂浆组2.36 mm～1 mm，试件尺寸按h/d相同的原则，缩小为Φ50 mm×100 mm。首先，生成高100 mm，直径为50 mm的压力室。试件压力室由上下加载墙和控制围压的圆桶柔性墙构成，上下加载墙的刚度略大于颗粒刚度，柔性圆桶墙的刚度应远小于颗粒刚度。之后，按各颗粒组粒径的平均值，根据表1的沥青混凝土配合比，计算每个粒径组的颗粒个数，在压力室内以均匀随机分布生成各颗粒组。压力室及数值模拟试件见图1。

图1 压力室及数值模拟试件

2.3 三轴数值模拟试验的实现

三轴数值模拟试验通过命令流及内置FISH语言编写程序得以实现。首先，建立伺服系统，通过Callback命令，监测每个运算循环中作用在上下加载墙及圆桶墙上的颗粒的反作用力，并据此调节上下加载墙的轴向相对运动和圆筒墙的径向运动，进行试样的等压固结过程，待达到指定围压后，自动脱离伺服系统，要求围压的误差不超过1%；之后，以应变控制的方式，进行静力加载过程，即剪切过程，在这个过程中，使用History命令记录轴向应力、应变，径向应力、应变，体积应变等；最后，通过判断语句，当剪切破坏达到某一指定程度时，结束剪切过程。

2.4 参数设置

在0.2 MPa、0.4 MPa、0.6 MPa三种围压下，以室内三轴试验的应力-应变曲线为标准，反复对细观参数进行调整。具体步骤为：第一步，在某一围压下，根据图2中所示沥青混凝土室内三轴试验应力-应变曲线的特征参数，按照控制变量的原则，逐一调整细观参数：弹性模量E*、刚度比K*、法向连接强度σ、切向连接强度τ及摩擦系数μ，使三轴数值模拟试验得到的应力-应变曲线的特征参数与该围压下室内三轴试验的对应结果近似相等。第二步，将得到的细观参数应用在其它两种围压的数值模拟中，根据这两种围压得到的结果，对细观参数进行微调；第三步，将微调后的细观参数重新应用在先前围压的模拟中，根据结果对细观参数进行再次调整；重复若干次二、三步的工作，并以偏应力峰值、黏聚力c及摩擦角φ为主要指标，判断选择的细观参数是否合适。室内三轴试验与三轴数值模拟试验偏应力峰值结果见表2。

表2 不同围压下最大偏应力对比

室内三轴试验与三轴数值模拟试验强度参数结果对比，见表3。

图2 三轴试验应力-应变曲线特征参数示意图

由表2、表3所示对比结果，发现在三种围压下，三轴数值模拟试验得到的偏应力值与室内三轴试验的偏应力值相近，且黏聚力和内摩擦角误差均在允许范围之内，由此获得一组较为理想的细观参数，见表4。

表4 三轴数值模拟试验细观参数

3 数值试验结果及分析

3.1 结果对比

三轴数值模拟试验与室内三轴试验的应力(σ1-σ3)-应变ε1曲线对比见图3，体应变εv-轴向应变ε1曲线对比见图4。

3.2 应力-应变及变形特性分析

由图3、图4对比三种围压下的三轴数值模拟试验与室内三轴试验结果可看出，三轴数值模拟试验得到的应力-应变关系与室内三轴试验得到的应力-应变关系基本一致，最大偏应力都随围压的增高而增大，且均在达到峰值应力后发生应变软化，只是在应力-应变曲线的应变软化部分略有偏离。这是由于在偏应力达到最大值之前，颗粒间的接触趋于紧密，接触力主要为压力，体积变形处于减缩状态；峰值应力后，部分黏结键受压断裂，颗粒所受限制和约束力变小，材料抗剪能力下降，试件体积变形处于减胀状态；随着加载的进一步推进，部分黏结键逐渐受拉，在受拉黏结键、颗粒间咬合以及围压约束的共同作用下，保留了试件的部分抗剪能力；但由于数值模拟试件的颗粒形状单一，排列方式简单，咬合力较弱，导致试件软化强度降低。

图3 室内三轴试验与三轴数值模拟试验应力-应变曲线

图4 室内三轴试验与三轴数值模拟试验εv-ε1曲线

三轴数值模拟试验与室内三轴试验试件的减缩、剪胀特性基本相同，但也略有差别。三轴数值模拟试验试件的减缩均略大，这是因为在生成数值模拟试件时，为提高软件运行效率，对沥青砂浆组颗粒进行了简化，导致数值模拟试件孔隙率略大于室内试验试件。其次，三轴数值模拟试验的剪胀都随加载的推进而减弱，这是因为在这一过程中，黏结键断裂加剧，失去黏结键约束的颗粒出现翻滚及跨越，致使试件孔隙增大，剪胀加剧，但由于围压的约束，颗粒间的接触力已不再以压力为主，部分颗粒间的拉力作用更加突显，这使得剪胀相对减弱。

3.3 剪切位移场的变化

颗粒流PFC3D软件可通过记录运算过程中的模型变量、状态，实现剪切过程的跟踪和可视化。通过记录三轴数值模拟试件中颗粒位移，可以再现轴向加载过程中试件剪切破坏的发展过程。

图5展示了围压为0.2 MPa时数值模拟试件剪切破坏的发展过程。其中，图5(a)为固结阶段的位移场，这时的颗粒位移场较为杂乱；图5(b)为强度峰值出现之前的位移场，这时的位移场关于中心截面对称，且位移沿主应力加载方向均匀递减；图5(c)至图5(d)为强度峰值出现之后的位移场，在这个过程中，试件逐渐从中心膨胀，且最后沿与轴线约成45°角的斜面剪断，这与室内三轴试验试件的破坏方式相同；图5(e)为试件破坏时的位移矢量图。

图5 数值模拟试件剪切破坏的发展过程

4 结 论

(1) 使用颗粒流PFC3D软件可近似模拟沥青混凝土的室内三轴试验。数值模拟试验得到的应力-应变曲线与室内试验结果基本一致。在选择细观参数的过程中，应根据室内三轴试验应力-应变曲线中的特征参数进行调整，这样可提高数值模拟的效率和成功率。

(2) 使用颗粒流PFC3D软件，可再现试件剪切破坏的发展过程，对沥青混凝土结构破坏的过程有更深入的发现和了解。在沥青混凝土试件的剪切破坏过程中，位移场关于中心截面对称，且位移沿主应力加载方向均匀递减，试件逐渐从中心膨胀，且最后沿与轴线约成45°角的斜面剪断。

(3) 对比三轴数值模拟试验与室内三轴试验得到的应力-应变曲线，前段拟合较好，只是在应变软化部分略有偏离。这是因为数值模拟试件是由刚性圆球颗粒组成，颗粒形状单一、咬合力差。若生成形状不规则的骨料颗粒，可能会与实际更相符，这值得进一步尝试和探讨。

[1] 郭鹏飞,何建新,刘 亮,等.浇筑式沥青混凝土配合比设计优选方法研究[J].水利与建筑工程学报,2012,10(4):42-46.

[2] 李志强,张鸿儒,侯永峰,等.土石坝沥青混凝土心墙三轴力学特性研究[J].岩石力学与工程学报,2006,25(5):997-1002.

[3] 陈 松,宁 聪,王 岩，等.某水库沥青混凝土心墙坝有限元计算分析[J].水利与建筑工程学报,2016,14(6):152-156.

[4] 王建祥,唐新军,凤 炜.土石坝沥青混凝土心墙材料参数对其工作性状的影响[J].南水北调与水利科技,2013,11(4):178-181.

[5] 朱 晟,张美英,戴会超.土石坝沥青混凝土心墙力学参数反演分析[J].岩土力学,2009,30(3):635-639.

[6] 余 林.沥青混凝土心墙高坝应力应变特性初探[J].水利与建筑工程学报,2016,14(6):107-112,125.

[7] 邓传贵,沈振中,凌春海.高沥青混凝土心墙坝动力反应特性[J].南水北调与水利科技,2012,10(5):147-152.

[8] 邵 磊,迟世春,张 勇,等.基于颗粒流的堆石料三轴剪切试验研究[J].岩土力学,2013,34(3):711-720.[9] 杨 贵,肖 杨,高德清.粗粒料三维颗粒流数值模拟及其破坏准则研究[J].岩土力学,2010,31(S2)：402-406.[10] 齐 阳,唐新军,李晓庆.粗粒料工程力学性质的细观模拟[J].南水北调与水利科技,2014,12(4):156-160.

[11] Deluzarche R, Cambou B. Discrete numerical modelling of rock-fill dams[J]. International Journal for Numerical and Analytical Methods in Geomechanics, 2006,30(11):1075-1096.

[12] You Z, Buttlar W G. Discrete element modeling to predict the modulus of asphalt concrete mixtures[J]. Journal of Materials in Civil Engineering, 2004,16(2):140-146.

[13] Cundall P A. A computer model for simulating progressive large scale movement in block rock system[J]. Proceedings of the International Symposium Rock Fracture, 1971,1(Ⅱ-8):129-136.

[14] 周里群,李 军,邢 国.基于PFC的沥青混凝土铣削仿真及刀具安装优化[J].湘潭大学自然科学学报,2012,34(4):95-98.

[15] 武利强,朱 晟,苏佩珍,等.成型方法对沥青混凝土力学性质影响的细观研究[J].水电能源科学,2012,30(2):106-110.

[16] 常明丰.多孔沥青路面微观力学特性与空隙衰变行为研究[D].西安:长安大学，2009:33-45.

[17] 裴建中.沥青路面细观结构特性与衰变行为[M].北京:科学出版社，2010:3-23.

Three Dimensional Mesoscopic Simulation of Hydraulic Asphalt Concrete Triaxial Test

LIU Xuan, TANG Xinjun, WANG Jianxiang

(SchoolofWaterConservancyandCivilEngineering,XinjiangAgriculturalUniversity,Urumqi,Xinjiang830052,China)

Based on the particle discrete element theory, three dimensional numerical simulation of the asphalt concrete is of great significance to the meso-performance of hydraulic asphalt concrete research. Through numerical simulation, the results show that the laboratory triaxial test of asphalt concrete can be simulated by particle flow software (PFC3D) if we select the appropriate contact constitutive model and mesoscopic parameters between elements according to the material composition and properties of asphalt concrete, and the development process of shear failure can be reproduced. The stress-strain curves and the volumetric strain-axial strain curves obtained by numerical simulation are basically consistent with the results of laboratory triaxial test.

particle flow code; asphalt concrete; mesoscopic simulation; stress-strain relation

10.3969/j.issn.1672-1144.2017.03.007

2017-02-20

2017-03-25

新疆高校科研计划重点资助项目(XJEDU2016I024)

刘 璇(1991—),女，新疆奇台人，硕士研究生，研究方向为水工结构工程数值分析。 E-mail：1379182998@qq.com

唐新军(1959—)，男，陕西蒲城人，博士，教授，博士生导师，主要从事当地材料坝设计理论教学与科研工作。 E-mail：tangxj59@163.com

TV431+.5

A

1672—1144(2017)03—0035—05