非平穩(wěn)隨機地震激勵在ANSYS中快速模擬

侯中學(xué)

(四川省交通勘察設(shè)計研究院有限公司，四川 成都 610017)

1 概述

地震動中的非平穩(wěn)性需要進(jìn)行著重考慮[1-2]。在幾乎所有的時程和隨機抗震分析中，地震動的非平穩(wěn)性用一個均勻調(diào)制隨機過程來模擬。但是，均勻調(diào)制方法得到的隨機過程只能表現(xiàn)出激勵功率譜與時間變化的對應(yīng)關(guān)系，而不能反映地震前期大量的高頻成分以及高頻成分迅速地衰減直至只剩低頻成分的過程[3-5]，因此，在實現(xiàn)對地震動非平穩(wěn)性模擬的方法中，用非均勻調(diào)制的隨機過程來模擬真實地震動在時域和頻率內(nèi)的變化是合適的。

當(dāng)前本領(lǐng)域研究人員的關(guān)注點主要放在尋找提高隨機振動計算效率的方法上。在復(fù)雜結(jié)構(gòu)的響應(yīng)分析時考慮多維多點平穩(wěn)隨機地震激勵，如果要考慮均勻調(diào)制的非平穩(wěn)隨機地震動激勵，則在編程建模、計算方法、計算硬件等方面都有很高要求，在這些前提下結(jié)構(gòu)的響應(yīng)計算已經(jīng)十分困難，如果再考慮非均勻調(diào)制模擬的非平穩(wěn)隨機地震動激勵，那么復(fù)雜結(jié)構(gòu)的響應(yīng)分析就會更加費力[6-10]。同時，由于隨機振動理論非常復(fù)雜，一般情況下，需要計算者自編程序來實現(xiàn)，這樣很難被工程師們所接受。不少學(xué)者對于隨機振動的計算效率和是否需要編程實現(xiàn)等問題進(jìn)行研究，也取得了不少的成果[11-18]。大跨復(fù)雜結(jié)構(gòu)遭受由非均勻調(diào)制方法得到的多維多點非平穩(wěn)隨機地震激勵下的動力響應(yīng)結(jié)果需要進(jìn)一步研究，下面將使用隨機振動法解決以上問題，通過精細(xì)積分法，可以提高計算效率的同時不需編制專門程序(完全通過有限元軟件，如ANSYS完成)，此方法對地震動特性有充分的考慮。

為了將非平穩(wěn)隨機振動理論應(yīng)用到實際工程抗震分析中去，本文在上述研究的基礎(chǔ)上，首先，多維多點非平穩(wěn)虛擬地震動激勵法與精細(xì)積分法進(jìn)行結(jié)合的計算方法，從計算效果上大大提高非均勻調(diào)制的非平穩(wěn)地震進(jìn)行激勵時結(jié)構(gòu)響應(yīng)的計算速度和精準(zhǔn)程度；給出單自由度體系由非平穩(wěn)激勵計算得到的理論解，同時檢驗該方法的精準(zhǔn)程度；其次，將本文方法在ANSYS中進(jìn)行快速模擬，節(jié)約了自編程序的時間，便于向?qū)嶋H工程應(yīng)用推廣；最后，以實際高墩橋梁為例，對該高墩橋梁在多維多點非均勻調(diào)制模擬方法得到的非平穩(wěn)地震動激勵下的結(jié)構(gòu)響應(yīng)進(jìn)行分析，完成計算效果的驗證。

2 非均勻調(diào)制非平穩(wěn)地震的虛擬激勵力的構(gòu)造

2.1 非平穩(wěn)激勵地震動場模擬

大跨復(fù)雜結(jié)構(gòu)m個支撐處受到3個方向激勵的多維多點非平穩(wěn)地震動場功率譜密度函數(shù)表示為：

(1)

其中任意子矩陣Skl(iω,t)是兩個水平分量(x,y)和一個豎向分量(z)的3×3矩陣,具體表示如下：

(2)

兩水平方向的功率譜密度函數(shù)被認(rèn)為是相同的,而豎向與水平向的相干系數(shù)取0.6,因此表示如下：

Sklxx(iω,t)=Sklyy(iω,t)=Sklxy(iω,t)=Sklyx(iω,t)

(3)

(4)

互功率譜密度矩陣中的元素可表示為:

(5)

其中,Skkxx(ω,t),Skkxx(iω)均為第k點上x分量的自功率譜密度函數(shù)；ρklxx(iω)為k,l兩點上x分量的相干函數(shù)。

2.2 直接求解虛擬激勵法虛擬力構(gòu)造

將式(1)中互功率譜矩陣進(jìn)行分解得到：

S0(iω,t)=P*PT=[A(ω,t)][V]*[SCP]{q0}{q0}T[SCP][V]T[A(ω,t)]T

(6)

其中,P為3m×r的虛擬力矩陣,r為矩陣S0(iω,t)的秩,上標(biāo)“*”和“T”分別為復(fù)數(shù)的共軛和轉(zhuǎn)置。

SCP為3 m×3 m維矩陣,表示如下：

(7)

