# A Multiobjective Particle Swarm Optimization Algorithm Based on Competition Mechanism and Gaussian Variation

Authors: Hongli Yu; Yuelin Gao; Jincheng Wang
Journal: Complexity
 (2020)

---

## Abstract

In order to solve the shortcomings of particle swarm optimization (PSO) in solving multiobjective optimization problems, an improved multiobjective particle swarm optimization (IMOPSO) algorithm is proposed. In this study, the competitive strategy was introduced into the construction process of Pareto external archives to speed up the search process of nondominated solutions, thereby increasing the speed of the establishment of Pareto external archives. In addition, the descending order of crowding distance method is used to limit the size of external archives and dynamically adjust particle parameters; in order to solve the problem of insufficient population diversity in the later stage of algorithm iteration, time-varying Gaussian mutation strategy is used to mutate the particles in external archives to improve diversity. The simulation experiment results show that the improved algorithm has better convergence and stability than the other compared algorithms.


---

## Body

## 1. Introduction

In many engineering problems, the problems are composed of multiple goals that influence and conflict with each other. In solving practical issues, people often encounter multiple objectives that need to obtain the best optimal solution at the same time, that is, multiobjective optimization problems. The optimization problem has more than one optimization objective and needs to be processed at the same time, which becomes a multiobjective optimization problem (MOP). Usually, the optimal solution of the multiobjective optimization problem obtained after analyzing the objective function is the Pareto optimal solution set [1]. Therefore, in solving the multiobjective optimization problem [2], the following three key elements need to be solved: (1) The solution set is as close to the Pareto front as possible. (2) Keep the good diversity of the population as much as possible. (3) Make the particles effectively and uniformly distributed in the solution space. In recent years, in order to effectively solve MOP, multiobjective optimization algorithms based on different optimization theories have been continuously proposed. Deb et al. [3] proposed a multiobjective evolutionary computing algorithm based on nondominated sorting, called NSGA-II, and introduced a new selection operator into the algorithm to reduce the complexity of the algorithm. Zitzler et al. [4] proposed SPEA-II using the idea of Pareto domination. Coello et al. [5] proposed an improved multiobjective particle swarm algorithm that uses the concept of Pareto dominance to determine the flight direction of a particle and it maintains previously found nondominated vectors in a global repository that is later used by other particles to guide their own flight. This algorithm improves the efficiency of solving multiobjective problems. Tsai et al. [6] proposed an improved multiobjective particle swarm optimizer with proportional distribution and jump improved operation, named PDJI-MOPSO, for dealing with multiobjective problems. PDJI-MOPSO maintains diversity of newly found nondominated solutions via proportional distribution and combines advantages of wide-ranged exploration and extensive exploitations of PSO in the external repository with the jump improved operation to enhance the solution searching abilities of particles. Introduction of cluster and disturbance allows the proposed method to sift through representative nondominated solutions from the external repository and prevent solutions from falling into local optimum. Mnif et al. [7] introduced a new approach called multiobjective firework algorithm (MFWA). Mellouli et al. [8], in order to solve the two-dimensional cutting stock problem, combined genetic algorithm with linear programming model to estimate the best Pareto frontier for these two goals. The Pareto front provided by this algorithm is very close to the optimal front. Ali et al. [9] proposed a modified variant of Differential Evolution (DE) algorithm for solving multiobjective optimization problems. The proposed algorithm, named Multiobjective Differential Evolution Algorithm (MODEA), utilizes the advantages of Opposition-Based Learning for generating an initial population of potential candidates and the concept of random localization in mutation step. Zhou et al. [10] proposed a generic transformation strategy that can be referred to as the Mucard strategy, which converts an MSCCOP into a low-dimensional multiobjective optimization problem (MOP) to simultaneously obtain all the (near-) optima of the constrained optimization problems in a single algorithmic run. Vargas et al. [11] studied the performance of the combination of adaptive penalty technology called APM and GDE3 algorithm. Sun et al. [12] proposed a novel multiobjective particle swarm optimization algorithm based on Gaussian mutation and an improved learning strategy. This method uses Gaussian mutation strategy to improve the consistency of external archives and current population. In order to improve the global optimal solution, different learning strategies are proposed for nondominant and dominant solutions. An indicator is proposed to measure the distribution width of nondominated solution sets generated by various algorithms. Coello et al. [13] proposed that the external set could be used to retain the nondominated solutions found in the iterative process, which improved the efficiency of the algorithm. Tsai et al. [6] used the specific global optimal value to replace the individual optimal value fully using the guiding role of the global optimal value. Zhang [14] proposed an improved MOPSO algorithm with a mutation operator that can maintain the diversity of optimal solutions and has good convergence. Zhang [15] proposed a MOPSO algorithm based on fuzzy dominance, and the experimental results showed the effectiveness of the proposed algorithm. Tao [16] proposed a multiobjective optimization algorithm combining PSO and difference algorithms. By generating common new particles and updating the particle velocity formula, the search efficiency of the algorithm was effectively improved. Li [17] proposed an improved MOPSO algorithm that updates the optimal position of all particles through the Pareto dominance relationship. The experiment shows that the proposed algorithm can obtain a better noninferior solution. Ni [18] proposed an adaptive dynamic recombinant PSO algorithm that adopts a high-level clustering algorithm. The experimental results show that this algorithm can improve the convergence speed and evolutionary ability of the algorithm.More and more improved algorithms and strategies are used to solve various multiobjective optimization problems. However, few researchers have improved the algorithm from the perspective of the balance between local search capability (exploitation) and global search capability (exploration) to solve multiobjective optimization problems. Related research has pointed out that the effective balance between the exploration and exploitation of intelligent algorithms will have a vital impact on the optimization performance of the algorithm [19]. Different from other study results, in this study, we proposed a multiobjective particle swarm optimization algorithm based on competition mechanism strategy and Gaussian mutation to balance the exploration and exploitation of the algorithm and enable the algorithm to search the optimal location and converge to the Pareto front more quickly. Finally, through the simulation test of the multiobjective test functions and the comparison with other multiobjective optimization algorithms, the results show that the algorithm proposed in this paper is superior to the compared algorithm in terms of convergence and population distribution.The rest of the paper is organized as follows. Section2 describes the multiobjective optimization problem and the basic concept. Section 3 presents the PSO algorithm. In Section 4, an improved multiobjective particle swarm optimization algorithm is introduced. In Section 5, the evaluation indicators and test functions used are introduced. The numerical experiment results and data analysis are described in Section 6. Finally, conclusions and respects are presented in Section 7.


## 2. Multiobjective Optimization Problem and the Basic Concept

### 2.1. Formal Definition of Multiobjective Optimization Problem

Generally, a multiobjective optimization problem (MOP) includes a set of objective functions and some constraints\enleadertwodots. Without loss of generality, an MOP withm objective functions and n decision variables can be described as follows [20, 21]:(1)minfx=f1x,f2x,…,fmxs.t.gix≤0,i=1,2,…,phjx=0,j=1,2,…,qx∈Rn,where x=x1,x2,…,xn∈Rn are the n-dimensional decision variables. f=f1,f2,…,fm is the objective function, which contains one or more objective functions; and f=f1,f2,…,fm is the multiobjective function m≥2; fkx:Rn⟶R,k=1,2,…,m is the k−th objective function; gi:Rn⟶R,i=1,2,…,p is the i−th inequality constraint; hj:Rn⟶R,j=1,2,…,q represents the j−th equality constraint condition. X=x∈Rn|gi≤0,i=1,2,…,p;hj=0,j=1,2,…,q is the feasible domain of decision variables.


