謝培,張亞輝,喬飛*
1.環境基準與風險評估國家重點實驗室,中國環境科學研究院
2.中國環境科學研究院環境分析技術測試中心
河流是陸源物質和污染物質進入海洋的重要通道,陸地進入海洋的物質85%由河流搬運[1]。河口區是水動力復雜的活躍水文生態系統和陸海相互作用的集中地帶,由于大多沿岸河口區人類活動頻繁,河口水體儲存了大量污染物質,在漲落潮不穩定的環境特性下,水體污染物質不斷振蕩難以擴散。
大亞灣是位于珠江口東側的亞熱帶半封閉性海灣,近幾十年隨著沿岸人口增長、經濟及相關產業的快速發展,大亞灣水體環境富營養化問題突出[2]。在大亞灣12條入海河流中,最大的河流輸入是西北部的淡澳河,已有研究表明,大亞灣入海河流中淡澳河入海污染物通量占比最高[3]。淡澳河為感潮河段,水流除受上游徑流影響外,還受海洋潮流的強烈影響,不同潮流特性對河口區影響差異大,河口內段受邊界約束潮波變形,因此明確該區域的潮流場及水動力特性,對于控制水環境污染具有極其重要的意義。
目前數值模擬是較為普遍且有效的環境建模方法,由于天然水體是三維的,三維數值模擬能全面反映水環境系統流場的基本特征。在三維數值模擬模型通用性上,國內鮮有模型能夠實現復雜邊界準確擬合以及水流分層可視化精細模擬,當前應用最為普遍的三維數值模擬模型包括美國的環境流體動力學(EFDC)模型和FVCOM模型、丹麥的Mike系列模型以及荷蘭的Delft3D模型等[4]。陳青毅等[5]利用三維數值模擬河口水庫的水流流動和水體交換,制定了合理的水庫引排水方案;姜恒志[6]應用EFDC模型建立了太湖三維水動力數學模型,分析了太湖流場特性;齊珺等[7]建立三維水動力模型考察了長江水系武漢段河床的沖淤變化。盡管國內外學者已在典型河流湖庫就水動力特征開展了一些研究[5-8],但鮮有對淡澳河下游感潮河段、河口及其毗鄰海區水體交換的整體探究,有關淡澳河流域水環境問題探究多集中于海水入侵帶來的災害,如土壤鹽漬化[9],以及陸地-海洋相互作用對地下水排泄的影響[10-11]等方面。為系統了解該水域復雜的水動力和咸淡水交互作用,本研究構建了淡澳河下游感潮河段、河口及其毗鄰海區為一體的三維水動力水質模型,旨在揭示徑流潮流作用下的水動力過程和內在影響規律。此外,考慮到充分認識不同徑流量與潮流交互作用引發河口水環境的變化,以系統全面地了解淡澳河口特性,筆者定量估算了不同水期鹽水入侵特征、咸淡水交互與水質的響應關系,以期為淡澳河流域水環境治理提供依據。
淡澳河為分洪渠與響水河匯合后至白壽灣出海口段,于 114°33'57"E、22°43'29"N 河口處匯入大亞灣,是大亞灣入海河流中流量最大且水質最差的河流。淡澳河河長約14 km,下游河口處寬度為500 m,上游河寬約為200 m,為較典型的上窄下寬喇叭狀河口。淡澳河河段由于海水的頂托作用,污染物不易向外擴散,河段水質較差。淡澳河上游的淡澳分洪渠為人工分洪渠,早期主要為緩解淡水河及其下游的洪水壓力,近年來隨著惠陽區的發展,大量工業廢水和生活污水通過淡澳分洪渠和淡澳河進入大亞灣。
1.2.1 模型選擇
模擬區域包含淡澳河感潮河段及整個河口區,河道與河口在橫向尺度上存在很大的差異,并且區域水動力過程復雜,對模型適用性的要求較高。EFDC模型垂向上采用σ坐標變換,能較好擬合近岸復雜的岸線和地形,在河流、水庫、湖泊、河口等多類水體都得到廣泛應用,并體現了很好的適用性。由于淡澳河口水動力條件復雜,采用長期應用和改進的河口海岸三維數值模型,該模型已長期應用于河口水動力過程和鹽水入侵等方向的研究,并取得了眾多成果[6-8,12-17]。
1.2.2 模型構建與驗證
1.2.2.1 網格設置
模型采用水平曲線非正交網格,范圍包括整個淡澳河河口、大亞灣和臨近海區,上游流量邊界設定為響水河與淡澳河交匯處上游約600 m處,下游為外海開邊界。由于模擬范圍較廣,河道、河口和海區各區域水體空間尺度變化較大,為盡可能滿足不同區域計算精度同時保證較高的計算效率,模型采用變尺度網格系統,河口內網格分辨率為15 m×15 m~40 m×80 m不等,河口外網格較稀疏,分辨率為80 m×80 m~125 m×200 m不等。垂向采用σ坐標,均勻分為5層。考慮到淡澳河河口區域淺灘較多,模型選用動邊界處理,運用干濕判別法處理潮灘移動邊界,臨界水深取0.2 m。
1.2.2.2 條件設置
模型上游設置為流量型邊界,下游設置為水位邊界,模型上下游邊界、姚田橋和虎爪斷面(國控點位)位置分布如圖1所示。本研究模擬2019年1月1日——12月31日淡澳河口——大亞灣區水動力水質過程,上邊界采用日流量過程,下邊界由潮位驅動,其潮流特性主要受來自太平洋的潮波制約,潮位數據資料(由各分潮調和常數合成得到)取自全球潮汐數值模式NAOTIDE計算結果(http://www.miz.nao.ac.jp/staffs/nao99/index_En.html)。模式岸線采用2019年海圖數字化資料及2019年水深實測資料。水質指標主要選擇COD、氨氮、TP等,水質濃度采用實測值。

圖1 網格邊界條件分布Fig.1 Grid division and boundary conditions of research region
此外,為提高模型計算效率,本次模擬采用動步長,基準步長設置為0.1 s,實際步長由模型根據實時克朗數自動調整,經統計平均時間步長約為1 s,平水期和枯水期稍長,豐水期稍短。
1.3.1 主要參數設置
模型設置水平渦黏系數為Ax=Ay=100 m2/s,擴散系數為1×10-5m2/s,粗糙高度是表征河流底板對水體阻力的變量,相當于曼寧公式中的糙率,采用單元劃分法[18]求解糙率,最終取值為0.001~0.01 m。
1.3.2 模型驗證
為保證模型的可靠性和適用性,水動力水質模型驗證采用3組不同水期資料:1月(枯水期)、6月(豐水期)和10月(平水期)。主要水動力水質參數驗證站位為虎爪斷面,計算結果表明模擬值與實測值均呈現較好的相關性(表1)。水位和鹽度模擬與實測值表現為強相關,相對誤差小于5%,主要水質指標中氨氮相對誤差小于10%,COD和TP相對誤差略高(在可接受范圍之內),主要原因可能在于COD和TP輸入量統計不夠,導致局部斷面模擬濃度出現一定的誤差。

表1 水動力水質參數模擬誤差統計Table 1 Error analysis of simulated hydrodynamic water quality parameters
淡澳河水體模擬縱向分為5層,分別提取底層(k=1)、中部(k=3)和表層(k=5)鹽度數據,通過不同水層含鹽量計算淡水與鹽水比例。淡澳河上游含鹽量約為0‰,中游含鹽量約為15‰,下游含鹽量約為35‰,河流入海區鹽度高于河流上中游鹽度,符合入海河流鹽度分布規律。以淡澳河虎爪斷面月內不同水層鹽度分布為例(圖2),計算得到該斷面表層、中部和底層的淡鹽水比例。水體底層含鹽量約為10.8‰,淡鹽水比例約為2∶1;水體中部含鹽量約為6‰,淡鹽水比例約為7∶1;水體表層含鹽量約為1.7‰,淡鹽水比例約為 16∶1。

圖2 水體垂向鹽度分配Fig.2 Distribution of vertical salinity of water body
大亞灣潮流性質以不規則半日混合潮型為主,由于受地形影響,外海潮波傳至大亞灣內變形較大,以致潮汐日現象顯著,潮波主導著水流運動,潮流動力作用較弱。灣內潮流基本是順水道(主航道)的往復流[19],漲潮平均流速為 0.065 m/s,落潮為0.073 m/s,落潮流速比漲潮流速大,與吳仁豪等[20]的研究結果一致。水動力模擬結果顯示,淡澳河虎爪斷面的表層平均流速為0.153 m/s,底層平均流速為0.094 m/s,表層流速約為底層流速的1.6倍。漲潮流向以NW(320°)為主,落潮流向以 SE(100°)為主,海浪以涌浪為主,風浪為輔,平均波高為0.2 m,常浪向為SE。港灣內缺少長期波浪觀測資料,短期觀測資料顯示,6級風時,港內水域波高僅0.6 m,12級臺風時,波高僅1.3 m。
大亞灣月內潮汐變化如圖3所示。從圖3可以看出,每月8~10 d為日潮,20~22 d為半日潮。不同潮型日潮(17日)和半日潮(11日)的潮流場分布如圖4所示。日潮漲急流速為0~0.963 m/s,其中表層漲急最大流速為1.019 m/s,底層漲急最大流速為0.878 m/s;半日潮漲急流速約為0~0.918 m/s,其中表層漲急最大流速為1.003 m/s,底層漲急最大流速為0.767 m/s,不同水層均表現為日潮漲急流速大于半日潮漲急流速。

圖3 大亞灣潮汐變化趨勢Fig.3 Trend of tide change in Daya Bay

圖4 日潮和半日潮漲急流域潮流場分布Fig.4 Distribution of tidal current at diurnal and semidiurnal tides during typical moment
考慮到不同潮汐類型對河流水質的影響,本研究分析了不同潮汐類型下的潮流場分布,由圖5可知,日潮高潮時,表層海水入侵到虎爪斷面下游約400 m處,底層海水入侵到虎爪斷面上游約1 100 m處;半日潮高潮時,表層海水基本停留灣內,底層海水入侵到虎爪斷面上游約700 m處。總體上表現為日潮和半日潮的高潮時刻底層海水可入侵到虎爪斷面上游700~1 100 m處,但表層海水并不能到達虎爪斷面上游。

圖5 日潮和半日潮流域表、底層潮流場分布Fig.5 Distribution of tidal current at diurnal and semidiurnal tides on the surface and bottom
通常河流淡水的鹽度在0.5‰以下,河水與陸架海水的混合水體鹽度為0.5‰~30‰,陸架海水的鹽度為30‰以上,因為底層鹽度梯度更大,界面混合影響小,可采用底層鹽度等值線作為河海水團分界線。模擬得到淡澳河上游姚田橋至虎爪斷面水體的底層鹽度分布如圖6所示。由圖6可知,0.5‰鹽度線分布在虎爪斷面上游,距離姚田橋約1 490~1 520 m處,在河口海岸水環境整治時可依據分界線對淡水和海水進行分區域治理。此外,由于鹽度分布受徑潮流等影響,季節差異大,分界線較為復雜,但水文要素往往有周期性規律,研究基于2019年水文資料初步劃分了河海水團分界線,考慮到分界線的穩定性和合理性可考慮取多年平均值來劃分。

圖6 淡澳河上游姚田橋至虎爪斷面底層水體鹽度分布Fig.6 Salinity distribution of bottom water from Yaotianqiao to Huzhao section in the upstream of Dan'ao River
淡澳河上、下游邊界處不同水層的含鹽量基本一致,其中上游邊界處含鹽量約為0‰,下游邊界處含鹽量約為35‰。中游地區由于海水入侵呈現明顯的鹽度分層流,研究以虎爪斷面為代表斷面,通過分析不同水期下水體表層和底層鹽度分布(圖7),確定水體淡鹽水分配比例。

圖7 中游河段不同水期水體鹽度分布Fig.7 Distribution of water salinity in the midstream in different water periods
枯水期時,表層和底層含鹽量約為10.4‰和20.7‰,經計算得到表層和底層咸淡水占比分別約為1∶2.5和10∶7;平水期時,表層和底層含鹽量約為1.7‰和10.8‰,計算得到表層和底層咸淡水占比分別約為1∶16和1∶2;豐水期時,表層和底層含鹽量約為0.3‰和5.0‰,計算得到表層和底層咸淡水占比分別約為1∶100和1∶6。由咸淡水占比可見,豐水期虎爪斷面水環境質量主要由淡澳河上游來水主導,枯水期時,虎爪斷面水質情況由淡澳河上游來水和海水共同作用,其中表層受上游來水影響大,底層受海水入侵影響大。
最新研究表明[21],大亞灣海域營養鹽絕大部分通過河流輸入,其中氮通量輸入9.4×105kg/a,磷通量輸入1.7×105kg/a,淡澳河是大亞灣入海污染物輸入的主要途徑,淡澳河輸入氮通量占76.7%,磷占82.4%。為進一步證明上游輸入對淡澳河口產生了明顯影響[22],摸清上游來水流量與河流水質的響應關系,以2019年8月實測水文水質數據為基準,即上游來水流量為6.0 m3/s,同時依據實際水文情況設定枯水期小流量為1.0 m3/s和豐水期大流量為10.0 m3/s。
由模擬結果可知(圖8),枯水期潮期間污染物濃度變化不明顯,虎爪斷面COD降幅為39%~42%,氨氮濃度降幅為27%~30%,TP濃度降幅為13%~15%;豐水期小潮期間污染物濃度增加幅度最為明顯,虎爪斷面COD增幅為1%~12%,氨氮濃度增幅為11%~14%,TP濃度增幅為6%~7%。由以上數據可知,上游來水流量變化直接影響虎爪斷面水質,流量越大對虎爪斷面的水質影響越大,以小潮期間最為顯著,建議進行相關工程優化調度,以減小小潮期間的下泄流量。

圖8 上游來水變化對虎爪斷面水質影響分析Fig.8 Analysis of the influence of upstream water on the water quality of Huzhao section
(1)大亞灣潮流性質以不規則半日混合潮型為主,每月約8~10 d為日潮,20~22 d為半日潮,落潮流速高于漲潮流速,日潮漲急流速略高于半日潮漲急流速。由不同潮汐類型對淡澳河潮流場的影響可知,日潮、半日潮的高潮時表層海水并不能到達虎爪斷面,但底層海水可入侵到虎爪斷面上游700~1 100 m處。
(2)基于2019年水文資料初步劃分河海水團交界線,位于虎爪斷面上游距姚田橋約1 490~1 520 m處,可為河口海岸帶淡水和海水分類整治提供依據。鹽水入侵表明,淡澳河水體含鹽量枯水期>平水期>豐水期;以虎爪斷面為例,枯水期含鹽量高達16‰,底層海水含量高于淡水,枯水期表層海水含鹽量較平水期增加明顯,相當于平水期的6倍。
(3)上游來水流量與水質響應模擬研究表明,來水流量直接影響虎爪斷面水質,其中,枯水期污染物濃度降低15%~40%,而豐水期污染物濃度增加10%左右,以小潮期間增加幅度最為顯著。基于此,可建議進行相關工程優化調度,以減小小潮期間的下泄流量。