鐘寄宸， 王 蕉*， 李翔宇,2， 楚錫華

(1.西南交通大學(xué) 四川省應(yīng)用力學(xué)與結(jié)構(gòu)安全重點(diǎn)實(shí)驗(yàn)室,力學(xué)與航空航天學(xué)院,成都 610031；2.西南交通大學(xué) 牽引動(dòng)力國(guó)家重點(diǎn)實(shí)驗(yàn)室應(yīng)用力學(xué)研究所,成都 610031；3.武漢大學(xué) 工程力學(xué)系，武漢 430072)

1 引 言

在宏觀尺度上，頁(yè)巖是典型的橫觀各向同性材料，其破壞通常表現(xiàn)為脆性破壞，難以監(jiān)測(cè)。研究者有必要對(duì)頁(yè)巖的力學(xué)特性進(jìn)行測(cè)試，以確定頁(yè)巖的各向異性對(duì)力學(xué)性能、斷裂模式和抗壓能力的影響因素。對(duì)頁(yè)巖或類(lèi)似頁(yè)巖的脆性材料的儲(chǔ)層結(jié)構(gòu)和力學(xué)性能的清晰認(rèn)識(shí)將為工程實(shí)踐提供啟示。

層狀巖石的單軸壓縮試驗(yàn)表明，巖石的抗壓強(qiáng)度隨著層理角的增加而下降[1-3]。層理面對(duì)沉積巖力學(xué)行為影響的巴西圓盤(pán)試驗(yàn)研究表明各向異性導(dǎo)致了復(fù)雜的斷裂力學(xué)行為[4，5]。單軸壓縮試驗(yàn)和巴西試驗(yàn)結(jié)合起來(lái)研究層狀巖石的力學(xué)行為,結(jié)果表明頁(yè)巖的單軸抗壓強(qiáng)度UCS、彈性模量和破壞模式受層理面影響顯著。衡帥等[8]總結(jié)了龍馬溪黑頁(yè)巖的四種典型破壞模式;Yang等[9]也總結(jié)了單軸壓縮條件下，巖石類(lèi)脆性材料的四種典型破壞模式。此外，也有研究探討了頁(yè)巖各向異性對(duì)井壁斷裂的影響[10,11],結(jié)果顯示，樣品的抗壓強(qiáng)度曲線(xiàn)呈U型。頁(yè)巖的力學(xué)行為很大程度上取決于層理面的粗糙度。為了進(jìn)一步研究頁(yè)巖的各向異性，Ge等[12]基于分形方法研究了巖石節(jié)理粗糙度的測(cè)量和量化。以上學(xué)者通過(guò)物理實(shí)驗(yàn)對(duì)頁(yè)巖的力學(xué)性能進(jìn)行了研究，并提供了直觀的物理現(xiàn)象和可靠的實(shí)驗(yàn)數(shù)據(jù)。但是物理實(shí)驗(yàn)受實(shí)驗(yàn)條件和實(shí)驗(yàn)設(shè)備的限制較多，且成本巨大。相比于傳統(tǒng)物理實(shí)驗(yàn)，數(shù)值方法研究成本較低且能收集到更多數(shù)據(jù)，所以很多學(xué)者采用數(shù)值方法對(duì)頁(yè)巖進(jìn)行研究。

