王秀菊， 張金龍， 戴 薇， 石 崇

（1.南京交通職業(yè)技術(shù)學(xué)院，南京 211188； 2.河海大學(xué)巖土工程科學(xué)研究所，南京 210098；3.中國(guó)電建集團(tuán)中南勘測(cè)設(shè)計(jì)研究院有限公司，長(zhǎng)沙 410014）

破碎圍巖地下洞室群工程中往往存在復(fù)雜的破碎帶及不連續(xù)面，對(duì)工程的建設(shè)有重要影響，但巖體地下工程圍巖具有非均勻性、非線性、非連續(xù)性特點(diǎn)，如何對(duì)圍巖力學(xué)參數(shù)取值是進(jìn)行變形預(yù)測(cè)與穩(wěn)定分析的關(guān)鍵之一，也是困擾巖土工程領(lǐng)域的一大難題. 近年來(lái)，隨著位移反分析方法的進(jìn)步，采用智能反演算法結(jié)合數(shù)值模擬方法，應(yīng)用計(jì)算機(jī)技術(shù)與巖土工程理論和工程實(shí)際相結(jié)合的方法來(lái)獲取圍巖變形參數(shù)，成為力學(xué)參數(shù)獲取的重要途徑［1-5］.

在圍巖變形反分析方面，許傳華等［6］結(jié)合索風(fēng)營(yíng)工程現(xiàn)場(chǎng)實(shí)際，建立了支持向量機(jī)和模擬退火位移反分析模型，并反分析其力學(xué)參數(shù). 聶衛(wèi)平等［7］綜合適當(dāng)?shù)母倪M(jìn)粒子群算法和支持向量機(jī)，進(jìn)行糯扎渡水電站的圍巖力學(xué)參數(shù)和初始地應(yīng)力場(chǎng)反演，獲得了良好的效果. 文輝輝等［8］利用BP神經(jīng)網(wǎng)絡(luò)和正交試驗(yàn)方法，反演了隧道圍巖物理力學(xué)參數(shù). 袁振華等［9］提出巖土力學(xué)參數(shù)GA算法與有限元耦合應(yīng)力反分析模型，該模型采用應(yīng)力反分析方法，克服了常規(guī)算法反演過(guò)程中易陷入局部極優(yōu)的難題. 周喻等［10］基于PFC3D和BP神經(jīng)網(wǎng)絡(luò)方法建立計(jì)算模型并反演出巖土體細(xì)觀力學(xué)參數(shù). 陳志超等［11］結(jié)合地下洞室裂隙化花崗巖實(shí)例，建立了BP神經(jīng)網(wǎng)絡(luò)和DEC離散元方法相結(jié)合的計(jì)算模型，反演出了區(qū)域應(yīng)力場(chǎng)值. 李澄清等［12］采用開(kāi)源的顆粒離散元程序LMGC開(kāi)展了土體雙軸壓縮數(shù)值試驗(yàn)，建立了BP人工神經(jīng)網(wǎng)絡(luò)反演系統(tǒng)，也獲得了較為精確的變形參數(shù).

本文依托某電站地下洞室群工程，通過(guò)螢火蟲(chóng)智能分析方法結(jié)合數(shù)值模擬模型進(jìn)行巖體參數(shù)反分析，并利用反演參數(shù)對(duì)典型地下洞室圍巖的穩(wěn)定性進(jìn)行分析與探討.

1 螢火蟲(chóng)算法-最小二乘支持向量機(jī)反分析方法研究

1.1 螢火蟲(chóng)算法

螢火蟲(chóng)算法（Firefly Algorithm，F(xiàn)A）［13-15］是一種新興的群智能優(yōu)化技術(shù)，它通過(guò)模擬自然界中螢火蟲(chóng)個(gè)體根據(jù)發(fā)光亮度尋找其他個(gè)體的特性，以發(fā)光亮度為依據(jù)，模擬螢火蟲(chóng)算法在模擬螢火蟲(chóng)個(gè)體之間相互吸引、相互靠近的過(guò)程，以達(dá)到尋優(yōu)的目的. 基于螢火蟲(chóng)算法的最小搜索模型有三個(gè)假設(shè)：

