控制系统仿真的应用

MATLAB 在控制系统仿真中的应用

从本章开始，正式进入到有关MATLAB 在控制系统仿真的领域中去。我们知道，对控制系统进行仿真，首先要建立系统的数学模型，这是计算机仿真的基础。我们下面先介绍控制系统数学模型的描述方法。

第一节 控制系统数学模型的基本描述方法

在控制系统仿真中，主要用4种形式的数学模型：传递函数.、零极点模型、结构图形式和状态方程模型。这些模型之间存在着内在的等效关系。 在不同的场合下可能使用的模型形式要求不同，需要了解模型之间的转换方法。

这一节主要介绍它们的MATLAB 实现。

一．控制系统的传递函数描述：（在MATLAB 中称tf 模型）

对系统的微分方程在零初始条件下做拉氏变换，则可得系统的传递函数(SISO 系统)： G(s)=1

1211121)()(+-+-++++++=n n n m m m a s a s a b s b s b s u s y 对线性时不变（线性定常）系统（LTI ）来说，a 、b 均为常数a1≠0。前面讲过多项式的表示方法。这里分子分母都为多项式，可将分子分母分别表示出来。即用分子分母的系数构成两个向量，唯一的确定出来： num=[b 1,b 2,…,b 1+m ]

den=[a 1,a 2,…,a 1+m ]

注意：构成分子，分母向量按降幂排列的顺序。

这只是一个简单的特例，很多时候，传递函数的分子、分母均为多项式相乘的形式，如：

)

