基于熵權變異系數融合算法和FLO-2D的朱家溝流域泥石流危險性評價

吳季寰，張春山，楊為民，孟華君，郭 涵，萬飛鵬

（1.中國地質科學院地質力學研究所，北京 100081；2.自然資源部活動構造與地質安全重點實驗室，北京 100081；3.自然資源部陜西寶雞地質災害野外科學觀測研究站，北京100081；4.中國地質大學（北京），北京 100083）

引言

青藏高原東北緣地處西秦嶺構造帶西延［1］，位于南北向地震帶與祁連山地震帶的交匯區，活動斷裂發育，地震活動性強，山地災害頻發。坐落在此的山地城鎮，如甘肅岷縣、武都等地區，多屬泥石流高發地［2，3］，泥石流災害已嚴重制約區域經濟社會發展。

甘肅岷縣朱家溝是青藏高原東北緣典型的小流域泥石流溝，歷史上多次發生泥石流災害，如1982年“8·8”朱家溝泥石流堵斷洮河中上游，形成長達15 km的回水區，河水越過防洪堤侵入岷縣城區，淹沒建筑物1 700余棟，毀壞房屋600余間，造成100余人傷亡，并導致下游約148.67 hm2農田絕收減產。近年來，受2008年汶川地震和2013年“7·22”岷縣-漳縣Ms6.6地震影響，朱家溝內物源劇增，導致該小流域泥石流危險性急劇升高，極有可能再次爆發大規模泥石流災害事件，極大威脅流域內及城區居民生命財產安全，因此亟需開展相應的泥石流災害危險性評估，為防災減災提供科學依據。

泥石流危險評價評價方法大致分為多因子疊加法［4-8］和數值模擬辦法［9-12］兩類，二者區別在于：1）評價方式不同：多因子疊加法以泥石流影響因素為數據基礎，對流域進行綜合評價分區，直接體現流域整體危險度，數值模擬法基于對流域內所爆發泥石流的泥深與流速的模擬，著重于對泥石流動態過程的反演，并以此為依據進行危險性分區；2）適用范圍不同：多因子疊加法適用于大范圍區域，而數值模擬辦法更適用于小流域或具體溝谷；3）評價目的不同：多因子疊加法旨在評價各流域于未來階段暴發泥石流的可能性，數值模擬法注重將泥石流強度與重現周期相結合［11］，模擬不同頻率下泥石流的運動過程，預測泥石流的運動特征及堆積區域。

本文結合以上兩類方法，通過詳細現場調查和資料分析，采用熵權與變異系數融合的多因子疊加法開展朱家溝泥石流易發性評價；在此基礎上，進行強降雨條件（P=1.6%，P=1%）下朱家溝及其高易發支溝的泥石流運動數值模擬，并與野外科考資料對比驗證，最后以流速、泥深、堆積范圍為基礎進行危險性區劃，該研究可為當地泥石流災害防治提供依據。

1 流域泥石流概況

朱家溝坐落于甘肅岷縣城區北部，是洮河一級支溝，流域面積約14 km2，主溝長度約6.7 km，與縣城一河之隔。流域位于青藏高原東北緣，東昆侖斷裂帶和西秦嶺北緣斷裂帶的應變傳遞和構造轉換的中間過渡區，屬中強震多發區，基本地震動峰值加速度為0.15 g，溝內哈樹—麻路灘斷裂為高角度逆沖兼右旋走滑斷裂，傾向北東，傾角約55～70°，斷層上盤以二疊系炭板巖夾石英砂巖、粉砂巖為主，表層中風化，風化殼厚度2～3 m，屬中等硬度塊狀碎屑巖組，下盤以三疊系板巖、砂質板巖、粉砂巖為主，表層強風化，風化殼厚度5～8 m，屬軟弱層狀碎裂巖。強表層巖土體受地震擾動，結構疏松，崩塌滑坡的觸發閾值降低，在地震觸發效應下，崩滑發育，溝道中上游淤積嚴重，具備豐富的泥石流固體物源條件。

