濱海礦山深部開采水力調查及穩定性分析

陳倩男 郭奇峰 萬思達 翟 濟 杜偉嘉 李 釗

(1.北京科技大學土木與環境工程學院，北京100083;2.金屬礦山高效開采與安全教育部重點實驗室，北京100083)

海洋蘊藏著豐富的海底礦產資源，包括海底表層沉積物和海底巖層中貯藏的礦產。目前開發利用最為廣泛的是海洋石油和天然氣，在整個海底礦產資源已開發產值中的占比超過90%。陸地礦產資源經過多年開采，淺部資源枯竭，國內外很多礦山進入深部開采階段［1-2］，并且諸多近海礦山開始進行海底礦產資源的開采［3-4］，其中日本的海底煤礦開采量達到總產量的30%，海底開采將成為未來礦產資源開發的發展方向。

海底基巖巖體存在大量節理裂隙，水和應力通過節理裂隙進行耦合作用，對海底資源開采過程的安全具有重要影響。因此，在海底資源開采必須深入調查水、節理裂隙和地應力的分布特征，研究水和應力的耦合作用機理。諸多學者開展了大量針對裂隙巖體中水和應力等因素耦合作用的研究，主要工作集中在理論模型研究，試驗分析和數值模擬應用3 個方面。王啟平［5］、馬海濤等［6］基于理論研究和試驗分析，采用數值模擬軟件，對強夯地基和露天采礦邊坡穩定進行了水力耦合分析，驗證了水力耦合在工程應用中的必要性和適用性。本課題以濱海礦山三山島金礦深井開采為背景開展水力研究，為水力耦合作用下的海底資源安全開采提供指導和保障。

1 水和地應力調查

1.1 地下水分布

礦體及上覆頂板中的基巖裂隙含水是三山島金礦深部開采區域直接充水層。根據含水層與隔水控礦構造F1的相對位置關系，直接影響礦坑充水的基巖裂隙含水帶可劃分為F3斷裂含水帶、F1上盤裂隙含水帶和F1下盤構造裂隙含水帶。

上述含水帶中，F1下盤構造裂隙含水帶分布廣泛，直接影響深部開采和井巷工程施工。根據構造發育特征和補給條件的不同，劃分為F1下盤構造裂隙富水區和F1下盤構造裂隙含水區。F1下盤構造裂隙富水區內，北西向導水裂隙發育，且裂隙規模和延伸較大，含水性較強;北西端有海水補給，東南方向與F3斷裂構造含水帶溝通，有區域鹵水補給。F1下盤構造裂隙含水區逐漸遠離F3，富水性較弱，北西向含水裂隙密度和規模逐漸變小，北西部有海水越流補給，并有少量的第四系孔隙水滲透補給。

1.2 地應力場分布

礦山深部開采面臨“三高一擾動”的問題，其中地應力是一個重要影響因素，其關系到井巷工程和開采工程的布置，直接影響深部開采的安全性。三山島金礦深部開采中，高應力顯現明顯，局部地區甚至出現巖爆現象。為保證現有開采工程的安全，為更深層次的開拓和開采設計提供安全指導和保障，用套孔應力解除法在深部－500 m 以下6 個水平開展9 個測點的地應力測量工作，結合2002 年在淺部3 水平4個測點的地應力測量結果［7］，共獲取三山島金礦13個測點的地應力測量數據，通過回歸分析，繪制三山島金礦深部主應力隨深度變化曲線如圖1 所示。

圖1 主應力隨深度變化曲線Fig.1 Variation curve of principal stresses with depth

礦區最大水平主應力、最小水平主應力和垂直主應力的回歸方程為

2 巖體滲透性

2.1 節理裂隙分布

對中段運輸平巷、采場、聯絡巷等進行現場節理裂隙調查，記錄節理傾向、傾角、走向、節理間距、充填物及滲水特征，共獲取854 組節理記錄，其傾向及傾角分布統計如表1 所示，繪制節理傾向玫瑰花圖和走向玫瑰花圖如圖2 和圖3 所示。

表1 節理分布特征Table 1 Joints distribution

圖2 節理傾向玫瑰花圖Fig.2 Rose diagram of joints＇strike

圖3 節理走向玫瑰花圖Fig.3 Rose diagram of joints＇dip

