姜云木，阮鑫鑫，劉章軍

(武漢工程大學(xué) 土木工程與建筑學(xué)院，武漢 430074)

地震災(zāi)害具有很強(qiáng)的破壞性與隨機(jī)性，同時(shí)在一次主震后往往會(huì)伴隨著多次余震的發(fā)生，且余震對(duì)結(jié)構(gòu)造成的二次破壞不容忽視。由于主余震序列地震動(dòng)記錄數(shù)量有限以及對(duì)具體場地條件的限制，無法滿足對(duì)結(jié)構(gòu)計(jì)算的需求，因此，主余震型地震動(dòng)的人工模擬備受關(guān)注。

建立合理的主余震地震動(dòng)參數(shù)之間的關(guān)系，是精細(xì)模擬主余震型地震動(dòng)的前提。相對(duì)于主余震復(fù)雜的震源機(jī)制等因素，工程師更關(guān)注主余震地震動(dòng)的幅值、持時(shí)與頻譜等工程特性的關(guān)聯(lián)。在早期的研究中，學(xué)者們主要對(duì)主余震之間的峰值加速度、持時(shí)等參數(shù)關(guān)系進(jìn)行了初步探索[1-3]。近年來，Zhang等[4]利用線性合成法擬合了主余震間震級(jí)、持時(shí)、能量等參數(shù)的關(guān)系。隨后，Muderrisoglu等[5]估計(jì)了主震后一段時(shí)間內(nèi)發(fā)生超越特定震級(jí)閾值的余震概率。此外，朱瑞廣等[6]通過Copula函數(shù)較為完整地分析了主余震強(qiáng)度參數(shù)的相關(guān)性。Ruiz-Garcia等[7]研究發(fā)現(xiàn)：主震的場地卓越圓頻率比余震的場地卓越圓頻率有偏小的趨勢，且具有中等的相關(guān)性。上述研究表明，主余震參數(shù)之間存在一定的相關(guān)性，但從工程實(shí)用角度來看仍需進(jìn)一步的研究。

主余震型地震動(dòng)的構(gòu)造方法可分為確定方法和隨機(jī)方法兩大類。在確定方法方面：Hatzigeorgiou等提出了重復(fù)法來構(gòu)造主余震序列，該方法假設(shè)主余震的特性一致，這顯然不符合實(shí)際情況；Li等[8]提出了隨機(jī)組合法構(gòu)造主余震序列，即從主余震記錄庫中分別隨機(jī)挑選出一條地震動(dòng)記錄組合成主余震序列，但該方法本質(zhì)上得到的主余震序列不具有隨機(jī)性。在隨機(jī)方法方面：Hu等[9]首先利用調(diào)幅過濾白噪聲方法生成主震，再利用分支序列法生成余震，但是該方法并沒有考慮主余震間的條件關(guān)系與衰減關(guān)系；Nithin等[10]基于Priestley理論，提出了一種基于主震的余震條件縮放模型，并結(jié)合特定的場景生成了主余震時(shí)程樣本。然而，該方法在本質(zhì)上屬于蒙特卡洛方法，為了保持模擬精度，需要成千上萬的隨機(jī)變量，導(dǎo)致了模擬時(shí)生成的樣本數(shù)量龐大以及概率信息不完備的問題。綜上所述，亟需建立一種能考慮主震與余震相互關(guān)聯(lián)且高效實(shí)用的模擬方法。

基于上述研究進(jìn)展，筆者首先在非平穩(wěn)地震動(dòng)過程的演變功率譜密度模型的基礎(chǔ)上，對(duì)實(shí)測主余震記錄的峰值加速度、場地土的卓越圓頻率、阻尼比以及調(diào)制函數(shù)參數(shù)進(jìn)行識(shí)別。然后，通過擬合優(yōu)度檢驗(yàn)，給出最優(yōu)邊緣分布和最優(yōu)Copula模型。在此基礎(chǔ)上，得到了主震參數(shù)條件下對(duì)應(yīng)余震參數(shù)的條件均值。最后，結(jié)合譜表示-隨機(jī)函數(shù)方法[11]，建立了主余震型地震動(dòng)的降維模型，實(shí)現(xiàn)了主余震型地震動(dòng)的高效模擬。此外，本文建立了主余震參數(shù)之間的實(shí)用計(jì)算公式，方便了工程應(yīng)用。

