CN111353634B - 一种基于模式初始场同化的二氧化硫排放源优化方法 - Google Patents

Abstract

本发明公开了一种基于模式初始场同化的二氧化硫排放源优化方法，属于大气污染数值模拟领域，包括如下步骤：步骤1：获取二氧化硫地面观测数据、气象观测数据、气象再分析格点数据和“自下而上”方法统计的排放源；步骤2：基于大气化学模式和三维变分(3Dvar)同化方法对二氧化硫地面观测资料进行同化，得到“最优的”二氧化硫浓度场和二氧化硫预报误差；步骤3：将二氧化硫预报误差转化为二氧化硫排放源误差，进而得到最优的二氧化硫排放源并进行验证。

Description

一种基于模式初始场同化的二氧化硫排放源优化方法

技术领域

本发明涉及空气质量模式预报技术领域，是一种基于模式初始场同化的二氧化硫排放源优化方法，属于大气污染数值模拟领域。

背景技术

二氧化硫是主要的大气污染物之一，也是空气质量模式预报的重要内容，但由于模式的初始场、沉降扩散等动力过程描述、气溶胶化学反应过程描述以及二氧化硫排放源准确性等因素的影响，二氧化硫预报仍存在较大的误差。

在上述因素中，排放源的不确定性被认为是影响二氧化硫预报准确性的主要因素。提高排放源的准确性，尤其是提高排放总量和时空分布的准确性，会对提高空气质量预报的准确度产生积极的影响。

二氧化硫排放源制作的传统方法主要是“自下而上”的方法，该方法需要首先统计大量企业和个体的活动水平、精准测算各单元的排放因子，估计出排放清单并计算二氧化硫排放强度的空间分布，然后针对空气质量模式的参数，在空间上将排放源进行三维网格化，在时间上将排放强度进一步细化，从而制作出空气质量模式运行需要的排放源。

利用传统方法制作排放源存在两个不足，一是制作过程需要耗费大量的人力物力资源，而随着经济社会的快速发展，实际排放源的变化速度很快，但难以进行快速有效地更新。二是将排放源进行三维网格化和时间区分的过程存在很大的不确定性，尤其在时间区分方面，“自下而上”方法统计的排放源只包含排放的日总量，并不包含小时变化信息，需要对排放源的小时因子(每小时排放量占日排放量的比重)进行假设才能实现时间区分，因此给排放源引入了更多的误差。

目前已有两种通过资料同化改进排放源的方法，一种是集合卡尔曼滤波方法，但该方法需要大量的集合成员，每个成员都是一次大规模的数值预报模拟，因此利用该方法进行排放源的反演对计算资源和计算机性能的要求非常高，难以满足大范围、长时间的模拟需求。第二种是四维变分法(4Dvar)，通过建立伴随模型 (inverse modelling)，可以实现对排放源的优化调整，然而4Dvar同样需要大量的计算资源，且伴随模式的编写复杂，维护和升级极其困难，应用起来会有很多限制。

本发明首次利用三维变分同化，通过改进模式的初始场，转化为对排放源的改进，从而大大减少了计算量。

发明内容

本发明的目的是提供一种基于模式初始场同化的二氧化硫排放源优化方法，能有效减小中国区域二氧化硫排放源的不确定性，得到逐小时、优化的二氧化硫排放源，且可以用于大气化学模式研究，提高模式的预报效果，克服以上各种传统方法的不足。

为达到上述目的，本发明采用的技术方案如下：

本发明涉及的基于三维变分同化(3Dvar)方法，对传统的“自下而上”方法统计的二氧化硫排放源进行优化的有效方法。具体来说是寻求一种通过3Dvar 同化地面观测资料得到二氧化硫浓度的最“优”初值场和预报误差，通过建立排放源误差和预报误差的关系，反演得到排放源的误差，从而降低排放源不确定性的一种方法。

本发明的基于模式初始场同化的二氧化硫排放源优化方法的设计思想是：

先获取二氧化硫地面观测数据、气象观测数据、气象再分析格点数据和“自下而上”方法统计的排放源；

再基于大气化学模式和三维变分(3Dvar)同化方法对二氧化硫地面观测资料进行同化，得到“最优的”二氧化硫浓度场和二氧化硫预报误差；

最后将二氧化硫预报误差转化为二氧化硫排放源误差，进而得到最优的二氧化硫排放源并进行验证。

具体来说，本发明的基于模式初始场同化的二氧化硫排放源优化方法，包括步骤：

步骤1：收集研究区域内的观测数据，这些观测数据包括：污染物观测数据、气象观测数据、气象再分析资料和先验排放源数据；

