生態(tài)清淤對于橋水庫水質(zhì)影響的數(shù)值模擬

程 揚， 賴錫軍

(1.中國科學(xué)院 南京地理與湖泊研究所 中國科學(xué)院流域地理學(xué)重點實驗室，江蘇 南京 210008； 2.中國科學(xué)院大學(xué)，北京 100049)

1 研究背景

湖泊水庫流域是人類社會經(jīng)濟發(fā)展的重要區(qū)域。流域內(nèi)頻繁的人類活動對湖泊水庫的生態(tài)環(huán)境形成了很大的壓力，水體日趨嚴(yán)重的富營養(yǎng)化現(xiàn)象就是其表現(xiàn)出來的突出問題之一[1]。水體富營養(yǎng)化威脅流域內(nèi)的供水安全和生態(tài)系統(tǒng)健康，嚴(yán)重制約了流域經(jīng)濟社會的可持續(xù)發(fā)展[2]。內(nèi)源污染是研究湖庫富營養(yǎng)化問題不可忽略的因素。即使外源污染逐漸得到有效控制，湖泊底泥的氮磷釋放仍有可能使水體保持在富營養(yǎng)化狀態(tài)[3]。生態(tài)清淤將污染底泥從湖底清除，從根本上控制內(nèi)源性污染，被廣泛應(yīng)用于國內(nèi)外的湖泊水庫治理中[4]。然而，不同湖泊實施清淤后對營養(yǎng)鹽釋放的控制效果不同，對生態(tài)清淤是否能從根本上改善富營養(yǎng)化狀態(tài)，國內(nèi)外仍存在爭議[5-6]。準(zhǔn)確認(rèn)識動態(tài)變化環(huán)境下沉積物-水界面營養(yǎng)鹽交換過程及其影響下水體水質(zhì)的變化是評估生態(tài)清淤工程能否實施的關(guān)鍵。目前，國內(nèi)對清淤效果的評價多為清淤后評估驗證[7-8]，而在清淤工程開展前缺乏必要的研究和論證，不僅為湖泊治理工作帶來隱患，同時也造成了資源資金的浪費。在生態(tài)清淤工程實施前利用模型模擬和預(yù)測清淤效果，可為生態(tài)清淤工程的決策提供重要支撐。

于橋水庫是天津市主要的飲用水源地。作為引灤入津工程的重點調(diào)蓄水庫，其在天津市供水系統(tǒng)中發(fā)揮著極其重要的作用。目前，于橋水庫整體水質(zhì)呈中營養(yǎng)-輕富營養(yǎng)狀態(tài)，且富營養(yǎng)化有加重的趨勢[9]。水庫水生態(tài)系統(tǒng)已遭到破壞，底泥對污染物的截留(自凈)能力已趨于飽和，底泥污染初顯。在不斷推進水庫上游污染治理且取得一定效果的形勢下，對于橋水庫進行清淤是進一步改善水質(zhì)的重要備選手段。

本研究以于橋水庫生態(tài)清淤為例，基于二維淺水方程和對流-擴散-反應(yīng)方程構(gòu)建了考慮底泥釋放的水動力-營養(yǎng)物質(zhì)輸移模型，分別模擬了不同工況、不同清淤范圍條件下清淤工程開展前后地表水環(huán)境的時間和空間變化情況，以期對清淤工程的設(shè)計和實施提供數(shù)據(jù)支撐和合理建議。

2 考慮沉積物-水界面物質(zhì)交換的二維水動力和營養(yǎng)物質(zhì)輸移模型

2.1 二維對流-擴散-反應(yīng)方程

采用二維對流-擴散-反應(yīng)方程模擬氮、磷營養(yǎng)物質(zhì)在湖泊水庫內(nèi)部的輸移降解過程[10-11]。其方程為：

(1)

式中：h為水深，m；C為水體中氮、磷的濃度，mg/L；u、v分別為x、y方向的流速分量，m/s；Dx、Dy分別為氮、磷在水體中x,y方向的擴散系數(shù)，m2/s；K為氮、磷在水體中的綜合降解系數(shù)，s-1；S為氮、磷的源匯項。

為了考慮底泥釋放和吸收等內(nèi)部源匯項，源項S需增加考慮沉積物底泥與表層水之間氮磷的交換。令Fex為底沉積物-水界面氮、磷交換通量，其紊動擴散理論[12]為：

(2)

