高堆石壩瞬變-流變參數(shù)三維全過程 聯(lián)合反演方法及變形預測

馬 剛 ，常曉林 ，周 偉 ，花俊杰

（1. 武漢大學 水資源與水電工程科學國家重點實驗室，武漢 430072；2. 武漢大學 水工巖石力學教育部重點實驗室，武漢 430072）

1 引 言

堆石體的變形是長期而復雜的過程，其變形包括瞬時變形和隨時間變化的流變變形，采用有限元等數(shù)值方法分析面板堆石壩的應力與變形已得到廣泛運用。在面板堆石壩的應力與變形分析中，選用合理的堆石體本構模型以及準確的模型參數(shù)是整個分析的關鍵，堆石體的參數(shù)一般由室內或現(xiàn)場試驗獲得，然而受試驗條件、縮尺效應的限制和堆石材料自身性質的離散性，使得測定的力學特性參數(shù)與實際值存在一定的差異，由此計算的堆石壩應力、變形與實測值差別較大，如室內流變試驗一般幾個月就完成了，而壩體的后期變形往往持續(xù)數(shù)年。因此，有必要利用壩體實測位移資料對堆石體的參數(shù)進行反演分析，并進行堆石壩后期變形預測。

以實測位移為基礎的位移反分析方法已成為巖土工程學科的一個重要組成部分。以往，巖土工程中的反演方法常采用直接法，將參數(shù)反演問題轉換成優(yōu)化問題，文獻[1]根據(jù)數(shù)座已建土石壩的原型觀測資料采用復合形法進行了流變參數(shù)的反演分析。但由于該類優(yōu)化問題的復雜性，傳統(tǒng)的優(yōu)化方法很難收斂到全局最優(yōu)解，隨著智能優(yōu)化算法的飛速發(fā)展，遺傳算法[2－4]、蟻群算法[5－6]、粒子群算法[7]、序列二次規(guī)劃法[8]等全局智能優(yōu)化算法被用于求解該問題，取得了較好的效果。在參數(shù)反演過程中需要反復進行有限元正分析，對于高堆石壩來說計算量非常大，而利用神經(jīng)網(wǎng)絡建立待反演參數(shù)與壩體位移之間的映射關系代替有限元正分析將大大提高計算效率[6－9]。徑向基函數(shù)(radial basis function，RBF)網(wǎng)絡是一種新型的神經(jīng)網(wǎng)絡模型[10]，與BP 網(wǎng)絡相比，具有生物學基礎和數(shù)學基礎，隱節(jié)點具有局部特性，逼近能力更強，具有較好的泛化能力和函數(shù)逼近能力，可實現(xiàn)性能優(yōu)越的非線性預測器。

在以往的研究中，堆石體參數(shù)反演分析主要是針對靜力本構模型參數(shù)[2－7]，也有對流變參數(shù)進行反演的[1，8]。反演流變參數(shù)時，通常選擇一個間歇期長的臨時斷面進行流變參數(shù)反演[8]，或者將滿蓄以前的變形看作是瞬時變形，將滿蓄以后的變形作為流變變形，以此為前提進行流變參數(shù)的反演分析。考慮到堆石壩實際施工過程和變形機制的復雜性，很難將流變變形與瞬時變形分開，因此，有必要對靜力本構模型和流變模型的參數(shù)進行綜合反演。本文在對靜力本構模型參數(shù)和流變參數(shù)進行敏感性分析的基礎上，選擇對壩體變形敏感的參數(shù)作為待反演參數(shù)，采用基于粒子遷徙的粒子群算法和徑向基函數(shù)神經(jīng)網(wǎng)絡構建反演平臺進行參數(shù)反演，對水布埡面板堆石壩進行了參數(shù)反演分析，并預測了壩體變形穩(wěn)定時的沉降值。

2 堆石體靜力本構模型和流變模型

2.1 南水模型

