【优化求解】基于风驱动算法WDO求解多目标最优matlab源码

106 阅读3分钟

1 简介

img

img

img

2 部分代码

clc
clear all
close all
%-------------------------------------------------------------------------
% Sample Matlab Code for the Adaptive Wind Driven Optimization.
% The coefficients of the WDO are optimized by the CMEAS algorithm.
% No additional cost function calls are needed.
%

tic; clear; close all; clc;
format long g;
delete('WDOoutput.txt');
delete('WDOpressure.txt');
delete('WDOposition.txt');
fid=fopen('WDOoutput.txt','a');
%--------------------------------------------------------------

% User defined parameters:
param.popsize = 20;     % population size.
param.npar = 10;% Dimension of the problem.
param.maxit = 100;% Maximum number of iterations.
maxV = 0.3;             % maximum allowed speed.
dimMin = -100;% Lower dimension boundary.
dimMax = 100;% Upper dimension boundary.

% Unlike the classical WDO, Adaptive WDO does not need the coefficients to
% be predetermined. Coefficients; alpha, RT, g, and c, will be selected by
% the CMAES algorithm:
rec.arx = rand(4,param.popsize);   %consistent with the CMAES indexing
%---------------------------------------------------------------

% Initialize WDO population, position and velocity:
% Randomize population in the range of [-1, 1]:
% Please note that WDO/AWDO always works in the range of [-1, 1], but is mapped
% to the upper/lower bounds of the problem before calling the pressure (cost,fitness) function.

pos = 2*(rand(param.popsize,param.npar)-0.5);
% Randomize velocity:
vel = maxV * 2 * (rand(param.popsize,param.npar)-0.5);
%---------------------------------------------------------------

% Call the pressure (cost,fitness) function
% Evaluate initial population one member at a time: (simple Sphere Function is utilized here)
for K=1:param.popsize,
  % map the AWDO 'pos' vector to problem itself using the upper/lower boundaries of each dimension.
  x = (dimMax - dimMin) * ((pos(K,:)+1)./2) + dimMin;
  % call the pressure (cost,fitness) function:
  pres(K,:) = fsphere(x);   %here insert your own cost function and make sure the mapping above carried out properly.
end
%----------------------------------------------------------------

% Finding best air parcel in the initial population :
[globalpres,indx] = min(pres);
globalpos = pos(indx,:);
minpres(1) = min(pres);% minimum pressure
%-----------------------------------------------------------------

% Rank the air parcels:
[sorted_pres rank_ind] = sort(pres);
% Do not sort the air parcels! Sorting mixes up CMAES indexing.
% pos = pos(rank_ind,:);
keepglob(1) = globalpres;
%-----------------------------------------------------------------

% Start iterations :
iter = 1;   % iteration counter
for ij = 2:param.maxit,
  % Update the velocity:
  for i=1:param.popsize
      % choose random dimensions:
      a = randperm(param.npar);
      % choose velocity based on random dimension:
      velot(i,:) = vel(i,a);
      vel(i,:) = (1-rec.arx(1,i))*vel(i,:)-(rec.arx(2,i)*pos(i,:))+ ...
          abs(1-1/rank_ind(i))*((globalpos-pos(i,:)).*rec.arx(3,i))+ ...
        (rec.arx(4,i)*velot(i,:)/rank_ind(i));
  end
  
  % Check velocity:
  vel = min(vel, maxV);
  vel = max(vel, -maxV);
  % Update air parcel positions:
  pos = pos + vel;
  pos = min(pos, 1.0);
  pos = max(pos, -1.0);
  
  % Call the pressure (cost,fitness) function
  % Evaluate initial population one member at a time: (simple Sphere Function is utilized here)
  for K=1:param.popsize,
      % map the AWDO 'pos' vector to problem itself using the upper/lower boundaries of each dimension.
      x = (dimMax - dimMin) * ((pos(K,:)+1)./2) + dimMin;
      % call the pressure (cost,fitness) function:
      pres(K,:) = fsphere(x);   %here insert your own cost function and make sure the mapping above carried out properly.
  end
  
  % call CMAES with the coefficients and corresponding pressure, so
  % that CMAES can return the new set of coefficient values for next
  % iteration:
[rec] = purecmaes_wdo(ij,rec,param.popsize,pres);
  
  %----------------------------------------------------
  % Finding best particle in population
[minpres,indx] = min(pres);
  minpos = pos(indx,:);         % min location for this iteration
  %----------------------------------------------------
  % Rank the air parcels:
[sorted_pres rank_ind] = sort(pres);
  % Do not sort the air parcels position, velocity and pressure!
  % Instead use "rank_ind" if needed.
  % pos = pos(rank_ind,:);
  % vel = vel(rank_ind,:);
  % pres = sorted_pres;
  
  % Updating the global best:
  better = minpres < globalpres;
  if better
      globalpres = minpres% global minimum pressure (cost,fitness) value
      globalpos = minpos;% global minimum position vector, note that it is in the range of [-1, 1].
  end
  % Keep a record of the progress:
  keepglob(ij) = globalpres;
  %     save WDOposition.txt pos -ascii -tabs;
end
%Save values to the final file.
pressure = transpose(keepglob);
filenamestr = ['WDOpressure.txt'];
save(filenamestr, 'pressure', '-ascii' , '-tabs');

% note that the 'globalpos' is the best solution vector found.
% 'globalpos' is in the range of [-1 1] and needs to be scaled with upper/lower bounds of the specific problem that you are trying to solve.

figure(1)
plot(keepglob)
xlabel('迭代次数')
ylabel('适应度函数')
% end-of-WDO.
%----------------------------------------------------------------------
%----------------------------------------------------------------------
%----------------------------------------------------------------------

% end-of-CMAES.
%----------------------------------------------------------------------
%----------------------------------------------------------------------
%----------------------------------------------------------------------

3 仿真结果

4 参考文献

[1]陈彬彬, 曹中清, 余胜威. 基于风驱动优化算法WDO的PID参数优化[J]. 计算机工程与应用, 2016, 52(014):250-253.

[2]寻琦, 赵惠玲, and 王肇平. "一种基于自适应风驱动优化算法的天线阵方向图综合方法.".

部分理论引用网络文献,若有侵权联系博主删除。

5 MATLAB代码与数据下载地址

见博客主页