非常規油氣藏體積壓裂數值模擬新進展

唐子春,王 朝,張子珂,胡 蝶,雷征東

(1.中國石油勘探開發研究院，北京 100083；2.中國石油大學(北京))

非常規油氣藏體積壓裂數值模擬新進展

唐子春1,2,王 朝1,2,張子珂1,2,胡 蝶1,2,雷征東1

(1.中國石油勘探開發研究院，北京 100083；2.中國石油大學(北京))

頁巖氣藏的高效開發與體積壓裂的合理設計密切相關，在調研大量相關文獻的基礎上，綜述了近年來國外非常規油氣藏水平井體積壓裂數值模擬技術。從天然裂縫與人工縫的相互作用、應力影以及支撐劑模擬等方面敘述了壓裂力學模型的發展歷程，介紹了非常規壓裂模型的原理，為該技術的深入研究提供參考;闡述了連續介質模型和離散裂縫模型在非常規油氣藏體積壓裂數值模擬中的應用，分析了其原理、數模思路以及適用性;總結了國外研究的主要方向，提出了針對長慶油田致密儲層體積壓裂數值模擬的發展方向。

長慶油田；非常規油氣藏；數值模擬；體積壓裂；離散裂縫模型

北美頁巖氣的勘探與開發至今已有百年的歷史，而真正進入商業化生產僅有20多年。水平井開發、體積壓裂等技術的突破促進了美國的頁巖氣革命，使非常規油氣在世界能源供給中占有越來越大的比例。從非常規油氣藏獲得工業性的經濟產量有賴于水平井近井地帶的改造體積，而在體積壓裂過程中產生的水力裂縫會與天然裂縫相互交錯，在水平井近井地帶形成復雜的裂縫網絡。因此，如何合理、準確地描述復雜縫網系統成為非常規油氣藏中數值模擬的關鍵。近年來，國外非常規油氣藏體積壓裂的數值模擬技術有了飛躍式的發展，而目前國內還沒有對此最新研究成果進行系統論述。本文針對該類問題，大量調研了近年來的相關研究，總結了雙重介質模型、離散裂縫模型的數值模擬方法，提出了相應的數模思路以及對未來技術的展望。

1 復雜縫網模型

不同于PKN模型、擬三維模型等常規壓裂模型，復雜縫網模型需要考慮很多因素，例如天然裂縫和水力裂縫的相互作用、應力影的影響等。Olson等人提出了一個縫網模型來預測水力裂縫的擴展，該模型雖然考慮了人工縫與天然縫的相互作用，但不能模擬壓裂液流動和支撐劑輸送；Rogers等人[1]提出了一個離散裂縫模型，簡化了天然裂縫與人工縫的關系，但沒有考慮壓裂縫之間的作用；Nagel等人[2,3]提出了一個擬三維縫網模擬器，可模擬天然裂縫與人工縫之間的相互作用，但該模型計算量大，不具有工程適用性；而最為經典的是Xu等人[4,5]提出的線網模型和Weng等人[6]提出的非常規壓裂模型。線網模型(Wire-Mesh)將復雜裂縫網絡簡化為一個橢球體，該橢球體由兩組正交、等距分布的垂直裂縫面組成。線網模型應用了巖石力學理論，計算了縫網隨著壓裂施工的擴展，考慮了天然裂縫與水力裂縫的力學作用，模擬了支撐劑在縫網中的輸送。非常規壓裂模型(UFM)是Kresse.O[7]于2011年提出的新模型，模型主要考慮以下幾個方面：

(1)基本的控制方程。①流動方程：一維泊肅葉圓管層流，假定流體為冪律流體；②質量守恒；③考慮液體濾失效應：限制條件是射孔點流速之和等于總排量。

(2)高度生長模型。主要基于擬三維壓裂模型計算，考慮流體壓力、地應力、斷裂韌性、彈性模量等因素，在裂縫尖端的上部、下部分別計算斷裂韌性，判斷裂縫高度延伸的可能方位。

