英曉明，嚴(yán)金輝，趙明利

（國家海洋局南海規(guī)劃與環(huán)境研究院，廣東廣州 510300）

0 引言

惠州大亞灣位于粵港澳大灣區(qū)的惠州市，西鄰大鵬灣，東接紅海灣，屬于半封閉海灣（見圖1）。大亞灣灣內(nèi)潮波因受狹長海灣和島嶼地形反射及變形影響，淺水分潮效應(yīng)較大[1]。根據(jù)惠州海洋站2006—2019 年的長期驗(yàn)潮資料，該站平均潮差為0.85 m，最大潮差為1.25 m，灣內(nèi)為不正規(guī)半日潮。大亞灣海域年平均氣壓為1 012.0 hPa，多年平均氣溫為22.9 ℃，累年平均降水量為1 734.0 mm，累年平均風(fēng)速為3.4 m/s，常風(fēng)向?yàn)闁|南向，強(qiáng)風(fēng)向?yàn)闁|北向。實(shí)測資料表明，2019 年3 月純洲島與灣頂之間海域的垂向平均最大流速為0.17 ～0.18 m/s；大亞灣中部海域垂向平均最大流速為0.18～0.34 m/s。該海域鵝洲站2018 年夏季和冬季浮標(biāo)觀測資料記錄，夏季有效波高平均值為0.41 m，平均波周期為4.1 s；冬季有效波高平均值為0.23 m，平均波周期為3.5 s。大亞灣內(nèi)島嶼眾多，化工企業(yè)聚集，港口碼頭密集，是粵港澳大灣區(qū)重要的海灣。大亞灣海域海洋生物資源豐富，屬于環(huán)境敏感區(qū)，有大亞灣水產(chǎn)資源自然保護(hù)區(qū)、珊瑚、重要產(chǎn)卵場、無居民海島、紅樹林、淡澳溪河口等環(huán)境敏感保護(hù)目標(biāo)。大亞灣是許多名貴魚類的繁殖、產(chǎn)卵和育肥場所[2]。廣東省海洋功能區(qū)劃顯示大亞灣還分布著多個旅游休閑娛樂區(qū)。

圖1 研究區(qū)域和水文觀測站位置示意圖Fig.1 Location of the study area and observation station

隨著中海殼牌、?？松梨诨葜菀蚁┑缺姸嗍?xiàng)目的落戶，大亞灣過往船舶密度增大，溢油事故風(fēng)險增加。大亞灣海域一旦發(fā)生船舶溢油事故，將嚴(yán)重影響其環(huán)境敏感保護(hù)目標(biāo)并破壞海洋生態(tài)環(huán)境。1970—2018 年登陸粵港澳大灣區(qū)的熱帶氣旋共有40 個，年均為0.8 個[3]，臺風(fēng)風(fēng)暴潮發(fā)生頻率較高。2017 年臺風(fēng)“苗柏”風(fēng)暴潮期間，大亞灣多艘船舶走錨；2018 年臺風(fēng)“山竹”風(fēng)暴潮期間，大亞灣大辣甲西側(cè)海域發(fā)生船舶走錨事故，所幸均沒有發(fā)生環(huán)境污染事故。風(fēng)暴潮期間因發(fā)生船舶事故導(dǎo)致破壞海洋環(huán)境安全的風(fēng)險是潛在的，應(yīng)用數(shù)學(xué)模型研究風(fēng)暴潮期間的船舶溢油事故能夠預(yù)測污染物擴(kuò)散范圍，為溢油事故應(yīng)急處置提供參考。已有學(xué)者對大亞灣溢油和?；沸孤妒鹿曙L(fēng)險進(jìn)行了模擬研究[1-5]。然而，已有研究是針對常規(guī)天氣下環(huán)境污染事故風(fēng)險進(jìn)行模擬，而極端天氣風(fēng)暴潮期間溢油事故的模擬研究較少。本文通過數(shù)值模擬方式，對大亞灣海域風(fēng)暴潮期間的溢油事故進(jìn)行情景模擬，為溢油事故應(yīng)急決策提供技術(shù)支撐。

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

1.1 水動力模型

