魯?shù)榈卣鸢l(fā)震斷層自發(fā)破裂模擬

廖力 李平恩 劉盼 楊建思

摘要：從彈性動力學(xué)方程出發(fā)模擬了魯?shù)榈卣鹪诎熔瘛『訑嗔训卣鹱园l(fā)破裂過程，探討影響魯?shù)榈卣鹌屏训囊蛩亍Ｑ芯拷Y(jié)果表明：魯?shù)榈卣鸬淖笮呋恼鹪礄C制以及震級主要是受背景應(yīng)力場的影響，斷層滑移分布受到斷層幾何結(jié)構(gòu)、水平應(yīng)力場方向及相對大小的影響，非平面復(fù)雜斷層幾何結(jié)構(gòu)是導(dǎo)致魯?shù)榈卣饛?fù)雜的滑動位移分布的原因。

關(guān)鍵詞：魯?shù)榈卣?自發(fā)破裂;曲線有限差分方法;包谷垴—小河斷裂

中圖分類號：P315.332文獻標(biāo)識碼：A文章編號：1000-0666（2019）01-0024-09

0引言

2014年8月3日云南魯?shù)榘l(fā)生6.5級地震，震源深度約12km。此次地震震中位于包谷垴—小河斷裂及昭通—魯?shù)閿嗔迅浇?7.1°N，103.3°E），是青藏高原巴顏喀拉塊體及其鄰近地區(qū)最新地震活躍期內(nèi)發(fā)生的一次中等強度地震（鄧起東等，2014）。張勇等（2014）利用地震波形數(shù)據(jù)反演了魯?shù)榈卣鹌屏堰^程，結(jié)果顯示是一次共軛破裂地震事件，約70%的地震能量釋放在近SN向的包谷垴—小河斷裂上;主震震源機制表明此次地震為一次高傾角左旋走滑地震。張振國等（2014）依據(jù)張勇等（2014）反演的破裂過程模擬了魯?shù)榈卣鸬膹姷孛孢\動并進行了烈度預(yù)測?，F(xiàn)場地質(zhì)調(diào)查發(fā)現(xiàn)魯?shù)榈卣鸬牡乇砥屏丫诎熔瘛『訑嗔褞Вɡ钗鞯龋?014），顯示此次地震的發(fā)震斷裂為包谷垴—小河斷裂，余震精定位（房立華等，2014）及無人機航測技術(shù)分析（于江等，2018）的結(jié)果也支持這一論斷。因此，本研究擬從彈性動力學(xué)方程出發(fā)模擬魯?shù)榈卣鹪诎熔瘛『訑嗔训牡卣鹱园l(fā)破裂過程，探討分析影響魯?shù)榈卣鸹瑒游灰品植嫉闹饕蛩亍?/p>

1研究方法與模型構(gòu)建

1.1數(shù)值方法

本文利用曲線有限差分方法模擬包谷垴—小河斷裂（圖1）上的魯?shù)榈卣鸬膭恿W(xué)自發(fā)破裂過程。該方法采用分裂節(jié)點描述斷層間斷面，網(wǎng)

格沿著斷層以及自由表面等不規(guī)則界面劃分。通過重新定義曲面斷層分裂節(jié)點物理量（Dayetal，2005;Dalguer，Day，2007），曲線坐標(biāo)系中的彈性動力學(xué)方程可以與摩擦準(zhǔn)則約束的物理變量牽引力結(jié)合起來。該方法既保持了傳統(tǒng)有限差分方法的特點，如計算效率高、容易施加等，又在模擬復(fù)雜斷層方面有較高靈活性，可以用來模擬現(xiàn)實中常見的復(fù)雜幾何斷層，如傾斜、彎曲、跳躍等斷層（張偉，2006;張振國，2014）。

1.2斷層結(jié)構(gòu)

包谷垴—小河斷裂東起昭通—魯?shù)閿嗔妖堫^山段，向西延伸至蓮峰斷裂，平均走向N30°W（李克昌等，1981）。震后地質(zhì)調(diào)查表明，該斷裂斷層產(chǎn)狀為60°∠60°（李西等，2014），斷層長度約40km，圖1給出了包谷垴—小河斷裂的不規(guī)則地表跡線示意圖。

現(xiàn)場地質(zhì)調(diào)查只能給出斷層在地表的出露，對于其深部構(gòu)造則不能給出更多約束信息。斷層的三維復(fù)雜幾何結(jié)構(gòu)一般由觀測資料反演獲得，特別是近場強震記錄、空間觀測技術(shù)等資料。綜合這些因素，本文構(gòu)造基礎(chǔ)斷層三維結(jié)構(gòu)，即斷層在地表的走向參考現(xiàn)場地質(zhì)調(diào)查給出的地表跡線進行約束，深部則按照傾角為90°的斷層向下延拓，整個斷層模型沿走向方向長約60km，沿傾向方向?qū)?0km。

