滑坡碎屑流顆粒分選效應的數值模擬

李天話 樊曉一 姜元俊

(1.西南科技大學土木工程與建筑學院 四川綿陽 621010； 2.中國科學院成都山地災害與環(huán)境研究所 四川成都 610041)

由不同形狀、大小、密度的顆粒構成的顆粒物質系統在運動過程中會出現顆粒分選效應，如受到制藥、運輸、冶金等工業(yè)廣泛關注的“巴西果”分離、“反巴西果”分離和“三明治”分離現象[1]。已有研究表明[2-3]，顆粒分選效應同樣出現在滑坡碎屑流的運動過程之中，具體表現為：粗大粒徑塊石易向滑體前緣及表層聚集，而細小碎屑主要分布于滑體中下部，該效應最終導致了滑坡碎屑流反粒序堆積結構[4]的形成。

滑坡碎屑流顆粒物質的分選效應作為滑坡碎屑流動力學的主要表現特征之一，對沿程溝道床面物質的“鏟刮”效應、滑體運動距離、影響范圍、動能及沖擊作用都有極重要的影響[3,5-7]，研究滑坡碎屑流顆粒分選效應和機理，對了解滑坡碎屑流運移堆積特征及滑坡沖擊災害的規(guī)劃防治工作具有重要的科學意義。

盡管早期學者利用模型試驗手段，通過慣性顆粒流理論、隨機振動篩機理及擠壓驅除機制等理論分析方法[7-9]，闡明了顆粒分選的形成機理，但各理論對顆粒流的密實度、離散度、黏性、級配等顆粒參數依賴較大，難以形成普適性的理論分析結果。

近年來，關于顆粒分選的研究多以現場調查[4,10]和模型試驗[2,6,11]為主。雖然通過現場調查對碎屑流流動堆積區(qū)的粒組組成進行定量分析能夠揭示滑體反粒序堆積特征，通過滑槽模型試驗手段能揭示大、小顆粒含量、粒徑差等因素對分選程度的影響，然而兩種方法均難以實現對顆粒分選形成過程的直觀分析。本文以模型試驗[12]的資料和數據為依據，運用三維顆粒離散元素法( Discrete Element Method, DEM) 分析軟件EDEM建立滑坡碎屑流顆粒流數值模型，分析顆粒流不同粒徑組分在運動過程中分布位置及動能的演化過程，初步揭示滑坡碎屑流顆粒分選效應的形成過程，進而深入分析顆粒分選對滑坡碎屑流動力學特征的影響。

1 研究方法

1.1 離散元素法

離散元素法是基于粒子間接觸的位置關系以及相互作用的力學關系來求得粒子的運動狀態(tài)的一種數值計算方法[13]，該方法由Peter Cundall于1971年首次提出并應用于巖土力學的研究，經過40多年的不斷深入與發(fā)展，DEM已被廣泛運用于滑坡碎屑流運移沖擊作用機制及承災體動力響應等問題的相關研究，對于模擬大變形或大位移條件下的滑坡碎屑流運動具有較好的適用性[14]。

離散元素法利用牛頓第二定律和力-位移法則確定顆粒的運動及受力情況，其核心是顆粒接觸特性，即接觸本構模型[15]。接觸模型的分析計算直接決定了顆粒所受的力和力矩的大小，對不同的仿真對象，需建立不同的接觸模型，以提高仿真結果的準確性[16]。由于室內模型試驗材料為無粘連的類球狀碎石，具有散粒物料特性，顆粒表面黏附力很小，因此本研究采用Hertz-Mindlin (no slip)接觸模型。

1.2 模型尺寸

如圖1所示，以滑槽模型試驗的模型尺寸為依據，建立數值模型?；鄣装彘L2.19 m，寬0.3 m，側板高0.35 m，用于重現碎屑流的流動環(huán)境?；鄣膬A角為35°；模型料槽長0.44 m，高0.2 m，寬0.3 m，用于裝填滑體物料，碎屑流的啟動通過一個觸發(fā)門來控制；斜滑槽底端為剛性擋墻，攔擋滑體顆粒，便于直觀獲取滑體在斜滑槽上的堆積體形態(tài)，分析不同粒組的顆粒分選及運移堆積過程。

圖1 滑槽模型示意圖Fig.1 Sketch map of the flume model

1.3 數值模型參數的確定

離散元數值模擬中，材料的泊松比、剪切模量和密度等本征參數表示各種材料的自身屬性，與外界條件無關?；谀Ｐ驮囼瀸崪y及計算結果，確定顆粒與幾何體的材料屬性如表1所示。

表1 材料屬性表Table 1 Material properties

