三沙灣余環(huán)流特征數(shù)值模擬研究

，，

(福建省水產(chǎn)研究所，福建 廈門 361013)

2017-07-19

福建省省屬公益類科研院所基本科研專項(2015R1003-4).

林建偉(1981-)，男，助理研究員，碩士研究生，主要從事海灣河口水動力環(huán)境研究.Tel：13696926085.E-mail：frilinjw@126.com

游遠新(1969-)，高級工程師，主要從事海域海島開發(fā)與保護研究.E-mail：13950066060@139.com

林建偉，游遠新，劉國昕.三沙灣余環(huán)流特征數(shù)值模擬研究[J].漁業(yè)研究，2017，39(5)：331-341.

三沙灣余環(huán)流特征數(shù)值模擬研究

林建偉，游遠新*，劉國昕

(福建省水產(chǎn)研究所，福建 廈門 361013)

本文基于ROMS模型，同時考慮溫度、鹽度和流場以及海表的凈熱通量、水氣通量、動量通量和地表徑流，建立三沙灣精細化水動力數(shù)值模型。結(jié)合實測水文資料，對基于ROMS的三沙灣水動力模型進行驗證，并研究了三沙灣的潮(余)流場分布特征，對潮(余)流場的空間分布進行了描述分析。結(jié)果表明：本文建立的三沙灣水動力模型能夠較好地反演三沙灣海域的潮(余)流場分布特征；由灣外至灣頂，余流呈減小的趨勢，灣外余流最大約為12 cm/s，灣內(nèi)較小，余流平均約為1 cm/s；垂向上，表層余流大于底層余流，表層余流指向灣外，而底層余流則指向灣內(nèi)；本模型較好地再現(xiàn)了三沙灣水體層化和垂向余環(huán)流結(jié)構(gòu)。

ROMS模型；三沙灣；余流；數(shù)值模擬

三沙灣形狀似伸展的右手掌，海灣被羅源、東沖半島環(huán)抱，僅在東南方向有一個狹口——東沖口與東海相通，口門寬約3 km，是個半封閉的海灣；四周為山環(huán)繞，海岸曲折復雜，主要由基巖、臺地和人工海岸組成，岸線總長度為553 km；水域開闊，海灣總面積為784 km2，西北側(cè)有賽江、霍童溪等中小河溪注入[1]。灣內(nèi)海底地形崎嶇不平，侵蝕和堆積地形都很發(fā)育，灣中有許多可航水道、暗礁、島嶼和淺灘；三都、東安、青山等島嶼是灣內(nèi)主島；東沖水道、青山水道和金梭門水道是灣內(nèi)主航道；灣內(nèi)各小灣頂及淺水航道兩側(cè)常有淺灘和干出灘發(fā)育。灣內(nèi)最大水深達90 m，是我國的天然良港之一[2]。

三沙灣獨特的地形和水域條件，為灣內(nèi)的海洋生物生長、繁殖提供了有利的條件。20世紀80年代，三沙灣官井洋被批準為大黃魚自然保護區(qū)之后，海上水產(chǎn)養(yǎng)殖迅速發(fā)展，成為全國最大的大黃魚網(wǎng)箱養(yǎng)殖地區(qū)[3]。同時，隨著環(huán)三沙灣地區(qū)的經(jīng)濟發(fā)展，該區(qū)域30余年來的圍填?；顒右草^頻繁，至2013年，圍填?？偯娣e達16 122.99 hm2[4]。根據(jù)該地區(qū)的相關(guān)規(guī)劃，環(huán)三沙灣地區(qū)的圍填海面積還會進一步增大。