### 2.2. Pareto Optimal Solution Concepts

The Pareto optimal solution was discovered by the Italian economist Pareto. It was originally limited to the field of economics; this rule has gradually extended to various fields of social life and was deeply recognized by people as well. The strict Pareto optimal solution can be described by multiobjective mathematical programming. Assuming that there are several objectives at the same time, these objectives are independent of each other and cannot be weighed and summed. New optimization theories are needed to solve such problems. Generally, if one subobjective is improved, other subobjectives will be sacrificed. So, it is impossible to improve all subobjectives at the same time. Pareto optimal solution is also called nondominated solution. In a multiobjective optimization problem, due to factors such as conflict and incomparability among various subobjectives, a solution is often the best in one subobjective and may be the worst in other subobjectives; if there is one solution, improving any subobjective function will inevitably weaken at least one other subobjective function, which is called nondominated solution or Pareto optimal solution. All Pareto optimal solutions constitute the Pareto optimal solution set, and these solutions are mapped by the objective function to form the Pareto optimal front. Pareto proposed the concept of nondominated set of multiobjective solutions in 1986, which is defined as follows: assuming any two solutionsS1 and S2 for all objectives, if S1 is better than S2, then we say that S1 dominates S2, and if S1 is not dominated by other solutions, then S1 is called a nondominated solution, also called a Pareto optimal solution.


### 2.3. Related Definitions of Multiobjective Optimization Problem

The following introduces some basic concepts in the multiobjective optimization problem [22, 23].Definition 1.
(feasible domain). In the decision space, the feasible domain is represented byX, and its expression is(2)X=x∈Rn|gi≤0,i=1,2,…,p;hj=0,j=1,2,…,q.Definition 2.
(feasible solution). For the pointx in the decision space, if x∈X, x is the feasible solution.Definition 3.
(Pareto dominates). For any vectorx1,x2, if and only if(3)∀i∈Rn,fix1≤fix2,∧∃j∈Rn,fjx1<fjx2,then x1 dominates x2, so x1≺x2.Definition 4.
(Pareto optimal solution). The Pareto optimal solution is also known as the nondominant solution. For a solutionx∗ in the feasible domain, if x∗ is not dominated by any other solution in the feasible domain, x∗ is called the Pareto optimal solution, and its definition is as follows:(4)¬∃x∈X:x≺x∗.Definition 5.
(Pareto optimal solution set, PS). The set of all nondominated solutions in the feasible domain is called the Pareto optimal solution set and can be defined as follows:(5)PS=x∗|¬∃x∈X:x≺x∗.Definition 6.
(Pareto optimal frontier, PF). The target vector set corresponding to the Pareto optimal solution set is the Pareto optimal frontier, also known as the Pareto optimal front end or Pareto equilibrium surface, which is defined as follows:(6)PF=Fx∗=f1x∗,f2x∗,…,fdx∗|x∗∈PS.


## 2.1. Formal Definition of Multiobjective Optimization Problem

Generally, a multiobjective optimization problem (MOP) includes a set of objective functions and some constraints\enleadertwodots. Without loss of generality, an MOP withm objective functions and n decision variables can be described as follows [20, 21]:(1)minfx=f1x,f2x,…,fmxs.t.gix≤0,i=1,2,…,phjx=0,j=1,2,…,qx∈Rn,where x=x1,x2,…,xn∈Rn are the n-dimensional decision variables. f=f1,f2,…,fm is the objective function, which contains one or more objective functions; and f=f1,f2,…,fm is the multiobjective function m≥2; fkx:Rn⟶R,k=1,2,…,m is the k−th objective function; gi:Rn⟶R,i=1,2,…,p is the i−th inequality constraint; hj:Rn⟶R,j=1,2,…,q represents the j−th equality constraint condition. X=x∈Rn|gi≤0,i=1,2,…,p;hj=0,j=1,2,…,q is the feasible domain of decision variables.


## 2.2. Pareto Optimal Solution Concepts

The Pareto optimal solution was discovered by the Italian economist Pareto. It was originally limited to the field of economics; this rule has gradually extended to various fields of social life and was deeply recognized by people as well. The strict Pareto optimal solution can be described by multiobjective mathematical programming. Assuming that there are several objectives at the same time, these objectives are independent of each other and cannot be weighed and summed. New optimization theories are needed to solve such problems. Generally, if one subobjective is improved, other subobjectives will be sacrificed. So, it is impossible to improve all subobjectives at the same time. Pareto optimal solution is also called nondominated solution. In a multiobjective optimization problem, due to factors such as conflict and incomparability among various subobjectives, a solution is often the best in one subobjective and may be the worst in other subobjectives; if there is one solution, improving any subobjective function will inevitably weaken at least one other subobjective function, which is called nondominated solution or Pareto optimal solution. All Pareto optimal solutions constitute the Pareto optimal solution set, and these solutions are mapped by the objective function to form the Pareto optimal front. Pareto proposed the concept of nondominated set of multiobjective solutions in 1986, which is defined as follows: assuming any two solutionsS1 and S2 for all objectives, if S1 is better than S2, then we say that S1 dominates S2, and if S1 is not dominated by other solutions, then S1 is called a nondominated solution, also called a Pareto optimal solution.


## 2.3. Related Definitions of Multiobjective Optimization Problem

The following introduces some basic concepts in the multiobjective optimization problem [22, 23].Definition 1.
(feasible domain). In the decision space, the feasible domain is represented byX, and its expression is(2)X=x∈Rn|gi≤0,i=1,2,…,p;hj=0,j=1,2,…,q.Definition 2.
(feasible solution). For the pointx in the decision space, if x∈X, x is the feasible solution.Definition 3.
(Pareto dominates). For any vectorx1,x2, if and only if(3)∀i∈Rn,fix1≤fix2,∧∃j∈Rn,fjx1<fjx2,then x1 dominates x2, so x1≺x2.Definition 4.
(Pareto optimal solution). The Pareto optimal solution is also known as the nondominant solution. For a solutionx∗ in the feasible domain, if x∗ is not dominated by any other solution in the feasible domain, x∗ is called the Pareto optimal solution, and its definition is as follows:(4)¬∃x∈X:x≺x∗.Definition 5.
(Pareto optimal solution set, PS). The set of all nondominated solutions in the feasible domain is called the Pareto optimal solution set and can be defined as follows:(5)PS=x∗|¬∃x∈X:x≺x∗.Definition 6.
(Pareto optimal frontier, PF). The target vector set corresponding to the Pareto optimal solution set is the Pareto optimal frontier, also known as the Pareto optimal front end or Pareto equilibrium surface, which is defined as follows:(6)PF=Fx∗=f1x∗,f2x∗,…,fdx∗|x∗∈PS.


## 3. Particle Swarm Optimization Algorithm

Particle swarm optimization (PSO) is a heuristic swarm intelligence algorithm that solves optimization problems by imitating the swarm behavior of birds [24]. The algorithm has high stability and good adaptability. It is an intelligent global optimization algorithm that has attracted the attention of many scholars in recent years. Each particle in the PSO algorithm is equivalent to a bird in the population, so each particle has its own speed and position. Through self-learning and social learning, the particles move in the solution space to obtain the global optimum solution. Assuming that the population size of the particle is N and the dimension of the space is D, the update formula for the velocity and position of the particle are as follows:(7)νijt+1=ωνijt+c1r1kpbestijt−xit+c2r2kgbestjt−xit,xijt+1=xijt+νijt+1,where i (i≤N) is the i−th particle, j (j≤D) is the velocity (position) of the j−th dimension, t is the current iteration number, velocity and position are limited in a certain range, ω is the inertia weight, c1 and c2 are two positive constants in 0,4 which represent the learning factors, r1 and r2 are two random numbers in 0,1, pbestij represents a coordinate that defines a particle pbesti in the j−th dimension of the individual best position, and gbestj represents a coordinate that defines the global best position in the j−th dimension of the population. ω is determined by the following formula:(8)ω=ωmax−t∗ωmax−ωminTmax,where ωmax,ωmin are the maximum and minimum values of inertia weight. Generally the value is set as 0.9 and 0.4, t represents the number of current iterations, and Tmax is the maximum number of iterations.


