Improved Seagull Optimization Algorithm Combined with an Unequal Division Method to Solve Dynamic Optimization Problems

: The numerical solution of the dynamic optimization problem is often sought for chemical processes, but the discretization of control variables is a difﬁcult problem. Firstly, based on the analysis of the seagull optimization algorithm, this paper introduces the cognitive part in the process of a seagull’s attack behavior to make the group approach the best position. Secondly, the algorithm adds the mechanism of natural selection, where the ﬁtness value is used to sort the population, and the best half is used to replace the worst half, so as to ﬁnd out the optimal solution. Finally, the improved seagull optimization algorithm (ISOA) is combined with the unequal division method to solve dynamic optimization problems. The feasibility of the method is veriﬁed by three practical examples of dynamic optimization in chemical industry.


Introduction
With the continuous development of society, the optimization of chemical processes, such as batch reaction process and tubular reaction process, has always been a research hotspot. Improving chemical efficiency and reducing energy consumption are the goals of chemical industry. However, as the chemical engineering model becomes more and more complex and the number of control variables increases, the stability of the model becomes lower, and inputs and outputs are related by nonlinear models. Therefore, we need to propose a new method to optimize it. Chemical dynamic optimization problems are usually complicated, which seeks one or more control variables to optimize the performance index of the objective function. At present, direct methods [1] and indirect methods [2] are often used to solve the optimal control problem. The direct method is to transform the continuous optimal control problem into a finite-dimensional nonlinear programming problem for solution through discretization. The indirect method mainly solves the optimality condition (necessary condition) of the original problem, thereby indirectly obtaining the optimal solution of the original problem. So far, the optimization methods for chemical engineering problems mainly include dynamic programming (DP) [3], control variable parameterization (CVP) [4], and also swarm intelligence algorithm [5,6].
Dynamic programming does not require gradient information, but it faces the difficult question of dimensionality. In 1989, Luus [7] proposed the iterative dynamic programming method, which mainly included the introduction of area reduction processing and the iterative application of dynamic programming. At the same time, the iterative dynamic programming algorithm needs to simultaneously discretize the two dimensions of time and space, and obtain the optimal solution through iterative calculation. However, if the dimension of time and space are relatively large, the implementation process of this method is very complicated.
In 1970, Sargent et al. [8] proposed the CVP algorithm for the first time. The main idea of the algorithm is the subsection control approximation strategy. At the same time, in 1979, Sargent et al. [9] used the CVP algorithm to solve dynamic optimization problems for the first time, and many people knew the algorithm. In recent years, the CVP algorithm has been improved by many people. Xu et al. [10] proposed a dynamic optimization problem grid reconstruction strategy based on pseudo Wigner-Ville analysis, using Wigner-Ville to perform grid reconstruction combined with variable time node CVP method to achieve switching between time nodes. Li et al. [11] proposed a variable time node control vector parameterization method to optimize the control variable and time node. Li et al. [12] proposed a control parameterization method with variable time nodes, using a multitime grid control vector parameterization method, so that each control variable had a corresponding time grid. Binder et al. [13] introduced wavelet analysis to further fine-tune the grid divided by the CVP method, so that the dynamic optimization trajectory was closer to the true value.
In recent years, intelligent algorithms have been widely used to solve dynamic optimization cases. The algorithm has the advantages of simplicity, ease of implementation, and strong search capabilities, so it is also a hot topic for scholars to study. Based on particle swarm algorithm (PSO) and control variable parameterization, Shi [14] et al. proposed to solve the dynamic optimization of the chemical process. First, the particle swarm was used for the first optimization, and then the optimized value was passed to the control variable parameter method for the second optimization. Xu et al. [15] proposed a biogeographic learning particle swarm algorithm to improve the efficiency of learning. Sarkar et al. [16] used genetic algorithms to optimize feed bioreactors: multiple control variables. Anand et al. [17] proposed a tabu search method to dynamically optimize the copolymerization reactor. Schluter et al. [18] proposed the design of an integrated process and control system based on the ant colony optimization algorithm. Egea et al. [19] proposed an enhanced scatter search method for the global dynamic optimization of nonlinear processes. Nikumbh [20] proposed biogeography-based optimization for dynamic optimization of chemical reactors.
In this study, we propose a new improvement strategy for the seagull optimization algorithm, which mainly includes the improvement strategy of the cognitive part and the natural selection mechanism. Finally, by combining the improved seagull optimization algorithm (ISOA) with the unequal segmentation method, three typical chemical examples are solved, and different methods are compared and analyzed. Each chemical case has a good performance index, and the results show the optimization ability of the algorithm.
The main contributions and innovations of this research are as follows: (1) The cognitive part strategy is introduced to prevent the algorithm from entering the local optimum.
(2) The use of the natural selection mechanism further enhances the optimization performance of the seagull optimization algorithm, so as to find the optimal solution faster.
(3) The improved seagull optimization algorithm combined with the unequal division method is used to solve dynamic optimization problems. (4) In order to verify the feasibility of the improved algorithm, this study introduces three chemical engineering cases for testing. At the same time, we analyze the performance indicators of the three cases and compare them with other methods. Finally, we also compare the execution efficiency of the equal division method and the unequal division method.

