王金梅, 張迎賓, 趙興權(quán),2, 余鵬程, 王 潘

(1. 西南交通大學(xué)土木工程學(xué)院, 四川 成都 610031; 2. 山東建筑大學(xué)土木工程學(xué)院, 山東 濟(jì)南 250101; 3. 四川省地質(zhì)工程勘察院, 四川 成都 610072)

DDA (discontinuous deformation analysis)方法[1]是以研究非連續(xù)塊體系統(tǒng)不連續(xù)位移和變形為目的的一種數(shù)值方法,能夠很好地模擬塊體間的滑動(dòng)、張開(kāi)和閉合,自提出以來(lái),已廣泛地應(yīng)用于滑坡[2-6]、落石[7-9]、隧洞坍塌[10-11]、開(kāi)挖爆破[12-13]和傾倒破壞[14]等許多工程領(lǐng)域.

許多學(xué)者對(duì)DDA進(jìn)行了更加深入的研究.作為一種非連續(xù)介質(zhì)力學(xué)數(shù)值分析方法,DDA方法的最大特點(diǎn)是可以反映巖石塊體之間接觸面的滑移、分離與傾覆旋轉(zhuǎn)等大變形大位移問(wèn)題,因而,塊體間的接觸相關(guān)問(wèn)題成為DDA方法的核心問(wèn)題之一.針對(duì)DDA方法在接觸上存在的一些缺陷,許多學(xué)者做了進(jìn)一步的修正和改進(jìn).焦玉勇等[15]在DDA方法中加入塊體自適應(yīng)剖分算法,并引入虛實(shí)節(jié)理轉(zhuǎn)換來(lái)模擬裂紋的擴(kuò)展、貫通到碎裂、坍塌的全過(guò)程.夏才初等[16]利用表征虛節(jié)理向?qū)嵐?jié)理轉(zhuǎn)化程度的虛節(jié)理連通率對(duì)虛節(jié)理力學(xué)性質(zhì)參數(shù)的弱化規(guī)律進(jìn)行研究,基于Jennings強(qiáng)度準(zhǔn)則和最大抗拉強(qiáng)度準(zhǔn)則導(dǎo)出了節(jié)理力學(xué)性質(zhì)參數(shù)弱化函數(shù),以模擬斷續(xù)節(jié)理的開(kāi)裂.初始的DDA程序由于其內(nèi)在的接觸判斷、接觸處理方式導(dǎo)致材料粘聚力不能準(zhǔn)確地被考慮,使DDA的穩(wěn)定性分析功能被大打折扣,對(duì)于此問(wèn)題,Wang等[17]采用接觸界面相對(duì)位移量作為界面粘聚力的失效判定,提出了與接觸界面相對(duì)位移相關(guān)的剪切強(qiáng)度準(zhǔn)則;Zhang等[18]針對(duì)邊邊接觸提出了一種用邊邊接觸上的兩個(gè)角邊接觸共同判定整個(gè)邊邊接觸上的接觸狀態(tài)的處理方式,以提高剪切作用下考慮粘聚力時(shí)節(jié)理上的計(jì)算精度;馬永政等[19]在DDA中引入剛塑性模型,補(bǔ)償接觸點(diǎn)在滑面上的黏聚力作用,即令粘聚力在整個(gè)滑移過(guò)程中均保持不變,使得用DDA方法分析邊坡的安全系數(shù)明顯提高.Jiao等[20]在接觸力的法向方向引入分段力-位移關(guān)系本構(gòu),切向采用摩爾庫(kù)倫準(zhǔn)則,以模擬巖石開(kāi)裂的整個(gè)過(guò)程.為了防止應(yīng)力波在模型邊界上的反射,Jiao[21]、張秀麗[22]等在DDA程序中加入了對(duì)黏性邊界的模擬,研究了節(jié)理面對(duì)應(yīng)力波傳播的影響.Bao等[23]針對(duì)DDA中用最短路徑法搜索角角接觸時(shí)接觸邊不唯一的缺陷,提出了兩種改進(jìn)的最短路徑搜索方法:通過(guò)施加一個(gè)臨時(shí)的角角接觸彈簧去決定接觸塊體的運(yùn)動(dòng)方向;當(dāng)接觸角發(fā)生嵌入時(shí),用接觸角在上一時(shí)步的軌跡來(lái)決定進(jìn)入邊.Fan等[24]對(duì)角角接觸的預(yù)處理提出了兩個(gè)臨時(shí)彈簧法,對(duì)角角接觸向角邊接觸的轉(zhuǎn)換提出了方向角法.

