裂隙介質(zhì)中地震波傳播特征近似

何現(xiàn)啟，彭凌星，朱自強(qiáng)，魯光銀

(1.湖南省交通規(guī)劃勘察設(shè)計(jì)院有限公司,長沙 410200；2.中南大學(xué)地球科學(xué)與信息物理學(xué)院,長沙 410083)

賀同江等[1]用迭積微分算子法實(shí)現(xiàn)了黏彈性介質(zhì)的地震波場正演模擬,對波傳播特征進(jìn)行了分析。宮猛等[2]從各向同性介質(zhì)中波場數(shù)值模擬的褶積微分算子法出發(fā),推導(dǎo)出了各向異性雙相介質(zhì)中波場傳播數(shù)值計(jì)算的褶積新算法。周銘等[3]應(yīng)用背景噪聲層析成像方法對新疆地區(qū)地殼S波速度結(jié)構(gòu)以及徑向各向異性進(jìn)行了研究,討論了各向異性產(chǎn)生的機(jī)理。Vigen等[4]討論奇異波相鄰方向上近似qS波群速度和極化的問題。Behura等[5]以各向同性黏彈性介質(zhì)為背景介質(zhì),利用攝動(dòng)法推導(dǎo)出黏彈性 VTI 介質(zhì)中P 波、PS 波的反射系數(shù)近似公式。鄧?yán)^新等[6]利用比奧特(Biot)理論研究了介觀尺度孔隙流體流動(dòng)介質(zhì)中地震波的傳播規(guī)律。司薌等[7]推導(dǎo)了傾斜裂隙等效TTI介質(zhì)PP波反射系數(shù)公式,建立了反射系數(shù)與裂隙參數(shù)之間的關(guān)系。Tariq[8]利用擾動(dòng)法求解了橢圓各向異性介質(zhì)中的聲波方程,并用Marmousi模型進(jìn)行了數(shù)值驗(yàn)證。張煜等[9]研究了地震勘探頻段內(nèi)針對含流體孔隙介質(zhì)邊界條件的面波的傳播特性。李雨生等[10]運(yùn)用交錯(cuò)網(wǎng)格高階有限差分法求解三維單斜各向異性介質(zhì)彈性波方程來模擬彈性波的傳播過程。Alexey等[11]用多參數(shù)擾動(dòng)法求解了正交各向異性介質(zhì)中的波動(dòng)方程,大大縮短了求解時(shí)間。何現(xiàn)啟等[12]用擾動(dòng)法推導(dǎo)了極端弱各向異性介質(zhì)中 qP、qS波群速度的一階表達(dá)式。何現(xiàn)啟等[13]研究黏彈性各向異性介質(zhì)中地震波的傳播特征,推導(dǎo)出黏彈性垂直定向裂隙的各向異性介質(zhì)中均勻、非均勻波的精確相速度、慢度和群速度公式。盧方超等[14]通過超聲實(shí)驗(yàn)和核磁共振(NMR)實(shí)驗(yàn),研究了煤孔、裂隙結(jié)構(gòu)應(yīng)力變形各向異性特征。何現(xiàn)啟等[15]通過數(shù)值模型計(jì)算,研究了黏彈性HTI(horizontal transverse isotropy)介質(zhì)中地震波的傳播特征。梁鍇等[16]推導(dǎo)了橢球各向異性介質(zhì)是橫向各向同性介質(zhì)中彈性波速度和偏振表達(dá)式。趙小龍等[17]用擾動(dòng)法推導(dǎo)了正交各向異性介質(zhì)巖石力學(xué)近似方程。Song等[18]用擾動(dòng)法多弱各向異性介質(zhì)中P、S波的反射、投射系數(shù)進(jìn)行了一、二階近似計(jì)算。

盡管對各向異性介質(zhì)中地震波傳播有不少報(bào)道,但對裂隙介質(zhì)中地震波傳播進(jìn)行系統(tǒng)研究的甚少,且由于裂隙介質(zhì)中地震波的精確傳播特征表達(dá)較復(fù)雜,非線性的隱性表達(dá)式不利于參數(shù)反演,因此有必要對其進(jìn)行線性化近似。

本文為了簡化Christoffel方程,建立了波向量坐標(biāo)系(x′，o′，z′),其中z′軸方向與波的傳播方向一致,y′為水平方向,其單位向量為e(1)、e(2)、e(3)表示。當(dāng)波的傳播方向給定時(shí),各向異性介質(zhì)中Christoffel方程的獨(dú)立彈性系數(shù)為六個(gè),而在各向同性介質(zhì)中Christoffel矩陣為對角陣,因此通過對各向同性介質(zhì)的擾動(dòng)便可對弱各向異性介質(zhì)進(jìn)行攝動(dòng)分析。

