2.1 运动方程与十八维相空间：把对手请上台

2.1 运动方程与十八维相空间：把对手请上台

本节摘要：三体问题的全部困难浓缩在一组方程里：三个天体各受另外两个的平方反比引力，九个二阶方程、十八个一阶方程，形式上与二体一样工整。本节把方程从向量式改写到哈密顿正则形式，解释为什么这一步是清点守恒量的会计准则，并用一段受力代码做一次完整的数值自检——它是第 4 章积分器的心脏零件。

方程到底哪里难？把三个天体的运动方程与两个天体的并排放在一起看，你几乎找不到形式上的差别——同样的平方反比、同样的光滑性、同样的确定性。既然差别不在"长相"，那就必须进到结构里去找。本节的任务是把方程写到最能暴露结构的形式为止。

图 2-1 三体的两两引力：六个力向量与成对抵消

为什么值得先把方程抄整齐

接上一节的二体推导，先补全三体的原貌。设三个质点位置为 \mathbf{r}_1, \mathbf{r}_2, \mathbf{r}_3，质量 m_1, m_2, m_3，第 i 个天体受到的引力是另外两体贡献之和：

三个向量方程，展开是九个二阶常微分方程；换成位置加速度变量后是一阶方程组，共十八个。与二体对比，真正的变化只有一处：每一个天体的"引力源"从一个变成两个，而且两个源头自己也在动。二体能化约成单体问题，靠的是引力源唯一且化成固定中心；三体这里，相对运动方程里的两项分别带着两个移动的源头，再也没法合并。

这处变化为什么是毁灭性的？因为二体的解法把"两体"化成"一体"用了三层杠杆——质心分离、角动量定面、径向方程凑成谐振子——而这三层杠杆全都不许引力源移动。方程只差一项，解的结构却整个坍塌，这是三体问题给后来者的第一课：非线性系统的困难常常不在方程的长度，而在方程里各项能不能被守恒律"对消"。

哈密顿形式：给守恒量立账的会计准则

十九世纪的力学改写给了这个问题一副更透亮的眼镜。定义广义动量 \mathbf{p}_i = m_i\dot{\mathbf{r}}_i，系统由哈密顿量统一描述：

运动方程变成正则形式：\dot{\mathbf{q}} = \partial H/\partial \mathbf{p}，\dot{\mathbf{p}} = -\partial H/\partial \mathbf{q}。十八个一阶方程被一个标量函数 H 收编。这个形式的好处在第 1 章已埋下伏笔，这里正式兑现：哈密顿框架下，每一条连续对称性严格对应一笔守恒账。哈密顿量不明显含时间，对应能量守恒；不明显含空间原点，对应总动量守恒；不明显含空间方向，对应总角动量守恒；不明显含伽利略相位，对应质心匀速。数下来恰好十笔——这就是"十笔经典账"的会计出处，不是运气，是对称性清单。

还有一个容易漏看的细节：势能项是负号。引力势能取负、动能取正，总能量为零的系统恰好"动能刚好够挣脱束缚"，这把尺子在第 5 章判三体结局时要反复用到。

动手自检：一段受力代码的完整实验

第 4 章的积分器要成千上万次调用受力计算，这里先把这枚心脏零件造出来并做体检。实验背景：验证我们写的引力代码忠实于牛顿第三定律与守恒律，因为积分器的能量误差一半来自受力函数的粗心。

操作如下：构造随机初始状态，计算九个受力向量，做三项体检——合力是否为零、每对内力是否等大反向、加速度是否只依赖相对位置。代码不长，值得逐行过一遍：

import math, random

G = 1.0
random.seed(7)

bodies = [  # 质量与位置随手取，量级故意拉开
    {"m": 1.0, "q": [0.0, 0.0]},
    {"m": 0.6, "q": [1.1, 0.4]},
    {"m": 0.2, "q": [-0.3, 0.9]},
]

