POJ 2451 Uyuw's Concert(半平面交nlgn)

2015-01-24 06:35:50 · 作者: · 浏览: 4
题目:给出一些向量,求出围成的多边形的核的面积
?本身没啥好说的,就是半平面交的模板题,听说 是Zzy为了他的那篇论文专门出的题,但是数据被POJ更改过,时限很宽,n^2随便水,而且我的nlgn比n^2耗时还多。
先引用一下nlgn算法的介绍,在zzy论文里也有详细介绍
step1. 将所有半平面按极角排序,对于极角相同的,选择性的保留一个。 O(nlogn)
step2. 使用一个双端队列(deque),加入最开始2个半平面。
step3. 每次考虑一个新的半平面:
? a.while deque顶端的两个半平面的交点在当前半平面外:删除deque顶端的半平面
? b.while deque底部的两个半平面的交点在当前半平面外:删除deque底部的半平面
? c.将新半平面加入deque顶端
step4.删除两端多余的半平面。
具体方法是:
a.while deque顶端的两个半平面的交点在底部半平面外:删除deque顶端的半平面
b.while deque底部的两个半平面的交点在顶端半平面外:删除deque底部的半平面
重复a,b直到不能删除为止。
step5:计算出deque顶端和底部的交点即可。
精度我是一个个试出来的。。。最终是1e-10
[cpp]?
#include?
#include?
#include?
#include?
#include?
#include?
#include?
#include?
#include?
#include?
#include?
#include?
#include?
#include?
#include?
#include?
#define LL long long?
#define eps 1e-10?
#define inf 10000?
#define zero(a) fabs(a) #define N 20005?
using namespace std;?
struct Point{?
??? double x,y;?
}p[N*2];?
struct Segment{?
??? Point s,e;?
??? double angle;?
??? void get_angle(){angle=atan2(e.y-s.y,e.x-s.x);}?
}seg[N];?
int m;?
//叉积为正说明,p2在p0-p1的左侧?
double xmul(Point p0,Point p1,Point p2){?
??? return (p1.x-p0.x)*(p2.y-p0.y)-(p2.x-p0.x)*(p1.y-p0.y);?
}?
Point Get_Intersect(Segment s1,Segment s2){?
??? double u=xmul(s1.s,s1.e,s2.s),v=xmul(s1.e,s1.s,s2.e);?
??? Point t;?
??? t.x=(s2.s.x*v+s2.e.x*u)/(u+v);t.y=(s2.s.y*v+s2.e.y*u)/(u+v);?
??? return t;?
}?
bool cmp(Segment s1,Segment s2){?
??? //先按极角排序?
??? if(s1.angle>s2.angle) return true;?
??? //极角相等,内侧的在前?
??? else if(zero(s1.angle-s2.angle)&&xmul(s2.s,s2.e,s1.e)>-eps) return true;?
??? return false;?
}?
void HalfPlaneIntersect(Segment seg[],int n){?
??? sort(seg,seg+n,cmp);?
??? int tmp=1;?
??? for(int i=1;i ??????? if(!zero(seg[i].angle-seg[tmp-1].angle))?
??????????? seg[tmp++]=seg[i];?
??? n=tmp;?
??? Segment deq[N];?
??? deq[0]=seg[0];deq[1]=seg[1];?
??? int head=0,tail=1;?
??? for(int i=2;i ??????? while(head ??????? while(head ??????? deq[++tail]=seg[i];?
??? }?
??? while(head ??? while(head ??? if(head==tail) return;?
??? m=0;?
??? for(int i=head;i ??????? p[m++]=Get_Intersect(deq[i],deq[i+1]);?
??? if(tail>head+1)?
??????? p[m++]=Get_Intersect(deq[head],deq[tail]);?
}?
double Get_area(Point p[],int &n){?
??? double area=0;?
??? for(int i=1;i ??????? area+=xmul(p[0],p[i],p[i+1]);?
??? return fabs(area)/2.0;?
}?
int main(){?
??? int n;?
??? while(scanf("%d",&n)!=EOF){?
??????? seg[0].s.x=0;seg[0].s.y=0;seg[0].e.x=10000;seg[0].e.y=0;seg[0].get_angle();?
??????? seg[1].s.x=10000;seg[1].s.y=0;seg[1].e.x=10000;seg[1].e.y=10000;seg[1].get_angle();?
??????? seg[2].s.x=10000;seg[2].s.y=10000;seg[2].e.x=0;seg[2].e.y=10000;seg[2].get_angle();?
??????? seg[3].s.x=0;seg[3].s.y=10000;seg[3].e.x=0;seg[3].e.y=0;seg[3].get_angle();?
??????? for(int i=0;i ??????????? scanf("%lf%lf%lf%lf",&seg[i+4].s.x,&seg[i+4].s.y,&seg[i+4].e.x,&seg[i+4].e.y);?
??????????? seg[i+4].get_angle();?
??????? }?
??????? HalfPlaneIntersect(seg,n+4);?
??????? printf("%.1f\n",Get_area(p,m));?
??? }?
??? return 0;?
}?