楊雅萍， 張文蓮， 孫曉云

(石家莊鐵道大學(xué)電氣與電子工程學(xué)院， 石家莊 050043)

巖土工程中，邊坡穩(wěn)定性預(yù)測是一個重要的領(lǐng)域。邊坡是一個動態(tài)開放和非線性的復(fù)雜系統(tǒng)，目前研究表明邊坡滑動的誘發(fā)因子之間具有不確定性、離散性、隨機性以及相互之間復(fù)雜的關(guān)聯(lián)作用[1]。在邊坡穩(wěn)定性確定性分析法中，極限平衡法和傳統(tǒng)的數(shù)值分析方法無法準確描述邊坡的非線性特征。近年來人工神經(jīng)網(wǎng)絡(luò)的發(fā)展，及其非線性動態(tài)處理的優(yōu)勢，神經(jīng)網(wǎng)絡(luò)也被逐漸應(yīng)用到邊坡工程實際中來預(yù)測邊坡穩(wěn)定性。李龍等[2]在傳統(tǒng)反向傳播(back propagation,BP)神經(jīng)網(wǎng)絡(luò)中引入附加動量、自適應(yīng)學(xué)習(xí)速率和模擬退火算法，提出了改進的BP神經(jīng)網(wǎng)絡(luò)，對網(wǎng)絡(luò)收斂速度慢和易陷入局部最小值的問題有較好的改進等。牛鵬飛等[3]引入核主成分算法(principal component analysis,KPCA)和Levenberg-Marquardt(L-M)算法建立了主成分分析(principal components analysis, PCA)-LMBP邊坡穩(wěn)定性網(wǎng)絡(luò)預(yù)測模型。臧焜巖等[4]在傳統(tǒng)BP神經(jīng)網(wǎng)絡(luò)中引入遺傳算法(genetic algorithm,GA)，用來更新網(wǎng)絡(luò)的權(quán)值和閾值，建立了GA-BP邊坡穩(wěn)定性預(yù)測模型。

以上3種神經(jīng)網(wǎng)絡(luò)模型輸入均采用Mohr-coulomb(M-C)準則中的內(nèi)聚力、內(nèi)摩擦角等巖土力學(xué)參數(shù)。雖然目前大多數(shù)學(xué)者均采用M-C準則中的力學(xué)參數(shù)，但M-C準則無法解釋巖質(zhì)邊坡的非線性變形破壞現(xiàn)象。改進Hoek-Brown(H-B)準則對巖體固有的非線性破壞特點可以很好地反映，而且該準則已在巖體強度預(yù)測及計算領(lǐng)域得到廣泛應(yīng)用[5]。目前，BP神經(jīng)網(wǎng)絡(luò)隱含層神經(jīng)單元數(shù)目在實際應(yīng)用中不存在一個合理的解析式，仍是以經(jīng)驗公式確定，而BP神經(jīng)網(wǎng)絡(luò)在工程應(yīng)用成敗的關(guān)鍵是其網(wǎng)絡(luò)的構(gòu)造[6]。對于網(wǎng)絡(luò)的構(gòu)造，廣義回歸神經(jīng)網(wǎng)絡(luò)(generalized regression neural network，GRNN)網(wǎng)絡(luò)結(jié)構(gòu)與輸入樣本維數(shù)和數(shù)量有關(guān)，網(wǎng)絡(luò)結(jié)構(gòu)無需根據(jù)經(jīng)驗設(shè)定。

針對上述兩個問題，提出了一種邊坡穩(wěn)定性預(yù)測新方法，即基于廣義H-B準則的巖體力學(xué)參數(shù)，采用交叉變異優(yōu)化的蝙蝠算法，構(gòu)建了改進的BA-GRNN神經(jīng)網(wǎng)絡(luò)模型，對邊坡的穩(wěn)定性和安全性進行預(yù)測。

1 算法原理

1.1 廣義Hoek-Brown準則

