景文龍, 李博涵, 楊守磊, 張 磊, 孫 海, 楊永飛, 李愛芬

(中國(guó)石油大學(xué)(華東)石油工程學(xué)院,山東青島 266580)

近年來,機(jī)器學(xué)習(xí)技術(shù)在石油行業(yè)中的應(yīng)用越來越多。Shirangi[1]通過人工神經(jīng)網(wǎng)絡(luò)、支持向量機(jī)等模型進(jìn)行油氣藏生產(chǎn)優(yōu)化;Deng等[2]基于原油價(jià)格影響因素,構(gòu)建了價(jià)格預(yù)測(cè)和交易回歸模型;Kidd[3]應(yīng)用機(jī)器學(xué)習(xí)技術(shù),判別海灣水域石油污染程度和污染位置; Chaki[4]在油藏描述當(dāng)中利用自適應(yīng)神經(jīng)模糊推理系統(tǒng)和模糊邏輯系統(tǒng)等機(jī)器學(xué)習(xí)工具提升儲(chǔ)層表征的精確性。機(jī)器學(xué)習(xí)適用于具有大量歷史數(shù)據(jù)并需要解決大規(guī)模變量問題、分析過程較為復(fù)雜的問題[5]。微觀數(shù)值模擬結(jié)果可以為機(jī)器學(xué)習(xí)提供充足、準(zhǔn)確的數(shù)據(jù)庫[6]。因此可以利用該技術(shù)進(jìn)行油氣滲流方面的理論研究。國(guó)內(nèi)外學(xué)者對(duì)基于數(shù)字巖心的滲透率計(jì)算已開展了大量的研究。Knackstedt等[7-8]研究了碳酸鹽巖溶洞大孔和微孔的特征,發(fā)現(xiàn)不同尺度的孔隙對(duì)碳酸鹽巖介質(zhì)的巖石滲透率以及其他物理特征和流體滲流性質(zhì)有重要影響;張磊等[9-12]模擬了真實(shí)頁巖巖心中的氣體流動(dòng)和有機(jī)質(zhì)中氣體的解析/擴(kuò)散現(xiàn)象,比較了3種滲透率計(jì)算方法,優(yōu)選了孔隙網(wǎng)絡(luò)模型法為滲透率的計(jì)算方法。筆者基于數(shù)字巖心和機(jī)器學(xué)習(xí)技術(shù),隨機(jī)生成大量數(shù)字巖心,通過孔隙網(wǎng)絡(luò)模型的方法對(duì)其進(jìn)行孔隙度和滲透率計(jì)算,將所得結(jié)果作為機(jī)器學(xué)習(xí)的數(shù)據(jù)庫,基于TensorFlow平臺(tái),利用BP人工神經(jīng)網(wǎng)絡(luò)方法對(duì)數(shù)據(jù)進(jìn)行提取和處理,對(duì)于給出的任意巖心孔隙度,通過機(jī)器學(xué)習(xí)可準(zhǔn)確預(yù)測(cè)出該巖心的滲透率。

1 基于TensorFlow的機(jī)器學(xué)習(xí)

1.1 TensorFlow系統(tǒng)

TensorFlow是一個(gè)可以將復(fù)雜的數(shù)據(jù)結(jié)構(gòu)體傳遞輸入至人工智能神經(jīng)網(wǎng)絡(luò)中并能夠進(jìn)行分析和處理過程的新一代人工智能學(xué)習(xí)系統(tǒng)。TensorFlow包含了兩層含義,首先Tensor,即張量,代表了N維的數(shù)組,而Flow,即流動(dòng),表示基于數(shù)據(jù)流圖的計(jì)算。TensorFlow可以理解為數(shù)組從數(shù)據(jù)流圖的一個(gè)端點(diǎn)流動(dòng)至另一個(gè)端點(diǎn)的計(jì)算過程。

