張鳳遠,鄒林君,崔 維,芮振華,李捷涵,潘 楊,張 偉,楊金智
(1.中國石油大學(北京)石油工程學院,北京 102249; 2.中國石油長慶油田分公司,陜西 西安 710018; 3.中國石油大學(北京)非常規油氣科學技術研究院,北京 102249; 4.中海油研究總院有限責任公司,北京 100028; 5.The Pennsylvania State University, State College, USA 16801)
中國常規油氣藏普遍進入開發中后期,頁巖油氣藏勘探[1-2]和開發技術受到重視[3-4],頁巖油氣藏開發進入快速發展時期[5]。水平井多段壓裂是實現頁巖油氣經濟開采的核心技術[6-9],壓裂效果評價和裂縫參數反演是頁巖油氣開發的重點和難點[10]。可以采用蘊含大量裂縫信息且十分經濟的壓裂液返排數據評價壓裂效果,美國Eagle Ford和Marcellus等頁巖油氣藏證明壓裂液返排數據在頁巖油氣開發中的重要作用。
人們建立多種非常規油氣藏滲流數學模型并反演裂縫參數。ABBASI M等采用單相水模型分析返排數據、反演裂縫存儲系數[11]。按照流體向裂縫竄流時間,CLARKSON C R等將致密油藏多段壓裂水平井壓裂液返排過程分為兩個階段,通過整合生產動態分析與歷史擬合方法反演裂縫半長和導流能力[12]。為表征壓裂液返排過程中的油水兩相流動,CLARKSON C R等建立半解析返排模型,提出分析兩相壓裂液返排數據的流動物質平衡方法和典型曲線方法[12]。文獻[13-14]將早期壓裂液返排看作油水兩相在封閉改造區域內的泄流過程,提出反演裂縫孔隙體積和滲透率的數學模型。KURTOGLU B等引入總流度和總產量概念,建立考慮油、氣、水三相壓裂液返排滲流模型[15]。這些模型與方法只考慮壓裂液返排初期,忽略儲層流體向裂縫竄流過程。
考慮壓裂液返排后期油藏中的油氣和濾失到儲層中的壓裂液向裂縫供給過程,CLARKSON C R等將試井分析中的動態動用區域(DDA)概念引入壓裂液返排動態分析過程,建立致密儲層半解析數學模型,同時考慮油水兩相在裂縫和基質內的復合流動,對壓裂液返排數據進行擬合,反演裂縫和儲層參數[16]。CLARKSON C R等將壓裂過程中的裂縫擴展模型與DDA模型相結合,對單孔介質[17]和雙孔介質[18]油藏分別建立改進的壓裂液返排數學模型,既可以分析壓裂液返排數據,還能計算壓裂液返排開始時裂縫和儲層中的壓力和飽和度分布,為壓裂液返排模型的初始化提供重要依據。隨體積壓裂技術的推廣,人們認識到兩翼縫壓裂液返排模型的局限性。CHEN Zhiming等[19]和JIA Pin 等[20]建立考慮復雜縫網的頁巖油藏壓裂液返排模型,實現油水兩相在不規則裂縫網絡內的流動模擬。ZHANG Fengyuan等對頁巖氣藏[21]和頁巖油藏[22]提出反演裂縫參數的壓裂液返排動態分析模型。這些壓裂液返排后期滲流數學模型大多用在給定的儲層和裂縫參數下進行正向模擬,不能反演裂縫參數,且缺少兼顧兩相流動和頁巖油氣藏復雜滲流機理的系統性裂縫參數反演方法。
以不穩定滲流理論為基礎,筆者建立一種適用于頁巖油藏和氣藏的壓裂液返排滲流數學模型。該模型不僅考慮裂縫—基質間油水或氣水兩相復合流動,而且考慮兩相流體在頁巖儲層內的復雜賦存和運移機理。引入擬壓力、擬時間函數,結合物質平衡方程和Duhamel原理,對模型進行半解析求解,提出分析壓裂液返排數據、反演裂縫初始孔隙體積和滲透率的直線分析法。
為模擬多段壓裂水平井壓裂液返排過程,建立裂縫內兩相流動與基質內兩相流相耦合的滲流物理模型(見圖1)。其中,km為基質滲透率;kf為壓力相關裂縫滲透率;xf為裂縫半長;wf為縫寬;ym為半縫間距;Rom為有機孔半徑;Rim為無機孔半徑;δw,im為無機孔內近壁區域的水層厚度;vnw為無機孔內近壁區域水相的流動速度;vbh為無機孔內體相區域烴類流體流動速度;vsh為體相區域和近壁區域交界面處的烴類流體滑脫速度。壓裂裂縫內充滿壓裂液(水相)及少量壓裂過程中殘留的油或氣,返排過程中裂縫和基質內為氣水或油水兩相流動。在壓裂液返排初期,壓裂液為井的主要產出流體,大部分產出水來自裂縫,少部分產出水來自壓裂過程中侵入到地層的壓裂液。隨返排進行,裂縫和基質之間壓差增大,產出流體主要來源于基質內油或氣向裂縫的竄流。