統計結果顯示，深部節理發育具有5 組優勢傾向。第1 組優勢傾向為20° ～30°，平均傾角62°;第2組優勢傾向為110° ～150°，平均傾角57°; 第3 組優勢傾向為190° ～200°，平均傾角66°; 第4 組優勢傾向為270° ～280°，平均傾角54°; 第5 組優勢傾向為290° ～320°，平均傾角55°。走向具有5°、20° ～40°、270° ～280°3 組優勢方位，傾角多分布在40° ～90°。深部節理裂隙的傾角均偏大，平均傾角超過55°。

2.2 巖體滲透特性

2.2.1 幾何測量計算

根據Snow［8］提出的滲透張量幾何計算方法，假定每個裂隙內沒有填充物，均無限延伸且相互連通，但各裂隙的滲流互不影響。以地理N 向為x 軸，地理E 向為y 軸，則各裂隙滲透系數張量Ki的表達式為

式中，ρ 為地下水密度，kg/m3; g 為重力加速度，m/s2; bi為第i 組裂隙寬度，m; μ 為地下水動力黏滯系數，kg/(m2·s); Kxxi= 1 － sin2αicos2βi; Kxyi= Kyxi=－sin2αicosβisinβi; Kzxi= Kxzi=－sinαicosαicosβi; Kyyi=1 －sin2αisin2βi; Kzzi= sin2αi;αi為第i 組裂隙傾向，0≤αi≤360°; βi為第i 組裂隙傾角，0≤βi≤90°。

作為二維張量，滲透系數張量的特征方程可表達為

滲透張量的3 個主值是上述特征方程( 5) 的3個特征根A1、A2和A3，則綜合滲透系數K［9-11］為

選取多個有代表性的上盤、礦體、下盤巖體典型表面裂隙特征調查窗口，詳細統計窗口內的節理裂隙產狀和分布特征。依據上述滲透系數幾何計算原理，采用MATLAB 軟件計算測量區域內上盤、礦體、下盤的綜合滲透系數，結果如表3 所示。

表2 滲透系數Table 3 Permeability coefficients m/s

2.2.2 滲透系數修正

Snow 幾何量測法在進行滲透系數計算時，假定所有節理裂隙在裂隙網絡交叉空間相互連通，但現實巖體中節理的延伸范圍有限，各組節理并未呈現相互交切的形態。為了更準確地掌握深部礦巖滲透特性，依據節理裂隙產狀、間距等分布特征，采用Monte －Carlo 隨機模擬法［12-13］對上盤、礦體、下盤3 類區域節理裂隙的滲透系數進行修正。模擬結果顯示，上盤巖體的連通率為84.76%，礦體連通率為88.49%，下盤巖體連通率為85.13%。張有天［14］提出采用連通系數即連通率對理想條件下計算出的滲透系數進行修正:

式中，k 為修正后的滲透系數; k0為計算值; η 為連通率。將幾何計算得出的滲透系數乘以節理裂隙連通率即獲得修正后的上盤滲透系數為6.236 1 ×10－7m/s，礦體滲透系數為8.349 2 ×10－7m/s，下盤滲透系數為6.743 3 ×10－7m/s。

3 深部開采水力耦合分析

3.1 計算模型

采用FLAC3D數值模擬軟件建立深部開采水力耦合數值模型，模型標高自－430 m 至－630 m。坐標系以礦體走向為x 軸，礦體厚度方向為y 軸，鉛垂方向為z 軸。計算模型在x 方向上的長度為400 m，y方向上的長度為300 m，z 方向上的高度為200 m，共劃分58 310 個單元，63 900 個節點。對中段F3斷層北側礦體進行房柱交替式盤區上向分層充填模擬開采，盤區長度選取80 m，盤區沿礦體走向共布置8 個采場進行礦房礦柱交替上升式開采，每個采場寬10 m，采場長度為礦體厚度，中段開采共分6 步進行，如圖4 所示的1 ～6。

圖4 開采步驟Fig.4 Mining steps

3.2 邊界條件

根據三山島金礦深部詳查報告中地下水位長期觀察和水位降低值的研究結果［15］，在模型上部和底部施加水頭邊界。依據深部地應力場分布規律計算深部最大主應力、最小主應力和垂直主應力值，對模型施加梯度應力邊界。

3.3 力學參數

