基于失穩模量理論的降雨條件下非飽和邊坡流滑與滑移失穩機理研究

馮學茂 ，劉敏 ，宋祉辰 ，趙煉恒 ，2，3?，林宇亮

［1.中南大學 土木工程學院，湖南 長沙 410075；2.重載鐵路工程結構教育部重點實驗室（中南大學），湖南 長沙 410075；3.湖南省軌道交通工程結構防災減災重點實驗室（中南大學），湖南 長沙 410075；4.廣西新發展交通集團有限公司，廣西 南寧 530029］

降雨入滲會導致土體抗剪強度降低［1-2］，引起山體滑坡，給國民經濟和人身財產安全帶來極大的威脅.近幾十年來，氣候變化的不確定性以及人類對大自然的開采使得災難性滑坡事件愈發頻繁地發生.許多經典滑坡案例表明，在滑坡事件中通常發生有限滑移破壞和快速流動破壞（流滑）共存的情況［3-4］，特別是后者，由于其破壞力強、發生速度快、波及范圍廣且常常難以預料，因而越來越受到人們的關注［5-8］.

邊坡失穩是其內部土體發生破壞后的最終表現.該過程通常與土層內滑動面的產生有關，從力學角度看，這意味著形成一個剪切應變集中的狹窄帶［9］，而這種局部化失穩的形成過程可以憑借確切的數學條件進行描述，并可通過適當的本構模型予以預測［10］.近幾十年來，涌現出了許多理論模型用來預測和解釋上述滑坡的發生及其觸發機理［11-16］.然而，絕大多數模型無法模擬土體液化導致的流體加壓現象以及非飽和土的水力耦合效應.對于快速流動滑坡破壞模式，其內部局部土體往往易發生由固態到液態的相轉變，即發生所謂的“靜態液化”現象［17］.無論是飽和土邊坡的剪切擾動還是非飽和土邊坡的降雨入滲，根據其內部土體所經受的邊界條件（即排水條件和不排水條件）不同，邊坡的失穩模式亦有所區別：對于存在潛在液化能力的土質邊坡，其在不排水條件下觸發液化所需的載荷擾動（Δτliq，Δswc）通常要遠小于在排水條件下所需的達到臨界狀態所需的荷載擾動（Δτsf，Δssf）.以上現象說明有必要借助合適的邊坡失穩準則來判斷和區分降雨導致的土體液化和傳統的剪切破壞.

本文針對上述可能出現的不同失穩模式開展了降雨條件下非飽和土邊坡的滲流穩定性研究.從經典塑性理論出發，引入考慮基質吸力的彈塑性本構模型，利用“失穩模量”的概念建立了相應于不同失穩模式的破壞準則，并詳細介紹了理論模型所需各項參數的標定方法；為了驗證理論預測模型的正確性，引入了火山灰土的室內水槽試驗，預測了火山灰土邊坡在不同初始條件下的失穩模式以及失穩時的體積含水率，并將模型預測結果與試驗結果進行對比.此外，將失穩模量法與傳統的安全系數法相聯系，給出了同等條件下基于安全系數法的非飽和土邊坡滲流穩定性預測.最后，通過一系列參數分析，探討了影響非飽和土邊坡滲流失穩模式的相關因素，以期為進一步理解非飽和土邊坡在降雨條件下失穩后演變成的不同運動模式提供理論依據.值得一提的是，由于非飽和土的降雨入滲過程并不是本文研究的重點，為了簡便起見，主要借非飽和土含水率的變化來反映降雨入滲的影響.

1 彈塑性材料的失穩準則

1.1 二階功理論

在小應變連續力學的背景下，Hill［18］提出了二階功的概念并由此建立了傳統的穩定性準則：

該理論不僅可以較好地描述材料不穩定性狀態的擴散效應，還能同應力應變增量響應解的存在性和唯一性聯系起來［19-20］.在任意給定的應力-應變路徑增量下，式（1）大于零，即d2W＞ 0，意味著材料具有足夠的穩定性，同時應力應變增量解具有唯一性.而當d2W＜ 0 時表明材料處于不穩定狀態.因此，二階功的喪失（即d2W=0）意味著材料將無法承受所施加的擾動，而在數學模型上則體現為本構控制矩陣為奇異矩陣.

1.2 失穩模量

考慮一般形式下的彈塑性材料經受塑性加載，根據經典塑性理論的一致性條件，有：

式中：f（σ'，χ）為屈服面函數；df為屈服面增量；σ'為廣義有效應力向量；χ則代表與硬化相關的變量.偏導數和分別量化了屈服面隨應力和硬化參數的變化.

根據彈塑性理論，彈塑性材料的應變可以分為彈性（可恢復）部分和塑性（不可恢復）部分，即：

