基于有限差分模型模擬地下水-湖泊相互作用的研究

龔嘉臨， 朱 松， 李非里， 金贊芳

(1．浙江工業大學環境學院， 浙江 杭州 310032；2.浙江省水利水電勘測設計院， 浙江 杭州 310002)

0 引言

目前，全球幾乎所有湖泊面臨著富營養化及其引起的藍藻水華問題，嚴重影響著湖泊的水質和生態功能[1-2]。在國內，湖泊富營養化問題也早已成為一個較為突出的環境問題，如太湖、滇池、鄱陽湖等大型湖泊均已處于富營養或重富營養狀態。而且，無錫飲用水供水危機等事件的發生也表明了我國湖泊水體富營養化問題的嚴重。然而，盡管發達國家對湖泊的富營養化研究和治理已經經歷了半個世紀，但污染湖泊的水質和生態恢復至今仍未能完全實現。這是因為湖泊的富營養化過程是一個不可逆且長期產生影響的過程[3]。首先，造成湖泊富營養化的外源有很多種，包括農業污染、工業污染、生活污水，大氣沉降、地下水入湖等。其次，湖泊底泥和腐爛的水生植物等內源物質也會將富集的營養物質釋放，從而影響湖泊。而在地表水源受到普遍污染的現在，人們進一步開發并利用地下水，導致淺層地下水逐漸成為氮元素的貯存庫。其中，硝酸鹽污染更是已成為地下水污染的重要問題[4-5]。由于地下水與湖泊緊密聯系，存在著水力交換，所以無論是哪種系統發生變化，都會對關聯的另一系統造成影響。

隨著對湖泊生態系統的不斷研究，湖泊周邊地下水對湖泊富營養化問題的影響受到了國內外研究人員的廣泛關注。自上世紀70年代，國外的研究人員就開始了對地下水與湖泊之間相互作用關系的研究。由于賦存于地面下的巖石空隙中的地下水在一般情況下不能被人們直接觀察，所以研究人員只有通過在現場的水文地質勘察與地下水動態的持續監測才能了解區域地下水系統的存在形式及運動規律。在研究初期，研究人員對地下水-湖泊系統的研究手段主要是利用滲流計和壓力計等現場監測儀器來估算淺層地下水和湖泊之間的交換量。

隨著計算機技術的發展和利用，越來越多的研究人員使用可視化模型來定量化模擬區域地下水和關聯湖泊之間的相互作用[6-8]。其中，有限差分模型MODFLOW作為應用最廣泛的模型軟件，受到了眾多研究人員的青睞。本文通過對地下水水流模型MODFLOW、溶質運移模型MT3DMS和湖泊模擬方法的介紹，總結其功能與應用特點。

1 有限差分模型MODFLOW簡介及特點

MODFLOW模型（Modular Three-dimensional Finitedifference Ground-water flow model）最初是由美國地質調查局（US Geological Survey）的MCDONALD等在上世紀80年代所開發出來的一款數值模擬程序，用于模擬孔隙介質中地下水流動。

MODFLOW基于離散化方法，將連續的區域空間和時間分為一系列離散的格點。在一定精度下，將區域空間離散為有限個規則單元的集合體，使得每個網格單元上的各種參數值近似為常數。在離散化之后，整個區域的計算問題就被概化為了有限個單元網格的計算。結合時間離散之后，利用有限差分方程構成線性方程組。對線性方程組聯立求解之后，獲得各個離散點上的水頭數值近似解。這種方法可以較好地反映區域地下水的狀態，且具有較高仿真度。

MODFLOW的一個最顯著特點是其模塊化的結構[9]。它在主程序的基礎上可以整合一系列相對獨立的子程序包，而這些子程序包是由具有類似功能的子程序所組成的。同時，研究人員可以根據實際情況選用其中合適的子程序包對地下水系統進行數值模擬。這種模塊化結構讓整個程序更容易理解、修改，并可對其添加新的子程序包，如河流子程序、模擬水平流動障礙的子程序包等。模塊化的結構使得MODFLOW的功能和應用范圍大大拓展。模塊化結構的MODFLOW不僅僅局限于最初的模擬地下水在孔隙介質中的流動，在合理使用的情況下可以解決裂隙水的流動等實際問題。

