金雯 談洪波 申重陽 申宇彤

1)中國地震局地震研究所,武漢 430071

2)湖北省地震局,武漢 430071

0 引言

據(jù)中國地震臺網(wǎng)測定,北京時間2022年9月5日12時52分四川甘孜藏族自治州瀘定縣(29.59°N,102.08°E)發(fā)生M6.8地震,此次地震是繼漾濞M6.4、瑪多M7.4、蘆山M6.1、馬爾康M6.0地震之后,近2年來在青藏高原中東部地區(qū)發(fā)生的第5次M≥6.0地震(李傳友等,2022),也是近40年來四川地區(qū)發(fā)生的M6.0～6.9地震中人員震亡最多的一次(許娟等,2023)。瀘定M6.8地震位于青藏高原東南緣鮮水河斷裂帶南東段的磨西斷裂帶附近,系鮮水河、安寧河和龍門山斷裂的交匯處(圖1)(單新建等,2023)。鮮水河斷裂帶是巴顏喀拉塊體與川滇菱形塊體的邊界,具有左旋走滑活動特征(郭長寶等,2015; 孫東等,2023)。研究表明,本次地震的發(fā)震斷層磨西段是全新世強烈活動的斷裂,其位于貢嘎山東側(cè)、大渡河西側(cè),斷裂沿線地質(zhì)條件和地貌地形較為復(fù)雜,歷史地震頻發(fā)(熊探宇等,2010; 李彩虹等,2022)。

注: 時間為1900—2022年,紅色五星表示瀘定M6.8地震,下同。

關(guān)于瀘定M6.8地震的破裂模型,已有一些地震學(xué)研究結(jié)果,如中國地震局地震預(yù)測研究所給出的結(jié)果[注]https://www.ief.ac.cn/kydts/info/2022/69528.html.顯示,瀘定地震為左旋走滑剪切型破裂,破裂面呈NW-SE走向,同時向南北兩側(cè)擴展,破裂方向與鮮水河斷裂帶南段磨西斷裂延伸方向基本一致,最大滑動量達1.84m; 中國地震局地球物理研究所[注]https://www.cea-igp.ac.cn/kydt/279423.html.利用全波形數(shù)據(jù)的迭代反褶積與疊加方法來反演地震破裂過程,得到斷層破裂模型,地震矩震級為6.7級; 王衛(wèi)民等[注]http://www.itpcas.cas.cn/new＿kycg/new＿kyjz/202209/t20220906＿6509485.html.通過遠場體波波形數(shù)據(jù)獲得了該地震的震源破裂過程模型,得出本次地震為高傾角走滑型事件,震源深度約為13km,最大滑動量為1.92m。這些地震學(xué)結(jié)果為更深入地理解瀘定地震機制和理論模擬研究提供了參考。

根據(jù)震前重力場變化圖像,中國地震局重力學(xué)科組對此次地震進行了較好的預(yù)測,瀘定地震正好位于重力學(xué)科2021年底劃定的年度重點危險區(qū)內(nèi)[注]湖北省地震局. 2021. 2022年度重力學(xué)科地震趨勢會商報告.,其主要判據(jù)之一為2019年9月—2020年9月重力測量結(jié)果(圖2):顯著正變化呈現(xiàn)NE向條帶狀,正負最大差異達 90×10-8m/s2; 震前重力變化大致以石棉—康定為中心呈現(xiàn)四象限分布特征。這種震前特征是否與同震破裂機制存在聯(lián)系,是否可用申重陽等(2011)提出的閉鎖剪力模式來解釋,是值得研究的問題。為此,本文基于彈性半無限空間矩形位錯理論,采用張勇②和王衛(wèi)民③發(fā)布的由遠場地震波反演的瀘定M6.8地震發(fā)震斷層模型,從理論上模擬了此次地震引起的地表同震重力和形變效應(yīng),并與實際觀測資料進行對比分析,以期理解瀘定地震的發(fā)震和孕育機制,為地震預(yù)測和地震危險性評估提供實例參考。

圖2 瀘定地震前實測重力場變化

1 同震位錯引起的地表重力和形變變化的理論公式

