Re: Better quadrature routine in octave

[Pedro Gonnet]

Ok, thanks for the pointers! I've re-worked the code to make it more
compact and to avoid most of the subsrefs. It is now 10-20% faster, but
still crawls in comparison to quadgk.

Any further suggestions?

Cheers, 
Pedro


On Thu, 2010-04-15 at 13:30 +0200, Jaroslav Hajek wrote:
> On Tue, Apr 13, 2010 at 4:54 PM, Pedro Gonnet <address@hidden> wrote:
> >
> > Hello again,
> >
> > Ok, I've checked Shampine's original paper on quadgk/quadva and there
> > are no real, i.e. systematic or parameterized, tests for the interval
> > transformation, so I wouldn't know what else to test my code against.
> >
> > How do you guys want to proceed on this? If anything needs
> > fixing/improving, just tell me what and I'll get to it.
> >
> > Cheers, Pedro
> >
> 
> I finally got into testing this, and for some functions of the test
> suite (e.g. floor(exp(x)) on 0..3) it is indeed a bit slow.
> I'm sure it can be improved, but I would start with simplifying the code.
> 
> Your code has a lot of repeated patterns. The quadruples of analogical
> vars like V_1, V_2, ... can be replaced by cell arrays, repeated
> similar code blocks can be moved into subfunctions. I started with
> some refactoring, the result is attached (this is *not* a working
> version) if you wish to continue.
> 
> Regarding the optimizations, there are some immediate targets. Stuff like
> 
>         % check for nans and clean them up
>         nans = [];
>         for i=1:length(fx), if ( ~isfinite (fx(i)) ), nans = [ i ,
> nans ]; fx(i) = 0.0; endif; endfor;
> 
>         ...
> 
>         % re-instate the nans
>         for i=nans, fx(i) = NaN; endfor;
> 
> is *much* better written in a vectorized way as
> 
>         % check for nans and clean them up
>         nans = find (! finite (fx));
>         fx(nans) = 0;
> 
>         % re-instate the nans
>         fx(nans) = NaN;
> 
> the repeated nested references are another potential source of
> sluggishness, I think lots of them can be eliminated.
> In general it's better to first simplify the code if possible, then
> start with optimizations.
> 
> regards
>

cquad.m
Description: Text Data