三明魯坑儲木場沖溝泥石流地質(zhì)災(zāi)害特征及發(fā)生過程模擬

簡 磊

(三明市三元區(qū)自然資源局， 三明，365000)

泥石流是一種災(zāi)源遠(yuǎn)、歷時短、強度高、破壞性強的地質(zhì)災(zāi)害[1]，對其開展預(yù)警和防治，目前較好的方法是通過危險性評價和預(yù)測，通過前期有效的識別規(guī)避泥石流災(zāi)害，將泥石流風(fēng)險降到最低[2]。傳統(tǒng)泥石流危險性評價難以達(dá)到預(yù)測泥石流發(fā)生過程的目的[3]，因此，研究擬引入數(shù)值方法對泥石流進行研究分析[4]。數(shù)值方法不僅計算精度高、客觀性強，更能分析動態(tài)識別泥石流的流速、沖擊狀態(tài)、堆積狀態(tài)等多元信息，因此有不少學(xué)者相繼開展了相關(guān)的研究工作。宋兵等[5]運用RAMMS模型模擬了白沙溝泥石流運動特征；Hungr[6]使用 DAN 3D模型對泥石流的流動過程進行了模擬；Christopher等[7]使用Laharz模型對(圣卡塔利納山脈7條)泥石流的運動距離和堆積范圍進行了模擬；林志遠(yuǎn)等[8]采用二維洪水與土石流數(shù)值模擬軟件(FLO-2D)對臺風(fēng)誘發(fā)泥石流進行了分析；常鳴等[9]運用FLO-2D對多種降雨頻率條件下的危險性進行評價。上述數(shù)值分析方法在泥石流過程模擬方面均取得了較好的效果，尤其是FLO-2D模型在泥石流危險性模擬方面潛力巨大。

研究基于FLO-2D模型從多維度分析和評價了三元區(qū)文筆山魯坑儲木場泥石流溝的危險性狀況，并以近期發(fā)生的山洪泥石流情況驗證泥石流評價結(jié)果的合理性，為相關(guān)研究工作提供依據(jù)。

1 泥石流場地概況

1.1 泥石流溝特征

魯坑儲木場沖溝(下稱LK沖溝)位于福建三明三元區(qū)境內(nèi)，該沖溝長約3 140 m，流域面積約1.65 km2，呈北西西-南東東向展布，后緣發(fā)育于東南側(cè)山體，前緣溝口，溝向為280°。沖溝總體處于丘陵地貌與中低山地貌交會處，頂?shù)讟?biāo)高為230～1 100 m，縱比降為0.10～0.40，總體呈上陡下緩。據(jù)此可將沖溝劃分為2個區(qū)域，即固體物源集中區(qū)和上游匯水區(qū)，各區(qū)的特征及分布情況(圖 1)。

圖1 魯坑儲木場沖溝流域分區(qū)圖Fig.1 Division drawing of the gully basin in Lukeng wood storage yard

上游匯水區(qū)高程為560～1 100 m，高差為540 m，匯水區(qū)面積為0.9 km2，該段主溝長約1 540 m，溝床平均縱坡降為351‰，溝谷呈“V”字形，溝寬為4～8 m。溝道底部主要出露侏羅世二長花崗巖，局部可見沖洪積漂石、碎(卵)石層，粗顆粒粒徑為10～50 cm。

固體物源集中區(qū)高程為237～560 m，面積約為0.75 km2，溝長約1 600 m，坡降為202‰，溝谷呈“V”字形，溝寬為5～15 m。漂石、碎(卵)石層，厚度為1～3 m。兩岸山坡基巖出露不良，坡面表層以殘坡積黏性土為主，局部可見崩坡積層，下伏花崗巖風(fēng)化層。黏土層下為基巖，巖性為花崗巖。溝道中有高大喬木。

1.2 泥石流溝發(fā)育情況

2019-05-15～2019-05-17期間，三明市普降暴雨，累計降雨量達(dá)到322.7 mm，最大降雨量為56.6 mm/h。強降雨期間，LK沖溝兩側(cè)山體發(fā)生大量的滑坡、泥石流災(zāi)害，初步統(tǒng)計規(guī)模較大的滑坡-泥石流災(zāi)害達(dá)13處，總規(guī)?？蛇_(dá)30 000 m3，溝道內(nèi)大部分滑坡體堆積于沖溝溝體中部及溝口，堵塞溝道。在強降雨情況下沖溝上游較大的匯水(1.3 km2)促使山洪形成，山洪攜帶溝內(nèi)松散沉積物，在溝口與淤塞的滑坡堆積體共同形成密度高、沖擊力強的泥石流，對溝口居民區(qū)造成沖擊，破壞部分民房及排洪溝。

此次調(diào)查在流域內(nèi)分別測量泥石流溝道(堆積物與上游洪水匯集后)中部(斷面1)和中下部(斷面2)2處洪痕斷面。

