移動荷載與土體中孔洞相互作用的2.5D邊界元法分析

摘要： 為了研究具有復雜表面的飽和土體與移動荷載相互作用機理，根據飽和土Biot理論，采用Fourier變換和勢函數分解法，推導了飽和土體頻域波數內的Green函數；建立了移動荷載作用下，飽和土體中孔洞動力響應的頻域波數域內邊界元法方程（2.5D boundary element methods）；利用快速Fourier逆變換法，得到了時間空間域內飽和土體的動力響應.研究結果表明：在頻域內，利用移動荷載方向的一致性建立的頻域波數域內邊界積分方程，可將3D空間問題轉化為頻域波數域的2D平面問題， 3D邊界元簡化為2D邊界元，使得計算面轉換為計算線，減小了計算規模.

關鍵詞： 格林函數；移動荷載；飽和土；Fourier積分變換

中圖分類號： TU 471文獻標志碼： AAnalysis of Dynamic Interaction between Hole Embedded in

Saturated Soil and Moving Loads

Using 2.5D Boundary Element MethodXU Bin1，2，GAO Liang2，LEI Xiaoyan3，XU Manqing1，LIU Linya3

（1. Department of Civil Engineering， Nanchang Institute of Technology， Nanchang 330029， China； 2. School of Civil Engineering and Architecture， Beijing Jiaotong University， Beijing 100044， China； 3. Engineering Research Center of Railway Environmental Vibration and Noise， Ministry of Education， East China Jiaotong University， Nanchang 330013， China）

Abstract：To study the mechanism of interaction between the saturated soil with complex surface and moving loads， a specific 2.5D boundary element method （BEM） for the problem of dynamic responses of the hole embedded in the saturated soil under a moving load was derived systematically. Based on Biots theory， the frequency wave number domain Greens function for saturated porous media was developed using the Fourier transform method and the potential decomposition approach. Then， dynamic responses in the timespace domain solution were further obtained by synthesizing the wave number solution via inverse fast Fourier transform （IFFT）. The results show that using the frequency wavenumber boundary integrated equations established by the consistency of the directions of moving load， a 3D spatial model for the dynamic interaction between the holes embedded in the saturated soil with moving loads is convertible to a 2D planar model in the frequency domain， and the 3D complex boundary can be simplified to a problem of 2D boundary element. Thus， the complicated boundary planes can be transformed into computed lines using the 2.5D wave number domain BEM formulation， consequently reducing the amount of computation to a manageable size.

Key words：Greens function； moving loads； saturated soil； Fourier integral transform

工程中穿越軟土地基或道路下的市政供水、天然氣等管道不可避免受路面交通荷載作用，當埋地管道一旦達到破壞強度就可能產生安全事故，造成水、天然氣等泄漏而引發燃燒與爆炸[12].上述工程問題均可歸結于移動荷載與軟土地基中孔洞結構相互作用的動力響應模型.

對于具有簡單幾何構形的半空間土體表面受移動載荷問題，可利用解析方法進行計算[35]，但對于移動荷載作用下，地下孔洞結構的動力響應，邊界較復雜，其幾何軸對稱性已破壞，難以采用解析法計算.特別是移動荷載作用下，土體中各點動力響應都隨時間和空間不同而改變，采用有限元法對土體內部單元網格劃分時，由于時間和空間的變化，不可避免造成單元數量急劇增加，計算困難.邊界元方法（BEM）具有自動滿足無窮遠處輻射條件及降維功能，已被證明是解決具有復雜邊界條件問題一種有效手段，并且在工程中廣泛應用.文獻[6]中以3個土骨架位移和孔隙水壓力為基本變量，建立了飽和土體邊界積分方程.采用相同西南交通大學學報第48卷第4期徐斌等：移動荷載與土體中孔洞相互作用的2.5D邊界元法分析的基本變量，文獻[7]中推導了3D飽和土體頻域內邊界積分方程.文獻[8]中利用3D頻域邊界積分方程，分析了飽和多孔介質中散射問題.文獻[9]中利用拉普拉斯域的Green函數和卷積積分法，推導了時域內邊界積分方程方法.采用3D邊界元法，文獻[10]中對飽和土體中群樁土共同作用的動力問題進行了分析.

