基于Massflow的達(dá)摩溝泥石流模擬與危險性分析

楊鑫 孟華君 路璐 王宇 吳季寰

摘要：

達(dá)摩溝小流域歷史上曾發(fā)生泥石流災(zāi)害事件，造成較大人員傷亡和財產(chǎn)損失。目前，該小流域主溝堆積大量煤矸石，主溝兩側(cè)坡面堆積有較多崩、坡積物。為分析泥石流再次發(fā)生的可能性與危險性，在無人機攝影測量和野外現(xiàn)場調(diào)查的基礎(chǔ)上，詳細(xì)分析達(dá)摩溝小流域泥石流發(fā)育條件，基于改進(jìn)MacCormack-TVD有限差分法，采用Massflow軟件對10，20，50，100 a一遇4種降雨概率條件下達(dá)摩溝暴發(fā)泥石流的工況進(jìn)行數(shù)值模擬，得到達(dá)摩溝小流域在不同降雨概率下的泥石流危險性。結(jié)果表明：①? 在此4種降雨概率下，大部分區(qū)域的泥石流流速維持在0～3 m/s，泥深保持在0～4 m的水平，而峰值流量最大可達(dá)19 m/s，泥深最大可達(dá)3.85 m；②? 在100 a一遇的降雨強度下，達(dá)摩溝高危險性地區(qū)占16.54%，中風(fēng)險地區(qū)占49.27%，低風(fēng)險地區(qū)占34.19%；③? 模擬結(jié)果較好地符合泥石流的時空發(fā)育特征，展現(xiàn)了降雨強度對泥石流運動的影響規(guī)律，可為泥石流的預(yù)警監(jiān)測提供參考。

關(guān) 鍵 詞：

泥石流； 危險性； Massflow； 數(shù)值模擬； 達(dá)摩溝

中圖法分類號： P642.23

文獻(xiàn)標(biāo)志碼： A

DOI：10.16232/j.cnki.1001-4179.2023.08.020

0 引 言

1867～1999年間，北京地區(qū)共發(fā)生超過50次泥石流災(zāi)害事件，造成2 400多人死亡，經(jīng)濟損失達(dá)數(shù)億元［1］。近年來，隨著極端天氣頻發(fā)和人類工程活動增強，北京地區(qū)先后發(fā)生2012年“7·21”和2018年“7·16”泥石流災(zāi)害群發(fā)事件［2］，雖處置得當(dāng)，避免了較大人員傷亡，但仍造成了較大財產(chǎn)損失和社會影響。

長期以來，大量專家學(xué)者基于野外調(diào)查、監(jiān)測預(yù)警等取得的一手資料，對北京市山區(qū)泥石流進(jìn)行了大量研究。在20世紀(jì)80年代，洪惜英、李容全、張仲德等就開始研究北京地區(qū)泥石流的發(fā)育特征、分布規(guī)律等基礎(chǔ)性工作［3-5］，旨在通過對泥石流災(zāi)害的科學(xué)認(rèn)識，有效降低地質(zhì)災(zāi)害造成的經(jīng)濟、財產(chǎn)損失，提出防治對策、保障社會安全。20世紀(jì)90年代，相關(guān)研究工作已經(jīng)從泥石流事件的現(xiàn)場調(diào)查和資料統(tǒng)計分析逐漸轉(zhuǎn)向成因機理分析及防治對策設(shè)計［6-7］和泥石流危險性區(qū)劃［8-10］，并嘗試開展北京地區(qū)泥石流預(yù)警預(yù)報工作［11-15］。21世紀(jì)后，隨著仿真模擬技術(shù)、信息傳輸技術(shù)的快速發(fā)展，基于泥石流啟動物理過程的泥石流災(zāi)害精細(xì)化風(fēng)險管理成為可能，不斷出現(xiàn)精度更高的預(yù)警模型。例如，白利平等［16-17］根據(jù)力學(xué)平衡條件，嘗試推出松散固體物質(zhì)起動判別式和雨強表達(dá)式。趙忠海［1］綜合考慮地形地貌、地質(zhì)條件、土壤類型以及降雨情況等因素，建立了泥石流是否發(fā)生的判斷公式。王海芝［18］通過整理1949年以來北京地區(qū)發(fā)生的泥石流與其發(fā)生時期的降水?dāng)?shù)據(jù)，建立了臨界雨量的定量估算I-D模型等。這些工作為北京市泥石流災(zāi)害管控提供了理論指導(dǎo)和技術(shù)支撐。然而，隨著京津冀協(xié)同一體化國家發(fā)展戰(zhàn)略的實施，北京市對泥石流災(zāi)害管控提出了更高的要求，逐漸由災(zāi)害點管控向“隱患點+風(fēng)險區(qū)”雙管控轉(zhuǎn)變。新的管控策略客觀上要求對隱患點泥石流形成、運動及承災(zāi)機理進(jìn)行深入研究［19］，而數(shù)值模擬方法是解決該問題的較好手段。通過數(shù)值模擬，可以再現(xiàn)泥石流運動演化過程，得到用于定量評價泥石流有關(guān)動力學(xué)、運動學(xué)要素和防治工程的參數(shù)［20］。