①螢火蟲(chóng)個(gè)體之間沒(méi)有任何區(qū)別，只考慮發(fā)光的亮度，亮度決定螢火蟲(chóng)個(gè)體的吸引力，亮度低的被亮度高的螢火蟲(chóng)吸引；②假如一定區(qū)域內(nèi)的某個(gè)螢火蟲(chóng)個(gè)體為局部最優(yōu)，則進(jìn)行隨機(jī)移動(dòng)；③螢火蟲(chóng)個(gè)體的適應(yīng)值函數(shù)決定其亮度.

上面三條件決定的螢火蟲(chóng)算法，螢火蟲(chóng)個(gè)體根據(jù)其他螢火蟲(chóng)亮度而發(fā)生追逐的行為即為最優(yōu)解的求解過(guò)程. 在FA算法中，螢火蟲(chóng)局部最優(yōu)搜索的關(guān)鍵是螢火蟲(chóng)的發(fā)光亮度和相互吸引度［16-17］. 假設(shè)所有個(gè)體的初始值相同，個(gè)體的亮度決定該個(gè)體對(duì)其他個(gè)體的吸引度，個(gè)體的目標(biāo)函數(shù)值決定個(gè)體的發(fā)光亮度，在固定的種群范圍內(nèi)，當(dāng)前最優(yōu)解即為種群中亮度最高的個(gè)體，為防止陷入局部最優(yōu)，處于最優(yōu)解的個(gè)體再在一定范圍內(nèi)隨機(jī)移動(dòng)，其他個(gè)體會(huì)根據(jù)個(gè)體間的吸引度向吸引度大的個(gè)體靠近，從而找出全局最優(yōu)值.

1.2 支持向量機(jī)算法原理

支持向量機(jī)算法［18］以統(tǒng)計(jì)學(xué)習(xí)理論為基礎(chǔ)，根據(jù)有限學(xué)習(xí)樣本構(gòu)建復(fù)雜模型，從而獲得更好的推廣能力. 最小二乘支持向量機(jī)（Least squares support vector machine，簡(jiǎn)稱LSSVM）是標(biāo)準(zhǔn)支持向量機(jī)的一種擴(kuò)展［19］，它將二次規(guī)劃問(wèn)題轉(zhuǎn)化為線性方程組求解，降低了計(jì)算的復(fù)雜性，提高了求解速度. LSSVM在訓(xùn)練的過(guò)程中，采取等式約束、經(jīng)驗(yàn)風(fēng)險(xiǎn)偏差的二次方和線性方程組求解的方式進(jìn)行，同時(shí)也需要遵循結(jié)構(gòu)風(fēng)險(xiǎn)最小化原則. 最小二乘支持向量機(jī)的標(biāo)準(zhǔn)輸出量可定義為：

式中：K(x,xi)為滿足Mercer條件的核函數(shù)；αi為回歸系數(shù)；b為回歸偏差；x為樣本變量；xi為第i個(gè)預(yù)測(cè)樣本值；l為樣本數(shù)目.

對(duì)于實(shí)際工程的設(shè)計(jì)和施工來(lái)說(shuō)，通過(guò)待反演參數(shù)獲得的計(jì)算值和實(shí)測(cè)值之間的差距越小，則說(shuō)明位移反分析獲得最有效果. 因此，目標(biāo)函數(shù)通常取以下形式：

其中：x={X1,X2,…,Xk,…,Xm}，為m個(gè)參數(shù)構(gòu)成的一個(gè)樣本，Xk為待反演的巖體參數(shù)，如彈性模量E、內(nèi)聚力c、內(nèi)摩擦角φ等；m為待反演的巖體參數(shù)個(gè)數(shù)；fi(X)為第i個(gè)位移計(jì)算值，ui為位移實(shí)測(cè)值；n為測(cè)點(diǎn)數(shù).

與傳統(tǒng)SVM 相比，LSSVM 計(jì)算復(fù)雜性更低，求解完成的解答中包含多數(shù)甚至全部的訓(xùn)練樣本，更加便于在線計(jì)算，能夠通過(guò)擴(kuò)展其他形式來(lái)處理動(dòng)態(tài)問(wèn)題，在求解過(guò)程中結(jié)合支持向量機(jī)、前向神經(jīng)網(wǎng)絡(luò)等進(jìn)行討論，與機(jī)器學(xué)習(xí)的本質(zhì)更加貼合.

