丁 甲 董東林 胡志強 張隴強 張百偉

(中國礦業(yè)大學(xué)(北京)地球科學(xué)與測繪工程學(xué)院， 北京 100083， 中國)

0 引 言

淮南位于安徽的中北部，該地區(qū)豐富的煤炭資源使其成為了華東地區(qū)乃至全國的能源供給保障。在長期的開采過程中，深部的煤炭資源逐漸枯竭，對該地區(qū)的煤炭開采提出了新的挑戰(zhàn)——從深部逐漸轉(zhuǎn)移至淺部?；茨厦禾餃\部煤層大多覆蓋于巨厚松散層之下，因此位于巨厚松散層底部的高承壓含水層成為了造成礦井水害的主要充水水源之一，根據(jù)以往的記錄，靠近松散層底部開采引起的突水事故也時有發(fā)生。處于巨厚松散層之下且基巖薄弱是淮南地區(qū)煤層淺部開采時面臨的現(xiàn)實條件，開采過程中上覆巖層原有應(yīng)力平衡遭到破壞，在一系列不平衡應(yīng)力的作用下發(fā)生變形破壞，重新平衡后的破碎巖體構(gòu)成了導(dǎo)水裂隙帶，若波及或溝通了松散層底部的承壓含水層，礦井生產(chǎn)將面臨突水或潰砂隱患的威脅(仝騰， 2019)。

口孜東礦位于阜東礦區(qū)東南，行政區(qū)劃隸屬于安徽省阜陽市管轄?？谧螙|煤礦新生界松散層厚度大，平均厚度591.60m，自上而下分為3個含水層4個隔水層。三隔厚度大，分布連續(xù)且穩(wěn)定與底部“紅層”可作為復(fù)合隔水層。基巖風(fēng)氧帶以泥質(zhì)巖類為主，風(fēng)化裂隙多充填泥質(zhì)，隔水性好。開采工作面主要采取綜放開采工藝，礦區(qū)伏在巨厚含水體松散層之下，并且該礦區(qū)松散層底部的含水層水體具有土層結(jié)構(gòu)不一、孔隙水壓高并且其中交替含有相當(dāng)數(shù)量低滲透性的薄夾層等特點(魏繼蓮， 2006)。

該礦區(qū)內(nèi)淺層地下水(120m以淺)受大氣降水和地表水體的補給。上部第二含水層(組)和第三含水層(組)上段是礦區(qū)的主要供水水源。中部第三隔水層(組)黏土厚度大和分布范圍廣，是區(qū)域內(nèi)的重要隔水層(組)。下部第四含水層(組)直接覆蓋基巖含水層之上，在開采淺部煤層時可沿基巖風(fēng)化帶垂直滲入補給，構(gòu)成礦井充水水源之一。因此，對于巨厚松散層水體下煤層開采造成的覆巖破壞規(guī)律認(rèn)識以及煤層頂板導(dǎo)水裂縫帶的發(fā)育高度預(yù)計、預(yù)留煤巖柱的安全性研究評價可以作為重要指標(biāo)來保障礦井的安全、高效生產(chǎn)。結(jié)合我國可持續(xù)發(fā)展的戰(zhàn)略要求，對于安全合理地解放處于松散含水層威脅下的煤炭資源，提高礦井生產(chǎn)帶來的巨大經(jīng)濟(jì)效益具有重要的現(xiàn)實意義(白斌， 2019)。

在具有我國特色的水文地質(zhì)條件下，如何能科學(xué)、合理地預(yù)留防水煤(巖)柱，已經(jīng)引起各個方面越來越多的關(guān)注和重視(張文斌等， 2020)。國內(nèi)對于覆巖穩(wěn)定性和安全煤柱的留設(shè)問題已經(jīng)有許多學(xué)者進(jìn)行了長期的討論，早在多年前錢鳴高就提出了覆巖破壞的關(guān)鍵層理論，認(rèn)為衡量巖層活動主要取決于起主要控制作用的堅硬巖層。因種種條件的限制，使得需要進(jìn)行大量工作的現(xiàn)場探測變得難以實施，運用數(shù)值模擬的方法漸漸走入人們的視野之中，運用數(shù)值模擬的方法能夠更加全面地了解采動過程中的應(yīng)力、位移、塑性破壞變化，并且對于判斷一些難以探測的斷層或者導(dǎo)水裂隙帶的發(fā)育高度也十分有必要。有學(xué)者運用人工神經(jīng)網(wǎng)絡(luò)的方法來對礦區(qū)生產(chǎn)工作面的變化情況進(jìn)行監(jiān)測，同時也根據(jù)變化分析規(guī)律，為預(yù)留煤柱的尺寸確定提供了依據(jù)； 也有采用相似材料模擬的方法對斷層附近防水煤柱的留設(shè)進(jìn)行了分析(施龍青， 2005； 李鵬軍， 2012； 顧大釗等， 2016)； 有學(xué)者采用基于FLAC3D的數(shù)值模擬方法分析了工作面回采時窄煤柱和寬煤柱的應(yīng)力場、位移場及塑性區(qū)的不同特征，從而判斷對煤柱的合理留設(shè)(孔德中等， 2014)，通過對滲流場的分析也可評價滲透破壞發(fā)生的可能性，有學(xué)者通過研究白鶴灘電站滲流場的變化規(guī)律及其水力梯度場的分布，并結(jié)合現(xiàn)場試驗來評價關(guān)鍵區(qū)域內(nèi)發(fā)生滲透破壞的可能性(周志芳等， 2020)； 也有學(xué)者建立了基于FEFLOW的三維地下水流模型并研究不同降雨強度下滲流場的變化，并以此來分析坡體可能發(fā)生的破壞變形規(guī)律(卓萬生， 2020)。

