An Adaptive Evolutionary Multi-Objective Approach Based on Simulated Annealing

A multi-objective optimization problem can be solved by decomposing it into one or more single objective subproblems in some multi-objective metaheuristic algorithms. Each subproblem corresponds to one weighted aggregation function. For example, MOEA/D is an evolutionary multi-objective optimization (EMO) algorithm that attempts to optimize multiple subproblems simultaneously by evolving a population of solutions. However, the performance of MOEA/D highly depends on the initial setting and diversity of the weight vectors. In this paper, we present an improved version of MOEA/D, called EMOSA, which incorporates an advanced local search technique (simulated annealing) and adapts the search directions (weight vectors) corresponding to various subproblems. In EMOSA, the weight vector of each subproblem is adaptively modified at the lowest temperature in order to diversify the search toward the unexplored parts of the Pareto-optimal front. Our computational results show that EMOSA outperforms six other well established multi-objective metaheuristic algorithms on both the (constrained) multi-objective knapsack problem and the (unconstrained) multi-objective traveling salesman problem. Moreover, the effects of the main algorithmic components and parameter sensitivities on the search performance of EMOSA are experimentally investigated.


Introduction
Many real-world problems can be modelled as combinatorial optimization problems, such as knapsack problem, traveling salesman problem, quadratic assignment problem, flowshop scheduling problem, vehicle routing problem, bin packing problem, and university timetabling problem (Papadimitriou and Steiglitz, 1998).Often, these problems are difficult to tackle due to their huge search space, many local optima, and complex constraints.Many of them are NP-hard, which means that no exact algorithms are known to solve these problems in polynomial computation time.In the last two decades, research on combinatorial optimization problems with multiple objectives has attracted growing interest from researchers (Ehrgott and Gandibleux, 2000;Landa-Silva et al., 2004;Ehrgott, 2005).Due to possible conflicting objectives, optimal solutions for a multi-objective combinatorial optimization (MOCO) problem represent trade-offs among objectives.Such solutions are known as Pareto-optimal solutions.Since the total number of Pareto-optimal solutions could be very large, many practical multi-objective search methods attempt to find a representative and diverse set of Pareto-optimal solutions so that decision-makers can choose a solution based on their preferences.
A number of metaheuristic algorithms including evolutionary algorithms (EA), simulated annealing (SA), tabu search (TS), memetic algorithms (MA) and others (Blum and Roli, 2003;Glover and Kochenberger, 2003;Burke and Kendall, 2005), have been proposed for solving single objective combinatorial optimization problems.Most metaheuristic algorithms try to find global optimal solutions by both diversifying and intensifying the search.Very naturally, these algorithms have also been extended to solve MOCO problems (Gandibleux et al., 2004).Among them, evolutionary multi-objective optimization (EMO) algorithms (Deb, 2001;Tan et al., 2005;Coello Coello et al., 2007) have received much attention due to their ability to find multiple Pareto-optimal solutions in a single run.Pareto dominance and decomposition (weighted aggregation or scalarization) are two major schemes for fitness assignment in EMO algorithms.
Since the 1990s, EMO algorithms based on Pareto dominance have been widely studied.Amongst the most popular methods are PAES (Knowles and Corne, 2000a), NSGA-II (Deb et al., 2002) and SPEA-II (Zitzler et al., 2002).These algorithms have been applied to many real-world and benchmark continuous multi-objective optimization problems.In contrast, EMO algorithms based on decomposition appear to be more successful on tackling MOCO problems.For example, IMMOGLS (Ishibuchi et al., 2003) was applied to the multi-objective flowshop scheduling problem, while MOGLS (Jaszkiewicz, 2002) and MOEA/D (Zhang and Li, 2007) dealt with the multi-objective knapsack problem.In both IMMOGLS and MOGLS, one subproblem with random weight vectors is considered in each generation or iteration.Instead of optimizing one subproblem each time, MOEA/D optimizes multiple subproblems with fixed but uniform weight vectors in parallel in each generation.In most EMO algorithms based on decomposition, local search or local heuristics can be directly used to improve offspring solutions along a certain search direction towards the Pareto-optimal front.These algorithms are also known as multi-objective memetic algorithms (MOMAs) (Knowles and Corne, 2004).
Multi-objective simulated annealing (MOSA) is another class of promising stochastic search techniques for multi-objective optimization (Suman and Kumar, 2006).Several early MOSA-like algorithms (Serafini, 1992;Czyzak and Jaszkiewicz, 1998;Ulungu et al., 1999) define the acceptance criteria for multi-objective local search by means of decomposition.To approximate the whole Pareto-optimal front, multiple weighted aggregation functions with different settings of weights (search directions) are required.Therefore, one of the main challenging tasks in MOSA-like algorithms is to choose appropriate weights for the independent simulated annealing runs or adaptively tune multiple weights in a single run.More recently, MOSA-like algorithms based on dominance have also attracted some attention from the research community (Smith et al., 2008;Sanghamitra et al., 2008).
In recent years, EMO algorithms and MOSA-like algorithms have been investigated and developed along different research lines.However, less attention has been given to the combination of EMO algorithms and simulated annealing.Recently, we investigated the idea of using simulated annealing within MOEA/D for combinatorial optimization (Li and Landa-Silva, 2008).Our preliminary results were very promising.Following that idea, we now make the following contributions in this paper: • We propose a hybrid between MOEA/D and simulated annealing, called EMOSA, which employs simulated annealing for the optimization of each subproblem and adapts search directions for diversifying non-dominated solutions.Moreover, new strategies for competition among individuals and for updating the external population are also incorporated in the proposed EMOSA algorithm.
• We compare EMOSA to three MOSA-like algorithms and three MOMA-like algorithms on both the multi-objective knapsack problem and the multi-objective traveling salesman problem.The test instances used in our experiments involve two and three objectives for both problems.
• We also investigate the effects of important algorithmic components in EMOSA, such as the strategies for the competition among individuals, the adaptation of weights, the choice of neighborhood structure, and the use of -dominance for updating the external population.
The remainder of this paper is organized as follows.Section 2 gives some basic definitions and also outlines traditional methods in multi-objective combinatorial optimization.Section 3 reviews related work, while the description of the proposed EMOSA algorithm is given in Section 4. Section 5 describes the two benchmark MOCO problems used here.Section 6 provides the experimental results of the algorithm comparison, while Section 7 experimentally analyzes the components of EMOSA and parameter sensitivities.The final Section 8 concludes the paper.

