Re: Better quadrature routine in octave

[David Bateman]

Pedro Gonnet wrote:
> Hi Octave-people,
>
> I recently finished a PhD on the topic of adaptive numerical quadrature
> in one dimension and also just recently got a paper accepted in ACM TOMS
> with two new quadrature algorithms that significantly out-perform
> QUADPACK's QAGS, which, as I understand, is the default integrator used
> in octave.
>
> I have attached a copy of the accepted paper and of one of the
> quadrature routines (the doubly-adaptive one, it's the one that does
> best) which should be call-compatible with quadl. If you think this
> would be a good addition to octave, I am willing to invest some time in
> making the code fit for whatever standards may be required for its
> inclusion, so just tell me what you need and I'll be happy to oblige.
>
> Cheers and many thanks -- for your time and for octave,
> Pedro
>   
The big problem with Matlab and Octave for quadrature is that its not a 
good idea to call the user function with scalar arguments but rather 
with the longest vector possible. So a lot of intelligent algorithms 
that reduce the number of function evaluations are outperformed by 
stupid ones that vectorize the call to the user function. I see your 
code tries to call the user function with vector arguments so thats a 
good point.

The are many quadrature algorithms in Octave itself and many others in 
octave-forge. Yes the quadrature used in the octave quad function is 
from qudpack, but I'd really suggest using the quadgk function instead. 
How does your code compare with that function in both matlab and octave 
for the number of function calls for a given precision and the speed?

In any case choice is a good thing and it would be great if you added 
your code to octave-forge in the "integration" package. See the page

http://octave.sourceforge.net/integration/index.html

Cheers
David