羅佳敏，姜云鵬，龐亮，馮鈺棟

( 1. 中國海洋大學 工程學院，山東 青島 266100；2. 交通運輸部天津水運工程科學研究所，天津 300456；3. 山東電力工程咨詢院有限公司，山東 濟南 250013)

1 引言

臺風是我國東南沿海常見的氣象災害，影響范圍廣、發(fā)生頻次高、變異性強、破壞力大，由臺風引發(fā)的風暴潮災害是破壞沿海地區(qū)人民生命財產(chǎn)最主要的方式之一[1]。浙江省是臺風影響頻率最高、受災最為嚴重的省份之一，2000-2016年共有46場臺風對浙江省造成影響[2]。在氣候變化的背景下，風暴潮發(fā)生的頻率不斷增加，精細化且能真實反映風暴增水變化的風暴潮數(shù)值模擬研究顯得尤為重要。

臺風是風暴增水最直接的驅(qū)動因素。風場的可靠性直接影響到風暴潮數(shù)值模擬的準確度。風暴潮數(shù)值模擬中多采用基于簡化的物理模型和數(shù)值模型的參數(shù)化臺風風場模式。1980年，Holland[3]根據(jù)澳大利亞低緯度海域臺風數(shù)據(jù)，提出了Holland模型。在最大風速半徑(Rmax)的基礎上，引入Holland參數(shù)，即徑向氣壓分布系數(shù)(B)，用來描述臺風徑向剖面結(jié)構(gòu)的變化，使模型具有更普遍的適用性[4-6]。研究發(fā)現(xiàn)，Rmax和B對臺風風場結(jié)構(gòu)及風速的影響較大[7-9]。Powell等[10]基于美國東南沿海地區(qū)的飛機探測數(shù)據(jù)對B進行了擬合研究，給出了B關于Rmax和緯度 ?的擬合關系式；Vickery和Wadhera[11]將B與Rmax關聯(lián)，提出了B與Rmax和中心氣壓差(ΔP)的擬合公式；在澳大利亞海域，Love和Murphy[12]、Hubbert等[13]、Harper和Holland[14]分別給出了B與ΔP的擬合公式；Jakobsen和Madsen[15]根據(jù)孟加拉灣海域臺風資料，提出了B與ΔP、1 min最大風速的擬合關系；在西北太平洋海域，江志輝等[8]結(jié)合《熱帶氣旋年鑒》數(shù)據(jù)，利用地球物理流體動力學實驗室（Laboratory of Geophysical Fluid Dynamics,LAGFD）的部分計算方法，得到ΔP關于Rmax的擬合曲線公式；李瑞龍[16]結(jié)合李小莉等[7]關于Rmax的計算方法，認為Rmax服從對數(shù)正態(tài)分布，同時擬合得到我國沿海臺風多發(fā)地區(qū)的Rmax計算公式；Fang等[17]利用中國氣象局（China Meteorological Administration, CMA） 最佳路徑數(shù)據(jù)集，根據(jù)1949-2016年的歷史臺風數(shù)據(jù)，給出B關于ΔP和緯度?的擬合關系式，構(gòu)建了適用于西北太平洋海域的參數(shù)化風場模型。

近年來，諸多學者針對我國沿海海域進行了臺風風場研究。唐建等[18]將CCMP（Cross Calibrated Multi-Platform）背景風場與多種經(jīng)驗風場模型結(jié)合，模擬了1105號臺風，發(fā)現(xiàn)Holland模式的模擬結(jié)果與實測數(shù)據(jù)較接近。程雨佳等[19]將聯(lián)合概率法與Holland參數(shù)風場模型結(jié)合，推算得到臺風影響下南海海域的波高重現(xiàn)期。潘冬冬等[6]發(fā)現(xiàn)在南海海域，將歐洲中期天氣預報中心（ECWMF）第五代再分析資料（ERA-5）背景風場與選用Hubbert等[13]提出的B計算方法構(gòu)建的Holland風場模型相疊加得到的臺風風場，能較好地重現(xiàn)臺風“山竹”的風速過程。方根深等[20]基于典型臺風過程建立風場模型并針對我國南部沿海進行了模型適用性研究，通過對比不同的B、Rmax計算公式，認為在西北太平洋海域，采用盧安平等[21]提出的B、Rmax組合方案進行“黑格比”臺風風場的模擬具有較高的精度。

