Soil groundwater pollution diffusion prediction method and system based on space-time diagram neural network

Abstract

The invention discloses a soil and groundwater pollution diffusion prediction method and system based on a time-space diagram neural network, and relates to the crossing field of environmental information modeling and artificial intelligence, and the prediction method comprises the following steps: S1, obtaining multi-source soil and groundwater monitoring data and external driving data; s2, generating a static hydraulic communication matrix and a dynamic convection dominant matrix, and constructing a multi-relation dynamic graph structure; s3, performing parallel processing on graph convolution operation of the static hydraulic connection matrix and the dynamic convection dominant adjacency matrix, capturing time dynamic characteristics, and obtaining a pollutant concentration prediction result; and S4, outputting structured prediction data and a visual product. By constructing a multi-relation dynamic graph structure fusing a hydrogeological physical mechanism and a dynamic flow field response, distortion of a pollution plume leading edge deduction trajectory caused by misjudgment of a spatial dependency relation can be avoided, so that the prediction precision of a non-uniform anisotropic pollution diffusion path is improved.

Description

TECHNICAL FIELD

[0001] The present application relates to the technical field of environmental information modeling and artificial intelligence, and particularly relates to a soil and groundwater pollution diffusion prediction method and system based on a spatio-temporal graph neural network. BACKGROUND

[0002] With the deepening of the concept of ecological environment protection and sustainable development, the prevention and control and governance of soil and groundwater pollution have become an important research direction in the field of environmental engineering. In recent years, with the rapid development of artificial intelligence technology, prediction methods based on deep learning, especially spatio-temporal graph neural networks that integrate spatial topology and temporal dynamic characteristics, have gradually been introduced into environmental system modeling due to their strong ability to model complex system dynamics. This kind of method abstracts monitoring wells as graph nodes and models the spatial relationship between nodes as graph edges. It uses graph convolution operations to aggregate neighborhood information and combines recurrent neural networks or spatio-temporal attention mechanisms to capture time series features, thereby achieving joint modeling and prediction of system state. In application scenarios such as traffic flow and urban crowd density, such methods have shown good prediction performance.

[0003] Traditional spatio-temporal graph neural networks generally use adjacency matrix definition based on geometric distance when constructing graph structures, such as using the negative exponential function of the distance between nodes as the edge weight, or setting a fixed distance threshold to determine the connection relationship. Its simplified processing has high computational efficiency and certain generalization ability in scenarios where spatial relationships are relatively uniform and dynamic changes are weak. Correspondingly, some improvement schemes attempt to introduce fixed attribute weighting, such as combining land use types or elevation differences to modify connection strength, to enhance the representation ability of the graph structure. On this basis, the model extracts spatial features through multi-layer graph convolution and works with the time encoding module to complete the deduction of the future state of the system.

[0004] However, the above technical solutions gradually show limitations in dealing with the highly complex, heterogeneous and dynamic evolving system of groundwater pollution diffusion. The migration process of underground pollutants does not follow the simple principle of spatial proximity, but is deeply coupled with the hydrogeological conditions and instantaneous flow field state of the underground medium. The migration of pollutants in porous media or fracture networks is essentially a transport process driven by multiple physical and chemical mechanisms such as convection, dispersion and adsorption, among which convection often plays a dominant role, and its direction and strength directly depend on the hydraulic gradient and spatial distribution of permeability of groundwater. Further, the permeability coefficient of the aquifer has significant spatial variability and anisotropy, and the water conductivity of different regions differs greatly, leading to the possibility of rapid migration of pollutants along high-permeability zones while bypassing low-permeability zones, forming a non-uniform and asymmetric diffusion front. In addition, external driving factors such as rainfall infiltration and artificial pumping will dynamically change the groundwater level and flow field direction, making the potential migration path of pollutants evolve over time. Under this background, if the graph structure based on fixed geometric distance or static attributes is still used, it will inevitably lead to misjudgment of the model's spatial dependence relationship. Nodes that are weakly hydrologically connected but have short geometric distances may be mistakenly assigned high connection weights, while nodes that actually have strong convection channels but are far apart may be ignored. The mismatch between the graph structure and the real physical process directly leads to the aggregation of incorrect neighbor information in the message passing process of the model, which distorts the deduced trajectory of the pollution plume front, and the prediction deviation is significantly amplified, especially in sparse monitoring wells or areas with sudden changes in flow field. SUMMARY

[0005] The present application overcomes the shortcomings of the prior art and provides a soil and groundwater pollution diffusion prediction method and system based on a spatio-temporal graph neural network.

