基于Flow-R和FLO-2D的溝谷型泥石流危險性評價

林文 周偉 李婧 楊莉 唐澤

摘要：泥石流是中國西南山區頻發的地質災害之一，具有突發性、流速快以及破壞力強等特點，往往造成巨大的人員傷亡和經濟損失。地質災害危險性評估是防災減災管理和防治環節中的有效措施之一，用于合理量化泥石流運動堆積線路沿程泥石流災害危險性的空間分布特征。以四川省都江堰市龍池鎮麻柳溝為研究對象，選取坡度、平面曲率和上坡匯水面積作為溝谷型泥石流啟動區識別因子，基于Flow-R模型識別出的泥石流啟動點，耦合FLO-2D模型模擬麻柳溝泥石流運動堆積過程，將泥石流堆積深度、運動速度和暴發頻率相結合分析泥石流的危險性，并采用混淆矩陣對模擬結果進行了精度評估。結果表明：基于Flow-R和FLO-2D耦合模型的麻柳溝泥石流的運動堆積模擬結果與實際情況吻合度較高。模擬泥石流堆積區面積2.19萬m2，堆積體積4.92萬m3，泥石流危險性由溝道中心向兩側逐漸降低，高危險區占比42.3%，中危險區占比23.1%，低危險區占比34.6%。麻柳溝形成區、流通區溝道中散落堆積物體積大約為30萬m3，麻柳溝泥石流仍處于高度活躍階段，危險性較高，今后在降雨和水動力充分的條件下大量松散物質可能將被重新激活，一旦發生泥石流將有可能造成重大地質災害。

關 鍵 詞：泥石流; 危險性評價; 泥石流啟動點; Flow-R; FLO-2D; 都江堰市

中圖法分類號： P642

文獻標志碼： A

DOI：10.16232/j.cnki.1001-4179.2022.05.023

0 引 言

泥石流是中國西南山區頻發的地質災害之一，具有突發性、流速快以及破壞力強等特點，往往會造成巨大的人員傷亡和經濟損失。泥石流危險性評價是災害學重點研究的課題之一。泥石流活動規模是泥石流危險性評價的一個重要指標[1]，它取決于控制流量觸發的運動特征、泥石流沿程侵蝕和泥石流峰值流量等因素。FLO-2D是一種被廣泛用于模擬洪水的軟件，近年來也被用來模擬溝谷型泥石流的運動和堆積，便于對溝谷型泥石流進行危險性評價。阮德修等[2]將FLO-2D與3DMine進行耦合，模擬了尾礦庫潰壩時的泥石流流動過程，并預測了尾礦庫潰壩和半潰壩情況下的致災程度。Chang等[3]通過實地調查，分析了汶川地震災區泥石流的形成條件，運用FLO-2D模型，模擬了泥石流運動及其危害范圍。張建石[4]通過統計分析，得出震后多年間的物源演變規律及特征，并利用FLO-2D模型模擬在不同降水頻率下肖家溝泥石流的沖出特征。Nocentini[5]將FLO-2D模型與DAN-W模型進行耦合，預測了研究區潛在泥石流的危險性。黃靖玲等[6]將HEC-HMS模型與FLO-2D模型耦合，對馬貴河流域進行了危險性評估。如果要獲得較好的FLO-2D模型模擬效果，關鍵在于能找出相對準確的泥石流啟動點。這就需要結合實地情況和專業經驗，還應考慮到一定的客觀因素。為了解決泥石流啟動點識別的問題，本文引入了一種識別溝谷型泥石流啟動點的軟件Flow-R。該軟件是一種基于GIS的重力災害評價模型，可高效快速識別溝谷型泥石流啟動點[7]。

四川省龍池鎮龍溪河流域自“5·12”地震后，流域內松散物源量增加，地質災害活動頻繁，大量研究表明，泥石流災害已成為汶川地震重災區龍池鎮最大的地質安全隱患[8]。2010年8月13日，流域在強降水條件下暴發泥石流災害，造成巨大損失。本文以麻柳溝為研究對象，將FLO-2D模型與Flow-R模型耦合使用，來識別泥石流大多啟動點，模擬泥石流淹沒區范圍，計算泥石流運動的最大動量和最大堆積深度[9];同時，結合流速、泥深與降水暴發頻率來探究麻柳溝的危險性。研究成果可為麻柳溝的泥石流災害的危險分析、預警及發展規劃提供科學的參考依據。

1 研究區概況