在統(tǒng)計分析大量巖體原位試驗成果后，得到了反映巖體產(chǎn)生破壞的半經(jīng)驗準則Hoek-Brown強度準則，該準則主要反映了巖體最大和最小主應(yīng)力的關(guān)系[7]，可表示為

(1)

式(1)中：σ1、σ3分別為最大、最小主應(yīng)力；σci為完整巖石材料單軸抗壓強度；mb、s、a為反映巖體特征的經(jīng)驗參數(shù)，可以由巖體地質(zhì)強指標、巖體的擾動因子和完整巖石經(jīng)驗常數(shù)計算確定。

1.2 廣義回歸神經(jīng)網(wǎng)絡(luò)

廣義回歸神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)以徑向基神經(jīng)網(wǎng)絡(luò)為基礎(chǔ)，具有很強的非線性映射能力和學(xué)習(xí)速度[8]。

GRNN網(wǎng)絡(luò)結(jié)構(gòu)有4層,其中輸入變量由輸入層直接傳遞給模式層，傳遞函數(shù)為

(2)

式(2)中：X=[x1,x2,…,xn]為輸入學(xué)習(xí)樣本；Xi為訓(xùn)練樣本；pi為指數(shù)分子為學(xué)習(xí)樣本X與訓(xùn)練樣本的Xi的歐氏距離平方和；徑向基函數(shù)為高斯函數(shù)，σ為高斯函數(shù)標準差，隱含層位置中基函數(shù)的形狀由σ決定。

輸出層輸出樣本為Y=[y1,y2,…,yn]，輸出層神經(jīng)元輸出為

(3)

式(3)中：Pi為每個模式層神經(jīng)元的輸出；SNj、SD分別為求和層來兩類神經(jīng)元的求和結(jié)果；yij為Yi的第j個元素。

由式(1)可知，σ為GRNN網(wǎng)絡(luò)需要優(yōu)化的唯一參數(shù)，σ的取值會影響網(wǎng)絡(luò)預(yù)測結(jié)果，合理的σ取值，可以增強GRNN網(wǎng)絡(luò)的泛化性，σ的較大取值，使得網(wǎng)絡(luò)預(yù)測值與所有樣本因變量的均值近似相等，預(yù)測效果不好；反之當(dāng)σ取值趨于0，預(yù)測值則會更加準確[8]。

GRNN具有柔性網(wǎng)絡(luò)結(jié)構(gòu)、高度容錯性及魯棒性[9]，針對樣本少和部分樣本數(shù)據(jù)精準度比較差的問題，能得到較好的預(yù)測效果。

1.3 蝙蝠算法

蝙蝠算法類似于粒子群算法，不同的是蝙蝠算法增加了頻率和響度因素，通過頻率變化影響蝙蝠速度，從而改變蝙蝠位置對蝙蝠個體進行尋優(yōu)。目前蝙蝠算法已成為搜索全局最優(yōu)解的有效方法[10]。算法具體闡述如下：隨機初始化虛擬蝙蝠的飛行速度vi、位置xi，脈沖頻率fi、脈沖響度Ai和脈沖率r，蝙蝠發(fā)現(xiàn)獵物時，改變頻率、減小響度和增加脈沖發(fā)射率來更新蝙蝠個體速度和位置，通過適應(yīng)度函數(shù)評價當(dāng)前個體位置的優(yōu)劣，對個體進行優(yōu)勝劣汰，選擇最優(yōu)解。

對進化過程中蝙蝠種群t時刻的搜索脈沖頻率fi、速度vi和位置xi進行更新，分別表示為

fi=fmin+(fmax-fmin)β

(4)

(5)

(6)

式中：β∈[0,1]為均勻分布的隨機數(shù)；X*為當(dāng)前所有蝙蝠的最優(yōu)解；蝙蝠的初始搜索脈沖頻率fi在[fmin,fmax]內(nèi)隨機確定。