洪峰斷面1：斷面整體呈矩形(圖2)，溝寬為3.7 m，上游至下游，沖溝洪痕斷面位置兩岸近于垂直，計算時取值為90°，測量洪痕的水位約為0.80 m，水力坡度為5°，河槽無草樹，河岸較陡，岸坡樹叢過洪時淹沒溝道植被發(fā)育，溝底有碎塊石，曼寧系數(shù)取0.05。

圖2 LK沖溝洪峰斷面1Fig.2 Section 1 of LK gulch flood peak

洪峰斷面2：斷面整體大致呈矩形(圖 3)，溝寬為4.1 m，上游至下游，沖溝洪痕斷面位置兩岸近于垂直，計算時取值為90°，當(dāng)時洪水的水位約為0.72 m，水力坡度為6°，河岸較陡，河底有卵石、大孤石，曼寧系數(shù)取0.05。

圖3 LK沖溝洪峰斷面2Fig.3 Section 2 of LK gulch flood peak

根據(jù)LK沖溝洪峰斷面計算結(jié)果(表1)，斷面1洪峰流量為11.63 m3/s，斷面2洪峰流量為12.47 m3/s。

表1 LK沖溝洪峰斷面流量計算Table 1 Peak peak flow of LK ditch section

2 泥石流數(shù)值模擬原理及參數(shù)率定

2.1 基本原理

研究主要利用FLO-2D程序開展區(qū)域內(nèi)流域的泥石流過程模擬。FLO-2D[10-11]基本原理是以某次降雨強度輸入量并假設(shè)泥石流運動過程的動量守恒。

連續(xù)性方程：

①

式中：t為時間，h為流深，vx、vy為速度，i為有效降雨強度。

在此基礎(chǔ)上，利用非牛頓流體模型與中央有限差分?jǐn)?shù)值方法求解運動控制方程。所求控制性方程包括連續(xù)性方程、運動方程和流變方程(公式①～③)。其中，連續(xù)性方程保證泥石流物質(zhì)單元間的質(zhì)量守恒，運動方程用于控制單元間的運動平衡，流變方程用于求解泥石流運動時顆粒之間的碰撞對泥石流流動阻力的影響。

運動方程：

②

式中：Sfx、Sfy為摩擦坡降；Sbx、Sby為底床坡降。

流變方程：

③

式中：Sy為屈服坡降，Sv為動力粘滯系數(shù)，Std為紊流-離散坡降；τy為屈服應(yīng)力，η為動力粘滯系數(shù)；γm為泥石流單位重，K為層流阻滯系數(shù)，n為曼寧系數(shù)。粘滯系數(shù)η和屈服應(yīng)力τy的大小主要取決于泥沙體積濃度，關(guān)系式如下。

τy=α1eβ1Cvη=α2eβ2Cv

④

通過FLO-2D軟件進行區(qū)域網(wǎng)格劃分，并將建立的方程在FLO-2D軟件上實現(xiàn)，即可形成泥石流運動模型架構(gòu)。在基本架構(gòu)下，根據(jù)高分辨率DEM數(shù)據(jù)及遙感影像數(shù)據(jù)獲取各網(wǎng)格區(qū)域的高程、坡降等地形數(shù)據(jù)；結(jié)合區(qū)域覆蓋層及植被發(fā)育情況，并根據(jù)經(jīng)驗數(shù)據(jù)和區(qū)域歷史數(shù)據(jù)校正，即可獲取包括屈服應(yīng)力、黏度等流變系數(shù)；泥石流流體的基本物理參數(shù)則可根據(jù)野外調(diào)查、室內(nèi)試驗及軟件提供的經(jīng)驗參數(shù)確定。獲取的參數(shù)利用FLO-2D軟件分析計算，通過運動方程和連續(xù)方程可以計算出每個網(wǎng)格中流體深度、流體流量，進而得知流體的運動范圍和物質(zhì)總量。通過動量方程計算出相鄰網(wǎng)格間流體的速度變化和動靜壓力，計算出泥石流的動量和實際沖擊力，從而預(yù)測泥石流流動和堆積范圍，最終得出泥石流的運動狀態(tài)和沖擊力。

2.2 模型建設(shè)及參數(shù)率定

根據(jù)泥石流災(zāi)害現(xiàn)場照片，判斷該泥石流為稀性泥石流，采集泥石流的沉積物質(zhì)，采用現(xiàn)場配漿計算泥石流的容重，運用FLO-2D模型分析計算，確定泥石流流體平均重度為1.37 g/m3，綜合選取體積濃度CV為30%。

2處沖溝泥石流流域為灌木植被區(qū)域，沖溝植被發(fā)育較好，利用FLO-2D模型分析計算，2處沖溝的層流阻滯系數(shù)(K)為15 000、曼寧粗糙系數(shù)(n)為0.05；結(jié)合當(dāng)?shù)匾酝?jīng)驗，2處沖溝的屈服應(yīng)力分別為0.821、0.006 22，動力粘滯系數(shù)分別為15.23、13.72。