目前，使用有限差分模型MODFLOW的軟件主要有GMS（楊百翰大學環境模型實驗室和美國軍工水道實驗室合作開發），Visual Modflow（加拿大Waterloo水文地質公司開發）和Visual Modflow flex。隨著地下水數值模擬軟件的廣泛應用，國內外研究人員主要通過MODFLOW模擬和預測區域地下水環境的動態變化，實現對實際問題的評價與動態預測，從而制定合理的措施解決問題[10]。SONG等[11]通過構建模型來確認開采計劃是否合理，從而確保巴拉素水源地的安全。馬馳等[12]根據模型結果預報西華水源地地下水位，獲得合理的水源地開采量，為西北地區地下水資源的評價開創了先例。尉鵬翔[13]利用Visual Modflow軟件對北京某污染場區進行建模（利用MODFLOW和MT3D），分析了污染物的遷移，為該污染場區地下水評價和污染治理提供了科學依據。ELKRAIL[14]結合GIS和MODFLOW對地下水的脆弱性進行評價，用于指導松花江所在區域地下水資源的開發利用和保護。

因為MODFLOW的最終目的是為了模擬實際的水文地質情況，故在使用該模型之前需要盡可能多地收集研究區的水文地質資料。在此基礎上，整個模型的構建步驟見圖1。

圖1 MODFLOW模型構建一般步驟

2 溶質運移模型MT3DMS

為了研究地下水和湖泊富營養化的關系，在利用MODFLOW獲得水流模型的基礎上，研究人員還需要根據溶質運移模型獲得相應的水質模型來進一步估算地下水對湖泊的污染負荷。

20世紀90年代，ZHENG開發了模塊化的三維溶質運移模型MT3D。這款模型開發完成之后就廣泛地應用于地下水中的溶質運移模擬和研究。隨后，可插入多重污染物組分生化反應模塊程序的MT3DMS由ZHENG等[15]于1998年開發。MT3DMS的功能在MT3D的基礎上得到了極大的擴散，其獨特之處在于囊括了有限差分法、基于粒子示蹤的歐拉-拉格朗日法以及高階有限差分TVD法3類主要的溶質運移解算法。之所以將3種解算法綜合起來是因為沒有一種單獨的方法可以使用所有的運移條件[16]。

與MT3D一致，MT3DMS也是以塊體中心法有限差分水流模型（如MODFLOW）為背景開發的。同時，MT3DMS保留了MT3D的模塊化結構。而且和MODFLOW一樣的是，MT3DMS同樣是模塊結構。其模塊化的結構可以單獨模擬運移的對流、彌散/擴散、源匯問題和化學反應問題，有利于合理分配內存；處理其它運移過程與化學反應的模塊可以隨時插入程序，而不需要改動現有代碼。因此，MT3DMS可以模擬含水層中的對流、彌散-擴散過程，也可以模擬地下水污染物質在運移過程工程的生物反應和化學反應。趙貝等[17]通過MODFLOW和MT3DMS構建了河北省東南部地區的一個垃圾填埋場的三維數學模型，并模擬預測了該區未來30 a內的特征污染物Cl-的遷移情況。模擬結果表明，污染物濃度隨著離污染源的距離增大而減少，主要危害地下水流方向上的地表水體和淺層地下水；湖泊對污染物具有持續捕獲作用；開采地下水會影響污染物的遷移速率和方向。GUSYEV等[18]利用MODFLOW，MT3DMS和同位素氚對新西蘭的西陶波湖捕集區進行地下水年齡模擬，用于構建流域地下水年齡分布。

3 湖泊模擬方法

自然環境中的地下水和湖泊存在著水力交換，從而導致地下水影響著湖泊的水平衡、營養化水平和酸緩沖能力。而湖泊周邊的地下水是作為補給源還是排泄匯，這取決于湖泊與周邊地下水的水力梯度。1979年，BORN等[19]依據湖泊流域的水文循環特征，對地下水-湖泊交換定義了3種不同的類型，見圖2。第1種類型是地下水穿流，湖泊的部分湖底接受地下水，而部分湖底發生滲漏；第2種類型是地下水進入湖泊，對湖泊進行補給；第3種是湖泊水對地下水進行補給。一般情況下，地下水-湖泊交換類型隨季節變化而發生變化。其中，第1種類型在山區流域中出現的頻率較多，且在具有恒定流向地表水體的平原地區也常出現，具體表現為湖泊在上游的入湖區接受地下水和下游的出湖區滲漏進入地下水；第2種類型常出現在沙漠地區或蒸發量遠遠大于降雨量的湖泊流域，當湖泊的周邊地勢較高時，湖泊位于流域的地勢較低處，湖泊水體的出流主要依靠蒸發的作用，入流則主要來自部分降水和地下水流；第3種類型通常在地表水位主要補給的湖泊流域出現，湖泊水體的出流受地勢影響，或湖底的滲透性能較好，從湖底補給的地下水能夠較快地排走，從而使含水層中的地下水位維持在比湖泊水位更低的水平。

