陳 偉

(廣東省環(huán)境地質(zhì)勘查院，廣東 廣州 510000)

0 引言

泥石流災(zāi)害屬于山區(qū)的一種突發(fā)性災(zāi)害，具有流速快、攜帶固體物質(zhì)能量強和破壞力大等特點，是我國主要地質(zhì)災(zāi)害類型之一(謝洪等，2006)。開展泥石流危險范圍的預(yù)測工作具有重要的現(xiàn)實意義。

利用新一代的地理信息系統(tǒng)(GIS)，結(jié)合適用于泥石流兩相流體分析的理論模型，綜合地質(zhì)、數(shù)值計算以及計算機科學等不同學科，對泥石流運動以及堆積過程進行模擬，進而為后期泥石流危險性評價區(qū)劃及防治提供關(guān)鍵參數(shù)，是當前泥石流理論研究以及工程防治設(shè)計的研究熱點和難點。從目前國內(nèi)外泥石流運動堆積過程研究現(xiàn)狀出發(fā)，探討Voellmy連續(xù)介質(zhì)模型原理以及GIS在泥石流研究中的應(yīng)用，并以典型泥石流溝為例，通過模擬泥石流運動、堆積的動力過程，獲取了泥石流體的最大速度、壓力、過流量、堆積高度以及堆積面積等關(guān)鍵參數(shù)。

根據(jù)對國內(nèi)外文獻的分析，目前對泥石流堆積扇運動及堆積預(yù)測通常采用3種方法：經(jīng)驗統(tǒng)計方法、物理模擬方法、動力模型分析。

日本是較早也是較多研究泥石流運動及堆積預(yù)測的國家之一。池谷浩等1979年根據(jù)流域面積推算泥石流沖出量，根據(jù)沖出量推算泥石流堆積長度和堆積寬度，率先從統(tǒng)計學角度探討了這一問題(池谷浩，1996);石川芳治(1999)研究了用模型實驗?zāi)M堆積區(qū)具有溝槽和隆起等微地貌時，泥石流堆積范圍應(yīng)如何修正以及泥石流體中細粒物質(zhì)對泥石流堆積擴散的影響等。

歐美國家也很重視泥石流危險范圍的研究。奧地利很早就進行了泥石流危險范圍的預(yù)測工作，并引用交通信號中的紅、黃、藍3色的特定含義，將泥石流危險區(qū)分為3類(Blaikie et al，1994);Brien等開發(fā)了水力模型，模擬泥石流和流域內(nèi)的堆積區(qū);Ghilardi等(2001)開發(fā)了1個數(shù)學模型，考慮了侵蝕和堆積過程，以及不同類型的堆積物;Pasuto等(2004)采用多學科方法，從歷史、地形地貌、地質(zhì)結(jié)構(gòu)、氣象、土壤和森林管理等角度，對泥石流進行了評價。

唐川等(1992)應(yīng)用泥石流堆積模型實驗，將實驗結(jié)果用多元回歸曲線擬合，得到泥石流堆積形態(tài)函數(shù)與泥石流參數(shù)關(guān)系表達式;吳積善等(1993)根據(jù)泥石流典型地區(qū)(主要是東川小江流域)的資料，建立了泥石流最大堆積面積、最大堆積長度、最大堆積寬度和泥石流堆積幅角與流域面積、松散固體物質(zhì)儲量、主溝長度和流域最大相對高差的相關(guān)式。

利用動力模型進行泥石流運動及堆積預(yù)測是目前研究的熱點，該模型通常采用數(shù)值方法，利用能量與運動轉(zhuǎn)換法則對泥石流運移及堆積過程進行模擬。動力模型包括離散元模型、質(zhì)量集中模型以及連續(xù)單元模型。其中，連續(xù)介質(zhì)模型是建立在流體力學之上，利用物質(zhì)－運動－能量消散的轉(zhuǎn)換方程來描述泥石流的動力過程。因此，連續(xù)介質(zhì)模型可以更為合理地刻畫出泥石流固體物質(zhì)的運動特性，從而被廣泛使用(Hungr et al，1984)。

1 連續(xù)介質(zhì)模型在泥石流中的應(yīng)用

在連續(xù)介質(zhì)模型中，泥石流流體被假設(shè)為一種非穩(wěn)定以及非均質(zhì)的流體，可以被2個主要流體參數(shù)來刻畫(Hungr et al，2009)，即流體高度H(x，y，t)，m;平均流速U(x，y，t)，m/s：

式(1)中，Ux為X方向的速度;Uy為Y方向的速度;T為平均速度的轉(zhuǎn)置矩陣。速度大小可以被定義為：

式(2)中，‖U‖表示對速度U取絕對平均值，這就確保U在向量空間中為嚴格的正速度，流體速度的方向可以用單位向量(nU)來定義：

Voellmy流變模型用下列質(zhì)量平衡公式來表達：

式(4)中，H為流體高度，m;Q x，y，()t，kg/(m2·s)為質(zhì)量源，當Q=0時，表示沒有物質(zhì)沉積。在X和Y方向上，流體平均深度平衡方程可以分別表示為：

