基于多點(diǎn)地質(zhì)統(tǒng)計(jì)學(xué)的數(shù)字巖心建模

吳玉其, 林承焰, 任麗華, 閆偉超,王 楊, 陳仕臻, 由春梅, 張 麗

(1.中國石油大學(xué)(華東)地球科學(xué)與技術(shù)學(xué)院,山東青島 266580; 2.山東省油藏地質(zhì)重點(diǎn)實(shí)驗(yàn)室,山東青島 266580;3.大慶油田有限責(zé)任公司勘探開發(fā)研究院,黑龍江大慶 163712; 4.山東科技大學(xué)資源與土木工程系,山東泰安 271019)

對于石油和天然氣行業(yè)而言,數(shù)字巖心技術(shù)作為新興的巖石物理試驗(yàn)數(shù)值模擬方法,已經(jīng)廣泛地應(yīng)用于地質(zhì)、地震、測井和開發(fā)以及提高采收率等各個(gè)領(lǐng)域。與常規(guī)巖石物理試驗(yàn)相比,數(shù)字巖心技術(shù)能夠同時(shí)進(jìn)行成巖作用、固體力學(xué)、聲學(xué)、電學(xué)、流體力學(xué)及流固耦合等試驗(yàn)的模擬[1-3]。開展該方向的研究對石油和天然氣行業(yè)勘探和開發(fā)的意義不言而喻[4]。一個(gè)準(zhǔn)確的三維數(shù)字巖心模型是進(jìn)行數(shù)字巖心分析技術(shù)的前提和基礎(chǔ)。隨著試驗(yàn)儀器的創(chuàng)新和新理論的突破,國內(nèi)外學(xué)者不斷提出新的構(gòu)建數(shù)字巖心模型的方法。到目前為止,數(shù)字巖心建模方法包括三大類,即物理試驗(yàn)法、數(shù)值重建法和混合法[5]。物理試驗(yàn)法指利用試驗(yàn)儀器(如高倍光學(xué)顯微鏡或X射線CT掃描儀等)對巖心樣品拍攝或掃描以獲取大量的巖心二維圖片,通過建模程序或軟件把二維圖片重構(gòu)成三維數(shù)字巖心的方法,主要有序列切片成像法[6]、激光掃描共聚焦顯微鏡法[7]和X射線CT掃描法[4]等。數(shù)值重建法是指以少量的二維薄片圖像為基礎(chǔ),利用二維圖片中包含的信息,通過隨機(jī)模擬法或沉積巖形成過程模擬法來重建三維數(shù)字巖心的方法,包括隨機(jī)模擬法和過程法。隨機(jī)模擬法大多是基于地質(zhì)統(tǒng)計(jì)學(xué)方法,利用少量的二維薄片信息模擬生成數(shù)字巖心模型，包括截?cái)喔咚闺S機(jī)域法[8]、模擬退火法[9-10]、馬爾科夫鏈蒙特卡洛法[11]、序貫指示模擬法[12]、多點(diǎn)統(tǒng)計(jì)學(xué)法[13]、相位恢復(fù)法[14-15]等。過程法是基于模擬沉積巖形成過程(沉積和成巖過程)的一種建模方法[16-17]?；旌戏ㄊ侵笖?shù)值重建法中兩種方法的混合,如Politis等[18]、劉學(xué)鋒和孫建孟[19]結(jié)合模擬退火法和過程法重建了三維數(shù)字巖心。上述方法中,目前常用的有3種方法:X射線CT掃描法、模擬退火法和過程法。雖然X射線CT掃描法是比較準(zhǔn)確的建模方法,但是價(jià)格貴試驗(yàn)周期長,而且不能解決樣品的分辨率和尺寸相互矛盾的問題。相比X射線CT掃描法,模擬退火法和過程法成本低效率高。不過模擬退火法建模是基于傳統(tǒng)的兩點(diǎn)地質(zhì)統(tǒng)計(jì)學(xué)理論,在傳統(tǒng)地質(zhì)統(tǒng)計(jì)學(xué)中計(jì)算兩點(diǎn)地質(zhì)變量空間之間相關(guān)性的重要工具是變差函數(shù),該函數(shù)一般只能表征空間上兩點(diǎn)之間的相關(guān)性,很難精確地模擬具有復(fù)雜幾何形態(tài)和拓?fù)浣Y(jié)構(gòu)的孔隙空間,這往往會(huì)使重構(gòu)出的孔隙結(jié)構(gòu)缺少長距離連通性[20-21]。過程法模擬沉積巖形成過程時(shí)假設(shè)顆粒全為球體,這與實(shí)際顆粒形態(tài)相差很大,而且成巖作用模擬過程比較簡單,不適于成巖作用復(fù)雜的巖石。為了解決上述問題,筆者依據(jù)多點(diǎn)地質(zhì)統(tǒng)計(jì)學(xué)法能夠復(fù)制地質(zhì)變量復(fù)雜空間結(jié)構(gòu)特征的優(yōu)點(diǎn),對多點(diǎn)地質(zhì)統(tǒng)計(jì)學(xué)法構(gòu)建三維數(shù)字巖心模型及重構(gòu)模型的準(zhǔn)確性評(píng)價(jià)進(jìn)行研究。

