吳永祺,張海明
北京大學地球與空間科學學院地球物理系,北京 100871
通過觀測儀器記錄的地面震動蘊含了有關地球介質和震源過程的信息,因此通過對比理論地震圖與實際觀測結果,可以了解地球內部結構和震源破裂過程.理論地震圖的計算需要設定震源的運動情況.早期的研究大多考慮均勻破裂的平面斷層模型,都做出了震源時間函數在斷層上各點一致(Hartzell,1978;Bouchon,1980)的假設.但是這種簡單的震源模型在很多情況下不能很好地描述震源的復雜性.例如實際的地震過程中,在斷層面上的滑動分布是不均勻的,它們可能由復雜的斷層幾何形狀、初始應力和摩擦系數等復雜性導致(Zhang et al.,2019).近二十年來,有很多基于復雜的震源模型的實際震例研究,例如1992年的Landers地震(Aochi and Fukuyama,2002;Wollherr et al.,2019),2008年的汶川地震(Wen et al.,2012;Zhang et al.,2019)和2010年的海地地震(Douilly et al.,2015).基于符合物理規律的震源破裂過程計算得到理論地震圖,并通過分析掌握波場的特征,是進一步進行反演獲得發震區域的應力場等信息的基礎.
地震波場與發震斷層上的破裂過程密切相關,而影響斷層破裂的因素之一是斷層的幾何復雜性.地震斷層不同的幾何形態,甚至平面斷層上微小的幾何擾動,都會導致不同的動態破裂過程,進而產生不同的地震波場(Madariaga,1983).斷層的幾何復雜性包括階躍(袁杰和朱守彪,2014;Hu et al.,2016),彎折(張麗芬等,2016),分叉(Kame et al.,2003;Bhat et al.,2007)和斷層面上的粗糙度(Shi and Day,2013;Luo and Duan,2018).多數破壞性的地震往往都發生在幾何形態復雜的斷層系統上,由于斷層空間延展尺度較大,斷層的幾何形狀會對破裂過程產生重要的影響(馬瑾等,1996).例如2008年的汶川地震,Zhang等(2019)考慮了斷層為平面斷層和非平面斷層的情形,非平面斷層的破裂會造成映秀縣和北川縣附近有明顯的地面運動,但平面斷層的破裂對北川縣附近沒有明顯的影響,汶川地震斷層的幾何復雜性極大地影響了震源的動態破裂和地震波場的分布.又例如1992年的Landers地震,該地震的一個特征是在幾個斷層之間的破裂轉移,基于動力學破裂的相關研究都能很好的重現該特征,并且都能給出與觀測結果擬合較好的地震波場結果(Aochi and Fukuyama,2002;Wollherr et al.,2019).
除了斷層的幾何形狀以外,初始破裂區位置、滑動弱化位移和初始應力等動力學參數也是影響地震波場的重要因素.初始成核區的不同位置,會導致不同的破裂過程,進而激發不同的地震波場.對于彎折斷層,在足夠的距離以積累能量的情況下,當破裂到達拐角處后,破裂方向的變化會輻射出強烈的地震波(Zhang et al.,2017).改變斷層上的滑動弱化位移大小,初始應力大小也會改變斷層的破裂過程.當滑動弱化距離較小,初始應力較大的情況下,往往容易產生超剪切破裂.超剪切破裂是研究人員重點關注的課題之一,在超剪切破裂過程中,斷層不同時刻輻射的剪切波互相干涉、疊加,形成馬赫波,使得地震波振幅更大,輻射范圍更廣,往往造成更嚴重的地震災害(Bernard and Baumont,2005;Dunham and Archuleta,2005;Dunham and Bhat,2008).超剪切地震的能量可以大幅度地傳播到更遠的地方,但是它在斷層附近的地震動要比亞剪切地震引起的地震動要弱(Huang et al.,2018).
復雜斷層系統的自發破裂過程模擬對計算方法提出較高的要求.目前常用的有兩類方法.一類為域方法,包括有限差分法(Madariaga et al.,1998;Day et al.,2005;Zhang et al.,2006,2014)、有限元方法(Aagaard et al.,2001a,b;Oglesby et al.,2003;Barall,2009).域方法在計算過程中會同時獲得破裂過程和地震波場.由于方法本身對于處理復雜斷層的情況有技術困難,因此采用這兩種方法研究三維的彎折和分叉斷層的研究相對較少.另一類方法是半解析的邊界積分方程方法(BIEM,Aochi et al.,2000;Chen and Zhang,2006;Zhang and Chen,2006a,b),將問題歸結到斷層平面上,因此相比于域方法,在處理復雜幾何斷層問題時(如彎折、分叉)更為簡便.
本文將基于三維空間非結構化的BIEM計算的自發破裂結果作為震源輸入,利用分層介質中的廣義反透射系數方法(Chen,1999)計算地震波場.針對彎折和分叉的情況,分別研究斷層幾何形狀的改變、初始成核區位置的不同以及超剪切破裂對地震波場的影響,最后討論介質模型對地震波場的影響.
本文通過BIEM計算斷層上的自發破裂過程,斷層應力與滑動速率的關系(Das and Aki,1977)為:
(x∈Γ)
(1)

在斷層自發破裂模擬的過程中,需要通過摩擦準則控制其破裂行為.而在同震破裂過程的模擬中,通常采用滑動弱化摩擦準則(Ida,1972):
(2)
其中,τ是當前時刻單元的剪切應力值,τu和τf分別是破裂強度和殘余應力,Dc為臨界滑動弱化位移,D為當前時刻的單元累積滑動量.
本文以自發破裂的結果作為輸入,運用廣義反透射系數法(Chen,1999)計算每個單元的地震波場,通過點源疊加的方式計算整個斷層系統的地震波場.當斷層距離地表幾公里以上時,地表的效應基本上可以忽略(Zhang and Chen,2006b),全空間介質模型與半空間模型的破裂過程基本一致.因此,本文的數值模擬中自發破裂過程采用無限空間模型,而計算地震波場采用半空間和分層半空間模型.
Bouchon(1980)利用離散波數法計算了半空間走滑斷層的位移場.本文與該結果進行比較,來驗證計算方法的正確性.計算模型如圖1所示,其中直立的斷層長30 km,寬5 km,埋深2 km,斷層破裂速度2.2 km·s-1,初始破裂位置在斷層最左側,從左向右均勻破裂.震源時間函數為斜坡函數,上升時間0.5 s.P波和S波波速分別為6.0 km·s-1和3.5 km·s-1,介質密度2.8×103kg·m-3.一共5條測線,在每條測線上各取10個觀測臺站,分別距離斷層跡線0,1,2,3,4,6,8,10,15,20 km.計算得到的位移場結果如圖2所示.通過比較可以發現,本文計算的結果和Bouchon的結果基本一致,在一定程度上驗證了本文用于計算地震波場的方法的可靠性.

圖1 Bouchon(1980)計算模型示意圖粗線表示接收點所在觀測線,觀測線共5條,測線1距斷層左側10 km,測線2穿過斷層左側,測線3距斷層左側15 km,測線4距斷層左側28 km,測線5距斷層左側50 km.每條測線上分布10個臺站,分別距離斷層跡線0,1,2,3,4,6,8,10,15,20 km.Fig.1 Geometry of the model in Bouchon (1980)The bold line indicates the observation line where the receiving point is located.There are 5 observation lines.Lines 1 is 10 km away from the left side of the fault,line 2 passes through the left side of the fault,line 3 is 15 km from the left side of the fault,line 4 is 28 km from the left side of the fault,and line 5 is 50 km from the left side of the fault.There are 10 stations on each survey line,which are 0,1,2,3,4,6,8,10,15,20 km away from the fault trace.
實際地震的幾何形狀非常復雜,如1992年蘭德斯MW7.3地震,2016年凱庫拉MW7.9地震.為了研究復雜的斷層模型,可以將斷層復雜的幾何形態抽象不同幾何模型的組合,如彎折、分叉、階躍.這樣便可以針對每一類幾何模型,研究其如何影響地震波場.其中彎折斷層是一種最為典型的斷層幾何結構,因此研究不同形態的彎折斷層,可以進一步理解彎折斷層會對地震波場產生何種影響.相較于彎折斷層,分叉斷層具有更加復雜的幾何形態,其產生的地震波場也將具有更加復雜的特征,同時考慮到邊界積分方法在處理復雜幾何模型的優勢,我們也將對分叉斷層模型進行考察.因此,本文將基于彎折斷層和分叉斷層進行地震波場性質的研究.本文以下的數值計算,均采用相同的半空間介質模型,P波和S波波速分別為6.0 km·s-1和3.5 km·s-1,介質密度2.8×103kg·m-3.斷層為走向90°的走滑斷層.
對于彎折斷層,不同的彎折角會影響破裂的傳播和停止.通常,過大的彎折角將抑制破裂的傳播,因此我們考察在其他動力學參數不變的情況下,彎折斷層不同彎折角對于地表地震波場分布的影響.斷層寬10 km,主平面長30 km,分支平面長15 km,如圖3所示.考慮彎折角分別為0°、20°和30°的情況,初始破裂位置在距離主平面左側20 km處.斷層上破裂強度為60 MPa,成核區內初始應力為72 MPa,成核區外初始應力為48 MPa,靜摩擦系數為0.6,動摩擦系數為0.1,成核區半徑2.5 km.圖4為不同彎折角的彎折斷層的破裂速度快照,隨著彎折角的增大,分支平面的破裂速度變小.圖5和圖6分別為地表峰值位移(PGD)分布和地表峰值速度(PGV)分布.如圖所示,對于東西(EW)分量和豎直(UD)分量,當彎折角為0°時(即平面斷層),PGD和 PGV都關于斷層對稱分布;當斷層彎折角增加為20°后,由于斷層角度的改變,PGD和PGV最大值出現在彎折平面斷層的走向方向(即東南方向)上,而東北方向上的地震波場減弱;當斷層彎折角增加到30°后,由于彎折斷層彎折角的增大會抑制破裂的傳播,對應方向上地震波場會減弱,能觀察到斷層右側地表地震波場明顯被削弱.對于南北(NS)分量,當彎折角為0°時,如圖5a2和圖6a2 所示,地震波場的NS分量表現出明顯的方向性效應,即在沿著破裂的方向上,垂直于主斷層平面走向方向的地震波場分量表現出最大值.而隨著彎折角的增大,由于彎折斷層的阻礙效應,如圖5c2和圖6c2 所示,在斷層右側的地震波南北分量會被削弱.由于初始破裂位置位于主平面斷層上,破裂向兩側傳播,因此過初始破裂位置且垂直于主平面斷層的位置,地震波場的PGD和PGV值最小.在固定其他計算參數的情況下,改變彎折平面斷層的彎折角,會使得PGD(或PGV)的最大值沿著彎折方向移動;同時,由于彎折斷層的阻礙效應,當彎折角過大時,會對該方向上的地震波場起到抑制作用(張麗芬等,2016).

圖3 彎折斷層模型斷層寬10 km,主平面(白色)長30 km,彎折平面(灰色)長15 km,斷層上沿距地表15 km.五角星表示初始破裂位置,距斷層左側20 km.Fig.3 Geometry of the bending faultThe width of the fault is 10 km,the length of main plane (white)is 30 km,and the length of bend plane (gray)is 15 km,the upper edge of the fault is 15 km from the surface.The star represents the initial rupture position and is 20 km from the left side of the fault.

圖4 不同彎折角的彎折斷層滑動速率快照每個子圖左側數字表示斷層破裂的不同時刻.(a)彎折角0°;(b)彎折角20°;(c)彎折角30°.Fig.4 Snapshots of the slip rate of bending faults with different bending anglesThe numbers on the left indicate different moment in the fault rupture process.(a)Bending angle 0°;(b)Bending angle 20°;(c)Bending angle 30°.

圖5 不同彎折角的彎折斷層引起的地表峰值位移場(PGD)分布圖從上到下分別為東西(EW)分量,南北(NS)分量和豎直(UD)分量的PGD分布圖.五角星為初始破裂位置,黑色實線為斷層在地表投影.(a)彎折角0°;(b)彎折角20°;(c)彎折角30°.Fig.5 Distribution map of the peak ground displacement (PGD)caused by bending faults with different bending anglesFrom top to bottom are the PGD distribution diagrams of the east-west (EW)component,the north-south (NS)component and the vertical (UD)component.The star is the initial rupture position,and the solid black line is the projection of the fault on the surface.(a)Bending angle 0°;(b)Bending angle 20°;(c)Bending angle 30°.

