Поиск собственных векторов и собственных значений разреженной матрицы с помощью ARPACK (вызывается из PYTHON, MATLAB или как подпрограмма FORTRAN)

Несколько дней назад я задал вопрос, как найти собственные значения большой разреженной матрицы. Ответов не получил, поэтому решил описать потенциальное решение.

One question remains: 
Can I use the python implementation of ARPACK 
to compute the eigenvalues of a asymmetric sparse matrix.  

В начале хотелось бы сказать, что совсем не обязательно вызывать подпрограммы ARPACK напрямую из программы-драйвера FOTRAN. Это довольно сложно, и у меня никогда не получалось. Но можно сделать следующее:

#

ВАРИАНТ 1: Питон

#

Можно установить numpy и scipy и запустить следующий код:

import numpy as np
from scipy.linalg import eigh
from scipy.sparse.linalg import eigsh
from scipy.sparse import *
from scipy import * 

# coordinate format storage of the matrix

# rows
ii = array([0, 0, 1, 1, 1, 2, 2, 2, 3, 3, 4])
# cols.
jj = array([0, 1, 0, 1, 2, 1, 2, 3, 2, 3, 4])
# and the data
data=array([1.,-1.,-1., 2.,-2.,-2., 1., 1., 1., 1., 1.])
# now put this into sparse storage (CSR-format)
m=csr_matrix( (data,(ii,jj)), shape=(5,5) )
# you can check what you did 
matrix([[ 1, -1,  0,  0,  0],
        [-1,  2, -2,  0,  0],
        [ 0, -2,  1,  1,  0],
        [ 0,  0,  1,  1,  0],
        [ 0,  0,  0,  0,  1]])

# the real part starts here
evals_large, evecs_large = eigsh(m, 4, which='LM')
# print the largest 4 eigenvalues
print evals_all
# and the values are
[-1.04948118  1.          1.48792836  3.90570354]

Что ж, все это очень мило, особенно потому, что это доставляет нам радость от чтения очень "хорошо написанного" руководства по ARPACK.

У меня проблема с этим, я думаю, что это не работает с асимметричными матрицами. По крайней мере сравнение результатов с матлабом было не очень убедительным.

#

ВАРИАНТ 2: МАТЛАБ

#
       % put your data in a file "matrix.dat"
       % row col. data
       % note that indexing starts at "1"
       1 1 1.
       1 2 -1.
       ......

       load matrix.dat
       M = spconvert(matrix)

       [v,d] = eig(M)

       % v - contains the eigenvectors
       % d - contains the eigenvalues 

Я думаю, что использование Matlab намного проще и работает для асимметричных матриц. Ну у меня разреженная матрица 500000х500000, так будет ли это работать в матлабе.... Еще чашка чая! Я должен отметить, что с помощью python я смог загрузить матрицу такого размера и вычислить ее собственные значения без особых проблем.

Ваше здоровье,


person Alexander Cska    schedule 07.07.2014    source источник
comment
есть ли вариант только для Fortran? Я имею в виду, как избежать ручного кошмара ARPACK   -  person Garini    schedule 08.08.2017
comment
Это было давно, так как я сделал это. В основном для симметричных матриц вы можете использовать реализацию Python, а для общих случаев — Matlab. ARPACK — это Фортран. Чтобы запустить его, вам нужно написать свою собственную функцию умножения векторных матриц и вызвать ее в сочетании с ARPACK. К сожалению, у меня всегда были некоторые сомнения по поводу моей реализации, и я просто использовал Matlab. В общем, инструкция довольно сложная. В свой вопрос я включил несколько примеров того, как работать с Python и Matlab.   -  person Alexander Cska    schedule 09.08.2017