Method for simulating environmental exposure concentration of chemical substances in lake meeting watershed scale

Abstract

The invention provides a lake chemical substance environmental exposure concentration simulation method meeting the watershed scale, which comprises the following steps: carrying out regional division on a lake according to a division rule to obtain a plurality of subareas, and recording geometric parameters of each subarea, the division rules comprise a first rule for dividing a unique water outlet of the lake and a water inlet of the lake into different zones and a second rule for dividing the zones according with chemical substance concentration management requirements; according to hydrological data of the current day from the watershed hydrological model, the estimated value of the water level of each subarea and the geometric parameters, calculating a water flow calculation value between the subareas and a water outlet flow calculation value of the water outlet subarea on the current day; based on optimization indexes, the partitions are optimized and adjusted by adjusting the water levels or volumes of the partitions, and the optimization indexes comprise a first index related to the calculated value of the water outlet flow of the water outlet partition on the current day and a second index related to the water flow of the partitions on the current day; and according to chemical substance data from a watershed hydrological model, calculating the intraday chemical substance concentration of the adjusted subarea.

Description

Technical Field

[0001] The present disclosure relates to the fields of environmental risk analysis, environmental risk assessment and environmental risk management of chemical substances, and more specifically, to a method for simulating the environmental exposure concentration of chemical substances in lakes that meets the watershed scale. Background Art

[0002] High-precision simulation of toxic and hazardous chemicals is an essential approach for exposure assessment and risk assessment of new chemical pollutants in new pollutant control. In existing basin-scale chemical simulation models, most lakes in the basin are simulated as a whole, that is, the lakes are not divided into regions, and the behavior of chemicals in the lakes is assumed to be uniformly distributed in horizontal space. This is common in basin multi-media chemical simulation. User needs are mainly focused on the concentration distribution of chemicals in the soil of different land space units, different river sections, and air of different space units. The concentration differences of chemicals at different locations in the lake are not paid much attention to. This is often sufficient to meet user needs under conditions such as small lake area, fast inflow and outflow flow, small lake water demand, and relatively simple lake ecosystem service functions. For example, the SWAT-KM model follows the design of the SWAT (Soil and Water Assessment Tool) basin hydrological model, takes a holistic view of lakes in the basin (including natural lakes, artificial lakes, and reservoirs), and realizes the daily coordinated simulation of chemicals in the water phase and sediment phase of the lake and basin units (sub-basins, river sections), etc.

[0003] There are many lake models for simulating lakes, such as the EFDC (Environmental FluidDynamics Code) model, which can spatially grid lakes and construct a hydrodynamic model of lake water flow by combining mass equations and motion equations. It can achieve high-precision simulation of time (the step size can be very low, such as seconds and minutes) and spatial points (grid accuracy can be as fine as meters). Obviously, if the lake model that simulates the lake alone cannot be linked with the basin, it is difficult to express the complex migration and transformation process of chemical substances in multiple spatial units and multiple media within the basin. However, by combining the lake model with the basin hydrological model, the daily coordinated simulation of chemical substances in the water phase (some include sediment phase) in the lake and the basin units (sub-basins, river sections) can be achieved. At this time, the technical bottleneck is how to achieve the spatiotemporal connection between the lake model and the basin hydrological model.

[0004] The hydrodynamic lake model has higher temporal accuracy. It is ideal to select a watershed hydrological model with the same high temporal resolution to interface and couple with it. However, the practical difficulty is that lakes in relatively large watersheds often have the need for high-precision chemical exposure simulation. Using a very high temporal resolution in a large watershed is likely to be inefficient or even a waste of resources - because a high temporal resolution requires a relatively high spatial unit resolution to have its simulation value. In large watershed simulations with both high spatial and temporal resolutions, the computational load of the simulation operation will increase sharply, forming a computational resource barrier and making it impossible to meet the actual needs of a large number of chemical exposure assessments at the watershed scale. Moreover, in terms of watershed hydrological models, high spatial and temporal resolutions often refer to grid models, such as the typical TOPKAPI-ETH model. This model has a large amount of calculation and is suitable for short-term precipitation process simulation, but not for long-term (e.g., more than ten years, several decades) large watershed simulations. And in many cases, its advantage of high temporal resolution cannot be fully utilized. With the same level of basic data support, it does not necessarily achieve significantly better results than a daily-precision hydrological model (such as the SWAT model).

[0005] The high-precision hydrodynamic lake model is coupled with a daily-precision hydrological model. For example, the EFDC model is coupled with the SWAT model. The SWAT provides daily-precision simulation values (such as the inflow into the lake, sediment inflow into the lake, etc.), and then it is assumed that the daily value is evenly distributed within a day or a time distribution within a day is assumed as the boundary condition for the EFDC. However, it is precisely this practice of distributing the daily value input to a higher temporal resolution that brings great uncertainties and, to a considerable extent, masks the simulation ability and utility of the high-precision hydrodynamic lake model. That is, the seemingly high precision actually has a limited reliability and is equivalent to a lower-precision simulation result. Moreover, in this coupling mode, the lake model has a slow calculation speed and often becomes the bottleneck step in the simulation process, which is not conducive to carrying out long-term large watershed simulations of chemical substances.

[0006] In summary, the lake models in the related technologies can only meet the basic simulation needs, and there are technical bottlenecks in how the lake model is spatio-temporally interfaced with the watershed hydrological model. If the temporal resolution of the watershed hydrological model is improved for interfacing, there are problems such as waste of computing resources and it is not suitable for long-term simulation. Summary of the Invention

[0007] The purpose of the present disclosure is to provide a method for simulating the environmental exposure concentration of chemical substances in lakes at the watershed scale to solve at least one of the problems existing in the prior art.

[0008] To achieve the above object, the present disclosure adopts the following technical solutions: The present disclosure provides a method for simulating the environmental exposure concentration of chemical substances in lakes at the watershed scale, including: The lake is regionally divided according to the division rules to obtain multiple sub-regions, and the geometric parameters of each sub-region are recorded. The division rules include a first rule and a second rule. The first rule is that the only water outlet of the lake and the water inlet of the lake are divided into different sub-regions. The second rule is that the regional division meets the management requirements of chemical substance concentration; According to the daily hydrological data from the watershed hydrological model, the estimated values of the water levels of each sub-region, and the geometric parameters of each sub-region, calculate the calculated value of the daily water flow between each sub-region and the calculated value of the daily water discharge of the water discharge sub-region containing the water outlet; Based on the optimization index, the sub-regions are optimized and adjusted by adjusting the water level of the sub-region or adjusting the volume of the sub-region to obtain the adjusted sub-regions. The optimization index includes a first index related to the calculated value of the daily water discharge of the water discharge sub-region and a second index related to the daily water flow between each sub-region; According to the chemical substance data from the watershed hydrological model, calculate the daily chemical substance concentration of each adjusted sub-region.

[0009] Optionally, The first index is the ratio of a first value to a second value. The first value is the difference between the calculated value of the daily water discharge of the water discharge sub-region and the actual value of the daily water discharge of the lake. The second value is the product of a first preset ratio value and the actual value of the daily lake flow profit and loss; and / or, The second index is the sum of the first weighted ratio values of each sub-region. The first weighted ratio value of a sub-region is the product of the first ratio value of the sub-region and the weight. The first ratio value of the sub-region is the absolute value of the ratio of the daily flow profit and loss value of the sub-region to the daily water inflow value of the sub-region. The weight is the ratio of the daily water inflow value of the sub-region to the sum of the daily water inflow values of each sub-region.

[0010] Optionally, the optimizing and adjusting the sub-regions by adjusting the water level of the sub-region or adjusting the volume of the sub-region based on the optimization index includes: Judge whether the first index is less than a first threshold and whether the second index is less than a second threshold: If the first index is greater than or equal to the first threshold and the second index is less than the second threshold, then optimize and adjust the sub-regions by adjusting the daily water level of the sub-region. After the optimization and adjustment, recalculate the calculated value of the daily water flow between each sub-region and the calculated value of the daily water discharge of the water discharge sub-region and re-judge; If the first index is less than the first threshold and the second index is greater than or equal to the second threshold, then optimize and adjust the sub-regions by adjusting the volume of the sub-region. After the optimization and adjustment, recalculate the calculated value of the daily water flow between each sub-region and the calculated value of the daily water discharge of the water discharge sub-region and re-judge; If the first index is greater than or equal to the first threshold and the second index is greater than or equal to the second threshold, then when the first index is greater than or equal to the second index, optimize and adjust the partition by adjusting the partition water level, and when the first index is less than the second index, optimize and adjust the partition by adjusting the partition volume. After the optimization and adjustment, recalculate the calculated value of the daily water flow between each partition and the calculated value of the daily water output of the water output partition, and rejudge; If the first index is less than the first threshold and the second index is less than the second threshold, then end the optimization and adjustment.

[0011] Optionally, the adjusting of the partition water level includes: Adjust the water level of the partition with the largest absolute value of the daily flow profit and loss value in each partition to reduce the first index; or, Batch-adjust the water levels of each partition to reduce the first index.

[0012] Optionally, the adjusting of the partition volume includes: Take the partition with the largest first weighted ratio value as the first partition, select a partition that has a water flow exchange with the first partition as the second partition, and divide a first set ratio of the volume of the first partition into the second partition or divide a first set ratio of the volume of the second partition into the first partition.