南水模型[11]克服了鄧肯-張E-B 模型不能考慮堆石體剪脹(縮)的缺點，能反映堆石體的壓硬性和各向異性。南水模型假定應力空間中存在雙屈服面，分別采用橢圓和冪函數(shù)表示，其屈服函數(shù)為

式中：p、q 分別為平均正應力和八面體剪應力；對于堆石體r、s 可取2。

南水模型中切線模量tE 的計算和鄧肯-張模型相同。

式中：iE 為初始切線模量；c 為材料凝聚力；φ 為內摩擦角；lS 為應力水平，它反映材料強度發(fā)揮程度；fR 為破壞比；k、n 為模型參數(shù)；ap 為一個標準大氣壓；1σ 、3σ 分別為大小主應力。

南水模型采用拋物線描述體積應變與軸向應變的關系，可以統(tǒng)一考慮堆石體的剪脹(縮)性。采用切線體積比代替E-B 模型中的切線體積模量tB ：

式中： Rs= RfS； cd、 nd、 Rd為模型參數(shù)。

2.2 長科院冪函數(shù)流變模型

目前，在高堆石壩流變分析中常采用的是能模擬高圍圧條件的冪函數(shù)流變模型[12]：

式中：sfε 、vfε 為某個應力狀態(tài)下最終軸向流變量和最終體積流變量；τ 為時間；s( )ε τ 、v( )ε τ 分別為0 τ→ 時段內累積的軸向和體積流變量；sλ 、vλ為反映軸向和體積流變速率的參數(shù)。

最終軸向流變量sfε 和應力水平lS 、圍圧3σ 的關系如下：

最終體積流變量vfε 和應力水平lS 、圍圧3σ 可用線性函數(shù)擬合：

式(4)～(6)及參數(shù)c、d、η 、m、ac 、ad 、cβ、dβ、vλ 完整地表達了堆石體的流變特性。

堆石壩都是采用分層分區(qū)填筑，其加載過程復雜，要準確確定某一澆筑層具體的初始流變發(fā)生時間以及應力狀態(tài)變化后的后續(xù)流變發(fā)生時間十分困難，因此，采用文獻[1]的方法，采用相對時間代替絕對時間。

2.3 參數(shù)敏感性分析

南水模型有8 個參數(shù)，長科院流變模型有9 個參數(shù)，若都進行反演，則工作量很大，因此，有必要進行參數(shù)敏感性分析。修正的Morris 法采用靈敏度指標S 來反映參數(shù)的敏感性[6]：

式中：S 為靈敏度指標；iy 為模型第i 次計算輸出值；0y 為初始參數(shù)對應的模型輸出值；ix 為第i 次計算時參數(shù)；0x 為初始參數(shù)；n 為計算次數(shù)。

文獻[6]對南水模型的參數(shù)進行了敏感性分析，認為k、φ、fR 、dc 、dn 、dR 對壩體的沉降比較敏感。本文采用一均質壩對長科院冪函數(shù)流變模型進行了參數(shù)敏感性分析，結果如圖1 所示，可以看出，參數(shù)c、η 、cα、cβ、vλ 對壩體沉降較為敏感。

圖1 長科院冪函數(shù)流變模型參數(shù)敏感性分析 Fig.1 Parameter sensitivity analysis of creep model

3 基于MPSO-RBF 網(wǎng)絡的參數(shù)反演

3.1 基于粒子遷徙的粒群算法

在粒群算法[13]中，每個粒子代表待優(yōu)化問題在多維空間中的一個潛在解。每個粒子具有位置和速度兩個特征，粒子位置對應的目標函數(shù)值即可作為該粒子的適應度值。每個粒子根據(jù)它自身的“經(jīng)驗”和同伴的“經(jīng)驗”在搜索空間中向更好的位置“飛行”，直到在整個搜索空間中找到最優(yōu)解或達到最大迭代次數(shù)為止。

