高面板壩地震動三維非一致輸入的波函數組合法①

姚 虞， 劉天云， 張建民

(清華大學水沙科學與水利水電工程國家重點實驗室，北京 100084)

高面板壩地震動三維非一致輸入的波函數組合法①

姚 虞， 劉天云， 張建民

(清華大學水沙科學與水利水電工程國家重點實驗室，北京 100084)

提出三維問題的波函數組合法求解不規則河谷的散射問題，應用到高面板壩動力有限元計算的地震動輸入中，可以滿足非一致輸入的要求，并回避人工邊界對外行波的反射問題。重點介紹三維問題的波函數組合法的理論推導和程序的算例驗證，并給出將該方法應用到高面板壩地震動輸入的基本步驟。

高面板壩； 三維； 非一致； 位移輸入

0 引言

高面板壩的抗震設計常常要進行動力有限元計算，在地震動輸入方面會面臨兩個問題：一是如何反映地震波的行波效應[1]，另一個是如何避免人工邊界對外行波的反射[2]。在目前常用的堆石壩動力分析程序中這兩個問題往往被忽略，代之以均勻輸入和固定邊界。對于300 m級的高面板壩，地震動作用的非一致性更加明顯，采用簡單的均勻輸入與實際有很大差別。

實現非一致輸入需要求解計算邊界即不規則河谷各點的振動過程。由于基巖的波阻抗(介質密度與波速的乘積)遠大于壩體，可以通過求解不規則空心河谷的散射問題得到邊界各點的振動過程。

圖1表示瑞利波或任意角度的壓縮波或剪切波入射到半空間不規則河谷的示意圖。要求解的問題是在任意入射波的作用下河谷邊界任意一點的震動過程。若能求解不規則河谷任意一點的震動過程，則該解可作為河谷內土石壩地震動響應有限元求解的輸入。這種輸入方式既考慮了地震動作用的非一致性，又回避了人工邊界對外行波的反射問題。

圖1 半空間不規則河谷入射波示意圖Fig.1 Incident wave at an irregular valley in a half-space

彈性波的散射問題已有很長的研究歷史。Sezawa在1927年完成了對球、圓柱和橢圓柱體在P波入射情況下的散射問題的研究，利用特殊波函數構造了問題的解。之后對散射問題的研究陸續展開，并形成了一些方法。Nishimura[3]解決了各向同性介質中球體空洞的動應力集中問題。同一時期，Pao和Mow等[4]還探討了圓柱體內含物的動應力集中問題。Baron等[5]率先對彈性介質中圓柱形空腔引起的脈沖散射問題進行了分析，使用積分變換和波函數展開法給出了壓縮波散射的解析解；對于球形內含物和空腔的瞬態效應，Mow和Miklowitz進行了求解[6-8]。Datta等[9-11]提出采用匹配漸進展開的方法分析了半空間內含柱形空洞對P波、SV波和SH波的散射問題，給出了遠場的漸近表達式及相應的數值解。Bostrom和Kristensson等[12]研究了半空間內球形空洞對Rayleigh波的散射問題。上世紀80年代，Pao等[13]運用廣義射線法討論了柱體的瞬態波散射。黎在良和劉殿魁[14]還建立了一個射線理論用來解各向異性介質中圓柱體對SH波的散射問題。

近20多年來，由于計算機技術的廣泛應用，用數值化方法求解彈性動力學問題取得了進展。實踐證明這些方法已成為解決波動問題的重要工具。如廖振鵬[15-16]的結構動力反應的數值模擬以及近場波動的數值模擬等。另外，對不規則的地形，Sanchez-Sesma等[17-19]先后采用邊界元法和波函數展開法對二維、三維問題進行了分析，Ray 等[20]、劉國利等[21]分別采用復變函數方法、積分變換和波函數展開相結合的方法給予分析。

1 方法原理介紹

由彈性力學中線彈性體位移表示的位移運動方程(Navier方程，不計體力)：

(1)

其中：λ、μ為拉梅常數；u為位移矢量。再考慮Helmholtz定理[22]：

(2)

其中：φ為速度場的標量勢；ψ為矢量勢。即任意矢量場可分解為一個無散場和一個無旋場之和。

