Artificial Electric Field Algorithm with Greedy State Transition Strategy for Spherical Multiple Traveling Salesmen Problem

The multiple traveling salesman problem (MTSP) is an extension of the traveling salesman problem (TSP). It is found that the MTSP problem on a three-dimensional sphere has more research value. In a spherical space, each city is located on the surface of the Earth. To solve this problem, an integer-serialized coding and decoding scheme was adopted, and artificial electric field algorithm (AEFA) was mixed with greedy strategy and state transition strategy, and an artificial electric field algorithm based on greedy state transition strategy (GSTAEFA) was proposed. Greedy state transition strategy provides state transition interference for AEFA, increases the diversity of population, and effectively improves the accuracy of the algorithm. Finally, we test the performance of GSTAEFA by optimizing examples with different numbers of cities. Experimental results show that GSTAEFA has better performance in solving SMTSP problems than other swarm intelligence algorithms.


Introduction
Combinatorial optimization problem is a kind of optimization problem. Optimization problems can be divided into two categories: one is continuous variable problems, the other is discrete variable problems. Problems with discrete variables are called combinatorial. In the case of continuous variables, it is generally a set of real numbers, or a function; In combinatorial problems, an object is found from an infinite or countable infinite set-typically an integer, a set, a permutation, or a graph. In general, these two types of problems have quite different characteristics, and the methods of solving them are quite different. Metaheuristic algorithms have been widely studied in recent years because of their advantages in multi-modal, large-scale, and highly constrained problems.
The traveling salesman problem (TSP) is a classical combinatorial optimization problem [1]. Many practical problems can be modeled as TSP, such as route planning [2], production scheduling [3], and emergency management [4]. In recent years, there are many scholars have made continuous exploration on solving TSP in two-dimensional space. There are some precise algorithms and heuristic algorithms that can effectively solve TSP. Gouveia, Leitner, and Ruthmair used a branch cutting algorithm [5]. Kinable, Smeulders, Delcour, and Spieksma used an accurate algorithm [6]. Asadpour et al. used Approximation Algorithm [7]. Both Genetic Algorithm [8] and Memetic Algorithm [9] are applied to asymmetric traveling salesman problem. Algorithms to solve the TSP problem also include Evolutionary Algorithm [10], Artificial Bee Colony Algorithm [11], Ant Colony Algorithm [12], Water Cycle Algorithm [13], Discrete Bat Algorithm [14], Random-key Cuckoo Search [15], Black Hole Algorithm [16], Discrete Comprehensive Learning Particle Swarm Optimization Algorithm [17], Simulated Annealing Algorithm [18], Fruit-fly Optimization Algorithm [19], and so on. However, the classic TSP cannot be modeled in some special cases. For example, suppose that a company has multiple salespeople living in different cities, and the company expects them to visit the city no less than a certain number of times to meet their minimum wage [20]. Therefore, the MTSP is introduced on the basis of the current problem [21]. ACO has been widely used in the research field of MTSP and achieved good results. In particular, Ghafurian and Javadian [22] designed an ant-based system to solve the MTSP problem at a fixed destination. Yousefkhoshbakht et al. [23] solved the MTSP problem by mixing the insertion, switching, and 2-opt operators with ant colony algorithm. Changdar et al. [24] proposed a genetic ACO algorithm for solving MTSP. However, it is difficult to optimize the parameters of ant colony optimization system because of the large number of inherent parameters. This further leads to its disadvantages such as slow computation speed and easy local convergence. When more complex problems are solved, unacceptable performance occurs. Other metaheuristic algorithms have been proposed to solve MTSP. Malmborg [25] proposed a genetic algorithm based on two chromosomes. Carter and Ragsdale [26] designed a genetic algorithm for the two-part chromosome technique. Venkatesh et al. [27] mixed artificial bee colony algorithm with weed invasion algorithm. Chen et al. [28] proposed an improved two-part Wolf pack search algorithm. Zhou et al. [29] proposed a improved PSO to solve MTSP.
Compared to the standard TSP, the spherical space TSP is more practical. In recent years, many scholars have studied TSP in three-dimensional space. Literature [30] solves 3D-TSP with multi-dimensional city location. Literature [31] adopts genetic algorithm to solve the TSP problem on the sphere. Literature [32] describes the distance distribution between random points on a sphere. The ant colony algorithm [33], firefly algorithm [34], an improved flower pollination algorithm based on greedy strategy [35], and the discrete improved cuckoo algorithm [36] are used to solve STSP. As we know, the surface of the earth on which we live is very close to a sphere. In many fields of study, such as chemistry, biology, and physics, structures such as atoms, molecules, and proteins are represented by spheres. Therefore, the solution of spherical MTSP can be directly applied to micro and macro systems such as path planning, robot path planning, picking, and placing parts. Therefore, it is of great significance to study spherical MTSP problem.
AEFA is a physics-based optimization algorithm designed by Yadav et al. in 2019 to solve continuous optimization problems. Its main characteristics include less optimization parameters, rapid convergence in the process of optimization, high accuracy, and low complexity of the algorithm [37]. Since its inception, AEFA has been widely used in continuous optimization problems [38], practical engineering problems [39], image matching problems [40], multi-objective optimization [41], nine parameters triode PV model best estimate [42], cancer detection [43], white blood cells [44], the optimal purchasing problem [45], feature classification [46], assembly line balancing problem [47], and wind turbines active loss in distribution network and the distribution of the voltage deviation problems [48] made a wide range of applications. In this paper, an improved AEFA is proposed. The greedy state transition strategy is introduced on the basis of AEFA, and the improved algorithm is used to solve SMTSP. Similar to STSP, the city coordinates of the SMTSP are extended to the surface of a sphere, and the coordinates of the points and the minimum spherical distance between the points are known. The improved AEFA is used to solve SMTSP. Experiments show that the proposed algorithm is feasible and effective in solving SMTSP.
The rest of this paper is organized as follows: "Spherical Geometry and Definition of the Problem" gives a brief introduction to spherical geometry and the definition of the problem. "Related Works" introduces the related work, including solution encoding and decoding, greedy state transfer strategy, AEFA, and GSTAEFA. In "Experimental Setup and Discussion of Results," the simulation experiment and numerical results are analyzed. Finally, "Conclusions and Further Works" gives summary and future prospects.

