風(fēng)城油田齊古組油砂細觀結(jié)構(gòu)和滲流特征

龐惠文，金衍，高彥芳，王琪琪

（1.中國石油大學(xué)（北京）a.油氣資源與探測國家重點實驗室；b.理學(xué)院；c.石油工程學(xué)院，北京 102249；2.西北大學(xué) 地質(zhì)學(xué)系，西安 710069）

稠油作為非常規(guī)油氣資源的重要組成部分，是中國未來常規(guī)油氣資源的有效接替[1]。中國西部和東北部地區(qū)稠油資源豐富，其中新疆準噶爾盆地風(fēng)城油田稠油儲量約為3.6×1012t[2]，開采潛力大，但平均黏度超過10 000 mPa·s（50 ℃），給高效開發(fā)帶來巨大挑戰(zhàn)[3]。蒸汽輔助重力泄油（SAGD）技術(shù)是中國目前常用的稠油開發(fā)技術(shù)，已廣泛應(yīng)用于包括風(fēng)城油田在內(nèi)的多個油田。風(fēng)城油田上侏羅統(tǒng)齊古組砂巖儲集層屬于陸相河道沉積，具有非均質(zhì)性強、泥質(zhì)含量高、原油黏度大等特點，在應(yīng)用SAGD 技術(shù)開采時，預(yù)熱周期長，蒸汽能耗大。目前形成了以快速啟動技術(shù)為代表的注水微壓裂改造技術(shù)，即在循環(huán)注蒸汽前，先對儲集層進行注水微壓裂改造，以形成溝通注采井間均勻的高滲擴容帶，為后期蒸汽腔的形成和發(fā)育提供有效通道[4]。但由于風(fēng)城油田儲集層孔隙結(jié)構(gòu)和滲流規(guī)律復(fù)雜，導(dǎo)致部分井改造效果不理想，亟需針對油砂細觀結(jié)構(gòu)和滲流特征開展精細化研究，為微壓裂改造方案的優(yōu)化和稠油油藏的高效開發(fā)提供理論依據(jù)。

針對儲集層細觀結(jié)構(gòu)和滲流特征的研究方法主要包括間接法、直接法和數(shù)字巖心重構(gòu)3 類。研究區(qū)齊古組砂巖儲集層結(jié)構(gòu)疏松，易受外部擾動，若采用以光學(xué)顯微鏡和掃描電鏡為主的定性或半定量細觀結(jié)構(gòu)分析方法，以及穩(wěn)態(tài)和瞬態(tài)滲透率實驗室測試方法，在制樣和壓力保持方面均存在很大困難[5]。因此，無損CT 配套重構(gòu)分析技術(shù)成為定量表征油砂細觀結(jié)構(gòu)和滲流特征的有效途徑，但準確區(qū)分巖心不同組分的顆粒接觸方式、孔隙網(wǎng)絡(luò)建模以及孔隙尺度流動模擬，是數(shù)字巖心技術(shù)面臨的主要難題[6]。

