A numerical simulation method of coupling vegetation ecology and groundwater

Abstract

The present application relates to the technical field of groundwater simulation, in particular to a groundwater numerical simulation method coupled with vegetation ecology. The method realizes the two-way coupling between groundwater and vegetation ecology through iteration, specifically including: based on the groundwater depth, the vegetation coverage is inverted; then based on the groundwater depth and the vegetation coverage, the actual evapotranspiration is inverted; finally, based on the evapotranspiration, the net recharge of the groundwater model is updated until convergence. The method replaces the complex physical model by establishing a statistical function relationship, and corrects the water level of the drainage and overflow cells in the iteration to ensure the calculation stability.

Description

Technical Field

[0001] This invention relates to the field of groundwater simulation technology, and more specifically to a groundwater numerical simulation method coupled with vegetation ecology. Background Technology

[0002] Numerical simulation of groundwater is a technique that uses mathematical models and numerical methods to simulate groundwater flow and its interaction with the surrounding environment. Its core lies in using numerical methods (such as the finite difference method or the finite element method) to solve the governing equations of groundwater flow based on Darcy's law, thereby simulating and analyzing the dynamic behavior of the groundwater system. MODFLOW, with its flexibility and wide range of applications, is one of the most commonly used numerical simulation software programs for groundwater flow. MODFLOW can effectively simulate groundwater flow processes under steady-state and unsteady-state conditions and can incorporate multiple factors such as atmospheric precipitation, soil, and rivers. However, there is still considerable room for improvement in the coupled simulation of groundwater and vegetation. In arid regions, vegetation plays a crucial role in the groundwater hydrological cycle. On the one hand, groundwater during droughts is a stable water source for deep-rooted plants; on the other hand, changes in vegetation cover affect land surface evapotranspiration, thus indirectly affecting the net recharge of groundwater. Therefore, understanding the interaction between groundwater and vegetation is of great significance for the rational management of water resources and ecological restoration in arid regions.

[0003] Although the interaction between groundwater and vegetation in arid regions is significant, most current groundwater-vegetation coupling models primarily focus on unidirectional coupling, considering only the impact of groundwater depth on vegetation. For example, they employ exponential or Gamma function models to describe the response of vegetation cover and community diversity to groundwater level changes. These unidirectional coupling models neglect the feedback effect of vegetation on groundwater, potentially leading to predictive biases and limiting their accuracy and applicability in eco-hydrological systems. Summary of the Invention

[0004] In view of this, the purpose of this invention is to provide a groundwater numerical simulation method coupled with vegetation ecology to overcome the problems existing in the current technology.