自Steketee(1958)將位錯理論引入地震學(xué)以來,很多學(xué)者研究了半無限空間均勻介質(zhì)地球模型的同震變形問題。Okada(1985)總結(jié)并整理了前人的工作,給出一套完整簡潔實用的同震變形計算公式,適用于計算任何剪切與張性斷層引起的位移、應(yīng)變和傾斜變形; Okubo(1991、1992)研究了半無限空間介質(zhì)內(nèi)剪切和張裂斷層的重力變化問題,利用類似于Okada(1985)的方法導(dǎo)出點源和有限斷層的同震重力位和重力變化的解析表達式; Sun(1992)和Sun等(1993)基于地球曲率和層狀結(jié)構(gòu),提出球形地球模型位錯理論; Fu(2007)研究了地球內(nèi)部橫向不均勻結(jié)構(gòu)的影響,對球形位錯理論進行了補充?？紤]到地震破裂的局域性和地震孕育過程的能量積累主要以彈性應(yīng)力應(yīng)變形式進行,地震位錯引起的重力及形變變化往往具有以發(fā)震斷層為中心、隨距離增大而較快衰減的特性,因此本文選擇彈性平面位錯理論進行研究。

1.1 單條斷層錯動

Steketee(1958)指出,在各向同性的彈性半無限空間介質(zhì)中,一個沿斷層面的位錯Σ在空間中任意一點(x1,x2,x3)產(chǎn)生的位移場為

(1)

在O-x1x2x3(左旋)直角坐標(biāo)系下,x1x2為水平地面,x3垂直向下; 設(shè)矩形斷層(斷層下盤)的長和寬分別為L和W,δ為斷層面傾角,U1、U2、U3分別對應(yīng)于任意位錯的走滑、傾滑和引張位錯分量,d為斷層底界深度(圖3),在自由地表某點(x1,x2,0)引起的重力變化Δg和地表變形Δu可以分別表示為(Okubo,1992)

Δg(x1,x2)={ρG[U1Sg(ξ,η)+U2Dg(ξ,η)+U3Tg(ξ,η)]+

ΔρGU3Cg(ξ,η)}‖-βΔh(x1,x2)

(2)

(3)

其中,G為萬有引力常數(shù);ρ為介質(zhì)密度; Δρ為張裂紋內(nèi)密度與介質(zhì)密度之差;Sg(ξ,η)、Dg(ξ,η)、Tg(ξ,η)、Cg(ξ,η)和S(ξ,η)、D(ξ,η)、T(ξ,η)為系數(shù)(Okubo,1992;Okada,1985); 自由空氣重力梯度β=0.3086×10-5/s2; Δh為地表高程變化(Okada,1985); 雙豎線符號“‖”可表示為(Okubo,1992)

f(ξ,η)‖=f(x1,p)-f(x1,p-W)-f(x1-L,p)+f(x1-L,p-W)

(4)

其中,p與斷層傾角δ關(guān)系式為

p=x2cosδ+dsinδ

(5)

1.2 多斷層錯動

設(shè)一地表坐標(biāo)系O-XYZ,X軸指向正東,Y軸指向正北,Z軸垂直向下,S(sx,sy,0)為斷層頂部地表跡線(或投影線)中點(圖4)。設(shè)斷層方位角為α,則由N條斷層位錯引起的地表任意一點(x,y,0)的地表形變u和重力變化ΔG可表示為多條單斷層模型的疊加,即

(6)

圖4 斷層運動地表坐標(biāo)系示意圖

(7)

由公式(6)、(7)可以得到,斷層尺寸(長、寬)、斷層產(chǎn)狀(傾向、傾角、走向)、斷層位置(埋深、方位、中心投影坐標(biāo))以及錯動性質(zhì)(錯動量的大小、錯動方式)是影響地表形變效應(yīng)和重力變化的關(guān)鍵因素。

2 瀘定M6.8地震同震破裂模型

選取張勇②給出的同震破裂模型(簡稱張勇模型,下同)和王衛(wèi)民③給出的同震破裂模型(簡稱王衛(wèi)民模型,下同)分別進行計算分析(圖5)。2個破裂模型包含多個子斷層,每個子斷層的參數(shù)包括中心位置(經(jīng)度、緯度)、中心深度、滑動量、滑動角、方位角和傾角。

圖5 瀘定M6.8地震滑動分布模型

張勇模型斷層面上的靜態(tài)滑動分布如圖5(a)所示,結(jié)果顯示震中29.59°N、102.08°E,震源深度16km,越靠近最佳矩心震源,滑動量越大,最大靜態(tài)位錯約為1.6m。其發(fā)震斷層的參數(shù)為:走向163°、傾角80°、滑動角8°,該斷層分別沿走向和傾向方向均勻地分成13×7塊子斷層,每個子斷層的尺度為3km×3km,并給出了每個子斷層的滑動角和滑動量。

王衛(wèi)民模型斷層面上的靜態(tài)滑動分布如圖5(b)所示,結(jié)果顯示震中29.59°N、102.08°E,震源深度13km,最大靜態(tài)位錯約2m。發(fā)震斷層參數(shù)為:走向166°、傾角78°、滑動角約0°,該斷層分別沿走向和傾向方向均勻地分成17×9塊子斷層,每個子斷層的尺度為3km×2km。