1 時(shí)-頻全非平穩(wěn)模型參數(shù)的識(shí)別

1.1 非平穩(wěn)地震動(dòng)過程的演變功率譜密度模型

相對(duì)于震級(jí)和震源機(jī)制的地震學(xué)因素和“震源-傳播途徑-局部場地”的物理模型，地震工程更加關(guān)注地震動(dòng)本身的幅值、持時(shí)與頻譜等工程特性。因此，本文選用工程中經(jīng)典的演變功率譜密度模型，該模型從地震動(dòng)的頻譜特性出發(fā)，能夠反映地震動(dòng)過程的非平穩(wěn)性和統(tǒng)計(jì)特性。同時(shí)，通過建立主余震演變功率譜參數(shù)之間關(guān)系，可以方便地體現(xiàn)地震動(dòng)的序列性。

根據(jù)非平穩(wěn)隨機(jī)過程的Priestley演變譜理論，非平穩(wěn)地震動(dòng)加速度過程的演變功率譜密度函數(shù)可表示為[12]

(1)

式中:SUg(t,ω)為非平穩(wěn)地震動(dòng)加速度過程Ug(t)的單邊演變功率譜密度函數(shù)；A(t,ω)為時(shí)-頻調(diào)制函數(shù)；S(ω)為平穩(wěn)地震動(dòng)加速度過程的單邊功率譜密度函數(shù)。

對(duì)于時(shí)-頻調(diào)制函數(shù)，采用Deodatis等[13]提出，劉章軍等[14]改進(jìn)的模型

其中,

(3)

式中:參數(shù)b=a+0.001;c=0.005；a為控制地震動(dòng)過程衰減快慢的參數(shù)，s-1。一般地，a越大，地震動(dòng)過程衰減越快。

對(duì)于平穩(wěn)地震動(dòng)加速度過程的單邊功率譜，采用Clough-Penzien模型[15]

(4)

式中:ωg和ξg分別為場地土的卓越圓頻率和阻尼比；ωf和ξf分別為基巖的卓越圓頻率和阻尼比。一般地，可取ωf=0.1ωg，ξf=ξg。

在式(4)中，譜強(qiáng)度因子S0可表示為

(5)

式中:Amax為地震動(dòng)峰值加速度的均值；r為峰值因子，在進(jìn)行參數(shù)識(shí)別時(shí)取r=3。

1.2 演變功率譜密度模型的參數(shù)識(shí)別

由1.1節(jié)提出的時(shí)-頻全非平穩(wěn)模型可知，演變功率譜密度模型的參數(shù)向量可以表示為λ，即λ=[ωg,ξg,a]。需要指出的是，在進(jìn)行場地土的卓越圓頻率和阻尼比識(shí)別時(shí)，對(duì)實(shí)測強(qiáng)震加速度記錄進(jìn)行了調(diào)幅，因此參數(shù)向量λ中不包含地震動(dòng)峰值加速度均值A(chǔ)max和峰值因子r。于是，可定義演變功率譜密度模型的頻域能量分布函數(shù)E(ω;λ)如下[16]

(6)

將式(1)和式(2)代入式(6)，積分可得本文模型對(duì)應(yīng)的頻域能量分布函數(shù)

E(ω;λ)=

(7)

記第i條實(shí)測記錄的演變功率譜為Sug,i(t,ω)，如果放松演變功率譜的能量隨時(shí)間分布再取其平均，那么對(duì)于第i條實(shí)測記錄，可以得到其等價(jià)平穩(wěn)過程的功率譜的函數(shù)Si(ω)如下[17]

