CN114441397A - 基于高斯烟羽和质心定位的pm2.5污染源定位方法 - Google Patents
基于高斯烟羽和质心定位的pm2.5污染源定位方法 Download PDFInfo
- Publication number
- CN114441397A CN114441397A CN202111304272.7A CN202111304272A CN114441397A CN 114441397 A CN114441397 A CN 114441397A CN 202111304272 A CN202111304272 A CN 202111304272A CN 114441397 A CN114441397 A CN 114441397A
- Authority
- CN
- China
- Prior art keywords
- coordinates
- calculate
- positioning
- above step
- calculation
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
- 238000000034 method Methods 0.000 title claims abstract description 45
- 238000009792 diffusion process Methods 0.000 claims abstract description 63
- 238000012544 monitoring process Methods 0.000 claims abstract description 54
- 238000004422 calculation algorithm Methods 0.000 claims abstract description 35
- 238000010845 search algorithm Methods 0.000 claims abstract description 11
- 238000003921 particle size analysis Methods 0.000 claims abstract description 10
- 238000007781 pre-processing Methods 0.000 claims abstract description 8
- 238000004364 calculation method Methods 0.000 claims description 42
- 230000009466 transformation Effects 0.000 claims description 9
- 238000006243 chemical reaction Methods 0.000 claims description 8
- 230000006872 improvement Effects 0.000 claims description 4
- NAWXUBYGYWOOIX-SFHVURJKSA-N (2s)-2-[[4-[2-(2,4-diaminoquinazolin-6-yl)ethyl]benzoyl]amino]-4-methylidenepentanedioic acid Chemical compound C1=CC2=NC(N)=NC(N)=C2C=C1CCC1=CC=C(C(=O)N[C@@H](CC(=C)C(O)=O)C(O)=O)C=C1 NAWXUBYGYWOOIX-SFHVURJKSA-N 0.000 claims description 3
- 230000008569 process Effects 0.000 claims description 3
- 239000000126 substance Substances 0.000 claims description 3
- 238000002474 experimental method Methods 0.000 description 6
- 238000010586 diagram Methods 0.000 description 5
- 235000019580 granularity Nutrition 0.000 description 5
- 230000004807 localization Effects 0.000 description 4
- 239000002245 particle Substances 0.000 description 4
- 238000004458 analytical method Methods 0.000 description 3
- 238000001514 detection method Methods 0.000 description 3
- 230000015572 biosynthetic process Effects 0.000 description 2
- 230000000694 effects Effects 0.000 description 2
- 239000013618 particulate matter Substances 0.000 description 2
- 230000009286 beneficial effect Effects 0.000 description 1
- 230000008859 change Effects 0.000 description 1
- 238000004891 communication Methods 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 238000005516 engineering process Methods 0.000 description 1
- 230000002068 genetic effect Effects 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 238000004088 simulation Methods 0.000 description 1
- 239000007787 solid Substances 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01N—INVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
- G01N15/00—Investigating characteristics of particles; Investigating permeability, pore-volume or surface-area of porous materials
- G01N15/06—Investigating concentration of particle suspensions
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01C—MEASURING DISTANCES, LEVELS OR BEARINGS; SURVEYING; NAVIGATION; GYROSCOPIC INSTRUMENTS; PHOTOGRAMMETRY OR VIDEOGRAMMETRY
- G01C21/00—Navigation; Navigational instruments not provided for in groups G01C1/00 - G01C19/00
- G01C21/20—Instruments for performing navigational calculations
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/20—Design optimisation, verification or simulation
- G06F30/27—Design optimisation, verification or simulation using machine learning, e.g. artificial intelligence, neural networks, support vector machines [SVM] or training a model
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06N—COMPUTING ARRANGEMENTS BASED ON SPECIFIC COMPUTATIONAL MODELS
- G06N3/00—Computing arrangements based on biological models
- G06N3/004—Artificial life, i.e. computing arrangements simulating life
- G06N3/006—Artificial life, i.e. computing arrangements simulating life based on simulated virtual individual or collective life forms, e.g. social simulations or particle swarm optimisation [PSO]
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Remote Sensing (AREA)
- Radar, Positioning & Navigation (AREA)
- Evolutionary Computation (AREA)
- Health & Medical Sciences (AREA)
- General Health & Medical Sciences (AREA)
- Software Systems (AREA)
- Chemical & Material Sciences (AREA)
- General Engineering & Computer Science (AREA)
- Life Sciences & Earth Sciences (AREA)
- Artificial Intelligence (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Dispersion Chemistry (AREA)
- Pathology (AREA)
- Immunology (AREA)
- Medical Informatics (AREA)
- Biochemistry (AREA)
- Computer Hardware Design (AREA)
- Geometry (AREA)
- Analytical Chemistry (AREA)
- Automation & Control Theory (AREA)
- Biomedical Technology (AREA)
- Biophysics (AREA)
- Computational Linguistics (AREA)
- Data Mining & Analysis (AREA)
- Molecular Biology (AREA)
- Computing Systems (AREA)
- Mathematical Physics (AREA)
- Management, Administration, Business Operations System, And Electronic Commerce (AREA)
Abstract
本发明提供了基于高斯烟羽和质心定位的PM2.5污染源定位方法，包括步骤S1空气质量监测站点位置数据预处理；S2基于气体扩散模型和加权质心定位算法，计算加权质心定位结果；S3基于改进的差分进化搜索算法，计算精细定位结果；S4以H小时为粒度分析定位结果，计算最终定位结果；S5利用高斯烟羽模型描述PM2.5污染源气体扩散。本发明方法结合了加权质心定位算法性能高和差分进化搜索算法准确度高的优点，提高了算法整体的性能和准确度。
Description
技术领域
本发明的技术方案涉及高斯烟羽模型和加权质心定位算法技术领域，特 别是设计到基于高斯烟羽和质心定位的PM2.5污染源定位方法。
背景技术
应用计算机技术来进行PM2.5污染源定位的算法可以很好的实时监控 和发现超标排放PM2.5污染源。现今，国内外有很多的团队对污染源定位算 法进行研究。比如：无风情况下的质心定位算法结合遗传算法、粒子滤波算 法等。其中在无风情况下，快速的进行污染源定位以及污染气体扩散的模拟， 质心定位算法性能最好且定位点数据和真实值也较为接近，高斯烟羽模型应 用最为广泛。目前也有有风情况下的加权质心定位算法，但其并不能利用多 条不同时段的数据来进行更为精准的定位，仅通过一个时刻的数据进行计算， 并不能排除数据的偶然性和定位点的精确性。
现有定位算法的应用，都在模拟环境中进行实验，没有很好的利用现有 大量真实的PM2.5浓度数据，没有实际的应用到PM2.5污染源定位当中去， 本方法解决了以上不足之处。
发明内容
有鉴于此，本发明旨在提出基于高斯烟羽和质心定位的PM2.5污染源定 位方法，以提供一种基于气体湍流扩散模型的加权质心定位的PM2.5污染源 定位方法，气体湍流模型是一种简化的高斯烟羽模型，克服了现有高斯烟羽 模型的计算复杂性的缺点；提供一种基于改进的差分进化搜索方法进行精确 定位以及多粒度分析的定位方法，克服了质心定位的缺陷；提供一种基于高 斯烟羽扩散模型的方法描述PM2.5的污染源扩散分布。
为达到上述目的，本发明的技术方案是这样实现的：
在基于气体湍流扩散模型的加权质心定位的基础上，基于改进的差分进 化搜索方法进行精确定位并从小时、天等多个粒度分析定位结果，最后基于 高斯烟羽扩散模型描述气体扩散，具体步骤如下：
第一步，空气质量监测站点位置数据预处理：
(1.1)空气质量监测站点经纬度坐标转换为笛卡尔平面坐标：
使用米勒投影法利用公式(1)和公式(2)实现地球经纬度坐标转换为笛卡 尔平面坐标(x，y)，
x＝(W/2)+(W/2π)×L1×π/180 (1)，
y＝(H/2)-(H/2m)×1.25×ln tan(0.25π+0.4L2×π/180) (2)，
其中L1和L2分别为经度和维度，W和H分别为地球的周长和地球周 长的一半，m为米勒投影常数，大约在正负2.3之间，(x，y)为笛卡尔平面坐 标。由此完成对空气质量监测站点经纬度坐标到笛卡尔平面坐标的转换。
(1.2)空气质量监测站点笛卡尔平面坐标转换：
得到上述(1.1)步的空气质量监测站点笛卡尔平面坐标，然后根据风向数 据以正北方向逆时针旋转的角度，对空气质量监测站点的笛卡尔平面坐标做 变换，使x轴的正方向为正北方向。由此完成对空气质量监测站点笛卡尔平 面坐标的转换。
至此完成空气质量监测站点位置数据预处理，得到预处理后的空气监测 站点笛卡尔平面坐标，以下简称为站点坐标；
第二步，基于气体扩散模型和加权质心定位算法，计算加权质心定位结 果：
(2.1)坐标转换，建立以风向为x轴正方向的新坐标系：
利用公式(3)对站点坐标进行坐标转换，
其中θ是风向与正北方向所形成的夹角。
由此完成建立以风向为x轴正方向的新坐标系。
(2.2)基于气体扩散模型计算加权质心定位算法的权值因子：
在上述(2.1)步所建立的新坐标系下，根据气体湍流扩散模型，利用公式 (4)计算加权质心定位算法的权值因子qx(i)和qy(i)，
其中C(i)是第i个监测站点的监测浓度，C(l)是监测站点中的最大监测 浓度，v是风速，k是扩散参数，α和β是指定参数。
由此完成加权质心定位算法的权值因子qx(i)和qy(i)的计算。
(2.3)基于加权质心定位算法计算初始定位坐标：
在上述(2.1)步所建立的新的坐标系下，根据加权质心定位算法，利用公 式(5)计算初始定位坐标x0’和y0’，
由此完成在上述(2.1)步所建立的新坐标系的初始定位坐标(x0’,y0’)的计 算。
(2.4)坐标转换，将新坐标系的初始定位坐标转换为原坐标系的初始定位 坐标(x0,y0)：
利用公式(6)将上述(2.3)的初始定位坐标转换回上述(1.2)步的坐标系，
由此完成在上述(1.2)步的坐标系的初始定位坐标(x0,y0)的计算。
至此完成加权质心定位结果的计算，以下简称为初始定位坐标；
第三步，基于改进的差分进化搜索算法，计算精细定位结果
(3.1)初始化，计算第一代群体G1：
以初始定位结果为中心，取十分之一站点坐标最大横纵坐标值为边长， 生成N个均匀分布的第一代个体xi(1)组成第一代群体G1。
(3.2)变异，通过差分策略计算变异群体G1’：
对上述(3.1)步计算的第一代群体G1中随机的选取两个个体，将其向量 差进行缩放后与随机选择的第三个个体求和产生变异个体，利用公式(7)计算，
vi(g+1)＝xr1(g)+F×(xr2(g)-xr3(g)) (7),
其中r1，r2，r3是三个互不相等的随机数，取值范围是[1,N]。F是缩放 因子，g是进化代数。在进化过程中，保证变异个体不超出站点坐标的边界。
由此完成由变异个体vi(1)组成的变异群体G1’的计算。
(3.3)交叉，将个体与变异个体进行交叉操作，计算交叉群体G1”:
对上述(3.1)步计算的第一代群体G1与上述(3.2)步计算的交叉群体G1’进 行交叉操作，利用公式(8)计算，
其中CR是交叉概率，取值为[0，1]；jrand随机取1或2，确保ui(g+1) 中至少从上述(3.2)步的变异个体vi(g+1)中获取一个参数。
由此完成由交叉个体ui(1)组成的交叉群体G1”的计算。
(3.4)改进，基于气体扩散模型计算搜索群体S1：
对原差分进化搜索算法进行改进，对上述(3.3)步计算的交叉群体G1”进 行进一步的计算，基于气体扩散模型使交叉群体G1”向预估的污染源点方向 进行偏移，加快收敛速度。以上述(3.3)步的交叉个体ui为气体源位置，计算 监测站点理论浓度与监测站点实际浓度的差异，将该差异相加求平均值，随 机取一个小于1的权值λ，利用公式(9)计算搜索个体si，si的坐标为(xs,ys)，
其中xk和yk是搜索移动轨迹。
由此完成由搜索个体si组成的搜索群体Si的计算。
(3.5)选择，计算下一代群体G2：
结合气体源扩散模型，利用公式(10)和公式(11)计算上述(3.3)步的搜索群 体Si与上述(3.1)步的第一代群体G1的适应值，采用贪婪选择的策略，利用 公式(11)选择较优的个体，
由此完成由第二代个体xi(2)组成第二代群体G2。
(3.6)进化，计算精细定位结果Mi：
使进化代数g＝g+1，判断进化代数g是否达到最大进化代数G。若判 断通过，则终止进化，使最佳的个体Mi作为精细定位结果输出；若判断不通 过，则返回执行上述(3.2)步。
由此完成精细定位结果Mi的计算。
至此完成精细定位结果的计算，以下称为精细定位坐标；
第四步，以H小时为粒度分析定位结果，计算最终定位结果。
选择当前时刻前H小时的数据输入第二步，计算出上述(3.2)步的H个 精细化定位坐标的平均值作为最终定位结果。
至此完成最终定位结果的计算，以下称为最终定位坐标；
第五步，利用高斯烟羽模型描述PM2.5污染源气体扩散：
计算出上述第四步所述的最终定位坐标作为污染源点，将污染源气体扩 散假设为连续点源扩散，基于高斯烟羽模型利用公式(12)进行污染气体扩散 的描述：
其中C(x,y,z,H)为任意点泄漏气体的浓度，单位为mg/m3；Q为泄漏源强 度，指单位时间泄露点源排放的气体量,单位为mg/s；C(x,y,z)为泄漏物质 在(x,y,z)点处的质量浓度，单位为mg/m3；x为下风向距离，单位为m；y 为横风向距离，单位为m；z为到地面的距离，单位为m；σy和σz,分别为 y轴和z轴方向的扩散参数；u为平均风速，单位为m/s；H为泄漏点源的有 效高度，单位为m。
其中扩散参数可由表1得出。
其中对于区域浓度的描述采用梯度划分，描述相对于污染源的百分比。
根据高斯烟羽模型，对所测区域进行网格划分，计算每一个网格内的浓 度数据，将所得结果利用等高线绘制。
至此完成污染源定位坐标计算和PM2.5的污染源扩散的描述。
相对于现有技术，本发明所述的基于高斯烟羽和质心定位的PM2.5污染 源定位方法具有以下有益效果：
(1)本发明方法结合了加权质心定位算法性能高和差分进化搜索算法 准确度高的优点，提高了算法整体的性能和准确度。
(2)本发明方法从小时、天等多粒度分析污染源定位结果，充分利用了 局部定位结果与整体定位结果，实现了动态定位，进一步提高了定位的准确 度。
(3)本发明方法利用高斯烟羽模型可视化的描述污染源的气体扩散情 况，结合现有能够良好拟合PM2.5扩散规律的算法，实际的应用到PM2.5 污染源定位当中去，使得现有的定位方法进入可应用的领域。
附图说明
构成本发明的一部分的附图用来提供对本发明的进一步理解，本发明的 示意性实施例及其说明用于解释本发明，并不构成对本发明的不当限定；在 附图中：
图1为本发明方法的流程图；
图2为本发明方法预处理部分中的坐标转换原理示意图；
图3为本发明方法中加权质心定位原理图；
图4为本发明方法中差分进化搜索算法原理示意图；
图5为本发明方法中根据高斯烟羽模型描述的气体扩散示意图；
图6为本发明方法中各代群体最小适应值曲线图；
图7为本发明方法中四种粒度定位实例一；
图8为本发明方法中四种粒度定位相对误差率变化一；
图9为本发明方法中四种粒度定位实例二；
图10为本发明方法中四种粒度定位相对误差率变化二。
具体实施方式
需要说明的是，在不冲突的情况下，本发明中的实施例及实施例中的特 征可以相互组合；
在本发明的描述中，需要理解的是，术语“中心”、“纵向”、“横向”、 “上”、“下”、“前”、“后”、“左”、“右”、“竖直”、“水平”、 “顶”、“底”、“内”、“外”等指示的方位或位置关系为基于附图所示 的方位或位置关系，仅是为了便于描述本发明和简化描述，而不是指示或暗 示所指的装置或元件必须具有特定的方位、以特定的方位构造和操作，因此 不能理解为对本发明的限制；此外，术语“第一”、“第二”等仅用于描述 目的，而不能理解为指示或暗示相对重要性或者隐含指明所指示的技术特征 的数量；由此，限定有“第一”、“第二”等的特征可以明示或者隐含地包 括一个或者更多个该特征；在本发明的描述中，除非另有说明，“多个”的 含义是两个或两个以上；
在本发明的描述中，需要说明的是，除非另有明确的规定和限定，术语 “安装”、“相连”、“连接”应做广义理解，例如，可以是固定连接，也 可以是可拆卸连接，或一体地连接；可以是机械连接，也可以是电连接；可 以是直接相连，也可以通过中间媒介间接相连，可以是两个元件内部的连通； 对于本领域的普通技术人员而言，可以通过具体情况理解上述术语在本发明 中的具体含义；
下面将参考附图并结合实施例来详细说明本发明；
图1所示实施例表明，本发明方法的流程是：空气质量监测点位置数据 预处理→根据气体扩散模型和质心定位算法，计算加权质心定位结果→根据 差分进化搜索算法，计算精细定位结果→多粒度分析定位结果，计算最终定 位结果→利用高斯烟羽模型描述污染源的扩散。
图2所示实施例表明，本发明方法空气质量监测站点经纬度坐标转换为 笛卡尔平面坐标部分中，使用米勒投影法实现地球经纬度坐标转换为笛卡尔 平面坐标(x,y)，其中地球的周长是W，W为12762744πm，米勒投影常数为 2.3。
图3所示实施例表明，本发明方法基于加权质心定位算法计算初始定位 坐标中，A点、B点、C点、D点和E点是空气质量监测站点，(x,y)是初始 定位点，根据气体扩散模型确定空气质量监测站点的影响力，就确定(x,y)为 这些空气监测站点所组成的多边形的质心，即初始定位坐标。
图4所示实施例表明，本发明方法基于气体扩散模型和改进的差分进化 搜索方法计算搜索群体中，图中以星型点为气体源位置，以空心点为空气监 测站点位置，实心点为预测气体源位置，利用气体扩散模型计算监测站点理 论浓度与监测站点实际浓度的差值，将该差值相加求平均值。
图5所示实施例表明，本发明方法利用高斯烟羽模型进行PM2.5污染区 域的绘制中，x轴正方向是自然风的风向，坐标原点O是泄露源点，气体在 横向y轴方向，纵向z轴方向的扩散浓度分布是高斯分布。
实施例1
本实例的基于高斯烟羽和质心定位的PM2.5污染源定位方法，是一种在 基于气体湍流扩散模型的加权质心定位的基础上，基于改进的差分进化搜索 方法进行精确定位并多粒度的分析定位结果，最后基于高斯烟羽扩散模型描 述气体扩散的方法，具体步骤如下：
第一步，空气质量监测点位置数据预处理：
(1.1)空气质量监测站点经纬度坐标转换为笛卡尔平面坐标：
使用米勒投影法利用公式(1)和公式(2)实现地球经纬度坐标转换为笛卡 尔平面坐标(x，y)，
x＝(W/2)+(W/2π)×L1×π/180 (1)，
y＝(H/2)-(H/2m)×1.25×ln tan(0.25π+0.4L2×π/180) (2)，
其中L1和L2分别为经度和维度，W和H分别为地球周长和地球周长一 半，m为米勒投影常数，(x，y)为笛卡尔平面坐标，由此完成对空气监测站 点经纬度坐标到笛卡尔平面坐标的转换，其中W＝12762744πm，(L1，L2)∈{(116.417,39.929),(116.407,39.886),(116.339,39.929),(116.352,39.878), (116.397,39.937),(116.461,39.987)}，m＝2.3。
(1.2)空气监测站点笛卡尔平面坐标转换：
得到上述(1.1)步的空气质量监测站点笛卡尔平面坐标，然后根据风向数 据是以正北方向逆时针旋转的角度，对空气质量监测站点的笛卡尔平面坐标 做变换，使x轴的正方向为正北方向。由此完成对空气质量监测站点笛卡尔 平面坐标的转换。
至此完成空气质量监测站点位置数据预处理，得到预处理后的空气监测 站点笛卡尔平面坐标，以下简称为站点坐标；第二步，基于气体扩散模型和 质心定位算法，计算质心定位结果：
(2.1)坐标转换，建立以风向为x轴正方向的新坐标系：
利用公式(3)对站点坐标进行坐标转换，
其中θ是风向与正北方向所形成的夹角。
由此完成建立以风向为x轴正方向的新的坐标系。
(2.2)基于气体扩散模型计算加权质心定位算法的权值因子：
在上述(2.1)步所建立的新的坐标系下，根据气体湍流扩散模型，利用公式 (4)计算加权质心定位算法的权值因子qx(i)和qy(i)，
其中C(i)是第i个监测站点的监测浓度，C(l)是监测站点中的最大监测浓 度，v是风速，k是扩散参数，α和β是指定参数，其中v＝8m/s，k＝4000， α＝1，β＝1。
由此完成加权质心定位算法的权值因子qx(i)和qy(i)的计算。
(2.3)基于加权质心定位算法计算初始定位坐标：
在上述(2.1)步所建立的新的坐标系下，根据加权质心定位算法，利用公式 (5)计算初始定位坐标x0’和y0’，
由此完成在上述(2.1)步所建立的新坐标系的初始定位坐标(x0’,y0’)的计算。
其中v＝8m/s，α＝1，β＝1。
由此完成在上述(2.1)步所建立的新的坐标系的初始定位坐标(x0’,y0’)的计 算。
(2.4)坐标转换，将新坐标系的初始定位坐标转换为原坐标系的初始定位 坐标(x0,y0)：
利用公式(6)将上述(2.3)的初始定位坐标转换回上述(1.2)步的坐标系，
由此完成在上述(1.2)步的坐标系的初始定位坐标(x0,y0)的计算。
至此完成加权质心定位结果的计算，以下简称为初始定位坐标；
第三步，基于改进的差分进化搜索算法，计算精细定位结果
(3.1)初始化，计算第一代群体G1：
以初始定位结果为中心，取十分之一站点坐标最大横纵坐标值为边长， 生成N个均匀分布的第一代个体xi(1)组成第一代群体G1。
(3.2)变异，通过差分策略计算变异群体G1’：
对上述(3.1)步计算的第一代群体G1中随机的选取两个个体，将其向量差 进行缩放后与随机选择的第三个个体求和产生变异个体，利用公式(7)计算，
vi(g+1)＝xr1(g)+F×(xr2(g)-xr3(g)) (7),
其中r1，r2，r3是三个互不相等的随机数，取值范围是[1,N]。F是缩放 因子，g是进化代数。在进化过程中，保证变异个体不超出站点坐标的边界， 其中F＝0.1。
由此完成由变异个体vi(1)组成的变异群体G1’的计算。
(3.3)交叉，将个体与变异个体进行交叉操作，计算交叉群体G1”:
对上述(3.1)步计算的第一代群体G1与上述(3.2)步计算的交叉群体G1’进 行交叉操作，利用公式(8)计算，
其中CR是交叉概率，取值为[0，1]；jrand随机取1或2，确保ui(g+1) 中至少从上述(3.2)步的变异个体vi(g+1)中获取一个参数。
由此完成由交叉个体ui(1)组成的交叉群体G1”的计算。
(3.4)改进，基于气体扩散模型计算搜索群体S1：
对原差分进化搜索算法进行改进，对上述(3.3)步计算的交叉群体G1”进 行进一步的计算，基于气体扩散模型使交叉群体G1”向预估的污染源点方向 进行偏移，加快收敛速度。以上述(3.3)步的交叉个体ui为气体源位置，计算 监测站点理论浓度与监测站点实际浓度的差异，将该差异相加求平均值，随 机取一个小于1的权值λ，利用公式(9)计算搜索个体si，si的坐标为(xs,ys)，
其中xk和yk是搜索移动轨迹，λ＝0.01。
由此完成由搜索个体si组成的搜索群体Si的计算。
(3.5)选择，计算下一代群体G2：
结合气体源扩散模型，利用公式(10)和公式(11)计算上述(3.3)步的搜索群 体Si与上述(3.1)步的第一代群体G1的适应值，采用贪婪选择的策略，利用 公式(11)选择较优的个体，
由此完成由第二代个体xi(2)组成第二代群体G2。
(3.6)进化，计算精细定位结果Mi：
使进化代数g＝g+1，判断进化代数g是否达到最大进化代数G。若判断 通过，则终止进化，使最佳的个体Mi作为精细定位结果输出；若判断不通 过，则返回执行上述(3.2)步。
由此完成精细定位结果Mi的计算。
至此完成精细定位结果的计算，以下称为精细定位坐标；
第四步，以H小时为粒度分析定位结果，计算最终定位结果。
选择当前时刻前H小时的数据输入第二步，计算出上述(3.2)步的H个精 细化定位坐标的平均值作为最终定位结果，其中H＝3*24。
至此完成最终定位结果的计算，以下称为最终定位坐标；
第五步，利用高斯烟羽模型描述PM2.5污染源气体扩散：
计算出上述第四步所述的最终定位坐标作为污染源点，将污染源气体扩 散假设为连续点源扩散，基于高斯烟羽模型利用公式(12)进行污染气体扩散 的描述：
其中C(x,y,z,H)为任意点泄漏气体的浓度，单位为mg/m3；Q为泄漏源强 度，指单位时间泄露点源排放的气体量,单位为mg/s；C(x,y,z)为泄漏物质 在(x,y,z)点处的质量浓度，单位为mg/m3；x为下风向距离，单位为m；y 为横风向距离，单位为m；z为到地面的距离，单位为m；σy和σz,分别为 y轴和z轴方向的扩散参数，m；u为平均风速，单位为m/s；H为泄漏点源 的有效高度，单位为m。
其中扩散参数可由表1得出。
表1.扩散系数计算公式
其中Q＝100mg/s，采用梯度划分，得到的区域浓度并不是真正的污染浓 度，而是相对于污染源的百分比，其中C＝50mg/m3，是相当于源强Q的50％。
根据高斯烟羽模型，对所测区域进行网格划分，计算每一个网格内的浓 度数据，将所得结果利用等高线绘制。
至此完成污染源定位坐标计算和PM2.5的污染源扩散的描述。
本实施例在全国监测点历史空气质量和气象数据集上进行了实验。其中 全国监测点历史空气质量包含7296条数据，共有6个月的数据，包括每小 时平均PM2.5浓度和24h平均PM2.5浓度；全国气象数据集包含3648条数 据，共6个月的数据，包括风向和风速。在Windows10系统下的Python和 MATLAB平台上进行实验。在以H小时为粒度分析定位结果，计算最终定 位结果中，选取H＝3*24，本实例的准确率为94.13％。如图6、图7、图8和 表1所示列出了本实施例的实验的结果。
如图6所示，在选择并计算下一代群体和进化并计算精细定位结果中， 设置最大进化代数G＝200时，适应值趋于收敛。
如图7所示，在以H小时为粒度分析定位结果，计算最终定位结果中， H＝1和H＝24的定位距离污染源较远，H＝3*24和H＝7*24的定位相近，距离 污染源较近，说明用于污染源定位的数据越多，定位结果越准确。
如图8所示，在以H小时为粒度分析定位结果，计算最终定位结果中， H＝1、H＝24和H＝3*24的定位相对误差率波动较大，H＝7*24的定位相对误 差率比较平稳。
表2. 3月份四种不同粒污染源定位相对误差率及平均相对误差率
表2表明，在以H小时为粒度分析定位结果，计算最终定位结果中， H＝1、H＝24的平均相对误差率较高，H＝3*24和H＝7*24的平均相对误差率 较低，污染源定位效果最好。
实施例2
为了验证本发明方法在检测点位分布不均匀的情况下，污染源定位的有 效性，本实施例在某城市监测点历史空气质量和气象数据集上进行了实验。 其中某城市监测点历史空气质量包含10986条数据，共有3个月的数据，包 括每小时平均PM2.5浓度和24h平均PM2.5浓度；全国气象数据集包含1824 条数据，共3个月的数据，包括风向和风速。在Windows10系统下的Python 和MATLAB平台上进行实验。在以H小时为粒度分析定位结果，计算最终 定位结果中，选取H＝3*24，本实例的准确率为92.04％。图4、图5和表3 列出了本实施例的实验的结果。
如图9所示，在以H小时为粒度分析定位结果，计算最终定位结果中， 在检测点位分布不均匀的情况下H＝3*24和H＝7*24的定位结果依然比较相 近，并且距离污染源较近。
如图10所示，在以H小时为粒度分析定位结果，计算最终定位结果中， 在检测点位分布不均匀的情况下H＝1、H＝24和H＝3*24的定位相对误差率 波动仍然较大，H＝7*24的定位相对误差率波动最小。
表3. 3月份四种不同粒污染源定位相对误差率及平均相对误差率
表3表明，在以H小时为粒度分析定位结果，计算最终定位结果中，相 对于实例1，H＝3*24的准确率超过H＝7*24的准确率，污染源定位效果最 好。
以上所述仅为本发明的较佳实施例而已，并不用以限制本发明，凡在本 发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在 本发明的保护范围之内。
Claims (1)
1.基于高斯烟羽和质心定位的PM2.5污染源定位方法，其特征在于：包括以下步骤：
S1、空气质量监测站点位置数据预处理：
S101、空气质量监测站点经纬度坐标转换为笛卡尔平面坐标：
使用米勒投影法利用公式(1)和公式(2)实现地球经纬度坐标转换为笛卡尔平面坐标(x，y)，
x＝(W/2)+(W/2π)×L1×π/180 (1)，
y＝(H/2)-(H/2m)×1.25×ln tan(0.25π+0.4L2×π/180) (2)，
其中L1和L2分别为经度和维度，W和H分别为地球的周长和地球周长的一半，m为米勒投影常数，大约在正负2.3之间，(x，y)为笛卡尔平面坐标；由此完成对空气质量监测站点经纬度坐标到笛卡尔平面坐标的转换；
S102、空气质量监测站点笛卡尔平面坐标转换：
得到上述S101步的空气质量监测站点笛卡尔平面坐标，然后根据风向数据是以正北方向逆时针旋转的角度，对空气质量监测站点的笛卡尔平面坐标做变换，使x轴的正方向为正北方向；由此完成对空气质量监测站点笛卡尔平面坐标的转换；
至此完成空气质量监测站点位置数据预处理，得到预处理后的空气监测站点笛卡尔平面坐标，以下简称为站点坐标；
S2、基于气体扩散模型和加权质心定位算法，计算加权质心定位结果：
S201、坐标转换，建立以风向为x轴正方向的新坐标系：
利用公式(3)对站点坐标进行坐标转换，
其中θ是风向与正北方向所形成的夹角；
由此完成建立以风向为x轴正方向的新坐标系；
S202、基于气体扩散模型计算加权质心定位算法的权值因子：
在上述S201步所建立的新坐标系下，根据气体湍流扩散模型，利用公式(4)计算加权质心定位算法的权值因子qx(i)和qy(i)，
其中C(i)是第i个监测站点的监测浓度，C(l)是监测站点中的最大监测浓度，v是风速，k是扩散参数，α和β是指定参数；
由此完成加权质心定位算法的权值因子qx(i)和qy(i)的计算；
S203、基于加权质心定位算法计算初始定位坐标：
在上述S201步所建立的新的坐标系下，根据加权质心定位算法，利用公式(5)计算初始定位坐标x0’和y0’，
由此完成在上述S201步所建立的新坐标系的初始定位坐标(x0’，y0’)的计算；
S204、坐标转换，将新坐标系的初始定位坐标转换为原坐标系的初始定位坐标(x0，y0)：
利用公式(6)将上述S203的初始定位坐标转换回上述S102步的坐标系，
由此完成在上述S102步的坐标系的初始定位坐标(x0，y0)的计算；
至此完成加权质心定位结果的计算，以下简称为初始定位坐标；
S3、基于改进的差分进化搜索算法，计算精细定位结果：
S301、初始化，计算第一代群体G1：
以初始定位结果为中心，取十分之一站点坐标最大横纵坐标值为边长，生成N个均匀分布的第一代个体xi(1)组成第一代群体G1；
S302、变异，通过差分策略计算变异群体G1’：
对上述S301步计算的第一代群体G1中随机的选取两个个体，将其向量差进行缩放后与随机选择的第三个个体求和产生变异个体，利用公式(7)计算，
vi(g+1)＝xr1(g)+F×(xr2(g)-xr3(g)) (7)，
其中r1，r2，r3是三个互不相等的随机数，取值范围是[1，N]；F是缩放因子，g是进化代数；在进化过程中，保证变异个体不超出站点坐标的边界；
由此完成由变异个体vi(1)组成的变异群体G1’的计算；
S303、交叉，将个体与变异个体进行交叉操作，计算交叉群体G1”：
对上述S301步计算的第一代群体G1与上述S302步计算的交叉群体G1’进行交叉操作，利用公式(8)计算，
其中CR是交叉概率，取值为[0，1]；jrand随机取1或2，确保ui(g+1)中至少从上述S302步的变异个体vi(g+1)中获取一个参数；
由此完成由交叉个体ui(1)组成的交叉群体G1”的计算；
S304、改进，基于气体扩散模型计算搜索群体S1：
对原差分进化搜索算法进行改进，对上述S303步计算的交叉群体G1”进行进一步的计算，基于气体扩散模型使交叉群体G1”向预估的污染源点方向进行偏移，加快收敛速度；以上述S303步的交叉个体ui为气体源位置，计算监测站点理论浓度与监测站点实际浓度的差异，将该差异相加求平均值，随机取一个小于1的权值λ，利用公式(9)计算搜索个体si，si的坐标为(xs，ys)，
其中xk和yk是搜索移动轨迹；
由此完成由搜索个体si组成的搜索群体Si的计算；
S305、选择，计算下一代群体G2：
结合气体源扩散模型，利用公式(10)和公式(11)计算上述S303步的搜索群体Si与上述S301步的第一代群体G1的适应值，采用贪婪选择的策略，利用公式(11)选择较优的个体，
由此完成由第二代个体xi(2)组成第二代群体G2；
S306、进化，计算精细定位结果Mi：
使进化代数g＝g+1，判断进化代数g是否达到最大进化代数G；若判断通过，则终止进化，使最佳的个体Mi作为精细定位结果输出；若判断不通过，则返回执行上述S302步；
由此完成精细定位结果Mi的计算；
至此完成精细定位结果的计算，以下称为精细定位坐标；
S4、以H小时为粒度分析定位结果，计算最终定位结果；
选择当前时刻前H小时的数据输入第二步，计算出上述S302步的H个精细化定位坐标的平均值作为最终定位结果；
至此完成最终定位结果的计算，以下称为最终定位坐标；
S5、利用高斯烟羽模型描述PM2.5污染源气体扩散：
计算出上述第四步所述的最终定位坐标作为污染源点，将污染源气体扩散假设为连续点源扩散，基于高斯烟羽模型利用公式(12)进行污染气体扩散的描述：
其中C(x，y，z，H)为任意点泄漏气体的浓度，单位为mg/m3；Q为泄漏源强度，指单位时间泄露点源排放的气体量，单位为mg/s；C(x，y，z)为泄漏物质在(x，y，z)点处的质量浓度，单位为mg/m3；x为下风向距离，单位为m；y为横风向距离，单位为m；z为到地面的距离，单位为m；σy和σz，分别为y轴和z轴方向的扩散参数；u为平均风速，单位为m/s；H为泄漏点源的有效高度，单位为m；
其中对于区域浓度的描述采用梯度划分，描述相对于污染源的百分比；
根据高斯烟羽模型，对所测区域进行网格划分，计算每一个网格内的浓度数据，将所得结果利用等高线绘制，完成污染源定位坐标计算和PM2.5的污染源扩散。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202111304272.7A CN114441397B (zh) | 2021-11-05 | 2021-11-05 | 基于高斯烟羽和质心定位的pm2.5污染源定位方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202111304272.7A CN114441397B (zh) | 2021-11-05 | 2021-11-05 | 基于高斯烟羽和质心定位的pm2.5污染源定位方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN114441397A true CN114441397A (zh) | 2022-05-06 |
| CN114441397B CN114441397B (zh) | 2023-07-11 |
Family
ID=81362366
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202111304272.7A Active CN114441397B (zh) | 2021-11-05 | 2021-11-05 | 基于高斯烟羽和质心定位的pm2.5污染源定位方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN114441397B (zh) |
Cited By (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN115792137A (zh) * | 2023-01-17 | 2023-03-14 | 河北先河环保科技股份有限公司 | 大气污染溯源方法及装置、终端 |
Citations (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20170091350A1 (en) * | 2015-09-24 | 2017-03-30 | Alexander Bauer | Near real-time modeling of pollution dispersion |
| CN107917987A (zh) * | 2017-11-09 | 2018-04-17 | 北京伟瑞迪科技有限公司 | 一种面向城市空气污染物溯源分析方法 |
| CN113065098A (zh) * | 2021-03-23 | 2021-07-02 | 北京航天创智科技有限公司 | 基于卫星遥感和地气系统数据的pm2.5浓度检测方法及装置 |
| CN113484204A (zh) * | 2021-06-02 | 2021-10-08 | 济南东之林智能软件有限公司 | 大气污染溯源方法、装置、电子设备和计算机可读介质 |
-
2021
- 2021-11-05 CN CN202111304272.7A patent/CN114441397B/zh active Active
Patent Citations (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20170091350A1 (en) * | 2015-09-24 | 2017-03-30 | Alexander Bauer | Near real-time modeling of pollution dispersion |
| CN107917987A (zh) * | 2017-11-09 | 2018-04-17 | 北京伟瑞迪科技有限公司 | 一种面向城市空气污染物溯源分析方法 |
| CN113065098A (zh) * | 2021-03-23 | 2021-07-02 | 北京航天创智科技有限公司 | 基于卫星遥感和地气系统数据的pm2.5浓度检测方法及装置 |
| CN113484204A (zh) * | 2021-06-02 | 2021-10-08 | 济南东之林智能软件有限公司 | 大气污染溯源方法、装置、电子设备和计算机可读介质 |
Non-Patent Citations (1)
| Title |
|---|
| 张鑫鹏;贺庆;刘佳;祝连庆;: "光度计与颗粒计数相结合的颗粒质量浓度检测方法研究", 半导体光电, no. 03 * |
Cited By (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN115792137A (zh) * | 2023-01-17 | 2023-03-14 | 河北先河环保科技股份有限公司 | 大气污染溯源方法及装置、终端 |
Also Published As
| Publication number | Publication date |
|---|---|
| CN114441397B (zh) | 2023-07-11 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN103234883B (zh) | 一种基于道路交通流量实时估算中心城区pm2.5浓度的方法 | |
| Lu et al. | Wind power evaluation and utilization over a reference high-rise building in urban area | |
| CN109086534B (zh) | 一种基于cfd流体力学模型的风电场尾流订正方法及系统 | |
| CN102663251B (zh) | 基于计算流体力学模型的风电场功率物理预测方法 | |
| WO2018161626A1 (zh) | 用于计算风电场发电量的方法和设备 | |
| US20190370418A1 (en) | Method and apparatus for arranging wind turbines based on rapid accessment fluid model and wake model | |
| Song et al. | Optimization of wind farm micro-siting for complex terrain using greedy algorithm | |
| CN108763825B (zh) | 一种模拟复杂地形的风场的数值模拟方法 | |
| Yan et al. | An improved brain storming optimization algorithm for estimating parameters of photovoltaic models | |
| Han et al. | Large eddy simulation for atmospheric boundary layer flow over flat and complex terrains | |
| CN109931903A (zh) | 一种基于改进鲸鱼优化算法的圆柱度评定方法 | |
| CN114580310A (zh) | 一种基于palm实现wrf模拟风场降尺度处理的方法 | |
| CN115329691A (zh) | 一种基于cfd与gis的超大城市风环境模拟方法 | |
| CN104112167A (zh) | 可发电风资源分布的获取方法 | |
| CN107919983A (zh) | 一种基于数据挖掘的天基信息网络效能评估系统及方法 | |
| CN114441397A (zh) | 基于高斯烟羽和质心定位的pm2.5污染源定位方法 | |
| Liu et al. | Cost and capacity optimization of regional wind-hydrogen integrated energy system | |
| CN109426901A (zh) | 一种中长期用电预测方法及装置 | |
| CN111985691B (zh) | 一种风电场升压站选址方法 | |
| CN115906487A (zh) | 基于Lasso回归分析的城区PM2.5污染扩散建模方法 | |
| CN113312836B (zh) | 一种短期风速预测方法 | |
| Zhao et al. | Three-dimensional Interval Identification of Permanent Magnet Spherical Motor Based on Improved Deep Neural Network | |
| CN116307018A (zh) | 一种基于wrf模式敏感性调参的风速预测方法及系统 | |
| CN104008305B (zh) | 用于千万千瓦风电基地的可发电风资源分布估计方法 | |
| Ghermandi et al. | Model comparison in simulating the atmospheric dispersion of a pollutant plume in low wind conditions |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant | ||
| OL01 | Intention to license declared |