層理特征對油頁巖水力壓裂裂縫擴展規(guī)律影響的數(shù)值模擬研究

李彥偉，朱超凡，曾壹堅，水浩澈，范存翰，郭 威

(1.吉林大學 建設工程學院，吉林 長春 130026；2.油頁巖地下原位轉(zhuǎn)化與鉆采技術(shù)國家地方聯(lián)合工程實驗室，吉林 長春 130026；3.頁巖油氣資源勘探開發(fā)省部共建協(xié)同創(chuàng)新中心，吉林 長春 130026；4.自然資源部復雜條件鉆采技術(shù)重點實驗室，吉林 長春 130026)

油頁巖是一種高灰分的固體可燃有機質(zhì)沉積巖，低溫干餾可獲得油頁巖油，我國油頁巖地質(zhì)資源量巨大，開采潛力巨大，是保障我國能源安全的重要戰(zhàn)略資源[1]。油頁巖地下原位轉(zhuǎn)化技術(shù)是實現(xiàn)油頁巖綠色無污染，高效開發(fā)，商業(yè)化發(fā)展的核心技術(shù)[2]。但是，油頁巖的低滲透性、低導熱性和低含油率，給其地下原位轉(zhuǎn)化開采帶來了巨大的困難和挑戰(zhàn)[3]。吉林大學通過水力壓裂技術(shù)和高效加熱技術(shù)在吉林農(nóng)安和扶余開展油頁巖原位轉(zhuǎn)化先導試驗工程[4]。水力壓裂是油頁巖實現(xiàn)高效開發(fā)核心技術(shù)，油頁巖水力壓裂技術(shù)不同于頁巖砂巖等儲層壓裂，油頁巖水力壓裂旨在形成一種適合于油頁巖原位轉(zhuǎn)化的“高效加熱裂解縫網(wǎng)”，為油頁巖提供高效加熱條件以及滲流傳熱通道。

針對層理面對水力壓裂裂縫擴展規(guī)律影響問題，國內(nèi)外眾多學者主要從層理面對水力裂縫擴展規(guī)律以及數(shù)值模擬方法進行研究[5-7]。水力壓裂裂縫擴展規(guī)律由原地應力狀態(tài)和層理面屬性共同決定，層理方向是水力壓裂裂紋擴展方向的主控因素[8]。層理面的弱膠結(jié)作用，在平行層理方向阻止裂紋擴展能力較弱，在垂直層理方向阻止裂紋擴展能力較強[9]。較弱膠結(jié)強度、較大結(jié)構(gòu)面張開度、較低水地應力差、較低排量能夠激活層理面，易形成縫網(wǎng)[10-11]。層理間距、井筒方位角、壓裂參數(shù)等也影響對水力裂縫擴展行為[12]。綜上所述，前人的研究主要集中在層理面膠結(jié)強度、層理方向、層間間距、水平應力差以及壓裂液排量和壓裂液黏度等對水力裂縫擴展影響。主要考慮頁巖砂巖的壓裂，其壓裂目標是形成體積改造，增大泄油面積。

水力壓裂數(shù)值模擬方法方面，有限元法可用于水力裂縫層內(nèi)[13]和層間的[14-16]擴展規(guī)律，認為水平地應力差和界面強度增加會增大層理面開啟難度，增加壓裂液黏度和垂向應力，提高壓裂液注入速度，有利于裂縫擴展高度延伸。黏聚力單元法主要是在有限元模型基礎上建立黏聚力模型[17]，可用于研究層理面，隨機天然裂縫等黏聚力材料屬性[18]對水力壓裂參數(shù)影響。認為影響天然裂縫網(wǎng)絡控制的主要因素是天然裂縫初始水力寬度和地應力場[13]，離散元法用于描述微觀顆粒之間的相互作用[19]，可以建立考慮力學各向異性水力壓裂裂縫擴展模型，可分析層理密度、層理強度、壓裂工程參數(shù)對水力裂縫擴展規(guī)律的影響[20]，采用前置高黏度/后置低黏度壓裂液交替注入壓裂工藝可最大限度提高壓裂裂縫復雜程度[21]。此外還有相法場[22]和邊界元[23]等方法。綜上所述，數(shù)值模擬技術(shù)主要基于有限元、邊界元、黏聚力單元法、離散元和相法場等，前人針對層理面研究采取黏聚力單元cohesive的方法主要是在模擬過程中在儲層巖石中插入一層cohesive 單元來表征層理，并不能真實反映水力裂縫沿任意路徑擴展。同時一層cohesive 單元忽略了層理面本身厚度問題。