1.3 參數(shù)反分析模型

根據(jù)以上原理，建立螢火蟲(chóng)-支持向量機(jī)參數(shù)反分析計(jì)算步驟如圖1所示，利用該流程可將數(shù)值計(jì)算結(jié)果作為樣本進(jìn)行學(xué)習(xí)，從而得到與監(jiān)測(cè)位移相互吻合的巖體力學(xué)參數(shù).

圖1 參數(shù)反分析模型Fig.1 Parameter back analysis model

2 工程概況及數(shù)值模型建立

江蘇溧陽(yáng)抽水蓄能電站屬破碎巖層大型地下廠房洞室群工程，該工程地質(zhì)條件十分復(fù)雜，總體上巖體以完整性較差為主，巖體結(jié)構(gòu)主要為鑲嵌碎裂結(jié)構(gòu)，部分屬于層狀結(jié)構(gòu)，圍巖穩(wěn)定問(wèn)題突出，施工安全壓力大. 圍巖類別：Ⅲ類占45%，Ⅳ類占52%，Ⅴ類約占3%.

根據(jù)溧陽(yáng)抽水蓄能電站地下洞室群工程地質(zhì)條件，建立地下洞室群塊體離散元數(shù)值模型［18］，如圖2 所示，地下洞室群工程數(shù)值計(jì)算模型共計(jì)約3 萬(wàn)個(gè)塊體，采取對(duì)邊界進(jìn)行法向位移約束的形式，其中主廠房區(qū)域整體計(jì)算模型如圖3所示.

圖2 溧陽(yáng)抽水蓄能電站地下洞室群布置圖Fig.2 Underground cave group layout of Liyang pumped storage power station

圖3 主廠房整體計(jì)算模型三維示意圖Fig.3 Three-dimensional diagram of the overall calculation model of the main workshop

3 主廠房圍巖力學(xué)參數(shù)反分析研究

地下洞室工程中，洞室的變形往往與側(cè)向約束條件有很大關(guān)系，因此必須考慮初始應(yīng)力場(chǎng). 在考慮主廠房區(qū)域地應(yīng)力影響時(shí)需要考慮兩個(gè)水平向側(cè)壓力系數(shù)，并轉(zhuǎn)化到模型坐標(biāo)系下，因此，在反演研究區(qū)域初始地應(yīng)力時(shí)，實(shí)際反演兩個(gè)側(cè)壓力系數(shù). 圍巖變形的主要影響因素為巖體彈性模量E和泊松比μ. 地下洞室變形側(cè)向變形主要取決于應(yīng)力場(chǎng)的變化，而泊松比影響較小，但是彈性模量對(duì)于圍巖變形有著較大影響. 因此本文在反演參數(shù)的選取中選擇彈性模量E，而采用確定的泊松比μ. 巖體力學(xué)建議取值如表1、表2所示.

表2 結(jié)構(gòu)面物理力學(xué)參數(shù)Tab.2 Physical and mechanical parameters of structural surfaces

根據(jù)正交設(shè)計(jì)原理設(shè)計(jì)試驗(yàn)方案，利用塊體離散元數(shù)值模型計(jì)算得到各個(gè)方案下的計(jì)算位移值，可得LSSVM算法的訓(xùn)練樣本如表3所示. 采用螢火蟲(chóng)算法優(yōu)化模型參數(shù)，得到反映巖體參數(shù)和圍巖位移之間映射關(guān)系的最優(yōu)LSSVM模型.

表3 正交試驗(yàn)獲得的正演訓(xùn)練樣本Tab.3 Forward training samples obtained by orthogonal experiment

將方案1～16作為訓(xùn)練樣本，基于螢火蟲(chóng)算法的最小二乘支持向量機(jī)算法，以實(shí)際監(jiān)測(cè)位移作為樣本輸入到已經(jīng)訓(xùn)練完成的LSSVM模型中，利用圖1流程進(jìn)行學(xué)習(xí)預(yù)測(cè). 得出LSSVM最佳參數(shù)組合如表4所示.