目前,學(xué)者們對(duì)DDA方法在其接觸上的改進(jìn)工作主要集中在接觸的搜索、虛實(shí)節(jié)理的轉(zhuǎn)換以及節(jié)理強(qiáng)度的處理等方面,對(duì)DDA方法及其程序計(jì)算過(guò)程中其節(jié)理上的參數(shù)賦值卻少有提及.在原DDA方法中,節(jié)理參數(shù)的定義是根據(jù)分析初始時(shí)刻塊體系統(tǒng)中塊體的分布和塊體間的節(jié)理性質(zhì)將其定義到塊體邊上的,計(jì)算時(shí)按接觸算法選擇某一接觸邊上的參數(shù)進(jìn)行計(jì)算.這種節(jié)理強(qiáng)度由接觸對(duì)中單邊決定的節(jié)理參數(shù)定義,對(duì)于含有多種節(jié)理的塊體系統(tǒng),在計(jì)算過(guò)程中任意塊體間發(fā)生接觸形成節(jié)理時(shí)的節(jié)理參數(shù)選取是不適用的,會(huì)導(dǎo)致賦有不同節(jié)理參數(shù)的兩條塊體邊發(fā)生接觸時(shí)節(jié)理參數(shù)選取的非一致性.對(duì)此,本文提出一種由節(jié)理兩側(cè)接觸邊共同決定的節(jié)理參數(shù)定義方式以期對(duì)其進(jìn)行改進(jìn).

1 DDA方法概述

1.1 DDA基本原理

DDA將研究對(duì)象看成不同塊體組成的塊體系統(tǒng),塊體之間的相互作用通過(guò)接觸彈簧來(lái)實(shí)現(xiàn),荷載按增量方法分時(shí)步施加到塊體,每一時(shí)步塊體的變形和位移滿足小變形和小位移假定,塊體的大變形和大位移是由眾多時(shí)步的小變形和小位移疊加的結(jié)果.設(shè)每一塊體通體有常應(yīng)力和常應(yīng)變,每一時(shí)步內(nèi),塊體i內(nèi)部任意一點(diǎn)(x,y)的總位移可用6個(gè)位移不變量表示,即

Di=(u0,v0,r0,εx,εy,γxy)T,

(1)

式中：(u0,v0)是塊體內(nèi)特殊點(diǎn)(x0,y0)的剛體位移;r0是塊體繞轉(zhuǎn)動(dòng)中心(x0,y0)的轉(zhuǎn)動(dòng)角;(εx,εy,γxy)是該塊體的正應(yīng)變和剪應(yīng)變.

塊體內(nèi)任意點(diǎn)(x,y)的位移(u，v)等于平動(dòng)、轉(zhuǎn)動(dòng)和變形共同作用的總和,可表示為

(2)

式中:Ti為塊體i的位移轉(zhuǎn)換矩陣，

各個(gè)塊體是連接的,并靠塊體間的接觸和對(duì)單個(gè)塊體的位移約束形成一個(gè)塊體系統(tǒng).通過(guò)勢(shì)能最小原理求得整個(gè)塊體系統(tǒng)的基本方程為

(3)

式中:

k為系統(tǒng)塊體總數(shù),因每個(gè)塊體有6個(gè)自由度(u0,v0,r0,εx,εy,γxy),式(3)給出的系數(shù)矩陣中每個(gè)元素Kij(i,j=1,2,…,k)是一個(gè)6×6的子矩陣,Di、Fi是6×1子矩陣；

子矩陣Kii與塊體材料特性和幾何尺寸有關(guān)；

Kij(i≠j)由塊體i和塊體j間的接觸條件決定；

Fi是在塊體i上分配給6個(gè)變形變量的荷載.

1.2 DDA接觸理論

