岳 漢
北京大學地球與空間科學學院,北京 100871
大地震依然是造成最多人類傷亡的自然災害.一次大地震發生后,解析其破裂過程是指導震后救災的重要內容,也是實現長期地震危險性研究的重要基礎研究支撐(Yue et al., 2020). Yue 等(2020)討論了震后快速響應以及破裂過程聯合反演在震后應急中的重要性,集中分析了多數據聯合反演對于破裂過程解析度的互補以及各個方法在大地震快速響應中的應用,而對于“復雜破裂過程”的破裂過程反演以及應急響應,現有研究框架尚未完全包括. 本文針對復雜破裂過程研究中的靈活性需求進行簡要概述,并提出結合反投影技術和多點源反演算法研究復雜地震破裂過程的技術路線.
大地震破裂過程研究的主要思路是:
(1)假設破裂斷層面為矩形,通過震源機制解獲得斷層面幾何參數,并通過震級估計其尺度.
(2)參考震源機制以及區域構造特征選擇合適的斷層面.
(3)搜集相關數據反演其破裂過程.
這種平面斷層的假設在很多情況下可以滿足破裂過程快速響應的需求,但是在破裂過程復雜的很多大陸內部強震中,平面斷層假設無法滿足解析破裂過程的需求. 現代地震學和大地測量測數據顯示,很多大地震地表破裂往往不局限于單一斷層面而包括若干條斷層的先后或者同時破裂(如 Sun et al.2018; Wang et al., 2018). 這些現象并不是因為地震變得更加復雜,而是因為觀測讓我們能夠解析到更多的地震破裂細節. 例如2012 年印度洋底地震(Meng et al., 2012; Yue et al., 2012)、2016 年新西蘭Kaikoura 地震(Duputel and Rivera, 2017; Wang et al., 2018)、2013 年巴基斯坦地震(Avouac et al., 2014)、2016 年日本熊本地震(Yue et al.,2017)、2019 年美國加州Ridgecrest 地震(Wang et al., 2020; Yue et al., 2021)等,這些地震在斷層幾何形態以及分段特征上都發生了明顯改變,破裂面都不限于單一平面斷層. 對于若干中等規模地震,如6~7 級地震,科研人員也發現了地震的復雜性,如2014 年云南魯甸地震(劉成利,2013,2014)、2017 年四川九寨溝地震(Sun et al., 2018)、2018年臺灣花蓮地震等(Lo et al., 2019). 應用了點源解的單一平面來參數化斷層面,往往讓這些地震的有限斷層反演無法準確表達地震的破裂機制,為地震危險性的快速估計帶來困難.
對于大陸內部復雜地震而言,研究人員通??梢詤⒖嫉刭|考察、衛星影響以及余震分布刻畫非平面斷層. 地表勘察以及衛星影像數據往往可以得到比較精細的地表斷層形態. 而精定位后的余震目錄往往可以在深度分布上刻畫斷層形態. 利用地表和余震數據,科研人員可以刻畫出精細的斷層形態,作為破裂過程反演的參考模型. 例如,2016 年日本熊本地震(Yue et al., 2017)、2019 年美國加州Ridgecrest 地震(Ross et al., 2019)等,都有研究人員通過上述方法獲得了多段非平面斷層面. 然而,衛星數據的獲得以及余震目錄的獲取往往需要若干天的時間,無法滿足快速應急的需求. 此外,發生在海中的地震,往往沒有足夠精細的表層觀測以及余震目錄實現上述斷層刻畫過程. 因此對于復雜地震的破裂過程,我們尚需要一套可以在全球范圍內普遍適用的、可以反演破裂在多斷層破裂上的發展過程的、靈活度較高的反演方法. 下文結合大地震研究中四種常用的方法進行簡要介紹.
對于復雜破裂過程,基于遠震地震波的多種反演技術為非單一平面的破裂過程反演提供了依據. 本文列舉了四種基于遠震地震波的重要方法,分別為單點源震源機制解、多點源震源機制解、反投影技術以及有限斷層模型. 對于一次大地震的復雜破裂過程,上述四種方法都可以發揮作用. 如圖1 所示是針對2017 年的Kaikoura 地震,應用幾個不同研究方法的成果. 這些方法都反映了地震東北向的傳播特征以及主要的兩次能量輻射. 這些研究結果都可以通過遠震地震波得到. 遠震地震波是震源快速響應的重要數據. 因為其具有全球可得的特點,可以形成一套全球地震普遍適用的方法體系,為全球范圍內復雜地震的快速響應提供了可能. 這些基于遠震地震波的方法中,有限斷層模型反映了最為準確的破裂細節,但是由于有限斷層模型需要先將斷層面參數化,對于復雜破裂過程而言,無法通過單點源解的最佳破裂面得到有效的斷層結構. 多點源方法和反投影技術(上游技術)為有限斷層模型(下游技術)的參數化提供了依據,但是這些方法之間的銜接目前依然依賴于研究人員的經驗和手動調試,這一過程通常比較耗時,且無法在快速響應中應用. 因此實現復雜地震快速響應的瓶頸不在數據層面,而在技術層面.復雜地震破裂過程的研究迫切需要一種可以將上游結果作為下游輸入參數的自動參數化手段,實現不同反演技術之間的銜接,實現整個流程自動化.

