胥瑞晨，逄 勇，胡祉冰，朱天依

(1.河海大學(xué)淺水湖泊綜合治理與資源開(kāi)發(fā)教育部重點(diǎn)實(shí)驗(yàn)室，江蘇 南京 210098；2.河海大學(xué)環(huán)境學(xué)院，江蘇 南京 210098)

太湖2007年暴發(fā)的大面積藍(lán)藻水華[1]，給太湖流域的人民生命財(cái)產(chǎn)安全造成了巨大的危害。藻類的運(yùn)動(dòng)和生長(zhǎng)與水動(dòng)力之間存在極其密切的關(guān)系，水齡作為水動(dòng)力的主要表現(xiàn)形式之一，一直以來(lái)被廣泛關(guān)注[2]。水齡的概念被定義為一個(gè)水包(一種不具有擴(kuò)散降解能力的理想水團(tuán)或者水質(zhì)點(diǎn))離開(kāi)其規(guī)定年齡為0的區(qū)域后到達(dá)某一位置所需要的時(shí)間[3-4]。水齡模型最早是由Bolin等[5]提出的，主要體現(xiàn)了水質(zhì)點(diǎn)由一點(diǎn)運(yùn)動(dòng)到另一指定位置所需的時(shí)間。此后根據(jù)不同計(jì)算精度及應(yīng)用范圍的相關(guān)需求，發(fā)展了很多關(guān)于水齡計(jì)算的水環(huán)境數(shù)學(xué)模型，其中箱式模型最早由Kashiwai[6]提出，具有計(jì)算簡(jiǎn)單的優(yōu)點(diǎn)，但同時(shí)存在假設(shè)條件和計(jì)算區(qū)域的限制，應(yīng)用范圍受到了很大的限制；質(zhì)點(diǎn)跟蹤模型[7]雖然可以實(shí)時(shí)跟蹤水質(zhì)點(diǎn)的移動(dòng)軌跡，但缺乏擴(kuò)散功能；隨機(jī)游動(dòng)模型[8]是質(zhì)點(diǎn)模型的優(yōu)化，加入了湍流擴(kuò)散模塊，但計(jì)算結(jié)果不夠準(zhǔn)確；時(shí)間尺度和年齡模型效果很好，缺點(diǎn)是計(jì)算復(fù)雜，很多模型方法不具備求解功能，限制了模型方法之間的參考與應(yīng)用[5]；對(duì)流擴(kuò)散模型[9]的應(yīng)用范圍較為廣泛，是一種基礎(chǔ)的水齡研究模型，可以統(tǒng)計(jì)出水質(zhì)點(diǎn)的駐留時(shí)間和交換矩陣等，但由于假設(shè)條件的限制，存在一定的數(shù)值震蕩和耗散情況，需要根據(jù)特殊情況進(jìn)行后期改進(jìn)。

Takeoka[10]早在1984年便對(duì)水力停留時(shí)間、水齡、過(guò)境時(shí)間和周轉(zhuǎn)時(shí)間等概念進(jìn)行了討論和研究。對(duì)于穩(wěn)定的大型湖泊而言，湖泊的水力停留時(shí)間通常被定義為一種描述水動(dòng)力能力的參數(shù)，通過(guò)水力停留時(shí)間的概念，可以準(zhǔn)確地描述湖體內(nèi)物質(zhì)的交換時(shí)間，其計(jì)算方法為湖體容積除以入湖流量[11]，但它一般只能表現(xiàn)一個(gè)湖體中水體交換的平均能力和周期，在湖體各個(gè)區(qū)域或者點(diǎn)位上，無(wú)法給出具體的變化規(guī)律。對(duì)于大型湖泊來(lái)說(shuō)，各個(gè)區(qū)域的水力停留時(shí)間存在很大的差異，F(xiàn)eng等[12]采用Matlab編程的方法對(duì)太湖水體的水力停留時(shí)間進(jìn)行了初步的推算，得到了湖體具體點(diǎn)位的水力停留時(shí)間。