對t時刻脈沖音量Ai和脈沖發(fā)射率ri進行更新，可表示為

(7)

(8)

式中：搜索脈沖頻率范圍為[fmin,fmax]；音量的衰減系數(shù)數(shù)為α；最大脈沖音量為A0；最大脈沖率r0；搜索頻率的增強系數(shù)為γ。

合理的音量衰減系數(shù)和搜索頻率可以有效平衡算法尋優(yōu)精度和收斂速度，在算法仿真中需要不斷調(diào)整音量衰減系數(shù)和搜索頻率。與其他粒子群算法相比，模型簡單，參數(shù)較少，收斂速度快已成為蝙蝠算法的突出優(yōu)點。

在算法尋優(yōu)過程中,由于蝙蝠種群個體缺乏變異機制，超級蝙蝠對種群中其他蝙蝠的吸引會使得種群多樣性顯著下降[11]。采用交叉變異算子來豐富蝙蝠種群的多樣性。

2 建立改進BA-GRNN網(wǎng)絡(luò)模型

2.1 輸入輸出變量的確定

研究影響邊坡穩(wěn)定性的因素分為邊坡巖土力學(xué)參數(shù)和邊坡幾何形態(tài)兩方面。對于巖土力學(xué)參數(shù)的影響，基于廣義Hoek-Brown經(jīng)驗準則[12-17]，采用巖體重度rG、地質(zhì)強度指標GSI、擾動因子D、單軸抗壓強度σci、巖體軟硬程度參數(shù)mi表征，對于邊坡幾何形態(tài)的影響，采用坡高H，坡角αs表征，將以上7個因素作為BA-GRNN網(wǎng)絡(luò)模型的輸入，輸出為邊坡安全系數(shù)和邊坡狀態(tài)[18-23]。搜集120個經(jīng)典實例作為邊坡樣本[24-28],最后10個作為預(yù)測樣本。對于網(wǎng)絡(luò)不同的輸入?yún)?shù)，其數(shù)值上相差很大，訓(xùn)練前對輸入數(shù)據(jù)進行歸一化處理可以提升網(wǎng)絡(luò)收斂速度，最后對預(yù)測結(jié)果網(wǎng)絡(luò)輸出安全系數(shù)和邊坡狀態(tài)進行反歸一化處理，得到與實際值相同量綱的值。

2.2 改進BA-GRNN模型建立

采用廣義H-B準則中的巖體力學(xué)參數(shù)，采用交叉變異算子改進的BA算法優(yōu)化的GRNN神經(jīng)網(wǎng)絡(luò)的光滑因子，建立改進BA-GRNN網(wǎng)絡(luò)，對邊坡的穩(wěn)定性和安全性進行預(yù)測，改進后網(wǎng)絡(luò)拓撲結(jié)構(gòu)如圖1所示。

圖1 BA-GRNN網(wǎng)絡(luò)拓撲圖

輸入層的7個神經(jīng)元分別對應(yīng)于2.1節(jié)中的7個影響因素，模式層神經(jīng)元數(shù)目對應(yīng)于2.1節(jié)中的110組訓(xùn)練樣本。在增加蝙蝠種群多樣性方面，分別采用最常用和最基本的單點交叉算子和基本位變異算子來優(yōu)化蝙蝠算法。單點交叉(one-point crossover)算子是指在個體兩兩配對的父代染色體編碼串中隨機設(shè)置交叉點A，然后在A處相互交換兩個父代的Ⅱ段染色體,交叉運算過程如圖2所示。

由圖2知，相互配對的父代染色體每一次交叉運算產(chǎn)生數(shù)量相同的子代染色體，而且根據(jù)交叉點選取的不同，兩個子代染色體來自兩個父代染色體基因座長度不同。

圖2 交叉算子示意圖

