NSGA3理解(NSGA3算法及其MATLAB版本实现)
GitHub上的C++版本:https://github.com/adanjoga/nsga3
晓风从platEMO中摘出来的NSGA-III算法代码:
NSGAIII_main.m
1 clc,clear 2 global N M D PopCon name gen 3 4 N = 400; % 种群个数 5 M = 3; % 目标个数 6 name = 'DTLZ1'; % 测试函数选择,目前只有:DTLZ1、DTLZ2、DTLZ3 7 gen = 500; %迭代次数 8 %% Generate the reference points and random population 9 [Z,N] = UniformPoint(N,M); % 生成一致性参考解 10 [res,Population,PF] = funfun(); % 生成初始种群与目标值 11 Pop_objs = CalObj(Population); % 计算适应度函数值 12 Zmin = min(Pop_objs(all(PopCon<=0,2),:),[],1); %求理想点,其实PopCon是处理有约束问题的,这里并没有用到 13 14 %% Optimization 15 for i = 1:gen 16 MatingPool = TournamentSelection(2,N,sum(max(0,PopCon),2)); 17 Offspring = GA(Population(MatingPool,:)); 18 Offspring_objs = CalObj(Offspring); 19 Zmin = min([Zmin;Offspring_objs],[],1); 20 Population = EnvironmentalSelection([Population;Offspring],N,Z,Zmin); 21 Popobj = CalObj(Population); 22 if(M<=3) 23 plot3(Popobj(:,1),Popobj(:,2),Popobj(:,3),'ro') 24 title(num2str(i)); 25 drawnow 26 end 27 end 28 if(M<=3) 29 hold on 30 plot3(PF(:,1),PF(:,2),PF(:,3),'g*') 31 else 32 for i=1:size(Popobj,1) 33 plot(Popobj(i,:)) 34 hold on 35 end 36 end 37 %%IGD 38 score = IGD(Popobj,PF)
UniformPoint.m
1 function [W,N] = UniformPoint(N,M) 2 %UniformPoint - Generate a set of uniformly distributed points on the unit 3 %hyperplane. 4 % 5 % [W,N] = UniformPoint(N,M) returns approximately N uniformly distributed 6 % points with M objectives on the unit hyperplane. 7 % 8 % Due to the requirement of uniform distribution, the number of points 9 % cannot be arbitrary, and the number of points in W may be slightly 10 % smaller than the predefined size N. 11 % 12 % Example: 13 % [W,N] = UniformPoint(275,10) 14 H1 = 1; 15 while nchoosek(H1+M-1,M-1) <= N 16 H1 = H1 + 1; 17 end 18 H1=H1 - 1; 19 W = nchoosek(1:H1+M-1,M-1) - repmat(0:M-2,nchoosek(H1+M-1,M-1),1) - 1;%nchoosek(v,M),v是一个向量,返回一个矩阵,该矩阵的行由每次从v中的M个元素取出k个取值的所有可能组合构成。比如:v=[1,2,3],n=2,输出[1,2;1,3;2,3] 20 W = ([W,zeros(size(W,1),1)+H1]-[zeros(size(W,1),1),W])/H1; 21 if H1 < M 22 H2 = 0; 23 while nchoosek(H1+M-1,M-1)+nchoosek(H2+M-1,M-1) <= N 24 H2 = H2 + 1; 25 end 26 H2 = H2 - 1; 27 if H2 > 0 28 W2 = nchoosek(1:H2+M-1,M-1) - repmat(0:M-2,nchoosek(H2+M-1,M-1),1) - 1; 29 W2 = ([W2,zeros(size(W2,1),1)+H2]-[zeros(size(W2,1),1),W2])/H2; 30 W = [W;W2/2+1/(2*M)]; 31 end 32 end 33 W = max(W,1e-6); 34 N = size(W,1); 35 end
funfun.m
1 function [PopObj,PopDec,P] = funfun() 2 3 global M D lower upper encoding N PopCon cons cons_flag name 4 5 %% Initialization 6 D = M + 4; 7 lower = zeros(1,D); 8 upper = ones(1,D); 9 encoding = 'real'; 10 switch encoding 11 case 'binary' 12 PopDec = randi([0,1],N,D); 13 case 'permutation' 14 [~,PopDec] = sort(rand(N,D),2); 15 otherwise 16 PopDec = unifrnd(repmat(lower,N,1),repmat(upper,N,1)); 17 end 18 cons = zeros(size(PopDec,1),1); 19 cons_flag = 1; 20 PopCon = cons; 21 %% Calculate objective values 22 switch name 23 case 'DTLZ1' 24 g = 100*(D-M+1+sum((PopDec(:,M:end)-0.5).^2-cos(20.*pi.*(PopDec(:,M:end)-0.5)),2)); 25 PopObj = 0.5*repmat(1+g,1,M).*fliplr(cumprod([ones(N,1),PopDec(:,1:M-1)],2)).*[ones(N,1),1-PopDec(:,M-1:-1:1)]; 26 P = UniformPoint(N,M)/2; 27 case 'DTLZ2' 28 g = sum((PopDec(:,M:end)-0.5).^2,2); 29 PopObj = repmat(1+g,1,M).*fliplr(cumprod([ones(size(g,1),1),cos(PopDec(:,1:M-1)*pi/2)],2)).*[ones(size(g,1),1),sin(PopDec(:,M-1:-1:1)*pi/2)]; 30 P = UniformPoint(N,M); 31 P = P./repmat(sqrt(sum(P.^2,2)),1,M); 32 case 'DTLZ3' 33 g = 100*(D-M+1+sum((PopDec(:,M:end)-0.5).^2-cos(20.*pi.*(PopDec(:,M:end)-0.5)),2)); 34 PopObj = repmat(1+g,1,M).*fliplr(cumprod([ones(N,1),cos(PopDec(:,1:M-1)*pi/2)],2)).*[ones(N,1),sin(PopDec(:,M-1:-1:1)*pi/2)]; 35 P = UniformPoint(N,M); 36 P = P./repmat(sqrt(sum(P.^2,2)),1,M); 37 end 38 %% Sample reference points on Pareto front 39 %P = UniformPoint(N,obj.Global.M)/2; 40 end
CalObj.m
1 function PopObj = CalObj(PopDec) 2 global M D name 3 NN = size(PopDec,1); 4 switch name 5 case 'DTLZ1' 6 g = 100*(D-M+1+sum((PopDec(:,M:end)-0.5).^2-cos(20.*pi.*(PopDec(:,M:end)-0.5)),2)); 7 PopObj = 0.5*repmat(1+g,1,M).*fliplr(cumprod([ones(NN,1),PopDec(:,1:M-1)],2)).*[ones(NN,1),1-PopDec(:,M-1:-1:1)]; 8 case 'DTLZ2' 9 g = sum((PopDec(:,M:end)-0.5).^2,2); 10 PopObj = repmat(1+g,1,M).*fliplr(cumprod([ones(size(g,1),1),cos(PopDec(:,1:M-1)*pi/2)],2)).*[ones(size(g,1),1),sin(PopDec(:,M-1:-1:1)*pi/2)]; 11 case 'DTLZ3' 12 g = 100*(D-M+1+sum((PopDec(:,M:end)-0.5).^2-cos(20.*pi.*(PopDec(:,M:end)-0.5)),2)); 13 PopObj = repmat(1+g,1,M).*fliplr(cumprod([ones(NN,1),cos(PopDec(:,1:M-1)*pi/2)],2)).*[ones(NN,1),sin(PopDec(:,M-1:-1:1)*pi/2)]; 14 end 15 end
TournamentSelection.m
1 function index = TournamentSelection(K,N,varargin) 2 %TournamentSelection - Tournament selection. 3 % 4 % P = TournamentSelection(K,N,fitness1,fitness2,...) returns the indices 5 % of N solutions by K-tournament selection based on their fitness values. 6 % In each selection, the candidate having the MINIMUM fitness1 value will 7 % be selected; if more than one candidates have the same minimum value of 8 % fitness1, then compare their fitness2 values, and so on. 9 % 10 % Example: 11 % P = TournamentSelection(2,100,FrontNo) 12 13 %------------------------------- Copyright -------------------------------- 14 % Copyright (c) 2018-2019 BIMK Group. You are free to use the PlatEMO for 15 % research purposes. All publications which use this platform or any code 16 % in the platform should acknowledge the use of "PlatEMO" and reference "Ye 17 % Tian, Ran Cheng, Xingyi Zhang, and Yaochu Jin, PlatEMO: A MATLAB platform 18 % for evolutionary multi-objective optimization [educational forum], IEEE 19 % Computational Intelligence Magazine, 2017, 12(4): 73-87". 20 %-------------------------------------------------------------------------- 21 22 varargin = cellfun(@(S) reshape(S,[],1),varargin,'UniformOutput',false); 23 [~,rank] = sortrows([varargin{:}]); 24 [~,rank] = sort(rank); 25 Parents = randi(length(varargin{1}),K,N); 26 [~,best] = min(rank(Parents),[],1); 27 index = Parents(best+(0:N-1)*K); 28 end
GA.m
1 function Offspring = GA(Parent,Parameter) 2 global encoding lower upper 3 %GA - Genetic operators for real, binary, and permutation based encodings. 4 % 5 % Off = GA(P) returns the offsprings generated by genetic operators, 6 % where P1 is a set of parents. If P is an array of INDIVIDUAL objects, 7 % then Off is also an array of INDIVIDUAL objects; while if P is a matrix 8 % of decision variables, then Off is also a matrix of decision variables, 9 % i.e., the offsprings will not be evaluated. P is split into two subsets 10 % P1 and P2 with the same size, and each object/row of P1 and P2 is used 11 % to generate TWO offsprings. Different operators are used for real, 12 % binary, and permutation based encodings, respectively. 13 % 14 % Off = GA(P,{proC,disC,proM,disM}) specifies the parameters of 15 % operators, where proC is the probabilities of doing crossover, disC is 16 % the distribution index of simulated binary crossover, proM is the 17 % expectation of number of bits doing mutation, and disM is the 18 % distribution index of polynomial mutation. 19 % 20 % Example: 21 % Off = GA(Parent) 22 % Off = GA(Parent.decs,{1,20,1,20}) 23 % 24 % See also GAhalf 25 26 %------------------------------- Reference -------------------------------- 27 % [1] K. Deb, K. Sindhya, and T. Okabe, Self-adaptive simulated binary 28 % crossover for real-parameter optimization, Proceedings of the 9th Annual 29 % Conference on Genetic and Evolutionary Computation, 2007, 1187-1194. 30 % [2] K. Deb and M. Goyal, A combined genetic adaptive search (GeneAS) for 31 % engineering design, Computer Science and informatics, 1996, 26: 30-45. 32 % [3] L. Davis, Applying adaptive algorithms to epistatic domains, 33 % Proceedings of the International Joint Conference on Artificial 34 % Intelligence, 1985, 162-164. 35 % [4] D. B. Fogel, An evolutionary approach to the traveling salesman 36 % problem, Biological Cybernetics, 1988, 60(2): 139-144. 37 %------------------------------- Copyright -------------------------------- 38 % Copyright (c) 2018-2019 BIMK Group. You are free to use the PlatEMO for 39 % research purposes. All publications which use this platform or any code 40 % in the platform should acknowledge the use of "PlatEMO" and reference "Ye 41 % Tian, Ran Cheng, Xingyi Zhang, and Yaochu Jin, PlatEMO: A MATLAB platform 42 % for evolutionary multi-objective optimization [educational forum], IEEE 43 % Computational Intelligence Magazine, 2017, 12(4): 73-87". 44 %-------------------------------------------------------------------------- 45 46 %% Parameter setting 47 if nargin > 1 48 [proC,disC,proM,disM] = deal(Parameter{:}); 49 else 50 [proC,disC,proM,disM] = deal(1,20,1,20); 51 end 52 if isa(Parent(1),'INDIVIDUAL') 53 calObj = true; 54 Parent = Parent.decs; 55 else 56 calObj = false; 57 end 58 Parent1 = Parent(1:floor(end/2),:); 59 Parent2 = Parent(floor(end/2)+1:floor(end/2)*2,:); 60 [N,D] = size(Parent1); 61 62 switch encoding 63 case 'binary' 64 %% Genetic operators for binary encoding 65 % One point crossover 66 k = repmat(1:D,N,1) > repmat(randi(D,N,1),1,D); 67 k(repmat(rand(N,1)>proC,1,D)) = false; 68 Offspring1 = Parent1; 69 Offspring2 = Parent2; 70 Offspring1(k) = Parent2(k); 71 Offspring2(k) = Parent1(k); 72 Offspring = [Offspring1;Offspring2]; 73 % Bitwise mutation 74 Site = rand(2*N,D) < proM/D; 75 Offspring(Site) = ~Offspring(Site); 76 case 'permutation' 77 %% Genetic operators for permutation based encoding 78 % Order crossover 79 Offspring = [Parent1;Parent2]; 80 k = randi(D,1,2*N); 81 for i = 1 : N 82 Offspring(i,k(i)+1:end) = setdiff(Parent2(i,:),Parent1(i,1:k(i)),'stable'); 83 Offspring(i+N,k(i)+1:end) = setdiff(Parent1(i,:),Parent2(i,1:k(i)),'stable'); 84 end 85 % Slight mutation 86 k = randi(D,1,2*N); 87 s = randi(D,1,2*N); 88 for i = 1 : 2*N 89 if s(i) < k(i) 90 Offspring(i,:) = Offspring(i,[1:s(i)-1,k(i),s(i):k(i)-1,k(i)+1:end]); 91 elseif s(i) > k(i) 92 Offspring(i,:) = Offspring(i,[1:k(i)-1,k(i)+1:s(i)-1,k(i),s(i):end]); 93 end 94 end 95 otherwise 96 %% Genetic operators for real encoding 97 % Simulated binary crossover 98 beta = zeros(N,D); 99 mu = rand(N,D); 100 beta(mu<=0.5) = (2*mu(mu<=0.5)).^(1/(disC+1)); 101 beta(mu>0.5) = (2-2*mu(mu>0.5)).^(-1/(disC+1)); 102 beta = beta.*(-1).^randi([0,1],N,D); 103 beta(rand(N,D)<0.5) = 1; 104 beta(repmat(rand(N,1)>proC,1,D)) = 1; 105 Offspring = [(Parent1+Parent2)/2+beta.*(Parent1-Parent2)/2 106 (Parent1+Parent2)/2-beta.*(Parent1-Parent2)/2]; 107 % Polynomial mutation 108 Lower = repmat(lower,2*N,1); 109 Upper = repmat(upper,2*N,1); 110 Site = rand(2*N,D) < proM/D; 111 mu = rand(2*N,D); 112 temp = Site & mu<=0.5; 113 Offspring = min(max(Offspring,Lower),Upper); 114 Offspring(temp) = Offspring(temp)+(Upper(temp)-Lower(temp)).*((2.*mu(temp)+(1-2.*mu(temp)).*... 115 (1-(Offspring(temp)-Lower(temp))./(Upper(temp)-Lower(temp))).^(disM+1)).^(1/(disM+1))-1); 116 temp = Site & mu>0.5; 117 Offspring(temp) = Offspring(temp)+(Upper(temp)-Lower(temp)).*(1-(2.*(1-mu(temp))+2.*(mu(temp)-0.5).*... 118 (1-(Upper(temp)-Offspring(temp))./(Upper(temp)-Lower(temp))).^(disM+1)).^(1/(disM+1))); 119 end 120 if calObj 121 Offspring = INDIVIDUAL(Offspring); 122 end 123 end
EnvironmentalSelection.m
1 function Population = EnvironmentalSelection(Population,N,Z,Zmin) 2 global PopCon 3 % The environmental selection of NSGA-III 4 5 %------------------------------- Copyright -------------------------------- 6 % Copyright (c) 2018-2019 BIMK Group. You are free to use the PlatEMO for 7 % research purposes. All publications which use this platform or any code 8 % in the platform should acknowledge the use of "PlatEMO" and reference "Ye 9 % Tian, Ran Cheng, Xingyi Zhang, and Yaochu Jin, PlatEMO: A MATLAB platform 10 % for evolutionary multi-objective optimization [educational forum], IEEE 11 % Computational Intelligence Magazine, 2017, 12(4): 73-87". 12 %-------------------------------------------------------------------------- 13 14 if isempty(Zmin) 15 Zmin = ones(1,size(Z,2)); 16 end 17 18 %% Non-dominated sorting 19 Population_objs = CalObj(Population); 20 [FrontNo,MaxFNo] = NDSort(Population_objs,PopCon,N); 21 Next = FrontNo < MaxFNo; 22 23 %% Select the solutions in the last front 24 Last = find(FrontNo==MaxFNo); 25 Choose = LastSelection(Population_objs(Next,:),Population_objs(Last,:),N-sum(Next),Z,Zmin); 26 Next(Last(Choose)) = true; 27 % Population for next generation 28 Population = Population(Next,:); 29 end 30 31 function Choose = LastSelection(PopObj1,PopObj2,K,Z,Zmin) 32 % Select part of the solutions in the last front 33 34 PopObj = [PopObj1;PopObj2] - repmat(Zmin,size(PopObj1,1)+size(PopObj2,1),1); 35 [N,M] = size(PopObj); 36 N1 = size(PopObj1,1); 37 N2 = size(PopObj2,1); 38 NZ = size(Z,1); 39 40 %% Normalization 41 % Detect the extreme points 42 Extreme = zeros(1,M); 43 w = zeros(M)+1e-6+eye(M); 44 for i = 1 : M 45 [~,Extreme(i)] = min(max(PopObj./repmat(w(i,:),N,1),[],2)); 46 end 47 % Calculate the intercepts of the hyperplane constructed by the extreme 48 % points and the axes 49 Hyperplane = PopObj(Extreme,:)\ones(M,1);%A的逆矩阵是 A1,则 B/A = B*A1,A\B = A1*B 50 a = 1./Hyperplane; 51 if any(isnan(a))%若a的元素为NaN(非数值),在对应位置上返回逻辑1(真),否则返回逻辑0(假) 52 a = max(PopObj,[],1)'; 53 end 54 % Normalization 55 %a = a-Zmin'; 56 PopObj = PopObj./repmat(a',N,1);%如果a、b是矩阵,a./b就是a、b中对应的每个元素相除,得到一个新的矩阵; 57 58 %% Associate each solution with one reference point 59 % Calculate the distance of each solution to each reference vector 60 Cosine = 1 - pdist2(PopObj,Z,'cosine');%cosine’ 针对向量,夹角余弦距离Cosine distance(‘cosine’)余弦相似度= 1-余弦距离 61 62 %杰卡德距离Jaccard distance(‘jaccard’)Jaccard距离常用来处理仅包含非对称的二元(0-1)属性的对象。很显然,Jaccard距离不关心0-0匹配[1]。 63 %夹角余弦距离Cosine distance(‘cosine’)与Jaccard距离相比,Cosine距离不仅忽略0-0匹配,而且能够处理非二元向量,即考虑到变量值的大小。 64 %对这两者,距离与相似度和为一。 65 %https://www.cnblogs.com/chaosimple/archive/2013/06/28/3160839.html 66 67 Distance = repmat(sqrt(sum(PopObj.^2,2)),1,NZ).*sqrt(1-Cosine.^2); 68 %Distance = pdist2(PopObj,Z,'euclidean');%上面注释的两个语句也可以哦 69 % Associate each solution with its nearest reference point 70 [d,pi] = min(Distance',[],1); 71 72 %% Calculate the number of associated solutions except for the last front of each reference point 73 rho = hist(pi(1:N1),1:NZ);%除了最后front的每一个参考点,计算关联解的个数 74 75 %% Environmental selection 76 Choose = false(1,N2); 77 Zchoose = true(1,NZ); 78 % Select K solutions one by one 79 while sum(Choose) < K 80 % Select the least crowded reference point 81 Temp = find(Zchoose); 82 Jmin = find(rho(Temp)==min(rho(Temp)));%所有参考点中与之关联个数最少的那个参考点 83 j = Temp(Jmin(randi(length(Jmin)))); 84 I = find(Choose==0 & pi(N1+1:end)==j); 85 % Then select one solution associated with this reference point 86 if ~isempty(I) 87 if rho(j) == 0 88 [~,s] = min(d(N1+I)); 89 else 90 s = randi(length(I)); 91 end 92 Choose(I(s)) = true; 93 rho(j) = rho(j) + 1; 94 else 95 Zchoose(j) = false; 96 end 97 end 98 end
NDSort.m
1 function [FrontNo,MaxFNo] = NDSort(varargin) 2 %NDSort - Do non-dominated sorting by efficient non-dominated sort. 3 % 4 % FrontNo = NDSort(F,s) does non-dominated sorting on F, where F is the 5 % matrix of objective values of a set of individuals, and s is the number 6 % of individuals to be sorted at least. FrontNo(i) denotes the front 7 % number of the i-th individual. The individuals have not been sorted are 8 % assigned a front number of inf. 9 % 10 % FrontNo = NDSort(F,C,s) does non-dominated sorting based on constrained 11 % domination, where C is the matrix of constraint values of the 12 % individuals. In this case, feasible solutions always dominate 13 % infeasible solutions, and one infeasible solution dominates another 14 % infeasible solution if the former has a smaller overall constraint 15 % violation than the latter. 16 % 17 % In particular, s = 1 indicates finding only the first non-dominated 18 % front, s = size(F,1)/2 indicates sorting only half the population 19 % (which is often used in the algorithm), and s = inf indicates sorting 20 % the whole population. 21 % 22 % [FrontNo,K] = NDSort(...) also returns the maximum front number besides 23 % inf. 24 % 25 % Example: 26 % [FrontNo,MaxFNo] = NDSort(PopObj,1) 27 % [FrontNo,MaxFNo] = NDSort(PopObj,PopCon,inf) 28 29 %------------------------------- Reference -------------------------------- 30 % [1] X. Zhang, Y. Tian, R. Cheng, and Y. Jin, An efficient approach to 31 % nondominated sorting for evolutionary multiobjective optimization, IEEE 32 % Transactions on Evolutionary Computation, 2015, 19(2): 201-213. 33 % [2] X. Zhang, Y. Tian, R. Cheng, and Y. Jin, A decision variable 34 % clustering based evolutionary algorithm for large-scale many-objective 35 % optimization, IEEE Transactions on Evolutionary Computation, 2018, 22(1): 36 % 97-112. 37 %------------------------------- Copyright -------------------------------- 38 % Copyright (c) 2018-2019 BIMK Group. You are free to use the PlatEMO for 39 % research purposes. All publications which use this platform or any code 40 % in the platform should acknowledge the use of "PlatEMO" and reference "Ye 41 % Tian, Ran Cheng, Xingyi Zhang, and Yaochu Jin, PlatEMO: A MATLAB platform 42 % for evolutionary multi-objective optimization [educational forum], IEEE 43 % Computational Intelligence Magazine, 2017, 12(4): 73-87". 44 %-------------------------------------------------------------------------- 45 46 PopObj = varargin{1}; 47 [N,M] = size(PopObj); 48 if nargin == 2 49 nSort = varargin{2}; 50 else 51 PopCon = varargin{2}; 52 nSort = varargin{3}; 53 Infeasible = any(PopCon>0,2); 54 PopObj(Infeasible,:) = repmat(max(PopObj,[],1),sum(Infeasible),1) + repmat(sum(max(0,PopCon(Infeasible,:)),2),1,M); 55 end 56 if M < 5 || N < 500 57 % Use efficient non-dominated sort with sequential search (ENS-SS) 58 [FrontNo,MaxFNo] = ENS_SS(PopObj,nSort); 59 else 60 % Use tree-based efficient non-dominated sort (T-ENS) 61 [FrontNo,MaxFNo] = T_ENS(PopObj,nSort); 62 end 63 end 64 65 function [FrontNo,MaxFNo] = ENS_SS(PopObj,nSort) 66 [PopObj,~,Loc] = unique(PopObj,'rows'); 67 Table = hist(Loc,1:max(Loc)); 68 [N,M] = size(PopObj); 69 [PopObj,rank] = sortrows(PopObj); 70 % if ~(all(PopObj(:) == PopObj1(:))) 71 % disp(2) 72 % end 73 FrontNo = inf(1,N); 74 MaxFNo = 0; 75 while sum(Table(FrontNo<inf)) < min(nSort,length(Loc)) 76 MaxFNo = MaxFNo + 1; 77 for i = 1 : N 78 if FrontNo(i) == inf 79 Dominated = false; 80 for j = i-1 : -1 : 1 81 if FrontNo(j) == MaxFNo 82 m = 2; 83 while m <= M && PopObj(i,m) >= PopObj(j,m) 84 m = m + 1; 85 end 86 Dominated = m > M; 87 if Dominated || M == 2 88 break; 89 end 90 end 91 end 92 if ~Dominated 93 FrontNo(i) = MaxFNo; 94 end 95 end 96 end 97 end 98 FrontNo(rank) = FrontNo; 99 FrontNo = FrontNo(:,Loc); 100 end 101 102 function [FrontNo,MaxFNo] = T_ENS(PopObj,nSort) 103 [PopObj,~,Loc] = unique(PopObj,'rows'); 104 Table = hist(Loc,1:max(Loc)); 105 [N,M] = size(PopObj); 106 [PopObj,rank] = sortrows(PopObj); 107 FrontNo = inf(1,N); 108 MaxFNo = 0; 109 Forest = zeros(1,N); 110 Children = zeros(N,M-1); 111 LeftChild = zeros(1,N) + M; 112 Father = zeros(1,N); 113 Brother = zeros(1,N) + M; 114 [~,ORank] = sort(PopObj(:,2:M),2,'descend'); 115 ORank = ORank + 1; 116 while sum(Table(FrontNo<inf)) < min(nSort,length(Loc)) 117 MaxFNo = MaxFNo + 1; 118 root = find(FrontNo==inf,1); 119 Forest(MaxFNo) = root; 120 FrontNo(root) = MaxFNo; 121 for p = 1 : N 122 if FrontNo(p) == inf 123 Pruning = zeros(1,N); 124 q = Forest(MaxFNo); 125 while true 126 m = 1; 127 while m < M && PopObj(p,ORank(q,m)) >= PopObj(q,ORank(q,m)) 128 m = m + 1; 129 end 130 if m == M 131 break; 132 else 133 Pruning(q) = m; 134 if LeftChild(q) <= Pruning(q) 135 q = Children(q,LeftChild(q)); 136 else 137 while Father(q) && Brother(q) > Pruning(Father(q)) 138 q = Father(q); 139 end 140 if Father(q) 141 q = Children(Father(q),Brother(q)); 142 else 143 break; 144 end 145 end 146 end 147 end 148 if m < M 149 FrontNo(p) = MaxFNo; 150 q = Forest(MaxFNo); 151 while Children(q,Pruning(q)) 152 q = Children(q,Pruning(q)); 153 end 154 Children(q,Pruning(q)) = p; 155 Father(p) = q; 156 if LeftChild(q) > Pruning(q) 157 Brother(p) = LeftChild(q); 158 LeftChild(q) = Pruning(q); 159 else 160 bro = Children(q,LeftChild(q)); 161 while Brother(bro) < Pruning(q) 162 bro = Children(q,Brother(bro)); 163 end 164 Brother(p) = Brother(bro); 165 Brother(bro) = Pruning(q); 166 end 167 end 168 end 169 end 170 end 171 FrontNo(rank) = FrontNo; 172 FrontNo = FrontNo(:,Loc); 173 end
IGD.m
1 function Score = IGD(PopObj,PF) 2 % <metric> <min> 3 % Inverted generational distance 4 5 %------------------------------- Reference -------------------------------- 6 % C. A. Coello Coello and N. C. Cortes, Solving multiobjective optimization 7 % problems using an artificial immune system, Genetic Programming and 8 % Evolvable Machines, 2005, 6(2): 163-190. 9 %------------------------------- Copyright -------------------------------- 10 % Copyright (c) 2018-2019 BIMK Group. You are free to use the PlatEMO for 11 % research purposes. All publications which use this platform or any code 12 % in the platform should acknowledge the use of "PlatEMO" and reference "Ye 13 % Tian, Ran Cheng, Xingyi Zhang, and Yaochu Jin, PlatEMO: A MATLAB platform 14 % for evolutionary multi-objective optimization [educational forum], IEEE 15 % Computational Intelligence Magazine, 2017, 12(4): 73-87". 16 %-------------------------------------------------------------------------- 17 18 Distance = min(pdist2(PF,PopObj),[],2); 19 Score = mean(Distance); 20 end
NSGAIII.m
1 function NSGAIII(Global) 2 3 %% Generate the reference points and random population 4 [Z,Global.N] = UniformPoint(Global.N,Global.M); 5 Population = Global.Initialization(); 6 Zmin = min(Population(all(Population.cons<=0,2)).objs,[],1); 7 8 %% Optimization 9 while Global.NotTermination(Population) 10 MatingPool = TournamentSelection(2,Global.N,sum(max(0,Population.cons),2)); 11 Offspring = GA(Population(MatingPool)); 12 Zmin = min([Zmin;Offspring(all(Offspring.cons<=0,2)).objs],[],1); 13 Population = EnvironmentalSelection([Population,Offspring],Global.N,Z,Zmin); 14 end 15 end