udec數(shù)值方法

1、一、概述一、概述二、基本原理二、基本原理四、四、UDEC（Universal Discrete Element Code）應(yīng)用簡(jiǎn)介）應(yīng)用簡(jiǎn)介三、運(yùn)動(dòng)方程的解法三、運(yùn)動(dòng)方程的解法推薦參考書推薦參考書：王泳嘉：王泳嘉.離散單元法及其在巖土力學(xué)中的應(yīng)用離散單元法及其在巖土力學(xué)中的應(yīng)用. 東北工學(xué)院出版社，東北工學(xué)院出版社，1991第六章第六章 離散單元法離散單元法 離散元法是專門用來解決不連續(xù)介質(zhì)問題的數(shù)值模擬方離散元法是專門用來解決不連續(xù)介質(zhì)問題的數(shù)值模擬方法。該方法把節(jié)理巖體視為由離散的巖塊和巖塊間的節(jié)理面法。該方法把節(jié)理巖體視為由離散的巖塊和巖塊間的節(jié)理面所組成所組成,允許巖塊平移、轉(zhuǎn)動(dòng)和變形

2、允許巖塊平移、轉(zhuǎn)動(dòng)和變形,而節(jié)理面可被壓縮、分離而節(jié)理面可被壓縮、分離或滑動(dòng)。因此或滑動(dòng)。因此,巖體被看作一種不連續(xù)的離散介質(zhì)。其內(nèi)部可巖體被看作一種不連續(xù)的離散介質(zhì)。其內(nèi)部可存在大位移、旋轉(zhuǎn)和滑動(dòng)乃至塊體的分離存在大位移、旋轉(zhuǎn)和滑動(dòng)乃至塊體的分離,從而可以較真實(shí)地從而可以較真實(shí)地模擬節(jié)理巖體中的非線性大變形特征。模擬節(jié)理巖體中的非線性大變形特征。一、概述一、概述 離散元法的一般求解過程為離散元法的一般求解過程為:將求解空間離散為離散元單將求解空間離散為離散元單元陣元陣,并根據(jù)實(shí)際問題用合理的連接元件將相鄰兩單元連接起并根據(jù)實(shí)際問題用合理的連接元件將相鄰兩單元連接起來來;單元間相對(duì)位移是基本變

3、量單元間相對(duì)位移是基本變量,由力與相對(duì)位移的關(guān)系可得到由力與相對(duì)位移的關(guān)系可得到兩單元間法向和切向的作用力兩單元間法向和切向的作用力;對(duì)單元在各個(gè)方向上與其它單對(duì)單元在各個(gè)方向上與其它單元間的作用力以及其它物理場(chǎng)對(duì)單元作用所引起的外力求合力元間的作用力以及其它物理場(chǎng)對(duì)單元作用所引起的外力求合力和合力矩和合力矩,根據(jù)牛頓運(yùn)動(dòng)第二定律可以求得單元的加速度根據(jù)牛頓運(yùn)動(dòng)第二定律可以求得單元的加速度;對(duì)其對(duì)其進(jìn)行時(shí)間積分進(jìn)行時(shí)間積分,進(jìn)而得到單元的速度和位移。進(jìn)而得到單元的速度和位移。 離散單元法所用的求解方法有靜態(tài)松弛法和動(dòng)態(tài)松弛法。離散單元法所用的求解方法有靜態(tài)松弛法和動(dòng)態(tài)松弛法。產(chǎn)生與發(fā)展：產(chǎn)生與

4、發(fā)展：該方法是該方法是20世紀(jì)世紀(jì)70年代發(fā)展起來的，目前已在數(shù)年代發(fā)展起來的，目前已在數(shù)值模擬理論與工程應(yīng)用方面取得了很大進(jìn)展。該方法的基本特值模擬理論與工程應(yīng)用方面取得了很大進(jìn)展。該方法的基本特征在于允許各離散塊體發(fā)生平動(dòng)和轉(zhuǎn)動(dòng)，甚至發(fā)生分離，征在于允許各離散塊體發(fā)生平動(dòng)和轉(zhuǎn)動(dòng)，甚至發(fā)生分離，彌補(bǔ)彌補(bǔ)了有限元法或邊界元法的介質(zhì)連續(xù)和小變形的限制，了有限元法或邊界元法的介質(zhì)連續(xù)和小變形的限制，因而特別因而特別適合塊狀及裂隙介質(zhì)的大變形及破壞問題的分析。適合塊狀及裂隙介質(zhì)的大變形及破壞問題的分析。 該方法首先由該方法首先由P.A.Cundall首先提出，于首先提出，于20世紀(jì)世紀(jì)80年代中期年

