请教一个分岔程序

relou

1. function s=rk45(t,q,flag,w)
   
2. 
   
3. %global m c k ks e g mb cb kc d f
   
4. 
   
5. m=0.5;
   
6. c=120;
   
7. k=4e4;
   
8. ks=40e4;
   
9. e=0.5e-3;
   
10. g=9.8;
    
11. mb=1;
    
12. cb=108;
    
13. kc=80e4;
    
14. d=1e-3;
    
15. f=0.05;
    
16. 
    
17. r=sqrt((q(1)-q(5))^2+(q(3)-q(7))^2);
    
18. %r=sqrt(q(1)^2+q(3)^2);
    
19. if r<d
    
20.     fx=0;
    
21.     fy=0;
    
22. else
    
23.     %fx=-(q(1)-q(5))*kc*(1-d/r)+(q(3)-q(7))*f*kc*(1-d/r);
    
24.     fx=-kc*(1-d/r)*(q(1)-f*q(3));
    
25.     %fy=-(q(3)-q(7))*kc*(1-d/r)-(q(1)-q(5))*f*kc*(1-d/r);
    
26.     fy=-kc*(1-d/r)*(f*q(1)+q(3));
    
27. end
    
28. 
    
29. s=zeros(8,1)
    
30. s(1)=q(2);
    
31. 
    
32. s(2)=fx/(m*w^2)-c*q(2)/(m*w)-k*q(1)/(m*w^2)-ks*(q(1)^2+q(3)^2)*q(1)/(m*w^2)+e*cos(t);
    
33. 
    
34. s(3)=q(4);
    
35. 
    
36. s(4)=fy/(m*w^2)-c*q(4)/(m*w)-k*q(3)/(m*w^2)-ks*(q(1)^2+q(3)^2)*q(3)/(m*w^2)+e*sin(t)-g/(w^2);
    
37. 
    
38. s(5)=q(6);
    
39. 
    
40. s(6)=-fx/(mb*w^2)-cb*q(6)/(mb*w)-kc*q(5)/(mb*w^2);
    
41. 
    
42. s(7)=q(8);
    
43. 
    
44. s(8)=-fy/(mb*w^2)-cb*q(8)/(mb*w)-kc*q(7)/(mb*w^2)-g/(w^2);
    
45. 
    
46. 
    
47. end
    

复制代码

1. function nonlinear
   
2. q0=[0;0;0;0;0;0;0;0]
   
3. w=linspace(2200,3001,100)
   
4. options=odeset('RelTol',1e-6,'AbsTol',1e-6)
   
5. y=[]
   
6. for n=1:length(w)
   
7.     T=2*pi
   
8.     ts=[0:T/200:500*T]
   
9.     [tt,z]=ode45('rk45',ts,q0,options,w(n))
   
10.     q0=z(end,:)
    
11.     %figure(1)
    
12.     subplot(2,1,1)
    
13.     plot(w(n),z(90000:200:end,1),'r.','markersize',1)
    
14.     hold on
    
15.     subplot(2,1,2)
    
16.     plot(w(n),z(90000:200:end,3),'k.','markersize',1)
    
17.     hold on
    
18. end
    
19.     %figure(2)
    
20.     %subplot(2,1,1)
    
21.     %plot(tt(50000:end),z(50000:end,1),'k-')
    
22. %%hold on
    
23. y=z(end,:)
    
24. save 'end_data.txt'y-ascii
    
25. save 'time.txt'tt-ascii
    
26. 
    
27. 我参考刘长春、吴敬东的《非线性转子-机匣系统的动力学特性》的文章写的一个程序，但是结果就和原文章的对不上，请各位前辈指点迷津，不胜感激

复制代码

meiyongyuandeze

给几章对比下吧

relou

补发文章中的一些参数和图片。

dj15175771160

你这分岔图采用的是什么方法？我也在做，跟参考的文献也是出入很大，能大家讨论下吗

huizi_nice

权限有限，看不到图，估计是初值有变吧，很多文献上给的程序和文中的图是对不上的