基于FK方法和GP14.3震源模型的2023年土耳其Mw7.8級地震寬頻地震動合成

巴振寧,韓書娟,趙靖軒,劉 悅,蘆 燕,陳三紅

(1. 中國地震局地震工程綜合模擬與城鄉抗震韌性重點實驗室(天津大學),天津 300350; 2. 天津大學 土木工程系,天津 300350;3. 中國地震災害防御中心,北京 100029)

0 引言

北京時間2023年2月6日9時17分土耳其發生Mw7.8地震,震源深度17.9 km,最高烈度為Ⅸ度,震中(北緯37.17°、東經37.03°)位于東南部的東安納托利亞斷裂帶(east anatolian fault zone,EAFZ)[1]。東安納托利亞斷裂帶位于安納托利亞板塊(亞歐板塊內部的次級板塊)和阿拉伯板塊的交界處如圖1所示,是土耳其境內的一條主干左旋走滑斷裂帶,起于西南的哈塔伊(Hatay),終止于東北的卡利若瓦(Karliova)三聯結點位置[2]。東安納托利亞斷裂帶地震活動頻繁,20世紀以來多次發生6級以上強震[3]。同日下午18時24分,再次發生Mw7.5地震,震源深度10 km,最高烈度亦為Ⅸ度,震中位于首次地震震中東北部96 km處[1]。這兩次地震構成了典型且罕見的“雙震”地震序列[4]。

圖1 土耳其地震的震中及其相關的構造環境Fig. 1 Epicenter of the Turkey earthquake and its associated tectonic environment

此次地震震級大、震源淺和地面運動強度大,造成了慘重的人員傷亡和嚴重的建筑結構破壞:41156人死亡、105505人受傷、超過8.4萬棟建筑倒塌和嚴重受損或急需拆除[5]。其中:位于土耳其東南部的人口密集的工業城市-加濟安泰普[6],距離此次Mw7.8地震震中僅37 km,是此次地震中受影響最為嚴重的城市之一。針對該Mw7.8地震,本文采用作者提出的FK方法開展0～10Hz的近斷層寬頻地震動合成[7],并重點分析加濟安泰普地區的地震動時空場分布特征,以期為該地區的震害原因分析和地震區劃等工作提供參考。

1 地下結構與震源模型建立

確定性物理寬頻地震動模擬依賴于能激發寬頻地震波的合理震源模型和能模擬寬頻地震波傳播的精確高效方法[8]。震源模型采用的是GP14.3運動學混合震源模型[9],該模型由Graves和Pitarka(2010)寬頻帶地面運動模擬方法(GP2010)改進后得到(GP14.3為在當前南加州地震中心寬帶平臺上的更新版本名稱)[10],其能夠同時考慮破裂面上低波數的確定性錯動和高波數的隨機性錯動,具備激發寬頻地震波場的能力。地震波傳播模擬方面采用作者提出的FK法,該方法是一種具有嚴密理論基礎的地震動合成方法,其適用于跨尺度一維速度結構寬頻帶地震波的傳播模擬。本節詳細介紹一維波速結構和運動學混合震源的構建。

1.1 加濟安泰普地區地下一維速度結構模型的構建

考慮到后續分析主要針對加濟安泰普市,地下一維速度結構基于該市進行構建,相應參數表1,對應的示意圖如圖2所示。

表1 加濟安泰普地區的地下一維速度結構Table 1 1D velocity structure in Gaziantep region

圖2 加濟安泰普地區的一維速度結構模型示意圖

