题意
给n个圆,求这些圆的面积的并
输入格式和上次用几何法求的格式一样
解法:自适应simpson积分
首先还是去掉被完全覆盖的圆,相当于计算所有圆在和某条x=t的直线相交的长度,然后把所有t拿出来求和.这个可以用simpson(比几何法好写很多)
#include<bits/stdc++.h>
using namespace std;
const int maxn=2e3+5;
const double eps=1e-6;
inline int read(){
char c=getchar();int t=0,f=1;
while(!isdigit(c)){if(c=='-')f=-1;c=getchar();}
while(isdigit(c)){t=(t<<3)+(t<<1)+(c^48);c=getchar();}
return t*f;
}
double rs,rb;
int t,ns,nb,n;
struct node{
double x,y;
double r;
}a[maxn];
struct line{
double l,r;
}p[maxn];
double xl[maxn],xr[maxn],ans;
inline double dist(int x,int y){
return sqrt((a[x].x-a[y].x)*(a[x].x-a[y].x)+(a[x].y-a[y].y)*(a[x].y-a[y].y));
}
int del[maxn],st,ed;
double l,r;
bool cmp2(node a,node b){
return a.x-a.r<b.x-b.r;
}
bool cmp3(line a,line b){
return a.l<b.l;
}
inline double getf(double x){
int cnt=0;double r,len=0,dis;
for(int i=st;i<=ed;i++){
if(x<=xl[i]||x>=xr[i])continue;
dis=sqrt(a[i].r*a[i].r-(x-a[i].x)*(x-a[i].x));
p[++cnt].l=a[i].y-dis;p[cnt].r=a[i].y+dis;
}
sort(p+1,p+1+cnt,cmp3);
for(int i=1;i<=cnt;i++){
r=p[i].r;int j;
for(j=i+1;j<=cnt;j++){
if(p[j].l>r)break;
if(r<p[j].r)r=p[j].r;
}
len+=r-p[i].l;i=j-1;
}
return len;
}
double val(double fl,double fmid,double fr,double len){
return (fl+4*fmid+fr)*len/6;
}
double simpson(double l,double mid,double r,double fl,double fmid,double fr,double s){
double m1=(l+mid)/2,m2=(mid+r)/2;
double f1=getf(m1),f2=getf(m2);
double g1=val(fl,f1,fmid,mid-l),g2=val(fmid,f2,fr,r-mid);
if(fabs(g1+g2-s)<eps)return g1+g2;
else return simpson(l,m1,mid,fl,f1,fmid,g1)+simpson(mid,m2,r,fmid,f2,fr,g2);
}
int main(){
//freopen("simpson.in","r",stdin);
//freopen("simpson.out","w",stdout);
scanf("%lf%lf",&rs,&rb);
t=read();
while(t--){n=0;ans=0;
memset(del,0,sizeof(del));
ns=read(),nb=read();
for(int i=1;i<=ns;i++){n++;
a[n].x=read(),a[n].y=read();
a[n].r=rs;
}
for(int i=1;i<=nb;i++){n++;
a[n].x=read(),a[n].y=read();
a[n].r=rb;
}
for(int i=1;i<=n;i++){
for(int j=i+1;j<=n;j++){
if(dist(i,j)<=a[j].r-a[i].r){
del[i]=1;break;
}
}
}
int tmp=0;
for(int i=1;i<=n;i++){
if(!del[i]){tmp++;
a[tmp]=a[i];
}
}
n=tmp;
/*for(int i=1;i<=n;i++){
printf("%.3lf %.3lf %.3lf\n",a[i].x,a[i].y,a[i].r);
}*/
sort(a+1,a+1+n,cmp2);
for(int i=1;i<=n;i++){
xl[i]=a[i].x-a[i].r;xr[i]=a[i].x+a[i].r;
}
for(int i=1;i<=n;i++){
double l=xl[i];r=xr[i];int j;
for(j=i+1;j<=n;j++){
if(xl[j]>r)break;
if(xr[j]>r)r=xr[j];
}
st=i,ed=j-1;i=j-1;
double mid=(l+r)/2;
double fl=getf(l),fmid=getf(mid),fr=getf(r);
ans+=simpson(l,mid,r,fl,fmid,fr,val(fl,fmid,fr,r-l));
}
printf("%.5lf\n",ans);
}
return 0;
}
/*
0.58 1.31
1
2 2
2 2
4 4
1 3
3 3
*/