在離散元素法中，材料宏觀的運動行為取決于顆粒間接觸的微觀力學參數，然而，物理模型試驗中滑體顆粒的運動是多因素、多參數共同作用的復雜力學過程，數值模擬不能完全復原實際模型試驗中的所有影響因素，因此，將數值模擬軟件中有限的微觀力學參數進行實測標定會導致得到的模擬運動過程與試驗原型產生較大差異。對數值模擬材料參數的選擇尚未形成一個完善成熟的理論方法，較多的方式是通過物理模型和數值模型結果的對比驗證來確定相關力學參數。

圖 2給出了不同時刻模型試驗與數值模擬所得滑體形態(tài)的對比結果。圖2表明，離散元數值模擬結果能夠較好地與物理模型試驗相吻合，進而表明表2所示的材料間的接觸屬性參數可用于進行下一步數值模擬研究，以彌補物理模型試驗在儀器設備以及技術手段等方面的局限性，實現對顆粒分選形成過程的直觀監(jiān)測及分析。

1.4 試驗設計

混合粒徑顆粒流在運移堆積過程中會產生摻混和分選效應，導致滑體流動性增強[3,6-7]，運動堆積體呈現“上粗下細”的反粒序分布特征[2,4,11]。為了分析滑坡碎屑流顆粒分選過程中不同粒徑滑體顆粒的動能及沿垂直于滑動方向(垂向)的分布位置演化等運移堆積特性，建立了如圖3所示的M1，M2，M33組初始滑體。

2.1.3 使用者評價。在交流過程中，李主任表示，截至目前，該產品每天都在運行，情況總體良好，煤種和水分適應性不錯，縮分精度、樣品重量、樣品各級粒度等各方面指標均達標。他特別指出，三德科技的該自動制樣系統有幾項關鍵技術解決了制樣過程中的難題。①全通破碎機，該設備在破碎過程中采用了傘旋自清掃技術和給料皮帶聯動技術，解決了6mm濕煤破碎沾堵的難題，達到濕煤全通的效果；②風透低溫快速在線干燥技術，有效解決了水分對制樣環(huán)節(jié)帶來的影響，保證制樣過程不粘不堵；③高性能制粉技術，該自動制樣系統制粉環(huán)節(jié)采用德國技術，制粉收集率達到95%以上，過篩率達到100%，遠高于同類產品。

圖3 3組滑體模型初始堆積狀態(tài)Fig.3 Initial deposition forms of 3 samples

如圖3及表3所示，初始滑體M1，M2和M3均由KL1至KL44個粒組組成，其中，以紅色標記的KL1表示粒徑為2～5 mm的滑體顆粒，以綠色標記的KL2表示粒徑為5～10 mm的顆粒，以黃色標記的KL3表示粒徑為10～20 mm的顆粒，以藍色標記的KL4表示粒徑為20～40 mm的顆粒。3組初始滑體M1,M2和M3中，4個粒組的質量百分比均相同，而各粒組的初始堆積位置不同。M1初始堆積體呈現反粒序堆積特征，具體表現為：堆積體沿垂直于滑槽方向，自下向上粒徑依次增大(圖3(a))；M3初始堆積體呈正粒序分布[17]特征：堆積體顆粒從下向上粒徑依次減小(圖3(c))；M2初始滑體中，4個粒組呈混雜排列特征(圖3(b)及圖4)，然而粗細顆粒在料槽中沉淀靜置的過程中，受重力作用影響，小粒徑顆粒會透過大顆粒之間的空隙，滲透到滑體中下部，導致初始滑體難以保持均勻的隨機混雜狀態(tài)。盡管如此，作為M1和M3之間的過渡模型，M2能夠體現各粒組隨機摻雜分布的特性。

表3 初始滑體級配組成表Table 3 Particle size distributions of the samples

圖4 滑體M2各粒組初始堆積狀態(tài)Fig.4 Initial particle deposition of M2

2 滑坡碎屑流顆粒分選效應

2.1 滑體顆粒垂向分布位置演化特征

國內外高速滑坡堆積體中已發(fā)現了大量的反粒序地質證據[2,4,11]。為分析顆粒分選效應影響下滑體反粒序分布結構的形成過程，研究滑體顆粒沿垂直于滑槽方向的分布位置演化特征，引入特征指標

(1)

式中，垂向分布位置指標HKLi(i=1,2,3,4)表示KL1-KL44個粒組距滑槽底板的平均高度，hj表示第j個顆粒距離滑槽底板的垂向高度(如圖5)，n表示該粒組滑體顆粒的總數量。

圖5 顆粒垂向高度示意圖Fig.5 Sketch map of the vertical distribution of particles

