張 凱, 馬小鵬, 王增飛, 劉 凡, 馬 瑋, 姚 軍

(1.中國石油大學(華東)石油工程學院,山東青島 266580; 2.北京中油瑞飛信息技術(shù)有限公司,北京 100007;3.中海油研究總院，北京 100010; 4.海洋石油高效開發(fā)國家重點實驗室，北京 100028;5.勝利油田油氣井下作業(yè)中心，山東東營 257001)

目前中國大多數(shù)油田經(jīng)過長期的勘探開發(fā),已經(jīng)進入注水開發(fā)的中后期階段,尤其是砂巖油氣藏儲層之間的非均質(zhì)特征非常強[1-3],嚴重影響了油藏注水開發(fā)效果。開展強非均質(zhì)性油藏中高滲通道識別技術(shù)的研究,對于油氣資源開發(fā)具有重要意義。Brigham等[4-9]利用油藏動靜態(tài)資料等識別復雜油氣藏中的高滲通道,方法成本高,適用性局限。自動歷史擬合技術(shù)依靠油藏數(shù)值模擬,根據(jù)油藏動態(tài)數(shù)據(jù)反演調(diào)整地質(zhì)參數(shù),具有很強的可靠性。目前油藏自動歷史擬合方法眾多[10-15],但是缺乏專門針對強非均質(zhì)性油藏的方法。筆者對強非均質(zhì)性油藏自動歷史擬合方法展開研究,提出混合求解方法:PCA方法[16-19]與DCT方法[20-23]混合;SPSA算法[24-27]與ABC算法[28-29]混合。

1 強非均質(zhì)性油藏自動歷史擬合方法

1.1 歷史擬合數(shù)學模型的建立

自動歷史擬合問題是一個典型的反問題,通常結(jié)合先驗地質(zhì)信息進行求解能夠得到接近準確解的結(jié)果[11]。為了得到準確的油藏地質(zhì)特征,需要利用已知的觀測數(shù)據(jù)與油藏模型參數(shù)m之間的關(guān)系建立數(shù)學模型,表達式為

dobs=g(m)+εr.

(1)

式中,dobs為生產(chǎn)動態(tài)數(shù)據(jù);g(·)為油藏數(shù)值模擬器;εr為觀測數(shù)據(jù)誤差。

為了求解反問題,通常對原函數(shù)進行正則化處理,利用先驗信息作為約束,對所求問題作合適的形式轉(zhuǎn)換,使反問題的解更加準確。其中貝葉斯理論就是正則化方法之一,基于貝葉斯理論,觀測數(shù)據(jù)與油藏模型參數(shù)之間關(guān)系式[30]為

P(m|dobs)∝P(dobs|m)P(m).

(2)

式中,m為油藏模型參數(shù);P(m|dobs)為貝葉斯理論中的后驗概率;P(dobs|m)為貝葉斯理論中的似然函數(shù);P(m)為貝葉斯理論中的先驗概率。

在地質(zhì)參數(shù)m給定的條件下,獲得的觀測數(shù)據(jù)存在觀測誤差εr,一般認為誤差符合均值為0,協(xié)方差矩陣為CD的高斯型概率分布,即εr～(0,CD)。因此可得到似然函數(shù)P(dobs|m)為

(3)

為了求得先驗概率,須利用已知的油藏信息建立先驗模型。本文中采用隨機建模[21]方法中的序貫高斯模擬方法生成多個可選的、等概率的先驗模型。目前隨機建模的方法包括序貫高斯模擬、截斷高斯模擬以及分形模擬等方法[31-32]。其中序貫高斯模擬主要運用高斯模型,根據(jù)離散的觀測點生成一個連續(xù)的隨機面,地質(zhì)參數(shù)符合高斯分布規(guī)律,能夠很好地反映實際地質(zhì)情況,在構(gòu)建地質(zhì)模型領(lǐng)域得到了廣泛的應用。因此可以得到先驗概率為

(4)

根據(jù)貝葉斯理論可以求得后驗概率為

P(m|d)∝exp[-O(m)].

(5)

其中

式中,O(m)為歷史擬合數(shù)學模型的目標函數(shù);mpr為先驗模型參數(shù)。

后驗概率P(m|dobs)中不僅包含了觀測數(shù)據(jù),還包含了先驗地質(zhì)信息。則歷史擬合反演問題轉(zhuǎn)化為使后驗概率P(m|dobs)最大化問題,此過程稱為MAP(Maximum A Posteriori)。后驗概率越大,自動歷史擬合所獲得模型參數(shù)越能反映真實的地質(zhì)情況。以O(shè)(m)作為目標函數(shù),則問題就轉(zhuǎn)換成求解目標函數(shù)最小值的問題。

1.2 PCA特征提取方法

主成分分析方法(principal component analysis,PCA)是目前應用很廣泛的一種特征提取方法,在數(shù)據(jù)降維和特征提取方面的有效性使得它在人臉識別領(lǐng)域獲得了廣泛的應用[16-17,19]。它主要利用K-L變換從原始的數(shù)據(jù)集中提取主要特征并構(gòu)成特征向量空間,而后將原有的數(shù)據(jù)影射到特征向量空間上并得到一組投影系數(shù),能夠很好地保留原始數(shù)據(jù)的主要特征信息。利用PCA方法對先驗模型的非均質(zhì)特征進行提取,能夠很好地提取出先驗模型中高(低)滲區(qū)域的信息,從而準確地反演出強非均質(zhì)油藏中高(低)滲區(qū)域。

首先根據(jù)先驗信息生成Ne個先驗模型,記為mj(j=1,2,…,Ne),平均值為mpr。對于實際油藏的歷史擬合問題,由于先驗模型滿足序列高斯分布,則其平均模型為

(6)

(7)

(8)

直接求CM的特征值和特征向量計算量較大,因此求解矩陣V=ATA的特征值λi及對應的正交化特征向量vi,減小計算量的同時更加準確地反映出先驗模型參數(shù)場的特征,選取前p(p≤Ne?Nm)個最大特征值及對應的特征向量;這些特征向量包含了先驗模型中的主要特征信息。利用式(8)求取協(xié)方差矩陣CM的正交歸一化特征向量,獲得特征向量提取矩陣Φ即為原參數(shù)場的特征提取場,

Φ=[φ1,φ2,…,φe],

(9)

(10)

對圖1表示的對數(shù)滲透率場進行特征提取,可見最大特征值λ1對應的特征提取場φ1中存在1條明顯的低(高)滲通道,特征值為λ5對應的特征提取場φ5中的高(低)滲通道條數(shù)增加到3條。隨著特征值的減小,對應的特征提取場中的高(低)滲通道數(shù)越來越多,因此選擇對應特征值大的特征提取場進行歷史擬合能夠減少先驗模型中無關(guān)信息的干擾,更加準確地反演出高(低)滲通道。

圖1 PCA方法特征提取場Fig.1 PCA method to extract feature field

根據(jù)PCA方法的K-L變換原理,目標函數(shù)中協(xié)方差矩陣的逆矩陣可以被近似表示為

(11)

實現(xiàn)特征提取場Φ代替先驗模型參數(shù)場M進行自動歷史擬合。此時,定義新的參數(shù)s(s∈Rp)對油藏模型參數(shù)進行參數(shù)化變換,

(12)

則目標函數(shù)變?yōu)?/p>

(13)

其中

使用PCA特征提取方法可將油藏模型的歷史擬合問題由Nm維降低到p維,并且有效避免了協(xié)方差矩陣求逆的過程,大規(guī)模油藏進行自動歷史擬合時由于擬合參數(shù)的維數(shù)過于龐大,CM的逆矩陣計算十分耗時,避免目標函數(shù)中協(xié)方差矩陣CM的逆矩陣計算能夠大幅度提高計算速度。結(jié)合優(yōu)化求解方法,通過更新參數(shù)s降低目標函數(shù)值,并基于式(13)反求油藏參數(shù),實現(xiàn)強非均質(zhì)性油藏歷史擬合問題的求解。

1.3 基于離散余弦變換的PCA降維方法

PCA方法對于高滲通道區(qū)分明顯的油藏模型歷史擬合,能夠很好地反演出接近真實地質(zhì)情況的油藏模型。但對于地質(zhì)參數(shù)如孔隙度、滲透率等分布復雜的油氣藏,PCA方法并不能夠很好地反演出準確的儲層參數(shù),因此在PCA方法的基礎(chǔ)上提出PCA方法與DCT方法混合進行特征提取。離散余弦變換方法(discrete cosine transform,DCT)是將一系列的離散點表示為頻率不同的余弦函數(shù)的疊加。在圖像處理領(lǐng)域,利用DCT方法良好的能量壓縮性能,通過保留大系數(shù)來保存大部分的圖像信息。在進行自動歷史擬合時選取低頻系數(shù)對先驗模型進行降維,能夠保留先驗模型的大部分信息,從而更加精準地描述油藏地質(zhì)信息。

