Pareto-Ranking Based Quantum-Behaved Particle Swarm Optimization for Multiobjective Optimization

A study on pareto-ranking based quantum-behaved particle swarm optimization (QPSO) for multiobjective optimization problems is presented in this paper. During the iteration, an external repository is maintained to remember the nondominated solutions, from which the global best position is chosen. The comparison between different elitist selection strategies (preference order, sigma value, and random selection) is performed on four benchmark functions and two metrics. The results demonstrate that QPSO with preference order has comparative performance with sigma value according to different number of objectives. Finally, QPSO with sigma value is applied to solve multiobjective flexible job-shop scheduling problems.


Introduction
Most real-world optimization problems have more than one objective, with at least two objectives conflicting with each other.The conflicting objectives lead to a problem where a single solution does not exist.Instead, a set of optimal trade-off solutions exists, which are referred to as the paretooptimal front or pareto front.This kind of optimization problems is referred to as multiobjective optimization problems.
In the past ten years, a wide variety of algorithms have been proposed to address such problems.Deb et al. [1] gave rise to a fast and elitist multiobjective genetic algorithm (NSGA-II) to reduce the computational complexity and adopted an elitism strategy.In [2], Konak et al. gave a tutorial about several variants of genetic algorithms and compared their performance.A particle swarm optimization (PSO) with pareto dominance was proposed in [3], in which a secondary repository and a special mutation operator were incorporated.A survey of the state-of-the-art multiobjective PSOs is presented in [4].From the above literatures, the great difference between genetic algorithm (GA) and PSO exists for solving the multiobjective optimization problems.In GA, a set of particles are needed to be selected into next generation, while, in QPSO, one particle has to be chosen from the nondominated optima as the global best position ().Therefore, different pareto-ranking strategies are required.
Inspired by the quantum mechanics, QPSO was proposed by Sun et al. in 2004 [5], which has been proved to outperform PSO in both convergence rate and global search ability [6,7].Similar to other population-based algorithms (GA, PSO), maintaining a diverse population is an important consideration in multiobjective optimization to ensure solutions uniformly distributed over the pareto front.Therefore, in this paper, QPSO with three pareto-ranking strategies is compared on some benchmark functions and two metrics are used to test the performance of each strategy.
The research on the multiobjective FJSP is not as widely as that on the monoobjective FJSP.Brandimarte [8] was first to apply the decomposition method to FJSP and solved the routing subproblem by using some existing dispatching rules and addressed the scheduling subproblem by using a tabu search algorithm.Davarzani et al. [9] used artificial immune system to solve the multiobjective FJSP and compared it with Approach by Localization [10], AL + CGA [11].Considering the advantages and disadvantages of stochastic algorithms and local search methods, some hybrid approaches were also proposed [12,13].In this paper, QPSO with sigma value is applied to solve multiobjective flexible job-shop scheduling.

Mathematical Problems in Engineering
The remainder of this paper is organized as follows.In Section 2, some basic concepts for multiobjective optimization are provided.Section 3 describes the three variants of pareto-ranking based QPSO in detail.Numerical tests on benchmark functions are presented in Section 4. In Section 5, QPSO with sigma value is applied to solve multiobjective flexible job-shop scheduling.Finally, conclusion is given in Section 6.
Definition 2 (pareto optimal set).One has Definition 3 (pareto front).For a given pareto optimal set  * , the pareto front is defined as Generally, the analytical expression of the line or surface which contain those points does not exist.It is only possible for us to determine the nondominated points and to produce the pareto front.

Pareto-Ranking Based Quantum-Behaved
Particle Swarm Optimization  [14].A PSO system simulates the knowledge evolvement of a social organism, in which particles represent candidate solutions of an optimization problem.The position of each particle is evaluated according to the objective function, and particles share memories of their "best" positions.These memories are used to adjust the particles' own velocities and their subsequent positions.It has already been shown that PSO is comparable in performance with and may be considered as an alternative to GA [15].Ina PSO system with  particles in the -dimensional space, the position, and velocity vector of particle  at the kth iteration is represented as   () = ( 1 (),  2 (), . . .,   ()) and   () = ( 1 (),  2 (), . . .,   ()).The particle moves according to the following equations: where,  = 1, 2, . . ., ,  = 1, 2, . . ., ,  1 and  2 are acceleration coefficient. 1 and  2 are random numbers distributed uniformly in (0, 1). is inertia weight.Vector   = ( 1 ,  2 , . . .,   ) is the previous best position of particle  (), and vector   = ( 1 ,  2 , . . .,   ) is the position of the best particle  in the swarm ().