根據(jù)魯坑儲木場周邊地質(zhì)環(huán)境要素以及泥石流過程模擬所需數(shù)據(jù)，研究選擇50a一遇的洪水流量(i)為11.93 m3/s(1)福建省地質(zhì)工程勘察院，三元區(qū)群發(fā)性泥石流調(diào)查報告，2019。作為模擬輸入數(shù)據(jù)，每個沖溝選擇一個輸入點，模擬網(wǎng)格設(shè)置為30 m。模擬洪峰流量時間約為1 h，總共模擬實際時間約為5 h。

3 模擬成果及評價

3.1 泥石流運動狀態(tài)

3.1.1 深度分布

通過模擬結(jié)果看出(圖4a)，泥石流主要沿溝道運移，溝道內(nèi)影響面積為1×104m2。流體在溝口處擴散，最終形成的堆積扇面積為20×104m2。溝道內(nèi)泥石流流體深度多數(shù)為0.5～4.1 m，在溝道轉(zhuǎn)交處附近，由于溝道的阻力和流體的彎道超高作用，使流體流深突然增大。當(dāng)泥石流沖出溝口后迅速擴散，流體深度一般在0.5 m以下。

3.1.2 流速分布

在泥石流最大速度模擬中(圖4b)，速度較大的區(qū)域主要集中在溝道內(nèi)，速度為1.3～6 m/s，尤其在溝道中上部,由于溝道縱坡降較大，局部流速超過7.5 m/s。當(dāng)泥石流到達(dá)堆積扇邊緣(下部城區(qū))后，速度降為1.3 m/s以下，這也與流體到達(dá)平坦地面的能量擴散作用相吻合。

3.1.3 壓力分布

壓力分析主要結(jié)合泥石流的重度，通過泥石流壓力判別泥石流的沖擊力(圖4c、4d)。泥石流在中上游的溝道底部附近壓力較大，受溝道的限制，在該區(qū)域內(nèi)具有較高的流速和流深，重度相對較大，沖擊力較強。而當(dāng)泥石流沖出溝口后，隨著側(cè)向限制的消失，泥石流迅速擴散，流速也隨之下降，沖擊力則明顯減小。

圖4 泥石流模擬成果圖Fig.4 Debris flow simulation results diagram

3.2 泥石流模擬成果評價

對比2處斷面實測泥石流數(shù)據(jù)及模擬數(shù)據(jù)，模擬得出的泥石流在溝道內(nèi)的流速與實測值接近，模擬效果較好(表2)。在流深方面均要略大于實測數(shù)據(jù)，尤其在靠近泥石流流量輸入端(固體物源集中區(qū))，流深明顯較大。其主要原因在于受軟件限制，流量輸入端輸入流量為局部的集中輸入，與實際的面狀流體存在差異，流體在輸入端附近未完全擴散，導(dǎo)致模擬成果失真。隨著模擬流體向下游運動，其流深、流速逐漸接近真實值(下游斷面2逐漸逼近真實值)，最終模擬得出流量也逐漸接近真實值。因此，在泥石流模擬過程中應(yīng)注意流量輸入口的選取，在合理設(shè)計輸入口后，可以較好地模擬出泥石流實際運動過程和狀態(tài)。泥石流流域范圍方面，模擬的泥石流范圍較實際范圍略小，主要原因在于下游樓房溝道較多，溝道粗糙系數(shù)在此急劇增加，極大的影響泥石流的流通，但模擬整體擴散圖形與實際圖形相似，也在一定程度上說明模擬模型的科學(xué)性(圖5)。

表2 泥石流模擬結(jié)果與實際情況對比Table 2 The comparison between simulation results and actual situation of the debris flow

圖5 實際泥石流范圍與模擬成果對比圖Fig.5 Comparison diagram of the actual debris flow area and the simulation results

4 結(jié)論

(1)實際調(diào)查顯示，魯坑儲木場泥石流總規(guī)模為30 000 m3，成災(zāi)機理主要是在強降雨情況下，溝道內(nèi)大部分滑坡體堆積于沖溝中部及溝口，堆積物堵塞溝道。在沖溝上游較大的匯水(1.3 km2)誘發(fā)下，形成泥石流。

(2)模擬成果顯示，魯坑儲木場泥石流流速為1.3～7.5 m/s，泥石流流體深度為0.5～4.1 m，在溝道轉(zhuǎn)交附近，由于溝道的阻力和流體的彎道超高作用，使流體流深突然增大，當(dāng)泥石流沖出溝口后迅速擴散，流體深度一般在0.5 m以下。最終形成的堆積扇面積約為20×104m2。

(3)模擬得出的泥石流在溝道流速與實測值接近。在流深方面，受軟件設(shè)置流量出口的限制，模擬成果略大于實測數(shù)據(jù)。

(4)FLO-2D軟件可以利用連續(xù)性方程、運動方程和流變方程有效構(gòu)建泥石流運動模型，并計算出泥石流的流速、流深和沖擊力等信息，在合理設(shè)計模型和參數(shù)情況下，可較好重現(xiàn)或預(yù)測泥石流發(fā)生情況。