韓 峰，徐 磊，張?zhí)?/p>

(1.浙江省水利水電勘測設(shè)計院，浙江杭州 310002;2.河海大學(xué)水利水電學(xué)院，江蘇 南京 210098)

各類壩基巖體都是在漫長地質(zhì)歷史中形成的天然地質(zhì)體，在給定力學(xué)模型的前提下，合理確定其力學(xué)參數(shù)是開展壩體-壩基系統(tǒng)應(yīng)力變形分析及穩(wěn)定性評價并取得合理成果的關(guān)鍵之一［1］。壩基巖體力學(xué)參數(shù)的確定一般采用試驗與工程類比相結(jié)合的方法，但由于壩基巖體固有的復(fù)雜性和不確定性，基于這一方法所確定的力學(xué)參數(shù)往往與實際情況不符，導(dǎo)致了以巖體為基礎(chǔ)的壩工設(shè)計分析成果與實際情況不符。隨著計算機、現(xiàn)代監(jiān)測儀器和監(jiān)測技術(shù)、巖石力學(xué)理論和數(shù)值計算方法的不斷發(fā)展、應(yīng)用，基于實測變形資料的力學(xué)參數(shù)反演分析方法逐漸成為解決這一問題的有效手段［2-4］。壩基巖體力學(xué)參數(shù)的反演本質(zhì)上是一個相關(guān)參數(shù)的尋優(yōu)過程，因而如何選擇一個高效的優(yōu)化算法是其核心問題之一［5］。自20世紀80年代以來，相關(guān)學(xué)者提出并發(fā)展了一些區(qū)別于傳統(tǒng)優(yōu)化方法的智能優(yōu)化算法，如人工神經(jīng)網(wǎng)絡(luò)(ANN)法、蟻群算法(ACO)、遺傳算法(GA)及粒子群優(yōu)化算法(PSO)等，其中，PSO算法因其收斂快、魯棒性好、隱含并行性及易于實現(xiàn)［6］等特點，一經(jīng)出現(xiàn)即得到了很快的發(fā)展。目前，已有相關(guān)學(xué)者將PSO算法應(yīng)用于初始地應(yīng)力場反演等諸多領(lǐng)域，效果良好［7-12］。

筆者將PSO算法引入壩基巖體力學(xué)參數(shù)反演領(lǐng)域，并提出基于PSO算法［13］和ABAQUS的壩基巖體力學(xué)參數(shù)聯(lián)合反演法。在介紹粒子群優(yōu)化算法基本理論的基礎(chǔ)上，構(gòu)建了壩基巖體力學(xué)參數(shù)的數(shù)學(xué)反演模型，給出了具體的實施流程并開發(fā)了相應(yīng)的程序，實現(xiàn)了基于施工和運行期實測變形資料的壩基巖體力學(xué)參數(shù)PSO-ABAQUS聯(lián)合反演。算例分析表明，所提出的聯(lián)合反演法不僅簡便易行，而且反演效率較高，是一種合理確定壩基巖體力學(xué)參數(shù)的通用方法。

1 壩基巖體力學(xué)參數(shù)的反演模型

反演分析按照所采用方法的不同，總體上可以分為正反分析法和逆反分析法兩大類。在壩基巖體力學(xué)參數(shù)反演領(lǐng)域，由于所涉及的研究對象為具有復(fù)雜工程特性的天然地質(zhì)體，因此通常不具備采用逆反分析法進行參數(shù)反演的前提條件，故筆者采用正反分析法開展參數(shù)反演。

對于正反分析法，基于實測變形數(shù)據(jù)(通常為實測位移)的壩基巖體力學(xué)參數(shù)反演本質(zhì)上是一個數(shù)學(xué)優(yōu)化問題，待反演的力學(xué)參數(shù)即為此優(yōu)化問題的設(shè)計變量，尋優(yōu)過程即為尋找滿足目標(biāo)函數(shù)極小的一組參數(shù)來使得正演分析的計算結(jié)果(通常采用有限元法)與監(jiān)測數(shù)據(jù)一致。

