基于FLO-2D的西藏若如村泥石流危險性分析

張奮翔， 張路青， 周 劍， 王 頌

(1.中國科學院 地質與地球物理研究所 頁巖氣與地質工程重點實驗室，北京 100029； 2.中國科學院 地球科學研究院， 北京 100029； 3.中國科學院大學， 北京 100049)

1 研究背景

泥石流是山區溝谷或坡面在暴雨、洪水等自然或人為因素作用下發生的一種挾帶大量泥砂、石塊或巨礫等固體物質的特殊洪流，會對人類生命財產或生存環境造成較大危害[1]。青藏高原地區地貌多變，峽谷、河谷、寬谷等交替出現，在喜馬拉雅塊體、拉薩塊體和羌塘塊體的相互作用下，新構造運動活躍[2]。其中，青藏鐵路、拉日鐵路、G318國道、G109國道等交通工程沿線風化劇烈，巖土體完整性較差，發育著大量泥石流，對交通安全及線路運維造成了極大威脅。因此，對此區域泥石流動力過程的研究十分必要。

隨著泥石流動力模型研究的深入以及計算機技術的發展和完善，數值模擬已逐漸成為泥石流動力過程研究的重要手段[3]。利用已有的或新建動力模型，選擇合適的數值方法進行計算分析是泥石流動力過程數值模擬的基礎。在數值算法和本構模型的研究相對成熟后，出現了描述泥石流運動的數值模擬軟件[4]。Debris-2D，FLO-2D和Geoflow等是國際上應用比較成熟的數值模擬軟件[5]，能夠開展較大尺度的工程問題分析。其中，FLO-2D軟件通過建立微分形式的質量守恒和動量守恒方程，采用顯式中心差分法進行求解，可以計算流深、流速和影響范圍隨時間的變化。作為二維洪災模型，FLO-2D的計算時間可以由用戶控制，操作簡單且實用，近年來應用廣泛。Bertolo等[6]運用FLO-2D模擬了小流域泥石流運動淤積過程；阮德修等[7]將FLO-2D與3DMine耦合模擬了尾礦庫潰壩災害下的泥石流流動過程；龔柯等[8]利用FLO-2D對四川省汶川縣綿虒鎮簇頭溝進行了泥石流危險性評價。

西藏地區由于地理位置偏遠、人口居住密度小、人類活動工程少等原因，泥石流災害的研究程度較低，且僅限于拉薩、日喀則、林芝等地區。戚國慶等[9]通過建立降雨型泥石流預測模型，對貢覺縣的3條泥石流溝進行了災害分析；韋方強等[10]進行了西藏古鄉溝泥石流的數值模擬；劉偉朋等[11]研究了拉日鐵路沿線(年木鄉-日喀則段)泥石流的分布特征。

若如村泥石流位于拉薩市堆龍德慶區德慶鄉，青藏公路、青藏鐵路及在建的青藏高速公路位于該泥石流堆積扇的前緣，繞避或下穿泥石流扇。2016年7月2日，德慶鄉暴發降雨頻率為5%的持續性暴雨，當日20時至次日20時累計降雨量為26.6 mm，計算最大雨強為14.09 mm/h，此次降雨使得原若如村泥石流的堆積范圍進一步擴大，對上述工程區域產生極大威脅。在詳細野外調查的基礎上，采用FLO-2D軟件對青藏鐵路、青藏高速公路及若如村的泥石流危險性進行了研究。

2 若如村泥石流概況

若如村泥石流距離拉薩市中心60 km，地處當雄至拉薩峽谷出口段。此地段屬高原寒溫帶濕潤季風氣候，年平均降水量近460 mm；堆龍曲從區內經過，流域面積約5 100 km2；處于雅魯藏布江縫合線與班公錯-怒江縫合線之間，巖性主要是花崗閃長巖和以砂巖為代表的堅硬陸源碎屑巖。G109公路沿線南側邊坡較高，平均坡度大于30°，在風化侵蝕和降雨的綜合作用下，坡面巖體破碎，有大量松散體堆積，因此，泥石流的發育具備良好的動力和物源條件。除了若如村泥石流，區段內還發育有近15條規模較大的泥石流。