5、代中期由王泳嘉等人引入到我國。目前，離散元法發(fā)展迅速，美國有由王泳嘉等人引入到我國。目前，離散元法發(fā)展迅速，美國有二維和三維的離散元程序二維和三維的離散元程序UDEC和和3DEC。我國有。我國有2D-BLOCK和和3D-BLOCK。應(yīng)用領(lǐng)域應(yīng)用領(lǐng)域：邊坡、巷道與采場(chǎng)、地下開采、地震、爆炸、核廢：邊坡、巷道與采場(chǎng)、地下開采、地震、爆炸、核廢料儲(chǔ)存、散體介質(zhì)運(yùn)動(dòng)、斷裂、地下水滲流、熱傳導(dǎo)等。料儲(chǔ)存、散體介質(zhì)運(yùn)動(dòng)、斷裂、地下水滲流、熱傳導(dǎo)等。二、基本原理二、基本原理（一）塊體接觸模型（一）塊體接觸模型（二）塊體接觸本構(gòu)關(guān)系（二）塊體接觸本構(gòu)關(guān)系 塊體接觸的本構(gòu)關(guān)系是指塊體間接觸的力與位移關(guān)系。最塊體

6、接觸的本構(gòu)關(guān)系是指塊體間接觸的力與位移關(guān)系。最簡(jiǎn)單的是邊簡(jiǎn)單的是邊-角接觸的庫侖角接觸的庫侖-莫爾關(guān)系，但最符合實(shí)際的是各種莫爾關(guān)系，但最符合實(shí)際的是各種節(jié)理模型。這里僅討論邊節(jié)理模型。這里僅討論邊-角接觸的庫侖角接觸的庫侖-莫爾關(guān)系。莫爾關(guān)系。FnFnun(a)ut(b)FtFt時(shí)刻t的位置時(shí)刻t+t的位置 塊體之間的相互作用，可以假想為通過邊塊體之間的相互作用，可以假想為通過邊-角間的角間的“疊合疊合”而發(fā)生，即一個(gè)塊體的邊與另一個(gè)塊體的角之間的疊合。作用而發(fā)生，即一個(gè)塊體的邊與另一個(gè)塊體的角之間的疊合。作用力的大小可以用力的大小可以用“疊合量疊合量”來表示。來表示。FnFnun(a)ut

7、(b)FtFt時(shí)刻t的位置時(shí)刻t+t的位置 設(shè)塊體間的相互作用力與其相對(duì)位移成正比，則在塊體接設(shè)塊體間的相互作用力與其相對(duì)位移成正比，則在塊體接觸處，由塊體相對(duì)法向疊合量觸處，由塊體相對(duì)法向疊合量un和剪切量和剪切量ut引起的力的法引起的力的法向增量向增量Fn和切向增量和切向增量Ft分別為：分別為：tttnnnKFKF)()()()()()(tFtFttFtFtFttFtttnnntt時(shí)刻時(shí)刻增量理論：應(yīng)變?cè)隽坷碚摚簯?yīng)變?cè)隽颗c應(yīng)力增量增量與應(yīng)力增量一一對(duì)應(yīng)關(guān)系一一對(duì)應(yīng)關(guān)系0nFjjnttcFFFtanmaxjjc對(duì)于塊體間不允許出現(xiàn)拉力，故對(duì)于塊體間不允許出現(xiàn)拉力，故 對(duì)于剪切力，其穩(wěn)定狀態(tài)有

8、庫侖對(duì)于剪切力，其穩(wěn)定狀態(tài)有庫侖-莫爾定量：莫爾定量： 式中式中 、為接觸處塊體的內(nèi)摩擦角和粘聚力。為接觸處塊體的內(nèi)摩擦角和粘聚力。 當(dāng)當(dāng)Ft趨近于趨近于 時(shí)，滑動(dòng)即趨于發(fā)生。時(shí)，滑動(dòng)即趨于發(fā)生。 maxtF 塊體間的相互作用還包括阻尼力，在接觸點(diǎn)其法向和切向分塊體間的相互作用還包括阻尼力，在接觸點(diǎn)其法向和切向分量量Dn、Dt與接觸點(diǎn)位移增量成正比，即：與接觸點(diǎn)位移增量成正比，即： tttnnnCDCDnCtC 式中式中 、為接觸點(diǎn)法向和切向阻尼系數(shù)，與剛度為接觸點(diǎn)法向和切向阻尼系數(shù)，與剛度系數(shù)系數(shù) 有關(guān)。可由阻尼參數(shù)有關(guān)?？捎勺枘釁?shù)與彈性剛度與彈性剛度的乘積得到：的乘積得到：ttnnKCK