(8)

由式(8)可得，第i條實(shí)測記錄的頻域能量分布函數(shù)Ei(ω)為

(9)

式中：Ti為第i條地震動(dòng)記錄的總持續(xù)時(shí)間；T0,i為第i條地震動(dòng)記錄的強(qiáng)震持時(shí)，具體可表示為[18]

(10)

(11)

(12)

式中:ai(t)為第i條實(shí)測地震記錄時(shí)程；ωu,i為第i條實(shí)測地震記錄時(shí)程的上限頻率。

當(dāng)取得一條地震動(dòng)實(shí)測記錄時(shí)，可以通過信號(hào)處理的方法得到實(shí)測記錄的平穩(wěn)功率譜Si(ω)，再根據(jù)式(9)與式(10)即可得到第i條實(shí)測記錄的頻域能量分布函數(shù)Ei(ω)。

根據(jù)主余震實(shí)測記錄得到的頻域能量分布函數(shù)曲線，采用最佳平方逼近準(zhǔn)則，針對(duì)每條曲線分別識(shí)別參數(shù)向量λ

(13)

本識(shí)別方法原理簡單，利用簡單的頻域能量相等原理，避免了在時(shí)頻域上參數(shù)識(shí)別的困難，減少了計(jì)算量。且編程簡單，僅需要3個(gè)主要的MATLAB工具箱函數(shù)便可以實(shí)現(xiàn)。

1.3 實(shí)測地震動(dòng)記錄的識(shí)別結(jié)果

本文在太平洋地震動(dòng)工程研究中心的NGA-West2地震動(dòng)數(shù)據(jù)庫中嚴(yán)格按照以下原則[19]篩選地震波：

(1)主震記錄與其對(duì)應(yīng)的余震記錄必須來自同一臺(tái)站。

(2)余震只選取與主震對(duì)應(yīng)且震級(jí)最大的余震作為研究對(duì)象。

(3)斷層距離應(yīng)大于10 km，以減少近斷層效應(yīng)的影響。

(4)主余震的震級(jí)應(yīng)該均大于5，排除對(duì)結(jié)構(gòu)影響不大的地震動(dòng)。

經(jīng)過嚴(yán)格篩選得到了17次地震的468組地震動(dòng)記錄。表1給出了本文所采用17次地震的地震名稱、臺(tái)站數(shù)量、震級(jí)、斷層距范圍和VS,30范圍等基本信息。

表1 本文選取地震動(dòng)記錄的詳細(xì)信息

本文根據(jù)GB 18306—2015《中國地震動(dòng)參數(shù)區(qū)劃圖》[20]中建議的5種場地類別對(duì)實(shí)測強(qiáng)震記錄進(jìn)行篩選。為了將國外的實(shí)測記錄與國內(nèi)的場地類別相結(jié)合，眾多學(xué)者研究建立了《中國地震動(dòng)參數(shù)區(qū)劃圖》中場地類別與剪切波速VS,30的對(duì)應(yīng)關(guān)系[21]。在本文中，采用文獻(xiàn)[22]的場地劃分標(biāo)準(zhǔn)。表2給出了不同場地類別對(duì)應(yīng)的剪切波速VS,30范圍與篩選得到的實(shí)測強(qiáng)震記錄的數(shù)量。

表2 場地類別與VS,30的對(duì)應(yīng)關(guān)系

對(duì)468組主余震加速度記錄進(jìn)行基線校正，并將峰值加速度調(diào)幅到200 cm/s2，利用1.2節(jié)的識(shí)別方法按照?qǐng)龅仡悇e依次對(duì)演變功率譜密度模型的參數(shù)向量λ進(jìn)行識(shí)別。

