线性系统理论MATLAB大作业

兰州理工大学2015级线性系统理论大作业

线性系统理论Matlab实验报告

1、在造纸流程中，投料箱应该把纸浆流变成2cm的射流，并均匀喷泗在网状传送

带上。为此，要精确控制喷射速度和传送速度之间的比例关系。投料箱内的压力是

需要控制的主要变量，它决定了纸浆的喷射速度。投料箱内的总压力是纸浆液压和

另外灌注的气压之和。由压力控制的投料箱是个耦合系统，因此，我们很难用手工

方法保证纸张的质量。

在特定的工作点上，将投料箱线性化，可以得到下面的状态空间模型：

---0.8+0.02]「().05I-

x=戈+u

-0.020JLo.ooi0

"卜，为]

其中，系统的状态变量液面高度，x2二压力，系统的控制变量ul二纸浆流量u2二

气压阀门的开启量。在上述条件下，试设计合适的状态变量反馈控制器，使系统具

有实特征根，且有一个根大于5

解：本题目是在已知状态空间描述的情况卜一要求设计一个状态反馈控制器，从而使

得系统具有实数特征根，并要求要有一个根的模值要大于5,而特征根是正数时系

统不稳定，这样的设计是无意义的，故而不妨采用状态反馈后的两个期望特征根为

-7,-6,这样满足题目中所需的要求。要对系统进行状态反馈的设计首先要判断其

是否能控，即求出该系统的能控性判别矩阵，然后判断其秩，从而得出其是否可控。

Matlab判断该系统可控性和求取状态反馈矩阵K的程序,如图1所示，同时求

得加入状态反馈后的特征根并与原系统的特征根进行了对比。

1-

2%设置的状态变量反馈控制器

3-A=[-0.80.02-0.020];B=[0.051;0.0010]:C=[10:01];D=0:

4-t=0:0.05:10;

5-Qc=ctrb(A,B);%首先判断系统的能控性

6—ifrank(Qc)==size(A)

