基于離散元的節(jié)理巖體水工隧洞地震動(dòng)力特性分析

包 磊，陳春武，潘 昆，羅曉琴，張 金，雷進(jìn)生

(1.國(guó)家林業(yè)局 昆明勘察設(shè)計(jì)院，昆明 650216;2.三峽大學(xué)土木與建筑學(xué)院，湖北宜昌 443002)

1 研究背景

天然巖體通?？煽醋魇怯赏暾麕r石和不連續(xù)面如層面、節(jié)理斷層等組成，不連續(xù)面的存在使得巖體的完整性和均勻性受到不同程度的破壞，巖體中的這些不連續(xù)面對(duì)巖體在靜力和動(dòng)力荷載作用下的力學(xué)行為均起到主導(dǎo)作用［1］。顯然，基于連續(xù)介質(zhì)理論的有限元方法難以克服因建立的分析模型與巖體介質(zhì)客觀屬性之間存在的差異對(duì)模擬結(jié)果產(chǎn)生的不利影響。目前關(guān)于隧洞圍巖穩(wěn)定性的研究主要集中于靜力分析，對(duì)一般隧洞在動(dòng)荷載作用下的研究也取得了一系列成果，但對(duì)于大跨度、高邊墻的節(jié)理巖體隧洞在動(dòng)力作用下的圍巖力學(xué)響應(yīng)研究較少;加之已有的研究大都基于小變形的假設(shè)，事實(shí)上，在高應(yīng)力或動(dòng)荷載的作用下，節(jié)理巖體要產(chǎn)生剪切滑移、離層脫落甚至失穩(wěn)垮塌，這些現(xiàn)象已非連續(xù)介質(zhì)所能描述。因此，本文主要討論地震作用下大跨度、高邊墻且含節(jié)理裂隙等軟弱結(jié)構(gòu)面條件下的水工隧洞圍巖的力學(xué)響應(yīng)，獲得大型地下隧洞不同部位圍巖在地震中的響應(yīng)規(guī)律，以期為同類大型地下工程的設(shè)計(jì)和施工及進(jìn)一步分析地震作用下圍巖的變形破壞過程提供參考。

2 動(dòng)力離散元求解原理

UDEC是一個(gè)處理不連續(xù)介質(zhì)的二維離散元程序，可用于模擬非連續(xù)介質(zhì)承受靜載或動(dòng)載作用下的響應(yīng)。非連續(xù)介質(zhì)被認(rèn)為是由離散的塊體集合體表示，塊體間的邊界接觸面以不連續(xù)面處理，允許各個(gè)離散塊體沿不連續(xù)面發(fā)生較大的位移和轉(zhuǎn)動(dòng)，而這恰好體現(xiàn)了節(jié)理巖體的重要變形和破壞機(jī)制。在離散單元法中，不連續(xù)面等結(jié)構(gòu)面在賦予相關(guān)節(jié)理參數(shù)后，主要采用Jset命令生成節(jié)理組。

在離散元法中，塊體的運(yùn)動(dòng)方程都可歸結(jié)為如下形式［2－3］，即

設(shè)an，vn和Fn分別為tn時(shí)刻的加速度、速度和力，開始時(shí)刻 t0=0，時(shí)步為Δt，則tn=t0+nΔt，可將式(1)改寫為

用顯式的Newmark法對(duì)上式積分，其過程如下:

則在時(shí)刻tn的加速度可由式(3)和式(4)求出，即

類似地可求得tn時(shí)刻的速度為

則易求塊體在tn+1時(shí)刻的位置，即

如此往復(fù)，利用上述數(shù)值法求解非線性微分方程(1)，中間時(shí)刻t(2n+1)/2時(shí)的速度是一個(gè)重要的中間變量，通過它就很容易求得塊體的新位置和當(dāng)前時(shí)刻的加速度。

根據(jù)上述的顯式積分方法，依次進(jìn)行3個(gè)主要的循環(huán)計(jì)算，即所有塊體→所有接觸→所有塊體間空洞。在輸入的動(dòng)載荷引入之前，要一直重復(fù)上述的循環(huán)計(jì)算。