1 方法原理

1.1 多點(diǎn)地質(zhì)統(tǒng)計(jì)學(xué)

多點(diǎn)統(tǒng)計(jì)學(xué)(multiple-point statistics,MPS)最早被Deutsch提出[22]?，F(xiàn)今MPS已經(jīng)廣泛應(yīng)用于儲(chǔ)層建模中[23-26],部分建模理論已經(jīng)比較成熟。MPS分為迭代法和非迭代法,其中迭代法包括模擬退火法、基于吉布斯取樣迭代法和馬爾科夫鏈蒙特卡洛法。1993年,Guardiano和 Srivastava[27]提出了一種非迭代算法,直接從“密集訓(xùn)練圖像”提取局部條件概率,并利用序貫指示模擬方法進(jìn)行模擬實(shí)現(xiàn)。但該方法在計(jì)算過程中需要對每個(gè)數(shù)據(jù)事件進(jìn)行設(shè)置,而且每模擬一個(gè)節(jié)點(diǎn)需要重新掃描一遍訓(xùn)練圖像,這嚴(yán)重制約了模擬的效率,直到2000年Strebelle[27]提出了單一標(biāo)準(zhǔn)方程模擬(single normal equation simulation,SNESIM)算法,應(yīng)用搜索樹(search tree)一次性存儲(chǔ)訓(xùn)練圖像中所有數(shù)據(jù)事件的條件概率分布,并保證了在模擬過程中能快速地提取條件概率分布函數(shù),進(jìn)而大大降低了模擬所需時(shí)間。

2004年,Okabe 和 Blunt[13]把MPS法用于數(shù)字巖心建模中,他們以巖心的二維圖片為訓(xùn)練圖像,通過旋轉(zhuǎn)該圖片獲得x、y和z方向上的條件數(shù)據(jù),然后利用優(yōu)選的搜索模板掃描訓(xùn)練圖像,產(chǎn)生一個(gè)隨機(jī)的三維數(shù)字巖心模型。將模擬結(jié)果與模擬退火法的重構(gòu)模型對比,MPS有效地解決了重構(gòu)模型孔隙空間長距離連通性的問題。但是Okabe和Blunt建模的不足之處在于他們假定模型內(nèi)部孔隙結(jié)構(gòu)在水平方向上和垂直方向上為各向同性。

1.2 SNESIM算法

SNESIM算法是MPS方法中進(jìn)行離散型變量模擬最常用的方法。通過多重?cái)?shù)據(jù)模板掃描訓(xùn)練圖像來獲取地質(zhì)變量不同尺寸的結(jié)構(gòu)特征,將掃描得到的所有地質(zhì)模式存儲(chǔ)于搜索樹中;在多個(gè)已知節(jié)點(diǎn)為條件數(shù)據(jù)的約束下,求取未知節(jié)點(diǎn)的條件概率分布。再基于序貫?zāi)M樣式,模擬每個(gè)未知節(jié)點(diǎn)。已經(jīng)模擬節(jié)點(diǎn)的值會(huì)變成下一個(gè)待模擬節(jié)點(diǎn)的條件數(shù)據(jù)。該算法通過使用多重?cái)?shù)據(jù)模板和訓(xùn)練圖像能夠考慮多尺度和多節(jié)點(diǎn)之間的關(guān)聯(lián)性,更加適用于具有長距離連通性的數(shù)字巖心模型的重構(gòu)。為了便于理解該算法構(gòu)建三維數(shù)字巖心的步驟,下面對算法中的幾個(gè)重要概念進(jìn)行簡要介紹。

1.2.1 數(shù)據(jù)模板與數(shù)據(jù)事件

數(shù)據(jù)模板(data template)是指由中心節(jié)點(diǎn)u和n個(gè)向量{u+ha,a=1,2,…,n}組成的數(shù)據(jù)幾何體,通常也稱為搜索模板(search template),常用τn表示。數(shù)據(jù)模板中心處u的狀態(tài)是未知值S(u),其他網(wǎng)格節(jié)點(diǎn)處的狀態(tài)一般為已知值S(u+ha)=S(ua),a=1,2,…,n。相對應(yīng)地,以u(píng)為中心,尺寸大小為n的數(shù)據(jù)事件(data event)dn是由搜索模板τn和n個(gè)向量處的數(shù)據(jù)值{S(uα),α=1,2,…,n}共同構(gòu)成。

1.2.2 訓(xùn)練圖像

影響數(shù)字巖心重構(gòu)模型準(zhǔn)確性的關(guān)鍵因素之一在于訓(xùn)練圖像的選擇。訓(xùn)練圖像本質(zhì)上就是一個(gè)概念模型,它要盡可能地包含要模擬的地質(zhì)變量所有的組合模式,如進(jìn)行數(shù)字巖心建模時(shí)的孔隙結(jié)構(gòu)模式。理論上,選擇較大尺寸的訓(xùn)練圖像會(huì)使建模效果要更好,因?yàn)橛?xùn)練圖像尺寸越大,它所包含的孔喉組合模式也就越多;但由于受計(jì)算機(jī)性能的限制,尺寸較大的訓(xùn)練圖像會(huì)大大降低建模效率,因此選擇訓(xùn)練圖像時(shí)要綜合考慮多種因素。

