[SciPy-dev] changes to scipy.sparse
Nathan Bell
wnbell@gmail....
Fri Oct 26 15:02:35 CDT 2007
There have been some recent changes to scipy.sparse that may
potentially affect others.
http://projects.scipy.org/scipy/scipy/changeset/3465
Sparse matrices (csr_matrix, csc_matrix, and coo_matrix) now support
64-bit indices and the following dtypes: int8, uint8, int16, int32,
int64, float32, float64, complex64, complex128
Unsupported dtypes will be accepted by the sparse objects, but they'll
be upcasted (by SWIG) to one of the above types when sent to one of
the sparsetools functions. Operations between different data types
should follow the expected upcasting procedure (e.g. complex64 *
float64 -> complex128).
The only issue I found when using non-floating point dtypes occurred
in linsolve.spsolve() which assume a FP dtype. I created a method
spmatrix.asfptype() that returns a sparse matrix withan appropriate FP
type when necessary. In addition, the member variable spmatrix.fptype
was removed since it pertained only to SuperLU. A few small edits to
linsolve.spsolve were made to account for these changes. I have not
tested the UMFPACK bindings yet, so let me know if there are any
problems.
I encountered an unexpected bug related to io.mio (matlab IO). The
mio methods were creating sparse matrices with non-native byte order.
The SWIG magic in sparsetools accounts for non-contiguous arrays, but
not endianness. As a temporary fix I cast all input arrays to native
byte order using:
44 def to_native(A):
45 if not A.dtype.isnative:
46 return A.astype(A.dtype.newbyteorder('native'))
47 else:
48 return A
This seems to work, however I'd like to improve my SWIG to account for
such issues instead. I believe the current SWIG bindings in the numpy
SVN check for, but do not change byteorder. Does anyone have SWIG
that accounts for this?
With the additional data types, the SWIG output for sparsetools now
weighs in a nearly 3MB and takes about 2.5 minutes to compile on my
1.8 GHz Athlon 64. While this is a rather large amount of C++ I don't
think we can assume that end-users will have SWIG, let alone the a
recent SVN version required by sparsetools for some systems. The size
and compile time are due to the large number of functions (23*9*2 =
414) generated by the templates.
A while ago I wrote a small C++ wrapper for numpy complex types so
templated code (like sparsetools) can make sense of T + T and T == 0
etc. The numpy SWIG examples seem to lack complex support, so I
thought this may be valuable to others.
http://projects.scipy.org/scipy/scipy/browser/trunk/scipy/sparse/sparsetools/complex_ops.h
Also, I plan to add support for block CSR and CSC formats to sparse in
the near future.
As always, comments and questions are welcome.
--
Nathan Bell wnbell@gmail.com
More information about the Scipy-dev
mailing list