寶 鑫 劉晶波， 李述濤 王 菲

1） 中國北京 100084 清華大學土木工程系

2） 中國北京 100036 軍事科學院國防工程研究院

引言

隨著技術進步和經(jīng)濟發(fā)展以及“海洋強國”戰(zhàn)略的提出，我國對于海洋資源的開發(fā)和利用已步入快速發(fā)展階段，大量近海交通工程、海上風機和海洋平臺的規(guī)劃和建設已提上日程.我國地處環(huán)太平洋地震帶西側，海域地質(zhì)構造復雜，地震頻發(fā)，海洋工程的地震安全性問題不容忽視（李小軍，2006），海域場地的地震反應分析作為海洋工程抗震設計和地震安全性評價的重要基礎應予以高度關注.

早期對于海域場地地震動的研究主要采用解析方法，將海域地震動的產(chǎn)生和傳播簡化為數(shù)學模型，通過理論推導求解海域場地的應力場或波場分布（Lindsay，1939；Brekhovskikh，1980；鄭天愉等，1985；朱鏡清，1988；馮啟民，1990）在揭示流固耦合開放系統(tǒng)的波動傳播規(guī)律及其影響因素等方面取得了具有理論指導意義的研究成果.近年來，美國、日本等國家陸續(xù)在其臨近海域布設地震觀測系統(tǒng)，利用近海地震臺站記錄可分析海域地震動的時頻特征和工程特性，Chen 等（2017）、陳蘇等（2018）和楊銘等（2020）分別利用美國加州和日本相模灣等地的海域地震記錄，揭示了近海地震動具有長周期成分豐富，豎向分量強度低，反應譜峰值大等特點.同時，隨著計算機性能水平和數(shù)值模擬方法的成熟與發(fā)展，數(shù)值模擬方法憑借其在求解復雜模型和邊界條件等方面的優(yōu)勢，被更多地應用于海域場地的地震反應分析.采用譜元法（Linket al，2009）、有限差分法（Okamotoet al，2017）、傳遞函數(shù)模型（Liet al，2017）和有限元法（陳少林等，2019）建立了海水-海床流固耦合分析模型與方法，為復雜海域場地的地震動場計算提供了實用分析工具.胡進軍等（2017）通過數(shù)值模擬，比較了上覆海水層對海底地震動的放大效應在豎直方向和水平方向的差異，揭示了海底豎向地震動，尤其是P 波引起的豎向震動更易受海水層影響.榮棉水等（2013）分析了海底沉積層對海域地震動場和地震動參數(shù)的影響表明，當考慮海底沉積層時，海域場地的地震動場存在放大現(xiàn)象，且地震動反應譜周期也有所增大.

然而，目前對海域場地地震反應的理論分析與數(shù)值模擬研究主要針對平坦的海底場地，而近海工程場地存在向海洋傾斜的單體邊坡，且邊坡一側直接與海水耦合.此類橫向不均勻的特殊場地類型導致其自由場的求解存在天然的困難，對現(xiàn)有的基于自由波場的地震波動輸入方法提出了挑戰(zhàn)；同時，場地與無限海水層和半無限海床的動力相互作用也是開展近海場地地震反應分析時需重點考慮的問題，其涉及的關鍵因素包括大范圍流-固耦合系統(tǒng)的數(shù)值模擬、無限固體域和流體域波動輻射效應的數(shù)值模擬等.

因此，綜合上述因素，本文擬采用基于聲流體單元的流固耦合算法（寶鑫，劉晶波，2017）模擬場地-海水動力相互作用，利用流體介質(zhì)動力人工邊界（劉晶波等，2017）和一致黏彈性人工邊界單元（劉晶波等，2006）模擬無限流、固介質(zhì)的輻射阻尼，提出一種基于人工邊界子結構，利用混合波場實現(xiàn)近海場地地震動輸入的方法，在此基礎上，建立近海場地地震反應分析模型，并對其適用性及計算精度進行驗證.

1 近海場地地震反應分析模型

近海場地與下部基巖、外部海水域和海床介質(zhì)共同構成半無限的開放系統(tǒng)，在地震反應分析中，需綜合考慮局部地形影響和場地與海水的耦合作用，截取有限的近場計算域進行建模分析，并對流體和固體介質(zhì)截斷邊界進行處理，以模擬其波動輻射效應.示意模型如圖1 所示，下文將對該數(shù)值模型所涉及的關鍵技術環(huán)節(jié)進行簡要介紹.