若如村泥石流流域面積約3.2 km2，主溝道長2.5 km，最高海拔5 440 m，相對高差1 330 m，兩側山坡的平均坡度為30°，主溝平均坡降30.6%。坡體上部有大量松散堆積體，沖溝內也有碎屑物堆積，為泥石流的發育提供了充足的物源。先前泥石流活動在若如村口形成了一個較大的泥石流堆積扇，堆積扇頂角約40°，扇長350 m，扇寬410 m，面積約0.07 km2。若如村泥石流概況見圖1。

圖1 若如村泥石流概況

在泥石流溝口處，青藏鐵路設置涵洞穿過堆積體，遠離堆積體后鐵路建立高架橋；堆積體前緣和青藏鐵路之間建有三級重力式擋土墻來保障行車安全，墻體高15 m左右，基底寬度為6～7 m；青藏高速公路目前正在修建中，堆積體前方為明洞開挖區(以下簡稱明挖區)，開挖深度近5 m；沿著泥石流溝口至高速路施工區設置主排導槽，槽寬1.5 m，深約0.5 m。該泥石流的主要威脅對象為青藏鐵路、青藏高速公路和若如村，此外，主排導槽在泥石流流動過程中可能面臨淤積失效的風險。

3 數值模型構建

3.1 高精度DEM網格劃分

根據FLO-2D泥石流運動數值模擬的方法和理論[8]建立模型。為了進行本次研究，在若如村泥石流堆積區前緣飛行無人機獲得了泥石流流域的航拍照片。將航拍照片的高程和GPS信息提取后，獲得流域的DEM數據進而轉化為FLO-2D軟件能夠識別的ASCII文件。主要步驟可分為：

(1)獲取DEM。將航拍照片利用Agisoft PhotoScan Professional軟件進行處理，得到流域的DEM數據。

(2)ASCII文件轉化。運用Uedit64中工具箱的轉換功能進行文件格式的轉換，得到FLO-2D建模所需的DTM數據。

(3)在FLO-2D中新建項目，將帶有高程數據的DTM文件加載進去作為基礎地形并進行網格劃分。選取的網格尺寸為2 m，網格數量146 972個。計算網格如圖2所示，劃分出計算區域，最后對選定的計算區域進行高程插值。模型的入水口設置在泥石流溝上游的物源區，出水邊界設置為堆龍曲河道。

圖2 FLO-2D計算網格

3.2 模型參數取值

FLO-2D數值模擬選取的主要參數有[12]：降雨強度i，等效曼寧系數n，泥石流土石材料比重Gs和體積濃度CV，層流阻滯系數K，屈服應力τ及黏滯系數η，其中τ、η用相關系數α1、β1、α2、β2表征。

η=α1eβ1CV

(1)

τ=α2eβ2CV

(2)

根據現場情況，將泥石流物源分為坡體上部和沖溝內兩部分。實地調查得到坡體上部松散體平均厚度約0.25 m，沖溝內碎屑物的平均厚度約0.15 m，結合衛星圖進一步得出：坡體上部松散體面積約1.7×105m2，沖溝內碎屑物面積約3×105m2。日降雨量按照實際工況26.6 mm計算，該泥石流流域內的集水量為8.5×104m3。按照以上數據得到泥石流土石材料的體積濃度為50.6%，取0.51。根據 FLO-2D 手冊結合現場情況，若如村泥石流土石材料的比重取值為2。

由土石材料重度和體積濃度，劃分為中阻泥石流，則泥砂比Rns=0.75，現場調查若如村泥石流的堆積扇厚度約5 m。根據王裕宜等[13]提出的不同泥砂比條件下的統一公式：

(3)

τ=0.021Rns-0.3exp(19.64Rns0.25·CV)

(4)

η=1.39×10-4Rns-0.47exp(18.07Rns0.13·CV)

(5)

式中：h為泥深，m。將h=5，Rns=0.75代入公式(3)，得到流動區域曼寧系數n=0.09；將公式(1)、(2)、(4)、(5)聯立，得：

α1=1.39×10-4Rns-0.47,β1=18.07Rns0.13

(6)

α2=0.021Rns-0.3,β2=19.64Rns0.25

(7)

將Rns=0.75代入公式(6)、(7)，得到α1=0.000159，β1=17.41，α2=0.0229，β2=18.28。

層流阻滯系數K的取值范圍為 24～50 000，鑒于土石材料的體積濃度值處于0.48～0.55，且具有黏性泥石流的特點，參考FLO-2D手冊中對于層流阻滯系數K的建議，取K=2 285。模擬參數取值見表1。

