李金洋 臧明東 徐能雄 秦 嚴 陳 宇

(中國地質大學(北京)工程技術學院， 北京 100083， 中國)

0 引 言

煤炭資源是經濟發(fā)展與社會進步的重要物質基礎。煤炭開采使得地下礦體形成空洞，圍巖的原始應力平衡被打破，導致圍巖應力重分布(馬鳳山等， 2018)，造成上覆巖層及地表劇烈破壞，巖層發(fā)生變形、開裂、破碎與冒落，形成大量采空區(qū)。采空區(qū)不但破壞范圍廣、變形量大，而且變形持續(xù)時間長。采空區(qū)巖層的殘余變形持續(xù)時間可達數年甚至數十年(郭廣禮， 2001)，對上方構筑物危害巨大。以往遇到采空區(qū)通常采取避讓的方式，盡可能不去觸碰涉及多層采空區(qū)的建設場地。然而，隨著我國公路鐵路網的不斷加密，尤其在山西等煤炭生產大省，越來越多的交通隧道不得不穿越采空區(qū)(張永波等， 2007； 段麗軍等， 2017)，采空區(qū)造成隧道襯砌結構變形破壞以及路面和路基下沉成為公路建設中的常見問題(葉懿尉等， 2018)。如何有效評價采空區(qū)殘余變形對隧道的影響成為穿越采空區(qū)隧道工程建設的關鍵。

1956年，波蘭學者Litwiniszyn(1956)提出隨機介質理論來研究巖層移動規(guī)律。隨機介質理論將礦山巖體視為不連續(xù)介質，由剛性顆粒組成，則巖層的移動可視為剛性顆粒在重力作用下向下方開采空穴的移動過程。20世紀60年代初期，劉寶琛等(1965)在隨機介質理論(Litwiniszyn， 1956)的基礎上解算出其簡化解，提出了概率積分法，實現了從理論方法到應用的發(fā)展，成為目前我國預計地表沉陷的主要方法(姜岳等， 2019)，后來也被廣泛用于采空區(qū)地表殘余沉降的預計。劉天泉等(1995)通過研究給出了基于概率積分法的裂隙帶高度計算方法。余學義等(1997)編譯了評價采空區(qū)對公路影響程度的概率積分法程序。李仁民等(2002)對概率積分法進行改進，提出了適用于多層采空區(qū)的地表沉陷預測方法。陳曉斌等(2007)建立了基于隨機介質理論的采空區(qū)地面變形計算模型，并編制了相關計算程序。王磊等(2011)和朱廣軼等(2014)提出了概率積分法等價變采厚計算模型。蔡音飛等(2016)提出了概率積分法的優(yōu)化方法，將概率積分法應用于起伏地表的采空區(qū)殘余沉降計算中。李亞蘭等(2018)基于隨機介質理論，推導出了半圓拱形斷面巷道開挖引起地表下沉和水平位移的計算公式，為深部工程的規(guī)劃設計提供了參考依據。

概率積分法在采空區(qū)地表及巖層殘余變形預測方面已經得到廣泛應用，但在預測穿越采空區(qū)隧道的殘余變形方面鮮見報道。本文以穿越山西省陽泉二礦采空區(qū)的新建桑掌隧道為案例，引入玻茲曼函數對等價變采厚概率積分法中的邊界空洞殘余可活化厚度進行修正，并引入時間函數對該地區(qū)采空區(qū)殘余變形進行了預測，分析了殘余變形對桑掌隧道建設的影響。研究成果可為隧道工程建設提供方法支撐。

表1 項目區(qū)地層巖性描述表

1 研究區(qū)地質背景及工程概況

1.1 區(qū)域自然地理及工程地質條件

陽泉二礦礦區(qū)位于山西省東部的陽泉市西南部郊區(qū)，所選研究區(qū)位于陽泉二礦西北部(圖1)，該區(qū)受內陸高原、暖溫帶季風氣候帶的大陸性氣候影響，四季分明，干燥多風，十年九旱，年平均氣溫為10.8℃，多年平均降雨量為572mm，雨量集中在7、8兩個月，降水量達340mm，占年平均降水量的66.4%。區(qū)內植被具有一定的垂直分布特征，沖洪積平原和山間河谷區(qū)以草甸植物群落為主，低中山區(qū)以木本、草灌植物群落為主。

圖1 研究區(qū)位置(據GoogleEarth衛(wèi)星影像修改)

