SPOJ CIRU The area of the union of circles

You are given N circles and expected to calculate the area of the union of the circles !

Input

The first line is one integer n indicates the number of the circles. (1 <= n <= 1000)

Then follows n lines every line has three integers

Xi Yi Ri

indicates the coordinate of the center of the circle, and the radius. (|Xi|. |Yi| <= 1000, Ri <= 1000)

Note that in this problem Ri may be 0 and it just means one point !

Output

The total area that these N circles with 3 digits after decimal point

Example

Input:
3

0 0 1
0 0 1
100 100 1


Output:
6.283

simpson自适应积分法

精度只需要1e-6，十分友好

调试语句懒得删

  1 #include<iostream>
  2 #include<cstdio>
  3 #include<algorithm>
  4 #include<cstring>
  5 #include<cmath>
  6 using namespace std;
  7 const double eps=1e-6;
  8 const int INF=1e9;
  9 const int mxn=1010;
 10 int read(){
 11     int x=0,f=1;char ch=getchar();
 12     while(ch<'0' || ch>'9'){if(ch=='-')f=-1;ch=getchar();}
 13     while(ch>='0' && ch<='9'){x=x*10-'0'+ch;ch=getchar();}
 14     return x*f;
 15 }
 16 //
 17 struct cir{
 18     double x,y,r;
 19     friend bool operator < (const cir a,const cir b){return a.r<b.r;}
 20 }c[mxn];int cnt=0;
 21 inline double dist(cir a,cir b){return sqrt((a.x-b.x)*(a.x-b.x)+(a.y-b.y)*(a.y-b.y));}
 22 //圆
 23 struct line{
 24     double l,r;
 25     friend bool operator <(const line a,const line b){return a.l<b.l;}
 26 }a[mxn],b[mxn];int lct=0;
 27 double f(double x){
 28     int i,j;
 29     lct=0;
 30     for(i=1;i<=cnt;i++){//计算直线截得圆弧长度 
 31         if(fabs(c[i].x-x)>=c[i].r)continue;
 32         double h= sqrt(c[i].r*c[i].r-(c[i].x-x)*(c[i].x-x));
 33         a[++lct].l=c[i].y-h;
 34         a[lct].r=c[i].y+h;
 35     }
 36     if(!lct)return 0;
 37     double len=0,last=-INF;
 38     sort(a+1,a+lct+1);
 39     for(i=1;i<=lct;i++){//线段长度并 
 40         if(a[i].l>last){len+=a[i].r-a[i].l;last=a[i].r;}
 41         else if(a[i].r>last){len+=a[i].r-last;last=a[i].r;}
 42     }
 43 //    printf("x:%.3f  len:%.3f
",x,len);
 44     return len;
 45 }
 46 inline double sim(double l,double r){
 47     return (f(l)+4*f((l+r)/2)+f(r))*(r-l)/6;
 48 }
 49 double solve(double l,double r,double S){
 50     double mid=(l+r)/2;
 51     double ls=sim(l,mid);
 52     double rs=sim(mid,r);
 53     if(fabs(rs+ls-S)<eps)return ls+rs;
 54     return solve(l,mid,ls)+solve(mid,r,rs);
 55 }
 56 int n;
 57 double ans=0;
 58 bool del[mxn];
 59 int main(){
 60     n=read();
 61     int i,j;
 62     double L=INF,R=-INF;
 63     for(i=1;i<=n;i++){
 64         c[i].x=read();    c[i].y=read();    c[i].r=read();
 65 //        L=min(L,c[i].x-c[i].r);
 66 //        R=max(R,c[i].x+c[i].r);
 67     }
 68     //
 69     sort(c+1,c+n+1);
 70     for(i=1;i<n;i++)
 71         for(j=i+1;j<=n;j++){
 72 //            printf("%d %.3f %.3f %.3f %.3f
",j,c[j].x,c[i].r,c[j].r,dist(c[i],c[j]));
 73             if(c[j].r-c[i].r>=dist(c[i],c[j]))
 74                 {del[i]=1;break;}
 75         }
 76     for(i=1;i<=n;i++)
 77         if(!del[i])c[++cnt]=c[i];
 78     //删去被包含的圆
 79 //    printf("cnt:%d
",cnt);
 80     double tmp=-INF;int blct=0;
 81     for(i=1;i<=cnt;i++){
 82         b[++blct].l=c[i].x-c[i].r;
 83         b[blct].r=c[i].x+c[i].r;
 84     }
 85     sort(b+1,b+blct+1);
 86 //    printf("lct:%d
",blct);
 87 //    int tlct=t;
 88     for(i=1;i<=blct;i++){
 89 //        printf("%.3f %.3f
",b[i].l,b[i].r);
 90 //        printf("tmp:%.3f
",tmp);
 91         if(b[i].r<=tmp)continue;
 92         L=max(tmp,b[i].l);
 93 //        printf("%d: %.3f %.3f
",i,L,a[i].r);
 94         ans+=solve(L,b[i].r,sim(L,b[i].r));
 95 //        printf("ANS:%.3f
",ans);
 96 //        printf("nlct:%d
",lct);
 97         tmp=b[i].r;
 98     }
 99     
100     
101 //    ans=solve(L,R,f((L+R)/2));
102     printf("%.3f
",ans);
103     return 0;
104 }