表1 FLO-2D輸入參數值

為全面分析流域內的危險性，以德慶鄉1998-2018年的降雨數據為基礎，分別對降雨頻率為5%、2%和1%下的若如村泥石流進行模擬計算，并以上述頻率設計泥石流流量過程線。

根據水文地質手冊，對降雨數據分析計算得到2%和1%頻率下的日降雨強度分別為48.1mm和57.8mm，采用雨洪法計算出該泥石流在不同降雨頻率下集水點的流量，結果見表 2。

由于 FLO-2D數值模擬時無法考慮泥石流對溝道的侵蝕效應，因此在清水流量的基礎上乘以流量放大系數，得到泥石流流量再生成流量過程線輸入到FLO-2D中：

(8)

式中：BF為流量放大系數。

表2 若如村泥石流不同降雨頻率下的集水點流量

根據表2數據，采用典型洪水過程放大法求取泥石流流量過程線[14]如圖3。

4 模擬結果驗證

由于實際24 h降雨量達到20年一遇的級別，將5%降雨頻率下的模擬結果與實際測量結果進行對比，對模擬結果的準確性進行評價。采用堆積范

圍、溝口最大沖出距離與橫向最大堆積寬度3個指標，將現場調查與遙感影像結合得到的實際堆積數據與模擬情況進行對比，見表3。根據表3數據，計算泥石流數值模擬精度因子A，公式如下：

(9)

式中：S1為實際堆積面積，m2；S2為模擬堆積面積，m2；S0為模擬與實際堆積范圍重疊的面積，m2。

圖3 不同降雨頻率下泥石流流量概化過程線

表3 模擬堆積與實際堆積數據對比

將表3數據代入公式(9)進行計算，得到此次模擬的精度因子為68.7%。通過野外調查發現，高速公路明挖坑周邊及隧道口、青藏鐵路路基下方都有殘留的泥石流松散體，若如村西北角房屋前方仍可見0.5 m厚的土石堆積體，與模擬結果基本相符。由于溝口最大沖出距離的模擬誤差率為19.4%，橫向最大堆積寬度的模擬誤差率為8.4%，與精度因子結合來看，此次模擬結果比較符合實際。

5 模擬結果分析

該泥石流的危險性分析主要包括兩方面：(1)主排導槽在不同降雨強度下的導流情況分析；(2)堆積區下穿交通工程和若如村在不同降雨頻率下的危險性分析。主排導槽從泥石流溝口向下延伸，經過鐵路涵洞上方，將沖溝上游的洪水引導至高速公路隧道明挖區。主排導槽深度約0.5 m，當其邊界的泥石流流深大于0.5 m時，流體會從導槽內溢出，此時主排導槽失效。下面將結合不同日降雨強度下泥石流流深、流動速度以及沖擊力，對主要威脅區域的危險性進行分析。

分別選取各威脅區域內最易產生堆積的網格記錄不同降雨工況下泥石流的流動情況，各威脅區域代表性點的具體位置如表4所示，取點位置如圖4。

表4 各威脅區域危險性分析的選點網格單元及位置

5.1 泥石流的流動深度分析

圖5給出了3種降雨頻率下模擬得到的泥石流最大流深。從堆積范圍來看，5%降雨頻率下的堆積范圍約15.8×104m2，2%頻率下約17×104m2，1%頻率下為18×104m2。在5%、2%和1% 3種降雨頻率下，青藏鐵路涵洞西出口及其以西鐵路段、鐵路涵洞、高速公路隧道明挖工程區以及若如村的大部分區域都被泥石流所覆蓋。鐵路涵洞東出口、高速路東側工程區以及河漫灘在5%降雨頻率下未見泥石流覆蓋，在2%降雨頻率下開始有流體堆積；1%頻率下的泥石流堆積范圍更廣，鐵路涵洞東出口完全被覆蓋。

圖6給出了上述各代表性點在3種不同降雨頻率下泥石流流深隨時間的變化曲線。結合圖5、6分析如下：

(1)在5%、2%和1%降雨頻率下，鐵路涵洞西出口的流深分別為1.4、1.9、2.1 m；與涵洞西出口相比，涵洞東出口受泥石流影響較小，在5%和2%降雨頻率下幾乎無流體堆積，在1%降雨頻率下的流深小于0.1 m。根據鐵路行車規范，當流體覆蓋道床時鐵路將無法通行，由此判斷5%降雨頻率下的泥石流流深就會影響鐵路正常通車。

