復雜采空區充填開采對地面建筑的影響研究

劉煥春 于建新 焦華喆 趙武勝

（1.河南理工大學安全科學與工程學院，河南焦作454000；2.河南省地下工程與災變防控重點實驗室，河南焦作454000；3.中國科學院武漢巖土力學研究所，湖北武漢 430071）

深部采動覆巖移動規律與控制是采礦學科的基本研究內容之一，眾多研究學者對巖層移動與地表沉陷規律進行了研究，取得了大量的研究成果[1-2]。金屬礦山由于礦體的賦存狀態及開采方法的特殊性，與煤礦地下開采存在較大不同[3-5]。任偉中等[6]通過地質力學模型試驗，分析了在厚覆蓋層條件下地表陷落特征隨不同開采水平的變化規律。Villegas T等[7]將Kiirunavaara礦山開采引起的地表變形隨時間的變化劃分為3個階段，并指出地表出現裂縫之前的應變擴展可能損害土木結構；并采用有限元和離散元軟件分析了開采對上盤破裂位置、破裂角以及上下盤的穩定性[8]。王運敏等[3]研究發現地表沉降與水平移動具有較強的一致性，整體變化趨勢具有階段性的快速增加與趨緩的特點。崩落法上覆巖層冒落對于緩解深部開采高應力、保證足夠覆蓋巖層厚度、控制巖層移動具有重要作用。Xia等[9]通過實際監測得出程潮鐵礦地下采礦引起的地表變形規律，與深部巖層移動密切相關。李庶林等[10]結合多通道微震監測、雙力偶點源理論、地表深鉆孔探測3種手段，得到了大紅山鐵礦采用無底柱分段崩落法開采不同時期高應力區、開裂帶和崩落帶的發展過程。宋許根等[4]依據實際監測結果，結合巖層移動機制，對塌陷發生以后不同時期、不同區域內的地表變形與擴展規律進行了分析，探討了降雨對地表變形與巖層移動的影響。陸玉根等[5]認為金屬礦山采用崩落法大規模深部開采后在地表形成了一個以裂隙為引導的塌陷區，明顯存在崩落帶和裂隙帶，彎曲變形帶不明顯，但其巖層變形過程仍具有“三帶”破壞特征。金屬礦山地下開采引起的地表塌陷主要采用現場實測、數值分析和相似模擬等方法進行研究，取得了一定的研究成果，保證了礦山開采的安全。但對于地表存在需要保護建筑物的情況，是不允許具有較大塌陷和變形的。而在開展現場監測時，首先采用數值方法對巖體變形及沉降規律進行把握是極有必要和指導意義的。

江西省內某金屬礦山采空區懸空面積大，地下空間復雜，礦柱分布不均，地表充填站位于采空區上方，一旦發生較大規模的礦柱崩塌、頂板大面積冒頂或上覆巖層的錯動，將直接威脅地面建筑的安全，甚至造成重大災害事故、產生嚴重的經濟損失和惡劣的社會后果。因此，為保證礦山安全生產，本項目利用大型有限元分析軟件ABAQUS對復雜采空區充填開采影響下的地表沉降規律以及充填站的穩定性進行研究。

1 工程概況

該金屬礦山位于江西省境內，礦體賦存于燕山晚期細粒黑云母花崗巖體與楊柳崗組含炭硅泥質灰巖的接觸帶上，呈似層狀、倒轉置鍋型產出，并隨接觸帶產狀的變化而變化，主要分布在0～24線間的1 200 m范圍內。礦體埋深40～300 m，礦體及頂底盤圍巖較穩定，水文地質條件簡單，礦體巖層分布如圖1所示。礦山采用地下方式開采，采用平硐開拓，鏟運機裝載，機車運輸，采礦方法為留不規則點柱、條柱的全面法，類似于房柱法，采場靜空跨度10～30 m不等，井下礦柱分布不均，部分頂板暴露面積過大，礦體回采多沿底板掘進，許多采空區頂部礦體殘留，開采難度極大，為避免資源的浪費，合理安全有效地開采殘余礦體，采用充填法對頂部殘礦進行開采。地表充填站標高682.0 m，在10～12線之間，位于采空區上方，占地面積3 726 m2，充填站位置、采空區范圍及礦柱分布如圖2所示。

2 數值模型

2.1 修正Mohr-Coulomb理論模型

在描述巖土材料的強度特性和破壞行為方面，Mohr-Coulomb強度理論得到了廣泛應用[11]。本項目采用文獻[12]的修正Mohr-Coulomb模型，可以避免數值計算的發散和收斂的緩慢。

