邊坡動力破壞機理的振動臺試驗研究

劉振平，劉福海，趙顯波

（1.黑龍江工程學院 土木與建筑工程學院，哈爾濱 150050；

2.東北電力設計院，長春130021；3.黑龍江省水利科學研究院，哈爾濱 150080）

在土木工程和水利工程中，有眾多的自然邊坡、人工邊坡如堤壩、路基路塹等邊坡工程。中國的許多重大工程如大型水電站、南水北調西線工程、西部山區高速公路等均處于強震區，特別是2008年汶川“5·12”地震以來，邊坡動力破壞機理引起學者們的廣泛關注，成為巖土地震工程領域中重要的研究課題之一。有關邊坡地震的原型觀測資料很少，所以室內振動臺試驗就成為研究地震作用下邊坡動力問題的重要手段之一。Lin等［1］采用大型振動臺試驗研究土坡在地震作用下的響應，認為砂土邊坡破壞面較淺，只是坡體表面的破壞。徐光興等［2－3］通過振動臺模型試驗研究了混合土邊坡的動力特性，并結合數值分析結果認為在強震作用下單一均質土層邊坡的破壞模式仍然是沿著某一弧形潛在滑動面失穩。許強等［4］通過多組振動臺試驗認為均質土坡和塊狀巖質斜坡的動力破壞模式都是坡頂拉裂，中下部剪切滑移破壞。陳新民等［5－6］采用振動臺試驗研究了粘土邊坡動力特性、動力反應和宏觀變形，發現在動力作用下邊坡坡頂和坡腳出現裂縫。葉海林等［7］巖質邊坡大型振動臺模型試驗，結果表明地震滑動面為上部拉裂縫和下部剪切滑移面形成貫通的破裂面。劉君等［8］把PIV技術應用在土質邊坡振動臺模型試驗中，獲取了邊坡破壞的完整過程。Wang等［9］在砂土邊坡振動臺試驗中采用PIV技術獲取了邊坡表面的位移場，探討了邊坡動力失穩機理。Wartman等［10－11］通過多組粘土邊坡振動臺模型試驗研究了邊坡產生的永久位移，并與基于峰值和殘余強度的Newmark公式計算結果進行了比較。Katz等［12］采用振動沙箱研究了砂土邊坡破壞的類型和頻率的關系。以上針對邊坡的動力特性、宏觀破壞現象及個別點的位移進行的研究，而沒有涉及邊坡的位移場，劉君、Wang等雖用PIV技術得到了整個邊坡的位移場，但沒有對應變場進行研究。

應變是直接對含有噪聲的離散位移進行差分計算得到的，那么即使微小的位移測量誤差也會被急劇地放大，使計算的應變不可靠［13］。一般的數據擬合平滑方法也不能滿足應變計算的精度要求，潘兵等［13－14］提出了局部最小二乘擬合的應變計算方法，其原理同Savitzky-Golay平滑方法。該方法能夠在一定程度上去除原始位移場的噪聲，缺點是由于每個位移點都要一個計算窗口平滑，計算量巨大，并且對于非均勻應變場，應變計算窗口的大小顯著影響應變的計算精度，目前還沒有合適的方法來選擇最佳的窗口大小。本文采用數字圖像位移和有限元數據平滑技術，得到了整個邊坡模型的位移場和應變場，同時探討了土質邊坡的地震破壞機理。

1 基于有限元數據平滑的應變計算

數字圖像位移測量技術現在已經成熟，采用互相關圖像匹配技術獲取位移場后，再用有限元數據平滑方法對位移場進行處理后，得到應變場。

1.1 基本原理

基本思路是引入泛函f（F）［15－19］，使該泛函取得極小值的自變函數F就是要求的光滑位移函數。

根據試驗數據點的位置分布，把區域A劃分為若干單元，在每一個單元內用形函數逼近光滑函數F，這個過程與有限元方法是一樣的。

