基于地震屬性建模的數(shù)值模擬及其應(yīng)用

李良泉, 熊曉軍, 侯秋平, 楊瑞齊, 陳　琴

(1.地球探測與信息技術(shù)教育部重點(diǎn)實(shí)驗(yàn)室(成都理工大學(xué)),成都 610059;2.中國石油川慶鉆探工程有限公司 地質(zhì)勘探開發(fā)研究院，成都 610051)

?

基于地震屬性建模的數(shù)值模擬及其應(yīng)用

李良泉1, 熊曉軍1, 侯秋平1, 楊瑞齊1, 陳琴2

(1.地球探測與信息技術(shù)教育部重點(diǎn)實(shí)驗(yàn)室(成都理工大學(xué)),成都 610059;2.中國石油川慶鉆探工程有限公司 地質(zhì)勘探開發(fā)研究院，成都 610051)

為了解決花崗巖潛山儲(chǔ)層內(nèi)部結(jié)構(gòu)復(fù)雜、非均質(zhì)性強(qiáng)、在數(shù)值模擬建模過程中模型參數(shù)填充誤差大的問題，本文提出了一種基于地震屬性建模的數(shù)值模擬技術(shù)。該技術(shù)首先優(yōu)選地震屬性建立儲(chǔ)層模型；然后采用波動(dòng)方程的數(shù)值模擬方法對(duì)模型進(jìn)行“正演+偏移”計(jì)算，并對(duì)比分析模擬得到的偏移剖面與原始地震剖面，逐步修正觀測系統(tǒng)參數(shù)，直至兩者的主要特征一致；最后，結(jié)合其他地震屬性(如裂縫檢測屬性)，綜合分析目的層的地震響應(yīng)特征。對(duì)渤海某區(qū)花崗巖儲(chǔ)層的數(shù)值模擬結(jié)果表明，該方法有效地刻畫了花崗巖“溶蝕孔隙+高角度裂縫”所形成的“串珠狀”反射特征，為該區(qū)后續(xù)的勘探研究奠定了基礎(chǔ)。利用本文提出的方法綜合分析儲(chǔ)層結(jié)構(gòu)的地震波場特征是一種可行的方法，該方法除了應(yīng)用于花崗巖儲(chǔ)層的數(shù)值模擬外，還可以應(yīng)用于其他復(fù)雜儲(chǔ)層的數(shù)值模擬。

數(shù)值模擬；地震屬性；花崗巖；儲(chǔ)層模型

隨著勘探的不斷深入，勘探目標(biāo)逐漸擴(kuò)大到非沉積巖儲(chǔ)層，其中基巖油氣藏正是一個(gè)潛在勘探目標(biāo)。目前世界范圍內(nèi)的基巖潛山油氣田中，約40%產(chǎn)于花崗巖潛山儲(chǔ)層。中國發(fā)現(xiàn)的花崗巖儲(chǔ)層油氣田大都屬于多期變質(zhì)的混合花崗巖，巖性復(fù)雜，地質(zhì)年代不太準(zhǔn)確[1]，對(duì)該類儲(chǔ)層的評(píng)價(jià)存在較大困難，屬于一個(gè)比較難的勘探領(lǐng)域。

本文的研究區(qū)是渤海某區(qū)的花崗巖儲(chǔ)層，目前已有很多學(xué)者對(duì)該區(qū)的花崗巖儲(chǔ)層做了大量的基礎(chǔ)研究：王昕等[2]對(duì)渤海海域PL花崗巖潛山儲(chǔ)層的發(fā)育機(jī)制和儲(chǔ)層特征進(jìn)行了仔細(xì)研究，并將潛山內(nèi)部細(xì)分為土壤帶、砂-礫質(zhì)分化帶、裂縫帶和基質(zhì)帶；王霄[3]對(duì)潛山的構(gòu)造演化、花崗巖儲(chǔ)層特征、油源對(duì)比及成藏條件進(jìn)行了剖析，分析潛山內(nèi)部4個(gè)帶的儲(chǔ)層特性，建立了油氣成藏模式；陳國成等[4]綜合巖礦和測井資料，分析了花崗巖儲(chǔ)層的巖性特征、儲(chǔ)集空間類型及物性特征；等。

