充填粗糙節理直剪數值模擬宏細觀分析

許萬忠, 林 杭, 曹日紅

(1. 昆明理工大學國土資源工程學院, 云南 昆明 650093; 2. 中南大學資源與安全工程學院, 湖南 長沙 410083; 3. 中國礦業大學深部巖土力學與地下工程國家重點實驗室, 江蘇 徐州 221116)

經歷億萬年的地質運動,巖體被各種分界面切割呈破碎狀態,這些軟弱界面包括節理、裂隙、褶皺等[1].軟弱界面的存在不但使巖體呈各向異性特征,也在很大程度上削弱了巖體的力學性質,同時也成為控制各類巖體工程穩定性的關鍵因素.在采礦工程及邊坡工程當中,節理經常受到剪切荷載作用而導致巖體滑移及崩塌的事故.因此,深入研究節理的的力學性質及破壞行為極為重要.截至目前,已有學者針對節理巖體的力學性質及破壞過程進行了相關研究[2-5]:蔣宇靜等[3]對巖石節理在不同接觸狀態下的破壞行為進行了研究;曹平等[6]利用劈裂所得巖石節理試樣進行直剪實驗并對形貌演化情況進行了分析;Ghazvinian等[7]采用水泥制作鋸齒形節理試樣并進行了不同法向荷載下的直剪實驗,對剪切峰值強度及殘余強度的演化特征進行了研究.以上研究促進了對節理剪切強度及形貌參數演化的理解,但是巖石實驗主要從宏觀角度進行分析,而不同粗糙節理內部顆粒的相互作用及其細觀參數的影響還需進一步研究.近年來,隨著計算機技術的不斷發展,涌現出了多種對裂隙巖石的力學性質和破壞行為進行模擬的數值計算方法[8-10],其中離散元數值方法尤為適合裂隙巖體壓剪狀態下的破裂及損傷演化研究.目前,PFC離散元數值模擬方法被國內外學者廣泛接受并用來模擬脆性材料的裂紋擴展及損傷演化過程[11-13],其在巖石力學及損傷的研究成果也尤為突出,主要體現在裂隙巖體加載條件下的裂紋起裂、剪切破壞等[4,14-15].基于PFC二維離散元程序建立不同粗糙程度的巖石充填節理模型;通過自編FISH程序語言制定巖石節理間充填層的力學參數與性質,對充填節理的力學性質與破壞過程進行模擬研究;分析節理JRC(粗糙度系數)、接觸面強度、充填體強度與抗剪強度間的關系.

1 節理直剪數值試驗

1.1 數值模型的建立

(1)

ΔU=viΔt,

Δθ3=(ω2-ω1)Δt,

(2)

(3)

(4)

(a) 顆粒接觸示意

(b) 粘結力鏈分解圖1 平行粘結示意Fig.1 Parallel bond model in PFC

新的力和力矩分別由式(5)、(6)計算得出.

(5)

(6)

作用在粘結上的最大法向應力和切向應力為

(7)

在粗糙節理模型建立時,首先生成由顆粒組成的矩形.模型的整體尺寸為:100 mm(長)×100 mm(寬).矩形中的剛性顆粒由程序隨機生成,強度分布符合威布爾函數,其中顆粒最小半徑為0.2 mm,顆粒半徑放大系數為1.66.模型中其余細觀參數如表1中所示.PFC模型由顆粒與顆粒粘結組成.顆粒通過粘結而形成塊體,顆粒類似于天然巖石中的礦物顆粒或混凝土中的谷粒物質.

表1 顆粒細觀參數Tab.1 Microscopic parameters for rock mass

平行粘結模型中,顆粒接觸關系如圖1(a)所示,顆粒是模型的主要框架,屬于剛體,接觸模量、接觸剛度比屬于顆粒固有性質.而平行粘結類似于天然巖石中的黏性礦物及混凝土中的水泥一樣充斥于顆粒物質的空隙中.其中粘結半徑系數、粘結法向強度、粘結切向強度均是平行粘結內部參數.粘結半徑為1表示顆粒兩側的粘結鍵長度等于顆粒直徑.