目前，國內(nèi)外開展泥石流數(shù)值模擬主要采用FLO-2D、Fluent、OpenLISEM及PFC等成熟商業(yè)軟件，研究結(jié)果顯示數(shù)值模擬方法在模擬泥石流運動過程、沖出量、堆積等方面具有較好的效果［21-24］。鑒于此，本文采用數(shù)值模擬方法，選擇北京市達(dá)摩溝作為研究區(qū)，利用更專業(yè)的泥石流模擬軟件Massflow模擬達(dá)摩溝泥石流運動過程，通過研究泥石流隱患點在不同降雨工況下的泥石流運動及沉積過程，厘清達(dá)摩溝泥石流危險性，為該小流域泥石流隱患點風(fēng)險管理提供依據(jù)。

1 研究區(qū)概況

達(dá)摩溝位于北京市門頭溝區(qū)清水鎮(zhèn)東南方向 3 km 處，是清水河支溝，位于東經(jīng)115°37′～115°39′，北緯39°52′～39°56′（見圖1）。該泥石流溝溝域面積為3.96 km2，相對高差647 m，主溝溝向325°，斷面大多呈“V”形，長2 680 m，溝寬30～100 m，平均15 m，溝床平均坡降230%。陡坎、卡口不多，河段寬窄均勻，主支溝交角大，溝道較順直平整。溝床兩岸坡面較陡，坡度30°～60°，平均40°。植被覆蓋率為70%。

達(dá)摩溝內(nèi)出露地層有侏羅系下統(tǒng)窯坡組（J1y）、龍門組（J1l），侏羅系中統(tǒng)九龍山組（J2j）、髫髻山組（J2t）及第四系全新統(tǒng)（Q4）。巖性主要為砂巖、粉砂巖、安山巖、玄武巖、凝灰質(zhì)砂巖、粉砂巖、頁巖及煤系地層等，第四系全新統(tǒng)主要由礫石、巖屑、黏土質(zhì)等殘、坡積物組成［25］。

研究區(qū)地處北京西部山區(qū)，屬中緯度大陸性季風(fēng)氣候，年平均氣溫10.2 ℃，極端最高氣溫37.6℃，極端最低氣溫西部-22.9 ℃。降水量自東向西逐漸減少，受中緯度大氣環(huán)流的不穩(wěn)定性和季風(fēng)影響，降水量年際變化大。

達(dá)摩溝流域歷史上曾多次暴發(fā)泥石流，最近一次發(fā)生于1950年8月4日，當(dāng)日最大降雨量達(dá)190.1 mm，最大雨強為56.5 mm/h。持續(xù)的暴雨造成達(dá)摩莊后港溝和王家港溝泥石流爆發(fā)，造成84人死亡、24人重傷，基礎(chǔ)設(shè)施損毀嚴(yán)重，給當(dāng)?shù)厝嗣裨斐闪舜罅康呢敭a(chǎn)損失［26］。

2 泥石流形成條件

2.1 物源條件

（1） 受達(dá)摩溝歷史煤礦開采影響，主溝及支溝溝道堆積大量煤矸石（見圖2（a）～（b）），其中僅少部分煤矸石堆有漿砌石護坡（圖2（a）），絕大部分煤矸石堆沒有任何防護措施，穩(wěn)定性較差，是達(dá)摩溝未來暴發(fā)泥石流的主要固體物源。

（2） 溝道兩側(cè)岸坡淺表層堆積大量坡積物。這些松散堆積物主要分布于坡面中下部淺表層，厚度不等，約30～80 cm（見圖2（c）～（d）），植被覆蓋較好。然而，根系延伸同時也使堆積更加松散，空隙發(fā)育，加之下伏基巖具有傾向主溝的優(yōu)勢結(jié)構(gòu)面（圖2（d）），當(dāng)降雨量較大時，雨水在緩慢入滲的同時極易沿大空隙灌入坡積物，抵達(dá)土巖接觸面，在該面形成積水，從而引起崩坡積物失穩(wěn)補給泥石流。圖2（c）～（d）可見，坡面下部已經(jīng)缺失了一部分淺表堆積物，下伏基巖面直接出露。

