基于ALE算法的隧道開挖爆破振動特性數值分析

王 崢 崢， 張 楊 生

( 大連理工大學 土木工程學院, 遼寧 大連 116024 )

基于ALE算法的隧道開挖爆破振動特性數值分析

王 崢 崢*， 張 楊 生

( 大連理工大學 土木工程學院, 遼寧 大連 116024 )

基于ANSYS/LS-DYNA軟件，分析隧道開挖過程中爆破振動對圍巖及初期支護的影響．為了使數值模擬能夠真正反映實際情況，采用更精確合理的爆炸數值計算方法：利用軟件內置炸藥模塊和狀態方程模擬爆破荷載的作用，并采用ALE算法模擬炸藥與巖石之間的接觸關系．在ALE算法中，為防止爆炸過程中網格的過分畸變給結果帶來不利影響，將炸藥定義成流體．分析結果表明：應力、速度均在爆炸發生的極短時間內達到峰值，而后迅速衰減，10 ms 后達到穩定狀態．上臺階爆破在圍巖拱頂處產生的水平振速峰值為下臺階爆破的6倍左右，在拱腳位置約為0.88倍；上臺階爆破在圍巖拱頂處產生的豎直振速峰值為下臺階爆破的8倍左右．下臺階爆破在圍巖拱頂處產生的應力峰值是上臺階爆破的1/5，在拱腳處相差不大；同一位置，初期支護結構質點振速峰值與單元應力峰值均比圍巖大．

狀態方程；數值模擬；ALE算法

0 引 言

城市地鐵區間隧道開挖中多采用鉆爆法進行硬巖地區的施工，其中，臺階法具有適應性廣、造價低等優點，被廣泛應用于城市淺埋隧道爆破施工中[1]．地鐵線路經過城市繁華地帶，會穿過商業區和居民區，地面環境比較復雜，爆破產生的振動、噪聲會給周圍居民的日常生活帶來不便．另外，城市隧道埋深較淺，爆破產生的次生影響最主要的就是振動效應[2]．因此控制爆破振動效應是施工順利進行的有效保證[3-4]．

迄今，國內外學者對隧道圍巖的物理力學特性進行了大量的研究．譬如，Tiwari等建立了鋼筋混凝土襯砌隧道內部爆炸的三維非線性有限元分析模型，研究受到內部爆破荷載作用的風化巖體隧道的動態響應[5]；Wang等建立了能夠模擬爆炸和沖擊波在土中傳播的三相土壤模型[6-7]；運用該模型，Wang等進行了地下結構在地下爆炸作用下的結構響應的完全耦合數值模擬[8]；黃達等通過FLAC數值模擬與監測資料相結合，從不同位置、抗剪強度和地應力場等因素方面，較詳細地分析了地下洞室圍巖應力、變形及穩定性分布影響規律[9]；耿萍等將LS-DYNA軟件與振動臺試驗結合，在地震作用下，研究了穿越斷層破碎帶隧道襯砌結構的內力分布規律及其與斷層間的距離關系，驗證了數值模擬結論的準確性[10]．

本文以大連地鐵學苑廣場站—海事大學站區間隧道為背景，采用數值模擬方法，借助有限元軟件ANSYS/LS-DYNA建立隧道上、下臺階爆破開挖的三維有限元模型，分析爆破荷載作用下圍巖的振動特性以及下臺階爆破對上部初期支護的影響．

1 工程概況

2 數值模型

2.1 模型說明

采用ANSYS/LS-DYNA有限元軟件建立三維模型，將隧道掘進方向作為Z軸，豎直方向作為Y軸，水平方向作為X軸．結合圣維南原理以及實際工程情況，利用ANSYS/LS-DYNA軟件建立橫向長×高×縱向長=70 m×40 m×8 m的數值計算模型．計算模型如圖1所示．其中隧道區域網格進行了局部加密．

(a) 尺寸

(b) 網格

圖1 隧道計算模型

Fig.1 Calculation model of the tunnel