對于數(shù)字巖心建模,訓(xùn)練圖像可以來源于X射線CT掃描的真實(shí)重構(gòu)模型,也可以是二維鏡下薄片,如鑄體薄片或者掃描電鏡圖像等。圖1(a)為二維訓(xùn)練圖像,它是由普通薄片二值化后得到的;其中黑色為顆粒,白色為孔隙。圖1(b)為用二維7×7的數(shù)據(jù)模板掃描訓(xùn)練圖像的過程。

1.2.3 條件概率分布函數(shù)

對于數(shù)字巖心中的某一節(jié)點(diǎn)或體素點(diǎn),其可能的狀態(tài)值為{sk,k=1,2},s1為顆粒,s2為孔隙。假設(shè)在一個(gè)隨機(jī)模式下,當(dāng)搜索模板掃描完訓(xùn)練圖像后,某一數(shù)據(jù)事件dn={S(uα)=skα,α=1,2,…,n}重復(fù)次數(shù)為c(dn),在此條件上,中心節(jié)點(diǎn)u取值{S(u)=sk,k=1,2}的重復(fù)次數(shù)為ck(dn),則中心節(jié)點(diǎn)取值S(u)為孔隙或顆粒的概率,用條件概率分布函數(shù)(conditional probability distribution function,CPDF)表示為

(1)

當(dāng)搜索模板掃描訓(xùn)練圖像后,得到的每個(gè)數(shù)據(jù)事件對應(yīng)的CPDF都會(huì)被存儲(chǔ)到一個(gè)樹狀數(shù)據(jù)結(jié)構(gòu)體(搜索樹)中。搜索樹的構(gòu)建保證了在模擬過程中搜索模板只需掃描一次訓(xùn)練圖像即可把所以數(shù)據(jù)事件記錄下來,這大大提高了建模的效率。

圖1 數(shù)據(jù)模板掃描訓(xùn)練圖像Fig.1 Search template scanning training image

1.2.4 多重網(wǎng)格

(2)

圖2 三重網(wǎng)格和三重搜索模板Fig.2 Three multiple grids and three multi-gridsearch templates

2 重構(gòu)數(shù)字巖心模型

2.1 數(shù)據(jù)預(yù)處理

以密西根盆地Waverly組的河流相Berea砂巖為例。Berea砂巖均質(zhì)性較好,顆粒以石英為主,包括少量的長石、白云石和黏土礦物,顆粒分選較好。Berea砂巖的數(shù)據(jù)通過X射線微米CT掃描儀采集[30],大小為4003體素,分辨率為5.345 μm,孔隙度為19.6%。從全巖樣中提取兩個(gè)1503體素的體積元為試驗(yàn)對象,令其中一個(gè)體積元為訓(xùn)練圖像(圖3(a)),孔隙度為20.10%;另一個(gè)體積元為真實(shí)模型(圖3(b)),孔隙度為19.95%,圖3中紅色為顆粒,藍(lán)色為孔隙。從真實(shí)模型中選取兩張正交的二維切片(圖3(c))作為已知信息,從中提取一部分像素點(diǎn)作為條件數(shù)據(jù)。

圖3 訓(xùn)練圖像、真實(shí)模型、條件數(shù)據(jù)和重構(gòu)模型Fig.3 Training image, real model, conditioning data and reconstructed models

2.2 建模步驟

以X射線CT掃描法建立的數(shù)字巖心為訓(xùn)練圖像,以兩張正交的二維切片的部分像素點(diǎn)為條件數(shù)據(jù),使用多點(diǎn)地質(zhì)統(tǒng)計(jì)學(xué)法重構(gòu)三維數(shù)字巖心的具體步驟為:

(1)選擇訓(xùn)練圖像。選擇X射線CT掃描獲得的三維圖像作為訓(xùn)練圖像(圖3(a)),三維訓(xùn)練圖像包含更加真實(shí)的立體的孔隙結(jié)構(gòu)模式。為了捕獲較大尺度的孔隙結(jié)構(gòu),設(shè)定搜索模板的搜索半徑為40×40×40,搜索角度為0°。

(2)使用四重搜索模板掃描訓(xùn)練圖像建立搜索樹。四重搜索模板不僅能保證捕獲不同尺度的孔隙結(jié)構(gòu)模式,而且能提高建模效率。

(3)將條件數(shù)據(jù)分配到相應(yīng)的網(wǎng)格,選定一條隨機(jī)路徑,訪問每一個(gè)待模擬節(jié)點(diǎn)。

