渤海灣風浪預報數值模式的研究

李大鳴，潘 番，羅 浩，解以揚

（1.天津大學水利工程仿真與安全國家重點實驗室，天津300072；2.天津市氣象科學研究所，天津300074）

渤海灣風浪預報數值模式的研究

李大鳴1，潘 番1，羅 浩1，解以揚2

（1.天津大學水利工程仿真與安全國家重點實驗室，天津300072；2.天津市氣象科學研究所，天津300074）

為提高渤海海域海浪數值預報工作的精度，從風浪模式的角度，借助SWAN模型的物理過程，針對風輸入項中線性系數的改進進行了論證。與榮城站觀測結果對比表明，默認線性系數0.001 5較為合理。結合Komen、Janssen和Westhuysen 3種不同風指數增長表達式和對應的白浪破碎表達式，設計和研究不同組合方案下，風浪模式在渤海灣海數值模擬中的適用性。3個方案模擬的波浪特征值與實測值進行比較，結果顯示，有效波高值與實測值變化趨勢相同，模擬平均周期均偏小，且Westhuysen方案模擬效果最好。在計算結果統計分析中，模擬結果與實測值的相關系數相差不大，模擬效果較理想。

渤海灣；風浪模式；預報數值模式；SWAN模型；WRF模式；風能輸入；白浪耗散

隨著社會的不斷發展與進步，災害性海浪這種不可忽視的自然災害，愈發得到人們的關注，其中與災害性近岸波浪的接觸最為直接與頻繁。災害性近岸波浪對于近海的海岸基礎建設、工業生產作業、漁業養殖捕撈等相關行業有著深遠的影響，與我國海洋社會經濟領域的健康發展密不可分。因此，研究海浪的形成發展傳播，進行海浪數值模擬，對波浪要素的實時預報具有重要意義。文中采用淺海、近岸海浪模式—SWAN模式，用于海岸、湖泊和河口的波浪模擬?；跉W拉近似的動譜平衡方程，考慮了影響海浪生成、發展的較多物理過程，還包含了當前海浪預報的較新成果。

國外對于SWAN模型的應用較早。Guan Changlone等［1］利用NDBC國家數據浮標中心的浮標數據，模擬了發生在卡羅萊納州海岸的兩次颶風過程。擬合結果與浮標數據資料對比表明，Komen形式比Janssen形式模擬效果更好。Ge Yijun等［2］從Komen、Janssen和Westhuysen 3種不同風輸入方案的角度，對強臺風Winnie引起的波浪進行模擬對比試驗。與其他2種方案相比，Westhuysen方案的模擬誤差最小。W.C.Dragani［3］等在對SWAN模型評估中，得出模擬波高低于實測值波高，在低風速情況下尤為顯著的結論。

國內SWAN模型以往主要應用于東海、南海和臺灣海峽3個海域的臺風浪研究中，而現今也逐漸應用于渤海海域［4］。賈曉等應用SWAN模式，模擬了恒定風速和臺風風場下的渤海海域波浪分布，通過模擬波高和海上測站波高比對，發現大風速條件下模擬波高偏大，于是提出修改模式中的參數，并驗證了修改風能輸入項后的模擬結果更接近于實測值［5］。楊德周等人應用SWAN模式，模擬渤海海域在風速小于30 m／s情況下的波浪場，并與實測值比較，結果顯示，模擬的最大有效波高偏?。?］。于是提出改進模式中風輸入的線性增長系數，修改后的模擬結果更接近于實測值。本文借助SWAN模型中的3種風能輸入指數表達式和對應的3種白浪破碎表達式，進行試驗、選擇，以便獲得適用于渤海灣的風浪數值模式。

1 風能輸入與白浪耗散

1.1 風能輸入

SWAN模式內部，同時引入了基于風-波相互作用的Miles切流不穩定機制和Phillips共振機制，描述海面風向波浪的能量轉移，使能量轉移模擬更加合理。Miles切流不穩定機制考慮的是波浪隨時間的指數增長，而Phillips共振機制考慮的是波浪隨時間的線性增長［7］。綜合這2種機制，SWAN模式內部風輸入項為線性增長和指數增長之和：

線性增長A計算公式如下：

式中：α為Pillips線性增長系數，H是估值過濾因子，σ為平均頻率，是全發展海況的波譜峰值頻率，θw是風向，θ為平均波向，U?為摩阻風速。

在SWAN40.85版本的模式中，選取了3種風輸入能量轉移的指數增長形式［8］，即常用的Komen和Janssen形式，以及較新的Van der Westhuysen形式。

1.1.1 Komen形式

該形式下的表達式［9］為U?／cph的函數：

式中：cph為相位速度；ρa，ρw分別為空氣和海水的密度；θw為平均風向。

1.1.2 Janssen形式

Janssen形式［10］基于擬線性波理論，表達式為

式中：β為Miles常數，可由無量綱臨界高度λ確定：

式中：κ是Karman常數，取0.41；ze是海面有效粗糙度；U?可以由給定的海面10 m風速U10和風能量譜密度E（σ，θ）確定：

式中：z0是海面粗糙度，α＝0.01，海面有效粗糙度ze依賴于z0和海面狀況，τw為波引發的壓力，τ為總表面壓力。

1.1.3 Van der Westhuysen形式

Van der Westhuysen形式［11］修改于Yan（1987年）的指數增長形式，提出的第3種指數增長計算方式表達式如下：

式中：D、E、F、H為常數有效系數。為了提高波浪在成熟期的增長進度與Snyder曲線的擬合度，Van der Westhuysen等將有效系數取值修改為D＝4.0×10-2，E＝5.52×10-3，F＝5.2×10-5，H＝-3.02×10-4。

1.2 白浪耗散

在SWAN海浪模式以及其他的第3代海浪模式中，最早都使用Hasselmann脈沖原理，計算白浪破碎引起的能量耗散，該表達式由波陡控制［12-13］。WAMDI Group重新確定了白浪耗散計算表達式，使之可以應用于有限水深情況，具體計算公式如下：

式中：Br為常數，表示飽和門限參數，取值1.75× 10-3；指數p控制耗散項的能量譜頻率縮放，取值范圍為2～4；Cg表示波群速度，k表示波數。

2 模式設置與檢驗

目前SWAN模式主要針對臺風過程這類大風速條件下的風浪模擬，而渤海地區在太平洋海洋氣團和東亞北部大陸氣團影響下，形成季節性風環境，故需檢驗SWAN模式中物理過程參數是否適用于季風條件下波浪發展，因此選擇黃、渤海區域進行模式試驗。

2.1 計算海區

本次試驗計算范圍為黃、渤海區域，即東經117.35°～127.38°，北緯35.07°～41.119°。模型矩形網格分辨率為0.1°×0.1°，網格節點101×61。計算區域如圖1所示。

圖1 計算區域的網格劃分Fig.1 Grid division of the computational region

2.2 地形資料

計算區域的海底地形數據來源于美國地球物理中心（U.S.National Geophysical Data Center）發布的ETOPO地形高程數據，精確到1·′。

圖2 黃、渤海地形等值線Fig.2 Contour for the terrain of Huang-Bo sea

圖3 黃、渤海三維地形Fig.3 Three-dimension graph for the terrain of Huang-Bo sea

經過數據轉換，插值處理后，繪制區域內插值后的地形等值線圖和三維地形圖，分別如圖2、3所示。

2.3 風場資料

本次風浪模式計算中，僅以風場作為海浪的驅動場，采用以完全可壓的非靜力預報模式可進行二次開發的開源軟件—WRF為基礎模型的TJ-WRF中尺度風場數值預報系統?；赩3.4.1版本WRF模型的TJ-WRF風場預報系統，小區域針對天津全境，時間分辨率設定為1 h，垂直分辨率51層，預報時效為3 d。采用粗、細網格兩重嵌套，粗、細網格水平分辨率分別為5 km和1 km。中心點為（N40°，E115°），格點數分別為441×369和151×169，采用的投影方式為蘭伯特投影。

利用WRF模式中的WRF-DA數據同化模塊，即同化常規觀測數據、衛星遙感輻射數據和多普勒雷達數據等多種觀測數據［16-17］，對各種地面、高空觀測資料進行預處理。由于大范圍區域初值很難通過觀測獲得，但是為獲得較好模擬效果需要一個準確初值，這種獲取初始條件的方法稱為數據同化。TJ-WRF模式中包含三維變分同化3D-Var和四維變分同化4D-Var兩種數據同化方法，該變分數據同化的計算方法是通過迭代算法，求解指定的目標函數，以獲取一個對真實狀態的最優估計。

對WRF模式預報的6級以上大風效果，采用國家氣象局規定用的檢驗降水預報的TS評分法和分區評分法進行檢驗，計算公式如下：

準確率：

式中：NA為預報正確次數、NB為空報次數、NC為漏報次數。以實況觀測為基準，如果預報風力級別與實況相同為正確，預報風力級別大于實況為空報，預報風力級別小于實況為漏報。

