許曉君，胡卸文,2，許澤鵬，韓雨菲

(1.西南交通大學(xué)地球科學(xué)與環(huán)境工程學(xué)院，四川 成都 610031；2.西南交通大學(xué)高速鐵路運(yùn)營安全空間信息技術(shù)國家地方聯(lián)合工程實驗室，四川 成都 610031)

0 引言

我國山區(qū)地質(zhì)災(zāi)害頻發(fā)，近年來，泥石流災(zāi)害成為制約我國山區(qū)經(jīng)濟(jì)發(fā)展的重要因素之一。泥石流造成的堵河常常也造成巨大災(zāi)害，如1985年西藏培龍溝泥石流堵斷帕隆藏布致使上游水位迅速升高約20 m，沖毀房屋；2010年舟曲泥石流堵斷白龍江形成堰塞湖，淹沒半個縣城[1]。目前國內(nèi)眾多學(xué)者對泥石流堵河影響因素、堵河判據(jù)進(jìn)行了研究，并取得了一定成果。

韋方強(qiáng)等[2]認(rèn)為主支流交匯關(guān)系是導(dǎo)致堵江的主要原因；陳春光等[3]研究表明入?yún)R角和主支流速度比影響著泥沙堆積位置和堆積形態(tài)；吳積善等[4]認(rèn)為泥石流性質(zhì)、泥石流規(guī)模、泥石流顆粒粒度、主河來水、河床形態(tài)、主支流交匯角以及交匯處地形等是形成堵塞壩的重要因素。周必凡等[5]提出了依據(jù)一次泥石流規(guī)模、主河流量、扇形地泥石流溝床條件等判斷是否堵塞主河的計算公式；唐川等[6]對岷江上游和云南小江多條堵河型泥石流溝進(jìn)行分析總結(jié)，推導(dǎo)出泥石流堵河判斷經(jīng)驗公式；張金山等[7]對岷江上游典型堵河型泥石流基本特征進(jìn)行分析總結(jié)，得出計算堵斷主河可能性的經(jīng)驗公式；何易平[8]利用非線性回歸方法對短時泥石流堵塞主河的條件進(jìn)行研究，得到了泥石流堵河的動量和流量判別式指標(biāo)；許向?qū)幍萚9]通過野外調(diào)查泥石流堵河特征，依據(jù)泥石流沖出距離、扇頂至主河距離以及主河寬度，提出了泥石流堵河判斷公式；王志剛等[10]運(yùn)用層次分析法選取泥石流堵河評價因子，構(gòu)造判斷矩陣得出權(quán)重，得到了泥石流堵河判斷標(biāo)準(zhǔn)。

圖1 2014年泥石流堆積范圍Fig.1 Accumulation extent of debris flow in 2014

位于四川省阿壩州壤塘縣宗科鄉(xiāng)臥龍村2組的石吾隴溝為典型“窄-陡型”泥石流溝，其溝口緊鄰主河宗科河。每年雨季石吾隴溝均有不同規(guī)模松散物沖出淤積于主河河床中，而2014年爆發(fā)的一次較大規(guī)模泥石流堵塞2/3的主河道，致使主河水位上漲，嚴(yán)重威脅溝口對岸10戶50余居民生命財產(chǎn)安全(圖1)。因此，本文在詳細(xì)調(diào)查石吾隴溝泥石流發(fā)育及堆積特征的基礎(chǔ)上，采用流體動力學(xué)計算軟件CFX對其堆積過程及范圍進(jìn)行模擬，主要分析泥石流堆積物堵塞主河的動力過程，對山區(qū)泥石流防治提供理論依據(jù)。

1 石吾隴溝基本特征

1.1 石吾隴溝發(fā)育特征

