本發(fā)明涉及油氣田開發(fā)工程數(shù)值模擬的技術(shù)領(lǐng)域，特別是一種頁巖基質(zhì)儲層孔隙空間表征方法。

背景技術(shù)：

為實現(xiàn)頁巖氣藏的高效開發(fā)，研究氣體在微納米孔隙中運移機制是關(guān)鍵力學(xué)問題之一，而數(shù)字巖心是進行微納尺度流動模擬的基礎(chǔ)。因此，進行頁巖微納米孔隙的識別和定量描述，建立準確的頁巖儲層孔隙空間表征方法成為研究的重點和難點。

頁巖孔隙結(jié)構(gòu)復(fù)雜，既有納米級的有機質(zhì)孔隙、還發(fā)育有無機礦物中納米-微米級孔隙，以及微米級-毫米級天然裂縫。泥頁巖熱演化程度不同，生成的孔隙類型不同。有機孔隙與無機孔隙中流體的相態(tài)具有較大區(qū)別，目前學(xué)者普遍認為頁巖有機質(zhì)中的有機孔隙是主要的存儲空間，因此研究有機孔隙的比例和分布等對頁巖氣后續(xù)開發(fā)具有重要的理論和實踐價值。由于頁巖孔隙極小，目前常用來研究頁巖孔隙結(jié)構(gòu)的方法主要有，掃描電鏡(SEM)、聚焦離子束掃描電鏡(FIB-SEM)、原子力顯微鏡(AFM)、納米CT、透射電子顯微鏡等掃描手段，結(jié)合能譜(ESD)或背散射圖像(BEI)還可以實現(xiàn)不同礦物成分的三維分布圖像。

基于上述設(shè)備，可以得到數(shù)字巖心，構(gòu)建方法分為兩大類，聚焦離子束掃描電鏡、納米CT等物理構(gòu)建法可以直接給出三維真實數(shù)字巖心，而掃描電鏡、原子力顯微鏡、透射電子顯微鏡等只能給出二維圖像資料，需要利用隨機算法重構(gòu)數(shù)字巖心。然而，不管是物理法直接構(gòu)建還是利用二維圖像資料隨機重構(gòu)，只能構(gòu)建單一尺度或單一成份的數(shù)字巖心，這對于非均質(zhì)性很強的頁巖來說都不夠準確，有學(xué)者 建議進行工業(yè)CT-微米CT-納米CT/FIB系列掃描分析的表征手段。但是即便如此，在同一尺度下，也很難同時考慮有機孔隙和無機孔隙。

技術(shù)實現(xiàn)要素：

本發(fā)明的目的在于提供一種頁巖基質(zhì)儲層孔隙空間表征方法，以滿足頁巖基質(zhì)儲層孔隙結(jié)構(gòu)分析及后續(xù)流動模擬的需求。

為實現(xiàn)上述技術(shù)目的，達到上述技術(shù)效果，本發(fā)明公開了一種頁巖基質(zhì)儲層孔隙空間表征方法，表征方法包括了巖心二維圖像采集、模擬退火法建立無機數(shù)字巖心、模擬馬爾可夫鏈蒙特卡洛法(MCMC)構(gòu)建有機孔隙數(shù)字巖心、疊加建立有機、無機孔隙數(shù)字巖心，具體的操作方法如下：

(1)巖心二維圖像采集：基于掃描電子顯微鏡(SEM)，獲取真實頁巖巖心在不同平面上的低分辨率無機孔隙和高分辨率有機孔隙掃描圖像，并將圖像分割轉(zhuǎn)換得到頁巖巖心二值圖像；

(2)模擬退火法建立無機數(shù)字巖心：在重建數(shù)字巖心過程中，除了保證孔隙度不變，還選取了兩個統(tǒng)計函數(shù)作為建模過程中的待擬合量，即：自相關(guān)函數(shù)、線性路徑函數(shù)；

a單點概率函數(shù)：

假設(shè)多相系統(tǒng)中第j相所占區(qū)域為v_j，其在整個系統(tǒng)中的體積分數(shù)為則定義第j相的相函數(shù)如下：

式中，表示系統(tǒng)中的任意一點，對于僅考慮孔隙和巖石骨架的兩相系統(tǒng)，可將相函數(shù)簡化為：

