基于等效群速度的非平穩地震動擬合方法研究

馬完君　盧建旗　唐麗華　謝志南

摘要：基于地震動能量傳播速度為群速度這一概念，回顧了相位微分和窄帶濾波兩種頻率相依波包時延計算方法及其在非平穩地震動相位擬合中的應用，歸納給出了一種非平穩地震動生成方法，通過模擬2條典型的海底與陸地地震動驗證了方法的有效性，并將該方法應用于區域非平穩地震動場模擬中。結果表明：利用該方法所得模擬地震動記錄平均加速度反應譜與原記錄加速度反應譜有較高的一致性，而且能很好地再現地震的非平穩特性。

關鍵詞：非平穩地震動;波包時遲;等效群速度;地震動模擬

中圖分類號：P315.914?文獻標志碼：A?文章編號：1000-0666（2019）04-0523-08

0?引言

加速度地震動是非平穩時間過程，不同頻率成分地震動組合得到窄帶波包所攜帶主要能量的到達時間是不一致的。一般而言，高頻成分的主要能量先到，低頻成分的主要能量后到。結構彈塑性地震反應分析結果表明，地震動非平穩特性對結構反應具有重要影響（Wang?et?al，2012;張郁山，趙鳳新，2014）。《超限高層建筑工程抗震設防專項審查技術要點》（中華人民共和國住房和城鄉建設部，2015）要求高度超過200?m的各類建筑結構需要進行大震彈塑性時程分析（王亞勇，2017）。因此，有必要研究可計入其非平穩特征的地震動模擬方法，為重要工程結構抗震設計以及大中型城市地震災害情景構建提供合理的地震動輸入。

地震動模擬方法有2類：地震學方法和工程方法。地震學方法一般應用于模擬區域地震動場，可大致區分為：確定性方法、隨機方法和混合方法。確定性方法受限于地下介質資料不足，目前大多用于1?Hz或以下頻帶范圍地震動的模擬（Roten?et?al，2011）。隨機方法產生的相位受限于高斯白噪聲是平穩隨機過程，模擬結果難以反映地震動的非平穩特性（Boore，2003，2009;Dabaghi?et?al，2018）。混合方法在如何實現高低頻地震動的合理銜接以保證地震動非平穩特征的一致性上仍需進一步研究（孫曉丹等，2011;方小丹等，2014）。工程方法一般應用于模擬工程結構場點位置處的地震動，但它一般不考慮震源以及地震動傳播的物理過程。為此，廖振鵬和魏穎（1988）提出了一種耦合工程方法與地震學方法的非平穩地震動模擬方法，基于能量傳播速度為群速度這一概念，建立了可用于統計區域范圍強地震動相位規律的等效頻散模型，進而提供了一個比相位差譜更容易理解的物理框架，以給出統計平均意義上符合區域地下介質波傳播特征的強地震動相位生成技術。

本文是廖振鵬和魏穎（1988）研究工作的延拓。首先回顧了窄帶波包能量行進速度——群速度與波包時延這2個概念，然后對比分析了已有的波包時延計算方法，采用Boore（2003）給出的相位微分法計算波包時延并建立等效群速度模型。基于該模型闡述了地震動非平穩相位的計算及地震動模擬方法。以典型非平穩地震動為例，通過分析模擬地震動和實際地震動對應的反應譜和RotD50殘差分布，驗證了該方法的合理性。在此基礎上，探討了其在區域非平穩地震動場模擬中的應用。

1?窄帶波包時延與群速度

1.1?波包時延的定義

基于傅里葉變換，可將地震動記錄分解為若干簡諧波的疊加。考慮震源附近地震動記錄的傅氏變換如下：

式中：Re[z]表示復數z的實部;exp表示自然對數為底的指數函數;ts（ω）表示震源附近處頻率為ω的地震分量起始時間。相應地，可將從震源附近輻射至距離R處地震動記錄表示為：

式中：HR（ω）涵蓋了幾何擴散與傳播介質對幅值的耗散等效應;R/cphase（ω）表示頻率為ω的地震動分量以相速度cphase（ω）傳播至R處所需要的時間;cphase（ω）=ω/k（ω），k（ω）表示空間波數。由于實際介質存在阻尼，介質中不同頻率成分波的相位傳播速度不一致。進一步考慮以ω0為中心頻率的窄頻帶（ω0-Δω，ω0+Δω）地震動分量的疊加，即：