圍填?；顒訉p小灣內(nèi)的水域面積，減小海灣納潮量，影響海灣的水交換能力和海灣的水動力環(huán)境。潮流作為海洋特質(zhì)輸運的動力基礎，將對海灣內(nèi)的溶解物質(zhì)和懸浮物質(zhì)(營養(yǎng)物、泥沙、污染物等)的輸移起重要作用。余流的量值雖小，但它指示了水體的輸移和交換情況，對灣內(nèi)物質(zhì)的長期輸移、擴散、沉積等有著密切的關(guān)系。所以開展三沙灣水動力及其余流特征基礎研究，充分認識三沙灣物質(zhì)輸運的動力基礎，對于今后研究圍填?；顒訉θ碁乘畡恿八w交換能力的影響，以及對航道回淤和水質(zhì)的影響具有重要意義。目前，海洋工程通用的模型有MIKE21、FVCON、DELFT3D和SMS等，各個模型有不同的算法和適用范圍。根據(jù)本次研究的目的，本文選用目前廣泛應用的開源模式ROMS來研究三沙灣的潮余流結(jié)構(gòu)特點。

1 模型介紹

ROMS(Regional Ocean Modeling System)是近年來由Rutger University與University of California at Los Angeles共同研究開發(fā)的三維非線性的斜壓原始方程模式，該模型在海洋科學中有著廣泛的應用[5-7]。ROMS在水平方向上可采用曲線正交坐標系，并使用“Arakawa C”交錯網(wǎng)格，可以增大關(guān)注區(qū)域的分辨率；而在垂向上采用S坐標系(Stretched Terrain-following Coordinates)[8]，可以加密表、底層，提高模型在溫躍層和底邊界層的解析度，并且可以客觀反映出流場受海底地形變化的影響，減少Sigma坐標帶來的計算誤差[9]。在時間算法上采用時間分裂技術(shù)，對內(nèi)外模態(tài)獨立求解，可以提高計算效率。ROMS在邊界條件的選擇、對流項的計算，提供了多種高精度的數(shù)值算法，尤其在湍流閉合模式上提供了諸如MY2.5[10]、GLS[11]、LMD[12]等經(jīng)典算法。此外，ROMS還耦合了泥沙模塊[13]、生態(tài)模塊[14]等子模塊，并通過C語言預處理技術(shù)(C-preprocessing)激活上述不同的物理及數(shù)值選項。

1.1控制方程

ROMS模式方程是采用Boussinesq近似和流體靜力近似的NavierStokes方程，在笛卡爾坐標系下的方程形式為：

(1)

(2)

(3)

ρ=ρ(T，S，P)

(4)

(5)

(6)

(7)

(8)

上述公式中，公式(1)代表連續(xù)性方程；公式(2)和公式(3)分別代表X方向和Y方向的動量方程；公式(4)代表靜力平衡方程；公式(5)是狀態(tài)方程；公式(6)和公式(7)分別是溫度以及鹽度的物質(zhì)輸運方程；公式(8)代表了除溫鹽以外如泥沙、營養(yǎng)鹽、污染物等的物質(zhì)輸運擴散方程。其中u、v為水平方向(X和Y方向)的速度分量，w為垂直方向的速度分量。f表示科氏力參數(shù)，g表示重力加速度，ρ0表示參考密度，p表示正壓和斜壓造成的壓強。ρ表示海水的局地密度。T為溫度，S為鹽度，C為其他物質(zhì)。KM和KH是垂向粘性和擴散系數(shù)，通過求解湍流閉合模式得到。Dx、Dy為X和Y方向的動量水平粘性項，Dc是物質(zhì)水平擴散項。Fx、Fy分別是X和Y方向的外部應力項，F(xiàn)c表示物質(zhì)的源匯項。

1.2邊界條件

海表面z=δ(x，y，t)處，考慮海表風應力與熱鹽通量：

(9)

(10)

(11)

τsx、τsy為海表風應力；H為海表凈熱通量；S為海表凈鹽通量。

海底z=-h(x，y)處，受底摩擦作用，海底熱、鹽通量均為0；

(12)

(13)

(14)

τbx、τby為海底摩擦力，與海底流速的平方成正比。

