作为地球系统中的重要组成部分，沙尘气溶胶通过吸收和散射大气中的长波和短波辐射，改变地气系统能量收支，还通过改变云的微物理过程进而改变云的光学特性、生命期和云量等，从而与气候产生复杂的相互作用；此外，沉降至海洋的沙尘还带来了丰富的矿物质而影响海洋生态系统，从而影响生物地球化学循环过程(Shao et al., 2011)。另外，大气中的沙尘会带来严重的环境污染，直接影响人类的身体健康以及生存居住环境，形成严重的环境问题(Bell et al., 2008)。

数值模式是模拟和预报沙尘过程、研究沙尘 循环及其影响的有力工具。自从20世纪80年代末开始，许多人就致力于研究和发展可用于模拟全球、区域和局地沙尘过程的模式。从20世纪80年代末开始，已经发展一系列的模式，包括全球沙尘输送模式(Huneeus et al., 2011)，以及模拟区域沙尘过程的模式系统如IAPS(孙建华等，2003；雷航等，2005；Lin et al., 2012)、CUACE/Dust(Zhou et al., 2008)、CEMSYS(Shao et al., 2003)等。总体来讲，这些模式大都能模拟出全球或区域沙尘循环过程的基本特征，但模式模拟结果仍存在较大的不确定性(Uno et al., 2006; Huneeus et al., 2011)。针对东亚区域沙尘天气开展的沙尘模式比较计划的结果表明，不同模式系统模拟的沙尘起沙量存在显著的差异，最大和最小的起沙量相差可达12倍(Uno et al., 2006)；而基于气溶胶模式比较计划(AeroCom)的多模式比较结果表明，不同模式对全球沙尘收支的模拟存在显著的差异，所有模式模拟的沙尘总起沙量变化范围为514～4313 Tg a-1，而不同模式模拟的大气沙尘总量的变化范围为6.8～29.5 Tg之间(Huneeus et al., 2011)。

模拟结果不确定性的关键原因之一就在于起沙过程的模拟，这主要由于沙尘起沙过程不仅依赖于近地面的大气状况，还受到植被覆盖、土壤特性等陆面过程参数的影响，因而需要起沙参数化方案能较好地考虑这些要素(Uno et al., 2006；Zhao et al., 2006；Kang et al., 2011)。目前，国内外已经发展了一系列的沙尘起沙方案，根据不同的假设和简化方式，Shao and Dong(2006)将其分为3类：经验性方案、基于简化起沙物理过程的方案以及基于详细起沙微物理过程的方案。由于起沙物理过程的简化方式以及所需的输入参数均有不同，因而不同起沙方案的模拟结果有所差别。Zhao et al.(2006)、Kang et al.(2011)等研究发现，不同起沙方案对沙尘过程模拟的影响显著，包括起沙量、沙尘浓度的空间分布及量值均存在显著差异。尽管上述研究取得了一些进展，但起沙方案对沙尘过程模拟的影响仍需要更多的模拟研究，并采用观测资料进行验证，以利于进一步改进模拟性能。

本文利用耦合了GOCART和Shao04两种起沙参数化方案的WRF/Chem(3.2版本)模式系统，对2002年3月19～22日东亚地区一次典型的强沙尘暴过程进行模拟，重点关注两种不同起沙方案 下模式对起沙过程的模拟及其差异。2002年3月19～22日的强沙尘暴过程是近年来强度强、影响范围广的典型强沙尘暴过程之一，主要影响了我国长江以北的广大地区，其中以我国西北、内蒙古和华北北部的沙尘天气强度最强。此次沙尘暴过程的前期，在蒙古国西南部、南部和东南部的大部分地区也观测到沙尘暴过程；在其后期，伴随着天气系统的移动，沙尘还传输到下游的韩国和日本地区。针对此次沙尘暴过程的发生发展、沙尘的分布及其影响，很多学者利用不同的沙尘天气数值模拟和预测系统，进行了一系列的模拟与预报研究(Shao et al., 2003；孙建华等，2003；Uno et al., 2006；Zhao et al., 2006)。值得指出的是，以往对沙尘天气过程的模拟及比较研究，大都采用不同的气象模式，并耦合了不同的起沙参数化方案进行数值模拟，难以区分不同气象场和沙尘起沙方案影响的相对贡献。本文研究的特点则在于采用了同一个中尺度气象模式框架，但分别耦合了不同的沙尘起沙方案，从而保证了在沙尘天气事件的模拟过程中，模式所用的气象驱动和环流背景是相同的，进而就可以分析两种不同的起沙方案对该沙尘天气过程模拟结果的异同，及其可能的原因。

WRF(Weather Research and Forecasting)模式是新一代的中尺度模式，是一个研究和预报通用的模式系统，现有ARW(Advanced Research WRF)和NMM(Nonhydrostatic Mesoscale Model)两个版本，其中本研究所用的版本是ARW版本(以下简称WRF)。WRF/Chem是在WRF的基础上引进了Chem模块，并实现了气象－化学的全面在线耦合(Grell et al., 2005)。

