李 博 楊 昊 歐素英 蔡華陽①劉 鋒 楊清書
(1. 中山大學(xué)海洋工程與技術(shù)學(xué)院 河口海岸研究所 廣東廣州 510275; 2. 河口水利技術(shù)國家地方聯(lián)合工程實驗室 廣東廣州510275; 3. 廣東省海岸與島礁工程技術(shù)研究中心 廣東廣州 510275; 4. 南方海洋科學(xué)與工程廣東省實驗室(珠海) 廣東珠海519000)
河口是陸海徑潮動力耦合的獨特區(qū)域, 既有上游徑流匯入, 又有外海潮波上溯, 還受人類活動疊加效應(yīng)的影響, 導(dǎo)致其動力結(jié)構(gòu)復(fù)雜多變。河口動力的典型特征包括徑潮相互作用、斜壓效應(yīng)、波流耦合等,其中, 徑潮相互作用對洪水下泄、鹽淡水混合、最大渾濁帶形成、營養(yǎng)鹽和污染物輸移等物理、化學(xué)和生物過程均有直接影響。潮波振幅梯度, 即潮波振幅沿程的平均變化率, 綜合反映了徑流、潮汐和底床摩擦等多因子影響在河口空間上的差異, 是探討河口復(fù)雜徑潮動力結(jié)構(gòu)的有效切入點。因此, 揭示潮波振幅梯度與上下游動力邊界(即上游徑流和外海潮汐)的耦合過程及機制是河口動力學(xué)研究的基礎(chǔ)科學(xué)前沿問題, 可為河口海岸區(qū)域自然災(zāi)害(如咸潮上溯、洪澇、赤潮等)的防治和水資源的高效開發(fā)利用等提供重要科學(xué)依據(jù)(Schuttelaarset al, 2013; Caoet al, 2020; 楊昊等, 2020)。
國內(nèi)外針對河口徑潮相互作用的研究由來已久。研究表明, 當潮波受到徑流和地形的調(diào)制作用時呈現(xiàn)出明顯的非線性變化(Hoitinket al, 2017; Jiet al,2019), 河口潮波衰減主要通過流量和底床的摩擦效應(yīng)實現(xiàn)(Godin, 1985; Savenije, 1992; Horrevoetset al,2004)。當流量增大時, 全日分潮和半日分潮的振幅以不同幅度衰減, 而河道底床下切, 水深增大可導(dǎo)致底床摩擦效應(yīng)減弱, 引起沿程潮波振幅不同程度的增加(Caiet al, 2012)。對河口地形進行概化并考慮下泄徑流影響, 采用一維潮波傳播解析模型可通過潮波振幅梯度和余水位梯度等參數(shù)揭示潮波傳播的時空變化特征及其動力學(xué)機制(Caiet al, 2016)。當河口動力以徑流作用為主時, 傳統(tǒng)調(diào)和分析(Pawlowiczet al,2002)的水位預(yù)報誤差明顯增大。針對流量影響下河口潮波傳播的非定常和非線性問題, 對水位序列采用連續(xù)小波變換(Guoet al, 2015; Moftakhariet al,2016)、經(jīng)驗?zāi)B(tài)分解(Panet al, 2018a)、經(jīng)驗正交函數(shù)(Panet al, 2019)、非平穩(wěn)潮汐調(diào)和分析(Matteet al,2013, 2014; Panet al, 2018b)等處理方法探討流量對潮波傳播的影響過程及機制均取得較好效果。
近年來, 流量對潮波傳播的衰減及其閾值效應(yīng)研究得到廣泛關(guān)注。研究表明, 在長江河口感潮河段,流量與潮波振幅梯度之間存在閾值效應(yīng), 即流量小于臨界閾值時, 隨著流量增大潮波衰減效應(yīng)增強, 但當潮波振幅梯度達到極小值后, 隨著流量增大衰減效應(yīng)反而減弱(Caiet al, 2019; 張先毅等, 2019)。在珠江磨刀門河口也觀察到類似現(xiàn)象, 即在甘竹-馬口河段潮波振幅梯度隨流量變化有先減小后上升的特點,在人類活動干預(yù)較弱的自然演變階段, 潮波振幅梯度在流量為5 000~10 000 m3/s 時達到極小值, 而在人類活動強烈影響后的恢復(fù)調(diào)整階段, 該流量閾值上升至15 000~20 000 m3/s, 且其對應(yīng)的潮波振幅梯度極小值增大(Yanget al, 2020; 張先毅等, 2020), 反映因采砂導(dǎo)致的河槽容積增大對潮波振幅梯度閾值的增大效應(yīng)。綜上所述, 對于河口徑潮相互作用及潮波振幅梯度與流量動力邊界關(guān)系的研究已有較多成果,但是對于定量辨識流量驅(qū)動下潮汐調(diào)和分析結(jié)果的階段性演變以及潮波振幅梯度與下游動力邊界關(guān)系的研究等相關(guān)科學(xué)問題還尚待深入。
位于我國粵港澳大灣區(qū)國家戰(zhàn)略核心經(jīng)濟圈的珠江三角洲是全球受到人類活動影響程度最大的區(qū)域之一。高強度的人類活動(如水庫建設(shè)、口門圍墾、無序采砂、航道疏浚等)導(dǎo)致來水來沙和地形地貌發(fā)生明顯變化, 進而顯著改變珠江河網(wǎng)的徑潮動力格局。水庫建設(shè)使河流挾沙量銳減, 三角洲河道下切加劇; 地形變化導(dǎo)致西江、北江的水沙分配格局發(fā)生變化, 河網(wǎng)中上部地區(qū)余水位梯度變緩(Liuet al, 2019);口門圍墾、口門整治疊加海平面上升的影響則導(dǎo)致磨刀門河口的形態(tài)向窄深化趨勢演變(Tanet al, 2019)。以上變化導(dǎo)致珠江主要出??? 即磨刀門河口的潮流界和潮區(qū)界顯著向上游移動, 咸潮上溯災(zāi)害風(fēng)險增大, 對粵港澳大灣區(qū)城市的居民生活、工農(nóng)業(yè)生產(chǎn)用水的水質(zhì)安全形成嚴峻挑戰(zhàn)。因此, 本文針對強人類活動驅(qū)動下的徑潮動力異變問題, 采用數(shù)據(jù)驅(qū)動模型探究流量驅(qū)動下珠江磨刀門河口潮波傳播特征的異變過程, 并分析潮波調(diào)和常數(shù)的階段性演變, 揭示潮波振幅梯度對上下游動力邊界的響應(yīng)過程及機制, 研究成果可為受到高強度人類活動影響的河口三角洲的治理和保護等提供科學(xué)依據(jù)。
珠江河網(wǎng)(圖1)流域面積達4.5×105km2, 具有三江(東江、西江、北江)匯流, 八口(自東向西依次為虎門、蕉門、洪奇門、橫門、磨刀門、雞啼門、虎跳門、崖門)入海的特點, 河網(wǎng)縱橫交錯導(dǎo)致徑潮動力相互作用過程復(fù)雜。珠江每年約有2 823×108m3淡水流量及72.4×106t 懸沙量流入南海(Liuet al, 2018), 平均流量約為8 952 m3/s, 口門附近潮差為1.0~1.7 m。磨刀門河口為西江入??? 是八大口門中輸水輸沙量最大的河口, 該河口屬典型的河控型河口。據(jù)三角洲頂端馬口水文站1959~2016 年月均流量資料, 其多年平均流量為7 078 m3/s, 月均流量最大值為29 000 m3/s(6 月), 最小值為1 210 m3/s (1 月), 洪枯季差異明顯(楊昊等, 2020)。潮汐以不規(guī)則半日潮為主, 日不等現(xiàn)象顯著, 燈籠山站平均漲落潮歷時分別為5.37 h 和7.25 h, 口門處三灶站潮型系數(shù)為1.30。
圖1 研究區(qū)域及所用潮位站或水文站位置Fig.1 The study area and the deployment of the tidal gauging stations and the hydrological station
隨著經(jīng)濟社會高速發(fā)展, 強人類活動已成為河口三角洲系統(tǒng)演變的第三驅(qū)動力(陳吉余等, 2008)。口門圍墾、聯(lián)圍筑閘、河道采砂及上游流域水庫建設(shè)等高強度人類活動極大影響磨刀門河口的河道形態(tài)和徑潮動力, 其“動力-沉積-地貌”體系發(fā)生顯著變化(吳超羽, 2018)。20 世紀80 年代前, 磨刀門地區(qū)主要人類活動為大規(guī)模灘涂圍墾, 至 1 9 9 7 年已有16.5 km2的淺灘被圍填成陸地(Tanet al, 2019)。圍墾、疏浚、口門整治等一系列工程措施使河口口門區(qū)水域面積減小, 入海水道水深增大, 至20 世紀90 年代初入??陂T向海延伸約13 km (何用等, 2018)。其次, 20世紀50~70 年代中期, 大規(guī)模的聯(lián)圍筑閘導(dǎo)致三角洲腹部水位壅高, 中上游水面比降變緩, 流速減小, 水流挾沙能力降低, 河床產(chǎn)生淤積(劉鋒等, 2011)。20世紀80 年代中期開始, 由于用沙需求激增, 在珠江河網(wǎng)區(qū)出現(xiàn)大規(guī)模河道采砂活動。磨刀門采砂活動主要集中在珠海大橋以北河段, 1991~2000 年共挖沙420×104m3(韓志遠等, 2010)。截至2006 年, 西江流域大壩總庫容達5.38×1010m3, 占其入海淡水流量的24.7% (Liuet al, 2018)。大壩建成后流域內(nèi)大部分沉積物被截留在上游地區(qū), 西江干支流輸沙量均呈下降趨勢, 下游泥沙供給量銳減(Wanget al, 2011)。在上述強人類活動的綜合影響下, 珠江磨刀門河口沿程河床大幅下切, 河槽容積增大, 導(dǎo)致近年來潮汐動力顯著增強, 咸潮上溯加劇(韓志遠等, 2010; 張先毅等, 2020)。
本文所用數(shù)據(jù)取自《廣東省水文年鑒》第八卷和廣東省水文局, 潮位原始數(shù)據(jù)高程為凍結(jié)基面, 已統(tǒng)一轉(zhuǎn)換至珠江基面。數(shù)據(jù)包括磨刀門河口沿程4 個潮位站(三灶、燈籠山、竹銀、甘竹站)逐日高、低潮位數(shù)據(jù)和相應(yīng)時段馬口水文站日均流量及逐日高、低潮位數(shù)據(jù), 并對潮位和流量數(shù)據(jù)進行分段三次Hermite插值處理后獲得的逐時數(shù)據(jù)。站點具體信息見表1。
表1 磨刀門河口各站點數(shù)據(jù)信息Tab.1 Information of the adopted tidal gauging stations and hydrological station in the Modaomen estuary
本文采用徑潮耦合的數(shù)據(jù)驅(qū)動模型(river-tidal harmonic analysis, R_TIDE)對逐時潮位數(shù)據(jù)進行調(diào)和分析(歐素英等, 2017), 模型可輸出重構(gòu)的逐時水位、特定分潮的振幅和遲角等。在徑流動力占主導(dǎo)的河口,潮波上溯受到徑流季節(jié)性周期變化的調(diào)制影響, 傳統(tǒng)的潮汐調(diào)和分析模型對河口區(qū)潮波特征的分析和預(yù)報誤差較大, 不能滿足徑潮信號分離的要求。為有效辨識潮波在流量影響下的傳播特征, 基于Matte 等(2013, 2014)提出的非平穩(wěn)調(diào)和分析(nonstationary tidal harmonic analysis, NS_TIDE)思路, 假設(shè)河口任意位置潮波振幅和遲角主要受上游流量和地形變化的非線性調(diào)制影響, 在此基礎(chǔ)上進一步假定分析時段內(nèi)河口地形邊界不變, 將信號分成由流量Q引起的水位變化和徑流調(diào)制影響下的潮水位變化, 通過信號分離定量分析感潮河段徑潮動力相互作用。
無徑流調(diào)制時, 對于月球和太陽引起的周期性潮汐現(xiàn)象可看作是許多假想天體引起的簡單的潮汐簡諧波動的總和(Pawlowiczet al, 2002), 即潮位曲線可近似看作是由許多振幅、周期、相位不同的分潮組成。潮位可表達成以下多個余弦函數(shù)的疊加形式:
其中,z(t)為站點的實測潮位;t為時間;z0為平均海面高度;N為分潮個數(shù);fi為分潮振幅的訂正因子, 為時間的函數(shù), 常取一年的中值;Hi為分潮平均振幅;σi為分潮角頻率;Vi為分潮初相角;ui為天文相角的交角訂正角,gi為由于海底摩擦、海水慣性引起的遲角。其中Hi和gi合稱為分潮調(diào)和常數(shù)。
Kukulka 等(2003)基于一維圣維南方程組推導(dǎo)得出河口三角洲內(nèi)任意位置x的潮波振幅變化:
其中,j=1,2,3,…,m,m為逐時序列的數(shù)據(jù)個數(shù);r為衰減系數(shù);ε0(tj)為河口口門位置x=0 處的潮波振幅;a'為常數(shù)。衰減系數(shù)r為潮波傳播速度、水深、流量和河道摩擦系數(shù)的函數(shù)(Jayet al, 1997)。若摩擦系數(shù)采用曼寧公式, 且摩擦系數(shù)與潮波傳播速度均為水深的函數(shù), 則衰減系數(shù)r可簡寫為水深h(tj)和流量Q(tj)的函數(shù):
或:
其中,h(tj)和Q(tj)表示隨時間變化的水深和流量;p0、p1、p2、β、γ均為待定參數(shù)。式(3)未考慮外海非潮動力邊界比如風(fēng)暴等因素的影響, 這是與Matte 等(2013,2014)提出的非定常調(diào)和分析模型的不同之處。假定研究時間段河道底床高程不隨時間變化, 則式(4)中水深項可歸到常數(shù)項, 衰減系數(shù)r可簡寫為
即河口潮汐非線性變化僅受流量驅(qū)動影響。結(jié)合式(1)、(2)、(5)可得流量影響下任意位置x的潮水位變化:
式(6)的待定系數(shù)bi、ci、γ可通過粒子群優(yōu)化算法(Kennedyet al, 1995)或求解帶約束的非線性多變量方程組來確定該參數(shù)的最優(yōu)值。由于流量對不同分潮的影響程度存在差異(Guoet al, 2015), 半日分潮的衰減效應(yīng)大于全日分潮, 將流量對每一個分潮族的調(diào)制影響均采用不同的待定系數(shù), 則式(6)可表示為:
其中,
式(7)即為流量驅(qū)動下的R_TIDE 模型,S(tj)為底床高程或水深變化、海平面、流量等引起的平均水面變化, 稱為余水位項;F(tj)代表不同流量條件下徑潮耦合引起的水位變化, 稱為徑潮耦合項。d0、d1、v0、v1、r0、r1為模型的待定參數(shù), 通過求解帶約束的非線性多變量方程組來確定。同時, 由于磨刀門河口徑流動力占優(yōu)勢導(dǎo)致徑潮動力非線性因素作用較強, 洪水流量時導(dǎo)致三角洲上段為徑流控制, 無潮汐波動。因此, 在R_TIDE 模型中引入臨界流量QC, 當連續(xù)2 d 潮差小于某個值(如0.001 m)時所對應(yīng)的最小流量即為臨界流量, 當流量大于QC時, 表示該站點及其以上河段不存在潮汐信號。
假設(shè)河口三角洲任意位置采用時間間隔為Δt的觀測潮位y(tj)(j=1,2,3,…,m)序列, 根據(jù)實測資料序列長度mΔt、分潮頻率差和 Rayleigh 準則判據(jù) Δσ=max[(mΔt)-1,σi-σj], 選擇待提取的分潮, 將參數(shù)tj、σj、Q(tj)代入式(9), 要求模型計算結(jié)果y(tj)和實測水位z(tj)的誤差平方和達到最小, 即:
選擇對應(yīng)資料長度的全日分潮、半日分潮等分潮族, 對實測潮位進行回歸擬合, 回歸模型效果以均方根誤差ERMS和相關(guān)指數(shù)R2表示, 即:
前期研究表明, 基于燈籠山-馬口河段月均余水位梯度(或潮波振幅梯度絕對值)與流量的雙累積曲線變化趨勢, 可知曲線斜率在1990 年與2000 年發(fā)生明顯改變(張先毅等, 2020)。根據(jù)上述結(jié)果, 可將磨刀門河口徑潮動力格局的演變分為三個階段: 1990 年前為自然演變和人類活動共同作用的階段, 但此階段人類活動尚未對河口徑潮動力格局產(chǎn)生明顯影響, 河口處在動態(tài)平衡期, 以自然演變?yōu)橹? 稱為“平衡期”(簡稱P1); 2000 年后, 河道中下游采砂活動基本已停止, 灘涂圍墾面積銳減, 河口在人類活動影響下已發(fā)生徑潮動力格局轉(zhuǎn)換, 為強人類活動干預(yù)后河口的自適應(yīng)調(diào)整階段, 稱為“調(diào)整期”(簡稱P2); 1991~2000年為人類活動影響最強烈的時期, 河口從自然演變?yōu)橹鞯钠胶鈶B(tài)向強人類活動干預(yù)后的平衡態(tài)轉(zhuǎn)變,為兩個平衡態(tài)的過渡階段, 稱為“過渡期”。
圍繞珠江磨刀門河口徑潮動力異變問題, 采用流量驅(qū)動的R_TIDE 數(shù)據(jù)驅(qū)動模型對徑潮信號進行分離。表2 為磨刀門河口兩個階段沿程站點的臨界流量QC、模型主要率定參數(shù)及相關(guān)的模型評價指標(ERMS和R2)結(jié)果。C1~C9為不同分潮族所使用的參數(shù), 其中C1代表長周期氣象分潮, 其周期約為15 d, 反映一個大小潮周期變化;C2、C3分別代表全日和半日分潮族,為主導(dǎo)該區(qū)域潮汐動力的分潮族;C4~C9分別代表短周期波動, 由分潮淺水變形效應(yīng)(如M4、MS4、M6、M8)產(chǎn)生, 反映地形、底床摩擦等邊界條件的影響。由表2 可知, 越往下游ERMS值越大, 表明模型對潮流優(yōu)勢區(qū)域水位的重構(gòu)效果略低于上游區(qū)域, 但是其R2均大于0.79, 模型結(jié)果合理。對于馬口站, 由于其為流量驅(qū)動邊界, 除受洪峰影響外還受邊界效應(yīng)影響, 因此其ERMS相對較大。調(diào)整期階段的ERMS比平衡期普遍增大,R2減小, 表明強人類活動對地形的改造導(dǎo)致水面線為適應(yīng)地形變化, 水位對流量的響應(yīng)敏感度有所降低, 由此一定程度上影響了模型的準確性。
圖2 為R_TIDE 模型在平衡期和調(diào)整期重構(gòu)得到的不同站點日均水位與實測值的對比結(jié)果, 其中黑色虛線表示實測值與重構(gòu)值完全一致。當流量大于臨界流量QC時潮波基本消失, 率定過程中流量大于QC時間段所對應(yīng)的潮波信號為虛假信號, 結(jié)果分析部分已剔除這部分信息。由圖2 可知, 馬口站下游站點模型擬合的效果較好(ERMS介于0.14~0.22 m,R2介于0.79~0.99)。但由于夏季洪峰出現(xiàn)次數(shù)多且持續(xù)時間長, 對徑潮信號分離產(chǎn)生較大影響, 模型誤差較大,因此當流量接近臨界流量QC時, 馬口站重構(gòu)水位比實測值偏大。
圖2 R_TIDE 數(shù)據(jù)驅(qū)動模型率定結(jié)果Fig.2 Calibration of the adopted R_TIDE data-driven model
由表2 可見, 調(diào)整期的QC和γ均大于平衡期, 主要原因為調(diào)整期階段潮汐動力增強, 潮區(qū)界向上游移動。流量指數(shù)γ反映地形邊界條件變化, 其值增大表明河床下切, 過水斷面增大。其中三灶站變化最大,反映入海水道水深劇增, 而上游各站點水深均有不同程度增大。竹銀、燈籠山站主要受采砂影響, 而馬口、甘竹站點除受采砂影響外, 還受上游流域水庫調(diào)蓄等影響, 導(dǎo)致馬口水文站的分流、分沙比發(fā)生變化, 進而顯著改變底床形態(tài)。竹銀及其下游站點為潮流優(yōu)勢河段,QC均一致, 而甘竹、馬口則受強烈徑流影響, 流量越大, 潮區(qū)界越往下游移動, 因此,QC逐漸減小。
表2 R_TIDE 數(shù)據(jù)驅(qū)動模型的部分參數(shù)及效果評價Tab.2 Calibrated parameters derived from R_TIDE data-driven model and its model performance
由于磨刀門水位及其梯度發(fā)生階段性異變, 導(dǎo)致潮波振幅及其梯度亦發(fā)生顯著變化。磨刀門河口潮汐動力主要由M2分潮主導(dǎo), 因此選用其作為代表性分潮進行分析。M2分潮振幅及其梯度的季節(jié)及年均變化如表3 所示。潮波向上游傳播過程中振幅沿程衰減, 洪季振幅小于枯季(三灶站除外), 表明徑流對潮波傳播存在明顯的衰減作用, 尤其是在磨刀門河口上游。平衡期燈籠山站振幅的洪枯季差異僅為0.03 m,甘竹站則為0.07 m, 調(diào)整期潮波振幅有明顯增大(三灶站除外), 且上游馬口站振幅的洪枯季差異增大(約為0.02 m), 其下游的甘竹、竹銀站則減小, 燈籠山、三灶站振幅的洪枯季差異基本不變。
表3 磨刀門河口各潮位站M2 分潮振幅及其梯度的季節(jié)變化及年均變化Tab.3 Seasonal and annual alterations in amplitude and its gradient of M2 tidal constituent observed at each tidal gauging station along the Modaomen estuary
為探究潮波傳播過程的時空不均勻性, 對比分析強人類活動前后各站點日均M2分潮振幅及相鄰站點間日均M2分潮潮波振幅梯度的階段性演變, 繪制時間過程線如圖3 所示, 圖中黑色實線表示階段性演變?yōu)?。圖3a 為強人類活動前后各站點日均M2分潮振幅變化值的年內(nèi)分布, 其中日均M2分潮振幅變化值定義為調(diào)整期的M2分潮振幅的日平均值減去平衡期所對應(yīng)的值, 記為Δη(以下日均M2分潮潮波振幅梯度變化值的定義類似, 記為Δδ)。由圖3a 可見, 站點距離口門越遠其值越大, 其中三灶站振幅年均減小0.01 m, 而馬口站則變化0.11 m。季節(jié)變化也有類似規(guī)律(馬口站夏季和甘竹站冬季除外), 表明振幅變化量的季節(jié)性差異有向上游累積的效應(yīng)。日均潮波振幅梯度的時空變化如圖3b 所示。潮波振幅梯度最小值出現(xiàn)在夏季, 最大值出現(xiàn)在冬季。三灶-甘竹河段的潮波振幅梯度在調(diào)整期階段的洪枯季差異減小(各河段洪枯季差異的階段性變化分別為 7.14×10-7、4.86×10-6、5.11×10-6m-1), 且越往上游減小量越大, 表明強人類活動對潮波振幅梯度的增大效應(yīng)在甘竹以下河段都有累積效應(yīng), 但在甘竹-馬口河段, 洪枯季差異雖減小但幅度小于竹銀-甘竹河段(洪枯季差異的階段性變化為4.24×10-6m-1), 表明強人類活動干預(yù)后甘竹-馬口河段潮波振幅梯度對流量響應(yīng)減弱。從潮波振幅梯度變化量來看, 三灶-燈籠山河段在秋季變化最大,燈籠山-甘竹河段則為夏季變化最大, 甘竹-馬口河段則是夏季。上述結(jié)果表明三灶-燈籠山河段水動力主要受海平面變化控制(秋季三灶余水位最大); 而甘竹-馬口河段受徑流控制為主, 由于該區(qū)域在調(diào)整期明顯受水庫調(diào)蓄影響(即洪季蓄水、枯季放水), 雖然潮波振幅梯度夏季階段性變化最大, 冬季最小, 但是其洪枯季差異減小, 是流量季節(jié)性差異減小的直接反映。
圖3 強人類活動前后磨刀門河口日均M2 分潮潮波振幅及其梯度的階段性差異Fig.3 Stepwise alteration in tidal amplitude and gradient under intensive human activities in the Modaomen estuary
圖4 為三灶-馬口河段M2分潮日均潮波振幅梯度與馬口站日均流量的關(guān)系曲線(圖4a)及其季節(jié)變化(圖4b~4e)。平衡期階段三灶-馬口河段潮波振幅梯度極小值(記為δM)為-1.41×10-5m-1, 對應(yīng)流量閾值為10 500 m3/s, 而在調(diào)整期階段則分別變?yōu)?1.38×10-5m-1和19 500 m3/s (其中流量閾值定義為當潮波振幅梯度達到極小值所對應(yīng)的流量值, 記為QM)。日均尺度下, 閾值效應(yīng)在調(diào)整期的冬季基本消失。從季節(jié)變化來看, 平衡期的秋季潮波振幅梯度極小值最小, 調(diào)整期則為夏季最小, 冬季最大, 大流量下潮波衰減效應(yīng)更明顯。調(diào)整期階段潮波振幅梯度極小值(約為-1.38×10-5m-1)和對應(yīng)流量閾值(約為19 500 m3/s)均大于平衡期。
圖4 三灶-馬口河段M2 分潮日均潮波振幅梯度(δ)與馬口站日均流量(Q)的閾值關(guān)系(a)及其季節(jié)變化(b~e)Fig.4 Relationship between the daily averaged M2 tidal amplitude gradient over SZ-MK reach (δ) and the river discharge (Q) observed at the MK hydrological station to illustrate the threshold effect (a) and its seasonal alteration (b~e)
潮波振幅梯度是徑潮動力相互作用的直接結(jié)果,因此類似流量閾值效應(yīng), 外海潮汐動力變化(即周期約為15 d 的大小潮)亦會導(dǎo)致潮波振幅梯度產(chǎn)生振幅閾值效應(yīng)(記為ηM)。圖5 為磨刀門河口三灶-馬口河段日均M2分潮潮波振幅梯度與三灶站日均M2分潮振幅的關(guān)系曲線。由圖5 可見, 日均潮波振幅梯度極小值在調(diào)整期略有增大, 約為2.49×10-7m-1, 對應(yīng)的振幅閾值基本不變(約增大0.14 cm)。圖5b~5e 為日均M2分潮潮波振幅梯度與三灶站日均M2分潮振幅閾值關(guān)系的季節(jié)變化。與圖4 類似, 除調(diào)整期的冬季振幅閾值效應(yīng)消失外(圖5e), 其余季節(jié)均存在較明顯的振幅閾值效應(yīng), 其中夏季振幅梯度極小值最小, 其次為春、秋兩季, 最大為冬季。而平衡期階段秋季的振幅閾值最小, 冬季最大。調(diào)整期冬季閾值效應(yīng)消失的主要原因是下泄流量較小, 即便經(jīng)過流域水庫調(diào)蓄,調(diào)整期口門振幅依舊未能達到使衰減效應(yīng)達到最大的臨界值(約0.43 m)。達到振幅閾值前, 由于大潮流速振幅較大導(dǎo)致其底床摩擦較大, 衰減效應(yīng)增強, 潮波振幅梯度減小。而達到振幅閾值后, 流速振幅增大導(dǎo)致的潮能耗散不足以平衡地形輻聚增強(大潮時水深增大)導(dǎo)致的潮能增加, 因此, 潮波振幅梯度略微增大。圖6 為相鄰站點間沿程潮波振幅梯度隨口門三灶站振幅的變化曲線。結(jié)果表明, 振幅閾值效應(yīng)僅存在于徑流動力占優(yōu)勢的河段, 主要存在于甘竹-馬口河段(僅調(diào)整期的冬季消失), 另外, 在竹銀-甘竹河段(春、夏季和平衡期的秋季)及燈籠山-竹銀河段(調(diào)整期的夏季)也有類似現(xiàn)象, 而在三灶-燈籠山河段, 該現(xiàn)象則不存在。對于甘竹-馬口河段, 平衡期階段, 其秋季潮波振幅梯度極小值最小, 為-3.55×10-5m-1,冬季值最大, 為-3.52×10-5m-1。雖然冬季余水位減小、底床摩擦增大, 但由于流量減小, 對潮能的耗散作用亦減弱, 潮波振幅梯度反而增大, 表明在年內(nèi)變化中甘竹-馬口河段主要受徑流影響。而在調(diào)整期階段, 冬季振幅閾值效應(yīng)消失, 這是由于冬季衰減效應(yīng)減小, 尚未達到閾值(約為-3.23×10-5m-1); 其余季節(jié)均存在閾值效應(yīng), 表明要在更大的外海潮汐動力驅(qū)動下才能使該河段發(fā)生閾值效應(yīng), 這與強人類活動驅(qū)動下各站點(三灶站除外)振幅均有不同程度的抬升有關(guān), 即兩個站點間振幅差值需要達到臨界值, 潮波振幅梯度極小值才會出現(xiàn)。
圖5 三灶-馬口河段M2 分潮日均潮波振幅梯度(δ)與三灶站日均M2 分潮振幅(ηSZ)的閾值關(guān)系(a)及其季節(jié)變化(b~e)Fig.5 Relationship between the daily averaged M2 tidal amplitude gradient (δ) in SZ-MK sector and the M2 tidal amplitude (ηSZ)observed at the SZ tidal gauging station to illustrate the threshold effect (a) and its seasonal alteration (b~e)
圖6 相鄰站點間河段M2 分潮日均潮波振幅梯度(δ)與三灶站日均振幅(ηSZ)關(guān)系曲線的季節(jié)變化Fig.6 Seasonal alteration in the relationship between the daily averaged M2 tidal amplitude gradient of the sub-sector along neighboring tidal gauging stations and tidal amplitude observed in SZ tidal gauging station
對于下游河段(竹銀以下), 該區(qū)域主要受口門圍墾、河道整治工程等地形邊界影響, 對徑流非線性作用的響應(yīng)不明顯, 潮波振幅梯度隨口門振幅的變化基本呈線性關(guān)系, 而燈籠山-竹銀河段出現(xiàn)閾值效應(yīng)(平衡期的夏季)主要是由于夏季衰減效應(yīng)已達到調(diào)整期階段的閾值(約為-3.92×10-5m-1), 而平衡期振幅梯度極值更小, 潮波衰減尚未達到該臨界值。竹銀-甘竹河段閾值效應(yīng)出現(xiàn)的原因與甘竹-馬口河段類似。
表 4 為三灶-馬口河段及相鄰站點間河段的潮波振幅梯度閾值及其對應(yīng)的振幅閾值。沿程潮波振幅梯度極小值變化規(guī)律與該河段M2分潮潮波振幅梯度的空間變化規(guī)律相似, 竹銀-甘竹河段衰減效應(yīng)減弱,到甘竹-馬口河段衰減效應(yīng)又增大, 表明甘竹站以上河段受強烈徑流調(diào)制, 潮波衰減效應(yīng)增強。在上游河段振幅沿程減小, 若要達到潮波振幅梯度閾值已不需要更大的口門振幅, 因此, 對應(yīng)振幅閾值沿程略有減小。
表4 三灶-馬口河段M2 分潮潮波振幅梯度極小值(δM)及其對應(yīng)的振幅閾值(ηM)的季節(jié)變化Tab.4 Seasonal alterations in the minimum of daily averaged M2 tidal amplitude gradient (δM) in SZ-MK sector and the critical tidal amplitude (ηM) observed at SZ tidal gauging station
上述結(jié)果表明, 磨刀門河口流量-余水位及其梯度關(guān)系在調(diào)整期階段發(fā)生顯著變化, 同流量下沿程各站點余水位及其梯度明顯下降(張蔚等, 2008; Yanget al, 2020)。然而三角洲頂端的馬口水文站年均流量僅減小816 m3/s (約為11%), 下游口門處三灶站M2分潮振幅年均值僅變化0.01 m, 表明上下游動力邊界的階段性變化不是導(dǎo)致潮波振幅梯度與上下游關(guān)系異變的主要原因, 地形異變及其導(dǎo)致的底床摩擦效應(yīng)變化才是導(dǎo)致徑潮動力格局轉(zhuǎn)換的主導(dǎo)因素。
余水位梯度作為表征徑潮動力的重要參數(shù), 其形成演變是探究地形異變對潮波傳播影響機制的重要依據(jù)和有效切入點, 其解析表達式可由一維動量守恒方程得到(Savenijeet al, 2005)
式中,x為距離口門處的距離;ρ為水體密度;U為斷面平均流速;Z為自由水面高程;D為水深;g為重力加速度;K為曼寧摩擦系數(shù)的倒數(shù)。忽略密度效應(yīng)和對流加速度項并假定流速具有周期性, 式(15)在一個潮周期內(nèi)積分可得式(16)(Caiet al, 2014)。
式(16)為余水位梯度的解析表達式, 其中Z為一個潮周期內(nèi)的余水位。表明余水位梯度主要與有效摩擦項相平衡(Buschmanet al, 2009; Sassiet al, 2013), 即余水位梯度越大, 潮波所受有效摩擦越大。圖7 為磨刀門河口全河段及各河段M2分潮日均潮波振幅梯度和日均余水位梯度的關(guān)系曲線。三灶-馬口河段的曲線均為下凹形, 表明余水位梯度增大過程中, 潮波振幅梯度減小速度減緩并趨近于其極小值, 達到極小值后變大, 表明潮波振幅梯度的變化已不受底床摩擦主導(dǎo)。竹銀以下河段曲線為上凸形, 表明余水位梯度增大時有效摩擦增大, 在該區(qū)域中, 靠近口門的三灶-燈籠山河段(圖7b)在調(diào)整期變化速率更快, 該變化主要由口門圍墾導(dǎo)致的水深變淺主導(dǎo), 而燈籠山-竹銀河段(圖7c)則受采砂和圍墾的共同影響。竹銀-甘竹河段(圖7d)曲線由平衡期的下凹轉(zhuǎn)變?yōu)檎{(diào)整期的上凸形, 根據(jù)該河段余水位梯度和潮波振幅梯度曲線,相同的余水位梯度條件下, 潮波振幅梯度增大, 表明有效摩擦對其影響減弱, 這與航道疏浚后水深增大有關(guān)。甘竹-馬口河段(圖7e)曲線為下凹形, 這表明影響該河段潮波傳播的多種因素耦合作用增強, 調(diào)整期余水位梯度顯著減小而水深增大, 當兩種效應(yīng)的耦合作用達到平衡時即達到潮波振幅梯度閾值, 因此調(diào)整期階段潮波振幅梯度閾值大于平衡期。
圖7 相鄰站點間河段M2 分潮日均潮波振幅梯度(δ)與日均余水位梯度(S)的變化Fig.7 Stepwise alteration in the relationship between the daily averaged M2 tidal amplitude gradient (δ) in SZ-MK sector and the corresponding residual water level gradient (S) observed in each sub-sector
余水位梯度的階段性變化是水面線適應(yīng)地形變化的直接反映。在人類活動干預(yù)前, 河床緩慢抬升,底坡降變化不大。1962~1977 年為緩慢淤積狀態(tài), 磨刀門河口竹洲頭-大排沙河段、大排沙-燈籠山河段平均水深分別減小0.49, 0.28 m; 1977~1999 年, 水道受到明顯沖刷, 河道深泓下切深度為0.59~1.70 m 不等,平均水深增加1.13 m (劉鋒等, 2011); 磨刀門河口沿程平均寬深比從1960 年的6.25 減小到1999 年的4.73(河道趨向窄深化趨勢發(fā)展), 同期河槽容積增大34.15×106m3(Liuet al, 2019)。表5 統(tǒng)計了三灶-竹銀河段1964、1977、1999、2016 年各水深區(qū)間的容積和面積。按水深分成0~2、2~5、5~10 和>10 m 四個區(qū)間, 分別計算各區(qū)間的河槽容積和面積。1964~1977 年, 各水深段容積和面積均顯著減小, 反映以圍墾為主導(dǎo)的淤積作用。1977~1999 年, 口門圍墾和采砂共存, 前者導(dǎo)致口門大部分水域被圍填成陸地, 0~2 m 水深區(qū)間容積和面積均銳減, 然而水深5 m 以下容積和面積顯著增大, 表明河床大幅下切,兩岸侵蝕效應(yīng)增強, 河道被沖刷, 容積和面積增加量分別為每年3.78×106m3和0.90×106m2。1999~2016年, 強人類活動逐漸減弱, 但潮波傳播過程已發(fā)生異變, 加上上游水庫建設(shè)導(dǎo)致泥沙減少, 導(dǎo)致河道水沙分配格局發(fā)生變化, 潮汐動力顯著增強, 沖刷效應(yīng)進一步加劇, 5 m 以下水深容積持續(xù)增大, 增加量為每年1.06×106m3, 面積僅增大2.93×106m2。當沖淤態(tài)勢由淤積轉(zhuǎn)變?yōu)闆_刷時, 河道容積增大, 底床摩擦減小, 潮波振幅梯度增大。
表5 三灶-竹銀河段1964、1977、1999、2016 年各水深區(qū)間的河槽容積和面積對比Tab.5 Comparison in water volume and surface area of the channel in different water depths in SZ-ZY sector observed in 1964, 1977,1999, 2016
本文基于珠江磨刀門河口沿程4 個潮位站的高、低潮數(shù)據(jù)及相應(yīng)時段馬口水文站高、低潮數(shù)據(jù)和流量數(shù)據(jù), 采用基于流量驅(qū)動的R_TIDE 數(shù)據(jù)驅(qū)動模型對河口潮波傳播的異變過程進行分析, 根據(jù)余水位及其梯度、M2分潮振幅及其梯度的階段性變化重點探討潮波振幅梯度與上下游動力邊界的關(guān)系異變過程及機制, 主要結(jié)論如下:
(1) 與平衡期相比, 在強人類活動影響后的調(diào)整期階段, 各站點日均M2分潮振幅及其梯度增大, 其中振幅在秋季變化最大(燈籠山-馬口分別為0.06、0.10、0.10、0.13 m), 冬季最小(馬口除外, 燈籠山、竹銀、甘竹分別為0.05、0.08、0.07 m), 年內(nèi)差異不大; 潮波振幅梯度變化則較為復(fù)雜, 三灶-甘竹河段沿程減小, 甘竹-馬口河段階段性變化增大。
(2) 日均M2分潮潮波振幅梯度與下游動力邊界之間存在閾值效應(yīng), 主要出現(xiàn)在甘竹-馬口河段。強人類活動驅(qū)動下潮波振幅梯度閾值增大(三灶-馬口河段增量約為2.50×10-7m-1), 振幅閾值基本不變(僅增大約0.14 cm)。由于底床有效摩擦減小導(dǎo)致沿程潮能耗散減弱, 因此當口門振幅相同時, 調(diào)整期的潮波振幅梯度大于平衡期, 隨著振幅增大, 底床摩擦效應(yīng)越明顯, 導(dǎo)致調(diào)整期潮波衰減速度減緩, 當兩個階段的口門振幅均達到振幅閾值時(此處假定振幅閾值相同),則調(diào)整期對應(yīng)的潮波振幅梯度大于平衡期。
(3) 對于潮波振幅梯度與大小潮的關(guān)系異變分析需綜合考慮地形、摩擦、流量等多因素影響。由于上下游動力邊界的階段性變化不足以導(dǎo)致上述關(guān)系的異變, 因此地形及其引起的底床摩擦變化則是徑潮動力異變的根本原因。從1964~2016 年, 10 m 以下水深的容積增大15.26×106m3, 0~2 m 水深容積則減小366.94×106m3, 反映河道窄深化的演變趨勢。由于地形的異變, 在口門振幅增大過程中, 達到振幅閾值前,由于底床摩擦增大, 小潮時的潮波振幅梯度大于大潮, 而達到振幅閾值后, 此時以地形輻聚為主導(dǎo)影響因子導(dǎo)致潮波振幅梯度增大。