整個數值計算模型包含圍巖、炸藥和初期支護．本次計算只考慮隧道掏槽爆破，炸藥采取集中裝藥的方式，裝藥位置如圖1(a)所示．

2.2 邊界處理

整個數值模擬計算過程分為靜力分析和動力分析兩部分．其中，靜力分析時，模型前、后、左、右四個面均設置為水平約束，模型的底部設置為固定約束，模型上表面設為自由邊界，而且考慮重力作用．進行動力分析時，模型的上部、洞口所在面施加自由邊界條件，且不計上部荷載，其余4個面邊界施加固定約束，為防止應力波反射對計算結果的影響，全部設為無反射邊界[11]．

2.3 理論算法和材料選擇

2.3.1 ALE算法介紹 ALE算法的核心是方程的求解不受模型材料的變形和移動的影響，在其獨立的網格上進行計算．這種方法解決了材料變形過大對計算結果的影響，同時解決了移動邊界產生的一些難以解決的問題．

以網格中包含物質種類的多少為依據，ALE算法分為單物質算法和多物質算法兩種．多物質ALE算法中，物質的交換可以在相應單元中進行[12]．隧道爆破開挖下結構的動力響應涉及巖石、空氣和炸藥等多種物質，計算過程中物質發生較大變形，網格也相對變化較大，因此必須重視解決大變形帶來的影響．

ALE算法允許炸藥、巖石和空氣自由穿梭于彼此間，通過物質的流動可以分析爆炸發生的過程．其能夠很好地避免計算中網格變形過大導致的計算中斷問題．因此，ALE算法能夠很好地應用于爆炸分析計算中．

2.3.2 炸藥燃燒模型 LS-DYNA軟件本身就含有模擬炸藥爆炸的相應單元，通過狀態方程控制炸藥材料，在模擬炸藥爆炸時，壓力和比容的關系可以通過如下JWL狀態方程來描述：

(1)

式中：p為壓力；V為相對體積；E0為初始比熱力學能；參數A、B、R1、R2、ω為試驗確定的常數．本文選取的炸藥參數如表1所示．

表1 炸藥參數

2.3.3 巖石材料彈塑性模型 炸藥爆炸時靠近爆源的巖體有明顯的應變率效應，極易發生屈服而導致巖體破碎，采用包含應變率的塑性隨動模型(MAT_PLASTIC_KINEMATIC)經證明是比較合適的[13]．表2給出了巖石材料參數．

2.3.4 初期支護材料模型 初期支護應采用噴射混凝土的支護形式，模型采用軟件自帶的損傷本構模型(*MAT_JOHNSON_HOLMQUIST_CONCRETE)，材料參數見表3．

表2 巖石參數

表3 混凝土H-J-C模型材料參數

3 數值計算結果及其分析

3.1 控制節點與單元的選取

數值模擬的重點是對隧道開挖面拱頂及拱腳節點振速和單元應力的數據分析．關鍵測點及單元Z向坐標值為掌子面處，布置在拱頂及拱腳位置，測點位置標號1、2、3，如圖2所示．

圖2 控制點位置布置圖

3.2 模擬結果及分析

3.2.1 振動速度分析

(1)由于2、3號測點關于炸藥中心對稱，上臺階爆破后只提取1、2號兩個測點處圍巖節點振速時程曲線，如圖3所示．可以看出，在爆破應力波的作用下，節點振速在極短時間內達到峰值，而后迅速衰減．其中，拱頂處水平方向振速峰值為 0.17 m/s，到達時刻1 ms，豎直方向振速峰值為1.87 m/s，到達時刻2 ms；拱腳處水平方向振速峰值為0.35 m/s，到達時刻2 ms，豎直方向振速峰值為0.25 m/s，整個過程速度變化較?。茟Σóa生的振動主要集中在爆破發生后10 ms內．

(a) 水平方向

(b) 豎直方向

圖3 上臺階爆破后振速時程曲線

Fig.3 The time history curve of velocity under blasting on the upper-step

