Дж. Деммель - Вычислительная линейная алгебра, страница 85

Теорема 7.1. Минимум величины йАЯь — Я»В()г по всем симметричным йхй-матрицам В достигается при В = Ть, при этом вАЯ» — Я»Ввг =!!Тьявг. Пусть Ть = »'Л»'т — спектральное разложение матрицы Ть. Минимум величины ~ЙАР» — Р»ПЦйг, когда Рь пробегает множество и х 1«-матриц с ортонормированными столбцами, таких, что арап(Р») = зрапЯ»), а Р пробегает множество диагональных Л х 1«-матриц, также равен ЙТь„()г и достигается для Рь = Я»»' и.0 = Л.

379 7.2. Метод Редея — Ритце Доказательство. Чтобы упростить записи, временно опустим индекс (е у матриц Ть и дь и будем писать Т = дтАд. Положим В = Т + Е. Мы хотим показать, что минимум величины [[Ад — дВЯ достигается при Е = О. Это будет достигнуто применением теоремы Пифагора в замаскированной форме: [~Ад — дВ[[~ ~= Л „[(Ад — дВ)т(Ад — дВ)] согласно утверждению 7 леммы 1.7 = Л [(Ад - д(Т+ Е))т(Ад - д(Т+ Е))] = Л „„[(Ад дТ)т(Ад дТ) (Ад дТ)т(дЕ) — (дЕ)т(Ад — дт) + (дЕ)'(дЕ)] Л [(Ад дТ)т (Ад дТ) (дтАд Т) Е Ет(дтАд Т) + ЕтЕ] = Л [(Ад дТ)т(Ад — дТ) + ЕтЕ] поскольку дтАд = Т > Л [(Ад дТ)т(Ад — дТ)] поскольку ŠŠ— симметричная положительно т полуопределенная матрица (см.

вопрос 5.5) = [[Ад — дТ[[~ ~согласно утверждению 7 леммы 1.7. Легко определить само значение минимума. Восстанавливая индексы, имеем [[Аде — дьТь[[э = ЦдьТь + д„Тье) — (дьТьЦз = [[д„Ть,[[з = [)Тье[[з. Если заменить дь произведением дь(7, где У вЂ” произвольная ортогональная матрица порядка й, то столбцы матриц дь и дел порождают одно и то же подпространство и справедливы равенства [)Аде — дьВ[[з = [)АдьУ вЂ” дьВЦ)т = ][А(дьП) — (дел)((У~ВП)[[т. Эти величины по-прежнему достигают минимума при В = Ть. Выбирая У = Ъ' так, чтобы матрица ПтТьП стала диагональной, мы решим вторую задачу минимизации в формулировке теоремы. П Эта теорема обосновывает использование чисел Ритца как приближений к собственным значениям. Предположим, что матрица дь вычислена алгоритмом Ланцоша.

В этом случае (см. (6.31)) 380 Глава 7. Итерационные методы для задач на собственные значения и все величины, упомянутые в теореме 7.1, легко находятся. Действительно, существуют хорошие алгоритмы для вычисления собственных значений и собственных векторов симметричной трехдиагональной матрицы Ть (см. разд. 5.3), а норма невязки ОТа 'Ог есть просто число зуг. (Из описания алгоритма Ланцоша известно, что число )7ь неотрицательно.) Все зто облегчает использование оценок для погрешностей в приближенных собственных значениях и собственных векторах, содержащихся в следующем утверждении. Теорема 7.2.

Пусть Тг, Ть и Дь — те же матрицы, что и в формуле (7.1). Пусть Тг = УЛУ'т — спектральное разложение матрицы Ть, где Ъ' = [ьз,..., ьг] — ортогональная матрица, а Л = йа8(Вз,..., Вг). Тогда: 1. Найдутся к (необязательно наибольших) собственных значений сгз, ..., аь матрицы А, такие, что 1Вз — сзз) < ЙТь Йг для г = 1,...,к. Если матрица Ць вычислена алгоритмом Лапцоша, то )Вз — сг,! < ~ОТьЯг —— ,Вь, где рг — единственный элемент блока Тг„, который (может бить) отличен от нулл; он находится в правом верхнем углу этого блока. 2 ~~АЯаиз) — Яььг)бйг = Ц~ТгьоД~г. Таким образом, разность между числом Ритце В; и некоторым собственным значением и матрицы А не превышает величины ЦОТьащ'йг, которол может быть много меньше, чем ~~ТьЯг. Если матрица Яь вычислена алгоритмом Ланцоша, то 'гТьаиДг = зуь~из(а)), где ьз(а) — последняя (зс-я) компонента вектора оз.

Эта формула дает дезаевый способ вычисления нормы невязки 3А(с,гщ)— (Яьоз)ВДг, не требующий умножения какого-либо вектора на матрицы Яг или А. 3. Не имея информации о спектре матрицы Т„, нельзя дать какую-либо содержательную оценку длл погрешности в векторе Ритца азиз. Если известно, что В; отделено расстоянием, не меньшим д, от прочих собственных значений матриц Тг и Т„, то угол В между Яавз и точным собственным вектором матрицы А можно оценить так: — з1п 2В < (7.2) 2 Если матрица Щ вычислена алгоритмом Ланцоша, то эта оценка упрощается к виду 1 . 1зг — з1п 2В < —.

2 Доказательство. 1. Числа Вз, ..., Вг входят в множество собственных значений матрицы Т = О" Т Так как ~~т — Т~~ = О Т"" = ~!Т „1Ь, то из теоремы Вейля (теоремы 5.1) следует, что собственные значения матриц Т и Т разнятся не более чем на ОТг„'Ог. Нужный результат следует из того, что собственные значения матриц Т и А одинаковы. 381 7.3. Алгоритм Лаицоша в точной арифметике 2. Имеем /)А(5)ьи;) — (Ц;,и;)В,/!г = //Я Адьо,) — С~ Яьиг)дг//г так как Тьи; = В,и; = 00тьаи,'Ог. По теореме 5.5, матрица А должна иметь собственное значение а, такое, что )сг — В;) < ЙТг„и;~(г. Если матрица Яь вычислена алгоритмом Ланцоша, то ЙТь о;~(г = Д|о;(а)~, потому что в Тьа может быть отличен от нуля только элемент Д, стоящий в правом верхнем углу.

3. Мы обратимся еще раз к примеру 5.4, чтобы показать, что нельзя получить содержательной оценки для погрешности в векторе Ритца, не имея информации о спектре матрицы Т„. Пусть т= '+д где 0 < е < д. Положим и = 1 и Щ = ~е~], тогда Т~ = 1+д, и приближенным собственным вектором является сам вектор ем В примере 5.4 было показано, что собственные векторы матрицы Т близки к [1, е/д)т и ( — е/д, 1)т.

Поэтому без нижней границы для д, т. е. для отделенности данного собственного значения матрицы Ть от всех прочих собственных значений, включая и собственные значения матрицы Т„, нельзя оценить погрешность в вычисленном собственном векторе. Если же такая нижняя граница имеется, то к матрицам Т и Т + Е = с)1а8(Ты Т„) можно применить вторую оценку теоремы 5.4; в результате получится неравенство (7.2). О 7.3. Алгоритм Ланцоша в точной арифметике Алгоритм Лаицоша для вычисления собственных значений симметричной матрицы А соединяет метод Ланцоша для построения крыловского надпространства (алгоритм 6.10) с процедурой Рэлея — Ритца, описанной в предыдущем разделе. Иными словами, из ортонормированных векторов Ланцоша строится матрица Яь = ~дш...,дг), и в качестве приближенных собственных значений матрицы А принимаются числа Ритца, т.

е. (см. (7.1)) собственные значения симметричной трехдиагональной матрицы Ть = ьД АЯ». Алгоритм 7.1. Алгоритпм Ланцоиьа для вычисления собственных значений и собственных векторов матрицы А = Ат в точной арифметике: д, = ЬД5! ~„6, = О, до = О 7от г'=1 1ой г = Адт из=а г г = г — сказ — Дз 1 до А =!!4!2 если Рз = О, то прекратить выполнение алгоритма 382 Глава 7. Итерационные методы для задач ла собственные значения 7з+ =г/Вз Вычислить собственные значения и собственные векторы матрицы Т и оценки погрешностей в них ецио /ог В этом разделе мы изучим характер сходимости алгоритма Ланцоша с помощью подробного разбора численного примера. Пример был сконструирован таким образом, чтобы проиллюстрировать как типичную картину сходимости, так и некоторые отклонения от нее, которые мы называем ложкой сходимостью.

Ложная сходимость может иметь место вследствие того, что начальный вектор й~ почти ортогонален собственному вектору для желаемого собственного значения, или потому, что среди собственных значений матрицы есть кратные (или очень близкие). Название данного раздела указывает, что влияние ошибок округления в нашем примере было (почти полностью) исключено. Разумеется, Ма11аЬ-программа, использованная для примера, выполнялась в арифметике с плавающей точкой (эта программа помещена на НОМЕРАСЕ/Ма11аЬ/ ЬапсзозРп11КеогФЬой.ш). Однако алгоритм Ланцоша в ней (в частности, внутренний цикл алгоритма 7.1) был реализован особо тщательным (и дорогостоящим) образом, чтобы как можно ближе воспроизвести точные результаты.

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

Тем не менее, мы укажем способ надежного вычисления собственных значений с помощью более дешевого алгоритма. Пример 7.1. Мы проиллюстрируем алгоритм Ланцоша и оценивание погрешностей в нем на примере большой диагональной матрицы А порядка 1000. Большинство ее собственньгх значений выбрано случайно из нормального гауссовского распределения. На рис.

7.1 дана картина собственных значений. Чтобы упростить понимание последующих графиков, диагональные элементы матрицы А были упорядочены так, чтобы с ростом индекса происходило их убывание. Элементы аи суть собственные значения Л;(А), а соответствующими собственными векторами являются столбцы е; единичной матрицы. Имеется небольшое число собственных значений с наибольшими абсолютными величинами; все остальные группируются вблизи центра спектра.

В начальном векторе Ланцоша д~ все компоненты, кроме одной, равны; об этом мы скажем ниже. Экспериментируя с диагональной матрицей, мы не теряем общности: проведение алгоритма Ланцоша для матрицы А при начальном векторе д~ эквивалентно проведению алгоритма для 1„1тАЯ при начальном векторе Ятд~ (см. вопрос 7.1). 383 7.3. Алгоритм Ланцоша а точной арифметике Собстаенные значения матрицы А 0 100 600 600 400 600 600 700 600 600 1МЮ Индекс Рис. 7.1. Собственные значения диагональной матрицы А. Чтобы лучше выявить характер сходимости, мы приводим несколько диаграмм типа тех, что показаны на рис.

7.2. На атом рисунке собственные значения каждой матрицы Та помещены в соответствующий столбец и. Верхняя диаграмма соответствует значениям )с от 1 до 9, а нижняя — значениям )с от 1 до 29. Столбец и содержит и плюсов, по одному для каждого собственного значения. В дополнительном столбце в правой части диаграммы помещены собственные значения матрицы А. Мы используем такой символьный код: наибольшее и наименьшее собственные значения каждой матрицы Та указаны крестиками, вторые по величине собственные значения (с обоих концов спектра)— квадратиками, третьи по величине в звездочками и, наконец, четвертые по величине — кружочками.