（3） 在坡腳分布有零星崩滑堆積體（圖2（e））和歷史泥石流堆積體（圖2（f））。其中，泥石流堆積體可見明顯的定向排列，大顆粒石塊長軸方向基本一致。這些堆積體規(guī)模大小不等，固結(jié)程度較高，不易失穩(wěn)，但在極端降雨工況下，受上游準(zhǔn)泥石流體沖擊侵蝕，會受到強烈側(cè)蝕作用，可能作為泥石流補給。

2.2 水文條件

達(dá)摩溝泥石流一般由強降雨激發(fā)。該地區(qū)降雨分布具有年際間不均衡、年內(nèi)變化大、降雨集中、降雨強度大等特點。該溝一般年降水量為550～650 mm，平均為575 mm，汛期平均雨量為400～500 mm；歷史年最大降水量932 mm，年最小降水量329 mm。最大日降雨量為224.8 mm，最大小時降雨量達(dá)56.5 mm。每年雨季降水量占全年降水量的70%左右，且多以大暴雨或特大暴雨形式出現(xiàn)。

3 研究方法

數(shù)值模擬是解決單溝泥石流隱患點危險性評價的有效方法。本文采用專業(yè)泥石流模擬Massflow軟件對研究區(qū)泥石流進(jìn)行運動過程模擬。

Massflow是中國科學(xué)院成都山地災(zāi)害與環(huán)境研究所歐陽朝軍研究員研發(fā)的一款數(shù)值模擬軟件。該軟件基于深度積分連續(xù)介質(zhì)力學(xué)理論，利用改進(jìn)MacCormack-TVD有限差分方法，兼顧復(fù)雜地形地貌，是具有二階精度和自適應(yīng)求解域特征的山地災(zāi)害動力學(xué)高效計算模擬軟件。與國內(nèi)外成熟軟件相比，該軟件不僅提供了豐富的物理模型以提高模擬準(zhǔn)確性，還能夠極大提升計算效率，縮短模擬時間。在甘肅岷縣二馬溝泥石流模擬［27］、金沙江白格滑坡模擬［28］以及泥石流防治工程防治效果模擬［29］中，Massflow均取得了良好的效果。

3.1 Massflow模型控制方程

4 模擬結(jié)果與分析

4.1 計算模擬結(jié)果

通過對10%，5%，2%，1%降雨頻率下達(dá)摩溝爆發(fā)泥石流的過程進(jìn)行模擬，可以得到不同降雨條件下達(dá)摩溝泥石流的運動過程和不同時刻的泥深以及流速分布，并在此基礎(chǔ)上開展達(dá)摩溝泥石流隱患點在不同降雨頻率下的危險性分析評價。部分模擬結(jié)果如圖3，4所示。

4.2 模擬結(jié)果分析

4.2.1 流速變化

根據(jù)計算結(jié)果（見表2），不同降雨工況下達(dá)摩溝泥石流流速變化具有以下3個特點。

（1） 總體流速較為穩(wěn)定。在10，20，50，100 a一遇降雨概率下，整個主溝道流速總體普遍維持在0～4 m/s的水平；只有下游的階地以及陡峭區(qū)域流速急劇增大，局部平均流速達(dá)到8 m/s以上（見圖3）。

（2） 峰值流速隨降雨強度增加，呈先增大后逐漸趨于穩(wěn)定的趨勢，見圖5。其中，10 a一遇和20 a一遇降雨頻率下的峰值流速差異較大，而20 a，50 a，100 a一遇降雨的峰值流速差異較小，總體維持在18～20 m/s的水平。

（3） 隨著降雨強度增大，泥石流高速流動區(qū)的距離明顯變長。對比50 a一遇和100 a一遇降雨工況下高速流動區(qū)的溝道長度，發(fā)現(xiàn)100 a一遇降雨條件下的峰值流速略小于50 a一遇降雨條件下的峰值流速，但其高速流動區(qū)更長。