9、C本構(gòu)模型本構(gòu)模型與速度無關(guān)的與速度無關(guān)的接觸模型接觸模型與速度有關(guān)的與速度有關(guān)的接觸模型接觸模型角角- -邊接觸模型邊接觸模型邊邊- -邊接觸模型邊接觸模型單狀態(tài)量摩擦模型單狀態(tài)量摩擦模型雙狀態(tài)量摩擦模型雙狀態(tài)量摩擦模型法向力：法向力：假定塊體之間的法向力假定塊體之間的法向力Fn正比于它們之間沿法向正比于它們之間沿法向“疊合疊合”n的大小，即：的大小，即：Fn=Kn n, Kn為法向剛度系數(shù)。為法向剛度系數(shù)。 “疊合疊合”量量n是在數(shù)值計(jì)算時(shí)人為假定的一個(gè)量，其值的選是在數(shù)值計(jì)算時(shí)人為假定的一個(gè)量，其值的選取與計(jì)算精度要求等因素有關(guān)。取與計(jì)算精度要求等因素有關(guān)。剪切力：剪切力：由于塊體所受的

10、剪切力與塊體運(yùn)動(dòng)和加載的歷史或由于塊體所受的剪切力與塊體運(yùn)動(dòng)和加載的歷史或路徑有關(guān)，所以對(duì)于剪切力要用增量路徑有關(guān)，所以對(duì)于剪切力要用增量F Ft t來表示。設(shè)兩塊體來表示。設(shè)兩塊體之間的相對(duì)位移為之間的相對(duì)位移為t t，則：，則：Ft=Ktt, , Kt為切向剛度系數(shù)，為切向剛度系數(shù)， t為兩塊體之間的相對(duì)位移。為兩塊體之間的相對(duì)位移。破壞條件：破壞條件：法向力和切向力所表示的力與位移關(guān)系法向力和切向力所表示的力與位移關(guān)系為彈性，但在某些情況下彈性關(guān)系是不成立的，需為彈性，但在某些情況下彈性關(guān)系是不成立的，需要考慮破壞條件。如當(dāng)巖塊受到張力分離時(shí)，作用要考慮破壞條件。如當(dāng)巖塊受到張力分離時(shí)，

11、作用在巖塊表面上的法向力和剪切力隨即消失。對(duì)于塑在巖塊表面上的法向力和剪切力隨即消失。對(duì)于塑性剪切破壞的情況，需要在每次迭代時(shí)檢查剪切力性剪切破壞的情況，需要在每次迭代時(shí)檢查剪切力是否超過是否超過c+Fntan，其中，其中，c為粘結(jié)力，為粘結(jié)力，為內(nèi)摩擦為內(nèi)摩擦角。如果剪切力超過該值，此時(shí)剪切力就取為該值。角。如果剪切力超過該值，此時(shí)剪切力就取為該值。角角- -邊接觸模型：邊接觸模型：該模型認(rèn)為塊體之間不存在拉力，且當(dāng)切該模型認(rèn)為塊體之間不存在拉力，且當(dāng)切向力向力Ft達(dá)到某一最大值時(shí)，就會(huì)發(fā)生塑性剪切滑移，并由達(dá)到某一最大值時(shí)，就會(huì)發(fā)生塑性剪切滑移，并由下式確定：下式確定： maxtantnt