TensorFlow系統(tǒng)在語音以及圖像識(shí)別等機(jī)器學(xué)習(xí)與深度學(xué)習(xí)領(lǐng)域應(yīng)用較廣。本文中所采用的機(jī)器學(xué)習(xí)方法——BP人工神經(jīng)網(wǎng)絡(luò)就是在TensorFlow平臺(tái)上搭建的。BP人工神經(jīng)網(wǎng)絡(luò)跟其他機(jī)器學(xué)習(xí)方法相比,其最大的優(yōu)點(diǎn)是采用多層網(wǎng)絡(luò),可以解決線性不可分的分類問題,而例如感知學(xué)習(xí)規(guī)則和梯度算法只能解決線性可分的分類問題,有一定的局限性;并且BP人工神經(jīng)網(wǎng)絡(luò)采用鏈條法來求解偏導(dǎo)數(shù),可以處理數(shù)字巖心滲透率的多個(gè)影響因素,故BP人工神經(jīng)網(wǎng)絡(luò)是最適合本次研究的機(jī)器學(xué)習(xí)方法。

1.2 基于BP人工神經(jīng)網(wǎng)絡(luò)的深度學(xué)習(xí)

人工神經(jīng)網(wǎng)絡(luò)的工作原理是按照信息處理的方式對(duì)人腦神經(jīng)元網(wǎng)絡(luò)進(jìn)行抽象化處理,從而建立某種簡(jiǎn)單的模型,并且按照不同的連接方式構(gòu)成不同的網(wǎng)絡(luò)。一般的人工神經(jīng)網(wǎng)絡(luò)建立的過程包含3個(gè)部分,這3個(gè)部分分別為單元/神經(jīng)元、連接/權(quán)重/參數(shù)、偏置項(xiàng)[13-16]。

神經(jīng)元作為神經(jīng)網(wǎng)絡(luò)架構(gòu)3個(gè)部分中最重要的部分,包含權(quán)重和偏置項(xiàng)的函數(shù),等待數(shù)據(jù)傳遞。接收數(shù)據(jù)后,執(zhí)行一些計(jì)算,然后使用激活函數(shù)將數(shù)據(jù)限制在一個(gè)范圍內(nèi)。將這些單元想象成一個(gè)包含權(quán)重和偏置項(xiàng)的盒子。盒子從兩端打開,一端接收數(shù)據(jù),另一端輸出修改后的數(shù)據(jù)。數(shù)據(jù)首先進(jìn)入盒子中,將權(quán)重與數(shù)據(jù)相乘,再向相乘的數(shù)據(jù)添加偏置項(xiàng)[17]。這是一個(gè)單元,也可以被認(rèn)為是一個(gè)函數(shù),如圖1所示。x表示神經(jīng)元i的輸入信號(hào),w表示連接權(quán)重,是利用神經(jīng)網(wǎng)絡(luò)解決問題時(shí)必須學(xué)習(xí)的數(shù)值;b為外部刺激,代表神經(jīng)網(wǎng)絡(luò)認(rèn)為其在將權(quán)重與數(shù)據(jù)相乘之后應(yīng)該添加的內(nèi)容;f為激勵(lì)函數(shù),也稱為映射函數(shù),在x軸上輸入數(shù)據(jù),并在有限的范圍內(nèi)輸出一個(gè)值,大多數(shù)情況下,它們被用于將單元中較大的值輸出轉(zhuǎn)換成較小的值;y為神經(jīng)元的輸出,其輸入輸出關(guān)系式為

圖1 典型人工神經(jīng)網(wǎng)絡(luò)的工作結(jié)構(gòu)

(1)

建立神經(jīng)元模型后,將各個(gè)神經(jīng)元進(jìn)行全連接便可建立人工神經(jīng)網(wǎng)絡(luò)模型。隨著機(jī)器學(xué)習(xí)和人工神經(jīng)網(wǎng)絡(luò)研究的不斷發(fā)展,目前國(guó)內(nèi)外學(xué)者已經(jīng)提出了多種人工神經(jīng)網(wǎng)絡(luò),其中誤差反向傳播神經(jīng)網(wǎng)絡(luò)是較常見的神經(jīng)網(wǎng)絡(luò)。本文中應(yīng)用BP人工神經(jīng)網(wǎng)絡(luò),其拓?fù)浣Y(jié)構(gòu)如圖2所示。

