基于PFC2D模擬的礦物粒徑非均質(zhì)效應(yīng)研究*

韓振華 張路青 周 劍

(①中國(guó)科學(xué)院地質(zhì)與地球物理研究所， 中國(guó)科學(xué)院頁(yè)巖氣與地質(zhì)工程重點(diǎn)實(shí)驗(yàn)室 北京 100029)(②中國(guó)科學(xué)院大學(xué) 北京 100049)

0 引 言

巖體作為一種天然的非均質(zhì)體，由結(jié)構(gòu)面和被結(jié)構(gòu)面包圍的巖石塊體共同組成。其中，巖石是由一種或多種強(qiáng)度和變形特征不同的礦物構(gòu)成，礦物顆粒往往會(huì)呈現(xiàn)出不同的晶形，并被不同的膠接材料黏結(jié)在一起(劉廣等， 2013)。巖石與巖體的重要區(qū)別就是巖體包含若干不連續(xù)面(如：裂隙、節(jié)理、層理等)，不連續(xù)面的存在使得巖體的強(qiáng)度遠(yuǎn)低于巖石強(qiáng)度。由于細(xì)觀層次上的非均質(zhì)性，不同巖石之間強(qiáng)度上也存在很大差異。巖石的非均質(zhì)性可定義為：巖石在形成過程中受成巖環(huán)境、成巖作用和構(gòu)造作用的影響，物質(zhì)組成及其結(jié)構(gòu)在空間分布、物理、力學(xué)屬性等方面的不均勻變化(王香增等， 2016)。巖石物質(zhì)組成的非均質(zhì)性主要體現(xiàn)在：礦物顆粒的形態(tài)及大小、強(qiáng)度差異及礦物顆粒的接觸方式(Lan et al.，2010)。非均質(zhì)性是導(dǎo)致巖石強(qiáng)度不同的根本原因之一(尤明慶等， 2000； 尹小濤等， 2011)，一直以來都是一個(gè)重要的研究課題。

巖石作為礦物顆粒的集合體，礦物粒徑非均質(zhì)性對(duì)巖石力學(xué)特性有很大的影響。一些學(xué)者很早就認(rèn)識(shí)到粒徑的重要作用，并開展了相應(yīng)的實(shí)驗(yàn)研究。Robina et al.(1996)通過對(duì)大理巖單軸壓縮實(shí)驗(yàn)研究表明巖石峰值強(qiáng)度隨粒徑增加而減小，并與粒徑平方呈反比。Fredrich et al.(1990)對(duì)4種不同粒徑的大理巖進(jìn)行三軸壓縮實(shí)驗(yàn)，發(fā)現(xiàn)巖石由脆性向塑性轉(zhuǎn)換所需的圍壓與粒徑成反比。Eberhardt et al.(1999)通過單軸壓縮實(shí)驗(yàn)分析了粒徑對(duì)巖石內(nèi)部裂紋啟動(dòng)和擴(kuò)展的影響，并指出巖石強(qiáng)度隨粒徑增大而減小。Sabri et al.(2016)基于巴西圓盤實(shí)驗(yàn)研究了粒徑對(duì)巖體斷裂韌度的影響，發(fā)現(xiàn)斷裂韌度隨粒徑增大呈現(xiàn)先增大后降低的規(guī)律。

還有部分學(xué)者利用數(shù)值模擬對(duì)巖石力學(xué)特性進(jìn)行研究，常用的方法有基于連續(xù)介質(zhì)力學(xué)的有限元(FEM)、有限差分(FEM)等和基于非連續(xù)介質(zhì)力學(xué)的離散元(DEM)等(魏懷鵬等， 2006)。但巖石非均質(zhì)性、不連續(xù)性及各向異性的特點(diǎn)使得連續(xù)介質(zhì)力學(xué)理論在模擬巖石破壞方面有很大局限性。離散元作為解決非連續(xù)介質(zhì)力學(xué)問題的重要數(shù)值模擬方法更能真實(shí)地還原巖石破壞過程，基于離散元理論的顆粒流程序PFC也成為研究巖石類材料的重要工具(周健等， 2000)。一些學(xué)者通過PFC數(shù)值模擬方法研究了粒徑對(duì)巖石力學(xué)特性影響。Yang et al.(2006)基于PFC2D 模型指出，彈性模量隨L/R(L為接觸面寬度，R為顆粒平均半徑)的增加呈增大的趨勢(shì)，當(dāng)L/R大于一定值(L/R>130)時(shí)彈性模量增長(zhǎng)緩慢并趨于穩(wěn)定，泊松比隨L/R的增大而減小。Potyondy et al.(2004)研究發(fā)現(xiàn)，在PFC3D模型中，彈性模量、單軸抗壓強(qiáng)度及抗拉強(qiáng)度隨粒徑的減小均呈增大趨勢(shì)。還有部分學(xué)者研究了模型中粒徑分布特征對(duì)宏觀特性的影響。Ding et al.(2014)研究發(fā)現(xiàn)，在模型平均粒徑不變的前提下，隨粒徑比Rmax/Rmin的增加，彈性模量和峰值強(qiáng)度降低，泊松比升高。Koyama et al.(2007)研究表明，即使粒徑參數(shù)相同，顆粒在模型內(nèi)隨機(jī)分布的位置不同，所得到的宏觀參數(shù)也不同，因顆粒分布狀態(tài)改變導(dǎo)致的模型宏觀參數(shù)變異系數(shù)隨L/R的增加而減小，當(dāng)L/R>50 時(shí)趨于穩(wěn)定。