Quantum-Behaved Particle Swarm Optimization .
The main disadvantage of PSO is that it is not guaranteed to be global convergent [16].Concepts of QPSO were developed to address this disadvantage and first reported by Sun et al. in 2004 [5].Trajectory analysis in [17] demonstrated the fact that the convergence of PSO may be achieved if each particle converges to its local attractor,   = ( 1 ,  2 , . . .,   ), defined as where  ∈ (0, 1).It can be seen that   is a stochastic attractor of particle  that lies in a hyperrectangle with   and   being two ends of its diagonal and moves following   and   .
In quantum world, the velocity of the particle is meaningless, so, in QPSO system, position is the only state to depict the particles, which moves according to the following equation [5]: where  is a random number uniformly distributed in (0, 1) and (), called mean best position, is defined as the mean value of personal best positions of the swarm: The parameter  in (4) is named as contraction-expansion (CE) coefficient, which can be adjusted to control convergence rate.The most commonly used method to control CE is linearly decreasing from  max to  min : where  is the current iteration step,  max is the predefined maximum iteration steps, and  max and  min are the maximum and minimum value of CE.
QPSO does not require velocity vectors for the particles and has fewer parameters to control, making the algorithm easier to implement.Experimental results performed on some well-known benchmark functions show that QPSO has better performance than PSO [5][6][7].
A great deal of efforts has been made to PSO for solving multiobjective optimization problems.A survey of state-ofthe-art works is presented in [4].The most contributive work is proposed by Ceollo et al. [3], in which the objective space is divided into hypercubes using an adaptive grid.The nondominated solutions are distributed into grids and then the grid is chosen according to a roulette wheel selection strategy.However, the grid is problems dependent and computation complicated.While in [18] a new density measure strategy named as sigma method was introduced to multiobjective PSO, so as to choose a  from the nondominated archive, the sigma method can guide particles to the pareto optimal front.However, if the initialized solutions in nondominated archive are bad, the algorithm will fall into local optimum due to the fast convergence.Yang et al. [19] proposed a hybrid method which combined the density information and sigma method so as to achieve a balance between global search and local search.Considering the calculation of the density value and crowding distance is limited by two objectives.When the number of objectives becomes large, the measure of density and crowding distance will become invalid.Therefore, the idea of preference order is introduced [20].

Preference Order.
As we all know, when the optimization problem has two objectives, the pareto optimal solutions can be plotted on a curve.Though the trade-off surface can also be visualized for three objectives, it is not easy to pick the final point.For problems with more than three objectives, it becomes extremely difficult to find the optimized solution through visualization.
The idea of preference order was firstly proposed in 1999 [21].A point  * ∈ Ω is said to be efficiency of order  if ( * ) is not dominated by any of the -element subsets of  Ω .
Consider a minimization example with three objectives in Table 1; combination of all subsets of the objectives is The dominant relation of all points in 2-element subsets is shown in Table 2, from which only Point 3 is nondominated by any points on all 2-element subsets.So Point 3 is said to be efficiency of order 2.
About efficiency of order, [22] gave three claims.
Claim 1.In a space with  objectives, there is no more than one point with efficiency of order  − 1.
Claim 2.  * is efficient of order  + 1, if it is efficient of order .
Claim 3.  * is not efficient of order  − 1, if it is not efficient of order .
The pseudocode of QPSO with preference order is shown below.(3) Repeat (2) until the maximum number of iterations is reached.

Sigma Method.
The sigma method was first proposed in [18] for finding the best local guide for each particle.Each nondominated solution in the external repository is assigned a value   ,  = 1, 2, . . ., ND, which is defined as (for two objectives shown in Figure 1.) Archive member Particle In the general case,  is a vector of (  2 ) elements, and each element is the combination of two coordinates.Take a problem of three objectives; for example, The pseudocode of QPSO with sigma value is shown below.

QPSO with Sigma Method
(1) Initialization: initialize the swarm size , the external archive size ND, and the maximum number of iterations  max , CE .Position of each particle   ,  = 1, 2, . . ., , is randomly initialized.Personal best position of each particle   is equal to   .Construct external nondominated set nds.
(2) Update the following: (a) Identify local best position using sigma method.
(i) Assign the value   to each particle in the external archive.(ii) Calculate the value   for each particle in the swarm.(iii) Calculate the distance between   and   .(iv) The particle  in the archive with   has the minimum distance to   selected as the best local guide for the particle .
(3) Repeat (2) until the maximum number of iterations is reached.

Numerical Experiments
4.1.Test Functions.Seven common test functions [20][21][22][23] are used to test the proposed algorithm (DTLZ1, DTLZ2, DTLZ3, and DTLZ4), as listed in Table 3.We call  = | ⃗   | difficulty factor.Table 4 gives the difficulty factor, number of runs, number of iterations, number of particles, and size of nondominated archive used in our test.The reason why these functions are chosen is due to the following features: (1) The relatively small implementation effort.
(2) The ability to be scaled to any number of objectives and decision variable.
(3) The global pareto front being known analytically.
(4) Convergence and diversity difficulties that can be easily controlled.

