气象学报  2019, Vol. 77 Issue (3): 552-562   PDF    
http://dx.doi.org/10.11676/qxxb2019.020
中国气象学会主办。
0

文章信息

李晓莉, 刘永柱. 2019.
LI Xiaoli, LIU Yongzhu. 2019.
GRAPES全球奇异向量方法改进及试验分析
The improvement of GRAPES global extratropical singular vectors and experimental study
气象学报, 77(3): 552-562.
Acta Meteorologica Sinica, 77(3): 552-562.
http://dx.doi.org/10.11676/qxxb2019.020

文章历史

2018-03-19 收稿
2018-09-30改回
GRAPES全球奇异向量方法改进及试验分析
李晓莉1,2 , 刘永柱1,2     
1. 中国气象局数值预报中心, 北京, 100081;
2. 国家气象中心, 北京, 100081
摘要: 基于总能量模的奇异向量扰动常用于构造集合预报的初始条件。以建立GRAPES(Global and Regional Assimilation PrEdiction System)全球集合预报系统为目的,基于前期研发的GRAPES全球模式奇异向量方法,在GRAPES全球切线性模式和伴随模式2.0版的框架下,开展了引入线性化边界层方案来改善奇异向量结构,并提高奇异向量计算效率的研究。通过连续试验,从奇异向量的扰动能量结构、扰动能量谱及扰动空间分布等方面,综合分析改进GRAPES全球奇异向量的结构及演变特征。试验结果表明,改进后的GRAPES奇异向量方法有效抑制了之前扰动能量在近地面层不合理的快速增长,同时,奇异向量最优扰动的结构更客观地体现了中高纬度区域大气初始条件中的斜压不稳定扰动及其演变,如在初始时刻奇异向量扰动能量主要位于对流层中层,并呈现出随高度向西倾斜的大气斜压特征;经过线性化演变,扰动能量向较大水平尺度转移,并在垂直结构上表现出向对流层高层上传及向对流层低层下传的特征等。针对GRAPES奇异向量迭代求解中伴随模式计算耗时为主的情况,改进伴随模式中广义共轭余差方案的调用方式,并采用大内存存储法来提高其计算效率,进而将奇异向量总计算时间缩短了25%。总之,改进后的GRAPES奇异向量方法,可应用于构建面向业务应用的GRAPES全球集合预报系统。
关键词: 奇异向量    线性化边界层方案    GRAPES模式    切线性模式    
The improvement of GRAPES global extratropical singular vectors and experimental study
LI Xiaoli1,2 , LIU Yongzhu1,2     
1. Numerical Weather Prediction Center of CMA, Beijing 100081, China;
2. National Meteorological Center, Beijing 100081, China
Abstract: The singular vectors (SVs) based on total-energy norm are generally used to represent initial uncertainty in ensemble forecasting. The GRAEPE SV method based on the total-energy norm using the GRAPES dynamical core model has been developed. Impacts and importance of linearized physical processes on the SVs have been widely studied in the literature. To improve the GRAPES SVs, based on the recently developed GRAPES tangent linear model (TLM) and adjoint model (ADM) version 2.0, the implementation of linearized planetary boundary layer (PBL) parameterization on the calculation of extratropical SVs is conducted. The characteristics of SVs and their linear evolutions are measured by energy partition, energy spectra and spatial distribution through continuous experiments over 1-month period. The unreasonable quick energy growth near the surface that was observed previously in the structure of SVs without linearized physics has been greatly improved. Furthermore, the structures of the upgraded SVs are more consistent with those of previous studies, which exhibit the characteristics of typical total-energy based SVs, i.e., the energy maximum is located in the middle troposphere with obvious westward tilt with height in the spatial structure of SVs at initial time; during their growth, there are upward energy transfer to the upper troposphere and downward energy transfer to the surface, and an upscale energy transfer is shown in the energy spectra. The results show that the upgraded SVs are capable of capturing the baroclinic instability in the troposphere. In order to improve the computation efficiency of GRAPES SVs, the ADM that is most computationally expensive is optimized by reducing the use of GCR (Generalized Conjugate-Residual) in the ADM and increasing computational memory. As a result, the total computation time of GRAPES SVs can be reduced by up to 25%. Overall, the upgraded GRAPES SVs are able to meet the expectation for constructing the initial perturbation for GRAPES global ensemble prediction.
Key words: Singular vectors    Linearized planetary boundary layer parameterization    GRAPES model    Tangent linear model    
1 引言

在集合数值天气预报领域,奇异向量方法主要是用于产生集合预报的扰动初始条件。通常来说,如果不考虑数值模式本身的误差,模式预报结果的不确定性主要来自于其初始条件的不确定性。基于大气切线性模式及伴随模式的奇异向量最优扰动是描述数值模式在线性相空间中动力最不稳定、增长最快的一组正交扰动,它体现了初始条件中不确定性的发展方向(Mureau, et al, 1993Molteni, et al, 1996),因此,奇异向量方法作为一种有效描述大气数值模式初始条件中误差概率密度函数的方法,可应用于构建集合预报的初始扰动(Buizza, et al, 1995)。奇异向量方法与其他初始扰动方法如增长模繁殖法(Toth, et al, 1997)、集合转换卡尔曼滤波(Bishop, et al, 2001)、集合卡尔曼滤波(Houtekamer, et al, 2005)以及非线性方法如条件非线性最优扰动法(Mu, et al, 2003)等成为目前集合预报业务应用及研究的主要初始扰动技术。自1992年起欧洲中期天气预报中心(ECMWF)率先建立了基于奇异向量初始扰动技术的全球中期集合预报业务系统(Palmer, et al, 1993),之后其他几个数值预报中心(如日本气象厅、法国气象局等)也相继发展了基于奇异向量技术的全球集合预报系统。