以CHY057臺(tái)站分別觀測的Chi-Chi地震主震和余震為典型實(shí)例。圖1(a)和圖2(a)分別給出了典型實(shí)例的主余震的實(shí)際能量曲線與模型能量曲線的對(duì)比，由圖1(a)和圖2(a)可知，模型能量分布曲線與實(shí)測記錄的能量分布曲線擬合效果較好。同時(shí)，圖1(b)、圖1(c)、圖2(b)、圖2(c)給出了本文選取所有的II類場地主余震實(shí)測能量曲線與模型能量曲線的均值和標(biāo)準(zhǔn)差對(duì)比?？梢钥闯觯饔嗾鹉Ｐ湍芰壳€的均值及標(biāo)準(zhǔn)差與實(shí)測記錄對(duì)比效果良好，充分說明了本文參數(shù)識(shí)別方法的有效性。

圖2 余震能量擬合曲線

2 主余震參數(shù)的最優(yōu)邊緣概率分布

在本文中，除了對(duì)主余震記錄的場地參數(shù)，即參數(shù)向量λ進(jìn)行識(shí)別外，還提取主余震地震動(dòng)實(shí)測記錄的峰值加速度PGA。這樣，考慮主震和余震，一共需要估計(jì)8個(gè)地震動(dòng)參數(shù)的邊緣概率密度分布。需要指出的是，提取主余震的峰值加速度未考慮場地類別，且均來自468組地震動(dòng)記錄調(diào)幅之前的統(tǒng)計(jì)結(jié)果。

圖3 主余震地震動(dòng)參數(shù)概率密度函數(shù)

在主余震參數(shù)的邊緣分布參數(shù)識(shí)別完成后，需要采用AIC信息準(zhǔn)則(akaike information criterion)判斷最優(yōu)的邊緣分布模型。AIC值的定義如下[24]

AIC=2k-2lnL

(14)

式中:k為分布模型中參數(shù)的個(gè)數(shù);L為模型最大似然函數(shù)。一般地，AIC值越小，代表模型擬合度越好。

通過AIC信息準(zhǔn)則判斷，并可得到主余震參數(shù)的最優(yōu)邊緣分布以及模型參數(shù),如表3所示。當(dāng)概率模型為對(duì)數(shù)正態(tài)分布時(shí)，Par1為均值，Par2為標(biāo)準(zhǔn)差；為Gumbel分布時(shí)，Par1為位置參數(shù)，Par2為形狀參數(shù)；為廣義極值分布時(shí)，Par1為形狀參數(shù)，Par2為尺度參數(shù)，Par3為位置參數(shù)；為Weibull分布時(shí)，Par1為比例參數(shù)，Par2為形狀參數(shù)。

表3 主余震時(shí)域參數(shù)與場地參數(shù)概率模型

3 基于Copula理論的主余震參數(shù)相關(guān)性分析

3.1 基本理論

Copula函數(shù)可以看作是一類聯(lián)結(jié)聯(lián)合分布函數(shù)和邊緣分布函數(shù)的“紐帶”，它描述了多維隨機(jī)變量之間的相關(guān)性。根據(jù)Sklar定理[25]，二維隨機(jī)變量(X1,X2)的聯(lián)合分布函數(shù)F(x1,x2)及其邊緣分布可寫成如下形式

F(x1,x2)=C(F1(x1),F2(x2))=C(u1,u2)

(15)

式中:FXi(xi)(i=1,2)為隨機(jī)變量Xi的邊緣概率分布函數(shù);C(·)為Copula分布函數(shù)。于是，二維隨機(jī)變量(X1,X2)的聯(lián)合概率密度函數(shù)可以表達(dá)為

f(x1,x2)=c(u1,u2)fX1(x1)fX2(x2)

(16)

式中：fXi(xi)(i=1,2)為隨機(jī)變量Xi的邊緣概率密度函數(shù)；c(·)為Copula密度函數(shù)，即

(17)

