Surface water-underground water coupling simulation method considering high-river-loss river characteristics

Abstract

The invention discloses a surface water-underground water coupling simulation method considering high-river-loss river characteristics, belongs to the technical field of surface water-underground water coupling simulation calculation, and can solve the problem of low simulation accuracy caused by large deviation between a simulation process and an actual river evolution process in the prior art. The method comprises the following steps: S1, determining a corresponding relation between the phreatic water evaporation capacity and the underground water burial depth of each sub-region in a plurality of sub-regions with different land surface coverage types in the area where the diversion is located according to hydrological geographic information of the area where the diversion of the main stream is located; s2, constructing an underground water numerical simulation model of the area where the diversion is located according to the corresponding relation and the hydrological geographic information of the area where the diversion is located; s3, according to the water diversion amount of the main stream and the diversion conditions of the diversion, the underground water level of the area where the diversion is located is simulated through the underground water numerical simulation model. The device is used for simulating the drainage basin groundwater level.

Description

Technical Field

[0001] The invention relates to a surface water-groundwater coupling simulation method taking into account high river loss river characteristics, and belongs to the technical field of surface water-groundwater coupling simulation calculation. Background Art

[0002] Simulating the water exchange process between surface water and groundwater, as well as the changes in groundwater level or groundwater depth under different water inflow conditions, is of great significance for studying the hydrological conditions and hydrological cycle processes in the basin, rationally planning and managing water resources, and realizing the conservation and intensive utilization of water resources.

[0003] Especially for inland river basins in extremely arid areas, the meteorological conditions of scarce precipitation and strong evaporation have created a unique desert riverside forest vegetation ecosystem and a fragile ecological environment in the basin area. The protection and restoration of the ecological environment in this area depends on the groundwater depth under certain conditions so that the roots of vegetation can absorb groundwater to maintain its growth, and the maintenance of groundwater depth in the region depends on the leakage and recharge of river water in the region. Therefore, it is particularly important to simulate the groundwater depth in the region under different water conditions in the river for ecological protection in extremely arid areas.

[0004] However, existing simulation methods often do not consider the impact of changes in surface cover types and groundwater depth within the basin on groundwater evaporation, nor do they consider the various water losses of incoming water. This causes a large deviation between the simulation process and the actual evolution of the basin, resulting in low simulation accuracy of existing technologies. Summary of the invention

[0005] The present invention provides a surface water-groundwater coupling simulation method taking into account the characteristics of high river loss rivers, which can solve the problem of low simulation accuracy caused by the large deviation between the simulation process of the prior art and the actual evolution process of the river.

[0006] The present invention provides a surface water-groundwater coupling simulation method considering the characteristics of a high river loss river, the method comprising: S1. According to the hydrological and geographical information of the area where the diversion of the main stream is located, determining the corresponding relationship between the phreatic evaporation amount and the groundwater depth of each sub-area in a plurality of sub-areas with different surface cover types in the area where the diversion is located; S2. constructing a groundwater numerical simulation model of the area where the diversion is located according to the corresponding relationship and the hydrological and geographical information of the area where the diversion is located; S3. According to the divertable water volume of the main stream and the diversion conditions of the diversion, the groundwater level in the diversion area is simulated using the groundwater numerical simulation model.

[0007] Optionally, before S3, the method further includes: The amount of water that can be diverted from the main stream is determined based on the hydrological and geographical information of the main stream.

[0008] Optionally, the determining the divertable water volume of the main stream according to the hydrological and geographical information of the main stream is specifically: Determine the water loss of the river channel of the main stream according to the hydrological and geographical information of the main stream; The amount of water that can be diverted from the main stream is determined based on the amount of water loss in the river channel and the hydrological and geographical information of the main stream.

[0009] Optionally, the water loss in the river channel includes water leakage from the river channel, water evaporation from the river channel and water overflow from the river channel.

[0010] Optionally, S3 is specifically: S31, determining the river water levels at a plurality of control sections on the diversion according to the diversion amount of the main stream and the diversion conditions of the diversion; S32. According to the river water levels at the multiple control sections on the diversion, the groundwater level in the area where the diversion is located is simulated using the groundwater numerical simulation model.

[0011] Optionally, the S31 is specifically: S311, determining the flow rates at multiple control sections on the diversion according to the divertable water volume of the main stream and the diversion conditions of the diversion; S312, determining the river channel water levels at the multiple control sections on the diversion according to the flow rates at the multiple control sections on the diversion and the shapes of the multiple control sections.

[0012] Optionally, the branch flow includes a main flow and multiple branch flows; S312 is specifically: Determine the flow rates at multiple control sections on the main stream according to the divertable water volume of the main stream and the diversion conditions of the diversion; The flow splitting ratio of each branch flow from the main flow is determined, and the flow rates at the multiple control sections on each branch flow are determined according to the flow splitting ratio and the flow rates at the multiple control sections on the main flow.

[0013] Optionally, after S3, the method further includes: S4. Calibrate the groundwater numerical simulation model using the historical groundwater level in the area where the diversion is located.

[0014] Optionally, the S4 is specifically: S41, determining a simulation error of the groundwater numerical simulation model according to the simulated groundwater level and the corresponding historical groundwater level; S42: When the simulation error does not satisfy a preset error range, adjusting the calculation parameters of the groundwater level until the simulation error satisfies the preset error range.

[0015] Optionally, the ground cover type is tree forest, shrub forest, grassland or bare land.