已有研究在一定程度上揭示了粒徑與材料力學(xué)特性之間的關(guān)系，但在粒徑非均質(zhì)性方面還缺乏系統(tǒng)的認(rèn)識(shí)。對(duì)粒徑的研究通常是通過改變平均粒徑R或粒徑比Rmax/Rmin值來實(shí)現(xiàn)，粒徑在Rmin～Rmax范圍內(nèi)連續(xù)分布。然而實(shí)際巖樣是由多種礦物組成，不同礦物顆粒的粒徑相差很大，如礫巖中粒徑大于2imm的顆粒碎屑含量占30%以上，而砂、粉砂等填隙物粒徑很??； 花崗巖中鉀長(zhǎng)石、斜長(zhǎng)石粒徑可達(dá)3imm以上，而云母粒徑小于1imm，這說明巖石中礦物顆粒粒徑并非連續(xù)分布。通過設(shè)置不同粒徑的組合，本次研究考慮了粒徑的非連續(xù)分布情況。研究中設(shè)計(jì)了6種不同的粒徑組合方案，以粒徑比Rmax/Rmin值作為非均質(zhì)性的量化指標(biāo)，旨在研究礦物粒徑非均質(zhì)性對(duì)宏觀力學(xué)特性的影響。由于Lac du Bonnet(LDB)花崗巖是一種很有代表性的脆性巖石，其巖石力學(xué)參數(shù)全面(Potyondy et al.，2004)，常被用來對(duì)模型進(jìn)行對(duì)比和校驗(yàn)。以LDB花崗巖為例，作者采用PFC2D 顆粒流程序來模擬巖石單軸壓縮破壞，分析了不同粒徑組合下的宏觀響應(yīng)，為掌握礦物粒徑對(duì)巖石強(qiáng)度及變形特性的力學(xué)效應(yīng)提供依據(jù)。

1 顆粒流分析法的基本原理

顆粒流分析法通常把離散介質(zhì)看作是無數(shù)顆粒單元的集合體，每個(gè)顆粒單元獨(dú)立運(yùn)動(dòng)，單元之間相互接觸、相互作用，對(duì)非連續(xù)介質(zhì)力學(xué)行為的研究就是通過模擬顆粒單元間的這些關(guān)系來實(shí)現(xiàn)的。如圖1所示，顆粒流分析法的基本計(jì)算原理為(Cundall， 1971； Cho et al.，2007)：(1)根據(jù)顆粒單元的位置和半徑確定單元相互間的接觸關(guān)系，采用接觸本構(gòu)模型(力-位移定律)計(jì)算接觸顆粒單元間的接觸力； (2)根據(jù)計(jì)算得到的接觸力，運(yùn)用牛頓第二定律，計(jì)算顆粒單元運(yùn)動(dòng)的加速度、速度和位移，更新顆粒單元的位置與相互的接觸關(guān)系重新計(jì)算接觸力。如此循環(huán)計(jì)算，直到試件破壞或滿足預(yù)定條件，計(jì)算結(jié)束。

圖1 顆粒-黏結(jié)系統(tǒng)的受力-位移關(guān)系Fig. 1 Force-displacement behavior of grain-bonding system

顆粒流分析法的基本接觸模型有3種：接觸剛度模型、滑動(dòng)模型和平行黏結(jié)模型。由于平行黏結(jié)(圖1)可以傳遞顆粒之間的力和力矩、黏結(jié)強(qiáng)度和巖石材料強(qiáng)度有較好的對(duì)應(yīng)關(guān)系，因此巖石類材料的細(xì)觀模擬多采用平行黏結(jié)模型。平行黏結(jié)模型可以假想為膠結(jié)聚合物，顆粒在接觸點(diǎn)處膠結(jié)在一起，同時(shí)黏結(jié)物自身受力時(shí)，允許產(chǎn)生一定的變形，當(dāng)接觸力超過黏結(jié)強(qiáng)度時(shí)，黏結(jié)破壞。

