自控课程设计报告

成绩

课程设计报告

题目控制系统的设计与校正

课程名称自动控制原理课程设计

院部名称机电工程学院

专业电气工程及其自动化

班级 10电气(1) 学生姓名董天宠

学号 1004103037 课程设计地点 C306 课程设计学时 1周

指导教师陈丽换

金陵科技学院教务处制

目录

一、设计目的 (3)

二、设计任务与要求 (3)

三、设计方案 (4)

四、校正函数的设计 (4)

4.1、校正前系统特性 (4)

4.2、利用MATLAB语言计算出超前校正器的传递函数 (6)

4.3校验系统校正后系统是否满足题目要求 (7)

五、函数特征根的计算 (8)

5.1校正前 (8)

5.2校正后 (9)

六、系统动态性能分析 (10)

6.1 校正前单位阶跃响应 (10)

6.2校正前单位脉冲响应 (11)

6.3校正前单位斜坡信号 (14)

七、校正后动态性能分析 (14)

7.1 校正后单位阶跃响应 (15)

7.2 校正后单位冲击响应 (15)

7.3 校正后单位斜坡响应 (16)

八、系统的根轨迹分析 (17)

8.1、校正前根轨迹分析 (17)

8.2、校正后根轨迹分析 (19)

九、系统的奈奎斯特曲线分析 (21)

9.1校正前奈奎斯特曲线分析 (21)

9.2 校正后奈奎斯特曲线分析 (22)

设计小结 (23)

参考文献 (24)

1.设计目的

1）掌握自动控制原理的时域分析法，根轨迹法，频域分析法，以及各种补偿（校正）装置的作用及用法，能够利用不同的分析法对给定系统进行性能分析，能根据不同的系统性

能指标要求进行合理的系统设计，并调试满足系统的指标。 2）学会使用MATLAB 语言及Simulink 动态仿真工具进行系统仿真与调试。

2.设计任务与要求

已知单位负反馈系统的开环传递函数0

K G(S)S(0.1S 1)(0.001S 1)

=

++，

试用频率法设计串联超前校正装置，使系统的相位裕度

045γ≥，静态速度误差系数1v K 1000s -=

1）首先， 根据给定的性能指标选择合适的校正方式对原系统进行校正，使其满足工作要求。要求程序执行的结果中有校正装置传递函数和校正后系统开环传递函数，校正装置的参数T ，α等的值。

2)利用MATLAB 函数求出校正前与校正后系统的特征根，并

判断其系统是否稳定，为什么? 3)利用MATLAB 作出系统校正前与校正后的单位脉冲响应曲线，单位阶跃响应曲线，单位斜坡响应曲线，分析这三种曲线的关系？求出系统校正前与校正后的动态性能指标

σ%、tr 、tp 、ts 以及稳态误差的值，并分析其有何变化？

4)绘制系统校正前与校正后的根轨迹图，并求其分离点、汇合点及与虚轴交点的坐标和相应点的增益K *值，得出系统稳定时增益K *的变化范围。绘制系统校正前与校正后的Nyquist 图，判断系统的稳定性，并说明理由？ 5)绘制系统校正前与校正后的Bode 图，计算系统的幅值裕量，相位裕量，幅值穿越频率和相位穿越频率。判断系统的稳定性，并说明理由？ 3设计方案

串联超前校正装置的传递函数为

1()(0)1Ts

Gc s Ts αα+=

>+

1根据对稳态误差的要求确定开环增益。

2按已确定的K ，绘制未校正系统的伯德图，并计算相角裕度0γ。

3确定需要补偿的相位超前角0~m ?γγ=-+?（510）。

4利用Matlab 语言计算出超前校正器的传递函数。 5校验校正后的系统是否满足题目要求。

4.校正函数的设计

.1 校正前系统特性 原传递函数0

K G(S)S(0.1S 1)(0.001S 1)

=

++具有一个积分环节，所以原

系统为I 型系统，则其速度误差系数Kv K =,使系统的速度误差系数

,即K=1000。

原传递函数0

K G(S)S(0.1S 1)(0.001S 1)

=

++

根据原系统的开环对数幅频特性的剪切频率c=100rad/s ，

求出原系统的相角裕度

约为0度，这说明原系统在

K=1000/s 时处于临界稳定状态，不能满足45度的要求。

为满足

45度的要求，串联校正装置提供的最大超前相角

m 必须大于等于45度。考虑到校正后系统的剪切频率c'

会稍大于校正前的剪切频率c ，因此，校正时应给校正装

置的最大超前相角

m

增加一个补偿角度

。

取15度

1sin 1sin m m

a ??

+=- m ω=

1aT

程序如下：

>> k=1000;

>> d1=conv(conv([1 0],[0.1 1]),[0.001 1]); >> scop=tf(k,d1);

>> w=logspace(0,4,50); >> bode(scop,w);

>> [Gm,Pm,Wcg,Wcp]=margin(scop)

结果为：

Gm =1.0100 Pm = 0.0584 Wcg =100.0000 Wcp = 99.4863

