MM5模拟垂直风剖面的准确性探讨及在垂向气象场的分析应用*

唐 灵，詹杰民，李毓湘，林嘉仕

(1．中山大学应用力学与工程系，广东 广州 510275；2．香港理工大学土木及结构工程系，香港)

MM5[1](第五代的NCAR/Penn State中尺度模型)是由MM4经过大量的改进发展而来。MM5 在处理模式大气降水物理过程时有多种方案可选择，这些方案通常分为两类，分别称作显式方案和隐式方案。其中显式方案处理模式可分辨尺度(或称为大尺度)降水物理过程，隐式方案处理模式不可分辨尺度(或称积云对流)降水过程[2]。MM5中共有5种显式云物理方案, 从简单的暖雨方案、简单冰相方案, 到物理过程比较详尽的Goddard、Reisner霰方案[3]。积云与大气的相互作用非常重要，常用的对流性积云参数化方法有Anthes-Kuo[4]、Grell[5-6]、Arakawa和Schubert[7]、Fritsch和Chappell[8]、 Kain和Fritsch[9]等方案。MM5的行星边界层的参数化有几种不同的处理,主要有：把行星边界层看作是由近地层和混合层组成的整体的总体边界层方案、不同情况下采用k理论或根据自由对流特征建立传输模式的高分辨边界层参数化方案、加入了反梯度传输项的MRF[10]方案及K与湍能有关的能量闭合方案等[11]。

本文着重讨论了对MM5垂直方向的σ层进行调整及加密方案以及各种方案对模拟的垂直风剖面的影响，并在此基础上讨论和分析了土地利用的更新对垂直风场造成的变化，同时分析了珠三角地区2001年9月15-18日空气污染时段物理气流输运过程的一种影响因素。

1 控制方程及边界条件

MM5所用到的基于沿地形的坐标(x,y,σ)的非静力模式基本方程主要有：

1)压力方程

ρ0gw+γp▽·V=

(1)

其中p′为压力扰动，在模式中，实际忽略了右边括号里的最后一项，它表示由于热膨胀效应而导致气压的增加。

2)动量方程

x方向分量

·▽u+

(2)

y方向分量：

·▽v-

(3)

σ方向分量：

·▽

(4)

3)热力学方程

·▽T+

(5)

另外还有水汽守恒方程和针对微物理变量(如云和降雨等)的预报方程，它们包括平流项和不同的源/汇项。在空间差分水平方向上MM5采用Arakawa B-grid上的跳点格式进行有限差分；在时间差分上采用二阶蛙跳式时间步长方案，对某些项采用时间分裂方式，这样可以提高计算效率。

MM5模式系统设计了固定侧边界、时变侧边界和流入/流出松驰侧边界等几种侧边界处理方案。对于非静力平衡系统中，选用流入/流出松驰侧边界条件能取得较好的模拟效果。模式上边界采用Klemp等设计的一种上边界处理方法，这种方法能使波能穿过上边界而不产生反射[2]，具体表达为

(6)

2 模拟的参数化设置及模拟时段描述

2.1 MM5模式及参数化设置

模拟采用的是MM5 v3.6版本，非静力模式,使用双向反馈的三重嵌套，分辨率分别为40.5，13.5和4.5 km,垂直方向上采用σ坐标。最粗的网格积分步长为60 s，1 h输出一次预报结果。本次模拟可分辨尺度降水采用加入了霰和冰的数量浓度预报方程的Reisner 2方案，边界层参数化方案为高分辨率的MRF方案和适用于高分辨率、采用4种稳定度方式、包括自由对流混合层的Blackadar方案；辐射采用云辐射方案；粗网格侧边界条件采用松弛流入流出方案，侧边界缓冲区采用了4个格点；地面温度取5层土壤模式。植被/土地类型采用USGS(U.S. Geological Survey)的25类全球资料。气象场采用NCEP再分析6 h间隔1°×1°分辨率的数据作为初估场。在MM5积分过程中，进行表面和高空场的牛顿松弛格点数据同化分析。

2.2 模拟时段天气描述

