象山港余流特征數(shù)值研究

萬猛，姚炎明，陳琴，李丹

（浙江大學(xué) 港口海岸與近海工程研究所，浙江 杭州 310058）

狹長(zhǎng)型海灣（形態(tài)系數(shù)＜0.5） 占我國海灣總量的1/4 以上（吳桑云，2000）。象山港地處浙江北部沿海，背面緊靠杭州灣，南鄰三門灣，東側(cè)為舟山群島，是一個(gè)呈東北-西南走向的狹長(zhǎng)型半封閉海灣（中國海灣志編纂委員會(huì)，1993）。由于象山港獨(dú)特的地形和水動(dòng)力特征，水體受到強(qiáng)潮、復(fù)雜地形和不規(guī)則岸線的相互影響，港內(nèi)水體運(yùn)動(dòng)具有較強(qiáng)的三維結(jié)構(gòu)。潮流等作為海洋物質(zhì)輸運(yùn)的動(dòng)力基礎(chǔ)，將對(duì)象山港內(nèi)的溶解物質(zhì)和懸浮物質(zhì)（營養(yǎng)物、泥沙、污染物等） 的輸移起重要作用。余流的量值雖小，但它指示了水體的輸移和交換情況，對(duì)灣內(nèi)物質(zhì)的長(zhǎng)期輸移、擴(kuò)散、沉積等有著密切的關(guān)系。因此，研究該海域的余流特征，將對(duì)象山港的物質(zhì)輸移規(guī)律和機(jī)制有一個(gè)整體的認(rèn)識(shí)，也將有助于解決該海域日益嚴(yán)重的海洋環(huán)境問題。

董禮先等（2000） 根據(jù)實(shí)測(cè)水文資料對(duì)象山港的鹽度分布、余環(huán)流結(jié)構(gòu)等進(jìn)行了分析，發(fā)現(xiàn)象山港余流結(jié)構(gòu)包括重力環(huán)流和水平環(huán)流，不同區(qū)域余環(huán)流的斷面結(jié)構(gòu)取決于灣內(nèi)重力環(huán)流和狹灣外水平環(huán)流二者的強(qiáng)弱對(duì)比。朱軍政（2009） 應(yīng)用FVCOM 模式建立了象山港海域三維正壓潮流數(shù)值模型，根據(jù)數(shù)值計(jì)算結(jié)果簡(jiǎn)要分析了潮流主導(dǎo)下象山港余流場(chǎng)的三維特性，結(jié)果表明，灣內(nèi)主槽為上部余流凈流向?yàn)惩?、下部余流凈流向?yàn)硟?nèi)的重力循環(huán)水流結(jié)構(gòu)，在島嶼周圍和地形變化劇烈處有多處渦流，且表層余流大于底層余流。解靜等（2012）基于FVCOM 建立了象山港海域的水動(dòng)力數(shù)值模型，在此基礎(chǔ)上建立了污染物擴(kuò)散模型，分析了潮余流對(duì)象山港海域污染物輸移擴(kuò)散的影響。數(shù)值模擬是研究余流及物質(zhì)輸運(yùn)的重要手段之一，但由于考慮因素比較單一，且通常缺乏實(shí)測(cè)資料驗(yàn)證，數(shù)值模擬精度明顯不足。隨著對(duì)河口、海岸地區(qū)的研究不斷深入，考慮因素隨之增加且更接近實(shí)際情況，依次出現(xiàn)了正壓模型（El-Shabrawy et al，2012）、斜壓模型（Guarnieri et al，2013）、環(huán)流模型（楊陽等，2007）、紊流閉合模型（張卓等，2007） 等。象山港東臨舟山群島，口門外島嶼眾多，岸線曲折，流態(tài)復(fù)雜，同時(shí)受長(zhǎng)江口、錢塘江徑流影響較大。因此，本文選取了包括象山港和舟山海域在內(nèi)的大范圍計(jì)算區(qū)域，基于Delft3D-Flow模塊建立了象山港海域的三維斜壓潮流模型，開展了該海域的潮（余） 流特征研究，對(duì)歐拉余流場(chǎng)的時(shí)空分布進(jìn)行了定量的比較和分析，并與實(shí)測(cè)資料計(jì)算分析得到的歐拉余流分布對(duì)比，很好的模擬了象山港的余流特征。此外，本文還比較了斜壓模式與正壓模式下歐拉余流場(chǎng)的差異，并簡(jiǎn)單分析了余流結(jié)構(gòu)的控制機(jī)制。