目前在含水層下開采煤層時防水(砂)煤柱的留設(shè)依據(jù)主要是根據(jù)煤層頂板的巖性和物理力學(xué)性質(zhì)進(jìn)行判斷分析，回采方法、采高、煤層傾角(李永明， 2012)等因素也會產(chǎn)生一定影響，但目前的研究大多都是對于覆巖性質(zhì)和回采方法的研究，而未考慮各個含水層產(chǎn)生的變化，亦未有對整個井田進(jìn)行抗?jié)B透破壞的評價，本文分3個部分，通過對巨厚松散層下開采的覆巖破壞規(guī)律、含水層滲流場變化及抗?jié)B透性破壞的綜合研究，為評價煤柱的安全性提供了科學(xué)的方法和依據(jù)。

1 研究區(qū)概況

1.1 區(qū)域構(gòu)造

淮南煤田位于中國華北陸塊東南緣淮南斷褶帶華北聚煤區(qū)南側(cè)安徽省中北部，其以淮南市為主體。煤田北部以劉府?dāng)嗔褳檫吔缗c蚌埠隆起相鄰，煤田南部以阜鳳逆沖斷層為邊界與合肥坳陷相靠，東部以新城口長豐斷裂為邊界，西部以阜陽斷層為邊界。煤田主體構(gòu)造形態(tài)是一個呈北西西向展布的大型復(fù)式向斜結(jié)構(gòu)，在平面上略有彎曲，褶皺軸的西部略有昂起。向斜的兩翼因受到逆沖推覆構(gòu)造作用(圖 1)，有新太古界五河雜巖、中元古界鳳陽群、古生界寒武系-奧陶系地層出露。而在復(fù)式向斜的內(nèi)部，地層平坦開闊，以石炭、二疊系含煤地層為主，是主要的含煤地層。煤層掩埋在新生界松散沉積層之下，傾角都在10°～20°，產(chǎn)狀平緩，南翼推覆斷塊內(nèi)的局部地層傾角陡立、偶呈倒轉(zhuǎn)。煤層由一系列次一級的褶曲組成，形態(tài)寬緩，自北向南有朱集—唐集背斜、尚塘—耿村向斜、陳橋—潘集背斜、謝橋—古溝向斜等。

圖 1 淮南煤田逆沖推覆構(gòu)造剖面示意圖Fig. 1 Schematic diagram of the Huainan coalfield reverse thrust overburden structure section

1.2 水文地質(zhì)概況

口孜東礦為全隱蔽礦床，位于新生界松散層覆蓋下，井田水文地質(zhì)條件分述如下：

1.2.1 含水層及隔水層劃分

本文沿用以往對松散層含隔水層劃分方法，四含(主要為礫石層夾黏土層、局部含砂層)參照淮南煤田松散層劃分原則劃為“紅層”。因此，口孜東礦新生界松散層可分別劃分為第一含水層、第一隔水層、第二含水層、第二隔水層、第三含水層、第三隔水層、第四含水層和“紅層”等8個含、隔水層，分別相當(dāng)于淮南煤田的上含上段、上段隔、上含下段、上隔、中含、中隔和紅層(N1)。

口孜東礦區(qū)內(nèi)松散層的厚度在426.18(88-1)～691.05(28-5)m之間，平均厚度為591.60m，主要受古地形控制，松散層厚度的整體變化趨勢由西向東，由中部向南、北部逐漸減小。東南殘丘處最薄處厚度僅有426.18m，而西邊界中部最厚處可達(dá)691.05m。

1.2.2 礦井充水水源

口孜東礦位于淮南煤田中部含水系統(tǒng)亞區(qū)，南部的阜-鳳逆斷層、北部尚塘集斷層切斷了與周邊的補排關(guān)系，加之近南北向的斜切斷層分割，使得礦井處在一個較為封閉水文地質(zhì)單元中。井田煤系地層上覆厚達(dá)426.18～691.05m的新生界松散層，其中：三隔厚度大，隔水性能良好，可以有效地阻止上覆水體向下入滲補給。

本井田礦井充水水源主要由新生界松散層底部四含水(該含水層不穩(wěn)定，局部可能富水并在底部紅層缺失時可能影響淺部煤層開采)、二疊系煤系砂巖裂隙水、石炭系太灰?guī)r溶裂隙水以及老空(塘)水、斷層水等組成。

