翟 瑋,張 蔚,季小梅,陳 婷

(1.河海大學(xué) 江蘇省海岸海洋資源開發(fā)與環(huán)境安全重點(diǎn)實(shí)驗(yàn)室,南京 210024; 2.上海灘涂海岸工程技術(shù)研究中心,上海 200061; 3.河海大學(xué) 疏浚技術(shù)教育部工程研究中心,南京 210024)

1 研究背景

河口是陸海交匯的過渡地帶,受到流域徑流和海洋潮汐的共同作用,水動(dòng)力條件復(fù)雜。余水位即潮周期內(nèi)的平均水位,濾去了周期性變化的潮汐,其時(shí)空變化反映了徑潮動(dòng)力間的非線性作用[1]。根據(jù)一維動(dòng)量守恒方程,余水位梯度項(xiàng)基本與摩阻項(xiàng)相平衡[2]。為探究徑流、潮汐及其相互作用對余水位的影響,多從動(dòng)量方程的摩阻項(xiàng)入手,將其分解為徑流作用、徑潮相互作用和潮不對稱作用產(chǎn)生的摩阻項(xiàng)[3-6]。近年來隨著技術(shù)手段的發(fā)展,學(xué)者們提出了新的方法來研究徑潮間的非線性作用,包括連續(xù)小波變換[7-8]、短期調(diào)和分析[9]、非穩(wěn)態(tài)調(diào)和分析[10]等。這些分析方法結(jié)合摩阻項(xiàng)的分解,有助于解釋評價(jià)徑流、潮汐及其相互作用對余水位的相對重要性。如Buschman等[4]基于連續(xù)小波變換方法和對摩阻項(xiàng)中非線性流速乘積的一階及三階展開,建立線性回歸模型,通過潮平均摩阻的變化預(yù)測余水位的變化。

長江口是世界第三大河口,徑流量巨大且季節(jié)性變化特征明顯,潮汐強(qiáng)度中等,徑潮相互作用顯著,深刻影響著此處的余水位。過去對長江口余水位的研究,大多采用實(shí)測資料的統(tǒng)計(jì)分析方法。劉曉婉等[11]發(fā)現(xiàn)崇明島南部近岸水位主要由潮汐驅(qū)動(dòng)形成。郭磊城等[12]基于實(shí)測徑流和潮位資料發(fā)現(xiàn),長江口洪季徑流量的增大使得余水位顯著抬高,洪枯季差異從上游大通到下游牛皮礁逐漸減小。楊正東等[13]發(fā)現(xiàn)長江口余水位的變化主要與徑流量有關(guān)。此外,宋永港等[14]基于數(shù)值模擬,對長江口北支潮位潮差時(shí)空變化展開研究,發(fā)現(xiàn)1—7月份的月平均潮位逐漸增大,而8—12月份逐漸減小,該變化趨勢由徑流產(chǎn)生的余水位決定。綜上而言,針對徑潮動(dòng)力下余水位時(shí)空變化的研究仍有待深入。余水位影響著水位與水深,同時(shí)作為徑潮動(dòng)力非線性的典型結(jié)果,在長江河口動(dòng)力格局變化的當(dāng)下也成為分析徑潮相互作用的有效切入點(diǎn),對防洪安全、水位預(yù)報(bào)、水資源管理等科學(xué)問題具有重要意義。

本文基于二維水動(dòng)力數(shù)值模型,獲取長江口主要站點(diǎn)高頻、長期的潮位和流速數(shù)據(jù),分析余水位時(shí)空變化特征;利用連續(xù)小波變換分解流速信號,通過摩阻公式計(jì)算潮平均摩阻的各分量,結(jié)合簡單線性回歸模型分析徑潮動(dòng)力相互作用下余水位的變化機(jī)制。

2 研究區(qū)域及研究方法

2.1 研究區(qū)域