[0005] To achieve the above objectives, the present invention adopts the following technical solution: This application provides a numerical simulation method for groundwater coupled with vegetation ecology, including: Step S1: Obtain geographic information, groundwater hydrological characteristics, and environmental data of the target area, and construct a steady-flow MODFLOW groundwater model for the target area based on the obtained data; Step S2: Obtain remote sensing products of vegetation cover, groundwater depth raster data and actual evapotranspiration raster data of the target area at the same time and spatial resolution, and perform spatial matching to form a standardized dataset. Step S3: Based on the standardized dataset, with groundwater depth as the independent variable and vegetation cover as the dependent variable, construct and calibrate a vegetation cover inversion model. Step S4: Based on the standardized dataset, with groundwater depth and vegetation cover as independent variables and actual evapotranspiration as dependent variable, construct and calibrate an actual evapotranspiration inversion model. Step S5: Based on the simulation results of the MODFLOW groundwater model, calculate the initial groundwater level field and initial net recharge of the iterative model. Step S6: Copy the MODFLOW groundwater model and delete the unconfined evaporation module in the MODFLOW groundwater model to obtain an iterative model for subsequent iterative calculations, and input the initial groundwater level field and initial net recharge into the iterative model. Step S7: Start iterative calculation, run the current iterative model, and obtain the groundwater level field output by the current iterative model; Step S8: Based on the groundwater level field, correct the water level of all overflow cells in the unconfined aquifer level field to the aquifer top elevation of the cell, and fill the water level values ​​of all drained cells in the unconfined aquifer level field using an interpolation method to obtain the corrected unconfined aquifer level field. Step S9: Subtract the corrected unconfined aquifer water level field from the surface elevation to obtain the groundwater depth field, and set all negative values ​​in the groundwater depth field to 0 to ensure that the groundwater level does not exceed the surface, thus obtaining the groundwater depth field of the current iteration. Step S10: Based on the vegetation cover inversion model and the groundwater depth field of the current iteration, calculate the vegetation cover distribution field of the current iteration. Step S11: Input the groundwater depth field and the vegetation cover distribution field of the current iteration into the actual evapotranspiration inversion model, and calculate the actual evapotranspiration field of the current iteration. Step S12: Calculate the net recharge amount for the next iteration of the iterative model based on the actual evapotranspiration field of the current iteration. Use the corrected water level field of the unconfined aquifer calculated in the current iteration of the iterative model as the initial water level field for the next iteration. Then, based on the net recharge amount for the next iteration, return to step S7 and perform iterative calculation with the updated parameters until the iterative model reaches the preset iteration termination condition. Step S13: After the iterative model reaches the preset iteration termination condition, output the groundwater level field, the inverted vegetation cover field, and the actual evapotranspiration field calculated by the iterative model after the iteration ends.

[0006] Furthermore, in the method described above, step S1 includes: Construct the MODFLOW model; Obtain the boundary conditions, mesh generation, aquifer structure, and permeability coefficient of the target region, and input them into the MODFLOW model; Set the groundwater flow process, initial water level conditions, solver, and convergence criterion parameters of the MODFLOW model to complete the construction of the steady-flow MODFLOW groundwater model for the target area.

[0007] Furthermore, in the method described above, step S5 includes: The water level simulated by the MODFLOW groundwater model is corrected to obtain the initial water level of the iterative model. The correction includes: interpolating and filling the water level of the drained cells, and correcting the water level of the overflow cells to the top elevation of the aquifer. By subtracting the groundwater evaporation from the net recharge obtained from the MODFLOW groundwater model simulation using a raster calculation tool, the initial net recharge of the iterative model is obtained.

[0008] Furthermore, in the method described above, step S12 includes: Based on the inverted actual evapotranspiration field, the current net replenishment, and the relaxation iteration coefficients, calculate the net replenishment for the next iteration of the iterative model; The corrected water level field of the unconfined aquifer calculated by the iterative model is used as the initial water level field for the next iteration; Return to step S7 and perform iterative calculations with the updated parameters until the preset iteration termination condition is met.

[0009] Furthermore, in the method described above, the step S7 of returning to the previous step and performing iterative calculations with the updated parameters includes: Before the iterative model performs the next iteration, it is determined whether the current iteration has reached the preset number of iterations or whether the calculation result has met the convergence condition. If the current iteration has reached the preset number of iterations, or the calculation result has met the convergence condition, then the iteration calculation will stop.

[0010] Furthermore, in the method described above, the convergence condition includes: In all non-drained cells, is the maximum absolute value of the water level change between two adjacent iterations less than the first set threshold? The relative change rates of the number of drained cells, the number of overflowed cells, the number of drained cells that became overflowed, and the number of overflowed cells that became drained are all less than the second set threshold.

[0011] The beneficial effects of this invention are as follows: This application establishes a functional relationship between vegetation cover and groundwater depth, and between evapotranspiration and both vegetation cover and groundwater depth, thus achieving bidirectional coupling between vegetation and groundwater. This statistical model-based functional relationship is universal and can flexibly adapt to different eco-hydrological systems and environmental conditions. For example, under different vegetation types and groundwater depths, only the parameters of the empirical function need to be adjusted, without the need to reconstruct a complex physical model. This scalability allows the invention to be widely applied to various eco-hydrological scenarios.