12、FFF 該模型的力與位移關(guān)系分別如下圖所示：該模型的力與位移關(guān)系分別如下圖所示：（a a）法向力與法向位移）法向力與法向位移（b b）切向力與切向位移）切向力與切向位移剛度系數(shù)的確定：剛度系數(shù)的確定：對(duì)于如圖所示的兩個(gè)接觸塊體，其對(duì)于如圖所示的兩個(gè)接觸塊體，其長度和寬度分別為長度和寬度分別為a、b，彈性常數(shù)為，彈性常數(shù)為E、?？傻闷浞ā？傻闷浞ㄏ騽偠认禂?shù)為：向剛度系數(shù)為： bEakabESKEnnnn2S12個(gè)單位，則設(shè)塊體厚度為塊體接觸計(jì)算模型塊體接觸計(jì)算模型切向剛度和法向剛度之間的關(guān)系為：切向剛度和法向剛度之間的關(guān)系為：為剪切模量所以，又Gbakkktnt2G)2(1EG12塊體接觸計(jì)算模

13、型塊體接觸計(jì)算模型（三）運(yùn)動(dòng)方程（三）運(yùn)動(dòng)方程- -牛頓第二運(yùn)動(dòng)定律牛頓第二運(yùn)動(dòng)定律 假設(shè)某一塊體與周圍假設(shè)某一塊體與周圍n n個(gè)塊體接觸，則其受個(gè)塊體接觸，則其受到到n n個(gè)力作用，將各個(gè)力在個(gè)力作用，將各個(gè)力在x x和和y y方向上分解，則方向上分解，則在兩個(gè)方向上的合力與合力矩分別為：在兩個(gè)方向上的合力與合力矩分別為：nixixff1niyiyff1niixiiyiyyfxxfM100塊體集合及作用于個(gè)別塊體上的力塊體集合及作用于個(gè)別塊體上的力 假設(shè)塊體的質(zhì)量為假設(shè)塊體的質(zhì)量為m，轉(zhuǎn)動(dòng)慣量為，轉(zhuǎn)動(dòng)慣量為I，考慮重力，考慮重力，則塊體質(zhì)心的運(yùn)動(dòng)方程為：則塊體質(zhì)心的運(yùn)動(dòng)方程為：mfuxx m

14、fuyy IM niiniiyiyyyniixixxxMcImgDFucumDFucum1.1.1.)()(（5-72）三、運(yùn)動(dòng)方程的解法三、運(yùn)動(dòng)方程的解法采用計(jì)算機(jī)實(shí)施數(shù)值計(jì)算的過程中通常會(huì)涉及采用計(jì)算機(jī)實(shí)施數(shù)值計(jì)算的過程中通常會(huì)涉及到動(dòng)態(tài)松弛法，所謂動(dòng)態(tài)松弛法是把非線性靜力學(xué)到動(dòng)態(tài)松弛法，所謂動(dòng)態(tài)松弛法是把非線性靜力學(xué)問題轉(zhuǎn)化為動(dòng)力學(xué)問題求解的一種數(shù)值方法。該方問題轉(zhuǎn)化為動(dòng)力學(xué)問題求解的一種數(shù)值方法。該方法的實(shí)質(zhì)是對(duì)臨界阻尼的振動(dòng)方程進(jìn)行逐步積分。法的實(shí)質(zhì)是對(duì)臨界阻尼的振動(dòng)方程進(jìn)行逐步積分。對(duì)于這種帶有阻尼項(xiàng)的動(dòng)態(tài)平衡方程，利用有限差對(duì)于這種帶有阻尼項(xiàng)的動(dòng)態(tài)平衡方程，利用有限差分法按時(shí)步在

15、計(jì)算機(jī)上迭代求解，分法按時(shí)步在計(jì)算機(jī)上迭代求解，即所謂的動(dòng)態(tài)松即所謂的動(dòng)態(tài)松弛法弛法。由于被求解的方程是時(shí)間的線性函數(shù)，整個(gè)由于被求解的方程是時(shí)間的線性函數(shù)，整個(gè)計(jì)算過程只需要直接代換，因此，對(duì)于非線性問題計(jì)算過程只需要直接代換，因此，對(duì)于非線性問題也能加以考慮，也能加以考慮，這是動(dòng)態(tài)松弛法的最大優(yōu)點(diǎn)這是動(dòng)態(tài)松弛法的最大優(yōu)點(diǎn)。 用動(dòng)態(tài)松弛法時(shí)，計(jì)算循環(huán)是以時(shí)步用動(dòng)態(tài)松弛法時(shí)，計(jì)算循環(huán)是以時(shí)步t t向前差向前差分進(jìn)行的。由于時(shí)步很小，每個(gè)單元在一個(gè)時(shí)步內(nèi)分進(jìn)行的。由于時(shí)步很小，每個(gè)單元在一個(gè)時(shí)步內(nèi)只能以很小的位移與其相鄰單元作用，而與其他單只能以很小的位移與其相鄰單元作用，而與其他單元無關(guān)，所以

16、在一個(gè)時(shí)步內(nèi)只能傳遞一個(gè)單元。元無關(guān)，所以在一個(gè)時(shí)步內(nèi)只能傳遞一個(gè)單元。計(jì)算循環(huán)示意圖：計(jì)算循環(huán)示意圖：計(jì)算循環(huán)示意圖計(jì)算循環(huán)示意圖力力- -位移的關(guān)系：位移的關(guān)系：K力力F F運(yùn)動(dòng)定律運(yùn)動(dòng)定律F=maF=ma位移位移u u力的邊界條件力的邊界條件位移邊界條件位移邊界條件離散元法采用動(dòng)態(tài)松弛法求解，其基本運(yùn)動(dòng)方程為：離散元法采用動(dòng)態(tài)松弛法求解，其基本運(yùn)動(dòng)方程為：該方法是采用顯式中心差分的動(dòng)態(tài)松弛法進(jìn)行求解，缺該方法是采用顯式中心差分的動(dòng)態(tài)松弛法進(jìn)行求解，缺點(diǎn)是計(jì)算時(shí)步要很小，且需要合理確定阻尼系數(shù)等。點(diǎn)是計(jì)算時(shí)步要很小，且需要合理確定阻尼系數(shù)等。該算法既可以由已知力求位移，也可以由已知位移求力

17、。該算法既可以由已知力求位移，也可以由已知位移求力。其總的算法為循環(huán)交錯(cuò)求解，首先根據(jù)邊界條件已知位移其總的算法為循環(huán)交錯(cuò)求解，首先根據(jù)邊界條件已知位移u(t)u(t)，則可以求出，則可以求出F(t)F(t)，然后可確定，然后可確定 , ,再加上前一輪再加上前一輪的的 ，對(duì)其積分可得，對(duì)其積分可得u(t+u(t+t/2),t/2),如此往復(fù)循環(huán)，如此往復(fù)循環(huán)，直到位移和力收斂，達(dá)到平衡狀態(tài)，或者達(dá)到某一循環(huán)數(shù)為直到位移和力收斂，達(dá)到平衡狀態(tài)，或者達(dá)到某一循環(huán)數(shù)為止。止。 tftkutuctum tu 2ttu計(jì)算循環(huán)示意圖：計(jì)算循環(huán)示意圖：離散元法交錯(cuò)循環(huán)求解特性圖離散元法交錯(cuò)循環(huán)求解特性圖力

18、力位移位移加速度加速度速度速度UDEC1.1. 非連續(xù)介質(zhì)被作為是凸塊或凹?jí)K的集合體，塊體可以是剛體非連續(xù)介質(zhì)被作為是凸塊或凹?jí)K的集合體，塊體可以是剛體或變形體?；蜃冃误w。2.2. 不連續(xù)面被看作是塊體之間的邊界條件。不連續(xù)面被看作是塊體之間的邊界條件。3.3. 沿著不連續(xù)面的運(yùn)動(dòng)是由切向和法向的線性及非線性力沿著不連續(xù)面的運(yùn)動(dòng)是由切向和法向的線性及非線性力- -位位移關(guān)系控制的。移關(guān)系控制的。4.4. 地質(zhì)體或類似材料可以由很多內(nèi)置的塊體和節(jié)理本構(gòu)模型來地質(zhì)體或類似材料可以由很多內(nèi)置的塊體和節(jié)理本構(gòu)模型來描述；還有可供開發(fā)的自定義模型。描述；還有可供開發(fā)的自定義模型。是一種非常適合于用來模擬