假定在對Δω取值較小情形下，有：

式中：HR（ω0）表示HR（ω0）的絕對值;ζ=（t-R/cg）Δω;cg=dω/dkω=ω0表示群速度。取HR（ω0）=π，Δω=π/10，ω0=4π，k0=π，cg=3?km/s，R分別取為45?km和60?km，圖1給出上述窄頻帶地震動分量組合的傳播示意圖，從圖中可知，窄頻帶地震動分量組合表現為波包形式，其幅值包絡由（sinζ/ζ）函數所控制，由于這一地震動的能量主要包含在波包內部，可以看出波包傳播的速度為群速度即能量傳播速度。進一步分析式（5）可知，包絡函數峰值所在位置對應于t-R/cg=0處，給定地震動起始時間，將包絡函數峰值所在時間與地震動起始時間差值定義為包絡時延時間。

1.2?波包時延計算方法

從上述定義可知，最簡單的波包時延計算方法可通過對地震動記錄進行窄帶濾波得到（廖振鵬，魏穎，1988），濾波后記錄的峰值所在位置即為波包時延（方法1）。以1990年加利福尼亞州高地5.6級地震中的海底地震臺S3E水平加速度記錄（Boore，2003）為例，選取中心頻率f=ω/2?π=0.5，1.0，2.5，5.0，10.0?Hz對該記錄在頻域內用寬度為0.2?Hz的矩形窗濾波，然后返回至時域。圖2分別給出了S3E原記錄與窄帶濾波后所得記錄，結果表明：窄帶濾波后的記錄包含若干幅值差異較小甚至相當的波包，直接從濾波記錄中提取包絡時延存在不唯一的問題。另外這一方法需要使用多次傅里葉逆變換，極大地增加了計算量。針對這一問題，廖振鵬和金星（1995）首先建立了波包時延與記錄傅里葉變換后所得相位導數之間的關系（方法2），即：

式（8）相當于在局部范圍內對相位進行了解卷繞，然而相位解卷繞缺乏嚴格的數學基礎（Shatilo，1992）。Boore（2003）建議了一種基于相位導數的波包時延計算方法（方法3）：

式中：Im[z]表示復數z的虛部;H（f）為原記錄h（t）的傅里葉變換;G（f）為原記錄與時間t的乘積所得記錄t×h（t）的傅里葉變換。

針對S3E記錄，圖3給出了上述3種方法的計算結果，從圖中可以看出，3種方法計算所得結果均反映了S3E記錄中呈現的高頻地震動分量先到、低頻地震動分量后到的特征。第一種方法由于選取中心頻率點較少，計算結果較為稀疏，后2種方法所得計算結果大致相同。但由于Boore（2003）的方法避免了相位解卷繞以及多次傅里葉逆變換，因此本文利用這一方法提取記錄波包時延。

2?等效群速度模型與地震動擬合

基于波包時延的概念，在震源輻射出來的波包中，將以最快速度到達記錄臺站的波包所經歷的時間t0作為地震動起始時間，得到不同中心頻率成分地震動波包行進的等效群速度為：

式中：R是震中距;t0=R/Um，Um是等效群速度的最大值;t（f）為相對于t0的波包時延。考慮到地震動受震源物理過程、波在深部地下介質和淺層場地中傳播過程的綜合影響，且一般而言傳播距離不等于震中距，稱U（f）為等效群速度。因此基于地震動記錄提取等效群速度，做如下處理：首先對地震動記錄進行基線校正和高通濾波處理（Boore，2003），因為極低頻地震動記錄有可能是失真的。然后根據99%能量持時提取記錄的強震段td=t2-t1，即根據I（t1）=0.001，I（t2）=0.991來確定，其中能量積分為：

式中：t1=t0。再由式（9）計算波包時延，這樣保證波包時延計算用到的記錄起始時間為最快波包到達臺站的時間，然后利用式（10）計算等效群速度。考慮到地震動傳播具有一定的隨機性，可將等效群速度表示為：

式中：U（f）表示等效群速度隨頻率變化的平均趨勢;η是一個具有零均值的隨機過程。考慮如下U（f）擬合公式（廖振鵬，金星，1995）：

式中：a0，a1，a2系數可通過最小二乘或其他優化方法求解。一旦求得U（f），η的樣本η便可以確定，但是得到的樣本均值可能不是零，所以要對η進行微調：