新生界松散層底部主要為含礫黏土、礫石層組成的“紅層”，根據(jù)抽水試驗資料，該層具有富水性弱、滲透性差的特點，其上方常有一不穩(wěn)定的薄砂層(口孜東礦劃為“四含”)，介于三隔和紅層之間，補給差，富水性弱，區(qū)域上可合并為相對隔水層。煤礦在淺部開采中，礦井接受補給量大小取決于四含的富水性(底部紅層變薄或缺失時)、基巖風(fēng)化帶的滲透性和巖性特征等因素。因此，在留中上限煤(巖)柱后，對淺部煤層開采一般不會造成威脅。

2 基于FLAC3D 的采動破壞數(shù)值模擬

2.1 設(shè)計模型方案2.1.1 111306 工作面概況

111306工作面走向長度為1965.5～1973.5m，平均走向長度為1969.5m； 可采走向長度為1909～1954m，平均可采走向長度為1936m； 傾斜長度為320m； 標(biāo)高為-642.5～739.3m； 工作面面積為 657264m2，工作面可采面積為 582541m2。

2.1.2 111306 工作面地質(zhì)構(gòu)造

111306 工作面總體構(gòu)造形態(tài)為單斜構(gòu)造，煤層走向近東西，傾向東南，走向上起伏較??； 工作面內(nèi)煤層傾角約為 8°～14°，平均傾角為 11°。

2.1.3 水文地質(zhì)條件

工作面上覆二疊系上石盒子組巖層主要為泥巖、砂質(zhì)泥巖、粉-細(xì)砂砂巖與薄煤層，向上至新生界松散層。上覆新生界厚 559.8～611.3m，平均厚 584.1m。主要由松散砂層組成，含潛水-承壓水，有4個含水層，即一含、二含、三含、四含。

選取靠近工作面附近的補-3鉆孔柱狀圖作為本次模擬的依據(jù)，按比例對其力學(xué)性質(zhì)進(jìn)行加權(quán)分配，鉆孔所取巖樣如圖 2所示。

圖 2 口孜東礦鉆孔巖樣Fig. 2 Borehole rock sample

用巖石力學(xué)特征試驗機對所取巖樣進(jìn)行各種力學(xué)參數(shù)的測試，并記錄數(shù)據(jù)(泊松比、抗拉強度、密度、抗剪強度等)，建立的具體模型方案如表 1所示。

表 1 模型方案分配表Table 1 Model program allocationTable

為了能準(zhǔn)確細(xì)致地分析了土-巖接觸面附近的情況，對其劃分了更小的網(wǎng)格，以便能更準(zhǔn)確地進(jìn)行分析。

2.2 建立初始模型

煤層采厚為4.3m，在采空區(qū)四周留設(shè)足夠的煤柱，防止因為采空區(qū)影響范圍太大而產(chǎn)生的誤差，模型尺寸為500m×500m×688.2m。設(shè)置當(dāng)前條件下的重力加速度為9.8m·s-2，側(cè)向地應(yīng)力與垂向應(yīng)力之比為0.5。添加含水層相應(yīng)的水位面，給出流體密度，并對水位面以上、下單元設(shè)置干密度，模型會自動計算含水層水壓產(chǎn)生的影響。對整體模型單元分配彈塑性莫爾-庫侖本構(gòu)模型，并對每層賦予模型方案中相應(yīng)的物理力學(xué)參數(shù)； 對模型底面邊界進(jìn)行全約束，即水平、垂直方向都不允許產(chǎn)生位移； 頂部為自由邊界，各個方向都可產(chǎn)生位移； 并且對模型剩余4個面進(jìn)行水平方向的速度約束，即不允許產(chǎn)生水平方向的位移； 原始地層模型如圖 3所示。

圖 3 原始地層模型圖Fig. 3 Original stratigraphic model

2.3 模型運行結(jié)果分析

對運行平衡之后的位移云圖進(jìn)行分析，采空區(qū)上方的直接頂在自身重力和周圍巖體的擠壓下第一時間產(chǎn)生了垂向位移且距離較大，隨著模型進(jìn)一步平衡，進(jìn)行垂向位移的巖體單元范圍越來越大，當(dāng)巖體單元失去下方的應(yīng)力支撐時將會在上方不平衡應(yīng)力的作用下產(chǎn)生向下的位移，直到重新建立應(yīng)力平衡，而這種不平衡將會從采空區(qū)開始從下向上傳遞直至地表。

觀察圖 4發(fā)現(xiàn)采空區(qū)范圍上方的中心區(qū)域整體分帶性較強，但隨著影響范圍的擴(kuò)大，更上方的單元分帶性已經(jīng)漸漸變得不明顯，并且中心區(qū)域的巖體單元位移變化不大，證明此區(qū)域已經(jīng)形成了新的應(yīng)力平衡，這為區(qū)分垮落帶和裂縫帶的高度提供了依據(jù)(楊達(dá)明等， 2019)。

圖 4 工作面推進(jìn)200m垂向位移云圖Fig. 4 Contour map of 200m vertical displacement

圖 5 最大主應(yīng)力云圖Fig. 5 Contour map of maximum principal stress