若將主震的某一參數(shù)(例如，場地土的卓越圓頻率)視為隨機(jī)變量X1，與其對(duì)應(yīng)的余震參數(shù)視為隨機(jī)變量X2，則在給定主震參數(shù)的條件下對(duì)應(yīng)余震參數(shù)的概率密度函數(shù)可表達(dá)為

鼓粒成熟期是大豆積累干物質(zhì)最多的時(shí)期，也是產(chǎn)量形成的重要時(shí)期。促進(jìn)養(yǎng)分向子粒中轉(zhuǎn)移，促粒飽增粒重，適期早熟則是這個(gè)時(shí)期管理的中心。這個(gè)時(shí)期缺水會(huì)使秕莢、秕粒增多，百粒重下降。秋季遇旱無雨，應(yīng)及時(shí)澆水，以水攻粒對(duì)提高產(chǎn)量和品質(zhì)有明顯影響。大豆黃熟末期為適收期。當(dāng)種子含水量達(dá)到13%時(shí)可以入庫。

(18)

進(jìn)一步地，條件均值可以表達(dá)為

(19)

從式(19)可知，在確定了主余震參數(shù)的邊緣分布以及Copula函數(shù)之后，通過給定主震參數(shù)的取值，就可以得到該條件下的對(duì)應(yīng)余震參數(shù)的取值。

3.2 最優(yōu)Copula函數(shù)的選擇

在選擇主余震參數(shù)的最優(yōu)Copula函數(shù)之前，需采用Kendall秩相關(guān)系數(shù)判斷參數(shù)相關(guān)性的大小。Kendall秩相關(guān)系數(shù)τ可由式(20)計(jì)算[26]

(20)

式中：x1，i和x2,i分別為主震和余震地震動(dòng)參數(shù)的第i個(gè)觀測值；N為樣本容量;sign[·]為符號(hào)函數(shù)，當(dāng)(x1,i-x1,j)(x2,i-x2,j)>0時(shí)，sign=1，否則sign=0。

主余震參數(shù)的Kendall相關(guān)系數(shù),如表4所示。需要說明的是，由于Ⅰ0和Ⅳ類場地的實(shí)測記錄數(shù)量較少，故沒有分析其相關(guān)性。從表4中可以看出主余震的PGA之間保持著中等相關(guān)性，場地參數(shù)除了場地阻尼比之外保持著中等相關(guān)性。盡管主余震場地阻尼比ξg之間趨近于獨(dú)立，但是筆者同樣發(fā)現(xiàn)不同場地余震的場地土卓越圓頻率ωg,aft和阻尼比ξg,aft之間均保持著中等負(fù)相關(guān)性(見表4)，因此可以采用Copula函數(shù)分析ωg,aft與ξg,aft的相關(guān)性。

表4 主余震參數(shù)之間的Kendall秩相關(guān)系數(shù)

根據(jù)主余震的參數(shù)之間的相關(guān)特性，本文采用6種Copula函數(shù)對(duì)主余震對(duì)應(yīng)的峰值加速度、場地土卓越圓頻率和調(diào)制函數(shù)參數(shù)以及余震的場地土卓越圓頻率與阻尼比之間的相關(guān)結(jié)構(gòu)進(jìn)行擬合，對(duì)參數(shù)間的正負(fù)相關(guān)性和各種尾部相關(guān)性要求至少有一種Copula函數(shù)能描述。這6種Copula函數(shù)及其特征分別為：Gumbel Copula(具有上尾相關(guān)性，適合描述兩個(gè)變量同時(shí)上漲的情況)、Clayton Copula(具有下尾相關(guān)性，適合描述兩個(gè)變量同時(shí)下降的情況)、Frank Copula(可以同時(shí)擬合上、下尾相關(guān)且可以同時(shí)描述變量正相關(guān)性與負(fù)相關(guān)性)、Plackett Copula(對(duì)變量間相關(guān)性的正負(fù)無要求)、AMH Copula(對(duì)正負(fù)相關(guān)性無要求，但是不適合高相關(guān)性情況分析)和FGM Copula(為Plackett Copula關(guān)于(η-1)的泰勒展開式的前兩項(xiàng))。Copula函數(shù)的具體分布函數(shù)及其參數(shù)η,如表5所示。

