lozi.m

Lorenz系统的李雅谱诺夫指数

function newx=Lozi(x);
% Lozi方程[差分方程]
 %newx=Lozi([x;y;p;q]);
% 方程如下:
   x(k+1)=-p*x(k)+y(k)+1
  y(k+1)=q*x(k)
  X=[1;0;1.7;.5];Y=[];
 for i=1:1000
     X=feval(@Lozi,X);
      Y(i,:)=X(1:2,1);
  end
  plot(Y(:,1),Y(:,2),'.');
  title('Lozi映射图'),xlabel('x'),ylabel('y')%
% Example(分岔图):
Z=[];
  for p=linspace(0,1.7,300);
     x=[1;0;p;.5];
     for k=1:300;
          x=Lozi(x);
          if k>60
               Z=[Z,p+x(1)*i];
          end
       end
   end
   plot(Z,'.','markersize',1)
   title('Lozi映射分岔图'),xlabel('p'),ylabel('x')
%% Example(最大Lyapunov指数谱图):
%
   d0=1e-8;
   Z=[];
   for p=linspace(0,1.7,300)
       le=0;
       lsum=0;
       x=[1;0;p;.5];
       x1=[1;d0;p;.5];
       for k=1:800
                    x=Lozi(x);
           x1=Lozi(x1);
           d1=sqrt((x(1)-x1(1))^2+(x(2)-x1(2))^2);
           x1=x+(d0/d1)*(x1-x);
           if k>100
               lsum=lsum+log(d1/d0);
           end
       end
       le=lsum/(k-100);
       Z=[Z,p+le*i];
   end
   plot(Z,'-')
   title('Lozi最大Lyapunov指数图'),xlabel('p'),ylabel('lyapunov')
   grid on
%% Author：yujunjie
% Author's email: [e-mail]
newx(1,1)=-x(3)*abs(x(1))+x(2)+1;
newx(2,1)=x(4)*x(1);
newx(3,1)=x(3);
newx(4,1)=x(4);