該一維地下結構模型包括3部分:1)依據USGS網站提供的Vs30數據[11],確定0～0.03 km深度范圍內的一維速度結構參數,其中剪切波速Vs取為加濟安泰普市界范圍內Vs30的平均值。2)依據INFANTINO等[12]給出的伊斯坦布爾地區波速模型確定方法,確定0.03～2 km及2～5 km深度范圍內的一維速度結構參數:在0.03～2 km之間,Vs隨深度z的增加以(z)0.5的速率增加;在2～5 km之間,假定Vs呈線性增加,通過對更密集的深度點進行線性插值、分層取均值得到每層土體的Vs。3)依據全球地殼模型Crust1.0,確定大于5 km深度范圍內的一維速度結構參數:包括上地殼(5～17.64 km)、中地殼(17.64～26.54 km)、下地殼(26.54～36.18 km)以及地幔層(>36.18 km)的Vp、Vs和密度。在整個地下一維速度結構模型中,除Crust1.0直接給出的大于5 km深度范圍的Vp和密度之外,Vp、密度和品質因子的取值方式如下:Vp的值為1.6Vs[12],土層密度的取值參考文獻[13],S波和 P波的品質因子分別由Qs=Vs/10和Qp=Vs/5(Vs的單位為m/s)確定[14]。

1.2 GP14.3運動學混合震源模型的構建

本文采用能同時考慮破裂面上低波數確定性錯動和高波數隨機性錯動的GP14.3運動學混合震源模型,該模型的建立過程如下:

1)基于USGS發布的震源機制(見表2)和斷層滑動量分布結果(見圖3)[15],建立低波數確定性的凹凸體震源模型。震源參數包括全局震源參數、局部震源參數和其他震源參數[16],其中:全局震源參數表示斷層的形狀、大小和位置等,如斷層走向、傾角和滑動角,斷層面沿走向長度、沿傾向寬度,斷層頂面深度,斷層面平均滑動量,斷層面劃分的矩形子源尺寸等;局部震源參數表示震源內部的滑動不均勻性,如凹凸體的面積、凹凸體的滑動量等,這些參數由斷層的滑動量和分布結果來確定;其他震源參數表示破裂的形成和終止,如破裂起始點和破裂傳播方式,這些參數由破裂前端分布結果來確定,其中:破裂傳播方式為由破裂起始點向四周傳播。以上參數詳見表3。

表2 USGS發布的震源機制參數Table 2 Focal mechanism parameters given by USGS

表3 土耳其Mw7.8地震的震源參數Table 3 Local source parameters of the Turkey Mw7.8 earthquake

圖3 USGS所給的斷層破裂面上的滑動量及破裂前端分布(色譜表示滑動量幅值,五角星表示破裂起始點,灰色箭頭表示滑動矢量,等值虛線表示破裂起始時間(以秒為單位))Fig. 3 Distribution of the slip and rupture-front on the fault rupture surface released by USGS (Slip amplitude is shown in color. The point of rupture initiation is denoted by a star. The gray arrow represents the slip vector. The equivalent dotted lines show the rupture initiation time in seconds)

2)在低波數確定性的凹凸體震源模型中引入高波數隨機性錯動,確保合成寬頻地面運動的帶寬有效性,進而構建GP14.3運動學混合震源模型[17]。①首先將已有的凹凸體斷層面上滑動量的空間分布經二維傅里葉變換至波數域,得到斷層面上確定性的滑動波數譜。②其次采用波數衰減滿足von Karman自相關函數的波數譜并引入隨機數φ=[-π, π]表達斷層破裂面上小尺度的隨機變量。③然后在波數域中結合確定性的低波數譜和隨機高波數譜,并利用二維逆傅里葉變換至空間域中。

以上引入高波數隨機性錯動的過程中所涉及的Hurst指數H、控制結合的銳度N、沿斷層走向和傾向的拐角波數kcx與kcd均參考文獻[9]進行取值,進而確定斷層面上的滑動量分布,如圖4(a)所示。關于破裂前端的分布,首先確定破裂速度的初始分布,根據子源位置、初始破裂速度分布以及各子源滑動量對破裂起始時刻的時間擾動,得到破裂起始時刻的空間分布(圖4(a)等值線表示間隔5 s的破裂前端分布)。上升時間表征子源完成破裂過程所需時間,本文采用Hartzell矩率函數,得到各子源的上升時間后根據斷層面上的平均上升時間按比例調整各子源上升時間,可得其在斷層面上的空間分布如圖4(b)所示。

圖4 運動學混合震源Fig. 4 Hybrid kinematic source model

2 模型的適用性驗證分析