為此，筆者以鄂爾多斯盆地旬邑地區(qū)油頁巖為研究對象，構(gòu)建應力-損傷-滲流的水力壓裂裂縫擴展模型，基于全局FEM-CZM(Finite Element Method-Cohesive Zone Method)方法開展數(shù)值模擬實驗研究，探究層理厚度、層理間距、地應力場對水力壓裂裂縫擴展的影響，以期為油頁巖水力壓裂工藝設計提供依據(jù)。

1 研究區(qū)選取與地質(zhì)特征

鄂爾多斯盆地南部旬邑地區(qū)油頁巖層理發(fā)育，地應力場分布不規(guī)律，油頁巖整體呈NE-SW 向分布，東西長約74 km，南北寬約27.5 km，面積約為1 496 km2[24]。鄂爾多斯盆地油頁巖地層為上三疊統(tǒng)延長組和侏羅系延安組、安定組，是鄂爾多斯盆地內(nèi)陸湖盆形成后的第一套自生自儲油頁巖層系，也是油頁巖發(fā)育的主要層系[25-26]。旬邑油頁巖主要賦存于三疊系延長組第二巖性段下部，埋深0～500 m，厚度5.74～12.49 m，平均厚度10.07 m[27]?；谑覂?nèi)實驗，油頁巖礦物組成主要為陸源碎屑礦物、黏土礦物、少量碳酸鹽礦物和其他礦物。陸源碎屑礦物主要為石英和長石，碳酸鹽礦物主要為方解石，其他礦物主要為黃鐵礦，見表1。目標地區(qū)的油頁巖含油率為6.0%，屬于中含油率油頁巖，成熟度為0.66 %，滲透率為0.013×10-3μm2，孔隙率為2.149 %，目標區(qū)域油頁巖層理發(fā)育，層理厚度平均在0～3 mm，層間距為0～25 mm，如圖1 所示。

表1 旬邑油頁巖地質(zhì)特征和礦物組成Table 1 Geological characteristics and mineral composition of oil shale in the Xunyi area

2 黏聚力單元法模擬水力壓裂理論

黏聚力單元法可模擬水力壓裂過程，水力壓裂過程中，壓裂液進入膠結(jié)面如果超過巖石的膠結(jié)強度，膠結(jié)面逐漸被打開，膠結(jié)面打開過程遵循牽引-分離損傷準則和損傷演化理論。在此過程中膠結(jié)面中一部分壓裂液沿著膠結(jié)面進行流動，一部分濾失到地層中。全局黏聚力單元法主要是基于ABAQUS 有限元軟件，通過二次開發(fā)技術(shù)在每個基礎單元外邊界批量嵌入黏聚力單元，其優(yōu)點是可實現(xiàn)水力裂縫在網(wǎng)格邊界上沿任意路徑擴展。同時多層黏聚力單元在平面內(nèi)所構(gòu)成的一個封閉面體系，可定量表征層理厚度屬性。因此，本文采取全局黏聚力單元FEM-CZM[28-31]數(shù)值模擬法。

2.1 線彈性牽引-分離準則

研究過程中，假定cohesive 黏結(jié)單元開始出現(xiàn)損傷時和損傷演化前，其承受的應力和相應的應變滿足線彈性關(guān)系，如下式所示，所對應的牽引-分離過程曲線如圖2 所示。

圖2 典型的牽引-分離曲線Fig.2 A typical traction-separation curve

式中：σ為應力矢量，MPa；σn、σs、σt分別為法向應力、第一切向應力、第二切向應力，MPa；K為彈性剛度矩陣；GPa；?為應變矢量；?n、?s、?t分別為法向應變、第一切向應變、第二切向應變；E為彈性模量，GPa。

