， ，，，

(1.南京大學(xué)地球科學(xué)與工程學(xué)院，江蘇 南京 210023；2.南京水利科學(xué)研究院， 江蘇 南京 210029)

隨著工業(yè)化的發(fā)展加速，地下水污染問題日益突出，極大地破壞著生態(tài)環(huán)境并威脅著民眾的健康[1-3]。根據(jù)污染物種類和含水層類型的不同，有多種不同的地下水污染治理與修復(fù)方法[4]，其中基于水力控制的污染物抽出-處理(pump and treat, PAT)技術(shù)最為典型。其主要優(yōu)點(diǎn)在于簡(jiǎn)便易行，可對(duì)地下水污染事件做出快速反應(yīng)，尤其適用于期望短期內(nèi)迅速降低污染水平的地下水水源地含水層系統(tǒng)，缺點(diǎn)在于修復(fù)工程復(fù)雜、且運(yùn)行成本高昂。因此，有效降低PAT系統(tǒng)的運(yùn)行成本對(duì)于該技術(shù)的推廣應(yīng)用具有重要意義。自20世紀(jì)80年代以來，地下水模擬-優(yōu)化模型就應(yīng)用于確定地下水PAT系統(tǒng)的最優(yōu)策略[5]。然而，與國內(nèi)外廣泛存在的模塊化地下水模擬軟件相比，地下水模擬-優(yōu)化模型方面的研究雖然很多[6-9]，但模塊化的模擬-優(yōu)化軟件相對(duì)很少，這種狀況極大地限制了地下水模擬-優(yōu)化模型在實(shí)際場(chǎng)地中地下水污染治理的優(yōu)化實(shí)踐。

MGO(Modular Groundwater Optimizer)[10]是一款模塊化的地下水模擬優(yōu)化軟件，它通過耦合地下水水流模型MODFLOW[11]和污染物運(yùn)移模型MT3DMS[12]結(jié)合通用管理目標(biāo)建立地下水PAT系統(tǒng)的模擬-優(yōu)化模型，并可根據(jù)需要選擇不同優(yōu)化算法對(duì)模型進(jìn)行求解。本文基于MGO軟件，運(yùn)用遺傳算法對(duì)美國猶他州的某地下水污染場(chǎng)地進(jìn)行抽出-處理系統(tǒng)的優(yōu)化研究，確定抽水井群的單井最優(yōu)動(dòng)態(tài)抽水量，實(shí)現(xiàn)以最低的系統(tǒng)運(yùn)行成本達(dá)到修復(fù)含水層的最優(yōu)目標(biāo)。

1 MGO軟件簡(jiǎn)介

MGO(Modular Groundwater Optimizer)是美國阿拉巴馬大學(xué)教授Zheng和Wang基于地下水水流模型MODFLOW和污染物運(yùn)移模型MT3DMS 開發(fā)出來的一款針對(duì)地下水PAT系統(tǒng)優(yōu)化的通用模塊化優(yōu)化管理軟件。利用MGO不僅能建立完全耦合地下水模擬模型(包括水流模型和運(yùn)移模型)的嵌入式模擬-優(yōu)化模型，而且對(duì)于特別復(fù)雜的地下水優(yōu)化問題還能結(jié)合實(shí)際情況建立新的基于替代模型(或響應(yīng)函數(shù))的地下水模擬-優(yōu)化模型。無論是嵌入式優(yōu)化模型還是基于響應(yīng)函數(shù)的優(yōu)化模型，MGO都可選擇遺傳算法(GA)、模擬退火(SA)和禁忌搜索(TS)等多種進(jìn)化算法中的任一算法來求解，從而保證能找到復(fù)雜地下水非線性模擬-優(yōu)化問題的最優(yōu)解。

MGO既能應(yīng)用于地下水系統(tǒng)的單一水量/水質(zhì)優(yōu)化管理，也能應(yīng)用于水量與水質(zhì)的聯(lián)合優(yōu)化管理，尤其適用于基于水力控制的地下水污染治理的PAT系統(tǒng)中的抽注水方案優(yōu)化管理。同時(shí)，MGO采用模塊化設(shè)計(jì)(見圖1)，能直接利用地下水模擬模型MODFLOW和MT3DMS中的輸入文件，大大簡(jiǎn)化了軟件的使用，便于推廣應(yīng)用。

圖1 MGO代碼設(shè)計(jì)流程圖

2 基于MGO的地下水PAT系統(tǒng)模擬-優(yōu)化模型