在第j個單元內，把光滑函數F表示成式中：N j為第j個單元的形函數；uej為第j個單元的節點自由度。

式（9）是有限元的標準格式，其中K＝K1＋λK2，λ是平滑系數，U表示所有單元的節點自由度。平滑系數的確定方法在下文專門論述。

只要確定了平滑系數，就可從式（9）中解出U，再從式（2）中得到平滑位移函數u，其應變的計算方法如下：

有限變形條形下的格林（Green）應變［20］（以受壓為正）為

式中：u和v分別為上述所求的x和y方向的位移函數；εx和εy分別為x和y方向的正應變；γxy為工程剪應變；γmax為最大剪應變。

小變形條件下柯西（Cauchy）應變的計算可不考慮上述應變表達式中的二次項。

1.2 形函數的確定

采用三角形單元，以便適應任意形狀的數據區域。同時為提高插值函數的光滑可導性，采用的是五次多項式的單元插值函數［15］，其表達式為式中：X為變量項；C為系數列向量。由于沿三角形單元三條邊滿足法向導數連續的條件，21個系數可縮減為18個。

1.3 平滑系數的確定

為了能夠自動得到最佳的平滑系數，Golub等［21］、Craven等［22］及 Bates等［23］提 出 了 廣 義 交 互驗證GCV方法，建立了關于λ的GCV函數X是n階滿秩陣，K2是奇異陣，其秩為r。擬總體剛度陣K1是對稱正定陣，因此式（15）成立，可以證明X是滿秩陣，因此y也能夠由式（16）式唯一確定。

求取廣義交互驗證函數GCV的最小值，就可以得到最佳平滑系數λ，具體計算方法詳見文獻［24］。

2 試驗設備及模型設計

2.1 試驗設備

試驗在大連理工大學工程抗震試驗室的水平與豎直雙向水下振動臺上進行，其主要性能參數如下：數字式控制方式，臺面尺寸3 m×4 m，最大載重100 k N，工作頻率0.1～50 Hz；滿載工況下最大加速度±1.0g（水 平 向）、±0.7g（豎 向）、最 大 速 度±50 cm／s（水平向）、±35 cm／s（豎向）、最大位移±75 mm（水平向）、±50 mm（豎向）。

2.2 試驗模型

在鋼制模型箱內堆制模型進行試驗，為便于觀察和圖像采集，一側為有機玻璃。模型壩壩高1.3 m，邊坡坡率1∶1.6，其模型如圖1所示。模型填土基本性能參數見表1。

表1 試驗填土基本性能參數

圖1 邊坡振動臺模型示意圖（單位：cm）

試驗加載的是頻率10 Hz的正弦增幅波（喇叭波），如圖2所示。從0開始逐漸增加，直至邊坡完全失穩。

采用的圖像采集設備是CANON EOS 450D高清數碼相機，最大分辨率4 272×2 848像素，配備10～22 mm廣角鏡頭，可以在距模型2 m處拍攝整個模型的圖片。存儲為JPG格式的圖片文件，最大連拍速度約3.3幀／s。振動過程中連續的采集圖像，用于數字圖像位移和應變識別。

圖2 輸入的加速度時程

3 試驗結果分析

3.1 模型的位移應變計算流程

開發了一套數字圖像測量分析程序，可以批處理多種格式的圖片，可以得到一系列的位移場和應變場，還能以等值線圖、云圖和矢量圖等形式輸出各種位移場、應變場和速度場。其計算流程和步驟如下：

位移計算流程：圖像像素塊劃分→計算總位移→計算模型箱位移→計算相對位移（總位移減去模型箱剛體位移）→刪除計算的位移壞點→輸出位移場等后處理。

得到位移場后，用有限元平滑位移場再微分計算應變場。其計算流程如下：生成單元和節點信息→計算單元剛度陣和單元節點荷載列陣→集成總體剛度陣和總體節點荷載列陣→計算平滑系數→求解節點未知參數→平滑位移場微分求應變→輸出應變場等后處理。