數(shù)值模擬作為勘探開發(fā)中一個(gè)必不可少的環(huán)節(jié)，尤其對(duì)于花崗巖這種復(fù)雜儲(chǔ)層，更加顯得重要。張華等[5]針對(duì)花崗巖裂隙儲(chǔ)層各向異性進(jìn)行了高階有限差分?jǐn)?shù)值模擬，并認(rèn)為在實(shí)際勘探中可以將花崗巖裂隙介質(zhì)等效為各向同性介質(zhì)；黃繼新等[6]探討了一種綜合運(yùn)用地震、測井、地質(zhì)資料的方法，以大港油田某區(qū)塊為例建立了基于地質(zhì)統(tǒng)計(jì)和地震反演的定量儲(chǔ)層地質(zhì)模型；龔再升[1]闡述了花崗巖儲(chǔ)層油氣藏的特點(diǎn)以及花崗巖儲(chǔ)層的巖性與分類，總結(jié)了花崗巖潛山地震響應(yīng)特征；等。

考慮到渤海某區(qū)花崗巖潛山儲(chǔ)層結(jié)構(gòu)復(fù)雜，內(nèi)部物性不均勻，儲(chǔ)集空間類型多為孔隙-裂縫型，儲(chǔ)層在縱橫向均具有較強(qiáng)的非均質(zhì)性[7]，測井參數(shù)統(tǒng)計(jì)時(shí)精確度不夠，且井間信息不足等因素，采用常規(guī)的儲(chǔ)層建模方法難以準(zhǔn)確地反映花崗巖儲(chǔ)層的地質(zhì)特征。因此，本文提出了一種基于地震屬性建模的數(shù)值模擬技術(shù)，利用該技術(shù)得到研究區(qū)的數(shù)值模擬剖面，并總結(jié)花崗巖儲(chǔ)層的地震反射特征。

1　方法原理

本文提出的基于地震屬性建模的數(shù)值模擬技術(shù)，首先建立基于地震屬性的儲(chǔ)層模型；然后對(duì)儲(chǔ)層模型采用波動(dòng)方程方法進(jìn)行數(shù)值模擬，得到模擬結(jié)果；最后綜合其他地震屬性進(jìn)行分析。圖1為基于地震屬性建模的數(shù)值模擬技術(shù)流程圖。

圖1　技術(shù)流程圖Fig.1　Technique flow chart for seismic attribution

1.1基于地震屬性的儲(chǔ)層建模

對(duì)于數(shù)值模擬，儲(chǔ)層模型的建立至關(guān)重要，也是制約數(shù)值模擬結(jié)果的重要因素；此外，模型的精度也直接影響數(shù)值模擬的準(zhǔn)確性?；诘卣饘傩缘慕＜夹g(shù)包含下述3個(gè)步驟。

a.首先優(yōu)選過井剖面，綜合利用目的層位、斷層發(fā)育，以及地質(zhì)縱向分帶等信息建立初始模型框架；然后整理和統(tǒng)計(jì)測井?dāng)?shù)據(jù)，對(duì)模型框架進(jìn)行物性參數(shù)填充，從而建立參數(shù)模型。

b.利用疊前EI反演的縱波阻抗數(shù)據(jù)換算出縱波速度。

常用的密度、速度、波阻抗之間轉(zhuǎn)換的經(jīng)驗(yàn)公式有Gardner公式[8]、Lindseth公式[9]等，本文根據(jù)Gardner經(jīng)驗(yàn)公式形式，利用最小二乘法進(jìn)行系數(shù)擬合，得到波阻抗(Z)與速度(v)之間的換算關(guān)系式。Gardner經(jīng)驗(yàn)公式如下[8]

ρ=a vm

(1)