在壩工建設(shè)和運行過程中，通過布置各類監(jiān)測儀器，可以獲取包括壩體應(yīng)力、變形在內(nèi)的多種監(jiān)測資料。理論上，監(jiān)測資料都可以作為壩基巖體力學(xué)參數(shù)反演分析的基礎(chǔ)信息，但考慮到監(jiān)測資料的精度以及待反演參數(shù)對監(jiān)測資料的敏感性，目前在壩基巖體力學(xué)參數(shù)反演分析中，大多選取實測位移作為反演分析的基礎(chǔ)信息［14］，鑒于此，本文所建立的壩基巖體力學(xué)參數(shù)反演模型中亦選擇實測位移資料作為基礎(chǔ)反演信息。

基于數(shù)學(xué)優(yōu)化基本理論，建立壩基巖體力學(xué)參數(shù)反演數(shù)學(xué)模型:

其中

式中:φ(X)——目標(biāo)函數(shù);X——m維反演參數(shù)向量，具體的反演參數(shù)需在地質(zhì)調(diào)查分析的基礎(chǔ)上確定;hi(X)=0——第i個等式約束;gj(X)≤0——第j個不等式約束;N——位移實測點數(shù)目;M——位移分量數(shù)目;U'pq(X)——與反演參數(shù)向量X對應(yīng)的第p測點的計算位移值;Upq——第p個測點的實測位移值。

對于待反演參數(shù)的一組特定取值，一般首先采用有限單元法來求解與之相對應(yīng)的位移場，進而提取測點處的位移值。

在基于PSO算法的壩基巖體力學(xué)參數(shù)反演中，計算量最大的部分是每一組試算力學(xué)參數(shù)所對應(yīng)的正演計算以及后續(xù)目標(biāo)函數(shù)值(各粒子適應(yīng)值)的計算，正演計算模型通常都很復(fù)雜，且計算規(guī)模很大，因而選擇一個功能強大且通用的正演分析計算工具十分必要?？紤]到通用有限元軟件平臺ABAQUS的適用性和計算效率［15］，本文采用ABAQUS作為正演分析的計算工具。

2 PSO算法

PSO算法是1995年由Kennedy和Eberhart提出的一種群體智能優(yōu)化算法［13］。在PSO算法中，優(yōu)化問題的每一個潛在解都是搜索空間中的一個粒子，若干個粒子構(gòu)成一個粒子群。粒子的適應(yīng)度即為被優(yōu)化問題的目標(biāo)函數(shù)值。每個粒子都具有位置和速度2個特征，位置決定粒子的適應(yīng)度，而速度決定下一代粒子運動的方向和距離。在迭代過程中，粒子通過跟蹤個體極值pBest和全局極值gBest實現(xiàn)代際更新，直到在整個搜索空間中找到最優(yōu)解或達到最大迭代次數(shù)為止。

假設(shè)在一個m維的搜索空間中，由n個粒子構(gòu)成一個粒子群X'={x1，…，xa，…，xn}，其中，第a個粒子的位置為 xa={xa1，xa2，…，xam}T，速度為 va={va1，va2，…，vam}T，其個體極值為 pa={pa1，pa2，…，pam}T，全局極值為pg={pg1，pg2，…，pgm}T，則第a個粒子的速度和位置按式(2)、式(3)完成代際更新。

式中:k——粒子更新代數(shù);wk——慣性權(quán)重，控制PSO算法的搜索能力，wk取值較大，則全局尋優(yōu)能力強，局部尋優(yōu)能力弱，wk通常依代數(shù)的增加而線性遞減;rand()——均勻分布于(0，1)之間的隨機數(shù);c1，c2——加速常數(shù)，均可取為2.0。式(2)右邊第1部分為粒子先前行為的慣性;第 2部分為“認知(cognition)”部分，表示粒子本身的思考;第3部分為“社會(social)”部分，表示粒子間的信息共享與相互合作。

圖1 PSO-ABAQUS聯(lián)合反演程序流程Fig.1 Flow chart of PSO and ABAQUS combined inversion

3 程序開發(fā)及流程圖