而對于外海開邊界條件，由10個潮汐調(diào)和常數(shù)推算開邊界的水位和流速，運用Chapman水位邊界條件和Flather二維流速邊界條件進行計算，外海三維流速及溫鹽用強制邊界條件。河流采用點源邊界條件。

1.3數(shù)值模擬及驗證

三沙灣島嶼眾多，岸線曲折復雜，流態(tài)復雜，因此模型選取了包括三沙灣、羅源灣及其周邊海域在內(nèi)的計算區(qū)域(圖1)。對計算區(qū)域采取正交曲線網(wǎng)格離散，網(wǎng)格數(shù)為307×217，網(wǎng)格空間尺度分辨率最小為270 m×270 m，最大為390 m×250 m。垂向分25層，根據(jù)水深平均等分。圖2為模型計算網(wǎng)格圖。模型計算時間步長取30 s。

1.3.1 模型區(qū)域邊界條件

模型開邊界采取水位控制。

模型初始場及開邊界的溫、鹽、流等參數(shù)均參考廈門大學臺灣海峽業(yè)務化模式結(jié)果插值得到，同時開邊界水位加入由潮汐調(diào)和常數(shù)計算得到的潮位。

陸域邊界同時考慮現(xiàn)行的法定岸線和三沙灣現(xiàn)狀圍填海區(qū)域。

閉邊界采取自由滑移邊界條件，與閉邊界垂直法向流速為零。

海表的熱通理，水通量及動量通量均來自福建省海洋預報臺的WRF大氣模式結(jié)果。

三沙灣海域水深地形數(shù)據(jù)采用海軍航海保證部2015年版三沙灣海圖(圖號：13971)和海軍航海保證處2014年版鹽田港及白馬港海圖(圖號：13975)數(shù)字化后得到。

在三沙灣項相應的網(wǎng)格上加入了“源”，以反映賽江的年平均徑流量69.69×108m3。

1.3.2 模型驗證