在經(jīng)典粒群算法(classic particle swarm optimization，CPSO)中，慣性權重ω 和加速因子1c ，2c 在優(yōu)化過程中保持不變。在PSO 算法中，合理地調整全局搜索和局部開發(fā)的關系是提高算法性能的關鍵之一，較大的ω 值有利于粒子充分的在搜索空間中探索，較小的ω 值有助于粒子在當前位置的附近搜索。根據(jù)文獻[14]提出的線性遞減策略，在線性遞減權重粒群算法中，在演化的早期階段，采用較大的ω 值，在演化后期逐漸減小ω 值，慣性權重在優(yōu)化過程中采用下式描述：

式中：iter 為當前迭代次數(shù)；i termax為允許最大迭代次數(shù)； ωmax、 ωmin分別為最大和最小慣性權重。

在群體智能算法中，提高算法性能的另一個關鍵是在演化的早期鼓勵粒子在整個搜索空間“漫游”，不至于聚集在超級粒子周圍而導致早熟收斂；在演化的后期，鼓勵粒子“飛向”搜索到的最優(yōu)值，以加速收斂。基于以上認識，文獻[15]提出了線性變化的加速因子，在該策略中代表“認知”部分1c 隨迭代次數(shù)逐漸減小，代表“社會”部分的2c 隨迭代次數(shù)逐漸增大：

式中：1ic ，1fc ，2ic 和2fc 分別是初始和最終“認知”、“社會”加速因子。

受自然界物種遷徙能保持種群多樣性的啟示，本文采用一種新的改進的粒群算法(MPSO)。算法初始化為一群隨機粒子，然后粒子被隨機劃分為若干個子粒群。每個子粒群獨立演化，演化策略采用考慮了線性遞減的慣性權重、線性變化的加速因子和自適應的變異算子。在演化的過程中，每隔若干迭代次數(shù)，進行一次粒子遷徙。粒子遷徙時，不僅將粒子目前的位置代入新的粒群中，還將粒子的個體極值pBest 引入新的粒群，以此加強子粒群間的信息交流并提高粒群的多樣性。

3.2 RBF 神經(jīng)網(wǎng)絡

徑向基函數(shù)(RBF)網(wǎng)絡是一種兩層前向神經(jīng)網(wǎng)絡，包括一個具有徑向基函數(shù)的隱層和一個具有線性神經(jīng)元的輸出層，它能以任意精度逼近任意函數(shù)，具有較強的逼近能力，特別適合解決函數(shù)逼近問題，網(wǎng)絡結構如圖2 所示。

圖2 RBF 神經(jīng)網(wǎng)絡結構圖 Fig.2 Structure of RBF ANN

選用高斯函數(shù)作為徑向基函數(shù)：

式中：σ 決定徑向基函數(shù)的形狀，σ 越大則基函數(shù)越平滑。

3.3 確定待反演參數(shù)

根據(jù)參數(shù)敏感性分析可知，k、φ、fR 、dc 、dn 、dR 、c、η 、cα、cβ、vλ 對壩體的變形比較敏感。由于堆石壩有多種材料分區(qū)，如果將上述所有敏感的參數(shù)均作為待反演的變量，工作量仍然較大，而且無法保證反演結果收斂到正確值，所以應該減少待反演參數(shù)，去掉測試方法較成熟的φ 和工程經(jīng)驗豐富的fR[2－3]。計算中的待反演參數(shù)共13個，分別為1k 、d1c 、d1n 、d1R 、2k 、d2c 、d2n 、d2R 、c、η 、cα、cβ、vλ ，下標1、2 分別對應主堆石和次堆石，主、次堆石采用相同的流變參數(shù)。

3.4 反演分析目標函數(shù)

堆石壩的參數(shù)反演就是尋找一組參數(shù)使計算位移值與實測位移值最佳逼近，由于堆石壩的測點眾多，因此，上面所說的最佳逼近是指總體上和平均意義上的最好近似。目標函數(shù)可取為監(jiān)測點的計算位移值與實測位移值差的二范數(shù)式，由于每個監(jiān)測點的監(jiān)測數(shù)據(jù)都是一個時間序列，故目標函數(shù)取為