Multi-objective Optimization
This section gives some basic definitions in multi-objective optimization and outlines two traditional multi-objective methods (weighted sum approach and weighted Tchebycheff approach) from mathematical programming.

Pareto Optimality
A multi-objective optimization problem (MOP) with m objectives for minimization1 can be formulated as: where x is the vector of decision variables in the feasible set X, F : X → Y ⊂ R m is a vector of m objective functions, and Y is the objective space.The MOP in ( 1) is called a multi-objective combinatorial optimization (MOCO) problem if X has a finite number of discrete solutions.
For any two objective vectors u = (u 1 , . . ., u m ) and v = (v 1 , . . ., v m ), u is said to dominate v, denoted by u ≺ v, if and only if u i ≤ v i for all i ∈ {1, . . ., m} and there exists at least one index i satisfying u i < v i .For any two solutions x and y, x is said to dominate y if F (x) dominates F (y).A solution x * is said to be Pareto-optimal if no solution in X dominates x * .
The set of all Pareto-optimal solutions in X is called Pareto-optimal set.Correspondingly, the set of objective vectors of solutions in the Pareto-optimal set is called Pareto-optimal front (POF).The lower and upper bounds of the POF are called the ideal point z * and the nadir point z nad respectively, that is, z * i = min u∈POF u i and A solution x is said to -dominate y if F (x) − dominates F (y), where = ( 1 , . . ., m ) with i > 0, i = 1, . . ., m.If x dominates y, x also -dominates y.However, the reverse is not necessarily true, see (Deb et al., 2005) for more on -dominance.

Traditional Multi-objective Methods
In mathematical programming, a multi-objective optimization problem is often converted into one or multiple single objective optimization subproblems by using linear or nonlinear aggregation of objectives with some weights.This is called decomposition or scalarization.In this section, we describe two commonly-used traditional methodsweighted sum approach and weighted Tchebycheff approach (Miettinen, 1999).

• Weighted Sum Approach
This method considers the following single objective optimization problem: where x ∈ X, λ = (λ 1 , . . ., λ m ) is the weight vector with λ i ∈ [0, 1], i = 1, . . ., m and m i=1 λ i = 1.Each component of λ can be regarded as the preference of the associated objective.A solution x * of the MOP in (1) is a supported optimal solution if it is the unique global minimum of the function in (2).

• Weighted Tchebycheff Approach
The aggregation function of this method has the following form: where x ∈ X, λ is the same as above, and z * is the ideal point.
Under some mild conditions, the global minimum of the function in (2) or ( 3) is also a Pareto-optimal solution of the MOP in (1).To find a diverse set of Pareto-optimal solutions, a number of weight vectors with evenly spread distribution across the tradeoff surface should be provided.The weighted Tchebycheff approach has the ability to deal with non-convex POF but the weighted sum approach lacks this ability.Normalization of objectives is often needed when the objectives are incommensurable, i.e. have different scales.1990s use weighted aggregation for their acceptance criteria in local search.The main differences between the various MOSA-like algorithms lie in the choice of weights during the search.The MOSA proposed by Serafini (1992) is a single-point method that optimizes one weighted aggregation function in each iteration (see Figure 1(a)).In this method, the current weight vector is smoothly changed as the search progresses.

Related Previous Work
Similarly, IMMOGLS and MOGLS optimize a single weighted aggregation function in each iteration (see Figure 1(d)).In contrast, two MOSA-like algorithms (Ulungu et al., 1999;Czyzak and Jaszkiewicz, 1998) are population-based search methods and optimize multiple weighted aggregation functions (see Figure 1(b) and (c)).In Ulungu's MOSA, a set of evenly-spread fixed weights is needed, while Czyzak's MOSA adaptively tunes the weight vector of each individual according to its location with respect to the nearest non-dominated neighbor.
In Zhang and Li (2007), a multi-objective evolutionary algorithm based on decomposition, called MOEA/D, was proposed.It optimizes multiple subproblems simultaneously by evolving a population of individuals.In this algorithm, each individual corresponds to one subproblem for which a weighted aggregation function is defined.Like Ulungu's MOSA, MOEA/D uses fixed weights during the whole search.Since subproblems with similar weight vectors have optimal solutions close to each other, MOEA/D optimizes each subproblem by using the information from the optimization of neighboring subproblems.For this purpose, recombination and competition Evolutionary Computation Volume x, Number x between neighboring individuals are taken into consideration for mating restriction.Unlike other EMO algorithms, MOEA/D provides a general framework which allows the application of any single objective optimization technique to its subproblems.
In Li and Landa-Silva (2008), we proposed a preliminary approach hybridizing MOEA/D with simulated annealing.In the present work, competition and adaptation of search directions are incorporated to develop an effective hybrid evolutionary approach, called EMOSA.In the proposed algorithm, the current solution of each subproblem is improved by simulated annealing with different temperature levels.After certain low temperature levels, the weight vectors of subproblems are modified in the same way as Czyzak's MOSA.Contrary to the original MOEA/D, no crossover is performed in our hybrid approach.Instead, diversity is promoted by allowing uphill moves following the simulated annealing rationale.
The hybridization of population-based global search coupled with local search is the main idea of many memetic algorithms (also called genetic local search), which have shown considerable success in single objective combinatorial optimization (Merz, 2000;Krasnogor, 2002;Hart et al., 2004).In these algorithms, genetic search is used to explore promising areas of the search space (diversification) while local search is applied to examine high-quality solutions in a specific local area (intensification).There have also been some extensions of memetic algorithms for multi-objective combinatorial optimization (Knowles and Corne, 2004).
In some well-known multi-objective memetic algorithms based on decomposition, such as MOGLS and MOEA/D, two basic strategies (random weights and fixed weights) have been commonly used to maintain diversity of search directions towards the POF.However, both strategies have their disadvantages.On the one hand, random weights might provoke the population to get stuck in local POF or approximate the same parts of the POF repeatedly.On the other hand, fixed weights might not be able to cover the whole POF well.In Czyzak and Jaszkiewicz (1998), the dispersion of non-dominated solutions over the POF is controlled by tuning the weights adaptively.
In our previous work (Li and Landa-Silva, 2008), the weight vector of each individual is modified on the basis of its location with respect to the closest non-dominated neighbor when the temperature goes below certain level.The value of each component in the weight vector is increased or decreased by a small value.However, the diversity of all weight vectors was not considered.The modified weight vector could be very close to others if the change is too big.In this case, the aggregation function with the modified weights might not be able to guide the search towards unexplored parts of the POF.To overcome this weakness, the adaptation of each weight vector should consider the diversity of all weight vectors in the current population.
As shown in Li and Zhang (2009), competition between solutions of neighboring subproblems in MOEA/D could affect the diversity of the population.Particularly, when the population is close to the POF, the competition between neighboring solutions provokes duplicated individuals in the population.Then, in some cases, competition should not be encouraged.For this reason, we need an adaptive mechanism to control competition at different stages during the search.
Like in many multi-objective metaheuristics, MOEA/D uses an external population to store non-dominated solutions.When the POF is very large, maintaining such population is computationally expensive.No diversity strategy is adopted to control Evolutionary Computation Volume x, Number x the size of the external population in MOEA/D.Popular diversity strategies are crowding distance in NSGA-II and the nearest neighbor method in SPEA-II.Both strategies have computational complexity of at least O(L 2 ), where L is the current size of the external population.The -dominance concept is a commonly-used technique for fitness assignment and diversity maintenance (Grosan, 2004;Deb et al., 2005).We use -dominance within EMOSA to maintain the diversity of the external population of non-dominated solutions.

