含空洞地層隧道動力響應解析解

臧萬軍，賈晨浩，王崢崢

(1.福建工程學院 土木工程學院，福建 福州 350118；2.地下工程福建省高校重點實驗室，福建 福州 350118；3.大連理工大學 土木工程學院，遼寧 大連 116024)

近年來，隧道地震災害頻發(fā)，地下結構抗減震問題已成為工程界熱點[1]，特別是鄰近地層含有空洞時，在地震波作用下空洞與隧道產生動力相互作用，將嚴重影響隧道襯砌結構安全。目前，對該問題的研究方法主要有模型試驗、解析法和數(shù)值模擬[2]。數(shù)值模擬的顯著特點是可適用于任意形狀隧道，應用廣泛，而解析法有助于探究問題的本質，克服軟件計算中由于網格劃分尺寸等問題造成計算結果不精確、計算慢等問題，具有較高的理論價值和工程實踐意義。由于入射波作用于地下圓形洞室時面臨著復雜的散射問題，目前為止解析解的數(shù)量非常有限。LEE等[3]首先給出了SV波入射下地下圓形非襯砌洞室的散射解析解；梁建文等[4?6]采用波函數(shù)展開法給出了半空間中洞室群在平面P波和SV波散射問題的級數(shù)解；王帥帥[7]推導出含減震層隧道襯砌在平面SV波入射下動應力集中系數(shù)級數(shù)解，并與數(shù)值模擬結果對比；臧萬軍等[8?9]采用ABAQUS有限元軟件，改變地震波傳入方向、溶洞大小、埋深等因素，得到溶洞隧道動力響應規(guī)律；PELLI等[10]給出了SV波作用下線彈性半平面空間中不同形狀洞室周圍的應力和變形場的解析解；楊公標等[11]給出了含空洞地層淺埋隧道開挖引起的位移解析解；RABETI等[12]通過改變輸入運動頻率、隧道埋深等參數(shù)，分析了地鐵盾構隧道在地震波激勵下對地表位移放大作用的規(guī)律；尹崇林等[13]通過冪級數(shù)解法得到了半圓拱形隧洞各邊界的應力和位移分布；朱翔宇等[14]采用FLAC3D軟件分析了溶洞和隧道之間的相互作用，得到溶洞對分離式隧道圍巖應力與變形影響的變化規(guī)律。以上學者研究了洞室或洞室群的動力響應位移和應力集中級數(shù)解，以及多種因素對隧道動力響應特性的影響規(guī)律，但未考慮地層空洞對隧道地震動力響應的影響，有待開展進一步研究。本文首先建立含空洞地層隧道理論計算模型，得到隧道動力響應位移解析解，并采用MIDAS/GTS有限元軟件建立三維模型，驗證解析解的合理性和準確性，隨后采用正交試驗方法開展空洞對隧道地震動力響應影響的敏感性因素分析，揭示隧道埋深、空洞與隧道中心距、空洞半徑和地震波入射角變化對含空洞地層隧道動力響應的影響規(guī)律。

1 理論解析解

1.1 計算模型建立

以一條半徑足夠大的圓弧替代水平自由地表[3]，圓形空洞和圓形隧道被鑲嵌在彈性、均勻和各向同性半空間中，兩者水平距離為d1；空洞半徑為a1，埋深為h1；隧道襯砌的內半徑為a2，外半徑為a3，埋深為h2。半無限空間屬性由拉梅常數(shù)λs，μs和質量密度ρs確定，縱波和剪切波波速分別為αs和βs。襯砌的材料屬性由拉梅常數(shù)λ2，μ2和質量密度ρ2確定，縱波和剪切波波速分別為α2和β2，理論簡化模型如圖1所示。

圖1 空洞和隧道理論模型Fig.1 Theoretical model of cavities and tunnels

1.2 波函數(shù)分析

一頻率為ω的平面SV波以入射角度θβ入射在O1極坐標系下的入射、反射S V波和反射P波勢函數(shù)可分別展開為Fourier-Bessel級數(shù)形式：

式中：ksβ=ω/βs為橫波波數(shù)；kαs=ω/αs為縱波波數(shù)；Jn(x)為第1類Bessel函數(shù)。

式中：θα為P波反射角，依據(jù)Snell定律，ksβsinθβ=ksαsinθα；當ni=0時，εni=1，ni≥1時，εni=2，k11和k12為反射系數(shù)，即

同理，O2極坐標系下的入射、反射SV波和反射P波勢函數(shù)的Fourier-Bessel級數(shù)形式為：

