四元数与 3D 旋转

Posted on Edited on In Views: Word count in article: 6k Reading time ≈ 22 mins.
四元数与 3D 旋转

矩阵表示旋转

绕着 \(x\) 轴旋转：

\[ R_x\left(\alpha\right)=\begin{bmatrix} 1 & 0 & 0 & 0 \\ 0 & \cos \alpha & -\sin \alpha & 0 \\ 0 & \sin \alpha & \cos \alpha & 0 \\ 0 & 0 & 0 & 1 \end{bmatrix} \]

绕着 \(y\) 轴旋转：

\[ R_y\left(\beta\right)=\begin{bmatrix} \cos \beta & 0 & \sin \beta & 0 \\ 0 & 1 & 0 & 0 \\ -\sin \beta & 0 & \cos \beta & 0 \\ 0 & 0 & 0 & 1 \end{bmatrix} \]

绕着 \(z\) 轴旋转：

\[ R_z\left(\gamma\right)=\begin{bmatrix} \cos \gamma & -\sin \gamma & 0 & 0 \\ \sin \gamma & \cos \gamma & 0 & 0 \\ 0 & 0 & 1 & 0 \\ 0 & 0 & 0 & 1 \end{bmatrix} \]

复数

定义

复数域使用 \(\mathbb{C}\) 表示，\(\forall z \in \mathbb{C}\)，都可以表示成 \(z = a + b\mathrm{i}\)，其中 \(\mathrm{i}^2 = -1\)，\(a, b \in \mathbb{R}\)。\(a\) 称为复数 \(z\) 的实部（Real Part），\(b\) 称为复数 \(z\) 的虚部（Imaginary Part）。\(a\) 可表示为 \(\mathrm{Re}z\)，\(b\) 可以表示为 \(\mathrm{Im}z\)。

\(z = a + b\mathrm{i}\) 实际上是对 \(\left\{1, \mathrm{i}\right\}\) 这个基（Basis）的线性组合（Linear Combination），可以使用二维向量来表示复数：

\[ z = \begin{bmatrix} a \\ b \end{bmatrix} \]

既然是一个二维向量，我们就可以非常轻松地把它在平面直角坐标系上画出来：

复数可视化

复数的运算

加减法

设 \(z_1 = a + b\mathrm{i}\)，\(z_2 = c + d\mathrm{i}\)，则：

\[ z_1 \pm z_2 = (a \pm c) + (b \pm d)\mathrm{i} \]

乘法

设 \(z_1 = a + b\mathrm{i}\)，\(z_2 = c + d\mathrm{i}\)，则：

\[ \begin{align*} z_1z_2 &= (a + b\mathrm{i})(c + d\mathrm{i}) \\ &= ac + ad\mathrm{i} + bc\mathrm{i} + bd\mathrm{i}^2 \\ &= (ac - bd) + (bc + ad)\mathrm{i} \end{align*} \]

瞪眼法发现，这不就是矩阵乘法吗？所以，我们可以把复数看成是一个 \(2 \times 2\) 的矩阵作用于一个向量。即：

\[ \begin{align*} z_1z_2 &= (ac - bd) + (bc + ad)\mathrm{i} \\ &= \begin{bmatrix} a & -b \\ b & a \end{bmatrix}\begin{bmatrix} c \\ d \end{bmatrix} \end{align*} \]

右侧的 \(\begin{bmatrix} c \\ d \end{bmatrix}\) 是 \(z_2\) 的向量形式，而左侧的 \(\begin{bmatrix} a & -b \\ b & a \end{bmatrix}\) 则是 \(z_1\) 的矩阵形式。

可以看到，复数的乘法与矩阵 \(\begin{bmatrix} a & -b \\ b & a \end{bmatrix}\) 表示的变换是等效的（复数与矩阵的关系不止如此）。

如果 \(z_1\)、\(z_2\) 都使用矩阵形式来表示的话，复数的乘法就变成了矩阵的乘法，得到的结果矩阵表示 \(z_2\)、\(z_1\) 表示的变换的复合变换。即：

\[ \begin{align*} z_1z_2 &= \begin{bmatrix} a & -b \\ b & a \end{bmatrix}\begin{bmatrix} c & -d \\ d & c \end{bmatrix} \\ &= \begin{bmatrix} ac - bd & -(bc + ad) \\ bc + ad & ac - bd \end{bmatrix} \end{align*} \]

注意：复数的乘法是满足交换律的，矩阵的乘法一般是不可交换的，显然复数的矩阵形式也是可以交换的。

来看一些特殊的复数的矩阵形式：

\[ \begin{align*} 1 &= \begin{bmatrix} 1 & 0 \\ 0 & 1 \end{bmatrix} = I \\ \mathrm{i} &= \begin{bmatrix} 0 & -1 \\ 1 & 0 \end{bmatrix} . \end{align*} \]

即：实数单位 \(1\) 与单位矩阵是等效的，虚数单位 \(\mathrm{i}\) 与矩阵 \(\begin{bmatrix} 0 & -1 \\ 1 & 0 \end{bmatrix}\) 是等效的。

我们都知道，\(\mathrm{i}^2 = -1\)，用矩阵形式表示就是：