Description of EMOSA
In this section, we propose EMOSA, a hybrid adaptive evolutionary algorithm that combines MOEA/D and simulated annealing.Like in MOEA/D, multiple single objective subproblems are optimized in parallel.For each subproblem, a simulated annealing based local search is applied to improve the current solution.To diversify the search towards different parts of the POF, the weight vectors of subproblems are adaptively modified during the search.The main features of EMOSA include: 1) maintenance of adaptive weight vectors; 2) competition between individuals with respect to both weighted aggregation and Pareto dominance; 3) use of -dominance for updating the external population.In this section, we first present the main algorithmic framework of EMOSA and then discuss its algorithmic component in more detail.

6:
for all s ∈ {1, . . ., Q} do 7: calculate the neighborhood Λ(s) for the s-th weight vector λ (s)   8: end for 9: repeat 10: estimate the nadir point z nad and the ideal point z * using EP 11: for all s ∈ {1, . . ., Q} do 12: apply SimulatedAnnealing(x (s) , λ (s) , T ) to x (s) and obtain y 13: compete y with the current population POP by UpdatePopulation(y, s) the current set of Q weight vectors Λ = {λ (1) , . . ., λ (Q) } and the candidate set Θ of weight vectors; 2) the current population POP = {x (1) , . . ., x (Q) } and the external population EP.The candidate weight vectors in Θ are generated by using uniform lattice points (see details in Section 4.2).The size of Θ is much larger than Q.The total number of weight vectors in the current weight set Λ equals the population size Q since each individual is associated to one weight vector.
The framework of EMOSA is shown in Procedure 1.Note that Θ, Λ, POP, and EP are global memory structures, which can be accessed from any subprocedure of EMOSA.In the following, we explain the main steps in EMOSA.
1.The initial weight vectors and population are generated in lines 1 and 2. The external population EP is formed with the non-dominated solutions in POP (line 3).
2. For each weight vector λ (s) , the associated neighborhood Λ(s) is computed based on the distance between weight vectors in lines 6-8.Once the current temperature level T goes below the final temperature level T min , the neighborhood should be updated.
3. The nadir and ideal points are estimated using all non-dominated solutions in EP (line 10).More precisely, These two points are used in the setting of -dominance.
4. For each individual in the population, simulated annealing local search is applied for improvement (line 12).The resulting solution y is then used to update other individuals in the population (line 13).
5. If the current temperature T is below the final temperature T min , the weight vector of each individual is adaptively modified (lines 17-19).
6.In line 4, the current temperature level T is set to T max (starting temperature level), and the current cooling rate is set to be α 1 (starting cooling rate).The current temperature is decreased in a geometric manner (line 15).To help the search escape from the local optimum, the current temperature is re-heated (line 20) to T reheat (< T max ) and a faster annealing scheme is also applied with α = α 2 (< α 1 ).
More details of the main steps in EMOSA are given in the following sections.

Generation of Diverse Weight Vectors
Here, Θ is the set of all normalized weight vectors with components chosen from the set: {0, 1/H, . . ., (H − 1)/H, 1}, with H a positive integer number.The size of Θ is C m−1 H+m−1 , with m the number of objectives.Figure 2 shows a set of 990 uniformlydistributed weight vectors produced using this method for H = 43 and m = 3.
In the initialization stage of EMOSA, Q initial weight vectors evenly spread are selected from the candidate weight set Θ. Note that the size of Θ is much larger than Q.This procedure SelectWeightVectors is shown in Procedure 2. In lines 1-3, the extreme weight vectors are generated.Each of these extreme weight vectors corresponds to the optimization of a single objective.In line 4, c is the number of weight vectors selected so far and all selected weight vectors are copied into a temporary set Φ.In line 6, the set A of all weight vectors in Θ with the same maximal distance to Φ are identified.In , the set B is formed with half of the recently selected weight vectors.In line 8, the weight vector in A with the maximal distance to B is chosen as the next weight vector in Λ.In this way, all weight vectors selected so far are well-distributed.For example, Figure 3 shows the distribution of 50 weight vectors selected from the initial 990 weight vectors shown in Figure 2.
The procedure SelectWeightVectors involves O(Q 2 × |Θ|) basic operations in lines 6 and 7 as well as O(|Θ| 2 ) distance calculations between the weight vectors in Θ. Ideally, the large size of Θ leads to a good approximation to the POF since each Pareto optimal solution is the optimal solution of a certain scalarizing function.However, as the size of Θ increases, the computational complexity of SelectWeightVectors and also AdaptWeightVectors (see details in the next subsection) increase.But, a small size of Θ might not be enough to approximate the whole Pareto front.To guarantee efficiency, H should be set properly.
In EMOSA, the neighborhood of each weight vector is determined using the same method as in MOEA/D.More precisely, the neighborhood Λ (s) of λ (s) for s = 1, . . ., Q, is formed by the indexes of its K closest neighbors in Λ.Note that the neighborhoods of weight vectors in MOEA/D remain unchanged during the whole search process while those in EMOSA are adaptively changed.This procedure involves O(Q 2 ) distance cal- In Czyzak's MOSA, the basic idea is to move each non-dominated solution away from its nearest neighbor.To implement this, the weight vector λ of each individual x is modified at each local move step according to the location of the nearest nondominated neighbor of x, lets say y.The components of λ are changed as follows: where δ > 0 is the variation level.Note that λ should be normalized after the above change.However, there are two weaknesses in this method.First, the setting of δ should be provided empirically.On the one hand, a large value of δ could push the current solution towards the neighborhood of other solutions.On the other hand, a very small value of δ might not be able to guide the search towards unexplored areas of the POF.Second, the diversity of all weight vectors in the current population is not considered in setting the value of λ.This will reduce the chance to approximate different parts of the POF with the same possibility.
Procedure 3 AdaptWeightVectors(s) To overcome the above weaknesses, we introduce a new strategy in EMOSA, illustrated in Figure 4, for adapting weight vectors when T < T min , i.e. when the simulated annealing enters the only improving phase.The corresponding procedure AdaptWeightVectors is shown in Procedure 3. Instead of changing the components of the weight vectors, our strategy picks one weight vector λ from the candidate set Θ to replace the current one λ (s) .For the current solution x (s) , we need to find its closest non-dominated neighbor x (t) , which corresponds to λ (t) .The selected weight vector λ must satisfy two conditions: 1) dist(λ (s) , λ (t) ) < dist(λ, λ (t) ), and 2) dist(λ, λ (s) ) ≤ dist(λ, Λ).Condition 1) indicates that the selected weight vector should increase the distance between the weight vectors of x (s) and x (t) .This could cause an increase in the distance between x (s) and x (t) in the objective space.Condition 2) guarantees that the selected weight vector should not be too close to other weight vectors in the current population.In line 2 of Procedure 3, all weight vectors in Θ satisfying these two conditions are stored in A. From this set, the weight vector with the maximal distance to λ (s) is selected.The procedure AdaptWeightVectors involves O(Q × |Θ|) distance calculations when the current temperature is below the value of T min .

