首页 > 代码库 > 关于半平面交

关于半平面交

嗯,这是一个很屌的东西。可以把他想象成数学中的线性规划问题,然后自然而然得想到就可以求最优解啦。

如何求解半平面交?????

做法一:暴力枚举点,用该点切割现有的凸多边形,这样的复杂度是O(n^2)

做法二:神奇的分治O(nlogn),然而我并不会。。。。

做法三:参见2006年朱泽园大神发明的排序增量法:这里就暂时不给详细介绍了。

复杂度O(nlogn),如果把快速排序改成基数排序,复杂度可以进一步降为:O(n)!!!!!!

半平面交的另外一个应用是求多边形的核:可以想象成站在一个多边形的点上,

该点能看到多边形的每一个部位的点集。这个直接求出半平面的面积就可以了。

#include<bits/stdc++.h>
#define N 2000
#define il inline
#define LL long long
#define db double
using namespace std;
struct point{
  db x,y;
  point() {}
  point(db X,db Y):x(X),y(Y) {}
  db operator *(const point & a)const{
    return x*a.y-y*a.x;
  }
  point operator +(const point & a)const{
    return point(x+a.x,y+a.y);
  }
  point operator -(const point & a)const{
    return point(x-a.x,y-a.y);
  }
  point operator *(const db & k)const{
    return point(x*k,y*k);
  }
}p[N];
struct line{
  point a,b;db Ang;
  line() {}
  line(point A,point B,db Aa):a(A),b(B),Ang(Aa) {}
}l[N];
int n;
#define eps 1e-8
il int CC(double kk){
  if(fabs(kk)<eps)return 0;
  return kk>0?1:-1;
}
il bool comp(const line & l1,const line & l2){
  int ff=CC(l1.Ang-l2.Ang);
  if(!ff)return CC((l2.a-l1.a)*(l2.b-l1.a))<0;
  return l1.Ang<l2.Ang;
}
void getnode(line l1,line l2,point & p){
  point a=l1.a,b=l1.b,c=l2.a,d=l2.b;
  db k=((c-a)*(d-c))/((b-a)*(d-c));
  p=a+(b-a)*k;
}
int de[N];
bool check(line l0,line l1,line l2){
  point p;getnode(l1,l2,p);
  return CC((l0.b-l0.a)*(p-l0.a))<0;
}
void HalfplaneIntersect(){
  int j=0;
  sort(l+1,l+n+1,comp);
  l[++j]=l[1];
  for(int i=2;i<=n;++i)
    if(CC(l[i].Ang-l[j].Ang)>0)l[++j]=l[i];
  n=j;de[1]=1,de[2]=2;
  int tail(2),head(1);
  for(int i=3;i<=n;++i){
    while(tail>head&&(!check(l[i],l[de[tail]],l[de[tail-1]])))tail--;
    while(tail>head&&(!check(l[i],l[de[head]],l[de[head+1]])))head++;
    de[++tail]=i;
  }while(tail>head&&(!check(l[de[head]],l[de[tail]],l[de[tail-1]])))tail--;
  while(tail>head&&(!check(l[de[tail]],l[de[head]],l[de[head+1]])))head++;
  int cnt=0;de[++tail]=de[head];
  for(int i=head;i<tail;++i)cnt++,getnode(l[de[i]],l[de[i+1]],p[cnt]);
  n=cnt;
}
db Ans(){
  db cnt=0.00;
  if(n<3)return 0;
  for(int i=2;i<n;++i)cnt+=((p[i]-p[1])*(p[i+1]-p[1]));
  cnt=fabs(cnt);cnt/=2;return cnt;
}
int main(){
  int T;scanf("%d",&T);
  while(T--){
    scanf("%d",&n);
    for(int i=1;i<=n;++i)scanf("%lf%lf",&p[i].x,&p[i].y);
    for(int i=1;i<n;++i)l[i]=line(p[i],p[i+1],atan2(p[i+1].y-p[i].y,p[i+1].x-p[i].x));
    l[n]=line(p[n],p[1],atan2(p[1].y-p[n].y,p[1].x-p[n].x));
    HalfplaneIntersect();printf("%.2lf\n",Ans());
  }
  return 0;
}

  

关于半平面交