一種考慮震源與傳播路徑的重力壩地震 響應(yīng)分析的新方法

賀春暉， 陳 林

(中國電建集團成都勘測設(shè)計研究院有限公司，四川 成都 610072)

0 前 言

影響大壩動力響應(yīng)分析準確性的因素主要包含三個方面的內(nèi)容：地震波輸入、壩-地基-庫水系統(tǒng)、壩體材料的強度評價準則。壩-地基-庫水系統(tǒng)在理論和數(shù)值模型方面都已經(jīng)相對成熟[1-3]，而集中于材料動力特性以及強度評價標準的研究也取得了一定的成果[4]，但是同樣對于大壩動力反應(yīng)影響很大的地震波輸入問題卻并沒有得到很好的解決。依據(jù)水工建筑物抗震規(guī)范，目前水工建筑抗震設(shè)計中確定基巖加速度峰值的方法一般采用《中國地震烈度區(qū)劃圖》確定基本烈度，對于大型工程設(shè)防依據(jù)專門的地震危險性分析提供的基巖峰值加速度成果評定。上述兩種方法都是基于已有的地震實測記錄，利用概率和統(tǒng)計的宏觀分析方法得出需要的設(shè)計地震動參數(shù)。現(xiàn)行的水工建筑抗震設(shè)計中，獲取地震動時程的方法有如下幾點不足：

(1)工程場址所在的局部小范圍區(qū)域內(nèi)，大多缺乏足夠的、可以用作統(tǒng)計樣本的實測記錄，導(dǎo)致基于這些樣本進行的外推本身準確性不高；

(2)外推時依據(jù)統(tǒng)計規(guī)律得到的衰減關(guān)系，不能反映場地特殊地形地質(zhì)條件的影響；

(3)確定峰值加速度后，根據(jù)規(guī)范給出的設(shè)計反應(yīng)譜求得時間歷程，無法反映由不同地震特異性的影響。

另外，不同的輸入機制對結(jié)果也會造成影響。對于高拱壩或重力壩這樣的空間大跨度結(jié)構(gòu)，采用均勻輸入顯然是不合適的。對于極少數(shù)有實測地震記錄的壩址[2]，可以考慮非均勻輸入，但是對于其他大多數(shù)情況，沒有足夠的依據(jù)來考慮由于局部場地條件導(dǎo)致的幅差和相差。目前已有部分學(xué)者開始將震源和水工結(jié)構(gòu)作為一個整體進行分析：張翠然等[5]通過隨機方法模擬了大崗山水電站場址的最大可信地震；Bayraktar[6]通過選取對應(yīng)于不同震源機制的實測地震記錄研究了斷層破裂與大壩動力響應(yīng)之間的關(guān)系。Keaton[7]通過選取合適的震源參數(shù)得到加速度時間歷程，以此為基礎(chǔ)計算胡佛大壩公路橋的地震非線性響應(yīng)。Zarafani[8]采用高低頻結(jié)合的方法計算近斷層大壩的反應(yīng)譜。

本文綜合考慮地震震源、傳播路徑以及場地反應(yīng)得到適合于大壩地震分析計算的空間非均勻地震動，即：作為大壩動力分析的輸入地震波能夠包含特定震源機制，傳播路徑和場地地形、地質(zhì)條件等信息，以實現(xiàn)合理的非均勻輸入，使大壩地震反應(yīng)分析的輸入環(huán)節(jié)更具可靠性和說服力。本文提出的方法通過局部場地在地震學(xué)和工程學(xué)之間建立聯(lián)系，大壩結(jié)構(gòu)動力分析基于現(xiàn)有的、較為成熟的模型進行，而局部場地地震波場模擬的上游工作則基于地震學(xué)方法展開。在地震學(xué)中，目前模擬工程適用的寬頻帶(0～10 Hz)地震波采用較多的是混合方法(Hybrid Method)[9-11]，即分別用隨機方法(Stochastic Method)[12]模擬高頻部分，確定性方法(Deterministic Method)模擬低頻部分，選取適當(dāng)?shù)慕徊骖l率，濾波后在頻域或者時域內(nèi)疊加得到寬頻帶地震波波形。隨機方法利用Boore的ω-2加速度幅值譜經(jīng)過傅里葉逆變換得到。確定性方法多采用有限差分法或有限元法等數(shù)值方法。有限差分法適應(yīng)復(fù)雜地形的能力弱，需要較多節(jié)點才能達到一定的精度，且有限元法則可能出現(xiàn)頻散或者偽波[13]。譜元法是一種高階有限元法，最早是Patera于1984年計算流體力學(xué)的時候提出，該方法融合了傳統(tǒng)有限元方法的靈活性和偽譜法的精度[14-15]。1998年Dimitri Komatitsch首次采用基于Legendre多項式和GLL積分的譜元法求解彈性波動方程，表明該方法具有收斂快、精度高等優(yōu)勢。本文利用譜元法進行探討性分析，選取模型范圍較小，通過假定震源頻率和細化網(wǎng)格，可以使波的頻率范圍包含重力壩結(jié)構(gòu)的主要振型。

