崩落法開采上覆崩落體降雨入滲及突水模擬研究

王 平 曾文旭 程愛平 許夢國 張玉山

（1.武漢科技大學資源與環(huán)境工程學院，湖北武漢430081；2.冶金礦產(chǎn)資源高效利用與造塊湖北省重點實驗室，湖北武漢430081）

強降雨期間塌陷坑內(nèi)匯集的雨水通過上覆崩落體和圍巖裂隙通道向采場入滲，滲流水導致碎石間的摩擦系數(shù)下降從而使其接觸力降低，增加了突水通道形成的可能性；雨水或碎屑物質(zhì)還可能在采場上方的較大縫隙或空區(qū)內(nèi)滯留，使局部壓力增大，同時在重力的共同作用下形成突水通道的可能性極高。例如程潮鐵礦使用無底柱分段崩落法開采多年來，在地表形成約93.4 hm2的塌陷區(qū)，在強降雨時共引發(fā)過4次大規(guī)模井下突水突泥事故，水裹挾大量泥土和礦巖碎屑從采場進路涌出，不僅使大量巷道和生產(chǎn)設備毀壞，還造成人員傷亡，給礦山造成巨大損失。

許多學者對井下突水突泥進行了相關研究，水從源頭到達井下引發(fā)突水突泥是一個復雜的過程。文獻［1—2］研究了具體礦山突水災害的形成機理，認為井下突水的形成需要具備流動通道和充足的水源這2個必要條件；文獻［3—4］通過數(shù)值模擬和相似實驗研究了由上覆崩落體導致采場突水突泥時關鍵部位的突變過程；文獻［5］進行了地表黃土層降雨入滲機理研究；文獻［6］建立了崩落法礦山井下采場突水突泥形成影響因素的指標體系。

綜上所述，對于降水在上覆崩落體中滲流的研究還較少，本研究結合松散礦巖的滲流特性及無底柱分段崩落法礦山采場突水機理，構建數(shù)值計算模型，分析不同工況下上覆崩落體中滲流速度和壓力分布狀況，研究降水在上覆崩落體中的入滲規(guī)律，對于預測和預防采場突水突泥具有重要意義。

1 入滲突水機制

上覆巖層的崩落方式有自然崩落或通過人工爆破崩落，崩落過程中形成的巖塊塊度較大，隨著井下后續(xù)的開采和出礦，崩落范圍越來越大，崩落范圍內(nèi)的巖塊也逐漸下降，致使上覆崩落體得以充分松散，可以將其看成是一種具有較強透水性的多孔介質(zhì)。流體在上覆崩落體內(nèi)滲流，具有流速快、雷諾數(shù)高等特點，流體在其中不僅要受到滲透壓力和孔隙水壓力的作用，還需要考慮黏性流體的剪切應力作用，這是以黏性力為主導、適合滲流速度較慢的Darcy定律所無法描述的［7］。Darcy方程描述的是牛頓流體在某種多孔介質(zhì)中滲流，滲流速度與水勢梯度正相關的滲流基本定律，該方法適應于地下水在含水層中的滲流計算；Navier-Stokes方程是不可壓縮黏性流體基于牛頓第二定律的表達式，它是在動能、勢能和靜壓能平衡同時考慮流體在重力、壓力和黏性阻力作用下的運動規(guī)律，該方程適用于管道內(nèi)流體的計算。Darcy方程和Navier-Stokes方程顯然不適合上覆崩落體內(nèi)的滲流計算。H C Brinkman在Darcy方程的基礎上考慮Navier-Stokes方程中的黏性流體剪切應力項，提出了Brinkman方程，此方程是基于牛頓第二定律，描述了在剪切應力和滲透壓力作用下快速滲流體的運動規(guī)律，對流體在礦巖松散體內(nèi)的滲流提供了較好的表達［8］。本研究將采用此方程來描述上覆崩落體內(nèi)的滲流。

Darcy定律是對飽和砂土中滲流流速和壓力之間線性關系的描述，通過飽和砂的滲流實驗得出：

式中，V為流體流速，m/s；μ為動黏系數(shù)，Pa·s；K為滲透率，m2；p為流體壓力，Pa；Z為位置高度，m；ρ為流體密度，kg/m3；?p為廣義壓力，Pa。

根據(jù)牛頓第二定律，流體微元所受重力與壓力之和等于其加速度和質(zhì)量之乘積，這是對歐拉方程的描述，Navier-Stokes方程是在歐拉方程的基礎上考慮了黏性力，式（2）為Navier-Stokes方程：

在Darcy方程的基礎上考慮Navier-Stokes方程中的黏性剪切應力項，其中壓力?p為壓力項和重力項合并后的廣義壓力，得到Brinkman方程：