模型劃分上盤、礦體、下盤以及F3斷層4 大組，F3斷層傾角近90°，垂直礦體走向，礦體傾角45°，開采后直接充填。依據物理力學試驗數據和滲透性計算數據，選取水力耦合數值模擬采用的力學參數和滲透系數如表3 所示。

表3 物理力學參數Table 3 Physico-mechanical parameters

3.4 模擬結果分析

3.4.1 采場位移變化

盤區內最大豎向位移均發生在新開采區域的中心位置，豎向位移朝開采區域集中，由圖5 所示開采至第二步的采場礦體和充填體豎向位移等值線圖可見，各開采步驟豎向位移出現明顯分區，即1 到4 號采場為1 組，5 ～8 號為1 組，在4 和5 采場交界帶豎向位移值相對較小，呈現出對稱性。開采過程中，采場頂板沉降逐漸增大，變化較為穩定，無明顯突變發生，在各采場開采接近頂柱時，頂板沉降達到最大值。

圖5 開采第二步的豎向位移等值線Fig.5 Vertical displacement isoline in 2nd mining step

為研究開采接近頂柱時的采場安全性，在中段每個采場上方對應的頂柱設置監測點共計8 個。各監測點在中段開采過程中的豎向位移因房柱交替開采呈現階段跳躍式增長，當開采高度超過中段高度一半之后，其沉降增幅明顯，開采結束后2 號監測點的最大沉降達到9.276 cm，存在安全隱患。塑性區云圖也顯示各采場回采接近頂柱時出現塑性區貫通，應加大支護密度和管理力度，關注采場頂板的安全穩定性。

3.4.2 采場應力變化

在6 個開采步驟過程中最大主應力值分別為50.073、60.024、52.398、67.959，52.603，43.048 MPa，具有代表性的開采第2 步的最大主應力等值線如圖6 所示。模擬結果顯示在整個盤區開采過程中，二步采引起的應力集中要大于一步采，中間和兩翼的最大主應力值較高，二步采階段采場穩定性相對較差。最大主應力的極值主要出現在1 號采場圍巖部位，5 號采場頂板應力集中顯現也較為明顯，開采過程中需要對上述2 個區域進行重點支護和加固，確保生產安全。

圖6 開采第2 步的最大主應力等值線Fig.6 Maximum principal stress isoline in 2nd mining step

各個開采步驟過程中均存在拉應力區，極值分別為2.656、1.532、1.484、1.515、1.362、1.225 MPa，拉應力主要集中在一步采和二步采的交替部位，其分布區域具有一定的時空變化規律，即一步采階段主要存在于2 和6 號采場，二步采階段主要存在于3 和7 號采場。針對上述拉應力區的分布，回采過程中應加大支護力度，保證充填體的強度和充填接頂效果。

4 結 語

(1) 針對海底資源深井開采面臨的富水和高應力問題，對三山島金礦地下水和地應力進行調查和測量，掌握了地下水的分布，建立了礦區地應力場模型，為礦山開采和優化提供了基礎資料。

(2) 開展深部采場和巷道圍巖節理裂隙調查，統計深部巖體節理裂隙分布特征，采用幾何測量計算和結構面隨機模擬修正的方法研究了深部礦巖體的滲透特性。

(3) 建立水力耦合數值模型，對深部礦體房柱交替式上向分層充填開采進行動態模擬。分析了分步開采過程中采場位移、應力隨開采推進的時空變化規律，對應力集中、拉應力破壞等危險區域進行了時空分析，確定了開采過程中需要重點支護的區域和時間，為海底資源水力耦合作用下的安全開采提供了保障。

［1］ 蔡美峰. 金屬礦山采礦設計優化與地壓控制——理論與實踐［M］北京:科學出版社，2001.

Cai Meifeng.Optimization of Mining Design and Control of Ground Pressure in Metal Mines: Theory and Practice［M］. Beijing: Science Press，2001.

［2］ 王春來，吳愛祥，劉曉輝，等. 深井開采微震活動容量維Df 變化特征［J］.北京科技大學學報，2010，32(11) :1379-1382.

Wang Chunlai，Wu Aixiang，Liu Xiaohui，et al.Variation characteristics of capacity dimension Df with microseismicity in deep mining［J］.Journal of University of Science and Technology Beijing，2010，32(11) :1379-1382.

［3］ 喬衛國，呂言新，魏烈昌，等. 三山島金礦海底高效開采工業試驗［J］.金屬礦山，2011(7) :43-46.