令 u(1)=φ，u(2)=×ψ，代入Navier方程，得

(3)

故只需滿足

(4)

在頻域內，對于單一頻率的波，式(4)的解可寫為分離變量的形式，即

(5)

將式(5)帶入式(4)，可得

(6)

式(6)稱為亥姆霍茲方程，是頻域內的波動方程。在給定邊界條件下求解式(6)得到位移勢函數之后便可求解位移場。

對于不規則河谷散射問題，要求滿足河谷邊界和半空間自由表面面力為零的邊界條件。為得到滿足邊界條件的解，本文提出一種滿足全局邊界條件的波函數組合法，可以考慮包括瑞利波在內的不同地震波的不同角度入射情況。該方法通過組合體波函數和Lamb問題的解以滿足半空間全局的邊界條件，來求解不規則河谷的散射問題，是一種半解析的方法。目前已用該方法編寫了二維和三維河谷的計算程序，并分別利用半圓形和半球形河谷與理論解作對比，驗證了程序的正確性。

波函數組合法的基本思路是：

(1) 確定入射波自由場應力、位移；

(2) 位移、應力波函數展開；

(3) 確定自由表面波函數產生的面力；

(4) 在半無限空間自由表面反向疊加波函數產生的面力，求解Lamb問題計算散射邊界次生面力；

(5) 匹配散射邊界處自由場面力、波函數產生的面力與Lamb解產生的面力，形成散射邊界無面力狀態；

(6) 疊加散射邊界處自由場、波函數和Lamb解產生的位移，得到散射邊界的總位移。

對于二維波函數組合法的理論推導及程序驗證，作者已在文獻[23]中有所介紹。本文主要介紹三維波函數組合法。

2 三維波函數組合法理論推導

球坐標下，式(6)有分離變量形式的通解：

(7)

用波函數表達三維半空間的應力場和位移場的方法與二維的方法基本一致，其表達式可參見文獻[24]。三維問題的難點在于三維非軸對稱Lamb問題的求解。

由于波函數展開法不能滿足水平自由表面的邊界條件，存在殘余正應力σ與參與剪應力τ。為滿足自由表面條件，提出在表面處施加與殘余面力相反方向的虛擬面力-σ與-τ(圖2)。

圖2 自由表面殘余面力的消除Fig.2 Elimination of residual force on a free surface

波數域內，三維蘭姆問題的控制方程為：

(8)

其中φ和Ψ滿足

(9)

滿足以上兩式的三個勢函數：

(10)

(11)

將三個邊界應力進行組合，有

(12)

(13)

解得

(14)

3 算例驗證

三維程序利用半球形河谷的散射問題進行驗證。驗證結果參考文獻[18]，該文獻使用波函數展開法，將三倍半徑范圍的半空間自由表面納入計算域邊界中，忽略三倍半徑以外自由表面的邊界條件。文獻中提供了不同頻率P波垂直向上入射的結果，筆者用本方法計算出相應的結果以進行對比(圖3)。

圖3 P波垂直入射半球形河谷邊界位移放大系數Fig.3 Displacement amplification factor of hemispherical valley boundary under vertical incidence of P wave

從圖3中可以看出，采用本文方法計算的結果與文獻[18]計算的結果接近，從而驗證了三維程序的正確性。

4 應用方法

地震從遠處傳到壩址建基面，輸入壩體內部，引起壩體的變形和破壞。為模擬地震輸入情況，需要以下條件：(1)入射方位與入射角；(2)入射波形式；(3)地震時程。

首先，將地震時程S(t)(一般是間隔0.02 s的離散記錄)通過FFT獲得功率譜曲線W(ω)；

第二，將功率譜曲線W(ω)等間隔離散，如：Δω=0.5，獲得一系列諧波Wi；

第四，對所有離散頻點進行第三步的計算，可獲得壩址建基面散射場位移譜Ws(x，ω)；

第五，應用逆FFT，獲得壩址建基面散射場的位移時程ws(x，t)。

(15)

式中下標s為結構結點量；b為壩基結點量；us為結構結點位移向量。

5 結論

