三維橢球形盆地對球面波的散射

馮 亭, 劉中憲,2, 黃 磊, 劉明珍

(1. 天津城建大學土木工程學院, 天津 300384; 2. 天津市軟土特性與工程環(huán)境重點實驗室, 天津 300384;3. 天津大學建筑工程學院, 天津 300372; 4. 青島上流遠大住宅工業(yè)有限公司, 山東 青島 266041)

0 引言

震害研究表明,局部場地對強地震動有明顯的放大效應[1-3],其本質是地震波在傳播過程中發(fā)生復雜的散射、衍射、相干作用以及波型之間轉換等。國內外眾多大中城市坐落在沉積盆地之上,或者許多城市內分布著古河道,由于沉積盆地卓越周期長,對于自振周期長的大跨橋梁、高層建筑等結構物在地震中造成更為嚴重的破壞,所以災害性地震中往往這些城市的經濟損失、人員傷亡嚴重。

該問題的定量分析方法有解析法[4-6]和數值法(有限單元法[7]、邊界單元法[8-9]和離散波數法[10])。關于沉積盆地對彈性波的散射問題,Trifunac[11]采用波函數展開法模擬分析了半空間中二維沉積盆地SH波的散射影響;劉殿魁等[12-13]利用復變函數法對平面SH波入射凹陷盆地的散射問題進行了解析求解;Lee[14]利用球波函數展開法對彈性波入射半球形沉積盆地的散射問題進行了模擬分析,但該模擬僅限于低頻波散射的解答;Bouchon[15],Kamalian[16],Zhao and Valliappan[17],模擬分析了二維沉積盆地對彈性波散射的影響,根據所得結果對影響局部場地(沉積盆地)放大效應的參數進行了詳細分析;Mossessian and Dravinsk[18]、Sanchez-Sesma and Luzon[19]利用間接邊界積分方程方法對彈性波入射三維沉積盆地半空間的散射情況進行了模擬分析,結果顯示三維散射情況與二維散射情況差異顯著,進行三維盆地場地效應分析,更接近實際情況。

值得指出的是,以上文獻對沉積盆地的地震響應研究主要針對平面波散射情況,沉積盆地對球面波散射問題目前研究較少,當波源距離較近時,入射波曲率對地震動的影響不可忽視,因此對球面波散射問題的研究顯得尤為重要。另外,目前對于三維沉積盆地對彈性波散射問題的研究多針對低頻波,然而高頻波對結構物的抗震性能影響不能忽略,且三維沉積盆地對高頻彈性波散射的研究還不夠系統(tǒng)深入,因此利用高精度數值模擬方法對三維沉積盆地對球面波散射問題在寬頻段情況下進行系統(tǒng)的參數分析和規(guī)律探討是有必要的。

本文針對三維沉積盆地對球面波散射的寬頻問題,建立了一種新的快速寬頻邊界元法。在GMRES迭代求解過程中,只有計算殘值向量和矩陣向量乘運算時才需要存儲系數大矩陣,通過修改GMRES迭代計算格式求解線性方程組,把存儲量降低到與計算自由度成正比[O(N)]。另外,應用OpenMP+MPI混合并行編譯模型,實現多核SMP集群的并行運算,提高了計算效率。基于該方法著重探討了入射波頻率、波源與盆地距離等參數對沉積盆地地震動放大效應的影響規(guī)律,研究結論對于沉積盆地區(qū)域城市的規(guī)劃建設、大中型工程抗震設計等具有一定的參考價值。

1 計算模型

均質半空間中存在一半橢球形沉積盆地,如圖1(a)所示,以盆地中心地表點為坐標原點o建立直角坐標系,三維橢球形沉積盆地的幾何表達式如下:

x=axcos(β)cos(φ)

(1)

y=aycos(β)sin(φ)

(2)

z=hcos(β)

(3)

式中:β表示交界面S0上任一點與z軸的夾角;φ表示交界面S0上任一點在xoy面的投影點與x軸正方向的夾角,0≤β≤π/2,0≤φ≤2π ;ax表示沉積盆地沿x方向的半寬;ay表示沿y方向的半寬;h表示沉積最大深度。