數(shù)值模擬是巖石力學(xué)和巖土工程的一個(gè)重要研究方法。在巖土力學(xué)中，廣泛使用Cundall等[13]提出的離散元法DEM為基礎(chǔ)的離散元模擬。采用嵌入光滑節(jié)理接觸的粘結(jié)顆粒離散元模型模擬橫觀各向同性巖石的力學(xué)行為[14-17]。在此基礎(chǔ)上，Xia等[18]研究了光滑節(jié)理的參數(shù)和校準(zhǔn)方法，本文采用了同樣的方法。Xu等[19]提出了一個(gè)新的模型，用以模擬各向異性巖石的力學(xué)行為。其中，巖石基體和層理面分別使用平節(jié)理接觸和光滑節(jié)理接觸進(jìn)行模擬。Lin等[20]用離散元法研究了層理面參數(shù)對(duì)雙孔層狀巖石力學(xué)性能的影響。實(shí)驗(yàn)結(jié)果表明巖石的抗壓強(qiáng)度和彈性模量隨層理面的變化呈中間低兩端高的形態(tài)。楊宇江等[21]引入多重分形理論，計(jì)算了單軸加載過(guò)程中峰前某一階段單元體積應(yīng)變的多重分形譜。除了壓縮試驗(yàn)，也有少量關(guān)于頁(yè)巖拉剪行為的研究[23]，發(fā)現(xiàn)試樣在混合形式的荷載條件下，裂縫通常以拉剪I/II型混合模式起裂，而后以I型為控制方式穩(wěn)定擴(kuò)展。上述學(xué)者通過(guò)數(shù)值方法研究了巖石的力學(xué)特性，極大地豐富了各向異性巖石力學(xué)響應(yīng)的數(shù)據(jù)庫(kù)。但是所做的貢獻(xiàn)主要集中在完整巖石的破壞表征和力學(xué)性能等方面，對(duì)于細(xì)觀損傷的研究尚不充分。

為了更好地理解和表征巖石內(nèi)部的細(xì)觀損傷，本文采用顆粒流模型分析不同層理角的頁(yè)巖樣品在單軸壓縮下的破壞情況，引入裂紋張量和裂紋密度，研究微裂紋隨加載過(guò)程的演化以及層理角對(duì)裂紋密度的影響。此外，使用平均配位數(shù)，分析配位數(shù)與試樣應(yīng)力狀態(tài)及細(xì)觀損傷的聯(lián)系。本文在細(xì)觀尺度上定義頁(yè)巖類(lèi)脆性材料的內(nèi)部損傷,并將其與宏觀應(yīng)力應(yīng)變相結(jié)合,對(duì)橫觀各向同性脆性巖石的破壞過(guò)程和內(nèi)在機(jī)理的理解提供參考。

2 PFC2D數(shù)值模型

2.1 接觸模型概述

由于實(shí)驗(yàn)手段的限制，在物理實(shí)驗(yàn)過(guò)程中，很難直接觀察材料內(nèi)部的細(xì)觀運(yùn)動(dòng)。為進(jìn)一步研究巖石破壞過(guò)程中材料的內(nèi)部損傷演化過(guò)程，本文采用基于DEM的顆粒流模型，每個(gè)樣本使用三種接觸模型，其中，線(xiàn)性平行粘結(jié)模型用于模擬巖石基體，光滑節(jié)理模型用于模擬層理面。當(dāng)顆粒間的力超過(guò)粘結(jié)強(qiáng)度時(shí)，粘結(jié)發(fā)生斷裂，顆粒之間的接觸行為由線(xiàn)性接觸模型描述。此外，摩擦力在這三種接觸模型中均存在，三種模型描述如下。

線(xiàn)性模型中接觸力的計(jì)算如下。

法向力Fn和切向力Fs為

Fn=Fn+knAΔδn,Fs=Fs-ksAΔδs

(1,2)

式中A為等效面積(與接觸顆粒半徑相關(guān))，kn和ks分別為線(xiàn)性模型法向和切向剛度，Δδn和Δδs分別為相對(duì)法向和相對(duì)切向位移增量。

平行粘結(jié)模型中接觸力的計(jì)算如下。

(3，4)

光滑節(jié)理模型中接觸力的計(jì)算如下。

(5，6)

2.2 頁(yè)巖數(shù)值模型描述

將層理角β從0調(diào)整到90°(β定義為層理面與x軸的夾角)，建立7個(gè)層理面厚度恒定的數(shù)值模型(層理面厚度db=1.4×rh i)，如圖4所示?？紤]到真實(shí)巖石不同層的厚度不盡相同，因此頁(yè)巖各層的厚度設(shè)置略有不同。