考慮沉積物-水界面的氮磷交換通量，則由公式(1)可得垂向平均平面二維氮、磷營養(yǎng)物質(zhì)輸移方程，其表達式為：

(3)

2.2 水動力-水質(zhì)耦合的守恒型方程

將二維淺水方程與二維氮、磷營養(yǎng)物質(zhì)輸運方程耦合，形成的方程組可以寫成以下守恒形式[13]：

(4)

其中：

(5)

(6)

(7)

(8)

式中：U為守恒物理向量；Fx，F(xiàn)y分別為x，y方向的通量向量；S為源匯項向量；S0x，S0y分別為x，y方向的床底坡降；Sfx，Sfy分別為x，y方向的摩阻坡降；g為重力加速度，m/s2；cω為風(fēng)的阻力系數(shù)；ρa為空氣的密度，kg/m3；ω為風(fēng)速，m/s；α為風(fēng)速與y軸的夾角。

2.3 數(shù)值方法

采用非結(jié)構(gòu)網(wǎng)格有限體積法求解該二維水動力-氮磷營養(yǎng)鹽輸移的控制方程組[13-14]。首先對方程在任意控制體積做體積分，后根據(jù)高斯定理將體積分化成面積分，再進行空間離散，時間導(dǎo)數(shù)項采用一階顯格式離散，離散后方程的法向通量采用HLLC算法[13]計算。

2.4 邊界條件

以上數(shù)值方法可用于計算域內(nèi)部單元界面的求解，而對于計算域的邊界單元，其左狀態(tài)為已知狀態(tài)，而右狀態(tài)為未知狀態(tài)，此時的數(shù)值通量求解就變成了邊界黎曼問題[15]。一般可以通過根據(jù)局部流態(tài)適當(dāng)選定的輸出特征的相容關(guān)系和指定邊界條件來確定未知狀態(tài)，常見的邊界條件可以為水位、單寬流量和水位-流量關(guān)系。

對于污染物的擴散輸移，存在三類邊界：水流入口處可給定邊界濃度過程或污染負(fù)荷過程；水流出口處可給定自由輸出的Neumann邊界條件；面域上，可給定污染負(fù)荷，加入源項中計算。

3 于橋水庫模型構(gòu)建

3.1 于橋水庫庫區(qū)及計算網(wǎng)格

于橋水庫位于天津市薊州區(qū)北部的州河上，是國家重點大型水庫(圖1)。正常蓄水水位21.16 m，平均水深 4.3 m，東西長約30 km，南北寬約8 km，總庫容15.59×108m3，控制流域面積2 060 km2。于橋水庫流域四面環(huán)山，流域內(nèi)小河流眾多，匯集成較大的3個水系——黎河、沙河與淋河。庫區(qū)水下地形東高西低，在低水位時有灘地出露。為了保障于橋水庫的飲用水供水安全，天津市在東側(cè)河口建設(shè)了于橋水庫入庫河口濕地工程。該工程匯集了黎河、沙河和引灤入津的來水，經(jīng)過凈化處理后進入于橋水庫。因此，本次計算中僅考慮水庫的主庫區(qū)，不納入濕地區(qū)塊。

圖1 于橋水庫流域及水系

為了保證較高的計算效率，又能較好地反映于橋水庫地形特征，采用非結(jié)構(gòu)網(wǎng)格將計算區(qū)域共劃分為2 280個計算單元(圖2)。其中以四邊形網(wǎng)格為主，三角形網(wǎng)格為輔，各單元的邊長從40m到250m不等，入庫河道主槽單元較為精細(xì)，而在地形變化不大的灘地和地形變化較大卻常年處于淹水狀態(tài)的區(qū)域則采用了較大尺度的單元。

圖2 于橋水庫計算網(wǎng)格

3.2 初始條件和邊界條件

將在初始時刻的入庫流量和出庫水位驅(qū)動下計算24 h得到準(zhǔn)恒定條件下的水流和水質(zhì)條件作為計算的初始場。

于橋水庫的主要入庫水量和營養(yǎng)物質(zhì)來自于果河(沙河、黎河)和淋河。計算時，以果河和淋河入庫流量和營養(yǎng)物質(zhì)濃度作為邊界條件。在下游于橋水庫的出口，給定水庫出口的水位過程作為邊界條件，水質(zhì)的邊界條件設(shè)成自由輸出的無反射邊界條件。

3.3 參數(shù)率定