表5 備選Copula函數(shù)模型[27]

對(duì)以上6種備選Copula函數(shù)同樣進(jìn)行AIC信息準(zhǔn)則檢驗(yàn)，即可選擇出最優(yōu)模型，檢驗(yàn)的具體結(jié)果和最優(yōu)Copula參數(shù)η，如表6所示?？梢钥闯?，PGA、場地參數(shù)ωg和時(shí)-頻調(diào)制函數(shù)參數(shù)a的最優(yōu)Copula模型以Plackett Copula為主；ωg,aft與ξg,aft的最優(yōu)Copula模型同樣以Plackett Copula為主。

表6 最優(yōu)Copula模型

3.3 余震參數(shù)的條件概率密度函數(shù)與條件均值

在確定了主余震參數(shù)的邊緣分布以及主余震參數(shù)之間的最優(yōu)Copula函數(shù)之后，由式(16)～式(19)可得在主震參數(shù)條件下余震的條件概率密度函數(shù)和條件均值。不同主震峰值加速度條件的余震峰值加速度的概率密度函數(shù)和邊緣概率密度的對(duì)比，如圖4所示。由圖4可知，在不同主震參數(shù)影響下，對(duì)應(yīng)余震參數(shù)的概率密度函數(shù)區(qū)別顯著，因此可以證明，考察余震參數(shù)時(shí)有必要考慮主余震參數(shù)之間的相關(guān)性。

圖4 不同主震峰值加速度下的余震峰值加速度的概率密度函數(shù)對(duì)比

條件均值與觀測值的對(duì)比圖，如圖5所示。由圖5可知，Copula條件均值能夠反映出給定主震參數(shù)條件下，余震參數(shù)的總體趨勢。因此，可以利用Copula條件均值預(yù)測給定主震參數(shù)條件下余震參數(shù)的取值。

圖5 Copula條件均值與觀測散點(diǎn)圖對(duì)比

在本文中，對(duì)主余震參數(shù)的條件均值曲線進(jìn)行了最高8次的多項(xiàng)式擬合，得到了主余震參數(shù)的衰減關(guān)系式。表7中給出了Copula條件均值多項(xiàng)式擬合系數(shù)，其中，PGAmain和PGAaft的單位為g；ωg,main和ωg,aft的單位為rad/s；amain和aaft的單位為s-1。由于Ⅰ0類場地與Ⅳ類場地?cái)?shù)據(jù)較少，因此不對(duì)這兩類場地進(jìn)行Copula均值計(jì)算，但是本文對(duì)總體記錄的場地參數(shù)進(jìn)行了分析并進(jìn)行了多項(xiàng)式擬合來彌補(bǔ)這個(gè)問題。由此便建立了主余震參數(shù)之間的衰減關(guān)系，通過給定一套主震參數(shù)，便可以得到一套對(duì)應(yīng)余震參數(shù)，這為工程應(yīng)用提供了方便。

表7 Copula條件均值的多項(xiàng)式擬合系數(shù)

4 主余震型地震動(dòng)過程的降維表達(dá)

4.1 非平穩(wěn)地震動(dòng)過程的降維表達(dá)

對(duì)于一個(gè)演變功率譜密度函數(shù)為SUg(t,ω)的零均值實(shí)非平穩(wěn)地震動(dòng)過程Ug(t)，其源譜表示為

(21a)

(21b)

式中:Δω為頻率步長，Δω=ωu/N，ωu為截?cái)囝l率，N為截?cái)囗?xiàng)數(shù);{Xk,Yk}為一組標(biāo)準(zhǔn)正交隨機(jī)變量，滿足如下基本條件

