大作业：洛伦兹系统参数空间探究

源材料包括一份实验报告、一个 Notebook 和四张输出图。本篇保留能够由方程或代码验证的内容，重新绘制关键图，并把原报告与 Notebook 不一致的地方逐项校正。个人信息和联系方式不进入博客。

研究对象

洛伦兹系统为

\begin{cases} \dot x=\sigma(y-x),\\ \dot y=x(r-z)-y,\\ \dot z=xy-\beta z, \end{cases} \qquad \sigma,r,\beta>0.

三个状态量最初来自简化的大气对流模型：

$x$ 表示对流强度；
$y$ 表示上升流与下降流的温差；
$z$ 表示垂直温度分布相对平衡态的偏差。

$\sigma$ 控制两种扩散时间尺度的比值， $r$ 与加热强度有关， $\beta$ 与几何尺度有关。大作业的目标不是证明“某个参数越大就一定越混沌”，而是回答：

不动点何时存在、何时稳定？
哪些参数区域必须靠数值实验继续判断？
数值图怎样做才不把暂态、采样不足或算法误差当成混沌？

先做完能严格推导的部分

不动点

原点

P_0=(0,0,0)

始终存在。当 $r>1$ 时，还有一对对称不动点

P_\pm= \left( \pm\sqrt{\beta(r-1)}, \pm\sqrt{\beta(r-1)}, r-1 \right).

系统在变换

(x,y,z)\longmapsto(-x,-y,z)

下不变，所以非原点平衡和吸引子的两翼总是成对出现。

原点的稳定性

雅可比矩阵为

J(x,y,z)= \begin{pmatrix} -\sigma&\sigma&0\\ r-z&-1&-x\\ y&x&-\beta \end{pmatrix}.

在原点，

\det(\lambda I-J) =(\lambda+\beta) \left[(\lambda+\sigma)(\lambda+1)-\sigma r\right].

其中一个特征值为 $-\beta$ ，另外两个满足

\lambda^2+(\sigma+1)\lambda+\sigma(1-r)=0.

因此：

$0<r<1$ ：三个特征值实部都为负，原点稳定；
$r=1$ ：一个特征值过零；
$r>1$ ：原点有正特征值，变成鞍点，同时 $P_\pm$ 产生。

这就是 $r=1$ 的叉式分岔。

对称不动点的稳定边界

在 $P_\pm$ ，特征多项式为

\lambda^3 +(\sigma+\beta+1)\lambda^2 +\beta(\sigma+r)\lambda +2\sigma\beta(r-1)=0.

对三次多项式

\lambda^3+a_1\lambda^2+a_2\lambda+a_3,

Routh–Hurwitz 条件是

a_1>0,\quad a_2>0,\quad a_3>0,\quad a_1a_2>a_3.

当 $r>1$ 时前三项自动满足，最后一项化为

\sigma(\sigma+\beta+3) +r(\beta+1-\sigma)>0.

若 $\sigma>\beta+1$ ，稳定性边界是

\boxed{ r_H= \frac{\sigma(\sigma+\beta+3)} {\sigma-\beta-1} }.

在经典参数 $\sigma=10,\ \beta=8/3$ 下，

r_H\approx24.7368.

$1<r<r_H$ 时 $P_\pm$ 稳定；越过 $r_H$ 后它们失稳。经典参数下这里是亚临界 Hopf 分岔。若 $\sigma\le\beta+1$ ，上面的分母非正，不能继续套“ $r<r_H$ ”这句话；此时 Routh–Hurwitz 不等式对所有 $r>1$ 都成立。

耗散性

散度恒为

\nabla\cdot\boldsymbol f =-\sigma-1-\beta<0.

一个无穷小体积元 $V(t)$ 满足

V(t)=V(0)e^{-(\sigma+1+\beta)t}.

所以相空间体积指数收缩。要注意：耗散只说明体积收缩，不等于自动证明混沌，也不能单独证明某条轨道有界。

经典参数下的吸引子