[0016] The beneficial effects that the present invention can produce include: The present invention can effectively improve the calculation accuracy of the phreatic water evaporation by determining the corresponding relationship between the phreatic water evaporation and the groundwater depth under different surface cover types in the diversion area. At the same time, the present invention constructs a groundwater numerical simulation model based on the above-mentioned corresponding relationship, and calculates the groundwater level in the diversion area according to the divertable water volume of the main stream and the diversion conditions of the diversion, thereby realizing the simulation of the groundwater level under different main stream water inflow conditions and different diversion water diversion conditions, which can effectively reduce the deviation between the simulation process and the actual evolution process of the river and improve the accuracy of the simulation.

[0017] The present invention determines the channel water loss of the mainstream by calculating the channel seepage water volume, channel evaporation water volume and channel overflow water volume of the mainstream, and determines the divertable water volume of the mainstream according to the channel loss water volume and hydrological and geographical information of the mainstream, which can improve the calculation accuracy of the divertable water volume of the mainstream, thereby further improving the accuracy of the simulation. BRIEF DESCRIPTION OF THE DRAWINGS

[0018] Figure 1 A flow chart of a surface water-groundwater coupling simulation method considering high river loss river characteristics provided by an embodiment of the present invention; Figure 2 A schematic diagram of the surface water-groundwater system recharge process provided by an embodiment of the present invention; Figure 3 A schematic diagram of the surface water-groundwater system drainage process provided by an embodiment of the present invention; Figure 4 A schematic diagram of the operation of a surface water-groundwater system coupling simulation system provided by an embodiment of the present invention; Figure 5 A schematic diagram of the branch flow splitting ratio provided in an embodiment of the present invention; Figure 6 A schematic diagram of the distribution of the control section behind the Xiahelike diversion sluice provided in an embodiment of the present invention; Figure 7 It is a schematic diagram of control sections Ⅰ-Ⅲ; Figure 8 It is a schematic diagram of the control section IV-VI; Fig. 9 This is a schematic diagram of the control section VII-VIII; Fig.10 The simulated regional submergence depth provided by the embodiment of the present invention With parameters Schematic diagram of the fitting; Fig.11A schematic diagram of fitting of groundwater evaporation and groundwater burial depth for multiple sub-regions with different surface cover types provided in an embodiment of the present invention, wherein (a) is a schematic diagram of fitting of tree forest, (b) is a schematic diagram of fitting of shrub forest, and (c) is a schematic diagram of fitting of grassland; Fig.12 A schematic diagram of the initial flow field setting of the model provided in an embodiment of the present invention; Fig.13 A schematic diagram of the initial parameter partitioning of the model provided by an embodiment of the present invention; Fig.14 A schematic diagram showing the comparison between the measured values ​​and the simulated values ​​of the groundwater levels of Alar Well 1, Alar Well 6 and Shazihekou Well 1 during the model calibration period provided by an embodiment of the present invention; Fig.15 A schematic diagram showing the comparison between the measured values ​​and the simulated values ​​of the groundwater levels of Shazihekou Well 6, Wusman Well 1 and Wusman Well 6 during the model calibration period provided in an embodiment of the present invention; Fig.16 A schematic diagram showing a comparison between the measured values ​​and the simulated values ​​of the groundwater levels of Aqike Well 1 and Aqike Well 6 during the model calibration period provided by an embodiment of the present invention; Fig.17 A schematic diagram showing the comparison between the measured values ​​and the simulated values ​​of the groundwater levels of Alar Well 1, Alar Well 6 and Shazihekou Well 1 during the model validation period provided by an embodiment of the present invention; Fig.18 A schematic diagram showing a comparison between the measured values ​​and the simulated values ​​of the groundwater levels of Shazihekou Well 6, Wusman Well 1 and Wusman Well 6 during the model validation period provided in an embodiment of the present invention; Fig.19 A schematic diagram showing a comparison between the measured values ​​and the simulated values ​​of the groundwater levels of Aqike Well 1 and Aqike Well 6 during the model validation period provided by an embodiment of the present invention; Fig. 20 A schematic diagram of the overall groundwater recharge and discharge characteristics and groundwater storage volume changes in the upper and middle reaches of the Tarim River mainstream provided by an embodiment of the present invention; Fig.21 A schematic diagram of groundwater recharge and discharge characteristics and groundwater storage volume changes in the Aral-Xinqiman River section provided by an embodiment of the present invention; Fig. 22 A schematic diagram of groundwater recharge and discharge characteristics and groundwater storage volume changes in the Xinqiman-Yingbaza River section provided in an embodiment of the present invention; Fig.23 A schematic diagram of groundwater recharge and discharge characteristics and groundwater storage volume changes in the Yingbaza-Usman River section provided by an embodiment of the present invention; Fig.24 A schematic diagram of groundwater recharge and discharge characteristics and groundwater storage volume changes in the Usman-Qala River section provided in an embodiment of the present invention. DETAILED DESCRIPTION

[0019] The present invention is described in detail below in conjunction with embodiments, but the present invention is not limited to these embodiments.

[0020] The embodiment of the present invention provides a surface water-groundwater coupling simulation method considering the characteristics of a high river loss river, such as Figure 1 As shown, the method includes: S1. According to the hydrological and geographical information of the area where the diversion of the mainstream is located, determine the corresponding relationship between the phreatic evaporation and the groundwater depth of each sub-area in the multiple sub-areas with different surface cover types in the area where the diversion is located; S2. Constructing a groundwater numerical simulation model for the diversion area based on the corresponding relationship and the hydrological and geographical information of the diversion area; S3. Based on the divertable water volume of the main stream and the diversion conditions of the diversion, the groundwater level in the diversion area is simulated using a groundwater numerical simulation model.

[0021] Specifically, the main stream can be connected to multiple branch streams, and each branch stream can draw water from the main stream.

[0022] Specifically, hydrogeographic information may include river channel information, land cover type, land use, groundwater level, river channel water level, hydrogeological conditions, etc. in the diversion area. River channel information may include information such as the length, width, depth and cross-sectional shape of each part of the river channel; land cover type may be tree forest, shrub forest, grassland or bare land.

