关键字：专家PID控制；位置式PID；增量式PID；MATLAB。
文章目录
前言
前言为作者的碎碎念，想看知识点的请转跳：
在作为电气工程师的今天对自己的前途与工作感到迷茫，于是想看看如果我们控制专业的人和人工智能蓬勃发展的今天是怎样接轨的。正在苦恼从何下手时忽然想到研究生的一门课程——智能控制，于是便重新阅读了起来，开篇第一个智能控制方法——专家PID让我既熟悉又陌生，故写下这篇长文记录学习过程。
提示：本文对一些复杂的流程与概念进行了简化与省略，如果想详细了解的小伙伴可以参考刘金琨老师的《智能控制 -第4版》的第二章（需要电子版请私戳）
一、专家PID控制是什么
1.专家系统
专家PID是专家控制的一种方式，而专家控制是将人工智能中专家系统引入控制得到的，故我们需要简单了解一下专家系统。
定义：专家系统是一类包含知识和推理的智能计算机程序，其内部包含某领域专家水平的知识和经验，具有解决专门问题的能力。
专家系统主要由知识库和推理机组成：
2.专家控制
通过以上的内容，我们可以简单的了解到，专家系统就是把专家的推理过程交给计算机，当我们有问题时，专家系统就会基于知识库与推理机交给我们一个答案。
思考：专家控制与专家系统有什么关联呢？
专家系统也分可以归类为为许多不同领域的专用系统，比如医疗、电气、教育…如果专家系统中的知识库和推理机都是用于控制相关领域，帮助我们解决控制上的问题，就和我们说的专家控制很类似了。
定义:将专家系统的李虎和技术控制论、方法与技术相结合，在未知环境下，效仿专家的经验，实现对系统的控制。
其基本结构如图：
或许有小伙伴会有疑惑，那专家控制和专家系统有区别吗？
答案是肯定的：
①专家系统能完成专门领域的功能，辅助用户决策；但专家控制能进行独立的、实时的自动决策；
②专家控制对比专家系统有更高的可靠性和抗干扰要求；
③专家系统处于离线的工作方式，而专家控制要求在线获取反馈信息，即要求在线工作方式。
逐条分析：
①专家系统类似于搜索引擎，输入问题获得答案，专家控制其实就是控制器或者组成控制器的一部分，需要根据输入进行输出控制系统；
②专家控制用于对系统的控制，对于运动控制或过程控制来说，控制器的可靠与抗干扰是必要的；
③专家系统可以将数据库直接部署在本地，随后所有的决策都可以离线处理，而专家控制这里的在线可以理解为与系统接入，实时获取数据从而进行控制。
2.1 直接型专家控制
直接型专家控制器是作为控制器本身在系统中存在的。
特点：
①取代传统控制器，直接控制被控对象；
②具有模拟（或延伸、扩展）操作工人智能的功能；
③任务和功能相对比较简单；
④在线、实时控制；
⑤通常由几十条规则构成，便于增删与修改。
2.2 间接型专家控制
间接型专家控制器是作为控制器的一部分在系统中存在的。
特点：
①与传统控制器结合，间接控制被控对象；
②具有模拟（或延伸、扩展）控制工程师智能的功能；
③实现优化适应、协调、组织等高层决策；
④在线或离线控制。
基于不同的控制领域专家间接专家控制器可分为：优化型专家控制器、适应型专家控制器、协调型专家控制器和组织型专家控制器。
3.专家PID控制
专家PID控制是一种直接型专家 控制器，如果直白一点就是把图3中的控制器替换为专家PID控制器。
专业点的定义是：基于受控对象和控制规律的各种知识，无须知道被控对象的精确模型，利用专家经验来设计PID参数。
3.1 离散化的PID输出
在详细介绍专家PID控制的算法之前需要对经典的具有PID控制器的系统结构图：
提示：如果有小伙伴对系统结构图比较陌生建议重新回顾《自动控制原理（第六版）》——胡松涛（需要电子版请私戳）
由图5可知PID的传递函数满足：
U
(
s
)
E
(
s
)
=
K
p
+
K
i
s
+
K
d
s
(1)
\frac{U \left( s \right)}{E \left( s \right)}=K_p+\frac{K_i}{s}+K_ds \tag{1}
E(s)U(s)=Kp+sKi+Kds(1)其中
U
(
s
)
U \left( s \right)
U(s)和
E
(
s
)
E \left( s \right)
E(s)分别是PID控制器的输出和系统误差，
K
p
K_p
Kp、
K
i
K_i
Ki与
K
d
K_d
Kd分别是比例系数、积分系数和微分系数。
根据公式(1)写出PID控制器的输出频域公式为：
U
(
s
)
=
K
p
E
(
s
)
+
K
i
E
(
s
)
s
+
K
d
s
E
(
s
)
(2)
U \left( s \right)={K_pE \left( s \right)}+{K_i\frac{E \left( s \right)}{s}}+{K_dsE \left( s \right)} \tag{2}
U(s)=KpE(s)+KisE(s)+KdsE(s)(2)根据公式(2)得到PID控制器的输出时域公式为：
u
(
t
)
=
K
p
e
(
t
)
+
K
i
∫
e
(
t
)
d
t
+
K
d
d
e
(
t
)
d
t
(3)
u \left( t \right)={K_pe \left( t \right)}+{K_i{\int_{}^{} e \left( t \right)\, dt}}+{K_d \frac{de \left( t \right)}{dt}} \tag{3}
u(t)=Kpe(t)+Ki∫e(t)dt+Kddtde(t)(3)频域到时域是通过拉普拉斯逆变换得到的：
u
(
t
)
=
L
−
1
[
U
(
s
)
]
=
K
p
L
−
1
[
E
(
s
)
]
+
K
i
L
−
1
[
E
(
s
)
s
]
+
L
−
1
[
s
E
(
s
)
]
=
K
p
e
(
t
)
+
K
i
∫
e
(
t
)
d
t
+
K
d
d
e
(
t
)
d
t
(4)
\begin{align} u \left( t \right) &= \mathscr{L}^{-1} \left[ U \left( s \right) \right] \\ &=K_p \mathscr{L}^{-1} \left[ E \left( s \right) \right]+K_i \mathscr{L}^{-1} \left[ \frac{E \left( s \right)}{s} \right]+\mathscr{L}^{-1} \left[ sE \left( s \right) \right] \\ &={K_pe \left( t \right)}+{K_i{\int_{}^{} e \left( t \right)\, dt}}+{K_d \frac{de \left( t \right)}{dt}} \\ \end{align} \tag{4}
u(t)=L−1[U(s)]=KpL−1[E(s)]+KiL−1[sE(s)]+L−1[sE(s)]=Kpe(t)+Ki∫e(t)dt+Kddtde(t)(4)其中运用到了拉普拉斯逆变换的性质：
f
(
t
)
=
L
−
1
[
F
(
s
)
]
∫
f
(
t
)
d
t
=
L
−
1
[
F
(
s
)
s
]
d
f
(
t
)
d
t
=
L
−
1
[
s
F
(
s
)
]
(5)
\begin{align} f \left( t \right) &= \mathscr{L}^{-1} \left[ F \left( s \right) \right] \\ \int_{}^{} f \left( t \right)\, dt &= \mathscr{L}^{-1} \left[ \frac{F \left( s \right)}{s} \right] \\ \frac{df \left( t \right)}{dt} &= \mathscr{L}^{-1} \left[ s F \left( s \right) \right] \\ \end{align} \tag{5}
f(t)∫f(t)dtdtdf(t)=L−1[F(s)]=L−1[sF(s)]=L−1[sF(s)](5)基于专家的PID控制需要将系统进行离散化，基于时域的公式可以帮助我们更好的理解离散过程，这也是博主在这里介绍上述时域化的原因。
对系统的离散化，需要我们对微分和积分的原理有一定的理解，博主记得在高中，老师为了帮助那时无知的我们理解抽象的微分和积分概念，将函数曲线无限细分为了小单元：
微分其实就是求导，求解函数当前点的斜率，以图6最小单元为例：
d
f
(
t
)
d
t
=
f
i
−
f
i
−
1
t
i
−
t
i
−
1
(6)
\frac{df \left( t \right)}{dt}=\frac{f_i-f_{i-1}}{t_i-t_{i-1}} \tag{6}
dtdf(t)=ti−ti−1fi−fi−1(6)其中
f
i
f_i
fi、
f
i
−
1
f_{i-1}
fi−1表示当前时刻函数的值、上一个时刻函数的值，
t
i
t_i
ti、
t
i
−
1
t_{i-1}
ti−1表示当前时刻、上一个时刻，如果引入采样时间的概念则公式(6)可以改写为：
d
f
(
t
)
d
t
=
f
i
−
f
i
−
1
T
s
(7)
\frac{df \left( t \right)}{dt}=\frac{f_i-f_{i-1}}{T_s} \tag{7}
dtdf(t)=Tsfi−fi−1(7)其中
T
s
T_s
Ts表示采样时间。
积分其实就是求面积，求解函数与坐标轴形成的封闭图形的面积，以图6最小单元为例，我们认为阴影部分的面积与实际面积相比可以忽略不计，则积分公式：
∫
t
i
t
i
−
1
f
(
t
)
d
t
=
T
s
f
i
(8)
\int_{t_i}^{t_{i-1}} f \left( t \right)\, dt=Tsf_i \tag{8}
∫titi−1f(t)dt=Tsfi(8)
基于上述积分与微分的离散化处理，基于时域的公式(3)的离散化结果表示为：
u
k
=
K
p
e
k
+
K
i
T
s
∑
n
=
1
k
e
n
+
K
d
e
k
−
e
k
−
1
T
s
(9)
u_k={K_pe_k}+{K_i T_s \sum_{n=1}^k e_n}+{K_d \frac{e_k-e_{k-1}}{T_s}} \tag{9}
uk=Kpek+KiTsn=1∑ken+KdTsek−ek−1(9)其中
k
k
k表示第
k
k
k次采样。公式(9)就是我们熟知的位置式离散PID输出公式，注意公式(3)也是位置式PID的时域输出。对熟悉PID控制的小伙伴来说应该对位置式PID并不陌生，与之对应的是增量式PID：
Δ
u
(
t
)
=
K
p
Δ
e
(
t
)
+
K
i
∫
Δ
e
(
t
)
d
t
+
K
d
d
Δ
e
(
t
)
d
t
(10)
\Delta u \left( t \right)={K_p \Delta e \left( t \right)}+{K_i{\int_{}^{} \Delta e \left( t \right)\, dt}}+{K_d \frac{d \Delta e \left( t \right)}{dt}} \tag{10}
Δu(t)=KpΔe(t)+Ki∫Δe(t)dt+KddtdΔe(t)(10)将上式进行离散化：
u
k
−
u
k
−
1
=
K
p
(
e
k
−
e
k
−
1
)
+
K
i
T
s
∑
n
=
1
k
(
e
n
−
e
n
−
1
)
+
K
d
(
e
k
−
e
k
−
1
)
−
(
e
k
−
1
−
e
k
−
2
)
T
s
=
K
p
(
e
k
−
e
k
−
1
)
+
K
i
T
s
e
k
+
K
d
e
k
−
2
e
k
−
1
+
e
k
−
2
T
s
(11)
\begin{align} u_k-u_{k-1} &={K_p \left( {e_k-e_{k-1}} \right) }+{K_i T_s \sum_{n=1}^k \left( {e_n-e_{n-1}} \right)}+{K_d \frac{\left( {e_k-e_{k-1}} \right)-\left( {e_{k-1}-e_{k-2}} \right)}{T_s}} \\ &= {K_p \left( {e_k-e_{k-1}} \right) }+{K_i T_s e_k}+{K_d \frac{e_k-2e_{k-1}+e_{k-2}}{T_s}} \end{align} \tag{11}
uk−uk−1=Kp(ek−ek−1)+KiTsn=1∑k(en−en−1)+KdTs(ek−ek−1)−(ek−1−ek−2)=Kp(ek−ek−1)+KiTsek+KdTsek−2ek−1+ek−2(11)进而可得增量式PID的离散输出：
u
k
=
u
k
−
1
+
K
p
(
e
k
−
e
k
−
1
)
+
K
i
T
s
e
k
+
K
d
e
k
−
2
e
k
−
1
+
e
k
−
2
T
s
(12)
u_k=u_{k-1}+{K_p \left( {e_k-e_{k-1}} \right) }+{K_i T_s e_k}+{K_d \frac{e_k-2e_{k-1}+e_{k-2}}{T_s}} \tag{12}
uk=uk−1+Kp(ek−ek−1)+KiTsek+KdTsek−2ek−1+ek−2(12)
根据图5和公式(1)不难理解PID控制器根据系统误差控制被控对象，主要目的是为了让误差尽可能的小甚至等于零。以典型的二阶系统单位阶跃响应误差曲线为例：
其中，误差减小的曲线部分为透明表示的Ⅰ、Ⅲ、Ⅴ、Ⅶ，误差增大的曲线部分为阴影表示的Ⅱ、Ⅳ、Ⅵ、Ⅷ。首先为了解释的方便，定义：
Δ
e
k
=
e
k
−
e
k
−
1
Δ
e
k
−
1
=
e
k
−
1
−
e
k
−
2
(12)
\begin{align} \Delta e_k &= e_k-e_{k-1} \\ \Delta e_{k-1}&= e_{k-1}-e_{k-2} \end{align} \tag{12}
ΔekΔek−1=ek−ek−1=ek−1−ek−2(12)根据控制对误差减小的需求专家经验制定了以下规则：
①当
∣
e
(
k
)
∣
>
M
1
|e \left( k \right)|>M_1
∣e(k)∣>M1时，说明当前采样时刻误差值很大。此时不论误差有何变化趋势，都直接令PID控制器的输出为定值，以迅速调整误差，使误差以最大速度减小，同时也要注意避免超调。此时的控制效果相当于开环控制。
u
k
=
C
(13)
u_k=C \tag{13}
uk=C(13)②当
e
k
Δ
e
k
>
0
{e_k}{\Delta e_k}>0
ekΔek>0或
Δ
e
(
k
)
=
0
\Delta e \left( k \right)=0
Δe(k)=0时，说明误差在朝误差增大的方向变化或误差为一个定值未发生变化
a.
∣
e
(
k
)
∣
≥
M
2
|e \left( k \right)| \geq M_2
∣e(k)∣≥M2说明误差较大，需要施加较强的控制作用，使误差绝对值向减小的方向变化，从而迅速减小误差绝对值，控制器输出为：
u
k
=
u
k
−
1
+
k
1
{
k
p
(
e
k
−
e
k
−
1
)
+
k
i
T
s
e
k
+
k
d
(
e
k
−
2
e
k
−
1
+
e
k
−
2
)
T
s
}
(14)
u_k=u_{k-1}+k_1\{ {k_p \left( e_k-e_{k-1} \right)}+{k_i T_s e_k}+{k_d \frac{\left( e_k-2e_{k-1}+e_{k-2} \right)}{T_s}} \} \tag{14}
uk=uk−1+k1{kp(ek−ek−1)+kiTsek+kdTs(ek−2ek−1+ek−2)}(14)b.
∣
e
(
k
)
∣
<
M
2
|e \left( k \right)| < M_2
∣e(k)∣<M2说明尽管误差朝绝对值增大的方向变化，但误差绝对值本身并不是很大，考虑施加一般的控制作用，扭转误差的变化趋势，使其朝着误差绝对值减小的方向变化，控制器的输出：
u
k
=
u
k
−
1
+
k
p
(
e
k
−
e
k
−
1
)
+
k
i
T
s
e
k
+
k
d
(
e
k
−
2
e
k
−
1
+
e
k
−
2
)
T
s
(15)
u_k=u_{k-1}+ {k_p \left( e_k-e_{k-1} \right)}+{k_i T_s e_k}+{k_d \frac{\left( e_k-2e_{k-1}+e_{k-2} \right)}{T_s}} \tag{15}
uk=uk−1+kp(ek−ek−1)+kiTsek+kdTs(ek−2ek−1+ek−2)(15)③当
e
k
Δ
e
k
<
0
{e_k}{\Delta e_k}<0
ekΔek<0且
Δ
e
k
Δ
e
k
−
1
>
0
{\Delta e_k}{\Delta e_{k-1}}>0
ΔekΔek−1>0，或者
e
k
=
0
e_k=0
ek=0时，说明误差的绝对值朝减小的方向变化，或已经达到了平衡状态。此时考虑保持控制器输出不变：
u
k
=
u
k
−
1
(16)
u_k=u_{k-1} \tag{16}
uk=uk−1(16)④当
e
k
Δ
e
k
<
0
{e_k}{\Delta e_k}<0
ekΔek<0且
Δ
e
k
Δ
e
k
−
1
<
0
{\Delta e_k}{\Delta e_{k-1}}<0
ΔekΔek−1<0时，说明误差处于极值状态
a.
∣
e
(
k
)
∣
≥
M
2
|e \left( k \right)| \geq M_2
∣e(k)∣≥M2误差朝绝对值差较大，需要施加较强的控制作用，控制器输出为：
u
k
=
u
k
−
1
+
k
1
k
p
e
m
k
(17)
u_k=u_{k-1}+k_1k_pe_{mk} \tag{17}
uk=uk−1+k1kpemk(17)b.
∣
e
(
k
)
∣
<
M
2
|e \left( k \right)| < M_2
∣e(k)∣<M2误差朝绝对值差较小，需要施加较弱的控制作用，控制器输出为：
u
k
=
u
k
−
1
+
k
2
k
p
e
m
k
(18)
u_k=u_{k-1}+k_2k_pe_{mk} \tag{18}
uk=uk−1+k2kpemk(18)⑤当
∣
e
k
∣
<
ϵ
|e_k|<\epsilon
∣ek∣<ϵ时，说明误差的绝对值很小，此时加入积分环节，减小稳态无差，控制器输出为：
u
k
=
k
3
k
p
e
k
+
k
4
k
i
T
s
e
k
(19)
u_k=k_3k_pe_k+k_4k_iT_se_k \tag{19}
uk=k3kpek+k4kiTsek(19)以上参数中
k
1
k_1
k1是增益放大系数要求
k
1
>
1
k_1>1
k1>1，
k
2
k_2
k2是增益抑制系数要求
0
<
k
1
<
1
0<k_1<1
0<k1<1，
M
1
M_1
M1和
M
2
M_2
M2是设定的误差界限要求
M
1
>
M
2
>
0
M_1>M_2>0
M1>M2>0，
ϵ
\epsilon
ϵ是误差进精度要求
ϵ
<
<
1
\epsilon<<1
ϵ<<1
注意：上述的规则并非是定死的适用于所有系统，对于不同的系统使用者也可以根据自己的经验在此规则基础上进行改进。
二、专家PID控制的代码实现
1.二阶传递函数
已知传递函数
G
p
(
s
)
=
133
s
2
+
25
s
(20)
G_p \left( s \right)=\frac{133}{s^2+25s} \tag{20}
Gp(s)=s2+25s133(20)以采样时间1ms进行离散化，设计专家PID控制器，并进行Matlab仿真，基于位置式PID的程序如下：
% Expert PID Control
% 二阶传递函数的阶跃响应
% 位置式
clc;
clear all;
close all;
n=500; % 设置离散点的个数
Ts=0.001; % 设置离散的采样时间
epsilon=0.001; % 设置加入积分换进的精度ε
% 对连续时间的传递函数进行离散化
% 两种方式可以同样设置传递函数（连续时间）
% s=tf('s');
% Gs=133/(s^2+25*s);
Gs=tf(133,[1,25,0]);
Gz=c2d(Gs,Ts,'z');
[num,den]=tfdata(Gz,'v');
time=zeros(1,n ); % 离散时间
r=ones(1,n ); % 离散输入R(s) r(t)
u1=zeros(1,n ); % 离散专家PID控制器输出U(s) u(t)
u2=zeros(1,n ); % 离散普通PID控制器输出U(s) u(t)
y1=zeros(1,n ); % 专家PID控制下的系统输出
y2=zeros(1,n ); % 普通PID控制下的系统输出
error1=zeros(1,n); % 专家PID控制下的系统误差
error2=zeros(1,n); % 普通PID控制下的系统误差
% 根据z变换后的公式可知需要
u1_1=0;u1_2=0; % u(k-1) u(k-2)
y1_1=0;y1_2=0; % y(k-1) y(k-2)
u2_1=0;u2_2=0; % u(k-1) u(k-2)
y2_1=0;y2_2=0; % y(k-1) y(k-2)
x1=[0,0,0]'; % 位置式PID：x(1)=error(k);x(2)=(error(k)-error(k-1))/Ts;x(3)=Ts*Σen (n=1~k)
x2=[0,0,0]';
x1_2_1=0;
x2_2_1=0;
kp1=2.0;ki1=0.20;kd1=1.00; % 专家PID原有PID参数
% kp1=3.5;ki1=0.05;kd1=0.25; % 专家PID原有PID参数
kp2=2.0;ki2=0.20;kd2=1.00; % 专家PID原有PID参数
% kp2=3.5;ki2=0.05;kd2=0.25; % 普通PID原有PID参数
% kp2=3.0;ki2=0.05;kd2=0.1; % 普通PID原有PID参数
error1_1=0; % 误差的第k个极值
error2_1=0; % 误差的第k个极值
% k1=5.0; % 增益放大系数
k1=2.0;
k2=0.6; % 抑制系数
% 循环更新数值(专家PID)
for k=1:1:n
time(k)=k*Ts; %离散的时间
u1(k)=kp1*x1(1)+kd1*x1(2)+ki1*x1(3); %PID控制器
%专家控制的规则
%规则1
%当误差的
if abs(x1(1))>0.80
u1(k)=(3/4)*kp1;
elseif abs(x1(1))>0.40
u1(k)=(2/3)*kp1;
elseif abs(x1(1))>0.20
u1(k)=(1/5)*kp1;
elseif abs(x1(1))>0.01
u1(k)=(1/10)*kp1;
end
%规则2
if (x1(1)*x1(2)>0)||(x1(2)==0)
if abs(x1(1))>=0.05
u1(k)=k1*(kp1*x1(1)+kd1*x1(2)+ki1*x1(3));
else
u1(k)=k2*(kp1*x1(1)+kd1*x1(2)+ki1*x1(3));
end
end
%规则3
if ((x1(1)*x1(2)<0)&&(x1(2)*x1_2_1>0))||(x1(1)==0)
u1(k)=u1(k);
end
%规则4
if (x1(1)*x1(2)<0)&&(x1(2)*x1_2_1<0)
if abs(x1(1))>=0.05
u1(k)=k1*kp1*error1_1;
else
u1(k)=k2*kp1*error1_1;
end
end
%规则5
if abs(x1(1))<=epsilon
u1(k)=((5/6)*kp1)*x1(1)+((1/3)*ki1)*x1(3);
% u1(k)=0.5*x1(1)+0.01*x1(3);
end
% 限制控制器的范围
if u1(k)>=10
u1(k)=10;
end
if u1(k)<=-10
u1(k)=-10;
end
%线性化的传递函数
y1(k)=-den(2)*y1_1-den(3)*y1_2+num(1)*u1(k)+num(2)*u1_1+num(3)*u1_2;
error1(k)=r(k)-y1(k);
%更新参数
u1_2=u1_1;u1_1=u1(k);
y1_2=y1_1;y1_1=y1(k);
x1(1)=error1(k);
x1_2_1=x1(2);
x1(2)=(error1(k)-error1_1)/Ts;
x1(3)=x1(3)+error1(k)*Ts;
error1_1=error1(k);
end
% 循环更新数值(普通PID)
for k=1:1:n
time(k)=k*Ts; %离散的时间
u2(k)=kp2*x2(1)+kd2*x2(2)+ki2*x2(3); %PID控制器
%线性化的传递函数
y2(k)=-den(2)*y2_1-den(3)*y2_2+num(1)*u2(k)+num(2)*u2_1+num(3)*u2_2;
error2(k)=r(k)-y2(k);
%更新参数
u2_2=u2_1;u2_1=u2(k);
y2_2=y2_1;y2_1=y2(k);
x2(1)=error2(k);
x2_2_1=x2(2);
x2(2)=(error2(k)-error2_1)/Ts;
x2(3)=x2(3)+error2(k)*Ts;
error2_1=error2(k);
end
figure(1);
plot(time,r,'b',time,y1,'r');
xlabel('time(s)');
ylabel('r,y');
legend('输入','专家PID输出');
figure(2);
plot(time,r-y1,'r');
xlabel('time(s)');
ylabel('error');
legend('专家PID误差');
figure(3);
plot(time,r,'b',time,y1,'r',time,y2,'k');
xlabel('time(s)');
ylabel('r,y1,y2');
legend('输入','专家PID输出','普通PID输出');
经过运行后的结果：
基于增量式PID的代码如下：
% Expert PID Control
% 二阶传递函数的阶跃响应
% 增量式
clc;
clear;
close all;
N=600; % 设置离散点的个数
Ts=0.001; % 设置离散的采样时间
epsilon=0.001; % 设置加入积分换进的精度ε
% 对连续时间的传递函数进行离散化
% 两种方式可以同样设置传递函数（连续时间）
% s=tf('s');
% Gs=133/(s^2+25*s);
sys_c=tf(133,[1,25,0]);
sys_d=c2d(sys_c,Ts,'z');
[num_d,den_d]=tfdata(sys_d,'v');
time=zeros(1,N ); % 离散时间
r=ones(1,N ); % 离散输入R(s) r(t)
u1=zeros(1,N ); % 离散专家PID控制器输出U(s) u(t)
u2=zeros(1,N ); % 离散普通PID控制器输出U(s) u(t)
c1=zeros(1,N ); % 专家PID控制下的系统输出
c2=zeros(1,N ); % 普通PID控制下的系统输出
error1=zeros(1,N); % 专家PID控制下的系统误差
error2=zeros(1,N); % 普通PID控制下的系统误差
u1_1=0;u1_2=0; % y(k-1) y(k-2)
c1_1=0;c1_2=0; % u(k-1) u(k-2)
u2_1=0;u2_2=0; % y(k-1) y(k-2)
c2_1=0;c2_2=0; % u(k-1) u(k-2)
x1=[0,0,0]'; % 增量式PID：x(1)=error(k)-error(k-1);x(2)=(error(k)-2*error(k-1)+error(k-2))/Ts;x(3)=error(k)
x2=[0,0,0]'; %
x1_2_1=0;
x2_2_1=0;
kp1=221.9756;ki1=15855;kd1=0.808; % 专家PID原有PID参数
kp2=221.9756;ki2=15855;kd2=0.808; % 普通PID原有PID参数
error1_1=0; % error1(k-1)
error1_2=0; % error1(k-2)
error1_3=0; % error1(k-3)
error2_1=0; % error2(k-1)
error2_2=0; % error2(k-2)
error2_3=0; % error2(k-3)
sign=0;
k1=1.2; % 规则2增益放大系数
k2=0.6; % 规则2增益抑制系数
k3=0.01; % 规则4增益抑制系数
% 循环更新数值（专家PID）
for n=1:1:N
time(n)=n*Ts; % 离散的时间
u1(n)=u1_1+kp1*x1(1)+kd1*x1(2)+ki1*x1(3); % PID控制器
%专家控制的规则
%规则1
%
if(abs(error1_1)>0.1)
if abs(error1_1)>0.80
u1(n)=10;
elseif abs(error1_1)>0.40
u1(n)=5;
elseif abs(error1_1)>0.20
u1(n)=1;
elseif abs(error1_1)>0.10
u1(n)=0.1;
end
else
%规则2
if (error1_1*(error1_1-error1_2)>0)||(error1_1-error1_2==0)
if abs(error1_1)>=0.05
u1(n)=u1_1+k1*(kp1*x1(1)+kd1*x1(2)+ki1*x1(3));
else
u1(n)=u1_1+k2*(kp1*x1(1)+kd1*x1(2)+ki1*x1(3));
end
end
%规则3
if ((error1_1*(error1_1-error1_2)<0)&&((error1_1-error1_2)*(error1_2-error1_3)>0))||(error1_1==0)
u1(n)=u1(n);
end
%规则4
if (error1_1*(error1_1-error1_2)<0)&&((error1_1-error1_2)*(error1_2-error1_3)<0)
if abs(error1_1)>=0.05
u1(n)=u1_1+k1*kp1*error1_1;
else
u1(n)=u1_1+k3*kp1*error1_1;
end
end
%规则5
if abs(error1_1)<=epsilon
u1(n)=500*x1(1)+5000*x1(3);
end
end
% 限制控制器的范围
if u1(n)>=10
u1(n)=10;
end
if u1(n)<=-10
u1(n)=-10;
end
%线性化的传递函数
c1(n)=-den_d(2)*c1_1-den_d(3)*c1_2+num_d(2)*u1_1+num_d(3)*u1_2;
error1(n)=r(n)-c1(n);
%更新参数
u1_2=u1_1;u1_1=u1(n);
c1_2=c1_1;c1_1=c1(n);
if n>1
x1(1)=error1(n)-error1(n-1);
else
x1(1)=error1(n);
end
if n>2
x1(2)=(error1(n)-2*error1(n-1)+error1(n-2))/Ts;
elseif n>1
x1(2)=(error1(n)-2*error1(n-1))/Ts;
else
x1(2)=error1(n)/Ts;
end
x1(3)=Ts*error1(n);
error1_1=error1(n);
if n>1
error1_2=error1(n-1);
end
if n>2
error1_3=error1(n-2);
end
end
% 循环更新数值（普通PID）
for n=1:1:N
time(n)=n*Ts; %离散的时间
u2(n)=u2_1+kp2*x2(1)+kd2*x2(2)+ki2*x2(3); %PID控制器
%线性化的传递函数
c2(n)=-den_d(2)*c2_1-den_d(3)*c2_2+num_d(2)*u2_1+num_d(3)*u2_2;
error2(n)=r(n)-c2(n);
%更新参数
u2_2=u2_1;u2_1=u2(n);
c2_2=c2_1;c2_1=c2(n);
if n>1
x2(1)=error2(n)-error2(n-1);
else
x2(1)=error2(n);
end
x2_2_1=x2(2);
if n>2
x2(2)=(error2(n)-2*error2(n-1)+error2(n-2))/Ts;
elseif n>1
x2(2)=(error2(n)-2*error2(n-1))/Ts;
else
x2(2)=error2(n)/Ts;
end
x2(3)=Ts*error2(n);
end
figure(1);
plot(time,r,'b',time,c1,'r');
xlabel('time(s)');
ylabel('r,c');
legend('输入','专家PID输出');
figure(2);
plot(time,r-c1,'r');
xlabel('time(s)');
ylabel('error');
legend('专家PID误差');
figure(3);
plot(time,r,'b',time,c1,'r',time,c2,'k');
legend('输入','专家PID输出','普通PID输出');
xlabel('time(s)');
ylabel('r,c1,c2');
经过运行后的结果：
有小伙伴在阅读代码的时候可能也发现了，博主在写代码的时候对于选择PID原有参数经过了一些测试与实验，主要使用的是临界比例度法，可以参考这篇文章先进PID控制——基于临界比例度法的PID整定。
还需要注意的是，专家经验中有许多参数如 M 1 M_1 M1、 M 2 M_2 M2、 k 1 k_1 k1、 k 2 k_2 k2和 ϵ \epsilon ϵ都是待确定的参数，而博主在这里设定的参数也是自己定的，具体的数值是怎样得到的还没有一个可以形成文字的公式或规则，如果有大佬也可以分享自己的调参经验。
2.三阶传递函数
有三阶传递：
G
p
(
s
)
=
523500
s
3
+
87.35
s
2
+
10470
s
G_p \left( s \right)=\frac{523500}{s^3+87.35s^2+10470s}
Gp(s)=s3+87.35s2+10470s523500同样以采样时间1ms进行离散化，设计专家PID控制器，并进行Matlab仿真，基于位置式PID的程序如下：
% Expert PID Control
% 三阶传递函数的阶跃响应
% 位置式PID
clc;
clear all;
close all;
n=500; % 设置离散点的个数
Ts=0.001; % 设置离散的采样时间
epsilon=0.001; % 设置加入积分换进的精度ε
% 对连续时间的传递函数进行离散化
% 两种方式可以同样设置传递函数（连续时间）
% s=tf('s');
% Gs=5.235e5/(s^3+87.35*s^2+10470*s);
Gs=tf(5.235e5,[1,87.35,1.047e4,0]);
Gz=c2d(Gs,Ts,'z');
[num,den]=tfdata(Gz,'v');
time=zeros(1,n ); % 离散时间
r=ones(1,n ); % 离散输入R(s) r(t)
u1=zeros(1,n ); % 离散专家PID控制器输出U(s) u(t)
u2=zeros(1,n ); % 离散普通PID控制器输出U(s) u(t)
y1=zeros(1,n ); % 专家PID控制下的系统输出
y2=zeros(1,n ); % 普通PID控制下的系统输出
error1=zeros(1,n); % 专家PID控制下的系统误差
error2=zeros(1,n); % 普通PID控制下的系统误差
% 根据z变换后的公式可知需要
u1_1=0;u1_2=0;u1_3=0; % u(k-1) u(k-2) u(k-3)
y1_1=0;y1_2=0;y1_3=0; % y(k-1) y(k-2) y(k-3)
u2_1=0;u2_2=0;u2_3=0; % u(k-1) u(k-2) u(k-3)
y2_1=0;y2_2=0;y2_3=0; % y(k-1) y(k-2) y(k-3)
%
x1=[0,0,0]'; % 位置式PID：x(1)=error(k);x(2)=(error(k)-error(k-1))/Ts;x(3)=Σen (n=1~k)
x2=[0,0,0]';
x1_2_1=0;
x2_2_1=0;
kp1=0.6;ki1=0.03;kd1=0.01; % 专家PID原有PID参数
kp2=0.6;ki2=0.03;kd2=0.01; % 普通PID原有PID参数
error1_1=0; % 误差的第k个极值
error2_1=0; % 误差的第k个极值
k1=2; % 增益放大系数
k2=0.6; % 抑制系数
% 循环更新数值(专家PID)
for k=1:1:n
time(k)=k*Ts; %离散的时间
u1(k)=kp1*x1(1)+kd1*x1(2)+ki1*x1(3); %PID控制器
%专家控制的规则
%规则1
%当误差的
if abs(x1(1))>0.8
u1(k)=0.45;
elseif abs(x1(1))>0.40
u1(k)=0.40;
elseif abs(x1(1))>0.20
u1(k)=0.12;
elseif abs(x1(1))>0.01
u1(k)=0.10;
end
%规则2
if (x1(1)*x1(2)>0)||(x1(2)==0)
if abs(x1(1))>=0.05
u1(k)=u1_1+k1*kp1*x1(1);
else
u1(k)=u1_1+k2*kp1*x1(1);
end
end
%规则3
if ((x1(1)*x1(2)<0)&&(x1(2)*x1_2_1>0))||(x1(1)==0)
u1(k)=u1(k);
end
%规则4
if (x1(1)*x1(2)<0)&&(x1(2)*x1_2_1<0)
if abs(x1(1))>=0.05
u1(k)=u1_1+k1*kp1*error1_1;
else
u1(k)=u1_1+k2*kp1*error1_1;
end
end
%规则5
if abs(x1(1))<=epsilon
u1(k)=0.5*x1(1)+0.01*x1(3);
end
% 限制控制器的范围
if u1(k)>=10
u1(k)=10;
end
if u1(k)<=-10
u1(k)=-10;
end
%线性化的传递函数
y1(k)=-den(2)*y1_1-den(3)*y1_2-den(4)*y1_3+num(1)*u1(k)+num(2)*u1_1+num(3)*u1_2+num(4)*u1_3;
error1(k)=r(k)-y1(k);
%更新参数
u1_3=u1_2;u1_2=u1_1;u1_1=u1(k);
y1_3=y1_2;y1_2=y1_1;y1_1=y1(k);
x1(1)=error1(k);
x1_2_1=x1(2);
x1(2)=(error1(k)-error1_1)/Ts;
x1(3)=x1(3)+error1(k)*Ts;
error1_1=error1(k);
end
for k=1:1:n
time(k)=k*Ts; %离散的时间
u2(k)=kp2*x2(1)+kd2*x2(2)+ki2*x2(3); %PID控制器
%线性化的传递函数
y2(k)=-den(2)*y2_1-den(3)*y2_2-den(4)*y2_3+num(1)*u2(k)+num(2)*u2_1+num(3)*u2_2+num(4)*u2_3;
error2(k)=r(k)-y2(k);
%更新参数
u2_3=u2_2;u2_2=u2_1;u2_1=u2(k);
y2_3=y2_2;y2_2=y2_1;y2_1=y2(k);
x2(1)=error2(k);
x2_2_1=x2(2);
x2(2)=(error2(k)-error2_1)/Ts;
x2(3)=x2(3)+error2(k)*Ts;
error2_1=error2(k);
end
figure(1);
plot(time,r,'b',time,y1,'r');
xlabel('time(s)');
ylabel('r,y');
legend('输入','专家PID输出');
figure(2);
plot(time,r-y1,'r');
xlabel('time(s)');
ylabel('error');
legend('专家PID误差');
figure(3);
plot(time,r,'b',time,y1,'r',time,y2,'k');
xlabel('time(s)');
ylabel('r,y1,y2');
legend('输入','专家PID输出','普通PID输出');
结果：
基于增量式PID的程序如下：
% Expert PID Control
% 三阶传递函数的阶跃响应
% 增量式PID
clc;
clear;
close all;
n=500; % 设置离散点的个数
Ts=0.001; % 设置离散的采样时间
epsilon=0.001; % 设置加入积分换进的精度ε
% 对连续时间的传递函数进行离散化
% 两种方式可以同样设置传递函数（连续时间）
% s=tf('s');
% Gs=5.235e5/(s^3+87.35*s^2+10470*s);
sys_c=tf(5.235e5,[1,87.35,1.047e4,0]);
sys_d=c2d(sys_c,Ts,'z');
[num_d,den_d]=tfdata(sys_d,'v');
time=zeros(1,n); % 离散时间
r=ones(1,n); % 离散输入R(s) r(t)
u1=zeros(1,n); % 离散专家PID控制器输出U(s) u(t)
u2=zeros(1,n); % 离散普通PID控制器输出U(s) u(t)
c1=zeros(1,n ); % 专家PID控制下的系统输出
c2=zeros(1,n ); % 普通PID控制下的系统输出
error1=zeros(1,n); % 专家PID控制下的系统误差
error2=zeros(1,n); % 普通PID控制下的系统误差
% 根据z变换后的公式可知需要
u1_1=0;u1_2=0;u1_3=0; % u(k-1) u(k-2) u(k-3)
c1_1=0;c1_2=0;c1_3=0; % c(k-1) c(k-2) c(k-3)
u2_1=0;u2_2=0;u2_3=0; % u(k-1) u(k-2) u(k-3)
c2_1=0;c2_2=0;c2_3=0; % c(k-1) c(k-2) c(k-3)
%
x1=[0,0,0]'; % 增量式PID：x(1)=error(k)-error(k-1);x(2)=(error(k)-2*error(k-1)+error(k-2))/Ts;x(3)=error(k)
x2=[0,0,0]';
x1_2_1=0; % Δe(k-1)
x2_2_1=0;
kp1=0.9804;ki1=26.4971;kd1=0.0094; % 专家PID参数
kp2=0.9804;ki2=26.4971;kd2=0.0094; % 普通PID参数
error1_1=0; % error1(k-1)
error1_2=0; % error1(k-2)
error1_3=0; % error1(k-3)
error2_1=0; % error2(k-1)
error2_2=0; % error2(k-2)
error2_3=0; % error2(k-3)
k1=2; % 增益放大系数
k2=0.4; % 抑制系数
% 循环更新数值（专家PID）
for n=1:1:n
time(n)=n*Ts; %离散的时间
u1(n)=u1_1+kp1*x1(1)+kd1*x1(2)+ki1*x1(3);%PID控制器
%专家控制的规则
%规则1
%
if abs(error1_1)>0.10
if abs(error1_1)>0.80
u1(n)=0.6;
elseif abs(error1_1)>0.40
u1(n)=0.4;
elseif abs(error1_1)>0.20
u1(n)=0.1;
elseif abs(error1_1)>0.10
u1(n)=0.05;
end
else
%规则2
if (error1_1*(error1_1-error1_2)>0)||(error1_1-error1_2==0)
if abs(x1(1))>=0.05
u1(n)=u1_1+k1*(kp1*x1(1)+kd1*x1(2)+ki1*x1(3));
else
u1(n)=u1_1+k2*(kp1*x1(1)+kd1*x1(2)+ki1*x1(3));
end
end
%规则3
if ((error1_1*(error1_1-error1_2)<0)&&((error1_1-error1_2)*(error1_2-error1_3)>0))||(error1_1==0)
u1(n)=u1(n);
end
%规则4
if (error1_1*(error1_1-error1_2)<0)&&((error1_1-error1_2)*(error1_2-error1_3)<0)
if abs(x1(1))>=0.05
u1(n)=u1_1+k1*kp1*error1_1;
else
u1(n)=u1_1+k2*kp1*error1_1;
end
end
%规则5
if abs(error1_1)<=epsilon
u1(n)=0.8*x1(1)+25*x1(3);
end
end
% 限制控制器的范围
if u1(n)>=10
u1(n)=10;
end
if u1(n)<=-10
u1(n)=-10;
end
%线性化的传递函数
c1(n)=-den_d(2)*c1_1-den_d(3)*c1_2-den_d(4)*c1_3+num_d(1)*u1(n)+num_d(2)*u1_1+num_d(3)*u1_2+num_d(4)*u1_3;
error1(n)=r(n)-c1(n);
%更新参数
u1_3=u1_2;u1_2=u1_1;u1_1=u1(n);
c1_3=c1_2;c1_2=c1_1;c1_1=c1(n);
if n>1
x1(1)=error1(n)-error1(n-1);
else
x1(1)=error1(n);
end
x1_2_1=x1(2);
if n>2
x1(2)=(error1(n)-2*error1(n-1)+error1(n-2))/Ts;
elseif n>1
x1(2)=(error1(n)-2*error1(n-1))/Ts;
else
x1(2)=error1(n)/Ts;
end
x1(3)=Ts*error1(n);
error1_1=error1(n);
if n>1
error1_2=error1(n-1);
end
if n>2
error1_3=error1(n-2);
end
end
% 循环更新数值（普通PID）
for n=1:1:n
time(n)=n*Ts; %离散的时间
u2(n)=u2_1+kp2*x2(1)+kd2*x2(2)+ki2*x2(3); %PID控制器
%线性化的传递函数
c2(n)=-den_d(2)*c2_1-den_d(3)*c2_2-den_d(4)*c2_3+num_d(1)*u2(n)+num_d(2)*u2_1+num_d(3)*u2_2+num_d(4)*u2_3;
error2(n)=r(n)-c2(n);
%更新参数
u2_3=u2_2;u2_2=u2_1;u2_1=u2(n);
c2_3=c2_2;c2_2=c2_1;c2_1=c2(n);
if n>1
x2(1)=error2(n)-error2(n-1);
else
x2(1)=error2(n);
end
x2_2_1=x2(2);
if n>2
x2(2)=(error2(n)-2*error2(n-1)+error2(n-2))/Ts;
elseif n>1
x2(2)=(error2(n)-2*error2(n-1))/Ts;
else
x2(2)=error2(n)/Ts;
end
x2(3)=Ts*error2(n);
error2_1=error2(n);
end
figure(1);
plot(time,r,'b',time,c1,'r');
xlabel('time(s)');
ylabel('r,c');
legend('输入','专家PID输出');
figure(2);
plot(time,r-c1,'r');
xlabel('time(s)');
ylabel('error');
legend('专家PID误差');
figure(3);
plot(time,r,'b',time,c1,'r',time,c2,'k');
xlabel('time(s)');
ylabel('r,c1,c2');
legend('输入','专家PID输出','普通PID输出');
figure(4);
plot(time,u1,'r',time,u2,'k');
xlabel('time(s)');
ylabel('u1,u2');
legend('专家PID控制器输出','普通PID控制器输出');
结果：
总结
博主在本文中着重介绍了专家PID控制的设计规则与原理， 也给出了基于二阶与三阶系统的增量式、位置式专家PID控制代码。不过上述的规则与代码也有很多可改进之处，博主受限于能力水平暂未找到一个更好的规则与实现方式，如果有大佬也可以在评论区指点一下，本文的错误也欢迎大家指正。
扫一扫
2029
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
