On May 4, 10:04 pm, Thierry Dumont <[email protected]> wrote:
> Do you mean that it is possible to define the RHS as a Cython *callback*
> function? or is there an other trick ? Can you give me a pointer to that ?
The code can for instance be found in the file sage/gsl/ode.pyx (the
gsl directory has other classes where this construction is also used),
and relies on the following trick: the ode solver does not call the
RHS defined by the user, but instead calls a C function called c_f
cdef int c_f(double t,double* y, double* dydt,void *params):
The array "params" contains as its first argument a pointer (suitably
wrapped) to the user-defined RHS (which could be written in Python).
If you check the default implementation of c_f, you will see that this
function does nothing else but unpacking the RHS and calling it, so
that the user never even has to know about c_f.
Now, if you do want to provide a compiled RHS, all you have to do is
supply your own c_f implementation. This is what the docs have to say
(ellipses for readability):
sage: ode_solver?
Unfortunately because Python functions are used, this solver is
slow on systems that require many function evaluations. It is
possible to pass a compiled function by deriving from the class
``ode_sysem`` and overloading ``c_f`` and ``c_j`` with C
functions
that specify the system. The following will work in the
notebook:
(...)
cdef class van_der_pol(sage.gsl.ode.ode_system):
cdef int c_f(self,double t, double *y,double *dydt):
dydt[0]=y[1]
dydt[1]=-y[0]-1000*y[1]*(y[0]*y[0]-1)
return GSL_SUCCESS
cdef int c_j(self, double t,double *y,double
*dfdy,double *dfdt):
(...)
I hope this helps. When I first saw this, I thought it was a very
elegant trick to ensure that you have both speed and readability.
All the best,
Joris
--
To post to this group, send an email to [email protected]
To unsubscribe from this group, send an email to
[email protected]
For more options, visit this group at http://groups.google.com/group/sage-devel
URL: http://www.sagemath.org