石吾隴溝地形陡峻，屬典型構(gòu)造侵蝕高山地貌。溝流域整體呈葉脈狀，流域面積約3.97 km2，流域最高海拔4 410 m，最低點位于溝口與宗科河交匯處，海拔3 545 m，相對高差865 m，溝域內(nèi)發(fā)育有一條主溝及一條較大支溝，主溝長2.64 km，支溝長2.30 km，主溝平均縱坡降248.9‰，支溝平均縱坡降236.8‰，泥石流流域地質(zhì)平面圖及溝道縱斷面圖見圖2、圖3。

圖2 石吾隴溝流域地質(zhì)平面圖Fig.2 The geologic plane of Shiwulong Gully

清水區(qū)集水面積2.18 km2，約占總流域面積55%，此段縱坡較陡，平均縱坡降373.5‰。此區(qū)域內(nèi)坡面上局部基巖出露明顯，溝谷形態(tài)不明顯、切割弱，未見有崩塌、滑坡等不良地質(zhì)現(xiàn)象發(fā)生，岸坡穩(wěn)定性較好，為泥石流提供了重要水利條件。

圖3 石吾隴溝主溝縱斷面圖Fig.3 The longitudinal cross-section of main gully

形成流通區(qū)流域面積約1.73 km2，約占總流域面積的43%。主溝溝道較寬，縱比降238.5‰，呈寬“U”型谷地貌，谷寬一般20～30 m，下切程度較弱，主流不穩(wěn)定，在高水位時可發(fā)生偏移；支溝溝道狹窄，縱比降216.7‰，呈“V”型溝谷地貌，其溝床拉槽下切及兩岸側(cè)蝕嚴(yán)重，下切深度一般8～10 m，側(cè)蝕寬度一般3～5 m，參與泥石流活動的物源絕大部分分布于此區(qū)段內(nèi)。

堆積區(qū)流域面積約0.05 km2，分布于溝道出口與主河道交匯處，此區(qū)段溝道平均縱坡降約97.9‰。泥石流堆積物主要為塊碎石土，塊碎石呈次棱角狀，粒徑30～60 cm，個別粒徑可達(dá)100 cm，分選性差。

1.2 物源特征

1.2.1溝道物源

形成流通區(qū)溝道整體縱坡比降較大，下蝕及側(cè)蝕嚴(yán)重，尤其是支溝。溝床堆積的松散固體物源以坡殘積含塊碎石土粉質(zhì)黏土為主。溝道物源總量5.153×105m3，可能參與泥石流活動量6.71×104m3，典型溝道物源照片及剖面見圖4、圖5。

圖4 溝道物源特征Fig.4 The characteristics of channel source

圖5 溝道典型物源橫斷面圖Fig.5 The cross-section of representative gully resource

在降雨及泥石流的沖刷作用下，溝道兩側(cè)岸坡以下蝕和側(cè)蝕的方式補(bǔ)給泥石流，下切深度一般8～10 m，側(cè)蝕寬度一般3～5 m。

下蝕產(chǎn)生的松散堆積物堆積于溝道內(nèi)成為泥石流的物源，其主要過程包括四個階段：下切侵蝕→形成溝壑→崩滑堵塞→沖刷潰決，動儲量按下切侵蝕型泥石流動儲量圖解法(圖6(a))依據(jù)以下公式進(jìn)行計算[11]：

式中：V01——下切侵蝕型動儲量/m3；

L1——溝床堆積體長度/m；

h——原溝床深度/m；

θ——自然休止角/(°)。

圖6 溝道物源啟動模型(喬建平)Fig.6 Statistic model of downward and lateral erosion induced material

側(cè)蝕造成兩側(cè)岸坡沿潛在塌滑面發(fā)生塌滑，塌滑后的松散物堆積于溝道內(nèi)，形成物源參與泥石流活動，其過程包括兩個階段：側(cè)緣侵蝕→失穩(wěn)滑坍，其動儲量按(圖6(b))以下公式計算[11]：

式中：V02——下切侵蝕型動儲量/m3；

L2——溝道堆積體長度/m；

l——實測坡面長度/m。

