Underground water pollution migration simulation method and device, electronic equipment and storage medium

Abstract

The invention belongs to the technical field of underground water numerical simulation, and provides an underground water pollution migration simulation method and device, electronic equipment and a storage medium, and the method comprises the steps: constructing a three-dimensional octree grid, a control equation and a definite solution condition of a hydrogeological structure of a target region; carrying out discrete solution on the control equation in the three-dimensional octree grid by adopting a finite volume method taking a node as a center; constructing a post-processing approximation function and a local posterior error indicator on octree leaf subunits of the three-dimensional octree grid; determining a unit set contributed by a global error by adopting a Drfler marking strategy, then performing splitting or merging operation on leaf units, and performing grid updating by adopting a high-order interpolation or volume weighted average method; and according to a preset simulation termination condition, repeatedly executing underground water pollution migration simulation and finally performing visual display. According to the method, the simulation accuracy of water pollution migration is improved in a posterior error estimation mode.

Description

Groundwater pollution migration simulation methods, devices, electronic equipment and storage media Technical Field

[0001] This invention relates to the field of groundwater numerical simulation technology, and in particular to a method, apparatus, electronic device, and storage medium for simulating groundwater pollution migration. Background Technology

[0002] Accurate and efficient prediction of pollutant migration patterns in groundwater systems is crucial for environmental protection, water resource management, and contaminated site remediation decisions. Numerical simulation is a core technology for achieving this goal, with mainstream methods including the finite difference method, the finite element method, and the finite volume method. Among these, the finite volume method is widely used in groundwater pollution migration simulation due to its ability to satisfy local and global mass conservation and its flexibility in handling complex unstructured meshes.

[0003] Traditional finite volume methods often employ fixed grids for spatial discretization. However, real-world groundwater pollution problems typically exhibit highly localized characteristics: concentration gradients are steep in areas such as pollutant concentration fronts, pollution source regions, and abrupt permeability interfaces, while most other areas show relatively gentle changes. Using a globally uniform, fine grid, while ensuring accuracy, leads to a dramatic increase in the number of grid cells, resulting in an exponential rise in computational costs and memory requirements. Conversely, using a coarse grid makes it difficult to capture fronts and local high gradient phenomena, producing significant numerical dispersion or oscillations, severely impacting the reliability of the simulation results. Summary of the Invention

[0004] In order to at least solve one of the technical problems existing in the prior art, the present invention provides a groundwater pollution migration simulation method, apparatus, electronic device and storage medium.

[0005] One aspect of the present invention provides a method for simulating the migration of groundwater pollution, comprising:

[0006] Acquire hydrogeological data of the target area, construct a three-dimensional octree mesh of the hydrogeological structure of the target area based on the hydrogeological data, and construct the governing equations and boundary conditions for groundwater seepage and pollutant migration.

[0007] The node-centric finite volume method is used to discretize and solve the governing equations in a three-dimensional octree mesh to obtain the head field and pollutant concentration field at the current time step.

[0008] Based on the head field and pollutant concentration field at the current time step, a post-processing approximation function is constructed in the octagonal leaf element of the three-dimensional octagonal grid. The residual flux is calculated through the post-processing approximation function to obtain the local posterior error indicator, which includes the residual within the element, the interface incompatibility error, the convection term error, and the upwind pattern error.

[0009] The global error estimate is determined based on the local posterior error indicator, and the Dörfler labeling strategy is used to determine the set of cells with the largest global error contribution and the set of cells with the smallest global error contribution based on the global error estimate.

[0010] The octree leaf cells of the set of cells with the largest global error contribution are split, the set of cells with the smallest global error contribution is merged, and the mesh is updated using a high-order interpolation or volume-weighted averaging method to obtain a new three-dimensional octree mesh.

[0011] Based on the preset simulation termination conditions, the groundwater pollution migration simulation is repeated to obtain the head field constant and pollutant concentration field data for all time steps under the preset simulation termination conditions, and then visualized using three-dimensional visualization.

[0012] According to the groundwater pollution migration simulation method, hydrogeological data of the target area are acquired; based on the hydrogeological data, a three-dimensional octree mesh of the hydrogeological structure of the target area is constructed; and the governing equations and boundary value conditions for groundwater seepage and pollutant migration are constructed, including:

[0013] Based on the hydrogeological data of the target area, a three-dimensional geological model including the distribution of aquifers and impermeable layers is constructed.

[0014] The three-dimensional geological model is initially discretized using a hexahedral mesh to generate a regular spatial partition covering the entire computational domain;

[0015] The regular spatial partitioning is processed using an adaptive grid data structure based on an octree to obtain a three-dimensional octree grid;

[0016] Based on hydrogeological data, Darcy's law and the law of conservation of mass are used to determine the governing equations for groundwater seepage, and the convection-dispersion equation is used as the governing equations for pollutant migration. The boundary conditions include using a constant head boundary as the boundary condition type and using the spatial distribution of the simulated head field and pollutant concentration field at the start time as the initial conditions.

[0017] According to the groundwater pollution migration simulation method described above, a node-centered finite volume method is used to discretize and solve the governing equations in a three-dimensional octree mesh, obtaining the head field and pollutant concentration field at the current time step, including:

[0018] The computational domain of the 3D octree mesh is divided into multiple non-overlapping control units;

[0019] By integrating the conservation equation for each control unit and converting the volume integral into the surface integral using the Gaussian divergence theorem, the flux relationship of physical quantities between adjacent control units is obtained.

[0020] Based on the physical flux relationship between adjacent control units, a three-dimensional octree grid is constructed to form a linear equation system for all control units. The linear equation system is then discretized using the finite volume method to obtain the head field and pollutant concentration field at the current time step.

[0021] According to the groundwater pollution migration simulation method, a post-processing approximation function is constructed in the octagonal leaf cells of a three-dimensional octagonal grid based on the head field and pollutant concentration field at the current time step. The residual flux is then calculated using the post-processing approximation function to obtain a local posterior error indicator, including:

[0022] For the leaf cells of the 3D octree mesh, a piecewise quadratic polynomial approximation solution is constructed based on the head field and pollutant concentration field at the current time step. Local error indicators are then calculated based on this piecewise quadratic polynomial approximation solution. for:

[0023]

[0024] in, The residuals within a cell are used to measure the degree to which the numerical solution violates the governing equations. This is the interface non-coordination error, used to quantify the discontinuity of the approximate solution at the interface between adjacent elements; This is the convection term error, used to quantify the error caused by the discontinuities in the solution under the influence of the flow field; The upwind scheme error is used to quantify the additional numerical dispersion introduced by the upwind scheme used to ensure numerical stability.

[0025] According to the groundwater pollution migration simulation method described above, the global error estimate is determined based on the local posterior error indicator, and the Dörfler labeling strategy is used to determine the set of cells with the largest global error contribution and the set of cells with the smallest global error contribution based on the global error estimate, including:

[0026] According to the local error indicator Calculate the global error estimator for:

[0027]

[0028] in, Represents the set of all leaf units. Leaf unit;

[0029] According to refined parameters Through constraint formula Find the set of cells with the largest global error contribution from the set of all leaf cells. ;

[0030] Based on coarsening parameters Through constraint formula Find the set of cells with the minimum global error contribution from the set of all leaf cells. .

[0031] According to the groundwater pollution migration simulation method, the octree leaf cells of the cell set with the largest global error contribution are split, and the cell set with the smallest global error contribution is merged to obtain a new three-dimensional octree mesh. The method also includes:

[0032] Obtain the octagonal leaf element and perform a splitting operation to obtain the sub-mesh element. The initial values ​​of the physical quantities of the sub-mesh element are obtained using a quadratic polynomial interpolation method. Based on the numerical solution information of the parent element and its neighboring elements, local smooth reconstruction is performed through the physical quantity spatial distribution function to obtain the initial values ​​of the physical quantities of the sub-element with spatial gradient.

[0033] Before performing the merging operation on the set of units that contributes the least to the global error, it will be confirmed that all child nodes of the same parent node satisfy the coarsening parameter. When merging child units, the volume-weighted average of the physical quantities of all child units will be used as the physical quantity of the new parent node after merging.

[0034] The three-dimensional octree mesh obtained through splitting and merging operations is subjected to balance constraints so that the octree level difference between adjacent leaf cells is no greater than 1, resulting in a new three-dimensional octree mesh.

[0035] According to the groundwater pollution migration simulation method, the groundwater pollution migration simulation is repeatedly executed according to a preset simulation termination condition to obtain the head field constant and pollutant concentration field data for all time steps under the preset simulation termination condition. These data are then visualized using three-dimensional visualization, including:

[0036] Iterative discrete solution, error and mesh adjustment are performed. The three-dimensional octree mesh generated in each iteration is used as the physical field of the previous iteration as the initial condition for the next iteration, until the preset simulation termination condition is reached to obtain the output data. The output data includes the head field constant, pollutant concentration field data, spatial coordinates and node level of all time steps.

[0037] The output data is visualized using a 3D graphics rendering engine, and the output data is stored in a structured manner.

[0038] Another aspect of the present invention provides a groundwater pollution migration simulation device, comprising:

[0039] The first module is used to acquire hydrogeological data of the target area, construct a three-dimensional octree mesh of the hydrogeological structure of the target area based on the hydrogeological data, and construct the control equations and boundary conditions for groundwater seepage and pollutant migration.

[0040] The second module is used to discretize and solve the governing equations in a three-dimensional octree mesh using the node-centered finite volume method to obtain the head field and pollutant concentration field at the current time step.

[0041] The third module is used to construct a post-processing approximation function in the octagonal leaf cells of the three-dimensional octagonal grid based on the head field and pollutant concentration field of the current time step. The residual flux is calculated through the post-processing approximation function to obtain the local posterior error indicator, which includes the cell internal residual, interface incompatibility error, convection term error and upwind pattern error.

[0042] The fourth module is used to determine the global error estimate based on the local posterior error indicator, and to determine the set of units with the largest global error contribution and the set of units with the smallest global error contribution using the Dörfler labeling strategy based on the global error estimate.

[0043] The fifth module is used to perform splitting operations on the octree leaf cells of the cell set with the largest global error contribution, merge operations on the cell set with the smallest global error contribution, and update the mesh using high-order interpolation or volume-weighted averaging methods to obtain a new three-dimensional octree mesh.

[0044] The sixth module is used to repeatedly execute the groundwater pollution migration simulation according to the preset simulation termination conditions, obtain the head field constant and pollutant concentration field data of all time steps under the preset simulation termination conditions, and then use three-dimensional visualization for visualization display.

[0045] Another aspect of the present invention provides an electronic device, including a processor and a memory;

[0046] The memory is used to store programs;

[0047] The processor executes the program to implement the method as described above.

[0048] This invention also discloses a computer program product or computer program, which includes computer instructions stored in a computer-readable storage medium. A processor of a computer device can read the computer instructions from the computer-readable storage medium and execute the computer instructions, causing the computer device to perform the methods described above.

[0049] The beneficial effects of this invention are as follows: Based on a comprehensive posterior error indicator including internal residuals, interface compatibility errors, convection errors, and upwind pattern errors, the local refinement and merging of the octree mesh are dynamically driven by the error contribution ratio. This automatically densifies the mesh in high-gradient regions such as pollutant migration fronts and source regions to improve resolution and simulation accuracy, and reasonably coarsens the mesh in regions with gentle physical fields to reduce computational scale, thus achieving self-optimization of computational resources. Through integrated 3D visualization, a highly efficient and interpretable system solution is provided for the problem of groundwater pollution migration. Attached Figure Description

[0050] Figure 1 is a schematic diagram of the groundwater pollution migration simulation process according to an embodiment of the present invention.

[0051] Figure 2 is a schematic diagram of the target area in an embodiment of the present invention.