對于S V波入射，存在臨界角θcr=sin-1(βs/αs)，當θβ>θcr時，反射P波以不均勻的表面波形式存在，需將式(2)和式(8)中的系數(shù)變?yōu)閇15]：

由大圓弧表面產生的散射S V波和P波可分別表示為：

式中：k2β和k2α分別為襯砌材料中橫波波數(shù)和縱波波數(shù)；Hn(x)為第1類Hankel函數(shù)。其中上標1，2和3分別表示為坐標系O1，O2和O3；下標i，r和s分別表示為入射波，反射波和散射波；式中A，B，C和D為待定系數(shù)。

為了引入邊界條件進行求解，通過Graf加法公式進行坐標變換[16]，將上述波函數(shù)置于同一坐標系中表示。

1.3 邊界條件

引入邊界條件：

1)半空間表面零應力

3)隧道襯砌與半空間交界面處應力與位移連續(xù)

平面SV波作用下，彈性介質中平面應變問題的位移和應力勢函數(shù)表達式分別為[17]：

由模型邊界條件，可聯(lián)立求解線性方程組，得到全部待定系數(shù)。將半無限空間和襯砌中的位移勢函數(shù)代入表達式(31)中即可得到空洞位移和隧道襯砌位移的級數(shù)解析解。

1.4 位移解析解求解

空洞位移為：

隧道襯砌位移為

在直角坐標系(x,y)下的位移值可通過下式變換得到：

采用MATLAB軟件編譯含空洞地層隧道動力響應位移解析公式，其求解步驟主要有：1)在對輸入的地震波加速度時程數(shù)據(jù)進行傅里葉變換后，將得到的頻譜數(shù)據(jù)代入位移解析解中。2)經程序運算后，對每個頻率對應的振幅幅值乘以頻譜放大系數(shù)[18]。3)對計算結果采用傅里葉逆變換，得到空洞和隧道的動力響應位移時程曲線。

2 有限元數(shù)值模擬

2.1 數(shù)值模型建立

依托某山嶺隧道工程，采用MIDAS/GTS有限元軟件開展三維數(shù)值模擬分析，巖土體和隧道襯砌部分均采用實體單元，并選用彈性本構關系，模型邊界采用軟件提供的曲面彈簧功能設置為自由場，以此吸收和過濾地震波的能量，解決模型與邊界相連處波的反射問題，實現(xiàn)地震波在巖體中的無限傳播。含空洞地層隧道數(shù)值計算模型如圖2所示：圖2(a)為隧道襯砌斷面，曲墻式三心圓拱，跨度為19.49 m，洞口高度13.06 m，初期支護厚度26 cm；圖2(b)為三維有限元模型，二次襯砌厚度55 cm；三維模型橫截面寬140 m，隧道長100 m；圖2(c)和2(d)為空洞位置剖面圖，空洞位于巖體縱軸方向50 m處，并簡化為空心球體模型，且空洞中心與隧道中心在同一水平面上，其他參數(shù)根據(jù)不同工況設定。

圖2 含空洞地層隧道三維數(shù)值模型Fig.2 Three-dimensional numerical model of tunneling with cavitation

2.2 隧道主要參數(shù)

表1表示各介質物理力學參數(shù)，根據(jù)工程地質勘察報告確定。

表1 各介質物理力學參數(shù)Table 1 Physical and mechanical parameters of each medium

2.3 選擇輸入地震波

選取具有代表性的阪神地震波作為輸入地震波，該地震波的總時長為46 s，波峰發(fā)生在4.52 s，加速度最大幅值為0.83g，最大頻率為18 Hz，卓越頻率大約在1.5 Hz左右。采用SeismoSignal地震波處理軟件對阪神地震波原始加速度時程曲線進行濾波分析和基線校正，經分析可知在0～10 Hz之間的傅里葉幅值波動變化較為全面，此頻率范圍內的信號可以代表該地震波的主要能量，對應總時長中2～12 s的部分。為了減少工作量，采用該時程內的地震波數(shù)據(jù)為動力荷載輸入?yún)?shù)，截取的地震波加速度時程曲線如圖3所示。

圖3 阪神地震波加速度時程曲線Fig.3 Time-history curve of Kobe seismic wave acceleration

3 解析解的驗證

將解析方法引入上述工程實例，計算結果與數(shù)值解進行對比分析。隧道埋深和空洞與隧道中心距設為20 m，地震波入射角度分別為γ＝0°，γ＝30°和γ＝60°，相同條件下數(shù)值模擬解與解析解的對比結果如圖4和圖5所示。

