賈宏宇 ，藍(lán)先林 ，陳 航 ，鄭史雄 ，李 晰 ，張永水

（1.西南交通大學(xué)土木工程學(xué)院，四川 成都 610031；2.北京工業(yè)大學(xué)城市與工程安全減災(zāi)省部共建教育部重點(diǎn)實(shí)驗(yàn)室，北京 100124；3.中鐵二院工程集團(tuán)有限責(zé)任公司，四川 成都 610031；4.重慶交通大學(xué)土木工程學(xué)院，重慶400074）

大跨度結(jié)構(gòu)物的多點(diǎn)地震激勵(lì)分析能較為真實(shí)反映地震作用下各支撐處不同的地面運(yùn)動(dòng)，能更加真實(shí)地考慮實(shí)際地震對(duì)長(zhǎng)大結(jié)構(gòu)物的作用，從而得到地震作用下更為真實(shí)的結(jié)構(gòu)動(dòng)力響應(yīng)[1-4].但是，有限的實(shí)測(cè)地震波根本不可能滿足抗震分析的需要，同時(shí)也不可能在擬建橋梁的各橋墩處剛好就有實(shí)測(cè)的地震波，且這些地震波正好能為擬建橋梁的抗震分析所用，這使得人工合成地震波成為天然地震波必要的補(bǔ)充[5-6].

對(duì)人工地震波合成，許多學(xué)者提出了各種方法，如三角級(jí)數(shù)法[7-9]、ARMA 方法[10]、小波變換法[11]等，其中三角級(jí)數(shù)法應(yīng)用最為廣泛[12].三角級(jí)數(shù)法假設(shè)相位角在[0，2π]內(nèi)均勻分布，從而使得三角級(jí)數(shù)法只能依靠包絡(luò)函數(shù)來(lái)模擬地震動(dòng)的強(qiáng)度非平穩(wěn)性[13]，縱所周知，人為選擇強(qiáng)度包絡(luò)函數(shù)具有任意性，故所生成的地震波表征的強(qiáng)度非平穩(wěn)性也具有任意性.地震動(dòng)的非平穩(wěn)性除了強(qiáng)度非平穩(wěn)性以外，還具有頻率非平穩(wěn)性.但是，頻率非平穩(wěn)性在三角級(jí)數(shù)法中不能得到體現(xiàn)[14].基于上述原因，提出關(guān)于相位差譜的地震波合成方法，而且這種方法能很好地考慮地震動(dòng)的強(qiáng)度非平穩(wěn)和頻率非平穩(wěn)性[15-16].此外，天然地震動(dòng)又具有空間性和多維性.就目前關(guān)于相位差譜合成地震波的研究而言，幾乎沒有相關(guān)研究考慮地震動(dòng)的多維性.為了能更好地模擬地震動(dòng)的真實(shí)特性，本文基于相位差譜合成了多維多點(diǎn)非平穩(wěn)地震波，能彌補(bǔ)天然地震波的不足，也為抗震規(guī)范人工波合成提供理論參考.

本文基于隨機(jī)振動(dòng)理論，模擬多維多點(diǎn)地震動(dòng)場(chǎng)，基于三角級(jí)數(shù)法，推導(dǎo)了基于相位差譜合成非平穩(wěn)地震波的理論，并分析了相位差譜的分布規(guī)律，最后以某剛構(gòu)橋?yàn)橐劳?，運(yùn)用合成的非平穩(wěn)地震波對(duì)剛構(gòu)橋進(jìn)行了地震響應(yīng)分析，并與反應(yīng)譜法進(jìn)行對(duì)比，驗(yàn)證其正確性和精度.

1 多維多點(diǎn)人工地震動(dòng)模擬

1.1 多維多點(diǎn)地震動(dòng)場(chǎng)

為了表示結(jié)構(gòu)m個(gè)地面支撐點(diǎn)處的地震激勵(lì)的相互影響，用加速度表示的互功率譜矩陣為

式中：矩陣主對(duì)角元素Skk（iω）表示任意點(diǎn)k處的自功率譜密函數(shù)；非主對(duì)角元素Skl（iω）表示任意兩點(diǎn)k、l的互功率譜密度函數(shù)；pkl（iω）表示相干函數(shù)；|pkl（iω）|表示相干函數(shù)的模，反映部分相干效應(yīng)；相干函數(shù)的指數(shù)部分- ωdkl/υapp表示相干函數(shù)的幅角，相位體現(xiàn)地震空間性中的行波效應(yīng)；dkl表示支撐點(diǎn)k和l之間的水平距離沿波的傳播方向的投影；υapp表示視波速，在實(shí)際地震中，視波速隨頻率變化，離散性較大.

