2016年阿克陶MW6.6地震的同震斷層滑動分布及庫倫應力變化

解 柳,吳紹宇,康明哲,曾慶偉,朱毅誠

(1.西安科技大學,陜西 西安 710054;2.陜西鐵路工程職業技術學院,陜西 渭南 714000;3.西安市政設計研究院有限公司,陜西 西安 710068)

0 引言

北京時間2016年11月25日22時24分,在新疆維吾爾自治區阿克陶縣境內發生了MW6.6地震。此次地震震中位于阿克陶縣木吉鄉和庫斯拉甫鄉之間的高寒山區,發震斷層為NWW走向的右旋走滑斷裂,破裂長度超過77 km[1]。不同機構和學者對此次地震斷層運動進行了研究,美國地質調查局(USGS)發布的破裂起始點坐標為(39.273°N,73.978°E),震源深度為17 km,認為它是具有拉張分量的右旋走滑事件。中國地震臺網(http://www.cea.gov.cn/publish/dizhenj/464/515/20161125231007712890916/index.Html)測定的地震震中位置為(39.27°N,74.04°E),震源深度10 km,為右旋走滑兼正斷層性質(表1)。張勇(http://www.ceaigp.ac.cn/tpxw/275080.html)研究的震源破裂結果表明該地震為單側破裂,從震中開始向東傳播。由于該地區地質、地形也相對復雜,有學者認為該地震為單一破裂。利用雷達數據提取LOS向形變圖,參考各機構的震源機制解的參數,進行同震的滑動分布反演[2]。也有學者認為該地震為雙斷層破裂[3],并利用升降軌數據獲取該地震的LOS向形變,以中國地震臺網確定的震源機制解的參數進行滑動分布反演,發現其破裂分別位于震中以東約7 km、約33 km處,主要的滑動量集中在地下深度0～20 km處,最大滑動量為0.84m,并基于彈性半空間形變模型采用兩段非均勻斷層滑動模型進行精細滑動分布反演,結果顯示斷層西段的第一次事件以右旋走滑為主,走向103°,傾角76°,斷層東段的第二次事件為右旋走滑兼具正斷作用,走向109°,傾角略緩約55°[4]。

表1 不同機構給出的阿克陶地震震源機制解Table 1 The focal mechanism solutions ofAketao earthquake given by different organizations

本文在采用貝葉斯算法獲取斷層幾何參數的基礎上,反演了單一斷層的滑動分布;采用粒子群算法獲取雙斷層的最優傾角值,并進行精細滑動分布反演,并以滑動分布反演結果對該地震的庫倫應力變化進行了分析。

1 研究區域及數據

1.1 構造背景

2016年發生的阿克陶MW6.6地震位于印度板塊沿北東向與歐亞板塊碰撞的西觸角地區[5],該區域處于南天山西段和帕米爾高原交會地區、木吉斷陷盆地西端地區附近(圖1)。不同機構給出的震源機制解如下表1。

從2000年以來距震中100 km范圍內共發生了2次6級地震、13次5級地震[6]。自阿克陶主震發生后,截止到2016年11月27日24時,記錄到5.0～5.9級余震 1次、4.0～4.9級余震5次,最大余震震級為 5.0級[1]。

注:黑白色沙灘球表示不同機構發布的震源機制結果,紅色與綠色圓點為余震分布圖1 阿克陶地震震中區域及鄰區地形、斷層、余震Fig.1 Topography,faults,and aftershocks in the epicenter area and adjacent areas of the Aketao earthquake

1.2 數據說明

1.2.1 數據介紹

本文采用歐空局(European Space Agency,ESA)哥白尼計劃中的地球觀測衛星分別于2014年4月和2016年6月發布的衛星數據Sentinel-lA升軌單視復數影像和Sentinel-1B降軌單視復數影像。其具體信息見表2。

表2 InSAR數據信息Table 2 InSAR data information

1.2.2 同震形變解算

利用SLC影像數據生成強度圖并結合DEM影像進行地理編碼、差分干涉處理,生成差分干涉影像。為提高干涉圖的質量,對干涉圖進行自適應濾波處理、相位解纏、地理編碼,最終提取LOS向形變信息(圖2)。在圖2中,(a)為降軌的LOS向形變場分布,地震斷層南側是木吉盆地,其干涉條紋明顯,其西部形變值約為14 cm,東部形變值約為23 cm。而北部山區地形起伏,影像失相干,故干涉條紋模糊,但不影響獲取形變范圍。(b)為升軌的LOS向形變場分布,東西兩部分形變量值相當,約為12 cm。相對于西部條紋分布而言,東部條紋分布范圍比較小,條紋緊密。

圖2 D-InSAR LOS向形變場Fig.2 Deformation field in LOS direction from D-InSAR

1.2.3 剖面分析

為更清晰地了解升降軌LOS向形變場信息,以兩個子破裂為中心取同一經度的剖面,對其數據進行四次多項式擬合(圖3)。淺藍色為A1-A2剖面的形變值,最大沉降值約為25 cm,淺藍色為A1-A2剖面的形變值,最大抬升值約16 cm,東部子破裂區的沉降量明顯大于西部[圖3(a)]。由剖面C1-C2,D1-D2可看出升軌的形變主要以抬升為主,最大抬升值約為13 cm,而最大沉降值僅約6 cm[圖3(b)]。對比圖3,由于升降軌的觀測模式不同,造成升降軌LOS向形變場的沉降、抬升趨勢相反,升軌的形變量值相對降軌的要小一些,但剖面圖的形變趨勢與整體形變趨勢一致。

圖3(a)表示降軌的剖面圖,淺藍色表示以東南部破裂為中心,經度為74.399°,緯度為39.05°至39.39°之間的取值(降軌解纏圖白線A1-A2),黑色為其擬合曲線;深藍色表示以西南部破裂為中心,經度為74.197°,緯度為39.05°至39.39°之間的取值(降軌解纏圖白線B1-B2),紅色為其擬合曲線;圖3(b)表示升軌的剖面圖,淺藍色表示以東南部破裂為中心,經度為74.302°,緯度為39°至39.47°之間的取值(升軌解纏圖白線C1-C2),黑色為其擬合曲線;深藍色以西南部破裂為中心,經度為74°,緯度為39°至39.47°之間的取值(升軌解纏圖白線D1-D2),紅色為其擬合曲線圖3 同震差分剖面圖Fig.3 Co-seismic difference section

2 斷層幾何與滑動分布

本文依據升降軌提取的形變場估計斷層的幾何參數并反演斷層的滑動分布。假設為均勻滑動模型,輸入相應的斷層幾何參數。將斷層面劃分成若干子斷層片,使用線性反演方法來估計斷層面的滑動分布。

2.1 斷層幾何與滑動參數

本文結合矩形位錯模型[7]采用貝葉斯算法(GBIS)將非線性反演應用于斷層幾何[8],對于反演,使用結合了Metropolis-Hastings算法的馬爾可夫鏈蒙特卡羅方法(MCMC)估計最佳參數及其不確定性[9-10]。為了提高反演效率和計算速度,使反演結果得到更好的收斂,采用四叉樹方法對D-InSAR數據進行二次采樣[11],設置降軌的采樣閾值為0.024 m,采樣后數據點降為1 266個;升軌的采樣閾值為0.014 m,采樣后數據點降為1 382個。以降軌和升軌數據同時作為約束條件,設置斷層走向的限制范圍為-180°至180°,傾角的限制范圍設置為0°至90°,迭代次數設置為106次,使得所有模型參數都根據其不確定性得到了很好的約束。反演結果顯示深度最佳值為17.35 km,與USGS發出的17 km相近,斷層的傾角最佳值為73.8°與USGS發布的76°相差不大,斷層的走向最佳值為106.9°,與USGS發布的107°基本吻合。震級與USGS也基本吻合,為MW6.64(表3)。

表3 斷層幾何與滑動參數Table 3 Fault geometry and sliding parameters

圖4(a)為計算斷層參數時,模擬的降軌同震形變場,存在5個條紋周期,對應的形變量值為14 cm,與圖2(b)中西部的形變量值一致,而與圖2(a)中東部的形變量值最大差異為6 cm。圖4(b)顯示出東部的殘差值明顯大于西部,其最大殘差值為～5.6 cm。圖4(c)為計算斷層參數時,模擬的升軌同震形變場,也存在5個條紋周期,對應的形變量值為14 cm,圖4(d)顯示出西南部條紋較密集,具有較大殘留信號,殘差值為～14 cm。根據表3斷層參數的模擬LOS方向的形變場結果來看,降軌模式優于升軌模式。

圖4為采用表3中的斷層參數計算的降、升軌模式的形變場模擬與殘差結果,長箭頭表示衛星飛行方向,短箭頭表示視線方向,每個條紋周期代表2.8 cm圖4 同震干涉圖的模擬與殘差Fig.4 Simulation and residual of coseismic interferogram

2.2 單一斷層模型的斷層滑動分布反演

根據表3中的斷層幾何參數結果建立單一斷層模型[12],起點坐標設置為(39.342°N,73.732°E),沿走向和下傾方向分別將斷層擴展到70 km和20 km,將斷層劃分為2 km×2 km的子斷層塊,共計350個。SDM反演斷層面上的滑動分布結果如圖5所示,顯示斷層運動以右旋走滑為主,斷層滑動分布主要集中在長度方向上30～70 km 之間,深度方向上5～20 km的區域。其中斷層面上的兩個滑動峰值對應的滑動量分別為 0.83 m、0.66 m,整個斷層面的平均滑動量為0.31 m,對應的地震矩震級約為MW6.65。

圖5 單一斷層模型的滑動分布Fig.5 Slip distribution of the single fault model

圖6(c)可以看出,單一滑動分布仍有較大的殘差,其殘差最大達到5.03 cm,說明單一斷層模型并不能很好地解釋觀測數據,由此判斷發震斷層可能為兩段斷層共同破裂而成,所以將斷層分為兩段斷層模型再次進行精細滑動分布反演。

黑色線為木吉斷裂圖6 單一斷層模型均勻滑動分布模擬的形變場Fig.6 Deformation field simulated by uniform sliding distribution of the single fault model

2.3 雙斷層模型的斷層滑動分布反演

以降軌LOS向形變結果為約束,參考單斷層反演結果,將斷層一分為二。斷層面Ⅰ起點設置為(39.308°N,73.78°E),長度為40 km,寬度為20 km,斷層面II起點設置為(39.191°N,74.332°E),長度設置為30 km,寬度為20 km,走向仍為106.9°,傾角采用粒子群算法[13],通過大量的反演實驗,確定斷層面Ⅰ最優傾角為70.79°,斷層面Ⅱ最優傾角為55.33°,圖7為雙斷層滑動量分布的反演結果,圖7(a)表示斷層面Ⅰ的結果,斷層面滑動峰值的深度11.61 km,最大值為0.68 m,(b)表示斷層面Ⅱ的結果,斷層面滑動峰值的深度為7.46 km,最大值為0.77 m,對應的地震矩震級約為MW6.62。

圖7 雙斷層模型的滑動分布Fig.7 Slip distribution of the double fault model

對比雙斷層滑動分布的結果計算的形變場[圖8(b)]與dLOS觀測值[圖8(a)]發現兩者整體分布形態具有一致性,數量級相同,形變量在斷層東部存在一定差異,其最大值為3.09 cm[圖8(c)],相比單一斷層模型圖8(c)而言,殘差值減小了1.94 cm。

2.4 庫侖應力

地震的同震位錯會使周圍斷裂或者臨近區域的應力狀態發生改變,從而引起庫倫應力的變化[14],并可能影響臨近斷層狀態,從而改變地震發生的概率[3]。本文在單一、雙斷層模型的基礎上,計算了震源深度17.35 km處的最優破裂面庫倫應力變化(圖9)。圖9中的黑色點表示2016年11月26日0點至11月29日24點發生的余震,可以看出地震的余震基本發生在斷層的南部,與應力的增強區相對吻合;綠色點表示震后自2017年1月1日到2020年10月28日發生的地震,圖9(a)、(b)分別表示為單一斷層和雙斷層的庫倫應力結果。對比可以看出,雙斷層在北部和南部的正庫倫應力變化增加了,且在雙斷層滑動分布的情況下計算的庫倫應力變化更符合震后地震的分布。

藍色區域表示庫倫應力的減小區,紅色區域表示庫倫應力的增強區;深藍色的線表示該區域附近的斷層,其中F1為帕米爾主斷裂帶,F2為木吉斷裂,F3為公格爾拉張系北段;黑色的點表示為余震的分布區域,綠色的點表示從2017年至今發生在該區域的大于3級的地震圖9 單一斷層、雙斷層最優破裂面的庫倫應力Fig.9 Coulomb stress of optimal fracture surface for single fault and double fault

3 討論

本次地震發生在新疆西南緣高寒地區,綜合前人對帕米爾地質構造研究、InSAR數據解算結果,反演阿克陶地震的滑動分布,在認識余震分布的基礎上進行庫倫應力變化分析,結合上述結果,展開如下討論:

雙斷層能否更好地反演斷層滑動分布？根據單一斷層反演得到的同震破裂長度為64.2 km,寬度為19 km,走向為105°,傾角為73.8°,與USGS第Ⅱ節面近似,但反演結果殘差約5 cm。從反演的單一滑動分布結果來看,兩個子破裂中間存在一個較小的滑動,可判斷出該地震存在復雜震源破裂,可能在短時間內發生了兩次破裂。繼而進行雙斷層反演,采用粒子群算法反演出西段斷層長度為48.4 km,寬度為19 km,走向為105.68°,最優傾角為70.79°,深度為5.94～22.94 km,東段破裂長度為14.99 km,寬度為19 km,走向為104.08°,最優傾角為55.36°,深度范圍為5.82～15.69 km,最大擬合殘差值約3 cm。對比擬合殘差結果,可見雙斷層的模型能更好地模擬斷層結構。

庫倫應力變化的計算,重點是創建地震的震源模型。本文使用反演的有限元斷層模型進行計算。輸入350個的小斷層塊的滑動量去計算庫倫應力的變化,后期將余震和后期發生的地震與之對應,可以看出,余震發生在正庫倫應力變化區的概率均相對較高,從2017年至今發生的地震在基于雙斷層反演結果計算的正庫倫應力變化區達到57%,可見,庫倫應力變化為正值時,很有可能會引發后一次地震事件的發生。阿克陶地震前后有多次地震,可見該區域的活動性是相當活躍,推斷該地區未來仍有發生強震的可能。

4 結論

本文通過多次對比實驗,綜合前人的研究結果,對滑動分布再一次進行了研究,給出了該地震的構造模型,得到了以下幾點結論:

(1)采用貝葉斯方法反演了阿克陶地震的單一斷層幾何,得到其最優走向為106.9°,傾角為73.8°。反演滑動分布結果顯示存在兩個滑動峰值,其中西部最大滑動量為0.66 m ,東部最大滑動量為0.83 m,對應的矩震級為MW6.65。

(2)采用粒子群算法計算阿克陶地震同震斷層西段(斷層面Ⅰ)與東段(斷層面Ⅱ)的傾角分別為70.79°、55.33°。反演結果顯示斷層面Ⅰ為走滑斷層,滑動峰值的深度11.61 km,最大滑動量為0.68 m;斷層面Ⅱ是以走滑為主兼正斷性質,滑動峰值的深度為7.46 km,最大滑動量為0.77 m;總的矩震級為MW6.62,整體擬合殘差比單一斷層模型的結果降低了約2 cm。

(3)最優破裂面的庫倫應力變化,其應力增強區與余震分布相對吻合。從震后至今該區域發生的3級以上地震,其發生在應力增強區的概率達到50%。

致謝:感謝中國地震局地球物理研究所房立華研究員提供的余震精定位數據;歐空局提供的Sentinel-1A與Sentinel-1B衛星數據;德國GFZ的汪榮江老師提供的SDM反演程序。部分圖由GMT繪制。感謝兩位評審專家對本文提出的寶貴意見！