粗粒料顆粒破碎數(shù)值模擬研究

楊 貴 ，許建寶 ，劉昆林

（1.河海大學(xué) 安全與防災(zāi)工程研究所，江蘇 南京 210098；2.河海大學(xué) 巖土力學(xué)與堤壩工程教育部重點(diǎn)試驗(yàn)室，江蘇 南京 210098）

1 引 言

堆石料作為典型的顆粒摩擦材料，具有顯著的顆粒破碎特性，且隨著堆石料顆粒尺寸的增大，顆粒破碎特征越明顯。顆粒破碎使得堆石料體積收縮并成為影響面板壩應(yīng)力、變形的決定性因素[1]。忽略顆粒破碎特性，將難以透徹理解和描述堆石料力學(xué)特性，是當(dāng)前高土石壩本構(gòu)模型中的研究熱點(diǎn)。

顆粒破碎的機(jī)制非常復(fù)雜，眾多學(xué)者從室內(nèi)試驗(yàn)，數(shù)值模擬和理論分析等手段進(jìn)行了廣泛研究[2-6]。一般根據(jù)試驗(yàn)前后顆粒級(jí)配的變化進(jìn)行顆粒破碎分析，忽略了試樣制備和固結(jié)過程對(duì)顆粒破碎的影響。高玉峰等[3]研究結(jié)果表明，試樣在制備過程中發(fā)生的顆粒破碎不可忽視，為了在固結(jié)開始前保持初始設(shè)計(jì)的級(jí)配，應(yīng)將設(shè)計(jì)級(jí)配曲線相應(yīng)下調(diào)6%來制備試樣。對(duì)試樣在剪切過程中顆粒破碎特征未能進(jìn)行廣泛研究，主要是由于受試驗(yàn)條件的限制。

顆粒流方法作為一種能夠有效模擬粗粒料力學(xué)特性的方法，在巖土工程的眾多領(lǐng)域得到廣泛應(yīng)用[7-9]。本文基于二維顆粒流程序，建立粗粒料顆粒破碎模型，開展堆石料雙軸剪切試驗(yàn)數(shù)值模擬研究，分析粗粒料在剪切過程中顆粒破碎的變化特征。

2 考慮粗粒料顆粒破碎的數(shù)值模型

數(shù)值模擬過程中關(guān)于顆粒破碎的模擬，目前主要有兩種方法：一種是采用組合顆粒來替代單一的模型顆粒，顆粒與顆粒之間通過接觸粘結(jié)和并行粘結(jié)進(jìn)行聯(lián)系，在外力荷載作用下當(dāng)顆粒之間的粘結(jié)力小于破碎力時(shí)，顆粒發(fā)生破碎[10]；另一種是采用離散方法，即當(dāng)顆粒的受力狀態(tài)滿足一定的破壞標(biāo)準(zhǔn)時(shí)顆粒發(fā)生破碎，采用按照一定形狀排列組合的顆粒代替原來的顆粒[9,11]。

兩種方法各有優(yōu)缺點(diǎn)，前一種方法能夠方便地模擬粗顆粒材料的形狀，缺點(diǎn)是顆粒之間的粘結(jié)力確定比較復(fù)雜，在外力作用下，顆粒不能進(jìn)行二次破碎，破碎后顆粒級(jí)配曲線很難獲得，無法與實(shí)際試驗(yàn)結(jié)果進(jìn)行對(duì)比驗(yàn)證，計(jì)算過程中因預(yù)先進(jìn)行了顆粒替換，試樣顆粒數(shù)目較多，計(jì)算效率較低。后一種方法的缺點(diǎn)在于不能模擬顆粒材料的復(fù)雜形狀，其優(yōu)點(diǎn)在于能夠?qū)崿F(xiàn)顆粒的二次破碎，只有當(dāng)顆粒的應(yīng)力達(dá)到其破壞標(biāo)準(zhǔn)時(shí)才發(fā)生破壞，試樣顆粒數(shù)目較少，計(jì)算效率較高，能夠方便地統(tǒng)計(jì)試驗(yàn)過程中顆粒級(jí)配，分析其變化規(guī)律并與室內(nèi)試驗(yàn)結(jié)果進(jìn)行對(duì)比。