在橋梁抗震分析中，取多個(gè)確定的視波速分別計(jì)算，其結(jié)果能包絡(luò)視波速離散性的影響，所以在考慮行波效應(yīng)的抗震分析中，取多個(gè)視波速是必要的，也是合理的.

為了研究地震動(dòng)的多維性，將一維功率譜密度函數(shù)矩陣（式（1））中每一個(gè)元素按式（4）擴(kuò)充為三維互功率譜譜密度函數(shù)矩陣.

式中：x、y和z分別表示兩水平向分量和豎直向分量.

各分量之間有一定比例關(guān)系，根據(jù)文獻(xiàn)[17]：水平分量各個(gè)自功率譜密度相同且相關(guān)，則可得到

豎直分量的互功率譜也可根據(jù)式（6）獲得.

將式（4）、（5）和（6）代入式（1），將只考慮一維（m×m）地震激勵(lì)的互功率譜密度函數(shù)矩陣擴(kuò)展到了能同時(shí)考慮3 個(gè)方向地震激勵(lì)相互影響的三維（3m× 3m）的互功率譜矩陣.

1.2 三角級(jí)數(shù)法合成多維多點(diǎn)地震動(dòng)

三維的地震加速度互功率密度函數(shù)矩陣是Herrmite 矩陣，而且是正定矩陣.根據(jù)式（7），對(duì)三維互功率譜矩陣進(jìn)行LDLT 分解，即

式中：P為3m×r的矩陣，r為矩陣S(ω)的秩.

P*為3m× 3m維矩陣，用lm'（iω）表示x、y和z方向中任意方向在P*矩陣中非主對(duì)角的元素，那么地震動(dòng)的幅值am'和相位角θm'計(jì)算方法為

在合成第m′點(diǎn)任意方向dx、dy和dz的平穩(wěn)地震動(dòng)時(shí)，本文考慮了第m′點(diǎn)任意方向地震動(dòng)與另外任意點(diǎn)各方向的相關(guān)性，那么，第m′點(diǎn)dx或dy、dz方向平穩(wěn)地震動(dòng)表示為

式中：am′b和 θm′b分 別為點(diǎn)m′與b在第k個(gè)頻率時(shí)的幅值和相位角；φbk為點(diǎn)m′與b在第k個(gè)頻率成分中的相位變化值；t為時(shí)間.

采用三角級(jí)數(shù)法合成人工地震波時(shí)，φbk為[0，2π]區(qū)間內(nèi)均勻分布的隨機(jī)變量.運(yùn)用包絡(luò)函數(shù)模擬地震動(dòng)非平穩(wěn)性中的強(qiáng)度非平穩(wěn)，但是強(qiáng)度包絡(luò)函數(shù)的選擇具有很大的任意性，也不能考慮地震動(dòng)的頻率非平穩(wěn)性.真實(shí)地震動(dòng)具有強(qiáng)度和頻率雙重非平穩(wěn)性，雙重非平穩(wěn)性與相位差有著密切的關(guān)系，因此，采用相位差譜合成人工地震波更能表征地震動(dòng)的雙重非平穩(wěn)性.

2 地震動(dòng)相位差譜

2.1 相位差譜和脈動(dòng)相位差譜

Ohsaki[15]認(rèn)為實(shí)際地震動(dòng)的相位譜是隨機(jī)變量，且具有一定相關(guān)性，并用相位差譜的形式來(lái)表示其相關(guān)性.相位差可表示為

式中：Δφ(f)為相位差，其定義域?yàn)閇-2π，0]，若其值大于0，則減去2π 轉(zhuǎn)換至定義域內(nèi)；φk和φk+1為相鄰兩個(gè)頻率對(duì)應(yīng)的相位角；f為頻率.

文獻(xiàn)[18]對(duì)真實(shí)強(qiáng)震加速度記錄的相位差譜進(jìn)行研究，回歸得到均值函數(shù)Δ(f)為

式中：M為震級(jí)；R為震源距；R0為常數(shù)，其值取15；a1（f）、a2（f）和a3（f）為與頻率有關(guān)的回歸系數(shù)，其取值參考文獻(xiàn)[18].

ε(f)的標(biāo)準(zhǔn)差σε根據(jù)式（15）計(jì)算獲得.

