萬艷
(福建省水產設計院,福建 福州 350003)
由于半封閉海灣的水體滯留時間較長,水交換和凈化能力較弱,人類活動頻繁,大部分海灣出現了嚴重的富營養化、水質惡化和生態系統退化等問題。盡管許多國家和地方政府已經開展了相應的環保行動,但評估其對水環境變化的有效性需要長期的監測技術和數據支持[1]。
從不同角度可以對海水交換率進行不同定義。PARKER等[2]定義了漲潮時流入灣內的海水中含第一次進入灣內的外海水的比率為海水交換率。柏井誠[3]則從流出考慮,即落潮時流出水含第一次流出灣外的灣內水的比率為海水交換率。早期人們通過實測數據計算海水交換率。近年來由于數值計算的興起,帶動國內外學者對一些海灣作了潮流和污染物擴散過程的數值模擬,研究方法主要包括箱式模型、標識質點數值跟蹤法、對流-擴散水交換模型和三維潮流場及物質輸送模型等。潘偉然[4]利用現場觀測數據,采用單箱模型和二維數值模型分別計算了湄洲灣海水的平均交換率、平均半更換期和各區段海水的半更換期;LIN等[5]利用二維淺水動力有限元模式(Shallow water Hydrodynamic Finite Element Model,SHYFEM)計算了三沙灣的海水交換率;YUAN等[6]利用有限體積海岸海洋模式(Finite-Volume Coastal Ocean Model,FVCOM)計算了土地復墾和跨灣大橋建設等人為干預行為如何影響膠州灣的水交換能力,結果表明土地復墾會顯著增加灣內水的停留時間,而跨灣大橋建設對水交換能力的影響十分有限。
不同海域的流場、地形和邊界條件的差異,導致不同海域水交換能力不同。福建東山八尺門海堤的建成,切斷了東山灣和詔安灣之間的水體交換,使海灣水體交換受阻,淤積嚴重,加上各類污染物的排放,使八尺門及其兩側水域的水環境惡化,嚴重影響了當地人民的生產和生活。八尺門海堤的打通,能否使東山灣和詔安灣的海水互相暢通,加快東山灣的水體交換頻率?使用數值模擬評估八尺門海堤開口后對周圍海域水動力和沖淤的影響是評判八尺門海堤是否應開口的前期研究工作之一。
東山灣是我國東南沿海一個典型的亞熱帶潮汐河口灣,水域面積為230 km2,大部分海域的水深小于5 m,灣口處東門嶼附近海域的水深達到或超過15 m[7]。灣內入海河流主要為漳江,流域面積為820 km2,是東山灣泥沙的主要來源[8]。東山歷年平均風速為7.1 m/s,冬季以東北風為主,夏季以西南風為主,最大風速達40 m/s。詔安灣水域面積為153 km2,與東山灣僅隔一座八尺門海堤,2/3海域的水深小于5 m,灣口水深為5~10 m[9]。前人通過數據統計分析和數值模擬等方式對東山灣進行了研究。鄭斌鑫等[10]獲得了東山灣深槽潮流和余流特征,漲潮流速最大可達0.83 m/s,落潮流速最大可達1.52 m/s,表層余流朝向灣口,底層余流朝向灣頂。李振云等[11]利用位于東山灣及附近海域5個潮位觀測站的實測數據,分析了東山灣的潮汐性質,結果表明,東山灣為正規半日潮,最大潮差為4.32 m。秦曉等[12]通過數值模擬計算了東山灣的自凈能力,結果表明,東山灣大潮和小潮的半交換周期分別為8.13 d和13.0 d,灣口水交換能力最強且向灣頂逐漸變弱。梁群峰等[8]結合東山灣的海圖和實測水深資料,分析了東山灣內近50 a的海底沖淤變化特征。
區域海洋模式(Regional Ocean Modeling System,ROMS)[13]近年來被廣泛應用于海灣、近岸及大洋等區域的研究。此前,郭民權等[14]通過ROMS建立了九龍江河口的泥沙輸運模型,很好地模擬了洪水期間九龍江河口的泥沙沖淤變化;張世民等[15]利用ROMS模擬并預測了東碇傾倒區疏浚懸沙輸移規律及海床沖淤變化;邢傳璽等[16]利用ROMS建立了普蘭店灣的潮流模型,并計算了普蘭店灣的水交換時間。
本文通過ROMS模式分析了八尺門開口前后的水動力環境和泥沙沖淤變化情況,模式方案選擇海堤開口423 m且在海堤附近疏浚1.39 km2。本文主要對水域的流速、流態、余流、海水半交換周期和泥沙沖淤的變化進行對比分析,探究八尺門海堤對東山灣和詔安灣的水動力和泥沙沖淤的影響。
本文采用的ROMS模型由美國羅格斯大學海洋與海岸科學研究所和加利福尼亞大學洛杉磯分校共同研發。該模型使用海洋原始方程模型進行并行計算,在近年的海洋科學研究中得到廣泛應用。模型在水平方向上采用曲線正交的坐標系,垂向上采用地形跟隨的S坐標[17],在淺海和復雜地形區域有更密集的垂直分層,能較好地模擬不規則的大陸架地形。此外,ROMS內含生態、海冰、泥沙和波浪等模塊,可將多模塊結合用于海洋的預報[18]。
ROMS模式可以定義粘性和非粘性兩種類型的泥沙,通過泥沙的中值粒徑、密度和沉降速率等參數可以區分砂、粉砂和粘土等多種泥沙種類。模式還需要初始化河床底沙的特性,包括底沙的層數、每層的厚度、空隙度以及各種泥沙的分配比例[19]。
ROMS中的懸沙輸運模塊可以表示為:

式中:u、v和Ω分別表示水平方向(x和y)和垂直方向(s)某一層的流速,x、y是水平坐標系,s是垂向地形擬合坐標;Hz是計算層的厚度;------c′w′表示湍流擴散引起的泥沙在垂向上的平均通量;vθ表示背景泥沙垂向擴散系數,取5×10-6m2/s;C代表泥沙濃度,Csource,m表示泥沙的源匯項。對于懸沙,需要考慮懸沙的沉降過程(匯)以及河床泥沙的再懸浮過程(源),源匯項表示為:

式中:ws,m表示泥沙的沉降速度;Es,m表示底部的侵蝕通量;m表示不同的泥沙種類。侵蝕通量通過ARIATHURAI等[20]的參數化方法實現,即:

式中:E0,m是河床侵蝕強度(單位:kg/(m2·s));φ是表層泥沙的孔隙度;τsf和τce,m分別表示床面剪切應力和臨界侵蝕應力。模型同時考慮了泥沙濃度對水動力的影響,考慮泥沙濃度后的水體密度為:

式中:ρwater表示水體的密度,它是溫度、鹽度和壓強的函數;Nsed表示所定義的泥沙種類;ρs,m和Cm分別表示第m種泥沙的密度和濃度。
本文的計算區域為福建省詔安灣、東山灣及相鄰海域(23.39°~24.03°N,117.09°~117.83°E,見圖1),網格數為1 200×1 200。本次計算采用曲線正并網格,模型在各海域的水平分辨率不一,八尺門海堤附近約為30 m,東山灣和詔安灣海域約為60 m,外海約為80 m。模型垂向分20層。表層風應力、短波輻射和蒸發降雨等驅動場由歐洲中期天氣預報中心(European Centre for Medium-range Weather Forecasts,ECMWF)大氣再分析產品的6 h間隔數據插值得到,初始場和開邊界的水位、溫度、鹽度和流速采用法國墨卡托海洋中心(Mercator Ocean)的業務化產品插值得到。在漳江河口地區需考慮相應月份漳江的流量以及含沙量,漳江流量取5月的多年平均值,約為40 m3/s。邊界條件的水平混合采用SMAGORINSKY[21]的混合方案,垂向混合采用MELLOR等[22]提出的湍流封閉模型。開邊界由業務化運行的臺灣海峽三維數值模型[17-18]所產生的環流水位疊加波賽冬全球潮汐模型(TOPEX/Poseidon global tidal model,TPXO7)[23]的16個分潮(2N2、J1、K1、K2、L2、M1、M2、MU2、N2、NU2、O1、OO1、

圖1 模型網格區域及地形圖Fig.1 Model grid domain and topography
P1、Q1、S2和T2)的潮汐調和常數計算得出。模型考慮了干濕網格,臨界水深設置為0.2 m。泥沙模塊的參數配置方法見文獻[14]。
根據福建省水產研究所提供的《福建省東山八尺門海堤工程水文動力測驗報告》粒度分析所得的粒徑譜分布特征,本文將實際泥沙看作是由兩種特征粒徑的泥沙組份構成。這兩種組份分別為:特征粒徑為3 μm的組份,代表粒徑小于4 μm的粘土,體積占比約為40%;特征粒徑為8 μm的組份,代表粒徑在4~16 μm之間的粉沙,體積占比約為60%。無論是泥沙濃度分布還是底質沖淤變化,模型最終的輸出結果都是對上述兩種粒徑組份進行加和運算后的結果。
模型初始化于2020年5月1日00時(世界時,下同),計算至2020年6月1日00時,分別計算了八尺門海堤開口前后情景。
為了驗證模式的準確性,選取了2020年5月24—25日東山灣和詔安灣不同站點的潮流數據以及含沙量數據,將其與模式結果進行驗證,結果見圖2—5。D1和D2位于東山灣,Z1和Z2位于詔安灣(見圖1),圖中表層、中層和底層分別為0.2、0.6和0.8倍的當地水深。
從整體來看,東山灣和沼安灣潮流流速和流向的觀測值和計算值有良好的一致性(見圖2和圖3)。計算流速和觀測流速總體的均方根誤差不超過0.36 m/s,表層、中層和底層誤差最小的都是D1,誤差分別為0.14 m/s、0.14 m/s和0.19 m/s,誤差最大的都是Z2,誤差分別為0.36 m/s、0.30 m/s和0.32 m/s。D1和Z1表現較優,究其原因是D1和Z1處于深槽航道中,流速較大。盡管在潮流驗證結果上,個別站點還存在一定的偏差,但從整體來看,流速和流向的觀測值和計算值的幅值相差很小,形態也基本吻合,模型能較為準確地模擬出東山灣和詔安灣的潮流場特征,具有較高的模擬精度。

