基于離散單元法的粗粒土蠕變級配演化規律

摘要：為探究粗粒土蠕變過程中顆粒級配演化規律，開展了粗粒土蠕變的離散單元法分析?；趹Ωg理論建立了粗粒土雙軸蠕變離散元數值模型；通過粗粒土室內三軸蠕變試驗結果標定了細觀模型參數；開展了不同應力水平下的粗粒土蠕變數值試驗，分析了顆粒破碎的模式及其演化規律、級配曲線及其分維值的演化規律。結果表明：基于應力腐蝕理論，采用離散元研究粗粒土在蠕變過程中的級配演化規律是可行的；粗粒土在蠕變過程中，顆粒破碎前期劇烈，后期微弱；隨著應力水平的增大，顆粒破碎率增大且達到穩定值所需時間更長；蠕變結束后，低應力水平下的試樣僅在端部存在少量顆粒破碎，高應力水平下的顆粒同時存在劈裂和棱角破碎兩種破碎模式；隨著應力水平的增大，級配曲線逐漸向最終級配曲線趨近，級配曲線在蠕變后期變化較小；粗粒土的級配在蠕變過程中始終具有分形特性，隨應力水平的增大分維值增大，分維值與顆粒破碎率呈良好的線性關系。研究成果有助于深入分析粗粒土蠕變規律和機理。

關 鍵 詞：粗粒土；顆粒級配；蠕變；腐蝕理論；離散元數值模型；三軸蠕變試驗；顆粒破碎；分維值

中圖法分類號：TU444 文獻標志碼：ADOI：10.16232/j.cnki.1001-4179.2024.S2.041

0 引言

粗粒土廣泛應用于水利、交通和建筑等領域，其顆粒間常為點接觸，在外力作用下易發生顆粒破碎，進而影響粗粒土的物理力學性質^［1］。其中，伴隨顆粒破碎及其顆粒位置重新調整而產生的粗粒土蠕變現象是工程重點關注的問題^［2-5］，如土耳其的Ataturk堆石壩竣工后壩頂沉降達到2.5 m^［6］；羅馬尼亞Lesu壩，在水庫建成運行2 a后，左壩肩面板與趾板間產生明顯的相對位移，導致周邊縫止水破壞^［7］。現有研究表明粗粒土蠕變機理主要為顆粒破碎，蠕變規律亦與顆粒破碎有直接關系^［8］，顆粒破碎的演化規律直接決定了蠕變規律。粗粒土顆粒破碎的度量一般可由級配來表示，如Marsal等提出的顆粒破碎指數就是依據粗粒土受力前后及最終級配計算得出的^［9-10］。故開展粗粒土的級配演化規律研究有助于深入分析粗粒土蠕變機理與規律。

目前級配演化規律研究方法主要分為室內試驗和基于離散單元法的數值試驗。在室內試驗方面，王瑞甫^［11］利用單軸壓縮蠕變試驗，探討粗粒土中細顆粒含量與蠕變參數之間的關系；孔德志等^［12］對堆石料的顆粒破碎進行了試驗研究，結果表明顆粒破碎與顆粒的級配、大小和受力情況等因素有關；程展林等^［13］通過CT三軸試驗，從細觀角度研究了粗粒土顆粒之間的相互作用以及顆粒的相對位置對其強度及變形性質的影響?？偟膩碚f，上述研究大多基于室內試驗前后的級配分析，而試驗過程中的級配演化規律未曾涉及，原因在于試驗進行過程中級配測試較為困難^［14］。而離散單元法將粗粒土視為離散介質，這與實際材料更加接近，且可以實現連續監測，能獲取試樣蠕變全過程相關參數，經濟迅速并可重復試驗^［15-16］。陳鋮等^［17］采用PFC^3D對粗粒土級配進行模擬，研究了顆粒級配對粗粒土強度與變形特性的影響；徐肖峰等^［18］通過大型直剪試驗和離散元方法研究了粗顆粒含量對粗粒土宏觀及細觀性質的影響機理；王瑤等^［19］基于應力腐蝕程序的離散單元法，討論了接觸模型和應力腐蝕模型細觀參數對蠕變的影響規律。近年來，越來越多的學者采用離散單元法研究粗粒土的變形、強度與顆粒破碎等問題，但對蠕變過程中粗粒土級配演化規律的研究不多^［20］。

