張德育 黃曉明 高英

(東南大學交通學院，江蘇南京210096)

隨著公路交通量的增加和重載車輛的增多，車轍已成為瀝青路面的主要病害之一，嚴重影響路面行車安全和舒適性，公路瀝青路面的車轍已經(jīng)成為人們普遍關注的問題.對于瀝青混合料的永久變形性能，各國道路工作者提出了多種研究方法，其中以現(xiàn)象試驗法最為普遍，但這種方法不僅需要消耗大量時間和成本，重復性和再現(xiàn)性也較差［1-2］.單軸蠕變試驗作為常用的瀝青混合料高溫永久變形性能試驗亦存在這種缺點.

瀝青混合料由集料、瀝青膠結(jié)料及空隙構成，集料在瀝青混合料中呈現(xiàn)不連續(xù)的顆粒特性，空隙結(jié)構亦導致瀝青混合料呈現(xiàn)不連續(xù)特性，因此，瀝青混合料是一種多相復合材料，其材料性質(zhì)較為復雜，實際的應力、應變場并不連續(xù)，集料、瀝青膠結(jié)料及空隙在混合料中扮演著不同的力學角色，呈現(xiàn)出不均勻和離散的力學特性.長期以來，道路工程領域一直沿用連續(xù)均質(zhì)的力學方法分析瀝青混合料的永久變形性能并指導材料設計，這與瀝青混合料的不連續(xù)特征不符，難以真實地反映瀝青混合料的材料力學特性.

近年來，離散元方法逐漸被引入到瀝青混合料性能分析中.離散元方法可以反映瀝青混合料的不連續(xù)特性，從細觀角度研究瀝青混合料的永久變形性能，建立的虛擬試驗的重復性和再現(xiàn)性好，經(jīng)濟方便，能夠克服傳統(tǒng)試驗的不足.但在以往的研究中，幾乎都是采用二維(2D)離散元模型對瀝青混合料進行分析［3-8］，而二維結(jié)構與三維(3D)結(jié)構相差較大，二維離散元模型具有一定的局限性.因此，文中采用離散元軟件PFC3D生成了包含集料、瀝青砂漿及空隙的瀝青混合料三維離散元虛擬試件，建立了瀝青混合料三維離散元虛擬單軸蠕變試驗，并與二維虛擬試驗及室內(nèi)試驗結(jié)果進行了比較，驗證了三維離散元虛擬試驗方法的正確性.

1 三維離散元虛擬試件的建立

為了兼顧離散元模型的計算效率和瀝青砂漿粘彈性參數(shù)的獲取，將粒徑大于1.18 mm的集料視為粗集料，粒徑小于1.18 mm的細集料和瀝青膠結(jié)料視為均質(zhì)的瀝青砂漿，因此，文中建立的瀝青混合料離散元虛擬試件由集料、瀝青砂漿和空隙構成.首先根據(jù)瀝青混合料的油石比、空隙率及集料的級配(見表1)等計算出各檔集料的體積，然后將各檔集料逐次投放到由墻構成的直徑為100 mm，高為150 mm的瀝青混合料虛擬試件尺寸空間中，并通過循環(huán)保證生成的集料球單元不重疊.生成的具有級配特征的集料如圖1所示.

表1 AC-16瀝青混合料級配Table 1 Gradation of AC-16 asphalt mixture

圖1 具有級配特征的集料Fig.1 Graded aggregates

待具有級配特征的集料生成后，在瀝青混合料試件尺寸空間內(nèi)生成148200個規(guī)則排列(每個單元與四周的6個單元相接觸)且半徑為1mm的離散單元(如圖2所示).遍歷新生成的離散單元，逐一判斷新單元與原集料球單元的位置是否重疊，如重疊則將新單元視為集料單元，否則視為瀝青砂漿單元，并將每個集料范圍內(nèi)的離散單元定義為剛性“聚?！?，以節(jié)省計算時間.判斷結(jié)束后刪除原有集料單元，生成如圖3所示的三維瀝青混合料虛擬試件.

圖2 規(guī)則排列的離散單元Fig.2 Rugularly-packed discrete elements

圖3 三維瀝青混合料虛擬試件Fig.3 3D virtual specimen of asphalt mixture

由于空隙對瀝青混合料的永久變形性能影響較大，為了在瀝青混合料離散元虛擬試件中生成空隙，在瀝青砂漿單元中隨機刪除5928個單元作為4%的空隙.三維離散元虛擬試件中的空隙分布如圖4所示.由于室內(nèi)試件及現(xiàn)場取芯試件內(nèi)部空隙分布的復雜性，文中將虛擬試件內(nèi)部的空隙分布視為隨機分布.

