林子揚(yáng), 王衛(wèi)紅, 楊宏珺

(1.西南科技大學(xué)環(huán)境與資源學(xué)院, 四川 綿陽 621010; 2.國家遙感中心綿陽科技城分部, 四川 綿陽 621010;3.西南科大四川天府新區(qū)創(chuàng)新研究院, 四川 成都 610299)

1 前言

我國境內(nèi)存在1.2萬余座尾礦庫,其中大多數(shù)為上游式筑壩法構(gòu)建的尾礦庫。為了使每年尾礦排量能低成本處理,同時(shí)延長尾礦庫的使用壽命,需要在現(xiàn)有尾礦庫的基礎(chǔ)上加高現(xiàn)有堆積壩[1]。從技術(shù)管理的問題來看,壩高越高,潰壩后影響范圍越大,尾礦庫的安全存在很大的潛在風(fēng)險(xiǎn)。

在尾礦庫潰壩方面,國內(nèi)外學(xué)者已開展系列研究。李全明等[2]建立尾礦庫潰決泥漿運(yùn)動(dòng)模型,提出尾礦庫潰壩范圍計(jì)算的評估方法;劉嘉欣等[3]研究漫頂潰壩機(jī)理,通過沿深度取積分平均值建立潰壩的二維數(shù)學(xué)模型;崔旭陽等[4]構(gòu)建了一種包含動(dòng)態(tài)指標(biāo)和靜態(tài)指標(biāo)的評價(jià)指標(biāo)體系,利用動(dòng)態(tài)加權(quán)貝葉斯網(wǎng)絡(luò)和時(shí)間權(quán)重、指標(biāo)權(quán)重來評估尾礦庫潰壩風(fēng)險(xiǎn);胡航等[5]強(qiáng)降雨洪水條件下下泄尾砂總量與潰壩后果之間的關(guān)系,并結(jié)合FLO-2D模擬軟件對潰壩災(zāi)害進(jìn)行了預(yù)測分析;呂松峰等[6]模擬了尾礦庫在瞬時(shí)局部潰壩條件下尾砂的運(yùn)動(dòng)情況。Rico M 等對事故尾礦庫的壩高、庫容與潰壩量等資料進(jìn)行回歸分析,得到尾礦庫潰壩量約占尾礦庫總庫容的1/3(包括水)。Concha Larrauri P等提出尾礦庫在潰壩后不會(huì)全部泄完。

2001—2018年我國共發(fā)生95起尾礦庫事故,事故中約380人喪生,其中潰壩事故占比高達(dá)56.9%[7-8]。因此,需要在加強(qiáng)尾礦庫日常監(jiān)管,降低事故隱患的同時(shí),對尾礦庫進(jìn)行潰壩模擬實(shí)驗(yàn),對實(shí)驗(yàn)結(jié)果和潰壩影響范圍進(jìn)行數(shù)據(jù)處理,得到應(yīng)對尾礦庫潰壩災(zāi)害時(shí)能提前對影響范圍內(nèi)的居民進(jìn)行預(yù)警,保障尾礦庫下游居民人身財(cái)產(chǎn)安全[9-10]。

綜上所述,目前研究僅限于單一洪水潰壩,而實(shí)際上,在遭遇地震或其他自然災(zāi)害時(shí),有可能發(fā)生單點(diǎn)或多點(diǎn)瞬時(shí)潰壩。因此,結(jié)合尾礦庫洪水潰壩數(shù)值模擬,局部尾礦庫瞬時(shí)潰壩數(shù)值模擬,提出尾礦庫多點(diǎn)瞬時(shí)潰壩數(shù)值模擬,研究采用尾礦庫安全評價(jià)項(xiàng)目,通過FLO-2D模擬尾礦庫洪水漫頂,單點(diǎn)瞬時(shí)潰壩和多點(diǎn)瞬時(shí)潰壩,來研究尾礦庫不同形式的潰壩在強(qiáng)度和影響范圍的差異,得到單點(diǎn)瞬時(shí)潰壩破壞力強(qiáng),影響范圍略小于多點(diǎn)瞬時(shí)潰壩,但兩種潰壩形式在強(qiáng)度和影響范圍上都大過于洪水漫頂潰壩成果,彌補(bǔ)尾礦庫潰壩實(shí)驗(yàn)形式單一的情況[11]。

2 研究區(qū)現(xiàn)狀

尾礦庫位于米易縣埡口鎮(zhèn)安寧河西岸爛壩山下的沖溝中,與選廠隔河相望,直線距離1 km,距撒蓮鎮(zhèn)政府5 km,距米易縣政府20 km,具體位置如圖1所示。