圖1 近海場地地震反應分析模型Fig. 1 Seismic response analysis model of offshore site

1.1 流體介質(zhì)動力人工邊界

劉晶波等（2017）基于波動理論，將流體介質(zhì)中的單側波動方程轉(zhuǎn)化為截斷邊界處的等效力學系統(tǒng)，從而提出一種離散化的流體介質(zhì)動力人工邊界條件（圖2），可用于吸收流體域截斷邊界處的外行波動.該流體介質(zhì)人工邊界由阻尼器與集中質(zhì)量構成.二維模型中對應的物理參數(shù)質(zhì)量M和阻尼C分別為

圖2 流體介質(zhì)動力人工邊界示意圖Fig. 2 Schematic diagram of dynamic artificial boundaries of fluid medium

1.2 一致黏彈性人工邊界單元

式中，G，ρ，cS和cP分別為固體介質(zhì)的剪切模量、密度、橫波和縱波波速；h為人工邊界單元的厚度；R為波源至人工邊界點的距離；αT和αN為人工邊界參數(shù)，劉晶波等（2006）給出的推薦值分別為αT＝0.5 和αN＝1，α＝αN/αT.

1.3 近海場地-海水動力相互作用有限元模型

基于通用有限元計算平臺Ansys 建立近海場地-海水動力相互作用有限元模型.其中，固體介質(zhì)和流體介質(zhì)分別采用Solid45 和Fluid80 單元建模，并通過耦合流-固界面節(jié)點的法向自由度，模擬海水與場地的動力相互作用（寶鑫，劉晶波，2017）.由于Solid45 和Fluid80 均為三維單元，應用于二維模型時，可僅在平面內(nèi)建立一層單元，并約束全部節(jié)點的出平面自由度.分別在流、固介質(zhì)截斷邊界處添加前兩節(jié)中介紹的流體介質(zhì)動力人工邊界和固體介質(zhì)一致黏彈性人工邊界單元.此外，為保證波動問題的數(shù)值模擬精度，離散化網(wǎng)格的尺寸 Δx應滿足（杜修力，2009）

式中，λmin為離散網(wǎng)格模型中波動傳播的最短波長，cmin為介質(zhì)中的最小波速，fmax為波動問題數(shù)值模擬的截止頻率.

2 基于混合波場的近海場地地震動輸入方法

在場地地震反應分析中，近場計算域的截取和人工邊界條件的施加將直接影響地震動的輸入過程.以域縮減法（Bielaket al，2003；Yoshimuraet al，2003）和波動法（劉晶波，呂彥東，1998）為代表，目前常用的地震動輸入方法一般以自由場地震動為輸入地震動場，通過將其轉(zhuǎn)化為等效地震荷載，在不影響人工邊界對于外行波動吸收的前提下，將地震動輸入至計算模型內(nèi)部.在此基礎上，為進一步簡化計算流程，劉晶波等（2018）根據(jù)波動法原理和有限元理論，提出一種土-結構相互作用系統(tǒng)地震波動輸入的人工邊界子結構法.該方法不依賴于人工邊界條件的具體形式，將等效地震荷載的計算轉(zhuǎn)化為人工邊界子結構模型的動力分析，即保證了計算精度又有效地提高了分析效率.但上述地震動輸入方法為了便于通過Shake91 和EERV 等一維化場地分析軟件計算得到自由場分布，主要分析橫向均勻、縱向成層的規(guī)則場地形式，而對于以圖1 中近海場地為代表的不規(guī)則場地，由于模型不同位置處的地形特征和地層分布存在明顯的差異，地震波入射下整體模型的自由波場難以預先求解，使得傳統(tǒng)的基于自由場地地震動的地震波動輸入方法難以實現(xiàn).為解決這一問題，本文提出一種基于人工邊界子結構，利用混合波場實現(xiàn)近海場地中地震P 波和SV 波垂直輸入的方法，實施思路如下.

將近海場地模型中的輸入波場分成左側、右側和底部輸入波場三部分，如圖3所示.對于模型兩側的波動輸入，可根據(jù)截斷位置處的地層特性建立等效的縱向成層、橫向均勻的場地模型，并根據(jù)輸入地震動計算相應的自由場分布，進而利用人工邊界子結構法，將該自由場分布轉(zhuǎn)化為等效地震荷載；而對于近場模型的底部截斷邊界，則僅考慮入射波場，并利用底部的人工邊界子結構將其轉(zhuǎn)化為等效地震荷載.將計算得到的不同截斷邊界處的等效地震荷載施加于近海場地模型，可完成地震動的輸入過程.