(3)裂縫相互作用模型[8]。該方法是最早實現水力裂縫與天然裂縫相互作用的模型之一，模型本身基于線彈性斷裂力學，沒有考慮液體影響。水力裂縫默認垂直于最小水平主應力方向起裂，單一裂縫時，非常規壓裂模型退化為擬三維壓裂模型。Warpinski[9]提出，當水力裂縫遇到天然裂縫時，會出現穿過和滑移兩種情況，而Blaton[10]的觀點是會出現張開、穿過、捕捉三種情況，但二者本質相同，分別由圖1b和c所示。

圖1 水力裂縫遇到天然裂縫的可能情況

(4)應力影模型。模型使用邊界元(BEM)方法，該方法考慮了三維應力修正項，模擬結果與CSIRO十分接近[11]。式(1)和式(2)分別代表第i個裂縫單元在受到應力影作用后的法向與切向應力值。

(1)

(2)

式中n，s分別代表張開、剪切應力；Dn，Ds代表張開、剪切位移；N代表裂縫體內總單元數；Cij是二維平面應變影響因子，該因子是用于作用于計算來自相鄰裂縫單元的額外影響應力值，具體表達式參見Crouch[12](1976)；Aij是三維修正系數，修正因子Cij。當距離增加時，三維條件下有限的裂縫高度會削弱裂縫單元間的影響，需要引入該系數進行修正。具體表達式參見Olson[13](2004)。

非常規壓裂模型在考慮應力干擾時有兩大假定：①裂縫寬度等于法向位移不連續度；②裂縫表面剪切應力等于0。基于這兩個假設才能計算切向位移不連續度及法向的應力干擾值。應力影模型對裂縫傳遞過程的影響體現在兩個方面：一是影響井場所有裂縫單元的初始地應力大小與方向，改變了裂縫的壓力與寬度；二是影響了裂縫尖端的應力場，可能導致裂縫沿其他方向擴展。多條裂縫的相互干擾模擬結果：當裂縫間距低于10 m時，水力裂縫之間表示出強烈的排斥現象，而間距40 m以上排斥現象變得比較微弱。

(5)支撐劑傳輸模型。為提高計算效率，水平方向(沿縫長)設定為一維模型，垂向(沿縫高)設定為三層模型，從底部到上部分別為砂堤、混砂液、無砂液，同時考慮支撐劑的沉降與侵蝕作用。Cohen等人[14]對非常規壓裂模型的各個參量進行了1800多次的擬合試驗，結果表明，支撐劑尺寸為40～70目時，可流動面積比30～50目大了3倍以上，分布也更為均勻。模型計算過程考慮了各個模型的相互作用，對每個時步，該模型基本計算流程如圖2所示。

圖2 UFM模型計算流程

非常規壓裂模型最大的優勢在于，除了考慮常規壓裂軟件的基本功能外，還考慮了多縫間應力干擾，包括單井多段間應力干擾與井間應力干擾，在現有的壓裂模型中有獨特的優勢。該模型還可以模擬與天然裂縫的相互作用，從而進行多井壓裂等復雜的壓裂設計，為壓裂領域較為前沿的壓裂設計模型。但是，模型本身的天然裂縫與水力裂縫的擴展準則是基于線彈性斷裂力學理論，沒有考慮流體的影響，可能會造成較大的偏差；地層天然裂縫識別困難，給建立天然裂縫網絡帶來一定的不確定性。目前，模型還處于進一步完善之中。

2 體積壓裂數值模擬

2.1 連續介質模型

連續介質模型包括雙孔模型、雙孔雙滲模型以及多重孔隙介質模型。傳統的雙孔模型假設基質到裂縫的流動為擬穩定流，而在非常規油氣藏中，流動并不是瞬間發生，需要多級基質網格來模擬致密儲層中的微孔隙、裂縫和有機質區域，獲取過渡流特征。Pruess[15]提出的MINC模型解決了這個問題，該模型將基質單元剖分成一系列次級單元，這些次級單元形成獨立的一維流動系統，并與對應的裂縫單元相連，后續研究又將模型中次級單元與裂縫表面的距離改進為對數變化。Zhang[16]使用ECLIPSE中的多重孔隙結構來描述體積壓裂區域基質和裂縫的關系，提出了一個綜合模擬方法，研究油藏和水力壓裂參數對頁巖氣藏生產動態的影響，模擬結果表明，該方法很好地獲取了基質到裂縫之間的過渡流動特征。