式中：δε為總應變增量；δεe和δεp分別為彈性應變增量和塑性應變增量.而當有效應力發生變化時總會相應地產生彈性變形，因此有：

式中：De為彈性剛度矩陣.

結合式（2）、式（3）、式（4）以及非關聯流動法則δεp=對于塑性硬化材料在標準應變控制加載條件下的情形，即χ=f（εp），有：

式中：Λ（＞0）為塑性乘子，量化塑性應變隨荷載增量的變化幅度而變化；g為塑性勢函數.其中，硬化模量H可表示為：

于是，式（5）可改寫為：

將上式重新調整并合并成含有Λ的項，即可得到：

至此即得到反映總應力與總應變關系的顯示公式.從式（8）中不難發現，當HIN→0 時，Λ→∞，此時本構矩陣為奇異矩陣，從試驗可控性的角度意味著達到材料的失穩狀態，因此也將HIN稱為“失穩模量”.

1.3 本構方程

為了推導出失穩模量HIN的顯式表達式，現考慮基于廣義塑性力學的Cam-Clay本構模型［21］：

式中：τ為剪應力；Mg為達到極限剪切破壞時的摩擦系數（Mg=tanφ'）；Mf為飽和條件下的抗液化強度（Mf=tanφLIQ，φLIQ為液化摩擦角）；為當前屈服面與有效平均主應力軸的交點，即前期固結應力（壓力）；為使塑性勢適應當前應力狀態的虛擬變量；σ*為平均土骨架應力，由下式給出：

式中：σnet為法向凈應力；s為基質吸力；k為吸力強度系數.

屈服面的硬化法則由其增量形式dσ*Y給出：

式中：rw為控制屈服面在受到浸濕后發生收縮的參數；λ為壓縮系數；ξ為破壞時的剪脹系數.

除了力學行為的描述外，該模型還需要反映非飽和土的水文變化規律，即土-水特征模型，用來描述土體基質吸力、滲透率與飽和度的關系.本文考慮Gardner水文模型［22］，有：

式中：θres、θsat分別為殘余體積含水率和飽和體積含水率；αw為反映進氣值的模型參數.

1.4 失穩模式

在降雨條件下的非飽和土邊坡中，外界的擾動通常較為復雜，一旦形成土體的不排水條件，則可能在特定的條件下引起土體的液化，進而造成邊坡發生流滑［6，23］.因此，不同的水力邊界條件對于邊坡的破壞模式有一定影響.現考慮兩種極端情況下的失穩模式：1）模式A，基質吸力控制下的滑移失穩，即排水剪切下的土體失穩；2）模式B，含水率控制下的流滑失穩，即不排水剪切失穩.如圖1所示.

圖1 無限邊坡失穩區周圍邊界條件示意圖Fig.1 Boundary conditions around the unstable zone of an infinite slope

在模式A 中，排水加載假設滑動層邊界具有無限大的透水性，因此基質吸力不再受局部應力應變的影響，而僅取決于當前的水力邊界條件.在模式B中，假設滑動層邊界的滲透性為零，導致剪切過程中產生的超孔隙水壓無法消散，基質吸力的變化來自多孔介質的飽和度和變形耦合的共同作用.特別地，由于在不排水失穩（模式B）中，很難直接測得基質吸力的變化率.為此，Buscarnera［24］引入了比水量νw的概念來反映基質吸力的影響：

式中：Fr為土水特征曲線的反函數，即s=Fr（Sr）.與失穩模式A不同的是，基質吸力不再僅僅受含水率的變化（dvw）的影響，同時也取決于土體的塑性變形（dεp）.

由失穩模量的定義，通過鏈式法則求導可得上述兩種失穩模式下的以及的顯式表達式，分別為：

其中，應力比η*由下式定義：

2 理論模型驗證——以水槽試驗為例

Greco等［25］針對松散火山灰土（臨界孔隙度為ncr=0.70）進行了一系列的室內水槽試驗，指出土在非飽和狀態下的體積塌陷是引起流滑的主要因素之一.本文以上述水槽試驗為參考，將模型預測結果與試驗結果進行對比來驗證理論模型的正確性.

2.1 模型參數標定

參考前人針對當地該類型土的室內試驗結果對模型所需的各項參數值予以標定.吸力強度系數k量化了基質吸力對排水抗剪強度的影響，可由直剪試驗［26］結果標定得到，如圖2（a）所示，通過對數據點進行擬合得到k=0.61；屈服面函數f和塑性勢函數g由式（10）給出，其中參數Mf與Mg分別與各自函數的頂點位置有關.通過對排水與不排水三軸試驗［7］數據擬合可以標定Mg=0.75（即φ'=37°），Mf=0.47（即φLIQ=25°）；硬化參數λ反映屈服面隨飽和狀態下的塑性應變的變化情況，而rw控制屈服面隨吸力值變化，由側限壓縮試驗［27］可得λ=0.12，rw=0.04；根據室內水槽試驗［25］的土體監測數據，用Gardner 水文模型擬合土-水特征曲線可得模型參數α=0.6.標定后的模型參數匯總見表1.

