北印度洋風浪流數值預報系統：I-設計與實現

蔡夕方，張志遠,2，樓偉，尹朝暉，黃司保（.海軍海洋水文氣象中心,北京006; 2.清華大學計算機科學與技術系,北京0008; .國家海洋環境預報中心,北京0008; . 958部隊氣象臺,海南三亞57202）

?

北印度洋風浪流數值預報系統：I-設計與實現

蔡夕方1，張志遠1,2，樓偉1，尹朝暉3，黃司保4
（1.海軍海洋水文氣象中心,北京100161; 2.清華大學計算機科學與技術系,北京100084; 3.國家海洋環境預報中心,北京100081; 4. 91458部隊氣象臺,海南三亞572021）

摘要：以5°S以北印度洋海域為目標區域，基于區域大氣模式WRF、海浪模式SWAN和海流模式ROMS，建立了北印度洋風浪流數值預報系統。該系統具有運行穩定、針對性強和產品定制靈活等特點，為海洋環境預報、艦艇航行保障、軍事行動準備等提供準確及時的常規和定制風浪流數值預報產品。

關鍵詞：北印度洋；風浪流；數值預報；定制產品

1　引言

印度洋位于亞洲、大洋洲、非洲和南極洲之間，北部封閉，海岸線曲折，主要有紅海、亞丁灣、阿曼灣、阿拉伯海、波斯灣、孟加拉灣等。印度洋地理位置重要，是連接太平洋和大西洋，貫通亞洲、歐洲、非洲、大洋洲的交通和石油輸送紐帶，也是中國海軍亞丁灣護航編隊必經之路，具有重要的經濟和軍事意義。

由于北印度洋受亞洲、歐洲、非洲大陸環繞，受陸地影響極大。隨著季節更替，海陸熱力差異造成較大的氣壓梯度變化，對北印度洋的風、浪影響比較大，同時受季風控制，北印度洋環流為季節洋流，冬季受東北季風影響，方向自東向西，夏季受西南季風影響，方向自西向東。劉金芳等[1]利用1980—1990年海上氣象船舶報資料對北印度洋海域風、浪要素進行了系統分析，指出北印度洋每年11月至翌年3月為東北季風時期，5—9月為西南季風時期，而4、10月為季風的轉換季節。李培等[2]較詳細的分析了北印度洋風等海洋氣象要素的分布特點及其年變化規律。楊永增等[3]基于MASNUM海浪數值預報系統的全球10天后報資料，分析了北印度洋區域波浪分布特征。由于該地區受季風控制顯著，夏季波浪大于冬季；在空間分布上，西部比東部風大、浪大，在亞丁灣、索馬里外海波浪最大。吳方華等[4]基于LICOM模式分析指出，冬季受東北季風影響，斯里蘭卡南段自東向西的北赤道流到非洲東岸匯入向南的索馬里流；夏季受西南季風影響，自西向東的西南季風流與赤道南側馬達加斯加流連成一體，形成一支向北的強大越赤道流，涌升作用將冷水帶到海表，從而使索馬里沿岸海溫達到全年最低。

海洋環境是一個中尺度天氣系統中大氣、海流、海浪、海冰、陸面徑流等各要素相互作用、相互影響的結果，近年來，國內外不同保障單位根據需求逐漸建立了針對不同海域的大氣—海浪—海流（簡稱為風浪流）數值預報系統，并展開了業務化運行及結果檢驗工作。美國海軍艦隊數值氣象與海洋中心運行全球和區域氣象與海洋模式，同化全球氣象與海洋資料，提供0—180 h常規和定制數值預報產品、艦隊和聯合作戰特定數值預報產品。全球風場采用全球環境模式（Navy Global Environmental Model，NAVGEM），區域部分采用三維的中尺度海洋大氣耦合模式（Coupled Ocean/Atmosphere Mesoscale Prediction System，COAMPS),由NAVGEM提供驅動場，嵌套多個區域。

我國從事海浪業務化預報的主要是國家海洋環境預報中心和海軍海洋水文氣象中心。國家海洋環境預報中心2013年建立了全球業務化海洋學預報系統（Chinese Global operational Oceanography Forecasting System，CGOFS），它是我國首個涵蓋全球大洋到中國近海的業務化海洋環境數值預報系統。海軍海洋水文氣象中心建立了針對西北太平洋尤其是臺灣周邊海域的風浪流數值預報系統[5]，并對預報結果進行了檢驗。

從2008年12月起，為了保護索馬里海域中國商船的安全，中國海軍實行的索馬里護航，使艦艇水文氣象保障的需求，由我國近海、西北太平洋等海域，擴展到北印度洋區域。本文面向保障需求，建立了北印度洋區域風浪流數值預報系統，并進行準業務化運行，從數值預報系統結構、預報結果檢驗、產品表現形式等方面進行設計與實現，不僅預報的準確性上達到相應的業務化保障要求，而且在根據保障需求制作不同區域、不同時間段的溫壓風濕氣象預報產品圖、不同海域精細化海浪預報、海流預報產品圖等方面進行了精心設計，達到了較好的表現效果，為今后數值預報業務化系統設計實現和應用服務，具有一定的借鑒意義。

2　模式和運行數據

2.1模式介紹

本系統的構建是基于大氣模式WRF（Weather Research and Forecasting Model）、海流模式ROMS （Regional Ocean Modeling System）和海浪模式SWAN（Simulating WAveNearshore）的，這3個區域模式比較成熟，在國內外得到廣泛應用。

（1）大氣模式WRF

WRF模式采用可壓縮非靜力原始方程組，垂直坐標采用地形追隨靜力氣壓垂直坐標，變量水平配置為Arakawa-C格式，采用3階Runge-Kuta時間分裂積分方案。模式運用了高分辨率的地形和下墊面分類資料，物理過程包括有：云微物理過程，積云參數化，長波輻射，短波輻射，邊界層參數化，表層和陸地面參數化，次網格擴散等物理過程。模式的水平分辨率、垂直方向層次、積分區域及各種物理過程可根據用戶需求調整。

（2）海流模式ROMS

ROMS[6-9]是一個三維非線性、自由表面的斜壓原始方程模式，在垂直方向上采用S坐標，該坐標保持了σ坐標隨地形漸變的特點，同時可以根據實際需求，在垂向上進行局部加密（溫躍層，海表面或者海底）。隨地坐標在垂向離散上很容易造成壓強梯度力計算誤差過大，對此ROMS中設計了專門的算法減小該誤差。ROMS在水平方向上采用C網格，其水平網格既可以是笛卡爾坐標形式，也可以設置成球面坐標。ROMS在平流項離散上，有二階中差、四階中差等差分方法。笛卡爾坐標系下，ROMS的基本控制方程在原始方程的基礎上采用了包辛內斯克近似和靜力近似。

（3）海浪模式SWAN

海浪模式SWAN[10-12]主要包含以下主要特征：海浪傳播、淺化、由定常水深或海流引起的折射、非定常水深或海流引起的頻移；風浪生成；三波和四波相互作用；白冠破碎、海底摩擦、水深變淺引起的海浪破碎；海浪增水；海浪遇到水下障礙物的透射和反射；海浪繞射。SWAN模式支持規則正交網格、曲線網格、三角網格，笛卡爾坐標，球面坐標網格剖分方案，提供了與多種模式的嵌套接口。

該模式采用的數值算法是三階和一階全隱式，優點是無條件穩定，高分辨率計算的時間步長不受空間分辨率的限制（在空間網格數相同的情況下，隨著空間分辨率增加不會增加計算量）；缺點是SWAN模式為了滿足全隱式差分格式這一計算限制，并行計算只能沿單一方向剖分，并行規模受到限制最大，大規模擴展潛力較小。該模式支持串行、MPI并行計算。

2.2數據準備

（1）模式地形選取和預處理

模式地形的好壞對模式的穩定性和可靠性至關重要。對于高分辨率的模式來說水深數據要求水平分辨率高、水深精度高。本系統利用ETOPO1[13]的水深數據來設置模式的地形和水深。該水深數據水平分辨為1′，集合了全球10多種不同的水深數據源。

海流模式分辨率為1/12°，沒有水深數據分辨率高，因此本文利用Blackman徑向濾波器[14]對模式水深進行預處理，將原始水深數據預處理以免影響海流模式模擬。

（2）海流啟動過程驅動場和強迫場

海流模式ROMS在正式預報之前，需要有較長時間的啟動過程，為此本文采用如下驅動場和強迫場對海流模式進行spin-up[15]。

驅動場方面，利用氣候態驅動場對模式進行20年的驅動，再利用2003—2013年的驅動場對模式進行驅動。使模式盡快達到相對穩定的狀態，此處主要選取Levitus[16]多年平均的溫、鹽場結合SODA、GDEM、WOA2009等三維溫、鹽、流數據，為驅動場進行積分。

強迫場方面，為使模式盡快達到與現實接近的狀態，要求海面強迫場的空間和時間分辨率越高越好，模式需要的驅動場主要包括海表面風場、海表面大氣壓、海表面氣溫、長波輻射、短波輻射等。本系統采用ECMWF再分析風場，該強迫場時空分辨率都較高，水平分辨率為0.5°，時間分辨率為6 h。海洋表面由海表風應力、熱通量和淡水通量驅動，其中表面熱通量由大氣模式輸出，海表熱通量的具體計算方法如下：

式中，Qwe為通過海氣邊界層下界面的熱通量。Qae為通過海氣邊界層上界面的熱通量，Qae是蒸發和降水引起的熱通量、短波長波輻射引起的熱通量和感熱通量三者之和，其中ρa為大氣密度，cap為大氣的定壓比熱，CT為Stantan數，θva為大氣溫度，θvs為海表溫度。Qlate為潛熱通量，通過公式得到，L是單位質量的淡水汽化熱。為蒸發量，通過公式得到，其中ha和hs分別為大氣比濕和飽和比濕。

海表風應力的計算利用公式：

3　預報區域與計算方案

3.1預報區域

整個預報系統的流程如圖1所示，核心由風浪流3個預報模式和產品制作應用組成，邊界場主要依托可實時接收到的全球大氣數值預報系統（GFS）、全球海流數值預報系統（POM）和大區域海浪系統（SWAN）提供。計算區域分為二重嵌套，大區（亞非—印度洋區域）、小區（北印度洋區域），大區域為小區域提供邊界和強迫場，是個過渡區域，不做預報和檢驗評估。小區域為本系統關注的北印度洋區域，預報范圍為（5oS—27oN，39o—105oE），預報時效為0—72 h。系統根據能實際接收資料情況，確定的用于風浪流3個模式的同化數據包括常規觀測資料GTS，非常規的HY-2A高度計有效波高資料、融合海表溫度資料MMGSST和海面高度異常資料MSLA等。系統同化模塊設計與實現過程，我們將另文闡述。

因渦旋、海洋鋒等精細化保障需要，為了較好地辨識中尺度渦旋（模式水平格距均小于第一斜壓Rossby變形半徑），本系統海浪和海流數值預報模式的水平分辨率為1/12°，根據大氣模式與海洋模式水平分辨率配比關系，大氣模式水平分辨率可比海洋模式約粗一倍，因此本系統的大氣模式水平分辨率設為20 km，垂向45層。考慮到海洋溫度、鹽度和密度在海洋上層變化較大而且存在明顯的躍層，海洋數值模式在垂向分層時，上層分層較密、層厚較小。深層海洋海水特征相對穩定，模式分層時在深層海洋相對稀疏、層厚較大，因此本系統海流模式的垂向分為50層，其中200 m以淺為間隔10 m一層，200 m以深，間隔逐步增大。

3.2數值計算

北印度洋風浪流數值預報系統以神威3000 A高性能計算機（計算峰值為23萬億次/秒）作為運行平臺，所用CPU為Intel 4核Xeon E5472處理器，主頻為3.0G Hz。每個節點由2個CPU共享32 GB內存組成，二級高速緩存為12 MB，計算節點由Infiniband高速網絡互聯。在系統業務運行中，每個CPU核均只運行一個MPI進程，每個進程只對應一個線程。

本系統的3個模式本身均已提供并行版本，本文均采用模式本身的MPI并行計算和通信方式。首先進行大氣模式WRF的計算，積分步長為120 s，并行規模為384 CPU核；其次，大氣模式后處理完成后，分別向海浪和海流模式提供強迫場數據；最后，海浪模式SWAN和海流模式ROMS同時運行，SWAN并行規模為128 CPU核，ROMS的并行規模為256個CPU核。從同化數據下載到位，啟動大氣模式到三個模式后處理結束，可在2 h內完成一次預報時長為72 h的預報。

圖1　系統流程

4　產品制作應用

本系統的產品制作過程，分為兩個層面，一方面按照數值預報區域制作常規風、浪、流各類要素的等值線、羽矢圖、填色圖等產品；另一方面，根據不斷變化的實時保障需求提供定制產品應用和快速推送服務，定制化的輸出特定區域、特定要素和特定時效的產品。本系統有3個模式輸出，統一利用NCL（NCAR Command Language）進行產品制作應用，NCL較其他繪圖軟件的優越性表現在，其除了圖形顯示功能外，還有完整的數據處理模塊，常用的數據處理方法如插值、經驗正交函數展開法（EOF），濾波（Filters），小波分析（Wavelets）都可以用它進行處理，同時，由于支持C語言和Fortran語言外部調用，使得其程序簡單易懂，且靈活多變。

