2021年青海瑪多MS7.4地震隨機有限斷層三維地震動模擬

李春果 王宏偉 溫瑞智 強生銀 任葉飛

(中國地震局工程力學研究所， 中國地震局地震工程與工程振動重點實驗室， 哈爾濱 150080)

0 引言

據中國地震臺網測定， 2021年5月22日02:04:11青海省果洛藏族自治州瑪多縣發生MS7.4地震， 震中位于(34.59°N， 98.34°E)， 震源深度17km。 此次地震是2008年汶川MS8.0地震后中國發生的震級最大的地震。 據多個機構提供的震源機制解(Global CMT、 USGS和中國地震局地球物理研究所等)， 此次地震為走滑型地震， 矩震級為MW7.4。 據中國地震臺網記錄(1)http: ∥www.ceic.ac.cn/history。可知， 截至2021年6月3日， 共發生MS≥3.0余震54次， 其中MS4～4.9地震15次、MS≥5.0地震1次， 余震的震源深度集中分布在8～10km， 呈NWW向分布， 分布區長約170km(詹艷等， 2021)， 主震大體位于余震分布區域的中間位置， 震級最大的余震為5月22日10:29:34發生的MS5.1地震， 位于此次地震序列的最NWW端， 距離主震震中約85km(圖 1)。

圖 1 青海瑪多MS7.4地震及其MS≥3.0余震(截至2021年6月3日)Fig. 1 Epicenters of the Maduo MS7.4 earthquake and its aftershocks with MS≥3.0(up to June 3， 2021).

根據震后實地震害調查， 并參考震區斷裂構造、 余震分布、 震源機制、 儀器烈度分布等， 結合強震動觀測記錄， 應急管理部正式發布了青海瑪多MS7.4地震宏觀烈度圖， 最高烈度為Ⅹ度(2)http:∥www.mem.gov.cn/xw/yjglbgzdt/202105/t20210528_386251.shtml。， 烈度圈長軸呈NWW走向(圖 1)。 瑪多地區的建筑抗震設防烈度為Ⅷ度， 抗震能力高， 倒塌房屋較少， 但道路、 橋梁等基礎設施破壞嚴重， 例如野馬灘大橋發生 “落梁”破壞， 橋梁梁板塌落、 橋墩受損(3)http:∥www.guoluo.gov.cn/html/106/332857.html。。 此外極震區人口密度小， 人員傷亡輕。

根據地震序列精定位結果， 尹欣欣等(2021)、 王未來等(2021)均認為瑪多MS7.4地震的發震構造為昆侖山口-江錯斷裂， 且在瑪多-甘德斷裂的東段觸發了一系列小地震(圖 1)。 昆侖山口-江錯斷裂位于巴顏喀拉塊體與東昆侖-柴達木塊體之間的東昆侖斷裂帶以南70km附近， 東昆侖斷裂帶庫賽湖段于2001年發生過昆侖山口西MS8.1地震。 印度板塊向N長期推擠亞歐大陸板塊， 青藏高原物質向E “逃逸”， 導致巴顏喀拉塊體與川滇塊體向E呈現不均勻擠出的差異運動， 使巴顏喀拉塊體東部NW和NWW向邊界斷裂表現出左旋拉張的特點， 并引發此地區大部分的左旋走滑型地震(陳立春等， 2010； 孫鑫喆等， 2010； 陳長云等， 2013)。 巴顏喀拉塊體內不同規模的斷裂帶具有孕育和發生強震的構造條件(梁明劍等， 2020)， 已發生多次破壞性大地震。

強震動觀測記錄是工程地震學研究的重要基礎數據， 能夠直觀地反映地面運動的強弱， 有助于地震宏觀烈度的評定、 震后應急救援及損失評估。 據中國地震局工程力學研究所強震動觀測組提供的強震動數據可知， 此次地震中共獲得自由場三分量強震動觀測記錄16組， 均為震中東北部或北部的遠場記錄， 震中距為175～450km， 距離震中最近的063DAW臺站東西(EW)、 南北(NS)、 垂直(UD)向峰值地面加速度(PGA)分別為46.0cm/s2、 40.6cm/s2、 19.1cm/s2， 由于震中附近缺少強震動觀測臺站， 在此次地震中并未獲取有效的近場強震動記錄。

