Matlab实验(五)
1.用三(五)种方法绘制单位阶跃响应曲线。
%%%利用实验四的系统%%%
n =[ 0 0.1606 0.8575 0.0661];
d =[1.0000 1.1016 0.8119 0.0661];
sys = tf(n,d);
sys = feedback(sys,1)
(1)step(无返回值)直接画出图像
figure(1);
step(sys)
(2)step(有返回值)+plot
figure(2);
[y,t]=step(sys);
plot(t,y,'-r')
grid on;
xlabel('Time(second)');
ylabel('Amplitude');
title('Step Reponse');
(3)使用lsim(注意u和t需要自己设定)
t=0:0.001:40;
u=ones(1,length(t));
figure(3);
lsim(sys,u,t);
(4)使用有返回值的lsim+plot
figure(4);
[y,t]=lsim(sys,u,t);
plot(t,y,'--b')
grid on;
xlabel('Time(second)');
ylabel('Amplitude');
title('Step Reponse');
(5)使用ltiview
ltiview({'step'},sys)
总结:类似的impluse(单位脉冲响应)也有响应的做法。有返回值与无返回值的step是有很大的区别的,在无返回值的step中,它可以自动帮你画出图像,在GUI中还可以显示响应的时域指标(tp、tr、ts等等)。有返回值要求的话得花费一番功夫,也有一个缺点,见下面第二题。
实验图像
Step版本(无返回值)
|
Step+plot
|
Lsim无返回值
Lsim+plot
Ltiview
2.对以上系统的单位负反馈系统进行稳态性能分析
function [pos,tr,ts,tp]=stepchar(g0,delta)
[y,t]=step(g0);
[mp,ind]=max(y);
dimt=length(t);
yss=y(dimt);
pos=100*(mp-yss)/yss;
tp=t(ind);
for i=1:dimt
if y(i)>=yss
tr=t(i);
break;
end
end
for i=1:length(t)
if y(i)<=(1-delta)*yss|y(i)>=(1+delta)*yss
ts=t(i);
end
end
%或者如下求取ts
% while i<length(t)
% if y(i)<=(1-delta)*yss|y(i)>=(1+delta)*yss
% ts=t(i);
% end
% i=i+1;
% end
end
%调用
n =[ 0 0.1606 0.8575 0.0661];
d =[1.0000 1.1016 0.8119 0.0661];
sys = tf(n,d);
sys = feedback(sys,1)
isstable(sys);
[pos,tr,ts,tp]=stepchar(sys,0.05)
step(sys)
%结果
sys =
0.1606 s^2 + 0.8575 s + 0.0661
----------------------------------
s^3 + 1.262 s^2 + 1.669 s + 0.1322
Continuous-time transfer function.
pos =
25.3768
tr =
1.5392
ts =
4.3682
tp =
2.6100
这个求取时域指标的代码来自某位xx,对比函数所求结果,发现不太对劲。结果差很多,这里给出个人的一个理解:
1、仿真时间与仿真精度会影响结果。采样时间越紧密,结果越准确,采样时间越长越准确。
2、函数版本可能最终结果不是稳态值(还有些震荡,不过很小)。
故,个人改进了一下,用dcgain函数求取终值。
修改部分:
%%%%修正部分
%yss=y(dimt);
yss=dcgain(g0);
最终结果:
pos =
27.3281
tr =
1.5637
ts =
4.5348
tp =
2.5801
基本吻合,所以改正没问题。
2.对以上系统的单位负反馈系统进行稳定性分析。
n =[ 0 0.1606 0.8575 0.0661];
d =[1.0000 1.1016 0.8119 0.0661];
sys = tf(n,d);
sys = feedback(sys,1)
%方法一
isstable(sys)
%方法二
eig(sys)
%方法三
pzmap(sys)
又或者如下自己编写一个函数
function y=ISstable(sys)
%输入系统,输出y
%1:稳定 0:不稳定
sys=tf(sys);
%sys=feedback(sys,1,-1);
[num,den]=tfdata(sys,'v');
r=roots(den);
a=find(real(r)>0);
if a>0
disp("系统不稳定");
else
disp("系统稳定");
end
%%%另外一种 原理一样
% temp=0;
% for i = 1:length(r)
% if real(r(i))<0
% temp=temp+1;
% else
% break;
% end
% end
% if temp==length(r)
% disp("系统稳定");
% y=1;
% else
% disp("系统不稳定");
% y=0;
% end
end