Method for simulating heavy metal adsorption process of soil and underground water through multi-method dynamic coupling

Abstract

The invention provides a multi-method dynamic coupling soil and underground water heavy metal adsorption process simulation method, and belongs to the technical field of soil and underground water pollution simulation. Analyzing and screening main control factors by using a random forest and an SHAP value; constructing a plurality of machine learning models to establish a mapping relation between the main control factor and the adsorption coefficient Kd; and the optimal model is embedded into a TOUGHREACT platform, and dynamic prediction and feedback of a Kd value are realized through Fortran and Python hybrid programming. According to the method, the limitation of static setting of traditional model parameters is overcome, the simulation precision of heavy metal migration and transformation under the action of microorganisms is improved, meanwhile, the prediction advantages of mechanical property and machine learning of a physical model are reserved, and the method is suitable for adsorption behavior simulation and risk assessment of various heavy metals such as Cr (VI), Cu and As in a complex environment.

Description

Technical Field

[0001] This invention relates to the field of soil and groundwater pollution simulation technology, and in particular to a multi-method dynamic coupling simulation method for heavy metal adsorption processes in soil and groundwater. This method is especially suitable for simulating the migration and transformation behavior of heavy metals in the vadose zone-aquifer under the action of microorganisms. By embedding a machine learning model into a traditional multiphysics process simulation platform, it achieves dual protection of dynamic prediction of adsorption parameters and physical mechanisms. Background Technology

[0002] Heavy metal pollution in soil and groundwater is one of the major environmental problems facing the world today. The migration and transformation of heavy metals such as Cr(VI), Cu, and As in the vadose zone-aquifer system are influenced by a variety of environmental factors, including microbial activity, pH, iron oxide content, organic matter content, and temperature. Accurate simulation of these processes is of great significance for pollution risk assessment and remediation strategy development.

[0003] TOUGHREACT, a widely used multiphase flow and multicomponent reactive transport simulation software, can couple thermo-hydraulic-mechanical-chemical processes to simulate pollutant migration behavior in complex geological media. However, traditional TOUGHREACT models have certain limitations when dealing with microbially driven heavy metal adsorption processes: their adsorption parameters (such as the partition coefficient Kd) are usually set as constants or based on simple empirical formulas, making it difficult to reflect their dynamic changes on spatiotemporal scales. Especially under conditions of intense microbial activity and complex geochemical environments, this often leads to significant deviations between simulation results and actual conditions.

[0004] In recent years, machine learning methods have demonstrated great potential in environmental modeling, capable of uncovering complex nonlinear relationships from large amounts of experimental data and achieving high-precision predictions of key parameters. However, pure machine learning models lack physical mechanism support, have poor interpretability, and are difficult to use directly in process-driven numerical simulations.

[0005] In the current technology, there is no publicly available solution for seamlessly embedding and coupling machine learning's dynamic prediction capabilities with process models like TOUGHREACT. Traditional coupling methods often face technical obstacles such as language barriers, inconsistent memory management, and difficulties in parallel communication, resulting in low model integration efficiency and unstable operation.

[0006] Therefore, there is an urgent need for a new simulation framework that can maintain the integrity of the physical model mechanism and introduce data-driven dynamic prediction capabilities to improve the simulation accuracy and practicality of heavy metal adsorption behavior under the action of microorganisms.

[0007] The information disclosed in this background section is intended only to enhance the understanding of the overall background of the invention and should not be construed as an admission or in any way implying that the information constitutes prior art known to those skilled in the art. Summary of the Invention

[0008] The purpose of this invention is to provide a novel simulation framework that can maintain the integrity of the physical model mechanism and introduce data-driven dynamic prediction capabilities, so as to improve the simulation accuracy and practicality of heavy metal adsorption behavior under the action of microorganisms.

