尹曉菲 王 芃, 張 偉 邵志剛曹文忠 熊仁偉
1) 中國(guó)北京 100036 中國(guó)地震局地震預(yù)測(cè)研究所
2) 中國(guó)廣東深圳 518055 南方科技大學(xué)地球與空間科學(xué)系
2022年1月8日01時(shí)45分27秒青海省海北藏族自治州門(mén)源回族自治縣發(fā)生MS6.9地震,震中位置為(37.77°N,101.26°E),震源深度為10 km.本次地震造成了17 069人受災(zāi),由于震中距離人口稠密地區(qū)較遠(yuǎn),尚無(wú)人員傷亡的相關(guān)報(bào)道(青海日?qǐng)?bào),2022),蘭新高鐵浩門(mén)至山丹軍馬場(chǎng)區(qū)間隧道群由于此次地震而發(fā)生局部塌方(央視網(wǎng),2022).門(mén)源MS6.9地震發(fā)生在青藏高原東北緣冷龍嶺斷裂、托萊山斷裂和肅南—祁連斷裂的階區(qū)部位.現(xiàn)場(chǎng)工作隊(duì)在冷龍嶺斷裂帶西段探測(cè)到長(zhǎng)約22 km的地表破裂帶(青海省地震局,2022).托萊山斷裂和冷龍嶺斷裂現(xiàn)為左旋兼擠壓斷裂(李強(qiáng)等,2013),歷史上無(wú)M7以上地震記載(姜文亮,2018),最近一次強(qiáng)震為2016年1月21日門(mén)源MS6.4地震,該地震的發(fā)生反映了青藏高原地塊向NE向不斷推擠生長(zhǎng)的過(guò)程(胡朝忠等,2016).門(mén)源MS6.9地震發(fā)生后,許英才等(2022)對(duì)門(mén)源地震早期序列(2022年1月8日至12日)進(jìn)行了重定位和震源機(jī)制研究,其結(jié)果表明目前門(mén)源地區(qū)還存在一定的應(yīng)力積累且應(yīng)力尚未得到充分釋放,該地區(qū)仍有發(fā)生強(qiáng)震的危險(xiǎn).由于地震災(zāi)害主要是地震在地表產(chǎn)生的強(qiáng)地面運(yùn)動(dòng)造成的,因此為了獲得地震波場(chǎng)傳播過(guò)程及其引起的地表響應(yīng),對(duì)門(mén)源MS6.9地震開(kāi)展強(qiáng)地面運(yùn)動(dòng)初步模擬及烈度估計(jì)具有重要意義.
地震波場(chǎng)的正演模擬需要考慮復(fù)雜地表.曲線(xiàn)網(wǎng)格有限差分方法(Zhang,Chen,2006;Zhanget al,2012)適用于含起伏地形的地震波場(chǎng)模擬,該方法中地表的形狀用任意曲線(xiàn)網(wǎng)格近似,并采用牽引力鏡像法處理自由表面條件.曲線(xiàn)網(wǎng)格有限差分的強(qiáng)地面運(yùn)動(dòng)模擬方法在地震后的災(zāi)害評(píng)估中得到了非常廣泛的應(yīng)用.Zhang等(2008)結(jié)合三維介質(zhì)模型、震源破裂模型和地表地形數(shù)據(jù)模擬了2008年汶川MS8.0地震的強(qiáng)地面運(yùn)動(dòng),研究表明斷層破裂方式和盆地構(gòu)造控制了地表峰值速度(peak ground velocity,縮寫(xiě)為PGV)的分布,地表起伏劇烈地區(qū)往往對(duì)應(yīng)較大的PGV數(shù)值,應(yīng)關(guān)注其震害問(wèn)題.張振國(guó)等(2014a,b)對(duì)2014年2月12日新疆于田MS7.3地震和2014年8月3日云南魯?shù)镸S6.5地震引起的強(qiáng)地面運(yùn)動(dòng)作了初步模擬和烈度預(yù)測(cè),研究顯示地震動(dòng)在山峰、山脊處具有較大幅值,該結(jié)果可指導(dǎo)震后的災(zāi)區(qū)重建工作.趙宏陽(yáng)和陳曉非(2017)利用1975年2月4日遼寧海城MS7.3地震的地震地質(zhì)資料,模擬計(jì)算了海城地震的波場(chǎng)傳播過(guò)程,分析了強(qiáng)地面運(yùn)動(dòng)的方向性效應(yīng)、盆地效應(yīng)和近斷層效應(yīng),得出理論烈度分布同震后調(diào)查烈度分布基本一致,驗(yàn)證了研究所用的震源模型和速度結(jié)構(gòu)的合理性.
本文擬根據(jù)張勇①?gòu)堄?022年1月10日與作者的個(gè)人交流提供的門(mén)源地震震源破裂過(guò)程的初步結(jié)果,利用曲線(xiàn)網(wǎng)格有限差分方法模擬門(mén)源MS6.9地震的強(qiáng)地面運(yùn)動(dòng),再結(jié)合地表峰值速度和烈度間的關(guān)系,計(jì)算地震烈度,并在此基礎(chǔ)上評(píng)估地震災(zāi)害分布特征,以期為門(mén)源地區(qū)的防震減災(zāi)提供科學(xué)依據(jù).
研究區(qū)范圍如圖1所示,可以看出,該區(qū)域的地貌特征變化大,地勢(shì)西南高東北低,高程介于1.3—5.0 km之間.門(mén)源MS6.9地震發(fā)生后至2022年1月20日,該地區(qū)已發(fā)生MS≥4.0余震23次,其中MS≥5.0余震2次.
圖1 研究區(qū)及周邊的地質(zhì)概況和2022年1月8—20日震源區(qū)MS≥3.0 地震分布Fig. 1 Geological overview of the study area and its surrounding regions,spatial distribution of MS≥3.0 earthquakes in source zone during the period of 8 to 20 January,2022
采用曲線(xiàn)網(wǎng)格有限差分方法(Zhang,Chen,2006;Zhanget al,2012)對(duì)門(mén)源MS6.9地震的強(qiáng)地面運(yùn)動(dòng)進(jìn)行模擬.計(jì)算中需要設(shè)置描述地形起伏的網(wǎng)格模型、反映地下物質(zhì)屬性的介質(zhì)模型,以及表示地震破裂過(guò)程的震源模型.
地形選取GTOPO30地形數(shù)據(jù),其水平分辨率大約為1 km.整個(gè)計(jì)算區(qū)域尺度為350 km×220 km,深度為60 km.將研究區(qū)域離散成700×440×120個(gè)網(wǎng)格,垂直方向采用等間距排列,單位網(wǎng)格為邊長(zhǎng)500 m的立方體.
考慮到面波對(duì)縱波速度的靈敏度和起伏地形的影響,Han等(2022)提出一種改進(jìn)的體波和面波數(shù)據(jù)聯(lián)合反演方法獲得中國(guó)大陸地殼和上地幔水平分辨率為0.5°的速度結(jié)構(gòu)(USTClitho2.0模型).該模型對(duì)于青藏高原地塊內(nèi)體波射線(xiàn)覆蓋相對(duì)稀疏的區(qū)域,加入面波頻散數(shù)據(jù)可更好地提高縱波速度和橫波速度的準(zhǔn)確度,因此本文以USTClitho2.0模型為基礎(chǔ)建立研究區(qū)域的網(wǎng)格化橫波速度和縱波速度模型.每個(gè)網(wǎng)格的東向、北向和垂直方向的長(zhǎng)度分別為43.75 ,55 和5 km,同時(shí)依據(jù)Brocher (2005)提出的密度和縱波速度的經(jīng)驗(yàn)關(guān)系確定各網(wǎng)格的密度.模擬采用的橫波速度結(jié)構(gòu)如圖2所示,可以看出,2—12 km深度存在高速層,沿北東方向地殼厚度逐漸減薄.雖然青海湖處于研究區(qū)域內(nèi)(圖1),但是由于其位于模擬計(jì)算區(qū)域的邊緣,且水深較淺,因此在建立介質(zhì)模型時(shí)忽略其影響.
圖2 橫波速度結(jié)構(gòu)Fig. 2 S-wave velocity structure
張勇①?gòu)堄?022年1月10日與作者的個(gè)人交流反演了該地震的震源破裂過(guò)程,結(jié)果顯示斷層滑動(dòng)出露地表,走向?yàn)?03°,傾角為88°,滑動(dòng)角為-6°,共有31×11=341個(gè)子斷層,子斷層的空間分辨率為2 km×2 km,采樣點(diǎn)為80個(gè),每間隔0.25 s給出一個(gè)滑動(dòng)速率,破裂過(guò)程持續(xù)時(shí)間約20 s,最大滑動(dòng)量為1.8 m.基于這一結(jié)果,本文設(shè)置了相應(yīng)的震源模型,如圖3所示.
圖3 本文所用的門(mén)源MS6.9地震震源模型Fig. 3 Source model of Menyuan MS6.9 earthquake used in this study
強(qiáng)地面運(yùn)動(dòng)模擬計(jì)算的時(shí)間步長(zhǎng)為0.01 s,總步數(shù)為1萬(wàn)步,模擬時(shí)長(zhǎng)共100 s,使用640個(gè)計(jì)算核心進(jìn)行計(jì)算.門(mén)源地震x分量的模擬速度場(chǎng)快照如圖4所示,可見(jiàn):地震的主要能量由位于震中附近下方的位錯(cuò)產(chǎn)生,且以水平走向錯(cuò)動(dòng)為主,因此在斷層垂向上產(chǎn)生了能量較強(qiáng)的S波;大約第2.36 s時(shí),地震波到達(dá)地表,在初 始破裂時(shí)刻,速度的最大值集中在斷層破裂的前鋒上;隨著遠(yuǎn)離發(fā)震斷層,地震波場(chǎng)能量逐漸減小.
圖4 門(mén)源MS6.9地震x分量粒子速度波場(chǎng)快照Fig. 4 Wavefield snapshots of the x component particle velocity of Menyuan MS6.9 earthquake
通過(guò)強(qiáng)地面運(yùn)動(dòng)模擬得到各網(wǎng)格點(diǎn)x,y,z三個(gè)方向不同時(shí)刻的速度,然后對(duì)三個(gè)方向的速度分量求矢量和獲得各網(wǎng)格點(diǎn)的運(yùn)動(dòng)速度隨時(shí)間的變化,根據(jù)各網(wǎng)格點(diǎn)的速度最大值確定研究區(qū)域的地表峰值速度(PGV)分布,并結(jié)合國(guó)家市場(chǎng)監(jiān)督管理總局和國(guó)家標(biāo)準(zhǔn)化管理委員會(huì)(2021)給出的PGV與烈度之間的關(guān)系得到對(duì)應(yīng)的地震烈度,如圖5所示,結(jié)果表明:沿平行斷層走向方向的地震動(dòng)衰減明顯小于垂直斷層走向方向;地震的最大烈度為Ⅷ度(PGV=37.63 cm/s),位于震源破裂起始點(diǎn)附近區(qū)域;Ⅷ度區(qū)主要涉及門(mén)源縣、祁連縣和肅南縣的部分區(qū)域;Ⅶ度區(qū)主要涉及大通縣、永昌縣、民樂(lè)縣等部分區(qū)域.中國(guó)地震局(2022)野外調(diào)查的烈度分布顯示等震線(xiàn)長(zhǎng)軸呈WNW走向,同理論模擬顯示的高烈度區(qū)主要沿WNW方向延伸的結(jié)果一致.理論烈度分布與中國(guó)地震局工程力學(xué)研究所強(qiáng)震動(dòng)觀(guān)測(cè)組(2022)給出的儀器烈度分布(圖6)較為接近,即高烈度區(qū)主要沿?cái)鄬幼呦蛘共?,且斷層南?cè)的地震烈度大于北側(cè).受模型分辨率和計(jì)算成本的限制,本文模擬的地震波場(chǎng)的最高有效頻率為0.4 Hz(2.5 s),缺少高頻成分,因此獲得的烈度較強(qiáng)震儀和烈度儀觀(guān)測(cè)的結(jié)果偏低,但相較于使用低分辨率介質(zhì)模型的模擬結(jié)果(徐劍俠等,2015),地震動(dòng)模擬的最高頻率得到了一定的提升.
圖5 門(mén)源MS6.9地震烈度分布圖Fig. 5 The simulated intensity of Menyuan MS6.9 earthquake
圖6 門(mén)源MS6.9地震儀器觀(guān)測(cè)的烈度分布(引自中國(guó)地震局工程力學(xué)研究所強(qiáng)震動(dòng)觀(guān)測(cè)組,2022)Fig. 6 Instrumental seismic intensity distribution of Menyuan MS6.9 earthquake (after Strong Motion Observation Group,Institute of Engineering Mechanics,China Earthquake Administration,2022)
地震災(zāi)害區(qū)域集中在發(fā)震斷層附近,并向WNW方向和ESE方向延伸,這主要由兩個(gè)原因造成:一方面,該地震為高傾角的走滑型地震,地震波能量主要沿走向方向傳播,表現(xiàn)出強(qiáng)地面運(yùn)動(dòng)的方向性效應(yīng);另一方面,起伏地表對(duì)地震波傳播具有重要影響(Zhanget al,2008),發(fā)震斷層的WNW方向和ESE方向分布有托萊山、大通山、達(dá)坂山和冷龍嶺等山脈(圖1),屬于山地地貌,因此該區(qū)域的地震動(dòng)高值可能與山脊地區(qū)地震波多次反射有關(guān).
強(qiáng)地面運(yùn)動(dòng)模擬結(jié)果由地形數(shù)據(jù)、介質(zhì)模型和震源模型共同決定,因此數(shù)據(jù)選擇對(duì)結(jié)果的可靠性有較大影響,為使結(jié)果更為可信,本文參考美國(guó)加州綜合地震破裂預(yù)測(cè)模型(Fieldet al,2014)中提出的“使用可獲得最優(yōu)解”原則對(duì)介質(zhì)模型進(jìn)行選擇.常用的介質(zhì)模型包括CRUST1.0模型(Laskeet al,2012)和CRUST2.0模型(Bassinet al,2000),二者均基于水平層狀介質(zhì)的假設(shè),水平分辨率分別為1.0°和2.0°,而Han等(2022)提出的USTClitho2.0模型的水平分辨率更高,因此研究中選擇USTClitho2.0模型作為介質(zhì)模型.本文基于門(mén)源地震震源區(qū)附近的地形數(shù)據(jù)、介質(zhì)模型和震源模型,使用曲線(xiàn)網(wǎng)格有限差分方法計(jì)算了該地震的近場(chǎng)地震波傳播過(guò)程,得到了理論的速度場(chǎng)快照和地震烈度分布.研究結(jié)果表明:總體上垂直斷層走向方向的地震動(dòng)衰減大于平行斷層走向方向的地震動(dòng)衰減,震中最大烈度為Ⅷ度;模擬得到的理論烈度同野外調(diào)查的地震烈度分布基本一致.受強(qiáng)地面運(yùn)動(dòng)方向性效應(yīng)和起伏地表的影響,地震災(zāi)害主要沿?cái)鄬拥腤NW方向和ESE方向分布.由于地震發(fā)生區(qū)域主要以山地地貌為主,該地震的發(fā)生造成了局部邊坡崩塌、滾石和凍土開(kāi)裂以及蘭新高鐵大橋橋面受損等次生災(zāi)害(頡滿(mǎn)斌,2022).考慮到研究區(qū)域仍有發(fā)生強(qiáng)震的危險(xiǎn)(許英才等,2022),今后的防震減災(zāi)工作中有必要加強(qiáng)地表起伏劇烈區(qū)域特別是山脊地貌的震害防御工作.
北京大學(xué)張勇教授為本文提供了震源破裂過(guò)程的初步結(jié)果,中國(guó)科學(xué)技術(shù)大學(xué)張海江教授為本文提供了三維速度結(jié)構(gòu),中國(guó)地震局工程力學(xué)研究所馬強(qiáng)研究員為本文提供了儀器觀(guān)測(cè)的地震烈度分布圖,中國(guó)地震局地震預(yù)測(cè)研究所徐岳仁研究員與作者就青藏高原東北緣的地質(zhì)概況進(jìn)行了討論,審稿專(zhuān)家為本文提出了寶貴意見(jiàn),作者在此一并表示感謝.