Sparse Matrix Multiplications

[Peter Otten]

dingo_1980 at yahoo.com wrote:

> I was trying to create a sparse matrix using scipy.sparse (100000 X
> 100000) with just the first row and first column filled with ones.
> Lets call this matrix Asp. This is how I created Asp
> 
> from scipy import sparse
> Asp = scipy.lil_matrix(100000,100000)
> for i in range(100000):
>    Asp[i,0] = i
>    Asp[0,i] = i
> Bsp = (Asp.tocsr())*(Asp.tocsr())
> 
> This multiplication to get Bsp is taking a long time (to the order of
> over an hour). Is there some other way of storing such huge sparse
> matrices and doing their multiplications or is this the most optimal
> way?

I think sparse matrices don't help here because the result isn't sparse.
Looking at your example I think you can see that B[i,k] == A[i,0]*A[0,k] for
i,k != 0,0, so for a matrix of that structure you might get away with
something like

# use at your own risk
import numpy

N = 10**4 # I get a MemoryError for a bigger exponent
b = numpy.array(range(N))
a = numpy.zeros((N, N)) + b
a *= a.transpose()
a[0,0] = (b*b).sum()
print a

But for expert advice you better ask on the numpy mailing list...

Peter