于曉東,馬永強(qiáng),王昱翔,謝 瑋,胡守旺
(中國石化 石油物探技術(shù)研究院,南京 211100)
我國裂縫性油氣藏分布廣泛,在各個地質(zhì)歷史時期都有分布,其在油氣勘探開發(fā)中有著重要地位[1-5]。在種類繁多的非常規(guī)儲層中,裂縫性儲層是影響產(chǎn)能的重要因素之一。裂縫性儲層勘探開發(fā)的難點在于如何預(yù)測裂縫的分布范圍、發(fā)育程度以及產(chǎn)狀[6-10]。
影響裂縫形成的因素很多,裂縫橫向、縱向變化大,物理性質(zhì)復(fù)雜,有很強(qiáng)的各向異性特征[11-14]。與其他油氣藏不同的是,裂縫多為后期生成,沒有相應(yīng)的沉積環(huán)境特征,使得裂縫性油氣藏更加難以勘探[15-16]。此外,由于裂縫的復(fù)雜性,井間裂縫方向和密度的預(yù)測難于依靠井中結(jié)果外推,當(dāng)探區(qū)內(nèi)缺乏測井?dāng)?shù)據(jù)時,就必須尋找其他解決方法[17-18]。
本次研究中,針對工區(qū)實際的地質(zhì)模型進(jìn)行了裂縫不同發(fā)育程度正演模擬,結(jié)合實際工區(qū)從地震剖面反射特征、優(yōu)選的4種地震屬性以及關(guān)聯(lián)維屬性分析對比裂縫發(fā)育區(qū)和裂縫欠發(fā)育區(qū)的特征,為今后對工區(qū)進(jìn)行的一系列裂縫預(yù)測工作提供借鑒。
有限差分法是一種最常見也是最成熟的一種正演模擬方法。它將波動方程中波場函數(shù)的空間導(dǎo)數(shù)和時間導(dǎo)數(shù)用相應(yīng)的空間、時間差分格式來替代,對介質(zhì)中彈性波場進(jìn)行數(shù)值模擬。該技術(shù)的出現(xiàn)和發(fā)展對于人們分析波動傳播規(guī)律,表征地下地質(zhì)構(gòu)造,分析實際地震資料以及油氣資源的開發(fā)有著至關(guān)重要的意義[9-10]。
一般有限差分?jǐn)?shù)值模擬是使用的笛卡爾坐標(biāo)系中的規(guī)則網(wǎng)格,但是在對復(fù)雜地質(zhì)體進(jìn)行模擬時,就會出現(xiàn)階梯狀邊界,這種階梯狀的邊界勢必會引起人為的繞射波,為了減弱這種現(xiàn)象就必須采用精細(xì)網(wǎng)格,這將大大增加數(shù)值模擬的計算量和計算儲存量。因此,我們需要開發(fā)可變網(wǎng)格和不規(guī)則網(wǎng)格的數(shù)值模擬方法。
這里使用三維速度-應(yīng)力方程交錯網(wǎng)格有限差分彈性波數(shù)值模擬方法。
從波動方程的三個基本方程—位移與應(yīng)變方程、位移與應(yīng)力方程以及應(yīng)力與應(yīng)變方程??梢酝茖?dǎo)出關(guān)于壓力P和體變系數(shù)K為變量的聲波方程:
(1)
其中:K為聲學(xué)介質(zhì)體積模量K=λ;V為縱波速度;P為聲壓;ρ為密度。
在式(1)中,為了計算的簡便,避免對彈性常數(shù)進(jìn)行空間微分,可以分別使用Vx、Vy、Vz表示彈性體中質(zhì)點在x、y、z方向上的速度分量,可以得到一階彈性波方程(式(2))。
(2)
對聲波方程(式(1))進(jìn)行降階處理,利用位移與質(zhì)點速度的關(guān)系可以得到式(3):
(3)
式(2)和式(3)可以表示為一階雙曲型標(biāo)量方程式(4)。
(4)
使用Taylor級數(shù)展開一元函數(shù)u(x):
u(x+dx)=u(x)+u′(x)Δx+
(5)
u(x-dx)=u(x)-u′(x)Δx+
(6)
可得:
(7)
一階向前差分:
(8)
一階向后差分:
(9)
一階中心差分:
(10)
差分算法就是將函數(shù)在對應(yīng)的點處用Taylor級數(shù)展開法來代替偏導(dǎo)數(shù),在對時間導(dǎo)數(shù)進(jìn)行正向和反向差分時產(chǎn)生了顯式和隱式兩種有限差分方法,一般情況下使用顯式法求解波動方程。差分算子是一種空間局部算子,其在空間域有著高分辨率,適用于復(fù)雜地質(zhì)構(gòu)造,其缺點在于在頻率域中分辨率低,其穩(wěn)定性以及收斂性受空間采樣率影響,但是該算法運(yùn)算速度較快。
交錯網(wǎng)格2L階精度有限差分系數(shù)計算公式:
u(i)(x0)+o(Δx2L+2)
m=1,2,…,L
(11)
一階導(dǎo)數(shù)2L階精度中心差分近似可用式(12)表示。
eLu(2L+1)(x0)Δx2L+1+o(Δx2(L+1)+1)
(12)
則有交錯網(wǎng)格2L階精度中心有限差分系數(shù)計算公式:
(13)
其中,差分系數(shù)為:
(14)
截斷誤差系數(shù)為式(15)。
(15)
當(dāng)L→0時:
(16)
在式(13)中:①L=1時,a1=1,e1=1/24;②L=2時,a1=9/8,a2=-1/24,e2=-3/640;③三維聲波波動方程差分格式。
壓力值p在(I,J,K)位置取值,νx在(I+1/2,J,K)處取值,vy在(I,J+1/2,K)處取值,vz在(I,J,K+1/2)處取值,如圖1所示。
三維聲波有限差分方程為:
圖1 三維聲波交錯網(wǎng)格Fig.1 3D acoustic staggered grid
表1 聲波波場分量和參數(shù)空間位置Tab.1 Components of acoustic wave field and spatial positions of parameters
(17)
(18)
(19)
(20)
從工區(qū)實際地震剖面和油藏剖面出發(fā),建立過排60井和排66井的地震地質(zhì)模型。工區(qū)石炭系地層含有不同程度的裂縫發(fā)育,模型中在排60、排66井之間布置了裂縫發(fā)育區(qū)和裂縫中等發(fā)育區(qū)。其次,為了對比裂縫發(fā)育區(qū)和無裂縫發(fā)育區(qū)的地震波場特征,在模型左側(cè)布置了無裂縫發(fā)育區(qū),這樣可以對比研究更好地指導(dǎo)裂縫預(yù)測。
對模型進(jìn)行波動方程正演模擬后可以得到無裂縫發(fā)育位置和裂縫發(fā)育位置的波場快照(圖5)。相比圖5(a)、圖5(b)可以看出,當(dāng)?shù)卣鸩▊鞑サ搅芽p時,裂縫的散射使得反射波雜亂,同時繞射波互相干涉疊加,使波前面不連續(xù),出現(xiàn)斷斷續(xù)續(xù)的現(xiàn)象。
圖2 油藏剖面Fig.2 3d acoustic staggered grid
圖3 地震剖面Fig.3 Seismic profile
圖4 地震地質(zhì)模型Fig.4 Seismic geological model
此外,還可以通過波動方程正演模擬得到無裂縫發(fā)育區(qū)和裂縫發(fā)育區(qū)的單炮記錄(圖6)。對比圖6(a)、圖6(b)發(fā)現(xiàn),裂縫發(fā)育區(qū)的單炮記錄,其波形更加不連續(xù),呈現(xiàn)出雜亂無章的現(xiàn)象。
對圖4所示的模型整體進(jìn)行波動方程正演模擬以及疊前深度偏移處理,從圖7中可以看到,裂縫越是發(fā)育的區(qū)域,其由于裂縫散射所導(dǎo)致的反射波雜亂使得地震波出現(xiàn)斷斷續(xù)續(xù)反射的現(xiàn)象,同時反射強(qiáng)度也相對較弱。
對模型進(jìn)行波動方程正演模擬及屬性提取后,可以得到圖8所示的地震屬性對比分析圖。圖8中可以看出,裂縫發(fā)育區(qū)域相關(guān)系數(shù)曲線平均值顯示為高值;地震弧長曲線、反射強(qiáng)度曲線以及均方根振幅曲線顯示為低值。
圖5 波動方程數(shù)值模擬波場快照Fig.5 Numerical simulation of wave field snapshot by wave equation(a)無裂縫發(fā)育區(qū);(b)裂縫發(fā)育區(qū)
圖6 波動方程數(shù)值模擬單炮記錄Fig.6 Numerical simulation of single gun recording by wave equation(a)無裂縫發(fā)育區(qū);(b)裂縫發(fā)育區(qū)
圖7 波動方程數(shù)值模擬疊前深度偏移剖面Fig.7 Numerical simulation of pre-stack depth migration profile by wave equation(a)深度域偏移剖面;(b)時間域偏移剖面
圖8 疊前深度偏移剖面中石炭系地震屬性分析Fig.8 Analysis of carboniferous seismic attributes in the prestack depth migration profile
圖9 疊前深度偏移剖面中石炭系地震信號關(guān)聯(lián)維分析Fig.9 Correlation dimension analysis of carboniferous seismic signals in the prestack depth migration profile
圖10 車排子地區(qū)地震反射特征Fig.10 seismic reflection characteristics in chepaizi region
此外,我們還對數(shù)值模擬的結(jié)果提取了關(guān)聯(lián)維屬性,從圖9可以看出,在裂縫欠發(fā)育和裂縫不發(fā)育的區(qū)域關(guān)聯(lián)維的值一般大于3,而在裂縫發(fā)育區(qū)關(guān)聯(lián)維的值一般小于3,可見關(guān)聯(lián)維屬性能夠較好地分辨裂縫發(fā)育區(qū)的空間位置以及發(fā)育程度。
對于車排子地區(qū)的石炭系火成巖來說,其巖性、巖相變化較快,因此地震記錄的信噪比低,圖10為排662、排66和排673井區(qū)的地震反射特征圖,其中排66、排673為高產(chǎn)井,其地震反射呈現(xiàn)中等-弱反射,同相軸錯斷,而排662井為裂縫欠發(fā)育區(qū),同相軸為強(qiáng)反射,連續(xù)性好。
針對上一部分優(yōu)選出的地震弧長屬性、相關(guān)系數(shù)屬性、反射強(qiáng)度屬性以及均方根振幅屬性,選取了車排子工區(qū)沿石炭系頂界面0 ms~15 ms和15 ms~200 ms時窗提取的各個屬性平面圖,如圖11所示。
從圖11可以看出,地震相關(guān)系數(shù)屬性普遍屬于高值,地震弧長、反射強(qiáng)度以及均方根振幅普遍屬于低值。
利用分形分維對工區(qū)進(jìn)行裂縫預(yù)測(圖12),圖12中藍(lán)色區(qū)域為裂縫發(fā)育區(qū),裂縫密度較大;紅黃色區(qū)域裂縫發(fā)育程度較低,相比藍(lán)色區(qū)域裂縫發(fā)育較少,為裂縫低密度區(qū)。由圖12可以發(fā)現(xiàn),裂縫發(fā)育區(qū)的關(guān)聯(lián)維值小于3,裂縫欠發(fā)育區(qū)的關(guān)聯(lián)維值大于3,裂縫預(yù)測發(fā)育區(qū)與車排子地區(qū)勘探成果圖有較好的對應(yīng)關(guān)系。
圖11 車排子工區(qū)石炭系頂界面以下0 ms~15 ms和15 ms~200 ms地震屬性平面圖Fig.11 Seismic attribute plans of 0 ms~15ms and 15 ms~200 ms below the top interface of carboniferous system in chepaizi working area(a)0 ms~15 ms地震弧長相關(guān)系數(shù);(b)0 ms~15 ms反射強(qiáng)度均方根振幅;(c)15 ms~200 ms地震弧長相關(guān)系數(shù);(d)15 ms~200 ms反射強(qiáng)度均方根振幅
圖12 分形分維裂縫密度預(yù)測圖Fig.12 Fractal fractal fracture density prediction diagram
圖13 車排子地區(qū)勘探成果圖Fig.13 Map of exploration results in chepaizi area
從實際工區(qū)需要利用交錯網(wǎng)格有限差分?jǐn)?shù)值模擬,對工區(qū)實際地質(zhì)情況進(jìn)行正演模擬研究,建立裂縫不同程度發(fā)育情況模型進(jìn)行正演模擬,最后進(jìn)行實例應(yīng)用分析得到以下結(jié)論:
1)在裂縫發(fā)育程度正演模擬中,可以看到裂縫發(fā)育強(qiáng)度越大,其反射波越雜亂,同相軸錯段以及地震波不連續(xù)的現(xiàn)象越嚴(yán)重。
2)裂縫發(fā)育區(qū)相較裂縫欠發(fā)育區(qū)其地震記錄會出現(xiàn)同相軸錯斷、不連續(xù)的現(xiàn)象。
3)裂縫發(fā)育區(qū)相關(guān)系數(shù)平均值為高值;地震弧長、反射強(qiáng)度以及均方根振幅屬性為低值。
4)車排子工區(qū)裂縫欠發(fā)育和裂縫不發(fā)育的區(qū)域關(guān)聯(lián)維的值一般大于3,而在裂縫發(fā)育區(qū)關(guān)聯(lián)維的值一般小于3,關(guān)聯(lián)維屬性能夠較好地分辨裂縫發(fā)育區(qū)的空間位置以及發(fā)育程度。
工區(qū)裂縫預(yù)測結(jié)果與之前的勘探成果吻合度較高,預(yù)測結(jié)果能夠較好地反映該地區(qū)裂縫發(fā)育特征,為該工區(qū)今后石炭系儲層的勘探開發(fā)提供了借鑒。