圖4 瀝青混合料空隙分布Fig.4 Distribution of air voids in asphalt mixture

2 微觀接觸模型及參數(shù)

2.1 微觀接觸模型

圖5 瀝青混合料微觀組成成分的接觸行為Fig.5 Interactions among microscale components in asphalt mixture

在瀝青混合料離散元模型中，瀝青混合料內(nèi)部微觀組成成分之間有4種接觸，分別為集料內(nèi)部單元、相鄰集料之間、瀝青砂漿內(nèi)部單元、瀝青砂漿和集料之間的接觸，如圖5所示.根據(jù)PFC3D中不同接觸模型的特點，采用接觸剛度模型描述集料內(nèi)部單元之間的接觸行為，采用接觸剛度模型和滑移模型描述相鄰集料之間的接觸行為，采用伯格斯接觸模型和接觸粘結(jié)模型描述瀝青砂漿內(nèi)部單元之間以及瀝青砂漿單元和集料之間的接觸行為.

2.2 模型參數(shù)

PFC3D中微觀接觸模型的參數(shù)很難直接從室內(nèi)試驗得到，需要建立微觀模型參數(shù)與宏觀試驗參數(shù)之間的關系，通過室內(nèi)試驗得到宏觀試驗參數(shù)，進而根據(jù)已建立的轉(zhuǎn)化關系得到微觀模型參數(shù).

考慮到瀝青砂漿的粘彈性，采用接觸粘結(jié)模型和微觀伯格斯接觸模型(如圖6所示)來描述瀝青砂漿內(nèi)部單元之間及瀝青砂漿與集料之間的接觸行為.圖6 中 E1、η1、E2、η2為宏觀伯格斯模型參數(shù)，Kmn、Cmn、Kkn、Ckn為微觀伯格斯接觸模型法向參數(shù)，Kms、Cms、Kks、Cks為微觀伯格斯接觸模型切向參數(shù)，m1、m2表示顆粒單元，fs為摩擦系數(shù).

圖6 伯格斯模型Fig.6 Burger’s model

由于文中主要研究瀝青混合料的蠕變特性，不考慮材料的破壞，接觸粘結(jié)強度取較大的值以避免單元粘結(jié)的破壞.瀝青砂漿的粘彈性參數(shù)可以由DSR動態(tài)剪切試驗獲取，文中采用的瀝青砂漿粘彈性參數(shù)分別為 E1=4.64 MPa，η1=83.34 MPa·s，E2=0.82MPa，η2=10.00 MPa·s.根據(jù)已有的研究成果［9］，微觀伯格斯接觸模型參數(shù)(見圖6(b)和圖6(c))與宏觀伯格斯模型參數(shù)(見圖6(a))之間的轉(zhuǎn)化關系為

式中，L為相鄰單元的球心距，υ為瀝青砂漿泊松比.

將集料視為彈性體，采用接觸剛度模型來描述集料內(nèi)部單元的接觸行為，采用接觸剛度模型和滑動模型來描述集料之間的接觸行為.由于離散單元之間的排列方式為規(guī)則的矩形排列，根據(jù)已有文獻的研究成果［10］，集料的接觸剛度與宏觀彈性模量的轉(zhuǎn)化關系為

式中:kn、ks分別為接觸剛度模型的法向剛度和切向剛度;E 為集料的宏觀彈性模量，文中取 55.5GPa［4，6，11］;R為集料單元的半徑;υ'為集料泊松比.

3 時溫等效原理

粘彈性離散元模型能夠更好地模擬瀝青混合料的粘彈性行為，而粘彈性模型的計算相當耗費時間［9］，以 CPU 2.8 GHz、RAM 2.0 GB 計算機進行真實時間下的虛擬試驗將會花費數(shù)年的時間.因此，如果計算中采用真實加載時間，將會對虛擬試驗的進行造成很大的困難，甚至無法進行.瀝青混合料作為粘彈性材料，其變形依賴于時間和溫度，為了減少計算時間，文中采用基于時溫等效原理的方法對虛擬試驗的計算效率進行優(yōu)化.

時溫等效原理是將材料對于外界荷載的響應進行時間和溫度的等效轉(zhuǎn)化，因此可以提高粘彈性材料的加載溫度進而減少加載時間，再將高溫下的結(jié)果通過時溫等效原理向低溫條件下進行轉(zhuǎn)化，轉(zhuǎn)化關系為

式中:ε(T，t)為真實溫度T和加載時間t下的蠕變應變;ε(Tr，tr)為參考溫度Tr和縮減時間tr下的蠕變應變;αT為移位因子.