下面的图不是报告截图，而是重新用固定步长 RK4、 $h=0.005$ 积分，并删去 $t<15$ 的暂态后绘制。这里故意使用 RK4，是为了让图题、代码方法和步长口径完全一致。

洛伦兹吸引子的 x-z 与 x-y 投影

轨迹在两个翼之间不规则切换。单看“蝴蝶形”只能说它像经典洛伦兹吸引子；要把“混沌”写成数值结论，还需要收敛测试和可信的最大李雅普诺夫指数。

原报告的积分方法与 Notebook 实际不一致

原报告写的是“四阶 Runge–Kutta，步长 $h=0.01$ ”。Notebook 实际调用：

solve_ivp(
    lorenz_system,
    [0, T],
    initial_state,
    t_eval=t_eval,
    method="RK45",
    rtol=1e-8,
    atol=1e-10,
)

这意味着：

实际方法是自适应 Dormand–Prince RK45，不是固定步长 RK4；
t_eval 相邻 $0.01$ 只是输出采样间隔，不是求解器内部步长；
rtol 与 atol 控制局部误差，内部步长会自动变化；
Notebook 删除前 $50\%$ 数据作为暂态，这一点与代码一致。

两种积分方法都能用，但报告必须如实写。最基本的数值自检是把容差收紧或把 RK4 步长减半，再比较统计量，而不是要求两条混沌轨迹逐点重合。

参数空间：先画解析边界，再谈混沌

下面这张图只依据 Routh–Hurwitz 条件，展示 $\sigma=10$ 时平衡点的局部稳定性。它不是李雅普诺夫热图，也没有把“平衡点失稳”直接涂成“混沌”。

洛伦兹系统解析平衡点稳定区域

从图中能严格读出的只有：

$r<1$ 时原点稳定；
$r>1$ 后看 $P_\pm$ ；
在橙色区域 $P_\pm$ 已失稳，但最终是周期轨、混沌吸引子还是其他行为，仍需数值或更深理论判断。

原 Notebook 的二维扫描只有

r\in[0,50]\ \text{的 10 个点}, \qquad \beta\in[1,10]\ \text{的 10 个点}.

也就是 $10\times10$ 网格。它适合检查代码流程，不足以支持“共振细带”“混沌与稳定区域交替”等精细结论。正式扫描至少应逐级加密网格，并验证分类边界是否收敛。

有限时间极值诊断：不能冒充渐近分岔图

Notebook 的原分岔程序只扫描 $r\in[0,50]$ ，共 51 个参数点，却又标出了 $99.5$ — $100.8$ 与 $145$ — $166$ 的窗口；这两段既在扫描范围外，也被横轴范围 $[0,50]$ 隐藏，不能作为该次实验的发现。

下面重新在 $[0,50]$ 内取 260 个 $r$ 值，对每条轨迹积分到 $t=120$ ，丢弃 $t<65$ 的部分，再保留 $x(t)$ 的严格局部极大值。图中只标记解析上确定的 $r=1$ 与 $r_H\approx24.74$ 。

洛伦兹系统在 r 从 0 到 50 时的有限时间局部极大值诊断图

读图方法：

定常平衡解没有严格局部极大值，因此不会在这种取点规则下留下“平衡支”；
$1<r<r_H$ 内的孤立点和负值支，是轨迹向 $P_\pm$ 衰减时尚未消失的有限时间暂态极值，也记录了落入哪个对称吸引域；它们不是渐近平衡分支，更不能解释成“一个极大值对应一个不动点”；
暂态之后仍出现有限几条稳定分支，才可能提示周期轨及倍周期结构；
出现密集点带，提示非周期行为，但仍不能仅凭“点很多”证明混沌。

这张图只是有限时间、有限步长的暂态诊断。只有在延长积分时间、增加丢弃暂态、加密参数网格并改变极值检测精度后仍保持的结构，才值得进一步作为渐近周期或混沌的证据。

原最大李雅普诺夫指数为什么不能直接采用