圖1 頁巖油氣藏二維兩相壓裂液返排滲流物理模型
頁巖儲層具有源儲一體的特點,儲層孔隙結構、巖石物性、流體賦存方式及運移機理與常規油藏的有明顯差別。為表征流體在頁巖基質內的復雜滲流機理,根據掃描電鏡圖像[23]將基質系統分為帶有納米級孔隙的有機質和無機質。不同流體在不同類型孔隙內部具有獨特的賦存和運移規律。
根據壓裂液返排物理模型,裂縫—基質系統內流體的流動方式由近及遠為:裂縫內壓裂液和少量殘余油或氣向井筒流動、儲層基質內流體受基質和裂縫之間壓差作用通過裂縫面進行供給、基質內流體向垂直裂縫面方向流動。壓裂液返排滲流模型可分為裂縫流動和基質流動,分別對兩種流動過程進行建模,在裂縫面進行壓力和流量耦合,建立壓裂液返排不穩定流動滲流數學模型。
基本假設條件:
(1)儲層均質等厚,流體在裂縫內流動遵循達西定律,忽略重力、毛細管力、井筒存儲效應和表皮效應。
(2)各條裂縫的幾何形態和裂縫參數相同,縱向穿透儲層,忽略裂縫端面以外的流體供給。
(3)裂縫和基質微可壓縮,油和水為微可壓縮流體,具有恒定的壓縮系數。氣體滿足真實氣體狀態方程。油、氣、水的黏度和體積因數隨壓力變化。
(4)有機孔和無機孔為平行排列的圓形微納米管,彼此之間沒有流體傳輸。有機孔和無機孔具有相同的孔隙壓力和應力敏感性。
(5)無機質具有親水性,壓裂液主要侵入基質的無機孔,無機孔內的水分子主要分布于近壁附近[24]。孔壁具有親水性,可以忽略液—固界面處的滑脫效應。油或氣賦存于無機孔的體相區域,油氣分子與水膜之間具有滑脫效應。
(6)有機質具有疏水性,可以忽略儲層改造過程中侵入到有機孔的壓裂液,有機孔內為單相氣或單相油。
(7)油相在有機孔內的賦存和流動狀態遵循頁巖油藏模型假設[25],即在連續流基礎上考慮滑脫效應和流體物性在孔壁區域和體相區域的不均勻分布特性。
(8)氣相在有機孔內的賦存和流動狀態遵循頁巖氣藏模型假設[21],即在連續流基礎上考慮滑脫效應、努森擴散、表面擴散和Langmuir單層吸附解吸規律。
首先建立定產量生產條件下裂縫流動模型,然后利用Duhamel原理獲得變產量生產條件下的模型解析解。基于頁巖氣藏[21]和頁巖油藏[22]壓裂液返排模型,建立適用于頁巖油氣藏的裂縫流動控制方程、初始條件和內外邊界條件:
(1)
pf(x,0)=pfi,
(2)
(3)
式中:下角標j為w時表示水相方程,j為o時表示油相方程,j為g時表示氣相方程;下角標i表示參數在初始條件下取值;x為沿裂縫方向的坐標;kf=kfiexp(-γf(pfi-pf)),γf為裂縫滲透率應力敏感模量;krj,f為裂縫中流體的相對滲透率;pf為裂縫壓力;φf為壓力相關的裂縫孔隙度,φf=φfiexp(-Cf(pfi-pf)),Cf為裂縫壓縮系數;μj,f為裂縫內流體黏度;Bj,f為裂縫內流體體積因數;t為時間;qj為從裂縫流入井筒的流體流量;Cej,f為裂縫有效壓縮系數,表達式為
(4)
式中:Sj,f為裂縫內流體飽和度;qsj為基質向裂縫的竄流量;Cj,f為裂縫內流體壓縮系數;h為縫高。
裂縫流動控制方程通過式(4)耦合基質流動過程,將基質向裂縫的竄流等效為裂縫有效壓縮系數的變化。式(1-3)中壓力相關的裂縫和流體物性參數使控制方程具有很強的非線性,通過引入擬壓力和擬時間函數(見表1),可將裂縫流動控制方程線性化,利用分離變量法得到解析解:

表1 裂縫和基質的擬壓力和擬時間定義
(5)

考慮式(5)的解析解過于復雜,無法直接用于參數反演,引入產量規整化擬壓力RNPj=(mj,f(pfi)-mj,f(pwf))/qj,其中pwf為井底流壓。利用文獻[26]的方法對式(5)進行簡化,得到定產量生產條件下的早期近似解和晚期近似解:
(6)
(7)
式中:Vfi為裂縫初始孔隙體積,Vfi=2xfwfhφfi;ε、α、β分別為單位轉換因數,ε=1.932,α=1,β=0.977。
實際返排過程中,井通常以變產量的方式進行生產,定義兩種形式的疊加擬時間tspj,1、tspj,2:
(8)
(9)
利用疊加原理,可以得到變產量生產條件下的早期近似解和晚期近似解:
(10)
(11)
2.2.1 滲流基本模型
類似于裂縫流動控制方程,基質流動控制方程為
(12)
根據裂縫面壓力連續性條件,將基質模型內邊界處的壓力設為平均裂縫壓力,耦合裂縫流動和基質流動模型。初始條件和邊界條件分別為
pm(y,0)=pmi,
(13)
(14)

(15)
式中:Sj,m為基質內流體飽和度;Cj,m為基質內流體壓縮系數;Cd為氣體解吸附貢獻的壓縮系數,對于液體,Cd=0。
式(12-14)中與壓力相關的基質和流體物性參數使得基質流動控制方程具有很強的非線性,引入擬壓力和擬時間函數(見表1),可將基質流動控制方程線性化。線性化的基質流動模型具有非齊次的內邊界條件,即裂縫平均壓力隨時間變化。利用分離變量法和Duhamel原理[21],得到基質方程的解析解為
(16)
(17)

2.2.2 微觀流動模型
2.2.2.1 有機孔單管
基于假設條件,有機孔壁具有疏水特性,孔內流體為單相氣或單相油。結合頁巖氣藏微觀流動模型[21],氣體在單根有機微納米管內的總質量流量Jg,om表示為
(18)
式中:Kn,om為有機質的努森系數;Jsl,om為氣體滑脫效應導致的質量流量;Jk,om為努森擴散導致的質量流量;Jsf,om為表面擴散導致的質量流量[21]。
若儲層為頁巖油藏,在連續流基礎上,考慮應力敏感、滑脫效應和流體物性在孔壁區域和體相區域的不均勻分布特性等,結合頁巖油藏微觀流動模型[25],油在單根有機微納米管內的總質量流量Jo,om表示為
(19)
式中:Rom=Rom,iexp[0.5(Cm-γm)(pmi-pm)];ρo,m和ρo,om分別為體相區域和近壁區域的油相密度;δo,om為有機孔內充滿油相時的近壁區域厚度;μo,m和μo,om分別為體相區域和近壁區域的油相黏度;λo,om為油相在有機孔壁的滑脫長度。
為與滲流基本模型結合,利用達西定律將油相或氣相的質量流量轉化為頁巖有機質的表觀滲透率kh,om:
(20)
式中:φom為有機質孔隙度;τom為迂曲度。
2.2.2.2 無機孔單管
基于假設條件,無機孔內既包含烴類流體又包含水。由于無機孔壁具有親水性,烴類流體占據體相區域,水相分布于近壁區域[27]。根據Hagen-Poiseulle方程,建立考慮滑脫效應的流動方程:
(21)
(22)
式(21-22)中:Rim=Rim,iexp[0.5(Cm-γm)(pmi-pm)];μh,m和μw,im分別為體相區域烴類流體和近壁區域水相流體的黏度;δw,im為近壁區域的水層厚度;c1和c2為待定系數。
考慮烴類流體和水在體相區域和近壁區域交界面處的滑脫效應,同時忽略水在孔壁液—固界面上的滑脫效應,流動方程的內外邊界條件為
(23)
vbh|r=Rim-δw,im=vnw|r=Rim-δw,im+vsh,
(24)
vnw|r=Rim=0。
(25)
當烴類流體為油時滑脫速度為vso,當烴類流體為氣時滑脫速度為vsg,分別用油或氣在固體界面上的滑脫速度[27]近似表示為
(26)
(27)
式(26-27)中:δv為切向動量調節因數,對于漫反射通常取為1;λg為氣體平均自由程;λo為油在體相區域和近壁區域交界面處的滑移長度,λo=10-9×[0.387exp(-2×109Rim/32.47)+0.74];vbg、vbo分別為無機孔內體相區域氣體和油的流動速度。
將式(21-22)代入式(23-25),可以得到無機孔內烴類流體和水在體相區域和近壁區域的流速剖面分別為
(28)
(29)
對體相區速度剖面從0到Rim-δw,im進行積分,對近壁區速度剖面從Rim-δw,im到Rim進行積分,可以得到烴類流體和水的體積流量,將質量流量表示為
(30)
(31)
式中:Jbh,im為烴類流體在無機孔內的質量流量;Jnw,im為水在無機孔內的質量流量;ρh,m和ρw,im分別為無機孔內體相區的烴類流體和近壁區域水的密度。
為與滲流基本模型結合,可以利用達西定律將油或氣的質量流量轉化為頁巖無機質的表觀滲透率kh,im:
(32)
式中:φim為無機質孔隙度;τim為無機質內孔隙的迂曲度。
基于有機孔和無機孔平行排列、無傳質的假設,基質內烴類流體的綜合表觀滲透率為
kh,m=ψkh,om+(1-ψ)kh,im,
(33)
式中:ψ為基質內的有機質體積分數。
基質內水相只存在于無機孔,水相的綜合表觀滲透率為無機質內的水相表觀滲透率:
(34)
式中表觀滲透率計算過程中涉及到的近壁區域水層厚度δw,im與基質含水飽和度相關。
2.3.1 平均壓力和飽和度
基于文中的裂縫模型和基質模型,提出考慮頁巖油氣藏復雜滲流機理的裂縫參數反演方法——直線分析法:利用式(10-11)構建劃分流動段的診斷曲線和求取裂縫參數的特征曲線,通過提取特征曲線上的斜率和截距計算裂縫初始孔隙體積和初始滲透率。在直線分析過程中,需要用到裂縫和基質的平均壓力、飽和度等參數。
為計算擬壓力和擬時間,需要計算裂縫和基質內平均壓力、平均飽和度隨時間的變化關系。分別對裂縫和基質建立物質平衡方程:
(35)
(36)