為了能同時(shí)考慮各種不同的水量/水質(zhì)優(yōu)化管理問題，基于MGO設(shè)計(jì)的地下水PAT系統(tǒng)優(yōu)化模型可表示為以下通用數(shù)學(xué)形式[10]：

目標(biāo)函數(shù)：

(1)

約束條件：

(2)

Qmin≤Qi≤Qmax

(3)

hmin≤hm≤hmax

(4)

(5)

cmin≤cm≤cmax

(6)

(7)

上述優(yōu)化模型中，式(1)為通用目標(biāo)函數(shù)，J是管理目標(biāo)，可以是總成本、總抽水量或含水層中污染物的總?cè)コ康?。用戶可以根?jù)實(shí)際情況選取目標(biāo)函數(shù)右端項(xiàng)中的一項(xiàng)或若干項(xiàng)(亦可表示為最大化)。Qi是某個(gè)管理期內(nèi)的第i個(gè)井在單位時(shí)間內(nèi)的抽/注水量，正值代表注水、負(fù)值代表抽水。對(duì)于單個(gè)管理期優(yōu)化問題，任何井的抽/注水量均為是常數(shù)，只需用一個(gè)參數(shù)表示；而對(duì)于多管理期的優(yōu)化問題，由于任何一個(gè)井在不同管理期的抽/注水量可以不同，則需要多個(gè)參數(shù)來表示井在不同管理期的抽/注水量。Mi第i井的溶質(zhì)移除量。F(q, h, c)用戶自定義成本函數(shù)，可以依賴于流量q,水頭h,濃度c。N是待優(yōu)化參數(shù)的總數(shù)目，下同。yi是個(gè)二進(jìn)制數(shù)，當(dāng)考慮參數(shù)i活躍時(shí)，yi=1；否則yi=0。di是第i參數(shù)相關(guān)的井徑深度。Δti是參數(shù)i的抽/注水持續(xù)時(shí)長(zhǎng)期。α1每個(gè)井的固定資產(chǎn)成本。α2是井筒單位深度的安裝、打鉆成本。α3是單位體積的抽水和處理成本。α4是單位質(zhì)量污染物的處理成本。α5是外部用戶指定的成本系數(shù)。

式(2)-(7)為通用約束條件。式(2)要求任何時(shí)間段內(nèi)的決策井?dāng)?shù)目不得超過一個(gè)給定的約束值(NW)。 式(3)是對(duì)任一管理期內(nèi)單位時(shí)間的抽水量約束，Qmin、Qmax分別是對(duì)應(yīng)井的最小、最大單位抽水量。式(4)是水頭約束，hm是監(jiān)測(cè)點(diǎn)位置的水頭，hmin、hmax分別是其最小、最大約束水頭。式(5)是水力約束，Δhmin是外部監(jiān)測(cè)點(diǎn)與內(nèi)部監(jiān)測(cè)點(diǎn)之間所允許的最小水頭差。式(6)是濃度約束，cm是任一監(jiān)測(cè)點(diǎn)的濃度，cmin、cmax是該監(jiān)測(cè)點(diǎn)的最小、最大濃度允許值。式(7)可作為均衡約束，Qm是與從井I1到井I2的抽/注水量之和成比例的抽水量，其中A、B為比例系數(shù)。同樣地，用戶可以根據(jù)優(yōu)化模型中目標(biāo)函數(shù)的組成來確定和選擇必要的約束條件。

優(yōu)化模型中涉及到兩種變量：決策變量和狀態(tài)變量。決策變量包括抽/注水井的數(shù)量、位置和流量；狀態(tài)變量就是指地下水的水位和污染物的濃度。水位和濃度狀態(tài)變量對(duì)應(yīng)于決策變量的響應(yīng)分別由水流模型MODFLOW和溶質(zhì)運(yùn)移模型MT3DMS來獲得更新。因此，對(duì)于嵌入式模擬-優(yōu)化模型來說，采用進(jìn)化算法求解時(shí)，都要針對(duì)決策變量(井流量、井位等)不斷重復(fù)調(diào)用水流模型MODFLOW與溶質(zhì)運(yùn)移模型MT3DMS來計(jì)算相應(yīng)的狀態(tài)變量(水位、濃度)，由此進(jìn)一步計(jì)算目標(biāo)函數(shù)值并判斷約束條件是否滿足，從而搜索得到模型的最優(yōu)策略。

