Homework9

李尧 2013301020048 天眷班
习题3.16

混沌 分形 非线性谐振子

摘要
虽然经典力学是确定性的逻辑，牛顿第二定律是确定性的方程，但是确定性的方程也有不确定的解，这些不确定的方程描述的物理现象被称之为混沌。进一步的研究发现，混沌的相空间的吸引子有分形结构。本文将以非线性谐振子为例详细地讨论混动现象。

背景
1963年美国气象学家洛伦兹在研究天气的过程中发现了混沌现象。这意味着由于解得不稳定性，即使描述运动的动力学方程是确定的，物理现象仍然不可预测。实际上，混沌现象在我们生活中并不少见，例如湍流。如下图所示

简单的烟雾从喷嘴喷出后，先是可预测的层流，然后是完全无法预测的湍流。人们在研究混沌的过程中发现了混沌与分形的联系，简单地说，混沌的吸引子具有分形的特征，即自相似。实际上，混沌可以说是时空间的分形，分形可以说是空间的混沌。如烟雾中的大涡旋里套着小涡旋，小涡旋套着更小的涡旋。

物理模型与计算方法
考虑一个受到驱动和摩擦的物理摆，他的运动方程由下式给出：
$\frac{d^2x}{dt^2}=-\frac{g}{l}sin(\theta)-g\frac{d\theta}{dt}+F_Dsin(\Omega_Dt)$
我们把这个一元二阶常微分方程变成等价的二元一阶常微分方程组：
$\frac{d\omega}{dt}=-\frac{g}{l}sin(\theta)-g\omega+F_Dsin(\Omega_Dt)$
$\frac{d\theta}{dt}=\omega$
这个方程组是不可解析求解的，但是可以用Euler-Cromer方法数值求解。这样我实际上使用的是以下的迭代关系式：
$\omega_{i+1}=\omega_i+\Delta t*（-\frac{g}{l}sin(\theta_{i+1})-g\omega_{i}+F_Dsin(\Omega_Dt)）$
$\theta_{i+1}=\theta_{i}+\Delta t*\omega_{i+1}$

根据以上的分析，模拟实验的代码如下：
program

解的稳定性分析
我们设满足上式的解为。
$\theta_{i}$
$\omega_i$
现在考虑在解附近有微小的$\delta$的偏离，由于$\delta$很小，我们可以只保留一阶项。方程改写为
$\delta\omega_{i+1}=\delta\omega_i+\Delta t*(-\frac{g}{l}cos(\theta)\delta\omega_i-g\delta\omega_{i})$
$\delta\theta_{i+1}=\delta\theta_i+\Delta t*\delta\omega_{i+1}$
我们假设$\Delta t$已经取得非常的小了，不必考虑算法本身的稳定性，只考虑数值解的稳定性。这是一个线性方程，可以改写成矩阵的形式（略），不妨设对应的迭代矩阵的特征根为$M_1$和$M_2$.那么解的稳定性充分条件可以写成
$-1<M_1,M_2<1$

结果与讨论
首先我们令g=9.8 l=9.8 $\Omega_D=\frac{2}{3}$ $\Delta t=0.0001$ 总的步数N=5000000 ，初始角度和初始角速度分别为0.2rad和0rad/s。调节$F_D$,我们发现混沌出现在$F_D=1.2$ 附近很小的区间类。

由图可知，非线性谐振子震荡的非常的没有规律，完全不可预测。可以画出对应的相图，如下：

相应的吸引子，如下：

上图的右下角有明显的分形的特征。经程序计算得其信息维数为0.9429.这是由于我们计算得到的图案不是数学意义上的分形，仅仅是起示意的作用，所以最终的结果接近其拓扑维数1.

接下来我们稍微改变一下初值条件，令$\theta_0=0.3$ 得到新的吸引子：

两张图的形状几乎不变，但是个别的点在曲线上的分布有明显的不同。为了直观的反应这一微小变化的影响，我们再画出$t-\theta$的图像。

相比第一幅图，一开始很接近，但是后来的行为完全不同。

结论
我们可以想象但无穷次的运行程序的话，两种初值条件下的吸引子完全相同。因此吸引子的形状与初值条件无关，即使不同初值条件的动力学行为完全不同，但是初值条件影响到后来的动力学的行为，也就是决定了吸引子上的点在时空间出现的顺序。

此外，出于兴趣。我还画了$F_D$从1.35到1.55的分插图：

致谢
王世兴同学为我的工作提供了必要的帮助。
《分形理论与应用》科学教育出版社