基于模糊數學和粒子群算法的邊坡參數反分析*

孫志彬 楊小禮 黃阜

(中南大學土木工程學院，湖南長沙410075)

邊坡是道路與鐵路工程中常見的構筑物形式之一，無論是對原有邊坡進行處治，還是利用現有材料構建邊坡以滿足人們生產與生活的需要，都必須進行邊坡的穩定性分析.而不管采用何種方法進行邊坡穩定性分析，邊坡的力學強度參數都非常重要.但由于邊坡內部巖土體結構復雜，內部強度隨空間變化較大，參數獲取較為困難.而模型試驗存在取樣困難、尺寸效應等問題.實測坡體位移比較簡便，因而采用位移反分析方法計算巖土體強度參數已經成為一種理想的巖土參數確定方法.

自Kavanagh等［1］提出反演彈性固體彈性模量的有限元方法以來，國內外學者對參數反演做了大量的研究.從最初的逆解法［2］、圖解法［3］逐漸發展為建立目標函數并求出最優解的優化法.優化法發展的初期，常常采用基于梯度的各種算法，如單純形法［4］、罰函數法［5］.隨著智能技術的逐漸發展，眾多學者將各種智能算法與優化法結合起來，用來尋找目標函數的最優解.常用的智能算法有小波向量機［6］、模擬退火算法［7］與遺傳算法［8］，或者將不同的智能算法結合起來，以彌補它們各自的不足［9］.

由優化法的求解過程可知，影響優化法求解的因素有兩種，即目標函數的構造和求解.大多數學者對優化法的研究集中在對目標函數的求解方法上，即如何利用更少的步驟求解出更精確的解，對目標函數的構造則研究較少.文中以浙江某滑坡為例，將模糊數學與粒子群算法結合起來，構造了一種新的參數反演的目標函數，并計算了該函數對邊坡土體強度參數的敏感度，從而為確定巖土邊坡強度參數開辟了新的思路.

1 目標函數的構造

1.1 傳統方法

利用優化法反演土體力學參數，需要構建目標函數，該函數通常為測點處位移實測值與計算值的平方和，即

式中:xi為第i個測點的實測位移值;x'i為第i個測點的計算位移值;J為待反算的邊坡強度參數矩陣.Ei、νi、ci、φi表示不同部分土體的強度參數，分別為彈性模量、泊松比、黏聚力、內摩擦角.

由公式(1)可知，隨著矩陣J中強度參數的變化，函數f也相應地變化.據此不斷調整J中的參數，使

此時矩陣J中的強度參數即為反算的結果.

傳統的觀測邊坡位移方式有兩種:一種是利用全站儀觀測邊坡表面若干測點的位移，如圖1(a)所示;另一種是利用測斜儀等儀器監測邊坡內部若干測點的位移，如圖1(b)所示.

圖1 邊坡測點位置Fig.1 Locations of the measuring points of slope

邊坡表面的位移影響因素較多，例如邊坡表面的土體崩落、雨水沖刷等.內部位移監測不僅可以有效減少干擾，還可對發生大規模滑坡前邊坡的內部滑動進行觀測.

數值分析結果表明，邊坡內部滑動面附近土體的剪應變最大.大量的觀測資料也表明，在滑坡發生前，滑動面附近的土體已經發生移動.雖然滑面土體范圍不大，但隨著時間的推移，其上部的土體將沿著滑面土體發生整體滑移，從而發生滑坡.因此滑動面附近土體對整個邊坡的穩定性起控制作用.但在發生大規模滑坡之前，該部分土體的位移從宏觀上看并不十分顯著.在邊坡參數反分析的過程中，該土體內部的測點的計算和實測位移差對式(1)中構造的優化函數f(J)的影響較小，甚至在誤差范圍內，很難對反算出的邊坡位移函數造成影響.由此可見，傳統方法將邊坡各測點位移差的平方和作為目標函數，忽略了在邊坡滑移過程中各部分土體性質的不同，無法精確地反演出邊坡的強度參數.

1.2 新的目標函數構造思路

由于邊坡內部土體變形與滑動性質的不同，將各測點位移以相同的權值進行計算存在不足之處.根據滑動特征，邊坡內部土體可分為3類:滑坡床、滑坡面和滑坡體.新的目標函數根據測點所在的位置，將其位移乘以加權系數代入式(1)，得式(2):

式中，λi為加權系數.為了確定某一測點位移的加權系數，需要確定該測點屬于哪部分的土體.但實際上，邊坡土體的精確劃分存在困難，因為各部分土體之間的界限是模糊的，在發生大規模滑動之前，滑坡體與滑坡面之間、滑坡面與滑坡床之間的轉換在空間上是漸變的.尤其對于“交界面”附近的土體，在性質上很難單獨劃分為某一類別的土體，往往“兼而有之”兩種土體的性質.可以利用模糊數學對土體的這種性質進行描述.據文獻［10］，可根據邊坡的內部位移監測數據尋找潛在滑動面，也可利用地質鉆孔資料，根據鉆出物的性質判斷潛在滑動面的位置.