表1中所示參數是根據實際巖石試樣(泥質灰巖)的宏觀強度參數標定后得到的.由于巖石在加載初期內部微裂隙受壓閉合而導致應力-應變曲線出現非線性增長,所以試樣的彈性模量為曲線直線段斜率.從圖2(a)不難看出,數值計算結果在峰值強度及彈性模量方面與實驗結果吻合良好.圖2(b)、2(c)中所示為數值模型最終破壞模式與實際巖石試樣單軸加載破壞狀態.

(a) 應力-應變曲線

(b) 數值模擬試件破壞形態(c) 真實試件破壞形態圖2 單軸加載下完整試樣應力-應變曲線對比Fig.2 Comparison of stress-strain curves and failure mode with numerical model and real sample under uniaxial compression

顯然,數值計算結果與實驗結果一致.PFC作為一款數值算法,其計算過程中材料及邊界條件均是理想化的,因此和所有數值算法一樣無法模擬出真實巖石中峰值前出現的軟化現象.在數值模擬中峰前曲線基本屬于線性或類線性,而真實試樣中存在泥質和較軟的伴生礦物,在加載初期會出現一定的軟化現象,從而導致圖2(a)中所示的情況.其次,在標定細觀參數時,所取的是切線彈性模量,即應力-應變曲線直線段的斜率,從而在圖2(a)中可以看到,數值曲線峰前段的斜率和實驗結果中直線段斜率相同,峰值強度相同,峰值應變不同.

利用Barton JRC標準輪廓線(圖3(a))為參照,生成不同粗糙程度的節理模型.PFC中需要利用FISH語言自編程序來實現該指令.具體過程為:先將不同粗糙度的巴頓輪廓線轉化為坐標,這一過程中為了保證輪廓線的準確度,將每條輪廓線分成100段;通過編程讀取每一段前后坐標值,之后導入PFC中,生成與輪廓線相似的彎曲墻體;利用墻與顆粒接觸的原理,識別接觸點并記錄下接觸點所屬顆粒的ID值;將與輪廓線墻體相接觸的顆粒找出并進行重新標記.從而,生成了不同粗糙度的充填節理模型.

如圖3(b)所示為PFC中生成的JRC為14～16的粗糙充填節理模型.模型中顏色較淺的顆粒代表巖石,中間顏色較深部分的顆粒代表充填物質.圖3(b)中右邊部分為放大區域,顆粒球心間的連線即為接觸鍵,而球外側的線條代表平行粘結.節理顆粒與圍巖間的虛線輪廓為充填物質與圍巖的接觸面.改變節理處顆粒的強度參數值與充填物質-圍巖接觸面粘結強度值來對充填節理剪切荷載下的力學性質進行研究與分析.雖然完整模型中的顆粒是隨機生成的,但是不同粗糙度的節理面是后期指定生成的,所以前后兩次生成的節理面一致,不對計算結果產生影響.

(a) Barton JRC標準輪廓線(b) JRC為14^16的充填節理模型圖3 Barton JRC標準輪廓線與充填節理模型Fig.3 Barton JRC standard profile and the filled joint model

1.2 直剪試驗模擬

PFC中通常采用移動控制墻來進行加載,圖4為本文中的直剪示意圖.

圖4 直剪示意圖Fig.4 Direct shear diagram for this study

墻體2為法向加載墻,模擬過程中該墻通過伺服機制對模型施加恒定的法向荷載.在離散元計算程序中,由伺服機制獲得的法向荷載并非完全恒定,而是在極小值的誤差范圍內跳動,相比于所施加的法向荷載而言,該誤差完全可以忽略不計.之后,程序通過加載墻5進行切向加載,此時墻體3運動方向和速率與墻5保持一致且剪切速率保持為0.03 mm/(106step),而墻1、4及6保持固定.

