人工蜂群算法MATLAB代码

%%%%%ARTIFICIAL BEE COLONY ALGORITHM%%%% %Artificial Bee Colony Algorithm was developed by Dervis Karaboga in 2005 %by simulating the foraging behaviour of bees. %Copyright ? 2008 Erciyes University, Intelligent Systems Research Group, The Dept. of Computer Engineering %Contact: %Dervis Karaboga (karaboga@erciyes.edu.tr ) %Bahriye Basturk Akay (bahriye@erciyes.edu.tr) clear all close all clc %% Set ABC Control Parameters ABCOpts = struct( 'ColonySize', 100, ... % Number of Employed Bees+ Number of Onlooker Bees % 雇佣蜂 和 观察蜂 的数量 'MaxCycles', 100,... % Maximum cycle number in order to terminate the algorithm % 算法终止条件 最大迭代次数 'ErrGoal', 1e-20, ... % Error goal in order to terminate the algorithm (not used in the code in current version) % 目标误差 'Dim', 5, ... % Number of parameters of the objective function % 目标函数自由变量的维度 （决策变量维度） 'Limit', 50, ... % Control paramter in order to abandone the food source % limit限制参数 'lb', -600, ... % Lower bound of the parameters to be optimized % 低阶 'ub', 600, ... %Upper bound of the parameters to be optimized % 高阶 'ObjFun' , 'griewank', ... %Write the name of the objective function you want to minimize 'RunTime',1); % Number of the runs % 程序运行次数 GlobalMins=zeros(ABCOpts.RunTime,ABCOpts.MaxCycles); % 记录 每次迭代的最小值 % 注意变量 Globa Mins 和 变量GlobalMin %% for r=1:ABCOpts.RunTime %% 初始化操作 % Initialise population Range = repmat((ABCOpts.ub-ABCOpts.lb),[ABCOpts.ColonySize ABCOpts.Dim]); % 维 度 区 间 Lower = repmat(ABCOpts.lb, [ABCOpts.ColonySize ABCOpts.Dim]); % 低 阶 Colony = rand(ABCOpts.ColonySize,ABCOpts.Dim) .* Range + Lower; % 初 始 化 Employed=Colony(1:(ABCOpts.ColonySize/2),:); %evaluate and calculate fitness ObjEmp=feval(ABCOpts.ObjFun,Employed); % 目标函数值 FitEmp=calculateFitness(ObjEmp); % 适应度值 %set initial values of Bas Bas=zeros(1,(ABCOpts.ColonySize/2)); % 记录邻域更新 limit相关量 Bas初始值 1 行 (ABCOpts.ColonySize/2) 列矩阵 %% GlobalMin=ObjEmp( find(ObjEmp==min(ObjEmp),1,'last') ); % 当前 最优值 GlobalParams=Employed( find(ObjEmp==min(ObjEmp),1,'last'),:); % 最优雇佣蜂 Cycle=1; while ((Cycle <= ABCOpts.MaxCycles)) %% Employed phase Employed2=Employed; % for循环 更新雇佣蜂 （更新邻域） for i=1:ABCOpts.ColonySize/2 Param2Change=fix(rand*ABCOpts.Dim)+1; % 确定 更新的维数 neighbour=fix(rand*(ABCOpts.ColonySize/2))+1; % 确定用于更新的另一个 蜜源 while(neighbour==i) % 使得用于更新的蜜源不等于 原蜜源 neighbour=fix(rand*(ABCOpts.ColonySize/2))+1; end; Employed2(i,Param2Change)=Employed(i,Param2Change)+( Employed(i,Param2Change)-Employed(neighbour,Param2Change))*(rand-0.5)*2; if (Employed2(i,Param2Change)<ABCOpts.lb) % 边界约束处理 Employed2(i,Param2Change)=ABCOpts.lb; end; if (Employed2(i,Param2Change)>ABCOpts.ub) Employed2(i,Param2Change)=ABCOpts.ub; end; end; %for end ObjEmp2=feval(ABCOpts.ObjFun,Employed2); FitEmp2=calculateFitness(ObjEmp2); [Employed ObjEmp FitEmp Bas]=GreedySelection(Employed,Employed2,ObjEmp,ObjEmp2,FitEmp,FitEmp2,Bas,ABCOpts); % 进行一次邻域搜索 % 更新蜜源 %% Normalize NormFit=FitEmp/sum(FitEmp); %基与轮盘赌的概率选择 计算得到各蜜源的概率值 %% Onlooker phase 观察蜂 操作 Employed2=Employed; i=1; t=0; while(t<ABCOpts.ColonySize/2) if(rand<NormFit(i)) t=t+1; Param2Change=fix(rand*ABCOpts.Dim)+1; neighbour=fix(rand*(ABCOpts.ColonySize/2))+1; while(neighbour==i) neighbour=fix(rand*(ABCOpts.ColonySize/2))+1; end; Employed2(i,:)=Employed(i,:); Employed2(i,Param2Change)=Employed(i,Param2Change)+(Employed(i,Param2Change)-Employed(neighbour,Param2Change))*(rand-0.5)*2;%find a new value in the neighborhood if (Employed2(i,Param2Change)<ABCOpts.lb) % 边界 约束 处理 Employed2(i,Param2Change)=ABCOpts.lb; end; if (Employed2(i,Param2Change)>ABCOpts.ub) Employed2(i,Param2Change)=ABCOpts.ub; end; ObjEmp2=feval(ABCOpts.ObjFun,Employed2); FitEmp2=calculateFitness(ObjEmp2); [Employed ObjEmp FitEmp Bas]=GreedySelection(Employed,Employed2,ObjEmp,ObjEmp2,FitEmp,FitEmp2,Bas,ABCOpts,i); % 贪婪选择 end; % end if i=i+1; if (i==(ABCOpts.ColonySize/2)+1) i=1; end; end; % end while %% Memorize Best CycleBestIndex=find(FitEmp==max(FitEmp)); CycleBestIndex=CycleBestIndex(end); CycleBestParams=Employed(CycleBestIndex,:); CycleMin=ObjEmp(CycleBestIndex); if CycleMin<GlobalMin GlobalMin=CycleMin; GlobalParams=CycleBestParams; end GlobalMins(r,Cycle)=GlobalMin; %% Scout phase ind=find(Bas==max(Bas)); ind=ind(end); if (Bas(ind)>ABCOpts.Limit) Bas(ind)=0; % 重新初始化 第 ind 个 雇佣蜂 所对应 的目标函数值 Employed(ind,:)=(ABCOpts.ub-ABCOpts.lb)*(0.5-rand(1,ABCOpts.Dim));%*2+ABCOpts.lb; %message=strcat('burada',num2str(ind)) end; % 重新计算 雇佣蜂所对应的蜜源的目标函数值及其适应度值 ObjEmp=feval(ABCOpts.ObjFun,Employed); FitEmp=calculateFitness(ObjEmp); % fprintf('Cycle=%d ObjVal=%g\n',Cycle,GlobalMin); Cycle=Cycle+1; end % End of ABC end; %end of runs figure(1); set(gcf,'Color',[1 1 1]); if ABCOpts.RunTime==1 semilogy(GlobalMins); else semilogy(mean(GlobalMins),'ab13456'); end title('Mean of Best function values'); xlabel('cycles'); ylabel('error'); fprintf('Mean =%g Std=%g\n',mean(GlobalMins(:,end)),std(GlobalMins(:,end)));