典型地震斷層錯動導致的重力變化及2022門源6.9級地震重力變化分析

鐘慧婷, 靳志同,2, 萬永革,2

(1. 防災科技學院, 河北 三河 065201;2. 河北省地震動力學重點實驗室, 河北 三河 065201)

0 引言

確定斷層的破裂方式對研究區的整個地震孕育過程研究及其震后地震趨勢的發展均起到非常重要的作用,而重力場是反映地球介質密度變化和各種環境下動力學特征的最基本、最直接的物理量,重力-地震的聯合反演具有顯著的互補效果和更完備額的地質學解釋[1],所以重力資料在反演震源參數和震源破裂中起著不可或缺的作用。如賈可可等[2]利用同震重力變化反演出日本“3·11”地震的震源破裂參數。鄭增記等[3]基于GRACE記錄的重力資料反演得到2012年蘇門答臘地震斷層參數。楊君妍等[4]通過GRACE重力衛星獲取的數據反演出了地震位錯Love數,并且得到了更準確的格林函數,反演出孕震區域局部構造信息。

除此之外,在地震預報領域,重力資料也發揮了一定的作用。祝意青等指出強震易發生在與構造活動有關聯的重力變化正、負異常區過渡的高梯度帶上,重力變化等值線的拐彎部位或四象限中心附近[5]。2001年昆侖山口西8.1級地震發生在重力變化高梯度帶上,重力負值變化低于-40 μgal[6];2021年云南漾濞6.4級地震震中處于重力變化的高梯度帶區域,判斷漾濞地區地下活動的能量已經完全釋放完,震后短期內無更大地震發生的可能[7]。這些研究表明重力資料可以為地球動力學過程、地震危險性判斷提供重要信息。由于重力監測資料中包含了地震斷層錯動產生的重力變化值,而地震產生的重力變化會對揭示地震孕育等長期動力學行為造成干擾,因此將地震導致的重力變化扣除更有利于揭示地震孕育等地球動力學過程。

2022年1月8日,青海門源發生了6.9級地震,震中位于101.28°E、37.82°N,震源深度為10 km,造成了房屋建筑受損、蘭新高鐵受損停運等問題。此次地震位于青藏高原東北部邊緣的祁連山地震帶,由一系列的深大活動斷裂組成。作為地震活動主體區域的祁連山次級地塊,四周被深大走滑活動斷裂包圍,形成一個相對獨立的較為活動的次級地塊[8]。此次6.9級地震的發震構造既不是冷龍嶺斷裂,也不是托萊山斷裂,而是發育在它們之間階區中的道溝斷裂[9-10],主要是由淺層的左旋走滑斷層造成的。

(紅色沙灘球表示本次地震主震的震源機制解,藍色三角號表示研究區域的前兆臺,紅色線代表研究區域內不同的斷層,左下角子圖表示研究地區所在的位置)圖1 研究區地形圖Fig.1 Topographic map in the study region

近100年來,該地區已發生高于6.0級的地震20余次,其中發生的距離此次地震最近的一次大地震是2016年門源縣6.4級地震。地震斷層錯動會導致地面重力變化,而中國地震前兆臺網的重力儀可以觀測到這種變化。這種變化信息對確定斷層的錯動方式、大小等都起到非常重要的作用。地震前后重力變化的原因一般有兩方面:一是震前地震孕育過程中的質量遷移和地表形變引起的重力變化,二是同震位錯引起周圍物理密度的變化;前一種原因往往是重力變化的主要因素,后一種原因導致的重力變化一般在距離震中幾十千米以外處只有約10 μgal量級且隨著距離的增加迅速衰減[11]。

2022年青海門源發生6.9級地震,該地震發生在周圍地震前兆臺布設比較密集的地區,為了甄別重力前兆數據中的該地震錯動引起的重力變化,從而為提取有意義的重力地震前兆提供基礎數據,并了解該地震對周圍物質運移和密度變化的影響,本文采用Okubo地震位錯模型,分析了各種類型的地震錯動引起的地表重力變化的特征,并計算了2022年青海門源縣6.9級地震在周圍地區引起的重力變化。

1 Okubo理論及計算公式