2.2 裂縫初始損傷準則

目前，軟件ABAQUS 常用的水力壓裂損傷準則有4 種：最大名義應力損傷準則、最大名義應變損傷準則、二次應力損傷準則和二次應變損傷準則，本模型采用的是最大名義應力損傷準則：

2.3 裂縫損傷演化準則

cohesive 黏結(jié)單元滿足初始損傷準則后，損傷演化準則是描述材料剛度退化速率。損傷變量D代表材料整體損傷，并考慮所有演化方式的綜合作用。初始損傷值為0，此時材料沒有發(fā)生損傷。定義了損傷演化模式之后形成損傷，隨著荷載的持續(xù)進行，損傷變量D按照相應損傷模式逐漸演化到1，此時材料剛度完全退化，形成破壞。

若考慮采用線性的位移損傷演化模式，ABAQUS軟件中使用的損傷演化可簡化為：

2.4 水力壓裂縫內(nèi)流體流動理論

cohesive 單元縫內(nèi)流體的流動主要分為沿著裂縫方向的徑向流動以及法向的流動，如圖3 所示。

圖3 Cohesive 黏結(jié)單元流體流動Fig.3 Schematic diagram showing fluid flow in cohesive zones

流體沿著cohesive 單元徑向流動的計算公式為：

流體沿著cohesive 單元法向流動的計算公式為：

式中：q為縫內(nèi)流量，m3/s；qt為流入頂部表面流量，m3/s；qb為流入底部表面流量，m3/s；d為cohesive 黏結(jié)單元的寬度，m；μ為壓裂液黏度，Pa·s；Δp為cohesive 黏結(jié)單元壓力梯度，MPa；ct、cb為cohesive 單元上下表面的濾失系數(shù)，m2/s；pf為cohesive 單元內(nèi)壓裂流體壓力，MPa；pt為cohesive 單元上表面孔隙壓力，MPa；pb為cohesive 單元下表面孔隙壓力，MPa。

3 裂縫擴展數(shù)值模型建立與驗證

3.1 數(shù)值模型建立

模型建立過程主要是基于FEM-CZM 方法，建立全局黏聚力單元，通過自主研究的網(wǎng)格節(jié)點在連續(xù)介質(zhì)網(wǎng)格內(nèi)添加一層cohesive 單元層，賦予cohesive 層理-基質(zhì)2 種不同面屬性。模型尺寸大小為0.1 m×0.1 m，在中心位置模擬注水孔，其半徑為2 mm，層理的傾角為0°，層理厚度的大小設置為0、1、2、3 mm。層理間距的大小設置為6、10、16、20 mm。如圖4 所示，其中X軸平行于σh(最小水平主應力)，Y軸平行于σH(最大水平主應力)，Z軸平行于σV(垂向主應力)，地應力場數(shù)據(jù)參考王永煒等[32]對鄂爾多斯盆地陸相中生界頁巖氣儲層地應力測試結(jié)果，應力場(σh、σH、σV)大小分別為(8、9、9 MPa)、(8、9、10 MPa)、(8、9、12 MPa)、(8、9、13 MPa)，模型整體采用有效應力和超靜水壓力原理，所對應的孔壓邊界為0，其中模型在X方向與Z方向采用固定約束自由度，約束自由度為0，XOZ平面內(nèi)的自由度不被約束。

圖4 油頁巖水力壓裂數(shù)值模型Fig.4 A numerical model for hydraulic fracturing of oil shale

數(shù)值模型巖石基質(zhì)和層理面的力學參數(shù)可參考鄂爾多斯盆地南部馬泉油頁巖巖石力學參數(shù)[33]，相關(guān)結(jié)果見表2。

表2 水力壓裂模型相關(guān)參數(shù)[33]Table 2 Related parameters of the hydraulic fracturing model[33]

3.2 數(shù)值模型的驗證

