Water network superstructure optimization method based on molecular group algorithm

Abstract

The invention discloses a water network superstructure optimization method based on a molecular swarm algorithm. The method mainly comprises two steps of establishing a multi-impurity water network superstructure optimization model and solving the optimization model based on an MGA. The multi-impurity water network is subjected to synchronous design by considering water reuse, regeneration reuse and regeneration cycle, and a water network synchronous design structure comprising various possible schemes is established. After the superstructure is analyzed, an optimization model of water network synchronous design is established. According to the model, the annual total cost of the water network is used as a target function, and material conservation in the water network, environmental limitation of impurity concentration in wastewater discharge and valuable substance recovery requirements are used as constraints. And then an optimal value is solved by adopting an MGA. The method is characterized in that for a multi-dimensional and multi-constraint optimization solution problem, the algorithm can better converge to a globally optimal solution, and the solution difficulty of the problem is reduced.

Description

Technical Field

[0001] The invention mainly involves the concepts of water resource circulation and reuse, multi-impurity water network superstructure modeling, intelligent optimization algorithms and other fields, and specifically involves an algorithm for solving multi-impurity water network superstructure optimization problems. Technical Background

[0002] Since the 1980s, the optimal allocation of water resources has gradually become a research hotspot in the industrial field. With the rapid development of industrialization, the demand for water resources has continued to rise, especially in water-intensive industries such as petrochemicals, metallurgy, papermaking and food processing. The importance of water resource management has become increasingly prominent. The large amount of wastewater containing impurities generated in the industrial production process poses a serious threat to the environment. At the same time, the excessive consumption of water resources has also aggravated the problem of water scarcity and pollution. Therefore, how to achieve efficient use and conservation of water resources while improving production efficiency through reasonable water resource management and wastewater treatment has become a key issue in the current industrial process optimization.

[0003] Modeling and solving water networks is a major challenge in realizing water-saving technologies in process industries. The core issue is to establish a suitable model for water network propositions and provide an effective algorithm. Water network design is divided into two methods: distributed design and synchronous design. Distributed design is to optimize the water use system and wastewater treatment system separately, but this method ignores the potential synergistic effects in water reuse and impurity treatment, resulting in waste of water resources and increased treatment costs. Especially in multi-impurity systems, wastewater often contains many different types of impurities. The interactions and treatment requirements between these impurities increase the complexity of optimization, making it difficult to obtain the global optimal solution. Synchronous design is to integrate the water use system and wastewater treatment system into a comprehensive network, namely the water network superstructure, by constructing an optimization framework covering water use, wastewater treatment and impurity recovery. This method can not only optimize water use, but also rationally plan the separation, treatment and recovery of impurities, thereby achieving the dual goals of energy conservation and emission reduction and economic benefits.

[0004] For optimization problems with constraints, the penalty function method is usually used when solving them with intelligent optimization algorithms. This method transforms the constrained optimization problem into an unconstrained optimization problem by incorporating the constraints into the objective function, thereby simplifying the problem-solving process. The basic idea of ​​the penalty function method is to introduce one or more penalty terms that penalize solutions that violate the constraints, so that the algorithm tends to find solutions that satisfy the constraints. However, for a multi-dimensional and multi-constrained optimization problem such as the water network superstructure, it is difficult for general intelligent optimization algorithms to converge to the global optimal solution even when combined with the penalty function method.

[0005] Glossary:

[0006] Molecular swarm algorithm: Based on the improvement of the gravitational search algorithm, during the evolution of the gravitational search algorithm, the material population will be updated according to the optimal position of the population. If a material finds a current optimal position, other materials will quickly gather to this optimal position. If this optimal position is a local optimal point, the algorithm cannot re-search in the solution space, resulting in stagnation or premature convergence. In order to solve this problem, a critical distance r is set. cd , when the distance between two individuals in the population is greater than r cd At the same time, set r cd is a quantity that decreases adaptively during the iteration process. In this way, at the beginning of the iteration, the individuals in the population can avoid premature aggregation and fall into the local optimum when searching for the area with the best solution through the interaction of attraction and repulsion. As the iteration proceeds, r cd Gradually decreases, the attraction between individuals increases, the repulsion decreases, and the algorithm converges faster. cd Approaching zero, it is equivalent to the gravitational search algorithm, and the population moves toward the optimal solution. Summary of the invention