2 計算結果與分析

2.1 充填節理剪切力學性質與破壞分析

本試驗中對節理模型施加1.0、1.5、2.0、2.5、3.0 MPa 5種不同的法向荷載.不同法向荷載下充填節理模型將表現出不同的力學性質與破壞程度.圖5為JRC為14～16的充填節理試樣在法向荷載為1 MPa和3 MPa下,剪應力、法向位移以及裂紋累計數的變化情況.

從圖5中可看出:兩種不同法向荷載下,充填節理模型力學性質存在一些差異,整體來看,隨著剪切位移不斷增加,節理剪切應力首先呈現出類線性增長趨勢;當剪切應力達到峰值時,繼續加載,剪切應力出現明顯的下降趨勢,并最終在達到一定值后基本保持穩定;雖然剪切應力保持在一定的范圍,但仍然存在一定的波動,這一現象也與實際剪切曲線相吻合[6];在剪切力未達到峰值之前,隨著剪切位移不斷增加,剪切應力增長明顯,而法向位移并未出現明顯的增長趨勢;剪切應力達到峰值后,隨著剪切位移的增加,剪切應力急劇下降,而法向位移一直處于增長趨勢,同時微裂紋的積累與法向位移呈現出類似的增長趨勢;隨著剪切應力的波動,法向位移依然存在一定的增長,而微裂紋增長在剪切殘余階段基本保持穩定.

(1) 剪切應力(2) 豎向位移(3) 微裂紋數量(a) 法向荷載1 MPa(1) 剪切應力(2) 豎向位移(3) 微裂紋數量(b) 法向荷載3 MPa圖5 不同法向荷載下剪切應力、法向位移、裂紋數量與剪切位移的關系Fig.5 Shear stress, normal displacement, and number of micro-cracks vs. shear displacement under different normal stresses

圖6為JRC為14～16和18～20的充填節理破壞情況.圖中節理面附近聚集了大量的微裂紋,其中黑色裂紋代表拉伸型裂紋,紅色細線代表剪切型裂紋.

從圖6中可看出:不論JRC為14～16還是18～20的充填節理,在法向荷載為3 MPa時節理面損傷程度較大;不論拉伸裂紋還是剪切型裂紋,較高法向荷載下的分布都更為密集.此現象主要是當法向荷載較低時,節理模型在剪切作用下出現較為明顯的剪脹現象,從而節理面附近區域在剪切過程中主要表現為摩擦破壞.充填物在剪切過程中產生變形導致顆粒粘結斷裂.同時,代表充填物質的顆粒在剛度和強度上都遠小于巖石.因此,剪切過程中節理面間的充填物也將在一定程度上起到“潤滑”的作用,節理面滑移也更為明顯.隨著法向荷載的增大,節理模型法向移動受到更大的限制,即使軟弱充填物質能增進節理剪切作用下的滑動程度,但較高法向荷載下節理面間的吻合系數較高,從而剪切荷載下節理面主要表現為剪斷與摩擦的復合式破壞.

圖7為JRC為14～16與18～20的充填節理在不同法向荷載下的細觀粘結力分布情況.圖中的黑色區域代表壓應力,紅色部分表示剪切應力.

從圖7中可以發現,不同粗糙度的充填節理在剪切過程中的接觸點位置明顯不同.JRC為14～16的充填節理接觸點主要集中在模型右側,該處的巖石節理面吻合度較高,此處最先產生破壞.相比于14～16的節理面,JRC為18～20的節理面接觸位置分布處于左側.當法向荷載為3 MPa時,JRC為14～16的節理接觸位置并未發生明顯變化,只是在法向荷載較大時,接觸面積和接觸程度明顯增加,此現象從圖7(b)所示的力鏈分布規律中也能看出.與14～16的節理類似,法向荷載較大時,JRC為18～20的節理面接觸程度更大,壓應力與剪應力分布區域明顯增加,且接觸位置與分布更為均勻,模型中節理面上存在多個接觸點.在剪切作用下,接觸點附近的巖石將被剪斷,之后隨著充填物一起產生摩擦作用.