本次进行了2个时间段的模拟工作，第一个时间段为2001年6月16日早上8∶00至18日早上8∶00(北京时，下同),在模拟时段内,华南地区主要受西南季风的影响,期间有零散的雨量,同时也受到向中国东南部伸展的副热带高压脊的影响出现高温晴朗天气。模拟时段内珠江三角洲地区主要盛行从海面吹向陆地的偏南向陆风。

第2个时间段为2001年9月15日8∶00至18日8∶00，根据香港天文台的气象资料，从9月15日开始，一股干燥大陆性气流开始影响华南地区，香港地区天晴及炎热。在轻微偏北风的情况下，香港大部分地区早上有烟霞，能见度下降至5 000 m以下，西部地区更下降至2 000 m以下。炎热及有烟霞的天气持续至9月19日。

2.3 垂直σ层调整和加密

MM5需要从NCEP或ECMWF等全球范围的气象场再分析数据读取模拟所需的数据作为第一猜场，通过INTERPF模块把气压层上的初始场数据插值到MM5所采用的σ层上，然后生成底、周边界条件和初始条件以供MM5进行模拟运算。通过调整INTERPF插值σ分层的值，可以对MM5原有的25个垂直层进行调整和加密，由于近地面行星边界层大气变化最活跃且模拟意义较大，因此加密的原则是越靠近地面越密，离地面越高越稀疏，并且尽量均匀分布。

3 MM5垂直层调整及加密模拟结果对比分析

3.1 原垂直方向25层及1 000 m以下近地面垂直层加密模拟结果

在原垂直方向25层的基础上，调整了1 000 m以下近地面层的密度。原25层方案1 000 m以下为11层，1 000 m以上为14层，调整后1 000 m以下为16层，1 000 m以上为9层。

模拟结果在2001年6月17日14∶00输出，模拟输出点gz1位于地势低平的城市市区(113°15′50″E，23°6′49″N),两种垂直层方案2 000 m以下中低空垂直风剖面模拟结果分析见图1和2。

图1 原25层方案垂直风剖面(gz1点，下同)

图2 调整后25层方案垂直风剖面

从gz1点两种方案垂直风速和垂直风向模拟的对比结果可以看到，不加密情况下风速从地面到2 000 m中低空的变化是呈一种类抛物线式的递增，其中越靠近地面增速很快，高度越高增速也越慢，风向角度从地面到高空也是呈现一种递增的趋势，反映出随着高度的增加，风向呈顺时针旋转，这是由于受到地转偏向力和地表摩擦力共同作用造成的。随着高度的增加，地表摩擦力对风向的影响越来越小，因而风速增加也变小，从800～2 000 m高度区间，风速几乎不再增加。而在加密的情况下，从地面到800 m高度，风速也是迅速增大，随高度增加风速的增加也越慢，这与原25层方案的变化趋势类似，然而从800～2 000 m的区间，风速继续小幅增加，从800 m高度的4 m/s增加到2 000 m高度的5 m/s，这与原25层的方案表现出了不同的特征，而且加密方案模拟出来的地表10 m高度风速比原方案高0.5 m/s,在2 000 m高度则比原方案约高1 m/s，而且因为新方案模拟的风速普遍比原方案高，风向偏转会比原方案小2°～5°。

我们知道，在中高纬西风带内或在低纬度地区都会出现高低空急流，急流是大气环流中在风场上的一个重要现象。在我国经常出现的是西南风低空急流[12]，一般出现在1 500～3 000 m的中低空。从天气角度来讲,一般指在850或700 hPa上最大风速中心超过一定数值(如12或15 m/s等)的气流，其平均长度1 000～2 000 km左右，宽度约数百km。对我国有较大影响的低空急流—西南风低空急流，常在华南前汛期和江淮梅雨期间出现。高空急流[13]是指出现在对流层顶附近或平流层中的一股强而窄的气流。一般指风速大于或等于30 m/s的强风带。因此在此次模拟的时段和地区范围内，正好受到了西南风低空急流影响，模拟时段内，珠江三角洲地区主要盛行从海面吹向陆地的海风，为偏南风，而到了1 000 m以上的高度，风场受到西南风低空急流的影响，风速将会加速增大，这个趋势在加密方案的情况下被模拟了出来，这说明加密方案使得模拟的准确性得到了提高。