本文提出的三維波函數組合法可解決三維半空間內任意空心體的散射問題。其基本思路是將地震波轉換到頻域之后離散，對每種頻率的波通過組合波函數和Lamb問題的解，使得邊界零面力條件得到滿足，再由波函數和Lamb問題的解得到邊界的頻域內的位移，最后轉換回時域，得到邊界的位移時程。由于堆石體的波阻抗相對于基巖很小，可將該邊界位移時程作為壩體動力有限元計算的輸入。

通過對半球形河谷在P波垂直入射條件下的散射問題進行編程計算，初步驗證了理論的正確性。將該方法應用到高面板壩地震動輸入中，可以考慮地震動作用的非一致性，并回避人工邊界對外行波的散射問題。

References)

[1] 王艷.非一致地震動場數值方法研究及在結構動力分析中的應用[D].北京：清華大學，2007.WANG Yan.Research on the Numerical Method for Asynchronous Seismic Wave Motions and its Application in Dynamic Analysis of Structures[D].Beijing：Tsinghua University，2007.(in Chinese)

[2] 李彬.地鐵地下結構抗震理論分析與應用研究[D].北京：清華大學，2005.LI Bin.Theoretical Analysis of Seismic Response of Underground Subway Structures and Its Application[D].Beijing：Tsinghua University，2005.(in Chinese)

[3] Nishimura G.Dynamic Problem of stress Concentration in the Vicinity of a Spherical Matter Included in an Elastic Solid Under Dynamical Force[J].Journal of the Faculty of Engineering，1955，24：101.

[4] Pao Y H，Mow C C.Dynamic Stress Concentration in an Elastic Plate with Rigid Circular Inclusion.No.SR-41[R].The MITRE Corporation，Bedford，Massachusetts，1961.

[5] Baron M L，Mattews A T.Diffaction of a Pressure Wave by a Cylindrical Cavity in an Elastic Medium[J].Journal of Applied Mechanics，1961，28(3)：347-354.

[6] Mow C C.On the Transient Motion of a Rigid Spherical Inclusion in an Elastic Medium and Its Inverse Problem[J].Journal of Applied Mechanics，1966，33(4)：807-813.

[7] Mow C C.Transient Response of a Rigid Spherical Inclusion in an Elastic Medium[J].Journal of Applied Mechanics，1965，32(3)：637-642.

[8] Miklowitz J.Scattering of a Plane Elastic Compressional Pulse by a Cylindrical Cavity[C]// Proceedings of the 11th International Congress of Applied Mechanics.Munich，1966：469-483.

[9] Datta S K，Akily N E.Diffraction of Elastic Waves by Cylindrical Cavity in a Half-space[J].The Journal of the Acoustical Society of America，64(6)：1978：1692-1699.

[10] Datta S K，Shah A H.Truell R，Scattering of SH Waves by Embedded Cavities[J].Wave Motion，1982，4(3)：265-283.

[11] Datta S K.Scattering of Elastic Waves[M]// Nematnasser S, edit.Mechanics Today.New York:Pergamon Press, 1978, 4:149-208.

[12] Bostrom A，Kristensson G.Scattering of a Pulsed Rayleigh Wave by a Spherical Cavity in an Elastic Half Space[J].Wave Motion，1983，5(2)：137-143.

[13] Pao Y H，Ku G C.Application of the Theory of Generalized Rays to Diffraction of Transient Waves by a Cylinder[J].Wave Motion，1983，5(4)：385-398.

[14] 黎在良，劉殿魁.各向異性介質中圓柱體對SH波散射的射線理論[J].地震工程與工程振動，1987，7(1)：1-8.LI Zai-liang，LIU Diam-kui.Ray Method of SH Wave Propagation in an Anisotropic Medium with a Circular Cavity[J].Earthquake Enginelring and Engineering Vibration，1987，7(1)：1-8.(in Chinese)

[15] 廖振鵬.近場波動的數值模擬[J].力學進展，1997，26(2)：193-216.LIAO Zhen-peng.Numerical Simulation of Near-field Wave Motion[J].Advances in Mechanics，1997，26(2)：193-216.(in Chinese)