式中：d1、d2和d3為回歸系數(shù)，當(dāng)處于巖石地基時(shí)，d1=0.089，d2= -1.124，d3= 0.316.

圖1給出了震距R= 100 km，不同震級(jí)時(shí)相位差譜均值變化規(guī)律.

圖1 相位差譜均值趨勢(shì)曲線Fig.1 Mean curve trend chart of phase difference spectrum

2.2 相位差譜的分布律

為了研究相位差譜的分布規(guī)律，對(duì)常用的實(shí)際地震波中的El Centro 波S00E 脈動(dòng)相位差譜進(jìn)行頻數(shù)分析[18]，得出隨機(jī)變量εb滿足對(duì)數(shù)正態(tài)分布規(guī)律，即隨機(jī)變量εb取對(duì)數(shù)有 εc=lnεb，那么εc服從正態(tài)分布：

則 εb滿足式（17）形式的對(duì)數(shù)正態(tài)分布：

隨機(jī)變量εc與εb的均值、均方差有如下關(guān)系：

式（16）～（18）中：參數(shù)mc、σc分別為隨機(jī)變量εc的均值和均方差；mb、σb為εb的均值和均方差.

2.3 基于相位差譜的人工地震動(dòng)合成

均值mb= 2π，σb可根據(jù)式（15）回歸得到，再根據(jù)式（18）獲得εb的均值和均方差.最后根據(jù)分布規(guī)律計(jì)算得到滿足相應(yīng)條件的隨機(jī)數(shù)εb，從而就得到脈動(dòng)相位差ε（f），即

再將ε（f）與式（13）得到的 Δ(f)相加，并結(jié)合式（10），就合成了非平穩(wěn)人工地震動(dòng).

3 實(shí)際算例

為了說(shuō)明本文方法的正確性和準(zhǔn)確性，采用傳統(tǒng)反應(yīng)譜法來(lái)進(jìn)行驗(yàn)證工作.以某高墩剛構(gòu)橋（89 m +160 m + 89 m）為對(duì)象.該橋位處6 級(jí)震區(qū)，震源距200 km.首先，選擇對(duì)應(yīng)的設(shè)計(jì)反應(yīng)譜，特征周期Tg= 0.3 s，

抗震重要系數(shù)Ci= 1.7，場(chǎng)地系數(shù)Cs= 0.9，阻尼調(diào)整系數(shù)Cd= 1.0，加速度峰值A(chǔ)= 0.15g，0.20g，反應(yīng)譜最大值Smax= 2.25CiCsCdA，場(chǎng)地類型為I 類場(chǎng)地，轉(zhuǎn)化為功率譜并模擬地震動(dòng)場(chǎng)（縱向、橫向和豎向地震波傳播速度分別為300、500、400 m/s）；然后基于相位差譜合成多維多點(diǎn)地震波（地震動(dòng)持續(xù)時(shí)間為50 s），再迭代修正合成地震波；最后得到滿足精度的人工地震波.本橋兩橋墩底處（1# 墩和2# 墩）設(shè)計(jì)反應(yīng)譜的最大值分別為0.15g和0.20g，見圖2.

3.1 非平穩(wěn)地震動(dòng)合成及修正

首先將設(shè)計(jì)反應(yīng)譜視為目標(biāo)反應(yīng)譜，且轉(zhuǎn)為功率譜，構(gòu)造具有多維性和地震動(dòng)空間性的互功率譜矩陣，基于式（10）和式（20）合成多維多點(diǎn)非平穩(wěn)地震動(dòng)時(shí)程.為了提高合成精度，對(duì)合成的地震動(dòng)進(jìn)行修正.在頻域內(nèi)，對(duì)合成地震動(dòng)的幅值進(jìn)行迭代修正.修正方法和步驟如下：

步驟1將合成地震動(dòng)時(shí)程轉(zhuǎn)化為反應(yīng)譜，并與目標(biāo)反應(yīng)譜做商，即

圖2 連續(xù)剛構(gòu)橋示意（單位：m）Fig.2 Schematic diagram of continuous rigid frame bridge(unit：m)

式中：SaT（ω，ζ）為目標(biāo)反應(yīng)譜；Sa（ω，ζ）為合成地震動(dòng)轉(zhuǎn)化的反應(yīng)譜；ζ為阻尼比.

步驟2精度不滿足要求，進(jìn)行幅值修正，再次合成地震動(dòng).