奇异向量方法是一种基于切线性和伴随技术的最优化问题,该方法计算出的奇异向量体现的是在一定最优化时间间隔内,相对于特定的权重算子,在切线性相空间增长最快的一组正交扰动。奇异向量扰动的结构与权重算子的选取、最优化时间间隔的长短、切线性模式及伴随模式中线性物理过程方案的使用情况等因素密切相关。Buizza(1994)Buizza等(1995)研究表明,基于包括干线性物理过程的切线性模式和伴随模式、采用总能量权重模(Total energy norm)和较长最优化时间间隔(如48 h)方法计算出的中高纬度地区奇异向量适用于构建全球集合预报系统的初始扰动,Molteni等(1996)Palmer等(1998)的研究也进一步证实了上述论点。

以发展基于奇异向量初始扰动技术的GRAPES全球模式集合预报系统为目的,刘永柱等(2013)利用GRAPES全球模式干动力框架下的切线性模式和伴随模式(任迪生等, 2011, 简称GRAPES全球切线性模式和伴随模式1.0版),开展了基于总能量模的奇异向量方法研究。研究结果显示,GRAPES奇异向量计算方案合理,求解正确,奇异向量结构和增长特征大体上能体现中高纬度地区对流层中斜压不稳定特征。然而,基于上述工作的GRAPES全球奇异向量在应用于构建集合预报初始扰动时还存在一些不足,主要表现在两个方面:(1)GRAPES奇异向量在线性演化过程中存在扰动能量在近地面层迅速增长的不合理特征,类似的特征在Buizza(1994)研究中也有论述,这种现象是切线性模式中线性化物理过程不完善或缺乏所致。GRAPES奇异向量扰动能量在近地面层的不合理快速增长,会影响扰动在线性演变过程中垂直方向的传播特征,也就是不能体现基于总能量模奇异向量在线性演变过程中扰动能量向对流层高层上传的特征(Hartmann, et al, 1995Hoskins, et al, 2000),与预期的GRAPES奇异向量所应体现的大气对流层中斜压不稳定扰动及其变化特征存在一定差距。(2)基于GRAPES全球切线性模式和伴随模式1.0版的奇异向量计算效率较低,在水平分辨率为2.5°、最优化时间间隔为48 h设置下,其计算时间在5 h左右,如应用于业务集合预报系统,计算效率将是重要瓶颈。

国际上相关研究(Buizza, 1994Buizza, et al, 1995Coutinho, et al, 2004)表明,在切线性模式及伴随模式中引入线性化物理过程(主要贡献来自线性化边界层垂直扩散方案),能有效抑制奇异向量扰动能量在近地面层不合理地快速增长。刘永柱等(2017, 2019)基于业务化的GRAPES全球模式2.0版(沈学顺等,2017),对GRAPES全球切线性模式和伴随模式的设计框架进行了重构,研发了GRAPES全球切线性模式和伴随模式2.0版,大幅度提高了计算效率和精度。基于GRAPES全球切线性模式和伴随模式2.0版框架,本研究进行了GRAPES奇异向量计算技术改进和计算效率优化,在GRAPES奇异向量计算中引入线性化边界层方案,通过连续1个月试验,从奇异向量扰动能量结构、扰动能量功率谱、奇异向量空间结构及其演变等方面,全面分析改进后GRAPES全球奇异向量结构及其合理性,为进一步应用于全球集合预报系统提供客观依据。同时,为了进一步提高GRAPES奇异向量计算效率,着重解决Lanczos迭代求解奇异向量过程中主要耗时部分——伴随模式的计算效率问题,进而缩短计算总耗时,满足业务实时应用的时效要求。

2 GRAPES全球奇异向量方法改进 2.1 基本原理

奇异向量扰动是在大气切线性模式中相对于特定的权重模,在最优化时间间隔内增长最快的一组正交扰动。对于一个初始扰动向量X0,经过GRAPES全球切线性模式(表示为L)一定时间的向前积分可以获得线性演化的扰动向量Xt=LX0,则GRAPES全球奇异向量求解可归结为演化扰动向量Xt模与初始扰动向量X0模比值的最大化问题

(1)

式中,λ为奇异向量值,E为衡量扰动大小的权重模,上式可转换为标准的奇异值分解问题

(2)

式中,Xi表示求解出的第i个奇异向量,λi为对应的第i个奇异值。为获得某个特定目标区域内增长最快的奇异向量扰动,Buizza(1994)应用了一个投影算子P,通过该投影算子将实际计算过程中目标区之外的奇异向量扰动设置为0。文中也采用以上局地投影算子方法,分别获得北半球中高纬度及南半球中高纬度这两个目标区的GRAPES奇异向量,则式(2)可转换为

