周紅， 李亞南,2， 常瑩

1 中國地震局地球物理研究所， 北京 100081 2 防災科技學院，河北三河 065201 3 礦冶科技集團有限公司礦山工程研究設計所， 北京 100160

0 引言

2021年5月21日21時48分，云南省大理白族自治州漾濞縣發(fā)生6.4級地震.據(jù)中國地震臺網(wǎng)中心測定：震源位置位于(北緯25.67°，東經(jīng)99.87°)，震源深度8 km.

震區(qū)所在的滇西地區(qū)為三江復合造山帶(賀傳松等，2004)的重要組成部分，是青藏高原東緣異?；钴S的陸內(nèi)變形區(qū).該區(qū)地質(zhì)結(jié)構(gòu)復雜，巖漿活動頻繁，構(gòu)造變形與變質(zhì)作用強烈(鐘大賚，1998)，是青藏高原東部地殼塊體旋轉(zhuǎn)擠出的重要通道和構(gòu)造轉(zhuǎn)折帶.震區(qū)主要有NNW向維西—喬后斷裂、NW向紅河斷裂和NE向程海斷裂、龍蟠—喬后斷裂.

漾濞地震發(fā)生在金沙江—元江地震帶，位于維西—喬后斷裂.維西—喬后斷裂在蘭坪—思茅中生界盆地內(nèi)，構(gòu)成蘭坪—思茅坳褶帶的東部邊緣，與鮮水河、安寧河、則木河、小江、紅河及金沙江斷裂帶一起，共同構(gòu)成川滇菱形塊體的邊界.

眾多學者對該地區(qū)的地質(zhì)構(gòu)造和地震活動進行了廣泛的研究.管燁等(2006)利用DSS地震測深、石油反射地震資料對蘭坪—思茅盆地內(nèi)的中軸構(gòu)造帶的構(gòu)造特征進行分析，研究結(jié)果表明中軸構(gòu)造帶在思茅盆地南部存在有深部地幔的上隆現(xiàn)象.王安建等(2013)對該地區(qū)的大地電磁測量、反射地震、折射地震綜合解釋分析，滇西蘭坪—思茅盆地中軸構(gòu)造帶是西南“三江”地區(qū)繼陸-陸碰撞、走滑、隆升之后正在發(fā)育的大陸裂谷系統(tǒng)，在地殼深部和上地幔形成顯著熱窿.李桂華和何家斌(2002)對2002年1月發(fā)生的漾濞地震序列位置分析，震源分布在N73°E走向的一個條槽內(nèi)，深度為11～14 km.常祖峰等(2014)通過對2013年洱源MS5.5級地震震源機制解結(jié)合地震烈度分布、滑坡崩塌體分布和微震分布特征推斷，發(fā)震構(gòu)造為維西—喬后斷裂.宋晉東等(2017)對2017年漾濞MS5.1級地震強震動記錄特征分析，發(fā)現(xiàn)獲取的強震記錄的PGA低于中國分區(qū)地震動參數(shù)衰減關系(國家質(zhì)量監(jiān)督檢驗檢疫總局和中國國家標準化管理委員會，2015；俞言祥等，2013)和美國2008年地震區(qū)劃圖用于西部地區(qū)的NGA-BA08衰減關系(Boore and Atkinson，2008；Wang et al.，2010).尹鳳玲等(2018)基于分層半無限空間黏彈性地球模型計算紅河斷裂帶上庫侖應力變化的演化過程，發(fā)現(xiàn)洱源附近、大理至大斗門以北地區(qū)、元江以南地區(qū)可能仍是地震潛在危險區(qū)段.潘睿等(2019)利用云南和部分四川地震臺網(wǎng)寬頻帶數(shù)字地震波形資料，對2017年漾濞MS5.1和MS4.8地震震源機制解和震源深度測定，研究認為發(fā)震斷層可能是維西—喬后斷裂中南段.受大陸側(cè)向碰撞影響，青藏高原東南緣三江地區(qū)構(gòu)造活動強烈，斷裂發(fā)育，存在中-長期大地震背景(鄧起東等，2014；王椿鏞等，2015；吳鵬等，2020).