[0007] In order to solve the above-mentioned modeling defects, the present invention discloses a water network superstructure optimization method based on molecular swarm algorithm.

[0008] To achieve the above purpose, the technical solution of the present invention is as follows:

[0009] A water network superstructure optimization method based on molecular swarm algorithm comprises the following steps:

[0010] Step 1: Obtain the original data of water network superstructure;

[0011] Step 2: establishing a multi-impurity water network superstructure model of the water network superstructure according to the original data;

[0012] Step three: using a molecular swarm algorithm to solve the multi-impurity water network superstructure model until a termination condition is reached to obtain an optimal multi-impurity water network superstructure.

[0013] As a further improvement, the termination condition is a preset number of iterations or convergence.

[0014] As a further improvement, the raw data of the water network superstructure includes:

[0015] Annual operation time H, fresh water cost coefficient C f , annual investment factor AR of the processing unit, investment cost coefficient IC of the processing unit t , Operation cost coefficient OC of processing unit t , cost function index α, 0<α≤1, cost coefficient C of waste treatment sw, the cost coefficient of the medicine C med 、Recycling price of valuable materials C rec , the mass load of impurity k in process i, Δm i，k，out , the residual rate β of impurity k in process i i，k , the maximum inlet concentration of impurity k in process i Minimum inlet concentration of impurity k in process i Concentration limit of impurity K discharged into the environment and the water flow requirement WR of water unit i i .

[0016] As a further improvement, the objective function of the multi-impurity water network superstructure model is as follows:

[0017]

[0018] in, is the fresh water cost, f i is the amount of fresh water entering the system, PU is the water use unit; is the annual investment cost of the processing unit, F i is the total water flow entering process i, TU is a processing unit; is the annual operating cost of the treatment unit, is the waste treatment cost of the treatment unit, W i is the wastewater discharge of process i; To process the unit drug costs, Economic benefits generated by the recovery of valuable substances, is the flow rate of impurity k entering process i.

[0019] For further improvement, the constraints of the objective function are as follows:

[0020] 1. Total water balance:

[0021] For water units:

[0022]

[0023] X ij is the amount of wastewater reused from process j to process i, X ji is the amount of wastewater reused from process i to process j;

[0024] For processing units:

[0025]

[0026] 2. Impurity mass balance:

[0027] For water units:

[0028]

[0029] C i,k,out is the outlet concentration of impurity k in process i, C j,k,out is the outlet concentration of impurity k in process j, Δm i,k,tot is the mass load of impurity k in process i;

[0030] For processing units:

[0031]

[0032] 3. Maximum import concentration limit:

[0033]

[0034] 4. Minimum import concentration limit:

[0035]

[0036] C i,k,in is the inlet concentration of impurity k in process i;

[0037] 5. Maximum outlet concentration limit:

[0038]

[0039] 6. Impurity recovery concentration limit:

[0040] C i,k,out -YC i,l,out ≤0,i∈{PU,TU},k,l∈C (9)

[0041] Y is an integer variable, indicating whether the material is recycled, Y=1 if recycled, Y=0 if not recycled;

[0042] 7. Water flow demand of water-using unit:

[0043]

[0044] For further improvement, the specific steps of step three are as follows:

[0045] Step 3.1: Initialize various hyperparameters in the molecular swarm algorithm MGA, including:

[0046] The number of molecules N, the number of variables n, the maximum number of iterations Maxgen, the gravitational constant G 0 , the upper and lower limits of the kth variable [L k , U k ], k = 1, 2...n, critical distance constant

[0047] Step 3.2: Initialize the initial multi-impurity water network superstructure scheme.