[0006] To achieve the above-mentioned purpose, the technical scheme adopted by the present application is as follows: In a first aspect, the present application provides a soil and groundwater pollution diffusion prediction method based on a spatio-temporal graph neural network, comprising the following steps:

[0007] S1, acquiring multi-source soil and groundwater monitoring data and external driving data, and performing spatio-temporal benchmark calibration and preprocessing;

[0008] S2, generating a static hydraulic connectivity matrix based on the spatial coordinates and hydrogeological parameters of the monitoring nodes, and generating a dynamic convection dominant matrix combined with the water head field data of numerical simulation, to construct a multi-relation dynamic graph structure;

[0009] S3, using a spatio-temporal graph neural network, performing graph convolution operations on the static hydraulic connectivity matrix and the dynamic convection dominant adjacency matrix through a multi-relation graph convolution layer, fusing the spatial features of the nodes, and then capturing the time dynamic features through a gated recurrent unit to obtain the pollutant concentration prediction results;

[0010] S4, performing a rolling multi-step prediction, periodically injecting actual observation data to optimize the prediction process, and outputting structured prediction data and visualized products.

[0011] In a preferred embodiment of the present application, in the step S1, the multi-source soil and groundwater monitoring data include pollutant concentration time series and groundwater level observation data; and the external driving data include daily precipitation, potential evapotranspiration and artificial pumping and injection data.

[0012] In a preferred embodiment of the present application, in the step S2, the static hydraulic connection matrix is generated by the following method: based on the permeability coefficient, porosity of the monitoring well node and the Euclidean distance between nodes, the initial weight is calculated, and the normalized static adjacency matrix is obtained after normalization processing.

[0013] In a preferred embodiment of the present application, in the step S2, the dynamic convection dominant matrix is generated by the following method: input the driving external input data into the groundwater flow numerical simulation engine, and solve the water head distribution field; based on the water head value of each monitoring well, the hydraulic gradient between nodes is calculated, the convection weight is determined combined with the equivalent permeability coefficient, and the normalized dynamic adjacency matrix is obtained after normalization processing.

[0014] In a preferred embodiment of the present application, the water head distribution field is obtained by solving the groundwater flow control equation, and the equation considers the permeability coefficient tensor component, specific yield and source and sink term.

[0015] In a preferred embodiment of the present application, in the step S3, the multi-relation graph convolution layer performs the following operations: performing static relation graph convolution and dynamic relation graph convolution on the input node hidden state respectively, and fusing the two branch outputs into a unified node representation through linear transformation after splicing.

[0016] In a preferred embodiment of the present application, in the step S3, the gating recurrent unit receives a node feature sequence of a plurality of consecutive time steps, calculates a current hidden state through an update gate, a reset gate and a candidate hidden state, and finally maps the current hidden state into a pollutant concentration prediction value through a fully connected layer.

[0017] In a preferred embodiment of the present application, in the step S3, the loss function of the spatio-temporal graph neural network includes a mean square error term and a physical consistency regularization term based on the mass conservation principle of pollutants, which is used to strengthen the model to follow the physical law.

[0018] In a preferred embodiment of the present application, in the step S4, the rolling multi-step prediction process injects actual observation data to reinitialize the feature sequence after every preset number of steps, which is used to reduce error accumulation.

[0019] In a second aspect, the present application provides a prediction system of a soil and groundwater pollution diffusion prediction method based on a spatio-temporal graph neural network, comprising:

[0020] A data acquisition and processing module is configured to acquire multi-source soil and groundwater monitoring data and external driving data, and perform preprocessing on the data.

[0021] A graph structure generation module is configured to generate a static hydraulic connectivity matrix and a dynamic convection dominant matrix, and construct a multi-relation dynamic graph structure that integrates the static medium conduction capacity and the dynamic flow field driving direction.

[0022] A spatio-temporal graph neural network prediction model is configured to perform graph convolution operations on the static hydraulic connectivity matrix and the dynamic convection dominant matrix in parallel, integrate node spatial features, and output pollution concentration prediction values at future time points.

[0023] A result output and visualization module is configured to perform a rolling multi-step prediction task, and generate structured data streams and geographic information system visualization products.

[0024] The present application solves the defects in the background art, and has the following beneficial effects:

[0025] (1) The present application provides a soil and groundwater pollution diffusion prediction method and system based on a spatio-temporal graph neural network, which constructs a multi-relation dynamic graph structure that integrates hydrogeological physical mechanisms and dynamic flow field responses. This structure includes static hydraulic connectivity relationships based on medium properties and dynamic convection dominant relationships based on instantaneous flow fields, enabling the graph neural network to accurately represent spatial dependency relationships between nodes and be highly consistent with the physical processes of groundwater pollution migration. This effectively improves the accuracy of spatial modeling and, compared to traditional graph structures based on static geometric distances, can avoid the distortion of pollution plume front trajectory caused by misjudgment of spatial dependency relationships, thereby improving the prediction accuracy of non-uniform anisotropic pollution diffusion paths, especially in areas with sudden flow field changes or sparse monitoring wells.

[0026] (2) In the present application, a multi-relation graph convolution layer is used to process static and dynamic graph convolution operations in parallel, and the outputs of the two branches are fused through linear transformation. This architecture enables the model to simultaneously learn the inherent conduction capacity of the medium and the driving direction of the instantaneous flow field, thereby enhancing the robustness and adaptability of spatial feature extraction. Compared to traditional methods that use a single static graph convolution and cannot effectively respond to flow field changes caused by external driving factors, the present application can improve the generalization ability of the model under complex hydrogeological conditions, thereby reducing the distortion of the pollution plume trajectory.

[0027] (3) In the model in the present application, a physically consistent regularization term based on the mass conservation principle of pollutants is introduced in the training process, a loss function is constructed by verifying the matching degree of node net flux and its own mass change, the regularization term forces the model to learn a dynamic mode that conforms to the physical law of groundwater solute transport, so that the prediction result not only approaches the measured value in value, but also follows the physical and chemical mechanisms such as convection and dispersion, thereby reducing the physically unreasonable prediction output, effectively enhancing the physical interpretability of the model, so that the prediction result is more in line with the actual pollution diffusion process, and provides a more reliable basis for repair decision. BRIEF DESCRIPTION OF DRAWINGS

[0028] In order to more clearly illustrate the technical solutions in the embodiments of the present application or the prior art, the drawings needed to be used in the embodiments or prior art description will be briefly introduced as follows. Obviously, the drawings in the following description are only some embodiments described in the present application, and other drawings can be obtained by those skilled in the art without creative labor.

[0029] Figure 1 is a timing diagram of the soil groundwater pollution diffusion prediction method of embodiment 1 of the present application;

[0030] Figure 2 is a multi-relation dynamic graph structure diagram of embodiment 1 of the present application, which shows the graph vertex set composed of monitoring well nodes, and the double edge structure constructed based on static hydraulic connection relationship and dynamic convection dominant relationship. DETAILED DESCRIPTION

[0031] The technical solutions in the embodiments of the present application will be described clearly and completely in the embodiments of the present application combined with the drawings. Obviously, the described embodiments are only some of the embodiments of the present application, not all. Based on the embodiments in the present application, all other embodiments obtained by those skilled in the art without creative labor are within the scope of protection of the present application.

[0032] In the following description, many specific details are set forth in order to provide a thorough understanding of the present application, but the present application can also be implemented in other ways different from those described herein, therefore, the scope of protection of the present application is not limited by the specific embodiments disclosed below.

[0033] Embodiment 1: as shown in the following table, a soil groundwater pollution diffusion prediction method based on spatio-temporal graph neural network, comprising the following steps: Figure 1

[0034] S1, acquire multi-source soil groundwater monitoring data and external driving data, and perform spatio-temporal benchmark calibration and pretreatment;

[0035] ​S2, generating a static hydraulic connection matrix based on the spatial coordinates of the monitoring nodes and the hydrogeological parameters, combining the head field data of the numerical simulation to generate a dynamic convection dominant matrix, and constructing a multi-relation dynamic graph structure;

[0036] S3, using a spatio-temporal graph neural network, performing graph convolution operations on the static hydraulic connection matrix and the dynamic convection dominant adjacency matrix in parallel through multi-relation graph convolution layers, capturing time dynamic characteristics through a gated recurrent unit after fusing node spatial features, and obtaining a pollutant concentration prediction result;

[0037] S4, performing a rolling multi-step prediction, periodically injecting actual observation data to optimize the prediction process, and outputting structured prediction data and visualization products.

