云南漾濞6.4級地震強地面運動的模擬和空間分布特征分析

周紅， 李亞南,2， 常瑩

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

0 引言

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

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

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

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

重建強地面運動對于確定地震引起的災害范圍、震害強度非常重要，而隨機有限斷層法是重要的研究工具之一.隨機有限斷層法從點源隨機法(Boore，1983，2003)發展起來：其基本思想是將一個發震斷層分成多個子斷層單元(Hartzell, 1978)，每個子斷層單元產生的地表地震動利用點源隨機法計算，全部點源地表地震動的疊加形成整個大震產生的地震動.這種方法成功的用于設定地震動的預測中(Satyam and Rao，2009；Zonno et al.，2009；Zonno and Carvalho，2006).隨機方法計算效率高，但是對于地震動低頻成分模擬效果欠佳，為了提高低頻地震動模擬效果，出現了混合模擬地震動方法，即：低頻地震動利用確定性方法模擬，而高頻地震動利用隨機方法模擬(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)，該方法不采用確定性方法模擬低頻，而是把低頻信息包括到震源里，即利用遠震記錄反演的震源低頻信息，如斷層應力降、位錯、破裂速度、破裂時間函數，結合斷層的幾何參數、地質產狀等形成震源信息，通過頻譜傳遞使得模擬結果包含了大量低頻信息.該方法在蘆山7.0級地震強地面運動模擬和九寨溝7.0級地震動模擬中均取得了有效應用.本文利用該方法模擬2021年云南漾濞6.4級地震的強地面運動，通過6個強震臺的模擬地震動與觀測記錄的對比，表明了模擬計算參數選取的有效性.在研究區域內，強地震動模擬獲得的PGV、PSV、PGA和PSA等值線分布近似為以震中為圓心的圓形，進一步討論了云南漾濞6.4級地震的震害強度分布和烈度分布范圍.

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

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

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

(1)

這里nl和nw分別是沿斷層的長和寬方向上子斷層單元的個數，Δtij是第ij子斷層單元到觀測點的相對延遲，aij(t)是利用隨機點源方法計算的第ij子斷層單元產生的地震動加速度，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子斷層的應力降和地震矩.實際的斷層破裂過程中，斷層面上各子斷層單元上的地震矩和應力降都不相同，根據公式(3)可以得到不同子斷層的拐角頻率.這種拐角頻率與動拐角頻率(Motazedian and Atkinson，2005)形式、大小、取值范圍都不相同.Zhou和Chang(2019)證明了非均勻拐角頻率(3)保留了動拐角頻率(Motazedian and Atkinson，2005)的優點，同樣避免斷層不同剖分尺度對于遠場輻射能量的影響.

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

2 漾濞6.4級地震模擬參數的確定

2.1 震源參數

圖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 場地系數κ

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

κ=a1R+a0,

(4)

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

圖3 研究區域和強震臺站分布白色線表示研究區域邊界，紅色三角形表示本文中模擬的臺站，白色實心點表示震中位置.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個臺站記錄振幅譜對數(a) 加速度EW分量; (b) 加速度NS分量； 藍線為加速度振幅譜對數，紅線為截止頻率段的線性擬合.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 衰減指數λ

(5)

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

2.4 附加時窗參數

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

(6)

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

dtij=tEW-tij，

(7)

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

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

2.5 模擬參數匯總

經過對實際地震記錄的分析、反演和擬合等處理工作，得到了模擬地震動需要的所有計算參數，匯總如表1：

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

3 模擬結果與實際結果的對比

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

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

圖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

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

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

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

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

圖11d是研究區域內PGA分布，分布形態近似圓形，中心點鄰近漾濞縣城，最大值超過400 cm·s-2，整個斷層位于300 cm·s-2等值線內，60 km以外地震動出現快速衰減，在很短的距離內由100 cm·s-2降為50 cm·s-2表現為等值線變密，震中距80 km以外，PGA減小到20 cm·s-2以下.圖12d是PGV分布，也呈圓形狀態，震中附近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左右.三個周期反應譜分布盡管幅值相差較大，但總體形態上具有相似之處，呈現圍繞震中圓形分布，震中60 km附近等值線變密，變化加大.對于圖12中PSV三個周期震中附近的值在40 cm·s-1以下，三個周期0.3 s、1 s和3 s變化幅度不大.

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

5 結論

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

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

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

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

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

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