## 4. An Improved Multiobjective Particle Swarm Optimization Algorithm

### 4.1. Competition Mechanism

The competition mechanism of this paper can quickly search the disposal solution set and build the external archives set in the MOPSO algorithm. First, select a particlex from the population s, and generally select the first particle in the population. Then, let s=s−x and compare each particle in the x population on the basis of the objective function value of Pareto dominance relations. If x≺y, particle y will be removed from the population s; otherwise, let x=y. Finally, let n=n∪x until s=ϕ , at which time n is the nondominant solution set to be solved. This method is also adopted when the nondominant solution set enters the external archive set. As more and more particles are removed, the number of times the algorithm runs is less, which can effectively reduce complexity of the algorithm and improve search speed of the algorithm. The pseudocode of the algorithm is shown in Table 1.Table 1
Constructing nondominated setn by competition mechanism.
Call function:n = find_nondominated_set (s)(1)n = ∅(2) while (|s| > 1){(3)x = first (s); //Assign the first particle(4)        //in populations to x(5)s = s − {x};(6) for (eachy∈s) {(7) if (x dominates y) s = s − {y};(8) else {x} = {y};(9) } //end for(10)n = n + {x};(11) }//end while


### 4.2. External Archives Maintenance Strategy

In solving a MOPSO algorithm, each iteration produces a set of Pareto solutions. Therefore, external archives are needed to store the nondominated solutions produced by each iteration, and the solution set forms the Pareto front [25]. After each iteration, the Pareto front is also updated. However, as the number of iterations increases, the size of external archives will increase, and the complexity of the algorithm will also greatly increase. Therefore, it is necessary to limit the size of external archives to reduce the complexity of the algorithm. This paper uses the descending order of crowded distance to limit the size of the external archives set.


### 4.3. Selecting pbest and gbest

In a MOPSO algorithm, selecting the gbest is very important, which directly affects the convergence speed and capabilities. In IMOPSO, the size of the particle population is fixed, and the particles will not be deleted from the population, but the position of the particles in the population needs to be adjusted to update pbest and gbest. In multiobjective conditions, gbest normally exists in a group of noninferior solutions and is not a single gbest position. When gbest and pbest are nondominated, each particle may have more than one pbest. Therefore, it is necessary to choose pbest and gbest by the appropriate method.


#### 4.3.1. Selecting pbest

The specific process is as follows:pbestit is used to record the individual position and save the nondominated solution of the particles in the evolution process. The updating formula of pbestit for the t−th generation particles is as follows:(9)pbestitpbestit,iffpbestit≺fxit,xit,iffpbestit≻fxit,randselectpbestit,xit,ifotherwise,where pbestit is the optimal position of the i−th particle in the previous t generation and xit is the position of the i−th particle in the t−th generation.


#### 4.3.2. Selecting gbest

The specific process is as follows: In the selection process of gbest, this paper adopts the Pareto principle, also known as the asymmetry principle or the 80/20 law; that is, 80% of the results in practical issues are produced by 20% of key factors [26]. Therefore, in this study, the global optimal value is randomly selected from the top 20% nondominated solutions of the external archive.


### 4.4. Parameter Improvement Strategy

Inertia weightω and the learning factors c1 and c2 in the PSO algorithm have a considerable influence on the searching ability of the population in the target region. However, the traditional linear adjustment strategies for these two types of parameters cannot effectively reflect the search process of the algorithm. Therefore, this article adopts a nonlinear dynamic adjustment strategy to more accurately reflect the search process of the algorithm and to more effectively balance the exploitation and exploration.


#### 4.4.1. Inertia Weight

The inertia weightω determines the influence of the velocity of the previous generation particle on the current velocity. The appropriate adjustment rule of ω can effectively balance the exploitation and exploration of the algorithm. If ω=0, the particle speed is determined by the current position, and the particle has no inheritance to the speed of the previous generation, so the algorithm is easy to fall into the local optimum. If ω≠0, the larger inertial weight can strengthen the global search capability, and the smaller inertial weight can strengthen the local search capability. In order to effectively balance the exploration and exploitation and improve the search performance of the algorithm, this paper adopts a new strategy for adjustment; the updated formula is as follows:(10)ωt=ωmax∗11+t/p1p2+ωmin,where ωt changes with the number of iterations. In this paper, p1 is one-third of the maximum number of iterations. p2 is 10. ωmax is 0.5. ωmin is 0.4. t is the number of current iterations. Figure 1(a) shows the change curve of the new inertia weight.Figure 1
The curves of (a) new inertia weight and (b) learning factor.
(a)(b)


#### 4.4.2. Learning Factor

The learning factor determines the influence of self-learning ability and social-learning ability on particle motion during the search process, which reflects the state of information exchange between particles. In recent years, many scholars have revised the learning factors. Some scholars have proposed asynchronous learning factors; that is, in the initial search phase of the algorithm, the particles have greater self-learning ability and smaller social-learning ability; in the later phase of the algorithm search, the asynchronous learning factors can enhance the ability of particles to move to the global optimal position and obtain high-quality particles, which makes the algorithm have a higher probability of converging to the global optimal solution. For the learning factors, the nonlinear function of the new inertia weight mentioned above [27] is adopted in this paper, and its updated formula is as follows:(11)c1t=1.167×ωt2−0.1167×ωt+0.66,c2t=3−c1t,where, through the above formula, it is known that the learning factor also adjusts dynamically with the adjustment of inertia weight. Figure 1(b) shows the change curve of the learning factor.


### 4.5. Time-Varying Gaussian Mutation

In solving the problem of multiobjective optimization, the fast convergence of the PSO may not be advantageous. Instead, it may make the population prematurely gather around certain particles or a certain area and lose diversity, which is easy to make the algorithm show premature phenomena. Therefore, to enhance the population diversity and avoid the algorithm falling into local optimum, based on literature [13], this paper designs a kind of time-varying particle mutation to produce new solutions through mutation operators to external archives. The perturbation formula is as follows:(12)xit+1xit+rg∗mutrange,Pm≥rand,xit,Pm<rand,Pm=1−ttmax5/Mr,mutrange=ubj−lbj∗Pm,where Pm is the mutation probability. rg is subject to a Gaussian distribution with a mean of 0 and a variance of 1. mutrange represents the scope of action of variation. ubj and lbj are the upper bound and lower bound of decision variables of the j−th dimension, respectively. Mr is the mutation parameter, and the value in this paper is 0.5.


### 4.6. The Specific Steps of the IMOPSO Algorithm

Step 1.
Initialize the group’s position and speed. Set the iteration times, population size, and algorithm parameters.Step 2.
Calculate the fitness value of each particle. Generate a nondominating solution set according to dominating relation.Step 3.
Update the external archive set.Step 4.
Arrange each particle of the external archival set in descending order according to the crowding distance and determine whether the set size number is exceeded. If it is exceeded, clear the nondominant solution beyond the size.Step 5.
Update the individual optimal position; if it is the first generation, directly select the initial position of each particle as the individual optimal value; otherwise, update according to the Pareto dominance relationship.Step 6.
Randomly select the global optimal location from the external archive set ranked in the top 20% nondominant solution.Step 7.
Update the speed formula. If the particle velocityνi>νmax, then let νi=νmax; if the particle velocity νi<νmin, then let νi=νmin.Step 8.
Update the position of the next generation of each particle and carry out time-varying Gaussian variation on the particles in each external archive according to the probabilityP to avoid premature precocity of the algorithm.Step 9.
Determine whether the condition (the maximum number of iterations) is met, and if so, end the loop. Otherwise, return to Step 2 to continue the iteration.


## 4.1. Competition Mechanism

The competition mechanism of this paper can quickly search the disposal solution set and build the external archives set in the MOPSO algorithm. First, select a particlex from the population s, and generally select the first particle in the population. Then, let s=s−x and compare each particle in the x population on the basis of the objective function value of Pareto dominance relations. If x≺y, particle y will be removed from the population s; otherwise, let x=y. Finally, let n=n∪x until s=ϕ , at which time n is the nondominant solution set to be solved. This method is also adopted when the nondominant solution set enters the external archive set. As more and more particles are removed, the number of times the algorithm runs is less, which can effectively reduce complexity of the algorithm and improve search speed of the algorithm. The pseudocode of the algorithm is shown in Table 1.Table 1
Constructing nondominated setn by competition mechanism.
Call function:n = find_nondominated_set (s)(1)n = ∅(2) while (|s| > 1){(3)x = first (s); //Assign the first particle(4)        //in populations to x(5)s = s − {x};(6) for (eachy∈s) {(7) if (x dominates y) s = s − {y};(8) else {x} = {y};(9) } //end for(10)n = n + {x};(11) }//end while


## 4.2. External Archives Maintenance Strategy

In solving a MOPSO algorithm, each iteration produces a set of Pareto solutions. Therefore, external archives are needed to store the nondominated solutions produced by each iteration, and the solution set forms the Pareto front [25]. After each iteration, the Pareto front is also updated. However, as the number of iterations increases, the size of external archives will increase, and the complexity of the algorithm will also greatly increase. Therefore, it is necessary to limit the size of external archives to reduce the complexity of the algorithm. This paper uses the descending order of crowded distance to limit the size of the external archives set.


## 4.3. Selecting pbest and gbest

In a MOPSO algorithm, selecting the gbest is very important, which directly affects the convergence speed and capabilities. In IMOPSO, the size of the particle population is fixed, and the particles will not be deleted from the population, but the position of the particles in the population needs to be adjusted to update pbest and gbest. In multiobjective conditions, gbest normally exists in a group of noninferior solutions and is not a single gbest position. When gbest and pbest are nondominated, each particle may have more than one pbest. Therefore, it is necessary to choose pbest and gbest by the appropriate method.


### 4.3.1. Selecting pbest

The specific process is as follows:pbestit is used to record the individual position and save the nondominated solution of the particles in the evolution process. The updating formula of pbestit for the t−th generation particles is as follows:(9)pbestitpbestit,iffpbestit≺fxit,xit,iffpbestit≻fxit,randselectpbestit,xit,ifotherwise,where pbestit is the optimal position of the i−th particle in the previous t generation and xit is the position of the i−th particle in the t−th generation.


### 4.3.2. Selecting gbest

The specific process is as follows: In the selection process of gbest, this paper adopts the Pareto principle, also known as the asymmetry principle or the 80/20 law; that is, 80% of the results in practical issues are produced by 20% of key factors [26]. Therefore, in this study, the global optimal value is randomly selected from the top 20% nondominated solutions of the external archive.


## 4.3.1. Selecting pbest

The specific process is as follows:pbestit is used to record the individual position and save the nondominated solution of the particles in the evolution process. The updating formula of pbestit for the t−th generation particles is as follows:(9)pbestitpbestit,iffpbestit≺fxit,xit,iffpbestit≻fxit,randselectpbestit,xit,ifotherwise,where pbestit is the optimal position of the i−th particle in the previous t generation and xit is the position of the i−th particle in the t−th generation.


## 4.3.2. Selecting gbest

The specific process is as follows: In the selection process of gbest, this paper adopts the Pareto principle, also known as the asymmetry principle or the 80/20 law; that is, 80% of the results in practical issues are produced by 20% of key factors [26]. Therefore, in this study, the global optimal value is randomly selected from the top 20% nondominated solutions of the external archive.


## 4.4. Parameter Improvement Strategy

Inertia weightω and the learning factors c1 and c2 in the PSO algorithm have a considerable influence on the searching ability of the population in the target region. However, the traditional linear adjustment strategies for these two types of parameters cannot effectively reflect the search process of the algorithm. Therefore, this article adopts a nonlinear dynamic adjustment strategy to more accurately reflect the search process of the algorithm and to more effectively balance the exploitation and exploration.


### 4.4.1. Inertia Weight

The inertia weightω determines the influence of the velocity of the previous generation particle on the current velocity. The appropriate adjustment rule of ω can effectively balance the exploitation and exploration of the algorithm. If ω=0, the particle speed is determined by the current position, and the particle has no inheritance to the speed of the previous generation, so the algorithm is easy to fall into the local optimum. If ω≠0, the larger inertial weight can strengthen the global search capability, and the smaller inertial weight can strengthen the local search capability. In order to effectively balance the exploration and exploitation and improve the search performance of the algorithm, this paper adopts a new strategy for adjustment; the updated formula is as follows:(10)ωt=ωmax∗11+t/p1p2+ωmin,where ωt changes with the number of iterations. In this paper, p1 is one-third of the maximum number of iterations. p2 is 10. ωmax is 0.5. ωmin is 0.4. t is the number of current iterations. Figure 1(a) shows the change curve of the new inertia weight.Figure 1
The curves of (a) new inertia weight and (b) learning factor.
(a)(b)


### 4.4.2. Learning Factor

The learning factor determines the influence of self-learning ability and social-learning ability on particle motion during the search process, which reflects the state of information exchange between particles. In recent years, many scholars have revised the learning factors. Some scholars have proposed asynchronous learning factors; that is, in the initial search phase of the algorithm, the particles have greater self-learning ability and smaller social-learning ability; in the later phase of the algorithm search, the asynchronous learning factors can enhance the ability of particles to move to the global optimal position and obtain high-quality particles, which makes the algorithm have a higher probability of converging to the global optimal solution. For the learning factors, the nonlinear function of the new inertia weight mentioned above [27] is adopted in this paper, and its updated formula is as follows:(11)c1t=1.167×ωt2−0.1167×ωt+0.66,c2t=3−c1t,where, through the above formula, it is known that the learning factor also adjusts dynamically with the adjustment of inertia weight. Figure 1(b) shows the change curve of the learning factor.


## 4.4.1. Inertia Weight

The inertia weightω determines the influence of the velocity of the previous generation particle on the current velocity. The appropriate adjustment rule of ω can effectively balance the exploitation and exploration of the algorithm. If ω=0, the particle speed is determined by the current position, and the particle has no inheritance to the speed of the previous generation, so the algorithm is easy to fall into the local optimum. If ω≠0, the larger inertial weight can strengthen the global search capability, and the smaller inertial weight can strengthen the local search capability. In order to effectively balance the exploration and exploitation and improve the search performance of the algorithm, this paper adopts a new strategy for adjustment; the updated formula is as follows:(10)ωt=ωmax∗11+t/p1p2+ωmin,where ωt changes with the number of iterations. In this paper, p1 is one-third of the maximum number of iterations. p2 is 10. ωmax is 0.5. ωmin is 0.4. t is the number of current iterations. Figure 1(a) shows the change curve of the new inertia weight.Figure 1
The curves of (a) new inertia weight and (b) learning factor.
(a)(b)


## 4.4.2. Learning Factor

The learning factor determines the influence of self-learning ability and social-learning ability on particle motion during the search process, which reflects the state of information exchange between particles. In recent years, many scholars have revised the learning factors. Some scholars have proposed asynchronous learning factors; that is, in the initial search phase of the algorithm, the particles have greater self-learning ability and smaller social-learning ability; in the later phase of the algorithm search, the asynchronous learning factors can enhance the ability of particles to move to the global optimal position and obtain high-quality particles, which makes the algorithm have a higher probability of converging to the global optimal solution. For the learning factors, the nonlinear function of the new inertia weight mentioned above [27] is adopted in this paper, and its updated formula is as follows:(11)c1t=1.167×ωt2−0.1167×ωt+0.66,c2t=3−c1t,where, through the above formula, it is known that the learning factor also adjusts dynamically with the adjustment of inertia weight. Figure 1(b) shows the change curve of the learning factor.


## 4.5. Time-Varying Gaussian Mutation

In solving the problem of multiobjective optimization, the fast convergence of the PSO may not be advantageous. Instead, it may make the population prematurely gather around certain particles or a certain area and lose diversity, which is easy to make the algorithm show premature phenomena. Therefore, to enhance the population diversity and avoid the algorithm falling into local optimum, based on literature [13], this paper designs a kind of time-varying particle mutation to produce new solutions through mutation operators to external archives. The perturbation formula is as follows:(12)xit+1xit+rg∗mutrange,Pm≥rand,xit,Pm<rand,Pm=1−ttmax5/Mr,mutrange=ubj−lbj∗Pm,where Pm is the mutation probability. rg is subject to a Gaussian distribution with a mean of 0 and a variance of 1. mutrange represents the scope of action of variation. ubj and lbj are the upper bound and lower bound of decision variables of the j−th dimension, respectively. Mr is the mutation parameter, and the value in this paper is 0.5.


## 4.6. The Specific Steps of the IMOPSO Algorithm

Step 1.
Initialize the group’s position and speed. Set the iteration times, population size, and algorithm parameters.Step 2.
Calculate the fitness value of each particle. Generate a nondominating solution set according to dominating relation.Step 3.
Update the external archive set.Step 4.
Arrange each particle of the external archival set in descending order according to the crowding distance and determine whether the set size number is exceeded. If it is exceeded, clear the nondominant solution beyond the size.Step 5.
Update the individual optimal position; if it is the first generation, directly select the initial position of each particle as the individual optimal value; otherwise, update according to the Pareto dominance relationship.Step 6.
Randomly select the global optimal location from the external archive set ranked in the top 20% nondominant solution.Step 7.
Update the speed formula. If the particle velocityνi>νmax, then let νi=νmax; if the particle velocity νi<νmin, then let νi=νmin.Step 8.
Update the position of the next generation of each particle and carry out time-varying Gaussian variation on the particles in each external archive according to the probabilityP to avoid premature precocity of the algorithm.Step 9.
Determine whether the condition (the maximum number of iterations) is met, and if so, end the loop. Otherwise, return to Step 2 to continue the iteration.


## 5. Algorithm Performance Evaluation Indexes and Test Functions

### 5.1. Performance Measurement

The quality evaluation mainly focuses on the distance between the solution produced by the algorithms and the Pareto optimal solution for MOPs and the extent covered by the solution produced by the algorithms. In this paper, two performance indicators are adopted.


#### 5.1.1. Generational Distance

Generational Distance (GD) is used to evaluate the convergence performance of multiobjective algorithms. It is used to calculate the average minimum Euclidean distance from each point in the solution setn to the reference set n∗. The calculation formula is defined as follows [28]:(13)GDn,n∗=∑y∈nminx∈n∗disx,y2n,where n is the solution set obtained by the algorithm, n∗ is a set of uniformly distributed reference points sampled from true Pareto front (PF), and disx,y is the Euclidean distance between the point y in the solution set n and the point x in the reference set n∗. The smaller GD is, the closer the Pareto optimal solution set obtained by the algorithm is to the true Pareto front, and the better convergence of the algorithm is. The ideal value of GD is 0; that is, the Pareto optimal solution obtained by the algorithm is on the true Pareto front.


#### 5.1.2. Spacing

Spacing (SP) is a parameter used to evaluate the performance of the Pareto optimal solution obtained by the algorithm. The mathematical expression is [29, 30](14)SP=1n−1∑i=1nd¯−di2,where n represents the number of nondominated solutions found; d¯ is the average of all di, and di is calculated as follows:(15)di=minj=1,…,n∑m=1kfmxi−fmxj,where k represents the number of objective functions. The smaller SP is, the better the distribution of the solution obtained by the algorithm is. The ideal value of SP is 0; that is, the Pareto optimal solution obtained by the algorithm is evenly distributed in the target space.


#### 5.1.3. Inverted Generational Distance

GD can only evaluate the convergence of the algorithm. In order to further evaluate the comprehensive performance of the algorithm, Inverted Generational Distance (IGD) was proposed. IGD represents the average value of the Euclidean distance between each reference point in the reference setn∗ and the closest solution to n. The closer the IGD value is to 0, the better the overall performance of the algorithm is. The calculation formula of IGD is as follows:(16)IGDn,n∗=∑x∈n∗miny∈n∗disx,y2n∗,where n is the solution set obtained by the algorithm, n∗ is a set of uniformly distributed reference points sampled from true PF, and disx,y is the Euclidean distance between the point x in the reference set n∗ and the point y in the solution set n.


### 5.2. Test Functions

To test the performance of the IMOPSO algorithm, this paper selects the classic multiobjective test functions for the simulation test [27, 31]. The test functions are expressed in Tables 2 and 3.Table 2
Test functions: ZDT.
NameObjective functionsDimensionVariable boundsZDT1minf1x=x1minf2x=gx1−x1/gxgx=1+9∑i=2nxi/n−130xi∈0,1ZDT2minf1x=x1minf2x=gx1−x1/gx2gx=1+9∑i=2nxi/n−130xi∈0,1ZDT3minf1x=x1minf2x=gx1−x1/gx−x1/gxsin10πx1gx=1+9∑i=2nxi/n−130xi∈0,1ZDT4minf1x=x1minf2x=gx1−x1/gxgx=1+10n−1+∑i=2nxi2−10cos10πxi10xi∈0,1xi∈−5,5ZDT6minf1x=x1minf2x=gx1−x1/gx2gx=1+10n−1+∑i=2nxi2−10cos10πxi10x1∈0,1Table 3
Test functions: DTLZ.
NameObjective functionsDimensionVariable boundsDTLZ1f1x=0.51+gxM∏i=1m−1xifm=2;M−1x=0.51+gxM1−xM−m+1∏i=1m−1xifMx=1+gxMsinxiπ/2gx=100xM+∑xi∈MX−0.52−cos20πx−0.530xi∈0,1DTLZ2minf1x=0.51+gxM∏i=1M−1cosxiπ/2minfm=2;M−1x=0.51+gxMsinxM−m+1π/2∏i=1M−1cosxiπ/2minfMx=1+gxMsinxiπ/2gxM=∑xi∈Mxi−0.5230xi=∈0,1DTLZ3f1=0.51+gxM∏i=1M−1cosfracxi2fm=2;M−1x=0.51+gxMsinπxM−m+1/2fMx=100xM+∑xi∈Mx−0.52−cos20πxi−0.5gx=∑xi∈Mx−0.5230xi∈0,1DTLZ4f1x=0.51+gxM∏i=1M−1cosxiπ/2fm=2;M−1x=0.51+gxMsinπxM−m+1/2∏i=1M−1cosπxi/2fMx=1+gxMsinx1aπgx=∑xi∈MX−0.52,a=10030xi∈0,1DTLZ5f1x=0.51+gxM∏i=1M−1cosπ1+2gxMx1/22+1+gxMfm=2;M−1x=0.51+gxM·sinπxM−m+1/2∏i=1M−mcosπ1+2gxMx1/22+1+gxMfMx=1+gxMsinπxi/2gx=∑xi∈Mx−0.5230xi∈0,1DTLZ6f1x=0.51+gxM∏i=1M−1cosπ1+2gxMxi/22+1+gxM,fm=2;M−1x=0.51+gxM·sinπxM−m+1/2∏i=1M−mcosπ1+2gxMxi/22+1+gxMfMx=1+gxMsinπxi/2gx=∑xi∈Mxi30xi∈0,1DTLZ7f1x=x1f2x=x2f3x=1+gxkhf1,f2,ggxk=1+9/xk∑xi∈xkxihf1,f2,g=k−∑i=1k−1fi/1+g1+sin3πfi30x∈0,1


## 5.1. Performance Measurement

The quality evaluation mainly focuses on the distance between the solution produced by the algorithms and the Pareto optimal solution for MOPs and the extent covered by the solution produced by the algorithms. In this paper, two performance indicators are adopted.


### 5.1.1. Generational Distance

Generational Distance (GD) is used to evaluate the convergence performance of multiobjective algorithms. It is used to calculate the average minimum Euclidean distance from each point in the solution setn to the reference set n∗. The calculation formula is defined as follows [28]:(13)GDn,n∗=∑y∈nminx∈n∗disx,y2n,where n is the solution set obtained by the algorithm, n∗ is a set of uniformly distributed reference points sampled from true Pareto front (PF), and disx,y is the Euclidean distance between the point y in the solution set n and the point x in the reference set n∗. The smaller GD is, the closer the Pareto optimal solution set obtained by the algorithm is to the true Pareto front, and the better convergence of the algorithm is. The ideal value of GD is 0; that is, the Pareto optimal solution obtained by the algorithm is on the true Pareto front.


### 5.1.2. Spacing

Spacing (SP) is a parameter used to evaluate the performance of the Pareto optimal solution obtained by the algorithm. The mathematical expression is [29, 30](14)SP=1n−1∑i=1nd¯−di2,where n represents the number of nondominated solutions found; d¯ is the average of all di, and di is calculated as follows:(15)di=minj=1,…,n∑m=1kfmxi−fmxj,where k represents the number of objective functions. The smaller SP is, the better the distribution of the solution obtained by the algorithm is. The ideal value of SP is 0; that is, the Pareto optimal solution obtained by the algorithm is evenly distributed in the target space.


### 5.1.3. Inverted Generational Distance

GD can only evaluate the convergence of the algorithm. In order to further evaluate the comprehensive performance of the algorithm, Inverted Generational Distance (IGD) was proposed. IGD represents the average value of the Euclidean distance between each reference point in the reference setn∗ and the closest solution to n. The closer the IGD value is to 0, the better the overall performance of the algorithm is. The calculation formula of IGD is as follows:(16)IGDn,n∗=∑x∈n∗miny∈n∗disx,y2n∗,where n is the solution set obtained by the algorithm, n∗ is a set of uniformly distributed reference points sampled from true PF, and disx,y is the Euclidean distance between the point x in the reference set n∗ and the point y in the solution set n.


## 5.1.1. Generational Distance

Generational Distance (GD) is used to evaluate the convergence performance of multiobjective algorithms. It is used to calculate the average minimum Euclidean distance from each point in the solution setn to the reference set n∗. The calculation formula is defined as follows [28]:(13)GDn,n∗=∑y∈nminx∈n∗disx,y2n,where n is the solution set obtained by the algorithm, n∗ is a set of uniformly distributed reference points sampled from true Pareto front (PF), and disx,y is the Euclidean distance between the point y in the solution set n and the point x in the reference set n∗. The smaller GD is, the closer the Pareto optimal solution set obtained by the algorithm is to the true Pareto front, and the better convergence of the algorithm is. The ideal value of GD is 0; that is, the Pareto optimal solution obtained by the algorithm is on the true Pareto front.


## 5.1.2. Spacing

Spacing (SP) is a parameter used to evaluate the performance of the Pareto optimal solution obtained by the algorithm. The mathematical expression is [29, 30](14)SP=1n−1∑i=1nd¯−di2,where n represents the number of nondominated solutions found; d¯ is the average of all di, and di is calculated as follows:(15)di=minj=1,…,n∑m=1kfmxi−fmxj,where k represents the number of objective functions. The smaller SP is, the better the distribution of the solution obtained by the algorithm is. The ideal value of SP is 0; that is, the Pareto optimal solution obtained by the algorithm is evenly distributed in the target space.


## 5.1.3. Inverted Generational Distance

GD can only evaluate the convergence of the algorithm. In order to further evaluate the comprehensive performance of the algorithm, Inverted Generational Distance (IGD) was proposed. IGD represents the average value of the Euclidean distance between each reference point in the reference setn∗ and the closest solution to n. The closer the IGD value is to 0, the better the overall performance of the algorithm is. The calculation formula of IGD is as follows:(16)IGDn,n∗=∑x∈n∗miny∈n∗disx,y2n∗,where n is the solution set obtained by the algorithm, n∗ is a set of uniformly distributed reference points sampled from true PF, and disx,y is the Euclidean distance between the point x in the reference set n∗ and the point y in the solution set n.


## 5.2. Test Functions

To test the performance of the IMOPSO algorithm, this paper selects the classic multiobjective test functions for the simulation test [27, 31]. The test functions are expressed in Tables 2 and 3.Table 2
Test functions: ZDT.
NameObjective functionsDimensionVariable boundsZDT1minf1x=x1minf2x=gx1−x1/gxgx=1+9∑i=2nxi/n−130xi∈0,1ZDT2minf1x=x1minf2x=gx1−x1/gx2gx=1+9∑i=2nxi/n−130xi∈0,1ZDT3minf1x=x1minf2x=gx1−x1/gx−x1/gxsin10πx1gx=1+9∑i=2nxi/n−130xi∈0,1ZDT4minf1x=x1minf2x=gx1−x1/gxgx=1+10n−1+∑i=2nxi2−10cos10πxi10xi∈0,1xi∈−5,5ZDT6minf1x=x1minf2x=gx1−x1/gx2gx=1+10n−1+∑i=2nxi2−10cos10πxi10x1∈0,1Table 3
Test functions: DTLZ.
NameObjective functionsDimensionVariable boundsDTLZ1f1x=0.51+gxM∏i=1m−1xifm=2;M−1x=0.51+gxM1−xM−m+1∏i=1m−1xifMx=1+gxMsinxiπ/2gx=100xM+∑xi∈MX−0.52−cos20πx−0.530xi∈0,1DTLZ2minf1x=0.51+gxM∏i=1M−1cosxiπ/2minfm=2;M−1x=0.51+gxMsinxM−m+1π/2∏i=1M−1cosxiπ/2minfMx=1+gxMsinxiπ/2gxM=∑xi∈Mxi−0.5230xi=∈0,1DTLZ3f1=0.51+gxM∏i=1M−1cosfracxi2fm=2;M−1x=0.51+gxMsinπxM−m+1/2fMx=100xM+∑xi∈Mx−0.52−cos20πxi−0.5gx=∑xi∈Mx−0.5230xi∈0,1DTLZ4f1x=0.51+gxM∏i=1M−1cosxiπ/2fm=2;M−1x=0.51+gxMsinπxM−m+1/2∏i=1M−1cosπxi/2fMx=1+gxMsinx1aπgx=∑xi∈MX−0.52,a=10030xi∈0,1DTLZ5f1x=0.51+gxM∏i=1M−1cosπ1+2gxMx1/22+1+gxMfm=2;M−1x=0.51+gxM·sinπxM−m+1/2∏i=1M−mcosπ1+2gxMx1/22+1+gxMfMx=1+gxMsinπxi/2gx=∑xi∈Mx−0.5230xi∈0,1DTLZ6f1x=0.51+gxM∏i=1M−1cosπ1+2gxMxi/22+1+gxM,fm=2;M−1x=0.51+gxM·sinπxM−m+1/2∏i=1M−mcosπ1+2gxMxi/22+1+gxMfMx=1+gxMsinπxi/2gx=∑xi∈Mxi30xi∈0,1DTLZ7f1x=x1f2x=x2f3x=1+gxkhf1,f2,ggxk=1+9/xk∑xi∈xkxihf1,f2,g=k−∑i=1k−1fi/1+g1+sin3πfi30x∈0,1


## 6. Experimental Analysis and Comparison

Through experiments, the IMOPSO algorithm in this paper was compared with the experimental results of the NSGA-II, SPEA-II, MOPSO, and NSGA-III [32] algorithms. The population size was set as 100, the maximum number of iterations was set as 10,000, and the size of the external archive was set as 100. The specific parameters of each algorithm were set as shown in Table 3. Among them, pc is the crossover probability, pm is the mutation probability, mu is the mutation rate, sep is the variable step length, ng is the grid control maximal particle number, cg is the cross probability, aph is the expansion probability, and t1 and t2 are mutation parameters. The NSGA-II, SPEA-II and MOPSO, and NSGA-III and IMOPSO independent algorithms were ran 30 times on each test function. Tables 4–6, respectively, show the convergence, distribution, and comprehensive performance of the algorithms. The test was completed in MATLAB 2018a with a Windows 10 system, and the computer was configured with an Intel Core i73.40GHz processor.Table 4
Parameters setting.
AlgorithmsParametersNSGA-IIpc=0.8, pm=0.3, mu=0.02, sep=0.1SPEA-IIpc=0.5, pm=0.5, ng=7, cg=0.15MOPSOc1=c2=2, w=0.5, wp=0.9, ng=7 ,aph=0.1NSGA-IIIt1=20, t2=20, pc=0.8, mu=0.3IMOPSOwmax=0.5, wmin=0.4Table 5
Statistical results of GD obtained by different methods for different test functions.
Test functionsGDNSGA-IISPEA-IIMOPSONSGA-IIIIMOPSOZDT1Average7.71E − 32.09E − 21.01E − 22.84E − 31.08E − 4Std. dev.7.28E − 42.31E − 29.59E − 45.95E − 41.37E − 5ZDT2Average7.69E − 32.32E − 22.53E − 29.60E − 44.15E − 3Std. dev.8.98E − 31.92E − 34.81E − 34.83E − 46.47E − 3ZDT3Average6.74E − 41.15E − 31.25E − 11.69E − 36.91E − 5Std. dev.1.59E − 45.85E − 42.91E − 24.51E − 43.25E − 5ZDT4Average2.34E − 21.22E + 06.41E − 36.19E − 31.24E − 3Std. dev.2.60E − 21.78E + 07.63E − 31.73E − 31.07E − 3ZDT6Average8.69E − 35.37E − 32.76E − 16.53E − 22.29E − 3Std. dev.3.69E − 35.43E − 33.43E − 13.29E − 21.79E − 3DTLZ1Average4.19E − 21.94E − 17.11E + 04.14E − 22.15E − 2Std. dev.5.18E − 22.78E − 11.92E + 03.36E − 22.09E − 2DTLZ2Average1.31E − 31.47E − 38.90E − 36.08E − 51.68E − 3Std. dev.6.54E − 32.00E − 42.70E − 34.65E − 53.73E − 4DTLZ3Average1.48E + 01.14E + 12.91E + 12.89E + 11.09E + 0Std. dev.3.74E + 08.16E + 01.03E + 15.39E + 13.87E − 1DTLZ4Average1.98E − 32.09E − 33.96E − 25.76E − 41.56E − 3Std. dev.2.30E − 41.36E − 42.59E − 21.39E − 42.13E − 4DTLZ5Average2.70E − 43.19E − 41.58E − 33.25E − 41.88E − 4Std. dev.7.01E − 58.19E − 51.06E − 38.32E − 56.59E − 5DTLZ6Average1.31E − 55.94E − 37.11E − 24.26E − 15.19E − 6Std. dev.4.53E − 59.37E − 24.28E − 31.72E − 22.97E − 7DTLZ7Average7.06E − 36.18E − 35.97E − 16.85E − 33.13E − 3Std. dev.1.56E − 31.60E − 34.03E − 12.24E − 31.06E − 3Table 6
Statistical results of SP obtained by different methods for different test functions.
Test functionsSPNSGA-IISPEA-IIMOPSONSGA-IIIIMOPSOZDT1Average7.33E − 38.97E − 31.28E − 21.08E − 26.47E − 3Std. dev.6.77E − 41.91E − 32.18E − 31.65E − 36.34E − 4ZDT2Average7.85E − 38.99E − 37.88E − 31.19E − 27.79E − 3Std. dev.2.34E − 34.61E − 35.74E − 35.76E − 31.75E − 3ZDT3Average8.16E − 31.05E − 21.37E − 21.06E − 27.101E − 3Std. dev.9.32E − 42.06E − 35.40E − 31.46E − 39.28E − 4ZDT4Average2.76E − 21.36E − 13.70E − 26.09E − 31.81E − 2Std. dev.2.79E − 23.62E − 13.93E − 22.83E − 25.27E − 2ZDT6Average1.73E − 21.61E − 28.82E − 24.23E − 21.19E − 2Std. dev.3.55E − 23.80E − 26.13E − 22.08E − 24.01E − 3DTLZ1Average5.81E − 21.44E + 05.51E + 06.13E − 23.93E − 1Std. dev.4.41E − 23.27E + 01.45E + 03.27E − 24.01E − 1DTLZ2Average6.17E − 25.20E − 26.15E − 25.78E − 25.95E − 2Std. dev.5.23E − 35.61E − 39.19E − 31.79E − 34.90E − 3DTLZ3Average7.71E + 02.56E + 06.63E + 11.70E + 11.32E + 0Std. dev.3.18E + 05.06E + 02.89E + 15.32E − 12.20E + 0DTLZ4Average5.44E − 24.426E − 24.97E − 25.17E − 26.12E − 2Std. dev.1.36E − 24.76E − 33.12E − 22.15E − 24.87E − 3DTLZ5Average1.04E − 21.32E − 21.38E − 21.52E − 29.79E − 3Std. dev.1.18E − 32.12E − 31.06E − 32.91E − 38.04E − 4DTLZ6Average1.13E − 22.793E − 11.38E − 22.19E − 21.22E − 2Std. dev.8.87E − 45.09E − 21.01E − 12.50E − 21.25E − 3DTLZ7Average7.03E − 24.55E − 25.97E − 16.48E − 23.56E − 2Std. dev.6.89E − 31.63E − 22.65E − 28.90E − 32.33E − 3The numerical experiment results are shown in Tables4–6, and the optimal results are marked in bold. Table 4 gives the numerical experiment result of the evaluation index GD, which is the convergence performance of the algorithm; Table 5 gives the numerical experiment result of the evaluation index SP, which is the algorithm distribution index; Table 6 gives the numerical experiment result of the evaluation index IGD, which represents the algorithm comprehensive performance. It can be seen from Table 4 that, compared with the other compared algorithms, in the test functions, ZDT, IMOPSO has obtained the best results 4 times, and NSGA-III has obtained the best results once; in the test functions, DTLZ, IMOPSO has obtained the best results 5 times and NSGA-III has obtained the best results 2 times. It can be concluded that the convergence of IMOPSO is superior to those of the other algorithms used for comparison, and, at the same time, it can be further concluded that the learning factor adjustment rules adopted in this paper have effectively improved the convergence of the IMOPSO. The numerical experiment results of the algorithms distribution are shown in Table 5. In the test functions, ZDT, IMOPSO obtained the best results 5 times; in the test functions, DTLZ, IMPSO obtained the best results 4 times, NSGA-II obtained the best results 2 times, and SPEA-II obtained the best result once. By analyzing the numerical results in Table 5, it can be seen that the distribution of the algorithm can be effectively improved by introducing time-varying Gaussian mutation strategy into IMOPSO. By analyzing Table 6, it can be seen that, in the test functions, ZDT, IMOPSO has obtained the best results 3 times, and NSGA-II and SPEA-II have obtained the best results once, respectively; in the test functions, DTLZ, IMOPSO has obtained the best results 4 times, and NSGA-II and MOPSO obtained the best results 2 times and once, respectively, so, on the whole, IMOPSO has better comprehensive performance (Table 7).Table 7
Statistical results of IGD obtained by different methods for different test functions.
Test functionsIGDNSGA-IISPEA-IIMOPSONSGA-IIIIMOPSOZDT1Average1.24E − 23.82E − 29.60E − 12.84E − 24.87E − 3Std. dev.2.03E − 34.08E − 22.39E − 16.75E − 32.15E − 4ZDT2Average2.02E − 29.61E − 21.87E + 05.26E − 24.54E − 2Std. dev.1.63E − 27.66E − 23.85E − 13.32E − 21.24E − 1ZDT3Average1.10E − 23.70E − 21.047E + 02.34E − 25.46E − 3Std. dev.5.72E − 32.70E − 22.05E − 15.57E − 32.19E − 4ZDT4Average1.85E − 12.21E − 11.90E + 16.04E − 12.66E + 0Std. dev.1.25E − 11.32E − 18.20E + 02.33E − 12.30E + 0ZDT6Average6.14E − 22.68E − 21.30E + 02.82E − 14.29E − 3Std. dev.2.44E − 29.81E − 31.91E + 01.20E − 11.97E − 4DTLZ1Average2.69E − 11.65E + 01.19E + 12.51E − 11.58E − 1Std. dev.3.14E − 11.33E + 04.63E + 02.10E − 11.44E − 1DTLZ2Average6.92E − 26.72E − 25.50E − 21.04E − 17.2E − 2Std. dev.1.76E − 33.72E − 32.46E − 41.02E − 27.2068E − 2DTLZ3Average8.33E + 07.76E + 11.76E + 2.11.25E + 17.42E + 0Std. dev.5.65E + 02.80E + 15.11E + 14.41E + 05.20E + 0DTLZ4Average9.97E − 27.14E − 23.75E − 12.01E − 16.55E − 2Std. dev.1.20E − 11.95E − 31.86E − 12.27E − 11.73E − 3DTLZ5Average5.98E − 31.23E − 21.31E − 21.52E − 26.62E − 3Std. dev.3.12E − 42.24E − 34.28E − 31.71E − 35.53E − 4DTLZ6Average6.39E − 31.76E − 23.178E + 05.56E − 26.87E − 3Std. dev.2.32E − 32.29E − 28.43E − 131.33E − 17.20E − 4DTLZ7Average2.06E − 14.43E + 01.23E − 11.22E − 19.78E − 2Std. dev.2.69E − 12.04E − 11.21E + 07.83E − 21.05E − 2In order to more intuitively show the characteristics of each optimization problem in this article and the status of each algorithm in solving the optimization problems, Figures2–13 show the corresponding relationship between the optimal solutions obtained by each algorithm and the true PF. In Figures 2–6, we can intuitively draw the following conclusions: In the test functions, ZDT, except that the convergence of IMOPSO in ZDT3 is worse than that of NSGA-III, IMOPSO is better than the other compared algorithms in terms of convergence and distribution. In addition, in the test functions, DTLZ, the convergence of IMOPSO in DTLZ2 and DTLZ4 is worse than that of NSGA-III, but it is better than the other compared algorithms in other test functions; in terms of SP, NSGA-II is better than IMOPSO in DTLZ1 and DTLZ6, and SPEA-II is better than IMOPSO in DTLZ6, but, in the remaining DTLZ test functions, the performance of IMOPSO is better than those of the other compared algorithms.Figure 2
ZDT1 test results.
(a)(b)(c)(d)(e)Figure 3
ZDT2 test results.
(a)(b)(c)(d)(e)Figure 4
ZDT3 test results.
(a)(b)(c)(d)(e)Figure 5
ZDT4 test results.
(a)(b)(c)(d)(e)Figure 6
ZDT6 test results.
(a)(b)(c)(d)(e)Figure 7
DTLZ1 test results.
(a)(b)(c)(d)(e)Figure 8
DTLZ2 test results.
(a)(b)(c)(d)(e)Figure 9
DTLZ3 test results.
(a)(b)(c)(d)(e)Figure 10
DTLZ4 test results.
(a)(b)(c)(d)(e)Figure 11
DTLZ5 test results.
(a)(b)(c)(d)(e)Figure 12
DTLZ6 test results.
(a)(b)(c)(d)(e)Figure 13
DTLZ7 test results.
(a)(b)(c)(d)(e)


## 7. Conclusions and Prospects

In order to improve MOPSO prone to premature convergence, poor distribution, and so forth, in solving multiobjective optimization problems, this paper proposes an improved multiobjective particle swarm optimization (IMOPSO) algorithm. In this paper, by introducing a competition mechanism strategy in IMOPSO and dynamically adjusting the inertia weight and learning factor of the algorithm, the convergence performance of the algorithm is effectively improved. In addition, in order to solve the shortcomings of insufficient distribution of the algorithm in the optimization process, time-varying Gaussian mutation is cited in the algorithm to increase the diversity of the algorithm and improve its distribution. Through the analysis of numerical experiment results, in the experimental results of the test functions, ZDT, and the test functions, DTLZ, as a whole, the convergence and distribution of the algorithm are better than those of the other compared algorithms. Finally, in the results of the comparison of the comprehensive performances of the algorithms, the comprehensive performance of the IMOPSO algorithm is also the best among the compared algorithms.It should be pointed out that this article only conducts the numerical experiment analysis of the benchmark function, and its performance in actual multiobjective optimization problems needs further verification and analysis, such as multiobjective job shop scheduling problems and redundancy allocation problem with cold-standby strategy. This will be our next research works.



---

*Source: 5980504-2020-12-01.xml*