黃光南,劉 洋,李紅星,張 華,李澤林

(1.東華理工大學(xué)放射性地質(zhì)與勘探技術(shù)國防重點學(xué)科實驗室,江西南昌330013;2.中國石油大學(xué)(北京)油氣資源與探測國家重點實驗室,北京102249;3.中國石油大學(xué)(北京)CNPC物探重點實驗室,北京102249)

各向異性現(xiàn)象在地下介質(zhì)中普遍存在,室內(nèi)研究和野外觀測表明大多數(shù)巖石存在明顯的各向異性性質(zhì)。引起巖石各向異性性質(zhì)的主要因素有：巖石的顆粒排列、晶格結(jié)構(gòu)、薄層疊加、裂隙裂縫以及方向應(yīng)力[1-4]。由于地震波在各向同性介質(zhì)與各向異性介質(zhì)中傳播的差異性,地震數(shù)據(jù)的處理、解釋和反演必須考慮到巖石的各向異性性質(zhì),如果仍然延用各向同性介質(zhì)的假設(shè),可能會導(dǎo)致一些錯誤的處理和解釋結(jié)果。從地球物理反演角度分析,敏感性是指模型參數(shù)產(chǎn)生一個微小擾動量時,對應(yīng)觀測數(shù)據(jù)的變化大小。速度反演通常利用雅克比矩陣對當前模型不斷更新,雅克比矩陣的元素一般是偏微分算子,這個偏微分算子就是速度敏感函數(shù)。如果某個彈性參數(shù)對應(yīng)的速度敏感函數(shù)很靈敏,則經(jīng)過幾次迭代就可以得到正確的反演結(jié)果;相反,如果該彈性參數(shù)對應(yīng)的速度敏感函數(shù)不靈敏,則無論經(jīng)過幾次迭代反演都很難得到正確的反演結(jié)果。

關(guān)于各向異性介質(zhì)速度敏感函數(shù)的計算方法,Every等[5]、Chapman等[6]和Zhou等[7]分別做了大量的研究工作。1992年,Every等在研究各向異性介質(zhì)的彈性參數(shù)時,將彈性參數(shù)剛度矩陣分解成各向同性介質(zhì)彈性參數(shù)矩陣與各向異性介質(zhì)擾動量矩陣之和,并且通過求解特征方程得到了縱波相速度表達式,縱波相速度表達式對彈性參數(shù)求偏微分得到了敏感性系數(shù)。1996年,Chapman等對各向異性介質(zhì)的速度敏感函數(shù)做了研究,利用“固定圖釘法”原理分析了qP波和qSV波在介質(zhì)對稱軸傾角為0,45°和90°方向的速度值。在此基礎(chǔ)之上,給定各向異性介質(zhì)參數(shù)一個擾動量,分析qP波和qSV波在0,45°和90°方向的相速度變化大小。2005年,Zhou等基于VTI介質(zhì)模型對速度敏感函數(shù)進行了詳細研究,分別利用特征值法與特征向量法求取了速度敏感函數(shù)表達式,其中特征值法包括相速度偏微分表達式和群速度偏微分表達式。

分析可知，上述學(xué)者們所做的研究工作尚存在一些不足之處。Every等推導(dǎo)出的敏感性系數(shù)表達式是從材料科學(xué)的角度提出的,表達式極其復(fù)雜,所以缺乏實用性。Chapman等推導(dǎo)的相慢度擾動量表達式非常有意義,但是該表達式只能針對特定的幾個角度,未能描述不同觀測角度時相慢度擾動量的變化。Zhou等推導(dǎo)的速度敏感函數(shù)是基于VTI介質(zhì)模型,其中速度敏感函數(shù)的群速度偏微分表達式非常復(fù)雜,在實現(xiàn)各向異性走時反演時比較困難。

參考前人基于VTI介質(zhì)模型的速度敏感函數(shù)計算方法,我們推導(dǎo)了基于TTI介質(zhì)模型以相速度偏微分表達的速度敏感函數(shù)計算方法。針對層間塊狀異常體模型,分別對背景介質(zhì)與異常體介質(zhì)計算qP,qSV和qSH波的速度敏感函數(shù),同時分析VTI介質(zhì)和TTI介質(zhì)情形下不同波模式的速度敏感函數(shù)特性。為了驗證TTI介質(zhì)速度敏感函數(shù)計算方法的正確性,利用TI介質(zhì)反射波非線性走時反演算法對層間塊狀異常體模型進行了走時反演。