圖6給出了M1,M2和M33組初始滑體各粒組垂向分布位置指標HKLi隨滑體運動時間的演化曲線，KL1-KL44條曲線均呈現先下降后上升的變化趨勢，表達了初始滑體從開始釋放沿斜滑槽下滑，到滑體長度不斷拉長而深度不斷減小，最后受到坡腳擋墻阻擋而壅高停積的運動過程。圖7給出了M1,M2和M33組初始滑體在不同時刻的運動快照。

圖6(a)及圖7(a)表明，初始滑體M1呈反粒序分布，隨著運動觸發(fā)，滑體沿斜滑槽向下運動，位于初始滑體最上部的KL4率先運動到滑體前緣部位，且粗大顆粒間碰撞作用顯著，導致滑體前緣密實度較低。隨著運動的持續(xù)發(fā)展，滑體長度沿運動方向延伸，滑體平均厚度逐漸減小，各粒組滑體垂向分布高度逐漸達到最小值，伴隨著擋墻處顆粒體的持續(xù)堆積壅高，滑體厚度逐漸增加，各粒組垂向分布高度增大，最終趨于穩(wěn)定，滑體運動停止。在整個運動過程中，圖6(a)所示4條曲線未出現相交，表明滑體反粒序分布結構顯著存在于整個運動過程之中。

圖6(c)及圖7(c)表明，初始滑體M3呈正粒序分布，隨著運動觸發(fā)，滑體沿斜滑槽向下運動，位于初始滑體最上部的KL1率先運動到滑體前緣部位，且細小顆粒比表面積大，顆粒間摩擦作用顯著，受表層及前緣細小顆粒的阻擋覆蓋，滑體整體密實度較高。隨著運動的持續(xù)發(fā)展，滑體長度沿運動方向延伸，滑體平均厚度逐漸減小，粗顆粒穿越表層細顆?；w，開始向滑體中上部運移，此時各粒組間摻雜均勻，滑體正粒序分布特征被打破，并開始向反粒序分布過渡。位于滑體前緣的細粒組分受約束相對較小，率先在擋墻處堆積，慣性作用更為顯著的大顆粒比表面積較小，摩擦耗能相對較少，在碰撞分離過程中保持了較大能量，翻越細粒堆積體后覆蓋停積，最終堆積體呈顯著的反序分布特征。圖6(c)所示4條曲線在0.79～1.03 s時間段內相交，表達了在顆粒分選效應的作用下，滑體正粒序分布向反粒序分布的過渡。

圖6(b)及圖7(b)表明,作為M1與M3的過渡模型，M2滑體初始狀態(tài)各粒組間摻雜相對較均勻，但粗顆粒之間孔隙大，接觸相對較松散，堆積高度相對較高。此外，沉淀靜置過程中，細顆粒難以避免地穿過粗顆粒間的孔隙向下運移，導致較粗粒組KL3和KL4的平均分布位置位于初始滑體中上部，較細粒組KL1和KL2平均分布位置位于滑體中下部。運動觸發(fā)后，反粒序分布特征進一步發(fā)展，隨著滑體在擋墻處的逐漸停積，堆積體高度持續(xù)增加，顯著的反粒序堆積特征形成。圖6(b)顯示，初始滑體M2中KL1與KL2間的平均分布高度差為2.00 mm，KL3與KL4間的平均分布高度差為26.46 mm。最終堆積體中，KL1與KL2間的平均分布高度差為27.48 mm，KL3與KL4間的平均分布高度差為46.25 mm。表明對于粗細顆粒摻雜的初始滑體，隨著滑坡碎屑流運動過程中的顆粒分選作用，將形成顯著的反粒序分布特征。

圖6 各粒組運動過程中的分布高度Fig.6 Location-time curves of vertical distribution of KL1-KL4

圖7 不同時刻滑體運動狀態(tài)快照Fig. 7 Motion state of different particle size components of a mobilized particle flow at different times

圖8 最終堆積狀態(tài)Fig.8 Forms of final depositions

2.2 滑體動能演化特征

滑坡發(fā)生時，滑體從高勢能位置向低勢能位置運動并伴隨著巨大的能量釋放。滑體顆粒間的碰撞及摩擦等動力特征可通過滑體動能的增減表達。受顆粒分選效應的影響，同樣規(guī)模的滑體沿相同坡度的溝道運動，其動能差別很大，對受災建筑結構的沖擊破壞作用亦有很大的差異[7]。圖9給出了M1-M33組初始滑體各粒組動能隨運動時間的變化曲線。

圖9 各粒徑組分運動過程中的動能時程曲線Fig.9 Kinetic energy-time curve of KL1-KL4

