请选择 进入手机版 | 继续访问电脑版

Hello Mat

 找回密码
 立即注册
查看: 5501|回复: 7

蚁群算法ACO优化的控制系统参数逼近算法

[复制链接]

1277

主题

1500

帖子

85

金钱

管理员

Rank: 9Rank: 9Rank: 9

积分
22528
发表于 2019-10-28 21:11:49 | 显示全部楼层 |阅读模式
蚁群算法ACO优化的控制系统参数逼近算法
参考链接:链接:https://pan.baidu.com/s/1ow0aKl3jcyvZrsFQAzghXA 提取码:ix7u
  1. clc,clear,close all
  2. warning off
  3. feature jit off
  4. load('U.mat')
  5. load('Y.mat')
  6. % zbest_tf = importfile1('zbest.mat');
  7. % 对Y进行抽样
  8. Y1 = Y(1:1000:end);
  9. % ACO 参数
  10. maxiter = 100;  % 迭代次数
  11. sizepop = 30;   % 种群数量
  12. popmin = [0.1,14,0.1]; % x_min
  13. popmax = [100,20,100]; % x_max
  14. Rou = 1.0;      % 信息素增量强度
  15. P0 = 0.01;      % 转移概率   
  16. nvar = 3;       % nvar个自变量
  17. % 初始化种群
  18. for i=1:sizepop
  19.     pop(i,:) = popmin + (popmax-popmin).*rand(1,nvar)*0.5;
  20.     fitness(i) = BJ_Fun4( pop(i,:), Y1 );
  21. end
  22. % 记录一组最优值
  23. [bestfitness,bestindex]=min(fitness);
  24. zbest=pop(bestindex,:);   % 全局最佳
  25. gbest=pop;                % 个体最佳
  26. fitnessgbest=fitness;     % 个体最佳适应度值
  27. fitnesszbest=bestfitness; % 全局最佳适应度值
  28. %% 迭代寻优
  29. for i=1:maxiter
  30.     disp(['当前迭代次数:', num2str(i)])
  31.     lamda = 1/i;       % 信息素挥发因子
  32.     [bestfit,flag] =min(fitness);
  33.     for j=1:sizepop
  34.         Pt(j) = (fitness(j)-bestfit)./bestfit;
  35.     end
  36.     for j=1:sizepop
  37.         % 转移概率值
  38.         if Pt(j)<P0
  39.             pop(j,:) = pop(j,:) + (2*rand(1, nvar)-1)*lamda/2;
  40.         else
  41.             pop(j,:) = pop(j,:) + (popmax-popmin).*(rand(1, nvar));
  42.         end   
  43.         
  44.         % 越界处理
  45.         for k=1:nvar
  46.             if pop(j,k)>popmax(k)
  47.                 pop(j,k)=popmax(k);
  48.             end
  49.             if pop(j,k)<popmin(k)
  50.                 pop(j,k)=popmin(k);
  51.             end
  52.         end
  53.         
  54.         % 适应度值
  55.         fitness1 = BJ_Fun4(pop(j,:), Y1);
  56.         fitness2 = BJ_Fun4(gbest(j,:), Y1);
  57.         if fitness1<fitness2
  58.             gbest(j,:) = pop(j,:);
  59.         else
  60.             pop(j,:) = gbest(j,:);
  61.         end
  62.         
  63.         fitness(j) = (1-Rou)*fitness(j) + BJ_Fun4(pop(j,:), Y1);
  64.         
  65.         if fitness(j) < fitnesszbest
  66.             fitnesszbest = fitness(j);
  67.             fitnessgbest(j) = BJ_Fun4(pop(j,:), Y1);
  68.         end
  69.     end
  70.     [fitness_iter(i),index]= min(fitnessgbest);
  71. end
  72. zbest = pop(index,:);  % 最优个体
  73. save zbest.mat zbest
  74. %% 显示
  75. disp('最优解')
  76. disp(['最优解: ', num2str( zbest )])
  77. fprintf('\n')
  78. disp(['最优解对应的目标函数值: ', num2str( min(fitness_iter) )])
  79. disp('传递函数');
  80. sys = tf( zbest(1), [1, zbest(2), zbest(1)] )
  81. %% 画图显示
  82. figure('color',[1,1,1])
  83. plot(fitness_iter,'ro-','linewidth',2)
  84. grid on
  85. xlabel('迭代次数');
  86. ylabel('适应度函数值');
  87. %% 响应曲线
  88. Ygap = 1000;
  89. [rint_BFO, yout_BFO] = BJ_Fun_response4(zbest, Y(1:Ygap:end));
  90. figure;hold on;
  91. plot(rint_BFO.*max(Y(1:Ygap:end)),'r-','linewidth',2);
  92. plot(yout_BFO.*max(Y(1:Ygap:end)),'b.-','linewidth',2);
  93. legend('阶跃控制信号','ACO算法优化的信号逼近')
  94. axis tight;
  95. hold off

  96. %% 响应曲线
  97. Amp= 5000;
  98. Ygap = 1;
  99. [rint_BFO2, yout_BFO2] = BJ_Fun_response5( Y(1:Ygap:end), zbest);
  100. figure;hold on;
  101. plot(rint_BFO2.*Amp,'k-','linewidth',2);
  102. plot(Y(1:Ygap:end),'r-','linewidth',2);
  103. plot(yout_BFO2*Amp,'b.-','linewidth',2);
  104. legend('阶跃控制信号','实际采集信号','ACO算法优化的信号逼近')
  105. axis tight;
  106. hold off