(3)

通过变量转换,将物理空间扰动向量Xi(t0)转换为适合于数学算法的欧拉空间上无量纲向量来求解,则式(3)可变化为

(4)

式中,LT为GRAPES全球切线性模式的转置,也就是伴随模式。对于式(4)定义的GRAPES奇异向量计算问题,采用Lanczos迭代算法来进行求解,在迭代过程中需多次积分切线性模式L及伴随模式LT。同时,由式(4)也可以看出,决定奇异向量结构的关键因素主要有3个方面:(1)衡量奇异向量扰动大小的权重模E的选择;(2)切线性和伴随模式的特点,主要为线性化物理过程参数化方案的使用;(3)切线性模式L向前及伴随模式LT向后的积分时间长度(最优化时间间隔)。

GRAPES奇异向量计算采用的权重模为总能量模(刘永柱等,2013),因需分析奇异向量扰动能量结构及演变,在此对总能量模的定义做个简介。GRAPES全球切线性模式和伴随模式的预报量为GRAPES模式预报变量中的水平风分量(u, v)、扰动位温(θ′)、扰动无量纲气压(Π′)对应的扰动量可表示为u′、v′、(θ′)′、(Π′)′。基于上述扰动量的总能量模E的计算为

(5)

式中,前面两项之和为扰动动能(KE)模,后两项之和表示扰动位能(PE)模,第3项和第4项分别表示位能中扰动位温(θ′)′和扰动无量纲气压(Π′)′分量的贡献。式(5)中,,为地形追随坐标,λφ分别代表模式球面坐标系中的经度和纬度,cp为干空气的定压比热。TrθrΠrρr分别表示模式参考温度、参考位温、参考无量纲气压及参考密度,计算式分别为

(6)

式中,Rd为干空气气体常数,pr为标准大气,g为重力加速度。

2.2 线性化边界层物理过程方案

基于GRAPES全球切线性模式和伴随模式2.0版,重构了GRAPES全球奇异向量计算模块,引入线性化边界层方案,拟解决前期工作中GRAPES奇异向量扰动能量在近地面层不合理迅速增长的问题。张林等(2008)对GRAPES全球模式所采用的MRF边界层方案的局地垂直扩散方案进行了线性化研究,初步构建了GRAPES全球切线性模式1.0版的线性化边界层方案。在GRAPES全球切线性模式2.0版中引入了GRAPES线性边界层方案,并进行了改进,包括:线性化边界层方案计算过程中轨迹基态采用非线性模式的全物理过程计算结果(刘永柱等,2017);增加了对MRF方案中混合层非局地垂直扩散过程的线性化处理。这里,对线性化边界层方案中自由大气和混合层的线性化处理进行简要介绍。

对于GRAPES全球模式的连续预报变量ψ(水平风分量(u, v)、扰动位温(θ′)和比湿(q)),由线性化边界层方案产生的对应扰动量ψ′的倾向为

(7)

式中,Kψ为湍流扩散系数:对于u′、v′扰动量,对应动量扩散系数KM,对于(θ′)′、q′扰动量,对应热量扩散系数KH。采用Mahfouf(1999)Laroche等(2002)的研究结果,在线性化边界层方案中不考虑湍流扩散系数扰动量(Kψ)的变化以避免噪声。文中所使用的线性化边界层方案,垂直扩散系数(Kψ)是按照MRF方案,分为混合层与自由大气分别处理。

在边界层混合层中,动量扩散系数(KM)和热量扩散系数(KH)分别为

(8)
(9)
(10)

式(8)—(10)中,k为卡曼常数,z为地面之上的垂直高度,h为边界层高度,wS=u*ϕM-1为混合层的特征速度,其中u*为地面摩擦速度,Pr为普朗特数,ϕMϕH分别表示混合层中动量和热量稳定度函数。

在边界层自由大气中,动量扩散系数及热量扩散系数为

(11)
(12)

式(11)—(12)中,l为混合长度,S为风垂直切变(S=|∂U/∂z|),fM(Ri)和fH(Ri)分别为自由大气部分基于局地理查孙数(Ri)的动量和热量稳定度函数。混合层和自由大气中动量和热量稳定度函数计算公式详见Hong等(1996)

2.3 线性化边界层方案的应用试验

基于2013年5月共31 d开展包括线性化边界层方案的奇异向量试验(简称边界层试验):切线性和伴随模式水平分辨率为2.5°,与刘永柱等(2013)设置相同,但将模式垂直层数提高至60层,与目前GRAPES_GFS全球业务系统的垂直层数相同;最优化时间间隔为48 h;Lanczos算法的最大迭代次数为50次,奇异值计算精度为0.01;分北半球中高纬度目标区域(30°—80°N)和南半球中高纬度目标区域(80°—30°S)分别进行计算。为分析引入线性化边界层方案对奇异向量结构的影响,基于上述试验设置,也开展了不包括线性化边界层方案的对照试验。

研究表明(Hartmann, et al, 1995Hoskins, et al, 2000),基于总能量模、采用干线性化物理过程(主要是线性化边界层方案)的奇异向量扰动能量垂直结构的主要特征为:在初始时刻,扰动能量主要分布在对流层中层;扰动在切线性模式的线性演化过程中,扰动能量具有向上传播至对流层高层,向下传播至近地面层的特点。具体表现为,在最后演化时刻,分别在对流层高层和低层出现奇异向量扰动能量大值区。