圖2 地下水-湖泊的交換形式

研究地下水和湖泊的相互作用有利于了解整個系統，從而為解決實際問題提供依據。而在數值模擬中，研究人員也可以通過以下幾種方式模擬地下水和湖泊的相互作用[20-22]。

3.1 湖泊區域作為高滲透系數的含水層

當進行MODFLOW模擬時，將湖泊作為構建的網格模型含水層的一部分，同時設置其所在區域的滲透系數為高滲透系數。李勇[23]和莫美仙[24]結合該方法和第2種方法，把滇池作為模型含水層中的一部分，設置湖泊所在區域單元具有高滲透系數，且將湖泊區域設置為邊界，從而來模擬滇池。前者利用ZONGBUDGET模塊計算出滇池東側地下水入湖水量為8 073.8 m3/d，硝態氮和總磷入湖量分別為32.4和0.91 t/a，北側地下水入湖水量為8 463.6 m3/d，硝態氮和總磷入湖量分別為54.72和0.832 t/a；后者計算出昆明盆地北部地下水入湖量在6.5×106～6.9×106m3/a之間，并根據地下水攜帶的氮、磷濃度平均值乘以入湖水量估算氮和磷的入湖量范圍分別為27.13～ 406.26，1.84～ 10.04 t/a。

這種方法模擬較為便利，可以計算湖泊水位，并模擬穩態和瞬態問題。但其缺陷在于只能模擬滲漏湖泊，無法包含其它的地表水流動，且容易出現質量不守恒問題。

3.2 湖泊作為邊界條件輸入模型

湖泊作為邊界條件輸入模型即把湖泊這一系統看成是MODFLOW模擬區域的其中一個已知邊界。在這種情況下，研究人員可以利用MODFLOW中的定水頭邊界或通用水頭邊界來表示湖泊。潘田[25]將研究區東北部連接的太湖定義為定水頭邊界，以太湖水位作為水頭條件，從而實現研究區地下水與太湖的水力聯系。同時，結合MT3DMS進行溶質運移模擬獲得太湖流域地下水硝態氮補給太湖。結果表明，太湖周邊淺層地下水補給太湖水量為0.2億m3/a，占總排泄量的6%；太湖周邊淺層地下水中的硝酸鹽氮入湖負荷約為14 245 t/a，占補給太湖氮污染負荷的27.8%。郭玉川[26]將博斯騰湖作為變水頭邊界，湖泊水頭值根據實際水位分時段賦值，通過構建模型，從而分析模擬了基于生態安全水位的不同調控方案的水均衡情況，對當地地下水位的調控方案進行優選。

該方法模擬的湖泊水位穩定且利于模擬滲漏和排水湖泊。但是，由于湖泊水位已經作為固定邊界輸入，在整個模型模擬過程中無法體現湖泊的水位變化且一般不包括其它地表水流。同時，當其作為定水頭邊界輸入模型的時候，湖泊邊界其實相當于無限的源匯項。

3.3 河流程序包代替湖泊模塊

河流模塊作為MODFLOW邊界條件中的常用模塊，也可以用于模擬其它地表水體。這種方法與通用水頭邊界條件一樣都是以外部應力形式導入模型。但與通用水頭邊界導入不同的是，使用河流模塊模擬湖泊等其它地表水體的時候可以定義水體范圍，同時在結果中可以顯示水域的變化情況及流場。而缺陷在于可能出現質量不守恒情況。

3.4 LAK3程序包

