Method for comprehensive assessment of regional groundwater environmental quality and risk of multi-source pollution factors

Abstract

The present application relates to the technical field of environmental monitoring, in particular to a multi-source pollution factor regional groundwater environmental quality and risk comprehensive evaluation method, comprising the following steps: quantifying groundwater pollution intensity based on adsorption saturation mechanism, using isotope fractionation correction to analyze multi-source weight, decomposing the item concentration and combining the pollution intensity to construct the risk structure tensor, extracting the dominant diffusion path and anisotropy index, and generating the groundwater environmental quality evaluation result, in the present application, the medium buffer state is determined by calculating the soil adsorption saturation degree index, the surface input total amount and the actual leaching penetration amount are dynamically distinguished, the real pollution intensity is quantified, the initial fingerprint is restored and the multi-source contribution weight is analyzed by using isotope fractionation reverse vector projection correction, the risk structure tensor matrix is constructed by combining the item concentration gradient and the toxicity factor, the dominant diffusion path and the anisotropy index are extracted by characteristic decomposition, and the transformation from scalar numerical evaluation to three-dimensional vector space risk structure is realized.

Description

Technical Field

[0001] This invention relates to the field of environmental monitoring technology, and in particular to a method for comprehensive assessment of regional groundwater environmental quality and risk from multiple sources of pollution. Background Technology

[0002] The field of environmental monitoring technology involves the determination of the concentration, distribution, and changing trends of pollutants in environmental media. Its core aspects encompass the monitoring of physical, chemical, and biological indicators of various environmental elements, including water, soil, and atmosphere. By employing sensor networks, laboratory analysis, and remote sensing, environmental data is systematically collected and analyzed to ascertain environmental quality, pollution source distribution, and environmental change patterns. This provides scientific basis and data support for environmental protection planning, pollution prevention and control decisions, and environmental management. Traditional regional groundwater environmental quality and risk comprehensive assessment methods refer to the evaluation process of groundwater quality and potential hazards within a specific area. These methods typically employ single-factor evaluation or Nemerow index methods. Groundwater samples are collected, and laboratory chemical analysis is used to obtain concentration data for each individual indicator. The measured concentration values ​​are directly compared with the limits in national groundwater quality standards. Based on the exceedance of standards, the water quality category is determined. Simultaneously, health risk assessment models are used to calculate carcinogenic or non-carcinogenic risk values ​​based on pollutant concentrations, human exposure parameters, and toxicity parameters to determine the safety of the groundwater environment.

[0003] Traditional methods rely on static monitoring of concentrations and comparison with standard limits, ignoring the adsorption and buffering effect of soil media on pollutants and the historical cumulative effect. This makes it difficult to distinguish the difference between the current actual leaching load and the total surface input load, leading to misjudgments of the intensity of pollution at the source. Risk calculations based on measured concentrations cannot effectively eliminate the interference of biodegradation on isotopic fingerprints, resulting in distorted multi-source apportionment results. Furthermore, the lack of vectorized analysis of the spatial diffusion directionality of pollutants and the anisotropic characteristics of multi-source superposition makes it difficult for the assessment results to characterize the actual migration trend and structural distribution pattern of risks in three-dimensional geological space. Summary of the Invention

[0004] To achieve the above objectives, the present invention adopts the following technical solution: a comprehensive assessment method for regional groundwater environmental quality and risk of multiple pollutants, comprising the following steps: S1: Obtain land use duration data and soil maximum adsorption capacity parameter, collect current soil pollutant background concentration data, calculate cumulative input total and deduct degradation loss, compare the deducted total with the soil maximum adsorption capacity parameter, and generate adsorption saturation index. S2: Compare the adsorption saturation index with the preset saturation limit standard, calculate the dissolved load as the unsaturated pollution production, select all load data as the saturated pollution production, analyze the degree of influence of groundwater on surface pollution source input, and generate groundwater pollution production intensity parameters. S3: Analyze the dual isotope ratio data to establish the distribution slope, match the distribution slope with the theoretical fractionation characteristics, perform reverse vector projection correction to restore the initial isotope characteristics, analyze the mixing relationship between the initial isotope characteristics and the end data, and generate the pollution source contribution ratio weight. S4: Based on the pollution source contribution ratio weight, decompose the total concentration data to obtain the component concentration data, calculate the spatial change rate of the component concentration data to establish a concentration gradient vector, combine the groundwater pollution generation intensity parameter and the toxicity impact factor corresponding to the pollutant type to weight the concentration gradient vector, perform the outer product operation, and construct the risk structure tensor matrix. S5: Perform eigenvalue decomposition on the risk structure tensor matrix, extract principal component vectors, identify the dominant diffusion path based on the principal component vectors and calculate the risk diffusion anisotropy index, integrate the dominant diffusion path and the risk diffusion anisotropy index, and generate groundwater environmental quality assessment results.

[0005] As a further embodiment of the present invention, the adsorption saturation index includes the cumulative load ratio and the remaining buffer capacity of the medium; the groundwater pollution generation intensity parameter includes the leaching flux, the penetration load, and the input response coefficient; the pollution source contribution ratio weight includes the contribution rate of industrial sources, the contribution rate of agricultural sources, and the contribution rate of domestic sources; the risk structure tensor matrix includes the spatial gradient cross-correlation component, the risk intensity weighted modulus, and the local diffusion potential energy tensor element; and the groundwater environmental quality assessment result includes the main axis of pollution plume migration, regional risk level zoning, and the coordinates of priority control areas.

[0006] As a further aspect of the present invention, the step of obtaining the adsorption saturation index specifically comprises: S101: Obtain the land use duration data and the current background concentration data of pollutants in the soil from the monitoring unit, call the unit time input load parameter of the corresponding monitoring unit attribute, perform multiplication operation on the land use duration data and the unit time input load parameter, sum with the background concentration data of pollutants, calculate the total pollutant load value received by the monitoring unit, and generate the theoretical cumulative input total of pollutants. S102: Based on the total theoretical cumulative input of pollutants, obtain the natural degradation loss coefficient of the soil medium, calculate the degradation loss value based on the land use duration data, calculate the net residual value of pollutants after environmental self-purification loss correction, determine the actual retention load data in the soil medium, and generate the cumulative soil retention. S103: Call the accumulated soil retention and the maximum soil adsorption capacity parameter of the monitoring unit, calculate the ratio of the current retention load value to the maximum adsorption capacity limit of the medium, analyze the retention status of pollutants in the soil medium, including the degree of filling of the soil medium with various pollutants and the remaining buffer space, and generate an adsorption saturation index.

