作为气溶胶主要组分之一，风蚀粉尘在全球气候变化中扮演着重要的角色：它的辐射强迫效应影响大气热量的收支^［1］；作为凝结核它影响大气的水文过程如液滴和冰晶的形成并进一步影响降水过程和大气热量收支状况^［2］；作为海洋生态系统的主要营养来源，它的沉降提高了海洋生态系统的生产力并成为影响全球气候变化的重要因素^［3］。风蚀粉尘释放、输送和沉降产生的环境问题也影响到人类的可持续发展。

粉尘释放模拟是粉尘气候效应评价和沙尘天气预报的基础。作为一种常用的粉尘释放模式，DPM模式^［4］已经被用于模拟中国北方及亚洲地区沙尘天气^［5-7］，并取得了较好的预报效果，但该模式在预报蒙古高原的戈壁及戈壁草原区域沙尘天气时^［7］，出现了漏报沙尘天气的现象。针对此类问题，前人通过优化地表的参数化方案如更新区域沙漠分布面积等数据来减小漏报的概率^［7］。但应该看到，沙尘天气预报的效果除与地表参数化方案及可靠性有关外，还与DPM模式的精度有关，其中起动摩阻风速模式是影响模拟精度的重要方面。

DPM模式推荐了两个用来模拟粗糙床面的起动摩阻风速Ut*的模式，分别是MB模式^［8-9］和SL模式^［10］。条件试验表明，前者对粗糙床面空气动力学粗糙度Z0非常敏感，仅适合模拟Z₀≤1 mm的粗糙床面的Ut*；后者预测的Ut*会随着Z0的增加而缓慢地增加，似乎可模拟更大的Z0情形下的Ut*^［11］。二者的模拟结果存在明显的差异，但由于缺乏实测数据的支持，目前还难以评价它们的可靠性。

在边界层湍流和粗糙床面相互作用的过程中，湍流总动量将在粗糙元和可蚀性地表之间分配，这被称为湍流切应力分解机理^［12］。目前模拟湍流切应力分解的方案可分为Raupach模式^［13-14］、MB模式^［8-9］和SL模式^［10］等。Raupach模式把起动摩阻风速比率（eff=Uts*Ut*， Uts*和Ut*分别代表光滑和粗糙床面的起动摩阻风速）表征为粗糙元密度、方向比率AR（粗糙元高度与直径的比值）、切应力空间分布系数和阻力系数的函数^［13-14］， 其被相关的风洞实验和野外实验不断地验证和修正^［15-19］。MB模式^［8-9］和SL模式^［10］基于粗糙床面的双边界层假定：粗糙元高度以上为代表整个床面气流特征的外边界层而粗糙元高度以下为代表粗糙元之间光滑床面气流特征的内边界层，且内外边界层风廓线皆服从对数廓线。双边界层成立需要满足低的粗糙元密度（常以小于0.1为界^［8-9］）和较大的粗糙元间距（常以间距与高度的比值大于3^［10］），而前人研究数据则暗示风蚀床面的粗糙元密度及Z0常常会超过上述条件^{［11，20-25］}，这意味着MB模式和SL模式可能不适合预测更高的Z0情形下的Ut*，但缺乏直接观测数据的支持。

从Raupach模式相关研究来看，粗糙元的方向比率AR和孔隙度P等因素对Ut*有显著的影响^［15］。对MB模式和SL模式而言，是否要把AR和P的因素整合到模式中还需要进一步的研究。

为了评价MB模式和SL模式的可靠性，拓展模式的适用范围，进一步理解粗糙元的方向比率和孔隙度对起动摩阻风速的影响，本文利用风洞实验来研究粗糙和光滑床面的空气动力学粗糙度、方向比率和孔隙度等因素对起动摩阻风速的影响，建立基于粗糙元结构参数和空气动力学粗糙度的新模式。

风洞实验在中国科学院寒区旱区环境工程研究所沙漠与沙漠化重点实验室的野外风洞中完成，该风洞几何尺寸和空气动力学性质见文献［21］。