圖1 三維沉積盆地計算模型Fig.1 Three-dimensional sedimentary basin calculation model

基巖半空間域設為E,三維沉積盆地域設為R,半空間域E與沉積盆地域R交界面設為S0,半空間地表面設為S2,沉積盆地地表面設為S1,區(qū)域R和E的介質均簡化為均勻、各向同性和線彈性。假設入射一球面波,在沉積盆地交界面將發(fā)生散射(衍射),波源埋深dz取10 km,波源與盆地正下方O1點的垂直距離設為ds。

本文利用全空間格林函數進行間接邊界元方法計算,采用ANSYS軟件建立模型,利用四邊形單元離散網格,如圖1(b)所示。為了滿足計算精度要求,考慮散射波的衰減效應,離散范圍一般情況下取5倍的波長就能滿足,本文計算范圍取為4ax。為保證彈性波散射的求解精度,一個波長至少離散10個單元,本文經過計算對比不同離散精度的模型,結果已經計算穩(wěn)定。值得指出的是本文方法適用于模擬任意形狀的沉積盆地,而本文暫以一半橢球形沉積盆地為例進行模擬分析。

2 頻域響應計算

針對三維地形對地震波的散射,本文利用全空間格林函數進行間接邊界元法模擬求解,計算區(qū)域介質假設為均勻各向同性、線彈性。不受體力作用情況下,彈性波在彈性各向同性固體介質中的位移場穩(wěn)態(tài)波動方程為:

日前，天津濱海新區(qū)通報了12起不作為不擔當的典型案例，集中表現出的問題主要有：黨的建設弱化，堅持以黨的政治建設為統(tǒng)領力度不夠，履行全面從嚴治黨責任不力，管黨治黨“寬松軟”；形式主義、官僚主義問題突出，貫徹上級決策部署打折扣，執(zhí)行政策“中梗阻”；群眾意識淡漠，服務發(fā)展服務群眾不積極，不思進取“躲推繞”；落實扶貧助困政策不到位，裝聾作啞當“木官”；推進生態(tài)環(huán)境治理不主動，履行安全生產管理服務職能不到位；執(zhí)行紀律制度規(guī)定不嚴格，作風漂浮“庸懶散”，等等。

(λ+μ)·u+μ2u=-ω2u

(4)

式中:λ、μ為介質的拉梅常數;為矢量微分算子;u=(u,v,w)為位移矢量;ω為簡諧波運動頻率,時間因子exp(iωt)已略去。

基于單層位勢理論,三維空間域內一點的散射波場位移和牽引力可表達為某邊界連續(xù)面S上的積分:

(5)

(6)

式中:Gij(x,y)是位移格林函數,表示由作用在y處j方向的一單位力引起的x處i方向的位移;Tij(x,y)是應力格林函數,表示由作用在y處j方向的一單位力引起的x處i方向的應力;φj(y)可看作在邊界面S上y處j方向上作用的虛擬荷載密度。

2.1 波場分析

根據間接邊界元法(IBEM)基本原理可知,基巖半空間域總波場是由自由場和散射場疊加構成,沉積盆地域總波場只由散射場構成。自由波場是指球面P波入射不包含沉積盆地域時的基巖半空間域的波場解答。根據惠更斯原理,散射波場是指在沉積盆地交界面和所取計算區(qū)域的基巖半空間自由地表面施加虛擬均布荷載構成。

則,基巖半空間域位移:

(7)

(8)

(9)

(10)

2.2 邊界條件及求解

(1) 沉積盆地域地表面S1與基巖半空間域地表面S2的應力自由邊界條件為:

(11)

(12)

把式(6)代入式(11)、(12)可得:

(13)

(14)

式(14):基巖半空間地表面S2上的自由波場應力為空。

(2) 沉積盆地交界面上的位移應力連續(xù)邊界條件為:

(15)

(16)

把(9)、(10)式分別代入(15)、(16)式可得:

(17)

(18)

