巨厚含水松散層下開采地表移動(dòng)變形規(guī)律研究

劉 輝，左建宇，蘇麗娟，程 樺，朱曉峻，張鵬飛，王保國

(1.安徽大學(xué) 資源與環(huán)境工程學(xué)院，安徽 合肥 230601；2.安徽省礦山生態(tài)修復(fù)工程實(shí)驗(yàn)室，安徽 合肥 230601；3.安徽大學(xué) 數(shù)學(xué)科學(xué)學(xué)院，安徽 合肥 230601；4.蚌埠市勘測(cè)設(shè)計(jì)研究院，安徽 蚌埠 233040)

0 引 言

儲(chǔ)量豐富的煤炭資源作為我國最重要的基礎(chǔ)能源，一直占據(jù)我國能源消費(fèi)總量一半以上，擔(dān)負(fù)著為經(jīng)濟(jì)發(fā)展提供主動(dòng)力的責(zé)任[1-4]。我國煤炭資源賦存狀態(tài)具有一定的地域特色，東部礦區(qū)較西部礦區(qū)普遍具有煤層埋藏深、松散層厚且潛水位高等地質(zhì)結(jié)構(gòu)特點(diǎn)，其復(fù)雜的地質(zhì)結(jié)構(gòu)下開采地表移動(dòng)變形較常規(guī)地質(zhì)條件更為劇烈，主要表現(xiàn)為：下沉量增大、下沉系數(shù)接近甚至超過1.0、下沉盆地影響范圍擴(kuò)大、沉陷后易形成積水區(qū)。因此對(duì)礦區(qū)地表生態(tài)環(huán)境和群眾居住生活產(chǎn)生了嚴(yán)重影響[5-8]。

針對(duì)巨厚松散層條件下開采的特殊地表移動(dòng)變形規(guī)律，國內(nèi)外學(xué)者從現(xiàn)場(chǎng)分析、理論分析和實(shí)驗(yàn)?zāi)M等多角度對(duì)其機(jī)理進(jìn)行了大量研究。研究表明巨厚松散層下開采地表移動(dòng)變形的特殊性不僅僅是由上覆巖層綜合巖性較軟所導(dǎo)致的，其不可忽略的因素是煤層開采引起上覆巖層內(nèi)含水層疏水滲流固結(jié)所帶來的附加沉降[9-12]。對(duì)于上覆巖層內(nèi)含水層疏水滲流固結(jié)引起的附加沉降，我國學(xué)者實(shí)驗(yàn)?zāi)M證明了含水層失水固結(jié)機(jī)理，量化分析了含水層固結(jié)沉降量，建立了松散層底部含水層失水固結(jié)模型，取得了一定的成果[13-18]。但是以往的研究手段主要為理論分析和相似材料模擬，難以全面考慮到應(yīng)力場(chǎng)和滲流場(chǎng)的耦合作用對(duì)巨厚含水松散層下開采地表移動(dòng)變形的影響，從而導(dǎo)致結(jié)果的偏差。

以菏澤礦區(qū)某礦典型的巨厚含水松散層地質(zhì)采礦條件為原型，在分析工作面實(shí)測(cè)數(shù)據(jù)的基礎(chǔ)上，采用FLAC3D數(shù)值模擬軟件，充分考慮應(yīng)力場(chǎng)和滲流場(chǎng)耦合作用的影響，對(duì)巨厚含水松散層下煤層開采地表移動(dòng)變形進(jìn)行了流固耦合模擬，總結(jié)了采動(dòng)過程中滲流固結(jié)沉降動(dòng)態(tài)規(guī)律，對(duì)比分析了有無含水層情況下煤層開采地表移動(dòng)變形差異，研究了松散層內(nèi)含水層位置對(duì)地表移動(dòng)變形的影響，為更準(zhǔn)確預(yù)測(cè)及控制巨厚含水松散層下煤層開采引起的地表移動(dòng)變形，提供一定的理論依據(jù)和技術(shù)參考。

1 研究區(qū)概況