1 數(shù)學(xué)模型

模型沿垂向采用σ 坐標(biāo)進(jìn)行分層。定義

模型控制方程如下：

連續(xù)性方程：

ξ 和η 方向的動(dòng)量方程：

鹽度方程：

式中：ζ 表示水位，m；d 表示水深，m；z 表示直角坐標(biāo)系下的垂向坐標(biāo)，m；H=ζ+d 表示總水深，m；表示直角坐標(biāo)系（x，y） 與正交曲線坐標(biāo)系（ξ，η） 的轉(zhuǎn)換系數(shù)；f 為科氏力參數(shù)；Fξ和Fη表示ξ 和η 方向上水平雷諾應(yīng)力的不平衡項(xiàng)；Mξ和Mη分別為ξ和η 方向上的動(dòng)量源或匯；ρ0為水體密度，kg/m3；Vv 為垂向渦動(dòng)系數(shù)，m2/s；DH，Dv 分別表示水平和垂向擴(kuò)散系數(shù)，m2/s；S 表示鹽度，PPT；λd表示一階降解系數(shù)；S0表示源匯項(xiàng)；σC0為普朗特-施密特?cái)?shù)；Pξ和Pη分別為ξ 和η 方向上的壓力梯度項(xiàng)。

為了分析河口分層的程度，F(xiàn)isher 等（2013） 提出了一個(gè)二階無量綱數(shù)——理查德森（Richarson）數(shù)，計(jì)算公式如下：

式中：Qf為河道流量，m3/s；b 為實(shí)測(cè)斷面寬度，m；△ρ 為表、底層的密度差，kg/m3；ρ 為平均密度，kg/m3；ut為根均方流速，m/s。若＞0.8，河口高度分層；0.08 ＜＜0.8，河口部分混合；＜0.08，則河口垂向高度混合。本文根據(jù)實(shí)測(cè)資料，選取口門處一實(shí)測(cè)斷面，經(jīng)過一個(gè)潮周期平均計(jì)算得到理查德森數(shù)=1.69。因此，本文選取斜壓模式來模擬象山港海域的水動(dòng)力場(chǎng)。

動(dòng)量方程中的壓力梯度分為水位梯度的正壓項(xiàng)和密度變化引起的斜壓項(xiàng)，Pξ和Pη計(jì)算公式如下：

式中：右側(cè)第1 項(xiàng)為水位梯度引起的正壓項(xiàng)；第2 項(xiàng)為密度變化引起的斜壓項(xiàng)；海水密度是鹽度和水溫的函數(shù)，由Eckart 經(jīng)驗(yàn)公式給定：

當(dāng)0 ＜t ＜40 ℃，0 ＜S ＜40 時(shí)，

2 模型應(yīng)用

象山港附近海域島嶼眾多，來流受到外海和長(zhǎng)江口、錢塘江等影響較大，流態(tài)復(fù)雜。因此，模型選取了包括象山港海域、舟山海域在內(nèi)的大范圍計(jì)算區(qū)域（圖1）。對(duì)計(jì)算區(qū)域采取正交曲線網(wǎng)格離散，并對(duì)象山港海域進(jìn)行了加密處理，網(wǎng)格數(shù)為601×369，ξ，η 方向水平分辨率約為90 m×90 m，外海的分辨率約為200 m×200 m。垂向分6 層，表層至底層所占比例分別為10%，20%，20%，20%，20%，10%。模型計(jì)算時(shí)間步長(zhǎng)取15 s。