塊體系統(tǒng)中,任一塊體可能會(huì)與其他塊體發(fā)生接觸.DDA方法中,塊體之間存在3種接觸方式:角對(duì)角接觸、角對(duì)邊接觸、邊對(duì)邊接觸.運(yùn)算中將一對(duì)邊邊接觸處理成兩對(duì)角邊接觸,相當(dāng)于塊體之間僅有角對(duì)角、角對(duì)邊兩種接觸,如圖1所示，圖中，節(jié)點(diǎn)編號(hào)即為角的編號(hào)，同時(shí)兩個(gè)角的編號(hào)代表一條邊，如，12代表角1、角2及邊12.若采用Mohr-Coulomb強(qiáng)度準(zhǔn)則,則涉及3個(gè)節(jié)理參數(shù):摩擦角φ、粘聚力c、抗拉強(qiáng)度σt,這些參數(shù)可根據(jù)相應(yīng)的材料試驗(yàn)得到.

對(duì)于角角接觸(點(diǎn)對(duì)點(diǎn)的接觸)、單純的角邊接觸(點(diǎn)線接觸),計(jì)算中只需用到φ這一個(gè)參數(shù).對(duì)于由邊邊接觸轉(zhuǎn)化而來(lái)的角邊接觸(線接觸),根據(jù)接觸的性質(zhì),計(jì)算中除涉及到φ外,還可能用到c、σt參數(shù).由于計(jì)算中將邊邊接觸處理成兩個(gè)角邊接觸,邊邊接觸上由c提供的抗剪能力和σt提供的抗拉能力則根據(jù)接觸長(zhǎng)度均分到兩個(gè)等價(jià)的角邊接觸上.塊體的運(yùn)動(dòng)過(guò)程中,不允許塊體之間受拉和嵌入.若塊體在某個(gè)接觸上發(fā)生嵌入,則在該接觸的法向和切向施加兩個(gè)給定變形狀態(tài)的剛硬彈簧,使得迭代計(jì)算后塊體系統(tǒng)在該部位不再發(fā)生嵌入.DDA方法中每一時(shí)步塊體的接觸狀態(tài)包括:閉合、滑動(dòng)、張開(kāi)3種狀態(tài),塊體接觸狀態(tài)由式(4)決定.

(4)

式中:

Kg、Ks分別為法向和切向彈簧剛度;

dg、ds分別為法向和切向嵌入距離;

l為1/2接觸長(zhǎng)度.

若接觸狀態(tài)為閉合,則在相應(yīng)接觸部位法向和切向施加兩個(gè)彈簧;若接觸狀態(tài)為滑動(dòng),則在接觸部位施加法向彈簧,切向加摩擦力;若接觸狀態(tài)為張開(kāi),則不施加彈簧,直至系統(tǒng)中所有的接觸部位都滿足不嵌入和無(wú)張拉準(zhǔn)則,則本時(shí)步計(jì)算完畢,進(jìn)入下一時(shí)步.

(a) 角-角接觸(b) 角-邊接觸(c) 邊-邊接觸圖1 DDA塊體接觸類(lèi)型Fig.1 Contact types of blocks in DDA

2 原DDA中接觸參數(shù)定義存在的問(wèn)題

巖體中由各種成因形成的有規(guī)律的破裂面稱(chēng)為節(jié)理,天然巖體是由節(jié)理和巖塊組成.從巖體力學(xué)屬性來(lái)看,可認(rèn)為完整巖體屬連續(xù)介質(zhì)力學(xué)范疇,節(jié)理巖體因受節(jié)理切割的影響,認(rèn)為其部分屬非連續(xù)介質(zhì)力學(xué)范疇.從巖體的力學(xué)強(qiáng)度來(lái)看,巖石的強(qiáng)度與組成巖體的巖塊和節(jié)理力學(xué)性質(zhì)有很大不同,通常,節(jié)理強(qiáng)度低于巖石強(qiáng)度,節(jié)理巖體的強(qiáng)度處在完整巖石強(qiáng)度和節(jié)理強(qiáng)度之間.因此,節(jié)理的存在造成了巖體介質(zhì)的不連續(xù)性,并使巖體與完整巖石的力學(xué)特性之間有很大差異,因而節(jié)理面的性質(zhì)將大大影響整個(gè)巖體的穩(wěn)定性.

節(jié)理面亦稱(chēng)結(jié)構(gòu)面,是由兩個(gè)表面組成.節(jié)理強(qiáng)度主要受到節(jié)理面的粗糙度和起伏度、面壁強(qiáng)度、開(kāi)度與填充物等因素影響,與構(gòu)成節(jié)理的兩個(gè)表面息息相關(guān).