[0052] Figure 3 is a schematic diagram of the adaptive mesh octree data structure according to an embodiment of the present invention, wherein (a) is the initial hexahedral mesh and the corresponding root node (level 0) of the octree structure; (b) is the eight second-level hexahedral meshes generated after the initial mesh is split once, corresponding to the eight second-level nodes generated in the octree; (c) is one of the second-level hexahedral meshes further split into eight third-level hexahedral meshes, corresponding to the eight third-level nodes generated in the octree; (d) is one of the third-level hexahedral meshes further split into eight fourth-level hexahedral meshes, corresponding to the eight fourth-level nodes generated in the octree.

[0053] Figure 4 is a three-dimensional schematic diagram of the initial data of the target area in an embodiment of the present invention, wherein (a) is the initial seepage head distribution map and (b) is the initial pollutant concentration field distribution map.

[0054] Figure 5 is a comparison of error distribution during the simulation process of an embodiment of the present invention.

[0055] Figure 6 shows the simulation results of the three-dimensional pollutant concentration distribution in an embodiment of the present invention, wherein (a) is the three-dimensional distribution of pollutant concentration for 90 days of simulation, (b) is the three-dimensional distribution of pollutant concentration for 180 days of simulation, and (c) is the three-dimensional distribution of pollutant concentration for 360 days of simulation.

[0056] Figure 7 is a partial schematic diagram of the concentration distribution and corresponding adaptive grid of a certain elevation horizontal section according to an embodiment of the present invention. (a) shows the concentration distribution and local grid state of the section at 90 days of simulation, (b) shows the concentration distribution and local grid state of the section at 180 days of simulation, and (c) shows the concentration distribution and local grid state of the section at 360 days of simulation.

[0057] Figure 8 is a schematic diagram of the groundwater pollution migration simulation device according to an embodiment of the present invention. Detailed Implementation

[0058] The embodiments of the present invention are described in detail below, examples of which are shown in the accompanying drawings. Throughout the description, the same or similar reference numerals denote the same or similar elements or elements having the same or similar functions. In the following description, suffixes such as "module," "part," or "unit" used to denote elements are used only for the purpose of illustrative purposes and have no specific meaning in themselves. Therefore, "module," "part," or "unit" can be used interchangeably. Terms such as "first," "second," etc., are used only to distinguish technical features and should not be construed as indicating or implying relative importance, or implicitly indicating the number of indicated technical features, or implicitly indicating the sequential relationship of the indicated technical features. In the following description, the consecutive reference numerals for method steps are for ease of review and understanding. Adjusting the implementation order of steps, in conjunction with the overall technical solution of the present invention and the logical relationship between the various steps, will not affect the technical effect achieved by the technical solution of the present invention. The embodiments described below with reference to the accompanying drawings are exemplary and are only used to explain the present invention, and should not be construed as limiting the present invention.

[0059] Referring to Figure 1, which is a schematic diagram of the groundwater pollution migration simulation process according to an embodiment of the present invention, including but not limited to steps S100-S600:

[0060] S100: Obtain hydrogeological data of the target area; based on the hydrogeological data, construct a three-dimensional octree mesh of the hydrogeological structure of the target area; and construct the governing equations and boundary conditions for groundwater seepage and pollutant migration.

[0061] In some embodiments, a three-dimensional geological model including the distribution of aquifers and impermeable layers is constructed based on hydrogeological data of the target area; the three-dimensional geological model is initially discretized using a hexahedral mesh to generate a regular spatial partition covering the entire computational domain; the regular spatial partition is processed using an octree-based adaptive mesh data structure to obtain a three-dimensional octree mesh; based on the hydrogeological data, Darcy's law and the law of conservation of mass are used to determine the governing equations for groundwater seepage, and the convection-dispersion equation is used as the governing equations for pollutant migration; the boundary conditions include using a constant head boundary as the boundary condition type, and using the spatial distribution of the simulated head field and pollutant concentration field at the start time as the initial conditions.

[0062] Specifically, as shown in Figure 2, the target area diagram illustrates the spatial layout and pollution diffusion path of a chromium-contaminated site. The factory buildings, slag heaps, and surrounding soil within the site are contaminated with Cr(VI). Hexavalent chromium ions enter the aquifer through processes such as adsorption and desorption, and rainwater leaching, forming a groundwater pollution source that migrates and spreads with the groundwater flow.

[0063] Based on the embodiment shown in Figure 2, the embodiments of the present invention simulate the target area and the types of geological bodies it encompasses, such as aquifers, impermeable layers, and other strata related to groundwater flow. Subsequently, based on constraint information such as borehole data and geological profiles, a three-dimensional model reflecting the actual geological structure is constructed using geological body modeling methods, and a corresponding geological attribute scalar field is formed through inversion, interpolation, and other means.

[0064] After completing the 3D geological modeling, a hexahedral mesh is used to spatially discretize the model for numerical calculations. The mesh size and density can be set according to the simulation accuracy requirements and the complexity of the geological structure.

[0065] To achieve adaptive dynamic adjustment of the mesh, this embodiment of the invention uses an octree data structure as the basis for mesh management. Specifically, referring to the adaptive mesh octree data structure diagram shown in Figure 3, it includes: (a) an initial hexahedral mesh and its corresponding octree root node (level 0); (b) eight second-order hexahedral meshes generated after the initial mesh undergoes one split, corresponding to eight second-order nodes in the octree; (c) one of the second-order hexahedral meshes is further split into eight third-order hexahedral meshes, corresponding to eight third-order nodes in the octree; and (d) one of the third-order hexahedral meshes is further split into eight fourth-order hexahedral meshes, corresponding to eight fourth-order nodes in the octree. For this purpose, an octree adaptive mesh class is constructed, which encapsulates data members such as 3D coordinates, physical properties, simulation parameters, parent and child node pointers, node center coordinates, and edge lengths, and provides functional methods such as node splitting, merging, and neighbor querying. Next, using the center coordinates of the initial hexahedral mesh as input, each initial mesh is used as the root node of the octree, and the initial octree mesh structure is constructed layer by layer, thereby providing support for subsequent dynamic mesh refinement and coarsening based on simulation errors.