現(xiàn)有的觀測資料表明，斷層兩側(cè)巖體的接觸面并非簡單平面，而是存在不同空間尺度的粗糙表面，數(shù)值計算表明粗糙表面會對自發(fā)破裂以及地震波輻射產(chǎn)生重要影響。因此筆者在基礎(chǔ)三維斷層面結(jié)構(gòu)基礎(chǔ)之上附加自相似擾動，擾動尺度約為1km。

地震，尤其是破裂面積較大的地震，其發(fā)生一般是由于一定區(qū)域達到破裂條件開始破裂并向外傳播。數(shù)值模擬地震的破裂過程也需要一個成核區(qū)，在成核區(qū)內(nèi)人為設(shè)置一個高于斷層強度的剪切應(yīng)力值，使之處于起始破裂狀態(tài)，之后引起斷層面上其他區(qū)域的自發(fā)破裂。根據(jù)地震精定位結(jié)果（房立華等，2014），魯?shù)榈卣鸬钠鹗计屏褏^(qū)域發(fā)生在斷層模型沿走向約20km處，參考Xie等（2015）的震源機制解及張勇等（2014）的震源破裂過程及精定位結(jié)果，震源斷層參數(shù)設(shè)定為斷層走向165°，傾角90°，震源矩心深度3km，地震自發(fā)破裂深度設(shè)定為6km，震源成核區(qū)在龍頭山附近。

1.3應(yīng)力狀態(tài)與動力學(xué)參數(shù)

地質(zhì)觀測資料以及研究區(qū)域的中小地震震源機制反演結(jié)果顯示研究區(qū)域的斷層以走滑斷層為主?；谠颇系貐^(qū)地震震源機制分析，孫業(yè)君等（2017）給出了該區(qū)域內(nèi)三軸構(gòu)造應(yīng)力分布與相對大小。根據(jù)其研究結(jié)果，本文中包谷垴—小河斷裂帶的應(yīng)力比設(shè)定為：

式中：σH，σh分別是最大、最小水平主壓應(yīng)力;σv為垂直主壓應(yīng)力，各應(yīng)力基本關(guān)系表達式為：

式中：ρ和ρw分別表示上覆巖體和水的密度;h表示質(zhì)點到地表的深度;g是重力常數(shù)。斷層及周圍巖體在構(gòu)造應(yīng)力的驅(qū)動下發(fā)生形變及移動，當(dāng)形變積累的應(yīng)力能力超過最大強度時，斷層發(fā)生相對錯動，釋放應(yīng)變能，產(chǎn)生地震波輻射。本文選取大尺度范圍內(nèi)的空間觀測結(jié)果作為數(shù)值模擬中的應(yīng)力場方向，并設(shè)定最大水平主壓應(yīng)力σH的方位角為330°，其他兩個主壓應(yīng)力軸與σH相互垂直，三者構(gòu)成右手系。圖1給出了2個水平應(yīng)力的方向和相對大小，其中箭頭對分別表示最大水平主壓應(yīng)力σH和最小水平主壓應(yīng)力σh，σv垂直于水平面。

本文用簡單速度弱化準(zhǔn)則（Slip-weakeninglaw）作為摩擦準(zhǔn)則（Ida，1972;Andrews，1976），臨界滑移距離Dc=0.4m，在小于臨界滑移距離時，靜摩擦系數(shù)逐漸減小直至臨界滑移距離轉(zhuǎn)變?yōu)閯幽Σ料禂?shù)，也反映了斷層面滑動過程中應(yīng)力釋放過程（圖2）。劉博研和史保平（2011）利用原地應(yīng)力測量數(shù)據(jù)（郭啟良等，2009）和簡單傾斜斷層模型分析汶川地震動力學(xué)參數(shù)，得到靜摩擦系數(shù)為0.5，動摩擦系數(shù)為0.2，由于鉆井深度在400m左右，該研究僅限于淺部地震有一定參考意義。本文設(shè)定靜摩擦系數(shù)μs為0.4，動摩擦系數(shù)μd為0.24。

1.4速度結(jié)構(gòu)

包谷垴—小河斷裂地處青藏高原東南緣（鄧起東等，2002），位于印度板塊與歐亞板塊碰撞帶的東部，速度結(jié)構(gòu)復(fù)雜，對地震破裂過程以及強地面運動產(chǎn)生重要調(diào)制作用。簡單的速度模型不能反映真實地震發(fā)生時對近場造成的地面運動和地震災(zāi)害。本文采用的是區(qū)域三維速度結(jié)果（鄭晨等，2016），對研究區(qū)域介質(zhì)可以提供較好的約束。

