A drifting trajectory model of dust particulate matter and its modeling method

Abstract

The present invention discloses a drifting trajectory model of dust particulate matter and its modeling method, including the following steps: 1) Data preprocessing stage: preprocess data; 2) Model pre-training stage: divide grids according to the geographical locations of site devices; 3) Model pre-training stage: perform model pre-training through simulated data and traditional physical methods; 4) Transfer learning fine-tuning stage: use the trained model to predict the real drifting trajectory of dust particulate matter. The present invention introduces the concept of physical field, replaces the calculation of velocity field with a network of encoder-decoder structure, and optimizes the calculation of convection term, diffusion term, and volume force term at the same time. Experiments show that the dust trajectory prediction model proposed by the present invention has a great improvement in calculation speed and calculation stability compared with traditional methods, and has better generality and interpretability compared with neural network methods.

Description

Technical Field

[0001] The present invention relates to the technical field of dust trajectory prediction, and particularly relates to a dust particulate matter drift trajectory model and a modeling method thereof. Background Art

[0002] With the continuous acceleration of the rate of urban and social development, more and more people have gradually begun to pay attention to environmental and living problems.

[0003] Most of the dust trajectory prediction models are based on traditional physical models or neural networks. Most traditional methods are based on experimental field research, numerical simulation research, meteorological field analysis research, and data mining methods. These traditional physical models have cumbersome calculation steps, a very large amount of calculation, and tend to be in an unstable state; methods based on neural networks, such as convolutional neural networks and recurrent neural networks, have their specific application fields and also lack interpretability.

[0004] The present invention introduces the concept of a field, replaces the calculation of the velocity field in traditional methods with a network of encoder-decoder structures, and simultaneously optimizes the calculation of the convection term, diffusion term, and body force term. Experiments show that the dust trajectory prediction model proposed by the present invention has greatly improved the accuracy and calculation speed compared with traditional methods, while improving the calculation stability, and has better generality and interpretability than methods based on neural networks. Summary of the Invention

[0005] The present invention mainly aims at the problems of heavy calculation tasks, lack of accuracy, and difficulty in calculating the dust movement velocity field in traditional methods, improves relevant traditional physical algorithms and combines them with neural network methods, and proposes a dust particulate matter drift trajectory model and a modeling method thereof.

[0006] The technical solution of the present invention is as follows:

[0007] The dust data used in the present invention includes the following information: the location information of the equipment site, equipment number, equipment longitude, latitude, equipment data status, temperature, humidity, wind speed, wind direction, PM2.5 concentration, and pressure.

[0008] A dust particulate matter drift trajectory model and a modeling method thereof, characterized by including the following steps:

[0009] 1) Data preprocessing stage: Preprocess the data;

[0010] 2) Model pre-training stage: Divide the grid according to the geographical location of the site equipment;

[0011] 3) Model pre-training stage: Perform model pre-training through simulated data and traditional physical methods;

[0012] 4) Migration learning fine-tuning stage: Use the trained model to predict the actual drifting trajectory of dust particles.

[0013] Further, 1) Data preprocessing: To improve the accuracy of model prediction, the following data preprocessing work needs to be carried out:

[0014] 1.1) Data deduplication: Remove duplicate data uploaded during device collection through a joined table query;

[0015] 1.2) Wind speed processing: Taking the geodetic longitude and latitude as the coordinate axes, calculate the longitudinal component u and the transverse component v of the wind speed according to the wind speed and wind direction parameters in the dust data;

[0016] 1.2.1) Calculate the angle θ between the wind direction and the geodetic plane coordinate axes;

[0017] 1.2.2) Calculate the transverse component u and the longitudinal component v of the wind speed V:

[0018] Transverse component u = Vcos(θ); Longitudinal component v = Vsin(θ);

[0019] 1.3) Geodetic longitude and latitude processing: To reduce the influence of site geodetic longitude and latitude errors and improve the accuracy of prediction results, the present invention uses the Gaussian projection algorithm to transform the ordinary geodetic coordinates (L, B) into Gaussian coordinates (x, y);

