石富強(qiáng) 朱琳 王瑩 丁曉光 邱玉榮 石軍

陜西省地震局，西安市碑林區(qū)水文巷4號 710068

0 引言

巴顏喀拉塊體東北緣位于南北地震帶中段，區(qū)域內(nèi)分布有龍門山斷裂、東昆侖斷裂、西秦嶺北緣斷裂等多條大型活動斷裂帶，為現(xiàn)今地殼變形顯著、地震活躍的區(qū)域。僅20世紀(jì)以來，該區(qū)域就發(fā)生了1933年疊溪7.5級、1976年松潘平武7.2級、2008年汶川8.0級、2017年九寨溝7.0級等多次7級以上強(qiáng)震。對該區(qū)域的地震危險性跟蹤一直是地震科技工作者關(guān)注的焦點(diǎn)，M7專項(xiàng)工作組根據(jù)歷史地震活動、大地測量及地球物理場數(shù)據(jù)對該區(qū)域給出了甘青川交界危險區(qū)、西秦嶺中西段危險區(qū)以及六盤山南-西秦嶺東危險區(qū)等強(qiáng)震高風(fēng)險區(qū)域。汶川 8.0級地震后，一些學(xué)者（Xiong et al，2010；Toda et al，2008；Wan et al，2010；單斌等，2009）基于歷史強(qiáng)震演化的斷層庫侖應(yīng)力變化的研究結(jié)果表明，東昆侖斷裂東段以及西秦嶺北緣斷裂西段的地震危險性值得關(guān)注。此外，該區(qū)域還存在歷史地震破裂空段、低b值、小震活動圖像、地震矩加速釋放等異常（陳為濤等，2013；M7專項(xiàng)工作組，2012）以及跨斷層水準(zhǔn)等地球物理場觀測數(shù)據(jù)異常①中國地震局第二監(jiān)測中心，2016，2017年度地震趨勢會商報告。。

庫侖應(yīng)力研究作為研究地震輪回中構(gòu)造應(yīng)力場時空變化的重要組成部分，已成為研究斷層應(yīng)力狀態(tài)及區(qū)域地震危險性的重要手段之一。King等（1994）研究指出，1992年美國Landers 7.3級地震在Big Bear產(chǎn)生的0.2MPa左右的應(yīng)力加載量加速了Big Bear 6.5級地震的發(fā)生。1999年Izm it 7.4級地震發(fā)生前，Stein等（1997）等基于庫侖應(yīng)力計算指出該區(qū)域地震危險性較高。而2004年P(guān)arkfield 6.0級地震的發(fā)生則受1983年Coalinga 6.5級地震及Nu?ez多個中等地震的應(yīng)力卸載的影響，出現(xiàn)了約11年的時間延遲（Bakun et al，1985）。此次九寨溝地震就發(fā)生在M7專項(xiàng)工作組（2012）給出的甘青川交界危險區(qū)內(nèi)，基于庫侖應(yīng)力的強(qiáng)震作用研究表明，區(qū)域內(nèi)1976年松潘平武7.2級地震、2008年汶川8.0級地震等都對此次九寨溝地震發(fā)震斷層有較強(qiáng)的庫侖應(yīng)力加載作用（徐晶等，2017；單斌等，2017）。

有限元數(shù)值模擬技術(shù)作為探索構(gòu)造變形等的虛擬實(shí)驗(yàn)手段，在地震危險性分析中得到了廣泛應(yīng)用（李玉江等，2013；Moreno et al，2010；Hergert et al，2010；陳連旺等，2008；祝愛玉等，2015）。為進(jìn)一步分析九寨溝地震對巴顏喀拉塊體東北緣主要活動斷裂的加卸載效應(yīng)以及該區(qū)域未來可能的強(qiáng)震風(fēng)險，本文構(gòu)建了巴顏喀拉塊體東北緣三維粘彈性有限元模型，模擬計算了九寨溝地震對周邊活動斷裂的加卸載效應(yīng)，分析了震后地震活動及其與九寨溝地震間的可能關(guān)系以及未來可能的強(qiáng)震風(fēng)險。

1 模型構(gòu)建

1.1 幾何模型及物性參數(shù)

