基于隨機場理論的水電工程渣料場邊坡穩定性分析

于曉天，曾亞武

(武漢大學土木建筑工程學院，武漢 430072)

0 引 言

隨著西部大開發的持續推進，近年來我國在西部山區興建了大量的水利水電工程。在這些水利水電工程的建設過程中不可避免地需要對原地貌進行相應的改造，從而產生大量的渣料，堆積后便形成渣料場。渣料場作為特殊的人工堆積邊坡，其均勻性和穩定性受到堆填歷史、堆填方式、降雨、工程擾動等諸多因素的影響，一旦失穩不僅有可能造成嚴重的生命財產損失，而且有可能造成嚴重的生態環境災難。在開展渣料場邊坡穩定性分析時，渣體抗剪強度參數的獲取比較困難，一方面通過試驗或經驗取值都難以準確反映渣體的真實狀態，即存在參數測不準的問題；另一方面，即使針對具體的渣體測得準確的抗剪強度參數，也只能代表測試試樣的結果，不能代表整個堆渣體的抗剪強度參數，即堆渣體的抗剪強度參數存在空間變異性的問題。對于前者，大多數的研究以可靠度分析來考慮抗剪強度參數的變異性，屬于隨機變量模型范疇，但其前提是假定研究范圍內抗剪強度參數是均勻分布的；而后者考慮的是空間不同點處局部與整體巖土物理力學性質之間的差異性，屬于隨機場理論的范疇。顯然采用隨機場理論來研究渣料場抗剪強度參數的空間變異性對其穩定性的影響，更加符合渣料場的工程實際，具有重要的工程意義。

目前，國內外學者利用隨機場理論對非均質邊坡穩定性問題進行了諸多有益的探究。如閆澍旺等[1]利用天津港現場的勘察資料首先建立了隨機場模型，再利用方差折減技術分析了邊坡的穩定性。祁小輝等[2]利用譜表現法建立隨機場模型，并采用非侵入式有限元分析方法，分析了邊坡可靠度。在國外，Griffiths和Fenton[3]基于隨機場理論和有限元強度折減法提出了邊坡可靠度分析的隨機有限元方法。Low[4]利用FORM方法研究了抗剪強度參數的二維各向異性空間變異性對挪威南部黏土邊坡失效概率的影響。Cho[5]利用Karhunen-Loeve(簡稱K-L)級數展開方法模擬土體參數空間變異性和MCS方法計算邊坡失效概率。

顯然，上述研究成果大多是圍繞自然形成的非均質邊坡穩定性進行的研究，相較于自然巖土體，渣料存在松散性明顯、抗剪強度較低、降雨易入滲等特點。本文將針對水電工程渣料場這一特殊的人工邊坡的穩定性問題，利用隨機場理論，采用喬列斯基分解技術和非侵入式隨機有限元法，研究渣料場抗剪強度參數的空間變異性對其穩定性的影響，為渣料場防護設計提供參考。

1 渣料參數空間變異性的模擬

1.1 自相關函數

在采用隨機場理論分析巖土體參數空間變異性時，常用自相關函數表征巖土體中任意兩點參數的自相關性，構建相應的隨機場。

由于實際工程中得到的實際測量數據往往十分有限，我們假設隨機場遵循平穩性假設，并且利用理論自相關函數來說明其空間相關性。蔣水華等[6]對5種理論自相關函數對邊坡穩定性可靠度的影響做出了詳細的說明，認為自相關函數類型對邊坡穩定性研究的影響不大。本文采用形式最為簡單，應用最為廣泛的指數型自相關函數來表征渣料場抗剪強度參數在空間任意兩點的相關性，其計算公式為：

(1)

式中：lx和ly分別為水平和垂直方向的相關距離，用于表征參數空間變異性的相關程度，相關距離越大，表示參數的空間相關性越強。

1.2 改進的喬列斯基分解法

采用隨機場理論模擬巖土體參數空間變異性時，最關鍵的一步是將隨機場離散為一組隨機變量。

喬列斯基分解是近些年發展起來的一種離散隨機場的有效方法[7]。相比較而言，利用喬列斯基分解技術離散相關非高斯隨機場，計算過程簡單、編程容易實現，且離散得到的隨機變量數量等于隨機場有限單元數量，物理意義直觀[8]。然而，喬列斯基分解技術只能對單參數隨機場進行離散，為此李典慶等[9]提出了改進喬列斯基分解技術用于相關多參數非高斯隨機場的離散，取得了一系列成果，下面簡單介紹改進喬列斯基分解技術的相關非高斯隨機場的模擬步驟。

