硬脆性頁巖斷裂韌性二維數值模擬研究

張景軒，范 曉，陳 波，劉向君，梁利喜，熊 健

(1.西南石油大學油氣藏地質及開發工程國家重點實驗室，四川 成都 610500；2.東方地球物理公司物探技術研究中心成都分中心，四川 成都 610000)

頁巖儲層作為頁巖氣賦存的場所，具有自生自儲的特點，儲層特征表現為極低孔、極低滲。頁巖氣資源的高效開發需要水平井鉆井、分段壓裂改造等一系列關鍵技術[1]。壓裂過程中裂縫起裂和延伸受頁巖斷裂韌性的影響，當裂縫尖端應力強度因子超過頁巖斷裂韌性值時，裂縫開始擴展[2]。因此，研究頁巖斷裂韌性及其影響因素對評價頁巖儲層可壓性有重要意義。

目前，針對巖石斷裂韌性研究，主要采用室內試驗測試和數值模擬方法[3-5]。衡帥[6]開展切口與層理呈不同方位的圓柱形試樣三點彎曲試驗，研究頁巖斷裂韌性的各向異性特征，并揭示其斷裂機制的各向異性。Atkinson等[7]通過直縫巴西圓盤實驗法(SNBD)測量斷裂韌性，并給出了I型和II型斷裂韌性計算公式。崔振東等[8]采用人字形切槽巴西圓盤法(CCNBD)，研究泥質砂巖的斷裂力學特性，討論了該試樣類型的有效尺寸和斷裂機制。張盛等[9]分別以5種不同形狀的圓盤試件測定了大理巖斷裂韌度，提出了用每種圓盤確定斷裂韌度的公式。由此可見，室內測試巖石斷裂韌性的方法較多，但眾多方法都有各自的適應性。國際巖石力學學會(ISRM)推薦人字形切槽巴西圓盤法為巖石斷裂韌性的測試方法[10-11]。因為該方法的制樣簡單，測試結果離散度較低，結果可靠。同時，ISRM也提出了巖石斷裂韌性測試樣品的最小有效試樣直徑為75 mm[12]。研究巖石斷裂韌性的數值方法，包括位移不連續法(DDM)、有限單元法(FEM)、擴展有限元法(XFEM)、離散元法(DEM)和邊界元法(BEM)等[13]。滿軻[14]在LS-DYNA程序的基礎上，進行了玄武巖動態斷裂韌性測試，得到了動態斷裂韌性隨著加載率增大而增大的規律。程遠方等[15]用擬三維模型描述煤層氣井壓裂裂縫擴展行為時提出，斷裂韌性是縫高的重要影響因素。由此可見，巖石斷裂韌性特征的研究已經取得了一定的成果。針對硬脆性頁巖斷裂韌性研究，陳建國等[16]基于SNBD法研究了渝東南酉陽地區的龍馬溪組頁巖斷裂韌性特性，并基于地球物理測井信息，建立了頁巖斷裂韌性的預測模型；梁利喜等[17]研究了四川長寧地區龍馬溪組頁巖斷裂韌性的預測模型，并討論了礦物組成、TOC含量對其影響。趙金洲等[18]基于頁巖儲層斷裂韌性、脆性等特征，提出了表征頁巖氣可壓性的評價方法；張明明[19]開展了三點彎曲半圓盤法測試頁巖Ⅰ-Ⅱ復合型斷裂韌性的實驗和模擬，建立了裂紋擴展角的預測方法。然而，這些研究并未考慮層理角度、層理密度及層理面力學性質等對頁巖斷裂韌性的影響，以及這些因素對頁巖裂縫前緣形態的影響。

在前人的研究基礎上，以龍馬溪組頁巖為研究對象，運用以統計損傷的本構關系與有限單元法結合的RFPA2D巖石破裂過程數值模擬軟件，構建人字形切槽巴西圓盤法數值模擬過程，從層理角度和層理面密度以及層理面力學性質等角度，研究頁巖斷裂韌性的變化規律；在此基礎上，討論頁巖裂縫前緣形態的影響因素。同時，也研究了頁巖斷裂韌性與脆性指數間的關系。