麻柳溝位于四川省都江堰市龍池鎮，溝口距離龍池鎮約0.65 km。溝口地理坐標為103°33′20″E，31°03′55.11″N（見圖1）。麻柳溝流域面積約為0.93 km2，流域相對高差為913 m，主溝長度為2.3 km，平均縱坡降401‰。巨大的相對高差對流域內的松散堆積物而言具有足夠的勢能條件，使得該溝易形成危害極大的高動能泥石流。

研究區地處龍門山斷裂帶。2008年5月12日汶川縣映秀鎮發生8.0級大地震，其主發震斷裂映秀-北川斷裂橫貫研究區，導致研究區內巖層破碎嚴重，成為典型的構造不穩定區域。震后流域內發生大量崩塌滑坡，為暴雨泥石流提供了充足的松散物質來源。

研究區屬于中亞熱帶濕潤季風氣候區，日照時間較少，降水充沛，陰雨天氣頻繁。多年平均降水量為1 134.8 mm（1987～2008年），5～9月降水量占全年降水量的80%。常年充分的降水使坡面松散堆積物直接運移到溝道堆積，形成泥石流溝道物源。2010年8月13日龍池鎮發生強降水事件，最大1 h降水量達75.0 mm，連續3 h降水量達到150.0 mm，導致泥石流暴發。

2 研究方法

Flow-R模型將坡度、平面曲率、上坡匯水面積作為溝谷型泥石流啟動點的識別因子，然后將其轉化為Flow-R模型能識別的ASCII文件，通過對不同的因子設置閾值，識別出研究區的啟動點。

基于DEM數據轉化為FLO-2D能識別的ASCII文件，結合現場實地調查進行網格劃分。調查發現研究區泥石流溝道流動寬度大于10 m，所以將本次數值模擬網格劃分為10 m的網格滿足本次研究的需要，對每個網格賦值唯一的高程值、粗糙系數等基礎地形因子，將Flow-R識別的啟動點導入網格，輸入啟動點泥石流流量及流變參數，模擬麻柳溝“8·13”特大泥石流運動堆積過程，構建混淆矩陣對運動堆積范圍進行精度驗證。

危險性評價建模是將耦合模型結果獲取的淹沒區流速、泥深和降水頻率相結合，將泥石流溝劃分為高危險區、中危險區和低危險區。工作流程如圖2所示。

2.1 泥石流啟動點識別

坡度（地形條件）、松散物源可用性（物源條件）、輸水量（水源條件）是影響泥石流形成的關鍵因素，當這些關鍵因素滿足一定條件時即可誘發泥石流[10]。Flow-R模型考慮了泥石流形成的三大影響因素，通過設置各個因子泥石流啟動識別閾值（詳見3.1節），將同時滿足閾值的柵格單元劃分為“啟動區”，將至少有一次不滿足閾值的柵格單元劃分為“不是啟動區”。

2.2 泥石流運動流通堆積模擬

FLO-2D模型用于模擬泥石流的運動過程和跳動特征[11]。模型以網格為單元計算泥石流流動，在軟件劃分出大小相同的規則網格后，將高程、粗糙系數等數值賦予給各網格。二維流動是通過運動方程的數值積分和流體體積守恒來實現的，即網格中固體物質和水的比值不斷變化，而總體積變化，且在模擬過程中，當低濃度流體進入網格遇到已停淤的高濃度流體時，匯流處流體發生混合并繼續流動，控制方程如下。

（1） 連續方程。

ht+（uh）x+（vh）y=l（1）

式中：l表示水力坡降;h表示泥深;v表示泥石流在y方向的速度;u表示泥石流在x方向的速度。

（2） 運動方程。

Sfx=Scx-hx-VxgVxx-VygVyy-1gVxt（2）

Sfy=Scy-hy-VygVyy-VygVyx-1gVyt（3）

式中：Sfx和Sfy分別表示x、y方向的摩擦坡降;Scx和Scy分別表示x、y方向的底床坡降。

（3） 流變方程。

τ=τc+τmc+τt+τv+τd（4）

式中：τc表示黏性屈服應力;τv表示黏性剪應力;τmc表示摩爾-庫倫剪應力;τt表示紊流剪應力;τd表示擴散剪應力。

2.3 泥石流危險性評價

對單溝泥石流危險性評價，根據唐川等的研究，結合泥石流運動速度和泥石流深度對泥石流強度進行劃分[12]（見表1）。在此基礎上，Chang等[3]以泥石流強度與暴發頻率相結合的方式進行泥石流危險性評價，將研究區危險性分為高、中、低3個等級（見圖3）。本文將耦合模型模擬結果結合泥石流運動深度、速度以及降水頻率，對研究區的泥石流危險性進行評價。

3 結果分析

3.1 泥石流啟動點識別

（1） 坡度。

