Раздел: Документация
0 ... 92 93 94 95 96 97 98 ... 177 □Пусть дана невырожденная квадратная матрица размерности N. Назовем ее А. Будем предполагать, что данная матрица может быть представлена в виде произведения нижней треугольной матрицы L и верхней треугольной матрицы U. □Перемножим матрицы L и U и поставим в соответствие полученным для элементов матрицы А выражениям их числовые значения. В результате получим систему из уравнений, каждое из которых описывается следующей обшей формулой (i™0, t...N-1rj-0, I...N-1): Так как н в матрице 1 и и матрице U главная диагональ заполнена, то в системе будет N2+N неизвестных. Чтобы уменьшить количество неизвестных до Ыг, Краут предложил считать, что все элементы главной диагонали матрицы L равны 1. Если принять это предположение, то систем» легко решается. □Краут обнаружил, что выведенная им система решается обычными подстановками в том случае, если неизвестные ищутся в правильном порядке. Детали следующие. •Начинаем поиск, перебирая столбцы матрицы А, а не строки (то есть верхним циклом должен быть никл по j, а не по i). Перейдя к очередному столб[[у, последовательно перебираем все элементы в нем посредством цикла по i. •Если i<j, то на основании А и определенных ранее элементов матриц L и U находим UpJ. Для этого используется следующая формула: И U . = А. . - V L, .и. , k = 0 В данной формуле сумма вычисляется, если i>0. Если же i-О, то сумму полагают равной нулю. •Если i>j, то вычисляется элемент матрицы U с индексами i j. Для этого применяется следующая зависимость: Н L. . =-- JU, . В этой формуле сумма п рос* i иты веется, если j>0. Если j=0, го L— A yiJ „ Проделав мысленно несколько итераций алгоритма Краута, вы обнаружите, что используемые в приведенных выше формулах элементы матриц L и U всегда создаются ло того, как будут востребованы. Также заметьте, что каждый отдельный элемент разлагаемой матрицы А участвует в вычислениях только один раз, давая элемент с такими же индексами п матрице L или в матрице U, В общем же, алгоритму Краута нужно проделать совсем немного вычислений, чтобы произвести LU-разложение. 1 In языке ирофаммироиаиия Mathcad алгоритм Краута реализуется следующим кодом: lu bfabi(A) л n +- rows( А) (L*~ kicntjty(n) U for j € 0..П - 1 for t € 0.. j L-L)
for i € j,. t) 4.j (L U) A. , -itl к =o J.J /v Программа lu brain довольно простая, но в пей есть несколько моментов, которые стоит пояснить. □В алгоритме Краута считается, что на диагонали матрицы L лежат единицы. Чтобы учесть это, в качестве основы матрицы L нужно взять единичную матрицу соответствующей размерности. Создать такую матрицу можно, задействовав встроенную функцию Identity. □В качестве основы матрицы U должна быть взята нулевая матрица той же размерности, что и матрица L Создать такую матрицу можно, отняв от матрицы L ее же □В рабочих формулах метода Краута суммы вычисляются лишь, если i>0 (формула для U ) или )Х> (формула для L). В противном случае суммы полагаются равными 0. Чтобы это учесть, используя минимальный объем кода, применяем функцию iffcond. true, false), где cond — условие, true - значение, которое должно быть возвращено, если условие выполнится, false — значение, возвращаемое, если условие не выполнится. Проверим программу lubrain: М := 3> 6 LUr=m brain(M)
Разложение было осуществлено верно. Однако любая ли невырожденная матрица может быть разложена метолом Краута? Оказывается, нет, К ггримеру, при попытке разложить следующую матрицу Mathcad выдаст сообщение об ошибке? Divide by zero in function evaluation (При вычислении функции происходит деление на нуль). ( I 2 6 4 8-1 ,-2 3 t\ j Проблема заключается в том. что если значением диагонального элемента матрицы U окажется 0, то при вычислении элемента матрицы L может возникнуть ошибка деления на 0. Что же делать, если это произойдет? То, что диагональный элемент оказался равным 0. это не более чем результат неудачного сочетания элементов в исходной матрице А. Но можно ли это сочетание изменить в более удачную сторону, не изменив саму систему? Конечно. С учетом того, что нет никакой разницы, в какой последовательности записаны уравнения в системе, мы можем поменять строки матрицы А местами. Если аналогичную модификацию произвести с вектором правых частей и вектором неизвестных, то система не изменится. Однако сочетание элементов в матрице коэффициентов, с точки зрения метода Краута, после перестановки может оказаться более удачным. Технически перестановка строк н матрице А осуществляется с помощью так называемой матрицы перестановок Матрица перестановок — это соразмерная матрице А матрица, в каждой строке и столбце которой имеется только один ненулевой элемент, равный 1. Несложно догадаться, что количество возможных матриц перестановок для матрицы размерности N составляет N1-1 (единичная матрица не является матрицей перестановки). Так, если N-.3, то матриц перестановок будет пять:
Исходя из правил матричного умножения очевидно, что, при умножении матрицы Л слева на матрицу перестановок, если в строке i ненулевой элемент находится в столбце j, то j-н строка матрицы А переместится на позицию i. Итак, нам необходимо написать программу, которая будет генерировать матрицы перестановок произвольной размерности. Сделать это в Mathcad не так уж и просто. Мы используем для этого следующий алгоритм. □Пусть нужно создать матрицу перестановок размерности N. Для начала сгенерируем вектор v, содержащий числа от 0 до N- 1. Каждое число будет адресовать одни столбец в матрице перестановок. □Запускаем цикл от 0 до N -1 и последовательно перебираем строки матрицы перестановок □Перейдя к очередной строке, случайным образом выбираем из вектора у элемент. Его значение укажет, в элемент, принадлежащий какому столбцу, должна быть помещена единица. □Так как в каждом столбце может быть только один ненулевой элемент, вырезаем из вектора v уже использованное значение индекса столбца. Это гарантирует, что данное значение не будет разыграно повторно. □После того, как матрица будет сгенерирована, щюверясм. не является ли она единичной. Если нет, то возвращаем ее как результат. Если же матрица получилась единичной, то генерируем матрицу перестановок заново, 0 ... 92 93 94 95 96 97 98 ... 177
|