從圖2中可以發(fā)現(xiàn),全連接層的每一個(gè)結(jié)點(diǎn)都與上一層的所有結(jié)點(diǎn)相連,可以把前邊提取到的特征綜合起來。由于其全相連的特性,一般全連接層的參數(shù)也是最多的,需要消耗較大內(nèi)存。

圖2 BP人工神經(jīng)網(wǎng)絡(luò)拓?fù)浣Y(jié)構(gòu)

以3層人工神經(jīng)網(wǎng)絡(luò)為例,其拓?fù)浣Y(jié)構(gòu)如圖3所示,其輸入表達(dá)式為

圖3 三層人工神經(jīng)網(wǎng)絡(luò)拓?fù)浣Y(jié)構(gòu)

(2)

把X、θ、a分別用矩陣表示,可以得到:θX=a。輸出層表達(dá)式為

(3)

從式(3)中可以發(fā)現(xiàn)對(duì)于每個(gè)ai,都依賴上層的X以及X對(duì)應(yīng)的參數(shù)θ,這樣從左至右逐級(jí)依賴的算法模型稱為前向傳播算法。

BP人工神經(jīng)網(wǎng)絡(luò)在人工智能領(lǐng)域應(yīng)用廣泛,但在石油工程領(lǐng)域的應(yīng)用相對(duì)較少。本文中將探討B(tài)P人工神經(jīng)網(wǎng)絡(luò)在數(shù)字巖心領(lǐng)域的應(yīng)用。

2 數(shù)字巖心計(jì)算及模型預(yù)測(cè)

2.1 數(shù)字巖心孔隙度和滲透率計(jì)算原理

對(duì)經(jīng)過X射線CT掃描生成的數(shù)字巖心圖像進(jìn)行濾波和分割處理。圖像濾波方法中最常用的為中值濾波法,是一種以排序統(tǒng)計(jì)為理論支撐的抑制圖像噪點(diǎn)的非線性信號(hào)處理技術(shù)。其原理是把數(shù)字圖像或數(shù)字序列中一點(diǎn)所代表的值用該點(diǎn)所在位置的鄰域中各點(diǎn)值的中值所替代,讓其周圍的像素值更加接近真實(shí)值,從而消除孤立的噪聲點(diǎn)。圖像分割方法中最常用的為分水嶺分割法,是一種基于形態(tài)學(xué)針對(duì)圖像梯度的分割算法。圖像分割是為了將巖心的CT掃描圖像分割成孔隙相和骨架相兩部分。數(shù)字巖心的孔隙度就是依據(jù)圖像分割后孔隙相的體積占巖心總體積的百分比。

數(shù)字巖心的滲透率計(jì)算方法為孔隙網(wǎng)絡(luò)模型法。在孔隙網(wǎng)絡(luò)模型中,影響滲透率的因素有形狀因子、傳導(dǎo)率以及閥壓。其中形狀因子G為

(4)

式中,A為孔隙孔喉截面面積,m2;P為孔隙孔喉截面形狀的周長(zhǎng),m。

在實(shí)際儲(chǔ)層中,巖心的孔隙和喉道的截面形狀是不斷改變的,其面積和周長(zhǎng)沿著孔隙和喉道中軸線不斷變化。在建立孔隙網(wǎng)絡(luò)模型時(shí),統(tǒng)計(jì)巖心的孔隙和喉道中軸線不同位置處的截面面積和周長(zhǎng),取平均值,之后用式(2)計(jì)算可得孔隙網(wǎng)絡(luò)模型中對(duì)應(yīng)的孔隙與喉道的形狀因子[18]。

