A method for evaluating the control of ground runoff pollutants by an ecological buffer zone

Abstract

The application relates to a method for evaluating the control of ground surface runoff pollutants by an ecological buffer zone, and belongs to the technical field of environmental engineering, and aims to solve the problems of staticity, inefficiency and difficulty in reflecting a dynamic process of an existing evaluation method. The method comprises the following steps: constructing a dynamic monitoring system fusing distributed sensing and remote sensing data; establishing a pollutant migration and transformation model based on slope water flow dynamics and coupling vegetation resistance and soil adsorption processes; quantifying the nonlinear relationship between buffer zone structure parameters and pollutant reduction efficiency; calculating a dynamic reduction rate through multi-scenario simulation; finally forming a comprehensive evaluation matrix containing process and result indexes, and outputting an optimized configuration scheme. Through high-frequency monitoring and physical-data double-driven modeling, minute-level dynamic evaluation is realized, the prediction accuracy of nitrogen and phosphorus and other pollutants is significantly improved, and cross-regional popularization and near-real-time intelligent decision-making are supported.

Description

Technical Field

[0001] This application belongs to the field of environmental engineering technology, specifically relating to a method for assessing the control of surface runoff pollutants by ecological buffer zones. Background Technology

[0002] Ecological buffer zones are key ecological engineering measures for controlling surface runoff pollutants and preventing agricultural non-point source pollution. They mainly purify runoff through processes such as vegetation interception, soil infiltration, and adsorption. However, current technologies for assessing their effectiveness are significantly inadequate, hindering their design optimization and precise management.

[0003] Existing assessment methods are mostly based on static empirical parameters or simplified models, which are difficult to reflect the complex dynamic coupling between hydrology, pollutant migration, and vegetation-soil systems during actual rainfall-runoff processes. The assessment results often deviate from the actual performance of buffer zones under variable weather conditions and long-term operation.

[0004] Furthermore, performance assessments rely heavily on intermittent manual sampling and laboratory analysis, resulting in highly sparse data in both time and space. This makes it impossible to obtain continuous changes in pollutant concentrations and hydrological parameters during runoff events, leading to delayed assessments and difficulty in supporting real-time management decisions.

[0005] Although some mechanistic models exist, they often fail to adequately characterize key processes such as the dynamic impact of vegetation resistance and the nonlinear adsorption of pollutants in the soil, or their parameters are rigid and lack the ability to adaptively correct using field data. Furthermore, existing methods have failed to systematically establish quantitative and nonlinear response relationships between buffer zone structural characteristics (such as vegetation composition, width, and soil structure) and their purification functions, making it difficult to use these models to guide structural optimization design in different scenarios.

[0006] Therefore, there is an urgent need in this field for a comprehensive assessment method that can achieve dynamic monitoring, couple key process mechanisms, and quantify structure-function relationships, so as to scientifically, in real time, and accurately evaluate and improve the pollutant control effectiveness of ecological buffer zones. Summary of the Invention

[0007] The purpose of this invention is to provide a method for assessing the control of surface runoff pollutants by ecological buffer zones, which can effectively solve the problems mentioned in the background art.

[0008] To achieve the above objectives, the technical solution adopted by the present invention is as follows:

[0009] A method for assessing the control of surface runoff pollutants by ecological buffer zones includes the following specific steps:

[0010] S1: Construct a dynamic monitoring system that integrates in-situ monitoring and remote sensing inversion data to obtain high spatiotemporal resolution hydrological and water quality, vegetation and soil parameter data for the ecological buffer zone area;

[0011] S2: Establish a dynamic model of pollutant migration and transformation based on the coupling of physical mechanisms and data-driven approaches. The dynamic model of pollutant migration and transformation is coupled with a hydrodynamic module describing shallow water flow on the slope, a resistance module that explicitly quantifies the effect of vegetation resistance, and a dynamic module describing the adsorption behavior of pollutants in the soil. The key parameters of the model are corrected online using the real-time monitoring data obtained in step S1.

[0012] S3: Quantify the nonlinear response relationship between the structural parameters of the ecological buffer zone and the pollutant reduction function. Based on the multi-dimensional structural parameters of vegetation, root system and soil obtained in step S1, construct a structural feature vector, and use a machine learning model to establish the mapping relationship between the structural feature vector and the pollutant reduction efficiency output in step S2. Then, identify the key structural factors affecting the reduction function and their effective thresholds through sensitivity analysis.

[0013] S4: Implement multi-scenario dynamic performance evaluation. Define multiple evaluation scenarios by combining different rainfall return periods, initial pollutant loads, and buffer zone degradation states. Drive the dynamic model of pollutant migration and transformation corrected in step S2 to simulate the scenario, calculate and output the dynamic reduction rate and key performance indicators of pollutants under each scenario.

[0014] S5: Generate a comprehensive evaluation index system and optimization configuration scheme, integrate the dynamic performance index obtained in step S4 and the structural adaptability parameters obtained in step S3, construct a comprehensive evaluation matrix that includes process performance, system stability and structural adaptability, and generate a targeted ecological buffer zone optimization configuration scheme based on the key structural factor thresholds identified in step S3.

[0015] Furthermore, in step S1, constructing a dynamic monitoring system includes:

[0016] At least three monitoring sections are set up along the main slope of the ecological buffer zone. Each section is equipped with monitoring points for the upstream inlet, the middle transition zone and the downstream outlet to form a three-dimensional monitoring grid.

[0017] Integrated water quality-hydrology-vegetation sensing units are deployed at each monitoring point to collect data on total nitrogen concentration, total phosphorus concentration, sediment concentration, flow velocity, water level, and rainfall intensity in surface runoff in situ and at high frequency.

[0018] Acquire high-resolution multispectral or radar remote sensing images covering the target area, and extract spatial distribution data of vegetation cover, leaf area index and soil volumetric water content through supervised classification and inversion models.

[0019] The in-situ monitoring data and remote sensing inversion data are spatiotemporally aligned and fused to form a multi-source heterogeneous dataset.

[0020] Furthermore, in step S2, establishing a dynamic model for pollutant migration and transformation includes:

[0021] Using the simplified Saint-Venant equations as the framework for the hydrodynamic module, the slope flow motion under the influence of rainfall, infiltration, and topography is simulated.

[0022] In the hydrodynamic module, the comprehensive resistance slope is decomposed into bed surface resistance and vegetation resistance. The vegetation resistance is dynamically calculated using a formula based on the vegetation stem drag coefficient, stem density, equivalent stem diameter, and hydrodynamic parameters.

[0023] Convection-diffusion-reaction transport equations for total nitrogen and total phosphorus were established, with second-order kinetic equations used for the source and sink terms to describe the adsorption and desorption processes of pollutants between the solid and liquid phases of the soil.

[0024] The hydrodynamics module, resistance module and adsorption kinetics module are coupled together and numerical solutions are obtained using the finite volume method.

[0025] Using the real-time upstream inflow water quality, water level, and flow rate data obtained in step S1, the Manning roughness coefficient and soil adsorption rate constant in the model are periodically corrected online using the extended Kalman filter algorithm.

[0026] Furthermore, in step S3, the nonlinear response relationship between the quantized structural parameters and the reduction function includes:

[0027] Based on the supervised classification results of remote sensing images, the proportion of vegetation functional types of herbaceous plants, shrubs, and trees was statistically analyzed;

[0028] Root length density, soil macropore ratio, pore connectivity and average pore size were obtained through root scanning analysis and soil CT scanning technology.

[0029] A multidimensional structural feature vector is constructed by integrating the buffer zone width, the proportion of the vegetation functional type, the root system and soil pore structure parameters, soil organic matter content, soil bulk density, leaf area index and vegetation coverage.

[0030] Using the structural feature vector as input and the pollutant reduction rate per unit width obtained through model simulation in step S2 as output, a gradient boosting decision tree model is trained to establish a structure-function nonlinear mapping relationship.

[0031] The Sobol global sensitivity analysis method was used to quantify the contribution of each structural feature to the uncertainty of pollutant reduction rate prediction, identify key control factors, and determine the functional threshold of the key control factors based on the model response curve.

[0032] Furthermore, in step S4, the multi-scenario dynamic performance evaluation includes:

[0033] Define multi-dimensional assessment boundary conditions that include different rainfall return periods, different initial pollutant load levels, and the health and degradation status of the buffer zone;

[0034] Based on the local rainstorm intensity formula, design rainstorm process lines under different return periods are generated, and initial concentrations of low, medium, and high levels of pollutants are set;

