CONTENT: 辛普森积分
DETAIL: 辛普森积分



BB

很早前在华东师范月赛看过，一个小点所以顺便补一下，屯模板用

Introduction

【模板】自适应辛普森法1 - Luogu P4525

Description
对函数f(x)=(cx+d)/(ax+b)对区间[l,r]积分
Sample Input

1 2 3 4 5 6

Sample Output

2.732937

Solution

#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N = (int)10000 * 220 + 15;

double a, b, c, d;

inline double f(double x) {
    return (c * x + d) / (a * x + b);
}

inline double simpson(double l, double r) {
    double mid = (l + r) / 2;
    return (r - l) * (f(l) + f(r) + 4 * f(mid)) / 6;
}

inline double calc(double l, double r, double eps, double ans) {
    double mid = (l + r) / 2;
    double lres = simpson(l, mid), rres = simpson(mid, r);
    if(fabs(lres + rres - ans) <= 15 * eps) {
        return lres + rres + (lres + rres - ans) / 15;
    }
    return calc(l, mid, eps / 2, lres) + calc(mid, r, eps / 2, rres);
}

int main() {
    double l, r;
    while(~scanf("%lf%lf%lf%lf%lf%lf", &a, &b, &c, &d, &l, &r)) {
        printf("%.6f\n", calc(l, r, 1e-6, simpson(l, r)));
    }
}


【模板】自适应辛普森法2 - Luogu P4526

Description
对函数f(x)=x^(a/x-x)对区间[0,inf)积分
Sample Input

2.33

Sample Output

1.51068

Solution
对a打表

#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N = (int)10000 * 220 + 15;

double a;

inline double f(double x) {
    return pow(x, a / x - x);
}

inline double simpson(double l, double r) {
    double mid = (l + r) / 2;
    return (r - l) * (f(l) + f(r) + 4 * f(mid)) / 6;
}

inline double calc(double l, double r, double eps, double ans) {
    double mid = (l + r) / 2;
    double lres = simpson(l, mid), rres = simpson(mid, r);
    if(fabs(lres + rres - ans) <= 15 * eps) {
        return lres + rres + (lres + rres - ans) / 15;
    }
    return calc(l, mid, eps / 2, lres) + calc(mid, r, eps / 2, rres);
}

int main() {
    while(~scanf("%lf", &a)) {
        if(a < 0) {
            puts("orz");
            continue;
        }
        double l = 1e-10, r = 16;
        printf("%.5f\n", calc(l, r, 1e-7, simpson(l, r)));
    }
}