王 青，葉榮輝，汪玉平，張夢(mèng)茹，孔 俊

(1.河海大學(xué)港口海岸與近海工程學(xué)院，江蘇 南京 210098； 2.珠江水利委員會(huì)珠江水利科學(xué)研究院，廣東 廣州 510611)

咸潮是河口地區(qū)常見的自然現(xiàn)象，造成嚴(yán)重的生態(tài)危害，直接影響到生活用水、工業(yè)用水以及農(nóng)作物灌溉，嚴(yán)重制約經(jīng)濟(jì)發(fā)展，因而如何做到準(zhǔn)確預(yù)測(cè)咸潮，特別是預(yù)測(cè)滿足取水指標(biāo)的低鹽過(guò)程一直是個(gè)科學(xué)難題。咸潮的早期研究主要采用實(shí)測(cè)分析、物理模型試驗(yàn)和數(shù)值解析方法。Savenije[1]基于大量實(shí)測(cè)河口咸潮上溯曲線，總結(jié)出“遞降”“鈴型”“拱型”和“駝背”4種類型的上溯曲線，用來(lái)模擬河口的咸潮；Grigg等[2]開展物理水槽實(shí)驗(yàn)，對(duì)咸潮的紊動(dòng)混合強(qiáng)度及鹽度時(shí)空分布對(duì)不同地形的響應(yīng)情況進(jìn)行了分析；Schijf等[3]在鹽水楔理論的基礎(chǔ)上，推算了咸潮距離的解析解。隨著科技的進(jìn)步，遙感反演、數(shù)值模擬、數(shù)理統(tǒng)計(jì)的方法也被廣泛用于咸潮研究中。Bowers等[4]根據(jù)水體中黃色物質(zhì)與鹽度呈負(fù)相關(guān)關(guān)系，提出了黃色物質(zhì)反演水體鹽度的遙感算法；數(shù)值模擬方法是目前研究咸潮問題最普遍的方法，使用較廣泛的河口海洋水動(dòng)力數(shù)值模型包括結(jié)構(gòu)網(wǎng)格模式，如Princeton Ocean Model(POM)[5]、Estuarine，Coastal and Ocean Model(ECOM)[6]、Tidal，Residual，Inter-tidal Mudflat(TRIM)[7]、Regional Oceanic Modeling System(ROMS)[8]等，以及無(wú)結(jié)構(gòu)網(wǎng)格模式，如Unstructured TRIM Model(UNTRIM)[9]、Finite Volume Community Ocean Model(FVCOM)[10-11]、Eulerian-Lagrangian Circulation Model(ELCIRC)[12]、Semi-implicit Eulerian-Lagrangian Finite-element Model(SELFE)[13]等；Qiu等[14]利用盧薩哈奇河近20年的測(cè)量資料進(jìn)行多元回歸分析，建立鹽度預(yù)報(bào)模型，通過(guò)對(duì)時(shí)間序列的分析，研究淡水徑流、降雨和潮汐水面高程導(dǎo)致的鹽度變化。然而已有研究主要關(guān)注咸潮上溯距離[15-18]，而對(duì)關(guān)注具體地點(diǎn)的鹽度變化的研究比較少。

本文以珠江口為研究對(duì)象，基于統(tǒng)計(jì)學(xué)方法，考慮歷史鹽度過(guò)程的記憶效應(yīng)，構(gòu)建預(yù)測(cè)日均鹽度的關(guān)系函數(shù)。通過(guò)分析上游徑流和外海潮差對(duì)鹽度的影響，結(jié)合統(tǒng)計(jì)回歸模型，建立日均流量、日最大潮差、記憶鹽度和日均鹽度的關(guān)系函數(shù)，以期為河口地區(qū)取水提供參考依據(jù)。

1 研究區(qū)概況

珠江口是水系結(jié)構(gòu)、動(dòng)力特性、人類活動(dòng)都較為復(fù)雜的河口之一，具有“三江匯流、八口入?！钡乃堤卣?，地理位置見圖1，其中上游流量測(cè)站為馬口站和三水站，外海潮位測(cè)站為三灶站，鹽度測(cè)站為平崗站。為了簡(jiǎn)化計(jì)算，對(duì)馬口站和三水站的流量值進(jìn)行相加，用相加得到的流量和進(jìn)行研究。上游流量、外海潮位和鹽度資料都為2011年9月至2013年12月的每小時(shí)實(shí)測(cè)數(shù)據(jù)。本文的預(yù)測(cè)模型研究的是日均鹽度值，因此對(duì)上游流量和鹽度數(shù)據(jù)處理得到日均流量和日均鹽度，對(duì)外海潮位處理得到日最大潮差。

