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

Hello Mat

 找回密码
 立即注册
查看: 6387|回复: 0

蚁群算法ACO的PID参数整定

[复制链接]

1294

主题

1520

帖子

110

金钱

管理员

Rank: 9Rank: 9Rank: 9

积分
22633
发表于 2019-10-28 21:12:37 | 显示全部楼层 |阅读模式
蚁群算法ACO的PID参数整定:
代码参考链接:链接:https://pan.baidu.com/s/15nqkVj9-WuwI-2HQ1gB7pg 提取码:howq
  1. clc,clear,close all
  2. warning off
  3. feature jit off
  4. load('U.mat')
  5. load('Y.mat')
  6. % 对Y进行抽样
  7. Y1 = Y(1:1000:end);
  8. % ACO 参数
  9. maxiter = 100;  % 迭代次数
  10. sizepop = 10;  % 种群数量
  11. popmin = [-1,-1,-1, 0.1, 0.1, 0.1 ]; % x_min
  12. popmax = [50,50,50, 15.0, 15.0, 15.0 ]; % x_max
  13. Rou = 0.8;     % 信息素增量强度
  14. P0 = 0.1;      % 转移概率   
  15. nvar = 6;      % nvar个自变量
  16. % 初始化种群
  17. for i=1:sizepop
  18.     pop(i,:) = popmin + (popmax-popmin).*rand(1,nvar)*0.5;
  19.     fitness(i) = BJ_Fun2( pop(i,:), Y1 );
  20. end
  21. % 记录一组最优值
  22. [bestfitness,bestindex]=min(fitness);
  23. zbest=pop(bestindex,:);   % 全局最佳
  24. gbest=pop;                % 个体最佳
  25. fitnessgbest=fitness;     % 个体最佳适应度值
  26. fitnesszbest=bestfitness; % 全局最佳适应度值
  27. %% 迭代寻优
  28. for i=1:maxiter
  29.     lamda = 1/i;       % 信息素挥发因子
  30.     [bestfit,flag] =min(fitness);
  31.     for j=1:sizepop
  32.         Pt(j) = (fitness(j)-bestfit)./bestfit;
  33.     end
  34.     for j=1:sizepop
  35.         % 转移概率值
  36.         if Pt(j)<P0
  37.             pop(j,:) = pop(j,:) + (2*rand(1, nvar)-1)*lamda/2;
  38.         else
  39.             pop(j,:) = pop(j,:) + (popmax-popmin).*(rand(1, nvar)-0.5);
  40.         end   
  41.         
  42.         % 越界处理
  43.         for k=1:nvar
  44.             if pop(j,k)>popmax(k)
  45.                 pop(j,k)=popmax(k);
  46.             end
  47.             if pop(j,k)<popmin(k)
  48.                 pop(j,k)=popmin(k);
  49.             end
  50.         end
  51.         
  52.         % 适应度值
  53.         fitness1 = BJ_Fun2(pop(j,:), Y1);
  54.         fitness2 = BJ_Fun2(gbest(j,:), Y1);
  55.         if fitness1<fitness2
  56.             gbest(j,:) = pop(j,:);
  57.         else
  58.             pop(j,:) = gbest(j,:);
  59.         end
  60.         
  61.         fitness(j) = (1-Rou)*fitness(j) + BJ_Fun2(pop(j,:), Y1);
  62.         
  63.         if fitness(j) < fitnesszbest
  64.             fitnesszbest = fitness(j);
  65.             fitnessgbest(j) = BJ_Fun2(pop(j,:), Y1);
  66.         end
  67.     end
  68.     [fitness_iter(i),index]= min(fitnessgbest);
  69. end
  70. zbest = pop(index,:);  % 最优个体
  71. save zbest.mat zbest
  72. %% 显示
  73. disp('最优解')
  74. disp(['最优解: ', num2str( zbest )])
  75. fprintf('\n')
  76. disp(['最优解对应的目标函数值: ', num2str( bestfitness )])
  77. disp('传递函数');
  78. sys = tf( zbest(4), [1, zbest(5), zbest(4)] )
  79. %% 画图显示
  80. figure('color',[1,1,1])
  81. plot(fitness_iter,'ro-','linewidth',2)
  82. grid on
  83. xlabel('迭代次数');
  84. ylabel('适应度函数值')
  85. %% 响应曲线
  86. Ygap = 100;
  87. [rint_BFO, yout_BFO] = BJ_Fun_response2(zbest, Y(1:Ygap:end));
  88. figure;hold on;
  89. plot(rint_BFO.*max(Y(1:Ygap:end)),'r-','linewidth',2);
  90. plot(yout_BFO.*max(Y(1:Ygap:end)),'b.-','linewidth',2);
  91. legend('阶跃控制信号','ACO算法优化的信号逼近')
  92. axis tight;
  93. hold off