礦區(qū)地處太行山中部，整體地勢西高東低，平均海拔1000m，相對高差約540m，呈典型的構造剝蝕地貌。在構造和長期剝蝕切割作用下，形成了山脈縱橫、峰起巒連、溝谷深切的地貌形態(tài)(圖1)。如圖2所示，區(qū)內出露地層主要為三疊系下統(tǒng)和尚溝組(T1h)的砂質泥巖夾長石砂巖和劉家溝組(T1l)的長石砂巖夾砂質泥巖以及二疊系上統(tǒng)上石盒子組(P2s)和二疊系下統(tǒng)下石盒子組(P1x)的砂巖、頁巖、泥質砂巖、砂質泥巖及泥巖(表1)。

圖2 研究區(qū)區(qū)域地質圖(據馬麗芳(2002)修改)

1.2 工程概況

陽泉二礦建礦于1951年5月，井田位于沁水煤田東北部，面積約60.06km2。陽泉二礦煤層穩(wěn)定，賦存豐富，剩余服務年限30年以上(山西省交通規(guī)劃勘察設計院， 2016)。煤礦開采深度463.3～713.5m，主要可采煤有3#、8#和15#煤層，平均厚度分別為1.96m， 2.26m以及6.42m，采煤方法為綜合機械化放頂煤開采，屬走向長壁后退式采煤法，頂板管理為全部垮落法。其中： 3#煤層工作面回采率98.6%， 15#煤層工作面回采率88.5%。新建桑掌隧道位于研究區(qū)中部(圖3)，全長925m，平均埋深80m。隧道從上石盒子組下段地層穿過，距離3#、8#和15#煤的平均垂直高差分別為275m、320m和400m(圖4)。由于3#和8#煤層的工作面停采時間較長，對隧道影響較小，因此本文主要研究15#煤層的工作面對隧道的影響。隧道主要受15#煤層的4個采煤工作面影響(圖3)。其中： ①號工作面位于隧道起點東南側，工作面長軸方向與隧道近似平行，該工作面于2019年9月停采； ②號和③號工作面沿長軸方向與隧道小角度相交，分別于2017年9月和2015年4月停采； ④號工作面位于隧道東側，沿長軸方向平行于②號和③號工作面，于2012年3月停采。

圖3 研究區(qū)工程地質平面圖(地層符號含義見表1)

圖4 I—I′剖面圖(剖面位置見圖3)

2 計算方法

2.1 計算模型

本文采用等價變采厚的概率積分法模型(王磊等， 2011； 朱廣軼等， 2014； 圖6)對隧道軸線所在平面的殘余變形進行計算，利用玻茲曼函數對邊界區(qū)殘余可活化厚度進行修正。同時，引入Knothe時間函數(Knothe， 1952)建立殘余下沉量與時間的定量關系，用于綜合考慮停采歷史不同的工作面在同一時刻的殘余變形量。

2.1.1 改進的等價變采厚概率積分法模型

如圖5所示，煤層開采后上覆巖層依次形成垮落帶、斷裂帶和彎曲帶(郭廣禮等， 2002)。在垮落帶內，巖體破碎，由于邊界處的巖梁作用，自老采空區(qū)邊界至采空中央具有明顯的分區(qū)性：空洞區(qū)、未充分壓實區(qū)和較充分壓實區(qū)(王磊等， 2011； 朱廣軼等， 2014)。針對長壁后退式開采形成的老采空區(qū)活化引起的地表殘余沉降，從老采空區(qū)的覆巖形態(tài)結構出發(fā)，采用等價變采厚的概率積分法模型(王磊等， 2011； 朱廣軼等， 2014)進行計算：

圖5 采空區(qū)覆巖結構示意圖(郭廣禮等， 2002)

(1)

m1=m(1-q)

(2)

(3)

式中：W0(x)為殘余沉降最大值；x為地表任意點的橫坐標；s為開采單元的橫坐標； AC、CD、DB為積分區(qū)間(圖6)；m′與m1分別表示邊界區(qū)與中部區(qū)殘余可活化厚度；We(x)為單元開采形成的單元下沉；m為采厚；q為殘余變形剛要出現時刻的地表下沉系數；r為影響半徑。