本模型利用本公益項目資助的外業(yè)調(diào)查潮流數(shù)據(jù)和收集到的同步潮流數(shù)據(jù)進行模型驗證。驗證點站位如圖1所示，包括了6個潮流觀測站(1#～6#)以及2個潮位同步觀測站(T1和T2)。圖3分別為1#～6#潮流觀測站在觀測期間的垂直平均流速、流向的驗證曲線和T1、T2潮位觀測站在觀測期間的驗證曲線。由圖3可見，1#～6#觀測點的計算流速、流向過程，以及T1、T2潮位變化過程曲線與實測過程變化趨勢較為一致；模型驗證情況詳見表1和表2。由表1可見，1#～6#觀測點的漲落潮平均流速的偏差值均較小，相對誤差均小于10%，其對應的流向誤差小于10°；由表2可見，T1、T2平均高、低潮位的偏差值在±10 cm以內(nèi)；均符合技術(shù)規(guī)范[15]對潮流場的驗證要求。表明模型采用的物理參數(shù)和計算參數(shù)基本合理，計算方法可靠，能夠模擬研究區(qū)域潮波運動特性。

表1 1#～6#潮流站流速、流向驗證情況表

站位Station漲潮平均流速/m·s-1Averagevelocityinfloodtide漲潮平均流向/°Averagedirectioninfloodtide落潮平均流速/m·s-1Averagevelocityinebbtide落潮平均流向/°Averagedirectioninebbtide觀測Observed計算Computed偏差Deviation觀測Observed計算Computed偏差Deviation觀測Observed計算Computed偏差Deviation觀測Observed計算Computed偏差Deviation1#0.410.459.76%354348-60.520.531.92%165157-82#0.510.545.88%1710-70.500.5510.00%18718813#0.390.427.69%6050-100.390.427.69%233228-54#0.590.648.47%4400.530.577.55%170180105#0.370.408.11%256248-80.340.378.82%9410176#0.590.613.39%3331-20.600.635.00%2082179

表2 T1和T2站潮位驗證情況表

2 結(jié)果與討論

2.1三沙灣潮流場數(shù)值模擬結(jié)果

圖4(a)和圖4(b)分別是三沙灣及其周邊海域在現(xiàn)狀邊界條件下的漲、落潮流場分布圖。由圖4(a)可見，在漲潮過程中，外海潮流以約0.7 m/s的流速由外海繞過浮鷹島和海島鄉(xiāng)，一小部分潮流朝W向由羅源灣口門漲進羅源灣，大部分潮流則朝WN向穿過三沙灣口門，沿著灣內(nèi)主潮流通道進入三沙灣海域；漲潮流在三沙灣內(nèi)主航道(深槽)的流場較強，如在三沙灣口門，由于口門過水斷面較小，其最大漲潮流速超過1 m/s，流向為NNW向，在鹽田港主航道上的最大漲潮流速可達0.8 m/s，流向為N向，在白馬河主航道上的最大漲潮流速為1.0 m/s左右，流向NW向，在三都鎮(zhèn)北側(cè)海域的漲潮流速約為0.6 m/s，流向為W向，在東吾洋東安島的西南側(cè)海域，流速約為0.8 m/s，流向為ENE向；而在水深較淺的海域流場較弱，流速小至0.1 m/s左右，漲潮歷時約6.3 h。由圖4(b)可見，落潮流沿著漲潮流的反方向流出灣外，在白馬河的落潮流向為SE向，在鹽田港內(nèi)的落潮流向為S向，在東吾洋東安島的西南側(cè)落潮流向為SW向，在三都鎮(zhèn)北側(cè)海域的落潮流向為ESE向；落潮時，在八都鎮(zhèn)、漳灣鎮(zhèn)、三都鎮(zhèn)、沙江鎮(zhèn)附近和長春鎮(zhèn)西側(cè)等澳內(nèi)淺灘地區(qū)有大面積的灘涂露出[圖4(b)中灣內(nèi)空白區(qū)域]；落潮流以較大的流速流出灣外；這基本上反映了三沙灣海域的地形特征，落潮歷時約6.0 h。

2.2三沙灣余流場數(shù)值模擬結(jié)果

余流是指海流中剔除潮流后的運動，是一種相對比較穩(wěn)定的海水流動，對污染物質(zhì)的長期輸運影響較大。本文基于ROMS模型，綜合考慮各方面環(huán)境要素，在潮流場模擬的基礎上，對三沙灣及周邊海域的冬季余流場進行了模擬。圖5(a)至圖5(c)分別為三沙灣及周邊海域冬季表、中、底層余流場分布圖。

由圖5(a)可見，冬季三沙灣口的表層余流基本態(tài)勢為南出北進，落潮余流與漲潮余流在灣口相互頂托，落潮余流場占主導地位；從矢量場可以清晰地看出，在北壁鄉(xiāng)的東南岸邊、鑒江鎮(zhèn)與北壁鄉(xiāng)之間(雞公山南側(cè)海域)，具有形成余流渦旋的趨勢；青山島東側(cè)余流整體上朝S向；青山島以北的三沙灣深槽，具有較強的余流場，其余海域余流場較弱；在內(nèi)、外灣口上，余流強度表現(xiàn)為灣內(nèi)、外的不對稱性，強流區(qū)發(fā)生在青山島至灣口及其東南側(cè)海域，最大余流流速為12 cm/s左右，弱流區(qū)發(fā)生在整個三沙灣內(nèi)，在鹽田港和東吾洋的余流流速約為1 cm/s，而在白馬河的余流流速約為3 cm/s，這與白馬河的賽江地表徑流有關(guān)。

由圖5(b)可見，中層余流流速較表層余流小，中層余流在灣口基本態(tài)勢仍為南出北進，但落潮余流場與漲潮余流場處于相對平衡的狀態(tài)，并在北壁鄉(xiāng)的東南岸邊、鑒江鎮(zhèn)與北壁鄉(xiāng)之間(雞公山南側(cè)海域)，已形成兩個較為明顯的余流渦旋結(jié)構(gòu)；不同于表層，青山島東側(cè)中層余流整體上朝N向，青山島以北的三沙灣深槽也沒有較明顯的強余流場，整個三沙灣海域的中層余流場較弱，但白馬河的中層余流場與表層一樣，由于地表徑流的輸入導致余流場強于鹽田港和東吾洋。在內(nèi)、外灣口仍存在余流強度的不對稱性，灣口較大一般約為9 cm/s，灣內(nèi)則明顯小很多。

由圖5(c)可見，底層余流流速相比中層余流流速差別不明顯，底層余流在灣口態(tài)勢仍為南出北進，但漲潮余流場處于主導地位，并在北壁鄉(xiāng)的東南岸邊、鑒江鎮(zhèn)與北壁鄉(xiāng)之間(雞公山南側(cè)海域)，形成兩個較為明顯的余流渦旋結(jié)構(gòu)；與中層余流場相同，青山島東側(cè)底層余流整體上朝N向，青山島以北的三沙灣深槽也沒有較明顯的強余流場，整個三沙灣海域的底層余流場較弱，白馬河的底層余流場仍強于鹽田港和東吾洋。在內(nèi)、外灣口仍存在余流強度的不對稱性，灣口較大一般約為8 cm/s，灣內(nèi)則明顯小很多。

由以上各圖可見，三沙灣的余流場呈明顯的層化現(xiàn)象，表層和底層余流在灣口、青山島東側(cè)、青山島北側(cè)深槽等區(qū)域均表現(xiàn)出不同的態(tài)勢。在灣外，苔菉鎮(zhèn)附近海域，表層余流表現(xiàn)為落潮，中層余流開始出現(xiàn)轉(zhuǎn)流，底層余流則表現(xiàn)為漲潮。在灣口，表層余流整體上表現(xiàn)為落潮，而底層余流則整體上表現(xiàn)為漲潮；表層余流表現(xiàn)渦旋結(jié)構(gòu)的趨勢，而底層則清晰地表現(xiàn)出兩個余流渦旋。在灣內(nèi)，表層余流較大，底層余流較小。

3 結(jié)論

本文基于ROMS模型，在同時考慮溫度、鹽度、流場和海表的凈熱通量、水氣通量、動量通量和地表徑流，以及三沙灣及周邊海域現(xiàn)狀圍填岸線和水深地形的基礎上，建立了三沙灣的潮流模型，結(jié)合實際測流資料，模擬了三沙灣潮流和余流特征，主要結(jié)論如下：

1)基于ROMS建立的三沙灣潮流模型計算結(jié)果與實測結(jié)果吻合得較好，能夠為下一步的水體交換模擬計算提供可靠的水動力數(shù)據(jù)。