在瀝青混合料的研究中，時溫等效原理常用來建立主曲線，進而可以將一定時間、溫度范圍內(nèi)的試驗測定結(jié)果拓延到更廣泛的時間溫度空間中去，大大減少材料研究中的試驗工作，并在相當程度上降低對試驗裝置的技術要求.同樣的，可以利用時溫等效原理減少粘彈性離散元模型的計算時間.瀝青混合料粘彈性伯格斯模型中真實溫度和時間下的蠕變應變可以表達為

式中，σ0為軸向應力.

當溫度升高到參考溫度時，粘彈性參數(shù)也會隨之改變.假定 E1、E2、η1r、η2r為參考溫度下的伯格斯模型參數(shù)，則在縮減時間下的蠕變應變?yōu)?/p>

瀝青混合料在真實溫度和參考溫度下的伯格斯模型參數(shù)是不同的，但不同溫度下的蠕變應變在不同的加載時間下可以進行等效.通過求解式(11)－(14)，可以得到參考溫度下的伯格斯模型參數(shù):

在以往的文獻中，離散元粘彈性模擬分析都是采用真實時間下的伯格斯模型參數(shù)，因此計算時間會很長.如果采用參考溫度和縮減時間下的伯格斯模型參數(shù)，計算時間就取決于縮減時間.當移位因子取較大值時，縮減時間較真實時間大幅度縮短，因此模型的計算時間將會大大減少.經(jīng)過試算，當移位因子不大于104時，計算結(jié)果較為穩(wěn)定(如圖7所示)，這與You等［12］的研究結(jié)論較為一致.因此，文中取104作為移位因子的取值，那么虛擬試驗的加載時間僅為真實時間的1/104，計算時間從數(shù)年減少到了數(shù)天.

圖7 不同移位因子下的瀝青混合料軸向應變Fig.7 Axial strain of asphalt mixture with different shifting factors

4 虛擬單軸蠕變試驗

4.1 三維虛擬試驗結(jié)果

在上述建立的虛擬試件基礎上，對虛擬試件施加0.7 MPa的靜載.在PFC3D中，通常是通過控制墻的移動速度來對試件進行加載，因此本虛擬試驗通過在試件的上下表面各設定一面墻，固定試件下表面的墻，通過PFC3D中的fish語言編寫伺服控制程序，不斷調(diào)整試件上表面墻的移動速度使軸向應力達到恒定的值.由于文中采用了基于時溫等效的計算優(yōu)化方法且移位因子取為104，當室內(nèi)試驗加載時間為100 s時，虛擬試驗僅需加載0.01 s，計算時間大大減少.

由于瀝青路面車轍通常發(fā)生在高溫條件下，尤其是接近瀝青軟化點的時候車轍量會迅速增加［13］.因此，文中的試驗溫度選擇具有代表性的60℃.瀝青混合料試件的直徑為100 mm，高為150 mm，最佳油石比為4.8%，空隙率為4%，混合料級配如表1所示.分別進行相同條件下的三維離散元虛擬試驗和室內(nèi)試驗，并進行比較分析，結(jié)果如圖8所示.

圖8 三維虛擬試驗與室內(nèi)試驗結(jié)果比較Fig.8 Comparing of results of 3D virtual test and laboratory measurement

從圖8可以看出，三維離散元虛擬試驗結(jié)果與真實試驗的結(jié)果較為吻合，但不完全一致，其原因可能是虛擬試件中集料的形狀與真實試件中集料的形狀有所差別，虛擬試件內(nèi)部空隙率的分布與真實試件也有所不同，而且室內(nèi)試驗的加載條件也沒有虛擬試驗理想.從總體上看，三維離散元虛擬試驗能夠較為準確地預估瀝青混合料的永久變形性能，可以作為瀝青混合料永久變形性能分析的輔助手段.

4.2 三維和二維虛擬試驗結(jié)果比較

為了比較三維和二維離散元模型預估瀝青混合料永久變形性能的準確性，文中還進行了相同試驗條件下的二維虛擬單軸蠕變試驗.二維離散元模型直接由三維離散元模型相互垂直的兩個切片獲得(如圖9所示).三維和二維離散元模型對比結(jié)果如圖10(a)所示，不同荷載下的三維和二維離散元模型對比結(jié)果如圖10(b)所示.

圖9 瀝青混合料二維離散元模型的獲取Fig.9 Illustration of obtaining 2D DE model of asphalt mixture

圖10 三維和二維離散元模型計算結(jié)果比較Fig.10 Comparing of results of 3D and 2D DE models