在第1節(jié)、第2節(jié)的基礎(chǔ)上，筆者提出了壩基巖體力學(xué)參數(shù)的PSO-ABAQUS聯(lián)合反演法，并編制了相應(yīng)的FORTRAN計算程序。該程序在獲取有限元計算模型及優(yōu)化參數(shù)等相關(guān)信息的基礎(chǔ)上，首先采用隨機方法形成初始粒子群，然后基于每個粒子的位置，自動生成用以計算應(yīng)力場的INP文件并調(diào)用ABAQUS求解器來完成壩體-壩基系統(tǒng)的應(yīng)力變形分析，進而計算粒子的適應(yīng)度(即目標(biāo)函數(shù)值)，最后基于PSO算法完成粒子更新，直至找到符合預(yù)定最小閾值的最優(yōu)粒子或達到最大迭代次數(shù)為止。若達到最大迭代次數(shù)仍未收斂，可調(diào)整反演模型或增加最大迭代次數(shù)，再次運行反演程序。圖1為PSO-ABAQUS聯(lián)合反演程序流程。

4 算例分析

為了驗證本文所提PSO-ABAQUS聯(lián)合反演法的有效性其相應(yīng)FORTRAN程序的正確性，進行如下算例分析。假定一建于均質(zhì)巖基上的重力壩，壩高為50 m，壩頂寬10 mm，壩底寬40 mm。建立此重力壩的有限元模型(圖2)，模型共劃分平面應(yīng)變四結(jié)點四變形等參數(shù)單元337個，結(jié)點380個。在有限元計算模型中選取壩頂上、下游，壩趾，壩踵以及下游壩面折坡處5個結(jié)點作為變形觀測點。

在算例分析中，取反演參數(shù)為壩基巖體的變形模量，依據(jù)工程經(jīng)驗，反演參數(shù)取值范圍為1.0～30.0 GPa。

除反演參數(shù)外，模型介質(zhì)的其他物理力學(xué)參數(shù)如下:壩體混凝土密度ρ1=2400 kg/m3，彈性模量E=20 GPa，泊松比v1=0.167;壩基巖體密度ρ2=2 000 kg/m3，泊松比v2=0.25，摩擦系數(shù)f=1.05，黏聚力c=0.95 MPa。計算中壩體混凝土采用線彈性本構(gòu)模型，壩基巖體采用基于Drucker-Prager屈服準則的理想彈塑性本構(gòu)模型。計算以壩基巖體的自重應(yīng)力場作為初始地應(yīng)力場，作用壩體的荷載包括壩體自重、上下游面水壓力以及揚壓力。

首先給定壩基巖體變形模量(取為9 GPa)進行有限元計算，并以5個測點處有限元計算所得的位移值(見表1)作為反演所需的實測變形數(shù)據(jù)，隨后應(yīng)用PSO-ABAQUS聯(lián)合反演程序反演壩基巖體力學(xué)參數(shù)(變形模量)，最后將反演所得的變形模量以及相應(yīng)的變形場與給定的變形模量及相應(yīng)的變形場進行對比分析，以驗證本文所提方法以及程序開發(fā)的正確性及有效性。

圖2 重力壩有限元模型及測點位置Fig.2 FEM mesh of gravity dam and location of measuring points

利用PSO-ABAQUS聯(lián)合反演程序，對算例開展了壩基巖體力學(xué)參數(shù)反演分析。反演中取粒子群組數(shù)為5組，目標(biāo)函數(shù)預(yù)定最小閾值為0.05。反演所得的壩基巖體彈性模量為8.797 GPa，誤差為2.3%。

圖3為在給定力學(xué)參數(shù)條件下、對壩基巖體賦予反演參數(shù)下的壩體水平向、垂向位移分布云圖。由圖3可知，給定參數(shù)和反演參數(shù)下模型計算位移場基本一致。由表1可知，5個測點的實測位移值與基于反演參數(shù)的計算位移值亦非常接近。

表1 各測點位移Table 1 Displacement of measuring points mm

圖3 壩體水平向和豎向位移云圖(單位:mm)Fig.3 Contour of horizontal and vertical displacement of dam(units:mm)