長江口是銜接了長江和東中國海的主要通道,全長約660 km,上至安徽大通,下至外海50 m等深線附近。長江口整體平面形態(tài)呈現(xiàn)為“三級分汊、四口入?！钡膹?fù)雜格局,河口呈典型喇叭型,上游寬度13 km,外?？陂T處寬度90 km,其地形及主要站點(diǎn)位置如圖1所示。由于屬亞熱帶季風(fēng)氣候,河口徑流季節(jié)性變化顯著,5—10月份徑流量占年徑流量71.7%,11月份至次年4月份徑流占28.3%[15]。全日分潮和半日分潮是進(jìn)入河口的主要分潮,其振幅比AD1/AD2約為0.24,河口內(nèi)潮汐類型屬于非正規(guī)半日潮。

圖1 長江口地形及主要站點(diǎn)位置Fig.1 Map of the Yangtze River Estuary and locations of hydrologic stations

2.2 研究方法

2.2.1 水動(dòng)力模型

本文采用課題組已建立的長江口整體平面二維水動(dòng)力數(shù)學(xué)模型,對徑潮動(dòng)力相互作用下余水位時(shí)空變化展開研究。模型采用非正交的三角形網(wǎng)格,覆蓋整個(gè)長江口區(qū)域,陸地邊界能較好地貼合實(shí)際岸線,模型網(wǎng)格總數(shù)為157 311,網(wǎng)格節(jié)點(diǎn)數(shù)為75 154,網(wǎng)格大小在50～10 000 m之間變化,精度滿足需要。模型地形資料采用長江水利委員會(huì)長江下游水文水資源勘測局測得的長江口水下實(shí)測地形,并通過14個(gè)潮位站點(diǎn)和5個(gè)流速站點(diǎn)的實(shí)測數(shù)據(jù)進(jìn)行率定與驗(yàn)證。大多數(shù)測站水位和流速數(shù)據(jù)誤差均在10%以內(nèi),結(jié)果良好,因此可認(rèn)為模型模擬結(jié)果能較好地反映長江口的水動(dòng)力條件。具體率定與驗(yàn)證過程見文獻(xiàn)[16]。本文模型參數(shù)設(shè)置如表1所示。

表1 模型參數(shù)設(shè)置Table 1 Setting of model parameters

為了充分消除潮汐的影響,模型上邊界設(shè)置在枯季的潮區(qū)界大通,并將2016年大通實(shí)測徑流量作為上邊界條件;外海邊界則設(shè)置在徑流作用完全消失的124.5°E處,分潮潮位由東中國海潮波數(shù)學(xué)模型[17]提供。模型采用冷啟動(dòng)方式,對2016年長江口的潮位和流速開展模擬。

2.2.2 連續(xù)小波變換及潮波分解

作為非穩(wěn)態(tài)調(diào)和分析的重點(diǎn)方法之一,連續(xù)小波變換(Continuous Wavelet Transform,CWT)通過對母小波進(jìn)行縮放和平移,對潮波信號進(jìn)行了多尺度的細(xì)化分析,利用小波變換在一個(gè)帶通濾波器上對潮波信號進(jìn)行處理,得到信號隨著時(shí)間和頻率的變化過程[18-19]。

CWT方法的特點(diǎn)在于它可以比較準(zhǔn)確地得到不同周期分潮簇的振幅與相位,非常適應(yīng)潮汐頻帶非幾何變化間隔的特性;其缺點(diǎn)是不能較好地區(qū)分同一分潮簇內(nèi)的不同分潮[8],且對小波的選擇具有較強(qiáng)的依賴性。

在徑潮動(dòng)力相互作用的背景下,CWT方法可以顯露潮波對變化的徑流量近乎瞬時(shí)的響應(yīng)[9]。利用CWT方法處理流速數(shù)據(jù),可將原始流速時(shí)間序列分解為不同時(shí)頻域信息,將同一周期的分潮簇結(jié)果歸類合并,從而得到主要的全日潮簇D1、半日潮簇D2以及四分之一日潮簇D4流速的振幅(U1、U2及U4)和初相位(φ1、φ2及φ4)。