圖2 東山灣潮流驗證圖Fig.2 Tidal current verification of Dongshan Bay

圖3 詔安灣潮流驗證圖Fig.3 Tidal current verification of Zhao'an Bay
圖4為D1和D2的懸沙濃度驗證效果,模型計算值整體比觀測值低,但基本能反映東山灣懸沙濃度隨時間的變化以及表層濃度低的趨勢。尤其在D1的中層和底層,模型能較好地模擬出懸沙濃度的峰值時刻。

圖4 東山灣懸沙濃度驗證圖Fig.4 The suspended sediment concentration verification of Dongshan Bay
根據模型計算結果,以5月24日13—18時作為落潮過程,以5月24日18時—25日01時作為漲潮過程,分別將海堤開口前后在落潮過程和漲潮過程時間內的平均流速進行對比,得到漲落潮流速變化百分比(見圖5)。

圖5 漲落潮平均流速變化百分比圖Fig.5 Percentage change of average velocity of high tide and low tide
落潮時(見圖5a、c),東山灣流速增加,深槽處流速小幅增加0.37%~1.92%;詔安灣流速減少,深槽處流速小幅減小1%~20%,兩岸部分流速增大了10%~20%;八尺門水道變化較大,西側水道流速大幅減小約53%,東側水道流速大幅增加100%以上。漲潮時(見圖5b、d)與落潮時總體相似,東山灣流速增加,但深槽處流速小幅減少1%~7%;詔安灣流速減小,深槽處流速小幅增加1%~17%;八尺門西側水道流速大幅減小約87%,東側水道大幅增加100%以上,流速增大區域相比落潮時更加西移。以上變化主要是八尺門海堤開口后,原在海堤處的匯潮線向西移動所致。
從八尺門水道的漲落急流態來看,落急時(見圖6a),水道內的分潮線從原來的海堤處移至東山島西北側處,水道東側潮流強度變大,西側潮流強度變小;漲急時(見圖6b),水道的東西側潮流匯合處從原來的海堤處西移約4~5 km至東山島西北側,水道東側的潮流強度變大,西側潮流強度變小。

圖6 八尺門水道漲落急垂向平均流態Fig.6 Vertical average flow pattern at high and low tide in Bachimen channel
總體來說,開口前后八尺門水道流速的變化非常劇烈,東側流速大幅度增加而西側流速大幅度減小,這是由于開口后原本位于八尺門海堤的匯潮線和分潮線西移約4~5 km。
余流場是采用模型計算出的5月11—31日的平均余流。由圖7可以看出,八尺門海堤對詔安灣和東山灣的余流場基本沒有影響,兩灣均表現出河口余流的主要特征,即表層流向灣外(見圖7a、b),底層流向灣頂(見圖7c、d),這和陳金泉等[24]的觀測結果一致。兩灣灣口及東山灣頂的漳江口區域的余流較強,達30 cm/s,分別由外海環流和河口徑流造成。兩灣深槽區域的余流強度可達10 cm/s,有利于深槽維持,有助于船舶航行。東山灣其余區域和詔安灣灣內的余流較弱。灣外外海環流的主要方向為從西到東,這與東山灣外海的北向環流一致,余流流速可達30 cm/s;底層海流略偏左(見圖7c、d),這與底部艾克曼層理論一致[25]。