步骤2：利用大气化学模式和3Dvar同化方法对污染物观测数据进行逐小时循环同化，得到“最优”的二氧化硫浓度分布场和二氧化硫浓度的预报误差；

二氧化硫预报包括二氧化硫浓度分布场和二氧化硫浓度。

步骤3：使用气象观测资料剔除不符合微风无降水条件的二氧化硫浓度预报误差数据样本，从而确保剩余样本中的二氧化硫浓度的预报误差是由排放源的不确定性引起；

步骤4：对研究区域内每个网格的任意时刻的二氧化硫的预报误差样本求中值；

步骤5：根据大气化学模式中二氧化硫排放源和环境浓度响应之间的转换关系：

式中，V_m＝22.4×10^-3m^-3为气体摩尔体积，ρ为实际空气密度，单位为kgm^-3，ρ_air＝1.29kgm^-3为标准状态下的空气密度，即气体摩尔体积，ΔS为单位面积，Δz为模式层高度，Δt＝1hr为时间；

将二氧化硫的预报误差转化为二氧化硫的排放源误差：

最后将二氧化硫的排放源误差叠加到先验排放源上，得到优化的二氧化硫排放源。

步骤1中为减小观测数据异常值对反演结果的影响，对所有观测数据进行极值检验、一致性检验，剔除异常的观测资料。

步骤3中通过筛选特定气象条件下符合要求的二氧化硫预报误差样本，确定由排放源导致模式污染预报误差贡献的方法。

步骤4中为提高二氧化硫浓度预报误差的代表性，对研究区域内每个网格任意时刻的二氧化硫的预报误差样本求中值。二氧化硫浓度预报误差与二氧化硫排放源误差转换关系的建立，基于模式中二氧化硫排放源和环境浓度响应之间的转换关系，将二氧化硫浓度预报误差反推回二氧化硫排放源误差。

有益效果

本发明的有益效果是：本发明提供的基于模式初始场同化的二氧化硫排放源优化方法同其它基于资料同化方法降低排放源不确定性的方法一样，均在“自下而上”排放清单的基础上，基于实际观测的污染物浓度资料，有效的降低了排放清单的不确定性，从而提高了模式预报水平。

但本发明相对于其他方法而言，具有计算量小，计算速度快，程序维护和升级相对简单的优点，可以进行长时间、大范围的排放源反演优化。此外，利用优化结果，也可以对排放源的变化进行评估，为减排政策的制定、政策落实情况和效果评估提供数据支撑。

附图说明

图1为本发明的基于模式初始场同化的二氧化硫排放源优化方法流程图。

图2为实施例中，2015年10月10日～11月10日采用MEIC-2010、优化的二氧化硫排放源和逐小时循环同化预报的二氧化硫浓度对比图，单位：μg m^-3。

具体实施方式

为了使本发明的目的、技术方案及优点更加清楚明白，以下结合实际例子，对本发明方法进行进一步详细说明。应当理解，此处所描述的具体实施例仅仅用以解释本发明，并不用于限定本发明。

一种基于模式初始场同化的二氧化硫排放源优化方法，它包括步骤如下：

