An atmospheric leakage pollution source identification algorithm combining Gaussian adjoint and regularization

Abstract

The present invention provides an algorithm for identifying single pollutant sources in urban environments using gas pollutant monitors, belonging to the technical field of atmospheric pollutant source identification. The method includes the following steps: obtaining local wind speed and direction from a meteorological station and collecting concentration and location data from the monitor; deriving Gaussian adjoint plume and puff equations and correcting coordinates; combining the Gaussian adjoint plume with a multi-monitor joint probability equation and substituting it into the wind field and monitor data to determine high-probability areas where pollution sources may be present; partitioning the high-probability areas and applying the Gaussian adjoint puff equation to obtain a response matrix between the potential source and the monitor; substituting the response matrix and the monitor concentration into a regularization algorithm to obtain the inverse source intensity of the potential source in each partition; multiplying the source intensity by the response matrix to calculate the forward simulated concentration at the monitor; and finally, substituting the forward simulated concentration into a Bayesian probability equation to calculate the posterior probability distribution of the potential source. The present invention introduces adjoint theory into a Gaussian diffusion model, using the adjoint form of the Gaussian equation to obtain high-probability areas for pollution sources and the response relationship between the pollution source and the monitor. The method then combines historical concentration records from the monitor to rapidly identify the location and intensity of atmospheric pollution sources. The calculation process is simplified, the calculation efficiency is improved, and it is convenient for implantation of miniaturized handheld devices.

Description

(1) Technical field

[0001] The present invention belongs to the technical field of atmospheric pollutant source identification, specifically an algorithm for identifying single pollutant sources in urban environments using gas pollutant monitors. This method can quickly and easily identify the location and release intensity of gaseous pollutant sources, providing support for source control of pollutants. (2) Background technology

[0002] Existing pollutant source identification methods are mainly divided into two categories according to their principles: forward calculation method and reverse calculation method. The forward calculation method realizes source identification by analyzing the matching degree between actual monitoring data and predicted pollutant response data, such as Russell et al. [1] The gradient optimization method is used to ensure that the robot's step length remains unchanged each time it moves, and the next movement direction is determined according to the concentration gradient of the pollutant, constantly approaching the source location. Lilienthal et al. [2] A Braitenberg vehicle source location identification algorithm is proposed. First, gaseous pollutant sensors are installed on both sides of the tool vehicle. Then, the steering of the vehicle is controlled based on the difference in concentration values ​​of the two sensors, so that it gradually moves to the high concentration position close to the source position. Such a method is often time-consuming and inefficient. The reverse calculation method is to reversely solve the parameters of the pollution source based on the relationship between the monitored concentration value and the concentration response at the corresponding moment. Compared with the forward calculation, the reverse calculation can quickly determine the location of the pollution source, the release intensity and even the pollutant release law through simulation calculation when there is enough information. For example, Skaggs and Kabala [3] Using the quasi-reversible method, the groundwater pollution source was identified by solving the convection-diffusion equation. Neupauer and Wilson [4] Introducing the adjoint method into the identification of groundwater pollutant source location and release time can identify the pollution source parameters with less prior information. [5] et al. applied the adjoint method to the air environment and combined multiple sensor data to identify the source of indoor air pollutants. [6] et al. applied the adjoint method to outdoor environments, using it to obtain the response relationship between pollution sources and sensors. Combining this with Bayesian methods, they were able to identify outdoor pollution source parameters. While inverse calculation methods are highly efficient for identifying pollution sources, they often require complex computational models and difficult to obtain the response relationship.

[0003] In summary, the existing pollution source identification methods have the following problems: (1) the forward calculation efficiency is low and is only applicable to specific occasions; (2) the reverse calculation method requires a complex calculation process, and the acquisition of the response relationship requires simulation calculations of all potential source areas and analysis of their posterior probabilities. In order to overcome the above shortcomings, the present invention introduces the adjoint theory into the Gaussian diffusion model, obtains the high-probability area of ​​the pollution source and the response relationship between the pollution source and the monitor through the adjoint form of the Gaussian equation, and quickly identifies the location and intensity of the atmospheric pollution source in combination with the historical concentration records of the monitor. This simplifies the calculation process, improves the calculation efficiency, and facilitates the implantation of miniaturized handheld devices.