7—disp('系统能控’)

8—else

9-disp('系统不可控’)

10-end

11-p=[-7-6]：%不放妨假定期望的特征根为-6分别满足给定的要求

12-ifrank(A)==size(A)

13-disp('极点个数与系统配置合理’)

14-else

15-dispC极点个数与系统配置不合理’)

16-end

17-K=place(A.B,p)%求取状态反馈u=-Kx的系数矩阵K

18-rootsnew=eig(Al)%加入状态反馈后新系统的特征根

19-rootsold=eig(A)%原系统的特征根

20-A1=A-B*K:B1=B:C1=C:D1=D：%加入状态反馈后新的状态矩阵

21

图1系统能控性、状态反馈矩阵和特征根的分析程序

上述程序的运行结果如图2所示:

图2系统能控性、反馈矩阵和特征根的运行结果

图2中为图Imatlab程序的运行结果，经过判断得知系统是可控的，同时极点

的配置个数与系统状态相符，求得了状态反馈矩阵K的值，并把原系统的特征根

(rootsold)和加入状态反馈后的特征根(rootsnew)进行对比。同时通过特征值

可以看出该系统是稳定的。

2、描述恒速制导导弹的运动方程为：

01ooo-T

-0.1-0.50001

x=().500000

0010000

0.510000

y=[o00Io]x

运用ctrb函数计算系统的能控型矩阵，并验证系统是不可控的；

计算从u到Y的传递函数，并消去传递函数中的分子和分母公因式，由此可以得到

能控的状态空间模型。在消去了公因子之后，请用tf2ss函数确定新的状态变量模

型；

证明(b)中得到的状态变量模型是能控的；说明恒速制导导弹是否稳定？

讨论状态变量模型的能控性和复杂性的关系(假设用状态变量的数目来度量复杂

性)。

解：该题是通过描述的恒速制导导弹的运动方程求解相应问题。

(a)运用ctrb函数计算系统的能控性矩阵,并判断该系统不可控，详细matlab

程序和判断结果如图3和图4所示。

图3是判断该系统能控性的matlab程序，通过求得能控性矩阵Qc,并通过秩

判据来判定该系统是否能控。

・

1-clc▲□

2%恒速制导导弹系统mat1ab分析

3%(a)

4-A=[0I000;-0.1-0.5000;0.50000;

5001000;0.51000];

6-B=[0I000]':

7一C=[00010]：

8一Qc=ctrb(A,B)%系统的能空性矩阵

9-ifrank(Qc)==5%判断系统的可控性(秩判据)

10-dispC系统能控,)

11-else

12-disp('系统不可控’)

13-end

14d

图3系统能控性的判别程序

判定的结果如图4所示:

CommandWindow-□X

FileEditDebugDesktopWindowHelp

T•NewtoMATLAB?WatchthisVideo,seeDemos,orreadGettinaStarted.X

jd

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

系统不可控

A»l

图4系统的能控性矩阵和能控性判定结果

通过matlab分析求得了系统的能控性矩阵Qc,同时通过秩判据判定该系统不可控。

(b)、(c)计算u到y的传递函数，并通过tf2ss函数确定新的状态变量模型，同

时判断该模型是能控的。具体程序如图5所示，判断的结果如图6示。

1-clc

2一A=[01000;-0.1-0.5000:

30.50000:001000:0.51000];

4-B=[010001'：

5-C=[00010]：

6-symss：

7一Gl=C*inv(s*eye(size(A))-A)*B:

8-num=[50];den=[105100];

9-sys=tf(num,den)

10-[Al,Bl,Cl,Dl]=tf2ss(num,den);

11-sys=ss(Al,Bl,Cl,DI)学新的状态空间模型

12-Qcl=ctrb(Al,Bl);%系统的能空性矩阵

13-ifrank(Qcl)==size(Al)%判断系统的可控性(秩判据)

14-disp('系统能控')

15-else

16-disp('系统不可控')

17-end|

图5确定新状态空间并判定能控性的程序

Winciow

BteEditDef2MgQesktopV^odowHelp~

<T■NewtoMATLAB?WatchthtsVideo,seeDemos,orreadGettmqStarted.X

Irans-fer-funcXion:

50

10s*4+5s*3+s*2

Corvtxnuous-*xmemode1.

fx系统能控

|OVR

图6系统的传递函数、新的状态空间模型和能控性判定结果

分析得知U到y的传递函数可通过状态空间描述的矩阵求得，同时通过tf2ss

函数确定了新的状态空间(Al,Bl,Cl,D1),运用函数ss求得新模型的状态方程,

再通过能控型矩阵判定系统的能控性。

显然得到系统是可控的，同时还要声明通过传涕函数求得空间描述和通过状态

矩阵求得结果不同，从而验证了传递函数对系统的内部描述不完整。

(d)判断恒速制导导弹系统稳定性

以下通过求得矩阵的特征值即传递函数的极值点来判断该系统是否稳定。图7

是求取极值点的程序，通过roots和eig函数来求取，目的进行必要的对比。图8

是通过两种途径获得的系统的极值点。

I1-clc■

2W恒速制导导弹系统matlab分析

3%(a)I

4-A=[01000;-0.1-0.5000;

50,50000:001000;0.51000];

6-B=[01000]'：

7-C=[00010];

8-symss;

9-Gl=C*iny(s*eye(size(A))-A)*B;■

10-sys=tf^50,[10,5,1,0,0])-

ii-TO：OJ：IOO『；)

12-y=step(sys,t);plot(t,y)

13-p=[10,5,1,0,0];

14-Xl=roots(p)•

15-X2=eig(A)

图7求取极值点的源程序

图8是图7程序的运行结果:

图8系统的传递函数和极值点

从求得的结果中可以看出其特征值的根的实部都不是正数，从而就说明了该系

统在李雅普洛夫意义下是稳定的。

图9系统的单位阶跃响应

通过程序给系统一个单位阶跃信号，从上图可以看出系统不是严格收敛的，而

是发散的。

(e)状态变量模型的能控性和复杂性的关系(用状态变量的数目来度量复杂性)。

讨论状杰变量模型的能控性与复杂性的关系。很直观地讲，一个系统要能控,

必须要其能控型判别矩阵的秩等于系统的阶数也即就是状态变量的数目，但是反过

来，系统越复杂，状态变量的个数越多，能控型判别矩阵要求满足的秩也就越大,

也即意味着越难达到要求，从而其能控性也就越不容易满足。从而可以得出结论,

即越复杂的系统越不容易达到完全可控。

3、垂直起降的飞机的线性化模型为：

=Ax+Blul+B2u2

其中

--0.03660.02710.01880.4555一

0.0482-1.01000.0024-4.0208

A=

0.10020.3681-.0.70701.4200

0010

--

■0.4422-■0.1761-

3.5446-7.5922

=B=

-5.520024.4900

00

系统的状态变量为水平速度为(节)、垂直速度工(节)、倾斜率不(度/秒)和倾

斜角/(度)；系统的控制输入为多和内，其中场用于控制垂直运动，的用于控制

水平运动。

(a)计算系统矩阵A的特征值，并由此判断系统是否稳定；

(b)利用poly函数确定A的特征多项式，计算特征根，并与(a)中得到的特征根相

比较；

(c)当只有叫发挥作用时，系统能控吗？当只有的发挥作用时，结果又如何？请比

较解释你的结论。

解：通过给定的垂直起降的飞机的线性化模型分析系统的属性

(a)计算系统矩阵A的特征值，并根据特征值判断系统是否稳定

1-clc

2-A=[-0.03660.02710.0188-0.4555;

30.04S2-1.01000.0024-4.0208;

40.10020.3681-0.70701.4200;0010];

5-rootsl=eig(A)瑞特征值的实部全小于卷则系统稳定(a)

6-Q=poly(A)%(b)

7-Ql=localpoly(Q)

8-roots2=roots(Ql)

9-Bl=[0.44223.5446-5.52000]'只有ul作用时判断系统的能控性

10-Qcl=ctrb(A,Bl);%只有ul作用时系统的能控性矩阵

11-dimA=size(A)

12-ifrank(Qcl)==4确断系统的可控性(秩判据)

13-dispC只有ul作用时系统可控')

14-else

15-disp('只有ul作用时系统不可控’)

end

B2=[0.1761-7.59224.49000]':%(c)只有u2作用时判断系统的能控性

18—Qc2=ctrb(A,B2);%只有u2作用时系统的能控性矩阵

19-ifrank(Qc2)==4%判断系统的可控性(秩判据)

20-dispC只有u2作用时系然可控’)

21-else

22-dispC只有作用时系统不可控’)

