基于Floyd算法的裂隙巖質(zhì)斜坡穩(wěn)定性分析

蘭志廣，陳潤(rùn)橋，聶細(xì)江，郭 猛，劉明亮

(天津市勘察設(shè)計(jì)院集團(tuán)有限公司,天津 300191)

1 概述

裂隙巖質(zhì)斜坡的穩(wěn)定性多由賦存于坡體內(nèi)的裂隙決定[1]。而由于裂隙的各向異性以及空間復(fù)雜性，針對(duì)裂隙巖質(zhì)斜坡的穩(wěn)定性研究仍十分困難。目前，針對(duì)裂隙巖質(zhì)斜坡穩(wěn)定性研究通常包括最危險(xiǎn)滑動(dòng)面的確定以及安全系數(shù)計(jì)算[2]。許多學(xué)者采用動(dòng)態(tài)規(guī)劃法或者智能算法來(lái)確定危險(xiǎn)滑動(dòng)面，并取得理想結(jié)果。如張文等[3-4]考慮裂隙的隨機(jī)性，綜合采用改進(jìn)的隨機(jī)動(dòng)力學(xué)法、裂隙頻率法以及Dijkstra最短路徑算法確定了裂隙巖質(zhì)斜坡的最危險(xiǎn)滑動(dòng)面；徐佩華等[5-6]考慮裂隙巖質(zhì)斜坡的應(yīng)力場(chǎng)特征，采用遺傳算法確定了最危險(xiǎn)滑動(dòng)面；段榮福[7]采用不平衡推力法計(jì)算出每個(gè)可能剪出口的剩余下滑力，結(jié)合最優(yōu)理論法確定了斜坡最危險(xiǎn)滑動(dòng)面。在確定最危險(xiǎn)滑動(dòng)面之后，就可以根據(jù)最危險(xiǎn)滑動(dòng)面上巖土體的抗剪強(qiáng)度參數(shù)來(lái)計(jì)算裂隙巖質(zhì)斜坡的安全系數(shù)。此外，隨計(jì)算機(jī)技術(shù)發(fā)展，裂隙巖質(zhì)斜坡穩(wěn)定性數(shù)值分析也取得迅速發(fā)展與廣泛應(yīng)用。如Grifiths，Nour等[8-9]利用有限元模型確定斜坡的最危險(xiǎn)滑動(dòng)面并計(jì)算出安全系數(shù)；韓建國(guó)[10]基于有限差分法，利用強(qiáng)度折減法對(duì)裂隙巖質(zhì)斜坡穩(wěn)定性及破壞模式進(jìn)行了評(píng)價(jià)；程謙恭等[11]利用離散單元法模擬評(píng)價(jià)裂隙巖質(zhì)斜坡穩(wěn)定性。

盡管針對(duì)裂隙巖質(zhì)斜坡的穩(wěn)定性研究取得了許多成果，但綜合考慮裂隙發(fā)育特征并合理確定裂隙巖質(zhì)斜坡安全系數(shù)的方法仍然較少。本文提出了一種基于裂隙網(wǎng)絡(luò)模擬技術(shù)及Floyd算法的裂隙巖質(zhì)斜坡穩(wěn)定性分析方法。該方法首先根據(jù)現(xiàn)場(chǎng)裂隙發(fā)育特征模擬生成主滑面二維裂隙網(wǎng)絡(luò)，然后采用Floyd算法搜索不同剪入-剪出口之間的危險(xiǎn)滑動(dòng)面，最后利用傳遞系數(shù)法計(jì)算所搜索的危險(xiǎn)滑動(dòng)面的安全系數(shù)，以最小的安全系數(shù)作為對(duì)裂隙巖質(zhì)斜坡整體穩(wěn)定性的評(píng)估。

2 裂隙網(wǎng)絡(luò)模擬

根據(jù)巖體結(jié)構(gòu)控制論，巖體穩(wěn)定性是由賦存于巖體內(nèi)部的結(jié)構(gòu)面控制的[12]。因此，在裂隙巖質(zhì)斜坡穩(wěn)定性分析前，首先要分析斜坡巖體內(nèi)裂隙的特征信息。本研究方法是在二維尺度上分析評(píng)價(jià)裂隙巖質(zhì)斜坡穩(wěn)定性，首先要獲得斜坡主滑面二維裂隙網(wǎng)絡(luò)。在實(shí)踐中，主滑面上的裂隙網(wǎng)絡(luò)通常難以調(diào)查到，因此需要根據(jù)露頭面裂隙信息模擬生成主滑面二維裂隙網(wǎng)絡(luò)(見(jiàn)圖1)。具體介紹如下。

