基于離散元方法的煤層壓裂裂縫模擬研究

段永強(qiáng)，楊兆中，梅永貴，王會(huì)來，李小剛，劉德強(qiáng)，路艷軍，彭 鵬

（ 1.中國石油華北油田公司，河北任丘 062550；2.西南石油大學(xué)油氣藏地質(zhì)和開發(fā)工程國家重點(diǎn)實(shí)驗(yàn)室，四川成都 610500）

1 引言

我國煤層氣藏具有低壓、低飽和、低滲透、非均質(zhì)性強(qiáng)、變質(zhì)程度高[1-4]以及儲(chǔ)層構(gòu)造復(fù)雜、臨界解吸壓力低[5-6]的特點(diǎn)。 水力壓裂是實(shí)現(xiàn)煤層氣工業(yè)開采必不可少的增產(chǎn)措施[8-9]。 水力裂縫的延伸規(guī)律是水力壓裂機(jī)理研究的關(guān)鍵，其主要受儲(chǔ)層地質(zhì)參數(shù)、壓裂液性能和施工參數(shù)的綜合影響。對(duì)于傳統(tǒng)壓裂下裂縫延伸模擬，國內(nèi)外學(xué)者做了大量的研究工作， 主要基于物理模擬和數(shù)值模擬兩種研究手段。

在物理模擬方面，Bell[10]、Abass 等[11-12]開展了中煤階煤巖的水力壓裂實(shí)驗(yàn)，發(fā)現(xiàn)裂縫面非常粗糙，裂縫的擴(kuò)展方向平行于最大主應(yīng)力方向， 且沿面割理方向擴(kuò)展。 水力壓裂的裂縫擴(kuò)展方向受水平地應(yīng)力和天然裂縫共同控制。 Warpinski 和Anderson[13]研究了水力裂縫穿越隔層的條件和行為， 研究認(rèn)為煤巖水力裂縫能否穿越界面主要取決于垂向壓應(yīng)力大小和界面性質(zhì)。 在國內(nèi)，2002 年靳鐘銘[14]利用壓裂模擬實(shí)驗(yàn)系統(tǒng)研究了大煤樣壓裂過程，并運(yùn)用宏觀損傷力學(xué)原理，建立了煤體壓裂本構(gòu)方程和破碎塊度分布公式。 2008 年杜春志[15]對(duì)型煤和原煤試件開展了水力壓裂真三軸模擬實(shí)驗(yàn)。 實(shí)驗(yàn)結(jié)果發(fā)現(xiàn)注入排量高，水平地應(yīng)力越大，則破裂壓力越大；反之，破裂壓力較小。2009 年，黃炳香等[16]根據(jù)相似準(zhǔn)則， 開展了型煤和水泥砂漿的材料配比實(shí)驗(yàn)和水力壓裂模擬實(shí)驗(yàn)。 結(jié)果表明，煤巖水力壓裂不僅產(chǎn)生了水力主裂縫，而且由于滲透水壓力對(duì)主裂縫兩側(cè)煤體形成滲透水楔作用，導(dǎo)致裂隙和層理面開裂擴(kuò)展。2011 年藺海曉等[17]通過對(duì)型煤試件和原煤試件擬三軸水壓致裂實(shí)驗(yàn)得出：水力裂縫總是沿最大地應(yīng)力的方向擴(kuò)展，同時(shí)，泵注流量越大，起裂壓力也越大。 2012 年楊焦生等[18]通過高煤階真三軸模擬實(shí)驗(yàn)發(fā)現(xiàn)：受煤巖割理的影響，煤巖壓裂施工壓力高，裂縫形態(tài)復(fù)雜。

通過現(xiàn)有的真三軸實(shí)驗(yàn)證明了帶有割理的煤巖在壓裂時(shí)產(chǎn)生了復(fù)雜的多裂縫。 傳統(tǒng)簡單的二維雙翼對(duì)稱裂縫模擬方法已不再適用煤層壓裂裂縫模擬。 物理模擬是認(rèn)識(shí)多裂縫擴(kuò)展規(guī)律最簡單最直接的方法，但是由于現(xiàn)有的裂縫監(jiān)測(cè)手段受限、 巖樣的尺度效應(yīng)已和較高的實(shí)驗(yàn)及取樣費(fèi)用致使物理模擬不能被廣泛的推廣。 而數(shù)值模擬方法被認(rèn)為是認(rèn)識(shí)裂縫擴(kuò)展規(guī)律的經(jīng)濟(jì)而有效方法。