選擇2015年全年的水量和水質(zhì)數(shù)據(jù)對建立的于橋水庫平面二維水動力-水質(zhì)模型進行率定和檢驗。本次于橋水庫清淤計算主要關(guān)注于總氮、總磷的指標(biāo)變化。模型率定的參數(shù)主要包括水動力參數(shù)、底泥參數(shù)和水質(zhì)參數(shù)。

水動力參數(shù)：于橋水庫水動力較弱，且?guī)靺^(qū)大部分區(qū)域多為常年淹沒狀態(tài)。計算時將曼寧糙率系數(shù)統(tǒng)一取為0.021，風(fēng)拖曳力系數(shù)為0.0026。

底泥參數(shù)：沉積物上表層平均孔隙度為0.82，沉積物底層平均孔隙度為0.8，總氮和總磷的理想擴散系數(shù)均為19.8×10-10m2/s；沉積物中的氮、磷釋放速率主要由沉積物間隙水和上覆水中氮、磷營養(yǎng)物質(zhì)的濃度梯度決定，根據(jù)中國科學(xué)院南京地理與湖泊研究所開展的于橋水庫底泥釋放調(diào)查和實驗數(shù)據(jù)[16]，設(shè)定在現(xiàn)狀條件下底泥間隙水總氮和總磷濃度分別為7 mg/L和0.8 mg/L，清淤后間隙水總氮和總磷濃度分別下降至3.5 mg/L和0.01 mg/L。

水質(zhì)參數(shù)：20℃時總氮和總磷的降解系數(shù)分別為0.04 d-1和0.02 d-1，考慮降解系數(shù)隨溫度變化，根據(jù)菲爾普斯公式[17]，不同水溫條件下的降解系數(shù)值與20℃值的關(guān)系為：

KT=K20β(T-20)

(9)

式中:β為溫度校正系數(shù)。對于總氮和總磷，分別取1.07和1.04;KT和K20分別為水溫T℃和20℃時的降解系數(shù)，d-1。

3.4 模型驗證

3.4.1 水庫水質(zhì) 基于模型參數(shù)和2015年的邊界條件，計算得到了于橋水庫全年的水質(zhì)動態(tài)變化過程，采用庫中心和出庫兩監(jiān)測點的總氮和總磷實測數(shù)據(jù)對模型計算結(jié)果進行驗證(圖3)。結(jié)果表明，模型較好地反演了2015年總氮、總磷濃度的上下波動變化過程。總氮在庫中心和出庫兩站的相對誤差分別為17%和24%；總磷的相對誤差分別為26%和32%。

圖3 2015年于橋水庫總氮和總磷濃度實測值與模擬值對比

3.4.2 底泥釋放通量 由于實驗數(shù)據(jù)的限制，選擇氨氮進行底泥釋放通量的驗證。對于模型參數(shù)，水溫為20℃時氨氮的降解系數(shù)取0.1 d-1，理想擴散系數(shù)和溫度校正系數(shù)與總氮取值相同，其他水動力參數(shù)、底泥參數(shù)及水質(zhì)參數(shù)與3.3節(jié)一致，沉積物間隙水中氨氮的濃度為6.0 mg/L。利用模型求得于橋水庫2015年8月庫區(qū)底泥的氨氮釋放通量155.60 mg/(m2·d)，與文帥龍等[16]對于橋水庫沉積物-水界面交換通量的實驗研究結(jié)果129.55 mg/(m2·d)基本一致，相對誤差20.11%，能夠反映底泥向水體釋放氨氮的水平。

3.4.3 實驗設(shè)計 水庫的水動力、水質(zhì)背景狀況和清淤方案的選擇都會影響分析結(jié)果。綜合考慮于橋水庫實施清淤工程的開展條件和環(huán)境背景，分別設(shè)計了兩種工況和兩種方案。

根據(jù)天津市水利科學(xué)研究院與中國科學(xué)院南京地理與湖泊研究所承擔(dān)的“于橋水庫環(huán)境現(xiàn)狀診斷及其對上游調(diào)水的響應(yīng)研究”成果，于橋水庫底泥總磷分布較為平均，但總氮含量在壩前、北岸峰山至宋家營大堤段、劉相營段兩側(cè)以及南岸東部等區(qū)域較高，因此在綜合考慮施工便利條件、資源資金的合理利用和清淤工程效果的前提下，設(shè)計了庫周清淤和全庫清淤兩個方案：(1)庫周清淤：清淤范圍以常水位以下、高程17.8 m以上的庫周區(qū)域(圖4)，具體包括北岸、東岸和南岸東部，水庫北岸富營養(yǎng)底泥清除厚度為30 cm，水庫東部及南岸為40 cm。富營養(yǎng)底泥清除總面積13.00 km2，共計清除底泥488×104m3。(2)庫區(qū)完全清淤：全庫清除表層底泥，清除深度為35 cm。