圖1 研究區(qū)地理位置

2 模型建立

2.1 咸潮影響因子分析

淡水徑流和潮汐動(dòng)力作用是影響咸潮的主要因素，同時(shí)風(fēng)力風(fēng)向、海平面上升以及地形等因素也對(duì)咸潮有影響。Neubauer等[19]提出在徑流占主導(dǎo)的河口，高徑流量可以抑制咸潮，反之，低徑流量下，咸潮上溯的距離更遠(yuǎn)。魏雪靜等[20]基于ECOMSED三維海洋水動(dòng)力模型，著重模擬和分析潮差對(duì)鹽度影響，指出隨著潮差的減小，表層和底層的鹽度差增大，鹽水分層現(xiàn)象更明顯，即小潮期更容易出現(xiàn)鹽水楔現(xiàn)象。潘存鴻等[21]基于長(zhǎng)系列的實(shí)測(cè)水文和地形資料，指出錢塘江口咸潮上溯受到徑流、潮汐和河床地形的影響。丁磊等[22]通過(guò)分析長(zhǎng)江口2013和2014年枯季數(shù)據(jù)指出長(zhǎng)江口咸潮上溯受地形地貌、徑流、潮流和風(fēng)因素共同影響。Ou等[23]指出徑流、近岸海流和風(fēng)對(duì)珠江口鹽度分布有明顯影響，不同形態(tài)的鹽度鋒主要是通過(guò)表層平流運(yùn)動(dòng)形成。Scully等[24]基于約克河實(shí)地觀測(cè)結(jié)果提出風(fēng)對(duì)鹽度分層有重要影響，順河口向上游的風(fēng)會(huì)減弱分層，一定強(qiáng)度的順河口向下游的風(fēng)會(huì)增加分層，而過(guò)強(qiáng)的順河口向下游的風(fēng)卻減弱分層，并提出了風(fēng)應(yīng)變理論解釋風(fēng)對(duì)河口鹽度分層的作用。近年來(lái)隨著全球氣候變暖，海平面上升伴隨著鹽水上涌，進(jìn)一步加劇了咸潮的影響態(tài)勢(shì)，Hilto等[25]在美國(guó)的切薩皮克灣觀測(cè)到了海平面上升導(dǎo)致鹽鋒上溯的現(xiàn)象。Zhang等[26]指出珠江口咸潮上溯除了受到徑流、風(fēng)和地形的影響，還受到海平面上升的影響。劉祖發(fā)等[11]基于FVCOM模型，發(fā)現(xiàn)珠江口海平面上升加劇咸潮上溯程度。

本文主要分析上游徑流和外海潮差兩個(gè)因素對(duì)咸潮的影響。采用上游日均流量代表上游徑流，根據(jù)相關(guān)文獻(xiàn)[14,27-28]，考慮用指數(shù)函數(shù)和冪函數(shù)來(lái)研究流量和鹽度之間的關(guān)系，其中指數(shù)函數(shù)形式對(duì)應(yīng)的擬合優(yōu)度指數(shù)為0.29，冪函數(shù)形式對(duì)應(yīng)的擬合優(yōu)度指數(shù)為0.27，可見采用指數(shù)函數(shù)形式較好。圖2和圖3分別為日均流量和日最大潮差與日均鹽度的關(guān)系。由圖2和圖3可見，當(dāng)日均流量較大時(shí)，日均鹽度基本被日均流量所控制，此時(shí)潮汐動(dòng)力基本不會(huì)對(duì)日均鹽度產(chǎn)生任何影響，但是影響鹽度的各種因素并不是孤立存在和單獨(dú)發(fā)揮作用的，而是相互耦合作用，考慮到日均流量和日最大潮差共同作用時(shí)可能出現(xiàn)特殊情況，將日最大潮差代替潮汐動(dòng)力作為一個(gè)影響因素。

