馬付紅，李瑩甄

(1.四川省地震局西昌地震中心站，四川 西昌 615000；2.防災科技學院，河北 燕郊 065201)

2017年8月8日21時19分，在四川省阿壩藏族羌族自治州九寨溝縣發(fā)生M7.0地震，地震震中位于(33.20°N,103.82°E),震源深度20 km(據中國地震臺網)。此次地震造成25人死亡，525人受傷，6人失聯， 176 492人受災，73 671間房屋受損；其最大烈度達Ⅸ度，Ⅵ度及以上總面積為18 295 km2(李志強等，2017)。大地震發(fā)生后人們除了關注其發(fā)震位置和破壞程度外，還關注其發(fā)震的物理過程，錯動方式，發(fā)震斷層面和破裂方向等，而這些物理量的確定可以用震源機制解進行描述。九寨溝M7.0地震發(fā)生后，多家研究機構和個人對其用不同的資料和方法得出了不同的震源機制，不同震源機制之間存在一定的差別，為此有必要進行多個震源機制中心解的確定，以便更準確的確定出此次地震的震源機制最優(yōu)解。最優(yōu)的震源機制解，將對后續(xù)的應力分析，應力場研究，地震應力前兆分析，地震靜態(tài)應力觸發(fā)(萬永革，2001)，大地震對周圍斷層影響(萬永革等，2009)等研究提供幫助。多個震源機制中心解是指對于幾個震源機制測量結果，確定一個震源機制，使所有震源機制與其旋轉角平方和最小(萬永革，2019)。目前計算同一地震多個震源機制中心解的方法除了有Levenberg─Marquadt算法(Levenberg，1944；Marquardt，1963；萬永革，2019)，還可以用旋切法、牛頓-拉普森法(Press et al,1992；萬永革，2019)、二分法等。求震源機制中心解首先要找到表達兩個震源機制差別的方法，本文用兩個震源機制空間旋轉角差別的方法來表示。表達兩個震源機制解空間旋轉角的差別是非線性問題，引入Levenberg─Marquadt算法將非線性問題轉化成線性問題進行迭代求解(Levenberg，1944；Marquardt，1963；萬永革，2019)。

1 發(fā)震背景

九寨溝M7.0地震發(fā)生在塔藏斷裂、虎牙斷裂、岷江斷裂和雪山梁子斷裂所圍限的區(qū)域(趙博等，2018)。塔藏斷裂為全新世左旋走滑斷裂，虎牙斷裂為全新世逆斷裂，岷江斷裂為西傾的全新世逆斷裂，岷江斷裂和塔藏斷裂分別為巴顏喀拉地塊的東、北邊界。本次地震發(fā)震斷層以左旋走滑運動為主，是巴顏喀拉塊體邊界斷裂持續(xù)活動的結果，主要表現為印支板塊向北擠壓，導致川西北地塊擠出(薄景山等，2018)。結合該區(qū)域地質構造背景和震源機制解，其發(fā)震斷層為塔藏斷裂南段，并與虎牙斷裂北端隱性趨于貫通(姚鑫等，2017；易桂喜等，2017)。

2 方法原理

獲取同一地震多個震源機制的中心解，要求所有測量震源機制與中心解的空間旋轉角平方和最小。在建立直角坐標系后最小空間旋轉角可以將兩個震源機制進行比較，把他們當中的一個震源機制的B、P、T軸通過坐標旋轉到另一個震源機制的B、P、T軸上，在空間上旋轉的角度值大小代表了兩個震源機制之間的差異(萬永革，2019)。所有測量震源機制與所求中心解的最小空間旋轉角平方和最小的表達式為：

(1)

在(1)式中(θci)min表示最小空間旋轉角，采用Levenberg-Marquardt方法進行迭代求解(Levenberg，1944；Marquardt，1963；萬永革，2019)。進行迭代求解時選定的初始解接近于估計的解，這樣估計解就可以表達為初始解附近的泰勒展開，忽略二階及以上導數項可以得：

(2)

在(2)式中，vi是第i個震源機制與所求中心解之間最小空間旋轉角的平方和。當震源機制總數為N時，將上式寫成矩陣表達式得：

(3)

將(3)式左邊的第一個雅克比矩陣設為J，第二個矩陣(在初始解基礎上改變的參數量)設為x，同時把右邊矩陣設為d，則變量x表示為：

x=(JTJ)-1JTd

(4)

對(4)式進行求解，則可得到改變量x(即Δφ，Δδ，Δλ)。令

φ0=φ0+Δφ，δ0=δ0+Δδ,λ0=λ0+Δλ

(5)

把(5)式代入(3)式，反復計算，直到解的改變量很小，然后終止迭代，這樣就可以得到非線性問題求解后的最優(yōu)解。為了讓求解快速收斂的同時，防止雅克比矩陣產生奇異，采用Levenberg-Marquardt算法時加上一個可變的阻尼，把(4)式改為：

x=(JTJ+κI)-1JTd

(6)

