[Numpy-discussion] Best way to broadcast a function from C

[James Crist]

I suspected as much. This is actually part of my work on numerical
evaluation in SymPy. In its current state compilation to C and autowrapping
*works*, but I think it could definitely be more versatile/efficient. Since
numpy seemed to have solved the broadcasting/datatype issues internally I
hoped I could reuse this.

I'll look into nditer (as this is for a library, it's better to do it
*right* than do it now), but right now it's looking like I may need to
implement this functionality myself...

Thanks,

-Jim


On Thu, Aug 21, 2014 at 4:50 PM, Nathaniel Smith <njs at pobox.com> wrote:

> On Thu, Aug 21, 2014 at 2:34 AM, James Crist <crist042 at umn.edu> wrote:
> > All,
> >
> > I have a C function func that takes in scalar arguments, and an array of
> > fixed dimension that is modified in place to provide the output. The
> > prototype is something like:
> >
> > `void func(double a, double b, double c, double *arr);`
> >
> > I've wrapped this in Cython and called it from python with no problem.
> What
> > I'd like to do now though is get it set up to broadcast over the input
> > arguments and return a 3 dimensional array of the results. By this I mean
> >
> > a = array([1, 2, 3])
> > b = array([2.0, 3.0, 4.0])
> > c = array([3, 4, 5])
> >
> > func(a, b, c) -> a 3d array containing the results of func for (a, b, c)
> =
> > (1, 2.0, 3), (2, 3.0, 4), (3, 4.0, 5)
> >
> > I'm not sure if this would qualify as a ufunc, as the result of one
> function
> > call isn't a scalar but an array, but the effect I'm looking for is
> similar.
> > Ideally it would handle datatype conversions (in the above `a` and `c`
> > aren't double, but `func` takes in double). It would also be awesome to
> > allow an argument to be a scalar and not an array, and have it be
> broadcast
> > as if it were.
> >
> > I'm just wondering what the best way for me to hook my code up to the
> > internals of numpy and get this kind of behavior in an efficient way.
> I've
> > read the "writing your own ufunc" part of the docs, but am unsure if what
> > I'm looking for qualifies. Note that I can change the inner workings of
> > `func` if this is required to achieve this behavior.
>
> I don't think it's currently possible to write a ufunc that maps
> scalars to fixed-length arrays.
>
> There's probably some not-too-terrible way to do this using the C
> nditer interface, but IIRC the docs aren't for that aren't very
> helpful.
>
> The simplest approach is just to do a bit of work at the Python level
> using the standard numpy API to do error-checking and set up your
> arrays in the way you want, and then pass them to the C
> implementation. This is not the best way to get a Right solution that
> will handle all the funky corner cases that ufuncs handle, but it's by
> far the fastest way to get something that's good enough for whatever
> you need Right Now.
>
> -n
>
> --
> Nathaniel J. Smith
> Postdoctoral researcher - Informatics - University of Edinburgh
> http://vorpus.org
> _______________________________________________
> NumPy-Discussion mailing list
> NumPy-Discussion at scipy.org
> http://mail.scipy.org/mailman/listinfo/numpy-discussion
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mail.python.org/pipermail/numpy-discussion/attachments/20140822/e1ac2db7/attachment.html>