重建強地面運動對于確定地震引起的災害范圍、震害強度非常重要，而隨機有限斷層法是重要的研究工具之一.隨機有限斷層法從點源隨機法(Boore，1983，2003)發(fā)展起來：其基本思想是將一個發(fā)震斷層分成多個子斷層單元(Hartzell, 1978)，每個子斷層單元產(chǎn)生的地表地震動利用點源隨機法計算，全部點源地表地震動的疊加形成整個大震產(chǎn)生的地震動.這種方法成功的用于設定地震動的預測中(Satyam and Rao，2009；Zonno et al.，2009；Zonno and Carvalho，2006).隨機方法計算效率高，但是對于地震動低頻成分模擬效果欠佳，為了提高低頻地震動模擬效果，出現(xiàn)了混合模擬地震動方法，即：低頻地震動利用確定性方法模擬，而高頻地震動利用隨機方法模擬(e.g., Pacor et al.,2005；Graves and Pitarka, 2010；Frankel, 2009; Atkinson et al., 2011；Shahjouei and Pezeshk, 2015).但是確定性方法模擬低頻地震動計算效率低，對計算條件又有很高的要求.周紅等(周紅，2018；Zhou and Chang，2019，2020)提出了斷層破裂過程控制的隨機有限斷層法(NNSIM)，該方法不采用確定性方法模擬低頻，而是把低頻信息包括到震源里，即利用遠震記錄反演的震源低頻信息，如斷層應力降、位錯、破裂速度、破裂時間函數(shù)，結(jié)合斷層的幾何參數(shù)、地質(zhì)產(chǎn)狀等形成震源信息，通過頻譜傳遞使得模擬結(jié)果包含了大量低頻信息.該方法在蘆山7.0級地震強地面運動模擬和九寨溝7.0級地震動模擬中均取得了有效應用.本文利用該方法模擬2021年云南漾濞6.4級地震的強地面運動，通過6個強震臺的模擬地震動與觀測記錄的對比，表明了模擬計算參數(shù)選取的有效性.在研究區(qū)域內(nèi)，強地震動模擬獲得的PGV、PSV、PGA和PSA等值線分布近似為以震中為圓心的圓形，進一步討論了云南漾濞6.4級地震的震害強度分布和烈度分布范圍.

1 NNSIM隨機有限斷層法方法介紹

基于遠震記錄反演得到的斷層破裂過程包含了斷層破裂的低頻信息，如破裂位錯、應力降、破裂速度和低頻破裂時間過程.周紅等(Zhou and Chang，2019,2020；周紅，2018)提出了隨機有限斷層法(NNSIM)，引入斷層的低頻破裂過程作為震源參數(shù)，打破了經(jīng)典隨機有限斷層法(Beresnev and Atkinson，1997，1998； Atkinson and Boore，2006；Atkinson et al.，2009)只限高頻段地震動模擬結(jié)果可靠的限制，提高了低頻段模擬結(jié)果的可靠性，得到了近場全頻段地震動模擬結(jié)果.

有限斷層方法的基本思想是：對于大地震，需要將斷層剖分成N個子斷層單元，每一個子斷層單元看作一個小的點震源，考慮到破裂延遲和傳播延遲，將每個點震源產(chǎn)生的地震動疊加，得到整個斷層產(chǎn)生的地震動.

(1)

這里nl和nw分別是沿斷層的長和寬方向上子斷層單元的個數(shù)，Δtij是第ij子斷層單元到觀測點的相對延遲，aij(t)是利用隨機點源方法計算的第ij子斷層單元產(chǎn)生的地震動加速度，a(t)是整個斷層在觀測點模擬的地震動加速度.而單個子震源的地震動加速度基于傅氏譜計算，單個點源的傅氏譜為

×exp(-πκf)SG(f),

(2)

公式(2)中的拐角頻率沿用周紅等(Zhou and Chang，2019,2020；周紅，2018)給出的子斷層非均勻的拐角頻率公式，即第ij子斷層拐角頻率為：

fij=4.9E(+6)β(Δσij/M0ij)1/3.

(3)

這里Δσij和M0ij分別是第ij子斷層的應力降和地震矩.實際的斷層破裂過程中，斷層面上各子斷層單元上的地震矩和應力降都不相同，根據(jù)公式(3)可以得到不同子斷層的拐角頻率.這種拐角頻率與動拐角頻率(Motazedian and Atkinson，2005)形式、大小、取值范圍都不相同.Zhou和Chang(2019)證明了非均勻拐角頻率(3)保留了動拐角頻率(Motazedian and Atkinson，2005)的優(yōu)點，同樣避免斷層不同剖分尺度對于遠場輻射能量的影響.

公式(2)模擬的地震動結(jié)果是否合理，關鍵在于計算參數(shù)的確定.下面詳細討論模擬漾濞6.4級地震的地震動所需要的計算參數(shù).

2 漾濞6.4級地震模擬參數(shù)的確定

2.1 震源參數(shù)

