计算几何学 | 凸包 | Convex Hull | C/C++实现
程序员文章站
2022-04-01 13:35:53
...
问题描述
求二维平面上的点集合P的凸包（Convex Hull）。凸包是指包含点集合P中所有点的最小凸多边形。请列举出该多边形的边及顶点上的所有点。
输入：
第1行输入点的数量n。接下来n行输入第 i 个点的坐标，坐标以2个整数、的形式给出。
输出：
第1行输出表示凸包的凸多边形的顶点数。接下来的几行输出凸多边形各个顶点的坐标(x, y)。输出时以凸多边形最下端左侧的顶点为起点，按逆时针方向依次输出坐标。
限制：
3 ≤ n ≤ 100000
-10000 ≤ ≤ 10000
所有点的坐标均不相同。
输入示例
7
2 1
0 0
1 2
2 2
4 2
1 3
3 3
输出示例
4
0 0
4 2
3 3
1 3
讲解
凸包可以理解为在木板上钉了许多钉子，然后用一根橡皮筋框住所有的钉子时所得到的多边形。在计算机几何学领域中，人们已经研究出了数个求凸包的算法。这里我们将要学习的是相对容易掌握的安德鲁算法（Andrew’s Algorithm）。
安德鲁算法：
1.将给定集合的点按x坐标升序排列。x相同的按y坐标升序排列
2.按下述流程创建凸包的上部：
将排序后的点按照x坐标从小到大顺序加入凸包U。如果新加入的点使得U不再是凸多边形，则逆序删除之前加入U的点，直至U重新成为凸多边形为止。
3.按照下述流程创建凸包的下部：
将排序后的点按照x坐标从大到小的顺序加入凸包L。如果新加入的点使得L不再是凸多边形，则逆序删除之前加入L的点，直至L重新成为凸多边形为止。
安德鲁算法借助栈操作构建凸包，这个操作的复杂度为。但实际上，整个算法中最消耗时间的是给点集合S排序，因此安德鲁算法整体的复杂度为。
安德鲁算法：
Polygon andrewScan(Polygon s) {
Polygon u, l;
if(s.size() < 3) return s;
sort(s.begin(), s.end());//以x、y为基准升序排序
//将x最小的2个点添加至u
u.push_back(s[0]);
u.push_back(s[1]);
//将x最大的2个点添加至l
l.push_back(s[s.size() - 1]);
l.push_back(s[s.size() - 2]);
//构建凸包上部
for(int i = 2; i < s.size(); i++) {
for(int n = u.size(); n >= 2 && ccw(u[n-2], u[n-1], s[i]) != CLOCKWISE; n--) {
u.pop_back();
}
u.push_back(s[i]);
}
//构建凸包下部
for( int i = s.size() - 3; i >= 0; i--) {
for(int n = l.size(); n >= 2 && ccw(l[n-2], l[n-1], s[i]) != CLOCKWISE; n--) {
l.pop_back();
}
l.push_back(s[i]);
}
//按顺时针方向生成凸包的点的序列
reverse(l.begin(), l.end());
for(int i = u.size() - 2; i >= 1; i--) l.push_back(u[i]);
return l;
}
AC代码如下
#include<iostream>
#include<cmath>
#include<vector>
#include<algorithm>
using namespace std;
#define EPS (1e-10)
#define equals(a, b) (fabs((a) - (b)) < EPS)
class Point {//Point类，点
public:
double x, y;
Point(double x = 0, double y = 0): x(x), y(y) {}
Point operator + (Point p) { return Point(x + p.x, y + p.y); }
Point operator - (Point p) { return Point(x - p.x, y - p.y); }
Point operator * (double a) { return Point(a * x, a * y); }
Point operator / (double a) { return Point(x / a, y / a); }
double abs() { return sqrt(norm()); }
double norm() { return x * x + y * y; }
bool operator < (const Point &p) const {
return x != p.x ? x < p.x : y < p.y;
}
bool operator == (const Point &p) const {
return fabs(x - p.x) < EPS && fabs(y - p.y) < EPS;
}
};
typedef Point Vector;//Vector类，向量
typedef vector<Point> Polygon;//Polygon 多边形
double dot(Vector a, Vector b) {//内积
return a.x * b.x + a.y * b.y;
}
double cross(Vector a, Vector b) {//外积
return a.x*b.y - a.y*b.x;
}
/*ccw模块*/
static const int COUNTER_CLOCKWISE = 1;
static const int CLOCKWISE = -1;
static const int ONLINE_BACK = 2;
static const int ONLINE_FRONT = -2;
static const int ON_SEGMENT = 0;
int ccw(Point p0, Point p1, Point p2) {
Vector a = p1 - p0;
Vector b = p2 - p0;
if( cross(a, b) > EPS ) return COUNTER_CLOCKWISE;
if( cross(a, b) < -EPS ) return CLOCKWISE;
if( dot(a, b) < -EPS ) return ONLINE_BACK;
if( a.norm() < b.norm() ) return ONLINE_FRONT;
return ON_SEGMENT;
}
Polygon andrewScan(Polygon s) {
Polygon u, l;
if(s.size() < 3) return s;
sort(s.begin(), s.end());//以x、y为基准升序排序
//将x最小的2个点添加至u
u.push_back(s[0]);
u.push_back(s[1]);
//将x最大的2个点添加至l
l.push_back(s[s.size() - 1]);
l.push_back(s[s.size() - 2]);
//构建凸包上部
for(int i = 2; i < s.size(); i++) {
for(int n = u.size(); n >= 2 && ccw(u[n-2], u[n-1], s[i]) != CLOCKWISE; n--) {
u.pop_back();
}
u.push_back(s[i]);
}
//构建凸包下部
for( int i = s.size() - 3; i >= 0; i--) {
for(int n = l.size(); n >= 2 && ccw(l[n-2], l[n-1], s[i]) != CLOCKWISE; n--) {
l.pop_back();
}
l.push_back(s[i]);
}
//按顺时针方向生成凸包的点的序列
reverse(l.begin(), l.end());
for(int i = u.size() - 2; i >= 1; i--) l.push_back(u[i]);
return l;
}
int main(){
int n;
cin>>n;
Polygon s;
Point p;
while(n--){
cin>>p.x>>p.y;
s.push_back(p);
}
cout<<s.size()<<endl;
s = andrewScan(s);
cout<<s.size()<<endl;
for(int i = 0; i < s.size(); i++){
cout<<s[i].x<<" "<<s[i].y<<endl;
}
}
注：以上本文未涉及代码的详细解释参见：计算几何学