圖2 研究區(qū)場(chǎng)地位置及模型第1層的初始TCE污染羽分布

3 實(shí)例應(yīng)用

3.1 污染場(chǎng)地概述與模擬-優(yōu)化模型

3.1.1 場(chǎng)地概述及模擬模型

研究場(chǎng)址位于美國猶他州大鹽湖以南約10 km的圖埃爾(Tooele)山谷。含水層主要由沖積層構(gòu)成，但在其中部有一向上抬升的基巖塊，使地下水從沖積沉積物中流向斷裂風(fēng)化的基巖，然后又回到?jīng)_積沉積物中。區(qū)域地下水流大致向東南往西北方向流動(dòng)。場(chǎng)地主要污染來自于其東南角的一個(gè)工業(yè)區(qū)，主要污染物為TCE。監(jiān)測(cè)結(jié)果表明，污染物主要分布在含水層的淺部，但由于場(chǎng)地地下水流的復(fù)雜性，淺部污染羽和深部污染羽的范圍并不一致，甚至東北角的污染羽已延伸到特征邊界之外(見圖2)。

現(xiàn)已確定采用基于水力技術(shù)的抽出-處理技術(shù)來治理和修復(fù)該含水層，該P(yáng)AT系統(tǒng)由美國某陸軍工程兵團(tuán)設(shè)計(jì)，包含16口抽水井、13口注水井和1口附加抽水井(E-fixed)。相應(yīng)的模擬模型已實(shí)現(xiàn)了抽水/不抽水條件下的多次校正，該模型是一個(gè)三維、穩(wěn)態(tài)的MODFLOW模型，有4層，165行，和99列。平面上剖分為正方形單元網(wǎng)格，其邊長(zhǎng)為60.96 m，垂向上自上而下4個(gè)模型層的平均厚度分別為45.72 m、30.48 m、45.72和91.44 m?；贛T3DMS構(gòu)建與之對(duì)應(yīng)的污染物運(yùn)移模型。模型時(shí)間離散化為有7個(gè)3年的應(yīng)力期來匹配優(yōu)化方案中定義的管理期。

現(xiàn)以校正后的模擬模型計(jì)算出的2003年1月的水頭分布與污染羽分布分別作為初始水頭、初始濃度分布。由于現(xiàn)存PAT系統(tǒng)模擬結(jié)果顯示，沒有可行解能滿足污染羽截?cái)噙吔绲募s束條件。經(jīng)過前期的專家論證，需新增4口注水井(IN1～I(xiàn)N4)，原有PAT系統(tǒng)中的16口抽水井、13口注水井和1口附加抽水井(E-fixed)，僅需保留抽水井E11和E15、注水井I4及E-fixed，但E11和E15均以最大抽水能力運(yùn)行，各井位置見圖2。需要說明的是，為了維持PAT系統(tǒng)的穩(wěn)定性，所有抽水井與注水井的水量總量相等。

3.1.2 優(yōu)化模型及其求解

本文在確定了井的位置基礎(chǔ)上，令每一個(gè)管理期的所有注/抽水井(E-fixed與I4除外)的流量為一個(gè)決策變量，以制定出一種最優(yōu)的動(dòng)態(tài)解。選擇使用GA進(jìn)行最優(yōu)解的搜索，相關(guān)的參數(shù)設(shè)置如下：POPSIZE=200，PCROSS=0.5， PMUTATE=0.05，NPOSIBL=32。

根據(jù)目標(biāo)函數(shù)的一般形式，由于井位、井?dāng)?shù)已定，總成本項(xiàng)主要來自抽水用電成本項(xiàng)和化學(xué)處理項(xiàng)，其中本文把抽水用電成本項(xiàng)當(dāng)作安裝成本，調(diào)整后的通用目標(biāo)函數(shù)可簡(jiǎn)化為如下PAT系統(tǒng)的最小運(yùn)行成本：

(8)

相應(yīng)的約束條件如下：

|Qsumin|≤|Q*|≤|Qsumax|

(9)

(10)

QI4=-(QEi+QNIi)-QE-fixed

(11)

Cm≤Cmax

(12)

0≤QEi≤│Li│

(13)

0≤QNIi≤Li

(14)