2研究結(jié)果

包谷垴—小河斷裂上魯?shù)榈卣鸬哪M主要是非平面上斷層的動力學(xué)模擬，即地震觸發(fā)之后的傳播、擴展到愈合等一系列的自發(fā)破裂過程。為了精細(xì)地反映出地震的破裂過程，動力學(xué)模擬采用60m網(wǎng)格。

圖3模擬了魯?shù)榈卣鸷?5s內(nèi)的自發(fā)破裂過程，本文截取的是從破裂開始到擴展，直至破裂接觸地表6個比較典型的位移分布結(jié)果。從圖中可以看出，破裂在斷層上的傳播具有較強的非對稱性，破裂從成核區(qū)（五角星）開始向地表傳播，在地表形成了破裂。圖4顯示地震矩能量釋放主要集中在第1～10s，矩震級MW=6.56。自發(fā)破裂模擬的震源時間函數(shù)和滑動位移分布與張勇等（2014）利用近震、遠震數(shù)據(jù)聯(lián)合反演的結(jié)果也較為接近。

從圖5可以看出，魯?shù)榈卣鹱园l(fā)破裂向深部傳播過程中，由于受到趨于靜水圍壓的三軸應(yīng)力的作用，斷層深部的有效剪切應(yīng)力趨于零，滑動在傳播到斷層深部邊界之前就逐漸停止而愈合。破裂在斷層走向兩側(cè)是正常傳播，地震的最終滑移量分布（圖5）顯示較強的非均勻性。最主要的是滑動位移從成核區(qū)向NW方向傳播直至地表，顯示出較強的左旋走滑性質(zhì)，這與區(qū)域背景應(yīng)力場的性質(zhì)關(guān)系較為密切。

3討論

3.1斷層幾何彎曲的影響

本文采用的斷層信息是從房立華等（2014）研究圖中取點，以斷層N30°W（李克昌等，1981）插值得到，斷層傾角以60°從地表向下拓展。因此，模型中的斷層面幾乎是一個平面，而實際斷層面常為不規(guī)則幾何彎曲結(jié)構(gòu)，從而導(dǎo)致斷層面上更為復(fù)雜的應(yīng)力狀態(tài)。為研究斷層幾何彎曲對自發(fā)破裂的影響，將模型中斷層向ES方向延拓50km，按照斷層的平均走向在包谷垴—小河斷裂段設(shè)置一些幾何彎曲，自發(fā)破裂的成核區(qū)設(shè)置在整個斷層的中部（龍頭山段），保持震源深度、模型的背景應(yīng)力場等參數(shù)設(shè)置不變模擬魯?shù)榈卣鹱园l(fā)破裂過程。從圖6結(jié)果可以看到，魯?shù)榈卣鹱园l(fā)破裂的方式和設(shè)置模型較為一致，都是從成核區(qū)向NW方向擴展，但破裂釋放的能量相對小一點，震級達到約MW6.35。相比于前述模型能量釋放比較均勻，在此模型中破裂能量釋放更快更集中，破裂的時間主要集中于0～5s，這種自發(fā)破裂的特點可能與成核區(qū)的斷層幾何彎曲結(jié)構(gòu)有關(guān)。圖7則顯示了魯?shù)榈卣鹱园l(fā)破裂在斷層走向和傾向上的最終滑移量，可以看到，斷層的最終滑移都是發(fā)生在走向上，這也顯示了魯?shù)榈卣鹱笮呋奶攸c。

3.2背景應(yīng)力場的影響

3.2.1最大主壓應(yīng)力方向的影響

孫業(yè)君等（2017）給出了該區(qū)域內(nèi)最大主壓應(yīng)力的方向，但因其采用近于水平主應(yīng)力中量值較大的主應(yīng)力方向近似，因此存在最大可能數(shù)十度的誤差。因此，本文模擬了將水平最大主壓應(yīng)力順時針和逆時針分別偏轉(zhuǎn)20°，即最大主壓應(yīng)力方向分別為310°，350°的情形，模擬結(jié)果如圖8所示。從模擬結(jié)果可以發(fā)現(xiàn)在最大主壓應(yīng)力向逆時針旋轉(zhuǎn)（從310°到350°）的過程中，自發(fā)破裂強度逐漸變大，地震震級逐漸變大，在地表的滑動位移也更大。這表明本文模型中逆時針旋轉(zhuǎn)的構(gòu)造應(yīng)力更有利于斷層上的整體破裂。