1 模型參數選取與驗證

1.1 數值模型介紹

RFPA2D(Rock Failure Process Analysis)用彈性損傷力學的本構關系來描述巖石細觀單元的力學性質，在考慮到細觀非均勻性的前提下，將巖石的非均勻性參數引入到基元中去，細觀基元的性質是簡單的彈塑性或者彈脆性，當基元應力達到破壞準則時則發生破壞，采用剛度退化處理，達到用連續介質力學方法處理物理非連續性的問題，從而實現巖石宏觀非線性[20]。

按照應變等價原理，應力作用在受損材料上引起的應變與有效應力作用在受損材料引起的應變等價，受損材料單元的本構關系可通過無損材料中的名義應力得到[21]，即：

(1)

式中：E和E0分別為損傷后單元的彈性模量和初始彈性模量，MPa；D為損傷變量，D=0時，表示材料無損傷，而當D=1時，表示材料完全喪失承載能力。

1.2 參數選取

為研究層理面對斷裂韌性以及裂縫前緣形態的影響，進一步探討頁巖斷裂韌性及誘導裂縫前緣形態影響因素，利用RFPA2D構建了層理性頁巖的巴西圓盤法數值模擬。該軟件用統計損傷的本構關系來描述巖石材料的非均勻性和缺陷分布的隨機性，再將細觀單元性質的統計分布與有限單元法結合，以彈性力學為應力分析工具，以彈性損傷理論及修正后的Coulomb破壞準則對單元進行變形及破壞處理，從而實現對非均勻材料破壞過程的模擬，設計了如圖1所示的技術路線。

圖1 參數確定流程

由于RFPA2D系統中所輸入的模型參數不等同于巖石力學實驗得到的巖石宏觀物性參數，而是細觀單元參數的平均值。巖石的孔隙度、滲透率、密度三個參數對裂縫的形態影響很小，在這里假設模型中這幾個參數均勻分布，故本研究以龍馬溪組頁巖基礎物性和力學參數值作為參考，主要確定細觀彈性模量、抗壓強度、泊松比三個力學參數。

首先，以物理實驗獲得的力學參數作為參考，對模型輸入的各細觀力學參數平均值，以及影響該力學參數大小的相關參數進行調節，初步確定力學模型參數。然后，對模型的參數進行驗證，以該套參數建立數值模型，并開展單軸壓縮實驗模擬，得到數值模型的宏觀力學參數，將獲得的力學參數與物理實驗力學參數進行對比。按照上述方法，結合測得的龍馬溪組巖心基本物性參數與力學參數的平均值，對軟件輸入的各參數進行調節，主要是力學參數的調節。最終確定模型參數見表1所示。

以表1中模型參數建立模型，模型尺寸為50 mm×25 mm的標準巖心，取層理面密度的平均值0.45條/mm，將模型共劃分網格400×200個。計算時采用位移加載，單步增量為0.005 mm。其模擬過程及計算結果如圖2所示。

表1 模擬參數

圖2 計算前后模型

1.3 模型驗證

根據應力-應變數據，繪制出單軸壓縮下數值模型的應力-應變曲線，取其中極大值為抗壓強度。取應力應變曲線中彈性階段進行擬合，通過擬合直線的斜率可得到彈性模量、泊松比的數據，如圖3所示。通過計算，單軸壓縮情況下，抗壓強度為98.872 MPa，彈性模量為243.77 MPa，泊松比為0.202 3。

圖3 力學參數計算

對比表2中數據，物理試驗與模擬計算所得結果較吻合，與第一組試驗的抗壓強度誤差為1.9%，彈性模量誤差為1.6%，泊松比誤差為1.8%；與第二組試驗的抗壓強度誤差為2.7%，彈性模量誤差為0.2%，泊松比誤差為1.2%，抗壓強度、彈性模量、泊松比誤差均在5%以內。因此，該模擬參數能較好地反應巖心的真實情況，適用于該龍馬溪頁巖的模擬計算。

表2 物理試驗力學參數與RFPA模擬宏觀力學參數對比

2 實驗結果討論與分析