對于水力裂縫的模擬，Kalantari Dahaghi[17]研究認為，對數間距局部網格加密對水平井體積壓裂區域有較好的模擬效果。在此基礎上，Cipolla和CMG公司的Rubin[18-19]基于雙重介質模型提出了一種叫做“雙滲模型-對數間距變化-局部網格加密(DK-LS-LGR)”的方法，結合微地震監測結果，對致密儲層中的復雜縫網進行整體建模。該方法在體積改造區域內部采用對數間距變化、局部網格加密的雙滲模型，在外部采用普通的雙滲網格。復雜縫網的建模首先通過微地震監測數據篩選體積改造區域的裂縫網絡，之后，通過對數間距局部網格加密模擬水平井的體積改造區域。結果表明，該方法歷史擬合效果好，可以較好反映頁巖的生產特征。Du等人[20-21]提出了一個綜合的工作流對Barnett頁巖氣藏進行模擬。在地質模型的建立過程中，天然裂縫建模基于井筒成像、測井和微地震監測的解釋結果，水力裂縫建模依賴于壓裂施工參數、支撐劑的分布以及微地震監測的結果，采用Oda方法將油藏地質模型和離散裂縫模型粗化為雙孔模型。現場應用結果表明，這個雙孔模型工作流可以較好地模擬非常規油氣藏的生產，使人們對致密儲層的生產有了更深的認識，同時，也為非常規油氣藏一體化模擬工作流的快速發展奠定了基礎。

2.2 離散裂縫模型

離散裂縫模型(DFM)是指用一系列指定網格代表裂縫，使裂縫離散分布在模型中，利用網格之間的傳導率控制基質與裂縫之間的流動。根據求解方法的不同，可將該方法分為非結構化離散裂縫模型(UDFM)和嵌入式離散列分模型(EDFM)[22-23]。

2.2.1 非結構化離散裂縫模型

常規的離散裂縫模型有賴于非結構化網格以及局部網格加密實現復雜縫網的精細描述，Karimi-Fard[24]、Fu[25]、Marcondes[26]等人應用有限元方法建立了離散裂縫模型，采用非結構化網格對裂縫進行了顯性剖分。

近年來，通過將地震、測井、微地震監測以及巖心分析結合到復雜縫網模型，進行不同方式的網格化，實現了離散裂縫模型到油藏數值模型的耦合。Cipolla等人[27-29]提出了地震-數模一體化工作流,該工作流包括以下步驟：①通過地震、測井解釋以及巖心分析手段，提供必要的巖石力學性質、地應力分布、斷層以及天然裂縫的分布，建立離散裂縫模型和地質力學模型，并運用微地震監測進行約束；②基于已建立的地質模型，使用非常規壓裂模型(UFM)輸入壓裂施工參數和巖石力學性質，通過計算耦合人工裂縫和天然裂縫，生成復雜縫網系統；③通過自動非結構網格化算法進行局部網格加密，生成油藏數值模型，由下一代數值模擬器INTERSECT進行計算。圖3為頁巖氣藏水平井衰竭式開發10年以后壓力分布模擬結果，從圖中可以看出，該工作流可以精細地描述水平井體積壓裂區域的開發特征。

圖3 體積壓裂水平井衰竭開采10年后的壓力分布

在此基礎上，Cipolla等人[30]提出了一個完井設計工具，將其耦合到了工作流中。該工具利用應力場、巖石力學特性和成像技術，對頁巖孔隙度、黏土含量、有機碳含量(TOC)等參數進行篩選，確定合理的壓裂段和射孔位置，進行嚴格高效的水平井分段壓裂完井設計，判別標準如表1所示。Cohen等人[31]也將工作流做了一些改進，包括壓裂段之間的應力干擾、壓裂液溫度和流變性的耦合等因素。目前，該流程已在斯倫貝謝公司的Mangrove軟件中實現，工作流如圖4所示。

地震-數模一體化工作流實現了壓裂模型與油藏數值模擬器之間的耦合，但是沒有有效的工具來模擬體積壓裂對非常規油氣藏單井產能的影響。Cohen等人[32]耦合了一個半解析的非常規產能模型(UPM)，進行多參數敏感性分析，研究壓裂施工參數設計和油藏性質與產能之間的關系。在此基礎上，Cohen等人[33]耦合了一個溫度模型，研究壓裂液體系和油藏溫度對產能的影響。

