滿 騰， 葛 轉， Herbert E.Huppert, Sergio Andres Galindo Torres*

(1.西湖大學 工學院 浙江省海岸帶環(huán)境與資源研究重點實驗室，杭州 310024；2.劍橋大學 國王學院，英國劍橋 CB2 1ST)

1 引 言

對于顆粒材料的瞬態(tài)流變性研究，顆粒材料柱體塌落行為常作為理想的研究對象，其運動規(guī)律雖體現(xiàn)瞬態(tài)動力學特性卻較為簡單，同時顆粒柱的塌落行為也可與更為復雜的滑坡和泥石流等自然顆粒流建立聯(lián)系，方便進一步將小尺度的柱體塌落問題轉化為大尺度的自然地球物理流的研究[12]。文獻[13-15]通過不同顆粒材料(沙粒、大米顆粒和小米顆粒等)塌落試驗研究以及量綱分析總結出顆粒材料柱體塌落擴散程度隨顆粒柱初始高徑比的變化趨勢，認為當柱體的初始高徑比α小于某個臨界值αc時，無量綱擴散比R*=(R∞-Ri)/Ri與初始高徑比α=Hi/Ri成線性關系;而當α>αc時，R*與初始高徑比α0.5成線性關系，其中R∞為柱體的最終塌落范圍的半徑，Ri為顆粒柱截面的平均初始半徑。此相關關系只針對二維顆粒柱塌落以及圓形截面柱體塌落的情況。

Staron等[16,17]在此基礎上研究了顆粒間接觸特性對顆粒柱整體塌落行為的影響，認為顆粒之間的摩擦系數(shù)對柱體的塌落擴散程度有較大的影響。前期研究通過改變顆粒材料的粒間摩擦系數(shù)、顆粒-底面摩擦系數(shù)、柱體系統(tǒng)相對尺寸以及柱體截面形狀等研究了相關參數(shù)對顆粒柱塌落的影響。首先，本文認為顆粒柱的塌落擴散比可以通過等效柱體高徑比αeff進行描述[18]，

(1)

式中μw為顆粒與底面的摩擦系數(shù)，μp為顆粒之間的摩擦系數(shù)，β為常數(shù),用以表征兩種摩擦之間的能耗比。其次，顆粒柱塌落的尺寸效應可以用有限尺寸標度模型(2)表述[19]，

R*=(Ri/d)- β1/νFr[(αeff-αc ∞)(Ri/d)1/ν]

(2)

式中Fr[·]為尺度標度函數(shù)，d為顆粒平均粒徑，β1=0.28±0.04和ν=1.39±0.14為常數(shù)，且與二維隨機網絡逾滲問題的有限尺寸標度模型參數(shù)相似。

本文在上述研究的基礎上總結顆粒材料塌落擴散距離的計算方法，并通過對顆粒材料柱體塌落過程的應力和應變率的分析，研究顆粒材料在瞬態(tài)過程的流變特性。

2 離散單元法及模型的建立

2.1 離散單元法和接觸模型

本文采用擴展多面體離散單元法(DEM)。使用DEM可以方便獲得顆粒材料柱體塌落過程中的整體動力學響應，特別是能夠獲得顆粒之間的相互作用力，便于進行細觀尺度上的力學分析，也有利于從顆粒尺度對問題進行深入研究。多面體單元由若干平面組成，其幾何構成主要為角點、棱邊和平面，在多面體單元的基礎上采用球體和任意多面體的Minkowski Sum來得到擴展多面體[20,21]，這樣顆粒之間相互作用的求解轉變?yōu)閮蓚€多面體的擴展球體之間的相互作用。

擴展多面體DEM的顆粒接觸有角點與角點的接觸、棱邊與棱邊的接觸以及角點與平面的接觸三種形式。在計算過程中，首先根據(jù)顆粒的位置以及其幾何形狀確定兩個顆粒的接觸變形δn，利用線性粘彈性模型求得顆粒之間的接觸面法向相互作用力。而顆粒之間的切向相互作用則需要在每一時間步計算兩個顆粒接觸點位置的切向位移，并確保利用線性模型計算顆粒之間的切向相互作用并且考慮顆粒之間的切向作用不得大于顆粒之間的滑動摩擦力。根據(jù)顆粒之間的相互作用力，可以實時利用牛頓第二定律計算每一個顆粒的加速度、速度和位移，獲取顆粒之間的相互作用力信息。本文DEM計算參數(shù)為正向剛度系數(shù)Kn=1×104N/cm,切向剛度系數(shù)Kt=5×103N/cm,顆粒碰撞回彈比e=0.1。本文采用多物理場仿真程序庫MechSys的DEM模塊進行計算和研究，該數(shù)值計算方法已在文獻[18-20]中得到多次驗證，可以進行顆粒材料柱體塌落的研究。

