陳 禎，陳俊智，馬毓婷，安 南，陳文倩

(昆明理工大學(xué) 國土資源工程學(xué)院，昆明 650032)

尾礦庫潰壩事故是影響較大的工程災(zāi)害事故，在礦山的安全生產(chǎn)中備受礦山企業(yè)關(guān)注與重視。在礦山生產(chǎn)中由于選礦后產(chǎn)生的尾砂材料需要長期積壓堆存。當(dāng)潰壩發(fā)生時(shí)庫內(nèi)的尾砂與水混合后會(huì)以類似“泥石流”的形式下涌，嚴(yán)重危害所涉及區(qū)域的人員以及財(cái)產(chǎn)安全[1]。此外，若災(zāi)害發(fā)生，尾砂中暫時(shí)不能處理的有害物質(zhì)泄露會(huì)造成周圍環(huán)境的污染[2]。且我國現(xiàn)存的上萬座尾礦庫大約有一成存在各種問題[3]。因此，加強(qiáng)對尾礦庫潰壩事故影響范圍的研究并建立相關(guān)數(shù)值模型，能夠有效減少潰壩災(zāi)害發(fā)生時(shí)造成的生命以及財(cái)產(chǎn)損失，對保障尾礦庫的安全運(yùn)行具有重大意義。

迄今為止，就潰壩事故的產(chǎn)生機(jī)理以及危險(xiǎn)系數(shù)和影響波及的范圍來看，國內(nèi)外學(xué)者已經(jīng)做了較多成果顯著的研究工作。BLIGHT[4]通過對環(huán)形尾礦壩的潰壩分析得出尾礦所處地表的濕度情況會(huì)影響尾礦漿的擴(kuò)散范圍，地表濕度較高的地方尾礦漿的移動(dòng)擴(kuò)散范圍較遠(yuǎn)；PASTER等[5]以彈塑性有限元法為原理建立相關(guān)模型，該模型可以用在滑坡、泥石流、潰壩等災(zāi)害的模擬中推算流體流動(dòng)隨時(shí)間堆積的程度；PIRULLI等[6]通過對比分析碎屑流與潰壩泥流，探究兩者共性，并在尾礦庫潰壩事故分析中運(yùn)用RASH 3D軟件提出一個(gè)流變模型。

國內(nèi)對潰壩等災(zāi)害現(xiàn)象的研究成果頗豐，陳青生、孫建華等[7]以建立數(shù)學(xué)模型的方式找出了潰壩現(xiàn)象發(fā)生后引起的泥漿流動(dòng)的表現(xiàn)規(guī)律；周彪[8]通過分析研究不同工況對滲流場的影響，并結(jié)合FLAC3D軟件分析壩體的滲流穩(wěn)定性；余樂文等[9]利用FLO 2D軟件對潰壩現(xiàn)象模擬分析，研究了潰壩后泥石流的堆積情況，并且預(yù)測了潰壩發(fā)生后的波及范圍；金佳旭等[10]利用ANSYS CFX軟件分析遼寧某礦山的尾礦庫潰壩事件，得出了潰壩所需時(shí)間、影響范圍和尾砂堆積范圍以及高度等情況。

據(jù)上，國內(nèi)外學(xué)者主要考慮尾礦潰壩泥漿動(dòng)力特征、滲流和潰壩后的影響范圍及尾砂堆積角度采用了相關(guān)理論分析和數(shù)值模擬方法分析尾礦壩穩(wěn)定性[11]，并利用基于巖土工程及邊坡工程成熟的理論和工程經(jīng)驗(yàn)與相應(yīng)可靠的分析軟件進(jìn)行數(shù)值模擬仿真。

因此，本文以云南某礦山磷石膏新建尾礦庫為例，根據(jù)現(xiàn)場測量參數(shù)，應(yīng)用相關(guān)公式設(shè)計(jì)計(jì)算三維模型大小，運(yùn)用Gambit軟件建立實(shí)際尺寸的三維模型并劃分輸入網(wǎng)格，運(yùn)用ANSYS FLUENT軟件進(jìn)行潰壩數(shù)值模擬。本文著重探討了潰壩后的影響范圍，為建庫安全評價(jià)和災(zāi)害預(yù)警防范提供有力依據(jù)及參考。

