Система линейных алгебраических уравнений.

Рассмотрим систему линейных уравнений с постоянными коэффициентами:

а₁₁х₁ + а₁₂х₂ + а₁₃х₃ + … + а_1тх_т = b₁,

а₂₁х₁ + а₂₂х₂ + а₂₃х₃ + … + а_2тх_т = b₂,

а₃₁х₁ + а₃₂х₂ + а₃₃х₃ + … + а_3тх_т = b₃,

………………………………………….

а_m1х₁ + а_m2х₂ + а_m3х₃ + … + а_mтх_т = b_m .

В матричной форме записи эта система принимает вид Ax = b, где

, , .

Будем предполагать, что матрица А задана и является невырожденной. Известно, что в этом случае решение системы существует, единственно и устойчиво к входным данным, т.е. рассматриваемая задача является корректной.

В системе Ах = b из п уравнений с п неизвестными можно было бы пытаться вычислять обратную матрицу А^-1 и умножать её на вектор b. Однако, на практике этот подход часто сталкивается с вычислительной неустойчивостью: вещественные числа хранятся в памяти лишь приближённо и ошибки округления могут накапливаться. Другой подход - так называемое LUP-разложение - в меньшей степени подвержен этой опасности и, к тому же, работает он раза в три быстрее.

Общая схема LUP-разложения Мы строим три (nn)-матрицы L, U и Р, образующие LUP-разложение матрицы А в следующем смысле: PA = LU причём

• L является нижне-треугольной матрицей, и ее главная диагональ заполнена единицами,

• U является верхне-треугольной матрицей,

• Р является матрицей перестановки.

Всякая невырожденная матриц обладает LUP-разложением. Зная его, мы сводим систему (2) к двум системам с треугольными матрицами, которые решаются просто. Умножим обе части уравнения Ах = b слева на матрицу Р, получаем уравнение

РАх = Рb (умножение на Р соответствует перестановке уравнений). Тогда LUх = Рb. Остаётся решить две системы с треугольными матрицами. Сначала из нижнетреугольной системы

Lу =Рb (4)

мы найдём у = Ux с помощью «прямой подстановки», а затем из нижнетреугольной системы Uх = у

найдём у= Uх с помощью «обратной подстановки».

Вектор х будет искомым решением: матрица Р обратима, поэтому

Ах = Р^-1LUx = Р^-1Ly= Р^-1Pb = b.

Вычисление lu-разложения (метод Гаусса).

Опишем рекурсивный алгоритм построения LU-разложения (nn)-матрицы А. Случай п = 1 тривиален: L = I₁ и U = А. При n > 1 разобьём матрицу А на четыре части:

где v—вектор-столбец из п - 1 элементов, w^T – вектор-строка из п-1 элементов, А’ - матрица размера (n - 1)  (n - 1). Paзложим А в произведение двух матриц:

.

Нули в первой и второй матрице обозначают вектор-строку и вектор-столбец из п - 1 элементов. Разность А' –vw^Т/a₁₁ называется дополнением Шура элемента а₁₁ в матрице А. Дополнение Шура - матрица размера (n-1)  (n - 1) и мы можем найти для него LU-разложение

A' – vw^T/a₁₁ = L’U’

с помощью рекурсивного вызова. Теперь A=LU и мы получили искомое LU-разложение. Заметим, что L’- нижнее-треугольная матрица с единицами на главной диагонали, а U’- верхнее-треугольная, и, следовательно, матрицы L и U также имеют такую структуру.

Конечно, при а₁₁ = 0 метод не годится (деление на ноль). По той же причине он не годится, если элемент в левом верхнем углу матрицы

А’-vw^Т/a₁₁ равен нулю (с деления на него начинается рекурсивный вызов).