% This function containts full information and implementations of the benchmark % functions in Table 1, Table 2, and Table 3 in the paper% lb is the lower bound: lb=[lb_1,lb_2,...,lb_d]% up is the uppper bound: ub=[ub_1,ub_2,...,ub_d]% dim is the number of variables (dimension of the problem)function[lb,ub,dim,fobj]=Get_Functions_details(F)switch Fcase'F1'fobj =@F1;lb=-100;ub=100;dim=30;case'F2'fobj =@F2;lb=-10;ub=10;dim=30;case'F3'fobj =@F3;lb=-100;ub=100;dim=30;case'F4'fobj =@F4;lb=-100;ub=100;dim=30;case'F5'fobj =@F5;lb=-30;ub=30;dim=30;case'F6'fobj =@F6;lb=-100;ub=100;dim=30;case'F7'fobj =@F7;lb=-1.28;ub=1.28;dim=30;case'F8'fobj =@F8;lb=-500;ub=500;dim=30;case'F9'fobj =@F9;lb=-5.12;ub=5.12;dim=30;case'F10'fobj =@F10;lb=-32;ub=32;dim=30;case'F11'fobj =@F11;lb=-600;ub=600;dim=30;case'F12'fobj =@F12;lb=-50;ub=50;dim=30;case'F13'fobj =@F13;lb=-50;ub=50;dim=30;case'F14'fobj =@F14;lb=-65.536;ub=65.536;dim=2;case'F15'fobj =@F15;lb=-5;ub=5;dim=4;case'F16'fobj =@F16;lb=-5;ub=5;dim=2;case'F17'fobj =@F17;lb=[-5,0];ub=[10,15];dim=2;case'F18'fobj =@F18;lb=-2;ub=2;dim=2;case'F19'fobj =@F19;lb=0;ub=1;dim=3;case'F20'fobj =@F20;lb=0;ub=1;dim=6;case'F21'fobj =@F21;lb=0;ub=10;dim=4;case'F22'fobj =@F22;lb=0;ub=10;dim=4;case'F23'fobj =@F23;lb=0;ub=10;dim=4;endend% F1function o =F1(x)
o=sum(x.^2);end% F2function o =F2(x)
o=sum(abs(x))+prod(abs(x));end% F3function o =F3(x)
dim=size(x,2);
o=0;fori=1:dimo=o+sum(x(1:i))^2;endend% F4function o =F4(x)
o=max(abs(x));end% F5function o =F5(x)
dim=size(x,2);
o=sum(100*(x(2:dim)-(x(1:dim-1).^2)).^2+(x(1:dim-1)-1).^2);end% F6function o =F6(x)
o=sum(abs((x+.5)).^2);end% F7function o =F7(x)
dim=size(x,2);
o=sum([1:dim].*(x.^4))+rand;end% F8function o =F8(x)
o=sum(-x.*sin(sqrt(abs(x))));end% F9function o =F9(x)
dim=size(x,2);
o=sum(x.^2-10*cos(2*pi.*x))+10*dim;end% F10function o =F10(x)
dim=size(x,2);
o=-20*exp(-.2*sqrt(sum(x.^2)/dim))-exp(sum(cos(2*pi.*x))/dim)+20+exp(1);end% F11function o =F11(x)
dim=size(x,2);
o=sum(x.^2)/4000-prod(cos(x./sqrt([1:dim])))+1;end% F12function o =F12(x)
dim=size(x,2);
o=(pi/dim)*(10*((sin(pi*(1+(x(1)+1)/4)))^2)+sum((((x(1:dim-1)+1)./4).^2).*...(1+10.*((sin(pi.*(1+(x(2:dim)+1)./4)))).^2))+((x(dim)+1)/4)^2)+sum(Ufun(x,10,100,4));end% F13function o =F13(x)
dim=size(x,2);
o=.1*((sin(3*pi*x(1)))^2+sum((x(1:dim-1)-1).^2.*(1+(sin(3.*pi.*x(2:dim))).^2))+...((x(dim)-1)^2)*(1+(sin(2*pi*x(dim)))^2))+sum(Ufun(x,5,100,4));end% F14function o =F14(x)
aS=[-32-1601632-32-1601632-32-1601632-32-1601632-32-1601632;,...-32-32-32-32-32-16-16-16-16-160000016161616163232323232];forj=1:25bS(j)=sum((x'-aS(:,j)).^6);end
o=(1/500+sum(1./([1:25]+bS))).^(-1);end% F15function o =F15(x)
aK=[.1957.1947.1735.16.0844.0627.0456.0342.0323.0235.0246];
bK=[.25.51246810121416];bK=1./bK;
o=sum((aK-((x(1).*(bK.^2+x(2).*bK))./(bK.^2+x(3).*bK+x(4)))).^2);end% F16function o =F16(x)
o=4*(x(1)^2)-2.1*(x(1)^4)+(x(1)^6)/3+x(1)*x(2)-4*(x(2)^2)+4*(x(2)^4);end% F17function o =F17(x)
o=(x(2)-(x(1)^2)*5.1/(4*(pi^2))+5/pi*x(1)-6)^2+10*(1-1/(8*pi))*cos(x(1))+10;end% F18function o =F18(x)
o=(1+(x(1)+x(2)+1)^2*(19-14*x(1)+3*(x(1)^2)-14*x(2)+6*x(1)*x(2)+3*x(2)^2))*...(30+(2*x(1)-3*x(2))^2*(18-32*x(1)+12*(x(1)^2)+48*x(2)-36*x(1)*x(2)+27*(x(2)^2)));end% F19function o =F19(x)
aH=[31030;.11035;31030;.11035];cH=[11.233.2];
pH=[.3689.117.2673;.4699.4387.747;.1091.8732.5547;.03815.5743.8828];
o=0;fori=1:4o=o-cH(i)*exp(-(sum(aH(i,:).*((x-pH(i,:)).^2))));endend% F20function o =F20(x)
aH=[103173.51.78;.051017.1814;33.51.710178;178.0510.114];
cH=[11.233.2];
pH=[.1312.1696.5569.0124.8283.5886;.2329.4135.8307.3736.1004.9991;....2348.1415.3522.2883.3047.6650;.4047.8828.8732.5743.1091.0381];
o=0;fori=1:4o=o-cH(i)*exp(-(sum(aH(i,:).*((x-pH(i,:)).^2))));endend% F21function o =F21(x)
aSH=[4444;1111;8888;6666;3737;2929;5533;8181;6262;73.673.6];
cSH=[.1.2.2.4.4.6.3.7.5.5];o=0;fori=1:5o=o-((x-aSH(i,:))*(x-aSH(i,:))'+cSH(i))^(-1);endend% F22function o =F22(x)
aSH=[4444;1111;8888;6666;3737;2929;5533;8181;6262;73.673.6];
cSH=[.1.2.2.4.4.6.3.7.5.5];o=0;fori=1:7o=o-((x-aSH(i,:))*(x-aSH(i,:))'+cSH(i))^(-1);endend% F23function o =F23(x)
aSH=[4444;1111;8888;6666;3737;2929;5533;8181;6262;73.673.6];
cSH=[.1.2.2.4.4.6.3.7.5.5];o=0;fori=1:10o=o-((x-aSH(i,:))*(x-aSH(i,:))'+cSH(i))^(-1);endendfunction o=Ufun(x,a,k,m)
o=k.*((x-a).^m).*(x>a)+k.*((-x-a).^m).*(x<(-a));end
initialization.m
% This function initialize the first population of search agentsfunction Positions=initialization(SearchAgents_no,dim,ub,lb)Boundary_no=size(ub,2);% numnber of boundaries% If the boundaries of all variables are equal and user enter a signle% number for both ub and lbif Boundary_no==1Positions=rand(SearchAgents_no,dim).*(ub-lb)+lb;end% If each variable has a different lb and ubif Boundary_no>1fori=1:dimub_i=ub(i);lb_i=lb(i);Positions(:,i)=rand(SearchAgents_no,1).*(ub_i-lb_i)+lb_i;endend
WOA
% The Whale Optimization Algorithmfunction[Leader_score,Leader_pos,Convergence_curve]=WOA(SearchAgents_no,Max_iter,lb,ub,dim,fobj)% initialize position vector and score for the leader
Leader_pos=zeros(1,dim);
Leader_score=inf;%change this to -inf for maximization problems%Initialize the positions of search agents
Positions=initialization(SearchAgents_no,dim,ub,lb);Convergence_curve=zeros(1,Max_iter);t=0;% Loop counter% Main loopwhile t<Max_iterfori=1:size(Positions,1)% Return back the search agents that go beyond the boundaries of the search spaceFlag4ub=Positions(i,:)>ub;Flag4lb=Positions(i,:)<lb;Positions(i,:)=(Positions(i,:).*(~(Flag4ub+Flag4lb)))+ub.*Flag4ub+lb.*Flag4lb;% Calculate objective function for each search agentfitness=fobj(Positions(i,:));% Update the leaderif fitness<Leader_score % Change this to > for maximization problemLeader_score=fitness;% Update alphaLeader_pos=Positions(i,:);endenda=2-t*((2)/Max_iter);% a decreases linearly fron 2 to 0 in Eq. (2.3)% a2 linearly dicreases from -1 to -2 to calculate t in Eq. (3.12)a2=-1+t*((-1)/Max_iter);% Update the Position of search agents fori=1:size(Positions,1)r1=rand();% r1 is a random number in [0,1]r2=rand();% r2 is a random number in [0,1]A=2*a*r1-a;% Eq. (2.3) in the paperC=2*r2;% Eq. (2.4) in the paperb=1;% parameters in Eq. (2.5)l=(a2-1)*rand+1;% parameters in Eq. (2.5)p =rand();% p in Eq. (2.6)forj=1:size(Positions,2)if p<0.5ifabs(A)>=1rand_leader_index =floor(SearchAgents_no*rand()+1);X_rand =Positions(rand_leader_index,:);D_X_rand=abs(C*X_rand(j)-Positions(i,j));% Eq. (2.7)Positions(i,j)=X_rand(j)-A*D_X_rand;% Eq. (2.8)elseifabs(A)<1D_Leader=abs(C*Leader_pos(j)-Positions(i,j));% Eq. (2.1)Positions(i,j)=Leader_pos(j)-A*D_Leader;% Eq. (2.2)endelseif p>=0.5distance2Leader=abs(Leader_pos(j)-Positions(i,j));% Eq. (2.5)Positions(i,j)=distance2Leader*exp(b.*l).*cos(l.*2*pi)+Leader_pos(j);endendendt=t+1;Convergence_curve(t)=Leader_score;[t Leader_score]end
main.m
% You can simply define your cost in a seperate file and load its handle to fobj % The initial parameters that you need are:%__________________________________________% fobj = @YourCostFunction% dim = number of your variables% Max_iteration = maximum number of generations% SearchAgents_no = number of search agents% lb=[lb1,lb2,...,lbn] where lbn is the lower bound of variable n% ub=[ub1,ub2,...,ubn] where ubn is the upper bound of variable n% If all the variables have equal lower bound you can just% define lb and ub as two single number numbers% To run WOA: [Best_score,Best_pos,WOA_cg_curve]=WOA(SearchAgents_no,Max_iteration,lb,ub,dim,fobj)%__________________________________________clear all
clcSearchAgents_no=30;% Number of search agentsFunction_name='F1';% Name of the test function that can be from F1 to F23 (Table 1,2,3 in the paper)Max_iteration=500;% Maximum numbef of iterations% Load details of the selected benchmark function[lb,ub,dim,fobj]=Get_Functions_details(Function_name);[Best_score,Best_pos,WOA_cg_curve]=WOA(SearchAgents_no,Max_iteration,lb,ub,dim,fobj);figure('Position',[269240660290])%Draw search spacesubplot(1,2,1);func_plot(Function_name);title('Parameter space')xlabel('x_1');ylabel('x_2');zlabel([Function_name,'( x_1 , x_2 )'])%Draw objective spacesubplot(1,2,2);semilogy(WOA_cg_curve,'Color','r')title('Objective space')xlabel('Iteration');ylabel('Best score obtained so far');axis tight
grid on
box on
legend('WOA')display(['The best solution obtained by WOA is : ',num2str(Best_pos)]);display(['The best optimal value of the objective funciton found by WOA is : ',num2str(Best_score)]);