23-end

图10矩阵A的特征值和Ul、U2分别作用的能控性判别程序

(b)利用poly函数确定A的特征多项式，计算特征值，并与(a)中的结果进行

对比

(c)当只有ul作用时，系统能控性；只有u2作用时，系统能控性。

针对以上二点问题，通过图10所示的matlab程序来判断所相问题，最终的结果在

图11中显示。

求取矩阵A的特征值和ul.u2分别作用时系统可控性的运行结果：

rootsi=

0.2758+0.2576i

0.2758-0.2576i

-0.2325

-2.0727

Q=

1.00001.7536-0.64720.06250.0686

QI=

Xl*4+1.75*X1*3-0.647*X1*2+0.0625*Xl+0.0686

roots2=

-2.0727

0.2758+0.2576i

0.2758-0.2576i

-0.2325

dimA=

44

只有ul作用时系统可控

只有u2作用时系统可控

fx»

图11特征值、特征多项式和Ul、u2分别作用的能控性结果

其中rootsl是通过eig函数求得的状态矩阵A的特征值，显然有两个特征值

具有正实部，故系统不稳定；Q1是通过poly函数确定的A的特征多项式，roots2

是通过roots函数求得的A矩阵的特征多项式的根,经过对比发现rootsl和roots2

的数值一样；只有ul或者u2作用是通过能控型矩阵Qc,用秩判据得到系统都是可

控的。dimA是通过size函数求得矩阵A的维数。对比的当ul与u2发挥作用时所