2 顆粒流模型及計(jì)算方案設(shè)計(jì)

2.1 顆粒流模型的建立

本次研究采用平行黏結(jié)模型，模型尺寸為100imm×50imm。首先生成試樣模型墻體，確定模型邊界，然后通過伺服調(diào)節(jié)法控制模型的“墻體”運(yùn)動(dòng)，模擬加載過程，從而實(shí)現(xiàn)單軸壓縮模擬試驗(yàn)。根據(jù)LDB花崗巖的礦物含量及顆粒特征尺寸(Potyondy et al.，2004)，計(jì)算出礦物的平均半徑為1.2imm(表1)。因此，模型顆粒最小半徑取0.8imm，粒徑比選為1.66，以保證模型顆粒平均半徑為1.2imm。采用“試錯(cuò)法”(徐金明等， 2010)反復(fù)改變細(xì)觀參數(shù)，使得模型的宏觀力學(xué)參數(shù)接近LDB花崗巖室內(nèi)試驗(yàn)力學(xué)參數(shù)(彈性模量70iGPa、峰值強(qiáng)度224iMPa、泊松比0.26)，最終確定的模型細(xì)觀物理力學(xué)性質(zhì)參數(shù)(表2)。校核后，花崗巖室內(nèi)試驗(yàn)與模型模擬宏觀力學(xué)參數(shù)吻合較好，彈性模量、泊松比和峰值強(qiáng)度分別為70iGPa、0.26iMPa和207iMPa。

表1 LDB花崗巖礦物含量及粒徑Table 1 Composition and average grain sizes of LDB granite

表 2 PFC軸壓實(shí)驗(yàn)?zāi)Ｐ图?xì)觀參數(shù)Table 2 Micro physico-mechanical parameters of particle

2.2 計(jì)算方案設(shè)計(jì)

本次研究?jī)?nèi)容為礦物粒徑分布對(duì)巖石宏觀力學(xué)特性的影響。模擬過程中，通過設(shè)置不同的粒徑組合來體現(xiàn)粒徑非均質(zhì)性。由于巖土材料內(nèi)部存在顆粒尺寸效應(yīng)，根據(jù)前人研究成果(尹小濤等， 2011)，當(dāng)粒徑非常小時(shí)，它對(duì)模型特性的影響很小，以致可以忽略粒徑的作用。但當(dāng)粒徑相對(duì)變大時(shí)，模型宏觀特性受粒徑的影響較大。因此，研究中首先采用均質(zhì)模型分析了顆粒的尺寸效應(yīng)，為后續(xù)非均質(zhì)模型方案設(shè)計(jì)提供參考。共設(shè)計(jì)了一種均質(zhì)模型及6種非均質(zhì)模型，具體如下：

圖2 非均質(zhì)方案1～6中模型粒徑分布(以平均粒徑0.6imm，粒徑比2.0模型為例)Fig. 2 Particle size distribution in scenario 1-6(case study of model with average grain size of 0.6imm and size ratio of 2.0)a. 連續(xù)粒徑組合-方案1； b. 10種粒徑組合-方案2； c. 8種粒徑組合-方案3； d. 5種粒徑組合-方案4； e. 3種粒徑組合-方案5； f. 兩種粒徑組合-方案6

(1)均質(zhì)模型：保持模型宏觀尺寸及其他細(xì)觀參數(shù)不變，固定粒徑比為1，使模型粒徑處于均勻分布狀態(tài)。改變粒徑，使其在0.30～2.00imm之間變化，以此來研究模型在均質(zhì)狀態(tài)下顆粒的尺寸效應(yīng)。

(2)非均質(zhì)模型：共設(shè)計(jì)6種方案(表3)，每種方案下設(shè)置5種平均粒徑模型，每種平均粒徑下設(shè)置4種粒徑比，最大與最小粒徑分布情況如表4所示。不同的方案中粒徑的數(shù)量不同，即最大與最小粒徑之間的過渡粒徑種類不同。

表3 6種非均質(zhì)模型方案Table 3 Six scenarios with heterogeneous grain sizes

(1)

其中，i為礦物的種類；Ri為第i種礦物的粒徑；Ci為第i種礦物的含量。模型的粒徑非均質(zhì)程度用粒徑比來表征，粒徑比越大，非均質(zhì)性越強(qiáng)。以平均粒徑0.6imm，粒徑比2.0的模型為例，在不同粒徑組合下模型內(nèi)部尺寸分布如圖2所示。

