poj3525 Most Distant Point from the Sea

发表于 | 更新于

Link

Solution

找一个半径最大的内切圆。 二分r，将凸包上的边向内移动，求半平面交看是否为空。 移动边的时候比较巧妙，用到了单位法向量。

Code

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
//Code by Lucida
#include<cstdio>
#include<cstring>
#include<algorithm>
#include<cmath>
#define red(x) scanf("%d",&x)
#define fred(x) scanf("%lf",&x)
template <class T> inline bool chkmx(T &a,const T &b){return a<b?a=b,1:0;}
template <class T> inline bool chkmn(T &a,const T &b){return a>b?a=b,1:0;}
const int MAXN=100+1;
typedef double ld;
using std::sort;
const ld eps=1e-7,undefined=1e250,pi=acos(-1.0);
int fcmp(ld x)
{
	if(-eps<x && x<eps) return 0;
	return x<0?-1:1;
}
template <class T> T sqr(T x){return x*x;}
template <class T> T abs(T x){return x>0?x:-x;}
struct vec
{
	ld x,y;
	vec(ld _x=0,ld _y=0):x(_x),y(_y){}
	vec operator -(){return vec(-x,-y);}
	ld norm(){return sqrt(x*x+y*y);}
	bool operator <(vec a)
	{
		if(x!=a.x) return x<a.x;
		return y<a.y;
	}
	void operator +=(vec a){x+=a.x,y+=a.y;}
	void operator -=(vec a){x-=a.x,y-=a.y;}
	void operator /=(ld a){x/=a,y/=a;}
	void operator *=(ld a){x*=a,y*=a;}
	bool operator ==(vec a){return x==a.x && y==a.y;}
};
typedef vec point;
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,ld t){return vec(a.x*t,a.y*t);}
vec operator /(vec a,ld t){return vec(a.x/t,a.y/t);}
ld inner(vec a,vec b){return a.x*b.x+a.y*b.y;}
ld outer(vec a,vec b){return a.x*b.y-a.y*b.x;}
vec rotate(vec a,ld rad){return vec(a.x*cos(rad)-a.y*sin(rad),a.x*sin(rad)+a.y*cos(rad));}
struct line
{
	point a;vec v;
	ld angle;
	line(){}
	line(point _a,vec _v):a(_a),v(_v){angle=atan2(v.y,v.x);}
	vec norm()
	{
		vec normal=rotate(v,pi/2);
		normal/=normal.norm();
		return normal;
	}
	bool operator<(const line &a)const {return angle<a.angle;}
}il[MAXN];
bool onleft(point a,line l){return fcmp(outer(l.v,a-l.a))>=0;}//在左或者在上
bool parl(line a,line b){return fcmp(outer(a.v,b.v))==0;}
point cross(line l1,line l2)
{
	vec v=l1.a-l2.a;
	ld t=outer(v,l2.v)/outer(l2.v,l1.v);
	return l1.a+l1.v*t;
}
int n;
int intersect(line *l,point *tio)
{
	sort(l+1,l+1+n);
	static point Qp[MAXN];
	static line Ql[MAXN];
	int head=1,tail=0;
	Ql[++tail]=l[1];
	#define count(x,y) (y-x+1)
	for(int i=2;i<=n;i++)
	{
		while(count(head,tail)>=2 && !onleft(Qp[tail-1],l[i])) tail--;
		while(count(head,tail)>=2 && !onleft(Qp[head],l[i])) head++;
		Ql[++tail]=l[i];
		if(parl(Ql[tail-1],Ql[tail]))
		{
			tail--;
			if(onleft(l[i].a,Ql[tail]))
				Ql[tail]=l[i];
		}
		if(count(head,tail)>=2) Qp[tail-1]=cross(Ql[tail-1],Ql[tail]);
	}
	while(count(head,tail)>=2 && !onleft(Qp[tail-1],Ql[head])) tail--;
	if(count(head,tail)<2) return 0;
	Qp[tail]=cross(Ql[head],Ql[tail]);
	int tol=0;
	for(int i=head;i<=tail;i++) tio[++tol]=Qp[i];
	return tol;
	#undef count
}
bool check(ld r)
{
	static line tr[MAXN];
	for(int i=1;i<=n;i++) tr[i]=il[i];
	for(int i=1;i<=n;i++)
		tr[i].a+=tr[i].norm()*r;
	static point tio[MAXN];
	int tol=intersect(tr,tio);
	return tol>2;
}
void WORK()
{
	static point p[MAXN];
	for(int i=1;i<=n;i++)
		fred(p[i].x),fred(p[i].y);
	for(int i=1;i<=n;i++)
		il[i]=line(p[i],p[i+1>n?1:i+1]-p[i]);
	ld L=0,R=1e10;
	while(fcmp(R-L))
	{
		ld MID=(L+R)/2;
		if(check(MID))
			L=MID;
		else
			R=MID;
	}
	printf("%.7lf\n",L);
}
int main()
{
	freopen("input","r",stdin);
	while(red(n),n) WORK();
	return 0;
}