布置了39组粗糙床面：25组细高、密实的粗糙元覆盖的床面，对应的方向比率分别为4、8、12、16、20，5组低矮的密实粗糙元，对应方向比率分别为0.25、0.5、0.75、1.0、1.25，9组孔隙粗糙元覆盖的床面（也是丛状分布的粗糙元），对应的孔隙度分别为0.15、0.55、0.75。具体的粗糙元结构参数定义及数值见文献［26］（未测量粗糙元高度为40 mm，个数密度为1 025床面的Ut*，实际共有38组床面的Ut*数据），孔隙度的定义为公式（1），ns，hs，ds是围成一个孔隙粗糙元的单个粗糙元的个数、高度、直径，在床面设计中令hs=H（图1）。

起动摩阻风速测试所用的沙样来自腾格里沙漠的沙漠沙，粒径100~450 μm，粒径分布符合双峰对数正态分布模式，可用两个中值粒径分别为170 μm和320 μm对数正态分布群体构成，后者占比为84%（表1）。

粗糙床面跃移起动的观测方法是在粗糙床面上均匀撒播约1 mm厚的沙粒，采用摄像头监测下风向沙粒起动（图2）。缓慢增加风机频率，当观测到胶带粘贴到5~10粒沙粒时，用风速廓线仪记录风速廓线数据。为了减小实验的误差，每种床面做5次平行实验，取平均值作为起动摩阻风速。这种判断沙粒起动的方法已被证明是有效的^［15］。裸露沙面的跃移起动观测方法与粗糙床面相同。

用风速廓线仪测定沙粒起动时的风廓线。风廓线仪由高度分别为3.0、6.0、10.0、15.0、30.0、60.0、120、200、350、500 mm的10个精细毕托管组成，使用前用标准毕托管标定。风速采集时长60 s，采样频率为1 Hz，用于分析的是60 s内采集的平均风速。

光滑床面的空气动力学粗糙度Z0s用中性层结下的对数风廓线方程（公式2）拟合计算，其中k为卡门常数，取0.4，U*为摩阻风速，Z0s为空气动力学粗糙度，UZ为高度为Z处的平均风速。利用风廓线仪分别测定自由风速为6、8、10、12、14、16、18、20 m·s^-1情况下的风廓线，利用公式（2）拟合可得到Z0s并取平均值。同样的，粗糙床面的空气动力学粗糙度Z0也是对公式（2）拟合并取平均值。起动摩阻风速Ut*也是根据对数风廓线公式（2）拟合并取平均值^［21］：

评估MB模式^［8-9］和SL模式^［10］的相对误差率E_r，采用公式（3）来表述，其中Utm*为模拟的起动摩阻风速，Ut0*为风洞实验观测的起动摩阻风速。

当Z0从0.02 mm增加到26.5 mm时，Ut*约从0.30 m·s^-1以对数函数形式增加到1.45 m·s^-1（图3），但Ut*明显的分散性意味着床面的某些结构参数可能对其有着重要的影响。鉴于粗糙元密度、高度对气流的阻滞效应可以通过Z0来表征^［21］，而粗糙元间距与起动摩阻风速成负相关关系，则意味着方向比率可能是造成摩阻风速数据分散性的主要原因。

在相近的或相同的AR情况下，Ut*与Z0具有更高的相关系数，则意味着方向比率对起动摩阻风速有重要的影响（图4）。当AR从0.25增加到8时，Ut*均随着Z0以对数函数Ut*=alnZ0+b增加，且经验系数a、b因AR的不同而有轻微的变化（图4A）。

当AR从12增加到25时，Ut*随着Z0以对数函数Ut*=alnZ0+b增加，且对应的经验系数a、b则随着方向比率具有明显的线性增加趋势（图4B）。

本文用丛状分布的粗糙元来模拟孔隙粗糙元对Ut*的影响。当孔隙粗糙元对应的Z0与密实粗糙元相近，且在AR相近的情况下，丛状分布的孔隙粗糙元（AR=15，20，25）与密实粗糙元（AR=16，20）服从同一对数函数（图4B），这表明孔隙率对跃移起动的影响可归结为Z0的影响。

MB模式预测的Ut*分别呈现出与风洞实验结果相对一致、大于实验结果0.2~1.2 m·s^-1及预测失效出现零值等3种情形（图4），分别对应Z0≤1mm、1≤Z0≤5mm、5<Z0<30mm，这表明MB模式仅适合预测Z0≤1mm的风蚀床面的起动摩阻风速。前人对蒙古高原戈壁草原地区沙尘天气的漏报^［7］，与所用的MB模式对Ut*的过高估计有关。

