基于河流與含水層相互作用的河流污染物運移模型研究

曹紅勝 張琳

（1. 濱州市政務服務中心，山東濱州 256600；2. 濱州市引黃灌溉服務中心，山東濱州 256600）

1 引言

隨著我國城鎮化的不斷推進，居住在城市周邊地區的人口數量越來越多［1-2］，如何有效保護水資源以及避免河流污染損害人類和生態健康成為一大挑戰［3-5］。目前，固體廢物和水污染是城市河流的主要污染源，并且人口的快速增長與基礎設施建設之間的不平衡，也使得河流中存在大量的病原體和有機物質。此外，大量文獻研究表明［6-8］，盡管城市地區的地表不易透水，但污染物仍會滲入地下污染地下水源，再加上過度抽水，各含水層間的相互作用增強，加劇了污染物在地下水系統不同層間的傳播和運移，使問題進一步惡化［9］。這種變化會影響河流的水利狀況，反過來也會加劇水質污染問題，因此，需要針對河流含水層系統的管理采取更加綜合有效的措施。

確定產生河流和含水層系統流動與運移動態的環境，對解決河流污染問題至關重要。然而，目前人們對此過程了解甚少，而且缺乏適當的模型來模擬河流與含水層間的相互作用，無法獲得污染物的運移情況。因此，本文提出了一種基于河流與地下水空間耦合分布的河流走廊污染物運移模型，該模型克服了現有單一系統建模方法和集水區尺度水文模型的局限性，真實模擬了河流走廊中河流和含水層之間的相互作用及污染路徑，有利于調查現有污染的影響，為提出相應的管理策略提供依據。

2 模型的建立

2.1 典型模型分析

河流和地下水交換流量與質量的速率是隨空間和時間變化的，取決于界面處的水力和水質特征。由于流量和質量的交換速率是動態的，因此會出現非線性交換模式。

為了描述該現象，需要建立耦合地表水（SW）和地下水（GW）成分的水文模型。已有的模型包括MODHMS，InHM，HydroGenoSphere，MIKE－SHE，tRIBS。但是，上述模型的比例尺雖然可以通過假設來簡化陸上水流路徑，但對于河段和走廊比例尺研究而言，此類簡化仍然是不合適的，尤其是不適合模擬洪水和湍流等現象以及相關的輸運過程，因為它們需要更詳細的信息。

隨著計算機的進步，耦合地表非穩定水動力和地下水流動解算器的模型逐漸成熟，包括SWIFT2D－SEWAT，MODFLOW LGR VSF newCFL，MODFLOW－OWHM 和2dMb MODFLOW。

但是，二維地表水流模型和三維地下水流模型的完全耦合仍處于起步階段，均沒有考慮洪水和濕地研究的移動邊界問題，仍未在輸運和水質耦合方面取得進展。

目前有多種典型的地表水和地下水模型，通常只適用于模擬河流或含水層的動態情況。但是，由于以下原因，它們不適用于涉及瞬態條件、不規則形狀河床的復雜情況下的河流含水層流量和質量估算：

（1）大多數是為解決典型河流或含水層問題而開發的單一系統模型。雖然針對影響系統的過程進行了相應的建模，與其他非建模系統相互作用產生的過程仍需要一組邊界條件，但是邊界條件的定義在多數情況下被簡化并且不符合實際情況。

（2）河道和漫灘中河流與含水層的相互作用是一個復雜的動態過程。相互作用的方向取決于地下水位和河流水位的相對位置。但是，河流與含水層之間的流量和質量交換通常被忽略或簡化為均勻或線性變化的預定義速率，這與實際情況是不符的。

（3）可能發生界面通量方向的反轉，從而導致每個系統交替充當源或匯。因此，河流－含水層界面的物質交換的模擬，受不斷變化的壓力梯度和作為源頭的整個系統的濃度分布的影響，需要并行求解河流與地下水流動和輸運方程，同時根據相互作用的變化不斷更新地表水和地下水系統狀態。對于外部耦合河流和地下水模型，一般無法計算由河流－含水層相互依存的動力學引起的流量振蕩。

2.2 新模型構建

2.2.1 模型架構

本文所建立的模型主要將兩種運移模型（QeS2用于地表水［10］，MT3DMS 用 于 地 下 水［11］）明 確 地耦合到一個整合的河流－含水層流動模型（2dMb＋MODFLOW）中。該模型不僅規定了流動和運移方程的聯立解，還計算了水和示蹤劑在河流與地下水系統連續體中的交換通量［12-15］。模型架構及模擬過程的簡化表示如圖1 所示。

圖1 模型架構及模擬過程的簡化

地表水動力與運移模型部分分別為2dMb 和QeS2 模型的修改版本，主要用于求解二維歐氏空間中的二維深度平均淺水和運移偏微分方程。由于地形的離散化用來獲取河流走廊區域的地表流和運移動態，因此無需指定河流的位置。2 個地下水模型部分分別為MODFLOW 和MT3DMS 的修改版本，主要用于求解二維/三維地下水流動和運移偏微分方程。每個流動模型組件中所使用的時間步長需要獨立計算，并且在整個模擬過程中要不斷更新，以符合Courant（CFL）條件。一般來說，河流模型的時間步長通常為秒或1/10 秒，具體取決于網格分辨率和流量條件。而地下水模型所需的時間步長通常為小時或天，比河流模型所需的時間步長高出104倍。

