半平面交

最新推荐文章于 2019-07-20 18:59:07 发布

穿袜子的流氓兔

最新推荐文章于 2019-07-20 18:59:07 发布

阅读量689

点赞数

分类专栏： ACM 几何文章标签：半平面交

ACM 几何专栏收录该内容

14 篇文章 0 订阅

订阅专栏

首先解决问题：什么是半平面？顾名思义，半平面就是指平面的一半，我们知道，一条直线可以将平面分为两个部分，那么这两个部分就叫做两个半平面。

然后，半平面怎么表示呢？二维坐标系下，直线可以表示为ax + by + c = 0，那么两个半平面则可以表示为ax + by + c >= 0 和ax + by + c < 0，这就是半平面的表示方法。

还有，半平面的交是神马玩意？其实就是一个方程组，让你画出满足若干个式子的坐标系上的区域（类似于线性规划的可行域），方程组就是由类似于上面的这些不等式组成的。

另外，半平面交可以干什么？半平面交虽然说是半平面的问题，但它其实就是关于直线的问题。一个一个的半平面其实就是一个一个有方向的直线而已。

半平面交的一个重要应用就是求多边形的核。多边形的核又是神马玩意？它是平面简单多边形的核是该多边形内部的一个点集，该点集中任意一点与多边形边界上一点的连线都处于这个多边形内部。就是一个在一个房子里面放一个摄像头，能将所有的地方监视到的放摄像头的地点的集合即为多边形的核。经常会遇到让你判定一个多边形是否有核的问题。

比如说：

POJ 3335 Rotating Scoreboard
http://acm.pku.edu.cn/JudgeOnline/problem?id=3335

POJ 1474 Video Surveillance
http://acm.pku.edu.cn/JudgeOnline/problem?id=1474

POJ 1279 Art Gallery
http://acm.pku.edu.cn/JudgeOnline/problem?id=1279

这三个题比较简单，裸的半平面交，下面是核心代码（暨半平面交的模板）：

/*半平面相交（直线切割多边形）（点标号从1开始）*/
Point points[MAXN],p[MAXN],q[MAXN];
int n;
double r;
int cCnt,curCnt;
inline void getline(Point x,Point y,double &a,double &b,double &c){<span style="white-space:pre">	//得到直线ax+by+c=0
    a = y.y - x.y;
    b = x.x - y.x;
    c = y.x * x.y - x.x * y.y;
}
inline void initial(){
    for(int i = 1; i <= n; ++i)p[i] = points[i];
    p[n+1] = p[1];
    p[0] = p[n];
    cCnt = n;
}
inline Point intersect(Point x,Point y,double a,double b,double c){<span style="white-space:pre">	//得到直线ax+by+c=0与点X和点Y所连直线的交点
    double u = fabs(a * x.x + b * x.y + c);
    double v = fabs(a * y.x + b * y.y + c);
    return Point( (x.x * v + y.x * u) / (u + v) , (x.y * v + y.y * u) / (u + v) );
}
inline void cut(double a,double b ,double c){<span style="white-space:pre">				//用直线ax+by+c=0切割多边形
    curCnt = 0;
    for(int i = 1; i <= cCnt; ++i){
        if(a*p[i].x + b*p[i].y + c >= EPS)q[++curCnt] = p[i];<span style="white-space:pre">		</span>//点集按顺时针排序，如果ax1+by1+c>0的话说明点（x1,y1）在直线右边，保留
        else {
            if(a*p[i-1].x + b*p[i-1].y + c > EPS){			//如果某一点不在多边形内核内，但是该点和其左右点连成的直线与直线ax+by+c=0有交点
                q[++curCnt] = intersect(p[i],p[i-1],a,b,c);             //那么保留交点
            }
            if(a*p[i+1].x + b*p[i+1].y + c > EPS){
                q[++curCnt] = intersect(p[i],p[i+1],a,b,c);
            }
        }
    }
    for(int i = 1; i <= curCnt; ++i)p[i] = q[i];
    p[curCnt+1] = q[1];p[0] = p[curCnt];
    cCnt = curCnt;
}
inline void solve(){
    //注意：默认点是顺时针，如果题目不是顺时针，规整化方向
    initial();
    for(int i = 1; i <= n; ++i){
        double a,b,c;
        getline(points[i],points[i+1],a,b,c);
        cut(a,b,c);
    }
    /*
    如果要向内推进r，用该部分代替上个函数
    for(int i = 1; i <= n; ++i){
        Point ta, tb, tt;
        tt.x = points[i+1].y - points[i].y;
        tt.y = points[i].x - points[i+1].x;
        double k = r / sqrt(tt.x * tt.x + tt.y * tt.y);
        tt.x = tt.x * k;
        tt.y = tt.y * k;
        ta.x = points[i].x + tt.x;
        ta.y = points[i].y + tt.y;
        tb.x = points[i+1].x + tt.x;
        tb.y = points[i+1].y + tt.y;
        double a,b,c;
        getline(ta,tb,a,b,c);
        cut(a,b,c);
    }
    */
    //多边形核的面积
    double area = 0;
    for(int i = 1; i <= curCnt; ++i)
        area += p[i].x * p[i + 1].y - p[i + 1].x * p[i].y;
    area = fabs(area / 2.0);
    //此时cCnt为最终切割得到的多边形的顶点数，p为存放顶点的数组
}
inline void GuiZhengHua(){
     //规整化方向，逆时针变顺时针，顺时针变逆时针
    for(int i = 1; i < (n+1)/2; i ++)
      swap(points[i], points[n-i]);//头文件加iostream
}
inline void init(){
     for(int i = 1; i <= n; ++i)points[i].input();
        points[n+1] = points[1];
}

