袁杰， 崔澤飛， 朱守彪*， 王進(jìn)廷

1 應(yīng)急管理部國(guó)家自然災(zāi)害防治研究院， 北京 100085 2 清華大學(xué)水利水電工程系， 北京 100084

0 引言

定量地模擬地震、孕育發(fā)生及其復(fù)發(fā)循環(huán)過(guò)程對(duì)于深入認(rèn)識(shí)地震及地震預(yù)測(cè)預(yù)報(bào)等具有非常重要的科學(xué)意義.但是，地震從孕育到發(fā)生是個(gè)十分復(fù)雜的力學(xué)過(guò)程.通常，強(qiáng)震孕育需要數(shù)十年、上百年甚至幾千年的時(shí)間(朱守彪和張培震，2009；Zhu and Zhang, 2010, 2013)，而地震的發(fā)生則只需數(shù)十秒；此外，由于地震孕育過(guò)程緩慢、變形速率小，該過(guò)程可視為準(zhǔn)靜態(tài)過(guò)程，可用靜態(tài)力學(xué)方程來(lái)求解；但地震發(fā)生時(shí)，斷層產(chǎn)生錯(cuò)動(dòng)，發(fā)射地震波，地面產(chǎn)生強(qiáng)烈震動(dòng)，只能用動(dòng)力學(xué)過(guò)程來(lái)描述.在地震、孕育發(fā)生的模擬過(guò)程中，不僅涉及力學(xué)狀態(tài)的改變，同時(shí)，如何來(lái)選擇計(jì)算中的時(shí)間步長(zhǎng)，也面臨著更大的困難：若按照漫長(zhǎng)的地震孕育過(guò)程來(lái)選擇時(shí)間步長(zhǎng)(如：步長(zhǎng)為半年或一個(gè)月)，則無(wú)法模擬短暫的地震發(fā)生過(guò)程(數(shù)十秒的范圍)；若以地震發(fā)生過(guò)程來(lái)選擇時(shí)間步長(zhǎng)(如：0.0001 s)，則無(wú)法模擬漫長(zhǎng)的地震孕育過(guò)程(數(shù)千年).因此，模擬過(guò)程中，還必須改變計(jì)算的時(shí)間步長(zhǎng).

為克服上述模擬中的困難，不少研究人員進(jìn)行了許多探索：利用靜態(tài)計(jì)算方法模擬斷層的孕震過(guò)程，利用動(dòng)態(tài)計(jì)算方法繼續(xù)模擬斷層的同震破裂過(guò)程.在斷層處于臨界狀態(tài)時(shí)，將靜態(tài)過(guò)程的計(jì)算結(jié)果作為初始條件，并改變計(jì)算方法(力學(xué)狀態(tài))從而繼續(xù)計(jì)算，進(jìn)而連續(xù)模擬孕震-同震過(guò)程 (Duan and Oglesby 2005a,b).然而，人為改變計(jì)算過(guò)程中的力學(xué)狀態(tài)會(huì)將應(yīng)力、位移、摩擦本構(gòu)等之間的復(fù)雜耦合關(guān)系解耦，進(jìn)而會(huì)干擾斷層失穩(wěn)演化的自然力學(xué)過(guò)程；同時(shí)，此方法也無(wú)法精確的判斷斷層的失穩(wěn)臨界點(diǎn).針對(duì)這種情況，Lapusta等(2000)在前人的研究之上(Tse and Rice, 1986; Rice and Ben-Zion, 1996; Ben-Zion and Rice, 1997),基于邊界積分方程方法，開發(fā)了一種可連續(xù)模擬斷層的孕震-同震及交替循環(huán)的計(jì)算方法.隨后又將此方法的應(yīng)用范圍由二維擴(kuò)展到三維 (Lapusta and Liu, 2009).特別值得指出的是，利用該方法， Barbot等(2012)模擬了San Andreas 斷層的Parkfield段的強(qiáng)震孕育發(fā)生過(guò)程，模擬結(jié)果再現(xiàn)了Parkfield斷層段上1966年至2004年間發(fā)生的MW6.0以上的所有地震事件，并且模擬給出的每次地震事件的發(fā)生時(shí)刻、成核位置、震級(jí)大小都與實(shí)際記錄的觀測(cè)資料很好的吻合.可以推測(cè)，若利用該方法繼續(xù)計(jì)算，則可以預(yù)測(cè)出San Andreas斷層的Parkfield段下一次大地震發(fā)生的時(shí)、空、強(qiáng)三要素；從而實(shí)現(xiàn)數(shù)值地震預(yù)測(cè)預(yù)報(bào).