表4 不同平均粒徑、粒徑比下模型中的最大與最小粒徑Table 4 Distribution of maximal and minimum size in the model with different mean grain size and size ratio

由于模型是由礦物顆粒的隨機(jī)分布組成的，PFC程序通過選取不同的隨機(jī)數(shù)種子來實(shí)現(xiàn)不同的顆粒分布。當(dāng)粒徑參數(shù)相同時(shí)，采用不同的顆粒分布狀態(tài)，所造成的模型內(nèi)部應(yīng)力分布不同，因此模型宏觀特性也是不同的。由顆粒分布狀態(tài)不同所造成的宏觀特性離散化是客觀存在的，但可以通過多次模擬取平均值的方式進(jìn)行宏觀趨勢(shì)化分析。本次研究針對(duì)每個(gè)模型試驗(yàn)均選取7個(gè)隨機(jī)分布狀態(tài)，因此，在5種平均粒徑、4種粒徑比下， 6種設(shè)計(jì)方案共包含算例840個(gè)，主要分析了粒徑非均質(zhì)性對(duì)宏觀彈性模量、峰值強(qiáng)度及泊松比的影響。彈性模量采用應(yīng)力為試樣峰值強(qiáng)度一半時(shí)應(yīng)力與應(yīng)變的比值，泊松比采用應(yīng)力為試樣峰值強(qiáng)度一半時(shí)側(cè)向應(yīng)變與軸向應(yīng)變的比值(尤明慶等， 2003)。

圖3 均質(zhì)模型中粒徑及顆粒分布對(duì)宏觀特性的影響Fig. 3 Effects of particle size and distribution on micro-mechanical properties in homogenous modela. 彈性模量； b. 峰值強(qiáng)度； c. 泊松比

3 模擬結(jié)果分析

3.1 粒徑均勻分布下的尺寸效應(yīng)研究

均質(zhì)狀態(tài)下，粒徑及其分布狀態(tài)對(duì)模型宏觀特性的影響如圖3所示。首先，模型的平均彈性模量和峰值強(qiáng)度均隨著粒徑的增大而減小，而泊松比隨粒徑的增大呈增大的規(guī)律，與前人研究成果一致(Potyondy et al.，2004； Ding et al.，2014)。粒徑由0.3imm增長(zhǎng)到2.0imm時(shí)，平均彈性模量降低了20.1%，峰值強(qiáng)度降低了13.4%，泊松比提高了9.4%，說明粒徑變化對(duì)彈性模量的影響比較明顯，其次是峰值強(qiáng)度，對(duì)泊松比影響較小。此外，由7組試樣模擬結(jié)果的標(biāo)準(zhǔn)差(圖3)可以看出，由顆粒分布狀態(tài)的不同導(dǎo)致的模型宏觀力學(xué)特性離散化程度隨粒徑的增大而增大。

由圖3不同粒徑所對(duì)應(yīng)的顆?？倲?shù)可以看出，粒徑小于0.40imm時(shí)，模型顆粒規(guī)模將過萬，這會(huì)導(dǎo)致計(jì)算效率低下，而粒徑過大則會(huì)導(dǎo)致顆?？倲?shù)太少，顆粒結(jié)構(gòu)在應(yīng)力作用下穩(wěn)定性變差，計(jì)算結(jié)果的波動(dòng)也較大。因此，雖然LDB花崗巖礦物顆粒尺寸較大(表1)，但在模擬過程綜合考慮計(jì)算規(guī)模和粒徑的敏感性，所采用的顆粒平均粒徑為0.3～0.7imm。

3.2 粒徑非均質(zhì)性對(duì)巖石宏觀力學(xué)特性的影響

巖石中礦物成分、粒徑、空間分布等都對(duì)其力學(xué)特性有影響，本次設(shè)計(jì)的6種非均質(zhì)方案綜合考慮了粒徑大小、數(shù)量及分布的影響。

圖4 非均質(zhì)方案1～6中模型粒徑非均質(zhì)性對(duì)彈性模量的影響Fig. 4 Effect of particle heterogeneity on elasticity modulus in scenario 1-6a. 連續(xù)粒徑組合-方案1； b. 10種粒徑組合-方案2； c. 8種粒徑組合-方案3； d. 5種粒徑組合-方案4； e. 3種粒徑組合-方案5； f. 兩種粒徑組合-方案6

3.2.1 粒徑非均質(zhì)性對(duì)彈性模量的影響

