四川綿竹走馬嶺特大泥石流動力特性的數(shù)值分析

陰 可，茍印祥，2，李承澤，吳漢輝，3

(1.重慶大學土木工程學院，重慶 400045;2.四川省交通運輸廳交通勘察設計研究院，四川成都 610017;3.重慶蜀通巖土工程有限公司，重慶 401147)

0 引言

我國是一個地質災害多發(fā)國家，其中泥石流這種地質災害對人民生命和財產(chǎn)造成的損失是十分慘痛的，關于泥石流的防治和研究也是長期的、持續(xù)的。目前，泥石流的防治工程研究主要有:泥石流的評估及預測預報［1-3］，泥石流流體模型研究［4-5］，泥石流的數(shù)值模擬分析［6-10］等，其中泥石流的數(shù)值模擬分析隨著計算機性能的大幅提高以及泥石流流體模型研究的進一步發(fā)展而逐步開始運用，三維泥石流的數(shù)值模擬不僅能直觀反映泥石流流動的整個過程，而且能計算出流動過程中的泥石流流場變化情況，并且可以與泥石流攔砂壩等防治工程設施進行流固耦合分析，為泥石流相關研究提供一些參考。

四川省綿竹市走馬嶺特大泥石流位于距綿竹市區(qū)約35km的清平鄉(xiāng)銀杏溝旅游區(qū)內鹽井村走馬嶺，為綿遠河左岸支流。該溝整個流域面積為5.7km2，最高點海拔高程2044m，最低點位于走馬嶺溝匯入綿遠河口，高程904m，相對高差1140m。主溝長度為3.5km，從溝口至溝源溝床坡降為85‰～361‰，溝源處坡降局部達620‰，平均溝床坡降為145‰，溝坡坡度一般為30°～45°。流域基巖地層依次為:寒武系下統(tǒng)清平組長石云母石英砂巖及頁巖、震旦系上統(tǒng)邱家河組白云巖和純灰?guī)r、砂頁巖夾泥質灰?guī)r、石英砂巖等。溝內覆蓋層為:坡積—洪積層卵石土，厚約4～15m，泥石流堆積層，溝槽內殘留的堆積層厚度為2～3m，溝口的堆積層厚度為3～12m。走馬嶺流域屬四川盆地西北部的龍門山推覆構造帶前緣，處于2008年“5.12”地震地震帶內。

2009年9月24日，清平鄉(xiāng)普降暴雨，引發(fā)“9.24”泥石流災害，2010年8月13日再降暴雨，最大每小時降雨量為70mm(略高于50a一遇)，日降雨量為227.5mm，走馬嶺再次爆發(fā)泥石流，一次性沖出量達8.3×105m3。泥石流將溝谷堵塞填高，堆積體厚度達8～15m，局部甚至達20m，河床平均淤積抬高7～10m。災害造成多幢房屋被掩埋，公路被沖毀，綿遠河改道，直接經(jīng)濟損失3000×104元以上。