于是，系統(tǒng)的單點概率函數(shù)(即孔隙度φ)可由統(tǒng)計平均值給出：

b自相關(guān)函數(shù)：

多相系統(tǒng)中第j相的自相關(guān)函數(shù)定義為：

式中，為系統(tǒng)中相距一定距離r的任意兩點，為計算方便，令r的單位為像素數(shù)，像素數(shù)與像素邊長的乘積即物理長度；對于三維情況，r的單位為體素數(shù)；

c線性路徑函數(shù)：

線性路徑函數(shù)是描述多孔介質(zhì)微觀結(jié)構(gòu)的一個重要函數(shù)，定義如下：

式中，為系統(tǒng)中相距一定距離r的任意兩點，為連接線段上的任意點，與自相關(guān)函數(shù)相同，在孔隙空間重建中以孔隙相為研究對象，L^(j)(r)可進一步簡化為L(r)；

利用頁巖巖心二值圖像，使用模擬退火法構(gòu)建數(shù)字巖心具體步驟如下：

a.初始時，在指定大小的三維空間內(nèi)隨機產(chǎn)生由兩種不同像素點組成的巖石孔隙和骨架體系，由這兩種體系組成的空間結(jié)構(gòu)的孔隙度要與由薄片分析得來的孔隙度相同；

b.以像素點為基礎(chǔ)，利用模擬退火算法對兩種體系的像素點進行演化，即通過交換兩種體系之間的像素點來改變整個系統(tǒng)的狀態(tài)，系統(tǒng)能否進行更新取決于新系統(tǒng)是否滿足Metropolis準則：

p≥RAN(1,0)

式中，p為在第k+1次搜索時狀態(tài)的接受概率，ΔE＝E^k+1-E^k，T為引進的控制參數(shù)，稱之為“溫度”，RAN(1,0)表示在[0.0，1.0]內(nèi)選取的隨機數(shù)；在演化過程中孔隙度保持恒定不變，擬合的目標函數(shù)為與孔隙空間結(jié)構(gòu)有關(guān)的自相關(guān)函數(shù)和線性路徑函；

(3)模擬馬爾可夫鏈蒙特卡洛法(MCMC)構(gòu)建有機孔隙數(shù)字巖心：

選擇雙體素組合方法，基于巖心二值圖像(狀態(tài)僅僅為0或者1)，利用xy、yz、xz三個平面來同時來構(gòu)建三維馬爾可夫鏈模型，具體模擬步驟如下：

a.模擬一維鏈開始的體素，采用薄片的孔隙度作為第一個體素狀態(tài)為孔隙的條件概率；

b.沿y方向模擬第一層第一行上的體素，第一行第二個體素使用2鄰居進行模擬，從第三個體素開始均使用3鄰居進行模擬，其條件概率通過xy平面的二維巖心薄片進行推導(dǎo)得出；

c.沿x方向模擬第一層各行體素進而構(gòu)建該層所有體素，對于雙體素(i,j)和(i,j+1)，在邊緣體素使用3、4鄰居進行模擬，內(nèi)部體素使用5、6鄰居進行模擬，其條件概率通過xy平面的二維巖心薄片進行推導(dǎo)得出；

d.沿z方向模擬各層體素進而構(gòu)建三維模型，對于雙體素(i,j,k)和(i,j+1,k)，第二層的第一行跟第一層的第二行構(gòu)造方法相似，在邊緣位置體素使用3、4鄰居進行模擬，內(nèi)部體素使用5、6鄰居進行模擬，不同的是其條件概率由其他方向二維薄片的條件概率依次重構(gòu)，從第二層第二行開始，邊緣體素使用9、10鄰居進行模擬，內(nèi)部體素14、15鄰居進行模擬，依次重構(gòu)最后得到基于馬爾科夫鏈蒙特卡洛方法的三維數(shù)字巖心，其中三維15鄰居可表示為：

(4)疊加建立有機、無機孔隙數(shù)字巖心：

三維數(shù)字巖心的數(shù)據(jù)體是用二進制0和1來進行表示的，基于物理尺寸相同的無機孔隙數(shù)字巖心和有機孔隙數(shù)字巖心，疊加數(shù)字巖心構(gòu)建步驟如下：