菏澤礦區(qū)某礦位于山東省菏澤地區(qū)鄆城縣城南約10 km處，處于巨野煤田的中北部，其范圍東起田橋斷層及田橋支斷層，西至煤系地層底界露頭，南起3925000緯線，北至25勘探線，南北長約14 km，東西寬約11 km，地表總面積達(dá)到了222 km2。礦區(qū)地面標(biāo)高+41.60～+45.38 m，淺層地下水埋深2～7 m，上覆新近系和第四系松散層平均厚度達(dá)到了590 m，局部超過700 m，屬于典型的巨厚含水松散層地質(zhì)條件，礦區(qū)綜合柱狀如圖1所示。

圖1 礦區(qū)綜合柱狀Fig.1 Stratigraphic column of mining area

通過對(duì)鉆孔柱狀圖以及水文地質(zhì)報(bào)告分析可將上覆巨厚松散層概括劃分為5層含水、隔水層，見表1。巨厚松散層內(nèi)第四系含水的砂層與隔水的黏土層相間沉積，屬于富水性中等的松散孔隙含水層，新近系含水的砂層單層呈犬牙交錯(cuò)狀相連，屬于富水性強(qiáng)的松散孔隙承壓含水層。

表1 礦區(qū)含水、隔水層劃分Table 1 Division of aquifers and aquiclude in mining area

某礦1308工作面走向長度約630 m，傾向長度約230 m，主采煤層為山西組3煤層，煤層平均厚度為3 m，平均傾角為3°，為近水平煤層，煤層采用長壁垮落法開采，采深為765 m，工作面上覆松散層平均厚度為582 m，煤層頂板主要由粗砂巖、中細(xì)砂巖、砂質(zhì)泥巖等組成，具有開采深度大、傾角小、松散層厚、松散層內(nèi)含有多層含水層的特點(diǎn)。1308工作面移動(dòng)盆地主斷面上方采用十字線布設(shè)了地表移動(dòng)觀測(cè)線，走向觀測(cè)線全長約1 550 m，累計(jì)布設(shè)47個(gè)觀測(cè)點(diǎn)，傾向觀測(cè)線全長約1 800 m，累計(jì)布設(shè)53個(gè)觀測(cè)點(diǎn)，采動(dòng)期間以平均35 d一次的頻率進(jìn)行觀測(cè)，至停采后共觀測(cè)12次。根據(jù)地表移動(dòng)觀測(cè)數(shù)據(jù)繪制了地表動(dòng)態(tài)下沉曲線(圖2)，地表移動(dòng)變形參數(shù)如下：

工作面1308下沉系數(shù)q1.092水平移動(dòng)系數(shù)b0.25主要影響角正切tan β1.60邊界角δ056.90°

圖2 工作面走向線動(dòng)態(tài)下沉曲線Fig.2 Dynamic sink curve of strike line of working face

可以看出，巨厚含水松散層下開采地表移動(dòng)變形與常規(guī)地質(zhì)條件下采動(dòng)存在一定的特殊性。相比較之下，巨厚含水松散層下開采具有移動(dòng)變形初始期時(shí)間短、起動(dòng)距短、初始期和活躍期下沉速度快、下沉系數(shù)大于1、主要影響角正切和邊界角小、下沉盆地影響范圍大、下沉盆地邊界處水平位移大于下沉量且下沉盆地衰退期長、邊緣收斂緩慢的特點(diǎn)。

2 模型設(shè)計(jì)與方案

2.1 模型設(shè)計(jì)

為研究巨厚含水松散層下開采地表移動(dòng)變形規(guī)律，以菏澤礦區(qū)某礦1308工作面地質(zhì)采礦條件為原型，采用FLAC3D數(shù)值模擬軟件建立了模型尺寸為3 000 m×1 500 m×824 m(長×寬×高)的三維數(shù)值模擬模型。模型單元尺寸在水平方向?yàn)?5 m，垂直方向按照巖層厚度不同進(jìn)行合理設(shè)置。模型采高為3 m，煤層設(shè)計(jì)為水平煤層，設(shè)計(jì)的工作面長度為1 200 m，寬度為400 m，工作面沿走向方向分為12次開采，開挖步長均為100 m，開采深度為765 m，其上覆巖層內(nèi)松散層厚度為582 m。為避免模型過小所導(dǎo)致的邊界應(yīng)力效應(yīng)，在工作面走向和傾向兩側(cè)各留900 m和550 m的邊界保護(hù)煤柱，三維數(shù)值模擬模型如圖3所示。