Spherical Geometry and Definition
of the Problem Figure 1 shows a sphere. The radius r of the sphere is the distance from the point on the sphere to the center. In threedimensional space, the coordinates of each point on the sphere need to satisfy the relationship as shown in Eq. (1): where r is the distance from a point on the sphere to the center, and it is the radius of the sphere, (x, y, z) are the coordinates of the points on the sphere in the three dimensions [34] (Fig. 2).

Representation of Points on a Sphere
Points on a sphere can be represented in cartesian coordinates so that the shortest distance between two points on the sphere can be calculated using vector methods. The Cartesian vector point function can represent that point [35]: In general, the standard equation for a sphere can be expressed in parametric form, and a spherical coordinate equation can be expressed with two parameters u and v (1) (2) P(u) = (x(u, v), y(u, v), z (u, v)).
between 0 and 1. A sphere of radius r with its center at the (0,0,0) can be expressed as [36] where the two parameters u and v represent the longitude and latitude lines, respectively. In this paper, we use a sphere with radius 1 for simulation. Therefore, it is easy to calculate the length of the path between two cities (Fig. 3).

Shortest Path Between Two Cities on Sphere
On the surface of the earth, a geodesic is the curve that shows the shortest distance between two cities. The point at which the tangent plane through the center intersects the earth forms a great circle. The length of the shortest line between two points on the sphere is the length of the lower arc between two points in the great circle that passes through them [34]. In Fig. 4, there are two cities City1 and City2, in order to obtain the shortest distance between them. Take the scalar product of two vectors [34]:    The smaller angle between v1 and v2 is θ. Scalar values are computed as follows [34]: The formula for the shortest distance between City1 and City2 on the sphere is [31] Combining formulas (6) to (8), the distance formula is deduced as follows [31]:

Definition of the Problem
In MTSP, there are multiple traveling agents visiting different cities, and each city can only be visited once and only by one person. Given a set of n cities that a salesperson will visit, TSP looks for the shortest possible trip for the salesperson, which visits each city only once. MTSP can be briefly described as follows: Give an undirected graph G = (V, A), It consists of a set of vertices V and a set of arcs A, where m represents the total number of salespeople. The objective function is to divide V into m nonempty subsets calculate the minimum cost sum of all sets. Then the objective function value of MTSP can be obtained [49] as follows: where the first sum represents the path through m salespeople, and the second sum represents the cycle of all the cities visited by the ith salesperson (the index of the first city visited by the ith salesperson is 1, and the index of the last city is n i ); x i j,j+1 represents the distance between the city j visited by the ith salesperson and j + 1 ; x i n i ,1 represents the distance between the last city visited by the ith salesperson and the first city; n i value should not be less than the minimum number of cities specified by each salesperson [49].