與傳統的有限元法相比，邊界元法只對邊界進行單元離散，因此，建模單元、計算時間均可節省.然而， 3D時間域內邊界元法仍然需占用計算機大量內存和計算時間.文獻[1112]在3D邊界元的基礎上，將結構作為2D分析或采用空間波數變換用來簡化問題，提出了波數域的2.5D邊界元法.與3D邊界元法相比， 2.5D邊界元法的邊界離散和計算要求大大減少.目前2.5D邊界元法已成功應用于彈性介質中的動力響應問題分析.文獻[13]根據具有恒定移動速度的點荷載Green函數，推導了任意角入射波下的無限長河谷3D響應的2.5D間接邊界積分方程.利用全空間Green函數，文獻[14]采用直接2.5D BEM，分析了粘彈性半空間中無限長隧道襯砌的3D地震荷載作用下的動力響應.采用離散波數的有限元與邊界元耦合，并結合傅里葉逆變換法.文獻[15]對軌道交通荷載引起地面3D空間的動力響應進行了評價.

上述文獻的研究表明：土體為單相彈性介質，忽略土體在移動載荷作用下土體中的水相、土顆粒之間的相互耦合作用，則無法考慮土體中孔壓的變化，以及孔壓變化對周圍土體動力響應的影響[3，16].多孔軟土中孔隙水壓的變化對土體強度有較大的影響，甚至引起震陷、液化等工程現象，這些都使得軟土地基具有非線性、非彈性的動力本質特征[1718].傳統的線彈性、粘彈塑性土體本構模型難以有效反映土體變形，特別是孔隙水壓力變化而導致的孔洞動態受力性狀.目前關于飽和土體中孔洞與移動荷載相互作用的2.5D BEM法尚未見報道.本文分別推導了由單位集中力作用在土骨架時所對應Green函數，以及單位散度源作用在孔隙介質流體時所對應Green函數.考慮荷載移動方向與飽和土體中孔洞軸線方向具有一致性，在三維邊界元法（3DBEM）的基礎上，采用沿飽和土體中孔洞縱向進行空間波數的積分變換，形成2.5D BEM，將三維移動荷載問題轉化為平面內邊界元積分問題，既可避免移動荷載引起的半無限空間有限元單元網格劃分隨時間而變化的難題，又可避免3D BEM在移動荷載方向的邊界劃分問題，能夠減少計算工作量.1飽和二相介質土Green函數為推導飽和地基土體模型的Green函數，定義時間t到頻域ω，空間z到波數域kz的Fourier變換對為

f-（ω）=∫+∞-∞f（t）e-iωtdt，

f（t）=12π∫+∞-∞f-（ω）eiωtdω，

f^（kz）=∫+∞-∞f（z）eikzzdz，

f（z）=12π∫+∞-∞f^（kz）e-ikzzdkz，（1）

式中：

“-”、“^”分別表示頻域、波數域量.

由于飽和多孔介質由土骨架及水相組成， Green函數張量G可分別由土骨架在單位集中力作用下所對應Green函數Gs、Gs以及飽和孔隙介質流體在單位散度源作用時所對應Green函數Gf、Gf 組成，表達式為

G=GsGf

GsGf.（2）

根據Biot理論[1718]，對于飽和二相介質，土骨架與水相分別在單位體積力Fi和fi（i=1，2，3）作用下，土體骨架位移ui（i=1，2，3）及孔隙水壓（pf）4個獨立變量表示的頻域內飽和土體運動方程為

μi，jj+（λ+μ）j，ji+ρgω2i-αgf，i=

-（i-β4f-i），

f，jj+β2ω2f-β3j，j=f-i，i，（3）