Vj,mi=yjxfhφmi,
(37)
(38)
式(37-38)中:yj為流體在基質內的調查距離。

2.3.2 直線分析法


圖2 二維兩相流動診斷曲線
(39)
由式(11)可知,壓裂液返排晚期RNPj和tspj,2呈線性關系,繪制在dRNPj/dlntspj與tspj的雙對數曲線(診斷曲線)上呈斜率為1的直線(見圖2[22]),與試井分析中的邊界控制流類似[28]。由于這一階段反映的是壓裂液返排后期裂縫內的壓力擾動傳播到裂縫邊界之后的流動狀態,稱為擬裂縫邊界控制流。若已知裂縫的其他基本參數,則可以利用解釋擬裂縫邊界控制流階段壓裂液返排晚期數據的直線分析法,獲得裂縫初始孔隙體積和滲透率,具體步驟:(1)將壓裂液返排晚期數據繪制在診斷曲線上,提取斜率為1直線對應的數據時間段;(2)將返排數據以tspj,2為橫坐標、RNPj為縱坐標繪制在直角坐標系上(擬裂縫邊界控制流特征曲線),提取特征曲線上對應步驟(1)中時間段的直線斜率(m2)和截距(b);(3)計算裂縫初始孔隙體積Vfi,j和滲透率kfi,j,2:
(40)
(41)
為了進行直線分析,需要將Vfi作為輸入參數代入物質平衡方程,求解平均壓力和平均飽和度,繪制診斷曲線和特征曲線。Vfi作為待求解的變量之一,需要采用迭代求取裂縫參數:
(1)根據現場經驗,給定一個裂縫初始孔隙體積Vfi作為迭代初始值。

(3)根據平均壓力和平均飽和度,計算擬壓力和擬時間(見表1)。
(4)繪制每一相流體對應的診斷曲線(見圖2[22]),提取斜率為1/2和1的直線段對應的數據時間段。
(5)繪制擬裂縫線性流和擬裂縫邊界控制流的特征曲線,對兩個不同的流動段提取直線段的斜率和截距,代入式(39-41)求取Vfi和kfi。
(6)將解釋的Vfi與迭代初始值進行對比。若相對誤差大于設定公差,則將新求出的Vfi作為迭代初始值代入步驟(2)。若相對誤差小于設定公差,則迭代收斂,輸出計算的Vfi和kfi作為解釋結果。
利用裂縫和基質控制方程的解析解實現對裂縫參數的反演,由于平均壓力和飽和度的計算過程涉及物質平衡方程的迭代求解,壓裂液返排滲流模型為半解析模型,利用數值模擬對滲流模型進行驗證。目前的數值模擬軟件無法考慮基質內存在兩相流條件下的微觀滲流機理,無法對烴類流體和水相賦不同的表觀滲透率函數。為驗證滲流模型準確性,利用只考慮基質內烴類流體流動(假定基質內水相不流動)的數值模型對采取相應處理的半解析模型進行驗證。

