林 皋

(大連理工大學海岸與近海工程國家重點實驗室，大連116024)(大連理工大學工程抗震研究所，大連116024)

地下結構地震響應的計算模型1)

林 皋2)

(大連理工大學海岸與近海工程國家重點實驗室，大連116024)(大連理工大學工程抗震研究所，大連116024)

地震時地下結構在圍巖的約束作用下發(fā)生變形，其動態(tài)特性與地面結構有很大不同.自二十世紀七八十年代以來，地下結構抗震設計與研究取得了很大進展.總的看來，工程設計中普遍采用的計算方法與設計導則大都建立在比較簡單假定的基礎上，實際的巖土介質條件都是十分復雜的.地下結構抗震研究的近期成果則表現在對地下結構動力分析中的波動散射問題提出了波函數展開法以及邊界積分方程方法等多種計算方法.但計算相對復雜，在工程設計中的推廣應用有一定困難.此文致力于地下結構計算模型的改進，使之具有良好的計算精度與效率，又便于工程應用.為此，提出了一種地下結構抗震響應分析的新的計算模型.模型具有較廣泛的適應性，可以進行河谷、孔洞、地下鐵道、隧洞等地下結構的散射與繞射分析.對于復雜層狀的地質條件，提出了格林函數求解簡便而有效的方法.數值算例論證了方法的精度和效率.

地下結構，地震響應，波的散射，格林函數，層狀半空間地基

引言

地下工程結構在城市生產和生活中發(fā)揮著重要作用.地下結構承擔著城市公共設施的許多重要功能，在通訊、能源供應(供電、供氣)、供水排水、交通運輸(地下鐵道、海底隧道)、國防和民防(防空、防爆)等許多方面，發(fā)揮著現代工業(yè)生產和城鎮(zhèn)生活中大動脈的作用.地下結構還常用作儲油庫、儲氣罐；地下結構也較常用來建設地下電站的廠房.隨著國民經濟和城市建設的發(fā)展，地下工程結構的高度和跨距也愈來愈大，這對地下結構的抗震設防提出了更高的要求.

歷史上多次大地震中地下結構的破壞造成城市供電、供水、通訊的中斷，引起國民經濟的重大損失.著名事例[1]如：1906年美國舊金山大地震(震級 8.2)，由于供水中斷造成火災蔓延，破壞區(qū)域達12.2km2(4.7平方哩)，大火燃燒數日，所造成的損失占地震總損失的80%；1923年日本關東大地震(震級8.2)，在東京80%的城區(qū)中引起火災，并有25條鐵路隧洞受有破壞，橫濱市由于排水總管破壞，除火災外又引起次生洪水災害.由于供電、供水等地下管線的破壞對民眾的健康和安全關系重大，因此，這些地下工程結構又稱為生命線工程.

20世紀80年代以前，國際上還缺乏地下結構或生命線工程抗擊地震荷載的設計規(guī)范化條款.對一些重要的地下管道，日本采用震度法(地震系數法)進行設計，這并不能反映地下結構所遭受地震作用的實際情況.美國舊金山大地震(1906)、日本關東大地震 (1923)、新潟地震 (1964)、宮城縣地震(1978)和日本海中部地震(1983)中管道等地下結構遭到很大程度的破壞，給城市的生產和生活造成了很大影響，促進了地下結構抗震研究的發(fā)展.這一時期(主要是二十世紀七八十年代)發(fā)表了大量研究論文[2]，涵蓋地下結構地震震害和地震響應觀測、影響因素、地下結構地震反應分析和計算模型的建立、設計要領和設計準則等各個方面.與此同時，日本進行了地下結構地震響應的觀測和模型試驗研究，其中比較有代表性的是田村和久保等的工作[37].

1977年 8月美國土木工程師學會生命線地震工程學科技委員會召開了地震工程學專業(yè)會議[2]，1981年第7屆世界地震工程會議期間舉行了生命線地震工程學最新進展委員會的會議.對地下結構或生命線工程抗震研究的情況進行了總結分析，深化了地下結構抗震特性的認識.

根據會議總結和會議相關論文的論述[27]，地下結構的地震響應具有以下特點：

(1)大多數現場觀測資料表明，地震時地下管線和水下隧道無論軸向和側向均隨圍巖發(fā)生著相同的運動.

(2)觀測到的地下管線和水下隧道的軸向應變遠較彎曲應變顯著.管道轉彎處的彎曲變形與直線段處的彎曲變形具有相同的量級.隧道的軸向剛度相對較大，對圍巖變形起一定的限制作用.

(3)地下管線和隧道地震時產生的慣性力，對結構自身的地震響應只產生非常小的影響.

(4)地下管線的存在對圍巖地震動的特性和擾動只有微弱的影響.

本文只討論由地震波傳播所引起的地下結構的地震響應分析.

1 地下結構地震響應計算模型的發(fā)展