式中： x1， x2， ，x13對應一組待反演的堆石體參數(shù)；m 為監(jiān)測斷面的個數(shù)；wi為第i 個斷面的權重系數(shù)；ni為第i 個斷面上監(jiān)測點的個數(shù)；為第i 個斷面上第j 個監(jiān)測點在第k 個時間點的沉降計算值；為相應的實測值。

4 水布埡面板堆石壩變形參數(shù)的綜合反演及變形預測

水布埡面板堆石壩，壩高為233 m，壩體和混 凝土面板的有限元網(wǎng)格分別見圖3、4，均采用8 節(jié)點等參實體單元離散，面板與墊層之間，設置接觸，接觸模型為Coulomb 摩擦，摩擦系數(shù)取0.8。堆石體南水模型參數(shù)見表1，流變模型參數(shù)由長江科學院根據(jù)室內流變試驗得出，見表2，反演參數(shù)取值范圍見表3。

4.1 訓練樣本的選取

為了保證訓練樣本具有足夠的代表性，用正交試驗設計方法生成27 組樣本，同時添加一些隨機樣本以保持樣本的多樣性與均勻性，共300 組樣本，訓練樣本占總樣本數(shù)的80%，測試樣本占20%，樣本集見表4。

圖3 堆石壩有限元網(wǎng)格 Fig.3 Finite element model of CFRD

圖4 混凝土面板有限元網(wǎng)格 Fig.4 Finite element model of concrete face

表1 堆石體南水模型材料參數(shù)表 Table 1 Parameters for double-yield-surface model of rockfill

表2 堆石料冪函數(shù)流變模型試驗參數(shù) Table 2 Experimental rheological parameters of rockfill

表3 待反演參數(shù)取值范圍 Table 3 Ranges of parameters to be inversed

4.2 神經(jīng)網(wǎng)絡的訓練

對樣本數(shù)據(jù)進行歸一化處理，將網(wǎng)絡的輸入數(shù)據(jù)限制在[0， 1]區(qū)間內。采用提前終止法和基于正則化的貝葉斯方法提高網(wǎng)絡的泛化能力[10]，圖5 是神經(jīng)網(wǎng)絡訓練過程誤差圖，可見RBF 神經(jīng)網(wǎng)絡學習速度很快，到第6 次迭代時測試誤差和訓練誤差均達到設計要求。

圖5 神經(jīng)網(wǎng)絡訓練過程 Fig.5 Neural network training process

4.3 參數(shù)反演

本文采用基于粒子遷徙的粒群算法對水布埡面 板堆石壩進行參數(shù)反演分析。算法參數(shù)如下：粒群規(guī)模為50，子粒群的個數(shù) p= 4，粒群的遷徙率r= 0.4，遷徙間隔等于5，最大迭代次數(shù)取500，ωmin= 0.4，ωmax= 0.9， c1i= c2f=2.0， c1f= c2i=0.5，變異概率 pm= 0.2，iterMax = 10。目標函數(shù)中最大監(jiān)測斷面的權重系數(shù)為0.6，其余2 個較小斷面的權重系數(shù)為0.2。圖6 為反演過程中目標函數(shù)收斂過程。反演分析得到的最優(yōu)參數(shù)組合見表5。

圖6 反演計算收斂過程 Fig.6 Iteration process of parametric inversion

表5 堆石體參數(shù)反演結果 Table 5 Inversion parameters of CFRD dam

由表5 可以看出，反演得到的參數(shù)與試驗參數(shù)差別較大，主堆石的k 值較試驗值小很多，次堆石的k 值略小于試驗值，反映最終流變量的參數(shù)c、cα、cβ較試驗值增大，表明最終流變變形量增加，反映流變速率的η 、vλ 較試驗值減小，表明流變變形速率變慢，壩體達到變形穩(wěn)定所需要的時間增加。