Performance Metrics.
Usually, pareto-based multiobjective algorithms consider two aspects: closeness to the global pareto front (distance) and spread along the pareto front (diversity) [2], which are defined as follows.

Generational Distance (GD).
It is a way of estimating how far the elements in the nondominated set are found so far from those in the pareto optimal set, which is defined as where ND is the number of points in nondominated set and   is the Euclidean distance between a point and the nearest member of Pareto optimal set.

Spacing (SP).
It is a metric measuring the distance variance of neighboring points in the nondominated set, which is defined as where , ,  = 1, 2, . . ., ND,  is the mean of all   , and  is the number of objectives.5 and 6 show the mean values of generational distance and spacing for three methods (QPSO + PO, QPSO + Rand, and QPSO + Sigma) with different number of objectives (2, 3, 4, 5, 6, 7, and 8).It is easy to note that QPSO with preference order has comparative performance with sigma value.When the number of objectives is small (e.g., 2, 3, 4, and 5), QPSO with sigma value performs better.When the number of objectives becomes larger (e.g., 6, 7, and 8), QPSO with preference order obtains solutions closer to the true pareto front, but it is very time-consuming.

Experiments Results. Tables
When  becomes larger, the number of nondominated particles in external archive ND and the number of subsets (   ) ( = 1, 2, . . .,  − 1) will accordingly become large.Therefore, the computing complexity of QPSO with preference order is ( 3 ).Furthermore, in QPSO with sigma value, each particle flies towards its nearest nondominated particle in external archive, which will easily result in a premature convergence and the particles will get trapped in a local optimum.

QPSO for Flexible Job-Shop Scheduling Problems
5.1.Flexible Job-Shop Scheduling Problems.The FJSP could be formulated as follows.There is a set of  jobs and a set of  machines. denotes the set of all machines.Each job  consists of a sequence of   operations.Each operation  , ( = 1, 2, . . ., ;  = 1, 2, . . .,   ) of job  has to be processed on one machine   out of a set of given compatible machines  , for (  ∈  , ,  , ⊆ ).The execution of each operation requires one machine selected from a set of available machines.In addition, the FJSP needs to set each operation's starting and ending time.The FJSP is needed to determine both an assignment and a sequence of the operations on the machines in order to satisfy the given criteria.If there is  , ⊂  for at least one operation, it is partial flexibility FJSP (P-FJSP); while there is  , =  for each operation, it is total flexibility FJSP (T-FJSP) [10,11].In this paper, the following criteria are to be minimized: (1) The maximal completion time of machines, that is, the makespan.
(2) The maximal machine workload, that is, the maximum working time spent on any machine.(3) The total workload of machines, which represents the total working time over all machines.
During the process of solving this problem, the following assumptions are made: (1) Each operation cannot be interrupted during its performance (nonpreemptive condition).
(2) Each machine can perform at most one operation at any time (resource constraint).
(3) The precedence constraints of the operations in a job can be defined for any pair of operations.
(4) Setting up time of machines and move time between operations are negligible.
(5) Machines are independent from each other.
(6) Jobs are independent from each other.
(7) All machines are available at time 0.
(8) Each job has its own release date.
(9) The order of operations for each job is predefined and invariant.The model of the multiobjective FJSP is then given as follows: min min min ∀ (, ) , (ℎ, ) , ; Inequation ( 14) ensures the operation precedence constraint.Inequity (15) ensures that each machine could process only  one operation at each time.Equation ( 16) states that one machine could be selected from the set of available machines for each operation.

Encoding Scheme.
Particle's position representation is the most important task for successful application of PSO or QPSO to solve the FJSP.In this paper, the particle's position is represented by a vector which consists of two parts (take Table 7, e.g.).The first part represents the sequence of the operations ( 1,1 , ).The second part represents the machines on which the operations in Table 7 are processed; take Figure 2 for example, 1000 represents  1,1 is processed on  1 and 0001 represents  1,2 is processed on  4 .The length of the vector is equal to ∑  =1   (1 + ).The advantage of this encoding method is that it is easy to apply crossover and mutation on operation part and machine part separately, while the disadvantage is the storage complexity.

Swarm Initialization.
The initial population has great influence on the performance of the algorithm.Here in this paper, a guided initialization method is used.In order to ensure the diversity of the population, 50% of the particles are initialized randomly, and the others are initialized according to the following rules.