朱家溝地勢北高南低，屬西秦嶺北部構造侵蝕剝蝕中低山地貌，海拔介于2 276～3 017 m，相對高差741 m，主溝縱比降約72～110‰，溝道上窄下寬，自北到南曲形延伸，溝道斷面自上而下由“V”型向“U”型過渡，兩岸分布沖溝33條，支溝17條，其中郎家溝、松樹溝、朱林溝和哈固溝是朱家溝的主要支溝。流域地處溫帶半濕潤與高寒濕潤氣候過渡帶，多年平均降水量560.8 mm，降水集中于5～9月，占全年降水總量的78%以上，且多為夏季短時強降雨，受其激發泥石流活動頻繁，嚴重威脅城區及流域內居民生命財產安全。

綜上，朱家溝一直是岷縣實施泥石流防災減災的重要地區。在岷縣特大冰雹山洪泥石流災后重建和岷縣-彰縣地震災后重建時，分別在朱家溝、松樹支溝、郎家支溝、王家支溝溝口修建排導槽，在松樹溝修建攔擋壩3道以及在王家支溝修建攔擋壩1道。

2 基于熵權與變異系數融合算法的泥石流易發性評價

2.1 融合算法原理

信息熵是用于度量事物的不確定性［13］，在地質災害領域表達為影響因子攜帶信息熵的大小與災害發生的概率成反比。通過信息熵計算熵權的過程如下：

（1）建立標準化評價矩陣：

式中m、n、x ij分別代表評價對象數、評價因子數以及第i個評價對象的第j個因子的標準值。

（2）確定第j個因子的信息熵權：

式中e j為第j個因子的信息熵，w1j為第j個因子的熵權。

變異系數是一個表達概率分布離散程度的歸一化量度，在評價體系中，某因子的變異系數愈大，則該因子愈能代表體系內部的差異性，對評價結果影響愈大。第j個因子變異系數權重如下：

式中Cv j為第j個因子的變異系數，w2j為第j個因子的變異系數權。

相對熵（式4），可用于衡量兩組離散向量之間差異。

根據相對熵的定義，引入融合向量W=(w1,w2,…,w n)與待評價向量組Wp×n，使該向量與熵權法權重向量以及變異系數法權重向量距離之和最小［13］。建立距離規劃模型以求解融合權重：

式中p、n的物理意義分別為災害評價模型數以及災害評價因子數，w i為第j個因子的融合權重。

2.2 泥石流易發程度的影響因子

從泥石流形成條件：降水、溝谷、物源三方面出發，選取影響泥石流發育的主要因素：形狀系數、流域面積、溝谷形心至斷層距離、溝谷密度、溝道長度、主溝縱比降、溝口單位寬度清水流量峰值、物源分布、相對高差、流域平均坡度、植被覆蓋率，淺表層工程巖組共計12個評價因子，各因子實際取值見表1。

表1 研究區各因子實際取值Table 1 The actual value of each factor in the study area

由于研究區范圍較小，同期次降雨下各支溝雨量接近，故降雨量不作為評價因子。但因一次性短歷時、強降雨形成泥石流過程中，流域內雨水短期入滲、蒸發較少，多以地表徑流形式激發物源啟動，并最終從溝口流出，所以溝口處單位寬度清水流量峰值可以直觀體現該流域產匯流能力。流域附近歷史最大日降雨量為61.5 mm（2010年5月10日），據岷縣氣象局2020年6～8月逐日降雨量資料，7月18日（以下簡稱7·18）59.3 mm的降雨量為2020年度當地雨季最大降水量。本文基于FLO-2D的Rain模塊，降雨分布如圖2，選擇適用于西部山區小流域的SKK模型［14］，進行朱家溝“7·18”強降雨-徑流模擬，以提取各支溝溝口清水流量過程峰值，見表1。

圖2 朱家溝“7.18”降雨曲線圖Fig.2 Zhujia gully'7.18'rainfall curve

依據行業規范［15］及相關研究成果［6.16］，將各因子分為：強（IV）、中（III）、弱（II）及不發育（I）四級，并將分級結果代入融合算法計算權重。分級標準與權重計算結果見表2，分級結果見圖3。