[0012] In this application, the evapotranspiration simulation only requires that the evapotranspiration inversion function be a function of groundwater depth and vegetation cover, a design that offers high flexibility. Both parametric and non-parametric models can be used to describe the relationship between evapotranspiration and groundwater depth and vegetation cover. This flexibility not only improves the model's adaptability but also allows for the selection of appropriate function forms based on the data characteristics of the specific study area, thus more accurately reflecting the actual situation.

[0013] This application simplifies the complex physical mechanisms relating evapotranspiration, vegetation cover, and groundwater depth, reducing reliance on a large number of complex parameters. This lower parameter requirement not only reduces the difficulty of data acquisition but also improves the operability and practicality of the model.

[0014] This application can seamlessly integrate with commonly used groundwater numerical simulation software (such as MODFLOW), fully leveraging the advantages of existing tools. For example, the initial model construction and iterative calculations can both be implemented based on MODFLOW, while the inversion of vegetation cover and evapotranspiration is completed through external statistical models. This compatibility not only reduces the difficulty of model development but also, with the help of MODFLOW's powerful groundwater simulation capabilities, further enhances the comprehensiveness and accuracy of the model of this invention. Furthermore, compatibility with existing software also makes this invention easy to promote and apply.

[0015] This application significantly improves the efficiency and stability of iterative calculations by designing a weighted update of net supply and a multi-index convergence criterion. This efficient iterative convergence mechanism enables the invention to significantly reduce computation time while maintaining accuracy. Attached Figure Description

[0016] To more clearly illustrate the technical solutions in the embodiments of the present invention or the prior art, the drawings used in the description of the embodiments or the prior art 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.

[0017] Figure 1This is a flowchart illustrating a groundwater numerical simulation method coupled with vegetation ecology, provided in an embodiment of the present invention. Figure 2 This is a schematic diagram showing the relationship between the convergence criterion index and the number of iterations in a groundwater numerical simulation method coupled with vegetation ecology provided by an embodiment of the present invention. Figure 3 This is a schematic diagram of the spatial distribution of groundwater depth and vegetation cover in a groundwater numerical simulation method coupled with vegetation ecology provided by an embodiment of the present invention. Detailed Implementation

[0018] To make the objectives, technical solutions, and advantages of this invention clearer, the technical solutions of this invention will be described in detail below. Obviously, the described embodiments are merely some embodiments of this invention, and not all embodiments. Based on the embodiments of this invention, all other implementation methods obtained by those skilled in the art without creative effort are within the scope of protection of this invention.