[0038] In some specific embodiments, in the step of S1, multi-source soil and groundwater monitoring data is obtained from a distributed environmental monitoring network, including pollutant concentration time series and groundwater level observation data. The pollutant concentration time series is sampled automatically every 72 hours, and the data fields include sampling timestamp, well location number, detected pollutant type, concentration value, detection method identifier, and confidence label. The groundwater level observation data is continuously collected by a pressure water level meter installed in each monitoring well, with an original sampling interval of every 15 minutes. After sliding average filtering, the mean value every 24 hours is extracted, and the fields include timestamp, well location number, water level elevation, and sensor calibration coefficient.

[0039] Among them, the external driving data includes daily precipitation recorded by the weather station, potential evapotranspiration, and daily pumping / injection volume, operation period and well location coordinates of artificial pumping / injection wells registered by the groundwater management platform. After all data are accessed into the system, unified space-time reference correction is performed: geographic coordinates are converted to UTM projected coordinate system, timestamp is unified to UTC+0, data are aligned to 24 hours as the basic time unit, and missing values are completed by cubic spline interpolation. The preprocessed data is written in column storage format into a time series database with ACID transaction characteristics.

[0040] In some specific embodiments, as shown in Figure 2 In the step of S2, the monitoring nodes constitute the vertex set V of the graph structure in space, defined as where n is the total number of monitoring wells, ranging from 32 to 128. Each node v i The position in the three-dimensional geographic coordinate system is determined by longitude λ i , latitude φ i and elevation z i . The coordinate data is obtained by differential GPS measurement, with a plane positioning accuracy better than ±0.1 m (95% confidence interval) and an elevation accuracy better than ±0.05 m (RMS error).

[0041] Further, each node v i is associated with two types of attribute vectors: static attribute vector Aiv static and dynamic attribute vector Aiv dynamic .

[0042] wherein the static attribute vector Aiv static contains hydrogeological parameters of the location where the well is located, specifically including permeability coefficient K i , porosity n i , aquifer type code (classification code: 1=sand layer, 2=silt layer, 3=gravel layer, 4=fractured bedrock), sediment thickness and median particle size d 50 of rock particles; the dynamic attribute vector Aiv dynamic contains the measured pollutant concentration ci(t), water level hi(t), water temperature Ti(t) and pH value pi(t) of the well at time t at each time step t. The static attribute vector is generated by geological drilling sampling, laboratory determination (according to GB / T 50123-2019 standard) and Kriging spatial interpolation method.

[0043] Preferably, the interpolation process adopts an exponential variation function model, nugget value C o =0.02, base value C=0.18, correlation distance a=80m, and interpolation grid resolution is 10m x 10m. Finally, the attribute values at the locations of each monitoring well are extracted by nearest neighbor search and bound to the corresponding nodes to form the node static feature matrix .

[0044] Further, the dynamic graph structure generation performs a dual construction mechanism of graph edge relationship, generating a multi-relationship adjacency structure containing static hydraulic connectivity relationship and dynamic convection dominant relationship. First, the static hydraulic connectivity adjacency matrix A is calculated wherein element a represents the potential hydraulic conduction capacity between node v i and v j based on medium attributes. The weight is calculated by the following deterministic function:

[0045] ;

[0046] wherein d is the Euclidean distance between node v i and v j , and the calculation formula is:

[0047] ;

[0048] wherein R is the average radius of the earth, taking the value 6371000m. The obtained The sum of each row element is 1 after normalization processing, and the normalization operation is defined as:

[0049] ;

[0050] Normalized matrix As a static graph convolution operation, the fixed adjacency weight remains unchanged throughout the life cycle of the model, and only triggers the recalculation process when the hydrogeological exploration data is updated. The matrix is stored in shared memory in sparse CSR format for direct call by the model inference stage.

[0051] Further, the dynamic graph structure generation module synchronously calculates the dynamic convection dominant adjacency matrix , which is updated according to the latest water head field data before each prediction period begins. The calculation of the matrix element uses the physical simulation-based scheme A, which is implemented as follows: the driving data obtained by the multi-source data acquisition module is input into the three-dimensional unsteady flow numerical simulation engine, which solves the groundwater flow control equation based on the finite difference method:

[0052] ;

[0053] wherein, denotes the partial derivative, denotes the rate of change of water head h with respect to x direction when y, z coordinates and time t are fixed, i.e. the hydraulic gradient in x direction; h(x, y, z, t) is the water head value at position (x, y, z) at time t, Kx, Ky, Kz are the permeability coefficient tensor components in x, y, z directions respectively (provided by the anisotropic geological model), Ss is the specific yield, with a value range of 0.01-0.15, N is the source-sink term, including rainfall infiltration recharge rate and pumping rate of pumping wells.