設某一測孔深度為h1，滑動面深為h2.隸屬函數是以測點深度為自變量的線性函數.土體的隸屬度函數U由子函數μ1、μ2、μ3組成，分別表示該土體對滑坡床、滑坡面與滑坡體的隸屬度，如下所示:

式中:h為該測點的深度.

由函數U可知，當測點分別位于滑坡表面、滑動面和基巖表面時，其對滑坡床、滑坡面和滑坡體3種土體的隸屬度為1.離上述3點的距離愈遠，其對相對土體的隸屬度愈小，隸屬度與距離呈線性關系，函數曲線如圖2所示.

由式(1)及隸屬度函數U可知，邊坡力學參數反分析的優化法的目標函數權重可表示為

目標函數為

圖2 隸屬度函數曲線Fig.2 Membership function curve

式中:P為滑坡床、滑坡體與滑坡面的權重系數矩陣，P=［η1η2η3］，η1、η2、η3分別表示3類土體的位移在反分析過程中的權重系數;Ui為某測點的隸屬度矩陣，Ui=［μi1μi2μi3］，μi1、μi2、μi3分別表示第i個測點對3類土體的隸屬度.

欲構造目標函數，除根據測點位置計算其隸屬度外，還需求取不同土體的位移權重系數，即矩陣P.由于滑動面附近的土體發生剪切破壞，對于一些尚未發生明顯滑移的邊坡，潛在滑動面附近土體的剪切位移可以表征滑坡的狀況.同樣，在相同尺寸的計算模型中，不同的強度參數取值也將影響該部分土體的剪切位移，所以理想的矩陣P滿足以下兩個條件:一是體現不同土體在滑坡過程中的作用，一般情況下為滑坡面＞滑坡體＞滑坡床;二是使各測點的位移乘以加權系數后對滑動面的位移有足夠大的敏感性.

現以浙江某滑坡為例，利用粒子群算法說明矩陣P的構造思路.該邊坡為土石混合體邊坡，邊坡形狀與測點布置如圖3所示，每個測斜孔長15.0m，每0.5m取一個測點，各個測點的位移如圖4所示.

圖3 邊坡形狀及測點位置Fig.3 Slope shape and location of measuring points

圖4 監測點的實際位移Fig.4 Actual displacements of monitoring points

式中，x'j(P)為加權后的測點位移，x'ij為數值計算得出的各測點位移，j為計算次數.為了保證邊坡各測點位移對滑動面的剪切位移的敏感性，要求

式中，x'1(P)、x'2(P)、x'3(P)分別為強度參數為時的測點位移值，maxs1m、maxs2m、maxs3m為相應強度條件下的邊坡最大剪切位移值.

若考慮不同類別土體在滑坡過程中的重要性:

上述問題即為在滿足式(10)的條件下求式(9)，可采用粒子群搜索算法.粒子群算法［11-12］是一種基于群體方法的演化計算技術，模仿一群飛鳥在一定區域中搜索食物的方法.鳥群的搜索空間對應于設計變量的變化范圍，食物對應于適應度函數的最優解.每個優化問題的潛在解對應于搜索空間中的一個粒子，每個粒子都有自己的位置和速度(決定飛行的方向和距離)，還有一個由被優化函數決定的適應值.各個粒子記憶、追隨當前的最優粒子，在解空間中搜索.每次迭代的過程不是完全隨機的，如果找到較好的解，將會以此為依據尋找下一個解.

在本問題中，粒子種群規模為10，最大迭代次數為15，其他參數取為慣性權重w=0.9，學習因子e1=e2=2.在計算過程中，若某個粒子的一個以上的分量的取值超過解空間的范圍，則不計算其適應度，并且失去同別的粒子競爭最優位置的機會，直至其被重新吸引至解空間的范圍內.

經過15次迭代以后，求得矩陣

2 目標函數對強度參數的敏感性分析

為了驗證新的目標函數對所要反演的邊坡強度參數的敏感性，文中特別計算了圖3所示的邊坡參數改變時的目標函數值(見表1和2)，并做出了目標函數對參數的敏感度曲線.為便于比較，同時計算出按傳統的優化方法反分析所構造的目標函數的函數值，并繪出敏感度曲線，結果如圖5所示.

表1 c值改變時目標函數的變化值Table 1 Objective function values varying withc

表2 φ值改變時目標函數的變化值Table 2 Objective function values varying withφ

圖5 目標函數對c值和φ值的敏感度Fig.5 Sensitivity of objective function to c and φ