本文在詳細(xì)調(diào)研地震地質(zhì)資料的基礎(chǔ)上，通過對所調(diào)研的地質(zhì)資料進(jìn)行假設(shè)構(gòu)建了巴顏喀拉塊體東北緣有限元模型：模型縱向深100km，地殼分層、速度結(jié)構(gòu)以及密度結(jié)構(gòu)根據(jù)Crust1.0速度結(jié)構(gòu)給出（Laske et al，2013）；模型中包含禮縣-羅家堡斷裂、虎牙斷裂、雪山斷裂、東昆侖斷裂、龍門山斷裂、青川斷裂、青川-平武斷裂、臨潭-宕昌斷裂、白龍江斷裂、塔藏斷裂、龍日壩斷裂、岷江斷裂、西秦嶺北緣斷裂以及基于有限斷層反演給出的此次九寨溝地震發(fā)震斷裂等等；區(qū)域孕震層深度（段星北，1997；張國民等，2002），除九寨溝地震發(fā)震斷層外，假定其他斷層面下邊界深度為20km；除九寨溝地震發(fā)震斷層施加位錯荷載外，其余斷層均為摩擦接觸斷層，斷層空間幾何結(jié)構(gòu)簡化為由斷層地表出露沿著傾角向下延伸的空間斜面；下地殼及地幔粘滯性結(jié)構(gòu)根據(jù)邵志剛等（2008）、徐晶等（2017）的研究給出。地殼單元尺寸控制在10km內(nèi)，地幔單元尺寸控制在20km內(nèi)，斷層面上單元尺寸控制在3km內(nèi)，模型共包含359676個單元、67796個節(jié)點(diǎn)，構(gòu)建的有限元模型及對應(yīng)的物性參數(shù)分別見圖1和表1。

1.2 邊界約束及模擬方案

地震是構(gòu)造應(yīng)力作用下斷層應(yīng)力狀態(tài)達(dá)到一定程度而斷層破壞形成的自然災(zāi)害。在數(shù)值模擬中，由于我們并不知道斷層當(dāng)前絕對應(yīng)力狀態(tài)以及斷層自身的承載極限，因此，很難做到對地震自發(fā)的模擬。為此，本文首先對模型施加以GPS約束的邊界荷載來模擬初始構(gòu)造應(yīng)力場；然后，再施加地震同震位錯，模擬計算給出九寨溝地震后區(qū)域的模擬構(gòu)造應(yīng)力場；最后，從震后模擬應(yīng)力場中扣除模擬初始應(yīng)力場，給出九寨溝地震的發(fā)生對區(qū)域主要活動斷裂的加卸載效應(yīng)。其中，九寨溝地震位錯模型參考波形反演結(jié)果（http：//www.cea-igp.ac.cn/tpxw275883.html）給出。所有模擬計算過程均在有限元開源代碼 Pylith（https：//geodynamics.org/cig/software/pylith/）平臺下進(jìn)行。

圖1 本文構(gòu)建的巴顏喀拉塊體東北緣有限元模型

表1 本文有限元模型力學(xué)參數(shù)

1.3 庫侖應(yīng)力計算

三維有限元模擬軟件并不能直接給出斷層面上的庫侖應(yīng)力，因此，需要通過計算給出的應(yīng)力張量來計算斷層面上的正應(yīng)力及剪應(yīng)力，進(jìn)而計算斷層面庫侖應(yīng)力變化。圖2為斷層面主方向（）與模型坐標(biāo)系（）相對關(guān)系，其中，分別為斷層走向（面向走向，右手在上盤）、傾向（逆沖為正）以及張拉方向（張為正）；分別為地理E、N向以及垂直地表向上的方向。α為斷層走向；β為斷層傾角。將斷層面應(yīng)力分量在模型坐標(biāo)系（、）下表示為

其中，T為矩陣轉(zhuǎn)置符號。

圖2 斷層變形模式示意圖（a）以及斷層面主方向與模型坐標(biāo)間的關(guān)系（b）

由圖2（b）可知，將模型坐標(biāo)系下的應(yīng)力張量繞z軸逆時針旋轉(zhuǎn)90°－α，然后，再繞主軸1逆時針旋轉(zhuǎn)β，便可得到主方向下的應(yīng)力張量。則