上文已說明MODFLOW是一種模塊化的結構，這為湖泊程序包的開發和使用提供了便利。LAK湖泊程序包（LAK1,LAK2,LAK3）從上世紀90年代開發至今，已被國內外研究人員廣泛應用。KROHELSKI等[27]利用 MODFLOW 和LAK3程序包來模擬威斯康星州西北部Dane鎮的3個淺水滲透湖泊（Fish湖，Mud湖和Crystal湖）及其周圍的淺層地下水。模擬結果表明，Fish湖水位在抽水條件下會降低；在沒有持續抽水的情況下，湖泊水位又會在幾年內恢復。FEINSTEIN等[28]模擬了威斯康星州西北部Waukesha鎮的湖泊、井和淺層地下水的相互作用。模擬結果表明，地下水是Beaver湖最重要的入流補給部分（等于總流入的59%）；對于Pine湖和North湖來說，地下水仍然是重要的來源（分別占總流入的16%和5%），但小于降水和地表水的貢獻。當在地表61 m以下以213.85 L/min速率抽取地下水5 a，對湖泊水位的影響也并不大，模擬的Pine湖的水位和出流量僅分別減少了0.91 cm和3%。但是，加入抽水速率增加到910 L/min，模擬的Pine湖的出流量明顯減少，相對基本流量條件下降了14%。在我國，應用LAK3模塊的情況較少。程春龍等[29]利用MODFLOW和LAK3模塊建立了天津擬建寧車沾水庫的耦合模型。模擬結果表明，水庫水位下降速率受到當地的蒸發強度影響較大，受滲漏影響極小，年水位降幅1.16 m。

LAK3程序包作為最新的模擬湖泊程序包，具有進行穩態和瞬態模擬的功能，且可以用于模擬地表水流和溶質運移。同時，和河流程序包一樣，LAK3程序包的湖泊模擬可以在結果中顯示。與前面3種方法相比，LAK3輸入較為復雜，包括湖泊水位、底泥層厚度、滲透系數、降雨、蒸散、徑流等參數。其缺陷主要有2點：①當滲透系數較高時，模擬具有潛在的不穩定性；②在沒有用戶干涉的情況下，可能出現錯誤結果。

4 存在的問題及發展趨勢

根據文獻檢索發現，近年來有關區域地下水-湖泊系統的三維可視化模型研究數并不少。這些研究涉及的湖泊及周邊地下水的水文地質條件復雜多樣。利用MODFLOW等模型的構建為區域地下水和湖泊的水資源預測規劃和管理提供了重要的依據。但是，綜合分析可以發現以下問題的存在。

（1）研究范圍常以人為邊界劃定。由于實際地質條件、資料收集、工作量和計算機容量所限，研究人員在進行模擬研究的過程中往往很難找到自然邊界作為邊界條件，這就導致模擬區域往往比實際的地下水含水層系統要小。這種情況下，MODFLOW模型常采用人為邊界。因為人為邊界的引入，也導致了模型預測的困難。

（2）水文數據資料收集和驗證困難。和國外相比，國內的野外實驗基地沒有形成規模，這也導致水文地質資料往往沒有長期的監測、考證資料，無法形成數據積累。這也導致了國內構建的模型難以進一步修正優化，達不到長期預測的效果。

（3）數據處理往往以多年平均、短期觀測資料為主。同樣是因為實際地質條件、資料收集和工作量所限，研究人員往往只能盡可能將收集到的水文地質信息進行合理概化，以求在模擬的時候接近實際情況。因此，這也造成了根據數據資料的差異，不同研究人員的研究成果往往存在一定差異。同時，這也必然導致研究成果與實際的偏差，以至于很多情況僅僅限于科學研究。

5 結論

MODFLOW這類用于模擬實際環境狀態的模型離不開各種大數據的支持。隨著我國社會經濟的不斷發展，對于城市及周邊地區地質數據的調查和完善提出了新的要求。進入21世紀以來，中國地質調查局和各地方政府合作開展的“城市地質調查項目”就是為了解決城市發展過程中的地質問題，并已在多個重要城市展開試點調查工作。這為將來的進一步地質研究提供了便利。而隨著互聯網的發展，大數據資料庫的建立使區域水文地質資料得到了不斷完善和擴充，加快了信息交流。同時，計算機技術的發展也使計算機容量不斷增大，為更大范圍的區域地下水模擬提供了可能。長遠來看，數據資料的擴充和共享交流對模型研究具有重要意義。