致密砂巖儲層微觀孔隙結構的分形研究

王 妍,馬 瑤,王若谷,陳永振,孫博亞,李文厚,吳 越

(1.西安石油大學 石油工程學院, 陜西 西安 710065;2.陜西延長石油(集團)有限責任公司 研究院,陜西 西安 710006;3.陜西省能源化工研究院,陜西 西安 710069;4.西安科技大學 地質與環境學院,陜西 西安 710054;5.西北大學 地質學系，陜西 西安 71069;6.長慶油田分公司 第十一采油廠,甘肅 慶陽 745000)

儲層非均質性決定著儲層物性的好壞,對儲層預測也有重要影響。儲層的孔隙分布于三維空間,但又不充滿三維空間,空間分布具有統計自相似性。分形理論提供了一種簡單的方法，將大尺度的變化與小尺度的變化聯系起來,反之亦然[1-2]。由于其在表征復雜、不規則特征方面的優勢,已逐漸被引入到分析孔隙結構的復雜性和評價儲層的非均質性中[3-4]。自Katz and Thompason利用掃描電鏡對砂巖的孔隙進行分形分析以來[5],眾多學者研究認為，孔隙結構在一定的孔徑范圍內表現出良好的分形特征,其特征可用分形維數來表征[6-11]。

分形維數求解方法多樣,如Friesen等[12]和Yang等[13]用導數法得到了分形維數與飽和度線性方程的關系,并計算了分形維數;賀承祖等[8]和張國輝等[14]從毛管壓力和濕相飽和度的線性方程出發,推導出巖石的分形維數;Krohn[4]采用掃描電子顯微鏡觀察巖石樣品的截面并求取分形維數;賈芬淑等[15]創建了根據圖像處理技術測定砂巖孔隙結構的分形維數的方法;彭瑞東等人[16]提出了基于CT圖像的灰度圖像分形維數的計算方法。

不同方法計算的巖石分形維數也不同,目前常用的計算方法主要有兩類,一類是基于壓汞資料的分形維數求解[8,14,17-22],另一類是基于圖像數據的分形維數求解[4,15-16, 23-29]。相對于圖像數據,壓汞資料便于獲取。現有研究大多根據毛細管壓力曲線來計算儲層分形維數,雖然有其優點,但也存在數據平滑等缺點,因此它們不能顯示細節和局部非均質性[2];而基于成像分析技術的分形與多重分析技術為描述儲層孔隙空間分布提供了新方法,從而為定量評價儲層的儲集性能提供新的認識[16, 25-26]。

目前，對儲層微觀結構的分形理論研究主要集中在中、高滲等常規儲層方面,針對致密儲層，尤其是致密砂巖的研究較少。為此，本研究以鄂爾多斯盆地隴東地區長7致密砂巖油藏的真實巖心作為研究對象,分別利用基于壓汞資料的MIFA法及基于CT掃描技術的分形盒維數研究，對孔隙空間的分形特征進行分析和描述,并對兩種方法計算的分形維數之間的差別進行分析與探討。

1 分形維數的表征方法

1.1 壓汞資料法

MIFA法是以毛管壓力曲線為資料,采用分形理論研究儲層孔隙結構的一種計算方法。非潤濕相注入多孔介質時,表面張力阻止非潤濕相進入孔隙,需額外施加壓力才能將非潤濕相液體注入巖石孔隙中,這個額外壓力即毛管壓力。1998年，賀承祖等根據毛細管壓力曲線建立了巖樣中的潤濕相流體飽和度與毛細管壓力之間的關系[8],即

(1)

式中:pc為任意孔隙半徑所對應的毛細管壓力,MPa;pmin為儲層巖石最大孔徑rmax對應的毛細管壓力,即突破毛細管壓力,MPa;Sw為毛細管壓力pc時儲層巖石中潤濕相的飽和度;D為分形維數。

由式(1)可知,如果巖石孔隙結構符合分形特征,則儲層中的潤濕相飽和度與毛管壓力之間應滿足冪函數關系。對式(1)兩端取對數,可得

lgSw=(3-D)lgpmin+(D-3)lgpc。