(4)模擬未知節(jié)點(diǎn)u時(shí),保留那些在最大搜索模板范圍內(nèi)的條件數(shù)據(jù),假設(shè)條件數(shù)據(jù)的數(shù)量為n,相應(yīng)的數(shù)據(jù)事件為dn,在搜索樹中檢索數(shù)據(jù)事件dn的CPDF。如果檢索過程中數(shù)據(jù)事件dn的重復(fù)數(shù)小于設(shè)置的最小重復(fù)數(shù),就把搜索模板中最遠(yuǎn)的條件數(shù)據(jù)去掉,在條件數(shù)據(jù)數(shù)量為n-1的條件下,再去檢索數(shù)據(jù)事件dn-1的CPDF;如果數(shù)據(jù)事件的重復(fù)數(shù)仍小于設(shè)定的最小重復(fù)數(shù),繼續(xù)重復(fù)上述操作。倘若條件數(shù)據(jù)的數(shù)量一直減小到n=1時(shí)仍不滿足要求,就把目標(biāo)邊緣概率賦值給CPDF。

(5)在一定的CPDF下,選擇u的一個(gè)狀態(tài)值S(u)作為下一個(gè)待模擬節(jié)點(diǎn)的條件數(shù)據(jù)。

(6)沿著隨機(jī)路徑重復(fù)(4)和(5)步驟,直到模擬完所有未知節(jié)點(diǎn)。

訓(xùn)練圖像選擇時(shí)應(yīng)該優(yōu)先考慮三維圖像,因?yàn)閹r心中的孔隙空間實(shí)質(zhì)上是三維展布的,具有立體的空間特征,這類似于沉積相中的分流河道具有立體的空間幾何特征。二維訓(xùn)練圖像往往只是捕獲了某一平面的孔隙空間模式特征,并不能完全反映孔隙空間三維幾何和拓?fù)浣Y(jié)構(gòu)特征如配位數(shù),而三維訓(xùn)練圖像則能提供更加真實(shí)的孔隙空間模型。目標(biāo)邊緣概率可以約束隨機(jī)生成模型的孔隙和顆粒比例,一般設(shè)置該值為目標(biāo)孔隙度和顆粒比例。還可以通過修改Servosystem系數(shù)控制孔隙所占的比例,Servosystem系數(shù)取值范圍為0～1,取值越大,模擬結(jié)果中孔隙和顆粒的比例與訓(xùn)練圖像中的兩者比例相差越小,但是該系數(shù)很高時(shí)建模效果反而會(huì)變差[31],建模時(shí)設(shè)置Servosystem系數(shù)為0.9。為了保證獲取概率較低的孔隙結(jié)構(gòu)模式,數(shù)據(jù)事件重復(fù)次數(shù)設(shè)置為1。按照以上步驟生成3個(gè)1503體素(網(wǎng)格)的隨機(jī)模型(圖3(d)～(f))。

3 重構(gòu)模型不確定性評(píng)價(jià)

對于數(shù)字巖心分析技術(shù),數(shù)字巖心模型的滲流特性等傳輸性質(zhì)本質(zhì)上取決于孔隙空間的幾何結(jié)構(gòu)和拓?fù)浣Y(jié)構(gòu)特征,重構(gòu)模型的可靠性不是指重構(gòu)模型孔隙空間如孔隙和喉道與真實(shí)模型中相同位置上的孔隙和喉道一定要相同,而是指整個(gè)巖心的孔隙空間幾何形態(tài)和拓?fù)浣Y(jié)構(gòu)具有等效性,相同概率的孔隙、喉道和孔喉配置關(guān)系出現(xiàn)即可。通過對比重構(gòu)模型和真實(shí)模型的變差函數(shù)、孔隙結(jié)構(gòu)參數(shù)和滲流特性可評(píng)價(jià)重構(gòu)模型的不確定性。

3.1 變差函數(shù)

變差函數(shù)可以反映地質(zhì)變量在空間兩點(diǎn)之間的相關(guān)性,可以評(píng)價(jià)孔隙結(jié)構(gòu)在空間中的變化。對真實(shí)模型和3個(gè)重構(gòu)模型在x、y和z方向上求取變差函數(shù),得到相應(yīng)的曲線(圖4(a)～(c))。計(jì)算結(jié)果表明重構(gòu)模型和真實(shí)模型在x、y和z方向上變差函數(shù)曲線十分接近,3個(gè)方向上的變差函數(shù)基臺(tái)值變化范圍均為0.15～0.17,滯后距變化范圍為18～20個(gè)體素,這表明重構(gòu)模型與真實(shí)模型在x、y和z某一方向上的孔隙空間變化具有很好的一致性。

圖4重構(gòu)模型與真實(shí)模型在x、y和z方向的變差函數(shù)曲線Fig.4 Variogram curves of reconstructed models and real model in x, y and z directions

3.2 孔隙結(jié)構(gòu)參數(shù)

僅依靠變差函數(shù)并不能表征孔隙結(jié)構(gòu)復(fù)雜的空間變化,本文中還利用能夠表征孔喉幾何形態(tài)和拓?fù)浣Y(jié)構(gòu)的孔隙結(jié)構(gòu)參數(shù)對重構(gòu)模型的準(zhǔn)確性進(jìn)一步評(píng)估。