對(duì)(1)式兩邊取對(duì)數(shù)可得

lnρ=mlnv+lna

(2)

根據(jù)測井曲線處密度(ρ)和速度(v)的值，采用最小二乘法對(duì)(2)式的系數(shù)a和m進(jìn)行擬合。

通過對(duì)研究區(qū)的花崗巖儲(chǔ)層的測井曲線統(tǒng)計(jì)分析，可得公式(2)中的系數(shù)為a=1.175 7，m=0.089 8。由此可得波阻(Z)抗與速度的轉(zhuǎn)換公式

v=0.859Z0.917

(3)

c.根據(jù)層位控制，利用換算出的縱波速度對(duì)初始模型進(jìn)行速度替換，建立最終模型。

1.2數(shù)值模擬方法原理

波動(dòng)方程數(shù)值模擬可以分為單程波動(dòng)方程數(shù)值模擬和雙程波動(dòng)方程數(shù)值模擬。為了模擬花崗巖儲(chǔ)層的主要反射特征，正確模擬它的一次反射波已滿足需要，因此本文采用單程波動(dòng)方程數(shù)值模擬方法[10]。

對(duì)于各向同性介質(zhì)，取二維標(biāo)量聲波方程作為延拓的基本方程

(4)

式中:P=P(x,z,t)為二維地震波場；t為時(shí)間；v(x,z)為縱、橫向都可變的地震波傳播速度。

對(duì)于波動(dòng)方程數(shù)值模擬而言，波場延拓至關(guān)重要。目前常用的頻率-波數(shù)域波場延拓方法有相移法(PS)[11]、相位移加插值方法(PSPI)[12]、裂步傅里葉方法(SSF)[13]、傅里葉有限差分方法(FFD)[14]、廣義屏方法(GSP)等。綜合考慮地震波場數(shù)值模擬的精度和計(jì)算效率等因素，本文選擇裂步傅里葉方法進(jìn)行頻率-波數(shù)域的波場延拓計(jì)算。

首先定義速度的倒數(shù)為慢度，設(shè)地下介質(zhì)每一點(diǎn)的慢度為s(x,z)，將慢度分解為橫向可變和橫向不變兩部分[10]

s(x,z)=s0(z)+Δs(x,z)

(5)

將式(5)代入式(4)，并將波動(dòng)方程變換到頻率域，得

(6)

略去Δs2得

(7)

(8)

使用參考慢度s0(z)作相移延拓

(9)

P1(x, zi+Δz; ω)=

(10)

在頻率-空間域中，對(duì)慢度的擾動(dòng)Δs(x,z)，使用第二次相移延拓

P1(x, zi+Δz; ω)·eiωΔs(x, zi)

(11)

在t=0時(shí)成像

(12)

式(9)、(11)構(gòu)成分步傅里葉波場延拓的基本公式，式(9)在頻率-波數(shù)域?qū)崿F(xiàn)，式(11)在頻率-空間域?qū)崿F(xiàn)。

2　實(shí)例分析

本文研究區(qū)的花崗巖儲(chǔ)層為中生界花崗巖，主要受表生喀斯特作用、有機(jī)酸溶蝕作用、構(gòu)造變形作用及微古地貌作用的共同控制，由表及里隨著風(fēng)化淋濾作用的逐漸減弱，可以劃分為土壤帶、砂-礫質(zhì)風(fēng)化帶、裂縫帶和基巖帶，儲(chǔ)集空間由孔隙型→裂縫-孔隙型→孔隙-裂縫型→裂縫型呈現(xiàn)規(guī)律性的變化[2]。

2.1實(shí)際剖面與典型屬性分析

圖2是研究區(qū)原始地震資料，t8層位為潛山頂界面，向下依次是潛山頂往下30 ms, 60 ms, 150 ms層位。t8層(潛山頂界面)位于波峰，呈連續(xù)的強(qiáng)振幅反射特征，該波峰下波谷能量在橫向上不連續(xù)；潛山內(nèi)部沒有出現(xiàn)連續(xù)反射波組，可見弱振幅、非連續(xù)或塊狀、空白狀雜亂反射特征。