圖3 數(shù)值模擬模型Fig.3 Numerical simulation model

模型的位移邊界條件設(shè)置為：約束模型底面各個(gè)方向的位移，模型頂面為自由面，前后左右四面約束x、y水平方向的位移，但可以發(fā)生z方向的移動(dòng)。滲流模型設(shè)置為各巖層均為各向同性且均勻等效的連續(xù)孔隙介質(zhì)，力學(xué)模型采用Mohr-Coulomb屈服準(zhǔn)則，數(shù)值模擬巖層參數(shù)根據(jù)地質(zhì)資料進(jìn)行概化合并見表2。

表2 數(shù)值模擬巖層參數(shù)Table 2 Numerical simulation of rock parameters

2.2 模擬方案

為研究巨厚含水松散層下開采地表移動(dòng)變形規(guī)律，明確巨厚含水松散層在下沉盆地形成中的作用以及含水層位置對(duì)地表移動(dòng)變形的影響，分別制定了動(dòng)態(tài)固結(jié)沉降、無含水層、含水層位置3種數(shù)值模擬方案。

1)動(dòng)態(tài)固結(jié)沉降方案。根據(jù)礦區(qū)含水、隔水層劃分表(表1)在松散層內(nèi)設(shè)置含水、隔水層，進(jìn)行流固耦合模擬，研究巨厚含水松散層下采動(dòng)過程中滲流固結(jié)沉降動(dòng)態(tài)規(guī)律。

2)無含水層方案。不考慮巨厚含水松散層在采動(dòng)中產(chǎn)生的流固耦合影響，進(jìn)行純力學(xué)開挖模擬，通過對(duì)比分析，研究有無含水層情況下煤層開采地表移動(dòng)變形差異。

3)含水層位置方案。在松散層內(nèi)設(shè)置厚度為40 m的單層含水層，通過調(diào)整含水層位置標(biāo)高，分別模擬含水層底部與基巖的距離為0、100、200、300、400 m的5種情況，研究巨厚松散層內(nèi)含水層位置標(biāo)高對(duì)地表移動(dòng)變形的影響。

3 巨厚含水松散層對(duì)下沉盆地影響

為明確巨厚含水松散層內(nèi)的多層含水層在煤層開采下沉盆地形成過程中的作用及影響，對(duì)動(dòng)態(tài)固結(jié)沉降方案模擬結(jié)果進(jìn)行歸納總結(jié)分析，并與無含水層方案模擬結(jié)果進(jìn)行對(duì)比分析，研究了巨厚含水松散層滲流固結(jié)沉降動(dòng)態(tài)規(guī)律，并對(duì)比分析了有無含水層情況下煤層開采地表移動(dòng)變形差異。

3.1 巨厚含水松散層滲流固結(jié)沉降動(dòng)態(tài)規(guī)律

動(dòng)態(tài)固結(jié)沉降方案通過循環(huán)開關(guān)流體滲流模塊從而實(shí)現(xiàn)流固耦合模擬：在煤層開挖后，先關(guān)閉FLAC3D中的流體滲流分析模塊，迭代計(jì)算至平衡，得到模型在單力學(xué)場(chǎng)中的土體不排水變形量(開挖變形沉降量)，然后開啟流體滲流模塊，流固耦合計(jì)算土體在該開挖階段內(nèi)的排水變形量(滲流固結(jié)沉降量)，流固耦合計(jì)算完成后再進(jìn)行下一步開挖，并如此往復(fù)循環(huán)計(jì)算，直至工作面開挖結(jié)束。

根據(jù)FLAC3D流固耦合模擬機(jī)制，對(duì)動(dòng)態(tài)固結(jié)沉降方案模擬計(jì)算結(jié)果進(jìn)行分析，繪制了隨工作面推進(jìn)受開挖變形和滲流固結(jié)耦合作用影響的地表最大下沉發(fā)展曲線圖，如圖4所示。從圖4可看出，當(dāng)工作面推進(jìn)至1 200 m開挖結(jié)束時(shí)，地表最大下沉為1.988 m，其中開挖變形沉降在地表總沉降中占主導(dǎo)地位，由開挖變形引起的沉降量為1.820 m，而滲流固結(jié)引起的沉降量為0.168 m，結(jié)合地表總沉降量可得出巨厚含水松散層滲流固結(jié)引起的沉降量占最終地表總沉降的8.5%。

