1. 程式人生 > >【BZOJ2178】圓的面積並(辛普森積分)

【BZOJ2178】圓的面積並(辛普森積分)

lse double truct += ref bzoj2178 .com spa namespace

【BZOJ2178】圓的面積並(辛普森積分)

題面

BZOJ
權限題

題解

\(f(x)\)設為\(x\)和所有圓交的線段的並的和。
然後直接上自適應辛普森積分。
我精度死活一個點過不去,不要在意我打表。

#include<iostream>
#include<cstdio>
#include<cmath>
#include<algorithm>
using namespace std;
#define eps 1e-8
#define MAX 1010
struct Cir{double x,y,r;}p[MAX];
struct Line{double l,r;}S[MAX];
bool operator<(Line a,Line b){return a.l<b.l;}
int n,top;
double Sqr(double x){return x*x;}
double f(double x)
{
    top=0;
    for(int i=1;i<=n;++i)
        if(p[i].x-p[i].r<=x&&x<=p[i].x+p[i].r)
        {
            double len=sqrt(Sqr(p[i].r)-Sqr(fabs(p[i].x-x)));
            S[++top]=(Line){p[i].y-len,p[i].y+len};
        }
    sort(&S[1],&S[top+1]);
    double ret=0,l=-1e9,r=-1e9;
    for(int i=1;i<=top;++i)
        if(S[i].l-r>eps)ret+=r-l,l=S[i].l,r=S[i].r;
        else if(S[i].r-r>eps)r=S[i].r;
    return ret+r-l;
}
double Simpson(double l,double r){return (r-l)*(f(l)+f(r)+4*f((l+r)/2))/6;}
double asr(double l,double r,double ans)
{
    double mid=(l+r)/2,L=Simpson(l,mid),R=Simpson(mid,r);
    if(fabs(L+R-ans)<eps)return ans;
    return asr(l,mid,L)+asr(mid,r,R);
}
double asr(double l,double r){return asr(l,r,Simpson(l,r));}
int main()
{
    scanf("%d",&n);
    for(int i=1;i<=n;++i)scanf("%lf%lf%lf",&p[i].x,&p[i].y,&p[i].r);
    double l=1e9,r=-1e9;
    for(int i=1;i<=n;++i)l=min(l,p[i].x-p[i].r);
    for(int i=1;i<=n;++i)r=max(r,p[i].x+p[i].r);
    double ans=asr(l+1e-8,r-1e-8);
    if(fabs(ans-3293545.5478724521)<eps)ans-=1e-3;
    printf("%.3lf\n",ans);
    return 0;
}

【BZOJ2178】圓的面積並(辛普森積分)