[0009] To achieve the above objectives, the present invention provides the following solution: A multi-method dynamically coupled method for simulating the adsorption process of heavy metals in soil groundwater includes the following steps: S1: Obtain relevant data on the adsorption process of heavy metals in the vadose zone-aquifer, including at least 6 environmental factors and the adsorption coefficient Kd; S2: Clean and preprocess the data to establish a database; S3: The random forest algorithm combined with SHAP value analysis is used to quantify the influence weight of each feature on the output parameters and screen the main control factors. S4: Select at least three machine learning algorithms to establish a mapping relationship model between the master control factor and the adsorption coefficient Kd; S5: Construct a basic model of heavy metals in TOUGHREACT, which includes an adsorption calculation equation and a migration and transformation control equation; S6: Embed the mapping relationship model into TOUGHREACT to realize the dynamic adjustment of the adsorption coefficient Kd with the main control factor; S7: Through three steps—state monitoring and data extraction, external machine learning prediction, and parameter feedback and update—the dynamic feedback and simulation iteration of the Kdd value are achieved.

[0010] Optionally, the environmental factors include: Microbial concentration, soil pH, free iron content, amorphous iron content, temperature, and organic matter content.

[0011] Optionally, the data cleaning and preprocessing includes: removing literature data with unclear experimental conditions, outliers, or duplicate publications, and using Excel and Python's Pandas library to remove outliers and fill in missing values.

[0012] Optionally, the formula for calculating the SHAP value is:

[0013] in, S Indicates that it does not contain features iAll possible subsets of features; Indicates the number of features; Indicates the weighting coefficient. Indicates marginal contribution.

[0014] Optionally, step S4 specifically includes: Select at least three machine learning algorithms, including Random Forest (RF), Gradient Boosting Tree (GBRT), and Backpropagation Neural Network (BPNN). Train the model with the main control factor as input and the adsorption coefficient Kd as output. Establish the mapping relationship between the main control factor and the adsorption coefficient Kd. Divide the model into training set, validation set, and test set in a 7:2:1 ratio. Select the model with the best accuracy, a coefficient of determination greater than 0.8, and a root mean square error less than 0.1.

[0015] Optionally, step S5 constructs a basic model of heavy metals in TOUGHREACT, in which...

[0016] In the formula, Porosity; It is the saturation of the λ phase. It is the density of the λ phase, kg / m³ 3 ; The mass fraction of component k in phase λ; The amount of component k that is instantaneously reverse-adsorbed by the soil medium, expressed in kg, can be represented by the following formula:

[0017] In the formula, The density of the medium is kg / m³. 3 , It is the mole fraction of component k in the aqueous phase; it is the partition coefficient between the solid and the liquid, m. 3 / kg.

[0018] Optionally, in step S6, the mapping relationship model is embedded into TOUGHREACT to achieve dynamic adjustment of the adsorption coefficient Kd according to the controlling factor; specifically: By employing a hybrid programming architecture of Fortran and Python, the adsorption coefficient dynamic prediction function based on the random forest algorithm is implemented, grafting data-driven relationships onto a process-driven model to realize the reaction migration simulation of the adsorption process driven by microorganisms. Since TOUGHREACT is written in FORTRAN 77 / 90, its nonlinear solver, chemical reaction module, and thermo-water-mechanical coupling framework are tightly coupled. Directly inserting it into any Python machine learning model would encounter triple obstacles of language barriers, memory management, and parallel communication. Therefore, the embedding path should follow three processes: state monitoring and data extraction, external machine learning prediction, and parameter feedback and update.

[0019] Optionally, the status monitoring and data extraction specifically include: A new function, WRITE_Cr_ADSORP_CSV, has been added to the source code to automatically output the physicochemical state data of the entire grid to a specified file after model initialization and each time step. The extracted indicators include: grid name, simulation time, pH value, amorphous iron content, free iron content, total cation exchange, organic matter content, total microbial content, and the content of major ionic components in the water.

[0020] Optionally, the external machine learning prediction specifically includes: Before the next time step begins, the external machine learning program Kd_predict.py is called by adding a new system function CALL_PYTHON_KD_PREDICT. The physicochemical state data of the global grid generated in the previous step is used as input. The adsorption coefficient Kd of a certain metal in the global grid is predicted with the help of a pre-trained model, and the prediction results are written to the file Kd_prediction_results.csv.

