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

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

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

0 引言

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

圖 1 青?，敹郙S7.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).

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

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

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

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

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

(1)

(2)

(3)

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

2 地震動模擬

隨機有限斷層法模擬地震動需首先給定震源破裂模型、 傳播路徑衰減及場地效應(yīng)模型等關(guān)鍵模擬輸入?yún)?shù)。

2.1 震源破裂模型

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

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

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

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

2.2 傳播路徑衰減和場地效應(yīng)模型

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

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

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

3 結(jié)果及討論

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

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

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

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

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

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

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

圖 6 應(yīng)力降提高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 結(jié)論

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

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

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

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

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