WRF/Chem模式(3.0版本以来)已经考虑了沙尘的起沙、传输和沉降过程，可用于沙尘的模拟研究。模式中起沙过程所采用的参数化方案是基于GOCART(Georgia Institute of Technology-Goddard Global Ozone Chemistry Aerosol Radiation and Transport model)的方案(Ginoux et al., 2004)，以下将其简称GOCART方案。模式的传输过程采用了质量守恒和标量守恒传输方案，保持大气中沙尘含量的守恒(Grell et al., 2005)。干沉降过程考虑了湍流扩散和分子扩散产生的沉降以及重力沉降，沉降速度分别采用阻力模拟(即空气动力学阻力，粘滞层阻力和边界层阻力)和Stokes公式进行计算(Grell et al., 2005)。由于湿沉降过程的相对复杂性，模式中未考虑沙尘的湿沉降过程。

在原有起沙方案(Ginoux et al., 2004)的基础上，Kang et al.(2011)在上述WRF/Chem(3.0版本)的基础上，引进了3种新起沙参数化方案，并比较了这3种方案的模拟结果差异(但并没有将新方案模拟结果与原起沙方案进行对比)。此后，我们在上述版本的基础上进行了发展，在WRF/Chem(3.2版本)中耦合进去了Shao(2004)起沙方案(以下简称Shao04方案)。以下将介绍GOCART和Shao04两种起沙方案。

GOCART方案根据垂直沙尘通量和风速的经验公式，直接计算各个粒径的起沙通量，第p档沙尘粒子的起沙通量

其中，C是常量，为0.8 μg s2 m−5；S是起沙概率函数(也称为风蚀度指数)，表示一个模式网格内能起沙的网格面积所占的比例；是各档沙尘的质量比例；是地面10 m的水平风速；和分别是第p档沙尘粒子的有效直径和其对应的地面10 m的临界风速。Ginoux et al.(2004)假设，其中为临界摩擦速度，采用Marticorena and Bergametti(1995)公式计算并考虑土壤湿度的影响。

Shao04方案则是基于起沙风蚀物理学理论，认为沙尘粒子的产生是因为跳跃起来的沙粒轰击土壤表面的直接释放过程以及轰击粒子自身的分裂作用，并给出由粒径为di的沙粒运动产生的粒径为dj的沙尘通量如下式(Shao，2004)：

其中，cy是比例系数，ηfj是质量比例，γ是与粒径分布有关的权重因子，g是重力加速度，pm(dj)和pf(dj)分别为粒径的全分布和最小分布，σm是有 效轰击率，是摩擦速度，ρb是土壤容重，p是 土壤的塑性压力。Q(i)为粒径为di的轰击沙粒的通量，采用White(1979)给出的方案计算如公式(3)所示

其中，c为常数，ρa为空气密度。u*t是临界摩擦速度，采用Shao and Lu(2000)方案计算且考虑了土壤湿度和地面粗糙元的影响(Kang et al., 2011)。

GOCART方案主要用于全球沙尘过程的模拟(Ginoux et al., 2004；Huneeus et al., 2011)，而Shao04方案已被应用到区域沙尘过程的模拟中，并展现了不错的模拟能力(Kang et al., 2011)。两方案的粒径范围均是20 μm以下，但粒径分档有所不 同，各粒径分布及相关参数如表 1和表 2所示。

Shao04方案计算起沙通量的过程需要输入的变量包括空气密度、土壤湿度、摩擦速度、叶面积指数(或植被覆盖)、土壤质地、粒径分布数据等。空气密度、土壤湿度、摩擦速度均由模式系统的其他模块计算得到；粒径分布采用Shao(2004)的数据。土壤质地、植被覆盖则利用了与实际情况更为接近的GIS数据集。所用的土壤质地则是张时煌等(2004a)借助于GIS技术融合几份土壤分类资料得到的，该资料符合美国农业部(USDA)土壤质地分类标准。由于WRF/Chem模式中所用的是多年统计得到的平均植被覆盖资料，考虑到植被覆盖存在逐年的变化以及其对于沙尘起沙的影响(Xu et al., 2006)，为了更好地实时模拟沙尘过程，我们利用张时煌等(2004b)借助于GIS等技术，从遥感资料反演得到的植被覆盖。需要说明的是Shao04方案所用的土壤质地和植被覆盖资料已被成功应用于东亚地区沙尘暴过程的模拟(Shao et al., 2003; 孙建华等，2003；雷航等，2005)。GOCART方案所需的空气密度、土壤湿度、10 m风速也同样由模式系统的其它模块得到。