1 磷石膏尾礦庫工程概況分析

云南某礦山磷石膏尾礦庫庫區(qū)位于山頂沖溝處，微地貌為山間溝谷、山腳坡地，場地東西兩面為山，溝谷近東西向，東高西低。從庫區(qū)到外界的交通比較方便，有建設(shè)好的公路。庫區(qū)地表分水嶺為近似圓弧形連綿山脈，南端最高山峰高程約為2 004 m；以溝谷為中心線，形成基本對稱的兩側(cè)岸，坡度角約為10°～20°，表面有植被覆蓋。但當(dāng)尾礦庫發(fā)生潰壩時(shí)，由于其位于山頂沖溝處，其庫漿泥沙有極大的重力勢能。所以，發(fā)生潰壩時(shí)庫漿將在極短的時(shí)間內(nèi)以極快的速度對下游造成毀滅性的破壞。

1.1 尾礦性質(zhì)

尾礦以粉粒為主，塑性指數(shù)IP為7.8，小于10；按照《尾礦堆積壩巖土工程技術(shù)規(guī)范》(GB50547-2010)中粒徑組類的劃分標(biāo)準(zhǔn)：擬堆積的磷石膏尾礦為尾粉土。并且磷石膏尾礦滲透系數(shù)為2.52×10-5cm/s，滲透性等級屬弱透水，尾礦pH值介于1.95～3.04，浸出液pH值介于1.46～6.4，屬酸性。

1.2 尾礦庫庫容設(shè)計(jì)概況

為增加庫容和方便防滲工作展開，云南某礦山企業(yè)將堆填區(qū)域內(nèi)的庫底及庫岸進(jìn)行清理平整、削坡。清理完成后場區(qū)與初期壩、庫尾擋水?dāng)r渣壩一同圍成一盆狀：庫底平緩，縱坡度不超過10%，庫岸坡度不陡于 1∶2。尾礦庫的計(jì)算容積見表1。

表1 尾礦庫容積計(jì)算表Table 1 Tailings pond volume calculation table

1.3 尾礦庫壩體設(shè)計(jì)概況

初期壩頂標(biāo)高1 945.0 m，壩軸處自然地面標(biāo)高1 922.0 m，溝谷段清基深度為4 m，最大壩高H初=23 m(不含清基深度)，壩頂寬B=4.0 m，壩體外坡比1∶0.5，內(nèi)坡比1∶0.5，壩軸長104.5 m。外坡設(shè)置一條1.5 m寬的道路，具體設(shè)置在標(biāo)高1 930.0 m處，內(nèi)坡在緊挨坡面處設(shè)置堆石排水棱體，排水棱體頂面標(biāo)高為1 935.0 m，頂部寬度為2.0 m，上游坡比1∶1.5，坡面上布置反濾層，可在一定程度上對上游流下來的水起引流作用，保護(hù)坡面。石壩體與岸坡的連接方式為斜面連接，這樣可以避免壩體與岸坡的連接面產(chǎn)生集中滲流現(xiàn)象，也可以減少軟弱夾層的形成。

由于堆石壩上下游坡比較緩，壩底縱向長度較大，清基工程量和筑壩石方量均較大，此外堆石壩方量較大，占用庫容空間也較大，使得設(shè)計(jì)庫容較小，無法滿足設(shè)計(jì)庫容要求。受制于尾礦庫徑深較短，因此不采用堆石壩設(shè)計(jì)。

2 潰壩尾礦漿流動(dòng)過程

尾礦庫潰壩后庫漿下泄引起的漿流本質(zhì)上屬于庫漿與水的混合物，類似泥石流，因此，本文的基本假定為庫漿砂是各向同性的連續(xù)均勻的介質(zhì)體，等同于高密度的水體。

2.1 尾礦庫潰壩下泄總量

