【BZOJ1038】【codevs1412】瞭望塔，半平面交/三分法

传送门1
传送门2
思路：
一开始看的时候感觉并没有什么直截了当的想法，想起在abclzr的博客里见到过这道题目，分类是半平面交，然后就往半平面交的方面上去想，发现塔顶的可行区域是一段下凸壳及其上面的区域，然后就有了一个想法，把可行区域求出来，然后最优选择一定是轮廓上的点往上到下凸壳或凸壳上的点往下到轮廓，当时的我只感觉这样是对的，造了一个数据，发现没有反例，然后就“大胆猜想，从不证明”地写了一发，而且很蠢地搞了一发$O(n^2)$枚举轮廓和凸壳上的点，交上去竟然A了，不过回头看看，发现轮廓到凸壳的距离实际上只是一个很蠢的由一坨直线组成的分段函数，轮廓和凸壳上的折点就是分段函数的端点，分段函数的最优值一定是在端点上取到的；而且求完半平面交后，存的凸壳上的点横坐标递增，所以可以$O(n)$扫一遍求出来的。
abclzr表示ljw(男)在去年的SD夏令营讲过这道题目，呃，蒟蒻表示没听过，也可能听过忘了吧。。。
后来看到BZOJ的discuss里有人说还可以三分，表示非常神奇，想了一会发现，对于轮廓上的每一段线段，它对应的分段函数区间一定是单峰的，因为它上方的下凸壳的线段连续且斜率递增，所以距离只能是递减，递增或先减后增，所以我们枚举每一段轮廓，在它上面三分一下位置，然后取较优值就可以了，复杂度$O(n^2\log S)$，S是指坐标范围。
这道题目还有模拟退火的做法，蒟蒻表示没学过，不过一题多解这种想法在平常练习中还是很重要的，有利于拓展思维，增长见识，从不同的角度考虑问题，或许能得到新的思路和方法
$O(n^2)$代码：

#include<cstdio>
#include<iostream>
#include<cmath>
#include<algorithm>
#define eps 1e-8
using namespace std; 
int n,be,ed;
int X[305],Y[305];
double ans=1e18;
struct Point{
    double x,y;
    Point(double X=0,double Y=0){x=X;y=Y;}
};
typedef Point Vec;
Vec operator +(Vec a,Vec b){return Vec(a.x+b.x,a.y+b.y);}
Vec operator -(Vec a,Vec b){return Vec(a.x-b.x,a.y-b.y);}
Vec operator *(Vec a,double p){return Vec(a.x*p,a.y*p);}
double cross(Vec a,Vec b){return a.x*b.y-a.y*b.x;}
struct Line{
    Point p;
    Vec v;
    double ang;
    Line(Point P=Point(0,0),Vec V=Vec(0,0)){p=P;v=V;ang=atan2(v.y,v.x);}
    bool operator <(const Line other)const
    {
        if (ang==other.ang) return v.x<other.v.x;
        return ang<other.ang;
    }
}line[305],q[305];
Point cp[305];
int dcmp(double x)
{
    if (fabs(x)<=eps) return 0;
    return x>0?1:-1;
}
bool onleft(Line a,Point b)
{
    return dcmp(cross(b-a.p,a.v))<=0;
}
Point inter(Line a,Line b)
{
    Point u=a.p-b.p;
    double t=cross(u,b.v)/cross(b.v,a.v);
    return a.p+a.v*t;
}
void SI(int n)
{
    sort(line+1,line+n+1);
    int head=1,tail=0;
    for (int i=1;i<=n;++i)
    {
        while (head<tail&&!onleft(line[i],cp[tail-1])) --tail;
        while (head<tail&&!onleft(line[i],cp[head])) ++head;
        if (head<=tail&&dcmp(cross(line[i].v,q[tail].v))==0)
            if (onleft(q[tail],line[i].p)) q[tail]=line[i];
            else; 
        else
            q[++tail]=line[i];
        if (head<tail) cp[tail-1]=inter(q[tail],q[tail-1]);
    }
    while (head<tail&&!onleft(q[head],cp[tail-1])) --tail;
    be=head,ed=tail-1;
    if (dcmp(cross(q[head].v,q[tail].v))==0);
    else ++ed,cp[tail]=inter(q[head],q[tail]);
}
main()
{
    scanf("%d",&n);
    for (int i=1;i<=n;++i) scanf("%d",X+i);
    for (int i=1;i<=n;++i) scanf("%d",Y+i);
    for (int i=2;i<=n;++i)
        line[i-1]=Line(Point(X[i-1],Y[i-1]),Vec(X[i]-X[i-1],Y[i]-Y[i-1]));
    line[n]=Line(Point(X[1],Y[1]),Vec(0,-1e6));
    line[n+1]=Line(Point(X[n],Y[n]),Vec(0,1e6));
    SI(n+1);
    double t;
    for (int i=1;i<=n;++i)
        for (int j=be;j<ed;++j)
            if (cp[j].x<=X[i]&&X[i]<=cp[j+1].x)
                ans=min(ans,t=fabs(inter(Line(cp[j],cp[j+1]-cp[j]),Line(Point(X[i],Y[i]),Vec(0,1))).y-Y[i]));

    for (int j=be;j<=ed;++j)
        for (int i=1;i<n;++i)
            if (X[i]<=cp[j].x&&cp[j].x<=X[i+1])
                ans=min(ans,fabs(cp[j].y-inter(Line(Point(X[i],Y[i]),Vec(X[i+1]-X[i],Y[i+1]-Y[i])),Line(cp[j].x,Vec(0,-1))).y));
    printf("%.3lf\n",ans);
}