此外，为了限制沙尘源地的范围，两方案均 采用了风蚀度指数，用以指示潜在沙尘源地分 布。Shao04方案应用了基于地形要素(不考虑高山)、多年统计的沙尘浓度(由能见度估算得到)和植被覆盖得到的风蚀度指数Se(Shao，2008)。当风蚀度指数大于或等于1，可作为沙尘源地。GOCART方案所考虑的风蚀度指数Sf是一个基于地形起伏变化并结合陆面植被覆盖得到的定量化因子(而Se为判断因子，不是定量化因子)，代表潜在起沙能力(Ginoux et al., 2004)。比较东亚地区这两种风蚀度指数的差异(图 1)，Se和Sf均基本抓住了东亚区域的主要沙尘源地，包括我国西北、内蒙古西部地区及其相邻的蒙古国地区、我国东北南部的科尔沁沙地，但Sf并没有考虑内蒙古中东部和华北北部(包括浑善达克沙地和河北北部的坝上草原)作为主要沙尘源地，而这可能是重要的沙尘源地(叶笃正等，2000；Gong et al., 2003；Zhang et al., 2003)。

采用上述介绍的沙尘模式系统，对2002年3月19～22日的强沙尘暴过程进行模拟。模拟的时间是2002年3月18日20:00至23日20:00，共5 d。模拟区域(如图 1)的中心点是(35°N，110°E)，空间有220×160个格点，空间分辨率是30 km。模式垂直层为28层。物理过程参数化方案选择如表 3所示。 每天4次、空间分辨率为1°(纬度)×1°(经度)的NCEP/FNL分析资料提供了积分所需的初始场和边界场。此外，积分时还采用grid-nudging四维资料同化方法与NCEP/FNL分析资料逼近。我们分别采用Shao04和GOCART起沙方案进行模拟，共两组试验。试验中未考虑沙尘的反馈作用，因而两组试验模拟结果的差异是由起沙方案的差异造成的。在下一部分，我们首先分析沙尘暴过程对应的环流系统及天气系统的演变，以及模式对它们的模拟能力，接着重点分析两组试验模拟的沙尘浓度的异同及与观测的偏差，然后通过分析这种差异和偏差的原因，给出了两种起沙方案的合理性及其不足。

由于缺少直接的垂直沙尘通量、沉降通量的直接观测资料，一般采用地面观测的颗粒物浓度 [如总悬浮颗粒物(TSP)、PM10浓度等] 来定量比较模式的模拟结果。由于沙尘的空间分布范围往往很大，在一些缺少直接沙尘浓度观测的地方，可以利用气象台站地面观测的能见度来间接反演地面沙尘浓度，而且气象台站地面观测的沙尘天气记录(即浮尘、扬沙、沙尘暴和强沙尘暴)和卫星观测的图像也可用来定性分析沙尘天气的强度和空间范围。尽管通过观测或反演的颗粒物浓度(如PM10或TSP)包含空气中的自然或人类活动释放的所有气溶胶，但在沙尘源地及邻近地区，它们还是在一定程度上反映了沙尘浓度的相对大小，特别是沙尘暴期间沙尘浓度急剧升高，使得PM10和TSP的浓度也急剧升高。

沙尘浓度与能见度的关系已有一些研究(Shao et al., 2003；Wang et al., 2008)。尽管通过能见度估算沙尘浓度的方法存在一些不确定性，但这仍不失为一种可行的途径(Shao，2008)。Shao et al.(2003)提出了一个通过能见度推算TSP浓度的经验公式：

其中，CTSP是沙尘浓度(单位：μg m-3)，Dv是能见度(单位：km)，该公式已被Uno et al.(2006)开展的DMIP计划所采用。

本文中所利用的是中国气象局提供的东亚地区地面气象台站观测数据，每3小时观测能见度等气象要素，并记录包括浮尘、扬沙、沙尘暴和强沙尘暴以及烟尘等在内的天气现象。需要说明的是，考虑到朝鲜观测的沙尘天气强度与上游的蒙古、中国可能存在不一致，即强度偏强(孙建华等，2003)，为了更好地进行对比，在朝鲜观测的沙尘暴和强沙尘暴天气均作为扬沙处理。

为了定量地进行对比，我们采用了中国环境保护部提供的银川、呼和浩特、太原、北京、青岛和沈阳的日平均PM10浓度，其中日平均的采样时间为前一天的12:00至当天的11:00。另外，本文还利用到了由能见度根据公式(4)推算得到的塔中、额济纳旗、乌拉特中旗和锡林浩特的沙尘浓度。除了青岛位于沙尘源地的下游地区，其余站点均位于源地及其附近地区。

在沙尘传输到的更遥远地区，本文利用了东亚酸雨网(Acid Deposition Monitoring Network in East Asia，EANET)观测的Oki、Sado、Tappi和Rishiri的逐小时平均PM10浓度。另外Oki和Rishiri观测的逐小时平均细颗粒物浓度也用于与模拟结果进行对比。这些用于验证模式结果的观测站点分布如图 1所示。