2.1 三維裂隙網(wǎng)絡(luò)模擬

三維裂隙網(wǎng)絡(luò)模擬技術(shù)是近半個(gè)世紀(jì)以來(lái)不斷發(fā)展和完善用以表征巖體內(nèi)部裂隙發(fā)育特征的最有效方法。該方法是基于現(xiàn)場(chǎng)實(shí)測(cè)裂隙數(shù)據(jù)，通過(guò)統(tǒng)計(jì)其跡長(zhǎng)、產(chǎn)狀、密度等參數(shù)的概率分布類(lèi)型，采用Monte Carlo模擬生成與現(xiàn)場(chǎng)裂隙統(tǒng)計(jì)相似的數(shù)學(xué)模型[13]。

裂隙巖質(zhì)斜坡巖體內(nèi)部三維裂隙網(wǎng)絡(luò)模擬主要步驟為：

1)裂隙位置模擬。

通常，采用泊松分布的方法在三維空間內(nèi)生成裂隙幾何中心的坐標(biāo)點(diǎn)。

2)裂隙尺寸模擬。

迄今為止，裂隙形狀尚未有統(tǒng)一定論，因此常假設(shè)裂隙形狀為圓盤(pán)、多邊形等。為簡(jiǎn)化建模程序，在工程實(shí)踐中，常將裂隙形狀假設(shè)為圓盤(pán)。基于裂隙形狀為圓盤(pán)的假設(shè)，裂隙尺寸用圓盤(pán)直徑表征。裂隙直徑一般服從正態(tài)分布、對(duì)數(shù)正態(tài)分布以及伽馬分布。

3)裂隙產(chǎn)狀模擬。

裂隙產(chǎn)狀通常服從Fisher分布、Bingham分布、雙變量正態(tài)分布以及經(jīng)驗(yàn)分布。

4)裂隙密度確定。

采用測(cè)線(xiàn)法可以獲得裂隙間距，基于此，可推導(dǎo)出三維空間內(nèi)裂隙密度。

5)Monte Carlo模擬。

采用Monte Carlo模擬對(duì)上述各項(xiàng)參數(shù)進(jìn)行融合，即可生成與現(xiàn)場(chǎng)實(shí)測(cè)裂隙具有相似統(tǒng)計(jì)特征的三維裂隙網(wǎng)絡(luò)模型。

2.2 主滑面二維裂隙網(wǎng)絡(luò)生成

基于裂隙為圓盤(pán)的假設(shè)，并忽略裂隙的空隙及其起伏，那么決定裂隙圓盤(pán)位置、大小與形態(tài)的參數(shù)有圓盤(pán)中心點(diǎn)的坐標(biāo)(xc，yc，zc)、圓盤(pán)直徑D以及裂隙的傾向α與傾角β，以這些參數(shù)即可表達(dá)三維裂隙網(wǎng)絡(luò)模型內(nèi)的裂隙圓盤(pán)。選擇主滑面與三維裂隙網(wǎng)絡(luò)模型相截切，即可得到主滑面上的二維裂隙網(wǎng)絡(luò)。具體地，主滑面二維裂隙網(wǎng)絡(luò)數(shù)學(xué)模型表達(dá)式為：

(1)

其中,A=sinβcosα;B=sinαcosβ;C=cosβ;θ為斜坡主滑面走向;(x0,y0)為斜坡主滑面任一點(diǎn)在三維裂隙網(wǎng)絡(luò)中的坐標(biāo)。

3 基于Floyd算法搜索危險(xiǎn)滑動(dòng)面