其中，SCPmx(ω)為第m個點x分量的功率譜密度函數(shù)。

V為3 m×3 m維行波因子矩陣,表示如下:

[V]=diag[e-iωT1x,e-iωT1y,e-iωT1z,e-iωT2x,e-iωT2y,
e-iωT2z,…,e-iωTmx,e-iωTmy,e-iωTmz]

(8)

其中，Tmx為空間地面運動x分量地震波傳播第m點需要的時間;V為行波效應(yīng)。

根據(jù)空間變化地面運動相干函數(shù)的分解方法,可將完全相干、部分相干和完全不相干虛擬力表示如下：

(9)

2.3 基于直接求解虛擬激勵法的結(jié)構(gòu)響應(yīng)計算

(10)

其中:

(11)

將每一個特征值下得到的響應(yīng)進(jìn)行疊加就得到了部分相干情況下的響應(yīng)功率譜密度如下：

(12)

完全相干:

(13)

式中:

(14)

完全不相干:

(15)

式中:

(16)

根據(jù)式(11)，式(14)，式(16),可以發(fā)現(xiàn)結(jié)構(gòu)的每一個自振周期都對應(yīng)著一個響應(yīng),那么就要求每次達(dá)到頻率點時,都進(jìn)行瞬態(tài)分析計算,雖然虛擬激勵法可以直接求解,相對于使用較多的隨機振動理論已經(jīng)簡化了計算過程,但是如果處理時離散的頻率點數(shù)量太多,那么意味著結(jié)構(gòu)需要計算的瞬態(tài)分析也就越多,計算的過程也就會變得復(fù)雜得多,采用將多維多點的非平穩(wěn)虛擬激勵法和精細(xì)化積分方法結(jié)合起來的結(jié)構(gòu)僅需要計算兩個瞬態(tài)分析就能確定結(jié)構(gòu)任意節(jié)點的響應(yīng)結(jié)果,可以有效簡化計算,減少計算時間。

3 精細(xì)積分理論

將每一個確定頻點下的動力方程式寫成狀態(tài)方程形式[10-11]，具體如下：

(17)

(18)

狀態(tài)方程(17)的一般解為:

(19)

其中，eHt為指數(shù)矩陣,對式(19)進(jìn)行離散積分,設(shè)時間步長為Δt,由遞推法可知,假如在ti時刻的響應(yīng)為V(ti)=Vi,那么ti-1時刻的響應(yīng)V(ti-1)=Vi-1表示如下：

(20)

其中，T=eHΔt為指數(shù)矩陣,據(jù)文獻(xiàn)[10-11]精細(xì)算法為:

T=eHΔt=[eHΔt/m]m

(21)

令τ=Δt/m,當(dāng)取m=2N很大時,Δτ非常小,由泰勒級數(shù)展開式有:

(22)

當(dāng)展開的項數(shù)N=20時,將泰勒級數(shù)節(jié)段后的截斷誤差與計算機進(jìn)行計算的舍入誤差相比,泰勒級數(shù)的截斷誤差要小很多,由以上對比可知,選取展開級數(shù)N=20不會因為截斷造成過大數(shù)值誤差,所以以上T(τ)的計算精度實際上已經(jīng)比計算機的精確解要高出許多,用泰勒計算展開計算是合理的。

[I+Tai]≡[I+Ta,i-1]2=[I+2Ta,i-1+Ta,i-1×Ta,i-1]
(i=1,2,…,N)

(23)

因此,依次類推：

[I+TaN]≡[I+Ta,N-1]2=…=[I+Ta0]m=T(τ)

(24)

每次運算單元矩陣式時[I]不參與計算,為了避免因為Tai的值過小,以至于在程序中被默認(rèn)忽略的情況,使用下列遞推方法進(jìn)行求解：

[Tai]=2[Ta,i-1]+
[Ta,i-1][Ta,i-1] (i=1,2,…,N)

(25)

假設(shè)在每個時間積分的區(qū)間(ti-1,ti)中,荷載值是線性變換的,那么非平穩(wěn)虛擬力從時間的角度上將F(ω,t)進(jìn)行等效離散,即分別在t0,t1,t2,…,tk處的隨機變量F0,F1,F2…并用函數(shù)表示出來，F(xiàn)k為非平穩(wěn)虛擬力按照以下公式表示為:

F(ω,t)=r0+r1(τ-ti-1)

(26)

其中,r0,r1均為不變向量。

將式(26)代入式(20)可得：

Vi=T·Vi-1-H-1[H-1(Fi-Fi-1)/Δt+Fi]+
TH-1[H-1(Fi-Fi-1)/Δt+Fi-1](i=1,2,…,k)

(27)

其中H-1的計算可根據(jù)式(18)可得:

(28)

根據(jù)在每一個固定頻率點上,V0=0,則可得到固定頻率點各時刻的響應(yīng)表達(dá)式為:

(29)

其中：

S3=TS2+S1

(30)

若用Ai,0,Ai,1,…,Ai,i分別表示F0,F1,…，F(xiàn)i前面的系數(shù),Ai,i下標(biāo)的第一個i為ti時刻,第二個i為第i個集中荷載點的位置,則式(30)可簡化表示為：

