机械系统动力学MATLAB仿真(上)
寒假期间预习李有堂编著、国防工业出版社出版的《机械系统动力学》,对其中第四章:单自由度系统的振动,第五章:两自由度系统的振动,的部分内容进行matlab编程仿真验证,并做简单研究。
目录
1 单自由度系统的振动
1.1 有阻尼系统自由振动
1.2 谐波激励下的有阻尼强迫振动1.3 强迫振动的应用1.3.1 隔振原理1.3.2 测振原理1.3.3 惯性振动筛
2 两自由度系统的振动
2.1 谐波激励下的强迫振动2.2 阻尼对强迫振动的影响2.3 阻尼减振器
1.1 有阻尼系统自由振动
有阻尼系统自由振动在欠阻尼、临界阻尼和过阻尼三种情况下的位移曲线和相图如下图所示。
实现代码如下。
1%frictional oscillation simulation
2%lack friciton, critical friction and over friction
3%2019.02.19
4%Author: Jing-Xuan Yang
5%page: 73
6
7zeta = [3.4641, 34.641, 173.2051];%three conditions of friction
8tspan = linspace(0,2,400);%time
9%first condition
10[t,x] = ode45('FreeOcillation',tspan,[0 1],[],zeta(1));%solve 2-ode
11subplot(2,1,1);
12plot(t,x(:,1),'-b');
13hold on
14subplot(2,1,2)
15plot(x(:,1),x(:,2),'-b');
16hold on
17%second condition
18[t,x] = ode45('FreeOcillation',tspan,[0 1],[],zeta(2));
19subplot(2,1,1);
20plot(t,x(:,1),'--g');
21hold on
22subplot(2,1,2)
23plot(x(:,1),x(:,2),'--g');
24hold on
25%third condition
26[t,x] = ode45('FreeOcillation',tspan,[0 1],[],zeta(3));
27subplot(2,1,1);
28plot(t,x(:,1),':r');
29hold on
30subplot(2,1,2)
31plot(x(:,1),x(:,2),':r');
32hold on
33%draw titles and legends
34subplot(2,1,1)
35xlabel('Time(\tau)');
36ylabel('Displacement x(\tau)');
37title('Displacement as a function of \tau');
38legend('\alpha = 0.1','\alpha = 1','\alpha = 5');
39subplot(2,1,2);
40xlabel('Displacement x(\tau)');
41ylabel('Velocity v(\tau)');
42title('Phase portrait');
43legend('\alpha = 0.1','\alpha = 1','\alpha = 5');
44
45function xdot = FreeOcillation(t,x,dummy,zeta)
46
47xdot = [x(2);-2.0*zeta*x(2)-1200*x(1)];
48
49end
1.2 谐波激励下的有阻尼强迫振动
谐波激励下有阻尼强迫振动的幅频特性曲线和相频特性曲线如下图所示。
实现代码如下。
1%Phase Diagram for Different Damping Ratio
2%Amplitude-Frequency Characteristic Curve
3%2019/01/17
4%Author: Jing-Xuan Yang
5%page: 84-85
6
7x = 0:0.01:3;
8phi = ones(size(x));
9H = ones(7,301);
10
11for i = 1:301
12 ksi = 1;
13 phi(i) = atan2(2*ksi*x(i),(1-(x(i))^2));
14 H(1,i) = 1/sqrt((1-(x(i))^2)^2+(2*ksi*x(i))^2);
15end
16%Phase Diagram for Different Damping Ratio
17figure(1)
18plot(x,phi,'b.-');
19hold on;
20
21for i = 1:301
22 ksi = 0.5;
23 phi(i) = atan2(2*ksi*x(i),(1-(x(i))^2));
24 H(2,i) = 1/sqrt((1-(x(i))^2)^2+(2*ksi*x(i))^2);
25end
26plot(x,phi,'ko-');