路川藤，錢明霞*，夏威夷，丁偉，丁佩

（1.南京水利科學研究院，江蘇南京 210029；2.水文水資源與水利工程科學國家重點實驗室，江蘇南京 210098）

0 引言

入海河口在徑流和潮流的共同作用下，形成了特有的物理現(xiàn)象—鹽淡水混合。鹽淡水混合具有時空分布不均的特點，垂向上，水體可按鹽淡水的混合程度分為高度分層型、緩混合型和強混合型[1]；平面上，水體受潮汐強度及海床地形等因素的影響，不同區(qū)域鹽水濃度不同[2]，同時，不同季節(jié)由于徑流存在差別，鹽水上溯距離也不同[3]。由于影響因素眾多，因此鹽淡水混合是極其復雜的科學問題。

目前，國內外關于河口鹽淡水混合的研究成果主要聚焦在機理和鹽水入侵及其工程應用上。在鹽淡水混合機理方面，徑潮流的此消彼長作用是影響鹽淡水平面分布、垂向混合的主要因素[4-5]，此外地形、工程及其氣象條件等因素對鹽淡水混合也有較大影響[6]；鹽水入侵研究表明，除外界因素外，鹽淡水混合后導致的水體密度的變化也會對水流流態(tài)產(chǎn)生影響[7]，不同河口的影響程度不同；鹽水入侵工程應用目前關注較多的是河口水源地飲用水安全問題，例如針對長江口青草沙水庫，眾多學者在分析研究海平面上升、南水北調、極端天氣、人類活動、河勢變化等對青草沙水庫取水產(chǎn)生影響的基礎上，提出了利用三峽水庫補淡壓咸、北支建閘等工程措施[8-9]，并建立了長江口鹽水入侵預報系統(tǒng)[10]，為水源地取水安全提供預警、預報。

本文在前人研究成果的基礎上，分析椒江河口鹽水入侵特征。椒江位于錢塘江南側約150 km 處（位置見圖1），是我國典型的強潮河口，河口口門處海門站多年平均潮差為4.0 m，最大潮差可達6.3 m[11]。通過本文研究，一方面可以明晰強潮河口的鹽水入侵特征，另一方面可以為相關科學研究以及工程建設提供技術支撐和依據(jù)。

圖1 模型范圍與網(wǎng)格示意圖Fig.1 Schematic diagram of model domain and grid

1 鹽水入侵數(shù)學模型構建

1.1 模型范圍

本文重點研究椒江河口徑流與潮流相互作用下的鹽水上溯規(guī)律，因此首先建立平面二維鹽水入侵數(shù)學模型。數(shù)值模擬軟件采用南京水利科學研究院自主研發(fā)的CJK3D-WEM 數(shù)值模擬系統(tǒng)[12]，該軟件適用于江河湖泊、河口海岸等涉水工程中的水動力、泥沙、水質、溫排、溢油的模擬預測。笛卡爾坐標系下的水流以及鹽水運動二維方程表達如下：

式中：H為水深，單位：m；z為水位，單位：m；u、v為速度分量，單位：m/s；t為時間，單位：s；f為科氏力，單位：s-1；g為重力加速度，單位：m/s2；Nx、Ny為x、y向水流紊動粘性系數(shù)，單位：m2/s；ρ0為水體密度，單位：kg/m3；τsx、τsy為x、y向的風切應力，τbx、τby為x、y向的底切應力，單位：kg/（m·s2）；?為鹽度，單位：‰；Kx、Ky為鹽水擴散系數(shù)，單位：m2/s。

采用三角形網(wǎng)格有限體積法求解式（1），對流項數(shù)值通量采用Roe 格式的近似Riemann 解，紊動項采用單元交界面的平均值來計算通過該界面紊動粘性項的數(shù)值通量，底坡項采用斜底方法處理，方程的具體求解過程見文獻[13]。