從圖10(a)可以看出:加載初期軸向應變的三維預測值大于二維預測值;隨著加載的進行，三維預測值則小于二維預測值，且差異隨之增大.從圖10(b)可以看出，軸向應變的二維預測值較三維預測值及室內(nèi)試驗值要大，且隨著荷載的增大差異也隨之增大.分析其原因可能是二維模型低估了集料的嵌鎖作用.三維離散元模型在預測瀝青混合料高溫永久變形性能方面較二維離散元模型更為理想，可以準確地預估瀝青混合料的高溫永久變形性能.

5 結(jié)語

文中建立了包括瀝青混合料三維離散元虛擬試件的生成、瀝青混合料微觀組成成分微觀接觸模型的選取、微觀接觸模型參數(shù)的獲取以及加載方式的實現(xiàn)等在內(nèi)的完整的瀝青混合料三維虛擬單軸蠕變試驗方法，提出了基于時溫等效原理的離散元粘彈性模型計算優(yōu)化方法，將瀝青混合料粘彈性離散元模型的計算時間從數(shù)年減少到了數(shù)天，從而解決了三維離散元粘彈性模型的計算時間問題.

二維虛擬試驗結(jié)果大于三維虛擬試驗及室內(nèi)試驗結(jié)果，且隨著荷載的增大差異也隨之增大，說明三維離散元模型在預測瀝青混合料高溫永久變形性能方面較二維離散元模型更為理想，所建立的瀝青混合料三維離散元虛擬單軸蠕變試驗能夠準確地預估瀝青混合料高溫永久變形性能，可以作為瀝青混合料永久變形性能分析的輔助手段.

［1］譚積青.粗集料形態(tài)與瀝青混合料級配組成虛擬力學試驗的基礎方法研究［D］.廣州:華南理工大學土木與交通學院，2006.

［2］胡霞光.瀝青混合料微觀力學分析綜述［J］.長安大學學報:自然科學版，2005，25(2):6-9.Hu Xia-guang.Review on asphalt mixture micromechanics analysis［J］.Journal of Chang＇an University:Natural Science Edition，2005，25(2):6-9.

［3］Chang G K，Meegoda J N.Micromechanical simulation of hot mix asphalt［J］.Journal of Engineering Mechanics，1997，123(5):495-503.

［4］Buttlar W G，You Z P.Discrete element modeling of asphalt concrete:microfabric approach ［J］.Transportation Research Record，2001，1757:111-118.

［5］王端宜，張肖寧，王紹懷.用虛擬試驗方法評價瀝青混合料的級配類型［J］.華南理工大學學報:自然科學版，2003，31(2):48-51.Wang Duan-yi，Zhang Xiao-ning，Wang Shao-huai.Evaluation on grading type of asphalt mixture with virtual test method［J］.Journal of South China University of Technology:Natural Science Edition，2003，31(2):48-51.

［6］You Z P，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.

［7］Abbas A，Masad E，Paapagiannakis T，et al.Micromechanical modeling of the viscoelastic behavior of asphalt mixtures using the discrete element method［J］.International Journal of Geomechanics，2007，7(2):131-139.

［8］Kim H，Wagoner W P，Buttlar W G.Simulation of fracture behavior in asphalt concrete using a heterogeneous cohesive zone discrete element model［J］.Journal of Materials in Civil Engineering，2008，20(8):552-563.

［9］Liu Y，Dai Q L，You Z P.Viscoelastic model for discrete element simulation of asphalt mixtures［J］.Journal of Engineering Mechanics，2009，135(4):324-333.

［10］Thornton C.The conditions for failure of a face-centered cubic array of uniform rigid spheres［J］.Geotechnique，1979，29(4):441-459.

［11］You Z P.Development of a micromechanical modeling approach to predict asphalt mixture stiffness using the discrete element method［D］.Urbana-Champaign:Department of Civil and Environmental Engineering，University of Illinois at Urbana-Champaign，2003.

［12］You Z P，Liu Y，Dai Q L.Three-dimensional microstructural-based discrete element viscoelastic modeling of creep compliance tests for asphalt mixtures［J］.Journal of Materials in Civil Engineering，2011，23(1):79-87.

［13］黃曉明，范躍武，趙永利，等.高速公路瀝青路面高溫車轍的調(diào)查與分析［J］.公路交通科技，2007，24(5):16-20.Huang Xiao-ming，F(xiàn)an Yue-wu，Zhao Yong-li，et al.Investigation and test of expressway asphalt pavement hightemperature performance［J］.Journal of Highway and Transportation Research and Development，2007，24(5):16-20.