圖3 基于混合波場的近海場地地震動輸入方法Fig. 3 Seismic wave input method of offshore site based on hybrid wave field

從波動傳播的角度分析該方法的實現(xiàn)過程：入射波場分別從側邊和底面邊界輸入至計算模型，與內(nèi)部局部地形相互作用產(chǎn)生反射波和散射波.對于側面邊界，由于輸入的自由波場同時包含了入射波和反射波成分，因此經(jīng)模型自由表面反射的外行波動被與其對應的等效地震荷載所抵消，而由局部地形產(chǎn)生的散射波則被人工邊界吸收；對于底面邊界，由于在等效地震荷載的計算中僅考慮了入射波場，因此反射波場和散射波場都將被人工邊界所吸收.此外，在兩種波場輸入位置的交界處，即模型兩側底部角點位置，可能存在輸入波場不協(xié)調(diào)的問題，由此引起的額外波動最終也將被人工邊界所吸收.

根據(jù)以上思路提出的基于混合波場的近海場地地震反應分析方法避免了預先求解不規(guī)則場地自由場分布的困難，通過構建由模型兩側的自由波場和模型底部的入射波場組成的混合波場，并利用人工邊界子結構法將其輸入至近場計算模型，實現(xiàn)此類不規(guī)則場地的地震反應分析.該方法的具體實施步驟如下：

1） 建立近海場地有限元模型，在此基礎上截取對應的人工邊界子結構模型，該模型僅由人工邊界及與其相鄰的一層內(nèi)部介質(zhì)單元構成，如圖4a 所示.

圖4 基于混合波場的近海場地地震波動輸入方法實現(xiàn)步驟（a） 利用子結構模型求解等效地震荷載；（b） 在近海場地模型中施加等效地震荷載Fig. 4 Implementation steps of seismic wave input method of offshore site based on hybrid wave field（a） Solving equivalent seismic loads through the substructure model；（b） Applying equivalent seismic loads on the offshore site model

2） 根據(jù)左右兩側邊界位置處的地層特性分別進行等效自由場分析，獲得自由波場，其中右側截斷邊界對應的自由場模型為上覆海水層的半空間場地模型，Bao 等（2020）給出了不同類型波動輸入下此類場地自由場分布的理論公式；固定人工邊界子結構模型的最外層節(jié)點，對所有與側邊人工邊界相連的子結構節(jié)點施加自由波場位移時程，而在其余的與底部人工邊界相連的子結構節(jié)點上輸入入射波場位移時程（圖4a）；對子結構模型進行動力分析，求得全部人工邊界節(jié)點上的反力，即為等效地震荷載.

3） 對近海場地有限元模型的人工邊界節(jié)點施加由上一步驟獲得的等效地震荷載（圖4b），進行動力時程計算，即可完成近海場地的地震反應分析.

3 算例與分析

建立近海場地有限元模型，對本文方法的有效性和計算精度進行驗證.假定近海場地、基巖和海床為均一介質(zhì)，材料參數(shù)列于表1.近場模型的整體尺寸為100 m×50 m，水深10 m，向海坡傾角為45° （圖5）.在截斷邊界處分別添加流、固介質(zhì)人工邊界.采用本文方法將圖6所示的持時為0.2 s 的脈沖波分別以SV 波和P 波的形式垂直輸入至計算模型.模型中有限元網(wǎng)格的最大尺寸為1 m，滿足式（3）給出的計算精度要求.

圖5 計算模型示意圖Fig. 5 Diagram of calculation model

圖6 脈沖波時程Fig. 6 Time history of pulse wave

表1 介質(zhì)材料參數(shù)Table 1 Material parameters of media

由于本文方法誤差的主要來源可能是底部人工邊界難以充分吸收下行反射波，導致該波動經(jīng)底部邊界再次反射后傳播至模型頂部的近海場地區(qū)域，對計算精度產(chǎn)生影響.可將模型底部邊界取至足夠遠，使反射波在計算時長內(nèi)不傳回頂部場地觀測點，從而獲得擴展網(wǎng)格解，作為本文方法的驗證.

3.1 SV 波垂直入射