在上述的數(shù)值方法中，要保證計(jì)算的收斂和穩(wěn)定，對(duì)于計(jì)算時(shí)步在小于一個(gè)限值的情況下才能滿足要求，其物理含義為:在非常短的時(shí)間內(nèi)，系統(tǒng)中某一點(diǎn)的擾動(dòng)只對(duì)與該點(diǎn)相鄰的單元有影響，而對(duì)系統(tǒng)中遠(yuǎn)離該擾動(dòng)點(diǎn)的其它點(diǎn)不會(huì)產(chǎn)生影響［4－5］。對(duì)于這樣的非線性問題的求解，Cundall給出了時(shí)步的相限值估算公式［5］，即

式中:Tmin為任一單元的最小固有振動(dòng)周期;mi，ki分別為塊體質(zhì)量和剛度;n為塊體數(shù)。

3 離散元?jiǎng)恿τ?jì)算模型的建立

以某大型水電工程尾水隧洞為工程背景，尾水主洞標(biāo)準(zhǔn)開挖斷面為22 m×36 m(寬×高)的城門洞型，為典型的大跨度、高邊墻工程，且地質(zhì)條件復(fù)雜。尾水隧洞位于，，緩傾角層狀巖體中，以厚至巨厚層砂巖為主，局部呈薄至中厚層狀，絕大部分為微風(fēng)化—新鮮巖體，以砂巖為主，圍巖類型以Ⅱ—Ⅲ類為主，巖石堅(jiān)硬，強(qiáng)度較高。廠房區(qū)無較大斷層通過，巖體中發(fā)育的結(jié)構(gòu)面主要有層面、層間軟弱夾層和節(jié)理裂隙，節(jié)理裂隙主要有2組，優(yōu)勢(shì)產(chǎn)狀分別為:275°～286°/NE∠55°～65°(J1)，70°～80°/NW∠62°～74°(J2)，節(jié)理裂隙發(fā)育段為Ⅲ類。地下廠房廠區(qū)實(shí)測(cè)地應(yīng)力的最大主應(yīng)力平均值為10 MPa左右，屬中等低應(yīng)力場(chǎng)，因此，巖體的動(dòng)力響應(yīng)仍以地質(zhì)結(jié)構(gòu)面控制型為主。

3.1 模型范圍

模型邊界必須足夠遠(yuǎn)，以致模型對(duì)邊界不產(chǎn)生影響。一般對(duì)單一地下開挖工程，邊界離開挖邊界的距離應(yīng)不小于開挖跨度的5倍［6］。為降低過小的計(jì)算尺寸對(duì)模擬分析結(jié)果的影響，模型x方向計(jì)算寬度取為300 m，y方向計(jì)算高度取250 m，其中洞室底板距模型下邊界為100 m。因地下隧洞埋深較淺，構(gòu)造應(yīng)力不明顯，按不利情況考慮，假設(shè)初始地應(yīng)力場(chǎng)僅涉及自重應(yīng)力場(chǎng)。對(duì)隧洞上方300 m以上至地表的巖體則計(jì)算出其重量作為外荷載施加于上邊界上。建立UDEC模型如圖1所示。

圖1 開挖后的計(jì)算模型Fig.1 Calculation model of tunnel after excavation

3.2 地震荷載輸入

考慮到剪切波對(duì)洞室響應(yīng)影響較大，本節(jié)采用速度形式的垂直向剪切正弦波作為動(dòng)力荷載輸入至模型底部，速度時(shí)程表征為v=λsin(2πft)，取速度峰值 vs=0.1 m/s，頻率 f=50 Hz，持時(shí) t=0.1 s，計(jì)算0.2 s。波形如圖2所示。

圖2 正弦剪切波速度時(shí)程曲線Fig.2 Time-history curve of the velocity of sine shear wave

3.3 邊界條件和阻尼

邊界條件采用UDEC自帶的自由場(chǎng)黏性邊界，左右兩側(cè)的邊界通過自由場(chǎng)網(wǎng)格與模型之間的黏性阻尼筒耦合來實(shí)現(xiàn)，以吸收來自模型內(nèi)部的入射波，底部采用自由場(chǎng)邊界并限制 x向和 y向位移［7－8］。本節(jié)動(dòng)力荷載是速度歷程形式的正弦剪切波，因無法直接將這種形式的動(dòng)力荷載施加在不反射邊界上，此時(shí)，可以按式(10)和式(11)將速度歷程形式的動(dòng)力荷載轉(zhuǎn)換為式(12)應(yīng)力歷程形式的動(dòng)力荷載［4］，轉(zhuǎn)換后相當(dāng)于在模型底部施加了一個(gè)幅值約為1 MPa的剪切應(yīng)力波。