[0023] Specifically, the hydrological and geographical information can be determined based on Google satellite images, water body remote sensing images, land use images, and other data of the diversion area. Among them, the water area images should be selected from the years with large water inflow in the mainstream and diversion in recent years to ensure that the river channel information of the diversion area is extracted to the maximum extent; other images should be selected as close to the year as possible to ensure the consistency of all hydrological and geographical information in time, thereby ensuring the timeliness of the construction of the groundwater numerical simulation model.

[0024] Specifically, the diversion conditions may include the diversion gate setting position, diversion flow rate, diversion water diversion demand, etc. In practice, the diversion flow rate can be controlled by adjusting the opening of the diversion gate. The diversion conditions can be specifically set according to the maximum opening of the diversion gate and the diversion water diversion demand.

[0025] Specifically, the operation process of the surface water-groundwater system mainly includes the recharge process and the discharge process, such as Figure 2 and Figure 3 As shown, the recharge process mainly includes seepage recharge of the mainstream river channel and seepage recharge of the diversion channel, and the discharge process mainly includes phreatic evaporation, etc.

[0026] Specifically, the existing groundwater numerical simulation model usually gives the phreatic evaporation rate when calculating the phreatic evaporation. Before the phreatic depth reaches the limit evaporation depth, the phreatic evaporation rate is constant. This calculation method obviously ignores the impact of the change of surface cover type and groundwater depth in the basin on phreatic evaporation, resulting in a large deviation between the calculation result of phreatic evaporation and the actual evaporation situation. However, this embodiment can effectively improve the calculation accuracy of phreatic evaporation by determining the corresponding relationship between phreatic evaporation and groundwater depth under different surface cover types in the diversion area.

[0027] In this embodiment, S1 may specifically be: S11, according to the hydrological and geographical information of the area where the diversion of the main stream is located, determining a plurality of phreatic water evaporation data and a plurality of groundwater depth data for each of a plurality of sub-areas with different surface cover types in the area where the diversion is located; the plurality of phreatic water evaporation data and the plurality of groundwater depth data for each sub-area correspond one to one; S12. According to the multiple phreatic water evaporation data and the multiple groundwater depth data of each sub-region, a polynomial fitting method is used to determine the corresponding relationship between the phreatic water evaporation of each sub-region and its groundwater depth.

[0028] In this embodiment, S2 may specifically be: According to the above correspondence and the hydrological and geographical information of the diversion area, a groundwater numerical simulation model of the diversion area was constructed based on Darcy's law and Flopy.

[0029] Specifically, the finite difference method is used to solve the groundwater numerical simulation model.

[0030] In this embodiment, constructing a groundwater numerical simulation model specifically includes establishing a groundwater movement control equation, and performing spatial discretization, boundary setting, time discretization, and initial parameter setting on the groundwater numerical simulation model according to the hydrological and geographical information of the diversion area.

[0031] In this embodiment, before S3, the method may further include: The amount of water that can be diverted from the main stream is determined based on the hydrological and geographical information of the main stream.

[0032] In this embodiment, the above-mentioned determination of the divertable water volume of the main stream according to the hydrological and geographical information of the main stream may specifically be: Determine the water loss of the main stream according to the hydrological and geographical information of the main stream; The amount of water that can be diverted from the main stream is determined based on the amount of water loss in the river channel and the hydrological and geographical information.

[0033] Specifically, the water loss in a river channel may include water leakage from the river channel, water evaporation from the river channel and water overflow from the river channel.

[0034] Specifically, the amount of water seepage in the mainstream river channel can be calculated based on Darcy's law, and the amount of water overflow in the mainstream river channel can be calculated using a polynomial fitting method.

[0035] Specifically, this embodiment constructs a mainstream river evolution model, which can calculate the channel seepage, channel evaporation and channel overflow of the mainstream according to the hydrological and geographical information of the mainstream, and then obtain the divertable water volume of the mainstream.

[0036] This embodiment determines the channel water loss of the mainstream by calculating the channel seepage water volume, channel evaporation water volume and channel overflow water volume of the mainstream, and determines the divertable water volume of the mainstream according to the channel loss water volume and hydrological and geographical information of the mainstream, which can improve the calculation accuracy of the divertable water volume of the mainstream, thereby further improving the simulation accuracy of the groundwater numerical simulation model.

[0037] In this embodiment, S3 may specifically be: S31, determining the river water levels at multiple control sections on the diversion according to the diversion amount of the main stream and the diversion conditions of the diversion; S32. According to the river water levels at multiple control sections on the diversion, the groundwater level in the diversion area is simulated using a groundwater numerical simulation model.

[0038] Further, S31 may specifically be: S311, determining the flow rates at multiple control sections on the diversion according to the diversion amount of the main stream and the diversion conditions of the diversion; S312. Determine the river channel water levels at the multiple control sections on the diversion according to the flow rates at the multiple control sections on the diversion and the shapes of the multiple control sections.

[0039] In this embodiment, the diversion includes a main flow and multiple branch flows; the flow rates at multiple control sections on the diversion include the flow rates at multiple control sections on the main flow and the flow rates at multiple control sections on each branch flow.

[0040] Therefore, S311 may specifically be: Determine the flow rates at multiple control sections on the main stream based on the divertable water volume of the main stream and the diversion conditions of the diversion; The diversion ratio of each branch flow from the mainstream is determined, and the flow rates at multiple control sections on each branch flow are determined according to the diversion ratio and the flow rates at multiple control sections on the mainstream.