Problem Description
A typical dynamic optimization problem (DOP) of a chemical engineering process can be stated as follows: where t 0 , t f are the initial time and terminal time, respectively, x(t) is the state vector; represents integral function; u max and u min are the lower and upper bounds of the control vector, respectively; x min and x max are the minimum and maximum values of the state vector, respectively; x(0) is the initial state, and f is a dynamic system equation.

Unequal Division Method
The piecewise constant method is the most common method in control variable parameterization (CVP). Within a time domain [t 0 , t f ], which is equally divided into n subregions using the equal division method, i.e., the length of each subregion is d = (t f − t 0 )/n, the Runge-Kutta method can be used to integrate the ordinary differential equations of the model. At the same time, the corresponding algorithm is used to optimize the control variables, so that the control variables produce a series of trajectories. Compared with the equal division method, the unequal division method may produce a more precise control trajectory, and thus obtain better performance indicators, therefore, we use the unequal division method. The steps for unequal division are as follows: 1 We assume that the time domain needs to be divided into n segments, a set of parameters are randomly initialized in [t 0 , t f ], namely τ 1 , τ 2 , . . . , τ n ∈ [t 0 , t f ]; 2 The final time parameter is as follows, t i represents time parameter.

Seagull Optimization Algorithm (SOA)
Inspired by the migration and attack behavior of seagulls, Dhiman et al. [21] proposed a random search algorithm based on swarm intelligence-seagull optimization algorithm (SOA). When using the seagull optimization algorithm, assume that the size of the seagull population is pop, the dimension of the problem space is D, and the position of the seagull is Finally, the migration and attack behavior of seagulls are the process of individual position updating. The seagull migration and attack process are shown in Figure 1.

Migration Behavior
The migration behavior belongs to the global search. Migration refers to th ment of seagulls from one position to another, but three conditions that seagulls meet during migration are avoiding collisions, direction of the best position, proaching the best position. ① Avoiding collisions: In order to prevent collisions between adjacent seagulls, tional variable A is introduced. At the same time, the algorithm uses the A to update the position of the seagull in the iterative process.

Migration Behavior
The migration behavior belongs to the global search. Migration refers to the movement of seagulls from one position to another, but three conditions that seagulls should meet during migration are avoiding collisions, direction of the best position, and approaching the best position. 1 Avoiding collisions: In order to prevent collisions between adjacent seagulls, an additional variable A is introduced. At the same time, the algorithm uses the variable A to update the position of the seagull in the iterative process.
where cs(iter) is the new position of the seagull after collision avoidance, iter is the current iteration number, ps(iter) is the initial position of the seagull, and A is the motion behavior of the seagull in a given search space. The equation of A is as follows: where the value of A is adjusted linearly, the value of f c decreases linearly to 0, and Maxiter is the maximum number of iterations. 2 Direction of the best position: After the seagull is satisfied that it will not collide with other individuals, the seagull will move in the direction of the best position. The equation is as follows: where B is a random number that balances global search and local search, zbest(iter) represents the best position in the population, and the value of rd is a random number of (0,1). 3 Approaching the best position: After the seagull moves to a position where it does not collide with other seagulls, it moves in the direction of the best position to reach a new position. ds(iter) = |cs(iter) + ms(iter)|, where ds(iter) represents the best fit search seagull.

Attack Behavior
The attack behavior belongs to the local search. When seagulls attack their prey, they move in a spiral shape in the air, using x, y, and z to represent their motion behavior. The equation for exercise behavior is as follows: where r is the spiral radius in the motion of the seagull, and µ, v are the correlation constants of the spiral shape. θ represents the angle, which is a random number in [0, 2π]. The equation for the attacking behavior of seagulls is as follows: where ps(iter) saves the best solution and updates the position of other search seagulls.