模型初始海流取為靜止，水位為零。初始鹽度場(chǎng)根據(jù)2011年夏季象山港海域60 個(gè)大面站實(shí)測(cè)鹽度資料插值給定。開邊界采取水位控制，考慮到模型計(jì)算范圍較大，將開邊界分為4 段，每條開邊界給定端點(diǎn)上的水位值，中間點(diǎn)采用線性插值的方法計(jì)算得出，各端點(diǎn)的水位值由歷史數(shù)據(jù)經(jīng)調(diào)和分析計(jì)算給定。開邊界鹽度則由2011年夏季象山港海域連續(xù)站大小潮期實(shí)測(cè)鹽度值給定。閉邊界采取自由滑移邊界條件，與閉邊界垂直方向流速為零。

圖1 模型范圍和計(jì)算網(wǎng)格

考慮到模型計(jì)算時(shí)間較短，且計(jì)算區(qū)域水深較淺，溫度隨時(shí)間的變化較小，水體溫度沿垂向變化不明顯，為計(jì)算方便，斜壓模型中假設(shè)計(jì)算區(qū)域內(nèi)水體溫度恒為15 ℃。由于缺少實(shí)測(cè)徑流資料，根據(jù)《海灣志》記載，象山港年平均徑流量約為41m3/s，徑流對(duì)象山港海域水動(dòng)力場(chǎng)影響相對(duì)較小。因此，本文選取象山港年平均徑流量來考慮徑流的影響。

3 計(jì)算結(jié)果驗(yàn)證及分析

模型計(jì)算時(shí)間為2011年7月24日0 ∶00-8月1日23 ∶00（包含一個(gè)大、小潮期），共216 h。驗(yàn)證資料為2011年夏季象山港海域?qū)崪y(cè)水文資料，本次測(cè)量范圍包括佛渡水道、牛鼻水道在內(nèi)的象山港海域，按照要求布設(shè)了15 個(gè)測(cè)流點(diǎn)和6 個(gè)潮位站，因篇幅所限，本文選取了其中的測(cè)流站、潮位站各2 個(gè)用來模型驗(yàn)證（圖2）。

3.1 模型驗(yàn)證

圖4 給出了C1、C2 兩個(gè)測(cè)站的潮位計(jì)算值與實(shí)測(cè)值的驗(yàn)證結(jié)果。計(jì)算潮位與實(shí)測(cè)潮位大小、相位符合較好。夏季小潮期和大潮期的流速、流向驗(yàn)證結(jié)果見圖5，夏季小潮期時(shí)間為2011年7月24日8 ∶00-7月25日9 ∶00，大潮期時(shí)間為2011年7月31日15 ∶00-8月1日16 ∶00。各站計(jì)算流速、流向與實(shí)測(cè)值基本符合，誤差基本控制在15 %以內(nèi)，滿足精度要求。

3.2 象山港的流場(chǎng)特征

象山港的漲、落潮流流向與岸線基本平行，呈明顯的往復(fù)流性質(zhì)。受到地形、邊界的影響，象山港漲落潮歷時(shí)明顯不對(duì)稱，漲潮歷時(shí)大于落潮歷時(shí)，漲潮流速則小于落潮流速。圖6 給出了象山港夏季大潮期漲、落急兩個(gè)特征時(shí)刻的垂向平均流場(chǎng)分布。可以看到，流速由灣外到灣內(nèi)逐漸減小，灣中部西滬港口附近由于過水?dāng)嗝婵s窄，流速較大。漲急時(shí)刻水體經(jīng)牛鼻水道和佛渡水道進(jìn)入灣內(nèi)，水體沿深槽向內(nèi)運(yùn)動(dòng)，至西滬港口處分出一支流進(jìn)入西滬港，主流則繼續(xù)沿深槽向?yàn)稠斖七M(jìn)。水體到達(dá)烏沙山上游時(shí)，受到兩側(cè)島嶼的阻擋，在強(qiáng)蛟附近分為兩支，分別進(jìn)入鐵港和黃墩港。落潮時(shí)則正好相反，水體由鐵港和黃墩港流出，在強(qiáng)蛟處交匯并沿深槽向下游運(yùn)動(dòng)，行經(jīng)西滬港口處與西滬港流出的支流匯合后沿牛鼻水道和佛渡水道流向外海。經(jīng)分析發(fā)現(xiàn)，大、小潮期漲落潮特性基本一致。