圖1 漾濞斷層面地震矩分布及破裂時間過程(a) 斷層地震矩分布(張旭等，2021)； (b) 地震中斷層面上各子斷層走向滑動速率時間序列； (c) 地震中斷層面上各子斷層傾向滑動速率時間序列.Fig.1 Distribution of seismic moments and rupture time on the fault of the Yangbi Earthquake (Zhang et al., 2021)(a) Distribution of seismic moment on the fault; (b) Slip-rate history at each sub-fault on the strike-profile across the hypocenter; (c) Slip-rate history at each sub-fault on the dip-profile across the hypocenter.

圖2 漾濞斷層應力降分布Fig.2 Distribution of stress drop on the fault of the Yangbi Earthquake

2.2 場地系數(shù)κ

公式(2)中場地項exp(-πκf)是一個高頻濾波器，它表現(xiàn)為每個場點地震加速度的高頻衰減(Anderson and Hough, 1984).我們選定近斷層200 km為研究區(qū)域(圖3白框)，研究區(qū)內(nèi)有六個強震臺(國家強震動臺網(wǎng)中心)，臺站的水平記錄提供了研究區(qū)域內(nèi)κ的變化特征.圖4藍線為六個臺站已完成零漂校正的NS(南北)和EW(東西)加速度記錄振幅譜對數(shù)，紅線為藍線的線性擬合，擬合斜率的變化反映了各個記錄的κ變化.將6個記錄的EW和NS分量κ繪制到圖5中，它們呈現(xiàn)出隨震中距線性變化的特征，即

κ=a1R+a0,

(4)

地震動模擬計算中我們按照該關系式確定各個模擬場點的κ取值.κ的物理含義代表臺站的場地條件，圖5顯示了6點κ的變化及其對應的最小二乘線性擬合直線.其中臺站53BCJ偏離擬合直線較遠，說明53BCJ臺站的場地介質(zhì)較為特殊.

圖3 研究區(qū)域和強震臺站分布白色線表示研究區(qū)域邊界，紅色三角形表示本文中模擬的臺站，白色實心點表示震中位置.Fig.3 The study area and the locations of strong ground motion stationsThe white square denotes the boundary of the study area, the red triangles mean simulated stations, the white solid point shows the epicenter of the Yangbi Earthquake.

圖4 6個臺站記錄振幅譜對數(shù)(a) 加速度EW分量; (b) 加速度NS分量； 藍線為加速度振幅譜對數(shù)，紅線為截止頻率段的線性擬合.Fig.4 The natural logarithm of the acceleration spectra recorded at 6 stations(a) The EW component of the acceleration; (b) The NS component of the acceleration； the blue line is the logarithm of the acceleration amplitude spectrum, and the red line is the linear fitting of the cut-off frequency.

圖5 κ隨震中距的變化關系Fig.5 The relationship between κ and epicenter distance

2.3 衰減指數(shù)λ

(5)

這樣需要反演的參數(shù)共六個：λ1，λ2，λ3，λ4，a0和a1.

2.4 附加時窗參數(shù)

公式(2)SG(f)中的時窗主體函數(shù)為反演的每個子斷層的破裂時間函數(shù)，而地震子斷層最長破裂時長小于4 s，斷層總破裂時長約為8 s.6個臺站記錄顯示距離震中越遠地震記錄持時越長，最遠的53LLP地震動持時超過30 s.把震中距與六個臺站記錄NS和EW加速度持時繪于圖6，顯然持時與震中距表現(xiàn)出線性相關，即：

(6)

為使得模擬記錄的持時能夠靠近實際記錄，針對每個斷層子單元增加了一個附加時窗，如選擇EW分量為模擬目標，附加時窗的時長為

dtij=tEW-tij，

(7)

這里tij為每個子斷層的破裂持時.

圖6 實際地震動加速度持時與震中距的關系Fig.6 The relationship between observed ground motion acceleration duration and epicentral distance

2.5 模擬參數(shù)匯總

經(jīng)過對實際地震記錄的分析、反演和擬合等處理工作，得到了模擬地震動需要的所有計算參數(shù)，匯總?cè)绫?：

表1 模擬地震動需要的計算參數(shù)Table 1 The computation parameters of ground motion simulation

3 模擬結(jié)果與實際結(jié)果的對比

地面運動的水平分量是造成地震破壞的主要成分，因此通常隨機有限斷層法主要用于模擬地表的地震動水平分量.方程(2)C里的V取常數(shù)0.71時，模擬的結(jié)果為地表水平分量(Atkinson et al.， 2011)，但不特定對應EW或者NS分量.為了使模擬參數(shù)帶有本地區(qū)的地質(zhì)特性(包括速度結(jié)構(gòu)、土層、地形等特征)，考慮到與觀測記錄的對比效果，本文選擇以單個水平分量EW分量觀測記錄作為目標函數(shù)，反演模擬參數(shù)，替代隨機有限斷層法經(jīng)常采用的統(tǒng)計參數(shù).當然選擇NS分量作為目標函數(shù)進行模擬和對比同樣可以.

