自适应控制.doc

10．自适应控制

严格地说，实际过程中的控制对象自身及能所处的环境都是十分复杂的，其参数会由于种种外部与内部的原因而发生变化。如，化学反应过程中的参数随环境温度和湿度的变化而变化（外部原因），化学反应速度随催化剂活性的衰减而变慢（内部原因），等等。如果实际控制对象客观存在着较强的不确定，那么，前面所述的一些基于确定性模型参数来设计控制系统的方法是不适用的。
所谓自适应控制是对于系统无法预知的变化，能自动地不断使系统保持所希望的状态。因此，一个自适应控制系统，应能在其运行过程中，通过不断地测取系统的输入、状态、输出或性能参数，逐渐地了解和掌握对象，然后根据所获得的过程信息，按一定的设计方法，作出控制决策去修正控制器的结构，参数或控制作用，以便在某种意义下，使控制效果达到最优或近似更优。目前比较成熟的自适应控制可分为两大类：模型参考自适应控制（Model Reference Adaptive Control）和自校正控制（Self-Turning）。
10.1模型参考自适应控制
10.1.1模型参考自适应控制原理
模型参考自适应控制系统的基本结构与图10.1所示：
10.1模型参考自适应控制系统
它由两个环路组成，由控制器和受控对象组成内环，这一部分称之为可调系统，由参考模型和自适应机构组成外环。实际上，该系统是在常规的反馈控制回路上再附加一个参考模型和控制器参数的自动调节回路而形成。
在该系统中，参考模型的输出或状态相当于给定一个动态性能指标，（通常，参考模型是一个响应比较好的模型），目标信号同时加在可调系统与参考模型上，通过比较受控对象与参考模型的输出或状态来得到两者之间的误差信息，按照一定的规律（自适应律）来修正控制器的参数（参数自适应）或产生一个辅助输入信号（信号综合自适应），从而使受控制对象的输出尽可能地跟随参考模型的输出。
在这个系统，当受控制对象由于外界或自身的原因系统的特性发生变化时，将导致受控对象输出与参考模型输出间误差的增大。于是，系统的自适应机构再次发生作用调整控制器的参数，使得受控对象的输出再一次趋近于参考模型的输出即与理想的希望输出相一致）。这就是参考模型自适应控制的基本工作原理。
模型参考自适应控制设计的核心问题是怎样决定和综合自适应律，有两类方法，一类为参数最优化方法，即利用优化方法寻找一组控制器的最优参数，使与系统有关的某个评价目标，如：J=e2(t)dt，达到最小。另一类方法是基于稳定性理论的方法，其基本思想是保证控制器参数自适应调节过程是稳定的。如基于Lyapunov稳定性理论的设计方法和基于Popov超稳定理论的方法。
系统设计举例
以下通过一个设计举例说明参数最优化设计方法的具体应用。
例10.1设一受控系统的开环传递函数为Wa(s)=，其中K可变，要求用一参考模型自适应控制使系统得到较好的输出。
解：对于该系统，我们选其控制器为PID控制器，而PID控制器的参数由自适应机构来调节，参考模型选性能综合指标良好的一个二阶系统：
Wm(d)=
自适应津决定的评价函数取
minJ=e2(t)dt ，e(t)为参考模型输出与对象输出的误差。
由于评价函数不能写成PID参数的解析函数形式，因此选用单纯形法做为寻优方法。（参见有关优化设计参考文献）。
在上述分析及考虑下，可将系统表示具体结构表示如下图10.2所示。

r(t) x0 x1 x2
1/s 1/s
- -
1.414
x7

x5 x3 x4 x6
PID 1/s

单纯形法寻优

图10.2 系统仿真结构图
为了进行数字仿真（连续系统离散相似法），图中的模型是采用典型环节（积分环节，惯性环节）型式表示。
对于图5-17所示系统，可写出各类环节的差数方程，其中，PID控制器的差分方程为：