針對影響浙江沿海的風暴潮災害，諸多學者采用不同的數(shù)值模型開展了相關研究。孫志林等[22]運用Delft-3D模型，探究了舟山漁港風暴潮位對不同臺風路徑的響應規(guī)律。丁駿等[23]采用MIKE 21水動力模型，建立了浙江沿海的二維風暴潮模型，準確模擬了1323號臺風風暴潮過程。陳潔等[24]利用MIKE 21模型，選用0608號臺風，進行海平面上升情景下的臺州市玉環(huán)市的風暴潮數(shù)值模擬，發(fā)現(xiàn)玉環(huán)市北側(cè)登陸的臺風風暴潮對海平面上升較敏感。何威等[25]選用Holland模型，對B、Rmax參數(shù)公式進行比選構(gòu)建9711號臺風模型，并使用MIKE 21模型研究了椒江河口不同形態(tài)變化對風暴潮的影響，認為河口收窄會使風暴高潮位抬升。

上述研究結(jié)果顯示，不同海域獲得的風場參數(shù)擬合方法差異較大，因此需綜合比較不同的計算方法以獲得適用于浙江海域的參數(shù)化風場模型，在此基礎上能夠?qū)Φ湫团_風誘發(fā)的風暴增水進行合理的數(shù)值分析。本文以該海域常用的Holland模型為切入點，比較風場模型中最大風速半徑Rmax和徑向氣壓分布系數(shù)B的計算方法，以0216號臺風“森拉克”、0414號臺風“云娜”為例，分析Rmax和B不同計算方法的組合對風速發(fā)展的影響，由此得到適用于浙江海域的參數(shù)化風場模型?；谏傻妮斎腼L場，利用MIKE 21模型進行該海域的臺風風暴潮數(shù)值模擬，與實測數(shù)據(jù)比較，進行天文潮位和風暴增水水位驗證與分析。

2 參數(shù)化風場模型及其關鍵參數(shù)

Holland模型[3]的徑向氣壓分布公式為

式中，Pr為距離臺風中心r處的氣壓；Pc為臺風中心氣壓；ΔP=Pa-Pc，為臺風中心氣壓差，Pa為臺風外圍氣壓， 取1 010 hPa；Rmax為臺風風場最大風速半徑；B為徑向氣壓分布系數(shù)。

梯度風場模型的計算公式為

式中，Vg(r)為距離臺風中心r處的梯度風速；ρA為空氣密度，取1.2 kg/m3；f=2Ωsin?，為科氏力系數(shù)，其中，?為臺風中心緯度。

式（2）Rmax中和B的表達式大多是基于實測或數(shù)值模擬的擬合結(jié)果[20]，應用較為廣泛的計算表達式見表1和表2。

由表1和表2可以看出，臺風風場模型中各個參數(shù)間相互影響，特別是Rmax和B與臺風中心氣壓差ΔP、中心緯度?、中心移動速度Vf等之間可能都存在一定的關系，并且不同海域獲得的風場參數(shù)統(tǒng)計擬合關系差異較大，需對計算方法進行綜合比較來確定參數(shù)組合計算模式。

表1 徑向氣壓分布系數(shù)B的主要計算表達式Table 1 Main calculation methods of radial pressure profile coefficient B

表2 最大風速半徑Rmax的 主要計算表達式Table 2 Main calculation methods of maximal wind velocity radius Rmax

為得到適用于浙江沿海的臺風風場，本文基于前人總結(jié)的參數(shù)經(jīng)驗公式，分別選取3個典型且適用于我國沿海海域的徑向氣壓分布系數(shù)B、最大風速半徑Rmax的計算方法，對應組合得到9種不同的參數(shù)組合方案，下文稱自匹配模式。B、Rmax的計算表達式及編號見表1和表2，9種自匹配模式的組合情況見表3。