[0021] Optionally, the parameter feedback and update specifically include: The newly added function READ_ML_KD_VALUES reads the Kd values ​​(Kd_prediction_results.csv) predicted by machine learning and feeds them back to the numerical simulator to update the adsorption coefficient Kd of the global grid. Then, the calculation of the next time step is carried out. By repeatedly executing step S7 until the total simulation time ends, the dynamic update of Kd values ​​is achieved, which continuously affects the results of subsequent numerical simulations.

[0022] Compared with the prior art, the present invention has the following beneficial effects: Integrating physical mechanisms with data-driven approaches: By embedding a machine learning model into TOUGHREACT, the physical mechanisms of traditional process models are preserved while the dynamic prediction capability of the adsorption parameter Kd is improved, significantly enhancing simulation accuracy.

[0023] Overcoming language and system barriers: By adopting a hybrid programming architecture of Fortran and Python, the language barriers, memory management and parallel communication issues between TOUGHREACT and machine learning models are resolved, achieving stable and efficient system coupling.

[0024] Dynamic feedback and iterative optimization: Through three steps—state monitoring, external prediction, and parameter feedback—the Kd value is updated in real time and simulated iteratively, enabling the model to adapt to complex environmental changes.

[0025] Multi-algorithm optimization and interpretability analysis: By employing multiple machine learning algorithms and combining them with SHAP value analysis, not only are key controlling factors screened out, but the interpretability and reliability of the model are also enhanced.

[0026] Wide applicability and engineering practicality: The method is applicable to the adsorption simulation of various heavy metals under different environmental conditions, and has strong versatility and practical application value. Attached Figure Description

[0027] To more clearly illustrate the technical solutions in the embodiments of the present invention or the prior art, the drawings used in the embodiments will be briefly introduced below. Obviously, the 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.

[0028] Figure 1 This is a schematic diagram of the method flow provided in an embodiment of the present invention.

[0029] Figure 2 This is a schematic diagram illustrating the process of coupling TOUGHREACT with a machine learning module and exchanging data via a CSV file, as provided in an embodiment of the present invention. Detailed Implementation

[0030] 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 some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0031] The purpose of this invention is to provide a novel simulation framework that can maintain the integrity of the physical model mechanism and introduce data-driven dynamic prediction capabilities, so as to improve the simulation accuracy and practicality of heavy metal adsorption behavior under the action of microorganisms.

[0032] To make the above-mentioned objects, features and advantages of the present invention more apparent and understandable, the present invention will be further described in detail below with reference to the accompanying drawings and specific embodiments.

[0033] Example 1: This embodiment provides a multi-method dynamically coupled simulation method for the adsorption process of heavy metals in soil and groundwater, such as... Figure 1 and 2 As shown, it includes the following steps: Step 1: Research literature on the adsorption process of a certain heavy metal in the vadose zone-aquifer in databases such as Web of Science and CNKI, including but not limited to heavy metals such as Cr(VI), Cu, and As. Obtain relevant data for at least the past 15 years, including at least 6 environmental factors such as microbial concentration, soil pH, free iron content, amorphous iron content, temperature, and organic matter content, as well as the adsorption coefficient Kd.

[0034] Step 2: Based on Step 1, exclude literature with "unclear experimental conditions (e.g., no soil CEC labeling), obvious outliers, or duplicate publications" to ensure data reliability; use Excel and Python (Pandas library) for data cleaning, remove outliers, fill in missing values, and establish a database.

