Rossler-混沌系统Lyapunov指数Matlab实现代码

Rossler混沌系统Lyapunov指数Matlab实现代码分类：转载2010-04-1110:49混沌系统的基本特点就是系统对初始值的极端敏感性，两个相差无几的初值所产生的轨迹，随着时间的推移按指数方式分离，lyapunov指数就是定量的描述这一现象的量。Lyapunov指数是衡量系统动力学特性的一个重要定量指标,它表征了系统在相空间中相邻轨道间收敛或发散的平均指数率。对于系统是否存在动力学混沌,可以从最大Lyapunov指数是否大于零非常直观的判断出来:一个正的Lyapunov指数,意味着在系统相空间中,无论初始两条轨线的间距多么小,其差别都会随着时间的演化而成指数率的增加以致达到无法预测,这就是混沌现象。Lyapunov指数的和表征了椭球体积的增长率或减小率,对Hamilton系统,Lyapunov指数的和为零;对耗散系统,Lyapunov指数的和为负。如果耗散系统的吸引子是一个不动点,那么所有的Lyapunov指数通常是负的。如果是一个简单的m维流形(m=1或m=2分别为一个曲线或一个面),那么,前m个Lyapunov指数是零,其余的Lyapunov指数为负。不管系统是不是耗散的,只要λ10就会出现混沌。微分动力系统Lyapunov指数的性质对于一维(单变量)情形,吸引子只可能是不动点(稳定定态)。此时λ是负的。对于二维情形,吸引子或者是不动点或者是极限环。对于不动点,任意方向的δxi,都要收缩,故这时两个Lyapunov指数都应该是负的,即对于不动点,(λ1,λ2)=(-,-)。至于极限环,如果取δxi始终是垂直于环线的方向,它一定要收缩,此时λ0;当取δxi沿轨道切线方向,它既不增大也不缩小,可以想像,这时λ=0。事实上,所有不终止于定点而又有界的轨道(或吸引子)都至少有一个Lyapunov指数等于零,它表示沿轨线的切线方向既无扩展又无收缩的趋势。所以极限环的Lyapunov指数是(λ1,λ2)=(0,-)。在三维情形下有(λ1,λ2,λ3)=(-,-,-):稳定不动点;(λ1,λ2,λ3)=(0,-,-):极限环;(λ1,λ2,λ3)=(0,0,-):二维环面;(λ1,λ2,λ3)=(+,+,0):不稳极限环;(λ1,λ2,λ3)=(+,0,0):不稳二维环面;(λ1,λ2,λ3)=(+,0,-):奇怪吸引子。李雅谱诺夫指数小于零，则意味着相邻点最终要靠拢合并成一点，这对应于稳定的不动点和周期运动;若指数大于零，则意味着相邻点最终要分离，这对应于轨道的局部不稳定，如果轨道还有整体的稳定因素（如整体有界、耗散、存在捕捉区域等），则在此作用下反复折叠并形成混沌吸引子。指数越大，说明混沌特性越明显，混沌程度越高.%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%Author:ThomasLeefunctiondX=Rossler_ly(t,X)%Rossler吸引子，用来计算Lyapunov指数%a=0.15,b=0.20,c=10.0%dx/dt=-y-z,%dy/dt=x+ay,%dz/dt=b+z(x-c),a=0.15;b=0.20;c=10.0;x=X(1);y=X(2);z=X(3);%Y的三个列向量为相互正交的单位向量Y=[X(4),X(7),X(10);X(5),X(8),X(11);X(6),X(9),X(12)];%输出向量的初始化，必不可少dX=zeros(12,1);%Rossler吸引子dX(1)=-y-z;dX(2)=x+a*y;dX(3)=b+z*(x-c);%Rossler吸引子的Jacobi矩阵Jaco=[0-1-1;1a0;z0x-c];dX(4:12)=Jaco*Y;求解LE代码：%计算Rossler吸引子的Lyapunov指数clear;yinit=[1,1,1];orthmatrix=[100;010;001];a=0.15;b=0.20;c=10.0;y=zeros(12,1);%初始化输入y(1:3)=yinit;y(4:12)=orthmatrix;tstart=0;%时间初始值tstep=1e-3;%时间步长wholetimes=1e5;%总的循环次数steps=10;%每次演化的步数iteratetimes=wholetimes/steps;%演化的次数mod=zeros(3,1);lp=zeros(3,1);%初始化三个Lyapunov指数Lyapunov1=zeros(iteratetimes,1);Lyapunov2=zeros(iteratetimes,1);Lyapunov3=zeros(iteratetimes,1);fori=1:iteratetimestspan=tstart:tstep:(tstart+tstep*steps);[T,Y]=ode45('Rossler_ly',tspan,y);%取积分得到的最后一个时刻的值y=Y(size(Y,1),:);%重新定义起始时刻tstart=tstart+tstep*steps;y0=[y(4)y(7)y(10);y(5)y(8)y(11);y(6)y(9)y(12)];%正交化y0=ThreeGS(y0);%取三个向量的模mod(1)=sqrt(y0(:,1)'*y0(:,1));mod(2)=sqrt(y0(:,2)'*y0(:,2));mod(3)=sqrt(y0(:,3)'*y0(:,3));y0(:,1)=y0(:,1)/mod(1);y0(:,2)=y0(:,2)/mod(2);y0(:,3)=y0(:,3)/mod(3);lp=lp+log(abs(mod));%三个Lyapunov指数Lyapunov1(i)=lp(1)/(tstart);Lyapunov2(i)=lp(2)/(tstart);Lyapunov3(i)=lp(3)/(tstart);y(4:12)=y0';end%作Lyapunov指数谱图i=1:iteratetimes;plot(i,Lyapunov1,i,Lyapunov2,i,Lyapunov3)程序中用到的ThreeGS程序如下：%G-S正交化functionA=ThreeGS(V)%V为3*3向量v1=V(:,1);v2=V(:,2);v3=V(:,3);a1=zeros(3,1);a2=zeros(3,1);a3=zeros(3,1);a1=v1;a2=v2-((a1'*v2)/(a1'*a1))*a1;a3=v3-((a1'*v3)/(a1'*a1))*a1-((a2'*v3)/(a2'*a2))*a2;A=[a1,a2,a3];计算得到的Rossler系统的LE为————0.0632310.092635-9.8924