張寶玉， 張昌鎖*， 王晨龍， 郝保欽

(1.太原理工大學(xué) 礦業(yè)工程學(xué)院，太原 030024；2.太原理工大學(xué) 應(yīng)用力學(xué)研究所，太原 030024)

1 引 言

顆粒流程序PFC屬于離散元法，可以解決非連續(xù)介質(zhì)力學(xué)問題，并能從細觀角度分析巖石的破壞機理。PFC用顆粒集合體來模擬巖石，通過定義顆粒間的黏結(jié)接觸模型使得顆粒產(chǎn)生相互作用來重現(xiàn)巖石的力學(xué)性質(zhì)。起初一般采用平行黏結(jié)模型(PBM)來模擬巖石，但研究發(fā)現(xiàn)PBM模擬巖石有如下缺陷[1-3]。壓拉比過低，不能滿足實際巖石的壓拉比要求；內(nèi)摩擦角偏小。為此，Potyondy[4]提出了平節(jié)理模型(FJM)，將圓形顆粒構(gòu)造成多邊形，可抑制接觸破壞后顆粒的旋轉(zhuǎn)，能夠模擬出大壓拉比。因此，本文采用FJM來模擬巖石。

2 平節(jié)理模型

平節(jié)理模型(FJM)是由晶粒及晶粒間賦予平節(jié)理接觸(FJC)構(gòu)成的，晶粒由圓盤顆粒和名義面(notional surfaces)組成，因此FJM相當于把顆粒構(gòu)造成了多邊形，可提供足夠的自鎖效應(yīng)，抑制接觸破壞后顆粒的旋轉(zhuǎn)，能夠模擬大壓拉比情況。晶粒間的接觸實際上是名義面間的接觸，接觸界面(interface)位于名義面之間。在PFC2D中，F(xiàn)JC是一條線段，該線段平均分成多段，每段代表一個單元，如圖1所示。每個單元的狀態(tài)是獨立的，可以是黏結(jié)的或是非黏結(jié)的。黏結(jié)單元作用機理為，當其所受法向拉應(yīng)力σ+>張拉強度σf時，單元破壞，生成張拉裂紋；當其所受切向應(yīng)力τ>τf時，單元破壞，生成剪切裂紋。非黏結(jié)單元作用機理為，不能承受拉力；當其所受切向應(yīng)力τ>τf時，沿接觸界面發(fā)生滑移。其中，黏結(jié)狀態(tài)和非黏結(jié)狀態(tài)下的剪切強度τf分別見式(1,2)，

圖1 平節(jié)理接觸模型[4]

(1，2)

3 正交試驗數(shù)值模擬

3.1 數(shù)值試驗

進行數(shù)值試驗包括建立模型和加載兩個部分。

(1) 建立模型

首先，生成上下左右4面墻，墻圍成的區(qū)域大小即為模型尺寸(50 mm×100 mm)；然后，按顆粒尺寸及分布形式在墻圍成區(qū)域內(nèi)生成重疊量較大的顆粒，計算彈開顆粒并使模型平衡；接著，通過伺服機制移動墻給模型施加一定的應(yīng)力，使得顆粒體系均勻；模型中會存在少量懸浮顆粒(顆粒上的接觸數(shù)少于3)，通過縮放這些懸浮顆粒的大小使其上的接觸數(shù)≥3來達到消除懸浮顆粒的目的；最后，給模型接觸賦予細觀參數(shù)并將模型狀態(tài)清零(包括顆粒速度和位移等)。建模完成后的數(shù)值模型如圖2所示。

圖2 數(shù)值模型

(2) 加載

單軸壓縮試驗。刪除模型兩側(cè)的墻，給上下兩面墻施加恒定速度來實現(xiàn)加載，通過上下墻監(jiān)測軸向應(yīng)力和應(yīng)變，設(shè)置測量圓監(jiān)測橫向應(yīng)變，據(jù)此可獲取單軸抗壓強度UCS、彈性模量E及泊松比ν。

直接拉伸試驗。刪除模型中所有的墻，在試樣上下端給一定厚度顆粒恒定速度來進行加載，設(shè)置測量圓監(jiān)測軸向應(yīng)力應(yīng)變，據(jù)此獲取抗拉強度TS。

雙軸壓縮試驗。根據(jù)伺服原理控制兩側(cè)的墻實現(xiàn)圍壓的施加，給上下兩面墻施加恒定速度實現(xiàn)加載，并通過上下墻監(jiān)測軸向應(yīng)力，圍壓σ3和軸向應(yīng)力σ1呈線性關(guān)系，

σ1=kσ3+b

(3)

(4)

(5)

3.2 試驗參數(shù)選取

3.2.1 宏觀參數(shù)選取

(6)

巖石的起裂強度σc i是判斷圍巖損傷范圍的重要指標。在PFC中可通過裂紋數(shù)來確定巖石的σc i,文獻[1，15]將裂紋數(shù)達到峰值強度對應(yīng)裂紋數(shù)的10%時所對應(yīng)的應(yīng)力視為σc i。

3.2.2 細觀參數(shù)選取

綜上，本文研究選取的宏細觀參數(shù)列入表1。

表1 選取的宏細觀參數(shù)

3.3 正交試驗設(shè)計

正交試驗是研究多因素多水平的試驗方法，主要是依據(jù)正交性在全面試驗中選擇具有代表性的點(試驗點分布均勻)進行試驗，用正交設(shè)計表來安排試驗并進行數(shù)據(jù)分析，具有高效性，結(jié)果可靠簡單易行[16]。