3.2.2三軸應(yīng)力場大小相對值的影響

魯?shù)榈卣鹫鹬袇^(qū)域位于以巧家為中心的三軸應(yīng)力場大小相對值R從高向低變化的過渡帶（孫業(yè)君等，2017），為研究三軸應(yīng)力場大小相對值對魯?shù)榈卣鹱园l(fā)破裂的影響，本文還模擬了R值設(shè)定為0.8，0.6和0.4的滑動分布。R值越大，在走滑斷層背景應(yīng)力場表示水平向2個主應(yīng)力大小越接近。從圖9可以看出：當(dāng)R=0.8時，無論在地表還是斷層深部，滑動位移相比其他兩種情況的模擬結(jié)果都更大，R=0.6時的滑動位移比R=0.4時稍大，震級也稍大。以上說明水平向2個主應(yīng)力大小越接近時，地震破裂釋放的能量越大。

3.2.3背景應(yīng)力場的影響

將模型設(shè)定的背景應(yīng)力場按以下方式進行設(shè)置，應(yīng)力比為：

各應(yīng)力基本關(guān)系表達式為：

在此應(yīng)力場設(shè)置下，斷層性質(zhì)表現(xiàn)為正斷層，其他參數(shù)設(shè)置保持不變，模擬魯?shù)榈卣鸷?～8s的自發(fā)破裂過程，斷層傾向滑動速率結(jié)果如圖10所示。從圖中可以看出，魯?shù)榈卣鹱园l(fā)破裂從震源成核區(qū)出發(fā)，向震源兩側(cè)擴展傳播，直到遇到人工強邊界而停止。因為斷層幾何結(jié)構(gòu)比較均勻，因此破裂顯示出比較明顯的均勻分布特性。斷層破裂在此應(yīng)力場設(shè)置下與前述2個模型的破裂方式完全不同，破裂主要是在斷層傾向上發(fā)生，顯示出比較明顯的正斷層性質(zhì)，斷層的最終滑移分布見圖11。由此可見背景應(yīng)力場是影響自發(fā)破裂方式的重要因素。

3.3滑動摩擦系數(shù)的影響

將模型的背景應(yīng)力場按3.2.3節(jié)的方式進行設(shè)置，但將動摩擦系數(shù)在斷層走向80km以南區(qū)域設(shè)為0.24，以北區(qū)域設(shè)為0.27，研究動摩擦系數(shù)對滑動位移分布的影響。由模型采用的滑動弱化準(zhǔn)則可知，高的動摩擦系數(shù)使斷層面上的剪應(yīng)力更高，使應(yīng)力降更小，從而使最終斷層滑移更小。將成核區(qū)設(shè)在沿斷層走向90km處，模擬結(jié)果（圖12）與3.2節(jié)成核區(qū)2側(cè)均勻的位移分布不同，成核區(qū)以南的斷層滑移更大，而成核區(qū)以北因為較大的動摩擦系數(shù)導(dǎo)致較小的應(yīng)力降，使最終滑動位移較小。

4結(jié)論

本文以魯?shù)榈卣馂檠芯繉ο?，考慮區(qū)域三維速度結(jié)構(gòu)差異、區(qū)域背景應(yīng)力場、地震自發(fā)破裂斷裂帶幾何結(jié)構(gòu)、建立能反映地表起伏的斷層自發(fā)破裂有限差分模型。再以現(xiàn)場地質(zhì)調(diào)查的滑動破裂結(jié)果為約束條件模擬魯?shù)榈卣鹱园l(fā)破裂過程，得到以下結(jié)論：

（1）本文模擬的魯?shù)榈卣鹱园l(fā)破裂過程是基于云南區(qū)域的背景應(yīng)力場，魯?shù)閿鄬拥膸缀钨Y料以及最新三維速度結(jié)構(gòu)建立模型進行模擬，以最新震源機制反演結(jié)果的矩心深度以及震源參數(shù)設(shè)定自發(fā)破裂的成核區(qū)，模擬的震源時間函數(shù)與張勇等（2014）的反演結(jié)果較為接近，自發(fā)破裂模擬的地表破裂約6km，與魯?shù)榈卣鸷蟮刭|(zhì)調(diào)查發(fā)現(xiàn)的地表破裂也較為吻合。

（2）研究區(qū)域的背景應(yīng)力場三軸主應(yīng)力相對大小決定了魯?shù)榈卣鸢l(fā)震震級、走滑型的震源機制，主應(yīng)力的方向變化會影響滑動位移分布。