沙尘过程的模拟一定程度上依赖于模式对于气象场的模拟。其中，地面风场直接影响起沙量的计算，大气环流场等会影响沙尘的传输路径和沉降过程。为此我们首先来考察模式对2002年3月19～22日沙尘暴过程对应的东亚地区环流形势及天气系统的模拟情况。

据孙建华等(2003)的分析，此次沙尘暴过程的主要影响天气系统为自中亚向东移动并加强的短波槽和低层的蒙古气旋及其伴随着的地面冷锋。从图 2中可以看到，3月19日08:00 500 hPa高度场上的低压槽位于新疆东部，至20日08:00，低压槽移动较快，移至蒙古东部，而且强度有所增加，20日14:00形成一闭合的低压中心，位势高度为5160 gpm(图略)，此后进一步加强并缓慢东移。21日20:00以后，500 hPa高度场上的低压槽有所减弱，并东移出海。伴随着低压槽的发展东移，蒙古气旋自19～20日有显著加强，20～21日气旋缓慢东移，且强度变化较小，22日以后东移出海。与气旋相伴随的，地面冷锋后大片区域风速很大(达6级风，甚至大风以上)，并观测到沙尘暴天气(图 3)。

从模拟的气象场来看，与分析资料相比，模式很好地模拟出500 hPa高空环流的演变，很好地再现了高空槽脊的移动及短波槽加深的特征(图略)。对于天气系统的模拟，如图 3e−3h所示，可以看出，模式也较好地模拟出了蒙古气旋的移动过程，中心的位置基本一致。高压和低压的强度虽然基本接近，但存在着偏差，即模拟的高压偏低2.5～5 hPa，低压偏低0～2.5 hPa。尽管如此，模式较好地模拟出地面风场的特征，沙尘源区的风速大值区的分布也比较接近。模式对于气象场(特别是环流场和地面风场)的合理再现，可以更好地模拟沙尘起沙和传输过程，从而有利于我们分析不同起沙方案对沙尘过程模拟的影响。

图 4给出了WRF/Chem模拟的地面沙尘浓度及其与地面观测沙尘天气事件的对比。总体上看，采用Shao04和GOCART方案，模式模拟的沙尘分布特征比较相似，均能基本抓住观测沙尘天气的时空演变特征，即沙尘自源地向下游地区传输的特征，其中蒙古南部、内蒙古及东北西部的沙尘向东南传输，新疆南部地区沙尘向西传输，而且模拟的地面沙尘浓度大值区也与观测的沙尘天气强度较强的区域相互吻合。这说明两方案均能抓住了主要的沙尘源地及其起沙过程。

分析两组试验模拟结果的差别，从空间分布范围来看，除了在蒙古东南部及内蒙古中东部地区范围偏小之外，采用GOCART方案模拟的沙尘分布范围更广。从量值上看，GOCART方案下模拟的 沙尘浓度均在5000 μg m−3以下，而Shao04方案下模拟的沙尘浓度在部分地区高于5000 μg m−3，在局地甚至可达20000 μg m−3以上。从不同区域的差别来看，在塔克拉玛干沙漠的绝 大部分地区，GOCART方案下模拟的沙尘浓度更大；在蒙古南部和内蒙古西部，19日08:00，GOCART方案下模拟的沙尘浓度较大，但19日14:00至20:00，采用Shao04方案模拟的沙尘浓度增加显著，在部分地区相比更大，此后20～22日，则以GOCART方案下的模拟值更大；在蒙古东南部和内蒙古中东部、华北北部等地区，Shao04方案下的结果则明显较大。两方案下模拟的沙尘浓度不仅在源区有较大差别，而且在下游也有明显差别，如渤海、黄海等(差异可达1000 μg m−3以上)。

对比台站观测资料，模拟偏差主要表现在以下5个方面：1)两组试验均未能模拟出蒙古西北部观测到的沙尘天气(其强度可达强沙尘暴，如图 4a1、4a2)，这可能与方案并不将此区域设为潜在沙尘源地而未模拟出起沙过程有关；2)采用GOCART方案模拟的沙尘浓度的空间分布范围比观测偏大，一方面，19日在观测的沙尘分布区的边缘也有分布，而采用Shao04方案模拟的沙尘分布范围在19日与观测更为一致；另一方面，采用GOCART方案还模拟出了19日观测未出现的东北西部和南部的沙尘分布，其原因可能是GOCART方案模拟出了东北西部的虚假起沙过程(如图 1b，此处为沙尘源区)；3)20～23日，GOCART方案下模拟的沙尘在内蒙古西北部及中国—蒙古交界地区有虚假分布，浓度较大，可达1000 μg m−3以上，而对应的，Shao04方案下模拟的沙尘分布范围则较小，与观测较为接近；4)采用Shao04方案模式较好地模拟出了20日14:00蒙古东部、内蒙古东部的沙尘分布带(图 4d1)，而采用GOCART方案时未能抓住该特 征(图 4d2)，这可能与GOCART方案未能模拟出该地区的起沙过程有关；5)采用Shao04方案模拟的塔里木盆地的沙尘分布范围比观测偏小，而GOCART方案下模拟的沙尘分布范围与观测相当接近。值得注意的是，在塔里木盆地，观测站点主要沿其周边进行分布，在塔里木盆地东部观测站点数目较少(图 1a)，因而尚不能验证GOCART方案下19日塔里木盆地东部的沙尘分布中心(浓度达1000 μg m−3)是否合理。