如圖6所示，由于采空區(qū)邊界空洞區(qū)殘余可活化厚度是線性變化的，等價變采厚的概率積分法在此處無法積分，前人在計算時常根據采厚m進行簡單估算(王磊等， 2011； 朱廣軼等， 2014)，但此假定與邊界空洞區(qū)實際情況不符，會造成邊界區(qū)殘余變形量比實際允許的殘余變形量大，導致地表殘余變形量的預測偏大。鑒于邊界空洞區(qū)實際上為非充分沉陷，本文采用玻茲曼函數對邊界區(qū)殘余可活化厚度進行修正，計算方法為：

圖6 等價變采厚的概率積分法模型示意圖(據朱廣軼等(2014)，郭慶彪(2017)修改)

m′=m[1-qρW(kL)]

(4)

其中

(5)

2.1.2 基于時序的殘余變形預測模型

采用上一小節(jié)所述方法，雖可獲得開采后地表殘余變形的最大值，但是殘余變形隨時間變化的中間過程并不清楚。Knothe時間函數(Knothe， 1952)作為開采沉陷預測中常用的時間函數被廣泛用于描述地表動態(tài)下沉過程(王文學等， 2017)。因此，本文引入Knothe時間函數(Knothe， 1952)來刻畫地表下沉量與時間的關系：

W(x，t)=W0(x)(1-e-λt)

(6)

式中：W0(x)為殘余沉降最大值，可由式(1)求得；λ為與上覆巖層力學性質相關的時間因素影響系數。

將式(1)～式(5)帶入式(6)中，可以得到基于時序的殘余變形預測模型：

(7)

2.2 參數選取

采用概率積分法對采空區(qū)上覆巖層與地表的殘余變形進行預測時，計算參數的選取是關鍵。根據《建筑物、水體、鐵路及主要井巷煤柱留設與壓煤開采規(guī)范》，計算參數可以根據采空區(qū)覆巖的巖性進行選取，如表2所示。

表2 覆巖參數表(國家安全監(jiān)管總局等， 2017)

根據現場勘察，新建桑掌隧道下伏采空區(qū)為3層采空區(qū)，采空區(qū)埋深463.3～713.5m。其中隧道下伏3#煤層和8#煤層與隧道的平均垂直高差分別為275m與320m，這兩層煤的停采時間均超過10年，因此，在計算時可將下伏多層采空區(qū)看作重復采動情況，進行相應的參數調整(國家安全監(jiān)管總局等， 2017)。下伏15#煤層與隧道的平均垂直高差為400m，煤層傾角較小，可近似看作0°。15#煤層覆巖為砂巖、泥巖和粉砂巖，根據表2，覆巖類型為中硬，經計算，非充分沉陷率ρW(kL)=0.81。據艾亞民等(1990)對陽泉地區(qū)煤礦的長期觀測與工程實踐，陽泉地區(qū)充分采動的下沉系數為0.85，重復采動情況下，下沉系數為q=0.85+0.1±0.05，本次計算下沉系數取q=0.92，其余計算參數如表3所示。

通過對時間函數的分析可知，參數λ的值越大，時間函數會在越短的時間內趨于一個定值。由于壓實區(qū)會先于空洞區(qū)完成殘余變形，所以壓實區(qū)的λ取值一般大于空洞區(qū)的取值。結合陽泉礦務局對二礦周邊礦區(qū)參數λ選取的經驗，本區(qū)內參數λ的取值如表4所示。

表3 計算參數表

表4 時間影響參數λ取值范圍

3 采空區(qū)殘余變形預測

首先，采用上節(jié)所述方法對陽泉二礦一工作面停采后293天的累積沉降量進行試算，預測的最大沉降量為140mm，與實測數據145mm相近(孫富強， 2018)。因此，本文提出的預測模型具有一定的可靠性，可用來對本研究區(qū)內采空區(qū)的殘余變形量進行預測。

將表3所示的計算參數代入式(1)，通過分區(qū)及積分疊加的方法計算各工作面停采后的最大殘余變形量。其中：工作面沿開采方向的最大殘余沉降量和最大殘余水平移動如圖7所示，垂直開采方向的最大殘余沉降量和最大殘余水平移動如圖8所示。工作面停采后地表最大殘余沉降量為647mm，最大殘余水平移動為255mm。

圖7 沿開采方向工作面殘余沉降與殘余水平移動

圖8 ①～④號工作面垂直開采方向殘余沉降與殘余水平移動