為了驗證數(shù)值模型的可靠性，以孫可明等[34]開展的頁巖層理方位及強度對水力壓裂影響的實驗為基礎，模擬在同一條件下的裂縫擴展路徑。數(shù)值模型的尺寸與原物理模擬試驗保持一致，為400 mm×400 mm，最大水平主應力為X方向(4 MPa)，最小水平主應力Y方向(2 MPa)，層理的傾角為20°。驗證模型所用的參數(shù)為：基質(zhì)的抗拉強度0.775 MPa，層理的抗拉強度為0.344 MPa。

圖5 為數(shù)值模擬結(jié)果與試驗結(jié)果裂縫擴展路徑對比。其中圖5a 為建立的幾何模型，圖5b 為數(shù)值模擬結(jié)果，圖5c 為引文試驗結(jié)果。從圖5 可以看出數(shù)值模擬和試驗結(jié)果裂縫擴展路徑基本一致，裂縫擴展規(guī)律基本吻合。說明采用全局FEM-CZM 方法模擬水力裂縫擴展路徑的可靠性。

圖5 數(shù)值模擬與試驗結(jié)果裂縫擴展路徑對比Fig.5 Comparison of fracture propagation paths derived from numerical simulations and experimental results

4 結(jié)果與討論

4.1 模擬結(jié)果分析

4.1.1 裂縫擴展過程形態(tài)

圖6 為層理厚度2 mm、層理間距10 mm、應力場(σh、σH、σV)大小為(8、9、10 MPa)不同模擬時刻油頁巖水力裂縫在基質(zhì)中的擴展情況和層理打開程度變化結(jié)果。根據(jù)水力裂縫在基質(zhì)-層理中的相互作用關(guān)系，可將水力裂縫的延伸過程分為3 個階段。第一階段：當注入時間為0～0.125 s 時，水力裂縫在兩層理面間的基質(zhì)中擴展，隨著壓裂液的持續(xù)注入，水力裂縫沿著垂向進行擴展，裂縫的延伸長度和寬度持續(xù)增加，如圖6a 所示。第二階段：當注入時間為0.125～0.300 s時，水力裂縫與層理面接觸，穿透層理面在基質(zhì)中擴展，裂縫的延伸長度和寬度緩慢增加，如圖6b 所示。第三階段：當注入時間為0.300～10.000 s 時，水力裂縫逐漸開啟層理面，沿著層理面進行擴展，伴隨著注入壓力升高，層理面完全被打開，且穿透層理面的水力裂縫寬度遠遠小于沿著層理面打開的寬度，如圖6c 所示。

圖6 不同階段下水力壓裂裂縫擴展形態(tài)Fig.6 Morphologies of hydraulic fracture propagation at different stages

4.1.2 裂縫擴展過程參數(shù)

進一步分析水力裂縫與層理面相互作用過程中參數(shù)，研究中提取了注入點的壓力、裂縫長度和水力裂縫擴展過程中層理溝通面積[35]曲線。如圖7 所示，分析注入點壓力曲線，當注入時間為0.05 s 時，儲層基質(zhì)巖石開始發(fā)生破壞，對應注入壓力為11.11 MPa；當注入時間為0.125 s 時，水力裂縫與層理面接觸，并未產(chǎn)生相互作用，形成局部憋壓，此時對應的壓力達到最大值15.71 MPa；當注入時間為0.17 s 時，水力裂縫與層理面之間發(fā)生相互作用，開始部分穿透層理面，并伴隨著層理面開啟；當注入時間為0.3 s 時，水力裂縫在層理面中延伸，對應注入點壓力降低到6.5 MPa；當注入時間為0.63 s 時，水力裂縫完全打開層理面，隨著時間的增長，注入壓力逐漸降至0 MPa。進一步分析裂縫長度和層理溝通面積，隨著注入時間的增加，裂縫長度和層理溝通面積逐漸增加，當注入時間段為0～0.3 s 時，裂縫長度增加趨勢大于層理溝通面積的，其主要原因是此階段對應開啟層理面的水力裂縫寬度較小。當注入時間大于0.3 s，水力裂縫逐漸開啟層理面，裂縫寬度持續(xù)增加，層理溝通面積增加趨勢大于裂縫長度增加趨勢。當注入時間為0.63 s 時，水力裂縫已經(jīng)完全打開層理面，裂縫長度和層理溝通面積趨于穩(wěn)定，不再增加。

