基于集合光滑的深度學(xué)習(xí)自動(dòng)歷史擬合方法

馬小鵬, 張 凱, 陳昕晟, 曹慶平, 姚傳進(jìn), 谷建偉

(1.中國(guó)石油大學(xué)(華東)石油工程學(xué)院,山東青島 266580; 2.中國(guó)石油大港油田分公司石油工程技術(shù)研究院，天津 300450)

油氣儲(chǔ)層精細(xì)建模對(duì)于指導(dǎo)油氣資源高效開(kāi)發(fā)具有重要意義,中國(guó)油氣資源以陸相沉積環(huán)境為主,陸相含油氣盆地占陸上油資源總量的90%[1],另外超過(guò)80%油氣藏采用注水開(kāi)發(fā),復(fù)雜的地質(zhì)條件和開(kāi)發(fā)動(dòng)態(tài)導(dǎo)致油氣藏精細(xì)描述十分困難。自動(dòng)歷史擬合技術(shù)[2]依據(jù)油氣田開(kāi)發(fā)生產(chǎn)動(dòng)態(tài)對(duì)測(cè)井、巖心等靜態(tài)數(shù)據(jù)建立的隨機(jī)油藏模型進(jìn)行動(dòng)態(tài)反演調(diào)整是目前油藏精細(xì)建模的唯一有效手段。盡管發(fā)展了許多高效的自動(dòng)歷史擬合求解方法[3-14],但根據(jù)動(dòng)態(tài)數(shù)據(jù)反演調(diào)整儲(chǔ)層模型仍然是一項(xiàng)艱巨的任務(wù),尤其是針對(duì)大規(guī)模復(fù)雜油藏建立高效的自動(dòng)歷史擬合方法依舊存在許多挑戰(zhàn)[15]。大規(guī)模油藏?cái)?shù)百萬(wàn)網(wǎng)格帶來(lái)的維數(shù)災(zāi)難、復(fù)雜地質(zhì)特征(如河流相分布)難以參數(shù)化表征等迫切需要行之有效的解決方法。近年來(lái),機(jī)器學(xué)習(xí)的興起為石油工程領(lǐng)域的發(fā)展帶來(lái)新的契機(jī),在地質(zhì)建模[16]、測(cè)井分析[17]、產(chǎn)量預(yù)測(cè)[18]等多個(gè)領(lǐng)域取得成功應(yīng)用。借助于機(jī)器學(xué)習(xí)理論方法,結(jié)合隨機(jī)奇異值分解(SVD)[19-21]的深度變分自編碼(DVAE)模型[22-26],筆者提出SVD-DVAE模型，將高維地質(zhì)特征(如滲透率映射到低維連續(xù)高斯空間)進(jìn)行重參數(shù)化表征,結(jié)合集合光滑多次數(shù)據(jù)同化(ES-MDA)方法[27-28]吸收生產(chǎn)歷史,實(shí)現(xiàn)復(fù)雜大規(guī)模油藏反演建模[29-31]高效求解。

1 方 法

1.1 油藏自動(dòng)歷史擬合數(shù)學(xué)模型

油藏模型參數(shù)m和油藏開(kāi)發(fā)生產(chǎn)動(dòng)態(tài)響應(yīng)dobs可以看作是油藏系統(tǒng)的輸入與輸出,表示為

dobs=g(m)+ε.

(1)

式中,g(·)為油藏系統(tǒng)的數(shù)值模型或油藏?cái)?shù)值模擬器;ε為生產(chǎn)動(dòng)態(tài)響應(yīng)的觀測(cè)誤差。

(2)

油藏生產(chǎn)開(kāi)發(fā)是由狀態(tài)方程、流動(dòng)方程等偏微分方程組成的復(fù)雜非線性系統(tǒng),難以直接建立油藏模型參數(shù)估計(jì)f,通常使用貝葉斯概率推理方法建立油藏模型參數(shù)估計(jì)數(shù)學(xué)模型。貝葉斯概率推理方法在水文地質(zhì)、石油工程、天氣預(yù)報(bào)等領(lǐng)域獲得廣泛使用,表示為

(3)

式中,P(m)為先驗(yàn)概率,表示模型參數(shù)的已知信息,如通過(guò)地震、測(cè)井等方式獲得的靜態(tài)數(shù)據(jù);P(dobs|m)稱為似然函數(shù),用于量化觀測(cè)數(shù)據(jù)與估計(jì)模型參數(shù)進(jìn)行數(shù)值模擬預(yù)測(cè)之間的擬合誤差;P(m|dobs)即后驗(yàn)概率,表示在已知觀測(cè)數(shù)據(jù)dobs的情形下模型參數(shù)m的概率分布。后驗(yàn)概率越大表明估計(jì)模型參數(shù)接近真實(shí)模型參數(shù)的可能性越大。