a.將低分辨數(shù)字巖心中每一個體素細化為i*i*i個更小體素，細 化后的無機孔隙數(shù)字巖心和有機孔隙數(shù)字巖心的物理尺寸和體素尺寸都完全相同，其中，i為無機孔隙數(shù)字巖心和有機孔隙數(shù)字巖心的分辨率之比；

b.為了構(gòu)建同時描述無機孔隙和有機孔隙特征的疊加數(shù)字巖心，按照布爾疊加算法，對二進制數(shù)據(jù)體進行操作，得到基于不同分辨率的無機孔隙數(shù)字巖心和有機孔隙數(shù)字巖心構(gòu)建的頁巖基質(zhì)孔隙數(shù)字巖心；

(5)氮氣吸附法:

通過氮氣吸附實驗可以獲得頁巖的吸附和解吸等溫線，并進一步通過計算可以得到孔隙結(jié)構(gòu)參數(shù)。通過氮氣吸附實驗獲得實際巖心的孔隙結(jié)構(gòu)參數(shù)，并與疊加的數(shù)字巖心對比檢驗?zāi)Ｐ偷臏蚀_性。

其中,掃描電子顯微鏡(SEM)獲得的圖像二值化方法為采用Otsu法、迭代法、手動法、最大熵法、自適應(yīng)閥值法、基本全局閾值法中的一種，優(yōu)選的為Otsu法。

其中,利用二維圖像資料重構(gòu)數(shù)字巖心采用高斯模擬法、模擬退火法、過程模擬法、多點統(tǒng)計法、馬爾可夫隨機重建法，優(yōu)選用模擬退火法構(gòu)建無機孔隙數(shù)字巖心，馬爾可夫隨機重建法構(gòu)建有機孔隙數(shù)字巖心。

其中,構(gòu)建同時描述無機孔隙和有機孔隙特征的疊加數(shù)字巖心，布爾疊加算法具體表達式如下所示：

Ω＝Ω₁+Ω₂

式中，Ω表示頁巖基質(zhì)孔隙數(shù)字巖心,Ω₁表示頁巖無機孔隙數(shù)字巖心，Ω₂表示頁巖有機孔隙數(shù)字巖心，上標s、op和ip分別表示巖 心骨架、有機孔隙和無機孔隙，如果用0代表骨架，1代表有機孔隙空間，2代表無機孔隙空間，則相應(yīng)的布爾疊加算法操作為：0+0＝0,0+1＝1,1+0＝2,1+1＝2。

本發(fā)明具有以下有益效果：

1.本發(fā)明通過疊加算法建立了同時考慮有機孔隙和無機孔隙的數(shù)字巖心，并通過實驗對建立的數(shù)字巖心進行了幾何拓撲結(jié)構(gòu)分析，孔隙半徑分布和配位數(shù)與實驗結(jié)果擬合程度較好。

附圖說明

圖1為本發(fā)明的頁巖無機孔隙和有機孔隙代表性掃描電鏡圖像。

圖2為本發(fā)明的頁巖無機孔隙和有機孔隙SEM二值圖像(白色表示巖石骨架，黑色表示巖石孔隙)。

圖3為本發(fā)明的模擬退火法建立數(shù)字巖心的程序流程。

圖4為本發(fā)明的馬爾可夫鏈蒙特卡洛法的三維數(shù)字巖心構(gòu)造圖(白色實線表示正在遍歷的體素，白色虛線表示還未遍歷的體素，深色表示已經(jīng)遍歷的體素，淺色表示正在遍歷體素的其他鄰居)。

圖5為本發(fā)明的有機孔隙數(shù)字巖心、無機孔隙數(shù)字巖心和二者疊加后的數(shù)字巖心。

具體實施方式

為了使本發(fā)明的目的、技術(shù)方案及優(yōu)點更加清楚明白，以下結(jié)合附圖及實施例，對本發(fā)明進行進一步詳細說明。

本發(fā)明公開了一種頁巖基質(zhì)儲層孔隙空間表征方法，表征方法包括了巖心二維圖像采集、模擬退火法建立無機數(shù)字巖心、模擬馬爾可夫鏈蒙特卡洛法(MCMC)構(gòu)建有機孔隙數(shù)字巖心、疊加建立有機、無機孔隙數(shù)字巖心，具體的操作方法如下：