圖6 不同彎折角的彎折斷層引起的地表峰值速度(PGV)分布圖(a)彎折角0°;(b)彎折角20°;(c)彎折角30°.其余圖例說明同圖5.Fig.6 Distribution of peak ground velocity (PGV)caused by bending faults with different bending angles(a)Bending angle 0°;(b)Bending angle 20°;(c)Bending angle 30°.The rest of the illustrations are the same as Fig.5.
在同一種幾何模型下,其他動力學參數相同,成核區位置的不同也會造成不同的破裂過程,進而導致不同的地震波場.因此,需要考慮不同成核區位置對地震波場的影響.計算模型參考圖3,斷層為彎折角30°的彎折斷層,主斷層長30 km,分支斷層長15 km.選取成核區位置位于距離主斷層左側5 km,10 km和15 km處.圖7和圖8分別為計算得到的PGD和PGV三分量結果.可以看到,無論是PGD還是PGV,成核區距斷層左側距離從5 km到15 km變化時,破裂從以從左向右破裂為主變化為雙側破裂,地震波場峰值分布也呈類似的變化.當初始破裂位置在斷層左側時,對地表地震波場的影響主要集中在右側,同時,在地震波場的南北分量上表現出明顯的方向性效應,即沿著斷層破裂方向上地震波場峰值最大,而逆破裂方向上地震波場峰值最小.隨著初始破裂位置的右移,斷層左側的地震波場逐漸變大,在斷層兩側的地表都會出現較大的破壞.值得注意的是,當初始破裂區距離彎折斷層的彎折處最遠,如圖7a 和圖8a 所示,斷層破裂在沿破裂方向上產生的地震波場造成的破壞也最大,這是因為距離越遠,破裂具有更大的能量,斷層將在地表產生更大的地震波場 (Xu et al.,2019).因此,在其他參數相同的情況下,對于彎折斷層,不同成核區位置的地震,會對地表造成截然不同的破壞,當破裂在斷層上的破裂距離越長,破裂具有更大的能量,地震在地表造成的破壞也越大.

圖7 彎折斷層不同初始成核區位置地表峰值位移場(PGD)分布圖(a)初始破裂位置距斷層左側5 km;(b)初始破裂位置距斷層左側10 km;(c)初始破裂位置距斷層左側15 km.其他圖例說明同圖5.Fig.7 Distribution of peak ground displacement (PGD)at different initial nucleation areas of bending faults(a)The initial rupture location is 5 km away the left side of the fault;(b)The initial rupture location is 10 km from the left side of the fault;(c)The initial rupture location is 15 km from the left side of the fault.Other legends are the same as Fig.5.

圖8 彎折斷層不同初始成核區位置地表峰值速度(PGV)分布圖(a)初始破裂位置距斷層左側5 km;(b)初始破裂位置距斷層左側10 km;(c)初始破裂位置距斷層左側15 km.其他圖例說明同圖5.Fig.8 Distribution of peak ground velocity (PGV)at different initial nucleation areas of bending faults(a)The initial rupture location is 5 km away the left side of the fault;(b)The initial rupture location is 10 km from the left side of the fault;(c)The initial rupture location is 15 km from the left side of the fault.Other legends are the same as Fig.5.
在同一種幾何模型下,其他動力學參數相同,改變初始應力大小或者滑動弱化距離會得到不同的破裂過程,其中,當初始應力較大或者滑動弱化距離較小的情況下可能產生超剪切現象,而斷層上超剪切破裂的發生往往會導致地表嚴重的破壞.因此,我們需要重點關注超剪切破裂對地表位移場的影響.我們使用無量綱化公式(Xu et al.,2015;鄭玲瓏等,2021)模擬平面斷層的動態破裂過程:
(3)


圖9 平面斷層自發破裂滑動速率快照每個子圖左側數字表示斷層破裂的不同時刻.(a)模型S1破裂自發停止;(b)模型S2亞剪切破裂;(c)模型S3超剪切破裂.Fig.9 Snapshot of the spontaneous rupture slip rate on the plane faultThe numbers on the left indicate different moment in the fault rupture process.(a)Model S1;(b)Model S2 with subshear rupture;(c)Model S3 with supershear rupture.