基于貝葉斯概率統(tǒng)計(jì),可以獲得油藏模型參數(shù)的后驗(yàn)概率估計(jì):

(4)

根據(jù)地震、測(cè)井、巖心取樣等手段取得油藏的靜態(tài)觀測(cè)數(shù)據(jù),使用隨機(jī)地質(zhì)統(tǒng)計(jì)建模方法可以建立大量的隨機(jī)初始油藏模型mi(i=1,2,…,Ne),結(jié)合高斯概率模型可建立油藏模型參數(shù)的先驗(yàn)概率,表示為

(5)

為便于計(jì)算,假設(shè)觀測(cè)誤差ε服從高斯分布,即ε～N(0,CD)。根據(jù)式(1)可以推得似然函數(shù)為

(6)

結(jié)合式(3)、(5)、(6)可得到油藏模型參數(shù)的后驗(yàn)概率為

(7)

由于式(7)涉及數(shù)值模擬g(m),難以獲得解析表達(dá)式,因此通常使用優(yōu)化算法、集合卡爾曼濾波、蒙特卡洛模擬等方法進(jìn)行數(shù)值求解。ES-MDA方法是Emerick和Reynolds 等提出的一種基于集合的數(shù)據(jù)同化歷史擬合求解方法。相比于傳統(tǒng)的集合卡爾曼濾波方法,ES-MDA方法具有收斂快速、不需要數(shù)值模擬多次重啟動(dòng)、模型不易崩潰等優(yōu)點(diǎn)。ES-MDA方法通過(guò)使用觀測(cè)數(shù)據(jù)的誤差協(xié)方差矩陣Cd乘以一個(gè)膨脹因子α多次迭代進(jìn)行數(shù)據(jù)同化,有效地避免了觀測(cè)數(shù)據(jù)過(guò)擬合,從而解決集合卡爾曼濾波方法模型崩潰的缺陷。ES-MDA方法中油藏模型參數(shù)的估計(jì)表示為

(8)

膨脹因子αt是ES-MDA方法中唯一需要確定的輔助參數(shù),對(duì)求解結(jié)果具有顯著影響,Emerick等提出了自適應(yīng)膨脹因子的ES-MDA方法。本文中使用標(biāo)準(zhǔn)ES-MDA方法中膨脹因子的設(shè)置,即αt=1/N,N為膨脹因子,通常取4～10。

1.2 深度自編碼變分模型

復(fù)雜大規(guī)模油藏的地質(zhì)特征(如滲透率、孔隙度等)通常需要數(shù)十萬(wàn)乃至上百萬(wàn)維的網(wǎng)格變量進(jìn)行地質(zhì)建模,而且非連續(xù)的地質(zhì)屬性參數(shù)化表征(如河流相等復(fù)雜相分布參數(shù))是自動(dòng)歷史擬合領(lǐng)域的關(guān)鍵難題,研究使用深度學(xué)習(xí)模型對(duì)復(fù)雜大規(guī)模油藏地質(zhì)特征進(jìn)行參數(shù)化表征。深度生成模型是近年來(lái)深度學(xué)習(xí)領(lǐng)域的研究熱點(diǎn)之一,在機(jī)器視覺(jué)、圖像生成、自動(dòng)駕駛汽車(chē)等領(lǐng)域得到廣泛關(guān)注。變分自編碼模型以及對(duì)抗性神經(jīng)網(wǎng)絡(luò)模型是深度生成模型領(lǐng)域最重要的兩種模型,其中變分自編碼模型假設(shè)數(shù)據(jù)是由隨機(jī)過(guò)程產(chǎn)生的,并涉及一個(gè)未觀察到的或不可觀察到的連續(xù)隨機(jī)變量z,稱為潛變量,這個(gè)隨機(jī)過(guò)程由兩步組成:①根據(jù)先驗(yàn)分布p(z)生成潛變量z;②給定z,根據(jù)條件分布p(x|z;θ)生成相應(yīng)的數(shù)據(jù)x。

根據(jù)貝葉斯理論,x的邊際似然為

(9)

根據(jù)式(9)的解析表達(dá)即可直接生成數(shù)據(jù)x。但是,P(x|z)通常涉及復(fù)雜的非線性方程,導(dǎo)致不可解析求解。對(duì)于這種不可解析求解的概率分布問(wèn)題最常用的方法有蒙特卡洛采樣和變分推理。蒙特卡洛采樣方法計(jì)算十分昂貴,變分自編碼模型從變分方法進(jìn)行推導(dǎo),得到邊際似然P(x)的變分下界:

logP(x)≥L(θ,φ;x)=Εqφ(z|x)[-logqφ(z|x)+

logPθ(x|z)].

(10)

最大化式(10)即可獲得邊際似然P(x)的最大似然估計(jì),變分自編碼模型采用隨機(jī)梯度技術(shù)以及自編碼神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)進(jìn)行優(yōu)化求解,對(duì)式(10)進(jìn)行隨機(jī)梯度變分貝葉斯估計(jì):

(11)

式(11)右側(cè)KL散度作為正則項(xiàng),第二項(xiàng)是負(fù)重構(gòu)誤差的期望。圖1為變分自編碼模型的結(jié)構(gòu)示意圖。編碼器(encoder)作為qφ(z|x),解碼器(decoder)作為Pθ(x|z),通過(guò)對(duì)數(shù)據(jù)集進(jìn)行訓(xùn)練即可獲得合適的模型參數(shù)θ和φ。

圖1 變分自編碼模型結(jié)構(gòu)示意圖

使用初始油藏模型集合進(jìn)行訓(xùn)練,潛變量z是油藏模型的低維表征,輸入解碼層可實(shí)現(xiàn)從低維表征重構(gòu)相應(yīng)的油藏模型。VAE模型需要確定的超參數(shù)包括隱藏層的層數(shù)、神經(jīng)元個(gè)數(shù)以及最關(guān)鍵的潛變量z的維數(shù)。

1.3 奇異值分解本征維數(shù)估計(jì)方法

奇異值分解方法(singular value decomposition,SVD)是一種無(wú)監(jiān)督學(xué)習(xí)方法,在數(shù)據(jù)降維、消除冗余、特征提取和噪音消除等領(lǐng)域都有廣泛應(yīng)用。另外,SVD方法還可用于高維數(shù)據(jù)的本征維數(shù)估計(jì),如確定線性模型中的變量數(shù)或確定動(dòng)力系統(tǒng)的自由度等。提出使用SVD方法對(duì)滲透率、孔隙度等高維模型參數(shù)進(jìn)行本征維數(shù)估計(jì),作為變分編碼模型潛變量層的維數(shù)。

SVD方法原理見(jiàn)圖2。任給M×N的矩陣A,可以分解為UM×M、ΣM×N和VN×N三個(gè)矩陣:

圖2 奇異值分解原理示意圖

A=U∑VT.

(12)

矩陣Σ是一個(gè)對(duì)角矩陣,對(duì)角線元素的值σi稱為奇異值,且奇異值具有沿著對(duì)角線逐漸減小的性質(zhì)。通過(guò)保留r?min(M,N)個(gè)奇異值可實(shí)現(xiàn)分解矩陣的壓縮,并用壓縮后的分解矩陣乘積近似原始矩陣A:

(13)

計(jì)算奇異值所能保留矩陣A的信息貢獻(xiàn)率,確定矩陣的本征維數(shù)為

(14)

該方法適用于小樣本情形,對(duì)于大量樣本情形,由于樣本數(shù)量增加導(dǎo)致次要奇異值部分被高估,無(wú)法真實(shí)反映奇異值的信息保留貢獻(xiàn)率。因此提出奇異值差分譜方法,確定矩陣本征維數(shù),奇異值差分譜定義為

λi=|σi+1-σi|,i=1,…,k-1.

(15)

(16)

通過(guò)式(15)差分操作,可以有效消除次要奇異值干擾。通過(guò)式(16)計(jì)算相應(yīng)的信息保留貢獻(xiàn)率,可確定本征維數(shù)。

1.4 結(jié)合ES-MDA和VAE-SVD方法的綜合自動(dòng)歷史擬合流程

整合奇異值分解方法、變分編碼模型以及ES-MDA數(shù)據(jù)同化方法形成綜合的自動(dòng)歷史擬合求解方法,如圖3所示。首先根據(jù)初始油藏模型計(jì)算地質(zhì)特征的本征維數(shù),并訓(xùn)練變分編碼模型。其次使用變分編碼模型的解碼器對(duì)油藏模型進(jìn)行降維參數(shù)化表征,ES-MDA迭代同化油田歷史觀測(cè)數(shù)據(jù),更新降維潛變量;使用變分編碼模型的解碼器重構(gòu)油藏模型并進(jìn)行數(shù)值模擬計(jì)算,用于下一次潛變量迭代更新,多次迭代求解,直至滿足收斂條件,如達(dá)到提前設(shè)定的迭代次數(shù),輸出最終擬合模型。

