什么是数值积分
数值积分可以用来求定积分的近似值。对于很多函数来说，我们是可以使用初等函数来表示出其积分的，对于这种函数，只需要求出不定积分然后代入值就能得到定积分了。
可是除此之外还有许多难求的函数和没法使用初等函数表示的函数。当我们想要求出它们的定积分的时候，需要使用数值积分来求解。
在ACM中一些题目需要使用数值积分来求解，以下列出一些求数值积分的方法，由简单到难，而对ACMer来说最重要的是复合Simpson，其精度较高，且可调精度，是乱搞积分几何的利器。
我从这学的：网易公开课 MIT 数值积分
法一·黎曼和
黎曼和是用将区间等长分为n段，然后用矩形去逼近函数，每段的长为Δx。
可以选择每段左侧的函数值作为矩形的高，也可以选择每段右侧的函数值作为矩形的高。
若设n+1个函数值从左至右为x0,x1?xn，可得如下公式：
Left Hand Riemann=Right Hand Riemann=Δxf(x0)+Δxf(x1)+?+Δxf(xn?1)=Δx∑i=0n?1f(xi)Δxf(x1)+Δxf(x2)+?+Δxf(xn)=Δx∑i=1nf(xi)
这种逼近比较粗糙，不过能比较好地传达数值积分的概念。
法二·梯形法
黎曼和虽然简单，但是精度堪忧，没法很好地模拟逼近函数，接下来介绍第二种方法。
可以看出用矩形逼近的时候有很多空缺，使用梯形去逼近，就能大大提高精度了，看上去很像了。
沿用之前的xi，由梯形公式我们可以得到如下公式：Trapezoid=Δx(f(x0)+f(x1))2+Δx(f(x1)+f(x2))2+?+Δx(f(xn?1)+f(xn))2=Δx(f(x0)2+f(x1)+?+f(xn?1)+f(xn)2)=Left Hand Riemann+Right Hand Riemann2
可以看出梯形法求出的就是左右黎曼和的平均值。
公式中第一个和最后一个需要除以2其他都能合出来1个。
法三·Simpson公式
前两种都比较简单，但是精度比较差，第三种方法是用抛物线去逼近。
使用Simpson公式首先需要n为偶数。
将整个区间分为n2段，每段的底为2Δx，高比较难搞，但是是有公式的：Simpson height=f(xl)+4f(xm)+f(xr)6
于是有公式：
Simpson=Δx3(f(x0)+4f(x1)+2f(x2)+?+2f(xn?2)+4f(xn?1)+f(xn))
法四·复合Simpson
Simpson公式的精度其实已经相当不错了，可是相对于ACM中所需的精度仍然有差距。
为了提高精度，我们需要多次重复利用Simpson公式。
我们先定义被积函数为f(x)，定义Simpson(l,r)=(r?l)f(l)+4f(l+r2)+f(r)6
这也就是常规的Simpson公式，再定义函数RSimpson为复合Simpson。
当我们要求RSimpson(l,r)，先令m=l+r2
当Simpson(l,r)≈Simpson(l,m)+Simpson(m,r)时我们就认为精度够了返回其中一个。
当不满足的时候我们就再分段去求RSimpson(l,m)+RSimpson(m,r)
这样得到的精度就比较高了，而且通过定义≈的范围可以调整精度。
总结公式如下：
RSimpson(l,r)= Simpson(l,r) approximate
RSimpson(l,m)+RSimpson(m,r) else
复合Simpson的实现
1 inline double getAppr(double fl, double fm, double fr, double l, double r) {
2 return (fl+4*fm+fr)*(r-l)/6.0;
3 }
4
5 double Simpson(double l, double r, double fl, double fr) {
6 double m = (l+r)/2, lm = (l+m)/2, rm = (r+m)/2;
7 double fm = f(m), flm = f(lm), frm = f(rm);
8 double vlr = getAppr(fl, fm, fr, l, r);
9 double vlm = getAppr(fl, flm, fm, l, m);
10 double vrm = getAppr(fm, frm, fr, m, r);
11 return fabs(vlr-vlm-vrm) < EPS ? vlr : Simpson(l, m, fl, fm)+Simpson(m, r, fm, fr);
12 }
复合Simpson的实现2(Natureal的代码)
1 inline double getAppr(double l,double r){
2 return (f(l) + 4.0*f((l+r)/2) + f(r)) * (r - l) / 6.0;
3 }
4
5 double Simpson(double l,double r){
6 double sum = getAppr(l,r);
7 double mid = (l+r)/2;
8 double suml = getAppr(l,mid);
9 double sumr = getAppr(mid,r);
10 return (fabs(sum - suml - sumr) < EPS) ? sum : Simpson(l, mid) + Simpson(mid, r);
11 }