词条 | 凸包 |
释义 | 定义1.对于一个集合D,D中任意有限个点的线性组合的全体称为D的凸包。 2.对于一个集合D,所有包含D的凸集之交称为D的凸包。 可以证明,上述两种定义是等价的 概念 1 点集Q的凸包(convex hull)是指一个最小凸多边形,满足Q中的点或者在多边形边上或者在其内。右图中由红色线段表示的多边形就是点集Q={p0,p1,...p12}的凸包。 2 一组平面上的点,求一个包含所有点的最小的凸多边形,这就是凸包问题了。这可以形象地想成这样:在地上放置一些不可移动的木桩,用一根绳子把他们尽量紧地圈起来,这就是凸包了。 平面凸包求法常见求法2.0 Graham's Scan法求解凸包问题概念 凸包(Convex Hull)是一个计算几何(图形学)中的概念。用不严谨的话来讲,给定二维平面上的点集,凸包就是将最外层的点连接起来构成的凸多边型,它能包含点集中所有点的。严谨的定义和相关概念参见维基百科:凸包。 这个算法是由数学大师葛立恒(Graham)发明的,他曾经是美国数学学会(AMS)主席、AT&T首席科学家以及国际杂技师协会(IJA)主席。(太汗了,这位大牛还会玩杂技~) 问题 给定平面上的二维点集,求解其凸包。 过程 1. 在所有点中选取y坐标最小的一点H,当作基点。如果存在多个点的y坐标都为最小值,则选取x坐标最小的一点。坐标相同的点应排除。然后按照其它各点p和基点构成的向量<H,p>与x轴的夹角进行排序,夹角由大至小进行顺时针扫描,反之则进行逆时针扫描。实现中无需求得夹角,只需根据向量的内积公式求出向量的模即可。以下图为例,基点为H,根据夹角由小至大排序后依次为H,K,C,D,L,F,G,E,I,B,A,J。下面进行逆时针扫描。 2. 线段<H, K>一定在凸包上,接着加入C。假设线段<K, C>也在凸包上,因为就H,K,C三点而言,它们的凸包就是由此三点所组成。但是接下来加入D时会发现,线段<K, D>才会在凸包上,所以将线段<K, C>排除,C点不可能是凸包。 3. 即当加入一点时,必须考虑到前面的线段是否会出现在凸包上。从基点开始,凸包上每条相临的线段的旋转方向应该一致,并与扫描的方向相反。如果发现新加的点使得新线段与上线段的旋转方向发生变化,则可判定上一点必然不在凸包上。实现时可用向量叉积进行判断,设新加入的点为pn + 1,上一点为pn,再上一点为pn - 1。顺时针扫描时,如果向量<pn - 1, pn>与<pn, pn + 1>的叉积为正(逆时针扫描判断是否为负),则将上一点删除。删除过程需要回溯,将之前所有叉积符号相反的点都删除,然后将新点加入凸包。 在上图中,加入K点时,由于线段<H,K>相对于<H,C>为顺时针旋转,所以C点不在凸包上,应该删除,保留K点。接着加入D点,由于线段<K, D>相对<H, K>为逆时针旋转,故D点保留。按照上述步骤进行扫描,直到点集中所有的点都遍例完成,即得到凸包。 复杂度 这个算法可以直接在原数据上进行运算,因此空间复杂度为O(1)。但如果将凸包的结果存储到另一数组中,则可能在代码级别进行优化。由于在扫描凸包前要进行排序,因此时间复杂度至少为快速排序的O(nlgn)。后面的扫描过程复杂度为O(n),因此整个算法的复杂度为O(nlgn)。 2.1 凸包最常用的凸包算法是Graham扫描法和Jarvis步进法。 对于一个有三个或以上点的点集Q,过程如下: 计算点集最右边的点为凸包的顶点的起点,如上图的P3点。 Do For i = 0 To 总顶点数 计算有向向量P3->Pi If 其余顶点全部在有向向量P3->Pi的左侧或右侧,则Pi点为凸包的下一顶点 Pi点加入凸包列表 GoTo 1 End If Next Exit Do 1: Loop 此过程执行后,点按极角自动顺时针或逆时针排序,只需要按任意两点的次序就可以了。而左侧或右侧的判断可以用前述的矢量点积性质实现。 特殊算法2.2 求凸包有很多方法,不过最适合OIer和ACMer的估计还是Graham's Scan这个方法了。它的大致方法是这样的: 首先,找到所有点中最左边的(y坐标最小的),如果y坐标相同,找x坐标最小的;以这个点为基准求所有点的极角(atan2(y-y0,x-x0)),并按照极角对这些点排序,前述基准点在最前面,设这些点为P[0]..P[n-1];建立一个栈,初始时P[0]、P[1]、P[2]进栈,对于P[3..n-1]的每个点,若栈顶的两个点与它不构成“向左转”的关系,则将栈顶的点出栈,直至没有点需要出栈以后将当前点进栈;所有点处理完之后栈中保存的点就是凸包了。 如何判断A、B、C构成的关系不是向左转呢?如果b-a与c-a的叉乘小于0就不是。a与b的叉乘就是a.x*b.y-a.y*b.x。 上面的这个Graham的实现比我原来按照USACO里的课文写得简单多了,主要是它通过简单的预处理保证了P[0]、P[1]以及P[n-1]肯定是凸包里的点,这样就可以避免在凸包“绕回来”的时候繁杂的处理。 中心法先构造一个中心点,然后将它与各点连接起来,按斜率递增的方法,求出凸包上部;再按斜率递减的方法,求出凸包下部。 水平法从最左边的点开始,按斜率递增的方法,求出凸包上部;再按斜率递减的方法,求出凸包下部。水平法较中心法减少了斜率无限大的可能,减少了代码的复杂度。 代码例代码一(在编辑器中将"_ "(下划线+空格)替换成两个空格即可编译; 注意要去掉开通的双字节中文空格,蛋疼的百科。) #include <iostream> #include <algorithm> using namespace std; struct point { _ _ int x; _ _ int y; } p[30005], res[30005]; //p标记图中所有的点,res标记凸包上的点 int cmp(point p1, point p2) { _ _ return p1.y < p2.y || (p1.y == p2.y && p1.x < p2.x); } bool ral(point p1, point p2, point p3) //用叉乘判断点的位置 { _ _ return (p2.x - p1.x)*(p3.y - p1.y) > (p3.x - p1.x)*(p2.y - p1.y); } int main() { _ _ int n, i; _ _ while(scanf("%d", &n) != EOF) //一共有n个点 _ _ { _ _ _ _ for(i = 0; i < n; i++) _ _ _ _ _ _ scanf("%d%d", &p[i].x, &p[i].y); __ _ _ if(n == 1) _ _ _ _ { _ _ _ _ _ _ printf("%d %d\", p[0].x, p[0].y); _ _ _ _ _ _ continue; _ _ _ _ } _ _ _ _ if(n == 2) _ _ _ _ { _ _ _ _ _ _ printf("%d %d\", p[0].x, p[0].y); _ _ _ _ _ _ printf("%d %d\", p[1].x, p[1].y); _ _ _ _ _ _ continue; _ _ _ _ } _ _ _ _ sort(p, p + n, cmp); _ _ _ _ res[0] = p[0]; _ _ _ _ res[1] = p[1]; _ _ _ _ int top = 1; _ _ _ _ for(i = 2; i < n; i++) _ _ _ _ { _ _ _ _ _ _ while(top && !ral(res[top], res[top - 1], p[i])) _ _ _ _ _ _ top--; _ _ _ _ _ _ res[++top] = p[i]; _ _ _ _ } _ _ _ _ int len = top; _ _ _ _ res[++top] = p[n - 2]; _ _ _ _ for(i = n - 3; i >= 0; i--) _ _ _ _ { _ _ _ _ _ _ while(top != len && !ral(res[top], res[top - 1], p[i])) _ _ _ _ _ _ top--; _ _ _ _ _ _ res[++top] = p[i]; _ _ _ _ } _ _ _ _ for(i = 0; i < top; i++) _ _ _ _ _ _ printf("%d %d\", res[i].x, res[i].y); //输出凸包上的点 _ _ } _ _ return 0; } 代码二#include <iostream> // 求点集合的凸包的gram算法。n是顶点个数,x,y是顶点 坐标。 #include <fstream> // order 是按照顶点和左下脚的角度的排序后数组。 #include <deque> // tu即是逆时针的凸包上的顶点。 #include <math.h> // using namespace std; //使用条件: 1。点可以任意给,可重复。 // 2。三个以及以上的点。 ifstream fin("input.txt"); // 3。已经考虑了边上有点的情况。 #define NN 1000 #define pi 3.1415827 typedef struct Cseg{ double x,y,tg; }Cseg; int n; double x[NN],y[NN]; deque <Cseg> order; deque <int> tu; Cseg seg1; deque <Cseg> ::iterator p1; deque <int> ::iterator p,q; void in(); void gram(); void makeorder(int s); double dist(double x1,double yy1,double x2,double yy2); double cross(double x1,double yy1,double x2,double yy2); void out(); int main() { in(); gram(); out(); return 0; } void out() { int i; for (i=0;i<tu.size();i++){ cout<<order[tu].x<<" "<<order[tu].y<<endl; } cout<<tu.size()<<" Edges Polydgon"<<endl; return; } void in() { int i; fin>>n; for (i=0;i<n;i++) fin>>x>>y; return; } void gram() { int i,mm; mm=0; for (i=1;i<n;i++) if (y[mm]>y+1e-9) mm=i; else if (fabs(y[mm]-y)<1e-9 && x[mm]>x+1e-9) mm=i; makeorder(mm); seg1.x=x[mm]; seg1.y=y[mm]; tu.clear(); tu.push_back(0); tu.push_back(1); tu.push_back(2); for (i=3;i<order.size();i++){ p=tu.end(); seg1.x=order.x; seg1.y=order.y; p--; q=p-1; if (cross(order[*p].x-order[*q].x,order[*p].y-order[*q].y,order.x-order[* q].x,order.y-order[*q].y)>1e-9) tu.push_back(i); else{ tu.pop_back(); i--; continue; //tu.push_back(i); } }//for return; } void makeorder(int s) { int i; double tg; order.clear(); for (i=0;i<n;i++){ if (i==s) continue; tg=atan2(y-y[s],x-x[s]); seg1.x=x; seg1.y=y; seg1.tg=tg; p1=order.begin(); while (p1!=order.end()){ if (fabs(tg-p1->tg)<1e-9){ if (dist(x[s],y[s],x,y)>dist(x[s],y[s],p1->x,p1->y)+1e-9) { p1->x=x; p1->y=y; } break; } else if (tg<p1->tg){ order.insert(p1,seg1); break; } p1++; }//while if (p1==order.end()) order.insert(p1,seg1); }//for seg1.x=x[s];seg1.y=y[s]; order.insert(order.begin(),seg1); //for (i=0;i<order.size();i++) // printf("i=%d %lf %lf %lf\",i,order.x,order.y,order.tg*180/pi); return; } double cross(double x1,double yy1,double x2,double yy2) { return (x1*yy2-x2*yy1); } double dist(double x1,double yy1,double x2,double yy2) { return pow((x1-x2)*(x1-x2)+(yy1-yy2)*(yy1-yy2),0.5); } 代码三P标程{pku 1113 } {$Q-,S-,R-} const pi=3.1415926575; zero=1e-6; maxn=1000; maxnum=100000000; var ans,temp :extended; n,tot :longint; x,y :array[0..maxn]of extended; zz,num :array[0..maxn]of longint; procedure swap(var ii,jj:extended); var t :extended; begin t:=ii;ii:=jj;jj:=t; end; procedure init; var i,j :longint; begin readln(n,temp); for i:=1 to n do readln(x[i],y[i]); end; function ok(x,midx,y,midy:extended):longint; begin if abs(x-midx)<=zero then begin if abs(midy-y)<=zero then exit(0); if midy>y then exit(1) else exit(2); end else begin if x<midx then exit(1) else exit(2); end; end; procedure qsort(head,tail:longint); var i,j :longint; midx,midy :extended; begin i:=head; j:=tail; midx:=x[(head+tail) div 2]; midy:=y[(head+tail) div 2]; repeat while ok(x[i],midx,y[i],midy)=1 do inc(i); while ok(x[j],midx,y[j],midy)=2 do dec(j); if i<=j then begin swap(x[i],x[j]); swap(y[i],y[j]); inc(i); dec(j); end; until i>j; if i<tail then qsort(i,tail); if j>head then qsort(head,j); end; function Plot(x1,y1,x2,y2:extended):extended; begin Plot:=x1*y2-x2*y1; end; function check(first,last,new:longint):boolean; var ax,ay,bx,by :extended; Pt :extended; begin ax:=x[last]-x[first];ay:=y[last]-y[first]; bx:=x[new]-x[first];by:=y[new]-y[first]; if Plot(ax,ay,bx,by)<-zero then exit(true) else exit(false); end; procedure Tbao; var i,j,tail :longint; begin tot:=0; zz[1]:=1;tail:=1; for i:=2 to n do begin while (zz[tail]<>1)and check(zz[tail-1],zz[tail],i) do dec(tail); inc(tail); zz[tail]:=i; end; inc(tot,tail-1); for i:=1 to tail-1 do num[i]:=zz[i]; zz[1]:=n;tail:=1; for i:=n-1 downto 1 do begin while (zz[tail]<>n)and check(zz[tail-1],zz[tail],i) do dec(tail); inc(tail); zz[tail]:=i; end; for i:=1 to tail-1 do num[tot+i]:=zz[i]; inc(tot,tail-1); end; function dist(a,b:longint):extended; begin dist:=sqrt( (x[a]-x[b])*(x[a]-x[b])+(y[a]-y[b])*(y[a]-y[b]) ); end; procedure main; var i,j :longint; begin qsort(1,n); Tbao; ans:=0; for i:=1 to tot-1 do ans:=ans+dist(num[i],num[i+1]); ans:=ans+dist(num[tot],num[1]); ans:=ans+temp*pi*2; writeln(ans:0:0); end; begin init; main; end. |
随便看 |
百科全书收录4421916条中文百科知识,基本涵盖了大多数领域的百科知识,是一部内容开放、自由的电子版百科全书。