[0054] It should be noted that the model space is discretized using structured hexahedral grids, with a horizontal grid size of 10m x 10m, 3 layers vertically, the total number of layers matches the aquifer structure, and the total number of grids is about 1.2 x 10 5 . The time is discretized using an implicit backward difference format, with a time step of 1 day. The initial water head field is obtained by inverse distance weighted interpolation from historical observation data, and the boundary conditions are set as follows: the west and east sides are constant head boundaries, with water head values determined by boundary monitoring well observations; the north and south sides are impermeable boundaries, with normal flow rate of zero; the upper boundary is a variable head boundary driven by source-sink terms, with net recharge = rainfall infiltration - evapotranspiration - surface runoff.

[0055] It can be understood that based on the simulated water head field, the water head value hi(t) at each monitoring well position is extracted and used to construct the dynamic convection adjacency matrix. For any node pair (i, j), the modulus and direction of the hydraulic gradient vector are calculated:

[0056] ;

[0057] like Then define the direction of water flow as from i to j, let ,in The equivalent permeability coefficient between i and j is calculated as the harmonic mean:

[0058] ;

[0059] like Then let A matrix consisting of all non-zero elements. Then, row normalization is performed:

[0060] ;

[0061] The resulting normalized matrix At each time step t, the variable adjacency weights of the dynamic graph convolution operation are input into the model.

[0062] In some specific implementations, in step S3, the prediction model of the spatiotemporal graph neural network adopts a hybrid architecture of graph convolution-gated recurrent units, which is a multi-relation graph convolution layer (MR-GCL). This layer performs two graph convolution operations in parallel within a single layer. For the l-th MR-GCL, its input is the hidden state representation of the nodes in the (l-1)-th layer. , where d is the hidden layer dimension, with a value of 128. This layer first performs a static relational graph convolution:

[0063] ;

[0064] in, The learnable parameter weight matrix is ​​initialized using a Xavier uniform distribution (range). ), The ReLU activation function is used. Then, dynamic graph convolution is performed:

[0065] ;

[0066] in, For another independent learnable parameter matrix, the initialization method is the same as above. The outputs of the two branches are concatenated and then fused through a linear transformation:

[0067] ;

[0068] Where || represents the tensor concatenation operation; The weight matrix is ​​initialized using a He normal distribution. The fused node representation... It serves as input to the next layer, or as input to the temporal module in the last layer. The model consists of three stacked MR-GCL layers, with the output dimension of each layer remaining constant at 128.

[0069] In this embodiment, the spatiotemporal graph neural network prediction model is connected to a time series modeling component after spatial feature extraction. This component adopts a two-layer gated recurrent unit (GRU) stacked structure. The first layer GRU receives the node feature sequence of n time steps from time tk to time t extracted by MR-GCL. Where k is the length of the input sequence, with a value of 14. The calculation of the update gate zt, reset gate rt, and candidate hidden state h̃t of the GRU unit follows the standard formula:

[0070] ;

[0071] ;

[0072] ;

[0073] ;

[0074] Where xt represents the input feature at time t; ht-1 represents the hidden state at the previous time step; ⊙ represents the Hadamard product; Wz, Wr, and Wh are learnable weight matrices (dimension 128×256); and bz, br, and bh are bias vectors. The hidden state sequence output by the first-layer GRU is subjected to max pooling to extract global temporal features, which are then input into the second-layer GRU for high-level semantic encoding. Finally, the output ht of the second-layer GRU at time t is mapped to an n-dimensional pollutant concentration prediction value through a fully connected layer. The weight matrix of the fully connected layer has a dimension of 128×n, and the activation function is a linear function to preserve the physical dimensions of concentration.

[0075] It should be noted that the model training process was performed on a training server equipped with four NVIDIA A100 GPUs (80GB of VRAM), using a supervised learning paradigm. The training dataset was constructed from historical observation data, with a sliding window extracting time series segments of length k+1=15. The input was the node feature matrix of the first 14 time steps. (m is the feature dimension of each node, with a value of 6, and the fields are concentration, water level, water temperature, pH, cumulative precipitation, and pumping volume), and the corresponding dynamic graph structure sequence. The label is the actual concentration vector at time step 15. The loss function takes a combined form:

[0076] ;

[0077] Where MSE is the mean squared error term; α is the concentration prediction weighting coefficient, with a value of 0.8; This is a physical consistency regularization term; β takes the value 0.2.