單顆粒的破碎過程一般采用單向抗壓試驗(yàn)和巴西圓盤劈裂試驗(yàn)來進(jìn)行模擬。認(rèn)為顆粒的破碎形式是拉裂破壞或剪切破壞。研究過程中，一般認(rèn)為，顆粒的實(shí)際受力狀態(tài)可以簡(jiǎn)化為靜水壓力和相互垂直的主應(yīng)力作用，如圖1 所示。

圖1 顆粒受力簡(jiǎn)化圖Fig.1 Schematic of particle stress state

部分學(xué)者認(rèn)為，顆粒在靜水壓力作用下不會(huì)發(fā)生破壞，只有在主應(yīng)力作用下才會(huì)發(fā)生破壞，且破壞形式為拉裂破壞[9]。在實(shí)際過程中，顆粒在外荷載作用下，周邊接觸顆粒的數(shù)目是有限的，不能完全等效為靜水壓力的作用，因此，有學(xué)者提出，顆粒在等向應(yīng)力荷載作用下也會(huì)發(fā)生破壞，且破壞類型為剪切破壞[11]。

對(duì)于顆粒的拉破壞，采用Oliver 等[9]根據(jù)圓盤劈裂原理建立的拉破壞準(zhǔn)則，即：

式中：σcritt為顆粒極限抗拉強(qiáng)度；σmax為顆粒最大主應(yīng)力；σmin為顆粒最小主應(yīng)力。

式中：σxx為顆粒x 向正應(yīng)力；σyy為顆粒y 向正應(yīng)力；σxy顆粒xy 方向的剪應(yīng)力。

式中：σt50為直徑為50 mm 顆粒的抗拉強(qiáng)度；d為顆粒直徑（mm）。

對(duì)于壓破碎的顆粒，當(dāng)顆粒所受的應(yīng)力狀態(tài)滿足下列公式時(shí)，則發(fā)生破壞。

式中：N為顆粒周邊的接觸數(shù)目；σi為顆粒周邊某一接觸的法向應(yīng)力。

考慮到不同直徑顆粒其抗壓強(qiáng)度不同，結(jié)合Hoek 等對(duì)不同尺寸巖石試樣的研究成果，以及Ben-Nun 等[11]數(shù)值模擬研究成果，定義顆粒的極限抗壓強(qiáng)度 σcritc為

式中：σc為顆粒單軸抗壓強(qiáng)度；fD為顆粒直徑影響系數(shù)；fCD為顆粒接觸影響系數(shù)；σc50為直徑為50 mm 顆粒的單軸抗壓強(qiáng)度；D為當(dāng)前顆粒周邊接觸顆粒平均直徑，當(dāng)接觸對(duì)象為墻體時(shí)，D為無窮大，此時(shí) fD=1 ；d為顆粒直徑；N為顆粒接觸數(shù)目。

數(shù)值模擬過程中對(duì)顆粒所受的應(yīng)力狀態(tài)進(jìn)行分析，當(dāng)顆粒間的接觸數(shù)目 ＜3 時(shí)，顆粒只發(fā)生拉破碎；當(dāng)顆粒數(shù)目 ≥3 時(shí)，顆粒發(fā)生壓破碎或拉破碎。當(dāng)顆粒滿足破碎強(qiáng)度準(zhǔn)則時(shí)，采用如圖2 所示的顆粒進(jìn)行替換，替換后顆粒的材料屬性與替換前一致。各顆粒的大小和位置關(guān)系如圖2 所示。組合顆粒的主軸方向φ 根據(jù)原顆粒所受的最大壓力方向確定。

圖2 單顆粒破碎組合關(guān)系Fig.2 Schematic of particle composition

從圖2 可以看出，采用組合顆粒后，會(huì)導(dǎo)致部分面積缺少。為了彌補(bǔ)缺失的面積，部分學(xué)者采用了面積不變的原則，也叫膨脹法，即在生成替換顆粒的過程中，將顆粒的半徑進(jìn)行放大。本文采用增加新增顆粒密度的方法進(jìn)行控制，即密度法，從而使得替換顆粒前后的總質(zhì)量不變。

對(duì)于拉破碎情況，設(shè)原顆粒的密度為ρ，顆粒半徑為R0，替換后的顆粒密度為ρ′，顆粒半徑分別為R1、R2和R3，則有：

則有：ρ′=1.2ρ

壓破碎的密度處理方法與拉破碎一致?？紤]到實(shí)際工程中，當(dāng)外荷載一定時(shí)，顆粒不能無限制破碎下去，且為了提高數(shù)值模擬的計(jì)算效率，定義最小不破碎顆粒直徑dmin，當(dāng)顆粒直徑大于該值時(shí)才會(huì)發(fā)生破碎，否則不發(fā)生破碎。