Olson[19]基于邊界元法對(duì)頁巖和致密砂巖提出了一個(gè)擬三維復(fù)雜裂縫網(wǎng)絡(luò)模型， 研究了多裂縫同時(shí)擴(kuò)展的裂縫形態(tài)。 Meyer 等[20-21]基于自相似原理及雙重介質(zhì)模型提出了離散化縫網(wǎng)模型， 該模型利用網(wǎng)格系統(tǒng)模擬解釋裂縫在3 個(gè)主平面上的擬三維離散化擴(kuò)展和支撐劑在縫網(wǎng)中的運(yùn)移及鋪砂方式，通過連續(xù)性原理及網(wǎng)格計(jì)算方法獲得壓裂后縫網(wǎng)幾何形態(tài)應(yīng)用離散化縫網(wǎng)模型。 2009 年，Xu 等[22]在對(duì)頁巖壓裂裂縫模擬時(shí)，將裂縫網(wǎng)絡(luò)等效成由正交的水平和垂直平面組成的橢球體，建立了線網(wǎng)模型[23-24]（ Wire-mesh），該模型能夠解釋在壓裂液同裂縫壁面和裂縫間的力學(xué)作用，可以明確的計(jì)算水力裂縫的數(shù)量和尺寸，實(shí)現(xiàn)了壓裂過程中裂縫的實(shí)時(shí)擴(kuò)展。 Dahi-Taleghani，Keshavarzi等[25-26]基于二維平面應(yīng)變彈性理論，應(yīng)用擴(kuò)展有限元法模擬研究了復(fù)雜裂縫擴(kuò)展問題。 Weng 等[27]在擬三維裂縫延伸模型的基礎(chǔ)上，考慮天然裂縫與水力誘導(dǎo)裂縫的相互作用及相鄰水力裂縫之間的應(yīng)力干擾，建立了非常規(guī)裂縫模型（ Unconventional Fracture Model，UFM）， 用以模擬體積改造中復(fù)雜裂縫網(wǎng)絡(luò)的延伸。2012 年Sesetty 等[28]模擬水力壓裂和預(yù)先存在的天然裂縫之間的相互作用， 提出了在水力壓裂過程建模這種相互作用的邊界元素為基礎(chǔ)的方法。Kresse 等[29-30]用有限差分法（ FDM）模擬了預(yù)存天然裂縫下的儲(chǔ)層的多裂縫起裂與延伸。

煤巖屬于割理分布的裂縫性儲(chǔ)層， 對(duì)于裂縫性儲(chǔ)層復(fù)雜多裂縫的模擬目前主要有離散裂縫網(wǎng)絡(luò)模型、線網(wǎng)模型、非常規(guī)裂縫模型、擴(kuò)展有限元模型以及邊界元模型等， 但是這些模擬方法以儲(chǔ)層連續(xù)介質(zhì)作為研究基礎(chǔ)，簡化了地層，使得模擬結(jié)果與實(shí)際情況存在偏離。離散單元法考慮了巖石的非連續(xù)結(jié)構(gòu)，在模擬割理分布下的水力裂縫延伸方面具有一定的優(yōu)勢(shì)。

2 離散單元法原理

離散元法的基本思想是把不連續(xù)體分離為剛性元素的集合，使各個(gè)剛性元素滿足運(yùn)動(dòng)方程，繼而求得不連續(xù)體的整體運(yùn)動(dòng)形態(tài)。在采用離散單元法求解時(shí)，每一迭代時(shí)步內(nèi)， 每個(gè)塊體單元僅對(duì)相鄰的塊體產(chǎn)生力的影響，這樣時(shí)步就需要取得足夠小，以使顯示法穩(wěn)定[31]。

離散單元法中塊體之間需要滿足平衡方程， 主要包括物理方程和運(yùn)動(dòng)方程。