建立的孔隙網(wǎng)絡(luò)模型為擬靜態(tài)網(wǎng)絡(luò)模型,特點(diǎn)是流體流動(dòng)時(shí)毛管力的作用遠(yuǎn)大于黏滯力的作用。具體體現(xiàn)在:確定毛管力后,界面的位置隨即確定,忽略壓力傳播的動(dòng)態(tài)變化,不考慮時(shí)間的影響,驅(qū)替順序由孔喉的毛管壓力決定。假設(shè)流體為牛頓流體,不可壓縮且黏度為常數(shù)。結(jié)合連續(xù)性方程和N-S方程可以得到橢圓型Possion方程描述不同截面形狀的毛細(xì)管中的流動(dòng),

(5)

式中,v為流體速度,m/s;μ為流體黏度,mPa·s;ρ為流體密度,kg/m3;f為單位質(zhì)量的體積力,為常數(shù)。

求解式(3)～(5)時(shí),首先建立坐標(biāo)系如圖4所示。其中z軸的方向?yàn)榇怪庇诩埫嫦蚶锊⑴c毛管中軸平行。

圖4 三角形孔隙截面坐標(biāo)系

假設(shè)沿z方向(垂直于紙面向里)施加一不變的壓力梯度,

(6)

式中,L為毛細(xì)管長(zhǎng)度,m;p0和pL分別為毛細(xì)管入口和出口的壓力,MPa;p為沿毛細(xì)管所產(chǎn)生的壓力降,MPa。

在該壓力梯度下,流體沿毛細(xì)管內(nèi)流動(dòng)。因此上述問題可轉(zhuǎn)化為一定邊界條件下的求解式,即

(7)

v(x,y)=0.

(8)

式中,v(x,y)為流體沿x和y方向的流動(dòng)速度,m/s。

對(duì)式(3)～(5)進(jìn)行求解后可得毛細(xì)管中體積流量與驅(qū)替壓力梯度的關(guān)系式為

(9)

利用孔隙網(wǎng)絡(luò)模型進(jìn)行流動(dòng)模擬時(shí),首先設(shè)定模型中所有孔隙和喉道都飽和水,之后進(jìn)行油驅(qū)水,孔隙和喉道中的驅(qū)替方式為活塞式驅(qū)替,其中驅(qū)替相和被驅(qū)替相之間形成彎液面,產(chǎn)生的毛細(xì)管力由Young-Laplace方程確定,

(10)

式中,pcow為閾壓,MPa;σow為油水界面張力,mN/m;r1和r2為界面曲面的兩個(gè)主半徑,m。

若孔隙形狀與油水接觸角已知,則可利用式(10)求毛細(xì)管壓力,并且只有當(dāng)驅(qū)替相壓力大于該毛細(xì)管力,驅(qū)替相才能夠驅(qū)替被驅(qū)替相進(jìn)入孔隙或喉道[18]。

建立孔隙網(wǎng)絡(luò)模型的一般方法是在一定的邊界條件下,對(duì)N個(gè)孔隙建立質(zhì)量守恒方程,并通過求解方程組得到孔隙壓力。對(duì)于可壓縮流體,必須滿足孔隙i的質(zhì)量守恒。對(duì)于牛頓流體,應(yīng)用連續(xù)性方程可以得到:

(11)

式中,qij是從孔隙j通過其連通喉道進(jìn)入相鄰孔隙i的體積流量,m3/s;pi和pj分別為孔隙i和j中的壓力,MPa;δij為連接孔隙i和j的孔喉的水力傳導(dǎo)系數(shù),m/s。

方程(11)表明孔隙i中流體的凈累積量為零。通過求解方程組,計(jì)算各孔隙的壓力和各孔喉處的流量,便可求得巖心滲透率。

2.2 數(shù)字巖心構(gòu)建過程

為了給機(jī)器學(xué)習(xí)提供原始數(shù)據(jù),首先需要進(jìn)行大量的隨機(jī)均質(zhì)數(shù)字巖心的構(gòu)建,并對(duì)隨機(jī)生成的均質(zhì)數(shù)字巖心進(jìn)行孔隙度和滲透率的計(jì)算,其中滲透率的計(jì)算是利用孔隙網(wǎng)絡(luò)模型的方法。

