(Preface)
自适应辛普森，一个非常有趣的算法。
感觉也不需要多么高深的数学知识，果然算法这种东西实践起来总比想象中容易许多。
辛普森积分
考虑对于任何一个函数图象，在(Delta x)极小的时候，都可以将这一范围内的函数视为一段抛物线。
假设这段抛物线拟合后的函数为(Ax^2+Bx+C)，然后推导一下：
[egin{align}int_l^r(Ax^2+Bx+C)&=frac A3(r^3-l^3)+frac B3(r^2-l^2)+frac C3(r-l)\&=frac{r-l}6(2A(l^2+lr+r^2)+3B(l+r)+6C)\&=frac{r-l}6((Al^2+Bl+C)+4(A(frac{l+r}2)^2+B(frac{l+r}2)+C)+(Ar^2+Br+C))\&=frac{r-l}6(f(l)+4f(frac{l+r}2)+f(r))end{align}
]
这就是所谓辛普森积分的公式了：
[int_l^rf(x)dxapproxfrac{r-l}6(f(l)+4f(frac{l+r}{2})+f(r))
]
考虑到这一公式的前提是(Delta x)极小，因此我们需要不断二分递归以保证精度，不难发现复杂度显然爆炸。
所以我们需要改进这一算法。
自适应辛普森
递归次数少，难以保持精度；递归次数多，复杂度又难以接受。
而所谓的自适应辛普森算法，就是在二分的过程中根据精度，自动控制区间分割的大小。
可以直接看代码。
/*友情提醒：I表示inline，可忽略；DB表示double；Con DB&表示const double&，可以直接视作double*/
namespace SimpsonInt
{
I DB f(Con DB& x) {return /*需要积分的函数*/;}
I DB Simpson(Con DB& l,Con DB& r) {return (r-l)*(f(l)+4*f((l+r)/2)+f(r))/6;}//辛普森积分
I DB ASR(Con DB& l,Con DB& r,Con DB& ans,Con DB& eps)//自适应辛普森算法
{
DB mid=(l+r)/2,L=Simpson(l,mid),R=Simpson(mid,r);
if(fabs(L+R-ans)<=15*eps) return L+R+(L+R-ans)/15;//控制精度（15据说是研究结果）
return ASR(l,mid,L,eps/2)+ASR(mid,r,R,eps/2);//二分递归，注意eps除以2
}
}
模板：自适应辛普森法1
- 求(int_L^Rfrac{cx+d}{ax+b}dx)。
- (-100le L<Rle 100)且(R-Lge 1)。
真正的模板题。
只要在上面给出代码中的函数(f(x))里填上(c*x+d)/(a*x+b)即可。
模板：自适应辛普森法2
- 求(int_0^{infty}x^{frac ax-x}dx)。
- 若积分发散，输出(orz)
- (ale50)。
可证明(a<0)时积分发散，(a>0)时积分收敛。（具体证明这里略去）
这道题看似是要求无限积分，实际上当(x)到(15)左右(x^{frac ax-x})就已经小得可以忽略不计了。
所以说，我们只要求(int_0^{15}x^{frac ax-x}dx)即可。
最后要做的就是在代码中的函数(f(x))里填上pow(x,a/x-x)。
(Postscript)
自适应辛普森算法是一个简单而又有用的算法。
具体来说，就是学了并不会造成多大负担，而万一碰到了这种题目又能发挥出巨大的作用。
所以，余认为，这还是非常值得一学的。