圖7 開口前后余流場(白色箭頭表示流向)Fig.7 Residual flow field before and after the opening of channel(white arrow represents the flow direction)
在八尺門水道內,開口后在原海堤處出現明顯的從水道東側到西側的余流(圖略),說明水道貫通后,水流總體從東山灣進入詔安灣,表層余流流速可達10 cm/s。該海域夏季盛行西南風,冬、秋季盛行東北風,5月為春末季風轉換期。余流的動力機制及其他因素(如風和外海環流)對余流的影響是我們未來要進一步研究的內容。
使用ROMS的被動示蹤物模塊計算海水的半交換周期,以AB、CD和EF為界定義東山灣和詔安灣的邊界,以GH和MN為界定義八尺門水道的邊界。2020年5月12日在各灣內設置被動示蹤物,初始化濃度為1,計算時間超過7 d,經歷了中潮和小潮的過程。本文將各灣示蹤物平均濃度經過12 h滑動平均后變為0.5所耗費的時間定義為海水半交換周期。
模式分別計算了東山灣、詔安灣、八尺門水道以及將東山灣和詔安灣當作一個整體的水體半交換周期。由表1可知,開口前,東山灣半交換周期為1.54 d,詔安灣半交換周期為9.60 d,東山灣和詔安灣總體半交換周期為3.00 d,八尺門水道半交換周期為1.92 d;開口后,東山灣半交換周期為1.50 d,詔安灣半交換周期為7.50 d,東山灣和詔安灣總體半交換周期為2.96 d,八尺門水道半交換周期為1.25 d。

表1 八尺門海堤開口前后海水半交換周期(單位:d)Tab.1 Seawater semi-exchange cycle before and after the opening of Bachimen seawall(unit:d)
開口后,由于兩個灣內的海水均多了一個交換口,兩個灣的半交換周期有不同程度的縮短,東山灣的半交換周期縮短了0.04 d,詔安灣的半交換周期縮短了2.1 d。東山灣的半交換周期較詔安灣明顯減小,是由于東山灣灣口的水深較深,且受外海環流以及漳江徑流影響,從東山灣內交換出去的海水更易被外海水以及徑流稀釋。八尺門水道的半交換周期縮短了0.67 d,說明兩個灣的海水可通過開口的八尺門水道進行水交換,從而改善八尺門附近的水動力環境。
泥沙沖淤采用2.2節的懸沙輸運模塊進行計算。東山灣和詔安灣的沖淤程度大體相似(見圖8),兩灣主要為淤積區域,淤積量為1~4 cm/a;兩灣的灣口以及漳江口深槽的沖刷程度較強,沖刷量為3~10 cm/a。開口前后沖淤程度的區別主要體現在八尺門水道內。開口前(見圖8c),越靠近八尺門堤壩的區域淤積量越大,東側水道和西側水道的大部分區域為淤積區域,淤積量為1~4 cm/a,在接近詔安灣頂的區域有一定的沖刷程度,沖刷量為1~5 cm/a。開口后(見圖8d),八尺門東側水道由淤積區域轉變為沖刷區域,沖刷量約為2~10 cm/a,西側水道由沖刷區域轉變為淤積區,淤積程度為1~3 cm/a。這是因為八尺門水道貫通后,水道東側的流速增大,沖刷變強,該處的泥沙被帶走,而水道西側流速減小,泥沙沉降于此。

圖8 沖淤厚度分布圖Fig.8 Souring and silting thickness distribution

圖8 (續)Fig.8(Continued)
本文通過ROMS模式建立了東山灣—詔安灣海域三維水動力-泥沙耦合模式,研究了八尺門海堤開口前后水動力環境和泥沙沖淤的變化。結果表明:當拆除八尺門海堤后,八尺門水道成為連通東山灣和詔安灣的水道,但因該水道比較窄,拆除后對水動力的影響主要在八尺門水道及鄰近海區,水道附近的水動力有明顯改善,而對東山灣和詔安灣整體的影響較小。海堤開口后,原本處于海堤位置的匯潮線西移,且有穩定的西向余流存在,兩個灣的海水可通過八尺門水道進行交換,海水自凈能力提高,但量值均較小。局部清淤對東山灣及詔安灣的改變影響較小。
總之,八尺門海堤已經建成60多年,其附近地形及岸線均有較大變動,但水動力已與地形達到了
新的平衡。單純的局部清淤和海堤開口雖然會對海堤兩側水域的水動力條件有較大改善,但對東山灣和詔安灣水動力的改進有限。
致謝:廈門大學東山太古海洋觀測與實驗站提供研究數據,林文欣、張思明幫助調試模型及畫圖。