[0019] Figure 1 This is a schematic flowchart illustrating a groundwater numerical simulation method coupled with vegetation ecology, provided in one embodiment of the present invention. Please refer to [link / reference]. Figure 1 This embodiment may include the following steps: Step S1: Obtain geographic information, groundwater hydrological characteristics, and environmental data of the target area, and construct a steady-flow MODFLOW groundwater model for the target area based on the obtained data; Step S2: Obtain remote sensing products of vegetation cover, groundwater depth raster data and actual evapotranspiration raster data of the target area at the same time and spatial resolution, and perform spatial matching to form a standardized dataset. Step S3: Based on the standardized dataset, with groundwater depth as the independent variable and vegetation cover as the dependent variable, construct and calibrate a vegetation cover inversion model. Step S4: Based on the standardized dataset, with groundwater depth and vegetation cover as independent variables and actual evapotranspiration as dependent variable, construct and calibrate the actual evapotranspiration inversion model. Step S5: Based on the simulation results of the MODFLOW groundwater model, calculate the initial groundwater level field and initial net recharge of the iterative model. Step S6: Copy the MODFLOW groundwater model and delete the unconfined evaporation module in the MODFLOW groundwater model to obtain the iterative model for subsequent iterative calculations, and input the initial groundwater level field and initial net recharge into the iterative model. Step S7: Start iterative calculation, run the current iterative model, and obtain the groundwater level field output by the current iterative model; Step S8: Based on the groundwater level field, correct the water level of all overflow cells in the unconfined aquifer level field to the aquifer top elevation of the cell, and fill the water level values ​​of all drain cells in the unconfined aquifer level field using an interpolation method to obtain the corrected unconfined aquifer level field. Step S9: Subtract the corrected unconfined aquifer water level field from the surface elevation to obtain the groundwater depth field, and set all negative values ​​in the groundwater depth field to 0 to ensure that the groundwater level does not exceed the surface, thus obtaining the groundwater depth field of the current iteration. Step S10: Based on the vegetation cover inversion model and the groundwater depth field of the current iteration, calculate the vegetation cover distribution field of the current iteration. Step S11: Input the groundwater depth field and vegetation cover distribution field of the current iteration into the actual evapotranspiration inversion model, and calculate the actual evapotranspiration field of the current iteration. Step S12: Calculate the net recharge amount for the next iteration of the iterative model based on the actual evapotranspiration field of the current iteration. Use the corrected unconfined aquifer water level field calculated in the current iteration of the iterative model as the initial water level field for the next iteration. Based on the net recharge amount for the next iteration, return to step S7 and perform iterative calculation with the updated parameters until the iterative model reaches the preset iteration termination condition. Step S13: After the iterative model reaches the preset iteration termination condition, output the groundwater level field, the inverted vegetation cover field, and the actual evapotranspiration field calculated by the iterative model after the iteration ends.

[0020] This embodiment achieves bidirectional coupling between groundwater and vegetation ecology through iterative calculation. The core process is as follows: first, vegetation cover is inverted based on groundwater depth; then, actual evapotranspiration is inverted based on depth and cover; finally, the net recharge of the groundwater model is updated based on actual evapotranspiration. This process is iterated until convergence. In this application, a statistical function relationship is established to replace the complex physical model, and the water level in the drainage and overflow cells is corrected during iteration to ensure computational stability.

[0021] It should be noted that the recommended values ​​for the maximum number of iterations N and the head change convergence threshold HCLOSE are 50~500 and 0.01~1m, respectively. Perform the k-th iteration (k=1,2,3…,N) and run the iterative model. The water level calculation results of the iterative model are read, and the water level in the overflow cells of the unconfined aquifer water level field is corrected to the elevation of the top plate of the aquifer. Then, the water level values ​​of the draining cells are filled using interpolation methods (such as linear interpolation, cubic polynomial interpolation, etc.) to obtain the corrected water level of the current iterative model. and the water level of the unconfined aquifer field Subtract the water level of the unconfined aquifer from the surface elevation. Obtain the groundwater depth field Negative values ​​are truncated to 0 to ensure that the groundwater level does not exceed the surface.

[0022] The specific steps for constructing and calibrating the vegetation cover inversion model are as follows: Based on the statistical relationship between vegetation cover and groundwater depth, an empirical function model (such as an exponential function model) is selected to invert vegetation cover. The calculation formula is:

[0023] in, This represents the average vegetation cover when the groundwater level is at an infinite depth. For the range of variation, D is the groundwater depth, which is a parameter reflecting the rate of change in vegetation cover.

[0024] The specific steps for constructing and calibrating the actual evapotranspiration inversion model are as follows: The actual evapotranspiration is inverted using functions relating actual evapotranspiration to groundwater depth and vegetation cover. For example, actual evapotranspiration can be decomposed into two parts: contributions from groundwater and non-groundwater factors, and the following formula can be used for inversion:

[0025] in, Actual evaporation (unit: mm). Where is the maximum actual evapotranspiration (in mm), and D is the groundwater depth (in m). λ represents the evaporation limit depth for groundwater (recommended value range is 1~15m), n is the shape parameter of the Avilyanov formula, generally 1~3, and λ is the proportion of evaporation from non-groundwater sources to total evaporation.

