tent映射的Lyapunov指数图
谁有tent映射的Lyapunov指数图和分岔图的程序啊？？？方程见附件tent映射序列程序
clear
x0=0.8;
y0=0.83;
u=1.7;
e=0;
t=50;
n=100;
x=zeros(1,n);
y=zeros(1,n);
for i=1:n
if (x0>=0)&(x0<=0.5)
x(i)=u*x0;
elseif (x0>0.5)&(x0<=1)
x(i)=u*(1-x0);
end
x0=x(i);
end
for i=1:n
if (y0>=0)&(y0<=0.5)
y(i)=u*y0;
elseif (y0>0.5)&(y0<=1)
y(i)=u*(1-y0);
end
y0=y(i);
end
for j=1:n
Ex=e+x(j);
e=Ex;
end
e=0;
for k=1:n-t
s(k)=x(k)*x(k+t);
r=e+s(k);
e=r;
end
EX=Ex/n
Rx=r/(n-t)
e=0;
for k=1:n-t
s(k)=x(k)*y(k+t);
r=e+s(k);
e=r;
end
Rxy=r/(n-t)
hold on
plot(x)
plot(y,'r')
tent映射图程序：
% Chapter 2 - Nonlinear Discrete Dynamical Systems.
% Program 2a - Graphical Iteration of the Tent Map.
% Symbolic Math toolbox required.
% Copyright Birkhauser 2013. Stephen Lynch.
% Plot a cobweb diagram (Figure 2.7(b)).
clear
% Initial condition 0.2001, must be symbolic.
nmax=200;
t=sym(zeros(1,nmax));t1=sym(zeros(1,nmax));t2=sym(zeros(1,nmax));
t(1)=sym(2001/10000);
mu=2;
halfm=nmax/2;
axis();
for n=2:nmax
if (double(t(n-1)))>0 && (double(t(n-1)))<=1/2
t(n)=sym(2*t(n-1));
else
if (double(t(n-1)))<1
t(n)=sym(2*(1-t(n-1)));
end
end
end
for n=1:halfm
t1(2*n-1)=t(n);
t1(2*n)=t(n);
end
t2(1)=0;t2(2)=double(t(2));
for n=2:halfm
t2(2*n-1)=double(t(n));
t2(2*n)=double(t(n+1));
end
hold on
fsize=20;
plot(double(t1),double(t2),'r');
x=;y=;
plot(x,y,'b');
x=;y=;
plot(x,y,'g');
title('Graphical iteration for the tent map','FontSize',fsize)
set(gca,'XTick',0:0.2:1,'FontSize',fsize)
set(gca,'YTick',0:0.2:1,'FontSize',fsize)
xlabel('x','FontSize',fsize)
ylabel('T','FontSize',fsize)
hold off
% End of Program 2a.
tent映射分岔图：
clear
x=0.3;
u=0.9:0.001:2;
n=length(u);
figure(2)
hold on;
for i=1:n
a=x;
for k=1:300
%f(i)=u(i)*a*(1-a);
if (a>=0)&(a<=0.5)
f(i)=u(i)*a;
a=f(i);
elseif (a>0.5)&(a<=1)
f(i)=u(i)*(1-a);
a=f(i);
end
end
z(1)=a;
for j=1:99
% z(j+1)=u(i)*z(j)*(1-z(j));
if (z(j)>=0)&(z(j)<=0.5)
z(j+1)=u(i)*z(j);
elseif (z(j)>0.5)&(z(j)<=1)
z(j+1)=u(i)*(1-z(j));
end
end
plot(u(i)*ones(1,100),z,'.','markersize',4);
end
hold off
Lyapunov指数图没有现成的，你可以在上述基础上自己做，Lyapunov指数的计算方法论坛有很多，自己搜索参考
页:
[1]