[0041] Specifically, this embodiment uses the Muskingum-Konji method to construct a diversion river evolution model. The diversion river evolution model can calculate the flow at multiple control sections on the main stream and the flow at multiple control sections on each diversion according to the diversion amount of the main stream and the diversion conditions of the diversion, combined with the diversion ratio of the diversion.

[0042] Specifically, this embodiment constructs a river surface water flow-water level conversion model in the diversion area based on the Manning formula. The river surface water flow-water level conversion model can convert the river flow at multiple control sections into river water levels according to the shapes of multiple control sections on the diversion.

[0043] In this embodiment, S312 may specifically be: According to the flow rates at multiple control sections on the diversion and the shapes of multiple control sections, and using the river surface water flow-water level conversion model in the diversion area, the river water levels at multiple control sections on the diversion are determined.

[0044] In summary, this embodiment constructs four models, namely, the mainstream river evolution model, the diversion river evolution model, the river surface water flow-water level conversion model in the diversion area, and the groundwater numerical simulation model in the diversion area. Among them, the mainstream river evolution model is used to calculate the river water loss (including river channel seepage water, river channel evaporation water and river channel overflow water) of each river section on the mainstream according to the hydrological and geographical information of the mainstream, and calculate the flow of water flowing to the diversion gate of the diversion, so as to determine the amount of water that can be diverted from the mainstream at the diversion gate. The diversion river evolution model is used to calculate the flow at multiple control sections on the diversion according to the diversion water of the mainstream and the diversion conditions of the diversion, combined with the diversion ratio of the branch. The river surface water flow-water level conversion model in the diversion area is used to convert the river flow at multiple control sections into river water levels according to the shapes of multiple control sections on the diversion. The groundwater numerical simulation model of the diversion area is used to simulate the groundwater level in the diversion area according to the river water level at multiple control sections on the diversion. It can be seen that the output of the previous model is successively used as the input of the next model, and the final output result is the simulation result of this embodiment.

[0045] In practice, if the models are independent of each other and data cannot be transmitted, it is necessary to frequently process the input and output data of each model manually or with the help of other tools, which reduces the model's operational efficiency and makes the model unsuitable for large-scale computing scenarios. Therefore, this embodiment couples the four models to form a complete surface water-groundwater system coupling simulation system, which can realize the automatic transmission of input and output data, greatly improve the operational efficiency of each model, and enable this embodiment to adapt to large-scale computing scenarios, which is conducive to the efficient simulation of groundwater levels under different main stream water inflow conditions and different diversion water diversion conditions. The operation schematic diagram of the surface water-groundwater system coupling simulation system is shown in the figure. Figure 4 shown.

[0046] Specifically, this embodiment uses the Python platform to couple four models.

[0047] Accordingly, all four models can be built on the Python platform, or imported into the Python platform after being built on other platforms to facilitate coupling on the Python platform.

[0048] In this embodiment, after S3, the method may further include: S4. Calibrate the groundwater numerical simulation model using the historical groundwater levels in the diversion area.

[0049] Furthermore, S4 may specifically be: S41, determining a simulation error of a groundwater numerical simulation model according to the simulated groundwater level and the corresponding historical groundwater level; S42. When the simulation error does not meet the preset error range, adjust the calculation parameters of the groundwater level until the simulation error meets the preset error range.

[0050] Among them, the preset error range can be specifically set according to simulation requirements.

[0051] Specifically, this embodiment adopts a parameter estimation (PEST) automatic calibration method, takes the groundwater level as a calibration indicator, and adjusts the calculation parameters of each model so that the simulation error meets the preset error range.

[0052] For example, the calculation parameters may include the permeability coefficient of the main stream and the branch stream, the hydraulic conductivity of the riverbed, the water supply degree, etc.

[0053] Specifically, this embodiment selects three error indicators, namely, Nash efficiency index (Nash efficiency index, referred to as NSE), determination coefficient ( ) and Mean Absolute Error (MAE). Using the three error indicators, three corresponding simulation errors can be calculated. All three simulation errors are dimensionless values. The value range of is [0, 1]. The closer the simulation error is to 1, the better the correction effect is. The smaller the value of MAE is, the better the correction effect is.

[0054] Next, this embodiment takes the Tarim River Basin as an example to specifically illustrate the surface water-groundwater coupling simulation method considering the characteristics of a high river loss river described in this embodiment.

[0055] The Tarim River Basin is a typical inland river basin in arid areas. The main stream of the Tarim River includes the Aral-Xinqiman River section, the Xinqiman-Yingbaza River section, the Yingbaza-Usman River section, and the Usman-Qala River section. Since the main stream of the Tarim River does not produce flow itself, its runoff depends on the discharge of the three upstream sources. Therefore, the Tarim River Basin is a runoff dissipation zone. Due to the unique climatic conditions and soil geological distribution of the basin, a unique desert oasis desert riverside forest ecosystem has been formed in the basin. The growth and recovery of vegetation depends on the groundwater depth under certain conditions. In order to protect the ecosystem of the basin, it is necessary to accurately simulate the changes in groundwater levels under different surface water inflow conditions.

[0056] Specifically, this embodiment takes the mainstream and diversion areas of the Tarim River Basin as simulation areas to construct a mainstream river evolution model, a diversion river evolution model, a river surface water flow-water level conversion model in the diversion area, and a groundwater numerical simulation model in the diversion area. The calibration period and verification period of the model are from January 2013 to December 2013 and from January 2014 to December 2014, respectively, and the time scale is month.

[0057] Specifically, the calculation process included in each model is as follows.

[0058] 1. Mainstream river evolution model 1. Calculate the amount of river leakage in the main stream according to Darcy's law.

[0059] Specifically, the calculation formula for the amount of river leakage in each river section on the mainstream is: (1) In formula (1), is the amount of water leakage in the river section, is the permeability coefficient of the river section (m / d), is the width of the river channel, is the length of the river section, is the river water level of the river section, is the groundwater level in the river section, The length of the river section that replenishes groundwater with river water.