圖3 基于深度學(xué)習(xí)與數(shù)據(jù)同化方法的綜合自動(dòng)歷史擬合流程

2 算 例

對(duì)一個(gè)二維河流相油藏進(jìn)行自動(dòng)歷史擬合研究,擬合參數(shù)為滲透率。該油藏真實(shí)滲透率分布如圖4(a)所示,低滲透背景相的滲透率為500×10-3μm2,高滲透河流相為5 000×10-3μm2。圖4(b)為初始模型集合中隨機(jī)抽取的28個(gè)模型,初始模型可根據(jù)訓(xùn)練圖像采用多點(diǎn)地質(zhì)統(tǒng)計(jì)隨機(jī)建模方法生成。不同于常規(guī)油藏,河流相油藏地質(zhì)相之間存在明顯的不連續(xù)邊界,因此滲透率的反演更具挑戰(zhàn)性。該油藏模型有45×45網(wǎng)格,每個(gè)網(wǎng)格30.5 m×30.5 m,厚度是8 m;含有3口注水井、4口生產(chǎn)井,采用Eclipse軟件作為油藏?cái)?shù)值模擬器。為增加歷史擬合預(yù)測(cè)模型的挑戰(zhàn)性,僅使用300 d共10個(gè)時(shí)間步的生產(chǎn)數(shù)據(jù)進(jìn)行自動(dòng)歷史擬合測(cè)試,以4口生產(chǎn)井的日產(chǎn)水量和日產(chǎn)油量數(shù)值模擬結(jié)果添加5%的高斯噪音作為觀測(cè)數(shù)據(jù)。

圖4 河流相油藏算例

2.1 SVD方法估計(jì)本征維數(shù)

利用SVD方法對(duì)2 000個(gè)隨機(jī)初始油藏模型進(jìn)行處理,計(jì)算初始油藏樣本的奇異值,如圖5(a)所示,少量的奇異值明顯突出。根據(jù)式(14)和式(16)計(jì)算奇異值貢獻(xiàn)率如圖5(b)和(c),可以看出,式(14)對(duì)樣本數(shù)量敏感,無(wú)法準(zhǔn)確估計(jì)奇異值對(duì)信息保留的貢獻(xiàn)率;式(16)中通過(guò)差分計(jì)算能夠顯著剔除次要奇異值的干擾,可以更加準(zhǔn)確地估計(jì)本征維數(shù)。

圖5 隨機(jī)初始油藏模型集合SVD計(jì)算奇異值及信息保留貢獻(xiàn)率

根據(jù)式(16)計(jì)算奇異值貢獻(xiàn)率D,并選取貢獻(xiàn)率為0.8、0.9及0.95時(shí)對(duì)應(yīng)的奇異值個(gè)數(shù),對(duì)應(yīng)的油藏模型滲透率場(chǎng)的本征維數(shù)估計(jì)分別為13、33和65。使用估計(jì)的本征維數(shù)作為變分編碼模型的潛變量層神經(jīng)元數(shù),建立深度變分編碼模型,模型結(jié)構(gòu)見(jiàn)表1。使用2 000個(gè)先驗(yàn)隨機(jī)油藏模型中1 800個(gè)模型作為訓(xùn)練樣本,另外200個(gè)模型作為驗(yàn)證樣本對(duì)模型進(jìn)行訓(xùn)練。

表1 變分自編碼模型結(jié)構(gòu)

隨機(jī)選取初始油藏模型中的一個(gè)滲透率場(chǎng),對(duì)變分編碼模型重構(gòu)結(jié)果進(jìn)行分析,見(jiàn)圖6(a)。圖6(b)～(d)為不同本征維數(shù)時(shí)變分自編碼模型對(duì)圖6(a)所示滲透率場(chǎng)的重構(gòu)結(jié)果。可看出高滲相的分布趨勢(shì)均能夠得到較好的重構(gòu)。主要差異在于相邊界,本征維數(shù)r為13時(shí)相邊界模糊,而r為65時(shí)相邊界重構(gòu)的更好。網(wǎng)格滲透率的分布頻率直方圖能更加直觀地說(shuō)明重構(gòu)差異,如圖6(e)所示,初始模型的滲透率顯著地集中在高滲相和低滲相。而重構(gòu)模型都存在中間過(guò)渡,圖6(f)中r=13時(shí),存在較多的網(wǎng)格滲透率分布在高滲和低滲值之間,r=65時(shí)對(duì)應(yīng)的網(wǎng)格滲透率分布在高、低滲兩端更加集中。

