尹 超 , 黃海軍 王道儒, 謝 琳
(1. 中國(guó)科學(xué)院海洋研究所, 山東 青島 266071; 2. 中國(guó)科學(xué)院大學(xué), 北京 100049; 3. 海南省海洋與漁業(yè)科學(xué)院, 海南 海口 571126)
臺(tái)風(fēng)是對(duì)海南影響較大的氣象災(zāi)害, 特別是水產(chǎn)養(yǎng)殖行業(yè)。海南臨高縣后水灣網(wǎng)箱養(yǎng)殖密度和產(chǎn)值均較高, 是臨高縣支柱產(chǎn)業(yè), 占全縣GDP一半以上。2011年, 遭受臺(tái)風(fēng)“納沙”(NESAT)和“尼格” (NALGAE)襲擊, 受損2050口網(wǎng)箱, 占全省的78.8%, 直接經(jīng)濟(jì)損失7.63億元; 2014年, 遭受臺(tái)風(fēng)“威馬遜”(RAMMASUN)及“海鷗”(KALMAEGI)襲擊, 受損1 383口網(wǎng)箱, 占全省的45%。到了2018年, 海南省開始推進(jìn)漁業(yè)產(chǎn)業(yè)結(jié)構(gòu)調(diào)整, 在陵水縣新建3個(gè)大型智能化養(yǎng)殖網(wǎng)箱, 從水深20米以淺近岸區(qū)域拓展到50米以深的深遠(yuǎn)海區(qū)域, 從養(yǎng)殖裝備直徑20~ 30米量級(jí)拓展到近100米量級(jí)。隨著養(yǎng)殖產(chǎn)業(yè)升級(jí), 為了維護(hù)生命和財(cái)產(chǎn)安全, 需要加強(qiáng)近岸海域的海浪預(yù)報(bào)系統(tǒng)的建設(shè)以提供科學(xué)依據(jù), 提高海洋防災(zāi)減災(zāi)的能力。
影響海浪預(yù)報(bào)的最重要因素是風(fēng)場(chǎng)預(yù)報(bào)的準(zhǔn)確性[1]。目前主要的兩大氣象業(yè)務(wù)預(yù)報(bào)系統(tǒng)分別是歐洲中期天氣預(yù)報(bào)中心的全球模式ECMWF(European Centre for Medium-range Weather Forecasts), 和美國(guó)全球預(yù)報(bào)系統(tǒng)GFS(Global Forecast System)。GFS是由美國(guó)國(guó)家海洋和大氣管理局NOAA(National Oceanic and Atmospheric Administration)推出的全球數(shù)值天氣預(yù)報(bào)計(jì)算模式, 該模式每天分別于世界時(shí)00時(shí)、06時(shí)、12時(shí)和18時(shí)發(fā)布四次預(yù)報(bào)結(jié)果, 預(yù)報(bào)時(shí)效長(zhǎng)達(dá)16天。從2015年開始空間分辨率已經(jīng)提高到了0.25°, 時(shí)間分辨率隨著預(yù)報(bào)時(shí)效增加而降低, 5天內(nèi)分辨率為1小時(shí), 5天至10天分辨率為3小時(shí), 10天至16天分辨率降為12小時(shí)。
目前國(guó)內(nèi)外的業(yè)務(wù)化海浪數(shù)值預(yù)報(bào)主要基于第3代海浪模式, 廣泛應(yīng)用的主要有WAM(Wave Modeling), WAVEWATCH Ⅲ和SWAN(Simulating Waves Nearshore)等。美國(guó)氣象局采用WAVEWATCH Ⅲ, 歐洲氣象中心采用WAM, 英國(guó)氣象局采用WAM和SWAN等。我國(guó)的國(guó)家海洋環(huán)境預(yù)報(bào)中心用WAM和WAVEWATCH Ⅲ模型進(jìn)行大洋的海浪預(yù)報(bào), 用WAVEWATCH Ⅲ和SWAN嵌套模型進(jìn)行中國(guó)海的海浪預(yù)報(bào)[2]。通常, WAVEWATCH Ⅲ模式適用于深水大洋波浪計(jì)算, SWAN模式更適和近岸波浪計(jì)算, WAVEWATCH Ⅲ和SWAN模式的嵌套常用于跨尺度的波浪計(jì)算[3]。
臺(tái)風(fēng)浪往往伴有風(fēng)暴潮, 風(fēng)暴潮會(huì)引起海岸增水[4,5], 水位的變化直接影響近岸波高的大小[6], 而波浪的輻射應(yīng)力的又會(huì)引起水位的改變[7]。潮汐引起的水位變化也會(huì)影響風(fēng)暴潮和海浪的增長(zhǎng)。網(wǎng)箱養(yǎng)殖區(qū)往往設(shè)置在半封閉海灣中, 以躲避強(qiáng)風(fēng)浪的直接作用。但其幾何特征有利于風(fēng)暴潮增水, 以至于增大了臺(tái)風(fēng)浪風(fēng)險(xiǎn)。所以, 近岸的風(fēng)浪潮流動(dòng)力過程存在著強(qiáng)相互作用, 而浪潮流耦合模型適用于近岸動(dòng)力環(huán)境問題的解決[8]。Dietrich耦合了ADCIRC(Advanced circulation model)和SWAN模式研究了颶風(fēng)Gustav (2008)過程的臺(tái)風(fēng)浪和風(fēng)暴潮特征[9]; 馮興如等基于ADCIRC+SWAN模型研究了浙江和福建海域臺(tái)風(fēng)浪變化特征[10]。
本研究通過海南島近岸養(yǎng)殖區(qū)臺(tái)風(fēng)浪精細(xì)化預(yù)報(bào)系統(tǒng)的建設(shè), 提高網(wǎng)箱養(yǎng)殖環(huán)境監(jiān)測(cè)預(yù)警水平。為漁業(yè)工作者投喂作業(yè)和維護(hù)設(shè)備提供合理的時(shí)間窗口, 并降低其勞動(dòng)風(fēng)險(xiǎn); 為養(yǎng)殖企業(yè)的生產(chǎn)決策和地方政府的應(yīng)急管理服務(wù)提供科學(xué)合理的依據(jù)。
ADCIRC模式由Luettich和Westerink于1992年開發(fā)并不斷完善發(fā)展[11]。該模式采用了有限元方法, 在二維條件下應(yīng)用了深度平均正壓和淺水方程, 通過求解廣義波連續(xù)性方程(Generalized Wave Continuity Equation, GWCE), 水深和動(dòng)量方程可以表述為:
式中: U和V是垂向積分平均流速; f為科氏力參數(shù); ζ為海表面自由高度; Ps為海水自由表面大氣壓; ρ0是海水密度; η為牛頓引潮勢(shì); τs,x, τs,y為海表面應(yīng)力; H=ζ+h為水深加自由表面高度; τb,x, τb,y為底切應(yīng)力; g為重力加速度; M, D和B分別代表了側(cè)向應(yīng)力梯度、動(dòng)量耗散和斜壓梯度。為了滿足CFL條件, 時(shí)間步長(zhǎng)取為1 s。ADCIRC模式帶有干濕網(wǎng)格判斷, 最小水深設(shè)為0.1 m。
SWAN模式包含了較多的近岸波浪物理過程, 例如折繞射、淺化、白冠破碎、底摩擦、波-波相互作用等。笛卡兒坐標(biāo)下控制方程表述為:
在耦合過程中, ADCIRC和SWAN使用同一套非結(jié)構(gòu)網(wǎng)格。所以在數(shù)據(jù)傳遞過程中不需要再調(diào)用計(jì)算資源進(jìn)行插值和數(shù)據(jù)分配, 保證了計(jì)算效率和穩(wěn)定性。因?yàn)镾WAN模式的有更大的時(shí)間步長(zhǎng), SWAN和ADCIRC模式的數(shù)據(jù)交換時(shí)間間隔設(shè)定為SWAN時(shí)間步長(zhǎng)[9]。由ADCIRC模式計(jì)算的風(fēng)速、海表面水位和流速每10分鐘傳遞給SWAN模式, 而SWAN模式將輻射應(yīng)力梯度傳回給ADCIRC模式。
水深地形采用美國(guó)國(guó)家海洋大氣管理局NOAA提供的ETOPO1數(shù)據(jù)集, 該數(shù)據(jù)為柵格化數(shù)據(jù), 分辨率為1′×1′。由于ETOPO1對(duì)近岸水深的觀測(cè)誤差較大, 近岸水深數(shù)據(jù)采用中華人民共和國(guó)海軍航海保證司令部出版的海圖, 分辨率從1︰150 000加密到1︰50 000再到重點(diǎn)海灣的1︰25 000(見圖1)。
近岸養(yǎng)殖區(qū)臺(tái)風(fēng)浪耦合數(shù)值預(yù)報(bào)模型計(jì)算區(qū)域范圍是105.5°—119.7°E和 10.8°—24.7°N, 包括海南島和南海大部分的海域。由于海南島海岸島嶼眾多, 為了精確描述岸線特征, 達(dá)到精細(xì)化預(yù)報(bào)要求, 我們從Google Earth提取了整個(gè)海南島的岸線, 使沿岸區(qū)域分辨率達(dá)到了100 m。網(wǎng)格分辨率從南海深海到陸架到近岸再到重點(diǎn)海灣逐漸提高(見圖2)。模型在開邊界處加入了俄勒岡州立大學(xué)開發(fā)的TPXO模型(http: //volkov.oce.orst.edu/tides/)的8個(gè)主要天文分潮(P1、Q1、K1、O1、K2、N2、M2和S2)。
圖1 海南島近岸地形分布圖 Fig. 1 Topography around Hainan Island
圖2 模式計(jì)算域網(wǎng)格 Fig. 2 Unstructured grid of the model domain
2014年09號(hào)臺(tái)風(fēng)“威馬遜”生成于西北太平洋, 2014年7月11日由熱帶低壓升級(jí)為熱帶風(fēng)暴, 7月14日升級(jí)為強(qiáng)熱帶風(fēng)暴, 并于下午5時(shí)升級(jí)為臺(tái)風(fēng)。7月15日“威馬遜”在菲律賓阿爾拜省拉普拉普登陸, 7月16日臺(tái)風(fēng)向西北偏西移動(dòng), 進(jìn)入南海海域。17日, 進(jìn)入南海北部, 形成為強(qiáng)臺(tái)風(fēng), 其路徑見圖3。臺(tái)風(fēng)“威馬遜”于7月18日15時(shí)30分在海南省文昌市翁田鎮(zhèn)登陸, 中心最大風(fēng)力達(dá)到了60 m/s, 中心氣壓為910 hPa, 10級(jí)風(fēng)圈半徑180 km, 7級(jí)風(fēng)圈半徑300 km。于19日7時(shí)10分在廣西防城港市再次登陸, 之后強(qiáng)度逐漸減弱, 于20日5時(shí)進(jìn)入云南境內(nèi), 并減弱為熱帶低壓[13]。超強(qiáng)臺(tái)風(fēng)“威馬遜”對(duì)海南、廣東和廣西沿海地區(qū)造成的直接經(jīng)濟(jì)損失達(dá)265.5億元。
圖3 海洋觀測(cè)站(藍(lán)點(diǎn))和臺(tái)風(fēng)“威馬遜”路徑(紅色實(shí)線) Fig. 3 Oceanographic stations (blue dots) and track of the RAMMASUN cyclone (red line)
潮位和波浪觀測(cè)數(shù)據(jù)由國(guó)家海洋局南海分局??谥行恼咎峁? 站位如圖3所示。潮位數(shù)據(jù)每小時(shí)一次, 波浪數(shù)據(jù)每3小時(shí)一次, 當(dāng)臺(tái)風(fēng)強(qiáng)度增強(qiáng)則加密到1小時(shí)一次。海口站波浪數(shù)據(jù)因?yàn)閮x器原因, 在接近波高最大值時(shí)停止工作, 缺失數(shù)據(jù)。
強(qiáng)迫場(chǎng)方面使用CFSv2(Climate Forecast System Version 2)再分析數(shù)據(jù), 包括海面10 m風(fēng)場(chǎng)和海表面氣壓場(chǎng)。CFSv2系統(tǒng)由美國(guó)環(huán)境預(yù)報(bào)中心NCEP(National Centers for Environmental Prediction)開發(fā), 它是全耦合模型包括了大氣、海洋、陸地和海冰模塊, 于2011年投入業(yè)務(wù)化運(yùn)行。該強(qiáng)迫場(chǎng)時(shí)間分辨率為1小時(shí), 空間分辨率為0.2°×0.2°。
模型計(jì)算時(shí)間從2014年7月9日18時(shí)至20日6時(shí)。模式的輸出變量設(shè)定含有效波高、譜峰周期、平均周期、波向, 輸出頻率為每小時(shí)一次。為檢驗(yàn)?zāi)J降慕Y(jié)果, 將計(jì)算的潮位和有效波高結(jié)果與臺(tái)風(fēng)“威馬遜”期間觀測(cè)波高進(jìn)行對(duì)比, 對(duì)比結(jié)果見圖4。由近岸觀測(cè)站潮位對(duì)比(圖4a和b)可知, 模擬的潮位在東方站和實(shí)測(cè)值基本一致, ??谡韭缘颓曳逯堤崆? 原因可能是再分析風(fēng)圈結(jié)構(gòu)和實(shí)際風(fēng)圈有差別。由近岸觀測(cè)站有效波高對(duì)比圖(圖4c和d) 可知, 有效波高增長(zhǎng)趨勢(shì)一致性較好, 最高值基本在同一時(shí)刻。海口站和東方站模擬有效波高值均略大于實(shí)測(cè)值。
2018年7月19日在海南島近岸生成了熱帶風(fēng)暴, 其過程持續(xù)到7月25日。圖5給出的是GFS預(yù)報(bào)的熱帶風(fēng)暴過程的風(fēng)速場(chǎng)和氣壓場(chǎng)分布圖, 預(yù)報(bào)時(shí)間均為2018年7月22日18時(shí)(UTC)。左圖為2018年7月21日0時(shí)(UTC)發(fā)布的預(yù)報(bào)數(shù)據(jù), 該組數(shù)據(jù)的預(yù)報(bào)時(shí)效約48小時(shí)。右圖為2018年7月22日0時(shí)(UTC)的發(fā)布數(shù)據(jù), 該組數(shù)據(jù)的預(yù)報(bào)時(shí)效約24小時(shí)。圖5a 和c顯示, 熱帶風(fēng)暴中心大體位于瓊州海峽中部, 海南島東部風(fēng)速較大, 能達(dá)到18 m/s。圖5b和d顯示, 預(yù)報(bào)的熱帶風(fēng)暴中心向西偏移, 大體位于海南島西北部海域, 海南島南部風(fēng)速較大, 達(dá)到了16 m/s, 比21日0時(shí)發(fā)布的預(yù)報(bào)數(shù)據(jù)極值偏小。
圖4 1409號(hào)臺(tái)風(fēng)“威馬遜”期間沿岸各站潮位及有效波高觀測(cè)值和模擬值 Fig. 4 Computed results of sea surface elevations and significant wave heights (red lines) vs. observations (black dots) during Typhoon RAMMASUN
圖5 GFS預(yù)報(bào)的2018年7月22日18時(shí)(UTC)海南島附近風(fēng)速場(chǎng)和氣壓場(chǎng) Fig. 5 Wind and surface air pressure fields predicted by GFS at 18: 00 UTC on July 22, 2018 near Hainan Island
圖6顯示的是模擬的2018年7月22日18時(shí)(UTC)海南島附近波浪場(chǎng), 圖6a采用2018年7月21日0時(shí)(UTC)發(fā)布的預(yù)報(bào)風(fēng)速場(chǎng)和氣壓場(chǎng), 而圖6b采用2018年7月22日0時(shí)(UTC)發(fā)布的預(yù)報(bào)風(fēng)速場(chǎng)和氣壓場(chǎng)。有效波高的分布與圖5風(fēng)速?gòu)?qiáng)度分布基本吻合。圖6a波浪高值區(qū)域在海南島東部海域, 最大值達(dá)到了5 m。圖6b高值區(qū)域在海南島南部海域, 相比于圖6a最大值減少了大概1 m。
將48小時(shí)、24小時(shí)預(yù)報(bào)值和海南陵水灣波浪觀測(cè)數(shù)據(jù)進(jìn)行對(duì)比(見圖7)。時(shí)間范圍選取7月20日16時(shí)-24日2時(shí), 其中紅線代表48小時(shí)預(yù)報(bào)值、藍(lán)線代表24小時(shí)預(yù)報(bào)值、黑色點(diǎn)代表實(shí)測(cè)值。48小時(shí)預(yù)報(bào)最大有效波高為3.2 m左右, 24小時(shí)預(yù)報(bào)最大有效波高為2.15 m左右, 二者在近岸的最大有效波高差值接近1.0 m。在7月23日4時(shí)之前, 48小時(shí)預(yù)報(bào)的有效波高較24小時(shí)預(yù)報(bào)值偏大, 之后較24小時(shí)預(yù)報(bào)值偏小且和實(shí)測(cè)值更為吻合。48小時(shí)預(yù)報(bào)結(jié)果存在極值偏高的現(xiàn)象, 且波高增長(zhǎng)時(shí)間及區(qū)間提前。24小時(shí)預(yù)報(bào)有效波高在波浪增長(zhǎng)、整體趨勢(shì)上和實(shí)測(cè)值更為吻合, 7月23日9時(shí)之后的數(shù)值較實(shí)測(cè)值偏大。對(duì)比圖5可以看出, 48小時(shí)預(yù)報(bào)最大風(fēng)速中心區(qū)域在海南島東海岸, 陵水波浪觀測(cè)站受其影響較大, 風(fēng)速比24小時(shí)預(yù)報(bào)值偏大, 導(dǎo)致波浪計(jì)算結(jié)果也較大。
圖6 模擬的2018年7月22日18時(shí)(UTC)海南島附近波浪場(chǎng) Fig. 6 Simulated significant wave heights at 18: 00 UTC on July 22, 2018 near Hainan Island
圖7 熱帶風(fēng)暴期間48小時(shí)和24小時(shí)預(yù)報(bào)的有效波高 Fig. 7 The 48 h and 24 h forecasts of significant wave heights during the tropical storm
最后對(duì)有效波高進(jìn)行定量檢驗(yàn)評(píng)估, 評(píng)估參數(shù)包括平均絕對(duì)誤差(MAE)、平均相對(duì)誤差(MRE)、均方根誤差(RMSE)和皮爾遜相關(guān)系數(shù)(R)[14-16]。其中, 平均絕對(duì)誤差、平均相對(duì)誤差和均方根誤差體現(xiàn)了觀測(cè)值同實(shí)測(cè)值的平均偏離情況, 皮爾遜相關(guān)系數(shù)則能反映觀測(cè)值同實(shí)測(cè)值變化趨勢(shì)的密切程度。檢驗(yàn)評(píng)估的計(jì)算公式如下所示:
有效波高的預(yù)報(bào)值與觀測(cè)值的檢驗(yàn)統(tǒng)計(jì)結(jié)果如表1所示。24小時(shí)和48小時(shí)預(yù)報(bào)的MAE均小于0.3 m, 24小時(shí)的MAE為0.22 m, 相對(duì)于48小時(shí)預(yù)報(bào)值提高了23%。48小時(shí)預(yù)報(bào)值MRE為20.75%, 而24小時(shí)預(yù)報(bào)值MRE為17.0%。48 小時(shí)預(yù)報(bào)值RMSE為0.439 m, 24小時(shí)預(yù)報(bào)值的RMSE提高到了0.272 m, 提高幅度較大。24小時(shí)和48小時(shí)的相關(guān)系數(shù)R均大于0.8, 24小時(shí)預(yù)報(bào)值略有提高達(dá)到了0.87。四項(xiàng)統(tǒng)計(jì)指標(biāo)24小時(shí)預(yù)報(bào)值均好于48小時(shí)的結(jié)果, 其中RMSE和MAE指標(biāo)提高較大。
表1 預(yù)報(bào)值與觀測(cè)值的統(tǒng)計(jì)檢驗(yàn)結(jié)果 Tab. 1 Statistical analysis of forecast results and observations
本文根據(jù)近岸網(wǎng)箱養(yǎng)殖區(qū)臺(tái)風(fēng)浪精細(xì)化預(yù)報(bào)的工作需求, 基于非結(jié)構(gòu)網(wǎng)格的ADCIRC+SWAN耦合模型, 采用GFS風(fēng)場(chǎng)和氣壓場(chǎng)數(shù)據(jù), 建立了適用于海南近岸的海浪預(yù)報(bào)系統(tǒng)。利用近岸觀測(cè)資料對(duì)2014年7月“威馬遜”臺(tái)風(fēng)過程進(jìn)行后報(bào)驗(yàn)證, 在此基礎(chǔ)上對(duì)2018年7月海南島近岸海域熱帶風(fēng)暴過程進(jìn)行海浪預(yù)報(bào)及結(jié)果統(tǒng)計(jì)分析檢驗(yàn)。主要結(jié)論如下:
(1) 采用高分辨率CFSv2風(fēng)場(chǎng)和氣壓場(chǎng)作為模式輸入, 在最大波高與波浪增長(zhǎng)衰減過程兩個(gè)方面, 計(jì)算的結(jié)果和近岸觀測(cè)一致性較好。模擬的有效波高值略大于實(shí)測(cè)值。
(2) GFS預(yù)報(bào)的風(fēng)場(chǎng)與氣壓場(chǎng)每6小時(shí)不斷更新, 熱帶氣旋的路徑、中心氣壓和中心強(qiáng)度均不斷變化。預(yù)報(bào)時(shí)效對(duì)海浪預(yù)報(bào)有直接影響。
(3) 熱帶風(fēng)暴過程的預(yù)報(bào)結(jié)果對(duì)比發(fā)現(xiàn), 48小時(shí)預(yù)報(bào)的最大波高比實(shí)測(cè)偏大較多, 24小時(shí)預(yù)報(bào)的波浪增長(zhǎng)衰減過程和實(shí)測(cè)值更為吻合。24小時(shí)預(yù)報(bào)的平均相對(duì)誤差達(dá)到了17%, 均方根誤差為0.272 m, 24小時(shí)和48小時(shí)預(yù)報(bào)的相關(guān)系數(shù)均大于0.8。
通過對(duì)預(yù)報(bào)結(jié)果的統(tǒng)計(jì)檢驗(yàn), 該模型的預(yù)報(bào)精度可以滿足近岸養(yǎng)殖區(qū)臺(tái)風(fēng)浪預(yù)報(bào)業(yè)務(wù)的需求。可以為養(yǎng)殖企業(yè)的生產(chǎn)決策和地方政府的應(yīng)急管理服務(wù)提供科學(xué)合理的依據(jù)。