式中，Kp为比例系数，TI为积分时间常数，Td为微分时间常数，它们是自适应津的调整参数，也是单纯形法的寻优参数。T是采样周期。
用MATLAB编写的仿真程序如下，（chp10-1 .m和suba.m）,程序中有关离散相似法仿真及单纯形法寻优的原理请参见有关参考书，如参考文献[ ]。
global R I9 L1 L2 L3 x p1 p2 p3
global Kp E1 F1 G1 E5 F5 G5 m
global K1 T4 T5
x=zeros(3,7);Q=zeros(1,7);p1=zeros(1000,1);k=zeros(3,1);
N=3;E=0.001;R=1;L=0.012;p2=zeros(1000,1);p3=zeros(1000,1);
L1=0.05;L2=5;L3=200;m=0;mn=0;
Kp1=zeros(100,1);Ti=zeros(100,1);Td=zeros(100,1);
Kp=2;
E1=1;F1=L1;G1=L1*L1/2;
E5=exp(-L1);F5=Kp*(1-E5);G5=Kp*(E5-1)+Kp*L1;
I9=0;
x(1,1)=1;x(2,1)=0.5;x(3,1)=0.9;
for j=2:N+1
for i=1:N
if i==j-1
x(i,j)=x(i,1)+L;
else
x(i,j)=x(i,1);
end
end
end
for j=1:N+1
Q(j)=suba(j);
end
Ex=0;
while Ex<1
Q0=0;
for j=1:N+1
if Q(j)>Q0
Q0=Q(j);
j0=j;
end
end
Q1=0;
for j=1:N+1
if j~=j0
if Q(j)>Q1
Q1=Q(j);
j1=j;
end
end
end
Q2=Q0;
for j=1:N+1
if Q(j)mn+1;
Kp1(mn,1)=x(1,j2);
Ti(mn,1)=x(2,j2);
Td(mn,1)=x(3,j2);
if Q2<1
if (Q0-Q2)end
x(i,N+2)=2*(x(i,N+2)-x(i,j0))/N-x(i,j0);
end
j=N+2;
Q(N+2)=suba(j);
if Q(N+2) for i=1:N
x(i,N+3)=x(i,N+4);
end
Q(N+3)=Q(N+4);
else
for i=1:N
x(i,N+3)=x(i,N+2);
end
Q(N+3)=Q(N+2)
end
else
for i=1:N
x(i,N+3)=x(i,N+2);
end
Q(N+3)=Q(N+2);
end
for i=1:N
x(i,j0)=x(i,N+3);
end
Q(j0)=Q(N+3);
else
if Q(N+2)>Q(j0)
for i=1:N
x(i,N+3)=x(i,N+2);
x(i,N+2)=x(i,j0);
x(i,j0)=x(i,N+3);
end
end
for i=1:N
x(i,N+3)=0.3*x(i,j0)+0.7*x(i,N+2);
end
j=N+3;
Q(N+3)=suba(j);
if Q(N+3)))/2;
end
end
for j=1:N+1
Q(j)=suba(j);
end
end
end
end
I9=1;
j=j2;
Fin=suba(j)
x(1,j2),x(2,j2),x(3,j2)
figure(1)
m1=1:1:m;m1=m1';
p1=p1(1:m,1);p2=p2(1:m,1);p3=p3(1:m,1);
plot(m1,p1,'r',m1,p2,'B',m1,p3,'g')
grid
pause
figure(2)

m2=1:1:mn;m2=m2';
Kp1=Kp1(1:mn,1);Ti=Ti(1:mn,1);Td=Td(1:mn,1);
plot(m2,Kp1,'r',m2,Ti,'B',m2,Td,'g')
grid

function resulta=suba(j)

global R I9 L1 L2 L3 x p1 p2 p3

global Kp E1 F1 G1 E5 F5 G5 m

global K1 T4 T5

K1=x(1,j);

T4=x(2,j);

T5=x(3,j);

x1=0;x2=0;x3=0;x4=0;x6=0;

w1=0;w2=0;w4=0;w5=0;u5=0;

t=0;F=0;

T6=L1*L2;

for i1=1:L3

x5=R-x6;

x3=x3+K1*(x5-w5+T6*x5/T4+T5*(x5-2*w5+u5)/T6);
u5=w5;w5=x5;
for i2=1:L2
x0=R-1.414*x1-x2;
u2=x1;u4=x4;
v1=(x0-w1)/L1;
x1=E1*x1+F1*x0+G1*v1;
w1=x0;
v2=(u2-w2)/L1;
x2=E1*x2+F1*u2+G1*v2;
w2=u2;
v3=(u4-w4)/L1;
x4=E1*x4+F1*x3;
x6=E5*x6+F5*u4+G5*v3;
w4=u4;
x7=x2-x6;
t=t+L1;
F=F+x7*x7*L1;
if I9==1
m=m+1;
p1(m,1)=x2;
p2(m,1)=x6;
p3(m,1)=x7;
end
end
end
resulta=F;