[0007] As a further aspect of the present invention, the steps for obtaining the groundwater pollution intensity parameter are as follows: S201: Call the soil medium adsorption characteristic parameters and pollutant chemical property data of the monitoring unit to establish the equilibrium distribution relationship of pollutant molecules in the dynamic mass exchange between the soil solid phase particle surface and the soil solution liquid phase. Calculate the ratio between the mass share of pollutants allocated to the solid phase medium and the mass share dissolved in the liquid phase water flow under the current environmental conditions. Quantify the migration activity and medium retention capacity of pollutants in the unsaturated state and generate the solid-liquid phase mass distribution coefficient. S202: Based on the solid-liquid phase mass distribution coefficient, obtain the surface pollutant input load data of the monitoring unit. For the transmission scenario where the soil medium has not yet reached the upper limit of adsorption capacity, perform calculations on the input load data and distribution coefficient, remove the components that are adsorbed and fixed by the soil medium, calculate the number of free components that infiltrate into the aquifer with gravity water, and generate dissolved leaching load data. S203: Obtain a preset saturation limit standard, compare the adsorption saturation index with the saturation limit standard, determine the buffering capacity of the soil medium, extract dissolved leaching load data as unsaturated pollution production when the index is lower than the standard, determine that the medium adsorption is ineffective and select all load data input from the surface as pollution production, determine the pollutant flux entering the groundwater, and generate groundwater pollution production intensity parameters.

[0008] As a further aspect of the present invention, the process of establishing the equilibrium distribution relationship of pollutant molecules on the surface of soil solid particles and the liquid phase of soil solution through dynamic mass exchange specifically includes: Soil media samples were collected within the monitoring unit area, and the percentage data of total organic carbon content, clay mineral component content data, and pH value of soil pore solution were measured. The system calls upon a pre-built pollution physicochemical fingerprint database to match the standard values ​​of organic carbon-water partition coefficients and molecular polar charge characteristic parameters corresponding to the types of pollutants to be evaluated. For organic pollutants with hydrophobic structures, based on the principle of linear distribution, the percentage data of total organic carbon content in the soil is multiplied by the standard value of the organic carbon-water partition coefficient to quantify the nonpolar adsorption potential energy parameter of soil organic matter relative to organic molecules. For inorganic pollutants in ionic form, a surface charge balance model is constructed based on ion exchange and surface complexation mechanisms, combined with the clay mineral component content data and pH values, to calculate the electrostatic adsorption intensity of ions on the soil colloid surface. Based on the molecular polar charge characteristic parameters, the dominant adsorption mechanism is identified. The nonpolar adsorption potential energy parameters and the electrostatic adsorption strength values ​​are integrated to construct an isothermal adsorption equation describing the bidirectional migration of solute molecules between the solid and liquid interfaces. Analyze the slope of the tangent line of the isothermal adsorption equation in the low concentration range to determine the constant ratio between the amount of solid phase adsorbed per unit mass and the concentration of liquid phase dissolved per unit volume.

[0009] As a further aspect of the present invention, the step of obtaining the pollution source contribution ratio weight specifically includes: S301: Acquire dual isotope ratio data of groundwater monitoring well samples, map them to a two-dimensional scatter coordinate system, calculate the spatial distribution fitting slope of the measured data points, call the preset theoretical fractionation slope interval of biochemical reaction, perform matching verification between the spatial distribution fitting slope and the theoretical fractionation slope interval, identify the drift trend of isotope signals under microbial action, determine the geometric correction direction and intensity of fingerprint backtracking, and generate fractionation correction vector parameters; S302: Based on the fractionation correction vector parameters, perform reverse geometric projection operation on isotope data points with offset characteristics, translate the coordinate position in the opposite direction of the fractionation evolution path until the projection point falls on the boundary line of the mixing area surrounded by the source metadata, extract the corrected coordinate values, and generate initial isotope fingerprint data. S303: Call the initial isotope fingerprint data and the end-member feature data of each pollution source in the region to construct a multi-source linear mixing equation system based on the principle of isotope mass conservation, analyze the mathematical composition relationship of each end-member component in the mixed sample, calculate the relative supply share of each pollution source to the groundwater sample, and generate the pollution source contribution ratio weight.

[0010] As a further aspect of the present invention, the process of calling a preset theoretical fractionation slope range of biochemical reactions, performing a matching verification between the spatial distribution fitting slope and the theoretical fractionation slope range, identifying the drift trend of isotope signals under microbial action, and determining the geometric correction direction and intensity of fingerprint tracing specifically involves: Query the dual-element isotope enrichment factor data of the pollutants to be evaluated under different electron acceptor environments for microbial degradation reactions, calculate the ratio of the two isotope enrichment factors, and construct a set of theoretical fractionation evolution slopes to characterize the degradation pathway. The least squares linear regression operation was performed on the measured dual isotope ratio data in a two-dimensional scatter coordinate system to obtain the statistical fitting slope and the corresponding slope confidence interval range that reflect the spatial distribution characteristics of the sample points. The statistical fitting slope is compared with the set of theoretical fractionation evolution slope values. When the statistical fitting slope is within the slope range of a specific degradation path, it is determined that the monitored water sample has undergone the corresponding type of biochemical fractionation. Based on the identified fractionation type, an evolution trend line is established in which the isotope ratio gradually increases with the reaction process. The vector opposite to the direction of the evolution trend line is defined as the geometric correction direction for restoring the initial characteristics. In a two-dimensional coordinate system, a polygonal mixed region is constructed, which is enclosed by metadata points of each potential pollution source. A projection ray is established starting from the measured data point and extending along the geometric correction direction. The intersection point of the projection ray and the boundary line of the polygonal mixed region is calculated. The Euclidean distance from the measured data point to the intersection point is measured and quantified into the geometric correction intensity for performing back projection correction.

[0011] As a further aspect of the present invention, the step of obtaining the risk structure tensor matrix specifically includes: S401: Obtain the measured total pollutant concentration data at the monitoring point, call the pollution source contribution ratio weight, perform a multiplication decomposition operation on the total concentration data of the mixed system based on the source weight, separate the concentration components belonging to each emission source, quantify the absolute contribution share of the pollution source at the current spatial location, and generate sub-pollutant concentration data. S402: Based on the pollutant concentration data, obtain the three-dimensional geographic coordinate information of the monitoring network, calculate the rate of change of the concentration field in each orthogonal direction by performing partial differential operations for spatial location, analyze the diffusion trend and migration front of each pollution plume in the aquifer, construct a vector set composed of directional components and modulus, establish a physical field description object of the transport characteristics of pollutants in the three-dimensional geological medium, and generate a three-dimensional concentration gradient vector. S403: Call the groundwater pollution intensity parameter and the toxicological factors corresponding to each pollutant type to construct a comprehensive weighted coefficient characterizing the source load and degree of harm, perform weighted operation on the three-dimensional concentration gradient vector, perform outer product expansion operation on the modulated vector group, construct a high-order matrix object of multi-source risk intensity information and spatial directionality characteristics, and generate a risk structure tensor matrix.