本文采用隨機有限斷層三維地震動模擬方法給出了青海瑪多MS7.4地震的地震動場， 對比少量觀測記錄峰值及地震動預測方程的預測中位值， 給出了部分虛擬觀測點的三分量地震動加速度時程、 速度時程及5%阻尼比擬加速度反應譜(PSA)， 基于模擬記錄給出了地震烈度分布， 評估了地震影響場。

1 隨機有限斷層三維地震動模擬方法

(1)

(2)

(3)

其中，φij表示第j個觀測點相對第i子斷層的方位角。

2 地震動模擬

隨機有限斷層法模擬地震動需首先給定震源破裂模型、 傳播路徑衰減及場地效應模型等關鍵模擬輸入參數。

2.1 震源破裂模型

震后多個機構給出了瑪多MS7.4地震的震源破裂過程反演結果。 USGS設定地震破裂面的走向和傾角分別為N106°E和76°， 震中位于(34.586°N， 98.255°E)， 震源深度10km， 破裂長度和寬度分別為182km和31.5km(據Somerville等(1999)裁剪反演破裂面)， 子斷層的尺寸為3.5km×3.5km， 基于遠震寬頻帶P波、 SH波和長周期面波數據的反演結果表明， 此次地震為均勻雙側破裂， 破裂至地表面， 最大滑動量可達3m(4)https:∥earthquake.usgs.gov/earthquakes/eventpage/us7000e54r/finite-fault。(圖2a)。 北京大學張勇教授設定地震破裂面的走向和傾角分別為N92°E和67°， 震中位于(34.586°N， 98.255°E)， 震源深度10km， 破裂面的長度和寬度分別為195km和20km(據Somerville等(1999)裁剪反演破裂面)， 子斷層的尺寸為5.0km×5.0km， 基于全球寬頻帶波形數據的反演結果表明， 此次地震為非均勻雙側破裂， 主要滑動量位于震中SEE向的破裂面上， 最大滑動量近4.0m， 地震破裂至地表(5)http:∥geophy.pku.edu.cn/tpxw/353689.htm。(圖2b)。 中國地震局地質研究所基于歐洲航空局(ESA)的Sentinel-1 SAR數據給出了此次地震的InSAR同震形變場和斷層滑動模型， 破裂面沿昆侖山口-江錯斷裂展布， 破裂面的形狀不規則， 傾向NNE， 傾角約為80°， 破裂面的長度和寬度約為160km和20km(據Somerville等(1999)裁剪反演破裂面)， 子斷層的尺寸約為5.0km×5.0km， 此次地震為非均勻雙側破裂， 滑動量主要分布于震中SEE方向的破裂面上， 最大滑動量可達5m(圖2c)， 以中國地震臺網測定的地震震源位置為起始破裂點， 破裂速度設定為0.8β(β為剪切波速)。

圖 2 青海瑪多MS7.4地震的震源破裂模型Fig. 2 The source rupture models of the Maduo MS7.4 earthquake in Qinghai Province.a USGS模型； b 北京大學張勇教授模型； c 中國地震局地質研究所模型； d 第1、 2個破裂隨機模型

同時， 本文基于Graves等(2010， 2015， 2016)提出的震源破裂隨機模型構建方法建立了破裂面上滑動量不均勻、 破裂速度與滑動量相關的震源運動學破裂隨機模型。 青海瑪多MS7.4地震的矩震級為MW7.4， 參考USGS的震源機制解， 將地震破裂面的走向、 傾角和滑動角分別設定為92°、 67°和-40°， 根據地震破裂面幾何(沿走向的長度L和沿傾向的寬度W)與MW、 破裂類型的經驗關系(Wellsetal.， 1994)及地殼孕震層厚度對破裂面寬度的限制(Pachecoetal.， 1992； Shawetal.， 2008)， 本研究將青海瑪多MS7.4地震震源破裂面的長度和寬度分別設定為L=104km、W=16.8km， 根據Somerville等(1999)建立的平均滑動量(D)與矩震級的經驗關系， 將此次地震破裂面上的平均滑動量設定為D=180.63cm， 將破裂面沿走向和傾向劃分為416個2.0km×2.1km的子斷層。 根據余震空間分布(圖 1)可知， 此次地震近似為雙側破裂事件， 考慮到地震極有可能破裂至地表面且震源深度略大于破裂面寬度， 設定此次地震的起始破裂點在破裂面上沿走向和傾向的比例分別為0.5和1.0， 即起始破裂點在破裂面上沿走向和傾向的位置分別為52km和16.8km， 破裂面上邊緣埋深1.54km。 本文生成了10個青海瑪多MS7.4地震的震源破裂隨機模型， 其中第1、 2個破裂隨機過程如圖2d 所示， 最大滑動量可達5.0m。