[0013] Optionally, the adjusting of the partition volume further includes: According to the volumes of the first partition and the second partition after dividing a first set ratio of the volume of the first partition into the second partition or dividing a first set ratio of the volume of the second partition into the first partition, update the water levels and geometric parameters including the adjacent side lengths and the planar area of the first partition and the second partition respectively.

[0014] Optionally, the selecting a partition that has a water flow exchange with the first partition as the second partition includes: Take the partition with the largest daily water flow exchange value with the first partition as the second partition; or, Take the partition that has the greatest impact on the daily flow profit and loss value of the first partition as the second partition.

[0015] Optionally, the adjusting of the partition volume includes: Take the partition with the largest daily water flow exchange value with the first partition as the second partition and perform iterative volume division until the second index is less than the third threshold or the number of iterations reaches the fourth threshold; Execute the volume division cycle for the first set number of times. One volume division cycle includes: taking the partition that has the greatest impact on the daily flow profit and loss value of the first partition as the second partition and performing iterative volume division for the second set number of times, and then taking the partition with the largest daily water flow exchange value with the first partition as the second partition and performing iterative volume division for the third set number of times.

[0016] Optionally, the volume value of a single volume division is set according to the maximum first weighted ratio value of each partition, and is reduced at a second set ratio when the first constraint condition is not satisfied. The first constraint condition includes: The volume value of a single volume division is less than or equal to the product of the second preset ratio value and the volume of the first partition and the product of the second preset ratio value and the volume of the second partition; The volume of the first partition after volume division is greater than or equal to the product of the third preset ratio value and the volume of the first partition; The water level change value of the first partition after volume division and the water level change value of the second partition after volume division are less than or equal to the fifth threshold value.

[0017] Optionally, calculating the chemical substance concentration of each adjusted partition on the same day includes: Judging whether the actual value of the lake effluent flow on the same day is zero: If not, calculate the chemical substance concentration of each adjusted partition according to the one-dimensional advection-diffusion equation; If so, calculate the chemical substance concentration of each adjusted partition according to the one-dimensional advection-diffusion equation combined with the return flow attenuation based on the set number of return times.

[0018] The beneficial effects of the present disclosure are as follows: The method for simulating the environmental exposure concentration of chemical substances in a lake that meets the basin scale provided by the present disclosure can scientifically divide the lake into regions. The obtained partitions are reasonable in terms of chemical substance migration and diffusion, meet the requirements of chemical substance concentration management (including monitoring requirements, etc.), and can also ensure the scientificity and rationality of the partition results by optimizing and adjusting the partitions in the case of missing lake information such as water level data, so as to accurately simulate the environmental exposure concentration of chemical substances in the lake. Based on the design of the present disclosure for the requirements of chemical substance concentration management, the lake model can be directly connected to the basin hydrological model, and the chemical substance concentration in the lake can be simulated without improving the time resolution of the basin hydrological model. The calculation speed can be improved while ensuring the simulation quality of the chemical substance concentration in the partition, and it is beneficial for long-term simulation. BRIEF DESCRIPTION OF THE DRAWINGS

[0019] The following further describes in detail the specific implementation manners of the present disclosure with reference to the accompanying drawings.

[0020] Figure 1 The flowchart showing the method for simulating the environmental exposure concentration of chemical substances in a lake that meets the basin scale provided by the embodiments of the present disclosure.

[0021] Figure 2 The overall flow schematic diagram showing the method for simulating the environmental exposure concentration of chemical substances in a lake that meets the basin scale provided by the embodiments of the present disclosure.

[0022] Figure 3 Schematic diagram showing the zonal planar distribution of a lake.

[0023] Figure 4 Show Figure 3 Schematic diagram showing the quasi-rectangular flow setting between the shown lake zones and the equal-area circles of the zones.

[0024] Figure 5 Schematic diagram showing the process of optimizing and adjusting the zones.

[0025] Figure 6 Schematic diagram showing the zonal planar distribution of the lake in the example.

[0026] Figure 7 Schematic diagram showing the dynamic change of the daily simulated concentration of a certain chemical substance in each zone of the lake in the example.

[0027] Figure 8 Schematic diagram showing the quarterly average concentration of a certain chemical substance in each zone of the lake, the inflow section, and the outflow section of the lake in the example.

[0028] Figure 9 Schematic diagram showing the structure of the computer system of the device for implementing the method for simulating the environmental exposure concentration of chemical substances in a lake that meets the basin scale provided by the embodiments of the present disclosure. Detailed implementation manners

[0029] To illustrate the present disclosure more clearly, the present disclosure will be further described below in conjunction with embodiments and the accompanying drawings. Similar components in the drawings are denoted by the same reference numerals. Those skilled in the art should understand that the content specifically described below is illustrative rather than restrictive, and should not be used to limit the protection scope of the present disclosure.

[0030] In the description of the present disclosure, it should be noted that the orientation or positional relationship indicated by terms such as "upper", "lower", etc. is based on the orientation or positional relationship shown in the drawings, and is only for the convenience of describing the present disclosure and simplifying the description, rather than indicating or implying that the device or element referred to must have a specific orientation, be constructed and operated in a specific orientation, and therefore should not be construed as a limitation of the present disclosure. Unless otherwise clearly specified and defined, the terms "installed", "connected", and "connected" should be understood in a broad sense. For example, it can be a fixed connection, a detachable connection, or an integral connection; it can be a mechanical connection, an electrical connection; it can be directly connected, or indirectly connected through an intermediate medium, and it can be the communication inside two elements. For those of ordinary skill in the art, the specific meanings of the above terms in the present disclosure can be understood according to specific circumstances.

[0031] It should also be noted that in the description of the present disclosure, relational terms such as first and second are only used to distinguish one entity or operation from another entity or operation, and do not necessarily require or imply any actual relationship or order between these entities or operations. Moreover, the terms "comprising", "including" or any other variant thereof are intended to cover non-exclusive inclusion, so that a process, method, article or device comprising a series of elements not only includes those elements but also includes other elements not expressly listed, or elements inherent to such process, method, article or device. Without further limitation, an element defined by the statement "comprising one..." does not exclude the presence of additional identical elements in the process, method, article or device comprising the said element.

[0032] An embodiment of the present disclosure provides a method for simulating the environmental exposure concentration of chemical substances in a lake at the basin scale, as Figure 1 shown, including steps S1, S2, S3 and S4. Among them, steps S1, S2 and S3 are lake area division steps, constituting a lake area division method, and step S4 is the simulation step of the daily chemical substance concentration in each sub-area of the lake after the lake area division. The introductions of steps S1, S2, S3 and S4 are as follows: S1. Divide the lake into multiple sub-areas according to the division rules, record the geometric parameters of each sub-area. The division rules include a first rule and a second rule. The first rule is to divide the only water outlet of the lake and the water inlet of the lake into different sub-areas, and the second rule is that the area division conforms to the management requirements of chemical substance concentration.

[0033] In a possible implementation manner, on the premise of meeting the first rule and the second rule, the division rules further include: A third rule: the planar area of the sub-area containing the water inlet of the lake corresponds to the water inflow value of the water inlet, that is, the planar area of the water inlet sub-area with a larger water inflow value of the water inlet should be larger; A fourth rule: data such as the water depth (lake water depth) and lake bottom roughness inside the sub-area are relatively uniform; A fifth rule: the planar shape of the sub-area is as close as possible to a circle or a square.

[0034] It can be understood that the first rule and the second rule are high-priority rules, and the third rule, the fourth rule and the fifth rule are low-priority rules. For the area division of the lake, the high-priority rules should be satisfied first, and on this basis, the low-priority rules are considered.

[0035] In a specific example, the inventors found that, from the analysis of the actual needs of chemical substance environmental exposure assessment and risk assessment, there is a need for the differential levels of chemical substance concentrations in several regions of a lake. In this embodiment, the management requirements for chemical substance concentrations in lake zoning, for example, include: monitoring points are set at different locations in a lake, and it is necessary to estimate and compare the daily changes in the chemical substance concentrations in the waters where the monitoring points are located. There are water intake demands (such as agricultural, ecological, and industrial water use) in some areas of the lake, and it is necessary to focus on the fluctuations in chemical substance concentrations in these areas. A certain area of the lake is a gathering place or habitat for some important species (such as aquatic organisms and waterfowl), and it is necessary to pay attention to the chemical substance concentration levels in this area, and so on.

[0036] In this embodiment, the management requirements for chemical substance concentrations in lake zoning have the following characteristics: The lake only needs to be divided into not too many partitions on the plane, and there is no rigid requirement for the planar areas of each partition to be close or even the same. The regional division is not a universal high spatial resolution.

[0037] For the chemical substance concentrations (or chemical substance environmental exposure concentrations) in lake zoning, it is sufficient to have daily value concentrations (or daily simulations), and there is no need to pay attention to hourly or finer time resolutions.

[0038] The chemical substance concentrations in lake zoning should co-vary with the chemical substance environmental behavior in the basin to achieve long-term simulations.