表1和2中，邊坡強度參數的設計值為c= 29kPa、φ=30°，設F(J)為傳統目標函數，H(J)為新的優化法目標函數，為了便于說明，只反演參數矩陣中的黏聚力與內摩擦角值.F(J)與H(J)可寫為F(c，φ)和H(c，φ).分別計算目標函數對c與φ敏感度變化.

目標函數對c的敏感度定義為

敏感性分析結果表明，邊坡位移對c的敏感度要大于對φ的敏感度.就兩個目標函數而言，無論是對于c還是φ，文中提出的目標函數的敏感度大多數情況下大于傳統方法的敏感度.這說明在利用優化法進行邊坡強度參數反分析時，文中的方法具有優勢.利用目標函數進行求解，得目標邊坡的強度參數為c=30.5kPa，φ=31.3°.

文中提出的算例為單一滑動面的邊坡，對于多層土體且僅有單個滑動面的邊坡，亦可以使用文中提出的方法.但若邊坡有若干次生滑動面，則需要重新考慮各個測點的隸屬度計算方法.

3 結語

傳統邊坡參數反分析的目標函數忽視各測點的位置與位移狀態，沒有反映邊坡的內部位移情況.因此，文中利用邊坡內部的監測資料，對邊坡位移進行加權后再構造目標函數.研究發現通過新的目標函數可以有效地求解邊坡的強度參數.層狀土質邊坡存在土層厚度不同造成的尺寸效應，利用文中提出的方法存在一定的局限性，在構建目標函數時，可以尋找一個與土層厚度無關的位移量.

［1］ Kavanagh K，Clough R.Finite element applications in the characterization of elastic solids［J］.International Journal of Solids and Structures，1971，7(1):11-23.

［2］ Sakurai S，Takeuchi K.Back analysis of measure displacement of tunnel［J］.Rock Mechanics and Rock Engeering，1983，16(3):173-180.

［3］ 楊林德.巖土工程問題的反演理論與工程實踐［M］.北京:科學出版社，1996:120-123.

［4］ Gioda G.Problem in elastic-plasticity identification of cohesion，friction angle and in-situ stress by pressure tunnel tests［J］.International Journal for Numerical Methods in Engineering，1980，15(5):1823-1834.

［5］ 王華寧，呂愛鐘.巷道裂隙巖體的損傷參數辨識［J］.巖土工程學報，2001，23(5):593-597.Wang Hua-ning，Lü Ai-zhong.Damage parameter identification of joint rock mass in underground opening［J］.Chinese Journal of Geotechnical Engineering，2001，23 (5):593-597.

［6］ 王復明，李曉龍，苗麗，等.基于小波支持向量機的圍巖力學參數辨識［J］.水力發電學報，2010，6(3):184-190.Wang Fu-ming，Li Xiao-long，Miao Li，et al.Mechanical parameters identification of surrounding rock based on wavelet SVM［J］.Journal of Hydroelectric Engineering，2010，6(3):184-190.

［7］ 魏連偉.基于模擬退火算法的水文地質參數識別［J］.天津大學學報，2003，36(5):618-621.Wei Lian-wei.Hydrogeological parameter identification on the simulated annealing genetic algorithm［J］.Journal of Tianjin University，2003，36(5):618-621.

［8］ 高瑋，鄭穎人.基于遺傳算法的巖土本構模型辨識［J］.巖石力學與工程學報，2002，21(1):9-12.Gao Wei，Zheng Ying-ren.Identification of the geomaterial constitutive model based on genetic algorithm［J］.Chinese Journal of Rock Mechanics and Engineering，2002，21(1):9-12.

［9］ 劉鵬程，紀晨.改進模擬退火－單純形綜合反演方法［J］.地球物理學報，1995，38(2):199-205.Liu Peng-cheng，Ji Chen.An improved simulated annealing-downhill simplex hybrid global inverse algorithm［J］.Chinese Journal of Geophysics，1995，38(2):199-205.

［10］ 孫志彬，楊小禮.基于深部位移的邊坡滑動特征分析［J］.長沙理工大學學報:自然科學版，2010，7(2):43-47.Sun Zhi-bin，Yang Xiao-li.Sliding characteristics analysis of land slide based on the deep displacement［J］.Journal of Changsha University of Science and Technology:Natural Science，2010，7(2):43-47.

［11］ Kennedy J，Eberhart R C.Proceedings of IEEE international conference on neural network［C］.New York: IEEE Press，1995:1942-1948.

［12］ 高瑋.基于粒子群優化的巖土工程反分析研究［J］.巖土力學，2006，27(5):795-798.Gao Wei.Back analysis algorithm in geotechnical engineering based on particle swarm optimization［J］.Rock and Soil Mechanics，2006，27(5):795-798.