實況觀測資料選取海上自動氣象站2011年8月1日0時～2012年2月29日23時的觀測資料，時間分辨率為1 h。根據自動氣象站所在位置，將黃、渤海海區分為渤海區域，渤海海峽，黃海北部和黃海中部4個區域，如圖4。對這4個區域的6級以上大風進行分區TS評分檢驗，得出24、48、72 h總體評分檢驗結果如表1所示。

通過對黃、渤海自動氣象站觀測資料6級以上大風的TS評分檢驗，可知對于海上6級以上大風的預報，WRF模式預報偏小的概率大于準確率，且各個分區的預報結果隨著時間的延長，預報準確率下降，漏報率增大，平均誤差增大。3個計算時段預報分區檢驗中，準確率最高的是渤海海峽，準確率最低的是渤海。

表1 總體TS評分檢驗結果Table1 Test result of the the whole TS method

圖4 自動氣象站位置和區域劃分Fig.4 Location of automatic weather station and area division

在黃、渤海區域（N30°～41°，E116°以東）內，有3個具有實測資料的海上浮標站可用，分別是54772（黃海北部，榮城，N37.5°，E122.5°）、54641（渤海西部，曹妃甸，N38.85°，E118.55°）和54558（渤海中部，綏中，N39.35°，E120.58°），測站分布如圖5。

圖5 海上浮標站分布Fig.5 Distribution of offshore buoys

圖6 72 h風場過程Fig.6 Process diagram of wind field in 72 hours

因此，為了比較風場驅動下的海浪嵌套模式模擬結果，選擇數據資料較全，參考價值較高的測站，即渤海海峽附近觀測站榮城站，進行參數化方案結果檢驗。

在利用天津市氣象科學研究所的TJ-WRF模式，預報風場模型結果中，選擇2012年1月18日20時至1月24日20時，2個連續3天風場過程進行SWAN海浪模式初步運行、矯正。每隔24 h繪制海面風場圖檢查風場過程是否合理，計算區域6個時刻風場過程如圖6所示。驗證風場模型，通過將風場圖與實測風場數據比較，風場的大小與方向與實測一致，結果合理可靠。

3 參數選取與試驗設計

在僅以風場作為驅動條件的海浪計算中，風能是海浪成長的唯一能量，風能與海浪間能量傳遞的過程決定海浪成長過程以及最終狀態，為盡可能準確模擬風浪，需要檢驗并調整適用于渤海灣風輸入形式及相應方案。根據風輸入項的組成，將風輸入項拆分為線性增長項和指數增長項分別進行驗證。在確定線性增長取舍的情況下，比較3種指數增長形式哪種更適用于渤海區域海浪模擬。

3.1 線性增長

在SWAN模式預設條件下，風輸入不考慮線性增長項，為比較線性增長項對波高計算的影響，設計以下3套比較方案：線性增長系數分別設α＝0、α＝0.001 5和α＝0.03，風輸入指數增項選擇Komen參數計算。

在模式輸出結果中，選取2012年1月20日18時，一個波高較大的時刻，繪制區域波高分布圖，結果分別如圖7所示。

圖7 不同線性增長系數下的波高分布Fig.7 Distribution of the wave height with different linear growth coefficients

根據天津氣象局劃分的波高等級，如表2。選取α＝0和α＝0.001 5兩種參數方案，計算所得波高值分別為H1和H2，比較上述兩種參數下同一點的波高等級，波級統計如表3所示。

觀察表3，在有無線性增長項情況下，波高在達到2 m時，對應風級屬于勁風，線性增長項的影響明顯，波高差值最大達到0.24 m。在微風狀態下，波高相差較小，僅有0.076 m。在區域范圍內波高分布圖幾乎完全相同。因此增長系數為0.001 5情況下，線性增長項對波高影響較小，即Phillips空氣湍動和波動效應對波浪發展影響較小，在波浪發展中基于風和波相互作用的Miles反饋機制起到決定作用。波高越高，湍流和波動效應即指數增長效應對波浪影響越來越大。

表2 波高分級表Table2 Grade of wave height

表3 線性增長系數α對波高影響Table3 Effect of wave height on linear growth coefficientα

以榮城測站為例，根據測站波高與模擬波高比較，如圖8。在加入風輸入線性增長項后，模擬波高最大值與實測值更接近。綜合考慮整體擬合程度，采用SWAN40.85版本中默認的線性增長系數0.001 5，加入風輸入線性增長項，有利于提高小波高模擬的精度，但是較大波高的情形需要依靠指數增長項調節。

圖8 榮城測站波高比較Fig.8 Comparison of wave heights in Rongcheng station

3.2 指數增長