4.1常規產品

以NetCDF（Network Common Data Form）[17]作為本系統標準的數據讀寫和存儲格式，將模式結果用線性插值方式插值到等經緯均勻分布的網格上。同時以NCL軟件為基礎實現海平面氣壓、10m高風速風向、氣溫、相對濕度、海平面氣壓、風和位勢高度等大氣要素和浪高、浪向、鹽度、溫度、流速流向等水文要素的自動化繪圖，用于可視化表達和顯示。常規產品通過目錄服務和網站發布的形式推送給保障人員。

4.2定制產品

除了常規產品，本系統針對經常變化的保障需求，以定制化產品提供服務。提供確定海域任意單點溫鹽密廓線、任意斷面溫鹽密分布圖和其他保障產品（中尺度渦、海洋鋒、躍層）等。

（1）任意區域要素預報

在產品制作前，根據保障需求確定特定區域的范圍，將新區域的左下角、右上角的經緯度寫在約定的文本文件中，然后利用shell腳本語句sed將相關經緯度替換畫圖的*.ncl文件，如圖2是根據保障任務的需要，定制了大區—中區—小區3個范圍的水深10 m層的流速和流向情況，并在大區中標注了冷渦（C）和暖渦（W）位置。

（2）任意單點剖面、任意連線斷面要素預報

單點和任意斷面的也是利用sed語句方式把單點經緯度或者斷面兩點的經緯度輸入到約定的文本文件中，替換*.ncl中的參數，然后利用NCL進行畫圖的，如圖3所示，是預報區域鹽度鋒的情況及任意單點（圖中紅點表示）鹽度廓線示意圖。

4.3并行制圖

傳統的數值預報數據可視化軟件NCL、Matlab等，都只能串行執行，在高性能計算機登錄節點提交任務，利用登錄節點畫圖，大量占用了登錄節點的計算和存儲資源，本系統利用NCL的并行版本ParNCL[18]提供并行設計方法，將畫圖任務和模式計算任務一樣，提交到計算節點上完成，這樣大大緩解了登錄節點的壓力，同時明顯縮短了整合系統運行的墻鐘時間。ParNCL v1.0完全繼承NCL的腳本和函數，主要是提供NCL語言一個并行運行環境。ParNCL利用并行NetCDF讀寫接口將數據讀入內存，并分發到多核上安裝NCL的執行方式運行畫圖任務。

圖2　定制區域海流預報圖

圖3　鹽度鋒和鹽度剖面預報圖

5　結論

隨著我國海洋戰略向印度洋挺進，包括海上石油命脈、海軍亞丁灣護航、護送敘利亞化學武器離敘等事件都促使海洋環境數值預報向北印度洋側重和發展。本文從建設一套穩定運行的業務預報系統出發，通過模式的不斷改進，豐富和完善模式參數化過程和物理過程，在保障產品制作等方面有新的突破。經過近一年的天氣學和海洋學統計檢驗驗證，本系統對海洋環境預報及時準確，各種統計檢驗指標正常，對海浪和海流尤其是有效波高的預報有較好的能力，具有一定的使用參考價值。該系統的準確性和穩定性都達到保障需求，可為面向北印度洋的水文氣象保障工作提供準確及時有效的數值預報產品支持。

接下來的工作中，將考慮進行多組重要天氣系統過境的個例檢驗來繼續驗證產品的可靠性，同時在模式的精細化預報和風浪流耦合物理過程方面開展深入的工作，使其逐步形成更具海洋環境預報保障特色的中尺度風浪流數值預報業務系統。

參考文獻:

[1]劉金芳,俞慕耕.北印度洋風浪場特點及最佳航線分析[J].熱帶海洋, 1998, 17(1): 17-25.

[2]李培,張弦,俞慕耕.北印度洋氣候特點分析[J].海洋預報, 2003, 20(3): 25-30.

[3]楊永增,孫玉娟,王關鎖,等.基于MASNUM海浪預報系統的北印度洋波浪特征模擬與預報分析[J].海洋科學進展, 2011, 29(1): 1-9.

[4]吳方華,李薇,劉海龍,等.一個海洋環流模式模擬的北印度洋經向環流及其熱輸送[J].大氣科學, 2008, 32(1): 14-26.

[5]張權,王燕,王關鎖,等.高分辨率浪-潮-流耦合數值預報系統的業務化運行與檢驗分析[J].海洋科學進展, 2011, 29(2): 129-136.