[0004] [1] Russell R A. Robotic location of underground chemical sources. Robotica, 2004, 22(1): 109-115.

[0005] [2] Lilienthal A, Duckett T. Experimental analysis of gas-sensitive Braitenberg vehicles. Advanced Robotics, 2004, 18(8): 817-834.

[0006] [3]Todd H. Skaggs, ZJ Kabala. Recovering the History of a GroundwaterContaminant Plume: Method of Quasi-Reversibility. Water Resources Research, 1995, 31(11): 2669-2673.

[0007] [4] RM Neupauer, JL Wilson. Backward probabilistic model of groundwater contamination in non-uniform and transient flow. Advances in Water Resources, 2002, 25(7): 733-746.

[0008] [5] Liu X, Zhai Z. Inverse modeling methods for indoor airborne pollutant tracking: literature review and fundamentals. Indoor Air, 2007, 17: 419-38.

[0009] [6]Keats, A., E.Yee, and F.-S.Lien.Bayesian inference for sourcedetermination with applications to a complex urban environment.AtmosphereEnvironment, 2007, 41(3), 465-479. (3) Summary of the invention

[0010] Technical problems solved

[0011] As cities continue to grow in size, rapidly locating air pollution sources at the city block scale is crucial for protecting the health and safety of residents. Existing pollution source identification methods are inefficient and computationally complex. To address these shortcomings, this paper proposes an atmospheric leakage source identification method that combines Gaussian adjoints with regularization, enabling rapid, convenient, and accurate identification of atmospheric pollution sources.

[0012] Technical Solution

[0013] 1. An atmospheric leakage pollution source identification algorithm combining Gaussian adjoint and regularization, characterized by comprising the following steps:

[0014] Step (1): Based on the data provided by the local meteorological station, obtain the wind direction α and wind speed v and the real-time atmospheric stability w of the measured area;

[0015] Step (2): When the pollutant concentration information detected by the pollutant monitor reaches a stable state, select N monitors to extract data, and try to make the connection line of the N sensors perpendicular to the wind speed direction; then, read the historical pollutant concentration information and location coordinates detected by the N monitors, and record the data as S1, S2, ..., S n ;

[0016] Step (3): Apply adjoint theory to the outdoor concentration control equation, and derive the adjoint form of the Gaussian plume and puff equation through adjoint flow field, variable conversion and Laplace transform:

[0017]

[0018]

[0019] C * is the accompanying concentration of pollutants, kg / m 3 ; Q is the release intensity of accompanying pollution sources, kg / s; Q S is the total mass released by the pollution source, kg; σ x , σ y , σ z are the diffusion coefficients in the x, y, and z directions, respectively, which are related to the atmospheric stability and downwind distance, m; x is the diffusion distance in the downwind direction, m; y and z are the distances perpendicular to the downwind direction, m; u is the wind speed in the flow field, m / s;

[0020] Step (4): As shown in step (3), the x used in the Gaussian adjoint equation is the diffusion distance relative to the monitor in the downwind direction of the adjoint wind field, and y and z are the distances perpendicular to the downwind direction. In actual situations, the spatial information provided is a Cartesian coordinate system with the north direction as the y-axis and the monitor as the origin. Therefore, the coordinate information needs to be corrected according to the wind direction α:

[0021]

[0022] x=d*cos(α2), y=d*sin(α2) (1-4)

[0023] α0 is the angle of wind speed relative to the north direction; α1 is the angle of the line between the pollution source and the monitor relative to the north direction; α2 is the angle between the line between the pollution source and the monitor and the wind speed direction; x0, y0 are the Cartesian coordinates of the monitor;

[0024] Step (5): Integrate the information of the N pollutant monitors in step (2) through the posterior probability equation of multi-sensor combination, calculate the posterior probability of the area to be measured, and take the area greater than the third quartile as the high probability area;

[0025] Step (6): Divide the high-probability area determined in step (5) into M intervals using the Thiessen polygon method. Assume that the pollution source exists at the center of each interval. The response matrix A of the N monitors corresponding to the center point of the M intervals can be obtained by equation (1-2). 11 、A 12 ,...,A 21 ,...,A nm ; Using the information S1 of monitor 1 and the response matrix A 11 、A 12 ,...,A 1m , the inverse source strength q1, q2, ..., q of M potential pollution sources corresponding to monitor 1 can be obtained by regularization method.m ;