走馬嶺特大泥石流治理工程擬采取“攔蓄+固源+排導+停淤+監(jiān)測預警+搬遷避讓”為主的工程措施。考慮在走馬嶺主溝上共修建7座攔砂壩(有效壩高8～12m，上游面坡坡度1∶0.70，下游面坡坡度1∶0.20，壩頂軸線長度 50～71m，1#～4#壩為天然地基，5#～7#壩采用樁基礎，壩身均采用C35塊石混凝土，密度為2400kg/m3)及其他防治措施，包括谷坊壩、兩岸防沖墻、泥石流停淤場和泥石流排導槽等。為了給該方案及設計提供一定的參考和優(yōu)化建議，特進行泥石流動力特性的數(shù)值模擬計算和分析。

1 泥石流動力特性的三維數(shù)值分析方法與思路

1.1 泥石流數(shù)值模擬的理論基礎

泥石流三維數(shù)值分析的理論基礎是計算流體動力學(Computational Fluid Dynamics，簡稱 CFD)，它可以看作是在流體三大基本方程(質量守恒方程、動量守恒方程、能量守恒方程)的控制下的流體數(shù)值模擬，通過模擬計算得到復雜流體問題的流場中各點處的速度、壓力、濃度等基本物理量的分布，以及這些物理量隨時間變化的規(guī)律等。本文采用基于有限體積法(又稱控制體積法)的軟件CFX進行流場計算和基于有限單元法的ANSYS軟件進行流固耦合計算。這里有限體積法的基本思路是:將計算區(qū)域劃分為網(wǎng)格，并使每個網(wǎng)格點周圍有一個互不重復的控制體積，將待解的控制方程對每一個控制體積積分，從而得出一組離散方程，其關鍵是在導出離散方程過程中，需要對界面上的被求函數(shù)本身及其導數(shù)的分布作某種形式的假定［11］。

1.2 泥石流漿體的流變模型

泥石流流體是非牛頓流體，近二十年來，泥石流流體模型研究有了較大的發(fā)展，常見的有:賓漢體模型、膨脹體模型、粘塑體模型、顆粒流模型、兩相流模型［4］等幾種。根據(jù)走馬嶺泥石流的特征和數(shù)值分析軟件的需要采用賓漢體模型，即簡化認為泥石流流體的流變特性符合以下應力應變關系:

1.3 泥石流數(shù)值模擬的流程

(1)建立泥石流流域幾何模型;(2)確定初始條件及邊界條件，建立計算域;(3)劃分網(wǎng)格，確定求解控制參數(shù);(4)有限體積法求解泥石流流場，分析結果;(5)將計算獲得的流場壓力等數(shù)據(jù)導入ANSYS軟件，加載到攔砂壩上;(6)有限元法求解攔砂壩的應力、應變及安全系數(shù)等;(7)整理計算結果，總結分析成果，提出相關意見和建議等。

2 泥石流數(shù)值分析結果

2.1 泥石流流場數(shù)值分析

走馬嶺泥石流數(shù)值分析采用的泥石流參數(shù)為:重度γ=1800kg/m3，進口泥石流流量依次為，主溝Q0=74.62 m3/s，①號支溝 Q1=725.83 m3/s，②號支溝Q2=40.90 m3/s，④號支溝 Q3=56.71m3/s，⑤號支溝 Q5=124.09 m3/s，⑥號支溝 Q6=64.35 m3/s，⑦號支溝Q7=57.00 m3/s。模擬湍流類型為k-ε模型，泥石流流體的流變特性參數(shù)采用CFX軟件自帶的非牛頓流體參數(shù)設置對應的選項。泥石流模擬計算的典型流動狀態(tài)見圖1，對河床壓力的等值云圖見圖2，流域的流速矢量圖見圖3。

模擬計算獲得的走馬嶺泥石流主溝斷面流場數(shù)據(jù)見表1。

圖1 泥石流239.99s和1634.15s時刻的流動狀態(tài)模擬圖Fig.1 Debris flow state simulation at 239.99s and 1634.15s

圖2 泥石流對河床壓力云圖(1634.15s)Fig.2 Debris flow pressure cloud on river bed(1634.15s)

圖3 泥石流流域流速矢量圖(1634.15s)Fig.3 Debris flow velocity vector diagram(1634.15s)

表1 走馬嶺泥石流主溝斷面流場數(shù)據(jù)表Table 1 Debris flow Cross-section field data at the master valley of Zoumaling

模擬計算獲得的走馬嶺泥石流5#、6#、7#支溝斷面流場數(shù)據(jù)見表2(各個支溝的橫斷面從各支溝上游到該支溝與主溝匯流處依次按200～300m間距布置)。

從以上走馬嶺泥石流數(shù)值計算結果分析表明:泥石流的流速從上游到下游呈逐漸增大的趨勢，同時還受河床坡降和河面寬度的影響，河床坡降越大或河床寬度越小的泥石流流速越大。泥石流的泥深隨流速的增加而逐漸減小，同時也受河床坡降和河床寬度影響，河床坡降越大或河床越寬泥深則越小。

表2 走馬嶺泥石流5#、6#、7#支溝斷面流場數(shù)據(jù)表Table 2 The slaves valleys of Zoumaling debris cross-section flow field data

2.2 泥石流流場與攔砂壩的流固耦合計算分析

采用CFX軟件模擬計算了走馬嶺泥石流的流場數(shù)據(jù)，將獲得的泥石流流場數(shù)據(jù)導入ANSYS軟件中進行流固耦合計算，即可得出按設計方案考慮設置的泥石流攔砂壩上作用的流體壓力和流體的運動特征，從而計算得到攔砂壩隨泥石流流動產(chǎn)生的應力、應變等結果。

選取1#泥石流攔砂壩作為代表，滿庫時泥石流流體壓力場分布見圖4。

圖4 泥石流攔砂壩滿庫時流體壓力形態(tài)圖(1#攔砂壩)Fig.4 Debris’s pressure and shape when the debris flow through the dam(1#dam)

可獲得各攔砂壩最大受力狀態(tài)見圖5。

圖5 1#～7#攔砂壩最大受力狀態(tài)統(tǒng)計圖Fig.5 The maximum stress state of 1#～7#dams

泥石流流場數(shù)值模擬顯示，泥石流在攔砂壩前的泥深隨著時間的增加逐漸加大，滿庫后漫過攔砂壩向下游流出，將泥石流瞬態(tài)模擬的各個時間點對攔砂壩的壓力和泥石流泥深數(shù)據(jù)導入ANSYS軟件有限元結構計算模塊，安全系數(shù)計算中，攔砂壩素混凝土采用材料的極限強度，通過計算得到攔砂壩安全系數(shù)與泥深關系見圖6。

圖6 1#～7#攔砂壩安全系數(shù)與泥深關系圖Fig.6 Safety factors and debris deep data’s relationship of 1#～7#dams

攔砂壩中1#～4#攔砂壩是采用的天然地基擴展基礎，基礎前部有榫齒狀構造，5#～7#攔砂壩主要采用的是樁基礎，基礎底部平整規(guī)則。從圖4～圖6可知，1#～4#攔砂壩受力狀態(tài)以及安全系數(shù)隨泥深、壩高的變化規(guī)律具有一致性。5#～7#攔砂壩受力狀態(tài)以及安全系數(shù)隨泥深、壩高的變化規(guī)律也具有一致性，1#～4#攔砂壩的安全系數(shù)總體上小于5#～7#攔砂壩的安全系數(shù)。

3 結論與建議

目前，還沒有關于泥石流動力特性數(shù)值計算完備的理論，以及專門針對泥石流數(shù)值分析的成熟軟件，因此泥石流的數(shù)值分析結果尚不能完全準確地定量反映泥石流流場特性，但是利用CFD軟件計算的結果可以反映泥石流流動的一般規(guī)律，對泥石流的學術研究和防治工作具有一定的指導意義和參考價值。

本文走馬嶺泥石流動力特性的數(shù)值計算和分析結果表明:

(1)泥石流對河床的最大壁面剪切應力為8.32×103Pa，位于⑤號支溝中上游;5#～7#攔砂壩下游300m段泥石流對河床壁面剪切應力也較大，達4.1×103～8.0×103Pa。

(2)泥石流的最大流速是13.57m/s，位于⑤號支溝中游;在5#～7#攔砂壩下游300m段，泥石流流速也相對較大，達到5～8m/s。

(3)泥石流對河床的最大壓力為1.49×105Pa，位于②號支溝與主溝交匯處;其次在⑥號支溝上游，5#～7#攔砂壩下游300m段泥石流對河床的壓力也較大，達2.2×104～1.0×105Pa。

(4)隨著攔砂壩壩體高度的增加，在泥石流的沖擊作用下，攔砂壩的安全系數(shù)減小;當攔砂壩過高時，壩體底部局部區(qū)域混凝土的應力在自重作用下就接近其材料的極限強度，壩體安全性較低。

(5)采用樁基礎的5#～7#攔砂壩安全系數(shù)高于采用天然地基或擴展基礎的1#～4#攔砂壩;隨著泥石流泥深的增加，壩體的安全系數(shù)逐步降低。當泥石流滿庫過流時，各攔砂壩的安全系數(shù)最低。

根據(jù)上述計算分析結果，建議走馬嶺泥石流治理工程應注意以下幾點:

(1)2#、4#、5#、6#支溝兩岸的河流岸坡應加強防沖刷及固岸措施。2#～7#支溝之間應當加強河床的清理疏通工作，這有利于形成相對穩(wěn)定的流動及降低泥石流下蝕與側蝕作用。

(2)泥石流出口段有較大面積的堆積區(qū)，該區(qū)內泥石流的速度逐漸降低，且泥石流流動呈面狀鋪開，可以開發(fā)出來做泥石流的停淤場，但需要將場地進行整理，加大停淤區(qū)的地形橫向坡度，以利于泥石流物質進一步橫向流動以及停淤。

(3)5#～7#攔砂壩壩體結構極限強度安全系數(shù)在4.5(壩體素混凝土采用極限強度)以上，也就是說，此時壩體部分的材料強度滿足安全性要求，并且有少量的富余。1#～4#攔砂壩壩體結構極限強度安全系數(shù)為2～3(壩體素混凝土采用極限強度)，受力位置在壩底轉角處，因此建議適當加強壩體強度，或采取降低壩高，分多級建壩，優(yōu)化壩體結構等措施以增強壩體安全儲備，保證攔砂壩的長期可靠性。

［1］韋方強，胡凱衡，J L Lopez.泥石流危險性分區(qū)及其在泥石流減災中的應用［J］.中國地質災害與防治學報，2007，18(1):23-27.WEI Fangqiang，HU Kaiheng，J L Lopez.Debris flow risk zoning and its application in disaster mitigation［J］.The Chinese Journal of Geological Hazard and Control，2007，18(1):23-27.

［2］唐川，朱大奎.基于GIS技術的泥石流風險評價研究［J］. 地理科學，2002，22(9):300-304.TANG Chuan，ZHU Daku.Assessment of debris flow risk ofYunnan province by using GIS ［J］.Scientia Geographica Sinica，2002，22(9):300-304.

［3］劉涌江，胡厚田，白志勇.泥石流危險度評價的神經(jīng)網(wǎng)絡法［J］.地質與勘探，2001，37(2):84-87.LIU Yongjiang，HU Houtian，BAI Zhiyon.Artificial neural network method for evaluating the dangerous degree of bebris flows［J］.Geology and Prospecting，2001，37(2):84-87.

［4］王光謙，倪晉仁，張軍，等.泥石流的顆粒流模型［J］. 山地研究，1992，10(1):1-10.WANG Guangqian，NI Jinren，ZHANG Jun，et al.Grainular flow model for debris flow［J］.Journal of Mountain Research，1992，10(1):1-10.

［5］倪晉仁，王光謙.泥石流的結構兩相流模型:I.理論［J］. 地理學報，1998，53(1):66-76.NI Jinren，WANG Guangqian.Conceptual two-phase flow model of debris flow:Ⅰ Theory［J］.Acta Geographica Sinica，1998，53(1):66-76.

［6］馬宗源，廖紅建，張駿.Bingham型黏性泥石流流體的三維數(shù)值模擬［J］.西安交通大學學報，2008，42(9):1146-1150.MA Zongyuan， LIAO Hongjian， ZHANG Jun.Three dimensional numerical simulation of Bingham viscous debris flow fluid［J］.JournalofXi’an Jiaotong University，2008，42(9):1146-1150.

［7］何菊明，彭云雄.液-固兩相流的運動描述與數(shù)值模擬［J］.湖北師范學院學報(自然科學版)，2006，26(3):26-29.HE Juming，PENG Yunxiong.Motion description and numerical simulation of fluid-solid two-phase flows［J］.Journal of Hubei Normal University(Natural Science)，2006，26(3):26-29.

［8］陳春光，姚令侃.泥石流與主河交匯區(qū)三維數(shù)值模擬［J］.重慶交通學院學報，2006，25(2):61-65.CHEN Chunguang，YAO Lingkan.Three dimensional numerical simulation of confluence of debris flow and main river［J］.Journal of Chongqing Jiaotong University，2006，25(2):61-65.

［9］劉學，王興奎，王光謙.基于GIS的泥石流過程模擬三維可視化［J］. 水科學進展，1999，10(4):388-392.LIU Xue， WANG Xingkui， WANG Guangqian.GISBased 3-D visualization of debris process simulation［J］.Advances in Water Science，1999，10(4):388-392.

［10］王純祥，白世偉，江崎哲郎，等.泥石流的二維數(shù)學模型［J］.巖土力學，2007，28(6):1237-1241.WANG Chunxiang，BAI Shiwei，ESAKI Tetsuro，et al.Two-dimensional mathematical model of debris flow［J］.Rock and Soil Mechanics，2007，28(6):1237-1241.

［11］王福軍.計算流體動力學分析 CFD軟件原理與應用［M］.北京:清華大學出版社，2004:272.WANG Fujun.Analysis of computational fuild dynamics analysis［M］.Beijing Tsinghua University press，2004:272.