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

[Nathan Woods]

On Thu, Oct 10, 2013 at 3:16 PM, Pauli Virtanen <pav at iki.fi> wrote:

> 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.
>
>
I agree fully with this thought. I think that we're essentially talking
about creating a new template for interaction between SciPy and C/Fortran
code. That's not something that you can just throw together.


> ***
>
> 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.
>

I'm not sure I follow you here. I think you're referring to only accepting
user functions that return a single argument, with a pointer to some
additional parameters? Can you give an example of such a parameter as it
would be used in, say, integrate?


>
> 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.
>
>
I guess I thought that we were talking about replacing most of the code in
quad(). As far as I can tell, the integrate library uses the C/Python API
to communicate with the underlying Quadpack code. One of the main
difficulties I see is that the C/Python API is fairly complex, and
therefore the wrapper code is difficult to modify. I'd be a lot happier if
this could be done with one of the simpler tools (probably Cython, since
that seems to be the way the community is leaning), so that the interface
doesn't become a magic black box. It helps that the Fortran community seems
to support this approach:
http://fortran90.org/src/best-practices.html#interfacing-with-python


> 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.
>

I hadn't noticed this. It's kind of a major difficulty.


>
> 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.
>
>
Yeah, one of the best things about SLATEC is that it's a standard. Don't
really want to mess with that.


> 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.
>

On-the-fly code generation? Yeah, I can see difficulties there. I ran into
the Python version of this when I wrote nquad, and I dealt with it by using
callable objects instead of functions to store the extra data I needed.
Perhaps something similar can be done here?


>
> 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
>
> _______________________________________________
> SciPy-Dev mailing list
> SciPy-Dev at scipy.org
> http://mail.scipy.org/mailman/listinfo/scipy-dev
>


Nathan Woods
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mail.python.org/pipermail/scipy-dev/attachments/20131010/c9e2e74d/attachment.html>