1 引言
强对流天气一般与雷暴大风、 短时强降水、 冰雹等相伴, 具有突发性强、 破坏力大等特点, 常常导致重大人员伤亡和财产损失（郑永光等, 2010）。目前, 强对流天气的研究主要采用两种手段: 一是运用大气再分析资料、 卫星和雷达等观测, 诊断分析其发生、 发展和演化（陈英英等, 2013; 朱平等, 2019）; 二是采用数值模拟方法, 探究初始场、 物理参数化方案、 下垫面以及环境背景等对模拟的影响（李娜等, 2013; 曹倩等, 2016; 况祥等, 2018; 张弛等, 2019; 陶局等, 2019; 周林帆等, 2019）。强对流天气与云微物理过程高度关联, 但云微物理的内部结构和演化高度复杂且尺度微小, 直接对其预报困难很大, 只有对微物理过程进行参数化处理。由于对关键物理过程理解还不全面和深刻, 加之问题高度复杂性, 采用简单参数化方案很难刻画问题的本质。目前, 无论是天气模式还是气候模式, 云微物理过程参数化方案还存在很大不确定性（Min et al, 2015）。
云微物理参数化方案(下称云微物理方案)一般可以分为两大类: 分档（Bin）方案和函数（Bulk）方案, 其中Bin方案根据水成物的微物理特征（如相态、 粒子大小、 形状、 密度等）对粒子分档, 建立每一类粒子演化方程进行预报, 而Bulk方案采用预先指定的经验函数（如Marshall-Palmer和Gamma函数）来描述水成物的总体分布特征, 预报量变为对应函数的一阶矩或多阶距。常用Bulk方案有单参和双参, 其中单参方案仅预报水成物混合比（Hong et al, 2004）, 双参方案同时预报水成物混合比和水成物数浓度（Morrison et al, 2009）。针对云微物理方案对强对流天气模拟影响, 国内外已有较多研究。大部分云微物理方案能较好地模拟出大尺度特征, 但模拟的对流系统形态、 演变及内部动力、 热力结构还存在差别（Luo et al, 2010; Yin et al, 2017; 张弛等, 2019）。例如, 针对一个中纬度飑线模拟, 用Bulk方案表现为多单体风暴, 发展迅速, 层云区狭窄, 而Bin方案呈现为单体风暴结构, 发展缓慢, 层云区宽广（Li et al, 2009a, 2009b）; 飑线的观测系统模拟试验显示, 双参方案的模拟效果更接近实况, 而单参方案的雨水蒸发过快, 后部层云降水模拟较差（Bryan et al, 2012）。
尽管已有许多学者探究不同云微物理方案对飑线模拟的影响, 但研究集中于强天气背景强迫。在弱天气尺度强迫下也会发生强对流, 但预报难度大, 是预报业务的难点之一（肖现等, 2015）。此外, 常用云微物理参数化方案很多, 复杂程度和不确定因素差别较大（如冰相的类型和数目）, 并且关键参数值（如粒子谱截距、 粒子间相互转换的阈值、 下落末速等）主要基于国外观测数据得出, 是否适合江淮地区的强对流预报也是值得研究的。因此, 本文将采用WRF模式中5种具有代表性云微物理参数化方案对2016年7月26日发生在江苏北部的一次弱天气尺度强迫下的强对流进行高分辨率模拟, 探究不同云微物理方案模拟差异的原因。
2 资料来源及方法介绍
2.1 资料来来源及试验设计
采用6 h一次1°×1°的FNL分析资料, 分析此次强对流过程的大尺度背景。同时利用江苏射阳站2016年7月26日08:00（北京时, 下同）和14:00探空资料, 以及位于盐城（33.42°N, 120.2°E, 海拔28.1 m）的新一代多普勒天气雷达资料, 对此次强对流过程演变进行分析。文中涉及的地图是基于国家测绘地理信息局标准地图服务网站下载的审图号为GS（2016）2923号的标准地图制作, 底图无修改。
模拟采用WRF3.8.1模式, 双向反馈两重嵌套, 中心点为（33.3°N, 117.8°E）, 外层格点数为190×154, 格距为9 km, 内层格点数为316×241, 格距为3 km, 垂直方向分为30层, 模式顶气压为50 hPa。积分时间为7月26日08:00至27日02:00, 积分时间步长为30 s, 每6 min输出一次结果。由6 h一次1°×1°的NCEP FNL分析资料提供初始条件和边界条件。主要物理参数化方案为: RRTM长波辐射方案, Dudhia短波辐射方案, Monin-Obukhov近地层方案, Noah陆面模式, YSU边界层方案, 在第一重嵌套中采用Kain-Fritsch积云对流方案, 第二重嵌套中关闭积云对流方案。
在其他物理过程参数化方案和初始场设置相同情况下, 通过改变云微物理方案进行敏感性试验, 5种方案分别为NSSL 1-momlfo（简称NSSL）、 Morrison、 Thompson、 WDM6和快速HUJI 方案, 其中前四种方案属于Bulk方案, 最后一种方案属于Bin方案, 具体细节见表1。
表1 五种云微物理方案的预报量Table 1 Forecast variables of five cloud microphysical schemes |
| 云微物理方案 | 质量变量/（kg·kg-1） | 数浓度变量/kg-1 |
|---|---|---|
| NSSL 1-momlfo | q c、 q r、 q i、 q s、 q g | 无 |
| Morrison | q c、 q r、 q i、 q s、 q g | N r、 N i、 N s、 N g |
| Thompson | q c、 q r、 q i、 q s、 q g | N i、 N r |
| WDM6 | q c、 q r、 q i、 q s、 q g | N n、 N c、 N r |
| HUJI fast | q c、 q r、 q i、 q s、 q g | N c、 N r、 N s、 N i、 N g、 N n |
q为混合比, N为数浓度, 下标c、 r、 i、 s、 g、 n分别代表云水、 雨水、 云冰、 雪、 霰和云凝结核 |
2.2 云微物理参数化方案和模拟雷达回波算子
NSSL方案基于3类冰相的单参Bulk方案（Lin et al, 1983）发展以来， 考虑了雨水蒸发、 云雨自动转换、 冰相粒子的升华、 凝华、 冻结、 融化、 相互碰并和冰核活化等34种云微物理过程, 并且采用Tao et al（1989）饱和修正方案, 从而保证在一个网格点上不存在过饱和、 半饱和状态。采用指数函数来描述雨水、 雪粒子、 霰粒子的粒子谱。此外, NSSL方案考虑了云冰和云水粒子的直径和下降末速度, 从而导致一些碰并过程的不同。该方案可以用于研究冰相和液相粒子较多的过程（Gilmore et al, 2004）。
Morrison方案考虑了核化、 自动转换、 收集、 繁生等40种云微物理过程。所有粒子的形状都假设为球形, 云滴和雨滴服从Gamma函数分布, 雨、 雪、 霰及云冰粒子采用指数函数分布。方案包括多种过程开关, 例如可根据情况选择预报霰粒子或雹粒子, 应用较为广泛（Morrison et al, 2009）。
Thompson方案是对Reinser et al（1998）方案的改进, 与其他Bulk方案的主要区别是: 一是假设雪粒子的分布同时依赖于冰水的含量和温度。二是雪粒子的分布由指数函数分布和Gamma函数分布共同决定的。三是假设雪粒子的形状是非球形的, 并且密度会随直径的变化而变化（Thompson et al, 2008）。
WSM6（WRF Single-Moment 6-class）方案大部分参数化过程与基于3类冰相的单参Bulk方案相似, 但在计算积冰及其他一些参数的设置上有差别。采用Dudhia（1989）和Hong et al（1998）的方案分开调整冰和水的饱和过程, 并在粒子下降过程中考虑了凝结、 融化过程, 提高了垂直廓线的精度。WDM6（WRF Double-Moment 6-class）方案是在WSM6方案的基础上发展而来, 采用双参数方案来处理暖雨过程, 并且增加了对云水和雨水的数浓度及云凝结核数的预报, 其他物理过程与WSM6方案类似。因此, 与其他双参数方案不同, WDM6方案可以用于研究气溶胶对云特性及降水过程的影响（Lim et al, 2010）。
快速HUJI方案是一种Bin方案。在保持原始HUJI方案优点前提下, 为节约计算时间, 快速HUJI方案（下称HUJI方案）把冰晶和雪粒子归为低密度冰, 把霰和冰雹归为高密度冰, 这样粒子分布函数由原来8类减小为4类（雨滴、 低密度冰、 高密度冰、 气溶胶）。与Bulk方案相比, HUJI方案中水凝物和气溶胶的尺度分布是在模式积分过程计算出来的, 并且采用随机运动方程来描述云中粒子间的碰撞聚合等过程（Khain et al, 2010）。
模拟的雷达反射率（Z e, 单位: mm6·m-3）是雨水混合比(q r)、 雪粒子混合比(q s)和霰粒子混合比(q g)的函数（Stoelinga, 2005）:
式中: N or、 N os、 N og分别是雨水、 雪粒子、 霰粒子的截距参数（单位: m-4）, 其对应的值分别为8×106、 2×107、 4×106。 、 、 、 分别是干空气、 液态水、 雪粒子、 霰粒子的密度（单位: kg·m-3）, 其中 、 、 的值分别为1000、 100、 400。可知, 雨水的反射能力要强于雪和霰。
3 天气形势
2016年7月26日08:00江淮上空200 hPa无急流（图略）。同一时刻500 hPa环流形势场[图1(a)]显示, 我国东北地区存在一个深厚的冷涡, 安徽和江苏地区受西太平洋副热带高压控制。700 hPa和850 hPa上[图1(b), (c)]均有来自北方的干冷气流和来自西南的暖湿气流在安徽和江苏的北部辐合, 但是冷空气来源不同, 其中700 hPa上来源于西西伯利亚和我国东北, 而850 hPa上冷空气主要来自我国东北。在辐合区出现了一个风速相对大值区, 但没有低空急流。在700 hPa气流辐合区的相对湿度为50%~70%, 而在850 hPa江苏和安徽北部的相对湿度达到了80%; 西南暖湿气流为江淮地区输送了充足的水汽, 且低层的水汽较中层更为充沛。14:00 500 hPa环流形势场（图略）无明显变化, 东北冷涡稳定少动, 安徽和江苏地区仍处于西太平洋副热带高压控制之下。
图1 2016年7月26日08:00 风场（a~c, 矢量, 单位: m·s-1）与500 hPa位势高度场（a, 黑色实线, 单位: dagpm）以及700 hPa（b）和850 hPa（c）相对湿度（阴影, 单位: %）Fig.1 The wind filed (a~c, vector, unit: m·s-1) and the geopotential height at 500 hPa (a, black solid lines, unit: dagpm) and the relative humidity (the shaded, unit: %) at 700 hPa (b) and 850 hPa (c) at 08:00 on 26 July 2016 |
2016年7月26日08:00射阳站的探空资料[图2(a)]显示, 对流有效位能（Convective Available Potential Energy, CAPE）达到2492.4 J·kg-1, 对流抑制能（Convective Inhibition, CIN）为112.2 J·kg-1, 925 hPa以下有很薄逆温层。由温度和露点温度的差值可以看出, 850 hPa以下大气水汽较为充足, 并存在弱暖平流, 而850~700 hPa存在弱冷平流, 形成了上干冷下暖湿的不稳定层结。但是, 风速和风向的垂直切变较弱。到14:00, CAPE值增加至3555.7 J·kg-1, 而CIN值减小至4.2 J·kg-1, 这为对流的发生提供足够大的不稳定能量, 此时自由对流高度降低到925 hPa, 0~3 km风的垂直切变大小为8 m·s-1, 大气处于极不稳定状态。以上分析表明, 此次强对流发生在强的热力不稳定和偏弱的垂直风切变环境下。
4 实况雷达回波分析
图3为6个时次盐城站1.5°仰角雷达反射率。2016年7月26日15:16在宿迁附近出现零散对流单体（图略）, 但强度较弱, 随后缓慢向东南方向移动。15:56[图3(a)], 对流单体在宿迁东北部合并, 形成有较强回波中心的块状回波区（记为A）, 中心强度>50 dBz, 并在宿迁西南部又出现分散对流单体（记为B）。此后, A和B继续向东南方向缓慢移动, 并在移动过程中继续发展加强。17:21 A、 B相连形成类似“人”形的回波区, 中心强度>55 dBz[图3(b)]。17:50回波带到达连云港、 淮安和盐城市境内, 并由人字型逐渐演变为弓形, 中心强度保持在55 dBz以上[图3(c)]。18:13回波带断裂为三部分（分别记为C、 D、 E）。在回波带继续向东南移动过程中, D开始减弱[图3(d)]。19:15 E移入黄海, D再次加强[图3(e)]。21:09线状回波消失[图3(f)], 此次强对流过程基本结束, 整个过程维持了约6 h。
5 结果与分析
5.1 雷达回波
图4、 图5和图6分别为2016年7月26日16:00、 16:30和17:00不同云微物理方案模拟的雷达组合反射率与观测值的对比。5种方案中, NSSL方案模拟的虚假强回波较少, 与实况较接近但触发位置靠前, 而其他4种方案都出现不同程度的虚假强回波。17:00 NSSL方案模拟出四个强对流单体但未形成线状对流带。Morrison方案模拟回波结构较松散, 强回波区与观测偏离较大, 如强对流中心错误地出现在海上, 并且在后续发展过程中, 位于江苏省境内的强对流单体迅速消亡。Thompson方案出现大范围虚假层云区, 且多个大于45 dBz的虚假强对流单体分布在层云区, 但随时间发展减弱, 它们未能在17:00形成强回波带。在整个强对流演化阶段, WDM6方案无论是虚假对流的强度还是范围在所有方案中都是最大, 尤其回波带宽度明显大于观测, 在海面以及安徽和江苏交接区域出现明显虚假强回波。同样, HUJI方案模拟的虚假回波也很多且触发位置靠前, 但强回波带的宽度和强度明显小于WDM6方案。
图4 2016年7月26日16:00不同云微物理方案模拟的雷达组合反射率（单位: dBz）Fig.4 The radar composite reflectivity with different cloud microphysics schemes at 16:00 on 26 July 2016.Unit: dBz |
图5 2016年7月26日16:30不同云微物理方案模拟的雷达组合反射率（单位: dBz）（a）中黑色线段AB为后文中提到剖面位置 Fig.5 The radar composite reflectivity with different cloud microphysics schemes at 16:30 on 26 July 2016.Unit: dBz.The black line AB in Fig.5(a) is the section position mentioned later |
为了进一步分析不同云微物理方案模拟强对流空间结构的差异, 图7给出了2016年7月26日16:30回波沿图5中线段“AB”垂直剖面图。观测的回波顶高约18 km, 强对流区位于对流前沿, 高达10 km, 并且后部有约10 km的层云区[图7(a)]。NSSL方案较好地反映了强对流区雷达反射率的结构, 再现了强回波对流顶的高度, 例如在对流区前沿模拟出一个发展至6 km的强对流中心、 中心强度达45 dBz以上, 但模拟回波位置整体较实际偏东南, 且后部层状云区水平范围很小、 高度也未能达到10 km。Morrison方案和Thompson方案均没有模拟出大于45 dBz的强对流区, 模拟效果较差。WDM6方案和HUJI方案均在对流前沿模拟出一个强对流区, 但模拟的回波顶高仅有14 km左右, 且强度和位置都与观测有较大差异。在强对流区的后部, WDM6方案和HUJI方案还出现一个远大于观测的层状云区。总之, 不同云微物理方案模拟的雷达反射率在位置和强度上存在明显差异, 相对来说, NSSL方案模拟的回波位置和强度与实际情况较吻合。
5.2 动力和热力特征
为了探究不同云微物理方案模拟的雷达回波间差异内在原因, 下面将从风场和温度场的时空分布两方面进行分析。飑线一般具有准二维结构, 因此在飑线上作垂直于飑线移动方向的剖面图就可以分析飑线的动力和热力特征（Thorpe et al, 1982）。图8为2016年7月26日16:30不同云微物理方案模拟的相对风暴风场（即水平风场减去系统水平移速以及垂直风速放大10倍后合成的风场）与扰动位温场的垂直剖面图。5种方案均在3~6 km高度模拟出强度相当的西北向气流。其中, NSSL方案的西北向气流与前部的东南暖湿气流辐合形成强烈的上升气流, 并一直延伸至18 km以上; 冷池出流与东南暖湿气流又在强对流前沿辐合, 新的对流单体正在触发, 从而保证强对流维持和传播。Morrison方案和Thompson方案中没有出现中低层西北气流与前侧的东南暖湿气流的辐合上升, 因此也没有模拟出此刻的强对流。WDM6方案和HUJI方案的后部入流在到达强对流区后受降水的拖曳作用迅速下沉, 加强了前侧东南暖湿气流的抬升, 使暖湿气流倾斜上升并穿过对流区, 最终在12 km左右流出回波区。在这两种方案中, 下沉的西北气流与冷池出流在地面合并, 并扩散至系统的前侧, 从而阻断了东南暖湿气流向强对流系统的输送。因此, WDM6方案和HUJI方案模拟的强对流在该时刻后开始衰亡。
图8 不同云微物理方案模拟的2016年7月26日16:30相对风暴风场[矢量, (u, v, 10w), 单位: m·s-1]和扰动位温（彩色区, 单位: K）沿图5(a)中线段AB的垂直剖面图Fig.8 The vertical cross section of storm-relative wind field [vector, (u, v, 10w), unit: m·s-1] and potential temperature perturbation (color area, unit: K) along the line AB in Fig.5 (a) with different cloud microphysics schemes at 16:30 on 26 July 2016 |
研究表明, 近地面冷池会与低层风切变相互作用, 在前沿产生垂直上升运动, 进而对飑线的维持和发展产生影响（Rotunno et al, 1988; 陈明轩等, 2012）。5种云微物理方案模拟的扰动位温场在0~18 km高度之间具有较相似的空间结构, 均表现为冷-热-冷的三层结构, 其中偏冷低层的厚度达6 km, 这是对流层中层干冷空气卷入以及降水蒸发冷却共同造成的。NSSL方案冷池范围较小, 其中心扰动位温约为-4 K, 对应上方约12 km处有一正扰动位温大值区, 中心值大于4 K, 这是由于强的水汽凝结潜热释放造成的。Morrison方案和Thompson方案由于没有模拟出强对流, 相应冷池范围很小, 尤其是Morrison方案几乎没有冷池。WDM6方案和HUJI方案模拟的冷池强度和范围均比其他方案大得多, 这与它们模拟的强对流回波带明显偏宽、 强度较强相呼应（图7）。在5种方案中, WDM6方案模拟的地面冷池范围最大、 强度最强, 其中心扰动温度甚至降至-6 K以下, 这也解释了为何WDM6方案模拟的层状云区较其他方案更宽广。综合前述, 不同云微物理方案模拟的风场与位温场等要素场存在差异, 从而造成对应强对流的位置和范围不同。
5.3 水成物及潜热释放的垂直分布
图9 不同云微物理方案模拟的2016年7月26日16:30雨（阴影区）和云水（等值线）（左）、 冰（阴影区）和雪（等值线）（中）以及霰（右）沿图5（a）中线段AB的垂直剖面（单位: g·kg-1）Fig.9 The vertical cross section of rain (shade) and cloud droplets (contour line) (left), ice (shade) and snow (contour line) (middle) and graupel (right) along the line AB in Fig.5 (a) with different cloud microphysics schemes at 16:30 on 26 July 2016.Unit: g·kg-1 |
从图9中可以看到, NSSL方案的雨水混合比主要分布在7 km以下, 大值在3~5 km, 有两个降水中心; 云水主要出现在1~11 km, 大值在4~6.5 km; 云冰在8 km以上, 大值位于12.5~17 km, 在云冰区还存在云水, 说明过冷水的存在; 雪主要存在6 km以上, 大值区为10~15 km, 其值比云冰小一个量级; 霰主要分布在3 km以上, 大值区在5.5~10.5 km, 三种固态水成物中霰含量最高。雨水、 雪、 霰的垂直分布与雷达反射率的特征一致, 其中低层雨水大值区与强雷达回波中心区相对应, 中高层分布的雪和霰与增高的回波顶相呼应。
Morrison方案中雨水、 云水、 云冰和霰的含量及分布高度均较NSSL方案小; 雪主要分布在5~16 km, 最大值为0.54 g·kg-1, 高于NSSL方案; 霰高层占优。Thompson方案中雨水、 云水、 云冰和霰的含量在五种方案中最低, 其中冰和霰几乎没有分布, 但雪含量较高, 最高值为0.79 g·kg-1, 分布在5~17 km。WDM6方案中水成物空间分布范围最广, 并出现多个极值中心, 这是WDM6方案模拟的强雷达回波范围比NSSL方案大的主要原因; 大范围分布的雨水产生较大范围的降水, 蒸发效应又导致该方案冷池范围比NSSL方案更大。HUJI方案雨水和云水的含量与NSSL方案相当, 但只有一个降水中心; 冰主要分布在4.5~16 km, 较NSSL方案范围宽但含量偏小; 雪与WDM6方案的分布比较接近, 无论大小还是空间范围都偏大; 霰集中分布在6~8 km, 强度和最大高度都小于NSSL方案。
对于几百千米级的中尺度持久对流活动, 对流活动的降水产生的相变潜热加热大气, 进一步促进对流活动（许焕斌, 2015）, 图10给出2016年7月26日16:30各方案模拟的相变潜热变化率的垂直剖面图。NSSL方案中相变潜热变化率正值区主要位于1~15 km, 且空间分布范围最广, 有多个极值中心, 其中在6~9 km达到了0.256 K·s-1, 这是因为NSSL方案在3 km以上模拟出大量冰, 而粒子间转换会产生相变潜热, 加热作用又会增大上升运动, 因此该方案模拟的上升气流更强劲。Morrison方案相变潜热变化率正值区主要位于5~14 km, 最大值约为0.004 K·s-1。Thompson方案相变潜热变化率的正值区零星分布于2~10 km, 最大值约为0.016 K·s-1。WDM6方案相变潜热变化率正值区主要位于6~12 km, 并出现多个极值中心, 最大值约为0.032 K·s-1。HUJI方案相变潜热变化率正值区主要位于2~12.5 km, 分布与NSSL方案相似, 但空间范围较NSSL方案要小, 最大值约为0.032 K·s-1。总之, 不同方案模拟的水成物的差异直接影响潜热的时空分布, 进一步影响不同方案的热动力过程。
6 结论和讨论
采用WRF3.8.1模式中5种不同云微物理参数化方案, 对2016年7月26日发生在江苏北部的一次弱天气尺度强迫下强对流进行模拟研究, 并通过对比分析热动力结构和水成物空间分布来探讨不同方案模拟效果差异的原因, 主要结论如下:
（1） 不同方案模拟的回波强度和范围差异很大, 只有NSSL方案能较好地再现强对流回波区的位置、 强度及层状云区, 而其余Morrison、 Thompson、 WDM6、 HUJI方案模拟回波与实况存在不同程度的差异。
（2） NSSL、 WDM6、 HUJI方案虽然能模拟出低层气流辐合, 但NSSL方案上升气流更强劲并与强雷达回波中心相对应; NSSL方案冷池范围集中而其余方案冷池范围偏大或偏小。水成物空间分布可以很好地解释不同方案模拟的雷达回波间差异: NSSL方案水成物聚集在较狭窄的水平区域, 中低层雨水含量较多, 中高层有较多的雪粒子和霰粒子; 与此对应, NSSL方案相变潜热释放最多, 上升运动也最强。
值得注意的是, 上述结论是在设定的试验方案并且5类云微物理方案中参数没有优化、 均采用默认值时得出的, 用它评价不同云物理参数化方案的优劣可能不客观, 因为即使不更改物理过程参数化方案内核, 优化参数取值也可能改善其模拟效果。此外, 初边界条件、 网格空间分辨率等因素也会影响模拟结果, 这些都需要进一步探究。