[省选前题目整理][SPOJ CIRU]The area of the union of circles(自适应Simpson积分求圆并面积)

题目链接

http://acm.hust.edu.cn/vjudge/problem/viewProblem.action?id=16754

题目大意

求圆的面积并

思路

http://blog.csdn.net/qpswwww/article/details/44201333

代码

#include <iostream>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <algorithm>
#include <cmath>

#define MAXN 111000
#define EPS 1e-6

using namespace std;

int n;

struct Point
{
    double x,y;
    Point(){}
    Point(double _x,double _y):x(_x),y(_y){}
    bool operator<(const Point &b)const{ return x==b.x?y<b.y:x<b.x; }
}points[MAXN];

struct Circle
{
    Point o;
    double r;
    Circle(){}
    Circle(Point _o,double _r):o(_o),r(_r){}
}circles[MAXN],tmp[MAXN];

bool cmp2(Circle a,Circle b)
{
    return a.r>b.r;
}

bool cmp(Circle a,Circle b) //按照圆的左边缘升序排序,第二关键字右边缘升序排序
{
    if(fabs(a.o.x-a.r-b.o.x+b.r)<=EPS) return a.o.x+a.r<b.o.x+b.r;
    return a.o.x-a.r<b.o.x-b.r; //!!!!
}

pair<double,double>intervals[MAXN];
int tot=0,st,ed;

bool check(Circle a,Circle b) //检查圆b是否被套在圆a里面
{
    return (a.o.x-b.o.x)*(a.o.x-b.o.x)+(a.o.y-b.o.y)*(a.o.y-b.o.y)<=(a.r-b.r)*(a.r-b.r); //!!!!
}

void prework()
{
    sort(circles+1,circles+n+1,cmp2);
    int cnt=0;
    for(int i=1;i<=n;i++)
    {
        bool flag=true;
        for(int j=1;j<=cnt;j++)
            if(check(tmp[j],circles[i]))
            {
                flag=false;
                break;
            }
        if(flag) tmp[++cnt]=circles[i];
    }
    n=cnt;
    for(int i=1;i<=n;i++)
        circles[i]=tmp[i];
}

void getCircleIntersec(Circle a,double x)
{
    intervals[++tot]=make_pair(a.o.y-sqrt(a.r*a.r-(x-a.o.x)*(x-a.o.x)),a.o.y+sqrt(a.r*a.r-(x-a.o.x)*(x-a.o.x)));
}

double f(double x)
{
    tot=0;
    for(int i=st;i<=ed;i++)
        if(x<circles[i].o.x+circles[i].r&&x>circles[i].o.x-circles[i].r) //!!!!!
            getCircleIntersec(circles[i],x);
    sort(intervals+1,intervals+tot+1);
    double ans=0,start=-1e12,end=-1e12;
    for(int i=1;i<=tot;i++)
    {
        if(intervals[i].first>=end)
        {
            ans+=end-start;
            start=intervals[i].first;
            end=intervals[i].second;
        }
        else end=max(end,intervals[i].second);
    }
    ans+=end-start;
    return ans;
}

double calc(double len,double fL,double fM,double fR)
{
    return (fL+4*fM+fR)*len/6;
}

double Simpson(double L,double M,double R,double fL,double fM,double fR,double sqr)
{
    double M1=(L+M)/2,M2=(M+R)/2;
    double fM1=f(M1),fM2=f(M2);
    double g1=calc(M-L,fL,fM1,fM),g2=calc(R-M,fM,fM2,fR);
    if(fabs(sqr-g1-g2)<EPS) return g1+g2;
    return Simpson(L,M1,M,fL,fM1,fM,g1)+Simpson(M,M2,R,fM,fM2,fR,g2);
}

int main()
{
    double ans=0;
    scanf("%d",&n);
    for(int i=1;i<=n;i++)
        scanf("%lf%lf%lf",&circles[i].o.x,&circles[i].o.y,&circles[i].r);
    prework();
    sort(circles+1,circles+n+1,cmp);
    for(int i=1,j;i<=n;i++)
    {
        double L=circles[i].o.x-circles[i].r,R=circles[i].o.x+circles[i].r;
        for(j=i+1;j<=n;j++)
        {
            if(circles[j].o.x-circles[j].r>R) break;
            else R=max(R,circles[j].o.x+circles[j].r);
        }
        double M=(L+R)/2;
        st=i,ed=j-1;
        i=j-1;
        double fL=f(L),fM=f(M),fR=f(R);
        ans+=Simpson(L,M,R,fL,fM,fR,calc(R-L,fL,fM,fR));
    }
    printf("%.3lf\n",ans);
    return 0;
}