Motivation
At present, for many practical application problems such as TSP problem, MST problem, MTSP problem, people cannot effectively obtain the optimal solution in a reasonable time range. Therefore, metaheuristic optimization algorithms are being applied by more and more researchers from various fields to obtain satisfactory solutions such as hybrid PSO-GA algorithm for constrained optimization problems [61], a hybrid GSA-GA algorithm for constrained optimization problems [59], GA-GSA hybrid algorithm for optimizing the performance of industrial systems using uncertain data [60], and a hybrid ITL-HHO algorithm for solving numerical and engineering optimization problems [62]. This paper mainly studies the improvement of artificial electric field algorithm and its application in practical optimization problems. The motive of the improvement is to improve the mining capacity and exploration capacity of the original algorithm, respectively, so as to further improve the overall optimization ability of the algorithm. Artificial electric field algorithm is a new intelligent optimization algorithm based on physical laws proposed in 2019, so its improvement and application is relatively small. Therefore, the purpose of improving the algorithm is to solve the practical optimization problem on one hand, on the other hand, to verify the optimization ability of the artificial electric field algorithm, overcome the shortcomings of the original algorithm, so as to broaden the application scope of the artificial electric field algorithm in practical life.

Solution Encoding and Decoding
Sequence coding is the simplest and most efficient way to express MTSP solutions, using a series of numbers to represent each solution. There are two common approaches to using the sequential encoding approach to represent an MTSP solution. One uses a path sequence and a breakpoint sequence, and the other uses a path sequence and a city number sequence. Here, we adopt the real encoding method, which is suitable for most continuous algorithms to solve the MTSP problem. The following example shows how the realcoded approach represents an MTSP solution. For example, suppose two salespeople have to visit 10 cities. We use a single sequence coding scheme to represent a feasible solution. Compared with double sequence coding, the calculation is simpler and the processing is more convenient.

Greedy State Transition Strategy
The state transition algorithm (STA) [50], a new heuristic search algorithm proposed by Zhou et al. [51], has outstanding performance in solving optimization problems and is widely used in many practical problems [52]. In this paper, we introduce three kinds of state transfer operators, and combine them to propose a state transfer strategy. The three operators are Swap, Reverse, and Insert of three basic transformations. In other words, instead of considering global optimality, the algorithm gets a local optimal solution in a certain sense. In the state transition strategy, we do not consider globally. After any operator produces a new solution, if the result is better, we accept the local change, so that the algorithm will keep getting closer to the global optimal solution. The MTSP problem is a combinatorial optimization problem, and the solution form is a sequence. Every state transfer will generate a feasible solution. Greedy thought is introduced, assuming that an optimal sequence can be obtained without considering the time complexity.
Swap Operator: Pick two positions i and j at random from the initial solution, and the corresponding values of i and j are exchanged.
where xnew is the immediate solution of x generated, i and j are two different positional indexes in x, respectively. Swap (x i,j ) swaps the values of the positions corresponding to i and j (Fig. 6).
Reverse Operator: Pick two positions i and j at random from the initial solution, and then reverse the values between i and j.
where xnew is the immediate solution of x generated, and i and j are two different positional indexes in x, respectively. Reverse (x i,j ) indicates the reverse order of all values from i to j (Fig. 7).
Insert Operator: Two positions i and j are randomly selected from the initial solution, and the corresponding value for position i is inserted after position j. where xnew is the immediate solution of x generated, and i and j are two different positional indexes in x, respectively. Insert (x i,j ) means that the value corresponding to position i is inserted after position j (Fig. 8). Figure 9 shows the proposed greedy state transition algorithm.