表1 理論預估模型參數標定匯總Tab.1 Summary of theoretical estimation model parameter calibration

圖2 模型參數試驗標定Fig.2 Model parameter test calibration

2.2 計算對比結果

表2 為Greco 等［25］對于不同初始狀態的火山灰土邊坡水槽試驗的結果.由表2 可知，對于不同初始孔隙度的樣本出現了不同的失穩模式.圖3（a）為試驗編號為D3 的樣本在體積含水率不斷增加的情況下其失穩模量的變化.該圖表明當考慮模式A，即φ'基質吸力保持恒定的排水條件下時，模型預測當體積含水率θf=0.694 時試驗邊坡破壞；而當考慮模式B，即當考慮土體的不排水失穩時，模型預測當θf=0.673 時邊坡破壞，更接近真實條件下土坡失穩時的體積含水率（θf=0.67），因此模型預測該破壞模式為流滑形式.圖3（b）為試驗中觀測到的D3累計濕潤鋒（h）以及平均坡表沉降（w）.由圖可知，當t=20 min時，壓力水頭發生陡增，與此同時坡表沉降增加，即內部土體發生體積塌陷，有靜態液化現象發生，最終演變為流滑失穩.

表2 水槽試驗初始條件及試驗結果Tab.2 Initial conditions and test results of the flume test

圖3 D3號樣本失穩模型預測與累計濕潤鋒及坡表沉降Fig.3 Model predictions of the failure mode，cumulated infiltration height h and mean soil surface settlement w measured during the D3 infiltration experiment

為了進一步說明模型的預測效率，圖4 為利用失穩模量法預測其余各組樣本出現的破壞形式以及相應失穩時的體積含水率.由計算結果可知，除D9-1 樣本外，其余各組均發生失穩.其中D4、D8 失穩后演變為流滑破壞，而D7 的后破壞模式為摩擦滑移破壞，與試驗觀測結果相符.而就失穩時的體積含水率來說，模型預測結果與試驗結果較為接近，其相對誤差小于1%，如表3所示.

表3 水槽試驗結果與模型預測結果對比統計Tab.3 Comparison and statistics of flume test results and model prediction results

圖4 基于失穩模量的模型預測失穩模式及其相應失穩時的體積含水率與實測數據對比Fig.4 Prediction of instability mode based on instability modulus model，comparison of volume moisture content at the corresponding instability with measured data

2.3 安全系數

失穩模量法從機理上解釋了非飽和土邊坡在不同加載條件下可能出現的不同失穩形式，而在巖土工程領域中，常通過安全系數來評估土體或邊坡的穩定性.對于無黏性土的無限邊坡模型而言，其安全系數可定義為當前應力比與臨界應力比的比值［28-29］，即：

式中：β為邊坡坡度為臨界應力比；rw為吸力硬化參數；λ為塑性壓縮性；Sr為有效飽和度；在式（15）中，當HIN=0 時，所對應的應力比即為臨界應力比ηf；G=G（Sr）是相對于St的水分特征曲線（WRC）的導數.式（18）、式（19）分別代表失穩模式A 和失穩模式B下的臨界安全系數.為了說明理論方法的正確性與適用性，下面用安全系數法針對以上算例給出相應的失穩預測如圖5 所示.由圖5 可知，使用安全系數法的計算結果與失穩模量法一致.特別地，由于D7試樣的初始孔隙度n0小于臨界孔隙度ncr，從而導致模式B 下的安全系數計算出現虛根，因此僅會發生模式A 情形下的失穩模式（滑移失穩）.由以上可知，無論是從安全系數的角度還是失穩模量的角度分析，都可以較好地預測非飽和土邊坡的滲流穩定性，二者能夠互為驗證.

圖5 利用FS預測破壞模式及其相應失穩時的體積含水率與實測數據對比Fig.5 Using FS to predict the failure mode and its corresponding volume water content during instability compared with measured data

3 參數分析

為了進一步探究在不同排水條件下（恒定含水率或恒定基質吸力），模型各參數對非飽和土滲流穩定性以及相應破壞模式的影響，選取上節分析中水槽試驗D4 樣本并以滲流過程中失穩模量HIN的演變趨勢進行說明.

3.1 初始孔隙度n的影響

