凸包

(9 mins to read)

计算几何中的凸包

给定n个点求凸包利用andrew算法(比graham好写一点而且常数小)将所有点按照x为第一关键字,y为第二关键字排序然后从前往后用栈维护出下凸壳,再从后往前用栈维护出上凸壳即可考虑当前栈中的最后两个点p,q,现在加入一个点o,那么只要看q-p和o-p这两条直线的左右关系即可,用叉积判断多边形用vector存储,凸包分为严格凸包和非严格凸包(一般都是建成严格的),区别就是如果叉积等于0的时候不再pop,那么就是非严格的,可能会存在多点共线的情况另外要注意多边形点数少于3的情况

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
Polygon ConvexHull(Polygon& s)
{
int n = s.size(), sz = -1;
sort(s.begin(), s.end());
Polygon res(2*n);
for(int i=0; i<n; i++)
{
while(sz>0 && crossop(res[sz-1], res[sz], s[i])<=0) --sz;
res[++sz] = s[i];
}
int pre = sz;
for(int i=n-2; i>=0; i--)
{
while(sz>pre && crossop(res[sz-1], res[sz], s[i])<=0) --sz;
res[++sz] = s[i];
}
res.resize(sz);
return res;
}

多边形面积(相邻点叉积求和,逆时针为正,顺时针为负)

1
2
3
4
5
6
db getarea(Polygon& s)
{
db ans = 0; int n = s.size();
for(int i=0; i<n; i++) ans += (s[i]^s[(i+1)%n]);
return ans/2;
}

凸包周长(相邻点距离求和即可)

1
2
3
4
5
6
db getcir(Polygon& s)
{
db ans = 0; int n = s.size();
for(int i=0; i<n; i++) ans += s[i].dis(s[(i+1)%n]);
return ans;
}

凸包直径(旋转卡壳,对于每条边(i,i+1),找到距离这条边最远的点j,容易发现具有单调性,可以用双指针实现)

1
2
3
4
5
6
7
8
9
10
11
12
db ConvexDiameter(Polygon& s)
{
int n = s.size();
if(n==2) return s[0].dis(s[1]);
int j = 2; db ans = 0;
for(int i=0; i<n; i++)
{
while(cross(s[i], s[i+1], s[j])<cross(s[i], s[i+1], s[(j+1)%n])) j = (j+1)%n;
ans = max(ans, max(s[i].dis(s[j]), s[i+1].dis(s[j])));
}
return ans;
}

凸包的闵可夫斯基和还是凸包(定义就是一个凸包中的点按照原点和另一个凸包中的点形成的向量移动后形成的图形)闵可夫斯基和后的凸包的边是由原来的两个凸包上的边构成的,所以归并排序一下即可

1
2
3
4
5
6
7
8
9
10
11
12
13
14
Polygon Minkowski(Polygon& a, Polygon& b) //可能有三点共线,可以再求一次凸包
{
Polygon c;
int n = a.size(), m = b.size();
for(int i=0, j=0; i<n||j<m; )
{
c.push_back(a[i%n]+b[j%m]);
if(i>=n) ++j;
else if(j>=m) ++i;
else if(sgn((a[(i+1)%n]-a[i])^(b[(j+1)%m]-b[j]))>=0) ++i;
else ++j;
}
return c;
}

动态规划中的凸包

一般用于斜率优化,考虑一个dp转移式$f_i = \max \limits_{k<i} \{ f_k + a_ib_k\}$可以把后面的东西看成$b_ka_i + f_k$​,即$kx+b$的直线形式,然后我们就是要对每个$a_i$​找到之前所有直线中在该点处的最大值(别人的斜率优化可能是按照斜率的,不过我觉得这样子更直观一点,也不需要移项,如果了解李超树,很容易发现这完全可以用李超树来做)当然用李超树来做似乎有点大材小用,容易发现答案一定是在所有直线形成的上凸壳上,而且每条直线会有一个贡献范围

  • 如果插入的直线斜率单调($b_k$​),且询问点也单调($a_i$​)这时候我们只需要用双端队列来维护上凸壳即可,在队尾不断插入直线并且维护凸壳,在队首及时排掉已经结束贡献了的点,对于当前询问的最优点就是队首

  • 如果插入的直线斜率单调($b_k$​),但询问点不单调($a_i$​)斜率单调,所以我们仍然可以线性地维护上凸壳,但是我们不能pop_front来维护最优点,因为最优点没有单调性了,所以我们需要在凸壳上二分,找到最优决策点

  • 如果插入的直线斜率不单调($b_k$​),且询问点也不单调($a_i$​)一种方法就是李超树了,另一种方法是用平衡树来维护,以下给出一份set实现的代码(很好用,如果插入的是直线就用这个,如果是线段的话就只能李超树了)k和b分别是斜率和截距,p类似于最优决策点(不用管)div函数是一个下取整函数,负数的时候的取整和正数会有区别isect是求两条直线的交点,然后进行判断(类似于李超树中的判断)使用的时候只需要add(k,b)来添加一条直线,qry(x)来询问在x点处的最大值当然如果要求最小值,只需要add(-k,-b),然后-qry(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
bool Cmp;
struct Line
{
mutable ll k, b, p;
bool operator < (const Line& oth) const { return Cmp ? p < oth.p : k < oth.k; }
};
struct ConvecHull : multiset<Line>
{
static const ll inf = LLONG_MAX;
ll div(ll a, ll b)
{
if(!b) return a>0 ? inf : -inf;
return a / b - ((a ^ b) < 0 && a % b);
}
bool isect(iterator x, iterator y)
{
if (y == end()) { x->p = inf; return false; }
if (x->k == y->k) x->p = x->b > y->b ? inf : -inf;
else x->p = div((y->b - x->b), (x->k - y->k));
return x->p >= y->p;
}
void add(ll k, ll b)
{
auto z = insert({k, b, 0}), y = z++, x = y;
while (isect(y, z)) z = erase(z);
if (x != begin() && isect(--x, y)) isect(x, y = erase(y));
while ((y = x) != begin() && (--x)->p >= y->p) isect(x, erase(y));
}
ll qry(ll x)
{
Cmp = 1; auto l = *lower_bound({0, 0, x}); Cmp = 0;
return l.k * x + l.b;
}
}

如何删除一条直线

删除操作不好做(类似并查集),我们可以离线后把每条直线的存活时间塞到线段树上的log个节点上(即线段树分治),对每个节点维护一个凸壳,然后对于一个询问就是线段树上包含该点的节点的凸壳上的最大值的最大值给出所有直线再建凸壳,可以利用andrew来实现(排序后,从前往后用栈维护),询问最值一次二分即可(没必要三分)