Local Search and Evolutionary Search
In EMOSA, each individual in the population is improved by a local search procedure, namely SimulatedAnnealing shown in Procedure 4.Then, the improved solutions are used to update other individuals in the neighborhood as shown in Procedure 5 UpdatePopulation.
Procedure 4 SimulatedAnnealing(x (s) , λ (s) , T ) update EP in terms of -dominance 6: calculate the acceptance probability P (y, x , λ (s) , T ) end if 11: until Stopping conditions are satisfied 12: return y In line 3 of Procedure 4, a neighboring solution x is generated from N (y), the neighborhood of the current solution y.If F (x ) is not dominated by F (y), then update EP in terms of -dominance.The components of are calculated by Here, β > 0 is a parameter to control the density of solutions in EP.
The update of EP using -dominance benefits two aspects: 1) maintain the diversity of EP and 2) truncate the solutions of EP in dense areas of the objective space.In this work, β is set to 0.002 for two-objective problems and 0.005 for three-objective problems.
In the simulated annealing component of EMOSA, the probability of accepting neighboring solutions is given by: where T ∈ (0, 1] is the normalized temperature level, τ is a problem-specific balance constant, and g is a weighted aggregation function of x with the weight vector λ, which can be the weighted sum function g (ws) or the weighted Tchebycheff approach g (te) or the combination of both.According to (5), more uphill (worse) moves are accepted with high probability when T is close to 1, but increasingly only downhill (better) moves are accepted as T goes to zero.The acceptance of worse neighboring solutions reduces the possibility of getting trapped in local optima.
The degree of uphill moves acceptance is determined by the annealing schedule, which is crucial to the performance of simulated annealing.In EMOSA, we apply the following two annealing schedules, shown in Figure 5, in two stages.
• Schedule 1: Before the current temperature T goes below the final temperature T min for the first time (early stage), local search seeks to improve all individuals using high temperature T max and slow cooling rate α 1 (line 4 in Procedure 1).This is the same as many other simulated annealing algorithms.The main task at this stage is to start from the initial population and then quickly find some representative solutions in the POF from the initial population.
• Schedule 2: Whenever the current temperature T goes below the final temperature T min (late stage), T is increased to a medium temperature level T reheat , and a faster cooling rate α 2 (< α 1 ) is applied (line 20 in Procedure 1).Since all solutions in the current population should by now be close to the POF, they should also be close to the optimal solutions of aggregation functions with modified weight vectors.For this reason, there is no need to start a local search from a high temperature level.Fast annealing could help the intensification of the search near the POF since not many uphill moves are allowed.
Moreover, the parameter τ in (5) should be set empirically.As suggested in Smith et al. (2004), about half of non-improving neighbors should be accepted at the beginning of simulated annealing.Based on this idea, we set τ to − log(0.5)Tmax / ∆ in EMOSA, where ∆ is the average of ∆g values in the first 1000 uphill moves.For these first 1000 uphill moves, the acceptance probability is set to 0.5.In EMOSA, the Sim-ulatedAnnealing procedure is terminated after examining #ls neighbors.Hence, at each temperature level, this procedure mainly involves Q × #ls function evaluations and the update of the external population EP.
Note that EMOSA optimizes different subproblems using the same temperature cooling scheme.Actually, each subproblem can also be solved under different cooling schemes.This idea has been used in previous work on simulated annealing algorithm for multi-objective optimization (Tekinalp and Karsli, 2007) and also for single objective optimization (Burke et al., 2001).
Procedure 5 UpdatePopulation(y, s) The key difference between EMOSA and other population-based MOSA-like algorithms, such as Ulungu and Czyzak's algorithms, is the competition between individuals.As shown in Li and Zhang (2009), the competition may affect the balance between diversity and convergence in MOEA/D.On the one hand, the competition among individuals could speed up the convergence towards the POF when the current population is distant to the POF.On the other hand, the competition between individuals close to the POF could cause loss of diversity.This is because one solution that is non-dominated with respect to its neighbors could be worse if the comparison is made using a weighted aggregation function.In this paper, we apply different criteria to update the current solution and its neighbors.The details of the process for updating the population are shown in Procedure 5.In lines 1-3, the current solution x (s) is compared to the solution y obtained by local search and the comparison is made using the weighted aggregation function.In lines 4-8, the neighbors of x (s) are replaced only if they are dominated by y.By doing this, the selection pressure among individuals is weak when the current population is close to the POF.Therefore, the diversity of individuals in the neighborhood can be preserved.

Two Benchmark MOCO Problems
In this paper, we consider two N P-hard MOCO test problems, the multi-objective knapsack problem and the multi-objective traveling salesman problem.These problems have been widely used when testing the performance of various multi-objective metaheuristics (Jaszkiewicz, 2002;Zhang and Li, 2007;Paquete and St ützle, 2003;Garcia-Martinez et al., 2007;Li et al., 2010).