1 TTI介質(zhì)速度敏感函數(shù)計算方法

1982年,Cerveny等[8]提出了各向異性介質(zhì)一階走時擾動方程,很多學(xué)者提出的地震走時反演表達式都是對該方程的演化[9-11]。Zhou等[12]推導(dǎo)出了走時擾動方程的3種表達形式,其中相對簡單的是一階走時擾動方程的相速度偏微分表達式：

(1)

式中：c和U分別為各向異性介質(zhì)的相速度和群速度;av∈{a11,a13,a33,a44,a66,θ0,φ0}為各向異性介質(zhì)模型的彈性參數(shù)、對稱軸傾角和方位角。由公式(1)可以看出,相速度的偏微分(?c/?av)是各向異性介質(zhì)走時擾動方程的重要組成成分,對各向異性介質(zhì)走時反演至關(guān)重要。

速度敏感函數(shù)是相速度(或群速度)對各向異性介質(zhì)參數(shù)的偏微分表達式。由于群速度是相速度的函數(shù),所以利用相速度偏微分表達速度敏感函數(shù)更加簡便。我們基于TTI介質(zhì)模型推導(dǎo)相速度對介質(zhì)彈性模量參數(shù)的偏微分表達的速度敏感函數(shù)計算方法。

TTI介質(zhì)模型由5個彈性參數(shù)(a11,a13,a33,a44,a66)或者Thomsen參數(shù)(α0,β0,γ,ε,δ*)來定義?；诿枋龈飨虍愋越橘|(zhì)中地震波傳播運動學(xué)特性的克里斯托夫方程,在TTI介質(zhì)中克里斯托夫矩陣的特征值,即相速度的表達式也可以寫成如下形式[13]：

(2)

其中,P和Q的表達式為

(3)

Q1,Q2和Q3的表達式為

(4)

式中：?為相慢度向量n和TTI介質(zhì)對稱軸ez′之間的夾角。n=(sinθcosφ,sinθsinφ,cosθ);ez′=(sinθ0cosφ0,sinθ0sinφ0,cosθ0);θ和φ分別為相慢度向量n的傾角和方位角,θ從垂直向下的z軸方向開始測量,φ從水平方向的x軸方向開始測量;θ0和φ0為TTI介質(zhì)對稱軸的傾角和方位角。將(3)式和(4)式代入(2)式并對av∈{a11,a13,a33,a44,a66,θ0,φ0}求偏導(dǎo)數(shù),得到3個相速度c1,c2和c3對TTI介質(zhì)彈性模量參數(shù)的偏微分表達式：

(5)

(6)

(7)

由TTI介質(zhì)的對稱軸向量與相慢度向量之間的關(guān)系,可以計算出兩個向量之間的夾角和夾角的余弦表達式：

(8)

(9)

將(6)式、(7)式和(9)式逐一代入(5)式,即可以得到相速度c1,c2和c3對TTI介質(zhì)彈性模量參數(shù)av∈{a11,a13,a33,a44,a66,θ0,φ0}的偏微分表達式,即速度敏感函數(shù)。利用彈性模量參數(shù)與Thomsen參數(shù)之間的關(guān)系式,可以進一步得到相速度對Thomsen參數(shù)的偏微分表達式。

2 速度敏感函數(shù)的數(shù)值模擬分析

利用數(shù)值模型對上述TTI介質(zhì)速度敏感函數(shù)計算方法進行測試,通過數(shù)值模擬來分析VTI介質(zhì)模型和TTI介質(zhì)模型的速度敏感函數(shù)特性。圖1 給出了測試采用的層間塊狀異常體模型,模型水平方向的長度為800m,深度為500m;模型的背景介質(zhì)為VTI介質(zhì),其彈性模量參數(shù)為a11=12,a13=5,a33=10,a44=4,a66=6,對稱軸傾角和方位角分別為θ0=0和φ0=0;模型內(nèi)的塊狀異常體為TTI介質(zhì),其彈性模量參數(shù)為a11=20,a13=10,a33=15,a44=8,a66=12,對稱軸傾角和方位角分別為θ0=450和φ0=0。

