Ортега Дж., Пул У. Введение в численные методы решения дифференциальных уравнений. Под ред. А.А.Абрамова (1986) (1095855), страница 20
Текст из файла (страница 20)
Определите, является ли матрица коэффициентов симметричной и диагонально доминирующей; в) выполните снова задания из пунктов а и б, аппроксимируя с' односторонними разностями (3.2.24); г) выполните снова задания из пунктов а и б, переписав уравнение в виде 1(1+ х')и'1 — с =-х' и воспользовавшисьзатем разностной схемой(32.29); д) выполните снова задания из пунктов а и бири граничных условиях е'(0) = 1, е(1) = = О, а затем при граничных условиях с'(0) + о(0) = 1. с'(1) + (1/2) е(1) = О. З.З. Решение систем линейных уравнений Одной из проблем, наиболее часто. встречающихся в задачах научного программирования, является нахождение решения системы линейных уравнений (см., например, предыдущий раздел, а также разделы 4.3, 5.2 и 8.2) .
В этом разделе мы ограничимся рассмотрением систем, в которых число уравнений равно числу неизвестных. Будем записывать такие системы в виде (3.3.1) Ах=Ь, где А — заданная матрица размера п Х и, Ь вЂ” заданныйвектор-столбенел координатами, а х — вектор с и координатами, который и подлежит определению.
В приложении 3 дается краткий обзор основных понятий и результатов, используемых в этом разделе. Если матрица А является невырожденной, то решение системы (3.3.1) можно представить в виде х=А 1Ь, (3.3.2) где А ' — матрица, обратная к А. Однако в большинстве случаев не рекомендуется вычислять обратную матрицу, поскольку затрачиваемое при этом машинное время значительно превосходит время, требуемое дпя непосредственного решения системы (3.3.1) . В этом разделе рассмотрим методы непосредственного решения системы (3.3.1) в предположении, что все элементы матрицы системы хранятся в оперативной памяти ЗВМ.
Мы уделим также внимание ленточным матрицам, т.е. матрицам, все ненулевые элементы которых располагаются лишь на нескольких диагоналях, которые и хранятся в памяти. Простейшим случаем нетривиальных ленточных матриц являются трехдиагональные матрицы, которые, например, встречались в разделе 3.2. Другие типы разреженньгх матриц рассматриваются в гл.
8. 80 Имеется один алгоритм решения систем (3.3.1), который почти всегда, если учитывать и все его модификации, используется для решения небольших и умеренно больших (скажем, и = 200) линейных систем с заполненными матрицами (т.е. с матрицами, почти все элементы которых отличны от нуля). Этот алгоритм называется гауссовым исключением.
Сначала продемонстрируем работу этого алгоритма на примере системы трех уравнений с тремя неизвестными 4х, — 9х, + 2хз =2, 2х, — 4х. + 4хз = 3, (3.3.3) — х~ +2х, +2х, =1, которая в матричной форме записывается в виде 4 — 9 2 — 4 4 х! (3.3.4) хз 2 хз Первый основной шаг гауссова исключения состоит в исключении первой переменной х; из второго и третьего уравнений. Если из второго уравнения системы вычесть первое, умноженное на 0,5, и из третьего уравнения вычесть первое, умноженное на — 0,25, то получим эквивалентную систему уравнений: х, 4 — 9 0 0,5 0 — 0,25 2 3 2,5 (3.3.5) хз 1,5 хз Второй основной шаг состоит в исключении х, из третьего уравнения. Это может быть сделано вычитанием из третьего уравнения второго, умножен- ного на — 0,5, что приводит к системе х, 2 2 4 — 9 2 0 0,5 3 (3.3.6) Х2 0 0 4 хз.
Этот алгоритм основывается на том факте (обычно доказываемом во ввод ных курсах по линейной алгебре), что при замене какого-либо уравнения системы линейной комбинацией этого уравнения и некоторого другого уравнения решение системы не меняется. Проделанные операции называются элементарными преобразованиями строк. К этому моменту мы завершили первую.
часть алгоритма гауссова исключения, которую обычно называют прямым исключением или приведением к треугольной форме. Эта часть завершается тогда, когда все элементы последней строки системы, кроме крайне правого, обращаются в нуль (см. 3.3.6) ) . Вторая часть алгоритма заключается в решении полученной верхней треугольной системы (3.3.6) . Это легко осуществляется с помощью процесса обратной подстановки. Последнее уравнение (3.3.6) имеет вид 4хз = 2,5. Следовательно, хз = 0,625.
Подставим теперь это значение во второе уравнение: 0,5 ха + (3) (0,625) = 2. В! б. Дж. Ортега Отсюда х, = 0,25. Подстановка этих значений х2 и хэ в первое уравнение, дает 4х1 — (9) (0,25) + (2) (0,625) = 2, или х, = 0,75. Чтобы проверить найденное решение, выполним умножение 4 — 9 2 075 2 — 4 4 025 — 1 2 2 0625 результат которого совпадает с правой частью (3.3.3) В случае общей системы размера и Х п алгоритм гауссова исключения состоит из тех же шагов, что н в примере 3 Х 3. Если система записана в виде а21х1 + ...
+а2лх„= Ь2, а,,х, + ... +аэлх„=Ь,, (3.3.7) ал2Х2 + ... +а„лх„=Ь„, то первый шаг состоит в исключении х2 из последних и — 1 уравнений. Это достигается вычитанием из второго уравнения первого, умноженного на а2,/а... из третьего уравнения первого, умноженного на а,2/а2 2, и т.д. Этот процесс приводит к преобразованной системе уравнений а11х1+ а12х2 + + а2лх~ =Ь1, а х +...+а х =Ь~) и22 х2 ... аэл хл— а~ )х + ... + а~ )х =Ь~ ) ал2 Х2 ° алл Хл л (3.3.8) где и22 ал =ан — и2, —, Ь,. =Ь,-Ь,,;/=2,..., . < а;, а21 Применяя теперь тот же самый процесс к последним п — 1 уравнениям системы (3.3.8), исключаем х2 из последних и — 2 уравнений и т.д., пока вся система не приведется к треугольной форме: Ь, Ь2' ) а~~ а~2 ...
а2л Х2 12) О) а22 ' ' а2 Х2 (3.3.9) 0 а х Ь~ ЛЛ Л Л где верхние индексы, вообще говоря, указывают, сколько раз изменялись соответствующие элементы. Этим завершается фаза прямого исключения (или приведения к треугольной форме) алгоритма гауссова исключения. Решение треугольной системы (3.3.9) теперь легко получается на фазе обратной подстановки, в ходе которой уравнения системы (3.3.9) решаются ет в обратном порядке: (л — 1), (л — 1) хл = Ьл Флл (л — 2) (л-2) ~1 (л — 2) Хл .1 =(Ь„, — ал 1 лХл)/ал (3.3.10) х, =(Ь, --а12х2 — ... — а1„х„)/а11. (1 — 1) Выше мы считали, что элемент а11 и все элементы аи (1 = 2,..., и), которые использовались в качестве делителей, отличны от нуля.
Однако дело не всегда обстоит именно так, и если какие-либо из этих элементов обращаются в нуль или становятся слишком маленькими, то алгоритм должен быть модифицирован. Мы вернемся к этому вопросу в следующем разделе. Процесс гауссова исключения можно очень компактно записать в виде алгоритма. Прямое исключение для 1( = 1,... „п — 1, для 1' = 1с + 1,..., п: 1и» -а;»/а»»., (3.3.11) для 1=В+ 1,..., и: а," '- а;1 — 1;»а»1, Ь; Ь; — 1;Ь„.
Обра тная подстанов ка для Й=п, п — 1,...,1: л х» (Ь» — 2' а»1 х,)/а»» . (3.3.12) 1=»+1 При составлении программы для ЭВМ, реализующей этот алгоритм, следу- ет обратить внимание на то, что последовательно щзеобраэуемые в ходе этого (») процесса элементы а" можно записывать в те же ячейки памяти,. где Ц располагались элементы исходной матрицы а,1, На это указывает пятая строка алгоритма.
Если это будет сделано, то исходная матрица, разумеет(») ся, будет испорчена. Аналогичным образом, новые элементы Ь( можно записывать на место исходных элементов Ь, . Вычисляемые в третьей строке алгоритма множители 1г» можно записывать на место элементов а;», так как после вычисления соответствующего значения 1г» элемент а1» больше не нужен. Покажем теперь, как алгоритм гауссова исключения связан с фактори- зацией (разложением в произведение) матрицы А. Введем нижнюю тре- угольную матрицу 2. с единичными диагональными элементами, положив расположенные ниже диагонали элементы 1 равными множителям, ис- пользованным при исключении у.й переменной из 1-го уравнения.
В случае примера (3.3.4) получим 1 О 0 С= 05 1 0 — 0,25 — 0,5 1 (3.3.13) (3.3.14) что и является треугольной факторизацией матрицы А. Это представление назьвают также АУ-разложением матрицы А, Если теперь, используя прямую подстановку, решить треугольную систему ~у =Ь относительно у, то получим вектор правой части в (3.3.6) . Таким образом, алгоритм гауссова исключения для решения системь! Ах = Ь эквивалентен следующему трехшаговому процессу. 1.
Факторизовать А = А(У. 2. Решить Ау =Ь. (3.3.! 6) З.Решить Ух=у. Этот взгляд на процесс гауссова исключения с точки зрения теории матриц оказывается весьма полезным для теоретических рассмотрений, а также служит основой для некоторых вычислительных вариантов проведения процесса исключения. Очень важным является вопрос об эффективности апгоритма гауссова исключения.
Чтобы ответить на него хотя бы частично, мы сейчас оценим число арифметических операций, необходимых для вычисления вектора х, являющегося решением системь1. Основная часть работы приходится на выполнение строки алгоритма (3.3.11), расположенной в самом внутреннем цикле. Эта строка включает одно сложение (операции сложения и вычитания мы считаем эквивалентными) и одно умножение. Для подсчета общего числа этих операций мы просто имитируем циклическу1о структуру алгоритма: и — 1 и и и — 1 и число сложений = Х Х Х 1 = Х Х (л — Ь) = ь=! 1=а+11=ь+! ь=! !=а+! (2н — 1) пз и-1 и — 1 (Л 1)(Л) Х (и — Ус)2 = Х 3с2 = К=! 6 (3.3.17) 3 и это же выражение определяет число умножений. Мы пока подсчитали только число операций, необходимое для приведения исходной матрицы к треугольной форме.
Определенная работа требуется также для модификации вектора правой части Ь, вычисления множителей (р, и выполнения обратной подстановки. Однако здесь в каждом случае требуется не более чем порядка и' операций (см. упражнение 3.3.3). Таким образом, при достаточно больших и основное время затрачивается на приведение матрицы к треугольному виду, причем соответствующее число операций оказывается пропорциональным и . Чтобы получить представление о количестве машинного времени, которое может потребоваться для решения с помощью гауссова исключения системы умеренного размера, предположим, что и = 100 и что операция сло- Обозначим теперь через У верхнюю .треугольную матрицу, которая полу» чается после завершения фазы исключения.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