國內外許多學者對位錯和斷層引起的重力變化進行過深入研究。Stekettee[12]最早將位錯理論引入地震學,Maruyama[13]給出了垂直和水平伸展斷層引起的地表位移場的完整解析解,陳運泰等[14-15]討論了結合半空間無限空間位錯理論和地表形變資料進行震源反演的一般方法,Okada[16]總結了前人的工作,給出了一套完整簡潔實用的同震變形計算公式。Okubo[11,17-18]研究了半無限空間介質內剪切和張裂斷層的重力變化問題,導出了點源和有限斷層的同震重力位和重力變化的解析表達式,Sun[19]解決了球形地球分層模型中地震位錯的同震重力變化計算問題。李振洪等利用InSAR技術獲取2022年青海門源地震的同震形變場,可以發現本次地震形變僅限局部地區,同時利用了彈性半空間的位錯模型確定了本次地震事件的斷層幾何參數[20],鑒于2022年青海門源地震產生的重力變化僅限于局部地區,本文介紹Okubo的半空間線彈性介質中地震位錯產生的重力變化計算方法。

Okubo理論主要計算了有限矩形斷層的重力變化,根據最后的推導結果,用雙豎約定符號標示為簡潔模式[21]:

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

(1)

重力變化可以表示為:

Δg(x1,x2)={ρG[U1Sg(ξ,η)+U2Dg(ξ,η)+
U3·Tg(ξ,η)]+ΔρGU6g(ξ,η)}‖-
βΔh(x1,x2)

(2)

式中:G為萬有引力常數;ρ為介質密度;Δρ=ρ′-ρ,ρ′為引張破裂后的填充物質的密度。S、D、T、C則表示走滑、傾滑、引張和引張破裂填充物的貢獻,(Sg,Dg,Tg,Cg)是對(S,D,T,C)進行微分得到的,β≈0.309 mgal/m為自由空氣重力梯度,Δh為高程變化,它們的具體表達式為:

(3)

(4)

(5)

Cg(ξ,η)=2I2cosδ-sinδ·lg(R+ξ)

(6)

q=x2sinδ-(d-x3)cosδ

(7)

(8)

同樣,高程變化可以表示為:

(9)

具體表達式為:

(10)

(11)

(12)

其中:

(13)

(14)

(15)

I5(ξ,η)=2(1-2v)I1secδ

(16)

當cosδ=0時:

(17)

(18)

式中:δ為傾角;v為泊松系數;d為有限矩形斷層的深度。

2 典型錯動斷層在其周圍產生的重力變化規律

為了解典型斷層引起的重力變化規律,本小節計算了單個左旋走滑斷層、正斷層、逆斷層的重力變化情況,并且進行了分析。走向上單個斷層長度為10 km,傾角方向上的斷層寬度為10 km,在滑動方向上的位錯為5 m,取介質密度為2 670 kg/m3,計算結果分析如下:

2.1 走滑斷層重力變化特征

圖2為單個走滑斷層引起的重力變化分布圖,走向為正東,傾角為90°,滑動角為0°,其中AA′為計算圖3的剖面位置。圖3為破裂中心位于不同深度的地表AA′剖面的重力變化曲線。

圖2 走滑斷層重力變化(AA′為圖3的剖面位置)Fig.2 Gravity changes induced by strike-slip fault (AA′is the position of profile in Fig.3)

圖3 同一剖面破裂中心位于不同深度的走滑斷層產生的重力變化(黑色水平線代表零重力值)Fig.3 Gravity changes along profile induced by strike-slip fault with fracture center at different depths (The black horizontal line represents the zero gravity value)

由圖2可以看出,垂直走滑斷層產生的重力變化總體呈現出四象限對稱分布;在近處,以斷層跡線逆時針算起,一三象限重力變化為正值,二四象限為負值;在遠處,一三象限重力變化為負值,二四象限為正值。圖3顯示,同一剖面重力變化曲線關于斷層所在位置對稱,重力變化大小隨著遠離斷層而逐漸縮小;對于相同破裂尺度和滑動量的垂直走滑斷層,破裂中心深度越深,震中附近的重力變化幅度越小,而自斷層位置向遠處的重力變化幅度衰減越慢。反之,震源破裂中心越淺,斷層產生的重力變化局限在斷層附近,峰值較大,但很快衰減。也就是說,相同破裂尺度和破裂量的淺表震源引起的重力變化值局限在震中的局部地區,而相同破裂尺度和破裂量的深部震源引起的重力變化會產生更大范圍的重力變化,但重力變化值較小。

2.2 傾滑斷層重力變化特征