圖4 地震-數模一體化工作流

油藏性質完井性質總有機碳>=3%低應力值天然氣含量>2.83m3/t電阻率>15Ω·m干酪根含量高黏土含量<40%泥巖孔隙度>4%楊氏模量>13.8GPa泊松比<0.2氣測滲透率>10-4μm2中子測井孔隙度<35%密度孔隙度>8%

2.2.2 嵌入式離散裂縫模型

雙重介質模型建立的縫網形態相對理想，但不能正確描述改造區域的復雜情況；而常規離散裂縫模型依賴于局部網格加密或非結構化網格，可以滿足高精度的要求，但網格數目增加，造成計算量呈數量級增長。為了解決這些問題，Li和Lee(2008)[34]首次提出了嵌入式離散裂縫模型(EDFM)，該模型使用結構化網格描述基質系統，通過計算裂縫與基質網格的交會部分，顯性地引入裂縫系統控制體作為單獨的一部分，將裂縫平面嵌入到基質網格中并進行邊界元離散化。

Moinfar[35-36]系統地發展了嵌入式離散裂縫模型，增加了傾斜裂縫的特性，并將其應用于全隱式組分模擬器中。由于兩個系統的物質平衡方程之間沒有直接關聯，需要通過非相鄰連接(NNC)實現相互之間的流體流動，即模型中需要計算的四種傳導率——裂縫系統中裂縫控制體之間的傳導率、基質網格與裂縫控制體之間的傳導率、相交裂縫之間的傳導率以及裂縫和井筒之間的傳導率。通過定義多種類型的傳導率，嵌入式離散裂縫模型可以動態地表征裂縫特征，在每個時間步末端更新裂縫的開度和滲透率，重新計算傳導率并應用到下一個時間步的流動方程中。Shakiba等人[37-38]利用微地震監測對嵌入式離散裂縫模型進行約束，同時應用組分模擬器UTCOMP更加精準而高效地對體積壓裂區域的縫網系統進行了數值模擬。

嵌入式離散裂縫模型繼承了雙重介質模型的優勢，不需要裂縫控制體與基質網格相鄰，也不要求基質和裂縫系統具有相同的網格尺寸。同時，該模型離散化表征了裂縫系統，不需要在模型中使用局部網格加密或非結構網格化。上述優勢使該模型可以在常規數值模擬器中應用，大大提高了計算效率和精確度。

3 發展趨勢

國外頁巖氣藏開發具有獨特的商業化模式，需要在衰竭式開采條件下將采收率最大化，其研究內容主要包括：①體積壓裂工藝參數對頁巖氣產能的影響，例如排量、支撐劑類型、壓裂液體系和溫度等；②完井工藝對頁巖氣藏生產動態的影響，例如分段多簇壓裂中的壓裂段數、段間距以及每段中的射孔簇數等因素，相應的數值模擬技術的突破與創新也著眼于這些方面。Mangrove的一體化工作流主要為單井模擬，而近年來對頁巖儲層中井間干擾現象的關注度越來越高；同時，隨著拉鏈式壓裂和井工廠作業的推廣，基于多井模型的非常規油藏數值模擬技術開始在國外發展。

不同于非常規頁巖儲層，長慶油田主要是致密砂巖儲層，在壓裂施工和射孔完井工藝上不需要進行有機碳含量、黏土含量等因素的篩選，即"甜點區"的優選；同時，長慶油田主要采用注水開發致密油，在開發方式上與國外存在根本性的差異。因此，體積壓裂數值模擬技術的發展，需要考慮超低滲、致密儲層中的應力敏感、啟動壓力梯度等特殊現象，并將研究重點放在注采井網的優化上來。

4 結論

(1)復雜縫網模型逐漸向復雜化、系統化的方向發展。非常規壓裂模型可以較好地模擬壓裂過程中縫網的形成及縫間應力干擾等問題，能夠并入油藏數值模型進行一體化分析，對現場問題提供了很好的解決方案。過于簡化的天然裂縫與水力裂縫的擴展準則，使裂縫擴展與實際有一定偏差，且對天然裂縫的描述仍然停留在統計意義上，很難準確描述地層裂縫系統，模擬的裂縫形態有很大不確定性。在以后的研究中，仍然需要大量的現場數據對非常規壓裂模型進行完善與修改，使之更加精確與高效。