[0026] Preferably, step S1 includes: Construct the MODFLOW model; Obtain the boundary conditions, mesh generation, aquifer structure, and permeability coefficient of the target region, and input them into the MODFLOW model; Set the groundwater flow process, initial water level conditions, solver, and convergence criterion parameters of the MODFLOW model to complete the construction of the steady-flow MODFLOW groundwater model for the target area.

[0027] The construction of the steady-flow MODFLOW groundwater model in this embodiment is based on the geographical information of the target area, the hydrological characteristics of groundwater, and existing environmental data. The specific implementation includes the following steps: (1) Determine the boundary conditions, grid division, aquifer structure, permeability coefficient, and other parameters of the study area, and input these data into the MODFLOW model. (2) Set the main processes of groundwater flow (such as recharge, evaporation, extraction, etc.). (3) Set the initial water level conditions, solver, and convergence criterion parameters, and run the model to obtain the initial water level results. The initial model will provide the starting values ​​for subsequent model coupling and iterative calculations. It should include a net recharge module (RCH) to simulate the groundwater recharge process, and an EVT or ETS module to simulate the groundwater evaporation process.

[0028] Preferably, step S5 includes: The water level simulated by the MODFLOW groundwater model is corrected to obtain the initial water level of the iterative model. The correction includes: interpolating and filling the water level of the drained cells, and correcting the water level of the overflow cells to the top elevation of the aquifer. By subtracting the groundwater evaporation from the net recharge obtained from the MODFLOW groundwater model simulation using a raster calculation tool, the initial net recharge of the iterative model is obtained.

[0029] The initial water level parameters of the iterative model are obtained by correcting the water level of the initial model. Correction includes correcting the water level calculation results for drained and overflow cells. For drained cells where the water level calculation result is missing, an interpolation method (such as linear interpolation, cubic polynomial interpolation, etc.) is selected to fill the missing value using the water level of adjacent valid cells. For overflow cells in the unconfined aquifer water level field, their hydraulic head value is truncated to the top elevation of the aquifer.

[0030] By subtracting the groundwater evaporation from the net recharge of the initial model using a grid-based computational tool, the initial net recharge of the iterative model is obtained. .

[0031] Preferably, step S12 includes: Based on the inverted actual evapotranspiration field, the current net replenishment, and the relaxation iteration coefficients, calculate the net replenishment for the next iteration of the iterative model; The corrected unconfined aquifer water level field calculated by the iterative model is used as the initial water level field for the next iteration; Return to step S7 and perform iterative calculations with the updated parameters until the preset iteration termination condition is met.

[0032] The net replenishment amount for the next iteration (k+1) is calculated based on the actual evapotranspiration field, and then weighted using the iteration update coefficients. The formula for the weighted calculation is as follows:

[0033] in, This is the net supply amount used by the model in the (k+1)th iteration. For relaxation iteration coefficients (recommended values ​​are 0.01 to 0.10), , These represent the net recharge and groundwater depth used in the k-th iteration, respectively. This refers to precipitation.

[0034] Based on the corrected unconfined aquifer water level field of the current iteration , which serves as the initial water level for the next iteration (k+1).

[0035] Preferably, the process returns to step S7, where iterative calculations are performed using the updated parameters, including: Before the iterative model performs the next iteration, it is determined whether the current iteration has reached the preset number of iterations or whether the calculation result has met the convergence condition. If the current iteration has reached the preset number of iterations, or the calculation result has met the convergence condition, then the iteration calculation will stop.

[0036] Preferred convergence conditions include: In all non-drained cells, is the maximum absolute value of the water level change between two adjacent iterations less than the first set threshold? The relative change rates of the number of drained cells, the number of overflowed cells, the number of drained cells that became overflowed, and the number of overflowed cells that became drained are all less than the second set threshold.