[0020] 1.3.1) Calculate the central meridian longitude, zone number and longitude difference (distance from the point to the central meridian) according to the longitude and latitude. The 3° central meridian is selected in this model, and its calculation process is as follows:

[0021] 1.3.1.1) Calculate the 3° zone n 3 : n 3 = int(L / 3 + 0.5);

[0022] 1.3.1.2) Calculate the central meridian L of the 3° zone n 3 : L 0 : L 0 = 3n 3 ;

[0023] 1.3.1.3) Calculate the longitude difference l (distance between the site and the central meridian): l = L - L 0 .

[0024] 1.3.2) Calculate the meridian arc length X according to the selection of the ellipsoid. The Krasovsky ellipsoid is selected in this model for calculation (it can be applied to the Beijing 54 coordinate system):

[0025] X = 111134.861B - 16036.480sin2B + 16.828sin4B - 0.022sin6B;

[0026] 1.3.3) Calculate the values of each parameter in the Gauss coordinate forward calculation formula, where B is the geodetic latitude, t is the tangent value of the geodetic latitude, e is the elliptic eccentricity, N is the radius of curvature of the prime vertical, and a is the long radius parameter of the ellipsoid:

[0027]

[0028] 1.3.4) Substitute the parameters calculated above into the Gauss coordinate forward calculation formula to obtain the coordinates (x, y) in the Gauss coordinate system:

[0029] x = X + NsinBcosBl 2 + NsinBcos 3 B(5 - t 2 - 9e 2 cos 2 B)l 4

[0030] y = NcosBl + Ncos 3 B(1 - t 2 + e 2 cos 2 B)l 3 + Ncos 5 B(5 - 18t 2 + t 4 )l 5

[0031] Furthermore, 2) Grid processing:

[0032] In the dust movement, the mutual influence between adjacent stations in space is very close. Therefore, the present invention no longer studies the dust law of a single station, but starts from an overall perspective and studies the distribution and variation law of dust in space through the velocity field, concentration field, temperature field, pressure field, and humidity field. In order to describe the physical field, it is necessary to further divide the grid according to the geographical location of the station equipment.

[0033] 2.1) Cluster analysis: In order to improve the prediction accuracy, before dividing the grid, the present invention performs DBScan cluster analysis according to the geographical location of the station equipment to obtain the cluster of station equipment in Gauss coordinates;

[0034] 2.2) Grid division: After obtaining a better cluster, divide a uniform square grid according to the geographical location of the station equipment within the cluster;

[0035] 2.3) Data interpolation: Linear interpolation is performed based on the distance between device sites and the dust emission data values, and finally the dust emission data is saved in the grid in the form of velocity field, concentration field, temperature field, pressure field, and humidity field.

[0036] Further, 3) Model pre-training:

[0037] Currently, the two mainstream methods applied to predict the drift trajectory of dust particles are the forward calculation method based on traditional physical methods and the neural network method based on numerical simulation. However, both of these methods have certain drawbacks: The calculation of traditional physical methods is too cumbersome and the steps are complicated; the neural network method lacks interpretability and it is generally difficult to obtain a real velocity field as a training label.

[0038] Based on this, the present invention provides a drift trajectory model of dust particles and its modeling method.

[0039] 3.1) Generate simulation data: Simulate and generate dust emission simulation data at time T 0 in the divided grid, including the velocity field, concentration field, temperature field, pressure field, and humidity field of the dust emission;

[0040] 3.2) Forward calculation based on traditional physical methods: Based on the simulation data, forward calculation is performed by updating the calculation of the body force term, convective term, and diffusion term according to traditional physical methods (mainly based on the N-S Navier-Stokes equations) to obtain the augmented dataset of the dust emission physical field at subsequent time nodes T 1 , T 2 , T 3 , …, T n , T n+1 :

[0041]