如圖9(a)所示，對于呈反粒序分布的初始滑體M1，各粒組最大動能與粒徑及始滑位置距底板距離呈正相關，粒徑越大，始滑位置距底板距離越大，動能峰值越大，且獲得最大動能所需時間越短。

如圖9(b)所示，對于粗細顆?；祀s的初始滑體M2，各粒組最大動能依然與粒徑呈正相關，粒徑越大，動能峰值越大。

如圖9(c)所示，對于呈正粒序分布的初始滑體M3，各粒組最大動能不再與其始滑高度或粒徑呈正相關。隨著運動的觸發(fā)，位于初始滑體最表層的KL1率先達到動能峰值，位于初始滑體最底層的KL4最后達到動能峰值，但底層初始滑體KL4獲得的最大動能(19.44 J)略大于表層初始滑體KL1(19.15 J)。位于初始滑體中下部的KL3的動能峰值為18.16 J，略大于中上部初始滑體KL2(17.72 J)。

圖9表明，滑體所得動能峰值受滑體垂向分布位置和粒徑共同影響。位于滑體表層的顆粒所受約束相對較小，能量損失相對較少，便于動能發(fā)育；大粒徑滑體顆粒慣性作用顯著，在頻繁的碰撞作用下，促進了顆粒間的能量傳遞[18]，且比表面積小，摩擦耗能小，動能得以充分發(fā)育。

2.3 滑坡碎屑流顆粒分選機制

(2)

式中，v表示顆粒的運動速度，h0表示該粒組初始堆積的質心距滑槽底板的距離，h表示顆粒沿垂直于底板方向的運動高度，g表示重力加速度，n表示顆粒的數量。

M1和M3兩組初始滑體均由4個粒組組成，每個粒組沿垂直于滑槽底板方向的分布層次不同，對應著不同的初始堆積高度。圖10給出了M1和M3兩組初始滑體對應分布層次的粒組的時程曲線對比情況。

圖10 各粒組時程曲線Fig.10 Comparisons of kinetic energy-time curves between different particle size components

自然條件下，受地質環(huán)境條件、坡體結構及滑體巖性等因素的影響[19]，不同滑坡碎屑流運動過程中，巖體破碎程度不同，粗大塊石與細小碎屑混雜。顆粒流數值模擬結果表明，滑坡碎屑流運動過程中，大粒徑滑體顆粒慣性作用顯著，顆粒間接觸以碰撞作用為主，碰撞過程中滑體密實度降低，小粒徑顆粒易穿過滑體間的縫隙向下部運移，此外，小粒徑顆粒比表面積大，摩擦耗能顯著，在上部滑體的擠壓作用下，動能發(fā)育受到進一步抑制，導致這些顆粒向滑體中后部聚集。而大顆粒比表面積小，碰撞作用顯著，動能得以充分發(fā)育[18]，運動過程中越過下部細小顆粒，向滑體前緣及表層聚集，前緣及表層大粒徑顆粒所受約束相對較小，在頻繁的碰撞過程中進一步促進了大顆粒的飛躍分離(如圖11)。因此，滑坡碎屑流顆粒分選效應具體表現為，粗大塊石運動速度較快并向滑體表層及前緣運移，細小碎屑運動速度較慢，并向滑體底層及后部聚集?；滤樾剂黝w粒分選效應是由振動篩分、小顆粒耗能較大率先停積、大顆粒碰撞分離等動力學機理共同作用的結果。

圖11 滑坡碎屑流大、小粒徑組分運動狀態(tài)示意圖Fig.11 Different trajectories of large and small particles

4 結論

本文以模型試驗的資料和數據為依據，運用三維顆粒離散元素法分析軟件EDEM建立滑坡碎屑流顆粒流數值模型，從滑體不同粒徑組分垂向分布位置和動能兩方面研究了滑坡碎屑流顆粒分選效應的形成過程，探討了顆粒分選效應的形成機制，得出以下主要結論：(1) 滑坡碎屑流顆粒分選效應會導致滑體顆粒向反粒序分布結構發(fā)育。呈反序分布的初始滑體在整個運動過程中會保持反序分布結構，粗細混雜或者呈正序分布的初始滑體會在顆粒分選效應的作用下形成反粒序分布結構。(2)滑體所得動能峰值受始滑深度和粒徑共同影響，始滑深度越淺，顆粒粒徑越大，動能峰值越大。(3) 初始堆積高度一定的情況下，粗大塊石在運動過程中易獲得較大的運動速度并運移至滑體表層及前緣，而細小碎屑運動速度較慢，并向滑體底層及后部聚集。(4) 滑坡碎屑流顆粒分選效應是由振動篩分、小顆粒耗能較大率先停積、大顆粒碰撞分離等動力學機理共同作用的結果。