孔隙度是表征孔隙空間的基本參數(shù)之一。在目標(biāo)邊緣概率和Servosystem 系數(shù)的控制下,重構(gòu)模型的孔隙度與真實(shí)模型的相差很小,最大差值僅為0.46%(表1)。根據(jù)孔隙是否連通,可把孔隙空間分為連通孔隙和孤立孔隙?？紫犊臻g的連通情況可用連通系數(shù)表示,連通系數(shù)越大,孔隙連續(xù)性越好。對總孔隙空間做連通性測試可以提取連通孔隙部分(圖5),再用減法運(yùn)算得到孤立孔隙(圖5)。真實(shí)模型中孤立孔隙所占孔隙比例不到3%(表1),連通性較好,隨機(jī)生成模型孤立孔隙比例略大,約為6%～8%,4個(gè)模型的總孔隙和連通孔隙體積相差很小(表1),表明重構(gòu)模型也具有較好的長距離連通性。迂曲度是表征孔隙空間復(fù)雜程度和彎曲程度的一個(gè)參數(shù),是指沿著某個(gè)方向流體質(zhì)點(diǎn)在孔隙中流經(jīng)的實(shí)際長度與孔隙兩端點(diǎn)間距離的比值,該值越大表示孔隙結(jié)構(gòu)越復(fù)雜。重構(gòu)模型的迂曲度與真實(shí)模型的迂曲度比較接近(表1),表明重構(gòu)模型與真實(shí)模型孔隙空間的復(fù)雜程度類似。

圖5 真實(shí)模型和重建模型的孔隙空間分布對比Fig.5 Comparison on pore space distribution of real model and reconstructed models

參數(shù)孔隙度/%總體積/μm3連通體積/μm3連通系數(shù)/%迂曲度 真實(shí)模型19.95102815787.40100091169.0097.351.91重構(gòu)模型119.90102558103.7096373850.0693.971.79重構(gòu)模型219.74101733515.9093248940.7291.661.97重構(gòu)模型319.49100445097.6093263273.0892.851.93

根據(jù)所建模型的拓?fù)湫再|(zhì),可以把孔隙網(wǎng)絡(luò)模型分為兩大類,即為規(guī)則拓?fù)浣Y(jié)構(gòu)孔隙網(wǎng)絡(luò)模型和真實(shí)拓?fù)浣Y(jié)構(gòu)孔隙網(wǎng)絡(luò)模型。真實(shí)拓?fù)浣Y(jié)構(gòu)孔隙網(wǎng)絡(luò)模型是基于數(shù)字巖心提取與孔隙空間幾何特征和拓?fù)涮卣鞯葍r(jià)的孔隙網(wǎng)絡(luò)。提取孔隙網(wǎng)絡(luò)模型的方法主要包括多向切片掃描法、孔喉居中軸線法、Voronoi多面體法和最大球法[30]。利用改進(jìn)的最大球法[32]提取真實(shí)模型和重構(gòu)模型的孔隙網(wǎng)絡(luò)模型并計(jì)算孔隙結(jié)構(gòu)參數(shù)包括孔喉半徑、孔喉比、孔喉形狀因子和配位數(shù)。在孔隙結(jié)構(gòu)參數(shù)中,形狀因子是描述多孔介質(zhì)模型孔隙喉道截面形態(tài)的重要參數(shù),也是影響微觀滲流的重要因素。形狀因子G的計(jì)算公式為

(3)

式中,A為孔隙或喉道的截面面積,μm2;P為孔隙或喉道截面的周長,μm。形狀因子G越大,孔喉截面越趨于圓滑,流體滲流的阻力越小。

真實(shí)模型和重構(gòu)模型的孔喉半徑、孔喉比、孔喉形狀因子和配位數(shù)的分布見圖6。盡管重構(gòu)模型中小孔隙和喉道(半徑小于10 μm)比例略高于真實(shí)模型中的小孔隙比例,但重構(gòu)模型和真實(shí)模型的孔喉半徑分布直方圖和孔喉比分布曲線總體上比較接近(圖6(a)～(c))。重構(gòu)模型孔隙中孤立孔隙比例較大(圖6(d)),約為13%;4個(gè)模型的配位數(shù)分布絕大多數(shù)集中在1～4,總體分布一致。重構(gòu)模型和真實(shí)模型孔喉形狀因子分布曲線集中區(qū)間均為0.02～0.04(圖6(e)～(f)),4個(gè)模型的孔隙和喉道的幾何形態(tài)都比較接近。

圖6 真實(shí)模型和重構(gòu)模型孔隙結(jié)構(gòu)參數(shù)分布Fig.6 Distribution of pore structure parameters for real model and reconstructed models

