基于MIKE3的豐枯水期防城灣水交換能力研究*

周爭橋，夏 維,邱文博,李冠宇

(1.交通運輸部南海航海保障中心廣州海事測繪中心，廣東廣州 510320; 2.國家海洋技術(shù)中心，天津 300111)

0 引言

防城灣地處北回歸線以南的低緯度地區(qū)，屬于亞熱帶海洋性季風氣候，冬季溫和，夏季多雨，季風明顯，受災害性天氣影響較明顯。防城灣分為東灣和西灣，兩灣被陸地隔開，通過狹長水道連通，灣區(qū)被主城環(huán)抱，人口密集，工業(yè)發(fā)達，每年有大量的污水排入海灣。因防城河的存在，目前污水口多放置在西灣。水交換能力是評價海灣環(huán)境容量的重要指標，表征了海灣自身的物理自凈能力，海灣的水質(zhì)情況與水交換能力的強弱直接相關(guān)[1]。在徑流和潮動力等條件下，開展防城港東、西灣海域水交換能力的研究，形成海灣交換能力的科學認識，以期制訂合理的對策，降低對海灣生態(tài)環(huán)境的影響，對防城港東、西灣未來的合理開發(fā)利用有著重要的科學意義。

箱式模型[2-3]是最早的水交換能力的研究方法，隨著數(shù)值計算模型的發(fā)展，對流-擴散模型[4]和Lagrange質(zhì)點追蹤[5]逐漸成為研究海灣水交換的主要方法。高抒等[2]基于狹長海灣多箱物理模型研究了象山港的水交換機制，該模型按垂向充分混合時間小于1 d和超過1 d兩種情況分別運行，結(jié)果均表明象山港是一個水體交換緩慢的海灣。趙亮等[6]基于ECOM水動力模型，采用Lagrange質(zhì)點追蹤法研究了膠州灣的整體和分區(qū)域水交換能力，指出深水區(qū)水交換時間為7 d，部分海區(qū)水交換時間最長可達2個月。李希彬等[7]應用三維水動力和水交換數(shù)值模型，得出渤海灣水體半交換周期為323 d，海灣西北部海域和南部海域的水交換率較低。陳振華等[1]基于 POM 水動力模式，考慮漫灘過程，得出欽州灣的整體平均水體半交換時間約為18 d。朱金龍等[8]利用二維對流-擴散模型，研究了岸線變遷影響下的芝罘灣水交換能力變化過程。曹雪峰等[9]基于無結(jié)構(gòu)三角網(wǎng)格的FVCOM海洋模式，在考慮徑流影響環(huán)境下，建立了欽州灣海域水動力數(shù)值模型，指出內(nèi)灣茅尾海的水體交換主要受制于潮汐和徑流，而外灣欽州灣的水體交換受制于潮汐和波浪。陶磊等[10]建立了基于歐拉觀點的渤海區(qū)域水齡模型，通過對比潮汐單獨作用下與潮汐和季風共同作用下渤海灣內(nèi)水齡的時空分布特征，研究了這兩個因素在渤海灣水交換中的作用，指出季風是促進渤海灣水交換的重要因素。目前，對流-擴散模型和Lagrange質(zhì)點追蹤是研究水交換時被廣泛接受和采用的方法，王聰?shù)萚11]利用Lagrange質(zhì)點追蹤和ECOM模式研究了大亞灣的水交換能力，分析了兩種方法結(jié)果差異產(chǎn)生的原因。

針對防城灣水交換的研究并不多，陳憲云等[12]采用傳統(tǒng)方法，利用形態(tài)因子面積及潮差因子計算防城港東灣2008年和2012年的納潮量，對比2008年，2012年防城港東灣納潮減少量占總納潮量的5.5%，導致該海域海水交換能力減弱。蔣磊明等[13]從納潮量、潮流場、余流和波浪4個因素分析防城港鋼鐵項目的圍填海工程對附近海域水動力環(huán)境的影響，認為圍填海工程所引起的海流、波浪、納潮量等因素的變化均較小，納潮量最大減少8%。上述研究為本研究提供了良好的研究基礎(chǔ)，為了對防城灣水交換能力有全面的認識，本文基于無結(jié)構(gòu)的三角網(wǎng)格建立了高分辨率的 MIKE3三維水動力數(shù)值模型，采用Lagrange質(zhì)點追蹤方法，選取水體半交換周期和水體交換律作為評價指標，開展防城灣水交換能力研究，以期為防城灣的科學管理決策提供支持。