圖1 層間塊狀異常體模型

2.1 VTI介質(zhì)的速度敏感函數(shù)特性

首先分析模型背景VTI介質(zhì)的速度敏感函數(shù)特性。圖2到圖4分別是VTI介質(zhì)qP,qSV和qSH波的速度敏感函數(shù),由這3張圖可以看到,VTI介質(zhì)中3種波的相速度對方位角φ的偏微分均等于0,這是因為VTI介質(zhì)的相速度不隨觀測方位角的變化而變化。

圖2是qP波的速度敏感函數(shù)圖。從圖2中可以看出?c1/?a11(圖2a),?c1/?a33(圖2c)和?c1/?a44(圖2d)的峰值比?c1/?a13(圖2b)的峰值大,這說明qP波的相速度對不同彈性參數(shù)的敏感性強度存在差別。另外,不同偏微分的峰值對應(yīng)的觀測角度也不一樣,例如,?c1/?a11的峰值在90°附近,?c1/?a33的峰值在0和180°附近,?c1/?a13和?c1/?a44的峰值在45°和135°附近。?c1/?θ0(圖2e)的數(shù)值從0到180°呈逐漸增大趨勢。

圖2 VTI介質(zhì)qP波的相速度(c1)對彈性模量參數(shù)以及介質(zhì)對稱軸傾角和方位角的偏微分a c1對a11的偏微分; b c1對a13的偏微分; c c1對a33的偏微分; d c1對a44的偏微分; e c1對θ0的偏微分; f c1對φ0的偏微分

圖3是qSV波的速度敏感函數(shù)圖。與qP波的速度敏感函數(shù)特性類似,qSV波的相速度對不同彈性參數(shù)的敏感性強度也存在差別,不同偏微分的峰值出現(xiàn)在不同的觀測角度范圍內(nèi)。

圖4是qSH波的速度敏感函數(shù)圖。qP波和qSV波的相速度與a11,a13,a33,a44,θ0,φ0有關(guān),而qSH波的相速度與a44,a66,θ0,φ0有關(guān)。可以看出,?c3/?a44(圖4a)和?c3/?a66(圖4b)的敏感性強度大致相同,但是它們的敏感性峰值存在于不同的觀測角度范圍內(nèi),?c3/?a44的峰值在0和180°附近,?c3/?a66的峰值在90°附近。

從圖2到圖4可知,VTI介質(zhì)不同彈性參數(shù)對應(yīng)的速度敏感性強度存在差別,而且敏感性峰值對應(yīng)的觀測角度也不同。

圖3 VTI介質(zhì)qSV波的相速度(c2)對彈性模量參數(shù)以及介質(zhì)對稱軸傾角和方位角的偏微分a c2對a11的偏微分; b c2對a13的偏微分; c c2對a33的偏微分; d c2對a44的偏微分; e c2對θ0的偏微分; f c2對φ0的偏微分

圖4 VTI介質(zhì)qSH波的相速度(c3)對彈性模量參數(shù)以及介質(zhì)對稱軸傾角和方位角的偏微分a c3對a44的偏微分; b c3對a66的偏微分; c c3對θ0的偏微分; d c3對φ0的偏微分

2.2 TTI介質(zhì)的速度敏感函數(shù)特性

繼續(xù)分析模型中異常體TTI介質(zhì)的速度敏感函數(shù)特性。圖5到圖9分別是qP波、qSV波和qSH波的速度敏感函數(shù)。由圖5到圖9可以看到,TTI介質(zhì)的速度敏感函數(shù)比VTI介質(zhì)的速度敏感函數(shù)復(fù)雜得多,這是因為巖石的對稱軸傾角θ0和方位角φ0,任何一個發(fā)生變化都會使速度敏感函數(shù)產(chǎn)生變化。TTI介質(zhì)的速度敏感函數(shù)粗看起來似乎沒有規(guī)律,實際上仍然可以看出各個彈性參數(shù)對應(yīng)的速度敏感函數(shù)強度大小,以及速度敏感函數(shù)峰值的分布區(qū)間。例如,qP波的速度敏感函數(shù)(圖5,圖6),?c1/?a13比其它5個參數(shù)的敏感性峰值更小,?c1/?a11的敏感性相對還較簡單,其它5個參數(shù)的敏感性變化更加劇烈。