式中:σn為施加的法向應(yīng)力荷載(MPa);σs為施加的切向應(yīng)力荷載(MPa);ρ為介質(zhì)的密度(kg/m3);CP為應(yīng)力波在介質(zhì)中傳播的P波(縱波)波速(m/s);CS為應(yīng)力波在介質(zhì)中傳播的S波(橫波)波速(m/s);vn為鉛直方向的質(zhì)點(diǎn)速度(m/s);vs為水平方向的質(zhì)點(diǎn)速度(m/s)。

對(duì)于巖土動(dòng)力非線性問題，阻尼參數(shù)的確定通常采用Rayleigh阻尼，圍巖的阻尼比一般為0.01～0.03［9］，本文阻尼比取0.01進(jìn)行計(jì)算。

3.4 本構(gòu)關(guān)系

數(shù)值計(jì)算中，對(duì)于巖石，采用理想彈塑性模型，以Mohr－Coulomb強(qiáng)度準(zhǔn)則作為其屈服準(zhǔn)則，屈服函數(shù)如下:

式中:σ1，σ3分別為最大、最小主應(yīng)力;φ為內(nèi)摩擦角;c為黏聚力;σt為巖石的抗拉強(qiáng)度。Nφ的表達(dá)式為

當(dāng)巖體內(nèi)某點(diǎn)應(yīng)力滿足fs＞0時(shí)，發(fā)生剪切破壞;當(dāng)滿足ft＞0時(shí)，產(chǎn)生拉伸破壞。

對(duì)于巖體節(jié)理，適用彈性模型，屈服準(zhǔn)則為接觸Coulomb 滑動(dòng)準(zhǔn)則［10］，即

式中σ，τ分別為節(jié)理面正應(yīng)力和剪應(yīng)力。

3.5 巖體物理力學(xué)指標(biāo)

一般而言，巖石節(jié)理和材料的動(dòng)態(tài)物理力學(xué)特性參數(shù)都是與應(yīng)變速率相關(guān)的，但對(duì)于巖體的變形模量E和泊松比ν等參數(shù)隨應(yīng)變率的變化幅度很小，通常用相應(yīng)的靜態(tài)參數(shù)表示［11－12］，故動(dòng)力分析巖體物理力學(xué)參數(shù)取值同靜力分析。結(jié)合工程勘察資料分別給出水工隧洞各巖層及結(jié)構(gòu)面的巖體力學(xué)參數(shù)，見表1和表2。

表1 尾水洞洞室圍巖的物理力學(xué)參數(shù)Table 1 Physical and mechanical parameters of surrounding rock of the tailrace tunnel

表2 尾水洞洞室結(jié)構(gòu)面物理力學(xué)參數(shù)Table 2 Physical and mechanical parameters of structural planes of the tailrace tunnel

4 計(jì)算結(jié)果與分析

計(jì)算過程包括靜力和動(dòng)力計(jì)算2部分，動(dòng)力分析以靜力計(jì)算為前提，計(jì)算步驟如下:

(1)首先靜力計(jì)算。施加自重和構(gòu)造荷載作用使模型固結(jié)，以模擬初始地應(yīng)力場(chǎng)的形成。

(2)其次開挖、支護(hù)。分層分步開挖、支護(hù)應(yīng)使模型達(dá)到平衡狀態(tài)，直至整個(gè)洞室開挖支護(hù)完成為止。

(3)最后動(dòng)力分析。施加地震作用力，研究大跨度、高邊墻且含節(jié)理裂隙等軟弱結(jié)構(gòu)面條件下的水工隧洞圍巖的動(dòng)力響應(yīng)。為分析洞室的開挖對(duì)地震波在介質(zhì)中傳播速度的影響，在模型中部水平方向和垂直方向各選取10個(gè)點(diǎn)進(jìn)行監(jiān)測(cè)。各監(jiān)測(cè)點(diǎn)布置如圖3所示。