圖1 尾礦庫位置圖

該尾礦庫修筑方式為上游式筑壩法,2022年之前該尾礦庫存在3座壩體,經(jīng)過不斷堆高修筑,現(xiàn)如今合并為一座壩體的尾礦庫,設(shè)計(jì)總壩高為246 m,總庫容設(shè)計(jì)為5 079.5萬m3,有效庫容為4 800萬m3,該尾礦庫由三等庫升級為二等庫,位于山谷中,初期壩下游是安寧河,有少量耕地,無工礦企業(yè)、水產(chǎn)基地及名勝古跡等重要設(shè)施。

庫區(qū)屬熱帶、亞熱帶森林氣候,季節(jié)干濕分明,夏天炎熱、多雨、潮濕,而冬天溫暖、多雨、干燥。年平均氣溫21～22 ℃,最高38.7 ℃,雨季期間,陡峭的河谷兩側(cè)頻繁發(fā)生山洪、河流暴漲和小型泥石流。

3 數(shù)據(jù)和方法

尾礦庫初期設(shè)計(jì)為溝口處修建3座壩體,隨著壩體升高,將3座壩體合并為1座壩體,因此選取該尾礦庫呈現(xiàn)出單點(diǎn)和多點(diǎn)潰壩模擬的成果會(huì)有明顯。

選取的不同潰壩條件分別為洪水漫頂潰壩數(shù)值模擬;瞬時(shí)單點(diǎn)潰壩數(shù)值模擬;瞬時(shí)多點(diǎn)潰壩數(shù)值模擬。其中洪水漫頂潰壩數(shù)值是以某斷面流體流量隨時(shí)間變化的過程來客觀反映泥石流災(zāi)害整體運(yùn)動(dòng)過程。瞬時(shí)潰壩數(shù)值模擬相較于洪水漫頂,主要在于其突發(fā)性和不可預(yù)見性,為了保障下游地區(qū)人民人身財(cái)產(chǎn)安全,需提前做好此類事故預(yù)測工作,為此選取瞬時(shí)單點(diǎn)潰壩。在此基礎(chǔ)上,由于尾礦庫壩頂寬度過長,在突發(fā)情況下不止一處會(huì)發(fā)生崩壞等現(xiàn)象,因此在實(shí)驗(yàn)中選取多點(diǎn)瞬時(shí)潰壩來模擬該情景下,下游受災(zāi)情況。

在實(shí)驗(yàn)軟件上,采用FLO-2D軟件實(shí)現(xiàn)尾礦庫潰壩數(shù)值模擬實(shí)驗(yàn),其具有對泥石流有較好的重現(xiàn)性,能表示泥深、流速,以及反映流體流變時(shí)空特征等優(yōu)點(diǎn)[12]。該軟件通過分析地形數(shù)字高程模型為實(shí)驗(yàn)藍(lán)本,將尾礦庫區(qū)域的地形制作成二維單元網(wǎng)格形式,計(jì)算得出每個(gè)網(wǎng)格單元在實(shí)驗(yàn)中對應(yīng)的流速,泥沙淹沒深度以及流量數(shù)據(jù)。

3.1 模擬流體參數(shù)設(shè)定

本次模擬是基于體積守恒和水力學(xué)原理構(gòu)建流體質(zhì)量守恒和動(dòng)量守恒方程的微分形式,采用顯式中心差分法對守恒方程進(jìn)行求解運(yùn)算,可以有效計(jì)算單元網(wǎng)格內(nèi)流體流深、流速和影響范圍的時(shí)間變化。

鑒于模擬模型原理,模擬流體以一種恒定的固液混合流體形式進(jìn)行運(yùn)動(dòng)堆積,其需要設(shè)定以下流體參數(shù):體積濃度折中取值0.5[13]。泥石流重度取值1.25 kg/m3。層流阻力系數(shù)本次模擬參照Woolhiser 實(shí)地測定的泥石流層流阻力系數(shù)設(shè)定該值為2 200[14]。屈服應(yīng)力和黏滯系數(shù):流體的屈服應(yīng)力選擇詹前登實(shí)驗(yàn)參數(shù)。

3.2 地表環(huán)境參數(shù)設(shè)定

在本次模擬中,根據(jù)高分辨率航空圖像,將數(shù)值模擬區(qū)域分類為地物識別,并根據(jù)地物類別估算滲透性和阻礙能力。通過曼寧系數(shù)描述地表阻礙能力,使用Green Ampt透射模型描述表面滲透能力[15]。

