算法笔记--辛普森积分公式

一图胜千言。

自适应辛普森模板：

 

namespace SimpsonIntegral {
    inline double simpson(double l, double r, function<double(double)> f) {
        return (f(l)+4.0*f((l+r)/2.0)+f(r))/6.0*(r-l);
    }
    inline double integral(double l, double r, double eps, function<double(double)> f) {
        double mid = (l+r)/2;
        double res = simpson(l, r, f), fl = simpson(l, mid, f), fr = simpson(mid, r, f);
        if(fabs(res - fl - fr) < eps) return fl + fr;
        else return integral(l, mid, eps/2, f) + integral(mid, r, eps/2, f);
    }
}

 

例题1：HDU 1724 Ellipse

代码：

 

#include<bits/stdc++.h>
using namespace std;
#define fi first
#define se second
#define LL long long
#define pli pair<LL, int>
#define fio ios::sync_with_stdio(false); cin.tie(0); cout.tie(0);

namespace SimpsonIntegral {
    inline double simpson(double l, double r, function<double(double)> f) {
        return (f(l)+4.0*f((l+r)/2.0)+f(r))/6.0*(r-l);
    }
    inline double integral(double l, double r, double eps, function<double(double)> f) {
        double mid = (l+r)/2;
        double res = simpson(l, r, f), fl = simpson(l, mid, f), fr = simpson(mid, r, f);
        if(fabs(res - fl - fr) < eps) return fl + fr;
        else return integral(l, mid, eps/2, f) + integral(mid, r, eps/2, f);
    }
}

int T;
double a, b, l, r;
int main() {
    fio;
    cin >> T;
    auto f = [&](double x){ return b*sqrt(1.0-x*x/a/a);};
    while(T--) {
        cin >> a >> b >> l >> r;
        cout << fixed << setprecision(3) << SimpsonIntegral::integral(l, r, 1e-5, f)*2 << "\n";
    }
    return 0;
}