然而，邊界積分方程方法使問(wèn)題的維數(shù)降低一維，例如：三維問(wèn)題降為二維問(wèn)題，二維降為一維問(wèn)題，使得求解的自由度下降.同時(shí)，邊界元僅適應(yīng)規(guī)則區(qū)域及邊界條件的模型，適應(yīng)于求解線性、均勻材料屬性的計(jì)算問(wèn)題.這樣，在處理實(shí)際的地質(zhì)問(wèn)題時(shí)(如：復(fù)雜的幾何邊界、復(fù)雜的物性結(jié)構(gòu)等)會(huì)遇到很大的困難.盡管對(duì)邊界積分方程方法進(jìn)行了不斷地優(yōu)化(如：Sathiakumar et al., 2020)，但目前進(jìn)展不大.因此，對(duì)于復(fù)雜情況的實(shí)際問(wèn)題，模擬斷層的孕震-同震及其循環(huán)過(guò)程，最好還是選擇能夠適應(yīng)復(fù)雜幾何形狀、復(fù)雜邊界條件、復(fù)雜的材料屬性的有限單元方法.雖然很多研究人員(如：Duan and Oglesby, 2005a, b; 朱守彪等，2008；朱守彪和張培震, 2009; Zhu and Zhang, 2010, 2013; Zhu, 2013; 馬林飛等，2018)利用有限單元法模擬了斷層的發(fā)震周期及斷層破裂過(guò)程，但沒(méi)有很好地連續(xù)模擬斷層孕震-同震(準(zhǔn)靜態(tài)-動(dòng)態(tài))及其循環(huán)的全過(guò)程.

1 計(jì)算原理及有限元模型

1.1 動(dòng)力學(xué)方程

利用有限元方法研究斷層自發(fā)破裂動(dòng)力學(xué)過(guò)程問(wèn)題，即整個(gè)區(qū)域內(nèi)求解二維波動(dòng)方程(王勖成，2003)：

σij,j+fi-ρui,tt-cui,t=0.

(1)

(1)式又稱為平衡方程.其中σij是應(yīng)力張量，fi是體力，ρ是彈性介質(zhì)的密度，c是阻尼系數(shù)，ui表示位移，ui,t和ui,tt分別是ui對(duì)t的一次導(dǎo)數(shù)和二次導(dǎo)數(shù)，即分別是速度和加速度；-ρui,tt和-cui,t分別代表慣性力和阻尼力.具體的動(dòng)力學(xué)基本方程和有限元方法見我們之前的文章(袁杰和朱守彪，2014)，以及相關(guān)有限元理論文獻(xiàn)(王勖成，2003)，此處不再贅述.

1.2 摩擦本構(gòu)關(guān)系

斷層上的摩擦關(guān)系非常復(fù)雜，研究人員通過(guò)擬合實(shí)驗(yàn)室實(shí)驗(yàn)數(shù)據(jù)，參照斷層破裂運(yùn)動(dòng)學(xué)反演結(jié)果，提出了數(shù)十種斷層面上的本構(gòu)關(guān)系假設(shè)(Bizzarri,2011).其中，運(yùn)用較為廣泛的是速率-狀態(tài)相關(guān)摩擦的關(guān)系和滑移弱化摩擦的關(guān)系.雖然利用這些摩擦本構(gòu)關(guān)系模擬斷層的破裂傳播過(guò)程，其模擬結(jié)果也能很好解釋一些地震現(xiàn)象，但使用不同的摩擦本構(gòu)關(guān)系，其模擬結(jié)果也不盡相同(Ryan, 2012; Ryan and Oglesby, 2014; Yuan et al., 2020)，到底哪一種摩擦本構(gòu)關(guān)系最貼合斷層的實(shí)際狀況，目前還沒(méi)有定論.本文的主要研究目的是模擬斷層孕震-同震準(zhǔn)周期循環(huán)過(guò)程，因此本文中所有模型擬使用形式簡(jiǎn)單且易于實(shí)現(xiàn)、得到非常廣泛的應(yīng)用的滑移弱化的摩擦本構(gòu)關(guān)系.