本节我们将利用站点观测的地面沙尘浓度的时间变化结果，考察耦合了不同起沙方案的WRF/ Chem模式的模拟性能。总体来看，两方案下模式均能抓住基本的沙尘过程，合理地再现了各个站点的浓度变化特征，但在个别站点也存在一定的偏差。

在位于沙尘源地的塔中、额济纳旗、乌拉特中旗及锡林浩特(图 5)，Shao04方案下模拟的沙尘浓度与观测的相关系数均在0.8以上，通过了99.9%的信度检验，说明模式很好地模拟出TSP浓度的时间变化特征。而GOCART方案下的模拟能力相对偏低，该方案对塔中、额济纳旗、乌拉特中旗模拟的沙尘浓度的相关系数分别为0.71、0.72、0.90，超过了99.9%的信度检验，但GOCART方案对锡林浩特站模拟的时间相关系数较低，仅为0.35。

在塔中，19日11:00至17:00观测到浮尘天气(表 4)，能见度为3～8 km，对应的沙尘浓度为800～1500 μg m−3。采用GOCART方案时模拟的沙尘浓度于19日00:00开始增加，19日12:00降至谷值，然后再次增加，20日02:00达到第二个峰值，其变化趋势与观测相反，而Shao04方案下模拟的沙尘浓度增加缓慢，也未模拟出浮尘天气下的沙尘浓度峰值，只模拟出20日00:00前后弱的沙尘过程(模拟的沙尘浓度最大达200 μg m−3)；对于塔中21日11:00至22日14:00的沙尘过程(最强为沙 尘暴天气)，GOCART方案下模拟的沙尘浓度的变化与观测比较一致，而Shao04方案下模拟的沙尘浓度虽然也显著增加，但减小更为迅速，至22日01:00降至50 μg m−3以下。

在额济纳旗和乌拉特中旗，采用GOCART方案时，尽管模式能抓住沙尘浓度的主要变化特征，但19～20日模拟的额济纳旗的浓度峰值出现得偏早2～3 h，且沙尘过程持续得偏久，另外21～22日模式还模拟出了额济纳旗相对较弱(浓度最大为1222 μg m-3)的沙尘天气事件(观测TSP浓度未明显增加)，即图 4所显示的20～22日内蒙古西北部的模拟偏差(见3.3节分析)。

在位于浑善达克沙地北部的锡林浩特，Shao04方案下模式很好地模拟出观测的沙尘浓度变化特征。但GOCART方案下，模拟的沙尘浓度只在沙尘暴过程的前期(19日20:00至20日12:00)有所增加，其强度偏小(模拟的沙尘浓度为823 μg m−3，而观测的最强沙尘天气为沙尘暴，能见度仅为500 m，对应可估算得TSP浓度为6806 μg m−3)，而且模式也未能模拟出沙尘暴后期的沙尘浓度变化；前期的沙尘浓度增加是因为有沙尘自上游地区(如蒙古南部地区)传输而来，而后期未模拟出沙尘天气则与模式未考虑蒙古国东南部和内蒙古中东部的起沙有关，这与图 4和3.3节的分析是相对应的。

需要指出的是，由于用能见度推算的TSP浓度包含大气中的所有粒子，所以模拟值应普遍低于推算的TSP浓度，但考虑到能见度观测存在一定的主观性，量值仍存在一些不确定性因素，因而我们更关注的是能见度(及由其推算的沙尘浓度)随时间的相对变化特征。

接下来，分析模拟与直接观测的沙尘浓度的对比，并辅以气象台站的观测记录(表 4)进行分析。考虑到背景空气中的颗粒物对观测的PM10有所贡献，而沙尘暴期间PM10主要由沙尘气溶胶组成，因而重点分析沙尘暴期间的对比结果。