为了进一步说明近地面层加密对气象场模拟准确度的影响，下面通过香港地面气象站KPS(京士柏站)的实测气象数据与数值模拟值的对比来进行分析。

从KPS站观测温度与模拟值的对比分析(图3(a))可以知道，近地面层加密后的模拟温度曲线变化趋势和峰值都比原25层方案的模拟温度曲线(虚线)更接近实测值(实线)，这说明近地面层加密后对温度模拟的准确性得到了提高。从KPS站观测风场与模拟风场的时间序列对比(图3(b))也可以看到，近地面层加密后模拟的风速大小和风向也更接近观测值，特别是在模拟的第15 h后到模拟结束这一阶段。

以上对比分析表明，在对近地面层进行加密之后，提高了2 000 m 以下低空气象场模拟的准确度，下面分析这2种方案下2 000～12 000 m区间垂直风场的模拟情况。

图3 (a) KPS站观测温度与模拟值的对比；(b) KPS站观测风场与模拟风场的时间序列对比

图4 (a) 原25层方案风速随高度分布；(b) 调整后25层方案风速随高度分布

从图4(a)和4(b)可以看到，由于加密后的方案增加了近地面层的密度，而总的层数还是25层没有变，因而高空层数比原方案更少，虽然模拟出来的2 000 m以上风速变化曲线大致变化趋势还是跟原方案接近，继西南风低空急流影响造成的2 000～4 000 m风速增大区间后，4 200 m以上风速开始减小，而到了12 000 m高空，受到高空气流的影响，风速将再次增大，显然由于高空垂向层数更密，原25层方案对高空气流的模拟比近地面层加密方案更准确，因此，有必要增加垂直方向总的σ层数。

3.2 增加垂向σ层总数对垂向风剖面模拟结果的影响

上面的分析结果表明，近地面加密后，得到的垂向风剖面模拟结果更符合实际情况，说明对近地面气象场的模拟结果有较大的改善作用，这对于研究近地面风剖面的变化有较大的帮助，因为通常我们对高空风剖面的关注会比较少。然而由于近地层的增加，高空层则变得更稀疏，那么高空垂向风剖面的变化则得不到很好的描述，因此，如果要研究中高空的气象场，那么就要增加总的垂直层数，当然这会带来更大的计算量。下面分别把垂直方向的总层数从25层增加到35层、40层和50层，并分析由此而带来的模拟垂直风剖面的变化。

图5 模拟垂直风速分布

图5(a)和5(b)分别是垂直方向25层和35层模拟输出的风速随高度变化曲线。25层方案的垂直层分配为2 000 m以下19层，2 000 m以上6层。35层方案的垂直层分配为2 000 m以下24层，2 000 m以上11层。从两者的模拟垂直风速结果对比可以看到，相比25层方案，35层方案在2 000 m以下的中低空和2 000 m以上的中高空都加密了5个垂直层，从而使得模拟风速的变化曲线发生了较为明显的变化。在大约500～2 000 m的高度区间，35层方案模拟的风速随高度变化比25层方案的结果要平缓，这说明在此高度区间，35层方案模拟结果反映出地表摩擦力对风速的影响已经较小，地表摩擦力对风速的影响主要集中在500 m以下的近地表层。而在中高空的5 000 m左右高度，35层方案模拟的风速达到一个极大值约12.5 m/s，25层方案模拟结果则是在4 000 m左右达到一个极大值约11 m/s。再看8 000 m以上的区间，25层方案由于在此区间只有一个垂直层了，因此模拟结果的细节严重不足，只能反映出一个风速增大的趋势，而到19 000 m以上的高空已经没有σ层了。35层方案的模拟结果显示，从8 000 m向上直到13 000 m左右，风速有一个小幅增加的过程，再往上风速回落，在16 000 m以上直到28 000 m区间，风速又会平稳增加。

下面再把垂直方向的σ层增加到40层和50层，并把35层、40层和50层的模拟结果进行对比分析。