如圖3所示，桑掌隧道受4個工作面影響。隧道起點距離②、③、④號工作面較遠，在一定范圍內只受到①號工作面的影響。隧道中段部分穿越了②號和③號工作面的上覆巖層，同時也受到①號和④號工作面的影響。4個工作面的停采時間各不相同，因此，采空區(qū)殘余變形對隧道的影響需要同時考慮工作面停采時間和工作面影響范圍兩個因素。

首先，考慮工作面停采時間這個因素。如前所述，本文已經通過式(1)計算出了4個工作面停采后潛在的最大殘余變形量。在此基礎上引入時間函數，并認為采空區(qū)變形總量為煤層厚度，即最終下沉系數為1。隧道開挖時， ①、②、③、④號工作面已分別停采1年、3年、5年和8年，通過式(6)計算得到各工作面在2020年之后在隧道軸線所在平面的殘余變形量。

然后，考慮工作面影響范圍這個因素。由于4個工作面距離較近，相互之間會有影響。針對區(qū)域內一點，通過對單個工作面影響下的殘余變形進行疊加獲取該點的累計殘余變形量，并據此做出隧道軸線所在平面的殘余沉降等值線(圖9)。

圖9 隧道軸線所在平面殘余沉降等值線圖

沿隧道軸線做剖面，可得隧道沿軸線方向上的殘余沉降，如圖10所示，以及隧道沿軸線方向上的殘余水平位移，如圖11所示。其中：殘余水平位移可分為垂直于隧道軸線方向和平行于隧道軸線方向兩個方向的水平位移(圖11)。

圖10 隧道殘余沉降(計算起點見圖9)

圖11 隧道殘余水平位移(計算起點見圖9)

4 采空區(qū)殘余變形對隧道的影響分析

通過對比圖7和圖8可以發(fā)現，沿開采方向，工作面殘余沉降呈“W”型，工作面殘余水平移動呈起伏波動的形態(tài)，兩者最大值均出現在開切眼與停采線附近的空洞區(qū)，而工作面中部的較充分壓實區(qū)殘余變形量相對較小。這主要是由于采用長壁垮落法開采時，開切眼與停采線附近出現了空洞區(qū)，致使其潛在可沉降空間較大。垂直開采方向，工作面殘余沉降呈“V”型，工作面殘余水平移動亦呈起伏波動形態(tài)，其中前者最大值出現在斷面中心位置，而后者最大值出現在靠近邊界的位置。這是由于垂直開采方向的懸空跨度較小，豎向變形集中發(fā)生在中部，而側向位移集中出現在邊界位置。

由圖9可知， ②、③和④號工作面沿開采方向平行展布，工作面殘余沉降較為集中地出現在②號工作面上方，在隧道軸線所在平面內潛在的最大沉降量為325mm。②、③和④號工作面停采較早，采空區(qū)沉降歷史相對較長，各工作面殘余沉降量較之①號工作面偏小。由于①號工作面停采時間只有一年，沉降不充分，殘余沉降量較大，在隧道軸線所在平面內潛在的最大沉降量為617mm。由于②、③和④號工作面距離較近，相互之間有影響，因此該區(qū)域內殘余沉降等值線圖的形態(tài)與①號工作面存在較大差異。而①號工作面距離其他3個工作面較遠，采空區(qū)變形相對獨立，同時，可以看出計算得出的①號工作面殘余沉降呈現分區(qū)特征(圖9)，邊界區(qū)殘余沉降較大而中部較充分壓實區(qū)殘余沉降較小，這與前人的研究結果是一致的(張永波等， 2010； 王磊等， 2011； 朱廣軼等， 2014； 郭慶彪， 2017； 陳懷玉等， 2019)。

由圖10可知，隧道的最大殘余沉降出現在隧道中段，殘余沉降量為325mm。由隧道中段向隧道兩端，殘余沉降逐漸變小，并在隧道起點和終點處達到最小，殘余沉降量均小于100mm。在隧道前段(740m附近)出現了一個殘余沉降的“局部峰值”(圖10)，產生這個“局部峰值”的原因在于，隧道在此處從②號工作面邊界的上方穿過，采空區(qū)的邊界空洞區(qū)相較于中部較充分壓實區(qū)具有較大的殘余沉降，而整個隧道的最大殘余沉降出現在②號工作面的另一邊界上方(1150m附近)，該處除了受②號工作面邊界空洞區(qū)的影響，還疊加了③號和④號工作面的影響，使得殘余沉降達到最大。

