Housner譜烈度的一種快速近似算法1

張立寶 王玉石 李小軍 蘭日清

1）中國地震局地球物理研究所，北京 100081

2）北京工業(yè)大學(xué)，建筑工程學(xué)院，北京 100124

引言

地震儀器烈度速報的關(guān)鍵科學(xué)問題之一，是尋找與震害相關(guān)性高的地震動強度表征物理參數(shù)指標(biāo)。由于地震引起的強地面運動是典型的非平穩(wěn)信號，幅值、頻率成分均隨時間變化顯著，地震動強度及其對建筑結(jié)構(gòu)、地質(zhì)災(zāi)害等的影響難以用地面加速度峰值PGA（Wald等，1999；李山有等，2002；Linkimer，2008）、地面速度峰值PGV（Wald等，1999；Wu等，2003；Atkinson等，2007）等單一參數(shù)指標(biāo)進行可靠描述。除地面運動峰值參數(shù)外，加速度反應(yīng)譜（Boatwright等，2001；Atkinson等，2007；李敏，2010）、傅立葉幅值譜（Sokolov等，1998；Chernov等，1999；Sokolov，2002）、累積絕對速度CAV（Fahjan等，2011）、阿里亞斯烈度Arias Intensity（Margottini等，1992）、日本氣象廳烈度IJMA（Karim等，2001，2002；Shabestari等，2001；金星等，2013）、持續(xù)時間（Trifunac等，1975，1977）等地震動參數(shù)也被用來描述地震動強度。利用在中國西部地震（如汶川地震）中獲得的強震動記錄進行檢驗，結(jié)果表明，地面運動峰值參數(shù)與地震烈度的相關(guān)性較低，與反應(yīng)譜相關(guān)的參數(shù)指標(biāo)可靠性相對較高（王玉石等，2008，2013；林淋等，2011；馬強等，2014；朱永莉等，2015）。

研究發(fā)現(xiàn)，著名地震工程學(xué)家Housner（1952）提出的地震動參數(shù)指標(biāo)“譜烈度”與地震破壞程度及地震烈度具有較高的相關(guān)性，在評估地震破壞時顯著優(yōu)于其它地震動參數(shù)（Nau等，1984；王玉石等，2008，2013；Bradley，2009；Ueong，2009），相關(guān)預(yù)測模型被廣泛應(yīng)用于地震動強度（地震烈度）快速評估中（Martínez-Rueda，1998；Bradley等，2009a，2009b，2010a，2010b；Ueong，2009；Bradley，2011）。日本科學(xué)家研制了譜烈度傳感器（Takubo等，1999；Koganemaru等，2000），并將譜烈度作為可靠的地震動強度指標(biāo)取代了地面加速度峰值，廣泛應(yīng)用于鐵路、燃氣管網(wǎng)與東京灣區(qū)的地震動強度監(jiān)測與地震緊急處置系統(tǒng)中（Katayama等，1986；Iwata，1991；Takubo等，1999）。中國臺灣大臺北區(qū)瓦斯股份有限公司也在臺北盆地內(nèi)安裝了譜烈度儀用于地震緊急處置（Shimizu等，2000；Ueong，2009）。王玉石等（2008）通過對美國加州地震數(shù)據(jù)的統(tǒng)計分析，綜合利用地面加速度峰值、地面速度峰值和譜烈度3個地震動參數(shù)，建立了1種地震動強度（烈度）的快速判別方法，顯著提高了地震儀器烈度的可靠性，并被“中國數(shù)字強震動臺網(wǎng)”地震動強度（烈度）速報軟件采用。李世成等（2011）通過對2009年姚安6.0級地震強震動記錄的統(tǒng)計分析，發(fā)現(xiàn)譜烈度是建立強震動峰值參數(shù)與地震烈度之間對應(yīng)關(guān)系的有效標(biāo)度。利用中國西部地區(qū)多次強震獲得的豐富強震動記錄，王玉石等提出可根據(jù)加速度反應(yīng)譜的卓越周期對譜烈度進行場地效應(yīng)修正，并建立了利用修正譜烈度確定中國西部地區(qū)儀器烈度的方法，進一步提高了地震儀器烈度的可靠性，并成功重現(xiàn)了汶川地震的地震動強度分布圖（王玉石等，2010；Wang等，2013）。研究結(jié)果表明，修正譜烈度與地震破壞程度具有很高的相關(guān)性，是目前已知的最可靠的地震儀器烈度物理參數(shù)指標(biāo)（Wang等，2012；王玉石等，2013；李世成等，2013）。