3 數(shù)值模擬結(jié)果驗(yàn)證

為進(jìn)行數(shù)值模型的驗(yàn)證，首先需要能夠模擬真實(shí)的室內(nèi)三軸試驗(yàn)成果。基于此，選取江蘇宜興抽水蓄能電站筑壩堆石料大型三軸剪切試驗(yàn)結(jié)果作為本次數(shù)值模擬的考題。試樣干密度為2.12 g/cm3，相對(duì)密度為0.95，最大粒徑為60 mm，不均勻系數(shù)為52.5，曲率系數(shù)為1.07。試驗(yàn)級(jí)配曲線如圖3 所示，試樣尺寸直徑為300 mm，高度為600 mm。

圖3 試驗(yàn)和數(shù)值模擬級(jí)配曲線Fig.3 Testing and simulation gradation curves

數(shù)值模擬過程中，采用如圖3 所示的試驗(yàn)級(jí)配進(jìn)行模擬，為提高計(jì)算效率限制試樣生成過程中最小顆粒直徑，并將試樣尺寸和試驗(yàn)級(jí)配同時(shí)放大10倍，生成后的試樣顆?？倲?shù)為2 473，已有研究表明，當(dāng)試樣尺寸的直徑和高度比保持不變，尺寸變化而帶來的差異可以通過細(xì)觀參數(shù)的變化來反映[12]。本文采用的顆粒接觸模型為非線性H-Z 模型，細(xì)觀參數(shù)如表1 所示。

表1 顆粒流數(shù)值模擬的參數(shù)Table 1 Parameters used in PFC2D simulation

表1為根據(jù)堆石料室內(nèi)試驗(yàn)600 kPa 圍壓下的結(jié)果進(jìn)行模型參數(shù)標(biāo)定后的結(jié)果，圖4為數(shù)值模擬結(jié)果與試驗(yàn)結(jié)果。在模型參數(shù)不變的條件下，開展了圍壓900 kPa 的三軸剪切試驗(yàn)?zāi)M，結(jié)果如圖4所示。

從圖中可以看出，無論是軸向應(yīng)力-應(yīng)變曲線，還是體積應(yīng)變-軸向應(yīng)變曲線，數(shù)值模擬結(jié)果都能與試驗(yàn)結(jié)果吻合較好；因此，可以認(rèn)為本文建立的粗粒料顆粒破碎數(shù)值模型和選取的模型參數(shù)能夠較好地模擬室內(nèi)試驗(yàn)開展的堆石料靜力特性。

圖4 堆石料三軸剪切試驗(yàn)結(jié)果Fig.4 Triaxial shear test results of rockfill material

4 數(shù)值模擬結(jié)果分析

4.1 圍壓影響分析

本文采用Marsal[13]提出的顆粒破碎指標(biāo)Bg定量描述顆粒破碎的程度。

式中：Wki為試驗(yàn)開始前某級(jí)配粒組的含量；Wkf為試驗(yàn)后對(duì)應(yīng)同一級(jí)配粒組的含量。

表2為數(shù)值模擬過程中，不同軸向應(yīng)變時(shí)刻對(duì)應(yīng)的顆粒破碎率。從表中可以看出，采用本文提出的方法在數(shù)值模擬過程中軸向應(yīng)變?yōu)?5%時(shí)刻的顆粒破碎率與試驗(yàn)結(jié)果基本一致，而采用膨脹法的數(shù)值模擬結(jié)果差異較大，結(jié)合不同圍壓下軸向應(yīng)變、體積應(yīng)變和軸向偏應(yīng)力條件下的數(shù)值模擬結(jié)果與試驗(yàn)結(jié)果的對(duì)比情況，表明本文建立的細(xì)觀數(shù)值模型基本能夠反映堆石料的在荷載作用下的顆粒破碎特征。從表2 還可以看出，試樣在制備和等向固結(jié)過程中存在明顯的顆粒破碎特征，隨著圍壓的增大，本次試驗(yàn)獲得的顆粒破碎率變化不明顯。不同圍壓下的試驗(yàn)制備和固結(jié)過程而引起的顆粒破碎均為1.13，約占整個(gè)破碎率的9.7%和8.1%。該結(jié)論與高玉峰等[3]的研究成果類似。其原因可能是由于試驗(yàn)初始孔隙率較小，相對(duì)密度較大，試驗(yàn)制備過程對(duì)顆粒破碎的影響明顯優(yōu)于試驗(yàn)圍壓對(duì)顆粒破碎的影響。