由于地震時地下結構受圍巖介質的約束作用而發(fā)生變形，引起震害，因此地下結構的地震響應分析應包括以下三方面的內容：

(1)圍巖介質，特別是復雜巖土介質中自由場的地震波動特性.由于地下結構的存在對總體地震波動場的擾動一般較小，因此，可首先研究地下結構建造以前自由場的波動特性.

(2)地震時地下結構周圍的地震波場特性.包括自由場和散射場兩部分.重點需要計算的是由于地下結構不同的幾何形狀和動力剛度特性引起的散射波場特性.

(3)地震時地下結構和周圍巖土介質間的動力相互作用特性.由于地下結構的質量和剛度與周圍巖土介質相比是一個微量，因此地下結構自身的慣性力對其地震變形只起很微小的作用.但是地下結構的剛度，尤其是軸向剛度仍然對圍巖的地震變形發(fā)生足夠的反作用，這使地震時地下結構和圍巖介質間產生有差異的地震變形.

二十世紀七八十年代以來，地下結構的抗震分析取得了很大的進展，但是，對上述三方面的研究深度，仍然不能使人滿意.下文將對計算模型、規(guī)范標準以及文獻中較廣泛采用的有限元數值計算方法等方面，根據筆者的體會，作簡要闡述.

自1972年Trifunac[8]發(fā)表關于半圓形河谷在SH平面波作用下的散射效應的論文以來，文獻中關于地下結構地震響應的計算模型方面發(fā)展了很多數值計算方法.主要針對不同形狀的河谷、沖積層、地下孔洞、地下管道、地下隧道等在SH波、SV波、P波和瑞利波作用下的散射效應提出了許多計算模型與方法.具體可參見有關文獻闡述[912].比較有代表性并獲得較廣泛應用的數值計算方法主要有波函數展開法和邊界積分方程方法兩種.

波函數展開(WFE)法主要適用于均質空間或半空間散射問題的求解[1316].大多數求解的為二維散射問題，三維問題的求解局限于軸對稱散射[17]，個別的也求解了非軸對稱散射問題[18].

邊界積分方程(BIE)方法又稱邊界元方法，其特點是只需進行散射體表面的離散，使問題降階一維；同時滿足無窮遠處的輻射條件.但需獲得計算域的基本解，增加了問題求解的困難，并增加了計算工作量.邊界積分方程方法又可區(qū)分為直接BIE方法和間接BIE方法兩種[1925].

將邊界積分方程方法與Haskell和Thomson[26-27]提出的傳遞矩陣方法或Kausel[28]提出的剛度矩陣方法相結合，可以求解層狀地基的格林函數，從而求解層狀半空間的散射問題.

文獻中也有的將邊界積分方程法與格林函數的離散波函數展開法相結合進行散射問題的求解[2930].一般情況下，邊界積分方程不適于求解復雜不均勻介質的散射問題，而有限元法則具有較廣泛的適應性，對計算域介質的不均質問題處理非常方便，將兩者相結合可以發(fā)揮各自的優(yōu)越性.近場采用有限元法分析，遠場輻射邊界則采用邊界積分方程求解[31].同理，也可以將有限元法與波函數展開法結合求解散射問題[32].

近年來，中國學者將前蘇聯(lián)學者Muskhelishvili所提出的求解彈性力學問題的復變函數法[33]應用于求解散射問題也取得了一定效果[3435].此外，文獻中還提出了一些其他方法，本文將不再贅述.

二十世紀七八十年代以后，一些國家的抗震設防標準，主要是土木工程的抗震規(guī)范和核電站的抗震規(guī)范,包含了地下結構的抗震設計準則.但相對于房屋和橋梁等地面結構來說，條文較簡略，難以滿足實際抗震設計的要求[3638].

現有規(guī)范中關于地下結構的抗震分析方法，主要針對地下埋設管道的計算.對于無限長管，假設管道的地震變形與均勻無限介質中波傳播產生的變形完全相同[3940]，據此計算管的軸向應變與彎曲應變.這一做法在核電結構等規(guī)范[41]中應用較為普遍.

一般規(guī)范標準中，很少有考慮復雜地基介質中的波傳播特性對地下結構地震響應影響的規(guī)定.只有日本沉埋隧道的抗震設計規(guī)程提出了考慮沿線地基地質條件變化的地震動輸入計算模型[38].其中用多質點彈簧阻尼體系來模擬圍巖介質中地震波的傳播.