根據(jù)式（18）對(duì)幅值進(jìn)行修正.

步驟3轉(zhuǎn)到第1 步進(jìn)行反應(yīng)譜對(duì)比，如滿足精度，則跳出此步驟，不滿足就進(jìn)入第2 步.

經(jīng)上述流程，人工合成橋梁墩底處初始（未迭代修正）多維多點(diǎn)非平穩(wěn)地震動(dòng)見圖3～5.同時(shí)，本文也給出初始合成的人工地震波轉(zhuǎn)化成的反應(yīng)譜與目標(biāo)反應(yīng)譜的對(duì)比，見圖6.

由圖3～6 可知：（1）從波形來(lái)看，基于相位差譜合成的多維多點(diǎn)非平穩(wěn)地震動(dòng)符合實(shí)際地震動(dòng)非平穩(wěn)性的3 個(gè)階段：增幅階段、平穩(wěn)階段和衰減階段；（2）從幅值來(lái)看，1# 墩底合成地震動(dòng)幅值比2#墩地震動(dòng)幅值大，這與兩墩底的目標(biāo)反應(yīng)譜峰值對(duì)應(yīng)，可從宏觀判斷地震波合成是否正確；（3）1# 墩和2# 墩底處初步合成地震動(dòng)轉(zhuǎn)化的反應(yīng)譜和目標(biāo)反應(yīng)譜基本一致，但1# 墩處幅值最大值相差約40%，2# 墩處相差約57%，且擬合反應(yīng)譜尖峰較多，不平滑，未達(dá)到目標(biāo)擬合精度15%.為了提高精度，需進(jìn)一步迭代修正.根據(jù)式（18）對(duì)初次合成的人工非平穩(wěn)地震動(dòng)時(shí)程進(jìn)行幅值修正，并再次對(duì)比擬合反應(yīng)譜和目標(biāo)反應(yīng)譜，其最終修正的人工非平穩(wěn)地震動(dòng)曲線見圖7，最終擬合反應(yīng)譜與目標(biāo)反應(yīng)譜的對(duì)比見圖8.基于MATLAB 采用Duhamel 積分將新生成地震波轉(zhuǎn)換為擬合反應(yīng)譜，相關(guān)參數(shù)：阻尼比為0.05；積分時(shí)間步長(zhǎng)為0.02 和周期為0～10 s.由于篇幅有限，此處僅列兩個(gè)橋墩處x方向人工合成的任意地震波對(duì)應(yīng)反應(yīng)譜和目標(biāo)反應(yīng)譜的對(duì)比.由圖7和圖8可知，經(jīng)過(guò)迭代合成的地震動(dòng)同樣具有前述非平穩(wěn)特性，且在1# 墩和2# 墩處幅值最大值相差分別變?yōu)榱思s6%和13%，且擬合反應(yīng)譜尖峰較少，曲線平滑.同時(shí)，經(jīng)過(guò)兩個(gè)墩處地震動(dòng)合成發(fā)現(xiàn)，目標(biāo)反應(yīng)譜峰值越大，擬合后峰值相差也越大，如2#墩處目標(biāo)反應(yīng)譜峰值0.2g，擬合反應(yīng)譜與目標(biāo)反應(yīng)譜差值比1# 墩處大.修正后1# 墩和2# 墩處的地震動(dòng)幅值縮減最大值分別為0.745 m/s2和1.580 m/s2.

圖3 x 方向非平穩(wěn)地震動(dòng)（未修正）Fig.3 Non-stationary ground motion at x direction（uncorrected）

圖4 y 方向非平穩(wěn)地震動(dòng)時(shí)程（未修正）Fig.4 Non-stationary ground motion at y direction（uncorrected）

圖5 z 方向非平穩(wěn)地震動(dòng)（未修正）Fig.5 Non-stationary ground motion at z direction（uncorrected）

圖6 擬合反應(yīng)譜和設(shè)計(jì)反應(yīng)譜（未修正）Fig.6 Fitting response spectrum and design response spectrum（uncorrected）

圖7 x 方向的非平穩(wěn)人工地震時(shí)程（修正）Fig.7 Non-stationary ground motion at x direction（corrected）

圖8 擬合反應(yīng)譜和目標(biāo)設(shè)計(jì)反應(yīng)譜（修正）Fig.8 Non-stationary fitting response spectrum and target design response spectrum（corrected）

3.2 橋梁地震響應(yīng)分析