(1)劃分隨機場有限單元，得出每個隨機單元中心點的坐標Qi=(xi,yi)，i=1, 2, …，n,n為隨機場有限單元的數目。

(2)拉丁超立方抽樣技術具有較好的一維投影和均勻分層分布的特性，樣本能均勻地覆蓋到變量概率分布的尾部，因此利用拉丁超立方抽樣配點法隨機產生一組每列含n個隨機樣本維度的獨立標準正態空間隨機樣本向量ζ。

(2)

式中：θ為隨機場的波動范圍。

最后利用概率變換可以得到相關非高斯隨機場。

由上述內容可知，相關非高斯隨機場特征值的每個分量不僅取決于參數的自相關函數、波動范圍、互相關系數，而且與有限單元的中心點坐標一一對應。

2 非侵入式隨機有限元法

非侵入式隨機有限元分析的最大特點在于確定性分析以及隨機分析過程相互獨立。確定性分析直接由現有的有限元計算軟件執行，無需修改軟件有限元源代碼，直接建立邊坡可靠度分析與通用有限元軟件的接口，使得計算效率大大提高，能更快、更好地揭示巖土體參數空間變異性對渣料場穩定性影響的規律。本文以Geo-studio軟件為例說明其主要的操作步驟如下：

(1)以渣料場抗剪強度參數均值在Geo-SIGMA/W模塊中建立渣料場有限元分析模型，劃分有限元網格，設置邊界條件，在SLOPE/W模塊中建立穩定性分析模型，并將SIGMA/W模塊設為SLOPE/W的母模塊。將邊坡穩定性有限元分析模型另存為擴展名為xml的計算源文件，并提取出有限單元中心點的坐標。

(2)利用拉丁超立方抽樣技術在獨立標準正態空間中進行一次抽樣，并利用改進喬列斯基分解方法模擬參數隨機場，得到隨機場樣本值的一次實現。將樣本值分別代替xml文件中相應的單元中心點上的參數均值，得到一個新的邊坡穩定性分析文件。

(3)利用批處理軟件直接調用Geo-studio軟件，對新生成的xml文件進行穩定性分析。

(4)計算完成后會在當前文件夾中自動生成一個相應的計算結果文件，從中可以提取出相應的臨界安全系數以及滑動面等信息。

(5)重復步驟(2)～(4)n次，計算得到相應的n個臨界安全系數，以及相應的滑動面等信息。

計算流程如圖1所示。

圖1 計算流程圖Fig.1 Flow diagram of calculatation

3 算 例

3.1 基本情況

為驗算本文方法的可行性和正確性，本文以文獻[6]中土質邊坡為算例進行分析。邊坡坡高H為10 m，坡度為1∶1，對應的土體抗剪強度參數及其空間變異系數如表1所示，抗剪強度參數符合對數正態分布，自相關函數采用指數型自相關函數。采用邊長為0.5 m的四邊形四節點二維實體結構等參單元，在坡面附近過度為三角形單元，將邊坡剖分為1 210個有限元單元，如圖2所示。采用理想彈塑性本構模型和摩爾庫倫屈服準則，利用Geo-slope中“進入和退出”的滑移面選擇方法計算該邊坡的臨界安全系數。

表1 邊坡土體統計參數取值表Tab.1 Statistical properties of soil parameters

圖2 算例邊坡有限單元剖分示意圖Fig.2 FEM model of slope

3.2 隨機場實現

波動范圍是表征巖土體參數空間自相關性的重要組成部分。顯而易見，空間任意兩處巖土體特性參數之間存在一定的差異性又存在一定的相關性，土體參數自身之間的這種自相關性一般隨著兩點間距離的增大而逐步減少，當距離大于某一定值后，自相關性可忽略不計，我們將這一定值稱為波動范圍。根據Vanmarcke[10]，波動范圍的計算公式為：

(3)

式中：ρ(τ)為自相關函數；τ為空間任意兩點間的相對距離；L為隨機場區域計算尺寸。

由此可見，波動范圍表示的是自相關函數與坐標軸所夾的面積。