本節通過將地震動合成結果與8個臺站的觀測記錄對比,有效檢驗了寬頻地震動模擬的可靠性和模型參數(包括地下一維速度結構模型和混合震源模型)的適用性,且FK法的精度及高頻穩定性本文作者在文獻[18-19]中已廣泛驗證。首先,基于建立的GP14.3運動學混合震源模型,采用 FK 法計算了距離此次Mw7.8地震的震中50 km范圍內地震動較強的8個臺站的加速度時程(臺站位置信息見圖5和表4);其次,將地震動合成結果和AFAD提供的強震記錄及其反應譜進行對比[20]。該次計算的時間步距取為 0.01 s,時間步的數量為8192。這里需要說明的是:本文所選取的臺站記錄來源于AFAD網站截止至2月21日的強震記錄數據。

表4 所驗證的8個臺站的信息

圖6對比了本文計算的8個臺站的三分量加速度時程和強震記錄,黑線代表強震記錄,紅線代表本文合成結果,每條加速度時程的時長均為80 s,PGA均標注于曲線的右上方,單位為cm/s2。強震記錄及本文計算結果均經4階Butterworth 低通濾波至 10 Hz。從地震動時程比較可知:前6個為記錄完整的臺站,合成地震動在4616 臺站偏小,在其它5個臺站符合較好;NAR和4615臺站記錄雖然存在部分缺失,但因其記錄到了主要波形成分而仍具有一定的參考意義,本文合成的兩個水平方向的地震動與觀測記錄前段呈現的幅值和波形特征基本上一致,但豎向的地震動比觀測記錄明顯偏小;除4630臺站的EW方向、4616臺站的NS方向、NAR和4615臺站的UD方向的上的計算結果和強震記錄的PGA相差較大之外,其他的計算結果和強震記錄符合較好;2718、4629和4632臺站處的計算結果和強震記錄在強震段持時符合較好,但其他臺站處有一定差別。總體上,本文計算的加速度時程與強震記錄有一定的一致性。

圖7對比了本文合成結果與強震記錄的5%阻尼比加速度反應譜。為便于比較和觀察,將所有EW分量縮小5倍以及UD分量縮小25倍。由圖7可知:NAR和4615臺站處的豎向合成結果的短周期段小于觀測記錄,NS和EW方向的合成結果在0.1～10 s內均與觀測記錄符合較好;除NAR和4615臺站外的其他6個臺站的三個方向的合成結果在短周期段內均與觀測記錄符合較好;而2703、4629、4630和4632這4個臺站的三個方向的合成結果在長周期段內高于觀測記錄。綜上所述,從平均意義上來看:計算的地震動的短周期段與觀測記錄符合較好,而長周期段有一定程度的高估現象。根據上述對地震動時程、峰值和反應譜的比較分析表明:本文所構建的一維速度結構模型和混合震源模型適用于2023年土耳其Mw7.8地震的寬頻地震動模擬,可進一步用于近場特性的分析。

圖7 模擬地震動反應譜與強震記錄對比Fig. 7 Response-spectrum comparison between the simulations and the records

3 加濟安泰普及附近區域地面運動場合成及分析

本節對土耳其Mw7.8地震進行了三分量寬頻地震動合成。本節設定研究目標區域為加濟安泰普及附近區域,并在上一節研究基礎上合成了該區域內的地面加速度峰值、速度峰值以及烈度分布。首先,設定目標區域:如圖8所示,設定東西向180 km×南北向150 km的目標矩形區域(紅色矩形,面積為27000 km2),以3 km間隔均勻設置61×51,共計3111個計算點;其次,采用本文方法及模型合成了目標區域地面運動場;最后,給出了目標區域的三分量加速度和速度峰值分布,分析了其空間分布特征及近斷層效應;給出了目標區域的烈度分布圖,并與USGS給出的土耳其Mw7.8地震的ShakeMap進行對照;給出了部分位置的時程結果,分析此次地震造成的地面永久位移及近斷層速度脈沖特征。

圖8 模擬區域布置示意圖(黑色矩形表示斷層投影面、星號表示破裂起始點)Fig. 8 Sketch of the simulation area (The black rectangle indicates the projection plane of the fault and the star indicates the rupture initiation point)