The 0-1 Multi-objective Knapsack Problem
Given n items and m knapsacks, the 0-1 multi-objective knapsack problem (MOKP) can be formulated as: subject to where x = (x 1 , . . ., x n ) is a binary vector.That is, x j = 1 if item j is selected to be included in the m knapsacks and x j = 0 otherwise, p ij and w ij are the profit and weight of item j in knapsack i, respectively, and c i is the capacity of knapsack i.
Due to the constraint on the capacity of each knapsack, the solutions generated by local search moves or genetic operators could be infeasible.One way to deal with this is to design a heuristic for repairing infeasible solutions.Some repair heuristics have been proposed in the literature (Zitzler and Thiele, 1999;Zhang and Li, 2007).In this paper, we use the same greedy heuristic adopted in Zhang and Li (2007) for constraint handling.The basic idea is to remove some items with heavy weights and little contribution to an objective function from the overfilled knapsacks.
In Knowles and Corne (2000b), neighboring solutions for the MOKP are produced by flipping each bit in the binary vector with probability k/n, called k-bit flipping here.For example, if the j-th bit of x is chosen for flipping, then x j = 1 − x j .After flipping, the modified solution x might be infeasible and a repair heuristic might be needed.In this paper, we suggest a new neighborhood structure, k-bit insertion, in which only the bits equal to zero may be flipped.This means that more items are added into the knapsacks but no items are removed.Consequently, the modified solution after flipping bits with value of zero is very likely to be infeasible but will have higher profits.Then, the repair heuristic is called to remove some items as described above until the solution becomes feasible again.

The Multi-objective Traveling Salesman Problem
The Traveling Salesman Problem (TSP) can be modeled as a graph G(V, E), where V = {1, . . ., n} is the set of vertices (cities) and E = {e i,j } n×n is the set of edges (connections between cities).The task is to find a Hamiltonian cycle of minimal length, which visits each vertex exactly once.In the case of the multi-objective TSP (MOTSP), each edge e i,j is associated with multiple values such as cost, length, traveling time, etc.Each of them corresponds to one criterion.Mathematically, the MOTSP can be formulated as: where π = (π(1), . . ., π(n)) is a permutation of cities in Π, the set of all permutations of I = {1, . . ., n} and c (i) s,t , s, t ∈ I is the cost of the edge between city s and city t regarding criterion i.Unlike the MOKP problem, no repair strategy is needed for solutions to the MOTSP problem because all permutations represent feasible solutions.
Like in the single objective TSP problem, the neighborhood move used here for the MOTSP is the 2-opt local search procedure constructed in two steps: 1) remove two non-adjacent edges from the tour and then 2) reconnect vertices associated with the removed edges to produce a different neighbor tour.

Performance Indicators
The performance assessment of various algorithms is one of the most important issues in the EMO research.It is well established that the quality (in the objective space) of approximation sets should be evaluated on the basis of two criteria: the closeness to the POF (convergence) and the wide and even coverage of the POF (diversity).In this paper, we use the following two popular indicators for comparing the performance of the multi-objective metaheuristics under consideration.
Ideally, the reference set should be the whole true POF.Unfortunately, the true POF of many MOCO problems is unknown in advance.Instead, a very good approximation to the true POF can be used.In this paper, the reference sets P * are formed by collecting all non-dominated solutions obtained by all algorithms in all runs.The smaller the IGD-metric value, the better the quality of the approximation set P .To have a lower value of the IGD-metric, the obtained solutions should be close to P * and should not miss any part of P * .Therefore, the IGD-metric can assess the quality of approximation sets with respect to both convergence and diversity.

Hypervolume (S-metric or Lebesgue measure)
This indicator (see Zitzler and Thiele (1999)) measures the volume enclosed by the approximation set P wrt a reference point r = (r 1 , . . ., r m ).In the case of minimization, this metric can be written as S(P, r) = VOL u∈P E(u, r) , where E(u, r) = {v|u i ≤ v i ≤ r i , i = 1, . . ., m} (i.e.E(u, r) is the hyperplane between two mdimensional points u and r)2 .Compared to the reference set required in the IGDmetric, establishing the reference point in the S-metric is relatively easier and more flexible.In our experiments, r can be set to be a point dominated by the nadir point: ), i = 1, . . ., m.The larger the S-metric value, the better the quality of the approximation set P , which means a smaller distance to the true POF and a better distribution of the approximation set.Similar to the IGD-metric, better S-metric values can be achieved when the obtained solutions are close to the POF and achieve a full coverage of the POF.Therefore, the S-metric also estimates the quality of approximation sets both in terms of convergence and diversity.

Experimental Settings
In this first set of experiments, we used the test instances with two and three objectives shown in Table 1.We used the same following parameters in the annealing schedule of the four algorithms: T max = 1.0, T min = 0.01, T reheat = 0.1, α 1 = 0.8, and α 2 = 0.5.Moreover, all four algorithms were allocated the same number of function evaluations  At every temperature level, the number #ls of neighbors examined by local search is set to 10 for each MOKP instance and 250 for each MOTSP instance.It should be noted that the total number of neighbors examined in SMOSA at each temperature level equals the total number of neighbors examined for the whole population in UMOSA and CMOSA at the same temperature level.For all test instances in both problems, the neighborhood size in EMOSA is set to K = 10.In both SMOSA and CMOSA, each component of the weight vectors is increased or decreased by 2.0/Q randomly or using the heuristic in eq.(4).