針對靜態或準靜態I型斷裂韌性測試，國內外學者作了大量研究，測試方法包括：V型切口的三點彎曲圓棒試樣(CB)、短棒試樣(SR)以及人字形切槽巴西圓盤試樣(CCNBD)等。CCNBD試樣在很多方面優于CB和SR試樣[8]：①試樣體積較??；②加工工藝較簡單，免除了加工尖銳直裂紋的困難，也能承受較高的臨界載荷；③不僅適合測試I型裂紋的斷裂韌度值，而且很容易實現純II型或I-II復合型裂紋斷裂模式。故選取人字形切槽巴西圓盤試樣進行斷裂韌性模擬，幾何參數如圖4所示。

圖4 CCNBD試樣示意

采用式(2)計算CCNBD試樣斷裂韌性值[17]：

(2)

由以上CCNBD測試巖石斷裂韌性值的原理可知，斷裂韌性與材料本身的物理力學參數無直接相關，僅與施加于試樣上的荷載以及試樣的幾何參數(D,B,a,a0,a1)有關。又由于模擬研究的各模型幾何參數保持一致，故可用RFPA2D程序計算的巖石破壞應力反映斷裂韌性的變化規律。首先基于巖石三軸實驗模擬，確定出頁巖的三軸力學參數，結合巴西圓盤的幾何參數，分別從層理的角度(即與加載方向的夾角，以10°為間隔依次從0°變化至90°)和層理的密度(以0.2條/mm為間隔依次從0.2條/mm變化至1.0條/mm)以及層理面力學性質(0.25～0.75倍基質強度折減)建立人字形切槽巴西圓盤模型，采用位移加載，單步增量取0.002 mm，模型運算條件如表1所示。

2.1 層理面對頁巖裂縫前緣形態影響

層理面是影響裂縫前緣形態的一個重要因素，由于層理面弱面的影響，裂縫在延伸至層理面時，可能出現轉向、分叉或穿過現象。在模型中劃分不同層理面，分別從層理面角度、層理面線密度、層理面力學性質三方面分析層理面對頁巖裂縫前緣形態的影響。

固定層理的面密度，以0.4條/mm為例，當頁巖層理與加載方向呈不同角度時，其裂縫前緣形態的變化如圖5所示。圖中可見，各個頁巖試樣均在預制裂縫處開始宏觀斷裂，但裂紋的擴展路徑卻差別較大：當層理夾角為0°時，在加載過程中裂縫沿著預制裂縫尖端起裂和延伸，未發生轉向，這可能是因為裂縫和層理都作為弱結構面且都相互平行；當層理角度小于40°時，在加載過程中，裂縫起裂時沿著預制裂縫延伸，當遇到層理后，裂縫并沒有穿過層理而是發生轉向，繼而沿著層理方向延伸，此時，裂縫轉向角近似等于層理角度，裂縫延伸過程中層理起到較重要的作用；當層理夾角在40～70°之間時，除了在延伸時前緣隨層理出現轉向外，圓盤端部出現多條沿層理方向的次生裂紋，這可能是由于在峰值載荷后軸向壓應力驟降，使卸荷后的拉張應力和剪切應力集中于試樣的端部區域，而層理面的膠結強度較小并且力學性質相對薄弱，最終導致圓盤端部沿層理面出現多條次生裂紋；當層理夾角大于70°時，層理在加載方向上的間隔變小，圓盤在軸向應力與拉張應力的作用下，層理面的剪切應力作用減弱，裂縫先沿預置裂縫方向擴展，后穿過層理發生分叉，并且出現次生裂紋的數量明顯減少。

裂紋自尖端起裂后，彈性應變能急劇釋放，表現有峰后載荷驟降、起裂偏離角度、次生裂紋產生等[8]，轉向后的裂縫沿層理面繼續延伸擴展，且在沿層理面繼續延伸的過程中會進一步溝通水力裂縫或天然裂縫而形成相對復雜的裂縫網絡[6]。以層理角度30°、60°、90°為例，分析密度0.2條/mm、0.4條/mm、0.6條/mm的裂縫前緣形態(見圖6)，層理密度對裂縫前緣形態影響規律較復雜：①層理密度較小時，裂縫在預置裂縫方向開裂后，受層理面影響轉向層理面方向，最終貫穿圓盤；若角度增大，伴隨有次生裂紋出現，增大至80°以上時，裂縫穿過層理面，會出現轉向和分叉；②隨著層理密度增大，裂縫不均勻擴展現象增多，圓盤端部伴隨產生的次生裂紋數量也逐漸增加，在角度較小時最明顯；在角度增大后，裂縫穿過層理面，并出現裂縫分叉。由此，進一步說明層理面對于裂縫前緣形態有著重要影響。