(2)主排導槽左邊界在5%、2%和1%降雨頻率下的流深分別為2.39、2.74和2.91 m，導槽深度僅有0.5 m，由此判斷在5%降雨頻率下主排導槽就會失效。從流深曲線來看，主排導槽在最初10 min里的流深小于0.5 m，因此，當泥石流流量較小時，主排導槽可以較好地發揮排導作用。

(3)由于地勢較低，高速公路明挖點成為各個方向的匯水區，其流深變化曲線的規律性不強。5%降雨頻率下的流深約5.4 m，而明挖坑深5 m，處于溢滿狀態，對施工人員安全和工程穩定有顯著影響；2%和1%降雨頻率下的流深超過6 m，對明挖區域的安全和穩定影響更大。

(4)在5%降雨頻率下，若如村西北角流深的最大值為0.67 m，村民的居住和通行將受到嚴重影響；2%和1%降雨頻率下的流深超過1 m，居民安全將面臨威脅。

(5)隨著時間的變化，各代表性點的流深都經歷了無堆積-突然增大-逐漸減小的過程，符合泥石流突然暴發和逐漸消退的特點。

5.2 泥石流流動速度與沖擊力分析

圖7給出了不同降雨頻率下各代表性點模擬過程中的泥石流流速變化，由圖7可以得出：

(1)5%降雨頻率下，泥石流在鐵路涵洞西出口的最大流速為0.52 m/s；在鐵路涵洞東出口的流速為0，這與流深的變化相對應；主排導槽左邊界的最大流速達到1.35 m/s，為工程區域內最大流速點，再次驗證了主排導槽的導流功能；與流深變化曲線類似，高速路明挖區的流速變化呈無規律曲線，可以得出其最大流速為0.65 m/s；若如村西北角房屋邊界的最大流速為0.1 m/s，威脅相對較小。

(2)與流深變化相同，各代表性點的模擬流速基本呈現不流動-突然增大-逐漸減小的過程；與流深變化不同的是，各點最大流速并非隨著日降雨強度的提高而增大，如主排導槽左邊界的最大流速在1%降雨頻率下為1.44 m/s，小于2%降雨頻率下的1.59 m/s，而高速路明挖區在5%降雨頻率下的最大流速也大于2%降雨頻率，說明流深增大產生的淤積效應可能引起泥石流流速的降低。

圖8為不同降雨頻率下各代表性點模擬過程中泥石流沖擊力變化。根據圖8可以得出，在5%降雨頻率下，鐵路涵洞西出口的最大泥石流沖擊力約1.1 kPa，東出口沖擊力為0；主排導槽左邊界在泥石流運移過程中的最大沖擊力為16.29 kPa；高速路明挖點的最大沖擊力為7.51 kPa；若如村西北角房屋承受的最大沖擊力僅25.06 Pa。

5.3 其他極端降雨條件下的泥石流危險性分析

為了更加全面地分析各威脅區域在極端降雨天氣下的泥石流危險性，以研究區域近20年的年最大24 h降雨量數據為基礎，選取其中的最小值19.5 mm作為極端降雨區間的下限，將1%降雨頻率下的24 h降雨量57.8 mm作為上限，對降雨量每隔5 mm進行插值，進行了19.5、24.5、29.5、34.5、39.5、44.5、49.5、54.5、57.8 mm共9種日降雨強度的模擬。將模擬結果進行處理，得到各威脅區域的最大流深和最終流深隨降雨量的變化曲線，如圖9所示。

由圖9可得，當日降雨強度在19.5～57.8 mm變化時：

(1)泥石流在鐵路涵洞西出口的最大瞬時流深大于1 m，在日降雨強度為54.5 mm時可達2.2 m，根據鐵路行車規范，列車在某個時間段內無法安全通行；然而，從最終流深變化曲線分析，涵洞西出口的泥石流流深會逐漸消退至0.4 m以下，對列車通行安全影響有所減小，19.5～34.5 mm日降雨強度下的最終流深在0.2m左右，對列車通行幾乎無影響。鐵路涵洞東出口的泥石流最大流深幾乎為零，對列車通行無影響。

(2)主排導槽左邊界的最大瞬時流深大于2m，在日降雨強度為54.5 mm時接近3 m；不同日降雨強度下主排導槽左邊界的最終流深在0.7～1.1 m之間變化，流深始終大于主排導槽深度，故判斷為失效。