用DDA方法進(jìn)行數(shù)值模擬時(shí),需要根據(jù)節(jié)理分布對(duì)系統(tǒng)進(jìn)行塊體切割,對(duì)完整巖石切割后的塊體,塊體之間以虛節(jié)理表示,計(jì)算時(shí)用完整巖石的強(qiáng)度進(jìn)行處理,虛節(jié)理破壞后即成為真實(shí)節(jié)理,節(jié)理強(qiáng)度隨之降低.原DDA方法程序中在處理節(jié)理強(qiáng)度參數(shù)時(shí)默認(rèn)將節(jié)理參數(shù)賦值到節(jié)理兩側(cè)塊體的塊體邊上,沿塊體外側(cè)將節(jié)理參數(shù)編號(hào)保存到塊體邊左側(cè)的節(jié)點(diǎn)中,如圖2所示.塊體i中每條邊上的節(jié)理參數(shù)保存到了與邊左側(cè)顏色一致的節(jié)點(diǎn)信息中(包括塊體的自由邊或臨空邊,程序中也需要對(duì)其賦值節(jié)理參數(shù)).計(jì)算時(shí),選擇接觸邊(一組邊邊接觸視為兩個(gè)角邊接觸)所賦的節(jié)理參數(shù)進(jìn)行接觸計(jì)算,如圖2中的塊體j的角7與塊體i的邊56形成的角邊接觸VC7-56(VC表示角邊接觸,其下標(biāo)分表為接觸角和接觸邊),則以塊體i的邊56上所賦的節(jié)理參數(shù)作為計(jì)算標(biāo)準(zhǔn).

圖2 節(jié)理參數(shù)對(duì)應(yīng)賦值點(diǎn)示意Fig.2 The corresponding vertex of the joint parameters diagram

計(jì)算初始時(shí),根據(jù)已知的節(jié)理參數(shù)和巖石參數(shù)對(duì)DDA塊體系統(tǒng)中節(jié)理接觸兩側(cè)的塊體邊賦上相同的節(jié)理參數(shù).在塊體系統(tǒng)的運(yùn)動(dòng)過(guò)程中,隨著初始接觸的破壞、新接觸的生成,當(dāng)任意兩個(gè)塊體由于接觸而形成新的接觸時(shí),若新接觸的兩條接觸邊上所賦的節(jié)理參數(shù)不一致時(shí),其接觸節(jié)理面上的計(jì)算參數(shù)的選取則存在多種情況,并給計(jì)算帶來(lái)極大的誤差.例如:n1、n2為二維塊體系統(tǒng)中的兩個(gè)任意塊體，如圖3所示.塊體n1的邊14和塊體n2的邊58 定義的節(jié)理參數(shù)編號(hào)分別為J1、J2,當(dāng)前時(shí)步下(第s1步),兩塊體尚未發(fā)生接觸.計(jì)算中某時(shí)刻(第s2步,s2>s1)假設(shè)塊體n1的邊14與塊體n2的邊58發(fā)生邊邊接觸EC14-58,(EC表示邊邊接觸,其下標(biāo)分別為接觸的兩條邊),如圖4所示.由于DDA方法中一個(gè)邊邊接觸等價(jià)于兩個(gè)角邊接觸,接觸EC14-58將被分別處理成圖4(a)VC1-58和VC5-14、圖4(b)VC8-14和VC5-14、圖4(c)VC1-58和VC4-58兩個(gè)角邊接觸.計(jì)算時(shí),圖4(a)中VC1-58接觸上的參數(shù)取J2節(jié)理參數(shù)編號(hào)下的參數(shù)值,VC5-14接觸上的參數(shù)取J1節(jié)理參數(shù)編號(hào)下的參數(shù)值;圖4(b)中接觸均取J2節(jié)理參數(shù)編號(hào)下的參數(shù)值;圖4(c)中接觸均取J1節(jié)理參數(shù)編號(hào)下的參數(shù).圖4(a)～(c)中,14、58兩種結(jié)構(gòu)面形成的同一性質(zhì)的節(jié)理應(yīng)該有相同的節(jié)理參數(shù).在原DDA方法模擬中,由于節(jié)理兩側(cè)接觸邊參數(shù)的不同,其在不同情況下接觸計(jì)算參數(shù)可能不一致,從而使得計(jì)算結(jié)果與實(shí)際情況不相符.

