單井抽出-回滲同步循環(huán)地下水水力控制模型研究

王新港，楊 昱，徐祥健，韓 旭，夏 甫，鄧 圣，肖 瀚*，姜永海

1. 中國(guó)環(huán)境科學(xué)研究院，國(guó)家環(huán)境保護(hù)地下水污染模擬與控制重點(diǎn)實(shí)驗(yàn)室，北京 100012

2. 河北工程大學(xué)能源與環(huán)境工程學(xué)院，河北 邯鄲 056038

單井抽出-回滲同步循環(huán)地下水水力控制技術(shù)是從地下水污染羽范圍內(nèi)采用單井抽出污染地下水，經(jīng)過(guò)處理后，在抽水井周邊同步等量回滲至包氣帶和含水層，通過(guò)抽出-回滲同步循環(huán)，人為改變污染羽及其周邊區(qū)域地下水流場(chǎng)和水力梯度，將受污染的地下水控制在一定范圍內(nèi)，防止污染物進(jìn)一步擴(kuò)散蔓延，進(jìn)而將受污染的地下水體與未受污染的清潔地下水體分隔開(kāi)，是地下水污染修復(fù)的重要手段[1-5]. 該技術(shù)的優(yōu)勢(shì)是操作簡(jiǎn)單、運(yùn)行成本較低、加速地下水循環(huán)流動(dòng)縮短地下水修復(fù)時(shí)間，影響該技術(shù)實(shí)施效果的重要因素是水動(dòng)力條件以及含水層水文地質(zhì)條件，水動(dòng)力條件主要由抽注水量控制，而含水層水文地質(zhì)條件主要指含水層滲透性[6-13].

在地下水污染水力控制技術(shù)領(lǐng)域內(nèi)，已有研究主要關(guān)注抽注技術(shù)參數(shù)的選取以及抽注水井布局和抽注水量的優(yōu)化對(duì)水力控制效果的影響，使抽注井的布局和抽注水量更為合理，從而縮短抽注時(shí)間、提高水力控制效率[14-19]. 如Bear等[20]定量地研究了注水井補(bǔ)給對(duì)含水層地下水運(yùn)動(dòng)規(guī)律的影響，為水力控制技術(shù)的發(fā)展奠定了理論基礎(chǔ)；Shar等[21]、Christ等[22]采用解析、半解析法分別研究了相同流量下任意雙抽水井及不同流量下多抽水井水力控制區(qū)域的邊界形態(tài)特征，為數(shù)值模擬提供了理論基礎(chǔ)；隨后，Wang等[23]用數(shù)值模擬和優(yōu)化耦合方法，采用拉普拉斯變換有限單元法對(duì)均勻等厚污染含水層的多個(gè)抽注井進(jìn)行位置布設(shè)，提高了捕獲污染的有效性；Matott等[24]采用分析單元法對(duì)井群系統(tǒng)的抽水井?dāng)?shù)量、位置和抽水量進(jìn)行了優(yōu)化；劉苑[25]和萬(wàn)鵬[26]均采用MGO模塊結(jié)合遺傳算法對(duì)水力控制系統(tǒng)進(jìn)行優(yōu)化，通過(guò)約束抽水井?dāng)?shù)量、抽水流量、水頭和抽水總量等參數(shù)，提高了水力控制能力. 目前地下水污染水力控制技術(shù)研究中數(shù)值模擬成為主要解決手段，常用的數(shù)值模擬軟件主要包括Groundwater Modeling Systems(GMS)、Visual MODFLOW、FEFLOW及Groundwater Vistas等[27-28]. 然而單井抽出-回滲技術(shù)特點(diǎn)是抽出-回滲循環(huán)同步進(jìn)行，通過(guò)人工干擾/控制實(shí)現(xiàn)水力控制下抽出-回滲同步進(jìn)行的研究相對(duì)較少[29].

研究采取物理模型和數(shù)值模型相結(jié)合的方法，搭建實(shí)驗(yàn)室尺度砂柱孔隙含水層物理模型，基于該物理模型構(gòu)建MODFLOW地下水流數(shù)值模型并校準(zhǔn)，使模擬地下水流可以較好地刻畫(huà)實(shí)際地下水流；利用校準(zhǔn)后的數(shù)值模型模擬不同抽注水量和含水層滲透性條件下實(shí)現(xiàn)水力控制時(shí)不同位置地下水水位的差異，基于潛水完整井穩(wěn)定流公式-裘布依公式的數(shù)據(jù)擬合，分析實(shí)現(xiàn)水力控制時(shí)含水層不同位置水位與抽注水量和含水層滲透系數(shù)的關(guān)系以及回滲半徑對(duì)水力控制的影響. 該研究論證了單井抽出-回滲同步循環(huán)地下水水力控制技術(shù)實(shí)現(xiàn)人工控制地下水流場(chǎng)和水力梯度的可行性，并量化了水動(dòng)力條件和水文地質(zhì)條件的影響，以期為該技術(shù)在實(shí)際場(chǎng)地的應(yīng)用提供重要參考.