利用以EW分量作為目標函數(shù)反演的計算參數(shù)(表1)模擬六個臺站的加速度和速度記錄的EW分量，圖7和圖8所示.因地震動三要素為峰值、持時、反應譜，其中反應譜最重要，它能比較全面的反應地震動全頻段的信息和對地表的破壞能力，對于模擬的結(jié)果繪制成圖，顯示這三方面的模擬效果，便于展示所使用的計算參數(shù)是否合理.

圖7 模擬加速度與實際臺站加速度對比紅線為觀測記錄(縮寫為Ob)，黑線為模擬加速度(縮寫為Sim).Fig.7 The simulated acceleration and the observed accelerationThe red line is the observed record (left column, abbreviated as Ob), and the black line is the simulated acceleration (right column, abbreviated as Sim).

圖8 模擬的速度與實際臺站速度對比Fig.8 The comparison between the simulated velocity and the observed velocity

由于模擬的時間函數(shù)由隨機白噪聲產(chǎn)生，對于每個場地計算5組時間函數(shù)，峰值PGA和PGV分別寫在圖7和圖8對應場點的模擬加速度子圖里(圖7和圖8右列黑線圖)，與實際記錄(圖7和圖8左列紅線圖)的PGA對比，六個臺站模擬結(jié)果PGA與實際記錄PGA差異不大，符合度較好，PGV除了臺站53YBX相差大一些，其他各臺站符合度較好.模擬結(jié)果的持時偏差相對于PGA大，可能是由于子斷層單元的附加時窗過長產(chǎn)生的，如何確定更合理的時窗，有待于在今后的研究中進一步提高.

圖9和圖10是模擬的加速度和速度的反應譜，即PSA和PSV.圖9和圖10 顯示：六個臺站模擬結(jié)果和實際記錄無論低頻(長周期)還是高頻(短周期)PSA和PSV結(jié)果都較為一致，符合度較好，并且震中距較短的圖(a—c)符合度優(yōu)于震中距較長的圖(d—f).對比圖7—10模擬結(jié)果和實際記錄的一致性說明，2.5節(jié)所給出的計算參數(shù)對于模擬研究區(qū)域的地震動有效，其模擬的地震動能夠代表該區(qū)域的地震動空間分布特征.

圖9 六個臺站模擬加速度反應譜與實際觀測記錄反應譜對比Fig.9 The response spectrum of the simulated acceleration and the observed acceleration at 6 stations

圖10 六個臺站模擬速度反應譜與實際觀測記錄反應譜對比Fig.10 The response spectrum of the simulated velocity and the observed velocity at 6 stations

4 地震動空間分布特征

第3節(jié)的結(jié)果對比顯示：2.5節(jié)給出的計算參數(shù)有效，完全可以用于模擬研究區(qū)域內(nèi)的地震動，以便展現(xiàn)漾濞地震的地震動空間分布特征.在計算地震動空間分布中，每個場點模擬10次隨機地震動時程，多個隨機地震動的使用使得計算的地震動空間分布不僅反應EW分量的特征，更可以反映水平分量的平均效果.以10條時程模擬的平均值作為該場點的PGA、PGV、PSV和PSA值，這樣的模擬結(jié)果能代表漾濞地震水平震動的空間分布，而非單個方向水平分量(如EW或者NS分量)的特征，每個場點模擬的時程越多，越能代表水平分量的平均效果.

在研究區(qū)域(24.6°N—27.8°N，98.5°E—101.5°E)內(nèi)，經(jīng)度和緯度均按照0.02°間隔進行取點，計算合成地震動加速度記錄.每條加速度的絕對值的最大值為對應的PGA值，同樣每條合成地震動的速度記錄絕對值的最大值為PGV值，加速度和速度記錄的反應譜的計算按胡聿賢(2006)介紹的方法計算.

圖11d是研究區(qū)域內(nèi)PGA分布，分布形態(tài)近似圓形，中心點鄰近漾濞縣城，最大值超過400 cm·s-2，整個斷層位于300 cm·s-2等值線內(nèi)，60 km以外地震動出現(xiàn)快速衰減，在很短的距離內(nèi)由100 cm·s-2降為50 cm·s-2表現(xiàn)為等值線變密，震中距80 km以外，PGA減小到20 cm·s-2以下.圖12d是PGV分布，也呈圓形狀態(tài)，震中附近PGV超過20 cm·s-1，震中距80 km以外，衰減到2 cm·s-1，最快變化在震中距60 km附近，漾濞縣城的PGV約為20 cm·s-1.