圖4為粒徑非均質(zhì)性對(duì)彈性模量的影響?？梢钥闯?， 6種方案中彈性模量隨模型平均粒徑的增加均呈減小趨勢(shì)，隨粒徑非均質(zhì)性(粒徑比)的增加也呈現(xiàn)減小的規(guī)律(方案5、方案6除外)。首先，粒徑在均質(zhì)狀態(tài)下，模型的彈性模量較高，在模型由均質(zhì)狀態(tài)向非均質(zhì)狀態(tài)轉(zhuǎn)化階段(粒徑比由1變化到1.2階段)，彈性模量呈現(xiàn)明顯減小的趨勢(shì)。隨著粒徑非均質(zhì)性的增大(粒徑比由1.2變化到2.0階段)彈性模量繼續(xù)降低(圖4a-4id)，但變化幅度較小。

需要指出的是，在3種粒徑組合(圖4e)和兩種粒徑組合(圖4f)模型中，當(dāng)粒徑比為2時(shí)，彈性模量并沒有達(dá)到最低，反而呈增高的趨勢(shì)。如圖4e所示， 3種粒徑組合下粒徑比為2時(shí)的彈性模量大于粒徑比為1.6時(shí)的彈性模量； 圖4f中， 兩種粒徑組合下粒徑比為2時(shí)的彈性模量大于粒徑比為1.2時(shí)的彈性模量。該變化規(guī)律與模型中細(xì)顆粒的含量有關(guān)，如表4所示，隨著粒徑比的增大模型中的最小粒徑減小，粒徑比為2的模型含最小的粒徑。在方案5和方案6中，最小粒徑顆粒含量分別為33.3%和50%，相比于其他方案，最小粒徑的顆粒含量及數(shù)量較高。根據(jù)對(duì)模型均質(zhì)狀態(tài)下尺寸效應(yīng)的研究(圖3)可知，粒徑越小，彈性模量越大。因此在這兩種方案中，粒徑比為2時(shí)，相比于粒徑非均質(zhì)性對(duì)彈性模量的降低作用，細(xì)顆粒對(duì)其提高作用更明顯。

為進(jìn)一步分析細(xì)顆粒對(duì)彈性模量的作用，以兩種粒徑組合為例，取0.4imm和0.6imm兩種平均粒徑的模型，以0.27imm和0.4imm分別作為兩模型中的最小粒徑(細(xì)顆粒粒徑)。模擬過程中，改變模型中細(xì)顆粒的含量，以研究細(xì)顆粒的作用，分析結(jié)果如圖5所示。可以發(fā)現(xiàn)，隨著模型中細(xì)顆粒含量的提高，細(xì)顆粒數(shù)量的增長(zhǎng)十分明顯，彈性模量呈現(xiàn)先減小后增大的規(guī)律。該規(guī)律解釋了圖4e和圖4f中彈性模量的變化，彈性模量在初期表現(xiàn)出降低的趨勢(shì)主要受粒徑非均質(zhì)性的影響，在后期呈現(xiàn)增大的規(guī)律主要受細(xì)顆粒的控制作用。

圖 5 細(xì)顆粒含量對(duì)彈性模量的影響(以兩種粒徑組合、平均粒徑0.4imm和0.6imm的模型為例)Fig. 5 Effects of fine particle content on elasticity modulus(case study of model with two kinds of grain size， mean grain size respectively 0.4imm and 0.6imm)

圖 6 非均質(zhì)方案1～6中模型粒徑非均質(zhì)性對(duì)峰值強(qiáng)度的影響Fig. 6 Effect of particle heterogeneity on peak strengtha. 連續(xù)粒徑組合-方案1； b. 10種粒徑組合-方案2； c. 8種粒徑組合-方案3； . 5種粒徑組合-方案4； e. 3種粒徑組合-方案5； f. 兩種粒徑組合-方案6

3.2.2 粒徑非均質(zhì)性對(duì)峰值強(qiáng)度的影響

圖6為粒徑的非均質(zhì)性對(duì)峰值強(qiáng)度的影響結(jié)果。可以看出，在6種方案下，峰值強(qiáng)度隨模型平均粒徑的增加均呈減小的趨勢(shì)。在非均質(zhì)狀態(tài)(粒徑比1.2、1.6、2.0)下，隨粒徑比的增大，峰值強(qiáng)度也呈減小的規(guī)律。此外，峰值強(qiáng)度隨粒徑比的變化規(guī)律與粒徑組合有關(guān)。在連續(xù)粒徑組合下(圖6a方案1)，均質(zhì)狀態(tài)時(shí)模型具有較高的峰值強(qiáng)度，隨著模型中粒徑種類的減少，粒徑比為1.2時(shí)的峰值強(qiáng)度逐漸提高(圖6b-6f)，并超過均質(zhì)狀態(tài)下模型的峰值強(qiáng)度。但當(dāng)粒徑比為1.6和2.0時(shí)，隨著模型中粒徑組合種類的減少，峰值強(qiáng)度沒有明顯增大或減小的規(guī)律。這說明，模型在粒徑非均質(zhì)程度不高的情況下，即粒徑比較小時(shí)，峰值強(qiáng)度受粒徑組合數(shù)影響較大。