(a) JRC為14^16(b) JRC為18^20圖6 法向荷載為1MPa和3MPa下的節理損傷情況Fig.6 Joint morphology damage at normal stresses of 1 MPa and 3 MPa

(a) JRC為14^16(b) JRC為18^20圖7 法向荷載為1 MPa和3 MPa下的細觀粘結力分布Fig.7 Distribution of bond stress at normal stresses of 1 MPa and 3 MPa

2.2 宏細觀參數分析

為了研究巖石節理間充填物性質及JRC對充填節理力學性質的影響,進行3種不同因素下的節理模型直剪數值實驗.影響因素包括:JRC、充填節理顆粒強度、充填物與巖石接觸面強度.其中JRC分為2～4、6～8、10～12、14～16和18～20共5個區間.在顆粒流中,顆粒的強度參數包括粘結強度、顆粒摩擦因數.而本文中的節理顆粒強度比是代表節理充填物顆粒的剛度、粘結強度、摩擦因數與代表完整巖樣的顆粒所對應的強度參數的比值.本文中5種節理顆粒強度摩擦因數比分別為:0.15、0.30、0.45、0.60、0.75.接觸面是指充填物質與圍巖接觸面,而接觸面粘結強度為節理顆粒與圍巖顆粒間的粘結強度與完整巖樣顆粒間粘結強度的比值.不同組內的相關細觀參數如表2所示.

表2 各工況組內的影響因素及數值Tab.2 Factors influencing the simulation cases

節理模型中原巖的細觀參數是根據實際巖石宏觀力學性質進行標定的.充填物質的強度則是根據原巖參數不斷弱化確定的,進行這一設置主要是考慮到實際工程當中節理間的充填物可能是低粘結的類砂質碎屑,也有可能是完整的侵入巖體,只是相比之下強度稍弱.所以進行數值模擬時，基于原巖強度參數對節理處顆粒參數進行弱化處理,且弱化處理過程中將粘結及摩擦因數進行同等弱化.

根據表2中A組所示的細觀參數對充填節理力學性質開展研究,主要分析節理粗糙程度對峰值剪切強度及剪切強度參數的影響.圖8為不同JRC影響下的節理剪切力學性質,由圖8(a)可以看出:隨法向荷載不斷增加,節理的峰值剪切強度呈線性增長,與室內實驗結果相似[6];在同等法向荷載下,隨著節理JRC的不斷增加,節理模型的峰值剪切強度呈現出明顯的增長趨勢;隨著節理粗糙度的不斷增大其抗剪強度也不斷增大.圖8(b)為節理抗剪強度參數隨著JRC的演化情況,隨節理JRC不斷增大,抗剪強度參數Cj(粘聚力)和φj(內摩擦角)均呈增大的趨勢,與黃飄[17]試驗結果的規律相同,相比而言Cj值的增長幅度較小,在JRC值超過12后,φj值出現明顯的增長.

根據表2中B組進行數值計算可以得出節理間充填物強度對其力學性質的影響,圖9為充填物強度對節理剪切力學性質的影響情況.與圖8類似,圖9(a)為不同法向荷載下的節理峰值剪切強度.不同的充填強度比下,隨著法向荷載增大,各組節理峰值剪切強度呈線性變化.當節理所處法向荷載一定時,隨充填物強度的增加其峰值剪切強度也出現明顯的增長趨勢.節理充填強度比為0.15時峰值剪切強度最小,0.75時最大.圖9(b)為不同充填強度對節理抗剪強度參數Cj和φj的影響.隨充填物強度比的增加,粘結力出現了明顯的增長,從0.15時的0.92 MPa增長至0.75時的4.3 MPa.由此可見,充填物強度比對節理面剪切粘結力存在較大的影響.隨著充填強度比的增加,節理的抗剪強度參數φj呈先下降后增加的趨勢.