筆者對(duì)地震屬性進(jìn)行了優(yōu)選，考慮到阻抗數(shù)據(jù)最能精細(xì)體現(xiàn)研究區(qū)縱橫向的地質(zhì)特征，且與速度關(guān)系密切，故選波阻抗作為建模過程中的約束條件。圖3是疊前反演的縱波阻抗剖面，紅色部分是波阻抗相對(duì)較低的區(qū)域，大多分布在潛山頂界向下60 ms范圍內(nèi)。

渤海A井完鉆深度為潛山頂界面向下255 m，鉆遇儲(chǔ)層厚度為125.4 m。其中, t8-(t8+30 ms)內(nèi)儲(chǔ)層厚度為37.8 m，t8-(t8+60 ms)內(nèi)儲(chǔ)層厚度為85.7 m。圖4是疊前裂縫檢測剖面，紅色代表裂縫發(fā)育體。根據(jù)測井?dāng)?shù)據(jù)統(tǒng)計(jì)：“土壤

圖4　疊前裂縫檢測剖面Fig.4　Fracture detection in pre-stack section

帶+風(fēng)化帶”主要分布在t8-(t8+30 ms)范圍內(nèi)，以孔隙發(fā)育為主，裂縫欠發(fā)育；“裂縫帶”主要分布在(t8+30 ms)-(t8+60 ms)范圍內(nèi)，大量發(fā)育高角度裂縫。圖中疊前裂縫檢測結(jié)果與測井信息相符。

2.2儲(chǔ)層數(shù)值模擬

2.2.1儲(chǔ)層建模

根據(jù)本文提出的基于地震屬性的建模方法建立渤海某工區(qū)的最終地質(zhì)模型(圖5)。

由于風(fēng)化帶和裂縫帶內(nèi)孔隙、裂縫的存在，導(dǎo)致t8-(t8+60 ms)范圍內(nèi)低速體大量發(fā)育。圖5中紅色代表低速體，可見低速體多發(fā)育在t8-(t8+60 ms)之間，與測井信息相符。該速度模型的準(zhǔn)確性為接下來的數(shù)值模擬提供了保障。

2.2.2數(shù)值模擬與分析

對(duì)建立的最終模型采用波動(dòng)方程數(shù)值模擬，得到偏移剖面(圖6)。

圖5　速度模型Fig.5　Velocity model

圖6　基于圖5的數(shù)值模擬偏移剖面Fig.6　Numerical simulation of shifting section based on Fig.5

由圖6可以看出，t8位于一個(gè)連續(xù)強(qiáng)波峰，其下方的波谷能量橫向上變換較大，潛山頂界(t8)以下波形雜亂，與實(shí)際剖面特征一致；數(shù)值模擬剖面中出現(xiàn)“串珠狀”特征，與典型的“串珠狀”特征(短段的波峰波谷交替出現(xiàn))相符合(圖7)。

圖7　塔中地區(qū)中古8井區(qū)典型“串珠狀”地震反射特征Fig.7　Typical paternoster reflection characteristics in Tazhong area資料來源：塔中地區(qū)中古8井區(qū)

此外，綜合地震屬性剖面(圖3、圖4和圖5)和數(shù)值模擬剖面(圖6)可見：(1)通過速度模型(圖5)和偏移剖面(圖6)的對(duì)比分析可以看出：風(fēng)化帶內(nèi)低速體不發(fā)育的區(qū)域?qū)?yīng)較連續(xù)的波谷(t8波峰下的第1個(gè)波谷)特征；低速體發(fā)育區(qū)域?qū)?yīng)的波能量呈中弱反射特征，并且橫向連續(xù)性變差。(2)通過疊前裂縫檢測剖面(圖4)和偏移剖面(圖6)的對(duì)比分析可以看出：裂縫大量發(fā)育區(qū)，正演記錄上出現(xiàn)“串珠狀”特征，波峰與波谷交替出現(xiàn)；隨著裂縫縱向展布增大，波峰和波谷能量均增強(qiáng)(串珠3與串珠4對(duì)應(yīng)的位置)。