[0048] Randomly generate a first-generation multi-impurity water network superstructure scheme X 0 , and record the sth generation molecular group X s for and generate the initial velocity of the kth dimension of the wth molecule is the wth molecule of the sth generation, representing a water network superstructure; w = 1, 2, ..., N, k = 1, 2, ..., n;

[0049] Step 3.3: Introducing the penalty function

[0050] The penalty function method is used to integrate the constraint conditions of the objective function into the objective function, thereby converting it into an unconstrained problem and obtaining the water network superstructure optimization model:

[0051]

[0052] Among them, g c (x) and h d (x) are the cth inequality constraint and the dth equality constraint respectively; m is the number of inequality constraints, x is a solution in the solution space, p is the number of equality constraints, T represents the matrix transpose, and n is the dimension of the solution space; R n is the solution space; the minimization optimization problem of the objective function is converted into a maximization optimization problem, and the following function is constructed after combining the penalty function:

[0053] fit(x)=1 / (Z(x)+G(x)+H(x)) (12)

[0054] Among them, fit(x) is the fitness function, G(x) and H(x) are the penalty functions constructed by the inequality constraints and equality constraints in the model respectively; Z(x) is the original objective function;

[0055]

[0056] r and b are numbers that decrease with the number of iterations; m is the number of inequality constraints, α is a constant, and β is a constant;

[0057] The expression of r in the formula is:

[0058]

[0059] a 1 is a constant, Maxgen is the maximum number of iterations;

[0060] r is a number that decreases with the number of iterations. In the early stages of iteration, r is a relatively large number and the penalty for violating the constraints is very large, which makes the mass of the inferior celestial bodies smaller and the gravitational effect smaller, making it easier for the inferior solutions to be eliminated.

[0061] The expression of b in the formula is:

[0062]

[0063] s is the current iteration number, a 2 is a constant

[0064] Convert the equality constraints to inequality constraints, that is:

[0065]

[0066] Step 3.4: Quality Normalization

[0067] The fitness value is normalized to obtain the quality, so that the quality varies in the range of [0, 1]. The calculation formula for quality normalization is as follows:

[0068]

[0069] Among them, fit w is the mass of the wth molecule. max and fit min are the maximum and minimum fitness values ​​in the population in a certain iteration, m w is the normalized mass of the w-th molecule;

[0070] Step 3.5: Calculate critical distance

[0071] First, calculate the distance between molecule w and molecule o in the kth dimension in the sth iteration;

[0072]

[0073] Then, according to the interaction force between molecules, a critical distance is defined When the distance between molecules w and o in the kth dimension is greater than When the distance between molecules w and o in the kth dimension is less than When , there is repulsion between molecules; is the force exerted by the sth generation molecule o on the molecule w in the kth dimension, s is the current iteration number, is the k-th dimension position of the s-th generation molecule w, is the position of the kth dimension of the sth generation molecule o;

[0074] At the same time, the critical distance r cdDynamically adjust as the number of iterations increases, the formula is as follows:

[0075]

[0076] is the critical distance of the kth dimension of the sth generation, is the initial critical distance, U k is the upper limit of the k-th dimension, L k is the lower limit of the kth dimension, γ is a constant, s is the current number of iterations, and e is the natural logarithm;

[0077] Step 3.6: Calculate the interaction force between different molecules based on the universal gravitation formula:

[0078] During the sth iteration, the wth molecule In dimension k, it is acted upon by molecule o as follows:

[0079]

[0080] Among them: ε is a constant used to prevent errors caused by the denominator being zero, M w (s) represents the mass of molecule w in the sth iteration, M o (s) represents the mass of molecule o in the sth iteration. wo (s) represents the Euclidean distance between molecules w and o in the s-th iteration; G(s) represents the gravitational constant in the s-th iteration, which will decrease with time, and γ is a constant; G 0 is the initial gravitational constant;

[0081] Based on critical distance The formula to determine whether the force at this time is gravitational or repulsive is as follows:

[0082]

[0083] In summary, the resultant force of the kth latitude in the sth generation acting on the celestial body w is for:

[0084]

[0085] Step 3.7: Calculate the acceleration

[0086] When a molecule is acted upon by other molecules, it will generate acceleration. Then the acceleration of molecule w in the kth dimension is The ratio of its force to mass is:

[0087]

[0088] Step 3.8: Update the position of the molecule

[0089] In each iteration, the molecule updates its velocity and position according to the calculated acceleration as follows:

[0090]

[0091] in, represents the velocity of molecule w in the kth dimension of the sth generation, represents the velocity of molecule w in the kth dimension of the s+1th generation, represents the position of molecule w in the kth dimension of the sth generation, represents the position of the molecule w in the kth dimension of the s+1th generation; rand w Represents a random number between 0 and 1;

[0092] The rule of cross-relocation is adopted to simulate the effect of infinite domain in finite domain, so as to reduce or eliminate the influence of boundary effect. The specific formula is as follows:

[0093]

[0094] Step 3.9: Repeat the iteration and output the results:

[0095] Repeat steps 3.3 to 3.8 until the termination condition is reached.

[0096] The advantages of the present invention are as follows:

[0097] The present invention takes the total annual cost of the water network as the objective function, and takes the conservation of materials in the water network, the environmental restrictions on the concentration of impurities in wastewater discharge, and the requirements for the recovery of valuable substances as constraints, and then uses MGA to solve the optimal value, which can better converge to the global optimal solution and reduce the difficulty of solving the problem. BRIEF DESCRIPTION OF THE DRAWINGS

[0098] Figure 1 is an optimization model of an embodiment of the present invention;

[0099] Figure 2 is a graph of the critical distance changing with the number of iterations;

[0100] Figure 3 is a flow chart of the MGA algorithm of the present invention;

[0101] Figure 4 This is a comparison chart between the molecular swarm algorithm and other algorithms. DETAILED DESCRIPTION

[0102] Taking a water network superstructure with two impurities, three water use units and three treatment units as an example, the optimization model is as follows Figure 1As shown. Among them, PU is a water use unit, TU is a processing unit, SU is a separation unit, and MU is a mixing unit. In this embodiment, the specific process of a multi-impurity water network superstructure optimization solution method based on MGA according to the present invention is described in detail. The drawings involved only illustrate certain implementation cases of the present invention and should not be regarded as limiting the scope. For technical personnel in this profession, in combination with specific examples and the actual situation of their own technical level, other relevant information can also be obtained according to the drawings.

[0103] The specific implementation includes the following steps:

[0104] Step 1: Get the original data, including:

[0105] Annual operation time H, fresh water cost coefficient C f , annual investment factor AR of the processing unit, investment cost coefficient IC of the processing unit t , Operation cost coefficient OC of processing unit t , α is the cost function index, 0<α≤1, the cost coefficient of waste treatment C sw , the cost coefficient of the medicine C med 、Recycling price of valuable materials C rec , the mass load of impurity k in process i, Δm i，k，out , the residual rate β of impurity k in process i i，k , the maximum inlet concentration of impurity k in process i Minimum inlet concentration of impurity k in process i Concentration limit of impurity K discharged into the environment Water flow requirement WR for treatment unit i i .

[0106] Step 2: Construct a multi-impurity water network superstructure model based on the objective function and constraints;

[0107] Step 2.1: Use formula (1) to construct the objective function Z:

[0108]

[0109] in, For fresh water costs, is the annual investment cost of the processing unit, is the annual operating cost of the treatment unit, The waste treatment costs of the treatment unit, To process the unit drug costs, Economic benefits generated by the recovery of valuable substances.

[0110] Step 2.2: Use equations (2) to (10) to construct constraints:

[0111] 8. Total water balance:

[0112] For water units:

[0113]

[0114] For processing units:

[0115]

[0116] 9. Impurity mass balance:

[0117] For water units:

[0118]

[0119] For processing units:

[0120]

[0121] 10. Maximum import concentration limit:

[0122]

[0123] 11. Minimum import concentration limit:

[0124]

[0125] 12. Maximum outlet concentration limit:

[0126]

[0127] 13. Impurity recovery concentration limit:

[0128] C i,k,out -YC i,l,out ≤0,i∈{PU,TU},k,l∈C(9)

[0129] 14. Water flow demand of water-using unit:

[0130]

[0131] Among them, f i is the amount of fresh water entering the system, F i is the total water flow entering process i, is the flow rate of impurity k entering process i, Y is an integer variable indicating whether the substance is recycled, and X ij is the amount of wastewater reused from process j to process i, W i is the wastewater discharge of process i, C i,k,out is the outlet concentration of impurity k in process i, C i,k,inis the inlet concentration of impurity k in process i.

[0132] Step 3: Use the MGA algorithm to iteratively solve the optimization problem.

[0133] Step 3.1: Initialize various hyperparameters in MGA, including:

[0134] The number of molecules N, the number of variables n, the maximum number of iterations Maxgen, the gravitational constant G 0 , the upper and lower limits of the kth variable [L k , U k ], k = 1, 2...n, critical distance constant

[0135] Step 3.2: Initialize the initial multi-impurity water network superstructure scheme.

[0136] Randomly generate a first-generation multi-impurity water network superstructure solution X within the required range 0 , and the sth generation molecular group is is the wth molecule of the sth generation, representing a water network superstructure. At the same time, the initial velocity of the kth dimension of the wth molecule is generated w=1,2,...,N,k=1,2,...,n.

[0137] Step 3.3: Introducing the penalty function

[0138] The water network superstructure optimization problem is a multi-constrained optimization problem. The penalty function method is used to integrate the constraints into the objective function, thus converting it into an unconstrained problem.

[0139] The water network superstructure optimization model can be expressed as:

[0140]

[0141] Among them, g i (x) and h j (x) are the inequality constraints and equality constraints respectively.

[0142] The original problem is a minimization optimization problem, which is converted into a maximization optimization problem, and the following function is constructed after combining the penalty function:

[0143] fit(x)=1 / (Z(x)+G(x)+H(x)) (12)

[0144] Among them, G(x) and H(x) are the penalty functions constructed by the inequality constraints and equality constraints in the model respectively:

[0145]

[0146]

[0147] The expression of r in the formula is:

[0148]

[0149] r is a number that decreases with the number of iterations. In the early stages of the iteration, r is a large number, and the penalty for violating the constraint is large, which makes the mass of the inferior celestial body smaller and the gravitational effect smaller, making it easier for the inferior solution to be eliminated.

[0150] The expression of b in the formula is:

[0151]

[0152] Since equality constraints are relatively strict constraints, they are sometimes difficult to satisfy for optimization algorithms. Therefore, we consider converting equality constraints into inequality constraints, namely:

[0153] |h j (x)|≤b(17)

[0154] b is also a number that decreases with the number of iterations. In the early stages of the iteration, b is a larger number, which is conducive to the convergence of the algorithm and global search. In the later stages of the iteration, b is smaller, which is conducive to improving the local optimization ability of the algorithm. This can relax the constraints and is conducive to the convergence of the algorithm.

[0155] Step 3.4: Quality Normalization

[0156] In the algorithm, the mass of an individual is associated with its fitness value. The better the fitness value, the greater the mass of the individual, and the greater the potential to find the global optimal solution around it. However, according to the universal gravitation calculation formula, gravity is proportional to mass. Using fitness to directly calculate may cause the calculation result of gravity to be too large or too small, thereby affecting the update of the individual's speed and acceleration, and further affecting the convergence performance of the algorithm. Therefore, it is necessary to normalize the fitness value to obtain the mass so that the mass varies in the range of [0, 1]. The calculation formula for mass normalization is as follows:

[0157]