分析圖 5可以看出，煤層回采時覆巖原本維持的應(yīng)力平衡遭到破壞，在回采初期，采空區(qū)上方出現(xiàn)紅色的卸壓區(qū)(王宏偉等， 2014)，因頂板巖石垮落導(dǎo)致此區(qū)域應(yīng)力顯著降低，但影響的區(qū)域十分有限。煤柱位置處高密度的集中應(yīng)力是導(dǎo)致煤柱塑性損傷的主要原因(張文斌等， 2020)，密集卸壓區(qū)的影響范圍和大小逐漸平衡之后，這種不平衡的應(yīng)力逐漸向上方傳遞，頂板卸壓區(qū)的影響范圍逐漸增大，并且此時的卸壓區(qū)呈現(xiàn)出一定的分帶性，在紅色的密集卸壓區(qū)造成的導(dǎo)水裂隙帶之上出現(xiàn)了黃色的應(yīng)力影響區(qū)域，這也合理地解釋了即使導(dǎo)水裂隙帶沒有直接溝通含水層，但是由于受到應(yīng)力的破壞，其上方的巖土體依然受到了影響，最明顯的變化就是在此作用下的巖土體垂向滲透系數(shù)會呈數(shù)倍增大，以此應(yīng)力影響的大小和范圍作為在GMS模型中開采造成的調(diào)整變量； 隨著模型的進(jìn)一步平衡，應(yīng)力影響的范圍和影響作用已經(jīng)逐漸趨于穩(wěn)定，這種最大主應(yīng)力分布的規(guī)律表明此時導(dǎo)水裂隙帶的發(fā)育已經(jīng)趨于平穩(wěn)，其最大高度不會隨著回采的推進(jìn)而繼續(xù)增大，用此時的導(dǎo)水裂隙帶高度來衡量煤柱的安全性將會十分可靠。

煤層開采形成采空區(qū)之后，在卸壓區(qū)范圍內(nèi)頂板失去下部支撐力，直接頂巖體單元在自身重力、頂部其他巖體單元的垂向壓應(yīng)力和側(cè)向巖體單元的側(cè)向擠壓應(yīng)力作用下發(fā)生崩塌垮落現(xiàn)象(王靜民， 2021)，巖體單元強度無法對抗拉應(yīng)力破壞和剪切破壞，當(dāng)垮落至采空區(qū)底面時重新獲得下部支撐應(yīng)力進(jìn)行平衡，且這種平衡會隨時間逐漸擴(kuò)大其影響范圍，直至巖體單元的抗拉強度可以支持巖體對抗上部傳遞的壓應(yīng)力且不發(fā)生垮落現(xiàn)象，此處即為垮落帶的發(fā)育高度。在垮落帶之上的巖體單元雖然可以承受上部壓力不發(fā)生垮落，但是會在垂向應(yīng)力和側(cè)向應(yīng)力的共同作用下發(fā)生剪切破壞，剪切破壞的高度也隨著時間推移而增加，但在向上發(fā)育的過程中最大主應(yīng)力逐漸減小，其受到剪切破壞的程度也逐漸減小，裂縫帶中的空隙會隨著高度增加越來越小，直至剪應(yīng)力無法對巖體單元造成剪切破壞，此時即為導(dǎo)水裂隙帶發(fā)育的最大高度，且此高度不會隨著工作面的繼續(xù)推進(jìn)而增大。如圖 6所示采空區(qū)上方出現(xiàn)的深紫色區(qū)域代表過去受拉和過去受剪切破壞并且現(xiàn)在依然受到剪切破壞，結(jié)合之前的垂向位移云圖和最大主應(yīng)力圖可以判斷深紫色區(qū)域頂部即為垮落帶發(fā)育的最大高度，經(jīng)測量其到采空區(qū)的高度為18.7m； 深紫色區(qū)域上方出現(xiàn)的影響范圍較大的綠色區(qū)域為過去受到剪切破壞狀態(tài)，判斷其符合裂縫帶的發(fā)育特征，測量其高度為56.6m。

表 2 含水層水文地質(zhì)參數(shù)表Table 2 Table of hydrogeological parameters of aquifer

圖 6 單元受力狀態(tài)圖Fig. 6 Unit force state image

本次模擬建立在礦區(qū)真實的條件之下，首先進(jìn)行模型方案設(shè)計，對初始地層進(jìn)行初始地應(yīng)力平衡，通過聯(lián)合垂向位移云圖、最大主應(yīng)力云圖、巖體單元受力狀態(tài)圖分析了開采后頂板覆巖破壞的規(guī)律，判斷了導(dǎo)水裂隙帶的發(fā)育高度，明確了采動后應(yīng)力的影響范圍和大小，為采動前后的滲流場變化和評價留設(shè)煤柱的安全性提供了參考依據(jù)。

3 松散層滲流場地下水?dāng)?shù)值模擬

為了比較深入研究口孜東煤礦的松散層底部含水層的孔隙水對工作面頂板的影響，根據(jù)口孜東礦的水文地質(zhì)資料及抽水鉆孔資料，采用GMS軟件，分析松散層的地下水滲流場分布，以及在不同滲透系數(shù)分區(qū)的條件下的水位情況。

3.1 水文地質(zhì)條件概化