圖11 模擬的加速度三個周期的反應譜及加速度峰值PGA的空間分布(YB表示漾濞縣城)Fig.11 The spatial distributions of the response spectrum amplitudes at three periods and the PGA of simulated acceleration (YB denotes Yangbi city)

圖12 模擬的速度三個周期的反應譜及峰值速度PGV的空間分布(YB表示漾濞縣城)Fig.12 The spatial distributions of the response spectrum amplitudes at three periods and the PGV of the simulated velocity (YB denotes Yangbi city)

圖11(a—c)對應模擬的加速度0.3 s、1 s和3 s周期的反應譜空間分布，圖11a顯示對于周期0.3 s的反應譜，斷層附近加速度值超過600 cm·s-2，震中距60 km處降到100 cm·s-2，之后以較為緩慢的速度降到50 cm·s-2以下；而1.0 s的反應譜震中的加速度值超過300 cm·s-2，斷層面上為200 cm·s-2，呈圓形分布，3 s代表的長周期PSA在震中不高，30 cm·s-2左右.三個周期反應譜分布盡管幅值相差較大，但總體形態(tài)上具有相似之處，呈現(xiàn)圍繞震中圓形分布，震中60 km附近等值線變密，變化加大.對于圖12中PSV三個周期震中附近的值在40 cm·s-1以下，三個周期0.3 s、1 s和3 s變化幅度不大.

從以上的分析可見漾濞地震斷層破裂震中的長10 km范圍內(nèi)，加速度PGA超過了400 cm·s-2(圖11-PGA),該范圍為本次地震的極震區(qū)，按照《中國地震烈度表》(GB/T17742-2020)(國家市場監(jiān)督管理總局和國家標準化管理委員會，2020)國標規(guī)定該值對應IX度(IX度對應加速度402～830 cm·s-2).而本次地震極震區(qū)達到IX度的原因主要是由于地震破裂能量集中在斷層面9 km×9 km的小范圍內(nèi)，且震中深度淺，大能量集中于小范圍使得地表極震區(qū)達到IX度.

5 結(jié)論

本文利用斷層破裂過程控制的隨機有限斷層法(NNSIM)計算了云南漾濞6.4級地震的強地面運動.隨機有限斷層法最大的優(yōu)勢是計算效率高(秒級)、計算頻帶寬(0～30 Hz 或者更高)；缺點是無法模擬復雜速度結(jié)構(gòu)的特殊作用，即不能詳細討論不同速度結(jié)構(gòu)對地震動的影響.為了在模擬結(jié)果中反應本地區(qū)介質(zhì)作用的特殊性，本文利用實際記錄反演地震動模擬所需要的計算參數(shù)：首先通過分析研究區(qū)內(nèi)強震臺站的觀測記錄的時程和頻譜，確定隨機有限斷層法模擬需要的參數(shù)κ和持時在研究區(qū)內(nèi)的變化特征；以觀測記錄反應譜為目標函數(shù)，采用反演手段獲得模擬地震動所需要的場地系數(shù)κ(代表場地特征)、傳播路徑衰減系數(shù)λ(代表路徑特征)的具體值；結(jié)合漾濞遠震記錄反演的震源破裂過程，模擬合成了研究區(qū)域內(nèi)各場點的加速度、速度和反應譜.由于模擬過程采用了反演獲得的計算參數(shù)，使得模擬結(jié)果一定程度上包含了本地區(qū)介質(zhì)變化的作用.基于模擬的地震動時程，計算研究區(qū)域的PGA、PSA、PGV和PSV，分析計算結(jié)果，可以得到如下的結(jié)論：

(1)研究區(qū)場地因子κ，隨著震中距線性變化；

(2)地震動持時隨著震中距線性變化；

(3)PGA、PSA、PGV和PSV空間等值線分布均呈現(xiàn)為以震中為圓心的圓形分布，沒有走向性，也沒有上下盤效應，這種空間特征說明此時地震表現(xiàn)為小震、點源性質(zhì)；

(4)PGA和PSA都說明漾濞縣城最大峰值加速度超過了400 cm·s-2，對應烈度IX度，因此漾濞縣城會出現(xiàn)不同程度的破壞.

致謝感謝國家自然基金面上項目(41774064)的資助；感謝國家強震動臺網(wǎng)中心提供的漾濞地震的強震資料；感謝國家超級計算廣東中心超算平臺.