[SciPy-user] sparse: slicing/indexing problems
Ed Schofield
schofield at ftw.at
Tue Jul 25 19:00:06 EDT 2006
On 25/07/2006, at 4:01 PM, William Hunter wrote:
> Dear Scipy-Users;
>
> I'm struggling (again/still) with sparse matrices, but now
> it has to do
> with slicing:
>
> 1) The following (inserting a sub-matrix into a bigger
> matrix at arbitrary
> positions and then updating it) works very well for dense
> arrays:
>
>
> from numpy allclose, import ix_, ones, zeros
> Ad = zeros((5,5),float)
> Adsub = 5*ones((3,3),float)
> indeks = [3,0,4]
> Ad[ix_(indeks,indeks)] += Adsub
>
>
> The last line above will typically be inside a <for> loop.
>
> 2) Doing this with sparse matrices is not 'possible',
> although I'm sure it
> can be done, but I've tried everything I know (slicing as
> per normal,
> over one axis only, etc.) This is the best I can do, but I
> don't like it:
>
> from scipy import sparse
> As = sparse.lil_matrix(5,5)
> col = 0
> for posr in indeks:
> row = 0
> for posc in indeks:
> As[posr,posc] += Adsub[row,col]
> row += 1
> col += 1
> # We want a True here...
> print allclose(Ad,As.todense())
>
>
> Any other suggestions, i.e., how do I get rid of the
> <for>'s?
I don't think you could gain much here by eliminating the loops. The
reason is that += will need to check each element anyway and take one
of four different code paths, depending on whether it's a zero being
converted to a nonzero, or a nonzero to a zero, etc. I'd probably
write it like this:
As = sparse.lil_matrix((5,5))
for row, posr in enumerate(indeks):
for col, posc in enumerate(indeks):
As[posr,posc] += Adsub[row,col]
This seems to work fine (I think your row and column increments were
swapped).
Using a single for loop and fancy indexing should work too (saving a
line of code), but doesn't at the moment, since support for fancy
indexing with in-place operations is incomplete (ticket #226). I'll
add this soon, but probably only after the 0.5.0 release.
> I'm using Enthon 2.4-1.0.0.beta4 (but my brand new Linux
> box is coming!!!)
Oooh, goody!
-- Ed
More information about the SciPy-User
mailing list