[0039] Therefore, the above division rules that take into account both scientificity and detection requirements are set in step S1, making the initial regional division of the lake area in step S1 reasonable and easy to achieve seamless connection with the multi-media daily value simulation results of chemical substances in the basin. Among them, the first rule can ensure seamless connection between partition simulation and basin simulation, and can avoid the chemical substances carried in individual inflow river sections flowing directly out of the outlet without passing through other partitions (that is, avoid "short-circuit" outflow), ensuring scientificity. The second rule limits that the regional division can serve actual management needs. The second rule is reflected, for example, in setting partitions at monitoring points, dividing the gathering areas of important species into partitions, etc. The third rule is more conducive to docking with basin simulation and is also conducive to the quick realization of subsequent optimization and adjustment of partitions. The fourth rule can ensure that the internal differences in partitions have no significant impact on the migration and diffusion of chemical substances, thus ensuring the rationality of the assumed conditions in subsequent calculations; the fifth rule is to ensure the rationality of the partition flow and the migration and diffusion of chemical substances based on the partition flow in this embodiment.

[0040] The first to fifth rules constitute rule set I. For example Figure 2 as shown, step S1 divides the lake according to rule set I.

[0041] For exampleFigure 3 The shown lake has a total of 4 water inlets (that is, the lake has 4 water inlet river sections or inflowing rivers, that is, the upstream river sections of the lake, including the main stream) and 1 water outlet (that is, the lake has 1 water outlet river section). Based on Rule Set I, the lake is divided into 7 sub-regions. The 4 water inlets are respectively divided into sub-regions ①, ②, ④, and ⑥, the water outlet is divided into sub-region ⑤, and there is no external water flow connection in sub-regions ③ and ⑦. It should be emphasized that, as Figure 3 shown, although, as Figure 3 shown in the lower part, the water inlet and the water outlet (the water outlet is the position of the weir in Figure 3 connected to the downstream river section) are close in position, they cannot be divided into the same sub-region, but as Figure 3 shown, they are respectively divided into sub-region ⑥ and sub-region ⑤. Even if the planar area of sub-region ⑥ is relatively small, it is still necessary to do so. In addition, as Figure 3 shown, the planar shapes of the 7 sub-regions are relatively square, without strange shapes or long and narrow shapes.

[0042] S2. According to the daily hydrological data from the basin hydrological model, the estimated values of the water levels in each sub-region, and the geometric parameters of each sub-region, calculate the calculated value of the daily water flow between each sub-region and the calculated value of the daily water discharge of the water discharge sub-region containing the water outlet.

[0043] Continuing the foregoing example, this embodiment is expected to be supported by simple and easily obtainable data to meet the requirements of simulating the concentration of lake sub-regions in the basin-scale chemical substance environmental exposure simulation. The application conditions and settings of this embodiment include: (1) being used in conjunction with a basin-scale chemical substance environmental exposure simulation model or a basin hydrological model, and providing the basic simulation data for lake sub-region simulation from the basin simulation results. For example, when the lake is regarded as a whole, the inflowing river section (upstream river section) and the outflowing river section (downstream river section) are clearly identifiable, and the water flow, sediment content, and the content of various forms of chemical substances of the inflow and outflow are simulated. (2) The lake has the characteristics of a shallow lake, that is, the "turning over of the lake" phenomenon is not obvious, and the vertical water flow movement is significantly weaker than the horizontal water flow movement. (3) The lake requires a relatively small number of sub-regions. (4) The time resolution of simulating the chemical substance exposure concentration in the lake sub-regions is the same as that of the basin hydrological model, both being daily values. (5) The lake has basic geometric terrain data, including water depth distribution, water level values, water discharge forms and their parameters (such as weir height, weir width, etc.), lake surface data (such as water surface map), etc. (6) The lake has only one water outlet (or outflowing river section, water outlet river section).

[0044] The specific steps of step S2 include: After obtaining multiple partitions in step S1, the geometric parameters and water levels of each partition are obtained. The basic assumptions include: the water depth is the same inside each partition and other hydrodynamic properties inside the partition are the same, and the central position of each partition in the plane is stable (even if the partition is optimized and adjusted subsequently, the central position of the partition in the plane remains unchanged). The geometric parameters of the partition mainly include: Partition water depth, denoted as d i represents the water depth of the i-th partition, and the set of water depths of each partition is denoted as {d}.

[0045] Partition planar area, denoted as As i represents the planar area of the i-th partition, and the set of planar areas of each partition is denoted as {As}. For example Figure 4 as shown, the partition is approximated by an equal-area circle to obtain an approximate circular shape of the partition. According to the planar area As i of the i-th partition, the radius R of the equal-area circle of the i-th partition can be calculated i , for example Figure 4 the radius R in is the radius of the equal-area circle of partition ⑦.

[0046] Partition volume (also known as water storage), denoted as V i represents the volume of the i-th partition, and the set of volumes of each partition is denoted as {V}. The volume V i of the i-th partition can be approximated, for example, by the product of the water depth d i of the i-th partition and the planar area As i of the i-th partition, that is (ignoring the influence brought by the shore connection slope, lake bottom inequality, etc.).

[0047] Distance between the central points of partitions, denoted as Ds i,j represents the distance between the central points of the i-th partition and the j-th partition. The distances between the central points of each partition can be expressed in the form of a matrix, denoted as {Ds i,j}, and the diagonal values of the matrix {Ds i,j} are all 0.

[0048] Adjacent side length between partitions, denoted as Al i,j represents the adjacent side length between the i-th partition and the j-th partition. Among them, Al i,j = 0 means that the i-th partition and the j-th partition are not adjacent. The adjacent side lengths between each partition can be expressed in the form of a matrix, denoted as {Al i,j}. The diagonal value Al i,j of the matrix {Al i,i} is the planar perimeter of the i-th partition.

[0049] The method for obtaining the above geometric parameters should be based on actual lake information. After directly determining the geometric centers of each partition in GIS software or planar graphics software based on remote sensing data or surveying and mapping data, the distance matrix {Ds i,j} between the center points of each partition and the adjacent side length matrix {Al i,j} between each partition can be obtained by the methods in the relevant technology. For example, taking partition ② in Figure 3 as an example, this partition is determined by the connecting lines of nodes ABCDA. AB, BC, and CD are the adjacent sides of partition ② and partitions ①, ⑦, and ③ respectively. Let their lengths be 15 km, 20 km, and 25 km respectively. The AD section is the lake shore adjacent to partition ②, and its length is 35 km. Then the values in the second row (of course, also the second column, as the adjacent side length matrix is a symmetric matrix) of the matrix {Al i,j} are: Al 2,1 = 15, Al 2,2 = 15 + 20 + 25 + 35 = 95, Al 2,3 = 25, Al 2,4 = 0, Al 2,5 = 0, Al 2,6 = 0, Al 2,7 = 20.

[0050] It should be noted that the water depth d of each partition can be reasonably determined or inferred with reference to actual data. If a water depth distribution map or lake bottom elevation map of the lake can be obtained, the graph can be read according to the partition and the spatial average value within the planar range of the partition can be recorded as the depth of this partition. Related to the water depth d is the water level h of each partition (expressed in elevation above sea level, or taking the same elevation reference point as the weir height). Let h i represent the water level of the i-th partition, and the set of water levels of each partition is represented by {h}. The user should provide an initial estimate of the water level geometry {h} of each partition.

[0051] In addition, the geometric dimensions of the water outlet of the water outlet partition (such as partition ⑤ in Figure 3) containing the water outlet can be obtained according to actual data. For example, when the weir is used as the water outlet, the weir width and the water head elevation of the water outlet partition (the water head elevation of the water outlet partition is the water level of the water outlet partition, represented by h o ; the water outlet elevation h 出 = h o - H, where H is the weir height) should be obtained to meet the requirement of being able to calculate the water flow rate.

[0052] It should be noted that the set of water depths {d} of each partition, the set of water levels {h} of each partition, and the weir height h oThere is a certain relevance, and these data are the key data for determining the water flow between partitions. When the subsequent step S3 iteratively optimizes and adjusts the partitions, adjusting the water level is to adjust the water level set {h} of each partition, which is a correction and adjustment of the water level estimated value; adjusting the volume is a systematic adjustment of multiple data including the water depth set {d} of each partition, the planar area set {As} of each partition, the volume set {V} of each partition, and the adjacent side length matrix {Al i,j} of each partition. Under the constraint that the total volume of the lake remains unchanged, adjusting the volume can also achieve the purpose of adjusting the water level set {h} of each partition. In order to reduce the pressure of the subsequent step S3 optimization and adjustment and accelerate the algorithm convergence, the data such as the estimated water depth and water level should be as accurate as possible when referring to and estimating the actual data.

[0053] After obtaining the geometric parameters of the above partitions, calculate the calculated value of the daily water flow between each partition and the calculated value of the daily water flow of the outlet partition including the outlet based on the Chezy Equation, as follows: The Chezy Equation is a widely used formula for calculating the flow velocity of open channels (unconfined). It ignores the differences in flow velocity distribution across the streamline cross-section and only represents it with a unified flow velocity u. When calculating the daily water flow velocity (flow rate) between each partition of the lake, the Chezy Equation in the following form is used:

[0054] Among them, u is the flow velocity between the partitions of the lake, with the unit of m / s; n m is the Manning coefficient of the water flow between the partitions of the lake, and the value of a relatively smooth section of the river basin can be adopted, such as about 0.05; Rh is the hydraulic radius, with the unit of m; Sl is the hydraulic gradient, dimensionless.