2.1 物理方程

塊體之間的法向力Fn等于法向“ 疊合”un與法向剛度系數(shù)kn的乘積。 所謂的“ 疊合”是計(jì)算時(shí)假設(shè)的一個(gè)量，即：

如果兩個(gè)離散單元的邊界相互“ 疊合”，則有兩個(gè)角點(diǎn)與界面接觸。此時(shí)，剪切力為節(jié)理的剪切剛度系數(shù)與塊體的相對(duì)位移us的乘積，即：

圖1 離散單元的作用力

當(dāng)巖體受到張力分離時(shí)， 作用在巖塊表面上的法向應(yīng)力和剪切力隨即消失。對(duì)于塑性剪切破壞的情況，需要在每次迭代計(jì)算時(shí)判斷剪切力是否滿足莫爾-庫倫準(zhǔn)則，當(dāng)滿足FS≥C+Fntanφ 時(shí)，塊體之間發(fā)生剪切滑動(dòng)。

2.2 運(yùn)動(dòng)方程

根據(jù)塊體的物理方程， 可計(jì)算出作用在某一巖體上的一組力以及它們的合力和合力矩。在此基礎(chǔ)上，根據(jù)牛頓第二運(yùn)動(dòng)定律確定塊體質(zhì)心的加速度和角加速度， 進(jìn)而可以確定在時(shí)步內(nèi)的速度和角速度以及位移和轉(zhuǎn)動(dòng)量。

在x 方向的加速度為：

式中：m-巖塊的質(zhì)量，F(xiàn)x-x 方向的合力。

通過差分計(jì)算， 可以得到巖塊質(zhì)心在x 方向的速度和位移：

式中：Δt-計(jì)算時(shí)步；t0-初始時(shí)間；t1=t0+Δt。 對(duì)y 方向的加速度、速度以及位移同樣由上式可以計(jì)算得出。

采用離散單元法模擬裂縫模擬的一般過程（ 見圖2）。

3 地層壓裂模擬

3.1 裂縫延伸幾何模型

簡化起見， 建立了300 m×300 m×80 m 的地質(zhì)模型，水平井鉆進(jìn)方向?yàn)樽钚∷街鲬?yīng)力方向，煤層氣井垂深為770 m，煤層厚度為6 m，其幾何模型（ 見圖3）。

3.2 參數(shù)輸入

圖2 3DEC 壓裂模擬的一般過程

表1 模擬參數(shù)

通過模擬結(jié)果發(fā)現(xiàn)， 煤層在壓裂后形成了復(fù)雜的裂縫網(wǎng)絡(luò)，同時(shí)，水裂縫的延伸方向沿著最大主應(yīng)力方向，同時(shí)割理的方向?qū)λα芽p的延伸影響明顯。當(dāng)水力裂縫遭遇割理時(shí)， 在低角度下水力裂縫會(huì)沿著割理延伸，在沿著割理延伸一定長度后，水力裂縫會(huì)繼續(xù)沿著原來方向延伸。 高角度的天然裂縫對(duì)水裂縫存在一定的遏制作用，限制水力裂縫的延伸。水力裂縫呈非對(duì)稱非均勻分布的形態(tài)。

圖3 動(dòng)態(tài)變化過程

圖4 裂縫破裂方式

多裂縫的產(chǎn)生過程主要是剪切和張開， 而天然裂縫發(fā)生破壞主要以剪切破壞為主， 在主裂縫上以張開為主。而在遠(yuǎn)端出現(xiàn)了裂縫張開點(diǎn)，這主要是在劃分網(wǎng)格時(shí)，塊體受到不均勻的力所產(chǎn)生的拉伸破壞造成的。這也意味著，在煤儲(chǔ)層由于構(gòu)造運(yùn)動(dòng)，部分區(qū)塊存在張開的天然裂縫。

在此基礎(chǔ)上， 對(duì)不同施工參數(shù)以及儲(chǔ)層條件下的煤層壓裂進(jìn)行模擬， 分析不同參數(shù)裂縫延伸的影響規(guī)律。

3.1 施工參數(shù)對(duì)裂縫形態(tài)的影響