式中： m=α∞ρf/；

αg=α-β4；

ρg=ρb-β4ρf；

ρb=（1-）ρs+ρf；

β1=M/[mω2-iω（η/k）]；

β2=1/（β1ω2）；

β3=ρfω2-α[mω2-iω（η/k）]；

β4=ρfω2β1/M，

其中， λ和μ為Lame常數， α、M分別為飽和土體中土骨架、水相的壓縮系數， 為飽和多孔介質的孔隙率， wi（i=1，2，3）為飽和土體中流體相對于土骨架的位移， ρs、ρf分別為飽和土體土骨架、孔隙水密度， α∞為孔隙介質彎曲系數， k、η分別為飽和土體水相的滲透系數及孔隙流體的黏質度， ω為移動荷載頻率.

飽和多孔介質中土骨架在集中力作用下，對應Green函數的Gs、Gs，可由式（3）求得，

μΔ 2Gs+（λ+μ）Δ （Δ Gs）+

ρgω2Gs-αgGs=-δ（R）I3×3，

Δ 2Gs+β2ω2Gs-β3Δ Gs=0，（4）

式中：

Δ為Laplace算子.

由文獻[19]可知，可采用3個標量勢函數f，Gs、s，Gs、f，Gs表示Gs、Gs，則有

Gs=Δ Δ f，Gs+Δ Δ s，Gs+

（Δ 2f，GsI3×3-Δ Δ f，Gs），

Gs=AfΔ Δ 2f，Gs+AsΔ Δ 2s，Gs.（5）

利用恒等式

-δ（R）I3×3=Δ Δ 14πR+

Δ 2I3×3-Δ Δ 14πR，（6）

式中：

R為兩點的空間距離.

將式（5）代入式（4）中，同時利用式（6），可得

[afΔ 2f，Gs+ρgω2Δ 2f，Gs]+

[asΔ 2s，Gs+ρgω2Δ 2s，Gs]=14πR，

[AfΔ 2f，s+dff，Gs]+

[AsΔ 2s，Gs+dss，Gs]=0，

Δ 2Gs+k2tGs=14πμR，（7）

式中：

af=λ+2μ-αgAf；

as=λ+2μ-αgAs；

df=β2Afω2-β3；

ds=β2Asω2-β3；

k2f=ρgω2/（λ+2μ-αgAf）；

k2s=ρgω2/（λ+2μ-αgAs）；

k2t=ρgω2/μ；

Re（kf）≤Re（ks）.

利用式（1），對式（7）進行z→kz的Fourier變換，整理后可得頻域、波數域內的勢函數為

f，Gs=iAs4k2f（afAs-asAf）[H0，（2）（γfr）-H0，（2）（-ikzr）]，

s，Gs=iAf4k2s（asAf-afAs）[H0，（2）（γsr）-H0，（2）（-ikzr）]，

Gs=i4ρgω2[H0，（2）（γtr）-H0，（2）（-ikzr）]，（8）

式中： “～”表示頻域波數內的量；

H0，（2）為第二類的0階Hankel函數；

r為兩點的平面距離.

根據式（1）進行z→kz的Fourier變換，由式（5）可得頻域、波數域內的Green函數Gs、Gs，再將式（8）代入后，可得在集中力作用下飽和多孔介質中土骨架所對應的頻域波數域內Green函數為

ij，Gs=[δf，（s）H0，（2）（γfr）+δs，（s）H0，（2）（γsr）-δt，（s）H0，（2）（γtr）]，ij-

[δt，（s）H0，（2）（γtr）]δij，

i，Gs=χp[-H0，（2）（γfr）+H0，（2）（γsr）]，ii，j=1，2，3，（9）

式中：

δf，（s）=iAs/[4k2f（afAs-asAf）]；

δs，（s）=iAf/[4k2s（asAf-afAs）]；

δt，（s）=i/（4ρgω2）；

χp=-Ask2fδs，（s）.

利用式（3），可得單位散度源作用在飽和土體中流體介質時所對應Green函數Gf、Gf為

μΔ 2Gf+（λ+μ）Δ （Δ Gf）+

ρgω2Gf-αgGf=0，

Δ 2Gf+β2ω2Gf-β3Δ Gf=-δ（R）.（10）

由于流體介質中不存在剪切波，因此可用2個標量函數f，Gf、s，Gf表示相應的Green函數Gf、Gf，

Gf=Δ f，Gf+Δ s，Gf，

Gf=AfΔ 2f，Gf+AsΔ 2s，Gf.（11）

同單位集中力作用在土骨架時所對應的頻域波數域內Green函數類似的步驟，可得單位散度源作用在飽和土體流體介質時，所對應的頻域波數域內Green函數為

i，Gf=δ（f）[H0，（2）（γfr）-H0，（2）（γsr）]，i，

Gf=δ（f）[-Afk2fH0，（2）（γfr）+

Ask2sH0，（2）（γsr）]，i=1，2，3，（12）

式中：

δf=ias/[4k2f（asAf-afAs）].2飽和土體中2.5D邊界積分方程2.1飽和土體頻域內的邊界積分方程根據互易定理及飽和土體的本構關系（式（3）），考慮飽和土體任意2個狀態（狀態1和狀態2），可得

ij，（1）ij，（2）+f，（1）ξ（2）=

ij，（2）ij，（1）+f，（2）ξ（1），（13）

式中：下標（1）、（2）分別表示飽和土體的2種運動狀態；

σij為飽和土體總應力分量；

ξ為飽和多孔介質中流體的單位體積增量，

ξ=wi，i.

將式（13）沿飽和土體任意區域Ω積分，并結合式（3），可得

∫S[j，（1）j，（2）-j，（2）j，（1）]dS-

∫S[f，（1）n，（2）-f，（2）n，（1）]dS=

∫Ω[f，（1）f，（2）-f，（2）f，（1）]dΩ-

∫Ω[gj，（1）j，（2）-gj，（2）j，（1）]dΩ，（14）

式中：

j為作用在區域Ω邊界S上的面力，

j=ijni；

ni為邊界S外法線方向余弦；

f=β1i，i/M；

n為沿邊界S的外法線方向流體位移，

n=β1（f，j-ρfω2j）nj/M.

將飽和土體的真實解作為狀態1，單位集中力或單位散度源所應于Green函數解為狀態2，則在不考慮體積力時，由式（14）可得在飽和土體土骨架任意點x上作用一集中力時，飽和土體的Somigliana等式為

i（x，ω）=∫S[ij，Gs（x，ζ，ω）j（ζ，ω）-

ij，Gs（x，ζ，ω）j（ζ，ω）]dS（ζ）-

∫S[ni，Gs（x，ζ，ω）f（ζ，ω）-

i，Gs（x，ζ，ω）n（ζ，ω）]dS（ζ），（15）

式中：

ni，Gs=β1[i，Gs（x，ζ，ω）/ζj-

ρfω2ij，Gs]nj（ζ）/M，

ij，Gs=ijk，Gs（x，ζ，ω）nk（ζ），

其中， Gs、Gs分別為單位集中力作用在土骨架時所對應Green函數.

利用相同方法，飽和土體中流體介質受單位散度源作用時的Somigliana等式為

i（x，ω）=∫S[j，Gf′（x，ζ，ω）j（ζ，ω）-

j，Gf′（x，ζ，ω）j（ζ，ω）]dS（ζ）-

∫S[n，Gf′（x，ζ，ω）f（ζ，ω）-

Gf′（x，ζ，ω）n（ζ，ω）]dS（ζ），（16）

式中：

j，Gf′=β1j，Gf/M；

j，Gf′=β1j，Gf/M；

n，Gf′=β1n，Gf/M；

Gf′=β1Gf/M.

根據式（2）及式（15）、（16）可得到飽和土體統一的Somigliana等式為

i（x，ω）=∫S[ij，G（x，ζ，ω）j（ζ，ω）-

ij，G（x，ζ，ω）j（ζ，ω）]dS（ζ）.（17）

利用文獻[20]的方法，考慮不同邊界的特性，可得頻域內飽和土體的邊界積分方程為

ciji（x，ω）=∫S[ij，Gs（x，ζ，ω）j（ζ，ω）-

ij，Gs（x，ζ，ω）j（ζ，ω）]dS（ζ），（18）

在光滑邊界時， cij=δij/2.2.2 頻域波數域的2.5D邊界積分方程移動荷載與飽和土體中的孔洞相互作用模型如圖1所示.

圖1移動點荷載與飽和土體相互作用示意圖

Fig.1Schematic diagram of the interaction between

the moving point loads and the saturated soil若荷載沿孔洞軸向以速度v移動，飽和土體及其中孔洞的動力響應均可表示為

ui（x，t）=ui（x⊥，z-vt），

ti（x，t）=ui（x⊥，z-vt）.（19）

對式（19）進行時間、空間Fourier變換，可得

i（x⊥，kz，ω）=

2πδ（ω+ω0-kzv）U^i（x⊥，kz），

i（x⊥，kz，ω）=

2πδ（ω+ω0-kzv）T^i（x⊥，kz），（20）

式中：

ω0為移動荷載初始頻率；

U^i（x⊥，kz）=∫∞-∞ui（x⊥，z-vt）eikz（z-vt）d（z-vt），

T^i（x⊥，kz）=∫∞-∞ti（x⊥，z-vt）eikz（z-vt）d（z-vt），

其中， x⊥為平面內坐標點.

將式（20）代入式（18）中，可得到頻率波數域的飽和土體中孔洞的2.5D邊界元積分方程為

cijU^j（x⊥，kz）=

∫Γij，G（ζ⊥-x⊥，-kz，kzv）T^j（x⊥，kz）dΓ（ζ⊥）-

∫ΓT～ij，G（ζ⊥-x⊥，-kz，kzv）U^j（x⊥，kz）dΓ（ζ⊥），（21）

式中：

Γ 為飽和土體及孔洞在 平面內的邊界；

ζ⊥為 平面內任意觀察點坐標.

對于飽和土體-孔洞體系內任一點，有

U^j（x⊥，kz）=

∫Γij，G（ζ⊥-x⊥，-kz，kzv）T^j（x⊥，kz）dΓ（ζ⊥）-

∫ΓT～ij，G（ζ⊥-x⊥，-kz，kzv）U^j（x⊥，kz）dΓ（ζ⊥）.（22）

當ζ⊥→x⊥時，式（12）中的Green函數T～ij，G（ζ⊥-x⊥，-kz，kzv）存在奇異性，采用文獻[21]的方法，借助輔助問題組合解，消除奇異性.3邊界元離散與算例3.1 邊界元離散移動荷載與孔洞相互作用的頻域波數域的2.5D邊界積分方程式（22），沿邊界Γ的積分可采用文獻[22]中常規的2D邊界元法計算.若采用等參單元離散，在局部坐標系-1≤η≤1下，對于節點x有

x（η）=∑Mi=1Nl（η）x（l），

i（η）=∑Mi=1Nl（η）i，（l），

i（η）=∑Mi=1Nl（η）i，（l），（23）

式中：

M為每個等參單元的節點數.

將飽和土體內任一點xm作為場點，若相連單元數為er （ er=1或2）個，且任一單元的節點數為Ne，則將式（23）代入式（22）中，可得離散的積分方程為

1-c∞U^i（xm）=∑Nek=1∑Ml=1T^j，（kl）∫Γk（ij，G-Uij，Ga）Nl（η）J（η）dη+

∑Nek=1∑Ml=1T^j，（kl）∫ΓkUij，GaNl（η）J（η）dη-