[0055] The planar shapes of the partitions of the lake are basically irregular shapes. In this embodiment, the water flow between the lakes in the source-sink partitions (when calculating the water flow between partitions, the upstream partition is called the source partition and the downstream partition is called the sink partition) is set to flow in a quasi-rectangular manner, that is For the i-th partition as the source partition and the j-th partition as the sink partition, the water flow between the source-sink partitions is as follows: Taking the distance Ds between the center points of the i-th partition and the j-th partition (i.e., the source-sink partitions) i,j as the length and the adjacent side length Al between the i-th partition and the j-th partition (i.e., the source-sink partitions) i,j as the width to form a "quasi-rectangular open channel". The geometric length and width of the quasi-rectangle can be obtained from the center point distance matrix {Ds i,j} between each partition and the adjacent side length matrix {Ali,j} Obtained directly.

[0056] Under the "quasi-rectangular" flow setting, the cross-sectional area Ax of water flow between the i-th partition and the j-th partition (i.e., the source-sink partition) is calculated as follows:

[0057] During the flow, the hydraulic radius Rh along the flow direction is calculated as follows:

[0058] Among them, different from other parameters, the hydraulic radius Rh needs to consider the influence of the lake shore on the partition adjacent to the lake shore, Ib i and Ib j are binary variables of 0-1 respectively: if the i-th partition is adjacent to the lake shore, then Ib i takes the value of 1, if the i-th partition is not adjacent to the lake shore, then Ib i takes the value of 0; if the j-th partition is adjacent to the lake shore, then Ib j takes the value of 1, if the j-th partition is not adjacent to the lake shore, then Ib j takes the value of 0.

[0059] When the center position of each partition plane is stable and unchanged, that is, the center position of each partition is fixed and the distance between each partition is fixed, and does not change with the optimization and adjustment of the partition, when calculating the water flow rate based on the Chezy formula, the hydraulic gradient between the i-th partition and the j-th partition (i.e., the source-sink partition) is determined by the water level difference of the source-sink partition:

[0060] Under the above quasi-rectangular flow setting, the daily water flow rate u i,j and the daily water flow Q i,j between the i-th partition and the j-th partition (i.e., the source-sink partition) can be calculated:

[0061] The daily water flow between each partition can be expressed in the form of a matrix, denoted as {Q i,j}, and the matrix {Q i,j} is an anti-symmetric matrix with a diagonal of 0, that is, for any i, j, Q i,j = Q j,i .

[0062] For example, according to the simplified formula for the discharge of a rectangular thin-walled weir , where b is the weir width, h 出 = h o - H, h ois the water outlet elevation of the water outlet zone, H is the weir height, m is the weir flow coefficient, and g is the acceleration due to gravity (unit: m / s 2 ), and the calculated value Q of the daily water outlet flow of the water outlet zone can be obtained o .

[0063] In addition, the hydraulic radius Rh i,j between the ith zone and the jth zone, the hydraulic gradient Sl i,j and the water flow rate u i,j can be combined to calculate the turbulent diffusion coefficient generated by the water flow between the ith zone and the jth zone for use in the calculation of the chemical substance concentration in each subsequent zone. The turbulent diffusion coefficient generated by the water flow between the ith zone and the jth zone is as follows:

[0064] where the turbulent diffusion coefficient D i,j has the unit of m 2 / s, k is the von Karman constant (for example, with a value of 0.40), dimensionless. The turbulent diffusion coefficients between each zone can be expressed in the form of a matrix, denoted as {D i,j}.

[0065] S3. Based on the optimization index, the zones are optimized and adjusted by adjusting the zone water level or the zone volume to obtain the adjusted zones. The optimization index includes a first index related to the calculated value of the daily water outlet flow of the water outlet zone and a second index related to the daily water flow between each zone.

[0066] In one possible implementation, in step S3, the first index is the ratio of a first value to a second value. The first value is the difference between the calculated value of the daily water outlet flow of the water outlet zone and the actual value of the daily lake water outlet flow, and the second value is the product of a first preset ratio value and the actual value of the daily lake flow profit and loss; and / or, the second index is the sum of the first weighted ratio values of each zone. The first weighted ratio value of a zone is the product of the first ratio value of the zone and the weight. The first ratio value of the zone is the absolute value of the ratio of the daily flow profit and loss value of the zone to the daily water inflow value of the zone, and the weight is the ratio of the daily water inflow value of the zone to the sum of the daily water inflow values of each zone.

[0067] The first index is an index that the daily flow profit and loss of the overall lake simulation should match the actual situation. Not only should the daily flow profit and loss of the overall lake simulation match the actual situation, but in a reasonable lake area division result, it should also ensure that the daily flow profit and loss of each zone is related to its geometric size (such as the water depth d i of the zone, the planar area As i, the sum of the adjacent side lengths Al of the partition and other partitions i etc.) and the water level h i match. This is essentially a process of solving the geometric shapes of each partition - if the shape information of a lake and its partitions is rich enough, in theory, all parameters of all partitions and the daily flow between partitions can be solved. However, this embodiment precisely aims at the situation where the lake information is insufficient (such as the lack of lake information like water level data), or when it is not necessary to provide all lake information, the calculation requirements can still be achieved. In a reasonable lake area division result, the change in the water storage volume of each partition in the sense of daily average (that is, the flow profit and loss amount calculated by partition) should match the incoming water volume of that partition, and there should be no situation where the water storage in some partitions is rapid while the water storage in other partitions is slow. This situation is a non-steady process and cannot be used as the representative state of each partition in the sense of daily average. Therefore, this embodiment also designs the above-mentioned second index related to the daily water flow between partitions.

[0068] In step S3, two means are provided for optimizing and adjusting the partition: adjusting the partition water level and adjusting the partition volume. The reason for proposing these two means is that the partitions obtained by dividing the lake area in step S1 may have certain irrationalities in either the planar division of the partition (referring to planar geometric properties such as area size, etc.) or the setting of the water level of each partition. As mentioned above, the water depth set {d} of each partition, the water level set {h} of each partition, and the outlet weir height h o have a certain correlation, and these data are the key data for determining the water flow between partitions. Adjusting the water level is to adjust the water level set {h} of each partition, which is to correct and adjust the water level estimated value, and the water level is preferentially adjusted when the first index does not meet the requirements. Adjusting the volume is a systematic adjustment of multiple data including the water depth set {d} of each partition, the planar area set {As} of each partition, the volume set {V} of each partition, and the adjacent side length matrix {Al i,j} of each partition. Under the constraint condition that the total volume of the lake remains unchanged, adjusting the volume will also achieve the purpose of adjusting the water level set {h} of each partition, and the volume is preferentially adjusted when the second index does not meet the requirements.

[0069] Continuing the previous example, the optimization and adjustment of the partition performed in step S3 is a constrained partition optimization and adjustment, and the optimization and adjustment should not violate the following basic constraints on water volume and water flow: Optimization and adjustment constraint 1: The inflow and outflow flows of each part of the lake are consistent with the basin simulation results, that is, the upstream and downstream flow constraints; Optimization and adjustment constraint 2: If there is a flow between the i-th partition and the j-th partition, then the water flow Q between these two partitions i,jThe water levels, positions (distance of the center point), geometric parameters (including depth, planar width, etc.) of these two partitions conform to the Chezy formula, i.e., the hydraulic constraints; Optimization and adjustment constraint 3: The difference between all the inflows (including the inflows from the inflow river section and the inflows from other partitions) and all the outflows (including the outflows from the outlet and the outflows to other partitions) of each partition is consistent with the water volume change of this partition, i.e., the mass constraint.

[0070] For example Figure 5 As shown, during the optimization and adjustment of the partitions, The exemplary description of the first index is as follows: After obtaining the geometric parameters of each partition in step S1, based on the Chezy formula, calculate the daily water flow of each partition, and then summarize and calculate the calculated value of the daily flow profit and loss of the overall lake simulation, and compare it with the actual daily flow profit and loss (i.e., , wQ k is the actual daily inflow of the k-th water inlet of the lake, InFSet is the number value of the water inlets of the lake, wQ o is the actual daily outflow of the lake), find their difference, and this difference should be within a certain threshold range. Because the daily water flow between partitions does not involve the increase or decrease of the overall water volume of the lake on the same day, in fact, for the daily flow profit and loss of the overall lake, it is only necessary to compare the calculated value Q o of the daily outflow of the outermost partition with the actual value wQ o The difference. Thus, this embodiment designs the first index , and sets the threshold according to 10% of the actual daily flow profit and loss of the overall lake (here, the first preset ratio value is exemplified by 10% and can be adjusted appropriately according to needs):

[0071] The exemplary description of the second index is as follows: Not only should the daily flow profit and loss of the overall lake simulation be consistent with the actual situation, but in a reasonable lake area division result, it should also ensure that the daily flow profit and loss of each partition is related to its geometric dimensions (such as the water depth d i of the partition, the planar area As i of the partition, the sum of the adjacent side lengths Al i of the partition and other partitions, etc.) and the water level h iMatch. This is essentially a process of solving the geometric shapes of each sub-region - if the shape information of a lake and its sub-regions is rich enough, in theory, all parameters of all sub-regions and the flow rates between sub-regions can be solved. However, this embodiment specifically targets the situation where the lake information is insufficient (such as the lack of lake information like water level data), or all lake information does not have to be provided, and the calculation requirements can still be met. In a reasonable result of lake region division, the change in the water storage volume of each sub-region in the sense of daily average (i.e., the flow profit and loss calculated by sub-region) should match the incoming water volume of that sub-region, and there should be no situation where the water storage in some sub-regions increases rapidly while the water storage in other sub-regions increases slowly. Such a situation is a non-steady process and cannot be used as the representative state of each sub-region in the sense of daily average. Therefore, a second index is designed :

