唐杰，劉英昌，韓盛元，孫成禹

中國石油大學(華東)地球科學與技術學院，青島 266580

0 引言

波形反演是一種數據驅動的高分辨率模型參數估計方法，其強大之處在于它能夠利用記錄的相位和振幅信息，通過伴隨狀態(tài)法獲得目標函數的梯度(Fichtner et al.，2006；劉立斌等，2020)，進而通過優(yōu)化算法獲得目標模型參數.Tarantola(1984)指出數據殘差在時間上的反向傳播波場與正向傳播波場的互相關有助于迭代地生成高分辨率的速度模型.結合伴隨狀態(tài)法，波形反演技術只需兩次模擬就可以估計不同類型參數(不僅僅是速度參數)的梯度.雖然波形反演主要用于速度模型構建(孫史磊等，2020)，但微地震波形中包含的信息也可以用于反演震源參數等信息(Virieux and Operto，2009).針對各向異性介質條件下的震源參數反演，Michel和Tsvankin(2014)對彈性VTI介質中微震事件的震源位置和矩張量進行波形反演梯度計算.Grechka(2015)分析了單井微震資料全矩張量反演的可行性，認為在VTI介質中的垂直井中由于方位各向異性導致慢度矢量偏離包含射線軌跡的平面，通過處理能夠得到理論上唯一的地震矩張量.Michel和Tsvankin(2019)的研究表明三維各向異性波形反演可以為微地震勘探提供高分辨率的速度模型和改進的定位結果，聯合矩張量反演能夠取得較好的效果.Li 等(2021)提出了一種考慮VTI介質中速度各向異性的微震震源機制反演方法，通過對合成波形和觀測波形的匹配，僅用P波進行地震矩張量反演提高魯棒性，研究表明各向異性參數的偏差會引起非DC成分的反演誤差.這些研究針對的是彈性各向異性介質中的波形反演，地層的各向異性吸收衰減作用對震源機制反演的影響并沒有得到充分的研究，地層的吸收衰減會影響微震信號的波形、頻譜和信噪比等，從而降低震源機制的反演精度，矩張量反演結果誤差會導致不同屬性參數解釋結果出現不同程度的偏差，因此需要對黏彈性各向異性介質中的震源機制反演開展進一步研究.

本文研究了各向異性黏彈性VTI介質中的微地震震源機制反演方法，首先基于分數階拉普拉斯算子解耦的波動方程實現了黏彈性VTI介質矩張量震源的正演模擬，分析了衰減各向異性與震源機制對波場的影響；然后利用伴隨態(tài)法計算了非均勻VTI介質中的矩張量源的目標函數梯度，采用梯度下降算法更新了震源參數，獲得了震源機制矩張量反演結果，并通過實驗測試驗證了該方法的有效性.

1 原理

1.1 各向異性衰減VTI介質中微地震波場模擬

微地震的震源機制可以采用矩張量M的形式表示(Jost and Herrmann，1989)：

(1)

pM=pISO+pCLVD+pDC,

(2)

各成分的矩張量分量如表1所示，ISO成分由歸一化單位矩陣表示；CLVD成分(Lay and Wallace，1995)被用于近似水力壓裂過程中發(fā)生的拉伸事件，與拉伸裂紋不同的是，CLVD源不涉及任何體積變化，因為沿一個軸的應變會通過沿其他軸的收縮或膨脹進行補償(Knopoff and Randall，1970).

表1 本文采用的震源機制及其分量Table 1 Focal mechanisms and components used in this paper

VTI介質是具有垂向對稱軸的橫觀各向同性介質，其彈性系數矩陣可由C11、C13、C33、C55與C66表示.衰減各向異性也會影響地震波的振幅與相位的變化規(guī)律，通過品質因子矩陣能夠表征各向異性衰減，矩陣元素可以通過復剛度系數的實部與虛部的比值來確定.黏彈性VTI介質的Q矩陣滿足：

(3)

黏彈性各向異性介質中的應力-應變關系可以表示為：

σ(x,t)=ψ(x,t)*?tε(x,t),

(4)