∑Nek=1（k≠e1，…，er）∑Ml=1[U^j，（kl）∫ΓkT～ij，GNl（η）J（η）dη]-

∑Nek=e1，…，er∑Ml=1U^j，（kl）∫ΓkT～ij，GNl（η）-δnklTij，Ga]J（η）dη×

U^j（xm）∑Nek=1（k≠e1，…，er）∫ΓkTij，GaJ（η）dη，（24）

式中： c∞=1/2.求解離散后的邊界積分方程可得到頻域解，采用快速Fourier逆變換法（IFFT），對波數到空間域進行Fourier逆變換，可得飽和土體中孔洞動力響應的時間空間域響應.

采用IFFT進行Fourier逆變換時，樣本點數N，空間域樣本間距Δz及波數域樣本間距Δkz應滿足

Δkz=2πNΔz.（25）

當Δz≤λmin/2時（λmin為最小波長），才能保證足夠大的截止頻率.3.2 與已知文獻結果對比為驗證本方計算方法正確性，分析不同荷載速度下，半空間飽和土體表面受移動點荷載作用下的動力響應.以圖1為例，圖中移動點荷載幅值為Fz，以不同速度v沿z軸正方向在飽和土體表面上移動.

飽和土體參數為：

μ=3.0 GPa， a∞=2.0，

λ=3.0 GPa， a=0.95，

bp=1.0×1010 （kg/m3）·s，

M=5.0 GPa， =0.3，

ρs=2 500 kg/m3，

ρf=1 000 kg/m3，

v=0.1vSH，0.5vSH，0.9vSH（vSH=μ/ρb）.

采用文中方法，邊界離散僅在與荷載運動方向垂直的x方向進行，本文采用2節點等參單元進行x方向的邊界離散，單元總數200個， x方向總計算長度100.0 m，觀察點（x=1.0 m， y=-1.0 m）的豎向位移在-2.0 m≤z*=z-vt≤2.0 m范圍內變化，如圖2，圖中： u*y=μuyaR/Fz ， aR=1.0 m.

由圖2可知，采用文中的2.5D BEM法與文獻[3]結果一致.

圖2本文解與文獻[3]解析解結果比較

Fig.2Comparison of the present results with

the analytical solutions of Ref. [3]

4結論（1） 不同于粗糙的近似線彈性模型的是，采用飽和半無限空間的軟土地基模型可以分析交通荷載作用下，土體中的水相、土顆粒之間的相互耦合作用以及考慮土體中孔壓的變化對周圍土體動力響應的影響.

（2） 對于具有與荷載運動方向具有一致性的飽和土體中孔洞結構的動力響應問題，采用2.5D BEM，可將移動荷載的3D問題轉化為平面內2D邊界問題，使邊界的離散只在與荷載移動方向相垂直的其它2個方向進行，避免了有限元法單元網格隨荷載移動而破壞的難題，減少了計算工作量.

（3） 利用文中2.5D BEM，結合有限元法，還可分析飽和土體中，隧道襯砌在地鐵交通荷載作用下的動力響應及地震波場散射問題.參考文獻：[1]王永強，牛星鋼，譚欽文. 重型車輛荷載下埋地天然氣管道的安全分析[J]. 中國安全生產科學技術，2011，7（8）： 109114.

WANG Yongqiang， NIU Xinggang， TAN Qinwen. Safety analysis for buried gas pipelines under heavy vehicle loads[J]. Journal of Safety Science and Technology， 2011， 7（8）： 109114.

[2]王文和，易俊，沈士明. 基于風險的城市埋地燃氣管道安全評價模型及應用[J]. 中國安全生產科學技術，2010，6（3）： 163166.

WANG Wenhe， YI Jun， SHEN Shiming. Safety evaluation model and application for the urban buried gas pipeline based on risk[J]. Journal of Safety Science and Technology， 2010， 6（3）： 163166.

[3]LU J F， JENG Dongsheng. A halfspace saturated poroelastic medium subjected to a moving point load[J]. International Journal of Solids and Structures， 2007， 44（2）： 573586.