本文針對算例邊坡，給出了4種不同水平與豎向波動范圍組合下的抗剪強度參數隨機場，如圖3所示。圖3中每一種隨機場模型都可視為一種可能的實際情況，每一個單元的中心參數視為本單元的參數，以不同的顏色區分參數的大小。很明顯，從圖3中可以看到所得計算范圍內的抗剪強度參數是不均勻的，且不同的波動范圍組合所得到的隨機場參數分布也不同，體現了抗剪強度參數的變異性及模擬方法的隨機性。從圖3中還可以看出，當δx=10δy時，抗剪強度參數高強度帶與低強度帶近似成條狀分布；當δy=δx時，高強度帶與低強度帶近似呈塊狀分布。此外，波動范圍越大，高強度帶(或低強度帶)越連續。一般來說，巖土體參數在空間上的各向異性，由于受地層沉積作用、土體堆積作用等影響，水平方向的相關距離大于垂直方向的相關距離[11]，如水平方向的相關距離變化范圍一般為10～40 m，而垂直方向的相關距離變化范圍為1～3 m。

圖3 不同的豎直向和水平向波動范圍組合下的隨機場實現Fig.3 Realizations of random fields for various com-binations of horizontal and vertical fluctuation scales

考慮到渣料場為人工逐漸堆填的，具有較明顯的成層特征，同時也考慮渣料本身構成由不均勻的渣體組成，因此本文的分析中取 ，水平方向波動范圍取40 m，垂直方向波動范圍取4 m。

3.3 穩定性分析結果

在土體參數為均值時，如圖2所示的邊坡穩定性安全系數為1.286。在考慮邊坡土體抗剪強度參數空間變異性的情況下，對該邊坡抗剪強度參數的空間分布進行了1 000次隨機模擬，相應地計算得出1 000個可能的安全系數及其滑動面位置。將該1 000個安全系數進行平均，均值為1.290，與邊坡土體抗剪強度參數取均值時的安全系數接近。

計算結果表明，在考慮抗剪強度參數空間變異性的情況下，算例邊坡的安全系數概率密度近似呈正態分布，其中最小安全系數為0.985，最大值為1.902，是抗剪強度參數空間隨機組合出現的極端情況，出現概率極低，大部分計算結果介于1.1～1.4之間，如圖4所示。

圖4 邊坡安全系數概率分布圖Fig.4 Probability distribution chart of safety factor

由圖4還可以看出，在考慮土體抗剪強度參數空間變異性時，53.3%的邊坡實際安全系數小于邊坡土體參數為均值時的安全系數，極端情況下邊坡甚至存在破壞的可能。如果以1.25作為穩定安全系數，則不考慮抗剪強度參數的空間變異性時，該邊坡是穩定性的(1.286)；而考慮抗剪強度參數的空間變異性后，安全系數小于1.25的概率為44%，即確保該邊坡穩定的概率只有56%。

本算例的分析結果說明采用本文方法能夠考慮土體抗剪強度參數的空間變異性對邊坡穩定性的影響，在實際工程中考慮土體抗剪強度參數空間變異性影響也是必要的。

4 工程實例

4.1 工程概況

某水電工程渣料場位于該工程下游河道左側，該渣料場在原設計渣料場基礎上擴容形成，總渣料方量達到60萬方。該渣料場共分為4～5個平臺，平臺寬10～30 m不等。由此可見該渣料場范圍大，堆填不規則。為確保該地區生態環境安全，需要對該渣料場進行安全防護和水土保持治理。從安全角度來看，該渣料場除應進行整體穩定性分析外，尚應進行分區穩定性分析。受文章篇幅限制，本文選取該渣料場某一分區的典型斷面進行穩定性分析，如圖5所示。

圖5 某水電工程渣料場局部現狀圖Fig.5 Current situation of abandoned dreg site in a hydropower project

4.2 計算參數

該渣料場位于河谷地帶，地勢相對平坦，避開了沖溝、泥石流等不良地質作用區。本項目渣料塊石含量較高，經過現場踏勘、測量并結合當地經驗，取飽和狀態下渣料的重度為21.2 kN/m ，黏聚力c的均值為10 kPa，內摩擦角φ的均值為33°。此外，與算例相同，黏聚力c的空間變異系數取0.3，內摩擦角φ的空間變異系數取0.2；水平的波動范圍δx=40 m，垂直波動范圍δy=4 m。該渣料場最終的物理力學參數及計算參數取值見表2所示。

表2 渣料計算參數取值表Tab.2 Properties of slag

4.3 穩定性分析