Vi=Ai,0F0+Ai,1F1+…Ai,i-1Fi-1+Ai,iFi(i=1,2,…,k)

(31)

其中,Ai,0,Ai,1,…,Ai,i是一種和結(jié)構(gòu)自身固有屬性相關(guān)和外界激勵無關(guān)的變量,但這種固有屬性是在結(jié)構(gòu)本身受到外界荷載激勵下所體現(xiàn)出來的。Ai,i為通過遞推方法進(jìn)行計算的,用Vi-1的系數(shù)Ai-1,0,Ai-1,1,…,Ai-1,i-1來依次推導(dǎo)下一步的系數(shù)為：

(32)

在系數(shù)矩陣中所有其他元素都是該系數(shù)矩陣的前兩列元素的重復(fù),同時具有明顯的規(guī)律性。所以只需要矩陣前兩列的元素,便可以確定整個系數(shù)矩陣。令F0=1,Fk=0(k>0),進(jìn)行一個瞬態(tài)分析便可以得到系數(shù)矩陣的第一列所有元素。同樣的,令F0=0,F1=1,Fk=1(k>2),在完成一次瞬態(tài)分析以后得到系數(shù)矩陣?yán)锏牡诙性亍Ｋ灾恍枰獌蓚€條件——兩個頻率點的瞬態(tài)分析就可以得到完整的系數(shù)矩陣,在得到完整的系數(shù)矩陣后,結(jié)構(gòu)中任意點在不同頻率下的非平穩(wěn)時程荷載作用后的響應(yīng)便可以求出。在運用過程中,可以在三個方向按照時間激勵先后順序系分別對結(jié)構(gòu)進(jìn)行激勵,可以求出結(jié)構(gòu)關(guān)鍵點響應(yīng)的系數(shù)矩陣,比如在運用絕對位移直接求解時,分別求出結(jié)構(gòu)關(guān)鍵點響應(yīng)的系數(shù)矩陣A,將施加在結(jié)構(gòu)上的各個方向的激勵時程與相對應(yīng)的系數(shù)矩陣A相乘后進(jìn)行疊加,求出以上固定頻率點的響應(yīng)結(jié)果。那么所有頻率點都按上述方法代入后,即可得到全部頻率下的結(jié)構(gòu)響應(yīng)。這種方法在一定程度上可以消除對每個頻率點都做多維多點的時程分析的必要性,只需通過兩個頻率點可計算出整個結(jié)構(gòu)的響應(yīng)結(jié)果,縮短了運算過程的步驟和計算時間。

4 工程算例

4.1 計算模型介紹及參數(shù)說明

以某高墩橋梁為實際工程算例,圖1為該算例中建立的三維有限元模型,采用梁墩固結(jié)的形式,墩底受到全部自由度上的約束,不過墩底全部平動方向的自由度在加載時會釋放,1號、2號、4號和5號墩的墩頂與梁底是完全耦合的,僅有3號墩的墩頂和梁底在縱向自由度不耦合,橋臺處約束了豎向位移和橫豎扭轉(zhuǎn)。

本文選用Clough-Penzien功率譜模型,參數(shù)譜強度因子S*=0.001 77 m2/s3,場地參數(shù)ωg=15.0 rad/s,ξg=0.6,ωf=1.5 rad/s,ξf=0.6,g(t)函數(shù)選用文獻(xiàn)[14]的形式及其參數(shù)。

4.2 結(jié)果分析

由于篇幅有限,且本文的目的在于介紹非均勻調(diào)制非平穩(wěn)隨機激勵的快速模擬計算理論,故這里僅給出2號橋墩受地震作用時的墩頂橫向剪力響應(yīng)時變功率譜和剪力響應(yīng)均方差如圖2，圖3所示。

從圖3可看出,響應(yīng)功率譜同時隨時間和頻率發(fā)生變化,說明地震動所具有的非平穩(wěn)特性,圖4從時變均方差的角度也可進(jìn)行佐證,因此說明該方法能夠有效模擬非均勻調(diào)制非平穩(wěn)隨機地震激勵。

5 結(jié)論

通過理論研究和計算驗證得出以下結(jié)論:

1)本文方法并未專門編制運算程序,節(jié)約了編程時間,可高效地分析多維多點非均勻調(diào)制非平穩(wěn)隨機地震激勵下結(jié)構(gòu)的響應(yīng)結(jié)果。

2)提出了絕對位移直接求解與精細(xì)積分相結(jié)合的ANSYS快速模擬方法,將N(離散頻率點個數(shù))個離散頻率點的時程分析轉(zhuǎn)換為2×3×m(m為支撐節(jié)點個數(shù))個時程分析,成百上千倍地提高了計算效率。

3)根據(jù)本文所述方法,既可以考慮完全非平穩(wěn)性的地震動,也能直接考慮地震動的場地效應(yīng)、相干效應(yīng)和行波效應(yīng)的影響,與傳統(tǒng)意義上的虛擬激勵法具有相同的效果。