凸包算法

把若干颗钉子钉在木板上，然后用一根橡皮筋套住所有钉子并拉紧——橡皮筋绷出的形状就是这些点的凸包（Convex Hull）。

凸包是包含所有给定点的最小凸多边形。它是计算几何的基础问题，也是地图计算、碰撞检测、模式识别、机器学习（支持向量机的几何基础）的核心工具。更重要的是，凸包问题揭示了计算几何与排序之间的深刻联系——它们有相同的计算复杂度下界。

为什么下界是 $\Omega(n \log n)$？#

凸包不是"普通"的几何问题，它有一个漂亮的理论下界：任何基于比较的凸包算法都需要 $\Omega(n \log n)$ 时间。

下界证明（规约论证）：可以将排序问题规约为凸包问题。给定 $n$ 个实数 $x_1, \ldots, x_n$ 的排序问题，构造点集 $\{(x_i, x_i^2) \mid i = 1, \ldots, n\}$——所有点都在抛物线 $y = x^2$ 上。由于抛物线是凸曲线，所有点都在凸包上，且凸包的顶点顺序对应 $x_i$ 的排序。

因此：会排序 $\Rightarrow$ 会算凸包；会算凸包 $\Rightarrow$ 会排序（从凸包顶点顺序读出）。排序的 $\Omega(n \log n)$ 下界直接传递给凸包，达到 $O(n \log n)$ 的算法即为渐近最优。

Graham 扫描算法#

由 Ronald Graham 在 1972 年提出，时间 $O(n \log n)$（主要耗在排序上），是最常用的精确凸包算法之一。

步骤：

1. 找最低点（$y$ 最小，若相同取 $x$ 最小的）$p_0$，称为极点（pivot）
2. 将其余点按相对 $p_0$ 的极角从小到大排序
3. 初始化栈，push $p_0$ 和排序后第一个点
4. 对剩余点依次处理：
5. - 若当前点使栈顶两点的连线向右转（非左转），弹出栈顶
6. - 重复检查，直到左转或栈中只剩一点
7. - 将当前点压栈
8. 栈中剩余点即为凸包顶点（按逆时针顺序）

左转/右转的判断：用叉积（Cross Product）：

$$\text{叉积}(O, A, B) = (A_x - O_x)(B_y - O_y) - (A_y - O_y)(B_x - O_x)$$

• $> 0$：$O \to A \to B$ 是逆时针（左转）
• $< 0$：是顺时针（右转）
• $= 0$：共线

Graham 扫描的关键：每个点最多被压栈一次、弹栈一次，整体压弹操作 $O(n)$，瓶颈在排序 $O(n \log n)$。

手算走查：8 个点上的 Graham 扫描#

规则说完了不等于看见了。取这 8 个点，把整个扫描过程逐步写出来——每个叉积都给具体数值，拿笔就能复算：

$$(0,0),\ (4,0),\ (5,2),\ (4,4),\ (2,3),\ (0,4),\ (1,2),\ (2,1)$$

第一步：找 $p_0$。 $y$ 最小的是 $(0,0)$ 和 $(4,0)$，同为 $y=0$，取 $x$ 更小的，$p_0 = (0,0)$。

第二步：其余 7 点按相对 $p_0$ 的极角排序。 因为 $p_0$ 在原点，极角就是 $\arctan(y/x)$：

注意 $(2,1)$ 和 $(1,2)$ 这两个明显在内部的点，在极角排序里混在中间——排序不知道谁在凸包上，甄别工作全部交给后面的栈。

第三步：扫描。 栈初始化为 $p_0$ 和极角最小的点。每步用 叉积(次栈顶, 栈顶, 待处理点) 判方向：正为左转（保留），负为右转（弹出栈顶后重判）。

结果：凸包 $= (0,0) \to (4,0) \to (5,2) \to (4,4) \to (0,4)$，逆时针，$h = 5$。三个内部点 $(2,1)$、$(1,2)$、$(2,3)$ 都曾被压进栈，又都被后来的点挤出去。