由于裂隙的強(qiáng)度參數(shù)相較巖石的強(qiáng)度參數(shù)小得多,因此在裂隙巖質(zhì)斜坡滑動(dòng)破壞時(shí),滑動(dòng)面往往沿著預(yù)先存在的裂隙擴(kuò)展。斜坡滑動(dòng)面經(jīng)過(guò)的裂隙越多,抗滑力越小。當(dāng)我們限定滑動(dòng)面的剪入口與剪出口時(shí),滑動(dòng)路徑的總長(zhǎng)度越小,經(jīng)過(guò)的裂隙長(zhǎng)度越大,則沿著此路徑滑動(dòng)時(shí)所受的抗剪力越小。當(dāng)滑動(dòng)路徑最短時(shí),則最大可能的經(jīng)過(guò)了巖體中的裂隙,沿此路徑滑動(dòng)時(shí)所受的抗剪力最小。由此分析可知,裂隙巖質(zhì)斜坡危險(xiǎn)滑動(dòng)面的滑動(dòng)主要沿剪入口與剪出口之間的最短路徑。搜索最短路徑的方法有很多種,由于Floyd算法容易理解,可以算出任意兩個(gè)節(jié)點(diǎn)之間的最短距離且代碼編寫(xiě)也很簡(jiǎn)單,因此我們采用Floyd算法搜索最短路徑(即危險(xiǎn)滑動(dòng)面)。Floyd算法具體介紹如下:

通過(guò)Floyd算法計(jì)算各個(gè)節(jié)點(diǎn)之間的最短路徑時(shí),需要引入兩個(gè)矩陣(矩陣D和矩陣P),矩陣D中的元素d[i][j]表示節(jié)點(diǎn)i(第i個(gè)節(jié)點(diǎn))到節(jié)點(diǎn)j(第j個(gè)節(jié)點(diǎn))之間的距離;矩陣P中的元素p[i][j]表示節(jié)點(diǎn)i到節(jié)點(diǎn)j經(jīng)過(guò)了p[i][j]記錄的值所表示的節(jié)點(diǎn)。假設(shè)有N個(gè)節(jié)點(diǎn),則需要對(duì)矩陣D和矩陣P進(jìn)行N次更新。初始時(shí),矩陣D中節(jié)點(diǎn)d[i][j]的距離為節(jié)點(diǎn)i到節(jié)點(diǎn)j的權(quán)值;如果i和j不相鄰,則d[i][j]=∞,矩陣P的值為節(jié)點(diǎn)p[i][j]的j值。接下來(lái)開(kāi)始對(duì)矩陣D進(jìn)行N次更新。第1次更新時(shí),如果“d[i][j]的距離”> “d[i][0]+d[0][j]” (d[i][0]+d[0][j]表示i與j之間經(jīng)過(guò)第1個(gè)節(jié)點(diǎn)的距離),則更新d[i][j]為“d[i][0]+d[0][j]”,更新p[i][j]=p[i][0]。同理,第k次更新時(shí),如果“d[i][j]的距離”> “d[i][k-1]+d[k-1][j]”,則更新d[i][j]為“d[i][k-1]+d[k-1][j]”,p[i][j]=p[i][k-1]。更新N次之后,算法完成,結(jié)合矩陣D和矩陣P即可獲得任意兩節(jié)點(diǎn)的最短路徑。

以圖2所示節(jié)點(diǎn)為例,應(yīng)用Floyd算法計(jì)算節(jié)點(diǎn)a與節(jié)點(diǎn)d之間的最短路徑,通過(guò)更新矩陣D和矩陣P,可得到a→d之間的最短路徑為a→c→d,最短距離為3+1=4。

由于Floyd算法針對(duì)節(jié)點(diǎn)進(jìn)行搜索,故在實(shí)際搜索時(shí),需將主滑面二維裂隙網(wǎng)絡(luò)中的裂隙跡線(xiàn)進(jìn)行離散化。裂隙跡線(xiàn)離散化程度越高,搜索結(jié)果越精確,但計(jì)算量也越大。此外,在實(shí)際搜索時(shí),也需要預(yù)先對(duì)滑動(dòng)面剪入口與剪出口進(jìn)行設(shè)定。在實(shí)際斜坡滑動(dòng)時(shí),危險(xiǎn)滑動(dòng)面會(huì)優(yōu)先選擇剪入口與剪出口處發(fā)育的裂隙,故在實(shí)際搜索過(guò)程中,將二維裂隙網(wǎng)絡(luò)的邊界假設(shè)為裂隙。離散化程度隨主滑面尺寸而定。

4 傳遞系數(shù)法計(jì)算安全系數(shù)