[0078] Furthermore, The calculation is based on the principle of conservation of pollutant mass, and the following operation is performed in each training batch: based on the predicted concentration changes Combined with aquifer porosity n i With effective thickness b i The value is 5.0-15.0m, and the mass increment of each node is calculated. ,in, To control the area, the value is π×(150). 2 m 2 ); based on the dynamic adjacency matrix The flow field direction is defined in the formula, and the mass flux between nodes is calculated. And verify the degree of matching between the net flux of each node and its own quality changes. Defined as the L2 norm of the mass imbalance error at all nodes. This regularization term is jointly optimized with the main loss during backpropagation.

[0079] Preferably, the AdamW algorithm is used for model parameter optimization, with an initial learning rate set to 3×10⁻⁶. -4 The weight decay coefficient is 1×10 -4 The batch size is 32. The training process lasts for 500 epochs, with performance evaluated on the validation set every 10 epochs. The validation set comprises 15% of the total data. An early stopping mechanism is used to prevent overfitting, with a tolerance threshold of no decrease in loss after 5 consecutive validations. Model parameters are stored in a distributed file system (Ceph cluster) in HDF5 format, including network weights, normalized statistics, and configuration metadata, namely version number, training timestamp, and hyperparameter records.

[0080] In some specific implementations, in step S4, the pollution diffusion prediction performs a rolling multi-step prediction task. At time t=T, the latest observation data from k=14 time steps and the corresponding dynamic graph structure are input into the trained model, and the predicted concentration value at time T+1 is output. Then, the system will The feature sequence is injected as new observations, and the dynamic graph structure generation module is invoked to recalculate based on the updated water level prediction or actual observations. This forms a closed-loop feedback mechanism.

[0081] It should be noted that this process is repeated to achieve rolling predictions for the next H=30 time steps. To reduce error accumulation, the system forcibly injects actual observation data (if any) after every 5 prediction steps and reinitializes the feature sequence.

[0082] In this implementation, the prediction results output includes two categories: structured data streams and visualization products. The structured data is published in JSON format via a RESTful API, conforming to the OpenAPI 3.0 specification. Response fields include: prediction timestamp, monitoring well number, pollutant type, predicted concentration value, 95% confidence interval, and data version number. The visualization product is generated by a geographic information system rendering engine, employing a dual-layer overlay mode: the base layer is the orthophoto of the study area and the distribution of monitoring well locations; the overlay layer is a pollutant concentration contour map generated based on an inverse distance weighted interpolation algorithm. The contour interval is dynamically set according to the pollutant toxicity threshold, and high-risk areas (concentrations exceeding the GB / T 14848-2017 Class III standard) are filled with a red gradient. The system supports a time slider control, which can dynamically play the evolution of the pollution plume over 30 days, with a frame rate set to 1 frame / s.

[0083] Example 2: As a preferred embodiment of the present invention, in areas where high-precision hydrogeological parameters are unavailable, a dynamic graph structure generation scheme is employed, using an alternative static attribute construction method. In this case, the static hydraulic connectivity weight... The calculation incorporates a terrain slope correction factor. This factor is extracted from a digital elevation model (DEM) with a spatial resolution of 5m × 5m, and the calculation formula is as follows:

[0084] ;

[0085] Among them, z i Indicates the elevation of monitoring well node i; z j Indicates the elevation of monitoring well node j;

[0086] And define: ;

[0087] Wherein, λ is the topographic sensitivity coefficient, determined to be 2.5 after calibration with historical data. This scheme effectively compensates for the uncertainty caused by the lack of permeability coefficient data in shallow groundwater systems in mountainous areas, and is applicable to elevation changes greater than 50 m / km. 2 The area.

[0088] In another preferred embodiment of the present invention, the spatiotemporal graph neural network prediction model uses a graph attention mechanism instead of fixed-weight convolution. Within a multi-relation framework, the static relation branches employ a self-attention mechanism to calculate dynamic weights.

[0089] ;

[0090] Among them, v i v j The input feature vectors (6 dimensions) for nodes i and j. It is a linear transformation matrix (d=128). This is the attention weight vector (initialized to a standard normal distribution). Dynamic relationship branches maintain fixed weights based on the physical flow field. This ensures the physical correctness of the convection direction. The attention outputs from the two branches are summed before being fed into the subsequent time module. This implementation demonstrates stronger local pattern capture capabilities in aquifers with complex fault zones and is suitable for fault densities greater than 0.5 faults / km. 2 Geological units.

