1 内容介绍
蝠鲼觅食优化器 (MRFO) 已显示出处理单目标现实世界问题的良好能力，这使其在解决多目标问题中的应用成为一个有趣的方向。因此，本文研究了 MRFO 优化器，以开发一种新的算法来处理多目标工程设计问题。为了实现这一目标，采用精英概念，通过将外部存档集成到标准 MRFO 中来保存 Pareto 解决方案集。该档案也被认为是一个存储库，根据其密度程度从中选择搜索代理以控制蝠鲼种群的收敛性和多样性。我们的算法的效率首先通过对十个测试函数的广泛实验进行验证，在几乎所有情况下，结果在收敛性和多样性方面都非常令人满意。然后，将其应用于四个多目标工程问题，并在解决具有多目标的现实世界问题方面显示出良好的前景。
2 仿真代码
function G=CreateHypercubes(costs,ngrid,alpha)
empty_grid.Lower=[]; empty_grid.Upper=[]; G=repmat(empty_grid,nobj,1); for j=1:nobj min_cj=min(costs(j,:)); max_cj=max(costs(j,:)); dcj=alpha*(max_cj-min_cj); min_cj=min_cj-dcj; max_cj=max_cj+dcj; gx=linspace(min_cj,max_cj,ngrid-1); G(j).Lower=[-inf gx]; G(j).Upper=[gx inf]; endnobj=size(costs,1);
end
% MOMRFO: Multi-Objective Manta Ray Foraging Optimizer for% % Handling engineering design problems% % % % Main paper: % % % % A. Got, D. Zouache, A. Moussaoui, MOMRFO: Multi-Objective % % Manta Ray Foraging Optimizer for handling engineering design % % problems, Knowledge Based-Systems, in press,% % in press, DOI: /10.1016/j.knosys..107880 % %%___________________________________________________________________%
global pop_size ub lb M D R %% Select the test functions fun='ZDT1'; mo_Xrange; nVar=D; %% Parameter settings pop_size=100; % Population size MaxIt=2000;% Maximum Number of Iterations Ar_size=100; % The maximum size of archive alpha=0.1;% Grid Inflation Parameter nGrid=30;% Number of Grids per each Dimension beta=4; % Leader Selection Pressure Parameter gamma=2;% Extra (to be deleted) archive Member Selection Pressure %% The number of independent runs Nbrun=1; %% START THE EXECUTION OF THE ALGORITHM for ru=1:Nbrun %% Initialize the population of Manta Rays MantaRay=CreateEmptyMantaRay(pop_size); for i=1:pop_size for j=1:nVar MantaRay(i).Position(1,j)=unifrnd(lb(j),ub(j),1); end MantaRay(i).Cost=mo_test_function(MantaRay(i).Position,fun); %% fitness function MantaRay(i).Best.Position=MantaRay(i).Position; MantaRay(i).Best.Cost=MantaRay(i).Cost; end %% Prepare the archive for the first iteration MantaRay=DetermineDomination(MantaRay); Archive=GetNonDominatedMantaRay(MantaRay); Archive_costs=GetCosts(Archive); G=CreateHypercubes(Archive_costs,nGrid,alpha); for i=1:numel(Archive) [Archive(i).GridIndex Archive(i).GridSubIndex]=GetGridIndex(Archive(i),G); end costs=GetCosts(MantaRay); Archive_costs=GetCosts(Archive); %% The main loop of MOMRFO algorithm for t=1:MaxIt Coef=t/MaxIt; Leader=SelectLeader(Archive,beta); % Select Gbesr if rand<0.5 r1=rand; Beta2=2*exp(r1*((MaxIt-t+1)/MaxIt))*(sin(2*pi*r1)); if Coef>rand MantaRay(1).Position=Leader.Position+rand(1,D).*(Leader.Position-MantaRay(1).Position)+Beta2*(Leader.Position-MantaRay(1).Position); %Equation (4) else IndivRand=rand(1,D).*(ub-lb)+lb; MantaRay(1).Position=IndivRand+rand(1,D).*(IndivRand-MantaRay(1).Position)+Beta2*(IndivRand-MantaRay(1).Position); %Equation (7) end else Alpha2=2*rand(1,D).*(-log(rand(1,D))).^0.5; MantaRay(1).Position=MantaRay(1).Position+rand(1,D).*(Leader.Position-MantaRay(1).Position)+Alpha2.*(Leader.Position-MantaRay(1).Position); %Equation (1) end for i=2:pop_size Leader=SelectLeader(Archive,beta); %% Select Gbest if rand<0.5 r1=rand; Beta2=2*exp(r1*((MaxIt-t+1)/MaxIt))*(sin(2*pi*r1)); if Coef>rand MantaRay(i).Position=Leader.Position+rand(1,D).*(MantaRay(i-1).Position-MantaRay(i-1).Position)+Beta2*(Leader.Position-MantaRay(i).Position); %Equation (4) else IndivRand=rand(1,D).*(ub-lb)+lb; MantaRay(i).Position=IndivRand+rand(1,D).*(MantaRay(i-1).Position-MantaRay(i).Position)+Beta2*(IndivRand-MantaRay(i).Position); %Equation (7) end else Alpha2=2*rand(1,D).*(-log(rand(1,D))).^0.5; MantaRay(i).Position=MantaRay(i).Position+rand(1,D).*(MantaRay(i-1).Position-MantaRay(i).Position)+Alpha2.*(Leader.Position-MantaRay(i).Position); %Equation (1) end end %% Adjust boundaries if necessary for i=1:pop_size MantaRay(i).Position=min(max(MantaRay(i).Position,lb),ub); end S=2; for i=1:pop_size Leader=SelectLeader(Archive,beta); %% Select Gbest MantaRay(i).Position=MantaRay(i).Position+S*(rand*Leader.Position-rand*MantaRay(i).Position); %Equation (8) MantaRay(i).Position=min(max(MantaRay(i).Position,lb),ub); MantaRay(i).Cost=mo_test_function(MantaRay(i).Position,fun); end %% ubdate the archive MantaRay=DetermineDomination(MantaRay); Archive=[Archive;MantaRay]; Archive=DetermineDomination(Archive); Archive=GetNonDominatedMantaRay(Archive); for i=1:numel(Archive) [Archive(i).GridIndex Archive(i).GridSubIndex]=GetGridIndex(Archive(i),G); end %% Remove EXTRA-solutions from the high crowded hybercubes in the archive if numel(Archive)>Ar_size EXTRA=numel(Archive)-Ar_size; Archive=DeleteFromRep(Archive,EXTRA,gamma); end Archive_costs=GetCosts(Archive); G=CreateHypercubes(Archive_costs,nGrid,alpha); costs=GetCosts(MantaRay); Archive_costs=GetCosts(Archive); %% Front represent the final Pareto front Front=Archive_costs'; disp(['Iteration ' num2str(t) ' Probleme test ' fun ' Run ' num2str(ru) ': Number of Pareto Solutions = ' num2str(numel(Front(:,1)))]); pause(0:0.1); figure(1); plot(Front(:,1),Front(:,2),'*') end end filename=strcat('MOMRFO_',fun); save(filename); disp('END OF EXECUTION, THANKS!');format short;
global M global Dfunction obj=mo_test_function(x,fun)
% SCH FON POL KUR % DT1 ZDT2 ZDT3 ZDT4 ZDT6 % DTLZ1 DTLZ2 DTLZ3 DTLZ4 DTLZ5 DTLZ6 % F1 F2 F3 F4 F5 F6 F7 F8 F9 D=length(x);% MOP test functions
if strcmp(fun,'SCH') obj(1) = x(1)^2; obj(2) = (x(1)-2)^2; end % FON if strcmp(fun,'FON') D=3; y1=(x-1/sqrt(3)).^2; g1=sum(y1); y2=(x+1/sqrt(3)).^2; g2=sum(y2); obj(1) = 1-exp(-g1); obj(2) = 1-exp(-g2); end% SCH
if strcmp(fun,'POL') A1=0.5*sin(1)-2*cos(1)+sin(2)-1.5*cos(2); A2=1.5*sin(1)-cos(1)+2*sin(2)-0.5*cos(2); B1=0.5*sin(x(1))-2*cos(x(1))+sin(x(2))-1.5*cos(x(2)); B2=1.5*sin(x(1))-cos(x(1))+2*sin(x(2))-0.5*cos(x(2)); obj(1) = 1+(A1-B1)^2+(A2-B2)^2; obj(2) = (x(1)+3)^2+(x(2)+1)^2; end% POL
if strcmp(fun,'KUR') x1=x(1:end-1); x2=x(2:end); y1=-10*exp(-0.2*sqrt(x1.^2+x2.^2)); obj(1) =sum(y1); y2=abs(x).^0.8+5*sin(x.^3); obj(2) =sum(y2); end% KUR
if strcmp(fun,'ZDT1') y=x(2:end); gx=1+9*sum(y)/(D-1); obj(1)=x(1); obj(2)=gx*(1-sqrt(x(1)/gx)); end% ZDT1
if strcmp(fun,'ZDT2') y=x(2:end); gx=1+9*sum(y)/(D-1); obj(1)=x(1); obj(2)=gx*(1-(x(1)/gx)^2); end% ZDT2
if strcmp(fun,'ZDT3') y=x(2:end); gx=1+9*sum(y)/(D-1); obj(1)=x(1); obj(2)=gx*(1-sqrt(x(1)/gx)-x(1)/gx*sin(10*pi*x(1))); end% ZDT3
if strcmp(fun,'ZDT4') y=x(2:end); gx=1+10*(D-1)+sum(y.*y-10*cos(4*pi*y)); obj(1) = x(1); obj(2) = gx*(1-sqrt(x(1)/gx)); end% ZDT4
if strcmp(fun,'ZDT6') y=x(2:end); gx=1+9*(sum(y)/(D-1))^0.25; obj(1)=1-exp(-4*x(1))*(sin(6*pi*x(1)))^6; obj(2)=gx*(1-(obj(1)/gx)^2); end %DTLZ1 if strcmp(fun,'DTLZ1') switch M case 2 xg=x(2:end); gx=100*(D-M+1+sum((xg-0.5).^2-cos(20*pi*(xg-0.5)))); obj(1)=0.5*x(1)*(1+gx); obj(2)=0.5*(1-x(1))*(1+gx); case 3 xg=x(3:end); gx=100*(D-M+1+sum((xg-0.5).^2-cos(20*pi*(xg-0.5)))); obj(1)=0.5*x(1)*x(2)*(1+gx); obj(2)=0.5*x(1)*(1-x(2))*(1+gx); obj(3)=0.5*(1-x(1))*(1+gx); end end% ZDT6
if strcmp(fun,'DTLZ2') switch M case 2 xg=x(2:end); gx=sum((xg-0.5).^2); obj(1)=(1+gx)*cos(x(1)*0.5*pi); obj(2)=(1+gx)*sin(x(1)*0.5*pi); case 3 xg=x(3:end); gx=sum((xg-0.5).^2); obj(1)=(1+gx)*cos(x(1)*0.5*pi)*cos(x(2)*0.5*pi); obj(2)=(1+gx)*cos(x(1)*0.5*pi)*sin(x(2)*0.5*pi); obj(3)=(1+gx)*sin(x(1)*0.5*pi); end end%DTLZ2
if strcmp(fun,'DTLZ3') switch M case 2 xg=x(2:end); gx=100*(D-M+1+sum((xg-0.5).^2-cos(20*pi*(xg-0.5)))); obj(1)=(1+gx)*cos(x(1)*0.5*pi); obj(2)=(1+gx)*sin(x(1)*0.5*pi); case 3 xg=x(3:end); gx=100*(D-M+1+sum((xg-0.5).^2-cos(20*pi*(xg-0.5)))); obj(1)=(1+gx)*cos(x(1)*0.5*pi)*cos(x(2)*0.5*pi); obj(2)=(1+gx)*cos(x(1)*0.5*pi)*sin(x(2)*0.5*pi); obj(3)=(1+gx)*sin(x(1)*0.5*pi); end end%DTLZ3
if strcmp(fun,'DTLZ4') switch M case 2 xg=x(2:end); gx=sum((xg-0.5).^2); obj(1)=(1+gx)*cos((x(1)^100)*0.5*pi); obj(2)=(1+gx)*sin((x(1)^100)*0.5*pi); case 3 xg=x(3:end); gx=sum((xg-0.5).^2); obj(1)=(1+gx)*cos((x(1)^100)*0.5*pi)*cos((x(2)^100)*0.5*pi); obj(2)=(1+gx)*cos((x(1)^100)*0.5*pi)*sin((x(2)^100)*0.5*pi); obj(3)=(1+gx)*sin((x(1)^100)*0.5*pi); end end %DTLZ 5 if strcmp(fun,'DTLZ5') xg=x(3:end); gx=sum((xg-0.5).^2); Q1 = x(1); % Q2=(pi/(4*(1+gx)))*(1+2*gx*x(2)); Q2= 0.5*(1+2*gx*x(2))/(1+gx); obj(1)=(1+gx)*cos(Q1*0.5*pi).*cos(Q2*0.5*pi); obj(2)=(1+gx)*cos(Q1*0.5*pi).*sin(Q2*0.5*pi); obj(3)=(1+gx)*sin(Q1*0.5*pi); end%DTLZ4
if strcmp(fun,'DTLZ6') xg=x(3:end); gx=sum(xg.^0.1); Q1 = x(1); % Q2=(pi/(4*(1+gx)))*(1+2*gx*x(2)); Q2= 0.5*(1+2*gx*x(2))/(1+gx); obj(1)=(1+gx)*cos(Q1*0.5*pi).*cos(Q2*0.5*pi); obj(2)=(1+gx)*cos(Q1*0.5*pi).*sin(Q2*0.5*pi); obj(3)=(1+gx)*sin(Q1*0.5*pi); end% DTLZ6
if strcmp(fun,'DTLZ7') xg=x(3:end); gx=1+(9/(D-M+1))*sum(xg); obj(1)=x(1); obj(2)=x(2); hf=3-(obj(1)/(1+gx))*(1+sin(3*pi*obj(1)))-(obj(2)/(1+gx))*(1+sin(3*pi*obj(2))); obj(3)=(1+gx)*hf; end% DTLZ7
if strcmp(fun,'F1') for j = 1:D y(j) = x(j)-x(1)^(0.5*(1+3*(j-2)/(D-2))); end obj(1)=x(1); obj(2)=1-sqrt(x(1)); for k= 1:14 obj(1) = obj(1)+2/14*y(2*k+1)^2; end for k= 1:15 obj(2) = obj(2)+2/15*y(2*k)^2; end end % F2 if strcmp(fun,'F2') for j = 1:D y(j) = x(j)-sin(6*pi*x(1)+j*pi/D); end obj(1)=x(1); obj(2)=1-sqrt(x(1)); for k= 1:14 obj(1) = obj(1)+2/14*y(2*k+1)^2; end for k= 1:15 obj(2) = obj(2)+2/15*y(2*k)^2; end end% F1
if strcmp(fun,'F3') for j = 1:D if mod(j,2)==1 y(j) = x(j)-0.8*x(1)*cos(6*pi*x(1)+j*pi/D); else y(j) = x(j)-0.8*x(1)*sin(6*pi*x(1)+j*pi/D); end end obj(1)=x(1); obj(2)=1-sqrt(x(1)); for k= 1:14 obj(1) = obj(1)+2/14*y(2*k+1)^2; end for k= 1:15 obj(2) = obj(2)+2/15*y(2*k)^2; end end% F3
if strcmp(fun,'F4') for j = 1:D if mod(j,2)==1 y(j) = x(j)-0.8*x(1)*cos((6*pi*x(1)+j*pi/D)/3); else y(j) = x(j)-0.8*x(1)*sin(6*pi*x(1)+j*pi/D); end end obj(1)=x(1); obj(2)=1-sqrt(x(1)); for k= 1:14 obj(1) = obj(1)+2/14*y(2*k+1)^2; end for k= 1:15 obj(2) = obj(2)+2/15*y(2*k)^2; end end% F4
if strcmp(fun,'F5') for j = 1:D if mod(j,2)==1 y(j) = x(j)-(0.3*x(1)^2*cos(24*pi*x(1)+4*j*pi/D)+0.6*x(1))*cos(6*pi*x(1)+j*pi/D); else y(j) = x(j)-(0.3*x(1)^2*cos(24*pi*x(1)+4*j*pi/D)+0.6*x(1))*sin(6*pi*x(1)+j*pi/D); end end obj(1)=x(1); obj(2)=1-sqrt(x(1)); for k= 1:14 obj(1) = obj(1)+2/14*y(2*k+1)^2; end for k= 1:15 obj(2) = obj(2)+2/15*y(2*k)^2; end end% F5
if strcmp(fun,'F6') for j = 1:D y(j) = x(j)-2*x(2)*sin(2*pi*x(1)+j*pi/D); end obj(1)=cos(0.5*x(1)*pi)*cos(0.5*x(2)*pi); obj(2)=cos(0.5*x(1)*pi)*sin(0.5*x(2)*pi); obj(3)=sin(0.5*x(1)*pi); for k= 1:3 obj(1) = obj(1)+2/3*y(3*k+1)^2; end for k= 1:2 obj(2) = obj(2)+2/2*y(3*k+2)^2; end for k= 1:3 obj(3) = obj(3)+2/3*y(3*k)^2; end end% F6
if strcmp(fun,'F7') for j = 1:D y(j) = x(j)-x(1)^(0.5*(1+3*(j-2)/(D-2))); end for j = 1:D z(j) = 4*y(j)^2-cos(8*y(j)*pi)+1.0 ; end obj(1)=x(1); obj(2)=1-sqrt(x(1)); for k= 1:4 obj(1) = obj(1)+2/4*z(2*k+1); end for k= 1:5 obj(2) = obj(2)+2/5*z(2*k); end end% F7
if strcmp(fun,'F8') for j = 1:D y(j) = x(j)-x(1)^(0.5*(1+3*(j-2)/(D-2))); end for j = 1:D z(j) = cos((20*y(j)*pi)/sqrt(j)) ; end a1 = y(3)^2+y(5)^2+y(7)^2+y(9)^2; a2 = y(2)^2+y(4)^2+y(6)^2+y(8)^2+y(10)^2; b1 = z(3)*z(5)*z(7)*z(9); b2 = z(2)*z(4)*z(6)*z(8)*z(10); obj(1) = x(1)+2/4*(4*a1-2*b1+2); obj(2) = 1-sqrt(x(1))+2/5*(4*a2-2*b2+2); end% F8
if strcmp(fun,'F9') for j = 1:D y(j) = x(j)-sin(6*pi*x(1)+j*pi/D); end obj(1)=x(1); obj(2)=1-x(1)^2; for k= 1:14 obj(1) = obj(1)+2/14*y(2*k+1)^2; end for k= 1:15 obj(2) = obj(2)+2/15*y(2*k)^2; end end if strcmp(fun,'MaF1') g = sum((x(3:end)-0.5).^2,2); obj = repmat(1+g,1,M) - repmat(1+g,1,M).*fliplr(cumprod([ones(size(g,1),1),x(1:M-1)],2)).*[ones(size(g,1),1),1-x(M-1:-1:1)]; end% F9
g= 100*(D-M+1+sum((x(M:end)-0.5).^2-cos(20.*pi.*(x(M:end)-0.5)),2)); obj = repmat(1+g,1,M).*fliplr(cumprod([ones(size(g,1),1),cos(x(1:M-1)*pi/2)],2)).*[ones(size(g,1),1),sin(x(M-1:-1:1)*pi/2)]; obj = [x(1:M-1).^4,obj(:,M).^2]; endif strcmp(fun,'MaF3')
J1 = 3 : 2 : D; J2 = 2 : 2 : D; Y = x - sin(6*pi*repmat(x(:,1),1,D)+repmat(1:D,size(x,1),1)*pi/D); obj(:,1) = x(:,1)+ 2*mean(Y(:,J1).^2,2); obj(:,2) = 1-sqrt(x(:,1)) + 2*mean(Y(:,J2).^2,2); endif strcmp(fun,'UF1')
if strcmp(fun,'UF2') J1 = 3 : 2 : D; J2 = 2 : 2 : D; Y = zeros(size(x)); X1 = repmat(x(:,1),1,length(J1)); Y(:,J1) = x(:,J1)-(0.3*X1.^2.*cos(24*pi*X1+4*repmat(J1,size(x,1),1)*pi/D)+0.6*X1).*cos(6*pi*X1+repmat(J1,size(x,1),1)*pi/D); X1 = repmat(x(:,1),1,length(J2)); Y(:,J2) = x(:,J2)-(0.3*X1.^2.*cos(24*pi*X1+4*repmat(J2,size(x,1),1)*pi/D)+0.6*X1).*sin(6*pi*X1+repmat(J2,size(x,1),1)*pi/D); obj(:,1) = x(:,1)+ 2*mean(Y(:,J1).^2,2); obj(:,2) = 1-sqrt(x(:,1)) + 2*mean(Y(:,J2).^2,2); end% UF2
if strcmp(fun,'UF3') J1 = 3 : 2 : D; J2 = 2 : 2 : D; Y = x - repmat(x(:,1),1,D).^(0.5*(1+3*(repmat(1:D,size(x,1),1)-2)/(D-2))); obj(1) = x(:,1)+ 2/length(J1)*(4*sum(Y(:,J1).^2,2)-2*prod(cos(20*Y(:,J1)*pi./sqrt(repmat(J1,size(x,1),1))),2)+2); obj(2) = 1-sqrt(x(:,1)) + 2/length(J2)*(4*sum(Y(:,J2).^2,2)-2*prod(cos(20*Y(:,J2)*pi./sqrt(repmat(J2,size(x,1),1))),2)+2); end% UF3
if strcmp(fun,'UF4') J1 = 3 : 2 : D; J2 = 2 : 2 : D; Y =x - sin(6*pi*repmat(x(:,1),1,D)+repmat(1:D,size(x,1),1)*pi/D); hY = abs(Y)./(1+exp(2*abs(Y))); obj(1) = x(:,1) + 2*mean(hY(:,J1),2); obj(2) = 1-x(:,1).^2 + 2*mean(hY(:,J2),2); end% UF4
if strcmp(fun,'UF5') J1 = 3 : 2 : D; J2 = 2 : 2 : D; Y = x - sin(6*pi*repmat(x(:,1),1,D)+repmat(1:D,size(x,1),1)*pi/D); hY = 2*Y.^2 - cos(4*pi*Y) + 1; obj(1) = x(:,1) + (1/20+0.1)*abs(sin(20*pi*x(:,1)))+2*mean(hY(:,J1),2); obj(2) = 1-x(:,1) + (1/20+0.1)*abs(sin(20*pi*x(:,1)))+2*mean(hY(:,J2),2); end% UF5
if strcmp(fun,'UF6') J1 = 3 : 2 : D; J2 = 2 : 2 : D; Y =x - sin(6*pi*repmat(x(:,1),1,D)+repmat(1:D,size(x,1),1)*pi/D); obj(1) = x(:,1) + max(0,2*(1/4+0.1)*sin(4*pi*x(:,1)))+2/length(J1)*(4*sum(Y(:,J1).^2,2)-2*prod(cos(20*Y(:,J1)*pi./sqrt(repmat(J1,size(x,1),1))),2)+2); obj(2) = 1-x(:,1) + max(0,2*(1/4+0.1)*sin(4*pi*x(:,1)))+2/length(J2)*(4*sum(Y(:,J2).^2,2)-2*prod(cos(20*Y(:,J2)*pi./sqrt(repmat(J2,size(x,1),1))),2)+2); end% UF6
if strcmp(fun,'UF7') J1 = 3 : 2 : D; J2 = 2 : 2 : D; Y = x - sin(6*pi*repmat(x(:,1),1,D)+repmat(1:D,size(x,1),1)*pi/D); obj(1) = x(:,1).^0.2 + 2*mean(Y(:,J1).^2,2); obj(2) = 1-x(:,1).^0.2 + 2*mean(Y(:,J2).^2,2); end% UF7
if strcmp(fun,'UF8') J1 = 4 : 3 : D; J2 = 5 : 3 : D; J3 = 3 : 3 : D; Y =x - 2*repmat(x(:,2),1,D).*sin(2*pi*repmat(x(:,1),1,D)+repmat(1:D,size(x,1),1)*pi/D); obj(1) = cos(0.5*x(:,1)*pi).*cos(0.5*x(:,2)*pi) + 2*mean(Y(:,J1).^2,2); obj(2) = cos(0.5*x(:,1)*pi).*sin(0.5*x(:,2)*pi) + 2*mean(Y(:,J2).^2,2); obj(3) = sin(0.5*x(:,1)*pi)+ 2*mean(Y(:,J3).^2,2); end% UF8
if strcmp(fun,'UF9') J1 = 4 : 3 : D; J2 = 5 : 3 : D; J3 = 3 : 3 : D; Y = x - 2*repmat(x(:,2),1,D).*sin(2*pi*repmat(x(:,1),1,D)+repmat(1:D,size(x,1),1)*pi/D); obj(1) = 0.5*(max(0,1.1*(1-4*(2*x(:,1)-1).^2))+2*x(:,1)).*x(:,2) + 2*mean(Y(:,J1).^2,2); obj(2) = 0.5*(max(0,1.1*(1-4*(2*x(:,1)-1).^2))-2*x(:,1)+2).*x(:,2) + 2*mean(Y(:,J2).^2,2); obj(3) = 1-x(:,2)+ 2*mean(Y(:,J3).^2,2); end% UF9
if strcmp(fun,'UF10') J1 = 4 : 3 : D; J2 = 5 : 3 : D; J3 = 3 : 3 : D; Y = x - 2*repmat(x(:,2),1,D).*sin(2*pi*repmat(x(:,1),1,D)+repmat(1:D,size(x,1),1)*pi/D); Y = 4*Y.^2 - cos(8*pi*Y) + 1; obj(1) = cos(0.5*x(:,1)*pi).*cos(0.5*x(:,2)*pi) + 2*mean(Y(:,J1),2); obj(2) = cos(0.5*x(:,1)*pi).*sin(0.5*x(:,2)*pi) + 2*mean(Y(:,J2),2); obj(3) = sin(0.5*x(:,1)*pi)+ 2*mean(Y(:,J3),2); end% UF10
if strcmp(fun,'Weld') obj(1)=1.10471*x(1)^2*x(2)+0.04811*x(3)*x(4)*(14.0+x(2)); obj(2)=65856000/(30*10^6*x(4)*x(3)^3); Z=0; % Penalty constant lam=10^15; Q=6000*(14+x(2)/2); d=sqrt(x(2)^2/4+(x(1)+x(3))^2/4); J=2*(x(1)*x(2)*sqrt(2)*(x(2)^2/12+(x(1)+x(3))^2/4)); alpha=6000/(sqrt(2)*x(1)*x(2)); beta=Q*d/J; tau=sqrt(alpha^2+2*alpha*beta*x(2)/(2*d)+beta^2); sigma=504000/(x(4)*x(3)^2); tmpf=4.013*(30*10^6)/196; P=tmpf*sqrt(x(3)^2*x(4)^6/36)*(1-x(3)*sqrt(30/48)/28); g(1)=tau-13600; g(2)=sigma-30000; g(3)=x(1)-x(4); g(4)=6000-P; % No equality constraint in this problem, so empty; geq=[];% Welded beam design problem
for k=1:length(g), Z=Z+ lam*g(k)^2*getH(g(k)); end % Apply equality constraints for k=1:length(geq), Z=Z+lam*geq(k)^2*getHeq(geq(k)); end obj=obj+Z; end% Apply inequality constraints
if strcmp(fun,'Speed') obj(1)=0.7854*x(1)*x(2)^2*(3.3333*x(3)^2+14.9334*x(3)-43.0934)-1.508*x(1)*(x(6)^2+x(7)^2)+7.4777*(x(6)^3+x(7)^3)+0.7854*(x(4)*x(6)^2+x(5)*x(7)^2); obj(2) =sqrt(((745*x(4))/(x(2)*x(3)))^2+1.69e7)/(0.1*x(6)^3); Z=0; % Penalty constant lam=10^15; g(1)=27/(x(1)*x(2)^2*x(3))-1; g(2)=397.5/(x(1)*x(2)^2*x(3)^2)-1; g(3)=(1.93*x(4)^3)/(x(2)*x(3)*x(6)^4)-1; g(4)=(1.93*x(5)^3)/(x(2)*x(3)*x(7)^4)-1; g(5)=((sqrt(((745*x(4))/(x(2)*x(3)))^2+16.9e6))/(110*x(6)^3))-1; g(6)=((sqrt(((745*x(5))/(x(2)*x(3)))^2+157.5e6))/(85*x(7)^3))-1; g(7)=((x(2)*x(3))/40)-1; g(8)=(5*x(2)/x(1))-1; g(9)=(x(1)/12*x(2))-1; g(10)=((1.5*x(6)+1.9)/x(4))-1; g(11)=((1.1*x(7)+1.9)/x(5))-1;% Speed reducer design problem
geq=[];% No equality constraint in this problem, so empty;
for k=1:length(g), Z=Z+ lam*g(k)^2*getH(g(k)); end % Apply equality constraints for k=1:length(geq), Z=Z+lam*geq(k)^2*getHeq(geq(k)); end obj=obj+Z; end% Apply inequality constraints
if strcmp(fun,'Bar') obj(1)=200*(2*x(1)+sqrt(2*x(2))+sqrt(x(3))+x(4)); obj(2)=0.01*((2/x(2))+2*(sqrt(2)/x(2))-2*(sqrt(2)/x(3))+(2/x(4))); end% Four-bar truss problem
if strcmp(fun,'Disk') obj(1)=4.9*(10^(-5))*(x(2)^2-x(1)^2)*(x(4)-1); obj(2)=9.82*10^6*(x(2)^2-x(1)^2)/((x(2)^3-x(1)^3)*x(4)*x(3)); Z=0; % Penalty constant lam=10^15; g(1) = 20+x(1)-x(2); g(2) = 2.5*(x(4)+1)-30; g(3) = (x(3))/(3.14*(x(2)^2-x(1)^2)^2)-0.4; g(4) = (2.22*10^(-3)*x(3)*(x(2)^3-x(1)^3))/((x(2)^2-x(1)^2)^2)-1; g(5) = 900-(2.66*10^(-2)*x(3)*x(4)*(x(2)^3-x(1)^3))/((x(2)^2-x(1)^2));% Disk brake design problem
geq=[];% No equality constraint in this problem, so empty;
for k=1:length(g), Z=Z+ lam*g(k)^2*getH(g(k)); end % Apply equality constraints for k=1:length(geq), Z=Z+lam*geq(k)^2*getHeq(geq(k)); end obj=obj+Z; end% Apply inequality constraints
end
if g<=0, H=0; else H=1; end end % Test if equalities hold function H=geteqH(g) if g==0, H=0; else H=1; end endfunction H=getH(g)
function Archive=DeleteFromRep(Archive,EXTRA,gamma)
if nargin<3
gamma=1;
end
for k=1:EXTRA
[occ_cell_index occ_cell_member_count]=GetOccupiedCells(Archive);
p=occ_cell_member_count.^gamma;
p=p/sum(p);
selected_cell_index=occ_cell_index(RouletteWheelSelection(p));
GridIndices=[Archive.GridIndex];
selected_cell_members=find(GridIndices==selected_cell_index);
n=numel(selected_cell_members);
selected_memebr_index=randi([1 n]);
j=selected_cell_members(selected_memebr_index);
Archive=[Archive(1:j-1); Archive(j+1:end)];
end
end
3 运行结果
4 参考文献
[1] Ag A , Dz A , Am B . MOMRFO: Multi-objective Manta ray foraging optimizer for handling engineering design problems. .
博主简介：擅长智能优化算法、神经网络预测、信号处理、元胞自动机、图像处理、路径规划、无人机等多种领域的Matlab仿真，相关matlab代码问题可私信交流。
部分理论引用网络文献，若有侵权联系博主删除。