圖2 測(cè)量海域及測(cè)站分布

圖3 潮位驗(yàn)證

圖4 垂向平均流速、流向驗(yàn)證

圖5 垂向平均流場(chǎng)分布

3.3 歐拉余流

余流是指從實(shí)測(cè)海流中除去周期性流（潮流等） 后的水體流動(dòng)，也即是由潮流調(diào)和分析得到的非周期性部分。實(shí)際測(cè)量的海流資料多為一個(gè)周日連續(xù)觀測(cè)資料，因此余流中實(shí)際上包含周期大于1 d的潮流。在流體力學(xué)中，描述流體運(yùn)動(dòng)有拉格朗日法和歐拉法兩種，因此有對(duì)應(yīng)的拉格朗日余流和歐拉余流。拉格朗日法以流體質(zhì)點(diǎn)為描述對(duì)象，描述位置隨時(shí)間的變化規(guī)律，如溢油擴(kuò)散、臺(tái)風(fēng)路徑預(yù)報(bào)等問題；歐拉法則描述了空間固定位置處的流體運(yùn)動(dòng)隨時(shí)間的變化情況，如計(jì)算某區(qū)域的流場(chǎng)、鹽度場(chǎng)、泥沙場(chǎng)分布。

基于上述理論，定點(diǎn)觀測(cè)的海流對(duì)應(yīng)的是歐拉流場(chǎng)，所對(duì)應(yīng)的即為歐拉余流。對(duì)夏季15 個(gè)測(cè)站漲、落潮歷時(shí)統(tǒng)計(jì)，發(fā)現(xiàn)大、小潮期平均潮流周期均為12 時(shí)14 分。因此，將固定測(cè)站連續(xù)兩個(gè)半日潮周期（24.5 h） 內(nèi)時(shí)間間隔10 min 的潮流進(jìn)行調(diào)和分析，去除周期性潮流后得到的即為固定測(cè)站處的余流。該余流包括咸淡水混合導(dǎo)致的密度流、風(fēng)海流、陸架環(huán)流和沿岸流等。

歐拉余流定義如下：

式中：T 為潮周期，文中取24.5 h。對(duì)大、小潮期各得一份數(shù)據(jù)，對(duì)二者取平均作為平均意義下的余流。對(duì)數(shù)值模擬得到的相同時(shí)間內(nèi)的潮流數(shù)據(jù)用最小二乘法做調(diào)和分析，得到的非周期流動(dòng)即為模擬的歐拉余流（圖6）。