[0072] For each sub-region, calculate the total incoming flow rate for the day by summarizing the flow rates of all sub-regions with which it has a flow connection and the total outgoing flow rate for the day , and calculate the index :

[0073] where, wQ i represents the incoming flow rate of the i-th sub-region at the water inlet on the day, represents the sum of the calculated values of the incoming flow rates from all sub-regions with which the i-th sub-region has a flow connection to the i-th sub-region on the day, represents the sum of the calculated values of the outgoing flow rates from the i-th sub-region to all sub-regions with which it has a flow connection on the day; Ie i is a 0-1 binary variable, which takes the value of 1 when the i-th sub-region is an outlet sub-region, and 0 otherwise

[0074] To accelerate the convergence rate of the second index , this embodiment adds a weight wt to the index i of each sub-region. This weight reflects the importance of this sub-region in terms of "receiving and sending" flow rates, and its setting method is as follows:

[0075] where, represents the sum of of each sub-region

[0076] In a possible implementation manner, in step S3, the optimizing and adjusting the sub-region based on the optimization index by adjusting the water level of the sub-region or adjusting the volume of the sub-region includes: Judge whether the first index is less than the first threshold and whether the second index is less than the second threshold: If the first index is greater than or equal to the first threshold and the second index is less than the second threshold, then optimize and adjust the partition by adjusting the water level of the partition on the current day. After the optimization and adjustment, recalculate the calculated value of the water flow on the current day between the partitions and the calculated value of the water flow out of the water outlet partition on the current day, and then re-judge; If the first index is less than the first threshold and the second index is greater than or equal to the second threshold, then optimize and adjust the partition by adjusting the volume of the partition. After the optimization and adjustment, recalculate the calculated value of the water flow on the current day between the partitions and the calculated value of the water flow out of the water outlet partition on the current day, and then re-judge; If the first index is greater than or equal to the first threshold and the second index is greater than or equal to the second threshold, then optimize and adjust the partition by adjusting the water level of the partition when the first index is greater than or equal to the second index, and optimize and adjust the partition by adjusting the volume of the partition when the first index is less than the second index. After the optimization and adjustment, recalculate the calculated value of the water flow on the current day between the partitions and the calculated value of the water flow out of the water outlet partition on the current day, and then re-judge; If the first index is less than the first threshold and the second index is less than the second threshold, then end the optimization and adjustment.

[0077] Continuing the foregoing example, as Figure 5 shown, in the process of iteratively adjusting the water level or volume to optimize and adjust the partition, the first index should be reduced and the second index as the goal. The user can set the first threshold of the first index and the second threshold of the second index according to the actual situation. Considering that this embodiment is for lake partition simulation in the case of insufficient information, where the prepared characteristics and required accuracy of each partition, for example, the first threshold and the second threshold are respectively set to 0.1. When the first index and the second index are respectively less than 0.1, it can be considered that the partition reaches a reasonable level, and the above-mentioned "constrained partition optimization and adjustment" link is completed.

[0078] In a possible implementation manner, in step S3, the adjusting the water level of the partition includes: Adjust the water level of the partition with the largest absolute value of the daily flow profit and loss value in each partition to reduce the first index; or, Batch-adjust the water levels of each partition to reduce the first index.

[0079] Continuing the foregoing example, if the setting or estimation of the partition water level is unreasonable, it will inevitably and systematically affect the simulated daily flow of the partition and will inevitably cause the lake's daily outflow (i.e., the daily water flow out of the outlet partition) to deviate from the actual situation. In this case, adjusting the geometric attributes of each partition has little effect. Therefore, to reduce the first index Taking this as the main goal, adjusting the water level is used as the main iterative optimization operation.

[0080] For example Figure 5 As shown, when determining that the first index is greater than or equal to the second index , adjust the water levels {h} of each partition. Specifically, the reason for the incorrect (or inaccurate) water level of each partition is uncertainty. If the user believes that the water levels of each partition are more accurate (for example, this determination can be made when part of them belongs to measured values or estimates based on historical measured values) and the water outlet elevation of the water outlet partition is inaccurate (such as only a rough estimate), then the water levels of individual partitions should be adjusted; conversely, the water levels of all partitions should be adjusted batchwise (either raised or lowered together). That is When determining that the water level data is relatively accurate, adjust the water levels of individual partitions; When determining that the water level data is relatively inaccurate, adjust the water levels of all partitions batchwise.

[0081] Among them, when adjusting the water levels of individual partitions, the partitions whose water levels should be adjusted can be identified according to the principle of the largest absolute value of the daily flow profit and loss, that is .

[0082] Whether adjusting the water levels of individual partitions or adjusting the water levels of all partitions batchwise, the adjustment (raising or lowering) of the water level should be in the direction of reducing the first index . The water level adjustment amount (or adjustment step) δh for each water level adjustment can be reasonably set based on actual experience on the basis of balancing the number of iterative steps and stability. For example, an empirical relationship between the water level adjustment amount δh and the first index can be constructed to adaptively and dynamically determine the water level adjustment amount δh for each water level adjustment.

[0083] In a possible implementation manner, in step S3, the adjusting the partition volume includes: Taking the partition with the largest first weighted ratio value as the first partition, selecting a partition that has a water flow exchange with the first partition as the second partition, and dividing a first set ratio of the volume of the first partition into the second partition or dividing a first set ratio of the volume of the second partition into the first partition.

[0084] Continuing with the previous example, adjusting the partition volume is an effective means to change the partition size and thus improve the daily flow between adjacent partitions.

[0085] For example Figure 5 As shown, when determining that the first index is less than the second index , adjust the partition volume (and its geometric properties) to reduce the second index is the main goal. Specifically, physical reasons for volume adjustment are required when the partition volume is incorrect (or inaccurate), or the partition shape deviates significantly from an equal-area circle. This embodiment proposes a practical strategy, that is, to select the partition with the largest value as the first partition mx, that is . Then, select a partition that has a water flow exchange (or a flow connection) with the first partition mx as the second partition mi, and divide a first set ratio of the volume δV of the first partition mx to the second partition mi or divide a first set ratio of the volume δV of the second partition mi to the first partition mx, so as to achieve the purpose of reducing the second index .

[0086] In a possible implementation manner, in step S3, the selecting a partition that has a water flow exchange with the first partition as the second partition includes: Select the partition with the largest daily water flow exchange value with the first partition as the second partition; or Select the partition that has the greatest impact on the daily flow profit and loss value of the first partition as the second partition.

[0087] In a possible implementation manner, in step S3, the adjusting the partition volume includes: Select the partition with the largest daily water flow exchange value with the first partition as the second partition and perform iterative volume division until the second index is less than the third threshold or the number of iterations reaches the fourth threshold; Execute the volume division cycle for the first set number of times. A volume division cycle includes: select the partition that has the greatest impact on the daily flow profit and loss value of the first partition as the second partition and perform iterative volume division for the second set number of times, and then select the partition with the largest daily water flow exchange value with the first partition as the second partition and perform iterative volume division for the third set number of times.

[0088] Continuing the previous example, as Figure 5 shown, this embodiment proposes two strategies for selecting the second partition mi: the FM strategy and the BM strategy.

[0089] Under the FM strategy, select the dual partition (the dual partition means a partition that has a flow exchange or a flow connection) with the largest daily flow exchange value with the first partition mx as the second partition mi. Among them, for example, if the first partition mx itself is an out-flow partition, then the second partition mi is selected as the largest daily inflow source partition of the first partition mx, otherwise the second partition mi is selected as the largest daily out-flow sink partition of the first partition mx (that is, the dual partition with the largest out-flow of the first partition mx).

[0090] Under the BM strategy, the dual partition that makes the greatest contribution to the daily flow profit and loss of the first partition mx is selected as the second partition mi. For example, if the outflow of the first partition mx itself is greater than the inflow, the second partition mi is selected as the partition with the largest daily outflow sink of the first partition mx; otherwise, the second partition mi is selected as the partition with the largest daily inflow source of the first partition mx. The inventors found through experiments that the effects of continuously using or simply alternating the FM strategy and the BM strategy in iterations are not good, the utility is not obvious, the convergence speed is slow, and it may also fall into local optimality. Therefore, this embodiment proposes a strategy of randomly using the FM strategy and the BM according to the random mechanism II, which can achieve good effects of improving the convergence speed and not easily falling into the local optimality trap. The design description of the random mechanism II is as follows: First, in the initial stage of iterative adjustment, the FM strategy is used to select the second partition mi for iterative adjustment of the volume until the second index <0.3 or the number of iterations reaches 300 times; Then, the BM strategy is used to select the second partition mi for iterative adjustment 3 times, and then the FM strategy is used to select the second partition mi for iterative adjustment 3 or 4 times. The above is an iteration cycle executed 6 - 7 times.