表4 圍巖參數(shù)反演結(jié)果Tab.4 Inversion results of surrounding rock parameters

將表4 中反演得到的參數(shù)代入數(shù)值計(jì)算模型，進(jìn)行正分析，選取典型（C1）斷面的6 個(gè)監(jiān)測(cè)點(diǎn)進(jìn)行現(xiàn)場(chǎng)實(shí)測(cè)變形和數(shù)值模擬監(jiān)測(cè)點(diǎn)變形進(jìn)行對(duì)比分析，如圖4、圖5 所示. 可見(jiàn)現(xiàn)場(chǎng)實(shí)測(cè)變形和數(shù)值模擬監(jiān)測(cè)點(diǎn)變形具有較好的對(duì)應(yīng)，各點(diǎn)變形與現(xiàn)場(chǎng)監(jiān)測(cè)數(shù)據(jù)基本一致. 開(kāi)挖過(guò)程中圍巖變形最大位于拱頂，現(xiàn)場(chǎng)實(shí)測(cè)值為15.1 mm，數(shù)值模擬值為12.2 mm；圍巖變形最小位于上游邊墻，現(xiàn)場(chǎng)實(shí)測(cè)值為-0.5 mm，數(shù)值模擬值為2.2 mm；由于受開(kāi)挖過(guò)程中母線洞開(kāi)挖的影響，明顯存在下游圍巖變形大于上游的現(xiàn)象，這些現(xiàn)象都與工程實(shí)際相符合，因此該模型獲得的反演參數(shù)較為準(zhǔn)確，可為相關(guān)研究提供參考（圖4 和圖5 中MC1-1，MC1-2，MC1-3，MC1-4，MC1-5，MC1-6都為監(jiān)測(cè)點(diǎn)編號(hào)）.

圖4 現(xiàn)場(chǎng)實(shí)測(cè)變形曲線Fig.4 Deformation curves measured on site

圖5 監(jiān)測(cè)點(diǎn)變形曲線的數(shù)值模擬Fig.5 Numerical simulations of deformation curves of monitoring points

4 多洞聯(lián)合作用下的變形協(xié)同分析

溧陽(yáng)抽水蓄能電站地下洞室群，共有3大洞室，分別是主廠房、主變洞和尾閘室，其中，主廠房共分為7層開(kāi)挖，主變洞共分為3層開(kāi)挖，整個(gè)地下洞室群共分8步進(jìn)行開(kāi)挖，其中，第4層開(kāi)挖為母線洞，第8層開(kāi)挖為交通洞.

4.1 分層開(kāi)挖下的圍巖位移場(chǎng)分布特征

分層開(kāi)挖步驟下的圍巖位移場(chǎng)分布如圖6所示. 由圖6（a）可知，第一層開(kāi)挖洞室頂拱部分后，洞室圍巖主要產(chǎn)生豎直方向的變形，其中，拱頂以下沉變形為主，拱底位置處由于產(chǎn)生向上的回彈變形而出現(xiàn)“底臌”現(xiàn)象. 隨著第二層進(jìn)行開(kāi)挖，洞室高邊墻逐漸形成，此時(shí)洞室圍巖變形逐漸轉(zhuǎn)變?yōu)橐运阶冃螢橹? 第二層開(kāi)挖后，開(kāi)挖面逐漸遠(yuǎn)離頂拱部分，頂拱圍巖變形受開(kāi)挖擾動(dòng)影響越來(lái)越小，圍巖變形速度減慢. 此時(shí)，洞室上下游邊墻的變形基本對(duì)稱分布. 隨著開(kāi)挖的進(jìn)行，主廠房下游側(cè)圍巖變形區(qū)和主變洞上游側(cè)變形區(qū)發(fā)生交接，兩洞室的圍巖變形開(kāi)始變?yōu)槭?次開(kāi)挖擾動(dòng)影響，在雙重開(kāi)挖擾動(dòng)作用下，主廠房邊墻開(kāi)始發(fā)生非對(duì)稱變形，下游側(cè)壁圍巖變形明顯增大. 洞室開(kāi)挖完成后，從圖6（h）中可以看出，主廠房變形最為嚴(yán)重，最大變形點(diǎn)位于下游側(cè)巖錨梁處，主廠房頂拱圍巖最大變形發(fā)生在上游側(cè)拱頂處，邊墻最大變形發(fā)生在下游側(cè)邊墻處；主變洞頂拱圍巖最大變形同樣發(fā)生在上游側(cè)拱頂位置處，邊墻最大變形發(fā)生在上游側(cè)邊墻.