[6] ShchepetkinAF, McWilliamsJ C. The regional oceanic modeling system(ROMS): a split-explicit, free-surface, topographyfollowing-coordinate oceanic model [J]. Ocean Modelling, 2005, 9 (4): 347-404.

[7] WangP, Song Y T, ChaoY, et al. Parallel computation of the regional ocean modeling system[J]. International Journal of High Performance Computing Applications, 2005, 19(4): 375-385.

[8] Moore AM, ArangoH G, Broquet G, et al. The Regional Ocean Modeling System(ROMS) 4-dimensional variational data assimilation systems: Part I-System overview and formulation[J]. Progress in Oceanography, 2011, 91(1): 34-49.

[9] Moore AM, ArangoH G, Broquet G, et al. The Regional Ocean Modeling System(ROMS) 4-dimensional variational data assimilation systems: Part II-Performance and application to the California Current System[J]. Progress In Oceanography, 2011, 91 (1): 50-73.

[10]王毅,余宙文.衛星高度計波高數據同化對西北太平洋海浪數值預報的影響評估[J].海洋學報(中文版), 2009, 31(6): 1-6.

[11] BooijN, Ris R C,HolthuijsenL A. A third-generation wave model for coastal regions: 1. Model description and validation[J]. Journal of Geophysical Research, 1999, 104(C4): 7649-7666.

[12] RisR C, Holthuijsen L A,BooijN. A third-generation wave model for coastal regions: 2. Verification[J]. Journal of Geophysical Research, 1999, 104(4): 7667-7681.

[13] Amante C, Eakins B W. ETOPO1 1 arc-minute global relief model: procedures, data sources and analysis[M]. Boulder, Colorado:US Department of Commerce, NOAA, National Environmental Satellite, Data, and Information Service, NGDC, Marine Geology and Geophysics Division, 2009.

[14] Arbic B K, Garner S T,Hallberg R, et al. The accuracy of surface elevations in forward global barotropic and baroclinic tide models [J]. Deep Sea Research Part II: Topical Studies in Oceanography, 2004, 51(25-26): 3069-3101.

[15] Bleck R, Boudra D. Wind-driven spin-up in eddy-resolving ocean models formulated in isopycnic and isobaric coordinates[J]. Journal of Geophysical Research: Oceans (1978-2012), 1986, 91 (C6): 7611-7621.

[16] Da Silva A M, Young CC, Levitus S. Atlas of Surface Marine Data 1994 Volumes 1-5 [M]//6-10 N N A. U. S. Department of Commerce, NOAA, NESDIS. 1994.

[17]張志遠,宋順強,劉利,等.浪流耦合模式數值模擬及檢驗分析[J].海洋技術, 2011, 30(4): 87-92.

[18] Rew R, DavisG..NetCDF: an interface for scientific data access [J]. IEEE Computer Graphics and Applications, 1990, 10(4): 76-82.

North Indian Ocean wind-wave-circulation numerical forecast system: I-design and implementation

CAI Xi-fang1, ZHANG Zhi-yuan1,2, LOU Wei1, YIN Zhao-hui3, HUANG Si-bao4
(1.Hydro-Meteorological Center of Navy, Beijing 100161 China; 2. Department of Computer Science and Technology, Tsinghua University, Beijing 100084 China; 3.National Marine Environmental Forecasting Center , Beijing 100081 China; 4. 91458 Army Troop Meteorological Center , Sanya 572021 China)

Abstract：A North Indian Ocean wind-wave-circulation numerical forecast system is designed based on Weather Research and Forecasting Model (WRF), Simulating Waves Nearshore Model (SWAN) and Regional Ocean Modeling System (ROMS). The target computational domain is in the north of 5°S of the Indian Ocean. This system has the advantages of stable operation, high pertinence and production customization agility. It can provide timely and accurate products, both conventional and customized, which are of atmosphere, wave and circulation. The production is for marine environment forecast, navigation and military action preparation.

Key words：North Indian Ocean; wind-wave-circulation; numerical forecast; customized product

通訊作者簡介：張志遠（1978-），男，工程師，博士研究生，主要從事海洋環境信息化和數值預報研究。E-mail：generalzzy@139.com

作者簡介：蔡夕方（1966-），男，高級工程師，本科，主要從事水文氣象預報保障研究。E-mail：bjcaixifang@163.com

基金項目：國家自然科學基金面上項目（41275098）；國家海洋局海洋公益性行業科研專項（201005033）

收稿日期：2014-07-13

DOI:10.11737/j.issn.1003-0239.2015.02.002

中圖分類號：P731

文獻標識碼：A

文章編號：1003-0239（2015）02-0007-07