多日連續(xù)降雨會(huì)導(dǎo)致尾礦庫內(nèi)積水，雨水與尾砂形成漿體，結(jié)合庫漿物理力學(xué)性質(zhì)，為保證設(shè)計(jì)安全，考慮潰壩發(fā)生時(shí)產(chǎn)生的潰口相應(yīng)高程處的全部庫容為該尾礦庫的下泄總量，具體庫容如表1。假設(shè)尾礦庫壩體全部容量泄出，取壩體全部容量為100萬m3。

2.2 尾礦庫潰口寬度

尾礦庫的潰口寬度參考黃河水利委員會(huì)提出的潰口寬度計(jì)算公式。

潰口寬度計(jì)算公式見式(1)。

(1)

式中，B0為尾礦庫潰口平均寬度，m；K1為壩體材料系數(shù)，尾礦庫壩體為石質(zhì)，取0.65；W為潰壩時(shí)的下泄總量，m3；B為潰壩發(fā)生時(shí)壩前的水面寬度，m；H為尾礦庫潰壩時(shí)的水深，m。

2.3 尾礦庫潰口的最大流量

假設(shè)潰壩時(shí)潰口在一瞬間完全形成，即瞬時(shí)全潰，不考慮潰口形成過程對最大流量的影響。本文采用謝任之提出的計(jì)算方法，這一方法在對尾礦庫類型的潰壩中運(yùn)用較為廣泛，潰口最大流量計(jì)算公式見式(2)[12]。

(2)

式中，QW為潰口最大流量，m3/s；λ為流量系數(shù)，由斷面形狀系數(shù)、沉溺系數(shù)、堰寬比等系數(shù)計(jì)算得出，瞬時(shí)全潰時(shí)取8/27；g為重力加速度，m/s2；H0為壩前平均水深，m。

2.4 庫漿下泄時(shí)間

潰壩發(fā)生時(shí)，潰口底部所在水平以上的所有庫容對應(yīng)的庫漿會(huì)全部泄出，下泄時(shí)間受制于庫容量和潰口高度。由于假設(shè)庫漿全部泄出因此潰口位于尾礦庫壩體底部。

庫漿下泄時(shí)間計(jì)算公式見式(3)。

(3)

式中，W為下泄總量，m3；Qm為尾礦壩潰口最大流量，m3/s；K2為修正系數(shù)，取3.5。

2.5 庫漿最大流速

潰壩發(fā)生前往往多降雨，尾礦庫內(nèi)積存大量雨水。雨水與尾礦混合后會(huì)形成飽和漿體。潰壩發(fā)生后，漿體從潰口流出，其流變性質(zhì)類似于稀性泥石流[13-14]。

庫漿最大流速計(jì)算公式見式(4)。

(4)

式中，Mm為河床粗糙系數(shù)；HC為尾礦漿的平均泥深，m；IC為尾礦庫的縱坡降，%。α為阻力系數(shù)，按式(5)計(jì)算。

(5)

式中，φC為泥石流修正系數(shù)；γg為尾礦堆積干，t/m3。

2.6 庫漿最大沖擊力

庫漿最大沖擊力計(jì)算公式見式(6)。

(6)

式中，γ為尾礦漿的飽和容重，kg/m3；Vmax為尾礦漿最大流速，m/s；g為重力加速度，m/s2；K3為泥漿不均勻系數(shù)，取3。

2.7 尾礦庫潰壩理論危害面積

尾礦庫潰壩后，庫漿從潰口流出并擴(kuò)散，擴(kuò)散后的潰壩庫漿危害面積公式見式(7)。

(7)

式中，h為尾礦漿的最大堆積厚度，m。

根據(jù)以上公式計(jì)算得出的僅是理論潰壩參數(shù)，但是在實(shí)際情況中，尾礦的堆積厚度受制于地形影響，無法達(dá)到設(shè)計(jì)厚度，因此在計(jì)算時(shí)，實(shí)際值取理論值的85%，結(jié)合理論計(jì)算公式與尾礦庫工程實(shí)際計(jì)算，得到所有參數(shù)的最終計(jì)算結(jié)果見表2。