1 震源與傳播路徑模擬

1.1 震源模擬

本文基于二維模型討論重力壩的動力響應(yīng)與行波效應(yīng)及震源相對位置的關(guān)系，暫且不討論斷層空間效應(yīng)的影響，故將震源簡化為點源，震源時間函數(shù)采用Ricker子波，時域表達式見式(1)，頻域表達式見式(2)，峰值頻率為4 Hz的Ricker子波的時域波形和幅值譜，如圖1所示。

(1)

(2)

1.2 譜元法模擬傳播路徑

譜元法的基本出發(fā)點和有限元一樣，也是由伽遼金加權(quán)余量法建立有限元的表達格式，不同之處在于它在每個單元上使用譜方法，即：將形函數(shù)取為無限可微函數(shù)(如Fourier級數(shù))，從而可以變換到譜空間求解。數(shù)值方法中，一般選取截斷的正交多項式(Chebyshev或Legendre多項式)。Dimitri Komatitsch首先基于Legendre多項式對彈性波方程求解。Legendre多項式遞推式如式(3)：

(3)

式中L0(x)=1，L1(x)=x。

譜元法和一般高階有限元方法差別就在于單元內(nèi)插值點和數(shù)值積分點的選取，Legendre譜元法將Legendre多項式的極值點作為單元內(nèi)插值點。以一維單元情形為例，取4階Legendre多項式，共有5個GLL點，其廣義坐標向量為(-1.000 0，-0.654 7，0.000 0，0.654 7，1.000 0)，相應(yīng)的5個拉格朗日插值多項式圖像如圖2所示。數(shù)值積分采用Gauss-Lobatto積分，積分點位置和插值點位置重合，各點積分權(quán)值為：

(4)

式中j=1,2,…,5。

對于5個GLL點的一維單元，則計算出各個點的權(quán)值分別為(0.100 0，0.544 4，0.711 1，0.544 4，0.100 0)，對于Lagrange插值多項式有Li(xj)=δij，質(zhì)量矩陣寫為：

(5)

矩陣中的元素數(shù)值積分后，由公式(6)：

(6)

得到嚴格對角化的質(zhì)量矩陣：

(7)

對于二維情形，可由一維插值函數(shù)的張量積表達，位移的展開形式如式(8)：

圖1 Ricker波的時域波形和幅值譜

(8)

φi(ξ)——一維插值形函數(shù)。

三維情形與二維類似。

圖2 4階譜單元的拉格朗日插值函數(shù)

2 重力壩地震分析模型

2.1 實現(xiàn)方法

傳統(tǒng)的重力壩動力計算中，采用的是沿截斷地基均勻輸入的方法，即從遠處傳來的地震波傳到近場邊界時，各個點的運動形態(tài)一致，對于大壩等大跨度結(jié)構(gòu)，截取的人工邊界的范圍也相應(yīng)較大，這種假設(shè)顯然是不合理的。采用圖3所示的方法將自由場加速度時程轉(zhuǎn)換成等效節(jié)點反力時程輸入建基面，并在地基截斷處考慮黏彈性人工邊界，即可模擬此非均勻性。而自由場加速度則是通過模擬震源和傳播路徑獲得，相比傳統(tǒng)方法更具有說服力。

如圖3所示，分三個步驟進行：

(1)用譜單元離散傳播路徑，通過選取適當(dāng)?shù)木W(wǎng)格尺寸和時間步長，確定有效的頻率范圍。設(shè)定震源，激發(fā)彈性波并利用譜元法求解波動方程，得到自由場建基面各點的加速度響應(yīng)。

(2)將近場地基作為子結(jié)構(gòu)，在邊界上施加黏彈性邊界條件以模擬輻射阻尼效應(yīng)，將第(1)步求得的加速度時程輸入建基面各點求得相應(yīng)的節(jié)點反力時程。

(3)最后加上重力壩結(jié)構(gòu)，并在建基面各點輸入等效節(jié)點反力時程，用有限元方法對重力壩進行動力響應(yīng)分析。