h1=16.5 mm,h2=21 mm,h3=19.5 mm
h4=21 mm,h5=18 mm,h6=21 mm
h7=24 mm

對(duì)層理角為0°～90°的矩形試樣(75 mm×150 mm)進(jìn)行單軸壓縮測(cè)試。顆粒與壁面的摩擦系數(shù)設(shè)為0.1。模擬中的每個(gè)數(shù)值模型大約由40300個(gè)顆粒組成。顆粒半徑遵循0.21 mm～0.35 mm 的均勻分布。在模擬過(guò)程中，上壁用于加載，下壁固定。根據(jù)Zhang等[24]的研究，加載速度設(shè)定為0.05 m/s(≤0.08 m/s)，時(shí)間步長(zhǎng)設(shè)定為2.5×10-8s，以確保整個(gè)加載過(guò)程是準(zhǔn)靜態(tài)的。

2.3 頁(yè)巖材料參數(shù)校正

表1 層理面的材料參數(shù)

表2 巖石基質(zhì)的材料參數(shù)

圖1 線(xiàn)性模型

2.3.1 破壞模式對(duì)比

將數(shù)值模擬結(jié)果與文獻(xiàn)[8]的物理實(shí)驗(yàn)的破壞模式進(jìn)行對(duì)比，如圖2所示。圖2(a1～d1)顯示了頁(yè)巖在單軸壓縮時(shí)的破壞模式，圖2(a2～d2)對(duì)應(yīng)為數(shù)值模擬加載至模型破壞時(shí)的應(yīng)變分布(應(yīng)變?cè)茍D計(jì)算過(guò)程參見(jiàn)文獻(xiàn)[25])，兩種方法得到的破壞模式具有相似性。圖2(a)的β=0；圖2(b)的β=30°，巖石中產(chǎn)生多個(gè)貫穿層理面和沿層理面的破壞面；圖2(c)的β=60°，巖石中產(chǎn)生了沿層理面的破壞面，其中，貫穿性破壞沿著60°的層理面，導(dǎo)致巖石試樣發(fā)生明顯的剪切滑移破壞；圖2(d)的β=90°，巖石中產(chǎn)生多個(gè)貫穿層理面和沿層理的破壞面，以后者為主。

圖2 平行粘結(jié)模型

除物理實(shí)驗(yàn)對(duì)應(yīng)的4組層理角外，對(duì)層理角為15°，45°和75°的情況進(jìn)行了數(shù)值模擬，破壞模式如 圖3 所示。結(jié)合圖5，對(duì)破壞模式進(jìn)行分析，隨著層理角從0°增加到90°，巖石的破壞模式完成了從基質(zhì)和層理面的混合破壞到沿層理面的單一破壞再到基質(zhì)和層理面的混合破壞的演化，如圖6所示。

圖3 光滑節(jié)理模型

圖4 不同層理角度頁(yè)巖的數(shù)值模型(β =0°和β =30°)

圖5 實(shí)驗(yàn)和數(shù)值模擬的破壞模式對(duì)照

圖6 不同層理角頁(yè)巖的破壞模式

2.3.2 宏觀力學(xué)性能對(duì)照

將物理實(shí)驗(yàn)和數(shù)值模擬所得的抗壓強(qiáng)度和彈性模量進(jìn)行比較，結(jié)果如圖7所示?？紤]到文獻(xiàn)[8，26]對(duì)于頁(yè)巖抗壓強(qiáng)度描述的差異，本文將文獻(xiàn)[8]的物理實(shí)驗(yàn)測(cè)得的抗壓強(qiáng)度和數(shù)值模擬的抗壓強(qiáng)度做歸一化處理以對(duì)比抗壓強(qiáng)度的變化規(guī)律。具體數(shù)值列入表3和表4。

表3 數(shù)值模型力學(xué)性能

表4 不同層理角度頁(yè)巖單軸壓縮試驗(yàn)結(jié)果

圖7 實(shí)驗(yàn)和數(shù)值模擬的對(duì)照