1 材料與方法

1.1 模式介紹

采用丹麥水力學研究所的三維水動力模型MIKE3開展本研究。該模型相比MIKE21更適合深水水體的模擬研究，被廣泛用于河口、海岸、海洋的研究中[14-17]。MIKE3 Flow Model采用標準Galerkin 有限元法進行水平空間離散，在時間上，采用顯式迎風差分格式離散動量方程與輸運方程。MIKE3模型是基于雷諾平均化的Navier-Stokes方程，滿足Boussinesq假設(shè)和靜水壓力分布。

1.2 數(shù)值模式配置

(1)計算區(qū)域設(shè)置

根據(jù)北部灣海域的潮波運動規(guī)律，在確保計算區(qū)域邊界效應對防城灣無影響的前提下，進行計算范圍的確定。開邊界為越南云屯縣-中國北海市大風江口，網(wǎng)格在防城灣、防城河口和東興灣進行了加密，計算區(qū)域及網(wǎng)格見圖1。計算區(qū)東西寬約150 km，南北長約75 km，總面積約為50 km2，外海的計算網(wǎng)格步長為1—2 km，在防城河流域周圍海域網(wǎng)格最小分辨率約為10 m，整個計算區(qū)域內(nèi)水平方向上共有節(jié)點數(shù)31 879個，網(wǎng)格數(shù)59 364個。垂向上采用sigma坐標，垂向共分為4個σ層。

圖1 計算區(qū)域網(wǎng)格及水深示意圖Fig.1 Diagram of grid and water depth in calculation area

(2)水深和岸界

外海水深數(shù)據(jù)采用美國國家地學測量中心提供的Digital Bathymetric Database Version 5.2(DBDB5)原始數(shù)據(jù)集，近岸水深采用電子海圖以及測量資料，在海圖水深利用過程中，通過沿海海洋站潮位與最低理論基準面關(guān)系對海圖水深進行訂正，確保其準確性。

采用海圖中岸界，并利用海岸線勘測資料和Google Earth最新岸線信息進行訂正。

(3)計算時間步長和底床糙率

模型計算時間步長根據(jù)CFL條件進行動態(tài)調(diào)整，確保模型計算穩(wěn)定進行，保證計算的CFL數(shù)小于0.8，最小時間步長0.1 s，最大時間步長120 s。

(4)模型初始條件

網(wǎng)格節(jié)點的水位與網(wǎng)格單元的流速都取0。潮汐在開邊界的振幅、遲角和潮流流速等要素使用TPXO9數(shù)據(jù)。TPXO9是2018年最新發(fā)布的版本，在中國近海有較高的準確度，共包括8個基本分潮(M2、S2、N2、K2、K1、O1、P1、Q1)，2個長周期分潮(Mf和Mm)，以及M4、Ms4、Mn4、2N2、S1等分潮。本文在開邊界通過潮位輸入了M2、S2、N2、K2、K1、O1、P1、Q1共8個分潮。

(5)外部強迫條件設(shè)置

風速和防城河流量設(shè)置根據(jù)多年統(tǒng)計資料進行確定。根據(jù)防城河中游的長歧水文站 1956—2004 年49 a資料分析，全流域多年平均徑流量為2.105×1010m3，多年平均徑流深2 353.0 mm。豐水期(4—10月)流量可達99.85 m3/s，枯水期(11月—次年3月)流量約為豐水期的33%，流量可達33.28 m3/s。根據(jù)防城港市氣象局提供的1994—2014年風速資料，豐水期多年平均風速為2.7 m/s，枯水期多年平均風速為3.7 m/s。

(6)工況設(shè)置

主要研究了兩種條件下防城港東灣和西灣水體交換情況：1)豐水期+ SSW向多年平均風(2.7 m/s)；2)枯水期+NNE向多年平均風(3.7 m/s)。

1.3 評價指標

