徐 曲,王玉平,曹華兵,杜建華,馮 杰
(1.湖北省應急救援中心,湖北 武漢 430070;2.武漢盛世金安安全科技有限公司,湖北 武漢430080;3.武鋼資源集團程潮礦業有限公司,湖北 鄂州 436050;4.中南財經政法大學信息與安全工程學院,湖北 武漢 430073)
地下礦產資源采出誘發的巖層沉陷、地表開裂是采掘區域巖體初始應力平衡被打破后,應力的轉移、釋放、再平衡的過程,隨開采范圍、深度加大圍巖發生漸進式破壞、變形,由采掘區域向地表逐漸發展,進而導致地表、井巷工程、地表構筑物的破壞,其發展規律受礦體空間形態、采掘工藝、工程與水文地質條件等多種因素的影響。巖層移動、破壞從廣義上來講屬于開采沉陷學的范疇,按照礦業行業劃分,煤礦開采領域相關的研究成果較豐富,理論性強,這主要是由于煤礦多屬于沉積礦層,礦體、巖層多呈層狀分布,巖層變形、破壞規律性較為一致。而金屬礦床多在地殼活動過程中巖漿沿地殼斷裂、薄弱地帶的侵入、噴出形成,其具有空間形態多變、圍巖種類繁多、強度多變、分布無序、巖體各級結構面錯綜復雜等特點。因此,概率積分法、砌體梁經典理論、關鍵層理論等方法均不能很好地應用于金屬礦山開采引起的地表沉陷、巖層破壞的預測。
金屬礦山開采引起的巖層移動、破壞規律是礦巖體空間分布、強度、結構面特征以及采掘工藝等一系列因素綜合影響的結果,對其研究是巨大的系統工程問題。由于金屬礦山開采引起的地表沉陷問題的復雜性,導致我國早先礦山設計在巖層移動范圍的測算上出現誤差,進而導致了眾多礦山關鍵工程破壞,例如金川集團二礦區14行風井破壞事故,程潮鐵礦地表運礦巷道、措施井開裂等。
針對金屬礦山開采過程中地表沉陷問題的研究方法主要包括經典理論分析、數值模擬、物理模型試驗等方法,其中三維數值模擬計算方法對于適合于復雜空間形態多因素協同作用誘發的地表沉陷、巖層破壞問題的研究最為有效。因此,本文以程潮鐵礦無底柱分段崩落法開采為背景,通過構建涵蓋礦體、圍巖體分布、采掘工藝、初始應力條件等信息的高精度三維空間模型,并通過準確的巖體力學實驗,選取適合的本構模型,且輔以地表沉降監測措施,可較為完善地解決復雜條件金屬礦山開采誘發的地表沉陷的測算問題。
本研究采用三維有限差分軟件(FLAC3D),根據礦巖體空間分布特征、強度特征、采掘順序建立數值計算模型,按照時間順序真實地模擬礦體開采,從而得出研究礦區內的巖層移動和破壞規律。
程潮鐵礦由多個分礦體組成,分礦體呈疊瓦狀分布,多呈透鏡狀、扁豆狀、脈狀、囊狀分布,礦體厚度變異大、夾巖復雜,礦體整體走向為北西西,傾向為南西,傾角為20°~30°,最大分礦體長度超過1 000 m,厚度超過100 m。
程潮鐵礦實際開采過程中分為東西兩個作業區域,西區礦體埋深深、厚度較大,東區礦體埋深淺、厚度相對較小。本文研究的地表構筑物位于西區,受東區礦體開采的影響較小,因此為了提升計算效率,本次建立的三維數值模型僅考慮了西區礦體和圍巖,主要針對西區開采對地表臨近構筑物的影響進行研究。根據地質勘探數據生成的程潮鐵礦區礦體空間分布模型,見圖1。