當節理強度比為0.3時,節理面的φj處于最低值.隨著充填強度比的不斷增加,與參數值Cj一樣,內摩擦角φj同樣呈明顯的增長趨勢.表2中C組所示即為改變接觸面粘結強度時各個參數的分布和變化情況.其中接觸面強度是變量,而JRC固定為2～4,節理顆粒強度比定為 0.3.

圖10中為充填物與巖石接觸面粘結情況對峰值剪切應力及抗剪強度參數的影響.隨著法向荷載的不斷增長,JRC為2～4的節理面峰值剪切應力明顯增長.反觀接觸面粘結強度對峰值剪切強度的影響,當法向荷載較小時,不同接觸面的粘結強度比下峰值剪切強度差距較大,可見較低法向荷載下接觸面粘結強度比對峰值剪切應力的影響較大.當法向荷載較大時,接觸面粘結強度對峰值剪切應力的影響程度明顯降低.主要是由于法向荷載較大時上下節理面被充分壓縮,接觸程度較大,剪切荷載作用下剪切應力主要受到粗糙度的控制和影響.

(a) 法向荷載對峰值強度的影響(b) JRC對粘聚力的影響(c) JRC對內摩擦角的影響圖8 JRC對節理力學性質的影響Fig.8 Influence of JRC value on the mechanical properties of the fill joint

(a) 法向荷載對峰值強度的影響(b) 充填體強度比對粘聚力的影響(c) 充填體強度比對內摩擦角的影響圖9 節理充填物強度比對節理力學性質的影響Fig.9 Influence of the filling strength ratio on the mechanical properties of the fill joint

(a) 法向荷載對峰值強度的影響(b) 接觸面粘結強度比對粘聚力的影響(c) 接觸面粘結強度比對內摩擦角的影響圖10 接觸面強度比對節理力學性質的影響Fig.10 Influence of the interface strength ratio on the mechanical properties of the fill joint

從圖10(b)可以看出,節理面抗剪強度參數Cj隨接觸面粘結的增加同樣呈增長的趨勢,而抗剪強度參數φj隨接觸面粘結強度增大而減小,與Cj的變化規律相反.由此可見,充填節理接觸面粘結強度對節理抗剪強度參數存在不同程度的影響.

除前述的幾種影響因素外,充填厚度也對充填節理剪切力學性質產生極為重要的影響.表2中的D組即為5種不同充填厚度下的直剪模擬方案.此組模擬方案當中,節理粗糙度選為JRC處于10～12的節理面.節理顆粒強度比為0.3,接觸面粘結強度比為1.0.充填厚度作為唯一變量,分別為 0.4、1.4、3.0、4.0 mm以及5.0 mm.

圖11中為不同充填厚度下的節理峰值剪切強度、抗剪強度參數以及粘結破壞數量分布規律.

從圖11(a)中看出:隨著法向荷載不斷增加,各充填厚度下的節理峰值抗剪強度依然呈線性變化;隨著充填厚度的增加,節理峰值抗剪強度越小;充填厚度為0.4 mm時,各法向荷載下的峰值抗剪強度處于最大值,而厚度為5 mm時均處于最小值;隨著充填厚度的增加,峰值抗剪強度并未線性遞減.

圖11(b)中所示為充填厚度的變化對抗剪強度參數的影響.隨著充填厚度的增加,節理面抗剪強度參數存在明顯的劣化趨勢.節理面抗剪強度參數Cj與φj均隨著充填厚度的增加而減小.節理面抗剪強度參數Cj與填厚度呈類線性關系,而φj呈非線性遞減的關系.充填節理的力學性質受充填厚度的影響較大,還要考慮到節理面起伏度的影響.對于JRC為10～12的模型而言,當充填厚度為5 mm 時,充填度大于100%,此時隨著充填厚度的增加節理強度主要受充填體本身強度的影響.因此,圖11中不論是峰值強度還是抗剪強度參數Cj與φj,當充填厚度從4 mm增加至5 mm時,強度參數變化程度不如前面幾種厚度變化顯著.

