I wrote a small library for creating Numpy ufuncs using Cython. Just copy the pxd and pxi linked at http://wiki.cython.org/MarkLodato/CreatingUfuncs and include the pxi at the top of your file. Then you can create 1-input/1-output or 2-input/1-output ufuncs operating of floats, doubles, or long doubles easily.
For those unfamiliar with ufuncs, these are functions that operate on NumPy arrays - for example, all of the mathematical operators of NumPy (add, subtract, multiply, exp, etc) are ufuncs. They are nice for the user since they handle broadcasting and all the other stuff you expect of NumPy. They are also nice for the writer of the function, since one needs only write a simple 1-D loop (or even a function that takes scalars as input and output) and NumPy takes care of the rest. That is, the author need not worry about multi-dimensional arrays, error handling, broadcasting rules, output buffers, etc. For more details, see http://docs.scipy.org/doc/numpy/reference/ufuncs.html. Although NumPy makes it easy to write a ufunc, there's still a fair amount of boilerplate. To make the process simpler, I created this library to do most of the work for you. Just define a function that takes a scalar (or two) as input and returns a scalar, then register the function using register_ufunc_d() (or similar). The wiki page shows how to create an "lgamma" ufunc. Here's another simple example: {{{ include "numpy_ufuncs.pxi" cdef double foo_dd(double a, double b): return (a-b) / (a+b) foo = register_ufunc_dd(foo_dd, "foo", "returns (x1-x2)/(x1+x2)", PyUFunc_Zero) }}} Now in Python you can do: {{{ In [1]: from foo import foo In [2]: import numpy as np In [3]: a = np.arange(16).reshape(4,4) In [4]: b = np.arange(30,34) In [5]: foo(a,b) Out[5]: array([[-1. , -0.9375 , -0.88235294, -0.83333333], [-0.76470588, -0.72222222, -0.68421053, -0.65 ], [-0.57894737, -0.55 , -0.52380952, -0.5 ], [-0.42857143, -0.40909091, -0.39130435, -0.375 ]]) In [6]: foo(32, 16) Out[6]: 0.33333333333333331 }}} The great thing is that it's easy to write, easy to use, and fast! I would like to add more common cases (particularly integer) and more generic cases. The NumPy C API only has 1-D loops for floating point (e.g. PyUFunc_d_d), so I would have to write this myself. For now this fits my needs. I'd be happy to hear any comments or suggestions you have. Perhaps this can be distributed with Cython? It would be really nice to be able to easily create ufuncs in Cython right out of the box. ~Mark Lodato _______________________________________________ Cython-dev mailing list [email protected] http://codespeak.net/mailman/listinfo/cython-dev
