武汉理工大学-计算机仿真实验作业答案

五、（10分）已知系统的传递函数为6

168682)(232+++++=s s s s s s G 。

语言建立系统传递函数模型，并求：

⑴ 该系统的单位阶跃响应；（2分）

⑵ 输入函数为u(t)时的响应；（3分）

（u(t)正弦信号，周期2秒，仿真时间8秒，采样周期0.1）；

(3) 输入函数为u(t)时的响应；（3分）

（u(t)方波输入信号，周期10秒，仿真时间20秒，采样周期0.05）

(4) 绘出系统的波德图（Bode ）。（2分）

解答：

num=[2 8 6];

den=[1 8 16 6];

sys=tf(num,den);

t=0:0.1:8;

y1=step(sys,t);

u=sin(t*pi);

y2=lsim(sys,u,t);

subplot(2,2,1);plot(t,y1);

grid;

title('阶跃响应曲线');

xlabel('响应时间');

ylabel('响应值');

hold on;

subplot(2,2,2);plot(t,y2);

grid on ;

title('对sin(t)的响应曲线');

xlabel('响应时间');

ylabel('响应值');

t=0:0.05:20

u=square(pi/5*t)

y3=lsim(sys,u,t);

subplot(2,2,3);plot(t,y3)

grid on ;

title('对方波信号的响应曲线');

xlabel('响应时间');

ylabel('响应值');

subplot(2,2,4);bode(sys);

grid ;title('bode 图');

运行结果：

六、（10分）设二阶动力学系统的传递函数如下，假设将无阻尼固有频

率固定为ωn ＝1 rad/s ，将阻尼比的值分别设置成ζ＝0，0.1，0.2，0.3，…，

MATLAB 语言编程，分析在这些阻尼比ζ的取值下该系统的阶跃响应。

2222)(n

n n s s s G ωςωω++= 解答：wn=1;

kesi=[0:0.1:1,2,3,4,5];

figure('color',[1 1 1]);

hold on

for i=kesi

num=wn.^2

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

step(num,den);

end

title('The Step Response of Two Order system')

运行结果:

七、（20分）对图示的车辆（汽车、摩托车等）悬架系统进行数字仿真。

m 1＝40 Kg 、

2、轮胎刚度k 1＝158 kN/m 、悬架刚度k 2＝

23 kN/m 、b 是悬架减振器阻尼系数，x 是簧下质量的振动位移、

y 是簧上质量的振动位移，u 是路面不平度函数。研究以u 为输入，

x 、y 为输出时，系统的动态响应。完成下述任务：

⑴建立系统数学模型，取车身垂直振动速度和位移、车轮垂直振动速度和

位移为四个状态变量，把系统数学模型转换为状态方程形式；（5分）

⑵绘制系统的函数方块图或状态变量图；（5分）

⑶取阻尼系数b 的值分别为600、1000、2000、4000，对系统

作阶跃输入动态仿真并绘制出系统的Bode 图。（7分）

⑷对仿真结果进行分析，给出分析结论。（3分）

解答：

取车身垂直振动速度和位移、车轮垂直振动速度和位移为四个状态变量 建立系统状态方程模型

函数方块图如下：

12122()()()()()m x k y x b y x k u x m y k y x b y x =-+-+-=---- 123142x x x y x x x x x y ====== [][]12122341121223421()1x k k x k x bx bx k u m y k x k x bx bx m =-++-++=-+- 1122122111113314422222200100000100x x x x k k k b b u k m m m m x x m x x k k b b m m m m ⎡⎤⎡⎤⎢⎥⎡⎤⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥+⎢⎥⎢⎥--⎢⎥=+⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦⎣⎦⎢⎥--⎢⎥⎣⎦⎣⎦

示波器显示结果：

源程序：

m1=40;

m2=400;

k1=158;

k2=23;

k1=1000*k1;k2=1000*k2;

p(1)=600;

p(2)=1000;

p(3)=2000;

p(4)=4000;

for n=1:4

b=p(n);

sys(n)=ss(A,B,C,D);

A=[0 0 1 0;0 0 0

1;-(k1+k2)/m1,k2/m1,-b/m1,b/m1;k2/m2,-k2/m2,b/m2,-b/m2]; B=[0; 0; k1/m1; 0]

disp('display:matrix A ,B');

disp('matrix A='),disp(A),

disp('matrix B='),disp(B),

ws=sqrt(k2/m2)/m1,disp(ws),

wt=sqrt(k2+k1/m1)/m1,disp(ws),

C=zeros(1,4);C(2)=1;D=0;

end

figure('name','悬架系统的阶跃响应','numbertitle','off'); t=0:0.01:25;

step(sys(1),sys(2),sys(3),sys(4),t);

grid;title('簧上质量位移响应曲线');

xlabel('响应时间');ylabel('响应值');

figure('name','悬架系统的bode图','numbertitle','off'); w=logspace(0.3,3.0);

bode(sys(1),sys(2),sys(3),sys(4),w);

grid;title('簧上质量的bode图');

运行结果：