2)三沙灣潮流模型計算結(jié)果再現(xiàn)了三沙灣的潮流分布特征：在灣內(nèi)深槽通道的流速大，流向與水道走向基本一致；三沙灣地形復雜，島嶼星羅棋布，水域多呈水道形式，潮流呈往復流，流向基本與岸線平行，潮流受地形的影響較為明顯；一般落潮流速大于漲潮流速。

3)根據(jù)余流場模擬計算結(jié)果，三沙灣的余流場總體表現(xiàn)為：三沙灣外的余流場強于灣內(nèi)，這指示著三沙灣外的水體輸移和交換能力強于灣內(nèi)；三沙灣及灣口的表層余流方向整體表現(xiàn)為東南，這種形態(tài)和已報道的三沙灣冬季余流方向相吻合[2]，而底層余流則無此特點，這說明了三沙灣表層余流受風場的影響較大。

4)三沙灣口的表、中、底層余流場，表現(xiàn)出清晰的余流渦旋結(jié)構(gòu)；表層的落潮余流占主導地位，而底層的漲潮余流占主導地位，這個結(jié)果充分體現(xiàn)三沙灣口存在垂直余環(huán)流結(jié)構(gòu)。

[1]福建省人民政府.福建省海洋功能區(qū)劃2011—2020年(文本)[R].福州，2012.