圖4與圖7為單個正斷層(滑動角為-90°)和逆斷層(滑動角為90°)引起的重力變化分布圖,走向為正東,傾角為45°。以正斷層為例,圖5為相同破裂中心位于不同深度的AA′剖面上的重力變化曲線,圖6為破裂中心位于6 km處,不同傾角的斷層破裂在AA′剖面產生的重力變化。逆斷層同理,圖像為7～9。

圖4 正斷層重力變化(AA′為圖5的剖面位置)Fig.4 Gravity variation induced by normal fault (AA′ is the position of profile in Fig.5)

圖5 同一剖面破裂中心位于不同深度的正斷層產生的重力變化Fig.5 Gravity changes along profile induced by normal fault with fracture center at different depths

圖6 同一剖面破裂中心位于6 km處不同傾角正斷層重力變化Fig.6 Gravity changes along profile induced by normal fault with fracture center at a depth of 6 km and different dip angles

圖7 逆斷層重力變化(AA′為計算8的剖面位置)Fig.7 Gravity variation induced by reverse fault (AA′ is the position of profile in Fig.8)

由于傾滑斷層在震源近處引起的重力變化值過大,為了更好地描述遠處的重力變化特征,本文將重力變化大于5 μgal的值不顯示,如圖4和7。可以發現,在遠處,正斷層在震源南北方向產生的重力變化為正值,在東西方向產生的重力變化為負值,逆斷層反之。圖5和圖8顯示,上盤重力變化的數值總體上比下盤重力變化的數值小;對于相同破裂尺度和滑動量的傾滑斷層,破裂中心深度越深,遠處的重力變化最大值越小。

圖8 同一剖面破裂中心位于不同深度的逆斷層的重力變化Fig.8 Gravity changes along profile induced by reverse fault with fracture center at different depths

從圖6和圖9可看出,在傾角小于45°時,斷層上盤產生的重力變化峰值較大,且隨著斷層的傾角的增大,在斷層上盤重力變化波峰越大;在傾角大于45°時,斷層下盤產生的重力變化峰值較大,且隨著傾角的增大,斷層下盤重力變化的峰值變小。當正斷層傾角較小時,上盤重力變化主要以負值為主,下盤重力變化主要以正值為主;當傾角較大時,上盤重力變化主要為正值,下盤為負值;逆斷層反之。

圖9 同一剖面同一深度不同傾角逆斷層重力變化Fig.9 Gravity variation along profile induced by reverse fault with fracture center at a depth of 6 km and different dip angles

3 2022年門源地震產生的重力場

本節采取Okubo理論,震中周圍介質密度選取為2 670 kg/m3,地震破裂模型采用李振洪等[20]根據InSAR資料聯合了升、降軌地表形變場反演的地震破裂模型;反演模型包含走向104°和109°兩個斷層段,走向104°的斷層段長度為10 km,寬度為16 km,傾角為80°,滑動角為0°,最大滑動量為2.5 m;走向為109°的斷層段長度為20 km,寬度為16 km,傾角為80°,滑動角為5°,最大滑動量為3 m[20]。將發震斷層離散成子斷層反演得到:該斷層主要受左旋滑動控制,屬于高傾角走滑型地震;子斷層最大滑動量為3.5 m,出現在地下約4 km處,滑動主要集中在地下2～7 km的區域[20]。本文利用多個子斷層疊加的方式,計算出重力變化情況,繪制出圖10。

(黑點表示城市,藍色三角號表示未能監測到重力變化的臺站,紅色三角號表示能夠檢測到重力變化的臺站,黑色花瓣輪廓線是重力變化為1μgal的分界線,沙灘球是本次地震的震源機制解)圖10 2022年1月8日門源地震后重力變化Fig.10 Gravity changes generated by Menyuan earthquake on January 8, 2022

圖10顯示,本次地震重力增加的區域范圍要比重力減小的區域范圍小。在數值上,重力變化負值區最大為-8.82 μgal,正值區最大變化為28.28 μgal。震中東北和西南方向的重力減小,震中西北和東南方向的重力增加;即在距離震中遠處斷層錯動前方重力變化為正,后方重力變化為負,與前文敘述單個走滑斷層引起的重力變化特征基本一致。