2.2 模型的建立

本文對兩種不同截面形狀(圓形截面和長方形截面)的顆粒材料柱進行研究，利用MechSys在設計區(qū)域內利用三維Voronoi圖進行分網，建立不規(guī)則擴展多面體顆粒模型，此方法生成的Voronoi顆粒的尺寸遵循高斯分布[20]，在此基礎上隨機選取20%的顆粒刪除，使得初始的堆積柱體具有80%的固體顆粒體積分數(shù)，且孔隙分布較為均勻，此方法可減少顆粒流動過程中的剪脹效應，同時增加系統(tǒng)的隨機性。此后，移除限制顆粒材料運動的邊界條件(如顆粒柱截面為圓形，則此邊界條件為圓管)使顆粒材料在自重條件下自然塌落。待顆粒停止運動，測量其最終擴散距離。對于圓形截面顆粒柱的塌落，本文測量四個方向的擴散距離并求平均，對于其他非圓形截面顆粒柱，則選取特定方向的擴散距離予以測量，如對于長方形截面顆粒柱，僅測量對角線方向的顆粒擴散距離以及垂直于截面棱邊中點方向的平均擴散距離。長方形截面長邊 6 cm，短邊3 cm，其初始高度從1 cm增加至 50 cm。為研究尺寸效應，圓形截面柱的直徑選取Di=1.0 cm,3.0 cm和7.0 cm。同時，顆粒之間的摩擦系數(shù)為0.2,0.4和0.6，而顆粒與底面的摩擦系數(shù)保持為0.4不變。

圖1 不同截面形狀顆粒柱的塌落情況以及塌落后的顆粒數(shù)量沿不同方向的分布情況

3 結果與討論

3.1 顆粒材料柱體塌落及其尺寸效應

首先考察顆粒之間的摩擦系數(shù)對顆粒柱塌落的影響。文獻[18]給出了不同摩擦系數(shù)條件下塌落擴散比R*與初始高徑比α之間的關系，增加顆粒間摩擦系數(shù)可以顯著降低其塌落擴散的距離。對于同一粒間摩擦系數(shù)的情況，雖然塌落擴散比R*隨初始高徑比的增加而單調遞增，但是其在雙對數(shù)坐標系下的斜率在某一臨界高徑比αc前后發(fā)生改變，一般認為這一臨界點表征著顆粒材料具備不同的動力學特性。

當顆粒材料柱的初始高徑比小于臨界高徑比時，顆粒材料處于擬靜態(tài);當大于臨界高徑比時,顆粒材料處于慣性流狀態(tài)或近似液態(tài)。更為重要的是，改變顆粒之間的摩擦系數(shù)不僅能夠改變其塌落擴散距離，也能改變臨界高徑比的大小。根據(jù)前期研究，當橫坐標為有效高徑比αeff時(如式(2))，本文得到的R*塌縮到一條曲線，從而在計算塌落擴散比時得以考慮摩擦系數(shù)的影響。本文進一步考察了兩種具備不同摩擦系數(shù)的顆粒混合在一起時的顆粒柱塌落擴散行為，此時顆粒與底面的摩擦系數(shù)為μw=0.4，顆粒1#之間的摩擦系數(shù)為μp 1=0.1，顆粒2#之間的摩擦系數(shù)為μp 2=0.6，顆粒1#和顆粒2#之間的摩擦系數(shù)通過μ12=2μp 1μp 2/(μp 1+μp 2)計算。圖2給出了兩種不同配比情況的顆粒柱塌落擴散比與初始高徑比的關系。當顆粒2#所占比例從10%增加到50%時，同樣的初始堆積柱體獲得的塌落擴散比顯著增加，因而，即使是顆粒混合的情況，增加高摩擦系數(shù)的顆粒也可以顯著提高顆粒材料整體的摩擦效應。

圖2 兩種不同顆粒配比條件下顆粒柱塌落擴散比R*與初始高徑比α的關系，其中顆粒2#所占比例分別為10%和50%