19、具有非線性力學(xué)行是一種非常適合于用來模擬具有非線性力學(xué)行為的非連續(xù)材料（包含多個(gè)相互切割的非連續(xù)為的非連續(xù)材料（包含多個(gè)相互切割的非連續(xù)面）。尤其是，它具有如下特征：面）。尤其是，它具有如下特征：四、四、UDEC（Universal Discrete Element Code）應(yīng)用簡(jiǎn)介）應(yīng)用簡(jiǎn)介5. 5. 有平面有平面- -應(yīng)變、平面應(yīng)變、平面- -應(yīng)力及軸對(duì)稱三種問題模型。應(yīng)力及軸對(duì)稱三種問題模型。6. 6. 具有用來描述巖石具有用來描述巖石- -結(jié)構(gòu)相互作用的結(jié)構(gòu)單元模型，如錨結(jié)構(gòu)相互作用的結(jié)構(gòu)單元模型，如錨桿、樁、梁、襯砌，噴射混凝土等。桿、樁、梁、襯砌，噴射混凝土等。7. 7. 能夠進(jìn)行靜力和動(dòng)力問題分析能夠進(jìn)行靜力和動(dòng)力問題分析8. 8. 沿著節(jié)理的暫態(tài)和穩(wěn)態(tài)流沿著節(jié)理的暫態(tài)和穩(wěn)態(tài)流. .9. 9.