圖 1 2016 年新西蘭Kaikoura 地震的多種震源研究方法比較. (a)GCMT 和Wphase 的點源反演結果;(b)多點源反演結果(修改自Duputel and Rivera, 2017);(c)反投影成像結果(修改自Hollingsworth et al., 2017);(d)有限斷層反演結果(修改自Wang et al., 2018)Fig. 1 Comparison between different inversion results for the 2016 Kaikoura earthquake. (a) GCMT and Wphase point source solution; (b) Multi-point-source solution (modified from Duputel and Rivera, 2017); (c) Back-projection (modified from Hollingsworth et al., 2017); (d) Finite fault model inversion (modified from Wang et al., 2018)
上述四種方法在學界內還處于大量單獨應用的“單兵作戰”模式. 而在單獨應用各個方法時通常因為方法本身的缺陷而無法了解地震的全面信息.下面對于各個方法的優缺點進行簡要的介紹. 由于單點源、反投影、有限斷層等方法的相關文章較多,本文只做簡要介紹. 而對于多點源方法,由于其滿足了震源研究靈活性的需求,本文將做重點介紹.
數學上可以證明,發生在斷層面上的點位錯源產生的位移場可以等效表示為雙力偶源產生的位移場(震源表示定理). 因此用空間中集中于一點的雙力偶源描述震源模型(震源機制解)成為描述地震破裂模型的基本形式. 震源機制解將震源描述為一個無限小的點源,描述地震斷層面的幾何特征,不包括尺度范圍,也可以稱為單點源解. 目前,基于遠震波形的單點源反演,如GCMT 點源解、(Ekstr?m et al., 2012)以 及Wphase(Kanamori and Rivera, 2008; Duputel et al., 2012b),已經自動成為監測全球大地震機制的重要手段. 不同于GCMT 和Wphase 的反演原理,Zhu 和Helmberger(1996)提出了用最大化區域震相時移后的相關函數的方法,建立了搜索震源機制的Cut-and-Paste(CAP)方法,對于中小型地震的區域數據反演震源機制解, 并在大量區域地震的研究中得到應用(如Zheng et al., 2009; Zheng et al., 2010; Zhao et al., 2013). 除了以波形為基礎的反演方法之外,利用P 波初動、P/S 振幅比等觀測方法也一直在大地震點源解的反演中發揮著重要的作用(如Dreger and Helmberger, 1993; 胡新亮等, 2004; 胡幸平等,2008). 在時間域單點源反演中通常不將破裂過程近似為脈沖源,而一般假設震源時間函數近似為一個三角形,通過搜索三角形時長的方法估計震源持續時間,如GCMT. 在有經驗格林函數的情況下,可以假設大地震和參考地震發生的位置相同,然后通過反卷積的方法提取比較精細的震源時間函數(Mueller, 1985).
N
方差(Xu Y et al., 2009)、Multi-taper(Meng et al., 2011)等技術的應用使得反投影技術可以快速自動化地成像大地震破裂過程. 不同于震源反演技術,反投影技術應用遠震地震波的高頻信息,不需要對破裂過程做出先驗假設,可以快速得到破裂過程的時空分布. 結合理論和經驗的到時校正,其對于破裂過程中能量位置的解析度可以達到10~20 km 的量級(Yue et al., 2020).反投影疊加得到的地震能量對于震源空間相關性比較敏感,并不與地震矩直接相關. 因此反投影技術具有快速成像、空間解析度高的優勢,但是對于地震矩、震源機制沒有解析度.N
個點源的目錄就比單點源多N
倍的信息.Kikuchi 和Kanamori(1982, 1986, 1991)自1980年代開始研發基于遠震體波的多點源反演算法,利用時間域迭代反演算法(剝皮法)實現了大地震的多點源反演. 本世紀的大地震記錄給多點源更多的發展空間,基于不同反演算法的多點源反演工作也在多個地震中得到了應用. 2004 年蘇門答臘地震后,Tsai 等(2005)借鑒GCMT 反演的方法研發出multi-CMT 反演,在限定地震位置的前提下反演其機制和地震矩. 對于這次9.0+級的地震,GCMT 使用的地震波頻寬已經達到飽和,無法正確獲得地震矩. GCMT 地震矩為9.0,明顯小于地球簡正模(Park et al., 2005)以及大地測量數據(Banerjee et al., 2005)給出 的 地震 矩(9.1~9.3). Tsai 等(2005)將蘇門答臘地震劃分成5 個點源,通過多次面波進行反演,反演得到的地震矩達到9.3,震源破裂時間為500 s,更符合區域數據得到的破裂過程. 2012 年印度洋底地震后, Duputel 等(2012)應用馬爾科夫鏈-蒙特卡羅(Markov Chain Monte Carlo, MCMC)方法,將這次地震分解為兩個點源,并且反演了點源的震級、位置和深度. 其得到的破裂模型與有限斷層和反投影技術得到的主要斷層破裂一致(如Meng et al., 2012; Yue et al., 2012).2016 年Kaikoura 地震后,Duputel 和Rivera (2017)利用相似技術,應用全球三維速度模型計算的波形格林函數,用四個點源模擬了破裂過程. 該模型表示,破裂后期的主要能量釋放階段,包含了兩個時間、空間和震級相近、而機制不同的逆沖和走滑分量的地震(M
7.6),分別對應了Kekerengu 斷層以及俯沖板塊上的破裂. 2018 年斐濟深震后,Jia 等(2019)通過MCMC 方法,解析了發生于毗連的兩個板塊上的兩個點源解. 對于中等規模的地震(M
6~7),多點源反演同樣能發揮作用. Shi 等(2018)研發了基于馬爾科夫鏈-蒙特卡羅(MCMC)的算法,對于2016 年6.2 級熊本地震前震進行了多點源反演. 對于中等規模的地震(6~7 級),作者應用時間域迭代反演算法,并且結合遠震、近震數據對于2017 年九寨溝地震(Sun et al., 2018)和2018 年花蓮地震(Lo et al., 2019)進行了多點源反演. 在九寨溝地震中,Sun 等(2018)通過波形初動分析以及多點源反演的方法,發現了破裂初期包括了與主斷層共軛的M
5.6 的逆沖型地震(圖2).在2018 年花蓮地震的研究中,Lo 等(2018)結合遠震體波以及近震波形,用6 個點源解釋這次地震的破裂過程(圖3),得出了從低角度走滑+逆沖轉化為高角度走滑、最終到垂直走滑的破裂機制,與通過大地測量學以及地震學手段得到的多斷層破裂模型符合較好. 該研究顯示花蓮地震涉及了外海板間斷層,以及板內走滑斷層(美倫斷層與嶺頂斷層)的先后破裂.
上述震例中,在斷層幾何形態不明朗的情況下,多點源方法作為具有較高自由度的方法,通過多個震源機制解反映破裂過程中斷層面幾何形態的變化. 在這些震例中,多點源算法的共同優勢是可以提高反演的自由度,解析復雜的破裂過程. 自由度高的優勢也帶來穩定性降低的不足. Yue 和 Lay(2020)針對反演穩定性提出了基于時間域反演算法的先驗信息的應用,通過限制點源時間函數、位置以及機制組合等策略增加對于單點源時空以及機制的約束,如圖4 所示.
Yue 和Lay(2020) 將改進的多點源反演方法應用于過去10 年內發生的7 個大地震中, 通過與有限斷層模型的對比,系統地討論了基于遠震數據的多點源反演的時空解析度以及反演結果穩定性問題. 并指出,基于遠震體波的多點源反演的空間解析度上限在20 km 左右,該誤差范圍對于參數化斷層模型依然比較粗糙. 因此,單純基于多點源反演的結果很難給出準確的破裂面幾何形態. Sun 等(2018)與Lo 等(2019)結合了大地測量數據、區域地震目錄以及多點源反演結果,構造了包含多個斷層面的斷層破裂模型,但是這些分析強烈依賴于大地測量數據,時效性不高,無法在地震快速響應方面取得突破.