表3 自匹配模式Table 3 Self-matching patterns

3 參數(shù)化風場的驗證

3.1 臺風選取

基于CMA最佳路徑數(shù)據(jù)集[30-31]提供的臺風軌跡、中心氣壓及臺風中心最大風速等數(shù)據(jù)，采用表3所示的自匹配模式，得到不同Rmax和B組合關系的9種參數(shù)風場模型，再現(xiàn)0216號臺風和0414號臺風在發(fā)展過程中的風速變化。為驗證風場的準確性，本文選取大陳、石浦、嵊泗氣象站的實測風速數(shù)據(jù)與數(shù)值模擬結(jié)果對比。圖1給出了兩場臺風的軌跡示意圖以及實測站點的位置。

3.2 風場驗證

根據(jù)現(xiàn)有資料，對海平面10 m高度風速進行驗證，圖2、圖3分別為0216號臺風、0414號臺風的風速對比結(jié)果。

從圖中可以看出，9種自匹配模式模擬結(jié)果與實測風速在變化趨勢上基本相似，但在不同臺風中，各模式在風速極值與臺風影響前期的風速模擬上差別較大。R3計算方法在極值風速的模擬上相對偏大，R1、R2方法模擬結(jié)果相近。B1、B2計算方法在臺風影響前期風速和極值風速的模擬上效果較好，B3計算方法整體風速偏小，不能很好地反映風速變化特征。

為更好地評估各種計算模式的計算偏差，選用相關系數(shù)C0與極大值相對誤差C1進行檢驗，公式如下：

式中，F(xiàn)i、Oi、Fmax、Omax分 別為模擬值、實測值、模擬值最大值和實測值最大值；N為樣本總?cè)萘俊?/p>

采用極大值相對誤差C1來評價模型在預測極值風速上的精度，采用相關系數(shù)C0來評價臺風整體過程中模擬風速的綜合效果。當C0=1.0時，表示兩變量完全線性相關；一般情況下，若C0＞0.8，就可認為兩變量相關性較好。相關系數(shù)和均方根誤差結(jié)果見表4。

由表可知，綜合不同臺風及不同站點的相關系數(shù)和相對誤差，自匹配模式B1R1在9個模式中是和實測結(jié)果符合度最好的，其次是B1R2。由此可知，選用江志輝等[8]提出的R1計算方法以及Jakobsen和Madsen[15]的B1方法進行臺風風場的構(gòu)造，能較準確地反映該區(qū)域臺風風場特征。

圖1 臺風軌跡及氣象站實測站點Fig. 1 Track of typhoon and location of weather observation sites

圖2 0216號臺風風速數(shù)值模擬結(jié)果與實測對比Fig. 2 Comparison of Typhoon 0216 wind speed between simulated and observed results

4 風暴潮數(shù)值模擬

4.1 模型設置

本文基于MIKE 21水動力模型，結(jié)合Holland 參數(shù)風場模型和全球潮汐模型TPXO 7.2，建立了天文潮-風暴潮耦合數(shù)值模型。該模型將二維水動力模塊與海面風強迫力相結(jié)合，在風暴潮模擬研究中具有良好的適用性[32-33]。

模型以浙江海域為研究對象，天文潮-風暴潮耦合模型計算海域為26.1°～31.0°N，120.1°～124.35°E，如圖4a所示?？紤]到浙江沿海島嶼眾多、岸線復雜，計算域采用不規(guī)則三角網(wǎng)格劃分，網(wǎng)格從外海至近岸逐漸加密，外海網(wǎng)格最大分辨率為9.5 km，近岸臺州沿岸網(wǎng)格分辨率為50～200 m，寧波、溫州海域網(wǎng)格分辨率為100～1 000 m。整個計算域共287 503個單元、146 733個節(jié)點，網(wǎng)格圖見圖4b。模擬區(qū)域的岸線數(shù)據(jù)采用美國國家海洋和大氣管理局（NOAA）的高精度岸線數(shù)據(jù)；寧波市石浦港至溫州市鰲江口一帶的近岸海域水深由海圖數(shù)字化得到，其他海域的水深數(shù)據(jù)引用NOAA的30″×30″的ETOPO1水深數(shù)據(jù)，模型最小時間步長為 0.1 s，最大時間步長為60 s。