巖心二維圖像采集：基于掃描電子顯微鏡(SEM)，獲取真實頁巖巖心在不同平面上的低分辨率無機孔隙和高分辨率有機孔隙掃描圖 像，掃描電子顯微鏡獲得的圖像二值化方法為采用Otsu法；

模擬退火法建立無機數(shù)字巖心：在重建數(shù)字巖心過程中，除了保證孔隙度不變，還選取了兩個統(tǒng)計函數(shù)作為建模過程中的待擬合量，即：自相關(guān)函數(shù)、線性路徑函數(shù)；

a單點概率函數(shù)：

假設(shè)多相系統(tǒng)中第j相所占區(qū)域為ν_j，其在整個系統(tǒng)中的體積分數(shù)為則定義第j相的相函數(shù)如下：

式中，表示系統(tǒng)中的任意一點，對于僅考慮孔隙和巖石骨架的兩相系統(tǒng)，可將相函數(shù)簡化為：

于是，系統(tǒng)的單點概率函數(shù)(即孔隙度φ)可由統(tǒng)計平均值給出：

b自相關(guān)函數(shù)：

多相系統(tǒng)中第j相的自相關(guān)函數(shù)定義為：

式中，為系統(tǒng)中相距一定距離r的任意兩點，為計算方便，令r的單位為像素數(shù)，像素數(shù)與像素邊長的乘積即物理長度；對于三維情況，r的單位為體素數(shù)；

c線性路徑函數(shù)：

線性路徑函數(shù)是描述多孔介質(zhì)微觀結(jié)構(gòu)的一個重要函數(shù)，定義如下：

式中，為系統(tǒng)中相距一定距離r的任意兩點，為連接線段上的任意點，與自相關(guān)函數(shù)相同，在孔隙空間重建中以孔隙相為研究對象，L^(j)(r)可進一步簡化為L(r)；

利用頁巖巖心二值圖像，使用模擬退火法構(gòu)建數(shù)字巖心具體步驟如下：

a.初始時，在指定大小的三維空間內(nèi)隨機產(chǎn)生由兩種不同像素點組成的巖石孔隙和骨架體系，由這兩種體系組成的空間結(jié)構(gòu)的孔隙度要與由薄片分析得來的孔隙度相同；

b.以像素點為基礎(chǔ)，利用模擬退火算法對兩種體系的像素點進行演化，即通過交換兩種體系之間的像素點來改變整個系統(tǒng)的狀態(tài)，系統(tǒng)能否進行更新取決于新系統(tǒng)是否滿足Metropolis準則：

p≥RAN(1,0)

式中，p為在第k+1次搜索時狀態(tài)的接受概率，ΔE＝E^k+1-E^k，T為引進的控制參數(shù)，稱之為“溫度”，RAN(1,0)表示在[0.0，1.0]內(nèi)選取的隨機數(shù)；在演化過程中孔隙度保持恒定不變，擬合的目標函數(shù)為與孔隙空間結(jié)構(gòu)有關(guān)的自相關(guān)函數(shù)和線性路徑函；

模擬馬爾可夫鏈蒙特卡洛法(MCMC)構(gòu)建有機孔隙數(shù)字巖心：

選擇雙體素組合方法，基于巖心二值圖像(狀態(tài)僅僅為0或者1)，利用xy、yz、xz三個平面來同時來構(gòu)建三維馬爾可夫鏈模型，具體模擬步驟如下：

a.模擬一維鏈開始的體素，采用薄片的孔隙度作為第一個體素狀態(tài)為孔隙的條件概率；

b.沿y方向模擬第一層第一行上的體素，第一行第二個體素使用2鄰居進行模擬，從第三個體素開始均使用3鄰居進行模擬，其條件概率通過xy平面的二維巖心薄片進行推導(dǎo)得出；

c.沿x方向模擬第一層各行體素進而構(gòu)建該層所有體素，對于雙體素(i,j)和(i,j+1)，在邊緣體素使用3、4鄰居進行模擬，內(nèi)部體素使用5、6鄰居進行模擬，其條件概率通過xy平面的二維巖心薄片進行推導(dǎo)得出；