在银川(图 6a)，模式能抓住20日和22日沙尘天气下沙尘浓度增加的特征，但GOCART方案下模拟的沙尘浓度分别为观测的2.3和4.7倍，而Shao04方案下模拟结果则与观测比较接近。在呼和浩特(图 6b)，两方案下模拟的20日沙尘浓度分别为902 μg m−3、1188 μg m−3，与观测的1354 μg m−3相比，量值相当，只是有所偏小。21日观测的沙尘浓度为206 μg m−3，GOCART方案下模拟的沙尘浓度略有偏大，而Shao04方案下的模拟值则偏小。在太原(图 6c)，21～22日两方案下模式均抓住了沙尘浓度显著增加的特征，只是量值有所偏小；在其它时段，气象台站观测到烟尘，这可能是由其它局地污染物(如煤的颗粒、生物质燃烧等)或(和)沙尘引起的。

此次沙尘天气在北京开始于20日09:15，此后迅速加剧，并持续至22日，其中20日10:50至15:30沙尘浓度高达12060 μg m−3(张仁健等，2002)(如表 4)。与观测的分析表明，两方案下模拟的日平均沙尘浓度的变化特征与观测均相当符合，但是在量值上整体偏小，观测值为模拟值的1～5倍，其中GOCART方案下的结果偏小更多(图 6d)。在青岛和沈阳，模拟的沙尘浓度在21～22日显著增加，这与观测比较接近(日平均PM10浓度达440 μg m−3以上)，但也存在着偏差，其中Shao04方案下模拟的21日沙尘浓度分别为观测的1.9和1.6倍，而22日两方案下模拟的沙尘浓度均偏小。值得一提的是，沈阳19日的平均沙尘浓度为297 μg m−3，由于采样时间为18日12:00至19日11:00，而18日23:00至19日11:00并未观测到沙尘天气，可见这是18日12:00至18日20:00这个时段内沙尘过程的影响所导致的，而根据地面气象台站资料，该时段内辽宁西部观测到了扬沙天气。虽然GOCART方案下模拟的沙尘浓度与观测比较接近，但由于模拟时间开始于18日20:00，可见这是虚假的沙尘浓度增加过程(如图 4a2)。

在下游的日本地区(图 7)，除了Rishiri站外，Shao04方案下模拟的沙尘浓度更大，与观测更为接近。在Oki和Sado两个站，两方案下模拟的沙尘浓度的主要时间变化特征基本相似，但量值上Shao04方案下模拟的沙尘浓度与观测比较接近，而采用GOCART方案模拟的沙尘浓度则比采用Shao04方案的结果小(与图 4对应)，显著低于观测结果。在Tappi站，观测的沙尘浓度最大值达2210 μg m−3(22日00:00)，这也是日本地区4个站点中观测到的最大值，而且在23日03：00浓度还存在另一峰值，达1156 μg m−3；与观测相比，Shao04方案下模式模拟出了21～23日的两个沙尘浓度峰值，但其中第一个峰值的量值偏小，模拟值只有739 μg m−3，而GOCART方案下只模拟出第一个浓度峰值，且量值只有323 μg m−3。

在Rishiri站，观测的沙尘浓度于21日20:00迅速增加，21日23:00达到峰值后迅速减少，但22日01:00达到谷值后又迅速增加，至22日04:00达到最大，此后才缓慢减少。相比而言，GOCART方案下模拟出了第一个浓度峰值，但浓度偏小2～3倍，同时在第二个峰值前后，沙尘浓度增加较为缓慢。另一方面，采用Shao04方案模拟的沙尘浓度增加过程与观测较为一致，但未出现观测中22日01:00的低值，而且峰值后浓度减小更为迅速，从而模拟的22～23日沙尘浓度偏小得多(图 7d和7f)。

上面比较了两方案下模拟的地面沙尘浓度的空间分布和时间变化的差异，并分析了存在的偏差。由于两组试验的模拟结果差别是由起沙过程的差别造成的(包括起沙总量和粒径分布的差异)。以下我们将分析模拟偏差的原因。

图 8给出了不同时刻两方案模拟的垂直沙尘通量的空间分布图，对比可以看出：从整体的空间分布范围来看，19日14:00两方案模拟的起沙区域均主要位于蒙古南部戈壁和内蒙古西部，差异较大的区域位于新疆塔里木盆地中东部，GOCART方 案模拟的起沙通量达10 μg m−2 s−1以上。对于20日的起沙过程，由于锋面系统经过东北、华北时，地面冷锋后风速较强(图 4b)，Shao04方案模拟的起沙区域主要在蒙古东南部、内蒙古中东部和河北北部等地区，而GOCART方案模拟的起沙区域则主要位于蒙古南部戈壁、内蒙古西部和东北西部。21日14:00，Shao04方案和GOCART方案计算的起沙通量分布与20日14:00比较近似，只是在塔克拉玛干沙漠地区均有所增强。22日14:00，两方案的模拟结果与20日14:00的分布比较近似，但Shao04方案的模拟结果在范围和强度上均有所减小，这与图 4d、4g中地面风速相比减弱的特征是一致的。进一步分析两方案模拟的垂直沙尘通量的空间分布特征和量值差异，可以看出：