這2種模型均是基于瀘定地震的遠場地震波形數(shù)據(jù)反演得到的,盡管其發(fā)震斷層均表現(xiàn)為左旋走滑特征,但在斷層破裂長度、寬度、震源深度、滑移量等方面具有明顯差異。

3 地表同震重力變化模擬結(jié)果

基于上述理論公式和2種斷層破裂模型,模擬計算瀘定M6.8地震引起的地表同震重力變化。根據(jù)同震效應(yīng)影響的范圍,選取計算區(qū)域為28°N～31°N、100.5°E～103.5°E,單元網(wǎng)格尺寸為0.025°×0.025°。

根據(jù)張勇模型的模擬結(jié)果(圖6(a)),地表重力變化圖像呈現(xiàn)明顯的四象限分布特征,以發(fā)震斷層為界,近震區(qū)(距發(fā)震斷層小于50km)以東區(qū)域呈現(xiàn)正的重力變化,最大正變化量約7×10-8m/s2,震中以西區(qū)域呈現(xiàn)負的重力變化,最大負變化量約-14 ×10-8m/s2; 遠震區(qū)(距發(fā)震斷層大于50km)的北東區(qū)域和西南區(qū)域呈現(xiàn)正重力變化,而西北區(qū)域和東南區(qū)域呈負重力變化。重力變化峰值出現(xiàn)在鮮水河斷裂南東端和安寧河斷裂北端交匯處。

圖6 瀘定M6.8地震斷層位錯引起的同震重力變化

根據(jù)王為民模型的模擬結(jié)果(圖6(b)),其遠震區(qū)的重力變化與張勇模型模擬結(jié)果基本特征一致; 近震區(qū)存在差異,主要體現(xiàn)在以發(fā)震斷層為界,震中EW向呈現(xiàn)正向變化,最大正變化量約6×10-8m/s2,SN向則呈現(xiàn)負向變化,最大負變化量約-8 ×10-8m/s2。

從圖6可以看出,2種模型計算的同震重力變化具有共性和差異性。共性主要表現(xiàn)在遠場,2種模型的模擬結(jié)果均具有相似的四象限特征,顯示出斷層源以左旋為主的走滑特征。差異性主要表現(xiàn)在震中近場細節(jié),張勇模型結(jié)果顯示震中以東為正變化,以西方向表現(xiàn)出負向變化; 王衛(wèi)民模型結(jié)果展示了震中SN向呈現(xiàn)負向變化,EW向則呈現(xiàn)正向變化。這種近場差異可以歸因于2種斷層模型的細節(jié)差異,首先,可能是因為張勇模型聚焦在單一破裂集中區(qū),而王衛(wèi)民模型涵蓋了2個破裂集中區(qū); 其次,張勇模型主要呈現(xiàn)走滑運動,并帶有逆沖分量(滑動角為8°),而王衛(wèi)民模型則主要表現(xiàn)為純走滑特征(滑動角為0°)。圖6 震中近場主要呈現(xiàn)局部“高頻”變化,目前還難以用重力觀測等手段予以驗證。同震重力變化峰值位置均出現(xiàn)在鮮水河斷裂南東端和安寧河斷裂交匯處; 同震造成震中東北和西南區(qū)域重力值繼續(xù)增大,表明這些區(qū)域受到的擠壓效應(yīng)在地震過程中不斷增強,處于物質(zhì)累積狀態(tài)。此外,發(fā)現(xiàn)這種擠壓效應(yīng)的增強與震前重力正變化區(qū)域的位置高度重合,這種疊加效應(yīng)可能會使相關(guān)區(qū)域地震危險性增強。

4 地表同震位移模擬結(jié)果

模擬計算瀘定M6.8地震引起的區(qū)域同震位移,結(jié)果如圖7所示。

圖7 瀘定M6.8地震引起的同震地表變形