首先分析SV 波入射下近海場地地震動場分布情況.以模型水平向總長度L＝100 m 對近海場地、向海坡以及海床位置的測點坐標進行歸一化處理，以該無量綱坐標為縱坐標，對比采用本文方法和通過擴展網(wǎng)格方法計算得到的位移波形，如圖7 所示.同時，為衡量近海場地的地震動放大效應，定義地震反應峰值比R為

圖7 SV 波入射下近海場地的水平（a）和垂直（b）位移波形Fig. 7 Displacement waveforms in horizontal （a） and vertical （b） directions on offshore site under incident SV wave

其中u為位移，下標x和z分別代表水平和豎直方向，下標0 代表輸入同一地震動下均勻半空間平坦場地的地震反應，其峰值為入射波峰值的二倍，x為測點水平坐標，t為時間.

分析SV 波入射下近海場地位移峰值放大系數(shù)的空間分布（圖8）可以看出，本文方法計算得到的近海場地峰值地震反應與擴展網(wǎng)格解吻合良好.在靠近向海坡的近海場地區(qū)域，水平地震峰值放大系數(shù)大于1；向海坡由上至下，峰值比逐漸減小，至向海坡底部角點處峰值比達到最小值；在遠離向海坡的區(qū)域，近海場地與海床處的峰值放大系數(shù)均接近于1，說明在此位置處的場地地震反應受向海坡局部地形的影響較小.進一步分析圖7 可知，由于采用混合波場進行波動輸入時，模型底面的人工邊界難以充分吸收下行的反射波，導致初至波峰過后的后續(xù)波動存在較小的計算誤差.但鑒于峰值地震反應是場地地震反應分析和結構抗震設計中最為重要的的地震動參數(shù)，本文方法對于此類不規(guī)則的近海場地地震反應分析具有較為良好的適用性.

圖8 SV 波入射下近海場地位移峰值放大系數(shù)的空間分布Fig. 8 Spatial distributions of peak displacement amplification coefficient on offshore site under incident SV wave

3.2 P 波垂直入射

進一步將圖6 所示的脈沖波以P 波的形式垂直輸入至計算模型中，比較近海場地、斜坡以及海床位置處的位移波形（圖9），同時計算P 波入射下近海場地位移峰值放大系數(shù)的空間分布（圖10）.

圖9 P 波入射下近海場地的水平（a）和垂直（b）位移波形Fig. 9 Displacement waveforms in horizontal （a） and vertical （b） directions on offshore site under incident P wave

圖10 P 波入射下近海場地位移峰值放大系數(shù)的空間分布Fig. 10 Spatial distributions of peak displacement amplification coefficient on offshore site under incident P wave

與SV 波入射時的計算結果（圖8）類似，P 波垂直入射下，采用本文方法與擴展網(wǎng)格方法計算得到的峰值地震反應也具有良好的一致性，且地震動場的分布規(guī)律與SV 波入射時的結果相似，均體現(xiàn)為靠近向海坡的近海場地區(qū)域地震動場被放大，而向海坡由上至下峰值地震反應逐漸減小.觀察圖9 可以看出，采用本文方法計算得到的初至波峰過后的后續(xù)波動（0.6—1.0 s）與擴展網(wǎng)格解相比也存在一定的誤差，這同樣是由于底部邊界對下行反射波吸收不充分導致的.整體而言，采用本文方法能夠較為合理準確的實現(xiàn)P 波入射下近海場地的地震動輸入和地震反應分析.

4 討論與結論

針對橫向不規(guī)則的近海場地地震反應問題，本文采用由場地近場截斷模型兩側的自由波場和模型底部的入射波場構成的混合波場作為輸入波場，改進基于人工邊界子結構的地震波從輸入方法以實現(xiàn)近海場地的地震動輸入.同時，綜合考慮場地-海水動力相互作用以及遠場流、固介質(zhì)的輻射阻尼，構建了一種近海場地地震反應分析模型與時域整體分析方法.

數(shù)值算例的結果表明，采用本文模型與方法能較為合理可靠地計算得到的近海場地的峰值地震反應.鑒于峰值地震反應是場地地震反應分析和結構抗震設計中最為研究人員所關心的地震動參數(shù)，本文方法對于此類不規(guī)則的近海場地地震反應分析具有較為良好的適用性.

采用本文提出的混合波場方法進行波動輸入時，由于模型底部的人工邊界難以充分吸收下行反射波，導致計算得到的場地地震動的后續(xù)波動存在一定的誤差.進一步合理優(yōu)化近海場地模型的底部輸入波場，是該領域的后續(xù)研究中有待解決的關鍵問題.