水文條件的概化就是通過編輯各巖層滲透系數(shù)等水文地質(zhì)參數(shù)(表 2)實現(xiàn)研究區(qū)含水層的概化。對流動系統(tǒng)地各項地質(zhì)條件，如地層的標(biāo)高、邊界條件、水動力聯(lián)系特征、各個地層的特定水文地質(zhì)參數(shù)等進(jìn)行抽象和概化。水文地質(zhì)條件的概化是數(shù)值模擬的基礎(chǔ)，因此建立地下水滲流場模型應(yīng)該充分地收集并且分析研究區(qū)的各種資料和信息。

綜合分析研究區(qū)域的含水層為松散層的底部的含水層，考慮到研究區(qū)域的埋藏條件、水力聯(lián)系、含水層的厚度、以及介質(zhì)的物理性質(zhì)，將研究區(qū)垂向自上而下劃分為8層，共4層含水層和4個隔水層。每一層的埋深由礦區(qū)的鉆孔水文資料可以得出。

3.2 地下水運動系統(tǒng)的概化

自然狀況下，整個區(qū)域的水流呈東西流向，南北兩邊的采區(qū)邊界都有逆斷層作為隔水邊界。該地區(qū)的鉆孔資料顯示，第三含水層為由東向西流，而第四含水層水位則為東北高，西南低。且第一、二、三含水層透水性能較好，分布均勻； 隔水層黏土分布穩(wěn)定，厚度較厚，隔水層能很好地阻隔上層的水流下滲。因此將研究區(qū)的地下水含水系統(tǒng)流動性概括為：多層位的非均質(zhì)各向同性準(zhǔn)三維非穩(wěn)定流系統(tǒng)。

3.3 地下水含水系統(tǒng)的概化

研究區(qū)以口孜東礦區(qū)為主要研究區(qū)域，根據(jù)地下空間的分布，地下水流以及鉆孔資料將本次研究模型的側(cè)向邊界設(shè)定為不同特性的一般流量邊界以及定流量邊界。詳述如下：

根據(jù)探查報告可知，在平面上，口孜東礦區(qū)南邊有阜-鳳逆斷層，北部的尚塘集斷層很好地切斷了與周邊地補給關(guān)系，設(shè)阜-鳳斷層為隔水邊界。第四含水層中，由于部分孔壁縫隙下滲水流較多設(shè)北部為補給邊界，南部斷層附近有部分抽水井，且抽水量不小，作為排泄邊界。東西方向的水位數(shù)值相差較小，水量交換少，因此將東西方向的邊界確定為弱透水邊界。

垂向上，頂部含水層接受大氣降水跟地表水體的入滲補給，概化為有水量交換的開放邊界； 底部紅層與基巖接觸，且以半固結(jié)狀態(tài)泥質(zhì)充填，中間有棕紅色黏土、砂質(zhì)黏土分布，可視為無水量交換的隔水邊界。

3.4 水文地質(zhì)結(jié)構(gòu)可視化模型的構(gòu)建

鑒于口孜東礦區(qū)基礎(chǔ)數(shù)據(jù)資料的收集情況，采用Borehole模塊并通過以下3個步驟來構(gòu)建三維地質(zhì)模型。

(1)輸入鉆孔數(shù)據(jù)。

(2)編輯鉆孔橫截面。

(3)構(gòu)建三維實體地質(zhì)模型。

所構(gòu)建的實體模型如圖 7所示。

圖 7 三維地層實體模型圖Fig. 7 Three-dimensional stratigraphic solid model diagram

圖 8 開采前三、四含水層地下水流場分布圖Fig. 8 Distribution of groundwater flow field in the third and fourth aquifers before mining

3.5 水文地質(zhì)概念模型的建立

定義第一含水層為非承壓含水層，其他3個含水層為承壓含水層，并且為每個含水層參數(shù)賦值，包括滲透系數(shù)或?qū)禂?shù)、給水度等； 確認(rèn)areal屬性，為降水入滲補給量和蒸發(fā)排泄量賦值，本次模擬涉及的源匯項主要有河流、入滲補給、蒸發(fā)以及通用水頭等； 確定邊界的補給排泄及賦值，對于定水頭邊界，通過定義弧段端點水頭高程來定義定水頭值，弧段其他地方水頭值以線性插值自動給出，然后對模擬計算區(qū)域進(jìn)行網(wǎng)格剖分。

導(dǎo)入第一含水層的頂板高程及各含水層的底板高程數(shù)據(jù)，再導(dǎo)入初始流場的水位點及水位值，對導(dǎo)入的點進(jìn)行插值計算，然后分別賦值給各含水層和初始流場，由GMS自動賦值到每個節(jié)點，初步建立起模型模擬的初始條件。運行MODFLOW進(jìn)行地下水模擬，模擬運算的結(jié)果即以各含水層地下水位等值線的形式如圖 8所示(由于對第一、二含水層幾乎不造成影響，僅用第三、四含水層水流場作為說明)。

在模型中確定工作面的坐標(biāo)點，根據(jù)FLAC3D的應(yīng)力影響范圍邊緣點在GMS模型中找到相對應(yīng)的影響范圍，根據(jù)應(yīng)力改變模型中相應(yīng)影響范圍內(nèi)單元的垂向滲透系數(shù)，如圖 9所示。

圖 9 彎曲破壞覆巖滲透系數(shù)調(diào)整圖Fig. 9 Vertical permeability coefficient adjustment chart