[0042] 3.2.1) N-S equations: The Navier-Stokes equations are fluid motion equations that describe the conservation of fluid momentum. The specific formula is as follows:

[0043]

[0044] Among them, is the material derivative of the fluid with respect to the velocity field. The Euler perspective mainly focuses on the grid points in space, so the material derivative is equal to the sum of the change rate of the fluid velocity at the grid point with respect to time and the change rate u·▽u of the fluid under the action of the velocity field. In terms of other parameters: ρ is the fluid density; u is the velocity vector; p is the pressure; f is the external force per unit volume of the fluid, and μ is the fluid viscosity.

[0045]

[0046] The above equation is the update process of the fluid state, and W 0 represents the initial state of the particle. After calculation by the body force term (f), the state is W 1 , and after the convective term is calculated, the state W 2 is obtained. Finally, after the diffusion term is calculated, the state W 3 is obtained.

[0047] 3.2.2) Body force term: The calculation formula of the body force term is f = ρg (where g is the acceleration due to gravity and ρ is the fluid density).

[0048] 3.2.3) Convective term: In fluid simulation, the convective term is automatically executed. The essence of convection is that under the action of different velocity fields, fluid micro-elements are transferred between adjacent grid points and the physical quantity values at the grid points are updated. Given a grid point position (x, y), the physical quantity field value at this point is The transverse component of the velocity field corresponding to the particle at the grid point is u, and the longitudinal component is v. The steps to solve the physical quantity field value at the grid point after one time step Δt are as follows:

[0049] 3.2.3.1) Backward tracking for one time step to obtain the grid position of the particle (x - uΔt, y - vΔt);

[0050] 3.2.3.2) Calculate the average change of the particle physical quantity field according to the physical field grid to construct a Catmull-Rom spline interpolation function

[0051] 3.2.3.3) Finally The calculation formula of is:

[0052] 3.2.4) Diffusion term: The calculation of the diffusion term is another very important part in fluid simulation. The calculation formula of the diffusion term of the fluid is where a v is the acceleration of the particle under the action of fluid viscous force, and μ is the viscosity coefficient of the fluid. Assume that at the (n + 1)-th moment, the physical quantity at the grid coordinates (x, y) is known, and now it is necessary to solve the physical quantity

[0053] at the n-th moment.

[0054] 3.2.4.1) Establish the backward Euler equation where Δt is the time step;

[0055] 3.2.4.2) Calculate the Laplace operator using the central difference method:

[0056]

[0057] where is the physical quantities around the grid, Δx is the grid width, Δy is the grid height. In this patent, square grids are used, so Δx = Δy;

[0058] 3.2.4.3) Let Factor out the constant term to get:

[0059]

[0060] 3.2.4.4) Transform it into matrix form and use the conjugate gradient method to solve this large-scale positive definite sparse matrix.

[0061]

[0062] (3.2.5) Update and calculate to obtain the dust physical field data at subsequent time nodes T 1 , T 2 , …, T n , T n+1 moment.

[0063] Furthermore, pre-train with the augmented dataset . (3.3) Pre-training: Use the dust physical information at time nodes T 1 , T 2 , …, T n , T n+1 calculated by the traditional physical method in section (3.2) as the dataset to input into the neural network model U for pre-training:

[0064] (3.3.1) Model selection: In the study of dust data in the present invention, the concepts of grids and fields are introduced. The dataset D is essentially a set of spatial distribution data expressing the variation of various physical variable fields over time, so it can be analogized to the pixel distribution data of an "image". Therefore, the U-net model, which performs excellently in the field of image segmentation prediction, is more suitable for dust physical field data compared to other networks. Thus, the present invention selects the U-net model as the pre-training model.

[0065] ​(3.3.2) Training details: The input channel size is 32×32×6×2, where 32×32 is the grid size, and 6×2 represents six physical field data (the transverse component of the velocity field, the longitudinal component of the velocity field, the concentration field, the temperature field, the pressure field, and the humidity field) and two adjacent time fields.