2.2.3 摩阻效應(yīng)

根據(jù)Buschman等[4]提出的針對余水位變化的新的分析方法,可將動(dòng)量方程中的摩阻項(xiàng)按照其來源分解為徑流作用、徑潮相互作用和潮不對稱作用產(chǎn)生的3個(gè)不同分量Sr、Srt及St。該方法的核心在于利用切比雪夫多項(xiàng)式對摩阻項(xiàng)中的非線性乘積U|U|進(jìn)行一階和三階展開,并將流速時(shí)間序列分解后經(jīng)無量綱化處理[5]代入,即:

(1)

(2)

對摩阻項(xiàng)在一個(gè)全日潮周期內(nèi)進(jìn)行積分,整理后得到公式如下:

(3)

(4)

(5)

日均流速和分潮簇流速振幅在經(jīng)過無量綱化處理后均<1,結(jié)合式(3)—式(5)分析可知,Sr、Srt和St大致范圍在0～1范圍內(nèi),反映各項(xiàng)作用對潮平均摩阻的貢獻(xiàn)。其中,Sr只與日均流速有關(guān),因而表達(dá)的是徑流作用對潮平均摩阻的貢獻(xiàn);Srt與平均流速和各分潮簇流速振幅的乘積相關(guān),反映徑潮相互作用對潮平均摩阻的貢獻(xiàn);而St則與分潮簇流速振幅和相位有關(guān),反映潮不對稱作用對潮平均摩阻的貢獻(xiàn)。

3 余水位時(shí)空變化特征

根據(jù)二維水動(dòng)力模型模擬結(jié)果,對長江口主要站點(diǎn)的模擬水位取日均值,從而濾去全日潮簇及更高頻分潮簇潮位對水位周期性的影響,得到余水位。圖2為2016年大通站的日均徑流量及主要站點(diǎn)的余水位。

圖2 2016年大通站日均徑流量及主要站點(diǎn)余水位Fig.2 Average daily runoff at Datong station and residual water level at major hydrological stations in 2016

從圖2(a)可以看出,長江口的徑流量季節(jié)性變化明顯,洪季(5—10月份)流量明顯大于枯季(1—4月份及11—12月份)流量,日均徑流量最大可達(dá)70 000 m3/s。根據(jù)圖2(b)可知,長江口余水位從下游吳淞向上游南京逐漸抬升,具有穩(wěn)定的梯度。在枯季,余水位沿程增幅較小,以兩周為周期的大小潮變化始終明顯,大潮時(shí)余水位要高于小潮時(shí);在洪季,余水位沿程增幅較大,在上游站點(diǎn)其變化與徑流量的變化十分相似,往下游余水位大小潮變化逐漸明顯。上述長江口余水位的時(shí)空變化體現(xiàn)了徑潮動(dòng)力的共同作用,但其變化原因還需通過理論方法作出解釋,具體見本文5.2節(jié)末尾。

4 分潮簇流速振幅時(shí)空變化

基于二維水動(dòng)力模型的模擬結(jié)果和連續(xù)小波變換方法,對各站點(diǎn)流速時(shí)間序列進(jìn)行調(diào)和分析,得到全日潮簇(D1)、半日潮簇(D2)和四分之一日潮簇(D4)流速的調(diào)和常數(shù)。圖3為各分潮簇流速振幅的洪枯季沿程分布。

圖3 各周期分潮簇流速振幅洪枯季沿程分布Fig.3 Distribution of velocity amplitude of tidal clusters in different periods during wet season and dry season