文獻中應用有限元方法進行地下結構的抗震分析較普遍.將地下結構與圍巖介質劃定一范圍作為計算域，用有限元離散進行分析.這種方法可以考慮地下結構的各種幾何形狀和圍巖介質的不同特性.但是由于缺乏有關地下結構抗震設計的導則和標準，因此，在網格的劃分和布置方面比較自由，不同研究者的計算結果離散性較大.此外，在計算域范圍、邊界條件的選擇等方面都含有一定的任意性.例如，側面邊界有的采用能量傳遞邊界，有的采用黏性邊界，有的則采用水平向可自由滑行的邊界；底面邊界有的采用黏性邊界，有的則采用一定深度處的固定邊界.地震動輸入的模型也有一定的任意性.總體上，對自由場地震動的分析，散射場的計算以及地下結構與圍巖介質的動力相互作用分析方面是否與實際情況相符都缺乏必要的檢驗，計算精度及其可靠性缺乏論證.將這樣的計算結果用來進行地下結構的抗震設計與抗震安全評價是難以令人滿意的.

地下空間的利用發(fā)展很快，適應各種需要的地下結構的高度、跨距等都逐漸增加，這對地下結構的抗震設計提出了新的要求.提高地下結構抗震研究的水平勢在必行.現有的地下結構地震觀測的資料還比較缺乏，并具有一定的局限性，多偏重于可能發(fā)生的地震破損形態(tài)，例如，沿隧道全線或斷面內若干測點的變形和加速度等，對地震作用下不同地下結構的變形特點尚缺乏深入認識.用以進行新形勢下地下結構的地震響應分析還有待進一步深入研究.

2 建議的地下結構地震響應分析的基本計算模型

本文的目的在于尋求一種理論上容易理解，計算不復雜，具有較好的計算精度，同時又便于在實際工程中應用的新方法，以便為地下結構抗震設計導則和規(guī)程的制訂提供一定的參考.

下面從土與結構動力相互作用 (dynamic soilstructure interaction)的基本理論出發(fā)來建立地下結構地震響應分析的方程[42].應用子結構分析方法，土與結構動力相互作用的計算方程可表示如下(圖1)

式中,ω為激勵頻率,ζ表示滯回阻尼系數.如需考慮高階阻尼影響時，S的表達式可參見文獻[43],進行地震作用分析時Ps=0，結構和無限地基的相互作用力Rb可表示為

圖1 結構與無限地基系統(tǒng)Fig.1 Structure-unbounded soil system

將式(2)和式(3)代入式(1)可得

在圖 1(b)中，將地基被挖去的部分表示為結構，則式(5)中矩陣Sbb-SbsSSsb表示地基被挖去部分的剛度(指被挖去部分凝聚于邊界b上的剛度).這樣式(5)左方括號內的矩陣+即表示地基未開挖前b邊界的動力剛度；而即表示地基未開挖前b邊界的自由場地震動[44].

圖2 散射體的計算模型示意圖Fig.2 Computational model for surface and subsurface irregularities and inhomogeneities

需要指出，式(8)的應用需要根據情況作出相應的變化.例如對于沖擊層填充的河谷(圖3)，計算公式相應地轉換為

圖3 填充河谷散射體Fig.3 Scattering by alluvial valley

圖4 散射場計算邊界的選擇Fig.4 Selection of the outer boundaries for scattering analysis

可以指出，這里所建立的計算公式，既可適用于二維結構，也可適用于三維結構；既可適用于均質半

3 數值驗算

本文將以均質半無限空間中的散射問題為主，通過算例檢驗計算模型的有效性.為計算均質半無限空間的動力剛度,比例邊界有限元方法(SBFEM)[45]是有效的方法.SBFEM為半解析半數值性的求解方法，其徑向可獲得解析解，環(huán)向則達到有限元的計算精度.由于只需要進行邊界面的離散，問題的維數降低一階；并可自動滿足無限遠處的輻射條件.因此SBFEM兼有FEM和BEM的優(yōu)點，也不需要基本解，計算簡便.

3.1 半圓形河谷的散射

如圖5所示,河谷邊界ABC,計算邊界DEFG，相應計算域大小為4a×2a(a為河谷半徑).設圍巖介質材料特性參數為：質量密度ρ0=1.5，剪切波速cs0=1，泊松比ν0=1/3；河谷中填充材料參數為：質量密度ρ1=1，剪切波速cs1=0.5，泊松比ν0=1/3.無量綱頻率 η =(ωa)/(πcs)，波長 λ =2/η.計算采用式(9).網格離散如圖5.各子域邊界采用了3節(jié)點二次單元離散，共144單元，288節(jié)點.相應DEFG邊界設48單元，97節(jié)點.入射波從底部垂直輸入，計算結果精度高，與Luco和de Barros的解[46]相符性良好 (圖 6).Luco和 de Barros，以及 Dravinski和Mossessian[47]采用的都是間接邊界積分方程(IBIE)方法.

圖5 半圓形河谷的計算圖形Fig.5 Schematic view of semi-circular valley

圖6 半圓形沉積河谷散射問題的檢驗Fig.6 Verificatio of the scattering by semi-circular valley

圖6 半圓形沉積河谷散射問題的檢驗(續(xù))Fig.6 Verificatio of the scattering by semi-circular valley(continued)

3.2 地下洞室襯砌的地震響應

