某復雜古滑坡堆積體在蓄水條件下的穩定性評價
2016-12-14 22:35:39闞露張發明
安徽理工大學學報·自然科學版 2016年1期
闞露+張發明
摘要：文章以位于某水電站壩址區右岸古滑坡Ⅱ區堆積體為研究對象，建立了一個較為復雜的三維模型，并從主應力、塑性區、剪應力和位移四個方面對古滑坡堆積體在蓄水條件下的變形和穩定性情況進行了詳細的模擬和分析。計算得到的四個方面的結果都相互吻合相互印證。文章以FLAC3D中的單元界面能自動依附于指定范圍內模型表面生成的特性為思路來模擬地下水位。研究結果表明，堆積體雖然未出現整體性滑動破壞，但要注意坡腳局部地區坍塌和土層內部局部滑弧滑動，具體直觀地在模型上指出堆積體上的不穩定區域，為后期針對性的加固治理方案提供參考，為工程優化設計提供借鑒意義。
關鍵詞：堆積體；三維數值模擬； 穩定性
中圖分類號：文獻標志碼：A文章編號：1672-1098（2013）03-0000-00
Abstract：This paper studies the deformation and stability of the Ⅱregion of the ancient landslide body on the right bank of a hydropower station by three-dimensional numerical simulation. Detailed simulation and analysis about deformation and stability is realized by the calculation of main stress， plastic zone， shear stress and displacement under water-reserving condition. Results from the four respects confirm each other. The unit interfaces in FLAC3D can attach on the model surfaces at specified areas automatically， which becomes the method to simulate the groundwater level in this article. The finding includes the following. Although there is no whole slip in the calculation， collapse on local areas of the toe of the slope and local circular slide inside the soil layer should be paid attention to. The unstable regions on the slope are pointed out on the slope clearly. References can be presented for the reinforcement method and the optimization design of the project.
Key words：landslide；numerical simulation； stability
崩塌、滑坡發生后形成的堆積體變形及穩定性問題一直是工程地質領域研究的重點問題之一，滑坡形成的堆積體結構松散，強度低，在受到擾動后容易失穩，對水電站的水庫正常運營會造成很大的影響，尤其是大型堆積體如果發生破壞，很可能會堵塞河道，后果不堪設想。迄今為止，研究者們已提出了如極限平衡分析法、塊體理論法[1]、數值方法[2-3]等計算條件相對簡單的傳統方法，也提出了人工神經網絡理論[4]、模糊數學理論[5]、灰色系統理論[6]、可靠性理論[7]等基于未確知性理論的計算方法。
本文采用FLAC3D數值模擬的方法對某水電站壩址區右岸古滑坡Ⅱ區堆積體進行穩定性分析， 由于FLAC3D采用的是顯式差分法求解微分方程，不形成剛度矩陣，每一步計算僅需要很小的內存。在求解過程中通過疊加每一時步的小變形獲得大變形求解，從而僅占用計算機很少的內存就可以模擬大變形，并且在模擬材料的塑性破壞和塑性流動方面具有優勢。但是由于FLAC前處理功能較弱，建立復雜的三維模型非常困難，故本文選用ANSYS進行堆積體前期建模，并劃分計算網格，再導入FLAC3D程序中計算，充分利用了兩種軟件的優點，優化分析過程。
2計算原理
21有限差分原理
有限差分法就是在利用數值計算方法求解偏微分方程時，用有限差分近似公式代替每一處導數，從而將求解偏微分方程的問題轉化成求解代數方程的問題。
FLAC有限差分公式：由高斯散度定理有
∫snids=∫Afxi （1）
式中：∫s為封閉曲面上沿邊界的積分； ni為曲面s的單位法向量；f為標量、向量或張量； xi為位置向量； ds為增量弧長；∫A為對表面積A積分。
定義f在面A上的梯度平均值為：
式中：<>表示求平均值。
將（2）代入（1）后得到
對一個三角形子單元
式中：s為三角形某一條邊的邊長，等式右邊的求和在三角形三條邊上加總。
應力應變：用每一邊速度矢量均值ui代替（4）式中的f，ui取各條邊兩端點的結點（即差分網格的角點）a和b的速度平均值，則
同理可以求出
eij=12uixj+ujxi （6）
由材料的本構方程和相應的邊界條件，就可以求得應力增量。對各向同性的材料，有
σij=λδijθ+3μeij （7）
式中：λ，μ為拉梅常數；θ為體積應變，當時i=j，δij為1，否則，δij為零。
這樣通過上述各式的迭代求解，就可以得出每一迭代時步對應的單元的應力應變值[9]。
22強度折減法
本文通過強度折減法來計算堆積體的穩定系數。強度折減法的原理是通過對巖土體的強度指標C和Φ值不斷地折減，反復計算，直到塑性區貫通整個坡體，即邊坡達到了臨界的破壞狀態，此時的折減系數就作為坡體的穩定系數Fs。
Cf=CFtrial （8）
Φf=tan-1（tanΦFtrial） （9）
式中：Cf為折減后的粘結力；Φf為折減后的摩擦角；Ftrial為折減系數。
3堆積體基本特征
古滑坡堆積體位于貴州某水電站選壩河段，分為四個區域，其中Ⅱ區巖土體結構最為復雜，控制底界以上物質主要以崩坡積塊碎石土組成，故以Ⅱ區作為本文研究對象。堆積體所在河段河谷為不對稱的“V”型河谷，總體為順向坡，呈明顯圈椅狀地貌，順河長600～700 m，橫河寬約730 m，前緣臨河，后緣高程為740～790 m。順河方向上呈上升臺階狀，上游側地勢低于下游側40～70 m。總體上地形坡度約25°。后緣緊靠灰巖陡壁，區域內發育兩條淺沖溝。覆蓋層物質組成為崩坡積碎塊石土為主。下伏基巖主要為J1夾層以下的二疊系下統棲霞組（P1q）中厚層夾薄層含炭質灰巖。圖1Ⅱ區地質分區平面示意圖及13條剖面線位置
J1夾層為該區早期層狀巖體失穩破壞的控制底界，夾層以上層狀座滑體已經完全剝蝕搬運走，目前該區上部覆蓋物質來源應為后緣陡壁后期崩塌堆積形成。平面示意圖和地質橫剖面示意圖分別見圖1和圖2。
4Flac3D數值模擬
41幾何模型
模型邊坡的主滑方向為S14°E，走向N76°E。模型以主滑方向的反方向為X軸正方向，豎直向上為Y軸正方向，走向NE方向為Z軸正方向。豎直方向坐標采用實際高程坐標，建立直角坐標系，如圖3 所示。根據現場提供的地質資料，Ⅱ區堆積體模型選取了13個剖面，分別是Z=0，Z=57，Z=128，Z=215，Z=272，Z=325，Z=380，Z=438，Z=501，Z=561，Z=627，Z=697，Z=786（剖面線位置見圖1），左邊界略包括了Ⅰ區部分巖體。模型邊坡底面高程190 m，坡頂高程850 m，沿走向方向長度為786 m，沿傾向方向長度為1400 m。模型在X軸上的范圍是-200～1200 m，在Y軸上的范圍是190～850 m，在Z軸上的范圍是0～786 m。
模型分為5個部分，如圖3，對應圖上顏色，1，5部分為基巖，2為j1夾層，3為Ⅰ區座滑層狀巖體，4為Ⅱ區表面堆積體，其中1，3，4，5部分選用solid45單元類型，夾層選用shell99單元類型。Ⅱ區堆積體模型共有15569個單元，4969個節點。計算模型邊界約束形式為：側邊界只對水平方向進行約束，底邊界在水平和豎直方向都進行約束，模型的上部邊界取為自由面。FLAC程序中，節點速度是主要變量，所以選取模型的邊界條件是通過約束模型邊界的節點速度實現的，即模型底部邊界的水平、豎直方向的速度約束和四周邊界水平方向的速度約束。在程序中表達為邊界上xvel，yvel設置0. 根據實際地勘報告，計算采用的各層物理力學
地下水面的生成：由于在FLAC3D中直接生成符合勘察資料的地下水面比較困難，所以本文中建立模型時充分利用FLAC3D中的單元界面能自動依附于指定范圍內模型表面生成的特性，生成水面。正常河水位為383 m，水庫正常蓄水位為400 m。根據勘察報告，本文采用地下水位從坡腳處高程400 m開始，以5°向上傾斜，來近似模擬實際地下水位形態，如圖4。沿著水位線，對三維模型采用切割命令，得到地下水面及其下部的幾何實體，并剖分網格，導入FLAC3D，并以之為輔助單元，生成與地下水面空間形態一致的單元界面。采用FISH語言[10]遍歷界面單元節點，生成水面，同時生成靜水壓力。
42計算結果分析
421主應力場規律分析計算表明，坡體的最大主應力及最小主應力基本為層狀分布，并在坡腳高程約400～460 m處出現應力集中，坡體大主應力σ1最大值為18139 MPa，小主應力值σ3最大值為10585 MPa。具體應力分布見圖5及圖6。
從剖面Z=650 m的應力分布圖（如圖7）可以看出，剖面附近的最大主應力（壓應力），基本順著坡面方向，并一直延伸到坡腳。而往邊坡內部，最大主應力方向與水平軸的夾角逐步變大，直至鉛直；由于巖層分界面的存在，使得其附近區域的最大主應力方向要比其他區域最大主應力方向變化大而且迅速，但并未影響主應力分布的總體走勢。這些都表明邊坡深部土體主要受鉛垂方向的壓應力作用，體現為受壓屈服。
來看，坡體的屈服區域較集中，范圍較小。坡體后緣（高程約780～850 m）分布有小塊拉張塑性區，坡腳處（高程約380～450 m）分布零散的剪切塑性區。堆積體發生張拉剪切破壞的可能性很小，即使發生，也僅是局部區域，不會對堆積體整體穩定性造成重大影響。從剖面塑性區分布圖（圖9，圖10）來看，僅僅在坡體的后緣和坡腳處出現了零星的塑性區，這表明堆積體處于正常的工作狀態。
m剖面塑性區分布圖坡腳處的剪切塑性區位于表面堆積體內，厚度大約5～10 m。需要強調的是，計算結果顯示的是以Mohr-Coulomb屈服準則為依據的塑性區分布情況。該屈服準則則認為材料進入屈服即破壞。
423剪應力（剪應變增量）規律分析判斷堆積體的（潛在）滑動面（帶），可根據其剪應力（應變增量）來判斷；剪應力較為集中或剪應變增量較大（絕對值）的部位，則為其（潛在）滑動面（帶）。變形破壞也都沿此處發生；剪應力分布較為分散（均勻）或者剪應變增量較小或基本上沒發生變化的部位，一般不會有潛在滑動面產生，因此，這些部位也不會發生較大的變形和破壞。應力集中帶主要出現在3個范圍內，坡體上下游側坡腳處和j1夾層附近，它們是堆積體最有可能發生破壞的部位，見圖11。J1夾層附近出現剪應變增量增高帶，說明坡體有可能的破壞模式是沿著以j1夾層附近為滑面進行滑動，但由于計算得到的堆積體整體穩定系數Fs為110，所以發生整體性滑動的可能性不大；坡體上下游側坡腳處的兩個剪應變增量增高帶，說明坡體可能沿土層內部局部滑弧滑動。Z=730剪應變增量剖面見圖12。
1） 根據強度折減法，得到的堆積體整體穩定系數Fs為110，所以發生整體性滑動的可能性不大。文章對滑坡主應力，塑性區，剪應力和位移四個方面的計算結果進行了詳細分析，并且四個方面的計算結果都相互吻合相互印證。
2） 可以直觀地在模型上看出在蓄水條件下，上下游坡腳處的碎石土堆積層是Ⅱ區的主要不穩區域，具體在兩個位置：①下游側坡腳處，高程（Y坐標）范圍大約位于420～425 m之間，走向方向（Z坐標）范圍大約處于130～150 m之間，②上游側坡腳處，高程（Y坐標）范圍大約位于410～415 m之間，走向方向范圍（Z坐標）范圍大約處于760～770 m之間，對后期堆積體的加固治理方案提供了直觀詳細的位置。3） 本文在模擬蓄水條件下的地下水位時利用了FLAC3D中的單元界面能自動依附于指定范圍內模型表面生成的特性，通過建立輔助單元生成水面并同時生成靜水壓力，思路清晰明了，適用于比較復雜的坡體模型。
4） 利用FLAC3D軟件對堆積體進行穩定性分析具有方便，經濟，直觀的優勢，但是由于建立的模型不可能完全模擬出現實堆積體的所有特性，選取的計算參數也不可能完全與實際相符，所以需要結合其他技術手段綜合評價，才能得出更為全面正確的結果。
參考文獻：
[1]李素梅. 塊體理論在巖石路塹邊坡穩定分析中的應用[J]. 云南交通科技， 2003， 19（4）：14-16.
[2]張均峰， 丁樺.邊坡穩定性分析的三維極限平衡法及應用[J]. 巖石力學與工程學報， 2005，24（3）：365-370.
[3]楊強， 朱玲， 薛利軍. 基于三維多重網格法的極限平衡法在錦屏高邊坡穩定性分析中的應用[J]. 巖石力學與工程學報， 2005， 24（增2）：5 313-5 318.
[4]陳昌彥， 王思敬，沈小克. 邊坡巖體穩定性的人工神經網絡預測模型[J]. 巖土工程學報， 2001， 23（2）：157-161.
[5]賈厚華， 賀懷建. 邊坡穩定模糊隨機可靠度分析[J]. 巖土力學， 2003， 24（4）：657-660.
[6]陳新民， 羅國煜. 基于經驗的邊坡穩定性灰色系統分析與評價[J]. 巖土工程學報， 1999， 21（5）：638-641.
[7]祝玉學. 邊坡可靠性分析[M]. 北京：冶金工業出版社， 1993.
[8]黃潤秋，許強. 顯式拉格朗日差分分析在巖石邊坡工程中的應用[J]. 巖石力學與工程學報， 1995， 14（4）：346-354.
[9]鄒棟， 鄭宏. 快速拉格朗日法及其在邊坡穩定性分析中的應用[J]. 礦業研究與開發， 2005， 25， （5）：80-83.
[10]陳育民， 徐鼎平. FLAC/FLAC3D基礎與工程實例[M]. 北京： 中國水利水電出版社， 2008.