对应的能控型判别矩阵的秩都为4,即其秩等于系统的阶数也就是矩阵A的维数，

从而说明在这两种情况下，系统均为能控。

4、为了探究月球背面(远离地球的一面)的奥秘，人们付出了不懈的努力。例如,

在地球-太阳-月球系统中，人们希望通信卫星能定点在不受月球遮挡的轨道上，并

为此开展了广泛的论证研究工作。图中给出了预期卫星轨道的示意图，从地球上看

上去，卫星轨道的光影恰似环绕月球的外层光晕，因此这种轨道又称为光晕轨道。

轨道控制的目的是，使通信卫星在地球可见的光晕轨道上运行，从而保证通信链路

的畅通，所需的通信链路包括从地球到卫星和从卫星到月球背面共两段线路。

卫星绕定点位置运动时，经过标准化和线性化的漂移运动方程为:

000100-000

000010000

■000001000

X—X+«1+u+

7.3809000201020

0-2.19040-200010

00-3.1904000001

其中，状态变量x是卫星在三个方向上的位置和速度漂移，输入〃g=L2,3)分别是

轨控发动机在4、〃和4方向上产生的加速度。

(a)卫星的定点位置是否稳定？(b)如果只有/发挥作用，卫星是否能控？

(c)如果只有归发挥作用，卫星是否能控？(d)如果只有内发挥作用，卫星是否能

控？

(e)如果能够测得〃方向的位置漂移，请确定由电到该位置漂移量的传递函数。(提

示：可以令观测输出为J=[o10000]x)

(f)用tf2ss函数，计算(。)中得到的传递函数的状态变量模型，并验证该轨迹子

系统是能控系统；

(g)采用状态反馈〃设计合适的反馈控制器，使(f)中得到的系统的闭环

极点为52=T土/和邑4=-10。

解：在给定的卫星绕定点位置运动时的标准化和线性化的漂移运动方程，通过

matlab分析一下几点问题。

1-clc

2一A=[000100;000010:000001:7.380900020：

30-2.19040-200;00-3.1904000]：

4-num=1:den=poly(A)

5-sys=tf(num,den)：

6—pzmap(sys)

7—rootsl=eig(A)

8—Bl=[000100]'：

9-Qcl=ctrb(A,Bl);

10-dimA=size(A)

11-ifrank(Qcl)==6

12-dispC只有ul作用时系统可控’)

else

dispC只有ul作用时系统不可控')

15-end

16—B2=[000010],：

17-Qc2=ctrb(A,B2);

18-ifrank(Qc2)=6

19-dispC只有u?作用时系统可控’)

20-else

21-dispC只有u2作用时系统不可控’)

22-end

23-B3=[000001]';

24-Qc3=ctrb(A,B3)：

25-ifrank(Qc3)==6

26-dispC只有u3作用时系统可控')

27-else|

28-dispC只有u3作用时系统不可控')

29-end

图12系统稳定性和ul、u2、u3分别作用时的能控性

(1)关于卫星的定点位置的稳定性和分别只有ul或者u2或者u3作用时，卫星的

能控性通过图12的程序来判断，判断结果在图13中显示。

卫星定位系统的稳定性和ul、u2、u3分别作用时的能控性判别结果如图13所示

rootsl

-2.1587

2.1587

0.0000+1.8626i

0.0000-1.8626i

0+1.7862i

0-1.7862i

66

只有ul作用时系统不可控

只有u2作用时系统不可控

只有u3作用时系统不可控

图13系统特征根和ul、u2、u3分别作用的能控性判别结果

Pole-Zeror.!ap

W依

A

J

e

u

_

6

e

E

-

1

.5-2-1.5-1-OS00.511.522.5

RealAxis

图14系统极值点分布图

通过图13可以看出系统的极值点(rootsl)中有大于零的点，直观的从图14

的系统极值点分布图中看出在虚轴的右半平面上有一个极值点，所以系统是不稳定

的；从图13中还可以发现系统在只有ul或者u2或者u3作用时，均不可控。

(2)确定由u2到漂移量的传递函数并确定传递函数所对应的状态变量模型，然后