圖 2 (a)2017 年九寨溝地震的發震位置,其南側分別發生了1973 年黃龍地震、1976 年的松潘震群. (b)2017 年九寨溝地震多點源反演. GCMT 點源解表示為大的震源球,左側4 個遠震體波波形以及單點源和多點源擬合的波形表示為黑色、藍色和紅色的波形,多點源反演的結果呈現于右下方,多點源疊加的震源機制解和有限斷層疊加的震源機制解畫在右上方(修改自Sun et al., 2018)Fig. 2 (a) 2017 Jiuzhaigou earthquake location. 1973 Huanglong earthquake and 1976 Songpan earthquake occurred on the south Huya fault. (b) Multi-point source inversion results of the 2017 Jiuzhaigou earthquake. Example waveforms of four southern teleseismic stations are plotted. Waveforms fitted by single point and multi-point solutions are plotted in blue and red, respectively. Solutions of 4 point sources of the MPS inversion result are plotted in the bottom (modified from Sun et al., 2018)

圖 3 (a)2018 年花蓮地震的發生位置以及區域構造. 花蓮地震發生于歐亞和菲律賓板塊的交接處,是造山和俯沖構造的交界點,破裂過程復雜. (b)花蓮地震多點源反演結果,不同顏色的震源球表示了遠震體波和近場三分量地震波形反演的多點源解,震源球的顏色代表發震斷層:紫色為板間斷層,藍色為美倫斷層,綠色為嶺頂斷層. 每個點源發生的時間、震級以及深度在右下角呈現(修改自Lo et al., 2019)Fig. 3 (a) 2018 Hualien earthquake rupture. (b) MPS inversion results of the 2018 Hualien earthquake. Point sources on the off-shore fault, Meilun fault and Linding faults are plotted in magenta, blue and green focal mechanisms, respectively (modified from Lo et al., 2019)
不同于震源機制解,有限斷層模型注重于描述斷層面的空間尺度,斷層面的幾何形態也不限制于單一平面. 如果有限斷層的反演中包含動態信息,如遠震體波、面波、區域震相等,有限斷層模型通常也可以描述地震破裂過程. 相對于震源機制解而言,有限斷層模型更加接近于真實的地震破裂模型,與野外地質以及大地測量觀測更加相關,同時也更加有助于強地面震動的模擬以及地震災害的估計.自從有了有限斷層模型的先驅性工作(Trifunac,1974),經過地震學家、大地測量學家數十年的努力,有限斷層的反演方法已經日趨成熟,并逐漸成為大家接受的并且廣泛使用的震源表示方法. 目前快速有限斷層反演在我國以及美國都已投入使用.有限斷層的快速反演通??梢栽谌舾尚r之內對于大地震的破裂過程給出第一手的模型,為救災工作提供參考(張培震等, 2009; Hayes, 2017).
快速有限斷層反演一般依賴于地震的點源解預設斷層面,然后通過遠震體波和面波進行有限斷層的反演. 對于大多數地震,單一斷層面的有限斷層反演可以囊括地震的主要破裂特征. 然而,現在地震學以及大地測量學發現在一些大地震破裂過程中,可以發生多個斷層的先后破裂,而破裂在多個斷層之間跳躍的過程中可能發生震源機制的變化. 過去幾十年的觀測已經積累了一些地震在多斷層之間跳躍的破裂過程的案例. 比如,2008 年汶川地震破裂了龍門山破裂帶的多條斷層,震源機制經歷了逆沖到走滑斷層的變化(Wang et al., 2008; Shen et al.,2009; Xu X et al., 2009; Zhang et al., 2009; Zheng et al., 2009);2012 年印度洋底8.6 級地震在空間上包含了5~6 個走滑型斷層的先后破裂,破裂范圍超過400 km(Meng et al., 2012; Yue et al., 2012);2016 年新西蘭Kaikoura 地震產生了至少10 個斷層的地表破裂,同時也造成了俯沖板塊表面的破裂(Wang et al., 2018; Xu et al., 2018). 這些地震的復雜破裂過程無法用單一斷層面的破裂模型進行解釋,用單斷層的有限斷層反演算法將會對破裂過程造成系統偏差.