3.1 地面運動峰值分布

圖9給出了目標區域地面加速度峰值(PGA)分布和速度峰值(PGV)分布圖,其中從左到右依次是NS、EW和豎向(UD)分量,從上到下依次是PGA和PGV。圖中矩形和星號分別表示斷層面的地表投影范圍和震中位置,藍色虛線區域為加濟安泰普市,綠色虛線區域分別為加濟安泰普市區區域和伊斯拉希耶鎮區域,黑色虛線表示主要道路。

圖9 地震動場PGA和PGV分布Fig. 9 Distribute of PGA and PGV

從圖9中可以看:加濟安泰普的市區范圍內的NS向PGA接近180 cm/s2,PGV接近100 cm/s,說明加濟安泰普市區范圍內的地震動強度較小。分析其原因為:走滑型斷層破裂釋放的能量主要集中在破裂前方區域,垂直斷層走向兩側的區域聚集能量較小,且遠離斷層投影區域后地震動幅值衰減較快,該市區雖距震中十分近,但其位置在這條西南-東北走向的東安納托利亞斷裂帶的東南方向,不屬于破裂前方區域,故該市區的地震動強度較小,這也與報道中提及的“加濟安泰普市區極少建筑出現崩塌或明顯損毀,大部分樓房從外部都看不出明顯損傷”的現象相符合[21]。由圖9還可知:位于加濟安泰普的西南部、直接處于地震斷層帶上的伊斯拉希耶鎮(37.03°N,36.63°E)地區內的EW向PGA高達560 cm/s2,PGV高達150 cm/s,這與附近2718臺站的觀測記錄(EW向PGA=584.2 cm/s2,EW向PGV=113.16 cm/s)較為接近,說明該地區的地震動強度較大。此外,從圖9的空間分布樣式來看:峰值較大區域集中在斷層投影面附近,這與震源模型的凹凸體位置相對應,體現了凹凸體對地面運動峰值的主導控制作用以及近斷層地震動的集中性效應。

3.2 烈度分布

在地震發生的初期,通常可利用烈度與強地面運動參數之間的轉換關系來快速獲取地震烈度分布特征,以快速可靠地估計烈度分布。本節基于強地面運動合成結果,應用WALD等[22]給出的加速峰值(PGA)、速度峰值(PGV)和烈度(IMM)的經驗關系式(1),給出震區的烈度分布圖。

(1)

式中:PGA單位為cm/s2,PGV的單位為cm/s,IMM為場地的烈度。基于上式即可將得到的地表水平PGA/PGV轉換為相應的地震烈度值。

圖10(a)和圖10(b)分別為USGS給出的土耳其Mw7.8地震的ShakeMap和本文基于強地面運動合成結果得到的地震烈度分布圖[23]。將圖10(a)范圍縮小至本文目標區域的大小,以便于本文結果與其對比。從圖10(b)中給出的烈度分布結果可看出:加濟安泰普市區范圍內的地震烈度為VII度,與圖10(a)中加濟安泰普市區范圍的烈度一致;處于西南部的伊斯拉希耶鎮地區的地震烈度高達IX度,與圖10(a)中伊斯拉希耶鎮地區的烈度一致,且與圖9所反映的伊斯拉希耶鎮地區的地震動強度較大的結論相符。

圖10 烈度分布Fig. 10 The intensity distribution

這里需要說明的是:本文所構建的GP14.3運動學混合震源模型沒有考慮彎折,為簡單連續的平面斷層,而USGS最終更新的震源模型是包含三段的彎折斷層,這也是造成本文一些地震動合成結果與臺站觀測記錄存在一定差距的主要原因。

3.3 時程結果分析