由于巖體當(dāng)前所處的應(yīng)力狀態(tài)將直接影響到巖體的孔隙率從而進(jìn)一步影響巖體的滲透率和滲透系數(shù)，根據(jù)公式可以得到所受應(yīng)力狀態(tài)對滲透系數(shù)大小的影響；

(1)

式中：k0為0應(yīng)力狀態(tài)的滲透率；K為滲透系數(shù);γ為流體的容重;μ為流體的動力黏度系數(shù);φ0為0應(yīng)力狀態(tài)下的孔隙。

3.6 模型校核

模型校核是地下水?dāng)?shù)值模擬中重要的環(huán)節(jié)，任何一個用于預(yù)測的地下水模擬模型，都必須證明其對地下水系統(tǒng)的模擬是正確的。模型校核的目的是不斷系統(tǒng)地調(diào)整各項參數(shù)，使得計算水頭值與觀測值盡可能接近。將本次模擬的置信度設(shè)置為95%，本次模擬得出三含、四含的觀測孔校核情況如圖 10示。

從圖中可以看到，在置信度為95%的前提下，本次模擬的觀測值和模擬值的條形顯示基本都顯示為綠色，即兩個數(shù)值都在置信區(qū)間的范圍之內(nèi)，認(rèn)為本次模擬結(jié)果是可靠的。

圖 10 開采后第三、四含水層觀測孔校核圖Fig. 10 Third and fourth aquifer observation hole calibration after mining

圖 11 第三、四含水層計算值和觀測值擬合圖Fig. 11 Fitted graphs of calculated and observed values of the third and fourth aquifers

3.7 模型參數(shù)反演

本次模擬我們在前面模型校核的基礎(chǔ)上選擇自動參數(shù)反演的方法自動調(diào)參。對反演參數(shù)區(qū)域進(jìn)行分區(qū)，賦以滲透系數(shù)、大氣降水入滲率速率等參數(shù)不同的指針值。在相應(yīng)模塊中反演各個區(qū)的名稱、類型、指針值、初值、最大值和最小值。軟件通過調(diào)整迭代次數(shù)，不斷調(diào)整參數(shù)使得觀測值和計算殘差最小，最終實現(xiàn)迭代優(yōu)化，如果不能實現(xiàn)收斂，則改變初始值再次進(jìn)行運算，直到計算值與觀測值達(dá)到殘差最小為止。

最后運行MODFLOW，查看觀測井的擬合情況，如果擬合程度高，觀測值和計算值會散落在置信斜線周圍，如果誤差大，散點圖會很分散。

從圖中可以發(fā)現(xiàn)觀測值和計算值的散點均分布在置信斜線附近，通過模型校核和模型參數(shù)反演兩個步驟的檢驗，認(rèn)為本次數(shù)值模擬的誤差較小，模擬結(jié)果可靠。

3.8 滲流場地下水?dāng)?shù)值模擬結(jié)果分析

綜合分析各個不同含水層的滲流場分布圖及觀測數(shù)據(jù)擬合圖可以得出礦區(qū)的地下水水頭和水位整體埋深情況。結(jié)果表明，模擬結(jié)果跟觀測數(shù)據(jù)高度擬合，相關(guān)系數(shù)為0.975，其中地下水流動方向為由北向南與資料中保持一致，在礦區(qū)北部及中部地勢較高處地下水水頭較高，水位埋深相對較深，接受大氣降雨入滲補給后由高處向地形低處徑流； 礦區(qū)東南殘丘處地勢相對較低，地下水水頭較低，地下水埋深較淺，局部地方有水體出露。

三維水文地質(zhì)模型再現(xiàn)了礦區(qū)新生界松散層不同巖性的空間分布情況和砂、黏土相互交替的沉積序列。通過模型可以觀察出松散層總體變化趨勢為由東向西，由南、北向中部逐漸增厚，頂部含水層以粉細(xì)砂巖為主，砂層和黏土層發(fā)育厚度極不穩(wěn)定，其中第一隔水層以灰黃-褐黃色砂質(zhì)黏土為主，組厚8.7～22.6m，平均16.4m，其中土層厚4.50～20.90m，平均11.77m，局部地段夾薄層粉細(xì)砂，分布較穩(wěn)定，能起隔水作用。上部第二含水層和第三含水層上段是礦區(qū)的主要供水水源。中部第三隔水層黏土厚度大和分布范圍廣，是區(qū)域內(nèi)的重要隔水層。下部第四含水層直接覆蓋基巖含水層之上，在開采淺部煤層時可沿基巖風(fēng)化帶垂直滲入補給，構(gòu)成礦井充水水源之一。

對比開采前后的各含水層的水流場分布圖，可以發(fā)現(xiàn)流場的變化非常微小，由此說明即使導(dǎo)水裂隙帶上方巖體的彎曲破壞使其垂向滲透系數(shù)明顯增大，但是其影響范圍依然無法波及到含水層位置，生產(chǎn)的安全性也可以保證。

4 抗?jié)B透性破壞評價研究

4.1 涌水潰砂機理

煤礦工作面在臨近含水松散層下生產(chǎn)推進(jìn)的過程中，砂混入水中形成混合流體，大量的潰入井下工作面，從而造成人員傷亡以及財產(chǎn)損失的一種礦井地質(zhì)災(zāi)害現(xiàn)象被稱為突水潰砂。

