何志亮

（中鐵水利水電規(guī)劃設(shè)計集團有限公司，江西南昌，330029）

由于暴雨或地震影響，山區(qū)通常會發(fā)生臨時或永久性河流阻塞。由這一過程形成的滑坡堆積壩可能會因漫頂、壩體突然坍塌或漸進性破壞引起的侵蝕而失效。大壩的滑坡是一個非穩(wěn)定滲流過程，邊坡穩(wěn)定分析與非穩(wěn)定滲流分析相結(jié)合有助于滑坡堆積壩的穩(wěn)定性研究。

許多學(xué)者采用不同的邊坡穩(wěn)定性分析方法，重點是結(jié)合水文模型與滲流穩(wěn)定分析。肖耀廷等人[1]采用了VIC模型，對滑坡體穩(wěn)定性進行了有限元評估與分析；鮮丹等人[2]以甘肅某村的滑坡體為研究對象，研究了初始含水率對邊坡穩(wěn)定性的影響原理；辛長虹等人[3]基于非飽和滲流分析理論，結(jié)合有限元數(shù)值方法，分析了大壩位移和應(yīng)力規(guī)律；付建新等[4]基于交替隱式有限差分法，將不同含水率作為初始條件進行二維入滲數(shù)值計算與分析；陳璽等人[5]對潰壩時的壩坡穩(wěn)定進行分析；周建云[6]分析了降雨入滲對水庫壩坡失穩(wěn)的影響，結(jié)合2D模型對水庫潰壩進行演算；孫蔚等人[7]考慮了滲流和巖土體的耦合作用，進行了滑坡堆積壩的穩(wěn)定性分析。但這些研究大多局限于二維分析，僅適用于三維長邊坡。然而，天然斜坡和滑坡堆積壩在狹窄的U形或V形山谷中的破壞是三維的，因此三維方法更適合于分析此類穩(wěn)定性問題。

1 數(shù)值模型

采用極限平衡法對邊坡穩(wěn)定性進行了評價。根據(jù)壩體內(nèi)部水的入滲情況，計算安全系數(shù)，尋找安全系數(shù)最低的臨界滑動面。

由兩個模型組成，瞬態(tài)滲流模型計算壩內(nèi)孔隙水壓力和含水量的變化。圖1為邊坡瞬態(tài)穩(wěn)定分析耦合模型通用流程。由于上游水庫水位逐漸升高，邊坡穩(wěn)定模型根據(jù)壩體孔隙水壓力和水分運動的變化計算邊坡的安全系數(shù)和臨界滑動面的幾何形狀。上游水庫水位（WL）在每個時間步由水量平衡方程確定，考慮流入流量（Qin）、壩體內(nèi)入滲率（Qsh）和水槽幾何形狀（寬度和坡度），耦合模型的總體輪廓如圖1所示。

圖1 邊坡暫態(tài)穩(wěn)定分析耦合模型通用流程Fig.1 General flow chart of coupling model for transient stabili-ty analysis of slope

1.1 滲流模型

在狹窄的山谷中形成的滑坡堆積壩滲流可能產(chǎn)生三維效應(yīng)。為考慮非飽和滑坡堆積壩三維滲流問題，建立三維滲流模型。為了評價飽和土的孔隙水壓力變化，采用了基于壓力的修正方程：

式中：h為水壓頭；K x(h)、K y(h)、K z(h)為x、y和z方向的導(dǎo)水率；C(h)=?θ/?h，為比濕容量，θ為土壤體積含水量；t為時間；S w為飽和比；S s為比庫容；x、y為水平空間坐標(biāo)；z為垂直空間坐標(biāo)，取向上為正。公式（1）表示非飽和域和飽和域中的流動。對于飽和域，K x(h)=K y(h)=K z(h)=K s，θ=θs，C(h)=0，S w=1，其中K s和θs分別為飽和導(dǎo)水率和飽和含水量。S s取決于固體基質(zhì)和流體的可壓縮性，因此在非飽和無側(cè)壓多孔介質(zhì)中，S s接近于0。

為求解式（1）中的方程，結(jié)合式（2）和式（3）：

式中：α和n v為方程參數(shù)；S e為有效飽和度；θs和θr分別為飽和含水量和殘余含水量。

1.2 邊坡穩(wěn)定性模型

用極限平衡分析法評價滑坡堆積壩的瞬態(tài)邊坡穩(wěn)定性，包括計算安全系數(shù)和尋找安全系數(shù)最低的臨界滑動面?；趯θ我饣浦臉O限平衡分析，本研究中考慮的邊坡是平緩的（小于45°）。一般情況下，安全系數(shù)只能由一個公式而不是聯(lián)立方程確定，因此迭代程序被簡化。邊坡穩(wěn)定性模型結(jié)合了基于動態(tài)規(guī)劃的最小化方法和隨機數(shù)生成方法，將該模型與滲流模型耦合，進行邊坡穩(wěn)定分析。