对一下账：一共压栈 8 次、弹栈 3 次，$8 - 3 = 5$ 正是最终顶点数；叉积一共算了 9 次。这就是"每点最多进出一次，扫描阶段 $O(n)$"的具体形状——$n = 8$ 时 9 次叉积，瓶颈确实在第二步那次排序上。

这张表里最该盯住的是第 6 步：叉积等于 1，是全表离 0 最近的一次判断。$(4,4) \to (2,3) \to (1,2)$ 三点几乎在一条直线上（斜率分别是 $1/2$ 和 $1$，差得不多），符号只勉强为正。如果这些坐标不是小整数而是浮点数——比如来自 GPS、来自传感器、来自上一步几何运算的输出——这个"1"完全可能算成 $-0.0000001$。下一节就是这件事的后果。

失败现场：退化输入与那个几乎为零的叉积#

计算几何的算法设计在"实数算术精确"的机器模型上，而真实机器只有浮点。这道裂缝是几何代码最大的 bug 来源，而且它有两种独立的表现形式。

表现一：退化输入让"正确"的实现给出不同答案。

这四行不是想象出来的。Kettner、Mehlhorn、Pion、Schirra 与 Yap 在《几何计算中鲁棒性问题的课堂例子》（Computational Geometry 40(1), 2008；会议版 ESA 2004）里，专门挑了"平面凸包"和"三维 Delaunay 三角化"这两个最常被教的算法，构造出让它们输出根本不是凸的多边形、漏掉本应在凸包上的极点、乃至陷入死循环的浮点输入，并且给出了系统化构造这类输入的方法。他们选这两个算法的理由很直白：正因为它们简单到能写在黑板上，所以"简单算法不会有鲁棒性问题"这个错觉才特别值得打破。

表现二：叉积本身算不准。

方向谓词 $\text{orient2d}$ 本质是一个 $2\times2$ 行列式：

$$\text{orient2d}(A, B, C) = \begin{vmatrix} A_x - C_x & A_y - C_y \\ B_x - C_x & B_y - C_y \end{vmatrix}$$

用 double 直接算，误差来自两处叠加：坐标相减时如果两数接近，会发生灾难性抵消（catastrophic cancellation）——有效位数大量丢失；然后乘法把剩下的误差放大，最后两个乘积相减又抵消一次。当真值接近 0 时，算出来的符号可以完全是错的。而扫描算法用这个符号决定"压还是弹"，一次符号错误产生的不是"精度略差的凸包"，而是拓扑非法的输出：自交的多边形、顺逆时针混乱的顶点序、无法闭合的环。

Shewchuk 在 1997 年的《自适应精度浮点算术与快速鲁棒几何谓词》（Discrete & Computational Geometry 18(3), 305–363）给出了现在的标准解法，思路是分级：

1. 先用最快的普通浮点算一个近似值，同时算出这次计算的误差上界；
2. 如果 $|近似值| >$ 误差上界，符号必然可靠，直接返回——绝大多数输入走到这里就结束了；
3. 否则升到更高精度重算，必要时一路升到精确算术。

关键在于"运行时间取决于结果的不确定程度"：绝大多数点根本不接近共线，所以平均代价和普通浮点几乎一样，只有真正危险的那少数几次才付高精度的钱。

CGAL 的默认内核名字把这套哲学写在了脸上：Exact_predicates_inexact_constructions_kernel——谓词精确，构造不精确。为什么这个奇怪的组合是对的？因为算法的控制流由谓词决定，谓词错一次就产生非法拓扑，而且错误无法在后续步骤中被修正；构造出的坐标只影响数值精度，误差不会传染给分支判断。该记住的是这条区分：几何代码里，决定"走哪条路"的计算必须精确，决定"值是多少"的计算可以近似。

Jarvis 步进算法（Gift Wrapping）#

由 R.A. Jarvis 在 1973 年提出，另一种直觉更简单的算法：像用礼品纸包裹礼物一样，每次找"最右转"的下一个点。

1. 从最左点开始
2. 当前点 $p$，找使所有其他点都在 $p \to q$ 左边的点 $q$（即极角最小的点）
3. $q$ 成为下一个凸包顶点
4. 重复直到回到起点

时间复杂度：$O(nh)$，$h$ 是凸包顶点数。

当 $h$ 很小（如 $h = O(\log n)$）时，Jarvis 优于 Graham；当 $h = O(n)$ 时，退化到 $O(n^2)$。

Chan 算法：最优输出敏感算法#

Timothy Chan 在 1996 年提出了一个精妙的算法，时间 $O(n \log h)$（$h$ 为凸包顶点数），这是凸包问题理论上的最优输出敏感算法（Output-Sensitive Algorithm）。

Chan 算法结合了 Graham 和 Jarvis 的优点： 1. 猜测 $h$ 的值（从 $m = 2$ 开始，猜错就把猜测值平方：$m \leftarrow m^2$，最多 $O(\log \log h)$ 轮） 2. 将点集分成大小约 $m$ 的子组，对每组做 Graham 扫描（$O(n \log m)$） 3. 用 Jarvis 步进进行全局扫描，但每步仅在各子组的凸包上二分查找极角最小点（$O(\log m)$ 每步），最多只走 $m$ 步——走满 $m$ 步还没绕回起点，就判定本轮猜小了，整轮作废

Chan 算法展示了"猜测 + 验证"的自举（Bootstrap）技术。它的关键不在"怎么猜"，而在"猜错的代价能被下一轮吃掉"——下一节把这笔账算清楚。

为什么 Chan 算法"猜错也不亏"#

第 $t$ 轮（$t = 0, 1, 2, \ldots$）的猜测值是 $m_t = 2^{2^t}$，也就是 $2, 4, 16, 256, 65536, \ldots$。于是 $\log_2 m_t = 2^t$，本轮总代价：

• 分组 Graham：$\lceil n/m\rceil$ 组，每组 $O(m\log m)$ → 合计 $O(n\log m)$
• Jarvis 步进：最多 $m$ 步，每步在 $\lceil n/m\rceil$ 个子凸包上各二分一次，$O(\frac{n}{m}\log m)$ → 合计也是 $O(n\log m)$

两项同阶，所以第 $t$ 轮代价 $= O(n\log m_t) = O(n \cdot 2^t)$——代价随轮次指数增长。这是全部诀窍所在。

假设 $n = 10^6$、真实凸包顶点数 $h = 100$，把每一轮列出来（代价按 $n$ 的倍数记）：

前三轮全白干，一共花掉 $1n + 2n + 4n = 7n$；而成功的那一轮花 $8n$。所有失败轮加起来还不如最后一轮贵——这就是等比数列被最后一项支配的性质：$\sum_{t=0}^{T} 2^t = 2^{T+1} - 1 < 2\cdot 2^T$。

成功发生在 $m_T \ge h$，即 $2^{2^T}\ge h$，即 $2^T \ge \log_2 h$，所以 $T = \lceil \log_2\log_2 h\rceil$，总代价

$$\sum_{t=0}^{T} O(n\cdot 2^t) = O(n\cdot 2^{T+1}) = O(n\log h)$$

上例里 $\log_2 100 \approx 6.64$，总代价 $15n$ 约是 $n\log_2 h$ 的 2.3 倍——常数因子小到 2 出头，代价就是这么来的。

换个角度看这件事：如果每轮只把 $m$ 翻倍（$2, 4, 8, 16, \ldots$），要猜到 $h$ 需要 $\log_2 h$ 轮，第 $t$ 轮代价 $O(n\log m_t) = O(nt)$，总代价 $\sum_{t\le \log h} O(nt) = O(n\log^2 h)$——多出一个 $\log h$ 因子。平方而不是翻倍，正是为了让轮数从 $\log h$ 降到 $\log\log h$，同时让代价序列陡到"失败可以忽略"。 这个"用指数增长的试探把不知道的参数试出来"的手法在算法设计里反复出现：动态数组的容量倍增、TCP 的慢启动、无界二分查找的指数试探，账都是同一本。