本文基于應力腐蝕理論和離散單元法建立粗粒土雙軸蠕變試驗數值模型，通過室內三軸蠕變試驗確定細觀模型參數。在此基礎上，開展不同應力水平粗粒土蠕變數值試驗，分析顆粒破碎模式、級配曲線及其分維值的演化規律。

1 數值模型的建立

1.1 應力腐蝕理論簡介

Potyondy^［21］通過類比巖石力學中的亞臨界裂縫擴展理論，基于膠結鍵拉應力函數實現了膠結鍵膠結半徑隨時間縮減的過程，模擬了巖土材料的裂縫隨時間發展的過程，可用于分析巖土材料的蠕變和靜態疲勞破壞等問題。膠結鍵直徑計算公式為

dd_c/dt=0，/σlt;σ_a

-β₁e^β₂^（σ/σ_t^），/σ_a≤σlt;σ_t

-∞，/σ_t≤σ（1）

式中：d_c為膠結鍵直徑，σ_a為激活應力，σ_t為膠結鍵抗拉強度，σ為膠結鍵最大拉應力，β₁和β₂為材料參數。

基于應力腐蝕理論采用離散單元法模擬粗粒土蠕變的程序為：①對試樣施加伺服偏應力并判定試樣是否破壞，試樣未破壞則執行應力腐蝕程序，試樣破壞則結束蠕變模擬；②在應力腐蝕程序運行時先判斷當前蠕變時長是否達到模擬時長，若未達到則執行應力腐蝕函數，若已達到則結束模擬；③遍歷每一個膠結鍵，判斷是否發生應力腐蝕，若不發生應力腐蝕則結束蠕變模擬，反之計算應力腐蝕速度，縮減膠結鍵半徑并更新蠕變時長；④重復施加應力進行迭代計算，直至滿足蠕變模擬時長、試樣不發生應力腐蝕或試樣破壞3個判定條件之一，結束蠕變模擬。

1.2 試樣建立

1.2.1 顆粒生成

顆粒生成包括級配設計和cluster團粒生成。根據分形理論將土料設計為具有質量自相似分形結構的級配，有

Mdlt;d_i/M_t=d^3-D_i-d^3-D_min/d^3-D_max-d^3-D_min（2）

式中：M（dlt;d_i）為粒徑d小于給定粒徑d_i的累積質量；M_t為顆?？傎|量；d_max和d_min為顆粒最大粒徑和最小粒徑，D為分維值。考慮到尺寸效應影響^［22］，最大粒徑d_max=20 mm約為試樣直徑的1/5，同時考慮到數值模擬的局限性，最小粒徑d_min=1 mm。分維值參考堆石料縮尺結果取D=1.5^［23-24］?，F有研究表明粗粒土級配最終都會趨向某一分形級配，一般取該最終分形級配分維值為2.5^［25-26］。將試樣分為20～10，10～5，5～2，2～1 mm共4個粒組，顆粒級配曲線如圖1所示。

通過典型紅砂巖粗粒土輪廓特征分析，cluster團粒生成時使用了4種模板。為使顆粒破碎后其體積守恒，模板顆粒兩兩相切。每個粒組用一種模板的顆粒代替，每種模板的等效粒徑半徑用面積等效原理計算，即

R=∑n/i=1r²_i（3）

式中：R為大顆粒的面積等效粒徑半徑，n為組成模板的小顆粒個數，r_i為第i個小顆粒的半徑。20～10，10～5，5～2，2～1 mm模板顆粒等效粒徑分別為14.56，9.59，4.0，2.0 mm，顆粒個數分別為60，49，115，102個。

1.2.2 墻體生成