[0035] The degradation state of the buffer zone is parameterized as a reduction in the vegetation resistance coefficient and the proportion of macropores in the soil in the model.

[0036] The rainfall recurrence period scenario, the initial pollutant load level, and the buffer zone status are fully combined to construct a multi-scenario input set;

[0037] For each scenario, the corrected model in step S2 is run to simulate the spatiotemporal changes in pollutant concentration, and the dynamic reduction rate of pollutants during the entire runoff event is calculated based on the concentration process lines at the upstream inlet and downstream outlet.

[0038] The peak reduction rate, average reduction rate, and reduction duration are extracted from the dynamic reduction rate process and output as key performance indicators.

[0039] Furthermore, in step S5, generating the comprehensive evaluation index system and optimized configuration scheme includes:

[0040] Three primary evaluation indicators are defined: process efficiency, system stability, and structural adaptability. Process efficiency is directly characterized by the dynamic reduction rate calculated in step S4. System stability is quantified by calculating the reciprocal of the standard deviation of the reduction rate in multiple consecutive runoff events. Structural adaptability is quantified by calculating the weighted sum of plant Shannon diversity index, soil organic matter content, and root biomass.

[0041] The weights of each sub-indicator of structural adaptability were determined using the analytic hierarchy process (AHP), and the total nitrogen reduction rate, total phosphorus reduction rate, sediment reduction rate, system stability index, and structural adaptability score were integrated into a visual comprehensive evaluation matrix.

[0042] Based on the evaluation results of each indicator in the comprehensive evaluation matrix, and combined with the threshold of key structural factors identified in step S3, specific optimization configuration suggestions for vegetation configuration adjustment, plant replanting, soil improvement, or vegetation management are generated.

[0043] Furthermore, it also includes time-series clustering analysis and knowledge transfer of historical monitoring data, specifically:

[0044] Extract multivariate time series data of multiple complete historical runoff-pollution events;

[0045] The dynamic time warping algorithm is used to calculate the morphological similarity distance between events, and the K-means clustering algorithm is applied to divide historical events into multiple typical pattern clusters;

[0046] The cluster centers of each cluster are extracted as representative event sequences, which drive the model in step S2 to generate structural functional expansion samples corresponding to different event patterns, thereby enhancing the generalization ability of the machine learning model in step S3.

[0047] Furthermore, it also includes: establishing a regional-scale spatial database of ecological buffer zone effectiveness to store buffer zone structural parameters, environmental driving parameters, and effectiveness indicators under different geographical locations, climatic conditions, and land use types;

[0048] Based on the aforementioned spatial database, the migration and application of cross-regional buffer zone design parameters are supported. That is, between areas with similar climate and land use, the core parameters of the verified buffer zone structure-function relationship can be migrated, and only the local soil infiltration parameters are adjusted for adaptability.

[0049] Furthermore, the online calibration of the model in step S2 and the multi-scenario simulation evaluation in step S4 are achieved through near real-time computing via an edge computing gateway deployed on-site.

[0050] The edge computing gateway loads a simplified pollutant migration and transformation model, periodically receives monitoring data from step S1 and performs localized performance evaluation calculations, and uploads the core state variables obtained from the evaluation to the cloud platform.

[0051] Soil pore structure parameters are obtained through soil CT scanning technology, specifically including:

[0052] An X-ray computed tomography system was used to scan the undisturbed soil column to obtain a three-dimensional volumetric data model of the soil structure.

[0053] Based on the aforementioned three-dimensional volumetric data model, the maximum sphere algorithm is applied to extract the pore network;

[0054] The proportion of large pore volume with an equivalent pore diameter greater than 75 micrometers is calculated as the large pore ratio. The number of pore tunnel connection points per unit volume is counted as the pore connectivity. The arithmetic mean of the equivalent diameters of all identified pores is calculated as the average pore diameter.

[0055] In summary, this application includes at least one of the following beneficial technical effects:

[0056] 1. This application constructs a dynamic monitoring system by integrating high-frequency in-situ monitoring and remote sensing data, and establishes a mechanism-data coupling model that can be corrected in real time, so as to realize the dynamic simulation and evaluation of the pollutant control effectiveness of the ecological buffer zone at the minute level, overcoming the shortcomings of traditional methods that rely on static parameters and lagging data and cannot reflect the dynamic process of actual runoff events.

[0057] 2. This application constructs a dynamic model with a clear mechanism and self-optimization by coupling key physicochemical processes such as hydrodynamics, vegetation resistance, and pollutant adsorption and transformation, and by using real-time monitoring data to perform online adaptive correction of the model's core parameters. This effectively solves the problems of insufficient simulation accuracy and poor adaptability of traditional empirical models or fixed-parameter mechanistic models under complex real-world conditions.

[0058] 3. This invention, through quantitative analysis of the nonlinear response relationship between the multi-dimensional structural parameters of the buffer zone and its purification function, can accurately identify key structural factors affecting effectiveness and their effective thresholds. This allows the evaluation conclusions to be directly transformed into specific optimized configuration schemes such as vegetation configuration and soil improvement, providing a quantitative decision-making basis for the scientific design and precise management of ecological buffer zones. Attached Figure Description

[0059] Figure 1 This is a schematic diagram of the overall scheme for assessing the method of controlling surface runoff pollutants in ecological buffer zones;

[0060] Figure 2 This is a schematic diagram of a pollutant migration and transformation model based on the coupling of physical mechanisms and data-driven approaches.

[0061] Figure 3 It is a flowchart of the process for multi-scenario dynamic performance evaluation and structure-function response relationship quantification. Detailed Implementation

[0062] To make the objectives, technical solutions, and advantages of this invention clearer, the invention will be further described in detail below with reference to the accompanying drawings and specific embodiments.

[0063] In the method for assessing the control of surface runoff pollutants by ecological buffer zones, the first step is to construct a dynamic monitoring system that integrates multi-source heterogeneous data. The specific operation process is as follows:

[0064] Step S101: Design and deploy a three-dimensional monitoring network; within the target ecological buffer zone area, determine the main slope direction of surface runoff collection and flow based on digital elevation model analysis or on-site topographic survey.

[0065] Along the main slope, at least three monitoring sections perpendicular to the water flow direction are set up. The spacing between the sections is determined according to the length of the buffer zone and the topographical changes to ensure spatial representativeness.

[0066] Within each monitoring section, three types of functional monitoring points are clearly set up: the upstream inlet point located at the starting boundary of the buffer zone, the middle transition zone point located inside the buffer zone, and the downstream outlet point located at the ending boundary of the buffer zone.

[0067] The spatial distribution of all monitoring points should form a three-dimensional monitoring grid that covers the entire width of the buffer zone and extends along the slope.

[0068] Step S102: Deploy the integrated sensor acquisition unit; at each monitoring point, install the integrated water quality-hydrology-vegetation sensor unit, which is fixed on an adjustable mounting bracket to ensure that the sensor is in full contact with the surface runoff or soil.

[0069] The sensing unit includes an embedded microcontroller, a multi-channel analog signal acquisition module, a digital communication interface, and a power management module.

[0070] An embedded microcontroller coordinates the operation commands of various sensors, performs data preprocessing, and performs temporary storage. A multi-channel analog signal acquisition module receives and converts standard analog signals output from various sensors. A digital communication interface supports wired or wireless communication protocols for transmitting acquired data to data aggregation nodes or edge computing gateways. A power management module supplies power to the entire unit, employing a combination of solar panels and batteries to ensure continuous and stable operation in outdoor environments.

[0071] Step S103: Configure the sensor and define key parameters; water quality monitoring adopts an in-situ online electrochemical and optical composite probe to simultaneously measure the concentrations of total nitrogen, total phosphorus and dissolved organic pollutants in surface runoff in real time.

[0072] The detection limit for total nitrogen is not less than 0.05 mg / L, the detection limit for total phosphorus is not less than 0.01 mg / L, and the sampling frequency is not less than once per minute.

[0073] Sediment concentration is measured using a sensor based on the principle of optical scattering, with a range of 0 to 20,000 mg / L and a resolution better than 10 mg / L.

[0074] Hydrological parameters were measured using a combination of an ultrasonic Doppler current meter and a pressure level gauge, and the instantaneous flow rate was calculated accordingly.

[0075] Surface runoff is calculated using the Manning formula, which is expressed as:

[0076]

[0077] in, Flow rate, measured in cubic meters per second; The cross-sectional area of ​​the water passage is expressed in square meters. The roughness coefficient is Manning's coefficient, the value of which is obtained through on-site calibration or by consulting engineering manuals; Hydraulic radius, in meters; This refers to the hydraulic gradient.