1 材料與方法

1.1 物理模型

1.1.1裝置搭建

單井抽出-回滲同步循環(huán)地下水水力控制技術(shù)試驗(yàn)裝置如圖1所示. 試驗(yàn)裝置由砂柱、布設(shè)于砂柱中軸線的抽水管、砂柱表層的布水管以及砂柱含水層內(nèi)部不同位置、不同深度的測(cè)壓管組成.

圖1 單井抽出-回滲同步循環(huán)地下水水力控制技術(shù)試驗(yàn)裝置Fig.1 Schematic diagram of single well pumpingrecharge synchronous cyclical groundwater hydraulic control technique device

砂柱由1 cm厚有機(jī)玻璃制成，直徑為100 cm，高度為100 cm. 砂柱內(nèi)均勻填充石英砂〔采購(gòu)于天津某石料廠的80～120目(100目=0.150 mm)石英砂〕，通過(guò)下部進(jìn)水且每5 cm進(jìn)行飽和夯實(shí)的方式裝填，構(gòu)建厚度為70 cm的均質(zhì)潛水含水層. 抽水管(長(zhǎng)100 cm、直徑2 cm、開(kāi)篩70 cm)垂直固定于砂柱中軸線位置并包裹160目過(guò)濾篩網(wǎng)防止石英砂進(jìn)入，管底位于砂柱底部(用于模擬抽水完整井)，管頂與蠕動(dòng)泵相連接. 布水管以抽水管為圓心，半徑32 cm范圍內(nèi)螺旋狀均勻鋪設(shè)于砂層表面(布水管內(nèi)徑0.6 cm，每隔5 cm開(kāi)孔)，一端與蠕動(dòng)泵相連接，使抽出的水等量同步均勻回滲至砂柱含水層中. 蠕動(dòng)泵型號(hào)為保定蘭格BT100-1L，泵頭為YII25，可調(diào)整流量范圍為0.1～180 mL/min.

1.1.2測(cè)試方法

裝置監(jiān)測(cè)點(diǎn)位如圖2所示，采用上密下疏、內(nèi)密外疏的方式，通過(guò)砂柱側(cè)面48個(gè)監(jiān)測(cè)孔在砂柱內(nèi)不同位置、不同深度安裝48個(gè)測(cè)壓管，用于監(jiān)測(cè)含水層水位. 垂向共8層，測(cè)壓深度分別為1、3、5、7、14、28、48、70 cm，每層設(shè)置6個(gè)測(cè)壓點(diǎn)，分別與中心抽水井井軸相距2、4、8、14、22、32 cm. 測(cè)壓管讀數(shù)實(shí)時(shí)反映砂柱內(nèi)不同位置、不同深度含水層水位的變化. 砂柱物理模型的滲透系數(shù)和孔隙度分別由抽水試驗(yàn)、補(bǔ)水試驗(yàn)測(cè)算[30]. 抽水試驗(yàn)測(cè)定滲透系數(shù)的原理是基于潛水完整井穩(wěn)定流公式-裘布依公式[31]：

圖2 單井抽出-回滲同步循環(huán)地下水水力控制技術(shù)試驗(yàn)裝置概念Fig.2 Conceptual diagram of single well pumpingrecharge synchronous cyclical groundwater hydraulic control technique device

式中：Q為抽水井流量，m3/d；K為含水層滲透系數(shù)，m/d；h為含水層中任意位置地下水水位，m；hw為抽水井水位，m；r為含水層中某處與抽水井井軸距離，m；rw為抽水井半徑，m.

實(shí)驗(yàn)室尺度砂柱物理模型搭建完成后，給定抽水量為2.5 cm3/s，運(yùn)行開(kāi)展抽水試驗(yàn)，獲取抽水穩(wěn)定后不同位置、不同深度觀測(cè)水位，計(jì)算得到砂柱模型填土的滲透系數(shù)為0.01 cm/s，孔隙度為0.4.

1.2 數(shù)值模型

1.2.1模擬軟件選取