圖7(a)顯示了物理實(shí)驗(yàn)和數(shù)值模擬中不同層理角頁(yè)巖的抗壓強(qiáng)度。數(shù)值模擬結(jié)果顯示為虛線(xiàn)，文獻(xiàn)[8]物理實(shí)驗(yàn)的結(jié)果顯示為實(shí)線(xiàn)。很明顯，從0°～60°，頁(yè)巖的抗壓強(qiáng)度隨著層理角的增加而降低，這種下降在30°～90°區(qū)間內(nèi)更明顯。這一現(xiàn)象與層理面的弱化效應(yīng)有關(guān)。隨著層理角的增加，層理面兩側(cè)的基質(zhì)更容易沿著層理面產(chǎn)生相對(duì)運(yùn)動(dòng)。由于層理面內(nèi)顆粒的粘結(jié)力和抵抗變形的能力相對(duì)較弱，這些區(qū)域更容易受到破壞。因此，頁(yè)巖的抗壓強(qiáng)度下降。在60°～90°范圍內(nèi)，頁(yè)巖的抗壓強(qiáng)度隨著層理角的增加而明顯增加。

當(dāng)層理角為75°時(shí)，巖石的抗壓強(qiáng)度增加。因此，在60°～75°區(qū)間可能存在某個(gè)角度閾值。當(dāng)層理角超過(guò)這個(gè)閾值時(shí)，層間弱化效應(yīng)對(duì)巖石抗壓強(qiáng)度的影響減少，抗壓強(qiáng)度增加，頁(yè)巖的抗壓強(qiáng)度在層理面的影響下顯示出明顯的各向異性。閾值大小還與層理面粗糙度和顆粒排列方式等密切相關(guān)，其具體作用機(jī)理還需要進(jìn)一步研究。

圖7(b)顯示了不同層理角頁(yè)巖的彈性模量。物理實(shí)驗(yàn)和數(shù)值結(jié)果都表明，頁(yè)巖的彈性模量隨著層理角的增加而增加，在75°～90°區(qū)間內(nèi)變化更為明顯。與抗壓強(qiáng)度類(lèi)似，彈性模量在層理面的影響下也顯示出明顯的各向異性，數(shù)值模擬結(jié)果與文獻(xiàn)[8]物理實(shí)驗(yàn)有良好的一致性。

3 細(xì)觀分析

3.1 裂紋張量定義

(7)

(8)

式中d為裂紋兩側(cè)顆粒的平均直徑，ni和nj分別為裂紋法線(xiàn)的i和j方向的分量，N為裂紋總數(shù)量，A為顆粒的總面積。

3.2 層理面對(duì)裂紋密度的影響

圖8顯示了不同層理角的頁(yè)巖樣本中裂紋密度隨軸向應(yīng)變的關(guān)系曲線(xiàn)，可以看出，當(dāng)層理角不大于60°時(shí)，裂紋密度隨著層理角的增加而增加，這與層理面的弱化效應(yīng)有關(guān)。隨著層理角的增加，層理面兩側(cè)的基體更容易沿著層理面產(chǎn)生相對(duì)運(yùn)動(dòng)。因此，在相同的荷載條件下，裂紋密度隨著層理角的增加而增加。

圖8 不同層理角頁(yè)巖的裂紋密度

層理角等于60°時(shí)，當(dāng)應(yīng)變從0變化到1.33×10-3，裂紋密度迅速增加。繼續(xù)加載，當(dāng)應(yīng)變從 1.33×10-3變化到2.16×10-3，裂紋密度曲線(xiàn)的斜率沒(méi)有明顯變化?？梢杂^察到應(yīng)力集中的現(xiàn)象大多發(fā)生在層理面附近，大量裂紋沿著層理面產(chǎn)生，此時(shí)，基體中的應(yīng)力較小，所以在基體中產(chǎn)生的微裂紋數(shù)量較少。裂紋密度的增長(zhǎng)主要由沿層理面的粘結(jié)破壞引起，當(dāng)應(yīng)變達(dá)到2.16×10-3時(shí)，巖石完全破壞，失去原有的承載能力。

