Geological environment risk prediction method and system based on multi-source heterogeneous data
By constructing a three-dimensional spatial attribute matrix and a hysteresis tensor matrix, and combining it with the groundwater velocity vector field, the reverse migration trajectory of pollutants is dynamically and iteratively tracked, solving the problem of insufficient source tracing accuracy in existing technologies, and realizing the precise location and risk prediction of carbon tetrachloride leakage sources.
Patent Information
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- SICHUAN GEOLOGICAL ENVIRONMENT SURVEY & RES CENT
- Filing Date
- 2026-04-29
- Publication Date
- 2026-05-29
AI Technical Summary
Existing technologies fail to effectively consider the delayed hysteresis characteristics of pollutants in aquifers when tracing the source of carbon tetrachloride groundwater pollution, resulting in insufficient accuracy in tracing the source and making it impossible to accurately locate the source of historical leaks.
By constructing a three-dimensional spatial attribute matrix and a hysteresis tensor matrix, and combining it with the groundwater velocity vector field, virtual particles are tracked in reverse through step-by-step dynamic iteration to obtain the reverse migration trajectory of pollutants and accurately locate the source of the leak.
It achieves a high degree of alignment between the pollutant migration trajectory and the actual leak source, improving the accuracy of source tracing and supporting the accurate location and prediction of geological and environmental risks.
Smart Images