图6 模拟垂直风速分布

图6(a)、6(b)和6(c)分别是35层、40层和50层的垂直方向风速随高度变化曲线。35层的垂直层分配为2 000 m以下24层，2 000 m以上11层，40层的垂直层分配为2 000 m以下27层，2 000 m以上13层，50层的垂直层分配为2 000 m以下25层，2 000 m以下25层。从35层和40层的风速分布对比可以看到，这两种情况下各个高度层的风速变化趋势及风速大小都基本吻合，只是在5 000～8 000 m的高度区间，40层的模拟垂直风速曲线出现了小幅波动，不过风速变化趋势跟35层的结果基本一致。继续增加垂直层到50层，从50层的模拟结果来看，20 000 m以下的风速变化曲线跟35层和40层的结果也是基本一致的，20 000～30 000 m的区间，从50层的模拟结果可以看到风速有一个先减小、后增加的过程，而35层和40层的模拟结果只反映出风速的增大趋势。由于高空垂直层的进一步增加，50层方案模拟出了更细致的风速变化趋势，且30 000 m以上高空的风速变化曲线也反映了出来，而在这个区间35层和40层方案已经没有垂直层故而无法显示。

从35层、40层和50层的风速分布曲线对比分析我们知道，对于20 000 m以下区间的风速模拟，这3种情况的模拟结果基本上趋于一致和稳定，进一步增加垂直层数已经不会对模拟结果产生较大的影响，这说明了加密垂直层到35层以上得到的模拟结果是可信的，且加密后的模拟结果与原来25层的模拟结果差异较大，这也说明加密低空和高空的垂直层对模拟的结果能起到较大的改善和提高作用。MM5模式里考虑了不同高度层的对流和上升下沉气流的卷入卷出及补偿运动浮力能等，低中高空不同高度各层之间有一个相互影响的过程，改善高空气象场的模拟精度也会对近地面层的模拟结果产生影响。然而垂直层也并不是越密越好，通过前面的讨论和分析，在选取的模拟工况下，2 000 m以下采用20～25层，2 000 m以上采用15～25层已经可以得到较好的垂直风剖面模拟结果。在敏感性测试过程中，我们发现，如果低空的垂直层太密，反而会使得模拟的风速变化趋势出现扭曲的失真现象。另外，垂直层数越多，计算量也会大幅增加，且由于计算能力的限制，如果水平网格较多，增加垂直层数，在实际模拟计算过程中可能产生不稳定从而导致程序溢出中断的现象。一般来讲，对于20 000 m以上高空的气象场关注较少，较多关注中低空和近地面气象场的数值模拟，所以综合而言，在模拟的过程中，垂直方向采用35～40层是比较合适的。

为了对模拟的结果进行验证，图7给出了香港地区一月份的垂直风廓线实测结果。由于垂直风廓线观测数据较少，给出的观测数据时段与模拟结果不在同一时段，观测资料的风速比模拟结果整体上大，不过风速大小随高度的变化规律与模拟结果基本相似。图7中，600 m以下风速随高度变化迅速增大，700～3 500 m高度区间，风速变化较小，4 000～11 000 m风速随高速变化是一个先加速增大，然后减小的过程，11 000 m以上风速随高度再次增大然后减小。这种变化规律与采用35和40个垂直层模拟结果比较吻合，进一步验证了采用35～40个垂直层进行模拟运算是比较合适的。

图7 香港地区垂直风廓线观测资料

4 在垂直层加密的基础上分析土地利用更新对垂直剖面流场的影响

在MM5全球25类土地利用数据的基础上，根据2001年的遥感卫星土地利用数据和香港土地利用分布图，对模式的土地利用数据进行了更新，图8是珠三角地区土地利用更新情况。

4.1 1993年和2001年土地类型下垂直风速廓线模拟结果分析

根据垂直方向σ层加密调整的讨论结果，综合计算量和稳定性等方面的考虑，在垂直方向上采用35个σ层的模拟方案，分别采用新旧土地利用数据对2001年9月15-18日进行了数值模拟。

