基于等效群速度的非平穩(wěn)地震動(dòng)擬合方法研究

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

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

關(guān)鍵詞：非平穩(wěn)地震動(dòng);波包時(shí)遲;等效群速度;地震動(dòng)模擬

中圖分類號：P315.914?文獻(xiàn)標(biāo)志碼：A?文章編號：1000-0666（2019）04-0523-08

0?引言

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

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

本文是廖振鵬和魏穎（1988）研究工作的延拓。首先回顧了窄帶波包能量行進(jìn)速度——群速度與波包時(shí)延這2個(gè)概念，然后對比分析了已有的波包時(shí)延計(jì)算方法，采用Boore（2003）給出的相位微分法計(jì)算波包時(shí)延并建立等效群速度模型?；谠撃Ｐ完U述了地震動(dòng)非平穩(wěn)相位的計(jì)算及地震動(dòng)模擬方法。以典型非平穩(wěn)地震動(dòng)為例，通過分析模擬地震動(dòng)和實(shí)際地震動(dòng)對應(yīng)的反應(yīng)譜和RotD50殘差分布，驗(yàn)證了該方法的合理性。在此基礎(chǔ)上，探討了其在區(qū)域非平穩(wěn)地震動(dòng)場模擬中的應(yīng)用。

1?窄帶波包時(shí)延與群速度

1.1?波包時(shí)延的定義

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

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

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

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

式中：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給出上述窄頻帶地震動(dòng)分量組合的傳播示意圖，從圖中可知，窄頻帶地震動(dòng)分量組合表現(xiàn)為波包形式，其幅值包絡(luò)由（sinζ/ζ）函數(shù)所控制，由于這一地震動(dòng)的能量主要包含在波包內(nèi)部，可以看出波包傳播的速度為群速度即能量傳播速度。進(jìn)一步分析式（5）可知，包絡(luò)函數(shù)峰值所在位置對應(yīng)于t-R/cg=0處，給定地震動(dòng)起始時(shí)間，將包絡(luò)函數(shù)峰值所在時(shí)間與地震動(dòng)起始時(shí)間差值定義為包絡(luò)時(shí)延時(shí)間。

1.2?波包時(shí)延計(jì)算方法

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

式（8）相當(dāng)于在局部范圍內(nèi)對相位進(jìn)行了解卷繞，然而相位解卷繞缺乏嚴(yán)格的數(shù)學(xué)基礎(chǔ)（Shatilo，1992）。Boore（2003）建議了一種基于相位導(dǎo)數(shù)的波包時(shí)延計(jì)算方法（方法3）：

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

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

2?等效群速度模型與地震動(dòng)擬合

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

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

式中：t1=t0。再由式（9）計(jì)算波包時(shí)延，這樣保證波包時(shí)延計(jì)算用到的記錄起始時(shí)間為最快波包到達(dá)臺站的時(shí)間，然后利用式（10）計(jì)算等效群速度?？紤]到地震動(dòng)傳播具有一定的隨機(jī)性，可將等效群速度表示為：

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

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

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

-1.002?9，S[η（f）]=0.285?7?；谌核俣饶Ｐ?，建立地震動(dòng)記錄的相位信息，進(jìn)而給出模擬地震動(dòng)。具體步驟包括：①基于式（12）計(jì)算出等效群速度樣本;②利用cphase（ω）=f∫f0[1/U（f）]df，計(jì)算得到相速度，在此基礎(chǔ)上計(jì)算相位φ（f）=-2πfR/cphase（f）+2πfR/Um;③結(jié)合基于模擬或?qū)嶋H記錄得到的地震動(dòng)幅值譜，基于傅里葉逆變換給出模擬地震動(dòng)。

3?基于等效群速度模型的地震動(dòng)模擬

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

3.1?典型海域與陸域場點(diǎn)非平穩(wěn)地震動(dòng)模擬

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

RotD50的殘差（Boore，2010），結(jié)果如圖7所示，本文建立的地震動(dòng)模擬方法在保證反應(yīng)譜擬合精度的前提下，很好地再現(xiàn)了地震的非平穩(wěn)特性。

4?結(jié)論與討論

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

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

參考文獻(xiàn)：

方小丹，魏璉，周靖.2014.長周期結(jié)構(gòu)地震反應(yīng)的特點(diǎn)與反應(yīng)譜[J].建筑結(jié)構(gòu)學(xué)報(bào)，35（3）：6-23.

廖振鵬，金星.1995.強(qiáng)地震動(dòng)相位的一個(gè)隨機(jī)模型[J].地震學(xué)報(bào)，17（3）：353-361.

廖振鵬，魏穎.1988.設(shè)計(jì)地震加速度圖的合成[J].地震工程與工程振動(dòng)，8（1）：12-28.

孫曉丹，陶夏新，劉陶鈞.2011.近場地震動(dòng)估計(jì)中高，低頻地震動(dòng)的疊加[J].自然災(zāi)害學(xué)報(bào)，20（4）：84-89.

王亞勇.2017.結(jié)構(gòu)時(shí)程分析輸入地震動(dòng)準(zhǔn)則和輸出結(jié)果解讀[J].建筑結(jié)構(gòu)，47（11）：6-11.

張郁山，趙鳳新.2014.基于希爾伯特變換的非平穩(wěn)地震動(dòng)模擬方法的驗(yàn)證[J].地震學(xué)報(bào)，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，F(xiàn)an?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.

中華人民共和國住房和城鄉(xiāng)建設(shè)部.2015.超限高層建筑工程抗震設(shè)防專項(xiàng)審查技術(shù)要點(diǎn)[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