逐层平滑地形的平缓地形追随坐标在高分辨率GRAPES模式中的应用研究*

李 超 陈德辉 李兴良 胡江林

国家气象中心，北京，100081

1 引 言

地形对大气运动具有重要的动力和热力作用。为了提高预报准确率，数值预报模式需要正确地描述地形对天气系统的影响。模式的垂直坐标与地形直接相关：恰当的垂直坐标可以使复杂的控制方程组转化为较简单的形式，下边界条件简单易给，便于使用精度和效率较高的计算技术，有效地减 小 气压梯 度 力 （Pressure Gradient Force，简称PGF）以及其他与模式面水平差分计算有关项的计算误差，从而正确描述模式中地形对大气的作用。

为了提高数值预报模式的模拟与预报能力，提高模式水平分辨率是非常必要的。然而，随着数值预报模式水平分辨率的不断提高，模式地形越来越精细，模式能分辨的地形越来越陡峭，垂直坐标面的坡度越来越大。采用传统的地形追随垂直坐标（Gal-Chen，et al，1975；以下简写为 Gal.C.S 坐标）引起的气压梯度力计算误差会相应增大，模拟的地形重力波形状也被严重扭曲（Schär，et al，2002），对高分辨率数值模拟带来不利影响。

中国气象局自主研发的GRAPES模式目前采用传统的高度地形追随坐标其中，Zs表 示地形高度，z表示坐标面绝对高度，zˆ表示坐标面相对高度，ZT为模式层顶高度。这种坐标形式可以使模式计算在一个“矩形”的有界域里进行，便于计算机编程自动计算；可以使实际大气中不规则的下边界变成光滑的、随地形起伏的模式下边界面，易于给定模式的下边界条件，也更易于描述气流随地形起伏的运动。地形追随坐标使得模式光滑的下边界面成为地表（陆面、海面）边界层的上边界面，便于模式动力、热力过程和地表边界层过程参数化方案的耦合（Phillips，1957），地形的动力和热力作用可以更准确地传递给上层大气。但是，Gal.C.S坐标的问题也很突出。地形追随坐标是非严格意义上的正交坐标（Zdunkowski，et al，2003），是由非地形追随坐标向地形追随坐标转化时，气压梯度力计算由一项变成了量级相当的两项的差，带来不可避免的气压梯度力项计算误差。地形追随坐标的坐标面随地形起伏，低层“隆起”的坐标面会一直伸展至模式层顶部。虽然坐标面的坡度随高度升高会有所降低，模式大气沿高层起伏模式面的运动与实际大气高层的准水平运动仍存在很大的偏差。随着GRAPES模式分辨率的提高，Gal.C.S坐标可识别的地形坡度增大，地形重力波会被严重歪曲，出现波动扭曲、破碎的现象。

Schär等（2002）提出平缓坐标（smoothed-level coordinate，简称SLEVE坐标）的概念，在不改变高度地形追随坐标定义的前提下，对坐标形式做变形，通过随高度变化的地形衰减系数控制各个模式面上的地形影响作用的大小，并设计试验比较了单尺度SLEVE1坐标和双尺度SLEVE2坐标，发现两种坐标有效地减小了地形引起的各种模拟偏差，SLEVE2坐标对缓解小尺度的模拟噪音效果更佳。李超等（2012a）选取不同的地形衰减函数，进行了一系列理论分析和理想试验，试验结果展示了平缓-混合坐标在减小气压梯度力计算误差、平流耗散等方面的模拟优势。李超等（2012b）基于GRAPES模式选取Gal.C.S坐标和平缓-混合SLEVE1坐标进行了初步的模式设计并进行了简单的个例模拟检验分析，SLEVE1坐标对预报量的小尺度噪音有所衰减，各个预报量的检验误差也有所减小。Li等（2015）基于 GRAPES模式选取 Gal.C.S坐标、SLEVE1坐标、SLEVE2坐标和一种以余弦函数为衰减基函数的COS坐标进行了理想试验和实际批量模拟试验。各种形式平缓-混合坐标相比较发现，SLEVE2坐标对计算误差衰减最大，较好地提高了预报能力。张旭等（2015）将 Klemp（2011）提出的一种随高度逐层平滑地形作用的平缓-混合坐标形式（下文简称Klemp坐标）应用于GRAPES区域模式，通过理想试验和简单个例模拟试验对这种坐标的优势进行了初步的探讨。李超等（2019）设计了一种新的COS坐标，有效解决了Li等（2015）中COS坐标理想试验结果最优而实际模拟结果较差的问题，新的COS坐标在底层设计也更加灵活。李超等（2020）比较了几种平缓-混合地形追随坐标对典型个例模拟的影响，发现平缓-混合坐标可以有效缓解GRAPES模式南风偏大、虚假降水和虚假天气系统等问题。