在图9中，OLU代表采用旧的1993年土地利用的模拟结果，NLU代表采用新的2001年土地利用的模拟结果(下同)。可以看到，土地类型改变对垂直方向风速的影响在本次模拟时段主要为2 500 m以下的中低空，再往上则两种情况下的模拟风速逐渐趋向一致。就时间上而言，2001年9月16日下午14∶00的影响程度比凌晨2∶00大。由于新土地类型的城市范围大大增加，土地利用类型的差异导致热容量、反照率等不一样，从而产生地面层热量的变化，由热量不均而产生的压力梯度对风场的影响比较明显，城市范围扩大造成的白天午后城市与郊区温差增大会使得这种压力梯度加大。在中尺度模拟的情况下，城市建筑物高度及密集度的影响被弱化，相对于城市化带来的地表粗糙度增加对风速的影响而言，温度差及其带来的压力差对较大范围内风速的影响是占主导地位的，所以16日下午14∶00，在近地面层新土地类型比旧土地类型模拟的风速明显增大。凌晨2∶00这种热量不均所产生的影响作用有较大的下降，因而两种情况下模拟的风速比较接近，这个时候gz1点的风场主要是受到大背景风场的影响，吹偏北风。

图8 模式采用的土地利用示意图

图9 gz1点模拟垂直风剖面

4.2 1993年和2001年两种土地利用类型情况下珠江口附近垂直剖面的流场及分析

在海陆分布差异及城市热岛效应造成热量分布不均的情况下，珠江口两岸及水面之间形成了一种比较独特的海陆风结构。下面在垂直层加密的基础上，分析新旧两种土地利用情况下模拟的珠江口区域海陆风垂直剖面流场的特点。

1)选取的珠江口附近垂直剖面示意图

图10 珠江口垂直剖面位置示意图

从图10可以看到，选取的三个垂直剖面分别为AD，DB和DC，这3个垂直剖面交汇于位于珠江口中心的轴线D。A，B，C和D的精确经纬度分别为：113.4°E22.35°N，114.1°E22.7°N，114.15°E22.4°N和113.75°E22.4°N。

图11 16日14:00 DB垂直剖面流场

图12 16日14:00 AD垂直剖面流场

图13 16日14:00 DC垂直剖面流场

2)垂直剖面流场结果及分析

图11-13分别输出16日14:00新旧两种土地类型在DB、AD、DC三个垂直剖面的流场，各剖面垂直方向采用σ坐标，垂直方向高度范围大约为地面至3 000 m高空，流线上不同的颜色代表风速大小。通过分析我们可以发现气流在这三个垂直剖面的流动规律。在DB截面上可以看到一个大的逆时针方向的涡，底层的气流被卷向高空，然后向西南方向输运，AD剖面的右边也可以看到一个大的逆时针方向的涡，这样东方高空的来流接着被卷到底层，然后向东输运，DC剖面的左边也同样可以看到一个大的逆时针方向的涡，在这种情况下，从西方输运来的底层气流就被继续向东输运，抵达香港，然后气流又被卷向高空，从而形成一个大的循环。另外在AD剖面的左边可以看到有一个顺时针方向的涡旋，这表明从东面来的高空气流在这个区域先下沉，在底层则向两边扩散，而AD剖面最左边(AD剖面在陆地上的部分)的气流则是上升气流，所以，AD剖面顺时针涡的右侧下沉区域就是DB、AD和DC剖面高空气流的辐合下沉及底层气流的辐散区。

通过分析3个剖面的位置，知道DB剖面右边的底层气流来自于经济和工业高速发展的珠三角深圳、东莞等区域，AD剖面左边的底层气流同样来自于工业发展迅速的珠海，中山等地区，相对而言，这些气流会容纳较多的工农业和生活污染物，这些污染物有一部分会通过这种循环输运过程被带到香港，使得空气指数升高。