The Principle of AEFA
In AEFA systems, the individual is described as a point particle that moves under electrostatic attraction. AEFA is based on Coulomb's law [37]. The initialization of an individual in AEFA is shown in Eq. (14).
where N represents the total number of individuals, that is, the total number of charged particles in the system, and D represents the dimension of the problem. In the algorithm, it is necessary to obtain the individual optimal particle in each iteration process as shown in Eq. (15).
In the process of iteration t, the force of any individual J acting on particle I is expressed by Eq. (16). where Q i (t) and Q j (t) are the amount of charge in different individuals. K is Coulomb's constant, ε is used to prevent the denominator from being 0, and R is the Euclidian distance between two individuals. Formulas 17 and 18 give the calculation method of R.
where parameters α and K 0 are two important parameters to balance the exploration and development capabilities of the algorithm. The resultant force of particle i in the iterative process is shown in Formula (19).
where rand is a random number in the range of [0, 1] that conforms to the uniform distribution. Then the resultant force can be calculated through the above equation, and the electric field strength E of the current particle can be calculated next.
The current individual acceleration is deduced by calculating the power plant strength and the resultant force.
Combined with Newton's law of motion, given the current velocity and acceleration of the individual, calculate the velocity and position at the next moment.
After mastering the above principles of AEFA, we can describe AEFA in the form of pseudocode, as shown in the figure below (Fig. 10): Figure 11 is the flowchart of AEFA. The algorithm of AEFA is simple in logic and easy to implement, with less parameter to be optimized.

Artificial Electric Field Algorithm Based on Greedy State Transition
Although AEFA has many advantages, it still needs to be improved. For example, the algorithm can be improved to speed up convergence and reduce the probability of falling into local solutions. Therefore, we add the greedy state transition strategy in this paper to improve its performance after combining with the original algorithm. The local search mechanism improves the ability of AEFA to search deeper in the solution space to obtain better solutions. In the optimization process of AEFA, the optimal individuals in each generation are selected for local search using GST strategy to obtain better individuals. After each GST operation, a second-best solution is obtained, which is returned to AEFA and the next search begins. Figure 12 is the pseudocode for GSTAEFA, and Fig. 13 is the flowchart of GSTAEFA.

GSTAEFA Complexity Analysis
The complexity of any optimization algorithm is mainly derived from the algorithm itself and the evaluation of the objective function. In our algorithm, the complexity mainly

Experimental Setup and Discussion of Results
We chose a more suitable hardware environment for the simulation of the algorithm, and maintained the principle of fairness in evaluating the algorithm performance. All

Test Problems
Since the problem in this study is new, some test questions are needed to evaluate the proposed metaheuristic algorithm. For this purpose, seven size test questions were considered from the data collected in the relevant article [34]. The number of cities in each question is 25, 50, 100, 150, 200, 300, and 400. Finally, in order to improve the persuasiveness of the algorithm, we increased the number to 500 cities to test the performance of the algorithm.

Experiment Setting
In the experimental process, to ensure fairness, the number of individuals of all the comparison algorithms was set to 30, and the algorithm controls the end of the algorithm according to the value of Iterations, and its value is 1000. To carry out the final experiment on the test problem, some problems are considered to obtain better solutions, and a fair comparative study is also provided for the algorithm. Each algorithm uses standard parameters. To reduce the randomness of the algorithm, each algorithm was run 30 times for each test problem. The minimum value, maximum value, average value, and standard deviation of the obtained objective function were compared, and Friedman rank sum test and P value test were performed on the results of all algorithms. GSTAEFA: NP = 30, Coulomb's constant K 0 is 500, the number of local searches is 100.
FPA: NP = 30, the transfer probability of global and local pollination of flowers was 0.8 [54].