1 裂隙介質(zhì)中相速度精確表達(dá)式

由裂隙介質(zhì)的彈性矩陣CEDA,可得到裂隙介質(zhì)Christoffel方程，求解方程可得裂隙介質(zhì)中的qP波、qSV波和SH波的相速度為[15,19-20]

(1)

式(1)中:D=[(c11-c55)F-(c33-c55)E2]2+4(c13+c55)2FE2，E=(-cosφsinθsinφ0+sinφcosφ0)，F(xiàn)=(cosφsinθcosφ0+sinφsinφ0)2+cos2θcos2φ。其中，φ為HTI介質(zhì)對稱軸方位角，φ0為裂隙介質(zhì)對稱軸方位角，θ為介質(zhì)對稱軸傾角，c為裂隙介質(zhì)彈性系數(shù)。

同理，可得到偏振向量及群速度的精確表達(dá)式。

2 裂隙介質(zhì)中相速度一階近似

建立野外觀測坐標(biāo)系xoz，其基本單位向量為(x,y,z)，波的傳播方向n=(nx,ny,nz)=(sinθcosφ, sinθsinφ,cosθ)，其中θ為極化角，φ為方位角。通過擾動(dòng)法可推導(dǎo)出弱各向異性條件下，極端各向異性介質(zhì)中地震波傳播特征的近似表達(dá)式，進(jìn)一步將HTI介質(zhì)的彈性系數(shù)矩陣代入相關(guān)公式即可推導(dǎo)出HTI介質(zhì)中地震波傳播特征近似公式。通過將HTI介質(zhì)的彈性系數(shù)矩陣進(jìn)行旋轉(zhuǎn)，可得到裂隙介質(zhì)(extensive dilatancy anisotropy,EDA)的彈性系數(shù)矩陣：

(2)

將式(2)代入極端各向異性介質(zhì)計(jì)算公式，即可推導(dǎo)出裂隙弱各向異性介質(zhì)中qP、qS波的一階相下速度[13,15]:

(3)

3 裂隙介質(zhì)中一階偏振向量

將HTI介質(zhì)中的模型坐標(biāo)轉(zhuǎn)換到測量(勘探)坐標(biāo),并用裂隙弱各向異性參數(shù)表示,可得qP、qS波的一階偏振單位向量[13,15-16,20]。

(4)

式(4)中：偏振角

4 裂隙介質(zhì)中一階群速度

將裂隙介質(zhì)的彈性系數(shù)矩陣代入極端各向異性介質(zhì)的計(jì)算公式可得[13,16,19]:

(5)

式(5)中：w′1±=F1(θ,φ)±[G1(θ,φ)cos2α+H1(θ,φ)sin2α],w′2±=F2(θ,φ)±[G2(θ,φ)cos2α+H2(θ,φ)sin2α]；并且有:

P(θ,φ)=sin2θ[sin2θ(εxcos4φ+εysin4φ+

δzcos2φsin2φ)+cos2θ(δysin2φ+

δxcos2φ)]+ sin2φsin4θ(ε16cos2φ+

ε26sin2φ)+χzsin2φcos2θsin2θ

(6)

εysin2φ(sin2φcos2θ+cos2φ) -

cos2θ(δysin2φ+δxcos2φ)-

ε26sin2φsin2θ(2cos2θcos2φ+cos2φ)-

(7)

εysin2φ(sin2φcos2θ-cos2φ)-

cos2θ(δysin2φ+δxcos2φ)+

ε26sin2φsin2θ(2cos2θcos2φ-cos2φ)-

2cos2θcos2φ)

(8)

[ε26sin2φ(4cos2φ-1)+ε16cos2φ(1-

(9)

εysin4φ)+cos2θ(δysin2φ+δxcos2φ)]+

2sin2θsin2φsin2θ(ε16cos2φ+ε26sin2φ)

(10)

cos2θ(δy-δx)+δzcos2φsin2θ]+

ε16sin3θcos2φ(2cos2φ-1)+

ε26sin3θsin2φ(2cos2φ+1)]+

(11)

εysin2φ(2sin2φsin2θ-1)+cos2θ×

(12)

(13)

sin2θ)[ε16(cos2φcos2φsinθ-

(14)

cos2φ)(sin2φ+cos2φcos2θ)+

(15)

ε16sinθ(cos2φ-sin22φ)-

ε26sinθcos4φ]-

(16)

cos2φsin2θ)+2εysin2φ(2cos2φ-

sin2φsin2θ)]+δy(2sin2φsin2θ-1)+

