A New Threshold Based Penalty Function Embedded MOEA / D

Recently, we proposed a new threshold based penalty function. The threshold dynamically controls the penalty to infeasible solutions. This paper implants the two different forms of the proposed penalty function in the multiobjective evolutionary algorithm based on decomposition (MOEA/D) framework to solve constrained multiobjective optimization problems. This led to a new algorithm, denoted by CMOEA/D-DE-ATP. The performance of CMOEA/D-DE-ATP is tested on hard CFseries test instances in terms of the values of IGD-metric and SC-metric. The experimental results are compared with the three best performers of CEC 2009 MOEA competition. Experimental results show that the proposed penalty function is very promising, and it works well in the MOEA/D framework. Keywords—Constrained multiobjective optimization; decomposition; MOEA/D; penalty function; threshold. I. I NTRODUCTION In this paper, we consider the following constrained multiobjective optimization problem (CMOP) [1]: Minimize F (x) = (f1(x), f2(x), . . . , fm(x)) ; Subject to gj(x) ≥ 0, j = 1, . . . , p; lk ≤ xk ≤ uk, k = 1, . . . , n, (1) wherex = (x1, . . . , xn) ∈ R is an n dimensional vector of decision variables, F is the objective vector function that consists ofm real-valued objective functions, and gi(x) ≥ 0 are inequality constraints. The objective and constraint functions, fi’s andgj ’s, could be linear or non linear real-valued functions.lk anduk are the lower and upper bounds (called bound constraints) of xk, k = 1, . . . , n, respectively, which define the search region S = {x = (x1, . . . , xn) | lk ≤ xk ≤ uk, k = 1, . . . , n}. A solutionx ∈ S is called a feasible solution, if it satisfies all the inequality constraints in (1). The set of all feasible solutions is called the feasible region. Mathematically, we can write: F = {x ∈ S ⊂ R|gj(x) ≥ 0, j = 1, · · · , p}. However, If a solution is not feasible, we call it infeasible. The set of all infeasible solutions is called the infeasible region. The feasible attainable objective set (AOS) can be defined as {F (x)|x ∈ F}. Since the objectives in (1) more often contradict each other, so it is hard to find a single solution in F that could minimize all the objectives at the same time. Instead, one looks for a set of optimal compromising/tradeoff feasible solutions. The best tradeoffs among the objectives can be defined in terms of Pareto-optimality [2], [3]. A solutionx is said to Pareto-dominate or simply dominate another solutiony, mathematically denoted as x y, if fi(x) ≤ fi(y), ∀i = 1, . . . ,m and fj(x) < fj(y) for at least onej ∈ {1, . . . ,m}1. This definition of domination is sometimes referred to as a weak dominance relation. A solution x ∈ F is Pareto-optimal to (1) if there is no solution x ∈ F such thatF (x) F (x). F (x) is then called a Pareto-optimal (objective) vector. The set of all Paretooptimal solutions is called the Pareto Set (PS) in the decision space and Pareto Front (PF) in the objective space [2]. In the majority of constrained optimization problems, the optimal solutions lie on the constraints’ boundaries. Thus, to arrive at these solutions, some algorithms evolve some good infeasible solutions with less constraint violation along with their feasible counterparts during the evolutionary process (e.g., see [4]–[6]). The primary purpose of evolving infeasible solutions in the search procedure is to utilize the information they transport. As EAs are stochastic search and optimization methods, rejecting infeasible individuals might lead the EA being stuck in local optima, particularly in problems with disconnected search space [7], [8]. Moreover, in some highly constrained optimization problems, it could be a demanding problem to find a single feasible solution [9], [10]. Therefore, constraint handling techniques used in multiobjective optimization (MOO) can be mainly distinguished by knowing how infeasible solutions are mixed up and evolved in the evolutionary process. 1This definition of domination is for minimization. All the inequalities should be reversed if the goal is to maximize the objectives in (1). “dominate” means “be better than”. 647 | P a g e www.ijacsa.thesai.org (IJACSA) International Journal of Advanced Computer Science and Applications, Vol. 7, No. 2, 2016 In [1], we introduced a new threshold based penalty function in the replacement and update scheme of MOEA/DDE [11], an improved version of MOEA/D [12], to penalize infeasible solutions. The threshold is adaptively adjusted by using the minimum and maximum constraint violation in the neighborhood of a solution. The infeasible solutions with constraint violation less than the threshold are less penalized than the ones with constraint violation greater than the threshold. As a result, we expect that some good infeasible solutions with less constraint violation will have a chance to evolve in the evolutionary process. The some preliminary experimental results, presented in [1], have proven the capability of the proposed algorithm for solving CF-series [13] test instances. In this paper, we present detailed experimental results and comment on the pitfalls of the proposed algorithm. The rest of this paper is organized as follows. Section II presents some basic concepts and the two versions of the proposed penalty function. Section III briefly introduces MOEA/D and the modified algorithmic framework of MOEA/D-DE. Section IV discusses the experimental settings. Section V presents and discusses experimental results on CF-series [13] test instances. Section VI compares our experimental results with the three best performers [14]–[16] of CEC 2009 MOEA competition. Finally, Section VII concludes this paper with an outline of the work carried out. II. BASIC CONCEPTS ANDTHE PROPOSEDPENALTY FUNCTION A. Degree of Constraint Violation The degree of constraint violation of a solution x ∈ S can be defined as [1], [3]:


I. INTRODUCTION
In this paper, we consider the following constrained multiobjective optimization problem (CMOP) [1]: Minimize F (x) = (f 1 (x), f 2 (x), . . ., f m (x)) T ; Subject to g j (x) ≥ 0, j = 1, . . ., p; where x = (x 1 , . . ., x n ) T ∈ R n is an n dimensional vector of decision variables, F is the objective vector function that consists of m real-valued objective functions, and g i (x) ≥ 0 are inequality constraints.The objective and constraint functions, f i 's and g j 's, could be linear or non linear real-valued functions.l k and u k are the lower and upper bounds (called bound constraints) of x k , k = 1, . . ., n, respectively, which define the search region S = {x = (x 1 , . . ., x n ) T | l k ≤ x k ≤ u k , k = 1, . . ., n}.
A solution x ∈ S is called a feasible solution, if it satisfies all the inequality constraints in (1).The set of all feasible solutions is called the feasible region.Mathematically, we can write: However, If a solution is not feasible, we call it infeasible.The set of all infeasible solutions is called the infeasible region.
The feasible attainable objective set (AOS) can be defined as {F (x)|x ∈ F }.
Since the objectives in (1) more often contradict each other, so it is hard to find a single solution in F that could minimize all the objectives at the same time.Instead, one looks for a set of optimal compromising/tradeoff feasible solutions.The best tradeoffs among the objectives can be defined in terms of Pareto-optimality [2], [3].
A solution x is said to Pareto-dominate or simply dominate another solution y, mathematically denoted as x y, if f i (x) ≤ f i (y), ∀i = 1, . . ., m and f j (x) < f j (y) for at least one j ∈ {1, . . ., m} 1 .This definition of domination is sometimes referred to as a weak dominance relation.A solution x * ∈ F is Pareto-optimal to (1) if there is no solution x ∈ F such that F (x) F (x * ).F (x * ) is then called a Pareto-optimal (objective) vector.The set of all Paretooptimal solutions is called the Pareto Set (PS) in the decision space and Pareto Front (PF) in the objective space [2].
In the majority of constrained optimization problems, the optimal solutions lie on the constraints' boundaries.Thus, to arrive at these solutions, some algorithms evolve some good infeasible solutions with less constraint violation along with their feasible counterparts during the evolutionary process (e.g., see [4]- [6]).The primary purpose of evolving infeasible solutions in the search procedure is to utilize the information they transport.As EAs are stochastic search and optimization methods, rejecting infeasible individuals might lead the EA being stuck in local optima, particularly in problems with disconnected search space [7], [8].Moreover, in some highly constrained optimization problems, it could be a demanding problem to find a single feasible solution [9], [10].Therefore, constraint handling techniques used in multiobjective optimization (MOO) can be mainly distinguished by knowing how infeasible solutions are mixed up and evolved in the evolutionary process.
In [1], we introduced a new threshold based penalty function in the replacement and update scheme of MOEA/D-DE [11], an improved version of MOEA/D [12], to penalize infeasible solutions.The threshold is adaptively adjusted by using the minimum and maximum constraint violation in the neighborhood of a solution.The infeasible solutions with constraint violation less than the threshold are less penalized than the ones with constraint violation greater than the threshold.As a result, we expect that some good infeasible solutions with less constraint violation will have a chance to evolve in the evolutionary process.The some preliminary experimental results, presented in [1], have proven the capability of the proposed algorithm for solving CF-series [13] test instances.In this paper, we present detailed experimental results and comment on the pitfalls of the proposed algorithm.
The rest of this paper is organized as follows.Section II presents some basic concepts and the two versions of the proposed penalty function.Section III briefly introduces MOEA/D and the modified algorithmic framework of MOEA/D-DE.Section IV discusses the experimental settings.Section V presents and discusses experimental results on CF-series [13] test instances.Section VI compares our experimental results with the three best performers [14]- [16] of CEC 2009 MOEA competition.Finally, Section VII concludes this paper with an outline of the work carried out.