構(gòu)建數(shù)字巖心的基礎(chǔ)是對二維巖心圖像的處理和分析，針對CT 獲得的圖像，需要進行圖像預(yù)處理和圖像分割。圖像預(yù)處理方法主要包括中值濾波、各向異性擴散濾波、全變差濾波、非局部濾波等[7]。圖像分割方法主要包括全局閾值法、局部閾值自適應(yīng)法、區(qū)域生長法、可變形表面水平集法、概率模糊聚類法、貝葉斯法和分水嶺分割法[8]。數(shù)字巖心建模方法主要有物理實驗法和數(shù)值重建法2 種，皆可以用來表征孔隙和基質(zhì)的空間分布關(guān)系[9]。物理實驗法通過對有序的巖心二維斷層掃描，直接進行三維重構(gòu)；數(shù)值重建法以少量二維圖像為基礎(chǔ)，運用特定數(shù)值方法構(gòu)建數(shù)字巖心。在重構(gòu)的數(shù)字巖心上，可以建立孔隙網(wǎng)絡(luò)模型來分析孔隙連通情況，孔隙網(wǎng)絡(luò)構(gòu)建方法主要包括多向掃描法、居中軸線法、Voronoi 多面體法、最大球體積法等[10]，其中居中軸線法較為常用。開展孔隙尺度流動模擬的方法主要分為3 類：第一類是基于連續(xù)性假設(shè)的宏觀尺度求解法，即依據(jù)經(jīng)典力學(xué)理論，采用傳統(tǒng)計算流體動力學(xué)直接求解Navier-Stokes 方程[11]；第二類是基于非連續(xù)假設(shè)的介觀尺度求解法，即依據(jù)統(tǒng)計物理學(xué)理論，通過追蹤粒子微團運動行為來描述流體流動特征[12]；第三類是基于非連續(xù)假設(shè)的微觀尺度求解法，該方法依據(jù)統(tǒng)計物理學(xué)直接描述分子或原子的運動[13]。上述3種方法實質(zhì)是對同一種物理現(xiàn)象的不同表征形式，三者可以通過一定的方式進行轉(zhuǎn)化[14]?，F(xiàn)階段第一類方法已經(jīng)被證明具有潛在應(yīng)用價值，由于其與實驗觀測值的相似性獲得了初步驗證，是較為常用的方法。目前數(shù)字巖心技術(shù)已廣泛應(yīng)用于包括致密砂巖、頁巖、砂礫巖等巖心的孔隙結(jié)構(gòu)和滲流特征的研究中[15-17]，但采用該方法對油砂孔隙空間的研究還未見相關(guān)報道。

本文以風(fēng)城油田侏羅系齊古組油砂為研究對象，通過CT 物理實驗法重構(gòu)了不同尺寸油砂的數(shù)字巖心，提出油砂圖像中顆粒、瀝青質(zhì)和孔隙3 種介質(zhì)圖像分割方法，建立包含與不包含瀝青質(zhì)的2 種油砂骨架結(jié)構(gòu)模型、孔隙空間模型和三維孔隙網(wǎng)絡(luò)模型，進行油砂孔隙尺度流動模擬，分析不同數(shù)字巖心的流場和壓力場分布特征，揭示風(fēng)城油田天然油砂樣品中瀝青質(zhì)相變對滲透率的影響機制，為明確風(fēng)城油田油砂孔隙結(jié)構(gòu)和滲流特征提供理論依據(jù)。

1 樣品處理

實驗樣品取自風(fēng)城油田侏羅系齊古組井下油砂巖心。由于油砂結(jié)構(gòu)疏松、對溫度敏感、易擾動，在冬季溫度較低時取心，運輸過程采用干冰降溫，到達實驗室后在-20 ℃的冷柜中儲存，取心、運輸和樣品保存均在低溫下進行。微米CT 實驗的樣品采用直徑為25.0 mm、高度為50.0 mm 的標(biāo)準巖心柱，通過液氮切割的方式制備樣品，制樣過程盡可能避免摩擦和震動，以降低對油砂巖心結(jié)構(gòu)的影響[18]。

在掃描過程中采用了Xradia 公司生產(chǎn)的Microxct-200微米CT，分辨率為0.5～35.0 μm，通常粗粒砂巖要求圖像分辨率為15.0～20.0 μm，疏松砂巖要求圖像分辨率為5.0～15.0 μm[19]，因而本文實驗采用的設(shè)備能夠準確識別樣品中的孔隙，并計算其大小，分析其分布情況。為了更加全面地表征油砂樣品的細觀結(jié)構(gòu)特征，對樣品進行了不同精度的多次掃描。先掃描直徑為25.0 mm、高度為50.0 mm 的標(biāo)準巖心柱，獲得3 403張分辨率為12.0 μm的灰度圖像；再對樣品中心直徑為10.0 mm、高度為6.5 mm 的圓柱形區(qū)域精細掃描，獲得1 306 張分辨率為5.0 μm 的灰度圖像；最后針對2組圖像分別進行圖像處理與數(shù)字巖心重構(gòu)。

2 油砂孔隙結(jié)構(gòu)數(shù)字巖心表征

