声振论坛

 找回密码
 我要加入

QQ登录

只需一步,快速开始

查看: 3628|回复: 11

[分形与混沌] Lyapunov指数的计算

[复制链接]
发表于 2006-11-23 21:09 | 显示全部楼层 |阅读模式

马上注册,结交更多好友,享用更多功能,让你轻松玩转社区。

您需要 登录 才可以下载或查看,没有账号?我要加入

x
T=[1.49 2 1;-2 1.7 0;4 -4 2];
x=[x1;x2; x3];
f=-x+T*tanh(x);

这个系统的lyapunov指数怎么计算呢?
高手请指教一下,谢谢!

本帖被以下淘专辑推荐:

回复
分享到:

使用道具 举报

 楼主| 发表于 2006-11-23 21:36 | 显示全部楼层

不好意思,发错了

T=[1.49 2 1;-2 1.7 0;4 -4 2];
x=[x1;x2; x3];
x'=-x+T*tanh(x);
发表于 2006-11-24 03:02 | 显示全部楼层

LE结果

用程序计算的

[ 本帖最后由 xinyuxf 于 2006-12-23 09:50 编辑 ]
发表于 2006-11-24 10:34 | 显示全部楼层
楼上的大哥,有源码吗???
 楼主| 发表于 2006-11-24 21:51 | 显示全部楼层

回复 #3 tt1341 的帖子

能把源代码贴上来吗?
谢谢!
 楼主| 发表于 2006-11-24 22:11 | 显示全部楼层

哪位高手指点一下啊?

哪位高手把算法贴上来吧,急用!!!
谢谢!!!
发表于 2006-11-25 09:37 | 显示全部楼层
最讨厌那些张口就要源码的人了。
这些东西只有自己去亲自编写程序了,才能有更深的了解。
想当初我为了写一个程序,花费了两个月时间
 楼主| 发表于 2006-11-25 14:24 | 显示全部楼层

回复 #7 chinacer 的帖子

论坛是交流和学习的地方,不是水吧。
我是新手,看了论坛里的一些程序,不是很明白。
我只是希望高手能够帮忙解答一下。
 楼主| 发表于 2006-11-25 16:43 | 显示全部楼层

我的程序

function dX = chaos(t,X)
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);

% 吸引子
dX(1) = -X(1)+149/100*tanh(X(1))+2*tanh(X(2))+tanh(X(3))
dX(2) = -X(2)-2*tanh(X(1))+17/10*tanh(X(2))
dX(3) = -X(3)+4*tanh(X(1))-4*tanh(X(2))+2*tanh(X(3))

% 吸引子的Jacobi矩阵
Jaco = [ 49/100-149/100*tanh(X(1))^2  2-2*tanh(X(2))^2               1-tanh(X(3))^2;
        -2+2*tanh(X(1))^2             7/10-17/10*tanh(X(2))^2        0;
        4-4*tanh(X(1))^2              -4+4*tanh(X(2))^2              1-2*tanh(X(3))^2];
        
dX(4:12) = Jaco*Y;




yinit = [0,0,0.1];
orthmatrix = [1 0 0;
              0 1 0;
              0 0 1];

y = zeros(12,1);
% 初始化输入
y(1:3) = yinit;
y(4:12) = orthmatrix;

tstart = 0; % 时间初始值
tstep = 1e-3; % 时间步长
wholetimes = 1e4; % 总的循环次数
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);

for i=1:iteratetimes
    tspan = tstart:tstep:(tstart + tstep*steps);     
    [T,Y] = ode45('chaos', 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 = orth(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)
 楼主| 发表于 2006-11-25 16:45 | 显示全部楼层

哪里出错了?

算出来的结果是
lp =1.0e-013 *
   -0.3331
   -0.5040
   -0.4030
但正确结果应该是三楼给出的,哪里出错了呢?
发表于 2006-12-7 10:31 | 显示全部楼层
楼上的,你的问题解决了没?我最近也在做这个东东,你好像%正交化    y0 = orth(y0); 这个有点早了
发表于 2017-10-9 13:49 | 显示全部楼层
谢谢楼主分享
您需要登录后才可以回帖 登录 | 我要加入

本版积分规则

QQ|小黑屋|Archiver|手机版|联系我们|声振论坛

GMT+8, 2024-12-23 06:34 , Processed in 0.075705 second(s), 19 queries , Gzip On.

Powered by Discuz! X3.4

Copyright © 2001-2021, Tencent Cloud.

快速回复 返回顶部 返回列表