2.2 SPECFEM2D程序網(wǎng)格敏感分析

二維譜單元網(wǎng)格型如圖4所示，長度為10 km，深度為3 km。震源在距離左側(cè)邊界2 km處，接收點距離右側(cè)邊界4 km，震源到接收點距離為4 km。單元采用4階Legendre多項式插值，每個單元含有5×5=25個GLL點，分別計算粗細網(wǎng)格的有效頻率范圍。粗網(wǎng)格譜單元尺寸為100 m×100 m，包含3 000個單元，總自由度數(shù)為97 042；細網(wǎng)格譜單元尺寸為50 m×50 m，包含12 000個單元，總自由度數(shù)為386 082。

粗細網(wǎng)格分別計算得到的不同激震頻率下接收點x向位移時程如圖5和圖6所示，分別用粗細網(wǎng)格計算不同頻率的波，ricker波主頻分別為2 Hz、4 Hz、6 Hz、8 Hz、10 Hz?？梢钥闯觯捎跀?shù)值離散，粗網(wǎng)格模型在主頻為4 Hz的時候波形和頻率成分都發(fā)生了改變，而細網(wǎng)格模型在主頻為4 Hz的時候仍然能夠保持原有的波形和頻率成分。

2.3 重力壩數(shù)值算例

2.3.1 計算條件及參數(shù)

取典型重力壩剖面，體型參數(shù)為高200 m，底寬140 m， 頂寬20 m， 坡比1 ∶0.7。 地基分別向上游、下游和豎直向下的方向延伸1.5倍壩高(300 m)，如圖7所示。壩體混凝土材料密度為2 400 kg/m3，彈性模量2.2E+10Pa，泊松比0.167，地基材料密度2 650 kg/m3，彈性模量2.5E+10Pa。壩體網(wǎng)格尺寸約10 m。壩體單元采用4節(jié)點平面應(yīng)力單元，地基單元采用4節(jié)點平面應(yīng)變單元。

圖3 計算過程示意

圖4 網(wǎng)格模型示意

圖5 粗(左)細(右)網(wǎng)格下不同激震頻率下波的傳播對比

圖6 粗(左)細(右)網(wǎng)格下不同激震頻率下觀測點位移幅值譜

譜通過自振分析得到重力壩前10自振階頻率以及各階模態(tài)的振型參與系數(shù)見表1。

圖7 壩體及地基有限元網(wǎng)格

模態(tài)12345678910頻率/Hz1.4103.1883.3725.6318.7168.86912.50113.34814.23916.212參與系數(shù)(x-分量)參與系數(shù)(y-分量)2.1660.498-2.0610.7680.086-1.7291.2110.5450.313-0.476-0.400-0.3550.2060.0380.0590.051-0.054-0.091-0.060-0.064

根據(jù)第3.2節(jié)所述，選取細網(wǎng)格模型模擬波的傳播過程，設(shè)震源位移時間函數(shù)為主頻4 Hz的Ricker子波，計算的加速度(歸一化)結(jié)果作為地震荷載輸入重力壩。波形和傅里葉譜如圖8所示。

圖8 輸入加速度波形和傅里葉譜

2.3.2 傳播介質(zhì)的影響

本節(jié)考慮傳統(tǒng)的均勻輸入方式和考慮地震波傳播過程的非均勻輸入對重力壩動力響應(yīng)結(jié)果的影響。如圖9(a)所示，輸入地震波的有限元節(jié)點分別編號為①～。圖9(b)給出了均勻輸入時各個節(jié)點的加速度時程，可以看到各點運動是一致的。圖9(c)給出了非均勻輸入時各節(jié)點的加速度時程，各個節(jié)點之間存在幅差和相差。非均勻輸入時分別考慮S波傳播速度Vs等于2.0 km/s、1.5 km/s、1.0 km/s三種情形。計算結(jié)果如圖10以及表2所示?？紤]幅差相差的非均勻輸入得到的最大主應(yīng)力結(jié)果明顯區(qū)別于均勻輸入：

(1)壩體上部，非均勻輸入的應(yīng)力結(jié)果小于均勻輸入，均勻輸入條件下，壩頸部上游面和下游面分別出現(xiàn)0.6 MPa和1.4 MPa的拉應(yīng)力。而非均勻輸入且波速等于2.0 km/s時，相應(yīng)部位的拉應(yīng)力分別降至0.4 MPa和1.0 MPa。當(dāng)波速進一步減小，即非均勻性增強時，壩體動力響應(yīng)減弱變得更加明顯。