圖 4 (a)第一行的左右分別為生成模擬波形的點源位置以及時間函數. 第二和第三行分別為有先驗約束的點源反演以及無先驗約束的點源反演的多點源機制以及時間函數. (b)模擬波形以及每個點源產生的模擬波形. 圖(b)左右分別為有無先驗約束的波形擬合(修改自Yue and Lay, 2020)Fig. 4 (a) Mechanisms and source time functions of synthetic earthquakes are plotted in the left and right panels of the top row,respectively. MPS inversion results with and without a priori constraints are plotted in the second and third rows. (b) Waveforms fit by either MPS inversions are plotted in the bottom panels. Synthetic and modeled waveforms are plotted in black and red curves in the top panel. Synthetic waveforms of each point source subevent are plotted as color coded waveforms at the bottom (modified from Yue and Lay, 2020)
相對于日本、智利等國家的俯沖帶地震而言,我國的地震災害主要來自于大陸內部的板內地震,更適合現代地球物理學觀測的解析能力. 這一特征適應了現代地球物理學以及大地測量學的觀測趨勢,也給我國的地震學研究提出了挑戰. 對于大陸內部的復雜板內地震,依然需要研發快速的、可以廣泛推廣的斷層面劃分算法,在地震破裂模型的快速反演、提高震災估計以及救援的時效方面提供有效的工具.
在復雜地震破裂過程的反演工作中,大地測量數據、地震地質資料、地震波波形分析等方法通常是劃分復雜斷層幾何的依據,而地震的快速響應幾乎必須依賴波形數據. 上述三種反演方法(多點源、反投影、有限斷層)都可以基于遠震體波進行處理, 在解析破裂過程中各有優缺點,在功能上可以互補. 三種方法的優缺點如表1所示.
總體來說, 先驗約束是在反演穩定性不高的時候通過經驗或者其他信息縮小參數空間的策略,因此反演穩定性越高的算法往往需要的先驗條件越少,但是這些方法對于破裂的多方面特征(如地震矩、震源機機制等)往往解析度不完備. 反投影方法需要的預設信息最少,通常只需要給出震源區的范圍以及網格信息,就可以對任意復雜的破裂過程給出高時間和空間精度的成像. 但是反投影技術應用高頻信息疊加而成,疊加方法也通常不是線性方法,無法恢復地震矩以及震源機制信息. 多點源反演在強的先驗模型下才能得出比較穩定的破裂模型,但是對于點源的時間和空間解析度較差. 有限斷層模型在所有方法中對于先驗約束的依賴性最強,通常需要對斷層面幾何形態以及破裂運動學過程給出比較完備的描述,其反演結果的時空解析度以及地震矩的解析度也是最高的.
從表1 的方法中可以看出,從上游到下游(反投影>多點源>有限斷層)遵循先驗約束逐漸增加,進而模型解析度逐漸增加的特征. 在邏輯上比較合理的策略是將上游方法的結論作為先驗約束輸入下游方法,從而約束下游方法反演結果的穩定性. 在一些綜合分析中(如Yue et al., 2012; Sun et al., 2018; Lo et al., 2019),也依賴研究人員的主觀判斷,將上游的信息轉換成為下游方法輸入參數.因此,將上述方法串行需要更加客觀的標準以及算法,才能夠實現快速、普適的震源模型參數化策略.