GMS是一個(gè)界面友好、處理功能強(qiáng)大的地下水?dāng)?shù)值模擬軟件[28,32]. MODFLOW模塊離散方法簡(jiǎn)單、求解方式多樣、采用了模塊化結(jié)構(gòu)，廣泛應(yīng)用于地下水滲流和溶質(zhì)運(yùn)移[33]. 該研究采用GMS中的MODFLOW模塊構(gòu)建實(shí)驗(yàn)室尺度砂柱孔隙含水層三維地下水穩(wěn)定流模型. MODFLOW模塊的水流控制方程[31]如下：

式中：W為水流源匯項(xiàng)，m3/d；t為模擬預(yù)測(cè)的時(shí)間，d；h為地下水水位值，m；μs為含水層儲(chǔ)水率，m-1；Kx、Ky、Kz為不同軸向上的滲透系數(shù)值，m/d.

1.2.2模型構(gòu)建

數(shù)值模型參照上述物理模型等比例構(gòu)建，模擬單井抽出-回滲循環(huán)過(guò)程含水層流場(chǎng)，4種模型構(gòu)建參數(shù)條件如表1所示. 模擬區(qū)域?yàn)樵囼?yàn)砂柱的石英砂填裝區(qū)(直徑100 cm、高度70 cm)，水平方向剖分為1 cm ×1 cm的正方形網(wǎng)格(100行、100列)，垂向剖分為4層，從上至下各層厚度分別為20、16、24、10 cm.

表1 模型構(gòu)建參數(shù)Table 1 Parameters for model construction

模擬區(qū)域下部及側(cè)向邊界設(shè)定為零流量邊界；模擬區(qū)域中心有抽水井排泄，設(shè)定為定流量邊界(抽水量分別設(shè)定為1、2.5、5、10 cm3/s)；模擬區(qū)域上部邊界有回滲補(bǔ)給(補(bǔ)給范圍為以抽水井井軸為圓心，半徑32 cm的圓)，設(shè)定為定流量邊界(回滲量分別設(shè)定為1、2.5、5、10 cm3/s，對(duì)應(yīng)的補(bǔ)給速率分別為3.11×10-4、7.78×10-4、1.56×10-3、3.11×10-3cm/s).

1.3 數(shù)值模型校準(zhǔn)

為保證模型模擬結(jié)果可靠，對(duì)模型進(jìn)行校準(zhǔn)驗(yàn)證. 物理模型的48個(gè)觀測(cè)點(diǎn)分別在抽出-回滲量為1、2.5、5、10 cm3/s達(dá)到抽出-回滲平衡時(shí)記錄實(shí)測(cè)水位，對(duì)比抽出-回滲量分別為1、2.5、5、10 cm3/s時(shí)數(shù)值模型的模擬水位與物理模型的實(shí)測(cè)水位，微調(diào)滲透系數(shù)，重復(fù)多次試算，直至大部分觀測(cè)點(diǎn)模擬水位與實(shí)測(cè)水位的誤差達(dá)到模擬精度要求.

數(shù)值模擬結(jié)果的準(zhǔn)確度通過(guò)納什效率系數(shù)(Nash-Sutcliffe efficiency coefficient, NSE)驗(yàn)證[34].

式中：Qo為實(shí)測(cè)水位值，cm；Qm為模擬水位值，cm；為t時(shí)的實(shí)測(cè)水位值，為實(shí)測(cè)水位值的總平均值，cm.

NSE取值為負(fù)無(wú)窮至1，NSE接近1，表示模型質(zhì)量好，模型可信度高；NSE接近0，表示模擬結(jié)果接近觀測(cè)值的平均值水平，即總體結(jié)果可信，但誤差大；NSE遠(yuǎn)遠(yuǎn)小于0，則模型不可信.

數(shù)值模型校準(zhǔn)后，滲透系數(shù)調(diào)整為0.005 cm/s，模擬水位和實(shí)測(cè)水位(見(jiàn)表2)分散比較結(jié)果如圖3所示，納什效率系數(shù)達(dá)到0.88，說(shuō)明模型精度較為理想，可以較好地刻畫(huà)物理模型中的實(shí)際地下水流場(chǎng).

表2 不同抽出-回滲量條件下模擬水位和實(shí)測(cè)水位對(duì)比Table 2 Comparison between simulated and observed groundwater levels under different pumping-return seepage conditions

圖3 不同抽出-回滲量條件下模擬水位和實(shí)測(cè)水位的分散比較Fig.3 Comparison diagram between simulated and observed groundwater levels under different pumping-return seepage conditions

2 結(jié)果與分析

2.1 實(shí)現(xiàn)水力控制時(shí)的模擬流場(chǎng)