圖5 TTI介質(zhì)qP波的相速度(c1)對介質(zhì)對稱軸傾角θ0(a)和方位角φ0(b)的偏微分

圖6 TTI介質(zhì)qP波的相速度(c1)對彈性模量參數(shù)的偏微分a c1對a11的偏微分; b c1對a13的偏微分; c c1對a33的偏微分; d c1對a44的偏微分

圖7 TTI介質(zhì)qSV波的相速度(c2)對彈性模量參數(shù)的偏微分a c2對a11的偏微分; b c2對a13的偏微分; c c2對a33的偏微分; d c2對a44的偏微分

圖8 TTI介質(zhì)qSV波的相速度(c2)對介質(zhì)對稱軸傾角θ0(a)和方位角φ0(b)的偏微分

圖9 TTI介質(zhì)qSH波的相速度(c3)對彈性模量參數(shù)以及介質(zhì)對稱軸傾角和方位角的偏微分a c3對a44的偏微分; b c3對a66的偏微分; c c3對θ0的偏微分; d c3對φ0的偏微分

3 速度敏感函數(shù)計算方法驗證

為了驗證本文推導(dǎo)出的速度敏感函數(shù)計算方法的正確性,編寫了TI介質(zhì)反射波非線性走時反演算法程序,利用該算法對層間塊狀異常體模型進行qP波走時反演,同時結(jié)合速度敏感函數(shù)分析彈性參數(shù)反演結(jié)果。根據(jù)反演理論,多參數(shù)反演會發(fā)生耦合效應(yīng),這樣會導(dǎo)致反演結(jié)果的不確定性增強,所以這里只對5個彈性參數(shù)進行反演。將背景介質(zhì)作為反演的初始模型,塊狀異常體作為反演目標。在模型上表面(地表)進行激發(fā)和接收,道間距為16m,炮點與檢波點均為51個。模型的網(wǎng)格間距為50m,對于5個彈性參數(shù)剖面,需要求解的未知量一共有935個,而反射波走時有2601個,反問題處于超定形式。圖10為層間塊狀異常體模型正演的射線路徑分布圖。

圖10 層間塊狀異常體模型正演射線路徑分布

圖11為塊狀異常體的速度敏感函數(shù)圖(圖中紫線代表?c1/?a11;紅線代表?c1/?a13;綠線代表?c1/?a33;藍線代表?c1/?a44)。由圖11可以看出,a11的速度敏感函數(shù)(紫色曲線)的敏感性最弱,由此預(yù)測較難得到正確的反演結(jié)果;a33的速度敏感函數(shù)(綠色曲線)在射線照射角度為45°和135°附近時最強,由于地表觀測方式在45°和135°附近的射線較多,可以預(yù)測能夠比較容易地得到正確的反演結(jié)果;a13和a44的速度敏感函數(shù)(紅色和藍色曲線)形態(tài)相似,但是a44比a13的敏感性強,在射線照射角度為0,90°和180°附近時較為敏感,由于水平與垂直方向的射線較少,所以它們的反演難度較a33困難。因此,預(yù)測4個彈性參數(shù)的敏感性強度排序為a33>a44>a13>a11。

圖11 塊狀異常體的速度敏感函數(shù)曲線

圖12為真實彈性參數(shù)剖面與反演剖面的對比圖,左側(cè)為真實剖面,右側(cè)為反演剖面。由于很難反演得到a11參數(shù)的正確結(jié)果,這里只給出了a13,a33和a44的反演剖面。由圖12可以看出,a33的反演效果最好,a13和a44的反演效果相當,反演結(jié)果與速度敏感函數(shù)的預(yù)測結(jié)果一致。由此驗證了TI介質(zhì)速度敏感函數(shù)計算方法的正確性,也證實了速度敏感函數(shù)對反演結(jié)果具有預(yù)測和指示作用。