()())(()(5s 2s 3s 1s s 6s 6s 2s 4s G 2332

2+++++++= 不能直接写出，可借助多项式乘法运算函数conv( )来处理，以便获得分子、分母多项式向量。

conv( )函数的调用方式为：

c=conv(a,b)

a,b 各表示一个多项式，c 表示a 和b 的乘积多项式，此函数允许嵌套使用：c=conv(a,conv(b,e))，则上例的G(s)可用下面的语句来输入： num=4*conv([1,2],conv([1, 6, 6],[1, 6, 6]));

den=conv([1,0],conv([1,1],conv([1,1],conv([1,1],[1,3,2,5]))));

还可以进一步地编写一个convs( )函数来一次性地求出若干个（十个）多项式的连乘积：

function a=convs(a1,a2,a3,a4,a5,a6,a7,a,8,a9,a10)

a=a1;

for i=2:nargin %在MATLAB 中它也是一个固定变量，表示在函数调用时实际

输入变量的个数。

eval([‘a=conv(a,a’ int2str(i) ‘);’]) %eval( )函数执行MATLAB 语句构成的字

符串。

end

这里定义的函数，可以一次执行十个多项式的乘积。如上例可以写成： convs([1,0],[1,1],[1,1],[1,1],[1,3,2,5])

注意：写完函数后一定要按函数名来保存该文件。

相应地，对离散时间系统的动态模型，用脉冲传递函数描述，其输入方法与此类似，不多举例。

如果是MIMO 系统，则用传递矩阵描述，如：

????

??????++++=1s 12s 2s 1s s G 2)( 则可表示为：

num= {[1,1];1}

den={[1 2 2 ];[1 1]}

二．控制系统的零极点模型：（zpk 模型）

是传递函数的另一种表现形式，对SISO 系统： G(s)=)

())()(()())()((321321n m p s p s p s p s z s z s z s z s k ++++++++ 将零点、极点及K 值输入即可建立零极点模型。

z=[-z 1,-z 2 …,-z m ]

p=[-p 1,-p 2 …,-p n ]

k=k

对于给出的传函来说，分别对分子分母作因式分解，则可以得出系统的零极点模型。这可以通过求出分子，分母多项式的根来实现。MATLAB 中提供了多项式求根的函数，roots ( )

调用格式： z=roots(a)

其中：z — 各个根所构成的向量 a — 多项式系数向量

系统增益k 即为原传递函数分子的最高项系数与分母最高项系数的比值。 分别求出分子分母多项式的根，即可得到系统的零，极点模型。

对于多输入多输出系统，应分别对每个输入求出系统的零极点模型最后才可以获得整个系统的零极点模型（为矩阵形式）。

两种模型之间的转换函数：

[z ,p ,k]=tf2zp(num , den);

[num , den ]=zp2tf(z ,p ,k)

三. 控制系统状态方程模型：( ss 模型)

LTI 系统的状态方程： x

= Ax + Bu Y = Cx + Du

只要将A ，B ，C ，D 几个矩阵输入进去即可。

对于离散系统来说，也与上面类似。

MATLAB 还提供了由系统状态方程转换为系统零极点模型及传递函数模型的函数。

[num , den]=ss2tf(A ,B ,C ,D , iu ) %iu —表示输入的序号（对多输入系统）

[z, p ,k]=ss2zp(A ,B ,C ,D , iu ) %[z ,p ,k]— 表示对第iu 个输入信号的传递

函数的零极点。

[A , B ,C ,D]= tf2ss (num , den)

[A , B ,C ,D]=zp2ss(z ,p ,k)

四． 在MATLAB 控制系统工具箱中的LTI 对象：

为了避免对一个系统采用多个分离变量进行描述，新版本的控制系统工具箱，将LTI 系统的各种描述封装成一个对象，即用一个变量来描述。

在控制系统工具箱中，有以上讲述的三种对象，即ss 对象，tf 对象和zpk 对象。每种系统模型的生成和模型间的转换均可以通过一个函数来实现。下面介绍这些函数：

1．dss( )函数：生成系统的状态空间模型。

格式： sys= dss (a ,b ,c ,d ,e )

该调用生成连续系统的状态空间模型：

)(.

t x E =)()(t Bu t Ax +

)(t y =)()(t Du t Cx +

E 为非奇异阵（E 阵奇异，为奇异系统）。

sys = dss (a ,b , c, d ,e ,Ts )生成离散系统的状态空间模型。

)1(+n Ex =)()(n Bu n Ax +

)()()(n Dx n Cx n y +=

Ts 为采样周期

2．filt()函数：生成DSP （数字信号处理）形式的离散传递函数：

格式：sys=filt(num ,den) --- 生成离散传递函数模型

sys=filt(num ,den ，Ts)---定义模型的采样周期

注：DSP （数字信号处理）形式的离散传递函数如下面的形式： 211

1

z 3z 21z 1z H ----+++=)(

3．ss()函数：生成状态空间模型，或者将传递函数及零极点模型转换成状态

空间模型。

格式：1）sys=ss(a,b,c,d)---生成连续系统的状态空间模型。形式为：

Bu Ax x

+= y =Cx + Du

2）sys =ss(a,b,c,d,Ts)---生成离散状态空间模型

x(n+1)=Ax(n)+Bu(n)

y(n)=Cx(n)+Du(n)

3）sys_ss=ss(sys)---将任意的LTI 对象sys 转换成状态空间模型。

4．tf()函数：生成传递函数模型，或将零极点模型及状态空间模型转换成

传递函数模型。

格式为：sys=tf(num,den)---生成连续时间系统传函模型: G(s)=n

1n 2n 1n 1m 2m 1a s a s a b s b s b ++++-- sys=tf(num ,den, Ts)---生成离散时间系统传函。

tfsys=tf(sys)---将任意的LTI 对象转换成传递函数模型。

5．zpk()函数：生成零极点模型或者将其他模型转化成零极点模型。 格式：sys=zpk(z ,p ,k) --- 连续系统的零极点增益模型。

sys=zpk(z ,p ,k,Ts)--- 离散时间系统的零极点增益模型。

zsys=zpk(sys) --- 将任意LTI 对象转换成零极点增益模型。

五. 控制系统的结构图描述及转换：

由于复杂的结构图建模可以用SIMUKINK 。在这里我们仅介绍简单的结构图描述，也可称为系统建模。

典型连接：

1．串联连接：

21子系统1： x

1 =A 1x 1+B 1 u 1 子系统2： x 2=A 2x 2+B 2u

2 y 1= c 1x 1+D 1u 1 y 2= C 2x 2+D 2u 2

串联后系统模型为：

???

?????.2.1x x

=??????21210A C B A ??????21x x +??????121D B B 1u y 2=[]212C C D ??

????21x x +D 2D 1u 1

2．并联连接：

u u = y y

并联后系统模型：

??????21x x =??????2100A A ??????21x x +??

????21B B u 1

y=[]21C C ??

????21x x +(1D +2D )u

3．负反馈连接：

r

1u 1y y

为了简化,设021==D D ，则反馈后：

1x

=1A 1x +1B 1u =1A 1x +1B (r-2y ) =1A 1x +1B r-1B (2C 2x +2D 2u )

=1A 1x +1B r-1B 2C 2x -1B 2D 1y

=r B x C B x C D B A 122111211+--)(

1122222222x C B x A u B x A x

+=+= 所以，状态方程为：

r 0B x x A C B C B C D B A x x 12121221121121??

????+????????????--=??????

在MATLAB 工具箱中，提供了子系统的连接处理函数：

1）series()函数：系统串联实现

格式：sys=series(1sys ,2sys )

[A,B,C,D]=series(A 1,B 1,C 1,D 1, A 2,B 2,C 2,D 2)

2）parallel()函数：系统并联实现。

格式：sys=parallel(1sys ,2sys )

3） feedback()函数：系统反馈连接

格式：sys=feedback(1sys ,2sys )

sys=feedback(1sys ,2sys ,sign)

sign —定义反馈形式：正反馈，sign=+1；负反馈，sign=-1。

4）append 函数：多个LTI 系统的组合

格式:sys=append(1sys ,2sys ,…,n sys )

1u 1y 2u 2y

N u N y

设各个子系统的模型为：1H (s), 2H (s),… n H (s)则合成后子系统为

分块对角阵： ???????????

?)()()(21s H s H s H n

如子系统为状态空间模型，则组合后：

??????????=N A A A 1 ????

??????=N B B B 1 =C ????

??????N C C 1 ??????????=N D D D 1

第二节 控制系统的稳定性分析

在计算机仿真中，控制系统稳定性的分析方法有：

① 求闭环特征方程的根；

② 化为零极点模型，看极点是否在s 右半平面；

③ 对状态空间形式（闭环），求A 阵的特征值eig(A)；

④ 用Lyapunov 方程求解。

在自控原理中，我们均采用间接方法劳斯判据，奈氏判据等，但由于在MATLAB 中很容易地求解多项式方程。因此，我们可以直接求出特征方程的根，可用多项式求根函数roots()。

例1：判断如下系统的稳定性：

可编程如下：numg=1; deng=[1 1 2 23];

numf=1; denf=1;

(num,den)= feedback(numg,deng,numf,denf,-1);

roots(den)

例2：判断状态空间形式描述的系统的稳定性，两个输入。

A=[0 1 –1；-6 –11 6；-6 –11 5]；

B=[0 0 1]；

C=[1 0 0]；

D=0；

[1z ,1p ,1k ]=ss2zp(A,B,C,D,1)

[2z ,2p ,2k ]=ss2zp(A,B,C,D,2)

ii=find(real(1p )>0),1n =length(ii);

%find （ ）：找出1p 数组中满足实部值大于0的所有元素的下标，并将结

果返回到ii 数组中。

ii1= find(real(2p )>0), 2n =length(ii1);

if (1n +2n >0),disp(‘system is Unstable ’)

else

disp(‘system is Stable ’);

end

一般当系统中含有非线性环节时，即不能用上述方法，一般需求解李雅普诺夫方程：T A P+PA=-Q 求出P 为正定对称矩阵，则系统是稳定的。

李雅普诺夫（Lyapunov ）方程求解函数Lyap()（采用Schur 方法编写的），其调用格式：

x=lyap(A,B,C)

该函数可求解：AX+XB=-C 方程

如要求解：T A P+PA=-Q 则可写成：p=lyap(A,Q)

注：判断稳定性的方法还可确定是否为最小相位系统（没有右半平面零极点的系统）。

第三节 控制系统的根轨迹

在控制系统工具箱中提供了系统根轨迹绘制及分析函数，还有一个GUI 分析工具。下面介绍一下这些根轨迹的绘制及分析函数：

1 rlocfind():计算给定根的根轨迹增益

格式：①[k,poles]=rlocfind(sys)

②[k,poles]= rlocfind(sys,p)

① 在LTI 对象的根轨迹图中显示出十字光标，当用户选择其中一点时，其

相应的增益由k 记录，与增益相关的所有极点记录poles 中。若要使用该函数，必须首先在当前窗口上绘制系统的根轨迹。

② 定义要得到增益的根矢量P ，即事先给出极点。除了显示出该根对应的

增益以外，还显示出该增益对应的其它根。

2 rlocus() 函数：功能为求系统根轨迹。

格式：rlocus(sys)--- 计算SISO 开环LTI 对象的根轨迹，增益自动选

取。

rlocus(sys,k)---显式设置增益。

[r,k]= rlocus(sys)--- 返回系统的增益k （向量）和闭环极点r

（向量），即对应于增益的闭环极点为r(i)。

3 sgrid()函数：绘制连续时间系统根轨迹和零极点图中的阻尼系数和自然频

率栅格线。（在根轨迹图中加ξ和n ω网格）

格式：sgrid---在连续系统根轨迹或零极点图中绘制栅格线。

栅格线由等阻尼系数和等自然频率线构成，阻尼系数ξ步长

为0.1范围从0到1，自然频率n ω步长为1弧度/秒，范围从

0到10。在绘制前当前窗口必须包含根轨迹图。

sgrid(z,n ω)---自己定义想要绘制的阻尼系数向量z 和自然频率向量

n ω的范围。

注：在画出的根轨迹图上面点击，可以给出对应于所点击点的极点值、对

应增益、阻尼比和超调量等。

4 zgrid()函数：绘制离散时间系统根轨迹和零极点图中的栅格。

格式：zgrid--- 在离散系统的根轨迹图或零极点图上绘制出栅格线。 zgrid(z ，n ω)---显式定义想要绘制的z 向量和n ω向量。

第四节 系统的频率响应

一．频率响应的计算：

若系统的传递函数为： G(s)= 1

11211121+-+-++++++++n n n m m m m a s a s a s a b s b s b s b 则频率特性函数

G （j ω）=1

1211121)()()()()()(+-+-++++++++n n n n m m m m a j a j a j a b j b j b j b ωωωωωω 若系统由状态空间模型来描述，则系统的频率特性为：

G(j ω)=C[j ωI+A]1-B+D

频率响应的求法：

1）对传函形式可用：

GWI=polyval(num,sqrt(-1)*w)./polyval(den,sqrt(-1)*w)

Polyval （）：多项式求值函数

sqrt(-1)*w ：构造j ω向量

2）对状态空间模型形式可用：evalfr()和freqresp()函数来求

evalft()格式：frsp=evalfr(sys ，f)

f ：复数频率点，f 值以复数给出

功能：计算系统在单个频率点的频率响应

freqresp()功能：计算系统频率响应

格式：H=freqresp(sys,ω)

H ：是三维向量

ω：实数频率向量，可指定一些实数频率点

如：计算下述系统在频率1，10，100处的频率响应。

H(s)=????

??????+-+121110s s s num={0,1;[1,-1],1};

den={1,[1 1];[1 2],1};

p=tf(num,den)

w=[1 10 100];

H=freqresp(p,w)

MATLAB 将返回：

0000

.16000.02000.05000.05000.00i i +--

0000

.12885.09423.000990.00099.00i i +- 对第2个频点

0000

.10300.09994.00100.00001.00i i +-

二．频率响应曲线的绘制：

提供了3种频率响应曲线的绘制函数：Bode 图绘制，Nyquist 曲线及Nichols 曲线的绘制，它们可用于SISO 或MIMO 的连续时间和离散时间系统。

1. Bode 图绘制：bode()函数

调用格式：①bode(sys): bode(num,den)

bode(A,B,C,D,iu) ---可自动地选择一个合适

的频率范围。

②bode(sys ，ω)---给出频率范围，这里频率范围一般由

ω=logspace(a,b,n)给出。

logspace(a,b,n)：表示在a 10 到b 10之间的

n 个点。

这两种格式可直接画出规范化的图形（有网格，有说明）。

③[mag,phase,ω]=bode(sys)或[m,p]=bode(sys)

这种格式只计算Bode 图的幅值向量和相位向量，不画出图形。

要在此基础上画图，可用：

subplot(211);semilogx(ω,20*log10(m) %对数幅频曲线

subplot(212);semilogx(ω,p) %对数相频曲线

④bode(1sys ,2sys ,…,N sys )

⑤bode(1sys ,2sys ,…,N sys ,ω)

这两种格式可在一个图形窗口同时绘多个系统的bode 图。

2. Nyquist 曲线的绘制：nyquist()函数

格式：①nyquist(sys)

② nyquist(sys,ω)

③nyquist(1sys ,2sys ,…,N sys )

④nyquist(1sys ,2sys ,…,N sys ,ω)

⑤[re,im,ω]=nyquist(sys) re —频率响应实部 im —虚部

3. Nichols(sys)曲线：nichols()函数

格式：①nichols(sys).

②nichols(sys,ω).

③nichols(1sys ,2sys ,…,N sys )

④nichols(1sys ,2sys ,…,N sys ,ω)

⑤[mag,phase,ω]=nichols(sys,ω)

= nichols(sys)

mag ：幅值 phase ：相位

四．计算控制系统的稳定裕度：margin( )函数

利用该函数可以直接求出系统的幅值裕度和相角裕度。

调用格式为：

①[m G ，m P ，cg ω，CP ω]= margin(num,den)

= margin(A,B,C,D)

= margin(sys)

m G --- 幅值裕度；m P ---相位裕度；cg ω---幅值裕度处对应的频率g ω； cp ω---相位裕度处对应的频率c ω。

② [m G ，m P ，cg ω，CP ω]=margin(mag ,phase,ω)---根据幅值、相位数

据计算稳定裕度。

③ margin(sys)---在当前图形窗口中绘制出系统裕度的Bode 图。

例：对本节内容的综合使用程序：

A=[0 1 –1；-6 –11 6；-6 –11 5]；

B=[0 0 1]；

C=[1 0 0]；

D=0；

ω=logspace(-1,1)；

[m,p]=bode(A,B,C,D,1,ω)

subplot(211),semilogx(ω,20*log10(m))；

subplot(212),semilogx(ω,p)；

[X,Y]=nyquist(A,B,C,D,1,ω);

plot(X,Y)

[m G ，m P ，cg ω，CP ω]=margin(A,B,C,D)

第五节 线性系统时间响应分析

系统仿真的实质就是对描述系统的数学模型进行求解，一般的求解方法为微分方程的数值解法（龙格库塔法ode23( ),ode45( )等），这是系统仿真的基本方法。

但我们没有必要直接采用上面的方法，因为工具箱已提供了直接求各种输入下系统响应的函数。

一般控制系统的输入有：阶跃，斜坡，加速度及脉冲输入等。在工具箱主要提供了如下三种函数：

1．step( )函数---求系统阶跃响应

其调用格式为：1）)(sys step ),(d e n n u m s t e p

=> ),,,,(iu D C B A step =>

2）),(t sys step ------定义仿真时间

这两种格式直接画出响应曲线

3）[y,x]=step(sys)---只计算仿真数据，不画图。

=step(sys,t)

y ---各个仿真时刻的输出向量。

x ---自动选择的状态向量的时间响应数据。

4）step(1sys ,2sys …)—同时仿真多个系统。

2．impulse( )函数：求取系统的脉冲响应

格式与step 函数相同。

3．lsim( )函数：求系统的任意输入下的仿真

格式：1）lism(sys ,u ,t)

u---给定输入构成的列向量，它的元素个数应与t 的元素个数

相一致。

这种形式还可以写成：

lsim(num ,den ,u ,t)

lsim(A ,B ,C ,D ,iu ,u ,t)

其中：iu---输入变量序号

2）lsim(sys ,u ,t,0x )---定义初始状态0x

3）[y,x]=lsim(sys ,u ,t ,0x )---不画出图形只算数据

4）lsim(1sys ,2sys ,…,N sys ,μ,t)---同时仿真多个LTI 对象

例：求下面系统在阶跃信号为0.011(t)时系统的响应。

num=20；

den=[1 8 36 40 20];

t=0:0.1:20;

1t =length(t);

u=ones(1t ,1)*0.01;

[y ，x]=lsim(num ,den ,u ,t)

plot(t,y)

第六节 延迟系统的仿真

带有时间延迟的控制系统传函模型 G(s)= 1

11211121+-+-++++++++n n n m m m m a s a s a s a b s b s b s b Ts e -=Ts e s G -)(? T 为延迟时间常数。

纯时间滞后环节可以由有理函数来近似，1892年法国数字家Pade 曾提出了一种著名的有理近似方法，后人将它命名为Pade 近似法。其表达式为：

Ts e -=

+++++++-+-+-5443322154433221)()()()(2/1)()()()(2/1TS P TS P TS P TS P TS TS P TS P TS P TS P TS 求这样的近似有理函数，控制系统工具箱提供了一个函数pade( )。 pade( )：计算pade 近似的系数。

它的调用格式为：

[num ,den]= pade(T ,n) --------- 返回近似模型

T---延迟时间常数 n---要求拟和的阶数

或：

[A ，B ，C ，D]=pade(T ,n)---------返回等效状态方程模型

拟和的阶数越高，精度也越高。

另，建立延迟系统的数学模型还有两种方法：

1）用函数TF 、ZPK 或SS 建模并同时给模型时间延迟属性“IoDelayMatrix 、

InputDelay 、OutputDelay ”赋值。

2）先用函数TF 、ZPK 或SS 建模，再用set 命令给模型时间延迟属性赋值。 例：已知系统传递函数为s e s s H 1.01

2)(-+= 用MATLAB 生成系统模型程序及运行结果如下：

%Create system model

Mun=2;

Den=[1 1];

sys=tf(num,den,’InputDelay ’,0.1)

>>Transfer function: exp(-0.1*s)*1

2+s 或：

%Create system model

Mun=2; Den=[1 1];

sys=tf(num,den)

set(sys, ’InputDelay ’,0.1)

第七节 状态空间法仿真实现

在控制系统工具箱中，有各种状态空间法实现的函数：

1．ss2ss( )：功能是进行相似变换

如有状态方程：???Du

Cx y Bu Ax x +=+=

引入非奇异变换矩阵T ，使变换后的状态变量x (t)=Tx(t)

则变换后系统的状态方程： ???u

D x C y u B x A x +=+=

其中： A =1TAT - TB B = 1-=CT C D D =

该变换为相似变换。

该函数的格式为：

sysT=ss2ss(sys,T)---完成状态空间中LTI 对象sys 的相似变换,T 为

变换矩阵

2．canon( )函数：功能是状态空间的正则实现

该函数计算连续时间或者离散时间LTI 对象的正则实现，它支持两种正则形式：模态矩阵(modal matrix)和伴随矩阵(companion matrix)形式。

格式：1）csys=canon(sys ,‘type ’)

‘type ’可以为以下字符串：‘modal ’---计算模态矩阵正则实现 ‘c ompanion ’---计算伴随矩阵正则实现 sys---如果系统不是以状态空间模型给出的，要首先进行模型转换。

2）[csys ,T]=canon(sys ,‘type’)

T---返回相似变换矩阵T ，当给出的LTI 对象不是ss 对象时T 为

是空矩阵。

注：模态实现要求系统的A 矩阵为可对角化的。

3．ctrb()函数：功能是计算可控性矩阵 Co=[B AB 2A B …1-n A ]

格式：Co=ctrb(A,B)

Co=ctrb(sys)

这里sys 为状态空间对象

如果rank(Co)=n 则系统可控

4．ctrbf()函数：系统的可控与不可控分解

格式：[Abar ,Bbar ,Cbar ,T,k]=ctrbf(A,B,C)

T---相似变换阵

k---是长度为n 的一个矢量，其元素为各个块的矩阵

Abar ,Bbar ,Cbar---对应于转换后系统的（A ，B ，C ）,通过它

们看出可控与不可控部分。

5. obsv() 函数：计算可观性矩阵

格式：ob=obsv(A,C)

ob=obsv(sys)

如果rank(ob)=n ,则系统可观。

6. obsvf():系统的可观与不可观分解

格式：[Abar ,Bbar ,Cbar ,T,k]=obsvf(A,B,C)

与ctrbf()函数相似。

7. minreal():状态空间的最小实现

格式：sysr=minreal(sys)---实现LTI 对象的最小实现，即删除状态空

间模型中的不可观和不可控状态，或者对

消零，极点模型中相同的零极点对。

sysr 具有最小的阶数。

如可给零极点模型 )(s G =)

1)(1()1(2++--s s s s s 用该函数后可返回：1

)(2++=s s s s G

8．模型降阶函数：

格式：1）?

??==)''lim,,(mod lim),(mod mde e sys red rsys e sys red rsys 降阶后系统rsys 具有与原系统相匹配的直流增益(DC gain)

elim---给出需要删除的状态的索引

2）rsys=modred(sys ,elim ,’del ’)

仅仅简单删除要删除的状态，不能保证匹配的直流增益。

9．SISO 系统的极点配置函数：acker()

格式：k=acker(a ,b ,p)--- 利用Ackermann 公式计算反馈增益矩阵

k,使u=-kx 的单输入系统具有指定的闭环

极点p 。

a,b---p=eig(a-bk)

10．MIMO 极点配置函数：place()

格式：k=place(A,B,p)----p=eig(A-Bk)

[k ,prec ,message]=place(A,B,p)

prec---系统闭环实际极点与希望极点的接近程度。

message---如果系统闭环的实际极点偏离希望极点10%以上，将给

出警告信息。

place 函数还可用来计算估计器增益：L=place(A ’,C ’,p)

11.求系统状态估计器函数：estim()

格式：①est=estim(sys,L)---生成状态和输出估计器

L--- 估计器增益矩阵，由函数形成或由滤波函数生成。 sys--- 状态空间模型

②est=estim(sys ,L ,sensors ,known)---生成更为普遍的状态和

输入估计器

sensors---指定可观输出； known---指定已知输入。

12．系统调节器生成函数： reg( )

格式：),,(L k sys reg rsys =---生成给定状态估计增益L 及状态反馈增益

阵k 下的sys 系统的调节器或补偿器est 。

sys 为状态空间模型.

),,,,,(controls

known sensors L k sys reg rsys = ---生成更为普遍的LTI 系统sys 的状态和输出调节器

sensors ---可测输出；known ---已知输入；controls ---控制输入

最优控制仿真

LQ 和LQG ：

LQ---线性二次型最优控制问题

LQG---线性二次型Gauss 最优控制问题（对系统中含有系统噪声和

量测噪声的系统的最优控制问题）。

1． ()lqr 函数：连续系统的LQ 调节器设计

格式：],,[e s k =),,,(R Q B A lqr 或 ],,[e s k =),,,,(N R Q B A lqr k ---状态反馈阵

s ---Riccati 方程的解

e ---闭环系统特征值

Q ，R---性能指标dt Nu x 2Ru u Qx x J 0T T T )(?∞

++=

2．()dlqr 函数：离散系统的LQ 调节器设计

格式：[k,s,e]=slqr(a,b,Q,R)

[k,s,e]=slqr(a,b,Q,R,N)

3. ()lqry 函数：系统的LQ 调节器设计

格式：[k,s,e]=lqry(sys,Q,R)

[k,s,e]=lqry(sys ,Q ,R ，N)

Q,R---性能指标：?∞

++=0T t T dt Nu y 2Ru u Qy y J )(

4．lqrd()：连续时间系统的离散LQ 调节器设计

格式：[d K ，s ,e] =lqrd(A,B,Q,R,s T )

LQG 控制问题可以分解为两个问题求解：

1）利用Kalman 滤波理论，从状态)(t x 中得到最优估计。

2）利用估计状态x

?代替x ，设计满足二次型性能指标的LQ 控制律。