潘 伶，張晉銘，呂志田，林旭健

（1.福州大學機械工程及自動化學院，福建 福州 350108；2.福州大學土木工程學院，福建 福州 350108）

瀝青混合料的破壞通常包括瀝青-集料界面黏附破壞、黏聚破壞以及黏附和黏聚混合破壞.現(xiàn)有的研究主要集中于定性分析瀝青包裹集料的程度［1］，或通過拉拔剪切試驗［2］、表面能測試［3］研究界面的失效，較少從微觀角度解釋界面的破壞機理.

瀝青-集料界面黏附性微觀力學行為的研究方法主要分為有限元法（FEM）、原子力顯微鏡（AFM）試驗法和分子動力學（MD）模擬法.FEM可以在較大時間和空間上進行模擬［4］，便于施加動態(tài)車輪載荷，但結果的準確度依賴于材料的本構關系、內(nèi)聚力理論等數(shù)學模型.AFM試驗法精度高、重復性好，可在干燥或潮濕、加熱或冷卻等多種環(huán)境下測試［5］，但加載方式簡單，與實際工況相差較大.MD模擬法目前多用于研究瀝青混合料的各項性能，如瀝青的材料性質(zhì)［6］和自愈性能［7?8］等，少數(shù)用于探究瀝青-集料界面的黏附性和力學行為.Xu等［9］建立了瀝青-單層集料界面的MD模型，在界面之間增加水分子層研究水對黏附功的影響，并進行拉伸模擬得到其應力-位移響應.

實際工況中瀝青在兩集料間承受載荷，本文建立了兩集料間瀝青的模型，在不同溫度（T）和拉伸速率（υ）下進行MD模擬，考察應力-位移關系，從原子尺度研究瀝青-集料界面的破壞形式和黏附機理.

1 模型的建立

1.1 瀝青模型

瀝青是由不同分子量的碳氫化合物及其非金屬衍生物組成的一種黑褐色混合物.根據(jù)溶解度，可分為瀝青質(zhì)、飽和分和膠質(zhì)［10］，各組分分子模型見圖1.

圖1 瀝青各組分分子模型Fig.1 Molecular structures of asphalt three?components

在Material studio軟件中通過加入相應數(shù)量的各組分分子，在垂直方向（Z向）建立收縮邊界條件的無定型晶胞，經(jīng)過幾何優(yōu)化，先后在正則系綜NVT和等溫等壓系綜NPT條件下馳豫200 ps，時間步長為1 fs，得到受限瀝青分子模型見圖2.

圖2 受限瀝青分子模型Fig.2 Molecular model of restricted asphalt

1.2 集料模型

集料是一種廣泛應用于建筑工程的粗顆粒材料，包括沙、花崗巖、石灰?guī)r和石英等.模型中集料采用沙和花崗巖中的主要成分SiO2，其晶胞模型見圖3（a）.為建立兩集料間瀝青模型以及方便數(shù)據(jù)提取，將該晶胞模型改造為正交幾何模型，在X，Y，Z方向定量重復晶胞，建立與瀝青尺寸相符的集料層模型，見圖3（b）.

圖3 SiO2晶胞模型和集料層模型Fig.3 Unit cell model of SiO2 and SiO2 aggregate layer model（size：nm）

1.3 兩集料間瀝青模型

為了模擬瀝青在兩集料間拉壓和剪切等復雜工況及宏觀實物試驗［1?3］，建立兩集料間瀝青模型，見圖4.由圖4可見：模型在X、Y方向的長度均為6.917 nm，SiO2集料層高度均為5.406 nm，兩集料層中間填充的瀝青分子高度為13.718 nm；將X、Y方向設置成周期性邊界，Z方向設置成非周期性邊界，上下集料分為用于施加邊界條件的剛性層、提供環(huán)境影響因素的恒溫層、提取應力參數(shù)的自由變形層；所有集料原子之間相互作用采用Lennard?Jones（9?6）勢［11?13］，模型中共有原子99 696個，其中瀝青層有64 416個原子.

圖4 兩集料間瀝青模型Fig.4 Model of asphalt confined in two aggregates

2 MD模擬

MD模擬過程采用大規(guī)模原子/分子并行模擬器（LAMMPS）編程實現(xiàn)，分子力場為聚合物一致性力場（PCFF），該力場的參數(shù)是通過大量的試驗數(shù)據(jù)和精確的量子計算擬合得出，適用于模擬聚合物、有機物、生物分子體系以及部分金屬.PCFF力場，包含非鍵結作用項、鍵伸縮項、鍵角彎曲項、二面角扭曲項、離平面振動項和多種交叉作用項，其勢函數(shù)和各項參數(shù)的選取見文獻［11?13］.