地震應力降是影響高頻地震動模擬結果的關鍵性因素， 然而地震應力降的不確定性很大， 可在0.1～100MPa范圍內變化。 大量觀測與理論計算結果表明， 大地震(破裂面寬度達到最大限值而破裂面長度繼續擴展的地震)的應力降與震源破裂面的尺度有關， 當破裂面長度不超過一定值時， 破裂面的平均位錯量與破裂面長度相關， 否則與破裂面寬度相關(Bodinetal.， 1996； Shawetal.， 2001； Manighettietal.， 2007)。 根據地震平均應力降Δσ與破裂面縱橫比L/W的關系(Dalgueretal.， 2008)， 本文分別給出了USGS、 北京大學張勇教授、 中國地震局地質研究所反演的震源破裂模型及震源破裂隨機模型統計意義上的平均應力降Δσ， 分別為3.413MPa、 5.189MPa、 4.406MPa和3.597MPa， 經驗估計值與全球范圍內大地震的平均應力降(約4MPa)較為接近(Allmannetal.， 2009； Baltayetal.， 2019)。

2.2 傳播路徑衰減和場地效應模型

根據Crust1.0全球地殼模型給出了研究區域的地殼介質水平成層模型， 其中震源處介質的P波和S波波速分別為6.10km/s和3.52km/s， 介質密度為2.72g/cm3。 幾何擴散模型采用震源距的三段線性模型(Atkinsonetal.， 1992)， 3個距離段的幾何擴散分別為R-1、R0、R-0.5， 采用與地殼厚度相關的距離分段點， 據Crust1.0模型可知地殼厚度為65km， 則距離分段點分別為97.5km和162.5km。 近似采用四川西部地區的S波品質因子QS=191.8f0.56(趙翠萍等， 2011)， 根據已有研究給出的P波品質因子(QP)與QS的關系(王勤彩等， 2005)， 設定QS/QP=1.5， 即近似有QP=128f0.56。 地震動路徑持時采用活動地殼區(ACR)的地震動路徑持時經驗模型表示(Booreetal.， 2015)。

地殼介質經驗放大模型采用美國國家地震減災計劃(NEHRP)B、 C類場地分界對應的場地條件(即場地地表以下30m土層的等效剪切波速VS30=760m/s)的放大效應經驗模型(Atkinsonetal.， 2006)。 局部場地放大效應采用VS30、 工程基巖面PGA相關的線性和非線性場地放大效應經驗模型(Seyhanetal.， 2014)， 該模型以B、 C類場地分界作為工程基巖面， 其中場地VS30根據Heath等(2020)基于地形數據的全球VS30估計值確定(6)https:∥earthquake.usgs.gov/data/vs30/。。 與場地相關的地震動高頻衰減參數κ近似取為0.02s。 模擬參數詳見表1。

表1 地震動模擬參數Table1 Input parameters for the stochastic finite-fault simulation of the Maduo earthquake

3 結果及討論

本文采用隨機有限斷層三維地震動模擬方法， 基于不同的震源破裂模型分別給出了震中周邊區域(32.09°～37.09°N， 95.84°～100.84°E)內4461個虛擬觀測點的地震動加速度時程， 虛擬觀測點的VS30在180～850m/s范圍內近似服從對數正態分布， 其平均值和自然對數標準差分別為491.6m/s和0.30m/s， 多數虛擬觀測點的VS30在360～660m/s范圍內， 主要為Ⅰ、 Ⅱ類場地(呂紅山等， 2007)， 大多數近斷層虛擬觀測點的VS30≤500m/s， 局部場地可能對近斷層的水平向地震動有一定影響。