重力變化一般在距離震中幾十千米以外處只有約10 μgal的量級,并隨著距離的增加迅速減小[11]。本次門源地震斷層錯動導致的重力變化從震源處向四周衰減得很快,在數百公里之外幾乎不發生重力變化。當前CG-5高精度重力儀在重力觀測上分辨率可以達到1 μgal[22],圖10中,黑色花瓣輪廓線是重力變化為1 μgal的分界線,花瓣內是能夠采用重力儀監測到地震錯動產生重力變化的地區,震中東北和西南方向、西北和東南方向的花瓣輪廓最遠大約在100 km處,超出該花瓣范圍的區域將難以監測到重力變化。

圖10中有63個用三角號表示的前兆臺,紅色三角(門源臺)標志表示能夠監測到此次門源地震重力變化的臺站,藍色為不能夠監測到重力變化的臺站。需要指出的是,在已知臺站中,有一個S9(門源臺)前兆臺站能夠監測到本次門源地震的重力變化,其余臺站處重力變化低于1 μgal,超過了重力儀監測的范圍,所以未能監測到重力變化。即便在本次門源地震中只有一個前兆臺能監測到重力變化,在監測和分析該地區長期重力變化從而揭示構造活動時,也應該扣除該地震造成的重力變化。

重力場的變化,能較好地反映地殼厚度的差異、地殼密度的變化和深部物質遷移等構造活動信息,重力場隨時間變化與地震的形成和發展有著內在聯系[23],所以重力資料在反演震源參數中也起著重要的作用。當然,單靠重力資料也是遠遠不夠的,需要聯合其他多種資料共同約束,從而得到更加準確的震源參數。

4 討論與結論

為分析地震錯動在地表導致的重力變化,本研究根據Okubo位錯理論公式,計算分析了典型斷層(走滑、傾滑)在地表的重力變化特征:(1)垂直走滑斷層產生的重力變化總體呈現出四象限對稱分布;對于同一剖面,重力變化大小隨著遠離斷層而逐漸縮小;對于相同破裂尺度和滑動量的垂直走滑斷層,深度越深,震中附近的重力變化幅度越小,而自斷層位置向遠處的重力變化幅度的衰減越慢。反之,震源越淺,斷層產生的重力變化局限。(2)傾滑斷層上盤重力變化的數值總體上比下盤重力變化的數值小;對于相同破裂尺度和滑動量的傾滑斷層,破裂中心深度越深,遠處的重力變化最大值越小。同時本文對2022年1月8日門源6.9級地震的破裂分布進行重力變化的計算和分析,得到:本次地震重力增加的區域范圍要比重力減小的區域范圍小;在數值上,重力變化負值區最大為-8.82 μgal,正值區最大變化為28.28 μgal,震中東北和西南方向重力降低,而震中西北和東南方向的重力增加。位于震中附近的門源臺的重力增加3.40 μgal,能夠被重力儀器監測到。

本文計算采用是李振洪等[20]的破裂模型,隨著資料的公開和共享,更加精確的震源模型會改變本文的計算結果,但本文基于InSAR的數據反演模型抓住了主要的特征,即使后續有更為精確的破裂模型,對我們的計算結果不會較大影響。

余震也會對本文的計算結果也會有一定的影響。萬永革等在計算蘭德斯地震余震產生的位移場時發現,地震斷層面及其附近余震產生了厘米量級的位移變化[24],對計算地表重力變化的數值影響較小;并且目前沒有更小地震的震源機制和破裂模型數據,所以本文中沒有計算更小地震對地表重力變化的影響。

除了上文所說的震源破裂模型、余震會對計算結果有影響,本文采取震源周圍介質為均勻彈性半空間,密度統一為2 670 kg/m3,這也會對計算結果造成一定的影響。在計算門源地震后的重力變化時,同樣還采用了密度為2 630 kg/m3、2 650 kg/m3、2 690 kg/m3進行計算,但是最后的計算結果影響都比較小,重力變化特征也基本一致。在不同密度的計算結果中,也是只有一個S9(門源臺)前兆臺站能夠監測到本次門源地震的重力變化,其余臺站處重力變化比較小。同樣震后余滑或介質黏彈性松弛等產生的影響也較小,所以可以忽略不計。

雖然上述種種問題還有待于今后進一步解決,但本文在一級近似情形下給出了不同類型斷層錯動在地表產生的重力變化模式,并根據2022年門源地震的破裂模型計算得到了該地震產生的重力變化,這對采用重力資料反演該地震的破裂模型和研究長期重力變化從而揭示其他動力學過程是有益的。