POJ 3525 Most Distant Point from the Sea 二分+半平面交

链接：http://poj.org/problem?id=3525

题意：给一个凸多边形的海岛，寻找海岛之中距离海边距离最长的一个点的距离。

思路：求凸多边形的最大内切圆。做法是二分半径加半平面交，将凸多边形的每条边向内部（垂直方向）收缩半径r，看每条边的半平面是否还会交出凸多边形。

P.S. 无意间找到这道题，顺便检验一下自己敲的模板。找一个好的模板真的很重要。

代码：

#include <iostream>
#include <cstdio>
#include <cmath>
#include <algorithm>
#include <vector>
#include <complex>
#include <queue>
#define eps 1e-10
#define INF 10000
using namespace std;
typedef complex < double > Point;
typedef pair < Point , Point > Halfplane;
inline double dist(Halfplane hal)
{
return sqrt((hal.second.imag()-hal.first.imag())*(hal.second.imag()-hal.first.imag())+(hal.second.real()-hal.first.real())*(hal.second.real()-hal.first.real()));
}
inline int sgn(double n)
{
return fabs(n)<eps?0:(n<0?-1:1);
}
inline double cross(Point a,Point b)
{
return (conj(a)*b).imag();
}
inline double dot(Point a,Point b)
{
return(conj(a)*b).real();
}
inline double satisfy(Point a,Halfplane p)
{
return sgn(cross(a-p.first,p.second-p.first))<=0;
}
Point crosspoint (const Halfplane &a,const Halfplane &b)
{
double k=cross(b.first-b.second,a.first-b.second);
k=k/(k-cross(b.first-b.second,a.second-b.second));
return a.first+(a.second-a.first)*k;
}
bool cmp(const Halfplane &a,const Halfplane &b)
{
int res= sgn(arg(a.second-a.first)-arg(b.second-b.first));
return res==0?satisfy(a.first,b):res<0;
}
vector < Point > halfplaneIntersection(vector<Halfplane>v)
{
sort(v.begin(),v.end(),cmp);
deque < Halfplane > q;
deque < Point > ans;
q.push_back(v[0]);
for(int i=1; i<int(v.size()); i++)
{
if(sgn(arg(v[i].second - v[i].first)-arg(v[i-1].second - v[i-1].first))==0)
continue;
while(ans.size()>0 && !satisfy(ans.back(),v[i]))
{
ans.pop_back();
q.pop_back();
}
while(ans.size()>0 && !satisfy(ans.front(),v[i]))
{
ans.pop_front();
q.pop_front();
}
ans.push_back(crosspoint (q.back(),v[i]));
q.push_back(v[i]);
}
while(ans.size() > 0 && !satisfy(ans.back(),q.front()))
{
ans.pop_back();
q.pop_back();
}
while(ans.size() > 0 && !satisfy(ans.front(),q.back()))
{
ans.pop_front();
q.pop_front();
}
ans.push_back(crosspoint(q.back(),q.front()));
return vector < Point > (ans.begin(),ans.end());
}
int main()
{
int tot;
while(scanf("%d",&tot))
{
vector < Halfplane > pol;
if(tot==0)
return 0;
double fa,fb;
scanf("%lf%lf",&fa,&fb);
double la=fa,lb=fb;
for(int i=1; i<tot; i++)
{
double a,b;
scanf("%lf%lf",&a,&b);
pol.push_back (Halfplane(Point(la,lb),Point(a,b)));
la=a,lb=b;
}
pol.push_back(Halfplane(Point(la,lb),Point(fa,fb)));
double l=0,r=20000;
while(r-l>eps)
{
double mid=(l+r)/2.0;
vector < Halfplane >res;
vector < Halfplane > ::iterator it;
for(it=pol.begin(); it!=pol.end(); it++)
{
double xx=((*it).first.imag()-(*it).second.imag())/dist((*it))*mid;
double yy=((*it).second.real()-(*it).first.real())/dist((*it))*mid;
res.push_back(Halfplane(Point((*it).first.real()+xx,(*it).first.imag()+yy),Point((*it).second.real()+xx,(*it).second.imag()+yy)));
}
vector < Point > aa = halfplaneIntersection(res);
if(aa.size()<3)
r=mid;
else l=mid;
}
printf("%.6lf\n",(l+r)/2.0);
}
return 0;
}