目前，水力停留時(shí)間和水齡之間具體的關(guān)系尚不明確，且研究較少。當(dāng)前水齡的模擬主要是采用EFDC模型，模擬出來(lái)的水齡仍需要通過(guò)試驗(yàn)實(shí)測(cè)進(jìn)行驗(yàn)證，存在很大的不確定性。本研究基于水齡和水力停留時(shí)間的概念推導(dǎo)得到了兩者之間的具體關(guān)系，并探討了影響太湖水齡的決定性因素。

1 研究區(qū)域概況

太湖(119°8′E～122°55′E，30°5′N～32°8′N)為中國(guó)第三大淡水湖泊，流域總面積為3.69萬(wàn)km2，湖泊總面積2 427.8 km2，水域面積為 2 338 km2。太湖平均水深1.9 m，最大水深不超過(guò) 3 m，是典型的大型淺水湖泊，2008年太湖實(shí)際年均水齡為300 d左右[13]。太湖區(qū)域夏季以東南風(fēng)為主導(dǎo)風(fēng)向，冬季以西北風(fēng)為主導(dǎo)風(fēng)向，風(fēng)速平均在 3.5～5.0 m/s。平均年降水量1 100～1 500 mm[14]。貢灣(120°13′E～120°24′E，31°19′N～31°28′N)為太湖東北的一個(gè)湖灣，水域面積為163.8 km2，平均水深為1.82 m，近年來(lái)水生態(tài)環(huán)境由于周邊經(jīng)濟(jì)發(fā)展而惡化，作為蘇錫兩地的金墅港、石帆港、錫東、南泉四大水廠的水源地，保證貢灣的水質(zhì)穩(wěn)定達(dá)標(biāo)極為重要[15]。目前，貢灣是“引江濟(jì)太”的第一站，貢灣水齡的影響研究對(duì)后期湖體的系統(tǒng)治理和維護(hù)起著至關(guān)重要的作用。

水利部太湖流域管理局提供的數(shù)據(jù)資料(圖1)顯示，太湖入湖水量2007—2017年總體呈現(xiàn)上升的態(tài)勢(shì)，其中2017年比2007年多出22.25億m3，這和“引江濟(jì)太”工程有著密不可分的聯(lián)系。由于太湖水位主要仍為人工控制，為了確保防洪安全，水位常年主要保持在3.10～3.40 m之間[16]，與2007—2017年太湖實(shí)測(cè)水位數(shù)據(jù)基本一致(圖1)，這意味著出湖水量幾乎也是同比增加的。近年來(lái)太湖實(shí)施了退漁還湖等措施，使得太湖水面面積增大，總體蓄水量也有逐漸增大的趨勢(shì)[17]。

圖1 太湖2007—2017年入湖水量和水位變化

2 研究方法

2.1 水齡與水力停留時(shí)間的關(guān)系確定及模型構(gòu)建

對(duì)流擴(kuò)散模型計(jì)算水力停留時(shí)間的殘余方程組[18]為

(1)

式中：rt為某一時(shí)刻某點(diǎn)的物質(zhì)質(zhì)量與初始質(zhì)量之比或質(zhì)量濃度與初始質(zhì)量濃度之比；M0為初始物質(zhì)質(zhì)量；Mt為某一時(shí)刻物質(zhì)質(zhì)量；ρ0為初始物質(zhì)質(zhì)量濃度；ρt為某一時(shí)刻物質(zhì)質(zhì)量濃度；τt為湖體的平均水力停留時(shí)間。由于現(xiàn)實(shí)中物質(zhì)質(zhì)量濃度不可能因物質(zhì)降解而下降至0，因此按照模型方程設(shè)定當(dāng)質(zhì)量濃度降低到e-1時(shí)，記為水力停留時(shí)間，也可稱為 E-folding 時(shí)刻[19]。