半平面交的最终奥义就是求出一个满足条件的凸多边形，而解决这个问题的前提就是直线切割多边形，让直线不断的去切割当前多边形，然后得到新的多边形，然后继续。。。最终得到一个符合条件的多边形，如果最终的多边形顶点数目为0，那么就说明半平面交无解（半平面交的解可以是一个凸多边形，一条直线，一个点，我们用顶点来记录这些解）。关于直线切割多边形，流程是这样的：对凸多边形（指的是当前多边形）的每一个顶点，如果这个顶点在直线的指定的一侧（暨在该半平面上），那么就将该顶点直接加入到新的多边形中，否则，看与该点相邻的多边形上的两个点（判断线段是否和直线相交），如果和直线相交，则把交点加入到新的多边形中。这样，就可以得到一个新的凸多边形。关于最初的多边形，我们可以设一个四方形的区域，比如说（-1000，-1000），（-1000,1000），（1000,1000），（1000,-1000），然后开始切割。

POJ 3525 Most Distant Point from the Sea （推荐）http://acm.pku.edu.cn/JudgeOnline/problem?id=3525求在多边形中画一个圆的最大半径可以使用半平面交+ 二分法解。对每个距离进行判定（是否存在区域可以放置圆心，多边形的的每条边向内推进

#include <cstdio>
#include <cstring>
#include <cmath>
#include <iostream>
using namespace std;
#define EPS 1e-10
const int MAXN = 200;
struct Point{
    double x,y;
    Point(){}
    Point(double _x,double _y):x(_x),y(_y){}
    void input(){
        scanf("%lf%lf",&x,&y);
    }
};
Point points[MAXN],p[MAXN],q[MAXN];
int n;
int cCnt,curCnt;
inline void getline(Point x,Point y,double &a,double &b,double &c){
     a = y.y - x.y;
    b = x.x - y.x;
    c = y.x * x.y - x.x * y.y;
}
inline void initial(){
    for(int i = 1; i <= n; ++i)p[i] = points[i];
    p[n+1] = p[1];
    p[0] = p[n];
    cCnt = n;
}
inline Point intersect(Point x,Point y,double a,double b,double c){
     double u = fabs(a * x.x + b * x.y + c);
    double v = fabs(a * y.x + b * y.y + c);
    return Point( (x.x * v + y.x * u) / (u + v) , (x.y * v + y.y * u) / (u + v) );
}
inline void cut(double a,double b ,double c){
    curCnt = 0;
    for(int i = 1; i <= cCnt; ++i){
        if(a*p[i].x + b*p[i].y + c >= 0)q[++curCnt] = p[i];
        else {
            if(a*p[i-1].x + b*p[i-1].y + c > 0){
                q[++curCnt] = intersect(p[i],p[i-1],a,b,c);
            }
            if(a*p[i+1].x + b*p[i+1].y + c > 0){
                q[++curCnt] = intersect(p[i],p[i+1],a,b,c);
            }
        }
    }
    for(int i = 1; i <= curCnt; ++i)p[i] = q[i];
    p[curCnt+1] = q[1];p[0] = p[curCnt];
    cCnt = curCnt;
}
inline int solve(double r){
    initial();
    for(int i = 1; i <= n; ++i){
        Point ta, tb, tt;
        tt.x = points[i+1].y - points[i].y;
        tt.y = points[i].x - points[i+1].x;
        double k = r / sqrt(tt.x * tt.x + tt.y * tt.y);
        tt.x = tt.x * k;
        tt.y = tt.y * k;
        ta.x = points[i].x + tt.x;
        ta.y = points[i].y + tt.y;
        tb.x = points[i+1].x + tt.x;
        tb.y = points[i+1].y + tt.y;
        double a,b,c;
        getline(ta,tb,a,b,c);
        cut(a,b,c);
    }
    if(cCnt <= 0)return 0;
    return 1;
}
int main(){
    while(scanf("%d",&n) != EOF){
        if(n == 0)break;
        for(int i = 1; i <= n; ++i)points[i].input();
       for(int i = 1; i < (n+1)/2; i ++)
            swap(points[i], points[n-i]);
        points[n+1] = points[1];
        double high = 10000000;
        double low = 0,mid;
        while(low + EPS <= high){
            mid = (high + low)/2.0;
            if(solve(mid))low = mid;
            else high = mid;
        }
        printf("%lf/n",high);
    }
    return 0;
}