3.2 邊坡模型的破壞過程描述

邊坡模型破壞過程基本可分3個階段，第1階段是整體變形階段（圖3（a）），當輸入加速度較小時，整個模型以均勻沉降變形為主，同時伴有相對于模型箱的水平往復運動；第2階段是滑移變形階段（圖3（b）），當輸入加速度大約在0.1g時，模型開始沿邊坡一側向下滑移，隨著輸入加速度的增大，滑移量也逐漸增大，在輸入加速度在0.45g時，坡頂出現張拉裂縫，坡體內部也出現剪切裂隙，形成明顯的滑動帶，邊坡達到極限破壞狀態；第3階段是破壞階段，滑動帶上的裂隙使得滑動體與邊坡相對分離，滑動體加速下滑，最后整個坡體坍塌，邊坡完全破壞。

從圖3（c）可以看出，邊坡模型動力失穩存在深層滑動帶，滑動面大概呈圓弧狀。

3.3 位移分析

圖4、圖5和圖6分別是8.42 s（對應輸入加速度為0.45g）時即在極限破壞狀態時的水平位移、豎向位移和總位移的等值線圖，其水平位移出現在坡面中部，豎向位移在坡頂附近，水平位移和總位移等值線基本平行于坡面，均符合一般的邊坡破壞時的位移分布規律。

圖3 邊坡模型位移矢量

圖4 8.42 s時（輸入加速度0.45 g）時水平位移場（單位：mm）

圖5 8.42 s時（輸入加速度0.45 g）時豎向位移場（單位：mm）

圖6 8.42 s時（輸入加速度0.45 g）時總位移場（單位：mm）

3.4 應變分析

從圖7和圖8最大剪應變和水平方向應變分布上可以看出，坡體中部到坡腳是剪切破壞，坡頂一定深度是拉剪破壞。其變形破壞模式與汶川地震觀察到的邊坡破壞現象相吻合［25］：土坡主要是坡頂向下一定深度內的拉破壞，坡腳向上延伸形成剪切滑移帶，最終二者連通形成貫通的破裂面。這也與許強等［4］的振動臺試驗結果以及鄭穎人等［26－27］用 FLAC程序做的數值分析結果基本一致。

圖7 8.42 s時（輸入加速度0.45 g）時水平應變場（單位：%）

8 8.42 s時（輸入加速度0.45 g）時最大剪應變場（單位：%）

3.5 位移應變時程分析

A點位于滑動體上，B點位于滑動體外（圖1），從圖9的位移時程可以看出，位移隨時間是逐漸增加的，沒有明顯的突變點，說明邊坡在地震作用下是漸進式破壞的。

圖9 位移時程曲線

由于位移時程曲線沒有明顯的突變拐點，如何判斷模型破壞的時刻就成為一個難題。地震作用下邊坡形成滑移帶后，滑動體與坡體相對分離，由于滑移帶的阻隔，滑動體以外的壩體振動對其影響將大大減小，滑動體將比較平穩的運動，而不會隨輸入地震波而振動。所以可以把位移時程曲線的曲率做為判斷壩坡失穩滑移的一個物理量，這里所指的曲率是有正負之分的，正號表示是凹曲線，負號表示是凸曲線，可定義為廣義曲率。

式中：κ為廣義曲率；y′是位移的一階導數；y″是位移的二階導數。

從A點的位移時程曲線曲率（圖10）可以看出，在8.42 s之前，曲率以零線為中心上下波動，說明在此之前沒有發生連續的滑動，而在此之后，曲率迅速趨于零，不再發生波動，該點此時已經失穩破壞，曲率發生突變表明滑移帶已經貫通形成，邊坡整體失穩，達到了極限抗震狀態。位于滑動體以外的B點，曲率自始自終都處于振動狀態，沒有發生失穩破壞。由此可知，用位移時程曲線曲率做為判斷邊坡動力破壞的物理量是可行的。