當(dāng)抽出-回滲量分別為1、2.5、5、10 cm3/s時(shí)，達(dá)到抽出-回滲平衡后，模擬流場(chǎng)三維空間分布如圖4所示；距離抽水井軸0、2、4、8、14、22、32 cm處的模擬水位以及抽水影響半徑見(jiàn)表3. 抽出-回滲量對(duì)水位降落漏斗的深度影響較大，抽出-回滲量越大，降落漏斗的深度越大. 當(dāng)抽出-回滲量分別為1、2.5、5、10 cm3/s時(shí)，最大降深分別達(dá)到1.76、4.55、9.75、18.65 cm，分別為含水層厚度的2.5%、6.5%、13.9%、26.6%. 然而，抽出-回滲量對(duì)水位降落漏斗的范圍幾乎無(wú)影響，抽水影響半徑幾乎不隨抽出-回滲量的變化而變化.

圖4 不同抽出-回滲量條件下模擬流場(chǎng)三維空間分布特征Fig.4 Characteristics of simulated three-dimensional spatial distribution of the flow field under different pumping-return seepage conditions

表3 不同抽出-回滲量條件下模擬水位及抽水影響半徑Table 3 Simulated water level and pumping influence radius under different pumping-return seepage conditions

2.2 實(shí)現(xiàn)水力控制時(shí)水位與抽注水量和含水層滲透性的關(guān)系

已有研究[6]表明，實(shí)現(xiàn)水力控制時(shí)含水層不同位置處的水位h主要取決于抽出-回滲量Q(水動(dòng)力條件)和含水層滲透系數(shù)K(水文地質(zhì)條件). 基于抽出-回滲量和含水層滲透性設(shè)置不同情景，利用校準(zhǔn)后的數(shù)值模型模擬不同情景實(shí)現(xiàn)水力控制時(shí)地下水水位的分布特征，采用潛水完整井穩(wěn)定流公式-裘布依公式開(kāi)展數(shù)據(jù)分析擬合，分析實(shí)現(xiàn)水力控制時(shí)含水層不同位置水位與抽注水量和含水層滲透系數(shù)的關(guān)系，獲取單井抽出-回滲同步循環(huán)地下水水力控制技術(shù)關(guān)鍵參數(shù)之間的定量關(guān)系.

選取《工程地質(zhì)手冊(cè)》提供的不同巖土類型(細(xì)砂、中砂、粗砂、礫砂)滲透系數(shù)的經(jīng)驗(yàn)值[35]設(shè)定4種水文地質(zhì)條件情景，選取不同抽出-回滲量(1、2.5、5、10 cm3/s)設(shè)定4種水動(dòng)力條件情景，共計(jì)16種組合情景(見(jiàn)表4).

表4 不同情景水位模擬結(jié)果Table 4 Simulation results of water level under different scenarios

基于裘布依公式，對(duì)16種組合情景模擬得到的抽水井處水位以及距抽水井井軸4、8、14、22、32 cm處的模擬水位進(jìn)行數(shù)據(jù)分析擬合，獲取Q/K與(h2-hw2)/(lgr-lgrw)之間的定量關(guān)系(見(jiàn)圖5)和擬合方程〔見(jiàn)式(5)〕. 由圖5可見(jiàn)，Q/K與(h2-hw2)/(lgr-lgrw)之間呈線性關(guān)系，擬合系數(shù)(R2)達(dá)到0.99. 擬合后得到的線性方程斜率(α)為154.84～197.85，遠(yuǎn)大于式(2)中的斜率(1.366)，主要是因?yàn)轸貌家拦降倪m用條件是補(bǔ)給項(xiàng)僅為側(cè)向滲流補(bǔ)給，不含上層回滲補(bǔ)給，而單井抽出-回滲同步循環(huán)地下水水力控制技術(shù)在抽水的同時(shí)從上層同步等量回滲補(bǔ)水，導(dǎo)致砂柱邊界地下水水位明顯高于中部地下水水位，水力梯度明顯增大.

圖5 16種情景下模擬水位的數(shù)據(jù)分析擬合結(jié)果Fig.5 Simulated groundwater levels analyses and fitting results under 16 various scenarios