圖3 0414號臺風風速數(shù)值模擬結(jié)果與實測對比Fig. 3 Comparison of Typhoon 0414 wind speed between simulated and observed results

4.2 模型的驗證與分析

4.2.1 臺風資料

0216號 臺 風 于2002年8月29日6時（UTC）生成，9月7日10時30分（UTC）在浙江蒼南登陸, 登陸時中心氣壓為965 hPa，近中心最大風力達12 級以上。登陸時, 溫臺等地出現(xiàn)了狂風巨浪、暴雨高潮等現(xiàn)象，沿海最大風速超40 m/s，且臺風過程正值天文大潮期，實測資料顯示臺風過程中整個浙江沿岸最高潮位均超警戒水位，浙中、浙南岸段受災尤為嚴重[34]。0414號臺風于2004年8月8日12時（UTC）生成，12日12時（UTC）在浙江溫嶺市石塘鎮(zhèn)登陸，登陸時近中心風速與中心氣壓分別為45 m/s和950 hPa。0414號臺風強度高、風力大、影響范圍廣，是1956-2004年間在浙江沿海登陸的最強臺風。臺風登陸適逢天文大潮起潮期，引起登陸點附近潮位站2.30～3.58 m的最大增水。

表4 風速計算值與實測值的誤差Table 4 Error between simulated and observed results of wind speed

圖4 模型水深圖和網(wǎng)格圖Fig. 4 Water depth map and grid map of modle

選取上文的參數(shù)化風場模型生成的0216號臺風和0414號臺風風場作為天文潮-風暴潮耦合模型的驅(qū)動風場，對兩場臺風影響下的天文潮潮位及風暴潮增水水位進行驗證，臺風路徑見圖5。

4.2.2 天文潮驗證

潮汐和風暴潮的非線性作用在風暴潮增水模擬中尤其重要[35]。檢驗天文潮-風暴潮耦合模型能否較好反映沿岸的水位變化情況，首先要進行潮位驗證。進行潮位驗證時，模式只給予開邊界潮汐水位驅(qū)動，不加入風場。本文根據(jù)國家海洋信息中心提供的2002年和2004年潮汐表中浙江省的7個潮位站的潮位預報數(shù)據(jù)進行天文潮潮位驗證，驗證站點位置見圖5。

兩場臺風模型的模擬時間分別為2002年8月29日6時至9月8日6時（UTC）、2004年8月6日0時至8月14日0時（UTC），軟啟動時間均為72 h，軟啟動期間不進行潮位驗證，天文潮驗證結(jié)果如圖6、圖7所示。

圖5 臺風路徑及驗證站點Fig. 5 Track of typhoon and location of verification sites

圖6 0216號臺風期間天文潮位數(shù)值模擬結(jié)果與潮汐表對比Fig. 6 Comparison of tide level between simulated and tide table value during Typhoon 0216

由圖可知，潮位模擬值與實測值整體趨勢較為一致，在各站點符合較好。

0216號臺風期間，從潮時差來看，各站高低潮位出現(xiàn)時間與潮汐表數(shù)據(jù)的相位差在1 h左右，誤差原因是健跳站、坎門站在8月29日至9月2日的高潮潮時提前約1 h；東門村站、海門港站、鎮(zhèn)海站在9月4-8日的低潮潮時滯后約1 h。從潮位差來看，7個站點的潮位高潮絕對平均誤差最大為0.16 m（鎮(zhèn)海站），最小為0.03 m（黃大岙站）；低潮絕對平均誤差最大為0.24 m（鎮(zhèn)海站），最小為0.10 m（海門站、坎門站）。其中，坎門站模擬高潮位偏低，海門站、海門港站在高潮位模擬中的誤差主要源于9月1-4日3次低高潮的潮位偏低。