圖4 地表最大下沉發(fā)展曲線Fig.4 Maximum surface subsidence development curve

為研究巨厚含水松散層滲流固結(jié)沉降動(dòng)態(tài)規(guī)律，通過計(jì)算得到了每一步開挖過程中開挖變形沉降增量和滲流固結(jié)沉降增量，并繪制了隨工作面推進(jìn)的地表沉降增量變化曲線，如圖5所示。通過對(duì)地表沉降增量變化曲線進(jìn)行分析可將開挖變形和滲流固結(jié)沉降增量的變化隨工作面推進(jìn)劃分為同步增長期、動(dòng)態(tài)變化期、同步減緩期3個(gè)階段。

圖5 地表沉降增量變化曲線Fig.5 Incremental change curve of surface subsidence

1)同步增長期。在工作面開采初期，工作面推進(jìn)至200 m階段內(nèi)，開挖變形和滲流固結(jié)沉降增量呈現(xiàn)出同步增長的趨勢(shì)，地表開始出現(xiàn)明顯下沉，起動(dòng)距約為1/4煤層埋深。

2)動(dòng)態(tài)變化期。工作面推進(jìn)至200～600 m，開挖變形和滲流固結(jié)沉降增量在整體上共同出現(xiàn)明顯的增大和減小波動(dòng)，其中滲流固結(jié)沉降增量隨工作面推進(jìn)的變化波動(dòng)更為劇烈，與此同時(shí)地表出現(xiàn)劇烈下沉變形。

3)同步減緩期。隨著工作面的繼續(xù)推進(jìn)至終采線1 200 m處，逐步接近走向充分采動(dòng)時(shí)，開挖變形和滲流固結(jié)沉降增量呈現(xiàn)出同步減緩的趨勢(shì)，但其整體上仍存在較小的變化波動(dòng)。

3.2 地表移動(dòng)變形差異對(duì)比分析

為更好地研究巨厚含水松散層的多層含水層對(duì)下沉盆地的影響，對(duì)比分析了動(dòng)態(tài)固結(jié)沉降方案與無含水層方案的模擬結(jié)果，并繪制了工作面推進(jìn)至1 200 m時(shí)動(dòng)態(tài)固結(jié)沉降方案與無含水層方案的工作面走向下沉曲線，如圖6所示。下沉曲線整體關(guān)于采空區(qū)中心對(duì)稱且曲線保持連續(xù)漸變，由于工作面已推進(jìn)至1 200 m，走向達(dá)到充分采動(dòng)，下沉曲線底部呈現(xiàn)平底狀，兩者均符合厚松散層條件下采動(dòng)引起的地表下沉曲線分布，但仍存在一定的差異性。

圖6 工作面走向下沉曲線Fig.6 Strike subsidence curve of working face

通過對(duì)模擬結(jié)果進(jìn)一步對(duì)比分析，得到了地表最大變形，見表3。由表3可看出：考慮松散層內(nèi)含水層影響進(jìn)行流固耦合計(jì)算后，動(dòng)態(tài)固結(jié)沉降方案相較于無含水層方案，最大下沉值增大了0.220 m、最大傾斜增大了0.178 mm/m、最大水平移動(dòng)增大了0.063 m，而最大曲率和最大水平變形分別減小了0.003 mm/m2和0.587 mm/m。

表3 地表最大變形Table 3 Numerical simulation of rock mechanics parameters

表4給出了動(dòng)態(tài)固結(jié)沉降方案與無含水層方案的地表移動(dòng)變形參數(shù)，通過對(duì)比分析可得：

表4 地表移動(dòng)變形參數(shù)對(duì)比Table 4 Comparison table of surface movement deformation parameters

1)考慮松散層內(nèi)含水層影響進(jìn)行流固耦合計(jì)算后，動(dòng)態(tài)固結(jié)沉降方案的下沉系數(shù)達(dá)到了1.096，相較于無含水層方案的0.974，下沉系數(shù)增大了12.5%；相比較之下，水平移動(dòng)系數(shù)變化極小，動(dòng)態(tài)固結(jié)沉降方案的水平移動(dòng)系數(shù)僅比無含水層方案減小了1.0%。