得到斜坡危險(xiǎn)滑動(dòng)面后,我們可以根據(jù)GB 50330—2013建筑邊坡工程技術(shù)規(guī)范計(jì)算斜坡的安全系數(shù)。由于基于Floyd算法搜索所得的危險(xiǎn)滑動(dòng)面形態(tài)為折線(xiàn)形,故采用傳遞系數(shù)法進(jìn)行斜坡安全系數(shù)計(jì)算。

根據(jù)傳遞系數(shù)法(見(jiàn)圖3),斜坡安全系數(shù)(Fs)可按下列公式計(jì)算:

Ri=Wicosαitanφi+cili。

Ti=Wisinαi。

Ψj=cos(αi-αi+1)-sin(αi-αi+1)tanφi+1。

其中，Ri為第i計(jì)算條塊單位寬度重力引起的抗滑力;Ti為第i計(jì)算條塊單位寬度重力引起的下滑力;Ψj為第i計(jì)算條塊對(duì)第i+1計(jì)算條塊的傳遞系數(shù);ci為第i計(jì)算條塊滑面的黏聚力;φi為第i計(jì)算條塊滑面的內(nèi)摩擦角;li為第i計(jì)算條塊滑面的長(zhǎng)度;αi為第i計(jì)算條塊滑面的傾角。

5 舉例分析

應(yīng)用上述方法對(duì)某裂隙巖質(zhì)斜坡進(jìn)行穩(wěn)定性分析。該裂隙巖質(zhì)斜坡特征以陡峻地形為主,坡腳、坡頂?shù)匦蜗鄬?duì)高差86 m,斜坡坡面平均產(chǎn)狀為30°∠65°。坡面多為基巖出露,巖性主要為閃長(zhǎng)巖,表現(xiàn)為塊狀構(gòu)造,主要成分為白色的斜長(zhǎng)石和深色的角閃石。此外,長(zhǎng)期受內(nèi)外動(dòng)力地質(zhì)作用的影響,坡面巖體結(jié)構(gòu)破碎,節(jié)理裂隙十分發(fā)育。

通過(guò)分析測(cè)量所得420條裂隙的跡長(zhǎng)及產(chǎn)狀信息,并采用動(dòng)態(tài)聚類(lèi)法進(jìn)行分析,可將其劃分為三個(gè)優(yōu)勢(shì)組,產(chǎn)狀分別為150°∠71°，89°∠43°和339°∠42°。據(jù)此確定裂隙的概率分布類(lèi)型,見(jiàn)表1。

表1 裂隙優(yōu)勢(shì)組主要參數(shù)匯總

對(duì)現(xiàn)場(chǎng)取樣巖石進(jìn)行室內(nèi)實(shí)驗(yàn),確定巖石黏聚力為24.6 MPa、內(nèi)摩擦角為53°、巖石密度為2.8 g/cm3;對(duì)于裂隙強(qiáng)度參數(shù),由于現(xiàn)場(chǎng)裂隙基本閉合且無(wú)填充,故黏聚力為0,裂隙內(nèi)摩擦角采用基本內(nèi)摩擦角,參照Brideau[14]研究取值為27°。

采用三維裂隙網(wǎng)絡(luò)模擬技術(shù),對(duì)該裂隙巖質(zhì)斜坡巖體內(nèi)裂隙進(jìn)行模擬。最終,模擬一個(gè)尺寸為110 m(X軸)×170 m(Y軸)×90 m(Z軸)的三維裂隙網(wǎng)絡(luò),其可視化模型見(jiàn)圖4。現(xiàn)場(chǎng)調(diào)查該裂隙巖質(zhì)斜坡可能沿臨空方向滑出,因此我們選取垂直于斜坡露頭面走向的鉛直平面為斜坡主滑面,并以此面截切三維裂隙網(wǎng)絡(luò)得到該裂隙巖質(zhì)斜坡主滑面的二維裂隙網(wǎng)絡(luò)(見(jiàn)圖5)。