本文同時考慮系統(tǒng)尺寸的大小對顆粒材料柱塌落擴散距離的影響?？紤]三種不同的圓形截面顆粒柱尺寸，其初始直徑分別為1 cm,3 cm和 7 cm，顆粒尺寸約為0.2 cm，所以相對系統(tǒng)尺寸Ri/d=2.5,7.5和17.5。當Ri/d=2.5時，臨界有效高徑比出現(xiàn)在αeff≈7.5；當Ri/d=7.5時，臨界有效高徑比出現(xiàn)在αeff≈3.0；當Ri/d=17.5時，臨界有效高徑比出現(xiàn)在αeff≈1.7。因而可以認為改變系統(tǒng)的相對尺寸，其臨界高徑比也隨之變化。如圖3所示，增加系統(tǒng)尺寸可以增加顆粒材料柱的相對塌落擴散距離，然而根據(jù)Cabrera等[22]的研究，這種擴散距離隨系統(tǒng)尺寸增加而增加的趨勢隨著相對尺寸Ri/d>50而逐漸消失，前期研究指出，擴散距離隨系統(tǒng)尺寸增加趨勢的消失具備冪律衰減的特性。根據(jù)顆粒柱塌落中表現(xiàn)出的冪律衰減規(guī)律，本文在前期研究中利用有限尺度標度理論去分析模擬數(shù)據(jù)，從而得到如式(3)所示的尺寸效應規(guī)律。其中，ν=1.39±0.14為臨界指數(shù)，αc ∞=0.2為系統(tǒng)尺寸無窮大時的臨界有效初始高徑比。而縱坐標為R*(Ri/d)β / ν，其中，β=0.28±0.04為臨界指數(shù)且與常見材料的實驗測量臨界指數(shù)相似。在考慮有限尺度標度理論的情況下，所有計算結果均分布在一條規(guī)律的曲線上，從而在考慮顆粒間摩擦系數(shù)的基礎上將系統(tǒng)相對尺寸也考慮在內。

圖3 不同相對系統(tǒng)尺寸對R*和αeff關系的影響

3.2 柱體截面形狀的影響

上述研究均以圓形截面顆粒柱的塌落為切入點，然而，非圓形截面顆粒柱往往表現(xiàn)出更為復雜的動力學特性。本文選取截面尺寸為6 cm ×3 cm的長方形截面柱作為研究對象。重點考察垂直于兩種棱邊方向以及沿對角線方向的顆粒柱塌落擴散情況。對于長方形截面柱，垂直于長棱邊方向稱為短棱向，垂直于短棱邊方向稱為長棱向，沿對角線方向仍稱為對角向。

圖4(a)為長方形截面柱DEM計算結果，其中，實心圓代表短棱向，空心圓代表長棱向，而×代表對角向計算結果?？梢姸汤庀蛳鄬τ谄渌较蚓哂懈鼜姷乃鋽U散能力。需要明確的是，正方形截面的棱向以及長方形截面的短棱向均具有最小的初始半徑Ri，因而該結果顯示，較小的初始半徑對塌落擴散有利。

然而，較小初始半徑Ri對塌落擴散結果有利的結論與3.1節(jié)中尺寸效應的結果相反，因而重新定義無量綱的塌落擴散比為

(3)

圖4 長方形截面顆粒柱不同測量方向對和αeff關系的影響

在此基礎上，進一步考慮利用現(xiàn)有尺寸效應工具，即有限尺度標度理論對截面形狀影響進行分析，相關結果如圖5所示?？梢钥闯?，盡管截面形狀改變和測量方向不同，考慮有限尺度標度理論后，DEM計算結果均與式(3)的結果吻合，其尺寸效應的臨界指數(shù)也與式(3)給出的參數(shù)一致。

圖5 考慮尺寸效應后的不同截面形狀對等效塌落擴散比的影響

3.3 顆粒柱塌落與地球物理流

塌落擴散距離這一指標雖然可以表征顆粒材料在自重條件下的輸運能力，但只展現(xiàn)了塌落堆積結果的其中一個維度，難以用于全面描述其堆積形態(tài)，而堆積形態(tài)對研究和預測地球物理流等復雜流體的運動和堆積更為重要。在前期研究中發(fā)現(xiàn)，當顆粒柱的初始高徑比較大時，顆粒材料的塌落具有擬液態(tài)的特性，其塌落更多地體現(xiàn)為上覆顆粒推動下層顆粒塌落及擴散。塌落后，在初始高徑比較小的情況下其堆積形態(tài)呈現(xiàn)為圓臺或圓錐體，但當初始高徑比較大時，顆粒更多地堆積在中心以及靠近外圍的區(qū)域，而使0.5R∞處出現(xiàn)顆粒較少的情況，其總體堆積形態(tài)可總結為如圖6(a)所示的中心開花，環(huán)狀分布的情況。

圖6 正方形截面顆粒柱塌落后形成的高度云圖及網絡中獲取的法屬波利尼西亞Bora-Bora島衛(wèi)星圖