程序运行后，绘制出受控对象和参考模型的输出及两者间的误差曲线如图10.3，控制器的三个参数（Kp,TI ,Td）的变化曲线如图10.4所示，图中可以看到，经过一段时间的调整（kd,TI,Td自动调整）后，受控对象输出能很好地接近参考模型的输出。

图10.3 受控对象和参考模型的输出 图10.4 Kp,TI ,Td的变化曲线
及两者间的误差曲线
读者可以改变系统的某些参数（如让程序中的ka=1或ka=3等）后再运行程序观看仿真结果。
10.2自校正控制
自校正控制技术特别适宜用于结构已知而参数未知但恒定或缓慢变化的随机系统。由于大多数工业对象正好具有这种特征，因此，自校正控制技术在工业控制过程中得到了广泛的应用。

10.2.1自校正控制系统的基本结构

图10.5 自校正控制系统的典型结构
与模型参考自适应控制系统一样，自校正控制系统也由两个环路组成，典型结构如图10.5所示。内环与常规反馈系统类似，由对象和控制器组成。外环由参数估计器和控制器设计计算机构组成。参数估计器的功用是根据受控对象的输入及输出信息，连续不断地估计受控对象的参数，而控制器则根据参数估计器不断送来的参数估计值，通过一定的控制算法，按某一性能指标，不断形成最优控制作用，由于存在着多种参数估计和控制器设计算法，所以自校正控制的设计方法很多，其中，以用最小二乘法进行参数估计，按最小方差来形成控制作用的最小方差自校正控制器最为简单，并获得较多应用。
10.2.2最小方差自校正控制器
以一个单输入、单输出、线性时不变的系统为讨论对象。假如系统经常处在随机扰动作用之下，如不加以控制、随机扰动将使系统输出相对其设定值产生很大波动。或者说，系统输出的稳态方差将很大。最小方差自校正控制就是根据给定的对象数学模型，综合一个最优控制津，使系统输出的稳态方差最小 。
设被控制对象的数学模型可表示为
A(q-1)y(k)= (q-d)B(q-1)u(k)+C(q-1)e(k) (10.2-1)
式中A(q-1)=1+a1q-1+、、、+anq-n
B(q-1)=b0+b1q-1+、、、+ b0q-n b00
C(q-1)=1+C1q-1+、、、+ c1q-n
U(k),y(k)为系统的输入，e(k)为均值为0的随机噪声，d为延时时间次数（d1）。
对于系统（10.2-1），设输出所跟踪的目标值yr(k)已知，要求求解u(k)使方差E{[y(k)-yr(k)]2}最小。如果系统中参数是已知的，这个问题就是最小方差控制问题。而如果参数未知（或不断缓慢变化）,要在线递推估计和不断修正模型的参数值且使系统输出方差保持最小，这就成为最小方差自校正控制器问题。
为求解该问题，首先用A(q-1)=1+q.q-1+an b0q-n并除C（q-1）=Hc1q-1+、、、cnq-n，进行d次除法后，C（q-1）可写成。
C（q-1）=A（q-1）R （q-1）+ q-dS(q-1） (10.2-2)
式中R（q-1）=1+r1q-1+、、、+sn-1q-d+1
S（q-1）=so+s1q-1+、、、sn-1q-n+1
由(10.2-2).(10.2-1)式，有，
C（q-1）y(k)=A（q-1）R（q-1）y(k )+q-dS(q-1)y(K)
=q-dR(q-1)B(q-1)u(k)+q-dS(q-1)y(k)+R(q-1)C(q-1)e(k)
C(q-1)[y(k)-yr(k)]=R(q-1)B(q-1)u(k-d)
+S(q-1)y(k-d)-C(q-1)yr(k)+C(q-1)R(q-1)e(k) (10.2-3)
上式右边的前三项u( )、y( )、yr( )是能够测定的，此外因为
R(q-1)B(q-1)=[b0,(b1+r1)q-1,、、、bmrd-1]
=[B0 B1、、、Bn+d-1]
且令
=[s0,s1,、、、sn-1,B0,B1,Bn+d-1,c1,c2,、、、cn]
X（k-d）=
则（10.2-3）式可写成
C(q-1)[y(k)-yr(k)]=X(k-d)-yr(k)+C(q-1)R(q-1)e(k) (10.2-4)
y(k)-yr(k)=c-1(q-1)[ x(K-d)-yr(k)]+R(q-1)e(k) （10.2-5）
由(10.2-4)欲使E{[y(k)-yr(k)]2}最小，须使
X(k-d)-yr(k)=o (10.2-6)
这时y(k)-yr(k)=R(q-1)e(k) (10.2-7)
如果已知，从（10.2-6）中求u(k-d)可表示为
u(k-d)=[yr(k)-]
u(k)= [yr(k+d)- ] (10.2-8)
式中，0、X0分别表示从、X中抽去B0和u(k-d)后的矢量。但如果0、B0是未知的，便不能由（10.2-8）确定 u(k),而须回到（10.2-6），考虑由X和Yr来求的估计问题，这时可以用递推最小二乘方法来进行，其递推公式为

可求得的估计值，并由此求得o,,则
u(k)=[yr(k+d)-Xo(k)]
由此而获得自校正控制输入
例10.2 考虑一自校正控制问题，设受控对象为
y(k)-0.5y(k-1)-0.14y(k-2)=u(k-2)+0.5u(k-3)+0.2e(k)
解：用MATLAB编制设计程序（ chp10_2.m），程序中设置了一个变化的目标值yr，程序运行结果如图10.6所示，该仿真结果表示，所设计的最小方差自校正控制器可使系统在有随机干扰的情况下使系统输出很好地跟踪设定值的变化。
% --------- Chp10_2 ---------
global q1 pm1 id y1_v r1
y1=zeros(26,1);u1=zeros(26,1);v1=zeros(26,1);r1=zeros(26,1);
q1=zeros(26,1);mu1=zeros(10,1);my1=zeros(10,1);mn1=zeros(10,1);
pm1=zeros(26,26);sf=zeros(2,1);f=zeros(10,1);
py1=zeros(800,1);py2=zeros(800,1);
c=1;
% Init
nf=1;di=1;n=2;m1=1;Id=0.99;um=80;mu=2;md=2;
p0=1;
for i=1:25
pm1(i,i)=p0;
end
y1_v=0;contr=0;l=n+mu+md;
f(1:5)=[0.5 0.14 1.0 0.5 0.2]';
% Set point valve
R=zeros(800,1);
for i=1:200
R(i)=20;
end
for i=201:400
R(i)=40;
end
for i=401:650
R(i)=30;
end
for i=651:800
R(i)=20;
end

for ij=1:800
% System_Model
for j=9:-1:1
my1(j+1)=my1(j);
mu1(j+1)=mu1(j);
mn1(j+1)=mn1(j);
end
my1(1)=y1_v;
mu1(1)=contr;
mn1(1)=nf*randn(1);
y1_v=0;
for i=1:n
y1_v=y1_v+f(i)*my1(i);
end
for i=1:mu
y1_v=y1_v+f(i+n)*mu1(i+md);
end
for i=1:m1
y1_v=y1_v+f(i+n+mu)*mn1(i);
end
% X_vector
for i=25:-1:1
u1(i+1)=u1(i);
end
u1(1)=mu1(1);
for i=1:n
r1(i)=y1(i);
end
for i=1:mu+md
r1(i+n)=u1(i+md);
end
%Computes estimates recursively
F=subrls(l);
for i=25:-1:1
y1(i+1)=y1(i);
end
y1(1)=y1_v;
for i=1:n
r1(i)=y1(i);
end
for i=1:mu+md
r1(i+n)=u1(i);
end
% Compute Y(k)
Qi=0;
for i=1:n
Qi=Qi+q1(i)*r1(i);
end
for i=1:mu
Qi=Qi+q1(i+n+1)*r1(i+n);
end
%Compute u(k)
B0=q1(n+1);
if B0==0
B0=1;
end
contr=(R(ij)-Qi)/B0;
if abs(contr)>um
if contr>0
contr=um;
else
contr=-um;
end
end
py1(ij)=y1_v;
py2(ij)=contr;
end
t=1:800;t=t';
plot(t,R,'g',t,py1,'b',t,py2,'r')
grid

% --------- subrls ---------
function result=subrls(l)
global q1 pm1 id y1_v r1
v1=zeros(30,1);
z1=0;ww=0;
for i=1:l
v1(i)=0;
for j=1:l
v1(i)=v1(i)+pm1(i,j)*r1(j);
end
end
for i=1:l
z1=z1+r1(i)*q1(i);
ww=ww+r1(i)*v1(i);
end
ww=1+ww;
for i=1:l
q1(i)=q1(i)+v1(i)*(y1_v-z1)/ww;
for j=1:l
pm1(i,j)=(pm1(i,j)-v1(i)*v1(j)/ww);
end
end
result=l;

图10.6 例10.2程序运行结果

14