圖10 平面斷層地表位移場從上到下分別為東西(EW)分量,南北(NS)分量和豎直(UD)分量.黑線、藍線和紅線分別表示斷層模型S1,S2和S3的位移場三分量.Fig.10 Surface displacement field of the plane faultFrom top to bottom are the east-west (EW)component,the north-south (NS)component and the vertical (UD)component.The black line,blue line and red line represent three components of displacements of the fault models S1,S2 and S3,respectively.

表1 不同初始應力平面斷層模型Table 1 Plane fault models with different initial stresses

圖11 平面斷層地表峰值速度(PGV)分布從上到下分別為東西(EW)分量,南北(NS)分量.(a)發生超剪切破裂的模型S3 地表PGV分布;(b)發生亞剪切破裂的模型S2 地表PGV分布.其他圖例說明同圖5.Fig.11 Distribution of peak ground velocity (PGV)of the plane fault From top to bottom are the east-west (EW)component,the north-south (NS)component.(a)Model S3 PGV with supershear rupture;(b)Model S2 PGV with subshear rupture.Other legends are the same as Fig.5.
考察夾角為30°的不同形狀的分叉斷層對位移場的影響.如圖12所示,主斷層長30 km,兩個分支斷層長20 km.斷層模型B1,B2,B3,B4各分支的張角(即與主平面的夾角,逆時針為正)由表2給出.動力學參數如2.1節所述.圖13顯示了各個模型中過斷層中心線上的位置處的滑動速率隨時間的變化.可以觀察到模型B1和模型B2分叉面A破裂占優,而模型B3和模型B4分叉面B破裂占優.隨著分叉面A角度越來越大,破裂從分叉面A逐漸轉移到分叉面B.同時,無論是在哪個分叉面破裂,優勢分支的角度越小,滑動速率越大.計算的地表位移場如圖14所示,地表位移場的波形特征與斷層破裂結果有著明確的對應關系.在接收點a,位移場結果如圖14a 所示,由于模型B1和模型B2的破裂時間更長,因此,能觀測到更長的波形.同時對于模型B4,由于其在斷層分支上的破裂速度最大,破裂最強烈,故該斷層對應的位移場峰值在幾個模型中也是最大的.在接收點b的位移場如圖14b所示,能看到在南北和豎直分量上,幾個斷層模型計算得到的位移場十分接近.這是由于將接收點沿著南北方向做投影,投影點會落在斷層的兩個分支上,考慮上文所提到的方向性效應,因此分叉斷層的兩個分支平面在南北分量和豎直分量的位移場值貢獻很小,位移場主要由相同的主平面主導,幾個模型得到位移結果相近.增大震中距,讓觀測點的投影落在斷層外,不同斷層模型的位移場差別將變大,結果如圖14c 所示.

圖12 分叉斷層模型斷層寬10 km,主平面長30 km,兩個分支平面長20 km,斷層上沿距地表15 km.五角星表示初始破裂位置,三角形表示接收點.接收點a震中距71 km,方位角56.0°;接收點b震中距75 km,方位角25.5°;接收點c震中距145 km,方位角26.7°.Fig.12 Branching fault modelThe width of the fault is 10 km,the length of main plane is 30 km,and the length of two branching planes is 20 km.The upper edge of the fault is 15 km from the surface.The star represents the initial rupture position.The triangle represents the receiver.The epicenter distance of receiver a is 71 km and the azimuth is 56.0°;the epicenter distance of receiver b is 75 km and the azimuth is 25.5°;the epicenter distance of receiver c is 145 km and the azimuth is 26.7°.

圖13 分叉斷層分支滑動速率時空分布圖虛線左側為分叉面A,右側為分叉面B.(a)模型B1;(b)模型B2;(c)模型B3;(d)模型B4.Fig.13 Space-time diagrams of the slip rate along the branch planes of the branching faultsThe left side of the dashed line is surface A,and the right side is surface B.(a)Model B1;(b)Model B2;(c)Model B3;(d)Model B4.