擬合方程(5)的斜率區(qū)間為154.84～197.85，存在一定程度波動(dòng)，為探究斜率波動(dòng)原因，基于表4設(shè)置的不同情景和現(xiàn)有模擬結(jié)果分析，利用校準(zhǔn)后的數(shù)值模型模擬不同回滲半徑下實(shí)現(xiàn)水力控制時(shí)擬合方程(5)的斜率變化. 目前回滲半徑為32 cm，進(jìn)一步考慮回滲半徑為15、20、25、30、35、40、45、50 cm 8種情況，在此基礎(chǔ)上對(duì)16種組合情景模擬得到的抽水井處水位以及距井軸4、8、14、22、32 cm處的模擬水位進(jìn)行數(shù)據(jù)分析擬合，對(duì)式(5)中的斜率進(jìn)行校正，結(jié)果如表5所示. 由表5可見(jiàn)，回滲半徑與斜率之間存在一定關(guān)系，且擬合度比較高，R2均達(dá)到0.99. 隨著回滲半徑增加，擬合后得到的線性方程斜率在一定范圍內(nèi)逐漸增大，且擬合方程呈明顯的線性正相關(guān)關(guān)系；當(dāng)回滲半徑接近邊界半徑時(shí)斜率范圍逐漸減小，其中回滲半徑等于邊界半徑時(shí)，斜率范圍為178.18～183.9，而回滲半徑為15 cm時(shí)，斜率范圍為164.99～206.34，范圍較大；此外，當(dāng)回滲半徑固定時(shí)，隨著監(jiān)測(cè)點(diǎn)位距離的增加，方程斜率逐漸增加. 主要原因是對(duì)于不變的水文地質(zhì)條件，當(dāng)流量確定后，回滲半徑大小不同，直接影響回滲范圍內(nèi)地下水水位的變化，因此合理地調(diào)整回滲范圍大小，有利于實(shí)現(xiàn)水力控制.

表5 16種情景下的方程斜率Table 5 Equation of the slope of equations under 16 various scenarios

3 討論

通過(guò)單井-抽出回滲同步循環(huán)地下水水力控制技術(shù)的數(shù)值模擬結(jié)果，使用裘布依公式進(jìn)行了擬合，對(duì)該問(wèn)題中水動(dòng)力參數(shù)、水文地質(zhì)參數(shù)以及回滲范圍進(jìn)行了初步分析，得到了參數(shù)之間的對(duì)應(yīng)關(guān)系，當(dāng)回滲量取值為零時(shí)，方程中各參數(shù)之間變?yōu)榻?jīng)典裘布依公式中的參數(shù)關(guān)系. 對(duì)不同巖土類型均質(zhì)含水層的模擬旨在概化不同類型滲透性場(chǎng)地中的參數(shù)關(guān)系.

在實(shí)際應(yīng)用上，對(duì)場(chǎng)地回滲過(guò)程要求為均勻穩(wěn)定入滲，然而自然場(chǎng)地一般為非均質(zhì)性，入滲一般是不均勻、不穩(wěn)定的，對(duì)該問(wèn)題的較好解決方式仍是使用數(shù)值模擬，通過(guò)不同尺度下大量水動(dòng)力條件和非均質(zhì)的水文地質(zhì)條件的組合，同時(shí)調(diào)整抽出-回滲量和回滲范圍進(jìn)行模擬研究，以此得出參數(shù)之間的進(jìn)一步關(guān)系，使之更貼近于實(shí)際應(yīng)用條件.

4 結(jié)論

a) 校準(zhǔn)后的數(shù)值模型模擬結(jié)果與物理模型實(shí)際觀測(cè)結(jié)果相吻合，表明數(shù)值模型可以較好地刻畫(huà)實(shí)現(xiàn)水力控制時(shí)物理模型中的實(shí)際地下水流場(chǎng)，可利用數(shù)值模型進(jìn)行單井抽出-回滲同步循環(huán)地下水水力控制技術(shù)研究.

b) 抽出-回滲量的變化對(duì)水位降落漏斗的范圍幾乎無(wú)影響，而對(duì)水位降落漏斗的深度影響較大，抽出-回滲量越大，降落漏斗越深.

c) 基于潛水完整井穩(wěn)定流公式-裘布依公式，對(duì)不同抽出-回滲量(水動(dòng)力條件)和含水層滲透性(水文地質(zhì)條件)情景模擬水位的擬合，獲取了實(shí)現(xiàn)水力控制時(shí)含水層不同位置水位與抽出-回滲量和含水層滲透系數(shù)之間的定量關(guān)系，可以為該技術(shù)在實(shí)際場(chǎng)地的應(yīng)用提供重要參考.

d) 實(shí)現(xiàn)水力控制時(shí)含水層不同位置水位與抽出-回滲量和含水層滲透系數(shù)之間的定量關(guān)系受到回滲半徑的影響，對(duì)于特定水文地質(zhì)條件的含水層，合理調(diào)整抽出-回滲量以及回滲范圍大小，有利于實(shí)現(xiàn)水力控制.