自适应辛普森法

算法简介以及模板可以在大白书P166~P170页找到，本篇文章主要帮助理解。

首先上一个资料，我就是看这个弄懂的。

http://wenku.baidu.com/link?url=w3VyU2KW8DXFkEgipi-2hjruxvi-XXMXMEWhb8vnuslIQ-ycsa-ESt0S3sIeHQ5l_DUAq7PZ-MOVkwaW815iPsPfqNXfRgRlCJBajtwJR3u

书上例题的解题报告

http://blog.csdn.net/xl2015190026/article/details/53516642

看大白书P169上的图 2-37

上面有一个公式∫f(x)dx=(Δx/3)*(y0+4*y1+y2)+(Δx/3)*(y2+4*y3+y4)+...+(Δx/3)*(yn-2+4*yn-1+yn)

这就是辛普森公式，其中(Δx/3)*(y0+4*y1+y2)算的是[x0,x2]的曲线下面积，其他同理，全部求和以后就是[a,b]的曲线下面积，也就是函数的定积分。这种方法不需要求出原函数，因此可以用来求积不出来的函数的定积分。

那为什么[x0,x2]的曲线下面积可以近似等于(Δx/3)*(y0+4*y1+y2)呢？

其实不过是将[x0,x2]的曲线用一个过(x0,y0),(x1,y1),(x2,y2)的抛物线近似替代，然后再用最基础的定积分求出此抛物线的定积分，并将之当成近似解罢了。

推导：

设近似替代的抛物线方程为f(x)=ax^2+bx+c。

因为此抛物线过三点(x0,y0),(x1,y1),(x2,y2)，所以得到方程组：

y1=ax0^2+bx0+c                                                               |

y2=ax1^2+bx1+c=a((x0+x2)/2)^2+b(x0+x2)/2+c           |          ①

y3=ax2^2+bx2+c                                                               |

对抛物线f(x)=ax^2+bx+c进行定积分。

∫f(x)dx=(1/3)ax^3+(1/2)bx^2+cx|

          =(1/6)(2ax^3+3bx^2+6cx)|

          =(1/6)(2ax2^3+3bx2^2+6cx2)-(1/6)(2ax0^3+3bx0^2+6cx0)

          =(1/6)[2a(x2^3-x0^3)+3b(x2^2-x0^2)+6c(x2-x0)]

          =(1/6)[2a(x2-x0)(x2^2+x2x0+x0^2)+3b(x2-x0)(x2+x0)+6c(x2-x0)]

          =(1/6)(x2-x0)[2a(x2^2+x2x0+x0^2)+3b(x2+x0)+6c]

          =(1/6)2Δx(2ax2^2+2ax2x0+2ax0^2+3bx2+3bx0+6c)

          =(1/3)Δx{(ax0^2+bx0+c)+4[a((x0+x2)/2)^2+b(x0+x2)/2+c]+(ax2^2+bx2+c)}

将方程组①代入上式，得

∫f(x)dx=(1/3)Δx{(ax0^2+bx0+c)+4[a((x0+x2)/2)^2+b(x0+x2)/2+c]+(ax2^2+bx2+c)}

          =(Δx/3)*(y0+4*y1+y2)
推导完毕

下面讲一讲自己对模板的理解：

主要讲两个地方。

一、

if(fabs(L+R-A)<=15*eps) return L+R+(L+R-A)/15;

首先要说一下那个15是什么鬼，其实我也不知道，但是能明显感觉到是用来修正误差的。

你看，为什么最后返回的是L+R+(L+R-A)/15而不是L+R也不是L+R+(L+R-A)/1呢?可能就是说明需要修正一下而且除以15修正效果最好吧，具体为什么这么修正，为什么15最好，其实我自己也不知道= =。

其实不修正也能AC。

二、

return asr(a,c,eps/2,L)+asr(c,b,eps/2,R);

观察代码后发现每递归一层eps就要除以2。

大白书上说“容易近似的地方少划几份，不容易近似的地方多划几份”。
事实上就算我们不除以2，按这个算法跑下去肯定也能做到“容易近似的地方少划几份，不容易近似的地方多划几份”。

因为只要不够近似，我们就会一直递归下去。而递归的越深，划得就越细。

所以我觉得除以2只是为了再次强化这个要求吧。

其实不除以2也能AC。