当Z0从0.01 mm增加到26.5 mm时，SL模式预测的Ut*从0.30 m·s^-1增加到0.68 m·s^-1（图4），明显低于风洞实验观测结果，特别是在Z₀>1 mm的情形下，这种低估效应更为显著。

图5表明，当Z₀≤2 mm时，MB模式的相对误差率-20%~40%；当Z₀>2 mm时，MB模式的相对误差率80%~100%或出现零值（图5中出现-100的情形），这意味MB模式在此情形下失效。SL模式的相对误差率-50%~-20%，相对误差率对Z0的变化不敏感，这表明该模式在不同Z0情形下均会显著地低估起动摩阻风速。

MB模式对Ut*的高估、失效及SL模式对Ut*的低估与双边界层的假定有关，而双边界层的成立需要满足粗糙元密度不超过0.1^［8-9］、粗糙元的流向间距与粗糙元高度比值在3以上的床面等条件^［10］，本文约有34%的数据不满足此条件^［26］，因而导致MB模式高估、失效及SL模式的低估或失效。

另外，SL模式对Ut*的低估也与它所依据的风洞实验数据有关^［10］。该风洞实验的Z0约为0.04~1.37 mm，对应的阻力分解系数集中于0.5~1.0，缺乏更大Z0情形下Ut*的数据，因而当Z0>2mm时，该模式会出现低估Ut*的情形。

为了说明粗糙床面湍流切应力的分配趋势，这里引入起动摩阻风速比率 e_ff，它为裸露沙面的起动摩阻风速Uts*和粗糙床面的起动摩阻风速Ut*的比值，其中实验测量的Uts*为0.25 m·s^-1，裸露沙面的空气动力学粗糙度Z0s为0.05 mm。

图6表明，在AR（0.25~25）相近或相同的情形下，e_{ff均随着lnZ0Z0s增加而减小，且eff数据与拟合函数的差值均在±0.1 m·s^-1以内，这意味着把方向比率整合于起动摩阻风速模式之中，就可以使模拟的误差在±0.1 m·s^-1以内。}

另外，在粗糙元丛状分布的情形下（即对应着孔隙率0.15，0.55，0.75的粗糙床面），这类粗糙床面的e_ff数值与相近AR情形下密实粗糙元接近（图6），这表明丛状分布或孔隙度对Ut*的影响可用Z0表征。

e_ff可以表示为公式（4），其中a'和b'为经验参数：

图7表明a'和b'分别是AR的一次函数和二次函数。公式（4）综合地反映了方向比率、Z0与Z0s的比值对e_ff的影响。由公式（3）可得到Ut*的函数，它是AR、Z0与Z0s的比值以及裸露沙面的Uts*的函数（图6）。

与SL模式相比，公式（4）的e_ff所涵盖的范围为0.15~0.85，其可以更好地模拟低、中、高粗糙元密度情况下的起动摩阻风速，在模拟沙化的戈壁草原地区的跃移起动和粉尘释放强度有更好的潜力。

需要说明的是，为了反映沙化风蚀的草原中不同群落类型的空气动力学粗糙度，梅凡民^［27］曾采用数字化群落类型图中建群种平均的高度、覆盖度和单个建群种的方向比率等参数来计算空气动力学粗糙度。同样地，可采用公式（4）来预测此群落类型的起动摩阻风速，此工作将会在另外的论文中讨论。

风洞实验表明，在中、高密度情形下，MB模式会高估（相对误差率80%~100%）起动摩阻风速而失效，SL模式会低估起动摩阻风速（-50%~-20%）；

粗糙元的方向比率对起动摩阻风速有显著的影响，而粗糙元的孔隙率对起动摩阻风速的影响可归结为其对应的空气动力学粗糙度效应，故起动摩阻风速的模式可不考虑孔隙度的影响。

本文建立的基于粗糙元的方向比率、粗糙床面与光滑床面的空气动力学粗糙度比值及裸沙起动摩阻风速的经验模式可更好地预测低、中、高密度的粗糙元覆盖条件下的起动摩阻风速，误差在±0.1 m·s^-1以内。