4 數(shù)值模擬

4.1 洪水漫頂潰壩數(shù)值模擬流量曲線

參照防排洪驗(yàn)算和尾礦庫建庫標(biāo)準(zhǔn),以1000年一遇的洪水頻率為尾礦庫潰決計(jì)算頻率,依照水文計(jì)算洪水流量和模擬流體參數(shù)計(jì)算庫區(qū)斷面的泥石流流量過程線如圖2所示。

圖2 泥石流流量過程線示意圖

圖2中1#,2#,3#代表3座壩體,以最北方向?yàn)?#壩體依次向南方排列。洪水流量過程線表示尾礦庫由于承受不住洪水容量潰壩后,洪水向壩體下游流動(dòng)的過程;泥石流流量過程線表示尾礦庫由于洪水漫頂后產(chǎn)生潰壩,洪水夾帶著尾砂向壩體下游流動(dòng)的過程。

4.2 洪水漫頂潰壩數(shù)值模擬成果

基于相關(guān)資料與參數(shù),開展尾礦庫洪水漫頂潰壩數(shù)值模擬,其計(jì)算結(jié)果如圖3所示。

圖3 尾礦庫洪水漫頂數(shù)值模擬結(jié)果示意

根據(jù)數(shù)據(jù)模擬成果分析,在1000年一遇的洪水頻率下,尾礦庫一旦發(fā)生潰決對尾礦庫下游地表環(huán)境造成嚴(yán)重破壞。尾礦庫壩體涉及范圍廣,下游共有3處流域通往安寧河。依據(jù)本次數(shù)值模擬,模擬災(zāi)害流體分別從壩體上的3處位置運(yùn)動(dòng)至下游溝口,模擬流體平均速度約為0.3 m/s,模擬體過流的最大厚度集中在2 m以下,占災(zāi)害影響面積的82.6%,最大沉積厚度在4 m以上的區(qū)域占總沉積面積的3%,泥深最高沉積厚度為7.3 m。

4.3 瞬時(shí)單點(diǎn)潰壩數(shù)值模擬泥石流流量計(jì)算

1)潰口寬度計(jì)算

潰口寬度是尾礦庫潰壩災(zāi)害破壞程度的主要影響因子之一,同時(shí)也是計(jì)算下泄流量的重要參數(shù)。通過水庫潰決資料得到計(jì)算潰口寬度的計(jì)算公式[16]。

(1)

式中:b——潰口平均寬度,m;

W——尾礦庫潰壩總量,萬m3;

B——尾礦庫水面寬度;

Hw——尾礦庫水面高,m;

K1——筑壩材料相關(guān)系數(shù),如黏土材質(zhì)取0.65,土壤材質(zhì)取1.3。

2)潰壩總量計(jì)算

尾礦庫災(zāi)害的關(guān)鍵取決于潰壩總量的分析,對于計(jì)算尾礦庫潰壩總量在世界上沒有統(tǒng)一的標(biāo)準(zhǔn),在數(shù)值模擬中將尾礦庫全庫容作為潰壩總量數(shù)據(jù),但在實(shí)際潰壩案例中,潰壩總量總是小于尾礦庫庫容,只存在部分庫容通過潰壩泄出。不同尾礦庫潰壩方式存在諸多不同的影響因素,因此尾礦庫潰壩總量的計(jì)算存在差異性和不唯一性。根據(jù)Rico M等的研究,本次實(shí)驗(yàn)選取潰壩量為總庫容的1/3。對此次潰壩模擬采用的潰壩總量為5 034萬m3,選取有效庫容為1 676.4萬m3,占尾礦庫總庫容的33.3%。

3)壩址處潰壩流量過程線計(jì)算

通過對水庫潰壩計(jì)算方法比對,主要是采用圣維南法,經(jīng)驗(yàn)公式等方法。對于尾礦庫潰壩采用借鑒水庫潰壩計(jì)算方法[17],其公式為

(2)

式中:QM——潰壩流沙的最大流量,m3/s;

m——不同形式的溝谷斷面參數(shù);

H0——尾礦庫上有水深。

對于不同形式的溝谷斷面,可選取相對應(yīng)或組合模型進(jìn)行取值,具體參數(shù)見表1。由于尾礦庫壩體形狀介于二次拋物線和組合拋物線之間,m選取為5/2。

表1 不同溝谷斷面形狀的特征參數(shù)