[0012] As a further aspect of the present invention, the steps for obtaining the groundwater environmental quality assessment results are specifically as follows: S501: Call the risk structure tensor matrix, perform linear algebraic eigenvalue decomposition on the numerical array inside the matrix, deconstruct the high-dimensional risk superposition effect data into independent projection components on the orthogonal coordinate system, extract the feature vector group representing the main extended dimensions of the risk field in three-dimensional space, calculate the eigenvalues ​​of each feature vector and quantize the vector modulus, and generate principal component feature vector parameters. S502: Based on the principal component feature vector parameters, sort and compare the modulus values ​​corresponding to each feature vector, select the feature vector with the largest weight of the modulus value as the principal axis direction, identify the dominant diffusion path of pollutants in the aquifer medium, calculate the risk diffusion anisotropy index that characterizes the directional intensity of risk diffusion, and generate risk diffusion spatial feature data. S503: Call the aforementioned risk diffusion spatial feature data, extract the dominant diffusion path vector and anisotropy index values, perform spatial mapping and overlay on the feature data and the regional hydrogeological background layer, divide the spatial distribution level of risk level gradient according to the anisotropy intensity values, and generate groundwater environmental quality assessment results.

[0013] As a further aspect of the present invention, the main extended dimension refers to the physical direction indicated by the mutually orthogonal eigenvectors obtained after eigenvalue decomposition of the risk structure tensor matrix in three-dimensional space, representing the independent physical migration trajectory of risk in three-dimensional space along the dominant diffusion direction, the lateral discrete direction, and the vertical settlement direction.

[0014] Compared with the prior art, the advantages and positive effects of the present invention are as follows: In this invention, the buffer state of the medium is determined by calculating the soil adsorption saturation index, dynamically distinguishing between the total amount of surface input and the actual leaching penetration, quantifying the true pollution intensity, correcting with isotope fractionation and reverse vector projection, restoring the initial fingerprint and analyzing the contribution weights of multiple sources, constructing a risk structure tensor matrix by combining the concentration gradient of sub-items and toxic factors, and extracting the dominant diffusion path and anisotropic indicators by feature decomposition, thus realizing the transformation from scalar numerical evaluation to three-dimensional vector space risk structure. Attached Figure Description

[0015] To more clearly illustrate the technical solutions in the embodiments of the present invention, the accompanying drawings used in the description of the embodiments will be briefly introduced below. Obviously, the accompanying drawings described below are only some embodiments of the present invention. For those skilled in the art, other drawings can be obtained based on these drawings without creative effort.

[0016] Figure 1 This is a schematic diagram of the steps of the present invention; Figure 2 This is a detailed schematic diagram of S1 of the present invention; Figure 3 This is a detailed schematic diagram of S2 of the present invention; Figure 4 This is a detailed schematic diagram of S3 of the present invention; Figure 5 This is a detailed schematic diagram of S4 of the present invention; Figure 6 This is a detailed schematic diagram of S5 of the present invention. Detailed Implementation

[0017] The technical solution of the present invention will now be described with reference to the accompanying drawings.

[0018] In embodiments of the present invention, words such as "exemplarily," "for example," etc., are used to indicate that something is an example, illustration, or description. Any embodiment or design described as "exemplary" in the present invention should not be construed as being more preferred or advantageous than other embodiments or designs. Specifically, the use of the word "exemplary" is intended to present the concept in a concrete manner. Furthermore, in embodiments of the present invention, the meaning expressed by "and / or" can be both, or either one.

[0019] In the embodiments of this invention, the terms "image" and "picture" may sometimes be used interchangeably. It should be noted that, without emphasizing the distinction between them, they convey the same meaning. Similarly, the terms "of," "corresponding (relevant)," and "corresponding" may sometimes be used interchangeably. It should be noted that, without emphasizing the distinction between them, they convey the same meaning.

[0020] In this embodiment of the invention, sometimes a subscript such as W1 may be written in a non-subscript form such as W1. When the difference is not emphasized, the meaning they express is the same.

[0021] To make the technical problems, technical solutions and advantages of the present invention clearer, a detailed description will be given below in conjunction with the accompanying drawings and specific embodiments.

[0022] Please see Figure 1 This invention provides a method for comprehensive assessment of regional groundwater environmental quality and risk from multiple sources of pollution, comprising the following steps: S1: Obtain land use duration data and soil maximum adsorption capacity parameter, collect current soil pollutant background concentration data, calculate cumulative input total and deduct degradation loss, compare the deducted total with soil maximum adsorption capacity parameter, and generate adsorption saturation index. S2: Compare the adsorption saturation index with the preset saturation limit standard, calculate the dissolved load as the unsaturated pollution production, select all load data as the saturated pollution production, analyze the degree of influence of groundwater on surface pollution source input, and generate groundwater pollution production intensity parameters. S3: Analyze the dual isotope ratio data to establish the distribution slope, match the distribution slope with the theoretical fractionation characteristics, perform reverse vector projection correction to restore the initial isotope characteristics, analyze the mixing relationship between the initial isotope characteristics and the end data, and generate the pollution source contribution ratio weight. S4: Based on the contribution ratio weight of pollution sources, decompose the total concentration data to obtain the component concentration data, calculate the spatial change rate of the component concentration data to establish the concentration gradient vector, combine the groundwater pollution generation intensity parameter and the toxicity impact factor corresponding to the pollutant type to weight the concentration gradient vector, perform the outer product operation, and construct the risk structure tensor matrix. S5: Perform eigenvalue decomposition on the risk structure tensor matrix, extract principal component vectors, identify the dominant diffusion path based on the principal component vectors, calculate the risk diffusion anisotropy index, integrate the dominant diffusion path and the risk diffusion anisotropy index, and generate groundwater environmental quality assessment results.

[0023] The adsorption saturation index includes the cumulative load ratio and the remaining buffer capacity of the medium; the groundwater pollution generation intensity parameters include the leaching flux, the breakthrough load, and the input response coefficient; the pollution source contribution ratio weights include the contribution rates of industrial sources, agricultural sources, and domestic sources; the risk structure tensor matrix includes the spatial gradient cross-correlation components, the risk intensity weighted modulus, and the local diffusion potential energy tensor element; and the groundwater environmental quality assessment results include the main axis of pollution plume migration, regional risk level zoning, and the coordinates of priority control areas.

[0024] Please see Figure 2 The specific steps for obtaining the adsorption saturation index are as follows: S101: Obtain the land use duration data and the current background concentration data of pollutants in the soil from the monitoring unit, call the unit time input load parameter of the corresponding monitoring unit attribute, perform multiplication operation on the land use duration data and the unit time input load parameter, sum with the background concentration data of pollutants, calculate the total pollutant load value received by the monitoring unit, and generate the theoretical cumulative input total of pollutants. In the specific implementation process, the system retrieves land use duration data for the monitoring unit plot by accessing the regional historical archive database. This data is determined to be the specific number of years from the start of operation of the chemical industrial park to the assessment benchmark date, set to 25 years in this embodiment. Simultaneously, the system collects surface soil samples at a depth of 0-20 cm from soil sampling points deployed in the vadose zone of the monitoring unit. These samples are then sent to a laboratory for gas chromatography-mass spectrometry (GC-MS) to determine the current background concentration of pollutants in the soil. The initial background concentration of trichloroethylene was measured to be 0.5 mg / kg. The system then accesses a pre-set source inventory database and matches the corresponding unit-time input load parameter based on the plot's "chemical production area" attribute. This parameter represents the average mass of pollutants entering the surface soil per unit time due to production leaks or atmospheric deposition, set to 20 kg / year. Subsequently, the system performs a multiplication operation, multiplying the land use duration data (25 years) by the unit-time input load parameter (20 kg / year) to obtain a total new input of 500 kg during the operation period. Next, the system sums the newly added total input with the initial stock calculated based on background concentration data (assuming a total surface soil mass of 100,000 kg, then the initial stock is 0.5 mg / kg × 100,000 kg = 0.05 kg, a small value that can be ignored or included). Finally, the system calculates the total pollutant load received by the monitoring unit during the assessment period, i.e., the theoretical cumulative total input of pollutants is 500.05 kg. This value represents the theoretical upper limit of the total mass of pollutants that should exist in the environmental media of this area, without considering any migration or transformation.