圖1 程潮鐵礦區礦體實體空間分布模型圖
根據西區礦體賦存位置以及地表投影范圍,確定礦體覆蓋區域尺寸如下:走向方向為1 200 m,傾向方向為800 m,礦體最大埋深接近1 200 m。本文中僅對-500 m水平以上礦體開采進行研究,故根據-500 m水平以上礦體的三維幾何外輪廓,按照移動角40°(保守取值)計算模型邊界。數值模型取礦體走向為X
軸,傾向為Y
軸,垂直方向為Z
軸,最終確定數值模型三維尺寸為:X
方向3 000 m,Y
方向2 800 m,Z
方向1 200 m,見圖2。
圖2 程潮鐵礦區數值計算模型示意圖
數值模型單元尺寸的劃分以及與地質模型、開采工藝的對應關系分析如下。
1.1.1 單元尺寸劃分
(1) 精細網格劃分區域:該區域主要為礦石賦存區域,其形狀為平行六面體,依據礦體X
方向覆蓋范圍確定左、右邊界,根據礦體Z
方向覆蓋范圍確定上、下邊界,再根據礦體傾角與礦體在Y
方向的投影邊界確定前、后邊界。該區域Z
方向最小單元尺寸根據礦山實際生產中的開采分段高度(17.5 m)確定,礦體區域采用六面體單元劃分網格,單元網格X
、Y
、Z
三向尺寸比例控制在1∶2范圍內,即六面體單元最長邊長度不大于最短邊長度的2倍,單元網格X
、Y
、Z
三向尺寸依次為33 m、28 m、17.5 m。(2) 外部圍巖區域:整個模型扣除上述精細網格劃分區域之外均為外部圍巖區域。該區域按照六面體單元劃分網格,為了提高計算速度,單元網格尺寸適當放大,單元網格X
、Y
、Z
三向尺寸依次為50 m、50 m、17.5 m,模型共劃分357 840個單元、378 288個結點。1.1.2 數值模型與地質模型的對應關系
數值模型建立后,首先導出每個單元的節點三維坐標,進而求取單元格質心坐標;然后將每個單元格的質心坐標導入地質模型,從而獲得該區域礦巖類型,并以此為依據進行單元分組及屬性賦值。
1.1.3 開采工藝模擬
依據開采順序,由上至下模擬礦體開采,每次開采整個分段高度(17.5 m)的礦體,礦石采出后頂部圍巖冒落區域按照橢球體計算,冒落區域最大寬度l
與開采范圍的關系依據下面公式確定:L
<50 m時l
=(1.
1~1.
5)L
(1)
L
>50 m時l
=(1.
5~2.
0)L
(2)
式中:L
為礦體長度(m);l
為冒落區域最大寬度(m)。數值模型建模過程中已將模型邊界擴展至開采影響區域之外,所以計算采取位移邊界條件,即模型四周限制水平位移,底部限制垂直位移;模型頂面為地表,無需施加邊界條件。
本構模型選取摩爾-庫倫(Mohr-Coulomb)強度準則,巖體力學參數在巖體力學試驗的基礎上采用Hoek-Brown準則進行折減,模擬計算采用的巖體力學和變形參數見表1。

表1 程潮鐵礦西區巖體的力學和變形參數
冒落區域為散體,其強度屬性隨時間發生變化,其強度計算公式如下:
ρ
=1 600+800(1-e-1.25)(3)
E
=15+200(1-e-1.25)(4)
υ
=0.05+(1-e-1.25)(5)
式中:ρ
為材料的密度(kg/m);E
為材料的彈性模量(MPa);υ
為材料的泊松比;t
為時間(a)。程潮北斷層是區域穩定的重要影響因素,其力學參數見表2。

表2 程潮北斷層的力學參數
模型應力條件根據礦區實測施加,第一主應力σ
位于水平面內,方向與礦體走向一致;最小主應力σ
位于水平面內,方向垂直礦體走向;第二主應力σ
位于水平面垂直方向。埋深H
位置的主應力分量計算公式如下:
(6)
式中:γ
為巖層的平均容重(kN/m),取值為28 kN/m;H
為埋深(m)。針對程潮鐵礦礦體開采誘發的巖層移動,已有研究表明,礦區地表變形和地表構筑物的破壞規律與巖體力學性質、礦體空間形態、采空區范圍、開采順序密切相關,因此在計算過程中,依據實際開采范圍、開采順序、無底柱分段崩落法采礦作業流程對礦體開采誘發的地表變形和地表構筑物的破壞規律進行了數值模擬計算。
根據礦區巖體強度和初始應力條件擬合開采前的原巖應力狀態,在此基礎上按照時間順序模擬不同階段的礦石開采。首采水平為-290 m,本研究的最大開采深度至-500 m水平,在此區間內按照17.5 m標高設置開采分段,共包括13個開采階段。
S
=2.5 mm/m、水平拉伸率T
=1.5 mm/m作為地表臨界變形值用來判定地表錯動區域范圍。X
方向地表水平位移分量和Y
方向地表水平位移分量合成為地表水平位移,進而根據下面公式(7)、(8)計算地表變形誘發的豎向傾斜率S
和水平拉伸率T
。假設點1、點2兩點距離為D
,點1的地表豎向沉降量為D
,點2的地表豎向沉降量為D
,則地表變形誘發的豎向傾斜率S
的計算公式如下:
(7)
點1的地表水平位移為D
,點2的地表水平位移為D
,則地表變形誘發的水平拉伸率T
的計算公式如下:
(8)
限于篇幅,本文僅對關鍵開采階段的數據進行分析論述,-360 m水平礦體開采完畢后,該礦區地表變形的模擬計算結果如下:X
方向為水平面礦體水平走向方向,Y
方向為水平面礦體垂直走向方向,Z
方向為豎直方向;礦體開采至-360 m水平時地表X
、Y
方向位移云圖分別見圖3(a)和圖3(b),將X
、Y
方向的地表水平位移分量進行矢量合成,得到地表水平位移合成見圖3(c),圖3(d)為地表豎向沉降等值線云圖。
圖3 程潮鐵礦區地表變形云圖
通過提取不同開采階段措施井所在位置的水平位移、垂直位移,可根據公式(7)、(8)計算得到對應位置的地表變形誘發的豎向傾斜率S
和水平拉伸率T
。程潮鐵礦區措施井所在位置各開采階段地表變形數據匯總于表3,各開采階段地表變形曲線見圖4。

