Разреженные собственные значения: scipy.sparse.linalg.eigs медленнее, чем scipy.linalg.eigvals

У меня есть странное явление, хотя scipy.sparse.linalg.eigs должен быть быстрее для разреженных матриц, я понимаю, что он работает медленнее, чем обычный метод eigvals для 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

Таким образом, метод calc_pde_numerical_jacobian строит плотную матрицу якобиана моей системы уравнений, а calc_pde_analytic_jacobian строит разреженную матрицу якобиана аналитически (формат csc). Хотя аналитический метод работает быстрее при построении разреженной матрицы якобиана, при использовании методов поиска собственных значений из scipy метод собственных значений разреженных матриц работает медленнее. Функции, которые я использую для вычисления собственных значений, таковы:

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)

Кто-нибудь знает, как это могло произойти?


python numpy scipy sparse-matrix eigenvalue
person Ohm    schedule 22.03.2017    source источник
comment
Какой размер матрицы?   -  person Warren Weckesser    schedule 22.03.2017
comment
Размер текущих матриц 512x512   -  person Ohm    schedule 22.03.2017
comment
Для разреженных матриц с небольшой пропускной способностью ситуация еще хуже.   -  person Mikhail Genkin    schedule 25.04.2018

Ответы (1)


arrow_upward
2
arrow_downward

Для достаточно больших и разреженных матриц разреженный решатель должен работать быстрее. Я выполнил следующий фрагмент для N в диапазоне (150, 550, 50) и 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

и получил следующее время (измеряется в мс):

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

В этом случае размер, при котором разреженный решатель становится конкурентоспособным, составляет 250-300.

Время, вероятно, будет зависеть от разреженности (то есть, какой процент матрицы отличен от нуля) и от структуры или шаблона ненулевых элементов. Для вашей проблемы разреженный решатель может быть не лучше, пока матрицы не превысят 512x512.

person Warren Weckesser    schedule 22.03.2017