復雜形態碎石顆粒的三維離散元模擬及試驗驗證

程曉穎, 喬婷, 秦建敏, 季順迎

(1. 大連理工大學 工業裝備結構分析優化與CAE軟件全國重點實驗室, 遼寧 大連 116024;2. 大連大學 建筑工程學院, 遼寧 大連 116622)

0 引 言

碎石顆粒多由山石風化破碎而成,因其具有較好的強度和硬度,被廣泛應用在交通土建工程中,如壩體工程和碎石樁的填料,修筑碎石路基結構等。碎石顆粒有尖銳棱角,形狀、紋理呈不規則狀態,具有典型的離散性、非連續分布特征,在數值方法上大多采用離散元法進行研究。

離散元法是20世紀70年代CUNDALL等[1]為解決巖石力學問題而提出的。目前,采用離散元法對具有明顯離散特性碎石顆粒的研究已取得較多成果,在模擬碎石的非規則形狀時大部分學者采用基于球單元的黏結和鑲嵌模型。徐旸等[2]通過三維掃描技術構造具有碎石顆粒形態的球黏結模型,分析顆粒形狀對碎石力學性能的影響。ZHAO等[3]采用球鑲嵌模型建立碎石離散元數值模型,并進行休止角和直接剪切試驗的數值研究,但基于球體填充方式構造的顆粒簇模型只能近似模擬碎石顆粒的基本幾何特征,雖然模擬精度會隨著球單元數量的增加而提高,但也會使模型計算效率急劇下降,且顆粒簇模型較難模擬碎石顆粒的尖銳棱角特征。

碎石的棱角特征對其力學性能有重要影響,非常有必要構造更適合描述碎石顆粒非規則形態的離散單元[4-5]。QIAN等[6]采用三維多面體塊對碎石顆粒進行離散元模擬,預測土工格柵加固碎石試樣的強度特性。BIAN等[7]采用多面體單元建立有砟軌道模型,探究動靜載荷下單個顆粒的滾動和滑動。多面體單元可合理模擬碎石顆粒尖銳棱角,且能模擬顆粒間重要的互鎖現象,但以上學者構造的離散多面體單元均為凸多面體,真實碎石顆粒除有凸出的棱角特征外還有一定的凹陷[8],凸多面體單元不能完全模擬,而用凹多面體單元同時模擬碎石顆粒棱角特征和凹陷特性的文獻又較少。

由于散體碎石顆粒間沒有黏結作用,碎石構筑物的穩定性主要依靠顆粒間的咬合互鎖和摩擦作用。休止角試驗常用于探究散體碎石顆粒間的摩擦作用,對碎石顆粒開展休止角試驗,可探究不同粒徑、形狀和級配碎石顆粒對碎石邊坡穩定性的影響[9-10]。此外,休止角還與許多重要的現象有關,包括巖石和邊坡坍塌、雪崩、巖石沉積行為等[11]。用離散元法建立顆粒材料休止角模型時,通過對比休止角試驗與離散元模擬結果,可對碎石顆粒材料進行參數校準,得到離散元模擬時的近似計算參數[12]?？梢?休止角試驗在表達碎石顆粒材料的流動特性和評價碎石構筑物穩定性方面具有較好的優勢。

在碎石構筑物中,碎石集料必須具有足夠大的強度來承受外部載荷,而測定集料強度較簡便、快速的試驗方法就是直接剪切試驗(以下簡稱為直剪試驗)?？紤]到碎石粒徑較大,因此采用大型直剪試驗儀研究非連續、非均質碎石顆粒的剪切性能。通過大型直剪試驗,可探究含水率和碎石含量對土石混合體凍融交界面剪切強度的影響[13],并分析法向應力、剪切速率、碎石粒徑等因素對剪切界面剪切力學特性的影響和接觸面剪切破壞機理[14-15]。另外,在進行室內試驗的同時,還可對碎石顆粒材料進行離散元直剪模擬,通過不同垂直壓力下的直剪試驗模擬,分析試樣在剪切過程中的體積變化和力學行為[16],探究土石混合料骨架結構的力傳遞路徑和混合料變形破壞機制[17-18]?？梢?通過室內大型直剪試驗和對應的離散元模擬,可探究碎石顆粒材料的剪切破壞和細觀特性。

本文為模擬碎石顆粒的棱角特征和凹凸特性,并合理表達各項力學性能,利用三維激光掃描技術獲取碎石顆粒的表面幾何形貌,將基于能量守恒接觸理論的多邊形網格離散元法引入到碎石顆粒離散元模型的構造中,建立具有真實碎石顆粒尖銳棱角和凹凸特性的多面體單元,開展碎石顆粒的休止角和直剪試驗的離散元模擬,結合室內試驗,分析碎石顆粒的流動性和剪切特性。

1 碎石顆粒的離散元構造