求解時需要離散沉積盆地交界面和整個自由地表面,并且在每個離散單元上施加虛擬均布荷載。覆蓋在離散單元上的斜面圓盤中心位置設為xm,虛擬均布荷載中心表示為ym,基巖半空間域的表面離散單元數設為N1,沉積盆地域地表面離散單元數設為N2,沉積盆地交界面離散單元數設為N0,假設虛擬荷載密度φj(ξ)為常量,則公式(13)、(14)、(17)、(18)的離散化線性方程如下所示:

(19)

(20)

(21)

(22)

以上離散化方程最終可整理成為(6N0+3N1+3N2,6N0+3N1+3N2)階的線性方程組,求解該線性方程組得到虛擬波源密度,進而可得到各離散點處的位移,再通過(7)、(8)式計算可得總的位移場。

文中已經過多次不同離散程序的設計對比,結果均已收斂,且以半球形均質沉積盆地模型為例,對程序編制的精確性進行了檢驗,精度檢驗詳見第3部分。另外考慮到(5)、(6)式積分時有奇異點存在,本方法中采用了一個與離散單元面積等效的斜面圓盤,并在其上施加均布荷載,并利用格林函數展開式求解,來解決積分奇異性問題。

3 精度檢驗

為檢驗程序編制的精確性,本文以半球形均質沉積盆地模型為例。將球面波設于沉積盆地正下方無限遠處(取dz=100 km),退化為平面P波垂直入射情況,將模擬結果與Mossenssian和Dravinski[18]給出的結果對比。計算參數取值:沉積盆地和基巖半空間的泊松比v2=v=1/3,沉積盆地與基巖半空間介質剪切波速比v2/v=0.5,密度比ρ2/ρ=2/3,材料黏滯阻尼比ξ=0.01,定義無量綱頻率為散射體等效直徑與入射波波長的比值,即:η=2a/λ=ωa/πcs=0.75,a為沉積半徑,λ為半空間入射波波長,cs為基巖半空間剪切波速。從圖2可以看出本文方法的計算結果與文獻給出的結果兩者吻合良好,由此驗證了本文方法的計算精度。

圖2 本文計算結果與Mossessian and Dravinsk所得地表位移幅值對比Fig.2 Comparison between the surface displacement amplitudeobtained by Mossessian and Dravinsk and thecalculation results in this paper

4 算例與分析

基于快速寬頻間接邊界元法,針對基巖半空間內三維半橢球形沉積盆地對球面波的散射進行研究。波源埋深取dz=10 km,考慮波源與盆地垂直距離的變化ds=0、2 km、5 km、10 km、50 km,入射波頻率取η=0.5、2.0、5.0。沉積盆地和基巖半空間的泊松比為v2=v=1/3,沉積盆地與半空間介質剪切波速比V2/V=0.5,密度比ρ2/ρ=2/3,材料黏滯阻尼比分別為ξ=0.02、ξ=0.01。為保證彈性波散射的求解精度,一個波長至少離散10個單元,具體離散單元數為:①當η=0.5、2.0時單元數為9 351。②當η=5.0時單元數為36 765。符號表示與精度驗證部分相同。

4.1 三維沉積盆地地表位移幅值單頻分析

圖3到圖5分別給出了膨脹波源在波源埋深與無量綱頻率相同的情況下,不同水平位置處地震波的放大效應和空間分布狀態(tài),并將結果與相應的平面波結果[9]進行對比分析。其中無量綱頻率η=0.5、2.0、5.0。在圖中位移幅值是用相應自由場(無沉積盆地存在)最大豎向位移幅值標準化。