對于南海地區(qū)臺風(fēng)風(fēng)暴潮的數(shù)值模擬，已有學(xué)者采用數(shù)學(xué)模型MIKE21、高級環(huán)流模型（Advanced Circulation Model，ADCIRC）、Delft 3D 三維水動力-水質(zhì)模型系統(tǒng)和非結(jié)構(gòu)網(wǎng)格海洋環(huán)流與生態(tài)模型（Finite-Volume Coastal Ocean Model，F(xiàn)VCOM）等開展了風(fēng)暴增水和海浪等研究[6-11]。本文選擇MIKE21 數(shù)值模式進(jìn)行風(fēng)暴潮與天文潮耦合情況下的水動力預(yù)測模擬，模型網(wǎng)格及模擬范圍見圖2。模型共有35 835個節(jié)點(diǎn)，63 770個網(wǎng)格單元，最小網(wǎng)格邊長約為17 m。

圖2 模型網(wǎng)格和水深地形Fig.2 Model mesh and bathymetry

大亞灣海域具有顯著的潮位“雙峰”現(xiàn)象和漲潮流速“雙峰”現(xiàn)象。研究表明，當(dāng)潮波由灣口傳播至灣頂時，上述現(xiàn)象主要由M2、M4和M6分潮共同作用造成[12]，M6分潮的異常增大是引起大亞灣內(nèi)部潮位“雙峰”現(xiàn)象的主要原因，2MS6分潮對大亞灣內(nèi)部潮汐的扭曲也具有一定作用[13]。多站位實(shí)測潮位調(diào)和分析研究結(jié)果表明，淺水分潮（M4、M6、MS4和2MS6）振幅的增長是造成“雙峰”現(xiàn)象在灣內(nèi)增強(qiáng)的主要原因[14]。模型水位邊界條件采用的M2、S2、N2、K2、K1、O1、P1、Q1、M4、MS4分潮調(diào)和常數(shù)初始值引自全球潮汐模型[15]，此外增加了M6和2MS6兩個分潮，這兩個分潮的調(diào)和常數(shù)初始值取自南海潮流模型結(jié)果[16]。海面風(fēng)風(fēng)場及氣壓場采用Holland 模型。臺風(fēng)氣壓模型如下[9,17]：

式中：Pc為臺風(fēng)中心氣壓；Pn為臺風(fēng)外圍氣壓；RMW為最大風(fēng)速半徑；r為距離臺風(fēng)中心的距離；Vmax為臺風(fēng)最大風(fēng)速；lat 為臺風(fēng)中心緯度；B為Holland 參數(shù)。臺風(fēng)對水動力的影響通過海水表面摩擦系數(shù)來體現(xiàn)，摩擦系數(shù)取值范圍為0.000 5～0.001 7。

1.2 溢油模型

溢油事故發(fā)生后，油膜主要在風(fēng)和海流共同驅(qū)動下發(fā)生平移運(yùn)動、擴(kuò)展以及蒸發(fā)、乳化、溶解等風(fēng)化過程。本文采用MIKE21 HD 模塊構(gòu)建平面二維水動力模型，在此基礎(chǔ)上耦合MIKE21 OS模塊建立溢油模型。關(guān)于模型中的風(fēng)漂移系數(shù)，一些學(xué)者取值0.03～0.04[18-19]，另一些學(xué)者取值0.02～0.03[20-22]，也有學(xué)者取值0.025～0.044[23]。風(fēng)漂移系數(shù)和風(fēng)速有關(guān)，當(dāng)風(fēng)速較大時風(fēng)漂移系數(shù)取大值，因此本文取0.044。

2 模型驗(yàn)證和流場分析

本文選擇2019年3月20—21日項(xiàng)目附近海域2個潮位觀測站和3個潮流觀測站的觀測值進(jìn)行常規(guī)情況下的潮位和潮流驗(yàn)證，選擇2016年7月30日17時（北京時，下同）—8 月2 日23 時赤灣、珠海、深圳和惠州4 個驗(yàn)潮站的潮位值以及2017 年8 月22 日00 時—24 日05 時赤灣、大萬山、深圳和惠州4 個驗(yàn)潮站的潮位值進(jìn)行臺風(fēng)風(fēng)暴潮期間潮位和風(fēng)速、風(fēng)向的驗(yàn)證。觀測站位見圖1。風(fēng)暴潮期間潮位觀測值引自文獻(xiàn)[7]。

2.1 模型驗(yàn)證和流場特征分析

2019年3 月潮位和潮流模型驗(yàn)證結(jié)果分別見圖3 和圖4。由圖3 可見，潮位模擬值和觀測值基本一致，T1 站和惠州站的平均潮位誤差分別為10%和-3%，潮位驗(yàn)證結(jié)果較好，模型能夠較好地模擬大亞灣海域的潮位“雙峰”現(xiàn)象。由圖4 可見，潮流模擬值和觀測值變化趨勢大體一致，V1、V2 和V3 站的平均流速誤差分別為2%、13%和14%，平均流向誤差分別為3%、-6%和44%。潮流驗(yàn)證結(jié)果總體較好，模型能夠反映大亞灣海域的潮流運(yùn)動特征。