圖 7 粒徑比對(duì)峰值強(qiáng)度的影響(以兩種粒徑組合、平均粒徑0.4imm和0.6imm模型為例)Fig. 7 Effect of size ratio on peak strength(case study of model with two kinds of grain size， mean grain size respectively 0.4imm and 0.6imm)

為進(jìn)一步分析粒徑在非均質(zhì)程度不高時(shí)，峰值強(qiáng)度的變化規(guī)律。以兩種粒徑組合(方案6)為例，增加幾種粒徑比： 1.05、1.1、1.15、1.3、1.4，取0.4imm和0.6imm兩種平均粒徑的模型，進(jìn)行模擬結(jié)果分析。如圖7所示，在模型粒徑非均質(zhì)程度不高的情況下(1.05

圖8 方案1～6中模型粒徑非均質(zhì)性對(duì)泊松比的影響Fig. 8 Effect of particle heterogeneity on Poisson’s ratioa. 10種粒徑組合-方案4； b. 5種粒徑組合-方案5； c. 2種粒徑組合-方案6

3.2.3 粒徑非均質(zhì)性對(duì)泊松比的影響

圖8為粒徑非均質(zhì)性對(duì)泊松比的影響，可以看出，泊松比受粒徑大小及粒徑非均質(zhì)性的影響比較小，但也表現(xiàn)出了明顯的規(guī)律性。以粒徑種類為10、5、2的3種模型為例(圖8a～8c)，可以看出，不同平均粒徑下，泊松比隨模型平均粒徑的增加均呈增加趨勢(shì)，隨粒徑比的增加也呈現(xiàn)出增加的趨勢(shì)，但變化幅度比較小。當(dāng)顆粒平均粒徑由0.3imm增加到0.7imm、粒徑比由1.0增加到2.0，泊松比變化范圍在0.23～0.25之間，其他顆粒組合狀態(tài)下，泊松比也表現(xiàn)了相似的變化規(guī)律。

4 討 論

4.1 顆粒隨機(jī)分布的影響

PFC程序可提供兩種默認(rèn)的顆粒生成方式：連續(xù)均勻分布和高斯分布。連續(xù)均勻分布是指在最大最小粒徑之間，粒徑是連續(xù)的，即本次模擬采用的方案1中的顆粒分布形式，是一種廣泛被采用的顆粒分布方式。本次研究，除了考慮粒徑連續(xù)均勻分布模型外，還增加了多種粒徑組合模型，并分析了它們?cè)诓煌奖认碌暮暧^響應(yīng)。上述顆粒分布指的是一個(gè)模型中粒徑大小的分布情況，實(shí)際上采用相同的粒徑參數(shù)，當(dāng)顆粒的隨機(jī)分布位置不同時(shí)，模型的宏觀特性也不同(Koyama et al.，2007)。

本文在對(duì)粒徑均質(zhì)模型尺寸效應(yīng)分析的同時(shí)，針對(duì)每個(gè)模型，通過7次隨機(jī)試驗(yàn)，研究得出顆粒隨機(jī)位置的不同會(huì)導(dǎo)致宏觀力學(xué)特性的離散化，且該離散程度隨粒徑的增大而增大。PFC通過選取不同的隨機(jī)數(shù)種子(set random seed)來實(shí)現(xiàn)不同的顆粒分布。顆粒在模型內(nèi)部的分布實(shí)際上由一系列的隨機(jī)數(shù)決定，隨機(jī)數(shù)種子就是產(chǎn)生隨機(jī)數(shù)的第一次使用值。顆粒生成機(jī)制是通過蒙特·卡羅方法，將種子值轉(zhuǎn)化為隨機(jī)數(shù)空間中的某一個(gè)點(diǎn)上，即模型內(nèi)部的某一位置，以后產(chǎn)生的隨機(jī)數(shù)都與前一個(gè)隨機(jī)數(shù)有關(guān)。因此，一個(gè)隨機(jī)數(shù)種子對(duì)應(yīng)一種顆粒分布狀態(tài)。