由水平經(jīng)向位移等值線(圖7(a))可以看出,水平經(jīng)向位移場具有北正南負的四象限特征。以發(fā)震斷層為界,西北和東北區(qū)域為正變化,西南和東南區(qū)域為負變化。在近震區(qū),該位移場的幅值衰減速度較快(張勇模型結(jié)果:0.3～1.3cm,-0.1～-1.1cm; 王衛(wèi)民模型結(jié)果:0.4～1.6cm,-0.2～-1cm),而在遠震區(qū)位移場的幅值衰減速度則相對緩慢。震中WN向出現(xiàn)正峰值,張勇模型模擬結(jié)果約為1.3cm,王衛(wèi)民模型模擬結(jié)果約為1.6cm; 震中ES向則出現(xiàn)負峰值,張勇模型模擬結(jié)果約為-1.1cm,王衛(wèi)民模型模擬結(jié)果約為-1cm。在震中附近的地區(qū),漢源縣表現(xiàn)出西向位移(張勇模型結(jié)果:1cm; 王衛(wèi)民模型結(jié)果:0.9cm),石棉縣也表現(xiàn)出西向位移(張勇模型結(jié)果:0.5cm; 王衛(wèi)民模型結(jié)果:0.5cm),而瀘定縣則出現(xiàn)東向位移(張勇模型結(jié)果:0.4cm; 王衛(wèi)民模型結(jié)果:0.5cm)。

由水平緯向位移等值線(圖7(b))可以看出,水平緯向位移場具有東正西負的四象限特征。以發(fā)震破裂斷層為界,北東、南東區(qū)域為正變化,而西南、西北區(qū)域則為負變化。近震區(qū)的幅值也顯示出快速衰減的特點(張勇模型結(jié)果:0.1～0.9cm,-0.3～-1.7cm; 王衛(wèi)民模型結(jié)果:0.2～1cm,-0.4～-1.6cm)。震中EN向出現(xiàn)正峰值,張勇模型結(jié)果約為0.9cm,王衛(wèi)民模型結(jié)果約為1cm; 而在WS向則出現(xiàn)負峰值,張勇模型結(jié)果約為-1.7cm,王衛(wèi)民模型結(jié)果約為-1.6cm。同時,漢源縣出現(xiàn)北向位移(張勇模型結(jié)果:0.4cm; 王衛(wèi)民模型結(jié)果:0.5cm),石棉縣也出現(xiàn)北向位移(張勇模型結(jié)果:0.3cm; 王衛(wèi)民模型結(jié)果:0.3cm),而瀘定縣同樣出現(xiàn)北向位移(張勇模型結(jié)果:0.8cm; 王衛(wèi)民模型結(jié)果:0.9cm)。

由垂直位移等值線(圖7(c))可以看出,垂直位移場具有以震中為中心的正負四象限對稱分布特征,垂直位移主要分布在發(fā)震斷層的近震區(qū),東北、西南區(qū)域為正變化,西北、東南區(qū)域為負變化,近場幅值變化較大(張勇模型結(jié)果:0.2～1.8cm,-0.2～-1cm; 王衛(wèi)民模型結(jié)果:0.3～1.5cm,-0.1～-1.5cm)。其中,漢源縣出現(xiàn)沉降位移(張勇模型結(jié)果:0.3cm; 王衛(wèi)民模型結(jié)果:0.3cm),石棉縣也出現(xiàn)沉降位移(張勇模型結(jié)果:0.2cm; 王衛(wèi)民模型結(jié)果:0.2cm),而瀘定縣出現(xiàn)隆起位移(張勇模型結(jié)果:0.3cm; 王衛(wèi)民模型結(jié)果:0.5cm)。

水平位移以及垂直位移圖像的對稱性主要體現(xiàn)了同震位錯的走滑型特征,不同震源模型造成了近震區(qū)和遠震區(qū)移幅值的差異。地表同震重力變化與垂直位移圖像在遠場區(qū)域呈現(xiàn)正相關(guān),擠壓隆起區(qū)重力值增大,拉張凹陷區(qū)重力值減少,說明該區(qū)域同震重力變化主要受斷層位錯引起的地殼密度變化的影響; 地表同震重力變化與垂直位移圖像在近場區(qū)域呈現(xiàn)負相關(guān),表明地表垂直位移是該區(qū)重力變化的主要因素。

5 討論與結(jié)論

本文基于Okubo平面矩形彈性位錯理論,采用已有的由地震波反演獲得的同震破裂模型,模擬研究了瀘定M6.8地震產(chǎn)生的地表同震重力變化、垂直位移和水平位移。結(jié)果表明:

(1)瀘定M6.8地震同震重力變化和形變分布形態(tài)與走滑型斷層位錯引起的地表重力與位移變化的理論結(jié)果(談洪波等,2009)具有較好的一致性,表明瀘定地震發(fā)震斷層主要呈現(xiàn)左旋走滑特征。