圖6(a)為大初始高徑比正方形截面顆粒材料柱塌落堆積云圖，該算例初始高徑比為20，顆粒間摩擦系數(shù)為0.4，顆粒與地面摩擦系數(shù)為0.4，可以看出，大部分顆粒分布在兩處，一部分在中心處呈圓錐體分布，另一部分呈圓環(huán)狀分布。該分布與常見的火山噴發(fā)形成的島礁的形態(tài)有一定的相似性。如圖6(b)所示，法屬波利尼西亞的Bora-Bora島的形態(tài)與顆粒柱塌落堆積形態(tài)類似，該島嶼中心為一死火山，外圍有高出水面的環(huán)狀島礁，與文獻[19]指出的顆粒柱塌落與火山碎屑流等地球物理流之間的聯(lián)系相似。

3.4 塌落過程中應力和應變率的計算

(i,j=1,2,3) (4)

3.5 塌落過程中的瞬態(tài)流變性

圖7 顆粒柱塌落過程的流變特性

(5)

圖8 顆粒塌落后過程中顆粒速度分布情況以及代表單元分網情況 (紅色邊界區(qū)域的大于0.6(近底板固相區(qū)域))

同時，本文研究了不同相對系統(tǒng)尺寸(Ri/d=15和30)的顆粒柱塌落過程的流變特性，兩組顆粒柱的初始直徑分別為6 cm和12 cm，初始高度分別為12 cm和24 cm，其初始高徑比相同，顆粒平均粒徑均為0.2 cm。由于系統(tǒng)的整體尺寸較大，顆粒柱較高，因此在顆粒柱塌落過程中產生了兩部分固相區(qū)域，如圖9所示，一部分為近底板部分流動速度較小的區(qū)域，該區(qū)域形成圓錐體狀的低速區(qū)域，另一區(qū)域為靠近顆粒柱頂端的具有高運動速度部分，此部分雖然顆粒流速較大，但顆粒之間的相對速度較小，而在這兩部分中間形成高速剪切帶。大系統(tǒng)的這種特性為流變性研究提供了豐富的信息。本部分研究集中于相對系統(tǒng)尺寸對流動起來的顆粒特性的影響，因此選取遠離底板區(qū)域進行統(tǒng)計分析以減少底板邊界條件的影響。

圖9 顆粒塌落后過程中顆粒速度分布情況以及代表單元固相系數(shù)分布(黑色虛線以上為不受底板影響塌落區(qū)域)

4 結 論

本文系統(tǒng)地分析了顆粒材料柱在不同條件下的塌落堆積情況，特別是分析了顆粒柱的相對系統(tǒng)尺寸和截面形狀對塌落擴散距離的影響。研究表明，減少顆粒之間的摩擦系數(shù)以及顆粒與底板直接的摩擦系數(shù)可以有效增加顆粒材料柱塌落擴散的相對距離，借助于αeff可以得到較為普適的描述不同摩擦系數(shù)條件下顆粒材料柱塌落擴散距離的計算方法。同時，發(fā)現(xiàn)增加顆粒材料柱的相對系統(tǒng)尺寸可以增加其塌落堆積的相對擴散距離，因而具備顯著的尺寸效應，對尺寸效應的研究和理解有助于將實驗室尺度的顆粒材料柱塌落的研究推廣至更為宏觀尺度的工程應用。針對尺寸效應，本文利用統(tǒng)計物理中常用的有限尺度標度理論進行分析，不同系統(tǒng)尺寸條件下，因考慮有限尺度標度理論，所有計算結果均塌縮到同一曲線，具有較好的擬合效果。

另外，對不同截面形狀的顆粒材料柱的塌落堆積形態(tài)也做了相應的總結，指出垂直于截面棱邊方向顆粒更易塌落擴散至較遠的距離。同時考慮等效塌落擴散比，不同方向上顆粒的塌落堆積亦符合有限尺度標度理論。前期研究表明，顆粒材料柱的塌落往往可以與真實地球物理流(泥石流和火山碎屑流等)相關聯(lián)，大初始高徑比條件下顆粒材料柱的堆積形態(tài)與某些火山噴發(fā)后物質堆積形態(tài)具有一定程度的相似性，該領域仍需進一步的研究。

本文也從細觀力學的角度探究顆粒材料在瞬態(tài)條件下的流變性能。以往對顆粒材料流變性的研究往往集中于穩(wěn)態(tài)條件下，顆粒柱的塌落過程則與穩(wěn)態(tài)相去甚遠，本文求得顆粒材料柱塌落過程中不同時刻不同位置的計算單元的局部應力-應變率關系，發(fā)現(xiàn)塌落過程中時刻的選取對于應力應變關系的影響并不明顯，然而改變顆粒材料中的摩擦系數(shù)以及系統(tǒng)的尺寸則對于其瞬態(tài)流變性有較大的影響。下一步將重點研究顆粒材料瞬態(tài)流變性量化分析，以期得到可以應用的瞬態(tài)本構關系方程。