以粒徑為0.5imm的均質(zhì)模型為例。采用隨機(jī)數(shù)種子為51時(shí)，模型內(nèi)局部顆粒分布情況如圖9a； 采用隨機(jī)數(shù)種子為9時(shí)，在相同位置顆粒分布情況如圖9b，兩種顆粒分布狀態(tài)有明顯的區(qū)別，因此造成了模型內(nèi)部應(yīng)力分布的不同，最終導(dǎo)致宏觀特性的不同。此外，相比于粒徑非均質(zhì)模型，均質(zhì)模型具有明顯較高的彈性模量(圖4)，而峰值強(qiáng)度卻沒有表現(xiàn)出相同的規(guī)律(圖6)，反而低于方案4～方案6中粒徑比為1.2的非均質(zhì)模型峰值強(qiáng)度。該現(xiàn)象與均質(zhì)模型在二維狀態(tài)下的顆粒分布狀態(tài)密不可分，如圖9所示，在均質(zhì)模型中，由于顆粒生成時(shí)位置在模型內(nèi)部是隨機(jī)分配的。因此，在生成各向同性應(yīng)力狀態(tài)的集合體后，顆粒并非全部都規(guī)則排列，在顆粒排列不規(guī)則的區(qū)域容易形成應(yīng)力集中，力鏈分布的不均勻性造成了宏觀強(qiáng)度的降低。當(dāng)均質(zhì)模型內(nèi)部顆粒完全規(guī)則排列時(shí)(Lan et al.，2010)，其峰值強(qiáng)度會(huì)明顯高于非均質(zhì)模型。采用三維均質(zhì)模型，這種局部顆粒不均勻排列的現(xiàn)象也會(huì)明顯減小，使得模型峰值強(qiáng)度會(huì)有明顯的提高(Ding et al.，2014)。

圖9 不同隨機(jī)種子數(shù)下模型內(nèi)部應(yīng)力分布及局部顆粒排列(以粒徑為0.5imm的均質(zhì)模型為例)Fig. 9 Contact force chain in the model and local grain distribu￣tion(case study of homogenous model with grain size of 0.5imm)a. 隨機(jī)種子數(shù)51(彈性模量86.6GPa、峰值強(qiáng)度218iMPa)； b. 隨機(jī)數(shù)9(彈性模量84.6iGPa、峰值強(qiáng)度200iMPa)

表 5 方案4各計(jì)算模型的實(shí)際孔隙率Table 5 Final porosity of each model in scenario 4

4.2 孔隙率的影響

巖石是一種天然的多孔材料，其內(nèi)部包含著大量不規(guī)則、跨尺度的孔隙，巖石孔隙率的增大會(huì)造成宏觀彈性模量的降低。本次研究?jī)?nèi)容為粒徑非均質(zhì)性對(duì)宏觀的特性的影響，因此模型孔隙率均設(shè)為0.16以減少模擬過程的研究變量。但需要注意的是，初始孔隙率僅用來確定生成模型所需的顆?？倲?shù)(式(2))，模型在達(dá)到應(yīng)力平衡狀態(tài)，即生成各向同性應(yīng)力狀態(tài)集合體的同時(shí)，模型的孔隙率會(huì)有所改變。

N=A(1-n)/πR2

(2)

以本次模擬的方案6為例，模型最終孔隙率統(tǒng)計(jì)結(jié)果如表5所示，可以看出，孔隙率隨粒徑增大呈增大趨勢(shì)，同時(shí)彈性模量和峰值強(qiáng)度表現(xiàn)出降低的趨勢(shì)。與設(shè)定孔隙率0.16相比，均質(zhì)模型的孔隙率明顯減小，主要原因?yàn)榫|(zhì)模型中大部分顆粒處于一種較密實(shí)的排列狀態(tài)(圖9)，這也是均質(zhì)模型彈性模量相對(duì)較高的原因。不同非均質(zhì)模型的孔隙率變化幅度較小，因此，粒徑比為1.2～2.0時(shí)，模型的彈性模量也變化較小?？梢钥闯觯綄?duì)彈性模量的作用機(jī)制是通過影響孔隙率實(shí)現(xiàn)的。

為進(jìn)一步分析不同初始孔隙率下，實(shí)際孔隙率及宏觀特性的變化，以平均粒徑0.4imm、粒徑比1.2的連續(xù)粒徑模型為例。如表6所示，當(dāng)設(shè)定的初始孔隙率從0.1變化到0.35(增長(zhǎng)250%)時(shí)，實(shí)際模型在單軸壓縮前對(duì)應(yīng)的孔隙率僅從0.151變化到0.164(增長(zhǎng)8.6%)，隨著孔隙率的變化，彈模模量、峰值強(qiáng)度及泊松比的變化趨勢(shì)如圖10所示。