數(shù)學模型構建以椒靈江干流和臺州灣為核心水域，北側邊界至寧波象山縣，南側邊界至樂清灣附近，東側邊界至椒江口門以東85 km 處的開闊水域，上游邊界包含椒靈江一級支流永安溪和始豐溪（見圖1）。模型上游邊界采用流量控制，外海邊界采用潮位控制。模型網(wǎng)格總數(shù)為50 572 個，椒靈江水域網(wǎng)格加密，最小網(wǎng)格邊長約為23.5 m，最大網(wǎng)格邊長為3 508 m，椒靈江干流河床地形采用2016 年實測地形，臺州灣及外海水深數(shù)據(jù)采用最新海圖水深數(shù)據(jù)，統(tǒng)一到85 高程，坐標系統(tǒng)采用2000 國家大地坐標系。數(shù)學模型計算時間步長取1 s；科氏力f=2ωsinφ，其中ω為地球自轉角速度，φ為地理緯度，取28.7°；糙率n采用附加糙率方式進行處理，計算公式為：

式中：H為水深。水流紊動粘性系數(shù)取0.01HU*（U*為摩阻流速）；動邊界水深取0.02 m；鹽水擴散系數(shù)取30HU*。

1.2 模型驗證

數(shù)學模型驗證包括水動力和鹽度兩部分。采用2016 年椒靈江河道同步水文測驗資料，包括3 個潮位站和5條垂線（見圖2），驗證時間為2016年12月8日00 時（北京時，下同）—16 日00 時，包括大潮、中潮、小潮的連續(xù)過程。圖3—4為潮位以及潮流流速流向驗證。從圖中可以看出，數(shù)學模型的潮位計算值與實測值吻合較好，西門、西岙、海門站高低潮位偏差基本在0.1 m 以內，大潮、中潮、小潮的潮流流速及流向與實測值基本吻合，流速偏差基本在10%以內，流向偏差小于15°，數(shù)學模型水動力驗證總體滿足相關規(guī)程要求。圖5 為各垂線鹽度驗證，由于實測資料鹽度采樣時間跨度（3 h）較大，因此鹽度的數(shù)學模型計算值無法與實際鹽度逐時對比，但總體來看，兩者大小及過程均吻合良好，鹽度自下游向上游逐漸減小，5#點最高鹽度值在20‰左右，3#點最高約為12‰，1#點為2.5‰左右。從水動力與鹽度的驗證結果來看，本文建立的數(shù)學模型能夠較好地模擬椒江鹽水運動特征。

圖2 潮位站及潮流垂線位置示意圖Fig.2 Schematic diagram of tide level station and tide vertical line position

圖3 水位過程驗證Fig.3 Tide level verification

圖4 潮流過程驗證Fig.4 Tidal current verification

圖5 鹽水入侵過程驗證Fig.5 Salinity intrusion verification

2 椒江河口鹽水入侵特征研究

2.1 計算條件

椒江為典型的山溪性河流，流量年內分配不均，月度差異明顯，根據(jù)保證率為50%的流域來水條件分析（見圖6），11月平均流量最小，約為23 m3/s，6 月最大，約為684 m3/s，全年平均流量約為160 m3/s，洪枯季流量差異可達30 倍以上。為反映不同徑流下的鹽水運動特征，上游徑流分別取10 m3/s、100 m3/s、200 m3/s、500 m3/s、1 000 m3/s和3 000 m3/s，外海潮位采用1.2 節(jié)驗證大潮、中潮、小潮潮位過程的結果，海門站大潮平均潮差約為5.2 m，小潮平均潮差約為2.8 m，根據(jù)2020 年海門站潮差累積頻率統(tǒng)計進行計算，得到大潮累積頻率約為5%，小潮約為87%（見圖6）。

圖6 2001年椒江徑流月度分配（a）與2020年海門站潮差累積頻率統(tǒng)計（b）Fig.6 Monthly distribution of Jiaojiang River runoff in 2001（a）and statistics of tidal range cumulative frequency of Haimen station in 2020（b）

2.2 鹽水分布特征