瞬時(shí)潰壩條件下,在壩體發(fā)生突然破壞后,壩址處下泄尾砂流量在開始時(shí)達(dá)到最大值,隨后快速下降,當(dāng)?shù)竭_(dá)臨界值時(shí),緩慢下降,對潰壩流量采用4次拋物線來進(jìn)行概化。此次潰壩流量4次拋物線成果見表2。

表2 潰壩流量4次拋物線概化參數(shù)

表2中,T表示為潰壩砂流下泄時(shí)長,t表示為潰壩下泄時(shí)間刻度;QM表示為潰壩流砂的最大流量;Qt表示為t時(shí)刻砂流流量。

繪制潰壩尾砂流流量過程曲線,主要參考潰壩中水力學(xué)經(jīng)驗(yàn)公式,采用以下計(jì)算參數(shù):

(1)潰壩總量W的計(jì)算。采用總庫容潰壩1/3,潰壩總量取值為1 676.4萬m3。

(3)尾礦庫潰壩砂流最大流量QM。

4.4 瞬時(shí)單點(diǎn)潰壩數(shù)值模擬流量曲線

參照水文電力手冊和尾礦庫建庫標(biāo)準(zhǔn),以尾礦庫單點(diǎn)潰壩計(jì)算頻率[18],依照模擬流體參數(shù)計(jì)算庫區(qū)斷面的尾砂流量過程線如圖4所示。

圖4 瞬時(shí)單點(diǎn)潰壩流量過程曲線圖

4.5 瞬時(shí)單點(diǎn)潰壩數(shù)值模擬成果

基于相關(guān)資料與參數(shù),開展尾礦庫瞬時(shí)單點(diǎn)數(shù)值模擬,其計(jì)算結(jié)果如圖5所示。

圖5 瞬時(shí)單點(diǎn)潰壩模擬成果圖

依據(jù)本次數(shù)值模擬,模擬災(zāi)害流體分別從壩體上的3處位置運(yùn)動(dòng)至下游溝口,模擬流體平均速度約為17.7 m/s,由于本次模擬有大量尾砂沖入安寧河內(nèi),選取模擬最大值位于地面部分,模擬體過流的最大厚度集中在6.5 m以下,占災(zāi)害影響面積的76.8%,最大沉積厚度在6.5 m以上的區(qū)域占總沉積面積的23.2%,尾砂最高沉積厚度為23.21 m。

4.6 瞬時(shí)多點(diǎn)潰壩數(shù)值模擬泥石流流量過程線

瞬時(shí)潰壩條件下,在壩體發(fā)生突然破壞后,壩址處下泄尾砂流量在開始時(shí)達(dá)到最大值,開始快速下降,當(dāng)?shù)竭_(dá)臨界值時(shí),緩慢下降。

參照水文電力手冊和尾礦庫建庫標(biāo)準(zhǔn),以尾礦庫3條溝口壩體壩長比作為權(quán)重參數(shù),本次實(shí)驗(yàn)中權(quán)重比為1#∶2#∶3#= 0.351 314 85∶0.326 197 28∶0.322 487 871,將潰口設(shè)置在三座壩體上計(jì)算頻率,依照模擬流體參數(shù)計(jì)算庫區(qū)斷面的尾砂流量過程線如圖6所示。

圖6 瞬時(shí)多點(diǎn)潰壩流量過程曲線圖

4.7 瞬時(shí)多點(diǎn)潰壩數(shù)值模擬成果

基于相關(guān)資料與參數(shù),開展尾礦庫多點(diǎn)瞬時(shí)潰壩數(shù)值模擬,其計(jì)算結(jié)果如圖7所示。

圖7 瞬時(shí)多點(diǎn)潰壩模擬成果圖

依據(jù)本次數(shù)值模擬,模擬災(zāi)害流體分別從壩體上的3處位置運(yùn)動(dòng)至下游溝口,模擬流體1#壩平均速度約為2.18 m/s,模擬流體2#壩平均速度約為2.17 m/s,模擬流體3#壩平均速度約為2.14 m/s,由于本次模擬有大量尾砂沖入安寧河內(nèi),選取模擬最大值位于地面部分,模擬尾砂的最大厚度集中在5.4 m以下,占災(zāi)害影響面積的80.6%,6.5 m以上最大沖積厚度區(qū)域占整個(gè)災(zāi)害影響面積的19.4%,尾砂最高沉積厚度為23.56 m。

5 下游威脅分析對比

