2017年九寨溝M7.0地震強地面運動估計與方向效應研究

張嚴方 胡其志,2 張 潔 游 姍

1 武昌工學院城市建設學院,武漢市白沙洲大道110號,430065 2 湖北工業大學土木建筑與環境學院,武漢市南李路28號,430068

由于強震事件的不可重復性與我國強震觀測能力分布的不均勻性,以及近場強震觀測資料的缺失,導致當前強地震動場的分析存在實測波形不足的困難,前人利用格林函數法、有限斷層法及三角級數疊加法等開展強地面運動估計與波形的合成[1-5]。2017-08-08九寨溝M7.0地震是2008年汶川M8.0地震與2013年雅安M7.0地震后四川地區發生的又一次破壞性地震,中國地震局發布的地震烈度分布顯示,宏觀震中比芒村地震烈度為Ⅸ度,近場PGA達707 Gal。現場建筑結構震害與地質災害情況揭示,近場強地震動具有較為顯著的方向效應特征,即平行斷層方向的結構震害顯著大于垂直斷層方向[6]。

本文基于二維有限單元法使用接觸單元模擬斷裂帶的非連續破裂過程,建立2017-08-08四川九寨溝M7.0地震近場強地面運動估計模型,并對應力降和屈服應力等關鍵震源參數進行研究。在利用國家強震數據中心提供的實際觀測波形開展地震動參數信度檢驗基礎上,給出PGA及儀器地震烈度分布結果,探討近場強地面運動方向效應特征及其成因機理。通過給出一種具備明確震源物理意義且參數設置簡單的強地面運動模擬方法,可在一定程度上克服當前觀測能力不足所造成的困難,對近斷層強地震動特性研究具有一定的科學意義和理論增量。

1 基于有限元的地震動模擬方法

1.1 斷層面接觸模型建立方法

發震斷層連續破裂過程的運動學方程描述了破裂行為的不連續性,將其作為非線性問題進行處理,分別使用平面應變單元和接觸單元模擬孕震區域的基巖和發震斷層。使用由2個固體單元組成的Goodman接觸單元模擬發震斷層[7],具有平行、垂直或旋轉等3種不同的運動方式,如圖1所示,當引入拉力時,接觸單元將會分離,作用力通過接觸面的傳遞將會被切斷[8]。

圖1 Goodman接觸單元示意圖Fig.1 Schematic diagram of Goodman contact element

當動剪切應力達到接觸面的剪應力屈服值τy時,單元將會發生滑動。此時,粘聚力C、正應力σ和內摩擦角φ存在以下關系式:

τy=C+σtanφ

(1)

1.2 初始應力分布

首先求解由重力產生的應力場,其應力分布可表示為:

(2)

式中,N為形函數矩陣,{p}={0-ρg}T,ρ為密度,g為重力加速度。

假設斷層面上由重力產生的剪應力在之前的歷史地震事件中被釋放,現階段組成斷層的接觸單元的剪應力被減小至0。然后在模型兩端邊界施加構造應力場,通過逐漸增加構造應力,直至接觸單元的剪應力到達屈服強度,確定地震事件破裂的臨界應力狀態,并將其轉移到動力分析模型上[9]。

1.3 動力分析

為分析斷層面上的破裂情況,需考慮滑動和分離現象。動力分析為非線性分析[8],第n時步的運動方程為:

(3)

式中,M、C、K分別為質量、阻尼及剛度矩陣,{g(t,s)}為由動應力降Δτ計算得到的外力,t和s為滑動發生的時間和節點對。在初始時步(t=0)有:

(4)

式中,l為接觸單元的長度,在接觸單元滑動時剪應力τ為定值,且大小為τ=τy-Δτ,其中τy為屈服應力。將該節點的外力{g(0,s)}施加于破裂起始位置處的節點對上,同時確保外力方向必須與初始剪應力方向一致。施加的節點外力會使應力重新分布,并使相鄰接觸單元的剪應力增加。

2 九寨溝M7.0地震強地面運動模擬

2.1 九寨溝地震概況

2017-08-08九寨溝M7.0地震震中為33.20°N、103.82°E,震源深度約20 km,宏觀震中位于九寨溝核心景區西部5 km處的比芒村,震中烈度為Ⅸ度,近場PGA達707 Gal,為左旋走滑型地震。易桂喜等[6]通過研究余震分布認為,該地震發震斷層為以左旋走滑為主的NW-SE向樹正斷裂,長度約38 km。表1及圖2為九寨溝地區CSMN強震測站、活動斷裂及地震序列等相關信息。

表1 CSMN臺站信息Tab.1 Information of CSMN stations

F01:東昆侖斷裂;F02:白龍江斷裂;F03光蓋山-迭山斷裂;F04:哈南-青山灣-稻畦子斷裂;F05:塔藏斷裂;F06:雪山梁子斷裂;F07:虎牙斷裂;F08:樹正斷裂;F09:岷江斷裂;F10:龍日壩斷裂圖2 九寨溝區域活動斷裂、CSMN強震測站及余震序列分布Fig.2 The distribution of active faults,CSMN strong earthquake stations and aftershock events in Jiuzhaigou area

2.2 有限元模型

九寨溝地震震源參數為走向152°、傾角88°、滑動角-3°,平均位錯1.4 m,應力降8.5 MPa?？紤]到發震斷層近乎垂直且震源機制為左旋走滑,故將該地震事件簡化為二維平面問題。結合區域活動斷裂分布,利用接觸單元及研究區內活動斷裂,建立九寨溝地區強地面運動數值模型,網格劃分情況見圖3。

圖3 九寨溝地區二維有限元模型Fig.3 Two-dimensional finite element model of Jiuzhaigou area

巖體介質的物理力學參數對數值模擬結果具有重要影響,利用國家強震動觀測臺網中心(CSMN)提供的九寨溝地震觀測波形計算建模區P波與S波波速。由圖4可知,研究區測站震中距與彈性波到達時間線性相關,因此將基巖視為彈性均質介質,重度為28.0 kN/m3,泊松比為0.3,VP=6.1 km/s,VS=3.5 km/s,與徐晶等[10]的研究結果基本一致。

圖4 九寨溝地區P波及S波波速Fig.4 Velocity of P wave and S wave of Jiuzhaigou area

2.3 參數敏感度與信度檢驗

在假設實際應力降未知的情況下,將應力降范圍設為6.0～12.0 MPa,取值間隔為0.5 MPa;將斷層接觸單元的屈服應力范圍設為2.0～5.0 MPa,取值間隔為1.0 MPa。經處理后可獲得研究區內水平地表任意位置的地震動加速度波形,將EW向與NS向地震動分量進行矢量合成,獲得水平向PGA,則:

PGA=max(abs(ah))

(5)

ah=aEW+aNS

(6)

式中,ah為某場址處水平向地震動加速度時程;aEW和aNS分別為EW向及NS向地震動加速度分量。

為直觀分析不同震源參數對水平向PGA的影響,取九寨百河(51JZB)及平武木座(51PWM)2個臺站場址處的PGA計算結果,其中近場九寨百河臺(51JZB)震中距約為31.73 km,遠場平武木座臺(51PWM)震中距約為91.30 km。

圖5(a)為當屈服應力取值4 MPa時,不同應力降條件下近場及遠場PGA的變化情況。由圖可見,隨著應力降取值由2.5 MPa增大至5.5 MPa,PGA在近場由192.52 Gal增大至318.15 Gal,遠場則由20.45 Gal增大至34.77 Gal,說明PGA在整體上隨應力降增大的同時,其近場地震動受到的影響較遠場更為顯著。究其原因可能在于,高頻地震動隨震中距增大而產生衰減效應,這與李宗超等[11]使用隨機有限斷層法獲得的結論一致。

圖5 不同條件下的PGAFig.5 PGA under different conditions

圖5(b)為當應力降取值8.5 MPa時,不同屈服應力條件下近場及遠場PGA的變化情況。由圖可見,隨著屈服應力取值由2.5 MPa增大至5.5 MPa,PGA在近場由128.40 Gal增大至380.44 Gal,遠場則由7.71 Gal增大至61.63 Gal。說明隨著屈服應力增大,PGA呈現出整體持續增大但遠場衰減的變化特征,且其對屈服應力的響應較應力降更加明顯。

引入模型偏差k及模型標準差s以分析2類震源參數對水平向PGA的影響[12],表達式為:

(7)

(8)

式中,n為臺站數量,k越大表示模擬結果越偏離觀測結果,s越大則離散程度越高。

通過對比不同震源參數條件下水平向PGA模擬值與實際觀測值的偏差及標準差發現,當屈服應力變化時,PGA偏差約為6.73,標準差約為0.35;而應力降變化時,PGA偏差約為3.24,標準差約為0.46,說明斷層接觸單元屈服應力與應力降均對PGA具有重要影響,且屈服應力影響更為顯著。在震源參數敏感性分析基礎上,通過對不同震源參數組合進行試算發現,應力降為8.5 MPa且屈服應力為4.0 MPa時,獲得的模擬值與CSMN強震觀測結果具有最好的一致性。圖6為該最優震源參數組合下的PGA觀測值與模擬值,可以發現,51JZB、51JZY和51JZW 3個近場強震臺站場址處的PGA明顯高于其他8個遠場臺站。

圖6 最優震源參數組合下的PGA觀測值與模擬值Fig.6 PGA observations and simulation results under the optimal source parameter combination

利用相對誤差百分比R來評價模擬誤差隨震中距的變化情況:

R=|PGAh,simu-PGAh,ob|/PGAh,ob×100%

(9)

式中,PGAh,simu與PGAh,ob分別表示水平向PGA的模擬值和觀測值。

圖7為相對誤差隨震中距的分布情況,其中近場強震臺站51JZB、51JZW和51JZY的相對誤差均小于5%;隨著震中距增加,相對誤差呈現明顯增大的趨勢,可能與有限元空間分網的高頻濾波效應及地震動傳播過程中的介質復雜性有關。總體而言,強地面運動模擬結果的相對誤差百分比整體上均小于40%,與李明會等[13]使用復合震源模型及Dang等[14]使用隨機有限斷層法獲得的結果相當,可用于后續強地面運動場特征分析。圖8為九寨百合臺(51JZB)地震動加速度觀測波形與模擬波形。

圖7 不同臺站處模擬誤差隨震中距的變化情況Fig.7 Variation of simulation error with epicenter distance at different stations

圖8 51JZB臺地震動加速度觀測波形與模擬波形Fig.8 The observed and simulated waveforms of ground motion acceleration at 51JZB station

3 強地面運動特征分析

水平向地震動加速度是近場強地面運動最重要的地震動參數之一,在表征地震能量釋放的同時,也是抗震設防的重要參考依據?；诙S有限元方法,本文在參數研究和不確定性分析基礎上獲得九寨溝M7.0地震水平向PGA空間分布。由圖9可見,近斷層地區水平向PGA等值線整體呈橢圓狀分布,橢圓長軸與發震斷層樹正斷裂走向具有較好的一致性。震中附近PGA最大值約為720 Gal,且整體上隨著震中距增大而不斷衰減?？傮w而言,本文地震動峰值加速度場結果與周紅[15]使用隨機有限斷層方法獲得的結論大體一致,符合對地震動加速度分布的常規認識。值得指出的是,在PGA分布呈現出較為顯著的方向效應的同時,亦呈現出斷層破裂的方向效應及其與斷裂帶空間分布的一致性,值得進一步展開分析。

圖9 九寨溝地震水平向地震動峰值加速度場Fig.9 Horizontal peak ground acceleration field of Jiuzhaigou earthquake

3.1 斷層破裂的方向效應

沿斷裂帶破裂方向的地震動幅值明顯比逆破裂方向大,這種斷層破裂的方向效應會使近斷層地區地震動產生劇烈的空間變化。為研究九寨溝地震斷層破裂的方向效應,以震中破裂起始位置為起點,計算沿破裂方向及逆破裂方向的水平向PGA。由圖10可見,2個方向的PGA均呈現出隨震中距增大而衰減的特點,但隨著震中距持續增大,2個方向的PGA開始呈現出差異性,沿破裂方向的PGA顯著大于逆破裂方向。當震中距分別為5 km、20 km及40 km時,逆破裂方向的PGA約為沿破裂方向的96.93%、93.92%和54.24%,說明斷層破裂的方向效應在近斷層地區較微弱,但隨著震中距增大,在遠場區域更加明顯。

3.2 地震動傳播的方向效應

建立垂直于發震斷層走向的剖面,將該方向PGA與沿斷層破裂方向PGA進行對比,以分析九寨溝M7.0地震的地震動方向效應。由圖11可見,隨著震中距增大,2個方向的PGA均出現較為顯著的衰減現象,但垂直斷層方向的PGA衰減速度更快,呈現出地震動傳播的方向效應。當震中距分別為5 km、20 km及40 km時,垂直斷層方向的PGA約為平行斷層方向的94.63%、44.02%和38.07%,說明地震動傳播的方向效應在近斷層地區非常顯著,且隨著震中距增大,在遠場區域愈發明顯。

3.3 斷裂帶的隔震效應

由圖2和圖9可知,水平向PGA的等值線在一定程度上與區域斷裂帶的空間分布一致,如雪山梁子斷裂、龍日壩斷裂及白龍江斷裂等。為分析區域活動斷裂在近場強地震動傳播中的作用,建立跨雪山梁子斷裂近NS向剖面C-C′及跨龍日壩斷裂北段EW向剖面D-D′,圖12為2個剖面的PGA衰減情況。

圖12 PGA衰減Fig.12 PGA attenuation

由圖12(a)可見,水平向地震動在由近場向遠場傳播過程中,跨雪山梁子斷裂NS向PGA由110 Gal突降至85 Gal,降幅為22.72%;由圖12(b)可見,跨龍日壩斷裂北段EW向PGA由42 Gal降低到35 Gal,降幅為16.67%,呈現出斷裂帶減震消能作用。究其成因機制,可能是斷層使地震動沿斷層面的透射明顯降低,呈現出顯著的隔震特性[16]。

4 結 語

基于二維有限單元法利用接觸單元模擬斷裂帶的非連續破裂過程,建立2017-08-08四川九寨溝M7.0地震近場強地面運動估計模型,在參數研究及信度檢驗基礎上,探討斷層破裂的方向效應、地震動傳播的方向效應及其成因機理。研究結果表明:

1)隨著震源接觸單元屈服應力及應力降的增大,水平向PGA整體增大,屈服應力對于強地面運動的影響更為顯著;當應力降為8.5 MPa、屈服應力為4.0 MPa時,強地面運動模擬結果的相對誤差百分比整體小于40%,近場強震臺站的相對誤差小于5%。

2)九寨溝M7.0地震水平向PGA模擬結果等值線整體呈橢圓狀,橢圓長軸與發震斷層樹正斷裂走向具有較好的一致性;震中附近PGA最大值約為720 Gal,且整體隨震中距的增大不斷衰減,跨龍日壩斷裂北段及跨雪山梁子斷裂PGA分別降低16.67%及22.72%,這種斷裂帶的隔震效應或與地震動沿斷層面的透射降低有關。

3)沿斷層破裂方向的PGA幅值明顯大于逆斷層破裂方向,呈現出較為顯著的斷層破裂方向效應,震中距40 km處逆斷層破裂方向的PGA約為沿破裂方向的54.24%;震中距40 km處垂直斷層方向的PGA約為平行斷層方向的38.07%,呈現出地震動傳播的方向效應。