圖3 2019年3月20日10時—21日12時潮位驗(yàn)證結(jié)果Fig.3 Tidal level verification results from 10:00 on March 20 to 12:00 on March 21,2019

圖4 2019年3月20日10時—21日12時潮流驗(yàn)證結(jié)果Fig.4 Tidal current verification results from 10:00 on March 20 to 12:00 on March 21,2019

大亞灣海域漲急時刻和落急時刻的流場分布見圖5。無論是漲急還是落急時刻，大亞灣中部海域和范和灣灣口的流速較大。漲急時刻的平均流速約為0.23 m/s；落急時刻平均流速約為0.18 m/s。漲急時刻，潮流主要由東南向進(jìn)入大亞灣，灣內(nèi)中央主流向北；落急時刻，灣內(nèi)中央主流向南，潮流主要由東南向流出大亞灣。

圖5 大亞灣海域漲急和落急時刻流場分布Fig.5 Distribution of tide current velocity field in Daya Bay at flood tide and ebb tide

2.2 風(fēng)暴潮驗(yàn)證

本文采用1713號臺風(fēng)“天鴿”風(fēng)暴潮的風(fēng)速和風(fēng)向過程對模式進(jìn)行風(fēng)的參數(shù)驗(yàn)證，驗(yàn)證結(jié)果見圖6，由圖可見，模擬的風(fēng)速和風(fēng)向與實(shí)測值趨勢相同，整體吻合較好。1604 號臺風(fēng)“妮妲”和1713 號臺風(fēng)“天鴿”風(fēng)暴潮期間的潮位過程驗(yàn)證結(jié)果見圖7，由圖可見，在兩次臺風(fēng)風(fēng)暴潮期間，潮位模擬值和觀測值基本接近，模型能夠反映臺風(fēng)風(fēng)暴潮的潮位變化過程。

圖6 1713號臺風(fēng)“天鴿”風(fēng)暴潮期間風(fēng)速過程比較結(jié)果Fig.6 Comparison results of wind speed during storm surge of 1713 Typhoon"Hato"

圖7 1604號臺風(fēng)“妮妲”和1713號臺風(fēng)“天鴿”風(fēng)暴潮期間潮位過程比較結(jié)果Fig.7 Comparison results of tidal level during storm surge of 1604 Typhoon"Nida"and 1713 Typhoon"Hato"

3 結(jié)果分析

由于近年來大亞灣海域沒有發(fā)生典型等級以上的涉海溢油事故案例，這里僅進(jìn)行情景事故模擬。假設(shè)臺風(fēng)風(fēng)暴潮期間發(fā)生特別重大的船舶污染事故，參照《防治船舶污染海洋環(huán)境管理?xiàng)l例》，船舶發(fā)生溢油情景事故的溢油量為1 000 t。漏油時間為1 h，在模型中溢油由240 000個“油粒子”表示。溢油方式選擇固定點(diǎn)源連續(xù)性溢油，油品選擇柴油，密度取0.89×103kg/m3。

3.1 1604號臺風(fēng)“妮妲”風(fēng)暴潮溢油情景模擬結(jié)果

1604 號臺風(fēng)“妮妲”風(fēng)暴潮模擬時間為2016 年7 月29 日20 時—8 月3 日05 時。在臺風(fēng)行進(jìn)過程中，假設(shè)在8月2日凌晨（臺風(fēng)中心臨近大亞灣前3 h）大辣甲西側(cè)海域發(fā)生溢油事故。風(fēng)暴潮過程期間溢油掃海面積、殘油量和油膜中心位置見表1，風(fēng)暴潮中止后油膜掃海面積見圖8。

表1 1604號臺風(fēng)“妮妲”風(fēng)暴潮過程期間溢油掃海面積、殘油量和油膜中心位置Tab.1 Oil spill sweeping area,residual oil volume and oil film center position during the storm surge process of 1604 Typhoon "Nida"

圖8 1604號臺風(fēng)“妮妲”風(fēng)暴潮過程期間溢油29 h后油膜掃海范圍Fig.8 Scope of oil film sweeping after 29 h oil spill during the storm surge of 1604 Typhoon"Nida"