上節研究表明，相對松散土更容易發生靜態液化，進而引發邊坡流滑失穩.因此，土體的初始孔隙度對邊坡的失穩模式有較大影響.圖6 為不同初始孔隙度條件下D4 樣本在相同降雨條件作用下的滲流穩定性計算結果.由圖6 可知，當土的初始孔隙度小于臨界孔隙度（ncr=0.70）時，較難發生流動破壞［圖6（a）］，土體以摩擦滑移破壞模式為主［圖6（b）］.這是由于相對密實土在加濕過程中發生體積膨脹，產生負的孔隙水壓力，不會發生靜態液化.相反，當土體初始孔隙度大于臨界孔隙度時，在加濕過程中發生體積收縮，產生正的孔隙水壓力，因此易發生靜態液化.由預測結果可知觸發流滑破壞時的體積含水率要小于觸發滑移破壞時的體積含水率，意味著松散土以流滑破壞形式為主.

圖6 相同條件下初始孔隙度對土體滲流穩定性的影響Fig.6 The influence of the initial porosity on infiltration stability of soils under the same conditions

3.2 壓縮系數λ的影響

圖7為不同壓縮系數條件下D4樣本在相同降雨條件作用下的滲流穩定性計算結果.由于λ反映了土體的硬化規律，該值越小時屈服面隨塑性應變的變化越快.因此當壓縮系數越小時，失穩模量對體積含水率的敏感性越高.就破壞時的體積含水率而言，當邊界條件為不排水時（恒定含水率條件），λ越大，土體體積收縮性越強，越容易在正的孔隙水壓作用下發生靜態液化從而發生流滑，其滲流穩定性越低［圖7（a）］.當邊界條件為排水時（恒定基質吸力條件），由于不考慮非飽和土的變形耦合效應，λ并不影響土體的滲流穩定性［圖7（b）］.

圖7 相同條件下壓縮系數對土體滲流穩定性的影響Fig.7 The influence of the compressibility coefficient on infiltration stability of soils under the same conditions

3.3 吸力硬化系數rw的影響

圖8為不同吸力硬化系數條件下D4樣本在相同降雨條件作用下的滲流穩定性計算結果.由于rw越大表明非飽和土屈服面隨基質吸力的變化越大，具體表現為圖中失穩模量隨體積含水率變化越敏感，因此滲流穩定性越低.由于在排水加載條件下并未考慮該參數，因此rw對該條件下的不穩定模量/安全系數的計算沒有影響.

圖8 相同條件下吸力硬化系數對土體滲流穩定性的影響Fig.8 The influence of the suction hardening coefficient on infiltration stability of soils under the same conditions

3.4 吸力強度系數k的影響

圖9為不同吸力強度系數條件下D4樣本在相同降雨條件作用下的滲流穩定性計算結果.吸力強度系數k直接影響土骨架應力σ*，因此無論是恒定含水率條件還是恒定基質吸力條件，參數k的提高對土體的滲流穩定性都有著積極作用.也就是說k的增大使得土體越難發生流動破壞或滑移破壞.

圖9 相同條件下吸力強度系數對土體滲流穩定性的影響Fig.9 The influence of the suction strength coefficient on infiltration stability of soils under the same conditions

3.5 孔徑分布系數α的影響

上述參數分析主要集中表現力學參數的影響，而水文模型參數同樣影響著非飽和土的滲流穩定性.圖10 即為Gardner 模型中不同α條件下，D4 樣本在相同降雨條件作用下的滲流穩定性計算結果.由于α與土體的粒徑分布有關，是與進氣壓力值有關的參數，因此α越大，相同體積函數率條件下基質吸力越低，非飽和土的滲流穩定性越低.

圖10 相同條件下孔徑分布系數對土體滲流穩定性的影響Fig.10 The influence of the aperture distribution coefficient on infiltration stability of soils under the same conditions

4 結論

1）基于失穩模量理論建立的非飽和土邊坡滲流穩定性準則能較好地預測非飽和土邊坡在排水條件下所發生的摩擦滑移失穩模式以及在不排水條件下發生的快速流動破壞模式（流滑）.

2）初始孔隙度越大通常意味著需要更大的體積含水率才能引起土體的失穩，而相對松散土在不排水條件下易發生土體液化，進而形成邊坡的流滑失穩；而相對密實土由于體積膨脹產生負的孔隙水壓力，因此不易發生液化，其破壞模式以有限滑移為主.

3）壓縮系數λ和吸力硬化系數rw與非飽和土的硬化特性有關.對于不排水加載條件下的流動破壞，λ（或rw）越大，土體的失穩含水率越低，即越容易發生滲流失穩.而當考慮基質吸力控制下的滑移破壞模式時，二者對土體的滲流穩定性幾乎沒有影響.

4）吸力強度參數k對非飽和土的滲流穩定性起著積極作用，該值越大，則觸發土體失穩的體積含水率越大，即滲流穩定性越高.

5）Gardner水文模型中的孔徑分布參數α對非飽和土的滲流穩定性起著消極作用，即α越大，非飽和土的滲流穩定性越低.