以上3點表明：①? 達(dá)摩溝泥石流峰值流速受限于降雨強度和溝道地質(zhì)條件兩個方面。當(dāng)降雨強度小于20 a一遇時，峰值流速會隨降雨強度增加而增大；當(dāng)降雨強度達(dá)到20 a一遇時，降雨強度繼續(xù)增大，并不會繼續(xù)提高峰值流速，而會延長高速流動的距離。②? 降雨強度增大導(dǎo)致泥石流流速和高速流動距離增加，本質(zhì)上是提高了泥石流體在流動過程中的沖擊速度和沖擊距離，即提高了侵蝕效率，大大提高了泥石流體的破壞能力。這既有利于沿程物源進(jìn)一步補給泥石流體，又進(jìn)一步增加了泥石流體的重度，提高了破壞力。

4.2.2 運動過程泥深變化

泥石流堆積深度是評估泥石流危險性的重要指標(biāo)。圖4，6是不同降雨工況下泥石流堆積情況統(tǒng)計趨勢。總體來看，達(dá)摩溝泥石流堆積特征主要可以概括如下：①? 泥石流最大泥深隨降雨強度增加呈先增大后減小的趨勢；②? 達(dá)摩溝在不同降雨條件下均不會形成明顯的堆積扇；③? 泥石流在溝道狹窄處泥深將會增大，偏向流通區(qū)上游的部分堆積深度較小，而流通區(qū)中下游和堆積區(qū)等泥石流前緣部分堆積深度較大。

造成以上泥石流堆積特征的原因主要有：①? 在降雨強度較低的情況下，水流對堆積物的侵蝕能力和搬運能力均較弱，泥石流流量較小；隨著雨量的增加，產(chǎn)匯流增加，水流的侵蝕搬運能力逐漸增強，泥石流流量變大，導(dǎo)致泥石流的堆積厚度增大；當(dāng)降雨強度增加到一定程度后，泥石流流量增加、流速變大，其裹挾能力更強，更容易將固體物質(zhì)搬運至更遠(yuǎn)的地方，因此高強度降雨條件下溝道內(nèi)的最大堆積厚度反而有所減小。②? 達(dá)摩溝狹窄的溝道不具備形成堆積扇的地形條件。主溝下游地勢相對平坦，有利于泥石流運動減速開始堆積，但歷史泥石流在流通區(qū)下游和近溝口處有堆積，導(dǎo)致目前流通區(qū)末端和近溝口處地勢仍然相對較高，不利于泥石流淤積蠕動形成扇體，而有利于未來泥石流活動進(jìn)一步向清水河流動；此外，在流通區(qū)下游以及堆積區(qū)南側(cè)有部分居民區(qū)并修建有道路，導(dǎo)致達(dá)摩溝主要為帶狀堆積。③? 達(dá)摩溝與清水河交接處溝道比降仍較大，不利于泥石流減速沉積，更可能出現(xiàn)的情況是泥石流體高速沖入清水河，因此，達(dá)摩溝溝口沒有明顯的堆積扇。④? 由于溝道較窄處泥石流的運動受到阻礙，泥石流流速變緩，導(dǎo)致其更易在此處堆積下來。⑤? 在泥石流堆積過程中，隨著泥石流不斷向周圍擴散堆積，其流速、動能不斷降低，泥石流邊運動邊堆積，最終整體停止運動，所以堆積區(qū)中心厚度最大，而發(fā)散區(qū)域堆積厚度逐漸降低。

4.2.3 危險性評價

泥石流的強度可以依據(jù)泥石流的流速與泥深的強度進(jìn)行劃分。本文參照常鳴［32］提出的泥石流強度影響以及危險性分區(qū)標(biāo)準(zhǔn)對模擬結(jié)果進(jìn)行危險性評價。泥石流強度影響劃分標(biāo)準(zhǔn)如表3所列。

根據(jù)表3，利用Arcgis對泥石流流速與泥深進(jìn)行柵格疊加計算，將計算結(jié)果分類為高危險、中危險、低危險，最終得到100 a一遇降雨工況下的達(dá)摩溝泥石流危險性分區(qū)如圖7所示。

在圖7中，高危險性地區(qū)占16.54%，中危險性地區(qū)占49.27%，低危險性地區(qū)占34.19%。低中危險地區(qū)主要分布在主溝邊緣地帶，主溝絕大部分區(qū)域處于中高危險區(qū)。考慮到達(dá)摩村位于主溝溝口處，涉及到大量人、財、物的潛在損失，說明達(dá)摩溝是一條風(fēng)險性較大的泥石流溝，因此需要加強對達(dá)摩溝的預(yù)警監(jiān)測。

5 結(jié)論與建議