(a) (b) (c) 圖3 第s1步塊體接觸狀態(tài)Fig.3 Contact state of blocks in the s1th step

(a) (b) (c) 圖4 第s2步塊體狀態(tài)(邊邊接觸)Fig.4 Contact state of blocks in the s2 th step(edge to edge contact)

3 接觸參數(shù)定義的改進(jìn)

原DDA方法中將節(jié)理強(qiáng)度參數(shù)賦值到塊體的每一條邊上,計(jì)算時(shí)只選取接觸的某一邊上的強(qiáng)度參數(shù)進(jìn)行計(jì)算,即節(jié)理強(qiáng)度由單邊決定.當(dāng)某一對(duì)邊邊接觸的兩個(gè)接觸邊上的強(qiáng)度參數(shù)不一致時(shí),同一性質(zhì)的節(jié)理在不同情況下接觸計(jì)算參數(shù)就可能不一致.造成這種錯(cuò)誤的主要原因是原DDA方法中節(jié)理強(qiáng)度由單邊參數(shù)決定的節(jié)理參數(shù)賦值方式與實(shí)際物理情況是不相符的,實(shí)際物理情況下,節(jié)理強(qiáng)度是由組成節(jié)理的兩個(gè)結(jié)構(gòu)面共同決定的.對(duì)此,本文提出一種改進(jìn)方法,即在一個(gè)分析初始時(shí)具有m種實(shí)節(jié)理的塊體系統(tǒng)中,本文額外定義(m2-m)/2種節(jié)理參數(shù),來(lái)處理計(jì)算過(guò)程中當(dāng)某一對(duì)新形成的邊邊接觸的兩個(gè)接觸邊上的強(qiáng)度參數(shù)不一致時(shí)節(jié)理強(qiáng)度參數(shù)取值的問(wèn)題.在石根華老師編寫(xiě)的原DDA程序[1]基礎(chǔ)上改進(jìn)后的程序計(jì)算流程圖如圖5所示.

圖5 改進(jìn)后的DDA計(jì)算流程Fig.5 Flow chart of the improved DDA

圖5中：Jvi為虛節(jié)理的節(jié)理參數(shù)編號(hào)；φvi、cvi、σvti為虛節(jié)理相應(yīng)的節(jié)理力學(xué)參數(shù)；cn為當(dāng)前計(jì)算時(shí)步內(nèi)搜索到的接觸數(shù)；ng為當(dāng)前計(jì)算時(shí)步.

設(shè)節(jié)理參數(shù)編號(hào)為Ja、Jb(a,b=1,2,..,m,其相應(yīng)的節(jié)理參數(shù)分別為(φa、ca、σta)、(φb、cb、σtb))的兩條接觸邊形成的節(jié)理，其節(jié)理參數(shù)編號(hào)定義為Jab或Jba(相應(yīng)的節(jié)理參數(shù)為(φab、cab、σtab),Jab=Jba),節(jié)理參數(shù)編號(hào)均為Ja的兩條接觸邊形成的節(jié)理，其節(jié)理參數(shù)編號(hào)仍為Ja,如表1所示.(m2-m)/2種節(jié)理參數(shù)可通過(guò)試驗(yàn)獲得.通過(guò)節(jié)理參數(shù)編號(hào)即可取得該節(jié)理的參數(shù),若采用摩爾庫(kù)倫強(qiáng)度準(zhǔn)則和最大抗拉強(qiáng)度準(zhǔn)則,則節(jié)理強(qiáng)度參數(shù)包括:摩擦角φab、粘聚力cab、抗拉強(qiáng)度σtab.

表1 節(jié)理接觸參數(shù)定義Tab.1 The definition of joint parameters

4 算例驗(yàn)證