2.2.2 污染負荷方案

河流的污染源類型主要分為點式和分布式兩種［16］。其中，第一種類型主要包括模擬溢出型負荷的瞬時點源和模擬出口型負荷的連續點源。后者主要包含線性斜率算法方案，因為它可以通過逐次線性逼近而適用于幾乎任何負荷模式［17］。瞬時點源和連續點源分別由以下方程式描述：

式中，Km為負荷率，%，主要取決于移動河床特征或自由水面面積。

2.2.3 模型設置與參數化

用于測試（瞬態）傳輸部分的解析解需要穩態和均勻的流場以及給定的色散系數，然而，SW 和GW流量部分產生此類條件的邊界條件與參數的定義并不直接，需要進一步考慮。在GW 流量分量的情況下，縱向壓力梯度和水力傳導率被定義為產生0.33 m/s的滲流速度，該速度使用達西定律估算，并可直接用于解析解。將橫向壓力梯度和電導率分別設置為0和Kx×10－5，以使流場為單向且沿縱向方向。在基于電導的模型部分中，參數kc/b（電導層厚度上的電導）值范圍在1～0.001 m/d 之間，在本模型中，使用的值為0.043 2 m/d。為了確保受地表水體滲透影響的含水層中流場的準確估計，地下水網格單元分辨率不應超過特征長度，因此，其值設置為10 m。

初始和邊界條件以及所用污染源的特征見表1。為了便于驗證模型，定義了2 種不同的場景“S1”和“S2”，分別測試SW 和GW 污染物傳輸部分以及SW－GW 相互作用模型部分。場景“S1”包括2 個連續的污染負荷點，1 個位于河流中，另1 個位于地下水中，分別測試SW 和GW 運輸模型部分的性能。場景“S2”旨在測試SW－GW 交互模型，僅使用邊界和初始條件進行強制。

表1 初始和邊界條件以及所用污染源的特征

3 結果與分析

用于測試傳輸組件的分析解決方案要求在整個區域內具有穩態和均勻流動條件。為了確保基準問題滿足上述條件，首先測試模型是否能夠再現所需的流場。

3.1 地下水運移

通過模擬結果可知，形成了均勻、穩定和單向的壓力梯度，并如前面所設計的一樣形成了0.33 m/s的穩定滲流速度。將GW 傳輸部分的數值結果與2個不同時刻的相應解析解進行對比，結果表明，模型結果與解析解非常吻合，其中具有微小的差異，主要是由數值色散引起的，誤差在可接受的范圍之內。

3.2 地表水運移

基于分布場景“S1”中達到穩態流動條件后的水位、速度、剪切速度、湍流黏度和擴散率分布結果，穩態條件非常接近，只有微小的偏差。由于區域大小有限及摩擦，河岸附近的速度較小。此外，由于北部和南部邊界條件的影響，導致了渦流擴散率和渦流黏度在這些區域分布也不均勻。但是，基準問題的范圍涵蓋了較廣泛的區域，該區域不受非均勻行為的影響，可根據分析解決方案對SW 傳輸組件的性能進行測試。所選區域的穩態證據中，所有與流量相關的變量在該區域中都是恒定的，并用陰影交叉線突出顯示。只有渦流擴散率沿流動方向略有變化，但與平均值的解析解中使用的值相對應。將數值解與解析解進行比較，結果表明，數值模型很好地捕捉了解析解在空間和時間上預測的總體濃度分布。微小差異是由SW 流分量模擬的驅動流變量的可變性以及數值解中遇到河岸的羽流引起的，其影響隨著時間的推移變得更加明顯。與地下水模型部分的誤差類似，可比區域的絕對誤差小于1%。

3.3 河流－地下水相互作用：流動和運輸

基于電導率模型中的輸運分量沒有解析解，因此，當外部耦合時，采用耦合模型得到的場景“S2”的結果與MODFLOW 和MT3DMS 獨立版本預測的結果進行比較。由結果可知，2 種模型方法模擬的滲透通量相似，表明所建立的新型模型準確性較高，在邊界附近的預測不匹配主要是由于表面流和運移模型中邊界條件的影響。

4 結論

（1）基于河流與地下水空間耦合分布的河流走廊污染物運移特征，本文建立了相應的數學模型，將目前已有的多個獨立模型完全耦合到河流的運移過程中。建立該模型是為了解決現有單一系統建模方法和現有尺度水文模型的局限性，提供一種分析復雜現象的方法。

（2）模擬結果表明，與現有單一系統模型結果相比，新型河流污染物運移模型能夠更準確地模擬相連河流含水層系統中的污染路徑，并且在當前城市高速發展的背景下，有利于調查現有污染的影響，并為提出相應的管理策略提供依據。

（3）雖然該模型仍存在一些局限性，但對城市河流污染的預防和恢復有著重要的指導意義，可在進一步了解水文系統之間的復雜相互作用及發生在河流含水層界面上的水力和生物化學過程基礎上，不斷更新及優化該模型。