3　結(jié) 論

a.對(duì)于花崗巖潛山儲(chǔ)層，由于其內(nèi)部結(jié)構(gòu)復(fù)雜、非均質(zhì)性強(qiáng)、內(nèi)部物性不均勻，常規(guī)建模技術(shù)難以滿足要求；利用本文提出的屬性建模方法，有效避免了參數(shù)填充不準(zhǔn)確的問題，提高了模型的縱橫向精度。

b.基于波動(dòng)方程數(shù)值模擬得到的偏移剖面能很好地反映潛山內(nèi)部儲(chǔ)層結(jié)構(gòu)的主要波場特征。模擬結(jié)果表明渤海某區(qū)花崗巖潛山頂界面呈強(qiáng)反射特征，潛山內(nèi)部由于物性不均勻?qū)е路瓷涮卣麟s亂。

c.通過模擬結(jié)果與疊前裂縫檢測剖面的對(duì)比分析，總結(jié)了渤海某區(qū)花崗巖潛山儲(chǔ)層“溶蝕孔隙+高角度裂縫”形成的“串珠狀”地震響應(yīng)特征，為后續(xù)儲(chǔ)層預(yù)測奠定了基礎(chǔ)。

d.本文研究表明，結(jié)合地震屬性與數(shù)值模擬結(jié)果，綜合分析儲(chǔ)層結(jié)構(gòu)的地震波場特征是一種可行的方法。該方法除了應(yīng)用于花崗巖儲(chǔ)層的數(shù)值模擬外，還可以應(yīng)用于其他復(fù)雜儲(chǔ)層的數(shù)值模擬，具有較強(qiáng)的實(shí)際應(yīng)用價(jià)值。

[1] 龔再升.繼續(xù)勘探中國近海盆地花崗巖儲(chǔ)層油氣藏[J].中國海上油氣,2010,22(4):213-220.

Gong Z S. Continued exploration of granitic-reservoir hydrocarbon accumulations in China offshore basins [J]. China Offshore Oil and Gas, 2010, 22(4): 213-220. (In Chinese)

[2] 王昕,周心懷,徐國勝,等.渤海海域蓬萊9-1花崗巖潛山大型油氣田儲(chǔ)層發(fā)育特征與主控因素[J].石油與天然氣地質(zhì),2015,36(2):262-270.

Wang X, Zhou X H, Xu G S,etal. Characteristics and controlling factors of reservoirs in Penglai 9-1 large-scale oilfield in buried granite hills, Bohai Sea [J]. Oil & Gas Geology, 2015, 36(2): 262-270. (In Chinese)

[3] 王霄.蓬萊9-1構(gòu)造花崗巖古潛山油氣成藏條件與成藏模式[D].成都:成都理工大學(xué)檔案館,2015.

Wang X. The Hydrocarbon Accumulation Conditions and Accumulation Mode of Granite Buried Hill in Penglai 9-1 Structure [D]. Chengdu: The Archive of Chengdu University of Technology, 2015. (In Chinese)

[4] 陳國成,陳華靖,田曉平.渤海PL油田花崗巖潛山儲(chǔ)層發(fā)育特征及控制因素[J].海洋地質(zhì)前沿, 2016,32(2):14-19.

Chen G C, Chen H J, Tian X P. Characteristics and controlling factors of granitic buried-hill reservoir in PL oilfield of Bohai Sea [J]. Frontiers of Marine Geology, 2016, 32(2): 14-19. (In Chinese)

[5] 張華,龔育齡.花崗巖裂隙儲(chǔ)層各向異性高階有限差分?jǐn)?shù)值模擬[J].東華理工大學(xué)學(xué)報(bào)(自然科學(xué)版),2012, 35(4):416-421.

