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

Я написал линейный решатель, использующий отражения/преобразования Хаусхолдера в 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}, очевидно.


c ansi linear-algebra
person nick_name    schedule 25.10.2011    source источник
comment
Крайне маловероятно, что в моем коде существует ошибка, поскольку я могу выполнять операции обращения матриц (это влечет за собой несколько линейных решений подряд). Инверсия всегда верна. Я считаю более вероятным, что моя тактика не будет плодотворной и что есть лучший способ выполнить часть обратной замены, которая не всегда дает вектор 0.   -  person nick_name    schedule 25.10.2011
comment
Как правило, использовать линейный решатель для нахождения собственных векторов - плохая идея (одна из причин - некорректность: (A - lambda I)x = 0 по определению имеет несколько линейно независимых решений). Если вы ищете собственный вектор, связанный с конкретным собственным значением, в LAPACK есть подпрограммы для этой задачи.   -  person Alexandre C.    schedule 25.10.2011

Ответы (2)


arrow_upward
4
arrow_downward

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

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

... из которого вы можете заменить его, не достигая вырожденного 0 вектора. Обычно вы также захотите доставить x в нормализованной форме.

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

person Drew Hall    schedule 25.10.2011
comment
@nick_name: Не то чтобы мне это не понравилось, но я бы пока не принял этот ответ. Вы можете получить что-то лучше, если вы подождете день или около того. - person Drew Hall; 25.10.2011

arrow_upward
1
arrow_downward

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

Если А — матрица

1  2  0
2  1  4
0  4  1

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

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

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

#  #  # | 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] равным единице.

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

person Jitse Niesen    schedule 25.10.2011
comment
Спасибо за подробный ответ! Матрица была просто примером. В качестве примера вы можете представить себе любую эрмитову матрицу с линейно независимыми строками. Возможно, тот, что я перечислил, не был хорошим выбором. - person nick_name; 25.10.2011