位于煤系上部的自由水體是突水潰砂發(fā)生的源動力，頂板覆巖遭到應(yīng)力破壞后其中心區(qū)域水力梯度將明顯增大，當(dāng)水力梯度超過了砂土混入水中的臨界條件時，裂隙通道附近砂層會在較大的動水壓力下對裂縫進(jìn)行滲透破壞，從而使得水砂更易沿著較大裂縫或通道潰入采空區(qū)(李建文， 2016)。

4.2 水動力條件分析

發(fā)生潰水涌沙時，可將含水層隨裂隙通道充水的過程當(dāng)作是承壓含水層中單井進(jìn)行定流量抽水(張玉軍等， 2006)。將完整的充水過程看作是地下水流向井的非穩(wěn)定流，此時裂隙通道內(nèi)或者其影響范圍內(nèi)任一點處的水力坡度就可概括為單井定流量承壓井流(李建文， 2016)，可歸納為如下數(shù)學(xué)模型：

(2)

此時任意時刻的任意一點的水位降升s，只要在其抽水產(chǎn)生的影響范圍之內(nèi)即滿足下列公式：

(3)

其中：

(4)

在地下水動力學(xué)中，采用井函數(shù)W(u)代替式中的指數(shù)積分式：

(5)

則就可改寫為：

(6)

式中：Q為流量；r為測點到抽水井位置的直線距離(m)；t為時間；T為導(dǎo)水系數(shù)；s為貯水系數(shù)。

由于突水潰砂發(fā)生時水頭下降至含水層底板，故水位降深sw=Ho(以含水層底板起算)，此時

hw=Ho-sw=0

(7)

因此，突水潰砂的涌水量公式即為：

(8)

式中：Q為涌水量(L·s-1)；K為滲透系數(shù)(m·d-1)；M為含水層厚度(m)；H為含水層水頭高度(m)；rw為裂隙通道半徑(m)；R為發(fā)生潰砂時的影響半徑(m)。

當(dāng)r約等于rw時，便可推導(dǎo)出當(dāng)前水力坡度的計算公式(裂縫通道內(nèi)部)：

(9)

4.3 潰砂發(fā)生的判別依據(jù)

當(dāng)采動產(chǎn)生的裂隙聯(lián)通上覆砂土含水層時，原本承壓的水體通過所產(chǎn)生的導(dǎo)水裂隙潰入工作面之中，最低水位也隨之迅速降低至與含水層底板相同的高度，此時產(chǎn)生的水力坡度也達(dá)到峰值，此峰值若是超過了其相應(yīng)臨界水力坡度時，位于含水層之中的砂土小顆粒將伴隨水流一并進(jìn)入工作面，這時在滿足條件J≥Jcr的條件下會發(fā)生潰砂現(xiàn)象。

由此可推導(dǎo)出發(fā)生潰砂時的臨界水頭高度Hcr的表達(dá)式為：

(10)

由式(10)可知，當(dāng)前含水層水頭的高度若小于臨界水頭高度時，突水潰砂的威脅將不存在； 相反，工作面生產(chǎn)將面臨潰砂的危險。

4.4 涌水通道半徑的確定

由之前關(guān)于潰水涌砂的機理可知充水過程可看作是以裂隙通道作為近似半徑rw的承壓含水層單井定流量抽水，此時，假設(shè)裂隙通道為圓形截面，其半徑rw作為此承壓完整井的井徑?？衫昧严秾挾葋韺ζ溥M(jìn)行計算(孟召平等， 2016)：

D=Δhcosα

(11)

Δh=M-∑h(Kp-1)

(12)

式中：D為裂隙寬度(m)； Δh為頂板臺階下沉量(m)；α為頂板基巖垮落角/(°)；M為采高(m)； ∑h為直接頂厚度(m)；Kp為垮落后的碎脹系數(shù)。

4.5 數(shù)據(jù)及結(jié)果的計算

由抗?jié)B透性破壞評價的機理可知，對實際水力坡度進(jìn)行計算時，必須要首先得到目標(biāo)層位的影響半徑以及滲透系數(shù)。根據(jù)已有資料的四含抽水試驗成果表以及新生界紅層補勘抽水試驗成果表對四含以及紅層進(jìn)行計算。

由于抽水試驗數(shù)據(jù)中只有一次水位降升數(shù)據(jù)，故在對其參數(shù)進(jìn)行計算時常用承壓穩(wěn)定流單井抽水的裘布依公式與庫薩金公式聯(lián)立后，給定初始的經(jīng)驗影響半徑值R，然后聯(lián)立迭代多次求解出R與滲透系數(shù)K，當(dāng)其誤差小于百分之五時，即得到計算合理的影響半徑值與滲透系數(shù)K值。這種計算方法雖然比起有觀測井的多次降升抽水試驗數(shù)據(jù)計算出來的結(jié)果誤差更大，但由于影響半徑在裘布依公式中是以lnR/rw的形式出現(xiàn)的，所以，該誤差對于滲透系數(shù)的計算所產(chǎn)生的影響并不大。

(13)

(14)

