自适应性辛普森法通过取区间左右端点和区间中点三点的函数值拟合二次函数的方式近似计算一段区间内的积分值，若与之前计算的值差别较大则将区间分为两份后递归计算。
Simpson公式：$$\int_a^b f(x)dx \approx \frac{\left(b-a\right)\left[f(a)+f(b)+4f(\frac{a+b}{2})\right]}{6}$$
推导过程：设\(f(x)\)为被积函数，\(g(x)=Ax^2+Bx+C\)为拟合的二次函数，则有：
$$\int_a^b f(x)dx \approx \int_a^b Ax^2+Bx+C$$
$$=\frac{A}{3}(b^3-a^3)+\frac{B}{2}\left(b^2-a^2\right)+C\left(a-b\right)$$
$$=\frac{b-a}{6}\left[2A\left(b^2+ab+a^2\right)+3B\left(b+a\right)+6C\right]$$
$$=\frac{b-a}{6}\left(2Ab^2+2Aab+2Aa^2+3Bb+3Ba+6C\right)$$
$$=\frac{b-a}{6}\left\{\left(Aa^2+Ba+C\right)+\left(Ab^2+Bb+C\right)+\left[4A\left(\frac{a+b}{2}\right)^2+4B\left(\frac{a+b}{2}\right)+4C\right]\right\}$$
$$=\frac{b-a}{6}\left[g(a)+g(b)+4g(\frac{a+b}{2})\right]$$
$$=\frac{b-a}{6}\left[f(a)+f(b)+4f(\frac{a+b}{2})\right]$$
写成程序就是这样：
template<typename Function>
double simpson(double l, double r, Function func) {
double mid = (l + r) / 2.0;
return (r - l)*(func(l) + 4.0 * func(mid) + func(r)) / 6.0;
}
但是一次拟合计算所得结果误差太大，于是在每次计算时将区间分成两半分别代入公式，将所得结果之和与上次计算所得的整个区间的积分值相比较，若差别不大则停止计算，若差别较大则继续分别递归计算两区间的值。
template<typename Function>
double _asr(double l, double r, Function func, double lastans) {
double mid = (l + r) / 2.0;
double lans = simpson(l, mid, func), rans = simpson(mid, r, func);
if (fabs(lastans - lans - rans) < eps)
return lans + rans;
else
return _asr(l, mid, func, lans) + _asr(mid, r, func, rans);
}
这样就可以智能地得到精度一定的数值积分。
例题：洛谷P4525 【模板】自适应辛普森法1 P4526 【模板】自适应辛普森法2
const double eps = 1e-8; // 要求的精度，一般相应加上一两位
template<typename Function>
double simpson(double l, double r, Function func) {
double mid = (l + r) / 2.0;
return (func(l) + 4.0 * func(mid) + func(r))*(r - l) / 6.0;
}
template<typename Function>
double _asr(double l, double r, Function func, double ans) {
double mid = (l + r) / 2.0;
double lans = simpson(l, mid, func), rans = simpson(mid, r, func);
if (fabs(ans - lans - rans) < eps)
return lans + rans;
else
return _asr(l, mid, func, lans) + _asr(mid, r, func, rans);
}
template<typename Function>
double asr(double left, double right, Function func) {
return _asr(left, right, func, simpson(left, right, func));
}
时间复杂度：玄学，与\(-\log_2eps\)正相关