取余运算的Lyapunov指数

一、取余运算

1、 画出取余运算的运动轨迹

N=100;              %给定迭代次数
x=ones(1,N)*0.6;    %对x赋初值
for i=2:N
x(i)=mod(2*x(i-1),1);
end
plot(x(2:N))        %从第二个点开始画图
xlabel('\fontsize{16}n')
ylabel('\fontsize{16}x')

问题：在50多步的时候，输出的结果与我们想象的有出入，并且导致之后的结果全为0。

原因：下图是x变量值的变化：



从表中可以发现在第43步时，已经产生了错误的结果，继而导致最终为0.

错误结果的产生与matlab中小数的存储方式有关

解决方案：为了消除有存储方式带来的误差，可以通过两种方式解决这个问题：

1、 通过使用取整函数roundn，对x进行百分位取整

2、 将数据扩大一个数量级，即令x为1-9之间的整数，在取余时，对10取余。

N=100;              %给定迭代次数
x=ones(1,N)*0.6;    %对x赋初值
for i=2:N
x(i)=roundn(mod(2*x(i-1),1),-2);%利用roundn函数取整到x的百分位
end
plot(x(2:N))        %从第二个点开始画图
xlabel('\fontsize{16}n')
ylabel('\fontsize{16}x')

2、计算Lyapunov指数（仿真方式）

取N=100,x0=0.6,∆x=0.1,则x1=0.7,

N=100;              %给定迭代次数
delta=0.1;
x0=ones(1,N)*0.1;    %对x赋初值
x1=ones(1,N)*(0.1+delta);
for i=2:N
x0(i)=roundn(mod(2*x0(i-1),1),-2);%利用roundn函数取整到x的百分位
x1(i)=roundn(mod(2*x1(i-1),1),-2);
end
L=log(abs(x1(N)-x0(N))/delta)/(N-1)

L =
0.0070

当N取的很大时，Lyapunov指数越来越小。

二、用Lyapunov第一方法求局部稳定性

与第一题无关

syms x y
f=[y-x*(x^2+y^2);-x-y*(x^2+y^2)];
v=[x,y];
R=jacobian(f, v); 				 %求雅克比矩阵
R=subs(R,[x,y],[0,0]);
eig(R)						 %雅克比矩阵的特征值

ans =
-i
i

特征值是一对实部为0的共轭复根，所以系统平衡状态0在李雅普诺夫意义下稳定。