緩傾斜薄礦體開采空區群精確數字化及穩定性分析

代碧波 孫麗軍 趙興東

（1.東北大學采礦地壓與控制研究中心，遼寧沈陽，110004；2.金屬礦山安全與健康國家重點實驗室，安徽馬鞍山，243000；3.中鋼集團馬鞍山礦山研究總院股份有限公司，安徽馬鞍山，243000；4.武漢科技大學資源與環境工程學院，湖北武漢，430081）

緩傾斜薄礦體的非充填法回采會形成大量的采空區，這些空區群的穩定性直接關系到井下生產的安全［1-2］。故而對空區群實際三維形態的測定和穩定性評價至關重要［3］。

采空區的實際三維形態和準確的空間位置是影響空區穩定性數值模擬結果準確性的關鍵因素。過去，在對空區穩定性進行數值模擬時，往往將采空區的形態視為空區形成之前（即開采前）進行采場設計時確定的規則形態，部分則是依據傳統的測量方法測定的簡化空區形態進行模擬計算［3-5］。然而，采空區的實際形態往往在實際開采過程中，因測量、爆破等因素的影響而呈不規則狀，加之采空區形成后受爆破震動的影響，以及空區圍巖的冒落、片幫，致使空區的實際形態和開采設計中確定的相差甚遠、這就勢必嚴重影響空區穩定性數值模擬結果的精度和可靠性［6-8］。因此，如何精確地確定空區的實際形態和空間位置，并以此為基礎建立空區群穩定性數值模擬模型進行數值模擬計算分析，是空區群穩定性數值模擬研究中應加以解決的一個重要問題［4］。

針對上述問題，本項目開展以精密探測獲得的空區實際三維形態和空間位置為基礎的空區群穩定性數值模擬研究，形成一種新型的基于空區實際三維形態精密探測的空區群穩定性數值模擬技術，以提高數值模擬結果的精度和可靠性，并結合礦山井下開采現狀進行采空區的精密三維數字化及數值模擬技術研究，將所形成的技術應用于采空區群穩定性數值模擬分析中，并在此基礎上對其進行穩定性分析，為空區群的綜合治理提供技術支撐。

1 采空區探測技術

對于地下不明的采空區探測，我國目前以鉆探為主，物探為輔。常用的物探方法主要有以下幾類：①重力測量；②電法測量；③地震測量；④MT法、CSAMT法和EMAP法；⑤放射性測量。由于物探手段對采空區探測不能明確空區體積與形態，這些方法只能是大概探測地下有無空區以及空區的大致位置［9］。近年來發展起來的激光3D掃描法是利用激光測距原理對采空區進行精密探測的一種技術手段，只要能將激光探測設備伸入到空區內部，便可對空區進行360°掃描，經過后處理建立空區三維實體模型便可獲取空區的具體參數，包括空區邊界坐標，容積等［8，10］。

2 空區群三維數字化

2.1 工程概況

本研究依托的礦山為中鋼南非鉻業有限公司Dilokong鉻鐵礦，該礦山特點：①礦石品質優良、開發利用價值巨大；②礦體傾角緩、厚度薄，僅開采單一礦層，開采環境差、效率低；③礦山空區分布廣泛，僅采用木支護，井下開采安全條件差。Dilokong鉻鐵礦目前開采的資源為LG6鉻礦層。LG6礦層厚度約1.0～1.2 m，上覆有厚度約0.5 m厚的LG6a礦層，中間輝巖厚度大約0.5～0.6 m。區內礦體變化不大，LG6與LG6a礦層連續性較好，成礦后斷裂影響不大，似層狀產出。礦體整體走向南北，傾向西，傾角14°左右，淺部礦體較陡，約18°，深部傾角較緩，約7～11°。礦體遍布礦權范圍，長約7 800 m，寬約500～3 700 m，埋深約在300～400 m，礦體規模巨大。礦石Cr2O3品位為41.07%～49.3%，平均品位44.62%，品位變化均勻。礦體密度為4.25 t/m3。

2.2 礦區地表三維數字化

根據礦區地形地質圖，對礦區地形等高線進行三維坐標賦值，然后根據三維空間關系生成面域，并根據高程范圍進行著色渲染和水流方向的繪制。經三維數字化的地表模型如圖1所示。

2.3 采空區群三維數字化

根據礦山開采實際，將已有采空區群進行精密三維數字化，實現空區群的數值三維體態查看和分析。空區群位于名稱為Hanekom的老礦山范圍內。根據開采現狀，礦山共進行了4個中段的回采，形成了大量的采空區。根據采空區位置和開采現狀，利用激光3D掃描進行空區數字化，建立空區群的三維形態模型，見圖2。整個礦區空區群總體積為135.708萬m3。

3 空區群的數值模擬分析

應用三維拉格朗日差分分析程序模擬分析空區群的穩定性、頂板覆巖的變化規律以及留設礦柱參數對采空區穩定性所產生的影響［11，12］。

3.1 三維計算模型的建立

本次計算結合礦山工程地質條件、礦體賦存條件及所用采礦法的特點，以三維化空區群為依據建立模型，從采空區頂板的應力變化和位移，礦柱承載能力變化、采空區頂板應力及位移狀態和地表變形等方面進行空區群穩定性的數值模擬。最終確定的模型如下：模型長度3 100 m、寬度3 000 m、高度1 100 m，模型共有287 928個單元，434 768個節點。生成主體計算模型之后，其中空區群模型的構建采用實測三維數字化后的立體模型與本主體模型進行耦合，如圖3所示。