圖4為了優(yōu)化反演中目標(biāo)函數(shù)的變化過程。由圖4可知，在較少粒子組數(shù)和較小閾值條件下，反演經(jīng)過27次迭代就達到了較高的精度，可以預(yù)計，若設(shè)定較多的粒子組數(shù)并采用較大的閾值，優(yōu)化迭代的次數(shù)還可以進一步減少，說明本文所提出的PSO-ABAQUS聯(lián)合反演法具有較高的效率。

圖4 目標(biāo)函數(shù)值變化過程線Fig.4 Change process of object function value

5 結(jié) 語

合理確定壩基巖體力學(xué)參數(shù)是開展壩體-壩基系統(tǒng)的應(yīng)力變形分析及穩(wěn)定性評價并取得合理成果的關(guān)鍵之一，但基于通常的試驗和工程類比所確定的力學(xué)參數(shù)往往與實際情況不符。鑒于此，本文將PSO算法引入壩基巖體力學(xué)參數(shù)反演領(lǐng)域，并提出基于PSO算法和ABAQUS的壩基巖體力學(xué)參數(shù)聯(lián)合反演法。在介紹粒子群優(yōu)化算法基本理論的基礎(chǔ)上，構(gòu)建了壩基巖體力學(xué)參數(shù)的數(shù)學(xué)反演模型，給出了具體的實施流程并開發(fā)了相應(yīng)的程序，實現(xiàn)了基于施工和運行期實測變形資料的壩基巖體力學(xué)參數(shù)PSO-ABAQUS聯(lián)合反演。算例分析表明，本文所提的PSO-ABAQUS聯(lián)合反演法不僅簡便易行，而且效率較高，是一種用于合理確定壩基巖體力學(xué)參數(shù)的通用方法。

［1］賀里堯.一般性壩基巖體力學(xué)參數(shù)取值研究路線的探討［D］.蘭州:蘭州大學(xué)，2009.

［2］張強勇，向文，楊文東，等.壩區(qū)巖體蠕變參數(shù)反演與邊坡開挖流變計算分析［J］.武漢大學(xué)學(xué)報:工學(xué)版，2008，41(5):72-76.(ZHANG Qiangyong，XIANG Wen，YANG Wendong，et al.Creep parameters inversion for dam zone rockmass and rheological computational analysis of slope excavation［J］.Engineering Journal of Wuhan University，2008，41(5):72-76.(in Chinese))

［3］常曉林，喻勝春，馬剛，等.基于粒子遷徙的粒群優(yōu)化算法及其在巖土工程中的應(yīng)用［J］.巖土力學(xué)，2011，32(4):1077-1082(CHANG Xiaolin，YU Shengchun，MA Gang，et al.Particle swarm optimization based on particle migration and its application to geotechnical engineering［J］.RocKand Soil Mechanics，2011，32(4):1077-1082.(in Chinese))

［4］張強勇，陳芳，楊文東，等.大崗山壩區(qū)巖體現(xiàn)場剪切蠕變試驗及參數(shù)反演［J］.巖土力學(xué)，2011，32(9):2584-2590.(ZHANG Qiangyong，CHEN Fang，YANG Wendong，et al.Field shear creep test and creep parameters inversion for dam zone rocKmass of Dagangshan hydropower station［J］.RocKand Soil Mechanics，2011，32(9):2584-2590.(in Chinese))

［5］陳益峰，周創(chuàng)兵.隔河巖壩基巖體在運行期的彈塑性力學(xué)參數(shù)反演［J］.巖石力學(xué)與工程學(xué)報，2002，21(7):968-975.(CHEN Yifeng，ZHOU Chuangbing.BacKanalysis on elastic-plastic mechanic parameters of rocKfoundation of Geheyan dam during running period［J］.Chinese Journal of RocKMechanics and Engineering，2002，21(7):968-975.(in Chinese))

［6］倪紹虎，肖明，何世海，等.地下工程并行優(yōu)化反演分析及算例驗證［J］.巖石力學(xué)與工程學(xué)報，2013，32(3):501-511.(NI Shaohu，XIAO Ming，HE Shihai，et al.BacKanalysis in underground engineering based on parallel computing and optimization algorithm and its verification［J］.Chinese Journal of RocKMechanics and Engineering，2013，32(3):501-511.(in Chinese))