[0158] Among them, fit w is the mass of the wth molecule. max and fit min They are respectively the maximum and minimum fitness values ​​in the population in a certain iteration.

[0159] Step 3.5: Calculate critical distance

[0160] First, the distance between molecule w and molecule o in the kth dimension in the sth iteration is calculated.

[0161]

[0162] Then, according to the interaction force between molecules, a critical distance is defined When the distance between molecules w and o in the kth dimension is greater than When the distance between molecules w and o in the kth dimension is less than When , there is repulsion between molecules.

[0163] At the same time, the critical distance r cd It can be adjusted dynamically as the number of iterations increases, and the formula is as follows:

[0164]

[0165] In the early stage of iteration, the mixed effect of attraction and repulsion between different molecules makes the distribution of the population in the solution space more uniform, which is conducive to improving the global search ability of the algorithm and avoiding premature convergence. As the number of iterations increases, the critical distance The decrease of the critical distance indicates that the intermolecular attraction is increasing, which is conducive to the molecules to move closer to the better individuals; at the same time, there is still repulsion between the molecules, which prevents the molecules from aggregating too early and falling into the local optimum. In this way, the local optimization ability of the algorithm can be gradually improved while ensuring the global search. Finally, when the critical distance When it continues to decrease or even approaches zero, the force between molecules is almost gravitational, which accelerates the convergence of the algorithm.

[0166] Step 3.6: Calculate the interaction force between different molecules based on the universal gravitation formula

[0167] During the sth iteration, the wth molecule In the k-dimensional space, it is acted upon by the molecule o. as follows:

[0168]

[0169] Where: ε represents a very small constant, which is set to prevent errors caused by the denominator being zero. w (s) represents the mass of molecule w, M o (s) represents the mass of molecule o. wo (s) represents the Euclidean distance between molecules w and o. G(s) represents the gravitational constant of the sth generation, which will decrease with time, and γ is a constant.

[0170] Based on critical distance The formula to determine whether the force at this time is gravitational or repulsive is as follows:

[0171]

[0172] In summary, the resultant force of the kth latitude in the sth generation acting on the celestial body w is for:

[0173]

[0174] Step 3.7: Calculate the acceleration

[0175] When a molecule is acted upon by other molecules, it will generate acceleration. Then the acceleration of molecule w in the kth dimension is The ratio of its force to mass is:

[0176]

[0177] Step 3.8: Update the position of the molecule

[0178] In each iteration, the molecule updates its velocity and position according to the calculated acceleration as follows:

[0179]

[0180] in, represents the velocity of molecule w in the kth dimension of the sth generation, represents the velocity of molecule w in the kth dimension of the s+1th generation, represents the position of molecule w in the kth dimension of the sth generation, Indicates the position of the molecule w in the kth dimension of the s+1th generation. rand w Represents a random number between 0 and 1.

[0181] Since the positions of the generated new generation of molecules may not be within the allowed upper and lower limits, the rule of cross-relocation is adopted. This method reduces or eliminates the impact of boundary effects by simulating the effects of infinite domains in finite domains. The specific formula is as follows:

[0182]

[0183] Step 3.9: Repeat the iteration and output the results

[0184] Repeat steps 3.3 to 3.8 until the termination condition is reached (such as reaching the maximum number of iterations).

[0185] Finally, the molecular positions of the optimal solution, i.e., the optimal multi-impurity water network superstructure, are output.

[0186] Figure 4This is a comparison chart between the molecular swarm algorithm and other algorithms. As can be seen from the figure, the molecular swarm algorithm has a better convergence speed, and due to the effect of intermolecular repulsion, it can avoid premature convergence of the algorithm, and can also jump out of the local optimum during the iteration process and continue to explore the solution space to find the global optimal solution.

[0187] The above description is only a preferred embodiment of the present invention and does not limit the present invention in any form. Any simple modification, equivalent change and modification made to the above embodiment according to the technical essence of the present invention without departing from the technical solution of the present invention still falls within the scope of the technical solution of the present invention.