[0037] The iterative process dynamically updates the model parameters through repeated execution. During iteration, the changes in the convergence criterion indicators need to be continuously monitored. Iteration terminates when any of the following conditions are met: (1) Reaching the maximum number of iterations When the number of iterations k reaches the preset maximum number of iterations N.

[0038] (2) The convergence criterion index is less than the threshold. The model is considered to have converged when the convergence criterion index of two consecutive iterations meets any of the following conditions.

[0039] Condition 1: In all non-drained cells, the maximum absolute value of the water level change between two consecutive iterations is less than the first set threshold HCLOSE, and the suggested value of HCLOSE is 0.01m.

[0040] Condition 2: The relative rate of change of the number of the following key cells is less than the second set threshold (5% is recommended): Number of dewatered cells (dry_count): The number of dewatered cells where the water level in the cell is lower than the bottom elevation of the specified cell.

[0041] Number of overflow cells (wet_count): The number of overflow cells where the water head is higher than the ground surface.

[0042] Number of cells changing from dewatered to overflow (dry2wet_count): The number of cells changing from dewatered to overflow.

[0043] Number of cells changing from overflow to dewatered (wet2dry_count): The number of cells changing from overflow to dewatered.

[0044] The recommended values for different convergence indicators can use absolute values. For example, the convergence criterion for the maximum water head change is max_chg < HCLOSE, and the recommended value of HCLOSE is 0.01 m. The recommended values can also use the relative change rate, where the relative change rate of the number of dewatered / overflow cells and the number of cells changing from dewatered to overflow / overflow to dewatered between two adjacent iterations is less than 5%.

[0045] In some optional embodiments, the design of the convergence criterion can also: 1. Use the error between the simulation results and the observed data as the convergence criterion, such as the mean square error (MSE) or the mean absolute error (MAE).

[0046] 2. Use statistical test methods (such as t-test, F-test) to judge whether there is a significant difference between the results of two adjacent iterations. The advantage is that it can more scientifically judge the convergence state of the model.

[0047] In specific practice, taking the Hailiutu River Basin in the Ordos Basin, China as the research area, a vegetation-groundwater coupling model in a stable situation is established and calculated. When simulating the stable flow of the groundwater aquifer structure in this basin, the aquifer is divided into two layers. The first layer is the unconfined aquifer, mainly composed of Quaternary sediments; the second layer is the confined aquifer, mainly composed of Cretaceous sandstone.

[0048] The source-sink terms of the MODFLOW model and the corresponding modules are: precipitation infiltration (RCH), river discharge (DRN), phreatic evaporation (ETS). Among them, the precipitation is taken as 360 mm / a, the maximum phreatic evaporation is set to 1055 mm / a, the limit phreatic evaporation depth is set to 2.7 m, and the change of phreatic evaporation with the distance from the evaporation surface follows the Avyryanov formula with n = 2.

[0049] The water level calculation results of the MODFLOW model are corrected to obtain the initial water level of the iterative model Specifically, for missing water level values ​​in the drainage cells, the effective water level calculation results of the adjacent cells are used to fill the missing values ​​using cubic polynomial interpolation. For overflow cells in the unconfined aquifer water level field, the water level is cut off to the top elevation of the aquifer.

[0050] MODFLOW model of net groundwater recharge It is obtained by subtracting groundwater evaporation from precipitation infiltration recharge.

[0051] Creating an iterative model Copy the MODFLOW model, delete the ETS module, and set the initial head of the model to... Set the model's supply amount to .

[0052] The total number of iterations N for the coupled model is set to 300, and the convergence threshold HCLOSE for the maximum head is set to 1m.