數(shù)字巖心表征流程主要包括以下3 部分：使用高分辨率CT 獲取巖心灰度圖像，對圖像進行處理與分割，數(shù)字巖心三維重構(gòu)（圖1）。但油砂的特殊性在于其包括顆粒、瀝青質(zhì)和孔隙3 種介質(zhì)，這使得傳統(tǒng)采用單一閾值劃分孔隙和骨架的方法不再適用。因此，如何準確區(qū)分出顆粒、瀝青質(zhì)和孔隙3 種介質(zhì)，建立適用于油砂3 種介質(zhì)圖像分割方法，是油砂數(shù)字巖心構(gòu)建的關(guān)鍵。

圖1 基于CT的油砂數(shù)字巖心重構(gòu)流程Fig.1.Reconstructing process for digital oil sand cores based on CT scan

2.1 介質(zhì)圖像分割

以致密砂巖為代表的多孔介質(zhì)，通常只包括孔隙和顆粒2部分[20]，其孔隙和顆粒灰度值相差較大，只要選取合適的閾值，就能將二者準確區(qū)分。油砂包括孔隙、顆粒和瀝青質(zhì)3 種介質(zhì)，顆粒分選性較差，瀝青質(zhì)分布無規(guī)律，不同介質(zhì)之間接觸關(guān)系復(fù)雜，且瀝青質(zhì)和顆粒的灰度差異較小，區(qū)分度過低，給兩者的劃分帶來了挑戰(zhàn)。因此，采用現(xiàn)有閾值分割方法難以實現(xiàn)對油砂灰度圖像中3種介質(zhì)的劃分。

為了解決上述問題，考慮油砂灰度圖像特點，提出了適用于油砂的3 種介質(zhì)圖像分割方法。流程如圖2 所示，首先區(qū)分出巖石骨架和孔隙（一次劃分），并將計算獲得的孔隙度與實測值進行對比，作為閾值依據(jù)（一次閾值調(diào)整）；再從巖石骨架中區(qū)分出顆粒和瀝青質(zhì)（二次劃分），并將計算獲得的瀝青質(zhì)含量與實測瀝青質(zhì)含量進行對比（二次閾值調(diào)整）。閾值可由非局部均值濾波過程來確定，進而獲得具有不同峰值的曲線，將通過濾波獲得的孔隙度與巖心實際測試孔隙度進行對比，若與實測值一致，則將該閾值作為孔隙介質(zhì)灰度閾值的劃分標(biāo)準，由此可分割出孔隙與骨架。同理將實測瀝青質(zhì)含量作為標(biāo)準值來確定顆粒和瀝青質(zhì)的劃分閾值，可將骨架進一步分割，最終獲得孔隙、顆粒和瀝青質(zhì)3種介質(zhì)的分割圖像。

圖2 油砂3種介質(zhì)圖像分割方法Fig.2.Three-medium image segmenting method for oil sands

在閾值劃分的基礎(chǔ)上，結(jié)合CT 圖像觀察油砂細觀結(jié)構(gòu)，可見油砂顆粒與瀝青質(zhì)接觸關(guān)系主要包括3 種類型。第一種是接觸，結(jié)構(gòu)單元中只有顆粒相互接觸，不同粒徑的顆粒聚集在一起，顆粒間以點接觸為主，局部顆粒呈線接觸。第二種是膠結(jié)，結(jié)構(gòu)單元中同時具有顆粒和瀝青質(zhì)，但瀝青質(zhì)較少，顆粒之間由少量瀝青質(zhì)充填膠結(jié)，瀝青質(zhì)呈薄膜狀黏附在顆粒周圍，此時顆粒之間呈點接觸或線接觸。第三種是懸浮，顆粒由瀝青質(zhì)完全包裹，俗稱“油包砂”，顆粒間互不接觸，呈懸浮狀態(tài)，顆粒間完全靠瀝青質(zhì)支撐。

2.2 三維孔隙結(jié)構(gòu)重構(gòu)

基于現(xiàn)有的CT 圖像，可以運用直接立體成像的方法實現(xiàn)油砂的三維重構(gòu)。該方法所獲取的是無損圖像，只要準確劃分出不同介質(zhì)，便能較為準確地表征孔隙空間。圖3 展示了重構(gòu)后的油砂三維孔隙空間和骨架結(jié)構(gòu)，這些實體結(jié)構(gòu)為后期開展孔隙級流動模擬提供了幾何模型，通過建立孔隙網(wǎng)絡(luò)模型，并通過求解Navier-Stokes 方程的方式來研究流體在油砂中的流動。