孔隙網(wǎng)絡(luò)模型能直觀地展示孔隙空間展布和孔喉配置關(guān)系。采用改進(jìn)后的最大球法提取4個(gè)模型的孔隙網(wǎng)絡(luò)模型(圖7),圖中球體為孔隙,管柱體為喉道。球體或管柱體越大,表示此處的孔隙或喉道半徑越大。無管柱體連接的球體表示此處的孔隙為孤立孔隙,孔隙配位數(shù)為零。觀察孔隙網(wǎng)絡(luò)模型能夠看出4個(gè)模型大多數(shù)孔隙被2～4個(gè)喉道連接,總體的孔喉配置關(guān)系相似。重構(gòu)模型和真實(shí)模型的孔隙結(jié)構(gòu)參數(shù)的相似性進(jìn)一步驗(yàn)證了重構(gòu)模型的孔隙空間幾何特征和拓?fù)涮卣髋c真實(shí)模型的相近性。

圖7 真實(shí)模型和重構(gòu)模型孔隙網(wǎng)絡(luò)模型Fig.7 Pore network models of real model and reconstructed models

3.3 滲流特性

基于Navier-Stokes方程和達(dá)西定律,以真實(shí)和重構(gòu)數(shù)字巖心模型為平臺(tái),進(jìn)行單相滲流數(shù)值模擬試驗(yàn)。滲流過程為穩(wěn)態(tài)流動(dòng),流體為不可壓縮的牛頓流體,流固界面為無滑移壁面。模型的z軸方向一端設(shè)置為入口,另一端為出口,入口壓力為130 000 Pa,出口為100 000 Pa,流體黏度設(shè)置為0.001 Pa·s。

3.3.1 絕對滲透率對比

根據(jù)單相流滲流試驗(yàn),計(jì)算出真實(shí)模型和重構(gòu)模型z軸方向上的絕對滲透率,真實(shí)模型、重構(gòu)模型1、重構(gòu)模型2、重構(gòu)模型3測得的絕對滲透率分別為1 032.56×10-3、986.96×10-3、993.34×10-3、975.27×10-3μm2,該樣品試驗(yàn)滲透率為1 100×10-3μm2。重構(gòu)模型的絕對滲透率與真實(shí)模型模擬出的和試驗(yàn)測得的絕對滲透率值比較接近,進(jìn)一步證明了重構(gòu)模型在滲流特性上的可靠性。重構(gòu)模型的滲透率略小于真實(shí)模型的滲透率,原因是重構(gòu)模型中孔隙的平均配位數(shù)略低于真實(shí)模型的平均配位數(shù)。

3.3.2 流速場對比

為進(jìn)一步對比重構(gòu)模型和真實(shí)模型單相滲流試驗(yàn)滲流的特性,應(yīng)用Navier-Stokes方程計(jì)算單相流的滲流速度,并使其可視化。圖8(a)顯示了流線在4個(gè)模型孔隙空間內(nèi)部的三維分布,流線從z軸負(fù)半軸流向正半軸。孔隙空間用藍(lán)色表示,流線流經(jīng)的部位表示流體在孔隙內(nèi)部的滲流路徑,流線的顏色深淺代表滲流的速度,顏色越紅速度越大,顏色偏綠速度偏小。圖8(a)顯示4個(gè)模型的流線均集中在連通性較好的大孔喉處,邊界處連通性不好的孔喉流線稀疏,孤立的孔隙內(nèi)部無流線;在半徑較大配位數(shù)較多的孔隙處,滲流阻力小,流體流量大,流速會(huì)相對較大,流線越紅;4個(gè)模型大孔喉處的流線顏色均偏紅,表明滲流速度比較接近。

圖8 真實(shí)模型和重構(gòu)模型的孔隙空間內(nèi)部流速場和壓力場分布Fig.8 Velocity field and pressure field distribution in pore space

3.3.3 壓力場對比

壓力場分布范圍也能反映孔喉分布。在單相流試驗(yàn)中,從入口端到出口端連通孔喉內(nèi)的壓力逐漸降低,同一壓降下,壓力場分布范圍與孔喉中毛細(xì)管力有關(guān),即是與孔喉半徑有關(guān)。圖8(b)展示了真實(shí)模型和重構(gòu)模型從入口端(z=0處)開始的3個(gè)壓降變化過程。顏色越淺表示壓力越小。同一壓降梯度下,壓力場分布范圍越大代表孔隙結(jié)構(gòu)中的毛細(xì)管力越小,孔喉半徑越大;4個(gè)模型在相同的壓降梯度下壓力場范圍比較接近,這進(jìn)一步證明了重構(gòu)模型的可靠性。

4 結(jié) 論

(1)提出了基于多點(diǎn)地質(zhì)統(tǒng)計(jì)學(xué)法構(gòu)建三維數(shù)字巖心模型的建模方法,創(chuàng)建了一套以三維圖像為訓(xùn)練圖像,二維切片為條件數(shù)據(jù),利用多點(diǎn)地質(zhì)統(tǒng)計(jì)學(xué)法重構(gòu)數(shù)字巖心模型的流程。

(2)以變差函數(shù)曲線、孔隙結(jié)構(gòu)參數(shù)和單相流滲流特性為評(píng)價(jià)參數(shù),驗(yàn)證了重構(gòu)模型與真實(shí)模型的孔隙空間具有相似的幾何結(jié)構(gòu)和拓?fù)浣Y(jié)構(gòu)特征。