試樣尺寸為D×H=100 mm×200 mm。試樣墻體共包含四面facet墻體，其中上下兩面墻為伺服加載板，分別模擬試樣帽和試樣底座，用于施加伺服偏應力；左右兩面墻為柔性伺服邊界，模擬橡皮膜，用于對試樣施加伺服圍壓。為防止蠕變模擬過程中顆粒發生泄漏，每面墻體的長度放大為試樣尺寸的1.5倍。

1.3 模型參數與模型驗證

本次粗粒土蠕變數值模擬使用接觸模型和應力腐蝕模型。參考王瑤等^［19］對模型參數的敏感性分析，對應力腐蝕模型參數和膠結鍵抗拉強度進行調試，獲取了一套基本符合室內試驗結果的模型參數，見表1。

為驗證模型參數的合理性，進行了圍壓σ₃=300 kPa、偏應力水平S=0.75的室內三軸蠕變試驗，試驗土料為微風化紅砂巖粗粒土，土料級配、試樣尺寸與數模試樣一致。室內試驗與數值模擬所得的軸向蠕應變ε_c與時間t關系曲線如圖2所示。由圖2可知：數值模擬曲線的衰減蠕變時長比室內試驗稍短，后期變形較室內試驗偏小，但大體趨勢一致。基于應力腐蝕理論，采用離散元研究粗粒土級配演化規律的方法是可行的，本文采用的模型參數是合理的。

2 數值模擬結果分析

離散元方法可實現連續監測，獲取不同時間節點的模擬結果，還可用其內嵌語言編寫不同的監測函數得到室內試驗難以獲取的細觀參數變化過程，從而更加直觀地認識粗粒土蠕變過程的級配演化規律。由于不同圍壓下的模擬結果基本相似，本文僅對圍壓300 kPa下的數值試驗結果進行分析。

2.1 顆粒破碎

2.1.1 顆粒破碎模式

顆粒破碎是指顆粒在外部荷載作用下產生結構的破裂或破損，分裂成粒徑相等或不等的多個顆粒的現象^［27］。根據室內試驗顆粒分析結果，可以推斷粗粒土在蠕變過程中存在顆粒破碎，但顆粒破碎模式卻難以直接觀察到。為此，在離散元分析中對顆粒碎片進行監測，以探討粗粒土蠕變過程中的顆粒破碎模式。

不同應力水平數值試樣在蠕變結束后的顆粒破碎模式如圖3所示。在低應力水平下，試樣僅在端部存在少量破碎顆粒，大部分顆粒表面分布有微裂隙，但并未貫穿，顆粒破碎量少；在高應力水平下，試樣不同部位的顆粒同時存在劈裂和棱角破碎現象。

微裂隙的數量在一定程度上可反映顆粒破碎的程度，圖4給出了微裂隙數量隨時間變化的過程。由圖4可知：①微裂隙的數量增長及其表征的顆粒破碎是一個逐漸發展的過程；②高應力邊界條件與低應力相比，微裂隙數量特別是拉伸微裂隙數量，達到穩定值所需時間長。

一定應力條件下，土顆粒裂隙隨時間逐漸發展的過程宏觀上表現為粗粒土的蠕變現象，高應力水平下，達到蠕變穩定所需時間更長。室內試驗結果可得到類似的結論，部分試驗后粗粒土存在明顯宏觀裂縫但并未貫穿，如圖5所示。

2.1.2 顆粒破碎演化規律

在離散元分析中可用斷裂膠結鍵數量與總膠結鍵數量之比（顆粒破碎率B_b）來反映顆粒破碎，不同應力水平下顆粒破碎率隨時間的演變過程如圖6所示。由圖6可知：顆粒破碎率隨應力水平的增大而增大，前期增長迅速，后期增長變緩并逐漸趨于一個穩定值；顆粒破碎率曲線拐點隨應力水平的增大逐漸右移，表明應力越大，顆粒破碎率達到穩定值所需時間越長。

2.2 級配曲線演化規律