圖3 直線分析輔助數據和曲線
基于式(10-11),分別模擬水相和油相壓裂液返排數據,繪制不同流動階段下的直線分析法特征曲線(見圖4)。提取直線段的斜率和截距,得到擬裂縫線性流階段水相和油相數據解釋的kfi,1分別為3.212、3.274 μm2,擬裂縫邊界控制流階段水相和油相數據解釋的Vfi分別為15.404×103、15.546×103m3,水相數據解釋的kfi,2為3.287 μm2。解釋結果與模型的輸入數據3.000 μm2和14.158×103m3相比誤差不超過10%,滿足工程上對參數反演的精度需求。擬裂縫邊界控制流階段油相數據解釋的kfi,2為4.255 μm2,與輸入數據誤差較大,原因是基質竄流使得油相方程具有強非線性。參數反演過程中,以水相數據解釋的kfi結果為準。

圖4 早期和晚期流動階段直線分析特征曲線
以美國Marcellus頁巖氣藏一口多段壓裂水平井為例,利用文中半解析模型和裂縫參數反演方法對該井的壓裂液返排數據進行分析,評價壓裂效果。該井分為22段壓裂,每段5簇。壓裂后進行24 d的返排,隨后連續生產8個月。排采階段的井底流壓、產量歷史及參數反演基本輸入參數由ZHANG Fengyuan等[21]給出。ZHANG Fengyuan等[21,29]采用基質—裂縫兩相復合流動模型和裂縫兩相—基質單相模型,對該井壓裂液返排數據進行分析,前者忽略流體在儲層基質內的復雜運移和賦存機理;后者雖然考慮頁巖氣復雜滲流機理,但假設基質內為單相氣體流動,忽略侵入儲層中的壓裂液對裂縫的供給過程。
首先給定一個Vfi作為迭代初始值,按照分析步驟求取新的Vfi對初始值進行更新。Vfi收斂后,得到裂縫平均壓力和平均含水飽和度曲線(見圖5(a)),以及基質內氣相和水相向裂縫竄流量曲線(見圖5(b))。初期產量以裂縫內流體為主,后期產量以基質供給為主。考慮復雜運移和賦存機理條件下氣相和水相表觀滲透率隨壓力的變化關系曲線(見圖5(c)),氣相表觀滲透率隨壓力降低而逐漸升高,與基于單相氣體假設得到的規律[21]一致。水相表觀滲透率隨壓力降低而逐漸減小,與現場認識吻合,即壓裂液返排后期侵入地層的壓裂液越來越難產出。對于水相和氣相返排數據,分別繪制不同流動階段下的特征曲線(見圖5(d-f))。由圖5(d-f)可以看出,擬裂縫線性流階段水相和氣相數據計算的裂縫初始滲透率分別為1.027和1.079 μm2;擬裂縫邊界控制流階段水相和氣相數據計算的裂縫初始滲透率分別為1.011和1.456 μm2,裂縫初始孔隙體積分別為4.69×103和5.31×103m3。不同相態、不同流動階段數據解釋的裂縫參數較為吻合且在合理的范圍內,證明該方法現場應用的可靠性。

圖5 文中模型Marcellus頁巖氣井現場應用結果
(1)建立一種普遍適用于多段壓裂水平井油藏和氣藏的壓裂液返排半解析滲流數學模型,不僅考慮裂縫—基質間油水或氣水兩相復合流動,而且考慮兩相流體在頁巖儲層內的復雜滲流機理。
(2)提出一種劃分兩相流動階段的診斷曲線和反演裂縫參數的直線分析法,能夠對裂縫初始孔隙體積和裂縫初始滲透率進行準確反演。壓裂液返排過程可分為擬裂縫線性流和擬裂縫邊界控制流,在診斷曲線上分別呈斜率為和1的直線。
(3)受流體在儲層中的復雜滲流機理影響,氣相和水相表觀滲透率隨壓力降低呈不同的變化趨勢。受現場數據測量精度影響,擬裂縫邊界控制流的直線段持續時間更長且更明顯,該階段的數據更適宜進行裂縫參數反演。