（1） 達(dá)摩溝小流域在強降雨條件下有發(fā)生泥石流的可能，模擬結(jié)果顯示：泥石流流速主要分布在0～4 m/s，隨著降雨強度的增大，泥石流平均流速逐漸增大，但泥石流峰值流量隨著降雨強度的增大呈現(xiàn)先增大后減小的趨勢。達(dá)摩溝泥石流堆積泥深主要分布在0～3 m之間，最大泥深可達(dá)3.85 m，平均泥深隨著降雨強度的增大而增大，峰值泥深則呈現(xiàn)先增大后減小的趨勢。

（2） 依據(jù)泥石流流速與泥深，利用Arcgis對達(dá)摩溝泥石流危險性進(jìn)行分區(qū)，結(jié)果顯示：高危險區(qū)占16.54%，中危險區(qū)占49.27%，低危險地區(qū)占34.19%。然而，中高危險區(qū)仍有居民居住，因此風(fēng)險較高，建議加強監(jiān)測。

（3） 泥石流流速及泥深與降雨強度的統(tǒng)計關(guān)系表明：泥石流運動不是簡單受降雨強度影響，還受到地形條件、植被覆蓋等多因素約束，但是不同工況下各因素貢獻(xiàn)有差異。因此，在對泥石流進(jìn)行預(yù)警監(jiān)測時，需針對不同小流域進(jìn)行具體分析，建立不同的監(jiān)測方案。

參考文獻(xiàn)：

［1］ 趙忠海.北京泥石流災(zāi)害預(yù)警研究現(xiàn)狀與技術(shù)方法探討［J］.城市地質(zhì)，2015，10（增1）：78-84.

［2］ 北京市地質(zhì)研究所，中國地質(zhì)科學(xué)院地質(zhì)力學(xué)所，北京林業(yè)大學(xué).北京市溝道侵蝕泥石流多元預(yù)警閾值模型研究及示范工程［R］.北京：北京市地質(zhì)研究所，中國地質(zhì)科學(xué)院地質(zhì)力學(xué)所，北京林業(yè)大學(xué)，2021.

［3］ 洪惜英，王禮先，朱靖才，等.北京市密云縣汗峪溝泥石流［J］.北京林學(xué)院學(xué)報，1983（4）：38-45.

［4］ 李容全，徐振溥.北京山區(qū)泥石流分布特點［J］.山地研究，1983（3）：42-48.

［5］ 張仲德.北京山區(qū)泥石流分布和成因分析［J］.北京師院學(xué)報（自然科學(xué)版），1984（2）：77-89.

［6］ 魏永明.1991年北京北部山區(qū)泥石流成因及防災(zāi)對策［J］.山地研究，1993（3）：187-192.

［7］ 張世清，畢小剛，劉振國.北京山區(qū)泥石流防治對策探討［J］.水土保持通報，1992（3）：46-51.

［8］ 張春山.北京北山地區(qū)泥石流災(zāi)害危險性評價［J］.北京地質(zhì)，1996（2）：11-20.

［9］ 張春山.北京地區(qū)泥石流災(zāi)害危險性評價［J］.地質(zhì)災(zāi)害與環(huán)境保護，1995（3）：33-40.

［10］ 趙波，董桂芝.北京西山清水河流域泥石流及危害預(yù)測［J］.北京地質(zhì)，1996（3）：1-14.

［11］ 謝又予，崔之久.北京山區(qū)泥石流的防治與預(yù)報［J］.貴州科學(xué)，1992（3）：132.

［12］ 吳正華，儲鎖龍.北京泥石流暴雨基本特征［J］.大氣科學(xué)，1992（4）：476-481.

［13］ 韋京蓮，董桂芝，趙波.北京山區(qū)近代泥石流活動規(guī)律及暴發(fā)周期分析［J］.中國地質(zhì)災(zāi)害與防治學(xué)報，1994（4）：48-53.

［14］ 韋京蓮，趙波，董桂芝.北京山區(qū)泥石流降雨特征分析及降雨預(yù)報初探［J］.中國地質(zhì)災(zāi)害與防治學(xué)報，1994（增1）：45-52

［15］ 吳正華.北京泥石流災(zāi)害及其降水觸發(fā)條件［J］.水土保持研究，2001（1）：67-72.

［16］ 白利平.北京山區(qū)泥石流災(zāi)害臨界雨量研究［J］.地質(zhì)災(zāi)害與環(huán)境保護，2006（4）：101-104.

［17］ 白利平，孫佳麗，南赟.北京地區(qū)泥石流災(zāi)害臨界雨量閾值分析［J］.地質(zhì)通報，2008（5）：674-680.