[0066] The output channel is 32×32×3×1, where 32×32 is the grid size, and 3×1 represents three physical fields predicted and generated (the transverse component of the velocity field, the longitudinal component of the velocity field, and the concentration field).

[0067] and one time field.

[0068] Furthermore, 4) Transfer learning and fine-tuning:

[0069] Input the real dataset into the pre-trained model U, and adjust the number of iterations (epoch),

[0070] the training batch size (batch-size), and the loss function (loss) to obtain the final predicted trajectory of the dust particulate matter drift.

[0071] (4.1) Model transfer: Retain all the parameters in the pre-trained U-Net model and use it as the training model for the real dataset;

[0072] (4.2) Data training: The input channel size is 32×32×4×2, where 32×32 is the grid size,

[0073] 4×2 represents four real physical field data (the concentration field, the temperature field, the pressure field, and the humidity field) and two adjacent time fields. The output channel is 32×32×3×1, where 32×32 is the grid size, and 3×1 represents three physical fields predicted and generated (the transverse component of the velocity field, the longitudinal component of the velocity field, and the concentration field) and one time field;

[0074] (4.3) Model fine-tuning: Adjust the value of the parameter epoch to 240, adjust the value of the training batch parameter batchsize to 8, set L2 as the loss function, and train to obtain the final dust drift trajectory. Among them, the L2 loss function is the least square error function. The purpose of using LSE is to minimize the sum of the squares D of the differences between the true value y i and the training result f(x i ), and its function formula is as follows: L2

[0075]

[0076] The beneficial effects of the present invention are as follows: The model of the present invention is superior to the traditional physical method in terms of accuracy, calculation time, and stability. Moreover, it has better generality and stronger interpretability compared to the numerical estimation method using only neural networks. Therefore, this model can be more widely applied to most scenarios. BRIEF DESCRIPTION OF THE DRAWINGS

[0077] Figure 1 is the overall architecture diagram of the present invention;

[0078] Figure 2 is the schematic diagram of Gaussian coordinates of the present invention;

[0079] Figure 3 is the grid schematic diagram of the present invention; Figure 3 In, the white color scale of each small square represents the size of the concentration field. The higher the white color scale, the higher the dust concentration in the grid;

[0080] Figure 4 is the schematic diagram of the semi-Lagrangian method of the present invention;

[0081] Figure 5 is the dust velocity field calculated by the traditional method of the present invention;

[0082] Figure 6 is the dust concentration field calculated by the traditional method of the present invention; Figure 6 consists of 32×32 small squares. The white color scale of each small square represents the size of the concentration field. The higher the white color scale, the higher the dust concentration in the grid;

[0083] Figure 7 is the schematic diagram of the U-net network structure used in the present invention. DETAILED DESCRIPTION OF THE INVENTION

[0084] The following further describes the present invention in conjunction with the accompanying drawings of the specification.

[0085] Refer to Figure 1-7 , a dust particulate matter drift trajectory model and its modeling method, the specific steps are as follows:

[0086] 1) Data preprocessing:

[0087] (1.1) Duplicate data removal: Perform a joined table query to screen out duplicate dust data uploaded by the device;

[0088] (1.2) Wind speed processing:

[0089] (1.2.1) Calculate the angle θ between the wind direction and the earth's plane coordinate axis;

[0090] (1.2.2) Calculate the transverse component u and the longitudinal component v of the wind speed V:

[0091] The horizontal component u = Vcos(θ); the vertical component v = Vsin(θ).

[0092] (1.3) Perform the Gauss projection on the site coordinates: The goal of this step is to convert (L, B) in longitude and latitude coordinates into (x, y) in Gauss coordinates. The specific steps can be referred to in the appendix Figure 2 ;

[0093] (1.3.1) Calculate the central meridian longitude, zone number, and longitude difference (distance from the point to the central meridian) based on the longitude and latitude. In this model, the 3° central meridian is selected, and its calculation process is as follows:

[0094] (1) Calculate the 3° zone n 3 : n 3 = int(L / 3 + 0.5);

[0095] (2) Calculate the central meridian L 3 of the 3° zone n 0 : L 0 = 3n 3 ;

[0096] (3) Calculate the longitude difference l (spacing between the site and the central meridian): l = L - L 0 ;

[0097] (1.3.2) Calculate the meridian arc length X according to the selection of the ellipsoid. In this model, the Krasovsky ellipsoid is selected for calculation (it can be applied to the Beijing 54 coordinate system):

[0098] X = 111134.861B - 16036.480sin2B + 16.828sin4B - 0.022sin6B;

[0099] (1.3.3) Calculate the values of each parameter in the Gauss coordinate direct calculation formula. Among them, B is the geodetic latitude, t is the tangent value of the geodetic latitude, e is the elliptical eccentricity, N is the radius of curvature of the prime vertical, and a is the long radius parameter of the ellipsoid:

[0100]

[0101] (1.3.4) Substitute the parameters calculated above into the Gauss coordinate direct calculation formula to obtain the coordinates (x, y) in the Gauss coordinate system:

[0102] x = X + NsinBcosBl 2 + NsinBcos 3 B(5 - t 2 - 9e 2 cos 2 B)l 4

[0103] y = N cos Bl + N cos 3 B(1 - t 2 + e 2 cos 2 B)l 3 + N cos 5 B(5 - 18t 2 + t 4 )l 5

[0104] 2) Grid processing:

[0105] (2.1) Clustering analysis: The clustering clusters C are obtained through the DBscan clustering algorithm, and the steps are as follows:

[0106] (1) Mark all objects as unvisited; randomly select an unvisited object P;

[0107] (2) If there are at least MinPts objects in the ε neighborhood of P, create a new cluster C and add P to the cluster C;

[0108] (3) Let N be the set of objects in the ε neighborhood of P. If there are at least MinPts objects in the ε neighborhood of P, add these objects to N; if P is not yet a member of any cluster, add P to the cluster C.

[0109] (2.2) Linear interpolation to obtain a uniform grid (the grid schematic diagram can be referred to in the appendix Figure 3 ):

[0110] (2.2.1) According to the Delaunay triangulation algorithm, first find 3 points around the interpolation point to form a triangle list triangles (the pseudocode is as follows):

[0111]

[0112] (2.2.2) Call the griddata function to perform data interpolation to obtain the dust emission speed field, concentration field, temperature field, pressure field, and humidity field represented by a uniform grid.

[0113] 3) Model pre-training:

[0114] (3.1) Simulate initial data: Use the dust emission data at 23:00 on the evening of July 20, 2021 as the simulated initial data at time T 0 moment.

[0115] (3.2) Forward calculation by physical traditional methods: The algorithm mainly relied on by the present invention for forward calculation by traditional methods is the Navier-Stokes equations, abbreviated as the N-S equations.

[0116] (3.2.1) Body force term: Update the body force term f = ρg (where g is the acceleration due to gravity and ρ is the density).

[0117] (3.2.2) Convection term (the convection process can refer to the appendix Figure 4 ): Given a grid point at position (x, y), the physical quantity field value at this point is The transverse component of the velocity field corresponding to the particle at the grid point is u, and the longitudinal component is v. The steps to solve the physical quantity field value at the grid point after one time step Δt are as follows: are:

[0118] (1) Backtrack one time step to obtain the grid position of the particle (x - uΔt, y - vΔt);

[0119] (2) Construct a Catmull-Rom spline interpolation function based on the physical field grid to calculate the average change in the physical quantity field of the particle

[0120] (3) Finally The calculation formula is:

[0121] (3.2.3) Diffusion term: The main purpose of this term is to solve the viscosity term of the fluid where a v is the acceleration of the particle under the action of the fluid viscous force, and μ is the viscosity coefficient of the fluid. Assume that at the (n + 1)-th moment, the physical quantity at the grid coordinates (x, y) is known, and now it is necessary to solve the physical quantity at the n-th moment The specific calculation steps are as follows:

[0122] (1) Establish the backward Euler equation where Δt is the time step;

[0123] (2) Use the central difference method to calculate the Laplace operator:

[0124]

[0125] where is the physical quantities around the grid, Δx is the grid width, and Δy is the grid height. This patent uses a square grid, so Δx = Δy;

[0126] (3) Let Factor out the constant term to get:

[0127]

[0128] (4) Transform it into matrix form and use the conjugate gradient method to solve this large-scale positive definite sparse matrix.

[0129]

[0130] (3.2.4) According to the body force term, convective term, and diffusion term, the subsequent time nodes T 1 , T 2 , T 3 , …, T n , T n+1 are obtained by the traditional method, and the dust physical field data at these moments

[0131] Furthermore, from the augmented data set

[0132] a pre-trained model U is trained. The velocity field and concentration field in the augmented data set can be referred to in Appendix Figure 5 and Appendix Figure 6 .

[0133] (3.3) Pre-training process:

[0134] (3.3.1) Network selection:

[0135] Based on the characteristics of dust drift movement (spatially, the dust data at adjacent stations has strong correlation), the concept of grid and field is introduced in the study of dust data in the present invention. The data set D is essentially a set of spatial distribution data expressing the variation of various physical variable fields over time, so it can be analogized to the pixel distribution data of an "image". Therefore, the U-net model, which performs excellently in the field of image segmentation prediction, is more adaptable to the dust physical field data in terms of the overall structure compared to other networks. Thus, the U-net model is selected as the pre-trained model in the present invention (the specific U-net model structure diagram can be referred to in Appendix Figure 7 ).

[0136] (3.3.2) Training details:

[0137] (1) The U-Net model is selected for pre-training. The value of the parameter epoch is set to 240, the value of the parameter batchsize is set to 8, the size of the data set is 2000 groups, and the ratio of the training set, test set, and validation set is 12:3:1. The regression index selects the L2 loss function.

[0138] (2) The input channel size is 32×32×6×2: where 32×32 is the grid size, and 6×2 represents six physical field data (the transverse component of the velocity field, the longitudinal component of the velocity field, the concentration field, the temperature field, the pressure field, and the humidity field) and two adjacent time fields. The output channel is 32×32×3×1: where 32×32 is the grid size, and 3×1 represents three physical fields (the transverse component of the velocity field, the longitudinal component of the velocity field, and the concentration field) and one time field predicted and generated.

[0139] (4) Model migration and fine-tuning:

[0140] (4.1) Model migration: Retain all the parameters in the pre-trained U-Net model and use it as the training model for the real dataset;

[0141] (4.2) Data training: The input channel size is 32×32×4×2, where 32×32 is the grid size, and 4×2 represents four real physical field data (the concentration field, the temperature field, the pressure field, and the humidity field)

[0142] and two adjacent time fields. The output channel is 32×32×3×1, where 32×32 is the grid size, and 3×1 represents three physical fields (the transverse component of the velocity field, the longitudinal component of the velocity field, and the concentration field) and one time field predicted and generated;

[0143] (4.3) Model fine-tuning: Adjust the value of the parameter epoch to 240, adjust the value of the training batch parameter batchsize to 8, set L2 as the loss function, and train to obtain the final dust drift trajectory. Among them, the L2 loss function is the least square error function. The purpose of using LSE is to minimize the square sum D of the difference between the true value y i and the training result f(x i ), and its function formula is as follows: L2 That is,

[0144]

[0145] In view of the problems of excessive computational complexity, poor accuracy and stability in traditional methods, and the lack of interpretability in neural network methods, the present invention combines an encoder-decoder neural network and an optimized traditional method, and proposes a dust particle drift trajectory model and its modeling method.

