辛普森积分

辛普森积分

---

定义

利用积分的定义对不规则面积求和，∫baf(x)dx代表从a到b的面积和，我们利用自适应的方式对面积进行二分，如果左边右边的面积和等于总的面积和那么就自适应((r-l)(f(l)+4f(mid)+f(r))/6)成功。

#include<iostream>
#include<cstring>
#include<cmath>
#include<algorithm>
using namespace std;
const double eps = 1e-12;
double f(double x)
{
    return sin(x)/x;
}
double simpson(double l,double r)
{
    double mid=(l+r)/2;
    return (r-l)*(f(l)+4*f(mid)+f(r))/6;
}
double asr(double l,double r,double s)
{
    double mid=(l+r)/2;
    double left=simpson(l,mid),right=simpson(mid,r);
    if(fabs(left+right-s)<eps)  return left+right;
    return  asr(l,mid,left)+asr(mid,r,right);
}
int main()
{
    double l,r;
    scanf("%lf%lf",&l,&r);
    printf("%.6lf\n",asr(l,r,simpson(l,r)));
    return 0;
}

---

圆面积并

判断x0位置处与所有圆的交集进行区间和并，然后利用自适应辛普森积分。

#include<iostream>
#include<cstring>
#include<algorithm>
#include<cmath>
#define x first
#define y second
using namespace std;
typedef pair<double,double> PDD;
const int N = 1e3+10;
const double eps = 1e-8;
int n;
struct Circle
{
    PDD r;
    double R;
}c[N];
PDD q[N];
int dcmp(double a,double b)
{
    if(fabs(a-b)<eps)   return 0;
    return a<b?-1:1;
}
double f(double x)
{
    int cnt=0;
    for(int i=0;i<n;i++)
    {
        double X=fabs(x-c[i].r.x),R=c[i].R;
        if(dcmp(X,R)<0)
        {
            double Y=sqrt(R*R-X*X);
            q[cnt++]={c[i].r.y-Y,c[i].r.y+Y};
        }
    }
    if(!cnt)    return 0;
    sort(q,q+cnt);
    double  res=0,st=q[0].x,ed=q[0].y;
    for(int i=1;i<cnt;i++)
        if(q[i].x<=ed)   ed=max(ed,q[i].y);
        else{
            res+=ed-st;
            st=q[i].x,ed=q[i].y;
        }
    return res+ed-st;
}
double simpson(double l,double r)
{
    double mid=(l+r)/2;
    return  (r-l)*(f(l)+4*f(mid)+f(r))/6;
}
double asr(double l,double r,double s)
{
    double mid=(l+r)/2;
    double left=simpson(l,mid),right=simpson(mid,r);
    if(fabs(s-left-right)<eps)  return s;
    return asr(l,mid,left)+asr(mid,r,right);
}
int main()
{
    scanf("%d",&n);
    double l=2000,r=-2000;
    for(int i=0;i<n;i++)
    {
        scanf("%lf%lf%lf",&c[i].r.x,&c[i].r.y,&c[i].R);
        l=min(l,c[i].r.x-c[i].R),r=max(r,c[i].r.x+c[i].R);
    }
    printf("%.3lf\n",asr(l-100,r+100,simpson(l,r)));
    return 0;
}

---