首先考慮在同一水平位置(膨脹波源位于沉積盆地正下方ds=0)處,不同入射波頻率下,對地震波的放大效應和空間分布狀態(tài)的影響。從圖中可以看出,①較低頻率(η=0.5)膨脹波入射下,主方向(z向)位移呈現聚焦效應,位移最大值集中在盆地中心位置,且位移幅值從盆地中心到邊緣位置逐漸減小。②高頻波(η=2.0,5.0)入射時,主方向位移出現波峰波谷交錯現象,這主要是由于從盆地底部透射的體波與盆地內體波轉換的面波產生的干涉效應。③頻率η=0.5時,盆地內部主方向位移相對于基巖半空間的放大效應不明顯,位移幅值僅為2.15,與自由場位移接近。入射波頻率越大,沉積盆地相對半空間的放大效應越顯著,且干涉效應越強烈,例如,η=5.0時,盆地中心位置處z向位移峰值達到21.7。④沉積盆地次方向(x、y向)位移出現兩點聚焦效應,隨著頻率增大,聚焦區(qū)域減小。將以上結果與相應平面波結果[9]對比發(fā)現地表位移幅值一般小于平面波入射情況,例如η=2.0,平面波入射時地表位移幅值為9.8,而球面波入射時地表位移幅值為7.6。

考慮膨脹波源水平位置變化,①高頻波(η=2.0,5.0)入射時,當波源到盆地水平距離較小(ds小于5 km)時,盆地主方向位移仍出現波峰波谷交錯現象,η=5.0,ds=2 km時,z向位移幅值達到18.2。②靠近點源一側區(qū)域主方向位移放大效應高于遠離點源一側區(qū)域,如η=2.0,ds=5 km時盆地左側區(qū)域主方向位移放大效應明顯高于右側。③受膨脹波源曲率對波的散射影響,當波源水平距離較近時,主方向位移波峰波谷是圍繞盆地中心呈圓環(huán)狀交替出現,但當入射波頻率較高,波源水平距離較遠時,如:η=5.0,ds=50 km時,主方向位移波峰波谷呈縱向交替呈現。這種相鄰區(qū)域地震動放大效應的差異性已被多次地震觀測記錄所證實,對于較大規(guī)模沉積盆地地震區(qū)劃需要考慮其內部地震動的差異性。④整體來看,隨著波源與盆地水平距離增大,地表主方向最大位移幅值降低,例如:入射波頻率η=5.0,波源水平距離ds=0時,盆地中心位置處z向位移幅值達到21.7,波源水平距離ds=50 km時,z向位移僅為1.3。⑤點源位置對次方向位移幅值分布影響較大,隨著點源水平距離增大,次方向位移聚焦區(qū)域增加,且聚焦區(qū)域分布更分散。⑥高頻波入射時,次方向位移出現一定的邊緣效應,如η=5.0,ds=50 km時y向位移幅值為1.24,這是由于盆地底部透射的體波和盆地內部由體波轉換為的面波以及從盆地邊緣處透射的體波相互疊加而出現的現象。⑦高頻膨脹波入射時,點源水平距離在10 km以內時,主方向、次方向最大位移幅值均相差不大,但點源水平距離ds=50 km時,主、次方向最大位移幅值均降低很多,如入射波頻率η=5.0,波源水平距離ds=5 km時,z方向最大位移幅值為12.1,x方向最大位移幅值為5.46;波源水平距離ds=50 km時,z方向最大位移幅值為1.21,x方向最大位移幅值為0.87。

圖3 膨脹波源不同水平位置處半橢球形沉積盆地地表位移云圖(η=0.5)Fig.3 Surface displacement nephogram of semi-ellipsoid sedimentary basin at different horizontal location of the compressional wave source (η=0.5)

圖4 膨脹波源不同水平位置處半橢球形沉積盆地地表位移云圖(η=0.2)Fig.4 Surface displacement nephogram of semi-ellipsoid sedimentary basin at different horizontallocation of the compressional wave source (η=0.2)

圖5 膨脹波源不同水平位置處半橢球形沉積盆地地表位移云圖(η=5.0)Fig.5 Surface displacement nephogram of semi-ellipsoid sedimentary basin at different horizontallocation of the compressional wave source (η=5.0)

4.2 膨脹波源不同水平位置處沉積盆地地表位移幅值譜分析

