1 内容介绍
蝠鲼觅食优化器 (MRFO) 已显示出处理单目标现实世界问题的良好能力,这使其在解决多目标问题中的应用成为一个有趣的方向。因此,本文研究了 MRFO 优化器,以开发一种新的算法来处理多目标工程设计问题。为了实现这一目标,采用精英概念,通过将外部存档集成到标准 MRFO 中来保存 Pareto 解决方案集。该档案也被认为是一个存储库,根据其密度程度从中选择搜索代理以控制蝠鲼种群的收敛性和多样性。我们的算法的效率首先通过对十个测试函数的广泛实验进行验证,在几乎所有情况下,结果在收敛性和多样性方面都非常令人满意。然后,将其应用于四个多目标工程问题,并在解决具有多目标的现实世界问题方面显示出良好的前景。
2 仿真代码
function G=CreateHypercubes(costs,ngrid,alpha)
nobj=size(costs,1);
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];
end
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: https://doi.org/10.1016/j.knosys.2021.107880 % %
format short;
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!');
function obj=mo_test_function(x,fun)
global M
global D
% MOP test functions
% 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);
% SCH
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
% POL
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
% KUR
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
% ZDT1
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
% ZDT2
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
% ZDT3
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
% ZDT4
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
% ZDT6
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
%DTLZ2
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
%DTLZ3
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
%DTLZ4
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
% DTLZ6
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
% DTLZ7
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
% F1
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
% F3
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
% F4
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
% F5
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
% F6
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
% F7
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
% F8
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
% F9
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
if strcmp(fun,'MaF3')
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];
end
if strcmp(fun,'UF1')
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);
end
% UF2
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
% UF3
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
% UF4
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
% UF5
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
% UF6
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
% UF7
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
% UF8
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
% UF9
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
% UF10
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
% Welded beam design problem
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=[];
% Apply inequality constraints
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
% Speed reducer design problem
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;
% No equality constraint in this problem, so empty;
geq=[];
% Apply inequality constraints
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
% Four-bar truss problem
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
% Disk brake design 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));
% No equality constraint in this problem, so empty;
geq=[];
% Apply inequality constraints
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
end
function H=getH(g)
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
end
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. 2021.
博主简介:擅长智能优化算法、神经网络预测、信号处理、元胞自动机、图像处理、路径规划、无人机等多种领域的Matlab仿真,相关matlab代码问题可私信交流。
部分理论引用网络文献,若有侵权联系博主删除。