# Sparse Matrix Multiplications

Discussion in 'Python' started by dingo_1980, May 13, 2008.

1. ### dingo_1980Guest

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?

dingo_1980, May 13, 2008

2. ### Peter OttenGuest

I think sparse matrices don't help here because the result isn't sparse.
Looking at your example I think you can see that A[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

Peter Otten, May 13, 2008

3. ### Peter OttenGuest

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

Peter Otten, May 13, 2008
4. ### Peter OttenGuest

Sorry, this is nonsense.

Peter Otten, May 13, 2008
5. ### Peter OttenGuest

Maybe not as bad as I feared. I seems you can't transpose and
inplace-multiply, so

# use at your own risk
import numpy

N = 5
b = numpy.array(range(N))
a = numpy.zeros((N, N)) + b
a = a * a.transpose()
a[0,0] = (b*b).sum()
print a

might work if there is sufficient memory available...

Peter

Peter Otten, May 13, 2008