圖6 不同潛變量維數(shù)對(duì)應(yīng)的變分自編碼模型重構(gòu)滲透率場(chǎng)與參考樣本對(duì)比

圖7對(duì)VAE模型隨機(jī)生成滲透率場(chǎng)的能力進(jìn)行了對(duì)比。r取13和65時(shí)都能夠很好地重構(gòu)具有河流相特征的滲透率模型,但r取65時(shí)河流相邊界更加清晰。

圖7 不同潛變量維數(shù)變分自編碼模型生成滲透率場(chǎng)效果對(duì)比

2.2 歷史擬合測(cè)試

對(duì)提出的ES-MDA結(jié)合VAE-SVD方法同化生產(chǎn)觀測(cè)數(shù)據(jù)進(jìn)行的自動(dòng)歷史擬合進(jìn)行了測(cè)試,并與ES-MDA結(jié)合SVD的傳統(tǒng)降維方法進(jìn)行了對(duì)比。圖8(a)展示了保留前65個(gè)奇異值的SVD方法結(jié)合ES-MDA方法進(jìn)行歷史擬合的預(yù)測(cè)模型,可看出,SVD方法能夠粗略地捕捉河流相的分布,但是無(wú)法預(yù)測(cè)河流相的準(zhǔn)確分布。如圖8(b)和(c)所示,ES-MDA結(jié)合VAE-SVD方法均能夠很好地預(yù)測(cè)高滲相的大致分布及形狀。另外對(duì)比取不同的本征維數(shù)擬合效果,可以看出,r取65時(shí)高滲相的相邊界預(yù)測(cè)更加清晰,對(duì)高滲相的分布能夠做出準(zhǔn)確預(yù)測(cè)。

圖8 自動(dòng)歷史擬合滲透率場(chǎng)反演結(jié)果對(duì)比

圖9對(duì)4口生產(chǎn)井的日產(chǎn)油量歷史擬合結(jié)果進(jìn)行了說(shuō)明,其中綠色菱形點(diǎn)表示觀測(cè)數(shù)據(jù)點(diǎn),灰色曲線是初始油藏模型集合的數(shù)值模擬預(yù)測(cè)結(jié)果,紅色曲線是歷史擬合更新模型集合的數(shù)值模擬預(yù)測(cè)結(jié)果。與初始隨機(jī)油藏模型集合相比,ES-MDA結(jié)合VAE-SVD自動(dòng)歷史擬合方法預(yù)測(cè)的油藏模型集合能夠很好地反映真實(shí)觀測(cè)數(shù)據(jù)的變化。

圖9 ES-MDA結(jié)合VAE-SVD方法生產(chǎn)觀測(cè)數(shù)據(jù)擬合結(jié)果

3 結(jié) 論

(1)使用奇異值分解方法對(duì)高維油藏模型參數(shù)進(jìn)行本征維數(shù)估計(jì)可以指導(dǎo)深度變分自編碼模型潛變量維數(shù)的設(shè)定。直接基于奇異值計(jì)算貢獻(xiàn)率方法在大樣本數(shù)量情形下估計(jì)的本征維數(shù)偏高,而奇異值差分譜可以消除樣本數(shù)量帶來(lái)的偏差,計(jì)算的貢獻(xiàn)率與奇異值分布規(guī)律保持一致,有助于深度學(xué)習(xí)等大樣本情景下準(zhǔn)確估計(jì)本征維數(shù)。

(2)參數(shù)場(chǎng)保留信息貢獻(xiàn)率與相應(yīng)的估計(jì)本征維數(shù)成正比,保留信息越多深度變分自編碼模型重構(gòu)效果越好,但是相應(yīng)地潛變量維數(shù)也越高;考慮到數(shù)據(jù)同化算法對(duì)于參數(shù)維數(shù)敏感性較低,因此保持盡可能高的信息貢獻(xiàn)率有助于精確地反演油藏模型參數(shù)。

(3)基于神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)的深度變分自編碼器能夠?qū)⒏呔S網(wǎng)格變量壓縮至低維連續(xù)變量,為大規(guī)模油藏歷史擬合求解提供了思路。與基于傳統(tǒng)的奇異值分解等矩陣分解自動(dòng)歷史擬合方法相比,綜合數(shù)據(jù)同化算法與深度學(xué)習(xí)模型構(gòu)建的自動(dòng)歷史擬合綜合流程能夠更好地保持油藏地質(zhì)特征。