同時還需要調整式（13）的求得的系數，即ai=ai（1+ε/fm），最后基于修正后η（f）樣本確定η這一均勻分布的取值范圍，本文取為95%置信區間-1.96S[η（f）]，1.96S[η（f）]，S[η（f）]為η（f）樣本的標準差。圖4中給出了基于S3E記錄求取對應等效群速度模型的流程圖，最終得到的參數分別為a0=1.110?9，a1=3.332?1，a2=

-1.002?9，S[η（f）]=0.285?7。基于群速度模型，建立地震動記錄的相位信息，進而給出模擬地震動。具體步驟包括：①基于式（12）計算出等效群速度樣本;②利用cphase（ω）=f∫f0[1/U（f）]df，計算得到相速度，在此基礎上計算相位φ（f）=-2πfR/cphase（f）+2πfR/Um;③結合基于模擬或實際記錄得到的地震動幅值譜，基于傅里葉逆變換給出模擬地震動。

3?基于等效群速度模型的地震動模擬

以上述S3E海底地震動以及1986年11月14日7級地震為例，利用震中距為79?km的Smart-1中心臺站三分量觀測地震動，基于上述方法給出符合記錄非平穩特征的模擬地震動。結果表明：模擬地震動在統計平均意義上與實測記錄的非平穩特征基本一致。最后討論了該方法在區域地震動場模擬中的應用。

3.1?典型海域與陸域場點非平穩地震動模擬

在上文給出的S3E記錄的一個相位擬合樣本的基礎上，在-1.96S[η（f）]，1.96S[η（f）]區間對η（f）進行隨機采樣，得到50組不同的相位樣本，然后將相位與原記錄傅里葉幅值譜結合，基于傅里葉逆變換得到模擬地震動。圖5給出了S3E原記錄與50組模擬地震動中的一個樣本，以及50組模擬地震動記錄計算所得加速度反應譜的平均值與原記錄加速度反應譜的對比圖，結果表明：模擬地震動記錄不僅在時間過程上較好地反映了原記錄的非平穩特性，且模擬地震動記錄所得加速度反應譜的平均值與原記錄加速度反應譜基本吻合。與之類似，本文將上述方法用于擬合Smart-1中心臺地震動三分量觀測記錄，對每一個地震分量給出50次模擬結果。表1給出對應于三分量記錄的群速度模型，地震動模擬結果如圖6所示。此外本文還將水平兩方向的50次模擬地震動進行了合成，求出其RotD50與原始的記錄的

RotD50的殘差（Boore，2010），結果如圖7所示，本文建立的地震動模擬方法在保證反應譜擬合精度的前提下，很好地再現了地震的非平穩特性。

4?結論與討論

本文回顧了相位微分和窄帶濾波2種頻率相依波包時延計算方法，闡明了窄帶波包能量行進速度——群速度與波包時延二者之間的關系，基于波包時延建立等效群速度模型及其在非平穩地震動相位擬合中的應用，并建議根據Boore（2003）給出的相位微分法計算波包時延。在此基礎上，建立了一種非平穩地震動生成方法，討論了該方法在場點地震動模擬中的應用，利用2條典型的海底與陸地地震動記錄模擬驗證了這一方法的有效性，給出了該方法在區域非平穩地震動場模擬中的推廣措施。該方法可進一步應用于歷史地震與設定地震作用下非平穩地震動場的模擬，為重要結構彈塑性地震響應分析、大型城市地震災害評估與預測提供合理的地震動輸入。后續研究工作還需深入分析斷層破裂效應以及不同類別場地或特殊場地如盆地效應對地震動場的影響。

本文的研究工作與Kumari等（2018）、廖振鵬和魏穎（1988）的結論基本一致，不同的是本文僅考慮點源，并且等效群速度模型計入了距離的影響，今后的研究有必要將二者結合，考慮震源在斷層面上的破裂過程，還可進一步探討剔除震源與場地影響的等效群速度模型構建方法。

參考文獻：

方小丹，魏璉，周靖.2014.長周期結構地震反應的特點與反應譜[J].建筑結構學報，35（3）：6-23.

廖振鵬，金星.1995.強地震動相位的一個隨機模型[J].地震學報，17（3）：353-361.

廖振鵬，魏穎.1988.設計地震加速度圖的合成[J].地震工程與工程振動，8（1）：12-28.

孫曉丹，陶夏新，劉陶鈞.2011.近場地震動估計中高，低頻地震動的疊加[J].自然災害學報，20（4）：84-89.