圖5 層理密度0.4條/mm時的結果

圖6 不同層理密度時不同層理傾角的裂縫擴展結果

圖7 不同力學性質層理條件下各層理傾角的裂縫擴展結果

取層理面線密度為0.6條/mm進行研究，分別將層理面力學性質按表1中頁巖基質的0.25、0.5、0.75倍強度折減。從圖7中可看出，層理面力學性質從層內強度的0.25倍增大至0.75倍過程中，力學性質也對裂縫前緣形態有重要的影響；當層理面與頁巖基質的力學性質越接近，層理越不易受剪切應力影響，繼而產生的次生裂紋數量減少，這可能是由于層理面的力學性質較小時，在試樣端部越易形成應力集中，促使圓盤沿層理產生次生裂紋；在角度較小時，隨層理面力學性質增強，試樣的破壞形態由順層理破壞轉變為裂縫轉向、分叉；當層理角度較大時，圓盤的裂縫前緣形態基本不發生變化，都是起裂轉向后穿過較多層理，最終貫穿圓盤。

2.2 層理面對頁巖斷裂韌性影響

斷裂韌性對裂縫擴展壓力產生影響，進而影響到裂縫的擴展及裂縫前緣的形態。經過模擬測試的頁巖力學參數計算得出其對應的斷裂韌性變化規律如圖8和圖9所示。①在層理角度相同條件下，隨著層理密度增大，模型中層理數量增多，數值模擬試樣的斷裂韌性整體呈減小趨勢；②在層理密度相同條件下，試樣的斷裂韌性隨角度變化的先減小后增大，與層理角度呈波折變化的趨勢；③在角度為20°時，裂縫在預置裂縫方向開裂后，受層理面影響轉向層理面方向，此時的斷裂韌性值最小；④在角度為90°時，層理面與加載方向垂直，裂縫在延伸的過程中，不易穿切層理面，故此時的斷裂韌性值最大。

圖8 不同層理密度條件下斷裂韌性與層理角度的關系

將層理面力學性質按基質強度的1/4、1/2、3/4倍折減，取層理面線密度為0.6條/mm，層理面角度為0～90°時為例，分析斷裂韌性的情況，如圖10所示。在層理角度相同的條件下，隨著層理面強度由層內巖石強度的1/4增大至3/4,數值模擬試樣的斷裂韌性逐漸增大；隨著層理面強度逐漸減小，巖石呈現的斷裂韌性值減小得越快；在層理面線密度為0.6條/mm時，隨著層理面強度由層內巖石強度的1/4增大至3/4，同一層理面角度間斷裂韌性值差異隨層理面角度增大而先增后減。

圖9 不同層理角度下斷裂韌性與層理密度的關系

圖10 不同層理力學性質下斷裂韌性與層理角度的關系

2.3 頁巖斷裂韌性與脆性指數的關系

頁巖的脆性也能夠顯著影響井壁的穩定性，關系到設備效率和安全鉆進時泥漿密度的選擇；是評價儲層力學特性的關鍵指標，同時還對壓裂的效果產生顯著影響。對于頁巖氣藏，頁巖的脆性越高，壓裂時裂縫形態越復雜，越容易形成縫網，達到體積改造的目的；若地層的脆性較低，在壓裂時地層不易被壓開，或者是形成對稱雙翼縫，延伸時裂縫前緣形態較為單一，不利于裂縫縫網的形成[22]。