A. Degree of Constraint Violation
The degree of constraint violation of a solution x ∈ S can be defined as [1], [3]: Obviously, if V (x) = 0, x is feasible; otherwise, it is infeasible.

C. The Proposed Penalty Function
The proposed penalty function uses a threshold value, τ for dynamically controlling the amount of penalty.
Suppose MOEA/D [12] decomposes the MOP into N subproblems.At each generation, MOEA/D retains N solutions x 1 , . . ., x N , where x i is the current solution to subproblem i.Let P be the mating and update range set in MOEA/D.Then define [1]: where V (x i ) is the degree of constraint violation of solution x i .
The threshold value, τ is defined as [1]: where the parameter s controls the threshold value.In [1], we used s = 0.3.
Our suggested penalty function encourages the algorithm to search the feasible region and the infeasible region near the feasible region.It is defined in the following two different ways: For i = 1, . . ., m where s 1 and s 2 are two scaling parameters with s 1 << s 2 .In the penalty functions, the penalty increases sharply when V (x) exceeds the threshold.This is realized by scaling the degree of constraint violation, V (x) of an infeasible solution by relatively high value of parameter s 2 than parameter s 1 in our penalty function formulations.In Eq. 7, the penalty is added to individual objective function values of an infeasible solution, while in Eq. 8, it is added directly to Tchebycheff aggregation function value of an infeasible solution.Furthermore, in [1], we tested Eq. 8 only.

III. MULTIOBJECTIVE EVOLUTIONARY ALGORITHM BASED ON DECOMPOSITION
Zhang and Li [12] suggested a simple yet efficient MOEA, multiobjective evolutionary algorithm based on decomposition (MOEA/D).MOEA/D approximates the PF by explicitly decomposing an MOP into several single objective optimization subproblems.These subproblems are then optimized concurrently and collaboratively by evolving population of solutions.An EA is employed for this purpose.The Euclidean distances between the aggregation coefficient vectors of these subproblems are calculated to identify the neighborhood of each subproblem.The information gathered from the neighboring subproblems is then used to optimize a subproblem.
In this work, we employed the penalty functions defined by Eqs. 7, 8 in the update scheme of MOEA/D-DE [11], one of the efficient versions of MOEA/D to solve CF-series [13] test instances.This resulted in a new algorithm, denoted by CMOEA/D-DE-ATP (For details of CMOEA/D-DE, please see [18] ).The pseudo-code of the modified update scheme is given in Algorithm 1.
Algorithm 1 Pseudo-code of the update scheme of CMOEA/D-DE-ATP.
1: Each new child solution y updates n r solutions from the set P of its neighboring solutions as follows: 2: Set c = 0 and then do the following: 3: if c = n r or P is empty then Randomly pick an index j from P ; 7: Compute the Tchebycheff aggregation function values of y and x j with the new objective values of Eq. 7 (or the new aggregation function values of y and x j with Eq. 8); 8: if g te (y|λ j , z) ≤ g te (x j |λ j , z) (or g te p (y|λ j , z) ≤ g te p (x j |λ j , z) ) then 9: Remove j from P and go to step 3; 12: end if

IV. EXPERIMENTAL SETTINGS
In our experiments, we use the same parameters' settings and weight vectors' selection criteria as is used in [13].Further, we use statistics of the inverted generational distance metric (IGD-metric) [12], [19] for comparing results on CF-series test instances, CF1-CF10.Also, the set coverage metric (SCmetric) [12] is used to compare the nondominated solutions obtained by different algorithms.Unless otherwise stated, we will use Eq. 6 with s = 0.7 and Eqs. 7, 8 with s 1 = 0.01 and s 2 = 20 in all experiments.Table II presents the average set coverage between the nondominated solutions of CMOEA/D-DE-ATP with Eq. 7 and Eq. 8.The results of this table reveal that, in terms of the SC-metric, the nondominated solutions found by CMOEA/D-DE-ATP with Eq. 8 are better than those obtained with Eq. 7 for test instances CF1, CF4-CF6, but are worse for test instances CF8 and CF10 vice versa.The table also shows that the nondominated solutions acquired from CMOEA/D-DE-ATP with both Eqs. 7, 8 are same for test instances CF2, CF3, CF7, and CF9.However, looking at the results of this table, it can be inferred that the performance of CMOEA/D-DE-ATP is comparable with both Eqs. 7, 8, as there is no big difference in the SC-metric values.