[0078] By regularly conducting on-site comparative measurements using standard measurement methods, the overall flow error is controlled within ±3%. Rainfall intensity is recorded in real time by a tipping bucket rain gauge with a time resolution of 1 minute.

[0079] Step S104: Acquire and process remote sensing data; at the same time, acquire high-resolution multispectral or radar remote sensing image data covering the target area, with a spatial resolution of not less than 0.5 meters and a revisit period of not more than 5 days.

[0080] Supervised classification algorithms are used to interpret vegetation types in images. First, training samples are selected based on field survey results. Then, algorithms such as support vector machines or random forests are used to perform classification and extract key vegetation parameters such as vegetation coverage and leaf area index.

[0081] Soil volumetric water content parameters were obtained by inverting the backscattering coefficients from synthetic aperture radar (SAR) images. Sentinel-1 satellite VV+VH dual-polarization combined data were used, with a spatial resolution of 10 meters and a revisit period of 6 days. The inversion employed a water-cloud model based on a physical scattering mechanism, which separates the scattering contributions from the vegetation layer and the soil layer. The expression is as follows:

[0082]

[0083] in, The total backscattering coefficient of the radar is expressed in decibels (dB). Contributes to direct scattering from the vegetation layer. , An empirical coefficient related to vegetation type. For vegetation coverage, The radar incident angle; The bidirectional transmittance coefficient of the vegetation layer. B is an empirical attenuation coefficient related to vegetation structure, and like A, it is a model empirical coefficient. This refers to vegetation water content or leaf area index. The contribution of direct backscattering from the soil layer is described using the Oh model:

[0084]

[0085] in, Soil volumetric water content (target inversion parameter), in cubic meters per cubic meter; Radar wave number, , The radar wavelength; The root mean square height of the soil surface, in meters, is obtained through LiDAR data or field measurements.

[0086] Empirical coefficients in the model , and soil surface parameters Calibration was performed based on measured samples within the study area: no fewer than 20 ground sampling points were set up within the study area, and measured data of soil volumetric moisture content were collected simultaneously (using a time domain reflectometer or the drying and weighing method). The multi-temporal SAR backscattering coefficients at the corresponding locations were extracted, and the least squares method was used for fitting to determine the SAR backscattering coefficients. , and The optimal value was obtained. After calibration, the calibrated model was applied to the entire region to retrieve the spatial distribution of soil volumetric water content pixel by pixel.

[0087] The inversion results must be validated by measured surface soil moisture data within the study area using a 5-fold cross-validation method. The root mean square error (RMSE) of the final inversion results should not exceed 0.04 cubic meters per cubic meter, and the coefficient of determination should be [not specified]. Not less than 0.75.

[0088] Step S105: Construct a spatiotemporally aligned dataset; fuse the in-situ monitoring data stream obtained in step S103 with the remote sensing inversion parameter product obtained in step S104.

[0089] All data are standardized to UTC, and the coordinate system is aligned to the same geographic and projected coordinate systems through spatial registration.

[0090] Ultimately, a time-synchronized and spatially matched multi-source heterogeneous dataset is formed, providing complete data input for subsequent dynamic modeling.

[0091] In summary, step S1 systematically constructs a space-ground collaborative dynamic monitoring system. This system achieves high-frequency acquisition of key hydrological and water quality parameters through a ground-based sensor network, periodically extracts the spatial characteristics of vegetation and soil moisture through remote sensing technology, and forms a complete data base through spatiotemporal alignment and fusion. This system provides indispensable high spatiotemporal resolution driving data and verification basis for the subsequent establishment and real-time correction of pollutant migration and transformation models, and is the fundamental prerequisite for this method to achieve minute-level dynamic assessment and accurate prediction.

[0092] For step S2: Establish a pollutant migration and transformation model based on the coupling of physical mechanisms and data-driven approaches, which includes the following implementation steps;

[0093] Step S201: Construct the hydrodynamic model framework; use the simplified form of the Saint-Venant equations describing shallow water flow on a slope as the hydrodynamic basis of the model. This equation set includes the continuity equation and the momentum equation, specifically expressed as follows:

[0094] Continuity equation:

[0095]

[0096] Momentum equation:

[0097]

[0098] The meanings and units of each combination term in the above equations are as follows: ∂ represents the partial differential symbol. The rate of change of water depth over time, expressed in meters per second; , They are respectively , Rate of change of flow rate per unit width in a given direction, in meters per second; for The rate of change of directional single-width momentum with time, in Newtons per square meter (or kilograms per second squared). for The gradient of directional momentum flux, in Newtons per square meter; for Directional momentum flux at The gradient in the direction, expressed in Newtons per square meter.

[0099] More specifically, Rainfall intensity, measured in meters per second; This represents the soil infiltration rate, expressed in meters per second. Let gravitational acceleration be the acceleration due to gravity, and its value be [value missing]. meters per second squared; The slope of the terrain; The slope is for overall resistance.

[0100] Among them, soil infiltration rate The Green-Ampt model is used for dynamic calculation, and the expression is:

[0101]

[0102] In the formula, The saturated hydraulic conductivity, in meters per second, is estimated using the soil transfer function based on the soil texture and bulk density parameters obtained in step S1. This represents the saturated volumetric water content. The initial volumetric water content is obtained by remote sensing inversion in step S104 or in-situ monitoring by soil moisture sensor in step S103. The wetted frontal suction is measured in meters and can be determined by consulting literature or through experiments based on soil texture type. The cumulative infiltration rate, expressed in meters, is calculated through internal model integration.

[0103] In the minute-level dynamic simulation, the model in step S2 is based on the real-time updated initial moisture content at each time step. Based on the soil moisture monitoring data from step S1, dynamically calculate the infiltration rate at the current moment. And substitute it as a source term into the continuity equation to solve.

[0104] Step S202: Explicit vegetation resistance effect; comprehensive resistance slope It is decomposed into two parts: bed surface resistance and vegetation resistance, namely Among them, bed surface resistance The Manning formula was used for calculation. Vegetation resistance. The following formula is used for dynamic calculation:

[0105]

[0106] In this formula, The drag coefficient of the vegetation stem (dimensionless). Plant stem density, expressed in units of plants per square meter (plants are a counting unit, and the dimension is considered as...). ); The equivalent diameter of the stem, in meters (dimensions). ); The dimensions are , representing the vegetation projection length per unit area (i.e., the characteristic length of the vegetation drag area within a unit water volume). The velocity of water flow is expressed in meters per second (dimensions not specified). ); Let gravitational acceleration be the acceleration due to gravity, and its value be [value missing]. meters per second squared (dimensions) ).

[0107] The physical meaning of this formula is: the drag force exerted by vegetation on a unit weight of water. Dimensional verification: In the molecule on the right side, Dimensionless The dimensions are , The dimensions are , The dimensions are Therefore, the overall dimension of the molecule is ; denominator The dimensions are Therefore, the fraction is dimensionless, and it is related to the slope of the resistance on the left. Their dimensions (dimensionless) are completely consistent.

[0108] Drag coefficient Equivalent diameter of stem The value is determined according to methods known in the art. Among them, the drag coefficient... Values ​​are assigned based on vegetation type and water flow conditions, referring to empirical ranges provided in publicly available literature or manuals in the field of vegetation hydraulics (such as "Introduction to Vegetation Hydraulics" or relevant standards). Common herbaceous plants... Values ​​range from 0.1 to 0.5 for shrubs and trees. The value ranges from 0.5 to 1.2; the equivalent diameter of the stem. By setting up quadrats within the buffer zone, representative vegetation stems were measured in-situ, and the average stem diameter at 10 cm above the ground was obtained. During measurement, no fewer than 20 target vegetation plants were randomly selected from each quadrat, and their stem diameters were measured using calipers. The arithmetic mean was taken as the equivalent stem diameter for that vegetation type. If multiple vegetation types existed within the same buffer zone, a weighted average was calculated based on their respective coverage area percentages to determine the comprehensive equivalent diameter.

[0109] Through steps S201 and S202, a hydrodynamic module is established that can simulate the flow law of slope water under the combined effects of rainfall, infiltration, topography and vegetation, providing the basis for the velocity field and flow field of pollutant migration.

[0110] Step S203: Establish pollutant transport and adsorption kinetic equations; for the target pollutants total nitrogen and total phosphorus, establish transport equations based on the convection-diffusion-reaction mechanism, with the source and sink terms focusing on describing the distribution process of pollutants between the solid and liquid phases of the soil.