4.1 動(dòng)荷載對(duì)隧洞圍巖動(dòng)力響應(yīng)分析

地下洞室圍巖中各監(jiān)測(cè)點(diǎn)的速度時(shí)程曲線可反映出介質(zhì)對(duì)動(dòng)力荷載的實(shí)時(shí)響應(yīng)。通過對(duì)頂拱、邊墻及底板監(jiān)測(cè)點(diǎn)的速度監(jiān)測(cè)發(fā)現(xiàn)，剪切波在0.01 s后傳至底板部位，在0.02 s后傳至拱頂部位，在0.015 s后傳至左右邊墻部位;在垂直方向，頂拱速度幅值約為0.043 m/s，小于底板速度幅值為0.086 m/s;水平方向，左邊墻速度幅值約為0.026 5 m/s，右邊墻速度幅值約為0.025 0 m/s。洞壁監(jiān)測(cè)點(diǎn)在考慮垂直方向正弦剪切地震波的作用下，各監(jiān)測(cè)點(diǎn)速度歷程反映了輸入地震波的振動(dòng)特性，且豎向動(dòng)力特性明顯比水平方向大。

圖4(a)為水平方向10個(gè)監(jiān)測(cè)點(diǎn)H10—H19的速度歷程曲線，圖4(b)為豎直方向10個(gè)監(jiān)測(cè)點(diǎn)V20—V29的速度歷程曲線。

圖3 監(jiān)測(cè)點(diǎn)布置Fig.3 Arrangement of monitoring points

圖4 不同位置監(jiān)測(cè)點(diǎn)x和y方向的速度時(shí)程曲線Fig.4 Time-history curves of wave velocity in direction x and y of monitoring points

結(jié)合圖3及圖4可以看出，地下洞室洞周附近點(diǎn)速度幅值明顯大于遠(yuǎn)離開挖輪廓面附近介質(zhì)監(jiān)測(cè)點(diǎn)的速度幅值，反映了地下洞室開挖斷面的存在對(duì)剪切波具有放大效應(yīng)。由此可見，形成開挖面的地下洞室在動(dòng)荷載作用下的穩(wěn)定性程度降低。整體而言，地下洞室洞監(jiān)測(cè)點(diǎn)速度幅值變化反映了波在介質(zhì)中傳播的衰減規(guī)律特性，即離波源越遠(yuǎn)的監(jiān)測(cè)點(diǎn)由于所受波擾動(dòng)小，因而速度幅值也較小。

4.2 動(dòng)荷載對(duì)隧洞圍巖位移影響

地震荷載作用下，大跨度、高邊墻節(jié)理巖體中，洞室開挖斷面對(duì)圍巖位移分布的影響可通過對(duì)開挖輪廓線附近的監(jiān)測(cè)點(diǎn)的位移歷程曲線來反映。圖5為水平方向監(jiān)測(cè)點(diǎn)x方向位移歷程曲線(洞室左邊墻監(jiān)測(cè)點(diǎn):H10—H14;右邊墻監(jiān)測(cè)點(diǎn):H15—H19)，圖6為豎直方向監(jiān)測(cè)點(diǎn)y方向位移歷程曲線(洞室底板部位監(jiān)測(cè)點(diǎn):V20—H24;頂拱部位監(jiān)測(cè)點(diǎn):V25—V29)。

圖5 水平方向監(jiān)測(cè)點(diǎn)x方向位移時(shí)程曲線Fig.5 Time－h(huán)istory curves of displacement in direction x of monitoring points at horizontal direction

由圖5和圖6可以看出，在靠近開挖斷面，位移值明顯大于遠(yuǎn)離開挖面部位，頂拱、底板及邊墻部位均表現(xiàn)出相似的規(guī)律。分析認(rèn)為，洞室的開挖，使原處于一定地質(zhì)環(huán)境中的巖體應(yīng)力釋放，巖體卸荷回彈，越靠近開挖面，卸荷量越大，圍巖變形也越明顯。