步骤1：收集中国环境监测总站(http://www.cnemc.cn/)观测的逐小时二氧化硫地面浓度资料和中国气象数据共享网(http://data.cma.cn)提供的中国国家级地面自动观测站的逐小时气象观测数据，为减小观测数据异常值对反演结果的影响，并对所有观测数据进行极值检验、一致性检验等质量控制，剔除异常的观测资料；

收集NECP(美国气象环境预报中心)提供的每6小时一次、分辨率为 1.0°×1.0°的气象再分析资料(FNL)；

收集由清华大学发展的中国多尺度清单排放(MEIC)。

步骤2：利用区域大气化学模式进行模拟，得到初始时刻化学场的浓度分布。

步骤3：使用3Dvar同化系统，利用二氧化硫地面观测资料对步骤2得到的二氧化硫浓度初始场(x⁰)进行同化，得到最优的二氧化硫浓度分析场(x^0a)和二氧化硫浓度增量场(δx⁰＝x^0a-x⁰，即二氧化硫的预报误差)。

步骤4：利用步骤3中得到的二氧化硫浓度分析场(x^0a)和先验排放源进行1 小时预报，得到1小时后的二氧化硫浓度预报场(x¹)，类似于步骤3，利用二氧化硫地面观测资料对二氧化硫浓度预报场(x¹)进行同化，得到二氧化硫浓度分析场(x^1a)和二氧化硫浓度预报误差，其中δx¹＝x^1a-x¹，δx¹同时记为该时刻二氧化硫浓度的预报误差。

步骤5：重复步骤4中逐小时预报、同化、再预报的过程。

步骤6：使用气象观测和模式模拟的风速和降雨资料及模式本身模拟的气象资料对模拟时间段内的逐小时二氧化硫预报误差(即δx¹)进行筛选，剔除大风和降水区域的二氧化硫预报误差数据。

步骤7：对研究区域内每个网格任意时刻的二氧化硫的预报误差样本

求中值：

式中，i，j分别为模式内东西、南北网格数，t代表每天的整点时刻，取值范围为1～24，n为模拟的天数，即为样本数量。

步骤8：区域大气化学模式中二氧化硫排放源与二氧化硫浓度之间的响应关系得到预报误差δx¹与排放源误差δE⁰的关系：

式中，V_m＝22.4×10^-3m^-3为气体摩尔体积，ρ为实际空气密度，单位为kg m^-3，ρ_air＝1.29kg m^-3为标准状态下的空气密度，即气体摩尔体积，ΔS为单位面积，Δz为模式层高度，Δt＝1hr为时间。所以可以根据模式单位网格中浓度的误差反演出排放源的误差：

将公式(1)带入公式(3)，即可得到二氧化硫的排放源误差：

式中，V_m＝22.4×10^-3m^-3为气体摩尔体积；

ρ_(i，j.t)为每个网格时刻t的实际空气密度，单位为kg m^-3；

ρ_air＝1.29kgm^-3为标准状态下的空气密度，即气体摩尔体积；

δz_(i，j)为每个网格的模式层高度，δt＝1hr为时间；

最后将二氧化硫排放源误差

叠加到先验排放源上，实现对先验排放源的优化并得到最优的二氧化硫排放源。

经过上述8个步骤，能够实现二氧化硫排放源的反演估计，得到最优的二氧化硫排放源，有效降低排放清单的不确定性。

实施例

以优化2015年10月10日-11月10日一个月的二氧化硫排放源为例，结合附图1所示步骤，进行二氧化硫排放源的优化，研究区域包括中国地区，分辨率为27km，具体包括如下步骤：

步骤1：收集试验时间内的二氧化硫地面观测数据和气象观测数据，对数据进行极值检验、时间一致性检验和空间一致性检验等质量控制。收集研究时间内的FNL资料。利用人为排放源清单MEIC-2010制作模式所用的排放源，并作为先验排放源。根据WRF-Chem模式的研究区域及网格分辨率，将落在同一网格内的观测资料进行求均值处理。(WRF-Chem模式是由美国国家海洋大气局 (NOAA)预报系统实验室(FSL)开发的，气象模式(WRF)和化学模式(Chem) 在线完全耦合的区域空气质量模式。)

步骤2：从10月1日00时(世界时，下同)起，进行为期10天的预报，该预报过程的目的是为研究时段的初始时刻(10月10日00时)生成化学初始场，采用MEIC-2010作为先验排放源。

步骤3：使用3Dvar同化系统，利用二氧化硫地面观测资料对步骤2得到的二氧化硫浓度初始场(x⁰)进行同化，得到二氧化硫浓度分析场(x^0a)和二氧化硫浓度增量场(δx⁰)。

步骤4：利用步骤3中得到的二氧化硫浓度分析场(x^0a)和先验排放源进行1 小时预报，得到1小时后的二氧化硫浓度预报场(x¹)，类似于步骤3，利用二氧化硫地面观测资料对二氧化硫浓度预报场(x¹)进行同化，得到二氧化硫浓度分析场(x^1a)和二氧化硫浓度预报误差，其中δx¹＝x^1a-x¹，δx¹同时记为该时刻二氧化硫浓度的预报误差。

步骤5：重复步骤4中逐小时预报、同化、再预报的过程，直到研究时段结束(2015年的11月10日00时)，得到2015年10月的逐小时二氧化硫浓度预报误差。

步骤6：使用气象观测和模式模拟的风速和降雨资料，对模拟时间段内的逐小时二氧化硫预报误差(即δx¹)进行筛选，剔除大风和降水区域的预报误差数据。具体实施方法为:利用中国自动观测站的降水资料，将站点数据插值到网格中，剔除有降水的站点周围四个格点的数据。对于风速超过4m s^-1，散度超过 10^-4s^-1的格点数据进行剔除。对于不符合晴空小风条件的区域进行时间插值，弥补缺失资料的区域。

步骤7：根据公式(1)对研究区域内每个网格任意时刻的二氧化硫的预报误差样本

求中值。

步骤8：利用WRF-Chem模式中二氧化硫排放源与二氧化硫浓度之间的响应关系，得到预报误差与排放源误差的关系，并反演得到二氧化硫的排放源误差。将二氧化硫排放源误叠加到先验排放源上，得到最优的二氧化硫排放源。

最后分别采用先验排放源(MEIC-2010)和本发明优化的二氧化硫排放源进行模拟分析，验证本方法优化的二氧化硫排放源的准确性。

本发明利用WRF-Chem模式，通过同化地面二氧化硫观测数据，对二氧化硫排放源进行了优化，建立了一种基于模式初始场同化方法的二氧化硫排放源优化方法；对2015年10月二氧化硫排放源的优化，并式用模式验证了优化的二氧化硫排放源的准确性；通过采用本方法优化的二氧化硫排放源进行模拟预报，有效提高了二氧化硫预报的准确性，具有重要的科学意义和推广应用价值；对现有的排放源评估统计技术的有效补充，方法简单易行，计算量小，经费投入少，能有效支持中国地区二氧化硫排放源的统计，为空气质量评估、减排政策制定提供依据。

Claims (5)

1.一种基于模式初始场同化的二氧化硫排放源优化方法，其特征在于包括步骤：

步骤1：收集研究区域内的观测数据，这些观测数据包括：污染物观测数据、气象观测数据、气象再分析资料和先验排放源数据；

步骤2：利用大气化学模式和3Dvar同化方法对污染物观测数据进行逐小时循环同化，得到最优的二氧化硫浓度分布场，以及二氧化硫浓度的预报误差；

步骤3：利用气象观测资料剔除不符合气象条件的二氧化硫浓度预报的误差数据样本；

步骤4：对研究区域内每个网格的任意时刻的二氧化硫硫浓的预报误差样本求中值；

步骤5：根据大气化学模式中二氧化硫排放源和二氧化硫浓度之间的响应关系得到预报误差δx¹与排放源误差δE⁰的关系：

式中，V_m＝22.4×10^-3m^-3为气体摩尔体积；ρ为实际空气密度，单位为kg m^-3；ρ_air＝1.29kg m^-3为标准状态下的空气密度，即气体摩尔体积；ΔS为单位面积；Δz为模式层高度，Δt＝1hr为时间；

将二氧化硫浓度的预报误差转化为二氧化硫的排放源误差：

最后将二氧化硫的排放源误差叠加到先验排放源上，得到优化的二氧化硫排放源。

2.根据权利要求1所述的一种基于模式初始场同化的二氧化硫排放源优化方法，其特征在于步骤3中，所述气象条件是指微风无降水条件。

3.根据权利要求1所述的一种基于模式初始场同化的二氧化硫排放源优化方法，其特征在于所述步骤2包括：

2.1)利用区域大气化学模式进行模拟，得到初始时刻二氧化硫的浓度分布；