雖然譜烈度作為可靠的地震動強度表征參數(shù)獲得了廣泛應(yīng)用，但是由于地面運動在水平面內(nèi)不同方向的分量存在差異，造成譜烈度計算時間相對較長，限制了譜烈度這一可靠的儀器烈度物理參數(shù)指標(biāo)在地震烈度速報等方面的應(yīng)用，亟需對其計算方法予以改進。

1 譜烈度定義

Housner（1952）給出譜烈度的定義如下：

式中，SV是阻尼比為ξ的相對速度反應(yīng)譜，T為周期，阻尼比ξ一般取20%。譜烈度是地震動強度的1個定量指標(biāo)，其物理意義是對應(yīng)地震動作用下自振周期為0.1—2.5s（即自振頻率為0.4—10Hz）范圍內(nèi)單自由度阻尼體系地震響應(yīng)相對速度峰值的平均值，它由強震動加速度記錄時程完全確定，不涉及任何宏觀震害現(xiàn)象。

單一加速度分量對應(yīng)的譜烈度計算方法見圖1（a），圖中曲線為阻尼比為20%的相對速度反應(yīng)譜。但是，強震動加速度記錄一般包括2個互相正交的水平向分量和1個豎向分量，分別計算每1個分量對應(yīng)的譜烈度，則會出現(xiàn)同1個地點有3個不同譜烈度值的情況。為了避免這一情況，一般采用Karim等（2002）建議的方法，即考慮到建筑結(jié)構(gòu)破壞主要是受水平向地震作用的影響，計算地震動強度參數(shù)時僅采用水平向分量，譜烈度值也取為水平面內(nèi)不同方向上地面運動分量對應(yīng)譜烈度的最大值（圖1（b））。為了獲得更加準(zhǔn)確的觀測點譜烈度值，需要在水平面360°內(nèi)不同方向上對加速度時程進行矢量正交分解，獲得每個方向上的加速度分量，然后計算每個方向上的加速度分量對應(yīng)的譜烈度，取其最大值為觀測點強震動記錄對應(yīng)的譜烈度值。取角度計算步距Δ＝1°時（Karim等，2002）需要進行180次反應(yīng)譜計算，在普通微型計算機上需要3—7min（依賴于記錄數(shù)據(jù)長度），在嵌入式微處理器上運算所需時間更長。當(dāng)需要計算百余條強震加速度記錄對應(yīng)的觀測點譜烈度值時，計算時間將達小時量級，嚴重影響了譜烈度作為儀器烈度物理參數(shù)指標(biāo)的時效性。

圖1 譜烈度計算方法（以汶川地震綿竹清平臺記錄為例）Fig.1 Calculation method for Housner spectum intensity at station 051MZQ in the Wenchuan earthquake

2 譜烈度快速算法

分析圖1（b）可以發(fā)現(xiàn)，水平面內(nèi)不同方向上的譜烈度曲線為中心對稱圖形，且近似為2個圓形相交后去除重合的部分。為了驗證這一現(xiàn)象是否為普遍規(guī)律，利用在汶川地震主震及余震（圖2（a））、蘆山地震（賴敏等，2014）、魯?shù)榈卣穑ù藿ㄎ牡龋?014）等共計386次MS＞3.0地震中獲取的2701組強震動加速度記錄，取角度計算步距Δ＝1°，分別計算了加速度記錄在水平面內(nèi)的譜烈度跡線，并利用譜烈度最大值及其對應(yīng)方向?qū)ψV烈度跡線進行歸一化與旋轉(zhuǎn)處理（圖2（b））。由圖2（b）可以看出，在水平面內(nèi)譜烈度曲線近似為2個半徑相等的圓相交，這確實是強震動加速度記錄揭示的強地面運動一般規(guī)律。對于不同強震動記錄兩圓的圓心距離各異，對于一小部分強震動記錄兩圓圓心十分接近而退化為一個圓形。根據(jù)這一發(fā)現(xiàn)，提出了1種譜烈度值快速算法，步驟如下：

（1）分別計算EW、N45°E、NS和N45°W方向上加速度分量的譜烈度值（圖3中A、B、C、D點），由中心對稱可得到WE、S45°W、SN和S45°E方向上加速度分量的譜烈度值（圖3中E、F、G、H點）。