(2)下臺階爆破后提取2個測點處圍巖節點振速時程曲線，如圖4所示．由圖可以看出，在爆破應力波的作用下，節點振速在較短時間內達到峰值，而后迅速衰減至趨于0．其中，拱頂處水平振速峰值為0.03 m/s，整個過程振動幅度較小，豎直振速峰值為0.24 m/s；拱腳處水平振速峰值為0.40 m/s，豎直振速峰值為0.30 m/s．

(a) 水平方向

(b) 豎直方向

圖4 下臺階爆破后振速時程曲線

Fig.4 The time history curve of velocity under blasting on the lower-step

對比圖3和圖4，可以看出，下臺階爆破后拱頂處節點水平方向峰值振速約為上臺階爆破后產生的節點振速的0.18倍，豎直方向峰值振速約為上臺階爆破后產生的節點振速的0.13倍，爆破振動均主要集中在爆破后的10 ms內；拱腳處振動結果相比，下臺階爆破產生的峰值振速比上臺階爆破產生的峰值振速大一些．兩次爆破后，質點振速時程曲線變化趨勢相似，均是在短時間內達到振速峰值，后迅速衰減，最后趨近于0．

(3)下臺階爆破后提取2個測點對應初期支護處節點振速時程曲線，如圖5所示．由圖5可以看出，在1 ms左右，拱頂處水平方向振速峰值達到最大，約為0.05 m/s，整個過程振速接近于0，豎直方向振速峰值約為0.20 m/s；在3 ms左右，拱腳處水平振速峰值達到最大，約為0.62 m/s，15 ms以后趨于0，豎直方向振速峰值約為0.48 m/s，15 ms以后趨于0．

(a) 水平方向

(b) 豎直方向

圖5 下臺階爆破后支護振速時程曲線

Fig.5 The time history curve of the supporting velocity under blasting on the lower-step

對比圖4和圖5可以看出，下臺階爆破后，同一位置初期支護質點振速峰值大于圍巖質點振速峰值，為其1.5倍左右．

3.2.2 單元應力分析

(1)上臺階爆破后提取2個測點處圍巖單元應力時程曲線，如圖6(a)所示．由圖可以看出，在爆破應力波的作用下，單元應力在極短時間內達到峰值，而后迅速衰減，短時間內達到穩定狀態．其中，拱頂處單元應力峰值為14.4 MPa，到達時刻2 ms，拱腳處單元應力峰值為12.6 MPa，到達時刻2 ms，爆破應力波產生的應力變化主要集中在爆破發生后的10 ms內．

(2)下臺階爆破后提取2個測點處圍巖單元應力時程曲線，如圖6(b)所示．由圖可以看出，在爆破應力波的作用下，拱頂單元應力在4 ms達到峰值3.1 MPa，而后迅速衰減至1 MPa左右，達到穩定狀態．拱腳處單元應力在2 ms達到峰值11.0 MPa，是拱頂處的3.5倍左右，而后迅速衰減至5 MPa以下，達到穩定狀態．

(a) 上臺階爆破

(b) 下臺階爆破

圖6 上、下臺階爆破后應力時程曲線

Fig.6 The time history curve of stress under blasting on the upper-step and lower-step

對比圖6(a)和(b)可以看出，下臺階爆破后拱頂處單元應力峰值約為上臺階爆破后產生的單元應力峰值的0.22倍，到達峰值的時刻延遲2 ms；拱腳處單元應力結果相比，下臺階爆破產生的峰值應力是上臺階爆破產生的峰值應力的0.87倍左右，到達峰值時刻相同．單元應力變化均集中在爆破后的10 ms內．

(3)下臺階爆破后提取2個測點對應初期支護處單元應力時程曲線，如圖7所示．由圖可以看出，在爆破應力波的作用下，拱頂單元應力在7 ms 內達到峰值23 MPa，10 ms后衰減至20 MPa，達到穩定狀態．拱腳處單元應力在3 ms內達到峰值30 MPa，是拱頂處的1.3倍左右，最后穩定在10 MPa上下波動．

圖7 下臺階爆破后支護應力時程曲線