[0026] Step (7): Inversely calculate the source strength q1, q2, ..., q m and the response matrix A 21 、A 22 ,...,A 31 ,...,A nm Substituting into the equation C = Aq, we can get the simulated concentration of M potential pollution sources corresponding to N-1 monitors. Then Substitute into the equation Calculate the likelihood function of multiple monitors and finally substitute it back into the Bayesian equation to get the maximum probability source location of the leakage source. is the concentration detected by monitor i, is the potential source location M l The simulated concentration for monitor i, σ 2 is the variance of the sensor.

[0027] 2. In step (5), the posterior probability equation of the multi-sensor combination used is: Among them, x i , u, α are the position and monitoring concentration of the i-th monitor, the dominant wind speed and direction in the region, M0 is the release amount of the pollution source; C * (x; x i , u, α) is the accompanying probability of the pollution source location calculated by the i-th monitor through equation (1-1); Under the conditions that the release amount of the pollution source is M0 and the source location is x, the concentration of the monitor is The posterior probability of

[0028] 3. In step (7), the Bayesian probability equation used is Where l is the location number of the potential source, P(M l ) is the source position M l Since each location may be the real source release location without knowing the prior information, P(M l ) is equal to 1 / m, where m is the number of potential sources. (IV) Description of the accompanying drawings

[0029] Figure 1 Flowchart of an atmospheric leakage pollution source identification algorithm combining Gaussian adjoint and regularization

[0030] Figure 2 Experimental case site layout

[0031] Figure 3Zoning diagram of high-probability areas of pollution sources using the adjoint method in the experimental case

[0032] Figure 4 The potential source strength obtained based on the information from monitor 1 in the experimental case

[0033] Figure 5 Posterior probability distribution of potential sources in experimental cases (V) Specific implementation methods

[0034] The specific embodiments of the present invention are further described in detail below in conjunction with the accompanying drawings and examples. The following examples are used to illustrate the present invention but are not intended to limit the scope of the present invention. After reading the contents of the present invention, those skilled in the art can make various changes or modifications to the present invention, and these equivalent forms also fall within the scope defined by the claims attached hereto.

[0035] (1) According to the data provided by the local meteorological station, the wind direction α in the measured area is due north, the wind speed v = 3.98 m / s, and the real-time atmospheric stability w is determined to be level D in Pasquill;

[0036] (2) Test site such as Figure 2 The pollution source coordinates are (0, 0), and the source intensity is constant at 1 kg / s. Three monitors are selected, named S1, S2, and S3. The connecting line of the monitors is perpendicular to the wind speed direction; the coordinates of monitors 1, 2, and 3 are (-9.67, 39.51, 1.6), (0, 39.51, 1.6), and (10, 39.51, 1.6), respectively.

[0037] (3) The wind direction in this case is due north, so the regional coordinates are corrected according to formulas (1-3) and (1-4);

[0038] (4) Substitute the corrected regional coordinates, concentration records, wind speed and direction data into formula (1-1) to obtain the accompanying concentration fields of the three monitors;

[0039] (5) Substitute the adjoint concentration field and S1, S2, S3 data information in (4) into the posterior probability equation of multi-sensor joint Then calculate the probability distribution of regional pollution sources and determine the area where high probability sources exist, an equilateral triangle with coordinates (0, 0) as the center, such as Figure 3 ;

[0040] (6) The high probability area of ​​pollution source is divided into 6 intervals by using Thiessen polygon method, such as Figure 3 , assuming that the pollution source exists at the center of each interval, the response matrix A of the three monitors corresponding to the six interval center points can be obtained through equation (1-2): 11 、A 12 ,...,A21 ,...,A 36 ; The information S1 of monitor 1 and the response matrix A 11 、A 21 ,...,A 16 Substituting regularization, we can obtain the inverse source strength q1, q2, ..., q6 of the six potential pollution sources corresponding to monitor 1, as Figure 4 ;

[0041] (7) The source intensity vectors q1, q2, ..., q6 and the response matrix A 21 、A 22 ,...,A 31 ,...,A 36 Substituting into the equation C = M(A)q, we can get the simulated concentrations of the six potential sources corresponding to monitors 2 and 3. Then Substitute into the equation Calculate the likelihood function of multiple monitors and finally substitute it back into the Bayesian probability equation The posterior probability distribution of the leakage source can be obtained, such as Figure 5 The maximum probable source locations are (0, 0) and (0, 2), both with a maximum posterior probability of 100% and source intensities of 0.967 and 1.733, respectively. The former has the same location as the pollution source, with a source intensity deviation of 0.023 kg / s; the latter has a location deviation of 2 meters and a source intensity deviation of 0.733 kg / s.