表2 不同分叉斷層模型Table 2 Different branching fault models

圖14 不同的分叉斷層模型上自發破裂導致的地表位移隨時間的變化每個子圖中從上到下分別為東西(EW)分量,南北(NS)分量和豎直(UD)分量.B1,B2,B3,B4為四個斷層模型.(a)接收點a計算結果;(b)接收點b計算結果;(c)接收點c計算結果.Fig.14 Displacement caused by spontaneous rupture on different branching faultsFrom top to bottom are the east-west (EW)component,the north-south (NS)component and the vertical (UD)component in the subplot.B1,B2,B3 and B4 are four fault models.(a)The result of receiver a;(b)The result of receiver b;(c)The result of receiver c.
除了上述各種因素以外,介質模型也是一個可能影響地表位移場的因素.考慮如表3、表4、表5所示的介質模型半空間模型L1,兩層介質模型L2,含低速層的三層介質模型L3.考慮如上文所示彎折角20° 的彎折斷層,斷層深度15 km.位移場計算結果如圖15所示.相較于其他介質模型,由于低速層速度更小,L3模型對應各震相到達更遲.同時,由于低速層對位移場的放大效應,能觀察到L3模型下的位移場振幅更大.

圖15 不同介質模型下的位移場從上到下分別為東西(EW)分量,南北(NS)分量和豎直(UD)分量.藍線、紅線和綠線分別表示半空間介質模型L1,兩層介質模型L2和含低速層的三層介質模型L3.Fig.15 Displacement field under different media modelsFrom top to bottom are the east-west (EW)component,the north-south (NS)component and the vertical (UD)component.The black line,the blue line and the red line represent the half-space medium model L1,the two-layer medium model L2,and the three-layer medium model L3 with low velocity layers,respectively.

表3 半空間介質模型L1Table 3 Half-space medium model L1

表4 兩層介質模型L2Table 4 Two-layer medium model L2

表5 含低速層的三層介質模型L3Table 5 Three-layer medium model with low-velocity layer L3
本文基于BIEM模擬了斷層的自發破裂過程,并根據此結果計算斷層對應的地震波場.討論了不同斷層幾何形狀,成核區位置、超剪切破裂以及介質模型分別對地震波場的影響,得到如下結論:
(1)動力學參數不變,對于不同彎折角的彎折斷層,彎折角度變化能明顯改變地表地震波場的分布.對于小彎折角的情形,沿彎折分支走向的方向上地震波場具有最大值;但是,對于大彎折角的情形,由于彎折斷層的抑制效應,會使得沿彎折斷層走向方向上的地震波場都較小.
(2)對于彎折斷層,其他動力學參數不變,改變初始成核區的位置,當破裂在斷層上的傳播距離越長,破裂的能量越大,地震在地表造成的破壞越大.
(3)對于彎折斷層,其他動力學參數不變,增大斷層面上的初始應力,能增大斷層上破裂的滑動速率,當初始應力足夠大時,斷層上會產生超剪切破裂,而超剪切的破裂的產生會明顯增大地震波場的振幅,激發高頻成分,使得部分震相到時提前,并增大永久位移,同時,超剪切破裂產生的地震波場在空間中的衰減更為緩慢,從一定程度上能解釋為什么發生超剪切的地震會造成更大的破壞.
(4)對于分叉斷層,在斷層上破裂時間越長,對應的地震波場具有更寬的波形;斷層上破裂強度越大,對應的地震波場的峰值也越大.
(5)低速層會增大地表波場的振幅,同時會使得波形上各震相到達時間更遲.
本文探討了斷層動力學自發破裂過程與其產生的地震波場的聯系,計算了不同幾何形態和不同動力學參數下自發破裂斷層引起的地震波場,可以由此評估不同斷層破裂對地表造成的破壞情況.實際地震具有更高的復雜性,包括更復雜的斷層幾何特征、更復雜的動力學參數、更復雜的介質信息,因此,有必要更加細致的考察斷層的各種復雜性是如何影響地震波場的特征,為下一步的動力學參數反演提供基礎.同時將模擬的地震波場分布與實際的地震破壞比較,以便于為未來的災害評估提供參考意見.這些是進一步的研究中予以考慮的重要內容.