(3)針對非均質(zhì)性較弱的巖石,多點(diǎn)地質(zhì)統(tǒng)計(jì)學(xué)能夠構(gòu)建較理想的三維數(shù)字巖心模型,但是該方法建立的多孔介質(zhì)模型中會(huì)出現(xiàn)較大比例的微小孔隙。該方法并不太適用非均質(zhì)性較強(qiáng)的巖心如包含孔洞的碳酸鹽巖,原因是復(fù)雜巖心中某些特殊的孔隙結(jié)構(gòu)模式如碳酸鹽巖中的孔洞占總體孔隙結(jié)構(gòu)模式的比例很小,使得SNESIM算法往往很難把這些概率極低的孔隙結(jié)構(gòu)模式復(fù)制到重構(gòu)模型中。

(4)重構(gòu)模型評(píng)價(jià)中僅使用了變差函數(shù)、孔隙結(jié)構(gòu)參數(shù)和單相流滲流特性等評(píng)價(jià)參數(shù),還可使用其他函數(shù)如多點(diǎn)概率分布函數(shù)和相對滲透率曲線等進(jìn)一步驗(yàn)證重建模型的準(zhǔn)確性。如果計(jì)算機(jī)配置允許,利用多點(diǎn)地質(zhì)統(tǒng)計(jì)學(xué)可以重構(gòu)更大尺度的數(shù)字巖心,再以重構(gòu)的數(shù)字巖心為平臺(tái)進(jìn)行滲流試驗(yàn)數(shù)值模擬,這對研究油氣在不同尺度上的滲流機(jī)制及建立微納米孔喉尺度的滲流機(jī)制與厘米—米巖心或儲(chǔ)層尺度的宏觀滲流規(guī)律之間的聯(lián)系具有更加重要的意義。

致謝感謝帝國理工大學(xué)Branko Bijeljic高級(jí)研究員提供的Berea砂巖數(shù)據(jù)。

:

[1] 孫建孟,閆國亮,姜黎明,等.基于數(shù)字巖心研究流體性質(zhì)對裂縫性低滲透儲(chǔ)層彈性參數(shù)的影響規(guī)律[J].中國石油大學(xué)學(xué)報(bào)(自然科學(xué)版),2014,38(3):39-44.

SUN Jianmeng, YAN Guoliang, JIANG Liming, et al. Research of influence laws of fluid properties on elastic parameters of fractured low permeability reservoir rocks based on digital core[J]. Journal of China University of Petroleum (Edition of Natural Science), 2014,38(3):39-44.

[2] 楊永飛,姚軍,王晨晨.水濕油藏油氣水三相滲流模擬[J].中國石油大學(xué)學(xué)報(bào)(自然科學(xué)版),2010,34(1):79-83.

YANG Yongfei, YAO Jun, WANG Chenchen. Oil-gas-water three phase flow simulation in water-wet reservoir[J]. Journal of China University of Petroleum (Edition of Natural Science), 2010,34(1):79-83.

[3] 鄒友龍,謝然紅,郭江峰,等.致密儲(chǔ)層數(shù)字巖心重構(gòu)及核磁共振響應(yīng)模擬[J].中國石油大學(xué)學(xué)報(bào)(自然科學(xué)版),2015,39(6):63-71.

ZOU Youlong, XIE Ranhong, GUO Jiangfeng, et al. Reconstruction of digital core of tight reservoir and simulation of NMR response[J]. Journal of China University of Petroleum (Edition of Natural Science), 2015,39(6):63-71.

[4] BLUNT M J, BIJELJIC B, DONG H, et al. Pore-scale imaging and modelling[J]. Advances in Water Resources, 2013,51(3):197-216.

[5] 林承焰,吳玉其,任麗華,等.數(shù)字巖心建模方法研究現(xiàn)狀及展望[J].地球物理學(xué)進(jìn)展,2018,33(2):679-689.

LIN Chengyan, WU Yuqi, REN Lihua, et al. Review of digital core modeling methods[J]. Progress in Geophysics, 2018,33(2):679-689.

[6] VOGEL H J, ROTH K. Quantitative morphology and network representation of soil pore structure[J]. Advances in Water Resources, 2001,24(3/4):233-242.

[7] FREDRICH J T. 3D imaging of porous media using laser scanning confocal microscopy with application to microscale transport processes[J]. Physics and Chemistry of the Earth, Part A: Solid Earth and Geodesy, 1999,24(7):551-561.

[8] QUIBLIER J A. A new three-dimensional modeling technique for studying porous media[J]. Journal of Colloid and Interface Science, 1984,98(1):84-102.

[9] TALUKDAR M S, TORSAETER O, HOWARD J J. Stochastic reconstruction of chalk samples containing vuggy porosity using a conditional simulated annealing technique[J]. Transport in Porous Media, 2004,57(1):1-15.

[10] JU Y, HUANG Y, ZHENG J, et al. Multi-thread parallel algorithm for reconstructing 3D large-scale porous structures[J]. Computers & Geosciences, 2017,101:10-20.

