劉媛媛，陳 可，郭一辰，戴北冰

(1.華北科技學(xué)院，北京 東燕郊 065201；2. 太原理工大學(xué)，山西 太原 030024；3.中山大學(xué)，廣東 廣州 510275)

0 引言

顆粒材料在自然界中普遍存在，比如砂石巖土材料、谷物顆粒等等。顆粒材料集合體可以呈現(xiàn)不同于液體及固體材料的特性，比如顆粒堆底面壓力的最大值并非對應(yīng)于顆粒堆的形狀頂點(diǎn)，而是處于形狀頂點(diǎn)的周圍的環(huán)狀區(qū)域，即“應(yīng)力凹陷”現(xiàn)象。對“應(yīng)力凹陷”現(xiàn)象的試驗(yàn)探索和理論解釋具有重要的基礎(chǔ)理論認(rèn)識意義，也一直是顆粒材料相關(guān)的研究熱點(diǎn)之一[1-6]。離散元法(Discrete Element Method)是20世紀(jì)70年代發(fā)展起來的用于計(jì)算散體介質(zhì)系統(tǒng)力學(xué)行為的一種數(shù)值方法，可以合理科學(xué)的模擬顆粒介質(zhì)的微觀力學(xué)特征，近年來成為研究顆粒體系相關(guān)問題的常用方法之一。

本文利用離散元手段模擬顆粒堆的形成過程，對形成的模擬顆粒堆的底部壓力分布與實(shí)際物理模型試驗(yàn)結(jié)果進(jìn)行對比，從細(xì)觀力鏈分布及其它細(xì)觀指標(biāo)的分析來進(jìn)一步解釋壓力凹陷現(xiàn)象底層的原因，有助于更全面深入的理解認(rèn)識顆粒體系的基本力學(xué)特性。

1 離散元分析方法介紹及模型建立

1.1 特點(diǎn)

離散元方法基本假定是塊體單元為理想剛體，各塊體的運(yùn)動只是空間位置的平移和繞形心的轉(zhuǎn)動，其自身的形狀和大小保持不變；塊體之間的接觸視為角-邊接觸或邊-邊接觸；塊體之間的接觸作用力由節(jié)理面的剛度、接觸點(diǎn)的相對位移及有關(guān)的阻尼力確定。離散元方法基本原理是集合中每一個單元都是獨(dú)立的，每個單元都具有相應(yīng)的尺寸、質(zhì)量、轉(zhuǎn)動慣量和接觸參數(shù)等屬性。它以牛頓第二定律和力-位移定理為基礎(chǔ)，對每一個單元首先確定與之接觸的單元，根據(jù)單元之間的重疊量，運(yùn)用力-位移定理計(jì)算單元之間的接觸力，從而得到單元的合力和合力矩，之后用牛頓第二定律確定單元的運(yùn)動規(guī)律，如此循環(huán)計(jì)算，直到系統(tǒng)中所有顆粒都計(jì)算完畢。顆粒流(Particle Flow Code)是目前應(yīng)用比較廣泛的離散元商用軟件之一。本次研究采用PFC2D進(jìn)行顆粒堆試驗(yàn)?zāi)M分析。在PFC計(jì)算中不需要給材料定理宏觀本構(gòu)關(guān)系和對應(yīng)的參數(shù)，這些傳統(tǒng)的力學(xué)特性和參數(shù)通過程序自動獲得，而定義它們的是顆粒和水泥的幾何和力學(xué)參數(shù)，如顆粒級配、剛度、摩擦力、粘結(jié)介質(zhì)強(qiáng)度等細(xì)觀力學(xué)參數(shù)。顆粒流PFC2D有不同的接觸模型，分別用于模擬不同的接觸材料。一般常用的有線性接觸模型、線性接觸黏結(jié)模型、平行黏結(jié)模型等。本文研究模擬對象為砂土，砂土之間的接觸接近于完全剛性變形，故選取線性模型[7]。

1.2 參數(shù)選取

