【自适应辛普森】积分计算

Simpson 公式

[int_l^rf(x)dxapprox frac{r-l}{6}[f(l)+4f(frac{l+r}{2})+f(r)] ]

Simpson 公式将 ((l,f(l)),(r,f(r)),(frac{l+r}{2},f(frac{l+r}{2}))) 看作一个抛物线，求得近似解

具体过程如下：
令 (g(x)=ax^2+bx+c) 过 ((l,f(l)),(r,f(r)),(frac{l+r}{2},f(frac{l+r}{2}))) 三点
积分求得

[G(x)=int_0^xg(x)dx= frac{1}{3}ax^3+frac{1}{2}bx^2+cx+d ]

那么有

[egin{eqnarray*} && int_l^rg(x)dx=G(r)-G(l)\ && = frac{1}{3}a(r^3-l^3)+frac{1}{2}b(r^2-l^2)+c(r-l)\ && = frac{r-l}{6}[2a(l^2+lcdot r+r^2)+3b(l+r)+6c]\ && = frac{r-l}{6}{(acdot l^2+bcdot l+c)+(acdot r^2+bcdot r+c)+4[a(frac{l+r}{2})^2+b(frac{l+r}{2})+c]}\ && = frac{r-l}{6}[g(l)+4g(frac{l+r}{2})+g(r)] end{eqnarray*}]

对于一个抛物线能取到等号，而将 (f(x)) 在区间 ([l,r]) 的端点和中点近似看作一个抛物线，就能取到一个近似值

double f(double x) { return ...; }
double simpson(double x) { return (r-l)*(f(l)+4*f((l+r)/2)+f(r))/6; }

自适应辛普森法

求出对区间 ([l,r]) 使用 Simpson 公式的结果，再求出对区间 ([l,mid]) 和 ([mid,r]) 使用的结果之和，若两者误差在精度要求范围内，就直接返回结果，否则递归求解两个子区间

double asr(double l,double r,double ans) {
	double mid = (l+r)/2,L = simpson(l,mid),R = simpson(mid,r);
	if (fabs(L+R-ans) <= eps) return L+R;
	return asr(l,mid,L)+asr(mid,r,R);
}