三维及更高维凸包#

三维凸包：在三维空间中，凸包是凸多面体（Convex Polyhedron）。时间复杂度 $O(n \log n)$，输出大小 $O(n)$（三维凸包的面数线性于点数——Euler 公式保证）。

高维凸包：在 $d$ 维空间中，凸包的面数（faces）可达 $O(n^{\lfloor d/2 \rfloor})$（上界定理），计算复杂度急剧增加。维度 $d \geq 4$ 时，凸包计算实践上极为昂贵——这被称为"维度灾难"的几何版本。

应用领域#

机器人与运动规划：障碍物的凸包是碰撞检测中的"保守近似"——若机器人与凸包不相交，则一定不与原始形状相交（减少精确碰撞检测的频率）。GJK 算法（Gilbert-Johnson-Keerthi，1988）基于凸包做高效碰撞检测，是物理引擎（Unity、Bullet Physics）的标准工具。

支持向量机（SVM）：SVM 的二分类问题（找最大间隔超平面）在几何上等价于计算两类点的凸包的距离——两个凸包的最近点对应支持向量，间隔等于两凸包距离。

地图 GIS：计算城市或区域的边界轮廓，找一组 GPS 点的最小覆盖范围。

计算机图形学：确定物体的包围体（Bounding Volume），用于场景渲染中的视锥剔除（Frustum Culling），以及光线追踪的加速结构。

图像处理：手势识别中，手部轮廓的凸包和凸缺陷（Convexity Defects，凸包与实际轮廓之间的区域）用于识别手指数量和手势类型。

现场：碰撞检测其实用不到凸包本身#

物理引擎每帧要判断成千上万对物体是否相交。按教科书思路，应该先算出两个物体的凸包，再判断两个凸多面体是否相交。真实的做法不是这样。

GJK 算法（Gilbert、Johnson、Keerthi，IEEE Journal of Robotics and Automation 4(2), 199–203, 1988）从不显式构造凸包。它只要求物体提供一个支撑函数（support function）：

$$s_A(\mathbf{d}) = \arg\max_{\mathbf{p}\in A}\ \langle \mathbf{p}, \mathbf{d}\rangle$$

翻译成人话就是"沿着方向 $\mathbf{d}$ 看，你身上最远的那个点是谁"。有了它，GJK 在闵可夫斯基差 $A \ominus B = \{\mathbf{a} - \mathbf{b}\}$ 上迭代地挑几个点搭出一个单纯形（二维是三角形，三维是四面体），判断原点是否落在里面——原点在里面等价于两物体相交。整个过程通常只需三五次迭代，每次迭代只问一次支撑函数。

为什么这比"先算凸包"好得多？看支撑函数对各种形状的成本：

第一行和最后一行是关键：球没有有限的多面体凸包，但它的支撑函数最简单；旋转一个物体不需要重算任何东西，只需把查询方向反着转回去。支撑函数是凸集的一种表示，而且往往比顶点列表更廉价、更通用。

这一点在高维上更要紧。上一节说过 $d$ 维凸包的面数可达 $O(n^{\lfloor d/2\rfloor})$——显式构造在 $d \ge 4$ 时已不现实。而支撑函数的代价与维数几乎无关。凸优化里的"支撑函数 / 共轭函数"是同一个对象，$\text{SVM}$ 求两类点凸包最近距离的对偶形式，本质上也是在用支撑函数而不是顶点列表说话。

该记住的是：很多时候你需要的不是凸包这个物体，而是凸包能回答的那个问题——"沿这个方向最远的点是谁"。把需求写成查询而不是数据结构，常常能绕过一个昂贵的构造步骤。

动态凸包#

如果点集动态变化（支持插入/删除），维护凸包的代价更高。

动态凸包（Dynamic Convex Hull）支持： - 插入点：$O(\log^2 n)$ 均摊（使用某些平衡 BST 实现） - 删除点：$O(\log^2 n)$ 均摊 - 查询"点是否在凸包内"：$O(\log n)$