(2)對于體積壓裂水平井的數值模擬，首先需要地質建模，通過力學模型與微地震手段的結合對復雜縫網進行描述；對于流動模型的選擇，如果選用雙重介質模型，則需要一種合理的粗化技術；如果選擇離散裂縫模型，則需要一種高效的非結構化網格技術以及局部網格加密方法對裂縫網絡準確刻畫。

(3)雙重介質模型假設的縫網形態較為均勻，不能準確描述實際情況；常規離散裂縫模型依賴于非結構網格化技術，具有較高的精確度，但是大量的網格數目造成了巨大的計算量；嵌入式離散裂縫模型綜合了兩種模型的優勢，較好地解決了一些復雜的問題。近年來發展較快的一體化工作流，基于連續介質模型和離散裂縫模型，綜合應用地質研究、微地震監測和壓裂參數等手段對體積壓裂水平井中的復雜縫網進行了表征，表現出良好的效果。

[1] STEPHEN R, DAVIDE E, RORY D, et al. Understanding hydraulic fracture geometry and interactions in the horn river basin through DFN and numerical modeling[C].SPE 137488, 2010.

[2] NEAL B N, IVAN G, MARISELA S,et al.Simulating hydraulic fracturing in real fractured rocks-overcoming the limits of pseudo 3D models[C].SPE 140480, 2011.

[3] NEAL B N, IVAN G, MARISELA S, et al. Stress shadowing and microseismic events: A Numerical Evaluation[C].SPE 147363, 2011.

[4]WENYUE X,MARC J, UTPAL G,et al.Wiremesh:A novel shale fracturing simulator[C].SPE 132218, 2010.

[5] WENYUE X, JAMES L, MIKE D. Quick estimate of initial production from stimulated reservoirs with complex hydraulic fracture network[C].SPE 146753, 2011.

[6]XIAOWEI W, OLGA K, CHARLES E C, et al. Modeling of hydraulic fracture network propagation in a naturally fractured formation[C].SPE 140253, 2011.

[7] KRESSE O, COHEN C, WENG X, et al. Numerical modeling of hydraulic fracturing in naturally fractured formations[J].Rock Mechanics,2011,15(5)516-535.

[8] GU H, WENG X, LUND J B, et al. Hydraulic fracture crossing natural fracture at non-orthogonal angles, a criterion, its validation and applications[C].SPE 139984, 2011.

[9] WARPINSKI N R, TEUFEL LW. Influence of geologic discontinuities on hydraulic fracture propagation[C].SPE 13224,1987.

[10] BLANTON T L. An experimental study of interaction between hydraulically induced and pre-existing fractures[C].SPE 10847, 1982.

[11]KRESSE O, WENG X, WU R, et al. Numerical modeling of hydraulic fractures interaction in complex naturally fractured formations[J]. Rock Mechanics and Rock Engieering,2013,46(3):555-568.

[12] CROUCH S L. Solution of plane elasticity problems by the displacement discontinuity method[J]. International Journal for Numerical Methods in Engineering, 1976, 10(2):301-343.

[13] OISON J E. Predicting fracture swarms——the influence of subcritical crack growth and the crack-tip process zone on joint spacing in rock[J]. Geological Society London Special Publications, 2004, 231(1):73-88.

[14] COHEN C E, ABAD C, WENG X, et al. Analysis on the impact of fracturing treatment design and re-servoir properties on production from shale gas reservoirs[C].IPTC 16400, 2013.

[15] KARSTEN P. A practical method for modeling fluid and heat flow in fractured porous media[C].SPE 10509, 1985.

[16] XU Z, CHANGAN D, FRANZ D, et al. Sensitivity studies of horizontal wells with hydraulic fractures in shale gas reservoirs[C]. IPTC 13338, 2009.

[17]AMIRMASOUD K D. Numerical simulation and modeling of enhanced gas recovery and CO2sequestration in shale gas reservoirs: A Feasibility Study[C]. SPE 139701, 2010.