1.2.2滑坡堆積物源

滑坡堆積物源總量1 600 m3，可能參與泥石流活動量300 m3，典型溝道物源照片及剖面見圖7、圖8。

圖7 滑坡物源特征Fig.7 The characteristics of landsilde source

圖8 典型滑坡物源橫斷面圖Fig.8 The cross-section of representative landslide resource

HP01位于斜坡體上，位于泥石流過流斷面之上，其參與泥石流活動的方式為在暴雨沖刷下，以坡面侵蝕或拉槽下切的方式補(bǔ)給泥石流活動，因此其可能參與泥石流活動的比例相對較小。

1.3 爆發(fā)特征

根據(jù)調(diào)查，石吾隴溝溝口老泥石流堆積扇發(fā)育，表明該溝歷史上曾發(fā)生過大規(guī)模的泥石流。而新近發(fā)生的泥石流則下切老泥石流堆積扇，并在老泥石流堆積扇邊緣形成新的堆積(圖9)。

圖9 石吾隴溝泥石流堵塞Fig.9 The debris flow deposits blocked the valley of Shiwulong

2014年8月，石吾隴溝流域經(jīng)歷連續(xù)強(qiáng)降雨后爆發(fā)泥石流，泥石流沖出的物質(zhì)堆積在主河宗科河內(nèi)，呈扇形分布，軸向長約50 m，扇緣寬約180 m，沖出泥石流固體物質(zhì)約有6 000 m3，堵塞了約2/3的河道。致使主河宗科河河水位上漲至河邊的鄉(xiāng)道，河水漫過鄉(xiāng)道。

2 泥石流運(yùn)動特征參數(shù)

2.1 容重

通過訪問居民確定泥石流漿體性狀，測量體積與質(zhì)量，計算泥石流容重，經(jīng)計算得泥石流漿體容重為1.65 t/m3。

2.2 流速

因泥石流是固液兩相混合流，其運(yùn)動過程十分復(fù)雜，因此目前提出的公式均是經(jīng)驗性的。石吾隴溝泥石流流體性質(zhì)偏于稀性，因此其流速采用西南地區(qū)(鐵二院陳光曦)稀性泥石流計算公式[11]：

式中：VC——泥石流平均流速/(m·s-1)；

γH——泥石流固體物質(zhì)重度/(t·m-3)；

φ——泥石流泥沙修正系數(shù)；

R——泥石流水力半徑(m)，天然溝道可用平均泥深代替；

Ic——泥石流水力坡度，天然溝道可用溝床縱坡代替。

經(jīng)計算，泥石流溝口流速為4.46 m/s(表1)。

表1 泥石流流速計算表

2.3 泥石流流量

(1)暴雨洪峰流量

采用《四川省中小流域暴雨洪水計算手冊》推薦的推理公式由設(shè)計暴雨計算洪峰流量。

式中：QP——暴雨洪水設(shè)計流量/(m3·s-1)；

ψ——洪峰徑流系數(shù)；

A——流域面積/km2；

S——單位歷時的暴雨平均強(qiáng)度/(mm·h-1)；

t——匯流時間/h；

n——暴雨遞減指數(shù)。

計算參數(shù)由設(shè)計泥石流流量計算手冊中相關(guān)圖表查得，流域特征值從1∶50 000數(shù)字高程模型獲取，暴雨設(shè)計洪峰流量結(jié)果見表2。

表2 泥石流流量計算參數(shù)及結(jié)果

(2)泥石流峰值流量

1.雨洪法

依據(jù)雨洪法計算泥石流峰值流量：

QC=(1+φ)QPDC

式中：QC——泥石流斷面峰值流量/(m3·s-1)；

φ——泥石流泥沙修正系數(shù)；

QP——暴雨洪峰流量；

DC——堵塞系數(shù)。

計算參數(shù)及結(jié)果見表2。

2.泥痕法

