这是一个倍周期的分岔图代码
n0=100;n=200;
x0=0.5;s=[];xx=[];
for r=-1:0.0005:2;
x(1)=x0;
for j=2:n;
x(j)=-(1+r)*x(j-1)+x(j-1).^3;
end;
s=[r*ones(1,n+1-n0);s];
xx=[x(n0:n);xx];xs=max(x(n0:n));
end; figure;plot(s,xx,'b.','markersize',0.5);
hold on;
n0=100;n=200;
x0=-0.5;s1=[];xx1=[];
for r=-1:0.0005:2;
x(1)=x0;
for j=2:n;
x(j)=-(1+r)*x(j-1)+x(j-1).^3;
end;
s1=[r*ones(1,n+1-n0);s1];
xx1=[x(n0:n);xx1];xs1=max(x(n0:n));
end; plot(s1,xx1,'r.','markersize',0.5); |