在整個試驗(yàn)中，通過賦予細(xì)觀參數(shù)來影響整個宏觀的模擬合理性，模型中的各參數(shù)指標(biāo)及取值見表1。所模擬的對象為砂土材料，取顆粒半徑為0.6～0.75 cm范圍內(nèi)的隨機(jī)分布[7]。參考相關(guān)研究工作[7,8]，顆粒間滑動摩擦系數(shù)取值0.977，顆粒和墻體之間的滑動摩擦取值0.7，顆粒間切向剛度和法向剛度取值均為1×107N/m，顆粒與墻之間切向剛度和法向剛度取值均為1×108N/m。通過姚玉相等[7]和周國慶等[9]在顆粒模擬中分析，如果僅使用滑動摩擦很難使系統(tǒng)在合理的分析步內(nèi)達(dá)到平衡，故加入阻尼比為0.3來使系統(tǒng)耗能更快，減小迭代次數(shù)。

表1 模型細(xì)觀參數(shù)

1.3 模型建立

本次模型的建立主要分為以下步驟：①以墻構(gòu)件建立下料漏斗邊界，漏斗下落口暫時封閉。邊界多邊形寬為W為50 cm，底部寬W0為15 cm，傾斜角度θ為30°，側(cè)邊長度H為40 cm，下料口距離底板高度H0為40 cm；②在漏斗邊界區(qū)域上部，在0.6～0.75 cm的半徑范圍內(nèi)隨機(jī)生成1176個顆粒，施加重力加速度9.81 m/s2，顆粒落入漏斗內(nèi)并達(dá)到重力場下的初始平衡；③刪除漏斗下方邊界墻，顆粒由漏斗下部出口自由落體；④模擬固定中心點(diǎn)源方式，在底板墻上堆積形成顆粒堆并達(dá)到平衡穩(wěn)定狀態(tài)。

圖1 模型建立示意圖

2 模型計(jì)算結(jié)果分析及討論

2.1 顆粒堆形態(tài)

通過中心點(diǎn)源法堆積并達(dá)到最終體系內(nèi)部靜力平衡后的顆粒堆形態(tài)如圖2(a)，外部形態(tài)與本文作者開展的室內(nèi)物理模型試驗(yàn)[10]圖2(b)較為一致。其中，如圖2(a)所示，模擬顆粒堆下部的底板用Facet墻結(jié)構(gòu)模擬；在裝料漏斗中由低到高的六個不同高度范圍內(nèi)，生成相應(yīng)顆粒，分別命名為LevelOne到LevelSix，用不同的顏色表示?？梢钥闯?，底層區(qū)域生成的顆粒，也是最先自由落體至底面墻的顆粒，多數(shù)處在顆粒堆底部的邊緣部位。生成位置較高，較遲落至底面的顆粒則多數(shù)位于顆粒堆的中部，甚至占據(jù)了底面中心部位?？梢?，由中心固定點(diǎn)源法堆積形成的顆粒堆，后續(xù)降落的顆粒會將前期已落定的顆粒推向兩側(cè)，勢必會導(dǎo)致中心顆粒會受到兩側(cè)顆粒的指向中心方向的水平推力，而兩側(cè)顆粒會受到底面施加的指向中心方向的水平摩檫力。

2.2 內(nèi)部力鏈

顆粒堆內(nèi)部力鏈分布如圖3所示，顆粒之間及顆粒與墻間的接觸力的大小以不同的顏色表示。顆粒堆在自重荷載條件下，總體來說，淺表層區(qū)域顆粒間接觸力普遍較小，而內(nèi)部接觸力較大。但顆粒堆內(nèi)部接觸力分布不均勻，其中最大值出現(xiàn)在堆體中心的兩側(cè)區(qū)域。力鏈呈網(wǎng)狀分布，并明顯存在由中心向兩側(cè)傾向的主力鏈，構(gòu)架出類似于拱形構(gòu)造力鏈形態(tài)。

圖2 顆粒堆形態(tài)

圖3 顆粒堆內(nèi)部接觸力鏈分布圖

2.3 底面壓力分布

提取顆粒與底部墻所有接觸力的數(shù)值及相應(yīng)位置，得到顆粒堆底部法向壓力在徑向上的分布，如圖4。為了得到接觸壓力的相對大小及其分布的相對位置，并和室內(nèi)物理模型試驗(yàn)進(jìn)行對比，將接觸力及位置分別除以γh和堆體半徑進(jìn)行標(biāo)準(zhǔn)化處理，得到標(biāo)準(zhǔn)化后的顆粒堆底面壓力沿徑向的分布如圖5。由于散粒體系內(nèi)接觸力離散性較大，為了了解更有價值的數(shù)據(jù)分布趨勢，對模擬結(jié)果數(shù)據(jù)進(jìn)行相鄰數(shù)據(jù)的移動平均得到趨勢線，如圖5中紅色虛線。模擬數(shù)據(jù)趨勢線與室內(nèi)物理模型試驗(yàn)數(shù)據(jù)進(jìn)行對比，吻合性較好，呈現(xiàn)出顆粒堆中心壓力凹陷特征。

