首先将所有点按横坐标、纵坐标排序，并删掉重复的点。然后分别求下凸壳和上凸壳。初始状态是点$P_0$和$P_1$。本质是单调队列，每次添加一个连边，保证删除队尾在凸包内的点。队尾的点在凸包内外，可以被转化为向量叉积的符号。遍历第一次后一定添加了$P_{n-1}$，就完成了下凸壳，然后倒序遍历，相当于把整个图像旋转$180 ^ \circ$再重新进行刚才的操作，就可以完成上凸壳。注意除非只有一个点，最后必然也会重复添加$P_0$，因此凸包的大小会多1。
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
| int ConvelHull() {
int k, m=0;
sort(p, p+n);
n = distance(p, unique(p, p+n));
For(i,n){
while(m>1 && Cross(ch[m-1]-ch[m-2], p[i]-ch[m-2]) <= 0) m--;
ch[m++] = p[i];
}
k = m;
rev(i, n-2, 0){
while(m>k && Cross(ch[m-1]-ch[m-2], p[i]-ch[m-2]) <= 0) m--;
ch[m++] = p[i];
}
if (n > 2) m--;
return m;
}
|
旋转卡壳
对于任何凸包的任何一条边$P_uP_{u+1}$，考虑凸包上它的对踵点$P_v$，即到它距离最远的那个点。注意这个距离可以通过三角形$P_uP_{u+1}P_v$的面积来比较。若我们从$P_{u+1}$开始，逆时针枚举凸包上所有点$P_v$，三角形$P_uP_{u+1}P_v$的面积呈单峰函数，即先增再减，从而当$P_uP_{u+1}P_{v+1} \le P_uP_{u+1}P_{v}$ 时，$v$取到三角形面积的最大值，即$P_v$是对踵点；注意当$P_uP_{u+1}P_{v+1} = P_uP_{u+1}P_{v}$时，$P_{v+1}$也是对踵点。这里注意，我们在建凸包的时候需要提前保证凸包的边上没有输入点，这样才能保证一条边的对踵点至多有两个。当我们逆时针枚举凸包上的所有边，对应的对踵点也必然以逆时针的方向移动，即满足决策单调性。再考虑到刚才提及的单峰函数性质，每次枚举边$P_uP_{u+1}$，都在前一条边$P_{u-1}P_{u}$所取的对踵点$v’$的基础上，尝试逐步自增$v’$到使得三角形面积最大的$v$即可。（一般资料称凸包上某点是另一点的对踵点，但这里延伸了这个概念。一点的对踵点一定是包含这一点的两条边的对踵点的其中一个。）
凸包的直径
凸包的直径是凸包上任意两个点的距离最大值，也可以化为：凸包上任意一条边的任意一个点与这条边的对踵点的距离的最大值。这样，按前面的算法，只需$O(n)$时间复杂度。
最小覆盖矩形
若要求最小覆盖矩形，就要在逆时针枚举凸包上的所有边的时候，在考虑其对踵点的同时，考虑其“最左点”和“最右点”。可以直观地理解这两个概念的定义：试把凸包旋转，使某一边成为水平底边，再以这水平底边为$x$轴建立平面直角坐标系，这时凸包横坐标最小的点即为“最左点”、横坐标最大的点即为“最右点”，而纵坐标最大的点就是对踵点。当然实际操作并不需要旋转坐标系之类的操作，本质就是比较点在这条边上投影的有向线段的大小，这用点积可以实现。这里的决策单调性和单峰函数性质同样满足，因此求法和对踵点是一样的。但特别注意最左点的初始值需要从对踵点开始，才能保证其单峰函数的性质。求完了三个极点之后，利用平面几何进行演算即可，注意利用点到直线的投影。 例题：最小矩形覆盖
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
| Point Project(Point a1, Point a2, Point b) {
Vector v = a2-a1;
return a1 + v*(Dot(v, b-a1) / Dot(v, v));
}
Point ret[4];
double solve() {
int m = ConvelHull();
if (m == 1) {
For(j, 4) ret[j] = ch[0];
return 0;
}
if (m == 2) {
ret[0] = ret[3] = ch[0];
ret[1] = ret[2] = ch[1];
return 0;
}
double ans = 1e100;
ch[m] = ch[0];
int v = 1, l, r = 1;
#define nxt(i) i=(i+1)%m
For(u, m) {
while(Area(ch[u], ch[u+1], ch[v]) <= Area(ch[u], ch[u+1], ch[v+1])) nxt(v);
if (!u) l = v;
while(Dot(ch[u+1]-ch[u], ch[l+1]-ch[u]) <= Dot(ch[u+1]-ch[u], ch[l]-ch[u])) nxt(l);
while(Dot(ch[u+1]-ch[u], ch[r+1]-ch[u]) >= Dot(ch[u+1]-ch[u], ch[r]-ch[u])) nxt(r);
Point Vp = Project(ch[u], ch[u+1], ch[v]);
Vector hv = ch[v] - Vp;
double h = Length(hv);
Point A[4];
A[0] = Project(ch[u], ch[u+1], ch[l]);
A[1] = Project(ch[u], ch[u+1], ch[r]);
double w = Length(A[0]-A[1]);
if (h * w < ans) {
ans = h * w;
A[2] = A[1] + hv;
A[3] = A[0] + hv;
For(j, 4) ret[j] = A[j];
}
}
return ans;
}
|