基本位變異(simple mutation)算子是對個體編碼基因座上的每一位隨機生成[0,1]內(nèi)的隨機數(shù)r′，依變異概率pm與r′的大小，判斷該位是否進行變異運算，對判斷pm>r′的基因座上的基因值作為變異運算。變異過程如圖3所示。

圖3 變異算子示意圖

改進BA-GRNN神經(jīng)網(wǎng)絡(luò)算法基本實現(xiàn)步驟如下。

步驟1學(xué)習(xí)樣本和預(yù)測樣本數(shù)據(jù)歸一化，輸入到改進后的網(wǎng)絡(luò)。

步驟2用二分法確定光滑因子初始取值范圍[a,b]。

步驟3初始化蝙蝠算法參數(shù)A0、r0、α、γ和最大迭代次數(shù)Nmax;[a,b]內(nèi)初始化蝙蝠種群；初始化交叉概率pc和變異概率pm。

步驟4計算初始化種群適應(yīng)度，記錄初始最優(yōu)解適應(yīng)度函數(shù)值fmin和其位置I，將位置I下的蝙蝠個體值賦給Best。

步驟5選擇新群體進行兩兩配對，根據(jù)交叉概率pc在相互配對的個體交叉點A點處相互交換父代個體的部分染色體，產(chǎn)生與父代數(shù)量相同的子代染色體。

步驟6對個體的每個基因座，根據(jù)變異概率pm對判斷為變異位的基因值做取反運算，如原基因值為0的取為1，反之原基因值為0的取為1，以此產(chǎn)生新的子代個體。

步驟7按式(3)～式(5)更新上述產(chǎn)生的個體的位置和速度。

步驟8若生成的均勻分布隨機數(shù)Rrand>r′,則給當(dāng)前最優(yōu)解一個隨機擾動，并在擾動后對新解進行越界處理，即對于超出給定范圍的值將其限制在邊界值上。

步驟9若生成均勻分布隨機Rrand>Ai，且適應(yīng)度函數(shù)f(xi)>f(x*)，則接受式(8)產(chǎn)生的解, 然后按式(6)、式(7)更新種群的響度和脈沖發(fā)射率。

步驟10計算所有蝙蝠個體的適應(yīng)度值并對其進行排序,找到的適應(yīng)度值最小值對應(yīng)的個體即為當(dāng)前最優(yōu)解，同時記錄最優(yōu)解對應(yīng)的位置。

步驟11若上一步中產(chǎn)生的局部最優(yōu)解適應(yīng)度Fnew

步驟12重復(fù)步驟(4)～步驟(8)直至滿足設(shè)定條件，本文采用達到最大迭代次數(shù)作為算法結(jié)束的標志。

步驟13輸出全局最優(yōu)解Best和最優(yōu)解位置。

步驟14上一步輸出的最優(yōu)解即為GRNN網(wǎng)絡(luò)的光滑因子，將10組預(yù)測樣本輸入到優(yōu)化好的網(wǎng)絡(luò)，即可輸出預(yù)測結(jié)果。

改進BA-GRNN網(wǎng)絡(luò)算法流程圖如圖4所示。

圖4 改進BA-GRNN網(wǎng)絡(luò)算法流程圖

3 實例分析與結(jié)果對比

基于MATLAB2020平臺，網(wǎng)絡(luò)對所選110組樣本進行訓(xùn)練后，在迭代300次后能得到較好預(yù)測結(jié)果。改進BA算法的初始化參數(shù)取值如表1所示。蝙蝠算法最優(yōu)個體值為0.178，最優(yōu)個體適應(yīng)度函數(shù)值為1.207 4×10-3。圖5為適應(yīng)度函數(shù)值的變化曲線。

終止代數(shù)=300

表1 改進BA算法初始參數(shù)

GRNN網(wǎng)絡(luò)光滑因子σ取值為0.178。為了進一步分析改進后BA-GRNN網(wǎng)絡(luò)的預(yù)測效果，與BA-GRNN神經(jīng)網(wǎng)絡(luò)和K-交叉驗證選取GRNN光滑因子網(wǎng)絡(luò)(Kfold-GRNN)兩種模型預(yù)測結(jié)果進行比較，表2對比了3種網(wǎng)絡(luò)預(yù)測結(jié)果。