圖4 空洞位移時程曲線Fig.4 Time-displacement curves of cavity

由圖4和圖5可知，數(shù)值模擬解和解析解位移時程曲線的變化趨勢相似度較高，其中前1 s和后2 s的位移幅值變化存在一定差異，但對整體結果影響較小。對比分析空洞和隧道在入射角度γ＝0°，γ＝30°和γ＝60°時數(shù)值模擬解和解析解的位移幅值峰值，空洞的位移幅值峰值相差9.7%，隧道的位移幅值峰值相差9.2%。數(shù)值解與解析解吻合較好，驗證了解析解的合理性。采用該解析解求解均質、各向同性的含空洞地層隧道動力響應位移時，可有效避免有限元劃分網格尺寸問題造成的數(shù)值模擬計算慢、計算結果有誤差等問題，為地下結構抗震分析提供了一種簡便的計算方法。

圖5 隧道位移時程曲線Fig.5 Time-displacement curves of tunnel

4 隧道動力響應敏感性分析

4.1 正交試驗方案

在解析解的基礎上，采用正交試驗方法進一步分析地震波作用下空洞對隧道動力響應特性的影響規(guī)律。正交設計方案采用L9(34)正交表，考慮4種因素，以隧道拱頂?shù)降乇淼木嚯x表示因素A(埋深)；以空洞與隧道中心的距離表示因素B(空洞與隧道中心距)；以空洞大小表示因素C(空洞半徑)；以地震波的入射角度表示因素D(入射角)，每種因素設定3個水平。表2表示正交設計方案，表3表示正交設計工況。

表2 正交設計方案Table 2 Orthogonal design scheme

表3 正交設計工況Table 3 Orthogonal design condition

采用解析解計算上述9種工況，得到的位移時程曲線如圖6所示。由于中位數(shù)不受極端數(shù)值影響，采用中位數(shù)作為所研究數(shù)據(jù)的代表值，表示每種工況位移時程曲線的一般水平。表4表示不同工況下隧道豎向位移代表值。由表4可知，中位數(shù)位移變化分布在11.5～23.6 cm之間。

表4 隧道豎向位移正交試驗結果Table 4 Orthogonal test results of vertical displacement of tunnel arch roof

圖6 各工況隧道豎向位移時程曲線Fig.6 Time-displacement curves of vertical displacement of tunnel in various working conditions

4.2 正交工況數(shù)據(jù)分析

通過極差分析可獲得不同因素水平變化對試驗指標的影響，其中極差由正交試驗中各因素水平計算結果的最大平均值和最小平均值相減得到，表5表示極差分析結果。通過方差分析可以研究不同因素對隧道拱頂豎向位移影響力的大小，以此確定具有顯著影響的因素。在置信水平α=0.1的條件下，分別對上述試驗結果進行方差分析，表6表示方差分析結果。

表5 各因素極差分析Table 5 Extremum difference analysis of various factors

表6 各因素方差分析Table 6 Variance analysis of various factors

對各因素極差分析結果運用正交試驗中的直觀分析法進行處理，得到各因素影響隧道豎向位移變化的直觀分析圖，位移敏感性分析結果如圖7所示。其中因素D極差最大，因素C的極差最小。各因素的敏感性由大到小依次為D＞B＞A＞C，對應入射角＞空洞與隧道中心距＞埋深＞空洞半徑，表明入射角對隧道豎向位移起控制作用，入射角度為0°時，對隧道結構抗震影響較大。埋深和空洞與隧道中心距在不同水平情況下對位移的影響規(guī)律較接近，隨著埋深增加以及空洞與隧道中心距的減少，隧道拱頂位移逐漸增大；空洞半徑的改變對位移變化幅度影響較小，當空洞直徑為0.3倍的隧道跨度時對隧道震動起到緩沖作用。

圖7 位移變化的敏感性分析Fig.7 Sensitivity analysis of displacement

5 結論

1)基于工程波動理論，采用波函數(shù)展開法給出了含空洞地層隧道在平面SV波激勵下動力響應位移解析解。

2)數(shù)值解和解析解位移時程曲線變化規(guī)律基本相似，曲線重合度較高，數(shù)值解位移幅值峰值與解析解位移幅值峰值相差9%左右，驗證了解析解的合理性。

3)地震波入射角對隧道豎向位移起控制作用；隨著埋深增加以及空洞與隧道中心距的減小，隧道拱頂位移逐漸增大；空洞直徑為0.3倍的隧道跨度時，空洞對隧道起到減震消能作用。