對(duì)比水力停留時(shí)間和水齡的概念，水齡的含義即為某水團(tuán)到達(dá)并完全替代原有水團(tuán)所用的時(shí)間。因此，在對(duì)流擴(kuò)散模型中的水齡其實(shí)就是理想狀態(tài)下粒子質(zhì)量濃度降解到0或者物質(zhì)質(zhì)量濃度上升到粒子質(zhì)量濃度所需的時(shí)間，相當(dāng)于一個(gè)100%的轉(zhuǎn)換關(guān)系(圖2)。綜上可知水齡與水力停留時(shí)間的關(guān)系為

atx=τtx/(1-e-1)=τtx/0.63

(2)

式中：atx為t時(shí)刻x點(diǎn)的水齡；τtx為t時(shí)刻x點(diǎn)的水力停留時(shí)間。

圖2 水齡與水力停留時(shí)間的關(guān)系

為了驗(yàn)證兩者的關(guān)系，本文采用同種示蹤劑，輸入初始質(zhì)量濃度為360 mg/L，擴(kuò)散系數(shù)為1 m2/s，降解系數(shù)為1.57×10-7s-1，保證出入湖流量平衡，根據(jù)模型計(jì)算得到湖體容積。基于水齡與水力停留時(shí)間的關(guān)系，根據(jù)湖體容積與入湖流量計(jì)算出全湖的年平均水齡，使用Origin8.5的非線性擬合工具箱，率定出一階降解方程具體參數(shù)取值，湖體每個(gè)網(wǎng)格的瞬時(shí)水齡便可以通過(guò)不同網(wǎng)格內(nèi)物質(zhì)的質(zhì)量濃度反演確定，其中模擬的平均水齡以最后一天的網(wǎng)格水齡與網(wǎng)格體積的加權(quán)平均得到[20]。具體數(shù)學(xué)模型構(gòu)建如下：

(3)

(4)

(5)

2.2 敏感性分析方法

由于實(shí)際影響參數(shù)較少，本文采用多元線性逐步回歸方法進(jìn)行敏感性分析。為了更快更好地分析參數(shù)的敏感性，即參數(shù)對(duì)結(jié)果方差的貢獻(xiàn)率，對(duì)所有的數(shù)據(jù)提前進(jìn)行標(biāo)準(zhǔn)化處理。在供選擇的各自變量中，依各自變量對(duì)因變量作用的大小，即偏回歸平方和的大小，由大到小把自變量依次逐個(gè)引入。對(duì)每一個(gè)變量進(jìn)行逐步回歸檢驗(yàn)，得到具體的響應(yīng)關(guān)系[21]?；貧w結(jié)束后所得方程組即為所求得的最優(yōu)回歸方程[22]：

(6)

(7)

3 結(jié)果與討論

3.1 模型參數(shù)率定與驗(yàn)證

基于水齡與水力停留時(shí)間的概念，最終確定了兩者之間存在63%的比例關(guān)系，采用一階指數(shù)衰減的方式擬合，其R2可達(dá)到1.000的擬合精度。為了再次驗(yàn)證其正確性，結(jié)合2013年出入湖流量、風(fēng)場(chǎng)、平均水位與庫(kù)容基礎(chǔ)數(shù)據(jù)，計(jì)算出平均水齡為 308 d，按照本文模型算法，最終率定得到了模型關(guān)鍵參數(shù)K=4.399×10-3，A=361.587 05，其全年半降解時(shí)間為157.58 d。

采用望虞河引水的邊界流量100 m3/s為模型輸入條件進(jìn)行驗(yàn)證，對(duì)太湖總體的水齡分布特征進(jìn)行模擬，結(jié)果如圖3所示。由于太湖邊界和島嶼的基礎(chǔ)地形構(gòu)建略有不同而對(duì)水齡模擬計(jì)算產(chǎn)生了一定的影響，但總體時(shí)空分布特征和水齡等值線與郝文彬等[23-24]研究成果基本一致，平均水齡誤差不超過(guò)10%，可以認(rèn)為此模型參數(shù)用于太湖地區(qū)的模擬研究基本可行。根據(jù)設(shè)定的水齡率定點(diǎn)位，模擬結(jié)果分別為196 d(東南風(fēng))和308 d(北風(fēng))，與Li等[25]在2011年采用EFDC模型模擬的200 d(東南風(fēng))和 310 d(北風(fēng))、黃春琳等[26]在2014年采用EFDC模型模擬的194 d(東南風(fēng))和308 d(北風(fēng))研究結(jié)果基本一致。

