III.[CQOI2006]凸多边形 /【模板】半平面交
半平面交开始~~
这里介绍一种做法:随机增量法,其可以在\(O(n\log n)\)的时间内完成半平面交的求解。
首先,我们可以用一个向量来表示直线;之所以使用向量来表示,是因为我们将强制该向量的左方表示半平面。
接着,我们考虑将所有直线按照向量的倾斜角从小往大排序(此处介绍一个很有用的函数atan2(y,x)
,其值等于向量\(\vec{v}=(x,y)\)的倾斜角(弧度制)。运用该函数即可轻松排序)。对于倾角相同的向量,我们保留最左方的一个。
排序以后就是关键部分了:
我们令一个双端队列储存当前半平面交中所有的直线。我们再令另一个长度为之前双端队列长度减一的双端队列,储存之前双端队列中两两相邻直线的交点。
当插入一条新直线时:
假如此时双端队列中,最左端的两条直线平行(因为保证倾角不同,所以只有可能方向相反),则无解。同理,最右端平行也会有同样结论。
之后,假如点的双端队列两端的点离开了新插入的直线,就直接连点带边一起pop
出去即可。
在全部pop
完毕之后,将新直线插入队列尾,并将新直线与之前队列尾直线的交点插入点的双端队列。
当全部直线插入完毕后,算法并未结束——毕竟,此时队列首尾的直线还可能有交点。于是,我们用队尾直线判断队首的点,再用队首直线判断队尾的点即可。
最终,如果队列中剩余\(0\)条或\(1\)条直线,则半平面交无解;
否则,计算首尾直线的交点并将其插入点的队列,则点的队列中所有点即构成凸包边界。
时间复杂度瓶颈在于排序。
代码:
#include<bits/stdc++.h>
using namespace std;
const double eps=1e-9;
const double pi=acos(-1);
double cmp(double x){
if(x>eps)return 1;
if(x<-eps)return -1;
return 0;
}
struct Vector{
double x,y;
Vector(double X=0,double Y=0){x=X,y=Y;}
friend Vector operator +(const Vector &u,const Vector &v){return Vector(u.x+v.x,u.y+v.y);}
friend Vector operator -(const Vector &u,const Vector &v){return Vector(u.x-v.x,u.y-v.y);}
friend Vector operator *(const Vector &u,const double &v){return Vector(u.x*v,u.y*v);}
friend Vector operator /(const Vector &u,const double &v){return Vector(u.x/v,u.y/v);}
friend double operator &(const Vector &u,const Vector &v){return u.x*v.y-u.y*v.x;}//cross times
friend double operator |(const Vector &u,const Vector &v){return u.x*v.x+u.y*v.y;}//point times
double operator ~()const{return sqrt(x*x+y*y);}//the modulo of a vector
double operator !()const{return atan2(y,x);}//the angle of a vector
void read(){scanf("%lf%lf",&x,&y);}
};
typedef Vector Point;
double area(Point *a,int sz){//the area of a figure(read points in the clockwise order or vice versa)
double ret=0;
for(int i=1;i<sz;i++)ret+=a[i]&a[i+1];
ret+=a[sz]&a[1];
return fabs(ret/2);
}
struct Line{//the left side of a line stands for the semiplane
Point x,y,z;//z stands for the vector of the line
Line(Point X=Point(),Point Y=Point()){x=X,y=Y,z=Y-X;}
friend bool operator <(const Line &u,const Line &v){//sort by angle, then by leftwardness
int k=cmp(!u.z-!v.z);
if(k==-1)return true;
if(k==1)return false;
return cmp(u.z&(v.y-u.x))==-1;
}
friend bool operator /(const Line &u,const Line &v){return cmp(u.z&v.z)==0;}
//check if two lines are parallel(of the same or opposite direction)
friend bool operator |(const Line &u,const Line &v){return cmp(!u.z-!v.z)==0;}//check if two lines have the same angle
friend bool operator &(const Line &u,const Point &v){return cmp((v-u.x)&u.z)!=1;}//if a point lies inside a semiplane
friend Point operator &(const Line &u,const Line &v){//calculate the intersection (u parallel to v MUSTN'T HOLD)
return u.x+u.z*((v.z&(u.x-v.x))/(u.z&v.z));
}
};
typedef Line Segment;
int m,n,L,R;
Point p[100100],res[100100];
Line l[100100];
Line a[100100];
int Semiplane_Intersection(){
sort(a+1,a+n+1);
l[L=R=1]=a[1];
for(int i=2;;i++){
while(i<=n&&a[i]|a[i-1])i++;//eliminate all but one line which hold the same angle
if(i>n)break;
if(L<R&&(l[L]/l[L+1]||l[R]/l[R-1]))return -1;
//check if there're two parallel lines(as they don't have the same angle, they are of opposite directions, which implies an empty union)
while(L<R&&!(a[i]&p[R-1]))R--;//eliminate all the off-set points, both left and right
while(L<R&&!(a[i]&p[L]))L++;
l[++R]=a[i];
if(L<R)p[R-1]=l[R]&l[R-1];
}
//considering the line on the front and the line on the back, they have an intersection and thus elimination could take place
while(L<R&&!(l[L]&p[R-1]))R--;
while(L<R&&!(l[R]&p[L]))L++;
if(R-L<=1)return -1;//a union doesn't exist when there is only one or less point.
p[R]=l[L]&l[R];
for(int i=L;i<=R;i++)res[i-L+1]=p[i];
return R-L+1;//the size of the intersection
}
int main(){
scanf("%d",&m);
for(int i=1,q;i<=m;i++){
scanf("%d",&q);
for(int j=1;j<=q;j++)p[j].read();
for(int j=1;j<q;j++)a[++n]=Line(p[j],p[j+1]);
a[++n]=Line(p[q],p[1]);
}
int len=Semiplane_Intersection();
if(len==-1){puts("0.000");return 0;}
printf("%.3lf\n",area(res,len));
return 0;
}