針對風輸入項和白浪破碎項不同方案設計試驗，研究各個方案對波浪成長的影響，找出最適合渤海灣風浪模擬的計算形式，本次試驗組合方案如表4所示。

圖9 不同方案波高等值線圖Fig.9 Contour of wave height in different programs

表4 風輸入、白浪破碎計算方案表Table4 Computational scheme of wind input and whitecapping dissipation

在上述的方案中，各個物理過程的參數除線性增長系數α取0.001 5外，其他均采用預設值。SWAN模式計算時間為2012年1月18日20時至2012年1月24日20時，2個3天過程。計算結果輸出時間間隔設定為1h，在計算結果中，選取2012年1月20日16時輸出的有效波高值，繪制波高分布圖，進行參數方案比較。Komen參數化、Janssen參數化和Westhuysen參數化3套方案，得到計算區域海浪波高分布，分別如圖9所示。

4 模擬結果檢驗與分析

在第3代波浪模式中：波浪的特征量分別根據作用量譜的某階矩計算。有效波高采用作用量譜的零階矩計算，平均周期采用作用量譜的二階矩和作用量譜的零階矩進行計算。為準確比較3種參數方案適用性，選取研究區域內的海上浮標如圖4，對波浪觀測結果進行點-點檢驗。

4.1 測站點-點結果檢驗

以榮城站為例，3套參數化方案下的模擬與實測的比較圖，如圖10和11所示。

圖11 不用方案模擬與實測平均周期比較Fig.11 Comparison of average period between simulated and measured in different programs

4.2 統計分析檢驗

根據統計分析計算公式，進行輸出結果和觀測資料樣本統計分析，定量分析3種組合方案的適用性，找出最適用于渤海灣海浪的計算模式。

榮成站在3個參數化方案下的有效波高、平均周期計算結果統計誤差由表5所示。

表5 榮城站參數化方案統計結果比較Table5 Comparison of statistical results in parametric schemes in Rong Cheng station

5 結論

1）在3幅計算區域指數增長波高等值線圖中，波高最大值和波高較大值分布范圍有明顯區別，說明3個方案由于采用不同風能輸入和海浪耗散理論，對計算得出波高極值的影響較大。

2）在榮城測站不同方案的有效波高比較圖中，開始時刻，模擬波高總體偏小。隨著計算時間的延長，模擬波高變化趨勢能夠與實測符合較好。在多次的波高比較中發現，一般在計算開始24 h以后，初始條件對此時刻的輸出結果影響微乎其微。

3）在榮城測站不同方案波高比較圖和不同方案周期比較圖中，對榮城站有效波高和平均周期模擬結果，與實測值的相關系數相差不大，分別都達到0.87和0.86，模擬效果較理想。在3個參數方案計算結果統計分析中，Westhuysen方案波高、周期誤差最小，相關系數最大。因此，Westhuysen參數方案更適合于淺海海浪模擬。

研究中還存在一些問題，例如在測站模擬結果中，波峰出現與實測值比較相對延后；模擬平均周期與實測值相比，整體偏小。因此，仍然需要關注影響淺水區域波浪發展的水深誘導破碎物理過程，調整相應參數，以提高模式在淺水區域適用性。

［1］GUAN Changlong，XIE Lian.Investigation of source terms in the SWAN wave model with NDBC buoy data during hurricanes［C］／／Proceedings of the International offshore and Polar Engineering Conference.Toulon，France，2004：217-220.

［2］GE Yijun，ZHONG Zhong，ZHANG Jinshan，et al.Comparison study on wind inputand whitecapping dissipation expressions in numerical simulation of typhoon-generated waves［J］.China Ocean Engineering，2008，22（4）：635-647.

［3］DRAGANIW C，GARAVENTO E，SIMIONATO C G，et al.Campos wave simulation in the outer Río de la Plata estuar-y：evaluation of SWAN model［J］.J Waterway，Port，Coastal，Ocean Eng，2008，134：299-305.

［4］管長龍.我國海浪理論及預報研究的回顧與展望［J］.青島海洋大學學報，2000，10（30）：549-556.

GUAN Changlong.A review of history and prospect for study of sea wave theory and its forecast in China［J］.Journal of Ocean University of Qingdao，2000，10（30）：549-556.

［5］賈曉，潘軍，NICLASEN B.SWAN模型風能輸入項的改進與驗證［J］.河海大學學報：自然科學版，2010，38（5）：585-591.

JIA Xiao，PAN Jun，NICLASEN B.Improvement and validation of wind energy input in SWAN model［J］.Journal of Hohai University：Natural Sciences，2010，38（5）：585-591.