(a) 法向荷載對峰值強度的影響(b) 充填體厚度對粘聚力的影響(c) 充填體厚度對內摩擦角的影響圖11 充填體厚度對節理力學性質的影響Fig.11 Influence of filling thickness on the mechanical properties of the fill joint

3 結 論

利用二維離散元程序建立不同粗糙程度的巖石充填節理模型,分析JRC、充填強度、接觸面粘結強度、充填厚度等參數對節理剪切強度的影響,主要結論如下:

(1) 法向荷載較大時,上下節理面接觸點分布較為均勻,節理峰值抗剪應力更大,節理面損傷程度也愈發明顯.

(2) 隨著JRC的增加,節理峰值剪切強度明顯增長,而當法向荷載較大時,這種影響將被放大.節理面抗剪強度參數Cj隨著JRC的增加呈類線性增長,而φj隨著JRC的增大呈非線性增長.

(3) 節理峰值剪切應力隨節理充填物的強度增加而增大.隨充填物強度比的增加,粘結力Cj出現明顯的增長趨勢,隨著充填強度比的增加,節理的抗剪強度參數φj呈先下降后增加.

(4) 法向荷載較小時,接觸面粘結強度比對剪切強度的影響較大;當法向荷載較大時,接觸面粘結強度對剪切強度的影響程度明顯降低.抗剪強度參數Cj隨接觸面粘結強度的增加呈現增長的趨勢,而抗剪強度參數φj隨著接觸面粘結強度增大而減小.

(5) 充填厚度對節理剪切力學特性存在明顯的影響,隨著充填厚度的增加,峰值強度、抗剪強度參數Cj和φj逐漸減小.當充填厚度增加至一定值后,充填節理剪切力學性質主要受充填體本身的強度影響.

致謝:中南大學“創新驅動計劃”項目(2016CX019);深部巖土力學與地下工程國家重點實驗室開放基金項目(SKLGDUEK1405).

[1] 夏才初,孫宗頎. 工程巖體節理力學[M]. 同濟大學出版社,2002:1-2.

[2] YANG S Q, XU T, H E L, et al. Numerical study on failure behavior of brittle rock specimen containing pre-existing combined flaws under different confining pressure[J]. Archives of Civil and Mechanical Engineering, 2015, 15(4): 1085-1097.

[3] 蔣宇靜,王剛,李博,等. 巖石節理剪切滲流耦合試驗及分析[J]. 巖石力學與工程學報,2007,26(11): 2253-2259.

JIANG Yujing, WANG Gang, LI bo, et al. Experimental study and analysis of shear-flow coupling behaviors of rock joints[J]. Chinese Journal of Rock Mechanicss and Engineering, 2007, 26(11): 2253-2259.

[4] 周喻,MISAR A,吳順川,等. 巖石節理直剪試驗顆粒流宏細觀分析[J]. 巖石力學與工程學報,2012,31(6): 1245-1256.

ZHOU Yu, MISRA A, WU Shunchuan, et al. Macro-and meso-analysis of rock joint direct shear test using particle flow theory[J]. Chinese Journal of Rock Mechanics and Engineering, 2012, 31(6): 1245-1256.

[5] 唐雯鈺,林杭. 不同鋸齒高度對軟弱結構面剪切特性的影響[J]. 中南大學學報:自然科學版,2017,48(5): 1300-1307.

TANG Wenyu, LIN Hang. Influence of dentate discontinuity height on shear properties of soft structure plane[J]. Journal of Central South University: Science and Technology, 2017, 48(5): 1300-1307.

