HDU - 3982:Harry Potter and J.K.Rowling(半平面交+圆与多边形求交)(WA ing)
pro:给定一枚蛋糕,蛋糕上某个位置有个草莓,寿星在上面切了N刀,最后寿星会吃含有草莓的那一块蛋糕,问他的蛋糕占总蛋糕的面积比。
sol:显然需要半平面交求含有蛋糕的那一块,然后有圆弧,不太方便求交。 所以我们可以直线构成的边界,求出平面交; 然后用这个多边形去和圆求交。
(百度了一下很多人都没过,好像是这题很卡精度,反正我每个地方都改过,还是WA,大概wa了4个小时了,要不以后再回来改。
当然也不排除有其他问题。
#include<bits/stdc++.h> #define rep(i,a,b) for(int i=a;i<=b;i++) using namespace std; const int maxn=200010; const double eps=1e-12; const double pi=acos(-1.0); struct point{ double x,y; point(){} point(double xx,double yy):x(xx),y(yy){} }; struct Circle{ point c; double r; }; struct line{ point a,b;//起点 point p;//起点到终点的向量 double angle; }; double det(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;} point operator *(point a,double t){ return point(a.x*t,a.y*t);} point operator +(point a,point b){ return point(a.x+b.x,a.y+b.y);} point operator -(point a,point b){ return point(a.x-b.x,a.y-b.y);} double Length(point A){return sqrt(dot(A,A));} double getangle(point a){ return atan2(a.y,a.x);} double getangle(line a){ return getangle(a.p);} int dcmp(double x){ if(fabs(x)<eps) return 0; if(x<0) return -1; return 1; } double TriAngleCircleInsection(Circle C, point A, point B) { point OA=A-C.c,OB=B-C.c; point BA=A-B, BC=C.c-B; point AB=B-A, AC=C.c-A; double DOA=Length(OA),DOB=Length(OB),DAB=Length(AB),r=C.r; if(dcmp(det(OA,OB))==0) return 0; //,三点一线,不构成三角形 if(dcmp(DOA-C.r)<0&&dcmp(DOB-C.r)<0) return det(OA,OB)*0.5; //内部 else if(DOB<r&&DOA>=r) //一内一外 { double x=(dot(BA,BC)+sqrt(r*r*DAB*DAB-det(BA,BC)*det(BA,BC)))/DAB; double TS=det(OA,OB)*0.5; return asin(TS*(1-x/DAB)*2/r/DOA)*r*r*0.5+TS*x/DAB; } else if(DOB>=r&&DOA<r)// 一外一内 { double y=(dot(AB,AC)+sqrt(r*r*DAB*DAB-det(AB,AC)*det(AB,AC)))/DAB; double TS=det(OA,OB)*0.5; return asin(TS*(1-y/DAB)*2/r/DOB)*r*r*0.5+TS*y/DAB; } else if(fabs(det(OA,OB))>=r*DAB||dot(AB,AC)<=0||dot(BA,BC)<=0)// { if(dot(OA,OB)<0){ if(det(OA,OB)<0) return (-acos(-1.0)-asin(det(OA,OB)/DOA/DOB))*r*r*0.5; else return ( acos(-1.0)-asin(det(OA,OB)/DOA/DOB))*r*r*0.5; } else return asin(det(OA,OB)/DOA/DOB)*r*r*0.5; //小于90度,以为asin对应的区间是[-90度,90度] } else //弧+三角形 { double x=(dot(BA,BC)+sqrt(r*r*DAB*DAB-det(BA,BC)*det(BA,BC)))/DAB; double y=(dot(AB,AC)+sqrt(r*r*DAB*DAB-det(AB,AC)*det(AB,AC)))/DAB; double TS=det(OA,OB)*0.5; return (asin(TS*(1-x/DAB)*2/r/DOA)+asin(TS*(1-y/DAB)*2/r/DOB))*r*r*0.5 + TS*((x+y)/DAB-1); } } point llintersect(line A,line B) { point C=A.a-B.a; double t=det(C,B.p)/det(B.p,A.p); return A.a+A.p*t; } point s[maxn]; line t[maxn],q[maxn]; bool cmp(line a,line b){ double A=getangle(a),B=getangle(b); point t=(b.a+b.p)-a.a; if(fabs(A-B)<eps) return det(a.p,t)>=0.0; return A<B; } bool onright(line P,line a,line b) { point o=llintersect(a,b); point Q=o-P.a; return det(Q,P.p)>0; //如果同一直线上不能相互看到,则>=0 } int tail,head; void halfplaneintersect(int N) { sort(t+1,t+N+1,cmp); int tot=0; rep(i,1,N-1) { if(fabs(getangle(t[i])-getangle(t[i+1]))>eps) t[++tot]=t[i]; } t[++tot]=t[N]; head=tail=0; rep(i,1,tot){ while(tail>head+1&&onright(t[i],q[tail],q[tail-1])) tail--; while(tail>head+1&&onright(t[i],q[head+1],q[head+2])) head++; q[++tail]=t[i]; } while(tail>head+1&&onright(t[head+1],q[tail],q[tail-1])) tail--; } point a[maxn],A,B; int main() { int N,T,Ca=0; Circle C; double x1,y1,x2,y2; scanf("%d",&T); while(T--){ scanf("%lf%d",&C.r,&N); C.c.x=0, C.c.y=0; rep(i,1,N) scanf("%lf%lf%lf%lf",&t[i].a.x,&t[i].a.y,&t[i].b.x,&t[i].b.y); N++; t[N].a=point(-10000,-10000); t[N].b=point(10000,-10000); N++; t[N].a=point(10000,-10000); t[N].b=point(10000,10000); N++; t[N].a=point(10000,10000); t[N].b=point(-10000,10000); N++; t[N].a=point(-10000,10000); t[N].b=point(-10000,-10000); point cr; scanf("%lf%lf",&cr.x,&cr.y); rep(i,1,N){ //保证草莓在刀的左边。 point p=t[i].b-t[i].a; point k=cr-t[i].a; if(det(p,k)<0) swap(t[i].a,t[i].b); t[i].p=t[i].b-t[i].a; } halfplaneintersect(N); double ans=0,sum=pi*C.r*C.r; q[tail+1]=q[head+1]; q[tail+2]=q[head+2]; rep(i,head+1,tail){ //cout<<q[i].a.x<<" "<<q[i].a.y<<" "<<q[i].a.x+q[i].p.x<<" "<<q[i].a.y+q[i].p.y<<endl; ans+=TriAngleCircleInsection (C,llintersect(q[i],q[i+1]),llintersect(q[i+1],q[i+2])); } printf("Case %d: %.5lf%%\n",++Ca,ans*100/sum); } return 0; }
It is your time to fight!