Guided Initialization
Step 0: Initialize.It includes the following: (a) Create an array for a random job order J (e.g., J4 J2 J3 J1).(b) Assign the first operation  ,1 of each job to a suitable machine following the order  4,1  2,1  3,1  1,1 .This machine is selected by calling Step 1(b);  ,1 is assigned directly to the selected machine.If there is no operation assigned to the selected machine,  ,1 starts from time 0. Otherwise, it starts from the stop time of the last operation on the selected machine.
Step 1: Find a Suitable Machine to Process an Operation.It includes the following: Step 2: Find a Suitable Order of Operations on a Machine.from Step 1(a), Ms and  , are identified, and then assign  ,+1 to a machine.Now an operation is needed to find on the waiting list to be processed on Ms.If there is more than one operation on the waiting list of Ms, select one of them depending on one of the following dispatching rules: shortest processing time, longest processing time, or first in first out.

Crossover and Mutation.
For combinational optimization problems, effective information exchanging is able to help find the better solution.The crossover operator on the operation part is shown in Figure 3, in which the offspring C1 and C2 inherit some operations from the parents P1 ({1}) and P2 ({2, 3}) directly, and the other operations are copied from another parent.The crossover operator on the machine part is shown in Figure 4, in which two positions are selected randomly, the same as that in GA.
The mutation of a particle is also divided into two steps; one is on operation part (Figure 5), in which, two positions are randomly chosen, and the operations between them are reversed.
Another mutation is on machine part, but the operation to be mutated cannot be randomly chosen.Because  the makespan and workload are determined by the critical path [13] in the scheduling sequence, the feasible mutation method to reduce the makespan is to choose a critical operation on the critical path, and then try to move it to another machine.Taking Figure 6, for example, the critical path is constructed by { 4,1 ,  2,1 ,  2,2 ,  2,3 }.The operation  4,1 is moved from  1 to  4 , and then the makespan is reduced from 12 to 11.

Experimental Results.
As the benchmark results analysis in Section 4, QPSO with sigma value outperforms QPSO with preference order for smaller number of objectives.In this section, only three objectives are defined (makespan, max machine workload, and total workload), so QPSO with sigma value is adopted.
The parameters are set as population size  = 100; size of the external archive ND = 100; dimensional size equals the total number of operations ∑  =1   ; the maximum number of iterations  max = 2000.
From Table 8, it is easy to notice that, for the problems of 4 × 5, 8 × 8, and 15 × 10, the results obtained by QPSO with sigma dominate other algorithms.For problem Figure 6: The Gantt graph of moving operation.
of 10 × 10, the result obtained by QPSO with sigma is not dominated by other algorithms.Therefore, QPSO with sigma value is competitive for solving multiobjective flexible jobshop scheduling problems.Figures 7-10 show the Gantt graph of the optimized scheduling sequence by QPSO with sigma.

Conclusion and Future Works
In this paper, QPSO with three elitist selection strategies are used to solve multiobjective optimization problems.The results show that QPSO with preference order has comparative performance with sigma value.When the number of objectives is small (e.g., 2, 3, 4, and 5), QPSO with sigma value performs better.However, when the number of objectives becomes larger (e.g., 6, 7, and 8), QPSO with preference order obtains solutions closer to the true pareto front but is very time-consuming.Therefore, in the future work, we will try       to combine these two methods to achieve a balance between global search and local search.
Moreover, we apply QPSO with sigma value to solve multiobjective flexible job-shop scheduling problems, the results of which demonstrate the competitive performance of the algorithm.
(a) Identify global best position using preferenceorder method.(i) Get combination of all subsets for  objectives.(ii) Compute the efficiency order for each nondominated solution.(iii) Sort nondominated solutions in descending order according to their efficiency order.(iv) Get .(b) Update position of each particle according to (3), (4), and (5).(c) Update  according to dominance rule.(d) Update external archive and remove solution with lower efficiency order if archive size ND is exceeded.(e) Update  according to (6).

Figure 1 :
Figure 1: Sigma value for a two-objective space.

Figure 2 :
Figure 2: Encoding of a particle.

Table 3 :
Description of the tested functions.

Table 4 :
Difficulty factor, number of runs, iterations, particles, and external archive.

Table 5 :
Mean value of generational distance.

Table 6 :
Mean value of spacing.

Table 7 :
Example of 3 jobs in 4 machines.
Identify all the machines that can process  ,+1 : ( ,+1 ).(ii) Calculate the waiting time of each machine   ∈ ( ,+1 ): if   does not process any operation at time , it is equal to processing time of  ,+1 on   .Otherwise, if it is processing an operation  , , its waiting time is the total processing times of operations on   's waiting list + (remaining time of  , ) + processing time of  ,+1 on machine   .(iii) The machine with shortest waiting time is chosen to process  ,+1 .If this machine does not process any operation at time , assign  ,+1 to it to be processed.Otherwise, add  ,+1 to the waiting list of this machine.

Table 8 :
Comparison of results in four problems.