在溝口上游200 m(海拔3 582 m)處，發(fā)現(xiàn)明顯的泥痕，取該斷面計算泥石流流量，計算參數(shù)及結(jié)果見表3。

表3 泥石流斷面流量計算

根據(jù)調(diào)查及計算，2014年泥石流的流量為35.38 m3/s，屬接近20年一遇的泥石流流量。

3 堵河程度分析

泥石流堵河分析已經(jīng)成為地質(zhì)災(zāi)害領(lǐng)域研究的熱點問題，目前很多學(xué)者對泥石流堵河判別研究做了大量工作，也取得了較多成果。目前，對泥石流堵河可能性判別多是針對某一流域或某條泥石流溝的堵河可能性判別，具有一定的適用條件。

唐川等[6]、張金山等[7]、何易平[8]分別提出了以下三種計算堵斷主河可能性的經(jīng)驗公式：

式中：Z——泥石流堵河判別指標(biāo)；

K——修正系數(shù)；

QC——泥石流流量；

γC——泥石流重度/(t·m-3)；

θ——泥石流溝與主河順流方向的夾角；

QB——主河流量；

B——堆積區(qū)主河寬度；

R——泥石流堵河判別指標(biāo)；

P——泥石流爆發(fā)頻率；

JC——泥石流溝床坡降；

JB——主河比降；

CF——泥石流堵河判別指標(biāo)；

Qc——泥石流單寬流量；

Qb——主河單寬流量；

γB——主河流體重度/(t·m-3)。

R大于10，易堵斷主河；R介于5～10，會對主河造成嚴(yán)重堵塞；R小于5，會造成主河輕微堵塞。Z≥1.0，則泥石流可能堵斷主河；Z=1.0～0.5，則泥石流造成部分堵河；Z≤0.5，則泥石流不會造成堵塞。CF<-8.572時會造成堵河。

依據(jù)經(jīng)驗公式，對不同頻率泥石流的堵河判別指標(biāo)進(jìn)行計算，計算參數(shù)見表4，堵河指標(biāo)計算結(jié)果見表5。

表4 泥石流堵河計算參數(shù)

表5 泥石流堵河計算結(jié)果

計算結(jié)果顯示，發(fā)生20年一遇的泥石流，依據(jù)張金山模型判斷，會發(fā)生嚴(yán)重堵河，而依據(jù)唐川模型和何易平模型，判斷結(jié)果分別為堵斷和不堵。根據(jù)野外調(diào)查，2014年發(fā)生接近20年一遇的泥石流，沖出的固體物質(zhì)侵占三分之二主河道，造成了較嚴(yán)重的堵河，因此張金山模型預(yù)測結(jié)果與石吾隴溝實際情況較為符合?；诖?，本文選取張金山模型預(yù)測50年一遇泥石流堵河情況。經(jīng)計算分析得出，50年一遇的泥石流會造成主河堵斷。

4 泥石流堵河數(shù)值模擬

4.1 CFX軟件介紹

CFX軟件是目前應(yīng)用最廣的流體動力學(xué)軟件之一，是利用計算機(jī)與數(shù)值模擬的方法，通過求解流體力學(xué)方程，對流體力學(xué)問題進(jìn)行模擬分析。CFX功能的實現(xiàn)，主要包括以下幾個步驟：模型建立、網(wǎng)格劃分、前處理、求解、后處理。

4.2 模型建立及網(wǎng)格劃分

泥石流山體模型為復(fù)雜的三維實體模型，本次建模采用CAD、Surfer、ANSYS等軟件進(jìn)行建模。本文以石吾隴溝溝口堆積區(qū)為對象建立模型，模型平面范圍300 m×200 m，流體域高度20 m，采用四面體結(jié)構(gòu)進(jìn)行網(wǎng)格劃分，共生成27 100個節(jié)點，134 167個單元，網(wǎng)格模型見圖10。

圖10 模型及網(wǎng)格劃分Fig.10 Size and meshes of model