Analysis and Comparison of Experimental Results
Next, we explore how the number of cities affects the solution and compare the test results of nine algorithms.
The complexity of the problem is mainly determined by the number of cities, the more cities, the more difficult to solve. Each algorithm runs 30 times independently for different numbers of cities. Table 1 shows the test values of all algorithms at different points on the sphere, where "City" represents the number of cities. As shown in Table 1, GSTAEFA was significantly superior to other algorithms. In all cases, GSTAEFA obtained the highest accuracy of the results, and at the same time, the average value in all cases was also the best. In the case of fewer cities, the results obtained by different algorithms are very close, but GSTAEFA has the best stability. It can be analyzed from the results that when the number of points on the sphere increases, the searching ability of other algorithms will be greatly reduced, and when the number of cities is 200, 300, 400, and 500, the stability of SA, GA, and ABC is better than that of GSTAEFA, while the results of GSTAEFA algorithm are 52.1173, 91.4491, 133.3439, and 178.2932, respectively. Development performance is optimal. Figures 14,15,16,17,18,19,20, and 21 show a line plot of the 30 times result optimal values. When the dimension of the problem is low, the results of the algorithm proposed by us are not different from other algorithms, but the accuracy of the results is always higher than other algorithms. With the continuous increase of the number of cities, the accuracy of GSTAEFA algorithm gradually opened a gap with other     22, 23, 24, 25, 26, 27, 28, and 29, it can be seen that GSTAEFA has good stability, and the results of each time are around the range of the variance diagram. It can be seen from the scatter diagram that the points of AEFA and FA are relatively discrete and have poor stability.
To visually demonstrate the power of all algorithms when solving SMTSP, we use the convergence curve shown in Figs. 30, 31, 32, 33, 34, 35, 36, and 37. Before the execution of the algorithm is finished, the results of other algorithms have been stabilized in a local solution range. The convergence effect of SA was similar to that of GSTAEFA, but the convergence speed and accuracy of SA were lower than that of GSTAEFA. Although the optimal solution procedure may not be given, we can conclude that on most data sets, GSTAEFA can obtain the best results and converge faster than the compared algorithm. Figures 38,39,40,41,42,43,44, and 45 intuitively show the optimal path searched by GSTAEFA. All cities and routes can be seen simultaneously in transparent mode. It should be noted that the routes shown in the figure are all approximately globally optimal, especially in large-scale problems.
We present a bar chart of five path lengths with different number of cities, as shown in Figs. 46, 47, 48, 49, 50, 51, 52, and 53. It can be seen that the path length of each    salesperson is smaller than other algorithms, so the sum of all paths must be better than other algorithms. So far, our proposed algorithm has a great advantage in finding the optimal solution, regardless of the size of the city. As shown in Figs. 54,55,56,57,58,59,60, and 61, we show a bar chart of CPU elapsed time for different algorithms. In the case of lower dimension, the running time of GSTAEFA is relatively large, but the convergence accuracy is the highest. And as the problem dimension increases, GSTAEFA's time performance becomes more and more outstanding.

Statistical Analysis
In order to conduct statistical analysis on the experimental results, two nonparametric tests, Friedman test [57] and Wilcoxon signed-rank test [58], were used to evaluate the performance of the proposed algorithm. Friedman's test is a statistical test for the consistency of multiple correlated samples. The test calculates the ranking of the shortest distance achieved by each algorithm on each instance, with a best value of 1 and a worst value of 9. Table 2 shows the average ranking obtained by the nine algorithms using the Friedman test at the 95% confidence level. As you can see, our proposed algorithm ranks first. Experimental results show that the proposed GSTAEFA search strategy is effective in solving the SMTSP problem, and its capabilities are the best among the compared algorithms.
Wilcoxon signed-rank test is a nonparametric test. Statistical test of rank is based on sample observation. The test results were used to compare the differences between the algorithms. Less than 5% indicates that the proposed algorithm is good enough. Table 3 shows the detection results. P values less than 0.05 indicate that GSTAEFA has a good performance in statistics. In the table, the value is greater than or equal to 0.05. All detection results are less than 0.05, indicating that the performance of GSTAEFA is irreplaceable.
From the above two non-parametric test, we can clearly see that GSTAEFA is superior to other advanced methods in solving SMTSP, and the algorithm is obviously different from other algorithms. It shows that GSTAEFA is irreplaceable when solving SMTSP.

Conclusions and Further Works
In this paper, a hybrid algorithm GSTAEFA based on traditional AEFA hybrid greedy state transition strategy is proposed. AEFA is an effective naturally inspired algorithm which relies on Coulomb's electrostatic law and Newton's laws of motion. It has been successfully applied to various types of optimization problems and obtained efficient solutions. Although AEFA has a satisfactory exploration experience, its development capability is its    main disadvantage. In order to overcome these problems of AEFA and improve its development ability, greedy state transfer strategy is considered in the proposed algorithm to effectively improve its development ability, which will make it more efficient, more effective, and easier to solve various combinatorial optimization problems. We apply GSTAEFA to solve a new spherical multi-traveler problem, which is considered as the shortest path problem for multiple travel agencies around the world, and test the performance of the algorithm with different dimensional examples. The results show that the artificial electric field algorithm based on greedy state transition strategy is superior to the classical algorithm. In addition, GSTAEFA is weak in time complexity, so it still has room for further optimization. In the future, we will consider reducing the time cost of the algorithm while ensuring the solution accuracy. We can try to reduce the local search times according to probability to reduce the time complexity. Finally, we focus on single-target SMTSP, so it is interesting that our algorithm can be used to solve multi-target projects. As a future research, the formulas proposed in this study can be used to obtain Pareto optimal or nondominated solutions through multi-objective precise and multi-objective metaheuristic solutions.