圖3 因子易發程度量化分級Fig.3 Quantitative classification of factors’susceptibility

表2 評級因子發育等級及權重Table 2 Development grade and weight of rating factors

2.3 泥石流易發程度評價結果

經多因子疊加計算，朱家溝內微流域泥石流易發程度評價結果介于［0.424，0.624］，利用統計學中的自然間斷點法進行分類（圖4）：高易發泥石流溝（0.526，0.624］、中易發泥石流溝（0.464，0.526］、低易發泥石流溝［0.424，0.464］。結合2020年8月野外實地調查分析（見圖1），流域內主要支溝：哈固溝、朱林溝、松樹溝和郎家溝的均為泥石流高易發溝。

圖1 朱家溝流域地質簡圖Fig.1 Geological schematic map of Zhujiagully watershed

圖4 朱家溝泥石流易發程度評價結果圖Fig.4 Assessment results of debris flow susceptibility in Zhujia gully

3 基于FLO-2D的朱家溝泥石流危險性評價

綜合考慮朱家溝內各支溝泥石流易發程度和獲取的ALSO12.5 m數字高程數據（DEM）的精度，選取哈固溝、朱林溝、松樹溝、郎家溝等高易發泥石流支溝以及主溝開展評價。

3.1 工況設置及計算參數

FLO?2D以非牛頓流體和中央有限差分法為基礎，迭代計算流體運動控制方程，模擬洪水、泥石流運動過程［12］。計算條件設置如下：

（1）降雨條件：根據國家氣象局提供的1981-2010年岷縣逐日降水數據，利用指數外推法對岷縣地區年降雨量頻率曲線（圖5）進行擬合（式6），可知雨量59.3 mm的2020年度朱家溝流域最大單日降雨：“7·18”強降雨，年降雨頻率為1.6%，屬50年一遇暴雨事件，設置此次降雨作為模擬條件。此外，根據式6選取單日雨量99.38 mm作為百年一遇極端降水條件進行模擬。

圖5 岷縣地區年降雨量頻率擬合曲線Fig.5 Fitting curve of annual frequency of rainfall in Minxian County

式中：R為降雨量，P為年降雨頻率，擬合公式置信概率為95%。

（2）有無防治工程條件：歷經“5·10”特大冰雹和“7·22”岷縣-漳縣地震后，朱家溝泥石流危險性大幅增加，且緊鄰縣城風險較大，定西國土局為有效防治地質災害，在模擬區內修建工程措施見表3。

表3 朱家溝泥石流治理工程參數Table 3 Parameters of debris flow control project in Zhujia gully

組合降雨、防治工程條件，模擬3種工況下高易發溝道泥石流的危險性，分別為：①無防治工程時“7·18”強降雨（50年一遇）泥石流，②有防治工程時“7·18”強降雨（50年一遇）泥石流以及③有防治工程時百年一遇強降雨泥石流爆發。以工況①、②作為對照組評判工程措施效用，以工況②、③作為對照組分析極端氣象下泥石流運動趨勢，以極端工況③作為目標組進行危險性區劃。

計算所需參數設置如下：1）重度及體積濃度：據李國營等在朱家溝的工作成果，主溝、郎家溝、松樹溝、朱林溝和哈固溝的泥石流容重γ（kN/m3）分別為［17.46，17.59，17.53，16.16，17.16］，流域平均容重為17.2 kN/m3，平均比重Gs=2.6 g/cm3，經計算流域內平均孔隙比為e=1.228，體積濃度C V=0.463，對應泥砂比Rns=0.70；2）層流阻滯系數：本文參考侯圣山等［10］研究，用工程地質類比法，取層流阻滯系數K=2 280；3）屈服應力及粘滯系數：根據1）中參數，聯立粘滯系數η（Pas）和屈服應力τy（kpa）計算公式以及王裕宜等［17］統一的泥砂比—體積濃度—流變參數關系式（式7），可得［α1，α2，β1，β2］為［0.012 2，0.015 6，17.41，18.28］，以及［η，τy］為［30.146，56.992］。