［7］李金風(fēng)，楊啟貴，徐衛(wèi)亞.基于改進粒子群算法CHPSO-DS的面板壩堆石體力學(xué)參數(shù)反演［J］.巖石力學(xué)與工程學(xué)報，2008，27(6):1229-1235.(LI Jinfeng，YANG Qigui，XU Weiya.BacKanalysis mechanical parameters of rockfill based on modified particle swarm optimization CHPSO-DS［J］.Chinese Journal of RocKMechanics and Engineering，2008，27(6):1229-1235.(in Chinese))

［8］蘇國韶，張克實，呂海波.位移反分析的粒子群優(yōu)化-高斯過程協(xié)同優(yōu)化方法［J］.巖土力學(xué)，2011，32(2):510-515.(SU Guoshao，ZHANG Keshi，LYU Haibo.A cooperative optimization based on particle swarm optimization and Gaussian process for displacement bacKanalysis［J］.RocKand Soil Mechanics，2011，32(2):510-515.(in Chinese))

［9］羅潤林，阮懷寧，黃亞哲，等.巖體初始應(yīng)力場的粒子群優(yōu)化反演及在FLAC3D中的實現(xiàn)［J］.長江科學(xué)院學(xué)報，2008，25(4):73-76.(LUO Runlin，RUAN Huaining，HUANG Yazhe，et al.Particle swarm optimization inversion method of initial ground stress and implementation in FLAC3D［J］.Journal of Yangtze River Scientific Research Institute，2008，25(4):73-76.(in Chinese))

［10］張鵠志，郭建青.粒子群優(yōu)化算法在確定越流含水層參數(shù)中的應(yīng)用［J］.水利水電科技進展，2011，31(3):13-16.(ZHANG Huzhi，GUO Jianqing.Application of particle swarm optimization algorithm to determination of leakage aquifer parameters［J］.Advances in Science and Technology of Water Resource，2011，31(3):13-16.(in Chinese))

［11］李占超，侯會靜.基于改進粒子群優(yōu)化算法的施工期拱壩結(jié)構(gòu)性態(tài)反演分析［J］.水利水電科技進展，2011，31(4):24-28.(LI Zhanchao，HOU huijing.Inverse analysis of structural behaviors of arch dams during construction period based on improved particle swarm optimization［J］.Advances in Science and Technology of Water Resource，2011，31(4):24-28.(in Chinese))

［12］趙迪，張宗亮，陳建生.粒子群算法和ADINA在土石壩參數(shù)反演中的聯(lián)合應(yīng)用［J］.水利水電科技進展，2012，32(3):43-47.(ZHAO Di，ZHANG Zongliang，CHEN Jiansheng.A combined application of particle swarm optimization algorithm and ADINA for parametric inversion of earth rocKdam［J］.Advances in Science and Technology of Water Resources，2012，32(3):43-47.(in Chinese))

［13］張麗平.粒子群優(yōu)化算法的理論及實踐［D］.杭州:浙江大學(xué)，2005.

［14］賈善坡，伍國軍，陳衛(wèi)忠.基于粒子群算法與混合罰函數(shù)法的有限元優(yōu)化反演模型及應(yīng)用［J］.巖土力學(xué)，2011，32(增刊2):598-603.(JIA Shanpo，WUGuojun，CHENWeizhong.Application of finite element inverse model based on improved particle swarm optimization and mixed penalty function［J］.RocKand Soil Mechanics，2011，32(Sup2):598-603.(in Chinese))

［15］徐磊.基于ABAQUS的強度儲備安全系數(shù)精確自動求解［J］.三峽大學(xué)學(xué)報:自然科學(xué)版，2011，33(3):10-13.(XU Lei.Accurate auto-calculation of strength reserve safety coefficients based on ABAQUS［J］.Journal of China Three Gorges University:Natural Sciences，2011，33(3):10-13.(in Chinese))