Notebook 中的最大指数函数同时积分两条初始距离为 $d_0$ 的轨迹，然后每隔一段时间读取距离 $d(t)$ ，累计

\log\frac{d(t)}{d_0}.

问题在于它从未把扰动重新缩回 $d_0$ ：

初期指数分离后，两轨距离会达到吸引子尺度并饱和；
后续采样不再处于切空间的线性区；
代码把多个“从初始时刻累计到现在”的距离重复平均；
因而结果不是严格的 Wolf、Benettin 或 QR 李雅普诺夫算法。

可靠的最大指数算法应周期重归一化。设扰动向量为 $\boldsymbol v$ ：

同时积分轨迹 $\dot{\boldsymbol x}=\boldsymbol f(\boldsymbol x)$ 与变分方程
$\dot{\boldsymbol v}=J(\boldsymbol x)\boldsymbol v;$
每隔 $\Delta t$ 记录 $q_k=\|\boldsymbol v\|/d_0$ ；
令 $\boldsymbol v\leftarrow d_0\boldsymbol v/\|\boldsymbol v\|$ ；
最后计算
$\lambda_{\max} \approx \frac{\sum_k\ln q_k}{\sum_k\Delta t}.$

若要全部三个指数，就同时积分一个 $3\times3$ 扰动矩阵，并在每个时间段做 QR 分解：

M=QR,\qquad \lambda_i\approx\frac1T\sum_k\ln|R_{ii}^{(k)}|.

还要检查延长总时间、缩短重归一化间隔、收紧积分容差后指数是否收敛。

“李雅普诺夫指数谱”那一格没有完整可靠的执行记录

Notebook 最后一格的 execution_count 为 null，但文件仍残留一段截断的标准输出，最后停在 $r=52.6$ 。这说明保存下来的执行状态与输出并不一致，不能把它当作一次完整、可靠的运行记录。即使把该格完整执行，它也只调用了一个最大指数函数，然后手工指定：

混沌时令 $\lambda_2=0$ ；
“周期”时令 $\lambda_2=-5$ ；
“稳定”时令 $\lambda_2=-2$ ；
再用指数和约束反推 $\lambda_3$ 。

所以那三条曲线不是计算得到的李雅普诺夫谱。恒等式

\lambda_1+\lambda_2+\lambda_3 =-\sigma-1-\beta

只能约束总和，不能凭一个数唯一恢复另外两个数。手工填的 $-5$ 、 $-2$ 不能作为结论。

哪些原结论应保留，哪些应撤回

结论	状态	理由
不动点、雅可比矩阵、 $r=1$ 叉式分岔	保留	可解析推导
$r_H$ 的 Routh–Hurwitz 边界	保留	条件与代数可核验
对称性与负散度	保留	直接代入方程可证
经典参数出现蝴蝶形吸引子	作为数值现象保留	可由独立积分复现
10×10 热图中的细小“共振区域”	撤回	分辨率不足且指数算法有误
“最大指数达到 1.5 以上”	撤回	原算法不可靠
$r>313$ 再次有序	本实验不作结论	扫描只到 $50$
机器学习能够分类参数区域	撤回	Notebook 中没有训练、验证或分类代码
三条李雅普诺夫指数谱	撤回	两条由人工赋值，且保存的执行状态与截断输出不一致

一套可复现的重做流程

解析层：先计算不动点、Routh–Hurwitz 边界、对称性和散度；
单轨层：用两套步长或容差复现时间序列与相图；
分岔层：逐级加密 $r$ 网格，统一暂态和极值检测规则；
指数层：用变分方程加周期重归一化或 QR 算法；
分类层：至少综合“是否收敛到平衡点、极值个数、最大指数”三个证据；
收敛层：延长积分时间、加密参数网格，报告分类发生变化的点；
表述层：解析结论、数值证据和推测分开写，不把图像外观当定理。

这个大作业真正有价值的地方，不是画出一只熟悉的“蝴蝶”，而是学会给数值结论划边界：代码算了什么、没算什么，以及还要通过哪些检查才能把“看起来像”升级为“证据支持”。