M084 · 数值积分的几何

Hamilton 相空间与积分器实验室

模拟轨道时能量误差很小，轨迹也可能是错的吗？ 完全可能。辛积分器按 Hamilton 方程的几何结构设计，能量误差可以长时间保持很小、有界，算出的物体却仍可能沿轨道越跑越超前或落后，偏离真实轨迹。能量表现、精度、辛性、时间可逆性是四回事，要分别检查。

并排运行 Euler、RK4、辛 Euler 和 Verlet。比较谁长期更贴近真实相图或轨道，能量怎样漂移，每一步映射是否保住局部正则几何、能否倒着走回去。最后把步长减半，挑战运动变化很快的开普勒近日点。

回放无量纲时间 $\text{Pending}$

实际恒定步长 $\text{Pending}$

每种方法的总步数 $\text{Pending}$

选中方法 $\text{Pending}$

当前相对能量误差 $\text{Pending}$

最大绝对能量误差 $\text{Pending}$

当前状态误差 $\text{Pending}$

末个保留状态误差 $\text{Pending}$

局部完整映射行列式 $\text{Pending}$

局部辛性误差 $\text{Pending}$

单步反演误差 $\text{Pending}$

参考构造 $\text{Pending}$

参考细化差异 $\text{Pending}$

选中方法记录 $\text{Pending}$

保留记录结束时间 $\text{Pending}$

正则维数 $\text{Pending}$

改变设置，比较轨迹

模拟轨道时能量误差很小，轨迹也可能是错的吗？

系统：振子 / 摆 / 开普勒

\mathcal M

0.00\,

初始坐标

q_0

1.00\,

初始动量

p_0

0.00\,

初始开普勒偏心率

e

0.40\,

请求步长比例

h/T_0

0.02\,

标称周期数

T/T_0

20.00\,

方法：Euler / RK4 / 辛 Euler / Verlet

\mathcal I

3.00\,

正则斑块半径

\rho

0.03\,

回放时间

t/T

0.00\,\%

拖动橙色初态环；对应的坐标和动量滑块可以用键盘调。开普勒模式下，圆环改的是偏心率，振子的初态控件不起作用。所有记录都预先算好，播放只是回看保存的采样。右下图把每个斑块的中心移到原点。开普勒斑块只是一个投影；行列式和辛性误差都按完整的四维正则映射计算。记录一旦发散或接近奇点就停止，并保留标记。

三个系统都是无量纲、单位质量的保守系统：谐振子、理想摆，以及在平面上沿开普勒轨道运动的质点。四种方法步长相同且恒定。标称周期按谐振子或开普勒的参考周期计算，不是有限振幅下摆的周期。振子和开普勒的参考是解析解；摆的参考用 RK4 在两种更细的步长下计算，两者的差别只是估计，不是严格的误差界。不模拟自适应步长、阻尼或物理碰撞。

先抓重点

核心问题: 模拟轨道时能量误差很小，轨迹也可能是错的吗？
一句话直觉: 精度、能量表现、辛性和可逆性彼此相关，但各自测的是不同的东西。

核心数学模型

三个系统的 Hamilton 方程

\dot{\mathbf q}=\mathbf p,\qquad\dot{\mathbf p}=-\nabla V,\qquad H=\frac{|\mathbf p|^2}{2}+V

取单位质量和无量纲单位，三个保守系统可以直接比较；谐振频率和开普勒引力参数都取一。

先踢再漂移

\mathbf p_{j+1}=\mathbf p_j-h\nabla V(\mathbf q_j),\qquad\mathbf q_{j+1}=\mathbf q_j+h\mathbf p_{j+1}

辛 Euler 是一阶方法。它和显式 Euler 的区别在于，漂移时用的是刚更新过的动量。

对称 Verlet 步

\mathbf p_{j+1/2}=\mathbf p_j-\frac h2\nabla V(\mathbf q_j),\quad\mathbf q_{j+1}=\mathbf q_j+h\mathbf p_{j+1/2},\quad\mathbf p_{j+1}=\mathbf p_{j+1/2}-\frac h2\nabla V(\mathbf q_{j+1})

把一次踢拆成前后两个半踢，得到二阶、可逆的映射。但真实 Hamilton 量并不因此精确守恒；近日点那段步长太粗，照样算不准。

检验步进映射是否保辛

D\Phi_h^{\mathsf T}J D\Phi_h=J,\qquad J=\begin{pmatrix}0&I\\-I&0\end{pmatrix}

读数用中心差分估算，给出矩阵各元素中最大的偏差。在四维里，只看行列式等于一是更弱的检验。开普勒斑块只是二维投影，不能当作完整的体积测量。

能量误差和状态误差是两回事

\epsilon_H=\frac{H_j-H_0}{\max(|H_0|,10^{-6})},\qquad\epsilon_z=\|\mathbf z_j-\mathbf z_{\rm ref}(t_j)\|_2

最大能量误差统计全部步；绘图时每种方法保留约 600 个样本。能量曲线平稳，相位误差仍可能越积越大。记录一旦停止，就不再往后预测。

怎样读对数坐标

y_H=\operatorname{sgn}(\epsilon_H)\log_{10}\left(1+\frac{|\epsilon_H|}{10^{-8}}\right),\qquad y_z=\log_{10}\left(1+\frac{\epsilon_z}{10^{-8}}\right)

同一种对数变换让微小误差和巨大误差同时可见，还保留能量漂移的正负。摆的参考解把每个粗步分成 32 个和 64 个 RK4 子步，再比较两者。标称周期数按缩放时间里每周期二倍圆周率计算。

常见难点

能量准不等于轨道准

容易误解能量误差很小，轨道就一定对。

正确理解能量几乎守恒时，相位误差照样会积累。要和同一时刻的参考状态比。

保辛不等于能量精确守恒

容易误解辛方法每一步都精确守恒 Hamilton 量。

正确理解步长足够小、能分辨运动时，能量误差通常在有限范围内来回振荡，但真实能量并没有精确守恒。

近日点仍要小步长

容易误解只要保住几何结构，大步长也没关系。

正确理解靠近中心的那段运动很快，仍然需要短步长。把步长减半，比较轨道误差，不要只看行列式。

投影面积不是相空间体积

容易误解开普勒投影斑块的面积不变，就证明了辛性。

正确理解开普勒问题的相空间是四维的，完整矩阵检验才给出相关的局部证据。

动手实验

01
比较全部方法
回看完整的谐振记录，切换选中的方法。
你会看到： 各种颜色用的是同一步长，能量表现却不同：有的增长，有的漂移，有的在有限范围内振荡。
02
步长减半
用长期 RK4 记录，再把恒定步长减半。
你会看到： 状态误差和能量漂移都随步长变化；局部阶数高，并不是唯一的标准。
03
离开小角度极限
比较大振幅往返摆和旋转摆两个预设。
你会看到： 这里的参考解是细化后的数值解，没有用谐振近似代替。
04
挑战开普勒近日点
比较已解析和未解析的近日点，检查完整映射误差，再把步长减半。
你会看到： 辛映射可以保住自己的几何结构，算出的轨迹却仍然很差。