可以看出，模擬與實(shí)測(cè)得到的歐拉余流結(jié)構(gòu)基本一致。在表層，余流流速普遍較大，佛渡水道南側(cè)的2 個(gè)測(cè)站表層余流流向指向?yàn)硟?nèi)，北側(cè)1 個(gè)測(cè)站流向偏E，流速大小為10～16 cm/s，這3 個(gè)點(diǎn)的表層流場(chǎng)主要受夏季長(zhǎng)江口、錢塘江入海徑流控制；灣頂附近區(qū)域余流流向?yàn)镹E，指向?yàn)惩猓嗔飨鄬?duì)較弱，徑流是水體運(yùn)動(dòng)的主要?jiǎng)恿Γ粸持胁繛跎成?、西滬港附近區(qū)域余流指向SE，余流大小為7～16 cm/s，這是由于受到過水?dāng)嗝婵s窄、北岸島嶼眾多和西滬港匯流的影響，水體運(yùn)動(dòng)較為紊亂，島嶼附近容易形成渦流；西澤附近區(qū)域受到佛渡水道來流影響，流向偏東南向，與佛渡水道入流匯合沿牛鼻水道流出，牛鼻水道3 個(gè)測(cè)站余流均指向?yàn)惩猓笮?1～14 cm/s。在底層，余流大小明顯小于表層，佛渡水道3 個(gè)測(cè)站流向與表層基本一致，但余流值只有5～6 cm/s；口門至西滬港附近區(qū)域主槽的余流方向則與表層相反，流向偏SW，指向?yàn)硟?nèi)，大小為2～9 cm/s，產(chǎn)生這種現(xiàn)象的原因一是受到地形底摩擦的影響，水體運(yùn)動(dòng)受阻；二是垂向斜壓效應(yīng)所致，徑流、潮流和地形的共同作用形成了垂向余環(huán)流結(jié)構(gòu)；灣中部烏沙山、西滬港區(qū)域的余流方向偏NW，大小為2～4 cm/s；口門外西澤附近區(qū)域余流方向與表層一致，均指向SE，與佛渡水道來流匯合后沿牛鼻水道流向?yàn)惩?，大小?～13 cm/s，牛鼻水道左岸余流受附近島嶼影響，流態(tài)較為復(fù)雜，余流也較大，這是受到長(zhǎng)江口、錢塘江夏季徑流的影響，水體沿佛渡水道、牛鼻水道形成一個(gè)沿岸流水交換通道。

對(duì)于表、底層，象山港歐拉余流場(chǎng)在島嶼附近表現(xiàn)出明顯的“多渦”結(jié)構(gòu)，在梅山附近產(chǎn)生了一個(gè)順時(shí)針環(huán)島渦；牛鼻水道靠近六橫島一側(cè)則產(chǎn)生了一個(gè)逆時(shí)針渦，這不利于物質(zhì)向?yàn)惩廨斶\(yùn)；灣中部烏沙山附近區(qū)域島嶼眾多，余流較為復(fù)雜，存在著很多變動(dòng)的順時(shí)針或逆時(shí)針的渦流；西滬港內(nèi)水深較淺，灘涂面積大，口門狹窄，余流較弱，與外界水體交換所需時(shí)間較長(zhǎng)，不利于物質(zhì)輸運(yùn)。

3.4 正壓與斜壓模式下的歐拉余流

為探討斜壓效應(yīng)在象山港歐拉余流場(chǎng)中的作用，本文建立了定解條件與前面斜壓模式相同的正壓潮流數(shù)值模型，并根據(jù)計(jì)算結(jié)果得到了正壓模式下的歐拉余流場(chǎng)分布（圖8）。通過比較，兩種模式計(jì)算得到的余流大小相差不大，但方向有所變化。在表層，正壓模式下的歐拉余流由灣頂流出，與佛渡水道入流匯合后，沿牛鼻水道流向外海，余流大小為4～28 cm/s，與實(shí)測(cè)點(diǎn)得到的歐拉余流較為符合；但在底層，正壓模式下口門附近的歐拉余流方向與斜壓模式下口門附近的歐拉余流方向相反，流向與表層基本一致，指向?yàn)惩?，但流速有所減小。與實(shí)測(cè)點(diǎn)資料得到的余流相比，口門至西滬港附件區(qū)域差異較明顯。造成這種差異的主要原因是狹灣內(nèi)外水體控制機(jī)制的不同。與正壓模式相比，斜壓模式考慮了計(jì)算區(qū)域內(nèi)垂向鹽度不均勻性引起的斜壓密度流對(duì)流場(chǎng)的影響，反映了鹽度場(chǎng)和流場(chǎng)的耦合作用。在灣頂附近，水體運(yùn)動(dòng)主要受徑流控制，斜壓動(dòng)力對(duì)歐拉余流的影響較小，兩種模式下得到的歐拉余流相差不大，余流方向均指向?yàn)惩?，水體運(yùn)動(dòng)以水平環(huán)流為主。而在口門至西滬港附近區(qū)域，水體受到上游徑流、潮流、地形的共同影響，垂向分層明顯，斜壓效應(yīng)增強(qiáng)，表層余流指向?yàn)惩猓讓佑嗔鲃t指向?yàn)硟?nèi)，水體運(yùn)動(dòng)為水平環(huán)流和垂向環(huán)流共同作用。因而斜壓模式更準(zhǔn)確模擬出了水體層化和垂向余環(huán)流結(jié)構(gòu)。