[6] 曹平,范祥,蒲成志,等. 節理剪切試驗及其表面形貌特征變化分析[J]. 巖石力學與工程學報,2011,30(3): 480-485.

CAO Ping, FAN Xiang, PU Chengzhi, et al. Shear test of joint and analysis of morohology characteristic evolution of joint surface[J]. Chinese Journal of Rock Mechanics and Engineering, 2011, 30(3): 480-485.

[7] GHAZVINIAN A, AZINFAR MJ, NOROZI P. Mechanical response of discontinuities of different joint wall contact strengths[J]. Arabian Journal of Geosciences, 2013, 6(11): 4419-4426.

[8] 陳必光,宋二祥,程曉輝. 二維裂隙巖體滲流傳熱的離散裂隙網絡模型數值計算方法[J]. 巖石力學與工程學報,2014(1): 43-51.

CHEN Biguang, SONG Erxiang, CHENG Xiaohui. A numerical method for discrete fracture network model for flow and heat transfer in two-dimensional fractured rocks[J]. Chinese Journal of Rock Mechanics and Engineering, 2014(1): 43-51.

[9] 杜守繼,職洪濤,周枝華. 巖石節理剪切過程中應力與滲流特性的數值模擬[J]. 巖石力學與工程學報,2008,27(12): 2473-2481.

DU Shouji, ZHI Hongtao, ZHOU Zhihua. Numerical simulation of stress and seepage properties of rock joint during shearing processes[J]. Chinese Journal of Rock Mechanics and Engineering, 2008, 27(12): 2473-2481.

[10] BAHAADDINI M, SHARROCK G, HEBBLEWHITE B K. Numerical investigation of the effect of joint geometrical parameters on the mechanical properties of a non-persistent jointed rock mass under uniaxial compression[J]. Computers and Geotechnics, 2013, 49(1): 206-225.

[11] 吳順川,周喻,高斌. 卸載巖爆試驗及PFC3D數值模擬研究[J]. 巖石力學與工程學報,2010,29(S2): 4082-4088.

WU Shunchuan, ZHOU Yu, GAO Bin. Study of unloading tests of rock burst and PFC3D numerical simulation[J]. Chinese Journal of Rock Mechanicss and Engineering, 2010, 29(S2): 4082-4088.

[12] MANOUCHEHRIAN A, SHARIFZADEH M, MARJI M F, et al. A bonded particle model for analysis of the flaw orientation effect on crack propagation mechanism in brittle materials under compression[J]. Archives of Civil and Mechanical Engineering, 2014, 14(1): 40-52.

[13] 黃達,岑奪豐. 單軸靜-動相繼壓縮下單裂隙巖樣力學響應及能量耗散機制顆粒流模擬[J]. 巖石力學與工程學報,2013(9): 1926-1936.

HUANG Da, CEN Duofeng. Mechanical responses and energy dissipation mechanism of rock specimen with a single fissure under static and dynamic unaxial compression using particle flow code simulations[J]. Chinese Journal of Rock Mechanics and Engineering, 2013(9): 1926-1936.

[14] PARK J W, SONG J J. Numerical simulation of a direct shear test on a rock joint using a bonded-particle model[J]. International Journal of Rock Mechanics and Mining Sciences, 2009, 46(8): 1315-1328.

[15] SARFARAZI V, GHAZVINIAN A, SCHUBERT W. Numerical simulation of the process of fracture of echelon rock joints[J]. Rock Mechanics and Rock Engineering, 2014, 47(4): 1355-1371.

[16] Itasca Consulting Group. PFC2D (partical flow code in 2 dimensions) theory and background[M]. Minnesota: Itasca Consulting Group Inc, 2002: 77-78.

[17] 黃飄. 不同粗糙度結構面試件剪切試驗特性研究[J]. 建筑監督檢測與造價,2015,8(3): 1-5.

HUANG Piao. Study on shear test properties of different structural plane roughness[J]. Supervision Test and Cost of Construction, 2015, 8(3): 1-5.