在經(jīng)典的滑移弱化摩擦本構(gòu)關(guān)系中，當(dāng)斷層兩側(cè)的相對(duì)滑移距離達(dá)到特征滑移距離后，斷層上的摩擦系數(shù)將一直不變，維持動(dòng)摩擦系數(shù)狀態(tài)，導(dǎo)致一次地震事件之后，斷層無(wú)法繼續(xù)積累能量，也就無(wú)法產(chǎn)生下一次地震(Duan and Oglesby, 2005a).為了能連續(xù)模擬出斷層孕震-同震周期循環(huán)過(guò)程，我們采用改進(jìn)后的滑移弱化摩擦本構(gòu)關(guān)系(Olsen-Kettle et al., 2008; 袁杰和朱守彪,2014; Yuan et al., 2020):

(2)

利用改進(jìn)后的滑移弱化摩擦本構(gòu)關(guān)系模擬斷層破裂過(guò)程時(shí)，在一次地震事件完成之后，斷層的摩擦系數(shù)將會(huì)再次強(qiáng)化，由動(dòng)摩擦系數(shù)變?yōu)殪o摩擦系數(shù).具體在計(jì)算中實(shí)現(xiàn)的做法是：當(dāng)斷層面上對(duì)應(yīng)兩節(jié)點(diǎn)之間相對(duì)滑移速度降為零，即停止滑移時(shí)，將摩擦系數(shù)變?yōu)殪o摩擦系數(shù).為保證計(jì)算過(guò)程的穩(wěn)定、自洽，本文中的摩擦系數(shù)再次強(qiáng)化是一個(gè)緩慢平滑的過(guò)程，而不是直接由動(dòng)摩擦系數(shù)重置為靜摩擦系數(shù).

1.3 動(dòng)力學(xué)方程求解方法

顯式算法和隱式算法是利用有限單元法求解動(dòng)力學(xué)問(wèn)題的兩種主要算法，其主要對(duì)應(yīng)的中心差分法和Newmark法，相關(guān)的理論介紹見有限元理論文獻(xiàn)(王勖成，2003).表1簡(jiǎn)略介紹了兩種算法的優(yōu)缺點(diǎn).

模型中，i和t分別表示不同省份和年份；C1、C2、C4分別表示投資、技術(shù)創(chuàng)新和人口對(duì)產(chǎn)業(yè)結(jié)構(gòu)升級(jí)的作用效應(yīng)；C3表示非老齡人口比例對(duì)產(chǎn)業(yè)結(jié)構(gòu)升級(jí)的影響，當(dāng)C3＞0時(shí)，表示老齡化程度越高，越不利于產(chǎn)業(yè)結(jié)構(gòu)升級(jí)；C5表示分配給勞動(dòng)人口的資源比例對(duì)產(chǎn)業(yè)結(jié)構(gòu)升級(jí)的影響，當(dāng)C5＜0時(shí)，表示分配給老齡人口的資源越多，越能促進(jìn)產(chǎn)業(yè)結(jié)構(gòu)升級(jí)。

表1 顯式算法和隱式算法的優(yōu)缺點(diǎn)Table 1 Advantages and disadvantages of explicit and implicit algorithms