复制代码
适应度函数如下:
  1. function fitness = BJ_Fun2(xx, Y1)
  2. ts=0.001;
  3. KP=xx(1); KI=xx(2); KD=xx(3);
  4. aa=xx(4); bb=xx(5); cc=xx(6);
  5. % sys = tf( [1.6], [1,1.5,1.6] );
  6. sys = tf( aa, [1, bb, aa] );
  7. dsys = c2d(sys, ts, 'z');
  8. [num,den] = tfdata( dsys, 'v' );
  9. u_1 = 0; u_2=0;
  10. y_1 = 0; y_2=0;
  11. error_1 = 0;
  12. x = [0,0,0]';
  13. for k=1:length(Y1)
  14.     timef(k) = k*ts;
  15.     r(k) = Y1(k)./(max(Y1)+eps);
  16.     u(k) = KP*x(1)+KI*x(2)+KD*x(3);
  17.     yout(k) = -den(2)*y_1 -den(3)*y_2 + num(2)*u_1 + num(3)*u_2;
  18.     error(k) = r(k) - yout(k);
  19.     % PID参数返回值
  20.     u_2 = u_1;
  21.     u_1 = u(k);
  22.     y_2 = y_1;
  23.     y_1 = yout(k);
  24.    
  25.     x(1) = error(k);              % Kp
  26.     x(2) = x(2) + error(k)*ts;    % KI
  27.     x(3) = (error(k)-error_1)/ts; % KD
  28.    
  29.     error_1 = error(k);
  30. end
  31. fitness = sum( abs(error) ) + 300*min( abs( error ) );
复制代码
响应度曲线:
  1. function [r, yout] = BJ_Fun_response2(xx, Y1)
  2. KP=xx(1); KI=xx(2); KD=xx(3);
  3. aa=xx(4); bb=xx(5); cc=xx(6);
  4. ts=0.001;
  5. % sys = tf( [1.6], [1,1.5,1.6] );
  6. sys = tf( aa, [1, bb, aa] );
  7. dsys = c2d(sys, ts, 'z');
  8. [num,den] = tfdata( dsys, 'v' );
  9. u_1 = 0; u_2=0;
  10. y_1 = 0; y_2=0;
  11. error_1 = 0;
  12. x = [0,0,0]';
  13. for k=1:length(Y1)
  14.     timef(k) = k*ts;
  15.     r(k) = Y1(k)./(max(Y1)+eps);
  16.     u(k) = KP*x(1)+KI*x(2)+KD*x(3);
  17.     yout(k) = -den(2)*y_1 -den(3)*y_2 + num(2)*u_1 + num(3)*u_2;
  18.     error(k) = r(k) - yout(k);
  19.     % PID参数返回值
  20.     u_2 = u_1;
  21.     u_1 = u(k);
  22.     y_2 = y_1;
  23.     y_1 = yout(k);
  24.    
  25.     x(1) = error(k);              % Kp
  26.     x(2) = x(2) + error(k)*ts;    % KI
  27.     x(3) = (error(k)-error_1)/ts; % KD
  28.    
  29.     error_1 = error(k);
  30. end
  31. % 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-4-24 07:05 , Processed in 0.235196 second(s), 25 queries .

Powered by Discuz! X3.4

Copyright © 2001-2021, Tencent Cloud.

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