坡度是泥石流形成的主要因素。大多數地區發生泥石流的地形坡度大于15°[6]，Ortigao等[13]認為，泥石流啟動時的地形坡度介于20°～25°之間。麻柳溝坡度主要分布在20°～50°，占整個流域的86.41%。結合麻柳溝地形坡度特征，將泥石流啟動的地形坡度閾值定為20°。

（2） 平面曲率。

麻柳溝泥石流是典型的溝谷型泥石流，泥石流的物源一般分布在曲率為凹形的區域。平面曲率表示地形的凹凸變化，可用于識別泥石流溝道等負地形，間接體現泥石流的物源。現有研究表明，泥石流啟動時的平面曲率閾值大多在-0.5/100 m-1到-2/100 m-1之間[14]。Delmonaco等[15]在平面曲率為-1/100 m-1的凹地形中識別出泥石流物源。麻柳溝負地形平面曲率統計結果表明，該溝平面曲率大于-3/100 m-1的比例達到84%。本文將研究區平面曲率閾值定為-3/100 m-1。

（3） 匯水面積。

上坡匯水面積，即匯流累積量，表示DEM中流經每個柵格的流水累積量，旨在識別研究區內任何活躍的河流、溝渠或暗流。Park等通過觀測實際發生的泥石流物源，構建了泥石流源區坡度與匯水面積的閾值預測模型[16]。該模型如下：

θ=15e32.34A-160.85（5）

式中：A表示匯水面積，θ表示坡度。參考公式，通過設定100，500 m2和1 000 m2匯水面積閾值提取水系，結果表明1 000 m2閾值最符合實際情況。對麻柳溝流域進行了匯水面積統計，結果顯示，該溝78%的匯水面積在0～1 000 m2范圍內，因此，將研究區匯水面積率閾值定為1 000 m2。

由圖5啟動點模擬結果可知：麻柳溝泥石流啟動區主要分布在溝道及溝道兩側，啟動點坡度分布在20°～45°之間，其中，超過一半的啟動點分布于30°～40°之間，啟動點高程主要分布在1 100～1 300 m。溝頂細小匯水溝道為泥石流提供了豐富的水源，促使下部溝道沉積物轉化為泥石流。

3.2 泥石流運動流通堆積模擬

參考水文地質手冊，結合實地考察，泥石流啟動點主要分布在楊家溝和磯子溝。采用雨洪法分別計算了楊家溝和磯子溝處啟動點的洪峰流量。泥石流在運動過程中挾帶有大量的松散物質，參考FLO-2D使用手冊[11]，泥石流中水和沉積物總混合物的體積可以通過將水體積乘以放大系數（BF）來確定。放大系數 BF 與泥石流的體積濃度密切相關[10]，采用野外配漿實驗獲取體積濃度為0.58，計算出麻柳溝啟動點流量。對洪峰流量采用簡單概化的五邊形方法來求取泥石流峰值流量過程線[17]。

泥石流的屈服強度τy及黏滯參數η與堆積范圍相關[18]。τy、n和η可由下式計算：

n=0.33C-0.15VeCv-0.15lnh（6）

η=α1eβ1Cv（7）

τy=α2eβ2Cv（8）

式中：α1、α2、β1和β2為經驗系數;Cv為泥石流體積濃度;h為泥深;n為曼寧系數。曼寧系數常被用來確定水和泥石流固體物質之間的接觸關系。

對麻柳溝進行現場調查發現，溝道內植被稀疏，而且固體堆積物平均厚度在10 cm以上。曼寧系數可以間接表現流體的流變特征，根據公式（6）計算出的麻柳溝曼寧系數為0.12。層流阻滯系數K、α和β均為經驗系數，一般通過查閱 FLO-2D 使用手冊或進行專項試驗獲取。本次研究，根據前人對經驗系數的取值[19]以及FLO-2D使用手冊[11]中的參數值列表，獲得了麻柳溝FLO-2D模擬的主要參數（見表3）。

麻柳溝數值模擬結果如圖6所示。由圖6（a）可知：溝口的堆積面積為2.19×104 m2，占總流域面積的2.3%，最大堆積厚度為4.5 m，平均堆積深度為2.25 m，沖出總量為4.92×104 m3。研究區最大堆積深度達5.90 m，最小堆積深度為0.04 m，溝道中心泥深最大，向溝道兩側泥深減少，溝道中和堆積扇前段泥深較大，速度較快，溝道兩側和堆積扇翼端堆積深度較小，運動速度較慢。

利用混淆矩陣對模擬結果進行精度評價，并采用準確率（A）和敏感度（S）指標對模擬結果進行評價，準確率和敏感度越高，表明模擬結果越好[20]。

