MATLAB中quad的误差问题

Question

y = @(s)sin(0.5*s.^2);
quad(y, 0, 26.9);

结果误差远远超过了缺省值，这是为什么？
而如果将上限增减0.01，即：
quad(y, 0, 26.89);
或quad(y, 0, 26.91);
误差就突然很小了，这到底是什么原因啊，百思不得其解，望各位高手指点一下啊，谢谢！
即使只增加0.001，
即quad(y, 0, 26.901);
就没误差了，但是加0.0001误差就产生了~

Answer 1

MATLAB中 quad 函数是用Simpson's rule方法，即利用给定区间内的二次曲线进行拟合法求积分的，属于数值积分。

quad 的完整调用格式是 [Q,fcnt] = quad(funfcn,a,b,tol,trace,varargin)

第三个参数tol的缺省值是 1e-6（我用的是R2009a的版本） ，所以

quad(y, 0, 26.9）

相当于 

quad(y, 0, 26.9,1e-6) 

是以 1e-6 为求积误差的。

另外，这个quad计算定积分使用的Simpson's rule方法，在区间较大时，无法保证多项式拟合能够无穷逼近，特别地区间较大时该方法对于highly oscillatory(高度震动)的函数可能失效。

解决方法：把[0,26.9]分解成两个小区间分别积分——在matlab中输入：

quad(y,0,13)+quad(y,13,26.9)则会给出正确结果：0.9185

或者，你积分时设置更小的积分精度或使用 quadl 函数。quadl 函数是用自适应Lobatto正交法求积分的，求解精度比qaud 高。

附图是把积分一直算到了28时的结果，图中显示的尖点都是计算错误的地方，可见区间越大，积分越不可信。

Answer 2

你这个问题非常好，我测试了一下，发现quad(y, 0, 26.9)=-0.3410要计算1873次函数，而quad(y, 0, 26.89)=0.9222是1997次，quad(y, 0, 26.91)=0.9124是2017次。
不出意外的话quad(y, 0, 26.9)原本应该计算2007次左右才对。在这里出现了程序bug，要深入到quad的源代码中才能排除问题。
另，quad(y, 26.9, 26.89)=0.0038计算了函数13次，数值0.0038，基本正确；quad(y, 26.91, 26.9)=0.0061计算了函数13次，数值0.0061也基本正确。
因此有quad(y, 0, 26.9)=quad(y, 0, 26.89)-quad(y, 26.9, 26.89)=0.9184；或quad(y, 0, 26.9)=quad(y, 0, 26.91)+quad(y, 26.91, 26.9)=0.9185；或者取两个值得平均有quad(y, 0, 26.9)=0.91845。本回答被提问者和网友采纳

Answer 3

你试一下:
quad(y, 0, 26.9,1e-7)
最后那个加一个精度，可以搞定。
另外你可以用符号进行积分。