表 1 反投影、多點源以及有限斷層模型反演方法的先驗約束以及解析度分析Table 1 A priori constraints and resolution analysis of back projection, multi-point source and finite fault inversion methods
上文分析了算法從上游到下游在先驗約束、結果穩定性上的差別,因此通過先驗約束少的上游算法逐步增加對下游算法的約束,在實現算法靈活性的同時,增加其穩定性的可行性. 本文提出一種通過反投影結果約束多點源先驗參數,再通過多點源反演結果構造多斷層幾何模型,從而增加破裂過程反演自由度的串行方法. 下文就串行方法的若干思路進行簡要描述. 反投影算法可以有兩種呈現方式:(a)每一時刻的集中輻射點源的離散表示(如圖2 反投影結果下圖);(b)聚束能量場的時空表示(如圖2 反投影結果上圖). 而每一種表示方法與現有多點源算法有比較容易的串行方式. 本文提出兩種從反投影到多點源串行算法的思路,并針對其各自的優勢進行分析.
(1)離散表示結果與時間域多點源算法的串聯
反投影結果的離散表示與時間域迭代算法的輸入比較契合. 在離散表示中,反投影的能量呈現為空間離散的若干個點源,這與時間域迭代算法需要提前設定的空間離散的震源位置情形類似. 而每個離散點源的能量釋放時間函數(包括起始時間以及震源時長)可以用反投影聚束結果給出比較好的約束. 因此,通過反投影聚束結果,我們可以比較容易給出各個點源的位置、時間信息. 這一方式極大地減少了時間域多點源反演中對于預設斷層面位置、網格劃分、點源起始時間以及時間函數的需求. 而多點源反演可以在目標點源位置和時間獲得其震源機制信息. 在點源個數方面,也可以根據反投影點源的能量大小,從大到小篩選若干個點源進行反演.
(2)聚束能量場時空分布與MCMC 多點源算法串聯
反投影的能量聚束直接體現為聚束能量的時空分布. 如果直接將聚束能量場的時空函數轉化為具有時空分布的先驗概率分布,并通過MCMC 方法進行多點源反演. 可以通過擬合實際波形獲得電源位置以及發震時間的后驗概率. 該思路的難點體現在兩個方面:(a)如何構建聚束能量到先驗概率的數學表達. (b)如何通過有先驗概率的采樣以及震源機制的采樣提高反演效率. 如果不對MCMC 的參數采樣做技術改進,MCMC 應該會直接搜索全部參數空間,包括每個點源的走向、傾角、滑動角以及地震矩. 而如果震源時間以及位置的采樣已經通過先驗概率給出,則基于該時空采樣的震源機制可以通過線性算法獲得. 這種結合采樣方法和線性方法的聯合策略有望提升MCMC 的反演效率.
類似將反投影結果應用于多點源反演的思路,張喆和許力生(2020) 將其應用于2013 年南斯科舍海嶺M
7.8 的地震中,并且取得了較好的反演結果. 而形成自動化的算法以及具有嚴格概率含義的先驗概率轉化方法,尚需要進行更多的探討. 同時,反投影成像技術應用不同的臺網可能獲得不同的時空分布. 這與三維介質下遠震體波的射線矯正相關. 通過三維射線矯正可以較好地解決該問題(Liu et al., 2018). 如果反投影體現的高頻輻射與多點源的低頻輻射存在差別,也可以通過增加先驗概率函數時空范圍的方法獲得二者兼顧的反演效果. 同時多點源反演的不確定性往往高于單點源反演,如GCMT 的角度誤差約15°. 如果用MCMC 方法實現多點源反演,往往可以對其不確定性給出較好的估計. 在后續工作中,結合地質構造信息,有可能進一步縮小震源機制的不確定性.實現反投影結果與多點源數據的串行方法,不僅為我們穩定地解析大地震破裂過程提供了策略,也為震后快速響應提供了可行性方案. 基于遠震體波的反投影,多點源和破裂過程反演對于數據的依賴較弱. 通??梢栽谡鸷髱资昼妰全@得. 三種方法的計算效率都較高,在參數穩定的情況下,完成一次串行反演通常可以在1 小時之內完成. 這為地震破裂過程速報以及救援快速部署提供了重要依據.在震源破裂過程快速反演之后,通過數值模擬的方法進一步計算震源區的強地面運動,也是一套可能為震后救災提供重要依據的模型驅動思路. 數值模擬方法一般基于有限斷層模型,對于破裂過程以及震源機制準確度要求較高. 因此通過多點源方法獲得更加準確的斷層面,以及在多點源機制約束下的有限斷層模型也可以進一步增加用數值模擬計算震源區強地面運動的準確性.