桂金詠, 高建虎, 雍學(xué)善, 李勝軍
1 中國石油勘探開發(fā)研究院西北分院, 蘭州 730020 2 中國石油天然氣集團(tuán)公司油藏描述重點實驗室, 蘭州 730020
?
基于雙相介質(zhì)理論的儲層參數(shù)反演方法
桂金詠1,2, 高建虎1, 雍學(xué)善1, 李勝軍1
1 中國石油勘探開發(fā)研究院西北分院, 蘭州 730020 2 中國石油天然氣集團(tuán)公司油藏描述重點實驗室, 蘭州 730020
傳統(tǒng)基于單相介質(zhì)理論的儲層參數(shù)反演方法將孔隙流體與固體骨架等效為單一固體,弱化了孔隙流體的影響,反演結(jié)果精度不高. 本文提出根據(jù)雙相介質(zhì)理論反演儲層參數(shù)的方法. 首先,在前人研究的基礎(chǔ)上,利用巖石物理模型建立彈性參數(shù)與孔隙度、飽和度、泥質(zhì)含量等儲層參數(shù)間的關(guān)系,進(jìn)而將雙相介質(zhì)反射系數(shù)推導(dǎo)為儲層參數(shù)的函數(shù);其次,根據(jù)貝葉斯反演理論,在高斯噪聲假設(shè)的基礎(chǔ)上,采用更加符合實際情況的修正柯西分布函數(shù)描述反射系數(shù)的稀疏性,推導(dǎo)出儲層物性參數(shù)目標(biāo)反演函數(shù);最后, 應(yīng)用差分進(jìn)化非線性全局尋優(yōu)算法來求解目標(biāo)反演函數(shù),使得反演結(jié)果與實際資料間誤差最小. 新方法旨在突出流體對介質(zhì)反射系數(shù)的影響,以期得到較高的儲層參數(shù)反演精度. 模型與實際資料測試均表明該方法可行、有效且反演精度較高.
雙相介質(zhì); 儲層參數(shù); 反演; 巖石物理模型; 貝葉斯框架
地震波特征分析,尤其是地震波振幅隨偏移距變化(AVO)分析,在儲層預(yù)測中有著非常重要的地位. 不同儲層以及同一種儲層的孔隙流體特征不同時,對地震振幅的影響規(guī)律也不同. 儲層內(nèi)部性質(zhì)如孔隙發(fā)育程度、孔隙流體類型、流體飽和度的改變會對地震振幅造成一定的影響(Hilterman, 1989,2001;Russell et al., 2003). 反映儲層內(nèi)部性質(zhì)的儲層參數(shù),如孔隙度、泥質(zhì)含量、飽和度等參數(shù)是預(yù)測油氣分布、估算油氣儲量、確定開發(fā)井位的重要參數(shù).
傳統(tǒng)的儲層參數(shù)獲取方法將含流體儲層等效為單一固體,利用地震波反射振幅隨偏移距變化特征開展疊前AVO反演,其理論基礎(chǔ)是Zoeppritz方程(Avseth et al.,2005; Grana and Della Rossa, 2010). Zoeppritz方程是描述地震振幅在兩種單相固體介質(zhì)分界面處的反射與透射規(guī)律的基本方程 (Zoeppritz, 1919). 然而,含油氣儲層并非完全固體的彈性介質(zhì). 通常情況下,含油氣儲層為孔隙、裂縫型,表現(xiàn)為固體骨架與流體雙相介質(zhì)特征. 由于流體的存在以及固體與流體的相互作用會弱化巖石的力學(xué)性質(zhì),彈性波在雙相介質(zhì)中的傳播比在單相介質(zhì)中的傳播更具復(fù)雜性,雙相介質(zhì)中不僅有第一類縱波,還有第二類縱波,地震波的反射與透射機(jī)理與單相介質(zhì)存在著較大的差別(Biot, 1956, 1962; 雍學(xué)善等,2006). 因此,傳統(tǒng)基于單相介質(zhì)理論的Zoeppritz方程難以準(zhǔn)確描述地震波在雙相介質(zhì)中的傳播過程和規(guī)律.
Biot理論的提出奠定了雙相介質(zhì)傳播理論的基礎(chǔ). 雙相介質(zhì)理論充分考慮了介質(zhì)的巖石骨架結(jié)構(gòu)和孔隙流體性質(zhì)以及局部特性與整體效應(yīng)的關(guān)系,將含流體儲層表述為固體相和流體相的復(fù)合體,且分別考慮了固體和流體以及二者相互耦合對地震波傳播的影響(Biot, 1956, 1962).Geertsma和Smit(1961)研究了飽和流體孔隙介質(zhì)分界面上縱波垂直入射時的反射透射問題. Nur和Wang(1989, 1999)總結(jié)了早期學(xué)者關(guān)于雙相介質(zhì)的理論、模型與實驗的研究成果,較為詳細(xì)地論述了流體飽和度、裂隙密度、孔隙度、孔隙流體壓力與圍壓、裂隙與孔隙空間幾何形態(tài)等因素對地震波傳播的影響以及孔隙性巖石地震波的傳播特點,并指出了雙相介質(zhì)理論在測井評價、強(qiáng)化回收率、斷層檢測、油氣區(qū)域圈定等領(lǐng)域的應(yīng)用前景. 王尚旭(1990)研究了雙相介質(zhì)中地震波的傳播規(guī)律,實現(xiàn)了雙相介質(zhì)中地震波傳播的有限元解法,并對單相和雙相介質(zhì)反射系數(shù)進(jìn)行了比較,指出當(dāng)雙相介質(zhì)孔隙度較大且孔隙流體的彈性模量和密度與固體顆粒相同量級時,兩者有明顯差異. Wu等(1990)和喬文孝等(1992)對彈性波在雙相介質(zhì)分界面處任意入射角下的反射、透射特征進(jìn)行了研究. 牟永光(1996)通過數(shù)值模擬實驗,觀測到了雙相介質(zhì)中慢速縱波和慢速橫波的存在并給出了雙相介質(zhì)AVO方程的具體表達(dá)式. 楊頂輝(2002)從微觀流場的流體相對流動各向異性速度和Biot宏觀流體的雙相各向異性介質(zhì)中彈性波方程出發(fā),實現(xiàn)了雙相各向異性介質(zhì)中“相對流體位移”的有限元波場模擬. 雍學(xué)善等(2006),高建虎等(2007)對雙相介質(zhì)AVO方程進(jìn)行了參數(shù)簡化及反演研究,為雙相介質(zhì)地震波傳播理論的實際應(yīng)用開辟了新的途徑. 王華青和田家勇(2011)從巖土力學(xué)、地球物理勘探、聲學(xué)研究等領(lǐng)域方面總結(jié)了雙相介質(zhì)彈性波傳播問題在理論、動態(tài)參數(shù)的測量技術(shù)、動力響應(yīng)分析方面的研究成果.這些成果的獲得為我們進(jìn)一步研究雙相介質(zhì)儲層參數(shù)反演提供了堅實的基礎(chǔ).
本文直接從雙相介質(zhì)地震波反射理論出發(fā),將雙相介質(zhì)反射系數(shù)簡化為儲層參數(shù)的函數(shù),并建立了新的目標(biāo)反演函數(shù)及非線性求解算法.通過模型反射系數(shù)敏感性分析、目標(biāo)函數(shù)收斂性分析以及實際資料測試,表明了方法的可行性與有效性.
2.1 雙相介質(zhì)反射系數(shù)方程
在兩個半無限的雙相介質(zhì)分界面處,平面縱波Pi以入射角αi入射到分界面,將產(chǎn)生以下幾類波:第一類反射縱波P11,第二類反射縱波P12,反射橫波S1,第一類透射縱波P21,第二類透射縱波P22,透射橫波S2,如圖1所示.
圖1 兩種介質(zhì)分界面上的反射與透射Fig.1 Reflection and transmission on the surface of the boundary between two media
定義快縱波反射系數(shù)RP11=AP11/APi,慢縱波反射系數(shù)RP12=AP12/APi,橫波反射系數(shù)RS1=AS1/APi,快縱波透射系數(shù)TP21=AP21/APi,慢縱波透射系數(shù)TP22=AP22/APi,橫波透射系數(shù)TS2=AS2/APi. 在平面波傳播斯奈爾定律以及固體法向、切向位移連續(xù)、法向總應(yīng)力連續(xù)、固體切向應(yīng)力連續(xù)、流體壓力及流體法向流量連續(xù)等六個邊界條件的假設(shè)基礎(chǔ)上,根據(jù)牟永光(1996)的推導(dǎo)結(jié)果,有:
(1)
由方程(1)可以看到,雙相介質(zhì)反射系數(shù)和透射系數(shù)方程并不像Zoeppritz方程那樣簡單,用縱波速度、橫波速度和密度等3個彈性參數(shù)就可以表達(dá),而需用10多個參數(shù)才能表達(dá). 其中,包括4個雙相介質(zhì)彈性參數(shù)(A、N、Q、R)、3個質(zhì)量參數(shù)(ρ11、ρ22、ρ12)、2個流體與固體振幅比參數(shù)(m1、m2)以及孔隙度(φ). 如此多的參數(shù)很難由地震資料反演求取且與常用的描述儲層內(nèi)部性質(zhì)的含氣飽和度、孔隙度等儲層參數(shù)差異較大,在實際工作中很難使用. 需要進(jìn)一步推導(dǎo),將雙相介質(zhì)反射系數(shù)方程改寫為儲層參數(shù)的函數(shù).
(2)
其中Km、Kd、Kf分別為基質(zhì)、干巖、混合流體的體積模量;μd為干巖剪切模量;ρd、ρf分別為干巖、流體的體積密度;α為曲折度,可以利用Berryman關(guān)系式α=1-r(1-1/φ)求取,r介于0和1之間(Berryman et al.,1995).
重點在于未知彈性參數(shù)Km、Kd、Kf、μd、ρd、ρf的求取. 巖石物理模型是建立彈性參數(shù)與儲層參數(shù)間定量關(guān)系的橋梁. 以砂巖天然氣儲層為例,通過多種巖石物理模型可將彈性參數(shù)E描述為含氣飽和度Sg、孔隙度φ以及泥質(zhì)含量Vsh的函數(shù):E=f(Sg,φ,Vsh).因此,可根據(jù)實際研究區(qū)地層巖石物理特點選擇合適的巖石物理模型將雙相介質(zhì)參數(shù)進(jìn)一步轉(zhuǎn)化為儲層參數(shù)的函數(shù).
在砂巖天然氣儲層中,飽和巖石一般由粘土、石英、水、氣等成分組成. 記Kc、Ks、Kw、Kg分別表示粘土、石英、水、氣的體積模量,μc、μs分別表示粘土、石英的剪切模量;ρc、ρs、ρw、ρg分別表示粘土、石英、水、氣的體積密度. 一般采用Voight-Ruess-Hill平均求取巖石基質(zhì)模量(Mavko et al., 2009):
(3)
不同巖石物理模型對巖石孔隙結(jié)構(gòu)的描述不同,主要體現(xiàn)在對干巖模量的求取方式上,具有代表性的有Krief模型、Nur模型、Pride模型以及Hertz-Mindlin模型等(Nur et al.,1998; Mavko et al., 2009). 對于一般干燥巖石,可根據(jù)Nur的臨界孔隙度理論,利用臨界孔隙度φc將孔隙介質(zhì)的機(jī)械和聲波特征分成兩個截然不同的區(qū)域,大于臨界孔隙度為懸浮域,小于臨界孔隙度為承載域,可以得到在承載域下的干巖體積、剪切、密度模量(Nuretal.,1998):
(4)
多種流體共存的混合模式主要有均勻分布和斑塊飽和分布(White, 1975). 在天然氣儲層中一般傾向于斑塊飽和分布(MavkoandMukerji,1998). 通過合理調(diào)節(jié)斑塊參數(shù)e可以得到流體體積模量Kf與含氣飽和度間的關(guān)系(Brie et al., 1995):
Kf=(Kw-Kg)(1-Sg)e+Kg.
(5)
由Wood公式建立含氣飽和度Sg同流體密度ρf的關(guān)系(Mavko et al., 2009):
ρf=Sgρg+(1-Sg)ρw.
(6)
最后,根據(jù)Gassmann理論,可以得到飽和巖石彈性參數(shù)與儲層參數(shù)間關(guān)系(Mavkoetal., 2009):
(7)
式(3)—(6)中飽和巖石組分模量參數(shù)一般通過巖石物理實驗或經(jīng)驗統(tǒng)計得到. 表1給出了本文采用的飽和巖石組分模量參數(shù). 在通過(3)—(7)式獲得飽和巖石彈性參數(shù)后,代入(2)式得到雙相介質(zhì)參數(shù). 最后可由雍學(xué)善等(2006)的方法得到流體與固體振幅比參數(shù)(m1、m2).
表1 飽和巖石組分模量參數(shù)Table 1 Modulus value of saturated rock components
至此,通過合理的巖石物理模型假設(shè),可以將雙相介質(zhì)參數(shù)轉(zhuǎn)化為儲層參數(shù)的函數(shù),進(jìn)而代入(1)式,雙相介質(zhì)反射系數(shù)方程便可以轉(zhuǎn)化為飽和度Sg、孔隙度φ以及泥質(zhì)含量Vsh的函數(shù).雙相介質(zhì)反射系數(shù)方程簡化為三個未知參數(shù)的函數(shù),較原反射系數(shù)方程未知參數(shù)大大減小,更加有利于疊前AVO反演. 通過改寫后的雙相介質(zhì)反射系數(shù)方程,反射系數(shù)與飽和度、孔隙度及泥質(zhì)含量等儲層參數(shù)間建立起了直接的關(guān)系,使得雙相介質(zhì)儲層參數(shù)的反演成為可能. 然而,這種關(guān)系為高度非線性關(guān)系,需要進(jìn)一步建立一套非線性反演機(jī)制才能從地震振幅中得到準(zhǔn)確的儲層參數(shù).
2.2 目標(biāo)反演函數(shù)
基于貝葉斯反演理論框架,構(gòu)建目標(biāo)反演函數(shù). 按照地震褶積模型:
d=Gr+n,
(8)
其中d=[d1,d2,…,dN]T是觀測得到的地震數(shù)據(jù),r=[r1,r2,…,rM]T是反射系數(shù)序列,G是N×M維子波矩陣,n=[n1,n2,…,nN]T表示觀測噪聲.由貝葉斯公式可得以下近似式:
P(r|d)∝P(r)·P(d|r),
(9)
其中:P(r|d)表示反射系數(shù)的后驗概率密度信息,P(r)表示反射系數(shù)先驗概率密度信息,P(d|r)表示似然函數(shù).
一般地面地震觀測數(shù)據(jù)的噪聲服從零均值、協(xié)方差矩陣為CS的高斯分布(Downton, 2005).當(dāng)給定反射系數(shù)r之后,則觀測數(shù)據(jù)與模型參數(shù)之間的似然函數(shù)P(d|r)可以通過該觀測方式下噪聲的分布特征來表示,即P(d|r)=
(10)
(11)
根據(jù)貝葉斯反演框架,后驗概率分布P(r|d)表達(dá)了從觀察數(shù)據(jù)中獲得的對未知量的認(rèn)知程度:
(12)
由貝葉斯理論可知,使(12)式取最大值的解為反射系數(shù)的最優(yōu)解,即對模型參數(shù)的最大后驗概率密度解.
對(12)式兩邊同時取對數(shù),整理后即可得到目標(biāo)反演函數(shù):
(13)
若任給一組飽和度、孔隙度、泥質(zhì)含量參數(shù),即(Sg,φ,Vsh)i,就可以得到一條與之相應(yīng)的反射系數(shù)曲線ri. 因此,反演問題就變成了通過選擇合適的一組參數(shù)(Sg,φ,Vsh)i,使得目標(biāo)函數(shù)J(ri)為最小的過程.
2.3 非線性反演
考慮到求解精度與效率,本文采用差分進(jìn)化(Differential Evolution,DE)算法對目標(biāo)反演函數(shù)(13)式進(jìn)行非線性求解. 同遺傳算法一樣,DE算法也是基于群體智能理論的優(yōu)化算法,通過群體內(nèi)個體間的合作與競爭產(chǎn)生的群體智能指導(dǎo)優(yōu)化搜索. 但與傳統(tǒng)遺傳算法相比,其采用實數(shù)編碼、基于差分的簡單變異操作和一對一的競爭生存策略,降低了遺傳操作的復(fù)雜性且易于并行實現(xiàn),效率較高且不需要借助問題的特征信息,適于求解一些利用常規(guī)的數(shù)學(xué)規(guī)劃方法所無法求解的復(fù)雜優(yōu)化問題(Storn and Price, 1997).
設(shè)初始種群規(guī)模是由N個D維參數(shù)向量表示的個體xi=[xi,1,xi,2,…,xi,D]組成,則該種群通過變異、交叉、選擇等三個過程即可得到下一代種群 (StornandPrice, 1997; 印興耀等,2013):
(1)變異操作
(2)交叉過程
其中,RAN(j)∈[0,1]代表均勻分布隨機(jī)數(shù);C為交叉概率常數(shù);RANn(i)∈[1,2,…,D].
(3)選擇標(biāo)準(zhǔn)
3.1 反射系數(shù)敏感性分析
疊前地震資料振幅隨偏移距的變化實質(zhì)是反射系數(shù)隨入射角的變化. 不同儲層參數(shù)值對應(yīng)的反射系數(shù)值差異越明顯,表示反射系數(shù)對儲層參數(shù)的敏感性越強(qiáng),利用反射系數(shù)反演該儲層參數(shù)的可能性也就越大. 設(shè)斑塊飽和參數(shù)e為3,依據(jù)表1中巖石物理參數(shù)和臨界孔隙度為0.4的Nur巖石物理模型以及(1)—(7)式(如無特別說明,文中分析實例的巖石物理模型及參數(shù)均同此),建立一個兩層泥、砂巖模型用來分析在不同入射角情況下含氣飽和度、孔隙度、泥質(zhì)含量的變化對雙相介質(zhì)反射系數(shù)的影響. 模型儲層參數(shù)及變化范圍如表2所示. 固定孔隙度、泥質(zhì)含量及含氣飽和度三參數(shù)中的兩個,另一參數(shù)在某一范圍內(nèi)變化,相應(yīng)的反射系數(shù)隨入射角變化情況如圖2所示.
表2 敏感性分析模型參數(shù)Table 2 Model parameters for sensitivity analysis
由圖2可以看到,含氣飽和度、孔隙度及泥質(zhì)含量的改變都會對反射系數(shù)的變化造成一定的影響且影響程度與入射角有關(guān).
在泥質(zhì)含量與孔隙度固定的情況下,可以發(fā)現(xiàn)隨著含氣飽和度的增加,反射系數(shù)間的差異逐漸減小,這種趨勢隨著入射角的增大越發(fā)明顯. 在入射角小于35°的范圍內(nèi),不同含氣飽和度值對應(yīng)的反射系數(shù)值差異較大,易于區(qū)分.當(dāng)入射角大于35°時,隨著入射角的增大,高含氣飽和度值對應(yīng)的反射系數(shù)值差異越來越小.
在含氣飽和度與孔隙度固定的情況下,模型的泥質(zhì)含量設(shè)定在砂巖巖性范圍內(nèi)0~0.5. 可以發(fā)現(xiàn)低泥質(zhì)含量(小于0.3)對應(yīng)的反射系數(shù)值在入射角小于35°的情況下,差異明顯.入射角大于35°時,差異逐漸變小,直至重合.而高泥質(zhì)含量(大于0.3)對應(yīng)的反射系數(shù)值在入射角臨近50°時,仍有明顯的差異.
在含氣飽和度與泥質(zhì)含量固定的情況下,可以發(fā)現(xiàn)不同孔隙度對應(yīng)的反射系數(shù)值差異較含氣飽和度、泥質(zhì)含量尤為明顯,且在50°后仍有一定差異.
整體上在泥巖、砂巖分界面處,雙相介質(zhì)反射系數(shù)對含氣飽和度、泥質(zhì)含量、孔隙度變化的敏感程度依次增強(qiáng). 當(dāng)入射角大于臨界角時,反射系數(shù)會發(fā)生劇烈變化,對儲層參數(shù)的敏感性變得十分復(fù)雜.
圖2 反射系數(shù)隨儲層參數(shù)及入射角的變化(色標(biāo)表示反射系數(shù)值),(a)、(b)、(c)分別為反射系數(shù)隨含氣飽和度、泥質(zhì)含量、孔隙度及入射角的變化Fig.2 Reflection coefficient versus reservoir parameters and incident angles (the color codes represent the reflection coefficient values): (a) reflection coefficient versus gas saturation and incident angles, (b) reflection coefficient versus shale content and incident angles, (c) reflection coefficient versus porosity and incident angles
圖3 目標(biāo)函數(shù)收斂情況Fig.3 The convergence of objective inversion function
實際地震資料的入射角范圍一般小于35°,在臨界角以內(nèi). 因此,在實際應(yīng)用中,利用振幅隨入射角的變化特征來反演儲層參數(shù)是可行的.
3.2 目標(biāo)函數(shù)收斂性分析
通過反射系數(shù)敏感性分析可以確認(rèn)利用振幅隨入射角的變化反演儲層參數(shù)是可行的,但反演目標(biāo)函數(shù)的收斂性關(guān)系到能否具備全局最優(yōu)解. 建立一個兩層砂泥巖模型,通過給定砂巖儲層參數(shù)取值范圍-50%~50%,設(shè)定入射角度為0°~40°,利用目標(biāo)反演函數(shù)(13)式分析雙相介質(zhì)反演目標(biāo)函數(shù)能否收斂到極小值,模型參數(shù)如表3所示.
表3 收斂性分析模型參數(shù)Table 3 Model parameter for convergence analysis
由圖3a—3c可以看到,目標(biāo)函數(shù)收斂域形態(tài)呈橢球狀,橢球的短軸方向較長軸方向易收斂. 整體上看,目標(biāo)函數(shù)對含氣飽和度、泥質(zhì)含量、孔隙度的收斂程度雖有所差別,但均能收斂到極小值附近.
固定含氣飽和度,分析孔隙度、泥質(zhì)含量取值與目標(biāo)函數(shù)的關(guān)系,由圖3a可以看到,目標(biāo)函數(shù)沿縱向變化快,沿橫向變化慢,說明砂巖儲層中孔隙度較泥質(zhì)含量易收斂.
固定泥質(zhì)含量,分析孔隙度、含氣飽和度與目標(biāo)函數(shù)的關(guān)系,由圖3b可以看到,目標(biāo)函數(shù)沿縱向變化快,沿橫向變化慢,說明砂巖儲層中孔隙度較含氣飽和度易收斂.
固定孔隙度,分析泥質(zhì)含量、含氣飽和度與目標(biāo)函數(shù)的關(guān)系,由圖3c可以看到,目標(biāo)函數(shù)沿橫向變化快,沿縱向變化慢,說明砂巖儲層中泥質(zhì)含量較含氣飽和度易收斂.
通過反射系數(shù)敏感性分析以及目標(biāo)函數(shù)收斂性分析,可以得出這樣的結(jié)論:在實際疊前地震資料入射角范圍內(nèi),反射系數(shù)整體上對孔隙度、泥質(zhì)含量、含氣飽和度的變化較為敏感且目標(biāo)函數(shù)能收斂到極小值,利用疊前地震資料反演儲層參數(shù)是可行的且能得到最優(yōu)解.
4.1 模型測試
依據(jù)表1中巖石物理參數(shù)和臨界孔隙度為0.4的Nur巖石物理模型,利用A區(qū)測井曲線計算得到的雙相介質(zhì)反射系數(shù)與40 Hz零相位雷克子波褶積,合成角度范圍為0°~35°的疊前地震角度道集(每道間隔5°),如圖4所示. 以此為基礎(chǔ)數(shù)據(jù)對不同初始模型、信噪比條件下的模型數(shù)據(jù)進(jìn)行反演結(jié)果分析.
反演過程中,DE算法的縮放控制參數(shù)F設(shè)定為0.8,交叉概率C設(shè)定為0.4,迭代次數(shù)設(shè)定為200次. 以長度為10個樣點的小時窗對測井曲線進(jìn)行平滑作為初始模型,反演結(jié)果如圖5所示. 可見,初始模型與測井曲線十分接近,反演結(jié)果亦與測井記錄非常吻合. 以長度為50個樣點的大時窗對測井曲線進(jìn)行平滑作為初始模型,反演結(jié)果如圖6所示. 圖6可見初始模型雖與測井曲線差異較大,但仍能得到較高精度的儲層參數(shù)反演結(jié)果. 這得益于DE算法對初始模型的依賴程度較弱,能有效避免陷入局部極值的優(yōu)點(印興耀等,2013).
在實際地球物理勘探中,受采集環(huán)境及儀器影響,采集到的疊前地震數(shù)據(jù)不可避免的含有一定程度的噪聲. 為了使模型分析結(jié)果更具實際指導(dǎo)意義,在合成地震記錄中加入隨機(jī)噪聲,得到信噪比分別為10、5的含噪地震記錄來分析噪聲對本文方法的影響. 含噪合成地震記錄如圖7所示.
采用與圖6相同的初始模型,利用圖7中所示含噪合成地震角度道集反演孔隙度、泥質(zhì)含量及含氣飽和度,反演結(jié)果如圖8、圖9所示.
通過對比圖6、圖8、圖9可以看到,當(dāng)加入噪聲時,孔隙度、泥質(zhì)含量及含氣飽和度反演結(jié)果在整體上均受到一定的影響,信噪比越低,影響越大. 同時也注意到,即使信噪比為5,反演結(jié)果雖有所偏離,但整體趨勢仍然接近實際曲線,主要原因在于目標(biāo)反演函數(shù)在構(gòu)建的過程中考慮了噪聲的影響. 含氣飽和度受噪聲影響相對于孔隙度、泥質(zhì)含量稍大,其他兩個儲層參數(shù)在反演過程中對噪聲具有一定的壓制,因為通過前面的分析可以知道雙相介質(zhì)反射系數(shù)本身對含氣飽和度的敏感性相對于孔隙度、泥質(zhì)含量較低. 因此,為確保含氣飽和度反演結(jié)果更加可靠,反演前需對疊前道集進(jìn)行相應(yīng)的疊前去噪處理.
4.2 實際資料測試
為了進(jìn)一步說明雙相介質(zhì)儲層參數(shù)反演方法的實用性,以實際測井、井旁道疊前地震道集進(jìn)行反演效果測試. 在進(jìn)行反演前,疊前道集已經(jīng)過去噪、保幅等處理,信噪比較高. 該測試數(shù)據(jù)孔隙度、泥質(zhì)含量、含氣飽和度測井曲線以及疊前道集如圖10所示,在630 ms左右處發(fā)育一套砂巖氣藏.
圖4 測井曲線及合成地震記錄Fig.4 The well log curves and synthetic seismic record
選用恰當(dāng)?shù)膸r石物理模型能保障儲層參數(shù)轉(zhuǎn)化到反射系數(shù)的正確性.利用圖10中孔隙度、泥質(zhì)含量、含氣飽和度曲線,分別根據(jù)Nur與Krief的巖石物理模型計算彈性參數(shù)曲線,與實測彈性參數(shù)曲線對比結(jié)果如圖11所示. 兩種巖石物理模型的主要差異在于對干巖體積模量與剪切模量的計算方式不同,進(jìn)而導(dǎo)致縱、橫波速度的計算結(jié)果不同. 由圖11可見,兩種巖石物理模型計算的縱、橫波速度差異較大,而密度相同,整體上Nur的巖石物理模型建模結(jié)果與實測曲線對應(yīng)較好.
高信噪比的疊前地震資料以及恰當(dāng)?shù)膸r石物理模型,說明該區(qū)域適合進(jìn)行雙相介質(zhì)儲層參數(shù)反演,儲層參數(shù)反演結(jié)果如圖12所示.
由圖12可以看到,儲層參數(shù)反演結(jié)果與實測結(jié)果吻合較好,630 ms處的氣層特征得到了較好的反映:高孔隙、低泥質(zhì)含量、高含氣飽和度.
將反演過程中的地震子波與由反演成果曲線計算得到的雙相介質(zhì)反射系數(shù)進(jìn)行褶積運(yùn)算,得到合成角度道集,分析合成角度道集與實際道集的吻合程度,如圖13所示.由圖13可以看到,利用反演結(jié)果合成得到的角度道集與實際地震角度道集在整體上對應(yīng)較好,振幅趨勢一致,殘差較小. 這說明反演結(jié)果是準(zhǔn)確、可靠的.
受計算效率及成本的限制,非線性反演方法目前仍停留在二維實際資料應(yīng)用階段. 一般而言,商業(yè)氣藏與貧氣藏在地震資料上具有相似的反射特征,如“亮點”等,一般方法較難區(qū)分(Avseth et al., 2005). 然而,地震解釋工作者仍可以根據(jù)含氣飽和度的高低來確定氣藏是否具有商業(yè)價值. 因此,準(zhǔn)確的含氣飽和度參數(shù)顯得尤為重要. 圖14為含氣飽和度二維反演剖面,可以看到含氣飽和度反演剖面上的高值部分與實際含氣位置(紅色橢圓位置)對應(yīng)較好. 根據(jù)該反演剖面,地震解釋工作者不僅可以確定含氣儲層的分布,還可以根據(jù)含氣飽和度值的高低判斷該儲層是否具有商業(yè)價值.
圖5 小時窗平滑初始模型儲層參數(shù)反演結(jié)果黑色虛線、藍(lán)色實線、紅色實線分別代表初始模型、測井?dāng)?shù)據(jù)、反演結(jié)果曲線(下同).Fig.5 Inversion results of reservoir parameters with initial model of small time window Black dash line represents the initial model curve, blue line represents the actual value, red line represents the inverted value (the same below).
圖6 大時窗平滑初始模型儲層參數(shù)反演結(jié)果Fig.6 Inversion results of reservoir parameters with initial model of big time window
圖7 含噪合成地震角度道集(a) 信噪比為10; (b) 信噪比為5.Fig.7 Synthetic seismic record with noise in the angle domain(a) Signal-noise ratio is 10; (b) Signal-noise ratio is 5.
圖8 地震數(shù)據(jù)信噪比為10時的儲層參數(shù)反演結(jié)果實線、虛線分別為實際測井、反演曲線(下同).Fig.8 Inversion results of reservoir parameters with signal-noise ratio is 10 Solid line represents the actual value, dash line represents the inverted value (the same below).
圖9 地震數(shù)據(jù)信噪比為5時的儲層參數(shù)反演結(jié)果Fig.9 Inversion results of reservoir parameters with signal-noise ratio is 5
圖10 實際資料儲層參數(shù)測井曲線及井旁道道集Fig.10 Actual reservoir parameters well log curves and the borehole-side seismic gathers
圖11 巖石物理模型彈性參數(shù)曲線黑色、紅色、藍(lán)色線分別表示實際測井曲線、Nur模型計算曲線、Krief模型計算曲線.Fig.11 Elastic parameter curves calculated by rock physical model Black line,red line and blue line respectively represent the actual curve, curve calculated by Nur′s model and curve calculated by Krief′s model.
儲層參數(shù)能較好地表征儲層的物性特征,評價儲層的開發(fā)潛力. 當(dāng)巖石含有流體時,表現(xiàn)為雙相介質(zhì)特征,常規(guī)基于單相介質(zhì)理論的儲層參數(shù)反演方法弱化了流體的影響,預(yù)測油氣存在不足. 本文提出的儲層參數(shù)反演方法,直接基于雙相介質(zhì)理論,突出流體的影響,利用雙相介質(zhì)振幅隨偏移距變化規(guī)律反演儲層參數(shù). 與常規(guī)方法相比,該方法更加貼合儲層實際情況. 模型與實際資料測試結(jié)果均表明,該方法具有較好的反演效果. 在實際應(yīng)用過程中有幾個方面需要特別注意:
(1) 本文方法受巖石物理模型的影響因素較大,不恰當(dāng)?shù)慕=Y(jié)果會導(dǎo)致儲層參數(shù)與地震振幅間的聯(lián)系發(fā)生嚴(yán)重偏離. 巖石物理實驗以及大量的巖石物理模型統(tǒng)計分析能有效地保證建模結(jié)果的準(zhǔn)確性;
(2) 同常規(guī)疊前地震反演方法相同,高信噪比以及保幅的疊前地震道集能保證反演結(jié)果的精度;
(3) DE算法雖然對初始模型的依賴較弱,但仍建議利用常規(guī)方法的反演結(jié)果作為初始模型,能夠減少反演迭代次數(shù),節(jié)省計算時間;
(4) 目標(biāo)函數(shù)的求解具有高度非線性,涉及大量、反復(fù)的正演過程,計算量較大,高性能計算設(shè)備的投入以及高效的正、反演算法的發(fā)展能夠有力推動雙相介質(zhì)理論走向大規(guī)模三維實際應(yīng)用.
圖12 儲層參數(shù)反演結(jié)果實線、虛線分別為實際測井、反演曲線.Fig.12 Inversion results of reservoir parameters Solid line represents the actual value, dash line represents the inverted value.
圖13 反演結(jié)果合成記錄與實際地震記錄對比(a) 黑色、紅色曲線分別為實際地震道、合成地震道; (b) 為合成地震道與實際地震道間的殘差.Fig.13 The comparison between synthetic seismic and actual seismic record(a) Black curve and red curve respectively represent for the actual seismic and synthetic seismic record; (b) The errors between synthetic seismic and actual seismic record.
圖14 含氣飽和度反演剖面Fig.14 The inverted gas saturation profile
Avseth P, Mukerji T, Mavko G. 2005. Quantitative Seismic Interpretation: Applying Rock Physics Tools to Reduce Interpretation Risk. Cambridge: Cambridge University Press.
Berryman J G, Wang H F. 1995. The elastic coefficients of double-porosity models for fluid transport in jointed rock.JournalofGeophysicalResearchSolidEarth, 100(B12): 24611-24627.
Biot M A. 1956. Theory of propagation of elastic waves in a fluid-saturated porous solid. I. low-frequency range.J.Acoust.Soc.Am., 28(2): 168-178.Biot M A. 1962. Mechanics of deformation and acoustic propagation in porous media.JournalofAppliedPhysics, 33(4): 1482-1498.
Brie A, Pampuri F, Marsala A F, et al. 1995. Shear sonic interpretation in gas-bearing sands. ∥SPE Annual Technical Conference & Exhibition. Dallas, Texas: SPE, 701-710.
Downton J E. 2005. Seismic parameter estimation from AVO inversion [Ph. D. thesis]. Calgary: University of Calgary.
Gao J H, Liu Q X, Yong X S, et al. 2007. The method study of pre-stack reservoir parameter inversion in dual phase media.AdvancesinEarthScience(in Chinese), 22(10): 1048-1054.
Geertsma J, Smit D C. 1961. Some aspects of elastic wave propagation in fluid-saturated porous solids.Geophysics, 26(2): 169-181.
Grana D, Della R E. 2010. Probabilistic petrophysical-properties estimation integrating statistical rock physics with seismic inversion.Geophysics, 75(3): O21-O37.
Hilterman F J. 1989. Is AVO the seismic signature of rock properties. ∥ 59th Annual Meeting, SEG, Expanded Abstracts. Dallas, Texas, 59.Hilterman F J. 2001. Seismic amplitude interpretation: 2001 distinguished instructor short course (No.4). Houston: Geophysical Development Corporation.Mavko G, Mukerji T. 1998. Bounds on low-frequency seismic velocities in partially saturated rocks.Geophysics, 63(3): 918-924. Mavko G, Mukerji T, Dvorkin J. 2009. The Rock Physics Handbook: Tools for Seismic Analysis of Porous Media. Cambridge: Cambridge University Press. Mou Y G. 1996. Reservoir Geophysics (in Chinese). Beijing: Petroleum Industry Press.Nur A M, Wang Z J. 1989. Seismic and acoustic velocities in reservoir rocks: experimental studies. Geophysics reprint series. Society of Exploration Geophysicists.Nur A M, Mavko G, Dvorkin J, et al. 1998. Critical porosity: A key to relating physical properties to porosity in rocks.TheLeadingEdge, 17(3): 357-362.
Nur A M, Wang Z J. 1999. Seismic and Acoustic Velocities in Reservoir Rocks: Recent Developments (Vol.3). Soc of Exploration Geophysicists.Qiao W X, Wang Y, Yan Z P. 1992. Reflection and transmission of acoustic wave at a porous solid/porous solid interface.ChineseJ.Geophys. (ActaGeophysicaSinica) (in Chinese), 35(2): 242-248.
Russell B H, Hedlin K, Hilterman F J, et al. 2003. Fluid-property discrimination with AVO: A Biot-Gassmann perspective.Geophysics, 68(1): 29-39.Storn R, Price K. 1997. Differential evolution—a simple and efficient heuristic for global optimization over continuous spaces.JournalofGlobalOptimization, 11(4): 341-359.
Wang S X. 1990. Finite element numerical solution and AVO problem of elastic wave in two phase media [Ph. D. thesis] (in Chinese). Beijing: University of Petroleum.
Wang H Q, Tian J Y. 2011. Review of the investigation of elastic-wave propagation in a fluid-saturated porous solid. ∥ Bulletin of the Institute of Crustal Dynamics (in Chinese). Beijing: Seismological Press: 28-34.
White J E. 1975. Computed seismic speeds and attenuation in rocks with partial gas saturation.Geophysics, 40(2): 224-232.
Wu K, Xue Q, Adler L.1990. Reflection and transmission of elastic waves from a fluid-saturated porous solid boundary.JournaloftheAcousticalSocietyofAmerica, 87:2349-2358.
Yang D H. 2002. Finite element method of the elastic wave equation and wave-field simulation in two-phase anisotropic media.ChineseJ.Geophys. (in Chinese), 45(4): 575-583.
Yin X Y, Kong X X, Zhang F C, et al. 2013. Pre-stack AVO inversion based on the differential evolution algorithm.OilGeophysicalProspecting(in Chinese), 48(4): 591-596.
Yong X S, Ma H Z, Gao J H. 2006. A study of AVO equation in dual-phase medium and parameter simplification.AdvancesinEarthScience(in Chinese), 21(3): 242-249.
Youzwishen C F. 2001. Non-linear sparse and blocky constraints for seismic inverse problems [Ph. D. thesis]. Alberta: University of Alberta.
Zoeppritz K. 1919. On the reflection and penetration of seismic waves through unstable layers.G?ttingerNachr, 1(7B): 66-84.
附中文參考文獻(xiàn)
高建虎, 劉全新, 雍學(xué)善等. 2007. 雙相介質(zhì)疊前儲層參數(shù)反演方法研究. 地球科學(xué)進(jìn)展, 22(10): 1048-1054.
牟永光. 1996. 儲層地球物理學(xué). 北京: 石油工業(yè)出版社.
喬文孝, 王宇, 嚴(yán)熾培. 1992. 聲波在兩種多孔介質(zhì)界面上的反射和透射. 地球物理學(xué)報, 35(2): 242-248.
王尚旭. 1990. 雙相介質(zhì)中彈性波問題有限元數(shù)值解和AVO問題[博士論文]. 北京: 石油大學(xué).
王華青, 田家勇. 2011. 流體飽和多孔介質(zhì)中彈性波傳播問題的研究進(jìn)展. ∥ 地殼構(gòu)造與地殼應(yīng)力文集. 北京: 地震出版社, 28-34.楊頂輝. 2002. 雙相各向異性介質(zhì)中彈性波方程的有限元解法及波場模擬. 地球物理學(xué)報, 45(4): 575-583.
印興耀, 孔栓栓, 張繁昌等. 2013. 基于差分進(jìn)化算法的疊前AVO反演. 石油地球物理勘探, 48(4): 591-596.
雍學(xué)善, 馬海珍, 高建虎. 2006. 雙相介質(zhì)AVO方程及參數(shù)簡化研究. 地球科學(xué)進(jìn)展, 21(3): 242-249.
(本文編輯 何燕)
Inversion of reservoir parameters based on dual-phase media theory
GUI Jin-Yong1,2, GAO Jian-Hu1, YONG Xue-Shan1, LI Sheng-Jun1
1PetroChinaResearchInstituteofPetroleumExploration&Development-NorthWest,Lanzhou730020,China2CNPCReservoirDescriptionKeyLaboratory,Lanzhou730020,China
Traditional inversion methods of reservoir parameters based on the one-phase media theory regard the medium a single solid-phase body of solid skeleton and pore fluid, which may weaken the effect of pore fluid and lead to a lower credibility in inversion results. To make full use of fluid information in seismic amplitude and improve the precision of reservoir parameter prediction, this paper suggests a new method to invert reservoir parameters based on dual-phase media theory.Firstly, on the basis of previous work, a rock physics model is used to establish the relationship between the elastic parameter and reservoir parameter. Then a new dual-phase reflection coefficient equation is derived, which is the function of reservoir parameters. Secondly, according to Bayesian inversion theory and the Gaussian noise assumption, the modified Cauchy distribution is employed to describe the sparsity of reflection coefficients and deduce the final objective inversion function. Finally the differential evolution (DE) nonlinear global optimization algorithm is adopted to solve the objective inversion function and minimize the error between the inversion results and real data.From the reflection coefficient sensitivity analysis, the sensitivities to the reflection coefficient caused by change values of gas saturation, shale content and porosity are different. While fixing the values of shale content and porosity, we found that with the increase of gas saturation, the difference between reflection coefficient decreases. This trend is more obvious with the increasing incident angle. While fixing the values of gas saturation and porosity, the variation of reflection coefficients with low shale content values (less than 0.3) and under the condition of the incident angle less than 35° are obvious. With fixed values of shale content and gas saturation, the difference of reflection coefficients changing with porosity is more obvious than gas saturation and shale content, and there are still some difference even after 50°. From the objective function convergence analysis, the shape of convergence domain likes an ellipsoid. The porosity is easiest to converge to minimum and the shale content comes second. The convergence rate changes with reservoir parameters, but they all can converge to minimum. According to the analysis of the model, the DE algorithm used to solve the objective function is not sensitive to the initial model and is immune to noise to a certain degree. Even the initial model has great deviation and the signal-noise ratio is 5, our method still has a satisfactory result. We think that reservoir parameters are sensitive to the reflection coefficient and objective function can converge to a minimum in the range of actual incident angle. Using the dual-phase media theory to conduct reservoir parameter inversion is feasible and can achieve the optimal solution.Reservoir parameters can be quantitatively described as the reservoir characteristics. Results show that the inversion method based on dual-phase media theory can highlight the influence of fluid-phase and the method is feasible, effective and accurate. To get better actual application results, there are some problems need to be noted: (1) Inappropriate modeling results can cause a serious deviation to the relationship between reservoir parameters and seismic amplitude. Rock-physics experiments and statistical analysis can effectively ensure the accuracy of the modeling results. (2) High signal to noise ratio and amplitude-preserved seismic data can guarantee the accuracy of the inversion results. (3) The DE algorithm has weak dependence on the initial model, while using inversion results from traditional methods as the initial model can reduce iteration times and computing consumption of the DE algorithm. (4) High-performance computation equipment and high efficiency of forward and inversion algorithms can vigorously promote the application of the dual-phase media theory to large-scale 3D issues.
Dual-phase media; Reservoir parameter; Inversion; Rock physics model; Bayesian framework
桂金詠, 高建虎, 雍學(xué)善等. 2015. 基于雙相介質(zhì)理論的儲層參數(shù)反演方法.地球物理學(xué)報,58(9):3424-3438,
10.6038/cjg20150934.
Gui J Y, Gao J H, Yong X S, et al. 2015. Inversion of reservoir parameters based on dual-phase media theory.ChineseJ.Geophys. (in Chinese),58(9):3424-3438,doi:10.6038/cjg20150934.
10.6038/cjg20150934
P631
2014-10-29,2015-06-16收修定稿
國家科技重大專項課題(2011ZX05007-006,2011ZX05019-008)聯(lián)合資助.
桂金詠,男,1986年生,碩士,工程師,主要從事疊前地震反演及儲層預(yù)測等方面的研究工作.E-mail:guijy@petrochina.com.cn