表2 堆石料試驗(yàn)和數(shù)值模擬顆粒破碎率Table 2 Testing and numerical simulation results of particle breakage rate of rockfill

隨著軸向應(yīng)變的增大，顆粒破碎逐漸增多，破碎率增大，近似呈非線性增長(zhǎng)；當(dāng)去除試樣制備和固結(jié)過程引起的顆粒破碎，不同圍壓和不同模擬方法獲得的顆粒破碎率歸一化后具有較好的一致性（如圖5 所示），可近似采用雙曲線函數(shù)進(jìn)行擬合。

式中：εa為軸向應(yīng)變。

圖5 歸一化顆粒破碎與軸向應(yīng)變關(guān)系Fig.5 Relationships between axial strain and normalized particle breakage

圖6為不同圍壓條件下三軸剪切試驗(yàn)初始和最終級(jí)配曲線及根據(jù)Einav[14]基于分形理論提出的粗粒料顆粒破碎理論曲線，公式如下：

式中：dM為試樣中的最大顆粒直徑；d為顆粒直徑；α為分形維數(shù)，對(duì)于多數(shù)材料設(shè)定為常數(shù)2.6。

從圖6 可以看出，隨著圍壓的增大不同級(jí)配條件下的試樣最終的級(jí)配曲線趨于Einav[14]提出的顆粒破碎的最終級(jí)配（圖中α=2.6）。表明，堆石料在顆粒破碎過程中存在一定的分形特征。

圖6 不同圍壓下試驗(yàn)級(jí)配曲線Fig.6 Gradation curves under different confining pressures

4.2 數(shù)值模擬方法對(duì)比分析

數(shù)值模擬成功的關(guān)鍵在于建立合適的數(shù)值模型，不同的數(shù)值模型可能獲得的試驗(yàn)規(guī)律是相反的，模型參數(shù)的選擇對(duì)數(shù)值模擬結(jié)果影響較大，但規(guī)律基本不變。基于此，采用相同的細(xì)觀數(shù)值模擬參數(shù)，分析了不同模擬方法對(duì)試驗(yàn)結(jié)果的影響。

圖7為采用密度法和膨脹法兩種分析方法獲得的數(shù)值模擬結(jié)果。從圖中可以看出，對(duì)于軸向應(yīng)變與偏應(yīng)力的關(guān)系，兩種不同數(shù)值模擬方法的結(jié)果差異較小，與試驗(yàn)結(jié)果基本吻合（見圖7(a)）。但在軸向應(yīng)變與體積應(yīng)變方面，兩者差異較大，采用膨脹法獲得的結(jié)果存在明顯的剪脹現(xiàn)象，見圖7(b)。在軸向應(yīng)變小于4%的情況下，兩種不同數(shù)值模擬方法獲得的數(shù)值模擬結(jié)果基本一致。在應(yīng)變大于4%的情況下，兩種不同數(shù)值模擬方法獲得的數(shù)值模擬結(jié)果差異較大。對(duì)比顆粒破碎率可知，在低應(yīng)變條件下顆粒破碎較小，因此，兩種分析方法差別不大；隨著軸向應(yīng)變的增大，顆粒破碎越多，兩種分析方法之間的差異逐漸顯露出來。

圖7 不同模擬方法結(jié)果分析Fig.7 Results obtained by different simulation methods

5 結(jié) 論

（1）采用Hoek 等對(duì)不同尺寸巖石試樣的研究成果和顆粒破碎過程中保持總質(zhì)量不變，通過改變密度來實(shí)現(xiàn)，建立的數(shù)值模型能夠近似反映粗粒料顆粒破碎特性。

（2）試驗(yàn)過程中，隨著軸向應(yīng)變的增大，顆粒破碎逐漸增多，去除試樣制備和固結(jié)過程引起的顆粒破碎，不同圍壓條件下的顆粒破碎率歸一化后基本重合，且可以近似采用雙曲線函數(shù)進(jìn)行擬合。

（3）顆粒破碎隨著圍壓的增大，逐漸增大，試驗(yàn)級(jí)配趨于Einav 提出的顆粒破碎的最終級(jí)配。