碎石顆粒棱角分明且具有凹凸特性,是典型的非規則幾何形貌。非規則顆粒系統在單元排布、運動形態及動力過程等方面與球形顆粒差異很大。為此,采用多面體離散元法構造多面體單元,模擬碎石顆粒的真實幾何形貌,表達其凹凸特性。

1.1 碎石任意形態多面體模型的構造

為模擬碎石顆粒的尖銳棱角和凹凸特性,用三維激光掃描技術獲取碎石顆粒的三維表面幾何形態。圖1為三維激光掃描系統,3D激光掃描儀基于激光測距原理記錄碎石表面的三維坐標信息。碎石顆粒形態各異,形狀不盡相同,為滿足其多樣性,本文掃描200多個不同形狀、粒徑的碎石顆粒以備使用。

圖1 三維激光掃描系統

由真實的碎石顆粒掃描得到的表面三維坐標信息有幾十萬甚至幾百萬個(見圖2),基于該坐標點云直接構造多面體單元并開展離散元模擬,會極大降低數值計算效率,因此有必要對掃描的點云數據進行簡化處理。

圖2 真實碎石顆粒及其對應的掃描點云信息

本文采用曲率刪減的方法不斷減少碎石模型的點云數據量,并對同一碎石模型不同點云數量下所構造顆粒的直徑和球形度進行分析比較,發現點云數量在100個左右時所構造的碎石模型相對真實顆粒的直徑和球形度誤差較小。因此,使用由約100個點云和約200個三角面圍成的多面體模擬碎石,在提高離散元數值計算效率的同時,還保留碎石的主要棱角特征和凹凸特性。

圖3為同一碎石顆粒曲率刪減后得到的不同點云數量的碎石表面點云信息。將最終簡化得到的點云坐標和創建的三角面邊界矩陣導入多面體離散元程序中,構造由三角面圍成的封閉不規則多面體單元,見圖4。

(a)46 657個

圖4 多面體單元

1.2 多面體單元能量守恒接觸模型

多面體離散單元間的接觸模型符合能量守恒接觸理論[19],即假設2個顆粒當前的接觸狀態,包括法向接觸力Fn和接觸力矩Mn,可以由接觸能量函數進行描述。接觸能量函數是關于2個顆粒當前位置和方向的函數,

(1)

(2)

(3)

以法向接觸方向n作為Fn的單位矢量,

(4)

接觸點xc可由下式確定:

(xc-x)×Fn=Mn

(5)

(6)

式中:x為參考點坐標;λ為一個自由參數。

假設2個顆粒接觸時發生體積重疊,1個合適的接觸勢能函數應考慮接觸時重疊區域的基本特征。FENG[20-21]提出關于接觸體積單調遞增的勢能函數:

ω=ω(Vc)

(7)

式中:Vc為接觸區域ΩA的體積。

勢能函數在法線方向的1階偏導為法向接觸力Fn,

Fn=-?xω(Vc)=-ω′(Vc)Sn

(8)

式中:Sn為接觸面積S1的矢量面積,

(9)

接觸面積S1見圖5中左側區域。

圖 5 2個顆粒接觸時重疊區域示意

單位法向接觸向量n和接觸點xc分別定義為:

(10)

(11)

(12)

式中:Gn為接觸面S1的一個向量化“靜態力矩”。

基于接觸體積的線性接觸能量函數為

ω(Vc)=knVc

(13)

其1階導數

ω′(Vc)=kn

(14)

式中:懲罰因子kn表示法向接觸剛度。

Sn由FENG[21]的一個算例計算得到:

(15)

多邊形網格離散元程序中的臨界時間步長由下式計算得到:

(16)

(17)

式中:R為最小顆粒的半徑;G為剪切模量;ρ為顆粒密度;E為彈性模量;υ為泊松比。

2 碎石顆粒休止角的離散元分析及驗證

休止角作為碎石顆粒重要的物理力學特征參數之一,不僅能反映碎石的流動特性,而且是評價碎石構筑物穩定性的重要指標。為探究碎石顆粒材料的流動特性,并驗證結合三維掃描技術所構造的多面體單元模擬碎石顆粒的合理性,建立多面體單元的休止角離散元模型,并進行對應的室內試驗。

2.1 構建休止角離散元模型

碎石等顆粒材料受重力等作用形成的穩定而不坍塌的堆積體,其最大坡角稱為碎石休止角。目前,針對特定散體顆粒材料測量休止角的方法主要有固定漏斗法、旋轉圓筒法、無底圓筒法等。對于粒徑較大的碎石,大多采用無底圓筒法[22-23]測量其休止角。

本文采用無底圓筒法建立碎石顆粒的休止角離散元模型,其中碎石顆粒為結合三維掃描系統構造的多面體單元,其級配曲線見圖6,對應的碎石顆粒見圖7,圖中的數值為對應碎石顆粒直徑。