2)對(duì)2種情況的工作面主斷面走向線下沉數(shù)據(jù)分析計(jì)算可得出，動(dòng)態(tài)固結(jié)沉降方案和無含水層方案的主要影響半徑分別為677、641 m，相比較之下考慮含水層影響進(jìn)行流固耦合計(jì)算后，主要影響半徑增大了5.62%。因2種模擬計(jì)算方案的工作面平均采深不變，動(dòng)態(tài)固結(jié)沉降方案的主要影響角正切比無含水層方案減小了5.28%。

3)動(dòng)態(tài)固結(jié)沉降方案和無含水層方案的邊界角分別為48.86°和50.22°，相比較之下動(dòng)態(tài)固結(jié)沉降方案考慮含水層影響進(jìn)行流固耦合計(jì)算后，邊界角減小了2.71%。

綜上所述，巨厚松散層內(nèi)含水層受采動(dòng)影響導(dǎo)致松散層內(nèi)孔隙水壓力消散，有效應(yīng)力增大，產(chǎn)生滲流固結(jié)沉降，引發(fā)下沉盆地在開挖變形沉降的基礎(chǔ)上再平衡[19-20]。從而導(dǎo)致地表下沉和水平移動(dòng)變形增大，主要影響半徑增大，主要影響角正切和邊界角減小，下沉盆地的地表移動(dòng)影響范圍擴(kuò)大，但最大曲率和最大水平變形的減小說明下沉盆地的整體形狀在巨厚松散層內(nèi)含水層疏水滲流固結(jié)沉降的作用下更為平緩[21-22]。

4 地表變形參數(shù)與含水層位置的關(guān)系

上文中對(duì)巨厚含水松散層內(nèi)的多層含水層在下沉盆地形成過程中所產(chǎn)生的影響進(jìn)行了研究，得知巨厚松散層內(nèi)含水層受采動(dòng)影響產(chǎn)生的疏水滲流固結(jié)現(xiàn)象會(huì)導(dǎo)致地表移動(dòng)變形增大、地表移動(dòng)影響范圍擴(kuò)大。為進(jìn)一步探索巨厚松散層內(nèi)不同位置標(biāo)高的含水層對(duì)地表移動(dòng)變形的影響，通過對(duì)含水層位置方案模擬數(shù)據(jù)結(jié)果進(jìn)行分析，研究了地表移動(dòng)變形參數(shù)與巨厚松散層內(nèi)含水層位置之間的關(guān)系。

為方便表述含水層位置方案中所設(shè)置的含水層在松散層內(nèi)所處的位置，定義含水層位置s為含水層底部距基巖的距離與松散層厚度的比值。

4.1 下沉系數(shù)變化規(guī)律

根據(jù)模擬結(jié)果，繪制了下沉系數(shù)q與含水層位置s之間的關(guān)系，如圖7所示。

圖7 下沉系數(shù)和水平移動(dòng)系數(shù)與含水層位置關(guān)系Fig.7 Correlation curve of subsidence factor and displacement factor with aquifer location

從圖7中可得：含水層底部與基巖之間的距離為0時(shí)，下沉系數(shù)為1.098。隨著含水層底部與基巖之間距離的增加，下沉系數(shù)呈現(xiàn)出線性增長趨勢(shì)，當(dāng)含水層底部與基巖之間的距離增加至400 m時(shí)，下沉系數(shù)達(dá)到最大值1.158。對(duì)含水層底部與基巖間距為0～400 m的5種模擬結(jié)果的下沉系數(shù)與含水層位置進(jìn)行相關(guān)性分析，得出下沉系數(shù)q與含水層位置s之間的相關(guān)性公式為：

q=0.087 7s+1.101 2

(1)

綜上可得，下沉系數(shù)與松散層內(nèi)含水層位置之間存在正線性關(guān)系。

4.2 水平移動(dòng)系數(shù)變化規(guī)律

根據(jù)模擬結(jié)果，繪制了水平移動(dòng)系數(shù)b與含水層位置s之間的關(guān)系圖，如圖7所示。

