反演近地表物性參數的地震層析成像方法綜述

劉玉柱,吳世林,劉偉剛,黃鑫泉,伍 正

(1.同濟大學海洋地質國家重點實驗室,上海200092;2.同濟大學海洋與地球科學學院,上海200092)

1 層析算法

基于射線理論的地震走時層析成像在工業界獲得了成功應用,現在仍是商業軟件的核心建模模塊。其基本原理建立在對正問題的一階Taylor級數展開的基礎上:

t=f(s)≈f(s0)+f′(s0)Δs

(1)

式中:t為射線的走時;f(s)為表達走時與慢度s之間非線性關系的函數;s0為背景慢度。在射線理論框架下,正問題可以用矩陣公式表達為:

f(s)=L(s)·s

(2)

式中:L(s)表示慢度模型s下的射線路徑長度矩陣。根據Fermat原理,射線路徑是相對穩定的,因此:

f′(s0)=L(s0)

(3)

將(3)式代入(1)式即得傳統的射線走時層析方程組:

L(s0)Δs=Δt

(4)

式中:Δs為慢度擾動向量;Δt為真實走時與理論合成走時之間的走時殘差向量。傳統層析成像通過求解上述方程組實現局部慢度向量的更新。如此,利用局部線性化+迭代方式實現非線性走時層析反演。

受計算能力的限制,20世紀70—80年代以來,地震層析成像一直采用經驗性的求解大規模線性方程組(4)的方式進行模型的更新[1]。典型的算法包括反投影法(BPT)、代數重建法(ART)和聯合代數重建法(SIRT)等。這類方法基于加權平均的方式對模型進行更新,計算效率高但反演精度低。另一類反演類算法包括最小二乘QR分解(LSQR)和截斷奇異值分解(SVD)等,基于奇異矩陣的廣義逆算法,反演效果更好,但計算量比較大,算法不穩定[2]。但無論經驗性解法還是廣義逆解法都需要存儲層析矩陣,對計算機內存要求很高。

借鑒全波形反演的實現過程[3],在目標函數的局部最優化基礎上[4]建立了一種走時層析成像實現方式,該方式通過尋找某一慢度向量s*,使得如下目標函數極小: