區(qū)域注漿擾動(dòng)下滲流場?化學(xué)場演化及耦合作用

郭 艷 ，桂和榮 ，魏久傳 ，胡滿聰 ，郭祥東 ，聶 鋒 ，陳永青 ，解 建 ，葉 爽 ，李 俊

（1.國家煤礦水害防治工程技術(shù)研究中心(宿州學(xué)院), 安徽 宿州 234000；2.山東科技大學(xué) 地球科學(xué)與工程學(xué)院, 山東 青島 266590；3.礦井水資源化利用安徽普通高校重點(diǎn)實(shí)驗(yàn)室, 安徽 宿州 234000；4.淮北礦業(yè)股份有限公司, 安徽 淮北 235000；5.安徽恒源煤電股份有限公司, 安徽 宿州 234000）

0 引 言

我國東部煤系賦存區(qū)地質(zhì)構(gòu)造及水文地質(zhì)條件較為復(fù)雜，是世界上煤礦水害最為嚴(yán)重的區(qū)域之一。隨著淺部易采煤炭資源逐漸枯竭，深部煤炭資源將成為我國煤礦開采的主攻方向。深部資源開采面臨高地壓、高水壓、高地溫等問題，地質(zhì)及水文地質(zhì)條件復(fù)雜。隨著深部煤炭資源的開發(fā)強(qiáng)度逐漸增大，面臨的水害問題也愈加嚴(yán)重。據(jù)初步統(tǒng)計(jì)，全國重點(diǎn)煤礦中受水害威脅的礦井占47.5%[1-4]。

從煤礦水害事故統(tǒng)計(jì)分析看，華北煤田95%以上的特大型淹井事故均因底板巖溶水突水所致。如桃園礦2013 年的“2?3”奧灰突水事故，最大突水量29 000 m3/h，死亡1 人[5-6]。

長期以來，底板巖溶水害防治一直是華北煤田煤礦防治水的重點(diǎn)。特別是煤礦向深部延伸開采，高地壓、高水壓、高地溫條件下，巖溶水突水威脅更大。目前，在淮北煤田普遍采用地面定向鉆高壓注漿技術(shù)，對(duì)太原群上組第三層薄層灰?guī)r（簡稱“三灰”）進(jìn)行區(qū)域性加固改造，在增加底板隔水層厚度和強(qiáng)度的同時(shí)，對(duì)穿過三灰?guī)r溶裂隙及導(dǎo)水通道（如導(dǎo)水?dāng)鄬拥龋┻M(jìn)行有效封堵[7-9]。

底板區(qū)域注漿工程實(shí)施中，注漿高壓以及漿液析水會(huì)在一定時(shí)期內(nèi)改變地下水流場的分布特征。多采區(qū)及多工作面區(qū)域高壓注漿改造三灰，補(bǔ)徑排路徑遭到人為干擾，使太灰水及與之有水力聯(lián)系的含水層滲流場發(fā)生變化。滲流場的變化必然會(huì)造成水化學(xué)場的改變。同時(shí)漿液析出水水量大（室內(nèi)試驗(yàn)可知，煤礦常用漿液相對(duì)密度1.2～1.4，析水率為56%～25%）；pH 值通常為10～12，甚至更高，地下水經(jīng)強(qiáng)堿性水泥漿液堿性作用后，必然造成含水層的水化學(xué)組成變化。大規(guī)模區(qū)域注漿，對(duì)地下水化學(xué)場和滲流場如何影響？耦合作用怎樣？這些均是亟待解決的科學(xué)問題。

在礦區(qū)地下水化學(xué)場?滲流場（簡稱“雙場”）耦合研究方面，已有研究主要立足于特征離子和長觀孔水位變化[10-11]。近年來，有大量文獻(xiàn)報(bào)道了利用TDS、同位素等分析流場特征[12-14]。許多軟件（如GMS、Feflow、Visual MODFLOW 等）在模擬采動(dòng)影響下地下水流場演變方面發(fā)揮了重要作用[15-18]。目前，關(guān)于滲流場與化學(xué)場耦合作用研究，基本都是基于原生自然環(huán)境下進(jìn)行，關(guān)于人為防治水工程、突水災(zāi)害等擾動(dòng)下的相關(guān)研究少有涉獵，特別是關(guān)于目前常用的地面定向鉆高壓注漿擾動(dòng)下的“雙場”耦合相關(guān)研究較少。

Feflow 軟件較適用于裂隙介質(zhì)三維地下水模擬，且能對(duì)斷層加密剖分，更適用于礦區(qū)地下水模擬[19]。筆者基于Feflow 軟件，以桃園煤礦為研究對(duì)象，建立滲流場和溶質(zhì)運(yùn)移化學(xué)場模型，探討“雙場”耦合作用下溶質(zhì)擴(kuò)散影響因素并預(yù)測(cè)注漿結(jié)束后溶質(zhì)濃度演變趨勢(shì)，為煤礦開展區(qū)域注漿影響下礦井水害的水化學(xué)預(yù)警奠定基礎(chǔ)，為煤礦有效開展水害防治提供科學(xué)依據(jù)，具有重要的理論和實(shí)踐意義。

1 研究區(qū)概況

淮北煤田桃園煤礦位于安徽省北部，區(qū)內(nèi)地勢(shì)平坦，無基巖出露，均為厚層松散層所覆蓋。其構(gòu)造位置位于宿南向斜西翼、宿南背斜東翼，基底為奧陶系和寒武系地層，構(gòu)造地質(zhì)及水文地質(zhì)條件復(fù)雜（圖1a）。礦井被F2斷層切割分成南、北兩塊，并且以F2斷層為界，地層走向發(fā)生變化，F(xiàn)2斷層以北為北北西，以南為北北東。礦井總體為一走向近南北、傾向東的單斜構(gòu)造。其中太灰含水層在礦井南部廣大地區(qū)和北部接受補(bǔ)給，向中北部F2斷層附近Ⅱ4采區(qū)徑流排泄（圖1b）。

圖1 桃園煤礦概況Fig.1 Overview of Taoyuan Coal Mine

礦井主要充水水源為煤層頂?shù)装迳皫r裂隙水、底板太原組灰?guī)r水、奧陶系灰?guī)r水。其中，太原組灰?guī)r含水層發(fā)育11 層薄層灰?guī)r，三灰、四灰較厚，巖溶發(fā)育，富水性較強(qiáng)，是太灰水防控的重點(diǎn)。為消除水害威脅，桃園煤礦對(duì)物探、鉆探和水文地質(zhì)試驗(yàn)查明有垂向?qū)畼?gòu)造的區(qū)域，采取地面定向順層鉆孔對(duì)太原組三灰進(jìn)行區(qū)域注漿改造（圖1c），變含水層為隔水層，同時(shí)阻斷三灰與其下伏灰?guī)r尤其是奧灰之間的水力聯(lián)系[20]。

針對(duì)10 號(hào)煤層水害防治，桃園煤礦從2014 年底開始，分別針對(duì)Ⅱ1、Ⅱ2及Ⅱ4采區(qū)的多個(gè)工作面進(jìn)行注漿改造，至2020 年5 月底板區(qū)域注漿工程全部結(jié)束，治理區(qū)主要分布在礦井中部和南部，并相互銜接（圖1c）。區(qū)域注漿面積總計(jì)4 587 459 m2，注漿量共計(jì)798 052 t。

2 注漿擾動(dòng)下“雙場”演化分析

在全面分析桃園井田地質(zhì)與水文地質(zhì)條件的基礎(chǔ)上，對(duì)客觀地質(zhì)體科學(xué)概化，首先利用Feflow 軟件建立滲流場模型，此基礎(chǔ)上選定較為穩(wěn)定的Cl?為模擬因子，建立溶質(zhì)運(yùn)移化學(xué)場模型，經(jīng)識(shí)別驗(yàn)證后得到注漿擾動(dòng)下含水層相關(guān)參數(shù)信息，并對(duì)模型的可靠性進(jìn)行評(píng)價(jià)，為“雙場”耦合作用研究奠定基礎(chǔ)。

2.1 地質(zhì)模型概化

由于部分太灰水水位觀測(cè)孔分布在井田西側(cè)，超出井田范圍，為更好利用觀測(cè)孔水位進(jìn)行識(shí)別驗(yàn)證，在桃園礦井范圍的基礎(chǔ)上將模型范圍向西延伸約150 m，至虛線范圍處。

從空間上看，研究區(qū)太灰水流整體上以水平運(yùn)動(dòng)為主、垂向運(yùn)動(dòng)為輔，地下水運(yùn)動(dòng)符合達(dá)西定律；考慮含水層之間的流量交換，灰?guī)r水運(yùn)動(dòng)可以概化為空間三維流；灰?guī)r水系統(tǒng)的垂向運(yùn)動(dòng)是含水層間構(gòu)造裂隙導(dǎo)水造成水頭差異引起的；地下水系統(tǒng)的輸入輸出隨時(shí)間、空間變化，故地下水為非穩(wěn)定流；參數(shù)隨空間變化，體現(xiàn)了系統(tǒng)的非均質(zhì)性，所以含水介質(zhì)概化為非均質(zhì)各向異性介質(zhì)[21]。

桃園煤礦與西部錢營孜煤礦、南部祁南煤礦處于同一個(gè)水文地質(zhì)單元上，相互間有一定的水力聯(lián)系，則桃園井田邊界概化為流量邊界，側(cè)向流量隨時(shí)間和位置的不同有所變化。灰?guī)r含水層垂向邊界包括1 灰頂板和4 灰底板，概化為隔水邊界。

綜上，概化的地質(zhì)模型如圖2 所示。

圖2 研究區(qū)地質(zhì)模型Fig.2 Geological model of study area

2.2 滲流場演化模擬

太灰水滲流為三維非穩(wěn)定流，可用Feflow 軟件建立注漿擾動(dòng)下太灰水滲流場模型。該承壓含水層運(yùn)動(dòng)方程如下：

式中：?為灰?guī)r水滲流區(qū)域；h為含水層的水位，m；Kxx、Kyy、Kzz為含水層沿x、y、z方向的滲透系數(shù)，m/d；Kn為邊界面法向滲透系數(shù)，m/d；S為承壓含水層貯水率，m?1；ω為含水層的源匯項(xiàng)代數(shù)和，d?1；H0為含水層的初始水位分布，m；τ2為滲流區(qū)域的側(cè)向邊界；τ4為滲流區(qū)域的下邊界，即承壓含水層底部的隔水邊界；為邊界面法線方向；q（x，y，z，t）為定義為側(cè)向邊界的單寬流量，m2/（d·m），流入為正，流出為負(fù)，隔水邊界為0[22]。

模型網(wǎng)格剖分采用三角剖分，剖分?jǐn)?shù)量10 000個(gè)，且對(duì)F2斷層加密剖分。經(jīng)網(wǎng)格剖分后，構(gòu)建三維地質(zhì)結(jié)構(gòu)體，頂板標(biāo)高?604.3～?1 252.2 m，底板埋深?646.6～?1 271.4 m，礦井北部地層相對(duì)較厚（圖3，因概念模型向西延伸，圖中頂板埋深淺至?260.7 m）。

圖3 模擬區(qū)網(wǎng)格剖分及三維地質(zhì)結(jié)構(gòu)體Fig.3 Mesh generation and three-dimensional geological structure of simulation area

統(tǒng)籌考慮注漿擾動(dòng)下長觀孔水位數(shù)據(jù)、井下放水孔水壓、礦井涌水量、水文地質(zhì)參數(shù)獲取時(shí)間、后期水樣測(cè)試等信息，太灰水滲流模型模擬時(shí)間選定為2021 年1 月至2021 年12 月（注漿結(jié)束后）。經(jīng)多次調(diào)參，獲得驗(yàn)證后滲流場模型。調(diào)整后參數(shù)見表1，調(diào)參后模型如圖4 所示。

表1 注漿期間驗(yàn)證后水文地質(zhì)參數(shù)Table 1 Hydrogeological parameters after verification during grouting

圖4 注漿后太灰滲流模型Fig.4 Model of percolation after grouting

注漿后初始水位（圖4a）顯示，礦井南部水位較高，多處于?110～?140 m 之間；向北至礦井中部Ⅱ4、Ⅱ2采區(qū)水位變化較大，總體向Ⅱ4采區(qū)西北部徑流排泄；礦井北部，水位較低，但變化不大。模型運(yùn)行365 d 后，2021 年12 月實(shí)際水位特征（圖4b）與2021 年1 月的初始水位相近，期間太灰長觀孔（98?觀3、2019?探3、2017?觀1、2014?觀1）水位值波動(dòng)不大。

圖5 顯示，各太灰水位長觀孔模擬水位與實(shí)測(cè)水位變化趨勢(shì)基本一致，相對(duì)來講中北部的98?3 太灰長觀孔水位變化與實(shí)際水位差距最大（推測(cè)此處水文地質(zhì)條件復(fù)雜，模型刻畫不夠精細(xì)），2014?觀1和2019?探3 模擬效果最好。通過對(duì)檢測(cè)孔各月水位實(shí)測(cè)與模擬值統(tǒng)計(jì)結(jié)果顯示，差值范圍在|0?1|、|1?3|、>|3| m 的點(diǎn)分別占到38.6%、33.3%以及28.1%，總體來看，誤差小于3 m 的點(diǎn)占到了71.9%，表明模擬取得了良好的效果。

圖5 注漿后太灰長觀孔中實(shí)際水位與模擬水位歷時(shí)演化Fig.5 Temporal evolution of actual water level and simulated water level in long view hole after grouting

2.3 化學(xué)場演化模擬

基于線性彌散定律和溶質(zhì)質(zhì)量守恒定律，建立水動(dòng)力彌散方程，采用對(duì)流彌散方程計(jì)算溶質(zhì)在三維地下水流中的遷移，控制方程如下：

式中：c為選定溶質(zhì)的質(zhì)量濃度，mg/L；θ為灰?guī)r含水層孔隙率；t為時(shí)間；qi為達(dá)西流速，m/d；qs為源匯項(xiàng)或者流出項(xiàng)的體積通量值，d?1；cs為源匯項(xiàng)或者流出項(xiàng)的質(zhì)量濃度，mg/L；Rk為表示各種反應(yīng)的特定術(shù)語；Dij為水動(dòng)力彌散系數(shù)張量，m2/d；c0為初始濃度，mg/L；τ1為定濃度邊界條件；τ2為彌散通量邊界條件；τ3為定溶質(zhì)通量邊界[22]。

基于桃園煤礦區(qū)域注漿實(shí)況及水巖作用分析可知，地下水常規(guī)離子中Cl?和是地下水環(huán)境中比較穩(wěn)定的常規(guī)離子。原含水層中含量較低，可近似為0，漿液析出水中穩(wěn)定，如若注漿后地下水檢測(cè)到，則說明受到漿液析出水的混合影響，但是由于質(zhì)量濃度總體較低，實(shí)測(cè)濃度及模型結(jié)果等誤差可能較大，所以此處選擇地下水中含量較高且比較穩(wěn)定的Cl?，作為本次模擬使用的溶質(zhì)[23]。

基于已建立滲流場模型，參數(shù)設(shè)置主要涉及灰?guī)r含水層孔隙率、溶質(zhì)初始濃度、邊界條件和水動(dòng)力彌散系數(shù)，其中含水層孔隙率和邊界條件根據(jù)實(shí)際注漿工況，分別設(shè)為時(shí)間的函數(shù)。當(dāng)滲流場確定后，水動(dòng)力彌散系數(shù)取決于3 個(gè)介質(zhì)參數(shù)，分別是有效擴(kuò)散系數(shù)Dd、縱向彌散度αL、橫向彌散度αT[24]。溶質(zhì)運(yùn)移模擬期設(shè)為2015 年1 月至2021 年11 月（2021 年最后一次水樣采集時(shí)間是2021 年11 月，為便于模擬值與實(shí)測(cè)值的驗(yàn)證，此處模擬期截止至2021 年11 月），基于定期定點(diǎn)（圖2）取樣測(cè)試的Cl?含量，識(shí)別驗(yàn)證后主要參數(shù)見表2，設(shè)溶質(zhì)運(yùn)移模型分別運(yùn)行1 a（365 d，至2015 年12 月，此時(shí)僅Ⅱ2采區(qū)有注漿）、3 a（1 085 d，至2017 年12 月，此時(shí)Ⅱ2、Ⅱ4采區(qū)注漿已結(jié)束）、7 a（2 525 d，至2021 年11 月，所有注漿均已結(jié)束），運(yùn)行結(jié)果如圖6 所示。

表2 溶質(zhì)運(yùn)移模型相關(guān)參數(shù)Table 2 Relevant parameters of solute transport model

圖6 溶質(zhì)運(yùn)移模擬結(jié)果Fig.6 Simulation results of solute transport

運(yùn)行1 a 后（圖6a），Ⅱ2采區(qū)南部徑流相對(duì)強(qiáng)烈地區(qū)溶質(zhì)運(yùn)移速度較快，濃度變化較大；運(yùn)行3 a 后（圖6b），Ⅱ2、Ⅱ4注漿區(qū)周邊溶質(zhì)擴(kuò)散均較明顯，且Ⅱ4注漿區(qū)西側(cè)及東側(cè)均為徑流排泄區(qū)（東側(cè)是Ⅱ6大巷抽放水所致），溶質(zhì)擴(kuò)散更為明顯。運(yùn)行7 a 后（圖6c），Ⅱ2、Ⅱ4注漿區(qū)溶質(zhì)擴(kuò)散均較明顯，且西側(cè)溶質(zhì)運(yùn)移更為明顯，溶質(zhì)質(zhì)量濃度變化較大，主要是由于兩采區(qū)西部滲透系數(shù)較大，地下水流動(dòng)性較好。雖然Ⅱ4注漿區(qū)注漿時(shí)間較晚，但其溶質(zhì)擴(kuò)散更為明顯，首先由于兩采區(qū)東部滲透性差異，Ⅱ4采區(qū)東部滲透系數(shù)為Ⅱ2采區(qū)東部的15 倍（Ⅱ4采區(qū)東部滲透系數(shù)0.03，Ⅱ2采區(qū)東部為0.002）；然后Ⅱ4采區(qū)東部Ⅱ6大巷放水量較大，其它采區(qū)放水量相對(duì)較少。Ⅱ1采區(qū)相對(duì)注漿時(shí)間較晚，且滲透系數(shù)較小，則溶質(zhì)運(yùn)移較慢。

繼續(xù)運(yùn)行至50 a（圖6d），結(jié)果顯示，Ⅱ4采區(qū)Cl?已均勻運(yùn)移擴(kuò)散開來，Ⅱ2采區(qū)Cl?也已基本擴(kuò)散，與周邊色階差異不大，而Ⅱ1注漿區(qū)Cl?與周邊色差較明顯，主要由于Ⅱ4、Ⅱ2采區(qū)滲透系數(shù)較大，而Ⅱ1采區(qū)巖石致密，滲透性較差，漿液析出水基本以“滯水”狀態(tài)存在。

模型運(yùn)行期間，各監(jiān)測(cè)點(diǎn)水位如圖7 所示，其中太灰長觀孔監(jiān)測(cè)點(diǎn)水位，已通過前述滲流場模型驗(yàn)證；R 開頭圖例均表示水樣采集點(diǎn)處水位，無水位監(jiān)測(cè)數(shù)據(jù)。由圖7 可以看出，Ⅱ4、Ⅱ2采區(qū)注漿時(shí)，注漿工作面內(nèi)監(jiān)測(cè)點(diǎn)（如R27 檢41 孔、R28 Ⅱ1029老空水、R8-10Ⅱ 1042 軌道、R12 北二灰?guī)r、R29Ⅱ 1027 底抽巷）灰?guī)r水位上升明顯。注漿區(qū)周邊監(jiān)測(cè)點(diǎn)（如Ⅱ4采區(qū)東部R13-21Ⅱ 6 大巷檢50 孔、Ⅱ2采區(qū)東部R11Ⅱ1026泄水巷、Ⅱ1采區(qū)西部南大巷放水孔等）以及距離較遠(yuǎn)的太灰水位長觀孔水位（如995-觀1、98-觀3 等）受注漿擾動(dòng)不夠明顯，其中太灰98-觀3 水位波動(dòng)主要是由于開采期間，不定量抽水引起的水位正常波動(dòng)。注漿結(jié)束后即進(jìn)入開采期，隨巷道內(nèi)灰?guī)r放水量增加，太灰含水層下降致安全水位。

圖7 模型運(yùn)行期間部分監(jiān)測(cè)點(diǎn)太灰水位動(dòng)態(tài)變化Fig.7 Dynamic change of limestone water level in Taiyuan Formation at some monitoring points during model operation

通過2021 年每隔1 個(gè)月采集的太灰水樣中Cl?質(zhì)量濃度驗(yàn)證可知（圖8），5 次采集的Ⅱ4采區(qū)Ⅱ6大巷檢50 孔太灰水Cl?實(shí)測(cè)質(zhì)量濃度為254.7～271.5 mg/L，平均267.5 mg/L，而模擬所得Cl?質(zhì)量濃度為256.3～279.9 mg/L，平均269.6 mg/L，模擬期內(nèi)隨時(shí)間增加，Cl?質(zhì)量濃度呈降低趨勢(shì)；5 次采集的Ⅱ 4 采區(qū)Ⅱ 6 大巷補(bǔ)檢6 孔太灰水中Cl?實(shí)測(cè)質(zhì)量濃度為220.1～270.8 mg/L，平均256.3 mg/L，而模擬所得Cl?質(zhì)量濃度為243.6～286.2 mg/L，平均260.4 mg/L，模擬期內(nèi)隨時(shí)間增加，Cl?質(zhì)量濃度呈降低趨勢(shì)；Ⅱ 1采區(qū)南三煤眼太灰水中Cl?實(shí)測(cè)質(zhì)量濃度為209.8～231.4 mg/L，平均218.5 mg/L，該模擬期內(nèi)Cl?模擬質(zhì)量濃度為203.7～204.3 mg/L，平均204.1 mg/L，模擬期內(nèi)隨時(shí)間增加，質(zhì)量濃度呈降低趨勢(shì)，但降低較慢。綜上，Ⅱ 4 采區(qū)溶質(zhì)運(yùn)移模擬效果較好，Ⅱ 1 采區(qū)模擬濃度偏低約14.4 mg/L，誤差約為實(shí)際值的6.6%，仍低于10%，總體上模擬效果較好。其中，中北部的Ⅱ 4采區(qū)Cl?質(zhì)量濃度多高于南部的Ⅱ1采區(qū)，主要由于桃園煤礦南部太灰水水位較高，多處于補(bǔ)給徑流區(qū)，中北部水位較低，處于排泄區(qū)，而沿地下水補(bǔ)給、徑流、排泄路徑，地下水中溶解固體物質(zhì)總量呈增大趨勢(shì)。

圖8 溶質(zhì)運(yùn)移模擬值與實(shí)測(cè)值驗(yàn)證分析Fig.8 Verification analysis of simulated and measured values of solute transport

自2021 年每隔一個(gè)月取桃園煤礦太灰水樣，因井下采樣條件有限，每次取樣地點(diǎn)稍有不同。其中自2022 年以來涉及Ⅱ4采區(qū)Ⅱ6大巷補(bǔ)檢6 孔、檢50 孔及Ⅱ1采區(qū)南三煤眼太灰水中Cl?實(shí)測(cè)質(zhì)量濃度見表3，可知水質(zhì)監(jiān)測(cè)點(diǎn)實(shí)測(cè)值與模擬值基本相近，進(jìn)一步證明了構(gòu)建模型的正確性。

表3 2022 年以來相關(guān)水質(zhì)監(jiān)測(cè)點(diǎn)Cl–實(shí)測(cè)質(zhì)量濃度Table 3 Measured concentration of Cl– at relevant water quality monitoring sites since 2022

模型運(yùn)行期間，各監(jiān)測(cè)點(diǎn)Cl?質(zhì)量濃度動(dòng)態(tài)變化如圖9 所示，其中太灰長觀孔Cl?質(zhì)量濃度無實(shí)測(cè)數(shù)據(jù)；R 開頭的監(jiān)測(cè)點(diǎn)表示水樣采集處Cl?質(zhì)量濃度，部分時(shí)間節(jié)點(diǎn)處濃度已通過驗(yàn)證?？梢钥闯鲎{時(shí)，Ⅱ4、Ⅱ2采區(qū)注漿工作面內(nèi)監(jiān)測(cè)點(diǎn)Cl?質(zhì)量濃度瞬時(shí)增加（如R27、R8-10、R29、R12），變化明顯，隨時(shí)間增加Cl?質(zhì)量濃度逐漸降低。隨與注漿區(qū)距離的增大，監(jiān)測(cè)點(diǎn)Cl?質(zhì)量濃度瞬時(shí)動(dòng)態(tài)響應(yīng)減弱（如R7、R11、R30 等）。距離注漿區(qū)較遠(yuǎn)處，Cl?質(zhì)量濃度瞬時(shí)響應(yīng)不明顯（如R1-6、R8-10、98-觀3、2014-觀1等），且短期內(nèi)隨時(shí)間增加，質(zhì)量濃度變化不明顯。

圖9 模型運(yùn)行期間部分監(jiān)測(cè)點(diǎn)處Cl?質(zhì)量濃度變化Fig.9 Changes of Cl? concentration at some monitoring points during model operation

3 注漿擾動(dòng)下“雙場”耦合作用分析

注漿擾動(dòng)下地下水環(huán)境的演化是由滲流場與化學(xué)場耦合作用引起的。為此，在上述滲流場與化學(xué)場演化分析基礎(chǔ)上，進(jìn)一步研究地下水滲流路徑上溶質(zhì)演化及其滲透系數(shù)、孔隙率等水力參數(shù)的影響機(jī)制，進(jìn)而構(gòu)建區(qū)域注漿擾動(dòng)下基于“雙場”耦合作用的溶質(zhì)擴(kuò)散預(yù)測(cè)模型。

3.1 滲流路徑上溶質(zhì)演化

以注漿擾動(dòng)下Ⅱ 4 采區(qū)為例（圖10），模型運(yùn)行初期（圖10a），滲流較強(qiáng)區(qū)域位于Ⅱ4采區(qū)西部，向西北方向徑流排泄；隨著Ⅱ4采區(qū)Ⅱ6大巷排水增加，Ⅱ4采區(qū)太灰水徑流方向逐漸由向西變?yōu)橄驏|（圖10b），注漿期間（第635～639 天）太灰水由注漿區(qū)向周邊徑流排泄（圖10c），特別向東西兩側(cè)及西北部排泄較多（西部滲透系數(shù)大于東部）；注漿后短期內(nèi)，太灰水仍是由注漿區(qū)向東西兩側(cè)特別是西側(cè)徑流排泄，由于注漿區(qū)滲透系數(shù)減小，隨時(shí)間增加，徑流強(qiáng)度減弱（圖10d）；隨時(shí)間延長，Ⅱ 6 大巷排水增加，Ⅱ 4 采區(qū)西部太灰水徑流方向逐漸由注漿區(qū)向西部排泄，變?yōu)橛晌鞑肯驏|部排泄，但強(qiáng)度較弱（圖10e），只要Ⅱ6大巷持續(xù)排水，Ⅱ4采區(qū)西部接受補(bǔ)給向東部徑流排泄的局部流場特征便會(huì)保持下去（圖10f 為Ⅱ 4 采區(qū)注漿后約5 a 時(shí)間）。

圖10 模型運(yùn)行過程中Ⅱ 4 采區(qū)局部“雙場”演化（箭頭代表徑流的方向和強(qiáng)度）Fig.10 Local “dual field”evolution in Ⅱ 4 mining area during model operation (Arrows represent the direction and intensity of runoff)

圖10 中D1-D5 及X1-X6 分別為注漿區(qū)東部和西部徑流路徑，基于Ⅱ 4 采區(qū)局部流場演化分析，分別選擇注漿后1～5 a 兩路徑處Cl?質(zhì)量濃度進(jìn)行分析，具體如圖11 所示?？梢钥闯?，Ⅱ 4 采區(qū)東部D1-D5 路徑上，補(bǔ)給區(qū)D1-D3 Cl?質(zhì)量濃度較高，注漿結(jié)束約2 a 后，D1-D3 處Cl?質(zhì)量濃度開始明顯下降，且D1 下降最明顯，隨時(shí)間運(yùn)移，下降幅度減??；D4-D5 徑流排泄處Cl?質(zhì)量濃度逐漸增加，且位于Ⅱ6大巷徑流排泄處的D5 監(jiān)測(cè)點(diǎn)增加較快，隨時(shí)間延長，D4 增加幅度逐漸加大，即徑流排泄區(qū)Cl?質(zhì)量濃度均逐漸增大，且排泄區(qū)Cl?質(zhì)量濃度最大。

圖11 兩徑流路徑上監(jiān)測(cè)點(diǎn)Cl?質(zhì)量濃度Fig.11 Cl? concentration at monitoring points on two runoff paths

Ⅱ 4 采區(qū)西部X1～X6 路徑上，補(bǔ)給區(qū)X1、X2 Cl?質(zhì)量濃度較高，而X2 處Cl?質(zhì)量濃度下降較快，主要是由于其西部未注漿區(qū)滲透系數(shù)相對(duì)較大，滲流相對(duì)較強(qiáng)，隨時(shí)間延長，Cl?下降變緩，并逐漸穩(wěn)定；隨徑流路線延長，X3～X6 監(jiān)測(cè)點(diǎn)Cl?質(zhì)量濃度依次增加，但隨時(shí)間延長，其變化趨勢(shì)特點(diǎn)不同，主要是由于期間徑流排泄路線有所變化而致。

由此可知，隨滲流場補(bǔ)徑排路徑發(fā)生變化，溶質(zhì)濃度也會(huì)相應(yīng)發(fā)生變化，一般情況下沿補(bǔ)徑排路徑，溶質(zhì)濃度呈增大趨勢(shì)。

3.2 溶質(zhì)運(yùn)移與水力參數(shù)的關(guān)系

基于該溶質(zhì)運(yùn)移模型，將某一參數(shù)按一定倍數(shù)調(diào)整，其它參數(shù)不變，運(yùn)行模型，觀察監(jiān)測(cè)點(diǎn)處溶質(zhì)質(zhì)量濃度，進(jìn)而研究相關(guān)參數(shù)對(duì)溶質(zhì)運(yùn)移的影響機(jī)理。

根據(jù)前述分析，溶質(zhì)運(yùn)移影響因素主要包括地層滲透系數(shù)、孔隙率、彌散度、水力梯度等。下面以礦井中部的Ⅱ 4 采區(qū)Ⅱ 6 大巷檢50 孔（圖2 中R13-21）和礦井南部的Ⅱ 1 采區(qū)南三煤眼（圖2 中R1-6）為例，予以闡述。

1）滲透系數(shù)。圖12a 顯示，當(dāng)Ⅱ4采區(qū)東部和Ⅱ1采區(qū)西部滲透系數(shù)在0.003～0.075 m/d 變化時(shí)，R13-21 點(diǎn)Cl?質(zhì)量濃度為167.3～288.0 mg/L，R1-6點(diǎn)Cl?質(zhì)量濃度為197.4～206.7 mg/L。即隨滲透系數(shù)的增加，監(jiān)測(cè)點(diǎn)溶質(zhì)濃度均增加，但是不同監(jiān)測(cè)點(diǎn)增加幅度不一樣，其中R13-21 點(diǎn)濃度增加較為明顯。圖12b 為模型運(yùn)行至2021 年底多個(gè)監(jiān)測(cè)點(diǎn)多個(gè)時(shí)間節(jié)點(diǎn)處Cl?質(zhì)量濃度統(tǒng)計(jì)結(jié)果（下同），可以看出，滲透系數(shù)較小時(shí)（小于0.01 m/d），Cl?質(zhì)量濃度變化幅度較大，說明此時(shí)Cl?運(yùn)移受滲透系數(shù)影響較大；而當(dāng)滲透系數(shù)大于0.1 m/d 時(shí)，Cl?質(zhì)量濃度受滲透系數(shù)影響較小，但是隨滲透系數(shù)增加，Cl?質(zhì)量濃度變化幅度增大。即當(dāng)滲透系數(shù)大于0.1 m/d 時(shí)，隨滲透系數(shù)增加，Cl?運(yùn)移受滲透系數(shù)影響逐漸增大。

圖12 滲透系數(shù)對(duì)溶質(zhì)運(yùn)移的影響Fig.12 Influence of permeability coefficient on solute transport

2）孔隙率。圖13a 顯示，當(dāng)Ⅱ4采區(qū)東部和Ⅱ1采區(qū)西部滲透系數(shù)在0.03～0.75 m/d 變化時(shí)，R13-21 點(diǎn)Cl?質(zhì)量濃度為368.0～192.8 mg/L，R1-6 點(diǎn)Cl?質(zhì)量濃度為207.5～158.3 mg/L。即隨孔隙率的增加，監(jiān)測(cè)點(diǎn)溶質(zhì)濃度均減小，但是減小幅度不一樣，R13-21 點(diǎn)Cl?質(zhì)量濃度減小較快，R1-6 點(diǎn)Cl?質(zhì)量濃度減小較慢。即在其他條件相同時(shí)，隨孔隙率增加，注入的漿液儲(chǔ)存于流經(jīng)巖石的量增多，結(jié)果使向前運(yùn)移擴(kuò)散的量減少，造成距離注漿區(qū)稍遠(yuǎn)的監(jiān)測(cè)點(diǎn)Cl?濃度減小。圖13b 顯示，孔隙率較小時(shí)（小于0.15），Cl?質(zhì)量濃度變化范圍較大，此時(shí)Cl?運(yùn)移受孔隙率影響較大?？傮w上隨孔隙率增大，Cl?質(zhì)量濃度變化范圍越來越小，即對(duì)Cl?運(yùn)移影響越小。

3）彌散度。圖14a 顯示，當(dāng)Ⅱ4采區(qū)東部和Ⅱ1采區(qū)西部彌散度在5～125 m 之間變化時(shí)，R13-21點(diǎn)Cl?質(zhì)量濃度為189.1～294.1 mg/L，R1-6 點(diǎn)Cl?質(zhì)量濃度為195.3～208.7 mg/L。即隨彌散度的增加，監(jiān)測(cè)點(diǎn)溶質(zhì)濃度均增加，不同監(jiān)測(cè)點(diǎn)增加幅度不一樣，R13-21 點(diǎn)濃度增加較為明顯，R1-6 點(diǎn)濃度增加不夠明顯。圖14b 顯示，隨礦井范圍內(nèi)彌散度同步增大，其對(duì)Cl?質(zhì)量濃度的影響逐漸增大。

圖14 彌散度對(duì)溶質(zhì)運(yùn)移的影響Fig.14 Effect of dispersion on solute transport

4）達(dá)西流速和水力梯度。模型運(yùn)行期間地下水滲透流速可以監(jiān)測(cè)得到，而某監(jiān)測(cè)點(diǎn)滲透系數(shù)是一定的，根據(jù)達(dá)西定律公式v=ki，可獲得水力梯度[25]。以R13-21（Ⅱ4采區(qū)東部檢50 孔）為例，期間滲透流速為0.004 5～0.005 6 m/d（圖15a），滲透流速與水力梯度同步變化，經(jīng)換算水力梯度為0.145～0.190（圖15b）。圖15a、15b、15c 顯示，注漿后一定時(shí)間內(nèi)隨水力梯度和滲透流速的增加，Cl?質(zhì)量濃度增加（281.07～411.48 mg/L）。另外，隨時(shí)間運(yùn)移，不同監(jiān)測(cè)點(diǎn)Cl?質(zhì)量濃度變化不同，但總的變化趨勢(shì)基本一致，先增大后減小。不同地點(diǎn)增大幅度不同，濃度峰值出現(xiàn)的時(shí)間也不同，一般在注漿后18～22 a 的時(shí)間內(nèi)Cl?質(zhì)量濃度達(dá)到峰值，隨后Cl?質(zhì)量濃度開始降低，即該區(qū)注漿的影響作用開始減弱，約在40 a 后達(dá)到平衡狀態(tài)。就R13-21 監(jiān)測(cè)點(diǎn)而言，從開始注漿至2033 年間（約注漿后18 a 內(nèi)），區(qū)域注漿對(duì)桃園井田化學(xué)場的影響呈增加趨勢(shì)，2033 年后（區(qū)域注漿結(jié)束18 a后），其影響逐漸減弱。圖15c 中R23-26 監(jiān)測(cè)點(diǎn)Cl?質(zhì)量濃度先降低再增大，主要是由于前期抽放水影響大于注漿擾動(dòng)影響所致；注漿后約10 a 開始，注漿擾動(dòng)影響大于抽放水。圖15d 也可以看出，由于抽放水影響，R23-26 監(jiān)測(cè)點(diǎn)水力梯度變化較大。

圖15 達(dá)西流速和水力梯度對(duì)溶質(zhì)運(yùn)移的影響Fig.15 Effects of Darcy flow rate and hydraulic gradient on solute transport

漿液相對(duì)密度對(duì)溶質(zhì)運(yùn)移也有一定的影響。比重不同，漿液析出水Cl?質(zhì)量濃度也不同。室內(nèi)漿液析水試驗(yàn)可知，相對(duì)密度1.2、1.3、1.4、1.5、1.6 的漿液，其析出水Cl?質(zhì)量濃度分別為395.2、439.3、496.5、501.6、513.1 mg/L，模型中注入濃度分別輸入不同比重漿液Cl?質(zhì)量濃度，運(yùn)行至2021 年12 月（圖16，以Ⅱ4采區(qū)東部檢50 孔監(jiān)測(cè)點(diǎn)為例）?？梢钥闯?，隨漿液相對(duì)密度增加，監(jiān)測(cè)點(diǎn)Cl?質(zhì)量濃度增加。2021 年在Ⅱ4采區(qū)東部檢50 孔所取的6 個(gè)太灰水測(cè)試結(jié)果顯示，Cl?質(zhì)量濃度為224.8～290.7 mg/L，平均為253.4 mg/L，在相應(yīng)模擬應(yīng)力期（2021 年）的Cl?質(zhì)量濃度范圍內(nèi)（247.0～266.1 mg/L）。

圖16 漿液相對(duì)密度對(duì)溶質(zhì)運(yùn)移的影響Fig.16 Effect of slurry specific gravity on solute transport

3.3 “雙場”耦合作用下溶質(zhì)擴(kuò)散預(yù)測(cè)模型

通過上述分析可知，某監(jiān)測(cè)點(diǎn)溶質(zhì)質(zhì)量濃度除與滲透系數(shù)、孔隙率、彌散度、水力梯度、注入溶質(zhì)濃度等參數(shù)有關(guān)外，還與注漿時(shí)間、與注漿點(diǎn)距離等參數(shù)有關(guān)?；隍?yàn)證后模型多次設(shè)參運(yùn)行所獲取的約690 組數(shù)據(jù)，借助SPSS 軟件擬合分析，建立了溶質(zhì)擴(kuò)散預(yù)測(cè)模型，以期預(yù)測(cè)分析區(qū)域注漿擾動(dòng)情況：

式中：C(Cl?)為預(yù)測(cè)某監(jiān)測(cè)點(diǎn)Cl?質(zhì)量濃度，mg/L；Cz為注漿區(qū)所注入漿液析出水中Cl?質(zhì)量濃度，mg/L；k為監(jiān)測(cè)點(diǎn)處地層滲透系數(shù)，m/d；I為監(jiān)測(cè)點(diǎn)處水力梯度，無量綱；αT為監(jiān)測(cè)點(diǎn)處橫向彌散度，m；t為注漿結(jié)束后時(shí)間，a；n為監(jiān)測(cè)點(diǎn)處孔隙率，無量綱；d為監(jiān)測(cè)點(diǎn)距注漿區(qū)的距離，m；C0為監(jiān)測(cè)點(diǎn)處Cl?初始質(zhì)量濃度，mg/L。

該模型將滲流場模型涉及參數(shù)（滲透系數(shù)）、溶質(zhì)運(yùn)移化學(xué)場模型涉及參數(shù)（孔隙率、彌散度、水力梯度等）與注漿工況（注入濃度、注漿時(shí)間等）集合在一起，可以預(yù)測(cè)得出注漿擾動(dòng)下某監(jiān)測(cè)點(diǎn)溶質(zhì)濃度。該模型變異系數(shù)分析（表4）顯示，決定系數(shù)R2為0.807，接近于1，可知其擬合效果較好。

表4 變異系數(shù)分析Table 4 Analysis of variance

4 結(jié) 論

1）驗(yàn)證后的注漿擾動(dòng)下太灰含水層滲流場模型和溶質(zhì)運(yùn)移模型，其模擬值與實(shí)測(cè)值相差均較小，模擬效果均較好。

2）溶質(zhì)運(yùn)移主要受滲透系數(shù)、彌散度、水力梯度、滲透流速、注漿時(shí)間、漿液比重等參數(shù)控制，并發(fā)現(xiàn)在注漿后18～22 a 內(nèi)Cl?質(zhì)量濃度達(dá)到峰值，隨后Cl?質(zhì)量濃度降低，約在40 a 后達(dá)到平衡狀態(tài)。

3）建立的“雙場”耦合作用下溶質(zhì)擴(kuò)散預(yù)測(cè)模型，誤差率小于10%的數(shù)據(jù)達(dá)81.4%，說明所建立的預(yù)測(cè)模型基本可靠。

研究具有較好的理論與實(shí)踐意義，可為區(qū)域注漿擾動(dòng)下受注目標(biāo)含水層水環(huán)境演化研究及煤礦水害預(yù)測(cè)預(yù)警提供科學(xué)依據(jù)。