[SPOJ-CIRU]The area of the union of circles/[BZOJ2178]圆的面积并
[SPOJ-CIRU]The area of the union of circles/[BZOJ2178]圆的面积并
题目大意：
求\(n(n\le1000)\)个圆的面积并。
思路：
对于一个\(x\)，我们可以用线段覆盖的方法求出被圆覆盖的长度。用\(f(x)\)表示横坐标为\(x\)时覆盖的长度，则我们可以对\(f(x)\)积分来得到答案。注意面积不连续的部分要分开求。
源代码：
#include<cmath>
#include<cstdio>
#include<cctype>
#include<vector>
#include<algorithm>
inline int getint() {
register char ch;
register bool neg=false;
while(!isdigit(ch=getchar())) neg|=ch=='-';
register int x=ch^'0';
while(isdigit(ch=getchar())) x=(((x<<2)+x)<<1)+(ch^'0');
return neg?-x:x;
}
const int N=1e3+1;
const double eps=1e-4;
inline double sqr(const double &x) {
return x*x;
}
int n;
struct Circle {
double x,y,r;
};
Circle c[N];
struct Segment {
double l,r;
bool operator < (const Segment &rhs) const {
return r<rhs.r;
}
};
Segment t[N];
bool mark[N];
std::vector<std::pair<double,int> > s;
inline double F(const double &x) {
for(register int i=1;i<=n;i++) {
if(fabs(c[i].x-x)>c[i].r) continue;
const double d=sqrt(sqr(c[i].r)-sqr(c[i].x-x));
s.push_back(std::make_pair(c[i].y-d,1));
s.push_back(std::make_pair(c[i].y+d,-1));
}
std::sort(s.begin(),s.end());
int tmp=0;
double st,ans=0;
for(register unsigned i=0;i<s.size();i++) {
if(tmp==0) st=s[i].first;
tmp+=s[i].second;
if(tmp==0) ans+=s[i].first-st;
}
s.clear();
return ans;
}
inline double simpson(const double &a,const double &b) {
const double c=(a+b)/2;
return (F(a)+4*F(c)+F(b))*(b-a)/6;
}
inline double asr(const double &a,const double &b,const double &eps,const double &s) {
const double &c=(a+b)/2,ls=simpson(a,c),rs=simpson(c,b);
if(fabs(ls+rs-s)<=15*eps) return ls+rs+(ls+rs-s)/15;
return asr(a,c,eps/2,ls)+asr(c,b,eps/2,rs);
}
int main() {
n=getint();
for(register int i=1;i<=n;i++) {
c[i].x=getint();
c[i].y=getint();
c[i].r=getint();
}
for(register int i=1;i<=n;i++) {
for(register int j=1;j<=n;j++) {
if(c[i].r>c[j].r||i==j) continue;
if(sqr(c[i].x-c[j].x)+sqr(c[i].y-c[j].y)<sqr(c[i].r-c[j].r)) {
mark[i]=true;
}
}
}
int tot=0;
for(register int i=1;i<=n;i++) {
if(!mark[i]) c[++tot]=c[i];
}
n=tot;
for(register int i=1;i<=n;i++) {
t[i]=(Segment){c[i].x-c[i].r,c[i].x+c[i].r};
}
std::fill(&mark[1],&mark[n]+1,0);
for(register int i=1;i<=n;i++) {
for(register int j=1;j<=n;j++) {
if(i==j) continue;
if(t[j].l<t[i].l&&t[i].r<t[j].r) mark[i]=true;
}
}
tot=0;
for(register int i=1;i<=n;i++) {
if(!mark[i]) t[++tot]=t[i];
}
std::sort(&t[1],&t[tot]+1);
double st=t[1].l,en=t[1].r,ans=0;
for(register int i=2;i<=tot;i++) {
if(t[i].l>en) {
ans+=asr(st,en,eps,simpson(st,en));
st=t[i].l;
}
en=t[i].r;
}
ans+=asr(st,en,eps,simpson(st,en));
printf("%.3f\n",ans);
return 0;
}
但是BZOJ上比较卡精度，而eps开太小又会T，解决方法是将Lyness优化去掉，将没有用的圆去掉。
#include<cmath>
#include<cstdio>
#include<cctype>
#include<vector>
#include<algorithm>
inline int getint() {
register char ch;
register bool neg=false;
while(!isdigit(ch=getchar())) neg|=ch=='-';
register int x=ch^'0';
while(isdigit(ch=getchar())) x=(((x<<2)+x)<<1)+(ch^'0');
return neg?-x:x;
}
const int N=1e3+1;
const double eps=1e-13;
inline double sqr(const double &x) {
return x*x;
}
int n;
struct Circle {
double x,y,r;
};
Circle c[N];
bool mark[N];
std::vector<int> cir;
std::vector<std::pair<double,int> > s;
inline double F(const double &x) {
for(register unsigned j=0;j<cir.size();j++) {
const int &i=cir[j];
if(std::abs(c[i].x-x)>c[i].r) continue;
const double d=sqrt(sqr(c[i].r)-sqr(c[i].x-x));
s.push_back(std::make_pair(c[i].y-d,1));
s.push_back(std::make_pair(c[i].y+d,-1));
}
std::sort(s.begin(),s.end());
double st,ans=0;
for(register unsigned i=0,tmp=0;i<s.size();i++) {
if(tmp==0) st=s[i].first;
tmp+=s[i].second;
if(tmp==0) ans+=s[i].first-st;
}
s.clear();
return ans;
}
inline double simpson(const double &fl,const double &fm,const double &fr,const double &len) {
return (fl+4*fm+fr)*len/6;
}
inline double asr(const double &a,const double &b,const double &fl,const double &fm,const double &fr) {
const double c=(a+b)/2;
const double flm=F((a+c)/2),frm=F((c+b)/2);
const double ls=simpson(fl,flm,fm,c-a),rs=simpson(fm,frm,fr,b-c),s=simpson(fl,fm,fr,b-a);
if(std::abs(ls+rs-s)<eps) return ls+rs;
return asr(a,c,fl,flm,fm)+asr(c,b,fm,frm,fr);
}
struct Node {
double p;
int v,id;
bool operator < (const Node &rhs) const {
return p<rhs.p;
}
};
std::vector<Node> t;
inline double dist(const Circle &a,const Circle &b) {
return sqrt(sqr(a.x-b.x)+sqr(a.y-b.y));
}
int main() {
n=getint();
for(register int i=1;i<=n;i++) {
c[i].x=getint();
c[i].y=getint();
c[i].r=getint();
}
for(register int i=1;i<=n;i++) {
for(register int j=1;j<=n;j++) {
if(i==j) continue;
if(dist(c[i],c[j])<c[j].r-c[i].r) {
mark[i]=true;
break;
}
}
}
int tot=0;
for(register int i=1;i<=n;i++) {
if(!mark[i]) c[++tot]=c[i];
}
n=tot;
for(register int i=1;i<=n;i++) {
t.push_back((Node){c[i].x-c[i].r,1,i});
t.push_back((Node){c[i].x+c[i].r,-1,i});
}
std::sort(t.begin(),t.end());
int beg;
double ans=0;
for(register unsigned i=0,tmp=0;i<t.size();i++) {
if(tmp==0) beg=i;
tmp+=t[i].v;
if(tmp==0) {
for(register int j=beg;j<=(int)i;j++) {
cir.push_back(t[j].id);
}
std::sort(cir.begin(),cir.end());
cir.resize(std::unique(cir.begin(),cir.end())-cir.begin());
const double &st=t[beg].p,&en=t[i].p,mid=(st+en)/2;
ans+=asr(st,en,F(st),F(mid),F(en));
cir.clear();
}
}
printf("%.3f\n",ans);
return 0;
}