圖3 太湖水齡特征模擬結(jié)果(單位：d)

根據(jù)本文模型模擬得到太湖2007—2017年水齡變化如圖4所示，可以看出太湖水齡并不是一成不變的，2013年模擬值與實(shí)際計(jì)算值之間的差距為2%，總體最大誤差不超過(guò)5%，屬于模型誤差可接受范圍。模擬結(jié)果表明，年均水齡從2007年的 315 d 降低到2017年的257 d，下降了18.41%。由于近年來(lái)太湖周邊的水利工程與人工調(diào)度，使得太湖年均水齡發(fā)生了很大的變化，其中2017年太湖水齡實(shí)際計(jì)算值為255 d，模擬誤差為0.8%，再次驗(yàn)證了本文模型的正確性。同時(shí)證明了利用一階指數(shù)衰減函數(shù)換算出來(lái)的水齡模型是可行的，并且相較于原有的復(fù)雜水齡模型來(lái)說(shuō)，在以對(duì)流擴(kuò)散為主的淺水湖泊中，這種模型具有更廣泛的普適性。

圖4 太湖2007—2017年水齡變化

3.2 參數(shù)敏感性分析

為了確定出入湖流量、風(fēng)向、風(fēng)速和初始水位4個(gè)主要邊界條件對(duì)水齡影響的具體權(quán)重，基于均勻分布限制條件，根據(jù)實(shí)際邊界條件取值范圍，隨機(jī)取值組成了200組邊界數(shù)據(jù)，并對(duì)其進(jìn)行水齡模擬。標(biāo)準(zhǔn)化處理模擬結(jié)果后，采用多元線性逐步回歸方法進(jìn)行水齡敏感性分析，得到R2=0.951>0.7，可認(rèn)為回歸結(jié)果可信，驗(yàn)證了出入湖流量、風(fēng)向、風(fēng)速和初始水位均為影響水齡的主要因素，其影響權(quán)重分別為64.41%、8.56%、9.57%和17.46%。以年均水齡較大的2013年作為研究基礎(chǔ)年，發(fā)現(xiàn)不同風(fēng)向?qū)μg的影響較為穩(wěn)定(表1)，平均水齡在300 d左右，東風(fēng)對(duì)平均水齡的影響程度達(dá)到了5.3%。當(dāng)外加出入湖流量從0 m3/s上升到400 m3/s時(shí)(圖5(a))，西北風(fēng)時(shí)的入湖流量對(duì)平均水齡影響最為顯著，為35.3%，西南風(fēng)對(duì)平均水齡影響最小，為27.6%。

表1 不同風(fēng)向下的太湖平均水齡及不同水齡占比

(a) 太湖

為了進(jìn)一步對(duì)同等水位下的其他3個(gè)因素做具體的敏感性解釋，分別對(duì)太湖與貢灣做了模擬研究，發(fā)現(xiàn)貢灣的水齡隨風(fēng)向變化較整體湖區(qū)來(lái)說(shuō)更為明顯(圖5(b))，西北風(fēng)時(shí)的水齡在4種流量下的變化率最大，能達(dá)到45.5%，南風(fēng)時(shí)由于大流場(chǎng)強(qiáng)烈抑制了貢灣的水體流動(dòng)，降低了流量對(duì)貢灣的水齡影響程度。這說(shuō)明周邊地形對(duì)區(qū)域小流場(chǎng)的形成具有很大的影響，風(fēng)向?qū)τ趨^(qū)域性水體流動(dòng)的影響要大于整體湖區(qū)，后期對(duì)于貢灣的望虞河引水應(yīng)盡量選在刮西風(fēng)和西北風(fēng)為主的時(shí)間段，以增加引水的環(huán)境效益。通過(guò)出入湖流量的變化比較，也能明顯看出在出入湖流量與初始水位不變的情況下，風(fēng)向的改變對(duì)湖體水質(zhì)的影響會(huì)低于外部條件同時(shí)變化的時(shí)間段，證明風(fēng)場(chǎng)與出入湖流量對(duì)于水動(dòng)力的改變有一定的協(xié)同作用。