The Steps of the Seagull Optimization Algorithm
Step 1: Initialize the size of the seagull population pop, set the parameters mainly including f c , u, v, θ and Maxiter, and initialize the position of pop seagulls randomly.
Step 2: Calculate the fitness value of each seagull, compare the fitness values between individuals, and find global optimal values of the current population.
Step 3: Enter the main loop and use Equations (4)- (8) to calculate the new position of the seagull after migrating.
Step 4: Use Equations (9)-(13) to calculate the final position of the seagull.
Step 5: Compare the fitness between individuals in the current seagull population again to find the global optimal value.
Step 6: Judging that the algorithm reaches the termination condition during the execution process. If it is reached, it ends; otherwise, the calculation goes to Step 4-Step 6 to continue the position update.
Step 7: Output the global optimal position and fitness value of the SOA algorithm.

Cognitive Part
The SOA has the advantages of a good optimization effect, simple operation, and few parameter settings. However, the SOA has a premature phenomenon, and it is easy to fall into a local optimum. For Equation (13), it belongs to the local search stage, which is easy to fall into a local extremum, and may not find the global optimal value. Therefore, we introduce the cognitive part into Equation (13) to reduce or jump out of the local extremum. The improvement of the equation is as follows: where w represents the inertia factor, and the value in this paper is 0.95. The role of w is to enhance the ability of each individual seagull to learn from the global optimal individual, so that the algorithm avoids falling into local extrema.

The Mechanism of Natural Selection
Inspired by genetic algorithms, in order to further improve the optimization performance of the algorithm, this paper introduces the mechanism of natural selection. In order to find the optimal value, keep the good seagulls and eliminate the bad ones. In each iteration, the individuals of the seagull population are sorted by the calculated fitness value, and the best half of the individuals are used to replace the worst half of the individuals, while retaining the historical best value of the original seagull individuals.
The seagull optimization algorithm is easy to fall into local extrema and the solution accuracy is not high. However, the improved seagull optimization algorithm adds the cognitive part and natural selection mechanism, so that the algorithm has better solving ability, and avoids falling into local extrema.