圖7 水力壓裂過程參數(shù)隨注入時間變化Fig.7 Changes in parameters of hydraulic fracturing process with injection time

4.2 影響因素分析

為進一步研究層理面對水力裂縫擴展的影響，分別設計了3 種實驗方案，研究層理厚度、層理間距、地應力場對水力裂縫擴展的影響。

4.2.1 層理厚度

基于本文建立的水力壓裂模型，在保持其他參數(shù)不變的條件下，研究在地應力場(8、9、10 MPa)，層理間距10 mm，層理厚度0、1、2 和3 mm 條件下水力裂縫擴展規(guī)律，層理厚度0 mm 指的是層理面沒有厚度屬性，但是有膠結(jié)強度等屬性。模擬結(jié)果如圖8 所示。當層理厚度為0 mm 時，水力裂縫與層理面交匯后沿著層理面擴展，裂縫寬度達到最大值；當層理厚度為1 mm 時，水力裂縫與層理面交匯后沿著層理面擴展，裂縫的長度明顯增加；當層理厚度為2 mm 時，水力裂縫與層理面交匯后在垂向上有微小的擴展，并最終沿著層理面進行擴展；當層理厚度為3 mm 時，水力裂縫與層理面交匯后沿著層理面擴展，裂縫寬度明顯降低。因此，較大層理厚度導致水力裂縫沿層理面擴展傾向更強，而沿垂直方向擴展相對較弱。

圖8 層理厚度對水力裂縫擴展的影響Fig.8 Effects of bedding plane thickness on hydraulic fracture propagation

進一步研究層理厚度對水力裂縫的破壞類型、裂縫長度、層理溝通面積的影響。如圖9 所示，裂縫的破壞類型以張拉破壞為主，當層理厚度為0 mm 時張拉破壞比例最高，這是因為此時所對應巖石材料內(nèi)部顆粒分子之間相互連接較弱，更容易達到層理的抗拉強度，發(fā)生張拉破壞；當層理厚度為1、2 和3 mm 時，隨著層理厚度增加，張拉破壞比例增加，其主要原因是地應力場對層理弱面有壓實作用，層理越厚，壓實作用就越弱，層理的黏聚力越小，發(fā)生張拉破壞的比例越高。隨著層理厚度增加，裂縫長度和層理溝通面積增加。這是因為水力裂縫均未能突破層理面，垂向方向裂縫擴展受到限制，層理越厚，層理面對水力裂縫攔截作用越強，水力裂縫在層理面中擴展空間越大，水力裂縫擴展路徑更長，更有利于開啟層理面，所對應的層理溝通面積越大。

圖9 不同層理厚度下水力裂縫擴展結(jié)果Fig.9 Hydraulic fracture propagation under different bedding plane thicknesses

4.2.2 層理間距

進一步研究地應力場(8、9、10 MPa)，層理厚度2 mm，層理間距8、10、16 和20 mm 條件下水力裂縫擴展規(guī)律。結(jié)果如圖10 所示，當層理間距為8 mm 時，水力裂縫直接貫穿層理面沿著縱向擴展，層理面未被打開；當層理間距為10、16 和20 mm 時，水力裂縫均沿著層理面進行擴展，這是因為層理間距主要影響水力裂縫到達層理面的時間，較大的層理間距意味著壓裂液需要更長的距離從起裂位置擴展到層理面，從而增加了裂縫擴展的阻力。

圖10 層理間距對水力裂縫擴展的影響Fig.10 Effects of bedding plane spacing on hydraulic fracture propagation