随着GRAPES模式分辨率进一步提高，GRAPES-3km模式中地形可识别的最大高度接近7 km，而上述基于 GRAPES-10km和 GRAPES-15km的模式研究工作最大识别地形均在6 km以下，更高分辨率下的平滑-混合坐标模拟试验还需要进行重新设计并检验结果。另外，张旭等（2015）对Klemp（2011）坐标的研究工作显示，在较低分辨率的GRAPES模式中，该坐标的个例模拟试验有较好的预报改进效果，而在分辨率较高的GRAPES-3km模式中，该坐标的批量预报试验效果还需要进一步的研究。本研究选取前期研究效果较好的SLEVE1坐标与Klemp坐标基于GRAPES-3km模式进行一系列的比较研究工作，为高分辨率GRAPES模式动力框架垂直坐标的选择提供一定的参考。

2 坐标介绍

李超等（2012a）给出通用坐标形式

式中，L为地形尺度谱的数目，为地形临界衰减高度，bh∗i为地形衰减系数，ZSi为地形高度，z为坐标面绝对高度，为坐标面相对高度。

对于GRAPES-3km模式目前采用的Gal.C.S坐标有对于单尺度平缓SLEVE1坐标有地形作用在高度上衰减为地形的1 / e。

Klemp （2011）提出在 COS坐标（李超等，2012a）基础上，对地形进行逐层平滑，得到坐标形式（简称Klemp坐标）。

式中，

式中，Zc为 地形作用衰减为0的高度，ZSk为每层坐标面上地形作用的高度。Klemp坐标的地形滤波方案满足图1所列出的计算流程，其中m是每一层的平滑滤波次数，k是垂直层数。每一层都以上一层经过m次滤波得到的地形为基础进行滤波。滤波采用拉普拉斯的二阶滤波方法

图1 Klemp 坐标地形逐层滤波方案Fig.1 The terrain filter scheme of Klemp coordinate

式中，i和j为水平格点数。滤波系数 βk是一个分段系数：表示坐标面高度，在2倍的地形高度以下滤波系数随高度升高线性增大，2倍的地形高度以上则为常数。

Klemp坐标与其他平缓-混合地形追随坐标差别在于地形补偿项中的ZSk是一个三维变量。需要注意的是，并不是把地形变成了三维变量，而是坐标面上的地形作用随高度是变化的。

高分辨率模式中可识别的地形更加精细。例如，GRAPES-15km模式可识别的地形最大高度是5800 m左右，而GRAPES-3km模式可识别的最大高度达到6800 m（图2）。3 km分辨率模式垂直分为51层，比15 km模式的31分层更细，如果使用GRAPES-15km模式中平缓-混合坐标的地形临界衰减高度（Li，et al，2015），GRAPES-3km 模式中地形较高点上空的模式面厚度就会过薄，导致计算不稳定。GRAPES-3km模式中平缓-混合坐标的设计应用面临新的问题，Klemp坐标对坐标面的进一步平滑会使这种计算不稳定更容易出现。下面的实际模拟试验对坐标的参数设计做了试验分析。

3 基于 Klemp 坐标的动力框架方程推导

高度地形追随坐标中地形作用对模式动力框架的影响主要体现在气压梯度力、垂直速度、散度等物理量的计算上。下面就基于Klemp坐标对垂直速度、散度等物理量做进一步推导。

图2 地形高度分布 （a.GRAPES-15km 模式，b.GRAPES-3km 模式；单位：m）Fig.2 Distribution of terrain height （a.GRAPES-15km model，b.GRAPES-3km model；unit：m）

（1）垂直速度

则：

又，

式中，z=Zs和z=ZT时 ，wˆSk=0。

得出：

式中，

又，

（2）水平气压梯度力

式中， φSk(x,y,z)=gZSk(x,y,z)为坐标面的位势高度。

（3）球面坐标下的水平散度

对于三维散度中的垂直速度的垂直变化，有：整理得

式中，

则

于是，三维散度（D3）可表示为

综合上述公式推导可以看出，与Li等（2015）中动力框架方程推导相比，基于Klemp坐标的模式动力框架中地形以及地形偏导数有关的项有所改变，动力框架中涉及到垂直速度、散度等的线性项、非线性项以及亥姆赫兹方程中相关的系数也要做相应的改变（薛纪善等，2008），这里不再做详细介绍。

4 理想试验和实际试验

