Т. Кормен, Ч. Лейсерзон, Р. Риверст, К. Штайн - Алгоритмы. Построение и анализ (2013) (1162189), страница 190
Текст из файла (страница 190)
Полученная в результате матрица размером (и — 1) х (и — 1) А' свит/а11 (28.9) называется дополнением 15уа (Яспш сошр!ешеп1) матрицы А по отношению к элементу аы. Мы утверждаем, что если матрица А невырождена, то невырождено и дополнение Шура. Почему? Предположим, что дополнение Шура, представляющее собой матрицу размером (и — 1) х (и — 1), вырождено. Тогда по теореме Г.1 она имеет ранг, строго меньший и — 1. Поскольку нижние и — 1 элементов в первом столбце матрицы равны О, нижние и — 1 строк должны иметь ранг, строго меньший и — 1. Таким образом, ранг всей матрицы строго меньше и. Применив результат упр. Г2.8 к уравнению (28.8), находим, что ранг матрицы А строго меньше и, так что согласно теореме Г.1 получаем противоречие с начальным условием невырожденности А.
где с = (сг,сз,...,га) = (п21,а22,..., аа1) пРедставлЯет собой вектоР-столбец длин1и» (и 1) ви = (виз) и13» ° ° ~ ша) = (п12, с12, °, п1а) является вектором- строкой длиной (и — 1), а А' — матрицей размером (и — 1) х (и — 1). Используя матричную алгебру (проверить полученный результат можно с помощью умножения), разложим матрицу А следующим образом: ббб Часть ИЬ Иэбранные темы Поскольку дополнение Шура невырождено, можно рекурсивно найти его 1.()- разложение. Запишем А' еэээт/оы = Г,'( где Ь' — единичная нижнетреугольная матрица, а à — верхнетреугольная матри- ца. Тогда, прибегнув к матричной алгебре, получим 1 О ~(сам эн о/аы !н 1 / 'э, О А' — ги~г/аы е/аы Т„э О ЕЛ/' е/аы 1.' О У~ что дает искомое 1Л)-разложение. (Заметим, что поскольку матрица э" — единичная ннжнетреугольная, таковой же является и матрица Ь; а поскольку матрица à — верхнетреугольная, таковой же является и матрица (7.) Конечно, если аы = О, этот метод не работает, поскольку при этом должно быть выполнено деление на О.
Он также не работает, если верхний левый элемент дополнения Шура А' — сэээт/аы равен нулю, поскольку в этом случае деление на нуль возникнет на следующем шаге рекурсии. Элементы, на которые в процессе 1.1)-разложения выполняется деление, называются ведущими (риком) и располагаются на диагонали матрицы (7. Причина, по которой в 1.()р-разложение включается матрица перестановок Р, состоит в том, чтобы избежать деления на нулевые элементы. Использование матрицы перестановки для того, чтобы избежать деленна на нуль (нли на малые величины, что вносит вклад в численную неустойчивость), называется выбором ведущего элемента (р)нойпй).
Важным классом матриц, для которых 1.1)-разложение всегда работает корректно, является класс симметричных положительно определенных матриц. Такие матрицы не требуют выбора ведущего элемента, так что описанная стратегия может быть применена без опасения столкнуться с делением на нуль. Мы докажем это (а также некоторые другие утверждения) в разделе 28.3. Наш код для 1.()-разложения матрицы А следует рекурсивной стратегии, хотя рекурсия в нем и заменена итеративным циклом (такое преобразование является стандартной оптимизацией процедуры с оконечной рекурсией, когда последней операцией процедуры является рекурсивный вызов ее самой; см.
задачу 7А), В коде предполагается, что размерность матрицы А хранится в атрибуте А. гошж Мы инициализируем матрицу У нулями ниже диагонали, а матрицу Ь единицами на диагонали и нулями — выше нее. Каждая итерация работает с квадратной подматрицей, используя ее верхний левый элемент в качестве ведущего для вычисления векторов е и и и дополнения Шура, которое становится квадратной подматрицей, обрабатываемой на следующей итерации. В61 Глава 2В.
Работа с матрицами ?.??-?)ЕСОМРОБ1Т!ОХ(А) 1 и = А.говав 2 ? и ?/ являются новыми матрицами размером и х п 3 Инициализируем (/ нулями ниже диагонали 4 Инициализируем ? единицами на диагонали и нулями выше нее 5 Гог?с = 1!оп б иьь = аьь 7 Гог ! = ?с + 1 го и 8 ?ия = аця/аьь // аия хранит и, 9 иы = аы // аы хранит ы, 10 !ог 1 = ?с + 1 го и 11 !ог 2' = )с + 1 го и 12 ап —— аб — ?се ив 13 гегпгп Г и ?/ Внешний цикл Гог, начинающийся в строке 5, выполняет каждый шаг рекурсии ло одному разу. В теле этого цикла в строке 6 ведущим элементом полагается элемент иьь = аьы В цикле Гог в строках 7-9 (который не выполняется, когда ?с = и) для обновления матриц ?/ и Г, используются векторы и и ю.
В строке 8 определяются элементы Г„находящиеся ниже диагонали, и выполняется сохранение и,/аьь в элементе !м, а в строке 9 вычисляются элементы ??, расположенные нвд диагональю, и выполняется сохранение ю; в иы. И наконец в строках !0 — !2 вычисляются элементы дополнения Шур~ которые затем сохраняются в матрице А (в строке 12 не нужно выполнять деление на аьы так как оно уже выполнено в строке 8 при вычислении ?сь). В связи с тройной вложенностью в циклы строки 12 время работы процедуры ?.??-?ЗесомрОЕ171014 равно вэ(п~). На рис.
28.1 показана работа процедуры ?Л.?-?3есомгоз1т10?ч. Здесь проиллюстрирована стандартная оптимизация процедуры, при юторой мы храним значащие элементы ?, и ?/ без привлечения дополнительной памяти, непосредственно в матрице А. Иначе говоря, мы можем установить соответствие между каждым элементом а, и либо 1, (если 1 > 2), либо и,. (если ! < 2) и обновлять матрицу А так, чтобы по окончании работы процедуры в ней хранились матрицы ? и (/.
Чтобы получить псевдокод этой оптимизации из приведенного выше, необходимо просто заменить каждое обращение к ! и и на обращение к а; можно легко убедиться в том, что такое преобразование сохраняет корректность алгоритма. Вычисление ?Я?Р-рахложении В общем случае при решении системы линейных уравнений Ах = б может оказаться, что необходимо выбирать ведущие элементы среди недиагональных элементов матрицы А для того, чтобы избежать деления на О.
Причем нежелательно не толью деление на О, но и просто на малое число, даже если матрица А невырождена, поскольку в результате можно получить численную неустойчивость. В связи с этим в качестве ведущего элемента следует выбирать наибольший возможный элемент. Чаешь РИ. Избранные текм 862 2 3 ! 5 (9 г3:::-)(;::$,. 2 3 ! 5 6 13 5 !9:,'У! 4 2 4 Зф;:Д,:4..
2 191023 сто)1б 9 18 ! 44) ! 2 4 Ю11 31 Й 4 9 21 2:-'!!,' 7 17 2 3 ! 5 314 2 4 ! 46;. ! У~3 (в) (в) [г) (в) 2 19 10 23 1 4 1 О О О 1 2 (л) Математические основы ШР-разложения аналогичны Ш-разложению. Вспомним, что имеется невырожденная матрица А размером и х и и требуется найти матрицу перестановки Р, единичную нижнегреугольную матрицу Ь и верхнетреугольную матрицу (/, такие, что РА = Х, (/. Перед разделением матрицы А на части, как при вычислении Ш-разложения, мы перемешаем ненулевой элемент, скажем, оы, откуда-то из первого столбца матрицы в позицию (1, 1) (если первый столбец содержит только нулевые значения, то матрица вырождена, поскольку ее определитель равен О (см.
теоремы Г.4 и Г.5)). Для сохранения множества исходных линейных уравнений мы меняем местами строки 1 и /с, что эквивалентно умножению матрицы А на матрицу перестановки (;) слева (см. упр. Г 1.4). Тоща мы можем записать ЯА как оА-( ", ), где и = (аз),пз(,..., а„!), за исключением того, что пы заменяет аы, ш Т (паз, оьз,..., оьа) т; а А' является матрицей размером (н — 1) х (и — 1). Поскольку оы ~ О, можно выполнить почти те же преобразования, что и при Ш-разложе- Рнс. 28.1.
Работа процедуры с()-Пясомвоз(т(он. (а) Матрица А. (8) Элемент аы = 2 (в черном круге) является ведущим, заштрихованный спшбец представляет собой о/аы, а заштрихованная строка — ш~. Вычисленные к этому моменту элементы (Г находятся над горизонтальной линией, а элементы 2 — слева от вертикальной линии. Матрица дополнении Шура А' — ош /аы занимаег т нюкиюю правую часть. (в) Теперь работа продаюкастся с матрицей дополнения Шура из части (б). Элемент аэг = 4 в черном круге яшшется ведущим, а заштрихованные столбец и строка представляют собой с/агз и шт (в разбиении дополнения Шура) соответственно.
Линии делят матрицу на вычисленные к этому моменту элементы (г (вверху), элементы Ь (слева) и новое дополнение Шура (внизу справа). (г) После следующего шага матрица А разлшкена. (Элемент 3 в новом дополнении Шура становится частью (7 при завершении рекурсии.) (д) Разложение А = б(/. Глава 2В. Работа с матрицами ббз нин, но теперь с гарантией, что деления на нуль не будет: ЯА = е/аы 1„1 О А' — еют/аы Как мы видели в Ш-разложении, если матрица А невырождена, то невырождено н дополнение Шура А' — еют/аы. Таким образом, мы можем индуктивно найти его ШР-разложение, дающее единичную нижнетреугольную матрицу 1', верхнетреугольную матрицу Г и матрицу перестановки Р', такие, что Р'(А' — ею /аы) = Е'Г .
Определим матрицу которая является матрицей перестановки в силу того, что она представляет собой произведение двух матриц перестановки (см. упр. П) А). Мы имеем РА = что и дает искомое ШР-разложение. Поскольку 1/ — единичная нижнетреугольная матрица, такой же будет н 1, и поскольку à — верхнетреугольная матрица, такой же будет и У.
Заметим, по в отличие от Ш-разложения, на матрицу перестановки Р' должны умножаться и вектор-столбец е/аы, н дополнение Шура А' — ею г/аы, Вот как выглядит псевдокод ШР-разложения. ( )ал (1" )( А, ( 'а", )( 1 О Р'е/аы 1„1 1 О Р'е/аы 1 -1 ( ) с)( 1„1 О А' — еют/ацз аы ю О А' — т/а с ам ют О Р'(А' — ею г/аы ) О Ь'Г аы ю т х О Г ) Часть Их Избранаые атаы 1Л)Р-Овсомеоз1т1оы(А) 1 и = А. госсе 2 и[1 .. и] — вновь созданный массив 3 1ог ь' = 1 1о и 4 к[1] =1 5 Гог)с = 11ои б р= О 7 1ог г' = )с 1о и 8 Ы]ась[ > р 9 р = [ось[ 10 )с' = 1 И ар ==О 12 еггог "вырожденная матрица" 13 Обменять з [)с] с х[й'] 14 1ог1 = 1 1о и !5 Обменять аы с аьч 1б аког 1 = )с + 1 1о и 17 аьь = оьь/оьь 18 Гог 7' = к + 1 1о и 19 ао — — об — асяаху Как и в процедуре 1.1)-13псомгощтюм, в псевдокоде 1ЛЗР-разложения 1.1)Р13исомвоз1тюм рекурсия заменяется итерацией.
В качестве усовершенствования непосредственной реализации рассмотренной рекурсии мы динамически поддерживаем матрицу перестановки Р в виде массива х, где я[1] = 7' означает, что 1-я строка массива Р содержит 1 в столбце 7'. Мы также реализуем код, который вычисляет Ь и У "на месте", в матрице А, т.е. по окончании работы процедуры )с если 1 > 7', ап —— и; если с <7'. На рис. 28.2 показано, как процедура 1Л)Р-13иСомгощтюм выполняет разложение матрицы.