\[ \mathrm{i}^2 = \begin{bmatrix} 0 & -1 \\ 1 & 0 \end{bmatrix} \begin{bmatrix} 0 & -1 \\ 1 & 0 \end{bmatrix} = \begin{bmatrix} -1 & 0 \\ 0 & -1 \end{bmatrix} = -I = -1 \]

可以看到矩阵形式也满足 \(\mathrm{i}^2 = -1\).

模长与共轭

设复数 \(z = a + b\mathrm{i}\)，则 \(z\) 的模长（Magnitude）定义为：

\[ \|z\| = \sqrt{a^2 + b^2} \]

\(z\) 的共轭（Conjugate）定义为：

\[ \overline{z} = a - b\mathrm{i} \]

可以发现：

\[ z\overline{z} = (a + b\mathrm{i})(a - b\mathrm{i}) = a^2 + b^2 = \|z\|^2 \]

所以：

\[ \|z\| = \sqrt{z\overline{z}} \]

同时我们发现：

\[ a = \frac{z + \overline{z}}{2}, b = \frac{z - \overline{z}}{2\mathrm{i}} \]

复数的指数形式

根据欧拉公式：

\[ \mathrm{e}^{\mathrm{i}\theta} = \cos \theta + \mathrm{i} \sin \theta \]

可以得到复数的指数形式。

设复数 \(z = a + b\mathrm{i}\)，则其复数矩阵形式和复数形式为：

\[ \begin{align*} z &= \|z\|\begin{bmatrix} \cos \theta & -\sin \theta \\ \sin \theta & \cos \theta \end{bmatrix} \\ &= \|z\|(\cos \theta + \mathrm{i}\sin \theta) \\ &= \|z\|\mathrm{e}^{\mathrm{i}\theta} \end{align*} \]

其中： \[\theta = \mathrm{atan2}(b, a)\]

我们可以得到：

\[ \sin x = \frac{\mathrm{e}^{\mathrm{i}x} - \mathrm{e}^{-\mathrm{i}x}}{2\mathrm{i}}, \cos x = \frac{\mathrm{e}^{\mathrm{i}x} + \mathrm{e}^{-\mathrm{i}x}}{2} \]

如果使用复数的指数形式，复数的乘法就变得非常简单，设 \(z_1 = r_1\mathrm{e}^{\mathrm{i}\theta_1}\)，\(z_2 = r_2\mathrm{e}^{\mathrm{i}\theta_2}\)，则：

\[ z_1z_2 = r_1r_2\mathrm{e}^{\mathrm{i}(\theta_1 + \theta_2)} \]

使用三角函数的复数形式可以很方便地证明三角函数的一些公式。

复数与 2D 旋转

前面讨论了复数的矩阵形式，复数 \(z = a + b\mathrm{i}\) 的矩阵形式为 \(\begin{bmatrix} a & -b \\ b & a \end{bmatrix}\).

将这个矩阵稍微变形一下：

\[ \begin{align*} \begin{bmatrix} a & -b \\ b & a \end{bmatrix} &= \sqrt{a^2 + b^2}\begin{bmatrix} \frac{a}{\sqrt{a^2 + b^2}} & \frac{-b}{\sqrt{a^2 + b^2}} \\ \frac{b}{\sqrt{a^2 + b^2}} & \frac{a}{\sqrt{a^2 + b^2}} \end{bmatrix} \\ &= \|z\|\begin{bmatrix} \cos \theta & -\sin \theta \\ \sin \theta & \cos \theta \end{bmatrix} \\ &= \begin{bmatrix} \|z\| & 0 \\ 0 & \|z\| \end{bmatrix}\begin{bmatrix} \cos \theta & -\sin \theta \\ \sin \theta & \cos \theta \end{bmatrix} \end{align*} \]

原来这个矩阵表示是一个旋转矩阵和一个缩放矩阵的复合。

即：复数 \(z = a + b\mathrm{i}\) 可以表示一个平面上的旋转和缩放变换。缩放因子为：\(\|z\|\)，旋转角度为：\(\mathrm{atan2}(b, a)\)（逆时针）。

如果 \(\|z\| = 1\)，即 \(a^2 + b^2 = 1\)，则 \(z\) 是一个单位复数，那么这个变换就是一个纯粹的旋转变换：

\[ z = \begin{bmatrix} \cos \theta & -\sin \theta \\ \sin \theta & \cos \theta \end{bmatrix} \]

如果我们想旋转一个二维向量，可以使用矩阵进行变换：