[0091] In a possible implementation manner, in step S3, the volume value of a single volume division is set according to the maximum first weighted ratio value of each partition, and is reduced by the second set ratio when the first constraint condition is not satisfied. The first constraint conditions include: The volume value of a single volume division is less than or equal to the product of the second preset ratio value and the volume of the first partition and the product of the second preset ratio value and the volume of the second partition; The volume of the first partition after volume division is greater than or equal to the product of the third preset ratio value and the volume of the first partition; The water level change value of the first partition after volume division and the water level change value of the second partition after volume division are less than or equal to the fifth threshold.

[0092] In a possible implementation manner, in step S3, the adjustment of the partition volume further includes: According to the volume of the first set ratio of the first partition divided into the second partition or the volume of the first set ratio of the second partition divided into the first partition, the water levels and geometric parameters including the adjacent side length and the planar area of the first partition and the second partition are updated respectively.

[0093] Continuing the foregoing example, as Figure 5 shown, the constraints in the iterative adjustment of the volume are designed as follows: After each iterative adjustment of the divided volume, the new and old volume ratios θ of the first partition mx and the second partition mi are calculated respectively, and according to θ 1 / 3Scale the water depth d of the updated partition and the adjacent side length Al between the first partition mx and the second partition mi mi,mx , and scale the planar area As of the updated partition by θ 2 / 3 , specifically as follows: , , , ,

[0094] After the above scaling, when the second partition mi and the first partition mx are still fitted with equal-area circles, the sum of their volumes is inconsistent with the adjusted division volume. The entire lake area can be regarded as a whole and normalized based on the total volume of the lake area:

[0095] Among them, represents the sum of the volumes of all partitions before this adjusted division volume, represents the sum of the volumes of partitions other than the first partition mx and the second partition mi before this adjusted division volume.

[0096] In this way, the lake volume is ensured to be consistent before and after various fitting assumptions and volume adjustment operations.

[0097] Dividing the volume δV of the first set ratio of the first partition mx to the second partition mi or dividing the volume δV of the first set ratio of the second partition mi to the first partition mx will result in a systematic change of "volume division → volume change of the second partition mi and the first partition mx → water level change of the second partition mi and the first partition mx (water depth change will cause water level change) → flow change between the second partition mi and the first partition mx and their respective adjacent partitions → changes in each partition". Therefore, the volume division amount δV should be reasonably set and can be dynamically adjusted. The initial value of the volume division amount δV can be tentatively set according to the value of of the first partition mx, and the value of the volume division amount δV is lowered when the volume division constraint condition is not satisfied. For example, Figure 5 as shown in, the value of the volume division amount δV is adjusted to 2 / 3 of the current value (for example, rounded to Figure 5 0.67 as shown in) and then the volume division is performed again. The volume division constraint conditions proposed in this embodiment are, for example, the following conditional expressions: Conditional expression (1), δV ≤ 0.1min(V mx,old , V mi,old ); Conditional expression (2), V mx,new = V mx,old - δV ≥ 0.1V mx,original ; Conditional expression (3), max(|h mx,new - h mx,old, h mi,new -h mi,old |) ≤ 0.01 m.

[0098] That is, the single - volume division quantity δV shall not be greater than 10% of the second sub - region mi and the original volume of the first sub - region. After the volume division, the new volume of the first sub - region mx shall not be less than 10% of the volume of the first sub - region (V mx,original ) during the initial set - up of the sub - regions, and shall not cause the water levels of the first sub - region mx and the second sub - region mi to change by more than 1 cm. These three constraints reflect the requirements for the robustness of volume adjustment, the consideration of various factors (with a practical basis) in the sub - region demarcation in step S1, and the requirements for the robustness of the exchange flow between sub - regions.

[0099] S4. Calculate the chemical substance concentration of each adjusted sub - region on the same day according to the chemical substance data from the watershed hydrological model.

[0100] In a possible implementation, in step S4, the calculation of the chemical substance concentration of each adjusted sub - region on the same day includes: Judge whether the actual value of the lake effluent flow on the same day is zero: If not, calculate the chemical substance concentration of each adjusted sub - region according to the one - dimensional advection - diffusion equation; If so, calculate the chemical substance concentration of each adjusted sub - region according to the one - dimensional advection - diffusion equation combined with the return flow attenuation based on the set number of return times.

[0101] Continuing with the previous example, first, the principle of the one - dimensional advection - diffusion equation is explained: The differential equation of the one - dimensional advection - diffusion equation (ADE) is as follows:

[0102] Among them, C is the chemical substance concentration, which is a function of time t and spatial position x, with the unit of kg / m 3 ; u is the one - dimensional flow velocity, with the unit of m / s; D is the turbulent diffusion coefficient; k is the degradation rate constant, with the unit of 1 / s.

[0103] In the lake sub - region simulation, the following boundary conditions apply:

[0104]

[0105] Among them, C0 is the average chemical substance concentration along the flow path at time 0, with the unit of kg / m 3 ; C in is the average chemical substance concentration at the inlet cross - section, with the unit of kg / m 3 .

[0106] The analytical solution of the one-dimensional advection-diffusion equation (ADE) under the above boundary conditions is as follows:

[0107] where erfc is the complementary error function.

[0108] The following introduces the application of the one-dimensional advection-diffusion equation (ADE) in this embodiment.