验证其为能控系统。执行程序如图15所示，该程序用于求解传递函数和状态模型,

并验证该模型的能控性。运行的结果如图16所示。

1-clc

2-symss;

3%卫星的状态瓜atlab分析(e)

4-A=[000100;000010;000001;7.380900020;

50-2.19040-200:00-3.1904000];

6-B2=[00Q010]';

7-C=[010000];

8-G=collect(factor(C*iny(s*eye(6)-A)<B2));

9-num=[62500000-46130625];

10-den=[62500000-74406250-101044521];

11-sys=tf(num,den)

12-[Al,Bl,CljDl]=tf2ss(num,den)

13-sys=ss(ALBl,CLDl)

14-Qc=ctrb(Al,Bl);

15-ifrank(Qc)==size(Al)

16-disp('系统能控')

17-else

18-disp(：系统不可控')

19-end

图15传递函数、状态变量模型和能控性求解程序

以上程序中求得了新系统的传递函数以及状态空间模型，并通过求取系统的能

控性矩阵，根据秩判据判定系统的可控性，由if语句来选取，把最终结果显示在

命令窗口。

图15程序的运行结果如图16所示：

Transferfunction:

6.25e006s'2-4.613e007

6.25eOC6s*4-7.441e006s*2-1.OleOOS

xlx2x3x4

X101.19016.17

x21000

x30100

x40010

xlx2x3

010-7.381

d=

ul

yl0

Contmuous-timemodel.

系统能控

A»

图16传递函数、状态空间描述和能控性的结果

图16中显示了由u2到n方向的位置漂移量的传递函数，以及通过tf2ss函数

得到该传递函数的状态变量模型，最后验证得到该模型是能控的。

(3)在给定状态空间描述的基础上采用状态反馈L2=-KX,使得(1)中得到的模型

的闭环极点为7+j,-1-j,TO,70.具体程序如图17所示，运行的结果如图18,

图19所示。

图17状态反馈设计的程序:

1-clc

2-symss：

3-num=[62500000-46130625]：

4-den=[62500000-74406250-101044521]

5-[Al,Bl,Cl,Dl]=tf2ss(num,den);

6—p=[-1+li-1-li-10-10]：

7-ifrank(Al)==4

8-disp（'极点个数与系统配置合理’）

9-else

10-disp（'极点个数与系统配置不合理'）

11-end

12-K=acker(Al,Bl,p)

13-agk=Al;bgk=B1;cgk=K;dgk=0;

14-subplot(2,1,1);

15-rlocus(agk,bgk,cgk,dgk)

16—grid

17-title（'加入反馈后系统的根轨迹图’）；

18-At=Al-B1*K:Bt=B1;Ct=C1;Dt=Dl;

19-[zl,pl,gainT]=ss2zp(At,Bt,Ct,Dt)

20-deg=degain(At,Bt,Ct,Dt):

21-t=0:0.1:10;

22-yc=step(At,Bt,Ct,Dt,1,t);ycl=yc/deg;

23-subplot(2,1,2);

24-plot(t,ycl)

25-grid

26-title（'加入反馈后系统的阶跃响应’）；

图17状态反馈设计程序

该程序首先判定极值点是否配置合理，求得反馈矩阵，并画出加入状态反馈后

的根轨迹，再求得加入状态反馈后系统的零极点及增益，最后画出系统加入反馈后

的阶跃响应图。

图18是上述程序运行得到的反馈矩阵和系统的零极点以及增益值。

极点个数与系统配置合理

22.0000143.1905240.0000216.1671

zl=

-2.7168

2.7168

Pl=

-10.0000+O.OOOOi

-10.0000-0.0000i

-1.0000+l.OOOOi

-1.0000-l.OOOOi

gainl=

执1.0000

图18状态反馈矩阵、零极点和增益

其中K是状态反馈矩阵，zl是加入状态反馈后的零点，pl是加入状态反馈后

的极点，gainl是系统增益，透过极点可以看出系统是稳定的。

通过图18可以看出再加入状态反馈K时，求得的极点正好是期望的值，也验

证了求取的正确性。最后求得加入状态反馈时系统的增益。通过图19可以看出系

统是渐近稳定的。

图19是加入状态反馈后系统的根轨迹和阶跃响应：

加入反馈后系统的根轨迹图

&91」.

s

_.佟…

x

v,4x99---

A

J17.&'

E

U

5

B<0:99-二1

E

-之96-…,

。.广、…

-15-10-505

RealXkxis

加入反馈后系统的阶跃响应

2345678910

图19加入反馈K时系统的根轨迹和阶跃响应

从图19看到系统是稳定的正好验证了极值点实部小于冬系统稳定，同时可以

发现系统在加入状态反馈后调节时间较快。

5、在&2风力机的一除模型中，采用浆距角控制风力机的转速，风速的变化视为

扰动，设计风力机转速的闭环P1控制，使转速恒定。

解：给定风力机的一阶模型，采用浆距角控制风力机转速，风速变化视为扰动，设

计风力机的PI控制，使转速恒定。

P1控制器的模型为K1+K2/S,从而可以求出该系统闭环系统的特征方程为

S^2+0.3397K1*S+O.3397K2=0

要使得该风力机稳定运行，则需让特征方程的根具有负实部。

symssklk2;

»s=solve(s^2+0.339*kl*s+0.339*k2)

-(339*kl)/2000-((114921*kr2)/1000000-(339*k2)/250)"(l/2)/2

((114921*kr2)/1000000-(339*k2)/250)71/2)/2-(339*kl)/2000