根據現場測量資料繪制該渣料場邊坡計算斷面圖如圖6所示。該部分渣料場邊坡坡高為24 m，分為兩級，第一級邊坡高10 m，坡度i=0.78，第二級邊坡高14 m，坡度i=0.83，中間馬道寬度為4 m。考慮到原地形比較平緩，且抗剪強度參數均高于渣料參數，因此為簡化分析，計算斷面不做區分。

圖6 典型計算斷面示意圖Fig.6 Typical section of abandoned dreg site

由于渣料大部分由石渣構成，孔隙率較大，難以形成穩定的地下水，因此考慮暴雨或久雨工況時，只按照飽和參數進行分析，不考慮滲流作用。若不考慮該渣料場邊坡渣料抗剪強度參數的空間變異性，則其安全系數為1.354，大于允許的臨界安全系數1.25。考慮渣料抗剪強度參數的空間變異性時，隨機生成1 000種抗剪強度參數的空間分布，經過1 000次模擬計算，獲得1 000個安全系數及相應的滑動面。將獲得的安全系數進行平均，均值為1.367，與渣料抗剪強度參數為均值時的安全系數相近。在1 000次模擬分析中，該渣料場邊坡安全系數的最小值為1.104，最大值為1.724，但二者出現概率極低，大部分結果介于1.2～1.4之間，其中小于1.25的占23.3%，小于1.354的占61%，說明確保該邊坡穩定的概率為76.7%，均值參數計算結果的保證率僅為39%。

渣料抗剪強度為均值時的滑動面如圖7中的實線所示，考慮渣料抗剪強度空間變異性時最小的安全系數對應的滑動面如圖7中虛線所示。由圖中可以發現兩者的滑動位置基本相同，滑動深度(面積)比較接近，后者略大于前者，可能由于極端情況下，渣料抗剪強度參數均取小值，近似均勻分布，從而使其滑動位置和滑動深度均與均值參數計算結果近似，但安全系數則小得多。

圖7 均值滑面與最小安全系數滑面對比圖Fig.7 Comparison between slip surface of homogenous slope and critical slip surface

事實上，在考慮渣料抗剪強度參數空間變異性的情況下，對于相同的安全系數隨機場，也可能存在不同的滑動面。本文挑選了邊坡安全系數為1.211的3個典型的滑動面進行比較分析，3個滑動面的示意圖如圖8所示。從圖8中可以發現，3個滑面中最大的滑動體斷面面積為162.52 m2，最小的滑動體斷面面積為45.21 m2，二者相差明顯，說明考慮渣料抗剪強度參數空間變異性時，滑動面也不是唯一的、確定的。

圖8 同一安全系數對應的多種滑面形式Fig.6 Multiple critical slip surface with the same safety factor

本例計算結果說明，對于渣料場邊坡僅通過均值抗剪強度參數進行穩定性分析難以保證其安全性和準確的滑動面位置和規模等，需考慮渣料抗剪強度參數的空間變異性，進行綜合的穩定性分析和評價。

5 結 語

本文通過自相關函數構建巖土體抗剪強度參數的隨機場，利用改進的喬列斯基分解法離散抗剪強度參數隨機場，并利用非侵入式隨機有限元方法分析巖土體抗剪強度參數空間變異性對邊坡穩定性的影響。通過算例分析驗證了本文方法的合理性和正確性，最后將本文方法應用于實際的水電工程渣料場邊坡穩定性分析，分析結果可供該渣料場防護設計參考。主要結論如下：

(1)利用改進喬列斯基分解法離散隨機場，計算過程簡單，編程易于實現，離散出的隨機變量與有限單元一一對應，便于理解，是一種行之有效的模擬方法。

(2)非侵入式隨機有限元方法在考慮巖土體抗剪強度參數的空間變異性時，可利用確定性有限元軟件進行分析，無需重新編程，計算效率高。

(3)考慮巖土體抗剪強度參數的空間變異性時，邊坡的穩定性安全系數近似成正態分布，其平均值與采用均值參數的計算結果相近。統計結果表明采用均值參數計算所得的邊坡穩定性安全系數的保證率不足50%。

(4)考慮巖土體抗剪強度參數的空間變異性時，計算所得邊坡滑動面與參數隨機分布有關，不同的參數分布可能得到不同的滑動面，即使計算所得安全系數相同，對應的滑動體面積也可能相差若干倍。

(5)工程實例分析結果表明，單一的安全系數不足以評價渣料場邊坡的穩定性。對于水電工程渣料場這一特殊的人工邊坡，應考慮渣料抗剪強度參數的空間變異性，并結合工程實際進行穩定性及滑動面的綜合分析和評價。