圖7 正壓模式下歐拉余流場(chǎng)分布

4 結(jié)語

象山港是一個(gè)強(qiáng)潮半日潮海灣，灣內(nèi)水體流動(dòng)具有較強(qiáng)的三維結(jié)構(gòu)。水平方向采用正交曲線坐標(biāo)，垂向利用σ 坐標(biāo)變換，基于Delft3D-Flow 模塊建立了象山港三維斜壓潮流數(shù)值模型，較好的模擬了象山港海域潮流場(chǎng)的時(shí)空分布特征。

由灣外至灣頂，歐拉余流呈減小的趨勢(shì)，灣外余流最大為33 cm/s，灣頂最小，潮致余流平均僅5 cm/s。垂向上，表層余流大于底層余流，表層余流指向?yàn)惩?，而底層余流則指向?yàn)硟?nèi)。潮流、徑流、復(fù)雜地形是影響象山港歐拉余流場(chǎng)的重要因素。

對(duì)比斜壓、正壓兩種模式下的歐拉余流分布，可以看出灣頂附近主要受徑流控制，水體運(yùn)動(dòng)以水平環(huán)流為主，口門至西滬港區(qū)域受徑流和潮流的共同影響，斜壓模式更準(zhǔn)確模擬出了水體層化和垂向余環(huán)流結(jié)構(gòu)。

模型計(jì)算得到了斜壓模式下象山港歐拉余流場(chǎng)的分布，并與15 個(gè)測(cè)站的實(shí)測(cè)水文資料經(jīng)調(diào)和分析得到的實(shí)測(cè)余流對(duì)比，符合良好，因而可以比較精確的復(fù)演象山港海域的水動(dòng)力環(huán)境，為今后進(jìn)一步模擬海區(qū)物質(zhì)輸移、擴(kuò)散、沉積提供了良好的水動(dòng)力場(chǎng)基礎(chǔ)。

El-Shabrawy M M,Fasseih K M,Zaki M A,2012.A Barotropic Model of the Red Sea Circulation.ISRN Civil Engineering,2012.

Guarnieri A, Pinardi N, Oddo P, et al, 2013. Impact of tides in a baroclinic circulation model of the Adriatic Sea. Journal of Geophysical Research:Oceans,118 (1) :166-183.

Guarnieri A, Pinardi N, Oddo P, et al, 2013. Impact of tides in a baroclinic circulation model of the Adriatic Sea. Journal of Geophysical Research:Oceans,118 (1) :166-183.

董禮先，蘇紀(jì)蘭，2000.象山港鹽度分布和水體混合:I.鹽度分布和環(huán)流結(jié)構(gòu).海洋與湖沼，31（2）：151-158.

解靜，梁書秀，2012.象山港污染物漂移擴(kuò)散規(guī)律研究.水道港口，33（5）：429-435.

吳桑云，王文海，2000.海灣分類系統(tǒng)研究. 海洋學(xué)報(bào)，22 (4) :83-89.

楊陽，周偉東，2007.大洋環(huán)流的診斷計(jì)算.海洋通報(bào)，26 (6)：3-8.

張卓，宋志堯，2007.近海區(qū)域二階紊流封閉模型的比較研究.海洋通報(bào)，29（1）：12-21.

中國海灣志編纂委員會(huì)，1993.中國海灣志:第五分冊(cè).北京：海洋出版社，166-167.

朱軍政，2009.象山港三維潮流特性的數(shù)值模擬.水力發(fā)電學(xué)報(bào)，28（3）：145-151.