下面选择GRAPES-3km模式目前采用的Gal.C.S坐标、前期研究模拟结果较好的SLEVE1坐标以及改进的Klemp坐标进行一系列模拟试验，对比分析改进的Klemp坐标在GRAPES-3km模式中的模拟性能。

4.1 地形重力波理想试验

低层大气中不同尺度的山脉会激发出形状尺度不同、传播方向各异的重力波。山脉重力波在大气底层的结构十分复杂，这也是引发天气状况的重要原因。下面是层结稳定的干空气穿过二维地形激发出山脉重力波的理想试验。上游大气廓线设置为：布朗特数N= 0.01，x方向水平速度U=10 m/s，上游地表面参考位温为280 K。

钟形地形叠加上小尺度的波动形成理想地形，其函数形式满足

试验计算区域设置为：（−25000，25000）m×（0，21000）m，Δx=250m， Δ ξ =210m，下边界为无通量的边界条件，上边界为9 km的吸收边界，侧边界为开放边界。无物理过程和扩散的条件下积分10 h，比较各种坐标的垂直风速的分布。

从积分模拟结果（图3）可见，平缓坐标可以有效缓解地形追随坐标带来的重力波变形破碎的问题。Klemp坐标模拟的重力波（图3d）形状更接近解析解（图3a），重力波变形的情况在中高层基本消失。

4.2 月连续试验

理论分析和理想试验中改进的Klemp坐标可以有效减小坐标面地形影响，模拟出更接近解析值的重力波，而在实际应用中地形更加复杂，物理和动力耦合起来后，坐标参数设置考虑的因素也更加复杂。下面选用GRAPES-3km模式对2018年9月1—30日进行24 h滚动预报，检验Klemp坐标在长期的批量预报试验中的改进效果。该试验水平分辨率为0.03°×0.03°，垂直方向分为51层；积分时间步长为30 s；试验初始场和边界场由分辨率为1°×1°的NCEP再分析资料提供，不做资料同化；模拟范围为（10°—60°N，70°—145°E），选用一维大气参考廓线。微物理过程采用WMS6方案，长波辐射过程采用RRTM方案，短波积云辐射过程采用SWRAD方案，地表层物理方案采用SFC方案，陆面过程采用NOAH方案，边界层参数化采用NMRF方案，无积云参数化方案。

GRAPES-3km模式可识别的地形最大高度接近7 km，较小的临界衰减高度（h∗）易造成计算不稳定，所以 SLEVE1坐标h∗取值 15 km，这相对于GRAPES_15km中SLEVE1坐标h∗有所增大（Li，et al，2015）。改进的 Klemp 坐标相对于 Klemp（2011）以及张旭等（2015）的参数设计也有所不同。考虑到3 km模式可识别的地形坡度较大、中国区域东西地形高度跨度较大、同一经纬度坐标面之间厚度差异较大也会带来计算误差，实际试验中改进的Klemp坐标各层平滑系数取每层平滑14次，衰减系数（bh∗）中zc取值30 km，参数n取 值3，z0暂时取值为0。可以看出这里平滑系数较 Klemp（2011）以及张旭等（2015）有所提高，平滑次数有所减少，衰减参数相对李超等（2019）也有所改变。参数的设计目前还未找到定量的标准，在最大化减小坐标面上地形影响作用的前提下，经过多种尝试获得1组能够较好体现预报效果的参数。在上述的参数设置下，3种坐标的坐标面分布如图4所示。Klemp坐标在SLEVE1坐标的基础上对坐标面上的小尺度地形进一步衰减，由此带来的气压梯度力计算误差和平流输送误差也会进一步减小（李超等，2012a）。

图3 稳定气流过山积分10 h三种坐标的垂直速度分布与解析值比较（a.解析值，b.Gal.C.S 坐标，c.SLEVE1 坐标，d.Klemp 坐标；等值线间隔为 0.05 m/s；时间步长为 0.3 s；横坐标是偏离中心位置的距离，纵坐标是距离地面的高度）Fig.3 Cross-section of vertical velocity after 10 hours of integration for steady-state flow over belle-shaped mountain（a.analytic solution，b.with Gal.C.S，c.with SLEVE1，d.with Klemp；the contour interval is 0.05 m/s；the abscissa denotes the distance away from the domain center，the coordinate denotes the vertical level height above the ground surface）

