scipy обобщенная собственная проблема с положительной полуопределенной

Привет, ребята!!!

Я хочу вычислить обобщенное собственное разложение формы:

Lf = лямбда Af

с помощью функции scipy.sparse.linalg.eigs, но получите эту ошибку:

/usr/local/lib/python2.7/dist-packages/scipy/linalg/decomp_lu.py:61: RuntimeWarning: диагональное число 65 равно нулю. Сингулярная матрица. RuntimeWarning) ** При входе в DLASCL параметр номер 4 имел недопустимое значение

Я передаю три аргумента, диагональную матрицу, положительную полуопределенную (PSD) матрицу и числовое значение K (первые K собственных значений). Функция eigs Matlab хорошо работает с теми же входными параметрами, но в SciPy, как я понял, для вычислений с PSD мне также нужно указать параметр sigma.

Итак, мой вопрос: есть ли способ избежать установки параметра сигмы, как в MatLab, или, если нет, как подобрать значение сигмы?

С нетерпением жду советов или подсказок... Заранее спасибо!


person maksmara    schedule 24.10.2012    source источник


Ответы (1)


Ошибка, по-видимому, означает, что в вашей обобщенной собственной проблеме

L x = lambda A x

матрица A не является положительно определенной (проверьте строку документации eigs - в вашем случае матрица, вероятно, сингулярна). Это требование для режима ARPACK 2. Однако вы можете попробовать указать sigma=0 для переключения в режим 3 ARPACK (но обратите внимание, что значение параметра which в этом случае инвертируется!).

Теперь я не уверен, что делает Matlab, но есть вероятность, что он вычисляет псевдоинверсию, а не инверсию A. Чтобы подражать этому, сделайте

from scipy.sparse.linalg import LinearOperator
from scipy.linalg import lstsq

Ainv = LinearOperator(matvec=lambda x: lstsq(A, x)[0], shape=A.shape)
w, v = eigs(L, M=A, Minv=Ainv)

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

В качестве альтернативы вы можете попробовать указать ненулевое значение sigma. То, что вы должны выбрать, зависит от задействованных матриц. Это влияет на выбранные собственные значения --- например, с which='LM' это те, для которых lambda' = 1/(lambda - sigma) велико. В противном случае его, вероятно, можно выбрать произвольно, конечно, для прогресса Крылова, вероятно, будет лучше, если преобразованные собственные значения lambda', которые вас интересуют, будут хорошо отделены от других собственных значений.

person pv.    schedule 24.10.2012