［6］楊德周，伊寶樹，徐艷青，等.SWAN淺水波浪模式在渤海的應用研究——Phillips線性增長比例系數的改進［J］.水科學進展，2005，16（5）：710-714.

YANG Dezhou，YI Baoshu，XU Yanqing，et al.Application of the SWAN wave model to Bohai Sea：improvement of Phillips linear growth term［J］.Advances in Water Science，2005，16（5）：710-714.

［7］PHILLIPS O M.On the generation of waves by turbulent wind［J］.J Fluid Mech，1957，2：417-445.

［8］MULLIGAN R P，BOWEN A J，HAY A E.Whitecapping and wave field evolution in a coastal bay［J］.Journal of Geophysical Research，2008，113（C03008）：1-16.

［9］KOMEN G J，HASSELMANN S，HASSELMANNK.On the existence of a fully developed wind-sea spectrum［J］.Phys Oceanogr，1984，14：1271-1285.

［10］JANSSEN P A E M.Quasi-linear theory of wind-wave generation applied to wave forecasting［J］.Journal of Physical Oceanography，1991，21：1631-1642.

［11］Van der WESTHUYSEN A J.Nonlinear saturation-based whitecapping dissipation in SWAN for deep and shallow water［J］.Coastal Engineering，2007，54（2）：151-170.

［12］李玉成.波浪在淺水區的變形及破碎［J］.海洋學報，1997，19（3）：111-118.

LI Yucheng.Deformation and fragmentation of waves in shallow water［J］.Acta Oceanologica Sinica，1997，19（3）：111-118.

［13］李玉成，董國海.緩坡上不規則波浪的破碎指標［J］.水動力學研究與進展（A輯），l993，8（1）：21-27.

LI Yucheng，DONG Guohai.Breaker indices of irregular waves on gentle beach［J］.Journal of Hydrodynamics（A），l993，8（1）：21-27.

［14］BOUWS E，KOMEN G J.On the balance between growth and dissipation in an extreme depth-limited wind-sea in the southern North Sea［J］.Journal of Physical Oceanography，1983，13：1653-1658.

［15］ALVES J，HENRIQUEG M，BANNER M L.Performance of a saturation-based dissipation-rate source term in modeling the fetch-limited evolution of wind waves［J］.Journal of Physical Oceanography，2003，33：1274-1298.

［16］XIE L，LIU B，LIU H Q，et al.Numerical simulation of tropical cyclone intensity using an air-sea-wave coupled prediction system［J］.Advance in Geosciences，2008，18：19-43.

［17］LIU B，LIU H Q，XIE L.A coupled atmosphere-wave-ocean modeling system：simulation of the intensity of an idealized tropical cyclone［J］.Monthly Weather Review，2010，139：132-152.

A study of the numerical forecast model of wind waves in Bohai Bay

LI Daming1，PAN Fan1，LUO Hao1，XIE Yiyang2
（1.State Key Laboratory of Hydraulic Engineering Simulation and Safety，Tianjin University，Tianjin 300072，China；2.Tianjin Institute of Meteorological Sciences，Tianjin 300074，China）

In order to enhance the accuracy of the numerical prediction of waves in Bohai Bay，the improvement of linear coefficients for the wind input is demonstrated with the aid of the physical processes of the SWAN model，from the perspective of wave patterns.Ithas been shown that0.001 5 is a more reasonable default of the linear coefficient in contrast to the observations of the Rongcheng buoy.Combined with the three wind exponential expressions of Komen，Janssen and Westhuysen and the corresponding whitecapping dissipation expressions，the applicability of the wave mode for the numerical simulation of waves in Bohai Bay was designed and researched with different portfolio schemes.Consequently，the characteristic values and measured values of the three schemes were compared and the results showed that the effective wave height has the same change trend with the measured ones.Both of them have a smaller average cycle during the simulation.Furthermore，Westhuysen has been found to be the best one among the three schemes.With the statistical analysis of the calculation results，it has been found that there is little difference between the correlation coefficients with the simulated and measured values，so the simulation results are satisfactory.

Bohai bay；wind wave；numerical forecast model；SWAN model；WRF model；wind energy input；whitecapping dissipation

10.3969／j.issn.1006-7043.201301045

P722.4；P731.33

A

1006-7043（2014）01-0132-09

http：／／www.cnki.net／kcms／detail／23.1390.U.20131112.1005.014.html

2013-01-23.網絡出版時間：2013-11-12.

國家自然科學基金資助項目（51079095）；國家自然科學基金創新研究群體科學基金資助項目（51021004）.

李大鳴（1957-），男，教授，博士生導師.

李大鳴，E-mail：lidaming＠tju.edu.cn.