Как вычислить собственный вектор стохастической матрицы столбца в С++

У меня есть стохастическая матрица столбца A (n-by-n реальная, неотрицательная матрица) и я хочу решить следующее уравнение на С++: Ax = x

Я предполагаю, что мне нужно найти собственный вектор x, где собственное значение должно быть равно 1 (правильно?), но я не мог понять это на С++. До сих пор я проверил некоторые математические библиотеки, такие как Seldon, CPPScaLapack, Eigen... Среди них Eigen кажется хорошим вариантом, но я не мог понять, как использовать какую-либо из них для решения приведенного выше уравнения.

Можете ли вы дать мне несколько предложений/фрагментов кода или идей для решения уравнения? Любая помощь высоко ценится.

Спасибо.


person Aleyna    schedule 05.12.2010    source источник
comment
О какой большой матрице мы говорим? В любом случае, поскольку вы уже знаете собственное значение, обратная итерация звучит как то, что вам нужно, чтобы получить соответствующий собственный вектор (ну, после приведения к форме Хессенберга или какой-либо другой более простой форме, конечно).   -  person    schedule 05.12.2010


Ответы (2)


Поскольку наибольший собственный вектор стохастической матрицы $M$ равен единице, вы можете найти этот собственный вектор итерацией (если только вы не очень плохо угадываете начальные значения).

Начните с произвольно выбранного начального вектора $v_1$, сумма значений (вероятностей) которого равна единице. Примените $M$ к $v_1$, чтобы получить $Mv_1$. Теперь перенормализуйте этот новый вектор $Mv_1$, то есть разделите на сумму его элементов, чтобы получить $v_2$. Это новый вектор вероятностей, и он будет ближе к желаемому собственному вектору (если только ваше первоначальное предположение не оказалось ортогональным собственному вектору).

Повторяйте этот процесс, пока $v_k$ не станет стабильной. Это должен быть вектор с $Mv_k = v_k$, как и хотелось.

person Carl Brannen    schedule 26.03.2011
comment
Большое спасибо! Я собираюсь попробовать это. Проблема в том, что моя матрица имеет размер более 1M на 1M. Вам не кажется, что каждая итерация займет много времени? - person Aleyna; 30.03.2011
comment
@Алейна; 1 000 000x1 000 000? Ой! Это может занять некоторое время процессора. - person Carl Brannen; 30.03.2011

Другой метод заключается в вычислении ядра вашей матрицы за вычетом единичной матрицы. Это может быть или не быть быстрее, чем использование степенной итерации, как объяснил Карл, в зависимости от размера матрицы и других собственных значений. Итерация мощности лучше, когда матрица становится больше, а второе собственное значение отдаляется от единицы.

Идея состоит в том, чтобы переписать Ax = x как Ax - x = 0. Затем использовать Ix = x, где I обозначает единичную матрицу. Таким образом, Ax - x = 0 эквивалентно Ax - Ix = 0 или (A-I) x = 0. Таким образом, ваш собственный вектор x лежит в ядре (или нулевом пространстве) A-I.

На этой учебной странице объясняется, как вычислить ядро ​​матрицы с помощью Eigen. Какой-то непроверенный код:

MatrixXf M;
/* initialize M */
FullPivLU<MatrixXf> lu_decomp(M);
VectorXf x = lu_decomp.kernel().col(0);
/* x now contains the vector you want */

Вы можете обнаружить, что ядро ​​пусто. Это означает, что либо матрица на самом деле не является стохастической, либо вам необходимо адаптировать порог (см. страницу, указанную выше).

person Jitse Niesen    schedule 28.03.2011
comment
спасибо за руководство! Я проверю это! Я попробовал метод ядра, но он оказался нулевой матрицей, хотя входная матрица, безусловно, стохастическая. - person Aleyna; 30.03.2011