Experimental Results
The mean and standard deviation values of the IGD-metric and S-metric found by the four algorithms are shown in Tables 3 and 4, respectively.From these results, it is evident that EMOSA outperforms the three other MOSA-like algorithms in terms of both performance indicators on all six MOKP test instances.The superiority of EMOSA is due to the stronger selection pressure to move towards the POF by having competition among individuals.We believe that this competition allows the local search in EMOSA to start from a high-quality initial solution at every temperature level.
We can also observe from Tables 3 and 4 that CMOSA is clearly inferior to the three other MOSA-like algorithms.This is due to the fact that CMOSA modifies the weight vector of each individual while examining every neighboring solution.This change of weight vectors in the early stage of the multi-objective search does not affect the distribution of the non-dominated solutions because the current population is still distant to the POF.However, the frequent change of weight vectors could potentially guide each individual to move along different search directions in the objective space.Of course, this will slow down the convergence speed towards the POF.This is why our proposed EMOSA approach only adapts the weight vector of each individual once the temperature is below T min .Figure 6 shows the non-dominated solutions found by the four algorithms on three bi-objective MOKP instances for the run for which each algorithm achieved its best value for the IGD-metric.We can see that EMOSA clearly found better solutions than the other three MOSA-like algorithms on KS5002 and KS7502.For the smallest instance KS2502, EMOSA is only slightly better.Moreover, both UMOSA and SMOSA perform quite similarly on these three instances while CMOSA has the tendency to locate the non-dominated solutions close to the middle part of the POF.
In Figure 7, we plot the non-dominated solutions obtained by the four algorithms on the largest MOKP instance with three objectives (KS7503) in the run for which each algorithm achieved its best value for the IGD-metric.It is clear that EMOSA produced the non-dominated front with the best diversity while the non-dominated front by CMOSA has the worst distribution.This is because EMOSA tunes the weight vectors according to distances of both objective vectors and weight vectors to maintain the diversity of the set of current weight vectors, which are less likely to be very similar.In contrast, CMOSA is unable to preserve the diversity of weight vectors particularly in the case of test instances with more than two objectives.
From Figure 7, we can also note that the set of non-dominated solutions found by UMOSA consists of many isolated clusters.This is because this MOSA-like algorithm with fixed weights has no ability to approximate all parts of the POF when the size of this front is very large.However, it can still produce a good approximation to the POF.This is why EMOSA uses fixed weights only in the early stage of its search process.It should be pointed out that the population size used in UMOSA and CMOSA is six times that of EMOSA for all 3-objective instances (50 vs. 300).This also demonstrates the efficiency of the strategy for adapting weights in EMOSA which requires a considerably smaller population size.
Tables 5 and 6 show the mean and standard deviation values of the IGD-metric and S-metric obtained by the four MOSA-like algorithms on the six MOTSP instances.Again, EMOSA performs better than the other three MOSA-like algorithms on all six MOTSP test instances in terms of the IGD-metric and the S-metric.In contrast, CMOSA shows again the worst performance.
Figure 8 shows the distributions of the non-dominated solutions obtained by the four MOSA-like algorithms in the run for which each algorithm achieved its best value for the IGD-metric value.It is clear that EMOSA found better non-dominated solu-  tions in terms of both convergence and diversity.Figure 9 shows the non-dominated solutions obtained in the best run of each algorithm for the three-objective instance -KROABC100.It can be seen that the diversity of the solutions found by EMOSA is the best among the four algorithms.
From all above results, we can also conclude that the overall performance of UMOSA is better than that of CMOSA and SMOSA for both problems.These results are also consistent with those reported in previous studies (Banos et al., 2007).The poor performance of CMOSA and SMOSA is caused by the frequent modification of weight vectors, which might slow down the convergence speed.
The four MOSA-like algorithms use different strategies to maintain the diversity of weight vectors.The running times of these algorithms rely on the computational complexity of the related strategies.Table 7 shows the average computational time used by the four MOSA-like algorithms on two 100-city MOTSP instances.From this table, we can make the following observations: • For the smaller 2-objective instance KROAB100, UMOSA is the fastest algorithm while CMOSA is the slowest.This is because UMOSA maintains a set of fixed weight vectors which does not have a computational cost because it is not required to tune the weight vectors.In contrast, both CMOSA and SMOSA modify weight vectors at each single local search move.The computational complexity of the strategy in CMOSA is much higher than that in SMOSA, however EMOSA is still faster.Therefore, it seems that it is not reasonable to tune the weight vectors at every local search move.
• For the larger 3-objective instance KROABC100, both EMOSA and UMOSA are slower than SMOSA.This is contrary to the above observation.The reason for this is that the the number of non-dominated solutions found by the former two algorithms is much larger than those found by SMOSA.In this case, most of the computational cost is allocated to maintaining the external population.This is also the reason why CMOSA is faster than EMOSA.

Experimental Settings
In this section, we compare EMOSA to three MOMA-like algorithms, namely MOEA/D, MOGLS, and IMMOGLS.All parameter settings used in EMOSA are the same as in the previous section.The population sizes used in the other three algorithms are given in Table 8.Like EMOSA, the neighborhood size used in MOEA/D is also 10 for all test instances.In MOGLS, the size of the mating pool is set to 20 for every test instance.In IMMOGLS, 10 elite solutions are randomly selected from the external  population.For the MOKP problem, single-point crossover is performed, and each bit of the binary string is mutated by flipping with probability 0.01.The neighborhood structures for both the MOKP problem and the MOTSP problem are the same as in the previous section.For the MOTSP instances, offspring solutions are sampled by using cycle crossover (Larranaga et al., 1999) in these three MOMA-like algorithms.

Experimental Results
The mean and standard deviation of the IGD-metric and S-metric values found by EMOSA and the three MOMA-like algorithms on six MOKP test instances are shown in Table 9 and Table 10, respectively.From this set of results, it is evident that EMOSA performs better than the three MOMA-like algorithms on all test instances.We think that the superiority of EMOSA is due to i) the use of a diverse set of weight vectors and ii) the use of simulated annealing for local search.
Among the four algorithms, IMMOGLS has the worst performance on all test instances.This is because IMMOGLS uses a fixed population size and updates the current population with all offspring solutions without competition.In EMOSA, individuals with similar weight vectors interact with each other via competition, which helps for a good performance of EMOSA in terms of convergence.Zhang and Li (2007) showed that MOEA/D constantly performs better than MOGLS on the MOKP test instances when local search heuristic is applied.When simple hill-climber local search is used in both algorithms, MOEA/D performs slightly worse than MOGLS regarding the IGD-metric since the population size is not large    enough, but better in terms of the S-metric.In contrast, EMOSA clearly outperforms MOGLS, which gives an indication that the use of advanced local search based metaheuristics can enhance the performance of EMO algorithms.
The final solutions found by EMOSA and the three MOMA-like algorithms on three 2-objective MOKP instances are plotted in Figure 10.For KS2502, EMOSA performs similarly to MOEA/D and MOGLS.However, the solutions found by EMOSA are better than those obtained by MOEA/D and MOGLS on two larger test instances KS5002 and KS7502.It is clear that the solutions of IMMOGLS are dominated by those of the other three algorithms.
Figure 11 shows the solutions found by the four algorithms on instance KS7503 for the run in which each algorithm achieved the best IGD-metric value.It is evident that EMOSA has the best performance wrt diversity.In MOEA/D, the non-dominated solutions are located in a number of clusters.This is because MOEA/D uses fixed weights during the whole search process.Therefore, the performance of MOEA/D depends on the population size.For large test instances, a reasonable large number of weight vectors is required to achieve good diversity in the final non-dominated set of solutions.However, although the population size in EMOSA is small, this algorithm is still capable of solving large test instances due to the effective adaptation of weight vectors.
Table 11 and Table 12 give the mean and standard deviation of IGD-metric and S-metric values found by EMOSA and the three MOMA-like algorithms on six MOTSP test instances.From these experimental results, it is clear that EMOSA found the best solutions in terms of both indicators.Figure 12 shows the distribution of the solutions obtained by the four algorithms on four bi-objective MOTSP test instances.As we can see, both EMOSA and MOEA/D have good performance in terms of convergence.The solutions found by both algorithms are quite close to each other.As for diversity, the distribution of the solutions found by EMOSA is more uniform.Both MOGLS and IMMOGLS have worse performance than EMOSA in all test instances.This could be due to the lack of efficient genetic operators for sampling promising offspring solutions for multi-objective local search.Figure 13 shows the distribution of the non-dominated solutions found by the four algorithms on the 3-objective KROABC100 instance.Unlike the concave POF of instance KS7503, the POF of instance KROABC100 is convex.We can observe that EMOSA has the best performance in finding a good distribution of non-dominated solutions.
7 Effect of Algorithmic Components and Analysis of Parameter Sensitivity in EMOSA 7.1 Effect of Algorithmic Components

