example5_41.m

MATLAB7辅助控制系统设计与仿真 第5章

GGAP = .8;          %代间隔，即每代产生多少个新个体
INSR = .9;           % 插入率，即插入多少个子代
XOVR =  1;           %交叉率
SP = 2;              % 选择压力
MUTR = 1;            % 变异率
MIGR = 0.2;          % 子代之间的迁移率
MIGGEN = 20;         % 用于迁移的子代数目

TERMEXACT = 1e-4;    % 当达到最小值时的终止值

SEL_F = 'sus';       % 选择函数名
XOV_F = 'recdis';    % 个体组合函数名
MUT_F = 'mutbga';    % 变异函数名

OBJ_F = 'objlinq'; % 目标函数名
% 得到目标函数的边界
   FieldDR = feval(OBJ_F,[],1);

% 根据在目标函数中定义的变量的数目计算SUBPOP和NIND
   NVAR = size(FieldDR,2);          % 变量的数目
   SUBPOP = 2 * floor(sqrt(NVAR));   % 子代的数目
   NIND = 20 + 5 * floor(NVAR/50);   % 每个子代中个体的数目
   MAXGEN = 300 * floor(sqrt(NVAR)); % 最大代数
   MUTR = MUTR / NVAR;           %依赖于NVAR的变异率

% 得到目标函数的最小值
   GlobalMin = feval(OBJ_F,[],3);

% 目标函数图形输出的标题
   FigTitle = [feval(OBJ_F,[],2) '   (' int2str(SUBPOP) ':' int2str(MAXGEN) ') '];

% 清除Bes和存储矩阵内容
   %初始化存储最优结果的矩阵
      Best = NaN * ones(MAXGEN,3);
      Best(:,3) = zeros(size(Best,1),1);
   % 存储最优个体的矩阵
      IndAll = [];

%生成实值种群
Chrom = crtrp(SUBPOP*NIND,FieldDR);
% 计数变量归零
gen = 0;
termopt = 0;

% 计算给定种群的目标函数值
   ObjV = feval(OBJ_F,Chrom);
   %记录目标函数估计的次数count number of objective function evaluations
   Best(gen+1,3) = Best(gen+1,3) + NIND;

%迭代运算直至符合终止条件或达到了最大代数
   while ((gen < MAXGEN) & (termopt == 0)),

   % 保存最优和平均的目标函数值以及最优的个体
      [Best(gen+1,1),ix] = min(ObjV);
      Best(gen+1,2) = mean(ObjV);
      IndAll = [IndAll; Chrom(ix,:)];

   % 适应度分配
      FitnV = ranking(ObjV,[2 0],SUBPOP);
            
   % 从中群众选择个体
      SelCh = select(SEL_F, Chrom, FitnV, GGAP, SUBPOP);
      
   % 重新组合选择的个体
      SelCh=recombin(XOV_F, SelCh, XOVR, SUBPOP);

   % 变异
      SelCh=mutate(MUT_F, SelCh, FieldDR, [MUTR], SUBPOP);

   % 计算对应子代的目标函数值
      ObjVOff = feval(OBJ_F,SelCh);
      Best(gen+1,3) = Best(gen+1,3) + size(SelCh,1);

   %在种群仲插入最优的子代以替换最劣的父代
      [Chrom, ObjV] = reins(Chrom, SelCh, SUBPOP, [1 INSR], ObjV, ObjVOff);

      gen=gen+1;

   % 检验, 如果最优的目标函数值接近于终止值，
   % 计算两者的差
      ActualMin = abs(min(ObjV) - GlobalMin);
   %如果最优的目标函数值小于终止值，
      if ((ActualMin < (TERMEXACT * abs(GlobalMin))) | (ActualMin < TERMEXACT))
         termopt = 1;
end
% 子代迁移
if ((termopt ~= 1) & (rem(gen,MIGGEN) == 0))
[Chrom, ObjV] = migrate(Chrom, SUBPOP, [MIGR, 1, 0], ObjV);
end
end
% 结果
   % 叠加目标函数估计的次数
   Results = cumsum(Best(1:gen,3));
   % 函数估计的次数、均值和最优结果
   Results = [Results Best(1:gen,2) Best(1:gen,1)];
   
% 绘图
   figure('Name',['Results of ' FigTitle]);
   subplot(2,1,1), plot(Results(:,1),Results(:,2),'-',Results(:,1),Results(:,3),':');
   subplot(2,1,2), plot(IndAll(gen-4:gen,:)');
   xlabel('u(k)');
   ylabel('f(x,u)');