Re: [Help-gsl] Performance issue with ODE solvers

[Michael Stauffer]

Hi,

Has anyone come up with anything new about this issue that I may have missed
on the list?

I've noticed what might be the same issue with ODE solver performance. I'm
porting some Matlab code to C and in Matlab I'm using the "ode45" solver
which is described as "based on an explicit Runge-Kutta (4,5) formula, the
Dormand-Prince pair. It is a one-step solver.". In GSL, I'm using
gsl_odeiv_step_rk4: "4th order (classical) Runge-Kutta." 

My C code, using GSL, ends up calling my ODE function anywhere from 10-30x
more often than the Matlab code does, which kills most or all of the speed
gain from coding in C. I'm using the same error tolerances for both. The GSL
code produces fewer time steps in the end, but is calling the ODE function
many more times so it could be an issue with solving for too high a
tolerance.

I'm new to ODE solvers, and am coming them strictly as a programmer, so I'm
pretty clueless beyond the basic principles of their performance. Does
anyone know how close the Matlab "ode45" solver type is to GSL's
gsl_odeiv_step_rk4? Should they end up running with a similar number of
iterations, and yield a similar number of time steps?

I'm going to look at Sundial and see if I can try that method for
comparison.

Many thanks,
Michael



> Date: Sun, 06 Jan 2008 21:36:18 +0000
> From: Brian Gough <address@hidden>
> 
> At Wed, 2 Jan 2008 17:51:57 +0100,
> Jochen K|pper wrote:
> > 
> > Sebastien,
> > 
> > have you received any replies, comments?
> > Do you meanwhile have more information on the issue?
> > 
> > Have people made similar observations for other problems?
> 
> The main thing I'd need to know before looking at it more is 
> how the user-specified error tolerances are used by each 
> package -- otherwise it's not clear that the runtimes can be 
> compared on a like-for-like basis. If that is ok, it may be 
> that GSL routines are overestimating the error, and hence 
> taking smaller steps.
> 
> > Sundials is BSD licensed - are there ways to (legally) 
> incorporate it
> > into GSL? Brian, how would we proceed to ask the original 
> authors for  
> > a contribution to GSL? What's your opinion on that?
> 
> The license is compatible, so it can be used side-by-side 
> with GSL as an independent packages.  I think that is the 
> best way to use it.
> 
> > On the other hand, there are several articles in scientific 
> journals,
> > including a long one in ACM TOMS. Anybody up to implementing it in  
> > GSL based on these descriptions? Actually, if no-one is 
> taking up on  
> > it, maybe with a found background in ODEs could mentor this as a  
> > Google summer of code project?
> 
> If somebody wanted to implement some of the algorithms in the 
> GSL framework that would be cool.
> 
> -- 
> Brian Gough