[BZOJ4445][Scoi2015]小凸想跑步（半平面交）

题目描述

传送门

题目大意：一个凸n边形，N个顶点按照逆时针从0～n-l编号。现在小凸随机选择多边形中的某个位置，标记为P点。将P点与n个顶点各连一条边，形成N个三角形。如果这时P点，0号点，1号点形成的三角形的面
积是N个三角形中最小的一个，小凸则认为这是一次正确站位。现在小凸想知道他一次站位正确的概率是多少。

题解

设选的位置坐标为(x,y)，用叉积化一下式子
$(x_1-x)(y_2-y)-(y_1-y)(x_2-x)\le (x_k-x)(y_{k+1}-y)-(x_{k+1}-x)(y_k-y)$
$(y_1-y_2-y_k+y_{k+1})x+(x_2-x_1-x_{k+1}+x_k)y+(x_1y_2-x_2y_1-x_ky_{k+1}+x_{k+1}y_k)\le 0$
这就是一个$ax+by+c\le 0$的直线形式，根据abc的符号判断直线的方向，特判ab为0的情况
然后所有的直线交出来一个凸壳，求这个凸壳的面积/原先的面积就行了
这题有点卡精度，eps开1e-15

代码

#include<algorithm>
#include<iostream>
#include<cstring>
#include<cstdio>
#include<cmath>
using namespace std;
#define N 200005

const double eps=1e-15;
const double  inf=1e18;
int dcmp(double x)
{
    if (x<=eps&&x>=-eps) return 0;
    return (x>0)?1:-1;
}
struct Vector
{
    double x,y;
    Vector(double X=0,double Y=0)
    {
        x=X,y=Y;
    }
};
typedef Vector Point;
struct Line
{
    Point p;Vector v;
    double ang;
    Line(Point P=Point(0,0),Vector V=Vector(0,0))
    {
        p=P,v=V;
        ang=atan2(v.y,v.x);
    }
    bool operator < (const Line &a) const
    {
        return ang<a.ang;
    }
};
Vector operator + (Vector a,Vector b){return Vector(a.x+b.x,a.y+b.y);}
Vector operator - (Vector a,Vector b){return Vector(a.x-b.x,a.y-b.y);}
Vector operator * (Vector a,double b){return Vector(a.x*b,a.y*b);}

int n,line,cnt,l,r;
double area,ans;
Line L[N<<1],q[N<<1];
Point pt[N],p[N<<1],poly[N<<1];

double Cross(Vector a,Vector b){return a.x*b.y-a.y*b.x;}
bool Onleft(Line l,Point P)
{
    Vector v=P-l.p;
    return dcmp(Cross(l.v,v))>=0;
}
Point GLI(Point P,Vector v,Point Q,Vector w)
{
    Vector u=P-Q;
    double t=Cross(w,u)/Cross(v,w);
    return P+v*t;
}
void init()
{
    L[++line]=Line(Point(inf,inf),Vector(-2*inf,0));
    L[++line]=Line(Point(-inf,inf),Vector(0,-2*inf));
    L[++line]=Line(Point(-inf,-inf),Vector(2*inf,0));
    L[++line]=Line(Point(inf,-inf),Vector(0,2*inf));
}
void halfp()
{
    sort(L+1,L+line+1);
    q[l=r=1]=L[1];
    for (int i=2;i<=line;++i)
    {
        while (l<r&&!Onleft(L[i],p[r-1]))
            --r;
        while (l<r&&!Onleft(L[i],p[l]))
            ++l;
        q[++r]=L[i];
        if (!dcmp(Cross(q[r].v,q[r-1].v)))
        {
            --r;
            if (Onleft(q[r],L[i].p))
                q[r]=L[i];
        }
        if (l<r) p[r-1]=GLI(q[r].p,q[r].v,q[r-1].p,q[r-1].v);
    }
    while (l<r&&!Onleft(q[l],p[r-1]))
        --r;
    if (r-l<=1) return;
    p[r]=GLI(q[r].p,q[r].v,q[l].p,q[l].v);
    for (int i=l;i<=r;++i) poly[++cnt]=p[i];
}
double Area(int n,Point *P)
{
    double ans=0;
    for (int i=2;i<n;++i)
        ans+=Cross(P[i]-P[1],P[i+1]-P[1]);
    return fabs(ans/2);
}
int main()
{
    scanf("%d",&n);
    for (int i=1;i<=n;++i) scanf("%lf%lf",&pt[i].x,&pt[i].y);
    area=Area(n,pt);
    for (int i=2;i<=n;++i)
    {
        int id=i,jd=i%n+1;
        double a=pt[1].y-pt[2].y-pt[id].y+pt[jd].y;
        double b=pt[2].x-pt[1].x-pt[jd].x+pt[id].x;
        double c=pt[1].x*pt[2].y-pt[2].x*pt[1].y-pt[id].x*pt[jd].y+pt[jd].x*pt[id].y;
        if (!dcmp(a)&&!dcmp(b)) continue;
        if (!dcmp(a))
        {
            if (dcmp(b)>0)
                L[++line]=Line(Point(inf,-c/b),Vector(-2*inf,0));
            else
                L[++line]=Line(Point(-inf,-c/b),Vector(2*inf,0));
            continue;
        }
        if (!dcmp(b))
        {
            if (dcmp(a)>0)
                L[++line]=Line(Point(-c/a,-inf),Vector(0,2*inf));
            else
                L[++line]=Line(Point(-c/a,inf),Vector(0,-2*inf));
            continue;
        }
        double K=-a/b,B=-c/b;
        double x1=0,x2=1;
        double y1=x1*K+B,y2=x2*K+B;
        Point P=Point(x1,y1),Q=Point(x2,y2);
        if (dcmp(-b)>0)
            L[++line]=Line(P,Q-P);
        else
            L[++line]=Line(Q,P-Q);
    }
    for (int i=1;i<=n;++i)
        L[++line]=Line(pt[i],pt[i%n+1]-pt[i]);
    init();
    halfp();
    ans=Area(cnt,poly);
    printf("%.4lf\n",ans/area);
}