[0111] The adsorption and desorption behavior of nitrogen and phosphorus in soil is described by a second-order kinetic equation:

[0112]

[0113] in, The amount of pollutants adsorbed by a unit mass of soil is expressed in milligrams per gram. This represents the amount of soil adsorbed per unit mass at adsorption equilibrium, expressed in milligrams per gram. This represents the second-order adsorption rate constant, expressed in grams per milligram per minute. Represents time, in minutes.

[0114] equilibrium adsorption capacity With adsorption rate constant The specific values ​​were obtained through batch equilibrium adsorption experiments on typical soil samples within the buffer zone of the study area. The experimental method is as follows:

[0115] (1) Soil sample collection and pretreatment: Within the buffer zone of the study area, select no fewer than three representative sampling points based on vegetation type and topographic differences to collect undisturbed topsoil (0-20 cm). After the soil samples are air-dried, pass them through a 2 mm sieve to remove gravel and plant debris, and set them aside for later use.

[0116] (2) Adsorption kinetics experiment: Weigh 2.00 g of the sieved soil sample and place it in a 50 mL centrifuge tube. Add 40 mL of an initial solution containing the target pollutant (total nitrogen or total phosphorus), with an initial concentration of 5 mg / L (calculated as nitrogen or phosphorus). The solid-liquid ratio (soil:water) is 1:20. Place the centrifuge tube in a constant temperature shaker and shake at 150 rpm at 25±1℃. Samples are taken at 5, 10, 20, 30, 45, 60, 90, 120, 180, 240, 360, 480, 720, and 1440 minutes after the start of shaking. Immediately centrifuge at 4000 rpm for 10 min and take the supernatant to determine the pollutant concentration. Calculate the adsorption capacity per unit mass of soil based on the difference between the initial concentration and the residual concentration at different time points. The adsorption capacity as a function of time was plotted. The second-order adsorption rate constant was obtained by fitting the second-order kinetic equation described above. .

[0117] (3) Isothermal adsorption experiment: Weigh 2.00 g of the sieved soil sample and place it in a 50 mL centrifuge tube. Add 40 mL of the initial solution containing the target pollutant, with initial concentration gradients of 0.5, 1.0, 2.0, 4.0, 8.0, 12.0, 16.0, and 20.0 mg / L (calculated as nitrogen or phosphorus). The solid-liquid ratio is also 1:20. Place the centrifuge tube in a constant temperature shaker and shake at 150 rpm at 25±1℃ until adsorption equilibrium is reached (the equilibrium time is determined according to the kinetic experiment results, generally 24 hours). After shaking, centrifuge at 4000 rpm for 10 min, and take the supernatant to determine the equilibrium concentration. Calculate the corresponding equilibrium adsorption amount. .by right Plot the values ​​and fit them using the Langmuir or Freundlich model to obtain the equilibrium adsorption amounts under different concentration conditions. , used for model input.

[0118] All experiments were conducted at the same temperature (25±1℃), with three parallel samples for each experiment, and the results were taken as the arithmetic mean.

[0119] The adsorption kinetic equation is coupled with the water flow equation established in step S201 to solve the problem, thereby dynamically simulating the entire process of pollutants migrating with water in the aqueous phase and being trapped and adsorbed by soil in the solid phase.

[0120] Step S204: Numerical solution of the model and real-time parameter correction; Before solving the model, spatial discretization is required. The buffer zone study area is divided into an unstructured triangular mesh, and the spatial discretization step length of the model is set to 1 meter.

[0121] The time integration step was set to 10 seconds, and the finite volume method was used to numerically discretize and solve the above coupled equations.

[0122] After the model is initialized and running, it receives real-time data streams from the monitoring system built in step S1 through the edge computing gateway deployed on site.

[0123] An extended Kalman filter algorithm is used to automatically correct key dynamic parameters in the model online using real-time monitored upstream inflow water quality, water level and flow rate data. These parameters include, but are not limited to, Manning roughness coefficient and soil adsorption rate constant.

[0124] The calibration process is performed in cycles of 5 minutes to ensure that the model's computational state is synchronized with the feedback from the actual monitoring system, thereby continuously improving the model's simulation accuracy and reliability.

[0125] In summary, step S2 constructs a dynamic model driven by both mechanism and data. First, a basic framework characterizing water flow motion is established based on hydrodynamic principles, and the resistance effect of vegetation on water flow is explicitly quantified. Then, a kinetic equation describing soil adsorption behavior is introduced and coupled with the hydrodynamic module to simulate pollutant migration and transformation. Finally, by setting the numerical solution scheme for the model and using real-time monitoring data for online calibration of key parameters, the model is given the ability to dynamically update and adaptively adjust. This model serves as a computational engine for achieving minute-level, high-precision dynamic simulation and evaluation of the pollutant control effectiveness of buffer zones.

[0126] For step S3: quantizing the buffer band structure-functional response relationship, it is carried out according to the following steps:

[0127] Step S301: Extract and quantify vegetation structure parameters; First, based on the supervised classification results of vegetation types obtained in step S104, extract the category label layer. From this layer, count the number of pixels classified as "herbaceous", "shrub", and "tree". Calculate the ratio of the number of pixels for each vegetation functional type to the total number of pixels in the buffer zone area to obtain the herbaceous percentage, shrub percentage, and tree percentage, expressed as a percentage.

[0128] Step S302: Determine the physical structure parameters of the root system and soil. To obtain root distribution characteristics, at least three representative quadrats are evenly distributed within the buffer zone according to vegetation type and topography. Within each quadrat, undisturbed soil columns at depths of 0 to 60 cm are collected using a root drill, or soil profiles are excavated and sampled layer by layer according to depth. The obtained root-containing soil samples are placed in a sieve and rinsed with running water until the roots are completely separated from the soil. The cleaned, intact root samples are scanned using a flatbed scanner with a resolution of at least 600 dpi to obtain digital images. Using professional root image analysis software, a color threshold is first set to separate the roots from the background to complete image segmentation, and then a skeletonization algorithm is applied to extract the central line skeleton of the root system. Finally, the root length density is calculated based on the total length of the skeleton and the corresponding undisturbed volume of the soil sample, in centimeters per cubic centimeter.

[0129] To obtain soil pore structure characteristics, undisturbed soil columns were collected from the buffer zone and scanned using an X-ray computed tomography (CT) system at a spatial resolution of 50 micrometers, obtaining a series of continuous two-dimensional grayscale slice images. Using image reconstruction software, such as a plugin for Fiji / ImageJ or dedicated CT analysis software, these two-dimensional slices were reconstructed into a three-dimensional volumetric data model of the soil structure. Based on this three-dimensional model, the pore network was extracted using the maximum sphere algorithm. The ratio of the total volume of pores with an equivalent pore diameter greater than 75 micrometers to the total soil pore volume was calculated to obtain the macropore proportion. The number of connection points between pore tunnels per unit volume was used as the coordination number to characterize pore connectivity. The arithmetic mean of the equivalent diameters of all identified pores was calculated to obtain the average pore diameter.

[0130] Step S303: Integrate multi-source structural parameters to construct a feature vector; create a list of structural parameters, and integrate the following parameters in order to form a feature vector:

[0131] 1. Buffer zone width: The average width of the buffer zone perpendicular to the water flow direction is measured using a geographic information system, and the unit is meters.

[0132] 2. Herbal proportion: from step S301, in percentage.

[0133] 3. Shrub percentage: from step S301, in percentage.

[0134] 4. Tree percentage: from step S301, in percentage.

[0135] 5. Root length density: from step S302, in centimeters per cubic centimeter.

[0136] 6. Proportion of macropores in soil: from step S302, in percentage.

[0137] 7. Pore connectivity: Derived from step S302, expressed as coordination number, and is a dimensionless parameter.

[0138] 8. Average pore size: from step S302, in micrometers.

[0139] 9. Soil organic matter content: determined by the potassium dichromate external heating method, and the unit is percentage.

[0140] 10. Soil bulk density: determined by the ring sampler method, in grams per cubic centimeter.

[0141] 11. Leaf area index: a remote sensing inversion product from step S104.

[0142] 12. Vegetation coverage: Remote sensing inversion product from step S104.

[0143] The above 12 parameters are arranged in a fixed order to form a 12-dimensional feature vector, which is used to characterize the comprehensive structure of a single buffer band sample.

