根据您的功能(以及x的分布方式),您可以通过首先对数据进行spline插值来获得更高的准确度:
pp = spline(x,y);
quadgk(@(t) ppval(pp,t), [range])那是快节奏的方式。 Ther是一种更快速,更直接的方法,但这种方法很简单且透明度更低:
result = sum(sum(...
bsxfun(@times, pp.coefs, 1./(4:-1:1)) .*... % coefficients of primitive
bsxfun(@power, diff(pp.breaks).', 4:-1:1)... % all 4 powers of shifted x-values
));作为一个例子,为什么所有这些都有用,我从here借用了这个例子。确切的答案应该是
>> pi/2/sqrt(2)*(17-40^(3/4))
ans =
1.215778726893561e+00定义
>> x = [0 sort(3*rand(1,5)) 3];
>> y = (x.^3.*(3-x)).^(1/4)./(5-x);我们发现
>> trapz(x,y)
ans =
1.142392438652055e+00
>> pp = spline(x,y);
>> tic; quadgk(@(t) ppval(pp,t), 0, 3), toc
ans =
1.213866446458034e+00
Elapsed time is 0.017472 seconds.
>> tic; result = sum(sum(...
bsxfun(@times, pp.coefs, 1./(4:-1:1)) .*... % coefficients of primitive
bsxfun(@power, diff(pp.breaks).', 4:-1:1)... % all 4 powers of shifted x-values
)), toc
result =
1.213866467945575e+00
Elapsed time is 0.002887 seconds.因此,trapz低于0.07低估了该值。使用后两种方法,误差小一个数量级。此外,spline方法的可读性较低的版本要快一个数量级。
所以,掌握这些知识:明智地选择:)