基于反演得到的參數(shù)進行水布埡面板堆石壩的應力、變形分析，由圖7 可以看出，監(jiān)測點的沉降計算值和實測值在數(shù)值和發(fā)展規(guī)律上均吻合得較好。圖8 為壩體在基本穩(wěn)定期的沉降和水平位移等值線圖，圖9 為面板在基本穩(wěn)定期的擾度圖和軸向變形圖，基于反演參數(shù)的計算結果表明，水布埡面板壩壩體沉降速率在2012 年底降低到12 mm/a，認為壩體變形達到基本穩(wěn)定，此時壩體的最大沉降位移為2.63 m，流變變形約占壩體總變形的18%，流變變形顯著，面板的最大擾度為0.9 m。

圖7 壩體最大斷面測點沉降實測值與計算值對比 Fig.7 Comparison between measured and calculated settlements of points in the largest section of dam

圖8 變形穩(wěn)定期最大斷面位移圖 Fig.8 Predicted displacements at largest section of dam in the deformation stable phase

圖9 變形穩(wěn)定期面板變形圖 Fig.9 Predicted displacements of concrete face in the deformation stable phase

5 結 論

（1）本文采用基于粒子遷徙的粒子群算法和徑向基函數(shù)神經(jīng)網(wǎng)絡構建參數(shù)反演平臺，克服了粒子群算法易陷入局部最優(yōu)和早熟收斂的缺點，采用經(jīng)過訓練的神經(jīng)網(wǎng)絡來描述模型參數(shù)和位移之間的映射關系，節(jié)省了參數(shù)反演的計算時間，提高了反演效率。

（2）對長科院冪函數(shù)流變模型進行了參數(shù)敏感性分析，認為參數(shù)c、η 、cα、cβ、vλ 對壩體沉降較為敏感。

（3）由于堆石壩填筑過程和變形機制復雜，很難將瞬時變形與流變變形分開，因此，對靜力本構參數(shù)和流變參數(shù)進行綜合反演。采用本文構建的參數(shù)反演平臺對水布埡混凝土面板堆石壩進行了參數(shù)反演分析，基于反演參數(shù)的堆石壩應力變形分析結果表明，測點沉降計算值與實測值在數(shù)值和發(fā)展規(guī)律上均吻合地較好；在2012 年底壩體沉降速率降低到12 mm/a，認為壩體變形達到基本穩(wěn)定，在變形穩(wěn)定期壩體的最大沉降為2.63 m，約占壩高的1.12%，流變變形約占壩體總變形的18%，流變變形顯著，面板的最大擾度為0.9 m，均在合理的范圍以內。

[1] 沈珠江， 趙魁芝. 堆石壩流變變形的反饋分析[J]. 水利學報， 1998， 6： 1－6. SHEN Zhu-jiang， ZHAO Kui-zhi. Back analysis of creep deformation of rockfill dams[J]. Journal of Hydraulic Engineering， 1998， 6： 1－6.

[2] 郭雪莽， 田俊明， 秦理曼. 土石壩位移反分析的遺傳方法[J]. 華北水利水電學院學報， 2001， 22(3)： 94－98. GUO Xue-mang， TIAN Jun-ming， QIN Li-man. Displacement back analysis of embankment dam using genetic algorithm[J]. Journal of North China Institute of Water Conservancy and Hydroelectric Power， 2001， 22(3)： 94－98.

[3] 張社榮， 何輝. 改進的遺傳算法在堆石體參數(shù)反演中的應用[J]. 巖土力學， 2005， 26(2)： 182－186. ZHANG She-rong， HE Hui. Application of improved genetic algorithm to back analyzing parameters of rockfill[J]. Rock and Soil Mechanics， 2005， 26(2)： 182－186.