式(8)中，α1取值為34.5K$(1K$表示1 000美元)，α3取值為1.01e-05K$/(m3)；M表示管理期，M=7；N=6；yji是二進(jìn)制數(shù)，當(dāng)?shù)趈管理期新建或優(yōu)化第i個(gè)井時(shí)，yji=1，否則yji=0。式(9)、(10)都作為全局約束條件：其中式(9)是對(duì)總抽水量的約束；Q*是參與模擬的總抽水量(包括E-fixed的抽水量)，Qsumax為最大總抽水量，取值為-41 427.42 m3/d，Qsumin為最小總抽水量，取值為0；式(10)是對(duì)抽/注水量和的約束，Qi為參與模擬的第i個(gè)井的抽/注水量，i=1,2,…,8；Qc為最大均衡誤差值，取值為5.45 m3/d。NW為各管理期參與模擬的井?dāng)?shù)，本次優(yōu)化各管理期NW均取值為8。式(11)作為均衡約束條件：QI4是現(xiàn)用注水井I4的注水量；QEi是第i個(gè)現(xiàn)用抽水井的抽水量，i=1,2；QNIi是第i個(gè)新建注水井的注水量，i=1,2,3,4；QE-fixed是E-fixed的抽水量，值為-7 767.68 m3/d。式(12)為濃度約束條件：Cm是第m個(gè)濃度約束邊界單元上TCE的濃度值，Cmax為不同管理期末的約束濃度上限，取值為5 μg/L。式(13)、(14)分別為現(xiàn)用井的抽/注水量、新建注水井的注水量的約束條件：QEi為第i個(gè)現(xiàn)存井的抽/注水量，i=1, 2, 3，Li為第i個(gè)井的最大抽/注水能力。PAT系統(tǒng)中參與優(yōu)化的井的最大抽/注水能力如表1所示。

圖3 執(zhí)行最優(yōu)策略最終管理期末得到模型第1層中的污染羽分布

圖4 執(zhí)行最優(yōu)策略最終管理期末得到模型第2層中的污染羽分布表1 優(yōu)化前PAT系統(tǒng)中待優(yōu)化抽/注水井的抽/注水能力及其模擬層

井Li/m3/d模擬層E11-3 365.992E15-3 314.212, 3I44 163.472IN13 107.071IN23 107.072IN33 107.071IN43 107.071

如表1所示，抽水井E15抽采2、3兩層，因此只要以其第2層的抽水量作為決策變量，而將第3層的抽水量處理為依賴變量，抽水量的分配依據(jù)兩層的導(dǎo)水系數(shù)可由優(yōu)化程序自動(dòng)計(jì)算得到。

3.2 優(yōu)化結(jié)果與分析

最優(yōu)抽水策略滿足所有約束條件，包括在4個(gè)模擬層的污染約束邊界單元處，污染羽濃度在每個(gè)管理期末均接近或小于5 μg/L。第1、2含水層是污染物主要分布層，這兩層在第7個(gè)管理期的污染羽分布見圖3和圖4。動(dòng)態(tài)策略下各管理期的抽/注水量見圖5、表2。原穩(wěn)態(tài)策略下的抽/注水量同見圖5。

圖5 最優(yōu)策略下各管理期的抽水(負(fù)值)/注水(正值)量分布表2 動(dòng)態(tài)策略下各管理期內(nèi)優(yōu)化得到PAT系統(tǒng)的抽水(負(fù)值)和注水(正值)方案

井號(hào)優(yōu)化前優(yōu)化后管理期抽/注水量/m3/d1234567E11-3 365.99-3 365.99-3 245.85-3 365.99-3 365.99-3 365.99-3 257.14-3 149.81E15-3 314.210-2 876.33-3 314.21-3 314.21-3 314.21-3 314.21-2 886.52NI13 107.071 202.871 109.663 107.071 904.3072 004.552 606.012 706.20NI23 107.072 906.852 407.983 107.073 107.073 107.073 107.072 906.85NI33 107.073 005.743 107.073 107.073 107.073 107.073 107.072 405.42NI43 107.073 005.743 107.073 107.073 005.743 005.742 606.013 107.07I42 019.601 010.784 158.302 019.543 322.113 223.232 912.522 677.75

優(yōu)化后的動(dòng)態(tài)策略抽水井E15在第1管理期內(nèi)抽水量為零，并且在第2、第7管理期內(nèi)低于最大負(fù)荷抽水率運(yùn)行，在整個(gè)項(xiàng)目期內(nèi)抽水井E15抽水量相比優(yōu)化前減少了18.02%，抽水井E11在第2、6、7管理期也是低于最大抽水量運(yùn)行，項(xiàng)目期內(nèi)抽水井E11的抽水量減少了1.89%，管理期的總抽水量減小了9.89%，此外，抽注水量之和小于給定最大值5.45 m3/d，說明來自三口抽水井的抽水量都從5口注水井注入含水層中。因此動(dòng)態(tài)策略下整個(gè)項(xiàng)目期的總抽水量減小，使得抽水用電成本和化學(xué)處理成本相應(yīng)減少。