式中：η為粘滯系數，τy為屈服應力，CV為體積濃度，Rn s為泥砂比，α1、α2、β1、β2為實驗系數。

（3）曼寧粗糙系數：利用王裕宜等［17］統一阻力糙率系數表征公式，結合實際調查主溝、郎家溝、松樹溝、朱林溝和哈固溝的泥深h（m）為［1.7，1.3，1.3，1.5，1.2］，可得各溝曼寧系數nc為［0.29，0.14，0.14，0.22，0.10］；5）集水點及泥石流流量過程曲線：集水點取自各溝的泥石流形成區底部流通區頂部，根據“7·18”降雨-徑流模擬提取的集水點處清水流量線和放大系數B F=1/（1-CV）=1.862可得點位處泥石流流量過程曲線（圖6）。

圖6 溝谷集水點分布及泥石流流量過程曲線Fig.6 The distribution of catchment points of gullies and flow process curves of debris flow

3.2 模擬結果

圖7在實際降雨條件下朱家溝泥石流運動過程的模擬結果。由圖7可知：（1）“7·18”強降雨時（圖7（a）、（b）、（d）、（e）），無工程條件下最大泥速6.03 m/s，最大泥深4.88～8.07 m，位于溝口堆積扇中部以及溝道匯集處；有工程條件下最大泥速4.54 m/s，三道攔擋壩前最大泥深為7.42～8.58 m，郎家溝至主溝排導槽段受后方補給與前方運動受阻影響，淤積較深可達4.2～5.75 m，溝口堆積區最大深度為3.23 m，攔擋壩與排導槽起到攔截回淤、減速、導流控界的作用。（2）百年一遇強降雨時（圖7（c）、（f）），水動力條件增強，雨水潤濕溝道，朱家溝內泥速提升、泥深增加，泥石流物質沖出溝口，堵塞洮河，直接威脅岷縣縣城。

圖7 多工況下泥深、泥速分布Fig.7 Distribution of mud depth and mud velocity under different conditions

經驗證（表4），模擬精度高于74.12%，能反映“7·18”降雨條件下區內泥石流運動特征與堆積范圍，且說明百年一遇降雨條件下泥石流預測較為可信。

表4 模擬結果驗證表Table 4 Validation table of simulation results

3.3 泥石流危險性評價

實地調查表明，朱家溝泥石流直進性明顯，破壞形式以沖毀、淤埋農田、民居、公路和擠壓河道，泥速、泥深與危險性呈正相關關系，參考侯圣山等［10］在岷縣耳陽河流域的研究（表5），在百年一遇強降雨條件和2020年“7·18”強降雨條件下，對朱家溝泥石流進行危險性分級（圖8）。

表5 泥石流危險性分區標準Table 5 Criteria for hazard zoning of debris flow

泥石流危險性評價結果表明：“7·18”泥石流總危險區面積為371 875 m2，高危險區占危險區總面積的45.72%，主溝中下游均處于高危險區，溝口排導槽將高危險性區限制于槽中，抑制堆積扇成形，此外松樹溝內三道攔擋壩處由于起到攔擋回淤作用表現為高危險性；中危險區占危險區總面積的30.70%，分布在松樹溝、朱林溝以及主溝上游區域、溝口排導槽出口周圍；低危險區占危險區總面積的23.58%，分布于哈固溝和郎家溝。P=1%概率下，泥石流總危險區面積為500 625 m2，與“7·18”降雨條件相比，危險區面積增大34.62%。其中高危險區占總面積71.62%，增加了188 525 m2，具體表現為主溝溝頭處危險性上升，朱家溝主溝沿岸房屋建筑易受泥石流威脅，以及溝口沖積扇成型，主要威脅洮河左岸溝口處與右岸30 m范圍內城區建筑；中危險區面積約占18.34%，包括朱林溝、哈固溝和松樹溝，但由于松樹溝內攔擋工程淤滿，超過極限狀態，使其周邊危險性較高；低危險區面積約占10.04%，郎家溝整體危險性較低，溝口排導設施仍可持續倒流。

