胡 艷 海,周 林 飛
(沈陽農(nóng)業(yè)大學(xué) 水利學(xué)院,遼寧 沈陽 110161)
長期以來,水庫水質(zhì)狀況和庫區(qū)水生態(tài)環(huán)境受到廣大學(xué)者的關(guān)注,其中以水庫水動(dòng)力條件為基礎(chǔ)進(jìn)行水庫污染物質(zhì)運(yùn)移消散過程模擬是主要的研究方式。國外關(guān)于水質(zhì)模型的研究起步于20世紀(jì)60~80年代,計(jì)算機(jī)技術(shù)的迅速發(fā)展促進(jìn)了雙向耦合模型對(duì)水質(zhì)模型變化所產(chǎn)生的影響,與此同時(shí)水質(zhì)模型也伴隨著水動(dòng)力系數(shù)、空間變量、主要的外界影響因素(例如溫度等)的引進(jìn)而變得更加完善[1]。20世紀(jì)80年代后,多維數(shù)值模型求解伴隨二維水質(zhì)模型的成熟得以實(shí)現(xiàn)。國內(nèi)關(guān)于水質(zhì)模型的研究比較晚,20世紀(jì)80年代,由于官廳水庫遭受嚴(yán)重的重金屬污染事件,由此引起了多個(gè)領(lǐng)域的研究學(xué)者對(duì)水庫水質(zhì)問題的關(guān)注,并意識(shí)到利用水質(zhì)數(shù)值模型解決水環(huán)境相關(guān)問題的重要性[2]。因此,眾多學(xué)者針對(duì)不同類型的湖泊、湖泊中各種水質(zhì)問題、富營養(yǎng)化問題以及不同的水質(zhì)因素構(gòu)建水質(zhì)模型,最終實(shí)現(xiàn)水質(zhì)模型可行性驗(yàn)證,同時(shí)為進(jìn)一步實(shí)現(xiàn)模型的簡化、促進(jìn)數(shù)值模型的應(yīng)用進(jìn)展創(chuàng)造條件[3-5]。龔春生等[6]通過模擬玄武湖混合流的水質(zhì)動(dòng)態(tài)變化過程為底泥污染和淺水湖泊混合流水質(zhì)變化模擬提供新途徑。20世紀(jì)90年代,水質(zhì)模型逐漸被廣泛應(yīng)用到水質(zhì)優(yōu)化方案設(shè)計(jì)方面[7]。21世紀(jì)后,水質(zhì)模型逐漸向多維、多層次方向發(fā)展[8-11]。例如,楊具瑞[12]構(gòu)建的垂直分層的二維水質(zhì)模型,計(jì)勇[8]構(gòu)建的水動(dòng)力-水質(zhì)-底泥模型,廖臨毓[13]構(gòu)建MIKE21水動(dòng)力-水質(zhì)耦合模型對(duì)COD、TN、TP三種污染物在3種調(diào)水方案下的流場和濃度場時(shí)空變化情況進(jìn)行模擬,馬寧[14]構(gòu)建水動(dòng)力-水質(zhì)耦合模型模擬水庫水動(dòng)力和水質(zhì)狀況并獲得水質(zhì)改善方式。
石佛寺水庫位于遼河干流,承擔(dān)防洪、供水、生態(tài)等重任。因此,掌握庫區(qū)水動(dòng)力條件以及水質(zhì)狀況顯得尤為重要。本次研究以石佛寺水庫作為研究對(duì)象,通過總結(jié)前人經(jīng)驗(yàn)構(gòu)建水庫的二維水動(dòng)力-水質(zhì)耦合模型,將石佛寺水庫來水量作為水庫水動(dòng)力模型的主要影響因素,以水質(zhì)模型濃度場變化和水動(dòng)力模型流場變化作為運(yùn)算基礎(chǔ),模擬運(yùn)算了庫區(qū)汛期和非汛期兩種工況下溶解氧(DO)、氨氮(NH3-N)、總氮(TN)、五日生化需氧量(BOD5)、化學(xué)需氧量(CODCr)、總磷(TN)6種水質(zhì)指標(biāo)的運(yùn)移消散狀況。最后,根據(jù)分析模擬結(jié)果提出水環(huán)境綜合治理的具體措施,用以指導(dǎo)水庫運(yùn)行,使遼河水流經(jīng)水庫后能夠達(dá)到地表水Ⅲ類水質(zhì)標(biāo)準(zhǔn)。
石佛寺水庫位于遼寧省沈陽市沈北新區(qū)黃家鄉(xiāng)和法庫縣依牛堡鄉(xiāng),地理坐標(biāo)為東經(jīng)123.427°~123.520°,北緯42.144°~42.188°,距沈陽市區(qū)約47 km,如圖1所示。石佛寺水庫具有防洪、供水、生態(tài)三大功能,其中最主要的是防洪,為典型的河道型平原水庫,是遼河干流上僅有的控制性水利工程。為改善區(qū)域生態(tài)環(huán)境,凈化遼河上游來水,實(shí)現(xiàn)其生態(tài)功能,2009年起在石佛寺水庫實(shí)施了生態(tài)建設(shè)工程,并在庫區(qū)內(nèi)遼河左岸種植水生植物蘆葦、蒲草和荷花,同時(shí)進(jìn)行生態(tài)蓄水,蓄水面積16.13 km2,蓄水位46.2 m。
為了了解遼河水流經(jīng)石佛寺水庫后水質(zhì)凈化狀況,在水庫的入口、中間和出口布設(shè)水質(zhì)監(jiān)測(cè)點(diǎn)位,分別標(biāo)記為采樣點(diǎn)1、采樣點(diǎn)2、采樣點(diǎn)3,如圖1所示。各監(jiān)測(cè)點(diǎn)監(jiān)測(cè)的水質(zhì)指標(biāo)為溶解氧(DO)、氨氮(NH3-N)、總氮(TN)、五日生化需氧量(BOD5)、化學(xué)需氧量(CODCr)、總磷(TP)等11個(gè)指標(biāo)。DO利用多參數(shù)光譜水質(zhì)分析儀測(cè)定,NH3-N采用水楊酸鈉法測(cè)定,TN采用全自動(dòng)間斷化學(xué)分析儀測(cè)定,BOD5采用稀釋與接種法測(cè)定,CODCr采用重鉻酸鹽法測(cè)定,TP采用硫酸鉀—鉬銻抗法測(cè)定。監(jiān)測(cè)頻率為每月2次,監(jiān)測(cè)標(biāo)準(zhǔn)以GB 3838-2002《地表水環(huán)境質(zhì)量標(biāo)準(zhǔn)》為依據(jù)。
監(jiān)測(cè)時(shí)間為2009~2017年。石佛寺水庫作為水源地,水標(biāo)準(zhǔn)執(zhí)行地表水Ⅲ類。
MIKE21水動(dòng)力模塊是在二維數(shù)值求解淺水方程的基礎(chǔ)上建立的[15-18]。二維水動(dòng)力數(shù)值模擬分析基于以下4項(xiàng)基本假定展開計(jì)算:不可壓縮性假定、靜水壓力假定、Boussinesq假定、Reynolds值均布假定。水動(dòng)力模塊控制方程應(yīng)用沿水深積分的不可壓縮的Navier-Stokes方程,其主要包括二維平面流的連續(xù)性方程、動(dòng)量方程[19],如式(1)~(3) 所示。本文選用單元中心的有限體積法作為數(shù)值解法,并采用MIKE21FM的非結(jié)構(gòu)化三角形網(wǎng)格進(jìn)行模型求解。
(1)
(2)
(3)
圖1 石佛寺水庫位置及水質(zhì)監(jiān)測(cè)點(diǎn)位置Fig.1 Location and water quality monitoring points of Shifosi Reservoir
(4)
(5)
(6)
庫湖水質(zhì)控制方程是全面考慮污染物質(zhì)水動(dòng)力學(xué)、物理運(yùn)輸和對(duì)流擴(kuò)散等過程并將質(zhì)量守恒原理作為前提建立起來的。該控制方程的數(shù)學(xué)表達(dá)式如式(7)~(8) 所示[20]。
(7)
S=Qs(Cs-C)
(8)
二維水動(dòng)力模型作為水質(zhì)模型的基礎(chǔ),只有建立水動(dòng)力模型并滿足模型參數(shù)率定和模型驗(yàn)證,方能加載水質(zhì)模型,其主要原理在于水動(dòng)力模型能夠提供水質(zhì)模型水質(zhì)濃度輸出的流速等水動(dòng)力條件。
水庫水動(dòng)力影響因素包括:進(jìn)出水條件、風(fēng)場、地形、庫內(nèi)水工構(gòu)筑物以及水庫的岸線形態(tài)等。其中,水庫的地形和岸線形態(tài)作為相對(duì)固定的因素,來水量為主要因素。
3.1.1模擬范圍
本文選用1∶5 000石佛寺水庫電子版地形圖作為地形資料,利用Auto CAD和南方CASS提取高程點(diǎn)并生成.xyz格式陸地land和水深water數(shù)據(jù)文件用于網(wǎng)格生成。采用Beijing_1954_GK_Zone_21N投影坐標(biāo),控制最大三角形面積20 000 m2、最小三角形角度30°以及最大節(jié)點(diǎn)個(gè)數(shù)150 000個(gè),利用MIKE ZERO網(wǎng)格生成器生成網(wǎng)格,計(jì)算網(wǎng)格如圖2所示。完成水深water文件插值并利用網(wǎng)格分析器進(jìn)行網(wǎng)格分析,獲得合理準(zhǔn)確網(wǎng)格文件并由MIKE Animator生成模擬區(qū)域的三維地形,如圖3所示。
圖2 石佛寺水庫計(jì)算網(wǎng)格Fig.2 Computational grid of Shifosi Reservoir
圖3 石佛寺水庫三維地形Fig.3 Three-dimensional topography of Shifosi Reservoir
3.1.2定解條件
(1) 初始條件。模擬時(shí)長為61 d,即2017年6月1日08:00:00至2017年7月31日08:00:00,時(shí)間步長為3 600 s。由于模型在短時(shí)間內(nèi)計(jì)算趨于穩(wěn)定,故在初始時(shí)刻t=0時(shí),庫區(qū)初始流速u和v設(shè)置為0,初始水位為δ=46.200 m。其他數(shù)據(jù)為實(shí)測(cè)值。
(2) 邊界條件。根據(jù)MIKE軟件前處理方法,上游邊界為石佛寺水庫入口處實(shí)測(cè)流量數(shù)據(jù),下游邊界為石佛寺水庫出口處實(shí)測(cè)水位數(shù)據(jù);水質(zhì)模型利用岸壁法確定閉邊界,即沿法線方向上的濃度值和流速值均為零;主要根據(jù)底床摩擦應(yīng)力確定底床邊界條件;采用可避免模型發(fā)散的“干濕水深判別法”確定移動(dòng)邊界,即水深大于濕水深0.1 m時(shí)按水域處理,水深小于干水深0.005 m時(shí)按陸域處理。
3.1.3模型的率定與驗(yàn)證
(1) 模型率定。選擇2017年6,7月石佛寺水庫入口處實(shí)測(cè)水位率定模型參數(shù),率定結(jié)果為:曼寧系數(shù)45 m1/3/s(MIKE軟件的曼寧系數(shù)與糙率值互為倒數(shù))、渦黏系數(shù)0.28。
圖4 出口合成流速模擬值與實(shí)測(cè)值對(duì)比Fig.4 Comparison of the simulated value and measured value of the exit flow velocity
圖5 入口水位模擬值與實(shí)測(cè)值對(duì)比Fig.5 Comparison of the simulated value and measured value of the entrance water level
3.2.1水質(zhì)模板選擇與水質(zhì)指標(biāo)選取
本文選用改進(jìn)的Eco Lab內(nèi)置水質(zhì)模板WQ with nutrients對(duì)石佛寺水庫水質(zhì)進(jìn)行模擬。該水質(zhì)模板包含狀態(tài)變量7個(gè),常量39個(gè),作用力6個(gè),輔助變量9個(gè),過程21個(gè)。
綜合考慮實(shí)際水質(zhì)監(jiān)測(cè)數(shù)據(jù)的超標(biāo)情況和生態(tài)數(shù)值實(shí)驗(yàn)室(ECO Lab)內(nèi)置模板中狀態(tài)變量的可模擬情況確定水質(zhì)指標(biāo)為:DO,NH3-N,BOD5,TN,CODCr,TP,其中TN和TP作為衍生結(jié)果輸出。
3.2.2定解條件與模型參數(shù)
(1) 初始條件。模擬范圍、模擬網(wǎng)格系統(tǒng)及模擬時(shí)間同水動(dòng)力模型??紤]模型整體一致,故在初始時(shí)刻t=0時(shí),水質(zhì)指標(biāo)的初始濃度值設(shè)置為0。
(2) 邊界條件。閉邊界法線方向上設(shè)定流速和濃度值均為零;開邊界,將水庫進(jìn)出口處水質(zhì)指標(biāo)實(shí)測(cè)濃度生成的時(shí)間序列文件作為上下游邊界條件。
(3) 模型參數(shù)。建立生態(tài)模型時(shí)污染物的運(yùn)移消散過程必然會(huì)發(fā)生,因此水質(zhì)的擴(kuò)散系數(shù)以及衰減系數(shù)的設(shè)置顯得至關(guān)重要[21]。其中,擴(kuò)散系數(shù)在各個(gè)方向上的上下限通常根據(jù)其與各網(wǎng)格速度分量之間的線性比例關(guān)系而確定;降解系數(shù)則通常根據(jù)原型觀測(cè)法[22]獲得,然后利用環(huán)境部門檢測(cè)結(jié)果和相關(guān)研究成果進(jìn)行修正和再次修正,最后確定各水質(zhì)指標(biāo)的降解系數(shù)如下:DO為0.075 14 mg/d;NH3-N為0.113 1 mg/d;BOD5為0.140 53 mg/d;TN為0.014 95 mg/d;CODCr為0.102 05 mg/d;TP為0.042 12 mg/d。
3.2.3水質(zhì)模型率定與驗(yàn)證
(1) 模型率定。選擇2017年6~7月石佛寺水庫主槽區(qū)實(shí)測(cè)五日生化需氧量、化學(xué)需氧量、氨氮、總氮、總磷率定模型參數(shù),率定結(jié)果為:橫向擴(kuò)散系數(shù)5 m2/s、縱向擴(kuò)散系數(shù)5 m2/s,化學(xué)需氧量衰減系數(shù)0.078 5 mg/d,五日生化需氧量衰減系數(shù)0.108 1 mg/d,氨氮衰減系數(shù)0.087 mg/d,溶解氧衰減系數(shù)0.057 8 mg/d,總氮衰減系數(shù)0.011 5 mg/d,總磷衰減系數(shù)0.032 4 mg/d。
(2) 模型驗(yàn)證?;谀P蛥?shù)的成功率定,將對(duì)模型進(jìn)行驗(yàn)證,模型驗(yàn)證時(shí)選用2017年5,7月出口處實(shí)測(cè)生化需氧量、化學(xué)需氧量、氨氮、總氮、總磷數(shù)據(jù)驗(yàn)證模型,驗(yàn)證結(jié)果如圖6所示。采用統(tǒng)計(jì)學(xué)中的偏差統(tǒng)計(jì)法,得到五日生化需氧量、化學(xué)需氧量、氨氮、溶解氧、總磷、總氮實(shí)測(cè)數(shù)據(jù)與計(jì)算模擬數(shù)據(jù)之間平均誤差分別為4.04%,7.39%,5.10%,3.80%,3.05%,1.40%,均在可接受范圍內(nèi),表明實(shí)測(cè)值與模擬計(jì)算值擬合程度較好,完全能夠滿足精度要求,所以本文所構(gòu)建的水質(zhì)模型是有效可行的。
利用驗(yàn)證合理的模型進(jìn)行模擬,將1號(hào)、2號(hào)和3號(hào)監(jiān)測(cè)點(diǎn)作為特征點(diǎn),用于反映石佛寺水庫流速和流場特點(diǎn),監(jiān)測(cè)點(diǎn)布置如圖1所示。
4.1.1模擬方案
本次水動(dòng)力數(shù)值模擬共分為兩種工況:即汛期和非汛期。根據(jù)石佛寺水庫多年運(yùn)行經(jīng)驗(yàn),6~9月為汛期,非汛期不考慮冰凍期,故非汛期僅模擬為10~11月和4~5月。汛期和非汛期初始流速均設(shè)置為0,初始水位均設(shè)置為46.200 m。汛期和非汛期上游邊界為石佛寺水庫入口處實(shí)測(cè)流量數(shù)據(jù),下游邊界為石佛寺水庫出口處實(shí)測(cè)水位數(shù)據(jù)。模擬時(shí)長非汛期為2016年10月1日~11月30日、2017年4月1日~5月31日,共122 d;汛期為2017年6月1日~9月30日,共122 d?,F(xiàn)選取汛期和非汛期水動(dòng)力模型運(yùn)行末期時(shí)典型流場分布情況進(jìn)行分析。
4.1.2模擬結(jié)果分析
(1) 根據(jù)圖7可知,汛期水位大部分區(qū)域介于46.280~46.256 m之間,非汛期水位大部分區(qū)域介于46.190~46.185 m之間,即汛期庫區(qū)水位大部分僅略高于非汛期。原因是雖然汛期水量較大,但會(huì)通過下游閘門控制水庫水位,使其保持在46.200 m左右,保障生態(tài)需水。
圖6 各水質(zhì)指標(biāo)模擬值與實(shí)測(cè)值對(duì)比Fig.6 Comparison of the simulated value and measured value of water quality indicators
(2) 圖8為汛期和非汛期不同水量條件下庫區(qū)典型流場流速分布情況,結(jié)果表明兩種工況下水庫具有較好的流動(dòng)狀態(tài)且未出現(xiàn)死水區(qū)(流速≤0.3 cm/s時(shí)為死水區(qū))和回水區(qū)。對(duì)比可以看出汛期的流速整體大于非汛期,說明不同的來水量水庫內(nèi)不同區(qū)域內(nèi)的流速變化明顯。兩種工況下,流速較大區(qū)域均出現(xiàn)在上游入口區(qū)域和下游出口區(qū)域,這是因?yàn)檫^水?dāng)嗝婷娣e小從而加速了水體運(yùn)動(dòng);過水?dāng)嗝婷娣e較大的區(qū)域,流速比較緩慢,占庫區(qū)面積也大;流速較小區(qū)域出現(xiàn)在左右岸岸邊,與岸邊地形相關(guān)??傮w來看,庫區(qū)內(nèi)水體流動(dòng)性良好、水流不急,有利于水體的交換以及污染物的去除。
圖7 水位分布(單位:m)Fig.7 Water level distribution
圖8 流場流速分布(單位:m/s)Fig.8 Flow field velocity distribution
來水量是石佛寺水庫水動(dòng)力的主要影響因素,因此選擇非汛期和汛期兩種不同情景對(duì)其進(jìn)行水質(zhì)模擬分析。具體模擬時(shí)間段為:非汛期為2016年10月1日至11月30日,2017年4月1日至5月31日,共122 d;汛期2017年6月1日至9月30日,共122 d。初始時(shí)刻水質(zhì)指標(biāo)濃度設(shè)置為0,采用水質(zhì)指標(biāo)濃度值生成的時(shí)間序列文件作為上下邊界條件。
4.2.1模擬結(jié)果
本次數(shù)值模擬獲得兩種情景下水庫DO,NH3-N,BOD5,TN,CODCr和TP在不同時(shí)間段的濃度分布情況。汛期情景下各污染物在模擬末期時(shí)分布圖如圖9所示,非汛期情景下各污染物在模擬末期時(shí)分布如圖10所示。
圖9 汛期各水質(zhì)指標(biāo)濃度分布(單位: mg/L)Fig.9 Distribution of water quality indicators in flood season
圖10 非汛期各水質(zhì)指標(biāo)濃度分布(單位: mg/L)Fig.10 Distribution of water quality indicators in non-flood season
4.2.2分 析
汛期和非汛期除DO以外各污染物濃度整體分布情況為:入口>中間>出口,原因是石佛寺水庫具有凈化污染物質(zhì)的作用,其中在非汛期盡管水流流速減緩,但是水體中的某些污染物質(zhì)會(huì)伴隨泥沙發(fā)生沉降,污染物質(zhì)被存儲(chǔ)在底泥中,進(jìn)而實(shí)現(xiàn)凈化水質(zhì)的作用。DO為出口>中間>入口,流經(jīng)水庫復(fù)氧能力增強(qiáng),水質(zhì)得到改善。
在汛期NH3-N,TN,CODCr和BOD5濃度小于非汛期的濃度,原因一是汛期庫區(qū)內(nèi)水生植物長勢(shì)茂盛,微生物具有很強(qiáng)的活性,有很好的水質(zhì)凈化作用;原因二是汛期有較大的水量,流速也較大,自然增強(qiáng)了水體的自凈能力。在汛期DO的濃度大于非汛期的濃度,原因是水體復(fù)氧途徑由原來的的大氣復(fù)氧,增加了水生植物向水體中輸送氧氣和藻類等沉水植物的光合作用產(chǎn)生氧氣,而且汛期流速快增加了水體的復(fù)氧能力。TP的濃度汛期和非汛期變化不明顯。
水質(zhì)從入口到出口逐漸得到凈化,說明遼河水流經(jīng)水庫后水質(zhì)得到了改善。但分析出口水質(zhì)數(shù)據(jù),仍有指標(biāo)在一定時(shí)間段內(nèi)達(dá)不到地表水Ⅲ類水質(zhì)要求,因此,仍需加強(qiáng)庫區(qū)內(nèi)的污染綜合控制,可以從以下3個(gè)方面進(jìn)行考慮:
(1) 控制外源污染。進(jìn)水污染物濃度為導(dǎo)致水庫污染的主要因素,因此可通過控制面源污染和點(diǎn)源污染,盡量減少經(jīng)降雨、下滲、徑流過程形成的污染源匯入水庫。在外源污染超過水庫的凈化能力時(shí),可根據(jù)石佛寺水庫水動(dòng)力水質(zhì)模擬結(jié)果進(jìn)行調(diào)水,進(jìn)一步改善水質(zhì)。
(2) 減少內(nèi)源污染。水庫在非汛期水流流速會(huì)減緩,泥沙會(huì)攜帶一些污染物質(zhì)沉降而儲(chǔ)存在水庫底泥中,尤其在水流相對(duì)緩慢的非汛期,為避免再次污染可以通過對(duì)庫底底泥疏浚來控制內(nèi)源的匯入。
(3) 種植水生植物。汛期的流場和水質(zhì)指標(biāo)濃度分布圖疊加后得出,生長水生植物的區(qū)域水流受阻流速變緩,水生植物的凈化作用使污染物濃度降低;在無水生植物的區(qū)域,水流流速較大的區(qū)域水質(zhì)指標(biāo)的濃度值相對(duì)較小些。因此,應(yīng)該在水流流速較小的區(qū)域種植凈化能力較強(qiáng)且阻水較弱的水生植物,在高流速區(qū)域種植阻水能力較強(qiáng)的水生植物。
(1) 石佛寺水庫水動(dòng)力數(shù)值模擬結(jié)果表明:來水量是水庫水流狀態(tài)的主要影響因素。汛期和非汛期兩種工況下,汛期的流速整體大于非汛期,說明不同的來水量水庫內(nèi)不同區(qū)域內(nèi)的流速變化明顯;兩種工況下水庫均具有較好的流動(dòng)狀態(tài),未出現(xiàn)死水區(qū)和回水區(qū),且水流不急,有利于水體的交換以及污染物的去除。
(2) 石佛寺水庫水質(zhì)數(shù)值模擬結(jié)果表明:遼河水流經(jīng)水庫后逐級(jí)得到了凈化,而且凈化效果汛期優(yōu)于非汛期。但在一定時(shí)間段內(nèi)仍有指標(biāo)達(dá)不到地表水Ⅲ類水質(zhì)要求,需加強(qiáng)庫區(qū)內(nèi)的污染綜合控制,具體措施如下:控制外源污染和調(diào)水凈化;通過底泥疏浚減少內(nèi)源污染;種植具有凈化功能的水生植物,在水流流速較小的區(qū)域種植阻水較弱的水生植物,在高流速區(qū)域種植阻水能力較強(qiáng)的水生植物。