層理角大于60°時(shí)，當(dāng)應(yīng)變從0到1.33×10-3變化時(shí)，裂紋密度增長(zhǎng)緩慢，當(dāng)應(yīng)變從1.33×10-3到2.16×10-3變化時(shí)，裂紋密度迅速增加。這一現(xiàn)象表明，巖石的破壞機(jī)制發(fā)生了改變，在60°～75°之間存在某個(gè)層理角閾值。當(dāng)層理角超過(guò)閾值時(shí)，破壞模式將從主要由層理面控制變?yōu)閷永砻婧蛶r石基體共同控制，層理面弱化效應(yīng)的影響減弱。因此，在初始加載階段，β=75°的頁(yè)巖的抗破壞能力比β=60°的強(qiáng)。當(dāng)層理角增加到90°時(shí)，巖石在初始加載階段抵抗變形的能力進(jìn)一步增強(qiáng)。

3.3 平均配位數(shù)分析

對(duì)15°的頁(yè)巖樣本展開(kāi)進(jìn)一步分析。樣本的應(yīng)力和裂紋數(shù)量曲線(xiàn)如圖9所示。根據(jù)裂紋數(shù)量的變化，整個(gè)加載過(guò)程可以分為三個(gè)階段，第I階段，裂紋數(shù)量緩慢上升；第II階段，裂紋數(shù)量增長(zhǎng)加快；第III段，裂紋的數(shù)量迅速上升，試樣迅速發(fā)生破壞。對(duì)照試樣的應(yīng)力曲線(xiàn)，在I和II階段一直保持上升狀態(tài)，僅通過(guò)應(yīng)力-應(yīng)變曲線(xiàn)，無(wú)法反映出試樣內(nèi)部的細(xì)觀損傷。

圖9 應(yīng)力與裂紋數(shù)量曲線(xiàn)

為了更好地表征巖石試樣的細(xì)觀損傷，引入平均配位數(shù)，分析配位數(shù)為3～7的顆粒數(shù)量和平均配位數(shù)的變化。平均配位數(shù)可表示為

(9)

式中Cave為平均配位數(shù)，N為顆?？倲?shù)量，Ck是id為k的顆粒的配位數(shù)。

本文使用的數(shù)值模型平均配位數(shù)約為 4.6，大于平均配位數(shù)的配位數(shù)定義為大配位數(shù)(5,6,7)，反之為小配位數(shù)(3,4)。從加載初時(shí)刻起，每間隔1.25 ms統(tǒng)計(jì)一次配位數(shù)，共統(tǒng)計(jì)12次，不同配位數(shù)顆粒數(shù)量和平均配位數(shù)曲線(xiàn)經(jīng)過(guò)歸一化處理。0代表12組取樣點(diǎn)的最小值，1代表12組取樣點(diǎn)的最大值。從圖10可以看出，巖石試樣的小配位數(shù)顆粒數(shù)量先下降后上升，大配位數(shù)顆粒數(shù)量先上升后下降，平均配位數(shù)先上升然后伴隨較小幅度的下降，然后急劇下降。小配位數(shù)顆粒數(shù)量曲線(xiàn)及大配位數(shù)顆粒數(shù)量曲線(xiàn)的轉(zhuǎn)折點(diǎn)都在臨界狀態(tài)附近(黑色豎線(xiàn)代表臨界狀態(tài)出現(xiàn)位置，即應(yīng)力峰值點(diǎn)附近)。

圖10 不同配位數(shù)顆粒數(shù)量及平均配位數(shù)的變化曲線(xiàn)