鹽淡水混合是徑潮流共同作用的結果。圖7與表1 給出了不同徑流、不同潮差條件下鹽水的最大上溯距離分布（指距離頭門島的距離）。從中可以看到，大潮時鹽水上溯距離明顯大于小潮，隨著徑流的增大，鹽水逐漸向下游運動，當徑流為10 m3/s時，0.1‰鹽水可以上溯至始豐溪與永安溪的交匯處三江村附近，5‰鹽水可以上溯至八仙巖附近；當徑流為200 m3/s 時，0.1‰鹽水僅能上溯至西岙附近，距離口外頭門島（見圖2中C1處）約60 km，5‰鹽水上溯至海門上游13 km 處；當徑流為500 m3/s 時，0.1‰鹽水僅能上溯至石仙婦上游約6 km 處，距離口外頭門島約40 km，5‰鹽水上溯至海門上游4 km處；當徑流為3 000 m3/s時，海門以上河段的鹽水濃度均小于0.1‰。

表1 不同條件下的鹽水最大上溯距離與鹽水鋒面比降Tab.1 Maximum upstream tracking distance of Brine and the gradient of salt water front under different conditions

圖7 不同徑流條件下的鹽水最大上溯距離Fig.7 Maximum upstream tracking distance of brine under different runoff conditions

在椒江干流沿程布置間距為5～5.5 km 的17個采樣點（位置見圖2）用于分析不同潮差、不同徑流條件下的椒江鹽水鋒面分布，結果見圖8。由圖可知，無論是小潮還是大潮，低流量條件下由于鹽水上溯距離遠，鹽水鋒面坡度較緩，當徑流量為10 m3/s時，小潮的鹽水濃度比降約為0.52‰/km，大潮約為0.51‰/km；隨著徑流量的增大，鹽水濃度比降逐漸增大（見表1），當徑流量為500 m3/s，小潮的鹽水濃度比降約為0.93‰/km，大潮約為0.73‰/km；當徑流量為1 000 m3/s，小潮的鹽水濃度比降約為0.94‰/km，大潮約為0.76‰/km。綜上可見，隨著徑流量的增大，鹽水達到最大上溯距離后，鹽水鋒面比降變化梯度呈減小趨勢。

圖8 不同徑流潮差條件下的鹽水鋒面分布Fig.8 Salt water front distribution under different runoff and tidal ranges

2.3 徑流對鹽水上溯距離變化梯度影響分析

在特定的徑流條件下，鹽水在漲落潮的作用下最大上溯位置不斷發(fā)生變化，漲憩時鹽水上溯距離最遠，落憩時最近。圖9 統(tǒng)計了椒江河口鹽水在大潮和小潮不同潮型作用下的鹽水上溯距離變化梯度與徑流變化的關系，即：

圖9 鹽水上溯距離變化與徑流的關系Fig.9 Relationship between saltwater upstream tracking distance change and runoff

式中：ε 為鹽水上溯距離隨徑流的變化梯度，單位：m/（m3·s）；?Q為不同徑流Q1、Q2的差值，單位：m3/s；?d為徑流Q1、Q2下鹽水上溯距離d1、d2之間的差值，單位：m。

由圖9可以看到，無論在漲憩還是落憩時刻，隨著徑流量的增大，鹽水上溯距離呈減速降低趨勢，漲憩時刻鹽水上溯距離變化梯度略大于落憩，小潮鹽水上溯距離變化梯度大于大潮。在10～100 m3/s的流量區(qū)間，漲憩時刻大潮鹽水上溯距離變化梯度為222 m/（m3·s），小潮為238 m/（m3·s），落憩時刻兩者數(shù)值分別為207 m/（m3·s）和233 m/（m3·s），而在200～500 m3/s流量區(qū)間，鹽水上溯距離變化梯度急劇降低至15～35 m/（m3·s），這說明在徑流抑制鹽水上溯的過程中，鹽水鋒面距離河口越遠，其受徑流的影響越大，反之，受徑流影響相對較弱。

2.4 潮差對鹽水上溯變化區(qū)間影響分析