兩種算法都可以用來(lái)模擬斷層的自發(fā)破裂過(guò)程(同震過(guò)程)，但是目前研究中主要采用的是顯式算法.這是由于斷層的瞬時(shí)破裂過(guò)程是高度非線性的問(wèn)題，采用隱式算法是需要很小的時(shí)間步長(zhǎng)能收斂，這時(shí)就體現(xiàn)不出隱式算法的優(yōu)點(diǎn)(時(shí)間步長(zhǎng)Δt不影響解的穩(wěn)定性).實(shí)際對(duì)比兩種算法模擬同震過(guò)程的結(jié)果時(shí)，我們發(fā)現(xiàn)：只有當(dāng)隱式算法的時(shí)間步長(zhǎng)小于0.001 s時(shí)，其計(jì)算精度才能達(dá)到顯式算法模擬結(jié)果.然而時(shí)間步長(zhǎng)同樣采用0.001 s時(shí)，隱式算法計(jì)算同震過(guò)程的求解時(shí)間是顯式算法的數(shù)十倍，所以單純模擬同震過(guò)程時(shí)，主要還是采用顯式算法.

但是在模擬斷層孕震過(guò)程時(shí)，由于顯式算法需要滿足穩(wěn)定條件，即時(shí)間步長(zhǎng)Δt必須小于Δtcr，基本上不可能模擬幾十年甚至上千年的孕震過(guò)程.因此本文基于Newmark隱式時(shí)間積分法，發(fā)展出一種時(shí)間步長(zhǎng)可以隨著斷層狀態(tài)自動(dòng)平滑縮放的計(jì)算方法，進(jìn)而模擬同時(shí)包含孕震-同震兩種狀態(tài)的地震復(fù)發(fā)模型.

1.4 模型幾何及有限元網(wǎng)格剖分

為簡(jiǎn)單起見，本文將實(shí)際逆沖斷層的三維情況抽象為二維模型，并將斷層簡(jiǎn)化為直線.圖1為研究所用的模型幾何及邊界條件：模型空間尺度為100 km×40 km的矩形，斷層的傾向長(zhǎng)度為40 km，傾角為30°.模型底部邊界條件是水平方向自由運(yùn)動(dòng)，垂直方向固定；右側(cè)邊界條件是：垂直方向上自由運(yùn)動(dòng)，水平方向上固定；左側(cè)邊界條件為施加15 mm·a-1的水平構(gòu)造載荷；于自由地表至地下2 km處施加重力加速度.此外，為防止地震波通過(guò)邊界反射而影響結(jié)果，模型在四周設(shè)置了吸收邊界.

圖1 模型幾何及邊界條件Fig.1 Model geometry and boundary conditions

模型全部采用三角形單元來(lái)剖分.斷層附近為研究的重點(diǎn)區(qū)域，為了保證精度，在該區(qū)域?qū)W(wǎng)格進(jìn)行細(xì)化，單元邊長(zhǎng)為100 m，此外，離斷層越遠(yuǎn)地方的單元尺度越大，模型最外圍部分單元的邊長(zhǎng)為500 m.模型中單元的節(jié)點(diǎn)數(shù)和總數(shù)分別為99714和196955.模型中的材料參數(shù)、摩擦系數(shù)等見表2.

表2 模型參數(shù)Table 2 Model parameters

2 模擬結(jié)果

模擬孕震-同震準(zhǔn)周期循環(huán)過(guò)程最大的難點(diǎn)是精準(zhǔn)地找到孕震-同震之間的臨界點(diǎn)，即斷層的臨界失穩(wěn)狀態(tài).此外，模擬程序中時(shí)間步長(zhǎng)如何隨著斷層的滑動(dòng)狀態(tài)做出相應(yīng)的自動(dòng)調(diào)整也是此研究中的難點(diǎn)之一.本文利用新的有限單元計(jì)算方法可以很好的解決上述難點(diǎn)，進(jìn)而連續(xù)穩(wěn)定的模擬斷層孕震-同震及其循環(huán)過(guò)程，模擬結(jié)果如下.