由圖11可知，在隧道軸線平面內，沿隧道軸線和垂直隧道軸線兩個方向，殘余水平位移均呈現起伏波動的形態(tài)。在隧道起點和終點附近均存在較大的殘余水平位移，最大可至168mm。因隧道與②號工作面的長軸方向呈小角度相交，隧道側向位移空間較大，因此在垂直隧道軸線方向上隧道的位移更為明顯。隧道殘余水平位移最大值出現的位置與殘余沉降“局部峰值”出現的位置基本相同，可以推斷，這兩者形成的原因是一致的。對比圖10和圖11可以發(fā)現，在殘余沉降最大的位置，殘余水平位移并不明顯。這是因為該處的位移受到②、③和④號工作面的共同影響，殘余沉降由3個工作面引起的沉降疊加形成，而水平位移因其引起的位移方向不同而產生互抵。

根據隧道受采空區(qū)殘余變形影響后的潛在狀態(tài)將其分成AB、BC、CD和DE 4段(圖9)，各段的變形特征如圖12所示。在垂直隧道軸線方向的水平位移方面，由于AC段和CE段的潛在水平位移方向相反，隧道將產生一定的彎曲變形，呈“S”形態(tài)； 在平行隧道軸線方向的水平位移方面，AB段和CD段具有向北的位移趨勢，而BC段和DE段具有向南的位移趨勢，這導致隧道在B點和D點附近將產生較大的壓縮變形，而在C點附近將產生較大的拉伸變形； 在沉降變形方面，隧道整體呈下沉趨勢，但沉降并不均勻。沿AC方向隧道的潛在沉降量依次增大，沿CE方向則依次減小，最大沉降將出現在C點附近。同時，由于隧道與沉降等值線斜交(圖9)，沉降中心并不在隧道的正下方，因此隧道將產生一定的軸向扭轉變形。需要指出的是，圖12展示了隧道的變形趨勢，并不代表隧道的真實變形量。

圖12 隧道變形特征分段示意圖

綜上所述，雖然②、③和④號3個工作面停采時間較早，各工作面的殘余沉降空間有限，但由于三者距離較近，相互之間存在影響，加之隧道從它們上方穿越，所以三者對隧道變形和長期穩(wěn)定性的綜合影響較大； 同時由于隧道與上述3個工作面的長軸方向呈小角度相交，隧道與殘余沉降等值線斜交，因此隧道的變形發(fā)生在多個方向上，受力狀態(tài)復雜，呈現側向移動、豎向下沉與軸向扭轉的復合變形趨勢。相反地，盡管①號工作面的殘余沉降較大，但由于其與隧道距離較遠，工作面長軸與隧道近似平行，隧道與工作面殘余變形等值線近似平行，因此①號工作面對隧道變形的影響較為單一而且范圍有限。

5 結 論

現有的變采厚概率積分法模型在進行邊界空洞區(qū)參數選取時，將采厚m看作邊界區(qū)的殘余可活化厚度，導致地表殘余變形量的預測值偏大。本文通過煤層上覆巖層巖性和工作面參數，采用玻茲曼函數確定邊界空洞區(qū)參數，對邊界殘余可活化厚度進行了修正，同時引入時間函數，構建了基于時序的殘余變形預測模型，利用該模型計算了桑掌隧道穿越的4個停采時間不同的采空區(qū)在2020年之后的殘余變形量，并分析了其對桑掌隧道長期穩(wěn)定性的影響，取得如下認識：

(1)工作面殘余沉降沿開采方向呈“W”型，在垂直開采方向呈“V”型，而工作面殘余水平移動均呈起伏波動形態(tài)。

(2)新建桑掌隧道線路上潛在的最大殘余沉降量和最大殘余側向位移分別為325mm和168mm，分別位于隧道中段和隧道起點附近。

(3)隧道的長期穩(wěn)定性受②、③和④號工作面的影響較大，受①號工作面的影響相對較小。隧道的變形呈現出側向移動、豎向下沉與軸向扭轉的復合變形趨勢。

(4)將時間函數引入概率積分法可以使停采時間不同的采空區(qū)在時間尺度上統(tǒng)一，消除沉降歷史的差異性，以便于多采空區(qū)綜合影響下殘余沉降的預測。本文的計算方法可以為分析采空區(qū)殘余變形對隧道的影響提供參考。