基于不同震源破裂模型的模擬記錄峰值地面加速度(PGA)和峰值地面速度(PGV)隨距離的衰減如圖 3 所示， 其中基于破裂隨機模型的模擬記錄峰值為10個破裂模型模擬記錄峰值的均值， 圖中給出了對數空間下均勻分布的19個距離分段的PGA、PGV對數平均值及±1倍標準差的范圍。 震源破裂模型對近斷層的模擬記錄峰值有一定影響， 除此之外， 基于不同震源破裂模型的模擬記錄峰值整體上十分接近。 俞言祥等(2013)為中國第5代區劃圖編制建立了青藏地震區基巖場地(VS30≥500m/s)水平向地震動峰值預測方程， 圖中對比了該預測方程給出的所有虛擬觀測點的預測中位值， 預測方程的長軸方向與破裂面的走向一致(N92°E)。 經對比發現， 模擬記錄的PGA、PGV與預測中位值的大小、 隨距離的衰減特性均具有很好的一致性。 本文同時對比了10組斷層距≤400km的遠場強震動觀測記錄的PGA和PGV， 對觀測記錄依次進行了零線校正、 記錄波形首位加余弦窗并補零、 巴特沃斯非因果帶通濾波(0.08～30Hz)處理，PGA、PGV的觀測值普遍高于模擬值及預測中位值。 由于模擬值與預測值具有一致的距離衰減， 且模擬中已考慮了局部場地效應對水平向地震動的放大效應， 我們推斷普遍偏低的模擬值可能與模擬中采用的地震應力降偏小有關。

圖 3 基于不同震源破裂模型的模擬記錄PGA和PGV隨距離的衰減Fig. 3 Attenuation of PGA and PGV for simulated ground motions produced by various source rupture models.