(2)

由式(2)可知,潤濕相飽和度與毛管壓力在雙對數坐標系下為一條直線,如儲層巖石孔隙具有分形結構的性質,則根據毛細管壓力資料可知,lgSw與lgpc之間呈線性相關關系,直線的斜率為D-3,據此可以得到巖石孔隙的分形維數D。

1.2 CT圖像分析法

傳統的儲層地質室內實驗，如普通薄片、鑄體薄片、掃描電鏡等能夠揭示巖石微觀孔隙結構,但是具有一定局限性。如普通薄片和鑄體薄片得到的二維圖像只是對樣品單一層面的取樣,無法反映立體結構;掃描電鏡獲得的圖像空間延續性差[30-31]。近年來，CT技術在石油行業廣泛應用,可以精確檢測儲層樣品的內部結構,反映儲層微觀孔隙的三維高精度結構[32-34]。通過對CT圖像進行分析,也可計算出分形維數。CT圖像為灰度圖像,目前常用手段是將灰度CT圖像轉化為黑白二值圖像再計算分析。

在實際操作中,盒維數由于其簡易性和可操作性被頻繁使用,基于二值圖像的盒維數算法(box-counting method,BC算法)主要過程步驟如下。

1) 讀取一幅M×N像素的二值圖像,每個像素點非即1,通常白色像素表示孔隙;

2) 用大小s×s像素的網格(盒子)去完全覆蓋二值圖像,分割比例r=s/M;

3) 任何一個網格塊只要包含目標對象就被認為是被占用(可計數)的網格塊,進而統計總共被占用的網格塊數Nr;

4) 改變s的大小不斷重復上述過程,可得到一系列(Nr,r);

5) 在雙對數坐標中用最小二乘法擬合[logNr, log(1/r)],所得的斜率即為該圖像的分形維數。

2 結果與討論

2.1 孔隙描述

所選樣品碎屑顆粒以石英為主,長石次之,長石顆粒以斜長石居多,蝕變中等偏多,常可見長石溶孔。粒間填隙物以泥鐵質為主,可見硅質膠結和鈣質膠結。硅質膠結較少,以石英次生加大為主,鈣質膠結發育較多,以方解石膠結為主。黏土礦物主要以伊利石、綠泥石為主。樣品孔隙發育一般,多為孤立狀,連通性較差,儲集空間類型以非組構選擇性孔隙(其中，粒間溶孔體積分數為1.5%，長石溶孔體積分數為0.8%)及少數組構選擇性孔隙(如粒間孔體積分數為1%)為主(見圖1),樣品還發育微裂縫,平均縫寬約16 μm，樣品可見伊利石、綠泥石及自生石英充填孔隙。對樣品進行薄片觀察可知，粒間溶孔體積分數為1.5%,長石溶孔體積分數為0.8%,粒間孔體積分數為1%。樣品壓汞法測定孔隙度7.94%,滲透率0.102×10-3μm2,為致密儲層。

2.2 基于MIFA法的分形維數

在三維歐氏空間中，孔隙結構的分形維數介于2與3之間。分形維數越接近2,孔喉表面越光滑、孔隙結構越均勻、均質性越強。分形維數越接近3,孔喉表面越不光滑、孔隙結構越不均勻、非均質性越強。

圖1 主要孔隙類型Fig.1 Main pore types

根據孔隙結構分形維數的適用條件,驗證了該致密砂巖樣品符合分形特征,即樣品的含水飽和度與壓汞毛管壓力在雙對數坐標系上是一條顯著的直線(見圖2),擬合直線斜率2.098 24,直線判定系數0.969 29。這說明孔隙結構整體具有“單段型”分形特征,其孔喉大小分布較集中,呈似單峰形態(見圖2)。無論是壓汞法測量的該樣品微觀均質系數0.389,還是通過對研究區該層位22個樣品統計的均質系數0.368 3,均表明儲層非均質性變化不明顯,為中等非均質性,與分形維數反映的微觀孔隙結構的非均質程度結果相吻合。

圖2 MIFA法分形特征及孔喉半徑特征Fig.2 Fractal characteristics of MIFA method and pore& throat radius characteristics