d.沿z方向模擬各層體素進而構(gòu)建三維模型，對于雙體素(i,j,k)和(i,j+1,k)，第二層的第一行跟第一層的第二行構(gòu)造方法相似，在邊緣位置體素使用3、4鄰居進行模擬，內(nèi)部體素使用5、6鄰居進行模擬，不同的是其條件概率由其他方向二維薄片的條件概率依次重構(gòu)，從第二層第二行開始，邊緣體素使用9、10鄰居進行模擬，內(nèi)部體素14、15鄰居進行模擬，依次重構(gòu)最后得到基于馬爾科夫鏈蒙特卡洛方法的三維數(shù)字巖心，其中三維15鄰居可表示為：

疊加建立有機、無機孔隙數(shù)字巖心：

三維數(shù)字巖心的數(shù)據(jù)體是用二進制0和1來進行表示的，基于物理尺寸相同的無機孔隙數(shù)字巖心和有機孔隙數(shù)字巖心，疊加數(shù)字巖心構(gòu)建步驟如下：

a.將低分辨數(shù)字巖心中每一個體素細化為i*i*i個更小體素，細化后的無機孔隙數(shù)字巖心和有機孔隙數(shù)字巖心的物理尺寸和體素尺寸都完全相同，其中，i為無機孔隙數(shù)字巖心和有機孔隙數(shù)字巖心的分辨率之比；

b.為了構(gòu)建同時描述無機孔隙和有機孔隙特征的疊加數(shù)字巖心，按照布爾疊加算法，對二進制數(shù)據(jù)體進行操作，得到基于不同分辨率的無機孔隙數(shù)字巖心和有機孔隙數(shù)字巖心構(gòu)建的頁巖基質(zhì)孔隙數(shù)字巖心，具體表達式如下所示：

Ω＝Ω₁+Ω₂；

式中，Ω表示頁巖基質(zhì)孔隙數(shù)字巖心，Ω₁表示頁巖無機孔 隙數(shù)字巖心，Ω₂表示頁巖有機孔隙數(shù)字巖心，上標s、op和ip分別表示巖心骨架、有機孔隙和無機孔隙，如果用0代表骨架，1代表有機孔隙空間，2代表無機孔隙空間，則相應(yīng)的布爾疊加算法操作為：0+0＝0,0+1＝1,1+0＝2,1+1＝2。

實施例2

實驗?zāi)康募胺椒ǎ簽榱嗽u價本發(fā)明的可行性和有效性，本實施例以某頁巖樣品為例進行研究實例，并且對該頁巖樣品進行礦石成分分析、X射線衍射圖譜分析、氮氣吸脫附分析，具體的實驗步驟參照標準實驗操作。

實驗結(jié)果：

對某頁巖樣品進行礦石分析后得到該頁巖樣品是一套碎屑巖和碳酸鹽巖混合沉積，是一套區(qū)域性優(yōu)質(zhì)烴源巖。由X射線衍射圖譜分析可知該樣品礦物成份以石英和黃鐵礦為主，含有少量的斜長石和黏土礦物，有機碳含量(TOC)6.5％，孔隙度為5.3％，并通過氮氣吸脫附實驗得到孔隙結(jié)構(gòu)參數(shù)為：孔隙半徑大體分布在2-120納米，主要分布在2-20納米。

基于疊加的數(shù)字巖心，頁巖有機孔隙半徑基本分布在3-25納米，頁巖無機孔隙半徑基本分布在10-120納米。疊加后的數(shù)字巖心可以得到類似的孔隙半徑分布，大體分布在3-120納米，主要分布在3-20納米，與氮氣吸附法獲得的孔徑分布結(jié)果基本吻合，驗證了本發(fā)明的準確性。

以上所述，僅為本發(fā)明較佳的具體實施方式，但本發(fā)明的保護范圍并不局限于此，任何熟悉本技術(shù)領(lǐng)域的技術(shù)人員在本發(fā)明揭露的技術(shù)范圍內(nèi)，可輕易想到的變化或替換，都應(yīng)涵蓋在本發(fā)明的保護范圍之內(nèi)。