声振论坛

 找回密码
 我要加入

QQ登录

只需一步,快速开始

楼主: TNC

[共享资源] 单自由度振动动画

  [复制链接]
发表于 2013-5-29 19:13 | 显示全部楼层
太牛了
回复 支持 反对
分享到:

使用道具 举报

发表于 2013-12-11 22:52 | 显示全部楼层
虽然该帖子有点“久远”了,我还是顶一下,顺便上传一个本人在学习《汽车振动学》时,做的一个Matlab程序,和楼主的差不多,但是相对来说比较简单,对于初学matlab(编程)的人来说,比较容易理解,对学习matlab编程有用!以下是程序:
clear all;
clear;
clc;
rectangle('position',[12,1.2,2,0.2],'FaceColor',[0.5,0.3,0.4]);
axis([0,15,-5,15]);
%画顶板
hold on
%画直线
y=1.5:.35:5.5;
M=length(y);
x=12+mod(1:M,2)*2;
x(1)=13;
x(end-3:end)=13;
D=plot(x,y);
%弹簧 D
C=(1/8:1/4:7/8)*2*pi;r=0.35;
t1=2*r*sin(C);
F1=fill(13+3*r*cos(C),5.5+t1,'r');
% 球
% set(gca,'ytick',[0:2:9]);
% set(gca,'yticklabels',num2str([-1:3]'));
plot([0,15],[5.9,5.9],'black');
H1=plot([0,13],[5.9,5.9],'y');
% 句柄[黄线]
Q=plot(0,5.5,'color','r');
% 运动曲线;
td=[];yd=[];
text(2,14,'阻尼振动','fontsize',16);
set(gcf,'doublebuffer','on');
%有阻尼自由衰减振动分析计算程序
k=80;   %N/m
m=10;   %Kg
c=0.219*sqrt(k*m)*2;%N.s/m
x0=-2;
v0=-10;
X1=[x0,v0]';
p=sqrt(k/m);
n=c/(2*m);
f=p/(2*pi);
ksai=n/p;
T0=1/f;
T1=T0*(1/sqrt(1-ksai^2));
T1_T=T1/T0;
A2_A1=1/exp(n*T1);
s1=-n+sqrt(n^2-p^2);
s2=-n-sqrt(n^2-p^2);
s=[1,1;s1,s2];
Cd=inv(s)*X1;
T=0;
while T<12;
pause(0.2);
Dy=real(Cd(1)*exp(s1*T)+Cd(2)*exp(s2*T))/4+1.15;
Y=(y-1.5)*Dy+1.5;
Yf=Y(end)+t1;
td=[td,T];yd=[yd,Y(end)];
set(D,'ydata',Y);
set(F1,'ydata',Yf,'facecolor',rand(1,3));
set(H1,'xdata',[T,13],'ydata',[Y(end),Y(end)]);
set(Q,'xdata',td,'ydata',yd) ;
T=T+0.1;
end
搜狗截图20131211223941.png

点评

赞成: 5.0
赞成: 5
楼主真乃学霸也 谢谢无私分享经验啊 太感谢了!!!必须顶啊  发表于 2013-12-12 11:45
赞成: 5
楼上辛苦了,很不错,值得学习!谢谢了!!!  发表于 2013-12-12 10:52
赞成: 5
楼上辛苦了,很不错,值得学习!谢谢了!!!  发表于 2013-12-12 10:47
发表于 2013-12-12 14:53 | 显示全部楼层
发表于 2013-12-12 18:08 | 显示全部楼层
为方便大家复制代码,我重新整理一下程序!
  1. clear all;
  2. clear;
  3. clc;
  4. rectangle('position',[12,1.2,2,0.2],'FaceColor',[0.5,0.3,0.4]);
  5. axis([0,15,-5,15]);
  6. %画顶板
  7. hold on
  8. %画直线
  9. y=1.5:.35:5.5;
  10. M=length(y);
  11. x=12+mod(1:M,2)*2;
  12. x(1)=13;
  13. x(end-3:end)=13;
  14. D=plot(x,y);
  15. %弹簧 D
  16. C=(1/8:1/4:7/8)*2*pi;r=0.35;
  17. t1=2*r*sin(C);
  18. F1=fill(13+3*r*cos(C),5.5+t1,'r');
  19. % 球
  20. % set(gca,'ytick',[0:2:9]);
  21. % set(gca,'yticklabels',num2str([-1:3]'));
  22. plot([0,15],[5.9,5.9],'black');
  23. H1=plot([0,13],[5.9,5.9],'y');
  24. % 句柄[黄线]
  25. Q=plot(0,5.5,'color','r');
  26. % 运动曲线;
  27. td=[];yd=[];
  28. text(2,14,'阻尼振动','fontsize',16);
  29. set(gcf,'doublebuffer','on');
  30. %有阻尼自由衰减振动分析计算程序
  31. k=80;   %N/m
  32. m=10;   %Kg
  33. c=0.219*sqrt(k*m)*2;%N.s/m
  34. x0=-2;
  35. v0=-10;
  36. X1=[x0,v0]';
  37. p=sqrt(k/m);
  38. n=c/(2*m);
  39. f=p/(2*pi);
  40. ksai=n/p;
  41. T0=1/f;
  42. T1=T0*(1/sqrt(1-ksai^2));
  43. T1_T=T1/T0;
  44. A2_A1=1/exp(n*T1);
  45. s1=-n+sqrt(n^2-p^2);
  46. s2=-n-sqrt(n^2-p^2);
  47. s=[1,1;s1,s2];
  48. Cd=inv(s)*X1;
  49. T=0;
  50. while T<12;
  51. pause(0.2);
  52. Dy=real(Cd(1)*exp(s1*T)+Cd(2)*exp(s2*T))/4+1.15;
  53. Y=(y-1.5)*Dy+1.5;
  54. Yf=Y(end)+t1;
  55. td=[td,T];yd=[yd,Y(end)];
  56. set(D,'ydata',Y);
  57. set(F1,'ydata',Yf,'facecolor',rand(1,3));
  58. set(H1,'xdata',[T,13],'ydata',[Y(end),Y(end)]);
  59. set(Q,'xdata',td,'ydata',yd) ;
  60. T=T+0.1;
  61. end
复制代码
补充一个多自由度的振动程序(参考张亚辉结构动力学基础)
  1. function newmark(t_max,dt)
  2. %t_max为持续时间,dt为时间步长
  3. M=1*[1 0 0;0 1 0;0 0 1];
  4. K=1*[3 -1 0;-1 2 -1;0 -1 1];
  5. C=0.015*M+0.02*K;%c=0.015*[0 0 0;0 0 0;0 0 0];
  6. u=[0 0 0]';
  7. v=[0 0 0]';
  8. a=[0 0 0]';
  9. t(1)=0;               %时间
  10. x(:,1)=u;             %位移
  11. x1(:,1)=v;            %速度
  12. x2(:,1)=a;            %加速度
  13. %newmark参数
  14. gama=0.5;
  15. delta=0.25;
  16. a0=1/(delta*dt^2);
  17. a1=gama/(delta*dt);
  18. a2=1/(delta*dt);
  19. a3=1/(2*delta)-1;
  20. a4=gama/delta-1;
  21. a5=dt*(gama/(2*delta)-1);
  22. a6=dt*(1-gama);
  23. a7=gama*dt;
  24. %等效刚度矩阵
  25. K1=K+a0*M+a1*C;
  26. i=1;
  27. t(1)=0;
  28. while t(i)<t_max
  29.       Q=-1*[sin(t(i)) sin(t(i)) sin(t(i))]';
  30.       q(:,i+1)=Q+M*(a0*x(:,i)+a2*x1(:,i)+a3*x2(:,i))+C*(a1*x(:,i)+a4*x1(:,i)+a5*x2(:,i));
  31.       x(:,i+1)=inv(K1)*q(:,i+1);
  32.       x2(:,i+1)=a0*(x(:,i+1)-x(:,i))-a2*x1(:,i)-a3*x2(:,i);
  33.       x1(:,i+1)=a1*(x(:,i+1)-x(:,i))-a4*x1(:,i)-a5*x2(:,i);
  34.       i=i+1;
  35.       t(i)=t(i-1)+dt;
  36. end  
  37. x1=x(1,:)';
  38. x2=x(2,:)';
  39. x3=x(3,:)';
  40. xx1=x1(1:t_max/dt)';
  41. xx2=x2(1:t_max/dt)';
  42. xx3=x3(1:t_max/dt)';
  43. t=dt:dt:t_max;
  44. plot(t,xx1,'-',t,xx2,'-.',t,xx3,'--')
  45. xlabel('t[sec]')
  46. ylabel('x(m)')
  47. legend('x1(t)','x2(t)','x3(t)')
  48. grid
复制代码
在Matlab命令窗口输入:输入:newmark(30,0.01)
另外,大家还可以参考一下这个帖子:http://forum.vibunion.com/thread-30794-1-1.html(关于多自由度 自由振动 MATLAB 仿真 的,里面有程序,又仿真结果,很详细!!)

点评

赞成: 5.0
不足之处,还请“专家”们积极补充、指正!!  发表于 2014-3-19 12:49
赞成: 5
  发表于 2014-3-5 17:54
赞成: 5
  发表于 2013-12-15 20:36
回复 支持 1 反对 0

使用道具 举报

发表于 2013-12-13 12:44 | 显示全部楼层
沧海一声笑 发表于 2013-12-12 14:53
谢谢分享!!

客气了 ,应该的!!!楼主,是Matlab编程初学者吗?我这还有些资料,有时间给你发过去!!
发表于 2013-12-15 20:37 | 显示全部楼层
牛小贱 发表于 2013-12-13 12:44
客气了 ,应该的!!!楼主,是Matlab编程初学者吗?我这还有些资料,有时间给你发过去!!

发表于 2014-2-25 15:50 | 显示全部楼层
牛小贱 发表于 2013-12-13 12:44
客气了 ,应该的!!!楼主,是Matlab编程初学者吗?我这还有些资料,有时间给你发过去!!

楼主,话说能不能给我发点资料  小弟最近也在编一个类似的程序,但是没有思路,不知从何下手!
发表于 2014-2-25 16:00 | 显示全部楼层
奋斗小夏 发表于 2014-2-25 15:50
楼主,话说能不能给我发点资料  小弟最近也在编一个类似的程序,但是没有思路,不知从何下手!

我在本贴的49楼已经留下这个比较简单的单自由度振动动画的程序了:http://forum.vibunion.com/forum. ... 6579&fromuid=10
另外,我已经把【MATLAB编程(最好中文教程)】资料上传到新浪爱问资料共享上去了,你可以看一下,链接为http://ishare.iask.sina.com.cn/f/67354723.html(PS:你可以看一下Matlab在振动分析应用这方面的资料,对你会有帮助!但是,既然你也要编一个类似的程序,那我建议你还是把这个小程序吃透,关键是要理解里面每个语句的含义……)希望对你有所帮助!

点评

赞成: 5.0
赞成: 5
  发表于 2014-3-5 17:54
回复 支持 1 反对 0

使用道具 举报

发表于 2014-2-25 21:30 | 显示全部楼层
牛小贱 发表于 2014-2-25 16:00
我在本贴的49楼已经留下这个比较简单的单自由度振动动画的程序了:http://forum.vibunion.com/forum.php? ...

好的 你辛苦了   非常地感谢你!
发表于 2016-9-3 14:38 | 显示全部楼层
厉害
发表于 2018-4-2 11:28 | 显示全部楼层
忒厉害了也
发表于 2018-4-18 13:09 | 显示全部楼层
您需要登录后才可以回帖 登录 | 我要加入

本版积分规则

QQ|小黑屋|Archiver|手机版|联系我们|声振论坛

GMT+8, 2024-3-29 13:46 , Processed in 0.080295 second(s), 26 queries , Gzip On.

Powered by Discuz! X3.4

Copyright © 2001-2021, Tencent Cloud.

快速回复 返回顶部 返回列表