結(jié)合圖11 可知：對比裂縫破壞類型，層理間距為8 mm，張拉破壞比例較高，其主要原因水力裂縫直接穿透層理面在基質(zhì)中擴展，在基質(zhì)中擴展時間較長，以張拉破壞為主；對比層理間距10、16 和20 mm，隨著層理間距增加，張拉破壞的比例整體上增加，一方面原因是水力裂縫在基質(zhì)中擴展時間增加，另一方面原因是較大的層理間距意味著流體抵達層理面時能量降低。對比裂縫長度和層理溝通面積，層理間距為8 mm 時，裂縫長度最長，但并未能溝通層理面，所對應層理溝通面積為0。隨著層理間距增加，裂縫長度和層理溝通面積均增加，裂縫長度一方面受層理間距影響，另一方面層理間距增加導致水力裂縫沿著層理面有多條擴展路徑，造成裂縫長度增加，同時隨著層理間距的增加，造成的水力裂縫開啟層理面寬度也增加，所對應的層理溝通面積增加。這表明層理間距越大，水力裂縫溝通的層理面積越大，越有利于體積改造。

圖11 不同層理間距下水力裂縫擴展結(jié)果Fig.11 Hydraulic fracture propagation under different bedding plane spacings

4.2.3 地應力場

為分析地應力場對水力裂縫擴展的影響，開展了層理厚度2 mm、層理間距10 mm 條件下，應力場(σh、σH、σV)大小分別為(8、9、9 MPa)、(8、9、10 MPa)、(8、9、12 MPa)、(8、9、13 MPa)，對應垂向地應力差(Δσ=σV-σh)為1.0、2.0、4.0 和5.0 MPa 條件下水力裂縫擴展動態(tài)模擬，結(jié)果如圖12 所示。當垂向地應力差為1.0 MPa 時，水力裂縫未能穿透層理面，沿著層理面擴展；當垂向地應力差為2.0 和4.0 MPa 時，水力裂縫少部分穿透層理面，大部分沿著層理面擴展；當垂向地應力差為5.0 MPa 時，水力裂縫直接穿透層理面在基質(zhì)中擴展。這是因為地應力場決定裂縫擴展方向，垂向地應力差越大，垂向應力會對層理有壓實作用，導致層理面的黏聚力增大，進而增大水力裂縫沿著層理面擴展阻力，導致裂縫更容易穿透層理面擴展。

圖12 垂向地應力差對水力裂縫擴展影響Fig.12 Effects of vertical in-situ stress difference on hydraulic fracture propagation

結(jié)合圖13 可知，當垂向地應力差為1.0、2.0 和4.0 MPa 時，隨著垂向地應力差增加，張拉破壞比例下降，其主要原因是以上3 種情況水力裂縫均未穿透層理面，垂向地應力差越大，對層理面壓實作用越強，盡管壓裂液可以進入層理面，但是在層理面中擴展存在較大阻力，因此發(fā)生張拉破壞比例下降，剪切破壞比例增大。垂向地應力差為5.0 MPa 時，水力裂縫直接穿透層理面，發(fā)生張拉破壞。對比裂縫長度和層理溝通面積，當垂向地應力差為1.0、2.0 和4.0 MPa 時，隨著垂向地應力差的增加，裂縫長度和層理溝通面積增加，這是因為垂向地應力差越大，水力裂縫在層理面中擴展受阻，所造成裂縫擴展形態(tài)多發(fā)生彎曲、分支的情況，整體上呈現(xiàn)增加的趨勢；當垂向地應力差為5.0 MPa時，水力裂縫穿透層理面，所對應裂縫長度最小。

圖13 不同垂向地應力差下水力裂縫擴展結(jié)果Fig.13 Hydraulic fracture propagation under different vertical in-situ stress differences

5 結(jié)論

a.層理厚度、層理間距越大，裂縫破壞類型以張拉破壞為主，對應裂縫長度和層理溝通面積越大；垂向地應力差較大時，裂縫破壞類型以張拉破壞為主，對應裂縫長度和層理溝通面積較小。

b.基于上述認識，建議在壓裂施工選址方面選擇層理厚度較大，層理間距較大，垂向地應力場較小的區(qū)域，更有利于水力裂縫溝通層理面，增大改造體積形成復雜縫網(wǎng)。

c.本次實驗揭示了考慮層理厚度情況下油頁巖水力壓裂裂縫擴展規(guī)律，下一步開展考慮壓裂縫網(wǎng)下的油頁巖原位轉(zhuǎn)化開發(fā)效果模擬，為油頁巖原位轉(zhuǎn)化工藝提供指導。