事故發(fā)生后，溢油將立即影響到附近珊瑚區(qū)、大亞灣水產(chǎn)資源自然保護(hù)區(qū)中部核心區(qū)、中部緩沖區(qū)和大辣甲無居民海島等4 個敏感目標(biāo)，2 h 后影響到大亞灣水產(chǎn)資源自然保護(hù)區(qū)南部緩沖區(qū)，3 h后影響到大亞灣水產(chǎn)資源自然保護(hù)區(qū)南部核心區(qū)，8 h后影響到小辣甲附近的珊瑚，10 h 后影響到許洲無居民海島和重要產(chǎn)卵場，11 h 后影響到大亞灣水產(chǎn)資源自然保護(hù)區(qū)北部實(shí)驗(yàn)區(qū)和沙魚洲無居民海島，12 h后影響到沙魚洲東側(cè)的珊瑚，13 h 后影響到大亞灣灣頂岸線，14 h 后影響到鵝洲和赤洲附近的珊瑚，16 h后影響到馬鞭洲附近海域，18 h后影響到寶塔洲無居民海島。臺風(fēng)“妮妲”風(fēng)暴潮中止后（溢油后29 h）油膜掃海范圍約為256.38 km2，殘油量為640.10 t。

3.2 1713號臺風(fēng)“天鴿”風(fēng)暴潮溢油情景模擬結(jié)果

1713 號臺風(fēng)“天鴿”風(fēng)暴潮模擬時間為2017 年8 月21 日02 時—24 日17 時。在臺風(fēng)行進(jìn)過程中，假設(shè)在8 月23 日06 時（臺風(fēng)中心臨近大亞灣南側(cè)）平海電廠碼頭附近海域發(fā)生溢油事故。風(fēng)暴潮過程期間溢油掃海面積、殘油量和油膜中心位置見表2，風(fēng)暴潮中止后油膜掃海面積見圖9。

表2 1713號臺風(fēng)“天鴿”風(fēng)暴潮過程期間溢油掃海面積、殘油量和油膜中心位置Tab.2 Oil spill sweeping area,residual oil volume and oil film center position during the storm surge process of 1713 Typhoon"Hato"

圖9 1713號臺風(fēng)“天鴿”風(fēng)暴潮過程期間溢油35 h后油膜掃海范圍Fig.9 Scope of oil film sweeping after 35 h oil spill during the storm surge of 1713 Typhoon"Hato"

事故發(fā)生后，溢油將在1 h 后影響到大亞灣水產(chǎn)資源自然保護(hù)區(qū)北部實(shí)驗(yàn)區(qū)，2 h后影響到大亞灣水產(chǎn)資源自然保護(hù)區(qū)中部緩沖區(qū)和桑洲無居民海島，3 h后影響到大亞灣水產(chǎn)資源自然保護(hù)區(qū)中部核心區(qū)和大辣甲無居民海島，4 h后影響到大辣甲附近的珊瑚，5 h后影響到大亞灣水產(chǎn)資源自然保護(hù)區(qū)中部緩沖區(qū)，7 h 后影響到大亞灣核電站北側(cè)岸線，8 h后影響到小辣甲附近的珊瑚，9 h后影響到小辣甲無居民海島，13 h后影響到重要產(chǎn)卵場、赤洲附近的無居民海島和珊瑚，14 h 后再次影響到大亞灣水產(chǎn)資源自然保護(hù)區(qū)北部實(shí)驗(yàn)區(qū)，15 h 后影響到沙魚洲和許洲無居民海島，17 h后影響到馬鞭洲北側(cè)的珊瑚，18 h 后影響到淡澳溪河口和鵝洲附近的珊瑚，20 h后影響到大亞灣灣頂岸線。風(fēng)暴潮中止后（溢油后35 h）油膜掃海范圍約為115.62 km2，殘油量為626.94 t。

3.3 流速變化分析

本節(jié)通過分析臺風(fēng)風(fēng)暴潮期間和無臺風(fēng)時的流速變化探討臺風(fēng)風(fēng)暴潮的影響。由圖10 可見，1604 號臺風(fēng)“妮妲”風(fēng)暴潮導(dǎo)致大亞灣內(nèi)大部分區(qū)域流速增加，尤其在灣口區(qū)域、馬鞭洲西側(cè)區(qū)域和大亞灣灣頂區(qū)域，灣口區(qū)域流速最大增加值超過0.9 m/s，馬鞭洲西側(cè)區(qū)域流速最大增加值超過0.8 m/s，大亞灣灣頂區(qū)域流速最大增加值超過0.5 m/s。1713 號臺風(fēng)“天鴿”風(fēng)暴潮導(dǎo)致大亞灣內(nèi)大部分區(qū)域流速增加，尤其在灣口區(qū)域、大亞灣核電站附近區(qū)域和大亞灣灣頂區(qū)域，灣口區(qū)域流速最大增加值超過1.0 m/s，大亞灣核電站附近區(qū)域超過0.8 m/s，大亞灣灣頂區(qū)域超過0.5 m/s。臺風(fēng)風(fēng)暴潮期間流速的增大加速了油膜的漂移和擴(kuò)散。

