线性系统结课论文.docx
《线性系统结课论文.docx》由会员分享,可在线阅读,更多相关《线性系统结课论文.docx(16页珍藏版)》请在冰豆网上搜索。
线性系统结课论文
线性系统理论Matlab实践
学院:
电气工程与信息工程学院
专业:
电气工程
姓名:
学号:
1.本题目是在已知状态空间描述的情况下要求设计一个状态反馈控制器,从而使得系统具有实数特征根,并要求要有一个根要大于5,所以本设计中取两个特征根分别为8,9。
要对系统进行状态反馈的设计首先要判断其是否能控,即求出该系统的能控性判别矩阵,然后判断其秩,从而得出其是否可控;
判断能控程序设计如下:
>>A=[-0.80.02;-0.020];
>>B=[0.051;0.0010];
>>Qc=ctrb(A,B)
Qc=
0.05001.0000-0.0400-0.8000
0.00100-0.0010-0.0200%能控性判别矩阵
>>Rc=rank(Qc)
Rc=
2
>>P=[89];
>>K=place(A,B,P)
K=
1.0e+003*
-0.0200-9.0000
-0.00780.4500
得出结果能控型判别矩阵的秩为2,故而该系统是完全可控的,故可以对其进行状态反馈设计。
程序中所求出的K即为所求状态反馈控制器的状态反馈矩阵,即由该状态反馈矩阵所构成的状态反馈控制器能够满足题目要求。
2、(a)运用ctrb函数计算系统的能控性矩阵,并验证恒速制导导弹的运动系统是不可控的
设计程序如下:
>>A=[01000;-0.1-0.5000;0.50000;001000;0.51000];
>>B=[0;1;0;0;0];
>>C=[00010];
>>Qc=ctrb(A,B)
Qc=
01.0000-0.50000.1500-0.0250
1.0000-0.50000.1500-0.0250-0.0025
000.5000-0.25000.0750
0005.0000-2.5000
01.00000-0.10000.0500
>>Rc=rank(Qc)
Rc=
4
从程序运行的结果可得,系统能控型判别矩阵的秩为4,而系统为5阶系统,故而就验证了该系统为不可控的。
而该系统的能控型矩阵就为程序中的Qc矩阵。
(b)求出其传递函数,并建立新的状态变量模型。
设计程序如下:
A=[01000;-0.1-0.5000;0.50000;001000;0.51000]
B=[0;1;0;0;0]
C=[00010]
D=0
%[num,den]=ss2tf(A,B,C,D)%状态空间转传函
%s1=tf(num,den)
%s2=zpk(s1)%传函转零极点
symss
G=simple(simple(C*(inv(s*eye(5)-A))*B+D))%有疑问与ss2tf
%G=50/s^2/(10*s^2+5*s+1)
[a,b,c,d]=tf2ss([50],[105100])%传函转状态空间
结果如下:
A=
01.0000000
-0.1000-0.5000000
0.50000000
0010.000000
0.50001.0000000
B=
0
1
0
0
0
C=
00010
D=
0
G=
50/s^2/(10*s^2+5*s+1)
a=
-0.5000-0.100000
1.0000000
01.000000
001.00000
b=
1
0
0
0
c=
0005
d=
0
(c)要求证明其上所求的新系统为可控的,只需求其新系统能控型判别矩阵的秩,看是否为4即可。
设计程序:
A=[01000;-0.1-0.5000;0.50000;001000;0.51000]
B=[0;1;0;0;0]
C=[00010]
D=0
%[num,den]=ss2tf(A,B,C,D)%状态空间转传函
%s1=tf(num,den)
%s2=zpk(s1)%传函转零极点
symss
G=simple(simple(C*(inv(s*eye(5)-A))*B+D))%有疑问与ss2tf
%G=50/s^2/(10*s^2+5*s+1)
[a,b,c,d]=tf2ss([50],[105100])%传函转状态空间
qc=rank(ctrb(a,b))
n=4
ifqc==n
disp('能控')
elseifqc~=n
disp('不能控')
end
结果如下:
A=
01.0000000
-0.1000-0.5000000
0.50000000
0010.000000
0.50001.0000000
B=
0
1
0
0
0
C=
00010
D=
0
G=
50/s^2/(10*s^2+5*s+1)
a=
-0.5000-0.100000
1.0000000
01.000000
001.00000
b=
1
0
0
0
c=
0005
d=
0
qc=
4
n=
4
能控
(d)要求判断系统的稳定性,可以采用李雅普洛夫特征值法进行判定;
设计程序:
A=[01000;-0.1-0.5000;0.50000;001000;0.51000];
>>d=eig(A)
d=
0
0
0
-0.2500+0.1936i
-0.2500-0.1936i
从求得的结果中可以看出其特征值的根的实部都不是正数,从而就说明了该系统是李雅普洛夫意义下稳定的。
(e)讨论状态变量模型的能控性元复杂性的关系。
很直观地讲,一个系统要能控,必须要其能控型判别矩阵的秩等于系统的阶数也即就是状态变量的数目,但是反过来,系统越复杂,状态变量的个数越多,能控型判别矩阵要求满足的秩也就越大,也即意味着越难达到要求,从而其能控性也就越不容易满足。
从而可以得出结论,即越复杂的系统越不容易达到完全可控。
3.(a)计算系统矩阵A的特征值,并判断系统是否是稳定的
程序:
A=[-0.03660.02710.0188-0.4555;0.0482-1.01000.0024-4.0208;0.10020.3681-0.70701.4200;0010];
>>d=eig(A)
d=
0.2758+0.2576i
0.2758-0.2576i
-0.2325
-2.0727
由系统矩阵A的特征根可知:
系统是不稳定的。
(b)利用poly函数确定A的特征多项式。
程序:
>>A=[-0.03660.02710.0188-0.4555;0.0482-1.01000.0024-4.0208;0.10020.3681-0.70701.4200;0010];
>>poly(A)
ans=
1.00001.7536-0.64720.06250.0686
综上a)计算结果可知,eig()函数和poly()函数计算得到的特征根的值是不相等的,但是两者的相同点是都可以判断系统的稳定性。
(c)分别判断当u1与u2发挥作用时,系统的能控性。
求得u1与u2发挥作用时不同的能控型判别矩阵,从而可判断出其能控性。
设计程序:
>>A=[-0.03660.02710.0188-0.4555;0.0482-1.01000.0024-4.0208;0.10020.3681-0.7071.42;0010];
>>B1=[0.4422;3.5446;-5.52;0];
>>B2=[0.1761;-7.5922;4.49;0];
>>Qc1=ctrb(A,B1)
Qc1=
0.4422-0.02392.5172-2.0267
3.5446-3.572025.8140-47.0978
-5.52005.2517-12.868626.3099
0-5.52005.2517-12.8686
>>Rc1=rank(Qc1)
Rc1=
4
>>Qc2=ctrb(A,B2)
Qc2=
0.1761-0.1278-1.94412.3338
-7.59227.6874-25.838149.9646
4.4900-5.951513.4004-27.6310
04.4900-5.951513.4004
>>Rc2=rank(Qc2)
Rc2=
4
求得当u1与u2发挥作用时所对应的能控型判别矩阵的秩都为4,即其秩等于系统的阶数,从而说明在这两种情况下,系统均为能控。
4.(a)卫星的定点位置是否稳定
程序如下:
A=[000100;000010;000001;7.380900020;0-2.19040-200;00-3.1904000];
>>d=eig(A)
d=
-2.1587
2.1587
0.0000+1.8626i
0.0000-1.8626i
0+1.7862i
0-1.7862i
可以看到系统矩阵A的特征值中有一个为正数,从而可得出该系统为不稳定的。
b)、c)、d).只有u1/u2/u3发挥作用时,系统的可控性
程序:
>>A=[000100;000010;000001;7.380900020;0-2.19040-200;00-3.1904000];
>>B1=[0;0;0;1;0;0];
>>B2=[0;0;0;0;1;0];
>>B3=[0;0;0;0;0;1];
>>QC1=ctrb(A,B1);
>>QC2=ctrb(A,B2);
>>QC3=ctrb(A,B3);
>>ifrank(QC1)==size(A)
disp('u1作用下可控')
else
disp('u1作用下不可控')
end
u1作用下不可控
>>ifrank(QC2)==size(A)
disp('u2可控')
else
disp('u2不可控')
end
u2不可控
>>ifrank(QC3)==size(A)
disp('u3可控')
else
disp('u3不可控')
end
u3不可控
e)确定由u2到该位置的漂移量的传递函数
程序:
>>A=[000100;000010;000001;7.380900020;0-2.19040-200;00-3.1904000];
>>B2=[0;0;0;0;1;0];
>>C=[010000];
>>symss
>>I=eye(6);
>>F=inv(s*I-A);
>>G=C*F*B2
G=
625*(10000*s^2-73809)/(6250000*s^4-7440625*s^2-101044521)
f)计算由上式得到的传递函数的状态变量模型,并验证该子系统是能控系统。
程序:
num=[6250000046130625];
>>den=[6250000074406250-101044521];
>>[a,b,c,d]=tf2ss(num,den)
a=
0-1.1905016.1671
1.0000000
01.000000
001.00000
b=
1
0
0
0
c=
01.000007.3809
d=
0
能控性判断程序:
>>Q=ctrb(a,b);
>>ifrank(Q)==size(a)
disp('能控')
else
disp('不能控')
end
能控
g)程序:
>>P=[-1+j-1-j-10-10];
K=acker(a,b,P)
K=
22.0000140.8095240.0000216.1671
5.根据风电机一阶数学模型,设计风力机转速的闭环PI控制,首先绘制出带有PI控制器的风力机控制系统,如图所示:
图中PI控制器的模型为K1+K2/s,从而可以求出该系统闭环系统的特征方程为
S^2+0.3397K1*S+0.3397K2=0
要使得该风力机稳定运行,则需让特征方程的根具有负实部。
设计程序如下:
symssk1k2;
>>s=solve(s^2+0.339*k1*s+0.339*k2)
s=
-(339*k1)/2000-((114921*k1^2)/1000000-(339*k2)/250)^(1/2)/2
((114921*k1^2)/1000000-(339*k2)/250)^(1/2)/2-(339*k1)/2000
使其为负数,从而可得出K1<=180,K2<=86,得到的PI控制器就能满足要求。
6.风力机三阶模型,采用状态反馈和极点配置法,设计风力机转速的闭环控制系统。
对于单输入时间线性时不变系统,基于状态反馈课任意配置全部闭环极点的条件为:
系统完全能控。
指定一组闭环极点-7-11-9,设计满足该极点的状态反馈矩阵。
设计程序:
>>A=[-0.8468-8.1598e-0080.5017;8676370000-867637000;1.1636e+0040.0019-1.1636e+004];
>>B=[-0.0262;0;0];
>>R=rank(ctrb(A,B))
R=
3
>>ifR==size(A)
disp('完全能控')
else
disp('不完全能控')
end
完全能控
>>P=[-7-11-9];
>>K=place(A,B,P)
K=
1.0e+005*
4.43120.0000-4.4311
7.对风力机三阶模型系统进行LQR控制器设计,由于系统矩阵和输入矩阵均为已知,则只需定出两个加权矩阵Q,R,不妨定加权矩阵Q=diag([1,1,1]);R=1;然后就可以直接利用lqr函数来进行LQR控制器的设计了。
程序设计如下:
>>A=[-0.8468-8.1598e-0080.5017;8676370000-867637000;1.1636e+0040.0019-1.1636e+004];
>>B=[-0.0262;0;0];
>>Q=eye(3);
>>R=1;
>>K=lqr(A,B,Q,R)
K=
1.0e+004*
-6.4633-0.00016.4633
8.设计风力机三阶模型的状态反馈及其状态观测器,并且为降维观测器,由于输出矩阵C的秩为1,故而降维观测器的维数为3-1=2;不妨设状态反馈的期望极点为-3,-4,-5,状态观测器的期望极点为-7,-8
设计程序:
>>A=[-0.8468-8.1598e-0080.5071;8676370000-867637000;1.1636e+0040.0019-1.1636e+004];
>>B=[-0.0262;0;0];
>>C=[001];
>>Qc=ctrb(A,B)
Qc=
1.0e+011*
-0.00000.0000-0.0000
0-0.00022.6453
0-0.00000.0000
>>Rc=rank(Qc)
Rc=
3%系统是完全能控的
>>Qo=obsv(A,C)
Qo=
1.0e+008*
000.0000
0.00010.0000-0.0001
-1.3376-0.00001.3375
>>Ro=rank(Qo)
Ro=
3%系统是完全能观测的
>>P1=[-3-4-5];%任意配置极点
>>k1=place(A,B,P1)
k1=
1.0e+005*
4.43700.0000-4.4368
>>P=[001;010;100];
>>AA=P*A*(inv(P))%对系统矩阵A进行线性非奇异变换
AA=
1.0e+008*
-0.00010.00000.0001
-8.676408.6764
0.0000-0.0000-0.0000
>>AA11=[AA(1,1)]
AA11=
-11636
>>AA12=[AA(1,2:
3)]
AA12=
1.0e+004*
0.00001.1636
>>AA21=[AA(2:
3,1)]
AA21=
1.0e+008*
-8.6764
0.0000
>>AA22=[AA(2:
3,2:
3)]
AA22=
1.0e+008*
08.6764
-0.0000-0.0000
>>BB=P*B
BB=
0
0
-0.0262
>>BB1=[BB(1,1)]
BB1=
0
>>BB2=[BB(2:
3,1)]
BB2=
0
-0.0262
>>AX=AA22';
>>BX=AA12';
>>P2=[-7-8];
>>KX=place(AX,BX,P2);
>>L=KX'
L=
1.0e+003*
7.5225
-0.0000
>>R=AA22-L*AA12
R=
1.0e+008*
-0.00007.8011
-0.0000-0.0000
>>S=R*L
S=
1.0e+005*
-1.1687
-0.0000
>>T=AA21-L*AA11
T=
1.0e+008*
-7.8011
0.0000
>>U=BB2-L*BB1
U=
0
-0.0262
程序中所求出的k1即为状态反馈控制器的状态反馈矩阵,L为降维观测器所求矩阵,而R,S,T,U矩阵则为构成状态反馈器所需要得几个中间矩阵。
由程序所求的矩阵则可构成所求状态反馈控制器与降维观测器。