圖6 多面體單元級配曲線

圖7 級配碎石顆粒試樣,mm

休止角數值模型見圖8(a)。給定2個邊界,一個為水平底面,尺寸為長×寬=2 m×2 m,另一個是直徑250 mm、高1 000 mm的無底圓筒。多面體碎石單元在圓筒內自然堆積,堆積高度為750 mm,此時圓筒內有600個多面體單元。數值試驗開始時,給無底圓筒施加恒定的速度,使其沿豎直方向緩慢向上移動,見圖8(b)。試驗過程中,多面體碎石單元從圓筒底部落下并最終形成穩定的堆積體,見圖8(c)。

(a)t=0

顆粒材料休止角的測量方法一般是直接計算[24],即測量圓錐狀堆積體的底面直徑和高度,再用反正切函數計算其休止角。碎石顆粒粒徑較大,不易像細沙等細小顆粒能堆成標準的圓錐狀堆積體;同時,碎石堆積體的坡面并非直線,其輪廓不規則,沒有明顯的圓錐頂點,運用常規的直接計算方法不能準確獲得碎石休止角,需要用更合適的方法測量碎石休止角。

根據碎石堆積體主視圖(見圖9(b)),提取碎石堆積體外輪廓線進行線性擬合,由直線斜率求得對應位置的休止角,即圖9(a)中A、B、C、D等4點對應的休止角度分別為29.89°、33.32°、29.38°、29.00°,計算得到碎石最終的休止角為30.40°。堆積體輪廓線線性擬合的相關系數R2均大于0.9,證明擬合坡面輪廓線的方法合理,所測碎石休止角有效。碎石顆粒休止角模型的主要計算參數見表1。

(a)碎石堆積體俯視圖

2.2 碎石顆粒室內休止角試驗驗證

為驗證碎石顆粒休止角離散元模擬結果的準確性,開展對應的室內休止角試驗。如圖10所示為自制的測量碎石休止角的試驗裝置及其示意圖,無底圓筒尺寸與離散元模型一致,將圓筒放置在水平底面上,裝入同級配的碎石顆粒。由液壓手動堆高車牽引圓筒使其緩慢、均勻地豎直向上移動,碎石顆粒在重力作用下落在水平底面上,最終形成穩定而不坍塌的碎石堆積體。

表 1 碎石顆粒休止角模型的主要計算參數

(a)試驗裝置

拍攝碎石堆積體試樣的主視圖(見圖11(b)),提取圖像中碎石堆積體的外輪廓線,對其進行線性擬合,以擬合的直線斜率計算對應位置的休止角θ。

單次試驗的測量結果具有一定的隨機性,誤差較大,如手動堆高車的提升速度、碎石顆粒的排布、相機拍攝的角度等都會對休止角的結果產生誤差。為提高碎石休止角試驗結果的準確性,在相同試驗條件下對同一級配碎石進行7次重復試驗并獲得對應休止角,取7次休止角試驗結果的平均值作為最終的碎石休止角,見表2。

表 2 碎石顆粒休止角的試驗結果

由表2可知,重復試驗中同一位置的休止角有一定的波動變化,證明重復試驗的必要性,本文7次重復試驗結果的平均值31.76°作為碎石顆粒的休止角是可信的。

碎石顆粒休止角室內試驗結果與離散元模擬結果對比見圖12,其中離散元模擬的休止角結果為30.4°,與室內試驗結果的相對誤差為4.2%。誤差的原因一方面是圓筒提升速度,在室內試驗中使用液壓手動堆高車手動操作控制,而離散元模擬是給圓筒施加恒定速度;另外,碎石顆粒在圓筒內的堆積排布也存在一定差異。盡管如此,二者的差異仍處于適度范圍內,故由多面體單元建立的模型在適當置信度下可以預測碎石的休止角。

圖12 碎石顆粒休止角室內試驗與離散元模擬結果對比

通過休止角室內試驗,驗證本文所構造的多面體單元不僅可以合理地模擬碎石顆粒,準確表達其流動性,而且為采用多面體離散元法研究碎石等不規則顆粒材料的復雜力學性能提供一定的依據。

3 碎石顆粒直剪離散元模擬及試驗驗證

在巖土工程中,常用直剪試驗測定土體的剪切強度。為進一步驗證結合三維掃描技術所構造的多面體單元在探究碎石力學性能方面的可靠性,建立碎石直剪離散元模型。針對碎石顆粒尺寸大的特點,采用大型直剪儀對碎石顆粒進行直剪離散元模擬和對應的室內直剪試驗驗證,在驗證多面體單元模擬碎石顆粒合理性的同時,探究碎石顆粒的剪切特性。

3.1 碎石顆粒直剪離散元模擬