圖2展示了斷層上盤上的典型點(diǎn)(中間點(diǎn))的位移隨時(shí)間的變化曲線.隨著構(gòu)造載荷的不斷擠壓，孕震期間典型點(diǎn)的位移隨時(shí)間呈線性增加，在整個(gè)計(jì)算的2000年時(shí)間里，該點(diǎn)的位移發(fā)生了5次跳躍式的突然增加.位移的每次跳躍式增加，表示發(fā)生了一次地震事件，兩次地震事件之間的時(shí)間間隔即為孕震時(shí)間.圖3是斷層典型點(diǎn)位置上下盤之間的位錯(cuò)隨時(shí)間的變化曲線.可以看出：(1)在時(shí)間較長(zhǎng)的孕震期間斷層兩盤對(duì)應(yīng)點(diǎn)之間位錯(cuò)變化緩慢，而在時(shí)間很短同震期間位錯(cuò)突然增加.(2)孕震時(shí)間與同震位錯(cuò)近似呈正比關(guān)系，當(dāng)孕震時(shí)間較長(zhǎng)時(shí)，緊接著到來(lái)的地震也會(huì)相對(duì)較大.這與前人(馬林飛等,2018)利用有限單元法模擬斷層準(zhǔn)周期行為的結(jié)果類似.但模擬方法不同，本文采用方法中的時(shí)間步長(zhǎng)可以自動(dòng)隨著斷層狀態(tài)做出相應(yīng)的調(diào)整，還可以模擬斷層的瞬時(shí)破裂過(guò)程.

圖2 斷層上盤上的典型點(diǎn)(中間點(diǎn))的位移隨時(shí)間的變化曲線Fig.2 Displacements change with time at the typical point on the hanging wall of the fault

圖3 斷層典型點(diǎn)位置上下盤之間的位錯(cuò)隨時(shí)間的變化曲線Fig.3 Sliding distance at the typical point on the fault vary with time

實(shí)現(xiàn)時(shí)間步自動(dòng)調(diào)整的思路是：(1)用較大的時(shí)間步長(zhǎng)計(jì)算斷層孕震過(guò)程.當(dāng)在某一計(jì)算時(shí)間步內(nèi)(以時(shí)間步長(zhǎng)107s為例)，斷層兩盤的相對(duì)滑移速率發(fā)生突變時(shí)，即認(rèn)為在這一時(shí)間步內(nèi)發(fā)生了地震事件.但由于時(shí)間步長(zhǎng)較大，無(wú)法捕捉同震過(guò)程，因此放棄這一步的計(jì)算結(jié)果，并將時(shí)間步長(zhǎng)減小(以時(shí)間步長(zhǎng)106s為例)，再進(jìn)行計(jì)算.此時(shí)，若斷層兩盤的相對(duì)滑移速率未發(fā)生突變，則認(rèn)為是孕震過(guò)程，繼續(xù)計(jì)算；若發(fā)生了突變，繼續(xù)減小時(shí)間步長(zhǎng).以此類推，直至?xí)r間步長(zhǎng)降到10-3s以下，進(jìn)而計(jì)算斷層的同震過(guò)程.(2)用較小的時(shí)間步長(zhǎng)計(jì)算斷層同震過(guò)程.當(dāng)在較長(zhǎng)一段時(shí)間內(nèi)(以時(shí)間步長(zhǎng)10-3s為例)，斷層兩盤的相對(duì)滑移速率約為零，即認(rèn)為同震過(guò)程結(jié)束，此時(shí)逐漸放大時(shí)間步長(zhǎng)，直至?xí)r間步長(zhǎng)增至107s，進(jìn)而計(jì)算斷層的孕震過(guò)程.圖4展示的是地震事件一和地震事件二之間，孕震過(guò)程和同震過(guò)程的時(shí)間步長(zhǎng)變化.

圖4 孕震期間和同震期間的自動(dòng)增量步時(shí)間的變化Fig.4 The change of the automatic time increments during interseismic and coseismic periods