修正的Mohr-Coulomb屈服準則的表達式為

其中：

式中，σ1＞σ2＞σ3為主應力；c、φ分別表示黏聚力和內摩擦角；σm為平均應力為等效應力；J2、J3分別為應力偏量第二不變量和第三不變量；sx=σm-σy，sy=σy-σm，sz=σz-σm；θ為Lode角。

采用分段函數來描述K(θ)，具體表達式為

式中，

取θT=25°，|θ|≤θT時，在ρ平面屈服函數跡線不做處理，和經典的Mohr-Coulomb準則一致，而當|θ|＞θT時，對屈服函數的跡線進行光滑處理。

經過四十余年的發展，以氮肥為起點的云天化集團取得了令人矚目的成就，同時也面臨著國家宏觀政策、安全環保形勢、市場競爭環境、農業需求等諸多方面的變化與挑戰。為了適應這些新變化、新挑戰，云天化主動進行戰略布局的調整，圍繞“綠色科技服務現代農業”的理念，積極探索轉型升級的新思路、新道路。

取塑性勢函數與屈服函數的表達式一致，即：

式中，φ為膨脹角，其中K(θ)也與φ有關，表達式與屈服函數中的K(θ)類似。若φ=φ，則為關聯流動；若φ≠φ，稱為非關聯流動，而當φ=0時塑性變形時材料的體積不發生變化。

2.2 三維模型建立

根據山體等高線和采空區內礦柱分布，建立10～12線間灰巖—礦體、礦體—花崗巖、采空區、礦柱的三維山體模型，如圖3所示。計算模型長約400 m，充填站距左側邊界150 m，距右側邊界250 m；模型寬度范圍為10～12線約100 m，模型底部高程為500 m，山頂最高處高程約790 m。采用四面體單元對計算模型進行離散，生成計算網格，模型共生成單元約90 000個。

模型采用THM耦合4節點軸對稱減縮積分單元，由于礦體開挖引起的位移主要是豎直方向位移，因此，在計算過程中，除山體上表面外，其余5個面施加法向約束。計算區域內礦柱的具體參數如表1所示，礦山巖體的物理力學參數如表2所示。

?

?

礦區地應力測試得到的最大水平主應力、最小水平主應力和垂直主應力隨深度H（單位為m）的線性回歸方程為

由于山體為不規則地表，無法通過常規方法施加地應力。為此，采用ABAQUS中子程序SIGINI通過編程進行施加，具體過程[13]如下。

根據山體表面的形狀，擬合山頂標高方程，得到不同位置處山頂標高；根據山頂標高，得到不同位置處的埋深；根據地應力測試結果，計算得到該位置處的3個主應力；將主應力轉化到局部坐標系中，在局部坐標系下施加初始應力。最終，山體內部形成沿著山體坡面方向分布的初始應力場，其中豎直方向上的正應力分布如圖4所示。

根據礦體開挖及充填的順序，數值計算共分5步：①地應力平衡，施加原巖應力，在巖體中形成原巖應力場；②底部開挖，根據礦體開挖順序首先模擬底部礦體的開挖，預留礦柱；③底部充填，底部開挖完成后，對底部進行分層充填；④上部開挖，充填完成后，對上部礦體進行開挖；⑤采空區充填，上部礦體開挖完成后，對采空區進行充填。

3 模擬結果分析

計算結果中的符號說明：+表示應力為拉應力，位移與坐標軸同向，如在采場底板表示底鼓，頂板表示下沉，采場側幫表示位移向采空區移動。-表示應力為壓應力，位移與坐標軸反向。x方向代表平行于勘探線方向，y方向代表垂直勘探線方向，z向代表豎直方向。

3.1 巖體應力與位移變化

礦體開挖后，巖體主應力分布和位移變化如圖5與圖6所示。從圖5中可以看出，開挖引起了山體應力調整，開挖充填完成后巖體幾乎全部受壓，在采空區頂板與礦柱出現了明顯的應力集中，礦柱及頂板部位x方向最大拉應力均小于5.0 MPa。巖體中y方向應力較其他2個方向調整最小，但在礦柱部位出現了少量拉應力區。從豎直方向應力的分布可以看出，礦柱在豎直方向主要承受壓應力，礦體開挖導致上覆巖層的應力向礦柱轉移，礦柱承受了很大的豎向壓應力，最大壓應力達到了63 MPa。

