馬金城 石 崇 張一平 張文浩 王樂(lè)榮
(1.河海大學(xué) 巖土力學(xué)與堤壩工程教育部重點(diǎn)實(shí)驗(yàn)室,南京 210098;2.河海大學(xué) 巖土工程科學(xué)研究所,南京 210098)
邊坡穩(wěn)定性分析是巖土工程十分重要的研究?jī)?nèi)容,涉及水利水電、公路、鐵路、礦山工程等諸多領(lǐng)域,是判斷邊坡穩(wěn)定性,是否需要加固及采取何種防護(hù)措施的主要依據(jù).在進(jìn)行邊坡穩(wěn)定分析時(shí),關(guān)鍵是如何得出滑面并利用滑動(dòng)力、抗滑力計(jì)算安全系數(shù)、優(yōu)化得出最危險(xiǎn)滑面.因此,對(duì)于危險(xiǎn)滑面搜索方法進(jìn)行理論和應(yīng)用研究非常有價(jià)值.
在常規(guī)的條分法邊坡穩(wěn)定計(jì)算中,如瑞典條分法和畢肖普法假定滑坡體為不變形的剛體,利用自重計(jì)算滑面上滑動(dòng)力與抗滑力,再利用抗滑力矩與滑動(dòng)力矩之比來(lái)定義穩(wěn)定性系數(shù)[1].而自然巖土體處于彈塑性變形狀態(tài),達(dá)到屈服時(shí)應(yīng)力可向周圍傳遞,因此條分法得到的各條分底部抗滑力、滑動(dòng)力往往與實(shí)際偏差很大.而數(shù)值模擬方法綜合考慮物理方程、平衡方程、幾何邊界條件,因此可以逼真地得到巖土體應(yīng)力場(chǎng),但不能直接給出穩(wěn)定性系數(shù).如果將這兩者有效結(jié)合起來(lái),利用數(shù)值模擬得到真實(shí)的應(yīng)力場(chǎng),采用矢量和法可定義穩(wěn)定性系數(shù),對(duì)邊坡穩(wěn)定性計(jì)算有很好參考價(jià)值.因此近些年很多學(xué)者致力于探索危險(xiǎn)滑面搜索方法的工作,如遺傳算法[2],蟻群算法[3],粒子群算法[4],模擬退火法[5]等.這些算法在邊坡滑面搜索中得以運(yùn)用,與早期使用的窮舉法、二分法相比,上述方法的精度和收斂性提高了很多,對(duì)邊坡穩(wěn)定研究和實(shí)際工程應(yīng)用起到了積極的作用,但蟻群算法容易出現(xiàn)停滯等現(xiàn)象[6],粒子群算法會(huì)過(guò)早陷入局部最優(yōu)等缺點(diǎn)[7],模擬退火算法和遺傳算法對(duì)參數(shù)選擇要求較高,如若參數(shù)選擇不當(dāng)會(huì)使搜索陷入局部極值或搜索結(jié)果不令人滿意[8].
狼群算法[9]是2011年提出的一種新型的智能優(yōu)化算法,該算法在全局搜索和局部開發(fā)能力上具有良好的適用性,被廣泛地運(yùn)用到路徑規(guī)劃等實(shí)際問(wèn)題中,具有魯棒性好,全局尋優(yōu)能力強(qiáng),收斂速度快等優(yōu)點(diǎn),但目前尚未運(yùn)用到邊坡穩(wěn)定性分析之中.本文基于傳統(tǒng)的狼群算法,嘗試結(jié)合邊坡穩(wěn)定性計(jì)算特點(diǎn)進(jìn)行優(yōu)化改進(jìn),建立狼群滑面位置搜索及穩(wěn)定計(jì)算模型,并用已有成果驗(yàn)證其可行性和準(zhǔn)確性.
對(duì)一個(gè)二維邊坡中潛在的滑面,計(jì)算穩(wěn)定時(shí)一般如圖1所示將滑面劃分為若干段,不同段相交點(diǎn)稱為一個(gè)狀態(tài)點(diǎn),從起始點(diǎn)到終止點(diǎn)通過(guò)狀態(tài)點(diǎn)首尾相連可形成一條路徑,每條路徑都是一條潛在的滑動(dòng)面.
圖1 滑面受力示意圖及單元體圖
設(shè)滑動(dòng)方向向右,逆時(shí)針為下滑力正方向,則其中某一段上的抗剪強(qiáng)度τf和切向應(yīng)力τ的計(jì)算公式如下:
第i段的平均滑動(dòng)力:
第i段的平均抗滑力:
式中:σx、σy、τxy分別為單元形心沿x,y方向的正應(yīng)力和切應(yīng)力;φ和c分別為材料的內(nèi)摩擦角和內(nèi)聚力;σα為斜面的法向應(yīng)力;α為斜面與水平面的夾角.
利用數(shù)值模擬的方法獲取邊坡真實(shí)應(yīng)力場(chǎng)后,計(jì)算出每一段上的滑動(dòng)力和抗滑力,然后采用矢量和法定義安全系數(shù),安全系數(shù)定義為:
式中:A和B分別為邊坡起始點(diǎn)和終止點(diǎn);n為滑動(dòng)面所經(jīng)過(guò)的段的數(shù)量;τ為沿滑動(dòng)面滑動(dòng)方向的剪應(yīng)力;dL為沿滑動(dòng)面每一段的長(zhǎng)度增量;τi為從第i-1個(gè)節(jié)點(diǎn)到第i個(gè)節(jié)點(diǎn)路徑上的剪應(yīng)力;τfi為從第i-1段到第i段路徑上的抗剪強(qiáng)度;Δd i為從第i-1階段到第i階段路徑的長(zhǎng)度[10].
對(duì)于已知滑面邊坡,采用上述極限平衡穩(wěn)定性分析方法可直接求出安全系數(shù),但對(duì)于滑面未知的邊坡,如何確定滑動(dòng)面上點(diǎn)的位置,找到最危險(xiǎn)滑面是非常困難的[11].從數(shù)學(xué)問(wèn)題看,確定滑面實(shí)際上是圖形路徑的優(yōu)化問(wèn)題.近年來(lái),狼群算法作為一種生物仿生算法被提出,可用于圖形最優(yōu)化問(wèn)題的求解.因此,本文針對(duì)二維邊坡的網(wǎng)格,將滑出滑入點(diǎn)當(dāng)做頭狼與母狼的位置,某點(diǎn)的應(yīng)力當(dāng)作當(dāng)前位置的獵物氣味濃度,網(wǎng)格節(jié)點(diǎn)當(dāng)作猛狼,網(wǎng)格內(nèi)部的某些點(diǎn)為探狼,通過(guò)不同種類狼之間的信息交互、“勝者為王”的頭狼產(chǎn)生機(jī)制和“強(qiáng)者生存”的狼群更新機(jī)制來(lái)尋求全局路徑的最優(yōu)解,可以用于解決邊坡最危險(xiǎn)滑面的搜索.
經(jīng)典狼群算法分為如下幾個(gè)步驟:
1)頭狼生成:在某一待尋優(yōu)區(qū)域內(nèi),將位置最優(yōu)的狼視為頭狼.頭狼是狼群的核心,負(fù)責(zé)狼群的指揮工作[12].
2)游走行為:選取除頭狼外最佳的數(shù)匹狼作為探狼,執(zhí)行游走行為,根據(jù)某種計(jì)算規(guī)則得到探狼處氣味濃度,然后和當(dāng)前頭狼氣味濃度進(jìn)行比較,選擇二者中較優(yōu)者為頭狼.
3)召喚行為:當(dāng)頭狼發(fā)出召喚行為后,猛狼以較大步長(zhǎng)向頭狼位置逼近.如果在奔襲過(guò)程中,猛狼位置優(yōu)于頭狼,則猛狼代替頭狼發(fā)出召喚行為,否則猛狼繼續(xù)奔襲,直至距離頭狼位置小于圍攻距離,此時(shí)轉(zhuǎn)化為圍攻行為[13].
4)圍攻行為:當(dāng)猛狼奔襲到距頭狼較近位置時(shí),猛狼聯(lián)合探狼對(duì)獵物進(jìn)行圍攻.
5)“強(qiáng)者生存”的狼群更新機(jī)制:剔除掉較差的狼,只有適應(yīng)性較好的狼才能生存下來(lái).
1)模型的建立:根據(jù)邊坡形狀及地層分布狀況,結(jié)合邊界條件,對(duì)邊坡進(jìn)行合理的有限單元網(wǎng)格劃分,單元格數(shù)量由上述因素及所需結(jié)果精度共同確定.得到有限單元網(wǎng)格后,計(jì)算邊坡的應(yīng)力場(chǎng),得到每個(gè)節(jié)點(diǎn)處的應(yīng)力和坐標(biāo)[14].將所得應(yīng)力場(chǎng)視為草原平面,也即狼群活動(dòng)的平面,將節(jié)點(diǎn)處的應(yīng)力值視為當(dāng)前位置的獵物氣味濃度.
2)頭狼生成:在狼群中,頭狼和母狼處于支配地位,頭狼負(fù)責(zé)狼群的指揮,母狼負(fù)責(zé)配合頭狼工作.頭狼和母狼對(duì)整個(gè)狼群的優(yōu)劣起決定作用.如圖2 所示,頭狼和母狼在可能滑入點(diǎn)與滑出點(diǎn)的某一區(qū)域內(nèi),初始頭狼與母狼分別為該區(qū)域內(nèi)橫坐標(biāo)值最小的點(diǎn).算法中引入母狼的原因是為配合頭狼工作,兩者共同把握整個(gè)狼群的位置,避免陷入局部最優(yōu).
圖2 可能圓心位置與滑出滑入點(diǎn)切線斜率
3)狼群區(qū)域劃分:如圖3 所示,以滑入點(diǎn)(X0,Y0)為坐標(biāo)原點(diǎn)建立搜索區(qū)域坐標(biāo)系,用相同間距為x d=(X n-X0)/n的縱向線將待尋優(yōu)區(qū)域劃分成n部分,其中X n為滑出點(diǎn)橫坐標(biāo).將每條縱向線鄰域視為一個(gè)搜索區(qū)域,分別記為L(zhǎng)Q1,LQ2,…,LQn.將各縱向線橫坐標(biāo)依次記為X1=(X n-X0)/n,X2=2(X n-X0)/n,…,X n=X n.在給定滑出滑入點(diǎn)位置后,需確定下一個(gè)點(diǎn)(X ML1,Y ML1)的搜索位置,潘家錚等人通過(guò)研究表明,對(duì)于簡(jiǎn)單均質(zhì)土坡,作邊坡中點(diǎn)的垂線和法線,以坡面中點(diǎn)為圓心,以1/4坡長(zhǎng)和5/4坡長(zhǎng)為半徑畫同心圓,最危險(xiǎn)滑弧的圓心在如圖2所示陰影部分ABDC范圍內(nèi).據(jù)經(jīng)驗(yàn)可知,圓心縱坐標(biāo)越小,滑出點(diǎn)和滑入點(diǎn)與圓弧相切切線斜率越大(滑出點(diǎn)和滑入點(diǎn)切線的斜率僅用于確定滑出滑入點(diǎn)的下一個(gè)點(diǎn)搜索范圍,并不對(duì)滑動(dòng)面是否為圓弧面做限制).設(shè)滑出點(diǎn)處的切線最大斜率絕對(duì)值為k1,滑入點(diǎn)處最大斜率絕對(duì)值為k2.
圖3 狼群搜索區(qū)域示意圖
如圖4(a)所示,LQ1的搜索范圍為:
如圖4(b)所示,在LQ i中,橫坐標(biāo)在i(X n-X0)/n的左右鄰域內(nèi)波動(dòng).
圖4 狼群搜索區(qū)域及探索步長(zhǎng)
由于滑動(dòng)面點(diǎn)為一下凹曲線,因此點(diǎn)(X0,Y0)與點(diǎn)(X ML1,Y ML1)所在直線為L(zhǎng)Q2搜索區(qū)域的下限,X=X1-δ,X=X1+δ是X=X1左右鄰域.點(diǎn)(X0,Y0)與點(diǎn)(X n,Y n)所在直線為L(zhǎng)Q2搜索區(qū)域上限.也即LQ2的搜索范圍為:
根據(jù)上述規(guī)則,由前一搜索區(qū)域所得猛狼位置來(lái)推算下一狼群搜索區(qū)域,直至最后一匹猛狼產(chǎn)生.這樣做可以有效縮小搜索區(qū)域的范圍,提高搜尋效率,減少搜索時(shí)間.
4)游走行為:探狼從每個(gè)待尋優(yōu)區(qū)域最上側(cè)的點(diǎn)開始以固定步長(zhǎng)y d向下游走探索.處于頭狼和母狼附近的探狼分別以固定步長(zhǎng)x d水平向右游走探索.
5)召喚行為:當(dāng)頭狼和母狼及探狼已處于各自狩獵位置后,頭狼發(fā)出召喚行為,召喚探狼周圍某一范圍內(nèi)猛狼.視探狼所處單元格內(nèi)部為召喚區(qū)域.處于單元格節(jié)點(diǎn)處的猛狼在接收到召喚指令后,向附近區(qū)域探狼奔襲,通過(guò)召喚區(qū)域上猛狼處氣味濃度綜合確定探狼處氣味濃度,也即由已知多邊形節(jié)點(diǎn)的坐標(biāo)和應(yīng)力通過(guò)機(jī)器學(xué)習(xí)的方法確定待求點(diǎn)的應(yīng)力狀態(tài)(具體確定過(guò)程見(8)).
6)圍攻行為:可能最優(yōu)路徑上的頭狼和猛狼經(jīng)上述方式產(chǎn)生后,對(duì)獵物展開圍攻,形成一狼群.
7)“強(qiáng)者生存”“優(yōu)勝劣汰”的狼群更新機(jī)制:將適應(yīng)度D作為目標(biāo)函數(shù),則狼群算法可將求最小安全系數(shù)的問(wèn)題轉(zhuǎn)換為求D最大值問(wèn)題,目標(biāo)函數(shù)定義為:
適應(yīng)度D為狼群對(duì)草原適應(yīng)性的好壞,D越大,則狼群越優(yōu).在完成每次可能最優(yōu)路徑適應(yīng)度D i計(jì)算后,將其適應(yīng)度D i與當(dāng)前適應(yīng)度最好的Dlead進(jìn)行比較,直至最后一條路徑計(jì)算結(jié)束,將適應(yīng)度最好的狼群視為優(yōu)勝狼群,也即:
其中:Dmax為全局最優(yōu)解,其對(duì)應(yīng)的各種路徑為最適宜路徑.
8)機(jī)器學(xué)習(xí)方法確定猛狼處氣味濃度:在Python 中使用sklearn庫(kù)中的Gradient Boosting Re-gressor模塊,將給定數(shù)據(jù)(結(jié)點(diǎn)坐標(biāo)和應(yīng)力值)分為訓(xùn)練集和測(cè)試集[15](將部分已知節(jié)點(diǎn)應(yīng)力和坐標(biāo)作為訓(xùn)練集,其余節(jié)點(diǎn)坐標(biāo)和應(yīng)力作為測(cè)試集,在訓(xùn)練集中找尋規(guī)律,在測(cè)試集中驗(yàn)證規(guī)律[16].可采用L1范數(shù)損失函數(shù)(loss function)來(lái)估量模型的預(yù)測(cè)值f(x)與真實(shí)值Y的不一致程度,其值的大小可以用來(lái)評(píng)價(jià)預(yù)測(cè)的準(zhǔn)確性,該函數(shù)把目標(biāo)值(Y)與估計(jì)值(f(x))的絕對(duì)差值的總和最小化[17]:
相比于傳統(tǒng)的插值法求中間部分應(yīng)力,機(jī)器學(xué)習(xí)方法可以很直觀地通過(guò)S值的大小來(lái)判斷對(duì)某點(diǎn)預(yù)測(cè)的準(zhǔn)確性.
根據(jù)以上原理,形成改進(jìn)狼群算法邊坡滑面搜索方法,計(jì)算流程圖如圖5所示.
圖5 改進(jìn)狼群算法流程圖
本文采用python語(yǔ)言編制程序,引入一經(jīng)典算例來(lái)驗(yàn)證狼群算法在邊坡滑面搜索以及計(jì)算安全系數(shù)的準(zhǔn)確性.該算例為澳大利亞計(jì)算機(jī)協(xié)會(huì)于1987年委托BDonald和PGiam 設(shè)計(jì)ACADS邊坡穩(wěn)定性分析程序的考核題,該題的參考滑動(dòng)面如圖6(a)所示,推薦答案的安全系數(shù)為1.00[18],具體參數(shù)為:容重γ=20 k N/m3,泊松比μ=0.25,黏聚力c=3.0 kPa,內(nèi)摩擦角φ=19.6°.由于該算例已經(jīng)過(guò)大量學(xué)者的引用和對(duì)比分析驗(yàn)證,因此其結(jié)果視為標(biāo)準(zhǔn)值.
基于改進(jìn)狼群算法理論的邊坡臨界滑動(dòng)面搜索方法搜索出最危險(xiǎn)滑動(dòng)面如圖6(b)所示,對(duì)應(yīng)的安全系數(shù)為1.03675,計(jì)算所得安全系數(shù)與推薦值1.00誤差約3.6%,采用本文方法得到的臨界滑動(dòng)面與推薦滑動(dòng)面基本一致,這表明基于改進(jìn)狼群算法的邊坡穩(wěn)定性分析是可行的,其值偏高的原因是數(shù)值應(yīng)力結(jié)果比剛體極限平衡協(xié)調(diào),因此有更高的承載能力.
圖6 邊坡滑面(單位:m)
溪洛渡水電站左岸谷肩堆積體的典型地質(zhì)剖面如圖7(a)所示.本文采用flac-3D 軟件建立如圖7(b)所示數(shù)值模型,選用摩爾-庫(kù)侖的本構(gòu)模型,以及表1中巖土物理力學(xué)參數(shù)進(jìn)行數(shù)值計(jì)算,得到其應(yīng)力場(chǎng)后,采用Python語(yǔ)言編寫程序,運(yùn)用改進(jìn)狼群算法搜索其最危險(xiǎn)滑面,并計(jì)算其安全系數(shù).
表1 巖土物理力學(xué)參數(shù)
如圖7(a)所示,IN13鉆孔傾斜儀與剖面較為接近,圖7(c)表明,溪洛渡水電站左岸谷肩堆積體的最危險(xiǎn)滑面位于滑坡堆積體與宣威組上部交界處.這與歷史監(jiān)測(cè)資料所揭示的滑動(dòng)面位置一致.文獻(xiàn)[19]在自然工況下采用改進(jìn)動(dòng)態(tài)規(guī)劃算法,圓弧形滑面計(jì)算所得安全系數(shù)為1.100,而在相同工況下采用狼群算法得到的安全系數(shù)為1.044,如圖7(d)所示.經(jīng)對(duì)比滑面可知,狼群算法搜索得到的滑面可保證更多比例滑面位于宣威組,而現(xiàn)場(chǎng)勘查表明該邊坡處于蠕滑狀態(tài),自然工況安全系數(shù)為1.00~1.05之間,計(jì)算結(jié)果接近蠕滑狀態(tài)自然工況下安全系數(shù)的上限1.05,故采用本方法得到的結(jié)果更加接近于實(shí)際.
圖7 工程案例計(jì)算結(jié)果
基于數(shù)值計(jì)算得到的邊坡應(yīng)力場(chǎng),提出一種基于改進(jìn)狼群算法的邊坡臨界滑動(dòng)面搜索方法,為邊坡穩(wěn)定分析提供一條新的有效思路.得到主要結(jié)論如下:
1)改進(jìn)狼群算法是利用數(shù)值模擬手段獲得真實(shí)的應(yīng)力場(chǎng)后,采用矢量和法來(lái)定義安全系數(shù).在邊坡穩(wěn)定性分析過(guò)程中可得到較為準(zhǔn)確的邊坡最危險(xiǎn)滑動(dòng)面和最小安全系數(shù).
2)在運(yùn)用改進(jìn)狼群算法的過(guò)程中,用機(jī)器學(xué)習(xí)方法代替?zhèn)鹘y(tǒng)的插值法來(lái)計(jì)算指定坐標(biāo)點(diǎn)處的應(yīng)力,可以通過(guò)損失函數(shù)直觀地看出計(jì)算結(jié)果的準(zhǔn)確性.
3)通過(guò)對(duì)經(jīng)典邊坡算例的分析,本文所提方法的計(jì)算結(jié)果與典型算例參考答案具有較好的一致性,證明了該算法的合理性.