0414號臺風期間，各站高低潮位的出現(xiàn)時間與潮汐表數(shù)據(jù)的相位差在1 h左右，誤差主要是東門村站、海門港站、海門站、健跳站在8月6-8日以及12-14日的低潮潮時滯后約1 h導致。6個站點的潮位高潮絕對平均誤差最大為0.23 m（海門港站），最小為0.09 m（東門村站）；黃大岙站模擬高潮位偏高，健跳站、坎門站、海門站、海門港站在8月10-13日的3次低高潮的模擬潮位偏低。低潮絕對平均誤差最大為0.27 m（健跳站），最小為0.05 m（坎門站）；除健跳站在低潮位的模擬中誤差較大外，其余站點低潮絕對平均誤差均在0.1 m內(nèi)。

可以看出，模型模擬結(jié)果與潮汐表數(shù)據(jù)較為一致，能夠較好地反映研究區(qū)域的潮位變化特征。但也存在一定誤差，例如個別站點在模擬后期出現(xiàn)低高潮位、相位滯后的問題，如健跳站、坎門站、鎮(zhèn)海站整體相比其他站點不能很好地再現(xiàn)潮位的變化等。產(chǎn)生的主要原因如下：

（1）模型本身的誤差。建立模型時，為確保計算效率，個別站點分辨率較低，導致模擬潮位誤差，如健跳站、鎮(zhèn)海站；同時，個別站點水深較淺、地形特殊，非線性效應顯著，易導致潮波畸變；另外，模型中未考慮漫灘也是造成結(jié)果差異的原因之一。

（2）初始潮位及驗證資料的誤差。一方面，模型選取TPXO 7.2在外邊界的模擬結(jié)果來確定模型開邊界潮位，與實況存在一定誤差；另一方面，驗證資料選自潮汐表，忽略了實際水位是潮汐項和擾動項之和。

（3）模擬過程中的誤差。模型以開邊界潮汐水位為驅(qū)動，進行由外海向近岸的計算，模擬范圍較大，會增大近岸的累積誤差[36]；此外，不同天文潮位時實際海床糙率不同，難以在模型中準確反映[23]。

4.2.3 風暴潮增水驗證

模擬總水位與模擬天文潮水位相減可以得到風暴潮增水水位。其中，模擬總水位是在模擬中添加潮位和風場的共同驅(qū)動得到。

0216號臺風與0414號臺風增水的模擬結(jié)果與實測數(shù)據(jù)的對比分別如圖8、圖9所示。從圖中可以看出，風暴潮增水的模擬結(jié)果與實測值的過程吻合度較好，主峰契合。

圖7 0414號臺風期間天文潮位數(shù)值模擬結(jié)果與潮汐表對比Fig. 7 Comparison of tide level between simulated value and tide table value during Typhoon 0414

0216號臺風期間，坎門站距離臺風中心較近，其過程風暴潮增水較大，出現(xiàn)的增水是標準型增水過程；海門站、健跳站呈現(xiàn)波動型增水過程；鎮(zhèn)海站呈現(xiàn)混合型增水過程。從風暴潮主振階段來看，各站峰值時間差在1 h內(nèi)，風暴潮增水極大值最大誤差出現(xiàn)在海門站，為0.19 m；最小誤差出現(xiàn)在坎門站，為0.09 m；4個站的平均誤差為0.14 m。除海門站以外，其余站點誤差均在0.15 m以下。從先兆波增水情況來看，4個站在趨勢上與實測增水情況符合，出現(xiàn)先兆波的周期性增水且增水時間長；但先兆波的周期性增水水位較實測值偏低0.05～0.60 m，其中，海門站、坎門站模擬效果稍好，除個別周期的先兆波增水水位的模擬誤差大于0.5 m外，其余的增水水位誤差在0.05～0.30 m范圍內(nèi)，健跳站、鎮(zhèn)海站模擬水位誤差較大，偏低0.10～0.65 m。除模型中岸線和水深數(shù)據(jù)的可靠性因素外，誤差可能和站點與臺風中心的距離、網(wǎng)格精度有關。先兆波是由外海強風低壓所引起的自由長波，由臺風中心向外傳開，并比臺風中心超前傳至岸邊，站點離臺風中心較遠會導致摩擦損耗以及累計誤差，例如離臺風行進路徑相對較遠的健跳站、鎮(zhèn)海站，先兆波增水模擬誤差較其他站點偏大；同時，先兆波發(fā)生時期，站點附近風力較小，風暴潮與天文潮的耦合作用偏弱，也會導致模擬增水過程波動不明顯。另一方面，網(wǎng)格精度的問題也可能會導致難以準確刻畫先兆波的增水情況。