國內外學者提出了數十種表征巖石脆性的指標，主要有基于強度、硬度、應力應變、貫入實驗、破碎程度和應變能等六大類巖石脆性評價方法，各有其局限性，有的并不適合評價頁巖[22]。采用基于巖石峰后應力跌落速率及能量比的巖石脆性特征評價方法來評價頁巖脆性，能夠充分反應巖體本身力學參數，及外在加載條件變化對于脆性破壞特征的影響[23]。通過計算峰后應力降速率BPOST和巖石破壞時所釋放彈性能與峰前儲存總能量比值BE，對脆性指數Bi進行綜合評價[24]，其計算公式如式(3)：

(3)

式中：A為修正系數，取值0.001；σp與為σr峰值應力和殘余應力，MPa；εp與εr為峰值應變和殘余應變。

采用與之前相同的數值模擬理論及方法，以龍馬溪地層巖石物理實驗數據為基礎，利用RFPA軟件對該地層進行數值模擬實驗。分別從層理的角度(以10°為間隔依次從0°變化至90°)和層理的密度(0.4條/mm，0.6條/mm)兩方面來研究其對脆性指數的影響。根據圖3中的應力-應變曲線，結合公式(3)可計算出頁巖的脆性指數，其結果如圖11所示。在此基礎上，在相同層理參數條件下，統計其斷裂韌性值與脆性指數間的關系，其關系如圖12所示。

圖11 脆性指數與層理角度的關系

從圖11中可看出，頁巖脆性指數有隨層理角度先增大再減小的趨勢；同層理密度條件下，層理角度20°時巖石的脆性指數最大，原因是此時破壞形態主要為裂縫尖端起裂后隨即沿層理面的剪切滑移破壞，破壞形態較為單一；而角度90°時破壞形態的脆性指數最小，此時的層理角度與位移加載方向垂直，相比角度20°時較易沿層理直接貫穿圓盤試樣，其破壞模式為穿切層理面的劈裂型層狀剪切破壞，故其表現的脆性較小；層理面作為片理狀構造的弱結構，隨著層理的密度減小，巖石的脆性隨之減小。

同時，從圖12中可看出，頁巖斷裂韌性與其脆性指數之間呈較好的負線性關系，即隨著巖石的脆性指數增大，表現的斷裂韌性減小。說明了頁巖的脆性越強，其可壓性越高，故其斷裂韌性越低。

圖12 不同層理密度和角度下脆性指數與斷裂韌性的關系

3 結論

以四川盆地長寧地區龍馬溪組為研究對象，在已有巖石力學特性的基礎上，用RFPA建立人字形切槽巴西圓盤試樣模型，模擬計算分析層理面對頁巖斷裂韌性的影響以及斷裂韌性測試過程中的裂縫前緣形態，得到了以下結論：

(1)隨著層理密度增大，斷裂韌性整體呈減小趨勢；斷裂韌性隨角度的變化先減小后增大，在20°時出現最小值；在層理角度相同的條件下，隨著層理面力學性質增強，斷裂韌性逐漸增大；在層理面線密度為0.6條/mm時，隨著層理面力學性質增強，斷裂韌性值差異隨層理面角度增大而先增后減。

(2)當層理角度為0°時，試樣發生劈裂破壞，裂縫未發生轉向；當層理角度較小時，裂縫前緣隨層理面發生轉向，并沿層理發生脆性剪切破壞；當角度增大時，層理面作用減弱，試樣端部出現次生裂紋；當角度增至80°以上時，裂縫在延伸過程中穿過層理會出現分叉現象。

(3)層理密度較小時，裂縫在預置裂縫方向開裂后，受層理面影響轉向層理面方向，最終貫穿圓盤，若角度增大，伴隨有次生裂紋出現，若增大至80°以上時，裂縫穿過層理面，并出現轉向和分叉；層理密度較大時，裂縫不均勻擴展現象增多，圓盤端部伴隨產生的次生裂紋數量也逐漸增加，在角度較小時最明顯，而當角度增大后，裂縫穿過層理面，并出現裂縫分叉。

(4)脆性指數與斷裂韌性之間呈現負相關關系，即隨著巖石的脆性指數增大，表現的斷裂韌性減小，并且兩者之間有較好的相關性；脆性指數有隨角度先增大再減小的趨勢。