[4]EASON G. The stress produced in a semiinfinite solid by a moving surface force[J]. International Journal of Engineering Science， 1965， 2（6）： 581609.

[5]PAYTON R G. An application of the dynamical Betti Rayleigh reciprocal theorem of moving point loads in elastic medium[J]. The Quarterly of Applied Mathematics， 1964， 20（1）： 299313.

[6]DOMINGUEZ J. Boundary element approach for dynamic poroelastic problems[J]. International Journal for Numerical Methods in Engineering， 1992， 35（2）： 307324.

[7]CHENG A H D， BADMUS D E， BESKOS D E. Integral equations for dynamic poroelasticity in frequency domain with BEM solution[J]. Journal of Engineering Mechanics， ASCE， 1991， 117（5）： 11361157.

[8]ZIMMERMAN C， STERN M. Boundary element solution of 3D wave scatter problems in a poroelastic medium[J]. Engineering Analysis with Boundary Elements， 1993， 12（4）： 223240.

[9]SCHANZ M. Application of 3D time domain boundary element formulation to wave propagation in poroelastic solids[J]. Engineering Analysis with Boundary Elements， 2001， 25（4/5）： 363376.

[10]MAESO O， AZNAREZ J J， GARCIA F. Dynamic impedances of piles and groups of piles in saturated soils[J]. Computers Structures， 2005， 83（10/11）： 769782.

[11]PAPAGEORGIOU A S， PEI D. A discrete wavenumber boundary element method for study of the 3D response of 2D scatterers[J]. Earthquake Engineering Structural Dynamics， 1998， 27（6）： 619638.

[12]TADEU A J B， KAUSEL E. Greens functions for twoandahalfdimensional elastodynamic problems[J]. Journal of Engineering Mechanics， ASCE， 2000， 126（10）： 10931097.

[13]LUCO J E， WONG H L， BARROS F C P. Threedimensional response of a cylindrical canyon in a layered halfspace[J]. Earthquake Engineering Structural Dynamics， 1990， 19（6）： 799817.

[14]STAMOS A A， BESKOS D E. 3D seismic response analysis of long lined tunnels in halfspace[J]. Soil Dynamics and Earthquake Engineering， 1996， 15（2）： 111118.

[15]SHENG X， JONES C J C， THOMPSON D J. Prediction of ground vibration from trains using the wavenumber finite and boundary element methods[J]. Journal of Sound and Vibration， 2006， 293（3/5）： 575586.

[16]徐志英，徐新寧. 土的非線性對土地下結構相互作用的影響[J]. 河海大學學報，1993，21（1）： 1115.

XU Zhiying， XU Xinning. The influential factors on the soilunderground structure interaction due to the nonlinear mechanical behavior of soils[J]. Journal of Hohai University： Natural Sciences， 1993， 21（1）： 1115.

[17]BIOT M A. Theory of propagation of elastic waves in a fluidsaturated porous solid， low frequency range[J]. Journal Acoustic Social America， 1956， 28： 168178.

[18]BIOT M A. Generalized theory of acoustic propagation in porous dissipative media[J]. Journal Acoustic Social America， 1962， 34（9A）： 12541264.

[19]ANDERSEN L， NIELSEN S R K. Boundary element analysis of the steadystate response of an elastic halfspace to a moving force on its surface[J]. Engineering Analysis with Boundary Elements， 2003， 27（1）： 2338.

[20]RASMUSSEN K M， NIELSEN S R K， KIRKEGAARD P H. Boundary element method in time domain for a moving timedependent force[J]. Computer Structure， 2001， 79（7）： 69701.

[21]NORRIS A N. Radiation from a point source and scattering theory in a fluidsaturated porous medium[J]. Journal Acoustic Social America， 1985， 77（2）： 20122023.

[22]姜弘道. 彈性力學問題的邊界元法[M]. 北京：中國水利水電出版社，2008： 6070.