圖8 0216號臺風期間增水水位數(shù)值模擬結(jié)果與實測對比Fig. 8 Comparison of storm tide between simulated and observed results during Typhoon 0216

0414號臺風期間，各站點距離臺風中心較近，過程增水水位較大，呈現(xiàn)標準型增水。從風暴潮主振階段來看，各站峰值時間差均在1 h內(nèi)，風暴潮極值增水最大誤差出現(xiàn)在海門港站，為0.22 m；最小誤差在海門站, 為0.15 m；3個站的平均誤差為0.18 m。除海門站以外，其余站點誤差均在0.20 m以下。綜上，可以認為模型能夠較好地模擬臺風過程中的風暴潮增水變化情況。

圖9 0414號臺風期間增水水位數(shù)值模擬結(jié)果與實測對比圖Fig. 9 Comparison of storm tide between simulated and observed results during Typhoon 0414

5 結(jié)論

本文基于Holland參數(shù)化風場模型，研究了臺風風場模型中關鍵參數(shù)最大風速半徑Rmax和徑向氣壓分布系數(shù)B對風速的影響，將3種典型的Rmax、B參數(shù)計算公式進行對應組合得到9種不同的參數(shù)化方案，對其在兩場臺風下的各站點風速模擬值進行了驗證分析；采用選定的Holland參數(shù)風場模型構(gòu)成的0216號臺風和0414號臺風風場作為輸入風場，利用MIKE 21數(shù)值模型進行浙江海域風暴潮數(shù)值模擬，并對實測站點的天文潮位及風暴潮增水水位進行了驗證分析，得出如下結(jié)果。

（1） 采用參數(shù)計算方法比選的方式確定Holland參數(shù)化風場模型，能使風場模擬更為準確。針對參數(shù)化風場模型中不確定性較強的最大風速半徑Rmax和徑向氣壓分布系數(shù)B，不同的參數(shù)風場公式獲得的風速模擬結(jié)果表現(xiàn)出明顯的差異。選用大西洋海域Graham[29]的最大風速半徑Rmax方法所得的風場模型在極值風速的模擬上偏大，Powell等[10]的徑向氣壓分布系數(shù)B方法在臺風影響前期的模擬風速偏小。采用江志輝等[8]、Fang等[27]提出的Rmax計算方法及Jakobsen和Madsen[15]、Fang等[17]的B計算公式組合得到的參數(shù)風模型能較準確地重現(xiàn)研究區(qū)域臺風過程的風速變化情況。

（2）浙江沿海地區(qū)，采用江志輝等[8]的Rmax和Jakobsen及Madsen[15]的B方法確定的參數(shù)化風場模型可較準確地模擬臺風過程中整體與極值風速的變化特征。

（3）基于本文提出的風場參數(shù)選取方法構(gòu)建臺風場，輸入MIKE 21風暴潮數(shù)值模型，在風暴潮增水的模擬精度上已經(jīng)達到要求。從風暴潮主振階段來看，各站增水最大時間誤差在1 h內(nèi)，增水最大值誤差在0.25 m以內(nèi)，平均誤差在0.15 m左右。從先兆波增水情況看，0216號臺風模擬過程中，各站風暴潮增水趨勢與實測情況相符，但在增水水位上偏低0.05～0.60 m不等，這可能是由于地形、岸線數(shù)據(jù)及網(wǎng)格精度等問題導致。未來可采用更高分辨率或是選用臺風發(fā)生期間的岸線水深數(shù)據(jù)，提高風暴潮數(shù)值模型的準確性。