［18］ 王海芝.北京地區(qū)暴雨泥石流預(yù)警閾值研究［J］.第四紀(jì)研究，2020，40（5）：1371-1380.

［19］ 黃曉虎.泥石流降雨啟動機制及早期預(yù)警模型的研究：以達(dá)摩溝為例［D］.長春：吉林大學(xué)，2016.

［20］ 涂劍，馬超，楊海龍.北京山區(qū)暴雨泥石流激發(fā)雨量條件［J］.中國水土保持科學(xué)，2017，15（5）：103-10.

［21］ 方群生，陳志和，唐川，等.基于FLO-2D數(shù)值模擬震區(qū)急陡型泥石流沖出量：以瓦窯溝為例［J］.中山大學(xué)學(xué)報（自然科學(xué)版），2020，59（3）：23-31.

［22］ 趙崢，劉鐵驥，陳亮，等.基于Fluent強震區(qū)王家溝泥石流運動特征研究［J］.災(zāi)害學(xué)，2021，36（2）：208-12.

［23］ 陳明，唐川，王飛龍，等.基于OpenLISEM模型的泥石流啟動沖出預(yù)測［J］.泥沙研究，2020，45（4）：59-65.

［24］ 胡明鑒，汪稔，陳中學(xué)，等.泥石流啟動過程PFC數(shù)值模擬［J］.巖土力學(xué)，2010，31（增1）：394-397，434.

［25］ 慎乃齊，韋京蓮，趙波，等.京西達(dá)摩溝泥石流形成條件與防治［J］.探礦工程（巖土鉆掘工程），2003（增1）：96-98.

［26］ 齊干，張長敏.達(dá)摩溝泥石流形成的物質(zhì)條件分析及防治對策［J］.水文地質(zhì)工程地質(zhì)，2011，38（5）：102-109.

［27］ 喬淵，劉鐵驥，陳亮，等.基于Massflow模型的甘肅省岷縣二馬溝泥石流危險性評價［J］.水利水電技術(shù)，2020，51（4）：184-192.

［28］ 曹水和，吳新明，鐘東.基于Massflow的金沙江滑坡數(shù)值模擬及堵江風(fēng)險預(yù)測［J］.地質(zhì)災(zāi)害與環(huán)境保護，2021，32（4）：3-7.

［29］ OUYANG C，HE S，XU Q，et al.A MacCormack-TVD finite difference method to simulate the mass flow in mountainous terrain with variable computational domain ［J］.Computers & Geosciences，2013，52：1-10.

［30］ HUNGR O，MCDOUGALL S.Two numerical models for landslide dynamic analysis ［J］.Computers & Geosciences，2009，35（5）：978-92.

［31］ 周凱琦.基于單層單相流和多層多相流的災(zāi)害動力學(xué)數(shù)值模擬［D］.武漢：華中科技大學(xué)，2017.

［32］ 常鳴.基于遙感及數(shù)值模擬的強震區(qū)泥石流定量風(fēng)險評價研究［D］.成都：成都理工大學(xué)，2014.

（編輯：高小雲(yún)）

Abstract：

Mudslide disaster events had occurred in the Dharma Gully sub-basin area historically and caused a large casualties and property losses.At present，a large amount of coal gangue was piled up in the main ditch of the sub-basin，and more landslide deposits were accumulated on both slopes sides of the main ditch.In order to analyze the possibility and risk of mudslides reoccurring，this paper analyzed the development of mudslides in Dharma Gully sub-basin area in detail through UAV photogrammetry and field investigation.Based on the improved MacCormack-TVD finite difference method，the Massflow software was used to simulate the mudslide conditions and obtained the mudslide risk in Dharma Gully under 4 rainfall probabilities of 10years，20years，50years，and 100years.The results showed that： ① Under the above 4 rainfall probabilities，the flow velocity of debris flow in most areas was 0～3 m/s，the mud depth was 0～4 m.The maximum flow rate could reach 19 m/s，and the maximum mud depth could reach 385 m；② Under the rainfall of 100-year，the high-risk areas in Damogou accounted for 1654%，the medium-risk areas accounted for 4927%，and the low-risk areas accounted for 3419%；③ The simulation results reproduced the temporal and spatial development characteristics of debris flow well，showed the influence of rainfall intensity on the movement of debris flow，which could provide some references for the early warning and monitoring of debris flow.

Key words：

debris flow；risk；Massflow；numerical simulation；Dharma Gully