1082. 【GDOI2005】选址

Description

很久以前，在世界的某处有一个形状为凸多边形的小岛，岛上的居民们决定建一个祭坛，居民们任务祭坛的位置离岛的顶点处越远越好。 你的任务是求凸多边形内一点，使其与各顶点的距离中最短的距离最远，点在边上也可以。 这样的点可能有多个，你只需输出这些点与各顶点的最短距离。

Solution

先考虑最后的答案在多边形内部。

可以证明，这个点到至少 3 个点的距离相同。

• 假设答案点只到 2 个点距离相同，那么这个点在这两个点的垂直平分线上，并且可以向外移动只到与另外的 1 个点距离相同。

基于上面的结论，我们可以枚举三角形，求出外心，如果外心在多边形内就可以直接计算答案。如果在外面，就取三条垂直平分线与多边形边的交点计算答案。

判断是否在内部可以用向量叉乘。

时间复杂度 \(\mathcal O(n^4)\)。

Code

#include<cmath>
#include<cstdio>
#include<algorithm>
#define N 105
#define inf 2147483647.0
#define eps 1e-60
#define db long double
using namespace std;
struct Point
{
    db x,y;
}p[N];
int n,a,b;
db ans;
db sqr(db x) {return x*x;}
bool check(db k1,db b1,db k2,db b2,db k3,db b3,db x,db y)
{
    int bjsum=0,bj[N];
    db xl[N][3],cc[N];
    for (int i=1;i<=n;++i)
        xl[i][1]=x-p[i].x,xl[i][2]=y-p[i].y;
    xl[n+1][1]=xl[1][1];
    xl[n+1][2]=xl[1][2];
    for (int i=1;i<=n;++i)
        cc[i]=xl[i][1]*xl[i+1][2]-xl[i][2]*xl[i+1][1];
    for (int i=1;i<=n;++i)
        if (cc[i]>0)
        {
            bj[i]=1;
            bjsum++;
        }
        else bj[i]=0;
    if (bjsum==0||bjsum==n) return true;
    for (int i=1;i<=n;++i)
    {
        if (bjsum==1&&bj[i]==1)
        {
            a=i,b=i+1;
            return false;
        }
        if (bjsum!=n&&bj[i]==0)
        {
            a=i;b=i+1;
            return false;
        }
    }
    if (!a&&!b)
    {
        for (int i=1;i<=n;++i)
        {
            if (bjsum!=0&&bj[i]==1) 
            {
                a=i;b=i+1;
                if (b>n) b=1;
                db kk=0,bb=0,xx=0,yy=0;
                if (abs(p[a].x-p[b].x)<eps) kk=inf;
                else kk=(p[a].y-p[b].y)/(p[a].x-p[b].x);
                if (kk!=inf) bb=p[a].y-kk*p[a].x;
                if (abs(kk-inf)<eps) xx=p[a].x,yy=k1*xx+b1;
                else xx=(b1-bb)/(kk-k1),yy=kk*xx+bb;
                if (check(k1,b1,k2,b2,k3,b3,xx,yy)) return false;
            }
            if (bjsum!=n&&bj[i]==0)
            {
                a=i;b=i+1;
                if (b>n) b=1;
                db kk=0,bb=0,xx=0,yy=0;
                if (abs(p[a].x-p[b].x)<eps) kk=inf;
                else kk=(p[a].y-p[b].y)/(p[a].x-p[b].x);
                if (kk!=inf) bb=p[a].y-kk*p[a].x;
                if (abs(kk-inf)<eps) xx=p[a].x,yy=k1*xx+b1;
                else xx=(b1-bb)/(kk-k1),yy=kk*xx+bb;
                if (check(k1,b1,k2,b2,k3,b3,xx,yy)) return false;
            }
        }
    }
    return false;
}
void calc(db xx,db yy)
{
    db sum=inf;
    for (int i=1;i<=n;++i)
        sum=min(sum,sqrt(sqr(p[i].x-xx)+sqr(p[i].y-yy)));
    ans=max(ans,sum);
}
int main()
{
    freopen("position.in","r",stdin);
    freopen("position.out","w",stdout);
    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)
        for (int j=i+1;j<=n;++j)
            for (int k=j+1;k<=n;++k)
            {
                db midx1,midx2,midy1,midy2,midx3,midy3,k1,k2,k3,b1,b2,b3,x1,x2,x3,y1,y2,y3;
                db pointx,pointy;
                x1=p[i].x;y1=p[i].y;
                x2=p[j].x;y2=p[j].y;
                x3=p[k].x;y3=p[k].y;
                midx1=(x1+x2)/2;midy1=(y1+y2)/2;
                midx2=(x2+x3)/2;midy2=(y2+y3)/2;
                midx3=(x3+x1)/2;midy3=(y3+y1)/2;
                if (y2==y1) k1=inf;
                else k1=(x1-x2)/(y2-y1);
                if (y3==y2) k2=inf;
                else k2=(x2-x3)/(y3-y2);
                if (y3==y1) k3=inf;
                else k3=(x3-x1)/(y1-y3);
                if (k1!=inf) b1=midy1-k1*midx1;
                if (k2!=inf) b2=midy2-k2*midx2;
                if (k3!=inf) b3=midy3-k3*midx3;
                if (abs(k1-k2)<eps) continue;
                if (k1==inf) pointx=midx1,pointy=k2*midx1+b2;
                if (k2==inf) pointx=midx2,pointy=k1*midx2+b1;
                if (k1!=inf&&k2!=inf)
                {
                    pointx=(b2-b1)/(k1-k2);
                    pointy=(k1*b2-k2*b1)/(k1-k2);
                }
                a=b=0;
                if (check(k1,b1,k2,b2,k3,b3,pointx,pointy))
                {
                    calc(pointx,pointy);
                }
                else
                {
                    db kk,bb,xx,yy;
                    if (b>n) b=1;
                    if (abs(p[a].x-p[b].x)<eps) kk=inf;
                    else kk=(p[a].y-p[b].y)/(p[a].x-p[b].x);
                    if (kk!=inf) bb=p[a].y-kk*p[a].x;
                    if (abs(kk-inf)<eps) xx=p[a].x,yy=k1*xx+b1;
                    else xx=(b1-bb)/(kk-k1),yy=kk*xx+bb;
                    calc(xx,yy);
                    if (abs(kk-inf)<eps) xx=p[a].x,yy=k2*xx+b2;
                    else xx=(b2-bb)/(kk-k2),yy=kk*xx+bb;
                    calc(xx,yy);
                    if (abs(kk-inf)<eps) xx=p[a].x,yy=k3*xx+b3;
                    else xx=(b3-bb)/(kk-k3),yy=kk*xx+bb;
                    calc(xx,yy);
                }
            }
    printf("%.3Lf\n",ans);
    return 0;
}