[0060] in, It can be determined by referring to the relevant data of the basin. and All of them can be obtained from the hydrological and geographical information of the mainstream.

[0061] 2. Calculate the water area based on the water surface width and river length of each river section on the mainstream, and then determine the amount of evaporation water in each river section based on the observed evaporation pan evaporation and evaporation conversion coefficient.

[0062] Specifically, the calculation formula for the river width of each river section on the mainstream is as follows: (2) (3) (4) (5) In formula (2) to formula (5), They are the water surface widths of Alaer Section, Xinqiman Section, Yingbaza Section, and Wusman Section respectively; The water inflow of Alar hydrological station, Xinqiman hydrological station, Yingbaza hydrological station and Wusiman hydrological station (×10 8 m 3 ).

[0063] Specifically, the calculation formula for the evaporation of water in each river section on the mainstream is: (6) (7) (8) (9) The meanings of the parameters in equations (6) to (9) are as follows: is the amount of water evaporated from the Aral-Xinqiman River section. are the water surface widths of the Alar section and the Xinqiman section respectively, is the length of the Aral-Xinqiman River section, The pan evaporation observed in the Alaer-Xinqiman River section is is the evaporation conversion coefficient of the Alaer-Xinqiman River section, This is the precipitation in the area where the Aral-Xinqiman River section is located.

[0064] The amount of water evaporated from the Xinqiman-Yingbaza River section. are the water surface widths of the Xinqiman section and the Yingbaza section, is the length of the Xinqiman-Yingbaza River section, The pan evaporation observed in the Xinqiman-Yingbaza River section is is the evaporation conversion coefficient of the Xinqiman-Yingbaza River section, This is the precipitation in the area where the Xinqiman-Yingbaza River section is located.

[0065] is the amount of evaporation water in the Yingbaza-Usman River section, are the water surface widths of the Yingbaza section and the Wusman section, is the length of the Yingbaza-Usman River section, The pan evaporation observed in the Yingbaza-Usman River section is is the evaporation conversion coefficient of the Yingbaza-Usman River section, It is the precipitation in the area where the Yingbaza-Usman River section is located.

[0066] is the amount of water evaporated from the Usman-Chala River section. is the water surface width of the Wusman section, is the length of the Usman-Qala River section, is the pan evaporation observed in the area where the Usman-Qala River section is located, is the evaporation conversion coefficient of the area where the Usman-Qala River section is located, It is the precipitation in the area where the Usman-Chala River section is located.

[0067] in, The value of is 0.67, The value of is 0.72, The value of is 0.65, The value of is 0.57.

[0068] 3. Calculate the overflow volume of the mainstream river based on the amount of river seepage and evaporation.

[0069] Since the overflow water volume of the river channel cannot be directly measured, this embodiment first determines the river channel loss water volume of each river section under different water inflow conditions, that is, the total water loss volume, based on the water balance principle, combined with the water inflow volume of the upper and lower sections of each river section and the water diversion volume within the river section, so as to obtain the river channel overflow water volume of each river section under different water inflow conditions. Then, according to the river channel overflow water volume of each river section under different water inflow conditions, a polynomial fitting method is used to determine the calculation formula for the overflow water volume of each river section.

[0070] Specifically, the calculation formula for the water loss of each river section on the mainstream is: (10) In formula (10), is the water loss in the river section, is the water inflow of the upper section of the river, is the water inflow of the lower section of the river, It is the amount of water diverted within the river section.

[0071] Specifically, the formula for calculating the river overflow volume based on the river loss volume is: (11) In formula (11), is the overflow volume of the river section, is the water loss in the river section, is the amount of water evaporated from the river channel in the river section, It is the amount of water leakage in the river section.

[0072] Then, this embodiment adopts a polynomial fitting method to fit the river overflow water volume of each river section under different water inflow conditions, and the calculation formulas of the river overflow water volume of each river section after fitting are respectively: (12) (13) (14) (15) In formula (12) to formula (15), are the river overflow volumes of the Aral-Xinqiman River section, Xinqiman-Yingbaza River section, Yingbaza-Usman River section, and Usman-Qiala River section after fitting, They are the water inflow of Aral hydrological station, Xinqiman hydrological station, Yingbaza hydrological station and Usman hydrological station.

[0073] Then, the mainstream river evolution model can determine the amount of water that can be diverted based on the water loss and water inflow of the mainstream of the Tarim River.

[0074] 2. Diversion river evolution model 1. Use the momentum balance method to calculate the diversion ratio of the branch flow in the diversion area.

[0075] By spatially superimposing and analyzing Google satellite images, water remote sensing images, and land use images of the basin, a schematic diagram of the river distribution in the area where the diversion is located after the water diversion gate can be drawn. Based on the schematic diagram of the river distribution, the runoff direction of the diversion after the water diversion gate can be determined. Usually, after the diversion is diverted by the water diversion gate, its runoff can be further diverted to form a mainstream and multiple branches. The diversion river evolution model can determine the diversion ratio of the branch according to the river shape of the mainstream and branch in the diversion area.

[0076] For example, Figure 5 As shown in the figure, the main stream is divided into two branches, namely branch 1 and branch 2. The calculation formula of the branch ratio of branch 1 is: (16) In formula (16), is the branching ratio of branching flow 1, that is, the proportion of the mainstream entering branching flow 1, which is a dimensionless number and its value range is [0, 1]; are the water flow control cross-sectional areas of branch flow 1 and branch flow 2 respectively; They are the deflection angles of the river channel of branch flow 1 and the river channel of branch flow 2 relative to the mainstream river channel.

