省选专练SCOI2015小凸想跑步

终于快把SCOI练完了。。。

毒瘤题

这个题难点不在半平面交，在于推导:

其次精度问题

eps 和INF 取值需谨慎。

大了会丢精度，小了会炸精度。

接着得开long double

于是流输出。

#include<bits/stdc++.h>
using namespace std;
#define double long double
const double INF=1e16;
const double eps=1e-16;
const int N=4e5;
/*int cmp(double x){
	if(fabs(x)<eps){
		return 0;
	}
	else return x>0?1:-1;
}
bool equal (double x,double y){
	return (cmp(x-y)==0);
}
struct Point{
	double x,y;
	Point (double _x=0,double _y=0):x(_x),y(_y){}
	friend Point operator + (Point A,Point B){return Point(A.x+B.x,A.y+B.y);}
	friend Point operator - (Point A,Point B){return Point(A.x-B.x,A.y-B.y);}
	friend Point operator * (Point A,double k){return Point(A.x*k,A.x*k);}
	friend Point operator / (Point A,double k){return Point(A.x/k,A.y/k);}
	double ang(){return (atan2(y,x));}
	void read(){scanf("%Lf%Lf",&x,&y);}
}p[N],a[N];
typedef Point Vector;
struct Line{
	Point p;
	Vector v;
	Line (Point _P=Point(0,0),Point _v=Vector(0,0)):p(_P),v(_v){}
	double Ang;
	double ang(){return v.ang();}
}l[N],convex_hull[N],temp[N];
bool comp(Line A,Line B){
	return A.Ang<B.Ang;
}
double Cross(Vector A,Vector B){
	return A.x*B.y-A.y*B.x;
}
bool in(Point p,Line l){
	return cmp(Cross(l.v,p-l.p))>=0;
}
Point getpot(Line A,Line B){
	return B.p+B.v*(Cross(A.v,A.p-B.p)/Cross(A.v,B.v));
}*/
struct Point{
	double x,y;
	Point(double _x=0,double _y=0):x(_x),y(_y){}
	friend Point operator +(Point A,Point B){return Point(A.x+B.x,A.y+B.y);}
	friend Point operator -(Point A,Point B){return Point(A.x-B.x,A.y-B.y);}
	friend Point operator *(Point A,double k){return Point(A.x*k,A.y*k);}
	friend Point operator /(Point A,double k){return Point(A.x/k,A.y/k);}
	double ang(){return atan2(y,x);}
	void read(){scanf("%lf%lf",&x,&y);x*=1.0;y*=1.0;}
}p[N],a[N];
int cmp(double x){
	if(fabs(x)<eps)return 0;
	return x>eps?1:-1;
}
bool equal(double x,double y){
	return cmp(x-y)==0;
}
typedef Point Vector;
struct Line{
	Point p;
	Vector v;
	Line(Point _p=Point(0,0),Vector _v=Vector(0,0)):v(_v),p(_p){}
	double Ang;
	double ang(){return v.ang();}
}l[N],temp[N],convex_hull[N];
bool comp(Line A,Line B){
	return A.Ang<B.Ang;
}
double Cross(Vector A,Vector B){
	return(double)(A.x*B.y-A.y*B.x);
}
Point getpot(Line A,Line B){
	return B.p+B.v*((double)Cross(A.v,A.p-B.p)/(double)Cross(A.v,B.v));
}
bool in(Point p,Line l){
	return cmp(Cross(l.v,p-l.p))>=0;
}
int n;
int lcnt=0;
double Area;
double squaring(int n,Point *p){
	double ans=0;
	for(int i=2;i<n;i++){
		ans+=Cross(p[i]-p[1],p[i+1]-p[1]);
//		cout<<a[i].x<<" "<<a[i].y<<endl; 
	}
	return fabs(ans/2);
}
int tot=0;
void Shrink(int &num){
	int siz=0;
	for(int i=1;i<=lcnt;i++){
		if(!siz){
			siz++;
			temp[siz]=l[i];
		}
		else{
			if(!(equal(temp[siz].Ang,l[i].Ang))){
				siz++;
				temp[siz]=l[i];
			}
			else{
				if(in(l[i].p,temp[siz])){
					temp[siz]=l[i];
				}
			}
		}
	}
	num=siz;
	for(int i=1;i<=siz;i++){
		l[i]=temp[i];
	}
}
void Half_Plane_Intersection(){
	sort(l+1,l+1+lcnt,comp);
	Shrink(lcnt);
	int head=1;
	int tail=2;
	convex_hull[1]=l[1];
	convex_hull[2]=l[2];
	for(int i=3;i<=lcnt;i++){
		while(head<tail&&!in(getpot(convex_hull[tail],convex_hull[tail-1]),l[i]))tail--;
		while(head<tail&&!in(getpot(convex_hull[head],convex_hull[head+1]),l[i]))head++;
		tail++;
		convex_hull[tail]=l[i];
	}
	while(head<tail&&!in(getpot(convex_hull[tail],convex_hull[tail-1]),convex_hull[head]))tail--;
	convex_hull[tail+1]=convex_hull[head];
	for(int i=head;i<=tail;i++){
		a[++tot]=getpot(convex_hull[i+1],convex_hull[i]);
	}
//	cout<<tot<<endl;
}
int main(){
	scanf("%d",&n);
	for(int i=1;i<=n;i++){
		int x,y;
		scanf("%d%d",&x,&y);
		p[i].x=x;
		p[i].y=y;
	}
	Area=squaring(n,p);
//	cout<<fixed<<setprecision(4)<<Area<<'\n';
	for(int i=2;i<=n;i++){
		int k1=i;
		int k2=i%n+1;//(1,n) 
		double A=p[1].y-p[2].y-p[k1].y+p[k2].y;
		double B=p[2].x-p[1].x-p[k2].x+p[k1].x;
		double C=(p[1].x*p[2].y-p[2].x*p[1].y-p[k1].x*p[k2].y+p[k2].x*p[k1].y);
		if(!cmp(A)&&!cmp(B))
			continue;
		if(!cmp(A)){//B*Y+C<=0
			if(cmp(B)>0){//Y<=(-(C/B))
				l[++lcnt]=Line(Point(INF,-(C/B)),Vector(-2*INF,0));
			}
			else{//Y>=(-(C/B))
				l[++lcnt]=Line(Point(-INF,-(C/B)),Vector(2*INF,0));
			}
			continue;
		}
		if(!cmp(B)){//A*X+C<=0
			if(cmp(A)>=0){
				l[++lcnt]=Line(Point(-(C/A),-INF),Vector(0,2*INF));
			}
			else{
				l[++lcnt]=Line(Point(-(C/A),INF),Vector(0,-2*INF));
			}
			continue;
		}
		//A*X+B*Y+C<=0
		double k=-(A/B);
		double b=-(C/B);
		//Y?=kx+b
		double x1=0;
		double x2=1;
		double y1=k*x1+b;
		double y2=k*x2+b;
		Point P=Point(x1,y1);
		Point Q=Point(x2,y2);
		if(cmp(-B)>0){
			l[++lcnt]=Line(P,Q-P);//Y>=kX+b
		}
		else{
			l[++lcnt]=Line(Q,P-Q);//Y<=kX+b
		}
	}
	//in the polygon
	for(int i=1;i<=n;i++){
		l[++lcnt]=Line(p[i],p[i%n+1]-p[i]);
	}
	for(int i=1;i<=lcnt;i++){
		l[i].Ang=l[i].ang();
	}
	Half_Plane_Intersection();
	double ans=squaring(tot,a);
	cout<<fixed<<setprecision(4)<<ans/Area;
}