二维凸包,这篇博客已经说得够好了,介绍了斜率逼近法、Jarvis算法,Graham算法,还有Andrew算法。我这篇博客只会非常详细的介绍Andrew算法。
数论小白都能看懂的平面凸包详解 - ShineEternal的笔记小屋 - 洛谷博客 (luogu.com.cn)
我相信凭借着我6个粉丝其中5个都是老熟人的传播量,应该不会因为乱贴别人链接导致啥问题()。也许会有朋友要问了,人家都写的这么好了,你博客的创新点在哪里呢?(谁问你了)我主要解决的是这三个问题:
- 关于三点共线等等的特殊情况,网上的有些代码会被hack掉,我在这里给出一个相对比较靠谱的代码。
- 每个人对Andrew算法的理解可能都有点点不一样,也许我的博客会更适合你。
- 我会尽量探讨,在不同的情况中,板子应该怎么做一些小改动。
那我们开始吧。
-
Andrew算法
- 算法思想
- 代码思路
- 代码实现
- 小小的改动
Andrew算法
算法思想
要对一个点集求凸包,就是用一个凸多边形把这些点集全部围起来,这些点最多只能在凸多边形上,不能超出这个凸多边形,并且这个凸多边形是所有可以满足上述条件的凸多边形中,面积、周长均是最小的那一个。并且,这样的一个凸多边形是唯一的。我们从感性上就可以感觉到,这个凸多边形的顶点一定是从点集中选取的。
假如黑色的点是原本的点集,点P是形成的凸包中的某一个顶点。P不属于初始,并且凸包把所有的点围了起来。那么,一定存在初始点集中两个点,把它们相连之后,围成的凸多边形仍然是凸包,并且面积比原来的凸包更小。有了这样一个认知,我们要做的事情就变成了:从点集中选出若干点构成凸包。也可以说,这些点肯定是最“外面”的凸包框住的,对吧?
首先,我们把所有点从左到右的排个序,找到点集的最左边的点A和最右边的点B,点A和点B肯定是最“外面”的点,也就是说,肯定是构成凸包的点。我们想象现在有一块柔软的布从上到下垂下来,盖住了这些点,形成了凸包的上边界:
这个上边界很明显,起始于A,终于B。我们先不想一次找完整个凸包,我们先找上边界。
把点按照这个方法排序,可以得到\(P=\{p_1,...,p_7\}\)
bool operator < (Point& B) {
if (!sgn(x - B.x)) return y < B.y;
return x < B.x;
}
初始点集\(U=\{A\}\),然后,A先试探性的和\(p_1\)相连,然后把\(p_1\)插入点集U。嗯,很显然\(U=\{A,p_1\}\)就是\(A,p_1\)的上边界。再看\(p_2\),如果把\(p_1p_2\)相连,那\(U=\{A,p_1,p_2\}\)就是这三个点的上边界。
现在压力来到\(p_3\),把\(p_2p_3\)相连后,情况开始不对了,为了维护美丽的上边界形状,我们删除\(\overrightarrow{p_1p_2},\overrightarrow{Ap_1}\),然后把\(p_3\)加进来,现在\(U=\{A, p_3\}\)。
重复上述过程,最终我们得到了\(U=\{A, p_3, p_7, B\}\)
纵观整个寻找上边界的过程,其实可以描述为:按照排序从小到大依次遍历节点,当遍历到第\(i\)个节点时,看节点\(i\)是否当前的上边界\(U\)的末尾两点构成的向量\(\overrightarrow{v}\)的右侧,如果是,则将点\(i\)加入\(U\),否则,从\(U\)中不断弹出最后一个节点,直到满足条件或者弹到没点可弹为止。
维护下边界,就是反着再来一遍:按照排序从大到小依次遍历节点,当遍历到第\(i\)个节点时,……(不要指望着我再告诉你怎么做了!自己推一下)就像这样:
嗯嗯,最后我们可以得到,\(U=\{A, p_3, p_7, B\}; D=\{ B, p_5, p_2, A\}\)
合起来就是凸包:
代码思路
问题已经解决了,现在我们只要把我们的思想用代码表示出来即可。
回过头来看凸包的形成过程,我们为什么能判断\(p_7\)是新的上边界,而把\(p_6\)弹出了呢?通过仔细的观察,我们发现是因为\(p_7\)这个点在向量\(\overrightarrow{p_3p_6}\)的左侧。如果在维护上边界的时候,点在向量左侧,那就说明\(p_7\)在“更上面”,不是吗?维护下边界也是一样:如果当前的点\(p_1\)在下边界末尾两点构成向量\(\overrightarrow{p_5p_2}\)的右侧,就把当前点加入下边界;如果当前的点在向量左侧,就把下边界末尾的点弹出,直到没点可弹或满足点在向量右侧为止。
判断点和向量(其实这里是有方向的直线)的位置关系,如果你已经看过我的[计算几何] 1 二维基础运算与概念 - 跳岩 - 博客园 (cnblogs.com)的话,应该可以已经非常轻松的写出来代码了罢()
// 判断直线a和点b的关系
// 1: a在直线左边; 0: a在直线上; -1: a在直线右边
int relation(Line a, Point b){
return sgn(cross(a.v, b - a.p));
}
代码实现
[P2742 USACO5.1] 圈奶牛Fencing the Cows /【模板】二维凸包 - 洛谷 | 计算机科学教育新生态 (luogu.com.cn)
现在我们就以这道题为例,看看怎么做一个二维凸包。
andrewScan的伪代码如下,真的代码其实就多了一点点边界的判断。注意看,我们给点排序的时间是\(O(n\log n)\),找上下边界的时间是\(O(n)\),故总的时间复杂度为\(O(n\log n)\)。
Polygon andrewScan(Polygon p){
sort(p.begin(), p.end()); // 给点排序
Polygon u, d; // u: 上边界; d: 下边界
u.push_back(p[0], p[1]); // 把最左边的两个点push到上边界
d.push_back(p[end], p[end - 1]); // 把最右边的两个点push到下边界
// 找上边界
for (int i = 2; i < p.size(); ++i){
while (u.size() >= 2 && p[i]不在u末尾向量的右边) 弹出u末尾元素;
u.push_back(p[i]);
}
// 找下边界
for (int i = end - 2; i >= 0; --i){
while (d.size() >= 2 && p[i]不在u末尾向量的右边) 弹出d末尾元素;
d.push_back(p[i]);
}
// 现在找到了上下边界, 我们把结果合并起来
// 本人喜欢逆时针存储
result = merge(d, u, anticlockwise);
return result;
}
AC的代码如下:
#include <bits/stdc++.h>
using namespace std;
#define ll long long
#define ZERO 1e-8
#define xx first
#define yy second
#define RIGHT -1
#define ON 0
#define LEFT 1
int sgn(double x){
if (fabs(x) < ZERO) return 0;
return x > 0 ? 1 : -1;
}
struct Point{
double x, y;
Point (double _x = 0, double _y = 0) : x(_x), y(_y) {}
Point operator + (Point& B) { return Point(x + B.x, y + B.y); }
Point operator - (Point& B) { return Point(x - B.x, y - B.y); }
bool operator < (Point& B) {
if (!sgn(x - B.x)) return y < B.y;
return x < B.x;
}
friend ostream& operator <<(ostream& out, Point& p) { out << "(" << p.x << ", " << p.y << ")"; return out; }
friend void operator >>(istream& in, Point& p) { in >> p.x >> p.y; }
};
typedef Point Vector;
typedef pair<Point, Point> Line;
typedef vector<Point> Polygon;
double cross(Vector v1, Vector v2){
return v1.x * v2.y - v1.y * v2.x;
}
int relation(Line l, Point p){
return sgn(cross(l.yy - l.xx, p - l.xx));
}
double get_dist(Point p1, Point p2){
return sqrt((p1.x - p2.x) * (p1.x - p2.x) + (p1.y - p2.y) * (p1.y - p2.y));
}
Polygon andrewScan(Polygon p){
Polygon u, d;
// 如果p的点少于3, 可以直接返回了
if (p.size() < 3) return p;
sort(p.begin(), p.end());
u.push_back(p[0]), u.push_back(p[1]);
d.push_back(p.back()), d.push_back(p[p.size() - 2]);
for (int i = 2; i < p.size(); ++i){
int k = u.size();
while (k >= 2 && relation(Line(u[k - 2], u[k - 1]), p[i]) != RIGHT){
k--; u.pop_back();
}
u.push_back(p[i]);
}
for (int i = p.size() - 3; i >= 0; --i){
int k = d.size();
while (k >= 2 && relation(Line(d[k - 2], d[k - 1]), p[i]) != RIGHT){
k--; d.pop_back();
}
d.push_back(p[i]);
}
// 这里就是做了逆时针合并, 记得u.st=A, u.ed=B, d.st=B, d.ed=A
// 所以合并的时候还要注意不要重复存储两遍A和B
reverse(d.begin(), d.end());
for (int i = u.size() - 2; i > 0; --i){
d.push_back(u[i]);
}
return d;
}
int n;
const int N = 1e5 + 5;
int main(void){
cin >> n;
Polygon pg;
Point p;
for (int i = 0; i < n; ++i){
cin >> p.x >> p.y;
pg.push_back(p);
}
Polygon res = andrewScan(pg);
double dist = 0;
for (int i = 0; i < res.size(); ++i){
dist += get_dist(res[i], res[(i + 1) % res.size()]);
}
cout << fixed << setprecision(2) << dist << endl;
return 0;
}
小小的改动
- 注意这段代码:
for (int i = 2; i < p.size(); ++i) {
int k = u.size();
while (k >= 2 && relation(Line(u[k - 2], u[k - 1]), p[i]) != RIGHT) {
k--; u.pop_back();
}
u.push_back(p[i]);
}
我们的代码中,认为如果当前点不在向量的右边,就要把上边界末尾的点弹出。也就是说,如果\(p_1,p_2,p_3\)三点共线,我们最终会把\(p_2\)弹出,就像下面这样:
如果我们遇到需要保留\(p_2\)的情况,那只要把!=RIGHT
改为== LEFT
就好啦。文章来源:https://www.toymoban.com/news/detail-618062.html
for (int i = 2; i < p.size(); ++i) {
int k = u.size();
while (k >= 2 && relation(Line(u[k - 2], u[k - 1]), p[i]) == LEFT) {
k--; u.pop_back();
}
u.push_back(p[i]);
}
- 而且,我们的代码也是经得起一些hack的:
if (p.size() < 3) : 返回的是p;
if (p中所有点都是三点共线): 返回p的左端点和右端点;
好!水完了!有什么表达不清的地方请联系我哦!文章来源地址https://www.toymoban.com/news/detail-618062.html
到了这里,关于[计算几何] 2 二维凸包/笨蛋(我)也能看懂的二维凸包算法的文章就介绍完了。如果您还想了解更多内容,请在右上角搜索TOY模板网以前的文章或继续浏览下面的相关文章,希望大家以后多多支持TOY模板网!