E[Xk]=E[Yk]=0，E[XjYk]=0

(22a)

E[XjXk]=E[YjYk]=δjk

(22b)

式中：E[·]為數(shù)學(xué)期望；δjk為Kronecker-delta記號(hào)。

在式(21)中，由于隨機(jī)變量Xk和Yk的概率分布未給定，因此不能直接用于模擬。定義正交隨機(jī)變量集{Xk,Yk}(k=1,2,…,N)為如下形式

(23)

式中，φk(k=1,2,…,N)為一組在(0,2π)上均勻分布且相互獨(dú)立的隨機(jī)相位角。顯然，式(23)定義的隨機(jī)變量集{Xk,Yk}(k=1,2,…,N)滿足基本條件式(25)。

將式(23)代入式(21a)中，即可得到非平穩(wěn)地震動(dòng)隨機(jī)過程的傳統(tǒng)隨機(jī)相位模擬公式

(24)

式中,Ug1(t)即為采用傳統(tǒng)隨機(jī)相位角的譜表示模擬過程?？梢钥闯鰝鹘y(tǒng)模擬方法在本質(zhì)上屬于經(jīng)典的蒙特卡洛方法，面臨著高維隨機(jī)變量所帶來的巨大挑戰(zhàn)。為了克服傳統(tǒng)蒙特卡洛方法隨機(jī)變量數(shù)量過大，模擬效率低下的缺陷，引入隨機(jī)函數(shù)的降維思想[28]，將標(biāo)準(zhǔn)正交隨機(jī)變量集{Xk,Yk}(k=1,2,…,N)定義為基本隨機(jī)向量的正交函數(shù)形式，具體如下

(25)

4.2 主余震序列的構(gòu)造

在主余震型地震動(dòng)過程的降維模擬中，主震地震動(dòng)過程和余震地震動(dòng)過程均采用1.1節(jié)中的演變功率譜密度模型，其主要區(qū)別在于模型參數(shù)的取值。通過給定的場地條件、地震烈度，確定主震的峰值加速度和持時(shí)。對(duì)于主震的場地土卓越圓頻率、阻尼比和調(diào)制函數(shù)，則根據(jù)場地類別采用參數(shù)識(shí)別結(jié)果的均值。進(jìn)一步，利用表7的得到對(duì)應(yīng)的余震參數(shù)，再把主震與余震的兩套參數(shù)代入式(2)、式(4)與式(1)中便可得到主震和余震的演變功率譜密度模型。

這樣，就實(shí)現(xiàn)了僅需兩個(gè)基本隨機(jī)變量即可精細(xì)的模擬主余震型地震動(dòng)過程的目的。本文方法能夠有效避免傳統(tǒng)蒙特卡羅方法帶來的樣本數(shù)量龐大與樣本集合概率信息不完備的問題，僅僅需要數(shù)百條樣本即可在全概率層面上反映主余震型地震動(dòng)過程的概率特性，具有與概率密度演化方法[30]結(jié)合的天然優(yōu)勢。

5 算例分析與驗(yàn)證

5.1 算例及分析

在本算例中，式(21)的計(jì)算參數(shù)為：頻率截?cái)囗?xiàng)數(shù)N=1 600；截?cái)囝l率ωu=240 rad/s；頻率步長：Δω=0.15 rad/s；時(shí)間步長：Δt=0.01 s；主震和余震的模擬持時(shí)分別取為20 s和15 s。主震和余震其他的地震動(dòng)參數(shù)如表8所示。需要說明的是，表8中主震的峰值加速度考慮8度設(shè)防地震，主震場地卓越圓頻率、阻尼比和時(shí)頻調(diào)制函數(shù)參數(shù)取為II類場地主震參數(shù)識(shí)別結(jié)果的均值。對(duì)于表8中余震參數(shù)，是利用所得主震參數(shù)通過表7求得的。

表8 主余震地震動(dòng)參數(shù)取值