复制代码
适应度函数如下:
  1. function fitness = BJ_Fun4(xx, Y1)
  2. ts=1e-2;
  3. aa=xx(1); bb=xx(2); cc=xx(3);
  4. % sys = tf( [1.6], [1,1.5,1.6] );
  5. sys = tf( aa, [1, bb, aa] );
  6. % dsys = c2d(sys, ts, 'z');
  7. % [num,den] = tfdata( dsys, 'v' );
  8. for k=1:length(Y1)
  9.     timef(k) = k*ts;
  10.     r(k) = 1.0;
  11. %     r(k) = Y1(k)./(max(Y1)+eps);
  12. %     r(k) = stepfun(timef(k),0);
  13.     u(k) = Y1(k)./(max(Y1)+eps);
  14. end
  15. rout =  lsim(sys, r, timef);
  16. yout =  lsim(sys, u, timef);
  17. error1 = r - rout';
  18. error2 = u - yout';
  19. fitness1 = sum( abs(error1) ) + 100*max( abs( error1 ) );
  20. fitness2 = sum( abs(error2) ) + 0*max( abs( error2 ) );
  21. fitness = fitness1+fitness2;
复制代码
响应度曲线:
  1. function [r, yout] = BJ_Fun_response4(xx, Y1)
  2. ts=1e-0;
  3. aa=xx(1); bb=xx(2); cc=xx(3);
  4. % sys = tf( [1.6], [1,1.5,1.6] );
  5. sys = tf( aa, [1, bb, aa] );
  6. % dsys = c2d(sys, ts, 'z');
  7. % [num,den] = tfdata( dsys, 'v' );
  8. for k=1:length(Y1)
  9.     timef(k) = k*ts;
  10. %     r(k) = 1.0;
  11.     r(k) = Y1(k)./(max(Y1)+eps);
  12. %     r(k) = stepfun(timef(k),0);
  13. end
  14. yout =  lsim(sys, r, timef);
  15. % error = r - yout';
  16. % fitness = sum( abs(error) ) + 300*min( abs( error ) );
复制代码









算法QQ  3283892722
群智能算法链接http://halcom.cn/forum.php?mod=forumdisplay&fid=73
回复

使用道具 举报

您需要登录后才可以回帖 登录 | 立即注册

本版积分规则

Python|Opencv|MATLAB|Halcom.cn ( 蜀ICP备16027072号 )

GMT+8, 2024-3-5 02:57 , Processed in 0.228785 second(s), 24 queries .

Powered by Discuz! X3.4

Copyright © 2001-2021, Tencent Cloud.

快速回复 返回顶部 返回列表