以2013年5月5日00时(世界时,下同)的试验结果为例,对比对照试验和边界层试验中北半球第1奇异向量和第4个奇异向量扰动总能量模在初始时刻和最后演化时刻的垂直结构(图 1)。在初始时刻,奇异向量扰动总能量模经标准化处理,其值为1,为方便比较,图 1中初始时刻的能量模放大了10倍。可以看出,在初始时刻,两个试验中的第1奇异向量和第4奇异向量的总能量分布特征比较相似,扰动总能量主要分布于模式第13—25层(对应850—500 hPa),其中能量极大值出现层次略有不同,对照试验中第1奇异向量是在模式第18层附近出现峰值,而边界层试验中是在模式20层左右出现。然而,在最后演化时刻,对照试验和边界层试验中奇异向量扰动总能量分布表现出明显的差异。对照试验中两个奇异向量的扰动能量演化呈现出一个共同特点,即扰动能量在对流层低层迅速增长,能量最大值出现在模式低层第10层(约920 hPa)附近,特别是对照试验的第1奇异向量扰动能量无论是强度还是结构均集中在对流层低层,以上结构与扰动能量在线性演化过程中向下传播特征有关,但很大程度上是由于缺乏线性化垂直扩散物理过程所致。对比边界层试验中第1奇异向量的扰动能量的演变可以发现,虽然扰动能量下传也占主导,但也清楚地表现出了扰动能量向对流层中高层传播的上传结构。此外,对照试验中第4奇异向量的扰动能量虽也出现上传和下传特征,但与边界层试验相比,扰动能量分布特征存在明显差异,边界层试验中第4奇异向量扰动能量上传的特征更为明显,表现为扰动能量主要分布在对流层中上层,最大值出现在模式层第30层左右。

图 1  对照试验(a、b)和边界层试验(c、d)中北半球第1(a、c)和第4(b、d)奇异向量扰动总能量模在初始时刻(虚线)和最后演化时刻(实线)的垂直结构 Fig. 1  Vertical distributions of total energy norms of the 1st (a, c) and 4th (b, d) SV over the Northern Hemisphere (NH) in the control (a, b) and PBL (c, d) experiments at initial time (dashed lines) and final evolved time (solid lines)

综上所述,引入线性化边界层方案能有效改进GRAPES奇异向量扰动能量在近地面层快速增长问题,使GRAPES奇异向量扰动能量的结构与典型的基于总能量模奇异向量扰动能量垂直分布特征更为一致。边界层试验中奇异向量结构特征的综合分析将在本文第3部分详述。

2.4 奇异向量计算效率改进

GRAPES全球奇异向量计算迭代求解过程中需多次积分切线性模式和伴随模式,尤其当采用48 h的最优化时间间隔时,计算量大且耗时长,在前期研究工作中(刘永柱等, 2013),基于GRAPES全球切线性模式和伴随模式1.0版的奇异向量计算时间长达5 h。显而易见,为实现奇异向量在业务GRAPES集合预报系统中的应用,亟待提高其计算效率。

对于2.3节中奇异向量计算的试验设置,在中国气象局业务高性能计算机系统IBM Flex P460上,采用16节点(每个节点16 CPU)的计算资源,GRAPES奇异向量的计算时间约为75 min,相比前期工作的计算耗时,已大为缩短,但面向业务应用时,计算效率仍需进一步提高。对奇异向量计算时间消耗进行分析,发现约98%的计算总时间(约73 min)来自于Lanczos算法的迭代过程中切线性模式和伴随模式积分的消耗。在总迭代过程中,每步迭代中切线性模式和伴随模式积分的总耗时基本一致,约86 s,其中伴随模式积分耗时较长,约59 s,占每次迭代耗时的70%。因此,缩短伴随模式的计算时间为改进奇异向量计算效率的重点。

对GRAPES全球伴随模式计算中耗时最长的亥姆霍兹模块的广义共轭余差方案(Generalized Conjugate-Residual,GCR)的调用方式进行优化。旧方案中需要额外计算一次广义共轭余差获取基态给广义共轭余差的伴随使用,现在通过在GRAPES全球切线性模式计算广义共轭余差时,采用存储法将基态保存到内存,因此在伴随模式中只需调用一次广义共轭余差伴随,减少了一次广义共轭余差的基态计算,进而提高了伴随模式的计算效率。以采用相同的计算资源(16节点,每个节点16 CPU)为例,优化伴随模式中广义共轭余差调用方式后的总迭代过程中切线性模式和伴随模式积分耗时约降至67 min,其中,每步迭代中伴随模式计算耗时约降至53 s。进而,根据大内存计算资源下,伴随模式并行计算加速高的特点(刘永柱等, 2017, 2019),在优化广义共轭余差方法的基础上,采用增加计算节点的大内存方法,如采用384个CPU计算资源(24节点, 每个节点16 CPU), 切线性模式和伴随模式积分总计算耗时进一步降至约54 min,其中,每步迭代中伴随模式计算时间缩至约38 s。

表 1给出了优化前、后50次Lanczos迭代求解GRAPES全球奇异向量时切线性模式和伴随模式积分耗时分析。可见对应每步迭代,优化后伴随模式计算效率改进显著,采用优化广义共轭余差和增加计算节点的优化方案可使每步迭代耗时缩至约65 s,则50次的Lanczos迭代中切线性模式和伴随模式运行时间约为3250 s,相比优化前,计算时间缩短了25%。这样,改进计算效率后,GRAPES奇异向量总计算时间约为55 min,基本满足了业务实时运行需求。

表 1  GRAPES全球奇异向量并行计算效率改进(单位:s) Table 1  The improvement of parallel efficiency of GRAPES SVs (unit: s)
迭代耗时 优化前 优化1 优化2
16个节点 16个节点优化广义共轭余差 24个节点优化广义共轭余差
切线性模式(次) 28 27 27
伴随模式(次) 59 53 38
切线性+伴随模式(次) 87 80 65
50次总耗时 4350 4000 3250
3 奇异向量扰动能量结构及演变分析 3.1 扰动能量垂直结构及演变

为进一步分析改进后GRAPES奇异向量扰动能量模的垂直结构及演变,图 2给出了2013年5月5日00时边界层试验中北、南半球前20个奇异向量平均的总能量模及其分量动能模的垂直廓线。与图 1类似处理,这里初始时刻的能量模也放大了10倍。可见在初始时刻,两个半球平均后奇异向量的总能量模(黑实线)主要分布在对流层中低层(对应850—500 hPa),呈现出单极值的结构特征,北半球能量最大值约位于模式第20层(约700 hPa),南半球出现峰值的层次略高。分析能量模分量进一步发现,初始时刻,动能模(图 2中黑虚线)及位能模(总能量模减去动能模)分布特征相似,但位能模占主导,这与其他基于总能量模的奇异向量在初始时刻位能模占主导的特征一致(Hoskins, et al, 2000Coutinho, et al, 2004Zadra, et al, 2004Diaconescu, et al, 2012)。

图 2  北半球(a)和南半球(b)前20个特征向量平均的总能量模(TE)及动能模(KE)在初始时刻和最后演化时刻的垂直结构 Fig. 2  Vertical distributions of total energy norm (TE) and kinetic energy norm (KE) averaged for the first 20 SVs at initial time and final evolved time over NH (a) and SH (b)

对比分析初始时刻和最后演化时刻的奇异向量扰动能量模的垂直结构,可以发现GRAPES奇异向量扰动能量的线性演变呈现了上传和下传特征,如北半球奇异向量总能量模(图 2a中红实线)分别在模式第30层(约350 hPa)及模式第10层(约920 hPa)出现极大值的双峰结构。相对而言,南半球奇异向量扰动能量的上传特征更为显著,表现为模式第30层附近出现单峰最大值。此外,在线性演化过程中动能(图 2中红虚线)已成为能量传播特别是上传的主导分量,这与Hoskins等(2000)Coutinho等(2004)论述的ECMWF奇异向量扰动总能量及其动能分量在线性演化过程中的变化特征一致,进一步证实前文所分析的结论,引入线性化边界层方案后的GRAPES奇异向量,表现出了典型的、基于总能量模奇异向量扰动能量的垂直结构及演变特征,合理地体现了大气对流层中斜压不稳定及其能量传播特征。

3.2 奇异向量扰动能量模分量结构及演变

除奇异向量扰动总能量模的垂直结构外,本节从整层总能量模(所有模式层上总能量模之和)出发,分析不同GRAPES奇异向量扰动总能量中动能模和位能模的构成及变化。图 34分别给出了月平均的北半球和南半球前10个奇异向量的总能量模中动能模(KE)分量和位能模(PE)分量的比例及增长特征。可见在初始时刻,两个半球(图 3a4a)前10个奇异向量的总能量模构成中均是位能模大于动能模,约是1.5:1。在最后演化时刻(图 3b4b),虽然不同奇异向量的扰动能量增长速度不同(前面奇异向量增长更快),但各个奇异向量中动能模均成为总能量模中的主导成分。其中,北半球动能模与位能模的比约为2:1,南半球动能模增长更为显著,与位能模的比约为3:1。

图 3  月平均的北半球前10个特征向量在初始时刻(a)及最后演化时刻(b)的总能量模分量动能模(KE)和位能模(PE)构成 Fig. 3  Energy partition of kinetic energy (KE) and potential energy (PE) terms of total energy norm at initial time (a) and at final evolved time (b) for monthly average of the first 10 SVs over NH
图 4  同图 3,但为南半球 Fig. 4  Same as Fig. 3 but for SH

以上结果表明,GRAPES奇异向量扰动总能量模的结构及演变特征与其他研究工作中类似的分析结果相符(Zadra, et al, 2004Diaconescu, et al, 2012)。初始时刻,总能量模中位能模占主导,体现了大气的斜压不稳定特征,在线性演化中扰动位能向动能转化,使得在最后演化时刻,动能模分量在总能量模中占主导。

3.3 GRAPES奇异向量扰动能量谱特征

通常使用能量谱来分析奇异向量扰动的空间尺度特征。不少研究结果表明,对于基于总能量模的奇异向量,其扰动能量在线性演变中存在着升尺度能量转移的能量谱特征(Buizza,1994Coutinho, et al, 2004Lawrence, et al, 2009)。对边界层试验中前10个奇异向量平均的扰动总能量进行了能量谱展开,分析在初始时刻(放大50倍)和最后演化时刻扰动能量的空间尺度特征,图 5a为2013年5月5日00时个例的结果,图 5b为月平均的结果。

图 5  2013年5月试验个例(a)和月平均(b)的前10个奇异向量平均总能量的谱分布 Fig. 5  The total energy spectra of averages for the first 10 leading SVs for a single case (a) and for one-month average (b) in May 2013

可以看出,无论是个例还是月平均特征,在初始时刻,奇异向量扰动总能量的大值分布在10—40波,最大值出现在20波左右(对应波长约为2000 km),在最后演化时刻,扰动能量主要分布在10—30波,最大值出现在12—13波(约3000 km),且能量分布更集中在最大值附近。对比两个时刻的能量谱特征,可以清楚地发现GRAPES奇异向量扰动能量在线性演变中也呈现出升尺度的能量转移特征。此外,GRAPES奇异向量扰动能量谱中出现最大值的水平尺度,也与其他研究(Coutinho, et al, 2004Zadra, et al, 2004Buehner, et al, 2006Lawrence, et al, 2009)中奇异向量扰动能量谱最大值分布尺度基本一致:初始时刻在20—22波附近出现能量峰值;而最后演化时刻能量峰值升尺度至12—14波。以上结果说明,文中GRAPES全球奇异向量扰动可以体现天气尺度的大气斜压不稳定能量及演变。

4 奇异向量扰动水平及垂直结构

本节从空间分布角度来分析GRAPES奇异向量的结构及演变特征。需要指出的是,基于总能量模的奇异向量扰动在空间分布上具有较强的局地性特征,扰动主要出现在天气系统斜压性较强的区域(Hoskins, et al, 2000)。为清楚地分析GRAPES奇异向量扰动的空间结构,初始时刻和最后演化时刻的GRAPES奇异向量扰动值均放大1000倍。

以2013年5月8日00时的边界层试验结果为例,选取北半球第1奇异向量扰动分布中比较显著的局地奇异向量扰动进行分析(图 6)。就影响天气系统而言,该区域受高压系统和极地低压的影响(图略),在初始时刻(图 6a),奇异向量位温扰动(模式第20层)出现在高压脊前较强西北气流控制区及高压脊与极地低压之间的位置。沿着图 6a中位温扰动中心做径向垂直剖面(图 6c,箭头所示起点为65°N,170°W;终点为55°N;90°W),可以看出,位温扰动主要集中在模式第15—28层,扰动大值在第20—25层(约700—500 hPa),且表现出显著的随高度向西倾斜的特征。这种随高度西斜的扰动结构为典型的基于总能量模奇异向量在初始时刻的特征(Buizza, 1994; Montani, et al, 2002; Coutinho, et al, 2004; Hoskins, et al, 2000),体现了大气对流层中下层常存在的、有利于天气系统发展的斜压结构。此外,其他研究结果也表明,初始时刻随高度向西倾斜的奇异向量扰动结构,在线性化演变过程中会逐渐消失,演变为较为垂直的正压结构。为此,需要分析GRAPES全球奇异向量位温扰动结构的演变是否具有这种变化特征?

图 6  初始时刻(a、c)和最后演化时刻(b、d)北半球第1奇异向量位温扰动在模式第20层的水平结构(a、b)及沿着箭头所示的径向垂直剖面(c、d) Fig. 6  Horizontal structures (a, b) at model level 20 and vertical cross-sections (c, d) (along arrows shown in a and b) of potential temperature perturbations of the 1st SV over NH at initial time (a, c) and final evolved time (b, d)

在最后演化时刻的位温扰动水平结构中(图 6b),可以看出,位温扰动迅速增长,并在背景轨迹气流的作用下移至下游地区,且扰动分布尺度变大,与上述讨论的GRAPES奇异向量扰动能量经过线性演化后向更大尺度转移的特征一致。从沿着扰动中心的径向垂直剖面(图 6d)可以看出,经线性演化后,位温扰动垂直结构发生了变化,由随高度向西倾斜的结构演变成了预期的、随高度无明显倾斜的正压结构特征,同时,扰动量值分布上也呈现出上传和下传的特征。

类似地,也显示了第1奇异向量中经向风(v风)扰动的水平及垂直结构(图 7)。为表现v风扰动水平结构的代表性特征,所选层次在初始时刻和最后演化时刻有所不同(初始时刻为模式第20层,最后演化时刻为模式第30层)。在初始时刻(图 7ac),v风扰动水平分布与位温扰动分布类似,v风扰动也主要集中在对流层中下层,同时也呈现出和位温扰动类似的随高度西倾的特征,这种特征在其他研究(Montani, et al,2002Coutinho,et al, 2004Lawrence, et al, 2009)也有提及。通过线性演变,v风扰动不仅水平尺度变大、强度增加(图 7b),在垂直方向也发展加强,分别在对流层上层及对流层低层出现极大值区(图 7d)。对比两个时刻v风扰动的垂直结构,可以清楚地看出v风扰动在线性演化过程中的上传和下传特征。此外,在最后演化时刻,v风扰动要明显强于位温扰动,这也从较具体的角度体现了GRAPES奇异向量扰动动能在线性演化中逐渐占主导的特征。