表2 理論潰壩參數(shù)和實(shí)際潰壩參數(shù)Table 2 Theoretical and practical dam break parameters

3 潰壩數(shù)值模擬

由于各類型尾礦庫的地理環(huán)境差異巨大且復(fù)雜多變，庫容物質(zhì)也不盡相同，因此，想建立統(tǒng)一模型進(jìn)行理論分析的可能性微乎其微。現(xiàn)如今，隨著科技的發(fā)展，運(yùn)用ANSYS FLUENT數(shù)值模擬手段能夠廣泛而可靠的適用不同庫型、不同條件、各種工況下的潰壩泥沙演進(jìn)過程，進(jìn)而為理論研究提供更加行之有效的參考數(shù)據(jù)。

3.1 三維建模

由于堆積子壩的庫容最大，因此研究堆積子壩潰壩災(zāi)害對壩體下游空間造成的影響最具有代表性。使用Space Claim軟件建立500 m×500 m×200 m的尾礦庫下游空間簡化三維模型和簡化潰口三維模型105 m×100 m×35 m(其中100 m為壩體寬度)。將潰口放置在下游空間模型左側(cè)，壩體右側(cè)地勢較左側(cè)高，潰壩時(shí)，泥漿由于地勢的阻擋，會(huì)向左下側(cè)流動(dòng)。假設(shè)潰壩事故發(fā)生時(shí)潰口在一瞬間完全潰決，潰口高度為35 m，為整個(gè)堆積子壩高度；潰口長度為整個(gè)壩軸線長，即105 m；計(jì)算區(qū)域內(nèi)共有908 170 m3下泄庫漿。潰壩計(jì)算模型及尺寸示意圖見圖1。

圖1 潰壩計(jì)算模型及尺寸示意圖Fig.1 Schematic diagram of dam break calculation model and dimensions

3.2 網(wǎng)格劃分

使用ICEM CFD軟件對三維模型劃分網(wǎng)格。設(shè)置簡化潰口模型一側(cè)為入口，下游空間簡化模型一側(cè)為出口，其余模型邊界設(shè)置為壁面，如圖2所示。對整個(gè)模型進(jìn)行結(jié)構(gòu)化網(wǎng)格劃分，網(wǎng)格為正六面體，網(wǎng)格大小為：Δx×Δy×Δz=1 m×1 m×1 m，將網(wǎng)格導(dǎo)入ANSYS FLUENT后，通過模擬計(jì)算直觀得出庫漿體積量在不同時(shí)段對網(wǎng)格的渲染程度來判斷潰壩后的大致淹沒區(qū)域。劃分網(wǎng)格數(shù)量為50 367 500個(gè)，如圖2所示。

圖2 ICEM CFD網(wǎng)格劃分細(xì)節(jié)Fig.2 Details of ICEM CFD meshing

3.3 模擬方法及邊界條件設(shè)定

模擬求解器設(shè)置：

1)根據(jù)模型，設(shè)定重力作用的方向?yàn)閦軸負(fù)方向，重力加速度大小取值為-10 N/kg；

2)庫漿在運(yùn)動(dòng)過程中為不可壓，因此在求解中選擇Pressure-Based求解器；

3)庫漿是非定常流的流體，因此設(shè)定瞬態(tài)計(jì)算，時(shí)間步長為0.01 s，模擬10 000步，模擬時(shí)長100 s，大于下泄時(shí)間78 s，可以觀察整個(gè)潰壩事故的庫漿完整運(yùn)動(dòng)情況；

4)庫漿流態(tài)為湍流流動(dòng)模型設(shè)置為標(biāo)準(zhǔn)k-ε湍流模型，該模型適合模擬流體的擴(kuò)散情況，符合庫漿流出后擴(kuò)散的實(shí)際情況；

5)根據(jù)現(xiàn)場采集的磷石膏尾砂配置的庫漿，其密度在1 674.7 kg/m3，黏度為0.43 Pa·s，因此設(shè)定材料參數(shù)時(shí)，以液態(tài)水為基礎(chǔ)修改密度和黏度進(jìn)行模擬；

