python - solve - numpy linear algebra tutorial

La forma más rápida de calcular k los valores propios más grandes y los vectores propios correspondientes con números

Question

Tengo una gran matriz simétrica NxN densa y quiero los vectores propios correspondientes a los valores propios k más grandes. ¿Cuál es la mejor manera de encontrarlos (preferiblemente usando numpy pero quizás en general usando blas / atlas / lapack si esa es la única manera de ir)? En general, N es mucho más grande que k (digamos N> 5000, k <10).

Numpy parece tener solo funciones para encontrar los valores propios más grandes k si mi matriz de inicio es escasa.

Answer 1

En realidad, las rutinas dispersas también funcionan con matrices de números densos, creo que usan algún tipo de iteración subespacial de Krylov, por lo tanto, deben computar varios productos de vectores matriciales, lo que significa que si es k << N, ¿las rutinas dispersas podrían ser (marginalmente? ) Más rápido.

Echa un vistazo a los documentos http://docs.scipy.org/doc/scipy/reference/tutorial/arpack.html

y el siguiente código (ve a tomar un buen café con amigos hasta que termine)

import numpy as np
from time import clock
from scipy.linalg import eigh as largest_eigh
from scipy.sparse.linalg.eigen.arpack import eigsh as largest_eigsh

np.set_printoptions(suppress=True)
np.random.seed(0)
N=5000
k=10
X = np.random.random((N,N)) - 0.5
X = np.dot(X, X.T) #create a symmetric matrix

# Benchmark the dense routine
start = clock()
evals_large, evecs_large = largest_eigh(X, eigvals=(N-k,N-1))
elapsed = (clock() - start)
print "eigh elapsed time: ", elapsed

# Benchmark the sparse routine
start = clock()
evals_large_sparse, evecs_large_sparse = largest_eigsh(X, k, which=''LM'')
elapsed = (clock() - start)
print "eigsh elapsed time: ", elapsed

Answer 2

En SciPy, puedes usar la función linalg.eigh , con el parámetro eigvals .

eigvals: tuple (lo, hi) Índices de los valores propios más pequeños y más grandes (en orden ascendente) y los eigenvectores correspondientes que se devolverán: 0 <= lo <hi <= M-1. Si se omite, se devuelven todos los valores propios y vectores propios.

Que en su caso debe ajustarse a (Nk,N-1) .