Claims

1. A water network superstructure optimization method based on molecular swarm algorithm, characterized in that: The following steps are involved: Step 1: Obtain the original data of water network superstructure; Step 2: establishing a multi-impurity water network superstructure model of the water network superstructure according to the original data; Step three: using a molecular swarm algorithm to solve the multi-impurity water network superstructure model until a termination condition is reached to obtain an optimal multi-impurity water network superstructure.

2. The water network superstructure optimization method based on molecular swarm algorithm according to claim 1, characterized in that: The termination condition is a preset number of iterations or convergence.

3. The water network superstructure optimization method based on molecular swarm algorithm according to claim 1, characterized in that: The raw data of the water network superstructure include: Annual operation time H, fresh water cost coefficient C f , annual investment factor AR of the processing unit, investment cost coefficient IC of the processing unit t , Operation cost coefficient OC of processing unit t , cost function index α, 0<α≤1, cost coefficient C of waste treatment sw , the cost coefficient of the medicine C med 、Recycling price of valuable materials C rec , the mass load of impurity k in process i, Δm i，k，out , the residual rate β of impurity k in process i i，k , the maximum inlet concentration of impurity k in process i Minimum inlet concentration of impurity k in process i Concentration limit of impurity K discharged into the environment and water unit PU i Water flow requirement WR i .

4. The water network superstructure optimization method based on molecular swarm algorithm according to claim 2, characterized in that: The objective function of the multi-impurity water network superstructure model is as follows: in, is the fresh water cost, f i is the amount of fresh water entering the system, PU is the water use unit; is the annual investment cost of the processing unit, F i is the total water flow entering process i, TU is a processing unit; is the annual operating cost of the treatment unit, is the waste treatment cost of the treatment unit, W i is the wastewater discharge of process i; To process the unit drug costs, Economic benefits generated by the recovery of valuable substances, is the flow rate of impurity k entering process i.

5. The water network superstructure optimization method based on molecular swarm algorithm according to claim 4, characterized in that: The constraints of the objective function are as follows:

1. Total water balance: For water units: X ij is the amount of wastewater reused from process j to process i, X ji is the amount of wastewater reused from process i to process j; For processing units:

2. Impurity mass balance: For water units: C i,k,out is the outlet concentration of impurity k in process i, C j,k,out is the outlet concentration of impurity k in process j, Δm i,k,tot is the mass load of impurity k in process i; For processing units:

3. Maximum import concentration limit:

4. Minimum import concentration limit: C i,k,in is the inlet concentration of impurity k in process i; 5. Maximum outlet concentration limit:

6. Impurity recovery concentration limit: C i,k,out -YC i,l,out ≤0,i∈{PU,TU},k,l∈C (9) Y is an integer variable, indicating whether the material is recycled, Y=1 if recycled, Y=0 if not recycled; 7. Water flow requirements of water-using units:

6. The water network superstructure optimization method based on molecular swarm algorithm according to any one of claims 1 to 5, characterized in that: The specific steps of step three are as follows: Step 3.1: Initialize various hyperparameters in the molecular swarm algorithm MGA, including: The number of molecules N, the number of variables n, the maximum number of iterations Maxgen, the gravitational constant G0, the upper and lower limits of the kth variable [L k , U k ], k = 1, 2...n, critical distance constant Step 3.2: Initialize the initial multi-impurity water network superstructure scheme. Randomly generate a first-generation multi-impurity water network superstructure scheme X0, and record the s-th generation molecular group X s for and generate the initial velocity of the kth dimension of the wth molecule is the wth molecule of the sth generation, representing a water network superstructure; w = 1, 2, ..., N, k = 1, 2, ..., n; Step 3.3: Introducing the penalty function The penalty function method is used to integrate the constraint conditions of the objective function into the objective function, thereby converting it into an unconstrained problem and obtaining the water network superstructure optimization model: minZ(x) Among them, g c (x) and h d (x) are the cth inequality constraint and the dth equality constraint respectively; m is the number of inequality constraints, x is a solution in the solution space, p is the number of equality constraints, T represents the matrix transpose, and n is the dimension of the solution space; R n is the solution space; the minimization optimization problem of the objective function is converted into a maximization optimization problem, and the following function is constructed after combining the penalty function: fit(x)=1 / (Z(x)+G(x)+H(x)) (12) Among them, fit(x) is the fitness function, G(x) and H(x) are the penalty functions constructed by the inequality constraints and equality constraints in the model respectively; Z(x) is the original objective function; r and b are numbers that decrease with the number of iterations; m is the number of inequality constraints, α=1, β=1; The expression of r in the formula is: a1=0.00005, Maxgen is the maximum number of iterations; r is a number that decreases with the number of iterations. In the early stages of iteration, r is a relatively large number and the penalty for violating the constraints is very large, which makes the mass of the inferior celestial bodies smaller and the gravitational effect smaller, making it easier for the inferior solutions to be eliminated. The expression of b in the formula is: s is the current iteration number, a2=0.00005 Convert the equality constraints to inequality constraints, that is: |h j (x)|≤b (17) Step 3.4: Quality Normalization The fitness value is normalized to obtain the quality, so that the quality varies in the range of [0, 1]. The calculation formula for quality normalization is as follows: Among them, fit w is the mass of the wth molecule. max and fit min are the maximum and minimum fitness values ​​in the population in a certain iteration, m w is the normalized mass of the w-th molecule; Step 3.5: Calculate critical distance First, calculate the distance between molecule w and molecule o in the kth dimension in the sth iteration; Then, according to the interaction force between molecules, a critical distance is defined When the distance between molecules w and o in the kth dimension is greater than When the distance between molecules w and o in the kth dimension is less than When , there is repulsion between molecules; is the force exerted by the sth generation molecule o on the molecule w in the kth dimension, s is the current iteration number, is the k-th dimension position of the s-th generation molecule w, is the position of the kth dimension of the sth generation molecule o; At the same time, the critical distance r cd Dynamically adjust as the number of iterations increases, the formula is as follows: is the critical distance of the kth dimension of the sth generation, is the initial critical distance, U k is the upper limit of the k-th dimension, L k is the lower limit of the kth dimension, γ = 1, s is the current number of iterations, and e is the natural logarithm; Step 3.6: Calculate the interaction force between different molecules based on the universal gravitation formula: During the sth iteration, the wth molecule In the k-dimensional space, it is acted upon by the molecule o. as follows: Among them: ε is a constant used to prevent errors caused by the denominator being zero, M w (s) represents the mass of molecule w in the sth iteration, M o (s) represents the mass of molecule o in the sth iteration. wo (s) represents the Euclidean distance between molecules w and o in the s-th iteration; G(s) represents the gravitational constant in the S-th iteration, which will decrease with time, and γ is a constant; G0 is the initial gravitational constant; Based on critical distance The formula to determine whether the force at this time is gravitational or repulsive is as follows: In summary, the resultant force of the kth latitude in the sth generation acting on the celestial body w is for: Step 3.7: Calculate the acceleration When a molecule is acted upon by other molecules, it will generate acceleration. Then the acceleration of molecule w in the kth dimension is The ratio of its force to mass is: Step 3.8: Update the position of the molecule In each iteration, the molecule updates its velocity and position according to the calculated acceleration as follows: in, represents the velocity of molecule w in the kth dimension of the sth generation, represents the velocity of molecule w in the kth dimension of the s+1th generation, represents the position of molecule w in the kth dimension of the sth generation, represents the position of the molecule w in the kth dimension of the s+1th generation; rand w Represents a random number between 0 and 1; The rule of cross-relocation is adopted to simulate the effect of infinite domain in finite domain, so as to reduce or eliminate the influence of boundary effect. The specific formula is as follows: Step 3.9: Repeat the iteration and output the results: Repeat steps 3.3 to 3.8 until the termination condition is reached.

Images