A=TP+TNN（9）

S=TPTP+FN（10）

式中：TP為實際淹沒區與模擬淹沒區均覆蓋的柵格數;TN為實際淹沒區與模擬淹沒區均未覆蓋的柵格數;FN為實際淹沒區覆蓋而模擬淹沒區未覆蓋的柵格數;N為研究區柵格總數。根據麻柳溝2011年4月26日谷歌影像提取麻柳溝泥石流淹沒區范圍（見圖6（b）），根據公式（9）和公式（10）計算出的模擬結果準確率為91.0%，敏感度為72.3%。表明耦合模型應用于麻柳溝泥石流的運動堆積模擬的結果與實際情況吻合度較高。

3.3 泥石流危險性評價

模擬還原了麻柳溝“8·13”泥石流運動堆積過程。由麻柳溝泥石流危險性分區圖（見圖7）可知：高危險區主要分布在溝道中心，面積約為58 300 m2，占危險區總面積的42.3%。由溝道中心向兩側泥石流危險性逐漸降低，中危險區面積為31 200 m2，占總面積的23.1%;低危險區面積為45 400 m2，占總面積的34.6%。高危險區域分布范圍較大，低危險和中危險分布范圍相對較小。溝道兩側和堆積扇翼端主要為低危險區域;泥石流形成區和堆積區邊緣主要為中危險區;高危險區則主要集中在溝道中和堆積扇前段。在形成區和流通區的松散堆積體總量約30萬m3，在極端降水天氣下極易發生次生泥石流災害，應該加強防災減災工作和泥石流風險管理工作。

4 結 論

本文以麻柳溝為研究對象，基于Flow-R模型識別泥石流啟動點后耦合FLO-2D模型，用以模擬泥石流運動;采用速度、泥石流堆積深度和降水暴發頻率相結合的方法，對單溝泥石流危險性進行了評價，可以得出如下結論。

（1） 麻柳溝泥石流啟動區絕大部分布在溝道及溝道兩側，啟動點坡度分布在20°～45°之間，其中超過一半的啟動點分布于30°～40°之間，高程主要分布在1 100～1 300 m范圍內。

（2） 模擬泥石流堆積區面積為2.19萬m2，堆積體積為4.92萬m3。計算形成區、流通區溝道中散落堆積物體積大約為30萬m3，耦合模型應用于麻柳溝模擬的準確率為91.0%，敏感度為72.3%，模擬精度較高，可以為土地資源規劃、災害防治設計提供參考依據。

（3） 麻柳溝高危險區分布范圍較大，占危險區面積的42.3%，溝道中心危險性較高，沿溝道中心向兩側危險性降低，高危險區主要集中在溝道中和堆積扇前沿，低危險區主要分布于溝道兩側和堆積扇翼端，泥石流上溝段和堆積區邊緣主要為中危險區。

綜上所述，基于Flow-R和FLO-2D耦合模型將泥石流源區和泥石流運動相結合，可以得出泥石流運動堆積線路沿程災害危險性空間分布特征。根據危險分布特征圖可知，麻柳溝泥石流仍處于高度活躍狀態，今后在降水和水動力充分的條件下，大量的松散物質可能被重新激活，一旦發生泥石流將可造成重大的地質災害。本文研究成果可為泥石流災害監測預警提供建設性考慮。本文研究未考慮物源量、巖性等其他因子對泥石流啟動點識別的影響，以及未考慮泥石流在運動過程中側蝕和下蝕的影響，使得數值模擬的結果稍顯偏小，且較為保守，在今后的研究中需要對其進行補充和完善。

參考文獻：

[1] LUNA Q B，REMATRE A，VAN ASCH T W J，et al.Analysis of debris flow behavior with a one dimensional run-out model incorporating entrainment[J].Engineering Geology，2012，128：63-75.

[2] 阮德修，胡建華，周科平，等.基于FLO-2D與3DMine耦合的尾礦庫潰壩災害模擬[J].中國安全科學學報，2012，22（8）：150-156.

[3] CHANG M，TANG C，VAN ASCH T W J，et al.Hazard assessment of debris flows inthe Wenchuan earthquake-stricken area，southwest China[J].Landslides，2017，14（5）：1783-1792.

[4] 張建石.汶川縣肖家溝泥石流物源演變及沖出規模研究[J].人民長江，2020，51（8）：37-43.

[5] NOCENTINI M，TOFANI V，GIGLI G，et al.Modeling debris flows in volcanic terrains for hazard mapping：the case study of Ischia Island（Italy）[J].Landslides，2015，12（5）：831-846.