MD模擬過程分為馳豫階段和拉伸階段.在馳豫階段：恒溫層和瀝青層使用正則系綜NVT將溫度控制在預設值，控溫方法選擇Nose?Hoover法，設置溫度阻尼系數(shù)為100 fs；將整個系統(tǒng)馳豫200 ps，最終體系的總能量收斂于5.40×106kJ/mol.在拉伸階段：固定下剛性層，在上剛性層施加拉伸速率υ，取υ=0.01、0.10、1.00、10.00和30.00 m/s；瀝青層設置為微正則系綜NVE，使拉伸過程中的能量被恒溫層消散；為保證拉伸過程中瀝青和集料可以達到完全分離，取垂直相對運動距離為6.0 nm.MD模擬中牛頓運動方程的時間積分采用Velocity?Verlet算法，時間步長取1 fs.

3 結果與討論

3.1 模型驗證

用內(nèi)聚能密度（CED）衡量分子間的相互吸引力，以評估瀝青分子間的相互作用［14］；用溶解度參數(shù)δ評估瀝青中膠體的穩(wěn)定性［15］，量化溶劑-溶質(zhì)之間的相互作用，排除溶劑-溶劑和溶質(zhì)-溶質(zhì)分子間的相互作用，δ的計算公式見式（1）.瀝青中分子間相互作用是通過范德華力和靜電力實現(xiàn)，因此δ也可用式（2）進行表達.

式中：δvdw和δele分別為范德華力和靜電力產(chǎn)生的溶解度參數(shù)分量.

T=298 K，瀝青的內(nèi)聚能密度以及溶解度參數(shù)的MD計算值和試驗值見表1.由表1可見：MD模擬結果與試驗值［16］一致；在溶解度參數(shù)中，范德華力分量遠大于靜電力分量.

表1 瀝青的內(nèi)聚能密度以及溶解度參數(shù)的MD計算值和試驗值Table 1 MD simulation and test values of cohesive energy density and solubility parameters of asphalt

3.2 應力-位移關系

3.2.1 應力-位移關系曲線

記錄拉伸過程中SiO2集料變形層在Z方向來自瀝青層的牽引力Ftrac，采用式（3）計算界面的拉伸應力σ，其應力峰值為σc.

式中：A為瀝青-集料界面面積.

υ=10.00 m/s、T=298 K下，MD模擬過程中界面的拉伸應力（σ）-位移（sz）曲線見圖5.由圖5可見：隨著Z方向位移sz的增大，拉伸應力σ先增大后減小，最終歸0；起始階段σ隨著sz的增加呈線性增加；在sz=0.95 nm時，σ達到峰值，即σc=22 MPa；σ達到最大值后，隨著sz的進一步增大，σ快速減小，瀝青處于軟化階段；當sz=3.50 nm時，σ減小的速率變慢；在sz=4.80 nm時，σ減小至0，表明分離完成.

圖5 MD模擬過程中界面的拉伸應力-位移曲線Fig.5σ?sz curve of interface in MD simulation

3.2.2 拉伸過程動態(tài)分析

MD模擬比常規(guī)實物試驗的優(yōu)勢是可以從原子尺度動態(tài)地觀察拉伸過程中界面的破壞過程.υ=10.00 m/s、T=298 K時，σ?sz曲線中對應點模型的正視圖和截面俯視圖見圖6.由圖6可見：在A點（sz=0 nm），系統(tǒng)處于馳豫結束的平衡狀態(tài)，由截面俯視圖可以看出瀝青分子緊緊團聚在一起，幾乎觀察不到集料層；在B點（sz=0.50 nm），可以透過瀝青分子間隙觀察到小部分集料界面，隨著位移增大，拉伸應力快速升高，由于瀝青是黏彈性材料，在瀝青分子之間出現(xiàn)細小的孔隙；在C點（sz=0.95 nm），拉伸應力達到峰值，通過俯視圖可以觀察到，下方的瀝青-集料界面裸露出較大塊的集料原子，表明此處瀝青與集料開始分離；從C點到D點（sz=2.40 nm），拉伸應力由峰值迅速下降，此階段俯視圖中瀝青-集料界面分離面積不斷增大，從正視圖也可以觀察到界面附近出現(xiàn)孔洞，并不斷增大，而瀝青內(nèi)部的微小孔隙融合；在E點（sz=3.90 nm），觀察到σ?sz曲線的斜率突然變小，同時界面上的分離面積不再增大，表明破壞形式由以黏附破壞為主轉為以黏聚破壞為主，這是由于界面分離到E點后，瀝青連接處面積減小，其黏聚強度小于界面的黏附強度；在F點（sz=5.30 nm），可以觀察到應力在0周圍波動，拉伸應力波動是由于集料表面黏附有瀝青分子，從正視圖可以觀察到上下集料已經(jīng)完全分離，殘留的組分主要為飽和分和膠質(zhì)分子.

