[Numpy-discussion] Dealing with types in extension modules
Nathan Bell
wnbell@gmail....
Sat Sep 13 11:32:20 CDT 2008
On Wed, Sep 10, 2008 at 11:59 PM, Lane Brooks <lbrooks@mit.edu> wrote:
>
> When writing an numpy extension module, what is the preferred way to
> deal with the all the possible types an ndarray can have?
>
> I have some data processing functions I need to implement and they need
> to be generic and work for all the possible numerical dtypes. I do not
> want to have to re-implement the same C-code for all the possible types,
> so the way I approached it was to use a C++ template function to
> implement the processing. Then I have a dispatching function that
> checks the type of the input ndarray and calls the correct template. Is
> there a better way?
>
In scipy.sparse there is a 'sparsetools' extension module that uses
C++ templates to support 14 (or so) numpy types. In sparsetools, SWIG
typemaps are used to dispatch the correct template based on the numpy
data type. For instance, a typemap decides to call foo<float> if the
array has type PyArray_FLOAT.
SWIG can be a little overwhelming, but it automates dispatching and
can even do automatic upcasting for you. For instance, if your C++
function is foo(double, double) and you pass it (float,int) then SWIG
can upcast both arguments to double before passing them into foo().
Also, you can force arrays to be C-contiguous (as opposed to strided
or Fortran order) and have the correct endianness.
Basically, you can use SWIG to tame the input numpy arrays.
The sparsetools module is here:
http://projects.scipy.org/scipy/scipy/browser/trunk/scipy/sparse/sparsetools
The file numpy.i does most of the heavy lifiting. For each sparse
matrix format, for instance 'csr', there is a separate SWIG file
('csr.i') that instantiates each function in the corresponding header
file ('csr.h'). The output of 'swig -c++ -python csr.i' is the
actual extension module ('csr_wrap.cxx'). If you look at one of these
(massive!) files you can see how SWIG dispatches the appropriate
function.
Also, in order to support complex types (e.g. complex128), there is a
wrapper in complex_ops.h that overloads the standard arithmetic
operators for complex numbers. Then complex types are handled exactly
as above.
http://projects.scipy.org/scipy/scipy/browser/trunk/scipy/sparse/sparsetools/complex_ops.h
Depending on your problem, SWIG may be overkill. OTOH you may be able
to get everything to need from the sparsetools source code. Feel free
to pillage it as you require :)
Should you go the SWIG path, I can help explain some of the more cryptic parts.
--
Nathan Bell wnbell@gmail.com
http://graphics.cs.uiuc.edu/~wnbell/
More information about the Numpy-discussion
mailing list