A site grading design method based on mixed integer programming considering site drainage

Abstract

The application provides a photovoltaic power station optimal field design method suitable for a region with field drainage requirements; the method takes the minimum excavation and filling amount of the whole field as an objective function, takes the field slope ratio, excavation and filling balance and field drainage requirements as constraint conditions, integrates a local coefficient matrix and a vector into a global coefficient matrix and a vector, can obtain an optimal field square scheme by solving a mixed integer programming optimal solution, and obtains an optimal solution of the drainage field design problem.

Description

Technical Field

[0001] This invention relates to a site leveling design method based on mixed integer programming that takes into account site drainage, which is mainly applicable to the site leveling design of photovoltaic power plants in areas with site drainage requirements. Background Technology

[0002] Photovoltaic modules are typically arranged in a north-south orientation, with the panels able to rotate east-west to face the sun. To minimize mutual shading between adjacent east-west facing panels, the east-west slope of the photovoltaic site needs to be controlled within a certain range (generally not exceeding 2%). In desert areas, photovoltaic panels are often cleaned using small cleaning robots that crawl along the north-south direction of the panels. The maximum crawling slope of the cleaning robot determines the maximum north-south slope of the site (generally not exceeding 6%).

[0003] Large-scale photovoltaic power plants often occupy tens of square kilometers. In areas with significant topographical variations, such as deserts, the earthwork volume for photovoltaic site leveling can reach millions or even tens of millions of cubic meters, resulting in a massive engineering project. Traditional site leveling design methods require first dividing the photovoltaic site into blocks, then using site leveling software to design the leveling for each block, and finally considering the coordination of boundaries between blocks. Traditional methods are extremely time-consuming and labor-intensive, limiting the number of blocks to a maximum of a dozen or so. Furthermore, it's difficult to effectively connect the blocks and determine the merits of different leveling schemes; multiple schemes must be compared before a relatively optimal one is selected. Traditional site leveling methods are labor-intensive, inefficient, and involve large amounts of excavation and filling.

[0004] Furthermore, in non-desert areas, photovoltaic site leveling in general areas needs to consider site drainage conditions, meaning that the leveling plan should not contain any depressions, and rainfall from any point should be able to flow by gravity to the edge of the site. Since photovoltaic sites generally occupy a large area, traditional block leveling methods cannot simultaneously minimize excavation and filling while ensuring adequate site drainage. Summary of the Invention

[0005] The purpose of this invention is to address the problems existing in the background technology by providing a site leveling design method based on mixed integer programming that considers site drainage. This method divides the site into continuous blocks, uses the minimum cut-fill volume for the entire site as the objective function, and takes the site slope ratio, cut-fill balance (or a given difference), and site drainage requirements as constraints. It integrates the global coefficient matrix and vectors using local coefficient matrices and vectors, ultimately obtaining the optimal site leveling scheme by finding the optimal solution through mixed integer programming.

[0006] The technical solution adopted in this invention is:

[0007] A site leveling design method based on mixed integer programming that considers site drainage, characterized by the following steps:

[0008] S1. Raw terrain data processing: Extract site coordinate data and generate terrain interpolation function F;

[0009] S2. Site Block Division and Coding: Based on the photovoltaic matrix layout, the entire site is divided into blocks, and all nodes, edges, and blocks are coded to form a topological structure. Let n be the number of blocks, nodes, internal nodes, and edges respectively. el n nd n ndin n ed ;

[0010] S3. Divide the area into rectangular grids: Divide each block into rectangular grids and encode them, recording the total number of grids as n. ms ;

[0011] S4. Initialize the global initial coefficient matrix and vector, including:

[0012] The initial global coefficient vector f,

[0013] The global inequality constraint coefficient matrix A can be represented as a block matrix:

[0014]

[0015] The constant vector b constrained by the global inequality can be represented as a block vector:

[0016]

[0017] Global equality constraint coefficient vector A eq This can be represented as a block vector:

[0018] A eq =(A eq11 A eq12 )

[0019] Global equality constraint constant b eq ;

[0020] S5. Correct the global initial coefficient vector f;

[0021] S6. Generate the local equality constraint coefficient vector A in the order of grid numbers. eq11-i And integrate it into the global equality constraint coefficient vector submatrix A eq11 middle;

[0022] S7. Generate the local inequality constraint coefficient matrix A for each edge according to the edge numbering order. 11-i Local inequality constraints on constant vector b 11 -i; and then integrate them sequentially into the global matrix A. 11 And vector b 11 middle;