地震波作用后，隧洞圍巖位移場(chǎng)分布如圖7所示。震后圍巖最大位移為0.297 mm，圍巖在經(jīng)過最優(yōu)開挖支護(hù)方案的實(shí)施以后，整體效果明顯，位移值增量不大，位移矢量方向多指向開挖臨空面。由于地震的頻率較集中于低頻，對(duì)長(zhǎng)周期結(jié)構(gòu)物或構(gòu)筑物有較大的影響，而巖體的自振頻率較小，在介質(zhì)自振頻率和地震波頻率接近的情況下達(dá)到共振，此時(shí)洞室的響應(yīng)破壞也最為明顯，本節(jié)輸入地震波頻率超過介質(zhì)共振頻率，故其對(duì)圍巖位移場(chǎng)影響較小。

圖6 豎直方向監(jiān)測(cè)點(diǎn)y方向位移時(shí)程曲線Fig.6 Time－h(huán)istory curves of displacement in direction y of monitoring points at vertical direction

圖7 震后圍巖位移矢量圖Fig.7 Displacement vectors of surrounding rock after earthquake

4.3 動(dòng)荷載對(duì)隧洞圍巖應(yīng)力影響

表3給出了各監(jiān)測(cè)震前、震后點(diǎn)水平應(yīng)力σxx和豎向應(yīng)力σyy的應(yīng)力幅值(負(fù)號(hào)表示受壓)。

表3 各監(jiān)測(cè)點(diǎn)震前、震后水平應(yīng)力和豎向應(yīng)力幅值Table 3 Horizontal stress and vertical stress amplitudes of monitoring points before and after earthquake

由表3可以看出，豎向正弦剪切地震波作用后，洞室各監(jiān)測(cè)點(diǎn)應(yīng)力值相對(duì)地震波作用前均有不同程度的增大，由于施加的動(dòng)力荷載幅值量級(jí)不大，約為1 MPa，所以受動(dòng)載影響應(yīng)力增大幅度不大。不難發(fā)現(xiàn)，相對(duì)震前豎向地震波的作用對(duì)頂拱和底板的水平應(yīng)力值的影響程度明顯大于其對(duì)豎向應(yīng)力的影響;而在邊墻部位，地震波作用后，監(jiān)測(cè)點(diǎn)水平應(yīng)力卻有一定程度減小，其對(duì)監(jiān)測(cè)點(diǎn)豎向應(yīng)力值影響不大。因此，在隧洞受豎向地震作用下，頂拱和底板部位應(yīng)加強(qiáng)支護(hù)和監(jiān)測(cè)，以確保洞室的安全。

洞室底板和左右邊墻中點(diǎn)監(jiān)測(cè)點(diǎn)在地震波作用下，先后分別表現(xiàn)出壓應(yīng)力和拉應(yīng)力，然而出現(xiàn)拉應(yīng)力時(shí)的量值小于出現(xiàn)壓應(yīng)力時(shí)的量值，局部微小拉應(yīng)力的出現(xiàn)因范圍和量值均很小，不足以使圍巖和節(jié)理破壞，對(duì)圍巖整體穩(wěn)定影響不大。

應(yīng)力場(chǎng)計(jì)算顯示，地震作用前后圍巖拉應(yīng)力及壓應(yīng)力均有一定程度增加，但幅度不大，洞室處于穩(wěn)定狀態(tài)。

4.4 動(dòng)荷載對(duì)隧洞圍巖塑性區(qū)影響

地震波作用前后，隧洞圍巖塑性區(qū)分布如圖8所示。與震前相比，震后的圍巖塑性區(qū)有較大幅度增加，數(shù)量由震前的28增至115，增幅為4倍，對(duì)比發(fā)現(xiàn)，震后隧洞邊墻及底板巖體塑性區(qū)范圍擴(kuò)大較明顯，地震作用前后拉伸破壞數(shù)量均為0。

圖8 地震波作用前后圍巖塑性區(qū)分布Fig.8 Distribution of plastic zones of surrounding rock before and after earthquake

5 結(jié)論

(1)地震波在節(jié)理巖體中傳播符合波的衰減規(guī)律;在垂直方向正弦剪切地震波的作用下，各監(jiān)測(cè)點(diǎn)豎向動(dòng)力特性明顯比水平方向大;同時(shí)，地下洞室洞周附近點(diǎn)速度幅值明顯大于遠(yuǎn)離開挖面附近介質(zhì)監(jiān)測(cè)點(diǎn)的速度幅值，即地下洞室開挖斷面的存在對(duì)剪切波具有放大效應(yīng)。