[0144] Step S304: Train the structure-function nonlinear mapping model; collect multiple samples from different buffer zones or different scenarios within the same buffer zone. Each sample contains two parts of data: first, the 12-dimensional structural feature vector constructed in step S303; second, the corresponding pollutant reduction efficiency label. The label is obtained as follows: for each buffer zone scenario corresponding to a sample, run the pollutant migration and transformation model corrected in step S2 to simulate the pollutant concentration process line at the downstream outlet, and then calculate the average reduction rates of total nitrogen, total phosphorus, and sediment according to the dynamic reduction rate formula defined in step S4. Divide this average reduction rate by the buffer zone width to obtain the reduction rate per unit width as the label.

[0145] All collected samples are randomly divided into training and test sets, typically with the training set comprising 70% to 80%. Using the structural feature vector as input X and the reduction rate per unit width as output Y, a gradient boosting decision tree model is constructed and trained using Python's scikit-learn library or the GradientBoostingRegressor module. During training, 5-fold cross-validation is used for hyperparameter tuning and performance evaluation within the training set. Key hyperparameters for model training are set as follows: learning rate of 0.1, maximum decision tree depth of 8, minimum number of samples required for leaf nodes of 20, and other parameters using default values. The goal of model training is to minimize the mean squared error between the predicted reduction rate and the true label.

[0146] Step S305: Perform sensitivity analysis and identify key thresholds; based on the trained gradient boosting decision tree model, use the Sobol global sensitivity analysis method to quantify the contribution of each structural feature to the uncertainty of the model's predicted output. Monte Carlo sampling calculations are performed using open-source sensitivity analysis libraries such as SALib, with the number of samplings set to 10,000. First, a reasonable value distribution range is defined for each input feature, as follows:

[0147] Buffer band width Uniformly distributed, with a minimum value of 5 m and a maximum value of 50 m, in accordance with the design specifications for typical ecological buffer zones;

[0148] Herbal proportion Uniform distribution, minimum value 0%, maximum value 100%, based on the actual possible range of values;

[0149] Shrub percentage Uniform distribution, minimum value 0%, maximum value 100%, based on the actual possible range of values;

[0150] Tree proportion Uniform distribution, minimum value 0%, maximum value 100%, based on the actual possible range of values;

[0151] Root length density Normal distribution, mean 2.0 cm / cm³, standard deviation 1.0, cutoff to [0.5, 5.0], based on statistical data from measured data;

[0152] Soil macropore ratio Normal distribution, mean 5%, standard deviation 2%, cutoff to [1%, 15%], based on statistical data from actual measurements;

[0153] Pore ​​connectivity Uniform distribution, minimum value 2, maximum value 10, based on the range reported in the literature;

[0154] Average aperture Uniformly distributed, with a minimum value of 50 μm and a maximum value of 300 μm, based on the range reported in the literature;

[0155] Soil organic matter content : Uniformly distributed, with a minimum value of 0.5% and a maximum value of 8%, based on the typical agricultural soil range;

[0156] Soil bulk density Normal distribution, mean 1.3 g / cm³, standard deviation 0.2, cutoff to [1.0, 1.8], based on statistical analysis of measured data;

[0157] Leaf area index Uniform distribution, minimum value 0.5, maximum value 6.0, based on typical vegetation cover range;

[0158] vegetation coverage Uniform distribution, minimum value 10%, maximum value 95%, based on the actual possible range of values.

[0159] Based on the aforementioned distribution range, the Sobol sampling method is used to generate the input feature matrix, and the trained gradient boosting decision tree model is called for prediction to obtain the model output. The set of simulation results (such as total phosphorus reduction rate). The core of the Sobol index calculation is variance decomposition, a first-order sensitivity index. The calculation formula is:

[0160]

[0161] in, It is the first The first-order Sobol exponent of each input feature; It is the model output. The total variance; It is in the fixed features At that time, for all other features Find the expected value of the condition. Is this conditional expected value relative to The variance. The index. The larger the value, the more characteristic it indicates. The greater the contribution to output uncertainty, the more critical the control factor.

[0162] For the identified key control factors (such as buffer band width) Within its value range (5 m to 50 m), 46 points are uniformly selected with a step size of 1 m, while keeping other feature values ​​fixed as the mean of the training dataset. These points are then substituted into the trained model to predict the corresponding pollutant reduction rate. A response curve between this factor and the reduction rate is plotted. Functional thresholds are identified by analyzing changes in the curve slope: the slope between adjacent points on the curve is calculated. When the absolute value of the slope first falls below 0.01 (i.e., the reduction rate increase is less than 1% / m), the corresponding width value is determined as the saturation width threshold; when the reduction rate begins to fall below 80% of the baseline value, the corresponding parameter value is determined as the functional failure lower limit. For example, based on measured data, when the width increases from 15 m to 16 m, if the sediment reduction rate increase is less than 1%, then 15 m is determined as the saturation width threshold under the current conditions; when the root length density is less than 2.5 cm / cm³, and the total phosphorus reduction rate drops below 80% of the baseline value, then 2.5 cm / cm³ is determined as the functional failure lower limit.

[0163] In summary, step S3 establishes a quantitative response relationship between the buffer zone's "structure" and the pollutant reduction "function." It first systematically acquires multi-dimensional structural parameters from vegetation, roots, and soil, integrating them into a feature vector. Then, a nonlinear machine learning model is trained using a gradient boosting decision tree algorithm to predict reduction efficiency based on structural features. Finally, global sensitivity analysis and threshold identification reveal the key structural factors influencing the function and their effective range. The quantitative relationship established in this step provides a direct scientific basis for subsequent performance evaluation and structure-based optimization design.

[0164] For step S4: Implementing multi-scenario dynamic performance evaluation, the following implementation steps shall be followed;

[0165] Step S401: Define the boundary conditions for multi-scenario evaluation; First, identify the three types of key boundary condition variables that affect the evaluation of buffer zone effectiveness.

[0166] The first type of variable is the rainfall return period, with four typical levels selected: 2 years, 5 years, 10 years, and 20 years. The second type of variable is the initial pollutant load, with three levels: low, medium, and high. The third type of variable is the buffer zone degradation state, defined as healthy and degraded states.

[0167] Step S402: Generate initial conditions for rainfall and pollutants; for the rainfall recurrence period, consult the hydrological manual or design specifications of the target area to obtain the local rainstorm intensity formula.

[0168] Using the rainstorm intensity formula, the design rainstorm intensity corresponding to the 2-year, 5-year, 10-year, and 20-year return periods is calculated respectively, and the corresponding design rainstorm process line is generated, that is, the curve of rainfall intensity changing with time.

[0169] For the initial pollutant load, the typical emission intensity of agricultural non-point sources is used as the benchmark for the medium level.

[0170] For example, the initial concentration of total nitrogen is set at 8.0 mg / L, the initial concentration of total phosphorus is set at 1.2 mg / L, and the initial concentration of sediment is set at 1500 mg / L.

[0171] The initial concentrations for low and high levels can be adjusted according to a preset ratio based on the local actual monitoring data or relevant standards, on the basis of the intermediate level benchmark value.

[0172] Steps S401 and S402 provide quantifiable and reproducible external driving conditions for numerical simulation, including rainfall inputs of different intensities and runoff inputs of different pollution levels.

[0173] Step S403: Define the degradation state of the buffer band and associate it with model parameters; clearly define the quantitative criteria for the degradation state of the buffer band.

[0174] When the vegetation cover of the buffer zone is less than 60%, or the soil bulk density is greater than 1.45 grams per cubic centimeter, the buffer zone is determined to be in a degraded state. In the pollutant migration and transformation model established in step S2, the degraded state is parameterized. Specifically, when the buffer zone is determined to be degraded, the vegetation resistance coefficient, which characterizes vegetation resistance, is reduced by 20%, and the macropore ratio parameter, which characterizes soil infiltration performance, is reduced by 30%.

[0175] Step S404: Construct a multi-scenario combination input set; combine the four rainfall recurrence period scenarios, the three initial load levels of pollutants, and the two buffer zone states defined in step S403 generated in step S402.

[0176] Thus, a total of 24 different evaluation scenarios are formed, which constitute a multi-scenario input set covering a variety of external conditions and internal states.

[0177] Step S405: Run the model to calculate dynamic performance indicators; for each of the 24 evaluation scenarios constructed in step S404, drive the pollutant migration and transformation coupling model corrected by real-time data in step S2 to perform simulation calculations.

[0178] The model simulation outputs the two-dimensional spatial migration path of pollutants within the buffer zone, the residence time distribution in different regions, and the time-series curves of the concentration of various pollutants at the downstream outlet as a function of time under each scenario.

