Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

Sparse eigenvalues : scipy.sparse.linalg.eigs is slower than scipy.linalg.eigvals

I have a peculiar phenomenon, though scipy.sparse.linalg.eigs should be faster for sparse matrices, I get that it runs slower than the normal eigvals method of scipy:

In [4]: %timeit m.calc_pde_numerical_jacobian(m.initial_state)
10 loops, best of 3: 41.2 ms per loop

In [5]: %timeit m.calc_pde_analytic_jacobian(m.initial_state)
1000 loops, best of 3: 1.42 ms per loop

In [6]: %timeit m.calc_analytic_pde_eigs(m.initial_state)
1 loop, best of 3: 374 ms per loop

In [7]: %timeit m.calc_numeric_pde_eigs(m.initial_state)
1 loop, best of 3: 256 ms per loop 

So the method calc_pde_numerical_jacobian construct a dense matrix of the Jacobian of my system of equations, and calc_pde_analytic_jacobian is constructing a sparse matrix of the Jacobian analytically (csc format). Though the analytical method is working faster in constructing the sparse matrix of the Jacobian, when use the eigenvalue finding methods from scipy, the sparse matrics eigenvalue method is slower. The functions that I use to calculate the eigenvalues are as such:

def calc_numeric_pde_eigs(self,state):
    return linalg.eigvals(self.calc_pde_numerical_jacobian(state))
def calc_analytic_pde_eigs(self,state):
    return sparse.linalg.eigs(self.calc_pde_analytic_jacobian(state),k=6,which='LR',return_eigenvectors=False)

Anyone knows how this could happen?

like image 387
Ohm Avatar asked Oct 30 '22 10:10

Ohm


1 Answers

For sufficiently large and sparse matrices, the sparse solver should be faster. I ran the following snippet for N in range(150, 550, 50) and N = 1000:

In [150]: from scipy import sparse

In [151]: from scipy import linalg

[...]

In [186]: N = 150

In [187]: m = sparse.random(N, N, density=0.05).tocsc()

In [188]: a = m.A

In [189]: %timeit sparse.linalg.eigs(m, k=6, which='LR', return_eigenvectors=False)
10 loops, best of 3: 20.2 ms per loop

In [190]: %timeit linalg.eigvals(a)
100 loops, best of 3: 9.66 ms per loop

and got the following timing (measured in msec):

N                    150   200   250   300   350   400   450   500   1000
sparse.linalg.eig   20.2  22.2  28.9  29.4  48.5  38.6  75.2   57.9   152
linalg.eigvals       9.7  17.0  24.5  37.0  52.7  63.3  82.5  105     482

In this case, the size where the sparse solver becomes competitive is 250-300.

The timing will likely depend on the sparsity (i.e. what percent of the matrix is nonzero) and on the structure or pattern of the nonzero elements. For your problem, the sparse solver might not be better until the matrices are larger than 512x512.

like image 99
Warren Weckesser Avatar answered Nov 15 '22 06:11

Warren Weckesser