圖2 日均流量與日均鹽度關(guān)系

圖3 日最大潮差與日均鹽度關(guān)系

2.2 統(tǒng)計(jì)模型的構(gòu)建

首先建立傳統(tǒng)的鹽度預(yù)測(cè)模型，僅考慮日均流量和日最大潮差對(duì)日均鹽度的影響，命名模型1，采用如下的函數(shù)形式：

S1=a1exp(-b1Qt+1)+cHt+1+e1

(1)

式中：S1為模型1的第t+1天的日均鹽度預(yù)測(cè)值；Qt+1為第t+1天的日均流量實(shí)測(cè)值；Ht+1為第t+1天的日最大潮差實(shí)測(cè)值；a1、b1、c1為擬合系數(shù)；e1為常數(shù)項(xiàng)。a1、b1、c1和e1的率定值分別為17.478 9、0.001 1、0.449 0和0.381 9。

為改進(jìn)傳統(tǒng)鹽度預(yù)測(cè)模型，引入自回歸模型，考慮過(guò)去鹽度的記憶效應(yīng)的影響，自回歸模型為

St=F(St-1,St-2,St-3)

(2)

式中：St為第t天的日均鹽度實(shí)測(cè)值；St-i為第t-i天的日均鹽度實(shí)測(cè)值。

采用兩種不同的方案進(jìn)行對(duì)比分析。

方案a：采用動(dòng)態(tài)模型，日均鹽度預(yù)測(cè)公式為

St+1=φ1St+φ2St-1+φ3St-2

(3)

式中：St+1為方案a的第t+1天的日均鹽度預(yù)測(cè)值；φ1、φ2、φ3為擬合系數(shù)。

方案b：引入伽馬分布函數(shù)[29]表征鹽度的記憶效應(yīng)：

(4)

(5)

(6)

圖4為方案a和方案b的日均鹽度預(yù)測(cè)值與實(shí)測(cè)值對(duì)比，圖5為方案a和方案b的日均鹽度預(yù)測(cè)峰值與實(shí)測(cè)峰值對(duì)比。由圖4、5可以看出，方案b預(yù)測(cè)的日均鹽度值明顯小于方案a，更加接近實(shí)測(cè)值，預(yù)測(cè)值和實(shí)測(cè)值的擬合優(yōu)度指數(shù)由方案a的0.635 5提高到方案b的0.820 1，而且方案a中的預(yù)測(cè)峰值不能實(shí)現(xiàn)和實(shí)測(cè)峰值的相互對(duì)應(yīng)，存在明顯的滯后性，說(shuō)明引入伽馬分布函數(shù)后能夠很好地減弱滯后效應(yīng)，方案b能更準(zhǔn)確地預(yù)測(cè)日均鹽度。

圖4 方案a和方案b的日均鹽度預(yù)測(cè)值與實(shí)測(cè)值對(duì)比

圖5 方案a和方案b的日均鹽度預(yù)測(cè)峰值與實(shí)測(cè)峰值對(duì)比

綜合考慮各項(xiàng)因素，建立一種新的鹽度預(yù)測(cè)統(tǒng)計(jì)模型，命名模型2，預(yù)測(cè)方程為

(7)

式中：S2為模型2的第t+1天的日均鹽度預(yù)測(cè)值；a2、b2、c2、d2為擬合系數(shù)；e2為常數(shù)項(xiàng)。

3 敏感性分析

3.1 徑流強(qiáng)度敏感性分析

徑流是影響咸潮最為敏感的因素之一，也是一定程度上可以受人為調(diào)控的因素?？紤]不同徑流大小對(duì)咸潮的影響，在上一節(jié)建立的模型的基礎(chǔ)上，設(shè)置日均流量增大50%和減小50%兩種試驗(yàn)，其他條件與上一節(jié)建立的模型相同，對(duì)比不同情況下上游徑流對(duì)咸潮的影響。令

SQ=a2exp(-b2Qt+1)

(8)

式中SQ為模型2中流量貢獻(xiàn)的日均鹽度預(yù)測(cè)值。