研究地下隧道襯砌的地震響應.de Barros和Luco[48]求解了這一問題，采用間接邊界積分方程方法模擬外圍半空間土介質的作用，再結合Donnell的殼體理論計算隧洞襯砌外半徑的波動響應.子域括號內的數字代表該子域相似中心的坐標.隧道結構如圖7(a)所示.隧道埋深H=8.33r1=7.573a.襯砌材料質量密度 ρ0=2.24×103kg/m3，彈性模量E0=1.6×1010N/m2，泊松比 ν0=0.2.圍巖地基介質材料質量密度ρ1=7.665×103kg/m3，彈性模量 E1=6.9×108N/m2，泊松比 ν0=0.45.位移按r=a進行規(guī)格化.假設SV波垂直入射，無量綱頻率η=ωa/(πcs)=0.132.計算的襯砌響應如圖8所示.可見與文獻解的相符性較好.

圖7 均質半無限地基中的隧道襯砌Fig.7 Tunnel and its liner embedded in uniform half-space

圖8 SV波垂直入射作用下隧道襯砌外表面的動力響應Fig.8 Dynamic response at the outer radius of the tunnel liner subjected to vertically incident SV wave

4 復雜地基條件下地下結構的地震響應分析

其中的水平分層地基比較常見，下文將重點研究.對于分塊不均質地基，當各塊之間可以用同一相似中心加以描繪時(圖9(b))，提出采用比例邊界有限元的處理方法[49]，對近場含不規(guī)則不均質捕虜體等復雜地基(圖9(c))則提出采用分步阻尼影響抽取法[50]進行求解.所提出的方法雖然是結合結構與地基動力相互作用分析進行闡述，但可適用于地下結構的抗震分析.

圖9 不均質地基的計算模型Fig.9 Computational model for inhomogeneous unbounded soil

對各向同性和非各向同性層狀半空間不均質地基的動力剛度及格林函數的求解，作者研究組從2012年起提出了基于積分變換的方法，并逐步得到完善[43,5153]，并可適用于橫觀各向同性的層狀介質地基的求解.

下面介紹其基本思想和主要方程.

4.1 層狀介質的波動方程

在圓柱坐標系下，層狀介質(圖10)中任一層以位移所表示的波動方程具有如下形式

圖10 水平層狀半空間Fig.10 Horizontally layered strata overlying half-space

對式(10)進行Fourier-Bessel變換，得出波數域的表達式.變換式如下

相應的反變換式為

式中k代表波數，其余系數的表達式為

式(15)中的兩個子式，上面的子式對應于荷載與變形相對于X軸對稱的情況；下面的子式則為荷載與變形相對于X軸成反對稱的情況.其中Jn(kr)為柱貝塞爾函數.Green函數計算時，針對圓形單元面上作用的均布水平向荷載和豎向荷載，可選擇n=1和n=0兩種情況(圖11)與之對應.當圓形單元半徑趨于極限Δr→0時，均布荷載轉化為集中水平力和集中豎向力作用的情況.即，在實際應用時，主要可限于n=1和n=0的兩種情況.

圖11 Green函數計算荷載單元Fig.11 Loaded element for evaluation of Green’s function

按式(12)進行Fourier-Bessel變換后得出頻率--波數域內解耦的SV-P波動方程(變量為ur.uz)和SH波動方程(變量為uθ)[5152]，也可參見文獻[28].

對橫觀各向同性介質中的波動，解耦方程同樣成立，但相應系數改變如下

式中,上標m=1代表SV-P波動，m=2代表SH波動.Im為單位陣，m=1為2階單位陣，m=2為1階單位陣.

對各向同性介質(方程(17)與方程(18))

對于橫觀各向同性介質

為了進行微分方程(21)的求解，引入如下定義的應力矢量pm作為位移um的對偶變量

可以證明pm滿足關系式[54]

以下的求解過程對m=1和m=2都適用，所以在下文中可將m省略.

利用式(26)，微分方程(21)可化為如下對耦形式

矩陣A,B,C,D和矩陣Kij(i,j=1,2)的變換關系如下

在狀態(tài)空間，對耦微分方程(27)可以聯(lián)立求解

齊次一階線性常微分方程(29)的解為指數函數.

采用鐘萬勰提出的精細積分方法[55]可以將解進行精確計算.

對層狀半空間中的任意一層，或其中的一層，設層厚為η=zb-za，則層上、下兩端的位移和應力將成立如下關系

或寫為

式中

為達到高精度的計算效果，T按如下方式進行計算[55]

式中，τ=η/b，b可取為任意整數.當取b=2N,N=20時，這相當于將η分成2N=1048557個微細薄層用式(34)計算，從而可以達到任意希望的精度，而計算工作量不大.文獻[55]中給出了計算T的遞推公式，應用方便.

為了進行層的結合和計算Green函數，解也宜表示成如下對耦變量的形式

式中矩陣F,G,Q也由T求出

對于水平分層介質，采用積分變換方法，可以求得波動方程在波數域的簡便形式(式(17)和式(18)或式(19)和式(20)).進一步引入對偶變量，又可將其化為一階齊次線性常微分方程(29)，其解可采用精細積分方法達到任意希望的精度.最后將解表示成對耦形式(36)，可便于分層之間的結合和Green函數的求解，說明如下.

4.2 相鄰層的結合

當將相鄰的兩層：層1[za,zb]和層2[zb,zc]進行結合時，應用式(36)，可有

從式(38)中消去ub和pb，即可求得合并后新層c[za,zc]的力與位移關系

式中

4.3 Green函數的計算

層狀半空間內部點Green函數的計算是地下結構抗震分析的核心部分.文獻中廣泛采用的剛度矩陣方法(sti ff ness matrix method)需要進行大型矩陣求解[28].這里提出對偶變量的轉換方法，計算簡捷方便，同時又準確高效.計算圖形如圖12和圖13所示.圖中源點(source)表示力的作用點，接受點(receiver)表示位移的計算點.已知源點和接受點的位置，可將層狀半空間劃分為3個子集，應用4.2節(jié)中闡述的方法可以分別建立各個子集形如式(36)所示的對耦變量方程，并相應求出其Fi,Qi,Gi矩陣(i=u,m,l，其中u表示上部，m表示中部，l表示下部)方程.Green函數的計算，基本上可劃分為圖13所示的6種工況.應用式(36)，即可求得各種工況的Green函數計算的矩陣表達式[53].

圖12 層狀半空間的Green函數Fig.12 Green’s function for multilayered half-space

圖13 層狀半空間不同工況的Green函數計算Fig.13 Computation of Green’s function for various cases of boundary conditions of multilayered half-space

下面給出不同工況的Green函數的矩陣表達式.下標a,b表示源點和接受點所在位置.

(1)工況(a)和(b)層狀地基表面點的Green函數

工況(a)下部為剛性地基

工況(b)下部為彈性半空間

(2)工況(c)和(d)下部固定的層狀地基內部點的格林函數

工況(c)

(3)工況(e)和(f)下部為彈性半空間的層狀地基內部點的格林函數

工況(e)

所有矩陣的階數對SV-P波為2×2，對SH波為1×1，計算很簡便.最后可將波數域各種工況格林函數的計算，寫成如下形式的統(tǒng)一的表達式(參見式(17)和式 (18)或式(19)和式(20),ur,uz和uB是解耦的，可分別求解

對圖 11表示的格林函數計算的荷載形式，例如，半徑為Δr的圓形單元上分布的Z向和X向均布荷載和，其相應的波數域的表達式為

在上述式中，令 Δr→0，πΔr2z→z或πΔr2x→x就可得出集中力z或x作用的表達式.再將式(53)、式(54)代入式(51)、式(52)并利用(13)進行反變換，就可得出圓柱坐標系物理域中格林函數的表達式，可統(tǒng)一寫成

4.4 數值檢驗

研究均質半無限地基上表面彈性層中半圓形河谷的散射(圖14).散射計算需要求解層狀半空間中的格林函數.計算中介質材料特性參數假設為：均質半空間 ρ0=4/3，r0=0.333，cs0=2；表面彈性層 ρ1=1，r1=0.333，cs1=1；河谷中沖積層ρ2=0.667，r2=0.333，cs2=0.5.幾何尺寸：河谷半徑r，表面彈性層厚h=1.5r.無量綱頻率η=ωa/(πcs).

在SV波垂直入射情況下半圓形充填河谷地表位移幅值計算結果如圖15所示，可見與文獻[56]的相符性很好.

圖14 半無限地基上表面彈性層中的半圓形河谷的散射Fig.14 Scattering of a semi-circular valley embedded in the surface layer overlying elastic half-space

圖15 SV波垂直入射情況下層狀半空間中半圓形河谷地表位移幅值Fig.15 Amplitude of displacement response of the semi-circular valley embedded in layered half-space under vertically incident SV wave

5 地下結構地震響應分析的其他有關問題

5.1地下結構地震響應的時域分析

地下結構地震響應的時域計算應用運動方程(1)，關鍵是相互作用力 R=或R=的計算.由于無限域中R是激勵頻率的函數，需將其轉換成時域的表達式，如式(51)所示(假設初值條件

式中，S∞(ω)代表無限域的動剛度(ω)或ω)；S∞(t)代表單位位移脈沖響應函數

式(58)只具有形式上的意義.由于S∞(ω)非平方可積，含有奇異性，式(58)不能直接用于脈沖響應函數的計算[57]，需進行必要的處理，將另文闡述.

5.2 地下結構地震響應的簡化計算

地下結構也可采用簡化的數值方法進行初步的抗震分析.以圖16所示的地下結構為例來加以說明.圖中虛線范圍以內為計算域，域內可采用有限元等數值方法進行分析.根據第2節(jié)計算模型的推導，虛線所表示的邊界應滿足相互作用力的條件R=或R=.也可選擇較大的計算域，并設立能量傳遞邊界，通過邊界進行自由場的地震動輸入.文獻中一般的能量傳遞邊界主要建立在均質無限地基假定的基礎上，不適于復雜地基條件下地下結構的地震響應分析.如采用時間、空間非耦合的局部人工透射邊界，計算準確度并不理想，需要采用高階局部人工邊界，但高階局部邊界計算精度的提高是以計算自由度的增加為代價的.我們提出的層狀地基耦合性(或全域性)時域傳遞邊界[58]，不僅計算精度高，而且計算效率也高，可供參考.

圖16 地下結構抗震分析的簡化模型Fig.16 A simplifie model for seismic analysis of under-ground structures

6 結束語

提出一種基于結構與地基動力相互作用方程的地下結構抗震分析模型，可以很方便地進行地震波散射場與地下結構和圍巖動力相互作用的計算.模型具有廣泛的適用性.可以考慮層狀半空間等復雜地基的影響.計算準確、高效.本文以二維均勻半空間中地下結構的地震響應為例，闡述地下結構抗震計算模型的設想.但是對復雜層狀半空間地下結構的地震響應只提出了計算方法和相關方程還來不及展開，以后將陸續(xù)進行闡述.

致謝 博士生李志遠、韓澤軍和課題組李建波、胡志強等在本文完成中發(fā)揮了重要作用,在此表示感謝.

1 Ariman T,Gregory E Muleski.A review of the response of buried pipelines under seismic excitations.Earthquake Engineering&Structural Dynamics,1981,9(2)：133-152

2 Wang LRL,O’Rourke MJ.Overview of buried pipelines under seismic loading.Journal of Technical Councils of ASCE,1978,104(TC1)：121-130

3 Okamoto S,Tamura C.Behavior of subaqueous tunnel during earthquakes.International Journal of Earthquake Engineering&Structural Dynamics,1973,1(3)：253-266

4 Tamura C,Okamoto S,Kato K.Observations on dynamic strains of submerged tunnel during earthquakes.Bulletin of Earthquake Resistant Structure Research Center,University of Tokyo,1972,6：1-21

5 Tamura C.Dynamic behavior of submerged tunnel during earthquake.Transactions of the Japan Society of Civil Engineers,1975,7：197-200

6 Kubo K.Behavior of underground water pipes during an earthquake//Proceedings of the Fifth World Conference on Earthquake Engineering,Rome,Italy,1974,569-578

7 Kubo K,Katayama T,Ohashi A.Present state of lifeline earthquake engineering in Japan//The Current State of Knowledge of Lifeline Earthquake Engineering,Rome,Italy,1974,569-578

8 Trifunac MD.Scattering of plane SH waves by a semi-cylindrical canyon.Earthquake Engineering&Structural Dynamics,1972,1(3)：267-281

9 Luco JE,De Barros FCP.Seismic response of a cylindrical shell embedded in a layered viscoelastic half-space.I：Formulation.Earthquake Engineering&Structural Dynamics,1994,23(5)：553-567

10 Wong HL.E ff ects of surface topography and site conditions//Proc.of the EPRI Workshop on Strong Ground Motion Simulation and Earthquake Engineering Applications,Los Altos,California,April 30-May 3,1984

11 Sanchez-Sesma FJ.Site e ff ects on strong ground motion.Soil Dynamics and Earthquake Engineering,1987,6(2)：124-132

12 Von Thun JL.Earthquake Engineering and Soil Dynamics II—Recent Advances in Ground-Motion Evaluation.[s.l.]：ASCE,1988

13 El-Akily N,Datta SK.Response of a circular cylindrical shell to disturbances in a half-space.Earthquake Engineering&Structural Dynamics,1980,8(5)：469-477

14 Wong KC,Shah AH,Datta SK.Dynamic stresses and displacements in a buried tunnel.Journal of Engineering Mechanics,1985,111(2)：218-234

15 Eshraghi H,Dravinski M.Transient scattering of elastic waves by dipping layers of arbitrary shape part 1：Antiplane strain model.Earthquake Engineering&Structural Dynamics,1989,18(3)：397-415

16 Eshraghi H,Dravinski M.Transient scattering of elastic waves by dipping layers of arbitrary shape.Part 2：Plane strain model.Earthquake Engineering&Structural Dynamics,1989,18(3)：417-434

17 Lee VW.Three-dimensional di ff raction of plane P,SV&SH waves by a hemispherical alluvial valley.International Journal of Soil Dynamics and Earthquake Engineering,1984,3(3)：133-144

18 Eshraghi H,Dravinski M.Scattering of plane harmonic SH,SV,P and Rayleigh waves by non-axisymmetric three-dimensional canyons：A wave function expansion approach.Earthquake Engineering&Structural Dynamics,1989,18(7)：983-998

19 Wong HL,Jennings PC.E ff ects of canyon topography on strong ground motion.Bulletin of the Seismological Society of America,1975,65(5)：1239-1257

20 Zhang L,Chopra AK.Three-dimensional analysis of spatially varying ground motions around a uniform canyon in a homogeneous half-space.Earthquake Engineering&Structural Dynamics,1991,20(10)：911-926

21 Wong HL.E ff ect of surface topography on the di ff raction of P,SV,and Rayleigh waves.Bulletin of the Seismological Society of America,1982,72(4)：1167-1183

22 S′anchez-Sesma FJ,Campillo M.Di ff raction of P,SV,and Rayleigh waves by topographic features：A boundary integral formulation.Bulletin of the Seismological Society of America,1991,81(6)：2234-2253

23 Pedersen HA,Sanchez-Sesma FJ,Campillo M.Three-dimensional scattering by two-dimensional topographies.Bulletin of the Seismological Society of America,1994,84(4)：1169-1183

24 Luco JE,Wong HL,De Barros FCP.Three-dimensional response of a cylindrical canyon in a layered half-space.Earthquake Engineering&Structural Dynamics,1990,19(6)：799-817

25 Luco JE,De Barros FCP.Seismic response of a cylindrical shell embedded in a layered viscoelastic half-space.I：Formulation.Earthquake Engineering&Structural Dynamics,1994,23(5)：553-567

26 Haskell NA.The dispersion of surface waves on multilayered media.Bulletin of the Seismological Society of America,1953,43(1)：86-103

27 Thomson WT.Transmission of elastic waves through a stratifie solid medium.Journal of Applied Physics,1950,21(2)：89-93

28 Kausel E.Fundamental Solutions in Elastodynamics：A Compendium.New York,USA：Cambridge University Press,2006

29 Bouchon M.A simple,complete numerical solution to the problem of di ff raction of SH waves by an irregular surface.The Journal of the Acoustical Society of America,1985,77(1)：1-5

30 Kawase H.Time-domain response of a semi-circular canyon for incident SV,P,and Rayleigh waves calculated by the discrete wavenumber boundary element method.Bulletin of the Seismological Society of America,1988,78(4)：1415-1437

31 Mossessian TK,Dravinski M.Application of a hybrid method for scattering of P,SV,and Rayleigh waves by near-surface irregularities.Bulletin of the Seismological Society of America,1987,77(5)：1784-1803

32 Shah AH,Wong KC,Datta SK.Di ff raction of plane SH waves in a Half-space.Earthquake Engineering&Structural Dynamics,1982,10(4)：519-528

33 Muskhelishvili NI.Some Basic Problems of the Mathematical Theory of Elasticity.4th Edn.Netherlands,Noorho ff,New York,University of Groningen,1963(in Russian)

34 Liu D,Gai B,Tao G.Applications of the method of complex functions to dynamic stress concentrations.Wave Motion,1982,4(3)：293-304

35 Liu Q,Zhao M,Wang L.Scattering of plane P,SV or Rayleigh waves by a shallow lined tunnel in an elastic half space.Soil Dynamics and Earthquake Engineering,2013,49：52-63

36王汝梁.地下管線的震害、抗震驗算、設計與措施.王珞珈等譯.北京：地震出版社，2007(Wang Ruliang.Buried Lifeline Earthquake Engineering and Related Research.Wang LJ,et al.Trans.Beijing：Seismicity Press,2007(in Chinese))

37侯忠良.地下管線抗震.北京：學術書刊出版社，1990(Hou Zhongliang.Earthquake Resistance of Buried Pipelines.Beijing：Academic Press,1990(in Chinese))

38 Earthquake Engineering Committee.The Japan Society of Engineers.Earthquake Resistant Design for Civil Engineering Structures in Japan.1988

39 Sakurai A,Takahashi T.Dynamic stresses of underground pipelines during Earthquake//Proceedings of 4th World Conference on Earthquake Engineering,Chile,1969,811-895

40 Newmark NM.Problems in wave propagation in soil and rock//Proceedings of International Symposium in Wave Propagation and Dynamic Properties of Earth Materials,Albuquerque,New Mexico,1967,7-26

41 Yeh GCK.Seismic analysis of slender buried beams.Bulletin of the Seismological Society of America,1974,64(15)：1551-1562

42 Lin G,Li ZY,Li JB.A substructure replacement technique for the numerical solution of wave scattering problem(to be published)

43 Song C.The scaled boundary finit element method in structural dynamics.International Journal for Numerical Methods in Engineering,2009,77(8)：1139-1171

44 Lin G,Han Z,Li J.Soil–structure interaction analysis on anisotropic stratifie medium.G′eotechnique,2014,64(7)：570-580

45 Wolf JP.Dynamic Soil-Structure Interaction.[s.l.]：Prentice Hall Int.,1985

46 Luco JE,de Barros FCP.Three-dimensional response of a layered cylindrical valley embedded in a layered half-space.Earthquake Engineering&Structural Dynamics,1995,24(1)：109-125

47 Dravinski M,Mossessian TK.Scattering of plane harmonic P,SV,and Rayleigh waves by dipping layers of arbitrary shape.Bulletin of the Seismological Society of America,1987,77(1)：212-235

48 de Barros FCP,Luco JE.Seismic response of a cylindrical shell embedded in a layered viscoelastic half-space.II：Validation and numerical results.Earthquake Engineering&Structural Dynamics,1994,23(5)：569-580

49 Lin G,Du J,Hu Z.Earthquake analysis of arch and gravity dams including the e ff ects of foundation inhomogeneity.Frontiers of Architecture&Civil Engineering in China,2007,1(1)：41-50

50 YinX,LiJ,WuC,etal.ANSYSimplementationofdampingsolvent stepwise extraction method for nonlinear seismic analysis of large 3-D structures.Soil Dynamics&Earthquake Engineering,2013,44(1)：139-152

51林皋，韓澤軍，李建波.層狀地基任意形狀剛性基礎動力響應求解.力學學報，2012，44(6)：1016-1027(Lin Gao,Han Zejun,Li Jianbo.Solution of the dynamic response of rigid foundation of arbitrary shape on multi-layered soil.Chinese Journal of Theoretical and Applied Mechanics,2012,44(6)：1016-1027(in Chinese))

52 Lin G,Han Z,Li J.General formulation and solution procedure for harmonic response of rigid foundation on isotropic as well as anisotropic multilayered half-space.Soil Dynamics&Earthquake Engineering,2015,70：48-59

53 Lin G,Han Z,Lu S,et al.Wave motion equation and the Green’s function for transverse isotropic multi-layered half-space.(to be published)

54 Zhong W,Williams FW,Bennett PN.Extension of the wittrickwilliams algorithm to mixed variable systems.Journal of Vibration&Acoustics,1997,119(3)：334-340

55 Zhong WX.On precise integration method.Journal of Computational&Applied Mathematics,2004,163(1)：59-78

56 de Barros FCP,Luco JE.Amplificatio of obliquely incident waves by a cylindrical valley embedded in a layered half-space.Soil Dynamics&Earthquake Engineering,1995,14(3)：163-175

57 Wolf JP,Song C.Finite-element Modelling of Unbounded Media.Chichester：Wiley,1996

58 Lin G,Lu S,Liu J.Transmitting boundary for transient analysis of wavepropagationin layered media formulated based on acceleration unit-impulse response.Soil Dynamics&Earthquake Engineering,2016,90：494-509

A COMPUTATIONAL MODEL FOR SEISMIC RESPONSE ANALYSIS OF UNDERGROUND STAUCTURES1)

