Т. Кормен, Ч. Лейсерзон, Р. Риверст, К. Штайн - Алгоритмы. Построение и анализ (2013) (1162189), страница 189
Текст из файла (страница 189)
В этой главе делается особый упор на задачу умножения матриц и решения систем линейных уравнений. Основы теории матриц излагаются в приложении Г. В разделе 28.1 речь идет о решении систем линейных уравнений с использованием ШР-разложения. Затем в разделе 28.2 изучается тесная взаимосвязь между умножением и обращением матриц. Наконец в разделе 28.3 рассматриваются важный класс симметричных положительно определенных матриц и их использование для решения переопределенных систем линейных уравнений методом наименьших квадратов. Одним из важнейших вопросов, возникающих на практике, является численная устойчивость (пшпепса! з!аЬ(1(гу). Из-за ограниченной точности представления действительных чисел в реальном компьютере в процессе вычислений могут резко нарастать ошибки округления, что приводит к неверным результатам. Такие вычисления называются численно неуетойчивымо (пшпепсайу ипз1аЫе).
Несмотря на важность данного вопроса, мы лишь поверхностно коснемся его в данной главе, так что мы рекомендуем читателям обратиться к отличной книге Голуба (Оо!иЬ) и ван Лоана (тап 1.оап) [143), в которой детально рассматриваются вопросы численной устойчивости. 28.1. Решение систем линейных уравнений Во многих приложениях требуется решать системы линейных уравнений. Систему линейных уравнений можно сформулировать в виде матричного уравнения, в котором каждая матрица или вектор принадлежат полю, обычно полю действительных чисел К.
В этом разделе рассматривается решение систем линейных уравнений с помощью метода, называемого 1Л)р-разложением. Глава 2В. Работа с матрицами В53 НаЧНЕМ С СИСТЕМЫ ЛИНЕЙНЫХ УРаВНЕНИй С П НЕИЗВЕСТНЫМИ Х1, Хг,..., Хп: Й11Х1 + о12Х2 + ° ° + о1пХп = Ь1, о21Х1 + о22Х2 + ' + о2пхп = Ь2, (28. 1) оп1Х1 + опгхг + + а пхп = Ь„ Решением уравнений (28.1) является множество значений х1, хг,..., Хп, которое удовлетворяет всем уравнениям одновременно. В этом разделе мы рассматриваем только случай, когда имеется ровно и уравнений с и неизвестными. Можно удобно переписать уравнения (28.1) в виде матрично-векторного урав- нения Х1 Ь, хг Ьг а11 аю .
а1п О21 О22 ' ' ' О2п Оп1 апг апп или, полагая А = (а; ), х = (х,) и Ь = (61), в эквивалентной форме как Ах = Ь. Если матрица А невырождена, имеется обратная к ней матрица А ', и х=А 16 (28.2) (28.3) является вектором решения. Можно доказать, что х является единственным ре- шением уравнения (28.2), следующим образом. Если имеется два решения, х и х', то Ах = Ах' = 6 и, обозначая тождественную матрицу как (, = (А 1А)х = А 1(Ах) = А 1(Ах') = (А 'А)х' В этом разделе мы в основном будем иметь дело с невырожденной матрицей А, или, что эквивалентно (согласно теореме Г.1), с матрицей А, ранг которой равен количеству и неизвестных. Существуют, однако, и другие варианты, которые заслуживают краткого рассмотрения.
Если количество уравнений меньше количества и неизвестных (или, говоря более обобщенно, если ранг матрицы А меньше и), то система линейных уравнений является недоопределеииой (цп1(егдегепп(пей). Обычно недоопределенная система линейных уравнений имеет б54 Часть гу Избрааные темы бесконечно много решений, хотя может не иметь решений вовсе, если уравнения системы оказываются несовместимыми. Если количество уравнений превышает количество и неизвестных, система линейных уравнений является переопределенпой (очегс(е[епп)пей), и у нее может не существовать решений. В разделе 28.3 рассматривается важная задача поиска хороших приближенных решений переопределенных систем линейных уравнений.
Вернемся к нашей задаче решения системы линейных уравнений Ая = Ь из и уравнений от и неизвестных. Мы можем вычислить А ь, а затем, воспользовавшись уравнением (28.3), умножить Ь на А 1 и получить я = А 16. На практике этого подхода избегают из-за его численной неустойчивости. К счастью, другой подход — 1ЛЗР-разложение — численно устойчив и обладает тем преимуществом, что на практике оказывается более быстрым.
Обзор 1Л)Р-разложения Идея, лежащая в основе 1ЛЗР-разложения, состоит в поиске трех матриц, Е, (5 и Р, размером и х и, таких, что (28.4) РА=ИУ, где Š— единичная нижнетреугольная матрица; (5 — верхнетреугольная матрица; Р— матрица перестановки. Матрицы Т,, Ьг и Р, удовлетворяюшие уравнению (28.4), называются йЮР-разложением (1Л)Р десошрогйбоп) матрицы А. Мы покажем, что всякая невырождеиная матрица А допускает такое разложение.
Преимущество вычисления ШР-разложения матрицы А основано на том, что система линейных уравнений решается гораздо легче, если ее матрица треугольна, что и выполняется в случае матриц Т и сГ. Найдя ШР-разложение матрицы А, мы можем решить уравнение (28.2), Ая = Ь, путем решения треугольной системы линейных уравнений, как показано далее. Умножая обе части уравнения Ах = 6 на Р, мы получим эквивалентное уравнение РАх = РЬ, которое в соответствии с упражнением Г.!.4 эквивалентно перестановке уравнений из (28.1).
Используя разложение (28.4), получаем Т(5х = РЬ. Теперь можно решить полученное уравнение, решив две треугольные системы линейных уравнений. Обозначим у = сГя, где к — решение исходной системы линейных уравнений. Сначала решим нижнегреугольную систему линейных уравнений (28.5) Ьу = РЬ найдя неизвестный вектор у с помощью метода, который называется прямой подстановкой. После этого, имея вектор у, решим верхнетреугольную систему ли- Глава за Работа о матрияопи В55 нейных уравнений (7к= у, (28.б) найдя х с помощью обратной подстановки. Поскольку матрица перестановки Р обращаема (см.
упр. Г2,3), умножение обеих частей уравнения (28А) на Р ' дает Р 'РА = Р 1Ш, так что А =Р 1Ш (28. 7) Следовательно„вектор х является искомым решением системы линейных уравне- ний Ак = Ь: А =Р и7к (согласно (28.7)) (согласно (28.б)) (согласно (28.5)) = Р-1Т,у — Р 1РЬ Теперь рассмотрим, как работают прямая и обратная подстановки, а затем приступим к самой задаче вычисления Шр-разложения.
Ьи[1) = 6.~,), Ьи13[ г У1 (21У1 + Уг (з1У1 + (згУ2 + Уз (гг1У1 + (п2У2 + (пЗУЗ + ' ' ' + Уп = Ьа~п~ Из первого уравнения получаем, что у1 = Ь„г16 Зная значение У1, его можно подставить во второе уравнение и найти у, = 6.[,[ - („у, Теперь можно подставить в третье уравнение два найденных значения, у1 и уз, и получить УЗ Ьи~з! (131У1 + 132У2) . Праман н обратнаи подстановки Прямая подстановка (Тогагагд зпЬзпш11оп) позволяет решить нижнетреугольную систему линейных уравнений (28.5) для данных Е, Р и Ь за время гЭ(пз). Для удобства мы представим перестановку Р в компактной форме с помощью массива я[1 ..
п[. Элемент я[1] приз = 1, 2,..., и указывает, что Р, „~й = 1 и Р;, = О для 5' Ф к[1[. Таким образом, в матрице РА на пересечении 1-й строки и 5-го столбца находится элемент а [1)сь а 1-м элементом РЬ является 6 Ьр Поскольку Т, является единичной нижнетреугольной матрицей, уравнение (28.5) можно переписать следующим образом: а56 Часть е71 Избранные тены В общем случае для поиска у, мы подставляем найденные значения уы уг,..., у; 1 в ьие уравнение и находим е — 1 После того как мы нашли у, найти х из уравнения (28.6) можно, воспользовавшись обратной подстановкой (Ьаск яиЬябшбоп), которая аналогична прямой. В этом случае мы решаем первым и-е уравнение и работаем в обратном направлении, к первому уравнению. Как и прямая подстановка, этот процесс требует времени 6(нг).
Поскольку П является верхнетреутольной матрицей, систему линейных уравнений (28.6) можно переписать как иыхг + иггхг + . + и1н-гх„-г + иг„гх„з + из„х„= у1, иггхг + + иг, -гх„г + иг,н-1х„-з + игнхн = Уг, ин — г,н — гха — г + ин — г,н — 1х~ — 1 + ин — г,нхн = Ун — г ин ьв гх„ з + и„ 1„х„ = ун г, и~ах~ = ун. Таким образом, можно последовательно найти х„, х„ы..., х1 следующим обра- зом: хн — ун/ин,н х„1 = (ун 1 — ин з ах„)(ин гн г, * -г = ~у -г — (и -г, -гх -1+ и -г, х ))lи -г, -г или, в общем случае, х; = у; — ~~~ и;ух /ии .
у=е+1 Процедура ШР-Яоьчи для заданных Р, Т,, У и 6 находит х путем комбинирования прямой и обратной подстановок. В псевдокоде предполагается„что размерность и хранится в атрибуте Т, гоесв и что матрица перестановки Р представлена массивом я. 857 1лааа 2В. Работа с матрицами ШР-Яоьуе(з, У, я, Ь) 1 и = 1.тоии 2 Пусть х и у — вновь созданные векторы длиной и 3 хогг = 1топ 4 у; = Ьи)б .7 =110уб 5 аког г' = и йозцп1о 1 6 х; = (у; — с 2", 1ибх ) /ин 7 ге1пгп х Процедура 1ЛЗР-Боьун находит у в строках 3 и 4 с помощью прямой подстановки, а затем вычисляет х с помощью обратной подстановки в строках 5 и б.
Наличие внутри каждого цикла зог неявною цикла суммирования приводит ко времени работы данной процедуры, равному 9(пз). В качестве примера применения данных методов рассмотрим систему линейных уравнений 344 х= 7 где А= 344 (7) и которую необходимо решить, найдя неизвестный вектор х. ШР-разложение матрицы А имеет вид О.б О.б 1 ( ) а5а Часть е7э'. Игранные тены (Вы можете убедиться в корректности разложения, проверив, что РА = И7.) Используя прямую подстановку, находим у из Ау = Р(э: 0.2 1 0 д, = З получив у= 14 путем вычисления сначала уп затем — уз и наконец — уз. Затем с помощью обратной подстановки находим к из Ук = у: 0 0.8 — 0.6 лз = 1.4 тем самым получив искомое решение исходной системы линейных уравнений л = 2.2 путем вычисления сначала лз, затем — лз и наконец — хь Вычисление 1Л)-разложения Мы показали, что если можно вычислить ШР-разложение невырожденной матрицы А, то для решения системы линейных уравнений Ах = Ь можно воспользоваться простыми прямой и обратной подстановками.
Осталось показать, как эффективно найти 1Л.1р-разложение матрицы А. Начнем со случая невырожденной матрицы А размером п х и и отсутствия матрицы Р (или, что эквивалентно, Р = 1н). В этом случае необходимо найти разложение А = Ь(э'. Матрицы Ь и (э' называются й()рлзложением (Ш десошроыйоп) матрицы А. Мы применим для построения ).(э'-разложения процесс, известный как метод исключения Гаусса (бацая е!цтлпайоп). Начнем с удаления первой переменной из всех уравнений, кроме первого, путем вычитания из этих уравнений первого, умноженного на соответствующий коэффициент. Затем та же операция будет проделана со вторым уравнением, чтобы в третьем и последующих уравнениях отсутствовали первая и вторая переменные.
Процесс будет продолжаться до тех пор, пока мы не получим верхнетреугольную матрицу — это и будет искомая матрица (э'. Матрица А составляется из коэффициентов, участвовавших в процессе исключения переменных. Мы реализуем описанную стратегию с помощью рекурсивного алгоритма. Итак, мы хотим построить Ш-разложение невырожденной матрицы А разме- Глава 28. Рабата с матрицами 859 ром и х и. Если и = 1, задача решена, поскольку мы можем выбрать Т = 11 и У = А. При и ) 1 разобьем А на четыре части; А=( е/а11 /„1 ) ( О А' — влит/а11 (28.8) Нули в первой и второй матрицах уравнения (28.8) представляют собой вектор- строку и вектор-столбец длиной и — 1 с нулевыми элементами. Матрица свит/а11, полученная тензорным произведением векторов с и щ и делением каждого элемента полученной в результате умножения матрицы на аы, представляет собой матрицу размером (и — 1) х (и — 1) (тот же размер имеет и матрица А', из которой она вычитается).