[0023] S8. Generate the local inequality constraint coefficient matrix A according to the grid numbering order. 21-i And integrated into the global inequality constraint coefficient submatrix A 21 middle;

[0024] S9. Correcting the global inequality constraint coefficient submatrix A 22 Global inequality constraint constant vector b 21 ;

[0025] S10. Correct the global matrix A sequentially according to the internal node numbering order. 31 ;

[0026] S11. Correction matrix A 33 A 41 A 43 A 44 A 54 ;

[0027] S12. Correction matrix b 21 b 41 b 51 ;

[0028] S13. Correct the global equality constraint constant b eq ;

[0029] S14. Solving the optimal solution problem of mixed integer programming:

[0030] min f T H

[0031]

[0032] S15. The optimal solution H of the mixed integer programming obtained according to S14 forms the site leveling design surface. This can be formed by connecting all adjacent nodes with the design elevations of each node obtained from the mixed integer programming solution to form the final site leveling design surface.

[0033] Furthermore, the global first-order coefficient vector f is (n nd +n ms +5n ndin A 1-dimensional vector.

[0034] Furthermore, the global inequality constraint coefficient matrix A is (2n ed +2n ms +9n ndin )×(n nd +n ms +5n ndin ( ) order matrix.

[0035] Furthermore, the constant vector b of the global inequality constraint is (2n ed+2n ms +9n ndin A 1-dimensional vector.

[0036] Furthermore, the global equality constraint coefficient vector A eq 1×(n nd +n ms +5n ndin )-dimensional vector.

[0037] In a second aspect of the present invention, the present invention provides a non-transitory computer-readable storage medium having a computer program stored thereon, characterized in that, when executed by a processor, the computer program implements the steps of the above-described site leveling design method based on mixed integer programming considering site drainage.

[0038] In a third aspect of the present invention, the present invention provides an electronic device including a memory, a processor, and a computer program stored in the memory and executable on the processor, characterized in that the processor, when executing the program, implements the steps of the above-described site leveling design method based on mixed integer programming considering site drainage.

[0039] The beneficial effects of this invention are:

[0040] 1) Optimal solution: The optimal field leveling scheme can be obtained by finding the optimal solution of the field leveling problem through mixed integer programming. The optimal solution of the drainage field leveling problem is theoretically obtained.

[0041] 2) Good applicability: This method can be used to design drainage leveling schemes for different photovoltaic site leveling projects.

[0042] 3) Fast design speed: By programming the design method, the site leveling scheme design can be completed in a very short time. Attached Figure Description

[0043] Figure 1 This is a schematic diagram of the original site area division and numbering according to an embodiment of the present invention;

[0044] Figure 2 A schematic diagram of grid division for a certain area;

[0045] Figure 3 A site leveling diagram and original site layout diagram are provided.

[0046] Figure 4 Design drawings for a photovoltaic power station using a traditional partial planar site leveling design scheme;

[0047] Figure 5 for Figure 4 The design drawings of a certain photovoltaic power station adopting the site leveling design scheme of the present invention.

[0048] In the attached diagram: 100 - Original site; 200 - Block division line; 300 - Coordinate system; 400 - North arrow; 500 - Grid division line; 600 - Final design site plan. Detailed Implementation

[0049] This embodiment provides a site leveling design method based on mixed integer programming that considers site drainage, including the following steps:

[0050] S1 - Generates the terrain interpolation function F.

[0051] Using the original terrain elevation point data (x, y, z) as samples, a two-dimensional terrain interpolation function F is generated. F is used to calculate the ground elevation value z0 = F(x0, y0) for a given plane coordinate (x0, y0).

[0052] S2 - Divide the area into rectangular blocks, number the blocks, nodes, and edges, and denote the number of blocks, nodes, internal nodes, and edges as n. el n nd n ndin n ed .

[0053] The photovoltaic site is divided into several rectangular blocks using straight lines parallel to the east-west and north-south directions. The size of each block is determined based on the arrangement of the photovoltaic matrix and the dimensions of the photovoltaic modules. Nodes, edges, and blocks are numbered, and the coordinates of nodes, the node composition of edges and blocks, and other information are recorded to form a topological structure.