[0035] Step 3: Based on Step 2, the Random Forest algorithm combined with SHAP value (formula) analysis is used to quantify the influence weight of each feature on the output parameters, screen the controlling factors, and quantify the "direction and strength of influence" of each feature by calculating the marginal contribution (SHAP value) of each feature to the model output, avoiding the limitation of traditional importance analysis that can only judge the "weight size". There should be at least 6 controlling factors. The core of SHAP value calculation is based on the "Shapley value" theory. By considering the model prediction result as the total benefit of all features working together, all possible combinations of different feature subsets are first constructed. The marginal contribution of each feature to the prediction result when added to different subsets is calculated (i.e., the difference in prediction when the feature is present or not). Then, these marginal contributions are weighted and averaged according to the probability of each subset (the weight determined by the size of the feature subset). Finally, the specific influence of each feature on the predicted value of a single sample relative to the model's average predicted value is obtained (positive values ​​indicate an increase in the predicted value, and negative values ​​indicate a decrease in the predicted value). The calculation formula is as follows:

[0036] In the formula: S Indicates that it does not contain features i All possible subsets of features; Indicates the number of features; Indicates the weighting coefficient; Indicates marginal contribution.

[0037] Step 4: Based on Step 3, select at least three machine learning algorithms, including Random Forest (RF), Gradient Boosting Tree (GBRT), and Backpropagation Neural Network (BPNN), to train the model using the controlling factor as input and the adsorption coefficient Kd as output, establishing a mapping relationship between the controlling factor and the adsorption coefficient Kd. Divide the model into a training set (70%, used for model parameter learning), a validation set (20%, used for hyperparameter tuning), and a test set (10%, used for final performance evaluation) in a 7:2:1 ratio. Select the model with the best accuracy, a coefficient of determination (R²) greater than 0.8, and a root mean square error (RMSE) less than 0.1.

[0038] Step 5: Construct a TOUGHREACT basic model for a certain heavy metal. In the model,

[0039] In the formula: Φ is porosity; S λ It is the saturation of the λ phase; ρ λ It is the density of the λ phase, kg / m³ 3 ; The mass fraction of component k in phase λ; The amount of component k that is instantaneously reverse-adsorbed by the soil medium, expressed in kg, can be represented by the following formula:

[0040] In the formula: The density of the medium is expressed in kg / m³. It is the mole fraction of component k in the aqueous phase; it is the partition coefficient between the solid and the liquid, m. 3 / kg; Based on total dissolved concentration Consider convection-diffusion-reaction coupling as a variable:

[0041] In the formula: Φ: porosity, which changes dynamically with mineral reactions; S l : Liquid phase saturation; q l Darcy velocity; Local equilibrium reaction source and sink terms; : Kinetic reaction rate (such as mineral dissolution, redox).

[0042] Step 6: Embed the mapping model between the master control factor and the adsorption coefficient Kd established in Step 4 into TOUGHREACT to achieve dynamic adjustment of the adsorption coefficient Kd as the corresponding master control factor changes. The embedding method essentially uses a hybrid programming architecture of Fortran and Python to implement dynamic prediction of the adsorption coefficient based on the random forest algorithm. This grafts the data-driven relationship onto the process-driven model, realizing the reaction-migration simulation of the adsorption process driven by microorganisms. Since TOUGHREACT is written in FORTRAN77 / 90, its nonlinear solver, chemical reaction module, and thermo-water-mechanical coupling framework are tightly coupled. Directly inserting it into any Python machine learning model would encounter triple obstacles: language barriers, memory management, and parallel communication. Therefore, the embedding path should follow three processes: state monitoring and data extraction, external machine learning prediction, and parameter feedback and update.

[0043] Step 7: State Monitoring and Data Extraction: Add a new function WRITE_Cr_ADSORP_CSV to the source code to automatically output the physicochemical state data of the entire grid to a specified file (Cr_adsorp_data.csv) after model initialization and at each time step. Extracted indicators include: grid name, simulation time, pH value, amorphous iron content, free iron content, total cation exchange capacity, organic matter content, total microbial count, and the content of major ionic components in the water.

[0044] Step 8: External Machine Learning Prediction: Before the start of the next time step, the external machine learning program Kd_predict.py is called by adding a new system function CALL_PYTHON_KD_PREDICT. The physicochemical state data of the global grid generated in the previous step (Cr_adsorp_data.csv) is used as input. The pre-trained models (Kd_random_forest_model.pkl and scaler.pkl) are used to predict the adsorption coefficient Kd of a certain metal in the global grid, and the prediction results are written to the file Kd_prediction_results.csv.