對油藏先驗模型參數(shù)場M進行離散余弦變換,D是參數(shù)場離散余弦變換后的DCT系數(shù)矩陣。定義DCT變換矩陣W1和W2,其元素分別定義為

(14)

(15)

(16)

對圖2表示的滲透率場進行DCT變換得到降維后的特征提取場,該變換能夠更好的保留先驗數(shù)據(jù)的信息,增強了歷史擬合的魯棒性。

本文中利用PCA和DCT結(jié)合進行先驗模型的特征提取,PCA特征提取的參數(shù)場ΦP為(p/2)×Ne維,DCT變換特征提取場ΦD為(p/2)×Ne維。將兩者組合得到特征提取場Φ代入原來的PCA方法過程即可實現(xiàn)PCA方法與DCT方法的混合,

Φ=[ΦP;ΦD].

(17)

圖2 DCT方法特征提取Fig.2 Discrete cosine transform to extract feature

1.4 歷史擬合目標函數(shù)的優(yōu)化求解

隨機擾動梯度近似算法[24-26](simultaneous perturbation stochastic approximation,SPSA)是通過添加擾動和重復迭代方法尋找線性或非線性目標函數(shù)的局部最優(yōu)值,適用于多維變量系統(tǒng),同其他隨機優(yōu)化方法(如遺傳算法)比較,可以與任何數(shù)值模擬器耦合。在實際應用中算法的近似梯度與實際梯度方向十分接近,具有很強的收斂性。但在迭代的過程中,由于步長的選取不當,可能會導致目標函數(shù)陷入局部最優(yōu)。本文中提出人工蜂群算法(ABC)與SPSA算法混合對目標函數(shù)進行求解。ABC算法[28-29]是一種較為系統(tǒng)的群集智能隨機優(yōu)化算法,具有操作簡單、控制參數(shù)少、搜索精度較高和魯棒性較強的特點,全局搜索能力較好,但局部搜索能力較差。將ABC算法的全局搜索模式引入SPSA算法中,可以提高SPSA算法的全局搜索能力。

(18)

(19)

(20)

在迭代的過程中,由于步長的選取具有局限性,可能會導致目標函數(shù)陷入局部最優(yōu),此時需要引入ABC算法的搜索模式尋求全局最優(yōu)值。在ABC算法中如果在一個食物源(問題的解)的位置處不能在預先設(shè)定的循環(huán)次數(shù)limit內(nèi)找到更優(yōu)的解,則該食物源被放棄并且被偵察蜂找到的新食物源所代替,偵察蜂搜索新的食物源依據(jù)為

sk+1=smin+rand[0,1](smax-smin).

(21)

將該搜索模式引入SPSA算法中,改進SPSA算法流程如圖3所示。

圖3 SPSA與ABC算法混合優(yōu)化流程Fig.3 SPSA and ABC algorithm hybrid optimization flow chart

2 實例應用

2.1 PCA方法自動歷史擬合實例

對一個油氣水三相的油藏模型進行PCA方法歷史擬合研究,采用Eclipse軟件作為油藏數(shù)值模擬器。該油藏存在三條明顯的高滲通道,需要反演調(diào)整的地質(zhì)參數(shù)為滲透率。用序貫高斯模擬隨機生成100個先驗模型,取部分如圖4所示。利用PCA方法對目標函數(shù)中的先驗模型參數(shù)矩陣進行特征提取,利用ABC算法與SPSA算法混合進行優(yōu)化求解,滲透率場的反演結(jié)果如圖5所示?？梢钥闯?當?shù)郊s第5步時,注水井Inj-1、生產(chǎn)井Pro-2、Pro-5連線上出現(xiàn)了明顯的高滲條帶,當?shù)?5步時,生產(chǎn)井Pro-4、Pro-3連線上同樣出現(xiàn)了一條明顯的滲透率條帶。由于生產(chǎn)井Pro-6是單口井生產(chǎn),與其他井之間不存在連通性,所以滲透率值變化很慢,直到迭代到25步時,才獲得一條模糊的滲透率條帶,至此反演所獲得的滲透率場已經(jīng)非常接近真實模型的滲透率場。結(jié)果表明,利用PCA方法提取地質(zhì)參數(shù)中的非均質(zhì)特征能夠很好地反演出油藏中的高(低)滲通道;同時將高維地質(zhì)參數(shù)數(shù)據(jù)轉(zhuǎn)化為低維參數(shù)使收斂速度更快。

