SCOI2015 小凸想跑步
内存限制:256 MiB 时间限制:1000 ms
问题描述
小凸晚上喜欢到操场跑步,今天他跑完两圈之后,他玩起了这样一个游戏。
操场是个凸 n 边形,N 个顶点按照逆时针从 0∼n−1 编号。现在小凸随机站在操场中的某个位置,标记为 P 点。将 P 点与 n 个顶点各连一条边,形成 N 个三角形。如果这时 P 点, 0 号点,1 号点形成的三角形的面积是 N 个三角形中最小的一个,小凸则认为这是一次正确站位。
现在小凸想知道他一次站位正确的概率是多少。
输入格式
第一行包含 1 个整数 n ,表示操场的顶点数和游戏的次数。
接下来有 N 行,每行包含两个整数 Xi 、Yi表示顶点的坐标。
输入保证按逆时针顺序输入点,所有点保证构成一个 n 多边形。所有点保证不存在三点共线。
输出格式
输出一个数,正确站位的概率,保留 4 位小数。
样例输入 1
5
1 8
0 7
0 0
8 0
8 8
样例输出 1
0.6316
样例输入 2
4
0 0
7 0
5 5
-1 4
样例输出 2
0.2475
数据范围
3≤N≤10^5,−10^9≤X,Y≤10^9
按照题意说的去做,列出一些不等式即可。发现P点可以取的范围是若干半平面的交,答案就是两块面积之比。由于数据范围限制,采用O(nlogn)的半平面交。
注意两点:
1.题目给出了P点在凸多边形内,这个限制要加入到用来求答案的半平面里。
2.如果用直线的一般式,就会出现一些较繁琐的特判。为避免这个采用纯向量推导。
下面是推导过程,方法来自Sparrow:
由题意得,(A1−A0)×(P−A1)≤(Ai+1−Ai)×(P−Ai)
由
题
意
得
,
(
A
1
−
A
0
)
×
(
P
−
A
1
)
≤
(
A
i
+
1
−
A
i
)
×
(
P
−
A
i
)
打开括号,整理得P×(Ai+1−Ai+A0−A1)−(A1×A0−Ai+1×Ai)≤0
打
开
括
号
,
整
理
得
P
×
(
A
i
+
1
−
A
i
+
A
0
−
A
1
)
−
(
A
1
×
A
0
−
A
i
+
1
×
A
i
)
≤
0
注意到我们想要的“点P在半平面l内的限制”的形式是(P−l.p)×l.v≤0
注
意
到
我
们
想
要
的
“
点
P
在
半
平
面
l
内
的
限
制
”
的
形
式
是
(
P
−
l
.
p
)
×
l
.
v
≤
0
展开得P×l.v−l.p×l.v≤0,这与上面的式子形式相似
展
开
得
P
×
l
.
v
−
l
.
p
×
l
.
v
≤
0
,
这
与
上
面
的
式
子
形
式
相
似
因此考虑将上面的约束转化为点P在半平面l内。那么l.v=Ai+1−Ai+A0−A1
因
此
考
虑
将
上
面
的
约
束
转
化
为
点
P
在
半
平
面
l
内
。
那
么
l
.
v
=
A
i
+
1
−
A
i
+
A
0
−
A
1
考虑如何求l.p,设k=A1×A0−Ai+1×Ai,那么有l.p×l.v=k
考
虑
如
何
求
l
.
p
,
设
k
=
A
1
×
A
0
−
A
i
+
1
×
A
i
,
那
么
有
l
.
p
×
l
.
v
=
k
之前已经构造出了l.v,k也是已知的,现在也采取构造法弄出l.p
之
前
已
经
构
造
出
了
l
.
v
,
k
也
是
已
知
的
,
现
在
也
采
取
构
造
法
弄
出
l
.
p
设ver为与l.v垂直、模长与l.v相等、且在l.v的顺时针方向的向量
设
v
e
r
为
与
l
.
v
垂
直
、
模
长
与
l
.
v
相
等
、
且
在
l
.
v
的
顺
时
针
方
向
的
向
量
那么有k|ver|2ver×l.v=k恒成立
那
么
有
k
|
v
e
r
|
2
v
e
r
×
l
.
v
=
k
恒
成
立
所以,l.p=k|ver|2ver
所
以
,
l
.
p
=
k
|
v
e
r
|
2
v
e
r
这样做就避免了特判。
下面的代码给整体下标都加了一。
#include<stdio.h>
#include<algorithm>
#include<cstring>
#include<cmath>
#define MAXN 200005
#define db double
using namespace std;
int N;
db S1,S2;
struct Point{
db x,y;
}A[MAXN];
typedef Point Vector;
struct Line{
Vector v;
Point p;
db ang;
void GetAng(){
ang=atan2(v.y,v.x);
}
}L[MAXN];int tot;
Vector operator+(Vector a,Vector b){
a.x+=b.x;a.y+=b.y;
return a;
}
Vector operator-(Vector a,Vector b){
a.x-=b.x;a.y-=b.y;
return a;
}
Vector operator*(db k,Vector a){
a.x*=k;a.y*=k;
return a;
}
db operator*(Vector a,Vector b){
return a.x*b.y-a.y*b.x;
}
bool operator<(Line a,Line b){
if(a.ang==b.ang)
return (a.p-b.p)*b.v>0;
return a.ang<b.ang;
}
db sqlen(Vector a){
return a.x*a.x+a.y*a.y;
}
bool onleft(Line b,Point a){
return (a-b.p)*b.v<=0;
}
Point Insc(Line A,Line B){
db k=(A.p-B.p)*A.v/(B.v*A.v);
return k*B.v+B.p;
}
Line Q[MAXN];
Point P[MAXN];int cnt;
void HPI(){
int i,t=0,j,head,tail;
sort(L+1,L+tot+1);
for(i=1;i<=tot;i++){
if(i==1||L[i-1].ang!=L[i].ang)t++;
L[t]=L[i];
}
head=1;tail=1;Q[tail]=L[1];
for(i=2;i<=t;i++){
while(head<tail&&!onleft(L[i],Insc(Q[tail],Q[tail-1])))tail--;
while(head<tail&&!onleft(L[i],Insc(Q[head],Q[head+1])))head++;
Q[++tail]=L[i];
}
while(head<tail&&!onleft(Q[head],Insc(Q[tail],Q[tail-1])))tail--;
if(head==tail)return;
Q[tail+1]=Q[head];
for(i=head;i<=tail;i++)P[++cnt]=Insc(Q[i],Q[i+1]);
}
int main(){
int i;
scanf("%d",&N);
for(i=1;i<=N;i++)scanf("%lf%lf",&A[i].x,&A[i].y);A[N+1]=A[1];
for(i=1;i<=N;i++)S1+=A[i]*A[i+1];
for(i=1;i<=N;i++)L[++tot].v=A[i+1]-A[i],L[tot].p=A[i];
Vector ver;db k;
for(i=2;i<=N;i++){
tot++;
L[tot].v=A[i+1]-A[i]-A[2]+A[1];
ver.x=L[tot].v.y;ver.y=-L[tot].v.x;
k=A[i+1]*A[i]-A[2]*A[1];
L[tot].p=-k/sqlen(ver)*ver;
}
for(i=1;i<=tot;i++)L[i].GetAng();
HPI();
if(cnt)for(i=1,P[cnt+1]=P[1];i<=cnt;i++)S2+=P[i]*P[i+1];
printf("%.4lf\n",S2/S1);
}