采用Floyd算法對(duì)該裂隙巖質(zhì)斜坡進(jìn)行危險(xiǎn)滑動(dòng)面搜索時(shí),為保證搜索精度同時(shí)減少計(jì)算量,統(tǒng)一將主滑面二維裂隙網(wǎng)絡(luò)中的每一條裂隙跡線(xiàn)離散成5個(gè)節(jié)點(diǎn)(4等分);將主滑面二維裂隙網(wǎng)絡(luò)的邊界每隔10 m設(shè)置一個(gè)節(jié)點(diǎn)。最終,選取剪入口(60 m,86 m),(70 m,86 m),(80 m,86 m),(90 m,86 m),(100 m,86 m),(110 m,86 m)與剪出口(0 m,0 m),(4.66 m,10 m),(9.32 m,20 m),(13.99 m,30 m),(18.65 m,40 m),(23.32 m,50 m)進(jìn)行分析。組合任一個(gè)剪入與剪出口,采用Floyd算法即可獲得36條危險(xiǎn)滑動(dòng)面(見(jiàn)圖6)。

根據(jù)GB 50330—2013建筑邊坡工程技術(shù)規(guī)范,采用傳遞系數(shù)法對(duì)獲得的36條危險(xiǎn)滑動(dòng)面進(jìn)行天然狀態(tài)下安全系數(shù)計(jì)算,結(jié)果如表2所示。

表2 各危險(xiǎn)滑動(dòng)面安全系數(shù)計(jì)算結(jié)果統(tǒng)計(jì)

由表2可知,不同剪入-剪出口組合之間危險(xiǎn)滑動(dòng)面所對(duì)應(yīng)的安全系數(shù)不同。該裂隙巖質(zhì)斜坡各危險(xiǎn)滑動(dòng)面最小安全系數(shù)為15.51,最大為68.37,均值為39.33。為保證工程安全性,選擇最小的安全系數(shù)作為對(duì)整個(gè)裂隙巖質(zhì)斜坡穩(wěn)定性的評(píng)估,則該裂隙巖質(zhì)斜坡安全系數(shù)為15.51,故天然狀態(tài)下該裂隙巖質(zhì)斜坡整體穩(wěn)定性好。

6 結(jié)論

通過(guò)本研究提出的方法,針對(duì)某裂隙巖質(zhì)斜坡進(jìn)行研究,最終可得到以下結(jié)論:

1)裂隙巖質(zhì)斜坡危險(xiǎn)滑動(dòng)面往往沿著預(yù)先存在的裂隙擴(kuò)展,因此裂隙巖質(zhì)斜坡穩(wěn)定性分析前必須考慮斜坡主滑面裂隙的發(fā)育特征。由于斜坡主滑面裂隙網(wǎng)絡(luò)在發(fā)生破壞前不可見(jiàn),因而根據(jù)現(xiàn)場(chǎng)露頭面裂隙信息模擬生成與現(xiàn)場(chǎng)裂隙統(tǒng)計(jì)相似的數(shù)學(xué)地質(zhì)模型的裂隙網(wǎng)絡(luò)模擬技術(shù)是裂隙巖質(zhì)斜坡穩(wěn)定性分析必不可缺的。

2)經(jīng)分析,裂隙巖質(zhì)斜坡危險(xiǎn)滑動(dòng)面的滑動(dòng)主要沿剪入口與剪出口之間的最短路徑。因而可視剪入口與剪出口之間的最短路徑為裂隙巖質(zhì)斜坡的危險(xiǎn)滑動(dòng)面。Floyd算法可以簡(jiǎn)單快速計(jì)算出任意兩個(gè)節(jié)點(diǎn)之間的最短路徑,因而能很好地被應(yīng)用于裂隙巖質(zhì)斜坡危險(xiǎn)滑動(dòng)面的搜索。

3)Floyd算法搜索得到的危險(xiǎn)滑動(dòng)面為折線(xiàn)形,因而可以采用傳遞系數(shù)法進(jìn)行安全系數(shù)計(jì)算。不同危險(xiǎn)滑動(dòng)面對(duì)應(yīng)的安全系數(shù)不同,選擇最小的安全系數(shù)對(duì)整個(gè)裂隙巖質(zhì)斜坡穩(wěn)定性進(jìn)行評(píng)估。對(duì)于研究案例,天然狀態(tài)下安全系數(shù)為15.51,表明天然狀態(tài)下該裂隙巖質(zhì)斜坡整體穩(wěn)定性好。

4)本文研究案例僅針對(duì)選定剪入-剪出口組合進(jìn)行危險(xiǎn)滑動(dòng)面搜索,計(jì)算所得最小安全系數(shù)與實(shí)際可能存在一定誤差,這個(gè)誤差可通過(guò)增密剪入-剪出口組合來(lái)消除或緩解。