Qiao Weiguo，Lu Yanxin，Wei Liechang，et al.Industrial test of highefficient submarine exploitation in Sanshandao Gold Mine［J］.Metal Mine，2011(7) :43-46.

［4］ 彭 康，李夕兵，趙國彥，等. 巖層微擾上向分層充填法在海底開采中的應用［J］.金屬礦山，2010(4) :50-52.

Peng Kang，Li Xibing，Zhao Guoyan，et al.Application of the upward slicing and filling with perturbation in undersea exploration［J］.Metal Mine，2010(4) :50-52.

［5］ 王啟平，謝能剛，史小路.基于強夯大變形的地基流固動力耦合分析［J］.北京科技大學學報，2004，26(4) :345-348.

Wang Qiping，Xie Nenggang，Shi Xiaolu. Fluid-solid dynamic coupling analysis of foundation soil during large deformation by dynamic consolidation［J］. Journal of University of Science and Technology Beijing，2004，26(4) :345-348.

［6］ 馬海濤，李迎佳，王云海.港里鐵礦水下開采流固耦合數值模擬與風險分析［J］.金屬礦山，2014(2) :36-39.

Ma Haitao，Li Yingjia，Wang Yunhai.Fluid-solid coupling numerical simulation and risk analysis of underwater mining of Gangli Iron Mine［J］.Metal Mine，2014(2) :36-39.

［7］ 喬 蘭，歐陽振華，來興平，等. 三山島金礦采空區地應力測量及其結果分析［J］.北京科技大學學報，2004，26(6) :569-571.

Qiao Lan，Ouyang Zhenhua，Lai Xingping，et al.In-situ stress measurement and its result analysis in Sanshandao Gold Mine of China［J］.Journal of University of Science and Technology Beijing，2004，26(6) :569-571.

［8］ Snow D T. Anisotropic permeability of fractured media［J］. Water Resources Research，1969，5(6) :1273-1289.

［9］ 仵彥卿.巖體水力學基礎( 七) ——巖體水力學參數的確定方法［J］.水文地質工程地質，1998，25(2) :42-48.

Chu Yanqing. Basic of rock mass hydraulics ( 7) : Determination methods of hydraulic parameters of rock mass［J］.Hydrogeology and Engineering Geology，1998，25(2) :42-48.

［10］ 黃克智，薛明德，陸明萬. 張量分析［M］. 北京: 清華大學出版社，2005.

Huang Kezhi，Xue Mingde，Lu Mingwan.Tensor Analysis［M］.Beijing:Tsinghua University Press，2005.

［11］ 周志芳，王錦國.裂隙介質水動力學［M］.北京: 中國水利水電出版社，2004.

Zhou Zhifang，Wang Jinguo.Theory on Dynamics of Fluids in Fractured Medium［M］.Beijing:China Water ＆ Power Press，2004.

［12］ 方 濤，柴軍瑞，胡海浪，等.Monte Carlo 方法在巖體裂隙結構面模擬中的應用［J］.露天采礦技術，2007(1) :7-9.

Fang Tao，Chai Junrui，Hu Hailang，et al. Application of Monte Carlo method to simulating the fracture interconnectivity of rock mass［J］.Opencast Mining Technology，2007(1) :7-9.

［13］ 杜景燦，陳祖煜，彌宏亮，等. 三維條件下應用遺傳算法與Monte－Carlo 法確定節理巖體的綜合抗剪強度［J］. 巖石力學與工程學報2004，23(13) :2157-2163.

Du Jingcan，Chen Zuyi，Mi Hongliang，et al.Determination of comprehensive shear strength for jointed rock massed in 3D condition using genetic algorithm and Monte-Carlo method［J］. Chinese Journal of Rock Mechanics and Engineering，2004，23(13) :2157-2163.

［14］ 張有天.巖石水力學與工程［M］.北京: 中國水利水電出版社，2005.

Zhang Youtian. Rock Hydraulics and Engineering［M］. Beijing:China Water ＆ Power Press，2005.

［15］ 三山島金礦. 山東省萊州市三山島礦區深部金礦詳查報告［R］.濟南:山東黃金集團有限公司，2007.

Sanshandao Gold Mine.Detailed Deep Survey Report of Sanshandao Gold Mine in Laizhou，Shandong Province［R］. Ji＇nan: Shandong Gold Group Co.，Ltd.，2007.