Claims

1. An atmospheric leakage pollution source identification algorithm combining Gaussian adjoint and regularization, characterized by: The following steps are involved: Step (1): Based on the data provided by the local meteorological station, obtain the wind direction α and wind speed v and the real-time atmospheric stability w of the measured area; Step (2): When the pollutant concentration information detected by the pollutant monitor reaches a stable state, select N monitors to extract data, and try to make the connection line of the N sensors perpendicular to the wind speed direction; then, read the historical pollutant concentration information and location coordinates detected by the N monitors, and record the data as S1, S2, ..., S n ; Step (3): Apply adjoint theory to the outdoor concentration control equation, and derive the adjoint form of the Gaussian plume and puff equation through adjoint flow field, variable conversion and Laplace transform: C * is the accompanying concentration of pollutants, kg / m 3 ; Q is the release intensity of accompanying pollution sources, kg / s; Q S is the total mass released by the pollution source, kg; σ x ,σ y ,σ z are the diffusion coefficients in the x, y, and z directions, respectively, which are related to atmospheric stability and downwind distance, m; x is the diffusion distance in the downwind direction, m; y, z are the distances perpendicular to the downwind direction, m; u is the wind speed in the flow field, m / s; Step (4): As shown in step (3), the x used in the Gaussian adjoint equation is the diffusion distance relative to the monitor in the downwind direction of the adjoint wind field, and y and z are the distances perpendicular to the downwind direction. In actual situations, the spatial information provided is a Cartesian coordinate system with the north direction as the y-axis and the monitor as the origin. Therefore, the coordinate information needs to be corrected according to the wind direction α: x=d*coS(α2), y=d*sin(α2) (1-4) α0 is the angle of wind speed relative to the north direction; α1 is the angle of the line between the pollution source and the monitor relative to the north direction; α2 is the angle between the line between the pollution source and the monitor and the wind speed direction; x0, y0 are the Cartesian coordinates of the monitor; Step (5): Integrate the information of the N pollutant monitors in step (2) through the posterior probability equation of multi-sensor combination, calculate the posterior probability of the area to be measured, and take the area greater than the third quartile as the high probability area; Step (6): Divide the high-probability area determined in step (5) into M intervals using the Thiessen polygon method. Assume that the pollution source exists at the center of each interval. The response matrix A of each of the N monitors corresponding to the center point of the M intervals can be obtained by equation (1-2). 11 、A 12 ,...,A 21 ,...,A nm ; Using the information S1 of monitor 1 and the response matrix A 11 、A 12 ,...,A 1m , the inverse source strength q1, q2, ..., q of M potential pollution sources corresponding to monitor 1 can be obtained by regularization method. m ; Step (7): Inversely calculate the source strength q1, q2, ..., q m and the response matrix A 21 、A 22 ,...,A 31 ,...,A nm Substituting into the equation C = Aq, we can get the simulated concentration of M potential pollution sources corresponding to N-1 monitors. Then Substitute into the equation Calculate the likelihood function of multiple monitors and finally substitute it back into the Bayesian equation to obtain the maximum probability source location of the leakage source, where: is the concentration detected by monitor i, is the potential source location M l The predicted concentration for monitor i, σ 2 is the variance of the sensor.

2. The algorithm according to claim 1, wherein: In step (5), the posterior probability equation of the multi-sensor combination is: Among them, x i , u, α are the position and monitoring concentration of the i-th monitor, the dominant wind speed and direction in the region, M0 is the release amount of the pollution source; C * (x; x i , u, α) is the accompanying probability of the pollution source location calculated by the i-th monitor through equation (1-1); Under the conditions that the release amount of the pollution source is M0 and the source location is x, the concentration of the monitor is The posterior probability of 3. The algorithm according to claim 1, wherein: In step (7), the Bayesian probability equation used is Where l is the location number of the potential source, P(M l ) is the source position M l The prior probability of P(M l ) is equal to 1 / m, where m is the number of potential sources.

Images