王亞勇.2017.結構時程分析輸入地震動準則和輸出結果解讀[J].建筑結構，47（11）：6-11.

張郁山，趙鳳新.2014.基于希爾伯特變換的非平穩地震動模擬方法的驗證[J].地震學報，36（4）：686-697.

Boore?D?M.2003.Phase?derivatives?and?simulation?of?strong?ground?motions[J].Bulletin?of?the?Seismological?Society?of?America，93（3）：1132-1143.

Boore?D?M.2009.Comparing?stochastic?point-source?and?finite-source?ground-motion?simulations：SMSIM?and?EXSIM[J].Bulletin?of?the?Seismological?Society?of?America，99（6）：3202-3216.

Boore?D?M.2010.Orientation-independent，nongeometric-mean?measures?of?seismic?intensity?from?two?horizontal?components?of?motion[J].Bulletin?of?the?Seismological?Society?of?America，100（4）：1830-1835.

Dabaghi?M，Der?Kiureghian?A.2018.Simulation?of?orthogonal?horizontal?components?of?near‐fault?ground?motion?for?specified?earthquake?source?and?site?characteristics.Earthquake?Engineering?&?Structural?Dynamics，47（6）：1369-1393.

Kumari?N，Gupta?I?D，Sharma?M?L.2018.Synthesizing?nonstationary?earthquake?ground?motion?via?empirically?simulated?equivalent?group?velocity?dispersion?curves?for?Western?Himalayan?Region[J].Bulletin?of?the?Seismological?Society?of?America，108（6）：3469-3487.

Razafindrakoto?H?N，Bradley?B?A，Graves?R?W.2018.Broadband?Ground‐Motion?Simulation?of?the?2011?MW6.2?Christchurch，New?Zealand，Earthquake[J].Bulletin?of?the?Seismological?Society?of?America，108（4）：2130-2147.

Roten?D，Olsen?K?B，Pechmann?J?C，et?al.2011.3D?simulations?of?M7?earthquakes?on?the?Wasatch?fault，Utah，Part?I：Long-period（0-1?Hz）ground?motion[J].Bulletin?of?the?Seismological?Society?of?America，101（5）：2045-2063.

Shatilo?A?P.1992.Seismic?phase?unwrapping：methods，results，problems[J].Geophysical?prospecting，40（2）：211-225.

Wang?J，Fan?L，Qian?S，et?al.2002.Simulations?of?nonstationary?frequency?content?and?its?importance?to?seismic?assessment?of?structures[J].Earthquake?Engineering?&?Structural?Dynamics，31（4）：993-1005.

中華人民共和國住房和城鄉建設部.2015.超限高層建筑工程抗震設防專項審查技術要點[S].

A?Method?for?Non-stationary?Ground?Motion?Simulation?Basedon?the?Equivalent?Group?Velocity?Model

MA?Wanjun1，LU?Jianqi1，TANG?Lihua2，XIE?Zhinan1

（1.Key?Laboratory?of?Earthquake?Engineering?and?Engineering?Vibration，Institute?of?Engineering?Mechanics，China?Earthquake?Administration，Harbin?150080，Heilongjiang，?China）（2.Earthquake?Agency?of?Xinjiang?Uygur?Autonomous?Region，Urumqi，830011，Xinjiang，?China）

Abstract

Based?on?the?concept?that?the?energy?propagation?velocity?of?ground?motion?is?group?velocity，this?paper?briefly?reviews?the?phase?differential?and?narrow-band?filtering?methods?for?calculating?frequency-dependent?wave?packet?time?delay?and?their?applications?in?non-stationary?ground?motion?phase?simulation.Then?a?non-stationary?ground?motion?simulation?method?is?established.The?effectiveness?of?this?method?is?verified?by?simulating?two?typical?ocean-bottom?and?land?ground?motions.Finally，the?application?of?this?method?in?regional?non-stationary?ground?motion?field?simulation?is?discussed.The?results?show?that?the?average?acceleration?response?spectrum?of?the?simulated?ground?motion?recorded?by?this?method?is?consistent?with?the?original?recorded?acceleration?response?spectrum.In?addition，the?non-stationary?characteristics?of?the?earthquake?can?be?well?reproduced.

Keywords：nonstationary?ground?motion;envelope?delay;equivalent?group?velocity;ground?motion?simulation