matlab计算Rossler吸引子的李雅普诺夫指数

function main

%求解LE代码：

% 计算Rossler吸引子的Lyapunov指数

% Rossler吸引子，用来计算Lyapunov指数

% a=0.25,b=1,c=5.5

% dx/dt = -y-z,

% dy/dt = x+ay,

% dz/dt = b+z(x-c),

clear;

yinit = [1,1,1];

orthmatrix = [1 0 0;

0 1 0;

0 0 1];

a = 0.25;

b = 1;

c = 5.5;

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);

for i=1:iteratetimes

tspan = tstart:tstep:(tstart + tstep*steps);

[T,Y] = ode45('Rossler_ly', tspan, y);

% 取积分得到的最后一个时刻的值

function main

%求解LE代码：

% 计算Rossler吸引子的Lyapunov指数

% Rossler吸引子，用来计算Lyapunov指数

% a=0.25,b=1,c=5.5

% dx/dt = -y-z,

% dy/dt = x+ay,

% dz/dt = b+z(x-c),

clear;

yinit = [1,1,1];

orthmatrix = [1 0 0;

0 1 0;

0 0 1];

a = 0.25;

b = 1;

c = 5.5;

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);

for i=1:iteratetimes

tspan = 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)

Le1=Lyapunov1(end),Le2=Lyapunov2(end),Le3=Lyapunov3(end)

title('Lyapunov指数谱图')

end

%计算得到的Rossler系统的LE为———— 0.0992 0.0755 -5.2644

%计算得到的Rossler系统的LE为———— 0.0992 0.0755 -5.2644

function dX = Rossler_ly(t,X)

a = 0.25;

b = 1;

c = 5.5;

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;

1 a 0;

z 0 x-c];

dX(4:12) = Jaco*Y;

%程序中用到的ThreeGS程序如下：

%G-S正交化

function A = 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];