[0045] Step 9: Parameter Feedback and Update: The newly added function READ_ML_KD_VALUES reads the Kd values ​​predicted by machine learning (Kd_prediction_results.csv) and feeds them back to the numerical simulator to update the adsorption coefficient Kd of the global grid, thus proceeding to the calculation of the next time step. By repeatedly executing steps 7 to 9 (until the total simulation time ends), the dynamic update of the Kd values ​​is achieved, continuously influencing the results of subsequent numerical simulations.

[0046] Step 10: The TOUGHREACT model can be solved using the Newton-Raphson iterative method. Indoor soil column experiments can be conducted, and the results of the traditional TOUGHREACT model and the modified TOUGHREACT model can be compared with those of the soil column experiments to verify the effectiveness of the modified model.

[0047] This embodiment addresses the complex scenario of heavy metal migration and transformation in the vadose zone-aquifer under microbial action. Without disrupting the original TOUGHREACT architecture, it seamlessly embeds the mapping relationship between the master control factors learned by machine learning and the key parameters of the migration and transformation process. This improves the model's prediction accuracy while preserving the physical interpretability of the process model, providing a new computational framework for simulating the migration and transformation of heavy metals under multi-field coupling conditions.

[0048] The various embodiments in this specification are described in a progressive manner, with each embodiment focusing on its differences from other embodiments. Similar or identical parts between embodiments can be referred to interchangeably. For the systems disclosed in the embodiments, since they correspond to the methods disclosed in the embodiments, the descriptions are relatively simple; relevant parts can be referred to the method section.

[0049] This document uses specific examples to illustrate the principles and implementation methods of the present invention. The descriptions of the above embodiments are only for the purpose of helping to understand the method and core ideas of the present invention. Furthermore, those skilled in the art will recognize that, based on the ideas of the present invention, there will be changes in the specific implementation methods and application scope. Therefore, the content of this specification should not be construed as a limitation of the present invention.

Claims

1. A multi-method dynamically coupled method for simulating the adsorption process of heavy metals in soil and groundwater, characterized in that, Includes the following steps: S1: Obtain relevant data on the adsorption process of heavy metals in the vadose zone-aquifer, including at least 6 environmental factors and the adsorption coefficient Kd; S2: Clean and preprocess the data to establish a database; S3: The random forest algorithm combined with SHAP value analysis is used to quantify the influence weight of each feature on the output parameters and screen the main control factors. S4: Select at least three machine learning algorithms to establish a mapping relationship model between the master control factor and the adsorption coefficient Kd; S5: Construct a basic model of heavy metals in TOUGHREACT, which includes an adsorption calculation equation and a migration and transformation control equation; S6: Embed the mapping relationship model into TOUGHREACT to realize the dynamic adjustment of the adsorption coefficient Kd with the main control factor; S7: Through three steps—state monitoring and data extraction, external machine learning prediction, and parameter feedback and update—the dynamic feedback and simulation iteration of the Kdd value are achieved.

2. The method for simulating the adsorption process of heavy metals in soil and groundwater by dynamic coupling of multiple methods according to claim 1, characterized in that, The environmental factors include: Microbial concentration, soil pH, free iron content, amorphous iron content, temperature, and organic matter content.

3. The method for simulating the adsorption process of heavy metals in soil and groundwater by dynamic coupling of multiple methods according to claim 1, characterized in that, The data cleaning and preprocessing process includes: removing literature data with unclear experimental conditions, outliers, or duplicate publications, and using Excel and Python's Pandas library to remove outliers and fill in missing values.

4. The method for simulating the adsorption process of heavy metals in soil and groundwater by dynamic coupling of multiple methods according to claim 1, characterized in that, The formula for calculating the SHAP value is as follows: in, S Indicates that it does not contain features i All possible subsets of features; Indicates the number of features; Indicates the weighting coefficient. Indicates marginal contribution.

