山丘區河流灘地邊緣內凹形態對溶質運移及水質惡化風險的影響

周之悅，夏繼紅，蔡旺煒，葉繼兵，劉 瀚，竇傳彬

(河海大學農業科學與工程學院，南京 210098)

0 引 言

灘地作為水域與陸域的過渡區域，是河道中水文、環境和生態過程相互作用的敏感區，對洪水調蓄、水質凈化、生物棲息等具有重要價值[1,2]。我國廣大的山丘區中分布著大量的中小型河流，其中不乏形態多樣的灘地。然而人類活動、異常水文過程、地質活動等因素易造成灘地形態及格局的破壞，使得灘地功能退化，影響河流生態系統健康[3,4]。Mejia等[5]通過開發基于瞬時響應函數和流域地貌關系的分析框架，探討了灘地格局塑造的河道斷面形態對流域水文過程的影響；蔡建楠[6]等通過構建含河道斷面形態、河岸帶寬度等指標在內的城市河流形態評價體系，發現河流形態評分較高的河段水體自凈能力也越強。伊紫函等[7]對灘地平面幾何形態進行了分類以研究灘地整體形態演變；王金平等[8]應用灰色關聯方法分析灘地景觀格局時空演化特征，發現近年來灘地破碎化特征明顯且常出現內凹向灘地腹地的邊緣形態。然而目前關注灘地形態的研究多涉及流域尺度上的地貌格局探討，對中、小尺度相關問題的定量研究較少[9]，尤其缺乏對灘地邊緣形態的研究，如，灘地邊緣內凹形態變化對水體溶質運移過程及水質惡化風險影響機理等問題尚待深入研究。本文以龍游縣靈山港為典型，通過構建定量表征灘地邊緣內凹的邊緣形態指數(Ie)，探討灘地內凹程度等級劃分方法，深入研究不同內凹程度對溶質運移的影響及水質惡化的風險，以期為山丘區河流灘地修復提供數據支撐和理論依據。

1 材料與方法

1.1 研究區概況

研究區域(E119°07′～119°12′，N28°49′～29°03′)為浙江省西部的靈山港，是錢塘江上游衢江右岸的一條重要支流。流域總面積 726. 9 km2，屬亞熱帶季風氣候區，年平均降水量為1 815 mm，年內分布不均，主要集中在3-9月份。多年平均徑流量20.8 m3/s，徑流系數為0.62。流域內山地土壤以黃紅壤亞類、侵蝕性紅壤、黃壤為主，母質為凝灰巖、片麻巖、板巖等。河床及灘地以堆積和侵蝕作用為主，堆積物一般厚為3～8 m。靈山港主流長43.79 km，河道平均比降為0.245%，沿線設有步坑口、龍游水文站。近40年來，靈山港內的灘地呈現逐漸破碎和萎縮的變化特征，下游河段(姜席堰-蘭石大橋 段)最為明顯[8]。該河段(圖1)總長6.1 km，平均比降為0.181%，河道寬度為50～160 m不等，共有18處灘地，其分布格局復雜多樣，流速分布空間差異較大。

圖1 研究區灘地分布格局及計算網格剖分示意圖Fig.1 The patterns of floodplains in the study area and sketch map of computational mesh

1.2 數據采集

分別于2017年4月和10月、2018年4月和10月對靈山港下游段開展了實地調查，采集河道地形、灘地邊緣形態、流量與水位、溶質(總磷TP)濃度等基礎數據。方法：①河道地形數據由當地水利部門提供；②灘地邊緣形態參數采用手持GPS記錄儀(MAGELLAN，美國)實地踏勘而得；③流量與水位數據主要源于水文監測站，另于校核斷面(圖1(a))處設置了臨時監測站點，記錄逐日水位；④水質指標TP通過現場采集水樣帶回實驗室參照《流動注射-鉬酸銨分光光度法 HJ 671-2013》測定。

1.3 二維水動力-水質模型構建

基于Mike 21構建靈山港下游段水動力-水質耦合模型。水動力學控制方程為基于Boussinesq假定和靜水壓假設的二維不可壓雷諾平均N-S方程組(即淺水方程)：

連續性方程：

(1)

動量方程：

(2)

(3)

基于對流-擴散方程耦合水動力模塊模擬溶質運移過程，水質模塊控制方程如下：

(4)

式中：c為標量濃度；Dx和Dy分別為x和y方向上的擴散系數。

通過對計算區域進行網格劃分，采用基于單元中心的有限體積法對物理控制方程進行數值求解。初始條件、邊界條件等參數通過具體模塊設置輸入，其中水動力模塊為基礎核心模塊。由于山丘區河流河床地貌形態變化復雜，水流方向多變，灘地的出露形態隨河流水位的動態變化而改變。遂依據2016年河道及灘地實測地形圖，選擇三角形非結構化網格對計算區域進行剖分，共生成網格330 377個，允許的最大網格面積為5 m2，網格最小角度為32°。網格分布特征如圖1(b)所示。