圖6 σ?sz曲線中對應點模型的正視圖和截面俯視圖Fig.6 Forward view and section top view of the corresponding point model inσ?sz curve

3.3 溫度的影響

為了研究溫度T對σ?sz曲線的影響，在υ=10.00 m/s下，取223、273、323、372 K 4個溫度進行MD模擬.不同溫度下界面的σ?sz曲線見圖7.由圖7可見：溫度對σ?sz曲線的應力峰值σc及σc對應的位移szc影響很大；溫度越低，σc越大；溫度升高導致達到σc所需的位移增大，這是由于瀝青是黏彈性材料，溫度越高，流動性就越高，需要更大的位移達到σc；在T=223 K時，達到σc后的軟化過程明顯要快于其他3個溫度.

圖7 不同溫度下界面的σ?sz曲線Fig.7σ?sz curves of interface at different temperatures

3.4 拉伸速率的影響

T=298 K時，不同拉伸速率下界面的σ?sz曲線見圖8.由圖8可見：不同拉伸速率下，應力峰值σc前曲線的斜率幾乎一致；σc后，拉伸速率越大，分離所需的位移越大.

圖8 不同拉伸速率下界面的σ?sz曲線Fig.8σ?sz curves of interface under different tensile rates

3.5 內(nèi)聚力模型

內(nèi)聚力模型（CZM）［17?18］是一種描述界面力學行為的模型，其通過作用在裂紋面上的張開力與張開位移的等效關系描述裂紋尖端附近復雜變形區(qū)域的力學行為.對MD模擬得到的拉伸應力-位移曲線進行CZM模型擬合，可為有限元等更大尺度模擬提供精確的數(shù)學模型.基于最小二乘法擬合拉伸應力-位移曲線，通過對比多項式、二次、三次和指數(shù)函數(shù)等常見內(nèi)聚力模型，得到式（4）所示的指數(shù)型CZM模型.

T=273 K時，MD模擬界面σ?sz曲線及其CZM模擬曲線見圖9.內(nèi)聚力關系表明，界面損傷是在應力達到峰值時發(fā)生的.由圖9可見，應力峰值后的軟化曲線表明隨著界面損傷的增加，界面拉伸應力從峰值逐漸降為0，直至界面完全分離.

圖9 MD模擬界面σ?sz曲線及其CZM模擬曲線Fig.9 MD simulation interfaceσ?sz curves and its CZM simulation curves

3.6 卸載位置的影響

在υ=10.00 m/s，T=273 K下，研究了加載-卸載-加載對拉伸過程中界面σ?sz曲線的影響，試驗過程分為3個階段：第Ⅰ階段，對兩集料間瀝青進行拉伸，分別加載到應力峰值后（σ=32 MPa）、應力峰值前（σ=29 MPa）；第Ⅱ階段，卸載；第Ⅲ階段，用第Ⅰ階段的拉伸速率重新進行加載.圖10為“加載-卸載-加載”過程界面的σ?sz曲線.由圖10（a）可見：達到應力峰值后卸載，σ會迅速下降；第Ⅲ階段再次加載后，σ也無法恢復到卸載時的拉伸應力值，說明此時系統(tǒng)在第Ⅰ階段的拉伸過程中已經(jīng)產(chǎn)生了不可逆的損傷，失去了大部分彈性.由圖10（b）可見：應力峰值前卸載，σ略微下降，并保持在較高水平；第Ⅲ階段再次加載后，應力峰值仍為32 MPa，說明在達到應力峰值前卸載，瀝青的黏彈性不會受到損傷.綜上，達到應力峰值前卸載，瀝青保持彈性狀態(tài)，此時的變形是可逆的；達到應力峰值后卸載，瀝青的黏彈性受到不可逆的損傷.

圖10 “加載-卸載-加載”過程界面的σ?sz曲線Fig.10σ?sz curves of interface under“l(fā)oading?unloading?loading”process

4 結論

（1）通過分子動力學（MD）模擬兩集料間瀝青的拉伸力學行為，從原子尺度動態(tài)地觀察界面破壞的過程，得到不同工況下的拉伸應力-位移曲線，當拉伸應力超過應力峰值后，破壞由以黏附破壞為主轉變?yōu)橐责ぞ燮茐臑橹?

（2）溫度T越高，應力峰值σc越??；拉伸速率υ越大，應力峰值σc越大；不同拉伸速率下，應力峰值之前的拉伸應力-位移曲線的斜率幾乎一樣；達到應力峰值后，υ越大，完成分離所需的位移越大.

（3）通過內(nèi)聚力模型擬合MD模擬得到拉伸應力-位移曲線，提出了適用于宏觀數(shù)值分析的內(nèi)聚力模型.

（4）瀝青在拉伸應力達到應力峰值σc之前保持彈性狀態(tài)，此時變形可逆；當拉伸應力超過應力峰值σc后，受到的損傷不可逆.