这篇文章主要讲解的是使用粒子群算法对PID参数进行寻优,大家可以进行参考。
基于粒子群算法的PID控制器优化设计
1. 理论基础
PID控制器的一般形式为
u
(
t
)
=
K
p
e
(
t
)
+
K
i
∫
0
t
e
(
τ
)
d
τ
+
K
d
d
e
(
t
)
d
t
u(t) = {K_p}e(t) + {K_i}\int_0^t {e(\tau )} d\tau + {K_d}\frac{{de(t)}}{{dt}}
u(t)=Kpe(t)+Ki∫0te(τ)dτ+Kddtde(t)
其中,
e
(
t
)
e(t)
e(t)是系统误差;
K
p
K_p
Kp、
K
i
K_i
Ki和
K
d
K_d
Kd分别是对系统误差信号及其积分与微分量的加权,控制器通过这样的加权就可以计算出控制信号,驱动受控对象。
因此,PID控制器的三个参数的选择相当重要,优化这三个参数相当重要。传统参数的选取都是靠人工经验进行选取,此文章将使用PSO对PID控制器中的这三个参数进行优化。
2. 问题描述
PID控制器的系统结构图如下
PID控制器的优化问题就是确定一组合适的参数
K
p
K_p
Kp、
K
i
K_i
Ki、
K
d
K_d
Kd,使得指标达到最优。这里选用ITAE指标,其定义为
J
=
∫
0
∞
t
∣
e
(
t
)
∣
d
t
J = \int_0^\infty {t\left| {e(t)} \right|dt}
J=∫0∞t∣e(t)∣dt
选取的被控对象为为下面的不稳定系统:
G
(
s
)
=
s
+
2
s
4
+
8
s
3
+
4
s
2
−
s
+
0.4
G(s) = \frac{{s + 2}}{{{s^4} + 8{s^3} + 4{s^2} - s + 0.4}}
G(s)=s4+8s3+4s2−s+0.4s+2
在Simulink环境下搭建上述的模型如下图所示:
3. 思路及步骤
3.1 优化设计过程
李用粒子群算法对PID控制器的参数进行优化设计,过程图如下所示:
上图中,粒子群算法与Simulink模型之间连接的桥梁是粒子(即PID控制器参数)和该粒子对应的适应值(即控制系统的性能指标)。优化过程如下:PSO产生粒子群(可以是初始化的粒子群,也可以是更新后的粒子群),将该粒子群中的粒子依次赋值给PID控制器的参数
K
p
K_p
Kp、
K
i
K_i
Ki、
K
d
K_d
Kd,然后运行控制器系统的Simulink模型,得到该参数对应的性能指标,该性能指标传递到PSO中作为该粒子的适应值,最后判断是否可以退出算法。
3.2 粒子群算法实现
粒子群算法中速度和位置是根据下面两个公式进行更新
v
t
+
1
=
ω
v
t
+
c
1
r
1
(
P
t
−
x
t
)
+
c
2
r
2
(
G
t
−
x
t
)
{v_{t + 1}} = \omega {v_t} + {c_1}{r_1}({P_t} - {x_t}) + {c_2}{r_2}({G_t} - {x_t})
vt+1=ωvt+c1r1(Pt−xt)+c2r2(Gt−xt)
x
t
+
1
=
x
t
+
v
t
+
1
x_{t+1}=x_t+v_{t+1}
xt+1=xt+vt+1
其中,
x
x
x表示粒子的位置;
v
v
v表示粒子的速度;
ω
\omega
ω为惯性因子;
c
1
c_1
c1、
c
2
c_2
c2为加速度常数;
r
1
r_1
r1、
r
2
r_2
r2为
[
0
,
1
]
[0,1]
[0,1]区间的随机数;
P
t
P_t
Pt是粒子迄今为止搜索到的最优为止;
G
t
G_t
Gt是整个粒子群迄今为止搜索到的最优为止。
PSO的流程如下:
- 初始化粒子群,随机产生所有粒子的位置和速度,并确定粒子的 P t P_t Pt和 G t G_t Gt。
- 对每个粒子,将其适应值与该粒子所经历过的最优位置 P t P_t Pt的适应值进行比较,若较好,则将其作为当前的 P t P_t Pt。
- 对每个粒子,将其适应值与整个粒子群所经历过的最优位置 G t G_t Gt的适应值进行比较,若较好,则将其作为当前的 G t G_t Gt。
- 按照上面的公式进行速度和位置的更新。
- 如果没有满足终止条件(通常为预设的最大迭代次数和适应值下限值),则返回步骤(2);否则,推出算法,得到最优解。
4. MATLAB程序
下面是Simulink图
其中函数部分的程序如下:
function z=PSO_PID(x)
assignin('base','Kp',x(1)); %粒子依次赋值给Kp
assignin('base','Ki',x(2)); %粒子依次赋值给Ki
assignin('base','Kd',x(3)); %粒子依次赋值给Kd
[t_time,x_state,y_out]=sim('PID_Model',[0,20]); %使用命令行运行控制系统模型
z=y_out(end,1); %返回性能指标
主程序如下:
%% 清空环境
clear
clc
%% 参数设置
w = 0.6; % 惯性因子
c1 = 2; % 加速常数
c2 = 2; % 加速常数
Dim = 3; % 维数
SwarmSize = 100; % 粒子群规模
ObjFun = @PSO_PID; % 待优化函数句柄
MaxIter = 100; % 最大迭代次数
MinFit = 0.1; % 最小适应值
Vmax = 1;
Vmin = -1;
Ub = [300 300 300];
Lb = [0 0 0];
%% 粒子群初始化
Range = ones(SwarmSize,1)*(Ub-Lb);
Swarm = rand(SwarmSize,Dim).*Range + ones(SwarmSize,1)*Lb % 初始化粒子群
VStep = rand(SwarmSize,Dim)*(Vmax-Vmin) + Vmin % 初始化速度
fSwarm = zeros(SwarmSize,1);
for i=1:SwarmSize
fSwarm(i,:) = feval(ObjFun,Swarm(i,:)); % 粒子群的适应值
end
%% 个体极值和群体极值
[bestf bestindex]=min(fSwarm);
zbest=Swarm(bestindex,:); % 全局最佳
gbest=Swarm; % 个体最佳
fgbest=fSwarm; % 个体最佳适应值
fzbest=bestf; % 全局最佳适应值
%% 迭代寻优
iter = 0;
y_fitness = zeros(1,MaxIter); % 预先产生4个空矩阵
K_p = zeros(1,MaxIter);
K_i = zeros(1,MaxIter);
K_d = zeros(1,MaxIter);
while( (iter < MaxIter) && (fzbest > MinFit) )
for j=1:SwarmSize
% 速度更新
VStep(j,:) = w*VStep(j,:) + c1*rand*(gbest(j,:) - Swarm(j,:)) + c2*rand*(zbest - Swarm(j,:));
if VStep(j,:)>Vmax, VStep(j,:)=Vmax; end
if VStep(j,:)<Vmin, VStep(j,:)=Vmin; end
% 位置更新
Swarm(j,:)=Swarm(j,:)+VStep(j,:);
for k=1:Dim
if Swarm(j,k)>Ub(k), Swarm(j,k)=Ub(k); end
if Swarm(j,k)<Lb(k), Swarm(j,k)=Lb(k); end
end
% 适应值
fSwarm(j,:) = feval(ObjFun,Swarm(j,:));
% 个体最优更新
if fSwarm(j) < fgbest(j)
gbest(j,:) = Swarm(j,:);
fgbest(j) = fSwarm(j);
end
% 群体最优更新
if fSwarm(j) < fzbest
zbest = Swarm(j,:);
fzbest = fSwarm(j);
end
end
iter = iter+1; % 迭代次数更新
y_fitness(1,iter) = fzbest; % 为绘图做准备
K_p(1,iter) = zbest(1);
K_i(1,iter) = zbest(2);
K_d(1,iter) = zbest(3);
end
%% 绘图输出
figure(1) % 绘制性能指标ITAE的变化曲线
plot(y_fitness,'LineWidth',2)
title('最优个体适应值','fontsize',18);
xlabel('迭代次数','fontsize',18);ylabel('适应值','fontsize',18);
set(gca,'Fontsize',18);
figure(2) % 绘制PID控制器参数变化曲线
plot(K_p)
hold on
plot(K_i,'k','LineWidth',3)
plot(K_d,'--r')
title('Kp、Ki、Kd 优化曲线','fontsize',18);
xlabel('迭代次数','fontsize',18);ylabel('参数值','fontsize',18);
set(gca,'Fontsize',18);
legend('Kp','Ki','Kd',1);
5.仿真结果
经过优化后的PID控制器的最优参数以及性能指标为
K
p
=
33.6326
,
K
i
=
0.1662
K_p= 33.6326, K_i=0.1662
Kp=33.6326,Ki=0.1662
K
d
=
38.7892
,
I
T
A
E
=
1.0580
K_d= 38.7892 , ITAE= 1.0580
Kd=38.7892,ITAE=1.0580
性能指标曲线图:
单位阶跃响应曲线:
K
p
K_p
Kp、
K
i
K_i
Ki、
K
d
K_d
Kd优化曲线:
大家如果有使用PSO算法对LQR中的Q参数进行优化的程序可以进行交流一下,谢谢!