由式（1）、（2）給出的新的應(yīng)力分量σ33即為斷層面上的正應(yīng)力變化 Δσn，σ23、σ13分別為斷層面上沿傾向和走向的剪應(yīng)力分量。令斷層面優(yōu)勢破裂方向的滑動角為θ，則斷層面優(yōu)勢破裂滑動方向的剪應(yīng)力變化為

根據(jù)斷層面庫侖應(yīng)力變化的定義（King et al，1994），給出斷層面優(yōu)勢破裂方向的庫侖應(yīng)力變化為

其中，μ為斷層面有效摩擦系數(shù)。

2 模擬結(jié)果分析

2.1 模擬初始應(yīng)力場

以1999～2007年GPS觀測速度為約束，通過本文構(gòu)建的有限元動力學(xué)模型計算了巴顏喀拉塊體東北緣的主應(yīng)力場（圖3）。模擬結(jié)果與基于小震震源機(jī)制反演的區(qū)域構(gòu)造應(yīng)力場（王曉山等，2015）以及基于GPS推算的區(qū)域應(yīng)變率場（陳長云等，2013）基本一致。區(qū)域主要處于NW-SE拉張、NE-SW西向擠壓的動力環(huán)境。九寨溝地震震中南側(cè)至龍門山斷裂，主張應(yīng)力逐漸減小，主壓方向逐漸偏轉(zhuǎn)并呈現(xiàn)出近EW向擠壓態(tài)勢。九寨溝地震震中附近主應(yīng)力方向與 USGS給出的震源機(jī)制解（https：//earthquake.usgs.gov/earthquakes/eventpage/us2000a5x1＃moment-tensor）的 P軸（方位角 104°）和 T軸（方位角 204°）基本一致。

2.2 九寨溝地震及其對周邊斷裂的加卸載效應(yīng)

以同震位錯為荷載模擬給出的九寨溝地震震中附近的同震位移場如圖4所示。由圖4可見，震中最大同震位移約為9cm，九寨溝地震的發(fā)生對東昆侖斷裂、塔藏斷裂、西秦嶺北緣斷裂東段以及龍門山斷裂均為左旋作用，同震位移場的左旋作用有利于左旋活動為主的東昆侖斷裂、塔藏斷裂以及西秦嶺北緣斷裂東段斷層面的應(yīng)力積累；而對于具有右旋逆沖活動的龍門山斷裂則為卸載作用。同理，對西秦嶺北緣斷裂西段應(yīng)為卸載影響。為進(jìn)一步定量分析九寨溝地震對各斷裂的加卸載效應(yīng)，我們基于本文模擬結(jié)果，根據(jù)式（2）、（4）計算了九寨溝地震在這些斷裂斷層面上引起的庫侖應(yīng)力變化（圖5）。

圖3 本文模擬給出的巴顏喀拉塊體東北緣初始主應(yīng)力場

圖4 數(shù)值模擬給出的九寨溝地震震中附近同震位移場

由圖5可見，九寨溝地震對塔藏斷裂西段、虎牙斷裂加載效應(yīng)顯著，達(dá)7056Pa；東昆侖斷裂東段、龍日壩斷裂、迭部-白龍江斷裂西端和東端、臨潭-宕昌斷裂東段以及西秦嶺北緣斷裂東段也均受到不同程度的加載作用；塔藏斷裂東段以及岷江斷裂受九寨溝地震的影響，斷層面庫侖應(yīng)力變化表現(xiàn)為卸載，最大卸載量達(dá)73984Pa。位于此次九寨溝地震破裂面正前方的雪山斷裂斷層面庫侖應(yīng)力變化較為復(fù)雜：由圖4模擬同震位移場可知，九寨溝地震對雪山斷裂西段的影響表現(xiàn)為擠壓，對東段的影響表現(xiàn)為拉張，根據(jù)庫侖應(yīng)力計算公式可知，當(dāng)斷層面應(yīng)力變化主要表現(xiàn)為垂直于斷層面的擠壓/拉張時，有效摩擦系數(shù)的取值對庫侖應(yīng)力變化影響顯著，甚至可以影響其極性變化。Wan等（2010）在分析汶川地震對周邊斷裂影響時也給出同樣的結(jié)論。從本文模擬結(jié)果看，當(dāng)有效摩擦系數(shù)較小（μ=0.1）時，九寨溝地震對雪山斷裂西段表現(xiàn)為庫侖應(yīng)力加載，對東段表現(xiàn)為卸載；當(dāng)有效摩擦系數(shù)較大（μ=0.7）時，對西段表現(xiàn)為卸載作用，而對東段表現(xiàn)為加載（圖5）。