圖10 位移時程曲線的曲率

4 結 論

將有限元數據平滑方法引入到數字圖像位移應變測量中是可行的，得到了邊坡振動臺模型整個試驗過程的位移場和應變場。初步表明土質邊坡的變形是漸進式的，坡體中部到坡腳是剪切破壞，坡頂一定深度是拉剪破壞，破壞時有深層的圓弧狀滑動面，用位移時程的廣義曲率做為判斷邊坡動力破壞的物理量是可行的。

［1］Lin M L，Wang K L.Seismic slope behavior in a largescale shaking table model test ［J］.Engineering Geology，2006，86（2／3）：118-133.

［2］徐光興，姚令侃，高召寧，等.邊坡動力特性與動力響應的大型振動臺模型試驗研究［J］.巖石力學與工程學報，2008，27（3）：624-631.

Xu G X，Yao L K，Gao Z N，et al.Large-scale shaking table model test study on dynamic characteristics and dynamic responses of slope［J］.Chinese Journal of Rock Mechanics and Engineering，2008，27（3）：624-631.

［3］徐光興，姚令侃，李朝紅，等.邊坡地震動力響應規律及地震動參數研究［J］.巖土工程學報，2008，30（6）：918-923.

Xu G X，Yao L K，Li C G，et al.Dynamic response of slopes under earthquakes and influence of ground motion parameters［J］.Chinese Journal of Geotechnical Engineering，2008，30（6）：918-923.

［4］許強，陳建君，馮文凱，等.斜坡地震響應的物理模擬試驗研究［J］.四川大學學報：工程科學版，2009，41（3）：266-272.

Xu J，Chen J J，Feng W K，et al.Study of the seismic response of slopes by physical modeling ［J］.Journal of Sichuan University：Engineering Science Edition，2009，41（3）：266-272.

［5］陳新民，沈建，魏平，等.下蜀土邊坡地震穩定性的大型振動臺試驗研究（I）-模型試驗設計［J］.防災減災工程學報，2010，30（5）：492-502.

Chen X M，Shen J，Wei P，et al.Large-scale shaking table test of seismic stability of Xiashu loess slopes（I）：design of model test［J］.Journal of Disaster Prevention and Mitigation Engineering，2010，30（5）：492-502.

［6］陳新民，沈建，魏平，等.下蜀土邊坡地震穩定性的大型振動臺試驗研究（II）-試驗結果及分析［J］.防災減災工程學報，2010，30（6）：587-594.

Chen X M，Shen J，Wei P，et al.Large-scale shaking table test of seismic stability of Xiashu loess slopes（II）：analysis of test results ［J］.Journal of Disaster Prevention and Mitigation Engineering，2010，30（6）：587-594.

［7］葉海林，鄭穎人，杜修力，等.邊坡動力破壞特征的振動臺模型試驗與數值分析［J］.土木工程學報，2012，45（9）：128-135.

Ye H L，Zheng Y R，Du X L，et al.Shaking table model test and numerical analysis on dynamicfailure characteristics of slope［J］.China Civil Engineering Journal，2012，45（9）：128-135.

［8］劉君，劉福海，孔憲京，等.PIV技術在大型振動臺模型試驗中的應用［J］.巖土工程學報，2010，32（3）：368-374.

Liu J，Liu F H，Kong X J，et al.Application of PIV in largescale shaking table model tests ［J］.Chinese Journal of Geotechnical Engineering，2010，32（3）：368-374.

［9］Wang K L，Lin M L.Initiation and displacement of landslide induced by earthquake-a study of shakingtable model slope test［J］.Engineering Geology，2011，122：106-114.

［10］Wartman J.Physical model studies of seismically induced deformations in slopes ［D］.California：University of California-Berkeley，1999.