圖6 3種降雨頻率下各代表性點的泥石流流深變化曲線

圖7 3種降雨頻率下各代表性點的泥石流流動速度變化曲線

(3)當日降雨強度大于19.5 mm時，明挖坑的瞬時流深大于5 m，此時坑內已被泥石流淤滿，施工安全及工程穩定將會受到嚴重影響；從最終流深來看，坑內泥石流流深會逐漸消退至2～4 m，可能延誤工期。

(4)日降雨強度在19.5～29.5 mm變化時，若如村內的最大瞬時流深小于0.8 m，可能影響居民出行；當日降雨強度大于29.5 mm時，村內的最大瞬時流深大于1m，可能危及居民的生命財產安全。從最終流深變化曲線分析，若如村的泥石流流深會逐漸消退至0.4 m以下，此時仍可能影響居民交通，對生命財產安全幾乎無影響。

(5)隨著日降雨強度逐漸增大，各威脅區域的最大瞬時流深均呈現先增大后減小的趨勢，在日降雨強度為54.5 mm時達到最大值，說明此降雨強度下的若如村泥石流危險性最強。

圖8 3種降雨頻率下各代表性點的泥石流沖擊力變化曲線

圖9 各威脅區域流深隨日降雨量的變化曲線

6 人類活動區的泥石流危險性評價

基于上述模擬獲得的泥石流流深與流速，對人類活動區整體的危險性進行評價。為了使評價結果更具代表性，將1%、2%和5% 3種降雨頻率進行綜合分析。根據FLO-2D手冊，泥石流的危險性與災害程度及暴發頻率相關，評價標準如圖10。泥石流的災害程度受流體的流深和流動速度控制，采用汶川泥石流強度劃分標準[8]，如表5所示。

圖10 泥石流危險性評價標準

各顏色表示意義如下：透明色(無標注色)為無危險區，表征安全區域；黃色為低危險區，表示該區域內，泥石流可導致工程建筑物輕微損傷；橙色為中危險區，表示該區域內工程建筑物部分摧毀，可能出現人員傷亡；紅色為高危險區，表示泥石流會導致該區域內工程建筑摧毀、人員傷亡慘重。從透明色到紅色依次分級為0～3。根據上述標準，得到若如村泥石流危險性評價圖見圖11。

表5 泥石流災害程度評價標準

圖11 若如村泥石流危險性評價圖(人類活動區)

從圖11來看，高危險區主要分布在主排導槽、涵洞西出口以及高速路明挖區，此區域內包括居住、鐵路通行以及高速公路施工在內的一系列人類活動將面臨極大的安全威脅，這與前文分析一致；鐵路涵洞、高速公路隧道以北的工程區以及若如村的小部分區域處在中危險區，可能出現工程設施摧毀和人員傷亡；涵洞東出口及其以東的青藏鐵路段、高速公路隧道以東的工程區基本無危險。

7 結 論

(1)在實際降雨頻率5%的情況下，泥石流堆積面積約1.58×105m2，最大流深為5.7 m，最大流速13.1 m/s，人類活動區域會出現鐵路無法正常通車、主排導槽失效和高速路明挖區溢滿等問題，若如村村民的住行安全也將受到影響。

(2)當日降雨強度在19.5～57.8 mm變化時：青藏鐵路若如村段在降雨期間無法通車，雨后列車通行不受影響；主排導槽始終失效；高速公路隧道明挖坑在降雨期間達到淤滿狀態，雨后淤積深度2～4 m，可能延誤工期；日降雨強度小于29.5 mm時，若如村民出行不便，大于29.5 mm時，泥石流會危及居民的生命財產安全；54.5 mm日降雨強度下的泥石流危險性最高，在防治工作中應高度重視。

(3)根據危險性分區評價圖，高危險區主要分布在主排導槽、涵洞西出口以及高速路明挖區，鐵路涵洞、高速路隧道以北工程區以及若如村的小部分區域處在中危險區，涵洞東出口及其以東鐵路段、高速公路隧道以東工程區基本無危險。

(4)在中、高危險區，需要加強泥石流防治工程建設，建議在鐵路涵洞西出口或若如村西北角布置降雨觀察點，根據雨情采取相應的預防措施。此外，建議增加主排導槽坡降以增大流體流速，考慮采用速流結構。