(2)在地震波作用過程中，圍巖介質(zhì)質(zhì)點(diǎn)位移幅值逐漸增大，地震波作用完畢，在巖體介質(zhì)阻尼的作用下，各監(jiān)測(cè)點(diǎn)位移幅值發(fā)生了很小的永久變形后趨于穩(wěn)定。在垂直方向正弦剪切地震波的初期，邊墻產(chǎn)生的位移值大于頂拱，隨時(shí)間的推移，邊墻受地震力作用變小，從而產(chǎn)生的位移變化量小，而底板圍巖變形受地震作用較大，即底板的變形更為敏感;由洞室斷面附近的水平和豎直方向監(jiān)測(cè)點(diǎn)位移歷程可以看出，洞室的開挖，使原處于一定地質(zhì)環(huán)境中的巖體應(yīng)力釋放，巖體卸荷回彈，在靠近開挖斷面，位移值明顯大于遠(yuǎn)離開挖面部位，頂拱、底板及邊墻部位均表現(xiàn)出相似的規(guī)律;而圍巖整體位移由于介質(zhì)自振頻率小于輸入地震波頻率，故地震波對(duì)圍巖位移場(chǎng)影響較小。

(3)豎向正弦剪切地震波作用后，洞室各監(jiān)測(cè)點(diǎn)應(yīng)力值相對(duì)地震波作用前均有較小程度的增大。相對(duì)震前情況，豎向地震波對(duì)頂拱和底板的水平應(yīng)力值的影響程度明顯大于其對(duì)豎向應(yīng)力的影響;而在邊墻部位，地震波作用后，監(jiān)測(cè)點(diǎn)水平應(yīng)力卻有一定程度減小，其對(duì)監(jiān)測(cè)點(diǎn)豎向應(yīng)力值影響不大，頂拱和底板部位應(yīng)加強(qiáng)支護(hù)和監(jiān)測(cè);地震作用前后圍巖最大拉應(yīng)力及最大壓應(yīng)力均有一定程度增加，但幅度不大。

(4)與震前相比，震后的圍巖塑性區(qū)有較大幅度增加，震后隧洞邊墻及底板巖體是塑性區(qū)范圍擴(kuò)大較明顯的部位。

6 討論

UDEC是一個(gè)處理不連續(xù)介質(zhì)的二維離散元程序，可用于模擬非連續(xù)介質(zhì)承受靜載或動(dòng)載作用下的響應(yīng)，非連續(xù)介質(zhì)被認(rèn)為是由離散的塊體集合體表示，塊體間的邊界接觸面以不連續(xù)面處理，允許各個(gè)離散塊體沿不連續(xù)面發(fā)生較大的位移和轉(zhuǎn)動(dòng)，而這些恰好體現(xiàn)了節(jié)理巖體的重要變形和破壞機(jī)制。

本文研究獲得大型地下隧洞不同部位圍巖在地震中的響應(yīng)規(guī)律，為進(jìn)一步分析地震作用下圍巖的變形破壞提供參考，可為類似大型地下洞室的工程設(shè)計(jì)和施工提供借鑒和參考。建議下一步研究工作如將離散單元法計(jì)算結(jié)果與其他較成熟有限元方法所得結(jié)論對(duì)比論證，使計(jì)算結(jié)果更符合客觀實(shí)際。

［1］GOODMAN R E.Methods of Geological Engineering in Discontinuous Rocks［J］.West Publishing，1976，8(2):135－140.

［2］套煉金，常 春.節(jié)理巖體的動(dòng)力響應(yīng)分析及工程應(yīng)用［J］.黑龍江礦業(yè)學(xué)院學(xué)報(bào)，2000，10(4):6－10.(TAO Lian-jin，CHANG Chun.Dynamic Response A-nalysis of Jointed Rock Mass and Its Application in Engineering［J］.Journal of Heilongjiang Mining Institute，2000，10(4):6－10.(in Chinese))

［3］陶連金，張悼元，王泳嘉.復(fù)雜工程巖休穩(wěn)定性評(píng)價(jià)的方法與實(shí)踐［M］.成都:成都科技大學(xué)出版社，1998.(TAO Long-jin，ZHANG Dao-yuan，WANG Yong-jia.Stability Evaluation for Complex Engineering Rockmass:Methods and Practice［M］.Chengdu:Chengdu University of Science and Technology Press，1988.(in Chinese))