根據(jù)數(shù)學模型計算結果，統(tǒng)計了在特定徑流條件下的小潮、中潮和大潮（對應潮型25 h 內的最大潮差分別為3.2 m、5.1 m、6.2 m）鹽水上溯距離變化區(qū)間（漲憩鹽水上溯距離-落憩鹽水上溯距離），結果見圖10。由圖可知，隨著潮差的增大，鹽水上溯變化區(qū)間總體呈線性增長趨勢，隨著徑流量的增大，鹽水上溯變化區(qū)間逐漸減小。當徑流為10 m3/s時，小潮鹽水上溯變化區(qū)間約為12 km，大潮約為21 km，而徑流為500 m3/s 時，大小潮的鹽水上溯變化區(qū)間明顯下降，小潮約為8 km，大潮約為14 km，說明鹽水鋒面距離河口越遠，潮差對鹽水上溯變化區(qū)間影響越大，反之，鹽水鋒面距離河口越近，潮差對鹽水上溯變化區(qū)間影響越小。

圖10 潮差對鹽水上溯區(qū)間影響Fig.10 Influence of tidal range on salt water upstream tracking region

2.5 椒江河口鹽水上溯距離函數(shù)關系式

自然條件下，鹽水上溯距離主要與徑流和潮差有關。根據(jù)表1中徑流與鹽水上溯距離的關系分析可知，在潮差不變的條件下，鹽水上溯距離與徑流的自然對數(shù)（lnQ）具有一定線性關系，因此用lnQ代表徑流因子，同時在徑流不變的條件下，鹽水上溯距離與潮差的線性關系較好（見圖10）。由以上分析可以構建椒江河口鹽水最大上溯距離與徑流、潮差的函數(shù)關系式：

式中：d為鹽水上溯距離（即距離頭門島的距離，單位：km）；Q為徑流，單位：m3/s；?h為潮差，單位：m；a、b、c為系數(shù)，由于鹽水上溯距離與徑流呈反比，與潮差呈正比，因此a＜0，b＞0。

根據(jù)數(shù)學模型計算不同徑流、不同潮差下漲憩和落憩各18 組鹽水上溯的最大距離，利用式（4）進行擬合，求解a、b、c。經(jīng)擬合計算，漲憩時a= -9.95，b= 2.18，c= 89.67；落憩時a= -8.08，b= 0.34，c= 75.96。圖11為鹽水上溯距離函數(shù)關系式預測值與數(shù)學模擬值的對比結果，經(jīng)統(tǒng)計，漲憩時刻兩者偏差為2%～3%，落憩時刻為5%～6%。漲憩和落憩時刻數(shù)學模擬值與預測值的相關性高達0.99，說明本文建立的椒江河口鹽水上溯距離函數(shù)關系式預測精度相對較高，可為椒江鹽水相關研究提供借鑒。

圖11 鹽水上溯距離函數(shù)關系式預測值與數(shù)學模擬值對比Fig.11 Comparison between predicted value and digital analog value of salt water upstream tracking distance function

3 結語

以強潮河口——椒江河口為研究對象，通過數(shù)學模型研究了不同徑流、潮差條件下的鹽水入侵特征。主要結論如下：

①椒江河口是典型的強潮河口，當徑流為10 m3/s 時，0.1‰鹽水可以上溯至始豐溪與永安溪的交匯處三江村附近，5‰鹽水可以上溯至八仙巖附近，當徑流為500 m3/s 時，0.1‰鹽水僅能上溯至石仙婦上游約6 km 處，5‰鹽水上溯至海門上游4 km 處，當徑流為3 000 m3/s時，海門以上椒江河段的鹽水濃度均小于0.1‰。

②在徑流變化條件下，鹽水上溯距離變化梯度漲憩大于落憩，小潮大于大潮，隨著徑流量的增大，鹽水上溯距離變化梯度逐漸增大；在潮差變化條件下，鹽水上溯距離變化區(qū)間呈線性增長趨勢，隨著徑流量的增大，同等潮差下的鹽水上溯距離變化區(qū)間逐漸減小。

③通過數(shù)學模擬結果擬合出了椒江河口鹽水最大上溯距離與徑流、潮差的函數(shù)關系式，漲憩時刻關系式預測值偏差為2%～3%，落憩時刻為5%～6%，預測值精度較高，可為椒江鹽水相關研究提供借鑒。