為了更好的展示斷層兩側(cè)相對(duì)滑移距離隨時(shí)間的變化，圖5a中忽略孕震期間斷層兩盤位錯(cuò)變化，只展示2000年時(shí)間內(nèi)五次地震事件的斷層兩盤位錯(cuò)隨時(shí)間的變化，圖中相鄰曲線的時(shí)間間隔是0.5 s，五次地震事件共造成的最大累計(jì)位錯(cuò)約為7 m.此外，圖5b—f分別展示了每次地震事件斷層兩側(cè)相對(duì)滑移距離隨時(shí)間的變化，每次地震的成核位置、破裂持續(xù)時(shí)間、最大位錯(cuò)等都是不同的，這說(shuō)明這些地震事件不是簡(jiǎn)單的重復(fù).結(jié)合圖3和圖5可以看出，同一條斷層，當(dāng)孕震時(shí)間越長(zhǎng)，其下一次地震事件造成的最大位錯(cuò)就會(huì)越大，即震級(jí)越大.圖6給出地震事件五中斷層在破裂過(guò)程中產(chǎn)生的地震波在介質(zhì)中傳播時(shí)，質(zhì)點(diǎn)振動(dòng)速度云圖在不同時(shí)刻的快照.破裂在斷層上成核后，向地表方向傳播，約6 s后傳至地表.圖5和圖6都表明本文所使用的方法可以模擬斷層詳細(xì)的同震破裂瞬時(shí)過(guò)程.

從上述模擬結(jié)果看，本文所使用的方法可以連續(xù)模擬斷層孕震-同震準(zhǔn)周期循環(huán)過(guò)程，不僅可以實(shí)現(xiàn)強(qiáng)震復(fù)發(fā)周期，還可以模擬斷層的同震破裂過(guò)程.此外，相較于其他數(shù)值計(jì)算方法，有限單元法可以很容易處理復(fù)雜幾何邊界、復(fù)雜物性結(jié)構(gòu)、流固耦合等等高度非線性問(wèn)題，因此今后利用此方法可以將更加簡(jiǎn)單、方便的研究一些復(fù)雜條件下的地震序列問(wèn)題.

3 問(wèn)題討論

3.1 初始場(chǎng)的設(shè)定

與前人(Duan and Oglesby, 2005a, b; Lapusta and Liu, 2009; Noda and Lapusta, 2010)模擬地震序列問(wèn)題直接給定初始應(yīng)力場(chǎng)不同，本文中設(shè)定的應(yīng)力場(chǎng)是通過(guò)施加重力和構(gòu)造載荷共同作用的結(jié)果.具體的做法是：先在自由地表至地下2 km處施加重力加速度，通過(guò)地應(yīng)力平衡的方法消除模型在受到重力作用后產(chǎn)生的大變形；然后再施加構(gòu)造載荷，隨著構(gòu)造載荷的不斷累積，直至發(fā)生地震事件.

重力作用在地學(xué)數(shù)值模擬中是個(gè)非常重要的問(wèn)題，對(duì)比本文研究的地震序列問(wèn)題，以往的研究者(Duan and Oglesby, 2005a,b; Lapusta and Liu, 2009; Noda and Lapusta, 2010)都不予考慮.究其原因，如圖7a所示，當(dāng)對(duì)整個(gè)模型施加重力后，隨著深度的增加，其應(yīng)力大小也在不斷的增加.然而地球內(nèi)部的應(yīng)力場(chǎng)絕不僅僅是一種簡(jiǎn)單的、線性的力學(xué)過(guò)程，而是多圈層相互作用、受多種物理和化學(xué)因素控制的、非線性作用過(guò)程.通常認(rèn)為由于孔隙壓、構(gòu)造應(yīng)力的作用，達(dá)一定深度后，斷層面上的有效正應(yīng)力就不再隨著深度的增加而增加.然而若不考慮重力作用，整個(gè)模型特別是自由地表附近沒(méi)有了重力加速度的約束，顯然也是不合理的.因此參照前人在地下2 km以下給定均勻初始應(yīng)力場(chǎng)的做法(Lapusta and Liu, 2009)，本文在自由地表至地下2 km處施加重力加速度，如圖7b所示，兼顧深部應(yīng)力場(chǎng)大小和重力加速度的作用.

3.2 摩擦系數(shù)對(duì)地震準(zhǔn)周期的影響

為了考察摩擦系數(shù)對(duì)地震準(zhǔn)周期的影響，本文通過(guò)改變靜摩擦系數(shù)和動(dòng)摩擦系數(shù)計(jì)算多組模型進(jìn)行對(duì)比分析，圖8給出了其中4種典型情況下斷層典型點(diǎn)位置上下盤之間的位錯(cuò)隨時(shí)間的變化.根據(jù)圖8，比較4種情況下的計(jì)算結(jié)果：