表3 程潮鐵礦區措施井對應位置各開采階段地表變形統計表

圖4 程潮鐵礦區措施井對應位置各開采階段地表變形曲線圖
由表3和圖4可以看出:
(1) 按照開采中段計算,-360 m水平礦體開采完成后,措施井對應位置地表變形誘發的豎向傾斜率為0.80 mm/m、水平拉伸率為0.78 mm/m,上述兩值均小于規范要求的限值,即說明開采至-360 m水平,措施井處于安全狀態。
(2) -430 m水平礦體開采后,措施井對應位置地表變形誘發的豎向傾斜率為1.18 mm/m、水平拉伸率為1.50 mm/m,水平拉伸率達到臨界值,即說明礦體開采至-430 m水平時,地表錯動區域邊界達到措施井附近。
(3) 但由于措施井為貫穿至地下采礦作業區域的圓筒結構,而巖層變形是由下部逐漸向上部圍巖擴展,即由開采區域至地表的變形呈遞減規律,在地表變形達到臨近破壞限值前,措施井井筒地下部分變形已經超過限值。通過對措施井井筒地下部位不同開采階段的地表變形分析可知,-360 m至-375 m水平礦體開采后,措施井井筒在地表200 m以下區域巖層水平拉伸率超過了1.5 mm/m限值,措施井井筒破壞區域首先在地表200 m以下范圍產生,且隨著開采深度的增加,其破壞程度加劇,破壞范圍向措施井井筒上部擴展。
(4) -500 m水平礦體開采后,措施井對應位置地表變形誘發的豎向傾斜率為1.38 mm/m、水平拉伸率為1.60 mm/m,水平拉伸率超過了臨界值,即說明開采至-500 m水平,措施井進入地表錯動區域,其破壞區域嚴重。
本文采用三維仿真計算方法對程潮鐵礦區內由于礦體地下開采引起的地表變形規律進行了研究,并對不同開采水平地表范圍和典型剖面地表移動角進行了預測,得到的結論如下:
(1) -342 m水平以上礦體開采對于措施井所在位置地表變形的影響甚微,對措施井穩定性無影響。
(2) 由于-360 m至-412 m水平礦體向下盤方向擴展,-360 m水平以下礦體開采對措施井所在位置地表變形的影響加劇。開采至-375 m水平時,措施井所在位置地表變形誘發的水平拉伸率為0.9 mm/m,尚未超過允許極限值,但由于措施井為貫通至地下采礦作業區域的井筒結構,根據模擬計算結果得出:-360 m至-375 m水平礦體開采后措施井井筒在地表200 m以下區域巖層水平拉伸率超過1.5 mm/m限值,措施井井筒破壞區域首先在地表200 m以下范圍產生,且隨著開采深度的增加,其破壞程度加劇,破壞范圍向措施井井筒上部擴展。
(3) 開采至-430 m水平時,地表變形誘發的水平拉伸率達到1.5 mm/m,措施井進入地表錯動區域范圍,措施井井筒破壞范圍擴展至地面,由于-430 m至-500 m水平礦體寬度變窄,此范圍礦體開采對于措施井穩定性的影響逐漸減弱。
(4) 措施井所在位置的地表錯動角在開采至-500 m水平時為50°左右,遠低于設計值(約60°)。