[11] WU K, NUNAN N, CRAWFORD J W, et al. An efficient Markov chain model for the simulation of heterogeneous soil structure[J]. Soil Science Society of America Journal, 2004,68(2):346-351.

[12] KEEHM Y. Permeability prediction from thin sections: 3D reconstruction and Lattice-Boltzmann flow simulation[J]. Geophysical Research Letters, 2004,31(4):1-4.

[13] OKABE H, BLUNT M J. Prediction of permeability for porous media reconstructed using multiple-point statistics[J]. Physical Review E, Statistical, Nonlinear, and Soft Matter Physics, 2004,70(6 Pt 2):66135.

[14] HASANABADI A, BANIASSADI M, ABRINIA K, et al. 3D microstructural reconstruction of heterogeneous materials from 2D cross sections: a modified phase-recovery algorithm[J]. Computational Materials Science, 2016,111:107-115.

[15] HASANABADI A, BANIASSADI M, ABRINIA K, et al. Efficient three-phase reconstruction of heterogeneous material from 2D cross-sections via phase-recovery algorithm[J].Journal of Microscopy, 2016,264(3):384-393.

[16] ?REN P L, BAKKE S. Process based reconstruction of sandstones and prediction of transport properties[J]. Transport in Porous Media, 2002,46(2):311-343.

[17] ?REN P L, BAKKE S. Reconstruction of Berea sandstone and pore-scale modelling of wettability effects[J]. Journal of Petroleum Science and Engineering, 2003,39(3/4):177-199.

[18] POLITIS M G, KIKKINIDES E S, KAINOURGIAKIS M E, et al. A hybrid process-based and stochastic reconstruction method of porous media[J]. Microporous and Mesoporous Materials, 2008,110(1):92-99.

[19] LIU X F, SUN J M, WANG H T. Reconstruction of 3-D digital cores using a hybrid method[J]. Applied Geophysics, 2009,6(2):105-112.

[20] ZHAO X, YAO J, YI Y. A new stochastic method of reconstructing porous media[J]. Transport in Porous Media, 2007,69(1):1-11.

[21] HAJIZADEH A, SAFEKORDI A, FARHADPOUR F A. A multiple-point statistics algorithm for 3D pore space reconstruction from 2D images[J]. Advances in Water Resources, 2011,34(10):1256-1267.

[22] DEUTSCH, CLAYTONV. GSLIB geostatistical software library and users guide[M]. Oxford: Oxford University Press, 1992:126.

[23] 吳勝和.儲(chǔ)層表征與建模[M].北京:石油工業(yè)出版社,2010:372-378.

[24] TAHMASEBI P, HEZARKHANI A, SAHIMI M. Multiple-point geostatistical modeling based on the cross-correlation functions[J]. Computational Geosciences, 2012,16(3):779-797.

[25] STRAUBHAAR J, WALGENWITZ A, RENARD P. Parallel multiple-point statistics algorithm based on list and tree structures[J]. Mathematical Geosciences, 2013,45(2):131-147.

[26] 張偉,林承焰,董春梅.多點(diǎn)地質(zhì)統(tǒng)計(jì)學(xué)在秘魯D油田地質(zhì)建模中的應(yīng)用[J].中國石油大學(xué)學(xué)報(bào)(自然科學(xué)版),2008,32(4):24-28.

ZHANG Wei, LIN Chengyan, DONG Chunmei. Application of multiple-point geostatistics in geological modeling of D Oilfield in Peru[J]. Journal of China University of Petroleum (Edition of Natural Science), 2008,32(4):24-28.

[27] STREBELLE S. Conditional simulation of complex geological structures using multiple-point statistics[J]. Mathematical Geology, 2002,34(1):1-21.

[28] REMY N, BOUCHER A, WU J. Applied geostatistics with SGeMS: a users guide [M]. Cambridge: Cambridge University Press, 2009:168-205.

[29] 張挺,李道倫,盧德唐,等.基于多點(diǎn)地質(zhì)統(tǒng)計(jì)法的多孔介質(zhì)重構(gòu)研究[J].中國科學(xué)(G輯:物理學(xué) 力學(xué) 天文學(xué)),2009,39(9):1348-1360.

ZHANG Ting, LI Daolun, LU Detang, et al.Research on the reconstruction method of porous media using multiple-point geostatistics[J]. Science in China Series G: Physics, Mechanics & Astronomy, 2009,39(9):1348-1360.

[30] DONG H, BLUNT M J. Pore-network extraction from micro-computerized-tomography images[J]. Physical Review E, Statistical, Nonlinear, and Soft Matter Physics, 2009,80(3 Pt 2):36307.

[31] LIU Y. Using the snesim program for multiple-point statistical simulation[J]. Computers & Geosciences, 2006,32(10):1544-1563.

[32] 閆國亮. 基于數(shù)字巖心儲(chǔ)層滲透率模型研究[D]. 青島:中國石油大學(xué)(華東), 2013.

YAN Guoliang. Research of permeability models of reservoirs based on digital cores [D].Qingdao: China University of Petroleum (East China), 2013.