[0025] S102: Based on the total cumulative input of pollutants, obtain the natural degradation loss coefficient of soil medium, calculate the degradation loss value based on land use duration data, calculate the net residual value of pollutants after environmental self-purification loss correction, determine the actual retention load data in soil medium, and generate the cumulative soil retention. In the specific implementation process, based on the hydrogeological survey report of the monitoring unit, the natural degradation loss coefficient of the target pollutant (trichloroethylene) in the soil medium is obtained. This coefficient is usually determined based on a first-order kinetic decay model, characterizing the average annual decay rate of the pollutant under soil microbial action and hydrolysis. In the example, the coefficient was determined to be 0.03 based on field pilot-scale experiments (i.e., 3% natural loss per year). Based on the aforementioned land use duration data (25 years), the system calculates the degradation loss value using an exponential decay algorithm. Specifically, using the theoretical cumulative input total amount of the generated pollutant (500.05 kg) as the base, the decay formula is applied to calculate the remaining amount after 25 years, or the annual cumulative decay value of the input amount is calculated year by year. For simplification, a comprehensive decay calculation is used, and the total biodegradation and volatilization loss over 25 years is calculated to be approximately 260 kg. Subsequently, the system performs a subtraction operation, subtracting the degradation loss value from the theoretical cumulative input total amount of the pollutant to calculate the net remaining pollutant value after environmental self-purification loss correction, i.e., 500.05 kg - 260 kg = 240.05 kg. This value is determined to be the actual retention load data in the soil medium, generating the cumulative soil retention. This process eliminates the share of pollutants that disappear due to the environment's self-purification capacity, ensuring that subsequent assessments of soil adsorption pressure are based on the actual existing stock, rather than theoretical input, thereby improving the accuracy of saturation calculations.

[0026] S103: Call the soil cumulative retention and the soil maximum adsorption capacity parameter of the monitoring unit, calculate the ratio coefficient of the current retention load value to the maximum adsorption capacity limit of the medium, analyze the retention status of pollutants in the soil medium, including the filling degree of the soil medium for various pollutants and the remaining buffer space, and generate the adsorption saturation index. In the specific implementation process, the maximum adsorption capacity parameter of the soil in the monitoring unit is invoked. This parameter is the saturated adsorption capacity measured by Langmuir isotherm adsorption experiments on the collected soil column, characterizing the upper limit of the pollutant mass that a specific volume of soil in the monitoring unit can adsorb and fix. Assuming the total maximum adsorption capacity of the vadose zone soil in the monitoring unit is determined to be 1200 kg, the system uses the cumulative soil retention amount generated by S102 (240.05 kg) as the numerator and the maximum adsorption capacity parameter of the soil (1200 kg) as the denominator, performing a division operation to calculate the proportionality coefficient of the current retention load value relative to the maximum adsorption capacity limit of the medium. The calculation result is 240.05 / 1200≈0.20. This value directly reflects the retention state of pollutants in the soil medium, indicating that the current soil medium's filling degree for pollutants is 20%, meaning that the soil still retains 80% of the remaining buffer space for intercepting newly introduced pollutants. The system generates this proportionality coefficient of 0.20 as the adsorption saturation index. If the index is close to 1.0, it indicates that the soil adsorption sites are full and the medium has lost its anti-fouling barrier function; while 0.20 in this embodiment indicates that the soil is in a good unsaturated state and has a strong blocking ability.

[0027] Please see Figure 3 The specific steps for obtaining groundwater pollution intensity parameters are as follows: S201: Call the soil medium adsorption characteristic parameters and pollutant chemical property data of the monitoring unit to establish the equilibrium distribution relationship of pollutant molecules in the dynamic mass exchange between the soil solid phase particle surface and the soil solution liquid phase. Calculate the ratio between the mass share of pollutants allocated to the solid phase medium and the mass share dissolved in the liquid phase water flow under the current environmental conditions. Quantify the migration activity and medium retention capacity of pollutants in the unsaturated state and generate the solid-liquid phase mass distribution coefficient. First, soil samples were collected from the monitoring unit, and the laboratory test results are shown in Table 1. The system then invoked a pre-set pollution physicochemical fingerprint database and matched it with the standard value of the organic carbon-water partition coefficient of trichloroethylene (TCE). The concentration of total organic carbon (TCE) is 126 L / kg, and its molecular polar charge characteristic parameters indicate that it is a nonpolar molecule. For TCEs, organic pollutants with hydrophobic structures, the system performs a multiplication operation based on the linear partition principle, multiplying the percentage data of total organic carbon content in the soil (1.5%) with the standard value of the organic carbon-water partition coefficient (126 L / kg). L / kg is used to quantify the nonpolar adsorption potential energy parameter of soil organic matter relative to organic molecules. For potentially coexisting heavy metal ions (such as cadmium), the system constructs a surface charge balance model based on the principle of ion exchange, combined with clay mineral component content data (20%) and pH value (pH 6.8), to calculate the electrostatic adsorption strength value (e.g., calculated to be 0.5 L / kg). In this embodiment, based on the molecular polar charge characteristic parameters, the dominant adsorption mechanism of TCE is identified as hydrophobic partitioning. Therefore, the system mainly integrates nonpolar adsorption potential energy parameters to construct an isothermal adsorption equation describing the bidirectional migration of solute molecules between the solid and liquid interfaces. The system analyzes the tangent slope of this equation in the low concentration range, ultimately determining a constant ratio between the unit mass of solid phase adsorption and the unit volume of liquid phase dissolution concentration to be 1.89 L / kg, generating a solid-liquid phase mass partition coefficient ( ). ).

[0028] Table 1 Physicochemical parameters of soil media samples from the monitoring unit

[0029] As shown in Table 1, the measured soil organic carbon content and clay content are the basic input data for calculating the solid-liquid partition coefficient.