［11］Wartman J，Seed R B，Bray J D.Shaking table modeling of seismically induced deformations in slopes ［J］.Journal of Geotechnical and Geoenvironmental Engineering，2005，131（5）：610-622.

［12］Katz O，Aharonov E.Andslides in vibrating sand box：what controls types of slope failure and frequencymagnitude relations？［J］.Earth and Planetary Science Letters，2006，247（3）：280-294.

［13］潘兵，謝惠民.數字圖像相關中基于位移場局部最小二乘擬合的全場應變測量［J］.光學學報，2007，27（11）：1980-1986.

Pan B，Xie H M.Full-field strain measurement based on least-square fitting of local displacement for digital image correlation method［J］.Acta Optica Sinica，2007，27（11）：1980-1986.

［14］Pan B，Xie H M，Guo Z Q，et al.Full-field strain measurement using a two-dimensional Savitzky-Golay digital differentiator in digital image correlation ［J］.Optical Engineering，2007，46（3）：1-10.

［15］Segalman D J，Woyak D B，Rowlands R E.Smooth spline-like finite-element differentiation of full-field experimental data over arbitrary geometry ［J］.Experimental Mechanics，1979，19（12）：429-437.

［16］Engelstad M J，Chambless D A，Swinson W F，et al.Hybrid stress analysis of vibrating plates using holographic interferometry and finite elements ［J］.Experimental Mechanics，1987，27（1）：23-30.

［17］Feng Z，Rowlands R E.Continuous full-field representation and differentiation of three-dimensional experimental vector data［J］.Computers and Structures，1987，26（6）：979-990.

［18］Sutton M A，Turner J L，Burck H A，et al.Full-field repersentation of discretely smapled surface deformational for displacement and strain analysis ［J］.Experimental Mechnaics，1991，31（2）：168-177.

［19］Freese C E，Gee L.A multilevel treatment of moiré fringe data using finite elements ［J］.Experimental Mechanics，1999，39（4）：304-310.

［20］陸明萬，羅學富.彈性理論基礎［M］.北京：清華大學出版社，2001：38-43.

［21］Golub G H，Heath M，Wahba G.Generalized cross validation as a method for choosing a good ridge parameter［J］.Technometrics，1979，21（2）：215-223.

［22］Craven P，Wahba G.Smoothing noisy data with spline functions：estimating the correct degree of smoothing by the method of generalized cross-validation ［J］.Numerische Mathematik，1979，31：377-403.

［23］Bates D M，Lindstrom M J，Wahba G，et al.Gcvpak-routines for generalized cross validation ［J］.Communications in Statistics-Simulation and Computation，1987，16（1）：263-297.

［24］Kent J T，Mohammadzadeh M.Global optimization of the generalized cross-validation criterion ［J］.Statistics and Computing，2000，10（3）：231-236.

［25］許強，董秀軍.汶川地震大型滑坡成因模式［J］.中國地質大學學報：地球科學，2011，36（6）：1134-1142.

Xu J，Dong X J.Genetic types of large-scale landslides induced by Wenchuan Earthquake［J］.Journal of China University of Geosciences：Earth Science，2011，36（6）：1134-1142.

［26］鄭穎人，葉海林，黃潤秋.地震邊坡破壞機制及其破裂面的分析探討［J］.巖石力學與工程學報，2009，28（8）：1714-1723.

Zheng Y R，Ye H L，Huang R Q.Analysis and discussion of failure mechanism and fracture surface of slope under earthquake ［J］.Chinese Journal of Rock Mechanics and Engineering，2009，28（8）：1714-1723.

［27］鄭穎人，葉海林，黃潤秋，等.邊坡地震穩定性分析探討［J］.地震工程與工程振動，2010，30（2）：173-180.

Zheng Y R，Ye H L，Huang R Q，et al.Study on the seismic stability analysis of a slope ［J］.Journal of Earthquake Engineering and Engineering Vibration，2010，30（2）：173-180.

（編輯王秀玲）