[4] 朱晟， 張美英， 戴會超. 土石壩瀝青混凝土心墻力學參數(shù)反演分析[J]. 巖土力學， 2009， 30(3)： 635－644. ZHU Sheng， ZHANG Mei-ying， DAI Hui-chao. Back analysis of mechanical parameters for asphalt-concrete core earth-rock dam[J]. Rock and Soil Mechanics， 2009， 30(3)： 635－644.

[5] 田俊明， 周晶. 基于蟻群算法的土石壩土體參數(shù)反演[J]. 巖石力學與工程學報， 2005， 24(8)： 1411－1416. TIAN Jun-ming， ZHOU Jing. Inversing soil mechanical parameters of embankment dam using ant colony algorithm[J]. Chinese Journal of Rock Mechanics and Engineering， 2005， 24(8)： 1411－1416.

[6] 康飛， 李俊杰， 許青. 堆石壩參數(shù)反演的蟻群聚類RBF網(wǎng)絡模型[J]. 巖石力學與工程學報， 2009， 28(增刊2)： 3639－3640. KANG Fei， LI Jun-jie， XU Qing. Ant colony clustering radial basis function network model for inverse analysis of rockfill dam[J]. Chinese Journal of Rock Mechanics and Engineering， 2009， 28(Supp.2))： 3639－3640.

[7] 李金鳳， 楊啟貴， 徐衛(wèi)亞. 基于改進粒子群算法CHPSO-DS 的面板壩堆石體力學參數(shù)反演[J]. 巖石力學與工程學報， 2008， 27(6)： 1229－1235. LI Jin-feng， YANG Qi-gui， XU Wei-ya. Back analyzing mechanical parameters of rockfill based on modified particle swarm optimization CHPSO-DS[J]. Chinese Journal of Rock Mechanics and Engineering， 2008， 27(6)： 1229－1235.

[8] 周偉， 徐干， 常曉林， 等. 堆石體流變本構模型參數(shù)的智能反演[J]. 水利學報， 2007， 38(4)： 389－394. ZHOU Wei， XU Gan， CHANG Xiao-lin， et al. Intelligent back analysis on parameters of creep constitutive model[J]. Journal of Hydraulic Engineering， 2007， 38(4)： 389－394.

[9] 張丙印， 袁會娜， 李全明. 基于神經(jīng)網(wǎng)絡和演化算法的土石壩位移反演分析[J]. 巖土力學， 2005， 26(4)： 547－552. ZHANG Bing-yin， YUAN Hui-na， LI Quan-ming. Displacement back analysis of embankment dam based on neural network and evolutionary algorithm[J]. Rock and Soil Mechanics， 2005， 26(4)： 547－552.

[10] 張德豐. MATLAB 神經(jīng)網(wǎng)絡應用設計[M]. 北京： 機械工業(yè)出版社， 2005.

[11] 沈珠江. 土體應力應變分析中的一種新模型[C]//第五屆土力學及基礎工程學術討論會論文集. 北京： 中國建筑工業(yè)出版社， 1990： 101－105.

[12] 程展林， 丁紅順. 堆石料蠕變特性試驗研究[J]. 巖土工程學報， 2004， 26(4)： 473－476. CHENG Zhan-lin， DING Hong-shun. Creep test for rockfill[J]. Chinese Journal of Geotechnical Engineering， 2004， 26(4)： 473－476.

[13] KENNEDY J， EBERHART R C. Particle swarm optimization[C]//Proceedings of IEEE International Conference on Neural Networks. Piscataway： [s. n.]， 1995： 1942－1948.

[14] CLERC M， KENNEDY J. The particle swarm-explosion， stability， and convergence in a multidimensional complex space[C]//Proceedings of IEEE Transactions on Evolutionary Computation. [S. l.]： [s. n.]， 2002， 6： 58－73.

[15] YANG X， YUAN J， YUAN J， et al. A modified particle swarm optimizer with dynamic adaptation[J]. Applied Mathematics and Computation， 2007， 189(2)： 1205－1213.