圖5 每次地震事件斷層上各點(diǎn)滑移隨時(shí)間的變化Fig.5 Snapshots of the slip profiles at every node along the fault vary with time in each event

圖6 地震事件五中不同時(shí)刻斷層周邊介質(zhì)的振動(dòng)速度云圖Fig.6 Contour distributions of particle velocities at different times around the fault in event 5

圖7 重力作用下，垂直于地表方向應(yīng)力云圖(a) 對(duì)整個(gè)模型施加重力; (b) 自由地表至地下2 km處施加重力.Fig.7 The vertical stresses of the models under applying gravity(a) Gravity is applied to the whole model; (b) Gravity is applied to the region between free surface and 2 km underground.

(1) 比較圖8a、圖8b和圖8d，在動(dòng)摩擦系數(shù)相同時(shí)，當(dāng)斷層上的靜摩擦系數(shù)越大，斷層的強(qiáng)震復(fù)發(fā)周期越長(zhǎng).

(2) 比較圖8b、圖8c和圖3，在靜摩擦系數(shù)相同時(shí)，當(dāng)斷層上的動(dòng)摩擦系數(shù)越大，斷層的強(qiáng)震復(fù)發(fā)周期越短.

(3) 比較圖8c和圖8d，可以看出，當(dāng)斷層上的靜、動(dòng)摩擦系數(shù)差值相同時(shí)，斷層的強(qiáng)震復(fù)發(fā)周期也相當(dāng).

由上述分析可見：斷層上的靜、動(dòng)摩擦系數(shù)差值直接影響著斷層的強(qiáng)震復(fù)發(fā)周期，差值越小，復(fù)發(fā)周期越短；差值越大，復(fù)發(fā)周期越長(zhǎng).

此外，為了提高效率，研究本質(zhì)問(wèn)題，本文采用二維有限元模型進(jìn)行計(jì)算，但本文所采用的方法是完全可以用來(lái)計(jì)算三維模型的.在今后分析復(fù)雜的實(shí)際地震問(wèn)題時(shí)，會(huì)將此方法擴(kuò)展應(yīng)用到三維模型.

4 結(jié)論

本文基于有限單元?jiǎng)討B(tài)隱式算法發(fā)展出一種新的計(jì)算方法，連續(xù)穩(wěn)定的模擬了斷層孕震-同震準(zhǔn)周期循環(huán)過(guò)程.新的計(jì)算方法既可以模擬斷層的強(qiáng)震復(fù)發(fā)周期，同時(shí)也可以模擬斷層同震的詳細(xì)破裂過(guò)程.此外，我們研究了摩擦系數(shù)對(duì)地震準(zhǔn)周期的影響，發(fā)現(xiàn)斷層上的靜、動(dòng)摩擦系數(shù)差值直接影響著斷層的強(qiáng)震復(fù)發(fā)周期，差值越小，復(fù)發(fā)周期越短；差值越大，復(fù)發(fā)周期越長(zhǎng).

由于有限單元方法能夠適應(yīng)復(fù)雜的斷層幾何、復(fù)雜的邊界條件、復(fù)雜的材料屬性，因此相對(duì)于其他數(shù)值模擬方法在計(jì)算真實(shí)的地震、地質(zhì)問(wèn)題時(shí)更具有優(yōu)勢(shì).在今后定量分析復(fù)雜的實(shí)際地震問(wèn)題時(shí)，可以運(yùn)用本文所發(fā)展的方法開展深入研究.因此本研究對(duì)于探索強(qiáng)震孕育發(fā)生的動(dòng)力學(xué)機(jī)制及地震預(yù)測(cè)預(yù)報(bào)等有重要的科學(xué)意義及實(shí)際應(yīng)用價(jià)值.

圖8 不同摩擦系數(shù)下，斷層典型點(diǎn)位置上下盤之間的位錯(cuò)隨時(shí)間的變化Fig.8 Slip between hanging wall and footwall at the typical point on the fault vary with time in different models in which friction coefficients are different