[0030] S202: Based on the solid-liquid phase mass distribution coefficient, obtain the surface pollutant input load data of the monitoring unit. For the transmission scenario where the soil medium has not yet reached the upper limit of adsorption capacity, calculate the input load data and distribution coefficient, remove the components that are adsorbed and fixed by the soil medium, calculate the number of free components that infiltrate into the aquifer with gravity water, and generate dissolved leaching load data. In the specific implementation process, the system acquires the surface pollutant input load data of the monitoring unit for the year and sets it as the newly added input amount of 20 kg for that year. For transport scenarios where the soil medium has not yet reached its adsorption capacity limit (based on the judgment in S103), the system utilizes the solid-liquid phase mass distribution coefficient generated in S201 (…). L / kg) and soil physical parameters (bulk density) g / cm³, porosity The system performs the calculations. First, it calculates the retardation factor. This hindering factor indicates that pollutants move in soil approximately 9.91 times slower than in water, and also means that, under equilibrium conditions, most pollutants are adsorbed. The system uses a distribution formula to calculate the proportion dissolved in the liquid phase. ,Right now The system performs a multiplication operation between the input load data (20 kg) and the specified ratio, removes the components adsorbed and fixed by the soil medium (approximately 18 kg), and calculates the amount of free components that infiltrate into the aquifer with gravity water. kg. This value represents the actual penetration amount under conditions where the soil has adsorption capacity, and the system generates it as dissolved leaching load data.

[0031] S203: Obtain the preset saturation limit standard, compare the adsorption saturation index with the saturation limit standard, determine the buffering capacity of the soil medium, extract the dissolved leaching load data as the unsaturated pollution production when the index is lower than the standard, and determine the pollution production when the index is not lower than the standard. In the specific implementation process, the system obtains the preset saturation limit standard from the evaluation standard database. This standard is usually set at 0.80 (i.e., 80% saturation) as the warning threshold for soil adsorption failure. The system compares the generated adsorption saturation index (0.20) with this saturation limit standard (0.80). In this embodiment, 0.20 < 0.80, and the result is that the buffering capacity of the soil medium is good and in the "unsaturated range". Therefore, the system extracts the dissolved leaching load data (2.0 kg / year) generated by S202 as the pollution generation. If, under another high pollution scenario, the calculated adsorption saturation index is 0.95, then 0.95 > 0.80, and the system will determine that the medium adsorption has failed. At this time, the soil no longer retains pollutants, but may become a secondary pollution source due to desorption. The system will directly select all the load data input from the surface (20 kg / year) as the pollution generation. In this example, the system ultimately determines that the pollutant flux entering the groundwater is 2.0 kg / year, and normalizes it to a flux value per unit area to generate groundwater pollution intensity parameters.

[0032] Please see Figure 4 The specific steps for obtaining the pollution source contribution ratio weights are as follows: S301: Acquire dual isotope ratio data of groundwater monitoring well samples, map them to a two-dimensional scatter coordinate system, calculate the spatial distribution fitting slope of the measured data points, call the preset theoretical fractionation slope interval of biochemical reaction, perform matching verification between the spatial distribution fitting slope and the theoretical fractionation slope interval, identify the drift trend of isotope signals under microbial action, determine the geometric correction direction and intensity of fingerprint backtracking, and generate fractionation correction vector parameters; The system collects and monitors well water samples and measures carbon isotopes ( ) and chlorine isotopes ( The ratio data of ) . Actual measured data shows that the isotopic fingerprint of a certain monitoring well is ( ) , The system maps dual isotope ratio data from multiple monitoring wells in the area to a two-dimensional scattered coordinate system. The slope of the spatial distribution fitting of the measured data points is calculated using the least squares method to obtain the slope value. The system calls the preset theoretical fractionation slope range for biochemical reactions and finds that the theoretical fractionation slope range for trichloroethylene under anaerobic conditions is [2.2, 2.8], while the slope range under aerobic oxidation conditions is [0.5, 1.2]. The system performs a matching check between the spatial distribution fitting slope (2.5) and the theoretical fractionation slope range, finding that it falls within the reductive dechlorination range. This indicates that the isotope signal drift trend under microbial action is moving upwards and to the right along a straight line with a slope of 2.5 (the direction of heavy isotope enrichment). Based on this, the system determines that the geometric correction direction for fingerprint tracing is the opposite direction of this trend line (i.e., downwards and to the left), and the correction strength depends on the distance of the measured point from the source region. The system packages this reverse unit vector and slope characteristics to generate fractionation correction vector parameters.