剪切频率

ω=99.5rad/sec，相角裕度0γ=0.0584deg和幅值0

C

裕度Gm=0.0864dB

4.2利用MATLAB语言计算出超前校正器的传递函数。

要计算出校正后系统的传递函数，就编写求超前校正器的传递函数的MATLAB程序，其中调用了求超前校正器传递函数的函数leadc()，leadc.m保存在matlab7.0\work\文件夹下，其中key=1时，为var=gama，是根据要求校正后的相角稳定裕度计算超前校正器；当key=2时，为var=wc，则是根据要求校正后的剪切频率计算校正器。若已知系统的开环

传递函数与要求校正后的相角稳定裕度或剪切频率，求系统串联超前校正器传递函数时，就可以调用此函数。leadc.m 编制如下：

function [Gc]=leadc(key,sope,vars)

% MATLAB FUNCTION PROGRAM leadc.m

%

if key==1

gama=vars(1);gama1=gama+5;

[mag,phase,w]=bode(sope);

[mu,pu]=bode(sope,w);

gam=gama1*pi/180;

alpha=(1-sin(gam))/(1+sin(gam));

adb=20*log10(mu);

am=10*log10(alpha);

wc=spline(adb,w,am);

T=1/(wc*sqrt(alpha));

alphat=alpha*T;

Gc=tf([T 1],[alphat 1]);

elseif key==2

wc=vars(1);

num=sope.num{1};den=sope.den{1};

na=polyval(num,j*wc);

da=polyval(den,j*wc);

g=na/da;

g1=abs(g);

h=20*log10(g1);

a=10^(h/10);

wm=wc;

T=1/(wm*(a)^(1/2));

alphat=a*T;

Gc=tf([T 1],[alphat 1]);

elseif key==3

gama=vars(1);wc=vars(2);gama1=gama+15;

num=sope.num{1};den=sope.den{1}; ngv=polyval(num,j*wc); dgv=polyval(den,j*wc); g=ngv/dgv; thetag=angle(g); thetag_d=thetag*180/pi; mg=abs(g);

gama_rad=gama1*pi/180;

z=(1+mg*cos(gama_rad-thetag))/(-wc*mg*sin(gama_rad-thetag)); p=(cos(gama_rad-thetag)+mg)/(wc*sin(gama_rad-thetag)); nc=[z,1];dc=[p,1]; Gc=tf(nc,dc); end

matlab 程序如下：

>> k=1000;

>> d1=conv(conv([1 0],[0.1 1]),[0.001 1]); >> sope=tf(k,d1);

>> gama=60;[Gc]=leadc(1,sope,[gama])

Transfer function: 0.02148 s + 1 -------------- 0.001056 s + 1

所以传递函数如右边所示 4.3

校验系统校正后系统是否满足题目要求：

程序如下： k0=1000; n1=1;

d1=conv(conv([1 0],[0.1 1]),[0.001 1]);s1=tf(k0,d1); s1=tf(k0*n1,d1); n2=[0.02148 1]; d2=[0.001056 1]; s2=tf(n2,d2);

sys=s1*s2;figure(1); margin(sys)

1

001056.01

02148.0G(S)++=

s s

剪切频率

c

ω

=167rad/sec ，相角裕度γ=45deg 和幅值裕度Gm =17.3dB

五、函数特征根的计算 5.1 校正前 开环传递函数：

K G(S)S(0.1S 1)(0.001S 1)

=

++

程序如下：

>> clear

>> k=1000;num=1;

>> den=conv(conv([1 0],[0.1 1]),[0.001 1]); >> s1=tf(k*num,den)

结果为：

Transfer function: 1000 --------------------------

0.0001 s^3 + 0.101 s^2 + s

故该系统闭环特征方程为：320.00010.10110000s s s +++= Matlab 程序如下：

>> clear

>> p=[0.0001 0.101 1 1000]; >> roots(p)

ans =

1.0e+003 *

-1.0099 -0.0000 + 0.0995i -0.0000 - 0.0995i

由于校正前系统单位负反馈的特征方程没有右半平面的根，故校正前的闭环系统稳定

5.2 校正后开环传递函数：

)