按式（3）計算團粒等效粒徑，得到不同時間和不同應力水平下的級配曲線，如圖7所示。圖7中分維值D=2.5的級配曲線為理論上的最終級配曲線。由圖7（a）可知：隨著應力水平的增大，大顆粒含量逐漸減少，平均粒徑變小，顆粒破碎率增大。由圖7（b）可知隨著時間的增長，粗粒土級配曲線逐漸向最終級配曲線趨近，級配在蠕變后期變化較小。該演變規律表明，粗粒土在蠕變過程中，顆粒破碎前期劇烈，后期微弱，這與顆粒破碎率B_b的演變規律是一致的。

2.3 分維值演化規律

級配曲線分析是分析顆粒破碎的主要手段之一，但從級配曲線中難以定量得出其演變規律，故將圖7的數據按照式（2）在Origin中進行非線性曲線擬合，得到不同應力水平和不同時間下的分維值如表2～3所列。不同時間下的分維值演化規律如圖8所示。由表2、表3和圖8可知：擬合曲線的相關系數R²始終在0.99以上，表明粗粒土的級配在蠕變過程中始終具有分形特性，其分維值隨應力水平的增大而增大；分維值隨時間變化，前期增長迅速，但增長速度逐漸變小，1 d后分維值變化很小。

分維值D是粗粒土顆粒大小分布特征的量化指標，分維值的大小表征顆粒發生破碎的程度。相同圍壓相同應力水平下D與顆粒破碎率B_b之間的關系如圖9所示?？芍?，數值模擬中的粗粒土顆粒破碎后，D與B_b存在良好的線性擬合關系，相關系數R²為0.930 3，D與B_b存在一定的內在聯系，這是由于顆粒破碎率B_b的增大使顆粒數目增多，孔隙結構更加合理，因而分維值D隨之增大。

3 結論

本文基于離散單元法模擬了粗粒土雙軸蠕變試驗，并與室內三軸蠕變試驗結果進行比較，驗證了數值模型的合理性，通過數值計算研究了粗粒土蠕變條件下的級配演化規律，主要結論如下。

（1）粗粒土在蠕變過程中，顆粒破碎率隨應力水平的增大而增大，前期增長迅速，后期增長變緩并逐漸趨于一個穩定值；顆粒破碎率曲線拐點隨應力水平的增大逐漸右移，表明應力越大，顆粒破碎率達到穩定值所需時間越長；低應力水平下的試樣僅在端部存在少量破碎顆粒，大部分顆粒表面分布有微裂隙但并未貫穿，高應力水平下試樣不同部位的的顆粒同時存在劈裂和棱角破碎現象。

（2）隨著應力水平的增大，粗粒土大顆粒含量逐漸減少，平均粒徑變小，顆粒破碎增大；隨時間的增長，粗粒土顆粒破碎前期劇烈，后期微弱，級配曲線逐漸向最終級配曲線趨近，其變化規律與顆粒破碎率的演變規律一致。

（3）粗粒土的級配在蠕變過程中始終具有分形特性，其分維值隨應力水平的增大而增大；隨時間增長，分維值的增速逐漸變??；分維值與顆粒破碎率存在良好的線性擬合關系。

參考文獻：

［1］ VARADARAJAN A，SHARMA KG，VENKATACHALAM K，et al.Testing and modeling two rockfill materials［J］.Journal of Geotechnical and Geoenvironmental Engineering，2003，129（3）：206-218.

［2］ 王琛，詹傳妮，何鵬，等.堆石蠕變與顆粒破碎特性三軸試驗［J］.土木建筑與環境工程，2011，33（5）：57-62.

［3］ ZHAO ZH，SONG EX.Particle mechanics modeling of creep behavior of rockfill materials under dry and wet conditions［J］.Computers and Geotechnics，2015（7）：137-146.

［4］ KUWANO R，JARDINE RJ.On measuring creep behaviour in granular materialsthrough triaxial testing［J］.Canadian Geotechnical Journal，2002，39（5）：1061-1074.

［5］ JUSTO JL，DURAND P.Settlement-time behaviour of granular embankments［J］.International Journal for Numerical and Analytical Methods in Geomechanics，2000，24（3）：281-303.