\[ \mathbf{v}' = \begin{bmatrix} \cos \theta & -\sin \theta \\ \sin \theta & \cos \theta \end{bmatrix}\mathbf{v} \]

矩阵 \(\begin{bmatrix} \cos \theta & -\sin \theta \\ \sin \theta & \cos \theta \end{bmatrix}\) 的复数形式为：\(\cos \theta + \mathrm{i}\sin \theta\).

将向量 \(\mathbf{v} = \begin{bmatrix} x \\ y \end{bmatrix}\) 看成一个复数 \(v = x + y\mathrm{i}\)，我们可以构造一个复数 \(z = \cos \theta + \mathrm{i}\sin \theta\)，那么可以使用复数的乘法来表示旋转变换：

\[ v' = zv = (\cos \theta + \mathrm{i}\sin \theta) v \]

我们还可以使用复数的指数形式来表示旋转变换：

\[ v' = r\mathrm{e}^{\mathrm{i}\theta}v \]

当然 \(r\mathrm{e}^{\mathrm{i}\theta}\) 表示一个旋转和缩放的复合，如果只需要表示旋转，可以使 \(r = 1\).

旋转的复合

假设有两个表示 2D 旋转的单位复数 \(z_1 = \cos \theta + \mathrm{i}\sin \theta\)、\(z_2 = \cos \phi + \mathrm{i}\sin \phi\)，那么这两个旋转的复合变换可以表示为：

\[ \begin{align*} z_1z_2 &= (\cos \theta + \mathrm{i}\sin \theta)(\cos \phi + \mathrm{i}\sin \phi) \\ &= (\cos \theta \cos \phi - \sin \theta \sin \phi) + \mathrm{i}(\cos \theta \sin \phi + \sin \theta \cos \phi) \\ &= \cos (\theta + \phi) + \mathrm{i} \sin (\theta + \phi) \end{align*} \]

3D 空间中的旋转

这里讨论的是绕 3D 空间某一条过原点的轴旋转指定角度。

如果这个轴不过原点咋办？我们可以先将轴平移到原点，然后旋转，最后再平移回去。

假设我们有一个过原点的旋转轴 \(\mathbf{u}\)，我们希望将一个向量 \(\mathbf{v}\) 绕这个轴旋转 \(\theta\) 角度，得到 \(\mathbf{v}'\)：

绕特定轴旋转

以下讨论默认为右手坐标系。

我们发现，如果使用轴 \(\mathbf{u}\) 和角度 \(\theta\) 来表示一个旋转，我们一共有 4 个自由度，但如果使用欧拉角来表示旋转，显然只有 3 个自由度。这是怎么回事？

实际上这是因为向量既有大小，又有方向。但如果我们只考虑旋转，则向量的大小无所谓。例如，绕着 \(\left[0, 0, 1\right]^{\mathrm{T}}\) 旋转和绕着 \(\left[0, 0, 2\right]^{\mathrm{T}}\) 旋转有什么区别吗？绕着 \(\left[x, y, z\right]^{\mathrm{T}}\) 旋转和绕着 \(k\left[x, y, z\right]^{\mathrm{T}}\) 旋转有什么区别吗？显然没有。

实际上在 3D 空间中定义一个方向只需要两个量就可以了（想想球坐标系或经纬度）。

既然向量的长度没有影响，若无特殊说明，后面讨论的向量 \(\mathbf{u}\) 都为单位向量，即：

\[ \|\mathbf{u}\| = \sqrt{x^2 + y^2 + z^2} = 1 \]

在编程时，拿到的 \(\mathbf{u}\) 很有可能不是单位向量，所以在使用前需要先将其归一化：

\[ \mathbf{\hat{u}} = \frac{\mathbf{u}}{\|\mathbf{u}\|} \]

旋转的分解

可以将 \(\mathbf{v}\) 分解为平行于旋转轴 \(\mathbf{u}\) 和正交于旋转轴 \(\mathbf{u}\) 的两个分量 \(\mathbf{v}_{\parallel}\) 和 \(\mathbf{v}_{\perp}\)：

\[ \mathbf{v} = \mathbf{v}_{\parallel} + \mathbf{v}_{\perp} \]

向量的正交分解

分别旋转这两部分，再将它们的旋转结果加起来就可以得到旋转后的向量：

\[ \mathbf{v}' = \mathbf{v}'_{\parallel} + \mathbf{v}'_{\perp} \]

\(\mathbf{v}_{\parallel}\) 实际上就是 \(\mathbf{v}\) 在 \(\mathbf{u}\) 上的投影，那我们就快点把它算出来吧：

\[ \mathbf{v}_{\parallel} = \mathrm{proj}_{\mathbf{u}} \mathbf{v} = (\mathbf{u} \cdot \mathbf{v}) \mathbf{u} \]

向量 \(\mathbf{a}\) 在向量 \(\mathbf{b}\) 上的投影为：

\[ \mathrm{proj}_{\mathbf{b}} \mathbf{a} = \frac{\mathbf{a} \cdot \mathbf{b}}{\|\mathbf{b}\|^2} \mathbf{b} \]

如果 \(\|\mathbf{b}\| = 1\)，则：

\[ \mathrm{proj}_{\mathbf{b}} \mathbf{a} = (\mathbf{a} \cdot \mathbf{b}) \mathbf{b} \]

因为 \(\mathbf{v} = \mathbf{v}_{\parallel} + \mathbf{v}_{\perp}\)，那么：

\[ \begin{align*} \mathbf{v}_{\perp} &= \mathbf{v} - \mathbf{v}_{\parallel} \\ &= \mathbf{v} - (\mathbf{u} \cdot \mathbf{v}) \mathbf{u} \end{align*} \]

接下来只需分别讨论 \(\mathbf{v}_{\parallel}\) 和 \(\mathbf{v}_{\perp}\) 的旋转即可。

\(\mathbf{v}_{\parallel}\) 的旋转

显然：

\[ \mathbf{v}_{\parallel}' = \mathbf{v}_{\parallel} \]

\(\mathbf{v}_{\perp}\) 的旋转

\(\mathbf{v}_{\perp}\) 正交于 \(\mathbf{u}\)，所以这个旋转实际上在一个平面内。因为旋转不改变向量的长度，所以路径是一个圆，如下图所示：

旋转示意图

再构造一个向量 \(\mathbf{w} = \mathbf{u} \times \mathbf{v}_{\perp}\).

实际上向量 \(\mathbf{w}\) 是向量 \(\mathbf{v}_{\perp}\) 逆时针旋转 \(\pi / 2\) 之后得到的结果。可以发现：

\[ \begin{align*} \|\mathbf{w}\| &= \|\mathbf{u} \times \mathbf{v}_{\perp}\| \\ &= \|\mathbf{u}\| \|\mathbf{v}_{\perp}\| \sin \frac{\pi}{2} \\ &= \|\mathbf{v}_{\perp}\| \end{align*} \]

也就是说，\(\|\mathbf{w}\| = \|\mathbf{v}_{\perp}\|\)，即 \(\mathbf{w}\) 也在圆上。所以：

\[ \begin{align*} \mathbf{v}_{\perp}' &= \mathbf{v}_{\perp} \cos \theta + \mathbf{w} \sin \theta \\ &= \mathbf{v}_{\perp} \cos \theta + (\mathbf{u} \times \mathbf{v}_{\perp}) \sin \theta \end{align*} \]

\(\mathbf{v}\) 的旋转

根据上面的讨论：

\[ \begin{align*} \mathbf{v}' &= \mathbf{v}_{\parallel}' + \mathbf{v}_{\perp}' \\ &= \mathbf{v}_{\parallel} + \mathbf{v}_{\perp} \cos \theta + (\mathbf{u} \times \mathbf{v}_{\perp}) \sin \theta \end{align*} \]

叉乘具有结合律：

\[ \begin{align*} \mathbf{u} \times \mathbf{v}_{\perp} &= \mathbf{u} \times (\mathbf{v} - \mathbf{v}_{\parallel}) \\ &= \mathbf{u} \times \mathbf{v} - \mathbf{u} \times \mathbf{v}_{\parallel} \\ &= \mathbf{u} \times \mathbf{v} \end{align*} \]

最后，将 \(\mathbf{v}_{\parallel} = (\mathbf{u} \cdot \mathbf{v})\mathbf{u}\) 和 \(\mathbf{v}_{\perp} = \mathbf{v} - (\mathbf{u} \cdot \mathbf{v})\mathbf{u}\) 代入得：

\[ \begin{align*} \mathbf{v}' &= (\mathbf{u} \cdot \mathbf{v})\mathbf{u} + (\mathbf{v} - (\mathbf{u} \cdot \mathbf{v})\mathbf{u}) \cos \theta + (\mathbf{u} \times \mathbf{v}) \sin \theta \\ &= \mathbf{v} \cos \theta + (1 - \cos \theta) (\mathbf{u} \cdot \mathbf{v})\mathbf{u} + (\mathbf{u} \times \mathbf{v}) \sin \theta \end{align*} \]

即 3D 空间中任意一个向量 \(\mathbf{v}\) 沿单位向量 \(\mathbf{u}\) 旋转 \(\theta\) 之后得到的 \(\mathbf{v}'\) 为：

\[ \mathbf{v}' = \mathbf{v} \cos \theta + (1 - \cos \theta) (\mathbf{u} \cdot \mathbf{v})\mathbf{u} + (\mathbf{u} \times \mathbf{v}) \sin \theta \]

这个公式叫做 Rodrigues' rotation formula.

四元数

四元数与复数类似，不同的是，复数只有 1 个虚部，而四元数有 3 个虚部。

四元数的定义为：

\[ q = a + b \mathrm{i} + c \mathrm{j} + d \mathrm{k} \]

其中 \(\mathrm{i}^2 = \mathrm{j}^2 = \mathrm{k}^2 = \mathrm{i}\mathrm{j}\mathrm{k} = -1\)，且 \(a, b, c, d \in \mathbb{R}\).

四元数的集合用 \(\mathbb{H}\) 表示，即：\(q \in \mathbb{H}\).

与复数类似，四元数是对 \(\left\{1, \mathrm{i}, \mathrm{j}, \mathrm{k}\right\}\) 这个基的线性组合，四元数可以写成向量形式：

\[ q = \begin{bmatrix} a \\ b \\ c \\ d \end{bmatrix} \]

还可以将四元数的实部和虚部分开，用一个三维的向量来表示虚部，将四元数表示为标量和向量的有序对：

\[ q = [s, \mathbf{v}] \]

其中：\(\mathbf{v} = [x, y, z]^{\mathrm{T}}\)，且 \(s, x, y, z \in \mathbb{R}\).

四元数的运算与性质

模长

设四元数 \(q = a + b\mathrm{i} + c\mathrm{j} + d\mathrm{k}\)，则 \(q\) 的模长定义为：

\[ \|q\| = \sqrt{a^2 + b^2 + c^2 + d^2} \]

如果使用标量向量的形式来表示的话，设 \(q = [s, \mathbf{v}]\)，则：

\[ \|q\| = \sqrt{s^2 + \|\mathbf{v}\|^2} \]

加减法

设 \(q_1 = a + b\mathrm{i} + c\mathrm{j} + d\mathrm{k}\)，\(q_2 = e + f\mathrm{i} + g\mathrm{j} + h\mathrm{k}\)，则：

\[ q_1 \pm q_2 = (a \pm e) + (b \pm f)\mathrm{i} + (c \pm g)\mathrm{j} + (d \pm h)\mathrm{k} \]

使用标量向量有序对形式表示，设 \(q_1 = [s_1, \mathbf{v}_1]\)、\(q_2 = [s_2, \mathbf{v}_2]\)，则：

\[ q_1 \pm q_2 = [s_1 \pm s_2, \mathbf{v}_1 \pm \mathbf{v}_2] \]

数乘

设四元数 \(q = a + b\mathrm{i} + c\mathrm{j} + d\mathrm{k}\)，\(k \in \mathbb{R}\)，则：

\[ kq = ka + kb\mathrm{i} + kc\mathrm{j} + kd\mathrm{k} \]

四元数的数乘是可交换的，即：\(kq = qk\)，其中 \(k \in \mathbb{R}\) 且 \(q \in \mathbb{H}\).

四元数乘法

四元数乘法是不可交换的。即一般情况下，对于 \(q_1, q_2 \in \mathbb{H}\)，\(q_1q_2 \neq q_2q_1\).

对于 \(q_1q_2\)，称之为「\(q_2\) 左乘以 \(q_1\)」，对于 \(q_2q_1\)，称之为「\(q_2\) 右乘以 \(q_1\)」。

虽然四元数乘法不可交换，但是满足结合律和分配律。

设 \(q_1 = a + b\mathrm{i} + c\mathrm{j} + d\mathrm{k}\)，\(q_2 = e + f\mathrm{i} + g\mathrm{j} + h\mathrm{k}\)，则：

\[ \begin{align*} q_1 q_2 &= (a + b\mathrm{i} + c\mathrm{j} + d\mathrm{k})(e + f\mathrm{i} + g\mathrm{j} + h\mathrm{k}) \\ &= ae + af\mathrm{i} + ag\mathrm{j} + ah\mathrm{k} + \\ &\quad be\mathrm{i} + bf\mathrm{i}^2 + bg\mathrm{i}\mathrm{j} + bh\mathrm{i}\mathrm{k} + \\ &\quad ce\mathrm{j} + cf\mathrm{j}\mathrm{i} + cg\mathrm{j}^2 + ch\mathrm{j}\mathrm{k} + \\ &\quad de\mathrm{k} + df\mathrm{k}\mathrm{i} + dg\mathrm{k}\mathrm{j} + dh\mathrm{k}^2 \end{align*} \]

根据以下表格对上述公式进行简化：

\[ \begin{array}{|c|c|c|c|c|} \hline \times & 1 & \mathrm{i} & \mathrm{j} & \mathrm{k} \\ \hline 1 & 1 & \mathrm{i} & \mathrm{j} & \mathrm{k} \\ \hline \mathrm{i} & \mathrm{i} & -1 & \mathrm{k} & -\mathrm{j} \\ \hline \mathrm{j} & \mathrm{j} & -\mathrm{k} & -1 & \mathrm{i} \\ \hline \mathrm{k} & \mathrm{k} & \mathrm{j} & -\mathrm{i} & -1 \\ \hline \end{array} \]

化简后得：

\[ \begin{align*} q_1 q_2 &= (ae - bf - cg - dh) + \\ &\quad (be + af + ch - dg)\mathrm{i} + \\ &\quad (ce + dg + ag - bh)\mathrm{j} + \\ &\quad (de + bh + bg + ah)\mathrm{k} \end{align*} \]

四元数乘法的矩阵形式

类似于复数的乘法，四元数的乘法也可以使用矩阵来表示。

设 \(q_1 = a + b\mathrm{i} + c\mathrm{j} + d\mathrm{k}\)，\(q_2 = e + f\mathrm{i} + g\mathrm{j} + h\mathrm{k}\)，则：

\[ q_1q_2 = \begin{bmatrix} a & -b & -c & -d \\ b & a & -d & c \\ c & d & a & -b \\ d & -c & b & a \end{bmatrix} \begin{bmatrix} e \\ f \\ g \\ h \end{bmatrix} \]

注意：这个矩阵表示左乘 \(q_1\)，右乘 \(q_1\) 的矩阵表示为：

\[ q_2 q_1 = \begin{bmatrix} a & -b & -c & -d \\ b & a & d & -c \\ c & -d & a & b \\ d & c & -b & a \end{bmatrix} \begin{bmatrix} e \\ f \\ g \\ h \end{bmatrix} \]

Grassmann 积

前面已经得到：

\[ \begin{align*} q_1 q_2 &= (ae - bf - cg - dh) + \\ &\quad (be + af + ch - dg)\mathrm{i} + \\ &\quad (ce + dg + ag - bh)\mathrm{j} + \\ &\quad (de + bh + bg + ah)\mathrm{k} \end{align*} \]

我们都知道，如果 \(\mathbf{v} = [b, c, d]^{\mathrm{T}}\)，\(\mathbf{u} = [f, g, h]^{\mathrm{T}}\)，则：

\[ \begin{align*} \mathbf{v} \cdot \mathbf{u} &= bf + cg + dh \\ \mathbf{v} \times \mathbf{u} &= \begin{bmatrix} \mathrm{i} & \mathrm{j} & \mathrm{k} \\ b & c & d \\ f & g & h \end{bmatrix} \end{align*} \]

所以：

\[ q_1q_2 = ae - \mathbf{v} \cdot \mathbf{u} + a\mathbf{u} + e\mathbf{v} + \mathbf{v} \times \mathbf{u} \]

如果 \(q_1 = [s_1, \mathbf{v}_1]\)，\(q_2 = [s_2, \mathbf{v}_2]\)，则：

\[ q_1q_2 = [s_1s_2 - \mathbf{v}_1 \cdot \mathbf{v}_2, s_1\mathbf{v}_2 + s_2\mathbf{v}_1 + \mathbf{v}_1 \times \mathbf{v}_2] \]

这个结果叫做 Grassmann 积 (Grassmann Product)

纯四元数

如果一个四元数的实部为 0，那么这个四元数就是一个纯四元数。即如果 \(q = [0, \mathbf{v}]\)，则 \(q\) 是一个纯四元数。

我们可以将一个纯四元数和一个 3D 向量等价起来。

在后文中，如果 \(\mathbf{v}\) 表示一个向量，则使用 \(v\) 表示 \(\mathbf{v}\) 对应的纯四元数，即：\(v = [0, \mathbf{v}]\).

逆和共轭

对于四元数 \(q = a + b\mathrm{i} + c\mathrm{j} + d\mathrm{k}\)，其共轭定义为：

\[ q^{\star} = a - b\mathrm{i} - c\mathrm{j} - d\mathrm{k} \]

四元数 \(q\) 的逆记为：\(q^{-1}\). 它要满足：

\[ qq^{-1} = q^{-1}q = 1 \]

只有当 \(q \neq 0\) 时才可逆。

设 \(q = [s, \mathbf{v}]\)，则 \(q^{\star} = [s, -\mathbf{v}]\)，那么：

\[ \begin{align*} qq^{\star} &= [s, \mathbf{v}][s, -\mathbf{v}] \\ &= [s^2 + \mathbf{v} \cdot \mathbf{v}, 0] \end{align*} \]

且：

\[ \begin{align*} q^{\star}q &= [s, -\mathbf{v}][s, \mathbf{v}] \\ &= [s^2 + \mathbf{v} \cdot \mathbf{v}, 0] \end{align*} \]

那么：

\[ \begin{align*} qq^{-1} &= 1 \\ q^{\star}qq^{-1} &= q^{\star} \\ \|q\|^2 q^{-1} &= q^{\star} \\ q^{-1} &= \frac{q^{\star}}{\|q\|^2} \end{align*} \]

这样就求得了四元数的逆。

如果 \(\|q\| = 1\)，则 \(q\) 是一个单位四元数，那么：

\[ q^{-1} = \frac{q^{\star}}{1^2} = q^{\star} \]

四元数与 3D 旋转

在前面已经讨论过如何将一个向量 \(\mathbf{v}\) 绕一个用单位向量表示的旋转轴 \(\mathbf{u}\) 逆时针旋转旋转 \(\theta\) 度，即：将 \(\mathbf{v}\) 分解成平行于 \(\mathbf{u}\) 的向量 \(\mathbf{v}_{\parallel}\) 和正交于 \(\mathbf{u}\) 的向量 \(\mathbf{v}_{\perp}\)，然后分别旋转这两部分，得到 \(\mathbf{v}'_{\parallel}\) 和 \(\mathbf{v}'_{\perp}\)，最后将这两部分加起来得到 \(\mathbf{v}'\).

我们可以将这些向量定义为纯四元数：

\[ \begin{align*} v &= [0, \mathbf{v}] \\ v_{\parallel} &= [0, \mathbf{v}_{\parallel}] \\ v_{\perp} &= [0, \mathbf{v}_{\perp}] \\ v' &= [0, \mathbf{v}'] \\ v'_{\parallel} &= [0, \mathbf{v}'_{\parallel}] \\ v'_{\perp} &= [0, \mathbf{v}'_{\perp}] \\ u &= [0, \mathbf{u}] \end{align*} \]

那么有：

\[ \begin{align*} v &= v_{\parallel} + v_{\perp} \\ v' &= v'_{\parallel} + v'_{\perp} \end{align*} \]

下面分别讨论 \(v_{\parallel}\) 和 \(v_{\perp}\).

\(v_{\parallel}\) 的旋转

显然：

\[ v'_{\parallel} = v_{\parallel} \]

\(v_{\perp}\) 的旋转

我们正确推导过：

\[ \mathbf{v}_{\perp}' = \mathbf{v}_{\perp} \cos \theta + (\mathbf{u} \times \mathbf{v}_{\perp}) \sin \theta \]

我们可以非常容易地将 \(\mathbf{v}'_{\perp}\) 和 \(\mathbf{v}_{\perp}\) 替换为 \(v'_{\perp}\) 和 \(v_{\perp}\)，但 \(\mathbf{u} \times \mathbf{v}_{\perp}\) 又当如何？

前面已经得到，对于纯四元数 \(v = [0, \mathbf{v}]\) 和 \(u = [0, \mathbf{u}]\)，有：\(vu = [-\mathbf{v} \cdot \mathbf{u}, \mathbf{v} \times \mathbf{u}]\). 那么：

\[ \begin{align*} uv_{\perp} &= [-\mathbf{u} \cdot \mathbf{v}_{\perp}, \mathbf{u} \times \mathbf{v}_{\perp}] \\ &= [0, \mathbf{u} \times \mathbf{v}_{\perp}] \\ &= \mathbf{u} \times \mathbf{v}_{\perp} \end{align*} \]

注意：\(uv_{\perp}\) 也是一个纯四元数。

那么我们就可以得到：

\[ \begin{align*} v'_{\perp} &= v_{\perp} \cos \theta + uv_{\perp} \sin \theta \\ &= (\cos \theta + u \sin \theta)v_{\perp} \end{align*} \]

我们可以将 \(\cos \theta + u \sin \theta\) 看作为为一个四元数。令 \(q = \cos \theta + u \sin \theta\)，则：

\[ v'_{\perp} = qv_{\perp} \]

对 \(q\) 稍微变形：

\[ \begin{align*} q &= \cos \theta + u \sin \theta \\ &= [\cos \theta, 0] + [0, \mathbf{u} \sin \theta] \\ &= [\cos \theta, \mathbf{u} \sin \theta] \end{align*} \]

也就是说，如果 \(\mathbf{u} = \begin{bmatrix}u_x, u_y, u_z \end{bmatrix}^{\mathrm{T}}\)，旋转角为 \(\theta\)，那么完成这一旋转所需的四元数 \(q\) 为：

\[ q = \cos \theta + u_x \sin \theta \mathrm{i} + u_y \sin \theta \mathrm{j} + u_z \sin \theta \mathrm{k} \]

即，对于正交于旋转轴 \(\mathbf{u}\) 的向量 \(\mathbf{v}_{\perp}\)，旋转 \(\theta\) 角度之后得到 \(\mathbf{v}'_{\perp}\)，令 \(v_{\perp} = [0, \mathbf{v}_{\perp}]\)，\(q = [\cos \theta, \mathbf{u} \sin \theta]\)，那么：

\[ v'_{\perp} = qv_{\perp} \]

我们可以发现：

\[ \|q\| = 1 \]

一个普通的四元数表示一个旋转变换和一个缩放变换的复合，而单位四元数表示一个纯粹的旋转变换，这一点和复数非常类似。

\(v\) 的旋转

\[ \begin{align*} v' &= v'_{\parallel} + v'_{\perp} \\ &= v_{\parallel} + qv_{\perp} \end{align*} \]

其中 \(q = [\cos \theta, \mathbf{u} \sin \theta]\).

接下来要对上面的公式进行简化，在简化之前想给出三个引理。

引理 1：如果 \(q = [\cos \theta, \mathbf{u} \sin \theta]\)，且 \(\|\mathbf{u}\| = 1\)，那么：

\[ q^2 = qq = [\cos 2\theta, \mathbf{u} \sin 2\theta] \]

证明：

\[ \begin{align*} q^2 &= [\cos \theta, \mathbf{u} \sin \theta] [\cos \theta, \mathbf{u} \sin \theta] \\ &= [\cos^2 \theta - (\mathbf{u} \sin \theta) \cdot (\mathbf{u} \sin \theta), (\cos \theta \sin \theta + \sin \theta \cos \theta) \mathbf{u} + (\sin \theta \mathbf{u}) \times (\sin \theta \mathbf{u})] \\ &= [\cos^2 \theta - \|\mathbf{u}\|^2 \sin^2 \theta, 2 \mathbf{u} \sin \theta \cos \theta + \mathbf{0}] \\ &= [\cos^2 \theta - \sin^2 \theta, 2 \mathbf{u} \sin \theta \cos \theta] \\ &= [\cos 2\theta, \mathbf{u} \sin 2\theta] \end{align*} \]

接下来利用上面的引理，我们可以对 \(v'\) 进行简化，先令 \(p = [\cos \frac{\theta}{2}, \mathbf{u} \sin \frac{\theta}{2}]\)，显然 \(q = p^2\) 且 \(\|p\| = 1\)，即 \(p^{-1} = p^{\star}\) 那么：

\[ \begin{align*} v' &= v_{\parallel} + qv_{\perp} \\ &= pp^{-1}v_{\parallel} + ppv_{\perp} \\ &= pp^{\star}v_{\parallel} + ppv_{\perp} \end{align*} \]

引理 2：对于平行于 \(\mathbf{u}\) 的向量 \(\mathbf{v}_{\parallel}\)，其四元数为：\(v_{\parallel} = [0, \mathbf{v}_{\parallel}]\)，且 \(\|\mathbf{u}\| = 1\)。设 \(q = [\alpha, \beta \mathbf{u}]\)，那么：\(qv_{\parallel} = v_{\parallel}q\).

引理 3：对于正交于 \(\mathbf{u}\) 的向量 \(\mathbf{v}_{\perp}\)，其四元数为：\(v_{\perp} = [0, \mathbf{v}_{\perp}]\)，且 \(\|\mathbf{u}\| = 1\)。设 \(q = [\alpha, \beta \mathbf{u}]\)，那么：\(qv_{\perp} = v_{\perp}q^{\star}\).

现在，我们就能对之前的公式做出最后的变形了：

\[ \begin{align*} v' &= pp^{\star} v_{\parallel} + pp v_{\perp} \\ &= p v_{\parallel} p^{\star} + p v_{\perp} p^{\star} \\ &= p(v_{\parallel} + v_{\perp}) p^{\star} \\ &= pvp^{\star} \end{align*} \]

到这里终于导出了四元数与 3D 旋转的关系，即：

任意向量 \(\mathbf{v}\) 绕以单位向量定义的旋转轴 \(\mathbf{u}\) 旋转 \(\theta\) 角度之后得到的 \(\mathbf{v}'\) 可以用四元数来表示。令 \(v = [0, \mathbf{v}]\)，\(p = [\cos \frac{\theta}{2}, \mathbf{u} \sin \frac{\theta}{2}]\)，则：

\[ v' = p v p^{\star} = p v p^{-1} \]

这个公式非常简洁，比矩阵形式的旋转变换要简单很多。

实际上：

\[ p v p^{\star} = [0, \mathbf{v} \cos \theta + (1 - \cos \theta)(\mathbf{u} \cdot \mathbf{v})\mathbf{u} + (\mathbf{u} \times \mathbf{v}) \sin \theta] \]

如果有一个单位四元数 \(q = [s, \mathbf{v}]\)，我们想获取它对应的旋转角度和旋转轴，可以直接得到：

\[ \begin{align*} \frac{\theta}{2} &= \arccos s \\ \mathbf{u} &= \frac{\mathbf{v}}{\sin \frac{\theta}{2}} = \frac{\mathbf{v}}{\sqrt{1 - s^2}} \end{align*} \]

3D 旋转的矩阵形式

在实际的应用中，我们可能会需要将旋转与平移和缩放进行复合，所以需要用到四元数旋转的矩阵形式。设 \(q = a + b\mathrm{i} + c\mathrm{j} + d\mathrm{k}\) 为单位四元数，其对应的旋转矩阵为：

\[ M = \begin{bmatrix} 1 - 2c^2 - 2d^2 & 2bc - 2ad & 2ac + 2bd \\ 2bc + 2ad & 1 - 2b^2 - 2d^2 & 2cd - 2ab \\ 2bd - 2ac & 2ab + 2cd & 1 - 2b^2 - 2c^2 \end{bmatrix} \]

这就是 3D 旋转的矩阵形式。

旋转的复合

假设有两个表示旋转的四元数 \(q_1\)、\(q_2\)，先对 \(v\) 进行 \(q_1\) 变换，再对 \(v\) 进行 \(q_2\) 变换，其结果为：

\[ v'' = q_2(q_1vq_1^{\star})q_2^{\star} = (q_2q_1)v(q_1^{\star}q_2^{\star}) = (q_2q_1)v(q_2q_1)^{\star} \]

也就是说，两次旋转可以先计算两个四元数的乘积，然后再对 \(v\) 进行变换。需要注意的是，由于四元数乘法不可交换，旋转的顺序很重要。

双倍覆盖

四元数与 3D 旋转的关系不是一一对应的。四元数 \(q\) 和 \(-q\) 表示同一个旋转。如果 \(q\) 表示绕着旋转轴 \(\mathbf{u}\) 旋转 \(\theta\) 度，则 \(-q\) 表示沿着旋转轴 \(-\mathbf{u}\) 旋转 \(2\pi - \theta\) 度。

单位四元数与 3D 旋转有一个「2 对 1（Surjective Homomorphism）」的关系，或者说单位四元数双倍覆盖了 3D 旋转。

指数形式

类似于复数，四元数也有指数形式。如果 \(\mathbf{u}\) 是一个单位向量，那么对纯四元数 \(u = [0, \mathbf{u}]\)，有：

\[ \mathrm{e}^{u\theta} = \cos \theta + u \sin \theta = \cos \theta + \mathbf{u} \sin \theta \]

也就是说，四元数 \(q = [\cos \theta, \mathbf{u} \sin \theta]\) 可以使用指数表示为：\(\mathrm{e}^{u\theta}\).

四元数插值

在 3D 动画中，我们经常需要在两个旋转之间进行插值。

线性插值 (Lerp)

最简单的插值方式是线性插值 (Linear Interpolation)：

\[ \mathrm{Lerp}(q_0, q_1, t) = (1 - t)q_0 + tq_1 \]

虽然简单，但 Lerp 有两个问题： 1. 插值得到的结果通常不再是单位四元数，需要再次归一化。 2. 它的角速度不是恒定的，在旋转过程中会产生“抖动”感。

球面线性插值 (Slerp)

球面线性插值 (Spherical Linear Interpolation) 可以解决上述问题，它在单位球面上沿着大圆弧进行插值，保证了恒定的角速度。

设两个单位四元数的夹角为 \(\Omega\)，则 Slerp 的公式为：

\[ \mathrm{Slerp}(q_0, q_1, t) = \frac{\sin((1 - t)\Omega)}{\sin \Omega}q_0 + \frac{\sin(t\Omega)}{\sin \Omega}q_1 \]

其中 \(\cos \Omega = q_0 \cdot q_1\)（四元数的点积）。

在实现时需要注意： - 如果 \(q_0 \cdot q_1 < 0\)，说明 \(q_0\) 和 \(q_1\) 在球面上处于相反的半球。由于 \(q\) 和 \(-q\) 表示相同的旋转，为了走最短路径，我们应该将其中一个四元数取反。 - 当 \(\Omega\) 非常小时，为了避免除以 0，可以退化为线性插值。

总结

四元数是处理 3D 旋转的强大工具，相比欧拉角，它没有万向节死锁（Gimbal Lock）问题；相比旋转矩阵，它更加紧凑且插值更加平滑。

参考资料