ε16sin2θsin2φ(cos2φsin2θ-

cos2φ)ε26sin2θsin2φ(cos2φ-

sin2θsin2φ)]+

(17)

5 算法驗(yàn)證

下面主要通過模型計(jì)算來驗(yàn)證式(3)的準(zhǔn)確性,模型主要參數(shù)見表1,相關(guān)計(jì)算結(jié)果見圖1～圖9。

表1 模型參數(shù)

由表1中模型1參數(shù)可計(jì)算HTI彈性系數(shù),通過坐標(biāo)旋轉(zhuǎn)得到裂隙介質(zhì)的彈性系數(shù)矩陣如下:

(18)

由表1中模型2參數(shù)可計(jì)算HTI彈性系數(shù),通過坐標(biāo)旋轉(zhuǎn)可得到裂隙介質(zhì)的彈性系數(shù)矩陣如下:

(19)

圖1 模型1裂隙介質(zhì)精確相速度平視圖Fig.1 Flat view of exact phase velocity of fracture media for model 1

圖2 模型1裂隙介質(zhì)的精確相速度俯視圖Fig.2 Top view of exact phase velocity of fracture media for model 1

圖3 模型1 裂隙介質(zhì)的近似相速度平視圖Fig.3 Flat view of approximation velocity of model 1

圖4 模型1 裂隙介質(zhì)的近似相速度俯視圖Fig.4 Top view of approximation velocity of model 1

圖5 模型1近似相速度相對誤差Fig.5 Relative error of approximation phase velocity for model1

圖6 模型2精確相速度Fig.6 Exact phase velocity of model 2

圖7 模型2近似相速度Fig.7 Approximation velocity of model 2

圖8 模型2 近似相速度的相對誤差Fig.8 Relative error of approximation phase velocity for model 2

圖9 xoy平面內(nèi)P、SH、SV的精確與近似相速度比較Fig.9 Comparation of exact phase velocity and approximation velocity for P、SH、SV in xoy

以上精確相速度通過解christoffel方程得到,近似解由式(3)得到,分析圖1～圖9可知,當(dāng)各向異性系數(shù)為0.1時(shí),S波(包括SH、SV波)近似相速度相對誤差小于10%,而P波的近似相速度相對誤差小于2%；當(dāng)各向異性系數(shù)為0.3時(shí),P波相速度相對誤差在10%以內(nèi),而S波(包括SH、SV波)近似相速度相對誤差最大到20%。當(dāng)各向異性系數(shù)為0.5時(shí),而P波近似相速度相對誤差小于20%,而S波(包括SH、SV波)近似相速度相對誤差最大達(dá)50%,由此可見,當(dāng)各向異性系數(shù)達(dá)0.5時(shí),對S波(包括SH、SV波)而言,式(3)計(jì)算精度已滿足不了要求,而對P波仍然適用。

6 結(jié)論

為了縮短數(shù)值計(jì)算時(shí)間,便于進(jìn)行速度分析及參數(shù)反演,需將地震波的隱性、非線性精確表達(dá)式進(jìn)行近似研究,轉(zhuǎn)換為線性表達(dá)式。在擾動(dòng)法求解弱各向異性條件下qP、qS的相速度、群速度的近似表達(dá)式的基礎(chǔ)上,推導(dǎo)了裂隙介質(zhì)中地震波傳播特征的近似表達(dá)式。

(1)依據(jù)介質(zhì)模型坐標(biāo)與觀測坐標(biāo),各種介質(zhì)彈性系數(shù)矩陣關(guān)系式推導(dǎo)出裂隙介質(zhì)中地震波相速度、偏振向量及群速度的近似表達(dá)。

(2)通過數(shù)值計(jì)算,研究了近似計(jì)算誤誤差,結(jié)果表明,當(dāng)各向異性系數(shù)達(dá)在0.3以內(nèi)時(shí),對一般裂隙各向異性介質(zhì)都適用,且具有較高的精確,在各向異性系數(shù)為0.5時(shí),P波仍然具有較高的精度,但S波的誤差增大到50%。

(3)因?yàn)閿_動(dòng)解是在弱各向異性條件下推導(dǎo)的,相關(guān)近似公式在弱各向異性條件下,精度較高,當(dāng)各向異性系數(shù)大于0.3時(shí),對S波應(yīng)謹(jǐn)慎使用。

近似計(jì)算公式雖然是在弱各向異性條件下推導(dǎo)得出的,但在各向異性系數(shù)達(dá)0.3時(shí),對一般各向異性介質(zhì)都適用,且具有較高的精確,為地震波速度分析及參數(shù)反演奠定了基礎(chǔ),對強(qiáng)各向異性介質(zhì)中地震波傳播特征近似研究是下一需要開展的工作。