基于FLAC3D 的臨渙礦馬頭門優化方案數值模擬分析

李承濤 吳夢宇

（安徽理工大學土木建筑學院，安徽 淮南 232001）

0 引言

馬頭門是煤礦井筒與煤礦運輸巷道的連接部分，是礦井的重要組成部分，包括了馬頭口上井筒段、下井筒段及馬頭門兩側的連接大銅室，因為其位于深部地層，故其受力情況多變且復雜[1]。馬頭門具有斷面大、服務年限長、特別容易受到應力集中等特點[2]，因此在馬頭門與井壁交界部位就容易出現較嚴重的應力集中現象，造成馬門頭發生破壞[3]。馬頭門發生破壞的原因有很多，除了馬頭門位置不佳、斷面大、設計不合理、支護方案效果不好、施工過程質量問題等原因外，還包括馬頭門所處的地層特征，由于馬門頭多位于軟弱巖層中，加上馬頭門支護施工難度大，在地應力特別是在深部地層高地應力的作用下非常容易發生破壞［4-5］。所以，研究馬頭門支護結構的穩定性對礦井建設的安全非常有必要，能夠保障礦工的生命安全和減少煤炭資源的損失，這不僅會使我國在煤炭工程領域安全生產方面擁有更加先進的技術，也可以促使馬頭門在其他領域得到更廣闊的應用前景。

目前，對于馬頭門的研究成果交多，尤其是采用數值模擬的研究方法具有成本小、速度快的優勢，逐漸得到應用，其中應用FLAC3D軟件研究馬頭門的方法效果較好[6]。姚直書[7]等利用FLAC3D技術對開挖和支護全過程5個階段進行了模擬，探究了馬頭門在開挖和掘進時圍巖受所變化的特性，展現了開挖和掘進過程對馬頭門的影響；田柯[8]通過FLAC3D軟件對馬頭門容易產生應力集中和圍巖容易破壞的區域進行了定性分析，并且對應力集中系數和深度做出了定量的確定。但利用FLAC3D軟件進行馬門頭優化設計的研究仍然不足。

基于此，以臨渙礦中央風井中馬頭門施工支護過程為例，利用FLAC3D軟件開展馬門頭支護方案優化設計的數值分析，以供相關技術人員參考。

1 工程模型建立

1.1 FLAC3D軟件介紹

本研究采用FLAC3D軟件進行數值模擬分析。FLAC3D軟件是一款用于仿真計算的軟件，能夠對土質、巖石等多種不同性質材料的受力特性以及對具有塑性、流動特性的材料進行三維的數值分析。該軟件在三維網格中通過調整多面體單元從而實現擬合實際的結構。該軟件可采用線性或非線性的單元材料，當外力作用在材料上時，發生屈服流動的材料，其三維網格能夠發生相應的變形和移動，這就是大變形模式。該軟件的技術特點為：顯式拉格朗日算法以及混合-離散分區技術，正是采用了這些技術，FLAC3D才能非常精準地對材料的塑性破壞和流動進行模擬。正因為它不需要剛度矩陣，所以它在極少的內存下就可以求解極多的單元結構。

FLAC3D軟件擁有以下3種計算優勢：

（1）混合離散法：當模擬的物體是塑性破壞或者塑性流動時，該方法相比采用“離散集成法”的有限元法更佳；

（2）動態運動方程：即使某個模擬的系統狀態為靜態，數值模擬過程中仍然使用動態的運動方程，故在模擬物體具有不穩定的狀態和過程時，FLAC3D軟件不存在數值上的障礙；

（3）“顯式解”方案：“顯式解”方案在線性本構關系與非線性本構關系的應力-應變求解上所花費的時間是相同的，與之相對的是隱式求解方案，該方案對于非線性問題處理的工作量比“顯式解”方案多得多。因為它不需要存儲剛度矩陣，這就意味著，在極少的內存下就可以求解極多的單元結構，并且在處理大變形和處理小變形的模擬上所消耗的時間相差無幾。

1.2 工程概況

臨渙礦跨度為13km，寬度約為4～5km，面積49.6624km2。該礦開采的模式為：使用立井、采用分區通風以及集中出煤的方式，礦建工程包括中央風井、-650m回風車場及回風石門、-450m回風聯巷。臨渙礦中央風井中馬頭門位于約-450m的位置，為單側馬頭門，井筒凈直徑為6.5m，馬頭門及巷道凈半徑為2.6m。

地面與馬頭門之間包含3層巖層，自上至下依次為泥巖、粉砂巖及泥巖，厚度分別為20.04m、3.27m以及16.69m，巖層物理與力學參數如表1所示，馬頭門及井筒襯砌物理與力學參數如表2所示。

表1 巖層物理與力學參數

表2 襯砌物理與力學參數

1.3 模型建立

根據實際工程深立井、馬頭門及其連接巷道的平面布置、斷面尺寸和水平位置，考慮邊界效應的影響范圍，使用Rhino+Griddle軟件建立模型網格，通過將網格導入FLAC3D建立三維有限元計算模型。模型尺寸為長50m、寬40m、高50m，單元總量為1509175。在模型側面及底面都施加位移約束，在模型頂面施加11.2MPa的豎向應力，如圖1所示。