Experimental Setup
In order to verify the feasibility and effectiveness of the improved algorithm, this paper introduces three typical chemical engineering cases for testing. The experimentation and algorithms are implemented in MATLAB R2017a (version 9.

Algorithm Steps
First, we assume that the number of segments in the time domain is n, and use Equation (3) for unequal division. Then, we express the control variable u as a constant function in each interval, that is, u = [u 1 , u 2 , . . . , u n ] represents optimizing the solution of the objective function. Finally, this paper uses the Runge-Kutta method to integrate the dynamic model and the algorithm to solve each interval. Each case is run independently 20 times, and the best value of the case is recorded.
Step 1: Setting parameters and initializing the seagull population, use Equation (3) to divide the time domain unequally, and then apply the Runge-Kutta method to solve differential equations and find the best value based on fitness.
Step 2: Iterative optimization. Use Equations (4)- (12) and (14) to update the position of the seagull, and then apply the Runge-Kutta method to calculate each interval differential equation. Finally, find the best value based on fitness, and add the mechanism of natural selection before the next iteration of the algorithm.
Step 3: Determine whether the algorithm meets the end condition. If yes, jump out of the loop and output the best result; if not, return to Step 2 to continue optimization.
Step 4: End the program. Table 1 shows the setting of experimental parameters, and Algorithm 1 is the algorithm flow chart.

Case 1: Batch Reactor
The batch reactor is widely used in the process industries. The batch reactor is used for solid dissolution, dissolution, product mixing, chemical reaction, crystallization, etc. A typical batch reactor consists of a tank with an agitator and an integrated heating/cooling system. Liquids and solids are usually charged through a connection port on the top cover of the reactor, vapor and gas are discharged from the top, and liquids and solids are usually discharged from the bottom. In this paper, the reaction process of the batch reactor is A→B→C. A is the reaction product, B is the intermediate product, and C is the product. After the reaction is completed, the concentration of the intermediate product B reaches the optimal value. The mathematical model of batch reactor [22] is as follows: where C A is the concentration of reactants, mol/L; C B is the concentration of intermediate products, mol/L, and t f is the time at the end of the reaction, h; T is the control temperature, K, and J is the performance index, mol/L. In the experiment, we set the pop to 600, and the Maxiter is 200.
Analysis of the experimental results of Case 1: Table 2 presents results obtained by our ISOA and other methods. The "-" in Table 3 means that the reference did not provide a value. Shi et al. [14] used the PSO-CVP optimization strategy to find the optimal value of 0.6105359. Rajesh et al. [22] [30] used IDP to obtain a value of 0.610775 (N = 80). In this paper, the optimal value of using equal division method is 0.61059223 (N = 30), and the optimal value of unequal division method is 0.610794203 (N = 30). In terms of method execution efficiency, the solution speed of ISOA combined with the equal division method is 272 s, but the solution speed of ISOA combined with the unequal division method is 118 s, and the solution efficiency of the latter is higher. Figures 2-4 are the optimal temperature control sequences. Figure 5 is the optimal trajectory of state variables.

Case 2: Parallel Reaction Problem of Tubular Reactor
Parallel reaction means that reactants can carry out two or more different re in parallel to obtain different products. For example, ethanol can undergo two re of dehydration and dehydrogenation in parallel. In our paper, we introduce a d optimization problem with saturation characteristics of control variables, which by many scholars. The process of this problem is A → B and A → C . By deter

Case 2: Parallel Reaction Problem of Tubular Reactor
Parallel reaction means that reactants can carry out two or more different re in parallel to obtain different products. For example, ethanol can undergo two re of dehydration and dehydrogenation in parallel. In our paper, we introduce a d optimization problem with saturation characteristics of control variables, which by many scholars. The process of this problem is A → B and A → C . By deter

Case 2: Parallel Reaction Problem of Tubular Reactor
Parallel reaction means that reactants can carry out two or more different reactions in parallel to obtain different products. For example, ethanol can undergo two reactions of dehydration and dehydrogenation in parallel. In our paper, we introduce a dynamic optimization problem with saturation characteristics of control variables, which is cited by many scholars. The process of this problem is A→B and A→C. By determining the optimal value of the control variable, the concentration of by-product B is maximized after the reaction is completed. The mathematical model of this reaction is as follows [31]: where x 1 (t) is the concentration of reactants A, mol/L; x 2 (t) is the concentration of byproduct B, mol/L, and J is the performance index, mol/L; u(t) represents the control variable, 10 6 , and t f is terminal time, s. In the experiment, we set the pop to 800, and the Maxiter is 200.
Analysis of the experimental results of Case 2: Table 3 presents results obtained by our ISOA and other methods. The "-" in Table 3 means that the reference did not provide a value. Rajesh et al. [22] used the ant colony structure optimization (ACSO) to obtain the value of 0.57284 (N = 4). Zhou et al. [26] used a control parameterized particle swarm optimization (CP-PSO) method to obtain the value of 0.573543, and proposed a control parameterized adaptive particle swarm optimization (CP-APSO) to obtain the value of 0.573544. Dadebo [30] applied IDP to obtain the value of 0.57353. Tanartkit and Biegler [31] used a nested simultaneous strategy to obtain the value of 0.5738, which is the best reported value. Vassiliadis [32] solved the model with a value of 0.57353. Biegler [33] proposed a combination of continuous quadratic programming and orthogonal collocation method to obtain values of 0.56910 (CVP) and 0.57322 control vector iteration (CVI). In this paper, the optimal value of using equal division method is 0.573073 (N = 40), and the optimal value of unequal division method is 0.573435 (N = 40). In terms of method execution efficiency, the solution speed of ISOA combined with the equal division method is 326 s, but the solution speed of ISOA combined with the unequal division method is 302 s, and the optimization efficiency of the latter is about 10% higher than that of the equal division method. Figures 6-8 are the optimal control trajectories, and Figure 9 is the optimal state variable trajectory.
where 1 ( ) x t is the concentration of reactants A , mol/L; 2 ( ) x t is the concentration of by-product B , mol/L, and J is the performance index, mol/L; ( ) u t represents the control variable, 10 6 , and f t is terminal time, s. In the experiment, we set the pop to 800, and the Maxiter is 200.
Analysis of the experimental results of Case 2: Table 4 presents results obtained by our ISOA and other methods. The "--" in Table 4 means that the reference did not provide a value. Rajesh et al. [22] used the ant colony structure optimization (ACSO) to obtain the value of 0.57284 ( N = 4). Zhou et al. [26] used a control parameterized particle swarm optimization (CP-PSO) method to obtain the value of 0.573543, and proposed a control parameterized adaptive particle swarm optimization (CP-APSO) to obtain the value of 0.573544. Dadebo [30] applied IDP to obtain the value of 0.57353. Tanartkit and Biegler [31] used a nested simultaneous strategy to obtain the value of 0.5738, which is the best reported value. Vassiliadis [32] solved the model with a value of 0.57353. Biegler [33] proposed a combination of continuous quadratic programming and orthogonal collocation method to obtain values of 0.56910 (CVP) and 0.57322 control vector iteration (CVI). In this paper, the optimal value of using equal division method is 0.573073 ( N = 40), and the optimal value of unequal division method is 0.573435 ( N = 40). In terms of method execution efficiency, the solution speed of ISOA combined with the equal division method is 326 s, but the solution speed of ISOA combined with the unequal division method is 302 s, and the optimization efficiency of the latter is about 10% higher than that of the equal division method. Figures 6-8 are the optimal control trajectories, and Figure 9 is the optimal state variable trajectory.

Case 3: Tubular Reactor
The tubular reactor was proposed by some authors [34,35]. The reaction follows: in a certain length of pipeline, the catalysts A and B are mixed, a product C is produced after the reaction is completed. The goal of the p distribute the catalyst in the designated pipeline so that the concentration product C reaches the optimal value. The mathematical model of the tu problem is as follows:

Case 3: Tubular Reactor
The tubular reactor was proposed by some authors [34,35]. The reaction follows: in a certain length of pipeline, the catalysts A and B are mixed, product C is produced after the reaction is completed. The goal of the distribute the catalyst in the designated pipeline so that the concentration product C reaches the optimal value. The mathematical model of the tu problem is as follows:

Case 3: Tubular Reactor
The tubular reactor was proposed by some authors [34,35]. The reaction process is as follows: in a certain length of pipeline, the catalysts A and B are mixed, and the target product C is produced after the reaction is completed. The goal of the problem is to distribute the catalyst in the designated pipeline so that the concentration of the target product C reaches the optimal value. The mathematical model of the tubular reactor problem is as follows: where x A (z) is the concentration of catalyst A, mol/L; x B (z) is the concentration of catalyst B, mol/L; z f is the pipe length of the tubular reactor, m; u(z) is the mixing fraction of catalyst A, %; J is the performance index of the model, mol/L. In the experiment, we set the pop to 800, and the Maxiter is 200. Analysis of the experimental results of Case 3: Table 4 presents results obtained by our ISOA and other methods. Rajesh et al. [22] used the ant colony structure optimization (ACSO) to obtain values of 0.47615 (N = 5) and 0.47527 (N = 20). Peng et al. [29] proposed an improved knowledge evolution algorithm (IKEA) to obtain values of 0.4757 (N = 20) and 0.47761-0.47768 (N = 100). Dadebo [30] applied IDP to obtain values of 0.47527 (N = 20) and 0.47695 (N = 40). Chen et al. [36] used nonuniform discretization-based control vector parameterization to obtain the value of 0.47771 (N = 15). The reference of [37] used trigonometric differential evolution (TDE) to obtain values of 0.47527 (N = 20) and 0.47683 (N = 40). The reference of [38] used bees algorithm (BA) to obtain the value of 0.4744 (N = 20). In this paper, the optimal value of using the equal division method is 0.47721 (N = 40), and the optimal value of unequal division method is 0.47770 (N = 40). The value of the latter is better than the former. All in all, the value of the unequal division method is better than the value of other references [22,29,30,37,38]. In terms of method execution efficiency, the solution speed of ISOA combined with the equal division method is 575 s, but the solution speed of ISOA combined with the unequal division method is 285 s. Figures 10-12 are the optimal control trajectories, and Figure 13 is the optimal state variable trajectory.

Conclusions
In order to achieve the sustainable development of society and reduce e sumption, this paper proposes an improved seagull optimization algorithm (I Figure 13. Optimal state variable trajectory.

Conclusions
In order to achieve the sustainable development of society and reduce energy consumption, this paper proposes an improved seagull optimization algorithm (ISOA) combined with an unequal division method to solve dynamic optimization problems through the analysis of chemical processes. First of all, we analyzed that the seagull optimization algorithm (SOA) may easily fall into local extrema. In order to improve this shortcoming, this paper introduced the cognitive part and the mechanism of natural selection to improve SOA, thereby effectively avoiding the algorithm from falling into the local optimum, and enhancing the optimization performance of the algorithm. Secondly, chemical dynamic optimization problems usually need to discretize the time domain, use a piecewise function to approach the optimal solution, and convert the problem into NLP for processing. Therefore, this paper used the unequal division method to divide the time domain so that the calculation becomes simple.
Finally, in order to verify the feasibility and effectiveness of improving SOA, the ISOA is combined with an unequal division method to solve dynamic optimization problems. The feasibility of the method is verified by comparing with other algorithms. In this paper, many reference methods use basically the uniform discretization of control variables, and the time domain needs to be divided into multiple segments. At the same time, we use the unequal division method of the control variable, and the time domain is divided into fewer segments, which makes the calculation more simple. For example, the time efficiency of ISOA combined with the unequal division method is higher than that of ISOA combined with the equal division method.