关长亮最优控制实验报告册

专    业:控制理论与控制工程

姓    名:   关长亮        

学    号:    2011270      

实验一 线性二次型最优控制问题

1、设二阶系统为， 

性能指标(注性能指标中加u（t）的平方项)，求最优控制和最优性能指标，并进一步画出闭环系统的阶跃响应曲线。

syms x1 x2;

x=[x1;x2];

A=[-1 0;1 0];

B=[1;0];

Q=[1 0;1 0];

R=1/2;

[K,P,E]=lqr(A,B,Q,R)

figure('pos',[50,50,200,150],'color','w');

axes('pos',[0.2,0.2,0.72,0.72])

C=[1,0];

D=0;

k1=K(1); 

Ac=A-B*K;Bc=B*k1;Cc=C;Dc=D; 

Step(Ac,Bc,Cc,Dc)

u=-inv(R)*B'*P*x

x0=[0;1];

J=1/2*x0'*P*x0

d=eig(A-B*inv(R)*B'*P)

运行结果如下：

K =

    0.7321    0.0000

P =

    0.3660    0.0000

    0.0000   -0.5000

E =

   -1.7321

   -0.0000

u =

- 2^(1/2)*x1 - 2^(1/2)*x2

J =

    0.8536

d =

   -1.4142

   -1.0000

有仿真得，

最优控制率为：

最优性能指标J=0.8536

闭环系统的阶跃响应图像为

2、考虑系统，其中

性能指标为，其中。

设计最优控制器，并求出黎卡提方程的解及闭环系统的特征值。

程序如下：syms x1 x2 x3;

x=[x1;x2;x3];

A=[0 1 0;0 0 1;-35 -27 -7];

B=[0;0;1];

R=1;

Q=[1 0 0;0 1 0;0 0 1 ];

N=[0;0;0];

[K,P,E]=lqr(A,B,Q,R,N)

u=-inv(R)*B'*P*x

d=eig(A-B*inv(R)*B'*P)

运行结果如下：

K =

    0.0143    0.1396    0.0908

P =

    5.2732    3.2786    0.0143

    3.2786    3.4266    0.1396

    0.0143    0.1396    0.0908

E =

  -2.4940 + 3.2297i

  -2.4940 - 3.2297i

  -2.1029          

u =

- (4116736823184415*x1)/288230376151711744 - (6287410398885*x2)/4503599627370496 - (6541462181055331*x3)/72057594037927936

d =

  -2.1029          

  -2.4940 + 3.2297i

  -2.4940 - 3.2297i

选用状态反馈矩阵为K = [0.0143  0.1396  0.0908]；

黎卡提方程的解是：P =

    5.2732    3.2786    0.0143

    3.2786    3.4266    0.1396

    0.0143    0.1396    0.0908

闭环系统特征值是：

d =

  -2.1029          

  -2.4940 + 3.2297i

  -2.4940 - 3.2297i

3、已知动态系统，性能指标为，其中。试计算最优状态反馈矩阵，并求其闭环系统的单位阶跃响应曲线。

程序如下：

syms x1 x2 x3;

x=[x1;x2;x3];

A=[0 1 0;0 0 1;-6 -11 -6];

B=[0;0;1];

Q=[100 0 0;0 1 0;0 0 1];

R=1;

[K,P,E]=lqr(A,B,Q,R)

figure('pos',[50,50,200,150],'color','w');

axes('pos',[0.2,0.2,0.72,0.72])

C=[1,0,0];

D=0;

k1=K(1); 

Ac=A-B*K;Bc=B*k1;Cc=C;Dc=D; 

Step(Ac,Bc,Cc,Dc)

运行结果如下：

K =

   26.1870   12.61    1.81

P =

  694.2199  217.9258   26.1870

  217.9258   94.1446   12.61

   26.1870   12.61    1.81

E =

  -3.9926          

  -1.9483 + 2.0654i

  -1.9483 - 2.0654i

最优反馈矩阵K=[26.1870   12.61    1.81]

4、已知动态系统，性能指标为，其中。试

（1）计算最优状态反馈矩阵；

（2）求其闭环系统的单位阶跃响应曲线以及状态的单位阶跃响应曲线；

（3）通过选取Q、R的不同，研究Q、R的不同对系统动态特性的影响。

（1）在q11=1，q21=2，q33=3下

最优控制矩阵K =

   10.0000   22.0801   15.7926

（2）编写程序如下：

syms x1 x2 x3 ;

x=[x1;x2;x3];

A=[0 1 0;0 0 1;0 -2 -3];

B=[0;0;1];

R=0.01;

Q=[1 0 0;0 2 0;0 0 3];

[K,P,E]=lqr(A,B,Q,R)

u=-inv(R)*B'*P*x

ap=[A-B*K];

bp=B;

C=[1,0,0];

D=0;

[ap,bp,cp,dp]=augstate(ap,bp,C,D);

cp=[cp;-K];

dp=[dp;0];

G=ss(ap,bp,cp,dp);

[y,t,x]=step(G);

plotyy(t,y(:,2:3),t,y(:,4))

（3）选用在q11=1，q21=2，q33=3 R=0.1（R不同）下

选用q11=q12=q13=40 R=0.01（Q不同）时，

推断Q的变化对快速性没有多大影响。到达0.98稳定区域用时更长，可以推断，R越小对快速性越好。

5、 线性系统为：                        ，

其目标函数是：

确定最优控制。

程序如下：

syms x1 x2;

x=[x1;x2];

A=[0 1;-5 -3];

B=[0;1];

R=1.6667;

Q=[500 200;200 100];

[K,P,E]=lqr(A,B,Q,R)

u=-inv(R)*B'*P*x

运行结果：

K =

   13.0276    6.7496

P =

   67.9406   21.7131

   21.7131   11.2495

E =

   -7.2698

   -2.4798

U=(3666940583351395*x1)/281474976710656-(75993324213883*x2)/11259906842624

6、 无人飞行器的最优高度控制，飞行器的控制方程如下 

是飞行器的高度；    是油门输入；设计控制律使得如下指标最小

初始状态                          。绘制系统状态与控制输入，对如下给定的     矩阵进行仿真分析

（1）

（2）

（3）

（4）

程序：syms x1 x2 x3 ;

x=[x1;x2;x3];

A=[0 1 0;0 0 1;0 0 -0.5];

B=[0;0;0.5];

R=200;

Q=[1 0 0;0 0 0;0 0 0];

[K,P,E]=lqr(A,B,Q,R)

u=-inv(R)*B'*P*x

ap=[A-B*K];

bp=B;

C=[1,0,0];

D=0;

[ap,bp,cp,dp]=augstate(ap,bp,C,D);

cp=[cp;-K];

dp=[dp;0];

G=ss(ap,bp,cp,dp);

x0=[10;0;0]

[y,t,x]=initial(G,x0);

plot(t,x);

hold on;

plot(t,y);

(1)

(2)

(3)

(4)

实验二：单级倒立摆稳定控制实验

    系统开环特征值为

eig(A)=

         0

   -0.0830

   -5.2780

    5.2727

特征根在复数平面的右半面上有一个零点，因此该开环系统不稳定。下面对其进行最优控制算法设计。

实验程序：

M = 1.096; 

b=0.1;

m = 0.109; 

l = 0.25; 

I=0.0034;

g = 9.8; 

T = 0.005;   

p = m+4*M2; 

A = [ 0 1 0 0 ; 

0 -(I+m*l^2)*b/[I*(M+m)+M*m*l^2]  m^2*g*l^2/[I*(M+m)+M*m*l^2] 0; 

0 0 0 1; 

0 -m*l*b/[I*(M+m)+M*m*l^2] m*l*g*(M+m)/[I*(M+m)+M*m*l^2] 0]; 

B = [0;(I+m*l^2)/[I*(M+m)+M*m*l^2];0;m*l/[I*(M+m)+M*m*l^2]];

C= [1 0 0 0;0 0 1 0]

eig(A)

Q11=1;Q22=1;Q33=1;Q44=1; 

Q=[Q11 0 0 0; 

0 Q22 0 0; 

0 0 Q33 0; 

0 0 0 Q44]; 

R=2; 

K = lqr(A,B,Q,R) 

Ac = [(A-B*K)]; 

Bc = [B]; 

Cc = [C]; 

Dc = [D]; 

T = 0:0.005:10; 

U = 0.1*ones(size(T)); 

[Y,X] = lsim(Ac,Bc,Cc,Dc,U,T); 

plot(T,Y) 

legend('Cart1 Position','Cart2 Position','Pendulum Angle') 

grid

LQR控制的阶跃响应如图所示

其中，红线代表摆杆角度，蓝线代表主动小车位移，绿线代表从动小车位移。从图中可以看出，响应的超调量很小，但稳定时间和上升时间偏大，小车的位置没有跟踪输入，而是向相反方向移动。在这里，我们首先来缩短稳定时间和上升时间。

可以发现，矩阵中，增加和使稳定时间和上升时间变短，并且使摆杆的角度变化减小。这里取，，则响应曲线如下：

如果再增大，，系统的响应还会改善。但在保证，和足够小的情况下，系统响应已经满足要求了。