(2)瀘定M6.8地震孕震過程符合閉鎖剪力孕震模式。地震孕育發(fā)生過程會伴隨重力場的變化(申重陽等,2003),震前重力變化可反映孕震現(xiàn)象(申重陽等,2009、2011)。申重陽等(2011、2018)根據(jù)2009年姚安M6.0地震前的典型重力變化與震源破裂機制解的比較分析,提出了地震孕育的閉鎖剪力模式:孕震區(qū)在外圍動力作用下應(yīng)先存在一定大小、優(yōu)勢方向明顯的剪應(yīng)力以及其作用下形成與之關(guān)聯(lián)的剪應(yīng)變,由于該優(yōu)勢剪切力及其剪應(yīng)變未超越巖體破裂閾值,使震前孕震區(qū)處于“閉鎖”或相對平衡的狀態(tài)。通過對大量震例的總結(jié)研究(申重陽等,2009、2010、2011、2020; 祝意青等,2020),認(rèn)為強震前孕震區(qū)重力場變化具有典型的梯度帶或四象限分布特征。瀘定M6.8地震前1年,震中區(qū)存在明顯的重力場變化正負四象限分布特征:以瀘定為中心,北東區(qū)域和西南區(qū)域受力擠壓呈正重力變化,西北區(qū)域和東南區(qū)域為拉張區(qū)域呈負重力變化,且震前重力變化和同震重力變化均具有類似的四象限特征,說明震前孕震源與同震破裂具有一定聯(lián)系,這種聯(lián)系與姚安M6.0地震前重力變化四象限分布特征類似,可能說明引起地震破裂的剪切力在震前就一直存在,即可用閉鎖剪力模式來解釋,這為孕震源的研究提供了新的證據(jù)。重力的震前觀測結(jié)果與同震模擬結(jié)果存在量值差異,可能反映出孕震源與同震破裂源之間的差異。

(3)地表實測GNSS結(jié)果(單新建等,2023)與本文理論模擬結(jié)果的對比,如表1所示。模擬結(jié)果與GNSS實測結(jié)果顯示的變形特性一致,位移整體呈現(xiàn)四象限特征分布:震中西南側(cè)具有遠離震中的SW向同震形變,震中東南側(cè)具有指向震中的NW向同震形變,震中東北側(cè)具有遠離震中的NW向的同震形變,震中西北側(cè)具有指向震中的SE向的同震形變。整體來說,實測結(jié)果與理論模擬結(jié)果在遠場區(qū)域吻合較好,而在近場區(qū)域差異稍大,這可能與本文所采用的發(fā)震斷層模型是基于遠場地震波數(shù)據(jù)反演獲得,其對近場形變控制能力相對較弱有關(guān)。雖然本文選取的斷層模型比較簡單,理論模擬結(jié)果與實測結(jié)果存在一定的偏差,但其基本特征是一致的,本文結(jié)果可為瀘定M6.8地震機理研究提供參考。

表1 模擬結(jié)果與GNSS實測結(jié)果

(4)瀘定M6.8地震對貢嘎山的隆升具有抑制作用。瀘定地震發(fā)生于磨西斷裂帶,該斷裂帶位于貢嘎山東側(cè)。同時,貢嘎山位于青藏高原東南緣,其主峰距震中約56.7km,貢嘎山是青藏高原隆升作用的重要組成區(qū)域之一(劉方斌等,2021)。有研究認(rèn)為貢嘎山花崗巖體是鮮水河斷裂至安寧河斷裂間擠壓彎曲段吸收、轉(zhuǎn)換川滇地塊SE向水平運動導(dǎo)致局部快速隆升的產(chǎn)物,其最南端隆升速率超過(3.3±0.8)mm/a(譚錫斌等,2010)。通過模擬貢嘎山主峰(29°35′44″N, 101°52′44″E)的同震水平形變和垂直形變(表2),2個模型模擬結(jié)果均顯示貢嘎山有向東的位移量(張勇模型結(jié)果:0.73cm; 王衛(wèi)民模型結(jié)果: 0.47cm)和向南的位移量(張勇模型結(jié)果:0.41cm; 王衛(wèi)民模型結(jié)果:0.71cm),垂向位移為負值(張勇模型結(jié)果: -0.66cm; 王衛(wèi)民模型結(jié)果: -0.52cm),表明貢嘎山出現(xiàn)沉降。對比沉降量與隆升速率,瀘定地震對貢嘎山的隆升產(chǎn)生了2年左右的延遲效應(yīng),這表明該地震對貢嘎山隆升產(chǎn)生了一定的抑制作用。因此,此次瀘定地震可能對貢嘎山及其周邊地區(qū)的地表變形和地殼運動產(chǎn)生一定影響。

表2 貢嘎山的同震效應(yīng)模擬結(jié)果

致謝：張勇研究員和王衛(wèi)民研究員提供了發(fā)震斷層模型數(shù)據(jù),審稿專家提出了寶貴修改意見,在此一并表示感謝。