以及：

式(5)、(6)中，Cx以及Cy為剖面系數(shù)，gz為垂直方向的重力加速度。在Voellmy模型中，垂直方向的接觸關(guān)系可以定義為非均質(zhì)的摩爾－庫倫關(guān)系。ka/p為土壓力系數(shù)，可以用下列關(guān)系式來表示：

式(7)中，φ為泥石流流體的內(nèi)摩擦角。

通過以上各式，可以得出Voellmy流變公式：

2 實例研究

2.1 泥石流溝基本特征

選取地震災(zāi)區(qū)汶川縣某典型泥石流溝作為實例研究，該溝位于岷江右岸，流域面積5.8 km2，主溝長約5.2 km，以“V”型谷為主，流域海拔范圍1 320～3 340 m之間，高差 1 020 m，主溝縱比降為200‰。該溝流域內(nèi)出露基巖主要以花崗巖、板巖及千枚巖為主。溝道兩岸多基巖出露，坡度較陡，平均坡度約40°～60°。該溝主要威脅到下游的阿壩師專、村莊、公路及農(nóng)田等承災(zāi)體。通過對該溝流域的調(diào)查特征分析，可將該溝可分為清水、物源、流通及堆積4個區(qū)，各分區(qū)特征見表1。

表1 泥石流流域分區(qū)特征表Table 1 Zonation features of the debris flow basins

通過實地調(diào)查，該溝于1958年及1979年8月13日分別暴發(fā)中—大規(guī)模泥石流。其中1979年發(fā)生的泥石流持續(xù)時間約lh，速度較快，“龍頭”高度為2.0～3.1 m，泥漿飛濺，泥位高度約3.1 m(康景文等，2011)。在“5·12”汶川特大地震的強大作用力之下，溝內(nèi)及兩岸松散固體物質(zhì)激增，在強降雨或綿雨條件下，很容易暴發(fā)嚴重的泥石流災(zāi)害。

2.2 泥石流特征值計算

根據(jù)相關(guān)計算公式，分別計算出該溝全流域的暴雨歷時、設(shè)計暴雨量、徑流深度、設(shè)計洪水總量、概化矩形歷時等參數(shù)，該溝的洪水總量、泥石流峰值流量、1次泥石流總量、1次泥石流沖出的固體物質(zhì)總量等計算結(jié)果(表2、表3、表4、表5)。

表2 暴雨洪水流量計算結(jié)果Table 2 Calculation of rainstorm and flood flow

表3 泥石流峰值流量計算結(jié)果Table 3 Calculation of peak flow of debris flow

表4 1次泥石流總量計算結(jié)果Table 4 Calculation of the total amount of one-time debris flow

表5 1次泥石流固體物質(zhì)總量計算表Table 5 Calculation of the total solid material of one-time debris flow

根據(jù)計算結(jié)果，該溝20年、50年及100年一遇的泥石流總量約4 230～9 052 m3。

2.3 泥石流運動及堆積數(shù)值模型建立

數(shù)字地面模型是泥石流數(shù)值模擬最基礎(chǔ)，也是最重要的步驟，這直接決定了模擬結(jié)果的精度以及正確性。通過對工作區(qū)紙質(zhì)地圖進行矢量化后，利用ArcGIS軟件，建立TIN(Triangulated Irregular Network)最優(yōu)化不規(guī)則三角網(wǎng)模型，通過對TIN模型數(shù)字高程點的采集，獲取模型網(wǎng)格點數(shù)字高程，并進行內(nèi)插計算，最終獲得研究區(qū)DEM。通過GIS模擬系統(tǒng)內(nèi)置分析代碼，獲取流體的流速、高度、過流量等關(guān)鍵參數(shù)，并將計算所得的結(jié)果以GIS格式數(shù)據(jù)的形式以3D效果進行展示。

在上述連續(xù)介質(zhì)模型原理的基礎(chǔ)上，采用動力數(shù)值模擬計算軟件，對該泥石流溝進行實例研究。模型長度=7 300 m，寬度=5 040 m，共計50 034個節(jié)點，49 580個單元。對該溝百年一遇的泥石流影響范圍進行預(yù)測，1次泥石流固體物質(zhì)總量為9 052 m3，建立的三維模型如圖1所示。

圖1 泥石流數(shù)值模擬三維模型Fig.1 3D numerical simulation model of debris flow

3 結(jié)果與討論

采用常規(guī)的泥石流流量過程概化線和洪峰流量進行模擬，針對泥石流的堆積范圍進行模擬，并為后期的危險性評價及區(qū)劃提供技術(shù)支撐。通過模擬計算，泥石流最大速度、最大壓力、過流量以及流體堆積高度、堆積特征如圖2、圖3、圖4、圖5、表6所示。

圖2 泥石流最大速度特征圖Fig.2 Characteristic diagram of maximum velocity of debris flow

圖3 泥石流最大壓力特征圖Fig.3 Characteristic diagram of maximum pressure of debris flow

圖4 泥石流最大過流量特征圖Fig.4 Characteristic diagram of maximum momentum of debris flow