[0077] Furthermore, it is easy to obtain the branching ratio of branching flow 2 as: (17) In formula (17), is the branching ratio of branching 2, is the branching ratio of branching 1.

[0078] Taking the Xiahelike Diversion Sluice in the Tarim River Basin as an example, the distribution diagram of the mainstream control section and the branch flow control section in the diversion area behind the sluice is as follows: Figure 6 The schematic diagrams of each control section are shown in Figures 7 to 9 As shown in Table 1, the diversion ratio of the branch flow after the gate is shown in Table 1.

[0079] Table 1 Branch flow ratio behind Xiahelike diversion sluice

[0080] 2. The Muskingum-Congee method is used to calculate the flow at multiple control sections on the diversion.

[0081] The Muskingum-Konji method takes into account the evolution of diversion rivers under strong leakage conditions, and takes the leakage water as the interval outflow. The calculation formula for the flow at each control section of the diversion can be obtained as follows: (18) In formula (18), The control section for diversion, For the period, They are flow calculation related parameters respectively.

[0082] Specifically, the calculation formulas for the flow calculation related parameters are: (19) (20) (twenty one) (twenty two) In formula (19) to formula (22), To calculate the step size, is the flow weighting factor, is the storage constant.

[0083] It can be seen that each river section corresponds to a set of flow calculation related parameters, forming a parameter group.

[0084] Specifically, The calculation formula is: (twenty three) In formula (23), is the reference flow; are the width and wave speed of the river control section corresponding to the reference flow, respectively; is the length of the river section; The average slope.

[0085] Specifically, The calculation formula is: (twenty four) In formula (24), is the length of the river section; is the wave velocity of the river control section corresponding to the reference flow.

[0086] Specifically, the calculation formula for the interval outflow is: (25) In formula (25), for Control section within the time period To control section J The amount of water leakage per unit length of the river section between is the length of the river section.

[0087] Taking the Xiahelike diversion gate as an example, the calculation results of its parameter group and leakage water volume are shown in Table 2.

[0088] Table 2 Parameter group and leakage volume of the diversion after the Xiahelike diversion sluice

[0089] Among them, the It represents the amount of water leakage per unit length of a river section. It corresponds to different river sections and represents the amount of water leakage per unit length of the corresponding river section.

[0090] 3. Surface water flow-water level conversion model for the river channel in the diversion area According to the river distribution diagram of the area where the diversion is located after the water diversion gate, the key control section after the water diversion gate can be determined. The river surface water flow-water level conversion model can be based on the shape of the control section and combined with the Manning formula to determine the flow-water level conversion relationship at each control section after the water diversion gate.

[0091] Specifically, taking the Tarim River Basin as an example, the number of controlled sections, river sections and parameter groups after each water diversion gate in the basin is shown in Table 3. Among them, the number of river sections and parameter groups after each water diversion gate is consistent.

[0092] Table 3 The river sections, parameter groups and number of controlled sections after each diversion sluice in the Tarim River Basin

[0093] Furthermore, based on the data of the Digital Elevation Model (DEM) of the diversion area and related design parameters, the shape of the key control section behind the water diversion gate can be obtained.

[0094] Taking the Xiahelike Diversion Gate as an example, the river surface water flow-water level conversion model can determine the flow-water level conversion relationship at each key control section based on the shape of each key control section and combined with the Manning formula.

[0095] Specifically, there are eight key control sections behind the Xiahelike diversion gate, which are numbered Ⅰ, Ⅱ, Ⅲ, Ⅳ, Ⅴ, Ⅵ, Ⅶ, and Ⅷ. The flow-water level conversion relationship at each key control section is: (26) (27) (28) (29) (30) (31) (32) (33) In formula (26) to formula (33), are the water levels at each key control section (m), are the flow rates at each key control section (m³ / s).

[0096] (iv) Numerical simulation model of groundwater in the diversion area.

[0097] 1. Determine the corresponding relationship between groundwater evaporation and groundwater depth under different surface cover types in the diversion area.

[0098] Specifically, this embodiment determines multiple phreatic water evaporation data and multiple groundwater depth data for each of multiple sub-areas with different surface cover types in the area where the diversion of the main stream is located based on the hydrological and geographical information of the area where the diversion of the main stream is located; the multiple phreatic water evaporation data and multiple groundwater depth data for each sub-area correspond one to one.

[0099] For example, the phreatic evaporation data and groundwater depth data corresponding to different surface cover types in the Tarim River Basin are shown in Table 4.

[0100] Table 4 Phreatic evaporation and groundwater depth corresponding to different land cover types in the Tarim River Basin

[0101] Furthermore, this embodiment determines the corresponding relationship between the phreatic evaporation of each sub-region and its groundwater depth according to a plurality of phreatic evaporation data and a plurality of groundwater depth data of each sub-region and adopts a polynomial fitting method.

[0102] For example, by fitting the data in Table 4 with a polynomial, we can obtain the corresponding relationship between groundwater evaporation and groundwater depth under different surface cover types: (34) (35) (36) In formula (34) to formula (36), are the phreatic water evaporation of tree forest, shrub forest and grassland respectively (mm); is the groundwater depth (m).

[0103] Specifically, the schematic diagram of the polynomial fitting effect corresponding to the three land cover types is shown in the figure below: Fig.11 shown.

[0104] In addition to considering vegetation cover, it is usually necessary to consider groundwater evaporation under bare soil cover.

[0105] Specifically, the calculation formula for the phreatic water evaporation of the sub-area covered by bare soil is: (37) In formula (37), is the phreatic evaporation under bare soil cover (mm); It is the ultimate strength of diving evaporation; is the observed evaporation volume of the evaporating dish (mm); Calculate relevant parameters for diving evaporation.