3.1.1 排量對(duì)水力裂縫形態(tài)的影響 壓裂液排量是壓裂施工的重要參數(shù)。 現(xiàn)有的煤層壓裂大部分采用活性水壓裂液， 由于活性水壓裂液的黏度較低攜帶石英砂困難，故現(xiàn)場(chǎng)大部分采用大排量的注入。排量不僅影響著壓裂縫高、對(duì)壓裂后裂縫形態(tài)也有重要的影響，為此根據(jù)現(xiàn)場(chǎng)常見的壓裂排量分析體積縫網(wǎng)的形態(tài)對(duì)壓裂設(shè)計(jì)有重要的意義。對(duì)不同的排量進(jìn)行裂縫模擬，模擬結(jié)果（ 見圖5）。

圖5 不同排量模擬結(jié)果

圖6 不同壓裂液黏度模擬結(jié)果

根據(jù)模擬結(jié)果得出：（ 1）隨著排量的增加，注入壓力不斷增加，裂縫寬度也不斷增加，壓力的波及范圍和裂縫擴(kuò)展范圍也增大。意味著，隨著排量的增加壓裂改造范圍增大；（ 2）當(dāng)水力裂縫遭遇割理時(shí)，水力裂縫容易被割理抑制造成裂縫發(fā)生轉(zhuǎn)向。 而轉(zhuǎn)向壓力高于裂縫延伸壓力下，致使水力裂縫在另一側(cè)延伸，形成非對(duì)稱裂縫。

3.1.2 壓裂液黏度對(duì)縫網(wǎng)形態(tài)的影響 壓裂液黏度影響水力裂縫形態(tài)的又一重要因素。 較高的壓裂液黏度不僅有利于攜砂而且在控制縫高也有一定優(yōu)勢(shì)， 但高黏壓裂液不利于復(fù)雜縫網(wǎng)的形成。為此，分析不同的壓裂液黏度對(duì)煤巖體積縫網(wǎng)的影響同樣意義重大。 對(duì)采用不同壓裂液黏度的壓裂模擬結(jié)果（ 見圖6）。

通過對(duì)不同壓裂液黏度的壓裂模擬發(fā)現(xiàn)：（ 1）隨著壓裂液黏度的增加，注入壓力、裂縫寬度明顯增加；（ 2）隨著壓裂液黏度的增加，裂縫延伸的長度增加，分支裂縫張開的長度和張開的數(shù)量明顯減?。唬?3）黏度的增加壓力波及的長度得到了增加， 但是這個(gè)波及面積也不斷減小。 水力裂縫更加接近雙翼對(duì)稱的長裂縫。

3.2 不同的水平應(yīng)力差對(duì)水力裂縫形態(tài)的影響

地應(yīng)力是形成復(fù)雜縫網(wǎng)的地質(zhì)因素， 最大最小水平主應(yīng)力的差異在一定程度上決定著割理張開的難易程度，進(jìn)而直接影響著體積縫網(wǎng)的形成。 為此，對(duì)不同的差應(yīng)力系數(shù)下的儲(chǔ)層進(jìn)行3DEC 模擬， 分析差應(yīng)力系數(shù)對(duì)復(fù)雜縫網(wǎng)的影響規(guī)律，模擬結(jié)果（ 見圖7）。

對(duì)不同的差應(yīng)力系數(shù)下的煤巖儲(chǔ)層壓裂模擬得出：根據(jù)模擬結(jié)果發(fā)現(xiàn)，隨著差應(yīng)力系數(shù)的增加，水力裂縫在最大水平主應(yīng)力方向延伸長度增加， 在最小水平主應(yīng)力方向出現(xiàn)減小的趨勢(shì)。 同時(shí)注入壓力也隨著應(yīng)力比的增加而增加。說明在高應(yīng)力差異系數(shù)下，不利于割理的張開與多裂縫的產(chǎn)生， 而更有利于形成對(duì)稱雙翼對(duì)稱裂縫。

圖7 不同應(yīng)力比模擬結(jié)果

4 結(jié)論與建議

通過采用離散元方法對(duì)煤巖裂縫延伸模擬的研究得出以下主要結(jié)論：