圖6為流量不變、減小50%和增大50%情況下SQ和S2的變化情況。由圖6可見，當(dāng)日均流量減少50%時(shí)，所對(duì)應(yīng)的SQ顯著增大，增幅基本在3倍以上，絕大部分在4～6倍，少數(shù)達(dá)到10倍，所對(duì)應(yīng)的S2增大幅度也比較明顯；當(dāng)日均流量增大50%時(shí)，所對(duì)應(yīng)的SQ變化幅度不大，大多數(shù)不變化，少數(shù)減少4倍，所對(duì)應(yīng)的S2減少幅度不明顯，基本沒有變化。由此可知上游徑流是影響咸潮的關(guān)鍵因素，當(dāng)上游流量比較大時(shí)，鹽度對(duì)上游流量的變化不是很敏感，而當(dāng)上游流量比較小時(shí)，鹽度對(duì)上游流量的變化比較敏感。

(b)S2

3.2 潮汐強(qiáng)度敏感性分析

選取2011年9月至2013年12月的實(shí)測(cè)數(shù)據(jù)進(jìn)行擬合，發(fā)現(xiàn)日最大潮差的系數(shù)b2僅有0.006，說(shuō)明本文情況下日最大潮差對(duì)日均鹽度的影響基本可以忽略不計(jì)。因此對(duì)鹽度預(yù)測(cè)模型進(jìn)行進(jìn)一步的改進(jìn)，命名為模型3，改進(jìn)后的預(yù)測(cè)方程為

(9)

式中：S3為模型3的第t+1天的日均鹽度預(yù)測(cè)值；a3、b3、d3為擬合系數(shù)；e3為常數(shù)項(xiàng)。

4 模型的驗(yàn)證

選取2011年9月至2012年8月的實(shí)測(cè)數(shù)據(jù)利用模型3進(jìn)行擬合并率定系數(shù)，得到a3、b3、d3和e3的率定值分別為7.50，0.001 53，0.73和0.001，擬合結(jié)果見圖7(a)，預(yù)測(cè)值和實(shí)測(cè)值的擬合優(yōu)度指數(shù)為0.909 4。對(duì)比模型1、模型3的日均鹽度預(yù)測(cè)值和實(shí)測(cè)值可見，模型1的擬合優(yōu)度指數(shù)僅為0.303 1，遠(yuǎn)小于模型3。而模型3的日均鹽度預(yù)測(cè)值無(wú)論在豐水期還是枯水期都與實(shí)測(cè)值誤差很小，說(shuō)明引入伽馬分布函數(shù)表征鹽度的記憶效應(yīng)能夠顯著提高預(yù)測(cè)的準(zhǔn)確度，同時(shí)，利用重要性系數(shù)減弱了滯后效應(yīng)帶來(lái)的影響。利用模型3，選取2012年9月至2013年12月的實(shí)測(cè)數(shù)據(jù)進(jìn)行模型驗(yàn)證，結(jié)果見圖7(b)。由圖7(b)可見，除了個(gè)別數(shù)值相差較大，其余預(yù)測(cè)值和實(shí)測(cè)值的相差很小，擬合優(yōu)度指數(shù)為0.842 1，校驗(yàn)結(jié)果比較理想，說(shuō)明本文提出的鹽度預(yù)測(cè)統(tǒng)計(jì)模型能夠準(zhǔn)確地預(yù)測(cè)珠江口的鹽度。

(a)參數(shù)率定

(b)模型驗(yàn)證

5 結(jié) 論

a. 珠江河口上游徑流對(duì)鹽度的影響顯著，日均流量和日均鹽度之間可用指數(shù)函數(shù)進(jìn)行擬合；外海潮差對(duì)鹽度的影響不顯著，在綜合考慮各項(xiàng)因素后，略去了外海潮差的影響。

b. 比較兩種日均鹽度的動(dòng)態(tài)預(yù)測(cè)方案，發(fā)現(xiàn)引入伽馬分布函數(shù)表征鹽度的記憶效應(yīng)能夠顯著提高預(yù)測(cè)的準(zhǔn)確度。同時(shí)，利用重要性系數(shù)減弱了滯后效應(yīng)帶來(lái)的影響。

c. 模型擬合優(yōu)度指數(shù)為0.842 1，校驗(yàn)結(jié)果比較理想，說(shuō)明提出的鹽度預(yù)測(cè)統(tǒng)計(jì)模型能夠準(zhǔn)確地預(yù)測(cè)珠江口的鹽度。