scoiday1T3&&bzoj4445小凸想跑步

9 篇文章 0 订阅

几何题啊,据说是裸的半平面交,but我还是WA了好几发。


逆时针给定一些点,求多边形内到第一条边面积比到其它边小的点的概率(用范围/总面积)


然后可以设该点为(x,y),到所有边的面积表示出来,然后建立不等关系,化简约分后我们可以得到一系列关于x,y的不等式,形如ax+by<c,然后就搞一搞半平面交。


然而有两个bug我没考虑到,一个就是要保证(x,y)的点在凸多边形内部,怎么保证?先把所有的半平面表示成直线,直线左侧表示半平面,然后所有边加进来,也构成一串不等关系。

我没考虑这种情况,所出来概率偏大,为什么?    因为这个点的范围到了多边形外侧啊。


上面一种方案对应第一种,空间是第二种的两倍,时间差不多也是.


那第二种怎么搞,我们可以只把第一条边的限制加进去,其余边都不管,为什么这么做是可行的呢?  可以感性的理解一下,如果在其他边外侧,它一定不可能是到第一条边最小了,所以在不等关系的约束条件之外,根本就搞不动了。(可以适当画图验证)

这样的时间和空间都节省了。



第二个bug然我搞了很久,因为之前样例已经可以过了,然而交了几发都WA了(又贡献了这么多个WA,真是汗颜!)

后来不断生成数据(其实最后是自己手推了一组数据),发现自己的输出居然是"nan",恩,一定是哪里挂了,然后我把数据带进去调试,发现我计算不等关系ax+by+c<0的时候,b有可能为0了,这种情况下我用的起始点分母有个b,就除以了0,就非法了,然后就GG了。


调出第二个bug之后果断又来了一发,结果就A了.


然后无聊就来找找出bug的原因(又交了几发,又WA了几次)。

在第一个bug中,如果要所有边建立不等式需要把空间开成2倍。


搞出了一系列不等关系,直接套一下半平面交。

#include<cstdio>
#include<algorithm>
#include<iostream>
#include<cmath>
#define eps 1e-10
using namespace std;
const int maxn=200000+20;
int dcmp(double x)
{
	if(fabs(x)<eps)return 0;
	else return x<0?-1:1;
}
struct point
{
	double x,y;
	point(){}
	point(double _x,double _y){x=_x;y=_y;}
	point operator +(const point &b){return point(x+b.x,y+b.y);}
	point operator -(const point &b){return point(x-b.x,y-b.y);}
	point operator *(const double &b){return point(x*b,y*b);}
};
struct line
{
	point p,v;
	double ang;
	line(){}
	line(point _p,point _v){p=_p;v=_v;ang=atan2(v.y,v.x);}
	bool operator <(const line &l)const
	{
		return ang<l.ang;
	}
};
double cross(point a,point b){return a.x*b.y-a.y*b.x;}
double dot(point a,point b){return a.x*b.x+a.y*b.y;}
bool onleft(line l,point p){return dcmp(cross(l.v,p-l.p))>0;}
/*
(y1+y4-y2-y3)x+(x2+x3-x1-x4)y+x1y2-x2y1+x4y3-x3y4<0

*/
point getsec(point p1,point v1,point p2,point v2)
{
	double x=cross(p2-p1,v2)/cross(v1,v2);
	return p1+v1*x;
}
double pans;
line l[maxn];
point t[maxn];
line q[maxn];
void Hp(line *L,int n)
{
	sort(L+1,L+n+1);
	int first,last;
	point p[maxn];		
	q[last=first=0]=L[1];
	for(int i=2;i<=n;i++)
	{
		while(first<last&&!onleft(L[i],p[last-1]))last--;
		while(first<last&&!onleft(L[i],p[first]))first++;
		q[++last]=L[i];
		if(dcmp(cross(q[last].v,q[last-1].v))==0)
		{
			last--;
			if(!onleft(L[i],p[last-1]))q[last]=L[i];
			//等价于if(onleft(q[last],L[i].p))q[last]=L[i]; 
		}
		if(first<last)p[last-1]=getsec(q[last].p,q[last].v,q[last-1].p,q[last-1].v);
	}
	while(first<last&&!onleft(q[first],p[last-1]))last--;
	if(last-first<=1)
	{
		printf("0.0000\n");
		return ;
	}
	p[last]=getsec(q[last].p,q[last].v,q[first].p,q[first].v);
	double ans=0;
	for(int i=first+1;i<last;i++)
	{
		ans+=fabs(cross(p[i]-p[first],p[i+1]-p[first]));
	}
	ans*=0.5;
	printf("%.4lf\n",ans/pans);
}
int cnt;
void pre(int n)
{
	cnt=0;
	double x1=t[1].x,y1=t[1].y;
	double x2=t[2].x,y2=t[2].y;
	for(int i=2;i<=n;i++)
	{
		double x3=t[i].x,y3=t[i].y;
		double x4=t[i+1].x,y4=t[i+1].y;
		double a=y1+y4-y2-y3;
		double b=x2+x3-x1-x4;
		double c=x1*y2-x2*y1+x4*y3-x3*y4;
		//ax+by+c<0
		//方向向量(-b,a); 
		//使得ax+by+c=0左边表示ax+by+c	
		point p;if(dcmp(b)!=0)p=point(0,-c/b);
		else p=point(-c/a,0);
		point v;v=point(-b,a);
		l[++cnt]=line(p,v); 
	}
	//l[++cnt]=line(t[1],t[2]-t[1]);//
	for(int i=1;i<=n;i++)l[++cnt]=line(t[i],t[i+1]-t[i]);
}
int main()
{
	int n;
	scanf("%d",&n);
	for(int i=1;i<=n;i++)scanf("%lf%lf",&t[i].x,&t[i].y);
	pans=0;
	for(int i=2;i<n;i++)
	{
		pans+=fabs(cross(t[i]-t[1],t[i+1]-t[1]));
	}
	pans*=0.5;
	t[n+1]=t[1];
	pre(n);
	Hp(l,cnt);
	
	return 0;
}


  • 1
    点赞
  • 0
    收藏
    觉得还不错? 一键收藏
  • 1
    评论
评论 1
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值