CHEN'S系统的最大LE曲线

octopussheng

近期对最大LE随参数变化的曲线进行了一些研究，找到一些方法，算出CHEN'S系统的最大LE曲线，因为以前没有看过该系统的这种曲线，各位帮忙看看是否正确！

xiaoqiu810818

OCT师兄，附件是我从陆博士的混沌分析工具箱里找的求解L-指数的程序，希望对你有所帮助。弟水平实在差，只能帮到这里啦:loveliness:

xiaoqiu810818

另外。有了结果别忘了告诉弟，其它问题欢迎讨论。Email:  xiaoqiu810818@163.com

octopussheng

呵呵，好的，一定会的！

octopussheng

你文档里面的内容我看了，是基于时间序列求解最大LE的程序！和我这里做的有很大的区别呀！

shmily137

师哥,上面这个程序缺少函数啊!!!!你调试通了吗?你那还有求其他求lyapunov指数的程序吗?

yibanlixue

敢问这是用什么方法画的啊

liliangbiao

这个程序是一个博客上面的，我在这里引用一下！没有和作者沟通，吧程序贴上，希望原作者见谅！
Z=[];
d0=1e-8;
for a=linspace(32,40,80)
lsum=0;
x=1;y=1;z=1;
x1=1;y1=1;z1=1+d0;
for i=1:1000
   [T1,Y1]=ode45('Chen',1,[x;y;z;a;3;28]);
   [T2,Y2]=ode45('Chen',1,[x1;y1;z1;a;3;28]);
   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>500
       lsum=lsum+log(d1/d0);
   end
end
Z=[Z lsum/(i-500)];
end
a=linspace(32,40,80);
plot(a,Z,'-');
title('Chen 系统最大lyapunov指数')
xlabel('parameter a'),ylabel('lyapunov exponents')



Z=[];
d0=1e-8;
for a=linspace(32,40,80)
lsum=0;
x=1;y=1;z=1;
x1=1;y1=1;z1=1+d0;
for i=1:10000
   [T1,Y1]=ode45('Chen',1,[x;y;z;a;3;28]);
   [T2,Y2]=ode45('Chen',1,[x1;y1;z1;a;3;28]);
   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>5000
       lsum=lsum+log(d1/d0);
   end
end
Z=[Z lsum/(i-5000)];
end
a=linspace(32,40,80);
plot(a,Z,'-');
title('Chen 系统最大lyapunov指数')
xlabel('parameter a'),ylabel('lyapunov exponents')

octopussheng

循环的步长增大，对结果影响大不大？
譬如我取i为100，和取10000，差别究竟会多大呢？

liliangbiao

不大， 我已经试验过了！ 我不知道原因是什么 但是程序的设计思想是正确的！？ 但是对于系统来说 所作出的指数不都是正确的，不能正确的判断！

octopussheng

程序设计思想肯定是对的呀！就是按照最大Lyapunov指数的定义来做的！真实很奇怪的！

我用这种方法做过几个例子，在参考资料上的那些点，其实结果都蛮准确的！就是全局的不知如何计算啊！

baofeng

这是我用上述程序得到的最大指数图，和正确的不一样，程序如下，还请大家指点啊
global a b c
Z=[];
d0=1e-8;
for a=linspace(32,40,80)
lsum=0;
x(1)=1;x(2)=1;x(3)=1;
x1(1)=1;x1(2)=1;x1(3)=1+d0;
x0=[1,1,1];
x_0=[1,1,1+d0];
for i=1:100
   [T1,Y1]=ode45('Chen',1,x0);
   [T2,Y2]=ode45('Chen',1,x_0);
   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
end
Z=[Z lsum/(i-50)];
end
a=linspace(32,40,80);
plot(a,Z,'-');
title('Chen 系统最大lyapunov指数')
xlabel('parameter a'),ylabel('lyapunov exponents')

%Chen 系统
function xdot = chen(t,x)
%dx/dt=a*(y-x)
%dy/dt=(c-a)*x+c*y-x*z
%dz/dt=x*y-b*z
global a b c
b=3;c=28;
xdot=[ 0,0,0];
xdot(1)=a*(x(2)-x(1));
xdot(2)=(c-a)*x(1)+c*x(2)-x(1)*x(3);
xdot(3)=x(1)*x(2)-b*x(3);
%其中a,b,c是系统参数。当a=35, b=3, c=28 时，Lorenz 系统是混沌的。

lebusishu

各位高手
请帮助我
我现在再编写lyapunov指数程序，我参照各位的程序。
可我不明白（ [T1,Y1]=ode45('Chen',1,x0); ）在调用这句话中，'chen'后的‘1’是什么意思？？不应该是时间变量吗？？
请指教！！！

咕噜噜

你参照的谁的程序，这里1肯定是不正确的，应该是一个时间段

ivy_1031

原帖由 octopussheng 于 2008-4-20 09:17 发表

                               
登录/注册后可看大图

循环的步长增大，对结果影响大不大？
譬如我取i为100，和取10000，差别究竟会多大呢？

关系是否大不能一概而论，这个跟系统本身的性质有关!一般绘制最好跟系统的分岔图(与该指数图的参数一致)进行对比！来观察结果是否一致！如果不一致很有可能是跟步长或者总的步数取值有关！很难做到同一套参数适合所有系统的。

建议：１、关于总步数：如果可以绘制每个李指数随时间变化的图，然后观察李指数大概在何时趋于平稳，然后根据这个时间选择总的步数！
　　　２、关于步长step：这个跟你所需的效果有关！一般计算的参数区域较大的话，建议step可以稍微取的大点，这样看起来较为“柔和”（自己找个了词形容，其实就是抖动的不会太剧烈）；如果看局部具体的现象，建议取小步长，以免遗漏！

总之，李指数图的绘制不是一个方法那里都适用的，毕竟是数值方法。所以要小心！可以参考，不能完全依赖！