Численное решение задач метода наименьших квадратов - Лоуcон Ч.
Скачать (прямая ссылка):
7,<(1 +2,/2)"-,т,/2тах|д}'/>|
6.4. Лоусон
81
f2. Q/i -алгоритм для симметричных матриц
QR-алгоритм Фрэнсиса *) для случая симметричных матриц и при условии использования сдвигов является одним из самых успешных методов во всем численном анализе. Имеется ясное представление о его свойствах и скорости сходимости. В этом параграфе мы опишем алгоритм и заложим основы для доказательства его сходимости.
Вначале симметричную л X n-матрицу А можно преобразовать к трех-диагональному виду посредством ортогональных преобразований подобия: л - 2 преобразований Хаусхолдера или (л - 1) (л - 2)/2 преобразований Гивенса. Поэтому без потери общности мы можем считать, что симметричная матрица, для которой нужно вычислить собственные значения, трех-диагональная.
QR-алгоритм со сдвигами для вычисления спектрального разложения симметричной трехдиагональной л X п-матрицы А можно описать формулами
Ах =А,
Qk(Ak - okI„) = Rk,
RkQkT + okI„ =Aktl,
где для k = 1,2,...
a) Qk - ортогональная матрица; b)Rk — верхняя треугольная матрица; с) ок - к-й сдвиг.
Прежде чем описать метод вычисления сдвигов ок, мы введем обозначения для элементов матриц Ак. Можно проверить (см. упражнение 18.46), что для трехдиагональной матрицы A i все матрицы Ак, к = 2,3,..., также будут трехдиатональными. В этой главе и в приложении В мы будем обозначать диагональные элементы каждой трехдиагональной матрицы Ак (*) . ,
через at , / = 1,..., л, а наддиагональные и поддиагональные элементы
. (fc) _
через bi , i f 2.....л. Теперь мы можем определить сдвиг ок следующим образом.
Каждый сдвиг ок есть то собственное значение нижней угловой 2 X 2-лод-матрицы матрицы Ак
С?» н- к (ш)
которое находится ближе к дл*\
Из (18.1)-(18.3) следует т
Ak + \ =QkAkQk ,
так что все матрицы Ак имеют одни и те же собственные значения. Обозначим их через Xi,..., Х„.
*) ДО-алгоритм был изобретен советским алгебраистом В.Н. Кублановской независимо от Фрэнсиса и одновременно с ним. (Примеч. пер.)
(18.1) (18.2)
(18.3)
Ь2
Если для некоторого * среди элементов есть нулевые, то матрица Ак разлагается в прямую сумму подматриц, в каждой из которых все наддиатональные и поддиагональные элементы отличны от нуля. Поэтому
нам достаточно рассмотреть задачу, в которой все о, ненулевые.
Далее, согласно лемме В.1 приложения В, из того, что все поддиагональные элементы не равны нулю, вытекает, что все собственные значения простые. Таким образом, мы можем ограничиться вычислением спектрального разложения для симметричных трехдиагональных матриц с простыми собственными значениями.
Сходимость -ал го ритма устанавливается следующей теоремой [194,195].
Теорема 18.5 {глобальная квадратичная сходимость QR-алгоритма со сдвигами). Пусть А = At - симметричная трехдиагональная п X п-матрица с ненулевыми поддиагональными элементами. Пусть Ак - матрицы, ортогонально подобные матрице А. которые генерируются QR-алгоритмом со сдвигами в соответствии с формулами (18.1)—(18.3) и (18.4). Тогда
а) каждая матрица A k трехдиагональная и симметричная;
о) элемент Ьп матрицы Ак, стоящий на пересечении п-и строки и п - 1-го столбца, стремится к нулю при к -+«»;
с) сходимость Ьпк^ к нулю асимптотически квадратичная, т.е. существует е > 0, зависящее от А и такое, что для всех к
\ъяк*"\«\ъ™\*.
Доказательство этой теоремы дано в приложении В.
На практике сходимость алгоритма обычно бывает кубической. Однако квадратичная сходимость - это лучшее, что удается доказать. Другие замечания по этому поводу читатель найдет в приложении В.
§ 3. Вычисление сингулярного разложения
Рассмотрим теперь построение сингулярного разложения m X л-матрицы А в предположении, что т> п. Замечания по поводу случая т < п см. в § 1. Сингулярное разложение будет вычислено в два этапа. На первом этапе А
преобразуется к верхней двухдиагональнои матрице I J посредством последовательности (не более чем из 2л - 1) преобразований Хаусхолдера [о] =Qn(-((QlA)H2)...H„) = QTAH, я (18.6)
где
В
Ч\е2 Яг е3
Яп-t еп
-- Ч
(18.7)
6»
i wwmtmm.....ш^п ъшртягмг,.....Е0сьгтУлирова, ь
элементы / + 1, .. ., т столбца /; матрица Hi - так, чтобы аннулировать элементы i + 1,. . ., и строки i - 1.
Заметим, что Н„ - это попросту единичная матрица. Она включена в (18.6), чтобы упростить обозначения; Q„ также будет единичной матрицей при т = и, но при т > п она, вообще говоря, отличается от единичной
Второй этап процесса состоит в применении специальным образом адаптированного -алгоритма к вычислению » сингулярного разложения матрицы В
B = USVT;
(18.8)
здесь Он V - ортогональные матрицы, a S диагональная. Из (18.8) можно получить сингулярное разложение А:
A=(Qlf)[S0](HV)T = u[S0] У''
(18.9)
Обсудим теперь вычисление разложения (18.8). Заметим прежде всего, что если какой-либо элемент et равен нулю, то матрица В распадается на блоки, которые можно обрабатывать независимо друг от друга. Ниже мы покажем, что если какой-либо элемент q( равен нулю, то применение некоторых преобразований также позволяет привести матрицу В к распадающемуся виду.
