AOJ 1177

1100点の幾何。ひたすら場合分けをがんばろう！！

#include<stdio.h>
#include<algorithm>
#include<math.h>
using namespace std;
const double EPS=1e-10;
const double INF=1e+10;
const double PI=acos(-1);
int sig(double r){return (r<-EPS)?-1:(r>+EPS)?+1:0;}
double Abs(double a){return max(a,-a);}
struct Pt{
    double x,y;
    Pt(){}
    Pt(double x,double y):x(x),y(y){}
    Pt operator+(const Pt &a)const{return Pt(x+a.x,y+a.y);}
    Pt operator-(const Pt &a)const{return Pt(x-a.x,y-a.y);}
    Pt operator*(const Pt &a)const{return Pt(x*a.x-y*a.y,x*a.y+y*a.x);}
    Pt operator-()const{return Pt(-x,-y);}
    Pt operator*(const double &k)const{return Pt(x*k,y*k);}
    Pt operator/(const double &k)const{return Pt(x/k,y/k);}
    double Abs2() const{return x*x+y*y;}
    double dot(const Pt &a)const{return x*a.x+y*a.y;}
    double det(const Pt &a)const{return x*a.y-y*a.x;}
    double Abs()const{return sqrt(x*x+y*y);}
};
double tri(const Pt &a,const Pt&b,const Pt&c){return (b-a).det(c-a);}
int iSP(Pt a,Pt b,Pt c){
    int s=sig((b-a).det(c-a));
    if(s)return s;
    if(sig((b-a).dot(c-a))<0)return -2;
    if(sig((a-b).dot(c-b))<0)return +2;
    return 0;
}
int iSS(Pt a,Pt b,Pt c,Pt d){
    return (iSP(a,b,c)*iSP(a,b,d)<=0&&iSP(c,d,a)*(iSP(c,d,b)<=0));
}
Pt pLL(Pt a,Pt b,Pt c,Pt d){
    b=b-a;d=d-c;return a+b*(c-a).det(d)/b.det(d);
}
int iCC(Pt a,double r,Pt b,double s){
    double d=(b-a).Abs();
    if(sig(d)==0&&sig(r-s)==0)return -1;
    if(sig(r-s-d)>0)return +2;
    if(sig(s-r-d)>0)return -2;
    return (sig(r+s-d)>=0)?1:0;
}
pair<Pt,Pt>pCC(Pt a,double r,Pt b,double s){
    double d=(b-a).Abs();
    double x=(d*d+r*r-s*s)/(d*2);
    Pt e=(b-a)/d,w=e*Pt(0,1)*sqrt(max(r*r-x*x,0.0));
    return make_pair(a+e*x-w,a+e*x+w);
}
double dSP(Pt a,Pt b,Pt c){
    if(sig((b-a).dot(c-a))<=0)return (c-a).Abs();
    if(sig((a-b).dot(c-b))<=0)return (c-b).Abs();
    return Abs(tri(a,b,c))/(b-a).Abs();
}
pair<Pt,Pt> pCL(Pt a,double r,Pt b,Pt c){
    Pt h=b+(c-b)*(c-b).dot(a-b)/(c-b).Abs2();
    double d=(h-a).Abs();
    double y=sqrt(max(r*r-d*d,0.0));
    Pt e=(c-b)/(c-b).Abs();
    return make_pair(h-e*y,h+e*y);
}
int use[30];
int main(){
    double r,x1,y1,x2,y2;
    while(1){
		scanf("%lf%lf%lf%lf%lf",&r,&x1,&y1,&x2,&y2);
        if(r<0.5)break;
        if(y2<0){
            double c=y1;
            y1=-y2;
            y2=-c;
        }else if(y1<=0){
            double d=x1;
            double e=x2;
			x1=y1;
			x2=y2;
			
            if(e<0){
                y1=-e;
                y2=-d;
            }else{
                y1=d;
                y2=e;
            }
        }
        if(x2<0){
            double c=x2;
            x2=-x1;
            x1=-c;
        }
        double ret=0;
        double L[5];
        double R[5];
        for(int i=0;i<5;i++)L[i]=R[i]=0;
        L[0]=r;
        R[0]=r;
        int ls=1;
        int rs=1;
        ret+=r*PI;
        if(dSP(Pt(x1,y1),Pt(x2,y1),Pt(0,0))>r){
            use[0]++;
            ret+=r*PI;
            goto end;
        }
        if(x1>0){
            ls++;
            if(Pt(x1,y2).Abs()<r){
                L[ls++]=r-Pt(x1,y2).Abs();
                if(L[2]>x2-x1){
                    L[ls++]=L[2]-(x2-x1);
                    if(L[3]>y2-y1){
                        L[ls++]=L[3]-(y2-y1);
                    }
                }
            }
            if(Pt(x2,y1).Abs()<r){
                R[rs++]=r-Pt(x2,y1).Abs();
                if(R[1]>y2-y1){
                    R[rs++]=R[1]-(y2-y1);
                    if(R[2]>x2-x1){
                        R[rs++]=R[2]-(x2-x1);
                        if(R[3]>y2-y1){
                            R[rs++]=R[3]-(y2-y1);
                        }
                    }
                }
            }
            if(ls>2&&rs>2&&iCC(Pt(x1,y2),L[2],Pt(x2,y2),R[2])==1){
use[1]++;               pair<Pt,Pt> W=pCC(Pt(x1,y2),L[2],Pt(x2,y2),R[2]);
                if(W.first.x<x2+EPS&&atan2(W.first.y,W.first.x)-EPS<atan2(y2,x1)&&W.first.y+EPS>y2){
                    Pt q=W.first;
                    Pt pL=q-Pt(x1,y2);
                    Pt pR=q-Pt(x2,y2);
                    ret+=(PI-atan2(y2,x1))*r;
                    ret+=(atan2(y1,x2))*r;
                    ret+=(PI/2-atan2(y1,x2))*R[1];
                    ret+=(atan2(pR.y,pR.x)-PI/2)*R[2];
                    ret+=(atan2(y2,x1)-atan2(pL.y,pL.x))*L[2];
                    goto end;
                }
                if(W.second.x<x2+EPS&&atan2(W.second.y,W.second.x)-EPS<atan2(y2,x1)&&W.second.y+EPS>y2){
                    Pt q=W.second;
                    Pt pL=q-Pt(x1,y2);
                    Pt pR=q-Pt(x2,y2);
                    ret+=(PI-atan2(y2,x1))*r;
                    ret+=(atan2(y1,x2))*r;
                    ret+=(PI/2-atan2(y1,x2))*R[1];
                    ret+=(atan2(pR.y,pR.x)-PI/2)*R[2];
                    ret+=(atan2(y2,x1)-atan2(pL.y,pL.x))*L[2];
                    goto end;
                }
            }
            if(ls>2&&rs>1&&iCC(Pt(x1,y2),L[2],Pt(x2,y1),R[1])==1){
use[2]++;               pair<Pt,Pt> W=pCC(Pt(x1,y2),L[2],Pt(x2,y1),R[1]);
                if(W.first.x+EPS>=x2&&atan2(W.first.y,W.first.x)-EPS<atan2(y2,x1)&&W.first.y+EPS>=y2){
                    Pt q=W.first;
                    Pt pL=q-Pt(x1,y2);
                    Pt pR=q-Pt(x2,y1);
                   ret+=(PI-atan2(y2,x1))*r;
                    ret+=(atan2(y1,x2))*r;
                    ret+=(atan2(pR.y,pR.x)-atan2(y1,x2))*R[1];
                    ret+=(atan2(y2,x1)-atan2(pL.y,pL.x))*L[2];
                    goto end;
                }
                if(W.second.x+EPS>x2&&atan2(W.second.y,W.second.x)-EPS<atan2(y2,x1)&&W.second.y+EPS>y2){
                    Pt q=W.second;
				//	printf("%f %f\n",q.x,q.y);
                    Pt pL=q-Pt(x1,y2);
                    Pt pR=q-Pt(x2,y1);
			//		printf("%f\n",PI-atan2(y2,x1));
                    ret+=(PI-atan2(y2,x1))*r;
				//	printf("%f\n",atan2(y1,x2));
                    ret+=(atan2(y1,x2))*r;
				//	printf("%f\n",atan2(pR.y,pR.x)-atan2(y1,x2));
                    ret+=(atan2(pR.y,pR.x)-atan2(y1,x2))*R[1];
				//	printf("%f\n",atan2(y2,x1)-atan2(pL.y,pL.x));
                    ret+=(atan2(y2,x1)-atan2(pL.y,pL.x))*L[2];
				//	printf("%f\n",ret);
                    goto end;
                }
            }
            if(ls>3&&rs>1&&iCC(Pt(x2,y2),L[3],Pt(x2,y1),R[1])==1){
use[3]++;               pair<Pt,Pt> W=pCC(Pt(x2,y2),L[3],Pt(x2,y1),R[1]);
                if(W.first.x+EPS>x2&&W.first.y+EPS>y1&&W.first.y<y2+EPS){
                    Pt q=W.first;
                    Pt pL=q-Pt(x2,y2);
                    Pt pR=q-Pt(x2,y1);
                    ret+=(PI-atan2(y2,x1))*r;
                    ret+=(atan2(y1,x2))*r;
                    ret+=(-atan2(y1,x2)+atan2(pR.y,pR.x))*R[1];
                    ret+=(atan2(y2,x1))*L[2];
                    ret+=(-atan2(pL.y,pL.x))*L[3];
                    goto end;
                }
                if(W.second.x+EPS>x2&&W.second.y+EPS>y1&&W.second.y<y2+EPS){
                    Pt q=W.second;
                    Pt pL=q-Pt(x2,y2);
                    Pt pR=q-Pt(x2,y1);
                    ret+=(PI-atan2(y2,x1))*r;
                    ret+=(atan2(y1,x2))*r;
                    ret+=(-atan2(y1,x2)+atan2(pR.y,pR.x))*R[1];
                    ret+=(atan2(y2,x1))*L[2];
                    ret+=(-atan2(pL.y,pL.x))*L[3];
                    goto end;
                }
            }
 
            if(Pt(x1,y2).Abs()>r){
use[4]++;               pair<Pt,Pt> D=pCL(Pt(0,0),r,Pt(x1,y1),Pt(x1,y2));
                if(D.first.y>0){
                    ret+=(PI-atan2(D.first.y,D.first.x))*r;
                }else ret+=(PI-atan2(D.second.y,D.second.x))*r;
            }else{
use[5]++;               ret+=(PI-atan2(y2,x1))*r;
                ret+=(atan2(y2,x1))*L[2];
                ret+=PI/2*L[3];
            }
            if(Pt(x2,y1).Abs()>r){
use[6]++;               pair<Pt,Pt> D=pCL(Pt(0,0),r,Pt(x1,y1),Pt(x2,y1));
                if(D.first.x>D.second.x){
                    ret+=(atan2(D.first.y,D.first.x))*r;
                }else{
                    ret+=(atan2(D.second.y,D.second.x))*r;
                }
            }else{
use[7]++;               ret+=(atan2(y1,x2))*r;
                ret+=(PI/2-atan2(y1,x2))*R[1];
                ret+=PI/2*R[2];
                ret+=PI/2*R[3];
            }
            goto end;
        }
        if(Pt(x1,y1).Abs()<r){
            L[ls++]=r-Pt(x1,y1).Abs();
            if(L[1]>y2-y1){
                L[ls++]=L[1]-(y2-y1);
                if(L[2]>x2-x1){
                    L[ls++]=L[2]-(x2-x1);
                    if(L[3]>y2-y1){
                        L[ls++]=L[3]-(y2-y1);
                    }
                }
            }
        }
        if(Pt(x2,y1).Abs()<r){
            R[rs++]=r-Pt(x2,y1).Abs();
            if(R[1]>y2-y1){
                R[rs++]=R[1]-(y2-y1);
                if(R[2]>x2-x1){
                    R[rs++]=R[2]-(x2-x1);
                    if(R[3]>y2-y1){
                        R[rs++]=R[3]-(y2-y1);
                    }
                }
            }
        }
         
        if(ls>1&&rs>2&&iCC(Pt(x1,y1),L[1],Pt(x2,y2),R[2])==1){
use[8]++;           pair<Pt,Pt> W=pCC(Pt(x1,y1),L[1],Pt(x2,y2),R[2]);
            if(W.first.x<x1+EPS&&W.first.y+EPS>y2){
                Pt q=W.first;
                Pt pL=q-Pt(x1,y1);
                Pt pR=q-Pt(x2,y2);
                ret+=(PI-atan2(y1,x1))*r;
                ret+=(atan2(y1,x2))*r;
                ret+=(atan2(y1,x1)-atan2(pL.y,pL.x))*L[1];
                ret+=(PI/2-atan2(y1,x2))*R[1];
                ret+=(atan2(pR.y,pR.x)-PI/2)*R[2];
                goto end;
            }
            if(W.second.x<x1+EPS&&W.second.y+EPS>y2){
                Pt q=W.second;
                Pt pL=q-Pt(x1,y1);
                Pt pR=q-Pt(x2,y2);
                ret+=(PI-atan2(y1,x1))*r;
                ret+=(atan2(y1,x2))*r;
                ret+=(atan2(y1,x1)-atan2(pL.y,pL.x))*L[1];
                ret+=(PI/2-atan2(y1,x2))*R[1];
                ret+=(atan2(pR.y,pR.x)-PI/2)*R[2];
                goto end;
            }
        }
        if(ls>1&&rs>3&&iCC(Pt(x1,y1),L[1],Pt(x1,y2),R[3])==1){
use[9]++;           pair<Pt,Pt> W=pCC(Pt(x1,y1),L[1],Pt(x1,y2),R[3]);
            if(W.first.x<x1+EPS&&W.first.y+EPS>y1&&W.first.y<y2+EPS){
                Pt q=W.first;
                Pt pL=q-Pt(x1,y1);
                Pt pR=q-Pt(x1,y2);
                ret+=(PI-atan2(y1,x1))*r;
                ret+=(atan2(y1,x2))*r;
                ret+=(atan2(y1,x1)-atan2(pL.y,pL.x))*L[1];
                ret+=(PI/2-atan2(y1,x2))*R[1];
                ret+=PI/2*R[2];
                ret+=(atan2(pR.y,pR.x)-PI)*R[3];
                goto end;
            }
            if(W.second.x<x1+EPS&&W.second.y+EPS>y1&&W.second.y<y2+EPS){
                Pt q=W.second;
                Pt pL=q-Pt(x1,y1);
                Pt pR=q-Pt(x1,y2);
                ret+=(PI-atan2(y1,x1))*r;
                ret+=(atan2(y1,x2))*r;
                ret+=(atan2(y1,x1)-atan2(pL.y,pL.x))*L[1];
                ret+=(PI/2-atan2(y1,x2))*R[1];
                ret+=PI/2*R[2];
                ret+=(atan2(pR.y,pR.x)-PI)*R[3];
                goto end;
            }
        }
        if(ls>2&&rs>1&&iCC(Pt(x1,y2),L[2],Pt(x2,y1),R[1])==1){
use[10]++;          pair<Pt,Pt> W=pCC(Pt(x1,y2),L[2],Pt(x2,y1),R[1]);
            if(W.first.x+EPS>x2&&W.first.y+EPS>y2){
                Pt q=W.first;
                Pt pL=q-Pt(x1,y2);
                Pt pR=q-Pt(x2,y1);
                ret+=(PI-atan2(y1,x1))*r;
                ret+=(atan2(y1,x2))*r;
                ret+=(atan2(y1,x1)-PI/2)*L[1];
                ret+=(atan2(pR.y,pR.x)-atan2(y1,x2))*R[1];
                ret+=(PI/2-atan2(pL.y,pL.x))*L[2];
                goto end;
            }
            if(W.second.x+EPS>x2&&W.second.y+EPS>y2){
                Pt q=W.second;
                Pt pL=q-Pt(x1,y2);
                Pt pR=q-Pt(x2,y1);
                ret+=(PI-atan2(y1,x1))*r;
                ret+=(atan2(y1,x2))*r;
                ret+=(atan2(y1,x1)-PI/2)*L[1];
                ret+=(atan2(pR.y,pR.x)-atan2(y1,x2))*R[1];
                ret+=(PI/2-atan2(pL.y,pL.x))*L[2];
                goto end;
            }
        }
        if(ls>2&&rs>2&&iCC(Pt(x1,y2),L[2],Pt(x2,y2),R[2])==1){
use[11]++;          pair<Pt,Pt> W=pCC(Pt(x1,y2),L[2],Pt(x2,y2),R[2]);
            if(W.first.x<x2+EPS&&W.first.x+EPS>x1&&W.first.y+EPS>y2){
                Pt q=W.first;
                Pt pL=q-Pt(x1,y2);
                Pt pR=q-Pt(x2,y2);
                ret+=(PI-atan2(y1,x1))*r;
                ret+=(atan2(y1,x2))*r;
                ret+=(atan2(y1,x1)-PI/2)*L[1];
                ret+=(PI/2-atan2(y1,x2))*R[1];
                ret+=(atan2(pR.y,pR.x)-PI/2)*R[2];
                ret+=(PI/2-atan2(pL.y,pL.x))*L[2];
                goto end;
            }
            if(W.second.x<x2+EPS&&W.second.x+EPS>x1&&W.second.y+EPS>y2){
                Pt q=W.second;
                Pt pL=q-Pt(x1,y2);
                Pt pR=q-Pt(x2,y2);
                ret+=(PI-atan2(y1,x1))*r;
                ret+=(atan2(y1,x2))*r;
                ret+=(atan2(y1,x1)-PI/2)*L[1];
                ret+=(PI/2-atan2(y1,x2))*R[1];
                ret+=(atan2(pR.y,pR.x)-PI/2)*R[2];
                ret+=(PI/2-atan2(pL.y,pL.x))*L[2];
                goto end;
            }
        }
        if(ls>3&&rs>1&&iCC(Pt(x2,y2),L[3],Pt(x2,y1),R[1])==1){
use[12]++;          pair<Pt,Pt> W=pCC(Pt(x2,y2),L[3],Pt(x2,y1),R[1]);
            if(W.first.x+EPS>x2&&W.first.y+EPS>y1&&W.first.y<y2+EPS){
                Pt q=W.first;
                Pt pL=q-Pt(x2,y2);
                Pt pR=q-Pt(x2,y1);
                ret+=(PI-atan2(y1,x1))*r;
                ret+=(atan2(y1,x2))*r;
                ret+=(-atan2(y1,x2)+atan2(pR.y,pR.x))*R[1];
                ret+=(atan2(y1,x1)-PI/2)*L[1];
                ret+=PI/2*L[2];
                ret+=(-atan2(pL.y,pL.x))*L[3];
                goto end;
            }
            if(W.second.x+EPS>x2&&W.second.y+EPS>y1&&W.second.y<y2+EPS){
                Pt q=W.second;
                Pt pL=q-Pt(x2,y2);
                Pt pR=q-Pt(x2,y1);
                ret+=(PI-atan2(y1,x1))*r;
                ret+=(atan2(y1,x2))*r;
                ret+=(-atan2(y1,x2)+atan2(pR.y,pR.x))*R[1];
                ret+=(atan2(y1,x1)-PI/2)*L[1];
                ret+=PI/2*L[2];
                ret+=(-atan2(pL.y,pL.x))*L[3];
                goto end;
            }
        }
        if(Pt(x1,y1).Abs()>r){
use[13]++;          pair<Pt,Pt> D=pCL(Pt(0,0),r,Pt(x1,y1),Pt(x2,y1));
            if(D.first.x<D.second.x){
                ret+=(PI-atan2(D.first.y,D.first.x))*r;
            }else ret+=(PI-atan2(D.second.y,D.second.x))*r;
        }else{
use[14]++;          ret+=(PI-atan2(y1,x1))*r;
            ret+=(atan2(y1,x1)-PI/2)*L[1];
            ret+=PI/2*L[2];
            ret+=PI/2*L[3];
        }
        if(Pt(x2,y1).Abs()>r){
use[15]++;          pair<Pt,Pt> D=pCL(Pt(0,0),r,Pt(x1,y1),Pt(x2,y1));
            if(D.first.x>D.second.x){
                ret+=(atan2(D.first.y,D.first.x))*r;
            }else{
                ret+=(atan2(D.second.y,D.second.x))*r;
            }
        }else{
use[16]++;          ret+=(atan2(y1,x2))*r;
            ret+=(PI/2-atan2(y1,x2))*R[1];
            ret+=PI/2*R[2];
            ret+=PI/2*R[3];
        }
        end:;
        printf("%.8f\n",ret);
    }
//  for(int i=0;i<17;i++)printf("%d ",use[i]);
}

死因：場合分けのパターン数を減らすために長方形をy軸正と重なるか第一象限におくということをしていたが、そこで典型的なスワップミスをしていた。