針對巖體破碎帶，文獻［9］中對斷層活化導水進行模擬后發(fā)現(xiàn)發(fā)生突水時斷層導水性明顯增強，文獻［10］引入滲透率突跳系數(shù)概念，認為滲流模型中單元破壞后滲透系數(shù)在原來的基礎上會增大若干倍；文獻［11］中對峰后破碎巖石采用立方定律和strain partition技術建立滲流—體應變關系方程和數(shù)值模型，該模型能夠定量描述峰后破碎巖石滲透性的急劇增大，即破碎巖石散體的突水也符合滲透率突跳系數(shù)的概念。這里將引用滲透率突跳系數(shù)概念，用滲透率急劇增大來模擬大突水通道的形成。

2 入滲突水物理模型

根據(jù)無底柱分段崩落法采場結構及上覆崩落體形態(tài)，簡化后建立如圖1所示的降雨入滲及突水模型，程潮鐵礦塌陷坑內(nèi)黃土層長時間受塌陷擾動和雨水沖刷，不僅嚴重流失而且存在許多裂隙通道，所以該模型忽略黃土層的存在；模型主要由上覆崩落體和崩落法采場進路巷道構成，其中分段高度取8 m，采場進路巷道高度取3 m；上覆崩落體整體下移時，不同粒徑散體的移動速度不一致，細小顆粒相對粗大顆粒移動速度較快，發(fā)生上覆巖層的自然分級現(xiàn)象［12］，其中在礦體上方并與礦體接觸的覆蓋巖層一般由人工爆破放頂形成，它的主要作用是崩礦時形成擠壓爆破條件和對冒落圍巖起到緩沖，覆蓋巖層的形成方式和整個上覆崩落體的分級現(xiàn)象使這一區(qū)域相對上覆圍巖自然崩落區(qū)域具有散體塊度小和較密實的特點，同時滲透率也相對較小，模型中將上覆崩落體劃分為覆蓋巖層區(qū)域和自然崩落區(qū)域，覆蓋巖層的厚度取約等于2個分段高度。

圖1中水從上邊界入口開始向上覆崩落體內(nèi)入滲，在采場流出進入回采巷道或流向下一分層。模擬中不考慮兩側其他水源對研究區(qū)域的影響，邊界為滑移屬性的隔水邊界，與流體間不存在黏滯效應。

3 入滲及突水模擬

COMSOL Multiphysics系統(tǒng)也稱用FEMLAB系統(tǒng)，是以有限元為基礎，通過求解偏微分方程或偏微分方程組來實現(xiàn)對科學和工程領域的各種物理過程的模擬。構建計算模型，運用COMSOL Multiphysics系統(tǒng)中自由和多孔介質(zhì)流動模塊，對場變量不隨時間變化的穩(wěn)態(tài)滲流及突水進行模擬。

3.1 計算模型及模擬過程

根據(jù)圖1所示的物理模型，在COMSOL數(shù)值軟件中建立上覆崩落體滲流及突水計算模型（圖2）。模型長×高約等于80 m×100 m，頂部邊界為進水邊界。模擬時入口根據(jù)不同工況可以設置為壓力邊界條件和流量邊界條件，雨水在靜水壓力和重力的驅動下從地表塌陷坑開始依據(jù)Brinkman滲流方程向上覆崩落體內(nèi)入滲，出口為自由邊界，其壓力值設為0 Pa，位于模型下方，分別對應回采巷道和下個分層。

程潮鐵礦礦區(qū)某天降雨強度達到63.7 mm，塌陷區(qū)降雨入滲當量速度（單位面積的入滲水量）為1.63×10-6m/s，根據(jù)參數(shù)反演法得到程潮鐵礦上覆崩落體的平均滲透系數(shù)為8.03×10-4m/s（即平均滲透率為1.07×10-10m2），自然崩落區(qū)和覆蓋巖層區(qū)滲透率根據(jù)不同孔隙率來選取，具體參數(shù)見表1所示［13］。假設流體的密度和黏度都是常數(shù)；水流從采場出口流出后在回采進路內(nèi)自由流動；模型中流體在各區(qū)域邊界兩側的壓力和流速相等，即流速和壓力保持連續(xù)性。

首先模擬發(fā)生突水前降雨入滲階段，入口設置為恒定流量邊界，保持流速當量為1.63×10-6m/s，在此模型入口相當于約3.74 kPa的靜水壓力，研究突水前上覆崩落體內(nèi)滲流流速和壓力分布情況；第二步假設覆蓋巖層之上的自然崩落巖體中存在不規(guī)則空隙，黏性流體容易在其中滯留使壓力增大，在出礦和重力作用下從采場爆堆涌出，形成突水大通道，所以把覆蓋巖層上界面到出口這部分區(qū)域定為突水通道的形成區(qū)域，引用滲透率突跳系數(shù)概念，將此區(qū)域的滲透率提高10倍，來模擬突水通道形成后階段，入口邊界條件分別保持恒定流量和恒定水壓，分別研究采場發(fā)生突水后入滲流量不變時和保持初始水壓不變時滲流場流速和壓力的分布狀況。

