基于數值模擬技術的井筒位置優化

吉學文，劉 輝，潘張偉，王 志，扶海鷹

(1.大興安嶺金欣礦業有限公司， 黑龍江大興安嶺市 165100;2.馳宏科技工程股份有限公司， 湖南長沙 410000)

1 工程背景

大興安嶺岔路口鉬鉛鋅多金屬礦為斑巖型礦床，主要圍巖是花崗斑巖、石英斑巖、閃長玢巖等。鉬礦體總體呈拉長的穹窿狀，主體隱伏，現已控制礦脈長度2600 m，剖面上控制寬度360～1260 m，最大寬達1760 m;礦體延深大于1500 m，垂直厚度一般200～900 m，最大為1097.48 m。設計生產規模為1650萬t/a，礦體按上下分區同時開采，采用大直徑深孔空場嗣后充填法和分段空場嗣后充填法自下而上開采，上采區開采范圍為+20～+440 m之間的礦體，上采區首采中段為140 m中段，下采區首采中段為－460 m中段，上下采區形成獨立的回風系統，1#風井為上采區的專用回風井筒，設計服務年限為20 a。

2 設計方案

按照井巷布置有關規定，若將1#回風井位置設置在全礦區地表移動界線范圍之外，距離上采區礦體約為1500 m，同時每個中段回風井石門也長達1500 m。這對于1#回風井僅服務于上采區礦體開采的回風，極大增加了基建投資及通風運營成本。為優化工程布置，根據礦體賦存狀況、采礦工藝要求及地表地形，在地表初步確定了回風井的3個位置，3個方案的回風井距離礦體300～500 m，按照采礦設計規范及建構筑物設計要求，通過FLAC3D三維數值模擬方法，對初步確定的回風井3個優選方案在不同時期(15 a，20 a)的穩定性進行驗證比較分析。

通過FLAC3D分析研究井筒位置的可行性，并與地表移動界限外井筒位置比較，確定合理可行的井筒位置。

3 數值模擬分析

有限差分法在計算機數值模擬中已被廣泛運用。采用拉格朗日法的有限差分數值模擬計算分析軟件已廣泛應用于巖土工程中，彌補了類比法和解析法難以實現大范圍礦床開采對井巷工程穩定性的分析的不足，同時也更好地考慮介質的各向異性、非均質特性及其隨時間的變化、復雜邊界條件和介質不連續性等復雜地質條件。

3.1 基本假設及模型建立

由于礦山地質結構的復雜性，影響采場穩定因素的多樣性，在模型建立和計算過程中采用以下簡化和假設:

(1)階段運輸巷道、主井、通風井等大量工程對整個礦山開采穩定性有一定的影響，但這些工程對宏觀礦區開采的穩定性影響較小，因此，本次研究不考慮這些工程的影響;

(2)在建模中忽略模擬范圍內的節理、裂隙和斷層，在模擬過程中，選擇具有代表性的3個巖層(風化帶巖組，微風化帶巖組和未風化帶巖組)建立圍巖;

(3)在模擬計算中，忽略爆破震動、地震波等動力因素及地下水對礦巖穩定性的影響。

本次模擬過程嚴格依照設計開采進度計劃，先進行風井的開挖，然后以年為單位進行各個中段的開采和充填，記錄下開采15 a和20 a對1#回風井的穩定性影響。經過折減后的礦巖物理力學參數見表1。

3.2 本構模型與邊界條件

莫爾－庫侖模型是一種彈塑性本構模型，它能反應巖土材料剪切破壞特性，并廣泛應用于巖土實際破壞分析中，剪切屈服函數和拉應力屈服函數組成了莫爾－庫侖模型破壞包絡線。

表1 礦巖物理力學參數

采用3Dmine和FLAC3D結合的方式建立數值計算模型。先以位置1進行建模，根據礦區地表等高線模型圖及整個礦體的空間形態，模型的計算域長×寬取為3200 m×2300 m，高度取－940 m標高到地表。采用3Dmine軟件對礦體實體模型進行塊處理，同時加上地表、礦體和1#回風井的約束條件，最終生產模型的FLAC3D命令流文件，在塊體大小的劃分中，由于模擬范圍大，同時開挖步驟復雜、中段高度為60 m，最終選用的最大單元格大小為50 m×50 m ×60 m，次級單元格大小為12.5 m ×12.5 m ×20 m，FLAC3D軟件讀取命令流后生成計算模型，位置1的模型共235989個單元。

模型四周采用可滾動的軸支座約束側向變形，底面采用固定支座邊界類型約束豎直方向變形，頂面為自由面。