[0053] Run the iterative model to obtain the head calculation result of the kth iteration. Correct the water level in the overflow cells of the unconfined aquifer level field to the elevation of the top plate of the aquifer, and use cubic polynomial interpolation to fill the water level values ​​of the draining cells to obtain the corrected water level of the current iterative model. and the water level of the unconfined aquifer field Subtract the water level of the unconfined aquifer from the surface elevation. Obtain the groundwater depth field And negative values ​​are truncated to 0.

[0054] The function used to set the inverted vegetation cover is:

[0055] in, It is 0.277. It is 0.11. It is 0.30.

[0056] The function used to invert actual evapotranspiration is set as follows:

[0057] in, It is 0.0029 m / d. The value is 7 m, and n is 2. For about The polynomial function is:

[0058] Parameters a, b, and c are set to 0.16, 0.025, and 0.263, respectively.

[0059] Set the weighted update coefficients used for weighted updates. =0.05, and the weighted calculation yields the net supply for the next iteration as: .

[0060] Update the iterative model and iterate repeatedly. Analyze the convergence criteria of adjacent iterative models as the number of iterations increases, and analyze the changes of convergence criteria (such as maximum head change, number of drain / overflow cells) with the number of iterations.

[0061] Figure 2 This shows the change of convergence criteria indices with the number of iterations. As the number of iterations increases, the maximum head change is less than 0.3m, which is less than HCLOSE. After more than 100 iterations, the relative change rate of the number of drainage cells and the number of overflow cells is less than 5%, which can be regarded as convergence.

[0062] The groundwater depth and vegetation cover obtained from the final groundwater level calculation of the iterative model are the numerical solutions of the coupled model, and their spatial distribution map is shown below. Figure 3 .

[0063] The same or similar parts in the above embodiments can be referred to each other, and the contents not described in detail in some embodiments can be referred to the same or similar contents in other embodiments.

[0064] It should be noted that in the description of this invention, the terms "first," "second," etc., are used for descriptive purposes only and should not be construed as indicating or implying relative importance. Furthermore, in the description of this invention, unless otherwise stated, "a plurality of" means at least two.

[0065] Any process or method description in the flowchart or otherwise herein can be understood as representing a module, segment, or portion of code comprising one or more executable instructions for implementing a particular logical function or process, and the scope of the preferred embodiments of the invention includes additional implementations in which functions may be performed not in the order shown or discussed, including substantially simultaneously or in reverse order depending on the functions involved, as will be understood by those skilled in the art to which embodiments of the invention pertain.

[0066] In the description of this specification, references to terms such as "one embodiment," "some embodiments," "example," "specific example," or "some examples," etc., indicate that a specific feature, structure, material, or characteristic described in connection with that embodiment or example is included in at least one embodiment or example of the invention. In this specification, the illustrative expressions of the above terms do not necessarily refer to the same embodiment or example. Furthermore, the specific features, structures, materials, or characteristics described may be combined in any suitable manner in one or more embodiments or examples.

[0067] Although embodiments of the present invention have been shown and described above, these embodiments are exemplary and should not be construed as limiting the present invention. Those skilled in the art can make changes, modifications, substitutions and variations to the above embodiments within the scope of the present invention.

Claims