［6］ 邵磊，遲世春，王振興.基于裂縫擴展的堆石料流變細觀模型［J］.巖土工程學報，2013，35（1）：66-75.

［7］ 梁軍，劉漢龍.面板壩堆石料的蠕變試驗研究［J］.巖土工程學報，2002，24（2）：257-259.

［8］ 梁軍，劉漢龍，高玉峰.堆石蠕變機理分析與顆粒破碎特性研究［J］.巖土力學，2003，24（3）：479-483.

［9］ MARSAL RJ.Large-scale testing of rockfills materials［J］.Journal of the Soil Mechanics and Foundation Engineering，ASCE，1967，93（2）：27-44.

［10］HARDIN CS.Crushing of soil particles［J］.Journal of Geotechnical Engineering，1985，111（10）：1177-1192.

［11］王瑞甫.細粒含量對粗粒土蠕變特性影響的試驗研究［J］.公路交通技術，2014，10（5）：1-3.

［12］孔德志，張丙印，孫遜.人工模擬堆石料顆粒破碎應變的三軸試驗研究［J］.巖土工程學報，2009，31（3）：464-469.

［13］程展林，丁紅順，吳良平.粗粒土試驗研究［J］.巖土工程學報，2007，29（8）：1151-1158.

［14］劉萌成，孟鋒，王洋洋.粗粒料顆粒破碎變化規律大型三軸試驗研究［J］.巖土工程學報，2020，42（3）：561-567.

［15］XU YF，FENG XB，ZHU HG，et al.Fractal model for rockfill shear strength based on particle fragmentation［J］.Granular Matter，2015，17（6）：753-761.

［16］蔣明鏡，王富周，朱合華，等.密實散粒體宏微觀特性的直剪試驗離散元數值分析［J］.河海大學學報（自然科學版），2010，38（5）：538-543.

［17］陳鋮，劉小清，羅正東，等.顆粒級配對粗粒土強度與變形特性影響的研究［J］.工程地質學報，2016，24（6）：1191-1198.

［18］徐肖峰，魏厚振，孟慶山，等.粗粒含量對礫類土直剪過程中強度與變形特性影響的離散元模擬研究［J］.工程地質學報，2013，21（2）：311-316.

［19］王瑤，魏進兵，溫仁節，等.粗粒土應力腐蝕蠕變影響因素的離散元分析［J］.水利與建筑工程學報，2020，18（3）：87-92，117.

［20］SEBASTIAN LG，LUIS E，LUIS F.Visualization of crushing evolution in granular materials under compression using DEM［J］.Imternational Journal of Geomechanics，2006，6（3）：195-200.

［21］POTYONDY DO.Simulating stress corrosion with abonded-particle model for rock［J］.International Journal of Rock Mechanics and Mining Sciences，2006，44（5）：677-691.

［22］趙婷婷，周偉，常曉林，等.堆石料縮尺方法的分形特性及縮尺效應研究［J］.巖土力學，2015，36（4）：1093-1101.

［23］譚彩，王琛，吳勇，等.無黏性粗粒土直接剪切試驗的強度尺寸效應［J］.四川大學學報（工程科學版），2016，48（增1）：94-99.

［24］WEN RJ，TAN C，WU Y，et al.Grain size effect on the mechanical behavior of cohesionless coarse-grained soils with the discrete element method［J］.Advances in Civil Engineering，2018，1-6.

［25］TURCOTTE DL.Fractals and fragmentation［J］.Journal of Geophysical Research Solid Earth，1986，91（B2）：1921-1926.

［26］MCDOWELL GR，BOLTON MD，ROBERTSON D.The fractal crushing of granular materials［J］.Journal of the Mechanics and Physics of Solids，1996，44（12）：2079-2101.

［27］杜俊，侯克鵬，梁維，等.粗粒土壓實特性及顆粒破碎分形特征試驗研究［J］.巖土力學，2013，34（增1）：155-161.

（編輯：黃文晉）