切分段线性控制方法
% 连续时间系统混沌化切控制方法
% 《动力系统的混沌化》陈关荣 汪小帆
%
% Example(函数图象):
% [T,Y]=ode45('chaos1',40,[2;2;1]);
% plot3(Y(:,1),Y(:,2),Y(:,3));
%
function dx=chaos1(t,x);
a=3;
%global a;
b=20;c=-20;
k=4;d=10;
A= [ a b 0
-b a 0
0 0 c];
if x(3)+norm(x(1:2))>k
u=k*[-x(1) -x(2) d]';
else
u=0;
end
dx=A*x+u;
% 求最大Lyapunov指数谱
% 连续时间系统混沌化切控制方法
% 《动力系统的混沌化》陈关荣 汪小帆
%
Z=[];
global a;
d0=1e-8;
for a=linspace(0.1,6,150);
y0=[2;2;1];
y=[2+d0;2;1];
lsum=0;
p=0;
for k=1:60
[T,Y0] = ode45('chaos1',.4,y0);
[T,Y] = ode45('chaos1',.4,y);
y0 = Y0(size(Y0,1),:);
y = Y(size(Y,1),:);
d1=norm(y-y0);
if d1~=0
y=y0+(d0/d1)*(y-y0);
if k>10
lsum=lsum+log(d1/d0);
end
else
p=p+1;
y(1)=y0(1)+d0;
end
end
le=lsum/(k-10-p);
Z=[Z,a+le*i];
end
plot(Z,'-')
title('Lyapunov exponents')
xlabel('parameter a'),ylabel('Maximum lyapunov exponents')
grid on
5 comments:
今天在振动论坛看到,你做的真好,我看你的好多程序,谢谢你的帮助。不过在计算 Lorenz Lyapunov exponent的程序里相关连接 GS.m 总是连接不上,为什么?还有这个程序,可不可以用在其它连续系统上,来计算该系统的Lyapunov exponent
http://ustbchaos.blogspot.com/2006/10/gs.html
放在同一个文件夹下
我运行该程序,为什么lyaponov指数恒为0呢?
我运行时没问题,那我也不知道你那里是怎么回事. :(
Post a Comment