1. A numerical simulation method for groundwater coupled with vegetation ecology, characterized in that, include: Step S1: Obtain geographic information, groundwater hydrological characteristics, and environmental data of the target area, and construct a steady-flow MODFLOW groundwater model for the target area based on the obtained data; Step S2: Obtain remote sensing products of vegetation cover, groundwater depth raster data and actual evapotranspiration raster data of the target area at the same time and spatial resolution, and perform spatial matching to form a standardized dataset. Step S3: Based on the standardized dataset, with groundwater depth as the independent variable and vegetation cover as the dependent variable, construct and calibrate a vegetation cover inversion model. Step S4: Based on the standardized dataset, with groundwater depth and vegetation cover as independent variables and actual evapotranspiration as dependent variable, construct and calibrate an actual evapotranspiration inversion model. Step S5: Based on the simulation results of the MODFLOW groundwater model, calculate the initial groundwater level field and initial net recharge of the iterative model. Step S6: Copy the MODFLOW groundwater model and delete the unconfined evaporation module in the MODFLOW groundwater model to obtain an iterative model for subsequent iterative calculations, and input the initial groundwater level field and initial net recharge into the iterative model. Step S7: Start iterative calculation, run the current iterative model, and obtain the groundwater level field output by the current iterative model; Step S8: Based on the groundwater level field, correct the water level of all overflow cells in the unconfined aquifer level field to the aquifer top elevation of the cell, and fill the water level values ​​of all drained cells in the unconfined aquifer level field using an interpolation method to obtain the corrected unconfined aquifer level field. Step S9: Subtract the corrected unconfined aquifer water level field from the surface elevation to obtain the groundwater depth field, and set all negative values ​​in the groundwater depth field to 0 to ensure that the groundwater level does not exceed the surface, thus obtaining the groundwater depth field of the current iteration. Step S10: Based on the vegetation cover inversion model and the groundwater depth field of the current iteration, calculate the vegetation cover distribution field of the current iteration. Step S11: Input the groundwater depth field and the vegetation cover distribution field of the current iteration into the actual evapotranspiration inversion model, and calculate the actual evapotranspiration field of the current iteration. Step S12: Calculate the net recharge amount for the next iteration of the iterative model based on the actual evapotranspiration field of the current iteration. Use the corrected water level field of the unconfined aquifer calculated in the current iteration of the iterative model as the initial water level field for the next iteration. Then, based on the net recharge amount for the next iteration, return to step S7 and perform iterative calculation with the updated parameters until the iterative model reaches the preset iteration termination condition. Step S13: After the iterative model reaches the preset iteration termination condition, output the groundwater level field, the inverted vegetation cover field, and the actual evapotranspiration field calculated by the iterative model after the iteration ends.

2. The method according to claim 1, characterized in that, Step S1 includes: Construct the MODFLOW model; Obtain the boundary conditions, mesh generation, aquifer structure, and permeability coefficient of the target region, and input them into the MODFLOW model; Set the groundwater flow process, initial water level conditions, solver, and convergence criterion parameters of the MODFLOW model to complete the construction of the steady-flow MODFLOW groundwater model for the target area.

3. The method according to claim 2, characterized in that, Step S5 includes: The water level simulated by the MODFLOW groundwater model is corrected to obtain the initial water level of the iterative model. The correction includes: interpolating and filling the water level of the drained cells, and correcting the water level of the overflow cells to the top elevation of the aquifer. By subtracting the groundwater evaporation from the net recharge obtained from the MODFLOW groundwater model simulation using a raster calculation tool, the initial net recharge of the iterative model is obtained.

4. The method according to claim 3, characterized in that, Step S12 includes: Based on the inverted actual evapotranspiration field, the current net replenishment, and the relaxation iteration coefficients, calculate the net replenishment for the next iteration of the iterative model; The corrected water level field of the unconfined aquifer calculated by the iterative model is used as the initial water level field for the next iteration; Return to step S7 and perform iterative calculations with the updated parameters until the preset iteration termination condition is met.

5. The method according to claim 4, characterized in that, The return step S7, which involves iterative calculations using the updated parameters, includes: Before the iterative model performs the next iteration, it is determined whether the current iteration has reached the preset number of iterations or whether the calculation result has met the convergence condition. If the current iteration has reached the preset number of iterations, or the calculation result has met the convergence condition, then the iteration calculation will stop.

6. The method according to claim 5, characterized in that, The convergence conditions include: In all non-drained cells, is the maximum absolute value of the water level change between two adjacent iterations less than the first set threshold? The relative change rates of the number of drained cells, the number of overflowed cells, the number of drained cells that became overflowed, and the number of overflowed cells that became drained are all less than the second set threshold.

Images