[0033] S302: Based on the fractionation correction vector parameters, perform reverse geometric projection operation on isotope data points with offset characteristics, translate the coordinate position in the opposite direction of the fractionation evolution path until the projection point falls on the boundary line of the mixing area surrounded by the source metadata, extract the corrected coordinate values, and generate the initial isotope fingerprint data. The system constructs a polygonal hybrid region in a two-dimensional coordinate system, bounded by metadata points from potential pollution sources. It is assumed that two potential pollution sources exist within this region: Source A (…). , ) and source B ( , The line connecting these two points forms the boundary line of the mixing region. The current measured data point P (-25‰, +2‰) is clearly outside this line, indicating that fractionation has occurred. Based on the fractionation correction vector parameters (slope 2.5, direction downward to the left) generated by S301, the system establishes the equation of a straight line passing through point P with a slope of 2.5: The system calculates the coordinates of the intersection point of the straight line and the line connecting source A and source B. After inverse geometric projection, the coordinate position is translated in the opposite direction of the fractionation evolution path, and the coordinates of the projected intersection point P' are calculated to be (-29‰, -0.5‰). This point P' falls exactly on the line connecting source A and source B, representing the original mixed state of the pollutants before biodegradation. The system extracts this corrected coordinate value (-29‰, -0.5‰) to generate initial isotopic fingerprint data, eliminating fingerprint bias caused by biodegradation.

[0034] S303: Call the initial isotopic fingerprint data and the end-member feature data of each pollution source in the region to construct a multi-source linear mixing equation system based on the principle of isotopic mass conservation, analyze the mathematical composition relationship of each end-member component in the mixed sample, calculate the relative supply share of each pollution source to the groundwater sample, and generate the pollution source contribution ratio weight. The system calls upon the generated initial isotopic fingerprint data P' (-29‰, -0.5‰) and the endmember feature data of source A (-30‰, -3‰) and source B (-28‰, +1‰). A multi-source linear mixing equation system is constructed based on the principle of isotopic mass conservation. Let the contribution ratio of source A be... The contribution ratio of source B is Then it satisfies: 1. 2. (Conservation of carbon isotopes); 3. (Chlorine isotope conservation); the system uses a matrix solution method to analyze the mathematical compositional relationships of each endmember component in the mixed sample. Substituting into numerical calculations: from Equation 2, we can obtain... Substitute into equation 3 to verify. (Considering experimental error, the results are within acceptable limits). The final calculation shows that the relative contribution of each pollution source to the groundwater sample is: Source A contributes 50%, and Source B contributes 50%. The system will display this result ( This is used to generate a weighted proportion of the pollution source's contribution, which is then used for subsequent concentration decomposition.

[0035] Please see Figure 5 The specific steps for obtaining the risk structure tensor matrix are as follows: S401: Obtain the measured total pollutant concentration data at the monitoring point, call the pollution source contribution ratio weight, perform a multiplicative decomposition operation on the total concentration data of the mixed system based on the source weight, separate the concentration components belonging to each emission source, quantify the absolute contribution share of the pollution source at the current spatial location, and generate sub-pollutant concentration data. Obtain the measured total pollutant concentration data for the monitoring point during the current sampling period, such as the total trichloroethylene concentration at that point. mg / L. The pollution source contribution weights (source A weights) generated by system call S303. Source B weight The system performs a source-weighted multiplicative decomposition operation on the total concentration data of the mixed system to calculate the concentration of the component attributed to source A. mg / L; Calculate the concentration of the component attributed to source B. mg / L. Through this step, the system successfully separated the independent concentration components belonging to different potential emission sources, quantifying the absolute contribution share of a specific pollution source at the current spatial location. This process decouples a single mixed concentration field into two independent source-induced concentration fields, generating individual pollutant concentration data, providing fundamental data for subsequent assessments of the spatial diffusion risk of different pollution sources.

[0036] S402: Based on the concentration data of individual pollutants, obtain the three-dimensional geographic coordinate information of the monitoring network, calculate the rate of change of the concentration field in each orthogonal direction by performing partial differential operations for spatial location, analyze the diffusion trend and migration front of each pollution plume in the aquifer, construct a vector set composed of directional components and modulus, establish a physical field description object of the transport characteristics of pollutants in the three-dimensional geological medium, and generate a three-dimensional concentration gradient vector. Based on the three-dimensional geographic coordinates (x, y, z) of 15 monitoring wells in the monitoring network and the source A component pollutant concentration data generated by S401, a three-dimensional concentration field grid for source A is constructed using Kriging interpolation. At each grid node, the system performs partial differential operations with respect to spatial location to calculate the rate of change of the concentration field in the three orthogonal directions: x (east-west), y (north-south), and z (vertical). For example, at a certain coordinate point, the calculated... mg / L / m, mg / L / m, mg / L / m. This set of values ​​indicates that the pollutants mainly diffuse westward. The system analyzes this result and identifies that the diffusion trend of the pollution plume from source A in the aquifer points in the negative x-axis direction. The system constructs these three components as a vector set consisting of the directional component and the modulus. A physical field description object for the transport characteristics of pollutants in a three-dimensional geological medium is established, generating a three-dimensional concentration gradient vector. The same operation is repeated for source B to obtain the gradient vector of source B.

[0037] S403: Call the groundwater pollution generation intensity parameter and the toxicological factors corresponding to each type of pollutant to construct a comprehensive weighted coefficient characterizing the source load and degree of hazard, perform weighted operation on the three-dimensional concentration gradient vector, perform outer product expansion operation on the modulated vector group, construct a high-order matrix object of multi-source risk intensity information and spatial directionality characteristics, and generate a risk structure tensor matrix. In the specific implementation process, the system calls the generated groundwater pollution generation intensity parameters (assuming the intensity of source area A). kg / year / m²) and the toxicological factors of trichloroethylene (including carcinogenic slope factor). The system constructs a comprehensive weighted coefficient characterizing the source load and the degree of hazard. The system generates a three-dimensional concentration gradient vector for S402. Perform weighted calculations and adjust the vector. Subsequently, the system performs an outer product expansion operation on the modulated vector group, that is, calculates the tensor. The calculation process is as follows: The system superimposes the tensors calculated from source A and source B to construct a high-order matrix object capable of simultaneously carrying multi-source risk intensity information and spatial directionality characteristics, generating a risk structure tensor matrix. The diagonal elements of this matrix represent the risk intensity in each direction, while the off-diagonal elements represent risk shear correlation.