以塊體的滑動(dòng)為例,如圖6所示.塊體n2、n3固定,初始時(shí),塊體n1在塊體n3上,并以一定速度v從塊體n3上滑動(dòng)到塊體n2上.圖中,1、2、5、6、9、10為塊體的部分節(jié)點(diǎn)編號(hào).在原DDA方法中,由于初始時(shí)塊體n1只與塊體n3接觸,可預(yù)先給定塊體n1、n3間的節(jié)理強(qiáng)度參數(shù),并賦值到塊體n1、n3接觸邊相應(yīng)的節(jié)點(diǎn)上;而塊體n1、n2尚未發(fā)生接觸.計(jì)算過(guò)程中,若其發(fā)生接觸,則其接觸上的節(jié)理強(qiáng)度參數(shù)只能根據(jù)其接觸面上所賦的參數(shù)按接觸算法進(jìn)行給定.設(shè)塊體n1剛好完全滑到塊體n2上時(shí)的速度為v0(v0=5 m/s),如圖7所示.

圖7中分別給出了計(jì)算的兩種工況:工況1中,塊體n1、n2的材料編號(hào)分別為R1、R2,邊12、56上的節(jié)理參數(shù)編號(hào)分別為J1、J2;工況2中,塊體n1、n2的材料編號(hào)分別為R2、R1,邊12、56上的節(jié)理參數(shù)編號(hào)分別為J2、J1;其力學(xué)參數(shù)和界面力學(xué)參數(shù)如表2～3所示.兩種工況下,當(dāng)塊體n1剛好完全滑倒塊體n2上后,塊體n1的速度v、位移隨時(shí)間的變化曲線如圖8所示.

兩種工況下,形成塊體n1、n2間的接觸節(jié)理的兩個(gè)結(jié)構(gòu)面性質(zhì)相同,應(yīng)該有相同的節(jié)理參數(shù),塊體n1應(yīng)該有相同的速度和位移.然而,從圖8可見(jiàn),兩種情況下,塊體n1以相同的速度從塊體n2上開(kāi)始滑動(dòng),經(jīng)過(guò)同一性質(zhì)的節(jié)理,其速度、位移并不相同.這是因?yàn)樵谠璂DA方法中,當(dāng)邊12與邊56發(fā)生邊邊接觸時(shí),程序自動(dòng)將其處理為兩個(gè)角邊接觸VC1-65、VC2-65,計(jì)算中選擇角邊接觸中接觸邊上的參數(shù)進(jìn)行計(jì)算(這里兩個(gè)角邊接觸的接觸邊均為56),因而當(dāng)兩個(gè)塊體交換材料時(shí),原DDA方法均選擇接觸邊56上的參數(shù)進(jìn)行計(jì)算,這就導(dǎo)致同一性質(zhì)的節(jié)理間的接觸節(jié)理參數(shù)不唯一,從而引起計(jì)算結(jié)果的不一致.

圖6 塊體滑動(dòng)算例Fig.6 Example for block sliding

(a) 工況1

(b) 工況2圖7 兩種塊體滑動(dòng)Fig.7 Two working condition of the block sliding

表2 塊體力學(xué)參數(shù)Tab.2 The mechanical parameters of blocks

表3 界面力學(xué)參數(shù)Tab.3 The mechanical parameters of joints

實(shí)際上節(jié)理性質(zhì)由其兩個(gè)表面共同決定,即使塊體1、2交換接觸邊12、56上的節(jié)理參數(shù)和接觸EC12-56仍為同一性質(zhì)的節(jié)理,其節(jié)理參數(shù)應(yīng)該不變.根據(jù)本文提出的節(jié)理參數(shù)定義方法,設(shè)原DDA方法中由節(jié)理參數(shù)編號(hào)為J1、J2的兩條接觸邊構(gòu)成的接觸其節(jié)理參數(shù)編號(hào)為J12,由本文改進(jìn)的DDA程序算得兩種工況下塊體n1的速度、位移變化如圖9所示.

由圖9可見(jiàn),塊體n1的速度、位移變化一致,說(shuō)明了本文的節(jié)理參數(shù)定義方式的有效性,該定義方式能在計(jì)算中反映更加貼近實(shí)際的節(jié)理參數(shù).

(a) 位移(b) 速度圖8 原DDA計(jì)算的塊體n1速度、位移曲線Fig.8 The velocity and displacement time history curves of block n1 by original DDA

(a) 位移(b) 速度圖9 改進(jìn)DDA計(jì)算的塊體n1速度、位移曲線Fig.9 The velocity and displacement time history curves of block n1 by improved DDA

5 結(jié)束語(yǔ)