4.3 前處理

本次石吾隴溝泥石流模型邊界條件設(shè)置上游入口為進(jìn)口(inlet)，設(shè)置下游及兩側(cè)為出口(outlet)，泥石流溝道地形面為壁面(wall)，頂部為開放式(opening)邊界，邊界條件見圖11。

圖11 邊界條件示意圖Fig.11 Schematic diagram of boundary condition

在入口邊界處設(shè)置泥石流流量以質(zhì)量流量形式流入，對20年一遇泥石流持續(xù)900 s、50年一遇泥石流持續(xù)1 200 s的堆積范圍進(jìn)行預(yù)測。在t=0 s時刻，取洪水流量，在泥石流持續(xù)1/3過程時刻，取泥石流峰值流量，泥石流流量隨時間變化如圖12～圖13所示。

圖12 20年一遇泥石流流量時間關(guān)系曲線Fig.12 Relation between inlet flow and time (P=5%)

圖13 50年一遇泥石流流量時間關(guān)系曲線Fig.13 Relation between inlet flow and time (P=2%)

4.4 求解

本次對石吾隴溝泥石流模擬，20年一遇總時長為900 s，50年一遇總時長1 200 s，步長為1 s，每隔5 s保存一次計算文件。求解過程中監(jiān)測入口和出口的質(zhì)量流量，質(zhì)量流量曲線見圖14～圖15。

圖14 20年一遇進(jìn)出口質(zhì)量流量曲線Fig.14 Curves of bulk mass flow rate of inlet and outlet (P=5%)

圖15 50年一遇泥石流進(jìn)出口質(zhì)量流量曲線Fig.15 Curves of bulk mass flow rate of inlet and outlet (P=2%)

由圖12、圖13可知，20年一遇泥石流進(jìn)口質(zhì)量流量在300 s達(dá)到最大值，50年一遇泥石流在400 s達(dá)到最大，與前處理階段設(shè)置相符；20年一遇泥石流出口流量在360 s達(dá)到最大值，50年一遇在450 s達(dá)到最大，說明在此時泥石流峰值流量流過出口邊界。

4.5 后處理及分析

后處理通過對數(shù)據(jù)結(jié)果進(jìn)行圖形化顯示和統(tǒng)計。生成等值面，顯示云圖和矢量圖。本文通過顯示速度云圖對不同頻率泥石流的運(yùn)動過程進(jìn)行分析。

4.5.120年一遇泥石流

20年一遇泥石流速度云圖見圖16～圖20。

圖16 t=20 s速度云圖(P=5%)Fig.16 Contour of velocity at t=20 s (P=5%)

圖17 t=50 s速度云圖(P=5%)Fig.17 Contour of velocity at t=50 s (P=5%)

圖18 t=300 s速度云圖(P=5%)Fig.18 Contour of velocity at t=300 s (P=5%)

圖19 t=360 s速度云圖(P=5%)Fig.19 Contour of velocity at t=360 s (P=5%)

圖20 t=900 s速度云圖(P=5%)Fig.20 Contour of velocity at t=900 s (P=5%)

對20年一遇泥石流速度云圖進(jìn)行分析，在0～50 s時段內(nèi)，隨著入口流量增大，泥石流流體沿溝道向下運(yùn)動，隨地形變得開闊，流體也不斷往兩側(cè)擴(kuò)展，在t=50 s時，流體達(dá)到兩側(cè)出口邊界；50～300 s時段內(nèi)，入口流量繼續(xù)增大，流體流經(jīng)范圍仍不斷增大，在t=300 s時，峰值流量通過入口；300～360 s時段內(nèi)，隨著峰值流量的流入，流體泛濫范圍仍繼續(xù)擴(kuò)大，在t=360 s時，流體泛濫范圍達(dá)到最大，并且此刻峰值流量流過出口邊界；在360～900 s時段內(nèi)，隨著入口流量不斷減小，并且流體從出口邊界流出，使得流體泛濫范圍也不斷減小。