圖6 不同開(kāi)挖步驟下的位移圖Fig.6 Displacement diagrams under different excavation steps

總體上，洞室頂拱和底板以豎向變形為主，高邊墻以水平變形為主. 頂拱變形受開(kāi)挖面擾動(dòng)影響越來(lái)越弱，圍巖變形速度越來(lái)越慢，邊墻圍巖水平變形受開(kāi)挖擾動(dòng)影響較大，尤其是受母線洞的開(kāi)挖影響，開(kāi)挖母線洞時(shí)，主廠房、主變洞邊墻變形發(fā)生明顯的增大. 在實(shí)際開(kāi)挖過(guò)程中，應(yīng)加強(qiáng)支護(hù)工作，尤其是主廠房下游側(cè)巖錨梁處的支護(hù).

4.2 分層開(kāi)挖下的圍巖應(yīng)力場(chǎng)分布特征

圖7為地下洞室群在不同開(kāi)挖步驟下的最大主應(yīng)力云圖. 由圖可以明顯地看出，由于地下廠房的開(kāi)挖擾動(dòng)作用，在洞室開(kāi)挖面附近產(chǎn)生了明顯的應(yīng)力重分布現(xiàn)象，形成一定范圍的拉應(yīng)力區(qū)，在遠(yuǎn)離洞室開(kāi)挖面一定距離，圍巖處于原巖應(yīng)力狀態(tài)，處于壓應(yīng)力區(qū).

由圖7（a）可知，地下洞室群第一層頂拱部分開(kāi)挖后，在拱頂、拱底產(chǎn)生較大的拉應(yīng)力區(qū)，引起廠房頂拱發(fā)生拱頂下沉、拱底臌起現(xiàn)象，通過(guò)4.1節(jié)的洞室開(kāi)挖位移場(chǎng)可以得到印證. 同時(shí)，在頂拱與邊墻交接處形成了明顯的應(yīng)力集中現(xiàn)場(chǎng). 由圖7（b）可知，從第二層開(kāi)挖開(kāi)始，地下廠房的高邊墻逐漸形成，洞室邊墻臨空部分越來(lái)越大，在應(yīng)力場(chǎng)作用下，邊墻產(chǎn)生以水平方向?yàn)橹鞯淖冃?，從圖中可以看出在洞室邊墻位置形成較大范圍的拉應(yīng)力區(qū)，在拉應(yīng)力的作用下洞室圍壓產(chǎn)生向洞室方向的變形，產(chǎn)生“片幫”現(xiàn)場(chǎng)，其中，由圖7（d）可知，在洞室第4步開(kāi)挖母線洞后，在主廠房、主變洞與母線洞交接部位形成明顯的應(yīng)力集中現(xiàn)象. 隨著開(kāi)挖的進(jìn)行，廠房頂拱部分距離開(kāi)挖面越來(lái)越遠(yuǎn)，開(kāi)挖擾動(dòng)作用對(duì)頂拱的影響越來(lái)越弱，集中的應(yīng)力逐漸得到釋放，應(yīng)力水平逐漸降低.

由圖7（h）可知，整個(gè)地下洞室群開(kāi)挖完成后，在3大洞室的邊墻部位形成了較大范圍的拉應(yīng)力區(qū)域，應(yīng)力集中主要發(fā)生在頂拱與邊墻、邊墻與母線洞交界處.

4.3 分層開(kāi)挖下的圍巖塑性區(qū)分布特征

圖8為地下洞室群在不同開(kāi)挖步驟下的各剖面塑性區(qū)分布圖. 可以看出，地下廠房開(kāi)挖后，在洞室表層圍巖處產(chǎn)生一定范圍的剪拉應(yīng)力帶.