［4］Itasca Consulting Group.UDEC(Universal Distinct Element Code)User’s Manual(Version 3.0)［M］.USA:Itasca Consulting Group，Inc.，1996.

［5］王泳嘉，邢紀(jì)波.離散元法及其在巖土力學(xué)中的應(yīng)用［M］.沈陽:東北工學(xué)院出版社，1991.(WANG Yong-jia，XINGJi-bo.The Discrete Element Method and Its Application in Rock and Soil Mechanics［M］.Shenyang:Northeast Institute of Technology Press，1991.(in Chinese))

［6］王 濤，熊 將，郭武祥，等.地震荷載作用下地下廠房圍巖穩(wěn)定的離散元計(jì)算方法研究［J］.長(zhǎng)江科學(xué)院院報(bào)，2009，26(12):58－61.(WANG Tao，XIONG Jiang，GUO Wu-xiang，et al.Research on Stability of Rock Mass around Underground Powerhouse Cavern by Discrete Element Method under Earthquake［J］.Journal of Yangtze River Scientific Research Institute，2009，26(12):58－61.(in Chinese))

［7］王魯明，趙 堅(jiān)，華安增，等.節(jié)理巖體中應(yīng)力波傳播規(guī)律研究的進(jìn)展［J］.巖土力學(xué)，2003，24(增):602－610.(WANG Lu-ming，ZHAO Jian，HUA An-zeng，et al.The Progress in Study of Regularity of a Stress Wave Propagation in the Jointed Rock Mass［J］.Rock and Soil Mechanics，2003，24(Sup.):602－ 610.(in Chinese))

［8］陶連金，張悼元，傅小敏.在地震載荷作用下的節(jié)理巖體地下洞室圍巖穩(wěn)定性分析［J］.中國(guó)地質(zhì)災(zāi)害與防治學(xué)報(bào)，1998，9(2):33－42.(TAO Lian-jin，ZHANG Zhuo-yuan，F(xiàn)U Xiao-min.Stability Analysis of Surround-ing Rockmass of Underground Excavation Rock Mass Under Seismic Load［J］.The Chinese Journal of Geological Hazard and Control，1988，9(2):33－ 42.(in Chinese))

［9］鄭穎人，肖 強(qiáng)，葉海林，等.地震隧洞穩(wěn)定性分析探討［J］.巖石力學(xué)與工程學(xué)報(bào)，2010，29(6):1081－1088.(ZHENG Ying-ren，XIAO Qiang，YE Hai-lin，et al.Study of Tunnel Stability Analysis with Seismic Load［J］.Chinese Journal of Rock Mechanics and Engineering，2010，29(6):1081－1088.(in Chinese))

［10］王 濤，陳曉玲，于利宏.地下洞室群圍巖穩(wěn)定的離散元計(jì)算［J］.巖土力學(xué)，2005，26(12):1936－1940.(WANG Tao，CHEN Xiao-ling，YU Li-hong.Discrete Element Calculation of Surrounding Rock Mass Stability of Underground Cavern Group［J］.Rock and Soil Mechanics，2005，26(12):1936－1940.(in Chinese))

［11］夏 祥，李俊如，李海波，等.爆破荷載作用下巖體振動(dòng)特征的數(shù)值模擬［J］.巖土力學(xué)，2005，25(1):50－56.(XIA Xiang，LI Jun-ru，LI Hai-bo，et al.UDEC Modeling of Vibration Characteristics of Jointed Rock Mass Under Explosion［J］.Rock and Soil Mechanics，2005，25(1):50－56.(in Chinese))

［12］趙 堅(jiān)，陳壽根，蔡軍剛，等.用UDEC模擬爆炸波在節(jié)理巖體中的傳播［J］.中國(guó)礦業(yè)大學(xué)學(xué)報(bào)，2002，31(2):111－115.(ZHAO Jian，CHEN Shou-gen，CAI Jun-gang，et al.Simulation of Blast Wave Propagation in Joined Rock Mass Using UDEC［J］.Rock and Soil Mechanics，2002，31(2):111－115.(in Chinese))