图 7  北半球第1特征向量经向(v)风扰动在模式第20层(a.初始时刻)和第30层(b.最后演化时刻)的水平结构及沿着箭头所示径向的垂直剖面(c.初始时刻,d.最后演化时刻) Fig. 7  Horizontal structures of v-wind perturbation of the 1st SV at model level 20 at initial time (a) and at model level 30 at finial evolved time (b) over NH, and vertical cross-sections along arrows (c. initial time, d. final evolved time)
5 总结和讨论

为发展GRAPES全球集合预报系统,在前期基于总能量模的GRAPES全球奇异向量研究工作的基础上,在GRAPES全球切线性模式和伴随模式2.0版的框架下,从两个方面进行了GRAPES全球奇异向量方法的改进。(1)通过引入线性化边界层方案来完善GRAPES全球奇异向量计算方案;(2)针对GRAPES奇异向量计算中主要耗时环节——伴随模式的计算效率,通过优化广义共轭余差方案的调用方式及采用大内存计算方法来提高奇异向量的计算效率。开展了连续试验,从奇异向量的扰动能量结构、扰动能量谱及扰动空间分布等方面,综合分析改进GRAPES全球奇异向量的结构及演变特征。得到如下结论:

(1) 改进后的GRAPES奇异向量扰动的线性演化特征合理,解决了之前GRAPES奇异向量扰动能量在近地面层不合理快速增长的问题。

(2) 综合分析表明,改进后的GRAPES奇异向量结构和演变特征符合典型的基于总能量模的奇异向量结构及演变特征:GRAPES奇异向量扰动能量主要位于对流层中层,其中位能占主导,在垂直结构上呈现出随高度向西倾斜的斜压特征,从水平尺度上,扰动能量集中在天气尺度;经过线性化演变,扰动能量在垂直结构上表现出向对流层高层上传及向对流层低层下传的特征,同时在水平尺度上也具有向较大水平尺度转移的特征。

(3) 通过优化伴随模式中广义共轭余差方案的调用方式及采用大内存计算方法,显著提高了GRAPES奇异向量计算效率,将计算时间较优化前缩短了25%,使得采用最优化时间间隔为48 h(适合集合预报应用)的GRAPES奇异向量的计算时间能够满足业务实时运行需求。

(4) 基于改进的GRAPES奇异向量方法所获得的奇异向量可以有效体现中高纬度地区GRAPES模式大气初始条件的斜压不稳定扰动,且计算效率显著提高,可应用于构建面向业务应用的GRAPES全球集合预报系统(相关工作内容将另文分析)。

需要指出的是,文中仅涉及针对中高纬度目标区GRAPES奇异向量的改进。下一步,将开展热带台风目标区GRAPES奇异向量技术及应用研究,拟从湿线性化物理过程对奇异向量结构的影响及如何将热带台风目标区奇异向量和中高纬度目标区奇异向量有机结合,共同应用于GRAPES全球集合预报系统等方面开展工作。

