bzoj 2178: 圆的面积并 （辛普森积分）

2178: 圆的面积并

Time Limit: 20 Sec  Memory Limit: 259 MB
Submit: 1679  Solved: 433
[Submit][Status][Discuss]

Description

给出N个圆，求其面积并

Input

先给一个数字N ,N< = 1000 接下来是N行是圆的圆心，半径，其绝对值均为小于1000的整数

Output

面积并，保留三位小数

题解：辛普森积分

辛普森积分可以用来求解一些较平滑曲线的面积，可以自动调整精度但还是有较大的误差。

具体怎么用呢？就是每次找到a,(a+b)/2,b三个位置利用三个位置的f值来近似计算不规则图形或函数的面积。

对于不同的题有不同的f

对这道题来说，f(i)就是x=i直线上所有被覆盖区域的长度。

#include<iostream>#include<cstdio>#include<cstring>#include<algorithm>#include<cmath>#define N 2003#define eps 1e-13using namespace std;int n,m,top,st,ed;bool mark[N];double xl[N],xr[N],ans;struct data{    double x,y,r;}a[N],sk[N]; struct line{    double l,r;}p[N];int cmp(data a,data b){    return a.r<b.r;}int cmp1(data a,data b){    return a.x-a.r<b.x-b.r;}int cmp2(line a,line b){    return a.l<b.l;}double pow(double x){    return x*x;}double dis(data a,data b){    return sqrt(pow(a.x-b.x)+pow(a.y-b.y));}double getf(double x){    double r,dis,len=0;    int sz=0;    for (int i=st;i<=ed;i++){        if (xl[i]>=x||xr[i]<=x) continue;        dis=sqrt(sk[i].r-pow(x-sk[i].x));        p[++sz].l=sk[i].y-dis;  p[sz].r=sk[i].y+dis;    }    sort(p+1,p+sz+1,cmp2);    int i,j;    for (i=1;i<=sz;i++){        r=p[i].r;        for (j=i+1;j<=sz;j++){            if (p[j].l>r) break;            r=max(r,p[j].r);        }        len+=r-p[i].l; i=j-1;    }    return len;}double calc(double l,double fl,double fmid,double fr){    return l/6.0*(fl+4.0*fmid+fr);}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=calc(mid-l,fl,f1,fmid),g2=calc(r-mid,fmid,f2,fr);    if (fabs(g1+g2-s)<eps) return g1+g2;    return simpson(l,m1,mid,fl,f1,fmid,g1)+simpson(mid,m2,r,fmid,f2,fr,g2);}void work(){    for (int i=1;i<=m;i++) xl[i]=sk[i].x-sk[i].r,xr[i]=sk[i].x+sk[i].r,sk[i].r*=sk[i].r;    int i,j; double l,r;    double fl,fr,fmid;    for (i=1;i<=m;i++){        l=xl[i]; r=xr[i];        for (j=i+1;j<=m;j++){            if (xl[j]>r) break;            r=max(r,xr[j]);        }        st=i; ed=j-1; i=j-1;        double mid=(l+r)/2;        fl=getf(l); fmid=getf(mid); fr=getf(r);        //cout<<fl<<" "<<fr<<" "<<fmid<<endl;        ans+=simpson(l,mid,r,fl,fmid,fr,calc(r-l,fl,fmid,fr));    }}int main(){freopen("15.in","r",stdin);    scanf("%d",&n);    for (int i=1;i<=n;i++)     scanf("%lf%lf%lf",&a[i].x,&a[i].y,&a[i].r);    sort(a+1,a+n+1,cmp);    for (int i=1;i<=n-1;i++)     for (int j=i+1;j<=n;j++)      if (dis(a[i],a[j])<=a[j].r-a[i].r) {        mark[i]=1;        break;      }    for (int i=1;i<=n;i++)     if (!mark[i]) sk[++m]=a[i];    sort(sk+1,sk+m+1,cmp1);    work();    printf("%.3lf\n",ans);}