圖1 水平馬頭門數值模型

對于所研究區域，中央風井井筒斷面為圓形，巷道及馬頭門硐室斷面為直墻半圓拱形。初期支護時，馬頭門及井筒均采用錨、噴、網的支護形式，馬頭門及井筒的錨桿長度為2.5m，間距、排距為0.8m×0.8m，預應力為60kN；馬頭門處錨索長度為6.2m，間距、排距也為0.8m×0.8m，預應力為120kN，井筒及馬頭門襯砌澆筑采用C50混凝土。在數值計算模型中，巷道、回風聯巷及其連接處共分4個步驟進行模擬，具體描述如下：

（1）馬頭門上部井筒掘進開挖，設置錨桿等參數，模擬讓壓過程；

（2）馬頭門上部井筒二次襯砌，計算平衡；

（3）馬頭門掘進開挖，設置錨桿等參數，模擬讓壓過程；

（4）馬頭門二次襯砌，計算平衡。

2 馬頭門初期支護方案數值分析

模型建立完成后，對井筒進行開挖。首先對上井筒段進行開挖，開挖完畢后對其進行支護；然后對馬頭門進行單側開挖，開挖完畢后進行二次襯砌。利用FLAC3D軟件得到其豎向位移、豎向應力以及最大主應力，將云圖導入Tecplot軟件進行分析，得出位移及應力云圖如圖2~圖4所示。由圖2可知，馬頭門開挖后進行了錨網噴支護，距離巷道、馬頭門空間越近，對應的巖體位移越大,此時巷道頂板最大下沉降位移超2.5cm，底板底鼓位移已超4.5cm。根據圖3和圖4可知，豎向應力最大達到了38MPa，最大主應力達到了110MPa，支護失效風險極大。因此，需要對支護方案進行優化。

3 馬頭門優化支護方案數值分析

3.1 優化方案一數值分析

考慮到有限差分軟件模擬單元需連續，導致圍巖應力釋放產生的變形將作用到襯砌結構上，這與實際情況存在差異。方案一為在井筒及馬頭門襯砌施加前，將圍巖應力釋放產生的位移清零，其他支護參數都保持不變，導出位移及應力云圖見圖5～圖7。

對圍巖應力釋放產生的位移清零進行計算后，對比圖5和圖2不難發現：優化方案一處馬頭門處襯砌整體下沉，其中地板下沉約0.5cm，頂板下沉位移超過0.7cm，最大位移已經降低了許多，說明方案一的優化設計對馬頭門的穩定性有很大的提高；而通過對圖3與圖6，圖4與圖7的對比發現，豎向應力和最大主應力在優化前后應力值大小未發生較大改變，說明優化方案一并不能改變豎向應力的大小。

圖2 z向位移（單位：cm）

圖3 z向應力（單位：MPa）

圖4 最大主應力（單位：MPa）

3.2 優化方案二數值分析

優化方案一的模擬結果證實，現有的措施基本能滿足支護要求。為了安全考慮，優化方案而在優化方案一的基礎上對馬頭門15m范圍內的巖體進行注漿加固，用以改善巖性，依舊采用Tecplot軟件得到的位移云圖和應力云圖如圖8～圖10所示。根據圖8可以看出，相較于圖5中馬頭門整體位移變化不大，說明注漿加固對馬頭門附近圍巖的位移影響不大；對比圖6和圖9，可以清楚地看出豎向應力也未發生明顯變化，只是從原來的38MPa降低為34MPa，對馬頭門整體影響不大；對比圖7和圖10能明顯看到最大主應力的巨大變化，圖7中最大主應力由110MPa變為圖10中約40MPa，說明馬頭門15m范圍內的巖體進行注漿加固對降低圍巖最大主應力有很好的效果，綜合未優化、優化方案一、優化方案二這三個支護方案來看，雖然優化方案一已經基本滿足馬頭門支護的效果，但優化方案二的支護效果更佳、結構穩定性更好。

圖5 優化方案一z向位移（單位：cm）

圖6 優化方案一z向應力（單位：MPa）

圖7 優化方案一最大主應力（單位：MPa）

圖8 優化方案二z向位移（單位：cm）

圖9 優化方案二z向應力（單位：MPa）

圖10 優化方案二最大主應力（單位：MPa）

4 結束語

通過FLAC3D軟件建立了臨渙礦馬門頭數值計算模型，分析了單側馬頭門采用錨網噴支護以及二次襯砌下的位移場和應力場，由于產生的位移較大，馬頭門支護失效風險較大，故需要進行優化方案設計。針對初期設計，設計了兩種優化方案，一種為在井筒及馬頭門襯砌施加前，將圍巖應力釋放產生的位移清零，其他支護參數都保持不變；另一種為在優化方案一的基礎上對馬頭門15m范圍內的巖體進行注漿加固，并對比優化前后的位移場和應力場。通過未優化與優化方案一對比可知，馬頭門處位移大幅下降，豎向應力和最大主應力未發生變化；通過優化方案一與優化方案二的對比可知，馬頭門處位移、豎向應力未發生明顯變化，但最大主應力發生明顯變化，故本模擬采用優化方案二的支護效果更佳、結構穩定性更好。