設lgSw=y,lgpc=x,(3-D)=a,lgpmin=b,則公式(2)變為y=ab-ax。利用回歸擬合法可得到精確的分形維數D,D=2.507。

2.3 基于BC算法的分形盒維數

應用MATLAB編程對樣品CT圖像進行分形盒維數的計算。

1)圖像預處理。為了計算分形維數,需要對CT圖像進行一系列處理,包括亮度調節(見圖3,4a)、圖像剪裁(見圖4b)、閾值分割等,從而獲取可進行分形分析的二值圖像。

使用Matlab程序將圖像RGB值乘以系數n,可得到不同亮度增強的圖像。如圖3所示,n為1.5～2.5時，圖像亮度依然較暗(見圖3b～d),n為3.5～4.5時(見圖3f～h)圖像亮度已經過亮,因此選擇n=3時(見圖3e)亮度調節后的圖像。圖像亮度調節后的圖像剪裁不僅可剪裁掉無效部分,保留圖像中心圓形區域,還可消除CT掃描過程中圖像邊緣偽影的影響。最后，在閾值分割過程中,應用Otsu算法進行分割,所獲得的二值圖像如圖4c所示。

圖3 不同亮度的CT圖像Fig.3 CT images of different brightness

圖4 CT圖像的預處理過程Fig.4 Pre-processing of CT images

2)分形分析。在三維歐氏空間中，平面內的分形曲線的分形維數介于1與2之間。因此二值圖像的盒維數范圍為1～2。數字圖像盒維數的計算是基于對像素矩陣的分析,表象上看是反映了研究對象對平面的占有情況,分形維數越大,越趨向于充滿整個平面,本質上來說，是反映了關心區域趨向于沿整個圖像范圍分布還是趨向于沿某條直線分布[24]。

作為示例,利用Matlab編程語言研究了鄂爾多斯盆地隴東地區某長7油層組致密砂巖孔隙結構的分形特征。通過對巖石樣品斷層掃描,獲得了1 000張左右的32位切片圖像,圖像分辨率為0.724 μm,圖像間隔1個像素點。CT圖像中包括粒間孔、溶孔等關心區域,通過一系列圖像分析處理,將這些孔隙提取出來,并轉化為黑白位圖。

灰度圖像被轉化成黑白位圖后,可通過對圖像中黑白像素的統計來計算孔隙度(白色代表孔隙,黑色代表固體顆粒),計算公式如下:

(3)

其中：Pixelswhite代表孔隙總像素個數；Pixelwhote為圖像的總像素數。

實測孔隙度為9.86%,計算孔隙度為10.02%，誤差為1.6%，結果表明，所選用的閾值分割算法是合理有效的。

本次研究通過應用幾何序列選取盒子的尺寸序列。具體來說,剪裁后圖像尺寸為400×400,盒子序列s滿足1≤s≤M/2,即尺寸盒子范圍s={1,2,4,8,16,32,64,128}。計算所得CT圖像盒維數為1.643。可以認為,這一維數反映了孔隙在平面內擴展分布的規律。數字圖像盒維數越大,表明孔隙分布越不規則,也就表明儲層非均質越強。因此,本樣品孔隙結構非均質性一般,與基于MIFA法計算的分維數代表的物理意義相同。

3 結 論

1)通過壓汞實驗,擬合求取了鄂爾多斯盆地隴東地區某長7油層組致密砂巖樣品的分形維數。

2)通過CT掃描獲取該砂巖樣品的高分辨率孔隙結構的灰度圖像,繼而通過圖像分析技術對其二值分割,獲取該圖像的黑白位圖,進而計算出黑白位圖的計算孔隙度φ2,與實驗孔隙度φ1相比,誤差在5%以內,驗證了閾值分割的合理性。

3)無論是基于MIFA法計算的分形維數還是基于CT圖像計算的分形盒維數,都能反映孔徑的分布特征。這兩類分形維數的計算與同層位22塊樣品計算的均值系數結果相吻合,表明該樣品孔隙結構較均勻,非均質程度一般。