прямым ходом метода Гаусса. Вычисление неизвестных из (3.5) называется обратным ходом метода Гаусса. Неизвестные i x вычисляются так: из последнего уравнения находим , n x из предпоследнего − 1 − n x и т.д. Итак , 1 ∑ + = − = n i j j ij i i x c y x ; 1 ,..., 1 − = n i . n n y x = (3.7) Изложенный алгоритм называется схемой единственного деления метода исключения Гаусса. Недостатки данного метода. Если ведущий элемент ) 1 ( − k kk a на каком- либо шаге окажется равным нулю, то эта схема формально непригодна (из-за деления на нуль метод падает), хотя заданная система уравнений может иметь единственное решение. Кроме того, если определитель системы не равен нулю, но в процессе вычислений встречаются ведущие элементы ) 1 ( − k kk a , которые достаточно малы по сравнению с другими элементами соответствующей строки, то это обстоятельство способствует усилению 109 отрицательного влияния погрешностей округления на точность результата. Вычислительная погрешность может существенно нарастать. Число действий (трудоемкость метода). Для реализации метода Гаусса (схема единственного деления) требуется примерно 3 3 2 n арифметических операций, причем подавляющее число этих действий совершается на этапе прямого хода. Русские математики Клюев и Коковкин-Щербак доказали, что в общем случае система линейных алгебраических уравнений не может быть решена точными методами с помощью меньшего числа операций, чем требуется в Гауссовском исключении. 3.2. Метод Гаусса с выбором главного элемента по столбцу Если в схеме единственного деления на некотором шаге ведущий элемент ) 1 ( − kkka окажется равным нулю, то схема единственного деления не может быть реализована, несмотря на то, что система bxAr r = (3.1) имеет единственное решение. В этом случае может быть применен метод Гаусса с выбором главного элемента. Главный элемент можно выбирать по строчке, по столбцу (частичный выбор) или по всей матрице A. Изложим наиболее простой и широко применяемый метод исключения Гаусса с выбором главного элемента по столбцу. Выбираем максимальный по модулю элемент в первом столбце матрицы A . Пусть | | max | | 1 1 1 inimaa≤ ≤ = , т.е. максимальный по модулю элемент первого столбца стоит в m -й строке. Меняем местами первое и m -е уравнения системы с соответствующей перенумерацией элементов этих уравнений. Делим новое первое уравнение на коэффициент при 1 x , имеем , 1 1 2 12 1 yxcxcxnn= + + + 11 1 1 aacjj= , nj,..., 3 , 2 = , 11 1 1 aby= , (3.8) и, исключив 1 x из второго, третьего,…, n -го уравнений точно так же, как в пункте 3.1, получаем систему 110 = + + + = + + + = + + + = + + + + ... ... ... ... ... , , , ) 1 ( ) 1 ( 3 ) 1 ( 3 2 ) 1 ( 2 ) 1 ( 3 ) 1 ( 3 3 ) 1 ( 33 2 ) 1 ( 32 ) 1 ( 2 ) 1 ( 2 3 ) 1 ( 23 2 ) 1 ( 22 1 1 3 13 2 12 1 nnnnnnnnnnnnbxaxaxabxaxaxabxaxaxayxcxcxcx(3.9) Далее первое уравнение системы (3.9) оставляем без изменения. Выбираем максимальный по модулю элемент во втором столбце матрицы системы (3.9), начиная со второй строки. Пусть | | max | | ) 1 ( 2 2 ) 1 ( 2 inilaa≤ ≤ = , т.е. максимальный по модулю элемент второго укороченного столбца стоит в l - той строке. Меняем местами второе и l -тое уравнения системы (3.9). Исключаем 2 x из всех уравнений, начиная с третьего, и т.д. После n шагов получаем систему уравнений с верхней треугольной матрицей. Обратный ход метода Гаусса (т.е. вычисление неизвестных ix , ni,..., 1 = ) осуществляем так же, как и в схеме единственного деления. Заметим, что дополнительная работа по выбору главных элементов в схеме частичного выбора требует порядка 2 n действий, что практически не влияет на общую трудоемкость метода. Известно, что для некоторых классов матриц при использовании схемы единственного деления главные элементы гарантированно располагаются на главной диагонали и поэтому применять частичный выбор нет необходимости. Именно так, например, обстоит дело с системами, имеющими положительно определенные матрицы, а так же матрицы с диагональным преобладанием. Замечание. Метод Гаусса с частичным выбором главного элемента можно применять к любой СЛАУ с невырожденной матрицей. Однако этот метод не всегда вычислительно устойчив. Например, пусть методом Гаусса с выбором главного элемента по столбцу решается система bxAr r = с матрицей коэффициентов − − − − − − = 1 1 1 1 1 1 1 1 1 0 1 1 1 0 0 1 A 111 Приведение такой матрицы к треугольному виду прямым ходом метода Гаусса равносильно следующей последовательности эквивалентных преобразований матрицы A : ⇒ − ⇒ − − − ⇒ − − − − − − = 8 0 0 0 4 1 0 0 2 0 1 0 1 0 0 1 4 1 0 0 4 1 0 0 2 0 1 0 1 0 0 1 2 1 1 0 2 1 1 0 2 0 1 0 1 0 0 1 1 1 1 1 1 1 1 1 1 0 1 1 1 0 0 1 AВ случае матрицы такого типа размерности nn× прямой ход метода Гаусса допускает рост элементов матрицы до ) 1 ( 2 − n. При больших n это может привести если не к переполнению разрядной сетки ЭВМ, то к сильному влиянию погрешностей округлений на конечный результат. Этот недостаток метода исключения Гаусса побуждает к разработке других точных методов, может быть более сложных, но численно более устойчивых, таких, например, как метод вращений, который будет изложен ниже. 3.3. Компактная схема Гаусса Как было показано в параграфе 3.1., для того, чтобы решить систему (3.1), нам достаточно знать коэффициенты ijc и правые части iy системы (3.5). Явно их можно вычислить следующим образом. Пусть осуществлены первые ) 1 ( − k шагов, т.е. уже исключены 1 1 ,..., − kxx, следовательно, система (3.1) приведена к виду: = + + = + + = + + + = + + + + = + + + + + − − − − − − − − − − ... ... ... ... , , ... ... ... ... ... , ... , ... ) 1 ( ) 1 ( ) 1 ( ) 1 ( ) 1 ( ) 1 ( 1 , 1 , 1 1 2 2 2 2 1 1 1 2 12 1 knnknnkknkkknkknkkkkknnkkkkknnkknnkkbxaxabxaxayxcxcxyxcxcxyxcxcxcx(3.10) Рассмотрим k -е уравнение системы (3.10) и предположим, что 0 ) 1 ( ≠ − kkka. Разделив обе части этого уравнения на ) 1 ( − kkka, получаем knnkkkkkyxcxcx= + + + + + , 1 1 , , (3.11) 112 где ) 1 ( ) 1 ( − − = kkkkkjkjaac, nkj,..., 1 + = , ) 1 ( ) 1 ( − − = kkkkkkaby Теперь уравнение (3.11) умножаем последовательно на ) 1 ( − kika и вычитаем, соответственно, из ( i -го) уравнения системы (3.10). В результате последняя группа уравнений системы (3.10) принимает вид = + + = + + = + + + + + + + + + + + + , ... ... ... ... , , ) ( ) ( 1 ) ( 1 , ) ( 1 ) ( , 1 1 ) ( 1 , 1 , 1 1 , knnknnkkknkknknkkkkkknnkkkkkbxaxabxaxayxcxcx(3.12) где , ) 1 ( ) 1 ( ) ( kjkikkijkijcaaa− − − = nkji,..., 1 , + = ; , ) 1 ( ) 1 ( ) ( kkikkikiyabb− − − = . ,..., 1 nki+ = Таким образом, в прямом ходе метода Гаусса коэффициенты уравнений преобразуются по формулам , ) 0 ( kjkjaa= njk,..., 2 , 1 , = ; , ) 1 ( ) 1 ( − − = kkkkkjkjaac, ,..., 1 nkj+ = nk,..., 2 , 1 = ; (3.13) , ) 1 ( ) 1 ( ) ( kjkikkijkijcaaa− − − = nkji,..., 1 , + = , nk,..., 1 = (3.14) Вычисление правых частей системы (3.5) осуществляется по формулам , ) 0 ( kkbb= , ) 1 ( ) 1 ( − − = kkkkkkabynk,..., 1 = ; (3.15) , ) 1 ( ) 1 ( ) ( kkikkikiyabb− − − = nki,..., 1 + = (3.16) По формулам (3.13) – (3.16) вычисляем коэффициенты ijc , правые части iy , , ,..., 1 ni= , ,..., 1 nij+ = системы (3.5) и затем осуществляем обратный ход, как и в схеме единственного деления. Реализация прямого хода метода Гаусса по формулам (3.13)-(3.16) называется компактной схемой Гаусса. 113 2.4. Метод оптимального исключения Методом оптимального исключения называется следующий алгоритм решения СЛАУ. Дана система = + + + = + + + = + + + ... .......... .......... , , 2 2 1 1 2 2 2 22 1 21 1 1 2 12 1 11 nnnnnnnnnnbxaxaxabxaxaxabxaxaxa(3.17) Делим первое уравнение на 11 a . Получаем ) 1 ( 1 , 1 ) 1 ( 1 2 ) 1 ( 12 1 + = + + + nnnaxaxaxУмножаем данное уравнение на 21 a и вычитаем его из второго уравнения системы (3.17). Получившееся второе уравнение делим на коэффициент при 2 x и исключаем 2 x из преобразованного первого уравнения. После первого шага система принимает вид: = + + + + = + + + + = + + + = + + + ... ... ... ... ... ... , , , 3 3 2 2 1 1 3 3 3 33 2 32 1 31 ) 2 ( 2 ) 2 ( 2 3 ) 2 ( 23 2 ) 2 ( 1 ) 2 ( 1 3 ) 2 ( 13 1 nnnnnnnnnnnnnbxaxaxaxabxaxaxaxabxaxaxbxaxax Теперь исключаем 1 x и 2 x из третьего уравнения. Преобразованное третье уравнение делим на коэффициент при 3 x и исключаем 3 x из первых двух уравнений, и т.д. В итоге получаем iiyx= , ni,..., 2 , 1 = , т.е. мы вычислили решение системы (3.17). Последовательное введение уравнений экономит память ЭВМ. Метод позволяет решать системы, порядок которых в 2 раза больше чем при решении с помощью схемы единственного деления. 3.5. LU- разложение Доказано, что с алгебраической точки зрения метод исключения Гаусса (схема единственного деления) эквивалентен разложению матрицы A на произведение двух треугольных матриц LUA= , где L – нижняя треугольная матрица, а U – верхняя. Пусть имеется матрица njiaA1 , ) ( = . Если 114 известно разложение LUA= , то решение системы bxAr r = сводится к решению двух систем с треугольными матрицами byLr r = и yxUr r = Решение СЛАУ с треугольными матрицами не вызывает затруднений. Возникает вопрос, когда такое разложение возможно. Ответ на этот вопрос дает теорема. Теорема об LU-разложении. Пусть все главные миноры матрицы A отличны от нуля ( 0 ≠ ∆ j, nj,..., 1 = ), тогда матрицу A можно представить единственным образом в виде произведения LUA= , где L – нижняя треугольная матрица с ненулевыми диагональными элементами, U – верхняя треугольная матрица с единицами по главной диагонали. Доказательство Доказательство будем вести по индукции по размерности матрицы. 1 = n: A=( ) 11 a→ 11 aL= , ) 1 ( = U, 0 11 ≠ a2 = n: = 22 21 12 11 aaaaA Матрицы L и U будем искать в виде = 22 21 11 0 lllL, = 1 0 1 12 uUТак как должно выполняться равенство LUA= , то, перемножая матрицы, для определения неизвестных ijl и iju получаем систему уравнений 11 11 al= , 12 12 11 aul= , 21 21 al= , 22 22 12 21 alul= + Решая эту систему, находим 11 11 al= , 11 12 12 aau= , 21 21 al= , 11 12 21 22 11 22 aaaaal− = По условию 0 11 ≠ a, 0 12 21 22 11 ≠ − aaaa, поэтому система имеет единственное решение и 0 11 ≠ l и 0 22 ≠ lПусть утверждение теоремы справедливо для матриц порядка ) 1 ( − kДокажем справедливость утверждения для матрицы порядка k . Представим матрицу A следующим образом: = − − − kkkkkabaAA1 1 1 r r , где = − − kkkkaaa, 1 1 1 r , ( ) 1 , 1 1 − − = kkkkaabr 115 Матрицы k L и k U будем искать в виде = − − kk k k k l l L L 1 1 0 r , = − − 1 0 1 1 k k k u U U r , где ) ,..., ( 1 , 1 1 − − = k k k k l l l r , T k k k k k u u u u ) ,..., , ( , 1 2 1 1 − − = r – вектор- строка и вектор-столбец. По предположению индукции 1 1 1 − − − = k k k U L A . Исходя из матричного равенства = A = − − − − − − − 1 0 0 1 1 1 1 1 1 1 k k kk k k kk k k k u U l l L a b a A r r r r , (3.18) для определения неизвестных 1 − k l r , kk l и 1 − k ur получаем систему уравнений 1 1 1 − − − = k k k a u L r r , (3.19) 1 1 1 − − − = k k k b U l r r , (3.20) kk kk k k a l u l = + − − 1 1 r r (3.21) Так как из равенства 1 1 1 − − − = k k k U L A и условия 0 1 ≠ − k A следует, что 0 1 ≠ − k L и 0 1 ≠ − k U , то система уравнений (3.19)-(3.21) однозначно разрешима. Отсюда 1 1 1 1 − − − − = k k k a L u r r , 1 1 1 1 − − − − = k k k U b l r r и 1 1 − − − = k k kk kk u l a l r r Докажем, что 0 ≠ kk l Из разложения (3.18) следует, что ) det( ) det( ) det( 1 1 − − = k kk k k U l L A , отсюда, т.к. 0 ) det( ≠ k A , 0 ) det( 1 ≠ − k L , 0 ) det( 1 ≠ − k U , то и 0 ≠ kk l Справедливость разложения LU A = доказана. Осталось доказать единственность разложения. Будем вести доказательство от противного. Пусть 2 2 1 1 U L U L A = = (здесь нижний индекс означает номер, а не размерность), отсюда 1 2 1 2 1 1 − − = U U L L (3.22) Получили слева нижнюю треугольную матрицу, а справа − верхнюю, а это возможно только в том случае, если 2 1 1 L L − и 1 2 1 − U U − диагональные. Но на диагонали матрицы 1 2 1 − U U стоят единицы, следовательно, E L L U U = = − − 2 1 1 1 2 1 и 2 1 U U = , 2 1 L L = . Теорема доказана. Если все главные миноры матрицы A отличны от нуля, то справедливы рекуррентные формулы, позволяющие найти элементы матриц L и U :
116 11 11 au= , jjau1 1 = , 11 1 1 ualjj= , nj,..., 3 , 2 = , ∑ − = − = 1 1 ippiipiiiiulau, ni,..., 3 , 2 = , ∑ − = − = 1 1 ippjipijijulau, iiippijpjijiuulal∑ − = − = 1 1 , ni,..., 3 , 2 = , niij,..., 2 , 1 + + = . Здесь 1 = iil(3.23) Замечание. Так как схема единственного деления эквивалентна LU- разложению, то из доказанной теоремы следует, что если все главные миноры матрицы A отличны от нуля, то схема единственного деления работает (нет деления на нуль). Поделитесь с Вашими друзьями: |