對比圖6(b)和圖7可以看出，下臺階爆破后，同一位置圍巖單元應力峰值遠小于初期支護單元應力峰值．

提取了圍巖及初期支護結構關鍵位置處節點振速時程曲線和單元應力時程曲線，表4給出了所有情況下的峰值．由表可知，上臺階爆破時，圍巖拱頂節點水平方向振速遠小于豎直方向振速，而拱腳節點水平方向振速大于豎直方向振速；單元應力相差1.8 MPa，應重點監控拱頂豎直方向和拱腳水平方向振動速度；下臺階爆破時，圍巖、支護拱頂節點振速和單元應力均小于拱腳位置，同一位置圍巖應力遠小于支護應力，應重點監控支護結構上的應力變化．

表4 峰值結果對比表

4 結 論

(1)上臺階爆破產生的節點水平、豎直振速峰值在圍巖拱頂處遠大于下臺階爆破產生的節點振速峰值，約為其5～8倍．爆破后，圍巖振速在極短時間內達到峰值，而后迅速衰減，振動持續10 ms以后，振速趨于0．

(2)上臺階爆破產生的單元應力峰值在圍巖拱頂處大于下臺階爆破產生的單元應力峰值，約為其5倍；在拱腳處，下臺階爆破產生的單元應力峰值與上臺階爆破產生的單元應力峰值相差微小，僅差1.6 MPa．爆破后，圍巖單元應力在極短時間內達到峰值，而后迅速衰減，10 ms后達到穩定狀態．

(3)下臺階爆破后，同一位置初期支護產生的振速峰值比圍巖的大，約為其1.5倍；單元應力峰值比圍巖的大，約為其3～7倍．

(4)根據模擬結果，在上臺階爆破施工時，應加強拱頂、拱腳位置的監測和超前支護；下臺階爆破施工時，應加強支護結構的監測和保護．

[1] 閆鴻浩,王小紅. 城市淺埋隧道爆破原理及設計[M]. 北京:中國建筑工業出版社, 2013.

YAN Honghao, WANG Xiaohong. Urban Shallow Tunnel Blasting Theory and Design [M]. Beijing: China Architecture & Building Press, 2013. (in Chinese)

[2] 劉 琳,明 峰. 地鐵隧道掘進爆破振動效應研究[J]. 公路, 2012(8):262-266.

LIU Lin, MING Feng. The research on subway tunnel excavation blasting vibration effect [J]. Highway, 2012(8):262-266. (in Chinese)

[3] 林學義,王蘭民. 隧洞的爆破地震動效應問題[C]// 工程爆破文集第六集. 深圳:海天出版社, 1997.

LIN Xueyi, WANG Lanmin. The tunnel blasting seismic effect [C] // Proceedings of Engineering Explosive: vol.6. Shenzhen: Haitian Publishing House, 1997. (in Chinese)

[4] 錢七虎,陳士海. 爆破地震效應[J]. 爆破, 2004, 21(2):1-5.

QIAN Qihu, CHEN Shihai. Blasting vibration effect [J]. Blasting, 2004, 21(2):1-5. (in Chinese)

[5] TIWARI R, CHAKRABORTY T, MATSAGAR V. Dynamic analysis of a twin tunnel in soil subjected to internal blast loading [J]. Indian Geotechnical Journal, 2016, 46(4):369-380.

[6] WANG Zhongqi, LU Yong. Numerical analysis on dynamic deformation mechanism of soils under blast loading [J]. Soil Dynamics and Earthquake Engineering, 2003, 23:705-714.

[7] WANG Zhongqi, HAO Hong, LU Yong. A three-phase soil model for simulating stress wave propagation due to blast loading [J]. International Journal for Numerical and Analytical Methods in Geomechanics, 2004, 28(1):33-56.

[8] WANG Zhongqi, LU Yong, HAO Hong,etal. A full coupled numerical analysis approach for buried structures subjected to surface blast [J]. Computers & Structures, 2005, 83:339-356.