表 6 實(shí)際孔隙率、平均粒徑及顆?？倲?shù)隨設(shè)定孔隙率的變化(以平均粒徑0.4imm、粒徑比1.2的連續(xù)粒徑模型為例)Table 6 Change of final porosity， average grain size and total grain number with design porosity(case study of model with average grain size of 0.4imm and size ratio of 1.2)

圖10 孔隙率對(duì)模型宏觀特性的影響Fig. 10 Effect of porosity on micro-mechanical propertiesa. 彈性模量； b. 峰值強(qiáng)度； c. 泊松比

由圖10可知，彈性模型和峰值強(qiáng)度隨設(shè)定孔隙率的增大呈降低的趨勢(shì)，泊松比呈增大的趨勢(shì)，但變化量均較小。宏觀特性隨孔隙率的變化所表現(xiàn)出的規(guī)律，實(shí)際上與模型顆?？倲?shù)和平均粒徑增大有關(guān)(表6)，初始孔隙率越大，模型顆粒數(shù)越少，因此造成的模型宏觀強(qiáng)度上的降低。因此，在生成模型時(shí)，需要選擇一個(gè)合適的孔隙率值，使得模型的最終孔隙率與平均粒徑能與初始設(shè)定值之間的差距最小。本次模擬采用的平均粒徑為0.4imm。根據(jù)表1，當(dāng)初始孔隙率為0.16時(shí)，模型的最終孔隙率和平均粒徑與設(shè)定值最接近。因此，在PFC二維模型中，要在給定空間內(nèi)生成顆粒，既能保證孔隙比和粒徑符合要求，又能使組合達(dá)到平衡，建議選用的孔隙率為0.14～0.18。

5 結(jié) 論

本文基于顆粒流程序PFC2D研究了粒徑非均質(zhì)性對(duì)巖石宏觀力學(xué)特性(彈性模量、峰值強(qiáng)度、泊松比)的影響。通過設(shè)置不同種類粒徑組合及粒徑比來體現(xiàn)粒徑非均質(zhì)性，共設(shè)計(jì)了6種研究方案，粒徑種類分別為：連續(xù)、10種、8種、5種、3種、2種，每種方案下設(shè)置5種平均粒徑模型： 0.3imm、0.4imm、0.5imm、0.6imm和0.7imm，每種平均粒徑下設(shè)置4種最大最小粒徑比： 1.0、1.2、1.6、2.0，每個(gè)模型通過變換隨機(jī)數(shù)取7個(gè)隨機(jī)分布狀態(tài)，進(jìn)行單軸壓縮試驗(yàn)，并得到以下幾點(diǎn)認(rèn)識(shí)：

(1)彈性模量和峰值強(qiáng)度隨模型平均粒徑的增大而減小，泊松比隨平均粒徑的增大呈增大的趨勢(shì)，由顆粒分布狀態(tài)所引起的宏觀特性的離散程度也會(huì)隨粒徑的增大而增大。

(2)通常情況下，彈性模量隨粒徑非均質(zhì)性增大呈減小的規(guī)律，但不同非均質(zhì)模型之間變化范圍較小。當(dāng)模型內(nèi)部細(xì)顆粒含量足夠高時(shí)，細(xì)顆粒對(duì)模型彈性模量具有提高的作用，此時(shí)細(xì)顆粒效應(yīng)會(huì)大于粒徑非均質(zhì)性效應(yīng)。此外，粒徑對(duì)彈性模量的作用主要是通過影響孔隙率實(shí)現(xiàn)的。

(3)當(dāng)模型粒徑非均質(zhì)程度不高時(shí)，峰值強(qiáng)度隨粒徑非均質(zhì)性增大呈增大的規(guī)律，主要原因?yàn)槟Ｐ蛢?nèi)細(xì)顆粒含量也是增大的，細(xì)顆粒對(duì)模型峰值強(qiáng)度具有提高作用； 隨著尺寸非均質(zhì)性的進(jìn)一步增大，峰值強(qiáng)度呈現(xiàn)減小的規(guī)律，此時(shí)粒徑非均質(zhì)性效應(yīng)大于細(xì)顆粒效應(yīng)。

(4)泊松比隨模型顆粒平均粒徑的增加均呈增加趨勢(shì)，隨粒徑比的增加也呈現(xiàn)出增加的趨勢(shì)，但粒徑大小及粒徑非均質(zhì)性對(duì)其影響較小。

(5)利用PFC程序生成模型顆粒時(shí)，需要選擇合適的初始孔隙率，否則會(huì)導(dǎo)致實(shí)際生成模型的孔隙率和粒徑不符合設(shè)定要求。