使其为负数,从而可得出Kl<=180,K2C86,得到的PI控制器就能满足要求，不妨

取Kl=80,K2=60.simulink仿真如图20所示。

G»in1

图20加入P1控制的系统仿真图

其中输入采用阶跃信号，扰动信号为单位脉冲。

在未加入扰动信号时，只有阶跃信号输入的系统响应如图21所示:

图21未加脉冲扰动只有阶跃信号输入的输出响应

图22在加入脉冲扰动情况下系统输出响应

通过图21和22对比发现系统能很好的通过反馈调节使系统在扰动情况下趋于

稳定，同时可以看出调节速度较快，所以得到的PI控制器可以满足所要求的指数。

6、在8.2风力机的三阶模型中，采用浆距角控制风力机的转速，风速的变化视为

扰动，电磁转矩视为常数，采用状态反馈和极点配置算法，设计风力机转速的闭环

控制系统。

解：理论上一般选取观测器的期望极点为传递函数极点的2到5倍为佳，在以下程

序中选取观测器的极值点为-2.2986c004,-300,-0.0286c004o具体求解程序如图

23所示.同时通过求取加入状态反馈后的伯德图，根轨迹以及单位阶跃响应来判定

系统的稳定性和调节速度。

通过图23的分析程序，我们能得到图24能控性矩阵和状态反馈矩阵K,图25

的加入状态反馈后系统的伯德图以及图26的根轨迹，图27的反馈系统的阶跃响应

曲线。

图23是利用极点配置算法求取风力机转速的闭环控制系统的程序：

2A=[-0.8468-8.1598e-0080.5071;8676370003-867637000：

31.1636e+0040.0019-1.1636e+004]:

4B=[-0.0262:0:0]：C=[001];D=0;

5Qc=ctrb(A,B)

6ifrank(Qc)==size(A)

disp('系统可控

8else

9dispC系统不可控’)

10end

11pe=[-2.2986e+004-300-0.0286e+004];%状态反馈的期望极点(选择合适的极点)

12w=[0.0110];t=0:0.05:10：

13K=place(A,B,pe)%状态反馈矩阵

14bode(A,B,K,0)：

15titleC加入状态反馈后系统的伯德图’)

16grid

17[ginpmwpcwgc]=jnargin(A,B,K,0);

18figure：

19rlocus(A,B,K,0):grid

20titled加入状容反馈后系统的根钛迹’)

21At=A-B*K;Bt=B;Ct=C:Dt=D:%加入状态反馈后系统的状态空间模型

22[zTplgainT]=ss2zp(At,Bt,Ct,Dt);