山中景色隨四時而變。出自北宋著名的畫家兼山水畫理論家郭熙寫的《山水訓(xùn)》：“真山之煙嵐，四時不同。春山淡冶而如笑，夏山蒼翠而欲滴，秋山明凈而如妝，冬山慘淡而如睡?！闭龑懗隽松剿S季節(jié)變化而迥然不同的物象特征，但這特征正如描寫的美人的姿態(tài)一樣“淡冶而如笑，蒼翠而欲滴，明凈而如妝，慘淡而如睡”，以此來比喻山水的春夏秋冬。

圖2顯示了斜坡內(nèi)的三維滑動體和垂直分割柱。圖3顯示了作用在圖2中典型柱上的力：Wij為柱的重量；Pij為作用在柱頂?shù)呢Q向外力；T ij和Nij為作用在柱腳上的剪力和總法向力；Q ij為作用在柱側(cè)面的所有柱間力的合力。將Q ij分解為平行于x軸的Q1和以一定角度β=tan-1(ηtanαxzij)傾斜的分量Q2，其中η為常數(shù)。此外，圖3中的αxzij和αyzij為柱基在xz和yz平面上與水平方向的傾角。

圖2 三維斜坡的階段狀態(tài)系統(tǒng)和劃分方案Fig.2 System of stage state and classification scheme of 3D slope

圖3 作用在典型柱上的力Fig.3 The force acting on the typical column

在整個滑動體的水平和垂直方向上應(yīng)用平衡條件，引入兩個安全系數(shù)Fsh和Fsv，這兩個系數(shù)都取決于η。對于較小的η值，通常滿足方程F=Fsh=Fsv，計算出的Fsh值受η變化的影響不大?？梢杂搔?0得出Fsh并基于此近似值計算3D安全系數(shù)，安全系數(shù)F可通過式（4）表示：

式中：c ij和φ為土的有效強度參數(shù)；u ij為柱基孔 隙 水 壓 力 ；和為作用于柱上下游表面的水壓力合力。

莫爾-庫侖方程可用于計算非飽和土的抗剪強度，從而得出公式（4）。然而，在這項研究中，由于負孔隙水壓力，剪切強度的增加被忽略了。

在公式（4）中，W ij、P ij和uij（無黏性土的cij=0）是時間相關(guān)的，可以在每個時間步長通過與3D瞬態(tài)滲流模型耦合確定。滲流模型的輸出（每個單元的孔隙水壓力和含水量）和邊界條件（地表水的深度）可用于確定這些物理量。

2 試驗研究

測量水分分布和觀察破壞面的試驗總結(jié)如表1所示。使用平均直徑為1 mm的混合硅砂在水槽中制成三角壩。方程（2）和（3）中的參數(shù)（包括θr）通過對pF儀試驗獲得的土壤持水量數(shù)據(jù)進行非線性回歸分析進行估算。

表1 試驗總結(jié)Table 1 Test results

利用含水率反射計測量了滲流過程中含水率隨時間的變化。為了測量壩坡在滑動過程中的運動，在槽壁面壩體內(nèi)放置紅色的泥沙條帶。在渡槽側(cè)邊放置數(shù)碼攝像機，拍攝渡槽突然滑動時的滑動面形狀。通過對水槽兩側(cè)的視頻分析，得出了壩體滑動過程中的滑動面形狀。

使用長500 cm、寬30 cm、高50 cm的水槽。設(shè)定水槽的坡度為20°。壩體形狀及尺寸見圖4和圖5。含水率反射計的排列如圖5所示。在水槽的B側(cè)開孔，并從該側(cè)將12個儀器插入壩體內(nèi)。進行了兩組試驗來測量壩體不同位置的水分剖面。

圖4 壩體形狀及壩頂斷面Fig.4 Dam shape and dam crest section

圖5 含水率反射計（1-12）布置圖和B側(cè)視圖Fig.5 Layout diagram of moisture content reflectometry(1-12)and view of B side

3 結(jié)果分析

庫水位逐漸上升導(dǎo)致水滲入壩體，增加了孔隙水壓力和壩體重量。當(dāng)壩體自重增大引起的剪應(yīng)力大于孔隙水壓力增大引起的抗剪應(yīng)力時，壩體就會發(fā)生滑動。