2.2)使用3Dvar同化系统，利用二氧化硫地面观测资料对步骤2.1)得到的二氧化硫浓度初始场x⁰进行同化，得到最优的二氧化硫浓度分析场x^0a和二氧化硫浓度增量场δx⁰＝x^0a-x⁰，即二氧化硫的预报误差；

2.3)利用步骤2.2)中得到的二氧化硫浓度分析场x^0a和先验排放源进行1小时预报，得到1小时后的二氧化硫浓度预报场x¹；

类似于步骤2.2)利用二氧化硫地面观测资料对二氧化硫浓度预报场x¹进行同化，得到二氧化硫浓度分析场x^1a和二氧化硫浓度预报误差，其中δx¹＝x^1a-x¹，δx¹同时记为该时刻二氧化硫浓度的预报误差；

2.4)重复步骤2.3)中逐小时预报、同化、再预报的过程。

4.根据权利要求1所述的一种基于模式初始场同化的二氧化硫排放源优化方法，其特征在于所述步骤3中，使用气象观测和模式模拟的风速和降雨资料及模式本身模拟的气象资料对模拟时间段内的逐小时二氧化硫预报误差即δx¹进行筛选，剔除大风和降水区域的二氧化硫预报误差数据。

5.根据权利要求1所述的一种基于模式初始场同化的二氧化硫排放源优化方法，其特征在于所述步骤4中，对研究区域内每个网格任意时刻的二氧化硫浓度的预报误差样本

求中值：

式中，i，j分别为模式内东西、南北网格数，t代表每天的整点时刻，取值范围为1～24，n为模拟的天数，即为样本数量；

步骤5中，将公式(3)带入公式(2)，得到二氧化硫排放源误差：

式中，V_m＝22.4×10^-3m^-3为气体摩尔体积；

ρ_(i，j.t)为每个网格时刻t的实际空气密度，单位为kg m^-3；

ρ_air＝1.29kg m^-3为标准状态下的空气密度，即气体摩尔体积；

δz_(i，j)为每个网格的模式层高度，δt＝1hr为时间；

最后将二氧化硫排放源误差

叠加到先验排放源上，得到最优的二氧化硫排放源。

Images