Lin Gao2)
(State Key Laboratory of Coastal and O ff shore Engineering，Dalian University of Technology，Dalian 116024，China)(Institute of Earthquake Engineering，Dalian University of Technology，Dalian 116024，China)

Under earthquake excitation,the deformation of underground structures is restricted by the surrounding soil media.The dynamic behavior of it displays quite di ff erently from aboveground structures.Considerable progress has been made on the design and research of the underground structures since seventies-eighties of twentieth century.However,in the main,the widely used computational methods and guidelines in engineering design practice are based on rather simple assumptions,in reality,the actual soil conditions might be much more complex than ideally boundary conditions.Recent achievements of earthquake research on underground structures lie in the development of various computational methods for wave scattering problems of underground structures,such as the wave function expansion method,the boundary integral equation method etc.As the computation is somewhat complex,which impedes its application and dissemination in the engineering design practice.The author devotes himself to the improvement of the computational model for seismic analysis of underground structures,such that it achieves higher accuracy and efficiency,meanwhile it proves to be convenient for engineering design.To this end,a new model for seismic response analysis of underground structuresis proposed.The model is versatile to deal with wave scattering and di ff raction by canyons,subsurface cavities,subways and tunnels etc.In case of the presence of complex soil conditions like the layered half-space,a simple and e ff ective technique is developed for the evaluation of Green’s functions.Numerical examples are provided to validate the accuracy and efficiency of the proposed approach.

underground structure,seismic response,wave scattering,Green’s function,layered half-space.

TU93

：A

10.6052/0459-1879-16-301

2016–10–28 收稿，2017–01–25 錄用,2017–03–02 網絡版發(fā)表.

1)國家重點研發(fā)計劃資助項目(2016YFB0201001).

2)林皋，中國科學院院士，教授，主要研究方向：水利工程與核電工程抗震.E-mail：gaolin@dlut.edu.cn

林皋.地下結構地震響應的計算模型.力學學報,2017,49(3)：528-542

Lin Gao.A computational model for seismic response analysis of underground stauctures.Chinese Journal of Theoretical and Applied Mechanics,2017,49(3)：528-542