[0054] like Figure 1 This example illustrates how the site is divided into 15 blocks using 4 east-west straight lines and 6 north-south straight lines. In the numbering, 1, 2, ..., 15 are block numbers; ①, ②, ..., ㉔ are node numbers; ⑥, ⑦, ..., ⑱, ⑲ are internal node numbers; and (1), (2), ..., (38) are edge numbers. In this example, n el =15, n nd =24, n ndin =8, n ed =38.

[0055] S3 - Divide the grid into rectangular sections, and let the total number of sections be n. ms .

[0056] Each block is divided into several grids using straight lines parallel to the east-west and north-south directions, with each grid having a similar size and area.

[0057] For all blocks, such as Figure 2 The horizontal plane is divided into 4 rows (y-axis) and 5 columns (x-axis) of straight lines running east-west and north-south, totaling 20 grids. The center coordinates of each grid are denoted as (x...). j y j), the elevation value z of each grid center i The area of ​​each grid cell is calculated using the interpolation function F. j ,j=1,2…20. Total number of grids n ms =300, assign a uniform number to all grid cells:

[0058] S4 - Set the global initial coefficient matrix and vector, including:

[0059] (1) Global initial coefficient vector f, f is (n nd +n ms +5n ndin f is a zero vector of order 364×1;

[0060] (2) The global inequality constraint coefficient matrix A, A is (2n ed +2n ms +9n ndin )×(n nd +n ms +5n ndin A is a zero matrix of order 748×364. In this example, A is a zero matrix of order 748×364. Representing matrix A in block matrix form is as follows:

[0061]

[0062] Among them: A 11 For 2n ed ×n nd A zero matrix of order 76×24, which is a zero matrix of order 76×24 in this example;

[0063] A 12 For 2n ed ×n ms A zero matrix of order 76×300 in this example;

[0064] A 13 For 2n ed ×n ndin A zero matrix of order 76×8 in this example;

[0065] A 14 For 2n ed ×4n ndin A zero matrix of order 76×32 in this example;

[0066] A 21 For 2n ms ×n nd A zero matrix of order 600×24, which is a zero matrix of order 600×24 in this example;

[0067] A 22 For 2n ms×n ms A zero matrix of order 600×300, which is a zero matrix of order 600×300 in this example;

[0068] A 23 For 2n ms ×n ndin A zero matrix of order 600×8, which is a zero matrix of order 600×8 in this example;

[0069] A 24 For 2n ms ×4n ndin A zero matrix of order 600×32, which is a zero matrix of order 600×32 in this example;

[0070] A 31 A 41 Both are 4n ndin ×n nd The zero matrices are of order 32×24 in this example;

[0071] A 32 A 42 Both are 4n ndin ×n ms The zero matrices are of order 32×300 in this example;

[0072] A 33 A 43 Both are 4n ndin ×n ndin The zero matrices are of order 32×8 in this example;

[0073] A 34 A 44 Both are 4n ndin ×4n ndin A zero matrix of order 32×32, which is a zero matrix of order 32×32 in this example;

[0074] A 51 For n ndin ×n nd A zero matrix of order 8×24, which is an 8×24 zero matrix in this example;

[0075] A 52 For n ndin ×n ms A zero matrix of order 8×300, which is an 8×300 zero matrix in this example;

[0076] A 53 For n ndin ×n ndin The zero matrix of order 8×8, which is an 8×8 zero matrix in this example;

[0077] A 54 For n ndin ×4n ndin The zero matrix of order 8×32 is used in this example.

[0078] (3) Global inequality constraints constrain constant vector b, where b is (2n ed +2n ms +9n ndin b is a 748×1 order zero vector. In this example, b is a 748×1 order zero vector. Vector b can be expressed in block vector form as follows:

[0079]

[0080] Where: b 11 For 2n ed The zero vector of order ×1, in this example it is the zero vector of order 76×1;

[0081] b 21 For 2n ms The zero vector of order ×1, in this example it is the zero vector of order 600×1;

[0082] b 31 b 41 Both are 4n ndin ×1 order zero vector, in this example b 31 b 41 All are zero vectors of order 32×1;

[0083] b 51 For n ndin The zero vector is of order ×1, and in this example it is the zero vector of order 8×1.

[0084] (4) Global equality constraint coefficient vector A eq A eq 1×(n nd +n ms +5n ndin A )-dimensional zero vector, in this example, A eq Given a zero vector of order 1×364, let vector A... eq Represented as a block vector:

[0085] A eq =(A eq11 A eq12 )

[0086] Among them: A eq11 For 1×n nd The zero vector of order A, in this example eq11 It is a zero vector of order 1×24;

[0087] A eq12 1×(n ms +5n ndin A is a zero vector of order 1, in this example, A. eq12 It is a zero vector of order 1×340.

[0088] (5) Global equality constraint constant b eq b eq It is zero.

[0089] S5 - Correcting the global initial coefficient vector f

[0090] For the global initial coefficient vector f, let its nth... nd +1~n nd +n ms The terms are the areas of each grid, s1 to s2. nms .

[0091] S6 - Generate local equality constraint coefficient vector A in the order of grid numbers. eq11-i And integrate it into the global equality constraint coefficient vector submatrix A eq11 In the middle. Perform the following operations on all grids in sequence until completion.

[0092] like Figure 2 As shown, the four nodes of the block containing the grid numbered i are a, b, c, and d, and their corresponding coordinates are: (x... a ,y a ), (x b ,y b ), (x c ,y c ), (x d ,y d The corresponding ground elevation is: z a z b z c z d Calculated using the interpolation function F. The area of ​​grid i is s. i The center coordinates are (x i y i The ground elevation value z at the center point i Calculated using the interpolation function F.

[0093] Calculate the number of times A is used for each grid step. i B i C i D i

[0094]

[0095]

[0096]

[0097]

[0098] Local equality constraint coefficient vector A eq11-i for:

[0099] A eq11-i =s i [A i B i C i D i ]

[0100] The equality constraint coefficient vector A eq11-i Integrating into the global equality constraint coefficient vector submatrix A eq11 Chinese: Update A eq11 Add A to the values ​​of terms a, b, c, and d in the equation. eq11-i The values ​​of items 1 through 4 are used as new values.

[0101] S7 - Generate the local inequality constraint coefficient matrix A for each edge according to the edge numbering order. 11-i Local inequality constraints on constant vector b 11-i And then sequentially integrate them into the global matrix A. 11 And vector b 11 In the middle, perform the following operations on all edges in sequence until completion.

[0102] like Figure 2 As shown, the operation is illustrated using edges numbered p and q in different directions within a certain block. Edge p is a north-south edge with nodes numbered a and b at both ends; edge q is an east-west edge with nodes numbered a and c at both ends. The local inequality constraint coefficient matrix A for edge p is shown. 11-p and local inequality constraints on constant vector b 11-p They are respectively:

[0103]

[0104] b 11-p =(y b -y a )×[k sn k ns ] T

[0105] Where, k sn k is the maximum permissible slope ratio of the site from south to north. ns k represents the maximum permissible slope ratio of the site from north to south. ns k sn ≥0.

[0106] q-edge local inequality constraint coefficient matrix A 11-q and local inequality constraints on constant vector b 11-q They are respectively:

[0107]

[0108] b 11-q=(x c -x a )×[k we k ew ] T

[0109] Where, k we k is the maximum permissible slope ratio of the site from west to east. ew k is the maximum permissible slope ratio of the site from east to west. ew k we ≥0.

[0110] Let the local inequality constraint coefficient matrix A of edge p be... 11-p Local inequality constraints on constant vector b 11-p Integrate into global matrix A 11 And vector b 11 In the middle: Update matrix A 11 Change the values ​​corresponding to columns a and b in rows 2p-1 and 2p to A respectively. 11-p The values ​​corresponding to columns 1 and 2 in rows 1 and 2; update b 11 Change the values ​​of terms 2p-1 and 2p to b respectively. 11-p The first and second values. The same operation is used to change the local inequality constraint coefficient matrix A of edge q. 11-q Local inequality constraints on constant vector b 11-q Integrate into global matrix A 11 And vector b 11 middle.

[0111] S8 - Generate the local inequality constraint coefficient matrix A in the order of mesh numbers. 21-i And integrated into the global inequality constraint coefficient submatrix A 21 middle.

[0112] For example, in step S6, the first-order coefficient A of the grid numbered i is calculated. i B i C i D i Local inequality constraint coefficient matrix A 21-i for:

[0113]

[0114] Let the local inequality constraint coefficient matrix A 21-i Integrating into the global inequality constraint coefficient submatrix A 21 In the middle: Update matrix A 21 Change the values ​​corresponding to columns a, b, c, and d in rows 2i-1 and 2i to A respectively. 21-i The values ​​corresponding to columns 1 to 4 in rows 1 to 2.

[0115] S9 - Corrected Global Inequality Constraint Coefficient Submatrix A 22 Global inequality constraint constant vector b 21 .

[0116] For the global inequality constraint coefficient submatrix A 22 Let the value in the i-th column of the 2i-1 and 2i-th rows be -1, where i = 1, 2, ..., n ms .

[0117] For the constant vector b constrained by the global inequality 21 Let its 2i-1th term (odd-numbered term) be z. i Its 2ith term (even-numbered term) is -z i , z i Let i be the ground elevation of the center point of the i-th grid, i = 1, 2, ..., n ms .

[0118] S10 - Correct the global matrix A sequentially according to the internal node numbering order. 31 .

[0119] like Figure 1 As shown, the operation is explained using the first internal node ⑥. The nodes adjacent to internal node ⑥ are ②, ⑤, ⑦, and ⑩. Let matrix A... 31 The value in column ⑥ (the number of the i-th internal node) of rows 1 to 4 (4i-1 to 4i, where i is the sequence number of the internal node) is 1, A 31 The values ​​in columns ②, ⑤, ⑦, and ⑩ (the node numbers adjacent to the i-th internal node) of rows 1 to 4 (4i-1 to 4i, where i is the sequence number of the internal node) are -1. Repeat the above operation in the order of the internal node numbers until all are completed.

[0120] S11 - Correction Matrix A 33 A 41 A 43 A 44 A 54 .

[0121] For A 33 Let the value of the i-th column in the 4i-1 to 4i rows be 1, i = 1, 2, ..., n. ms ;

[0122] For A 41 Let A 41 =-A 31 ;

[0123] For A 43 Let A 43 =-A 33 ;

[0124] For A44 Let A 44 =ME; where: M is a constant, taken as M = [max(zi) - min(zi)] × 10; E is 4n ndin ×4n ndin An identity matrix of order 1.

[0125] For A 54 ,make

[0126] S12-correction matrix b 21 b 41 b 51 .

[0127] For b 21 Let its 2i-1th term (odd-numbered term) be z. i Its 2ith term (even-numbered term) is -z i , z i Let i be the ground elevation of the center point of the i-th grid, i = 1, 2, ..., n ms .

[0128] For b 41 Let each of its values ​​be M, and the value of M is the same as in step S11.

[0129] For b 51 Let its values ​​be -1.

[0130] S13 - Correction of global equality constraint constant b eq

[0131] make

[0132]

[0133] Where dV is the difference between the fill and the cut, dV = fill volume - cut volume.

[0134] S14 - Solving the optimal solution problem of mixed integer programming

[0135] After calculating the coefficient matrices, the site leveling problem with drainage requirements is transformed into a mixed-integer programming optimization problem:

[0136] min f T H

[0137]

[0138] Find the optimal solution H for the above mixed integer programming problem, where H is n. nd +n ms +5n ndin A dimensional vector, whose first n... nd The term represents the design elevation h of each node in the optimal site leveling scheme.

[0139] S15 - Forms the site plan design surface.

[0140] Connecting the design elevation h of each adjacent node forms the final site leveling design surface, such as... Figure 3 The number 600 in the middle refers to.

[0141] The application of this invention in a photovoltaic power station project illustrates the comparison between the effects of this invention and previous methods.

[0142] A photovoltaic power station covers an area of ​​1.5 square kilometers with an installed photovoltaic capacity of 52MW. The site leveling of a certain area within the project must meet drainage design requirements, meaning that rainwater from any location on the site can flow off the site. During the bidding phase, the project adopted a traditional partial planar approach for site leveling design, dividing the entire site into four zones, each requiring a planar site leveling design, such as... Figure 4 As shown, the lighter-colored areas represent excavation, and the darker-colored areas represent fill. The total excavation volume for this project is 69,000 m³. 3 The total fill volume is 67,000 m³. 3 The design took approximately 6 hours.

[0143] If the site leveling design is carried out using the method of this invention, the site is divided into 776 zones, such as... Figure 5 As shown, the lighter-colored areas represent excavation, and the darker-colored areas represent fill. The total excavation volume for this project is 42,000 m³. 3 The total fill volume was 41,000 m³. 3 The design process takes approximately 0.5 hours. Compared to traditional methods, this invention saves approximately 39% of the excavation and filling work and reduces design time by approximately 92%.

[0144] This method can be implemented as a software module for data processing. The invention can be implemented using software plus the necessary general-purpose hardware platform. For example, it can be implemented using existing processors, or by a dedicated processor used for this or other purposes in a suitable system, or by a hardwired system. Embodiments of the invention also include non-transitory computer-readable storage media, which include machine-readable media for carrying or having machine-executable instructions or data structures stored thereon; such machine-readable media can be any available medium accessible by a general-purpose or special-purpose computer or other machine with a processor. For example, such machine-readable media can include RAM, ROM, EPROM, EEPROM, CD-ROM or other optical disk storage, magnetic disk storage or other magnetic storage devices, or any other medium that can be used to carry or store the required program code in the form of machine-executable instructions or data structures and can be accessed by a general-purpose or special-purpose computer or other machine with a processor. When information is transmitted or provided to a machine via a network or other communication connection (hardwired, wireless, or a combination of hardwired and wireless), that connection is also considered a machine-readable medium.

Claims

1. A site leveling design method based on mixed integer programming considering site drainage, characterized in that, Includes the following steps: S1. Raw terrain data processing: Extract site coordinate data and generate terrain interpolation function F; S2. Site Block Division and Coding: Based on the photovoltaic matrix layout, the entire site is divided into blocks, and all nodes, edges, and blocks are coded to form a topology; let n be the number of blocks, nodes, internal nodes, and edges respectively. el n nd n ndin n ed ; S3. Divide the area into rectangular grids: Divide each block into rectangular grids and encode them, recording the total number of grids as n. ms ; S4. Initialize the global initial coefficient matrix and vector, including: The initial global coefficient vector f, The global inequality constraint coefficient matrix A can be represented as a block matrix: The constant vector b constrained by the global inequality can be represented as a block vector: Global equality constraint coefficient vector A eq This can be represented as a block vector: A eq =(A eq11 A eq12 ) Global equality constraint constant b eq ; S5. Correct the global initial coefficient vector f; S6. Generate the local equality constraint coefficient vector A in the order of grid numbers. eq11-i And integrate it into the global equality constraint coefficient vector submatrix A eq11 middle; S7. Generate the local inequality constraint coefficient matrix A for each edge according to the edge numbering order. 11-i Local inequality constraints on constant vector b 11-i And then sequentially integrate them into the global matrix A. 11 And vector b 11 middle; S8. Generate the local inequality constraint coefficient matrix A according to the grid numbering order. 21-i And integrated into the global inequality constraint coefficient submatrix A 21 middle; S9. Correcting the global inequality constraint coefficient submatrix A 22 Global inequality constraint constant vector b 21 ; S10. Correct the global matrix A sequentially according to the internal node numbering order. 31 ; S11. Correction matrix A 33 A 41 A 43 A 44 A 54 ; S12. Correction matrix b 21 b 41 b 51 ; S13. Correct the global equality constraint constant b eq ; S14. Solving the optimal solution problem of mixed integer programming: S15. The optimal solution H of the mixed integer programming obtained according to S14 is used to form the field level design surface.

2. The site leveling design method based on mixed integer programming considering site drainage according to claim 1, characterized in that: The global first-order coefficient vector f is (n nd +n ms +5n ndin A 1-dimensional vector.

3. The site leveling design method based on mixed integer programming considering site drainage according to claim 1, characterized in that: The global inequality constraint coefficient matrix A is (2n ed +2n ms +9n ndin )×(n nd +n ms +5n ndin ( ) order matrix.

4. The site leveling design method based on mixed integer programming considering site drainage according to claim 1, characterized in that: The constant vector b of the global inequality constraint is (2n) ed +2n ms +9n ndin A 1-dimensional vector.

5. The site leveling design method based on mixed integer programming considering site drainage according to claim 1, characterized in that: Global equality constraint coefficient vector A eq 1×(n nd +n ms +5n ndin )-dimensional vector.

6. A non-transitory computer-readable storage medium having a computer program stored thereon, characterized in that, When executed by a processor, the computer program implements the steps of the site leveling design method based on mixed integer programming considering site drainage as described in any one of claims 1 to 5.

7. An electronic device comprising a memory, a processor, and a computer program stored in the memory and executable on the processor, characterized in that, When the processor executes the program, it implements the steps of the site leveling design method based on mixed integer programming considering site drainage as described in any one of claims 1 to 5.

Images