Adaptation of Weight Vectors
In our early work on EMOSA (Li and Landa-Silva, 2008), we tried the same strategy as in CMOSA to adapt weight vectors.As we discussed in section 3, that strategy has some weaknesses.In this section, we demonstrate the improvement of EMOSA due to the new diversity strategy.We modified the EMOSA approach shown in Procedure 1 by replacing our weights adaptation strategy with the strategy used in CMOSA and applied the modified algorithm to instance KROABC100.All parameters in this modified EMOSA remain the same as those in section 6.2.1.Figure 14 shows the non-dominated solutions found by this modified EMOSA.Compared to the results in Figure 13, the modified EMOSA performs worse in terms of diversity.This indicates that our new di- versity strategy works better in guiding the multi-objective search towards unexplored parts of the POF by maintaining the diversity of weight vectors.

Competition Strategies
As we said before, EMOSA distinguishes itself from other algorithms by the adaptive competition scheme between individuals.We carried out three further experiments to study the effect of the competition among individuals on the performance of EMOSA.
To study the effect of weighted aggregation function on the convergence of EMOSA, we compared two versions of EMOSA, one using the weighted sum approach and the other one using the weighted Tchebycheff approach.Both versions use the  same parameter settings described in Section 6.2.1.Figure 15 shows the non-dominated solutions found by these two versions.From these results, we can observe that the weighted Tchebycheff approach performs very similar to the weighted sum approach on KROAB100.Zhang and Li (2007) made a different observation, that the former has worse performance on the MOKP problem.The reason for this is that the weighted sum approach has higher effectiveness in preserving the quality of infeasible solutions in the repair heuristic.Also, we tried the version of EMOSA with no competition between individuals.That is, the current solution of each weighted aggregation function does not compete with its neighboring solutions.To implement this, lines 4-8 in Procedure 5 were removed.We tested this version of EMOSA on the KROAB100 instance using the same parameters as in section 6.2.1.Figure 16 shows the non-dominated solutions found by EMOSA without competition between individuals.Clearly, the convergence of this version is worse than the EMOSA with competition.
Moreover, we tested another version of EMOSA with competition based on only scalarization.In Procedure 5, lines 5-7 were removed.Instead, each neighboring solution x (k) is replaced with y if g(y, λ (k) ) < g(x (k) , λ (k) ).This variant of EMOSA was tested on instance KROABC100.The final solutions found by this variant are plotted in Figure 17.It can be observed that the modified EMOSA failed to find a set of welldistributed non-dominated solutions despite having the adaptation of weight vectors.As mentioned before, this is caused by the loss of diversity in the population.When the current population is close to the POF, neighboring solutions are very likely to be mutually non-dominated.However, one solution can still be replaced by its neighbors in terms of a weighted sum function.In this case, EMOSA might miss the chance to approximate the part of POF nearby such a solution.This is why EMOSA does not work well on KROABC100 wrt diversity.

Choice of Neighborhood Structure
The choice of appropriate neighborhood structure plays a very important role in local search methods for combinatorial optimization.In this paper, we have proposed a k-bit insertion neighborhood structure for the MOKP problem.In this section, we compare two versions of EMOSA on the instance KS2502, one using the k-bit insertion neighborhood and another version using the k-bit flipping neighborhood.Figure 18 shows the final solutions found by both versions in the best of 20 runs regarding the IGD-metric.It can be observed that the non-dominated solutions obtained by EMOSA with the k-bit insertion neighborhood dominate all solutions produced by the version with the k-bit flipping neighborhood.Moreover, the mean IGD-metric value obtained when using the k-bit flipping neighborhood is 76.934, which is worse than the value obtained when using the k-bit insertion (16.245).

Updating the External Population with -dominance
In MOEA/D, no diversity strategy is applied to manage the external population.This is not a problem when the size of the POF is small.However, for large fronts, the external population might contain more and more non-dominated solutions as the search progresses.The computational time needed to update the external population may exceed substantially the time used by all other components of EMOSA.This is why we incorporate -dominance into EMOSA to help in keeping the diversity of EP automatically. Figure 19 shows the results obtained by EMOSA without using -dominance on instance KROABC100.Compared to the results in Figure 9, the distribution of nondominated solutions is more dense that the distribution obtained with EMOSA using -dominance.Up to 20,000 non-dominated solutions are saved in the external population while EMOSA with -dominance reports only about 3000 solutions at the end of each run.In terms of the computational time, we also noticed that EMOSA without -dominance used about 1600 seconds to complete one run, while EMOSA with -dominance only spent less than 170 seconds to complete a single run, which represents almost a ten-fold reduction.

Parameters in Two-stage Cooling Scheme
Like the simulated annealing algorithms for single objective optimization, the performance of EMOSA might also be influenced by the setting of parameters in the cooling scheme.To verify this, we studied two parameters T reheat (reheat temperature value) and α 2 (fast cooling rate).EMOSA with the combinations of T reheat = {0.1,0.3, 0.5, 0.8, 1.0} and α 2 = {0.3,0.5, 0.8} were tested on KROAB50 and KROAB100 in 20 runs.ent combinations of T reheat and α 2 .It can be observed that EMOSA with lower reheat temperature has better performance.For the smaller instance KROAB50, the best result was found by EMOSA using a faster cooling rate (0.5).However, such a strategy does not perform best for KROAB100.For this larger instance, EMOSA with a slower cooling rate (0.8) worked better.For the small instance, the current population of EMOSA is very likely to be well-converged towards the Pareto front.In this case, EMOSA needs to intensify the search near the Pareto front.Therefore, the temperature levels in EMOSA should be in lower values.In contrast, EMOSA with slow cooling rate is able to increase the diversity of the search for large instances.

