1 / 56

# 第 六 章 - PowerPoint PPT Presentation

I am the owner, or an agent authorized to act on behalf of the owner, of the copyrighted work described.

## PowerPoint Slideshow about ' 第 六 章' - odessa

An Image/Link below is provided (as is) to download presentation

Download Policy: Content on the Website is provided to you AS IS for your information and personal use and may not be sold / licensed / shared on other websites without getting consent from its author.While downloading, if for some reason you are not able to download a presentation, the publisher may have deleted the file from their server.

- - - - - - - - - - - - - - - - - - - - - - - - - - E N D - - - - - - - - - - - - - - - - - - - - - - - - - -
Presentation Transcript

MATLAB在控制系统中的应用

§6.1 控制系统的数学描述与建模

1.4 Ω

2H

0.32F

Vs=1V

ELECTSYS.M

function xdot=electsys(t,x)

v=1; R=1.4; L=2; C=0.32;

xdot=[x(2)/C;1/L*(V-x(1)-R*x(2))];

SIMU.M

t0=0;

tfinal=15; %时间间隔为0~15

x0=[0.5,0]; %初始化

tol=0.001; %精度

trace=0 %若非零则打印出每一步的计算值

[t,x]=ode23(\'electsys\',t0,tfinal,x0,tol,trace);

subplot(211)

plot(t,x)

title(\'Time Response of a RLC series ciruit\');

xlabel(\'Time-sec\');

text(8,1.05,\'Capatitor Voltage\');

text(8,0.05,\'current\');

vc=x(:,1);

i=x(:,2);

subplot(212);

plot(vc,i)

title(\'current versus capatitor voltage\');

xlabel(\'Capacitor Voltage\');

ylabel(\'current\');

subplot(111);

1. 多项式和特征多项式的根 root(P) poly(r) poly2sym(p)

P=[1 9 31.25 61.25 67.25 14.75 15]

r=roots(P)

r=[-1 –2 –3+j4 –3-j4]

p=poly(r)

poly2sym(p)

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

p=poly(A); r=roots(p); eig(A);

poly2sym(p)

2. 传递函数的零点和极点 tf2zp( )

num=[1 11 30 0];

den=[1 9 45 87 50];

[z,p,k]=tf2zp(unm,den)

z=[-6 –5 0]’; k=1; p=[-3+4*j,-3-4*j,-2,-1]’;

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

3. 部分分式展开 [r,p,k]=residue(b,a)

b=[2 0 9 1];

a=[1 1 4 4];

[r,p,k]=residue(b,a)

1. 将微分方程化为状态方程

2. 矩阵的对角化

1. 数学模型的转换（见下表： ）

(1) ss2tf—实现由状态空间描述到传递函数形式的转换。

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

a=[0 1 –1;-6 –11 6;-6 –11 5];

b=[0 0 1]’; c=[1 0 0]; d=0

[num,den]=ss2tf(a,b,c,d)

(2) ss2zp—由状态空间形式转化为零极点增益形式。

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

a=[0 1 –1;-6 –11 6;-6 –11 5];

b=[0 0 ;0 1;1 0]’; c=[1 0 0;0 1 0]; d=[2 0;0 2];

[num,den]=ss2zp(a,b,c,d,1);

[num,den]=ss2zp(a,b,c,d,2);

(3) tf2ss—由传递函数形式转化为状态空间形式。

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

num=[0 0 –2;0 –1 –5; 1 2 0]; den=[1 6 11 6];

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

2. 系统模型的连接（见下表： ）

（1）append( )—将两个状态空间描述的动态特性连接在一起。

[A,B,C,D]=append(A1,B1,C1,D1,A2,B2,C2,D2)

(2) series( )—将两个状态空间形式表示的系统进行串联。

[A,B,C,D]=series(A1,B1,C1,D1,A2,B2,C2,D2)

(3) parallel( )—两个状态空间系统的并联。

[A,B,C,D]=parallel(A1,B1,C1,D1,A2,B2,C2,D2)

[num,den]=parallel(num1,den1,num2,den2)

(4) feedback( )—用于两个系统的反馈连接。

[A,B,C,D]=feedback(A1,B1,C1,D1,A2,B2,C2,D2)

[A,B,C,D]=feedback(A1,B1,C1,D1,A2,B2,C2,D2,sign)

[num,den]=feedback(num1,den1,num2,den2)

[num,den]=parallel(num1,den1,num2,den2,sign)

u

y

g(s)

-

h(s)

numg=[2 5 1]; deng=[1 2 3];

numh=[5 10]; denh=[1 10];

[num,den]=feedback(numg,deng,numh,denh)

MATLAB提供了许多用来分析控制系统模型属性的函数，例如常

（1）gram( ),dgram( )—其可控、可观测的格来姆(Gram)矩阵。

Gc=gram(A,B,’c’) Gc=dgram(A,B,’c’)

Go= gram(A,B,’o’) Gc=dgram(A,B,’o’)

（2）ctrb( ),obsv( )—求系统的可控性和可观测性。

Ac=ctrb(A,B)

Ao=obsv(A,C)

CANDO.M

%判断系统在不同参数下的可控性和可观测性

for alpha=[-1 0 2]

alpha

num=[1 alpha];

den=[1 -15 14 0];

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

Ac=ctrb(A,B);

rAc=rank(Ac)

Ao=obsv(A,C);

rAo=rank(Ao)

end

(3) damp( ),ddamp( )—求系统的阻尼系数和自然频率。其用法如下：

[Wn,z]=damp(A)

mag=ddamp(A） %得到A的特征值向量。

[mag,Wn,z]=ddamp(A,Ts) %Ts为抽样时间间隔。

damp( ),ddamp( )函数用于计算自然频率和阻尼系数，当不带输出变量时，可在屏幕上显示出特征值表、阻尼系数和自然频率。A变量可以取三种形式，如下表：

A=[1 2 –1 2;2 6 3 0;4 7 –8 –5;7 2 1 6]

damp(A)

num=[3 –5.1 6];

den=[4 3.2 –1.9 8];

ddamp(den,0.01)

（4）printsys( )—显示或打印线性系统，其用法如下：

printsys(A,B,C,D)

printsys(A,B,C,D,ulabels,ylabels,xlabels)

printsys(unm,den,’s’)

printsys(unm,den,’z’)

B=[-1 0 0 1]’;

C=[-2 5 6 1];

D=8;

printsys(A,B,C,D)

or: printsys(A,B,C,D,’current’,’voltage’,’Vr I1 Vc I’)

or: [num,den]=ss2tf(A,B,C,D)

printsys(num,den)

(5) dcgain( ),ddcgain( )—计算控制系统的稳态增益，其用法如下：

k=dcgain(A,B,C,D); k= ddcgain(A,B,C,D);

k= dcgain(num,den); k= ddcgain(num,den);

B=[-1 0 0 1]’; C=[-2 5 6 1]; D=8;

k= dcgain(A,B,C,D)

1）are( )——求解里卡第(ATx+xB=-C)方程，其用法为：

x=are(A,B,C)

A=[1 2 –1 2;2 6 3 0;4 7 –8 –5;7 2 1 6];

B=eye(size(A));

C=[1 2 3 0;2 4 4 6;3 4 1 5;0 6 5 0];

x=are(A,B,C)

2）lyap( ),lyap2( ),dlyap( )——求解李亚普诺夫方程。其用法如下：

x=lyap(A,B,C) 可求解连续系统的一般形式的李亚普诺夫方程：

Ax+xB=-C (A,B,C应当有适当的维数)

x=lyap(A,C)可求解特殊形式的李亚普诺夫方程：

Ax+xAT=-C

x=dlyap(A,B,C)用于求解离散系统的李亚普诺夫方程：

lyap2( )与lyap( )的功能相同，但采用了特征值分解技术求解李亚普

A=[1 2 –1 2;2 6 3 0;4 7 –8 –5;7 2 1 6];

B=eye(size(A));

C=[1 2 3 0;2 4 4 6;3 4 1 5;0 6 5 0];

x=lyap(A,B,C)

§6.2 控制系统的分析方法

function sstable(A,B,C,D)

% 判断系统的稳定性以及系统是否为最小相位系统

[z,p,k]=ss2zp(A,B,C,D); %求出系统的零点、极点

ii=find(real(z)>0);

n1=length(ii);

jj=find(real(p)>0);

n2=length(jj);

if (n1>0)

disp(\'The System is Unstable.\');

else

disp(\'The System is stable.\');

if (n2>0)

disp(\'The System is a Nonminimal Phase One.\');

else

disp(\'The System is a Minimal Phase One.\');

end

end

function sstable(A,B,C,D)

% 判断系统(A,B,C,D)的稳定性以及系统是否为最小相位系统

if(nargin==2)

[z,p,k]=tf2zp(A,B) %求出系统的零点、极点

elseif(nargin==4)

[z,p,k]=ss2zp(A,B,C,D); %求出系统的零点、极点

end

ii=find(real(z)>0);

n1=length(ii);

jj=find(real(p)>0);

n2=length(jj);

if (n1>0)

disp(\'The System is Unstable.\');

else

disp(\'The System is stable.\');

if (n2>0)

disp(‘But the System is a Nonminimal Phase One.\');

else

disp(\'The System is a Minimal Phase One.\');

end

end

1. 时域分析的一般方法

step( )—求取系统的阶跃响应。其用法：

[y,x]=step(num,den,t)

[y,x]=step(A,B,C,D,iu,t)

t为选定的仿真时间，y为系统在仿真时刻各个输出所组成的矩阵，x为自动选择的状态变量的时间响应数据。

step(num,den,t) or step(num,den)

step(A,B,C,D,iu,t) or step(A,B,C,D,iu)

EXP4.M

num0=20;

den0=[1 8 36 40 0];

t=1:0.1:10;

numc=num0;

denc=[zeros(1,length(den0)-length(num0)),num0]+den0;

numc

denc

[y,T,x]=step(numc,denc,t);

plot(t,y,t,x);

title(\'The Step Response\')

xlabel(\'Time_Sce\');

text(4,1.5,\'The Output\');

text(5,6.5,\'The State\')

impulse( )—求取脉冲响应的函数，用法与step( )函数基本一致。

2. 常用的时域分析函数 （见下表）

(1) initial( )—求连续系统的零输入响应。用法如下：

[y,x,t]=iniial(A,B,C,D,x0)

[y,x,t]=iniial(A,B,C,D,x0,t)

A=[1 –1 0.5;2 –2 0.3;1 –4 –0.1];

B=[0 0 1]’; C=[0 0 1]; D=0;

x0=[1 0 0]’;

t=0:0.1:20;

initial(A,B,C,D,x0,t)

title(‘The Initial Condition Response’)

(2) dinitial( )—求离散系统的零输入响应。用法如下：

[y,x,t]=diniial(A,B,C,D,x0)

[y,x,t]=diniial(A,B,C,D,x0,n)

(3) lsim( )—对任意输入的连续系统进行仿真，用法如下：

[y,x]=lsim(A,B,C,D,u,t)

[y,x]=lsim(A,B,C,D,u,t,x0)

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

num=[1 6.8 13.85 8.05]

den=[1 11.2 46.4 88.4 77.4 25.2]

t=0:0.2*pi:2*pi

u=sin(t)

lsim(num,den,u,t)

3. 时域分析应用实例

ωn=4， ξ分别为0.1,0.2,…,1.0,2.0时系统的单位阶跃响应。

wn=4;

kosai=[0.1:0.1:1,2];

figure(1)

hold on

for i=kosai

num=wn.^2;

den=[1,2*i*wn,wn.^2]

step(num,den)

end

title(\'The Step Response of Two Order System\');

w=2:2:12;

kosai=0.6

figure(1)

hold on

for wn=w

num=wn.^2;

den=[1,2*kosai*wn,wn.^2]

step(num,den)

end

title(\'The Step Response of Two Order System\');

R(s)

C(s)

1+es

os=40;

tmax=0.8;

z=log(100/os)/sqrt(pi^2+(log(100/os))^2);

wn=pi/(tmax*sqrt(1-z^2));

num=wn^2;

den=[1 2*z*wn wn^2];

t=0:0.02:4;

c=step(num,den,t);

plot(t,c);

xlabel(\'Time-Sec\');

ylabel(\'y(t)\');

d=wn^2

e=(2*z*wn-1)/d

title(\'The Step Response of Two Order System\');

grid;

1. 频域分析的一般方法

1）对数频率特性图（波特图）

MATLAB提供了函数bode( )来绘制波特图，其用法如下：

[mag,phase]=bode(num,den,w)

wn=5;

kosai=[0.1:0.2:2];

w=logspace(-1,1,100) %生成对数横坐标，区间为0.1~10

num=[wn.^2];

for ii=kosai

den=[1 2*ii*wn wn.^2];

[mag,pha,w1]=bode(num,den,w);

subplot(2,1,1);

hold on

semilogx(w1,mag);

subplot(2,1,2);

hold on

semilogx(w1,pha);

end

subplot(2,1,1);

grid on

title(\'Bode Plot\');

ylabel(\'Gain dB\');

text(5.5,4.5,\'0.1\');

subplot(2,1,2);

grid on

ylabel(\'Phase deg\');

text(4,-20,\'0.1\');

text(2.5,-90,\'2.0\');

2） 极坐标图（奈奎斯特图）

MATLAB提供了函数nyquist来绘制系统的奈奎斯特图。

[re,im,w]=nyquist(num,den,w)

[re,im,w]=nyquist(A,B,C,D)

k1=1300;

k2=5200;

w=8:1:80;

num1=[k1];num2=[k2];

den=[1 52 100 0];

%subplot(2,1,1);

figure(1)

nyquist(num1,den,w);

axis([-1.0,0,-0.04,0.04]);

grid;

%subplot(2,1,2);

figure(2)

nyquist(num2,den,w)

grid;

3）频率响应

MATLAB提供了频率响应函数freqs( ),其用法如下：

y=freqs(num,den,w)

%Sysfreqs.M

num=4;

den=[1 2 4];

w=0:0.01:3

g=freqs(num,den,w);

mag=abs(g);

plot(w,mag);

ylabel(\'Magnitude\');

grid;

axis([0 3 0.5 1.2])

title(\'幅频特性\')

• 例1: 已知某开环系统如下所示：

%SysExp1.M

k=26;

z=[];p1=[1 -6];p2=[1 -6 2];

[num1,den1]=zp2tf(z,p1,k);

[num2,den2]=zp2tf(z,p2,k);

subplot(2,2,1);

nyquist(num1,den1);

subplot(2,2,2)

hold on

[numc1,denc1]=cloop(num1,den1);

impulse(numc1,denc1)

step(numc1,denc1)

subplot(2,2,3)

nyquist(num2,den2)

title(\'Nyquist Diagrams\');

subplot(2,2,4)

hold on

t=0:0.1:3

[numc2,denc2]=cloop(num2,den2);

impulse(numc2,denc2)

step(numc2,denc2,t);

axis([0,3,-100,20]);

%SysExp2.M

a=[-0.6 -1.044 0 0;1.044 0 0 0;0 0.9578 -0.7 -0.3162;0 0 0.3162 0];

b=[1;0;0;0];

c=[0 0 0 0.3162];d=0;

w=logspace(-1,1,100); %生成对数横坐标，区间为0.1~10

figure(1);bode(a,b,c,d);

figure(2);

nyquist(a,b,c,d);

figure(3);

[ac,bc,cc,dc]=cloop(a,b,c,d);

margin(a,b,c,d);

figure(4);

impulse(ac,bc,cc,dc);

deno=[1 36 205 750];

r=roots(deno)

-3.0000+4.0000i

-3.0000-4.0000i

%SysExp3.M

num1=25;den1=[1 6 25]; %近似的二阶系统

num2=750;den2=[1 36 205 750]; %原三阶系统

figure(1);

hold on;

bode(num1,den1,\'r\');

bode(num2,den2,\'g\'); %比较频率响应

figure(2);

hold on

step(num1,den1,\'r\');

step(num2,den2,\'g\'); %比较单位阶跃响应

MATLAB绘制轨迹的函数为rlocus(num,den,K)，num和den分别是系

%SysExp4.M

num=[1 5];den=[1 5 6 0];

rlocus(num,den);

[k,p]=rlocfind(num,den) %确定增益及相应的闭环极点

title(\'Root Locus\');

gtext(\'k=0.5\') %用鼠标标示文本

%SysExp5.M

num=[1 2];

den1=[1 4 3];

den=conv(den1,den1);

k=0:0.1:150;

figure(1)

rlocus(num,den,k);

title(\'Root Locus\');

[k,p]=rlocfind(num,den);

%检验系统的稳定性

figure(2)

k1=55;

num1=k1*num;

[numc1,denc1]=cloop(num1,den,-1);

impulse(numc1,denc1);

title(\'Impulse Response(K=55)\');

figure(3)

k2=56;

num2=k2*num;

[numc2,denc2]=cloop(num2,den,-1);

impulse(numc2,denc2);

title(\'Impulse Response(K=56)\');

• 学号未位数为0、5的同学做第1题，学号未位数为1、6的同学做第2题，学号未位数为2、7的同学做第3题，学号未位数为3、8的同学做第4题，学号未位数为4、9的同学做第5题。
• 1. 系统开环传递函数为：
• 1）试绘制系统的Bode, Nyquist, Nichols图，并确定系统的幅值裕量、相位裕量及它们相应的转折频率；
• 2）绘制 闭环系统的根轨迹，并确定临界稳定的K值；
• 3）绘制开环系统、闭环系统的阶跃响应和脉冲响应曲线；
• 4）结合上述要求进行理论分析。
• 单位反馈系统的开环传递函数为：
• 要求同上题。

• （1）
• （2）
• 的闭环系统的稳定性，对不稳定的系统通过什么措施使之稳定，并通过脉冲响应、阶跃响应加以验证。
• 已知系统开环传递函数为
• （1）求出闭环系统稳定的K值范围；
• （2）计算闭环主导极点具有ζ=0.5时的性能指标δ%，ts(0.05) ;
• （3）写出闭环传递函数近似为二阶系统的形式；
• （4）通过Bode图、阶跃响应、脉冲响应比较原系统与近似二阶系统的差异。
• 5. 已知单位反馈系统的开环传递函数分别为：

（1）

（2）

（3）