BZOJ 2178 圆的面积并

#include<cstdio>
#include<cmath>
#include<algorithm>
#define N 1005
using namespace std;
struct circle
{
double r, x, y;
}c
;
struct interval
{
double l, r;
friend bool operator < (interval a, interval b)
{
if(a.l!=b.l)return a.l<b.l;
else return a.r>b.r;
}
}inter
;
const double INF = 1000000000.0, eps = 1e-13;
int n;
bool ban
;
bool cover(circle a, circle b)
{
double dis = sqrt((a.x-b.x)*(a.x-b.x)+(a.y-b.y)*(a.y-b.y));
return dis + b.r <= a.r;
}
bool cmp(circle a, circle b)
{
if(a.x-a.r != b.x-b.r)return a.x-a.r < b.x-b.r;
else return a.x+a.r > b.x+b.r;
}
double f(double x)
{
int icnt=0;
double d, h, ret=0, l = -INF;
for(int i = 1; i <= n; i++)
{
d = abs(c[i].x-x);
if(d>=c[i].r)continue;
h = sqrt(c[i].r*c[i].r - d*d);
inter[++icnt]=(interval){c[i].y-h, c[i].y+h};
}
if(!icnt)return 0;
sort(inter+1,inter+1+icnt);
for(int i = 1; i <= icnt; i++)
{
if(inter[i].r <= l)continue;
l = max(l,inter[i].l);
ret += inter[i].r-l;
l = inter[i].r;
}
return ret;
}
double Simpson(double l, double r, double fl, double fm, double fr)
{
double mid = (l+r)/2, lmid = (l+mid)/2, rmid = (mid+r)/2;
double flm = f(lmid), frm = f(rmid);
if(abs ((fl+4*flm+fm)/6*(mid-l) + (fm+4*frm+fr)/6*(r-mid) - (fl+4*fm+fr)/6*(r-l) ) < eps)
return (fl+4*fm+fr)/6*(r-l);
else return Simpson(l,mid,fl,flm,fm) + Simpson(mid,r,fm,frm,fr);
}
int main()
{
scanf("%d",&n);
for(int i = 1; i <= n; i++)
scanf("%lf%lf%lf",&c[i].x,&c[i].y,&c[i].r);
for(int i = 1; i <= n; i++)
if(!ban[i])
for(int j = i+1; j <= n; j++)
if(!ban[j])
if(cover(c[i],c[j]))
ban[j]=1;
int pos = 0;
for(int i = 1; i <= n; i++)
if(!ban[i])
c[++pos]=c[i];
n=pos;
sort(c+1,c+1+n,cmp);
double l=-INF, ans=0;
for(int i = 1; i <= n; i++)
{
if(c[i].x + c[i].r <= l+eps)continue;
l = max(l,c[i].x-c[i].r);
ans += Simpson(l,c[i].x+c[i].r,f(l),f((l+c[i].x+c[i].r)/2),f(c[i].x+c[i].r));
l = c[i].x + c[i].r;
}
printf("%.3lf\n",ans);
}