Zhang H, Gong Y L. Anisotropic numerical simulation of fracture reservoir in granite with high-order finite difference[J]. Journal of East China Institute of Technology (Natural Science Edition), 2012, 35(4): 416-421. (In Chinese)

[6] 黃繼新,彭仕宓,黃述旺.河流相儲(chǔ)層建模中測井和地震資料整合研究[J].中國礦業(yè)大學(xué)學(xué)報(bào),2007,36(1):126-131.

Huang J X, Peng S M, Huang S W, Study of integrating seismic and logging data in fluvial reservoir modeling[J]. Journal of China University of Mining & Technology, 2007, 36(1): 126-131. (In Chinese)

[7] 譚偉雄,王俊瑞,鄧強(qiáng),等.花崗巖儲(chǔ)層儲(chǔ)集性能定量評(píng)價(jià)方法及應(yīng)用[J].中國海上油氣,2015,27(2):31-38.

Tan W X, Wang J R, Deng Q,etal. Quantitative evaluation method for granite reservoir properties and its applications[J]. China Offshore Oil and Gas, 2015, 27(2): 31-38. (In Chinese)

[8] Gardner G H F, Gardner L W, Gregory A R. Formation velocity and density-the diagnostic basics for stratigraphic traps[J]. Geophysics, 2012, 39(6): 770.

[9] Lindseth R O. Synthetic sonic logs — a process for stratigraphic interpretation[J]. Journal of Virology, 1994, 68(11): 7559-7565.

[10] 熊曉軍.單程波動(dòng)方程地震數(shù)值模擬新方法研究[D]. 成都:成都理工大學(xué)檔案館, 2007.

Xiong X J. Method Research of One-Way Wave-Equation Seismic Numerical Simulation [D]. Chengdu: The Archive of Chengdu University of technology, 2007. (In Chinese)

[11] Gazdag J. Wave equation migration with the phase-shift method[J]. Geophysics, 2012, 43(7): 1342-1351.

[12] Gazdag J. Migration of seismic data by phase shift plus interpolation[J]. Geophysics, 1984, 49(2):124-131.

[13] Stoffa P L, Fokkema J T, Freire R M D L,etal. Split-step Fourier migration[J]. Geophysics, 1990, 55(4): 410-421.

[14] Ristow D, Rühl T. Fourier finite-difference migration[J]. Geophysics, 1994, 59(12): 1882-1893.

Numerical simulation and its application based on modeling of seismic attribution

LI Liang-quan1, XIONG Xiao-jun1, HOU Qiu-ping1, YANG Rui-qi1, CHEN Qin2

1.KeyLaboratoryofEarthExploration&InformationTechniquesofChinaMinistryofEducation,ChengduUniversityofTechnology,Chengdu610059,China;2.GeophysicalProspectingCompany,ChuanqingDrillingEngineeringCompany,CNPC,Chengdu610051,China

A numerical simulation method based on modeling constrained by seismic attribution is proposed to cope with the complex structure, strong heterogeneity and inhomogeneity in exploration of granite buried hill reservoirs. Firstly, reservoir model is established by selecting preferred seismic attribution followed by simulation by wave equation to get a satisfied result; then, the simulated cross section and original cross section are correlated and the observation parameters are revised until the major characteristics of the two cross sections reach unanimous. The numerical simulation on granite reservoir in Bohai area indicates that this method can efficiently describe “the string of beads-like” reflection characteristics caused by high angle fractures and corroded fissures in granite, which is very important to subsequent exploration of oil and gas. The numerical simulation method can apply to granite reservoirs and other complex reservoirs.

numerical simulation; seismic attribution; granite; reservoir model

10.3969/j.issn.1671-9727.2016.04.09

1671-9727(2016)04-0454-06

2016-04-30。

國家自然科學(xué)基金項(xiàng)目(41274130)。

李良泉(1993-),男,碩士研究生,研究方向:復(fù)雜儲(chǔ)層地震波場數(shù)值模擬, E-mail:452484842@qq.com。

P631.4

A