Is there a higher-order version of cumtrapz()?

Viewed 1737

Introduction

Suppose I have N points x(1:N) at which I have function values f(1:N), for example:

x = [ 0.0795, 0.1327, 0.1395, 0.5133, 0.6470, 0.7358, 0.7640 ];
f = [ 0.0388, 0.4774, 0.4547, 0.0784, 0.3241, 0.2818, 0.9667 ];

I want to calculate the cumulative integral of f with respect to x using these data.

Low-Order Solution

In MATLAB, I can do this easily using cumtrapz():

>> result = cumtrapz( x, f )

result =

     0    0.0137    0.0169    0.1165    0.1434    0.1703    0.1879

Problem

Unfortunately, cumtrapz() uses the trapezoidal method for numerical integration, which is insufficient for my purposes.

Higher-order methods exist, like Simpson's Rule, but to my knowledge there isn't a function that performs the cumulative version of Simpson's Rule for nonuniform grids at the MATLAB File Exchange, or anywhere else.

Does a higher-order version of cumtrapz() already exist? If not, what would I have to do to implement it myself?

1 Answers

I don't know of another method however you could use an interpolation with a pchip, spline, or some other method to increase the resolution. Then use cumtrapz to get a closer approximation of the numeric integral.

It would be up to you to determine what method is applicable to your function.

Example using a sine function and spline

>> x = linspace(0,pi,5);
>> f = sin(x);
>> intF = cumtrapz(x,f);
error = 2-intF(end)
error =
    0.1039

>> x2 = linspace(x(1),x(end),numel(x)*10); %Up sample by 10x
>> f2 = interp1(x,f,x2,'spline');   %Interpolate with spline
>> intF2 = cumtrapz(x2,f2);
>> error = 2-intF2(end)  %MUCH LESS ERROR
error =
   -0.0038

enter image description here

Related