塊體間的接觸處理相關(guān)問(wèn)題是DDA方法計(jì)算中的核心問(wèn)題,在原DDA方法中,程序在處理節(jié)理強(qiáng)度參數(shù)時(shí)默認(rèn)將節(jié)理參數(shù)賦值到節(jié)理兩側(cè)塊體的塊體邊上,計(jì)算時(shí)選取某一側(cè)接觸邊上的節(jié)理參數(shù)進(jìn)行計(jì)算.這種處理方式,對(duì)于塊體系統(tǒng)運(yùn)動(dòng)過(guò)程中不同參數(shù)性質(zhì)的塊體間接觸時(shí)接觸節(jié)理參數(shù)的選取是不適宜的,將造成塊體間同一節(jié)理上的節(jié)理參數(shù)取值的不一致.

針對(duì)上述問(wèn)題,本文從節(jié)理的實(shí)際物理意義出發(fā),提出了一種由節(jié)理兩側(cè)接觸邊共同決定的節(jié)理參數(shù)定義方式,編寫(xiě)了相應(yīng)的程序整合到原DDA程序中,并給出了具體算例.本文的定義方式更符合物理實(shí)際的節(jié)理參數(shù)定義,可以在計(jì)算中反映更加貼近實(shí)際的節(jié)理參數(shù);克服了原DDA方法中節(jié)理參數(shù)選取時(shí)出現(xiàn)的非一致情況;程序簡(jiǎn)潔、思路清晰,可以很好的整合到原程序中去,能夠方便地實(shí)現(xiàn)了分析初期和分析過(guò)程中塊體間節(jié)理參數(shù)的賦值.

[1] SHI Genhua. Discontinuous deformation analysis: a new numerical model for the statics and dynamics of block system[D]. Berkeley: University of California, 1988.

[2] SITAR N, MACLAUGHLIN M M, DOOLIN D M. Influence of kinematics on landslide mobility and failure mode[J]. Journal of Geotechnical & Geoenvir onmental Engineering, 2005, 131(6): 716-728.

[3] WU Jianhong. Applying discontinuous deformation analysis to assess the constrained area of the unstable Chiu-fen-erh-shan landslide slope[J]. International Journal for Numerical & Analytic Methods in Geomechanics, 2010, 31(5): 649-666.

[4] KVELDSVIK V, EINSTEIN H H, NILSEN B, et al. Numerical analysis of the 650 000 m2?knes rock slope based on measured displacements and geotechnical data[J]. Rock Mechanics and Rock Engineering, 2009, 42(5): 689-728.

[5] WU Jianhong, CHEN Chunhwa. Application of DDA to simulate characteristics of the Tsaoling landslide[J]. Computers & Geotechnics, 2011, 38(5): 741-750.

[6] ZHANG Yingbin, CHEN Guangqi, ZHENG Lu, et al. Effects of near-fault seismic loadings on run-out of large-scale landslide: a case study[J]. Engineering Geology, 2013, 166(8): 216-236.

[7] CHEN Guangqi. Numerical modelling of rockfall using extended DDA[J]. Chinese Journal of Rock Mechanics and Engineering, 2003, 22(6): 926-931.

[8] MATSUYAMA H, NISHIYAMA S, OHNISHI Y. Practical studies on rockfall simulation by DDA[J]. Journal of Rock Mechanics and Geotechnical Engineering, 2011, 3(1): 57-63.

[9] CHEN Guangqi, ZHENG Lu, ZHANG Yingbin, et al. Numerical simulation in rockfall analysis: a close comparison of 2-D and 3-D DDA[J]. Rock Mechanics and Rock Engineering, 2013, 46(3): 527-41.

[10] WU Jianhong, OHNISHI Y, NISHIYAMA S. Simulation of the mechanical behavior of inclined jointed rock masses during tunnel construction using discontinuous deformation analysis (DDA)[J]. International Journal of Rock Mechanics and Mining Sciences, 2004, 41(5): 731-743.

[11] TSESARSKY M, HATZOR Y. Tunnel roof deflection in blocky rock masses as a function of joint spacing and friction-A parametric study using discontinuous deformation analysis (DDA)[J]. Tunnelling & Underground Space Technology, 2006, 21(1): 29-45.