模擬計算主要記錄3個不同風井位置條件下應力場、位移場和塑性區的變化過程，通過對比分析選取合理的風井位置。

3.3 模擬結果

以位置1進行分析，主要對模型內的位移變化、應力值變化規律和塑性區分布大小3個方面進行分析。

(1)變形規律分析。開采15 a對回風井擬選1#位置地表造成的沉降約為65 mm，開采20 a的最終沉降值約為85 mm。圖1中，地表最大位移區域在15 a后出現向右下角移動的趨勢，是由于在模擬開挖過程中，15 a后開挖的80 m和20 m中段的礦體主要集中在風井位置的右下角。距井筒腰部最近的礦體約為50 m，見圖2，由于該礦體為獨立礦體，開挖跨度小，開采對井筒的影響不大;下采區礦體開采20 a，由于中間還隔有近400 m厚度的未采動區域，由計算結果可知，此條件下下采區的采動對1#回風井的穩定性影響不大;開采15 a和20 a井筒底部的沉降值分別約為80 mm和120 mm。

圖1 20 a地表沉降

圖2 剖面豎向位移云圖

(2)應力值變化規律。井筒開挖，最大主應力出現在井筒的底部，為9.82 MPa，整個開采過程到20 a時，經過開挖和回填，模型達到新的力學平衡，內部應力重新進行分配，最大主應力出現在開采礦體的頂板和礦柱上，如圖3所示，最大主應力為15.32 MPa，但是井筒底部受到擾動影響，應力變化為8.77 MPa，出現了一定的應力釋放。

(3)塑性區變化規律分析。在整個開采過程中，礦體及圍巖的剪切破壞和拉伸破壞都主要集中在開采空區的頂板和礦柱區域，隨著開采范圍的不斷擴大，膠結充填體邊緣也出現了部分拉伸破壞，但是靠近1#回風井位置并未出現塑性區。

同理，對其它2個方案進行分析，得到以下數值模擬結論:

圖3 最大主應力云圖

(1)開采到20 a時，影響1#回風井穩定性的主要因素是上采區的開采擾動，由于上下采區之間存在厚大的未采動區域，極大地減弱了下部開采對回風井的影響;

(2)模擬結果對比見表2，同時根據水平位移，可計算得出3個位置在不同年份條件下的傾斜、曲率和水平變形(見表3)，均符合《有色金屬采礦設計規范(GB50771－2012)》要求，結合井筒與礦體的三維空間位置，位置3變形位移量較小，符合國家建構筑物規范，位置3為最佳選擇。

表2 3個擬選位置的數值模擬結果

表3 數值模擬計算移動變形值

4 井筒優化結果對比

按照有關礦山設計規定，通風井設計在地表移動線以外，風井距離礦體較遠，與風井優化位置3相比，每個中段增加巷道及通風路線1000 m，共有6各中段與通風井連接，增加了工程投資5280萬元;增加了通風機功率200 kW，年增加了電耗共1584000度。

5 結論

(1)根據模擬結果，推薦井筒位置3為最佳方案，減少了通風井巷工程量，節約了工程投資及生產成本。

(2)采用FLAC3D軟件對井筒位置進行模擬計算，得到了整個礦區的采動應力場及位移場的分布規律，說明了該方法的可行性。

(3)由模擬結果推薦井筒位置3的方案，但數值模擬的計算結果都是在現有的巖體力學參數、原巖應力場及簡化的開采工藝過程條件下進行的，由于巖石條件的多變性、巖石力學性質試驗取樣的局限性、巖體結構面的復雜性，計算所需的巖體力學參數難以準確確定。因此，建議在生產期，設置位移監測點，監測井筒位移的變化。

［1］王發芝.冬瓜山銅礦深部采場充填技術［J］.礦業研究與開發，2008，28(8):4-5.

［2］李克順，吳 姍，宋衛東.大冶鐵礦充填法采場結構參數校核與建議［J］.金屬礦山，2010(9):16-22.

［3］王金安，趙志宏，侯志鷹.淺埋堅硬覆巖下開采地表塌陷機理研究［J］.煤炭學報，2007，32(10):1051-1056.

［4］《采礦設計手冊》編輯部.采礦設計手冊(礦床開采卷)［M］.北京:中國建筑工業出版社，1987，895-899.

［5］GB50771－2012.有色金屬采礦設計規范［S］.

［6］吉學文，唐紹輝，李愛兵.某地下礦山采場與圍巖穩定性三維有限元模擬分析［J］.礦業研究與開發，2003，23(2):7-9.

［7］喬 蘭，王雙紅，蔡美峰.某地下礦山巖層及地表移動規律的有限元模擬研究［J］.金屬礦山，2000(4):23-25.