由于太湖藍(lán)藻暴發(fā)期主要集中在夏季，主導(dǎo)風(fēng)向?yàn)闁|南風(fēng)，通過(guò)圖6(a)可以發(fā)現(xiàn)當(dāng)風(fēng)向不變時(shí)，風(fēng)速對(duì)太湖水齡的影響隨速度增大而降低，根據(jù)不同風(fēng)向，風(fēng)速會(huì)在一定的大小之前起決定性作用，之后則是出入湖流量成為了主導(dǎo)因素，這是由于風(fēng)場(chǎng)前期對(duì)水流的作用較大，促進(jìn)了水循環(huán)導(dǎo)致的，當(dāng)風(fēng)速到達(dá)一定程度時(shí)，整個(gè)流場(chǎng)由于出入湖流量不變形成了一個(gè)固定的結(jié)構(gòu)，水動(dòng)力變化的主要因素則變成了出入湖流量。而對(duì)于貢灣來(lái)說(shuō)，風(fēng)速對(duì)局部區(qū)域產(chǎn)生的影響并不都是穩(wěn)定上升或者下降的(圖6(b)與圖7)，而是具有一定的波動(dòng)性，在東南風(fēng) 5 m/s 情況下，貢灣的小流場(chǎng)會(huì)與太湖大流場(chǎng)發(fā)生一定程度的抵消，導(dǎo)致偏南部分的流速減緩，同時(shí)使得局部形成小流場(chǎng)，抑制了水體交換，導(dǎo)致區(qū)域水齡升高，同理也能解釋梅梁灣在夏季水體交換能力不足的原因。

(a) 太湖

(b) 貢灣

圖7 區(qū)域水齡在東南風(fēng)不同風(fēng)速下的變化關(guān)系

4 結(jié) 論

a. 構(gòu)建了水齡與水力停留時(shí)間之間的轉(zhuǎn)換模型，率定驗(yàn)證了相關(guān)參數(shù)，確認(rèn)了水力停留時(shí)間與水齡之間存在63%的比例關(guān)系。采用一階指數(shù)衰減的方式擬合，其R2達(dá)到了1.000，進(jìn)而將質(zhì)量濃度與水齡進(jìn)行數(shù)值轉(zhuǎn)化，模擬得到每一個(gè)網(wǎng)格的水齡反演結(jié)果，能真實(shí)反映太湖具體點(diǎn)位的水齡變化情況。太湖每年的水齡不完全相同，由于“引江濟(jì)太”工程的實(shí)施，2017年太湖年均水齡降低至 257 d，2007—2017年的水齡模擬值與實(shí)際計(jì)算值之間的誤差不超過(guò)5%，表明本文改進(jìn)后的模型在太湖水齡研究方面具有可行性。

b. 多元線性逐步回歸法分析表明，出入湖流量、風(fēng)向、風(fēng)速和初始水位為影響水齡的主要因素，其權(quán)重分別占到了64.41%、8.56%、9.57%和17.46%。大型淺水湖泊的水動(dòng)力在地形不變的前提下，出入湖流量對(duì)水齡的影響最大，并且風(fēng)場(chǎng)對(duì)于具有小流場(chǎng)的貢灣和梅梁灣來(lái)說(shuō)，存在很大的不確定性，在后期工程引水期間應(yīng)結(jié)合具體氣象預(yù)報(bào)來(lái)確定引水的時(shí)間。