[0091] Furthermore, the system executes a periodic calibration process. At the beginning of each month, the system automatically compares the predicted concentrations with the measured concentrations of all monitoring wells over the past 30 days, calculating the root mean square error (RMSE) and the Nash efficiency coefficient (NSE). If the NSE is below 0.65, a model retraining process is triggered: the most recent 180 days of data are extracted from the database as a new training set, the static graph convolution parameters of the network's bottom layer (the first two layers of MR-GCL) are frozen, only the dynamic branches and time module parameters are fine-tuned, and the production model is updated after retraining for 100 epochs. This mechanism ensures that the model adapts to the slow changes in aquifer parameters and the updates of the monitoring network over the long term.

[0092] The physical interpretability of the method is quantitatively analyzed using gradient backpropagation. After each prediction, the system calculates the loss function for the elements of the dynamic adjacency matrix. partial derivatives This value reflects the intensity of the influence of specific flow field connections on the prediction results. Edges with high gradient values ​​are identified as critical migration channels, and their spatial trajectories are vector-superimposed to generate a main pollution transport path map, which serves as the basis for decision-making regarding the placement of interceptor wells in remediation projects.

[0093] Furthermore, the calculation of the dynamic convection-dominant adjacency matrix adopts Scheme B under data sparsity conditions, which directly calculates the hydraulic gradient based on observed water levels. In this case, the system requires at least one synchronized water level observation from the entire network every 72 hours, using a network of mobile pressure gauges to ensure time consistency error is less than 10 minutes. For time steps without synchronized observations, linear interpolation is used to supplement water level data, and confidence markers are introduced into the adjacency matrix. ,when At hour, The final weights were adjusted to This is used to reduce errors introduced by asynchronous data.

[0094] The method introduces lag order expansion into the model input features. In addition to the current time-time features, first-order difference features are also constructed. With second-order difference features These represent the rate of concentration change and acceleration, respectively. This significantly improves the model's ability to identify pollution breakthrough events early, increasing the lead time for breakthrough point prediction from an average of 3.2 days to 6.8 days in the test case.

[0095] Example 3: A soil and groundwater pollution diffusion prediction system based on spatiotemporal graph neural networks, comprising:

[0096] The data acquisition and processing module is used to acquire pollutant concentration time series, groundwater level observation data and external driving data of monitoring well nodes, and to perform unified spatiotemporal benchmark correction, missing value interpolation and completion and format standardization preprocessing on the data;

[0097] The graph structure generation module is used to generate a static hydraulic connectivity matrix based on the spatial coordinates and hydrogeological parameters of the monitoring well nodes, and to calculate the dynamic convection dominance matrix by combining the head field data output by the three-dimensional unsteady flow numerical simulation engine, and to construct a multi-relationship dynamic graph structure that integrates the static medium conduction capacity and the dynamic flow field driving direction.

[0098] A spatiotemporal graph neural network prediction model is used to process graph convolution operations of static hydraulic connectivity matrix and dynamic convection dominance matrix in parallel using multiple relation graph convolutional layers, integrate node spatial features, capture the temporal dynamic evolution of pollutant concentration, and output the predicted value of pollutant concentration at future time.

[0099] The results output and visualization module is used to perform rolling multi-step prediction tasks and generate structured data streams and geographic information system visualization products.

[0100] It should be noted that the soil and groundwater pollution diffusion prediction system can implement the steps in the soil and groundwater pollution diffusion prediction method in the above embodiments and achieve the same technical effect. Refer to the description in the above embodiments, which will not be elaborated here.

[0101] The above description is based on the preferred embodiments of the present invention. It will be apparent to those skilled in the art that the present invention is not limited to the details of the exemplary embodiments described above, and that the invention can be implemented in other specific forms without departing from its spirit or essential characteristics. Therefore, the embodiments should be considered exemplary and non-limiting in all respects. The scope of the invention is defined by the appended claims rather than the foregoing description, and all variations falling within the meaning and scope of equivalents of the claims are intended to be included within the present invention. No reference numerals in the claims should be construed as limiting the scope of the claims.

[0102] Furthermore, it should be understood that although this specification describes embodiments, not every embodiment contains only one independent technical solution. This narrative style is merely for clarity. Those skilled in the art should consider the specification as a whole, and the technical solutions in each embodiment can also be appropriately combined to form other embodiments that can be understood by those skilled in the art.

Claims