[0146] The content described in the embodiments of this specification is only an enumeration of the implementation forms of the inventive concept. The protection scope of the present invention should not be regarded as limited to the specific forms stated in the embodiments. The protection scope of the present invention also extends to equivalent technical means that those skilled in the art can think of according to the inventive concept.

Claims

1. A method for modeling the drifting trajectory of dust particles, characterized by including the following steps: 1) Data preprocessing stage: Preprocess the data; The data preprocessing stage includes the following steps: 1.1) Data duplicate checking stage: Through joined table query, filter out the dust data repeatedly uploaded by site equipment; 1.2) Wind speed processing stage: According to the wind speed and wind direction parameters in the dust data, calculate the lateral component and longitudinal component of the wind speed. The specific implementation steps are as follows: 1.2.1) Calculate the angle θ between the wind direction and the axes of the earth's plane coordinate system; 1.2.2) Calculate the wind speed V: u = V cos(θ), v = V sin(θ); where u represents the transverse component and v represents the longitudinal component; 1.3) Geographical coordinate Gauss projection stage: Convert the longitude and latitude coordinates (L, B) of the station into Gauss coordinates (x, y) through Gauss projection calculation, where L is the geodetic longitude and B is the geodetic latitude; The specific implementation steps are as follows: 1.3.1) Calculate the central meridian longitude, zone number, and longitude difference based on the longitude and latitude. The 3° zone central meridian is selected for the model. The steps are as follows: 1.3.1.1) Calculate the zone number n of the 3° zone 3 : n 3 = int(L / 3 + 0.5); 1.3.1.2) Calculate the 3° zone number n 3 of the central meridian L 0 : L 0 = 3n 3 ; 1.3.1.3) Calculate the longitude difference l: l = L - L 0 ; 1.3.2) Calculate the meridian arc length X according to the selection of the ellipsoid; The Krasovsky ellipsoid is selected for the model calculation: X = 111134.861B - 16036.480sin2B + 16.828sin4B - 0.022sin6B; 1.3.3) Calculate the values of each parameter in the Gauss coordinate direct calculation formula, where, B is the geodetic latitude, t is the tangent value of the geodetic latitude, e is the elliptical eccentricity, N is the radius of curvature of the prime vertical, and a is the long radius parameter of the ellipsoid: 1.3.4) Use the Gauss coordinate direct calculation formula to solve the coordinates (x, y) in the Gauss coordinate system: x = X + N sin B cos Bl 2 + N sin B cos 3 B(5 - t 2 - 9e 2 cos 2 B)l 4 y = N cos Bl + N cos 3 B(1 - t 2 + e 2 cos 2 B)l 3 + N cos 5 B(5 - 18t 2 + t 4 )l 5 ; 2) Network processing stage: Divide the grid according to the geographical location of the station equipment; The grid processing stage includes the following steps: 2.1) Clustering analysis stage: Before dividing the grid, perform DBScan clustering analysis according to the geographical location of the station equipment to obtain the clustering clusters of the station equipment in Gauss coordinates; 2.2) Grid division stage: After obtaining a better clustering cluster in step 2.1), divide a uniform square grid according to the geographical location of the station equipment within the clustering cluster; 2.3) Data interpolation stage: Perform linear interpolation according to the distance between equipment stations and the dust emission data value, and finally save the dust emission data in the grid in the form of a velocity field, concentration field, temperature field, pressure field, and humidity field; 3) Model pre-training stage: Perform model pre-training through simulated data and traditional physical methods; The model pre-training stage includes the following steps: 3.1) Simulation data stage: Simulate and generate the dust emission simulation data at time T in the divided grid, including the velocity field, concentration field, temperature field, pressure field, and humidity field of the dust emission; 0 ​ 3.2) Forward calculation stage of traditional physical methods: Based on the simulated data, according to the traditional physical methods, that is, updating the calculation of the volume force term f, the convection term, and the diffusion term according to the N-S Navier-Stokes equation and performing forward prediction to obtain the subsequent time nodes T 1 , T 2 , T 3 , …, T n , T n+1 of the dust physical field dataset The specific implementation steps are as follows: 3.2.1) Calculate the body force term f: f = ρg, where g is the acceleration due to gravity and ρ is the fluid density; 3.2.2) Convection term calculation: Given a grid point at position (x, y) with the physical quantity field value at this point being The transverse component of the velocity field corresponding to the particle at the grid point is u, and the longitudinal component is v. Solve for the physical quantity field value at the grid point after one time step Δt The relevant steps are as follows: 3.2.2.1) Trace backward for 1 time step to obtain the grid position of the particle (x - uΔt, y - vΔt); 3.2.2.2) Calculate the mean change of the particle physical quantity field by constructing a Catmull-Rom spline interpolation function based on the physical field grid 3.2.2.3) Final The calculation formula is: 3.2.3) Diffusion term calculation: The main purpose of this step is to solve the viscosity term of the fluid where a v is the acceleration of the particle under the action of the fluid viscous force, and μ is the viscosity coefficient of the fluid; assuming that at the (n + 1)-th moment, the physical quantity at the grid coordinates (x, y) is known, and now it is necessary to solve the physical quantity at the n-th moment The specific calculation steps are as follows: 3.2.3.1) Establish the backward Euler equation where Δt is the time step; 3.2.3.2) Calculate the Laplace operator using the central difference method: Among them is the physical quantity around the grid, Δx is the grid width, and Δy is the grid height; 3.2.3.3) Let Factor out the constant term to obtain: 3.2.3.4) Convert it into matrix form and use the conjugate gradient method to solve the large-scale positive definite sparse matrix: 3.3) Dataset generation stage: Update and calculate the subsequent time nodes T according to the body force term, convection term, and diffusion term 1 , T 2 , …, T n , T n+1 The dust physical field data at time Furthermore, pre-training is performed using the dataset ; 3.4) Pre-training stage: Use the dust physical information at the time nodes of T 1 , T 2 , T 3 , …, T n , T n+1 calculated by the traditional physical method in step 3.3) as the dataset for pre-training the neural network model U: 1 T 2 T 3 T n T n+1 dust physical information as the dataset for pre-training the neural network model U: 3.4.1) Model selection: Select U-Net as the training network; 3.4.2) Training details: The input channel size is 32×32×6×2, where 32×32 is the grid size, and 6×2 represents 6 physical field data and 2 adjacent time fields. The 6 physical field data are the transverse component of the velocity field, the longitudinal component of the velocity field, the concentration field, the temperature field, the pressure field, and the humidity field respectively; the output channel is 32×32×3×1, where 32×32 is the grid size, and 3×1 represents 3 physical fields and 1 time field generated by prediction. The 3 physical fields are the transverse component of the velocity field, the longitudinal component of the velocity field, and the concentration field respectively; 4) Transfer learning fine-tuning stage: Use the trained model to predict the real drifting trajectory of dust particles; The transfer learning fine-tuning stage includes the following steps: 4.1) Model transfer stage: Retain all the parameters in the pre-trained U-Net model in step 3.4) and use it as the training model for the real dataset; 4.2) Data training stage: The input channel size is 32×32×4×2, where 32×32 is the grid size, and 4×2 represents 4 real physical field data and 2 adjacent time fields. The 4 real physical field data are the concentration field, the temperature field, the pressure field, and the humidity field respectively; the output channel is 32×32×3×1, where 32×32 is the grid size, and 3×1 represents 3 physical fields and 1 time field generated by prediction. The 3 physical fields are the transverse component of the velocity field, the longitudinal component of the velocity field, and the concentration field respectively; 4.3) Model fine-tuning stage: Adjust the value of the parameter epoch to 240, set the value of the training batch parameter batchsize to 8, set L2 as the loss function, and finally obtain the prediction result; Among them, The loss function L2 uses the least squares error (LSE) function. The purpose of using LSE is to minimize the true value y i and the squared sum D of the differences between the training result f(x i ). Its function formula is L2 .

Images