[0038] Please see Figure 6 The specific steps for obtaining groundwater environmental quality assessment results are as follows: S501: Call the risk structure tensor matrix, perform linear algebraic eigenvalue decomposition on the numerical array inside the matrix, deconstruct the high-dimensional risk superposition effect data into independent projected components on the orthogonal coordinate system, extract the feature vector group representing the main extended dimensions of the risk field in three-dimensional space, calculate the eigenvalues ​​of each feature vector and quantize the vector modulus, and generate principal component feature vector parameters. The system calls the superimposed risk structure tensor matrix (for simplicity, the single-source matrix result from S403 is used as an example; the actual result is a multi-source superposition matrix). Symmetric matrices undergo linear algebraic eigenvalue decomposition. This is achieved by solving the characteristic equation. Three eigenvalues ​​were calculated. , , (This is a simplified example; in reality, multiple sources will have secondary eigenvalues ​​after stacking.) The system calculates the eigenvector corresponding to each eigenvalue to obtain the corresponding largest eigenvalue. eigenvectors Its direction roughly points to [-0.5, 0.2, 0], which is the original gradient direction. This process deconstructs the high-dimensional risk superposition effect data into independent projected components on an orthogonal coordinate system. The system extracts feature vector sets representing the main extended dimensions of the risk field in three-dimensional space ( ) and the corresponding eigenvalues ​​( The quantized vector modulus is used to generate principal component eigenvector parameters. As shown in Table 2, the magnitude of the eigenvalues ​​reflects the diffusion intensity of risk in the corresponding eigenvector direction.

[0039] Table 2. Example of Risk Structure Tensor Feature Decomposition Results

[0040] As shown in Table 2, the first principal component has an absolute advantage, indicating that the risk has a strong directionality.

[0041] S502: Based on the principal component feature vector parameters, sort and compare the modulus values ​​corresponding to each feature vector, select the feature vector with the largest weight of the modulus value and use it as the principal axis direction, identify the dominant diffusion path of pollutants in the aquifer medium, calculate the risk diffusion anisotropy index that characterizes the directional intensity of risk diffusion, and generate risk diffusion spatial feature data. Based on the principal component eigenvector parameters in Table 2, the eigenvalues Perform sorting and comparison to confirm. The maximum value is selected. The system selects the feature vector whose modulus value has the largest weight. Using the direction of the axis as the principal axis, the dominant diffusion path of pollutants in the aquifer was identified as west-northwest. Subsequently, the system calculated the Anisotropy Index (AI), which characterizes the intensity of the risk diffusion directionality. The calculation formula was set as follows: Substitute the values ​​into the calculation: The results indicate that the risk field exhibits significant linear anisotropy (values ​​close to 1 represent linearity, and values ​​close to 0 represent spherical isotropy). The index of 0.68 reveals that the pollution plume presents a narrow, strip-like distribution, rather than spreading uniformly in all directions. The system integrates the dominant path vector with the anisotropy index to generate spatial characteristic data of risk diffusion.

[0042] S503: Call the spatial feature data of risk diffusion, extract the dominant diffusion path vector and anisotropy index values, perform spatial mapping and overlay on the feature data and the regional hydrogeological background layer, divide the spatial distribution level of risk level gradient according to the anisotropy intensity value, and generate groundwater environmental quality assessment results. The system retrieves spatial feature data on risk diffusion and extracts the dominant diffusion path vector (west-northwest) and the anisotropy index value (0.68). Within the GIS platform, this feature data is spatially mapped and overlaid with a regional hydrogeological background layer (including aquifer thickness and permeability coefficient distribution maps). Based on the anisotropy intensity values, the system classifies the spatial distribution levels of the risk grade gradient and sets anisotropy indices. The region is designated as a "high-risk, high-migration zone," and the indicators are at... The area between these points is designated as a "medium-risk diffusion zone." In this example, the system delineates a narrow, elongated zone of strong migration along the dominant diffusion path vector and marks it as a key control area. This result indicates that although the average concentration in the surrounding area may not be high, the groundwater environmental quality faces extremely high risks in this specific direction, necessitating the priority deployment of interception wells downstream of this path. The system ultimately outputs a report containing a risk zoning map, migration trajectory predictions, and key control recommendations, generating a groundwater environmental quality assessment result.

[0043] The above are merely specific embodiments of the present invention, but the scope of protection of the present invention is not limited thereto. Any variations or substitutions that can be easily conceived by those skilled in the art within the technical scope disclosed in the present invention should be included within the scope of protection of the present invention. Therefore, the scope of protection of the present invention should be determined by the scope of the claims.

Claims

1. A comprehensive assessment method for regional groundwater environmental quality and risk from multiple pollution sources, characterized in that, Includes the following steps: S1: Obtain land use duration data and soil maximum adsorption capacity parameter, collect current soil pollutant background concentration data, calculate cumulative input total and deduct degradation loss, compare the deducted total with the soil maximum adsorption capacity parameter, and generate adsorption saturation index. S2: Compare the adsorption saturation index with the preset saturation limit standard, calculate the dissolved load as the unsaturated pollution production, select all load data as the saturated pollution production, analyze the degree of influence of groundwater on surface pollution source input, and generate groundwater pollution production intensity parameters. S3: Analyze the dual isotope ratio data to establish the distribution slope, match the distribution slope with the theoretical fractionation characteristics, perform reverse vector projection correction to restore the initial isotope characteristics, analyze the mixing relationship between the initial isotope characteristics and the end data, and generate the pollution source contribution ratio weight. S4: Based on the pollution source contribution ratio weight, decompose the total concentration data to obtain the component concentration data, calculate the spatial change rate of the component concentration data to establish a concentration gradient vector, combine the groundwater pollution generation intensity parameter and the toxicity impact factor corresponding to the pollutant type to weight the concentration gradient vector, perform the outer product operation, and construct the risk structure tensor matrix. S5: Perform eigenvalue decomposition on the risk structure tensor matrix, extract principal component vectors, identify the dominant diffusion path based on the principal component vectors and calculate the risk diffusion anisotropy index, integrate the dominant diffusion path and the risk diffusion anisotropy index, and generate groundwater environmental quality assessment results.

2. The method for comprehensive assessment of regional groundwater environmental quality and risk based on multiple pollutants according to claim 1, characterized in that, The adsorption saturation index includes the cumulative load percentage and the remaining buffer capacity of the medium; the groundwater pollution generation intensity parameter includes the leaching flux, the penetration load, and the input response coefficient; the pollution source contribution ratio weight includes the contribution rate of industrial sources, the contribution rate of agricultural sources, and the contribution rate of domestic sources; the risk structure tensor matrix includes the spatial gradient cross-correlation component, the risk intensity weighted modulus, and the local diffusion potential energy tensor element; and the groundwater environmental quality assessment results include the main axis of pollution plume migration, regional risk level zoning, and the coordinates of priority control areas.

3. The method for comprehensive assessment of regional groundwater environmental quality and risk based on multiple pollutants according to claim 1, characterized in that, The steps for obtaining the adsorption saturation index are as follows: S101: Obtain the land use duration data and the current background concentration data of pollutants in the soil from the monitoring unit, call the unit time input load parameter of the corresponding monitoring unit attribute, perform multiplication operation on the land use duration data and the unit time input load parameter, and sum them with the background concentration data of pollutants to calculate the total pollutant load value received by the monitoring unit and generate the theoretical cumulative input total of pollutants. S102: Based on the total theoretical cumulative input of pollutants, obtain the natural degradation loss coefficient of the soil medium, calculate the degradation loss value based on the land use duration data, calculate the net residual value of pollutants after environmental self-purification loss correction, determine the actual retention load data in the soil medium, and generate the cumulative soil retention. S103: Call the accumulated soil retention and the maximum soil adsorption capacity parameter of the monitoring unit, calculate the ratio of the current retention load value to the maximum adsorption capacity limit of the medium, analyze the retention status of pollutants in the soil medium, including the degree of filling of the soil medium with various pollutants and the remaining buffer space, and generate an adsorption saturation index.

4. The method for comprehensive assessment of regional groundwater environmental quality and risk of multiple pollutants according to claim 3, characterized in that, The specific steps for obtaining the groundwater pollution intensity parameters are as follows: S201: Call the soil medium adsorption characteristic parameters and pollutant chemical property data of the monitoring unit to establish the equilibrium distribution relationship of pollutant molecules in the dynamic mass exchange between the soil solid phase particle surface and the soil solution liquid phase. Calculate the ratio between the mass share of pollutants allocated to the solid phase medium and the mass share dissolved in the liquid phase water flow under the current environmental conditions. Quantify the migration activity and medium retention capacity of pollutants in the unsaturated state and generate the solid-liquid phase mass distribution coefficient. S202: Based on the solid-liquid phase mass distribution coefficient, obtain the surface pollutant input load data of the monitoring unit. For the transmission scenario where the soil medium has not yet reached the upper limit of adsorption capacity, perform calculations on the input load data and distribution coefficient, remove the components that are adsorbed and fixed by the soil medium, calculate the number of free components that infiltrate into the aquifer with gravity water, and generate dissolved leaching load data. S203: Obtain a preset saturation limit standard, compare the adsorption saturation index with the saturation limit standard, determine the buffering capacity of the soil medium, extract dissolved leaching load data as unsaturated pollution production when the index is lower than the standard, determine that the medium adsorption is ineffective and select all load data input from the surface as pollution production, determine the pollutant flux entering the groundwater, and generate groundwater pollution production intensity parameters.

5. The method for comprehensive assessment of regional groundwater environmental quality and risk based on multiple pollutants according to claim 4, characterized in that, The process of establishing the equilibrium distribution relationship of pollutant molecules on the surface of soil solid particles and the liquid phase of soil solution through dynamic mass exchange is as follows: Soil media samples were collected within the monitoring unit area, and the percentage data of total organic carbon content, clay mineral component content data, and pH value of soil pore solution were measured. The system calls upon a pre-built pollution physicochemical fingerprint database to match the standard values ​​of organic carbon-water partition coefficients and molecular polar charge characteristic parameters corresponding to the types of pollutants to be evaluated. For organic pollutants with hydrophobic structures, based on the principle of linear distribution, the percentage data of total organic carbon content in the soil is multiplied by the standard value of the organic carbon-water partition coefficient to quantify the nonpolar adsorption potential energy parameter of soil organic matter relative to organic molecules. For inorganic pollutants in ionic form, a surface charge balance model is constructed based on ion exchange and surface complexation mechanisms, combined with the clay mineral component content data and pH values, to calculate the electrostatic adsorption intensity of ions on the soil colloid surface. Based on the molecular polar charge characteristic parameters, the dominant adsorption mechanism is identified. The nonpolar adsorption potential energy parameters and the electrostatic adsorption strength values ​​are integrated to construct an isothermal adsorption equation describing the bidirectional migration of solute molecules between the solid and liquid interfaces. Analyze the slope of the tangent line of the isothermal adsorption equation in the low concentration range to determine the constant ratio between the amount of solid phase adsorbed per unit mass and the concentration of liquid phase dissolved per unit volume.

6. The method for comprehensive assessment of regional groundwater environmental quality and risk of multiple pollutants according to claim 4, characterized in that, The specific steps for obtaining the weight of the pollution source contribution ratio are as follows: S301: Acquire dual isotope ratio data of groundwater monitoring well samples, map them to a two-dimensional scatter coordinate system, calculate the spatial distribution fitting slope of the measured data points, call the preset theoretical fractionation slope interval of biochemical reaction, perform matching verification between the spatial distribution fitting slope and the theoretical fractionation slope interval, identify the drift trend of isotope signals under microbial action, determine the geometric correction direction and intensity of fingerprint backtracking, and generate fractionation correction vector parameters; S302: Based on the fractionation correction vector parameters, perform reverse geometric projection operation on isotope data points with offset characteristics, translate the coordinate position in the opposite direction of the fractionation evolution path until the projection point falls on the boundary line of the mixing area surrounded by the source metadata, extract the corrected coordinate values, and generate initial isotope fingerprint data. S303: Call the initial isotope fingerprint data and the end-member feature data of each pollution source in the region to construct a multi-source linear mixing equation system based on the principle of isotope mass conservation, analyze the mathematical composition relationship of each end-member component in the mixed sample, calculate the relative supply share of each pollution source to the groundwater sample, and generate the pollution source contribution ratio weight.

7. The method for comprehensive assessment of regional groundwater environmental quality and risk of multiple pollutants according to claim 6, characterized in that, The process of calling the preset theoretical fractionation slope range of biochemical reactions, performing matching verification between the spatial distribution fitting slope and the theoretical fractionation slope range, identifying the drift trend of isotope signals under microbial action, and determining the geometric correction direction and intensity of fingerprint tracing is specifically as follows: Query the dual-element isotope enrichment factor data of the pollutants to be evaluated under different electron acceptor environments for microbial degradation reactions, calculate the ratio of the two isotope enrichment factors, and construct a set of theoretical fractionation evolution slopes to characterize the degradation pathway. The least squares linear regression operation was performed on the measured dual isotope ratio data in a two-dimensional scatter coordinate system to obtain the statistical fitting slope and the corresponding slope confidence interval range that reflect the spatial distribution characteristics of the sample points. The statistical fitting slope is compared with the set of theoretical fractionation evolution slope values. When the statistical fitting slope is within the slope range of a specific degradation path, it is determined that the monitored water sample has undergone the corresponding type of biochemical fractionation. Based on the identified fractionation type, an evolution trend line is established in which the isotope ratio gradually increases with the reaction process. The vector opposite to the direction of the evolution trend line is defined as the geometric correction direction for restoring the initial characteristics. In a two-dimensional coordinate system, a polygonal mixed region is constructed, which is enclosed by metadata points of each potential pollution source. A projection ray is established starting from the measured data point and extending along the geometric correction direction. The intersection point of the projection ray and the boundary line of the polygonal mixed region is calculated. The Euclidean distance from the measured data point to the intersection point is measured and quantified into the geometric correction intensity for performing back projection correction.

8. The method for comprehensive assessment of regional groundwater environmental quality and risk of multiple pollutants according to claim 6, characterized in that, The specific steps for obtaining the risk structure tensor matrix are as follows: S401: Obtain the measured total pollutant concentration data at the monitoring point, call the pollution source contribution ratio weight, perform a multiplication decomposition operation on the total concentration data of the mixed system based on the source weight, separate the concentration components belonging to each emission source, quantify the absolute contribution share of the pollution source at the current spatial location, and generate sub-pollutant concentration data. S402: Based on the pollutant concentration data, obtain the three-dimensional geographic coordinate information of the monitoring network, calculate the rate of change of the concentration field in each orthogonal direction by performing partial differential operations for spatial location, analyze the diffusion trend and migration front of each pollution plume in the aquifer, construct a vector set composed of directional components and modulus, establish a physical field description object of the transport characteristics of pollutants in the three-dimensional geological medium, and generate a three-dimensional concentration gradient vector. S403: Call the groundwater pollution intensity parameter and the toxicological factors corresponding to each pollutant type to construct a comprehensive weighted coefficient characterizing the source load and degree of harm, perform weighted operation on the three-dimensional concentration gradient vector, perform outer product expansion operation on the modulated vector group, construct a high-order matrix object of multi-source risk intensity information and spatial directionality characteristics, and generate a risk structure tensor matrix.

9. The method for comprehensive assessment of regional groundwater environmental quality and risk of multiple pollutants according to claim 8, characterized in that, The specific steps for obtaining the groundwater environmental quality assessment results are as follows: S501: Call the risk structure tensor matrix, perform linear algebraic eigenvalue decomposition on the numerical array inside the matrix, deconstruct the high-dimensional risk superposition effect data into independent projection components on the orthogonal coordinate system, extract the feature vector group representing the main extended dimensions of the risk field in three-dimensional space, calculate the eigenvalues ​​of each feature vector and quantize the vector modulus, and generate principal component feature vector parameters. S502: Based on the principal component feature vector parameters, sort and compare the modulus values ​​corresponding to each feature vector, select the feature vector with the largest weight of the modulus value as the principal axis direction, identify the dominant diffusion path of pollutants in the aquifer medium, calculate the risk diffusion anisotropy index that characterizes the directional intensity of risk diffusion, and generate risk diffusion spatial feature data. S503: Call the aforementioned risk diffusion spatial feature data, extract the dominant diffusion path vector and anisotropy index values, perform spatial mapping and overlay on the feature data and the regional hydrogeological background layer, divide the spatial distribution level of risk level gradient according to the anisotropy intensity values, and generate groundwater environmental quality assessment results.

10. The method for comprehensive assessment of regional groundwater environmental quality and risk of multiple pollutants according to claim 9, characterized in that, The main extended dimension refers to the physical direction indicated by the mutually orthogonal eigenvectors obtained after eigenvalue decomposition of the risk structure tensor matrix in three-dimensional space, representing the independent physical migration trajectory of risk in three-dimensional space along the dominant diffusion direction, the lateral discrete direction, and the vertical settlement direction.

Images