[0106] Specifically, The calculation formula is: (38) In formula (38), is the ultimate strength of diving evaporation, Buried deep for diving.

[0107] Specifically, the parameters Submerged Impact, different diving depths Lower Parameters The values ​​of are shown in Table 5.

[0108] Table 5 Parameters at different burial depths Table of values

[0109] To facilitate calculation, this embodiment uses a linear fitting method to calculate different burial depths. Lower Parameters The value of is fitted.

[0110] Specifically, the submerged depth needs to be determined before fitting. The calculation formula for standardization is: (39) In formula (39), is the standardised groundwater depth (m), To dive deep, is the average value of the burial depth (m), is the standard deviation of the burial depth series. Combining Table 4 and Table 5, The value of is 1.25, The value of is 0.7906.

[0111] Specifically, the standardized diving depth With parameters Perform fitting and obtain parameters The calculation formula is: (40) In formula (40), is the standardized groundwater burial depth (m).

[0112] Diving depth With parameters The schematic diagram of the fitting effect is as follows Fig.10 As shown, it can be seen that the diving depth With parameters The fitting effect is better.

[0113] 2. Use Flopy to establish a three-dimensional simulation space for the groundwater numerical simulation model.

[0114] 2.1 Spatial discretization. A basin DEM map with a resolution of 30m×30m was selected as input to draw the geometric shape of the aquifer. Then, the gridding method was used to perform spatial discretization based on the size of the simulated area, the distribution of the water system, and the hydrogeological characteristics. The simulated area of ​​the Tarim River basin is 19031km 2 , the length of the simulation area is x Direction (east-west direction) is 525.67km. y The direction (north-south direction) is 116.60 km. After spatial discretization, the model has a total of 61,542 grids with a grid size of 1 km × 1 km, of which the number of active units is 20,081.

[0115] 2.2. Set the boundaries. Set the upper boundary of the vertical boundary to the water table, and the lower boundary of the vertical boundary to the top of the sandy loam. The lateral boundaries can be divided into two types: basin boundaries, water body boundaries, and boundaries with known flow rates.

[0116] 2.3. Set the initial flow field and parameters in the simulation area. According to the hydrogeological and soil conditions in the basin, the hydrogeological parameters of the three-dimensional simulation space are divided into zones, and the initial flow field of the three-dimensional simulation space is set. Fig.12 As shown in Figure 2, the hydrogeological parameters are divided into Fig.13 shown.

[0117] 3. Establish the groundwater movement control equation based on Darcy's law.

[0118] The governing equation for groundwater movement is shown below: (41) In formula (41), is the source and sink term (m / d), is the diving water level elevation (m), ) is the thickness of the aquifer, is the initial water level (m), A boundary moment t Water level (m), is the water supply degree of the unconfined aquifer, is the permeability coefficient of the unconfined aquifer (m / s), The second type of boundary Single-width flow at time (m 2 / d), is the elevation of the aquifer bottom, is the water level boundary, is the flow boundary, is the external normal vector of the flow boundary.

[0119] This embodiment can simulate the phreatic water level elevation of the diversion area according to the river water levels at multiple control sections on the diversion and use a groundwater numerical simulation model, and then the groundwater level of the diversion area can be obtained according to the phreatic water level elevation.

[0120] Specifically, after the above four models are constructed in this embodiment, the four models are coupled to obtain a surface water-groundwater system coupling simulation system.

[0121] Specifically, this embodiment adopts a parameter estimation (PEST) automatic calibration method, takes the groundwater level as a calibration index, and calibrates the surface water-groundwater system coupling simulation system by correcting relevant parameters of the coupling simulation system.

[0122] For example, the coupled simulation system can be calibrated by correcting relevant parameters such as permeability coefficient, riverbed hydraulic conductivity, water supply degree, etc. in the coupled simulation system to minimize the difference between the observed value and the simulated value of the groundwater level at the monitoring well.

[0123] Specifically, this embodiment selects three error indicators, namely Nash efficiency index (NSE), determination coefficient (R 2 ) and Mean Absolute Error (MAE), and use these three error indicators to quantitatively describe the simulation error of the model.

[0124] The calculation method of each error index is: (42) (43) (44) In formula (42) to formula (44), is the total number of monitoring sequences of monitoring wells, For the The measured groundwater level in the monitoring wells, For the Simulated groundwater levels in monitoring wells, For the The average groundwater level measured in the monitoring wells is For the Simulated groundwater depth in the monitoring wells.

[0125] Specifically, the simulation errors obtained according to the above three error indicators are dimensionless values, NSE and The value range of is [0, 1]. The closer the simulation error is to 1, the better the correction effect is. The smaller the value of MAE is, the better the correction effect is.

[0126] For example, taking the Tarim River Basin as an example, there are four monitoring stations in the basin, and each monitoring station corresponds to two monitoring wells. In this embodiment, the full calendar year from January to December 2013 is selected as the calibration period of the surface water-groundwater system coupling simulation system to determine the parameters and water balance of the coupling simulation system; the full calendar year from January to December 2014 is set as the verification period of the surface water-groundwater system coupling simulation system.

[0127] During the calibration and validation period, the simulation errors of each monitoring well calculated according to the above three error indicators are shown in Table 6.

[0128] Table 6 Simulation errors of each monitoring well during the calibration and validation period

[0129] According to Table 6, during the calibration period, the minimum NSE value of the eight monitoring wells is 0.717. The minimum value is 0.871, and the maximum value of MAE is 0.639. During the validation period, the minimum NSE value of the eight monitoring wells is 0.8. The minimum value is 0.878, and the maximum value of MAE is 0.168. This shows that the surface water-groundwater system coupling simulation system constructed by the present invention has good performance in both the calibration period and the validation period, and can accurately describe the surface water processes in the basin, including the evolution of the mainstream river and the evolution of the diversion river, as well as the discharge processes such as the recharge and evaporation of groundwater, so that it can accurately simulate and simulate the dynamic changes of the groundwater level under different water inflow conditions and different water diversion conditions in the basin.