[6] 黃靖玲，徐粒，周華真，等.基于歷史災害情景再現的小流域山洪風險研究[J].災害學，2020，35（3）：194-200.

[7] HORTON P，JABOYEDOFF M，RUDAZ B，et al.Flow-R，a model for susceptibility mapping ofdebris flows and other gravitational hazards at a regionalscale[J].Natural Hazards and Earth System Sciences，2013，13（4）：869-885.

[8] 許強.四川省8.13特大泥石流災害特點、成因與啟示[J].工程地質學報，2010，18（5）：596-608.

[9] CHEN H，DADSON S，CHI Y.Recent rainfall-induced landslides and debris flows in northern Taiwan[J].Geomorphology，2018，77：112-125.

[10] RICKENMANN D，ZIMMERMANN M.The 1987 debris flows in Switzerland：documentation and analysis[J].Geomorphology，1993，8（2-3）：175-189.

[11] FLO-2D Software Inc.FLO-2D，2-Dimensional Flood Routine Model Manual[Z].State of Arizona：FLD-2D Software，Inc.，2004.

[12] 唐川，周鉅乾，朱靜，等.泥石流堆積扇危險度分區評價的數值模擬研究[J].災害學，1994，9（4）：7-13.

[13] ORTIGAO J A R，KANJI M A.Landslide classification and risk management[M].Berlin：Springer，2004.

[14] BLAHUT J，WESTEN C，STERLACCHINI S.Analysis of landslide inventories for accurate prediction of debris-flow source areas[J].Geomorphology，2010，119（1-2）：36-51

[15] DELMONACO G，LEONI G，MARGOTTINI C，et al.Large scale debris-flow hazard assessment：a geotechnical approach and GIS modeling[J].Natural Hazards and Earth System Sciences，2003，3：443-455.

[16] PARK D W，LEE SR，VASU N N，et al.Coupled model for simulation of landslides and debris flows at local scale[J].Natural Hazards，2016，81：1653-1682.

[17] 夏添.震區泥石流危險性評價及預警減災系統研究[D].成都：成都理工大學，2013.

[18] O’BRIEN J S，JULIEN P，FULLERTON W.Two-dimensional water flood andmudflow simulation[J].Journal of Hydraulic Engineering，1993，119：244-261.

[19] 程思.都江堰市龍溪河流域震后多溝同發泥石流危險性及易損性研究[D].成都：成都理工大學，2015.

[20] CARRANZA E，CASTRO O T.Predicting lahar-inundation zones：case study in West Mount Pinatubo，Philippines[J].Natural Hazards，2006，37：331-372.

（編輯：趙秋云）

Gully debris flow hazard assessment based on Flow-R and FLO-2D coupling models

LIN Wen，ZHOU Wei，LI Jing，YANG Li，TANG Ze

（State Key Laboratory of Geohazard Prevention and Geoenvironment Protection，Chengdu University of Technology，Chengdu 610059，China）

Abstract：

Debris flow is one of the frequent geological disasters in southwest? mountainous areas of China.It is characterized by sudden occurrence，fast velocity and strong destructive force，which often causes huge casualties and economic losses.Risk assessment of geological hazards is one of the effective measures in disaster prevention and mitigation management and prevention，which can reasonably quantify the spatial distribution characteristics of debris flow hazard along the debris flow movement and accumulation route.Taking Maliu gully in Longchi Town，Dujiangyan City，Sichuan Province as the research object，slope，plane curvature and upslope catchment area were selected as the identification factors for the initiation area of gully debris flow.Based on the starting point of debris flow identified by Flow-R model，the FLO-2D model was coupled to simulate the movement and accumulation process of debris flow in Maliu gully.The debris flow hazard can be determined by combination of accumulation depth，velocity and outbreak frequency of debris flow.In addition，the accuracy of the simulation results was assessed by confusion matrix.The results showed that the simulation results of Maliu gully debris flow by Flow-R and FLO-2D coupling models was in good agreement with the actual conditions.The simulated debris flow accumulation area is 2.19×104m2，and the accumulation volume is 4.92×104m3.The debris flow risk gradually decreases from the center of the channel to both sides.High risk area accounts for 42.3%，medium risk area accounts for 23.1%，and low risk area accounts for 34.6%.In formation area and circulation area of Maliu gully，the volume of scattered accumulation materials in the channel is about 300000 m3.The Maliu gully debris flow is still in a highly active stage with high risk.In the future，under the condition of sufficient rainfall and hydrodynamic force，a large quantity of loose materials may be reactivated，which may cause major geological disasters once it occurs.

Key words：

debris flow;hazard assessment;starting point of debris flow;Flow-R;FLO-2D;Dujiangyan City