20年一遇泥石流在整個運(yùn)動程中，泥石流從入口流入，流經(jīng)溝道，在溝口處堆積，堆積范圍擴(kuò)大至溝口宗科河內(nèi)，侵占宗科河河道約3/4，造成河道嚴(yán)重堵塞。與實際發(fā)生的20年一遇泥石流沖出堵塞主河道范圍基本一致。

4.5.250年一遇泥石流

50年一遇泥石流速度云圖見圖21～圖25。

圖21 t=20 s速度云圖(P=2%)Fig.21 Contour of velocity at t=25 s (P=2%)

圖22 t=45 s速度云圖(P=2%)Fig.22 Contour of velocity at t=45 s (P=2%)

圖23 t=400 s速度云圖(P=2%)Fig.23 Contour of velocity at t=400 s (P=2%)

圖24 t=600 s速度云圖(P=2%)Fig.24 Contour of velocity at t=600 s (P=2%)

圖25 t=1 200 s速度云圖(P=2%)Fig.25 Contour of velocity at t=1 200 s (P=2%)

對50年一遇泥石流速度云圖進(jìn)行分析，其運(yùn)動過程與20年一遇泥石流運(yùn)動過程類似。與20年一遇泥石流相比，50年一遇泥石流規(guī)模更大、持續(xù)時間更長、速度更快。泥石流在0 s時刻從入口流入，由流動速度更快，故在t=45 s時刻前緣“龍頭”便運(yùn)動至主河道對岸；在45～400 s時段內(nèi)，隨入口流量增大，流體流經(jīng)范圍不斷增大，由于50年一遇泥石流持續(xù)時間更長，故在t=400 s時，達(dá)到峰值流量，峰值流量通過入口，沿溝道向下運(yùn)移，并沿兩側(cè)擴(kuò)散，在t=600 s時，流體泛濫范圍達(dá)到最大；在600～1 200 s時段內(nèi)，隨著入口流量不斷減小，并且流體從出口邊界流出，使得流體泛濫范圍也不斷減小。

50年一遇泥石流在整個運(yùn)動過程中，泥石流從入口流入，流經(jīng)主河道堆積淤埋，侵占全部宗科河河道，造成宗科河堵斷。

5 結(jié)論

(1)石吾隴溝溝域內(nèi)物源豐富，物源動儲量為6.71×104m3，溝域內(nèi)地形陡峭，平均縱坡降為248.9‰，所在地區(qū)雨量充沛，泥石流發(fā)育條件充分。

(2)2014年石吾隴溝爆發(fā)了較大規(guī)模泥石流，泥石流峰值流量35.38 m3/s，相當(dāng)于20年一遇降雨頻率下爆發(fā)的泥石流規(guī)模，沖出的固體物質(zhì)堆積于溝口宗科河內(nèi)，侵占三分之二主河道，造成較為嚴(yán)重的堵河。

(3)對已提出的泥石流堵河經(jīng)驗公式進(jìn)行綜合分析，張金山模型較適用于石吾隴溝，用其對石吾隴溝不同降雨頻率的泥石流堵河程度進(jìn)行預(yù)測，結(jié)果顯示20年一遇泥石流會造成主河嚴(yán)重堵塞，而50年一遇泥石流會造成主河堵斷。

(4)利用CFX軟件對不同頻率泥石流堆積過程進(jìn)行模擬，結(jié)果顯示20年一遇泥石流隨時間增長，泥石堆積濫范圍不斷增大，在360 s時，堆積范圍達(dá)到最大，泥石流蔓延至宗科河內(nèi)，侵占主河道，造成宗科河嚴(yán)重堵塞；50年一遇泥石流相較于20年一遇泥石流持續(xù)時間更長、堆積范圍更廣，在600 s時堆積范圍達(dá)到最大，蔓延至宗科河內(nèi)，造成宗科河堵斷。