為了進(jìn)一步研究修正前與修正后地震波對(duì)橋梁結(jié)構(gòu)動(dòng)力響應(yīng)的影響，也為驗(yàn)證人工合成地震波的正確性.在圖2剛構(gòu)橋?qū)?yīng)的有限元模型中，在1#墩和2# 墩底基礎(chǔ)處分別施加合成的非平穩(wěn)人工地震動(dòng)，對(duì)橋梁結(jié)構(gòu)進(jìn)行多點(diǎn)激勵(lì)分析.圖9是基于本文合成地震動(dòng)的時(shí)程分析法和傳統(tǒng)反應(yīng)譜法的工況說(shuō)明圖例.圖9中列出了5 種工況：工況1，1# 墩和2# 墩峰值加速度都為0.15g時(shí)程分析法（一致激勵(lì)）；工況2，1# 墩和2# 墩峰值加速度分別為0.15g和0.2g時(shí)程分析法（多點(diǎn)激勵(lì)）；工況3，1# 墩和2# 墩峰值加速度都為0.2g時(shí)程分析法（一致激勵(lì)）；工況4，1# 墩和2# 墩峰值加速度都為0.15g反應(yīng)譜法（一致激勵(lì)）；工況5，1# 墩和2# 墩峰值加速度都為0.2g反應(yīng)譜法（一致激勵(lì)）.圖中：UE 為一致激勵(lì)；PDS為相位差譜法；ME 為多點(diǎn)激勵(lì)；RS 為反應(yīng)譜法.

圖10～12 給出了各工況荷載作用下，剛構(gòu)橋墩底剪力和彎矩的對(duì)比情況，由圖10可知，加速度峰值越大，墩底的剪力和彎矩值也越大，初步判斷地震作用下結(jié)構(gòu)響應(yīng)規(guī)律符合實(shí)際情況.采用反應(yīng)譜法比基于相位差譜合成地震動(dòng)的時(shí)程方法（本文方法）計(jì)算得到的剪力和彎矩要大，當(dāng)加速度為0.2g時(shí)，1# 墩的剪力最大相差約3 973.5 kN，2# 墩的彎矩值最大相差889 146.2 kN·m；2# 墩比1# 墩剪力值小，但2# 墩比1# 墩彎矩值大，其中采用反應(yīng)譜法計(jì)算的剪力值相差在10%以內(nèi)，彎矩值相差在17%以內(nèi)，采用本文時(shí)程分析法計(jì)算的剪力值相差在9%，彎矩相差在12%以內(nèi).這說(shuō)明了本文方法和傳統(tǒng)反應(yīng)譜法計(jì)算結(jié)果基本一致，也說(shuō)明了本文方法正確性.

圖9 連續(xù)剛構(gòu)計(jì)算工況圖例Fig.9 Analysis case of continuous rigid frame bridge

圖10 一致激勵(lì)與多點(diǎn)激勵(lì)的墩底剪力彎矩（相位差譜法）Fig.10 Shearing and moment at the bottom of the piers under UE and ME（phase difference spectrum method）

圖11 墩底處剪力Fig.11 Shearing of the piers bottom

圖12 墩底處彎矩Fig.12 Moment of the piers bottom

4 結(jié) 論

本文基于相位差譜合成了多維多點(diǎn)非平穩(wěn)地震動(dòng)，采用頻域方法對(duì)合成地震動(dòng)進(jìn)行迭代修正.并研究了人工非平穩(wěn)地震動(dòng)對(duì)剛構(gòu)橋多點(diǎn)激勵(lì)下的響應(yīng)情況.得到如下結(jié)論：

（1）基于相位差譜合成的多維多點(diǎn)非平穩(wěn)地震動(dòng)不用包絡(luò)函數(shù)控制地震波的波形，較好地排除了人工選擇強(qiáng)度包絡(luò)函數(shù)的任意性，優(yōu)于傳統(tǒng)的三角級(jí)數(shù)法；

（2）基于相位差譜合成多維多點(diǎn)地震動(dòng)不但考慮了地震動(dòng)的空間相關(guān)性，還考慮了地震動(dòng)的多維性，更加符合真實(shí)地震動(dòng)特征；

（3）本文方法與傳統(tǒng)反應(yīng)譜法計(jì)算結(jié)果基本一致，采用反應(yīng)譜法計(jì)算的1# 和2# 墩剪力值相差在10%以內(nèi)，彎矩值相差在17%以內(nèi)，而本文時(shí)程分析法分別在9%和12%以內(nèi).