大型直剪試驗儀包含上、下剪切盒和上蓋板,其中上、下剪切盒的尺寸分別為0.4 m×0.4 m×0.2 m和0.5 m×0.4 m×0.2 m,且在上剪切盒底部增加長、寬分別為0.4 m×0.1 m的擋板,防止剪切時下剪切盒中外露碎石顆粒被擠出,見圖13(a)。剪切盒和上蓋板均為三角面組成的剛性邊界,碎石顆粒為多面體單元。

(a)離散元直剪模型

直剪離散元模型中碎石顆粒的級配與休止角試驗一致,將多面體碎石單元混合均勻后分層放入剪切盒中,使用靜壓法在上蓋板施加壓力,對剪切盒內的碎石進行預壓實。此時,剪切盒中孔隙率為0.384,有1 200個多面體碎石單元。

直剪試驗開始,在上蓋板上施加恒定法向載荷(20、50和100 kPa)的同時,在下剪切盒上施加恒定的速度,使其緩慢水平向左移動。剪切過程中上剪切盒始終保持固定不動,見圖13(b)。當剪切位移達到75 mm即剪應變為15%時,結束試驗。碎石顆粒直剪離散元模型的主要計算參數見表3。

表 3 直剪離散元模型的主要計算參數

3.2 碎石顆粒室內直剪試驗驗證

為驗證碎石直剪離散元模型的合理性,對碎石顆粒進行大型直剪試驗。本文使用的大型直剪儀示意見圖14,剪切盒的上下盒尺寸、碎石級配、施加的法向載荷均與離散元模型保持一致。將級配碎石混合均勻放入剪切盒內壓實,通過控制法向制動器在上蓋板上施加法向載荷,控制水平作動器使下剪切盒水平移動,完成直剪試驗。

圖14 大型直剪儀示意

法向載荷分別為20、50和100 kPa時,室內直剪試驗與離散元模擬的剪切應力-剪切應變對比曲線見圖15,可知在剪切過程中,基于多面體單元建立的直剪離散元模型的剪切應力-剪切應變曲線與室內試驗結果基本吻合。因此,當離散元模型參數選擇合適時,多面體離散元法能較準確地模擬碎石顆粒的剪切特性。

圖15 不同法向載荷條件下剪切應力-剪切應變曲線

當法向載荷分別為20、50和100 kPa時,取剪切應力-剪切應變曲線上的峰值為抗剪強度。其中室內直剪試驗的抗剪強度分別為41.15、73.66和114.14 kPa,離散元模擬的抗剪強度分別為42.81、68.21和121.44 kPa,對應的相對誤差分別為4.03%、7.40%和6.40%,表明多面體碎石單元可以合理模擬碎石顆粒,且在較小誤差下表達其抗剪特性。

由于碎石顆粒的離散特性,顆粒間幾乎無黏結,碎石剪切強度的變化規律不能簡單地用線性的摩爾庫倫強度準則表征。本文采用INDRARATNA等[25]提出的非線性強度準則對碎石強度進行分析,即:

(18)

式中:τp為碎石顆粒的抗剪強度;σn為法向載荷;σc為碎石母巖的單軸抗壓強度,本文取130 MPa[26];a、b均為無量綱的擬合參數。

采用式(18)對碎石顆粒的抗剪強度與法向應力的非線性關系進行擬合時,考慮3個有效散點確?；趦绾瘮禂M合的有效性。正則化后碎石剪切應力比與法向應力比關系的離散元模擬結果及試驗結果見圖16?？梢钥闯?采用冪函數擬合獲得的離散元模擬和試驗的擬合優度R2均大于0.98,說明此函數能較好反映碎石顆粒的非線性剪切強度。

圖16 正則化后碎石法向應力比與剪切應力比關系的離散元模擬結果及試驗結果

4 結 論

本文結合三維激光掃描技術獲取碎石顆粒真實的表面幾何形貌,構造碎石顆粒任意形態多面體單元。為探究該類多面體單元在模擬碎石顆粒力學性能上的有效性,進行了碎石顆粒的休止角和直剪離散元模擬,并將離散元模擬結果與室內試驗的結果進行對比,有以下結論:

(1)多面體單元可以合理地表達真實碎石顆粒的棱角特征和凹凸特性,能較準確地復刻碎石顆粒的自然形態。

(2)本文建立的碎石休止角數值模型能合理模擬碎石的流動特性,并獲得碎石顆粒的休止角,計算結果可以指導碎石的現場施工、運輸和堆積儲存。多面體碎石直剪數值模型的剪切應力-剪切應變曲線和剪切應力比與法向應力比關系同試驗結果能較好吻合,表明多面體離散元模型能有效表達碎石顆粒材料的剪切特性。

本文的研究結果為模擬碎石等非規則顆粒材料提供新的思路和參考,可通過多面體單元建立相關的離散元力學模型,探究非規則顆粒材料的復雜力學性能。