水交換是指水體通過對流和擴散等物理過程與周圍水體相互混合，在海灣、近岸海域的水交換研究中，一般是通過定義各種時間尺度來描述水體交換能力的強弱或交換速度快慢，并使用半交換時間、交換率和存留時間等描述其時間尺度。水體交換率是指灣外水體通過潮汐作用進出河口或海灣，每次進入或流出灣內(nèi)的外海水占整個進潮量或退潮量的比率。而水交換半周期是指在潮汐和徑流的共同作用下，灣內(nèi)水體交換出50% 到灣外時的周期數(shù)，在一定程度上表征了海灣的水交換能力，在欽州灣[1]、膠州灣[6]、丁字灣[18]等海灣均有所應用。

水體交換半周期在不同的潮周期下出現(xiàn)不同的變化與當時的潮汐現(xiàn)象、徑流作用、風的作用、溫度鹽度變化等因素有關(guān)，因此要確定防城港海域水體交換的周期不能僅根據(jù)某一兩次潮周期下的計算來確定，而應該是建立在多個潮周期的計算上，統(tǒng)計出各海域的平均水體交換率和平均水體交換半周期。

本文以海灣口門處作為分界線，選取大潮落潮中間時刻釋放保守物質(zhì)，保守物質(zhì)濃度海灣內(nèi)設(shè)置為1，灣外和徑流設(shè)置為0。由于潮流的動力作用，防城灣內(nèi)保守物質(zhì)的平均濃度必然會隨漲、落潮過程出現(xiàn)潮周期性的變化，也會出現(xiàn)大小潮的半月周期性變化，本文利用PL64TAP低通濾波器[19]濾掉潮流部分，分別計算豐枯水期防城灣東、西灣保守物質(zhì)濃度變化。本文對多個潮周期的水體交換率和水體交換半周期計算取平均值。

1.4 模型驗證

為了驗證模型的準確性，利用水動力模型模擬了該海域的潮位和潮流狀況，潮位和海流的實測站位布設(shè)見圖2，經(jīng)緯度坐標和觀測時間如表1。

圖2 潮位和潮流驗證點位置圖Fig.2 Location map of verification points of tide level and tidal current

表1 潮位和潮流驗證點坐標Table 1 Location of verification points of tide level and tidal current

2 結(jié)果與分析

2.1 模型驗證結(jié)果

(1)潮位驗證

模型輸出了仙人洞臨時潮位站位置2019年9月29日—2019年9月30日潮位，將潮位數(shù)據(jù)和模型輸出的潮位數(shù)據(jù)訂正至同一基面，繪制潮位過程曲線(圖3)。結(jié)果顯示，模型可以較好地模擬潮位過程。

圖3 仙人洞驗潮站潮位實測模擬過程曲線Fig.3 Simulation process curve of tidal level measurement at Xianrendong tidal station

(2)潮流驗證

采用4個潮流實測站位表層流速資料對模擬結(jié)果進行驗證，觀測時間為2019年9月29日—2019年9月30日。潮流驗證曲線見圖4。 驗證結(jié)果表明，對應觀測點潮流模擬結(jié)果與實測資料基本吻合，能夠較好地反映項目周邊海域潮流狀況。

圖4 各站位表層潮流驗證曲線Fig.4 Surface current verification curves in each stations

潮位和潮流模擬結(jié)果顯示，本文建立的水動力模型可以真實反映防城灣實際的水體流動特點，可作為該海域水交換研究的模型。

2.2 流場特征分析

模擬結(jié)果(圖5)表明海域潮流性質(zhì)為不正規(guī)日潮，漲潮時，海水從外海進入近海各個海灣，落潮時，海水向南退入外海，落潮流流速大于漲潮流流速，漲潮流流向為偏北向，落潮流流向則為偏南向，在岸邊的速度方向與岸線較為一致。在漲潮過程中，外海潮波傳入防城港灣，流向都指向口門，在潮波傳播過程中，受地形約束強度較大，到達漁萬島南端分為兩支，一支沿牛頭嶺向西北方向流去，另一支沿暗埠口江向東北方向流去。流速在牛頭嶺深槽附近流到最大，過了牛頭嶺和企沙半島后，潮流進入東、西灣內(nèi)灣，由于水深變淺，水流阻力增大，流速逐漸減小。無論在外灣還是內(nèi)灣，潮流都受岸線和深槽的影響，流向有順著深槽的趨勢。落潮是漲潮的逆過程，流場情況類似。在低潮時，漁萬島東北和西北面海域流速較小，較大面積的灘涂露出。西灣狹窄水道處漲落潮流速可達1 m/s以上，東灣張羅潮最大流速為0.8—0.9 m/s。