POJ 3384 Feng Shui （推荐）
http://acm.pku.edu.cn/JudgeOnline/problem?id=3384
半平面交实际应用，用两个圆覆盖一个多边形，问最多能覆盖多边形的面积。
解法：用半平面交将多边形的每条边一起向“内”推进R，得到新的多边形，然后求多边形的最远两点

#include <cstdio>
#include <cstring>
#include <cmath>
#include <iostream>
using namespace std;
#define EPS 1e-10
const int MAXN = 300;
struct Point{
    double x,y;
    Point(){}
    Point(double _x,double _y):x(_x),y(_y){}
    void input(){
        scanf("%lf%lf",&x,&y);
    }
};
Point points[MAXN],p[MAXN],q[MAXN];
int n;
double r;
int cCnt,curCnt;
void getline(Point x,Point y,double &a,double &b,double &c){
 a = y.y - x.y;
 b = x.x - y.x;
 c = y.x * x.y - x.x * y.y;
}
inline void initial(){
    for(int i = 1; i <= n; ++i)p[i] = points[i];
    p[n+1] = p[1];
    p[0] = p[n];
    cCnt = n;
}
inline Point intersect(Point x,Point y,double a,double b,double c){
     double u = fabs(a * x.x + b * x.y + c);
    double v = fabs(a * y.x + b * y.y + c);
    return Point( (x.x * v + y.x * u) / (u + v) , (x.y * v + y.y * u) / (u + v) );
}
inline void cut(double a,double b ,double c){
    curCnt = 0;
    for(int i = 1; i <= cCnt; ++i){
        if(a*p[i].x + b*p[i].y + c >= 0)q[++curCnt] = p[i];
        else {
            if(a*p[i-1].x + b*p[i-1].y + c > 0){
                q[++curCnt] = intersect(p[i],p[i-1],a,b,c);
            }
            if(a*p[i+1].x + b*p[i+1].y + c > 0){
                q[++curCnt] = intersect(p[i],p[i+1],a,b,c);
            }
        }
    }
    for(int i = 1; i <= curCnt; ++i)p[i] = q[i];
    p[curCnt+1] = q[1];p[0] = p[curCnt];
    cCnt = curCnt;
}
inline double pdis(Point a,Point b){
    return sqrt((a.x-b.x)*(a.x-b.x)+(a.y-b.y)*(a.y-b.y));
}
int main(){
    scanf("%d%lf",&n,&r);
    for(int i = 1; i <= n; ++i)points[i].input();
    points[n+1] = points[1];
    initial();
    for(int i = 1; i <= n; ++i){
        Point ta, tb, tt;
        tt.x = points[i+1].y - points[i].y;
        tt.y = points[i].x - points[i+1].x;
        double k = r / sqrt(tt.x * tt.x + tt.y * tt.y);
        tt.x = tt.x * k;
        tt.y = tt.y * k;
        ta.x = points[i].x + tt.x;
        ta.y = points[i].y + tt.y;
        tb.x = points[i+1].x + tt.x;
        tb.y = points[i+1].y + tt.y;
        double a,b,c;
        getline(ta,tb,a,b,c);
        cut(a,b,c);
    }
    int ansx = 0,ansy = 0;
    double res = 0;
    for(int i = 1; i <= cCnt; ++i){
        for(int j = i + 1; j <= cCnt; ++j){
            double tmp = pdis(p[i],p[j]);
            if(tmp > res){
                res = tmp;
                ansx = i;
                ansy = j;
            }
        }
    }
    printf("%.4lf %.4lf %.4lf %.4lf/n",p[ansx].x,p[ansx].y,p[ansy].x,p[ansy].y);
    return 0;
}

穿袜子的流氓兔

关注

0
点赞
踩
0

收藏

觉得还不错? 一键收藏
0
评论
半平面交

首先解决问题：什么是半平面？顾名思义，半平面就是指平面的一半，我们知道，一条直线可以将平面分为两个部分，那么这两个部分就叫做两个半平面。然后，半平面怎么表示呢？二维坐标系下，直线可以表示为ax + by + c = 0，那么两个半平面则可以表示为ax + by + c >= 0 和ax + by + c 还有，半平面的交是神马玩意？其实就是一个方程组，让你画出满足若干个式子的坐标系上
复制链接

扫一扫