lyapunov_exponent.m

来自「非线性动力学中lyapunov、duffing等方程的求解」· M 代码 · 共 21 行

clear;

%computes lyapunov exponent 
r = 4.0;

%%%%% first converge to attractor
x = rand;   %choose random initial data
for i=1:200;
   x = r*x*(1-x);     %apply logistic map
end;

%%%%% find (1/n) Sum ln(f'(x_i))  -  limit should be lyapunov exponent
sum = 0;
for i=1:500;
    x = r*x*(1-x);     %logistic map
    sum = sum + log(abs(r - 2*r*x));
    sum/i                                    %print out average
    average(i) = sum/i;             %collect into vector for later plotting
end;

plot(average);      %plots sucessive averages so a trend can be observed