2.4 細(xì)觀結(jié)構(gòu)特征

細(xì)觀結(jié)構(gòu)的特征描述，可采用組構(gòu)張量來定量分析顆粒間接觸方向、顆粒朝向等細(xì)觀指標(biāo)的分布情況[11,12]。由于本次模擬是采用一定尺寸范圍內(nèi)隨機(jī)分布的二維圓盤來模擬顆粒材料，并且研究焦點(diǎn)為顆粒堆內(nèi)部粒間力的傳遞方向，考慮到顆粒堆為左右對稱的形態(tài)，故僅對顆粒堆左側(cè)一半的區(qū)域內(nèi)顆粒接觸力法向方向分布進(jìn)行統(tǒng)計(jì)分析。根據(jù)Rothenburg等[13]提出的用來近似描述以上細(xì)觀指標(biāo)離散分布及各向異性的連續(xù)函數(shù)：

圖4 顆粒堆底部接觸壓力徑向分布圖

圖5 顆粒堆底部壓力標(biāo)準(zhǔn)化值徑向分布圖

Ec(φ)=E0[1+ancos2(φ-φn)]

(1)

式中，Ec(φ)為顆粒接觸法向的分布函數(shù)。其中an反映分布的各向異性程度，而φn對應(yīng)表示各向異性的主方向，E0是各向同性狀態(tài)時的分布函數(shù)，取值為常數(shù)1/2π。對模擬試驗(yàn)的顆粒堆體左側(cè)一半的區(qū)域內(nèi)粒間接觸法向數(shù)據(jù)進(jìn)行統(tǒng)計(jì)分析，并與上文中的連續(xù)函數(shù)對比，可反推得到粒間接觸法向這一指標(biāo)的各向異性特征：an為0.008；φn為30°。在模擬顆粒堆內(nèi)部粒間法向力方向的分布具有一定程度的主方向的定向性，主方向大約偏離豎向30°。這種分布定向性說明顆粒堆內(nèi)部對重力荷載的傳遞會偏離通常認(rèn)為的沿重力的豎向傳遞，從而可一定程度上解釋底部壓力的分布并非與顆粒堆剖面形狀一致，而產(chǎn)生堆頂投影底面部位的壓力凹陷現(xiàn)象。

2.5 討論

從顆粒堆體內(nèi)部力鏈形態(tài)和細(xì)觀結(jié)構(gòu)特征分析，可反映出由中心點(diǎn)源方式形成的顆粒堆，在重力作用下的力的傳遞具有偏離豎向的主力鏈定向特征。這種特征與表示顆粒堆內(nèi)力傳遞方向的Edwards拱圖[8]假定不謀而合。Edwards拱圖是較早的理論假定模型，如圖6所示。在此模型中，進(jìn)行了定性描述：當(dāng)以中心點(diǎn)源方式堆積顆粒時，顆粒所受重力假定為以一系列的層疊的拱式結(jié)構(gòu)沿斜向直線傳遞，會形成底面壓力最大值出現(xiàn)在堆體頂點(diǎn)兩側(cè)的位置，而中心點(diǎn)壓力為零。很顯然，這個Edwards拱圖假定夸大了壓力凹陷的程度，但可以從一定程度上解釋了顆粒堆體底面中心壓力凹陷現(xiàn)象的可能原因，得到了眾多后續(xù)研究的借鑒。

圖6 顆粒堆內(nèi)力傳遞的Edwards拱圖及對應(yīng)的底部壓力分布圖[8]

3 結(jié)論

(1) 利用離散元手段模擬顆粒堆的中心點(diǎn)源方式堆積過程，分析堆體底部壓力和堆體內(nèi)部力鏈的分布特征，得到與室內(nèi)物理模型試驗(yàn)較吻合的堆體底部中心部位壓力凹陷的分布形態(tài)。

(2) 基于堆體內(nèi)部主力鏈分布形態(tài)和粒間接觸力的法向角度各向異性特征分析，結(jié)合Edwards拱圖模型，解釋了顆粒堆體底部中心壓力凹陷的原因來自顆粒集合體內(nèi)部力的斜向傳遞特征，是顆粒堆體中土拱效應(yīng)的體現(xiàn)。