www.gusucode.com > 几种多目标优化算法集合，包括MOEAD，MOPSO，NNIA，NSGA2等 > NNIAToolbox/NNIA.m
function NNIA()
% NNIA.m
% NNIA performs experiment presented in our paper:
% Maoguo Gong, Licheng Jiao,Haifeng Du, Liefeng Bo. Multiobjective Immune Algorithm with Nondominated Neighbor-based Selection.
% Evolutionary Computation Journal, MIT Press, 2007, in press.
% avialable from http://see.xidian.edu.cn/iiip/mggong/publication.htm
%
% Authors: Maoguo Gong and Licheng Jiao
% April 7, 2006
% Copyright (C) 2005-2006 by Maoguo Gong (e-mail: gong@ieee.org)
%
% REFERENCE
% Please refer to the following paper if you use the toolbox NNIA
% Maoguo Gong, Licheng Jiao,Haifeng Du, Liefeng Bo. Multiobjective Immune Algorithm with Nondominated Neighbor-based Selection.
% Evolutionary Computation Journal, MIT Press, 2007, in press.
% NNIA is Copyrighted by the Authors.
%
% RELEASE
% version 1.0, 2006/04/10, tested on Matlab 7.0
% license granted for research use ONLY
% Copyright (C) 2005-2007 by Maoguo Gong (e-mail: gong@ieee.org)
% Please refer to the Readme.txt for a detailed description.
%--------------------------------------------------------------------------
clear all;
%--------------------------------------------------------------------------
disp(sprintf('Nondominated Neighbor Immune Algorithm (NNIA) min'));
disp(sprintf('Authour: Maoguo Gong and Licheng Jiao'));
disp(sprintf('Last Modified: Oct. 10, 2007'));
EMOinstruction;%% display the instruction for running the programming.
%--------------------------------------------------------------------------
TestNO=input('press the enter key after inputting the serial number of test problem:');
Trial=input('input the number of independent runs:');
NA=20; % size of active population
CS=100; % clonal scale
NM=100; % size of dominant population
gmax=500; % maximum number of iterations
[bu,bd,testfunction]=getbud(TestNO); % bu denotes the upper boundary of variable;bd denotes the nether boundary of variable
c=size(bu,2);
pm=1/c;
if bu==bd return;end
%--------------------------------------------------------------------------
Datime=date;
Datime(size(Datime,2)-4:size(Datime,2))=[];%%test date
TestTime=clock;%%test time
TestTime=[num2str(TestTime(4)),'-',num2str(TestTime(5))];
Method='NNIA';
paretof=[];
runtime=[];
for trial=1:Trial
timerbegin=clock;
%--------------------------------------------------------------------------
POP=rand(NM,c).*(ones(NM,1)*bu-ones(NM,1)*bd)+ones(NM,1)*bd;
ME=[];%Initialization
%--------------------------------------------------------------------------
pa=OVcom(POP,TestNO);
[DON,DONt]=IDAf(pa);
nodom=find(DON==1);
MEpa=pa(nodom,:);MEPOP=POP(nodom,:);
numnod=size(MEPOP,1);
[ClonePOP,Clonepa,RNDCDt]=UDPf(MEPOP,MEpa,NA);
it=0;Cloneti=0;DASti=0;PNmti=0;DONjudti=0;RNDCDti=0;AbAbfitcomti=0;
while it<gmax
%--------------------------------------------------------------------------
cloneover=[];
[cloneover,Clonet]=Clonef(ClonePOP,Clonepa,CS);
[cloneover,DASt]=Recombinationf(cloneover,bu,bd,ClonePOP);
[cloneover,PNmt]=Mutationf(cloneover,bu,bd,pm);
%--------------------------------------------------------------------------
clonepa=OVcom(cloneover,TestNO);
NPOP=[MEPOP;cloneover];
Npa=[MEpa;clonepa];
[NDON,DONjudt]=IDAf(Npa);
Nnodom=find(NDON==1);
NEpa=Npa(Nnodom,:);NEPOP=NPOP(Nnodom,:);
Nnumnod=size(NEPOP,1);
[MEPOP MEpa,RNDCDt]=UDPf(NEPOP,NEpa,NM);
numnod=size(MEPOP,1);
[ClonePOP,Clonepa,RNDCDt]=UDPf(MEPOP,MEpa,NA);
%Update Dominant Population
%--------------------------------------------------------------------------
it=it+1;
Cloneti=Cloneti+Clonet;DASti=DASti+DASt;PNmti=PNmti+PNmt;
DONjudti=DONjudti+DONjudt;RNDCDti=RNDCDti+RNDCDt;
disp(sprintf('time: %d generation: %d number of nodominate: %d',trial,it,numnod));
end %the end of iterations
%--------------------------------------------------------------------------
%Save the output solutions
[NS(trial),NF]=size(MEpa);Trials=trial;
paretof=[paretof;MEpa];
runtime(trial)=etime(clock,timerbegin);
Clonetime(trial)=Cloneti;DAStime(trial)=DASti;PNmtime(trial)=PNmti;
DONjudtime(trial)=DONjudti;RNDCDtime(trial)=RNDCDti;
eval(['save ', testfunction Method Datime TestTime ,' Method gmax NA CS paretof runtime Trials NS NF TestTime Datime testfunction ']) ;
end %the end of runs
%--------------------------------------------------------------------------
Frontshow(MEpa);% plot the Pareto fronts solved by the last run
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [DA,time]=IDAf(pa)
%--------------------------------------------------------------------------
tic;
[N,C]=size(pa);
DA=ones(N,1);
for i=1:N
temppa=pa;
temppa(i,:)=[];
LEsign=ones(N-1,1);
for j=1:C
LessEqual=find(temppa(:,j)<=pa(i,j));
tepa=[];tepa=temppa(LessEqual,:);
temppa=[];temppa=tepa;
end
if size(temppa,1)~=0
k=1;
while k<=C
Lessthan=[];
Lessthan=find(temppa(:,k)<pa(i,k));
if size(Lessthan,1)~=0
DA(i)=0;k=C+1;
else
k=k+1;
end
end
end
end
time=toc;
%%-------------------------------------------------------------------------
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [EPOP,Epa,time]=UDPf(POP,pa,NM);
%--------------------------------------------------------------------------
tic;
[Ns,C]=size(pa);padis=[];
i=1;
while i<Ns
deltf=pa-ones(Ns,1)*pa(i,:);
deltf(i,:)=inf;
aa=find(sum(abs(deltf),2)==0);
POP(aa,:)=[];pa(aa,:)=[];
[Ns,C]=size(pa);i=i+1;
end
if Ns>NM
for i=1:C
[k,l]=sort(pa(:,i));
N=[];M=[];N=pa(l,:);M=POP(l,:);
pa=[];POP=[];pa=N;POP=M;
pa(1,C+1)=Inf;pa(Ns,C+1)=Inf;
pai1=[];pad1=[];pai1=pa(3:Ns,i);pad1=pa(1:(Ns-2),i);
fimin=min(pa(:,i));fimax=max(pa(:,i));
pa(2:(Ns-1),C+1)=pa(2:(Ns-1),C+1)+(pai1-pad1)/(0.0001+fimax-fimin);
end
padis=pa(:,C+1);pa=pa(:,1:C);POP=POP;
[aa,ss]=sort(-padis);
EPOP=POP(ss(1:NM),:);Epa=pa(ss(1:NM),:);
else
EPOP=POP;Epa=pa;
end
time=toc;
%%-------------------------------------------------------------------------
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [NPOP,time]=Clonef(POP,pa,CS)
%--------------------------------------------------------------------------
tic
NC=[];
[N,C]=size(POP);
[POP,pa,padis]=CDAf(POP,pa);
aa=find(padis==inf);
bb=find(padis~=inf);
if length(bb)>0
padis(aa)=2*max(max(padis(bb)));
NC=ceil(CS*padis./sum(padis));
else
NC=ceil(CS/length(aa))+zeros(1,N);
end
NPOP=[];
for i=1:N
NiPOP=ones(NC(i),1)*POP(i,:);
NPOP=[NPOP;NiPOP];
end
time=toc;
%--------------------------------------------------------------------------
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [POP,pa,padis]=CDAf(tPOP,tpa)
%--------------------------------------------------------------------------
[Ns,C]=size(tpa);padis=[];
for i=1:C
[k,l]=sort(tpa(:,i));
N=[];M=[];N=tpa(l,:);M=tPOP(l,:);
tpa=[];tPOP=[];tpa=N;tPOP=M;
tpa(1,C+1)=Inf;tpa(Ns,C+1)=Inf;
tpai1=[];tpad1=[];tpai1=tpa(3:Ns,i);tpad1=tpa(1:(Ns-2),i);
fimin=min(tpa(:,i));fimax=max(tpa(:,i));
tpa(2:(Ns-1),C+1)=tpa(2:(Ns-1),C+1)+(tpai1-tpad1)/(0.0001+fimax-fimin);
end
pa=tpa(:,1:C);POP=tPOP;padis=tpa(:,C+1);
%--------------------------------------------------------------------------
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [NPOP,time]=Recombinationf(POP,bu,bd,CPOP);
%--------------------------------------------------------------------------
tic;
eta_c=15;
[N,C]=size(POP);
NPOP=POP;
for i=1:N
r1=rand;
if r1<=1.1
aa=randperm(size(CPOP,1));bb=aa(1);
for j=1:C
par1=POP(i,j);par2=CPOP(bb,j);
yd=bd(j);yu=bu(j);
r2=rand;
if r2<=0.5
if abs(par1-par2)>10^(-14)
y1=min(par1,par2);y2=max(par1,par2);
if (y1-yd)>(yu-y2)
beta=1+2*(yu-y2)/(y2-y1);
else
beta=1+2*(y1-yd)/(y2-y1);
end
expp=eta_c+1;beta=1/beta;alpha=2.0-beta^(expp);
r3=rand;
if r3<=1/alpha
alpha=alpha*r3;expp=1/(eta_c+1.0);
betaq=alpha^(expp);
else
alpha=1/(2.0-alpha*r3);expp=1/(eta_c+1);
betaq=alpha^(expp);
end
chld1=0.5*((y1+y2)-betaq*(y2-y1));
chld2=0.5*((y1+y2)+betaq*(y2-y1));
if rand<=0.5
aa=max(chld1,yd);NPOP(i,j)=min(aa,yu);
else
aa=max(chld2,yd);NPOP(i,j)=min(aa,yu);
end
end
end
end
end
end
time=toc;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [NPOP,time]=Mutationf(POP,bu,bd,pm)
%--------------------------------------------------------------------------
tic;
[N,C]=size(POP);
eta_m=20;
NPOP=POP;
for i=1:N
for j=1:C
r1=rand;
if r1<=pm
y=POP(i,j);
yd=bd(j);yu=bu(j);
if y>yd
if (y-yd)<(yu-y)
delta=(y-yd)/(yu-yd);
else
delta=(yu-y)/(yu-yd);
end
r2=rand;
indi=1/(eta_m+1);
if r2<=0.5
xy=1-delta;
val=2*r2+(1-2*r2)*(xy^(eta_m+1));
deltaq=val^indi-1;
else
xy=1-delta;
val=2*(1-r2)+2*(r2-0.5)*(xy^(eta_m+1));
deltaq=1-val^indi;
end
y=y+deltaq*(yu-yd);
NPOP(i,j)=min(y,yu);NPOP(i,j)=max(y,yd);
else
NPOP(i,j)=rand*(yu-yd)+yd;
end
end
end
end
time=toc;
%--------------------------------------------------------------------------
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function Frontshow(Epa)
%--------------------------------------------------------------------------
if size(Epa,2)==2
%%plot pareto fronts in 2-D
f1=Epa(:,1);f2=Epa(:,2);
plot(f1,f2,'r*'); grid on;
xlabel('Function 1');
ylabel('Function 2');
hold on
elseif size(Epa,2)>=3
%%plot pareto fronts in 3-D
f1=Epa(:,1);f2=Epa(:,2);f3=Epa(:,3);
plot3(f1,f2,f3,'kd'); grid on;
xlabel('Function 1');
ylabel('Function 2');
zlabel('Function 3');
hold on
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%