從圖7中可得：5種模擬結(jié)果的水平移動(dòng)系數(shù)為0.287～0.300，數(shù)值波動(dòng)較小。當(dāng)含水層底部與基巖之間的距離為0時(shí)，水平移動(dòng)系數(shù)為0.295，隨著含水層與基巖之間距離的增加，水平移動(dòng)系數(shù)呈現(xiàn)出先增大后減小的變化規(guī)律，在含水層底部與基巖間距為200 m時(shí)達(dá)到0.300的最大值。對(duì)含水層底部與基巖間距為0～400 m的5種模擬結(jié)果的水平移動(dòng)系數(shù)與含水層位置進(jìn)行相關(guān)性分析，得出水平移動(dòng)系數(shù)b與含水層位置s之間的相關(guān)性公式為：

b=-0.073 1s2+0.040 5s+0.293 8

(2)

綜上可得，水平移動(dòng)系數(shù)與松散層內(nèi)含水層位置之間存在先增大后減小的二次函數(shù)關(guān)系。

4.3 主要影響角正切變化規(guī)律

根據(jù)模擬結(jié)果，繪制了主要影響角正切tanβ與含水層位置s之間的關(guān)系，如圖8所示。

圖8 主要影響角正切和邊界角與含水層位置關(guān)系Fig.8 Correlation curve of tangent of major influence tan β and limit angle δ0 with aquifer location

由圖8可得：含水層底部與基巖之間的距離為0時(shí)，主要影響角正切為1.115。隨著含水層底部與基巖之間距離的增加，主要影響角正切呈現(xiàn)出先增大后減小的變化規(guī)律，在含水層底部與基巖間距為300 m時(shí)達(dá)到1.145的最大值。對(duì)含水層底部與基巖間距為0～400 m的5種模擬結(jié)果的主要影響角正切與含水層位置進(jìn)行相關(guān)性分析，得出主要影響角正切tanβ與含水層位置s之間的相關(guān)性公式為：

tanβ=-0.238 3s2+0.165 4s+1.112 3

(3)

綜上可得，主要影響角正切與松散層內(nèi)含水層位置之間存在先增大后減小的二次函數(shù)關(guān)系。

4.4 邊界角變化規(guī)律

根據(jù)模擬結(jié)果，繪制了邊界角δ0與含水層位置s之間的關(guān)系圖，如圖8所示。

從圖8可得：含水層底部與基巖之間的距離為0時(shí)，邊界角為63.77°。當(dāng)含水層位置與基巖不再相鄰時(shí)，邊界角急劇減小至58.54°，并隨著含水層與基巖之間距離的增加，呈現(xiàn)出明顯的線性增長趨勢(shì)，增大至59.37°。對(duì)含水層底部與基巖間距為100～400 m的4種模擬結(jié)果的邊界角與含水層位置進(jìn)行相關(guān)性分析，得出邊界角δ0與含水層位置s之間的相關(guān)性公式為：

δ0=1.600 2s+58.322

(4)

綜上可得，邊界角與松散層內(nèi)含水層位置之間存在先減小后增大的關(guān)系。

5 結(jié) 論

1)某礦巨厚含水松散層采礦條件下受采動(dòng)引起的滲流固結(jié)沉降占地表總下沉量的8.5%；隨工作面推進(jìn)，開挖變形和滲流固結(jié)沉降增量的變化可劃分為同步增長期、動(dòng)態(tài)變化期、同步減緩期3個(gè)階段。

2)巨厚含水松散層采礦條件下，松散層內(nèi)含水層受采動(dòng)影響產(chǎn)生的疏水滲流固結(jié)現(xiàn)象會(huì)導(dǎo)致地表移動(dòng)變形增大、地表移動(dòng)影響范圍擴(kuò)大，但最大曲率和最大水平變形的減小說明下沉盆地的整體形狀更為平緩。

3)地表變形參數(shù)與巨厚松散層內(nèi)含水層位置之間存在相關(guān)性：隨著松散層內(nèi)含水層位置的上移，下沉系數(shù)線性增大，水平移動(dòng)系數(shù)和主要影響角正切先增大后減小，而邊界角先減小后增大。