(1)選取CT掃描之后的巖心切片(200 μm×200 μm),該切片具有明顯的孔隙和骨架結(jié)構(gòu),利用馬爾可夫鏈蒙特卡羅法(MCMC)進(jìn)行數(shù)字巖心的構(gòu)建,如圖5所示。具體構(gòu)建步驟[19]:首先設(shè)定一個(gè)體素,作為所構(gòu)建的一維鏈的起點(diǎn),令其條件概率為孔隙;之后沿y方向模擬第一層第一行上的體素;然后沿x方向模擬第一層其他行的體素,繼而完成巖心第一層的構(gòu)建;最后,沿z方向模擬各層的體素進(jìn)而完成整個(gè)三維數(shù)字巖心的構(gòu)建。

圖5 隨機(jī)生成數(shù)字巖心

(2)采用分水嶺分割的方法對(duì)生成的數(shù)字巖心進(jìn)行圖像分割,將巖心的孔隙相和骨架相分別提取出來,如圖6所示。并計(jì)算其孔隙度。

圖6 數(shù)字巖心分割

(3)從數(shù)字巖心的孔隙相中篩選連通孔隙,依據(jù)孔隙半徑的差別區(qū)分出不同的連通孔隙集合,最終建立孔隙網(wǎng)絡(luò)模型,如圖7所示?；诳紫毒W(wǎng)絡(luò)模型模擬計(jì)算數(shù)字巖心滲透率。

圖7 孔隙網(wǎng)絡(luò)模型建立

2.3 機(jī)器學(xué)習(xí)訓(xùn)練及驗(yàn)證

隨機(jī)生成300塊均質(zhì)數(shù)字巖心,并對(duì)其進(jìn)行孔隙度和滲透率的計(jì)算,通過對(duì)300組孔隙度和滲透率數(shù)據(jù)的學(xué)習(xí),最終實(shí)現(xiàn)在給定孔隙度的前提下預(yù)測(cè)滲透率,某區(qū)塊的真實(shí)巖心如圖8所示。

圖8 某區(qū)塊真實(shí)巖心

洗油烘干后對(duì)其進(jìn)行孔隙度和滲透率的試驗(yàn)測(cè)試,測(cè)試結(jié)果為φ=14.64%,k=4.16×10-3μm2。將孔隙度加入到機(jī)器學(xué)習(xí)中,最終預(yù)測(cè)的滲透率為4.03×10-3μm2。可以發(fā)現(xiàn),通過機(jī)器學(xué)習(xí)所得的巖心滲透率與試驗(yàn)所測(cè)巖心滲透率基本相等,誤差僅為3.1%,能夠滿足油氣田開發(fā)現(xiàn)場(chǎng)工程應(yīng)用的要求。

3 結(jié)束語

通過隨機(jī)生成大量均質(zhì)數(shù)字巖心,并利用孔隙網(wǎng)絡(luò)模型的方法對(duì)其進(jìn)行孔隙度和滲透率的計(jì)算,將所得的結(jié)果作為機(jī)器學(xué)習(xí)的數(shù)據(jù)庫,之后基于BP人工神經(jīng)網(wǎng)絡(luò)的方法對(duì)巖心的孔隙度、滲透率數(shù)據(jù)進(jìn)行提取和處理,對(duì)于給出的任意巖心孔隙度,通過機(jī)器學(xué)習(xí)可準(zhǔn)確地預(yù)測(cè)出該巖心的滲透率,最后進(jìn)行實(shí)例驗(yàn)證。結(jié)果表明,通過機(jī)器學(xué)習(xí)技術(shù)預(yù)測(cè)滲透率的方法準(zhǔn)確高效,與巖心的實(shí)測(cè)滲透率誤差僅為3.1%,可在實(shí)際生產(chǎn)中進(jìn)行應(yīng)用,避免了大量的試驗(yàn)操作,提高了巖心滲透率的計(jì)算效率。