微分混沌方程组最大Lyapunov指数算法

求Lorenz系统
在参数a=16，b=4，r=46.92下的最大Lyapunov指数

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('Lorenz',[0,1],[x;y;z;16;4;46.92]); [T2,Y2]=ode45('Lorenz',[0,1],[x1;y1;z1;16;4;46.92]); 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>100 lsum=lsum+log(d1/d0); end end le=lsum/(i-100)

le =
 
    1.4954

Chen系统和Rossler系统也可以用以上的算法只不过将Lorenz换成Chen或Rossler即可。