[0130] Specifically, the comparison diagram between the measured and simulated groundwater levels of each monitoring well during the model calibration period is shown in the figure below. Figures 14 to 16 The comparison diagram of the measured and simulated values ​​of the groundwater level in each monitoring well during the model validation period is shown in Figures 17 to 19 shown.

[0131] Furthermore, after calibrating the surface water-groundwater system coupling simulation system, this embodiment uses the surface water-groundwater system coupling simulation system to conduct the following research and analysis on the surface water-groundwater conversion law in the basin.

[0132] 1. The historical water inflow and water diversion monitoring data of each diversion gate are used as model input. The time series is 2014-2022. The groundwater level and groundwater recharge and discharge relationship under different water inflow conditions of the main stream and different water diversion conditions of the diversion gate are simulated and analyzed.

[0133] Specifically, the groundwater storage variable is calculated based on the model simulation of groundwater recharge and discharge, and the calculation formula is: (45) In formula (45), is the groundwater storage variable (×10 8 m 3 ); The amount of water that is replenished to groundwater by surface water leakage (×10 8 m 3 ); is the lateral replenishment water volume (×10 8 m 3 ); is the amount of water discharged laterally (×10 8 m 3 ); is the amount of water discharged laterally (×10 8 m 3 ), since groundwater exploitation is not considered in this embodiment, this item is zero; is the diving evaporation (×10 8 m 3 ).

[0134] 2. Analyze the groundwater recharge, discharge characteristics and groundwater storage in the upper and middle reaches of the Tarim River from 2014 to 2022. Fig. 20 shown.

[0135] 3. The upper and middle reaches of the Tarim River are divided into four sub-areas according to the location of the hydrological stations. The groundwater recharge and discharge characteristics and groundwater storage variables of each sub-area are analyzed, and the lateral recharge and discharge characteristics of groundwater between the sub-areas are analyzed. The schematic diagram of the groundwater recharge and discharge characteristics and groundwater storage variable changes in each river section in the Tarim River basin is obtained as shown in the figure. Figure 21 to Figure 24 shown.

[0136] The above are only a few embodiments of the present application and do not constitute any form of limitation to the present application. Although the present application is disclosed as above with preferred embodiments, it is not intended to limit the present application. Any technician familiar with the profession, without departing from the scope of the technical solution of the present application, using the technical content disclosed above to make slight changes or modifications are equivalent to equivalent implementation cases and fall within the scope of the technical solution.

Claims

1. A surface water-groundwater coupling simulation method considering the characteristics of a high-loss river, characterized in that: The method comprises: S1. According to the hydrological and geographical information of the area where the diversion of the main stream is located, determining the corresponding relationship between the phreatic evaporation amount and the groundwater depth of each sub-area in a plurality of sub-areas with different surface cover types in the area where the diversion is located; S2. constructing a groundwater numerical simulation model of the area where the diversion is located according to the corresponding relationship and the hydrological and geographical information of the area where the diversion is located; S3. According to the divertable water volume of the main stream and the diversion conditions of the diversion, the groundwater level in the diversion area is simulated using the groundwater numerical simulation model.

2. The method according to claim 1, characterized in that Before S3, the method further includes: The amount of water that can be diverted from the main stream is determined based on the hydrological and geographical information of the main stream.

3. The method according to claim 2, characterized in that The method of determining the divertable water volume of the main stream according to the hydrological and geographical information of the main stream is specifically as follows: Determine the water loss of the river channel of the main stream according to the hydrological and geographical information of the main stream; The amount of water that can be diverted from the main stream is determined based on the amount of water loss in the river channel and the hydrological and geographical information of the main stream.

4. The method according to claim 3, characterized in that The water loss in the river channel includes water leakage, water evaporation and water overflow in the river channel.

5. The method according to claim 1, characterized in that The S3 is specifically: S31, determining the river water levels at a plurality of control sections on the diversion according to the diversion amount of the main stream and the diversion conditions of the diversion; S32. According to the river water levels at the multiple control sections on the diversion, the groundwater level in the area where the diversion is located is simulated using the groundwater numerical simulation model.

6. The method according to claim 5, characterized in that The S31 is specifically: S311, determining the flow rates at multiple control sections on the diversion according to the divertable water volume of the main stream and the diversion conditions of the diversion; S312, determining the river channel water levels at the multiple control sections on the diversion according to the flow rates at the multiple control sections on the diversion and the shapes of the multiple control sections.

7. The method according to claim 6, characterized in that The branch flow includes a main flow and multiple branch flows; the S312 is specifically as follows: Determine the flow rates at multiple control sections on the main stream according to the divertable water volume of the main stream and the diversion conditions of the diversion; The flow splitting ratio of each branch flow from the main flow is determined, and the flow rates at the multiple control sections on each branch flow are determined according to the flow splitting ratio and the flow rates at the multiple control sections on the main flow.

8. The method according to claim 1, characterized in that After S3, the method further includes: S4. Calibrate the groundwater numerical simulation model using the historical groundwater level in the area where the diversion is located.

9. The method according to claim 8, characterized in that The S4 is specifically: S41, determining a simulation error of the groundwater numerical simulation model according to the simulated groundwater level and the corresponding historical groundwater level; S42: When the simulation error does not satisfy a preset error range, adjusting the calculation parameters of the groundwater level until the simulation error satisfies the preset error range.

10. The method according to claim 1, characterized in that The land cover type is tree forest, shrub forest, grassland or bare land.

Images