圖4 庫周清淤范圍

2014年于橋水庫入庫河口濕地工程開工建設(shè)，建成后上游來水排入水庫前需先經(jīng)過庫前區(qū)域的濕地生態(tài)系統(tǒng)的凈化，入庫水體所含營養(yǎng)鹽等污染物得到一定程度的削減。基于此背景，本研究針對河口濕地啟用前后兩種工況分別進行模擬和分析：(1)2015年實際水流和水質(zhì)(下文簡稱“現(xiàn)狀”)：采用2015年實際觀測數(shù)據(jù)驗證計算后的水動力-物質(zhì)輸移參數(shù)，給定清淤后清淤位置的沉積物間隙水氮、磷濃度；(2)未來河口濕地正常運行(下文簡稱“未來”)：采用濕地設(shè)計的引灤入津水量54 m3/s作為入庫水量，庫水位保持在20.6m；入庫水質(zhì)總氮和總磷濃度分別為3.14 mg/L和0.11 mg/L。

4 生態(tài)清淤的影響分析

4.1 總體影響特征

在水庫外源不變的情況下，生態(tài)清淤可降低水庫氮、磷營養(yǎng)鹽濃度(表1)，改善庫區(qū)水環(huán)境富營養(yǎng)化狀況。

表1 清淤前后水體TN、TP年平均濃度

對于TN，若入庫水質(zhì)維持現(xiàn)狀，庫周清淤和全庫清淤兩種方案可使庫中心年平均濃度從1.273降至1.270和0.908 mg/L，分別下降0.23%和28.70%；若未來河口濕地正常運行，入庫水質(zhì)改善，兩種方案可分別使TN年平均濃度從1.144降至1.139和0.775 mg/L，分別下降0.48%和32.29%。無論是在現(xiàn)狀或未來的水質(zhì)背景條件下，庫周或全庫清淤方案下，庫區(qū)水體中TN濃度均有所下降。但是，庫周清淤方案效果有限，即使清除了庫岸區(qū)域的高濃度底泥，庫區(qū)水體TN濃度減少并不明顯，而全庫清淤治理效果顯著，現(xiàn)狀和未來背景下，TN濃度均大幅下降。

對于TP，現(xiàn)狀工況下，清淤前庫中心年平均濃度為0.066 mg/L，庫周和全庫清淤實施后，年平均濃度降至0.065和0.049 mg/L，分別下降0.52%和26.06%；未來工況下，TN年平均濃度為0.036 mg/L，庫周和全庫清淤實施后，年平均濃度降至0.035和0.019mg/L，分別下降1.64%和48.05%。總體變化規(guī)律與TN相似，除現(xiàn)狀背景下全庫清淤的效果比TN略低外，其他工況和清淤方案下的治理效果均優(yōu)于TN。

雖然清淤后于橋水庫水體中TN、TP均有所減少，但是值得注意的是，這種改善并不能簡單地認(rèn)為是長效的，沉積物底泥和水體中的營養(yǎng)鹽存在一個再平衡的過程，上覆水濃度變化后，底泥中的營養(yǎng)鹽含量也會隨之變化。底泥清除后，水庫來水營養(yǎng)鹽含量和本底若仍維持較高的水平，營養(yǎng)鹽仍然會在底泥累積使得底泥再次成為內(nèi)部釋放源[18]。如我國杭州西湖幾十年來經(jīng)歷了3次大規(guī)模的庫周和全湖清淤，清淤后富營養(yǎng)化狀態(tài)有所改善，但維持的時間不長[19-20]。因此，水庫的清淤工程僅是富營養(yǎng)化治理的一個方面，以改善水庫富營養(yǎng)化狀況為目的水污染防治應(yīng)該綜合外部營養(yǎng)輸入的控制、內(nèi)部生物的營養(yǎng)利用等措施共同開展。

從TN、TP濃度的年內(nèi)變化情況(圖5)看，生態(tài)清淤工程實施前后變化趨勢基本不變，清淤并不能改變水庫營養(yǎng)鹽的變化過程，庫區(qū)水體營養(yǎng)鹽的變化主要取決于外部的輸入。

4.2 水質(zhì)分布空間變化