由圖3可知,長江口內(nèi)半日潮簇D2流速振幅最大,枯季在吳淞處可達(dá)0.9 m/s,全日潮簇D1和四分之一日潮簇D4流速振幅大致相當(dāng),均<0.2 m/s。潮波從下游的吳淞向上游傳播過程中,受徑流頂托和摩阻耗散,各分潮簇流速振幅大體上沿程衰減。其中在吳淞至徐六涇段,斷面收縮明顯,潮能集中,因此流速振幅衰減較慢。與其他兩項(xiàng)不同,D4為淺水分潮,其能量也能通過淺水非線性作用從天文分潮轉(zhuǎn)移而來,因此在吳淞至徐六涇段流速振幅明顯增大。

對比洪枯季分潮簇流速振幅,各站點(diǎn)枯季流速振幅均大于洪季,可以得知洪季更大的徑流量加劇了分潮的衰減,洪季潮波上溯至鎮(zhèn)江處就基本消散。

5 摩阻作用及對余水位的控制

5.1 潮平均摩阻時(shí)空變化

根據(jù)Buschman等[4]推導(dǎo)的式(3)—式(5),計(jì)算洪枯季長江口主要站點(diǎn)徑流作用、徑潮相互作用、潮不對稱作用對潮平均摩阻的貢獻(xiàn),結(jié)果見圖4(其中S=Sr+Srt+St)。

圖4 洪枯季徑流作用、徑潮相互作用、潮不對稱作用對潮平均摩阻的貢獻(xiàn)Fig.4 Contributions of river flow, river-tide interactions and tidal asymmetry to tidally averaged friction in wet season and dry season

從圖4可以看出,長江口的潮平均摩阻以徑流作用及徑潮相互作用產(chǎn)生的摩阻為主,潮不對稱作用的貢獻(xiàn)微乎其微。究其原因,分潮簇流速振幅一般比日均流速小,加上算式內(nèi)還包含分潮簇相位差的影響,使得St較小。

分析各摩阻分量的沿程變化規(guī)律,無論洪季還是枯季,Sr都從吳淞向上游逐漸增加,在南京取得最大值。產(chǎn)生這種現(xiàn)象的原因是:自上游向下傳播過程中,由于和潮汐的相互作用以及河道斷面的放寬,徑流逐漸減弱,因此徑流作用對潮平均摩阻的貢獻(xiàn)逐漸減小。Srt從吳淞向上游大體上逐漸減小,在枯季沿程減小較慢,在洪季減小較快,這與圖3中分潮簇洪枯季流速振幅沿程變化相似,表明徑流增加會(huì)加劇潮波的衰減,從而顯著降低徑潮相互作用對潮平均摩阻的貢獻(xiàn)。

對比洪枯季各站點(diǎn)Sr和Srt兩項(xiàng),Sr均>Srt,說明在長江高徑流的背景下徑流作用對潮平均摩阻的貢獻(xiàn)始終更大。在江陰及以上站點(diǎn),Sr占據(jù)主導(dǎo),而在下游的徐六涇和吳淞站點(diǎn),Srt對潮平均摩阻的貢獻(xiàn)也較大,不可忽略。

5.2 摩阻分量對余水位的擬合

在一維動(dòng)量守恒方程中,余水位梯度項(xiàng)基本與潮平均摩阻項(xiàng)相平衡,摩阻分量與余水位梯度間大致呈線性關(guān)系,因此通過摩阻分量的線性回歸可以預(yù)測擬合余水位的變化。由5.1節(jié)可知St較小,因而選取潮平均摩阻中主要的Sr和Srt兩項(xiàng),建立二元線性回歸模型,在各站點(diǎn)根據(jù)二者的日均值預(yù)測擬合對應(yīng)的余水位,結(jié)果見圖5。

圖5 二元線性回歸模型預(yù)測的余水位和數(shù)值模型模擬的余水位對比Fig.5 Comparison of residual water level predicted by binary linear regression model and simulated by numerical model

在線性回歸方程中,p值反映檢驗(yàn)的顯著性水平,而決定系數(shù)R2反映自變量對于因變量變化的解釋程度,本文中即余水位變化被摩阻分量變化定性解釋的程度。