V. EXPERIMENTAL RESULTS
Figures 1 and 2 show, in the objective space, the distributions of the 100 and 150 nondominated population members for the seven 2-objective, CF1-CF7, and the three 3-objective, CF8-CF10, CF-series test instances.These solutions are selected based on the criteria as mentioned in [13] from the final population of the run with the best (i.e., lowest) IGD-metric value among the 30 independent runs.These figures also show Fig. 1: Plots of the nondominated front with the best IGD value and all the 30 final nondominated fronts found by CMOEA/D-DE-ATP when using Eq.7 (columns 1 and 3) and Eq. 8 (columns 2 and 4) for CF1-CF6.Real PF Obtained PFs Fig. 2: Plots of the nondominated front with the best IGD value and all the 30 final nondominated fronts found by CMOEA/D-DE-ATP when using Eq.7 (columns 1 and 3) and Eq. 8 (columns 2 and 4) for CF7-CF10.all the 30 final nondominated fronts of these selected 100 and 150 nondominated solutions.
It is very clear from these figures that CMOEA/D-DE-ATP with both Eqs. 7, 8 found good approximations for the four 2-objective, CF1, CF2, CF4, and CF6, and two 3-objective, CF8, CF9 test instances.However, it performed poorly on test instances CF3, CF5, CF7, and CF10.It is also apparent from the plots of 30 nondominated fronts of test instances CF4 and CF6 that CMOEA/D-DE-ATP fails to find the whole PF in some runs.
The PF of CF3 is concave and discontinuous.Therefore, it could be hard for the algorithm than all other 2-objective test instances.Although the PFs of CF4 and CF5, CF6 and CF7, and CF9 and CF10 are identical, the poor performance of CMOEA/D-DE-ATP on test instances CF5, CF7 and CF10 could be due to the presence of harder objective and constraint functions in these test instances than test instances CF4, CF6, and CF9.
Figure 3 shows the evolution of the average IGD-metric values versus function evaluations of the nondominated so-lutions in the current population.This figure shows that CMOEA/D-DE-ATP with both Eqs. 7, 8 converges at the same rate in terms of IGD-metric values for six CF-series test instances CF1, CF3, CF4, CF6, CF8 and CF9.However, it converges slightly faster in terms of IGD-metric values for the other four CF-series test instances CF2, CF5, CF7, and CF10 with Eq. 7 than with Eq. 8.
Figure 4 depicts the average generation feasibility versus generations' graphs.This figure demonstrates that CMOEA/D-DE-ATP with both Eqs. 7, 8 approaches to the feasible regions at the same rate for test instances CF1-CF6.It converges slower to the feasible regions with Eq. 7 than with Eq. 8 for the test instances CF7-CF10.This permits further exploration of the infeasible regions near the PF and could be one of the reasons for the better performance of CMOEA/D-DE-ATP with Eq. 7 on the three 3-objective test instances, CF8-CF10.
As it can be seen from Figure 4 that 50 % or more of the initial populations for test instances CF1, CF2, CF4 and CF5 are feasible.These feasible solutions are propagated in the subsequent generations by the replacement and update scheme of the algorithm and thus produce better feasible solutions due to the DE operator.Furthermore, the feasibility ratio becomes 1 after the initial 30 to 40 generations for these test instances.
The reason for the quick convergence to the feasible region is the higher adopted update number of neighboring parent solutions (as in our settings n r = 6 when T = 60 and n r = 10 when T = 100) that are replaced by a better child solution in the update scheme of the algorithm.This speedy convergence to the feasible region is good for test instances like CF1, CF2, and CF4, but it causes problems for harder test instance like CF5.The PF of CF5 is a piecewise continuous curve with three pieces like CF4, but its objective and constraint functions are quite different and harder than CF4.
On the other hand, about 25 % or below of the initial populations for test instances CF6-CF10 is feasible.Here, the proposed constraint handling technique has more chances to evolve better infeasible solutions during the evolutionary process.Particularly, in the two 3-objective test instances CF8 and CF9, the average feasibility ratio at the last generations of the algorithmic runs is 0.6 and about 0.9, respectively (see Figure 4).This way the infeasible regions near the feasibility boundaries in these two instances are well explored and could be a reason for the better performance of the algorithm on these two instances.
Moreover, in test instance CF6, the feasibility ratio of CMOEA/D-DE-ATP with both Eqs. 7, 8 becomes 1 after the initial 40 generations, while in test instance CF7, it takes 200 generations of CMOEA/D-DE-ATP with Eq. 7 and 50 gener- ations with Eq. 8 to become 1.Thus, the quick convergence to the feasible region is good in case of CF6, but could be a reason for the poor performance of CMOEA/D-DE-ATP with both Eqs. 7, 8 in case of CF7.
In test instance CF10, the feasibility ratio takes about 140 generations of CMOEA/D-DE-ATP with Eq. 7 and about 80 generations of CMOEA/D-DE-ATP with Eq. 8 to become 1.Again, the less exploration of the infeasible regions could be the reason for the poor performance of CMOEA/D-DE-ATP on test instance CF10.