（2）選取EW、N45°E、NS和N45°W方向上加速度分量譜烈度的最大值（圖3中D點），并分別選取最大值對應(yīng)方向逆時針旋轉(zhuǎn)45°的加速度分量譜烈度值（圖3中E點）和順時針旋轉(zhuǎn)45°的加速度分量譜烈度值（圖3中C點），三點作圓。

（3）取圓心的矢徑長度與圓的半徑之和為觀測點的譜烈度值快速計算結(jié)果FSI。這一算法的本質(zhì)是找到4個連續(xù)相隔45°方向上的加速度分量對應(yīng)譜烈度的最大值，并利用它與其±45°方向上的譜烈度點進行三點作圓，取圓心的矢徑長度與半徑之和為觀測點的譜烈度近似值。采用快速算法后，計算大幅簡化，獲得1個觀測點譜烈度近似值僅需4次反應(yīng)譜計算，使計算速度提高了45倍，計算時間僅為2—6s。

圖2 譜烈度跡線統(tǒng)計特征Fig.2 Statistical characteristics for the trace of Housner spectum intensity in the horizontal plane

圖3 加速度記錄譜烈度快速算法Fig.3 Fast calculating algorithm for Housner spectum intensity of a three-component acceleration record

3 可靠性驗證

利用在上述386次MS＞3.0地震中獲取的2701組強震動加速度記錄，分別計算得到了對應(yīng)的譜烈度SI與快速譜烈度FSI，并通過兩者的相關(guān)性與相對誤差的統(tǒng)計分析，對所提出的譜烈度快速算法進行了可靠性檢驗。所用強震動記錄的震級在MS3.0—8.0、PGA在0.1—1000cm/s2范圍內(nèi)均有分布，最大PGA為957.7cm/s2，可以認為具有代表性（表1）。

由圖4（a）可以看出，快速譜烈度FSI與譜烈度SI相關(guān)性很強，相關(guān)系數(shù)為0.9999；由圖4（b）可以看出，快速譜烈度FSI的相對誤差均在±5%以內(nèi)，相對誤差在±1%之間的比例高達62.31%。對于地震儀器烈度確定、地震動強度度量等應(yīng)用，此誤差在可以接受的范圍之內(nèi)。

表1 算法可靠性驗證所用強震動記錄的數(shù)量統(tǒng)計Table 1 Acceleration records used in reliability veritification of fast calculating algorithm for Housner spectum intensity

圖4 加速度記錄譜烈度快速算法可靠性Fig.4 Reliability of fast calculating algorithm for Housner spectum intensity

4 結(jié)論與展望

譜烈度作為可靠的地震動強度參數(shù)指標(biāo)，在地震儀器烈度確定與地震緊急處置中得到廣泛應(yīng)用。但由于地震引起的強震動在水平面內(nèi)的方向性差異，傳統(tǒng)算法需要在水平面內(nèi)計算不同方向上加速度分量對應(yīng)的譜烈度值并取其最大值，單組記錄在普通微型計算機上計算需要3—7min，在嵌入式微處理器上計算所需時間更長；對于破壞性大地震（如汶川地震），傳統(tǒng)算法的譜烈度單機計算時間將達小時量級，嚴重影響了譜烈度作為儀器烈度物理參數(shù)指標(biāo)的時效性。

通過對386次MS＞3.0地震中獲取的2701組強震動加速度記錄進行統(tǒng)計分析，發(fā)現(xiàn)水平面內(nèi)不同方向上的譜烈度跡線為中心對稱圖形，且近似為2個半徑相同的圓形相交后去除重合部分這一普遍規(guī)律，并根據(jù)這一發(fā)現(xiàn)，利用三點作圓提出了1種譜烈度快速算法。此算法在計算三分量加速度記錄對應(yīng)的水平面內(nèi)譜烈度最大值時，計算速度提高了45倍，且計算誤差在±4.5%以內(nèi)，在保證譜烈度作為地震儀器烈度參數(shù)指標(biāo)可靠性的同時提高了其時效性。期望依賴于其高可靠性與時效性，譜烈度及其快速算法可以在國家地震烈度速報與預(yù)警工程的地震儀器烈度速報中得到應(yīng)用。

致謝：感謝國家強震動臺網(wǎng)中心提供的強震動觀測數(shù)據(jù)下載服務(wù)以及匿名審稿專家的中肯意見與建議。