新的土地利用能更真实地反映珠江口附近城市密集区域与水面之间的海陆分布差异，因而在新旧两种土地类型的模拟垂直剖面流场对比中，我们看到，相对于旧土地类型，新土地类型模拟结果显示珠江口区域这种大气循环输运作用得到了加强，无论是DB、AD还是DC截面，新土地类型模拟的风速大小都比旧类型有较大的增加，AD剖面左边的下沉区也得到了明显地加强，如果这种加强得到继续，则可能通过循环输运过程把更多的珠江口两岸工业污染物带到香港地区，从而加重当地的空气污染程度，例如，位于香港西部靠近珠江口的东涌站测到的空气污染指数相比其他站点存在偏高的情况，这种输运过程可能是造成这种情况的其中一个因素。

5 结 论

1) 通过不同的垂直σ层的数值模拟结果的对比分析和讨论，我们可以知道，通过调整和增加垂直方向上的σ层可以提高近地面和中高空气象场的模拟准确度，而且可以更精细地反映出垂直风剖面随高度增加的变化，使得模拟的结果更接近实际的气象条件，综合而言，在模拟的过程中，垂直方向采用35～40层是比较合适的。

2) 为了控制适当的计算量和保证计算的稳定性，需要合理安排垂直层的间距和总的垂直层数，对于研究近地面的气象场而言可以较多地进行近地面地垂直层加密，对于研究高空气象场比如说高空急流而言，可以增加高空的垂直层以保证得到更好的模拟结果，在本次选取的模拟工况下，2 000 m以下采用20～25层，2 000 m以上采用15～25层已经可以得到较好的垂直风剖面模拟结果。

3) 新的土地利用类型更真实地反映了珠江口附近城市密集区域与水面之间的海陆分布差异，从而更准确模拟该区域海陆风垂直剖面的流场结构，更清晰地反映珠江口两岸地区在这种海陆风模式下的大气输运过程，从而提高对该区域大气模拟的准确度并对分析空气污染成因起到积极的作用。

参考文献：

[1] GRELL G A , DUDHIA J , STAUFFER D R. A description of the fifth-generation Penn State/NCAR mesoscale model(MM5)[R/OL]. NCAR/TN-398+STR. NCAR Technical Note, 1995.http://www.mmm.ucar.edu/mm5/documents/mm5-desc-doc.html.

[2] 张金善,钟中,黄瑾. 中尺度大气模式MM5简介[J]. 海洋预报，2005, 22(1):31-40.

[3] 陈静,薛纪善,颜红.物理过程参数化方案对中尺度暴雨数值模拟影响的研究[J].气象学报,2003,61(2):203-217.

[4] ANTHES R A. A cumulus parameterization scheme utilizing a one-dimensional cloud model[J]. Mon Wea Rev, 1977, 105:270-286.

[5] 王建捷,周斌,郭肖容.不同对流参数化方案试验中凝结加热的特征及对暴雨中尺度模拟结果的影响[J].气象学报,2005,63(4):405-417.

[6] GRELL G A. Prognostic evolution of assumption used by cumulus parameterizations[J]. Mon Wea Rev, 1993, 121: 764-787.

[7] ARAKAWA A, SCHUBERT W H. Interaction of a cumulus cloud ensemble with the large-scale environment, Part I [J].J Atmos Sci, 1974, 31: 674-701.

[8] FRITSCH J M , CHAPPELL C F. Numerical prediction of convectively driven mesoscale pressure system, Part I: convective parameterization[J].J Atmos Sci, 1980, 37: 1722-1733.

[9] KAIN J S , FRITSCH J M. A one-dimensional entraining/detraining plume model and its application in convective parameterization [J].J Atmos Sci, 1990, 47(23):2784-2802.

[10] HONG S Y, PAN H L. Nonlocal boundary layer vertical diffusion in a medium-range forecast model[J]. Mon Wea Rev, 1996, 124:2322-2339.

[11] 江勇,赵鸣,汤剑平,苏炳凯. MM5中新边界层方案的引入和对比试验[J]. 气象科学, 2002, 22(3):253-263.

[12] 袁信轩. 我国江南西南风低空急流形成的天气学分析[J]. 气象学报,1981,39(2):245-251.

[13] 况雪源，张耀存. 东亚副热带西风急流季节变化特征及其热力影响机制探讨[J]. 气象学报，2006，64(5):564-575.