[1]沈珠江.面板堆石壩應(yīng)力應(yīng)變分析的若干問題[C]//中國(guó)混凝土面板堆石壩十年學(xué)術(shù)研討會(huì)論文集(1985-1995).北京:中國(guó)水利水電出版社,1995:43－50.SHEN Zhu-jiang.Problems of concrete faced rockfill dam stress and deformation analysis[C]//Proceedings of the 10th Chinese Academic Conference on Concrete Faced Rockfill Dam(1985-1995).Beijing:China Water &Power Press,1995:43－50.

[2]劉漢龍,秦紅玉,高玉峰,等.堆石粗粒料顆粒破碎試驗(yàn)研究[J].巖土力學(xué),2005,26(4):562－566.LIU Han-long,QIN Hong-yu,GAO Yu-feng,et al.Experimental study on particle breakage of rockfill and coarse aggregates[J].Rock and Soil Mechanics,2005,26(4):526－566.

[3]高玉峰,張兵,劉偉,等.堆石料顆粒破碎特征的大型三軸試驗(yàn)研究[J].巖土力學(xué),2009,30(5):1237－1241.GAO Yu-feng,ZHANG Bing,LIU Wei,et al.Experimental study on particle breakage behavior of rockfills in large-scale triaxial tests[J].Rock and Soil Mechanics,2009 30(5):1237－1241.

[4]陳生水,傅中志,韓華強(qiáng),等.一個(gè)考慮顆粒破碎的堆石料彈塑性本構(gòu)模型[J].巖土工程學(xué)報(bào),2011,(10):1489－1495.CHEN Sheng-shui,FU Zhong-zhi,HAN Hua-qiang,et al.An elastoplastic model for rockfill materials considering particle breakage[J].Chinese Journal of Geotechnical Engineering,2011,(10):1489－1495.

[5]BAGHERZADEH-KHALKHALL A,MOHAMMADI S.Numerical simulation of particle breakage of angular particles using combined DEM and FEM[J].Powder Technology,2011,20(11):15－29.

[6]BAGHERZADEH-KHALKHALL A,MOHAMMADI S.Micromechanics of breakage in sharp-edge particles using combined DEM and FEM[J].Particuology,2008,(6):347－361.

[7]楊貴,肖楊,高德清.粗粒料三維顆粒流數(shù)值模擬及其破壞準(zhǔn)則研究[J].巖土力學(xué),2010,31(增刊2):402－406.YANG Gui,XIAO Yang,GAO De-qing.Study ofthree-dimensional particle flow numerical simulation and failure criterion for coarse-grained soils[J].Rock and Soil Mechanics,2010,31(Supp.2):402－406.

[8]周健,王家全,曾遠(yuǎn),等.顆粒流強(qiáng)度折減法和重力增加法的邊坡安全系數(shù)研究[J].巖土力學(xué),2009,30(6):1549－1554.ZHOU Jian,WANG Jia-quan,ZENG Yuan,et al.Slope safety factor by methods of particle flow code strength reduction and gravity increase[J].Rock and Soil Mechanics,2009,30(6):1549－1554.

[9]OLIVIER T D V,JEAN-CLAUDE C.Numerical model of crushing of grains inside two-dimensional granular materials[J].Powder Technology,1999,(105):190－198.

[10]劉君,劉福海,孔憲京.考慮破碎的堆石料顆粒流數(shù)值模擬[J].巖土力學(xué),2008,29(增刊):107－112.LIU Jun,LIU Fu-hai,KONG Xian-jing.Particle flow code numerical simulation of particle breakage of rockfill[J].Rock and Soil Mechanics,2008,29(Supp.):107－102.

[11]BEN-NUM O,EINAV I.The role of self-organization during confined comminution of granular materials[J].Philosophical Transactions of the Royal Society A,2010,36(8):231－247.

[12]YANG G,LIU K L,LIU H L.Dynamic deformation behaviour of coarse aggregate under cyclic loading by PFC3D[J].European Journal of Environmental and Civil Engineering,2013,17(Supp.1):282－293.

[13]MARSAL R J.Large-scale testing of rockfill materials[J].Journal of Soil Mechanics and Foundation Division,ASCE,1967,93(2):27－43.

[14]EINAV I.Breakage mechanics-Part I:Theory[J].Journal of the Mechanics and Physics of Solids,2007,55(6):1274－1297.