23dcg=dcgain(At,Bt,Ct,Dt):

24figure

25yc=step(At,Bt,Ct,Dt,1,t);yc=yc/deg:

26plot(t,yc);

27titleC加入状态反馈后系统的阶跃响应‘)

28grid

图23系统状态反馈矩阵的求取程序

图24能控性矩阵和状态反馈矩阵K

图25是加入反馈矩阵K时系统的伯德图:

加入状态反馈后系统的伯德图

80

60

40

m

p

)

名

n

j

u

e2O

s

-2O

-4O0

®

«

)

为

ue

d

10-210°102104106

Frequency(rad/sec)

图25加入状态反馈后系统的伯德图

加入状态反馈后系统的根轨迹

w

x

v

JA

UB

e&

-E

RealAxisx10*

图26加入状态反馈后系统的根轨迹图

图27是加入状态反馈后系统的单位阶跃响应:

加入状态反馈后系统的阶跃响应

1.4IIIIII

1.2

1

08

0.6

0.4

0.2

0II

C)12345678910

图27加入状态反馈后系统的单位阶跃响应

在求取系统状态反馈时，首先应该判断该系统的能控性，所以程序开始通过求

取能控型矩阵，运用秩判据来判断系统能控性，把结果显示在命令窗口。图24中

Qc是能控型矩阵，判断结果为系统可控，并求得状态反馈矩阵Ko通过图25,图

26,图27分析得知系统在加入反馈后，系统不仅稳定而且调节和响应速度较快。

7、给定风力机的三阶模型，采用浆距角控制风力机的转速，风速变化视为扰动，

设计风力机的LQR控制器。

解：在设计LQR控制器时，首先要确定加权矩阵Q(不妨选取

Q=diag([4e004,0,le004]),R=1),然后直接利用Iqr函数求取LQR控制器的状态反

馈矩阵k。求取程序如图28所示：

1

2-A=[-0.8468-8.1598—0080.5071:8676370000-867637000;

31.1636e+0040.0019-1.1636e+004];

4-B=[-0.0262;0：0];

5-Q=diag([4e+0040le+004])：%不妨定加权矩阵Q=d1agH2,2,2]):

6%R=1;然后就可以直接利用Iqr函数来进行LQR控制器的设计了。

i-R=l;

8-K=lqr(A,B,Q,R)%k即为LQR控制器中的状态反馈矩阵

9-bode(A,B?K,0);

10-title('加入状态反馈后系统的伯德图’)

11-grid

12-[gmpmwpcwgc]=margin(A,B,K,0):

13-figure;

14-rlocus(A,B,K,0);

15-title('加入状态反馈后系统的根轨迹')

16-grid

17-At=A-B*K;Bt=B;Ct=C;Dt=D;

IS—[zlplgainl]=ss2zp(At,Bt,Ct,Dt);

19-dcg=dcgain(At?Dt);

20一figure

21-yc=step(At,Bt,Ct,Dt?1,t);yc=yc/dcg;

22一plot(t,yc);

23-title('加入状态反馈后系统的阶跃响应')

24-grid

图28LQR控制器的求取程序

CommandWindow-□x

FileDebugDesktopWindowHelp青

i1NewtoMATLAB?WatchthisVideo,seeDemos,orreadGettncx

K=

-211.0180-0.00000.0013

A»l

OVR

图29LQR控制器的状态反馈矩阵k

同时在图28所示的程序中求取加入状态反馈后系统的伯德图，根轨迹和单位阶跃

响应,如图30,图31,图32所示。

加入状态反馈后系统的伯德图

!!!!!

1i•iiilliiiiinl•ei•iiml

11i•iieaitil1,aiiiiiia

aiia•iiaaait•iaan•tiit•taia•aaaii•iilliiaii•aataa•iiitiia

IiI•!!!!•■111Mill1•lI•eHI■(•••IlIiiiinai1••am・iiiiiII

aiiitii>>aat•itan(iai•iin••■■anii•••naii••ana•iaaiiia

1IliIi111•11I1IM1(iaaiiiiaiiiiti1•aIIMii11IIM