批量试验结果显示，与FNL再分析资料相比，SLEVE1和Klemp坐标的高层温度场24 h预报平均偏差有较明显减小。100和250 hPa温度场偏差（图5）相比GRAPES-3km模式中目前的Gal.C.S坐标在青藏高原和蒙古国上空有两个较集中的减小区域。沿28°N剖面的V风场平均偏差（图6）显示平缓坐标在高层的偏差有较明显的减小，特别是平滑地形的Klemp坐标在100 hPa以上几个偏差大值中心均消失。为分区域检验降水预报结果，将中国划分为8个区域（图7）分别检验24 h降水累加ETS评分（图8），结果显示：从中国区域来看，Klemp坐标对中雨、大雨、暴雨的评分较Gal C.S坐标有一定程度提高，也整体优于SLEVE1坐标。对受青藏高原大地形影响较明显的长江中下游、华南、西南东部地区，Klemp坐标小雨评分较高，可能与其对小尺度地形的平滑作用有关。其他降水量级的评分没有特别显著的规律。但是在西南东部地区，平缓混合坐标的评分几乎全部优于Gal.C.S坐标，特别是Klemp坐标的降水评分有较大程度地提高。这也说明平缓混合坐标带来的计算误差对西南东部地区即青藏高原下游的云贵高原和四川盆地等地的降水有较大影响。

图4 GRAPES-3km 模式中 3 种坐标沿 28°N 的模式面分布 （a.Gal.C.S 坐标，b.SLEVE1 坐标，c.Klemp 坐标）Fig.4 Coordinate surfaces of the GRAPES-3km model along 28°N （a.for Gal.C.S，b.for SLEVE1，c.for Klemp）

图5 3 种坐标 （a.Gal.C.S 坐标，b.SLEVE1 坐标，c.Klemp 坐标） 批量试验 24 h 预报100 hPa （a1—c1） 和 250 hPa （a2—c2） 温度平均偏差场 （单位：K）Fig.5 Monthly averaged 24 h forecast bias of temperature at 100 hPa （a1—c1） and 250 hPa （a2—c2） against analysis with （a） the Gal.C.S，（b） the SLEVE1 and （c） the Klemp coordinate （unit：K）

图6 3 种坐标 （a.Gal.C.S 坐标，b.SLEVE1 坐标，c.Klemp 坐标） 批量试验 24 h 预报沿 28°N 剖面 V 风速的平均偏差场 （单位：m/s）Fig.6 Longitude-height cross sections along 28°N of monthly averaged 24 h forecast bias of V wind against analysis with（a） the Gal.C.S，（b） the SLEVE1 and （c） the Klemp coordinate （unit：m/s）

图7 降水评分中国区域划分 （1.东北地区，2.新疆地区，3.西北地区，4.华北地区，5.西藏地区，6.西南地区东部，7.长江中下游地区，8.华南区）Fig.7 Subregions in China for precipitation ETS score（1.Northeast China，2.Xinjiang region，3.Northwest China，4.North China，5.Tibet region，6.Southwest China，7.the middle-lower reaches of Yangtze River，8.South China）

图8 2018 年 9 月 1—30 日 24 h 降水预报累加检验平均 ETS 评分（a.中国区域，b.长江中下游区域，c.华南区域，d.西南东部区域）Fig.8 Precipitation ETS scores for 24 h forecast（a.entire China，b.the middle-lower reaches of Yangtze River，c.South China，d.Southwest China）

5 结 论

平缓地形追随坐标可以通过调整地形衰减系数来约束坐标面上的地形作用，在这个基础上进一步对逐层坐标面进行平滑又可以将小尺度的地形影响剔除，得到更加平滑的坐标面，进一步减小各种计算误差。

文中将逐层平滑地形的平缓地形追随坐标应用到高分辨率的GRAPES-3km模式中进行理想试验和批量模拟试验，得到如下的初步结论：

（1）相对于其他平缓地形追随坐标，逐层平滑地形的平缓地形追随坐标地形重力波模拟更接近解析值，波形变形和破碎问题进一步缓解。

（2）月连续模拟试验中逐层平滑地形的平缓地形追随坐标对高层月平均的温度场、风场模拟误差有一定程度减小，中国范围内月平均的降水评分也有所提高，西南地区东部的评分提高更加明显。

（3）无论是普通的平缓地形追随坐标还是平滑地形的平缓地形追随坐标，参数的设置都至关重要。高分辨率的GRAPES-3km模式对中国范围的模拟对坐标参数的设置更加敏感。保证计算稳定、追求坐标面平滑、得到较好的模拟结果三者同时达到需要对参数反复试验调整。

平滑地形的平缓地形追随坐标是目前国际上最先进的平缓地形追随坐标形式，也在高分辨率的GRAPES模式中得到了较好的应用结果。该坐标对模式分辨率、模式分层等有较高的依赖，所以在对业务模式进行改进升级时，也需要同步地调整该坐标的参数设置。