ANSI C + Численная линейная алгебра - Использование линейного решателя для поиска собственного вектора по заданному собственному значению (выпуск)

Question

Я написал линейный решатель, использующий отражения / преобразования Домохозяина в ANSI C, который решает Ax = b с учетом A и b.Я хочу использовать его, чтобы найти собственный вектор, связанный с собственным значением, например:

(A-lambda*I)x = 0

Проблема в том, что вектор 0 всегда является решением, которое я получаю (до того, как кто-то скажет это, да, у меня естьправильное собственное значение с вероятностью 100%).

Вот пример, который довольно точно иллюстрирует проблему:

Учитывая A-lambda*I (пример просто эрмитов):

1 2 0 | 0
2 1 4 | 0
0 4 1 | 0

Размышления / преобразования домохозяина приведут к чему-то подобному

# # # | 0
0 # # | 0
0 0 # | 0

Обратная подстановка обнаружит, что решение {0,0,0}, очевидно.

Answer 1

Прошло много времени с тех пор, как я написал eigensolver, но я, кажется, вспоминаю, что трюк состоял в том, чтобы реорганизовать его с (A - lambda*I) * x = 0 до A*x = lambda*x.Тогда ваши шаги Домовладельца или Гивенса дадут вам что-то вроде:

# # # | #
0 # # | #
0 0 1 | 1

..., из которого вы можете выполнить обратную замену, не достигнув вырожденного 0-го вектора.Обычно вы также хотите поставить x в нормализованном виде.

Моя память здесь довольно ржавая, поэтому я рекомендую проверить Golub & Van Loan для окончательного ответа.Есть несколько хитростей, чтобы заставить это работать надежно, особенно для несимметричного случая.

Answer 2

Это в основном тот же ответ, что и @Drew, но объясняется немного по-другому.

Если A - матрица

1  2  0
2  1  4
0  4  1

тогда собственные значения лямбда = 1, 1 + sqrt (20), 1-sqrt (20). Для простоты возьмем lambda = 1. Тогда расширенная матрица для системы (A - lambda*I) * x = 0 равна

0  2  0 | 0
2  0  4 | 0
0  4  0 | 0

Теперь вы делаете Домовладелец / Гивенс, чтобы уменьшить его до верхней треугольной формы. Как вы говорите, вы получаете что-то вроде

#  #  # | 0
0  #  # | 0
0  0  # | 0

Однако последний # должен быть нулем (или почти нулем). То, что вы получите, зависит от деталей ваших преобразований, но если я сделаю это вручную, я получу

2  0  4 | 0
0  2  0 | 0
0  0  0 | 0

Теперь вы делаете обратную замену. На первом этапе вы решаете уравнение в последнем ряду. Однако это уравнение не дает никакой информации, поэтому вы можете установить x[2] (последний элемент вектора x) на любое значение, которое вы хотите. Если вы установите его на ноль и продолжите обратную замену этим значением, вы получите нулевой вектор. Если вы установите его на одно (или любое ненулевое значение), вы получите ненулевой вектор. Идея ответа Дрю состоит в том, чтобы заменить последний ряд на 0 0 1 | 1, который устанавливает x[2] в 1.

Ошибка округления означает, что последний #, который должен быть равен нулю, вероятно, не совсем нулевой, а какое-то небольшое значение, например 1e-16. Это можно игнорировать: просто возьмите его в ноль и установите x[2] в единицу.

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