完全动态凸包的高效实现是计算几何中的研究热点，涉及复杂的分层数据结构。

跨域连接#

• 排序算法：把 $n$ 个实数映射到抛物线上，所有点都落在包上且顶点顺序就是排序结果——排序的对数线性下界因此原样传递过来。这是规约论证最干净的一次演示：证明新问题的下界不必从头分析，只需把已知难的问题线性时间地塞进去。反方向也成立，两个问题在比较模型下同难。
• 数值方法：方向判断是一个二阶行列式，坐标相减时的灾难性抵消能让接近共线的三点算出错误符号。符号错一次产生的不是精度略差的结果，而是拓扑非法的输出：自交多边形、顺逆序混乱、无法闭合。所以工业库分级处理——先浮点算并估误差界，不可靠时才升到精确算术。决定分支的计算必须精确，决定数值的计算可以近似。
• 计算材料设计：多组分体系里哪些相稳定，由生成焓对组成作图后的下凸包决定：落在包上的相稳定，落在上方的会分解成包顶点的组合。这不是类比而是同一个几何判据，高通量筛选新材料算的正是这条包络线，以及候选相到包络线的距离——距离越大，越不可能被合成出来。
• 投资组合理论与分散化：可行组合在风险收益平面上构成一个凸集，有效前沿就是它的上边界，被支配的组合都在内部。支撑函数的语言在这里同样成立：给定一个偏好方向，最优组合就是该方向上的支撑点。这解释了切点组合为何随无风险利率移动——移动的是方向，边界本身没变。
• 动量守恒与碰撞：物理引擎判断两物体是否相交时通常不显式构造包，只要求物体回答"沿这个方向最远的点是谁"。球没有有限的多面体包，它的支撑函数却只有一行；旋转物体不必重算任何东西，把查询方向反着转回去即可。而用包近似真实形状会系统性高估碰撞体积，所以凹形物体要先做凸分解再逐块判定。

参考文献#

• Graham, R. "An Efficient Algorithm for Determining the Convex Hull of a Finite Planar Set." Information Processing Letters 1(4), 1972.
• Chan, T. "Optimal Output-Sensitive Convex Hull Algorithms in Two and Three Dimensions." Discrete & Computational Geometry 16(4), 1996.
• Kettner, L., Mehlhorn, K., Pion, S., Schirra, S. & Yap, C. "Classroom Examples of Robustness Problems in Geometric Computations." Computational Geometry 40(1), 61–78, 2008.（会议版 ESA 2004；平面凸包与三维 Delaunay 的浮点失效实例及其构造方法）
• Shewchuk, J. R. "Adaptive Precision Floating-Point Arithmetic and Fast Robust Geometric Predicates." Discrete & Computational Geometry 18(3), 305–363, 1997.
• Gilbert, E. G., Johnson, D. W. & Keerthi, S. S. "A Fast Procedure for Computing the Distance Between Complex Objects in Three-Dimensional Space." IEEE Journal of Robotics and Automation 4(2), 199–203, 1988.（GJK 算法原论文）
• de Berg, M. et al. Computational Geometry: Algorithms and Applications. 3rd ed. Springer, 2008. （第1章）
• O'Rourke, J. Computational Geometry in C. 2nd ed. Cambridge University Press, 1998.

延伸阅读#

• CGAL 官方手册 — https://doc.cgal.org/ （Exact_predicates_inexact_constructions_kernel 等内核的选择依据）
• Ericson, C. Real-Time Collision Detection. Morgan Kaufmann, 2004.（GJK、包围体层次与凸形状表示的工程细节）
• Boyd, S. & Vandenberghe, L. Convex Optimization. Cambridge University Press, 2004.（支撑函数、共轭函数与凸集的对偶描述；免费在线版）
• Preparata, F. P. & Shamos, M. I. Computational Geometry: An Introduction. Springer, 1985.（凸包算法的早期系统整理）