[0109] (1) The one-dimensional advection-diffusion equation (ADE) is used to estimate the aqueous phase concentration in the partition In actual operation, the total aqueous phase concentration of this partition simulated yesterday is used as the average concentration C0 at the time delay process of time 0, and the total concentration of the aqueous phase chemical substances in the sources flowing into this partition (including the inflow river section, and the partitions flowing into this partition determined according to the flow matrix Q mat is used as C in . Substituting the flow velocity u and the turbulent diffusion coefficient D calculated according to the Chezy formula into the above analytical formula (d), the chemical substance concentration C(x, t) at each time and location in this partition on the current day can be calculated. Of course, in actual calculations, if the vertical mixing diffusion coefficient D m caused by meteorological conditions is considered, then the turbulent diffusion coefficient D in the above analytical formula (d) should be substituted with the sum of the turbulent diffusion coefficient D caused by the water flow between partitions and the vertical mixing diffusion coefficient D m .

[0110] For example, the method of approximating the partition by an equal-area circle as shown, etc., can be continued. Taking the radius R of the equal-area circle of this partition as x, C(R, t) is obtained as the concentration of the representative concentration of the partition at time t; further, one day is decomposed into 24 h, and C(R, t) for each hour is calculated, and its average value is taken as the average concentration of the partition on the current day: Figure 4 where C

[0111] represents the average chemical substance concentration of the j-th partition on the current day, and R j represents the radius of the equal-area circle of the j-th partition. j

[0112] Assume that parameters related to hydrodynamic conditions such as the flow velocity u and the turbulent diffusion coefficient D and the biodegradation rate k are relatively stable within a day and are generally acceptable. Therefore, it is acceptable to use the formula (e) of the average concentration of the above partition on the current day averaged over 24 h as the daily value. However, if the chemical substance has obvious photodegradation properties, then it is not appropriate to set the average chemical substance concentration C j, instead, the biodegradation rate k should be made to vary with time t in the above analytical formula - specifically, the hourly incident light intensity on the water surface, the extinction coefficient of the lake body, etc. can be calculated based on the calendar and geographical location, which will not be further specified in this embodiment.

[0113] (2) Estimating the dissolved concentration in the water phase of each sub-region There is water flow between sub-regions. Not only can the chemical substances be approximated by the advection-diffusion equation, but the suspended solids (also known as suspended sediment, Suspended Solids) in the water phase (watertable) can also be set according to formulas (a)-(d). Here, the degradation rate constant is not for degradation, but for the settling removal rate (Settling Rate) of the suspended solids. In addition, the advection-diffusion equation approximation of the suspended solids between sub-regions is based on the actual inflow sediment volume wSed of the k-th water inlet of the lake in the daily hydrological data from the watershed hydrological model. k . Representing the suspended solid concentration by S and the suspended solid settling rate by ks, the concentration expression and the average daily concentration expression for the suspended solids can be listed:

[0114]

[0115] Among them, S j is the daily suspended solid concentration of the j-th sub-region, S0 is the average suspended solid concentration along the flow path at time 0, and S in is the average suspended solid concentration at the inlet cross-section.

[0116] With the support of the daily suspended solid concentration S j in the j-th sub-region, the daily chemical substance concentration C j in the j-th sub-region can be further calculated to obtain the dissolved concentration in the water phase of this sub-region through the equilibrium distribution method - using the organic carbon-water partition coefficient K oc of the chemical substance and the proportion f oc of the organic carbon in the suspended solids.

[0117] The following introduces the calculation of the chemical substance concentration in each sub-region by applying the one-dimensional advection-diffusion equation combined with the return flow attenuation under the water storage mode (i.e., the actual value of the lake's daily water outflow is zero).

[0118] In the case where the lake has an outflow, the one-dimensional advection-diffusion equation can be directly used to calculate the chemical substance concentration in each sub-region. However, if the lake has no outflow on a certain day (i.e., wQ o= 0 (the lake is in the "water storage mode"), there may be a large deviation in the result of directly solving the one-dimensional advection-diffusion equation. This is because the simulation in this embodiment assumes that the flow between each partition reaches a relatively stable state within a day, and the diffusion coefficient is estimated based on this state, and the diffusion result is the concentration distribution of chemical substances. However, in the case of the water storage mode (there is still inflow from the incoming river), the flow carried by the incoming river can easily reach some streamlines to the outlet partition within a short time (such as several hours). When it reaches the outlet partition but does not leak out, the diffusion of these chemical substances cannot fully reflect the result of one-dimensional advection-diffusion, but the partition near the inlet of the lake can still partially reflect the aforementioned one-dimensional advection-diffusion result.

[0119] Therefore, this embodiment proposes to approximate the diffusion estimation of chemical substances in the water storage mode with "virtual flow of reverse attenuation". The estimated concentration distribution of chemical substances in each partition is denoted as {C o=0}, where the subscript o = 0 is the abbreviation of out = 0, which means that the lake has no outflow or the outflow is zero. The specific steps and designs are as follows: Based on the water level {h} of the previous day in the "water storage mode", the optimized partition geometric parameters and the water flow matrix {Q i,j} and the turbulent diffusion coefficient matrix {D i,j} between each partition are calculated. The one-dimensional advection-diffusion equation is conventionally applied to obtain the estimated concentration distribution {C} of chemical substances in each partition. To show the characteristics in the water storage mode, the results of this round of calculation are marked with subscript 0: {Q i,j}0, {D i,j}0, {C}0.

[0120] The first round of reverse calculation. At this time, it is set that the flow between each partition is completely "opposite" to the flow of the previous round marked as 0, and the flow is weakened. The calculation results of this round are marked with subscript 1. According to experience and the principle of potential energy dissipation analysis, it is more reasonable to calculate the reverse flow value of this round as 40% of the flow of the previous round, but in the opposite direction, that is

[0121] At the same time, the turbulent diffusion coefficient of this round is also reduced by 40% of the previous round (there is no direction problem), that is

[0122] Then, with {C}0 as the initial state, the numerical values of the water flow matrix {Q i,j}1 between each partition and the turbulent diffusion coefficient matrix {D i,j}1 after the first round of reverse calculation are substituted into the one-dimensional advection-diffusion equation to calculate the new concentration distribution {C}1 of chemical substances in each partition.

[0123] Second round of return flow calculation. At this time, it is set that the flow rates in each sub-region are completely "opposite" to those in the previous round, and the flow rates are weakened. The return flow rate value in this round is calculated as 70% of the flow rate in the previous round, but in the opposite direction, that is

[0124] At the same time, the turbulent diffusion coefficient in this round is also converted to 70% of that in the previous round (there is no direction problem), that is

[0125] Then, taking {C}1 as the initial state, substituting the water flow rate matrix {Q i,j}2 between each sub-region after the second round of return flow calculation and the numerical value of the turbulent diffusion coefficient matrix {D i,j}2 after the second round of return flow calculation into the one-dimensional advection-diffusion equation to calculate the new chemical substance concentration distribution {C}2 in each sub-region.

[0126] After setting the return virtual flow and applying the one-dimensional advection-diffusion in the above first round and second round, good results can usually be obtained, that is, the influence of the sub-regions far from the lake inlet on the lake inlet gradually weakens and converges to an almost completely mixed concentration distribution.

[0127] This way of applying the one-dimensional advection-diffusion with the return virtual flow, because the water flow rate matrix {Q i,j} between each sub-region, the turbulent diffusion coefficient matrix {D i,j}, the estimated chemical substance concentration distribution {C} in each sub-region and its calculation logic are set based on the original variables, is relatively simple in program implementation and has good application value.

[0128] The inventor found that generally, two rounds of return (as described above) can meet the requirements. If further return calculations are carried out, it may instead lose the approximation ability for the phenomenon that "the influence of the sub-regions far from the lake inlet on the lake inlet gradually weakens". At this time, it may be better to directly simulate with a completely mixed simplified model on the "water storage day" (that is, give up the pursuit of the concentration distribution differences in the sub-regions and return to the level of simulating the lake as a single unit in the basin simulation). Of course, for specific simulation cases, the inventor believes that appropriately increasing the return may also be beneficial.

[0129] In summary, in the lake chemical substance concentration simulation method provided in this embodiment, the lake area division can scientifically divide the lake area, and the obtained sub-areas are reasonable in terms of chemical substance migration and diffusion, etc., meeting the requirements of chemical substance concentration management. Moreover, even in the case of missing lake information such as water level data, the scientificity and rationality of the sub-area results can be ensured by optimizing and adjusting the sub-areas, so as to accurately simulate the environmental exposure concentration of chemical substances in the lake. Based on the design of this embodiment for the requirements of chemical substance concentration management, the lake model can be directly connected to the watershed hydrological model, and the lake chemical substance concentration can be simulated without increasing the time resolution of the watershed hydrological model. The calculation speed can be improved while ensuring the simulation quality of the sub-area chemical substance concentration simulation, and it is beneficial for long-term simulation.

[0130] The lake chemical substance concentration simulation method provided in this embodiment will be further described below through a specific example.

[0131] Taking the simulation of a certain lake (denoted as DL) in a certain pilot watershed (denoted as YD) as an example. The method provided in this embodiment is implemented by writing a program in Fortran code to form a relatively independent module Module LKZ, which is connected and mixed with the watershed-scale environmental exposure simulation model SWAT-KM model. The SWAT-KM model containing the Module LKZ module was used to conduct sub-area simulation on the DL lake in the YD watershed, realizing all the contents of this embodiment, and realizing the systematic simulation of multiple media such as soil, river sections, air, and vegetation in the watershed.

[0132] Other watershed data (including meteorological data) and parameters required in SWAT-KM are omitted here. Only the inputs and outputs required for the sub-area simulation of the DL lake involved in this embodiment will be described.

[0133] As Figure 6 shown, according to the inflow and outflow conditions, the DL lake area is divided into 7 sub-areas, which correspond to the local monitoring points and reflect the management requirements of lake chemical substances. Sub-area 5 is the outlet sub-area where the single outlet of the lake is located. Sub-areas 2, 3, and 7 all have multiple inflowing rivers, while sub-area 6 has no inflow and outflow. The size of the sub-areas takes into account the magnitude of the average annual inflow as much as possible. Except for sub-areas 4 and 7, the other sub-areas are relatively square.

[0134] The water depths of each sub-region of DL Lake are estimated based on the publicly available lake bed elevation, average water level, and water depth sub-region estimation; data such as the elevation and width of the outlet weir are available in public information; the water level elevations of each sub-region are set according to the flow direction and the measured public data at some points and in some seasons, with a certain error. The estimated error level is approximately no more than 10 cm and is set at 0.05 m. Since the elevation data of the outlet weir is more accurate, when adjusting the water level strategy in this embodiment, the method of batch-adjusting the water level is mainly adopted.

[0135] Based on the above simple data setting, the SWAT-KM model including the LKZ Module can achieve the daily multi-media simulation of the YD Basin, including the simulation of the DL Lake sub-region. For example, from the beginning of 2008 to the end of September 2024, the simulation proceeds smoothly. The comparison curves of the simulated values and the actual values of the daily concentration changes of a certain harmful chemical substance in each sub-region during some sample periods are respectively as Figure 7 shown in the curve groups of Lake Zone#1, Lake Zone#2, Lake Zone#3, Lake Zone#4, LakeZone#5, Lake Zone#6, and Lake Zone#7. Among the seven curve groups: the colored curve is the simulated value of the daily concentration of a certain chemical substance in the sub-region, and the gray curve is the simulated value of the daily concentration of the corresponding chemical substance in other sub-regions (as a control, showing the concentration dynamics of cooperation but not equal concentrations in different sub-regions of the lake). Figure 7 As shown, the abscissa is the day, and the ordinate is the chemical substance concentration, with the unit of mg / m 3 .

[0136] After comparing with the actual monitoring results, the concentration on the same day and at the same location simulated in this embodiment (indicating the sampling point concentration with the concentration of the sub-region and river section) is in good agreement with the monitoring sampling concentration. According to the daily results of the simulation of the lake sub-region and the daily simulation results of other units in the basin that are synergistically simulated with it, the average value within a certain statistical period can be drawn. As an example, Figure 8 shows the quarterly average value of the water body concentration in the YD Basin (including DL Lake and its surrounding river sections), Figure 8 as shown, the colors of the sub-region and the river section represent the concentration values. The redder the color, the greater the concentration value. The digital numbers such as 16, 19, 20, 21, 43, 28, 39, 47, etc. are the unit numbers in the simulation.

[0137] Using the SWAT-KM model linked with the Module LKZ module, and calculating based on the start and end times of the code marked by the CPU of a personal computer (configuration: CPU Intel i7-10700 2.9 GHz, 4 cores, 16G of memory, 2T hard drive, windows 10), it saves 93% of the time compared to the SWAT-KM model linked with the EFDC model (time resolution set to 10 min). Moreover, when comparing the simulation results of the SWAT-KM model linked with the EFDC model with those of the SWAT-KM model containing the Module LKZ module, there is no obvious advantage, that is, there is no better compliance with the monitoring data. This shows that under the condition of limited input information, the lake zoning simulation method designed in this embodiment can, while ensuring the simulation quality, be much faster in calculation speed than the algorithm of the watershed model linked with the hydrodynamic model, providing a better basis for carrying out long-term continuous operation simulations, which is exactly what is needed in the environmental exposure assessment and environmental risk assessment of chemical substances.

[0138] As Figure 9 shown, a computer system suitable for implementing the method for simulating the environmental exposure concentration of chemical substances in lakes that meet the watershed scale provided by the above embodiment includes a central processing module (CPU), which can perform various appropriate actions and processes according to the program stored in the read-only memory (ROM) or the program loaded from the storage part into the random access memory (RAM). In the RAM, various programs and data required for the operation of the computer system are also stored. The CPU, ROM, and RAM are connected to each other through a bus. The input / output (I / O) interface is also connected to the bus.

[0139] The following components are connected to the I / O interface: an input part including a keyboard, a mouse, etc.; an output part including, for example, a liquid crystal display (LCD), etc. and a speaker, etc.; a storage part including a hard disk, etc.; and a communication part including a network interface card such as a LAN card, a modem, etc. The communication part performs communication processing via a network such as the Internet. A drive is also connected to the I / O interface as needed. Removable media, such as magnetic disks, optical disks, magneto-optical disks, semiconductor memories, etc., are installed on the drive as needed so that the computer programs read from them can be installed into the storage part as needed.

[0140] Particularly, according to this embodiment, the process described in the above flowchart can be implemented as a computer software program. For example, this embodiment includes a computer program product, which includes a computer program tangibly contained on a computer-readable medium, and the above computer program includes program codes for executing the method shown in the flowchart. In such an embodiment, the computer program can be downloaded and installed from the network through the communication part, and / or installed from the removable media.

[0141] The flowcharts and schematic diagrams in the accompanying drawings illustrate the possible architectures, functions, and operations of the system, method, and computer program product of this embodiment. In this regard, each block in the flowchart or schematic diagram can represent a module, a program segment, or a part of the code, and the above-mentioned module, program segment, or part of the code contains one or more executable instructions for implementing the specified logical function. It should also be noted that in some alternative implementations, the functions marked in the blocks may occur in a different order than that marked in the accompanying drawings. For example, two consecutive blocks shown may actually be executed substantially in parallel, and they may sometimes be executed in the reverse order, depending on the functions involved. It should also be noted that each block in the schematic diagram and / or flowchart, as well as the combination of blocks in the schematic and / or flowchart, can be implemented by a dedicated hardware-based system that performs the specified functions or operations, or can be implemented by a combination of dedicated hardware and computer instructions.

[0142] As another aspect, this embodiment also provides a non-volatile computer storage medium, which can be the non-volatile computer storage medium included in the above-mentioned device in the above-mentioned embodiment, or can exist separately and be unassembled into the terminal. The above non-volatile computer storage medium stores one or more programs, and when the above one or more programs are executed by a device, the device implements the method for simulating the environmental exposure concentration of chemical substances in lakes at the basin scale provided in the above-mentioned embodiment.

[0143] As another aspect, this embodiment also provides a computer program product, including a computer program, which when executed by a processor implements the method for simulating the environmental exposure concentration of chemical substances in lakes at the basin scale provided in the above-mentioned embodiment.

[0144] Obviously, the above-mentioned embodiments of the present disclosure are merely examples for clearly illustrating the present disclosure, rather than limitations on the implementation manners of the present disclosure. For those of ordinary skill in the art, other different forms of changes or variations can be made based on the above description. It is impossible to list all the implementation manners here. Any obvious changes or variations derived from the technical solutions of the present disclosure still fall within the protection scope of the present disclosure.

Claims

1. A method for simulating the environmental exposure concentration of chemical substances in lakes at the basin scale, characterized in that Including: Dividing the lake into multiple sub - regions according to the division rules, recording the geometric parameters of each sub - region. The division rules include a first rule and a second rule. The first rule is that the only water outlet of the lake and the water inlet of the lake are divided into different sub - regions. The second rule is that the regional division meets the requirements of chemical substance concentration management; Calculating the calculated value of the daily water flow between each sub - region and the calculated value of the daily water discharge of the water - discharging sub - region including the water outlet according to the daily hydrological data from the watershed hydrological model, the estimated value of the water level of each sub - region, and the geometric parameters of each sub - region; Based on the optimization index, optimizing and adjusting the sub - regions by adjusting the water level of the sub - regions or adjusting the volume of the sub - regions to obtain the adjusted sub - regions. The optimization index includes a first index related to the calculated value of the daily water discharge of the water - discharging sub - region and a second index related to the daily water flow between each sub - region; Calculating the daily chemical substance concentration of each adjusted sub - region according to the chemical substance data from the watershed hydrological model.

2. The method according to claim 1, wherein The first index is the ratio of a first value to a second value. The first value is the difference between the calculated value of the daily water discharge of the water - discharging sub - region and the actual value of the daily water discharge of the lake. The second value is the product of a first preset ratio value and the actual value of the daily flow profit and loss of the lake; and / or, The second index is the sum of the first weighted ratio values of each sub - region. The first weighted ratio value of a sub - region is the product of the first ratio value of the sub - region and the weight. The first ratio value of the sub - region is the absolute value of the ratio of the daily flow profit and loss value of the sub - region to the daily water inflow value of the sub - region. The weight is the ratio of the daily water inflow value of the sub - region to the sum of the daily water inflow values of each sub - region.

3. The method according to claim 2, wherein The optimizing and adjusting the sub - regions by adjusting the water level of the sub - regions or adjusting the volume of the sub - regions based on the optimization index includes: Judging whether the first index is less than a first threshold and whether the second index is less than a second threshold: If the first index is greater than or equal to the first threshold and the second index is less than the second threshold, then optimize and adjust the sub - regions by adjusting the daily water level of the sub - regions. After the optimization and adjustment, recalculate the calculated value of the daily water flow between each sub - region and the calculated value of the daily water discharge of the water - discharging sub - region and re - judge; If the first index is less than the first threshold and the second index is greater than or equal to the second threshold, then optimize and adjust the sub - regions by adjusting the volume of the sub - regions. After the optimization and adjustment, recalculate the calculated value of the daily water flow between each sub - region and the calculated value of the daily water discharge of the water - discharging sub - region and re - judge; If the first index is greater than or equal to the first threshold and the second index is greater than or equal to the second threshold, then optimize and adjust the sub - regions by adjusting the water level of the sub - regions when the first index is greater than or equal to the second index, and optimize and adjust the sub - regions by adjusting the volume of the sub - regions when the first index is less than the second index. After the optimization and adjustment, recalculate the calculated value of the daily water flow between each sub - region and the calculated value of the daily water discharge of the water - discharging sub - region and re - judge; If the first index is less than the first threshold and the second index is less than the second threshold, then end the optimization and adjustment.

4. The method according to any one of claims 1 to 3, characterized in that The adjusting the water level of the sub - regions includes: Adjust the water level of the partition with the largest absolute value of the daily flow profit and loss value in each partition to reduce the first index; or, Batch adjust the water levels of each partition to reduce the first index.

5. The method according to any one of claims 1-3, characterized in that, The adjusting partition volume includes: Take the partition with the largest first weighted ratio value as the first partition, select a partition that exchanges water flow with the first partition as the second partition, and divide a first set ratio of the volume of the first partition into the second partition or divide a first set ratio of the volume of the second partition into the first partition.

6. The method according to claim 5, wherein The adjusting partition volume further includes: According to the volumes of the first partition and the second partition after dividing a first set ratio of the volume of the first partition into the second partition or dividing a first set ratio of the volume of the second partition into the first partition, update the water levels and geometric parameters including the adjacent side lengths and planar areas of the first partition and the second partition respectively.

7. The method according to claim 5, characterized in that The selecting a partition that exchanges water flow with the first partition as the second partition includes: Take the partition with the largest daily water flow exchange value with the first partition as the second partition; or, Take the partition that has the greatest impact on the daily flow profit and loss value of the first partition as the second partition.

8. The method according to claim 7, wherein The adjusting partition volume includes: Take the partition with the largest daily water flow exchange value with the first partition as the second partition and iterate the volume division until the second index is less than the third threshold or the number of iterations reaches the fourth threshold; Execute the volume division cycle for the first set number of times. A volume division cycle includes: take the partition that has the greatest impact on the daily flow profit and loss value of the first partition as the second partition and iterate the volume division for the second set number of times, and then take the partition with the largest daily water flow exchange value with the first partition as the second partition and iterate the volume division for the third set number of times.

9. The method according to claim 7, characterized in that The volume value of a single volume division is set according to the largest first weighted ratio value of each partition, and is reduced by the second set ratio when the first constraint condition is not satisfied. The first constraint condition includes: The volume value of a single volume division is less than or equal to the product of the second preset ratio value and the volume of the first partition and the product of the second preset ratio value and the volume of the second partition; The volume of the first partition after volume division is greater than or equal to the product of the third preset ratio value and the volume of the first partition; The water level change value of the first partition after volume division and the water level change value of the second partition after volume division are less than or equal to the fifth threshold.

10. The method according to claim 1, characterized in that The calculating the daily chemical substance concentration of each adjusted partition includes: Judge whether the actual value of the daily lake effluent flow is zero: If not, calculate the chemical substance concentration of each adjusted partition according to the one-dimensional advection-diffusion equation; If so, calculate the chemical substance concentration of each adjusted partition according to the one-dimensional advection-diffusion equation combined with the return flow attenuation based on the set number of return times.

Images