圖3 油砂三維孔隙空間和骨架結(jié)構(gòu)重構(gòu)Fig.3.Reconstruction of three-dimensional pores and framework of oil sands

3 細觀尺度油砂滲流模擬

3.1 孔隙網(wǎng)絡(luò)建模

由于孔隙空間較為復(fù)雜，為開展孔隙尺度流動模擬，需要建立孔隙網(wǎng)絡(luò)模型，進而獲得孔隙和喉道的大小、數(shù)量等參數(shù)。中軸線算法的拓撲結(jié)構(gòu)準確，計算量小，適合進行復(fù)雜多重介質(zhì)樣品的研究。本文以中軸線算法作為孔隙網(wǎng)絡(luò)提取方法，主要步驟包括孔隙空間骨架提取、孔隙和喉道的識別。為了獲得能代表孔隙空間特征的簡化的球棍結(jié)構(gòu)，采用距離排序同倫細化算法來提取孔隙空間的中軸線骨架，通過分析骨架結(jié)構(gòu)來分割孔隙空間，進而識別孔隙和喉道。

3.2 細觀滲流模擬

3.2.1 控制方程及其求解

采用流體力學(xué)方法求解孔隙中流體流動的本質(zhì)為求解Navier-Stokes 方程，即將孔隙中流體的流動看作自由流，可以通過如下控制方程描述[21]。

動量守恒方程：

連續(xù)性方程（質(zhì)量守恒方程）：

多孔介質(zhì)中的流體流動非常緩慢，雷諾數(shù)非常小，此時的流動可視為蠕動流，即求解Navier-Stokes方程，整合后的控制方程為[22]

設(shè)定初始值和邊界條件后便可預(yù)測流體流動速度和壓力分布，通過在入口或出口邊界處對流動速度進行積分，便可獲得流過孔隙空間的流體流量，結(jié)合Darcy公式，便可獲得宏觀滲流參數(shù)，Darcy公式如下：

在用數(shù)值方法求解上式時，由于需要對相交網(wǎng)格進行離散以便更好地處理無滑移邊界條件，因此采用有限體積法對控制方程進行求解。與傳統(tǒng)劃分網(wǎng)格的方式不同，本文采用三維數(shù)字巖心的代表性體積單元作為網(wǎng)格，壓力未知數(shù)位于單元中心，速度未知數(shù)位于單元表面，并按照笛卡爾坐標(biāo)系進行分解。

3.2.2 代表性體積單元分析

雖然在裁剪原始油砂數(shù)字圖像時已經(jīng)盡可能保留巖心結(jié)構(gòu)信息，但基于觀察結(jié)構(gòu)特征定義代表單元的方法是定性的，要保證所建立的模型能夠表征整個樣品特征，需要進行代表性體積單元分析。本文采用的方法是以孔隙度作為約束，即通過對比分析隨著所選立方體邊長的增大，三維數(shù)字巖心中孔隙度的變化來實現(xiàn)。

用于CT的樣品是直徑為10.0 mm，高度為6.5 mm的圓柱體油砂，通過數(shù)字巖心重構(gòu)技術(shù)和提取技術(shù)，以該樣品圓柱中心點為中心，依次構(gòu)建邊長為0.5 mm，1.0 mm，1.5 mm，……，5.0 mm的立方體（圖4），對比不同尺寸下的孔隙度和滲透率。

圖4 代表性體積單元選取示意Fig.4.Schematic diagram of selection of a representative volume unit

應(yīng)用前文的方法，對不同尺寸樣品進行三維重構(gòu)及孔隙網(wǎng)絡(luò)建模。由于油砂中的瀝青質(zhì)溫度敏感性強，具有兩相性，需要分別分析2 種狀態(tài)下的孔隙網(wǎng)絡(luò)模型（圖5）。在常溫（25 ℃）下瀝青質(zhì)不能流動，被看作油砂骨架的一部分，此時模擬的是流體的有效滲流；當(dāng)溫度大于150 ℃時，瀝青質(zhì)成為可以流動的液態(tài)，此時油砂骨架部分只有骨架顆粒，模型計算的是油砂的絕對滲透率。這2 種孔隙網(wǎng)絡(luò)模型，分別對應(yīng)超稠油儲集層快速啟動階段和SAGD 生產(chǎn)階段的儲集層油砂滲流狀態(tài)。