根據(jù)資料，利用計算機建立迭代計算程序，計算得出的結(jié)果如表 3與表 4所示。

表 3 四含抽水試驗計算結(jié)果表Table 3 The fourth aquifer contains pumping test calculation resultsTable

表 4 新生界紅層補勘抽水試驗計算結(jié)果表Table 4 Cenozoic red beds replenishment survey pumping test calculation resultsTable

對比計算與補勘鉆孔計算結(jié)果，吻合度較高。

涌水通道半徑的計算由式(11)與式(12)確定，由于D為倒水通道寬度，在假設(shè)導(dǎo)水通道界面為圓形的情況下，D則為導(dǎo)水通道直徑。由FLAC3D覆巖破壞數(shù)值模擬結(jié)果給定垮落角α為60度，根據(jù)資料，取采高為4.3m。將表 3與表 4中所計算出來的單井抽水的各個參數(shù)帶入式中得出各點位的模擬突水潰砂時的涌水通道半徑范圍如表 5與表 6所示。

表 5 四含涌水通道半徑取值范圍Table 5 The fourth aquifer gushing water channel radius taking range

表 6 紅含涌水通道取值范圍Table 6 The red beds gushing water channel radius taking range

根據(jù)含水層巖性，四含選擇黏土土粒密度均值，取2740kg·m-3。紅層密度則根據(jù)鉆孔資料取平均值2688.2kg·m-3，由此分別計算出四含與紅層的臨界水力梯度分別為：

Jcr四含=1.66

Jcr紅層=1.62

其中：孔隙率n取值為鉆孔數(shù)據(jù)均值，其值較低，約為6%。

表 7 四含臨界水頭高度計算值Table 7 Calculated value of critical head height of the fourth aquifer

表 8 紅層臨界水頭高度計算值Table 8 Calculated value of critical head height of the red beds

5 結(jié) 論

(1)通過對工作面附近鉆孔柱狀圖的分析，建立模型地層層序，并對工作面覆巖的物理力學(xué)性質(zhì)進(jìn)行收集和實測，在此基礎(chǔ)上建立基于FLAC3D的開采覆巖破壞數(shù)值模擬，運行模型達(dá)到平衡后通過分析單元狀態(tài)圖和導(dǎo)水裂隙帶的形成機理判斷其發(fā)育高度，垮落帶的高度為18.7m，垮落帶和裂縫帶的總體發(fā)育高度為56.6m，當(dāng)前開采條件下留設(shè)的煤柱高度遠(yuǎn)高于此值，導(dǎo)水裂隙帶無法波及到主要充水含水層； 在應(yīng)力的影響下覆巖產(chǎn)生變形，使水力坡度和巖體垂向滲透系數(shù)發(fā)生了變化，通過分析最大主應(yīng)力圖確定開采造成的應(yīng)力破壞影響范圍和大小。

(2)根據(jù)口孜東礦的水文地質(zhì)資料及抽水鉆孔資料，建立了基于GMS的松散層滲流場地下水?dāng)?shù)值模型，并通過模型校核和模型參數(shù)反演確定了本次模擬的可靠性，通過模型可以觀察出松散層總體厚度變化趨勢為由東向西，由南、北向中部逐漸增厚，下部第四含水層水體直接覆蓋于基巖含水層上方，在開采淺部煤層時可沿基巖風(fēng)化帶垂直滲入補給，構(gòu)成礦井充水水源之一。調(diào)整應(yīng)力影響范圍內(nèi)網(wǎng)格單元的垂向滲透系數(shù)，對比開采前后的滲流場變化，可以發(fā)現(xiàn)雖然覆巖發(fā)生破壞形變后在一定范圍內(nèi)的垂向滲透系數(shù)突增，但是各個含水層的滲流場變化卻非常微小，說明了當(dāng)前條件下留設(shè)的煤柱高度不僅可以保證導(dǎo)水裂隙帶無法發(fā)育至含水層，且導(dǎo)水裂隙帶上方巖體的彎曲變形造成的影響也極少波及到含水層，認(rèn)為當(dāng)前煤柱的安全性可以保證。

(3)根據(jù)抗?jié)B透性破壞評價的機理，對各個影響因素進(jìn)行分析，并以此作為突水潰砂的判斷依據(jù)，根據(jù)已有資料的四含抽水試驗成果表以及新生界紅層補勘抽水試驗成果表對四含以及紅層進(jìn)行計算，通過計算機迭代計算得到了四含及紅層的抽水試驗計算結(jié)果，且結(jié)果吻合度較高。由FLAC3D覆巖破壞數(shù)值模擬結(jié)果給定垮落角α為60°，將計算出來的單井抽水的各個參數(shù)帶入式中得出各點位的模擬突水潰砂時的涌水通道半徑范圍，并確定了第四含水層和紅層的臨界水力坡度，并在此基礎(chǔ)上計算了以多處鉆孔為依據(jù)的第四含水層和紅層的臨界水頭高度的都介于10～25m之間，初始水位都介于10～20m之間，但是第四含水層與紅層的實際水頭值遠(yuǎn)遠(yuǎn)大于該處突水潰砂的臨界水頭值，從而也可以保證生產(chǎn)的安全。