马上注册,结交更多好友,享用更多功能,让你轻松玩转社区。
您需要 登录 才可以下载或查看,没有账号?我要加入
x
本帖最后由 tudexin 于 2010-11-20 02:34 编辑
论坛里的各位老师,毕业论文需要非线性这一块知识,刚接触非线性这一块时间不长,暑假才开始接触,以前不知道有这么好的论坛的存在。以前只是在MTALAB中文论坛里逛,那里面对于非线性这一块的东西还是比较少。对论坛看了2天后,对分岔图和庞加莱图的做法才有个初步的认识,对照论坛里的频闪法做了汽车半车模型的分岔图,可是发现问题很多。 我做出的分岔图不分岔,感觉像是幅频图,很是郁闷。不知道是程序编写错误,还是系统模型建立错误。恳请论坛里的各位专家指点迷津。 先贴上代码:function dy=biyelunwen21(t,y,flag,f);
mb=1180;
mf=50;
mr=45;
cf2=1385.4;
cr2=1385.4;
cf1=3012;
cr1=3012;
kf1=140000;
kr1=140000;
k1=12394;
k2=-73696;
k3=3170400;
A=0.05;
J=633.615;
a=1.123;
b=1.377;
xrd=A*sin(2*pi*f*t+pi/9);
dxrd=2*pi*f*A*cos(2*pi*f*t+pi/9);
xfd=A*sin(2*pi*f*t);
dxfd=2*pi*f*A*cos(2*pi*f*t);
dy(1,1)=y(2);
dy(2,1)=((k1+2*k2*(y(3)-y(1)+a*y(7))+3*k3*(y(3)-y(1)+a*y(7))^2)*(y(3)-y(1)+a*y(7))+cf2*(y(4)-y(2)+a*y(8))+(k1+2*k2*(y(5)-y(1)-b*y(7))+3*k3*(y(5)-y(1)-b*y(7))^2)*(y(5)-y(1)-b*y(7))+cr2*(y(6)-y(2)-b*y(8)))/mb;
dy(3,1)=y(4);
dy(4,1)=((k1+2*k2*(xfd-y(3))+3*k3*(xfd-y(3))^2)*(xfd-y(3))+cf1*(dxfd-y(4))-(k1+2*k2*(y(3)-y(1)+a*y(7))+3*k3*(y(3)-y(1)+a*y(7))^2)*(y(3)-y(1)+a*y(7))-cf2*(y(4)-y(2)+a*y(8)))/mf;
dy(5,1)=y(6);
dy(6,1)=((k1+2*k2*(xrd-y(5))+3*k3*(xrd-y(5))^2)*(xrd-y(5))+cf1*(dxrd-y(6))-(k1+2*k2*(y(5)-y(1)-b*y(7))+3*k3*(y(5)-y(1)-b*y(7))^2)*(y(5)-y(1)-b*y(7))-cr2*(y(6)-y(2)-b*y(8)))/mr;
dy(7,1)=y(8);
dy(8,1)=(-a*((k1+2*k2*(y(3)-y(1)+a*y(7))+3*k3*(y(3)-y(1)+a*y(7))^2)*(y(3)-y(1)+a*y(7))+cf2*(y(4)-y(2)+a*y(8)))+b*((k1+2*k2*(y(5)-y(1)-b*y(7))+3*k3*(y(5)-y(1)-b*y(7))^2)*(y(5)-y(1)-b*y(7))+cr2*(y(6)-y(2)-b*y(8))))/J;
分岔图程序:clear;
options=odeset('RelTol',1e-6);
for f=0.1:0.1:10
tt=1/f;
y0=[0,0,0,0,0,0,0,0];
[t,y]=ode45('biyelunwen21',[0:tt/100:1000*tt],y0,options,f) ;
hold on
plot(f,y(90000:100:end,1),'r.')
end
最后出来的随激励频率f的分岔图是这样的,不是连续的,而且问题的关键是不分岔,而且出现了跳跃,我把跳跃的部分又细作了下,还是跳跃,而且越看越像幅频图,郁闷。肯求无水,octopussheng,咕噜噜等前辈指定迷津!
|