GOCART方案模拟的起沙通量的空间分布更为平滑，而Shao04方案模拟的起沙通量的局地特征明显，存在局地极大值区。这种差别的原因在于GOCART中的风蚀度指数是基于1°(纬度)×1.25°(经度)的分布得到的，而Shao04方案利用了0.05°(纬度)×0.05°(经度)的高分辨率陆面数据(植被覆盖、土壤质地和潜在源地分布)来计算起沙过程。

从起沙通量的量值来看，GOCART方案模拟的起沙通量基本在500 μg m−2 s−1以下(如图 8a2，除了19日14:00，蒙古国的南部靠近中国－蒙古边境的一局地地区，达500～1000 μg m−2 s−1)，大部分区域则在100 μg m−2 s−1以下，而Shao04方案的起沙通量在部分地区可达1000 μg m−2 s−1以上，在局地甚至达5000 μg m−2 s−1以上。

总体来看，并结合3.3节和3.4节的分析，表 5总结了不同源区两种方案的模拟结果。首先，两方案均模拟出了19日蒙古南部和内蒙古西部的起沙过程，起沙中心也比较接近，从而模式比较合理地再现了19～20日该区域的沙尘分布，也均抓住了额济纳旗和乌拉特中旗的沙尘浓度变化；其次，GOCART方案未能模拟出蒙古东南部和内蒙古中东部的起沙过程，从而模式未能抓住20～21日该区域的沙尘分布特征，而Shao04方案模拟出了20～21日该区域强度达1000 μg m−2 s−1以上的起沙过程，因而模式很好地抓住该区域的沙尘分布，也较好地模拟出锡林浩特以及下游日本站点地面沙尘浓度的量值变化(见3.2节和3.3节分析)；再次，18日21:00至19日10:00，GOCART方案在东北 西部模拟有起沙过程，中心区域的起沙通量达10～100 μg m−2 s−1(图略)，导致19日东北西部和南部模拟有虚假沙尘分布；最后，在塔克拉玛干沙漠，Shao04方案模拟的起沙量偏小，而GOCART方案则较为合理。此外，两方案均未能模拟出19日蒙古西北部的起沙过程，因而模式未能抓住该区域观测到的沙尘天气过程。

另一方面，GOCART方案计算的起沙区域范围变化较小，但随着天气系统的移动，Shao04方案计算的起沙区域也发生移动。从图 9中满足起沙条件的区域可以看出，GOCART方案由于临界风速普遍比10 m风速小一个量级，起沙区域范围较大(虽然随着天气系统的移动也会移动，但变化不大)，而Shao04方案采用的临界摩擦速度与摩擦速度量级相当，只在摩擦速度相比较大(对应的10 m风速也较大)的区域才能起沙，因而起沙区域随着风速大值区的移动而发生移动，且起沙区域范围也相对较小。

需要指出的是，与Shao04方案相比，GOCART方案模拟的偏多的沙尘起沙区域，其垂直沙尘通量普遍较小，在100 μg m−2 s−1以下，如19日14:00位于塔里木盆地中西部以及20日14:00位于内蒙古西北部的起沙区域(图 8)。这主要是由于在这些区域，摩擦速度较小，不满足Shao04方案起沙条件，但依然满足GOCART方案中较小的临界风速条件，因而发生起沙，但由于风速小，所以垂直沙尘通量也较小。

进一步分析起沙过程风速相关因子的量值(图 9)。在沙尘源区附近，Shao04方案的临界摩擦速度为0.4～1.0 m s−1，源地边缘可到1 m s−1以上，GOCART方案的临界起沙风速在10 m高度上的值为0.4～0.8 m s−1。事实上，在东亚地区开展的站点观测研究显示出：塔克拉玛干沙漠2个站点3～4月的临界摩擦速度为0.48和0.56 m s−1(Ishizuka et al., 2009)，戈壁地区春季的临界摩擦速度为0.54～0.69 m s−1(Kimura and Shinoda， 2010)、0.34～0.42 m s−1(Li and Zhang， 2011)，内蒙古浑善达克沙地、科尔沁沙地的临界摩擦速度则分别为0.4 m s−1(张宏升等，2007)和0.62 m s−1(李晓岚和张宏升，2012)；而Kurosaki and Mikami(2007)基于气象台站观测资料得到塔克拉玛干沙漠、戈壁地区、黄土高原的临界起沙风速在塔克拉玛干沙漠为4.4±0.6 m s−1、8.9±2.2 m s−1、6.9±1.2 m s−1，野外观测的结果则显示出临界起沙风速在戈壁地区为8 m s−1(Park et al., 2010)，在浑善达克沙地、科尔沁沙地则分别为6 m s−1(张宏升等，2007)和9.5 m s−1(李晓岚和张宏升，2012)。尽管不同时期的临界摩擦速度和临界起沙风速可能有所差别，但很明显，GOCART方案的临界起沙风速确实偏小，而Shao04方案的临界摩擦速度则基本合理。

