北印度洋风浪流数值预报系统:I-设计与实现
2015-03-10蔡夕方张志远楼伟尹朝晖黄司保海军海洋水文气象中心北京006清华大学计算机科学与技术系北京0008国家海洋环境预报中心北京0008958部队气象台海南三亚57202
蔡夕方,张志远,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)提供。计算区域分为二重嵌套,大区(亚非—印度洋区域)、小区(北印度洋区域),大区域为小区域提供边界和强迫场,是个过渡区域,不做预报和检验评估。小区域为本系统关注的北印度洋区域,预报范围为(5ºS—27ºN,39º—105ºE),预报时效为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