[18] CRAIG L C, ELYEZER L, JAMES C E, et al. Re-servoir modeling in shale-gas reservoirs[C].SPE 125530, 2009.

[19] BARRY R. Accurate simulation of non-darcy flow in stimulated fractured shale reservoirs[C].SPE 132093, 2010.

[20] CHANGAN D, XU Z, BRAD M, et al. A workflow for integrated Barnett shale gas reservoir modeling and simulation[C].SPE 122934, 2009.

[21] CHANGAN M D, XU Z, LANG Z, et al. Modeling hydraulic fracturing induced fracture networks in shale gas reservoirs as a dual porosity system[C].SPE 132180, 2010.

[22] HUI M H， MALLISON B T. System and method for predicting fluid flow characteristics within fractured subsurface reservoirs:US[P],2010.

[23] HUI M H， MALLISON B T. System and method of grid generation for discrete fracture modeling: US[P],2011.

[24] MOHAMMAD K F, ABBAS F. Numerical simulation of water injection in fractured media using the discrete-fracture model and the Galerkin method[C].SPE 83633, 2003.

[25] FU Y, YANG Y K, DEO M. Three-dimensional, three-phase discrete-fracture reservoir simulator based on control volume finite element formulation[C].SPE 93292, 2005.

[26] MARCONDES F. An element-based finite-volume method approach for naturally fractured compositional reservoir simulation.13th Brazilian Thermal Sciences Meeting, Uberlandia, Brazil, December.

[27] CRAIG L C, XIAOWEI W, MARK G M, et al. Integrating microseismic mapping and complex fracture modeling to characterize hydraulic fracture complexity[C]. SPE 140185, 2011.

[28] CRAIG L C, TONY F, MICHAEL J W, et al. Seismic-to-simulation for unconventional reservoir development[C].SPE 146876, 2011.

[29] MAHSANAM M, CRAIG L C. A workflow for modeling and simulation of hydraulic fractures in unconventional gas reservoirs[C].SPE 153022, 2012.

[30] CRAIG L C, XIAOWEI W, HITOSHI O, et al. New algorithms and integrated workflow for tight gas and shale completions[C].SPE 146872, 2011.

[31]HASHEMI R, CAOTHIEN S, ROVIG J. Field measurement of viscosified fluids[C]. SPE 17462, 1988.

[32] COHEN C E, WENYUE X, XIAOWEI W, et al. Production forecast after hydraulic fracturing in naturally fractured reservoirs: coupling a complex fracturing simulator and a semi-analytical production model[C].SPE 152541, 2012.

[33] COHEN C E, ABAD C, WENG X, et al. Analysis on the impact of fracturing treatment design and re-servoir properties on production from shale gas reservoirs[C].IPTC 16400, 2013.

[34] LIYONG L, SEONG H L. Efficient field-scale si-mulation for black oil in a naturally fractured reservoir via discrete fracture networks and homogenized media[C].SPE 103901, 2006.

[35] ALI M, ABDOLJALIL V, KAMY S, et al. Deve-lopment of a novel and computationally-efficient discrete-fracture model to study IOR processes in na-turally fractured reservoirs[C].SPE 154246, 2012.

[36] ALI M, KAMY S, RUSSELL T J, et al. Coupled geomechanics and flow simulation for an embedded discrete fracture model[C].SPE 163666, 2013.

[37] JOSESA C F, MAHMOOD S, ALI M, et al. Implementation of a preprocessor for embedded discrete fracture modeling in an IMPEC compositional reservoir simulator[C].SPE 173289, 2015.

[38] MAHMOOD S, KAMY S. Using embedded discrete fracture model and microseismic monitoring data to characterize the complex hydraulic fracture networks[C]. SPE 175142, 2015.

編輯：王金旗

1673-8217(2017)03-0108-06

2016-10-22

唐子春，1991年出生，2014年畢業于中國石油大學(北京)石油工程專業，現為中國石油大學(北京)石油與天然氣工程專業在讀碩士研究生，主要從事超低滲透油藏工程與數值模擬研究。

國家科技重大專項“低滲、特低滲油藏水驅擴大波及體積方法與關鍵技術”(2017ZX05013-002)。

TE357.7

A