圖5 大潮漲急(左)和落急(右)時刻流場Fig.5 Flow field of high (left) and low (right) tide during spring period

2.3 水體交換率

豐枯水期防城灣東灣、西灣保守物質(zhì)濃度變化見圖6。對保守物質(zhì)隨時間變化分析可知，在潮汐徑流和季風的共同作用下，兩個海灣區(qū)域的保守物質(zhì)濃度均以較快的速度減小，西灣河口區(qū)域的保守物質(zhì)濃度迅速減小，保守物質(zhì)濃度減小速度明顯大于防城灣東灣。由于天文潮的動力作用，兩個海灣灣內(nèi)平均濃度必然會隨漲、落潮過程出現(xiàn)潮周期性的變化，也會出現(xiàn)大小潮的半月周期性變化。

圖6 東、西灣豐水期和枯水期保守物質(zhì)濃度變化曲線Fig.6 The variation curve of the concentration of the conservative substance in the East Bay and the West Bay during the wet and dry period

本文對豐枯水期不同時間段的保守粒子濃度進行了分析(圖7-8)，可知：東灣的保守物質(zhì)濃度的空間分布呈現(xiàn)由灣口向灣內(nèi)增加，主水道保守物質(zhì)濃度明顯大于兩側(cè)淺水，說明主水道的水體交換能力較強，隨著時間的增加，保守物質(zhì)濃度減小速率降低。西灣保守物質(zhì)濃度的空間分布呈現(xiàn)由河口向灣外增加，沖淡水主水道保守物質(zhì)濃度明顯大于兩側(cè)淺水，說明主水道的水體交換能力較強，隨著時間的增加，保守物質(zhì)濃度減小速率降低。

圖7 防城灣豐水期保守物質(zhì)濃度分布Fig.7 Distribution of conservative substance concentration in Fangcheng Bay during the wet period

2.4 水體交換時間

經(jīng)計算，防城灣西灣、東灣水體半交換時間如表2所示，防城灣西灣豐水期水體半交換周期約為1.7 d，枯水期水體半交換周期約為3.6 d，東灣豐水期水體半交換周期約為15.2 d，枯水期水體半交換周期約為27.8 d。從半交換周期來看，豐水期半交換周期明顯短于枯水期，徑流對污染物擴散起到加速作用，作用較季風明顯。

表2 防城灣海域水體交換情況Table 2 Water exchange situation of Fangcheng Bay

3 結(jié)論

防城灣海域水交換時間主要受防城河徑流和潮流的影響，防城河豐水期徑流量99.85 m3/s，枯水期徑流量為33.28 m3/s，數(shù)值模擬顯示西灣最大流速可達1 m/s,東灣最大流速為0.9 m/s。

因防城河的存在，西灣保守物質(zhì)濃度的空間分布呈現(xiàn)由河口向灣外增加，沖淡水主水道保守物質(zhì)濃度明顯大于兩側(cè)淺水，說明主水道的水體交換能力較強。東灣的保守物質(zhì)濃度的空間分布呈現(xiàn)由灣口向灣內(nèi)增加，主水道保守物質(zhì)濃度明顯大于兩側(cè)淺水。兩灣隨著時間的增加，保守物質(zhì)濃度減小速率降低。

西灣因潮流流速較大，加上防城河徑流的影響，其水交換時間明顯較東灣短，半交換期在豐水期為1.7 d，枯水期為3.6 d；東灣因與西灣通過狹長水道相連，在豐水期可明顯看到徑流的影響，半交換時間為15.2 d，枯水期其影響效果不明顯，半交換時間為27.8 d。

綜上，防城灣西灣水交換能力較好。在未來的防城灣使用規(guī)劃中，選擇水交換能力較強的西灣中部的東側(cè)沿岸，東灣應該選擇在中部海域。