三維沉積盆地對膨脹波的散射作用非常復雜,為了更全面反映波源水平位置對膨脹波散射的頻譜特性,圖6、7分別給出了膨脹波源入射下沉積盆地地表沿坐標軸上幾個典型點位的位移幅值譜。其中η∈(0,5.0),本文假設膨脹波源位置沿x軸負方向變化,則觀察點位取x軸上x/ax=-0.8,-0.5,0.0,0.5,0.8,y軸上y/ay=0.0,0.5,0.8。

從圖中容易看出,①當膨脹波源位置為ds=50 km時,沉積盆地地表不同點位位移幅值較小。②當膨脹波源與沉積盆地水平距離較近(ds≤10 km)時,在較低頻率(η=0.8)時,沉積盆地地表不同點位的主方向位移幅值差別不大且位移放大效應較小;在較高頻率域內,沉積盆地地表不同點位的位移頻譜特性差別很大,這主要是由于沉積盆地對膨脹波的散射作用,以及散射波之間復雜的相干效應。與平面波入射情況相比,球面波入射之下地表位移幅值譜曲線均有所降低。

圖6 膨脹波源不同水平位置處半橢球形沉積盆地地表x軸上位移幅值譜[Plane(y=0)]Fig.6 Displacement amplitude spectrum on x axis of the surface of semi-ellipsoid sedimentary basinat different horizontal position of the compressional wave source [Plane(y=0)]

圖7 膨脹波源不同水平位置處半橢球形沉積盆地地表y軸上位移幅值譜[Plane(x=0)]Fig.7 Displacement amplitude spectrum on y axis of the surface of semi-ellipsoid sedimentary basin at different horizontal position of the compressional wave source [Plane(x=0)]

考慮膨脹波源水平位置變化:①整體來看,沉積盆地地表主方向位移,隨著入射波頻率增大,位移頻譜曲線震蕩越劇烈,且隨著膨脹波源與沉積盆地水平距離越近,位移幅值頻譜值也基本越大。②不同點位處,第一階峰值頻率均在η=1.0附近,不同膨脹波源水平位置下,第一階峰值頻率基本一致,說明對于同一種入射波,不同波源位置下,位移“共振”頻率基本相同。③沉積盆地地表中心區(qū)域點位的位移幅值譜值一般大于其他點位。表明沉積盆地地表位移聚焦效應一般位于盆地中心區(qū)域。對該模型(沉積盆地與半空間介質剪切波速比V2/V=0.5),沉積盆地地表主方向位移譜峰值達到24.9。④膨脹波源與沉積盆地水平距離較近(ds≤5 km)時,沉積盆地中心區(qū)域點位位移幅值在高頻段內達到峰值,因此當膨脹波源與盆地水平距離較近時,在實際工程中應特別注意沉積盆地對高頻段波的放大效應。

5 結論

三維半橢球形沉積盆地對球面波散射的地震動特性依賴于入射波頻率、膨脹波源水平位置等因素。數值結果研究表明:

(1) 整體來看,隨著膨脹波源與盆地水平距離增大,地表主方向最大位移幅值降低,地表次方向位移聚焦區(qū)域增加且分布更分散。

(2) 膨脹波源水平位置為ds=0時,較低頻率膨脹波入射下,地表主方向位移呈現聚焦效應,且聚焦效應一般位于盆地中心區(qū)域,高頻波入射時,地表主方向位移出現波峰波谷交錯現象;次方向位移出現兩點聚焦效應。

(3) 位移頻譜分析表明,沉積盆地地表主方向位移,隨著入射波頻率增大,位移頻譜曲線震蕩越劇烈,且隨著水平距離越近,位移幅值頻譜值也基本越大。在較低頻率時,地表不同點位的主方向位移幅值基本相同;在較高頻率域內,地表不同點位的位移頻譜特性差別很大。地表不同點位處,第一階峰值頻率均在η=1.0附近,且不同膨脹波源水平位置下,第一階峰值頻率基本一致,即位移“共振”頻率基本相同。

(4) 本文結果與相應平面波結果對比發(fā)現地表位移幅值一般小于平面波入射情況,例如η=2.0,平面波入射時地表位移幅值為9.8,而球面波入射時地表位移幅值為7.6。