[12] NING Youjun, YANG Jun, AN Xinmei, et al. Simulation of blast induced crater in jointed rock mass by discontinuous deformation analysis method[J]. Frontiers of Structural and Civil Engineering, 2010, 4(2): 223-232.

[13] NING Youjun, YANG Jun, AN Xinmei, et al. Modelling rock fracturing and blast-induced rock mass failure via advanced discretisation within the discontinuous deformation analysis framework[J]. Computers & Geotechnics, 2011, 38(1): 40-49.

[14] 孫東亞,彭一江,王興珍. DDA數(shù)值方法在巖質(zhì)邊坡傾倒破壞分析中的應(yīng)用[J]. 巖石力學(xué)與工程學(xué)報(bào), 2002,121(1): 39-42.

SUN Dongya, PENG Yijiang, WANG Xingzhen. Application of DDA method in stability analysis of topple rock slope[J]. Chinese Journal of Rock Mechanics and Engineering, 2002, 121(1): 3942.

[15] 焦玉勇,張秀麗,劉泉聲,等. 用非連續(xù)變形分析方法模擬巖石裂紋擴(kuò)展[J]. 巖石力學(xué)與工程學(xué)報(bào),2007,26(4): 682-691.

JIAO Yuyong, ZHANG Xiuli, LIU Quansheng, et al. Simulation of rock crack propagation using discontinuous deformation analysis method[J]. Chinese Journal of Rock Mechanics and Engineering, 2007, 26(4): 682-691.

[16] 夏才初,許崇幫. 非連續(xù)變形分析(DDA)中斷續(xù)節(jié)理擴(kuò)展的模擬方法研究和試驗(yàn)驗(yàn)證[J]. 巖石力學(xué)與工程學(xué)報(bào),2010,29(10): 2027-2033.

XIA Caichu, XU Chongbang. Atudy of fracturing algorithm of intermitteng joint by DDA and experimental validation[J]. Chinese Journal of Rock Mechanics and Engineering, 2010, 29(10): 2027-2033.

[17] WANG Lizhong, JIANG Hongyi, YANG Zhongxuan, et al. Development of discontinuous deformation analysis with displacement-dependent interface shear strength[J]. Computers and Geotechnics. 2013, 47(1): 91-101.

[18] ZHANG Yingbin, XU Qiang, CHEN Guangqi, et al. Extension of discontinuous deformation analysis and application in cohesive-frictional slope analysis[J]. International Journal of Rock Mechanics & Mining Sciences, 2014, 70(9): 533-545.

[19] 馬永政,鄭宏,朱合華,等. DDA法計(jì)算邊坡安全系數(shù)的黏聚力影響分析[J]. 巖土工程學(xué)報(bào),2009,31(7): 1088-1093.

MA Yongzheng, ZHEN Hong, ZHU Hehua, et al. Effect of cohesion on evaluating slope stability factor of safety by DDA method[J]. Chinese Journal of Geotechnical Engineering, 2009, 31(7): 1088-1093.

[20] JIAO Yuyong, ZHANG Xiuli, ZHAO Jian. Two-dimensional DDA contact constitutive model for simulating rock fragmentation[J]. Journal of Engineering Mechanics, 2011, 138(2): 199-209.

[21] JIAO Yuyong, ZHANG Xiuli, ZHAO Jian. Viscous boundary of DDA for modeling stress wave propagation in jointed rock[J]. International Journal of Rock Mechanics and Mining Sciences, 2007, 44(7): 1070-1076.

[22] 張秀麗,焦玉勇,劉泉聲,等. 節(jié)理對(duì)爆炸波傳播影響的數(shù)值研究[J]. 巖土力學(xué),2008,29(3): 717-721.

ZHANG Xiuli, JIAO Yuyong, LIU Quansheng, et al. Numerical study on effect of joints on blasting wave progation in rock mass[J]. Rock and Soil Mechanics, 2008, 29(3): 717-721.

[23] BAO Huirong, ZHAO Zhiye. An alternative scheme for the corner-corner contact in the two-dimensional Discontinuous Deformation Analysis[J]. Advances in Engineering Software, 2010, 41(2): 206-212.

[24] FAN Huo, HE Siming. An angle-based method dealing with vertex-vertex contact in the two-dimensional discontinuous deformation analysis[J]. Rock Mechanics and Rock Engineering, 2015, 48(5): 2031-2043.