Population Size (Q) and the Number of Uniform Weight Vectors (H)
In EMOSA, the set of current weight vectors are chosen from a set of predefined uniform weight vectors.Such a set is controlled by a parameter H.In this section, we study the influence of H on the performance of EMOSA.We tested EMOSA with a smaller H=23 (300 uniform weight vectors) and a larger H=53 (1485 uniform weight vectors) for KROABC100.Figure 21 shows the distribution of non-dominated solutions found by EMOSA with H = 23 and H = 53.We can observe from this figure that EMOSA with the larger value of H performed better than that with the smaller value of H wrt diversity.It is easy to understand that more uniform weight vectors are needed in EMOSA when the size of Pareto front is very large.However, a very large value of H will increase the computational time of adapting weight vectors in EMOSA.
As we mentioned before, EMOSA produced better results than other MOMA-like algorithms although its population size is much smaller.The main reason for this is that EMOSA uses a strategy to adapt weight vectors.We also tested EMOSA with population sizes of 100, 200, 300, 400 and 500 on instance KROABC100.Correspondingly, the obtained mean IGD values were 4308.13, 4078.63, 4396.37, 5280.28 and 5792.63.It is clear that EMOSA with population sizes smaller than 300 found similar IGD-metric values.In contrast, EMOSA with larger population sizes performed worse.

Conclusions
Both evolutionary algorithms and simulated annealing are popular stochastic search techniques for tackling multi-objective combinatorial optimization.However, the hybridization of both methods has not yet been studied in the context of MOCO problems.Most of the existing research work on multi-objective memetic algorithms only focuses on the use of simple hill-climber local search within evolutionary algorithms.In this paper, we proposed an adaptive evolutionary multi-objective approach combined with simulated annealing.We call this approach Evolutionary Multi-objective Simulated Annealing (EMOSA).The algorithm seeks to optimize multiple weighted aggregation functions.To approximate various parts of the POF, we have suggested a new method to tune the weight vectors of these aggregation functions.Moreover, -dominance was used to update the external population of non-dominated solutions in EMOSA for enhancing computational efficiency.
We also compared EMOSA to six other algorithms using weighted aggregation (three MOMA-like algorithms and three MOSA-like algorithms) on both the multiobjective knapsack problems and the multi-objective traveling salesman problem.Our experimental results show that EMOSA performs better than all other algorithms considered in terms of both the IGD-metric (inverted generational distance) and the Smetric (hypervolume).In particular, EMOSA is capable of finding a very good distribution of non-dominated solutions for the 3-objective test instances of both problems.The main reason for this is that our new strategy for tuning weight vectors directs the search towards different portions of the POF in a more effective way than the old strategy.Moreover, we also studied here the effects of important components on the performance of the proposed EMOSA algorithm.We have shown that the use ofdominance significantly reduces the computational time for maintaining the external population of non-dominated solutions.
It should be noted that the POF shapes of all test instances in both MOCO problems considered in this paper do not present difficulties to EMOSA using the weighted sum approach.As future work, in order to deal with some hard POF shapes, such as discontinuity, we intend to design a robust version of EMOSA to tackle such MOCO problems.Moreover, we also intend to investigate the ability of EMOSA to solve other challenging MOPs, such as continuous MOPs with many local optima, MOCO problems with complex constraints, many-objective problems, and dynamic MOPs.The C++ source code of EMOSA is available from the authors at request.

Figure 1 :
Figure 1: Graphical illustration of different strategies used in multi-objective metaheuristics based on decomposition.
settings: T ← T reheat and α ← α 2 21: until stopping conditions are satisfied 22: return EP Two sets of weight vectors and two populations of individuals are maintained: 1) Evolutionary Computation Volume x, Number x
6.1.1Inverted Generational Distance (IGD-metric) This indicator (see Coello Coello and Cruz Cortés (2005)) measures the average distance from a set of reference points P * to the approximation set P .It can be formulated as Evolutionary Computation Volume x, Number x

Figure 7 :
Figure 7: Non-dominated solutions found by the four MOSA-like algorithms on the three-objective knapsack instance KS7503.

Figure 8 :Figure 9 :
Figure 8: The non-dominated solutions found by MOSA-like algorithms on bi-objective MOTSP instances.

Figure 11 :
Figure 11: Non-dominated solutions found by MOMA-like algorithms on the KS7503 instance.

Figure 12 :
Figure 12: Non-dominated solutions found by MOMA-like algorithms on bi-objective TSP instances.

Figure 13 :
Figure 13: Non-dominated solutions found by MOMA-like algorithms on the KROABC100 instance.

Figure 20 :Figure 21 :
Figure 20: The mean IGD-metric values found by EMOSA with 15 combinations of T reheat and α 2 .

Table 1 :
Test instances of MOKP and MOTSP and computational cost.

Table 2 :
The setting of population size in MOSA-like algorithms.

Table 3 :
Mean and standard deviation of IGD-metric values on MOKP instances.

Table 4 :
Mean and standard deviation of of S-metric values on MOKP instances.each of 20 runs for the same test instance as detailed in column # nfes of Table1.The population sizes for the three population-based MOSA-like algorithms are given in Table2.We use a smaller (than in UMOSA and CMOSA) population size in EMOSA for most of the test instances.The reason for this is that EMOSA has the ability to search different parts of the POF by adapting weight vectors.Note that SMOSA is a single point based search method. in

Table 5 :
Mean and standard deviation of IGD-metric values on MOTSP instances.

Table 6 :
Mean and standard deviation of S-metric values on MOTSP instances.

Table 7 :
The average running time (seconds) of four MOSA-like algorithms on two 100-city MOTSP instances.

Table 8 :
The population size used in the MOMA-like algorithms.

Table 9 :
Mean and standard deviation of IGD-metric values found by MOMA-like algorithms on MOKP instances.

Table 10 :
Mean and standard deviation of S-metric values found by MOMA-like algorithms on MOKP instances.

Table 11 :
Mean and standard deviation of IGD-metric values found by MOMA-like algorithms on MOTSP instances.

Table 12 :
Mean and standard deviation of S-metric values found by four MOMA algorithms on six MOTSP instances.