由圖8（a）可知，3大洞室頂拱部分開(kāi)挖后，在拱頂部位形成了剪切應(yīng)力帶，第2步開(kāi)挖開(kāi)始后，廠房頂拱的應(yīng)力集中現(xiàn)象更加明顯，拱頂產(chǎn)生少量剪切破壞，隨著開(kāi)挖面的不斷深入，頂拱部位圍壓逐漸趨于穩(wěn)定，后續(xù)洞室開(kāi)挖對(duì)頂拱塑性區(qū)分布影響不明顯；洞室高邊墻逐漸開(kāi)挖形成后，在洞室圍巖尤其是臨空面附近圍巖形成了較大區(qū)域的剪切和張拉塑性破壞區(qū)，其中，母線洞巖體開(kāi)挖后，在母線洞與主廠房、母線洞與主變洞交接部位形成了明顯的應(yīng)力集中現(xiàn)象.

圖8 不同開(kāi)挖步驟下的典型剖面塑性區(qū)圖Fig.8 Plastic zone diagrams of typical sections under different excavation steps

隨著廠房分層開(kāi)挖的不斷進(jìn)行，洞室圍巖中產(chǎn)生的塑性區(qū)范圍也在不斷增大，所產(chǎn)生的破壞形式主要為剪切破壞和剪切張拉破壞，引起的破壞區(qū)主要集中在洞室上下游邊墻部位、邊墻與母線洞交接部位. 對(duì)比分析主廠房上下游側(cè)的圍巖塑性區(qū)范圍可知，主廠房下游側(cè)產(chǎn)生的塑性區(qū)范圍要比上游側(cè)塑性區(qū)大，主變洞上游側(cè)產(chǎn)生的應(yīng)力集中要比下游側(cè)應(yīng)力集中明顯，主要原因是母線洞的開(kāi)挖擾動(dòng)對(duì)主廠房的下游側(cè)、主變洞的上游側(cè)影響比較大.

5 結(jié)論

將螢火蟲(chóng)算法應(yīng)用于最小二乘支持向量機(jī)的核函數(shù)及其參數(shù)的優(yōu)化，建立了螢火蟲(chóng)算法-最小二乘支持向量機(jī)反分析模型. 結(jié)合溧陽(yáng)抽水蓄能電站基本地質(zhì)條件和監(jiān)測(cè)數(shù)據(jù)，采用數(shù)值模擬和智能算法結(jié)合的方法進(jìn)行了參數(shù)反分析研究. 得到如下結(jié)論：

1）采用螢火蟲(chóng)優(yōu)化算法對(duì)LSSVM的兩個(gè)參數(shù)隨機(jī)搜索，具體思路為首先正交試驗(yàn)設(shè)計(jì)得到訓(xùn)練樣本和測(cè)試樣本，由螢火蟲(chóng)算法對(duì)LSSVM的參數(shù)進(jìn)行優(yōu)化，獲得最佳LSSVM預(yù)測(cè)模型；

2）根據(jù)最佳LSSVM預(yù)測(cè)模型，建立溧陽(yáng)抽水蓄能電站地下洞室群塊體離散元數(shù)值模型，并對(duì)比分析數(shù)值模擬結(jié)果與現(xiàn)場(chǎng)監(jiān)測(cè)數(shù)據(jù)，結(jié)果表明數(shù)值結(jié)果與現(xiàn)場(chǎng)監(jiān)測(cè)數(shù)據(jù)基本吻合，所反映的圍巖變形規(guī)律一致，驗(yàn)證了反演參數(shù)的合理性；

3）通過(guò)分析分層開(kāi)挖下的圍巖位移場(chǎng)、應(yīng)力場(chǎng)和塑性區(qū)，發(fā)現(xiàn)洞室頂拱和底板以豎向變形為主，高邊墻以水平變形為主，應(yīng)力集中主要發(fā)生在頂拱與邊墻、邊墻與母線洞交界處，引起的破壞區(qū)主要集中在洞室上下游邊墻部位、邊墻與母線洞交接部位.