3.2 結果分析

模擬結果如圖3～圖6所示，圖3和圖4分別是工況一中流體在上覆崩落體內(nèi)的流速和壓力分布云圖，圖5和圖6分別是地表塌陷坑到井下采場（沿圖3中折線A-B-C）的流速和壓力分布曲線。其中工況一是突水之前降雨入滲模擬，入口保持流量恒定；工況二表示入口邊界條件依舊為流量恒定時的突水工況；工況三表示維持入口邊界條件為恒定水壓時的突水工況。

從圖5中可以看出，工況一中滲流場流速在接近出口時發(fā)生大幅度變化，工況二、三中在突水通道區(qū)域開始發(fā)生大幅度變化，在出口方向上呈指數(shù)型增長，在此之前速度保持小幅度增長；工況一中流速突變之前平均值約為2×10-6m/s，出口處平均流速約1.5×10-5m/s，這是由于滲流面積總體縮小造成的；滲流壓力表現(xiàn)為自然崩落區(qū)域卸壓較慢，而靠近出口的卸壓速度快速提高，這是由于滲流場下半?yún)^(qū)域的滲透率較上半部分低，對上半?yún)^(qū)域流體的卸壓起到阻礙作用，而靠近出口時由于流體速度極大提高，使壓力得以迅速釋放。

工況二中采場上方突水通道打開后，保持入口邊界補給水量不變，入口邊界水壓從3.74 kPa降到了2.3 kPa；滲流速度和工況一相比，突水通道打開使流速增速加大，在采場出口處平均流速達到2.5×10-5m/s，這表明采場上方突水通道打開之后會起到大幅度的卸壓作用，流速也會極大提高，滲流場中流向突水采場的水量會增大。

工況三滲流場入口為恒定壓力的邊界條件，突水通道打開之后，對自然崩落區(qū)域的滲流失去阻礙作用，而且此時流體已獲得較快的速度，使上半?yún)^(qū)域流體壓力得以迅速釋放（如圖6）；滲流速度變化和工況二相似，但每個階段都明顯增大，出口處流速最高達到6×10-5m/s，約是工況一的4倍；這說明當突水水源保持恒定高水壓力時，一旦發(fā)生突水災害后所產(chǎn)生的破壞力也越大。

4 礦山記錄數(shù)據(jù)對比

程潮鐵礦2001年進行強含水帶疏干作業(yè)之后，地下水對礦山開采的影響很微弱，礦區(qū)在2005年之前平均年降雨量1 218.73 mm，4—7月雨季時日降雨量30～216.6 mm，礦區(qū)地形特點有利于周邊地表徑流向陷落區(qū)匯集，大氣降雨和地表徑流沿陷落區(qū)或移動區(qū)進入地下是開采期間最主要的充水因素。2005年時統(tǒng)計近幾年井下雨季平均涌水量，其中西區(qū)平均涌水量4 000 m3/d，當時西區(qū)在-307.5 m和-325 m水平同時采礦，2個分層共有約120個回采巷道，計算得單個回采巷道日平均涌水流量約33.3 m3/d；工況一為日降雨強度達到63.7 mm時上覆崩落體內(nèi)雨水入滲模擬，采場滲流出口處平均流速約1.5×10-5m/s，掌子面爆堆豎直投影面積約13.35 m2，計算得單個回采巷道涌水量為17.3 m3/d。

模擬結果與實際數(shù)據(jù)有一定偏差，主要是源于模型尺寸、邊界條件、參數(shù)分布等與現(xiàn)場存在差別，另外也沒有考慮崩落塌陷區(qū)之外的圍巖裂隙滲流。雖然有不可避免的誤差存在，但模擬結果與實際記錄數(shù)據(jù)在同一數(shù)量等級，模擬計算得到的滲流流速和壓力變化規(guī)律與實際情況符合，今后可以考慮更多條件，建立更完善的立體模型進行模擬。

5 結論

（1）采用Brinkman方程比較適合上覆崩落體內(nèi)滲流的表征，根據(jù)上覆礦巖散體滲透率上大下小的特點，水向下滲流過程中，卸壓速度總體上呈先慢后快的趨勢，崩落體上部分滲流流體壓力較高，是促發(fā)突水通道形成的一個重要因素；

（2）地表塌陷坑匯水深度越大，井下采場一旦發(fā)生突水所造成的損失程度越大，甚至形成井下泥石流，在降雨天氣做好塌陷坑內(nèi)的防排水工作能有效減少井下突水事故發(fā)生幾率或降低事故危害。