此外，此次九寨溝地震還對汶川地震破裂區(qū)東北端的青川斷裂以及青川-平武斷裂的西南段有一定的庫侖應(yīng)力加載，最大加載量可達(dá)1139Pa。邵志剛等（2010）研究汶川地震的震后影響發(fā)現(xiàn)，青川斷裂庫侖應(yīng)力加載非常顯著，實(shí)際上青川地區(qū)也是汶川地震余震活動非常活躍的區(qū)域（石軍，2008），但從2013年開始，該區(qū)域地震活動水平明顯減弱。此次九寨溝地震后不到2個月，9月30日青川（32.27°N，105°E）再次發(fā)生5.4級地震。分析認(rèn)為，此次青川地震可能是由九寨溝地震的觸發(fā)作用所致。

圖5 有限元模擬給出的不同有效摩擦系數(shù)下斷層面庫侖應(yīng)力變化

3 結(jié)論與討論

基于構(gòu)建的三維有限元模型，以同震位錯為荷載模擬了九寨溝地震的發(fā)生及其對巴顏喀拉塊體東北緣主要活動斷裂的加卸載效應(yīng)。數(shù)值模擬結(jié)果顯示，九寨溝地震的發(fā)生對處于甘青川交界危險區(qū)內(nèi)的東昆侖斷裂東段至塔藏斷裂西段庫侖應(yīng)力加載效應(yīng)顯著，最大可達(dá)7056Pa；對六盤山南-西秦嶺東危險區(qū)的西秦嶺北緣斷裂東段也有一定程度的加載效應(yīng)，約為346Pa。此外，計算結(jié)果還顯示，九寨溝地震對龍日壩斷裂、虎牙斷裂、青川-平武斷裂西段、迭部-白龍江斷裂西段和東段、臨潭-宕昌斷裂東段等也均有不同程度的庫侖應(yīng)力加載效應(yīng)；對岷江斷裂、塔藏斷裂東段庫侖應(yīng)力卸載效應(yīng)顯著，最大卸載量可達(dá)73984Pa。

地震是區(qū)域構(gòu)造應(yīng)力積累到一定程度且超越斷層面承載極限而引起斷層剪切破壞的動力過程。但實(shí)際上我們很難獲取斷層面的絕對應(yīng)力狀態(tài)和斷層面的承載極限。庫侖應(yīng)力作為由地震發(fā)生而引起的構(gòu)造應(yīng)力場變化在斷層面優(yōu)勢破裂方向的體現(xiàn)，在地震危險性定性分析判斷中有著廣泛的應(yīng)用。汶川地震后，許多學(xué)者通過庫侖應(yīng)力計算分析了汶川地震對區(qū)域活動斷層的加卸載效應(yīng)及地震危險性（Toda et al，2008；單斌等，2009；邵志剛等，2010；李玉江等，2013）。隨后的2013年蘆山7.0級地震、岷縣漳縣6.6級地震以及此次九寨溝7.0級地震均發(fā)生于文獻(xiàn)報道的庫侖應(yīng)力加載區(qū)域。這些強(qiáng)地震的時空叢集活動也表明，巴顏喀拉地塊東段在未來一段時間內(nèi)可能繼續(xù)保持這種活躍態(tài)勢。從庫侖應(yīng)力變化的角度來看，此次九寨溝地震對甘青川交界危險區(qū)內(nèi)的主要活動斷裂（東昆侖斷裂東段和塔藏斷裂）以及位于六盤山南-西秦嶺東危險區(qū)的主要活動斷裂（西秦嶺北緣斷裂東段），均具有不同程度的庫侖應(yīng)力加載效應(yīng)。