從圖6可以看出，礦體開挖后，在x方向（平行于勘探線方向）上部巖體向采空區移動，地表的變形表現為采空區左側向右側移動，x方向最大位移出現在采空區的左上方，最大位移量約6.2 mm。礦巖位移狀態中，垂直位移比水平位移重要的多，因為過量的垂直位移將直接導致頂板的冒落，因此，該部分主要對頂板的垂直位移進行分析。在z方向（豎直方向），礦體開挖導致頂板下沉與底鼓，礦體頂板最大下降量約18.2 mm，充填在一定程度上減小了底鼓量，最大底鼓量約3.9 mm。

3.2 礦柱受力及變形分析

礦體開挖后礦柱受力如圖7所示，礦柱受壓縮變形如圖8所示。在2個水平方向上，礦柱出現了拉應力，最大拉應力約1.2 MPa；在豎直方向上，礦柱主要承受壓應力，且表現為采空區中部的礦柱受力大，采空區邊界處的礦柱受力小，豎直方向最大壓應力約63 MPa。由于山體與開挖體的不對稱性，在礦體開挖過程中引起了礦柱的水平位移，但礦柱在2個水平方向上的位移量較小，最大水平位移量不足1.5 mm。由于上覆巖體向下移動、底板向上移動，礦柱豎向位移遠大于水平方向位移，越靠近采空區中部、礦柱尺寸越小，礦柱的壓縮量越大，最大豎向壓縮量約15.58 mm。

3.3 地表變形規律

為了分析礦體開挖及采空區充填作用下地表的變形情況，以充填站為中心，布置如圖9所示的A、B2類監測線，其中沿山坡方向布置 5 條（A-50、A-25、A0、A+25、A+50），相鄰監測線間距25 m；在垂直山坡方向布置 7 條（B-60、B-40、B-20、B0、B+20、B+40、B+60），相鄰監測線間距20.0 m。A類監測線上的位移變化如圖10所示，B類監測線上的位移變化如圖11所示。

從圖10可以看出，沿山坡方向，5條監測曲線上x方向的位移變化趨勢相同。在充填站左側（山腳側）120 m的位置處，地表x方向位移接近于零；隨著監測點靠近充填站，x方向的位移不斷增加，監測站位置處x方向位移約為4.0 mm；地表x方向位移最大值約6.0 mm，出現在充填站右側（山頂側）40 m的位置。由于計算模型在y方向尺寸的限制，地表y方向的位移較小，最大約0.6 mm。沿山坡方向，在充填站左側（山腳側）120 m的位置處，地表z方向位移接近于零；隨著監測點沿山坡向上，z方向的位移不斷增加，監測站位置處地表z方向位移約為2.6 mm；在充填站右則120 m以后，地表z方向的位移趨于穩定，最大值約15.0 mm。

從圖11可以看出，礦體開挖及采空區充填后，越靠近山腳，地表x方向的位移越小，越靠近山頂，地表x方向的位移越大。不同監測線之間，地表x方向的位移變化曲線基本一致。充填站處地表x方向最大位移約4.0 mm。由于計算模型在y方向尺寸的限制，充填站附近地表y方向的位移較小，最大約0.15 mm。礦體開挖及采空區充填后，地表z方向位移與x方向位移變化趨勢相同。即：越靠近山腳，地表z方向的位移越小，越靠近山頂，地表z方向的位移越大。充填站處地表x方向最大位移約2.8 mm。

4 結論

（1）礦體開采及采空區充填造成巖體應力的重分布，礦體頂板下沉，底板鼓起，頂板最大下沉量約18.2 mm，充填在一定程度上減小了底鼓量，最大底鼓量約3.9 mm。

（2）礦體開采造成礦柱承受較大的豎向壓力，在巖體壓力作用下，礦柱豎直方向最大壓應力約63 MPa；同時，礦柱局部出現拉應力，最大拉應力約1.2 MPa，可能出現少量拉剪破壞。在巖體壓力作用下，礦柱出現明顯的壓縮，且越靠近采空區中部、礦柱尺寸越小，礦柱的壓縮量越大，最大壓縮量約15.58 mm。

（3）充填站附近，最大水平位移約4.0 mm、最大沉降約2.6 mm。越靠近山腳（遠離采空區），地表的位移越小，越靠近山頂，地表的位移越大。地表x方向位移最大值約6.0 mm，出現在充填站右側（山頂側）40m的位置。地表最大沉降量約15.0 mm，出現在充填站右則（山頂側）120 m附近。

（4）研究成果為后期現場地表及建筑物沉降觀測方案的制定提供了參考，為礦山進行殘礦安全回采提供了依據。