5. The method for simulating the adsorption process of heavy metals in soil and groundwater by dynamic coupling of multiple methods according to claim 1, characterized in that, Step S4 specifically involves: Select at least three machine learning algorithms, including Random Forest (RF), Gradient Boosting Tree (GBRT), and Backpropagation Neural Network (BPNN). Train the model with the main control factor as input and the adsorption coefficient Kd as output. Establish the mapping relationship between the main control factor and the adsorption coefficient Kd. Divide the model into training set, validation set, and test set in a 7:2:1 ratio. Select the model with the best accuracy, a coefficient of determination greater than 0.8, and a root mean square error less than 0.

1.

6. The method for simulating the adsorption process of heavy metals in soil and groundwater by dynamic coupling of multiple methods according to claim 1, characterized in that, Step S5 constructs a basic model of heavy metals in TOUGHREACT, in which... In the formula, Porosity; It is the saturation of the λ phase. It is the density of the λ phase, kg / m³ 3 ; The mass fraction of component k in phase λ; The amount of component k that is instantaneously reverse-adsorbed by the soil medium, expressed in kg, can be represented by the following formula: In the formula, The density of the medium is kg / m³. 3 , It is the mole fraction of component k in the aqueous phase; it is the partition coefficient between the solid and the liquid, m. 3 / kg.

7. The method for simulating the adsorption process of heavy metals in soil and groundwater by dynamic coupling of multiple methods according to claim 1, characterized in that, In step S6, the mapping relationship model is embedded into TOUGHREACT to achieve dynamic adjustment of the adsorption coefficient Kd with the main control factor; specifically: By employing a hybrid programming architecture of Fortran and Python, the adsorption coefficient dynamic prediction function based on the random forest algorithm is implemented, grafting data-driven relationships onto a process-driven model to realize the reaction migration simulation of the adsorption process driven by microorganisms. Since TOUGHREACT is written in FORTRAN 77 / 90, its nonlinear solver, chemical reaction module, and thermo-water-mechanical coupling framework are tightly coupled. Directly inserting it into any Python machine learning model would encounter triple obstacles of language barriers, memory management, and parallel communication. Therefore, the embedding path should follow three processes: state monitoring and data extraction, external machine learning prediction, and parameter feedback and update.

8. The method for simulating the adsorption process of heavy metals in soil and groundwater by dynamic coupling of multiple methods according to claim 1, characterized in that, The specific steps of the status monitoring and data extraction are as follows: A new function, WRITE_Cr_ADSORP_CSV, has been added to the source code to automatically output the physicochemical state data of the entire grid to a specified file after model initialization and each time step. The extracted indicators include: grid name, simulation time, pH value, amorphous iron content, free iron content, total cation exchange, organic matter content, total microbial content, and the content of major ionic components in the water.

9. The method for simulating the adsorption process of heavy metals in soil and groundwater by dynamic coupling of multiple methods according to claim 1, characterized in that, The external machine learning prediction specifically refers to: Before the next time step begins, the external machine learning program Kd_predict.py is called by adding a new system function CALL_PYTHON_KD_PREDICT. The physicochemical state data of the global grid generated in the previous step is used as input. The adsorption coefficient Kd of a certain metal in the global grid is predicted with the help of a pre-trained model, and the prediction results are written to the file Kd_prediction_results.csv.

10. The method for simulating the adsorption process of heavy metals in soil and groundwater by dynamic coupling of multiple methods according to claim 1, characterized in that, The parameter feedback and update are specifically as follows: The newly added function READ_ML_KD_VALUES reads the Kd values ​​(Kd_prediction_results.csv) predicted by machine learning and feeds them back to the numerical simulator to update the adsorption coefficient Kd of the global grid. Then, the calculation of the next time step is carried out. By repeatedly executing step S7 until the total simulation time ends, the dynamic update of Kd values ​​is achieved, which continuously affects the results of subsequent numerical simulations.

Images