式中σ=(σxx,σyy,σzz,σyz,σxz,σxy)T表示應力張量，ε=(εxx,εyy,εzz,εyz,εxz,εxy)T表示應變張量，x表示空間位置矢量，ψ表示松弛函數矩陣，矩陣中的每個分量可以表示為：

(5)

式中Γ表示伽馬(Gamma)函數，t0表示參考頻率ω0對應的參考時間，γij=arctan(1/Qij)/π表示與Qij相關的無量綱參數，H(t)表示階躍函數.

將各元素分量代入式(4)可以得到：

(6)

(7)

參考Caputo和Mainardi(1971)對分數階導數的定義，黏彈性介質應力-應變關系可以表示為：

(8)

其中，δ是克羅內克函數.將VTI介質彈性系數矩陣代入式(8)可以獲得時間分數階形式：

(9)

時間分數階導數的求解需要大量內存與計算時間，而通過近似可以將其轉化為由分數階拉普拉斯算子表征的空間分數階導數，同時波場也可以解耦為衰減項與頻散項，提高計算效率并有利于后續(xù)工作開展.在弱衰減的假設下，時間分數階導數可以由分數階拉普拉斯算子代替(Zhu and Harris，2014)，即：

(10)

則：

(11)

結合動量守恒公式與應變-位移關系可得到黏彈性VTI介質分數階拉普拉斯算子波動方程的速度-應力形式：

(12)

(12′)

通常正應力和剪應力是直接施加在相應的網格點上，而在交錯網格離散的偽譜法中正應力和剪應力不是在同一位置計算的.假設矩張量源作用于正應力分量的網格位置，則正應力分量的位置將作為中心點.為了得到一個對稱的力矩張量源，將剪應力增量均勻分布在真實矩張量源位置附近的四個應力網格點上.因此，在矩張量源的真實位置周圍總共有12個相鄰的網格點，這些網格點用剪應力分量數值實現.作用于交錯網格節(jié)點i、j、k(即正應力分量的網格位置)的矩張量源的完整公式如式(13)(Shi et al.，2018)：

(13)

其中V=Δx·Δy·Δz是網格單元的有效體積，Δt是時間間隔.由于應力分量的時間導數用于黏彈性動力學方程，所以使用矩張量的時間導數.

圖1 運算效率對比(a)計算時間；(b)加速比.Fig.1 Comparison of calculation efficiency(a)Calculation time;(b)Acceleration ratio.

1.2 矩張量波形反演

微地震震源矩張量反演方法有很多，包括采用全波形反演或集中考慮特定的體波相位(如P波和S波)的反演方法等.本文使用全波形反演來估計黏彈性VTI介質中的微地震震源的矩張量，用伴隨狀態(tài)法計算目標函數的梯度.

目標函數采用L2范數：

(14)

其中dpre(m)是模擬波場記錄，dobs是觀測記錄.采用伴隨狀態(tài)法計算目標函數相對于模型參數的導數?F(m)/?m(Michel and Tsvankin，2014)，根據Kim等(2011)可獲得目標函數相對于矩張量元素的導數：

(15)

其中xs是震源位置，T是記錄時間，ε+是伴隨應變張量，S是震源時間函數.矩張量M的梯度只需計算在震源位置的伴隨應變張量.

(16)

其中：

(17)

由于式(16)可以將波場解耦為衰減項和頻散項，所以在波場反傳時，保持頻散項符號不變，同時將衰減項對應的符號改變可以補償波場傳播過程中損失的能量，即將式(17)修改為：

(18)

為了克服反傳過程中高頻噪聲的影響，本文使用巴特沃茲低通濾波器在波數域對波場衰減項進行濾波處理(唐杰等，2021).對波場進行濾波時不應對傳播算子信號進行濾波(Zhu and Carcione，2014)，需要將傳播算子從頻散項中分離，并加入低通濾波器，式(18)變?yōu)椋?/p>