[0066] In this embodiment of the invention, establishing a mathematical model describing the processes of groundwater seepage and pollutant migration is the core step in conducting numerical simulations.

[0067] First, a seepage equation is established based on Darcy's law and the law of conservation of mass. For isotropic and homogeneous saturated aquifers, steady-state or unsteady-state groundwater flow equations can be used:

[0068] Steady-state seepage equation: ,in Let be the hydraulic conductivity tensor. For water head;

[0069] Unsteady seepage equation: ,in For water storage rate, For time, For source and sink items.

[0070] Secondly, pollutant transport equations describe the diffusion, convection, and other processes of pollutants in groundwater. For example, a commonly used solute transport equation is expressed as:

[0071]

[0072] in, Porosity As a delay factor, For pollutant concentration, Let be the dispersion coefficient tensor. It is the groundwater flow velocity (according to Darcy's law) (Calculated) The volumetric flow rate per unit volume of fluid source / sink. The concentration of the source and sink terms.

[0073] While establishing the mathematical model, reasonable boundary conditions and initial conditions need to be set to ensure the uniqueness of the solution. The boundary conditions supported by this invention include:

[0074] Define the head boundary and specify the boundary head value.

[0075] A fixed flow boundary is defined, specifying the flow rate per unit area at the boundary. Initial conditions are defined as the simulation start time. The distribution of water head and pollutant concentration in the entire region was determined by obtaining measured data from the target area.

[0076] S200 employs a node-centric finite volume method to discretize and solve the governing equations in a three-dimensional octree mesh, obtaining the head field and pollutant concentration field at the current time step.

[0077] In some embodiments, the computational domain of the three-dimensional octree mesh is divided into multiple non-overlapping control units; the conservation equation is integrated for each control unit, and the volume integral is transformed into the surface integral using the Gaussian divergence theorem to obtain the physical flux relationship between adjacent control units; a linear equation system of all control units of the three-dimensional octree mesh is constructed based on the physical flux relationship between adjacent control units, and the linear equation system is discretized using the finite volume method to obtain the head field and pollutant concentration field at the current time step.

[0078] For example, based on the initial octree mesh and boundary conditions, the core of this step is to use the finite volume method to numerically discretize and solve the governing equations in order to obtain the seepage field and concentration field at the current time step.

[0079] Spatial discretization is performed using a node-centric finite volume method. This embodiment of the invention divides the computational domain into non-overlapping control units. The conservation equations are integrated over each control unit, and the Gaussian divergence theorem is used to transform the volume integral into a surface integral on the unit interface, thereby establishing flux relationships between adjacent units and naturally ensuring local mass conservation. For time discretization, an implicit scheme is employed to improve computational stability.

[0080] Specifically, the groundwater seepage control equations are discretized and solved to obtain the hydraulic head distribution. and the Darcy velocity field calculated therefrom .

[0081] In this embodiment, for a node Centered on, with a volume of The hexahedral control unit will incorporate the seepage equation within the control volume and at the time step. Inner integral:

[0082]

[0083] Implicit difference is used for the time term, and the Gaussian divergence theorem is applied to the diffusion term to transform it into a summation of fluxes over the interfaces of each element. After discretization and rearrangement, the following information about the nodes is obtained. and its adjacent nodes Linear algebraic equations:

[0084]

[0085] in, and For the conduction coefficient between adjacent units Coefficients related to geometric dimensions and water storage coefficient. For including source and sink items And the constant term of the information from the previous time step. The conduction coefficient between adjacent cells. Calculated from the harmonic mean, for example, for the eastward interface:

[0086]

[0087] In the formula, , They are nodes Its eastern neighbor node The permeability coefficient, The boundary area, The distance between the centers of the two nodes is denoted as .

[0088] Assemble such equations for all grid cells, thus forming a large system of linear equations. Solving this system of equations will yield the head values ​​for all nodes at the new time step (n+1). Furthermore, according to Darcy's Law Calculate the Darcy velocity at each unit interface. This provides input for the pollutant migration equation.

[0089] The velocity field obtained from seepage solution Substitute the pollutant migration equation and use the finite volume method for discrete solution.

[0090] Similarly, for nodes Integral pollutant migration equation centered on the control unit:

[0091]

[0092] During the discrete process:

[0093] diffusion term Discretization is performed using a central difference scheme to ensure accuracy.

[0094] Convection term To maintain numerical stability and reduce non-physical oscillations, an upwind discretization scheme is used. That is, discretization is performed through the element interface. Convection current depends on the direction of flow velocity:

[0095]

[0096] in, Let be the Darcy velocity at the interface. For the interface area, The outward normal unit vector, and These represent the concentration values ​​for the upstream and downstream units, respectively.

[0097] After discretization and rearrangement, we obtain information about the nodes. Algebraic equation for concentration:

[0098]

[0099] Where the coefficient , This includes contributions from diffusion, convection, and reaction terms. Including source and sink terms and concentration information from the previous time step, all unit equations are assembled to form a system of linear equations. Solving this gives the pollutant concentration field at the new time step. .

[0100] S300, based on the head field and pollutant concentration field of the current time step, constructs a post-processing approximation function in the octagonal leaf element of the three-dimensional octagonal grid, and calculates the residual flux through the post-processing approximation function to obtain the local posterior error indicator, which includes the internal residual of the element, the interface incompatibility error, the convection term error and the upwind pattern error.

[0101] In some embodiments, piecewise quadratic polynomial approximations are constructed for the leaf cells of the three-dimensional octree mesh based on the head field and pollutant concentration field at the current time step, and local error indicators are calculated based on the piecewise quadratic polynomial approximations. for:

[0102]

[0103] in, The residuals within a cell are used to measure the degree to which the numerical solution violates the governing equations. This is the interface non-coordination error, used to quantify the discontinuity of the approximate solution at the interface between adjacent elements; This is the convection term error, used to quantify the error caused by the discontinuities in the solution under the influence of the flow field; The upwind scheme error is used to quantify the additional numerical dispersion introduced by the upwind scheme used to ensure numerical stability.

[0104] For example, based on the numerical solution of step S200, a post-processing approximation function with higher regularity is constructed, and a local posterior error indicator is calculated based on this function, thereby providing an error basis for mesh adaptive adjustment.

[0105] In each octagonal leaf node unit Above, a higher-order post-processing approximation solution is constructed based on the current numerical solution. In this embodiment of the invention, Using a piecewise quadratic polynomial:

[0106]

[0107] coefficient The establishment of [the system] follows two basic principles:

[0108] The physical flux at each cell interface remains consistent with the discrete flux calculated by the original finite volume method; the average value within the cell or the value at its center point remains unchanged from the value calculated by the original finite volume method.

[0109] This embodiment transforms the piecewise constant solution, which cannot be directly calculated for gradients, into an approximate function that is smooth within the cell and provides continuous gradient information.

[0110] Based on post-processing approximate solution Residual analysis of the original governing equations and numerical schemes, for each leaf element. Construct a comprehensive local posterior error indicator .

[0111] S400 determines the global error estimate based on the local posterior error indicator, and uses the Dörfler labeling strategy based on the global error estimate to determine the set of cells with the largest global error contribution and the set of cells with the smallest global error contribution.

[0112] In some embodiments, based on local error indicators Calculate the global error estimator for:

[0113]

[0114] in, Represents the set of all leaf units. For leaf units, where the global error estimator It characterizes the overall discrete error level of the current numerical solution;

[0115] According to refined parameters Through constraint formula Find the set of cells with the largest global error contribution from the set of all leaf cells. In some embodiments, this embodiment sets a refinement ratio parameter. The goal is to find a set containing the minimum number of units such that the sum of the error contributions of the units in that set is at least 0.7 times the global error.

[0116] In some embodiments, all leaf units are sorted according to their local error indicators. The values ​​are sorted in descending order, starting with the unit with the largest error, and the squared error values ​​are accumulated sequentially. When the cumulative error contribution first meets or exceeds the threshold When the accumulation stops, all the accumulated units constitute the set to be refined. And mark it as needing further refinement.

[0117] Based on coarsening parameters Through constraint formula Find the set of cells with the minimum global error contribution from the set of all leaf cells. .

[0118] In some embodiments, the coarsening ratio parameter The goal is to find a set of units that contributes the least to the error, such that the sum of their errors does not exceed 0.1 times the global error.

[0119] For example, all leaf units are sorted according to their local error indicators. The values ​​are sorted in ascending order, starting from the unit with the smallest error, and the squared error values ​​of each unit are accumulated sequentially. When the cumulative error contribution is about to exceed the threshold When the accumulation stops, all the cells that have been accumulated constitute the set to be coarsened. And mark it as to be coarsened.

[0120] For example, referring to Figure 5, the left column of Figure 5 shows the error distribution results of using an initial static uniform grid at simulation times of 90 days, 180 days, and 360 days, and the right column shows the error distribution results of using the adaptive grid of the method of the present invention at the same time. It can be determined that the adaptive grid of the present invention has smaller error.

[0121] S500 performs a split operation on the octree leaf cells of the cell set with the largest global error contribution, performs a merge operation on the cell set with the smallest global error contribution, and updates the mesh using a high-order interpolation or volume-weighted averaging method to obtain a new three-dimensional octree mesh.

[0122] In some embodiments, octree leaf elements are split to obtain sub-mesh elements. The initial physical quantities of the sub-mesh elements are obtained using a quadratic polynomial interpolation method. Based on the numerical solution information of the parent element and its neighboring elements, a local smooth reconstruction is performed using the physical quantity spatial distribution function to obtain the initial physical quantity values ​​of the sub-elements with spatial gradients. For example, for each leaf node marked as to be refined, an octree splitting operation is performed to uniformly split it into eight sub-elements of equal geometric size, which become internal nodes. Using the numerical solution information of the parent element and its neighboring elements, a locally smooth physical quantity spatial distribution function is reconstructed, and then initial values ​​with spatial gradients are assigned to the center or integration point of each sub-element.

[0123] Before performing a merging operation on the set of elements that contributes the least to the global error, it is confirmed that all child nodes of the same parent node satisfy the coarsening parameter. When merging child elements, the volume-weighted average of the physical quantities of all child elements is used as the physical quantity of the new parent node after merging. The 3D octree mesh after splitting and merging operations is subjected to balance constraints to ensure that the octree level difference between adjacent leaf elements is no greater than 1, resulting in a new 3D octree mesh. For example, the balance constraint automatically enforces a 2:1 balance constraint, which requires that the octree level difference between any two adjacent leaf elements in the final mesh is no greater than 1, ensuring a smooth transition of mesh resolution. This is achieved by recursively checking and performing necessary temporary refinement on coarser elements that do not meet the constraints.

[0124] S600, based on the preset simulation termination conditions, repeatedly executes the groundwater pollution migration simulation to obtain the head field constant and pollutant concentration field data for all time steps under the preset simulation termination conditions, and then uses three-dimensional visualization for visualization display.

[0125] In some embodiments, discrete solutions, errors, and mesh adjustments are performed iteratively. The three-dimensional octree mesh generated in each iteration is used as the initial condition of the physical field of the previous iteration for the next iteration, until the preset simulation termination condition is reached to obtain output data. The output data includes the head field constant, pollutant concentration field data, spatial coordinates, and node levels for all time steps. The output data is visualized using a three-dimensional graphics rendering engine and stored in a structured manner.

[0126] It should be noted that the solution-error estimation-grid adjustment process is a closed-loop iterative process. The physical field obtained in the previous step through interpolation or averaging is used as the initial condition, and steps S200 to S500 are repeated. This iterative loop continues until the total simulation duration is met, at which point the system outputs the final simulation results. The output data is stored in a structured format, covering the spatial coordinates, node levels, head values, and pollutant concentration values ​​of all octagonal leaf nodes, providing a complete data foundation for subsequent analysis and visualization.

[0127] Furthermore, to present the simulation results intuitively, the output data is imported into a 3D graphics rendering engine for visualization. Color mapping is performed based on pollutant concentration values. Referring to the 3D pollutant concentration distribution simulation results shown in Figure 6, this embodiment of the invention uses volume rendering technology to 3D display the spatial distribution and migration path of pollutants, or uses slice rendering to present the concentration distribution of a specified cross-section. The results are shown in Figure 7, where (a) is the concentration distribution and local mesh state of the cross-section at 90 days of the simulation, (b) is the concentration distribution and local mesh state of the cross-section at 180 days of the simulation, and (c) is the concentration distribution and local mesh state of the cross-section at 360 days of the simulation. Specifically, this embodiment of the invention can simultaneously overlay and display the hierarchical structure or dynamic change process of the octree mesh, thereby achieving synchronous, dynamic, and visual display of the pollutant migration process and the adaptive evolution state of the mesh, providing intuitive technical support for groundwater pollution migration analysis.

[0128] Figure 8 is a schematic diagram of a groundwater pollution migration simulation device according to an embodiment of the present invention. The device includes a first module 810, a second module 820, a third module 830, a fourth module 840, a fifth module 850, and a sixth module 860.

[0129] The system comprises three modules: The first module acquires hydrogeological data of the target area, constructs a three-dimensional octree mesh of the hydrogeological structure based on this data, and establishes the governing equations and boundary conditions for groundwater seepage and pollutant migration. The second module uses a node-centric finite volume method to discretize and solve the governing equations within the three-dimensional octree mesh, obtaining the head field and pollutant concentration field at the current time step. The third module constructs a post-processing approximation function in the octagonal leaf cells of the three-dimensional octree mesh based on the head field and pollutant concentration field at the current time step. It then calculates the residual flux using the post-processing approximation function to obtain a local posterior error indicator, which includes the residual within the cell, interface incompatibility error, convection term error, and windward error. The fourth module is used to determine the global error estimate based on the local posterior error indicator, and to determine the set of cells with the largest global error contribution and the set of cells with the smallest global error contribution using the Dörfler labeling strategy based on the global error estimate; the fifth module is used to perform a splitting operation on the octree leaf cells of the set of cells with the largest global error contribution, a merging operation on the set of cells with the smallest global error contribution, and to update the mesh using a higher-order interpolation or volume-weighted averaging method to obtain a new three-dimensional octree mesh; the sixth module is used to repeatedly execute the groundwater pollution migration simulation according to the preset simulation termination conditions to obtain the head field constant and pollutant concentration field data for all time steps of the preset simulation termination conditions, and then use three-dimensional visualization for visualization.

[0130] Exemplarily, with the cooperation of the first, second, third, fourth, fifth, and sixth modules in the device, the embodiment device can implement any of the aforementioned groundwater pollution migration simulation methods, namely, acquiring hydrogeological data of the target area, constructing a three-dimensional octree mesh of the hydrogeological structure of the target area based on the hydrogeological data, and constructing the governing equations and boundary conditions for groundwater seepage and pollutant migration; using the node-centered finite volume method, discretizing and solving the governing equations in the three-dimensional octree mesh to obtain the head field and pollutant concentration field at the current time step; constructing a post-processing approximation function in the octree leaf cells of the three-dimensional octree mesh based on the head field and pollutant concentration field at the current time step; and calculating the residual flux through the post-processing approximation function to obtain the local posterior error indicator, wherein... The local posterior error indicators include intra-cell residuals, interface incoordination errors, convection term errors, and upwind pattern errors. Global error estimates are determined based on these indicators. The Dörfler labeling strategy is then used to identify the set of cells contributing the most and least to the global error. The octree leaves of the set contributing the most to the global error are split, while the set contributing the least is merged. A higher-order interpolation or volume-weighted averaging method is used to update the mesh, resulting in a new 3D octree mesh. The groundwater pollution migration simulation is repeated based on preset termination conditions to obtain head constants and pollutant concentration field data for all time steps under these conditions. These data are then visualized using 3D visualization. The beneficial effects of this invention are as follows: Based on a comprehensive posterior error indicator including internal residuals, interface compatibility errors, convection errors, and upwind pattern errors, the local refinement and merging of the octree mesh are dynamically driven by the error contribution ratio. This automatically densifies the mesh in high-gradient regions such as pollutant migration fronts and source regions to improve resolution and simulation accuracy, and reasonably coarsens the mesh in regions with gentle physical fields to reduce computational scale, thus achieving self-optimization of computational resources. Through integrated 3D visualization, a highly efficient and interpretable system solution is provided for the problem of groundwater pollution migration.

[0131] This invention also provides an electronic device, which includes a processor and a memory;

[0132] The memory stores the program;

[0133] The processor executes a program to perform the aforementioned groundwater pollution migration simulation method; the electronic device has the function of carrying and running the software system for groundwater pollution migration simulation provided in the embodiments of the present invention, such as a personal computer, minicomputer, mainframe, workstation, network or distributed computing environment, standalone or integrated computer platform, or communicating with charged particle tools or other imaging devices, etc.

[0134] This invention also provides a computer-readable storage medium storing a program that is executed by a processor to implement the groundwater pollution migration simulation method described above.

[0135] In some alternative embodiments, the functions / operations mentioned in the block diagrams may not occur in the order shown in the operation diagrams. For example, depending on the functions / operations involved, two consecutively shown blocks may actually be executed substantially simultaneously, or the blocks may sometimes be executed in reverse order. Furthermore, the embodiments presented and described in the flowcharts of this invention are provided by way of example to provide a more comprehensive understanding of the technology. The disclosed methods are not limited to the operations and logic flows presented in the embodiments of this invention. Alternative embodiments are contemplated, in which the order of various operations is changed and sub-operations described as part of a larger operation are executed independently.

[0136] This invention also discloses a computer program product or computer program, which includes computer instructions stored in a computer-readable storage medium. A processor of a computer device can read the computer instructions from the computer-readable storage medium and execute the computer instructions, causing the computer device to perform the aforementioned groundwater pollution migration simulation method.

[0137] Furthermore, although the invention has been described in the context of functional modules, it should be understood that, unless otherwise stated, one or more of the described functions and / or features may be integrated into a single physical device and / or software module, or one or more functions and / or features may be implemented in a separate physical device or software module. It is also understood that a detailed discussion of the actual implementation of each module is unnecessary for understanding the invention. Rather, considering the properties, functions, and internal relationships of the various functional modules in the apparatus disclosed in the embodiments of the invention, the actual implementation of the module will be understood within the scope of conventional skill of an engineer. Therefore, those skilled in the art can implement the invention as set forth in the claims using ordinary techniques without excessive experimentation. It is also understood that the specific concepts disclosed are merely illustrative and are not intended to limit the scope of the invention, which is determined by the full scope of the appended claims and their equivalents.

[0138] If the aforementioned functions are implemented as software functional units and sold or used as independent products, they can be stored in a computer-readable storage medium. Based on this understanding, the technical solution of this invention, essentially, or the part that contributes to the prior art, or a portion of the technical solution, can be embodied in the form of a software product. This computer software product is stored in a storage medium and includes several instructions to cause a computer device (which may be a personal computer, server, or network device, etc.) to execute all or part of the steps of the methods described in the various embodiments of this invention. The aforementioned storage medium includes various media capable of storing program code, such as USB flash drives, portable hard drives, read-only memory (ROM), random access memory (RAM), magnetic disks, or optical disks.

[0139] The logic and / or steps represented in the flowchart or otherwise described herein, for example, can be considered as a sequenced list of executable instructions for implementing logical functions, and can be embodied in any computer-readable medium for use by, or in conjunction with, an instruction execution system, apparatus, or device (such as a computer-based system, a processor-included system, or other system that can fetch and execute instructions from, an instruction execution system, apparatus, or device). For the purposes of this specification, "computer-readable medium" can be any means that can include, store, communicate, propagate, or transmit programs for use by, or in conjunction with, an instruction execution system, apparatus, or device.

[0140] More specific examples of computer-readable media (a non-exhaustive list) include: electrical connections (electronic devices) having one or more wires, portable computer disk drives (magnetic devices), random access memory (RAM), read-only memory (ROM), erasable and editable read-only memory (EPROM or flash memory), fiber optic devices, and portable optical disc read-only memory (CDROM). Furthermore, computer-readable media can even be paper or other suitable media on which the program can be printed, because the program can be obtained electronically, for example, by optically scanning the paper or other medium, followed by editing, interpreting, or otherwise processing as necessary, and then stored in computer memory.

[0141] It should be understood that various parts of the present invention can be implemented in hardware, software, firmware, or a combination thereof. In the above embodiments, multiple steps or methods can be implemented in software or firmware stored in memory and executed by a suitable instruction execution system. For example, if implemented in hardware, as in another embodiment, it can be implemented using any one or a combination of the following techniques known in the art: discrete logic circuits having logic gates for implementing logical functions on data signals, application-specific integrated circuits (ASICs) having suitable combinational logic gates, programmable gate arrays (PGAs), field-programmable gate arrays (FPGAs), etc.

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

[0143] Although embodiments of the invention have been shown and described, those skilled in the art will understand that various changes, modifications, substitutions and alterations can be made to these embodiments without departing from the principles and spirit of the invention, the scope of which is defined by the claims and their equivalents.

[0144] The above is a detailed description of the preferred embodiments of the present invention, but the present invention is not limited to the embodiments described. Those skilled in the art can make various equivalent modifications or substitutions without departing from the spirit of the present invention, and these equivalent modifications or substitutions are all included within the scope defined by the claims of this application.

Claims

1. A method for simulating the migration of groundwater pollution, characterized in that, include: Hydrogeological data of the target area is acquired. Based on this data, a three-dimensional octree mesh of the hydrogeological structure of the target area is constructed, along with the governing equations and boundary conditions for groundwater seepage and pollutant migration. The governing equations are discretized and solved using the node-centric finite volume method within the three-dimensional octree mesh to obtain the head field and pollutant concentration field at the current time step. Based on these fields, a post-processing approximation function is constructed for each octree leaf cell of the three-dimensional octree mesh. The residual flux is then calculated using this approximation function to obtain a local posterior error indicator, which includes the cell internal residual, interface incompatibility error, convection term error, and upwind pattern error. The global error estimate is determined based on the local posterior error indicator. Finally, the Dörfler labeling strategy is used to identify the set of cells with the largest and smallest global error contributions based on the global error estimate. The octree leaf cells of the set of cells with the largest global error contribution are split, the set of cells with the smallest global error contribution is merged, and the mesh is updated using a high-order interpolation or volume-weighted averaging method to obtain a new three-dimensional octree mesh. Based on the preset simulation termination conditions, the groundwater pollution migration simulation is repeated to obtain the head field constant and pollutant concentration field data for all time steps under the preset simulation termination conditions, and then visualized using three-dimensional visualization.

2. The groundwater pollution migration simulation method according to claim 1, characterized in that, The process of acquiring hydrogeological data of the target area, constructing a three-dimensional octree mesh of the hydrogeological structure of the target area based on the hydrogeological data, and constructing the governing equations and boundary conditions for groundwater seepage and pollutant migration includes: constructing a three-dimensional geological model including the distribution of aquifers and impermeable layers based on the hydrogeological data of the target area; initially discretizing the three-dimensional geological model using a hexahedral mesh to generate a regular spatial partition covering the entire computational domain; processing the regular spatial partition using an adaptive mesh data structure based on an octree to obtain a three-dimensional octree mesh; determining the governing equations for groundwater seepage using Darcy's law and the law of conservation of mass based on the hydrogeological data, and using the convection-dispersion equation as the governing equation for pollutant migration; the boundary conditions include using a constant head boundary as the boundary condition type, and using the spatial distribution of the simulated head field and pollutant concentration field at the initial time as the initial conditions.

3. The groundwater pollution migration simulation method according to claim 1, characterized in that, The method employs a node-centric finite volume approach to discretize and solve the governing equations within a three-dimensional octree mesh, yielding the head field and pollutant concentration field for the current time step. This includes: dividing the computational domain of the three-dimensional octree mesh into multiple non-overlapping control units; integrating the conservation equations for each control unit and converting the volume integral to the surface integral using the Gaussian divergence theorem to obtain the physical flux relationships between adjacent control units; constructing a system of linear equations for all control units within the three-dimensional octree mesh based on these flux relationships; and then discretizing and solving the linear equations using the finite volume approach to obtain the head field and pollutant concentration field for the current time step.

4. The groundwater pollution migration simulation method according to claim 1, characterized in that, The process of constructing a post-processing approximation function in the leaf cells of a 3D octree mesh based on the head field and pollutant concentration field at the current time step, and calculating the residual flux using the post-processing approximation function to obtain a local posterior error indicator, includes: constructing piecewise quadratic polynomial approximation solutions for the leaf cells of the 3D octree mesh based on the head field and pollutant concentration field at the current time step, and calculating the local error indicator based on the piecewise quadratic polynomial approximation solutions. for: in, The residuals within a cell are used to measure the degree to which the numerical solution violates the governing equations. This is the interface non-coordination error, used to quantify the discontinuity of the approximate solution at the interface between adjacent elements; This is the convection term error, used to quantify the error caused by the discontinuities in the solution under the influence of the flow field; The upwind scheme error is used to quantify the additional numerical dispersion introduced by the upwind scheme used to ensure numerical stability.

5. The groundwater pollution migration simulation method according to claim 4, characterized in that, The step of determining the global error estimate based on local posterior error indicators, and then using the Dörfler labeling strategy to determine the set of units with the largest global error contribution and the set of units with the smallest global error contribution based on the global error estimate, includes: determining the global error estimate based on local error indicators. Calculate the global error estimator for: in, Represents the set of all leaf units. Leaf element; based on refinement parameters Through constraint formula Find the set of cells with the largest global error contribution from the set of all leaf cells. According to the coarsening parameters Through constraint formula Find the set of cells with the minimum global error contribution from the set of all leaf cells. 。 6. The groundwater pollution migration simulation method according to claim 5, characterized in that, The process of splitting the octagonal leaf cells of the set of cells contributing the most to the global error and merging the set of cells contributing the least to the global error to obtain a new 3D octagonal mesh further includes: obtaining sub-mesh cells by splitting the octagonal leaf cells; obtaining the initial physical quantity values ​​of the sub-mesh cells using a quadratic polynomial interpolation method; performing local smooth reconstruction using a physical quantity spatial distribution function based on the numerical solution information of the parent cell and its neighboring cells to obtain the initial physical quantity values ​​of the sub-cells with spatial gradients; confirming that all child nodes of the same parent node satisfy the coarsening parameter before merging the sub-cells, and using the volume-weighted average of the physical quantities of all sub-cells as the physical quantity of the new parent node after merging when merging sub-cells; and applying balance constraints to the 3D octagonal mesh after the splitting and merging operations to ensure that the octagonal level difference between adjacent leaf cells is no greater than 1 to obtain a new 3D octagonal mesh.

7. The groundwater pollution migration simulation method according to claim 1, characterized in that, The process involves repeatedly executing a groundwater pollution migration simulation based on a preset simulation termination condition to obtain the head field constant and pollutant concentration field data for all time steps under the preset simulation termination condition. This data is then visualized using 3D visualization. The process includes: iteratively executing discrete solutions, adjusting errors and meshes; using the physical quantity field from the previous iteration as the initial condition for each iteration's generated 3D octree mesh, and performing the next iteration until the preset simulation termination condition is met, thus obtaining output data. The output data includes the head field constant, pollutant concentration field data, spatial coordinates, and node levels for all time steps; visualizing the output data using a 3D graphics rendering engine; and storing the output data in a structured manner.

8. A groundwater pollution migration simulation device, characterized in that, include: The first module is used to acquire hydrogeological data of the target area, construct a three-dimensional octree mesh of the hydrogeological structure of the target area based on the hydrogeological data, and construct the control equations and boundary conditions for groundwater seepage and pollutant migration. The second module is used to discretize and solve the governing equations in a three-dimensional octree mesh using the node-centric finite volume method to obtain the head field and pollutant concentration field at the current time step. The third module is used to construct a post-processing approximation function in the octree leaf cells of the three-dimensional octree mesh based on the head field and pollutant concentration field at the current time step, and calculate the residual flux through the post-processing approximation function to obtain the local posterior error indicator, which includes the cell internal residual, interface incompatibility error, convection term error, and upwind pattern error. The fourth module is used to determine the global error estimate based on the local posterior error indicator, and to determine the set of cells with the largest global error contribution and the set of cells with the smallest global error contribution based on the global error estimate using the Dörfler labeling strategy. The fifth module is used to perform splitting operations on the octree leaf cells of the cell set with the largest global error contribution, merge operations on the cell set with the smallest global error contribution, and update the mesh using high-order interpolation or volume-weighted averaging methods to obtain a new three-dimensional octree mesh. The sixth module is used to repeatedly execute the groundwater pollution migration simulation according to the preset simulation termination conditions, obtain the head field constant and pollutant concentration field data of all time steps under the preset simulation termination conditions, and then use three-dimensional visualization for visualization display.

9. An electronic device, characterized in that, It includes a processor and a memory; the memory is used to store a program; the processor executes the program to implement the groundwater pollution migration simulation method as described in any one of claims 1-7.

10. A computer-readable storage medium, characterized in that, The storage medium stores a program that is executed by a processor to implement the groundwater pollution migration simulation method as described in any one of claims 1-7.

Images