圖10 1604號臺風(fēng)“妮妲”和1713號臺風(fēng)“天鴿”風(fēng)暴潮過程期間和無臺風(fēng)的流場和流速變化Fig.10 Comparison of current fields and current velocity change among the storm surge process of 1604 Typhoon"Nida",1713 Typhoon"Hato"and the scenario without typhoon

4 結(jié)論和建議

大亞灣海域島嶼眾多，潮汐具有獨(dú)特的“雙峰”現(xiàn)象。本文通過MIKE21 數(shù)值模式，建立了高分辨率的惠州大亞灣附近海域的風(fēng)暴潮和溢油數(shù)學(xué)模型，通過引入M4、MS4、M6和2MS6淺水分潮，精確模擬了大亞灣潮位的“雙峰”特征。通過Holland 模型對1713號臺風(fēng)“天鴿”風(fēng)暴潮的風(fēng)場進(jìn)行了模擬，驗(yàn)證結(jié)果較好，運(yùn)用模型模擬分析了在1604 號臺風(fēng)“妮妲”和1713 號臺風(fēng)“天鴿”兩場典型臺風(fēng)風(fēng)暴潮期間，大亞灣海域發(fā)生1 000 t 特大溢油情景事故對敏感目標(biāo)的影響。

在1604號臺風(fēng)“妮妲”風(fēng)暴潮期間，假設(shè)大辣甲西側(cè)海域發(fā)生溢油事故，油膜最大擴(kuò)散范圍約為256.38 km2，溢油發(fā)生后，3 h 內(nèi)將影響附近珊瑚區(qū)、大亞灣水產(chǎn)資源自然保護(hù)區(qū)中部核心區(qū)、大辣甲無居民海島、中部緩沖區(qū)、大亞灣水產(chǎn)資源自然保護(hù)區(qū)南部緩沖區(qū)和大亞灣水產(chǎn)資源自然保護(hù)區(qū)南部核心區(qū)等敏感目標(biāo)，臺風(fēng)風(fēng)暴潮中止后，殘油為溢油量的64.0%；在1713 號臺風(fēng)“天鴿”風(fēng)暴潮期間，假設(shè)平海電廠碼頭附近海域發(fā)生溢油事故，油膜最大擴(kuò)散范圍約為115.62 km2，溢油發(fā)生后，3 h 內(nèi)將影響附近大亞灣水產(chǎn)資源自然保護(hù)區(qū)北部實(shí)驗(yàn)區(qū)、大亞灣水產(chǎn)資源自然保護(hù)區(qū)中部緩沖區(qū)、桑洲無居民海島、大亞灣水產(chǎn)資源自然保護(hù)區(qū)中部核心區(qū)和大辣甲無居民海島等敏感目標(biāo)，臺風(fēng)風(fēng)暴潮中止后，殘油為溢油量的62.7%。本文建立的溢油模型可以應(yīng)用于常規(guī)天氣和臺風(fēng)風(fēng)暴潮情況下的溢油擴(kuò)散模擬，能夠?yàn)橐缬褪鹿蕬?yīng)急決策提供技術(shù)支撐。

本文所采用的模型沒有考慮大洋環(huán)流對近岸區(qū)域的影響，同時由于大亞灣溢油事故較少，缺少驗(yàn)證資料，因此沒有深入分析溶解、擴(kuò)散等參數(shù)對溢油的影響，這是本文存在的不足，后續(xù)將根據(jù)資料情況進(jìn)一步開展研究。

對于溢油事故風(fēng)險防控，應(yīng)做好應(yīng)急預(yù)案，充分準(zhǔn)備溢油防控應(yīng)急物資，平時做好應(yīng)急演練；加強(qiáng)船舶管控，在臺風(fēng)風(fēng)暴潮來臨前做好防臺準(zhǔn)備；在臺風(fēng)風(fēng)暴潮期間若發(fā)生溢油事故，除借助衛(wèi)星遙感監(jiān)測技術(shù)外，還可以通過數(shù)值模擬及時預(yù)測溢油擴(kuò)散范圍和到達(dá)敏感目標(biāo)的時間；在臺風(fēng)風(fēng)暴潮過后，及時安排船舶、人員利用除油物資開展油膜打撈等應(yīng)急處置。