×[cos(πγij)(-?2)γij-(Cij/ρ)-0.5sin(πγij)

(19)

對于微地震震源機制反演的研究，評估和量化反演參數的不確定性也是非常重要的，震源機制反演標準化誤差由式(20)給出：

(20)

2 各向異性黏彈性VTI介質微地震正演

2.1 衰減各向異性對正演結果的影響

為了表征衰減的各向異性，Zhu和Tsvankin(2006)根據Thomsen參數(Thomsen，1986)的思路提出一組Thomsen型衰減參數(Thomsen-style notation for attenuation)，包括兩個參考各向同性參數(QP0、QS0)與三個無量綱參數(εQ、δQ和γQ).QP0與QS0為分別表征對稱(垂向)方向P波和SV波品質因子，參數εQ表征P波衰減各向異性強弱，作用與Thomsen參數中的ε類似，參數γQ表征SH波的衰減系數在水平方向與垂直方向的差異，作用與Thomsen參數中的γ類似，參數δQ在描述P波垂直方向衰減各向異性變化的快慢程度，與描述P波速度各向異性的δ相似.描述SV波衰減各向異性的σQ參數受Thomsen參數和Thomsen型衰減參數控制，即SV波衰減與ε、δ、εQ和δQ有關，σQ非獨立參數.εQ、δQ、γQ、σQ分別為：

(21)

設置一個均勻模型進行波場模擬，分析衰減各向異性對正演結果的影響，模型大小為1000 m×1000 m×1000 m，網格間距10 m，網格點數100×100×100，模型參數為：縱波速度為2.42 km·s-1，橫波速度為1.4 km·s-1；品質因子QP為35，QS為30；速度各向異性參數ε=0.4，δ=0.15，γ=0.1；衰減各向異性參數εQ=-0.125，δQ=0.94，γQ=-0.1.

隨角度變化的縱橫波速度如圖2a中的實線和虛線所示，P波衰減系數與SV波衰減系數隨角度變化分別如圖2b、c所示.VTI介質中縱橫波速度和衰減均呈現各向異性特征，圖2d為彈性波與黏彈性波場快照對比，震源采用位于模型中心的主頻為35 Hz的Ricker子波，時間步長為0.4 ms，震源采用爆炸源，圖中四部分分別為彈性波、僅考慮頻散項、僅考慮衰減項和黏彈性波，從圖2d中可以看出，相對于彈性波，僅考慮頻散項只有相位的延遲，僅考慮衰減項只有振幅的減小，而黏彈性介質中的波場既有相位的延遲又有振幅的減小.

圖2 各向異性介質參數分析(a)速度各向異性，其中實線為P波，虛線為SV波；(b)P波衰減各向異性；(c)SV波衰減各向異性；(d)波場快照.Fig.2 Analysis of media anisotropy parameters(a)Velocity anisotropy.Solid curve is P-wave,and dashed curve is SV wave;(b)P-wave attenuation anisotropy;(c)SV wave attenuation anisotropy;(d)Wave field snapshot.

分別進行無衰減、各向同性衰減和三個各向異性衰減模型(參數如圖3標注所示)正演模擬，獲得0.15 s時刻的垂直速度分量波場快照如圖3所示，其中圖3a為無衰減波場快照，圖3b為各向同性衰減波場快照，圖3c、d、e分別為三個各向異性衰減模型的波場快照，圖3f、g、h分別為各向異性衰減與圖3b各向同性衰減的差值.在距震源約280 m位置處觀測點獲得0.25 s波場記錄如圖4b，其對比曲線如圖4c.在各向同性介質中，爆炸源只產生P波，而從圖3a可以看出，在各向異性介質中爆炸源同樣會產生S波.除速度各向異性參數外，三組衰減各向異性參數均對波場具有一定影響.εQ表征P波在水平方向和垂直方向的衰減系數差異，0.12～0.17 s時間段曲線f存在起伏，即P波衰減受到各向異性影響，與各向同性衰減產生了差異；δQ為連接P波水平方向和垂直方向速度的過渡參數，也對P波衰減各向異性存在影響，在0.12～0.17 s的時間段曲線g存在起伏；而γQ表征SH波的衰減各向異性，在P波段曲線h基本無起伏，僅在0.18～0.25 s時間段的S波記錄中曲線h存在起伏.在各向異性衰減VTI介質的波場模擬中，應當綜合考慮多個各向異性參數的影響.

圖3 介質參數對波場的影響(虛線為截面位置)(a)Q =∞；(b)εQ=0，δQ=0，γQ=0；(c)εQ=-0.3，δQ=0，γQ=0；(d)εQ=0，δQ=0.4，γQ=0；(e)εQ=0，δQ=0，γQ=0.4；(f)(c)—(b)；(g)(d)—(b)；(h)(e)—(b).Fig.3 Influence of media parameters on wave field (the dashed lines are positions of sections)

圖4 介質參數對地震記錄的影響(a)觀測系統(tǒng)；(b)地震記錄；(c)地震記錄間差異.Fig.4 Influence of media parameters on seismic records(a)Observation system；(b)Seismic record；(c)Difference between seismic records.

2.2 不同類型震源機制正演結果

利用圖2所采用的模型參數，模擬不同類型震源在黏彈性VTI介質中的波場傳播，選用的震源機制如表1中所示，圖5a為MISO震源的結果，圖5b為MDC震源的結果，圖5c為MCLVD震源的結果，圖5d為混合類型震源MCOMBO的結果，不同類型震源計算得到的波場快照明顯不同.MISO源在各向同性介質當中只能激發(fā)出P波，但在各向異性介質當中既能激發(fā)出P波又能激發(fā)出S波.MISO源激發(fā)出的P波能量要明顯強于S波.非各向同性震源激發(fā)的信號中S波振幅和能量較大，考慮到微地震信號能量較弱，信噪比較低，P波可能會淹沒在噪聲中，在微地震監(jiān)測中應考慮S波的拾取與利用；不同機制震源激發(fā)的信號在不同方向強弱不同，極性也存在差異，在成本的限制下應盡量考慮多口不同方位井進行多分量信號采集.

圖5 均勻黏彈性VTI介質中不同類型源的波場快照(虛線為截面位置)(a)MISO；(b)MDC；(c)MCLVD；(d)MCOMBO.Fig.5 Wavefield snaps of different sources in a homogeneous viscoelastic VTI media(dashed lines are the positions of sections)

考慮斷層面與[x，z]坐標平面重合的位錯源，根據Aki和Richards(2002)提及的方法，若介質為VTI介質，則點位錯源的矩張量可由式(22)計算得到：

(22)

圖6給出了不同θ角度的純剪切源的震源機制沙灘球與鉆石圖，當θ=30°或者θ=60°時，斷層面不是位于VTI介質對稱平面內，由圖可知純剪切錯動在各向異性介質中仍然可能產生非DC成分，見表1的結果.當θ=0°或者θ=90°時，斷層面位于VTI介質對稱平面內，純剪切源不會產生非DC成分，仍然為雙力偶機制.

利用圖2所采用的模型參數，加載圖6中的位錯源，圖7給出了均勻黏彈性VTI介質中不同方向(由角度θ定義)位錯源的垂直速度分量波場快照，不同的θ角度會導致振幅大小和極性分布發(fā)生明顯的變化.背景介質各向異性強弱和斷層面的幾何位置等參數都會影響非DC成分，這種非DC成分也能夠反應震源區(qū)的各向異性信息.

圖6 VTI介質不同θ角度的位錯源震源機制沙灘球與鉆石圖(a)θ=0°；(b)θ=30°；(c)θ=60°；(d)θ=90°.Fig.6 Focal mechanism beach balls and diamond plots of dislocation sources with different θ angles in VTI media

圖7 均勻各向異性黏彈性VTI介質中不同方向(由角度θ定義)位錯源的vz分量(虛線為截面位置)(a)M0°；(b)M30°；(c)M60°；(d)M90°.Fig.7 Component vz of different direction (defined by angle θ)dislocation sources in homogeneous viscoelastic VTI media (dashed lines are positions of sections)

基于常Q模型的分數階拉普拉斯本構方程給出了一階速度-應力方程中相應的包含震源機制影響的時域黏彈性VTI介質波動方程，通過波場模擬可以發(fā)現介質的衰減各向異性以及震源機制將使記錄的波場更為復雜，速度各向異性會對記錄的不同方位的地震波走時產生影響，衰減各向異性會對記錄的地震波的行程時間(頻散的影響)和振幅(衰減的影響)等參數產生影響，從微震數據中分析衰減各向異性有助于表征水力壓裂中的裂縫特性；由速度和衰減各向異性引起的行程時間、幅度以及波形的變化可用于評估地下介質各向異性的強度.

3 各向異性黏彈性VTI介質震源機制反演

3.1 速度與衰減各向異性參數對反演的影響

為了驗證各向異性黏彈性介質中震源機制全波形反演方法的有效性，建立三層黏彈性VTI介質模型，如圖8a所示，模型尺寸為1000 m×1000 m×1000 m，各層的參數如表2所示，第一層和第三層為各向同性衰減的VTI介質，代表儲層的上覆和下伏地層，由于不含孔隙流體，其衰減呈現各向同性，第二層為各向異性衰減的VTI介質，代表含孔隙流體的頁巖地層，在該層進行水力壓裂.該模型采用10 m的網格間距進行離散，震源位于模型中間位置，源函數是主頻35 Hz的Ricker子波，觀測系統(tǒng)如圖8b，在震源周圍設置三口垂直井，井中每50 m放置一個檢波器，共60個檢波器.通過正演模擬可以獲得檢波器位置的三分量微地震數據，不含噪聲的記錄如圖9所示.

圖8 (a)模型；(b)觀測系統(tǒng)Fig.8 (a)Model；(b)Observation system

圖9 井中地震記錄(a)vz；(b)vy；(c)vx.Fig.9 Borehole seismic record

表2 模型參數Table 2 Model parameters

采用全波形反演方法可以獲得黏彈性各向異性介質中的矩張量，反演時對第二層介質采用的模型參數如表3所示，其中A8為對模型介質參數進行高斯平滑的情況，因此震源區(qū)參數未發(fā)生改變.采用的震源機制MCOMBO如圖10a所示，震源機制各分量見表1，其震源機制成分Lune圖如圖10b所示(Tape and Tape，2013)，文中矩張量反演初始值統(tǒng)一使用隨機數組Minitial，見表1.

表3 反演采用的模型參數Table 3 Model parameters used in inversion

圖10 本文采用的震源機制(a)沙灘球；(b)Lune圖.Fig.10 Focal mechanism used in the paper(a)Beach ball;(b)Lune graph.

圖11給出了采用不同介質參數反演目標函數隨著迭代次數的變化情況，隨著迭代次數的增加，目標函數都呈下降趨勢，其中A1采用黏彈性各向異性介質模型，目標函數迭代后值最小，說明模擬波形與觀測波形匹配得最好；A2采用彈性各向異性介質模型，忽略了衰減的影響，模擬波形的振幅與觀測波形的振幅會存在差異，導致目標函數值迭代后相對較大；存在其他參數誤差時，目標函數值迭代后相對無干擾時均有一定程度的升高.

圖11 采用不同介質參數條件下的目標函數隨迭代次數的變化Fig.11 Variation of the objective function with the number of iterations under different media parameters conditions

圖12a為采用不同介質參數得到的反演矩張量的各個分量，圖12b為分解得到的不同成分.由圖可知采用黏彈性各向異性介質參數能夠得到正確的反演結果，而速度或者Q值存在誤差時會導致反演結果出現較大的偏差，衰減各向異性誤差也會影響反演結果.圖12c為不同介質參數條件下反演得到的震源機制沙灘球以及震源機制偏差，橫軸為選取不同介質參數的情況，縱軸為反演結果的偏差.由圖可知采用黏彈性各向異性介質模型能夠獲得正確的反演結果，本文方法在介質參數正確的情況下能夠準確反演震源機制；采用彈性各向同性介質模型或者速度存在誤差時，震源機制的反演結果會存在誤差，結合圖2及圖3的分析可知ε、δ和γ控制水平縱波速度和橫波速度偏離垂直方向的變化，這些參數的誤差可能會改變某些接收器記錄信號的到達時間和振幅，速度以及各向異性參數都會影響反演震源機制各分量的百分比.此外模型誤差(如速度結構、衰減等相關的誤差)的影響會隨著距離的增加而增大，從而對反演結果產生不利影響，采用距離震源區(qū)較近的觀測井數據能夠減少此誤差.

圖12 不同介質參數時震源機制反演結果(a)矩張量；(b)成分；(c)震源機制沙灘球.Fig.12 Inversion results of focal mechanism with different media parameters(a)Moment tensor;(b)Component;(c)Focal mechanism beach ball.

3.2 噪聲的影響

觀測數據中的高斯噪聲會在地震記錄中造成嚴重的扭曲，噪聲主要影響振幅，進而會影響矩張量的反演結果.為了分析數據中的噪聲對微震震源機制反演的影響，在合成波形中加入不同信噪比的噪聲.圖13給出了添加不同噪聲水平的垂直分量地震記錄和單道三分量記錄，當信噪比低于1時，地震記錄受到大量高頻成分噪聲的干擾，圖14為不同信噪比下反演目標函數迭代變化，存在噪聲時合成記錄無法與含噪信號完好匹配，伴隨應變張量受到影響，目標函數值出現升高.

圖13 不同信噪比條件下垂直分量地震記錄及三分量單道記錄對比(a)SNR=5；(b)SNR=1；(c)SNR=0.5；(d)SNR=0.1；(e)vz；(f)vy；(g)vx.Fig.13 Comparison of vertical component and three-component single-trace records under different SNR conditions

圖14 不同信噪比條件下目標函數隨迭代次數的變化Fig.14 Variation of the objective function with the number of iterations under different SNR conditions

圖15a為采用不同噪聲水平得到的反演矩張量的各個分量，圖15b為分解得到的不同成分.圖15c給出了不同信噪比條件下的震源機制反演結果，橫軸為選取不同信噪比的情況，縱軸為反演結果的偏差.反演時采用了黏彈性各向異性介質模型，由圖可知全波形反演方法具有一定程度的抗噪性.信噪比較高時盡管振幅值中增加了噪聲，但結果仍具有較高的精度.如果使用較少的觀測數據或者信噪比過低時，誤差也會相應增大.

圖15 不同噪聲水平下震源機制反演結果(a)矩張量；(b)成分；(c)震源機制沙灘球.Fig.15 Inversion results of focal mechanism with different levels of noise(a)Moment tensor;(b)Component;(c)Focal mechanism beach ball.

3.3 觀測系統(tǒng)的影響

三維波形反演應用于微地震數據的一個主要挑戰(zhàn)是有限的源-接收器覆蓋范圍(數據孔徑)，檢波器分布對矩張量反演結果的可靠性有很大影響，因此微地震監(jiān)測項目的設計應考慮到這一點.表面陣列有一個更大的孔徑(即大偏移距深度比)，部署成本通常也更低，能夠采集到多方位的地震信號，但由于微地震信號能量弱、衰減快，地表監(jiān)測的主要問題是相對井監(jiān)測記錄信噪比較低，本節(jié)只討論井監(jiān)測觀測系統(tǒng)的影響.井下觀測通常使用三分量傳感器，能夠獲得更高信噪比的微地震信號，這將提高反演結果的可靠性.但是如果需要鉆新的井，則井下監(jiān)測可能比地表監(jiān)測昂貴得多，另外由于大多數井的一維特性，它們比設計良好的表面陣列采集的焦球信息更少，反演結果更容易受震源位置與井的關系的影響.本節(jié)采用不同數目、方位垂直井和不同檢波器進行反演，分析井下觀測系統(tǒng)對微地震震源機制反演的影響.

圖16 不同觀測系統(tǒng)下矩張量各分量反演結果(a)走滑斷層；(b)傾滑斷層；(c)正斷層；(d)逆斷層.Fig.16 Inversion results of moment tensor components with different observation systems(a)Strike-slip fault;(b)Dip-slip fault;(c)Normal fault;(d)Reverse fault.

圖17給出了不同檢波器分布條件下的震源機制反演結果，橫軸為選取不同震源機制的情況，縱軸為反演結果的偏差.由圖可知對于不同的震源機制單井三分量和單井單分量反演得到的偏差較大，震源機制沙灘球與真實值存在較大的差異.三井三分量記錄反演的偏差高于三井單分量，其原因為噪聲和Q值誤差綜合干擾，導致三分量記錄受到的影響高于單分量記錄，若僅考慮噪聲影響時，三分量記錄抗噪性高于單分量記錄.在垂直橫向各向同性地層中鉆取的垂直井中，利用單井微震觀測資料，對于地震矩張量的六個獨立分量不能唯一地進行反演.通過假設震源的某些物理特性使原本獨立的矩張量分量相互關聯(如拉伸斷裂假設等)，使矩張量反演正則化，有助于減少反演結果的模糊性(Grechka et al.，2016).

圖17 不同觀測系統(tǒng)下矩張量反演結果沙灘球(a)走滑斷層；(b)傾滑斷層；(c)正斷層；(d)逆斷層.Fig.17 Beach balls of moment tensor inversion results with different observation systems(a)Strike-slip fault;(b)Dip-slip fault;(c)Normal fault;(d)Reverse fault.

如果微地震事件靠近監(jiān)測井，如圖18a中1號井和2號井，由于能夠采集到的更多的焦球信息，如圖18b所示通過全波形反演可以從單井數據集中獲得完整的矩張量(如1號井和2號井)；當觀測井位于遠場范圍時(如3號井和4號井)，僅靠單井資料無法完成矩張量反演，因此獲得足夠的震源輻射花樣信息對于震源機制反演具有重要的意義(Song and Toks?z，2011).

圖18 (a)觀測系統(tǒng)；(b)不同方位井記錄下矩張量反演結果沙灘球Fig.18 (a)Observation system；(b)Beach balls of moment tensor inversion results with different azimuth borehole records

4 結論

本文提出了黏彈性VTI介質中震源矩張量反演方法.首先推導了黏彈性VTI介質波動方程并采用GPU并行偽譜法正演模擬微地震記錄，然后通過全波形反演算法獲得了微地震震源機制，并通過實驗測試驗證了該方法的有效性，得到的主要結論包括：

(1)頁巖油氣儲層通常具有各向異性或者由于水力壓裂而產生各向異性，基于常Q模型和分數階拉普拉斯算子推導的黏彈性VTI介質波動方程能夠描述速度各向異性和衰減各向異性，更符合地下介質的實際情況.正演結果表明介質的衰減各向異性以及震源機制將使記錄波形更為復雜，分析衰減各向異性有助于表征水力壓裂中的裂縫特性.

(2)震源機制全波形反演需要獲得伴隨波場計算矩張量梯度，由于介質具有各向異性和黏彈性，在計算伴隨波場時需要進行衰減補償，本文推導的波動方程能夠在反傳過程對波場進行各向異性衰減補償，提高矩張量梯度計算的準確性.

(3)全波形反演方法可以有效反演黏彈性VTI介質中的震源矩張量.數值測試證明采用黏彈性各向異性介質模型能夠獲得正確的反演結果，將黏彈性各向異性添加到矩張量反演中可以幫助詳細研究非DC成分是否真的存在或是由計算誤差引起的.反演過程中的模型參數誤差、觀測數據噪聲和觀測系統(tǒng)等也會對反演結果造成影響，獲取良好的地層介質參數和設計合理的觀測系統(tǒng)能夠更準確地反演震源機制以評價壓裂效果.

除黏彈性各向異性外，微地震震源機制的反演還涉及許多其他問題，如震源定位結果、震源時間函數和模型網格大小等因素的影響，本文不做進一步展開.由于波形反演效果依賴背景介質模型的可靠性，為提高抗干擾能力可以在反演過程中同時對震源機制和模型參數聯合反演，以減少模型參數不準確帶來的誤差，但這也會帶來參數串擾問題以及更大的計算量，有待進一步研究.