总体说来，Shao04方案模拟的起沙区域较为合理，而GOCART方案模拟的起沙区域则偏大，存在一些虚假起沙区域(表 5)。如GOCART方案模拟了19～20日东北西部持续的起沙过程，从而导致了模式模拟的东北西部和南部存在虚假沙尘分布；20日内蒙古西北部的虚假起沙过程也造成了 模拟的地面沙尘浓度在内蒙古西部和河套地区持续偏强(模拟的沙尘浓度可达500～1000 μg m−3之间，而观测无沙尘天气记录)。

与起沙通量的空间分布特征类似，两方案均模拟出了蒙古南部、中国内蒙古西部、河西走廊和河套地区的起沙过程，差别主要在于：一是GOCART方案未模拟出蒙古东南部和中国内蒙古中东部的起沙；二是在塔克拉玛干沙漠Shao04方案模拟的起沙区域更小；三是GOCART方案在东北西部模拟的起沙区域更大，强度更强(如图 10所示)。GOCART和Shao04方案模拟的5 d总起沙量分别为8.56、17.37 Tg(各粒径的起沙量及其占总起沙量的比例可见表 1和表 2)，粒径越大的粒径档，比例也越大。另外，Shao04模拟的10～20 μm粒径区间起沙量 占的比例超过了一半，达61.4%，而GOCART方案模拟的总起沙量在2.0～20.0 μm 4个粒径区间较为均匀地分布，其中最大的占24.9%，最小的占18.0%，这与质量比例的取值均为0.25有关。

本文利用耦合了GOCART和Shao04两种起沙参数化方案的WRF/Chem模式对2002年3月19～22日发生在东亚地区的强沙尘暴过程进行模拟，考察耦合了不同起沙方案的WRF/Chem对此次沙尘天气过程的模拟能力，并进一步分析了两种不同起沙方案对模式模拟结果的影响及其可能原因。总体来看，耦合了两种不同方案的WRF/Chem总体上均能较合理地模拟出主要的起沙区域、起沙强度的变化以及沙尘浓度的时空演变特征，模式对沙尘源地附近及下游地区地面沙尘浓度时间变化特征的模拟与站点观测结果也十分接近。但相对而言，不管是对沙尘浓度的空间分布还是地面站点沙尘浓度的时间变化，耦合了Shao04方案的WRF/Chem对此次强沙尘暴过程的模拟结果要优于耦合了GOCART方 案下的结果，这主要是由于Shao04方案对沙尘起沙的物理过程有相对更好的描述，从而使得Shao04方案模拟的起沙区域和垂直沙尘通量更为合理。

与实际观测的资料进一步比较表明，耦合了两种不同方案的模式模拟结果与实测相比均存在一定的偏差，而这与不同起沙方案的物理机制紧密相关。就GOCART方案而言，由于该方案中采用的临界风速偏小，从而导致模式在一些区域模拟出虚假的起沙过程，虽然对应的起沙通量并不大(100 μg m−2 s−1以下)。另一方面，风蚀度指数是GOCART方案表征潜在沙尘源地的重要指标，由于GOCART方案中未将蒙古东南部、中国内蒙古中东部等地区考虑作为潜在沙尘源地，因此模式也未能模拟出该区域的起沙过程，导致模式对该区域沙尘浓度的模拟存在偏差，尤其在下游的日本站点，模式模拟的沙尘浓度也偏小。

就Shao04方案而言，该方案下模拟的塔克拉玛干沙漠的沙尘范围与观测相比偏小，且沙尘浓度也偏小，与观测偏差较大。塔里木盆地地区在沙尘过程期间，风速相比其他区域更小，而在未发生轰击作用的情况，空气拖曳力的夹卷作用对起沙过程的影响是不可忽视的。已有研究指出，空气拖曳力的夹卷作用产生的塔里木盆地起沙通量能到30～50 μg m−2 s−1，相对应的底层沙尘浓度可到150～300 μg m−3(Klose and Shao， 2012)，因此Shao04方案模拟结果偏低的可能原因之一，是Shao04方案中未能考虑空气拖曳力的夹卷作用对沙尘起沙的可能贡献。此外，鉴于土壤特性资料的缺乏，Shao04方案中用到的土壤粒径分布的不准确，可能也是影响模式模拟性能的重要指标之一。因此广泛收集整理沙尘起沙参数化方案所需要的相关地表GIS信息，也是提高沙尘模拟与预测的重要途径之一。

致谢 感谢德国科隆大学邵亚平教授提供了沙尘起沙模式以及对本研究工作的帮助。感谢中国气象局提供了观测的地面台站数据，中国环保局、东亚酸雨网(EANET)提供了站点观测的PM10和细颗粒物浓度。