圖11顯示了加載過(guò)程中平均配位數(shù)與裂紋數(shù)量的對(duì)應(yīng)關(guān)系，可以看出，平均配位數(shù)和裂紋數(shù)量變化趨勢(shì)有很好的一致性。從初始狀態(tài)到臨界狀態(tài)，平均配位數(shù)和裂紋數(shù)量都緩慢上升。超過(guò)臨界狀態(tài)后，裂紋數(shù)量迅速上升，平均配位數(shù)迅速下降。在樣本發(fā)生破壞之前(即I和II階段)，平均配位數(shù)和裂紋數(shù)量的變化可以相互佐證，I階段，平均配位數(shù)穩(wěn)步上升，對(duì)應(yīng)裂紋數(shù)量從無(wú)到有，并線(xiàn)性增加；II階段，平均配位數(shù)有小幅下降，對(duì)應(yīng)裂紋數(shù)量繼續(xù)增加，但增速明顯大于I階段。最后樣本發(fā)生破壞，平均配位數(shù)急劇下降，對(duì)應(yīng)裂紋數(shù)量急劇上升。因此，配位數(shù)能夠在一定程度上反映試樣內(nèi)部細(xì)觀結(jié)構(gòu)的變化。

圖11 平均配位數(shù)與裂紋數(shù)量曲線(xiàn)

單軸壓縮過(guò)程中，巖石顆粒受到擠壓。在外力的作用下，部分小配位數(shù)顆粒的穩(wěn)定結(jié)構(gòu)破壞，顆粒之間變得緊密以維持內(nèi)部的穩(wěn)定狀態(tài)，配位數(shù)上升。此時(shí)顆粒之間的相互作用更加緊密，接觸力上升；當(dāng)接觸力進(jìn)一步增大，試樣達(dá)到臨界狀態(tài)，大配位數(shù)顆粒的接觸結(jié)構(gòu)發(fā)生脆性破壞，迅速分裂成多個(gè)小配位數(shù)顆粒的接觸結(jié)構(gòu)(對(duì)應(yīng)著大配位數(shù)顆粒的減少和小配位數(shù)顆粒的增多)，裂紋大量增加，試樣發(fā)生破壞。因此，顆粒集合體在單軸壓縮條件下的應(yīng)力狀態(tài)及內(nèi)部細(xì)觀損傷程度與平均配位數(shù)的變化關(guān)系是對(duì)應(yīng)的。

4 結(jié) 論

本文通過(guò)離散元模擬再現(xiàn)了不同層理角的頁(yè)巖試樣在單軸壓縮下的破壞模式。在此基礎(chǔ)上，系統(tǒng)地研究了層理角對(duì)巖石破壞過(guò)程的影響。為描述巖石內(nèi)部細(xì)觀損傷，引入裂紋張量和平均配位數(shù)，研究了單軸壓縮條件下，層理角對(duì)巖石內(nèi)部細(xì)觀損傷的影響以及平均配位數(shù)與巖石應(yīng)力狀態(tài)和細(xì)觀損傷的聯(lián)系。結(jié)論如下，隨著層理角從0°增加到90°，巖石破壞模式從基質(zhì)和層理面的混合破壞到沿層理面的單一破壞再到基質(zhì)和層理面的混合破壞。同等條件下，頁(yè)巖的裂紋密度會(huì)隨著層理角的增加而增加。相較于應(yīng)力應(yīng)變曲線(xiàn)，平均配位數(shù)曲線(xiàn)能夠更準(zhǔn)確地反映巖石的細(xì)觀損傷；應(yīng)力應(yīng)變曲線(xiàn)峰值滯后于裂紋數(shù)量曲線(xiàn)和平均配位數(shù)曲線(xiàn)的轉(zhuǎn)折點(diǎn)。盡管在第III階段初，試樣內(nèi)部已經(jīng)出現(xiàn)大量裂紋，但裂紋尚未聚結(jié)產(chǎn)生貫穿性裂紋，所以試樣仍有一定強(qiáng)度，繼續(xù)加載，裂紋貫通，試樣達(dá)到臨界狀態(tài)，迅速發(fā)生破壞。小配位數(shù)顆粒的數(shù)量先下降后上升，大配位數(shù)顆粒的顆粒數(shù)量先上升后下降；顆粒集合體在單軸壓縮條件下的應(yīng)力應(yīng)變及裂紋數(shù)量曲線(xiàn)與平均配位數(shù)曲線(xiàn)有良好的一致性。