基于dijcke算法的邊坡穩(wěn)定性分析

基于dijcke算法的邊坡穩(wěn)定性分析

0基于dijpstor的極限平衡分析方法的分析確定坡度的最危險滑動表面的位置以及坡度的抗滑系數(shù)對邊坡工程的實際應用非常重要。目前,研究邊坡穩(wěn)定性的主要方法有地質力學法(赤平投影法)、極限平衡法、有限單元法、塊體單元法、顯式拉格郎日差分法等。極限平衡分析并不考慮巖土體內部的應力應變關系,難以分析邊坡破壞的發(fā)生和發(fā)展過程,在求解安全系數(shù)時需要事先假定滑裂面形狀和進行條間力簡化等。利用有限元進行邊坡穩(wěn)定性分析考慮了巖土體的本構關系,可以合理的反映巖土體的真實工作狀態(tài),不需要事先假設滑動面和進行條間力的簡化等。但工程中的有限元計算雖然能反映邊坡的變形特征、應力狀態(tài)和破壞機理,但是尚不能給出一個直觀的安全系數(shù)。剛體極限平衡分析方法,雖然能夠得到直觀的安全系數(shù)和危險滑弧位置,但是在破壞機理研究上有一定的局限。邊坡工程的研究需要發(fā)展一種將有限元計算結果同極限平衡分析方法相結合的極限平衡有限元方法。本文從高邊坡穩(wěn)定性分析現(xiàn)狀入手,將圖論中解決最短路問題的Dijkstra算法引入到搜索邊坡最危險滑動面及其安全系數(shù)中去,將邊坡的穩(wěn)定問題轉化為圖論中尋找最短路問題,對Dijkstra算法作了相應的改進,編制了相應的計算程序。通過一個經(jīng)典邊坡考核算例對基于Dijkstra算法的極限平衡有限元方法進行驗證,計算結果表明,該方法應用于巖土邊坡穩(wěn)定性分析是合理可行的。將本文提出的理論方法應用于糯扎渡水電站右岸泄洪洞出口邊坡的穩(wěn)定性分析。為了比較分析,采用邊坡穩(wěn)定性分析的極限平衡法和強度折減法對糯扎渡右岸泄洪洞出口邊坡進行了計算分析。計算結果表明:用極限平衡法計算得到的邊坡安全系數(shù)比較小;強度折減法計算得到邊坡安全系數(shù)比極限平衡計算和基于Dijkstra算法的極限平衡有限元方法計算得到的邊坡安全系數(shù)要大;而基于Dijkstra算法的極限平衡有限元方法計算得到的安全系數(shù)介于兩者之間。3種方法計算得到的潛在危險滑動面位置基本一致,該方法應用于復雜的巖體邊坡是可行的。通過對洛古水電站右岸壩址邊坡的穩(wěn)定性分析,闡述了基于Dijkstra算法的極限平衡有限方法對于逆坡向發(fā)展的滑動面的安全系數(shù)計算處理方法。比較分析該方法和3DEC的計算結果,表明基于Dijkstra算法極限平衡有限元方法適用于逆坡向發(fā)展的滑動面的穩(wěn)定性分析。1極限有限元法高邊坡穩(wěn)定性分析目前使用較多的研究方法有以下5種:剛體極限平衡法、極限分析法(上、下限分析法)、有限元(有限差分)法、離散元法和非連續(xù)變形分析法等。1.1平衡方程的性質剛體極限平衡法是邊坡穩(wěn)定分析領域中最古老,也是目前工程應用較多的一種方法。它以Mohr-Coulomb抗剪強度理論為基礎,先設定滑體的滑動面,并將滑動面所劃分出來的滑體劃分成若干垂直條塊,建立作用在這些條塊上力及力矩的平衡方程,根據(jù)這些平衡方程求解出各個條塊受到的滑動力及其所能提供的抗滑力,再將這些滑動力和抗滑力疊加得到整個滑動體上產(chǎn)生的滑動力和所能承受的抗滑力,兩者相除得到安全系數(shù)sF,這通常稱為條分法。由于極限平衡法是通過一定的假定將本質上不靜定的問題簡化為靜定問題來求解,這一作法的合理性問題一直受到人們的普遍關注。總結下來,剛體極限平衡法有以下問題:(1)條塊之間不允許出現(xiàn)拉力;(2)條分面上的剪應力不超過Mohr-Coulomb準則提供的抗剪強度;(3)剛體極限平衡法不能搜索逆坡向發(fā)展的滑弧,它會出現(xiàn)如圖1所示的錯誤的條塊;(4)剛體極限平衡方法需要事先假定滑弧位置,所以其滑動面形狀就為規(guī)則的幾何形狀;(5)剛體極限平衡方法的破壞準則主要基于Mohr-Coulomb準則。1.2應力、應變分析有限單元法和有限差分法是近年來比較活躍的一種邊坡穩(wěn)定分析方法。它全面滿足靜力許可、應變相容和應力、應變之間的本構關系;不受邊坡幾何形狀的不規(guī)則和材料的不均勻性限制;破壞面的形狀和位置不需事先假定;分析邊坡破壞的發(fā)生和發(fā)展過程,模擬邊坡開挖及加固的施工過程,考慮巖土體與支擋結構的共同作用及其變形協(xié)調;其分析結果可以提供應力、應變的全部信息。應用于邊坡工程實踐主要分為間接法和直接法。(1)利用剛體極限平衡的方法搜索危險滑弧和安全系數(shù)間接法是采用剛體極限平衡法原理分析有限單元法或有限差分法的計算結果,應用剛體極限平衡的方法搜索危險滑弧和安全系數(shù)的一種方法。和剛體極限平衡法相比,現(xiàn)有的間接法除了其計算結果符合各種力學原理外,仍然要事先假設滑弧位置,其采用的材料模型亦局限于Mohr-Coulomb準則。(2)基于折減和容重增加有限元法的滑動面參數(shù)確定直接法直接使用有限單元法或有限差分法,通過不斷降低邊坡巖土體強度或增加巖土體的自重使邊坡巖土體達到臨界狀態(tài),然后得到邊坡的安全系數(shù)和潛在危險滑動面。采用這類方法無需事先假定滑動面的形狀和位置,通過不斷降低巖土體的強度或增加巖土體的自重,破壞將很“自然地”發(fā)生在邊坡巖土體抗剪強度不能抵抗剪應力的位置,從而得到最危險滑動面及相應的安全系數(shù)。這兩種方法分別稱之為強度折減有限元法和容重增加有限元法。對于這兩種方法,特別是強度折減法,以下問題還處于爭議之中:(1)破壞準則還存在疑議;(2)實際計算得到的安全系數(shù)會偏大;(3)很多情況下不能得到滑弧位置;(4)折減標準不能統(tǒng)一(等比例或是不等比例折減)。2邊坡穩(wěn)定與dijkt算法2.1x0和x0的最短路問題圖論的應用非常廣泛,它是應用數(shù)學的一個重要分支,起源于著名的柯尼斯堡七橋問題。圖論要解決的問題就是圖上任意兩點之間的最短路徑問題。最短路問題的定義為:設有一個交通系統(tǒng)連接著若干個目標,x0是該系統(tǒng)中的一個固定目標。在該系統(tǒng)中試求從x0到其他各目標的最短路線。圖2描述了一個從起點x0到終點xn之間的最短路問題。最短路問題有解的一個必要條件是D(G)中存在根在x0的支撐樹形圖;反之,若D(G)中存在根在x0的支撐樹形圖,則最短路問題一定有解。以下我們假定D(G)中存在根在x0的支撐樹形圖。(2)用于邊坡穩(wěn)定性問題的最短路問題在邊坡有限元計算分析過程中,通?？梢缘玫秸麄€邊坡的應力場和位移場,但是并不能給出一個直觀的安全系數(shù)和最危險滑動面的位置,而邊坡的安全系數(shù)對于工程實際卻有非常重要的意義。如何利用有限元計算得到的應力場來求解邊坡的最危險滑動面及其安全系數(shù)是一個非常有意義的研究課題。圖3為一簡單邊坡的有限元計算網(wǎng)格,如果將節(jié)點看作網(wǎng)絡圖中的節(jié)點,節(jié)點與節(jié)點間的單元的邊看作網(wǎng)絡圖中的邊,那么,這個剖分好的網(wǎng)格就可以表示成圖論中所描述的一個圖,并且,這個圖是一個簡單無向圖。在邊坡穩(wěn)定性問題中最小安全系數(shù)的定義為:滑動體沿一條滑動面滑動,其抗滑力除以滑動力(另外一種說法是其中每點的最大應力球半徑除以現(xiàn)實的應力球半徑)的值達到最小。由安全系數(shù)的定義可以發(fā)現(xiàn),其隱藏的是一個最短路問題。因為,滑動面是一個與基礎邊坡相獨立的面,其邊界就是邊坡的一個邊界到另外一個邊界的連線,如果將滑動面進入點和滑出點放在網(wǎng)絡圖的節(jié)點上,圖3中邊坡的最危險滑動面可以看作一個從起點x0到終點xn之間的最短路問題。這樣,邊坡的穩(wěn)定性問題就可以用最短路問題來表達,邊坡的最危險滑動面及其安全系數(shù)可以用求解最短路問題的算法解決。解決最短路問題的算法有兩種:(1)Dijkstra算法;(2)Floyd算法。由于在邊坡穩(wěn)定性問題分析過程中,雖然Floyd算法比Dijkstra算法計算效率更高,但由于Floyd算法在搜索到負的安全系數(shù)的時候,程序就會停止計算,容易得到錯誤的計算結果,所以最終選擇Dijkstra算法來求解最短路問題。(3)ra算法的計算步驟當所有道路的長wij≥0時,假設S是頂點集V的真子集,且vs∈S,并以表示V/S。若P=vs,…,vi,vi是圖4中從起點vs到的最短路,則顯然vi∈S,,且從vs到vi的最短路即為從vs沿P到vi的路。Dijkstra算法正是基于這種思想構造的,Dijkstra算法的計算步驟:(1)令P(vs)=0,λ(vs)=0,S0={vs},對vj≠vs,令T(vj)=∞,i=0,k=s;(2)對每個,取T(vj)=min{T(vj),P(vk)+wkj},若T(vj)=P(vk)+wkj,則令λ(vj)=k;(3)計算,其值記為T(vl),若T(vl)<∞,則令P(vl)=T(vl),Si+1=Si∪{lV},k=l,否則停止,此時對vj∈Si,d(vs,vj)=P(vj);而對,d(vs,vj)=T(vj);(4)若i=v-1,則停止,此時對,d(vs,vj)=P(vj),若i<v-1,令i=i+1,轉(2)。Dijkstra算法計算所需的時間與v2成正比。2.2樹形圖的生長對于邊坡穩(wěn)定性問題來說,其數(shù)學表達轉換為圖論中的賦雙權簡單網(wǎng)絡圖來表達。也就是說在G=(V(G),E(G),W(G))中,尋找一棵根在0v的最小樹形圖,其末梢位于lv或vm,…,vn。由G=(V(G),E(G),W(G))所定義賦雙權簡單網(wǎng)絡圖中,必定存在一棵根在v0的最小樹形圖,其末梢位于vl或vm,…,vn。因為圖是簡單圖且連通,那圖必存在樹,由于樹形圖的節(jié)點包含于樹的節(jié)點,所以簡單圖中必含有樹形圖,又因為圖的邊有限,所以其根在v0且末梢位于vl或vm,…,vn的樹形圖的棵數(shù)有限,而每棵樹形圖都帶有一個固定的權值,所以在所有的樹形圖中,必存在最小權值的樹形圖。這就保證了由邊坡穩(wěn)定性問題所演化得到的數(shù)學模型有解。(2)計算時約束滑動力vda)Dijkstra算法的本質是每做一步保證所得到的樹形圖的末梢節(jié)點是從根v0過來的最短路。所以,要保證在計算到vi時,路P(v0,vi)為最短路。b)Dijkstra算法中還有一個前提條件是權值為非負,但是在實際計算中,滑動力是可能為負值的,由安全系數(shù)的計算公式得到:Fs=τf/τ,Fs也可能為負值,這個也需要在計算中進行修改。c)計算中需要的數(shù)據(jù):滑動力與抗滑力。在以安全系數(shù)為目標函數(shù)的最短路求解過程中,由于Dijkstra算法的前提條件已經(jīng)不能用于雙賦權簡單圖最短路求解的過程,所以必須對Dijkstra算法進行修改。2.3修改dijksta算法(1)抗滑力安全問題的處理對于模型里的負權,在安全系數(shù)的定義中可以知道,負的滑動力,就是說是滑體沿滑動面反向滑動的力,這對安全系數(shù)是有利的,因為由于所有的抗滑力都是正值,如果有一個滑動力為負值,相當于使得分母變小,從而安全系數(shù)會變大,所以負的抗滑力對邊坡安全是有利的,本文對它的處理是,定義大小關系如下:在這個關系中,第一個關系要首先被滿足,在滿足第一個關系后,再進行第二個關系的判別。由這個大小的定義可得:0是最小的數(shù),所有數(shù)值都比0大,這與Dijkstra算法中的非負權就等價了。(2)最短路的出現(xiàn)—保證路P(v0,vi)為最短路對于這個問題,有一種最直接的方法,就是窮舉所有以v0為起點為終點的路,然后比較其大小,而不是Dijkstra算法中比較所有的路。但是這種方法是不可行的。這里,需要研究P(v0,vi)為最短路的幾個必要條件。下面是這個問題的數(shù)學表達:假設v0為起點,vi∈S為終點的所有已經(jīng)求得的路為最短路,其所占用的邊集為E(S),其沒有占用的邊集為,如果由Dijkstra算法得到的路P(v0,vj),為最短路,需要滿足以下兩個條件。a)如果在中的邊的安全系數(shù)sF(e)中存在小于路P(v0,vj)的安全系數(shù)值sF(P),那么,將不能保證由Dijkstra算法得到的路P(v0,vj),為最短路,因為,其安全系數(shù)的路將有可能通過那些安全系數(shù)比P(v0,vj),小的路,使得那些安全系數(shù)本來比P(v0,vj),大的路P(v0,vi),vi∈S的安全系數(shù)減少,從本質上說和Dijkstra算法里面要求路的權值是遞增的相違背。所以要求:b)當條件(1)被滿足了之后,這里還不能確定這條路P(v0,vj),為最短路,如圖4所示。當所有屬于的邊的安全系數(shù)大于等于Fs(P)時,我們注意到:假設S是頂點集V的真子集,且vs∈S,并以表示V/S。若是圖4中從起點vs到的最短路,則顯然vi∈S,,且P(vs,vi)為從vs到vi的最短路。我們用反證法得到其成立的條件,現(xiàn)在假設Q(vs,vi)為從vs到vi的最短路,則可令:則要使得下式成立:由于Fs(P′(vs,vj))比Fs(P(vs,vj))小,這與P是從起點vs到的最短路矛盾。這個時候,就能確定這條路P(v0,vj),為最短路。如式(9)成立,則通過反證法就可推導得到P(vs,vi)為最短路,將此問題轉化為數(shù)學表達式,得:將fs(vj)表示為將sF(P)表示為將Fs(Q)表示為那么存在以下關系:在這個前提下,式(9)可表示為如果式(9)得到滿足,則式(9)可以推出矛盾,所以假設條件就不成立,那么P(vs,vi)為從vs到vi的最短路成立。設d=k1f,c=k1(e+n1),b=k2f,a=k2(e+n2),由題設關系可得n2>n1>0,由于前提條件(1)可得,這里可認為a,b,c,d,e,f均大于0,則可以得因為f>0,因為k1>0,k2>0,不等式兩邊同除k1k2得那么因為k1>0,k2>0,如果等式右邊小于等于0,那么等式自然成立,此時因為n2>n1>0,所以若等式右邊大于0,即此時要滿足式(18),那么,現(xiàn)在就得到了保證式(16)成立的條件:由式(23)可知,當k1,k2滿足式(23)后,式(9)大于0,與假設矛盾,P(vs,vi)為從vs到iv的最短路成立。當滿足了式(4)、(5)的條件之后,Dijkstra方法引入邊坡穩(wěn)定性分析的前提條件就能滿足,就可以用數(shù)學方法嚴格的證明所求得的最危險滑動面及其安全系數(shù)為在所有危險滑動面及安全系數(shù)中的最小值。3在確定滑動面的限制條件下運行并進行優(yōu)化引入了Dijkstra算法后,是不能直接用于計算的,因為Dijkstra算法的退出條件是:若Sk包含所有樹形圖的末梢節(jié)點,而要將節(jié)點加入到Sk就必須滿足式(4)、(5)中規(guī)定的條件。對于式(4)中的條件,滿足并不困難,只需定義一個新的必要運算。對于式(23),由于d,f為最短路或假設的最短路的滑動力,而b為單元邊上的滑動力,所以d≈f>>b,可以推導得k1∈[1±Δ],0<k2<<1,n2>n1且n2≈n1,這個條件在計算中很容易得到滿足。而對于式(4),要滿足是非常困難的,因為:如果存在一條邊,其安全系數(shù)比較小,但其位置離最短路比較遠,最短路不會通過該邊,于是,式(5)就一直不能得到滿足,直到搜索到該邊,而把這條邊加入到P(vs,vj)后,條件才能得到滿足,但是,這個過程是非常漫長的,所以程序不能直接應用于邊坡穩(wěn)定性分析中去。所以,這里必須引入一些限制條件,但這些限制條件應該是符合邊坡的力學特性,這里引入了3個限制條件,如圖5所示的3種情況。圖5所示的情況中,說明了以下問題:(1)邊坡的滑動面不會出現(xiàn)尖角(二階導數(shù)很大或兩相鄰邊的夾角小于90°)。要出現(xiàn)這種情況,在尖角處必會出現(xiàn)非常大的應力矢量方向的改變,或者說出現(xiàn)非常大的應力集中,而且這個應力的方向與整體方向的夾角超過90°,這種情況在實際有限元計算中是不可能出現(xiàn)的。因為,出現(xiàn)很大應力集中后,邊坡已經(jīng)破壞了,如圖5中(b)、(c)所示。具體這個夾角的度數(shù),可以在程序中根據(jù)實際情況加以限制。(2)邊坡滑動面的發(fā)展不會出現(xiàn)鋸齒形,如圖(c)所示,這種情況在局部小范圍內可能會出現(xiàn),但是在整體的邊坡滑動中,是不符合力學性質的。(3)在通常情況下,邊坡滑動面不會向坡內發(fā)展(如圖5(a)中的邊坡為向“右”發(fā)展),如果這個條件被確立,將會大大的減少計算時間,但是這個條件在計算中也是可以被屏蔽掉的。對于第一點,其圖形學的表達為:在一條最短路上,連接同一個節(jié)點的兩條邊的夾角α要大于90°,如圖6所示,這個條件的判別的數(shù)學表達為將這兩個條件加入導Dijkstra算法中的判別式去,使得程序的運算更為有效。由于規(guī)定了滑動面的“生長”方向,那么式(4)的判別范圍將被限定在滑動面起點與所到達的節(jié)點之間,這將大大簡化判斷滑動面是否為最短路所需要的計算量,也能更快的找到該節(jié)點所對應的最短路。另外規(guī)定了滑動面“生長”方向后,使得搜索的滑動面數(shù)量大大的減少,優(yōu)化了程序的運行。(1)減少了潛在的最短路的數(shù)量,由于許多“最短路”的發(fā)展不符合上述兩個要求,使得這個“最短路”的搜索提前結束,節(jié)省了許多計算消耗。(2)增大了Dijkstra結束條件成立的概率。當有了上述兩個條件的限制之后,式(4)的比較區(qū)域就變?yōu)?(1)對于順坡向發(fā)展的滑動面,其式(4)的判斷區(qū)域變?yōu)閳D7中的III區(qū);(2)對于存在逆坡向發(fā)展滑動面,那么搜索區(qū)域在這里的處理方法是,給定一個與滑動面開展深度h成比例的發(fā)展最大寬度l,搜索區(qū)域就處于如圖8所示的III中。經(jīng)過優(yōu)化后,可得到一有效性滿足計算要求的程序流程。經(jīng)過優(yōu)化后的程序運行步驟如下所示。(2)對每個vj∈Si,如果滿足式(2)與式(3),則計算若當k=k0時,上式達到最小,則令λ(vj)=k0。(3)計算,其值記為T(vl),若T(vl)<∞,則令:式中,k=λ(vl)。如果不能滿足兩條件中的一個,檢查PT內是否存在從vs到vl的路,如果不存在或者存在但Fs(P(vs,vi))<PT(vl)Fs,記為否則,不保存這次結果。(5)若Sk包含所有樹形圖的末梢節(jié)點,則停止,此時P(vs,vEi)=PT(vEi)P,vEi∈SE,否則轉(2)。經(jīng)過優(yōu)化后,可得到一有效性滿足計算要求的程序流程。4輸入數(shù)據(jù)的計算采用dijksta算法Dijkstra算法所需要的輸入?yún)?shù)是每條邊的抗滑力與滑動力。4.1關聯(lián)節(jié)點的確定這里使用一個半徑為R圓確定與固定節(jié)點相關聯(lián)的節(jié)點,如圖9所示。r>R時,認為這兩個節(jié)點沒有關聯(lián),反之r<R時,兩節(jié)點相關聯(lián)。且使用自適應半徑,以有限元計算中單元的最大直徑作為R,當關聯(lián)節(jié)點大于規(guī)定數(shù)目m時,取半徑最短的前m個節(jié)點作為關聯(lián)節(jié)點。但是,當所選取的關聯(lián)節(jié)點有Dijkstra算法的輸出節(jié)點的時候,以輸出節(jié)點優(yōu)先,這里的認為是將兩點之間的半徑減少1/3后再使用前面的算法判定關聯(lián)節(jié)點。4.2vivj邊的滑動力在一環(huán)如圖10所示,滑動力的定義:滑動力是指滑動體下滑對邊坡“基座”所產(chǎn)生的力。對坡向向左的邊坡而言,逆時針為滑動力正方向。在用圖形表達這個抗滑力時,存在3種情況:邊與x軸夾角為銳角、鈍角、滑動面向坡內發(fā)展,如圖11~13所示。對于前面兩種情況,都將邊vivj以上的部分看作是邊坡滑動體的一部分,因為以滑動體上的滑動力為逆時針方向的切應力,所以vivj邊上如圖11,12所示的方向為正方向,而在有限元的定義中,正應力以面的外法線方向為正。根據(jù)上面定義τy以x軸正方向為正,但是在滑動體的切應力定義中,它是負(順時針)的。但是在應力圓的計算中,以順時針為正,如果以這個方法定義滑動體滑動力的正方向,那么滑動力是負的,而τy是正的,那么在基于y軸正平面用應力圓計算滑動力的時候,最終結果與計算結果是反號的。對于逆坡向發(fā)展的邊,用相同的分析方法得到相同的結果。如圖14所示,滑動力的計算公式為4.3o的x+y/2這里提出了基于不同材料模型的抗滑力計算公式。最大承受剪應力大小為的長度,計算公式為其中,σo為應力圓的中心,由圖15可知:σo=(σx+σy)/2。(30)R為最大抗滑力的應力圓的半徑,由圖14可以看出,R應該是點σo到屈服邊界(Mohr-Coulumb準則的屈服邊緣——直線τ=-σtan?+c)的距離,根據(jù)點到直線的距離的公式:要注意的一點是,在Mohr-Coulumb準則計算抗滑力的時候,所用的σx,σy應該是該單元的最大正應力和最小正應力。(2)xy型滑動面的抗滑力的形成根據(jù)Drucker–Prager準則的屈服面計算公式:對于式(32),可以寫成一個關于σx,σy,σz,τxy,τzx,τyz的函數(shù),在三維問題中,要求的是滑動面所在平面的τ,而其他5個參數(shù)可以當成已知條件代入,本文只對平面問題進行分析。在平面問題中τzx,τyz=0(在平面應力問題里面σz=0),t可以表示成一個關于τxy的函數(shù):令q=q(τxy),將r=r(τxy)表示成q的函數(shù):那么t可以表示成一個的q函數(shù)。而在式(32)中的其他3個參數(shù)p,d,β都是與τxy無關的。所以t(τxy)=ptanβ+d,需要注意的是,這個時候σx為滑動面所在的平面,而不是笛卡兒坐標系下的x正平面上的正應力,σy為與滑動面所在平面垂直的平面上的正應力,σz為平面外的正應力(平面應力問題中σz=0)。那么式(32)就可以表示為一個q的三次多項式(等式兩邊同乘q2)。其中,由F=0求解出q,再由q的定義就可以解出材料破壞時,滑動弧所在面上所能承受的剪切應力τxy,也就是滑動面上的抗滑力。本程序所需要輸入的是抗滑力與滑動力在單元邊上的值,所以,無論是經(jīng)典的取值法還是文中提出的取值法都適用于該程序。4.4材料界面的實現(xiàn)根據(jù)前面的理論,編制了相關的計算程序及前處理后處理程序,操作界面如圖16所示,后處理界面如圖17所示。程序需要的輸入數(shù)據(jù)有:(1)節(jié)點信息(每種材料作為一個輸入文件),信息包括:節(jié)點號、節(jié)點x坐標、節(jié)點y坐標、σx,σy,σz,τxy;(2)單元信息(計算網(wǎng)格數(shù)據(jù));(3)每種材料的材料參數(shù);(4)滑動面進入點和剪出點信息和上邊界信息。5基于dijpstor算法的算例分析以下通過一個經(jīng)典考題來驗證本文的邊坡穩(wěn)定性分析方法的合理性。澳大利亞計算機協(xié)會于1987年委托BDonald和PGiam設計ACADS邊坡穩(wěn)定性分析程序的考核題,第一道考核題是個二維均質土坡,邊坡基本形狀見圖18,材料參數(shù)見表1。28家單位提供了答案,陳祖煜與Fredlund,Baker一起參與了計算工作,考核題最終推薦答案的邊坡安全系數(shù)為1.00,最危險滑動面位置及形狀見圖18。由于該算例調查的規(guī)模較大,得到的分析成果比較可靠。因此,為了驗證基于有限元計算結果的改進Dijkstra算法的極限平衡有限元方法的正確性,這里將利用上述考核算例邊坡作為一個簡單算例。采用基于Dijkstra算法的極限平衡有限元方法對這個算例進行分析,算例邊坡的網(wǎng)格剖分見圖19,并且同時采用了強度折減法對該邊坡進行了分析,計算得到的潛在危險滑動面位置和邊坡安全系數(shù)見圖20。利用強度折減法計算得到的考核算例邊坡安全系數(shù)為1.02,利用基于Dijkstra算法的極限平衡有限元方法得到考核算例邊坡安全系數(shù)為1.0015673,與最終推薦答案的邊坡安全系數(shù)1.00比較接近;并且計算得到的潛在危險滑動面位置也基本一致,這說明采用基于Dijkstra算法的極限平衡有限元方法來分析邊坡的穩(wěn)定性是合理可行的。本文建立的邊坡穩(wěn)定性分析方法考慮了巖、土體材料的本構關系,可以反映真實的力學性態(tài),不需要事先假設滑動面和進行條間力的簡化,可以自然的形成滑動破壞面。6茂佐渡右岸排洪洞邊坡分析為了對比分析,對糯扎渡右岸泄洪洞出口邊坡還采用了極限平衡方法和強度折減法計算邊坡的最危險滑動面及其安全系數(shù)。(1)邊坡穩(wěn)定性分析極限平衡方法分別采用簡化Bishop法、Janbu法、Spencer法和Morgenstern-Price法4種方法進行了邊坡穩(wěn)定性分析,安全系數(shù)取4種方法計算得到的安全系數(shù)的平均值。(2)邊坡臨界強度指標基于強度折減法的有限元計算過程為:保持巖土體的重力加速度為常數(shù),通過逐步減小抗剪強度指標,將c,?值同時除以折減系數(shù)Fs,得到一組新的強度指標c′,?′再進行有限元,反復計算直至邊坡達到臨界破壞狀態(tài),此時采用的強度指標與巖土體原具有的強度指標之比即為該邊坡的安全系數(shù)Fs。當對巖體力學參數(shù)進行不斷折減的時候,判定一個邊坡是否達到臨界破壞的依據(jù)是:選取幾個具有代表意義的關鍵點(在選取關鍵點時,參照極限平衡的計算結果,使得關鍵點的位置盡可能處于潛在危險滑動面上),當巖體力學參數(shù)不斷折減時,關鍵點的位移會在某一時刻發(fā)生迅速增大變化,這一時刻即為邊坡的臨界狀態(tài),此時的折減系數(shù)作為邊坡的安全系數(shù)。6.1邊坡穩(wěn)定性分析糯扎渡水電樞紐右岸泄洪隧洞出口地段邊坡布置于工程地質F區(qū),工程地質條件較差。出口地段地形坡度30°~45°,開挖邊坡與#3導流洞出口邊坡相連,洞臉邊坡最大開挖坡高以#3導流洞出口處最大,約為126m,洞臉邊坡走向約N25E,傾向SE,設計開挖坡角為48。開挖邊坡部位地表有1~3m厚的坡積層,局部基巖裸露,基巖巖性為花崗巖,邊坡巖體主要為弱風化,坡頂有部分全、強風化巖體(厚度約16m左右),F28斷層(產(chǎn)狀為N15—27W,SW75—85)斜切坡面。F38、F39、F28斷層帶及全風化、部分強風化花崗巖為Ⅴ類邊坡巖體,強風化及部分弱風化上部花崗巖為Ⅳ類邊坡巖體。全、強風化巖體及F28斷層帶組成的邊坡,穩(wěn)定性差,弱風化花崗巖組成的邊坡,其穩(wěn)定性主要受結構面控制。為了準確估計泄洪洞出口邊坡的穩(wěn)定性,選取了3個具有代表意義的二維剖面進行邊坡穩(wěn)定性分析,剖面布置如圖21所示:其中1-1剖面沿泄洪洞軸線;2-2剖面沿左右坡拐角線附近;3-3剖面經(jīng)過第(5)組節(jié)理帶,為可能產(chǎn)生平面滑動剖面。限于篇幅,這里只對3個剖面中有代表意義的3-3剖面進行邊坡穩(wěn)定性分析。對于3-3剖面,首先采用極限平衡法和強度折減法進行邊坡穩(wěn)定性分析;然后應用基于Dijkstra算法的極限平衡有限元方法對邊坡進行穩(wěn)定性分析。3-3剖面的巖體材料分層如圖22所示。有限元計算時采用的巖體力學參數(shù)如表2所示。計算時,在極限平衡與強度折減計算中,節(jié)理和斷層采用弱化材料來模擬;而基于Dijkstra算法的極限平衡有限元計算中,節(jié)理和斷層采用結構面來模擬。6.2極端平衡法的計算極限平衡法計算結果如表3所示。圖23為3-3剖面用極限平衡法計算得到的潛在危險滑動面,工況分別為自然邊坡和開挖完成后。6.3位移對比分析極限平衡有限元邊坡穩(wěn)定性分析的強度折減法采用FLAC3D有限差分軟件進行計算。結合極限平衡的計算結果,關鍵點的位置盡可能處于潛在危險滑動面上,開挖前取2個關鍵點:a點(取豎向位移)和b點(取水平位移);開挖后取c點(取豎向位移)和d點(取水平位移),位置見圖24。圖25,26給出了邊坡兩種工況下位移與強度折減系數(shù)的關系曲線。圖中可看出,自然邊坡在巖體力學參數(shù)折減到1.8時a點的豎向位移和b點的水平位移開始發(fā)生急劇變化,認為此時的邊坡安全系數(shù)為1.8;邊坡開挖完成后,材料參數(shù)折減到2.0時c點的豎向位移和d點的水平位移開始發(fā)生急劇變化,認為此時的邊坡安全系數(shù)為2.0。6.4數(shù)值計算結果為了方便基于Dijkstra算法的極限平衡有限元方法的數(shù)據(jù)前處理,采用ABAQUS進行有限元計算,然后利用計算得到的各種工況下的應力–應變數(shù)據(jù)、網(wǎng)格信息和材料單元信息,通過上述程序進行邊坡穩(wěn)定性分析,計算結果如圖27,28所示。6.5邊坡變形面溫度比較分析3-3剖面極限平衡法、強度折減法和基于Dijkstra算法的極限平衡有限元方法的計算結果,如表4所示。從表4可以看出,采用剛體極限平衡計算得到的安全是最小的,邊坡在自然狀態(tài)下安全系數(shù)為1.552,開挖完成后安全系數(shù)為1.578;而強度折減法計算得到安全系數(shù)比極限平衡計算和基于Dijkstra算法的極限平衡有限元方法計算得到的安全系數(shù)要大,邊坡在自然狀態(tài)下安全系數(shù)為1.800,開挖完成后安全系數(shù)為2.000;基于Dijkstra算法的極限平衡有限元方法計算得到的安全系數(shù)在兩者之間,邊坡在自然狀態(tài)下安全系數(shù)為1.595,開挖完成后安全系數(shù)為1.683。在各種工況下,極限平衡和基于Dijkstra算法的極限平衡有限元方法計算得到的潛在危險滑動面位置基本一致。從圖23(b)和圖28可以看出,開挖完成后危險滑動面的位置以表層滑動為主,從邊坡的坡積層開始,沿著N25°—35°W和SW45°—70°順層節(jié)理,然后從強風化層滑出。7基于dijpstor算法的極限平衡有限元方法對于逆坡向滑動面來說,在極限平衡方法進行條分的時候,所劃分得到的計算條塊必須要和坡面相交,如果遇到滑動面逆坡向發(fā)展的時候,將無法正確的劃分計算條塊,而且在計算滑動面上滑動力的時候,逆坡向發(fā)展的滑動面將無法計算其滑動力與抗滑力,這將造成安全系數(shù)計算的不準確,所以一般的剛體極限平衡方法不容許滑動面逆坡向發(fā)展,也就沒法搜索得到正確的結果,如圖1所示?；贒ijkstra算法的極限平衡有限元方法對這個問題迎刃而解,因為該方法采用單元的切應力作為滑動力,只要計算得到單元邊上的切應力,乘以邊的長度,就可以求得危險滑動面上的滑動力,不會產(chǎn)生錯誤的條分的概念,這就是基于Dijkstra算法的極限平衡有限元方法與其他研究方法不同的地方,也是該方法的優(yōu)勢所在。在應用基于Dijkstra算法的極限平衡有限元方法來計算逆坡向發(fā)展的危險滑動面的安全系數(shù)時,需要對滑動面上的滑動力控制條件進行一定的修改,以使得計算得到的安全系數(shù)更為符合實際,具體的修改方法以洛古水電站右岸壩址邊坡為實例進行說明。7.1數(shù)值模型的建立洛古水電站是西溪河流域兩庫五級電站中的第二級,是一座以發(fā)電為主的混合式電站,其水庫為西溪河流域水電梯級近期開發(fā)的龍頭水庫,總庫容3730×104m3,裝機容量110MW。樞紐由碾壓混凝土攔河壩、引水系統(tǒng)、發(fā)電廠房等組成。壩體右岸邊坡坡度較陡,風化層深厚,卸荷帶和斷層發(fā)育,主要斷層為F107等,卸荷裂隙多填充巖塊及次生黃泥。右岸緩傾角結構面傾向坡外,其產(chǎn)狀與巖層流面產(chǎn)狀基本一致,多充填巖屑、方解石,部分夾有泥質,多具有擠壓特征。洛古水電站壩址邊坡在天然狀態(tài)下右岸巖體穩(wěn)定性即處于臨界狀態(tài),在進行工程開挖后右岸邊坡能否保持穩(wěn)定,是施工期安全的一個重要問題。由工程地質條件分析,邊坡的滑動面位置應在是沿F107斷層向邊坡內部發(fā)展,然后沿緩傾角節(jié)理滑出的,但因F107斷層的傾向是逆坡向,所以一般的極限平衡軟件無法搜索這樣的滑動面,無論是指定滑動面還是自動搜索。以下采用基于Dijkstra算法的極限平衡有限元方法和3DEC軟件對該邊坡進行分析。根據(jù)實測資料,建立壩址開挖邊坡的二維有限元模型,見圖29,模型以右岸壩軸線剖面為計算平面,模擬的結構面為F107斷層(傾角75°),通過壩軸線的五條緩傾角節(jié)理(傾角25°)。計算采用的直角坐標系:X軸從河谷指向右岸;Z軸為豎直方向,向上為正。計算范圍側向為河道中心線向右岸220m,豎直方向高程1848~2160m。該邊坡的計算參數(shù)見表5。7.2tf06斷層應力–應變耦合控制洛古水電站右岸壩址邊坡穩(wěn)定性問題的特殊性在于:逆坡向的F107斷層與順坡向的緩傾角節(jié)理在邊坡內部相交,使邊坡形成幾個軟弱的滑動面,危險滑動面的發(fā)展應該是沿著F107斷層向邊坡內部發(fā)展,然后順著某條緩傾角節(jié)理滑出。根據(jù)3DEC的塊體單元計算結果,表明滑動體是由F107斷層與緩傾角節(jié)理所圍成的塊體單元組成。通過有限元計算分析,可以得到整個邊坡的應力–應變狀態(tài),以下來分析一下F107斷層的應力狀態(tài)。圖30為F107斷層的最大主應力圖,可以看出,在F107上面存在拉應力區(qū)。對于本文的方法來說,各個連接上計算得到的滑動力將都是負值,而在緩傾角節(jié)理上的滑動力都是正值,這就使得滑動力沿F107斷層和緩傾角節(jié)理發(fā)展的時候滑動力相互抵消,這樣會造成計算得到的安全系數(shù)偏大。但是從力學角度