[9] 黃 達,黃潤秋,張永興. 斷層位置及強度對地下洞室圍巖穩定性影響[J]. 土木建筑與環境工程, 2009, 31(2):68-73.

HUANG Da, HUANG Runqiu, ZHANG Yongxing. Analysis on influence of fault location and strength on deformation and stress distribution of surrounding rocks of large underground openings [J]. Journal of Civil, Architectural & Environmental Engineering, 2009, 31(2):68-73. (in Chinese)

[10] 耿 萍,吳 川,唐金良,等. 穿越斷層破碎帶隧道動力響應特性分析[J]. 巖石力學與工程學報, 2012, 31(7):1406-1413.

GENG Ping, WU Chuan, TANG Jinliang,etal. Analysis of dynamic response properties for tunnel through fault fracture zone [J]. Chinese Journal of Rock Mechanics and Engineering, 2012, 31(7): 1406-1413. (in Chinese)

[11] 何忠明,蔡 軍,王利軍,等. 地鐵隧道爆破開挖引起地表位移沉降的數值模擬[J]. 中南大學學報(自然科學版), 2015, 46(9):3496-3502.

HE Zhongming, CAI Jun, WANG Lijun,etal. Numerical simulation of ground subsidence due to tunnel blasting excavation [J]. Journal of Central South University (Science and Technology), 2015, 46(9):3496-3502. ( in Chinese)

[12] 賈祖朋. 基于MOF界面重構的多物質ALE方法[J]. 計算物理, 2010, 27(3):353-360.

JIA Zupeng. A multi-material arbitrary Lagrangian-Eulerian method based on MOF interface reconstruction [J]. Chinese Journal of Computational Physics, 2010, 27(3):353-360. (in Chinese)

[13] 夏 祥. 爆炸荷載作用下巖體損傷特征及安全閾值研究[D]. 武漢:中國科學院武漢巖土力學研究所, 2006.

XIA Xiang. Study on damage characteristic and safety threshold of rock vibration by blast [D]. Wuhan:Institute of Rock and Soil Mechanics of the Chinese Academy of Sciences, 2006. (in Chinese)

Numerical analysis for blasting vibration characteristics in tunnel excavation based on ALE algorithm

WANG Zhengzheng*, ZHANG Yangsheng

( School of Civil Engineering, Dalian University of Technology, Dalian 116024, China )

Based on ANSYS/LS-DYNA software, the influences of blasting vibration in tunnel excavation on surrounding rock and initial support structure are analyzed. In order to make the numerical simulation truly reflect the actual situation, a more accurate method of explosion simulation is adopted. The built-in explosives module in software and equation of state are used to simulate the explosive loads. In the meantime, the ALE algorithm is used to simulate the contact relationship between the structure of explosives and rocks. In ALE algorithm, explosives are defined into the fluid in order to avoid excessive grid distortion caused by explosion, which adversely affects the result of the calculation. The analytical results show that: Stress and velocity all reach a maximum within a very short period of time after the blasting, and then decay to the stable state at 10 ms. The peak value of horizontal vibration velocity at the dome roof of the surrounding rock of the upper-step blasting is about 6 times that of the lower-step blasting, and about 0.88 times at the arch foot position. The vertical vibration velocity peak value at the dome roof of the surrounding rock of upper-step blasting is about 8 times that of the lower-step blasting. The peak value of stress generated by the lower-step blasting at the dome roof of the surrounding rock is 1/5 of the upper-step blasting, and the difference is small at the arch foot. At the same position, the peak value of mass point vibration velocity and element stress of primary support structure are larger than those of surrounding rock.

equation of state; numerical simulation; ALE algorithm

1000-8608(2017)03-0279-06

2016-08-15；

2017-04-10．

住建部科技項目(2014-K5-005)；大連市建委科技項目(052033).

王崢崢*(1982-)，男，博士，副教授，E-mail:wangzhengzheng@dlut.edu.cn；張楊生(1989-)，男，碩士生，E-mail:1120837032@qq.com.

TU443

A

10.7511/dllgxb201703010