由表2可知，3種網(wǎng)絡(luò)對邊坡狀態(tài)預(yù)測結(jié)果完全相同。由圖6可以看出，對于邊坡安全系數(shù)的預(yù)測，改進后的BA-GRNN網(wǎng)絡(luò)模型對安全系數(shù)真實值的擬合效果更好。由3種網(wǎng)絡(luò)邊坡安全系數(shù)預(yù)測值得平均絕對百分比誤差(mean absolute percentage error,MAPE)和均方誤差(mean square error，MSE)，如表3所示，進一步得到改進后的BA-GRNN網(wǎng)絡(luò)MAPE和MSE均小于前兩種網(wǎng)絡(luò)。由上述知，改進BA-GRNN網(wǎng)絡(luò)對邊坡穩(wěn)定性的預(yù)測效果較好。

表2 3種網(wǎng)絡(luò)模型預(yù)測結(jié)果對比

圖6 安全系數(shù)對比

表3 3種網(wǎng)絡(luò)模型預(yù)測結(jié)果及誤差對比

4 算例驗證

選用經(jīng)典邊坡算例，該邊坡為天然均質(zhì)邊坡，邊坡巖體強度參數(shù)及幾何參數(shù)如表4所示。將邊坡算例在FLAC3D中進行建模，模型只考慮重力作用，底部邊界為固定約束，水平約束施加于兩側(cè)邊界，坡面設(shè)為自由邊界[29-31]，由此所建模型如圖7所示。

表4 邊坡巖體強度參數(shù)及幾何參數(shù)

針對圖7模型，利用FLAC3D中model factor-of-safety命令計算邊坡安全系數(shù)，得安全系數(shù)為1.15，通過邊坡計臨界破壞狀態(tài)的最大剪應(yīng)變云圖可得邊坡臨界破壞面，如圖8所示。

圖7 邊坡有限元模型

圖8 邊坡臨界破壞最大剪應(yīng)變云圖

通過本文網(wǎng)絡(luò)對經(jīng)典算例計算，得到安全系數(shù)和邊坡狀態(tài)的值為1.158 0和1.000 0，通過與FLAC3D的計算結(jié)果進行對比，可得本文方法與FLAC3D建模計算結(jié)果相近，網(wǎng)絡(luò)預(yù)測安全系數(shù)MAPE為0.70%?？梢妼τ趯嶋H邊坡案例，改進后的BA-GRNN網(wǎng)絡(luò)對安全系數(shù)和狀態(tài)的預(yù)測是準確的。

5 結(jié)論

(1)結(jié)合廣義Hoek-Brown準則，將Hoek-Brown準則所涉及的巖體參數(shù)作為穩(wěn)定性影響因素，對巖質(zhì)邊坡進行穩(wěn)定性預(yù)測。分析結(jié)果顯示，所選擇影響因素對巖質(zhì)邊坡有更好的適用性。

(2)利用改進的BA算法優(yōu)化的GRNN網(wǎng)絡(luò)對邊坡的安全系數(shù)和狀態(tài)進行預(yù)測，通過與BA-GRNN和Kfold-GRNN神經(jīng)網(wǎng)絡(luò)邊坡預(yù)測模型對比，改進BA-GRNN網(wǎng)絡(luò)安全系數(shù)的平均絕對百分比誤差在8%以下，均方誤差只有0.17%，而且對實際邊坡的安全系數(shù)和是否失穩(wěn)可以進行準確的判斷。

(3)在邊坡穩(wěn)定性的確定分析中，傳統(tǒng)的數(shù)值分析法不僅計算煩瑣，計算量大，時間長，而且對經(jīng)驗的依賴性較強，所使用的方法可以很好地解決這些問題。