def accelerations(bodies):
    """返回每个天体的加速度向量"""
    acc = [[0.0, 0.0] for _ in bodies]
    for i, bi in enumerate(bodies):
        for j, bj in enumerate(bodies):
            if i == j:
                continue
            dx = bj["q"][0] - bi["q"][0]
            dy = bj["q"][1] - bi["q"][1]
            r = math.sqrt(dx * dx + dy * dy) + 1e-12   # 软化项，防零距离爆炸
            f = G * bj["m"] / (r ** 3)
            acc[i][0] += f * dx
            acc[i][1] += f * dy
    return acc

def total_mass_times_acc(bodies):
    acc = accelerations(bodies)
    return [sum(bodies[k]["m"] * acc[k][d] for k in range(3)) for d in range(2)]

# 体检一：合力为零（牛顿第三定律的整体表现）
Fx, Fy = total_mass_times_acc(bodies)
print("合力 x 分量:", Fx, " y 分量:", Fy)     # 均为 0 量级的小数

# 体检二：逐对内力等大反向
dx = bodies[1]["q"][0] - bodies[0]["q"][0]
dy = bodies[1]["q"][1] - bodies[0]["q"][1]
r = math.hypot(dx, dy)
f01 = G * bodies[0]["m"] * bodies[1]["m"] / r ** 2
dx2 = bodies[2]["q"][0] - bodies[0]["q"][0]
dy2 = bodies[2]["q"][1] - bodies[0]["q"][1]
r2 = math.hypot(dx2, dy2)
f02 = G * bodies[0]["m"] * bodies[2]["m"] / r2 ** 2
print("1-2 对内力大小:", f01, " 1-3 对内力大小:", f02)
​

结果：合力两个分量都输出在 10^{-17} 量级——浮点意义下严格为零，说明内力成对抵消；逐对内力大小与手算一致。解读这个结果有两层：工程层，这枚零件可以直接搬进第 4 章的积分器，合力为零保证总动量在数值演化中不被代码本身偷偷破坏；结构层，体检一实际验证的是"质心不加速"这笔账，它不依赖初值——这正是守恒量的本性，也是它比"算一遍对一遍"值钱的原因。

变式值得做两个。第一个，把某个质量设为零再跑：合力依然为零，第三个天体变成纯粹的"测试粒子"，这正是限制性三体的代码形态——近似的门槛在代码里就是这么一行。第二个，做标度实验：把所有位置放大 \lambda 倍、时间放大 \lambda^{3/2} 倍，方程的解依然成立（对应开普勒第三定律的标度对称性）。跑一遍前后两组轨道再对比，你会发现轨道形状逐点相似、只是时间轴被拉伸——这个"放大缩小不改解"的性质，就是等质量周期解扫描（第 4 章）只需搜两个参数而不必搜全部初值的深层原因。

两个易错点

软化项不是免费的。 代码里给距离加了个 10^{-12} 的小垫片，防止两体重合时除零。它让近距离交会的力略偏小——对大多数演示无伤大雅，但若你要研究碰撞轨道，这个垫片就是系统性偏差，得换成碰撞正则化这类专门技术。任何数值技巧都要像这样标明它的"有效范围"。

九个二阶方程与十八个一阶方程不是两道题。 初学者常把降阶当成"简化"。实际上两者完全等价，换成十八个一阶方程的动机纯粹是工程性的：数值积分算法的标准接口吃的是一阶方程组。概念上一阶二阶随意，工程上必须统一。

本节要点回顾

• 方程的工整是假象：与二体唯一的形式差别是"每个引力源变成两个移动源"，三层化约杠杆同时失效；
• 哈密顿形式是会计准则：一条对称性对应一笔守恒账，十笔经典账全部由此清点出来，势能取负是判束缚与否的尺子；
• 受力代码三重体检：合力为零、内力成对、只依赖相对位置——通过的零件才能进第 4 章的积分器；
• 标度对称性：位置与时间按固定幂次缩放解不变，这是后文参数扫描能砍掉维度的理论依据；
• 软化项有代价：每个数值技巧都要标明有效范围，碰撞问题的正则化是专门学问。

方程上台了，账本也开了户。下一节开始兑现：十笔账各自能换来多少维度的折扣，降维之后的相空间还剩多大、长什么样。