UOJ Logo zhouzixuan的博客

博客

半平面交

2016-03-28 22:17:39 By zhouzixuan
【写在前面】:
    很久之前就想填的坑
    既然市选都考挂了,就来填填坑好了
【参考】:

http://www.csie.ntnu.edu.tw/~u91029/Half-planeIntersection.html

http://blog.csdn.net/non_cease/article/details/7820361

【定义】:
    什么是半平面?顾名思义就是半个平面
    更准确一点的定义:二维平面上一条直线将平面分成两个部分,其中一个部分就是半平面
    如果有多条直线,那么这些直线对应半平面的交集就是半平面交
    注意,对于题目给出的直线,有可能不能形成有界的区域,但是我们可以根据题目数据范围认为的规定一个非常大的边界即可
    如何快速求出半平面交,以及有什么用处呢?这就是以下将涉及的内容

【增量法】:
    考虑我们认为加入边界后,使得半平面交一定是一个凸多边形(很好证明:因为你根本画不出凹多边形的情况)
    因此,我么考虑加入一个新的半平面,然后对于凸多边形上的每一个点判断是否在该半平面中,然后判断该直线和凸多边形的边是否有交点,并且把存在于半平面中的点和交点加入新的多边形集合中,最终构造出的凸多边形就是加入该半平面后的半平面交
    复杂度O(N^2) 算法是在线的
【随机增量法】:
    在增量法的基础上加上对初始半平面顺序的随机洗牌
    可以证明:每回合期望加入2个点,删除n/i个点,因此复杂度为O(N log N)
    然而这样就不是在线的算法了
【队列法】:
    其实我并不知道这种算法的名字,好像是朱泽园自创的?
    1.初始将两个半平面加入双端队列
    2.加入新的半平面时,从左端和右端分别弹出,直到两个半平面构成的交点在当前半平面中,然后加入新的半平面
    3.现在需要将双端队列两端的半平面拼起来,因此需要对两边分别判断,直到一段的前两个的交点在另一端队首半平面内
    复杂度O(N log N)
【直线方程求半平面交】:
    在给出一些直线求线性规划可行区域时往往需要注意Ax+By+C=0中A B C的符号问题
    因此我们不妨考虑取直线上的两个点作为向量,然后默认向量的左边为半平面
    因此需要根据A B C讨论一下
【例题】:
    poj1755 半平面交求线性规划可行域
poj 1755
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <algorithm>
#include <cmath>
using namespace std;
const int N=105;
const double zero=1e-16;
const double INF=1e9;
struct point
{
    double x,y;
    inline point operator - (const point &a) {return (point){x-a.x,y-a.y};}
    inline point operator + (const point &a) {return (point){x+a.x,y+a.y};}
    inline point operator * (double b) {return (point){x*b,y*b};}
    inline point operator / (double b) {return (point){x/b,y/b};}
    inline double operator * (const point &a) {return x*a.y-y*a.x;}
    inline double operator ^ (const point &a) {return x*a.x+y*a.y;}
} P[N];
struct Node {double x,y,z;} speed[N];
struct node
{
    point p1,p2; double angle; node() {}
    node (point _p1,point _p2) {p1=_p1; p2=_p2;}
    inline void calc() {angle=atan2(p2.y-p1.y,p2.x-p1.x);}
} p[N],q[N];
int n,tot; bool ans;
inline long long getint()
{
    long long x=0; char c=getchar(); bool flag=false;
    while ((c!='-')&&((c<'0')||(c>'9'))) c=getchar();
    if (c=='-') flag=true,c=getchar();
    while ((c>='0')&&(c<='9')) x=x*10+(long long)(c-'0'),c=getchar();
    if (flag) return -x; else return x;
}
void init()
{
    n=getint();
    for (int i=1; i<=n; i++) speed[i].x=getint(),speed[i].y=getint(),speed[i].z=getint();
}
inline int dcmp(double x)
{
    if (x>zero) return 1;
    if (x<-zero) return -1;
    return 0;
}
inline double cross(point a,point b,point c)
{
    double nowans=(b-a)*(c-a);
    return nowans;
}
inline point intersect(node a,node b)
{
    double s1=cross(b.p1,a.p2,a.p1);
    double s2=cross(a.p2,b.p2,a.p1);
    return (b.p1*s2+b.p2*s1)/(s1+s2);
}
inline bool judge(node a,node b,node c)
{
    point now=intersect(b,c);
    return dcmp(cross(a.p1,a.p2,now))>0;
}
inline bool parllel(node a,node b)
{
    double nowans=(a.p2-a.p1)*(b.p2-b.p1);
    return dcmp(nowans)==0;
}
inline bool MakeHalf(double a,double b,double c)
{
    int da=dcmp(a),db=dcmp(b),dc=dcmp(c);
    if ((da==0)&&(db==0)) return dc>0;
    if (db==0) {p[++tot]=(node){(point){-c/a,da},(point){-c/a,-da}}; return true;}
    if (db>0) p[++tot]=(node){(point){0,-c/b},(point){1,(-c-a)/b}};
    if (db<0) p[++tot]=(node){(point){1,(-c-a)/b},(point){0,-c/b}};
    return true;
}
inline bool cmp(const node &a,const node &b)
{
    int d=dcmp(a.angle-b.angle); if (d!=0) return d<0;
    return dcmp(cross(a.p1,a.p2,b.p2))<0;
}
bool HalfSection()
{
    int head=1,tail=0;
    for (int i=1; i<=tot; i++) p[i].calc();
    sort(p+1,p+tot+1,cmp); int cnt=1;
    for (int i=2; i<=tot; i++) if (dcmp(p[i].angle-p[i-1].angle)!=0) p[++cnt]=p[i];
    for (int i=1; i<=cnt; i++)
    {
        while ((head<tail)&&(!judge(p[i],q[tail],q[tail-1]))) tail--;
        while ((head<tail)&&(!judge(p[i],q[head],q[head+1]))) head++;
        if ((head<=tail)&&(parllel(p[i],q[tail]))) return false;
        q[++tail]=p[i];
    }
    while ((head+1<tail)&&(!judge(q[head],q[tail-1],q[tail]))) tail--;
    while ((head+1<tail)&&(!judge(q[tail],q[head+1],q[head]))) head++;
    return tail-head+1>2;
}
void solve()
{
    for (int i=1; i<=n; i++)
    {
        tot=0; ans=true;
        p[++tot]=(node){(point){0,0},(point){INF,0}};
        p[++tot]=(node){(point){INF,0},(point){INF,INF}};
        p[++tot]=(node){(point){INF,INF},(point){0,INF}};
        p[++tot]=(node){(point){0,INF},(point){0,0}};
        for (int j=1; j<=n; j++)
        {
            if (i==j) continue;
            //x*(1/a[i]-1/a[j]-1/c[i]+1/c[j])+y*(1/b[i]-1/b[j]-1/c[i]+1/c[j])+1/c[i]-1/c[j]>0
            double a=1/speed[j].x-1/speed[i].x;
            double b=1/speed[j].y-1/speed[i].y;
            double c=1/speed[j].z-1/speed[i].z;
            ans=ans&&MakeHalf(a,b,c);
        }
        ans=ans&&HalfSection();
        if (ans) printf("Yes\n"); else printf("No\n");
    }
}
int main()
{
    init();
    solve();
    return 0;
}

评论

暂无评论

发表评论

可以用@mike来提到mike这个用户,mike会被高亮显示。如果你真的想打“@”这个字符,请用“@@”。