(2)壩底處，非均勻輸入導(dǎo)致最大主應(yīng)力值有所增大，幅差相差在壩體引起的各點運動的不一致有增大壩底處動力響應(yīng)產(chǎn)生的拉應(yīng)力的趨勢。

另外，圖11給出了壩頂?shù)膭恿憫?yīng)過程，可以看到非均勻輸入條件下壩頂點的加速度幅值有所降低。

圖9 均勻輸入和非均勻輸入的加速度時程

(a)均勻輸入

(b)非均勻輸入，Vs=2.0 km/s

(c)非均勻輸入，Vs=1.5 km/s (d)非均勻輸入，Vs=1.0 km/s

圖11 不同輸入機制下壩頂?shù)募铀俣软樅酉蚍磻?yīng)時程

輸入機制最大主應(yīng)力極值/MPa壩頂加速度極值/m·s-2均勻輸入1.8910.40非均勻輸入Vs=2.0 km/s1.428.06Vs=1.5 km/s1.186.56Vs=1.0 km/s0.733.50

對于一般的工程場址，P波波速約在2～3 km/s，S波波速在1～2 km/s。由上述計算結(jié)果可知，實際工程中，地震波的行波效應(yīng)對壩體動力響應(yīng)會產(chǎn)生一定的影響。對于重力壩，表現(xiàn)為降低壩頂反應(yīng)，同時增大壩底附近的拉應(yīng)力(主要為壩底中部)。對于拱壩結(jié)構(gòu)，其跨度更大，加之峽谷兩側(cè)的地形影響，行波效應(yīng)更加明顯，對此需要做專門的三維分析，此處不予詳述。

2.3.3 震源位置的影響

前面的內(nèi)容都是假設(shè)震源位于重力壩上游側(cè)，波由上游方向入射。本節(jié)討論波從不同方向傳入的影響，設(shè)置位于重力壩下游側(cè)的震源，比較不同震源位置對重力壩動力響應(yīng)的影響。圖12表明，考慮靜動荷載疊加時震源位置對應(yīng)力分布結(jié)果影響不大，相同的波速下，上游側(cè)震源和下游側(cè)震源的應(yīng)力分布規(guī)律基本一致。但是從圖13給出的純動力反應(yīng)結(jié)果可以看出，震源位置對重力壩的純動力響應(yīng)有顯著的影響：當(dāng)震源位于上游側(cè)時，壩踵部位的拉應(yīng)力較大；而當(dāng)震源位于下游側(cè)時，壩趾部位的拉應(yīng)力較大。從整體應(yīng)力分布來看，整個受拉區(qū)有向震源一側(cè)偏移的趨勢。不同震源位置條件下，壩體上部應(yīng)力分布規(guī)律仍基本一致。

(a)上游側(cè)震源

(b)下游側(cè)震源

(a)上游側(cè)震源

(b)下游側(cè)震源

3 結(jié) 論

本文提出的基于地震學(xué)方法為大壩動力分析提供輸入的方法，比傳統(tǒng)輸入方法更具有科學(xué)性和說服力。通過此方法，可以考慮非均勻輸入以及震源和結(jié)構(gòu)相對位置的影響，并且通過重力壩算例分析得到了初步的結(jié)論：

(1)考慮非均勻輸入的重力壩整體動力響應(yīng)明顯低于均勻輸入，但是在壩底部中間小范圍區(qū)域會出現(xiàn)更高的拉應(yīng)力，且這種影響隨著結(jié)構(gòu)空間尺度的增加而增加，在大型水利工程的設(shè)計中應(yīng)該予以考慮。

(2)震源和結(jié)構(gòu)的相對位置使得重力壩拉應(yīng)力區(qū)趨向于靠著震源一側(cè)偏移，即震源側(cè)的動力響應(yīng)更高。對這一方面進一步的研究，可以為工程選址和大跨度結(jié)構(gòu)局部加固提供有意義的成果。

本文僅通過較小規(guī)模的二維算例，研究了重力壩的動力響應(yīng)，實現(xiàn)了從震源模擬到結(jié)構(gòu)響應(yīng)的過程，但是未考慮實際地形和速度結(jié)構(gòu)的影響，根據(jù)此思路的進一步的研究工作已經(jīng)取得成果：基于三維模型利用譜元法進行地震波傳播模擬，并且考慮地形地質(zhì)條件的影響，結(jié)合地震學(xué)高頻隨機模擬方法，合成寬頻帶地震波，以此作為地震動輸入用于計算峽谷地形中拱壩的地震響應(yīng)。