about trapz and quad integral

10 次查看（过去 30 天）
Adam
Adam 2014-8-4
回答： Yu Jiang 2014-8-5
Hi all;
I have the folowing data:
x_1=data(:,1);
y_1=data(:,2);
xx=linspace(x(1),x(end),10000);
pp=spline(x,y);
yy=ppval(pp,xx);
%%integral over [x(1),x(end)]
I=trapz(xx,yy)
Int1=quad(@ppval,x_1(1),x_1(end),[],[],pp)
the value of the integral calculated using quad and trapz are different!!!!!!
I=trapz(x,y) gives I=1.353085426296705e+09
and Int1=2.329171651204020e+10 and I obtained the following warning:
Warning: Maximum function count exceeded; singularity likely.
Why there is such differnce ?
how can I avoid this warning?
thank yoy for your help
cheers
  2 个评论
Yu Jiang
Yu Jiang 2014-8-4
Hi Adam
I think the question can be better answered if the following issues can be explained:
  • First, the code uses trapz (See Documentation) to integrate over the interval of [x(1),x(end)] but use quad (See Documentation) over the interval of [x_1(1),x_1(end)]. Can you confirm that x and x_1 are the same?
  • Can you share your data set, i.e., the matrix data? Or can you let me know how I can generate this data set to reproduce the issue?
- Yu
Adam
Adam 2014-8-5
Hi Yu,
thanks for your reply.
Find attached please the data: data=psf110
cheers

请先登录，再进行评论。

请先登录，再回答此问题。

回答（1 个）

Yu Jiang
Yu Jiang 2014-8-5
Let us first take a look at the differences between quad (See Documentation) and trapz (See Documentation). quad approximates the integral of a given function by using the quadrature method, i.e., breaking the area down into rectangles. On the other hand, trapz approximates integral by breaking the area down into trapezoids. Therefore, trapz generally provides more accurate results.
From the data set you provided, a huge spike around x=0 can be observed. This seems to be related with the warning message “singularity likely.” According to the diagnostic section in the documentation , 'Maximum function count exceeded' indicates that the integrand has been evaluated more than 10,000 times. A nonintegrable singularity is likely. There seems to be no way to change that maximum number.
If you would like to calculate this integral by using quadrature method, there are two workarounds as follows.
  • Note that instead of using trapz, which uses recursive adaptive Simpson quadrature, you may want to use the function integral (See Documentation) , which uses global adaptive quadrature. In this case, the last line in your code can be changed to
Int1 = integral(@(x) ppval(pp,x),xx(1),xx(end))
I have tested it. The result is 1.3307e+09.
  • If quad has to be used and you would like to improve the accuracy and avoid the warning, I would suggest you break down the problem into multiple integration problems over smaller intervals, calculate them individually, and sum up the results.

类别

在 Help Center 和 File Exchange 中查找有关 Numerical Integration and Differentiation 的更多信息

标签

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by