2.2 PCA結(jié)合DCT方法歷史擬合實例

以一個油氣水三相的復雜油藏模型,采用Eclipse軟件作為油藏數(shù)值模擬器,進行PCA與DCT混合特征提取方法歷史擬合研究,對滲透率場進行反演。采用序貫高斯模擬隨機生成100個先驗模型,選取部分如圖6所示。然后利用PCA與DCT方法混合對先驗模型進行特征提取,最后利用SPSA與ABC算法混合進行求解,結(jié)果如圖7所示,可以看出求得的滲透率場在生產(chǎn)井Pro-1、Pro-4和注水井Inj-5、Inj-9的連線上出現(xiàn)了高滲條帶;在生產(chǎn)井Pro-3和注水井Inj-4、Inj-7以及生產(chǎn)井Pro-2和注水井Inj-2、Inj-6連線上也出現(xiàn)了高滲條帶。在交叉井之間,利用該方法反演出的滲透率值也十分接近真實滲透率場。對比PCA與DCT方法混合與傳統(tǒng)的SVD方法進行歷史擬合所獲得的結(jié)果,可以看出相比PCA結(jié)合DCT方法,SVD方法反演的滲透率場中高滲區(qū)域并不明顯,與真實滲透率場的差距較大。與單一的PCA方法進行歷史擬合相比,PCA與DCT混合的方法能夠獲得更加精細的擬合結(jié)果,滲透率場更加接近真實地質(zhì)情況。

圖4 初始隨機滲透率場實現(xiàn)Fig.4 Initial random permeability field

圖5 滲透率場反演過程Fig.5 Permeability field inversion process

圖6 隨機滲透率場實現(xiàn)Fig.6 Random permeability field

圖7 滲透率場反演結(jié)果Fig.7 Results of history matching

圖8和9分別為累積指標和單井指標擬合結(jié)果。由圖8可以看出,PCA與DCT方法混合擬合得到的累積產(chǎn)水量與地層壓力隨時間的變化曲線能夠很好地反映真實生產(chǎn)數(shù)據(jù)的變化;由圖9可以看出,擬合得到的生產(chǎn)井Pro-1的日產(chǎn)油量以及生產(chǎn)井Pro-3的日產(chǎn)水量隨時間的變化曲線很好地反映了真實生產(chǎn)數(shù)據(jù)的變化。滲透率場與生產(chǎn)數(shù)據(jù)的擬合結(jié)果基本滿足實際需要,因此提出的PCA結(jié)合DCT方法可以作為一種有效的方法進行復雜強非均質(zhì)性油藏的自動歷史擬合。

圖8 累積指標的擬合結(jié)果Fig.8 History matching results of cumulative index

圖9 單井指標的擬合結(jié)果Fig.9 History matching results of single well index

3 結(jié) 論

(1)將SPSA算法與ABC算法混合對目標函數(shù)進行求解,SPSA算法能夠很好地尋找目標函數(shù)的局部最優(yōu)值,近似梯度與實際梯度方向十分接近使SPSA算法具有很強的收斂性,引入ABC算法的全局搜索模式提高SPSA算法的全局搜索能力。SPSA算法與ABC算法混合提高了模型求解速度與精度。

(2)利用PCA特征提取的方法基本能夠識別出油藏中的高滲通道,同時將高維地質(zhì)參數(shù)數(shù)據(jù)轉(zhuǎn)化為低維參數(shù)使優(yōu)化算法收斂速度更快;但是由于單一的PCA方法不能提取出較為全面的儲層信息,滲透率場反演結(jié)果與真實滲透率場存在較大偏差。

(3)PCA和DCT方法混合能夠提取更為全面的地質(zhì)特征,并且提高歷史擬合的魯棒性,很好地反演出油藏中的高滲通道區(qū)域并且反演結(jié)果更加接近真實場。提出的混合求解方法可以作為一種有效的方法進行該類油藏的自動歷史擬合。