 NunodeSousa -4 years ago 330
Python Question

# Numpy matrix product - sparse matrices

Let us consider a matrix A as a diagonal matrix and a matrix B a random matrix, both with size N x N.
We want to use the sparse properties of the matrix A to optimize the dot product, i.e. dot(B,A).

However if we compute the product using the sparcity properties of the matrix A, we cannot see any advantage (and it is much slower).

``````import numpy as np
from scipy.sparse import csr_matrix
# Matrix sizes
N = 1000

#-- matrices generation --
A = np.zeros((N,N), dtype=complex)
for i in range(N):
A[i][i] = np.random.rand()
B = np.random.rand(N,N)

#product
%time csr_matrix(B).dot(A)
%time np.dot(B,A)
``````

Results:

CPU times: user 3.51 s, sys: 8 ms, total: 3.52 s
Wall time: 3.74 s

CPU times: user 348 ms, sys: 0 ns, total: 348 ms
Wall time: 230 ms

How to do it correctly? jotasi

The difference stems from the fact that you convert `B` to a sparse matrix during the timing (minor effect) and even worse, that `dot` is not aware of the fact, that `A` is sparse. If you were to do the conversion before the dot product the sparse dot product is actually faster:

``````import numpy as np
from scipy.sparse import csr_matrix
# Matrix sizes
N = 1000

#-- matrices generation --
A = np.zeros((N,N), dtype=complex)
for i in range(N):
A[i][i] = np.random.rand()
B = np.random.rand(N,N)

Asparse = csr_matrix(A)
Bsparse = csr_matrix(B)

%timeit np.dot(B, A)
%timeit csr_matrix(B).dot(A)
%timeit Bsparse.dot(A)
%timeit csr_matrix.dot(B, Asparse)
%timeit csr_matrix.dot(Bsparse, Asparse)
``````

Gives:
`np.dot(B, A)`: 1 loop, best of 3: 414 ms per loop
`csr_matrix(B).dot(A)`: 1 loop, best of 3: 2.22 s per loop
`Bsparse.dot(A)`: 1 loop, best of 3: 2.17 s per loop
`csr_matrix.dot(B, Asparse)`: 10 loops, best of 3: 32.5 ms per loop
`csr_matrix.dot(Bsparse, Asparse)`: 10 loops, best of 3: 28 ms per loop

As you can see the sparse dot product is much faster in all cases where `A` is in a sparse matrix format which is making `dot` aware of the fact, that `A` is sparse. In your timing, the function is actually doing a conversion of `B` to a sparse format and then a dot product with a dense matrix `A`.

Recommended from our users: Dynamic Network Monitoring from WhatsUp Gold from IPSwitch. Free Download