3.2 巖體力學參數

根據礦山的工程地質特征，礦體和圍巖特性，本次數值模擬過程中選取了對模型計算有影響的3種礦巖體，即輝巖、斜長巖和礦體。數值模擬過程中所用巖體力學參數見表1。原巖應力按常規進行假設計算，忽略構造應力，按照自重應力計算。

3.3 數值模擬結果

對建立的空區群模型進行計算，空區群內留設10 m×6 m間柱，按照下行式和前進式的回采順序逐步形成空區群，并在每個回采步驟記錄礦巖體應力和位移變化特征，經過模擬分析，得到了次生應力平衡之后的最大主應力云圖、最小主應力云圖、位移云圖。以采空區中央為中心，創建剖面，得到剖面位置的最大主應力云圖、最小主應力云圖、位移云圖。其數值分析模型計算第1步驟回采和第12步驟回采結果如圖4和圖5所示。

3.4 模擬結果分析

由以上數值模擬結果可知，原始地應力呈現層狀分布，在進行回采過程中，明顯擾動了原巖應力場，采空區部分的上覆巖層中最大主應力和最小主應力均存有共性，呈現似層狀分布，在臨近空區的礦巖交界處應力變化明顯，應力變化沿礦體傾向逐漸過渡。

進行第1步驟回采時，最大主應力為8.44 MPa，出現在上盤與底板的交界處，并在角隅處伴隨有應力集中現象，最大位移出現在下盤中央部位，為12.81 mm，表現為底鼓現象。進行第2步驟回采時，次生應力場的分布規律與第1步驟的相似，最大主應力為15.44 MPa，出現在下盤與頂柱的交界處，最大位移為16.99 mm，出現在本步驟內的上盤中央位置。

后續回采過程中，最大主應力值和最大位移值隨著回采步驟逐漸增大。在回采包括第10步驟在內的礦體時，最大主應力均出現在本步驟上盤和底板的交界處，表現為應力集中現象，但進行第11和第12步驟回采時，最大主應力始終在第10步驟中的上盤和底板交界處，至第12步驟回采結束并達到次生應力平衡時，最大主應力為66.14 MPa。最大主應力值隨回采步驟的變化如圖6所示。

在第5回采步驟之前，最大位移出現在上個步驟內的上盤中央位置，最大值為33.14 mm，第6至第8回采步驟中，最大位置始終在第4步驟的上盤中央位置，最大值為40.85 mm，第9回采步驟至第12回采步驟中，最大位移始終在第9步驟上盤中央位置，至回采結束時，最大位移為116.06 mm。最大位移隨回采步驟的變化如圖7所示。

從圖7可以看出，最大位移在第10回采步驟后趨于穩定。

另一方面，從塑性區計算可以看出，回采至第9步驟時，在上盤中央位置出現拉伸屈服區域，在計算結束后，大部分退出屈服狀態。至回采結束，未見有塑性區的貫通現象。表明現有空區群表現為穩定狀態。

由上述模擬結果可知，在形成空區群過程中明顯擾動原巖應力場，在相同條件下，采空區最大壓應力和最大拉應力隨礦柱寬度增大而減小，隨采空區頂板暴露面積的增大而增大。在數值模擬過程中，應力分布狀態以及塑性區發生了比較明顯的變化。當礦柱達不到與采空區頂板暴露面積最佳耦合的尺寸時，模型在計算過程中出現了屈服區域，計算結束后仍有較大范圍的屈服區域，此時采空區處于不穩定狀態，有發生頂板冒落的危險，當礦柱達到與采空區暴露面積最佳耦合狀態的尺寸時，有部分采空區圍巖在開采過程中曾經進入到拉伸和剪切屈服狀態，但在次生應力平衡之后又退出了屈服狀態。

另外，當增大采空區頂板暴露面積時，采空區上覆巖層的壓力傳遞到礦柱上，使得頂板中的壓力拱跨度增大，在采空區頂板產生拉應力區域，當此壓力拱跨度大于合理拱軸線對應拱形跨度，拉應力超過巖體抗拉強度時，便會在頂板位置產生拉伸屈服區域，造成采空區的破壞。在礦柱尺寸一定的情況下增大頂板暴露面積時，采空區頂板中央的下沉量也隨之增大，礦柱中的位移也明顯增大，與理論分析結果相符。

4 結論

（1）按照現有空區群內留設10 m×6 m間柱的情況下，經過數值模擬分析，空區群表現為穩定狀態。

（2）數值計算模型中，礦柱靠近頂板的位置為應力集中部位，并且由于礦體緩傾斜，使得應力集中部位偏向于礦柱靠近上部空區一側。

（3）給空區群模型留設合理尺寸的礦柱時，采空區頂板和礦柱在模型計算結束后均沒有拉伸屈服和剪切屈服區域的貫通，并且采空區頂板中央的拉應力均沒有超過礦體抗拉強度，頂板和礦柱內的位移趨于定值，并處于工程可控范圍，因此處于穩定狀態。

（4）采空區頂板中央的下沉量及礦柱內的最大主應力和位移對礦柱尺寸表現敏感，并且頂板暴露面積越大，此敏感度越高，表明在空場條件下，礦柱尺寸對采空區群的穩定性影響很大。