進行模擬計算時以常水位流量(Q=15 m3/s)為上游邊界條件，以相應水位(h=47 m)為下游邊界條件。因研究河段較短，不考慮降雨、蒸發及物質降解。參照實際情況，設置河道全水域TP背景值為0.023 8 mg/L。分析模擬結果時，當TP濃度<0.040 mg/L即認為恢復至背景值水平。待水動力條件穩定后，上游邊界持續2 h輸入TP濃度為8 mg/L。以溶質投放時刻作為模擬計算的初始時刻，計算步長為30 s，共計3 600步(30 h)。

1.4 灘地邊緣內凹形態定量表征

借鑒河道蜿蜒度的概念[10]，本文將灘地邊緣曲線長度(S)與近灘水域開口寬度(L，與主河道水流方向平行)的比值定義為彎曲度指數Ie，即Ie=S/L，用以定量表征灘地邊緣內凹形態。Ie值越大表明灘地邊緣越彎曲，內凹形態向灘地腹地延伸的程度越高。

1.5 水質惡化風險的表征方法

溶質運移是指溶解在水中的物質輸運和轉移的空間和時間模式[11]。當某處水域的營養鹽或污染物的濃度在較長時間內偏高時，則該水域的水質惡化風險較高。在本研究中，①采用近灘水域TP與主河道TP均值的比值(即TP的相對濃度R)來反映溶質在近灘水域的聚積程度，R越高水質惡化風險越高；②采用近灘水域溶質濃度等值線圖反映TP濃度的空間梯度，等值線越密集表明水質惡化風險的空間變異性越劇烈，而等值線圖分布的空間范圍越大則說明內凹形態影響的水域面積越大；③根據近灘水域形心所在斷面的溶質平均濃度的時間變化曲線，及該濃度與相應主河道斷面的溶質平均濃度差值ΔC的時間變化曲線，來反映溶質在近灘水域中的累積與削減過程。

2 結果與分析

2.1 灘地邊緣內凹形態分布特征

野外調查結果顯示，共有14處灘地邊緣內凹形態(編號為F1～F14)分布于靈山港下游段的11個灘地上。由表1可知，這些內凹形態的邊緣曲線長度S、近灘水域開口寬度L(平行于河道主流方向)、近灘水域面積A等特征參數的分布范圍分別是43.1～191.0 m、16.8～59.6 m和206～1 460 m2，內凹形態邊緣彎曲度Ie的分布范圍則為1.51～6.70。根據下游段灘地Ie的分布頻率，將灘地邊緣內凹形態劃分為低彎(Ie<2)、中彎(2 ≤Ie< 4)和高彎(Ie≥ 4)。則F1～F6為低彎，F7～F9為中彎，F10～F14為高彎。

表1 靈山港灘地邊緣內凹形態特征參數Tab.1 Floodplain edge morphological characteristics in the downstream of Lingshan River

注：S為灘地內凹形態邊緣曲線長度；L為內凹形態開口寬度；A為內凹形態所圍近灘水域面積。

2.2 模型校核及可靠性分析

基于2016年的流量監測數據進行模型參數率定，確定主河道曼寧系數為28～30，灘地曼寧系數為17～25。基于水文監測數據得校核斷面2017年逐日水位的模擬與實測值吻合良好，Nash-Sutcliffe效率系數(NSE)為0.98，接近于1，模型質量較好，總體相對誤差(RE)為6.29%(見圖2)。表明所構建的水動力數值模型具有良好的穩定性及可靠性，適用于模擬溶質運移過程。

圖2 校核斷面2017年水位驗證曲線Fig. 2 Validation of daily surface water level in 2017

2.3 內凹形態對溶質濃度空間分布的影響

選取位于研究區中部的內凹形態F11作為典型分析灘地內凹形態對近灘水域溶質濃度空間分布的影響。通過調整局部邊緣形態使F11的Ie值從4.86(高彎)逐漸下降為2.43(中彎)和1.82(低彎)，形態特征如圖3所示。由近灘水域溶質濃度等值線圖(圖3)可知，在近灘水域中，TP濃度與主河道存在明顯差異的水域范圍隨著Ie的降低而減小，同時等值線的數量也隨之減少。這一結果表明，灘地邊緣內凹形態的彎曲度等級越高，其近灘水域中存在水質惡化風險的水域面積將會越大。

圖3 不同邊緣形態下近灘水域溶質濃度等值線圖Fig.3 Solute concentration map of nearshore water

2.4 內凹形態對溶質濃度運移時間的影響

模型輸出不同形態條件下近灘水域溶質濃度及其與主河道濃度差值的時間變化曲線。由圖4(a)知，當Ie等級從高彎(M0)變為低彎(M2)后，近灘水域溶質濃度從峰值下降至背景水平(<0.040 mg/L)所需的時間從超過22 h縮短至11.8 h。而由圖4(b)知，溶質濃度在近灘水域高于主河道所持續的時間從超過23 h(M0-高彎)下降至9.8 h(M2-低彎)，且勢差峰值明顯降低。上述結果表明，灘地邊緣形態的變化對近灘水域溶質變化歷程具有顯著的影響。表現為隨著Ie下降，近灘水域溶質濃度的變化過程與主河道的差異縮小。