参考文献
刘永柱, 沈学顺, 李晓莉. 2013. 基于总能量模的GRAPES全球模式奇异向量扰动研究. 气象学报, 71(3): 517–526. Liu Y Z, Shen X S, Li X L. 2013. Research on the singular vector perturbation of the GRAPES global model based on the total energy norm. Acta Meteor Sinica, 71(3): 517–526. DOI:10.3969/j.issn.1004-4965.2013.03.020 (in Chinese)
刘永柱, 张林, 金之雁. 2017. GRAPES全球切线性和伴随模式的调优. 应用气象学报, 28(1): 62–71. Liu Y Z, Zhang L, Jin Z Y. 2017. The optimization of GRAPES global tangent linear model and adjoint model. J Appl Meteor Sci, 28(1): 62–71. (in Chinese)
刘永柱, 龚建东, 张林, 等. 2019. 线性化物理过程对GRAPES 4DVAR同化的影响. 气象学报, 77(2): 196–209. Liu Y Z, Gong J D, Zhang L, et al. 2019. Influence of linearized physical processes on the GRAPES 4DVAR. Acta Meteor Sinica, 77(2): 196–209. (in Chinese)
任迪生, 沈学顺, 薛纪善, 等. 2011. GRAPES伴随模式底层数据栈优化. 应用气象学报, 22(3): 362–366. Ren D S, Shen X S, Xue J S, et al. 2011. Optimized design of stack for GRAPES's adjoint model. J Appl Meteor Sci, 22(3): 362–366. DOI:10.3969/j.issn.1001-7313.2011.03.013 (in Chinese)
沈学顺, 苏勇, 胡江林, 等. 2017. GRAPES_GFS全球中期预报系统的研发和业务化. 应用气象学报, 28(1): 1–10. Shen X S, Su Y, Hu J L, et al. 2017. Development and operation transformation of GRAPES global middle-range forecast system. J Appl Meteor Sci, 28(1): 1–10. (in Chinese)
张林, 朱宗申. 2008. GRAPES模式切线性垂直扩散方案的误差分析和改进. 应用气象学报, 19(2): 194–200. Zhang L, Zhu Z S. 2008. Estimation of linearized vertical diffusion scheme in GRAPES model. J Appl Meteor Sci, 19(2): 194–200. DOI:10.3969/j.issn.1001-7313.2008.02.009 (in Chinese)
Bishop C H, Etherton B J, Majumdar S J. 2001. Adaptive sampling with the ensemble transform Kalman filter. Part Ⅰ:Theoretical aspects. Mon Wea Rev, 129(3): 420–436. DOI:10.1175/1520-0493(2001)129<0420:ASWTET>2.0.CO;2
Buehner M, Zadra A. 2006. Impact of flow-dependent analysis-error covariance norms on extratropical singular vectors. Quart J Roy Meteor Soc, 132(615): 625–646. DOI:10.1256/qj.05.66
Buizza R. 1994. Sensitivity of optimal unstable structures. Quart J Roy Meteor Soc, 120(516): 429–451. DOI:10.1002/(ISSN)1477-870X
Buizza R, Palmer T N. 1995. The singular-vector structure of the atmospheric global circulation. J Atoms Sci, 52(9): 1434–1456. DOI:10.1175/1520-0469(1995)052<1434:TSVSOT>2.0.CO;2
Coutinho M M, Hoskins B J, Buizza R. 2004. The influence of physical processes on extratropical singular vectors. J Atmos Sci, 61(2): 195–209. DOI:10.1175/1520-0469(2004)061<0195:TIOPPO>2.0.CO;2
Diaconescu E P, Laprise R. 2012. Singular vectors in atmospheric sciences:A review. Earth-Sci Rev, 113(3-4): 161–175. DOI:10.1016/j.earscirev.2012.05.005
Hartmann D L, Buizza R, Palmer T N. 1995. Singular vectors:The effect of spatial scale on linear growth of disturbances. J Atmos Sci, 52(22): 3885–3894. DOI:10.1175/1520-0469(1995)052<3885:SVTEOS>2.0.CO;2
Hong S Y, Pan H K L. 1996. Nonlocal boundary layer vertical diffusion in a medium-range forecast model. Mon Wea Rev, 124(10): 2322–2339. DOI:10.1175/1520-0493(1996)124<2322:NBLVDI>2.0.CO;2
Hoskins B, Buizza R, Badger J. 2000. The nature of singular vector growth and structure. Quart J Roy Meteor Soc, 126(566): 1565–1580. DOI:10.1256/smsqj.56601
Houtekamer P L, Mitchell H L. 2005. Ensemble Kalman filtering. Quart J Roy Meteor Soc, 31(613): 3269–3289.
Laroche S, Tanguay M, Delage Y. 2002. Linearization of a simplified planetary boundary layer parameterization. Mon Wea Rev, 130(8): 2074–2087. DOI:10.1175/1520-0493(2002)130<2074:LOASPB>2.0.CO;2
Lawrence A R, Leutbecher M, Palmer T N. 2009. The characteristics of Hessian singular vectors using an advanced data assimilation scheme. Quart J Roy Meteor Soc, 135(642): 1117–1132. DOI:10.1002/qj.v135:642
Mahfouf J F. 1999. Influence of physical processes on the tangent-linear approximation. Tellus A, 51(2): 147–166. DOI:10.3402/tellusa.v51i2.12312
Molteni F, Buizza R, Palmer T N, et al. 1996. The ECMWF ensemble prediction system:Methodology and validation. Quart J Roy Meteor Soc, 122(529): 73–119. DOI:10.1002/(ISSN)1477-870X
Montani A, Thorpe A J. 2002. Mechanisms leading to singular-vector growth for FASTEX cyclones. Quart J Roy Meteor Soc, 128(579): 131–148. DOI:10.1256/00359000260498824
Mu M, Duan W S, Wand B. 2003. Conditional nonlinear optimal perturbation and its application. Nolin Processes Geophys, 10(6): 493–501. DOI:10.5194/npg-10-493-2003
Mureau R, Molteni F, Palmer T N. 1993. Ensemble prediction using dynamically conditioned perturbations. Quart J Roy Meteor Soc, 119(510): 299–323. DOI:10.1002/qj.v119:510
Palmer T N, Molteni F, Mureau R, et al. 1993. Ensemble Prediction//Proceedings of the ECMWF Seminar on Validation of Models over Europe. Reading, UK: Shinfield Park, 21-66
Palmer T N, Gelaro R, Barkmeijer J, et al. 1998. Singular vectors, metrics, and adaptive observations. J Atoms Sci, 55(4): 633–653. DOI:10.1175/1520-0469(1998)055<0633:SVMAAO>2.0.CO;2
Toth Z, Kalnay E. 1997. Ensemble forecasting at NCEP and the breeding method. Mon Wea Rev, 125(12): 3297–3319. DOI:10.1175/1520-0493(1997)125<3297:EFANAT>2.0.CO;2
Zadra A, Buehner M, Laroche S, et al. 2004. Impact of the GEM model simplified physics on extratropical singular vectors. Quart J Roy Meteor Soc, 130(602): 2541–2569. DOI:10.1256/qj.03.208