[SciPy-Dev] sparse.block_diag improvement

[Severin Holzer-Graf]

Hi,

i had to generate very large sparse block diagonal matrices.
I found out that the existing sparse.block_diag is way to slow for
this and took over 30 secs to generate the matrix i needed.
The reason why the current sparse.block_diag implementation does that
very slow is because it generates everything in a row-wise manner,
which involves looping and a lot of memory.

I wrote my own method that creates the indptr and indices directly.
This works for equally shaped matrices that should be arranged in a
sparse block diagonal matrix.
A benchmark can be found here: http://de.tinypic.com/r/30ue59h/8
It's about a x600 speedup.

My question is, is there interest in integrating this in Scipy, and if
yes, how can this be done best?
I tried extending the current block_diag in a fork, but the indexing
trick only works for equally sized matrices. So either add a check in
the beginning of the existing block_diag and than executing the fast
way, or writing another function, especially for equally sized
matrices.

Best
Severin