[16] 廖振鵬.結構動力反應的數值模擬[M]//中國科學技術前沿(第四卷).北京：高等教育出版社， 2001：427-458.LIAO Zhen-peng.Numerical Simulation of Structure Dynamic Response[M]//Science and Technology at the Frontier in China (Vol.4).Beijing：Higher Education Press，2001：427-458.(in Chinese)

[17] Sanchez-Sesma F J，Herrera I，Aviles J.A Boundary Method for Elastic Wave Diffraction：Application to Scattering of SH Waves by Surface Irregularities[J].Bulletin of The Seismological Society of America，1982，72(2)：473-490.

[18] Sanchez-Sesma F J.Diffraction of Elastic Waves by Three-dimensional Surface Irregularities[J].Bulletin of the Seismological Society of America，1983，73(6)：1621-1636.

[19] Sanchez-Sesma F J，Miguel A B，Ismael H.Surface Motion of Topographical Irregularities for Incident P，SV，and Rayleigh Waves[J].Bulletin of the Seismological Society of America，1985，75(1)：263-269.

[20] Han R PS，Yeh K Y，Liu G L，et al.Scattering of Plane SH-waves by a Cylindrical Canyon of Arbitrary Shape in Anisotropic Media[J].International Journal of Engineering Science,1992,30(12): 1773-1787.

[21] 劉國利，劉殿魁.位移階躍SH波對半圓形凹陷地形的散射[J].力學學報，1994：26(1)：70-80.LIU Guo-li，LIU Dian-kui.The Scattering of Displacement Step SH-waves by a Semi-cylindrical Canyon[J].Acta Mechanics Sinica，1994，26(1)：70-80.(in Chinese)

[22] [美]阿肯巴赫，編.徐植信，洪錦如，譯.彈性固體中波的傳播[M].上海：同濟大學出版社，1992.Achenbach J D，edit.Translated by XU Zhi-xin，HONG Jin-ru.Wave Propagation in Elastic Solids[M].Shanghai：Tongji Univerity Press，1992.(in Chinese)

[23] 姚虞，劉天云，張建民.求解二維彈性半空間散射問題的波函數組合法[J].地震工程學報，2013，35(1)：109-113.YAO Yu，LIU Tian-yun，ZHANG Jian-min.The Method by the Combination of Eave Functions to Solve the Problem of Scattering in 2D Elastic Half Space[J].China Earthquake Engineering Journal，2013，35(1)：109-113.(in Chinese)

[24] 鮑亦興，毛昭宙.彈性波的衍射與動應力集中[M].北京：科學出版社，1993.PAO Yih-hsing，MOW Chao-chow.The Diffraction of Elastic Waves and Dynamic Stress Concentrations[M].Beijing：Science Press，1993.(in Chinese)

Wave Function Combination Method Used for Three-dimensional Non-uniform Input of Earthquake Motion in High Concrete-faced Rockfill Dam

YAO Yu， LIU Tian-yun， ZHANG Jian-min

(StateKeyLaboratoryofHydroscienceandEngineering，TsinghuaUniversity，Beijing100084，China)

The scattering problem of elastic waves is of great engineering significance，especially in earthquake engineering.Frequently，the response of the scatterer is needed to evaluate the safety of structures during earthquakes.The ground motion can be used as dynamic input for the finite-element method calculation，which is a way to realize the non-uniform input.This article puts forward the wave function combination method for a three-dimensional (3D)problem to solve the irregular valley scattering problem.The application of the earthquake input to the dynamic finite-element calculation for a high-face rockfill dam can meet the requirement for non-uniform input and，at the same time，avoid the problem of the reflection of the outward wave by the artificial boundary.The theoretical derivation of the wave function combination method for the 3D problem and the validation of the program are introduced，as are the basic steps needed to apply this method to the earthquake input.

high concrete-faced rockfill dam； three dimension； non-uniform； displacement input

2014-08-20

國家自然科學基金項目(51038007,51079074,51129902)

姚 虞(1989-)，男，博士生，研究課題為高壩抗震計算.E-mail：yaoyutsinghua@126.com

P315.31

A

1000-0844(2015)02-0324-005

10.3969/j.issn.1000-0844.2015.02.0324