計(jì)算表明,所選取的5個(gè)站點(diǎn)二元線性回歸模型的p值均<0.01,說明檢驗(yàn)具有顯著統(tǒng)計(jì)學(xué)差異,擬合方程可靠。圖5(a)—圖5(e),從上游南京至下游吳淞各站點(diǎn)二元線性回歸各模型的R2分別為0.958、0.955、0.891、0.711、0.453,除最下游的吳淞站點(diǎn)外,均擬合良好,說明Sr和Srt的變化能較好地解釋余水位的變化。除Sr和Srt的二元線性回歸外,采用R2的大小作為依據(jù),對比Sr、Srt、S各自的一元線性回歸擬合結(jié)果,如表2所示。

表2 一元、二元線性回歸擬合對比Table 2 Comparison between univariate and binary linear regression fittings

從表2可看出,摩阻分量之和S的一元線性回歸即可較好地預(yù)測余水位的變化,體現(xiàn)出摩阻控制水位。由此可對圖5中洪季擬合結(jié)果較枯季更好作出解釋。因?yàn)楹榧境逼骄ψ韪?動(dòng)量方程中其他項(xiàng)的影響顯得相對較小,因此余水位變化能被摩阻變化更好地預(yù)測。同理,由于吳淞站點(diǎn)位置最靠下,在洪季徑流作用也較弱,因而擬合結(jié)果比其上游的站點(diǎn)要差一些。由于Sr始終>Srt,因此前者的一元線性擬合結(jié)果更優(yōu)。在江陰以下站點(diǎn),Srt對總摩阻的貢獻(xiàn)變大,因此考慮了Srt變化的二元線性回歸的擬合結(jié)果較一元線性回歸有一定提升。

根據(jù)上述分析對圖2所示長江口余水位時(shí)空變化進(jìn)行解釋。由于潮平均摩阻一直存在,平衡著余水位梯度,因此自吳淞向上游余水位一直增大。洪季徑流量巨大,徑流作用產(chǎn)生的摩阻占據(jù)主導(dǎo)且遠(yuǎn)大于枯季,使得洪季余水位梯度更大,余水位向上游增幅明顯?？菁緩搅髁孔兓淮?余水位的大小潮變化來源于潮平均摩阻的大小潮變化,歸根到底是外海潮差變化所引起。

6 結(jié) 論

本文基于長江口二維水動(dòng)力模型獲得主要站點(diǎn)一年的潮位和流速數(shù)據(jù),分析余水位時(shí)空變化特征;利用連續(xù)小波變換方法對流速進(jìn)行分解,得到主要分潮簇流速的調(diào)和常數(shù);通過摩阻公式計(jì)算徑流作用、徑潮相互作用和潮不對稱作用對潮平均摩阻的貢獻(xiàn),并結(jié)合簡單線性回歸模型分析徑潮動(dòng)力相互作用下余水位的時(shí)空變化。得到的主要結(jié)論如下:

(1)在徑潮動(dòng)力相互作用下,長江口余水位自下游向上游逐漸抬升,表現(xiàn)出明顯的洪枯季變化和大小潮變化;洪季余水位沿程增幅較大,枯季余水位沿程增幅較小,大潮時(shí)余水位比小潮時(shí)更高。

(2)長江口的潮平均摩阻以徑流作用和徑潮相互作用產(chǎn)生的摩阻為主,潮不對稱項(xiàng)的貢獻(xiàn)可忽略,且洪季的摩阻大于枯季;因徑流量巨大,各站點(diǎn)徑流作用對潮平均摩阻的貢獻(xiàn)始終超過徑潮相互作用。

(3)潮平均摩阻控制余水位,通過建立線性回歸模型用摩阻的變化較好地解釋說明了余水位的變化,所采用的回歸模型包括一元與二元,決定系數(shù)R2多達(dá)到0.7以上;長江余水位的變化主要由徑流作用控制,在下游站點(diǎn)考慮了徑潮相互作用的二元回歸模型較一元模型擬合效果有一定提升,徑潮相互作用的影響不可忽視。