穩(wěn)定的水流量來自水槽上游，測量出壩體中水分含量。測得的水分含量是距離水槽側(cè)面18 cm處的平均值（見圖5）。由于水在壩體內(nèi)運動的三維效應(yīng)，水分運動不是恒定的。盡管兩個試驗中使用的排放量幾乎相等（見圖6），但3D-1和3D-2測量的水分分布略有不同，這可能由兩個試驗中飽和導(dǎo)水率的變化所致。圖7顯示了不同水分剖面的模擬和試驗結(jié)果比較，具有良好的一致性。

圖6 水分剖面比較（示例3D-1和3D-2）Fig.6 Comparison of water profiles(experiment 3D-1 and 3D-2)

圖7 水分分布比較（模擬和試驗）Fig.7 Comparison of water distribution(simulation and experi-ment)

對比兩個試驗3D-3和3D-4以觀察邊坡的破壞情況。試驗3D-3（圖8）中，在930 s處觀察到斜坡失效。而在試驗3D-4的情況下，在1 030 s處觀察到邊坡的破壞，滑動面也比3D-3更深（圖9），這可能由兩個試驗在壓實和飽和導(dǎo)水率方面的不均勻性不同而導(dǎo)致。

圖8 930 s時的滑動面（試驗：3D-3）Fig.8 Sliding surface at 930 seconds(experiment:3D-3)

圖9 1 030 s時的滑動面（試驗：3D-4）Fig.9 Sliding surface at 1 030 seconds(experiment:3D-4)

在三維邊坡穩(wěn)定性初步分析的基礎(chǔ)上，在每個階段平面上生成了上千個狀態(tài)。滲流模型中用于模擬的其他水力條件/參數(shù)和網(wǎng)格系統(tǒng)為Qin=29.8 cm3/s，K s=0.000 3 m/s，θs=0.287，θr=0.045，α=5.5 m-1，n v=3.2，Δt=0.01 s，塊體尺寸為10 mm。邊坡穩(wěn)定模型采用的單元大小Δx=5 cm、Δy=3 cm。

在t=0 s和t=770 s時，模擬的臨界滑移面如圖10所示。模擬安全系數(shù)（F）一開始時大于1。770 s時，模擬的安全系數(shù)小于1，但試驗中觀察到的破壞時間約為930 s。

模擬過程中，考慮蒸發(fā)作用是以減少29 cm3/s的流量實現(xiàn)的，考慮水力傳導(dǎo)率的不確定性是通過減少0.000 28 m/s的飽和傳導(dǎo)率來進行的。這兩種情況下，大壩都在790 s時失效。Ks=0.000 25 m/s時，模擬失效時間為830 s。因此，失效時間也取決于飽和導(dǎo)水率。初步分析（圖11）表明，使用此3D方法計算的安全系數(shù)，當(dāng)滑坡堆積壩破壞時為1.05~1.06。因此，如果使用更嚴(yán)格的邊坡穩(wěn)定性分析方法，破壞時間將延長并接近觀察到的時間。此外，在計算中也忽略了水槽側(cè)壁的摩擦力。水槽兩面破壞面（圖10）與試驗吻合較好。

圖10 模擬臨界滑動面（3D）Fig.10 Simulation of criticalsliding surface(3D)

圖11 用3D方法模擬不同Ks和安全系數(shù)下的破壞時間Fig.11 Simulation of failure time by 3D method with different Ks and safety factors

滑動面的形狀取決于許多因素，如土體物理力學(xué)性質(zhì)、植被分布、孔隙水壓力等，可進一步完善模型，采用更嚴(yán)格的邊坡穩(wěn)定分析方法，是三維邊坡穩(wěn)定分析實際應(yīng)用的關(guān)鍵。

4 結(jié)語

通過水槽試驗和數(shù)值模擬，對滑坡堆積壩的三維瞬態(tài)邊坡穩(wěn)定性進行了分析，得出以下結(jié)論：

（1）在上游水庫穩(wěn)定泄流試驗水槽中，對滑坡堆積壩的突然破壞進行了研究。庫水位逐漸上升使水滲入壩體，增大了動剪應(yīng)力，當(dāng)動剪應(yīng)力大于抗剪應(yīng)力時，壩體突然垮塌。

（2）采用極限平衡法建立了邊坡穩(wěn)定與非穩(wěn)定滲流耦合模型。通過數(shù)值模擬和水槽試驗，研究了滑坡潰壩機理。結(jié)果表明，數(shù)值模擬和試驗結(jié)果在壩體水分運動和預(yù)測臨界滑動面方面非常接近。

（3）然而，在模擬中，壩體發(fā)生破壞的時間較早。今后應(yīng)通過引入更嚴(yán)格的邊坡穩(wěn)定性分析方法對模型進行進一步的改進，使其既可用于滑坡堆積壩的邊坡穩(wěn)定性分析，也可用于天然邊坡的穩(wěn)定性分析。