[SciPy-Dev] Beginning work on robustly allowing function pointers in integrate - quad

[Pauli Virtanen]

Hi,

10.10.2013 21:27, Brian Lee Newsom kirjoitti:
> Thanks for all your replies on this matter.  After consideration it seems
> that a reasonable place for work would be on allowing integrate.quad to
> robustly accept function pointers from a compiled language, so that they
> may bypass the underlying callbacks to python.  The current implementation
> does not support additional arguments and because of this cannot be used
> with higher dimensions of integration using dblquad/tplquad/nquad when it
> is needed most.

This would certainly be an useful addition. Other places where this
could come useful are in optimization, ODEs, f2py, and 3rd party
packages, and this should be kept in mind while designing the
implementation.


There are a couple of problems that need to be solved:

(i) Converting arguments from Python objects to something generally
usable for a C/Fortran/etc. routine and actually constructing the
C-level call.

(ii) How to pass the raw function + data pointers in to the low-level
computational routines?

(iii) Stashing and retrieving the arguments somewhere.

I didn't yet think this through well, so what's below is just some
first-order thoughts. What perhaps should first be done is to just write
down a proposal of how this should work, before starting coding.

***

In detail:

(i)

What sort of C/Fortran signatures should be accepted? The simplest
option is probably to just support user functions of the form

	return_t func(..., void *params)

where `...` are the arguments as they are passed in by the Fortran/C
algorithm code, and `*params` is then tacked on by indirection. This
removes the need for multiple trampoline routines.

As `params`, we could support e.g. the contents of any Python buffer
objects. The user would then be able to pack the parameters then as
necessary into the array, and interpret them as they see fit in the
objective function.


(ii)

This would need to involve some additions to the code in quad(); for
instance it should check the type of function pointer passed in. Ctypes
makes this possible, but it may be useful to have our own representation
for it (which can be constructed from the ctypes one) that's easier to
deal with.

There are also a couple of different types of ways to integrate with
lower-level code in Scipy: Cython, raw C, f2py.

f2py can take raw function pointers with some trickery.
The C code varies, possible to make to work (e.g. PyCapsules and so on).
Cython code is also reasonably easy to make to work.

Because of this, it may be best to keep the type checks etc. on the
Python side, and just pass the verified function pointer to the
low-level routine.


(iii)

Many of the algorithms we currently have is ancient Fortran code, which
does not accept any "additional arguments" for the functions it calls.

To work around this, there a couple of options. One option is to patch
the codes in Scipy so that an additional pointer argument is always
passed along. This is the simplest solution, but makes the computational
codes to diverge from the originals.

The second option is to generate a new low-level function on the fly
that calls the routine with the correct extra parameters. Difficult to
make portable, so doesn't seem feasible.

The third option is to store the extra data to a global variable / TLS
and be careful with re-entrancy when handling it. quad() currently does
this, so extra args could be added easily in the same mechanism. It may
be possible to refactor some of this a component that can be reused.


-- 
Pauli Virtanen