圖5 研究區(qū)不同尺寸代表性體積單元的孔隙網(wǎng)絡(luò)模型Fig.5.Pore network models of representative volume units with different sizes

對于孔隙而言，圓球大小表征體積大小，顏色表征半徑（紅色為大孔隙，藍色為小孔隙）；對于喉道而言，顏色表征半徑，管徑大小代表喉道粗細。有效滲流條件下，樣品A1—A3 的孔隙結(jié)構(gòu)差異明顯；樣品A5 未得到連通的孔隙，無法形成孔隙空間網(wǎng)絡(luò)結(jié)構(gòu)；樣品A4、A6—A10 的孔隙結(jié)構(gòu)較為相似，可以表征原有樣品的孔隙結(jié)構(gòu)（圖6a）。絕對滲流條件下，樣品A1—A4 的孔隙結(jié)構(gòu)有明顯差異，樣品A5—A10 的孔隙結(jié)構(gòu)較為相似，樣品A5—A10基本可以表征原有樣品的孔隙結(jié)構(gòu)（圖6b）。

圖6 研究區(qū)不同尺寸代表性體積單元的孔隙度和滲透率Fig.6.Porosities and permeabilities of representative volume units with different sizes

為對油砂樣品的滲流特征開展定量分析，分別計算不同滲流條件下樣品的孔隙度和滲透率，并進行對比。樣品A1—A3 的滲透率與樣品A4—A10 的滲透率差距較大，因此不能用樣品A1—A3 代表原有樣品的孔隙度。樣品A4—A10滲透率差距不大，與迂曲度和孔隙網(wǎng)絡(luò)參數(shù)有關(guān)，從樣品A4 到A10，孔隙度按照一定規(guī)律減小，有效孔隙度減小的速度更快（圖6a）。樣品尺寸隨著樣品編號增大逐漸增大，絕對滲流條件下，有效孔隙度與孔隙度基本一致，說明樣品的孔隙連通性很好，孔隙度隨著樣品尺寸的增大而減小，滲透率隨著樣品尺寸的增大而減小，與孔隙度呈正相關(guān)，樣品A4、A5 孔隙度和滲透率的變化較大，樣品A6 是可以替代原有樣品的最小尺寸樣品（圖6b）。綜上所述，為準確表征風(fēng)城油田侏羅系齊古組油砂特征，可以選取的代表性體積單元最小應(yīng)為3.0 mm×3.0 mm×3.0 mm 的樣品。

根據(jù)代表性體積單元選取原則，為了更加細致地研究油砂孔滲參數(shù)與滲流特征，對標(biāo)準巖心柱掃描圖像進行截取，獲得直徑為1.5 cm，高度為3.0 cm 的數(shù)字巖心，在該巖心中心截取6.0 mm×3.0 mm×9.0 mm的區(qū)域作為研究對象W，在W 的4個角和中心位置分別截取3.0 mm×3.0 mm×3.0 mm 的立方體樣品，并依次編號為C1、C2、C3、C4 和C5（圖7）。分別計算5 個子樣品的孔隙度和x、y和z方向的絕對滲透率。

圖7 油砂子樣品選取示意Fig.7.Schematic diagram of selection of oil sand subsamples

（1）孔滲相對誤差分析 經(jīng)測定，樣品W 的孔隙度為9.52%，子樣品C1、C2、C3、C4 和C5 的孔隙度分別為11.53%、12.34%、14.49%、6.26%和8.53%，子樣品的孔隙度與樣品W 整體孔隙度的相對誤差可以用以下公式計算：

計算結(jié)果為，子樣品的孔隙度與樣品W 整體孔隙度的相對誤差分別為21.11%、29.62%、52.21%、34.24%和10.40%。

計算得到子樣品3個方向的滲透率，分別用Kx、Ky和Kz表示，計算子樣品各方向滲透率與樣品W整體滲透率的相對誤差，公式如下：

經(jīng)計算可知，子樣品孔隙度與樣品W 整體孔隙度相對誤差為10.40%～52.21%。而子樣品各方向滲透率與樣品W 整體滲透率的相對誤差為0.44%～119.11%（表1），說明了油砂樣品非均質(zhì)性較強。

