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

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

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

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

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

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

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

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

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

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

表1 尾礦庫容積計算表Table 1 Tailings pond volume calculation table

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

初期壩頂標(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，坡面上布置反濾層，可在一定程度上對上游流下來的水起引流作用，保護坡面。石壩體與岸坡的連接方式為斜面連接，這樣可以避免壩體與岸坡的連接面產(chǎn)生集中滲流現(xiàn)象，也可以減少軟弱夾層的形成。

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

2 潰壩尾礦漿流動過程

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

2.1 尾礦庫潰壩下泄總量

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

2.2 尾礦庫潰口寬度

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

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

(1)

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

2.3 尾礦庫潰口的最大流量

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

(2)

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

2.4 庫漿下泄時間

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

庫漿下泄時間計算公式見式(3)。

(3)

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

2.5 庫漿最大流速

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

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

(4)

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

(5)

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

2.6 庫漿最大沖擊力

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

(6)

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

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

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

(7)

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

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

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

3 潰壩數(shù)值模擬

由于各類型尾礦庫的地理環(huán)境差異巨大且復(fù)雜多變，庫容物質(zhì)也不盡相同，因此，想建立統(tǒng)一模型進行理論分析的可能性微乎其微?，F(xiàn)如今，隨著科技的發(fā)展，運用ANSYS FLUENT數(shù)值模擬手段能夠廣泛而可靠的適用不同庫型、不同條件、各種工況下的潰壩泥沙演進過程，進而為理論研究提供更加行之有效的參考數(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è)高，潰壩時，泥漿由于地勢的阻擋，會向左下側(cè)流動。假設(shè)潰壩事故發(fā)生時潰口在一瞬間完全潰決，潰口高度為35 m，為整個堆積子壩高度；潰口長度為整個壩軸線長，即105 m；計算區(qū)域內(nèi)共有908 170 m3下泄庫漿。潰壩計算模型及尺寸示意圖見圖1。

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

3.2 網(wǎng)格劃分

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

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

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

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

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

2)庫漿在運動過程中為不可壓，因此在求解中選擇Pressure-Based求解器；

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

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

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

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

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

邊界條件設(shè)置：

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

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

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

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

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

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

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

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

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

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

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

4 結(jié)論

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

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

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

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

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