（ 1）對(duì)于煤巖壓裂裂縫擴(kuò)展的研究方法主要有物理模擬和數(shù)值模擬，物理模擬由于本身成本高、裂縫觀察受限等原因難以廣泛推廣。 數(shù)值模擬是裂縫擴(kuò)展模擬的主要方法，而常規(guī)的有限元、有限差分、邊界元等方法在處理非連續(xù)介質(zhì)時(shí)存在不足， 離散單元法是處理煤巖割理發(fā)育的非連續(xù)地層具有明顯優(yōu)勢(shì)， 使得模擬結(jié)果更加符合實(shí)際情況。

（ 2）通過對(duì)不同壓裂施工參數(shù)模擬結(jié)果可知，壓裂液的黏度與壓裂液的排量對(duì)多裂縫的形成有著重要的影響，隨著壓裂液排量的增加，煤巖中割理更容易張開形成復(fù)雜縫網(wǎng)，壓裂液黏度的增加不利于割理的張開，但在一定程度上能夠增加壓裂施工凈壓力， 增加主裂縫的寬度。

（ 3）通過對(duì)不同的最大最小水平主應(yīng)力比模擬可知，應(yīng)力比對(duì)割理的張開也存在著較大影響，隨著應(yīng)力比的增加，分支割理不易張開，壓裂時(shí)主要是對(duì)稱雙翼裂縫，在地應(yīng)力比下分支裂縫的數(shù)量較多，并且延伸長度較大。

［ 1］ 王紅巖，劉洪林，趙慶波，等. 煤層氣富集成藏規(guī)律［ M］.北京：石油工業(yè)出版社，2005.

［ 2］ 鄧澤，康永尚，劉洪林，等. 開發(fā)過程中煤儲(chǔ)層滲透率動(dòng)態(tài)變化特征［ J］.煤炭學(xué)報(bào)，2009，34（ 7）：947-951.

［ 3］ 張培河. 影響我國煤層氣可采性主要儲(chǔ)層參數(shù)特征［ J］.天然氣地球科學(xué)，2008，18（ 6）：880-884.

［ 4］ 葉建平，史保生，張春才. 中國煤儲(chǔ)層滲透性及其主要影響因素［ J］.煤炭學(xué)報(bào)，1999，24（ 2）：118-122.

［ 5］ 馮三利，葉建平.中國煤層氣勘探開發(fā)技術(shù)研究進(jìn)展［ J］.中國煤田地質(zhì)，2003，15（ 6）：19-23.

［ 6］ 劉貽軍，等. 中國煤層氣儲(chǔ)層特征及開發(fā)技術(shù)探討［ J］.天然氣工業(yè)，2004，24（ 1）：68-71.

［ 7］ Jeffrey R G，Vlahovic W，Doyle R P，et al. Propped Fracture Geometry of Three Hydraulic Fractures in Sydney Basin Coal Seams［ C］.Paper 50061 presented at the Asia Pacific Oil &Gas Conference and Exhibtion，Perth，Australia，1998.

［ 8］ 李道品.低滲透油田開發(fā)［ M］.北京：石油工業(yè)出版社，1999：2-16.

［ 9］ Bell G J，Jones A H，Morales R H，et al.Coal Seam Hydraulic Fracture Propagation on a Laboratory Scale［ C］.Proceedings of the 1989 Coalbed Methane Symposium，The University of Alabama/Tuscaloosa，1989：417-425.

［ 10］ Abass H H，Van Domelen M L，El Rabaa W M. Experimental observations of hydraulic fracture propagation through coal blocks［ C］.SPE Eastern Regional Meeting，31 October-2 November，Columbus.Ohio，USA，1990.

［ 11］ Abass H H，Hedayati S，Kim C M. Mathematical and experimental simulation of hydraulic fracturing in shallow coal seams［ J］.SPE Eastern Regional Meeting，22-25 October，Lexington，Kentucky，USA，1991.

［ 12］ Beugelsdijk L J L，De Pater C J，Sato K.Experimental hydraulic fracture propagation in a multi-fractured medium［ C］.SPE Asia Pacific Conference on Integrated Modelling for Asset Management.Society of Petroleum Engineers，2000.