[2]中國海灣志編纂委員會.中國海灣志：第七分冊[M].北京：海洋出版社，1994：54-55.

[3]蔡清海.福建三沙灣海洋生態(tài)環(huán)境研究[J].中國環(huán)境監(jiān)測，2007，23(6)：101-105.

[4]劉榮杰，張杰，馬毅.三沙灣30余年來圍填海遙感監(jiān)測與分析[J].海洋開發(fā)與管理，2014，(9)：17-21.

[5]Haidvogel D B，Arango H G，Hedstrom K，et al.Model evaluation experiments in the North Atlantic Basin：simulations in nonlinear terrain-following coordinates[J].Dynamics of Atmospheres and Oceans，2000，32(3-4)：239-281.

[6]Warner J C，Geyer W R，Lerczak J A.Numerical modeling of an estuary：A comprehensive skill assessment[J].Journal of Geophysical Research，2005，110：1-13.

[7]Wilkin J L，Arango H G，Haidvogel D B，et al.A regional ocean modeling system for the Long-term Ecosystem Observatory[J].Journal of Geophysical Research，2005，110(c6)：1-13.

[8]Ezet T，Arango H G，Schepetkin A F.Developments in terrain-following ocean models：intercomparisons of numerical aspects[J].Ocean Modeling，2002，4(3-4)：249-267.

[9]Shchepetkin A F，Williams J C.A method for computing horizontal pressure-gradient force in an oceanic model with an onaligned vertical coordinate[J].Journal of Geophysical Research，2003，108(C3)：3090.doi：10.1029/2001JC001047.

[10]Mellor G L，Yamada T.Development of a Turbulence Closure Model for Geophysical Fluid Problems[J].Reviews of Geophysics and Space Physics，1982，20(4)：851-875.

[11]Warner J C，Sherwood C R，Arango H G，et al.Performance of four turbulence closure models implemented using a generic length scale method[J].Ocean Modelling，2005，8(1-2)：81-113.

[12]Large W G，Mcwilliams J C，Doney S C.Oceanic Vertical Mixing：A preview and a model with a nonlocal boundary layer parameterization[J].Reviews of Geophysics，1994，32(4)：363-403.

[13]Warner J C，Sherwood C R，Signell R P，et al.Development of a three-dimensional，regional，coupled wave，current，and sediment-transport model[J].Computers & Geosciences，2008，34(10)：1284-1306.

[14]Fennel K，Wilkin J，Levin J，et al.Nitrogen cycling in the Middle Atlantic Bight：Results from a three-dimensional model and implications for the North Atlantic nitrogen budget[J].Global Biogeochemical Cycles，2006，20(3)：145-180.

[15]中華人民共和國交通運輸部.JTS/T 231-2—2010海岸與河口潮流泥沙模擬技術(shù)規(guī)程[S].北京：人民交通出版社，2010.

NumericalstudyofthetidalresidualcurrentsinSanshaBay

LINJianwei，YOUYuanxin*，LIUGuoxin

(Fisheries Research Institute of Fujian，Xiamen 361013，China)

Based on the ROMS model，considering the temperature，salinity，currents，net heat，moisture，momentum and runoff，this paper built the hydrodynamic numerical model in Sansha Bay.Combined with the observed data，the hydrodynamic numerical model was validated，the spatial distribution and distribution characteristics of tidal residual currents were researched and analyzed.The results showed that the model based on ROMS could reflect the distribution characteristics of tidal and residual currents in Sansha Bay well；the tidal residual currents decreased from outer to inner of the bay；on the vertical section，the surface residual currents were larger than the bottom residual currents；the surface currents flowed out of the bay，while the bottom currents flowed into the bay；this model clearly performed the water stratification and vertical structute of residual circulation in Sansha Bay.

ROMS model；Sansha Bay；tidal residual currents；numerical simulation

P731.2

A

1006-5601(2017)05-0331-11