山丘區(qū)河流灘地邊緣內(nèi)凹形態(tài)對溶質(zhì)運移及水質(zhì)惡化風(fēng)險的影響

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

(河海大學(xué)農(nóng)業(yè)科學(xué)與工程學(xué)院，南京 210098)

0 引 言

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

1 材料與方法

1.1 研究區(qū)概況

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

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

1.2 數(shù)據(jù)采集

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

1.3 二維水動力-水質(zhì)模型構(gòu)建

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

連續(xù)性方程：

(1)

動量方程：

(2)

(3)

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

(4)

式中：c為標(biāo)量濃度；Dx和Dy分別為x和y方向上的擴散系數(shù)。

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

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

1.4 灘地邊緣內(nèi)凹形態(tài)定量表征

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

1.5 水質(zhì)惡化風(fēng)險的表征方法

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

2 結(jié)果與分析

2.1 灘地邊緣內(nèi)凹形態(tài)分布特征

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

表1 靈山港灘地邊緣內(nèi)凹形態(tài)特征參數(shù)Tab.1 Floodplain edge morphological characteristics in the downstream of Lingshan River

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

2.2 模型校核及可靠性分析

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

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

2.3 內(nèi)凹形態(tài)對溶質(zhì)濃度空間分布的影響

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

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

2.4 內(nèi)凹形態(tài)對溶質(zhì)濃度運移時間的影響

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

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

2.5 內(nèi)凹形態(tài)對近灘水域水質(zhì)惡化風(fēng)險影響

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

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

3 討 論

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

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

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

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

4 結(jié) 論

論文在野外調(diào)查的基礎(chǔ)上，建立了灘地邊緣內(nèi)凹形態(tài)指數(shù)計算方法及等級劃分方法，運用數(shù)值模擬探究不同彎曲等級對溶質(zhì)運移及水質(zhì)惡化的影響，為灘地治理和水環(huán)境改善提供參考。

(1)山丘區(qū)灘地邊緣內(nèi)凹形態(tài)可根據(jù)邊緣彎曲度指數(shù)Ie劃分為低彎(Ie< 2)、中彎(2 ≤Ie< 4)和高彎(Ie≥ 4)。

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

(3)當(dāng)Ie為低彎時，內(nèi)凹形態(tài)所圍近灘水域幾乎不存在水質(zhì)惡化風(fēng)險；當(dāng)Ie為中彎時，近灘水域水質(zhì)惡化風(fēng)險較低；當(dāng)Ie為高彎時，水質(zhì)惡化風(fēng)險則較高。但是根據(jù)維持生境多樣性的需求，不宜將內(nèi)凹形態(tài)降至低彎。

□