2. 刚体运动学基础

在刚体转动数值模拟中，最先要解决的不是动力学方程，而是“如何把姿态和转动状态用数值表达清楚”。刚体姿态有三个转动自由度，但表达方式并不唯一：可以用旋转矩阵、欧拉角、四元数，也可以用轴角。不同表示之间的转换、坐标系约定、插值方式，直接决定了仿真是否会出现奇异、漂移或符号跳变。本节围绕这些基础工具展开，为后续角动量与转动惯量计算打底。

一、角速度矢量与角加速度

刚体的角速度是一个三维矢量 ω = [ωx, ωy, ωz]，方向沿瞬时转轴，大小等于转动快慢。它描述的是整个刚体的转动状态，而不是某个点的速度。刚体上位置为 r 的点，其线速度为 v = ω × r（相对转轴原点）。角加速度 α = dω/dt 同样是矢量，方向不一定与 ω 平行——只有当转轴方向不变时二者才共线。

这里必须强调参照系约定：ω 可以表示在惯性系（世界系）中，也可以表示在随体坐标系（本体系）中，两者通过旋转矩阵联系。数值模拟中常见做法是把 ω 存在本体系下，因为转动惯量张量在本体系中是常量，欧拉方程写起来最简洁。全文约定：世界系用上标 w，本体系用上标 b，旋转矩阵 R 表示“把本体系矢量变换到世界系”，即 v^w = R · v^b。

二、旋转矩阵与坐标系变换

旋转矩阵是 3×3 正交矩阵，满足 R^T R = I、det(R) = +1。它没有奇异点，复合旋转就是矩阵相乘，非常适合做坐标变换和连续积分。绕 x、y、z 轴的基本旋转分别为：

Rx(θ) = [1, 0, 0; 0, cosθ, -sinθ; 0, sinθ, cosθ]
Ry(θ) = [cosθ, 0, sinθ; 0, 1, 0; -sinθ, 0, cosθ]
Rz(θ) = [cosθ, -sinθ, 0; sinθ, cosθ, 0; 0, 0, 1]

坐标系变换的核心是“同一个矢量在不同基下坐标不同”。若已知本体系分量 v^b，则世界系分量为 v^w = R v^b；反过来 v^b = R^T v^w。角速度的变换稍特殊：ω^w = R ω^b 成立，但角速度本身不是普通矢量在旋转下的简单像，其时间导数关系为 dR/dt = [ω^w]× R = R [ω^b]×，其中 [·]× 是反对称矩阵。这个关系是旋转矩阵数值积分的出发点。

三、欧拉角：直观但有奇异

欧拉角用三个角度描述姿态，最常用的是 ZYX 顺序（先绕 z 偏航，再绕 y 俯仰，最后绕 x 滚转），对应航空中的 yaw-pitch-roll。其旋转矩阵为 R = Rz(ψ) Ry(θ) Rx(φ)。欧拉角直观、易读，适合人机交互和结果展示，但存在万向节锁：当俯仰角 θ = ±90° 时，偏航与滚转退化为同一自由度，雅可比矩阵奇异，数值求解会出现跳变。

因此工程上的常见策略是：内部积分用四元数或旋转矩阵，仅在输入输出和可视化时转成欧拉角。由旋转矩阵反解 ZYX 欧拉角的公式为 θ = -asin(R31)，ψ = atan2(R21, R11)，φ = atan2(R32, R33)，并需处理 |R31|≈1 的奇异分支。

四、数值表示与插值方法

四元数 q = [w, x, y, z]（单位四元数）是数值模拟中最常用的姿态表示：无奇异、只有 4 个分量、乘法复合高效。其导数关系为 dq/dt = 0.5 · q ⊗ [0, ω^b]，可用一阶或 RK4 积分。旋转矩阵与四元数可互相转换，四元数转矩阵时注意先归一化。

插值是姿态动画与轨迹生成的关键。对旋转矩阵直接做线性插值会破坏正交性，必须用球面线性插值（Slerp）：给定 q0、q1，先判断点积符号（若为负则取 -q1 保证走短弧），再按 q(t) = (sin((1-t)Ω) q0 + sin(tΩ) q1) / sinΩ 计算，其中 cosΩ = q0·q1。对欧拉角逐分量线性插值虽然简单，但在大角度或跨奇异区时会产生不自然的路径，只适合小角度近似。

五、一个可复现的数值示例

设某刚体初始姿态为单位四元数 q0 = [1,0,0,0]，本体系角速度恒为 ω^b = [0.5, 0.2, 0.1] rad/s，仿真步长 dt = 0.01 s，积分 1 秒。采用一阶四元数更新：

dq = 0.5 * q ⊗ [0, ωx, ωy, ωz] * dt
q_new = normalize(q + dq)

运行后应观察到四元数模长保持在 1 附近（归一化后严格为 1），姿态平滑变化；若把同一 ω^b 直接当作世界系角速度去更新欧拉角，则轨迹会明显偏离，因为忽略了本体系随刚体转动这一事实。判断实现是否正确，可以检查两点：一是 R^T R 与单位阵的最大偏差应小于 1e-6；二是用 R 反解出的欧拉角在非奇异区应与直接积分结果一致。这两个检查是刚体运动学数值实现中最基本、也最容易被忽略的验证手段。

最后提醒：所有约定（旋转顺序、矩阵作用方向、角速度所在坐标系、四元数分量顺序 wxyz 还是 xyzw）必须在代码和文档中统一写明。数值模拟中大量“看起来对但结果不对”的问题，根源往往不是算法，而是这些约定在模块之间不一致。