（3）背景應(yīng)力場對地震自發(fā)破裂起決定性作用，而不同的應(yīng)力場參數(shù)包括水平主應(yīng)力方向的改變以及相對大小的改變，會影響地震破裂的傳播、滑動位錯的分布以及地震波的輻射。

（4）斷層的幾何結(jié)構(gòu)會改變斷層面上應(yīng)力場的大小，從而控制斷層自發(fā)破裂，特別是發(fā)震斷層的幾何結(jié)構(gòu)、參數(shù)等存在較強的非均勻性時，不同的起始破裂區(qū)域往往能導(dǎo)致相差較大的地震。地震起始破裂區(qū)域的位置對地震破裂形態(tài)有重要作用，可能對相應(yīng)強地面運動規(guī)律與地震災(zāi)害分布也有重要影響。

（5）在滑動弱化準(zhǔn)則控制下的動摩擦系數(shù)會影響斷層破裂引起的滑動位移分布，動摩擦系數(shù)越小的斷層單元最終滑動位錯越大。

感謝南方科技大學(xué)張振國助理教授、中國科技大學(xué)張文強博士對自發(fā)破裂代碼修改的建議和意見。

參考文獻：

鄧起東，程紹平，馬冀，等.2014.青藏高原地震活動特征及當(dāng)前地震活動形勢[J].地球物理學(xué)報，57（7）：2025-2042.

鄧起東，張培震，冉勇康，等.2002.中國活動構(gòu)造基本特征[J].中國科學(xué)，32（12）：1020-1030.

房立華，吳建平，王未來，等.2014.云南魯?shù)镸S6.5地震余震重定位及其發(fā)震構(gòu)造[J].地震地質(zhì)，36（4）：1173-1185.

郭啟良，王成虎，馬洪生，等.2009.汶川MS8.0大震前后的水壓致裂原地應(yīng)力測量[J].地球物理學(xué)報，52（5）：1395-1401.

李克昌，候?qū)W英，趙維城，等.1981.滇東北地區(qū)地震地質(zhì)特征[J].地震研究，4（1）：57-63.

李西，張建國，謝英情，等.2014.魯?shù)镸S6.5地震地表破壞及其與構(gòu)造的關(guān)系[J].地震地質(zhì)，36（4）：1280-1291.

劉博研，史保平.2011.MS8.0汶川地震斷層的應(yīng)力狀態(tài)以及對余震危險性的影響[J].地球物理學(xué)報，54（4）：1002-1009.

孫業(yè)君，趙小艷，黃耘，等.2017.云南地區(qū)震源機制及應(yīng)力場特征[J].地震地質(zhì)，39（2）：390-407.

于江，張彥琪，李西，等.2018.無人機航測技術(shù)在2014年魯?shù)镸S6.5地震震區(qū)活動構(gòu)造調(diào)查中的應(yīng)用[J].地震研究，41（2）：166-172.

張偉.2006.含起伏地形的三維非均勻介質(zhì)中地震波傳播的有限差分算法以及在強地面運動模擬中的應(yīng)用[D].北京：北京大學(xué).

張勇，許力生，陳運泰，等.2014.2014年8月3日云南魯?shù)镸W6.1（MS6.5）地震破裂過程[J].地球物理學(xué)報，57（9）：3052-3059.

張振國，孫耀充，徐建寬，等.2014.2014年8月3日云南魯?shù)榈卣饛姷孛孢\動初步模擬及烈度預(yù)測[J].地球物理學(xué)報，57（9）：3038-3041.

張振國.2014.三維非平面斷層破裂動力學(xué)研究[D].合肥：中國科學(xué)技術(shù)大學(xué).

鄭晨，丁志峰，宋曉東.2016.利用面波頻散與接收函數(shù)聯(lián)合反演青藏高原東南緣地殼上地幔速度結(jié)構(gòu)[J].地球物理學(xué)報，59（9）：3223-3236.

AbstractInthisstudy，wesimulatetheprocessofspontaneousruptureoftheLudianearthquakeatBaogunaoXiaohefault，anddiscusstheinfluencingfactorsoftheLudianearthquakerupture.TheresultsshowthatthefocalmechanismoftheleftlateralstrikeslipandthemagnitudeoftheLudianearthquakearemainlyaffectedbythebackgroundstressfield.Thedistributionoffaultslipisaffectedbythechangeoffaultgeometry，directionandrelativesizeofhorizontalstressfield.NonplanarcomplexfaultgeometrystructurecausesthecomplexslidingdisplacementdistributionoftheLudianearthquake.

Keywords：Ludianearthquake;spontaneousrupture;curvefinitedifferencemethod;BaogunaoXiaohefault