基于不同震源破裂模型， 圖 4 給出了PGA最大的虛擬觀測點及最接近震中(#3740)、 瑪多縣城(#3683)、 果洛藏族自治州市區(#1966)的各4個虛擬觀測點的模擬記錄加速度時程、 速度時程和PSA。 近斷層模擬記錄的水平向PGA主要在300～500cm/s2范圍內， 水平向PGV主要在30～50cm/s范圍內。 模擬記錄垂直向的加速度、 速度時程的幅值遠低于水平向， 其PSA幅值也低于水平向， 其中近斷層模擬地震動記錄尤為明顯， 局部場地效應對水平向地震動的放大效應會導致水平向地震動高于垂直向。 其次， 采用的較大的P波震源譜衰減率(2.4)可能會降低高頻P波的幅值， 尤其是近斷層區域。 近斷層虛擬觀測點模擬記錄中在S波前抵達的P波比較突出， P波的幅值甚至高于S波。 我們還注意到模擬記錄的速度時程具有較為顯著的速度脈沖， 這說明較大的地震破裂面產生了更為豐富的長周期地震動。 由于震源破裂尺度較大， 破裂持續時間較長， 近斷層地震動的持續時間可達40s。

圖 4 模擬記錄加速度時程、 速度時程及PSAFig. 4 Acceleration and velocity time histories， and PSA.a USGS模型； b 北京大學張勇教授模型； c 中國地震局地質研究所模型； d 第1個破裂隨機模型； e 第2個破裂隨機模型

地震烈度是衡量地震對地表及工程結構破壞程度的重要指標。 根據《中國地震烈度表 GB/T17742-2020》， 利用虛擬觀測點的三分量模擬記錄給出了合成加速度時程和合成速度時程， 得到了三分量合成PGA和PGV， 最終計算了虛擬觀測點基于模擬記錄的地震烈度及其空間分布， 如圖 5 所示。 基于模擬記錄的地震烈度圖長軸方向與采用的震源破裂模型的走向一致， 極震區沿斷層破裂面分布， 破裂面上的滑動分布則控制著最高烈度區的具體位置， 基于中國地震局地質研究所的破裂模型模擬記錄的近斷層區烈度分布與地震宏觀烈度圖的極震區分布更接近。 受局部場地效應的影響， 烈度區呈不規則的橢圓形， 或出現零散的高烈度異常區。 基于模擬記錄的地震烈度區的范圍普遍小于應急管理部正式發布的結果(圖 1)， 其中破裂斷層附近的最高烈度為Ⅸ度或Ⅷ度， 低于正式發布的最高烈度Ⅹ度， 瑪多縣城、 果洛藏族自治州市區的烈度分別為Ⅵ～Ⅶ度、 Ⅲ～Ⅳ度， 均低于正式發布的Ⅶ度、 Ⅵ度， 基于模擬記錄的地震烈度對低烈度區的低估尤為明顯， 整體偏低1～2度。

圖 5 基于不同地震破裂模型模擬記錄給出的青海瑪多MS7.4地震的地震烈度分布Fig. 5 Seismic intensity estimated by the simulated ground motions based on various source rupture models.a USGS模型； b 北京大學張勇教授模型； c 中國地震局地質研究所模型； d 破裂隨機模型

模擬采用的地震應力降是基于破裂面縱橫比的統計意義上的平均值， 地震應力降近似服從對數正態分布， 其對數標準差可高達0.6(log10)， 考慮到地震應力降具有很大的不確定性， 模擬中采用的地震應力降偏低可能導致基于模擬記錄的地震烈度普遍偏低。 根據地震應力降不確定性與地震動不確定性的關系， 即τlog10(PSA)=0.702exp(-0.062T)τlog10(Δσ)(Wangetal.， 2019)， 假設應力降的對數標準差τlog10(Δσ)為0.6， 以短周期的PSA近似表示PGA， 2.0s周期的PSA近似表示PGV， 則PGA的對數標準差(τlog10(PGA))和PGV的對數標準差(τlog10(PGV))分別約為0.42和0.37， 因此當地震應力降的變化范圍為平均值±1倍標準差時， 相應的地震烈度應±1.12(烈度≥Ⅵ)或±1.23(烈度<Ⅵ)。 進一步基于USGS震源破裂模型模擬地震動場， 其中地震應力降為平均應力降+1倍標準差(13.59MPa)， 應力降提高1倍標準差后， 基于模擬地震動的地震烈度分布與宏觀地震烈度圖的Ⅶ及以上烈度區具有較好的一致性(圖 6)。

圖 6 應力降提高1倍標準差后， 基于USGS震源破裂模型的地震烈度估計值Fig. 6 Seismic intensity estimated by the simulated ground motions based on the USGS model with the stress drop increased to mean plus one standard deviation.

4 結論

本文利用青海瑪多MS7.4地震的多個震源破裂模型或斷層滑動模型以及震源破裂隨機模型， 基于隨機有限斷層的三維地震動模擬方法， 考慮VS30相關的局部場地放大效應， 給出了震中附近一定區域內4461個虛擬觀測點的地震動三分量加速度時程， 對比了少量遠場強震動臺站的觀測記錄峰值及青藏地震區的地震動峰值預測中位值， 提供了部分近斷層臺站的加速度、 速度時程及PSA， 最后基于模擬記錄給出了地震烈度空間分布， 結果表明：

(1)模擬記錄的峰值(PGA、PGV)與俞言祥等(2013)地震動預測方程的預測中位值具有很好的一致性， 但模擬記錄的峰值普遍低于遠場觀測記錄的峰值。

(2)基于隨機有限斷層三維地震動方法的模擬記錄可用來評估地震烈度及其影響場， 地震烈度的長軸方向與破裂面走向一致， 極震區沿破裂面分布， 破裂面上的滑動分布決定了最高烈度區的位置， 受局部場地效應的影響， 烈度區呈不規則的橢圓形。 與應急管理部正式發布的烈度圖相比， 基于模擬記錄的地震烈度整體上偏低1～2度。

(3)地震應力降是決定地震動強度、 地震烈度評估結果的重要參數， 模擬中采用的地震平均應力降可能低估了青海瑪多MS7.4地震的應力降， 進而導致模擬結果普遍偏小。 由于地震應力降的不確定性， 地震烈度的評估結果的不確定性約在±1度以上。

致謝中國地震局地質研究所InSAR工作組提供了青海瑪多MS7.4地震的斷層滑動模型； 北京大學張勇教授提供了青海瑪多MS7.4地震的震源破裂過程反演結果； 中國地震局工程力學研究所為本研究提供了強震動觀測數據。 在此一并表示感謝!