據報道,加濟安泰普一西部城鎮-努爾達伊附近出現幾百米長的地表斷裂,水平位移約3～4 m[24]。本節首先計算努爾達伊鎮中心(37.18°N,36.74°E)地表上一觀測點的時程結果,如圖11所示,其中從左到右依次是平行走向、垂直走向和豎向分量,從上到下依次是加速度、速度和位移時程。從圖11所合成的時程結果中可發現:在水平方向上有比較明顯的永久位移:平行于斷層走向的地表永久位移達到了105 cm,垂直于斷層走向的地表永久位移達到了100 cm,而豎向的地表永久位移約為10 cm。合成結果所反映的現象說明近斷層效應較為明顯,且這種現象可能導致建筑發生更加明顯的破壞。

圖11 Mw7.8土耳其地震中努爾達伊地區的三方向的加速度、速度和位移時程Fig. 11 Time histories of acceleration, velocity and displacement in three directions of Nurdagi region in Mw7.8 Turkey earthquake

圖12分別給出了在地表上沿斷層西南-東北走向以10 km的間隔均勻設置的15個觀測點的位置以及觀測點平行斷層走向、垂直斷層走向兩個方向的速度時程。如圖12(a)所示,R1觀測點位于破裂起始點在地表的投影處,R1～R15分別對應震中距為0、10、20、30、40、50、60、70、80、90、100、110、120、130和140 km的各個觀測點。圖12(b)為本文所合成的平行于斷層走向的速度時程結果,在速度時程中可看到單向速度脈沖(特別是R6～R9這四個觀測點);圖12(c)為本文所合成的垂直于斷層走向的速度時程結果,可看到明顯的雙向速度脈沖,這是因為此次地震的斷層為傾角較大的走滑斷層,方向性速度脈沖主要體現在垂直斷層走向的分量上[25],特別是位于破裂前端的R6～R8這三個觀測點,滿足方向性效應的條件,雙向速度脈沖現象更為明顯。由圖 12(b)和圖12(c)可以看出:在破裂起始點于地表的投影處(R1觀測點)速度脈沖現象不明顯,在距離斷層破裂起始點一定的距離后開始看到較為明顯的速度脈沖,比如距斷層破裂起始點30 km處的R4觀測點,速度脈沖已經很明顯。

圖12 震中距分別為0、10、20、30、40、50、60、70、80、90、100、110、120、130和140 km的各觀測點Fig. 12 Observation points with epicenter distance of 0、10、20、30、40、50、60、70、80、90、100、110、120、130 and 140 km

4 結論

針對2023年2月6日土耳其南部Mw7.8地震開展了近斷層寬頻地震動合成研究。綜合考慮Vs30數據及全球地殼模型Crust1.0,構建了土耳其東南部城市-加濟安泰普及附近地區的地下一維速度結構模型。依據USGS發布的信息建立了GP14.3運動學混合震源模型,進而采用作者提出的FK方法合成了8個臺站的三分量加速度時程,通過與臺站的時程記錄及其反應譜對比,驗證了本文方法的可靠性及模型的適用性。最后計算并分析了加濟安泰普附近27000 km2目標區域內的地震動場,得到如下結論:

1)本文合成的加濟安泰普市區范圍內的NS向PGA接近180 cm/s2,PGV接近100 cm/s;本文合成的位于加濟安泰普的西南部和直接處于地震斷層帶上的伊斯拉希耶鎮地區的EW向PGA高達560 cm/s2,PGV高達150 cm/s,與附近2718臺站的觀測記錄一致。

2)此次土耳其Mw7.8地震呈現明顯的近斷層地震動集中性效應、地面永久位移以及破裂方向性效應。

3)基于強地面運動合成結果得到的地震烈度分布顯示:加濟安泰普市區范圍內的地震烈度為VII度;處于西南部伊斯拉希耶鎮地區的地震烈度高達IX度。

4)從努爾達伊地區的地震動合成結果中可發現在水平方向上有比較明顯的永久位移。

由于本文震源模型未考慮實際彎折,同時采用的方法也不能考慮起伏地表的影響和地下結構的橫向非均勻性,這些因素可能是造成本文模擬結果與強震記錄在某些臺站差異還較大的主要原因。因此,在接下來的地震動模擬研究中,有必要對震源進行深入研究并考慮復雜地形的影響,以建立更為可靠的模型,進行更精細化的復雜場地地震動模擬,為地震區劃和抗震設防等工作提供參考。