［ 13］ 靳鐘銘，弓培林，靳文學(xué).煤體壓裂特性研究［ J］.巖石力學(xué)與工程學(xué)報(bào)，2002，21（ 1）：70-72.

［ 14］ 杜春志.煤層水壓致裂理論與應(yīng)用研究［ D］.北京： 中國礦業(yè)大學(xué)，2008.

［ 15］ 黃炳香. 煤巖體水力致裂弱化的理論與應(yīng)用研究［ D］.北京：中國礦業(yè)大學(xué)，2009.

［ 16］ 藺海曉，杜春志.煤巖擬三軸水力壓裂實(shí)驗(yàn)研究［ J］.煤炭學(xué)報(bào)，2011，36（ 11）：1801-1805.

［ 17］ 楊焦生，王一兵，李安啟，等.煤巖水力裂縫擴(kuò)展規(guī)律試驗(yàn)研究［ J］.煤炭學(xué)報(bào)，2012，37（ 1）：73-77.

［ 18］ Olson J E，Dahi -Taleghani A. Modeling simultaneous growth of multiple hydraulic fractures and their interaction with natural fractures［ C］.SPE Hydraulic Fracturing Technology Conference，19-21 January，The Woodlands，Texas，USA，2009.

［ 19］ B R Meyer，L W. Bazan，R H Jacot，et al. Optimization of Multiple Transverse Hydraulic Fractures in Horizontal Wellbores［ C］.SPE Unconventional Gas Conference，23-25 February，Pittsburgh，Pennsylvania，USA，2010.

［ 20］ B R Meyer，L W.Bazan.A Discrete Fracture Network Model for Hydraulically Induced Fractures - Theory，Parametric and Case Studies［ C］.SPE Hydraulic Fracturing Technology Conference，24 -26 January，The Woodlands，Texas，USA，2011.

［ 21］ Xu W，Thiercelin M J，Walton I C. Characterization of hydraulically-induced shale fracture network using an analytical/semi-analytical model［ C］.SPE Annual Technical Conference and Exhibition. Society of Petroleum Engineers，2009.

［ 22］ Wenyue Xu，Marc Thiercelin，et al. Wiremesh： A Novel Shale Fracturing Simulator［ C］.Paper SPE 132218 presented at CPS/SPE Intemational Oil & Gas Conference and Exhibition.Beijing China，2010.

［ 23］ 張士誠，牟松茹，崔勇. 頁巖氣壓裂數(shù)值模型分析［ J］.天然氣工藝，2011，31（ 12）：81-84.

［ 24］ Dahi-Taleghani A，Olson J E. Numerical modeling of multistranded-hydraulic-fracture propagation： Accounting for the interaction between induced and natural fractures［ J］.SPE J，2011，16（ 3）：575-581.

［ 25］ Wang S，Li Y，Liu H，et al. Hydraulic Fracture Propagation in Unconventional Reservoirs： The Role of Bedding Plane［ J］.SL： Structural Longevity，2012，8（ 3）：173-191.

［ 26］ Weng X，Kresse O，Cohen C，et al. Modeling of hydraulicfracture-network propagation in a naturally fractured formation［ J］.SPE Production and Operations，2011，26（ 4）：362-368.

［ 27］ Sesetty V，Ghassemi A. Modeling and analysis of stimulation for fracture network generation［ C］.Proceedings of the 37th workshop on geothermal reservoir engineering，Stanford University，Stanford，CA.2012.

［ 28］ Kresse O，Cohen C，Weng X，et al. Numerical modeling of hydraulic fracturing in naturally fractured formations［ C］.45th US Rock Mechanics/Geomechanics Symposium.American Rock Mechanics Association，2011：26-29.

［ 29］ Kresse O，Weng X，Gu H，et al. Numerical modeling of hydraulic fractures interaction in complex naturally fractured formations ［ J］.Rock mechanics and rock engineering，2013，46（ 3）：555-568.

［ 30］ 王泳嘉，等.離散單元法及其在巖土力學(xué)中的應(yīng)用［ M］.東北大學(xué)出版社，1991.