Метод ортогонализации и метод сопряженных градиентов Введение К решению систем линейных алгебраических уравнений приводятся многие задачи численного анализа. Известное из курса высшей алгебры правило Крамера для решения систем линейных алгебраических уравнений практически невыгодно, так как требует слишком большого количества арифметических операций и записей. Поэтому было предложено много различных способов, более пригодных для практики. Используемые практически методы решения систем линейных алгебраических уравнений можно разделить на две большие группы: так называемые точные методы и методы последовательных приближений. Точные методы характеризуются тем, что с их помощью принципиально возможно, проделав конечное число операций, получить точные значения неизвестных. При этом, конечно, предполагается, что коэффициенты и правые части системы известны точно, а все вычисления производятся без округлений. Чаще всего они осуществляются в два этапа. На первом этапе преобразуют систему к тому или иному простому виду. На втором этапе решают упрощенную систему и получают значения неизвестных. Методы последовательных приближений характеризуются тем, что с самого начала задаются какими-то приближенными значениями неизвестных. Из этих приближенных значений тем или иным способом получают новые «улучшенные» приближенные значения. С новыми приближенными значениями поступают точно также и т.д. При выполнении определенных условий можно придти, вообще говоря, после бесконечного числа шагов. Рассмотрим два точных метода: метод ортогонализации и метод сопряженных градиентов. 1. Метод ортогонализации 1.1 Метод ортогонализации в случае симметрической матрицы Пусть дана система (1) порядка n. Чтобы избежать в дальнейшем путаницы, над векторами поставим черточки. Решение системы будем разыскивать в виде , (2) где – n векторов, удовлетворяющих условиям при (3) Здесь рассматривается обычное скалярное произведение векторов в n-мерном векторном пространстве, т.е. если и , то . Пусть такие векторы найдены. Как это делается, будет показано ниже. Рассмотрим скалярное произведение обеих частей системы (1) с (4) Используя (2) получим: (5) или, в силу выбора векторов , . (6) Итак, для определения коэффициентов получили систему с треугольной матрицей. Определитель этой системы равен . (7) Следовательно, если , то возможно найти и находятся они без труда. Особенно легко определятся , если матрица А симметрическая. В этом случае, очевидно, (8) и, следовательно, =0 при . (9) Тогда система для определения примет вид (10) и . (11) Метод можно обобщить. Пусть каким-то образом удалось найти систему 2n векторов так, что =0 при . (12) Умножая обе части равенства (1) на и используя представление через , как и ранее, получим: . (13) Опять получилась система линейных алгебраических уравнений с треугольной матрицей для определения . Несколько усложнив вычисления можно получить систему диагонального вида. Для этого построим три системы векторов , так что имеют место равенства: (14) (15) (16) Тогда , (17) так как при i<r (18) и при i>r (19) Таким образом, (20) Остановимся подробнее на первом из описанных методов. Рассмотрим случай, когда матрица А симметрическая и положительно определенная. Последнее означает, что для любого вектора квадратичная форма его компонент больше или равна нулю, причем равенство нулю возможно в том и только том случае, если вектор нулевой. Как мы видели ранее, нужно построить систему векторов , удовлетворяющих условиям =0 . (21) Это построение можно осуществить следующим образом. Исходим из какой-то системы линейно независимых векторов , например из системы единичных векторов, направленных по координатным осям: (22) Далее проводим «ортогонализацию». Принимаем и ищем в виде . (23) Из условия находим: (24) Ищем в виде . (25) Условия влекут за собой (26) Далее поступаем также. Процесс будет осуществим, так как все . Это же обеспечит нам разрешимость системы для определения коэффициентов . Заметим, что в нашем случае это будет процесс настоящей ортогонализации, если в пространстве векторов ввести новое скалярное произведение при помощи соотношения . (26) Нетрудно проверить, что введенное таким способом скалярное произведение будет удовлетворять всем требованиям, которые к нему предъявляются. При решении системы n уравнений по настоящей схеме требуется произвести (28) операций умножения и деления. 1.2 Метод ортогонализации в случае несимметрической матрицы В случае несимметрической матрицы процесс ортогонализации проводится точно также. Пусть векторы уже построены. Тогда ищется в виде (29) Коэффициенты определяются из системы (30) Система в случае несимметрической матрицы будет треугольной. Аналогично строится система «биортогональных» векторов, т.е. система 2n векторов, удовлетворяющих условию (12). При этом – n произвольных линейно независимых векторов, а векторы строятся последовательно в виде (31) Коэффициенты находятся из системы (32) Также поступаем, отыскивая коэффициенты и , при построении систем векторов (14) и (15), удовлетворяющих условиям (16). При этом получим две системы: (33) из которых и определяем и . Остановимся еще на одном методе ортогонализации. Будем рассматривать строки матрицы А как векторы: (34) Ортонормируем эту систему векторов. Первое уравнение системы делим на . При этом получим (35) где (36) Второе уравнение системы заменится на (37) где (38) Аналогично поступаем дальше. Уравнение с номером i примет вид (39) где (40) Процесс будет осуществим, если система уравнений линейно независима. В результате мы придем к новой системе , где матрица С будет ортогональной, т.е. обладает свойством СС¢=I. Таким образом, решение системы можно записать в виде . (41) Практически, вследствие ошибок округления, СС¢ будет отлична от единичной матрицы и может оказаться целесообразным произвести несколько итераций для системы . 2. Метод сопряженных градиентов 2.1 Первый алгоритм метода Пусть требуется решить систему линейных алгебраических уравнений (1) с положительно определенной матрицей A порядка n. Рассмотрим функционал , (2) представляющий многочлен второго порядка относительно x1 , x2 , …, xn . Обозначим через решение системы (1), т.е. . В силу симметричности и положительной определенности матрицы, имеем: При этом знак равенства возможен лишь при . Таким образом, задача решения уравнения (1) сводится к задаче отыскания вектора , обращающего в минимум функционал (2). Для отыскания такого вектора применим следующий метод. Пусть – произвольный начальный вектор, а (4) – вектор невязок системы. Покажем, что вектор невязок имеет направление нормали к поверхности в точке . В самом деле, направление нормали совпадает с направлением быстрейшего изменения функции в точке . Это направление мы найдем, если найдем среди векторов , для которых , такой вектор, что имеет наибольшее значение. Но Но среди векторов постоянный длины достигает максимального значения, если имеет направление вектора или ему противоположное. Утверждение доказано. Будем двигаться из точки в направлении вектора до тех пор, пока функция достигает минимального значения. Это будет при , т.е. при . (5) Вектор (6) и принимаем за новое приближение к решению. В методе сопряженных градиентов следующее приближение находится так. Через точку проведем гиперплоскость (n-1) – го измерения (7) и через обозначим новую невязку системы . (8) Вектор направлен по нормали к поверхности в точке , а вектор параллелен касательной плоскости в этой точке. Поэтому . (9) Гиперплоскость (7) проходит через точку , так как . При любом вектор параллелен некоторой нормальной плоскости к поверхности в точке . Найдем среди них тот, который лежит в гиперплоскости (7), т.е. ортогонален к . Из условия ортогональности имеем: , или . (10) Вектор (11) имеет направление нормали к сечению поверхности гиперплоскости (7) в точке . Будем двигаться из точки в направлении вектора до тех пор, пока функция достигнет минимума. Это будет при . (12) Вектор примем за новое приближение к решению системы. Вектор невязок (13) имеет направление нормали к поверхности в точке . Покажем, что он будет ортогонален к и . В самом деле, используя (9), (11), (12), (13), имеем: Рассмотрим гиперплоскость (n-2) – х измерений , (14) проходящую через точку . Эта гиперплоскость содержит и , так как мы ранее видели, что , а . Вектор при любом параллелен гиперплоскости (7), так как . Подберем так, чтобы он был параллелен и гиперплоскости (14), т.е. потребуем ортогональности к вектору . Будем иметь: , или (15) Вектор (16) будет иметь направление нормали к сечению поверхности гиперплоскостью (14) в точке . Из точки сместимся в направлении этого вектора так, чтобы функция достигла минимального значения. Это будет при , (17) (18) примем за новое приближение к . Новый вектор невязок будет: . (19) Продолжая процесс, получим последовательности векторов , , , определяемые рекуррентными соотношениями: (20) Для этих векторов имеют место следующие соотношения: (21) (22) В самом деле, в силу самого построения при i¹j Далее, при i>j Если i=j+1, то правая часть равна нулю, в силу определения , если же i>j+1, то , по доказанному, и . Продолжая понижение индекса у вектора , через несколько шагов придем к скалярному произведению (по определению ). Таким образом, соотношения (21) доказаны. Для доказательства (22), в силу равноправия индексов i и j, предположим, что i>j. Тогда . Так как в n-мерном векторном пространства не может быть более n взаимно ортогональных векторов, то на некотором шаге получим , т.е. будет решением системы (1). На рис. 1 показана геометрическая картина нашего построения при n=3. Рис. 1 2.2 Второй алгоритм метода Приведем другой алгоритм метода. Будем обозначать последовательные приближения к решению через и введем обозначения: . (23) Первые два приближения и возьмем так, чтобы . (24) Предположим, что уже известно приближение (i³1), вычислены и справедливо равенство . (25) Будем искать минимум функционала (2) на множестве векторов . (26) Приравнивая к нулю частные производные от по и для определения и , получим систему: (27) или, учитывая (25), (28) Обозначим через решение этой системы: (29) и за (i+1) – е приближение к решению примем: (30) Из системы (27) следует, что , (31) а так как то из (31) следует: (32) Докажем, что если (33) то при всех i (34) что будет доказывать и сходимость, и конечность второго алгоритма. В самом деле, при условиях (33) и т.е. условие (24) выполнено. Предположим, что уже доказаны равенства (35) и докажем равенство При предположении (35) и, следовательно, Но из соотношений (20) имеем: т.е. Докажем коллинеарность векторов и (36) Из (20) и (29) имеем: а это и доказывает коллинеарность векторов (36). Вектор дает минимум функционала в плоскости, проходящей через и натянутой на векторы и , а мы показали, что этот минимум лежит на прямой, проходящей через в направлении вектора . Но на этой прямой минимум функционала достигается на векторе . Это и означает, что Это и доказывает справедливость (34) при всех i. На первый взгляд кажется, что первый алгоритм лучше, так как на каждом шаге он требует лишь одного умножения матрицы А на вектор , а во втором алгоритме требуется два умножения матрицы А на вектор и , но опыт показал, что применение первого алгоритма приводит к быстрому накоплению ошибок округления, так что для матриц большого порядка возможно существенное отклонение от точного решения. Второй алгоритм менее чувствителен к ошибкам округления и поэтому требует меньшего количество шагов для получения хорошего приближенного решения. Метод сопряженных градиентов целесообразно использовать для решения систем уравнений, в которых матрица А имеет много нулевых элементов. При решении системы по этому методу элементы матрицы участвуют в арифметических операциях лишь при умножении матрицы на вектор, а умножение матрицы на вектор можно организовать так, чтобы в арифметических операциях участвовали только ненулевые элементы. Заключение В данной работе были рассмотрены метод ортогонализации и метод сопряженных градиентов, а также представлена программа на языке программирования С++, реализующая метод ортогонализации на ЭВМ, и ее результаты работы. Список литературы 1. Березин И.С. и Жидков Н.П. Методы вычислений. т. 1. М.: «Наука», 1965. 633c. 2. Воеводин В.В. Численные методы алгебры (теория и алгоритмы). М.: «Наука», 1966. 3. Подбельский В.В. и Фомин С.С. Программирование на языке Си. М.: «Финансы и статистика», 2000. 599 с. 4. Калиткин Н.Н. Численные методы. М.: «Наука», 1978. 512 с.