6)求解方法設(shè)置為PISO，該求解方法適用于非定常計(jì)算；

7)設(shè)置殘差為1×105，已知潰口處的速度條件，因此選擇標(biāo)準(zhǔn)初始化中的入口初始化。

邊界條件設(shè)置：

1)設(shè)置簡化潰口模型一側(cè)為速度入口，速度取計(jì)算值；

2)下游空間簡化模型一側(cè)為唯一自由流出口，權(quán)重為1；

3)其余模型邊界設(shè)置為靜止壁面邊界，無滑移情況；

4)壁面邊界中，壁面的粗糙高度和粗糙系數(shù)根據(jù)壁面材料設(shè)定，由于模擬的壁面為粗糙的巖土質(zhì)地面，因此取粗糙高度為0.1，取粗糙高度為0.5；

5)下游計(jì)算區(qū)域左右兩側(cè)均為封閉邊界。

邊界條件示意圖如圖3所示。

圖3 模型邊界條件示意圖Fig.3 Schematic diagram of model boundary conditions

3.4 模擬結(jié)果及分析

模擬結(jié)束后，通過生成流動(dòng)動(dòng)畫，觀測尾礦壩潰決后不同時(shí)間的庫漿淹沒范圍圖(圖4)。

圖4 尾礦壩潰決后不同時(shí)間的庫漿淹沒范圍圖Fig.4 Map of slurry inundation range at different times after tailings dam break

由模擬結(jié)果的不同時(shí)刻截圖可知：1)在潰壩初期，庫漿主要呈柱狀涌出并開始向兩邊擴(kuò)散，此時(shí)潰口庫漿流速較快。2)隨著淹沒面積的形狀和趨勢可以看出，由于潰口靠近下游區(qū)域左側(cè)，尾礦漿從壩址處決堤后先接觸左側(cè)的封閉邊界，而后慢慢向下方和右側(cè)擴(kuò)散，但主要沿著左側(cè)邊界朝著潰口正對方向逐漸進(jìn)入下游區(qū)域，大約在50 s內(nèi)到達(dá)模擬下邊界，與此同時(shí)計(jì)算區(qū)域右半部分的淹沒面積顯著增大。3)結(jié)合模擬可以發(fā)現(xiàn)，自50 s開始，涌入右半部分的泥漿開始加速蔓延擴(kuò)散，然后泥漿從整個(gè)下游邊界出口流出。由于庫內(nèi)泥漿下泄量巨大，潰壩發(fā)生后50 s時(shí)，庫漿就能到達(dá)潰口下流500 m左右，這表明僅50 s庫漿就能縱穿整個(gè)下游區(qū)域。4)到80 s時(shí)，整個(gè)下游區(qū)域幾乎已經(jīng)完全被淹沒。

4 結(jié)論

1)依據(jù)計(jì)算得出的參數(shù)建立實(shí)體三維模型，使用ANSYS FLUENT進(jìn)行數(shù)值模擬，表明理論計(jì)算模型可以為潰壩影響預(yù)測提供相對可靠的分析方法并相互驗(yàn)證。

2)通過潰壩模擬及潰壩影響范圍的研究等為尾礦庫潰壩事故的早期預(yù)警人員疏散等提供極其重要的相關(guān)參數(shù)及信息。

3)本模擬中無論建立的計(jì)算模型是三維還是二維，實(shí)際的潰壩庫漿影響范圍面積是一致的，因此可以在分析模擬中，可使用二維模型替代三維模型進(jìn)行分析，能夠節(jié)省工作量。

4)尾礦庫潰壩，庫漿下泄是一個(gè)速度快，時(shí)間極短的過程，且該尾礦庫下流重點(diǎn)研究區(qū)域兩側(cè)封閉。因此，根據(jù)預(yù)測的潰壩后庫漿的影響范圍進(jìn)行必要的安全防控方案的編制和預(yù)防是必要的。

5)尾礦庫潰壩是一個(gè)復(fù)雜而影響巨大的事故，對庫漿下泄過程的計(jì)算及影響范圍的模擬分析對今后尾礦庫安全的管理也是極其重要的。