圖4 近灘水域溶質濃度變化曲線及 近灘水域-主河道濃度勢差變化曲線Fig.4 The temporal curve of solute concentration in nearshore water and the concentration difference between nearshore water and the main channel

2.5 內凹形態對近灘水域水質惡化風險影響

野外調查結果顯示14處灘地邊緣內凹形態的近灘水域TP濃度(0.152～1.825 mg/L)與對應主河道斷面的TP濃度(0.023～0.190 mg/L)差異較大，TP相對濃度R最大可達10倍。由圖5可知，模擬結果(當主河道TP濃度恢復至背景水平時)與實測結果均顯示，隨著Ie從低彎升至高彎，R呈極顯著(P<0.01，DUNCAN)升高。這一結果表明隨著灘地邊緣內凹形態的彎曲度等級升高，近灘水域水質惡化風險也隨之增高。另外，相同彎曲度等級下，R的實測值均略低于模擬值。如高彎等級下，R的模擬值和實測值的均值分別為7.01±1.18和5.53±2.71。產生這一差別的主要原因是模擬計算時不考慮溶質的降解和植物吸收。

圖5 不同彎曲度等級下近灘水域TP相對濃度RFig.5 Relative TP concentration under different Ie levels

3 討 論

上述結果顯示灘地邊緣內凹形態的彎曲度指數Ie越大則近灘水域內溶質的相對濃度越高，溶質濃度下降速度越慢。根據余根聽等[12]所開展的室內物理模型試驗結果，正弦型蜿蜒岸的背水坡振幅越大(等價于Ie越大)，溶質在潛流帶中的停留時間越長。因而，在近灘水域內，地表水與潛流可能會對溶質運移過程發生協同作用，進一步延長溶質在水體中的滯留時間，從而引發水質惡化風險，影響灘地生態系統的穩定性與功能性。

溶質動力學與水流的物理運動緊密相關。對于灘地錯落分布的山丘區河流來說，對流作用對溶質運移的驅動影響要大大強于溶質勢差主導的擴散作用。在暫態存儲理論中，溶質被認為能夠暫時存儲在河流中的緩流或靜水區，進而能夠延長溶質向下游水體的釋放時間[13]。而增強主流區與溶質暫態存儲區之間的水流交換作用對河流系統健康至關重要[14]。從這個意義上說，對F11的形態改造本質上是通過控制邊緣形態來控制主流區與溶質暫態存儲區的對流關系，進而影響了暫態存儲區對溶質的存儲能力。

自然河流系統中，灘地周邊廣泛存在的緩水區、回水區及靜水區能夠為水生生物提供多樣化的生境條件，尤其為魚類在山丘區河流中提供了重要的休憩、產卵及索餌的場所[15]。在對灘地邊緣形態進行優化改造時，一方面應當避免營養物質或污染物等溶質在近灘水域中停留時間過長而發生水質災害，另一方面還應滿足生境多樣化的需求。根據本文研究結果，Ie為中彎(2 ≤Ie< 4)時，近灘水域不會出現過于明顯的溶質積聚現象，水質惡化風險不高。實際上中彎條件還能提供一定面積比例的緩流區，能為水生生物提供良好的生境。對于高彎形態，通過改造局部形態就有可能使得其彎曲度從“高彎”下降為“中彎”。這類基于Ie指數進行設計的“微改造”既能避免對原有灘地整體形態及生態功能的強烈擾動，又能減少工程量。

此外，即便Ie處于相同等級，其近灘水域面積的大小可能差異較大(如F13和F14，見表1)。在實際情況中，當近灘水域的面積較大時，其TP相對濃度較高，反之，TP相對濃度較低(如F13的TP相對濃度R的實測值僅為1.84)。在后續研究中可將近灘水域面積大小與彎曲度指數Ie結合起來對灘地邊緣內凹形態進行更精細的劃分。

4 結 論

論文在野外調查的基礎上，建立了灘地邊緣內凹形態指數計算方法及等級劃分方法，運用數值模擬探究不同彎曲等級對溶質運移及水質惡化的影響，為灘地治理和水環境改善提供參考。

(1)山丘區灘地邊緣內凹形態可根據邊緣彎曲度指數Ie劃分為低彎(Ie< 2)、中彎(2 ≤Ie< 4)和高彎(Ie≥ 4)。

(2)灘地邊緣內凹形態的變化能顯著影響回水區的溶質分布和運移過程。Ie的變化對近灘水域溶質濃度的空間分布與時間變化過程均具有顯著的影響。隨著Ie降低，近灘水域溶質濃度的聚積程度、空間異質性、和削減耗時均明顯降低。

(3)當Ie為低彎時，內凹形態所圍近灘水域幾乎不存在水質惡化風險；當Ie為中彎時，近灘水域水質惡化風險較低；當Ie為高彎時，水質惡化風險則較高。但是根據維持生境多樣性的需求，不宜將內凹形態降至低彎。

□