圖5 泥石流最大流體堆積高度特征圖Fig.5 Characteristic diagram of maximum accumulation height of debris flow

表6 泥石流運動及堆積計算表Table 6 Calculation of the movement and accumulation of debris flow

從以上分析結(jié)果圖可以看出，該溝泥石流暴發(fā)后，在泥石流的流通區(qū)和物源區(qū)，泥石流流體的最大流速可以達到6.96 m/s，快速的泥石流流體可以產(chǎn)生很強的側(cè)蝕和底蝕作用，加大了泥石流流體的固體物質(zhì)攜帶能力，并大幅度地增加了固體物質(zhì)運移距離，表現(xiàn)出強烈的侵蝕作用。在溝谷的中下段，特別是出山口至主河之間的溝段，泥石流流速已經(jīng)有大幅度的降低，以沉積作用為主。根據(jù)數(shù)值模擬計算結(jié)果，在1次泥石流固體物質(zhì)總量為9 052 m3的情況下，泥石流流體攜帶固體物質(zhì)沖出溝口后，堆積扇面積可高達0.147 km2，將對溝口大量居民、民房、公路等造成嚴重危害。

另外，根據(jù)模擬，泥石流的堆積扇扇頂剛好沖入岷江，對岷江主河道造成擠迫，雖未完全阻斷，但會造成河流的回水，對位于溝口所處的岷江上游區(qū)段造成潛在威脅。

4 結(jié)論

(1)泥石流危險范圍預(yù)測是泥石流危險性評價及區(qū)劃的一項重要研究內(nèi)容，特別是針對地震災(zāi)區(qū)泥石流災(zāi)害防治有著重要的意義。

(2)在預(yù)測泥石流危險范圍時，動力模型分析方法有著經(jīng)驗統(tǒng)計方法以及物理模擬等方法不具備的優(yōu)勢，可以更為合理地刻畫出泥石流固體物質(zhì)的運動特性，從而被廣泛使用。

(3)選取地震災(zāi)區(qū)汶川縣某典型泥石流溝作為實例研究，該溝流域面積5.8 km2，主溝長約5.2 km，該溝20年、50年及100年一遇的泥石流總量約為4 230～9 052 m3。

(4)采用動力數(shù)值模擬計算軟件，建立泥石流溝的數(shù)值模型，對該溝百年一遇的泥石流影響范圍進行預(yù)測，計算結(jié)果顯示該溝在1次泥石流固體物質(zhì)總量為9 052 m3的情況下，最大速度為6.96 m/s，最大壓力99.96 kPa，最大過流量為5.55 m2/s，最大流體堆積高度為3.1 m，堆積面積為0.147 km2。不僅對下游居民造成嚴重威脅，還有可能阻塞岷江主河道。

(5)根據(jù)野外調(diào)查情況，并結(jié)合數(shù)值模擬研究，證明基于GIS技術(shù)平臺以及Voellmy連續(xù)介質(zhì)模型可以很好地模擬泥石流運動、堆積的動力過程，在泥石流危險范圍預(yù)測工作中具有較高的準確性。

池谷浩.1996.云仙水無川泥砂流失量的推算方法［J］.水土保持科技情報，(4)：30－34.

康景文，黃建波，周勇，等.2011.汶川縣威州鎮(zhèn)羊嶺溝泥石流特征及成因分析［J］.施工技術(shù)，40(增刊)：180－184.

唐川，劉希林.1992.泥石流動力堆積模擬和危險范圍預(yù)測模型［M］.北京：科學出版社.

謝洪，鐘敦倫，韋方強，等.2006.我國山區(qū)城鎮(zhèn)泥石流災(zāi)害及其成因［J］.山地學報，24(1)：79－87.

吳積善，田連權(quán)，康志誠，等.1993.泥石流及其綜合治理［M］.北京：科學出版社.

石川芳治.1999.地震にょる土石流の発生に係わる地形，地質(zhì)條件［J］.砂防學會誌，51(5)：35－42.

BLAIKIE P，CANNON T.1994.At Risk：Natural Hazards，Peoples Vulnerability and Disasters［M］.London，UK：Psychology Press.

GHILARDI P，NATALE L，SAVI F.2001.Modelling debris－flow propagation and deposition［J］.Physics and Chemistry of the Earth：Part C，26(9)：651 －656.

HUNGR O，MORGAN G C，KELLERHALS R.1984.Quantitative analysis of debris torrent hazards for design of remedial measures［J］.Canadian Geotechnical Journal，21(4)：663－677.

HUNGR O，MCDOUGALL S.2009.Two numerical models for landslide dynamic analysis［J］.Computers ＆ Geosciences，35(5)：978－992.

JULIEN P Y，O'BRIEN J S.1997.Slected notes on debris flow dynamics［J］.Lecture Notes in Earth Sciences，64：144 －162.

PASUTO A，SOLDATI M.2004.An integrated approach for hazard assessment and mitigation of debris lows in the Italian Dolomites［J］.Geomorphology，61(1/2)：59 － 70.