圖12 層間塊狀異常體模型彈性參數(shù)(a13,a33,a44)的真實剖面(a)與反演剖面(b)

4 結(jié)束語

速度敏感函數(shù)有群速度偏微分和相速度偏微分兩種,由于群速度是相速度的函數(shù),所以相速度偏微分表達速度敏感函數(shù)更加簡便。我們基于TTI介質(zhì)模型推導(dǎo)出了相速度偏微分形式的速度敏感函數(shù)計算方法,并對層間塊狀異常體模型的背景VTI介質(zhì)和塊狀TTI介質(zhì)異常體分別模擬計算了速度敏感函數(shù)。對速度敏感函數(shù)特性進行分析得出5點結(jié)論：①同一種波模式,不同彈性參數(shù)的速度敏感性強度不同;②速度敏感函數(shù)在不同觀測角度的敏感性強度不同;③VTI介質(zhì)的速度敏感函數(shù)在觀測傾角為0～90°和90°～180°呈對稱關(guān)系;④VTI介質(zhì)的速度敏感函數(shù)不隨觀測方位角發(fā)生變化;⑤TTI介質(zhì)的速度敏感函數(shù)比VTI介質(zhì)復(fù)雜得多,它的敏感性隨觀測角度的變化而變化。

由于速度敏感函數(shù)是各向異性介質(zhì)走時擾動方程的重要組成成分,文中利用TI介質(zhì)反射波非線性走時反演算法對層間塊狀異常體模型進行反演,得到了a13,a33和a44參數(shù)的反演結(jié)果,從而驗證了本文速度敏感函數(shù)計算方法的正確性,也證實了速度敏感函數(shù)對反演結(jié)果具有預(yù)測和指示作用。

致謝:感謝澳大利亞阿德萊德大學(xué)周兵高級研究員對本文給予的指導(dǎo)和幫助。

參 考 文 獻

[1] Crampin S.Effective anisotropic elastic constants for wave-propagation through cracked solids[J].Geophysical Journal of the Royal Astronomical Society,1984,76：135-145

[2] Helbig K.Systematic classification of layered-induced transverse isotropy[J].Geophysical Prospecting,1981,29：550-577

[3] Thomsen L.Weak elastic anisotropy[J].Geophysics,1986,51(10)：1954-1966

[4] 李芳,曹思遠,姚健.任意各向異性介質(zhì)相(群)速度的計算[J].地球物理學(xué)報,2012,55(10):3420-3426

Li F,Cao S Y,Yao J.Calculation of phase and group velocities in an arbitrary anisotropic medium[J].Chinese Journal of Geophysics,2012,55(10):3420-3426

[5] Every A G,Sachse W.Sensitivity of inversion algorithms for recovering elastic constants of anisotropic solids from longitudinal wavespeed data[J].Ultrasonics,1992,30(4):43-48

[6] Chapman C H,Miller D E.Velocity sensitivity in transversely isotropic media[J].Geophysical Prospecting,1996,44:525-549

[7] Zhou B,Greenhalgh S A.Analytic expressions for the velocity sensitivity to the elastic moduli for the most general anisotropic media[J].Geophysical Prospecting,2005,53:619-641

[8] Cerveny V,Jech J.Linearized solutions of kinematic problems of seismic body waves in inhomogeneous slightly anisotropic media[J].Journal of Geophysics,1982,51:96-104

[9] Jech J,Psencik I.First-order perturbation method for anisotropic media[J].Geophysical Journal International,1989,99:369-376

[10] Chapman C H,Pratt R G.Traveltime tomography in anisotropic media-I[J].Geophysical Journal International,1992,109：1-19

[11] Pratt R G,Chapman C H.Traveltime tomography in anisotropic media-II[J].Geophysical Journal International,1992,109：20-37

[12] Zhou B,Greenhalgh S A.Non-linear traveltime inversion for 3-D seismic tomography in strongly anisotropic media[J].Geophysical Journal International,2008,172:383-394

[13] Zhou B,Greenhalgh S A.Raypath and traveltime computations for 2D transversely isotropic media with dipping symmetry axes[J].Exploration Geophysics,2006,37(2):150-159