在穩(wěn)定策略下，總成本為1 258.85K$，其 中抽水用電成本、化學(xué)處理成本各占73.79%、26.21%。而優(yōu)化后動(dòng)態(tài)策略下抽水用電成本相比優(yōu)化前減少10.62%；化學(xué)處理成本相比優(yōu)化前減少13.08%?？偝杀緸? 117.02K$，相比優(yōu)化前總成本減少了11.27%。

4 結(jié)語

本文通過MGO軟件對(duì)美國猶他州某處TCE污染場(chǎng)地的PAT修復(fù)系統(tǒng)進(jìn)行模擬優(yōu)化管理，實(shí)現(xiàn)了對(duì)修復(fù)系統(tǒng)的動(dòng)態(tài)優(yōu)化管理。優(yōu)化結(jié)果表明優(yōu)化后的修復(fù)策略減少了11.27%的修復(fù)成本，同時(shí)保證了修復(fù)效率。實(shí)例研究的結(jié)果表明MGO可以有效地為地下水修復(fù)提高經(jīng)濟(jì)、高效的管理措施，為今后的地下水污染修復(fù)提供一種可靠的工具，具有較好的應(yīng)用前景。

[1]呂書君. 我國地下水污染分析[J]. 地下水.2009.31(1): 1-5.

[2]Soderstrom M, Boldemann C, Sahlin U, et al. The quality of the outdoor environment influences childrens health-across-sectional study of preschools. Acta Paediatrica, 2013.102: 83-91.

[3]Ozdemir S, Johnson FR. Estimating willingness to pay: do health and environmental researchers have different methodological standards Applied Economics, 2013.45: 2215-2229.

[4]李瑋,王明玉,韓占濤,等. 棕地地下水污染修復(fù)技術(shù)篩選方法研究—以某廢棄化工廠污染場(chǎng)地為例[J]. 水文地質(zhì)工程地質(zhì).2016.43(3): 131-141.

[5]Minsker BS editor. Long-Term Groundwater Monitoring: The State of the Art[M]. ASCE/EWRI, 0-7844-0678-2, Reston, VA; 2003.

[6]肖傳寧,盧文喜,安永凱,等. 基于兩種耦合方法的模擬-優(yōu)化模型在地下水污染源識(shí)別中的對(duì)比[J]. 中國環(huán)境科學(xué).2015.35(8): 2393-2399.

[7]吳劍鋒,彭偉, 錢家忠, 等. 基于INPGA的地下水污染治理多目標(biāo)優(yōu)化管理模型:I—理論方法與算例驗(yàn)證[J]. 地質(zhì)評(píng)論.2011.57(2): 277-283.

[8]吳劍鋒,彭偉, 錢家忠, 等. 基于INPGA的地下水污染治理多目標(biāo)優(yōu)化管理模型: II.實(shí)例應(yīng)用[J]. 地質(zhì)評(píng)論.2011.57(3): 1-7.

[9]Zhang JL, Li YP, Huang GH. A robust simulation-optimization modeling system for effluent trading—a case study of nonpoint source pollution control[J]. Environmental Science and Pollution Research.2014.21: 5036-5053.

[10]Zheng C, Wang PP. MGO-A Modular Groundwater Optimizer Incorporating MODFLOW/MT3DMS, Documentation and User’s Guide[R]. The University of Alabama and Groundwater Research Ltd., Tuscaloosa, AL, 2003.

[11]McDonald MG, Harbaugh AW. A modular three-dimensional finite-difference groundwater flow model. Techniques of Water Resources Investigations of the U.S. Geological Survey, Book 6[R], U.S. Government Printing Office, Washington, D.C., 1988.

[12]Zheng C, Wang PP. MT3DMS: A Modular Three-Dimensional Multispecies Transport Model for Simulation of Advection, Dispersion and Chemical Reactions of Contaminants in Groundwater Systems; Documentation and User’s Guide, Contract Report SERDP-99-1[R]. U.S. Army Engineer Research and Development Center, Vicksburg, MS. 1999.