式中κ為計算中可變的數，I為單位矩陣。把測量的各個震源機制與中心震源機制的最小空間旋轉角的標準差作為計算的最終誤差范圍，其標準差表達式為：

(7)

在(7)式中s為所有測量震源機制與中心解的最小空間旋轉角的平方和。

3 數據與處理

本文收集了九寨溝M7.0地震后美國地質勘探局(USGS)(https://earthquake.usgs.gov/earthquakes/eventpage/us2000a5x1/moment-tensor)、全球質心矩張量項目(gCMT)(https://www.globalcmt.org/)、中國地球物理研究所(IG-CEA)、中國地震臺網中心(CENC)(http://www.cenc.ac.cn/eportal/fileDir/cenc/resource/cms/2018/06/2018060811421545641.pdf)、中國地震預測研究所(IEF-CEA)、易桂喜等(2017)、謝祖軍等(2018)、楊宜海等(2017)、劉旭宙(2017)、皺立曄(2018)、梁姍姍等(2018)、王瑩等(2017)以及韓立波等(http://www.cea-igp.ac.cn/Uploads/image/20170811/2017年8月8日四川九寨溝7-v1.0.jpg)機構和研究人員的震源機制解(見表1)。根據求中心解的原理，選擇不同的初始解，將產生不同的震源機制和標準誤差。本研究均以各機構和個人的震源機制作為初始解，把每次計算所得的目標函數進行比較，選標準誤差最小的作為中心解。以研究機構和個人產出的震源機制解作為初始震源機制解進行計算得出的中心震源機制結果及標準差見表2。從表2可以看出不管以哪個機構或個人得到的震源機制為初始震源機制，得到的中心震源機制解的差別都很小，標準差到小數點后4位都是相同的。根據標準誤差的描述以韓立波等的震源機制作為初始解，計算得到的最小三維空間旋轉角的標準差最小(表2第4列14行)，其值為23.12492。通過求中心解的程序得出本研究的中心解節(jié)面Ⅰ的走向152.52°、傾角84.40°、滑動角-3.95°，節(jié)面Ⅱ的走向242.91°、傾角86.07°、滑動角-174.39°。壓應力軸P走向107.78°，不確定范圍為85.53°～130.53°，傾伏角6.75°，不確定范圍為-9.32°～22.31°。張應力軸T走向17.64°，不確定范圍為-4.61°～40.39°，傾伏角1.17°，不確定范圍為-14.78°～17.20°。中間應力軸B走向為277.83，不確定范圍為168.51°～409.34°，傾伏角83.15°，不確定范圍為61.90°～86.05°。該地震的中心解(a)和空間三維輻射花樣(b)見圖1。

表1 九寨溝7.0級地震多個震源機制解

表2 九寨溝7.0級地震的中心解、標準誤差和最小旋轉角結果

圖1 九寨溝地震的中心震源機制解(a)及空間三維輻射花樣(b)(a)中的粗黑色弧線表示中心震源機制的兩個節(jié)面，淡灰色細弧線覆蓋區(qū)域為其不確定范圍；黑色的圓點表示中心震源機制解的P軸、T軸和B軸，其周圍對應顏色的封閉曲線表示其不確定性的范圍，雪花、小圓圈和加號的點表示各個機構得到的震源機制解P軸、T軸和B軸，黑色細弧線表示各個機構和作者得到的震源機制節(jié)面。(b)中的壓縮區(qū)域和膨脹區(qū)域分別用黑色和灰色表示。

4 結論與討論

(1)本文收集了多家研究機構和個人用不同資料和方法得出的多個震源機制解，分別以他們的震源機制解作為初始解進行計算比較。以標準誤差最小的韓立波等的震源機制為初始解，得出九寨溝M7.0地震的中心解節(jié)面Ⅰ的走向152.52°、傾角84.40°、滑動角-3.95°，節(jié)面Ⅱ的走向242.91°、傾角86.07 °、滑動角-174.39°。壓應力軸P走向107.78°，不確定范圍為85.53°～130.53°。張應力軸T走向17.64°，不確定范圍為-4.61°～40.39°。中間應力軸B走向為277.83，不確定范圍為168.51°～409.34°。(2)九寨溝M7.0地震中心解的確定，首先要找出兩個震源機制解空間旋轉角差別的表達方法，再引入Levenberg─Marquadt算法使非線性問題轉化成線性問題進行迭代求解。該中心解不是多個震源機制的平均值，而是所有測量震源機制與中心解的最小空間旋轉角的平方和最小為標準解。(3)之前機構和個人得出九寨溝M7.0地震的多個震源機制解之間存在一定的差別，用求中心解的原理，求出該地震震源機制最優(yōu)解，為該區(qū)域進行應力及應力場、地震應力前兆分析、地震靜態(tài)應力觸發(fā)和大地震對周圍斷層影響等研究發(fā)揮作用。(4)后續(xù)可以用本文提及的其他求中心解的方法，求出該地震的中心解，并與本方法求得的中心解進行比較研究。

致謝：感謝防災科技學院萬永革教授提供了求中心解的程序。