根據(jù)上文分析，以清淤對水環(huán)境改善效果較好的2015年5月1日為例，對比分析不同工況、不同范圍的清淤工程實施前后庫區(qū)總氮和總磷的空間分布情況。

現(xiàn)狀水質(zhì)背景條件下(圖6)，清淤前TN空間分布總體表現(xiàn)為庫周較高而中心范圍較低，東部入庫口向西部出庫口逐漸降低，這與徐媛等[8]對于橋水庫富營養(yǎng)化空間分布的調(diào)查研究結(jié)果一致，表明庫區(qū)水體TN濃度受水流影響較大，上游來水導(dǎo)致入庫口TN濃度居高不下[21]。庫周清淤方案實施后，可以明顯看到庫周淺水區(qū)的TN濃度下降明顯，但庫中心區(qū)域底泥向庫區(qū)水體中釋放的營養(yǎng)鹽含量依然很高，庫區(qū)水體TN濃度下降不明顯。全庫清淤后，整個庫區(qū)水體中TN含量都明顯降低，其中，東部水域TN下降較西部地區(qū)更加明顯，可能的原因是庫區(qū)西部水域本底值較低，TN濃度下降值要比東部低。

圖5 清淤前后水體TN、TP濃度年變化過程

圖6 現(xiàn)狀水質(zhì)背景下清淤前后水體TN、TP空間分布的變化

對于TP，清淤前其分布狀況與TN總體一致，但濃度差距較小，整體水質(zhì)相對穩(wěn)定。庫周清淤后對沿岸高濃度地帶治理效果較好，全庫清淤后整個庫區(qū)水體東西分布均勻，較清淤前明顯降低，清淤工程實施效果顯著。

河口濕地正常運行后的未來水質(zhì)背景下(圖7)，上游來水污染物濃度大幅降低，外源污染的減少使水庫水質(zhì)較現(xiàn)狀背景得到有效改善。空間分布上，上游來水TN濃度的降低導(dǎo)致入庫口TN濃度大幅降低，上下游TN污染空間差異減小。TP空間分布與TN總體一致。

上述結(jié)果表明，庫周清淤仍只能使庫周淺水區(qū)和東部范圍的N、P濃度下降，全庫清淤效果明顯，東西部N、P濃度差進一步降低，庫周淺水區(qū)和中心水體的N、P濃度基本一致。

圖7 未來水質(zhì)背景下清淤前后水體TN、TP空間分布的變化

5 結(jié)論和討論

(1)本文基于二維淺水方程和對流-擴散-反應(yīng)方程，建立了考慮沉積物-水界面物質(zhì)交換的二維水動力和營養(yǎng)物質(zhì)輸移模型。結(jié)果表明，模型能夠有效地模擬水庫底泥氮、磷營養(yǎng)鹽的動態(tài)釋放及其與外界污染輸入共同作用下水庫的氮、磷營養(yǎng)鹽的時空分布特征，可用于生態(tài)清淤工程實施效果的預(yù)測。

(2)于橋水庫底泥是重要的污染來源，底泥的清除對水庫氮、磷營養(yǎng)鹽含量有較大影響。無論是在現(xiàn)狀或未來的入庫污染負(fù)荷條件下，庫周清淤或全庫清淤方案實施后，于橋水庫水體的氮、磷營養(yǎng)物質(zhì)含量均會有所下降，但年內(nèi)變化過程基本不變。庫周清淤對沿岸帶淺水區(qū)水體氮、磷有很好的改善作用，但對庫區(qū)的整體改善效果不佳(濃度降低率小于2%)；而全庫清淤后，整個庫區(qū)水體的氮、磷濃度均明顯下降(濃度降低率達到26%～48%)，治理效果相對明顯。清淤后底泥營養(yǎng)鹽仍有可能隨著上覆水中營養(yǎng)鹽含量的上升而重復(fù)積累，使得新生底泥在較短時期內(nèi)表現(xiàn)出污染“源”的特征。因此，全面的水庫污染防治應(yīng)該綜合內(nèi)源、外源污染控制、改善生態(tài)結(jié)構(gòu)等措施共同開展。

(3)該模型的建立為生態(tài)清淤工程的效益評估提供了科學(xué)預(yù)測工具。在應(yīng)用模型分析清淤的環(huán)境生態(tài)效應(yīng)時，應(yīng)充分開展沉積物-水界面物質(zhì)交換實驗獲取準(zhǔn)確的模型參數(shù)，提高模型預(yù)測分析的可靠性。