原文链接https://www.cnblogs.com/zhouzhendong/p/BZOJ2178.html
题目传送门 - BZOJ2178
题意
给出 $n(n\leq 1000)$ 个圆,求面积并。
所有圆的圆心坐标和半径都是绝对值不大于 1000 的整数。
题解
自适应辛普森积分模板题。注意先删掉被其他圆包含的圆。
但是 bzoj 大概是加过数据了,网上大部分直接自适应辛普森的代码都 TLE 了。
有一种卡常方法效果很好:
把 x 坐标按照整点划分成 $O(1000)$ 个区间,对于每一个区间,删除不涉及这个区间的所有圆。然后就过了。
代码
#include <bits/stdc++.h>
#define y1 __zzd12323
using namespace std;
const int N=1005;
const double Eps=1e-9,INF=1e9;
double sqr(double x){
return x*x;
}
int n;
struct Circle{
double x,y,r;
}c[N],_c[N];
bool cmpr(Circle a,Circle b){
return a.r>b.r;
}
struct Seg{
double L,R;
}s[N];
bool cmp(Seg a,Seg b){
return a.L<b.L;
}
double F(double x){
int m=0;
for (int i=1;i<=n;i++){
if (fabs(x-c[i].x)>c[i].r-Eps)
continue;
double M=c[i].y,a=sqrt(sqr(c[i].r)-sqr(x-c[i].x));
s[++m].L=M-a,s[m].R=M+a;
}
sort(s+1,s+m+1,cmp);
double res=0.0,mx=-INF;
for (int i=1;i<=m;i++)
if (s[i].L>mx)
res+=s[i].R-s[i].L,mx=s[i].R;
else if (s[i].R>mx)
res+=s[i].R-mx,mx=s[i].R;
return res;
}
double Simpson(double a,double b,double Fa,double Fb,double Fab){
return (b-a)/6*(Fa+Fab*4+Fb);
}
double Asr(double a,double b,double v,double Fa,double Fb,double Fab){
double ab=(a+b)/2;
double LFab=F((a+ab)/2),L=Simpson(a,ab,Fa,Fab,LFab);
double RFab=F((ab+b)/2),R=Simpson(ab,b,Fab,Fb,RFab);
if (fabs(v-L-R)<Eps)
return v;
return Asr(a,ab,L,Fa,Fab,LFab)
+Asr(ab,b,R,Fab,Fb,RFab);
}
int killed[N];
double Dis(double x1,double y1,double x2,double y2){
return sqrt(sqr(x1-x2)+sqr(y1-y2));
}
int main(){
scanf("%d",&n);
double a=INF,b=-INF;
for (int i=1;i<=n;i++){
scanf("%lf%lf%lf",&c[i].x,&c[i].y,&c[i].r);
a=min(a,c[i].x-c[i].r);
b=max(b,c[i].x+c[i].r);
}
sort(c+1,c+n+1,cmpr);
for (int i=1;i<=n;i++)
for (int j=i+1;j<=n;j++)
if (Dis(c[i].x,c[i].y,c[j].x,c[j].y)+c[j].r<c[i].r+Eps)
killed[j]=1;
int _n=0;
for (int i=1;i<=n;i++)
if (!killed[i])
c[++_n]=c[i];
n=_n;
double ans=0;
for (double x=a;x<b;x++){
for (int i=1;i<=n;i++)
_c[i]=c[i];
int _n=n;
n=0;
for (int i=1;i<=_n;i++)
if (c[i].x-c[i].r<i+1&&c[i].x+c[i].r>i)
c[++n]=_c[i];
n=_n;
for (int i=1;i<=n;i++)
c[i]=_c[i];
double Fa=F(x),Fb=F(x+1),Fab=F((x+x+1)/2);
ans+=Asr(x,x+1,Simpson(x,x+1,Fa,Fb,Fab),Fa,Fb,Fab);
}
printf("%.3lf",ans);
return 0;
}