1. A method for predicting the diffusion of soil and groundwater pollution based on spatiotemporal graph neural networks, characterized in that, Includes the following steps: S1. Acquire multi-source soil and groundwater monitoring data and external driving data, and perform spatiotemporal reference calibration and preprocessing; S2. Based on the spatial coordinates and hydrogeological parameters of the monitoring nodes, a static hydraulic connectivity matrix is ​​generated. Combined with the head field data from numerical simulation, a dynamic convection dominance matrix is ​​generated, and a multi-relationship dynamic graph structure is constructed. S3. Using a spatiotemporal graph neural network, the graph convolution operations of the static hydraulic connectivity matrix and the dynamic convection-dominant adjacency matrix are processed in parallel through multiple relational graph convolutional layers. After fusing the spatial features of nodes, the temporal dynamic features are captured through a gated recurrent unit to obtain the pollutant concentration prediction results. S4. Perform rolling multi-step forecasts, periodically inject actual observation data to optimize the forecast process, and output structured forecast data and visualization products.

2. The method for predicting the diffusion of soil and groundwater pollution based on spatiotemporal graph neural networks according to claim 1, characterized in that: In step S1, the multi-source soil groundwater monitoring data includes pollutant concentration time series and groundwater level observation data; the external driving data includes daily precipitation, potential evapotranspiration and artificial water injection data.

3. The method for predicting the diffusion of soil and groundwater pollution based on spatiotemporal graph neural networks according to claim 1, characterized in that: In step S2, the static hydraulic connectivity matrix is ​​generated as follows: the initial weights are calculated based on the permeability coefficient, porosity, and Euclidean distance between the monitoring well nodes, and then normalized to obtain a normalized static adjacency matrix.

4. The method for predicting the diffusion of soil and groundwater pollution based on spatiotemporal graph neural networks according to claim 1, characterized in that: In step S2, the dynamic convection dominance matrix is ​​generated by inputting driving external input data into a groundwater flow numerical simulation engine to solve for the head distribution field. The hydraulic gradient between nodes is calculated based on the head value of each monitoring well. The convection weight is determined by combining the equivalent permeability coefficient and then normalized to obtain the normalized dynamic adjacency matrix.

5. The method for predicting the diffusion of soil and groundwater pollution based on a spatiotemporal graph neural network according to claim 4, characterized in that: The head distribution field is obtained by solving the groundwater flow control equation, which takes into account the permeability coefficient tensor components, specific yield, and source and sink terms.

6. The method for predicting the diffusion of soil and groundwater pollution based on spatiotemporal graph neural networks according to claim 1, characterized in that: In step S3, the multi-relation graph convolutional layer performs the following operations: static relation graph convolution and dynamic relation graph convolution are performed on the hidden states of the input nodes respectively, and the outputs of the two branches are concatenated and fused into a unified node representation through linear transformation.

7. The method for predicting the diffusion of soil and groundwater pollution based on spatiotemporal graph neural networks according to claim 1, characterized in that: In step S3, the gated loop unit receives node feature sequences from multiple consecutive time steps, calculates the current hidden state by updating the gate, resetting the gate, and candidate hidden states, and finally maps it to the predicted pollutant concentration value through a fully connected layer.

8. The method for predicting the diffusion of soil and groundwater pollution based on spatiotemporal graph neural networks according to claim 1, characterized in that: In step S3, the loss function of the spatiotemporal graph neural network includes a mean squared error term and a physical consistency regularization term based on the principle of conservation of pollutant mass, which is used to strengthen the model's adherence to physical laws.

9. The method for predicting the diffusion of soil and groundwater pollution based on spatiotemporal graph neural networks according to claim 1, characterized in that: In step S4, the rolling multi-step prediction process injects actual observation data after each preset number of steps to reinitialize the feature sequence in order to reduce error accumulation.

10. A prediction system based on the method for predicting the diffusion of soil and groundwater pollution based on a spatiotemporal graph neural network according to any one of claims 9, characterized in that, include: The data acquisition and processing module is used to acquire multi-source soil and groundwater monitoring data and external driving data, and to perform preprocessing on the data; The graph structure generation module is used to generate static hydraulic connectivity matrices and dynamic convection dominance matrices, and to construct a multi-relationship dynamic graph structure that integrates the static medium conduction capacity and the dynamic flow field driving direction. A spatiotemporal graph neural network prediction model is used to process graph convolution operations of static hydraulic connectivity matrix and dynamic convection dominance matrix in parallel, fuse node spatial features, and output predicted pollutant concentration values ​​for future time moments. The results output and visualization module is used to perform rolling multi-step prediction tasks and generate structured data streams and geographic information system visualization products.

Images