VI. COMPARISON WITH THE THREE BEST PERFORMERS OF CEC 2009 MOEA COMPETITION
In this section, we compare the results of CMOEA/D-DE-ATP with Eqs. 7, 8 with the three best performers [14]- [16] in CEC 2009 MOEA competition on the CF-series test instances.

VII. CONCLUSIONS
A penalty function that penalizes infeasible solutions based on an adaptive threshold value has been introduced into the update and replacement scheme of MOEA/D-DE.This resulted in a new algorithm, CMOEA/D-DE-ATP for CMOO.The proposed penalty function is presented in two forms given by Eqs. 7, 8.The performance of CMOEA/D-DE-ATP is tested on CF-series test instances in terms of the values of IGD-metric and SC-metric.
From the experimental results in this paper, we can make the following conclusions.
• Overall, CMOEA/D-DE-ATP produced better results with the proposed penalty function defined by Eq. 7 than when it is defined by Eq. 8.That is, it is better to add the penalty to individual objective function values before calculating the aggregation function values than directly adding the penalty to aggregation function values of an infeasible solution for better performance achievement on CF-series test instances.

TABLE I :
THE IGD-METRIC STATISTICS OF CMOEA/D-DE-ATP USING EQS. 7, 8. THE RESULTS IN BOLDFACE INDICATE THE BETTER RESULTS; IF NOT, THEY ARE IDENTICAL.by CMOEA/D-DE-ATP with Eqs. 7, 8.These statistics are based on 30 independent runs.As it can be seen from this table that CMOEA/D-DE-ATP can find better best values with Eq. 7 for one 2-objective, CF6 and three 3-objective, CF8-CF10 and with Eq. 8 for four 2-objective, CF2, CF4, CF5, and CF7 test instances.The best values for test instance CF1 are identical.This table also shows that CMOEA/D-DE-ATP with both Eqs. 7, 8 performs similarly on test instance CF3.However, improved mean and st.dev.values can be found when CMOEA/D-DE-ATP employs Eq. 7 for most of the test instances.In particular, the improved performance can be seen for the three 3-objective test instances, CF8-CF10.This suggests that adding the penalty to individual objective function values as is done in Eq. 7 before calculating the aggregation function values is a good choice for better performance on CF-series test instances.It can also be seen from Table I that CMOEA/D-DE-ATP with both Eqs. 7, 8 finds small values for the mean of IGDmetric on CF1, CF2, CF4, CF6, CF8, CF9.Empirically, these results illustrate that the final nondominated solutions found by CMOEA/D-DE-ATP for these test instances approximate the PF very well in a sense.
Table I presents the best (i.e., lowest), mean, and standard deviation of the IGD-metric values for CF-series test instances found

TABLE II :
THE AVERAGE SET COVERAGE BETWEEN CMOEA/D-DE-ATP WITH EQ. 7 AND WITH EQ. 8 ON CF-SERIES TEST INSTANCES.THE RESULTS IN BOLD-FACE INDICATE THE BETTER RESULTS; IF NOT, THEY ARE IDENTICAL.

Table
[16]compares the best (i.e., lowest), mean, and standard deviation values of the IGD-metric obtained from our algorithm, CMOEA/D-DE-ATP with Eqs. 7, 8, and the three best performers[14]-[16]in CEC 2009 MOEA competition for the CF-series test instances.The table clearly shows that CMOEA/D-DE-ATP has found the best (i.e., lowest) IGDmetric values for four test instances CF1, CF6, CF8 and CF9 and the second best value for one test instance CF3 with Eq.

•
The comparison of CMOEA/D-DE-ATP with the three best performers in CEC 2009 special session and competition indicated that CMOEA/D-DE-ATP has found the best (i.e., lowest) IGD-metric values for five test instances CF1, CF4, CF6, CF8 and CF9 and the second best values for two test instances CF2 and CF3.In particular, our algorithm overall found better statistics for tests instances CF1, CF8, and CF9.