總體而言，朱家溝泥石流危險性大。其每年雨季向洮河瀉流，泥石流物質在溝河交匯處淤積，抬高河床，在極端的氣象條件下（圖8（b））有阻斷洮河，形成堰塞壩的風險，進而會迫使上游水位增漲，淹沒岷縣城區。并且朱家溝溝口下游2.5 km處為龍王臺水電站，溝內流體瀉入洮河后，存在部分流體順流向下匯入水電站庫區，擠占庫容的可能性。

圖8 朱家溝泥石流危險性分區圖Fig.8 Zhujia gully debris flow risk zoning map

4 結論與討論

本文在青藏高原東北部山地泥石流災害調查的基礎上，選取具有區域代表性的泥石流孕災小流域——朱家溝為研究區，用于探索適合青藏高原東北緣小流域泥石流溝的危險性評價方法。所得方法將依托熵權-變異系數融合算法的多因子疊加模型和數值模擬技術相結合，前者通過衡量泥石流形成的地形、地質條件評價小流域內各溝道泥石流發育程度，判別在降雨誘發下，區內易于暴發泥石流并匯入主溝，參與主溝泥石流流體物質補給，最終共同致災的支溝；后者則結合前者結果，引入泥石流觸發因素，通過模擬極端降雨條件下流域內泥石流的爆發過程和運動特征，進行泥石流危險性區劃，二者相互補充，充分考慮青藏高原東北緣小流域的泥石流孕災條件，用于泥石流危險性評價具有充分的合理性。經2020年朱家溝“7·18”泥石流實例驗證，該法具有較高的準確度，用其進行朱家溝泥石流危險性評價得到以下結論和認識：

（1）朱家溝內各支溝泥石流易發程度主要受溝口單寬清水流量峰值、主溝長度、淺表層巖組等條件控制，區內17條支溝中高易發溝6條，中易發溝5條，低易發溝6條，其中哈固溝、朱林溝、松樹溝和郎家溝等主要支溝均為高易發程度泥石流溝。

（2）泥石流流速和溝谷淤積數值模擬結果表明，在50年一遇降雨強度以下，朱家溝溝口的排導槽以及松樹溝內三座攔擋壩對溝內泥石流起到了導流控界和攔截減速的作用，但基本已達到其能力極限。在百年一遇的極端降雨條件下，對溝口居民和房屋、岷縣縣城與下游龍王臺水電站威脅較大，建議在主溝修建攔擋壩等防治工程，并及時清淤。

（3）“7·18”朱家溝泥石流危險區集中分布于主溝、哈固溝、朱林溝、松樹溝和郎家溝溝道，高危險區分布于主溝中下游區域。當遭遇百年一遇強降雨條件（P=1%），朱家溝泥石流危險區面積約500 625 m2，較“7·18”泥石流危險性區面積增加34%，其中高危險區分布于朱家溝主溝、溝口處與洮河右岸30 m范圍，對岷縣縣城構成重大威脅。評價結果可為岷縣縣城的地質災害防治、國土空間規劃提供支撐。

小流域泥石流危險性區劃的影響因素眾多，從泥石流高易發支溝的選取到模擬時所選定的極端降雨條件均直接影響評價結果。本文雖建立起包含了12個地形地質方面致災因子的泥石流易發程度評價體系，但忽略了因子相關性，未能將高度相關因子剔除或整合，一定程度上影響了高易發支溝的選擇，同時，僅重現了百年一遇降雨條件下的泥石流過程，難以全面反映爆發不同規模泥石流時流域的危險性分布。因此，為得到更準確的評價結果，在下一步的工作中應重點提煉泥石流致災因子，凝練評價體系，以及開展不同雨強條件下小流域泥石流危險性綜合區劃。

致謝

本文在撰寫過程中使用了岷縣氣象局提供的朱家溝流域降雨資料以及岷縣岷陽鎮朱家溝泥石流治理工程勘察資料，在此表示誠摯的感謝。