依據(jù)洪水漫頂潰壩數(shù)值模擬結(jié)果,模擬災(zāi)害流體沿著尾礦庫下游溝道向下運(yùn)動(dòng),并在流域溝口處進(jìn)行了擴(kuò)展堆積,形成扇形的破壞區(qū),從而對溝口的居民建筑、農(nóng)田耕地和交通干道造成了大面積的淤埋,甚至有流體向下運(yùn)動(dòng)堆積至安寧河流域造成河道堵塞,引發(fā)次生災(zāi)害?，F(xiàn)選取3處下游流域溝口區(qū)域進(jìn)行模擬破壞分析,下游影響區(qū):該溝口區(qū)域距離壩頂約1.2 km,溝口區(qū)域內(nèi)模擬流體橫向波及長度達(dá)220 m,約15 984 m2耕地被模擬流體淤埋,1棟居民建筑邊緣被流體波及,1棟疑似臨時(shí)建筑的房屋被流體淤埋破壞,但災(zāi)害流體并未繼續(xù)向下堆積至安寧河道之中。

依據(jù)尾礦庫瞬時(shí)潰壩數(shù)值模擬結(jié)果,潰壩的尾砂沖過安寧河達(dá)到對岸,瞬時(shí)潰壩最遠(yuǎn)處距離壩頂約為2.3 km,溝口區(qū)域內(nèi)模擬流體橫向波及長度達(dá)1 506 m,約有3.35×105m2農(nóng)田耕地被尾砂覆蓋,其中23棟居民建筑邊緣被流體波及,4棟建筑的房屋被流體淤埋破壞。其中單點(diǎn)瞬時(shí)潰壩對潰壩點(diǎn)處的壩體造成嚴(yán)重影響,對另外兩座壩體影響較小,而多點(diǎn)瞬時(shí)潰壩分別對壩體造成的影響相差不大,因此單點(diǎn)瞬時(shí)潰壩對2#壩體下的房屋耕地受災(zāi)最為嚴(yán)重。

6 結(jié)論

本文利用對四川省攀枝花市米易縣某尾礦庫DSM影像進(jìn)行了處理,再通過FLO-2D軟件對尾礦庫進(jìn)行多種潰壩數(shù)值模擬實(shí)驗(yàn),得到以下結(jié)論:

(1)不同潰壩條件下的潰壩流量不同,瞬時(shí)單點(diǎn)潰壩的潰口流量速度為17.7 m/s,短時(shí)間內(nèi)對下游地區(qū)的沖刷破壞性較強(qiáng),瞬時(shí)多點(diǎn)潰壩的潰口流量相對有所降低,速度為2.17 m/s,具有持續(xù)性破壞作用,洪水漫頂造成的潰壩災(zāi)害流量為0.3 m/s最小,造成的破壞性最小,但由于伴隨著大量降雨,潰壩災(zāi)害的同時(shí)伴隨著洪水等次生災(zāi)害。

(2)不同潰壩條件下的潰壩破壞強(qiáng)度不同,瞬時(shí)潰壩模擬結(jié)果中,單點(diǎn)瞬時(shí)潰壩對對應(yīng)壩體下游受災(zāi)地區(qū)破壞力最強(qiáng),其余壩體下游受災(zāi)相較于多點(diǎn)潰壩較弱,多點(diǎn)潰壩災(zāi)害受災(zāi)程度較均勻,受災(zāi)面積擴(kuò)大,但兩種潰壩形式下受災(zāi)區(qū)域大致相同。

(3)不同潰壩條件下的尾砂沉積方式也不同,洪水漫頂潰壩中尾砂堆積主要在初期壩下游地區(qū);單點(diǎn)瞬時(shí)潰壩尾砂主要堆積于主體受災(zāi)壩中后段,堵塞初期壩;多點(diǎn)瞬時(shí)潰壩尾砂堆積于潰壩口和初期壩位置,在河道處沉積的尾砂也是最多的。

(4)在三種潰壩模擬中,洪水漫頂潰壩破壞力最小,尾砂沉積量最少;單點(diǎn)潰壩模擬,破壞力最強(qiáng),尾砂主要集中于主壩體口;多點(diǎn)潰壩模擬,潰壩破壞力較分散,破壞力由兩邊向中間集中,尾砂同時(shí)也主要集中于壩體中部,對受災(zāi)區(qū)域沖擊面積最廣。

(5)尾礦庫潰壩后,尾砂大多流入安寧河中,使其破壞程度得到有效削減,減少了潰壩災(zāi)害對城鎮(zhèn)的直接危害。