FJM中涉及很多細觀參數(shù)，正交試驗設(shè)計之前，為降低參數(shù)標定的難度，確定一些參數(shù)的值如下。

表2 因素水平

4 結(jié)果分析

表3 正交數(shù)值試驗方案及結(jié)果Tab.3 Orthogonal numerical test scheme and results

4.1 多因素方差分析

方差分析(F檢驗)用于多個樣本均數(shù)差別的顯著性檢驗。根據(jù)方差分析計算過程[19]計算出試驗因素及誤差的自由度分別為3和7，查表可知顯著性水平α=0.05時，F(xiàn)0.95(3,7)=4.35,顯著性水平α=0.01時，F(xiàn)0.99(3,7)=8.45；計算各試驗因素的F值并繪成圖3。認為當4.35<因素F值<8.45時，因素對結(jié)果有顯著影響；當因素F值> 8.45 時，因素對結(jié)果有非常顯著的影響；因素的F值越大，其對結(jié)果的影響越大。據(jù)圖3分析如下。

圖3 多因素方差分析F值

4.2 回歸分析

根據(jù)4.1節(jié)的分析，用巖石宏觀參數(shù)的主要影響因素對其進行簡單的線性擬合，結(jié)果見式(7～12)，可見只有式(7,8,10)的R2大于0.85，擬合效果較好，說明宏細觀參數(shù)間的線性關(guān)系良好；其余公式的R2則較小，效果差，說明宏細觀參數(shù)間的關(guān)系復(fù)雜，需要進一步分析，因此將試驗結(jié)果平均值變化趨勢繪于圖4，分析如下。

圖4 試驗結(jié)果平均值變化趨勢

E=0.858Ef-4.677Kf+12.021，R2=0.976

(7)

(8)

UCS=-0.037Ef-10.397Kf+31.312μf+

3.333σf+4.378(cf/σf)-41.644Rs d+

R2=0.819

(9)

Ts=0.007Ef+0.25σf-6.145Rs d+

R2=0.867

(10)

R2=0.607

(11)

R2=0.798

(12)

(1) 圖4(a)顯示，Ef和kf對E的影響趨勢線變化很明顯，且Ef和kf與E均呈良好的線性關(guān)系，而E隨其余參數(shù)的變化趨勢線基本水平。因此，用E的主要影響因素Ef和kf對其進行簡單線性擬合的效果就很好，見式(7)，對應(yīng)表4中第1個公式。

表4 宏細觀參數(shù)間的擬合公式

5 細觀參數(shù)標定流程及實例驗證

5.1 FJM細觀參數(shù)標定流程

圖5 標定流程

5.2 實例驗證

劉翌晨[21]從齊熱哈塔爾引水隧洞中鉆取片麻花崗巖試樣進行物理試驗來獲取其宏觀力學(xué)參數(shù)，試驗結(jié)果平均值列入表5，本文以此為依據(jù)，進行片麻花崗巖的細觀參數(shù)標定。

表5 片麻花崗巖試樣宏觀參數(shù)試驗值和模擬值

表6 片麻花崗巖細觀參數(shù)標定結(jié)果

數(shù)值模擬及物理試驗得到的單軸壓縮應(yīng)力-應(yīng)變曲線如圖6所示，其中編號為Y-1,Y-2,Y-3和Y-4的曲線分別是物理試驗中4個試樣對應(yīng)的結(jié)果?？梢钥闯?，模擬得到的曲線表現(xiàn)的特征與物理試驗得到的曲線表現(xiàn)的特征接近，但數(shù)值模擬得到的曲線沒有體現(xiàn)壓密階段特征(應(yīng)力-應(yīng)變曲線輕微下凸)，這是由于在建模過程中為得到均勻的顆粒體系進行了伺服，壓密了顆粒體系，因此施加上應(yīng)力即出現(xiàn)彈性變形特征(應(yīng)力-應(yīng)變曲線近似直線)。

圖6 單軸壓縮應(yīng)力-應(yīng)變曲線

數(shù)值試驗和室內(nèi)物理試驗試樣破壞結(jié)果如圖7所示，可見模擬結(jié)果與物理試驗結(jié)果吻合較好。單軸壓縮情況下出現(xiàn)貫通試樣上下端的主裂紋，由巖石內(nèi)部拉伸破壞造成，體現(xiàn)了片麻花崗巖的劈裂張拉破壞特征；三軸壓縮下試樣出現(xiàn)了剪切破裂帶。

圖7 試樣破壞形態(tài)

6 結(jié) 論

(1) 正交數(shù)值試驗結(jié)果顯示采用FJM可以模擬出大壓拉比。對試驗結(jié)果進行多因素方差分析得到了各細觀參數(shù)對宏觀參數(shù)的影響程度，確定了影響各宏觀參數(shù)的主要細觀參數(shù)。

(3) 提出FJM細觀參數(shù)標定流程，然后以齊熱哈塔爾片麻花崗巖物理試驗結(jié)果平均值為依據(jù)進行細觀參數(shù)標定，用標定好的模型進行數(shù)值試驗，結(jié)果顯示模擬得到的宏觀參數(shù)值、應(yīng)力-應(yīng)變曲線及破壞形式接近物理試驗結(jié)果，驗證了本文PFC2D平節(jié)理模型細觀參數(shù)標定方法的可行性。