表1 研究區(qū)油砂子樣品滲透率計算結(jié)果Table 1.Calculated permeability of subsamples of oil sands

（2）滲透率與孔隙度的關(guān)系 對比油砂5 個子樣品3 個方向的滲透率與孔隙度可知，滲透率與孔隙度基本呈正相關(guān)，即孔隙度越大，滲透率越大。但在個別異常點二者也存在負相關(guān)的現(xiàn)象，原因在于該點上孔隙連通性差，導(dǎo)致滲透率隨著孔隙度的增加而降低（圖8）。

圖8 研究區(qū)油砂子樣品不同方向滲透率與孔隙度的關(guān)系Fig.8.Relationships between permeability and porosity in different directions

（3）流場和壓力場特點 圖9 為子樣品C1、C3 和C5 的流場和壓力場分布，不同樣品的流場和壓力場分布差異較大。x、y和z方向上的流場流動速度為C3＞C1＞C5，這與滲透率的大小一致，而壓力場的變化也符合滲透率的變化規(guī)律，說明流場和壓力場的分布可以更為形象直觀地體現(xiàn)油砂的非均質(zhì)性。

圖9 研究區(qū)油砂子樣品流場和壓力場分布Fig.9.Distribution of flow field and pressure field of oil sand subsamples in the study area

（4）滲透率各向異性 通過計算得出5 個子樣品不同方向滲透率的比值，z方向滲透率與x方向滲透率比值為0.83～1.33，z方向滲透率與y方向滲透率比值為0.53～1.13，x方向滲透率與y方向滲透率比值為為0.64～1.18?？傮w而言，油砂各個方向上的滲透率差異較大，這與油砂本身的非均質(zhì)性和其成巖過程有關(guān)，說明油砂具有一定程度的各向異性。

4 結(jié)論

（1）油砂溫度、應(yīng)力敏感性強，通過直接自然斷面圖像研究油砂結(jié)構(gòu)和進行滲透率參數(shù)測試均易破壞其本身結(jié)構(gòu)，采用無損、非接觸的CT 重構(gòu)的數(shù)字巖心技術(shù)可有效獲取未受擾動的油砂內(nèi)部結(jié)構(gòu)。

（2）風(fēng)城油田齊古組油砂中顆粒和瀝青質(zhì)存在接觸、膠結(jié)和懸浮3 種結(jié)構(gòu)關(guān)系，油砂圖像中顆粒、瀝青質(zhì)和孔隙3 種介質(zhì)的分割效果決定數(shù)字巖心三維重構(gòu)模型的準確程度。

（3）為了選取能代表風(fēng)城油田齊古組油砂復(fù)雜結(jié)構(gòu)的數(shù)字巖心，代表性體積單元應(yīng)至少是邊長為3.0 mm的立方體。

（4）風(fēng)城油田齊古組油砂具有較強的非均質(zhì)性和各向異性，瀝青質(zhì)的相態(tài)（成為骨架一部分或者成為孔隙流體）顯著影響油砂的孔隙度、孔隙迂曲度和滲透率，不包含瀝青質(zhì)的骨架模型滲透率比包含瀝青質(zhì)的骨架模型滲透率高約2個數(shù)量級。

符號注釋

Ex、Ey、Ez——分別為x方向、y方向和z方向滲透率相對誤差；

f——相對誤差函數(shù)；

F——流體受到的外力，Pa；

I——單位張量；

K——絕對滲透率，m2；

Kx、Ky、Kz——分別為x方向、y方向和z方向滲透率，mD；

L——流動方向樣品長度，m；

p——流體壓力，Pa；

△p——流體壓差，Pa；

▽p——流體壓力梯度，Pa/m

Q——通過樣品的流量，m3/s；

S——流體通過多孔介質(zhì)的截面積，m2；

t——時間，s；

T——轉(zhuǎn)置符號；

u——流體流動速度，m/s；

▽u——流體流動速度梯度，s-1；

μ——流體動力黏度，Pa·s；

▽——哈密頓算子；

ρ——流體密度，kg/m3；

φi——子樣品i的孔隙度，%；

φz——樣品整體孔隙度，%。