Форсайт Дж., Малькольм М., Моулер К. - Машинные методы математических вычислений (1040536), страница 41
Текст из файла (страница 41)
Можно организовать систематический поиск значения г, которое приносит минимум (надо надеяться, нулевой) функции о ,„ (А (г)); для этой цели можно использовать универсальные программы минимизации, такие, как ГМ(г(. Мы с большим успехом применялн эту технику для весьма разнообразных задач. В принципе можно изменить алгоритм Я70 таким образом, ЧтОбЫ ПрОИЗВЕдЕНИЕ бЕ1 ((7) бЕ1 ((7т) ВСЕГда ИМЕЛО фИКСнрОВаинмй знак + или —. Тогда станет возможно находить знак определителя матрицы А (г) и использовать программы вычисления нулей вроде ХЕК01й), а не программы минимизации. У нас не было большой практики в реализации этой идеи. Проблема собственных значений В этом пункте будет описано не реальное приложение Я70, а скорее указание на то, как сингулярное разложение связано с более широким классом важных матричных задач.
Пусть А — квадратная матрица. Собственные значения А— это числа Х, для которых Ах==Хх имеет ненулевое решение х. Поскольку это эквивалентно требованию, чтобы де1 (А — И) =О, то тсоретически собственные значения можно находить как корни такого полинома. Рассмотрим, однако, матрицу А ( где е — некоторое малое число, которым все же нельзя пренебречь, Правильные собственные значения суть Х, =-1+е и Х,=1 — е. Полипом с1е1(А — И) равен Х* — 21+ (1 — е'). Если е мало, то е' много меньше и в коэффициентах полинома необходимо иметь удвоенную точность по сравнению с точностью 9. НАИМВНЫНИЬ КВАДРАТЫ 238 матричных элементов или собственных значений.
Эта трудность становится еще более явственной для матриц большего порядка. Поэтому современные методы вычисления собственных значений избегают использования полиномов или алгоритмов для нахождения корней полиномов. Связь между сингулярным разложением и собственными значениями наиболее проста для матриц, являющихся симметричными и положительно полуопределенными. Симметрию легко определить и проверить: А симметрична, если аы — — ап для всех 9 и ~. Положительная полуопределенность — это гораздо более интимное свойство; она означает, что хгАх)0 для всех х. Грубо говоря, это значит, что диагональные элементы А довольно велики по сравнению с внедиагональными элементами. Действительно, достаточным условием будет выполнение для каждого ~ не авенства Р ап ~ ~ч' „)а,,), ! ~! однако эти неравенства не являются необходимыми. Нетрудно показать, что для действительной симметричной матрицы А все собственные значения — действительные числа, а если А положительно полуопределена, то все ее собственные значения неотрицательны.
Если А — -УХГ' и А=А"', то А9 АТА -РтгЦ9ЦЕ~/г = Ь'ВЧ~г. Отсюда Обозначим через п2 столбцы матрицы и; рассматривая зй столбец полученного уравнения, находим А'о, =- о';о,. Поскольку собственные значения матрицы А' суть квадраты соб- ~ ствеиных значений А, а эти последние неотрицательны, то мы заключаем, что собственные значения симметричной положитель- ' но полуопределенной матрицы равны ее сингулярным числам,, а собственными векторами будут столбцы Р. Если А симметричная, но не положительно полуопределениая, то некоторые из ее собственных значений отрицательны.
Оказывается, что модули собственных значений равны сингулярным числам, а знаки собственных значений можно найти, сравнивая столбцы У со столбцами к'; мы, однако, не будем вдаваться в детали. Если А не является симметричной, то между еа собственными значениями и сингулярными числами нет какой- либо простой связи. Все это показывает, что алгоритм ЯЛЭ можно использоват' для вычисления собственных значений симметричных матриц 9.с Вычисление Р1зложения 939 ФА. Вычнсление разлеження ЯЛ) — это также наименование подпрограммы для вычисления сингулярного разложения произвольной матрицы. Она основана на АЛГОЛ-процедуре Голуба и Райнша (Голуб, Райнш (!97!)). Прежде чем описывать, как работает ЯЧЭ, полезно рассмотреть более элементарный матричный алгоритм — гауссово исклю.
чение — с несколько иной точки зрения. Пусть 1 6 11 2 7 12 3 8 !3 4 9 14 5 10 15 Если бы Л, была матрицей коэффипиентов системы линейных уравнений, то первым шагом в решении системы было бы вычитание кратных первой строки из других строк с целью получить нули в первом столбце. (В этом простом примере нет необходимости выбирать ведущий элемент,) Мы получили бы ! 6 11 — !Π— 20 0 — 5 0 — !О А, .— 0 — 15 — 30 0 — 20 -40 Как можно описать преобразование Л, в Л„используя матрич- ные обозначенняр Пусть Однако такой способ ие особенно эффективен, и мы не рекомендуем его.
Нашей единственной целью было продемонстрировать, что сингулярное разложение тесно связано с матричными задачами иа собственные значения. Тщательное изучение алгоритма ЯЧО поможет также и в понимании алгоритмов для собственных значений. к нлименьшие квлдялты ее=(1 О О О О).
С помощью внешнего произведения т,е,' образуем матрицуМ,= т. =1 — т,е,: 10000 — 2 1 0 0 О -ЗО100 — 4 О 0 1,0 — 5 0 0,0 1 м,= Конечно, подпрограмма гауссова исключения в действительности не формирует матрицу М, и не выполняет матричного умножения М,А,. Она просто оперирует со строками А„получая из нее А,. Однако математическое описание этого процесса компактней всего в терминах матрицы М,. Кроме того, введение М, позволяет рассматривать преобразование с точки зрения операций над столбцами А„ а не над ее строками. Пусть а, есть 1зй столбец А. Тогда М,а.
=- (I — т,е,') а = а — (т,ет) аг (ег =а, — (а„) т1. В частности, М,а,=а,— 1 т,= Легко проверить, что А, = М,Аг. О 0 0 0 9.Е ВЫЧИСЛЕНИЕ РАЗЛОЖЕНИЯ М,а,=а,— 6 т,= М,а,=а,— 11 т,= Это показывает, что столбцы А,можно получить путем вычитания некоторых кратных вектора т, из столбцов А,. Следующий шаг исключения даст матрицу 1 б 11 0 -5 -10 0 0 0 О 0 О О О 0 е," = ~0 1 0 0 0) М, =! — т,е,', А, =М,А,. Алгоритм Я~Р основан на аналогичных, хотя и несколько более сложных идеях.
Реальная подпрограмма ныполняет различные операции как над строками, так и над столбцами матрицы. Наиболее прозрачное описание алгоритма использует ортогональные аналоги матриц М, и М,. Цель алгоритма состоит в том, чтобы найти ортогональные матрицы У и Г, такие, что матрица УРАР л Этот шаг можно описать с помощью гО 0 2 3 4 4 6 — 5 — 10 — 15 — 20 111 — 10 — 20 — зо — 40 э нлименьшие квхдглты 242 диагональная.
Обе эти матрицы получаются как произведения ортогональных матриц, называемых хаусхолдеровыми отраже- НИЯЛ4И. Алгоритм распадается на два этапа. Первый этап заключается в приведении А к дврхдиагональной форме, т. е. матрице, у которой ненулевыми могут быть лишь элементы диагонали и первой наддиагонали. Второй этап — это итерационный процесс, в котором наддиагопальные элементы уменьшаются до пренебрежимо малой величины, что дает желаемую диагональную матрицу. Снова рассмотрим введение нулей в а, — первьш столбец матрицы А,. Гауссово исключение использовать нельзя, поскольку М, неортогональная. Возьмем вектор 8.4! б -2.000 3.000 4.000 5.000 Он получается из а, прибавлением к первой компоненте числа 2а,~,:', Пусть и,-=1 нр „,„, А, .— 1/„А,, Преобразованная матрица А, действительно вычисляется стол- бец за столбцом; У,а, =- я-л а,— 14Г и,и,а, 1' ° а,— и, (так как Р игал =1) — 7.410 1 О О 0 0 Заметим, что в первый столбец были введены желаемые нули и что его длина не изменилась.
9.4. ВЫЧИСЛЕНИЕ РАЗЛОЖЕНИЯ 243 У,ЕА = их — ф,'и,и~~ах = а, — (2.796)и1 -17.529 1.409 -0.387 — 2.183 — 3.949 Гу,а, = а, — ф, 'и,и~~а, а, — (4.592)и, -27.642 2.817 — 0.774 — 4.366 -7.957 Следовательно, — 7.416 — 17.529 — 27.642 0 1.409 2.817 0 -0.387 — 0.774 0 -2.183 — 4.366 0 -3,979 — 7.957 Столбцы А, получаются вычитанием кратных вектора и, нз столбцов А,. Эти кратные порождаются скалярными произведениями, а не отдельными элементами, как в гауссовом исключении. Преобразующая матрица У, и есть хаусхолдерово отражение. В ее ортогональности, критически важной лля процесса, можно убедиться, проверив, что Уг(7,=(: 2Р-1п пг+ф-1н ит = 7.
Эта ортогональность гарантирует, что длина каждого столбца в А, равна длине соответствующего столбца в А,. Геометрически преобразование любого столбца ЕО в Уеаг состоит в отражении ат относительно плоскости, проходящей через 9. илимеиьшие квлдулты 244 начало координат и перпендикулярной к и,. Зта плоскость выбрана так, что а, отражается в вектор с четырьмя нулевыми компонентами.
Прежде чем вводить дальнейшие нули под диагональю, преобразованием вида А,= А,1', получают нули в первой строке. Нули, уже стоящие в первом столбце, не должны быть испорчены, а длину первой строки нуж- но сохранить, так что в данном примере, где а равно всего лишь 3, можно получить только один нуль. Для матрицы большего по- рядка на этом шаге было бы получено и — 2 нуля.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