Figure CN122109501A_ABST
Abstract
Description
Technical Field
[0001] This invention relates to the field of processing and detection technology, and in particular to a method and system for predicting geological environmental risks based on multi-source heterogeneous data. Background Technology
[0002] An investigation into groundwater pollution in chemical industrial parks is a necessary prerequisite for site remediation and risk management. Carbon tetrachloride, a common chemical raw material, often enters and remains in the groundwater system due to early pipeline leaks. To clarify the main responsibility for remediation and carry out targeted source control projects, it is essential to accurately trace and pinpoint the historical sources of carbon tetrachloride leaks.
[0003] In the current groundwater pollution source tracing project, the mainstream approach is to rely on numerical simulation technology to conduct reverse hydrodynamic tracing. Engineers establish a groundwater flow field model based on site survey data, release tracking particles starting from the point where carbon tetrachloride exceeds the standard in the current monitoring well, extract the water flow velocity field data in the model, perform inversion calculations in the opposite direction of water flow, and deduce the spatial coordinates of the particles at past time points, thereby inferring the location of the pollution source.
[0004] However, in actual engineering scenarios, existing solutions have a serious flaw that affects the accuracy of source tracing. When performing historical trajectory reversal, existing source tracing systems directly equate the reverse flow velocity of pure groundwater itself with the historical reverse reversal velocity of carbon tetrachloride. They do not introduce the hysteresis characteristic parameters of carbon tetrachloride during actual migration in the aquifer medium into the numerical inversion calculation system. This neglect of the difference in material migration velocity at the algorithm level will lead to a serious spatial deviation and misalignment between the final calculated source tracing trajectory and the actual pollution leakage location in history, which cannot meet the engineering requirements for accurate location and responsibility determination of historical pollution sources. Summary of the Invention
[0005] The main objective of this invention is to provide a geological environment risk prediction method and system based on multi-source heterogeneous data. This addresses the problem that existing technologies directly use the reverse flow velocity of pure water as the velocity for pollutant source tracing, which makes it difficult to obtain the difference in migration velocity between carbon tetrachloride and groundwater, resulting in a serious spatial misalignment between the predicted trajectory and the actual leakage source.
[0006] To achieve the above objectives, this invention provides a geological environment risk prediction method based on multi-source heterogeneous data, the method comprising the following steps: Acquire multi-source heterogeneous geological environment data of the target area, and construct a three-dimensional spatial attribute matrix and a three-dimensional groundwater velocity vector field based on the multi-source heterogeneous geological environment data. The three-dimensional spatial attribute matrix includes at least the effective porosity and soil organic carbon mass fraction of each spatial grid. Obtain the organic carbon partition coefficient of the target pollutant, and construct a three-dimensional spatial hysteresis tensor matrix based on the organic carbon partition coefficient and the three-dimensional spatial attribute matrix. The three-dimensional spatial hysteresis tensor matrix is used to characterize the heterogeneous interception and hindrance differences of different spatial grids on the migration velocity of the target pollutant. The current pollution monitoring location of the target area is obtained as the starting point for tracking. Based on the three-dimensional groundwater velocity vector field and the three-dimensional spatial hysteresis tensor matrix, the virtual particles located at the starting point are subjected to step-by-step dynamic iterative reverse tracking based on spatial grid attributes to obtain the reverse migration trajectory of the target pollutants at different historical time nodes. The actual source location of the leak that caused the geological environment was obtained based on the reverse migration trajectory of the target pollutant, and the geological environment risk was predicted based on the actual leak source location.
[0007] Optionally, the process of constructing the three-dimensional spatial attribute matrix includes the following steps: The effective porosity and organic carbon mass fraction of undisturbed soil samples at different depths within the target area were measured. Spatial three-dimensional interpolation is performed on the effective porosity measurement value and the organic carbon mass fraction measurement value to generate a three-dimensional effective porosity matrix and a three-dimensional soil organic carbon mass fraction matrix covering the target area; The three-dimensional effective porosity matrix and the three-dimensional soil organic carbon mass fraction matrix are fused using a data structure to generate the three-dimensional spatial attribute matrix.
[0008] Optionally, the process of constructing the three-dimensional groundwater velocity vector field includes the following steps: Obtain the aquifer permeability coefficient and hydrological boundary conditions of the target area; A three-dimensional unsteady groundwater flow model is established based on the aquifer permeability coefficient and the hydrological boundary conditions. Based on the three-dimensional unsteady groundwater flow model, the output includes the groundwater velocity vector of each spatial grid at each discrete time step. The groundwater velocity vectors of all spatial grids are combined to form the three-dimensional groundwater velocity vector field.
[0009] Optionally, constructing a three-dimensional spatial hysteresis tensor matrix based on the organic carbon partition coefficient and the three-dimensional spatial attribute matrix includes: Based on any target spatial grid in the three-dimensional spatial attribute matrix, obtain the soil bulk density corresponding to the target spatial grid, and call the effective porosity and soil organic carbon mass fraction corresponding to the target spatial grid; The hysteresis factor corresponding to the target spatial grid is obtained based on the soil bulk density, the effective porosity, and the soil organic carbon mass fraction. Obtain the hysteresis factors corresponding to all spatial grids, and arrange all hysteresis factors in a matrix to generate the three-dimensional spatial hysteresis tensor matrix.
[0010] Optionally, based on the three-dimensional groundwater velocity vector field and the three-dimensional spatial hysteresis tensor matrix, a step-by-step dynamic iterative reverse tracking based on spatial grid attributes is performed on the virtual particle located at the tracking starting point, including: Preset the discrete time step for reverse inference; Determine the current spatial grid position of the virtual particle at the current simulation moment; Based on the current groundwater velocity vector corresponding to the current spatial grid in the three-dimensional groundwater velocity vector field, and the hysteresis factor corresponding to the three-dimensional spatial hysteresis tensor matrix, the instantaneous equivalent inverse velocity vector of the virtual particle in the current spatial grid is obtained.
[0011] Optionally, after obtaining the instantaneous equivalent inverse velocity vector of the virtual particle in the current spatial grid, the method further includes: The inverse spatial displacement of the virtual particle within the current discrete time step is obtained by multiplying the instantaneous equivalent inverse velocity vector with the discrete time step. The virtual particle's spatial position vector at the current simulation moment is subtracted from the inverse spatial displacement to complete the iterative update and obtain the spatial position vector of the previous historical moment. Repeat the steps of spatial grid judgment, instantaneous equivalent reverse velocity vector acquisition and spatial position vector update until the set total backtracking time is reached, so as to form the reverse migration trajectory of the target pollutant composed of continuous historical spatial position vectors.
[0012] Optionally, the target pollutant includes carbon tetrachloride.
[0013] Optionally, the geological environmental risk prediction based on the actual leak source location includes the following steps: The actual leak source location is taken as the starting point for the forward diffusion evolution; Based on the three-dimensional groundwater velocity vector field, the three-dimensional spatial hysteresis tensor matrix, and the natural decay coefficient of the target pollutant, a positive geological evolution prediction model is constructed. Based on the positive geological evolution prediction model, the three-dimensional concentration distribution field of the target pollutant diffuses to the surrounding area within a set time period in the future is obtained; Based on the three-dimensional concentration distribution field and the preset environmental health risk assessment standards, the geological environmental risk prediction results of the target area are obtained and output.
[0014] Optionally, after performing geological environmental risk prediction based on the actual leak source location, the method further includes: Based on the geological environment risk prediction results, new geological environment monitoring points will be set up in the target area where the three-dimensional concentration distribution field exceeds the safety risk threshold. The actual feedback data of the newly added geological environment monitoring points is obtained, and the actual feedback data is fed back to the three-dimensional groundwater velocity vector field for closed-loop iterative correction.
[0015] To achieve the above objectives, the present invention also provides a geological environment risk prediction system based on multi-source heterogeneous data, the system comprising: A data flow field construction module is used to acquire multi-source heterogeneous geological environment data of the target area, and construct a three-dimensional spatial attribute matrix and a three-dimensional groundwater velocity vector field based on the multi-source heterogeneous geological environment data. The three-dimensional spatial attribute matrix includes at least the effective porosity and soil organic carbon mass fraction of each spatial grid. The hysteresis tensor analysis module is used to obtain the organic carbon partition coefficient of the target pollutant. Based on the organic carbon partition coefficient and the three-dimensional spatial attribute matrix, a three-dimensional spatial hysteresis tensor matrix is constructed. The three-dimensional spatial hysteresis tensor matrix is used to characterize the heterogeneous interception and hindrance differences of different spatial grids on the migration velocity of the target pollutant. The dynamic iterative source tracing module is used to obtain the current pollution monitoring location of the target area as the tracking starting point. Based on the three-dimensional groundwater flow velocity vector field and the three-dimensional spatial hysteresis tensor matrix, the module performs step-by-step dynamic iterative reverse tracking of the virtual particles located at the tracking starting point based on spatial grid attributes to obtain the reverse migration trajectory of the target pollutants at different historical time nodes. The risk prediction and locking module is used to obtain the actual source location of the geological environment leakage based on the reverse migration trajectory of the target pollutant, and to perform geological environment risk prediction based on the actual source location.
[0016] The beneficial effects that this invention can achieve are as follows: This invention acquires multi-source heterogeneous geological data of the target area to construct a three-dimensional spatial attribute matrix containing effective porosity and soil organic carbon mass fraction, as well as a three-dimensional groundwater velocity vector field. Then, it combines the organic carbon partition coefficient of the target pollutant to construct a three-dimensional spatial hysteresis tensor matrix to characterize the differences in interception and hindrance of the heterogeneous medium. Finally, based on the velocity field and the hysteresis tensor matrix, it conducts step-by-step dynamic iterative reverse tracking of virtual particles based on grid attributes. This solves the problem that existing source tracing technologies directly equate the pure groundwater dynamic velocity with the actual migration velocity of pollutants, making it difficult to decouple the velocity attenuation of pollutants due to adsorption in heterogeneous strata, thus causing serious spatial misalignment and off-target positioning between historical trajectory and the actual leakage source. Specifically, this invention utilizes a three-dimensional spatial hysteresis tensor matrix to assign a unique hysteresis feature to each microscopic spatial grid at the data layer. In historical source tracing and simulation, when the virtual particle performs a backtracking at a discrete time step, it accurately reads the flow velocity and unique hysteresis factor of the current grid and performs instantaneous dynamic deceleration. This rigorous step-by-step dynamic iteration mechanism successfully decouples the physical groundwater drive from the chemical medium interception, eliminating the engineering hazard of the virtual particle overtaking the real leakage source due to excessive speed during simulation. This ensures that the final generated reverse migration trajectory closely matches the real physical path of the pollutant's slow evolution underground. Attached Figure Description
[0017] To more clearly illustrate the specific embodiments of the present invention or the technical solutions in the prior art, the accompanying drawings used in the description of the specific embodiments or the prior art will be briefly introduced below. In all the drawings, similar elements or parts are generally identified by similar reference numerals. In the drawings, the elements or parts are not necessarily drawn to scale.
[0018] Figure 1 This is a flowchart illustrating the method in Embodiment 1 of the present invention; Figure 2 This is a structural block diagram of the system in Embodiment 2 of the present invention.
[0019] The realization of the objective, functional features and advantages of the present invention will be further explained in conjunction with the embodiments and with reference to the accompanying drawings. Detailed Implementation
[0020] The technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only a part of the embodiments of the present invention, and not all of the embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those of ordinary skill in the art without creative effort are within the scope of protection of the present invention.
[0021] It should be noted that all directional indications (such as up, down, left, right, front, back, etc.) in the embodiments of the present invention are only used to explain the relative positional relationship and movement of each component in a specific posture. If the specific posture changes, the directional indication will also change accordingly.
[0022] In this invention, unless otherwise explicitly specified and limited, the terms "connection," "fixed," etc., should be interpreted broadly. For example, "connection" can be a fixed connection, a detachable connection, or an integral part; it can be a mechanical connection or an electrical connection; it can be a direct connection or an indirect connection through an intermediate medium; it can be the internal communication of two components or the interaction between two components, unless otherwise explicitly limited. Those skilled in the art can understand the specific meaning of the above terms in this invention according to the specific circumstances.
[0023] Furthermore, if the embodiments of this invention involve descriptions such as "first" or "second," these descriptions are for descriptive purposes only and should not be construed as indicating or implying their relative importance or implicitly specifying the number of technical features indicated. Therefore, a feature defined with "first" or "second" may explicitly or implicitly include at least one of those features. Additionally, the meaning of "and / or" throughout the text includes three parallel solutions; for example, "A and / or B" includes solution A, solution B, or a solution where both A and B are satisfied simultaneously. Furthermore, the technical solutions of the various embodiments can be combined with each other, but this must be based on the ability of those skilled in the art to implement them. When the combination of technical solutions is contradictory or impossible to implement, it should be considered that such a combination of technical solutions does not exist and is not within the scope of protection claimed by this invention.
[0024] Example 1: Reference Figure 1 This embodiment provides a geological environment risk prediction method based on multi-source heterogeneous data, the method including the following steps: Acquire multi-source heterogeneous geological environment data of the target area, and construct a three-dimensional spatial attribute matrix and a three-dimensional groundwater velocity vector field based on the multi-source heterogeneous geological environment data. The three-dimensional spatial attribute matrix includes at least the effective porosity and soil organic carbon mass fraction of each spatial grid. Obtain the organic carbon partition coefficient of the target pollutant, and construct a three-dimensional spatial hysteresis tensor matrix based on the organic carbon partition coefficient and the three-dimensional spatial attribute matrix. The three-dimensional spatial hysteresis tensor matrix is used to characterize the heterogeneous interception and hindrance differences of different spatial grids on the migration velocity of the target pollutant. The current pollution monitoring location of the target area is obtained as the starting point for tracking. Based on the three-dimensional groundwater velocity vector field and the three-dimensional spatial hysteresis tensor matrix, the virtual particles located at the starting point are subjected to step-by-step dynamic iterative reverse tracking based on spatial grid attributes to obtain the reverse migration trajectory of the target pollutants at different historical time nodes. The actual source location of the leak that caused the geological environment was obtained based on the reverse migration trajectory of the target pollutant, and the geological environment risk was predicted based on the actual leak source location.
[0025] It should be noted that when dealing with complex industrial sites with a long history of chemical leaks and risks of geological evolution, the underground aquifers and vadose zones exhibit extremely strong heterogeneous characteristics, containing interwoven silty clay lenses, medium-coarse sand layers, and other diverse geological structures. Therefore, the system first needs to acquire multi-source heterogeneous geological environmental data of the target area. In actual engineering operations, this data mainly originates from high-density drilling and sampling in the field, in-situ hydrogeological tests, and laboratory physicochemical analysis of undisturbed soil samples.
[0026] After receiving these multi-dimensional discrete exploration data, the system uses spatial interpolation algorithms in three-dimensional geostatistics to transform the originally discrete data points into a continuous three-dimensional spatial attribute matrix covering the entire target area. This matrix is equivalent to establishing a digital entity for the invisible heterogeneous space underground. Each spatial grid inside corresponds to key geological attributes such as effective porosity and soil organic carbon mass fraction. At the same time, based on the obtained aquifer permeability coefficient and boundary conditions, the system drives the groundwater flow dynamics model to solve and output the three-dimensional groundwater velocity vector field of the entire field at different discrete time steps.
[0027] It should also be noted that existing detection models often directly equate the pure water flow velocity with the actual migration velocity of pollutants, leading to severe misalignment and distortion in the reverse trajectory calculation when traversing heterogeneous strata. To address this issue, the system first obtains the organic carbon distribution coefficient of the target pollutant. The inducing source in the geological environment is often a non-aqueous liquid with strong hydrophobicity. When such substances flow through porous media, they spontaneously and continuously undergo dynamic distribution and adsorption onto the organic matter on the surface of soil particles. Through its built-in underlying computational logic, the system couples the inherent organic carbon distribution coefficient of the target pollutant with physical properties such as the effective porosity, soil organic carbon mass fraction, and soil bulk density of the specific spatial grid. For each discrete grid in the three-dimensional spatial attribute matrix, its specific hysteresis drag multiple is calculated, and then integrated to generate a three-dimensional spatial hysteresis tensor matrix. This tensor matrix breaks the erroneous assumption of constant drag across the entire field in traditional models, enabling the model to accurately quantify and identify the differences in nonlinear retention and drag of pollutants by local micro-strata.
[0028] Subsequently, the system obtains the current pollution monitoring location coordinates of the target area as the starting point for tracking and begins to execute a step-by-step dynamic iterative reverse tracking based on spatial grid attributes. Considering that historical pollution tracing tasks typically span decades, the cumulative time can easily cause an exponential amplification of spatial misalignment errors. Therefore, a rigorous discrete iterative deduction logic is adopted. In the specific reverse deduction process, before each small time step of tracing, the virtual particle located at the tracking starting point is forcibly positioned by the system in the spatial grid. The system then calls the current hydrodynamic velocity and specific hysteresis resistance corresponding to that specific grid in real time from the three-dimensional groundwater velocity vector field and the three-dimensional spatial hysteresis tensor matrix. By dividing the called pure hydrodynamic velocity by the specific hysteresis resistance of that grid, the true instantaneous equivalent reverse velocity of the pollutant at that microscopic location is extracted. Then, combined with the set discrete time step, the actual reverse spatial displacement within this time period is calculated, and the historical spatial position of the particle is iteratively updated, effectively preventing the virtual particle from crossing the real historical leakage point due to excessive speed in the model.
[0029] Through the rigorous and realistic dynamic iterative simulation described above, the reverse migration trajectory of the virtual particles will eventually converge to a specific endpoint or a high-probability polygonal region in the three-dimensional model space. Based on the final convergence region of this trajectory, the system can accurately identify and pinpoint the actual source of leakage that historically caused the deterioration of the geological environment in the area.
[0030] After obtaining the actual location of the leak source, the system uses it as the initial input condition for the positive geological evolution prediction model. Combining the natural hydrological and meteorological evolution patterns of the target area, the natural decay properties of pollutants, and environmental health risk assessment standards, the system calculates the diffusion trend of pollutants to surrounding sensitive receptors within a set time period in the future, thereby outputting high-confidence geological environmental risk prediction results. This provides accurate data support and decision-making basis for risk management, responsibility identification, and subsequent engineering remediation planning for the entire site.
[0031] In this embodiment, the process of constructing the three-dimensional spatial attribute matrix includes the following steps: The effective porosity and organic carbon mass fraction of undisturbed soil samples at different depths within the target area were measured. Spatial three-dimensional interpolation is performed on the effective porosity measurement value and the organic carbon mass fraction measurement value to generate a three-dimensional effective porosity matrix and a three-dimensional soil organic carbon mass fraction matrix covering the target area; The three-dimensional effective porosity matrix and the three-dimensional soil organic carbon mass fraction matrix are fused using a data structure to generate the three-dimensional spatial attribute matrix.
[0032] In actual geological environment investigation and risk assessment projects, in order to ensure the accuracy of physical environment data and the high confidence of subsequent source tracing and inference, the specific implementation details and processing logic of the three-dimensional spatial attribute matrix construction process described in the embodiments are explained in detail here. Among them, the engineering execution end will conduct gridded high-density drilling in the site area and collect undisturbed soil samples from the strata at set depth intervals. After these undisturbed soil samples are processed by geotechnical testing instruments and carbon analysis equipment, the system can receive the real physicochemical parameters of each discrete sampling point. Effective porosity characterizes the proportion of interconnected physical spaces in underground soil or rock that allow groundwater to pass through; while the soil organic carbon mass fraction characterizes the chemical adsorption and retention potential of the formation medium for hydrophobic target pollutants such as heavy non-aqueous phase liquids.
[0033] Since the number and spatial location of field borehole samples are always finite and discrete, this embodiment prioritizes the use of a three-dimensional spatial Kriging interpolation formula based on geostatistics for the specific processing. If the system uses a simple arithmetic mean or linear transition algorithm, it will ignore the lenticular abrupt changes in the strata or the anisotropy in the direction of sedimentary bedding, leading to model distortion. Therefore, this embodiment introduces an unbiased optimal estimation interpolation logic, which assigns specific spatial weights to known sampling points at different distances and orientations by solving the spatial variogram at the bottom layer. The expression satisfies: ; In the formula, The estimated attribute value represents the unknown spatial grid of the target. This attribute value specifically points to the estimated effective porosity or organic carbon mass fraction. The data comes from the calculation output result after the system runs the formula in this interpolation step. The three-dimensional spatial coordinate vector representing the unknown spatial grid of the target in the three-dimensional geological model; This represents the total number of valid known sampling points participating in this interpolation calculation within the search range defined by the unknown target space grid; Representing the The spatial weight coefficients assigned to the unknown spatial grid of the target by each known sampling point are dimensionless. The data comes from the weight allocation values obtained in real time by the underlying system based on the spatial variogram (semivariance function) and the unbiased optimal estimation condition matrix. Representing the Actual measured values of a known sampling point; Representing the The three-dimensional spatial coordinate vector of a known sampling point in a three-dimensional geological model.
[0034] Based on the above expression, the estimated attribute value of the target unknown spatial grid in the three-dimensional geological model is equal to the sum of the actual measured values of all valid known sampling points within the search range set by the target unknown spatial grid and their corresponding spatial weight coefficients.
[0035] Furthermore, in regions with distinct horizontal sedimentary facies, the above expression can identify that the continuity of properties in the horizontal direction is much higher than that in the vertical direction, thus enabling calculations... Higher weights are assigned to adjacent points in the horizontal direction. This processing logic, which closely matches the anisotropy of real geology, enables the system to smoothly and physically accurately reconstruct the outline of deeply buried clay water-blocking layers or high-permeability sand layers in digital space. Finally, the system performs the above interpolation operation on the effective porosity measurement value and the organic carbon mass fraction measurement value, respectively, and successfully generates a three-dimensional effective porosity matrix and a three-dimensional soil organic carbon mass fraction matrix that cover the entire target area without any spatial discontinuities.
[0036] Furthermore, the system integrates the three-dimensional effective porosity matrix and the three-dimensional soil organic carbon mass fraction matrix into a unified data structure. Using a standardized three-dimensional mesh generation standard, the system co-encapsulates mesh cells with identical spatial coordinate systems from both matrices into a single data structure. During subsequent step-by-step dynamic iterative reverse tracing, when a virtual particle moves to any microscopic coordinate point, the system can instantly and synchronously extract the effective porosity and organic carbon mass fraction at that location using a single mesh index pointer. This not only completely eliminates the risk of mismatch in spatial matching of multi-source heterogeneous data but also significantly improves the underlying operational efficiency of computing devices when calculating the three-dimensional spatial hysteresis tensor matrix.
[0037] In this embodiment, the construction process of the three-dimensional groundwater velocity vector field includes the following steps: Obtain the aquifer permeability coefficient and hydrological boundary conditions of the target area; A three-dimensional unsteady groundwater flow model is established based on the aquifer permeability coefficient and the hydrological boundary conditions. Based on the three-dimensional unsteady groundwater flow model, the output includes the groundwater velocity vector of each spatial grid at each discrete time step. The groundwater velocity vectors of all spatial grids are combined to form the three-dimensional groundwater velocity vector field.
[0038] Understandably, the aquifer permeability coefficient is a core parameter characterizing the inherent physical ability of porous media to allow groundwater fluids to pass through. The engineering execution end obtains permeability coefficient data for different aquifers and weakly permeable layers through in-situ hydrogeological testing methods such as pumping tests, micro-water tests, and well logging curve interpretation. Meanwhile, hydrological boundary conditions constitute the external control factors for model calculation. The system receives external dynamic constraint data, including the water level of surface water bodies such as rivers around the site (given head boundary), the regional atmospheric rainfall infiltration recharge (flow boundary), and the dynamic operating parameters of production pumping wells or injection wells on site (source and sink terms).
[0039] It is also understandable that conventional groundwater flow models often only output Darcy velocity (i.e., apparent specific flow rate). However, Darcy velocity assumes that the water flow fills the entire cross-section of the aquifer. In the real micro-geological environment, groundwater can only flow within interconnected effective pore channels. Since the actual cross-sectional area of the water flow is much smaller than the total cross-sectional area of the aquifer, the actual water flow velocity that drives the migration of pollutants such as heavy non-aqueous liquids will be much greater than Darcy velocity. Therefore, it is necessary to deeply couple the hydraulic gradient and permeability coefficient that are solved in real time with the three-dimensional effective porosity matrix that has been accurately constructed in the previous steps.
[0040] Its specific processing logic satisfies: ; In the formula, The target space grid is represented by the time coordinate. The groundwater velocity vector at that time; The permeability tensor representing the target space grid; The effective porosity of the target spatial grid; The target space grid is represented by the time coordinate. The hydraulic gradient vector at time.
[0041] Based on the above expression, the Darcy velocity is obtained by multiplying the permeability coefficient tensor with the spatial partial derivative (hydraulic gradient), and then further divided by the effective porosity. This successfully reduces the dimensionality of the hydrodynamic results of the macroscopic continuous medium to the microscopic pore scale. More importantly, the calculation logic directly calls the effective porosity data in the three-dimensional spatial attribute matrix generated above across modules. This processing method enables the groundwater flow model to sense and respond to extremely small abrupt changes in the pore structure of the formation when calculating the flow velocity. For example, when facing a high-porosity microfracture developed in a dense rock mass, the system can not only calculate the water head convergence effect at that location, but also accurately reconstruct the extremely high actual water flow surge velocity inside the fracture through the local effective porosity.
[0042] Finally, based on the aforementioned three-dimensional unsteady groundwater flow model and velocity vector calculation logic, the system outputs groundwater velocity vectors for each spatial grid at each discrete time step (e.g., daily or hourly). In memory, the system performs ordered time-axis stacking and spatial combination of these massive velocity vector data spanning a vast three-dimensional space and a long time dimension, formally forming the aforementioned three-dimensional groundwater velocity vector field. At this point, the three-dimensional spatial attribute matrix representing the static physical environment of the target area and the three-dimensional groundwater velocity vector field representing the dynamic fluid driving force have achieved a bottom-level logical integration, providing a complete and continuous data support system for the subsequent seamless retrieval of the target pollutant organic carbon distribution coefficient, the solution of the three-dimensional spatial hysteresis tensor matrix, and the final step-by-step dynamic iterative reverse tracking.
[0043] In this embodiment, constructing a three-dimensional spatial hysteresis tensor matrix based on the organic carbon partition coefficient and the three-dimensional spatial attribute matrix includes: Based on any target spatial grid in the three-dimensional spatial attribute matrix, obtain the soil bulk density corresponding to the target spatial grid, and call the effective porosity and soil organic carbon mass fraction corresponding to the target spatial grid; The hysteresis factor corresponding to the target spatial grid is obtained based on the soil bulk density, the effective porosity, and the soil organic carbon mass fraction. Obtain the hysteresis factors corresponding to all spatial grids, and arrange all hysteresis factors in a matrix to generate the three-dimensional spatial hysteresis tensor matrix.
[0044] In the actual data processing, the system scans each discrete unit in the three-dimensional spatial attribute matrix one by one using a traversal algorithm. For any target spatial grid, the system obtains the soil bulk density corresponding to that grid. During this process, the soil bulk density data can originate from previous field borehole sampling, laboratory physical property measurements, and spatial interpolation calculations. Simultaneously, the system uses internal grid index pointers to synchronously access the effective porosity and soil organic carbon mass fraction already encapsulated in the target spatial grid data structure. It can be understood that effective porosity characterizes the actual physical connectivity within the microgrid that allows fluids (groundwater) to pass through, while the soil organic carbon mass fraction reflects the potential chemical retention and adsorption capacity of the grid medium for target pollutants (especially hydrophobic organic pollutants).
[0045] Based on the soil bulk density, effective porosity, and soil organic carbon mass fraction, combined with the pre-obtained organic carbon partition coefficient of the target pollutant, the hysteresis factor corresponding to the target spatial grid is obtained. Specifically, in the underlying calculation logic of the system: since the transport of the target pollutant in the aquifer is not synchronous with pure water, but is delayed due to the dynamic adsorption of organic matter on the surface of soil particles, the system first calculates the ratio of soil bulk density to effective porosity to characterize the mass of solid medium per unit effective pore volume. Subsequently, the system multiplies this ratio with the soil organic carbon mass fraction and the organic carbon partition coefficient of the target pollutant. This product term is used to accurately quantify the dynamic distribution ratio (i.e., adsorption potential energy) of the target pollutant between the solid phase (soil) and the liquid phase (groundwater). Finally, the system adds a value of one to the above product result to obtain the specific hysteresis factor of the target spatial grid. This calculation step deeply mathematically couples the physical parameters characterizing the geological structure with the chemical parameters characterizing the pollutant characteristics, so that the output hysteresis factor can truly reflect the nonlinear interception and hindrance effect of the specific grid on the target pollutant.
[0046] Finally, by using multi-threaded parallel computing or iterative loop mechanisms, the hysteresis factors corresponding to all spatial grids in the entire target area are obtained. After completing the hysteresis factor calculation for the entire field, all hysteresis factors are matrix-arranged and structurally encapsulated according to the original spatial topology (i.e., three-dimensional coordinate system) of each spatial grid in the three-dimensional geological model to generate the three-dimensional spatial hysteresis tensor matrix.
[0047] By generating this three-dimensional spatial hysteresis tensor matrix, the system upgrades the original static single empirical parameter into a three-dimensional tensor set that can reflect the heterogeneous characteristics of global space. When the source tracing algorithm scans a clay layer grid containing a large amount of organic matter, the hysteresis factor at the corresponding position of the matrix will increase significantly, thereby indicating that the system will greatly reduce the reverse inference speed of virtual particles at this point. When encountering a coarse sand layer grid with scarce organic matter, the hysteresis factor approaches one, indicating that the resistance is minimal. This processing method breaks the technical limitation of using global constant resistance in traditional models, and provides spatial resistance guidance for subsequent high-precision step-by-step dynamic iterative reverse tracing. It effectively avoids spatial misalignment and source misjudgment that occur when the reverse tracing trajectory passes through complex strata.
[0048] In this embodiment, based on the three-dimensional groundwater velocity vector field and the three-dimensional spatial hysteresis tensor matrix, a step-by-step dynamic iterative reverse tracking of the virtual particle located at the tracking starting point is performed, using spatial grid attributes. This includes: Preset the discrete time step for reverse inference; Determine the current spatial grid position of the virtual particle at the current simulation moment; Based on the current groundwater velocity vector corresponding to the current spatial grid in the three-dimensional groundwater velocity vector field, and the hysteresis factor corresponding to the three-dimensional spatial hysteresis tensor matrix, the instantaneous equivalent inverse velocity vector of the virtual particle in the current spatial grid is obtained.
[0049] Understandably, the setting of the discrete time step needs to balance computational accuracy and computational cost. The smaller the step size, the higher the trajectory resolution of the virtual particles in the heterogeneous flow field. Especially in regions where the flow field gradient changes drastically, a smaller step size can capture more detailed flow axis deflection. After the setting is completed, the system defines the spatial coordinates of the monitoring well where pollution exceeds the standard as the tracking starting point and releases virtual particles representing the target pollutant at this point.
[0050] Next, the loop iterative calculation begins, and at each iteration, the current spatial grid where the virtual particle is located at the current simulation moment is determined. Since the underground space has been digitized into a three-dimensional discrete grid system containing row, column, and layer indices, the system locates the grid index where the particle is currently located by mapping the instantaneous spatial coordinates of the virtual particle to the three-dimensional coordinate system.
[0051] Based on the current groundwater velocity vector corresponding to the current spatial grid in the three-dimensional groundwater velocity vector field, and the hysteresis factor corresponding to the three-dimensional spatial hysteresis tensor matrix, the instantaneous equivalent inverse velocity vector of the virtual particle in the current spatial grid is obtained.
[0052] The actual migration velocity of heavy non-aqueous liquids in groundwater is limited by the adsorption effect. Their macroscopic migration velocity always lags behind the actual groundwater flow velocity. If the groundwater flow velocity is directly used for reverse tracking, virtual particles will leap too fast in the reverse direction within the model, causing the ultimately located source to be far beyond the actual leak point. Therefore, it is necessary to perform real-time idling correction on the flow velocity at each microscopic time step, based on the specific resistance environment of the particle's location. The correction process can be understood as follows: the instantaneous equivalent reverse velocity vector of the virtual particle in the current spatial grid is equal to the groundwater flow velocity vector corresponding to that grid at the current moment divided by the grid's specific hysteresis factor.
[0053] When a virtual particle moves into a clay grid with high organic matter content and strong resistance, the large hysteresis factor will significantly reduce the instantaneous equivalent reverse velocity vector of the virtual particle in the current spatial grid, and the displacement of the virtual particle within that step will be shortened accordingly; while when it enters a coarse sand grid with very low resistance, the velocity will recover to a level close to that of water flow.
[0054] In this embodiment, after obtaining the instantaneous equivalent inverse velocity vector of the virtual particle in the current spatial grid, the method further includes: The inverse spatial displacement of the virtual particle within the current discrete time step is obtained by multiplying the instantaneous equivalent inverse velocity vector with the discrete time step. The virtual particle's spatial position vector at the current simulation moment is subtracted from the inverse spatial displacement to complete the iterative update and obtain the spatial position vector of the previous historical moment. Repeat the steps of spatial grid judgment, instantaneous equivalent reverse velocity vector acquisition and spatial position vector update until the set total backtracking time is reached, so as to form the reverse migration trajectory of the target pollutant composed of continuous historical spatial position vectors.
[0055] In the above steps, the reverse spatial displacement of the virtual particle within the current discrete time step is obtained by multiplying the instantaneous equivalent reverse velocity vector by the discrete time step. In actual model simulations, due to the extremely heterogeneous and continuously changing spatial properties of the groundwater flow field and the medium's resistance force, the system cannot directly calculate a long-distance cross-grid displacement. Based on the Taylor expansion and numerical approximation principles in calculus, within a sufficiently small discrete time step, the motion of the virtual particle can be considered as uniform linear motion. By multiplying the instantaneous true migration velocity (after removing the medium's hysteresis effect) by the time slice, the actual physical distance traversed by the particle in three-dimensional space within this extremely short time window can be obtained. The spatial position vector of the virtual particle at the current simulation moment is subtracted from the reverse spatial displacement to complete the iterative update and obtain the spatial position vector of the previous historical moment.
[0056] Understandably, the above steps are essentially an explicit Euler numerical integration algorithm. In highly heterogeneous underground media, if the system calculates a month's displacement at once, the particles are likely to directly penetrate multiple stratigraphic grids with completely different lithologies, leading to severe penetration errors and trajectory distortion. However, through subtraction iteration, the system forces the virtual particles to respond promptly to changes in resistance and velocity each time they cross the grid boundary, thereby suppressing numerical dissipation and truncation errors to the greatest extent.
[0057] Finally, the system repeatedly executes the steps of spatial grid judgment, instantaneous equivalent reverse velocity vector acquisition, and spatial position vector update in the background engine until the set total backtracking time is reached. During the tens of thousands of loop iterations, the system records the spatial position vector after each calculation. When the system clock is backtracked to the set historical reference year (i.e., the total backtracking time), the loop terminates. The system performs smooth spline fitting and graphical connection on all the discrete historical spatial position vectors recorded in memory in time sequence, and finally forms the reverse migration trajectory of the target pollutant composed of continuous historical spatial position vectors.
[0058] In this embodiment, the target pollutant includes carbon tetrachloride. In actual geological environment investigation and assessment projects of chemical industrial parks, old heavy industrial bases, and historically abandoned sites, carbon tetrachloride, as an industrial cleaning agent, extractant, and chemical raw material that has been widely used in history, is the most typical and destructive characteristic pollutant causing complex pollution of deep groundwater and soil.
[0059] Carbon tetrachloride is a typical heavy, non-aqueous liquid. Its physical density is significantly greater than that of groundwater, and it also has extremely strong hydrophobic properties and low solubility. When underground storage tanks or pipelines in a chemical plant area rupture, causing a large-scale leak of carbon tetrachloride at the surface or in shallow layers and infiltrating into the underground geological environment, driven by both its own gravity and capillary forces, it can not only rapidly penetrate the shallow unsaturated vadose zone, but also overcome buoyancy to continue migrating downwards into deeper aquifers until it encounters a dense clay impermeable layer at the bottom, at which point it stops infiltrating and forms a high-concentration pure-phase accumulation pool.
[0060] Traditional methods that rely solely on macroscopic groundwater seepage velocity for time reversal inevitably lead to virtual particle regression speeds far exceeding the actual infiltration speed of pollutants, resulting in severe distortions in the source tracing trajectory and significant misalignments in the spatial location of the responsible source. Therefore, using carbon tetrachloride as the target pollutant in this source prediction method maximizes the activation and validation of the technical advantages of the aforementioned constructed three-dimensional spatial attribute matrix and three-dimensional spatial hysteresis tensor matrix. In the actual data processing steps, the system's underlying algorithm precisely retrieves the inherent organic carbon distribution coefficient of carbon tetrachloride and performs grid-level deep mathematical coupling and tensor operations with the heterogeneous effective porosity and soil organic carbon mass fraction of each discrete grid in the three-dimensional spatial attribute matrix. This allows for extremely precise extraction and quantification of the actual retention resistance multiple of carbon tetrachloride when facing different strata lithologies.
[0061] In this embodiment, the geological environmental risk prediction based on the actual leak source location includes the following steps: The actual leak source location is taken as the starting point for the forward diffusion evolution; Based on the three-dimensional groundwater velocity vector field, the three-dimensional spatial hysteresis tensor matrix, and the natural decay coefficient of the target pollutant, a positive geological evolution prediction model is constructed. Based on the positive geological evolution prediction model, the three-dimensional concentration distribution field of the target pollutant diffuses to the surrounding area within a set time period in the future is obtained; Based on the three-dimensional concentration distribution field and the preset environmental health risk assessment standards, the geological environmental risk prediction results of the target area are obtained and output.
[0062] In this embodiment, after performing geological environmental risk prediction based on the actual leak source location, the method further includes: Based on the geological environment risk prediction results, new geological environment monitoring points will be set up in the target area where the three-dimensional concentration distribution field exceeds the safety risk threshold. The actual feedback data of the newly added geological environment monitoring points is obtained, and the actual feedback data is fed back to the three-dimensional groundwater velocity vector field for closed-loop iterative correction.
[0063] Understandably, in the first step of performing geological environmental risk prediction, the system takes the actual leakage source location obtained from the aforementioned simulation as the starting point for the evolution of forward diffusion, and sets the source coordinates as the boundary condition source term for the continuous or instantaneous release of pollutants into the underground aquifer; subsequently, the system constructs a forward geological evolution prediction model based on the aforementioned three-dimensional groundwater velocity vector field, the three-dimensional spatial hysteresis tensor matrix, and the natural decay coefficient of the target pollutant.
[0064] The model employs a forward evolution partial differential governing equation based on the convection-diffusion-reaction mechanism, the specific expression of which satisfies: ; In the formula, The transient mass concentration of the target pollutant at three-dimensional spatial coordinates (x, y, z) and time node t is represented by the data sourced from the real-time solution output of this forward evolution prediction model. This represents the hysteresis factor corresponding to this specific spatial grid. The tensor representing the groundwater dynamic dispersion coefficient; This represents the groundwater velocity vector of the spatial grid at the corresponding time. The natural decay coefficient represents the target pollutant.
[0065] In the above expression, it represents the rate of change of the concentration of the target pollutant over time at any spatial coordinate and time node. It is equal to the difference between the spatial dispersion term and the convective transport term at that location, divided by the specific hysteresis factor of that particular spatial grid, and then subtracted from the natural decay and degradation term of the target pollutant.
[0066] It has extremely high engineering adaptability in geological environmental risk prediction. At the algorithm level, it deeply decouples and recombines the physical migration (convection and dispersion) and chemical behavior (adsorption lag and natural decay) of pollutants. Dividing by the lag factor accurately simulates the dragging and deceleration effect of heterogeneous strata on the pollutant diffusion front; while subtracting the natural decay term truly restores the mass loss of carbon tetrachloride in the long underground evolution process. By using a numerical solver (such as the finite difference method) to iteratively solve the partial differential equation, the three-dimensional concentration distribution field of the target pollutant spreading to the surrounding area in the future within a set time period (such as the next 10 years or 30 years) can be obtained.
[0067] After obtaining a high-fidelity three-dimensional concentration distribution field, the geological environmental risk prediction results for the target area are further obtained and output based on the three-dimensional concentration distribution field and the preset environmental health risk assessment standards. In the geological environmental assessment system, the absolute value of pollutant concentration is not the final risk indicator; it must be converted into a substantial probability of harm to human health or the ecological environment. At this point, the system introduces a health risk assessment model to calculate the carcinogenic risk index of carbon tetrachloride, a typical carcinogen, whose expression satisfies: ; In the formula, Represents the probability of carcinogenic risk in a specific spatial grid region for the recipient population; This represents the maximum target pollutant concentration predicted for this spatial grid within a specified future time period; The average daily groundwater intake rate representing the recipient population; Represents exposure frequency; Represents the duration of exposure; The average weight of the recipient population; The average time representing the carcinogenic effect; Carcinogenicity slope factor representing the ingestion of target pollutants (such as carbon tetrachloride) via oral intake.
[0068] By traversing the three-dimensional concentration distribution field across the entire domain, the carcinogenic risk probability corresponding to each grid cell is calculated. When the risk level is greater than one in a million, it is usually considered to be an unacceptable geological environmental risk. Based on this, the system generates a three-dimensional risk zoning map containing chromatograms of different risk levels, thereby providing an intuitive and accurate output of the geological environmental risk prediction results for the target area.
[0069] However, geological environment prediction is not a one-time open-loop calculation. Due to the extreme complexity of underground systems, any mathematical model has inherent parameter uncertainties. Therefore, in this embodiment, after predicting the geological environment risk based on the actual leakage source location, the system introduces a dynamic correction mechanism. This mechanism includes: First, based on the geological environment risk prediction results, new geological environment monitoring points are deployed in the target area where the three-dimensional concentration distribution field exceeds the safety risk threshold. Engineers will, based on the three-dimensional risk zoning map output by the system, conduct physical drilling and install high-precision online groundwater level and water quality monitoring sensors at key nodes such as the predicted high-risk plume front and upstream of sensitive receptors (e.g., residential water wells). The actual feedback data (including the actual observed groundwater head and pollutant concentration) from the new geological environment monitoring points is acquired in real time through an IoT interface, and this actual feedback data is fed back to the three-dimensional groundwater velocity vector field for closed-loop iterative correction.
[0070] The objective function for closed-loop iterative correction is expressed as follows: ; In the formula, The total error of the objective function representing the closed-loop correction of the model; This represents the total number of newly added geological environment and water quality monitoring points; The data weighting coefficient for the m-th water quality monitoring point; The target pollutant concentration at the m-th monitoring point is predicted by the model. This represents the actual feedback concentration data transmitted back from the m-th newly added monitoring point; This represents the total number of newly added geological environmental water level monitoring points; The data weighting coefficient represents the nth water level monitoring point; The groundwater head elevation at the nth monitoring point, as predicted by the model; This represents the actual head elevation reported by the nth newly added monitoring point.
[0071] The above expression represents the total error of the objective function of the closed-loop correction, which is equal to the sum of the squares of the differences between the predicted concentration and the actual feedback concentration at all water quality monitoring points multiplied by the water quality weight coefficient, plus the sum of the squares of the differences between the predicted head and the actual feedback head at all water level monitoring points multiplied by the water level weight coefficient.
[0072] To minimize The total error is the optimization target. The aquifer permeability coefficient K and hydrological boundary conditions are continuously adjusted in the background to regenerate a more accurate three-dimensional groundwater velocity vector field. This upgrades the original static, unidirectional prediction model into a site digital twin system with self-learning and self-evolution capabilities. As time progresses and new monitoring data continues to be added, the underlying flow field parameters and hysteresis tensor matrix are continuously refined to better reflect the actual geological prediction system of the site.
[0073] Example 2: See attached document Figure 2 This embodiment provides a geological environment risk prediction system based on multi-source heterogeneous data. The system includes: A data flow field construction module is used to acquire multi-source heterogeneous geological environment data of the target area, and construct a three-dimensional spatial attribute matrix and a three-dimensional groundwater velocity vector field based on the multi-source heterogeneous geological environment data. The three-dimensional spatial attribute matrix includes at least the effective porosity and soil organic carbon mass fraction of each spatial grid. The hysteresis tensor analysis module is used to obtain the organic carbon partition coefficient of the target pollutant. Based on the organic carbon partition coefficient and the three-dimensional spatial attribute matrix, a three-dimensional spatial hysteresis tensor matrix is constructed. The three-dimensional spatial hysteresis tensor matrix is used to characterize the heterogeneous interception and hindrance differences of different spatial grids on the migration velocity of the target pollutant. The dynamic iterative source tracing module is used to obtain the current pollution monitoring location of the target area as the tracking starting point. Based on the three-dimensional groundwater flow velocity vector field and the three-dimensional spatial hysteresis tensor matrix, the module performs step-by-step dynamic iterative reverse tracking of the virtual particles located at the tracking starting point based on spatial grid attributes to obtain the reverse migration trajectory of the target pollutants at different historical time nodes. The risk prediction and locking module is used to obtain the actual source location of the geological environment leakage based on the reverse migration trajectory of the target pollutant, and to perform geological environment risk prediction based on the actual source location.
[0074] The above are merely preferred embodiments of the present invention and do not limit the scope of the patent. Any equivalent structural or procedural transformations made based on the description and drawings of the present invention, or direct or indirect applications in other related technical fields, are similarly included within the scope of patent protection of the present invention.
Claims
1. A geological environment risk prediction method based on multi-source heterogeneous data, characterized in that, The method includes the following steps: Acquire multi-source heterogeneous geological environment data of the target area, and construct a three-dimensional spatial attribute matrix and a three-dimensional groundwater velocity vector field based on the multi-source heterogeneous geological environment data. The three-dimensional spatial attribute matrix includes at least the effective porosity and soil organic carbon mass fraction of each spatial grid. Obtain the organic carbon partition coefficient of the target pollutant, and construct a three-dimensional spatial hysteresis tensor matrix based on the organic carbon partition coefficient and the three-dimensional spatial attribute matrix. The three-dimensional spatial hysteresis tensor matrix is used to characterize the heterogeneous interception and hindrance differences of different spatial grids on the migration velocity of the target pollutant. The current pollution monitoring location of the target area is obtained as the starting point for tracking. Based on the three-dimensional groundwater velocity vector field and the three-dimensional spatial hysteresis tensor matrix, the virtual particles located at the starting point are subjected to step-by-step dynamic iterative reverse tracking based on spatial grid attributes to obtain the reverse migration trajectory of the target pollutants at different historical time nodes. The actual source location of the leak that caused the geological environment was obtained based on the reverse migration trajectory of the target pollutant, and the geological environment risk was predicted based on the actual leak source location.
2. The geological environment risk prediction method based on multi-source heterogeneous data as described in claim 1, characterized in that, The process of constructing the three-dimensional spatial attribute matrix includes the following steps: The effective porosity and organic carbon mass fraction of undisturbed soil samples at different depths within the target area were measured. Spatial three-dimensional interpolation is performed on the effective porosity measurement value and the organic carbon mass fraction measurement value to generate a three-dimensional effective porosity matrix and a three-dimensional soil organic carbon mass fraction matrix covering the target area; The three-dimensional effective porosity matrix and the three-dimensional soil organic carbon mass fraction matrix are fused using a data structure to generate the three-dimensional spatial attribute matrix.
3. The geological environment risk prediction method based on multi-source heterogeneous data as described in claim 1, characterized in that, The process of constructing the three-dimensional groundwater velocity vector field includes the following steps: Obtain the aquifer permeability coefficient and hydrological boundary conditions of the target area; A three-dimensional unsteady groundwater flow model is established based on the aquifer permeability coefficient and the hydrological boundary conditions. Based on the three-dimensional unsteady groundwater flow model, the output includes the groundwater velocity vector of each spatial grid at each discrete time step. The groundwater velocity vectors of all spatial grids are combined to form the three-dimensional groundwater velocity vector field.
4. The geological environment risk prediction method based on multi-source heterogeneous data as described in claim 1, characterized in that, The step of constructing a three-dimensional spatial hysteresis tensor matrix based on the organic carbon partition coefficient and the three-dimensional spatial attribute matrix includes: Based on any target spatial grid in the three-dimensional spatial attribute matrix, obtain the soil bulk density corresponding to the target spatial grid, and call the effective porosity and soil organic carbon mass fraction corresponding to the target spatial grid; The hysteresis factor corresponding to the target spatial grid is obtained based on the soil bulk density, the effective porosity, and the soil organic carbon mass fraction. Obtain the hysteresis factors corresponding to all spatial grids, and arrange all hysteresis factors in a matrix to generate the three-dimensional spatial hysteresis tensor matrix.
5. The geological environment risk prediction method based on multi-source heterogeneous data as described in claim 4, characterized in that, Based on the three-dimensional groundwater velocity vector field and the three-dimensional spatial hysteresis tensor matrix, a step-by-step dynamic iterative reverse tracking method based on spatial grid attributes is performed on the virtual particle located at the tracking starting point, including: Preset the discrete time step for reverse inference; Determine the current spatial grid position of the virtual particle at the current simulation moment; Based on the current groundwater velocity vector corresponding to the current spatial grid in the three-dimensional groundwater velocity vector field, and the hysteresis factor corresponding to the three-dimensional spatial hysteresis tensor matrix, the instantaneous equivalent inverse velocity vector of the virtual particle in the current spatial grid is obtained.
6. The geological environment risk prediction method based on multi-source heterogeneous data as described in claim 5, characterized in that, After obtaining the instantaneous equivalent inverse velocity vector of the virtual particle in the current spatial grid, the method further includes: The inverse spatial displacement of the virtual particle within the current discrete time step is obtained by multiplying the instantaneous equivalent inverse velocity vector with the discrete time step. The virtual particle's spatial position vector at the current simulation moment is subtracted from the inverse spatial displacement to complete the iterative update and obtain the spatial position vector of the previous historical moment. Repeat the steps of spatial grid judgment, instantaneous equivalent reverse velocity vector acquisition and spatial position vector update until the set total backtracking time is reached, so as to form the reverse migration trajectory of the target pollutant composed of continuous historical spatial position vectors.
7. The geological environment risk prediction method based on multi-source heterogeneous data as described in claim 1, characterized in that, The target pollutant includes carbon tetrachloride.
8. The geological environment risk prediction method based on multi-source heterogeneous data as described in claim 1, characterized in that, The geological environmental risk prediction based on the actual leak source location includes the following steps: The actual leak source location is taken as the starting point for the forward diffusion evolution; Based on the three-dimensional groundwater velocity vector field, the three-dimensional spatial hysteresis tensor matrix, and the natural decay coefficient of the target pollutant, a positive geological evolution prediction model is constructed. Based on the positive geological evolution prediction model, the three-dimensional concentration distribution field of the target pollutant diffuses to the surrounding area within a set time period in the future is obtained; Based on the three-dimensional concentration distribution field and the preset environmental health risk assessment standards, the geological environmental risk prediction results of the target area are obtained and output.
9. The geological environment risk prediction method based on multi-source heterogeneous data as described in claim 8, characterized in that, After performing geological environmental risk prediction based on the actual leak source location, the method further includes: Based on the geological environment risk prediction results, new geological environment monitoring points will be set up in the target area where the three-dimensional concentration distribution field exceeds the safety risk threshold. The actual feedback data of the newly added geological environment monitoring points is obtained, and the actual feedback data is fed back to the three-dimensional groundwater velocity vector field for closed-loop iterative correction.
10. A geological environment risk prediction system based on multi-source heterogeneous data, characterized in that, The system includes: A data flow field construction module is used to acquire multi-source heterogeneous geological environment data of the target area, and construct a three-dimensional spatial attribute matrix and a three-dimensional groundwater velocity vector field based on the multi-source heterogeneous geological environment data. The three-dimensional spatial attribute matrix includes at least the effective porosity and soil organic carbon mass fraction of each spatial grid. The hysteresis tensor analysis module is used to obtain the organic carbon partition coefficient of the target pollutant. Based on the organic carbon partition coefficient and the three-dimensional spatial attribute matrix, a three-dimensional spatial hysteresis tensor matrix is constructed. The three-dimensional spatial hysteresis tensor matrix is used to characterize the heterogeneous interception and hindrance differences of different spatial grids on the migration velocity of the target pollutant. The dynamic iterative source tracing module is used to obtain the current pollution monitoring location of the target area as the tracking starting point. Based on the three-dimensional groundwater flow velocity vector field and the three-dimensional spatial hysteresis tensor matrix, the module performs step-by-step dynamic iterative reverse tracking of the virtual particles located at the tracking starting point based on spatial grid attributes to obtain the reverse migration trajectory of the target pollutants at different historical time nodes. The risk prediction and locking module is used to obtain the actual source location of the geological environment leakage based on the reverse migration trajectory of the target pollutant, and to perform geological environment risk prediction based on the actual source location.