[0179] Based on the simulated upstream inlet concentration process line and downstream outlet concentration process line Considering the transport time (i.e., hydraulic residence time) of water flowing from the upstream inlet to the downstream outlet, concentrations at the same moment cannot be directly compared. Therefore, the average residence time of the water flow within the buffer zone is first calculated using a hydrodynamic model. The calculation formula is as follows:

[0180]

[0181] in, The length of the buffer zone along the direction of water flow is expressed in meters. The average flow velocity within the simulated buffer zone is expressed in meters per second. Residence time. The value of is dynamically determined based on the simulated runoff process, and is usually taken as the weighted average value of the flow velocity within the simulation period.

[0182] Based on stay time By aligning the upstream inflow concentration with the downstream outflow concentration over time, the dynamic reduction rate of pollutants during the entire runoff event can be calculated. :

[0183]

[0184] in, Represents the upstream inlet in time The pollutant concentration, i.e., the residence time of the water entering the buffer zone. Later than time Reaching the downstream outlet; Represents the downstream outlet in time The formula ensures that the comparison involves the concentration changes of the same water body before and after purification through a buffer zone.

[0185] The calculation time window should cover the entire process of a single runoff event, starting from the start of runoff generation and continuing until the runoff volume decreases to 10% of its peak flow. In the early stages of the event (… During the period when the upstream water has not yet reached the downstream outlet, the reduction rate is not calculated or is marked as invalid.

[0186] Step S406: Output multi-scenario assessment results; calculate the dynamic reduction rate throughout the entire runoff event process obtained in step S405. The analysis was conducted to extract and output three types of key performance indicators.

[0187] First is the peak reduction rate, which is the percentage of peak reduction during the entire event. The maximum value reached. Secondly, the average reduction rate, which is the average reduction rate over the entire calculation window. Take the arithmetic mean. Thirdly, reduce the duration, i.e. The length of time that the value remains positive and above a certain set threshold (e.g., 5%).

[0188] Ultimately, for each of the 24 assessment scenarios, the above three performance indicators for the three types of pollutants—total nitrogen, total phosphorus, and sediment—are output respectively.

[0189] In summary, step S4 achieves a comprehensive and dynamic evaluation of the buffer zone's effectiveness by constructing and simulating a systematic set of multiple scenarios. First, three key boundary conditions—rainfall, pollution load, and system state—are identified, and specific quantification methods are provided. Then, these conditions are combined into various scenarios to drive the mechanistic model for simulation. Finally, based on the concentration process lines output by the model, key indicators reflecting the dynamic control performance of the buffer zone are calculated and extracted. This process ensures that the evaluation results are no longer single static values, but rather a spectrum reflecting the buffer zone's performance under different external pressures and its own state, providing rich decision-making information for subsequent comprehensive evaluation and optimization.

[0190] For step S5: generating a generalizable evaluation index system and optimization suggestions, the specific implementation steps include the following:

[0191] Step S501: Define three-level evaluation indicators; integrate the dynamic reduction rate output from step S4 and the structural parameters obtained from step S3 to construct a comprehensive evaluation indicator system containing three primary indicators.

[0192] The first primary indicator is defined as process effectiveness, directly characterized by the dynamic reduction rate calculated in step S4. The second primary indicator is defined as system stability, aiming to reflect the fluctuation of buffer zone effectiveness over time or events. The third primary indicator is defined as structural adaptability, aiming to reflect the potential of the buffer zone's own structure to support ecological functions.

[0193] Step S502: Calculate the system stability index; to quantify system stability, it is necessary to collect pollutant reduction rate data of the target buffer zone in five or more consecutive typical rainstorm events. The standard deviation of the reduction rate data series from these consecutive events is calculated and denoted as . System stability index It is then defined as the reciprocal of the standard deviation, and the calculation formula is:

[0194]

[0195] in, Represents the rate of pollutant reduction. This represents the standard deviation of the reduction rate series. A smaller standard deviation indicates more stable performance. The larger the value, the better.

[0196] Steps S501 and S502 complete the direct quantification of the buffer band's "process performance" and the indirect quantification of its "system stability." The former reflects instantaneous performance, while the latter reflects long-term reliability.

[0197] Step S503: Calculate the structural adaptability score; to quantify structural adaptability, three basic parameters need to be measured.

[0198] The first item is the Shannon diversity index for plants. By conducting quadrat surveys of all plant species within the buffer zone, recording the quantity of each species, and following ecological formulas... Perform calculations, where For the first The relative abundance of each species.

[0199] The second item is soil organic matter content. The percentage by mass is determined using the potassium dichromate external heating method. The third item is root biomass. The unit is the dry weight of the root system per unit area, which is obtained by setting up quadrats in the buffer zone, digging and washing the roots, drying them, and then weighing them.

[0200] because (Shannon Diversity Index) (Soil organic matter content) The three indicators (root biomass) have different dimensions and numerical magnitudes. Directly weighting and summing them would lead to the larger-magnitude indicator dominating the scoring result. Therefore, each indicator is first standardized, mapping them uniformly to the [0,1] interval to eliminate the influence of dimensions. The standardization formula is:

[0201]

[0202] in, These are the original indicator values. and These represent the theoretical minimum and maximum values ​​of the indicator within the evaluation area. The boundary values ​​for each indicator were determined based on statistical analysis of measured data and literature review. The value range is from 0 to 5. The value ranges from 0.5% to 8%. The value ranges from 50 g / m² to 2000 g / m².

[0203] After standardizing the three indicators respectively, we obtain the standardized values. All three are dimensionless numbers, and their values ​​are all in the interval [0,1]. Structural adaptability score Defined as the weighted sum of the three, the calculation formula is:

[0204]

[0205] in, Represent , , The weighting coefficients, the sum of the three is 1.

[0206] Step S504: Determine weights and construct a comprehensive evaluation matrix; structural adaptability scoring. The three sub-indicators (Shannon Diversity Index) Soil organic matter content Root biomass The weighting coefficients of the variables were determined using the analytic hierarchy process (AHP) combined with general knowledge in the field. Based on the consensus in ecological buffer zone function research, vegetation diversity, soil fertility, and root retention capacity contribute similarly to the structural adaptability of buffer zones. Among these, vegetation diversity reflects the system's stability potential, soil organic matter reflects nutrient supply capacity, and root biomass is directly related to soil erosion resistance and porosity development. Based on this, the following judgment matrix was constructed:

[0207] index Its importance relative to itself is 1, and its importance relative to itself is 1. Its importance is 2, relative to The importance level is 1;

[0208] index Compared to Its importance is 1 / 2, its importance relative to itself is 1, and its importance relative to... Its importance is 1 / 2;

[0209] index Compared to Its importance is 1, relative to Its importance is 2, and its importance relative to itself is 1.

[0210] The meaning of this matrix is: and Equal relative importance (value is 1) Compared to Slightly important (value is 2) Compared to Less important (value is 1 / 2). Calculate the largest eigenvalue of this matrix. The weight coefficients are obtained after normalizing the corresponding feature vectors: , , Consistency ratio This satisfies the consistency requirement. It should be noted that the determination of this weighting coefficient is based on the fact that each indicator has been standardized to the [0,1] interval, that is, after eliminating the influence of dimensions, the relative contribution of each indicator to structural adaptability is evaluated.

[0211] The total nitrogen reduction rate, total phosphorus reduction rate, sediment reduction rate, and system stability index were included. Structural adaptability score These five indicators are visualized in the form of a radar chart, forming a pentagonal comprehensive evaluation matrix.

[0212] Step S505: Generate a targeted optimization configuration scheme; based on the scores of each indicator in the comprehensive evaluation matrix and the threshold values ​​of key structural parameters identified in step S3, automatically match and generate an optimization configuration scheme.

[0213] For example, if the assessment shows that the sediment reduction rate is insufficient, but the width of the buffer zone has exceeded the 15-meter saturation threshold identified in step S3, the optimization suggestion is to adjust the vegetation configuration and increase the proportion of shrub vegetation area to more than 40%.

[0214] If the assessment shows that the total phosphorus reduction rate is low and the root length density measured in step S302 is below the functional failure limit of 2.5 cm per cubic centimeter, the optimization suggestion is to replant deep-rooted leguminous plants such as Amorpha fruticosa and consider applying biochar to improve the soil pore structure.

[0215] If the calculated system stability index If the value is below the preset experience threshold of 0.8, the optimization suggestion is to implement vegetation management measures such as segmented rotation mowing to maintain the continuous vitality and function of the vegetation.

[0216] Step S506: Generate a generalizable evaluation report; integrate all processes and results from steps S501 to S505, including raw data, calculation process, indicator scores, visualization matrix and optimization suggestions, to generate a structured evaluation report.

[0217] The standardized report format ensures that its conclusions and recommendations apply not only to the current buffer zone assessment, but also that its methodology and key thresholds can serve as a reference for the assessment and design of buffer zones in other similar areas.

[0218] In summary, step S5 integrates and elevates the dynamic simulation data, structural parameters, and mechanistic analysis results obtained in the preceding steps into an operable and scalable decision support tool. By constructing a multi-dimensional comprehensive evaluation index system, it not only quantifies the instantaneous purification efficiency of the buffer zone but also assesses its long-term operational stability and the adaptability of the structure itself. Finally, the evaluation conclusions are correlated with scientific thresholds in the structure-function relationship to automatically generate targeted engineering and management optimization suggestions, thereby completing a full technical closed loop from dynamic monitoring and mechanistic simulation to precise evaluation and optimization design.

[0219] In addition, this embodiment also includes time-series cluster analysis of historical monitoring data. The specific implementation process is as follows:

[0220] First, complete runoff-pollution event data recorded over a continuous period of three years or longer are extracted from the established monitoring database. For each individual event, the data is compiled into a multivariate time series object. This object contains five parallel, timestamp-aligned series representing the time-varying rainfall intensity, flow rate, total nitrogen concentration, total phosphorus concentration, and sediment concentration during that event.

[0221] Next, to measure the morphological similarity between different event time series, a dynamic time warping algorithm is used as the distance metric. For any two multivariate time series of events, their DTW distance is calculated. Specifically, multivariate DTW implementations provided by libraries such as tslearn in Python or the dtw-python package can be used to simultaneously align the five variable sequences and calculate the cumulative distance. The smaller the distance value, the more similar the overall processes of the two events are.

[0222] Based on the calculated pairwise DTW distance matrix, K-means clustering algorithm is applied to classify all historical events into patterns. Before implementing clustering, the optimal number of clusters K needs to be determined. The specific method is as follows: preset a reasonable range for K, for example, from 2 to 6; for each candidate K value, perform K-means clustering and calculate the average silhouette coefficient of all samples; finally, select the K value that maximizes the average silhouette coefficient as the optimal number of clusters. Based on extensive application experience, this optimal K value is often 4.

[0223] When K=4, the four event clusters generated by clustering can be assigned typical scenario labels based on the morphological characteristics of their cluster center sequences. For example, one cluster center exhibiting high rainfall intensity and short total duration can be interpreted as a "short-duration, high-intensity" event; another cluster center exhibiting moderate rainfall intensity but long duration can be interpreted as a "long-duration, medium-intensity" event; the third cluster may exhibit complex multi-peak pollution load characteristics, interpreted as a "compound" event; and the fourth cluster exhibiting low intensity and long duration characteristics can be interpreted as a "low-intensity, persistent" event. These labels provide semantic explanations for understanding event types.

[0224] Finally, to leverage the clustering results to improve the generalization ability of the model in step S3, the cluster centers need to be converted into usable inputs for the model. The cluster centers of the four clusters are extracted, each center being a representative multivariate time series. Then, for each cluster center sequence, it is simulated as an input to a runoff event, and the mechanistic model from step S2 is run to calculate the comprehensive performance indicators of the buffer zone under this typical event (such as the average reduction rate). Combined with the structural parameters of the buffer zone, expanded "structure-function" training samples are constructed according to the methods described in steps S303 and S304. These new samples, generated from the cluster centers and covering different typical event patterns, are added to the original training set to retrain or fine-tune the GBDT model. This process is equivalent to supplementing the model with knowledge of rainfall-pollution combination patterns that were not directly observed, thereby effectively improving its predictive generalization ability.

[0225] Simultaneously, a regional-scale buffer zone effectiveness database was established, using a PostgreSQL+PostGIS spatial database architecture. Each record contains the following fields: geographic coordinates (POINT type), average annual rainfall (mm), slope (°), buffer zone width (m), plant species (TEXT array), soil type (FAO classification), and pollutant reduction rate (JSONB format, including TN / TP / sediment).

[0226] The database initially contains over 10,000 entries and supports filtering and querying by climate zone (humid / semi-humid / arid) and land use (farmland / orchard / pasture), enabling cross-basin parameter migration. For example, when migrating the optimal width parameter of the farmland buffer zone in the North China Plain to a similar climate zone in Northeast China, only the soil infiltration parameter needs to be adjusted, while the remaining structure-function relationships remain unchanged. The migration accuracy has been verified to be over 85%.

[0227] This embodiment also achieves local intelligent evaluation through an edge computing gateway. The gateway hardware uses an ARM Cortex-A53 quad-core processor with a main frequency of 1.2 GHz, 2 GB of memory, and a built-in lightweight TensorFlow Lite inference engine. The coupled model in step S2 is simplified into a reduced-order model, retaining the core hydrodynamic and adsorption modules, reducing the spatial resolution to 5 meters, and increasing the time step to 30 seconds. The gateway performs a local evaluation every 10 minutes, uploading only key state variables (such as peak reduction rate, residence time, and system stability index) to the cloud platform via the LoRaWAN protocol. The response latency for a single evaluation is less than 30 seconds, meeting near real-time management requirements.

[0228] To verify the effectiveness of this embodiment, a 30-meter-wide grass-shrub composite buffer zone (2.3 km²) in a small agricultural watershed in the middle reaches of the Yangtze River was selected as the test object. During the 2023 flood season, seven effective runoff events were recorded. Monitoring data showed that under a 5-year return period rainstorm, the average reduction rate of total nitrogen was 62.3%, total phosphorus was 78.5%, and sediment was 89.1%. The correlation coefficient R² between the model simulation results and the measured outflow concentrations reached 0.91, with a peak concentration error of <8%. Sensitivity analysis showed that the Sobol first-order exponents for total phosphorus reduction were 0.38 and 0.42, respectively, indicating they were the main control factors. Based on this, the optimization suggestion was to replant *Salix integra* 10 meters downstream of the buffer zone, increasing its shrub proportion from the current 25% to 40%, which is expected to increase the total phosphorus reduction rate to over 85%.

[0229] The foregoing has shown and described the basic principles, main features, and advantages of the present invention. Those skilled in the art should understand that the present invention is not limited to the above embodiments. The embodiments and descriptions in the specification are merely illustrative of the principles of the invention. Various changes and modifications can be made to the invention without departing from its spirit and scope, and all such changes and modifications fall within the scope of the present invention as claimed. The scope of protection of this invention is defined by the appended claims and their equivalents.

[0230] It will be apparent to those skilled in the art that the present invention is not limited to the details of the exemplary embodiments described above, and that the present invention can be implemented in other specific forms without departing from the spirit or essential characteristics of the present invention. Therefore, the embodiments should be regarded as exemplary and non-limiting in all respects.

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

Claims

1. A method for assessing the control of surface runoff pollutants by ecological buffer zones, characterized in that, Includes the following steps: S1: Construct a dynamic monitoring system that integrates in-situ monitoring and remote sensing inversion data to obtain high spatiotemporal resolution hydrological and water quality, vegetation and soil parameter data for the ecological buffer zone area; S2: Establish a dynamic model of pollutant migration and transformation based on the coupling of physical mechanisms and data-driven approaches. The dynamic model of pollutant migration and transformation is coupled with a hydrodynamic module describing shallow water flow on slopes, a resistance module that explicitly quantifies the effect of vegetation resistance, and a dynamic module describing the adsorption behavior of pollutants in soil. Using the real-time upstream inflow water quality, water level, and flow rate data obtained in step S1, the Manning roughness coefficient and soil adsorption rate constant in the model are periodically corrected online using the extended Kalman filter algorithm. S3: Quantify the nonlinear response relationship between the structural parameters of the ecological buffer zone and the pollutant reduction function. Based on the multi-dimensional structural parameters of vegetation, root system and soil obtained in step S1, construct a structural feature vector. Using the structural feature vector as input and the pollutant reduction rate per unit width obtained by the model simulation in step S2 as output, train the gradient boosting decision tree model, establish the nonlinear mapping relationship between structure and function, use the Sobol global sensitivity analysis method to quantify the contribution of each structural feature to the uncertainty of the pollutant reduction rate prediction, identify key control factors, and determine the functional threshold of the key control factors based on the model response curve. S4: Implement multi-scenario dynamic performance evaluation. Define multiple evaluation scenarios by combining different rainfall return periods, initial pollutant loads, and buffer zone degradation states. Drive the pollutant migration and transformation dynamic model corrected in step S2 to simulate the scenario. Calculate and output the dynamic reduction rate of pollutants under each scenario. Extract the peak reduction rate, average reduction rate, and reduction duration from the dynamic reduction rate process as key performance indicators for output. S5: Generate a comprehensive evaluation index system and optimization configuration scheme, integrate the key performance indicators obtained in step S4 and the structural adaptability parameters obtained in step S3. The structural adaptability parameters include the Shannon diversity index of plants, soil organic matter content and root biomass, construct a comprehensive evaluation matrix that includes process performance, system stability and structural adaptability, and generate a targeted ecological buffer zone optimization configuration scheme based on the functional thresholds of the key control factors identified in step S3.

2. The method for assessing the control of surface runoff pollutants by ecological buffer zones according to claim 1, characterized in that, In step S1, constructing a dynamic monitoring system includes: At least three monitoring sections are set up along the main slope of the ecological buffer zone. Each section is equipped with monitoring points for the upstream inlet, the middle transition zone and the downstream outlet to form a three-dimensional monitoring grid. Integrated water quality-hydrology-vegetation sensing units are deployed at each monitoring point to collect data on total nitrogen concentration, total phosphorus concentration, sediment concentration, flow velocity, water level, and rainfall intensity in surface runoff in situ and at high frequency. Acquire high-resolution multispectral or radar remote sensing images covering the target area, and extract spatial distribution data of vegetation cover, leaf area index and soil volumetric water content through supervised classification and inversion models. The in-situ monitoring data and remote sensing inversion data are spatiotemporally aligned and fused to form a multi-source heterogeneous dataset.

3. The method for assessing the control of surface runoff pollutants by ecological buffer zones according to claim 1, characterized in that, In step S2, establishing a dynamic model for pollutant migration and transformation includes: Using the simplified Saint-Venant equations as the framework for the hydrodynamic module, the slope flow motion under the influence of rainfall, infiltration, and topography is simulated. In the hydrodynamic module, the comprehensive resistance slope is decomposed into bed surface resistance and vegetation resistance. The vegetation resistance is dynamically calculated using a formula based on the vegetation stem drag coefficient, stem density, equivalent stem diameter, and hydrodynamic parameters. Convection-diffusion-reaction transport equations for total nitrogen and total phosphorus were established, with second-order kinetic equations used for the source and sink terms to describe the adsorption and desorption processes of pollutants between the solid and liquid phases of the soil. The hydrodynamics module, resistance module and adsorption kinetics module are coupled together and numerical solutions are obtained using the finite volume method. Using the real-time upstream inflow water quality, water level, and flow rate data obtained in step S1, the Manning roughness coefficient and soil adsorption rate constant in the model are periodically corrected online using the extended Kalman filter algorithm.

4. The method for assessing the control of surface runoff pollutants by ecological buffer zones according to claim 1, characterized in that, In step S3, the nonlinear response relationship between the quantization structure parameters and the reduction function includes: Based on the supervised classification results of remote sensing images, the proportion of vegetation functional types of herbaceous plants, shrubs, and trees was statistically analyzed; Root length density, soil macropore ratio, pore connectivity and average pore size were obtained through root scanning analysis and soil CT scanning technology. A multidimensional structural feature vector is constructed by integrating the buffer zone width, the proportion of the vegetation functional type, the root system and soil pore structure parameters, soil organic matter content, soil bulk density, leaf area index and vegetation coverage. Using the structural feature vector as input and the pollutant reduction rate per unit width obtained through model simulation in step S2 as output, a gradient boosting decision tree model is trained to establish a structure-function nonlinear mapping relationship. The Sobol global sensitivity analysis method was used to quantify the contribution of each structural feature to the uncertainty of pollutant reduction rate prediction, identify key control factors, and determine the functional threshold of the key control factors based on the model response curve.

5. The method for assessing the control of surface runoff pollutants by ecological buffer zones according to claim 1, characterized in that, Step S4, implementing multi-scenario dynamic performance evaluation includes: Define multi-dimensional assessment boundary conditions that include different rainfall return periods, different initial pollutant load levels, and the health and degradation status of the buffer zone; Based on the local rainstorm intensity formula, design rainstorm process lines under different return periods are generated, and initial concentrations of low, medium, and high levels of pollutants are set; The degradation state of the buffer zone is parameterized as a reduction in the vegetation resistance coefficient and the proportion of macropores in the soil in the model. The rainfall recurrence period scenario, the initial pollutant load level, and the buffer zone status are fully combined to construct a multi-scenario input set; For each scenario, the corrected model in step S2 is run to simulate the spatiotemporal changes in pollutant concentration, and the dynamic reduction rate of pollutants during the entire runoff event is calculated based on the concentration process lines at the upstream inlet and downstream outlet. The peak reduction rate, average reduction rate, and reduction duration are extracted from the dynamic reduction rate process and output as key performance indicators.

6. The method for assessing the control of surface runoff pollutants by ecological buffer zones according to claim 1, characterized in that, Step S5 involves generating a comprehensive evaluation index system and an optimized configuration scheme, including: Three primary evaluation indicators are defined: process efficiency, system stability, and structural adaptability. Process efficiency is directly characterized by the dynamic reduction rate calculated in step S4. System stability is quantified by calculating the reciprocal of the standard deviation of the reduction rate in multiple consecutive runoff events. Structural adaptability is quantified by calculating the weighted sum of plant Shannon diversity index, soil organic matter content, and root biomass. The weights of each sub-indicator of structural adaptability were determined using the analytic hierarchy process (AHP), and the total nitrogen reduction rate, total phosphorus reduction rate, sediment reduction rate, system stability index, and structural adaptability score were integrated into a visual comprehensive evaluation matrix. Based on the evaluation results of each indicator in the comprehensive evaluation matrix, and combined with the functional thresholds of the key control factors identified in step S3, specific optimization configuration suggestions for vegetation configuration adjustment, plant replanting, soil improvement, or vegetation management are generated.

7. The method for assessing the control of surface runoff pollutants by ecological buffer zones according to any one of claims 1 to 6, characterized in that, This also includes time-series clustering analysis and knowledge transfer of historical monitoring data, specifically: Extract multivariate time series data of multiple complete historical runoff-pollution events; The dynamic time warping algorithm is used to calculate the morphological similarity distance between events, and the K-means clustering algorithm is applied to divide historical events into multiple typical pattern clusters; The cluster centers of each cluster are extracted as representative event sequences, which drive the model in step S2 to generate structural functional expansion samples corresponding to different event patterns, thereby enhancing the generalization ability of the machine learning model in step S3.

8. The method for assessing the control of surface runoff pollutants by ecological buffer zones according to claim 7, characterized in that, Also includes: Establish a regional-scale spatial database of ecological buffer zone effectiveness to store buffer zone structural parameters, environmental driving parameters, and effectiveness indicators under different geographical locations, climatic conditions, and land use types; Based on the aforementioned spatial database, the migration and application of cross-regional buffer zone design parameters are supported. That is, between areas with similar climate and land use, the core parameters of the verified buffer zone structure-function relationship can be migrated, and only the local soil infiltration parameters are adjusted for adaptability.

9. The method for assessing the control of surface runoff pollutants by ecological buffer zones according to claim 1, characterized in that, The online calibration of the model in step S2 and the multi-scenario simulation evaluation in step S4 are achieved through near real-time computing via an edge computing gateway deployed on-site. The edge computing gateway loads a simplified pollutant migration and transformation model, periodically receives monitoring data from step S1 and performs localized performance evaluation calculations, and uploads the core state variables obtained from the evaluation to the cloud platform.

10. The method for assessing the control of surface runoff pollutants by ecological buffer zones according to claim 4, characterized in that, Soil pore structure parameters are obtained through soil CT scanning technology, specifically including: An X-ray computed tomography system was used to scan the undisturbed soil column to obtain a three-dimensional volumetric data model of the soil structure. Based on the aforementioned three-dimensional volumetric data model, the maximum sphere algorithm is applied to extract the pore network; The proportion of large pore volume with an equivalent pore diameter greater than 75 micrometers is calculated as the large pore ratio. The number of pore tunnel connection points per unit volume is counted as the pore connectivity. The arithmetic mean of the equivalent diameters of all identified pores is calculated as the average pore diameter.

Images