圖6為用本文提出的方法模擬的Ⅱ類場地的主余震代表性時(shí)程曲線及模擬的144條主余震地震動(dòng)過程代表性樣本的均值及標(biāo)準(zhǔn)差與相應(yīng)目標(biāo)值的比較結(jié)果。從圖6可以看出，余震衰減比主震更快，包含的能量更小，且144條代表性樣本集合的模擬結(jié)果均與目標(biāo)值擬合較好，這充分說明了本文方法的有效性。

圖6 主余震代表性樣本以及樣本集合均值和標(biāo)準(zhǔn)差

5.2 主余震頻譜特征驗(yàn)證

在本節(jié)中，將本文模擬方法得到的144條Ⅱ類場地的主余震代表性樣本與篩選得到的主余震實(shí)測地震動(dòng)記錄對(duì)比。為方便起見，將模擬得到的余震時(shí)程與主余震實(shí)測記錄均調(diào)幅至200 cm/s2。計(jì)算Ⅱ類場地的實(shí)測地震動(dòng)記錄的加速度反應(yīng)譜與傅里葉幅值譜，并同本文的模擬結(jié)果進(jìn)行對(duì)比，如圖7所示。可以看出，不論是主震還是余震，實(shí)測強(qiáng)震記錄均值都在模擬均值加減一倍標(biāo)準(zhǔn)差的范圍內(nèi)，這表明本文的模擬方法可以準(zhǔn)確地反映主余震特征，具有良好的工程適用性。

圖7 模擬的主余震地震動(dòng)加速度的反應(yīng)譜和傅里葉幅值譜與實(shí)測記錄的比較

6 結(jié) 論

本文在非平穩(wěn)地震動(dòng)的演變功率譜密度模型的基礎(chǔ)上，建議了該模型參數(shù)的識(shí)別方法，并利用實(shí)測主余震地震記錄識(shí)別了模型參數(shù)。基于Copula理論，實(shí)現(xiàn)了主余震參數(shù)之間的相關(guān)性分析。結(jié)合譜表示-隨機(jī)函數(shù)方法，建立了主余震型地震動(dòng)的降維模型，生成了主余震型地震動(dòng)的代表性時(shí)程。文中實(shí)現(xiàn)了從主余震實(shí)測記錄出發(fā)、獲取參數(shù)的邊緣分布、到生成具有考慮參數(shù)相關(guān)性的主余震型地震動(dòng)過程的代表性時(shí)程，再到與實(shí)測主余震型地震動(dòng)記錄的頻譜特性進(jìn)行對(duì)比的全過程，得出以下結(jié)論：

(1)主震的峰值加速度比余震有偏大的趨勢；余震與主震相比，場地土卓越圓頻率較大，阻尼比較小，時(shí)頻調(diào)制函數(shù)參數(shù)a較大，這代表著余震能量比主震更小，衰減速度比主震更快。

(2)主余震參數(shù)之間大多保持著一定的相關(guān)性，其中峰值加速度PGA、場地土卓越圓頻率和調(diào)制函數(shù)參數(shù)a之間保持著中等的相關(guān)性，而阻尼比之間表現(xiàn)為低相關(guān)性。不同場地的場地土卓越圓頻率和阻尼比之間均保持著中等負(fù)相關(guān)性。此外，在不同主震參數(shù)取值的條件下，余震對(duì)應(yīng)參數(shù)的概率密度函數(shù)有顯著的區(qū)別。

(3)主余震型地震動(dòng)過程的降維模擬方法僅需數(shù)百條樣本即可在全概率信息上反應(yīng)主余震地震動(dòng)過程的概率特性，這為應(yīng)用概率密度演化理論進(jìn)行主余震型地震動(dòng)作用下工程結(jié)構(gòu)的隨機(jī)動(dòng)力反應(yīng)分析與可靠度評(píng)估奠定基礎(chǔ)。