分岔圖的總結(jié)

1、混沌動(dòng)力學(xué)行為研究的程序說(shuō)明1.分岔圖：一個(gè)耦合發(fā)電機(jī)系統(tǒng)：耦合系統(tǒng)函數(shù)：function dx=ouhe1(t,x)dx(1,1)=-x(4)*x(1)+x(2)*(x(3)+x(5);dx(2,1)=-x(4)*x(2)+x(1)*(x(3)-x(5);dx(3,1)=x(3)-x(1)*x(2);dx(4,1)=0;dx(5,1)=0;分岔圖程序：clear;clc;Z=;for a=linspace(0.5,10.5,500); %舍棄前面迭帶的結(jié)果，用后面的結(jié)果畫圖：即0.5-10.5分為500點(diǎn) T,Y=ode45('ouhe1',1,1;1;1;2;a); T,Y

2、=ode45('ouhe1',20,Y(length(Y),:); Y(:,1)=Y(:,2)-Y(:,1); for k=2:length(Y) f=k-1; if Y(k,1)<0 if Y(f,1)>0 y=Y(k,2)-Y(k,1)*(Y(f,2)-Y(k,2)/(Y(f,1)-Y(k,1); Z=Z a+abs(y)*i; end else if Y(f,1)<0 y=Y(k,2)-Y(k,1)*(Y(f,2)-Y(k,2)/(Y(f,1)-Y(k,1); Z=Z a+abs(y)*i; end end endend plot(Z,'.

3、9;,'markersize',1);title('ouhe映射分岔圖');xlabel('a'),ylabel('|y|')2.功率譜：對(duì)上面耦合系統(tǒng)的功率譜的研究clear allx=zeros(1,10001);y=zeros(1,10001);z=zeros(1,10001);%數(shù)組置零x(1)=1;y(1)=1;z(1)=1;h=0.001;k=10000;a=3;u=2;for i=1:kx(i+1)=x(i)+h*(-u*x(i)+y(i)*(z(i)+a); %歐拉離散y(i+1)=y(i)+h*(-u*y(i)+

4、x(i)*(z(i)-a);z(i+1)=z(i)+h*(z(i)-x(i)*y(i);endX1=fft(z,16384); %對(duì)x做傅里葉變換，取8192個(gè)點(diǎn)p=X1.*conj(X1)/16384; %求x的模及功率譜密度，單位：dB 同樣可求y或zc=100*0:8191/16384; %取雙邊，也可取單邊c=0:4095/0.8192; %title('題目') %頂端題目plot(c,log10(p(1:8192),'k')%畫出左半部分 axis(0 4 -3 6); %限制橫、縱坐標(biāo)范圍%plot(c,abs(X1(1:4096) 有時(shí)候也可用x

5、的絕對(duì)值表示功率大小，沒求對(duì)數(shù)xlabel('itfrm/HZ','fontsize',18,'fontName','times new Roman','fontweight','bold','color','k'); %加橫坐標(biāo),it表傾斜，rm表復(fù)正ylabel('power spectrum/dB','fontsize',18,'fontName','times new Roman','fo

6、ntweight','bold','color','k'); %縱坐標(biāo)標(biāo)示3.最大李雅譜指數(shù)程序：仍對(duì)上述耦合系統(tǒng)clear; %與ouhe系統(tǒng)的分岔圖（fotran）符合的很好clc;d0=1e-8;le=0;lsum=0;x=1;y=1;z=1;x1=1;y1=1;z1=1+d0;for i=1:500 T1,Y1=ode45('ouhe1',0,1,x;y;z;4.05;7.0); T2,Y2=ode45('ouhe1',0,1,x1;y1;z1;4.05;7.0); n1=length(Y1); n

7、2=length(Y2); x=Y1(n1,1); y=Y1(n1,2); z=Y1(n1,3); x1=Y2(n2,1); y1=Y2(n2,2); z1=Y2(n2,3); d1=sqrt(x-x1)2+(y-y1)2+(z-z1)2); x1=x+(d0/d1)*(x1-x); y1=y+(d0/d1)*(y1-y); z1=z+(d0/d1)*(z1-z); if i>100 lsum=lsum+log(d1/d0); endend le=lsum/(i-100)最大李雅譜指數(shù)譜程序：仍對(duì)上述耦合系統(tǒng)clear; %可調(diào)參數(shù)區(qū)間和步長(zhǎng)，如a的第5和36行clc;LE1=;d0=1

8、e-8;for a=linspace(0.5,10.5,300); le=0; lsum=0; x=1; y=1; z=1; x1=1; y1=1; z1=1+d0; for i=1:150 T1,Y1=ode45('ouhe1',0,1,x;y;z;2;a); T2,Y2=ode45('ouhe1',0,1,x1;y1;z1;2;a); n1=length(Y1); n2=length(Y2); x=Y1(n1,1); y=Y1(n1,2); z=Y1(n1,3); x1=Y2(n2,1); y1=Y2(n2,2); z1=Y2(n2,3); d1=sqrt(

9、x-x1)2+(y-y1)2+(z-z1)2); x1=x+(d0/d1)*(x1-x); y1=y+(d0/d1)*(y1-y); z1=z+(d0/d1)*(z1-z); if i>50 lsum=lsum+log(d1/d0); end end le=lsum/(i-50); LE1=LE1 le;end a=linspace(0.5,10.5,300); plot(a,LE1,'-'); title('largest Lyapunov exponents of ouhe1'); xlabel('parameter a'),ylabe

10、l('largest Lyapunov exponents'); grid 雙參數(shù)空間的最大李雅譜指數(shù)譜程序：仍對(duì)以上耦合系統(tǒng)clear; %可調(diào)參數(shù)區(qū)間和步長(zhǎng)clc;global u aN1=linspace(0,0,200);N2=linspace(0,0,400);for I=1:200 u=1.5+I*0.025; d0=1e-8; for L=1:400 a=0.5+L*0.025; le=0; lsum=0; x=1; y=1; z=1; x1=1; y1=1; z1=1+d0; for i=1:150 T1,Y1=ode45('ouhe1',0,1

11、,x;y;z;u;a); T2,Y2=ode45('ouhe1',0,1,x1;y1;z1;u;a); n1=length(Y1); n2=length(Y2); x=Y1(n1,1); y=Y1(n1,2); z=Y1(n1,3); x1=Y2(n2,1); y1=Y2(n2,2); z1=Y2(n2,3); d1=sqrt(x-x1)2+(y-y1)2+(z-z1)2); x1=x+(d0/d1)*(x1-x); y1=y+(d0/d1)*(y1-y); z1=z+(d0/d1)*(z1-z); if i>50 lsum=lsum+log(d1/d0); end en

12、d le=lsum/(i-50); LE1(I,L)=le; N2(L)=a; end N1(I)=u;end X,Y=meshgrid(N1,N2); Z=LE1; pcolor(X,Y,Z); %畫偽彩圖 %colormap jet,shading interp %連續(xù)變化的變異飽和色圖，表面畫偽彩圖 %contourf(X,Y,Z)%畫等高線 title('largest Lyapunov exponents of ouhe1') xlabel('parameter itu') ylabel('parameter ita') zlabel(

13、'最大李雅普指數(shù)delta','FontSize',12)4.李雅譜指數(shù)程序：仍對(duì)以上耦合系統(tǒng)clear;clc;x=1;y=1;z=1;h=0.003;a=3.5;u=2;V=eye(3);S=V;b1=0;k=50000;for i=1:k x1=x+h*(-u*x+y*(z+a); y1=y+h*(-u*y+x*(z-a); z1=z+h*(z-x*y); x=x1;y=y1;z=z1; J=-u a+z y z-a -u x -y -x 1; J=eye(3)+h*J; B=J*V*S; V,S,U=svd(B); a_max=max(diag(S);

14、S=(1/a_max)*S; b1=b1+log(a_max);endLyapunov=(log(diag(S)+b1)/(k*h)李雅譜指數(shù)譜程序：仍對(duì)以上耦合系統(tǒng)clear;%奇異值分解法計(jì)算ouhe系統(tǒng)的李雅普諾夫指數(shù)譜clc;Z1=;Z2=;Z3=;x=1;y=1;z=1;h=0.002;u=2;%a=3;k=10000;for a=linspace(0.5,10.5,1000); V=eye(3); S=V; b1=0; lp=0; for i=1:k x1=x+h*(-u*x+y*(z+a); y1=y+h*(-u*y+x*(z-a); z1=z+h*(z-x*y); x=x1;y

15、=y1;z=z1; J=-u a+z y z-a -u x -y -x 1; J=eye(3)+h*J; B=J*V*S; V,S,U=svd(B); a_max=max(diag(S); S=(1/a_max)*S; b1=b1+log(a_max); end lp=(log(diag(S)+b1)/(k*h); Z1=Z1 lp(1); Z2=Z2 lp(2); Z3=Z3 lp(3); end a=linspace(0.5,10.5,1000); plot(a,Z1,'-',a,Z2,'-',a,Z3,'-'); title('Ly

16、apunov exponents of ouhe'); xlabel('parameter a'),ylabel('lyapunov exponents');grid on6.三維、二維以及時(shí)間序列圖：耦合系統(tǒng)的函數(shù)程序：function dx=ouhe(t,x)dx=zeros(3,1);a=3;u=2;dx(1)=-u*x(1)+x(2)*(x(3)+a);dx(2)=-u*x(2)+x(1)*(x(3)-a);dx(3)=x(3)-x(1)*x(2);圖程序：clc;clear;t,x=ode45('ouhe',0 1000,1 -

17、1 1);subplot(2,2,1); plot(x(:,1),'k','markersize',0.5);xlabel('t/ms','fontsize',12,'fontName','times new Roman','fontweight','bold','color','k');ylabel('x','fontsize',12,'fontName','times new

18、Roman','fontweight','bold','color','k');subplot(2,2,2); plot(x(:,2),'k','markersize',0.5);xlabel('t/ms','fontsize',12,'fontName','times new Roman','fontweight','bold','color','k');ylab

19、el('y','fontsize',12,'fontName','times new Roman','fontweight','bold','color','k');subplot(2,2,3); plot(x(:,3),'k','markersize',0.5);xlabel('t/ms','fontsize',12,'fontName','times new Roman'

20、;,'fontweight','bold','color','k');ylabel('z','fontsize',12,'fontName','times new Roman','fontweight','bold','color','k');subplot(2,2,4); plot3(x(:,1),x(:,2),x(:,3),'k','markersize',0.5);g

21、rid on;xlabel('x','fontsize',12,'fontName','times new Roman','fontweight','bold','color','k');ylabel('y','fontsize',12,'fontName','times new Roman','fontweight','bold','color','k');zlab