Амосов А.А., Дубинский В.А., Копченова Н.В. Вычислительные методы для инженеров (1994) (1095853), страница 44
Текст из файла (страница 44)
Существует несколько способов избавления от уже вычисленных собственных чисел и соответствующих собственных векторов с целью избежать их повторного вычисления. Эти способы принято называть аенерпыаанием. Более подробную информацию о них можно найти, например, в [191, [62]. ~ 8.3. Метод обратных итераций 1 Вычисление собственных векторов методом обратных итераций. многие из методов решения проблем собственных значений лучше приспособлены для вычисления собственных значений, чем собственных векторов. Поэтому целесообразно рассмотреть задачу вычисления собственного вектора ед при условии, что уже найдено достаточно точное приближение Л '.
к собственному значению Л . Если исходить непосредственно из определения собственного вектора, то еу следует искать как нетривиальное решение однородной системы уравнений (А — Л Е)х= О (8.24) (А — Л'Ь)х = О. (8.25) Так как матрица А — Л*Е заведомо невырождена, то решением системы (8.25) является только х = О. Следовательно, непосредственное численное решение системы (8.25) не дает возможность вычислить собственный вектор.
Одним из эффективных методов вычисления собственных векторов является нетод обратных итераций. В этом методе приближения к собственному вектору определяют последовательным решением систем уравнений (А — Л *Е) у' ~"" = х' "' 1 (8.26) с последующей нормировкой решения: 227 с вырожденной матрицей А — Л~Е. Однако Л известно лишь прибли- женно и в действительности при таком подходе вместо системы (8.24)- придется решать систему >П х(о) = Е сеь (и1 »> у'" = Е а(е;. (и1 Так как (А — Л*.е)у(" = Е а;(Л, — Л *.) е,, то систему уравнений (8.26) .) (и1 при е = О можно записать в виде Еа,(Л,— Л~)е,= Ес;е,. (и1 ~ (и1 с; Приравнивая коэффициенты при е;, получим а; =, . Следова- $ тельно, с, 1 Л.
— Л*. Л.— (и1 (8.28) Если ~Л. — Л'.~ < ~ Л; — Л'.~ для всех 1 ~ 1'„то второе слагаемое в правой части формулы (8.28) мало по сравнению с первым. Поэтому у(') с. — — е и вектор р»> оквзыезетея близким по и прввле ию к Л. — Л*. 1 1 вектору е. Можно показать, что вектор х()б), вычисляемый на Й-й итерации, имеет вид сбе) +, с е, я(>б) — ф(>с) 228 В качестве начального приближения берут ненулевой вектор я(о) с произвольно выбираемыми или даже случайными компонентами. Часто удовлетворительным является выбор в(о) = (1, 1> ..., 1)т, Чтобы понять механизм действия метода, рассмотрим случай, когда А — матрица простой структуры, а Л вЂ” простое собственное значение. Представим векторы х(О) и у'1' в виде линейных комбинаций собственных векторов е), е2, ..., е»): где ~ф< "~ ~ -+ ~ с1 ~ при Й-+ оо.
Вектор х'"' сходится к е1 по направлению со скоростью геометрической прогрессии, знаменатель которой Л*~ =;~, Гл, — л,7 Если абсолютная погрешность значения Л*. много меньше рассто- 1 яния от Л до ближайшего из остальных собственных чисел (что эквивалентно выполнению условия ~Л~ — Л*.~ < ~Л, — Л*~ для всех з Ф Я, то метод обратных итераций сходится очень быстро. Чаще всего достаточно сделать 1 — 3 итерации. Пример 8.6. Используя метод обратных итераций, найдем на 6-разрядной десятичной ЭВМ собственный вектор матрицы (8.4), отвечающий собственному значению Л1 м Л * = -7.8728.
1 Возьмем х'о> = (1, 1, 1)т. Тогда система (8.26) при Ф = 0 примет вид 9 8728р~ Оу~ + 5дз = 1 1.2у~ + 2.4729уг + 6дз = 1, й — Уг + 0.3728рз = 1. Вычисляя ее решение методом Гаусса, получим у~~1~ = — 123909, у~~1~ = 99398.3, '~> = 65749. После нормировки имеем х'~~ = -0.720736, х'~' 0 578166 х~з~~ 0 382440 Подставляя в правую часть системы (8.26) вместо вектора х<о> вектор х<~> и вычисляя решение, находим у',~' = — 59671.5, у'~> = — 47867,7, у'~' = 31663.1. Полученные после нормировки значения х<~> = -0.720737, х<т> = — 0.578166, х< т> = 0.382440 с машинной точностью совпадают со значениями, полученными на первой итерации.
Итерации следует прервать и считать ре- зультатом вектор е~ и (-0.720737, -0.578166, 0.382440)т. 3 а м е ч а н и е 1. Так как Л '- почти совпадает со значением Л, .1 при котором йеС(А — Л Е) = О, то матрица А — Л*Е очень плохо ) обусловлена. Возникает естественное опасение, что большие погрешности, неизбежные при реализации вычислительного процесса (8.26), (8.27) на ЗВМ, могут существенно повлиять на свойства приближений. К счастью, это не так; ошибка, возникающая при 229 численном решении системы (8.26) (которая может быть по величи не сравнима с у<"'), оказывается почти пропорциональной вектору е .
В данном случае плохая обусловленность системы не ухудшает, а улучшает ситуацию. В справедливости сказанного легко убедиться, если повторить вычисления из примера 8.5, используя микрокалькулятор. Полученные значения у<1) и у<~) наверняка будут иметь мало общего с указанными в примере. Тем не менее приближения к<1) и в<~) практически совпадут с вычисленными выше, отличаясь, возможно, только знаком. 3 а м е ч а н и е 2. Записав равенство (8.26) в виде у<)< " = (А — Л*Ь~ 1х<"), замечаем, что метод обратных итераций — это 1 просто степенной метод, примененный к матрице (А — Л *.Ю) '.
2. Метод обратных итераций с использованием отношения Рэлея. Одной из проблем применения метода обратных итераций является необходимость получения хорошего приближения к собственному значению. В случае, когда А — симметричная матрица, можно попытаться использовать для оценки Л~ отношение Рэлея. Этот подход приводит к следующему комбинированному методу: Л<") = р(х<")) = (Ат«""), х<") ), >с)~0, (А — Л < "' Е) у< "'1' = х< "', (8.29) (8.30) ~< й+1) в< Ус+1) (8.31) Предполагается, что начальный вектор а<о) нормирован, т.е.
~х«)) $ = = 1. Если начальное приближение хорошо аппроксимирует по направлению вектор е, то метод (8.29) — (8.31) сходится очень быстро. Например, если Л вЂ” простое собственное значение, то сходимость оказывается кубической. Одним из способов получения удовлетворительного начального приближения к<о) является выполнение нескольких итераций степенного метода. 3 а м е ч а н и е. В случае, когда метод (8.29) — (8.31) сходится> он позволяет одновременно эффективно вычислять и собственное значение Л и соответствующий собственный вектор с . 230 $8 4. ЯЛ-алгоритм В настоящее время лучшим методом вычисления всех собственных значений квадратных заполненных матриц общего вида (умеренного порядка) является ЯВ-алгоритм>.
1. Основной >,>>В-алгоритм. Опишем итерационную процедуру, явля>ощуюся основой алгоритма. Она существенно использует возможность разложения произвольной матрицы в произведение ортогональной и верхней треугольной матриц, т.е. так называемое Щ-разложение (см. 1 5.10). На 1-й итерации с помощью метода отражений или метода вращений вычисляют ЧВ-разложение матрицы А'о> = А, имеющее вид А>о> = ц>Л (8.32) Затем строят матрицу А'>> = 1с>ф.
Заметим, что из равенства (8.32) следует, что Я> — — Я>А'о> и поэтому А''> = (~г'А>о> Ц,. Таким обра- Л> х х О Лз х ... х О 0 Лз (8.33) 0 0 0 ... Л„ Этот метод независимо был предложен в 1960 г. в России математиком В.Н.Кублановской и в 1961 г. в Англии системным программистом Дж.Фрэнсисом. 231 зом, матрицы А''> и А'с> подобны (см. ~ 8.1) и поэтому имеют общий набор собственных значений Л>, Л2, ..., Л На 2-й итерации находят ДЯ-разложение матрицы А~», имеющее вид А<>> = фА2, и вычисляют матрицу А'2' = Л~ф, подобную матрице А»>. На (1 + 1)-й итерации вычисляют разложение А> "> = ф„>Я~ > и строят матрицу А> >'"> = В~ >ч и Неограниченное продолжение этого процесса дает последовательность матриц А> > >, А< 2>, ..., А> и>, подобных матрице А. К сожалению, в общем случае, когда собственные значения матрицы могут быть кратными или комплексными, строгое изложение предельных свойств последовательности А < "> потребовало бы введения ряда довольно сложных математических понятий.
Поэтому. ограничимся указанием на то, что при определенных условиях в случае, когда собственные значения матрицы А вещественны и различны по модулю, причем >Л>~ > (Л~) > ... > ~Л,„~, последовательность А< "> сходится к верхней треугольной матрице Л вида Здесь крестиками помечены элементы, в общем случае не равные нулю. Известно, что в рассматриваемом случае элементы а' Ь~ матриц $г А<"', стоящие ниже главной диагонали, сходятся к нулю со скоростью геометрической прогрессии, причем Ь ~ а<.Ь.> ~ ~( С ~ — '), з > г, (8.34) сходится к нулю со скоростью геометрической прогрессии, знаменатель которой д = ~ ~ в 1. Л1 ~л,, Для того чтобы уменьшить число арифметических операций, матрицу А предварительно с помощью подобных преобразований отражения или вращения приводят к так называемой форме Хессенбсрза'.
Ьп Ь1г ... Ь1,з-1 Ь1в Ьг1 Ьгг ... Ьг щ1 Ьга О Ьзг ... Ьз,а-1 Ьзт (8.35) О О ... Ь , Ь „ Матрица Н, в которой равны нулю все элементы Ь,1, такие, что 1 > 7 + 1 Герхард Хессенберг (1874 — 19%) — немецкий математик. 232 т,е, скорость сходимости а<.".' к нулю определяется величиной отношеФг л;, ния Л; к Л (заметим, что ~ — ~ < 1 при 1 > Я. л,1 Следует подчеркнуть, что ДЛ-алгоритм не обеспечивает обычной поэлементной сходимости и речь идет о сзодилости к треугольной матрице (8.33) по форле. Она характеризуется сходимостью к нулю поддиагональных элементов, тогда как наддиагональные элементы от итерации к итерации могут существенно изменяться.