1001056.0)(1001.0)(11.0(G(S)GC(S))

102148.0(1000+++=

+S S S S S

>> num=[21.48 1000];

>> den=[0.0000001056 0.0003145 0.1031 1 0]; >> s=tf(num,den); >> s1=feedback(s,1)

21.48 s + 1000

------------------------------------------------------------

1.056e-007 s^4 + 0.0003145 s^3 + 0.1031 s^2 + 2

2.48 s + 1000

矫正后闭环传递函数

()1000

+ s 22.48 + s^2 0.1031 + s^3 0.0003145 + s^4 007-1.056e 1000

+ s 21.48=

Φs

程序如下:

>> clear

>> p=[0.0000001056 0.0003145 0.1031 22.48 1000]; >> roots(p)

ans =

1.0e+003 *

-2.6382 -0.1416 + 0.2078i -0.1416 - 0.2078i -0.0567

由于校正后系统单位负反馈的特征方程没有右半平面的根，故校正后的闭环系统稳定。六、系统动态性能特性

校正前开环传递函数：

K G(S)S(0.1S 1)(0.001S 1)

=

++

6.1 校正前单位阶跃响应

6.2 校正前单位脉冲响应>>K=1000;

>>den=[0.0001 0.1010 1 0]; >>G1=tf(k,den);

>>G0=feedback(G1,1);

>>t=0:1:100;

>>impulse(G0,t);

>>grid;

由阶跃响应求动态性能参数

要计算出阶跃响应动态性能参数，就编写求解阶跃响应动态性能参数的MATLAB程序，其中调用了函数perf()，perf.m保存在matlab7.0\work\文件夹下，其中key=1时，表示选择5%误差带，当key=2时表示选择2%误差带。y，t是对应系统阶跃响应的函数值与其对应的时间。函数返回的是阶跃响应超调量sigma(即σ)、峰值时间t p、调节时间t s。perf.m编制如下：

function [sigma,tp,ts]=perf(key,y,t)

%MATLAB FUNCTION PROGRAM perf.m

%

%Count sgma and tp

[mp,tf]=max(y);

cs=length(t);

yss=y(cs);

sigma= (mp-yss)/yss

tp=t(tf)

%Count ts

i=cs+1;

n=0;

while n==0,

i=i-1;

if key==1,

if i==1,

n=1;

elseif y(i)>1.05*yss,

n=1;

end;

elseif key==2,

if i==1,

n=1;

elseif y(i)>1.02*yss,

n=1;

end;

end

end;

t1=t(i);

cs=length(t);

j=cs+1;

n=0;

while n==0,

j=j-1;

if key==1,

if j==1,

n=1;

elseif y(j)<0.95*yss,

n=1;

end;

elseif key==2,

if j==1,

n=1;

elseif y(j)<0.98*yss,

n=1;

end;

end;

end;

t2=t(j);

if t2

if t1>t2;

ts=t1

end

elseif t2>tp,

if t2

ts=t2

else

ts=t1

end

end

程序如下：

>> clear

>> global y t;

>> s1=tf(1000,[0.0001 0.101 1 1000]); >> sys=feedback(s1,1); >> figure(1); >> step(sys);

>> [y,t]=step(sys); >> perf(1,y,t)

结果为：sigma =-2.0943 tp =0.4686 ts =0.4978 ans = -2.0943

6.3 校正前单位斜坡信号

在Simulink 窗口里菜单方式下的单位斜坡响应的动态结构图如下：

校正前单位斜坡响应曲线如图所示：

七、校正后动态性能分析

校正后开环传递函数：

)

1001056.0)(1001.0)(11.0(G(S)GC(S))

102148.0(1000+++=

+S S S S S

7.1 校正后单位阶跃响应

程序如下：

>> n1=[21.48 1000];

>> d1=conv(conv(conv([1 0],[0.1 1]),[0.001 1]),[0.001056 1]); >> s1=tf(n1,d1);

sys=feedback(s1,1);

step(sys)

结果为：

7.2 校正后单位冲击响应

程序如下：

>> n1=[21.48 1000];

>> d1=conv(conv(conv([1 0],[0.1 1]),[0.001 1]),[0.001056 1]); >> s1=tf(n1,d1);

>> sys=feedback(s1,1);

>> impulse(sys)

由阶跃响应求动态性能参数：

>> clear

global y t

>> num=[21.48 1000];

>> den=[0.0000001056 0.0003145 0.1031 1 0];

>> s=tf(num,den);

>> sys=feedback(s,1);

>> figure(1);

>> step(sys)

>> [y,t]=step(sys);

>> perf(1,y,t)

结果为：

sigma = 0.7102 tp = 0.04148 ts = 0.0520 ans =0.6013

校正前单位斜坡响应曲线如下所示：

单位脉冲、阶跃、斜坡响应曲线的相互对应关系是： 单位脉冲响应的积分是单位阶跃响应曲线。 单位阶跃响应的积分是单位斜坡响应。 八、系统的根轨迹分析

8.1、校正前根轨迹分析 校正前的开环传递函数为：

K G(S)S(0.1S 1)(0.001S 1)

=

++

程序如下：

>> n1=[1000];

>> d1=conv(conv([1 0],[0.1 1]),[0.001 1]); >> sys=tf(n1,d1); >> rlocus(sys)

确定分离点坐标：

分离点坐标d=-4.99增益k*=0.00249 与虚轴的交点：

与虚轴交点0+89.3j, 0-89.3j 增益k*=0.804

8.2、校正后根轨迹分析 校正后开环传递函数：)

1001056.0)(1001.0)(11.0(G(S)GC(S))

102148.0(1000+++=

+S S S S S

程序如下：

>> n1=[0.02148 1000];

>>d1=conv(conv(conv([1 0],[0.1 1]),[0.001 1]),[0.001056 1]); >>sys=tf(n1,d1);rlocus(sys)

确定分离点坐标：